COOPERATIVE OPTIMAL CONTROL METHOD AND SYSTEM FOR WASTEWATER TREATMENT PROCESS

Abstract

In a cooperative optimal control system, firstly, two-level models are established to capture the dynamic features of different time-scale performance indices. Secondly, a data-driven assisted model based cooperative optimization algorithm is developed to optimize the two-level models, so that the optimal set-points of dissolved oxygen and nitrate nitrogen can be acquired. Thirdly, a predictive control strategy is designed to trace the obtained optimal set-points of dissolved oxygen and nitrate nitrogen. This proposed cooperative optimal control system can effectively deal with the difficulties of formulating the dynamic features and acquiring the optimal set-points.

Description

CROSS REFERENCE TO RELATED APPLICATION

This application claims priority to Chinese Patent Application No. 201810499231.X, filed on May 23, 2018, entitled “Cooperative optimal control system for wastewater treatment process”, which is hereby incorporated by reference in its entirety.

TECHNICAL FIELD

In this present invention, a cooperative optimal control system (COCS) is designed to coordinate the different time-scale performance indices, realize the real-time optimal control of wastewater treatment process (WWTP), and improve the system performance. It is feasible to apply the proposed COCS to the dynamic WWTP, and realize the tracking control of the dissolved oxygen (S_O) and the nitrate nitrogen (S_NO). This system can not only address the cooperative problem of the different performance indices, but also conserve the operation cost, improve the effluent quality, and guarantee the efficient stable operation. The technology of this present invention is a part of advanced manufacturing technologies, which belongs to the field of control engineer and wastewater treatment engineer.

BACKGROUND

WWTP is conducted to effectively adsorb, decompose and oxidize the pollutants based on the physical, chemical and biological reaction, so that the pollutants can be degraded and separated from the wastewater, thereby realizing the purification. The efficient operation plays an important role in improving the optimal control performance and optimizing the system performance, which is an important strategy for alleviating scarcity of water resources. It has not only good economic benefits, but also significant environmental and social benefits. Therefore, the research results of the present invention have broad application prospects.

WWTP is a complex industrial system, and it has obvious nonlinearity and time-varying features. In addition, WWTP comprises different dynamic response time, multiple performance evaluation indices, which make it difficult to realize the optimal control. Since the performance indices are conflicted with each other, and influenced by different system response time, how to achieve the cooperation between the performance indices with different time scales is of great significance to improve the system performance and ensure its efficient and stable operation. It is significant to improve the operation performance based on the cooperative performance indices, which can not only improve the optimal operations, but also guarantee the control performance. Therefore, it is necessary to explore a novel optimal control method, where different time-scale performance indices can be formulated, and the optimal set-points of the control variables can be acquired, so that the requirements of the optimal control can be satisfied. Appropriate optimal control method can not only adjust the relationship of the performance indices, so as to guarantee the effluent quality and reduce the operation cost, but also provide suitable dynamic optimal set-points to reduce and cope with the occurrence of the abnormal conditions, and ensure the efficient and stable operation. In addition, by improving the automation level of WWTP, it can also effectively decrease the operation management and operation personnel, as well as minimize the operation cost.

SUMMARY

In this present invention, a COCS is designed for WWTP, where the different time-scale performance indices models are formulated based on the system response time, the different time-scale performance indices models are then optimized by cooperative optimization algorithm, so that the optimal set-points of S_O and S_NO can be obtained, and then a predictive control strategy is designed to trace the obtained optimal set-points of S_O and S_NO. In this present invention, a COCS is designed. This proposed COCS can address the different time-scale performance indices formulation problem by the two-level performance indices models. It can optimize the multiple time-scale performance indices by the cooperative optimization algorithm based on the data-driven assisted model, so that the problem about the optimal set-points of S_O and S_NO can be solved. Moreover, a predictive control strategy is designed to trace the obtained optimal set-points of S_O and S_NO. Based on COCS, not only the multiobjective optimal control problem can be approached, but also the performance indices can be improved, as well as the efficient and stable operation.

