Abstract
In this work, an efficient strategy was proposed for efficient solution of the dynamic model of SWRO system. Since the dynamic model is formulated by a set of differential-algebraic equations, simultaneous strategies based on collocations on finite element were used to transform the DAOP into large scale nonlinear programming problem named Opt2. Then, simulation of RO process and storage tanks was carried element by element and step by step with fixed control variables. All the obtained values of these variables then were used as the initial value for the optimal solution of SWRO system. Finally, in order to accelerate the computing efficiency and at the same time to keep enough accuracy for the solution of Opt2, a simple but efficient finite element refinement rule was used to reduce the scale of Opt2. The proposed strategy was applied to a large scale SWRO system with 8 RO plants and 4 storage tanks as case study. Computing result shows that the proposed strategy is quite effective for optimal operation of the large scale SWRO system; the optimal problem can be successfully solved within decades of iterations and several minutes when load and other operating parameters fluctuate.
1. Introduction
The shortage of freshwater resources is expected to worsen with the growth of population and industrialization, as well as climate change [1, 2]. Seawater desalination is one of the most promising approaches to get freshwater resources in the world and is considered as the most important strategy to develop new freshwater for the coastal counties [3, 4]. With the process of low-cost, high rejection membrane technique and high efficiency energy recovery device, seawater reverse osmosis (SWRO) technique is becoming the most popular and attractive seawater technique for its economy and convenience [5–7].
Recently, more attention was paid to optimal operation and energy management to further reduce the energy consumption of SWRO system [8]. Modeling methods based on first principle and data-driven are used for practical control and fault diagnosis [9]. And with development of advanced control techniques such as predictive control, sliding-mode control, and optimal control [10–13], more potential benefit can be expected.
Kim et al. comprehensively studied the RO process especially on the minimization of product cost with system engineering method [14]. Sassi and Mujtaba and Palacin et al. evaluated the optimal operation of SWRO system through minimizing the specific energy consumption [15, 16]; with the consideration of thermodynamic restriction, Zhu studied the energy cost optimization problem. Among which factors such as stages and energy recovery efficiency were discussed to get the optimal operation condition [17]. After a comprehensive first-principle based mathematical model has been developed and validated by plant data, Li studied the optimal plant operation of brackish water reverse osmosis (BWRO) desalination to reduce specific energy consumption (SEC). His computing results show that about 16% reduction of SEC can be achieved by optimizing operating condition [18]. Since all the cost spent on SWRO system includes not only operational cost but also capital cost, Geraldes studied the optimization problem through minimizing the objective function including all these costs; but to make solution of the problem with DAEs easier and simple, Geraldes used differential technology to discretize the equations and then obtained the optimal results with low accuracy [19]. Since the operating environments change frequently, Sassi and Mujtaba studied the optimal operation problem which is the variation of load and feed temperature. In his study, storage tank was used as the buffering unit between water production and water supply, and water level of storage tank can be freely adjusted to add the flexibility of operation [20].
The researches above are of great significance to improve the economical performance of desalination, but most of them put emphasis on the optimal operation of steady-state process. As we know, in the actual process many operating parameters such as feed temperature, electricity price, and water supply requirement change dynamically over time [21]; the performance of equipment and reverse osmosis membrane also constantly change over time, so dynamic optimization to improve the operation will be more meaningful. However, since the dynamical optimization problem includes a set of differential and algebraic equations (DAEs), its efficient and stable solution is fairly tough when inequality constraints enforced on the bounds. Direct methods such as variational principle based method often failed because of their disadvantage in dealing with active inequalities [22, 23].
Though with development of computing technique, simultaneous approach with collocation on finite element is more suitable for this kind of problems, and large scale solvers such as sparse SQP and IPOPT become more powerful to solve large scale nonlinear problems, there are still many works to do to make the solution more efficient and stable, especially when the optimal problem is not convex and with the characteristic of strong nonlinear. But the real time optimization need the optimal operation problem be solved efficiently. In this paper, to solve optimization problem of RO process system with network structure more efficiently, based on the well-developed dynamic models of RO process and water storage process, an efficient optimization strategy to optimize the dynamic operation of SWRO system is proposed, which will be helpful to the reduction of energy cost and to the realization of real time optimization.
2. Dynamic Model of SWRO System
According to the flowsheet of SWRO system, RO unit is the key part for freshwater production. After the pretreatment, the seawater was pumped into RO modules by high pressure pump, from which the freshwater and the salt were separated through the work of solution-diffusion; then with simple posttreatment, the freshwater was pumped to terminal user [19] to meet the needs of daily life. In the RO module, the pressure, flow rate, and concentration change along the channel (shown in Figure 1). These variables satisfy the following equations [19]: here, denotes solvent flux, and represent the axial velocity in feed channel. is the bulk concentration of RO module, and is the pressure drop along the RO module. The parameters and variables from the above equations can be obtained from the solution-diffusion mass transport relations and principle of energy conservation; all the equations involved are listed as follows [19, 24–27].
RO process model equations are
See [24]
See [19]
See [27]
Consider (Masaaki Sekino, 1995)
As to the large scale SWRO system with multiple RO units and multiple storage tanks, the storage tanks perform as the buffering units between water production of RO units and water supply. Since there are pipes to connect each RO unit with storage tanks, the feed water can come from each RO unit; the dynamic process can be formulated as and equations Bound constraint is as follows: Initial condition is as follows: Here denotes the number of RO units, denotes the number of storage tanks, denotes the flow rate from the th RO units to the th storage tanks, denotes the load requirement of freshwater of the th storage tank, denotes the total flow rate of freshwater supply, and the denotes the concentration of the supplied freshwater.
3. Formulation of Dynamic Optimization Problem of SWRO System
Energy cost is the largest portion of the SWRO system’s operation cost; generally it accounts for over 50% of the operation cost. So minimizing the energy cost can be set as the objective function. It is affected by the electricity price and the specific energy consumption (SEC) which is defined as [28]
With the consideration of daily fluctuation of feed temperature, electricity price, and the load requirement, the daily objective function can be denoted by
Since the control variables cannot be adjusted too frequently and the change of electricity price is based on an hour, it is reasonable to adjust the control variables one hour by one hour, and so (12) can be renewed in another form, and the optimization problem named Opt1 can be denoted by
Objective Function.
Equality Constraints. RO process model ((1)–(6)). Dynamic model of storage tanks ((7)-(8)). Specific energy cost (11).
Inequality Constraints. Water quality: . Scaling resistant: . Concentration polarization: . Pressure constraints: . Flowrate of RO units: . Water level . Initial conditions of RO at each hour are as follows: Initial condition of storage tanks , (10).
4. The Strategy to Solve Opt2
4.1. Simultaneous Approach Based on Collocation on Finite Element
Since the Opt1 includes differential equations as well as strongly nonlinear algebraic equations, it belongs to a kind of differential-algebraic optimization problems (DAOPs). The Opt1 problem of SWRO system with RO plants and storage tanks will have () differential equations and more nonlinear algebraic equations, since there are also constraints to satisfy the water quality and the equipment’s safety; it is really hard to get the solution efficiently through direct ways such as variational principle. Here we use the simultaneous approach to fully discretize all the control variables and status variables, which is quite suitable and has many advantages for this kind of problems. We use the following monomial basis representation for the differential profiles, which is popular for Runge-Kutta discretization [29]:
Here is the value of the differential variable at the beginning of element , is the length of element , denotes the value of its first derivative in element at the collocation point , and is a polynomial of order , satisfying here is the location of the th collocation point within each element. Continuity of the differential profile is enforced by
In addition, the control and algebraic profile are approximated using a Lagrange basis representation which takes the form: Here and represent the values of the algebraic and control variables, respectively, in element at collocation point , is the value satisfying , and is the Lagrange polynomial of degree .
Through the simultaneous method mentioned above, the Opt1 in the form of DAEs is discretized into the NLP problem denoted by (19a)–(19i); here the problem is named Opt2.
Here , .
Since all the differential variables and continuous variables are fully discretized into algebraic variables along the spacious or time horizon, and large number of finite elements is generally selected to ensure the accuracy of the discretization, all the factors led to the great increase of the size of discretized problem of Opt2. As the result, though with the development of computing technique such as spare matrix and automatic differentiation, large scale solvers such as IPOPT and CONOPT can solve problems with millions of variables, and its efficient and stable solution is a big challenge. Special work need to be paid attention to to reduce the computing time and to keep the stability of the problem.
There are still two key factors influencing the efficient solution of the discretized problem of Opt2: the first is to configure the number of finite element. Too large number of finite elements will lead to rapid increase in the size of variables and equations, which will greatly increase the computing time and memory occupation. The second problem is to obtain initial values for all the discretized variables, which is very critical to the stability of Opt2 formulated in (19a)–(19i). If the initial values are not carefully given, time-consuming and even failure of convergence will happen during the solution of Opt2. Based on the above reasons, techniques for assignment of initial values and suitable division of finite element are carefully designed as below.
4.2. Assignment of Initial Value through Simulation on Finite Element Individually
If all the control variables are fixed as constant values, all the variables of Opt2 can be obtained by simulation of the whole process model. According to the flowsheet of SWRO system, the permeate flow rate and concentration of the permeate water are the input parameters for the dynamic process of storage tanks, so variables of storage tanks can only be obtained after the simulation of RO process.
The principle of collocation on finite element is shown in Figure 2. Though all the variables are discretized, to keep the continuity, the final values of the first finite element are equal to the initial value of the second finite element. So to get all the values of RO process, the control variables are fixed as constant values, and then the simulation is carried out from one finite element to the next along the space or time horizon. After the simulation of RO process was finished, the values of , , and so on were used as input parameters and the simulations of storage tanks are then carried out in the same way of RO process.
After all the values of the discretized variables of Opt2 were obtained, we relax the control variables and enforce the constraints again. Then those obtained values were set as the initial value for Opt2, and nonlinear solver of IPOPT is used for the solution. The whole process of the technique to get good initial values is shown in Figure 3.
(a) Simulation of RO process
(b) Simulation of storage tanks
(c) The whole process of the strategy
4.3. Method for Suitable Division of Finite Element
To balance the relationship of the discrete accuracy and computing efficiency, suitable ways to determine the number of finite element are necessary. Truong et al. improved an adaptive mesh refinement method [30] to improve the computing efficiency; Lan and Taylor designed a moving grid computing strategy [31], though the computational efficiency is not high, the method is more stable. Binder proposed a mesh refinement strategy, making it superior to the average grid distribution strategy [32]. Biegler et al. proposed a dynamic optimization method [33–35] based on moving finite element technique; the finite element is automatically adjusted through error analysis. In this work, with the consideration of the special structure of SWRO system, and the relationship between the change of objective function value and the number of finite element which affect discrete state and control variables, a simple but effective method to get suitable division of finite element is designed as follows.(1)Set the accuracy of the objective function which is defined as ; is a small positive value such as .(2)Set the number of finite element as Nef1 and the corresponding collocation point by conventional way.(3)Solve the discretized problem of Opt2, and get the value of objective function .(4)Add the number of finite element to Nef1′.(5)Resolve the discretized problem of Opt2, and then obtain another value of objective function .(6)Compute the value of .(7)If the required accuracy is satisfied, then stop; otherwise, let , and go to step .
Since IPOPT is one of the best solvers to solve large scale dynamical problems, here it was used for the solution of Opt2 quickly.
5. Case Study
The large scale SWRO system has 8 RO plants and 4 storage tanks, and the RO membranes are the products of Dow Chemical Company named SW30HR. Each RO plant has decades of pressure vessels in parallel, and several RO membranes in the pressure vessel connected in series are used to separate salt and pure water. The storage tanks are used as the buffering units between water supply and water production. The electricity price and other information for the system are listed in Tables 1-2 and shown in Figure 4. Here NPT rays and NMOD represent the number of pressure vessel and the number of membranes in pressure vessels. FF denotes the fouling factor, and and denote the efficiency of high-pressure pumps and energy recovery unit.
The optimization problem of this system was discretized by simultaneous approach based on collocation on finite element. Since the control variables are not allowed to adjust frequently and the electricity price changes at different hour, 24 finite elements and 3 collocation points are selected to discretize the dynamic variables of storage tanks. For the RO process, 3 collocation points are used and the number of finite element is determined by the proposed technique. To get good initial values, the developed method in chapter 3.2 was used for the solution of Opt2. To comprehensively understand the effect of finite element on the solution, information of problem Opt2 with different finite element of RO process is listed in Table 3. Table 3 shows that when the number of finite element of RO increased to 40, the discretized variables and the equations of Opt2 reached 495371 and 494411, respectively; the size of the problem is quite large.
Nonlinear solver IPOPT [36, 37] which is based on interior point method was used for the solution of Opt2, and computing results with different finite element of RO process were listed in Table 4. The results show that the problem of Opt2 can be successfully solved when the finite element increases from 2 to 40, and the solution converged to the best point within decades of iterations. The computing time of simulation and optimization increased much more greatly when the finite element increases. In particular, it can be seen that the optimization time increases almost in the form of exponential. When 40 finite elements and 3 collocation points were selected, the optimization time will reach 3519.953 CPU seconds.
It can also be seen that the change of objective function became much smaller with the increase of finite elements. Compared with the objective function of 40 finite elements, the relative errors were obtained and showed in Figure 5. As the finite elements increased to 10, the relative errors become fairly small, and the required accuracy can be reached at this time. Compared with the computing time with 40 finite elements, more than 94% computing time was saved. To show the effectiveness of the above method, initial values for Opt2 were given by reasonable assumption near the set point, the optimal solution failed to converge when the finite element is 40; and even when the finite element is 10, more than 1187 iterations and more than 4680 seconds were spent to get the optimal solution.
As 10 finite elements of RO unit were selected, the solution of problem of Opt2 was resolved under default operating conditions. The solution led to the profiles of control variables-feed flow rates and feed pressures along time horizon, from which the best control trajectory as well as the various status variables are achieved. Some of these results can be seen from Figures 6, 7, 8, and 9.
Since the feed temperature changes frequently and has significant effect on performance of RO membranes, it is necessary to optimize the process quickly through adjusting the control variables to the aimed value. To validate the performance of the proposed strategy with large fluctuation of feed temperature, increase the feed temperature from 21°C to 29°C and keep other parameters as default, and then solve the Opt2 individually. Summary of the computing results can be seen from Table 5. From which it can be seen that temperature has significant effect on the energy cost. The optimal energy cost reduces quickly as the temperature increases, when feed temperature increases from 21°C to 29°C, the optimal energy cost reduces from 3.8654 Wan CNY to 3.4202 Wan CNY. It can also be seen that the change of feed temperature has little effect on the computing performance and convergence of the proposed method. All the solution can be successfully achieved within decades of iterations and within several minutes. The results show that the proposed strategy for the SWRO system not only is quite effective but also has wide temperature adaptability. All the results show the good potential to the on-line optimization of actual SWRO system.
6. Conclusions
It is of significance to reduce the energy cost of seawater reverse osmosis (SWRO) system through the system engineering approach, especially when the system is operated in dynamical environments due to the changes of load and other operating parameters.
In this work, a dynamic optimization strategy for optimizing operation of large scale SWRO system was proposed. Based on first principle, the dynamic models of RO units and storage tanks which are developed by were formulated by differential-algebraic equations, and the objective function is formulated to minimize the energy cost. The optimization problem in the form of DAOPs was fully discretized into large scale NLP through simultaneous method, with which all the dynamic variables were approximated by polynomial equations with collocation point on each finite element. To solve the large scale NLP problem named Opt2 efficiently, a method was developed to get the initial value of all variables through simulation element by element and step by step. Then a simple but effective rule to determine the number of finite element was used to reduce the scale and at the same time to keep enough accuracy of Opt2.
Case study of a large scale SWRO system which includes 8 RO plants and 4 storage tanks was studied to validate the proposed strategy. Computing results demonstrate that the method is quite efficient for the solution of dynamic optimization of SWRO system. The problem can be solved in decade’s iterations within several minutes. The minimized objective function and accordingly the best daily profiles of control variables can be achieved. The proposed strategy also performs well with large fluctuation of feed temperature, which can fully satisfy the requirement of on-line real time optimization.
Appendix
| : | Membrane water permeability (m·s−1·Pa−1) |
| : | Intrinsic membrane water permeability (m·s−1·Pa−1) |
| : | Membrane TDS permeability (m/s) |
| : | Intrinsic membrane TDS permeability (m/s) |
| : | Bulk concentration along feed channel (kg/m3) |
| : | Feed salt concentration (kg/m3) |
| : | Salt concentration of membrane surface (kg/m3) |
| : | Permeate concentration of th RO unit (kg/m3) |
| : | Brine concentration (kg/m3) |
| : | Salt concentration of output freshwater from th tanks (kg/m3) |
| : | Dynamic viscosity (m2/s) |
| : | Hydraulic diameter of the feed spacer channel (m) |
| : | Electricity price at different hours (CNY/kw·h) |
| : | Water level of th storage tank (m) |
| : | The length of the th finite element |
| : | Height of the feed spacer channel |
| : | Solvent flux |
| : | Solute flux |
| : | Empirical parameter |
| : | Mass transfer coefficient (m/s) |
| : | Length of the RO module (m) |
| : | The number of RO units |
| : | The number of storage tanks |
| : | Number of leaf in RO module |
| : | Pressure along feed channel (bar) |
| : | Pressure drop along RO spiral wound module (bar) |
| : | Feed pressure (bar) |
| : | Pressure in permeate side (bar) |
| : | Brine pressure (bar) |
| : | Feed flow rate (m3/h) |
| : | Output flow rate of the th tank (m3/h) |
| : | Permeate flow rate from th RO unit (m3/h) |
| : | Flow rate from th RO unit to the th tank (m3/h) |
| : | Brine flow rate (m3/h) |
| : | Total water supply flowrate (m3/h) |
| : | Gas law constant |
| Re: | Reynolds number (dimensionless) |
| : | Water recovery ratio (%) |
| Sc: | Schmidt number (dimensionless) |
| SEC: | Specific energy consumption (kw·h/m3) |
| Sh: | Sherwood number (dimensionless) |
| Sp: | Salt passage coefficient (%) |
| : | Area of storage tank (m2) |
| : | Operational temperature (K) |
| : | Axial velocity in feed channel (m/s) |
| : | Input axial velocity in membrane (m/s) |
| : | Output axial velocity in membrane (m/s) |
| : | Width of the RO module (m) |
| : | Salt rejection coefficient (%) |
| : | Algebraic variables in element at collocation point |
| : | Control variables in element at collocation point |
| : | Differential variable |
| : | Differential variable at the beginning of element . |
| : | Constant parameters |
| : | Pressure loss of osmosis pressure (bar) |
| : | Density of permeate water (kg/m3) |
| : | Kinematic viscosity (kg m−1 s−1) |
| : | Friction factor |
| : | Mechanical efficiency of high pressure pump |
| : | Energy recovery efficiency |
| : | Polynomial of order |
| : | Lagrange polynomial of degree . |
| : | Bulk |
| : | Module feed channel |
| : | Membrane surface |
| : | Permeate side |
| : | Salt. |
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This work was supported by National Natural Science foundation of China (nos. 60904058 and 61374142), National Key Technology Research and Development Program (2009BAB47B06), and Natural Science Foundation of Zhejiang (LY12F03001). The authors would like to extend their sincere thanks to Mr. Yang Chen from Carnegie Mellon University for many helpful discussions.