A Monthly Hydropower Scheduling Model of Cascaded Reservoirs with the Zoutendijk Method

Abstract

A monthly hydropower scheduling determines the monthly flows, storage, and power generation of each reservoir/hydropower plant over a planning horizon to maximize the total revenue or minimize the total operational cost. The problem is typically a complex and nonlinear optimization that involves equality and inequality constraints including the water balance, hydraulic coupling between cascaded hydropower plants, bounds on the reservoir storage, etc. This work applied the Zoutendijk algorithm for the first time to a medium/long-term hydropower scheduling of cascaded reservoirs, where the generating discharge capacity is handled with an iterative procedure, while the other head-related nonlinear constraints are represented with exponential functions fitting to discrete points. The procedure starts at an initial feasible solution, from which it finds a feasible improving direction, along which a better feasible solution is sought with a one-dimensional search. The results demonstrate that the Zoutendijk algorithm, when applied to six cascaded hydropower reservoirs on the Lancang River, worked very well in maximizing the hydropower production while ensuring the highest firm power output to be secured.

1. Introduction

Hydropower is a major green and renewable energy resource that can be massively developed and utilized at reasonable costs, and is also favorable in countries that are rich in it. An urgent deep decarbonation to ensure the sustainable development of economy and society places hydroelectricity in a critical role in energy structure [1,2]. In particular, the generation efficiency and water utilization rate can be substantially improved by a reasonable hydropower scheduling of cascaded reservoirs, which, however, is typically a complicated, high-dimensional, and multi-objective dynamic optimization problem that involves a nonlinear dynamic head, stochastic inflow, and hydraulic coupling between upstream and downstream hydropower plants [3,4]. Thus, it is essential to find an effective and efficient method or procedure to address this hydropower scheduling problem, which aims to maximize comprehensive benefits under some boundary conditions and other constraints.

The hydropower scheduling problem has long been studied in previous works, which experimented with many optimization algorithms including mathematical programming. The linear programming (LP), for instance, was employed by Feng et al. [5] who presented a weekly hydropower scheduling problem of cascaded hydropower plants that was formulated into a mixed integer linear programming (MILP) model. Formulating a nonlinear optimization into an LP problem, however, is a skillful task, requiring the nonlinear functions to be linearized, as illustrated by Zheng et al. [6], who presented a new three-triangle based method to linearly concave the hydropower output function (HOF), or in many cases, many integer variables to be introduced, causing the “combinatorial explosion” that makes it unlikely to derive a satisfactory solution in reasonable time. Nonlinear programming (NLP) is also a popular option, as applied by Barros et al. [7], who formulated hydropower generation into a nonlinear function of the water discharge and head, but it could not guarantee a global optimum. The dynamic programming (DP) is particularly applicable to a stage-by-stage decision making process, just like in the reservoir operation, but was only limited to small scale hydro systems due to the well-known “dimensional difficulty”, which needs to be alleviated by improving the solution efficiency, as shown by Ji et al. [8], who skillfully fitted the functional structure into the DP formulation to reduce the computation time.

Heuristic programming (HP), also known as intelligent algorithms, is another category of optimization algorithm that has been extensively applied in the hydropower scheduling of cascaded reservoirs. Wang et al. [9], for instance, applied particle swarm optimization (PSO) to a multi-reservoir system to optimize the distribution of water storage. The evolutionary algorithms (EA), which starts its search from a randomly generated initial population to attain the global optimal solution, were explored by Arunkumar and Jothiprakash [10], who introduced the chaotic technique to initiate a population randomly and also proposed a search strategy to enhance the performance of the EA in optimizing a multi-reservoir system. Bozorg et al. [11] investigated the efficacy of the gravity search algorithm (GSA) based on the gravity law and mass interactions, and applied it to some benchmark functional objectives and operation problems involving a single and four reservoirs. Shang et al. [12] made efforts to make the genetic algorithm (GA) more efficient when applied to maximize the hydropower production by proposing a proportional reproduction selection-based operator and a steady-state reproduction selection-based operator. HP is adaptable to a wide range of optimization problems, but very often with a volatile performance in securing a satisfactory solution in reasonable time.