This invention develops a COCS for WWTP. This proposed COCS can solve the problem of description of performance indexes at different time scales by using a two-layer performance index modeling method; and utilize a cooperative optimization algorithm based on a data driven assisted model to solve the problem. The different time-scale performance index models are then optimized by the cooperative optimization algorithm, where the set-points of S_O and the concentration of S_NO can be obtained. Then a predictive control strategy is developed to trace the obtained optimal set-points of S_O and S_NO. Based on the proposed COCS, not only the system performance can be improved, but also efficient and stable operation can be guaranteed.

A cooperative optimal control system is developed for WWTP, its features include the design of control system framework and the optimal controller, the specific steps are:

(1) the design of the control system framework, the steps are:

1) the control system includes the serve, the PLC control cabinet, the blower, an oxidation-reduction potential (ORP) electrode sensor, a S_O sensor, a suspend solid (SS) sensor, a temperature sensor, a pH sensor, an ammonia nitrogen (S_NH) sensor, a total nitrogen (TN) sensor, a chemical oxygen demand (COD) sensor, where

the serve is the operation center of control system to optimize the control procedure and realize the monitoring functions;

the PLC control cabinet is used to complete the automatic detection, monitor and control of the equipment status and operation conditions;

the blower acts as an executor of the control system to realize the aeration and ventilation of aerobic tank;

an ORP electrode sensor, located in the anaerobic tank, is configured to measure the values of ORP;

a S_O sensor, located in the aerobic tank, is configured to measure the values of S_O;

a SS sensor, located in the aerobic tank, is configured to measure the values of SS;

a temperature sensor, located in the effluent tank, is configured to measure the values of temperature;

a pH sensor, located in the in the effluent tank, is configured to measure the values of pH;

a S_NH sensor, located in the effluent tank, is configured to measure the values of S_NH;

a TN sensor, located in the effluent tank, is configured to measure the values of TN;

a COD sensor, located in the effluent tank, is configured to measure the values of COD;

2) data collection and transmission in WWTP, first, it is necessary to achieve the data acquisition through an ORP electrode sensor, a S_O sensor, a SS sensor, a temperature sensor, a pH sensor, a S_NH sensor, a TN sensor, a COD sensor; then it is essential to connect the PLC control cabinet and the sensors by using RS 485, realize the communication, complete the upload of the collected data, and transfer the collected data to the serve through the coordinated communication standards;

(2) the optimal controller, which comprises the establishment of performance indices, the development of optimization algorithm and the design of the controller, the steps are:

1) select the related process variables of the pumping energy (PE): S_NO, mixed liquor suspended solids (MLSS), and choose the related process variables of the aeration energy (AE) and the effluent quality (EQ): S_O, SS, S_NH, S_NO;

2) formulate the two-level models based on the different time scales, the upper-level model is designed for PE, and the lower-level models are designed for AE and EQ:

F ₁(t₁)=l₁(x_u(t₁)),  (1)

f ₁(t₂)=l₂(x_l(t₂),x_u*(t₁)),

f ₂(t₂)=l₃(x_l(t₂),x_u*(t₁)),  (2)

where F₁(t₁) is a model of the PE at time t₁, l₁(x_u(t₁)) is a mapping function of the PE model, f₁(t₂) is a model of the AE at time t₂, l₂(x_l(t₂), x*_u(t₁)) is a mapping function of the AE model, x*_u(t₁) is an optimal set-points of nitrate nitrogen S_NO* at time t₁, f₂(t₂) is a model of the EQ at time t₂, l₃(x_l(t₂), x*_u(t₁)) is a mapping function of the EQ model, x_u(t₁)=[S_NO(t₁), MLSS(t₁)] is the input variables vector of the pumping energy at time t₁, S_NO(t₁) is concentration of nitrate nitrogen S_NO at time t₁, MLSS(t₁) is the concentration of the MLSS at time t₁, x_l(t₂)=[S_O(t₂), SS(t₂), S_NH(t₂)], S_O(t₂) is the concentration of S_O at time t₂, SS(t₂) is the concentration of SS at time t₂, S_NH(t₂) is the concentration of S_NH at time t₂, and [S_O(t₂), SS(t₂), S_NH(t₂), S_NO*(t₁)] are input variables vector of AE and EQ at time t₂;

(3) design a cooperative optimization algorithm to optimize upper-level and lower-level optimization problems to obtain optimal set-points of control variables, an optimization period of upper level is two hours, and an optimization period of lower level is thirty minutes, which includes:

{circle around (1)} formulate the upper-level and lower-level problems:

Min F₁(S_NO(t₁),MLSS(t₁)),  (3)

Min[f₁(S_O(t₂),S_NH(t₂),SS(t₂),S_NO*(t₁)),

f₂(S_O(t₂),S_NH(t₂),SS(t₂),S_NO*(t₁))],  (4)

where Min F₁(S_NO(t₁), MLSS(t₁)) is the upper-level optimization problem, Min [f₁(S_O(t₂), S_NH(t₂), SS(t₂), S_NO*(t₁)), f₂(S_O(t₂), S_NH(t₂), SS(t₂), S_NO*(t₁))] is the lower-level optimization problem;

{circle around (2)} set the number of particle population in the upper level optimization I₁, the number of the particle population in the lower level optimization I₂, a maximum number of iterations in the upper level optimization N₁, and a maximum number of iterations in the lower level optimization N₂, I₁=50, I₂=50, N₁=20, N₂=50;

{circle around (3)} introduce a single particle swarm optimization (SPSO) algorithm to optimize the upper-level optimization problem, a position and a velocity of the ith particle are:

s _i(t₁)=[s_i,1(t₁),s_i,2(t₁)],  (5)

v _i(t₁)=[v_i,1(t₁),v_i,2(t₁)],  (6)

s_i(t₁) is the position of the ith particle at time t₁, s_i,1(t₁) is the value of S_NO at time t₁, s_i,2(t₁) is the value of MLSS at time t₁, v_i(t₁) is the velocity of the ith particle at time t₁, i is the number of particles, i=1, 2, . . . , 50, an update process of s_i(t₁) and v_i(t₁) are

v _i,d(t₁+1)=0.7v_i,d(t₁)+0.72(p_i,d(t₁)−s_i,d(t₁))+0.72₂(g_d(t₁)−s_i,d(t₁)),  (7)

s _i,d(t₁+1)=s_i,d(t₁)+v_i,d(t₁+1),  (8)

where d is the space dimension, d=1, 2, v_i,d(t₁) is the velocity of the ith particle in the dth dimension at time t₁, p_i,d(t₁) is the individual optimal solution of the ith particle in the dth dimension at time t₁, g_d(t₁) is the global optimal solutions of the ith particle at time t₁;

{circle around (4)} if SPSO reaches a preset maximum number of evolutions N₁, stop the iterative evolution process, transfer the value of S_NO* to the lower level; if SPSO does not reach the preset maximum number of evolutions N₁, return to step {circle around (3)};

{circle around (5)} introduce a multiobjective particle swarm optimization (MOPSO) algorithm to optimize the lower-level optimization problem, the position of the jth particle a_j(t₂) and the velocity of the jth particle b_j(t₂) can be represented as a_j(t₂)=[a_j,1(t₂), a_j,2(t₂), a_j,3(t₂), a_j,4(t₂)], a_i,1(t₂) represents the value of S_O at time t₂, a_i,2(t₂) represents the value of S_NH at time t₂, a_i,3(t₂) represents the value of SS at time t₂, a_i,4(t₂) represents the value of S_NO* at time t₂, b_j(t₂)=[b_j,1(t₂), b_j,2(t₂), b_j,3(t₂), b_j,4(t₂)], j is the number of particles, j=1, 2, . . . , 50; during the iterative evolution process, the obtained non-dominated solutions are conserved in the external archive Z(t₂), Z(t₂)=[z₁(t₂), z₂(t₂), . . . , z_j(t₂), . . . , z₅₀(t₂)], the update rule of the external archive is:

ž _j(t₂)=z_j(t₂)+0.09∇D(z_j(t₂)),  (9)

where z_j(t₂) is the jth non-dominated solution at time t₂ before the archive is updated, ž_j(t₂) is the jth non-dominated solution at time t₂ after the archive is updated, z_j(t₂)=[z_j,1(t₂), z_j,2(t₂)], ž_j(t₂)=[ž_j,1(t₂), ž_j,2(t₂)], z_j,1(t₂) and ž_j,1(t₂) are the values of S_O before and after the archive is updated, z_j,2(t₂) and ž_j,2(t₂) are the values of S_NO before and after the archive is updated, ∇D is the gradient descent direction;

{circle around (6)} establish a multi-input-multi-output radial basis assisted model (RBSM) based on the non-dominated solutions in Z(t₂):

$\begin{array}{cc}{B}_j\left(t_2\right)=\sum _{k=1}^8o_{j,k}\left(t_2\right)\times {\theta }_{j,k}\left(t_2\right),& \left(10\right)\end{array}$

where B_j(t₂) is an output vector of RBSM, B_j(t₂)=[B_j,1(t₂), B_j,2(t₂)]^T, B_j,1(t₂) is a predicted value of the aeration energy at time t₂, B_j,2(t₂) is a predicted value of the effluent quality at time t₂, o_j(t₂)=[o_j,1(t₂), o_j,2(t₂), . . . , o_j,8(t₂)]^T are connection weights, θ_j(t₂)=[θ_j,1(t₂), θ_j,2(t₂), . . . , θ_j,8(t₂)]^T is the output vector of the neurons in hidden layer, the sum of the squared errors between the output of RBSM and the actual system is expressed as

e(z_n(t₂))=min(B_n(t₂)−Q(t₂))^T(B_n(t₂)−Q(t₂)),  (11)

where e(z_n(t₂)) is the sum of the squared errors between the outputs of the nth non-dominated solution B_n(t₂) and the actual system Q(t₂), n∈[1, I₂], Q(t₂)=[Q₁(t₂), Q₂(t₂)] is the real outputs of the aeration energy and effluent quality in the actual system, select the solution corresponding to the minimal sum of the squared error as the global optimal solution;

{circle around (7)} if MOPSO reaches the preset maximum number of evolutions N₂, stop the iterative evolution process and output the optimal set-points of dissolved oxygen S_O*; if MOPSO does not reach the preset maximum number of evolutions N₂, return to step {circle around (5)};

(4) design the controller based on the predictive control strategy:

{circle around (1)} define the cost functions in the predictive control strategy:

$\begin{array}{cc}{J}_1\left(t\right)=\frac{1}{2}\sum _{q=1}^5{\left(z_1\left(t\right)-y_1\left(t\right)\right)}^2+\frac{1}{2}\sum _{m=1}^4\Delta {u\left(t\right)}^T\Delta u\left(t\right),J_2\left(t\right)=\frac{1}{2}\sum _{q=1}^5{\left(z_2\left(t\right)-y_2\left(t\right)\right)}^2+\frac{1}{2}\sum _{m=1}^4\Delta {u\left(t\right)}^T\Delta u\left(t\right),& \left(12\right)\end{array}$

where z₁(t) and z₂(t) are the optimal set-points of S_O* and S_NO*, y₁(t) and y₂(t) are the predicted values of S_O and S_NO;

{circle around (2)} update control laws based on the predictive control strategy, the updated rule is:

u(t+1)=u(t)+Δu(t),  (13)

where Δu(t) is the control law at time t, Δu(t) are control variations, whose expressions are shown as:

Δu₁(t)=½∂J₁(t)/∂u₁(t)+½∂J₂(t)/∂u₁(t),

Δu₂(t)=½∂J₁(t)/∂u₂(t)+½∂J₂(t)/∂u₂(t),  (14)