The reservoir scheduling involves a variety of objectives including the power generation, flood control benefit, and ecological benefits. Diego et al. [13], for instance, performed a comparative assessment on revenues from the short-term optimizations driven by two econometric models to enhance the expected revenue in hydropower systems, especially contributed by those with large storage capacity. Daniel et al. [14] formulated a deterministic and a stochastic mathematical model to maximize the profits of selling energy produced by cascaded hydropower plants in a deregulated market. Yang et al. [15] presented a multi-reservoir flood control model that incorporated operating rules to alleviate the risk incurred by a traditional flood control model. Wang et al. [16] established a multi-objective optimization model to maximize the power generation while considering an ecological base flow, solved with a multi group cooperation operation particle swarm optimization (MGCL-PSO) model based on population cooperation. Wei et al. [17] proposed a multi-objective optimization method and a multi-attribute decision making strategy to balance the objectives involving flood control, power generation, water supply, and ecology in reservoir operation. An optimal flood control operation model of cascaded reservoirs was proposed by Zhou et al. [18] to coordinate the upstream reservoirs to reduce the downstream inflows at some flood control locations.

Table 1. Summary of abbreviations.

Abbreviations	Corresponding Meaning
LP	Linear programming
MILP	Mixed integer linear programming
HOF	Hydropower output function
NLP	Nonlinear programming
DP	Dynamic programming
HP	Heuristic programming
PSO	Particle swarm optimization
EA	Evolutionary algorithms
GSA	Gravity search algorithm
GA	Genetic algorithm
MGCL-PSO	Multi group cooperation operation particle swarm optimization

shows the summary of abbreviations.

As one of the feasible direction methods, the Zoutendijk method, which is applicable to nonlinear programming problems with linear constraints, starts from an initial feasible solution, around which a descending feasible direction will be found to explore a better solution by searching along this direction [19]. The Zoutendijk method has been applied to numerous optimization problems, showing very good performance in securing the global optimum and ensuring a high speed of convergence [20]. The method, however, has not been reported yet, and was applied to a mid/long-term optimal scheduling of cascaded hydropower plants [21].

This work will demonstrate how the Zoutendijk method can be applied to a monthly hydropower scheduling problem of six cascaded reservoirs located on the Lancang River, which involves nonlinear functions including the hydropower output and the capacity to generate discharge, which will be handled differently since the capacity of generating discharge is much less sensitive to the variation in the water head.

The following work is organized as follows. Section 2 presents the monthly hydropower scheduling formulation with the objective to maximize the total hydropower production while ensuring the highest firm power output, Section 3 shows how to convert this scheduling model into one that Zoutendijk is ready to solve, Section 4 validates the model and its methods with a case study, and Section 5 concludes the work.

2. Problem Formulation

Traditionally, a monthly hydropower scheduling takes into account the nonlinear dynamic head, hydraulic coupling, generating capacity, etc., with the objective to maximize the firm power output first and then the energy production during a planning horizon. Since the two objectives are different in priority, they are weighted differently in magnitude to convert the two objectives into an objective, expressed mathematically as

$$\max {\mathrm{W}}_1\cdot F+{\mathrm{W}}_2\cdot {\displaystyle \sum _{i=1}^N{\displaystyle \sum _{t=0}^{T-1}{P}_{it}}}$$

(1)

where W₁ and W₂ are the weights with ${\mathrm{w}}_1\gg {\mathrm{w}}_2$ to prioritize the firm power output (F) in MW over the total power output of all the hydropower plants; i and t are the subscripts for the reservoir and time-step, respectively; P_it in MW is the power output in MW in time-step t.

Constraints include:

2

The water balance,

$$V_{i,t+1}=V_{it}+({\displaystyle \sum _{j\in \Omega (i)}{Q}_{jt}}+I_{it}-Q_{it})\cdot \frac{\Delta t\cdot 24\times 3600}{1000,000}$$

(2)

with