where Δu(t) are the variations of the manipulated variables oxygen transfer coefficient ΔK_La and internal recycle ΔQ_a, Δu(t)=[Δu₁(t), Δu₂(t)];

{circle around (3)} transfer the obtained ΔK_La and ΔQ_a to the actual system of WWTP;

(3) the inputs of the cooperative optimal control based wastewater treatment system is ΔK_La and ΔQ_a, the outputs of the wastewater treatment system is S_O and S_NO, the effect of the optimal control results is reflected by daily average of PE value, the daily average of AE value, the daily average of EQ value, and the tracking control results of S_O and S_NO.

The Novelties of this Present Disclosure Contain:

(1) Since WWTP is a complex and dynamic biochemical reaction process, conflicts with coupling performance indices, it is necessary to coordinate the performance indices in this present invention, and then improve the system performance, as well as guarantee the efficient and stable operation. However, due to the nonlinearity, strong coupling feature, and the conflicted multiple time-scale performance indices, it is difficult to realize the optimal control in WWTP. According to the operation features, a COCS method based on the cooperation optimization algorithm and predictive control strategy is proposed, which can balance the coupling performance indices and improve the operation efficiency.

(2) In this present invention, a COCS based on the cooperation optimization algorithm and predictive control strategy is studied to realize the optimal control in WWTP. This proposed optimal control method can efficiently capture the dynamic features of the performance indices, and acquire the optimal set-points of S_O and S_NO based on the dynamic working scenarios. Moreover, this proposed COCS can improve the optimization indices, and optimize the control performance.

Attention: for the convenient description, a data-driven assisted model based cooperative optimization algorithm and a predictive control strategy are used to describe and optimize the performance index, and trace the optimal set-points of S_O and S_No. The other optimal control methods based on different optimization algorithm and control strategy also belong to the scope of this present invention.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 shows the control system scheme of WWTP.

FIG. 2 shows the results of PE in COCS.

FIG. 3 shows the results of AE in COCS.

FIG. 4 shows the results of EQ in COCS.

FIG. 5 shows the results of S_O in COCS.

FIG. 6 shows the results of S_NO in COCS.

FIG. 7 is a chart showing the contribution degree index of process variables.

DETAILED DESCRIPTION

(1) The design of the control system framework for WWTP, the steps are:

1) the control system includes the serve, the PLC control cabinet, the blower, an ORP electrode sensor, a S_O sensor, a SS sensor, a temperature sensor, a pH sensor, a S_NH sensor, a TN sensor, a COD sensor, where

the serve is the operation center of control system to optimize the control procedure and realize the monitoring functions;

the PLC control cabinet is used to complete the automatic detection, monitor and control of the equipment status and operation conditions;

the blower acts as an executor of the control system to realize the aeration and ventilation of aerobic tank;

an ORP electrode sensor, located in the anaerobic tank, is configured to measure the values of ORP;

a S_O sensor, located in the aerobic tank, is configured to measure the values of S_O;

a SS sensor, located in the aerobic tank, is configured to measure the values of SS;

a temperature sensor, located in the effluent tank, is configured to measure the values of temperature;

a pH sensor, located in the in the effluent tank, is configured to measure the values of pH;

a S_NH sensor, located in the effluent tank, is configured to measure the values of S_NH;

a TN sensor, located in the effluent tank, is configured to measure the values of TN;

a COD sensor, located in the effluent tank, is configured to measure the values of COD;

2) data collection and transmission, at first, it is necessary to achieve the data acquisition through an ORP electrode sensor, a S_O sensor, a SS sensor, a temperature sensor, a pH sensor, a S_NH sensor, a TN sensor, a COD sensor; then it is essential to connect PLC and the sensors by using RS 485, realize the communication, complete the upload of the collected data, and transfer the collected data to servers through coordinated communication standards;

(2) the optimal controller, which comprises the establishment of the performance indices, the development of the optimization algorithm and the design of the controller, the steps are:

1) select the related process variables of PE: S_NO, MLSS, and choose the related process variables of AE and EQ: S_O, SS, S_NH, S_NO;