$$\{\begin{array}{l}{V}_{i,0}=V_i^{\mathrm{ini}}\\ V_{i,T}=V_i^{\mathrm{end}}\\ sp{l}_{it}+q_{it}=Q_{it}\end{array}$$

(3)

where V_it = the storage in hm³ at the beginning of time-step t of reservoir i; Q_it = the outflow in m³/s in time-step t from reservoir i; I_it = the local inflow in m³/s in time-step t into reservoir i; $\Omega (i)$ = the set of reservoirs immediately upstream of reservoir i; $\Delta t$ = the number of days in time-step t; V_iⁱⁿⁱ and V_i^end = the initial and target storages in hm³ at the beginning and end of the planning horizon, respectively; spl_it = the spillage in m³/s in time-step t from reservoir i; q_it = the generating discharge in m³/s in time-step t from hydropower plant i.

The reservoirs on a monthly reservoir scheduling will have some storage capacity allocated for flood control during flood seasons, while they will have their forebay water levels remaining as high as possible during dry seasons to gain a higher water head and allow the hydro plants to produce more energy. The initial storage can usually be observed or forecasted when monthly generations are scheduled shortly after or before the beginning of the scheduling horizon, while the final storage is targeted in an over-year perspective.

2

Upper and lower bounds on storage and release, respectively,

$$\{\begin{array}{l}{V}_i^{\mathrm{dead}}\le V_{it}\le V_{it}^{\max }\\ Q_i^{\min }\le Q_{it}\le Q_i^{\max }\end{array}$$

(4)

where $V_i^{\mathrm{dead}}$ = the dead storage in hm³ of reservoir i; $V_{it}^{\max }$ = the upper bound on the storage at the beginning of t of reservoir i, equal to the flood control limited storage during flooding seasons and the normal storage during dry seasons; Q_i^min and Q_i^max = the lower and upper bounds on the release from reservoir i in time-step t.

3

Firm hydropower output,

$${\displaystyle \sum _{i=1}^N{P}_{it}}\ge F$$

(5)

where N = the number of hydropower plants/reservoirs.

4

Hydropower output and the capacity of generating discharge,

$$P_{it}=A_i\cdot q_{it}\cdot h_{it}$$

(6)

$$q_{it}\le G_i^{\max }(h_{it})$$

(7)

with

$$h_{it}=Z_i^{\mathrm{u}}(\frac{V_{it}+V_{i,t+1}}{2})-Z_i^{\mathrm{d}}(Q_{it})$$

(8)

where A_i = the power generating efficiency in MW; s/m⁴; h_it = the water head in time-step t of hydropower plant i in m; G_i^max (.) = the capacity of generating discharge of i in m³/s, a function of water head; Z_i^u (.) and Z_i^d (.) = the forebay and tailwater elevations in meters, dependent of the water storage and release, respectively, of reservoir i.

3. Solution Techniques

The problem is typically a nonlinear programming, with the nonlinearity coming from constraints (6)–(8), among which the firm hydropower output, along with the capacity of generating discharge that generally has a small variety over the water heads, will be fixed as linear constraints, but updated during the solution with

$$q_{it}\le G_i^{\max }({\widehat{h}}_{it}^{[n]})$$

(9)

$${\displaystyle \sum _{i=1}^N(A_i\cdot {\widehat{h}}_{it}^{[n]}})\cdot q_{it}\ge F$$

(10)

with the water head estimated by

$${\widehat{h}}_{it}^{[n]}=0.5\cdot [{\widehat{h}}_{it}^{[n-1]}+Z_i^{\mathrm{u}}(\frac{V_{it}^{[n]}+V_{i,t+1}^{[n]}}{2})-Z_i^{\mathrm{d}}(Q_{it}^{[n]})]$$

(11)

where ${\widehat{h}}_{it}^{[n]}$ = the water head in meters in time-step t of hydropower plant i, estimated at the beginning of the n^th cycle of updating water heads; $V_{it}^{[n]}$ and $Q_{it}^{[n]}$ = the storage and release at the solution derived at the beginning of the n^th cycle of updating water heads.

Figure 1 illustrates how the problem can be solved in a high-level flowchart, where the original problem with nonlinear constraints will be reformulated into a nonlinear optimization that involves only linear constraints, to which the Zoutendijk method is applicable. A new objective will be used to derive a solution feasible to the nonlinear optimization problem without considering the generating capacity (9), which will then be updated at this solution, which will then be modified to also meet the constraint on the generating capacity. Starting with a feasible solution, the LDP will be solved with the Zoutendijk method to obtain solutions used to update the generating capacity until convergence in updating the water heads.

3.1. Reformulation into Lagrange Dual Problem

$$Z_i^{\mathrm{u}}(V)={\alpha }_i\cdot {(V-{\widehat{V}}_0)}^{{\beta }_i}+{\widehat{Z}}_{\mathrm{u}}$$

(12)

With the forebay and tailwater elevations analytically expressed by

$$Z_i^{\mathrm{d}}(Q)={\lambda }_i\cdot {(Q-{\widehat{Q}}_0)}^{{\theta }_i}+{\widehat{Z}}_{\mathrm{d}}$$

(13)

the water head can be determined as

$$\begin{array}{l}{h}_{it}=h_i(V_{it},V_{i,t+1},Q_{it})\\ =Z_i^{\mathrm{u}}(\frac{V_{it}+V_{i,t+1}}{2})-Z_i^{\mathrm{d}}(Q_{it})\\ ={\alpha }_i\cdot {\left(\frac{V_{it}+V_{i,t+1}}{2}-{\widehat{V}}_0\right)}^{{\beta }_i}+{\widehat{Z}}_i^{\mathrm{u}}-{\lambda }_i\cdot {\left(Q_{it}-{\widehat{Q}}_0\right)}^{{\theta }_i}-{\widehat{Z}}_i^{\mathrm{d}}\end{array}$$

(14)

and the hydropower output will be

$$\begin{array}{l}{P}_{it}=P_i(q_{it},V_{it},V_{i,t+1},Q_{it})\\ =A_i\cdot q_{it}\cdot h_{it}\\ =A_i\cdot q_{it}\cdot \left[{\alpha }_i\cdot {\left(\frac{V_{it}+V_{i,t+1}}{2}-{\widehat{V}}_0\right)}^{{\beta }_i}+{\widehat{Z}}_i^{\mathrm{u}}-{\lambda }_i\cdot {\left(Q_{it}-{\widehat{Q}}_0\right)}^{{\theta }_i}-{\widehat{Z}}_i^{\mathrm{d}}\right]\end{array}$$

(15)

where ${\alpha }_i$, ${\beta }_i$, ${\widehat{V}}_0$ and ${\widehat{Z}}_{\mathrm{u}}$ = the coefficients/parameters to be calibrated to fit the relationship curve of storage versus forebay elevation of reservoir i; ${\lambda }_i$, ${\theta }_i$, ${\widehat{Q}}_0$ and ${\widehat{Z}}_{\mathrm{d}}$ = the coefficients/parameters to be calibrated to fit the relationship curve of the release versus tailwater elevation of reservoir i.

Thus, the problem will have the objective (1), expressed as

$$\max {\mathrm{W}}_1\cdot F+{\mathrm{W}}_2\cdot {\displaystyle \sum _{i=1}^N{\displaystyle \sum _{t=0}^{T-1}{P}_i(q_{it},V_{it},V_{i,t+1},Q_{it})}}$$

(16)

subject to (2)–(4), (9), (10), and (15). Substituting (15) for “P_i (…)” in objective (16) gives a nonlinear programming problem with linear constraints, which allows its solution with methods of feasible directions.

3.2. Initial Feasible Solutions

(1) Initiate a feasible solution without generating capacity

The constraint of generating capacity needs to be estimated based on an initial solution, denoted as

$$\begin{array}{ll}{x}^{[0]}& =[q^{[0]},V^{[0]},Q^{[0]},sp{l}^{[0]},F^{[0]}]\\ & =[q_{it}^{[0]}:(i,t),V_{it}^{[0]}:(i,t),Q_{it}^{[0]}:(i,t),sp{l}_{it}^{[0]}:(i,t),F^{[0]}]\end{array}$$

(17)

where $x_{it}:(i,t)$ represents a vector that consists of variables (x_it) for all the possible combinations (i, t), and $q_{it}^{[0]}$, $V_{it}^{[0]}$,$Q_{it}^{[0]}$, and $sp{l}_{it}^{[0]}$ that determine the water head (${\widehat{h}}_{it}^{[0]}$),

$${\widehat{h}}_{it}^{[0]}=Z_i^{\mathrm{u}}(\frac{V_{it}^{[0]}+V_{i,t+1}^{[0]}}{2})-Z_i^{\mathrm{d}}(Q_{it}^{[0]})]$$

(18)

that will be determined by solving a problem that

$$\max {\mathrm{W}}_1\cdot Y-{\mathrm{W}}_2\cdot {\displaystyle \sum _{i=1}^N{\displaystyle \sum _{t=1}^{T-1}{\left(V_{it}-V_{it}^{\max }\right)}^2}}$$

(19)

Subject to: (2)–(4) and the generating yield (Y),

$${\displaystyle \sum _{i=1}^N{q}_{it}}\ge Y$$

(20)

aimed to approximately maximize the firm hydropower output and the generating efficiency by operating reservoirs at their maximum water levels. To have all decision variables initiated, let

$$F^{[0]}=\underset{t}{\min }{\displaystyle \sum _{i=1}^N[(A_i{\widehat{h}}_{it}^{[0]})\cdot q_{it}]}$$

(21)

(2) Modify the solution to meet generating capacity

The solution used to estimate the generating capacity, when included or updated, may not meet the constraint on the generating capacity with a violation,

$$vio{l}_{it}=\max [0,q_{it}-G_i^{\max }({\widehat{h}}_{it}^{[n]})]$$

(22)

which is used to modify the solution to be feasible with

$$\{\begin{array}{l}{q}_{it}^{[n]}:=q_{it}^{[n]}-vio{l}_{it}\\ sp{l}_{it}^{[n]}:=sp{l}_{it}^{[n]}+vio{l}_{it}\end{array}$$

(23)

3.3. Application of the Zoutendijk Method

(1) Procedure of the Zoutendijk Method

The Zoutendijk algorithm involves two major steps: the first is to find a search direction ($\overrightarrow{d}$) and the second to determine the best step ($\alpha$) moving along this direction.

A nonlinear programming problem with linear constraints can be expressed as:

$$\min f(x),s.t.Ax\ge b\&Ex=e$$

(24)

which has a total of $m$ inequality constraints and $n$ equation constraints. The solution steps of the Zoutendijk method can be summarized as follows:

Step 1: Initiate a feasible solution [$x^{(0)}$] and let $k=0$;

Step 2: Decompose the $A$ into submatrices: $A_1$ and $A_2$, and accordingly $b$ into: $b_1$ and $b_2$ so that $A_1{x}^{(k)}>b_1,A_2{x}^{(k)}=b_2$ according to the constraints on $x^{(k)}$;

Step 3: Solve a linear programming:

$$\begin{array}{l}\min \nabla f{(x^{(k)})}^{T}d\\ s.t.\{\begin{array}{l}{A}_{1}d\ge 0\\ A_{2}d=0\\ -1\le d_i\le 1,\mathrm{for}i=1,\cdots ,n\end{array}\end{array}$$

(25)

to derive the optimum: $d^{*}$, with $-1\le d_i\le 1$ enforced to avoid unbounded solutions and ensure an optimum to the linear programming.

Step 4: If precision is achieved with $|\nabla f{(x^{(k)})}^T{d}^{*}|<\epsilon$, then output the current $x^{(k)}$ as the optimum; otherwise, continue to Step 5.

Step 5: if $A_1{d}^{*}\ge 0$, then seek the best $x^{(k+1)}=x^{(k)}+t_k{d}^{*}$ by a one-dimensional search, and then let $k:=k+1$, then turn to step 2; otherwise, turn to Step 6.

Step 6: Let $u=A_1{x}^{(k)}-b_1,v=A_1{d}^{*}$, calculate $\overline{t}=\min \{-\frac{u_i}{v_i},v_i<0\}$, then update the solution $x^{(k+1)}=x^{(k)}+t_k{d}^{*}$ by solving $\underset{0\le t\le \overline{t}}{\min }f(x^{(k)}+t_k{d}^{*})$, then let $k:=k+1$; turn back to Step 2.

(2) Seek a feasible increasing direction

Start at the beginning of the n^th cycle of updating water heads, and initiate a solution for the Zoutendijk method to be applied,

$$x^{(0)}=x^{[n]}$$

(26)

where the superscript number in parentheses and square brackets represents the number of cycles of updating water heads and iterations in the Zoutendijk method, respectively. Start at the beginning of the k^th iteration,

$$x^{(k)}=[q^{(k)},V^{(k)},Q^{(k)},sp{l}^{(k)},F^{(k)}]$$

(27)

where the Zoutendijk method seeks a feasible increasing direction (d) by

$$\max {\displaystyle \sum _{i=1}^N{\displaystyle \sum _{t=0}^{T-1}(\frac{\partial f}{\partial q_{it}}}}{d}_{it}^{\mathrm{q}})+{\displaystyle \sum _{i=1}^N{\displaystyle \sum _{t=1}^{T-1}(\frac{\partial f}{\partial V_{it}}}}{d}_{it}^{\mathrm{V}})+{\displaystyle \sum _{i=1}^N{\displaystyle \sum _{t=0}^{T-1}(\frac{\partial f}{\partial Q_{it}}}}{d}_{it}^{\mathrm{Q}})+\frac{\partial f}{\partial F}{d}^F$$

(28)

where

$$\frac{\partial f}{\partial F}=W_1$$

(29)

$$\frac{\partial f}{\partial q_{it}}={{\displaystyle \sum _{i=1}^N{\displaystyle \sum _{t=0}^{T-1}{W}_2\frac{\partial P_i}{\partial q_{it}}|}}}_{x=x^{(k)}}$$

(30)

$$\frac{\partial f}{\partial Q_{it}}={{\displaystyle \sum _{i=1}^N{\displaystyle \sum _{t=0}^{T-1}{W}_2\frac{\partial P_i}{\partial Q_{it}}|}}}_{x=x^{(k)}}$$

(31)

$$\frac{\partial f}{\partial V_{it}}={{\displaystyle \sum _{i=1}^N{\displaystyle \sum _{t=0}^{T-1}2W_2\frac{\partial P_i}{\partial V_{it}}|}}}_{x=x^{(k)}}$$

(32)

subject to:

$$d_{i,t+1}^V=d_{it}^V+({\displaystyle \sum _{j\in \Omega (i)}{d}_{jt}^Q}-d_{it}^Q)\cdot \frac{\Delta t\cdot 24\times 3600}{1,000,000}$$

(33)

$$d_{it}^{spl}+d_{it}^q=d_{it}^Q$$

(34)

$${\displaystyle \sum _{i=1}^N{A}_i{\widehat{h}}_{it}^{[n]}{d}_{it}^q}\ge d^F\mathrm{if}{\displaystyle \sum _{i=1}^N{A}_i{\widehat{h}}_{it}^{[n]}{q}_{it}^{(k)}}=F$$

(35)

$$d_{it}^V\ge 0\mathrm{if}{V}_{it}^{(k)}=V_i^{\mathrm{dead}}$$

(36)

$$d_{it}^V\le 0\mathrm{if}{V}_{it}^{(k)}=V_{it}^{\max }$$

(37)

$$d_{it}^Q\ge 0\mathrm{if}{Q}_{it}^{(k)}=Q_i^{\min }$$

(38)

$$d_{it}^Q\le 0\mathrm{if}{Q}_{it}^{(k)}=Q_i^{\max }$$

(39)

$$d_{it}^q\ge 0\mathrm{if}{q}_{it}^{(k)}=0$$

(40)

$$d_{it}^q\le 0\mathrm{if}{q}_{it}^{(k)}=G_i^{\max }({\widehat{h}}_{it}^{[n]})$$

(41)

$$d_{it}^{spl}\ge 0\mathrm{if}sp{l}_{it}^{(k)}=0$$

(42)

$$d^F\ge 0\mathrm{if}{F}^{(k)}=0$$

(43)

$$-1\le d_{it}^V,d_{it}^Q,d_{it}^q,d_{it}^{spl},d_t^{\mu },d^F\le 1$$

(44)

to derive the optimum, denoted as d^(k).

(3) Line search with the golden section method

The gold section method is a line search procedure that starts from a feasible solution and seeks the optimum on a feasible improving direction by

$$\underset{\underset{\_}{\zeta }\le \zeta \le \overline{\zeta }}{\max }f(x^{(k)}+\zeta \cdot d^{(k)})$$

(45)

with the feasible range [$\underset{\_}{\zeta },\overline{\zeta }$] determined to enforce:

$$\{\begin{array}{l}{V}_i^{\mathrm{dead}}\le V_{it}^{(k)}+\zeta \cdot d_{it}^{\mathrm{V}}\le V_{it}^{\max }\\ Q_i^{\min }\le Q_{it}^{(k)}+\zeta \cdot d_{it}^{\mathrm{Q}}\le Q_i^{\max }\\ 0\le q_{it}^{(k)}+\zeta \cdot d_{it}^{\mathrm{q}}\le G_i^{\max }(h_{it}^{[n]})\\ sp{l}_{it}^{(k)}+\zeta \cdot d_{it}^{spl}\ge 0\\ {\displaystyle \sum _{i=1}^N{A}_i{\widehat{h}}_{it}^{[n]}\cdot [q_{it}^{(k)}+\zeta \cdot d_{it}^{\mathrm{q}}]\ge F^{(k)}+\zeta \cdot d^F}\\ F^{(k)}+\zeta \cdot d^F\ge 0\end{array}$$

(46)

4. Case Studies

4.1. Engineering Background

The models and solution techniques presented in this work were applied to the cascaded hydropower plants on the Lancang River, located in southwest China. The total water head drop of the main stream of the Upper Mekong River reaches 5500 m, of which 91% concentrates in the Lancang River, and impressively, 90% of the total drop is concentrated above the Nuozhadu Hydropower Plant. The exploitable hydropower resources of the Lancang River amount to about 30 million kilowatts. Up to now, 11 hydropower plants have been put into operation in the Lancang River. Figure 2 shows the hydraulic layout of the six reservoirs studied in this work, which are all located on the Lancang River in Yunnan Province, China, where the flood seasons are usually from July to September, with the dry seasons in the other months.

Table 2. The basic parameters of cascaded hydropower reservoirs.

Number	Name	Annual Inflow (m3/s)	Installed Capacity (MW)	Storage Capacity (108 m3)	Dam Height (m)	Operability
1	Wunonglong	743	990	2.65	137.5	Daily
2	Huangdeng	902	1900	1.418	203	Seasonal
3	Xiaowan	1210	4200	145.57	1245	Over-year
4	Manwan	1230	1670	3.716	1002	Seasonal
5	Dachaoshan	1340	1350	7.42	906	Seasonal
6	Nuozhadu	1740	5850	217.776	821.5	Over-year

summarizes the basic information of the six cascaded hydropower reservoirs to be studied in this work.

4.2. Results of Curve Fitting

As described in Section 2 and Section 3, the segmented curves Z_i^u (.) and Z_i^d (.) were fitted into two exponential functions to let them be partially derivable so that the Zoutendijk method could be applied to the hydropower scheduling model.

Table 3. Parameters in $Z_i^{\mathrm{u}}(V)={\alpha }_i\cdot {(V-{\widehat{V}}_0)}^{{\beta }_i}+{\widehat{Z}}_{\mathrm{u}}$.

$\mathit{i}$	Name	${\mathit{\alpha }}_{\mathit{i}}$	${\widehat{\mathit{V}}}_0$	${\mathit{\beta }}_{\mathit{i}}$	${\widehat{\mathit{Z}}}_{\mathbf{u}}$	Average MSE
1	Wunonglong	0.13	−5.43	1.01	1868.89	0.00%
2	Huangdeng	0.08	648.49	0.91	1582.38	0.02%
3	Xiaowan	0.04	4662	0.82	1165.98	0.05%
4	Manwan	0.11	136.59	0.89	980.42	0.05%
5	Dachaoshan	0.11	425.35	0.84	884.54	0.00%
6	Nuozhadu	0.02	9554	0.82	760	0.06%

and

Table 4. Parameters in $Z_i^{\mathrm{d}}(Q)={\lambda }_i\cdot {(Q-{\widehat{Q}}_0)}^{{\theta }_i}+{\widehat{Z}}_{\mathrm{d}}$.

$\mathit{i}$	Name	${\mathit{\lambda }}_{\mathit{i}}$	${\widehat{\mathit{Q}}}_0$	${\mathit{\theta }}_{\mathit{i}}$	${\widehat{\mathit{Z}}}_{\mathbf{d}}$	Average MSE
1	Wunonglong	0.002	199.96	0.97	1818.23	0.02%
2	Huangdeng	0.005	149.99	0.89	1473.8	0.03%
3	Xiaowan	0.002	−0.01	1.00	990.83	0.00%
4	Manwan	0.014	197	0.76	896.36	0.02%
5	Dachaoshan	0.009	300	0.9	811.63	0.05%
6	Nuozhadu	0.38	−31.98	0.45	591.49	0.00%

present the parameters calibrated and the fitting errors, measured by the mean squared error (MSE). The results show that the fitting accuracy was mostly above 99.5%, suggesting that Equations (12) and (13) are scientifically credible in representing the forebay and tailwater elevations, respectively.

Figure 2 shows the exponential functions fitted to the curves of the Nuozhadu Reservoir, which has the largest installed capacity in the Lancang River. As described in Table 3 and Table 4, the function of the fitting curve in Figure 3a,b were $Z_6^u(V)=0.02\cdot {(V-9554)}^{0.82}+760$ and $Z_6^{\mathrm{d}}(Q)=0.38\cdot {(Q+31.98)}^{0.45}+591.49$, respectively. Obviously, the fitting accuracy of the two functions was extremely high, with almost every point on the curve where the function was located.

4.3. Scheduling Results of the Cascaded Hydropower Plants

The monthly hydropower scheduling model of the six cascaded reservoirs was solved by coding the Zoutendijk method on Gurobi 9.1.1, which was run on an AMD Ryzen 7 3700X CPU at 3.60 GHz.

Table 5. The performances in solving the model.

Firm Output (MW)	Total Output (MW)	CPU Time (Second)	Number of Optimizations	Is Optimal
146,476.11	2,529,067	7.539	86	Yes

summarizes the performances of the Zoutendijk method in solving the model, and the results demonstrate its capability to secure the optimal solution in quite a short time (in 8 s) to the problem involving a large-scale hydropower system. The firm output, at the top priority, reached 146,476.11 MW, which is the maximum it can be, and the total of hydropower outputs of the cascade, at the second priority, reached 252,9067 MW, which is also the maximum it can be without any sacrifice on the firm output. The Zoutendijk method converged after the 86^th iteration. In general, the algorithm had an excellent performance in solving the scheduling problems of cascaded hydropower reservoirs.

Figure 4 shows that the objective functional value monotonically increased until reaching its convergence criterion after about 86 iterations. In the early stage of optimization, the improvement of the objective value was more obvious, while in the later stage, the improvement speed decreased significantly. It was observed that the spillages were the first decision variables that were optimized in the line search of the golden section method, with a very fast evolution, leading to substantial improvements on the objective in the early stage. When the spillage approached its optimum, the algorithm turned its optimization to the storage, generating discharge and outflow, which are closely coupled with each other, with a small space to explore in the line search, resulting in a slow convergence in the later stage. The increasing monotonicity is a very good property, which ensures the convergence of the method.

As shown in Figure 5, the hydropower scheduling model, when solved with the Zoutendijk method, could ensure a stable power output during dry seasons. As demonstrated in the results, the firm hydropower output from all cascaded hydropower plants could be maintained at a high level thanks to the cooperation of all reservoirs by coordinating their active storage capacities. At the same time, the power output during flood seasons was significant in ensuring the total energy production by utilizing abundant inflows into the reservoirs.

Figure 6 illustrates the results of the optimal scheduling process of the six cascaded hydropower reservoirs. Without help from large reservoirs upstream, the Wunonglong and Huangdeng hydropower reservoirs, with limited refilling capability themselves, had a small amount of spillage scheduled during the flood seasons. The change in water levels at all hydropower plants followed the same trajectory, that is, drawing down during the dry seasons and then rising during the flood seasons to have as little spillage occurring as possible and to generate hydroelectricity as much as possible. The reservoirs had some of their storage capacity set aside during flood seasons to meet the need of flood control. At Manwan, however, its water level remained unchanged throughout the year because it had no duty in flood control and its flood limit water level was equal to the normal water level. No spillage also occurred at Manwan thanks to its upstream reservoir, Xiaowan, which has an over-year operability.

5. Conclusions

This work demonstrates how the Zoutendijk method could be applied, for the first time, to a monthly hydropower scheduling problem of cascaded reservoirs, which was formulated into a nonlinear programming with linear constraints, where the nonlinear hydropower output was represented with an exponential function to ensure its partial derivability while the generating discharge capacity was enforced with an iterative strategy. The model and method were applied to the Lancang River, which involved six cascaded hydropower plants, and the case studies demonstrated the following:

(1)

The Zoutendijk method is well applicable to the monthly hydropower scheduling of cascaded reservoirs, deriving results that are reasonable and reliable;

(2)

The exponential functions used to fit the forebay and tailwater curves showed a very high fitting accuracy at more than 99.5%, showing a great prospect to make segmented curves derivable when formulating a hydropower scheduling problem.

(3)

The Zoutendijk method can significantly increase the total hydropower production while ensuring the highest firm power output of the whole cascade.

(4)

This solution procedure is very fast in securing the optimum to the problem.

This work, however, only conducted a preliminary study and verified the effectiveness and efficiency of the Zoutendijk method, which may need to be compared with other widely used heuristic algorithms to show its strength and to be tested in solving larger scale problems that have more complex constraints.