2) formulate the two-level models based on the different time scales, the upper-level model is for PE, and the lower-level models are for AE and EQ:

F ₁(t₁)=l₁(x_u(t₁)),  (1)

f ₁(t₂)=l₂(x_l(t₂),x_u*(t₁)),

f ₂(t₂)=l₃(x_l(t₂),x_u*(t₁)),  (2)

where F₁(t₁) is the model for PE at time t₁, l₁(x_u(t₁)) is the mapping function of PE model, f₁(t₂) is the model for AE at time t₂, l₂(x_l(t₂), x*_u(t₁)) is the mapping function of AE model, x*_u(t₁) is the optimal set-points of nitrate nitrogen S_NO* at time t₁, f₂(t₂) is the model for EQ at time t₂, l₃(x_l(t₂), x*_u(t₁)) is the mapping function of EQ model, x_u(t₁)=[S_NO(t₁), MLSS(t₁)] is the input variables vector of PE at time t₁, S_NO(t₁) is the concentration of S_NO at time t₁, MLSS(t₁) is the concentration of MLSS at time t₁, and the initial values of the two variables are [0.85, 1.56], x_l(t₂)=[S_O(t₂), SS(t₂), S_NH(t₂)], S_O(t₂) is the concentration of S_O at time t₂, SS(t₂) is the concentration of SS at time t₂, S_NH(t₂) is the concentration of S_NH at time t₂, and [S_O(t₂), SS(t₂), S_NH(t₂), S_NO*(t₁)] is the input variables vector of AE and EQ at time t₂, and the initial values are [1.9, 11.6, 3.8, 0.95];

(3) design a cooperative optimization algorithm to optimize upper-level and lower-level optimization problems to obtain the optimal set-points of the control variables, where the optimization period in the upper level is two hours, in the lower level is thirty minutes, the steps are:

{circle around (1)} formulate the upper-level and lower-level problems:

Min F₁(S_NO(t₁),MLSS(t₁)),  (3)

Min[f₁(S_O(t₂),S_NH(t₂),SS(t₂),S_NO*(t₁)),

f₂(S_O(t₂),S_NH(t₂),SS(t₂),S_NO*(t₁))],  (4)

where Min F₁(S_NO(t₁), MLSS(t₁)) is the upper-level optimization problem, Min [f₁(S_O(t₂), S_NH(t₂), SS(t₂), S_NO*(t₁)), f₂(S_O(t₂), S_NH(t₂), SS(t₂), S_NO*(t₁))] is the lower-level optimization problem;

{circle around (2)} set the number of the particle population in the upper level optimization I₁, the number of the particle population in the lower level optimization I₂, the maximum number of iterations in the upper level optimization N₁, and the maximum number of iterations in the lower level optimization N₂, I₁=50, I₂=50, N₁=20, N₂=50;

{circle around (3)} introduce the single particle swarm optimization (SPSO) algorithm to optimize the upper-level optimization problem, the position and the velocity of the ith particle can be shown as:

s _i(t₁)=[s_i,1(t₁),s_i,2(t₁)],  (5)

v _i(t₁)=[v_i,1(t₁),v_i,2(t₁)],  (6)

s_i(t₁) is the position of the ith particle at time t₁, s_i,1(t₁) is the value of S_NO at time t₁, s_i,2(t₁) is the value of MLSS at time t₁, v_i(t₁) is the velocity of the ith particle at time t₁, i is the number of particles, i=1, 2, . . . , 50, an update process of s_i(t₁) and v_i(t₁) are

v _i,d(t₁+1)=0.7v_i,d(t₁)+0.72(p_i,d(t₁)−s_i,d(t₁))+0.72₂(g_d(t₁)−s_i,d(t₁)),  (7)

s _i,d(t₁+1)=s_i,d(t₁)+v_i,d(t₁+1),  (8)

where d is the space dimension, d=1, 2, v_i,d(t₁) is the velocity of the ith particle in the dth dimension at time t₁, p_i,d(t₁) is the individual optimal solution of the ith particle in the dth dimension at time t₁, g_d(t₁) is the global optimal solutions of the ith particle at time t₁;

{circle around (4)} if SPSO reaches the preset maximum number of evolutions N₁, stop the iterative evolution process, transfer the value of S_NO* to the lower level; if SPSO does not reach the preset maximum number of evolutions N₁, return to step {circle around (3)};

{circle around (5)} introduce the multiobjective particle swarm optimization (MOPSO) algorithm to optimize the lower-level optimization problem, the position of the jth particle a_j(t₂) and the velocity of the jth particle b_j(t₂) can be represented as a_j(t₂)=[a_j,1(t₂), a_j,2(t₂), a_j,3(t₂), a_j,4(t₂)], a_i,1(t₂) represents the value of S_O at time t₂, a_i,2(t₂) represents the value of S_NH at time t₂, a_i,3(t₂) represents the value of SS at time t₂, a_i,4(t₂) represents the value of S_NO* at time t₂, b_j(t₂)=[b_j,1(t₂), b_j,2(t₂), b_j,3(t₂), b_j,4(t₂)], j is the number of particles, j=1, 2, . . . , 50; during the iterative evolution process, the obtained non-dominated solutions are conserved in the external archive Z(t₂), Z(t₂)=[z₁(t₂), z₂(t₂), . . . , z_j(t₂), . . . , z₅₀(t₂)], the update rule of the external archive is:

ž _j(t₂)=z_j(t₂)+0.09∇D(z_j(t₂)),  (9)

where z_j(t₂) is the jth non-dominated solution at time t₂ before the archive is updated, ž_j(t₂) is the jth non-dominated solution at time t₂ after the archive is updated, z_j(t₂)=[z_j,1(t₂), z_j,2(t₂)], ž_j(t₂)=[ž_j,1(t₂), ž_j,2(t₂)], z_j,1(t₂) and ž_j,1(t₂) are the values of S_O before and after the archive is updated, z_j,2(t₂) and ž_j,2(t₂) are the values of S_NO before and after the archive is updated, ∇D is the gradient descent direction;

{circle around (6)} establish the multi-input-multi-output radial basis assisted model (RBSM) based on the non-dominated solutions in Z(t₂):

$\begin{array}{cc}{B}_j\left(t_2\right)=\sum _{k=1}^8o_{j,k}\left(t_2\right)\times {\theta }_{j,k}\left(t_2\right),& \left(10\right)\end{array}$

where B_j(t₂) is the output vector of RBSM, B_j(t₂)=[B_j,1(t₂), B_j,2(t₂)]^T, B_j,1(t₂) is the predicted value of the aeration energy at time t₂, B_j,2(t₂) is the predicted value of the effluent quality at time t₂, o_j(t₂)=[o_j,1(t₂), o_j,2(t₂), . . . , o_j,8(t₂)]^T are the connection weights, θ_j(t₂)=[θ_j,1(t₂), θ_j,2(t₂), . . . , θ_j,8(t₂)]^T is the output vector of the neurons in hidden layer, the sum of the squared errors between the output of RBSM and the actual system is expressed as

e(z_n(t₂))=min(B_n(t₂)−Q(t₂))^T(B_n(t₂)−Q(t₂)),  (11)

where e(z_n(t₂)) is the sum of the squared errors between the outputs of the nth non-dominated solution B_n(t₂) and the actual system Q(t₂), n∈[1, 50], Q(t₂)=[Q₁(t₂), Q₂(t₂)] is the real outputs of AE and EQ in the actual system, select the solution corresponding to the minimal sum of the squared error as the global optimal solution;

{circle around (7)} if MOPSO reaches the preset maximum number of evolutions N₂, stop the iterative evolution process and output the optimal set-points of dissolved oxygen S_O*; if MOPSO does not reach the preset maximum number of evolutions N₂, return to step {circle around (5)};

(4) design the tracking control method based on the predictive control strategy:

{circle around (1)} define the cost functions in the predictive control strategy:

$\begin{array}{cc}{J}_1\left(t\right)=\frac{1}{2}\sum _{q=1}^5{\left(z_1\left(t\right)-y_1\left(t\right)\right)}^2+\frac{1}{2}\sum _{m=1}^4\Delta {u\left(t\right)}^T\Delta u\left(t\right),J_2\left(t\right)=\frac{1}{2}\sum _{q=1}^5{\left(z_2\left(t\right)-y_2\left(t\right)\right)}^2+\frac{1}{2}\sum _{m=1}^4\Delta {u\left(t\right)}^T\Delta u\left(t\right),& \left(12\right)\end{array}$

where z₁(t) and z₂(t) are the optimal set-points of S_O* and S_NO*, y₁(t) and y₂(t) are the predicted values of S_O and S_NO;

{circle around (2)} update the control laws based on the predictive control strategy, the updated rule is:

u(t+1)=u(t)+Δu(t),  (13)

where Δu(t) is the control law at time t, Δu(t) are the control variations, whose expressions are shown as:

$\begin{array}{cc}\Delta u_1\left(t\right)=\frac{1}{2}\partial J_1\left(t\right)/\partial u_1\left(t\right)+\frac{1}{2}\partial J_2\left(t\right)/\partial u_1\left(t\right),\Delta u_2\left(t\right)=\frac{1}{2}\partial J_1\left(t\right)/\partial u_2\left(t\right)+\frac{1}{2}\partial J_2\left(t\right)/\partial u_2\left(t\right),& \left(14\right)\end{array}$

where Δu(t) are the variations of the manipulated variables oxygen transfer coefficient ΔK_La and internal recycle ΔQ_a, Δu(t)=[Δu₁(t), Δu₂(t)];

{circle around (3)} transfer the obtained ΔK_La and ΔQ_a to the actual system of WWTP;

(3) the inputs of the cooperative optimal control based wastewater treatment system is ΔK_La and ΔQ_a, the outputs of the wastewater treatment system is S_O and S_NO, the effect of the optimal control results is reflected by the daily average of PE value, the daily average of AE value, the daily average of EQ value, and the tracking control results of S_O and S_NO.

The control system scheme based on COCS is shown in FIG. 1, and the optimal control results are presented in FIGS. 2-6. FIG. 2 gives the daily average PE values, X axis shows the time, and the unit is day, Y axis is the average PE value, and the unit is €/KWh. FIG. 3 gives the daily average AE values, X axis shows the time, and the unit is day, Y axis is the average AE value, and the unit is €/KWh. FIG. 4 gives the daily average EQ values, X axis shows the time, and the unit is day, Y axis is the average EQ value, and the unit is €/m³. FIG. 5(a) gives S_O values, X axis shows the time, and the unit is day, Y axis is control results of S_O, and the unit is mg/L. FIG. 5(b) gives the control errors of S_O, X axis shows the time, and the unit is day, Y axis is control errors of S_O, and the unit is mg/L. FIG. 6(a) gives the S_NO values, X axis shows the time, and the unit is day, Y axis is control results of S_NO, and the unit is mg/L. FIG. 6(b) gives the control errors of S_NO, X axis shows the time, and the unit is day, Y axis is control errors of S_NO, and the unit is mg/L.

Claims

1. A cooperative optimal control method for wastewater treatment process (WWTP), comprising the following steps: (1) select process variables of pumping energy, which include nitrate nitrogen (SNO), mixed liquor suspended solids (MLSS), and choose process variables of aeration energy and the effluent quality, which include dissolved oxygen (SO), suspended solids (SS), ammonia nitrogen (SNH), SNO;(2) formulate two-level models based on different time scales, wherein an upper-level model is designed for the pumping energy, and lower-level models are designed for the aeration energy and effluent quality: F1(t1)=l1(xu(t1)),  (1)f1(t2)=l2(xl(t2),xu*(t1))f2(t2)=l3(xl(t2),xu*(t1)),  (2)