# Transient Simulations in Hydropower Stations Based on a Novel Turbine Boundary

Authors: Yanna Liu; Jiandong Yang; Jiebin Yang; Chao Wang; Wei Zeng
Journal: Mathematical Problems in Engineering
 (2016)

---

## Abstract

Most accidents in hydropower stations happened during transient processes; thus, simulation of these processes is important for station design and safety operation. This study establishes a mathematical model of the transient process in hydropower stations and presents a new method to calculate the hydraulic turbine boundary based on an error function of the rotational speed. The mathematical derivation shows that the error function along the equal-opening characteristic curve is monotonic and has opposite signs at the two sides, which means that a unique solution exists to make the error function null. Thus, iteration of the transient simulation is unique and monotonous, which avoids iterative convergence or false solution and improves the solution efficiency compared with traditional methods. Simulation of an engineering case illustrates that the results obtained by the error function are reasonable. Then, the accuracy and feasibility of the mathematical model using the proposed solution are verified by comparison with model and field tests.


---

## Body

## 1. Introduction

Because hydropower stations play an important role in the peak regulation and valley filling of a power grid, the hydraulic turbine needs to frequently change its operating conditions and experiences many transient processes. The transient process in hydropower stations, including the interactions among hydraulics, mechanism, and electricity, is complicated. The closure of guide vanes and spherical valve induces a change in the flow inertia, which causes changes in the turbine rotational speed and hydraulic pressure in the piping system. When the working condition dramatically changes during transients, drastic changes in the water-hammer pressure and high rotational speed may lead to serious accidents that will endanger the safety of the hydraulic structure and turbine unit [1–3] and affect the power grid stability [4]. Therefore, simulating the transient process of hydropower stations is necessary. The calculation accuracy is directly related to the design of the water diversion system, safe operation of the hydropower plant, and power quality.The calculation methods of the water-hammer pressure are analytical [5, 6], graphical [7, 8], and computer-numerical [9]. The analytical method is based on the chain equation of Allievi and the simplifications of the hydraulic turbine as a valve. This method is suitable for simple tube Pelton turbines. On the basis of the analytical method, the graphical method combines the water-hammer wave reflection and superposition with the turbine characteristic curves and guide-vane closing scheme. The traditional numerical simulation can be classified into two types: time- and frequency-domain methods. The time-domain method includes the method of characteristics (MOC) [8–13], finite difference (FD) [14, 15], and finite volume (FV) [16, 17]. In these methods, the water delivery system is divided into a limited number of calculated cells. MOC is mostly applied owing to its advantages: high efficiency of computation, easy implementation for boundary conditions, easy parallelization with FD and FV for pipeline processing [18, 19], and so on. By using MOC, various transition processes can be simulated. In the frequency-domain method [9], the basic equations are linearized to obtain the transfer function of the pipeline. Combined with the transfer functions of the governor and turbine, obtaining the dynamic responses of the hydropower system becomes easy. Irrespective of the time or frequency domain, the transient simulation model includes the pipe model and the turbine boundary.The turbine boundary in the transient simulations can be divided into two categories [20]. One is based on the geometric size of the turbine [21–23]. The other is the characteristic curve based on the model tests, which is generally adopted in the transient simulations. The characteristic curves are usually transformed into curves with the unit speed as an independent variable and the unit discharge and unit torque as dependent variables. To enable transient simulation, the turbine characteristic curves are formulated into piecewise polygonal functions using an auxiliary grid [24–26]. Then, the operating point consisting of unit parameters is expressed as (1)
Q
1
′
=
f
1
n
1
′
,
y
=
A
1
+
A
2
n
1
′
(2)
M
1
′
=
f
2
n
1
′
,
y
=
B
1
+
B
2
n
1
′
.By combining the two equations and other boundary equations, a one-element cubic equation of the square root of the head can be obtained, and the equation can be solved using the Newton–Simpson iterative method [27]. The roots obtained by this method are correlated to the initial values and it can only obtain the root near the initial values. Meanwhile the first-order derivative of the function which is taken as the denominator should be neither too small nor null; otherwise, the iteration cannot be carried out. Generally, we take the root of last instant as the initial value. If the unit parameters obtained from the root of the iteration just lie in the line segment assumed previously, the root is considered as the correct solution and it will be taken as the initial values of the iteration in next instant. If the unit parameters are not in the assumed line segment, we should extend the search range and solve the cubic equation again until we get the right line segment or search all over the equal-opening curve. If all of the line segments of the equal-opening curve cannot get a right root, we will reduce the accuracy of the iteration properly and search again. If all the methods proposed above cannot get a root, then we will select the root whose operating point is the closest to that of last instant as the right solution. Obviously, this calculation method suffers from two shortcomings. First, the root search direction is ambiguous; when the desired operating point is beyond the line segment, it is not clear whether the next search should target the last or next segment of the current line segment. Second, the roots obtained by the Newton–Simpson iterative method are correlated to the initial values; if the initial values are incorrect, the iteration results may be beyond the assumed line segment. Therefore, the equation may either become unsolvable or be incorrectly solved.This study establishes a mathematical model of a transient process by MOC and a novel turbine boundary, as shown in Figure1. In terms of the novel turbine boundary, an error function on the rotational speed is constructed, and a new solution method based on the monotonicity of this error function is presented. The novel solution method is theoretically analyzed, and the superiority of the new solution is validated by transient simulations, model test, and field test.Figure 1
Characteristic grid.


## 2. Mathematical Model and Calculation Method

### 2.1. Mathematical Model of a Novel Turbine Boundary

The mathematical model of the transient process comprises the pipe, turbine boundary, generator, and speed governor models.Two sets of characteristic curves are generally used to represent the turbine characteristics. Here, the traditional characteristic curves are grouped into a nonuniformB-spline surface [28], as shown in Figure 2. Thus, the turbine characteristics can be expressed as(3)
n
1
′
u
,
v
=
f
1
u
,
v
,
Q
1
′
u
,
v
=
f
2
u
,
v
,
M
1
′
u
,
v
=
f
3
u
,
v
.Figure 2
NonuniformB-spline surface of the turbine characteristics.
(a)
Space-curved surface of the characteristic curves
(b)
Space-curved surface ofn
1
′
u
,
v
=
f
1
u
,
v
(c)
Parameterization planeThe measured data points(
n
1
′
,
Q
1
′
,
M
1
′
) of the turbine characteristics are first expressed as a grid B-spline space-curved surface, as shown in Figure 2(a). The space-curved surface is composed of three independent three-dimensional surfaces, as shown in Figure 2(b), or projections on a two-dimensional plane, as shown in Figure 2(c). Parameter u represents the position of an operating point in the opening curve and is evaluated by uniform parameterization u
i
=
i
/
m
(
i
=
0,1
,
…
,
m
). Parameter v stands for relative opening v
j
=
y
j
/
y
n
(
j
=
0,1
,
…
,
n). The ranges of the two parameters are both in [
0,1
], and each group of (
u
,
v
) values corresponds to an operating point in the space-curved surface; that is, when u and v are determined, the three unit parameter values of the operating point become known.According to the one-dimensional continuity and momentum equations of the pipe flow, the equations along the characteristic lines at the front and back of the turbine joint using MOC are as follows [14]:(4)
Q
P
=
Q
C
P
-
C
Q
P
·
H
P
,
(5)
Q
P
=
Q
C
M
+
C
Q
M
·
H
S
.
Q
C
P, C
Q
P, Q
C
M, and C
Q
M in (4) and (5) are all values of the previous instant and are known values in the current instant. Q
C
P and Q
C
M are determined by the flow and head of the previous instant, respectively, and C
Q
P and C
Q
M are both correlated to the wave velocity and pipe area, respectively, and are positive constants.The equations for the unit parameters are expressed as(6)
Q
P
=
Q
1
′
D
1
2
H
P
-
H
S
,
(7)
n
=
n
1
′
H
P
-
H
S
D
1
,
(8)
M
t
=
M
1
′
D
1
3
H
P
-
H
S
.The first derivative of the differential equation of the generator is expressed as(9)
J
d
ω
d
t
=
M
t
-
M
g
-
30
e
g
P
r
n
r
2
π
Δ
n
.Under an off-grid operation,M
g
=
0 and e
g
=
0. Thus, the generator equation can be simplified as(10)
n
=
n
0
+
0.1875
M
t
+
M
t
0
Δ
t
G
D
2
.The guide-vane opening is defined in advance; thus(11)
y
=
f
t
.Parameterv represents the relative opening; therefore(12)
v
=
y
y
m
a
x
.If the turbine is operating under a frequency or power control, (11) is eliminated, whereas the speed governor equation is added.The unknowns in these equations include the headsH
P and H
S of the inlet and outlet of the rotating wheel, respectively, discharge Q
P, rotating speed n, torque M
t, unit speed n
1
′, unit discharge Q
1
′, unit torque M
1
′, guide-vane opening y, and parameters u and v. In general, eleven equations with eleven unknowns exist.


### 2.2. Calculation Steps of the Transient Process

During the calculation of the load rejection in the transient process, the closure law of the guide vane is given, and the guide-vane opening is known at every instant. Thus, solving the abovementioned boundary conditions is simplified to seeking the operating points that satisfy the equations. With a known opening value, the equations can be solved as follows [27].(1) We assume a value for parameteru and obtain the three unit parameters from (3). Inserting X
=
H
p
-
H
S, C
1
=
Q
C
P
/
C
Q
P
+
Q
C
M
/
C
Q
M, C
2
=
1
/
C
Q
P
+
1
/
C
Q
M, and E
=
0.1875
Δ
t
/
G
D
2 into (4)–(6) and simplifying the result yield(13)
X
2
+
C
2
·
Q
1
′
·
D
1
2
X
-
C
1
=
0
.Then, the two roots are obtained as follows:(14a)
X
1
=
-
C
2
Q
1
′
D
1
2
+
C
2
D
1
2
Q
1
′
2
+
4
C
1
2
,
(14b)
X
2
=
-
C
2
Q
1
′
D
1
2
-
C
2
D
1
2
Q
1
′
2
+
4
C
1
2
.
X represents the square root of the head and is always a positive real number. According to (14a) and (14b), X
1
>
0 when C
1
>
0. When C
1
<
0 and Q
1
′
≤
-
-
4
C
1
/
C
2
2
D
1
4, both X
1 and X
2 are positive real numbers, and we need to select the right one. Under other conditions, when both X
1 and X
2 are nonpositive real numbers, Step (
1
) is repeated.(2) We insertX into (7) and (8) and obtain rotating speed n and torque M
t. We insert M
t into (10) and obtain rotating speed n
∗.(3) We defineΩ
=
n
-
n
∗. When Ω
=
0, parameter u is reasonably assumed. X and the three unit parameter values of the operating point are desired. Hence, we proceed to Step (4); otherwise, Steps (1)–(3) are repeated.(4) We insertX into (4)–(6) to obtain the values of H
P, H
S, and Q
P.The flowchart is shown in Figure3. When C
1
>
0, the search interval of Ω
=
0 is the whole equal-opening curve (i.e., from u
=
0 to u
=
1); when C
1
<
0, the search interval is the part in which Q
1
′
≤
-
-
4
C
1
/
C
2
2
D
1
4 is in the pump and reverse-pump regions. In general, this study needs to solve two problems: (
1
) searching for the operating point where Ω
=
0 in a known equal-opening curve and (
2
) distinguishing between X
1 and X
2 when C
1
<
0 and Q
1
′
≤
-
-
4
C
1
/
C
2
2
D
1
4.Figure 3
Flowchart of the solution process of the turbine boundary equations based on the space surface.Assumption 1.
Ω is a monotonic function in the search segment of the equal-opening curve with parameteru.Assumption 2.
Ω has opposite signs at the two boundary points of this segment.From Assumptions1 and 2, a unique point must exist where Ω
=
0. The search direction of the root can be determined according to the monotonicity, and the first problem in this study can be solved. For the second problem, that is, choosing between X
1 or X
2, if Ω is a monotonic function, by inserting X
1 or X
2 into Ω, the value that can make the signs of Ω at the search boundaries opposite is required, and the other value should be discarded. In summary, to solve problems (
1
) and (
2
), we need to validate Assumptions 1 and 2.


## 2.1. Mathematical Model of a Novel Turbine Boundary

The mathematical model of the transient process comprises the pipe, turbine boundary, generator, and speed governor models.Two sets of characteristic curves are generally used to represent the turbine characteristics. Here, the traditional characteristic curves are grouped into a nonuniformB-spline surface [28], as shown in Figure 2. Thus, the turbine characteristics can be expressed as(3)
n
1
′
u
,
v
=
f
1
u
,
v
,
Q
1
′
u
,
v
=
f
2
u
,
v
,
M
1
′
u
,
v
=
f
3
u
,
v
.Figure 2
NonuniformB-spline surface of the turbine characteristics.
(a)
Space-curved surface of the characteristic curves
(b)
Space-curved surface ofn
1
′
u
,
v
=
f
1
u
,
v
(c)
Parameterization planeThe measured data points(
n
1
′
,
Q
1
′
,
M
1
′
) of the turbine characteristics are first expressed as a grid B-spline space-curved surface, as shown in Figure 2(a). The space-curved surface is composed of three independent three-dimensional surfaces, as shown in Figure 2(b), or projections on a two-dimensional plane, as shown in Figure 2(c). Parameter u represents the position of an operating point in the opening curve and is evaluated by uniform parameterization u
i
=
i
/
m
(
i
=
0,1
,
…
,
m
). Parameter v stands for relative opening v
j
=
y
j
/
y
n
(
j
=
0,1
,
…
,
n). The ranges of the two parameters are both in [
0,1
], and each group of (
u
,
v
) values corresponds to an operating point in the space-curved surface; that is, when u and v are determined, the three unit parameter values of the operating point become known.According to the one-dimensional continuity and momentum equations of the pipe flow, the equations along the characteristic lines at the front and back of the turbine joint using MOC are as follows [14]:(4)
Q
P
=
Q
C
P
-
C
Q
P
·
H
P
,
(5)
Q
P
=
Q
C
M
+
C
Q
M
·
H
S
.
Q
C
P, C
Q
P, Q
C
M, and C
Q
M in (4) and (5) are all values of the previous instant and are known values in the current instant. Q
C
P and Q
C
M are determined by the flow and head of the previous instant, respectively, and C
Q
P and C
Q
M are both correlated to the wave velocity and pipe area, respectively, and are positive constants.The equations for the unit parameters are expressed as(6)
Q
P
=
Q
1
′
D
1
2
H
P
-
H
S
,
(7)
n
=
n
1
′
H
P
-
H
S
D
1
,
(8)
M
t
=
M
1
′
D
1
3
H
P
-
H
S
.The first derivative of the differential equation of the generator is expressed as(9)
J
d
ω
d
t
=
M
t
-
M
g
-
30
e
g
P
r
n
r
2
π
Δ
n
.Under an off-grid operation,M
g
=
0 and e
g
=
0. Thus, the generator equation can be simplified as(10)
n
=
n
0
+
0.1875
M
t
+
M
t
0
Δ
t
G
D
2
.The guide-vane opening is defined in advance; thus(11)
y
=
f
t
.Parameterv represents the relative opening; therefore(12)
v
=
y
y
m
a
x
.If the turbine is operating under a frequency or power control, (11) is eliminated, whereas the speed governor equation is added.The unknowns in these equations include the headsH
P and H
S of the inlet and outlet of the rotating wheel, respectively, discharge Q
P, rotating speed n, torque M
t, unit speed n
1
′, unit discharge Q
1
′, unit torque M
1
′, guide-vane opening y, and parameters u and v. In general, eleven equations with eleven unknowns exist.


## 2.2. Calculation Steps of the Transient Process

During the calculation of the load rejection in the transient process, the closure law of the guide vane is given, and the guide-vane opening is known at every instant. Thus, solving the abovementioned boundary conditions is simplified to seeking the operating points that satisfy the equations. With a known opening value, the equations can be solved as follows [27].(1) We assume a value for parameteru and obtain the three unit parameters from (3). Inserting X
=
H
p
-
H
S, C
1
=
Q
C
P
/
C
Q
P
+
Q
C
M
/
C
Q
M, C
2
=
1
/
C
Q
P
+
1
/
C
Q
M, and E
=
0.1875
Δ
t
/
G
D
2 into (4)–(6) and simplifying the result yield(13)
X
2
+
C
2
·
Q
1
′
·
D
1
2
X
-
C
1
=
0
.Then, the two roots are obtained as follows:(14a)
X
1
=
-
C
2
Q
1
′
D
1
2
+
C
2
D
1
2
Q
1
′
2
+
4
C
1
2
,
(14b)
X
2
=
-
C
2
Q
1
′
D
1
2
-
C
2
D
1
2
Q
1
′
2
+
4
C
1
2
.
X represents the square root of the head and is always a positive real number. According to (14a) and (14b), X
1
>
0 when C
1
>
0. When C
1
<
0 and Q
1
′
≤
-
-
4
C
1
/
C
2
2
D
1
4, both X
1 and X
2 are positive real numbers, and we need to select the right one. Under other conditions, when both X
1 and X
2 are nonpositive real numbers, Step (
1
) is repeated.(2) We insertX into (7) and (8) and obtain rotating speed n and torque M
t. We insert M
t into (10) and obtain rotating speed n
∗.(3) We defineΩ
=
n
-
n
∗. When Ω
=
0, parameter u is reasonably assumed. X and the three unit parameter values of the operating point are desired. Hence, we proceed to Step (4); otherwise, Steps (1)–(3) are repeated.(4) We insertX into (4)–(6) to obtain the values of H
P, H
S, and Q
P.The flowchart is shown in Figure3. When C
1
>
0, the search interval of Ω
=
0 is the whole equal-opening curve (i.e., from u
=
0 to u
=
1); when C
1
<
0, the search interval is the part in which Q
1
′
≤
-
-
4
C
1
/
C
2
2
D
1
4 is in the pump and reverse-pump regions. In general, this study needs to solve two problems: (
1
) searching for the operating point where Ω
=
0 in a known equal-opening curve and (
2
) distinguishing between X
1 and X
2 when C
1
<
0 and Q
1
′
≤
-
-
4
C
1
/
C
2
2
D
1
4.Figure 3
Flowchart of the solution process of the turbine boundary equations based on the space surface.Assumption 1.
Ω is a monotonic function in the search segment of the equal-opening curve with parameteru.Assumption 2.
Ω has opposite signs at the two boundary points of this segment.From Assumptions1 and 2, a unique point must exist where Ω
=
0. The search direction of the root can be determined according to the monotonicity, and the first problem in this study can be solved. For the second problem, that is, choosing between X
1 or X
2, if Ω is a monotonic function, by inserting X
1 or X
2 into Ω, the value that can make the signs of Ω at the search boundaries opposite is required, and the other value should be discarded. In summary, to solve problems (
1
) and (
2
), we need to validate Assumptions 1 and 2.


## 3. Search Direction ofX

### 3.1. Directional Derivative ofΩ

According to the directional derivative theorem, the directional derivative of Ω in line segmentP
i
P
i
+
1 shown in Figure 2 can be expressed as [29](15)
∂
Ω
∂
u
=
∂
Ω
∂
n
1
′
cos
⁡
α
+
∂
Ω
∂
Q
1
′
cos
⁡
β
+
∂
Ω
∂
M
1
′
cos
⁡
γ
.We defineΔ
L
=
Δ
n
1
′
2
+
Δ
Q
1
′
2
+
Δ
M
1
′
2 as the length of the segment, cos
⁡
α
=
Δ
n
1
′
/
Δ
L, cos
⁡
β
=
Δ
Q
1
′
/
Δ
L, and cos
⁡
γ
=
Δ
M
1
′
/
Δ
L as the cosines of the angles separately formed by segment P
i
P
i
+
1 and the three coordinate axes in space. Then(16)
∂
Ω
∂
u
=
1
Δ
L
∂
Ω
∂
n
1
′
Δ
n
1
′
+
∂
Ω
∂
Q
1
′
Δ
Q
1
′
+
∂
Ω
∂
M
1
′
Δ
M
1
′
.Therefore, the sign of∂
Ω
/
∂
u depends on the sign of the sum of Δ
n
1
′
∂
Ω
/
∂
n
1
′, Δ
Q
1
′
∂
Ω
/
∂
Q
1
′, and Δ
M
1
′
∂
Ω
/
∂
M
1
′. When ∂
Ω
/
∂
u
>
0, Ω is a monotonically increasing function within P
i
P
i
+
1; when ∂
Ω
/
∂
u
<
0, Ω is monotonically decreasing. From the calculation method presented in Section 2.2 and (7)–(10), Ω can be expressed as(17)
Ω
=
X
n
1
′
D
1
︸
n
-
E
M
1
′
D
1
3
X
2
+
n
0
+
E
M
t
0
︸
n
∗
.The following parts will considerX
=
X
1 as an example for analysis. X
=
X
2 can be analyzed by the same method. The partial derivatives of (17) are as follows:(18)
∂
Ω
∂
n
1
′
=
X
D
1
,
∂
Ω
∂
Q
1
′
=
-
2
C
2
D
1
2
X
C
2
D
1
2
Q
1
′
2
+
4
C
1
n
1
′
2
D
1
-
E
D
1
3
M
1
′
X
,
∂
Ω
∂
M
1
′
=
E
D
1
3
X
2
.Evidently,∂
Ω
/
∂
n
1
′ and ∂
Ω
/
∂
M
1
′ are constants that are greater than zero, and the sign of ∂
Ω
/
∂
Q
1
′ varies under the following scenarios. (
1
) When the operating point is located in the pump region or pump braking regionwhere n
1
′
<
0 and M
1
′
>
0, ∂
Ω
/
∂
Q
1
′ is positive. (
2
) When the operating point is located in the turbine braking region or reverse-pump region where n
1
′
>
0 and M
1
′
<
0, ∂
Ω
/
∂
Q
1
′ is negative. (
3
) When the operating point is located in the turbine region and n
1
′
<
2
E
D
1
4
M
1
′
X (the estimated magnitude of n
1
′ is not greater than 100), ∂
Ω
/
∂
Q
1
′ is positive; otherwise, ∂
Ω
/
∂
Q
1
′ is negative.The signs ofΔ
n
1
′, Δ
Q
1
′, and Δ
M
1
′ depend on the shape of the equal-opening curve. Because the characteristic curve of pump-turbine is the most complex, this study considers this curve as an example. In particular, the equal-opening curve of pump-turbine has the following characteristics [30] (as shown in Figure 4). Point A is the start point of the pump region. Point B is the inflection point of the pump operating region. Point C is the point where the unit discharge of the pump region is zero. Point D is the point where the unit rotating speed is zero. Point E is the point where n
1
′
=
2
E
D
1
4
M
1
′
X in the turbine region. Point F is the point where the unit torque is the largest in the operating region of the turbine. Point G is the point where the unit discharge is the largest in the operating region of the turbine. Point H is the point where the unit rotating speed is the highest in the operating region of the turbine. Point I is the runaway point. Point J is the point where the unit discharge is zero in the turbine region. Point K is the point where the unit rotating speed is the lowest in the braking region of the turbine. Point L is the end point of the reverse-pump region.Figure 4
Key points of the characteristic curves.Table1 lists the changes in the signs of Δ
n
1
′
∂
Ω
/
∂
n
1
′, Δ
Q
1
′
∂
Ω
/
∂
Q
1
′, and Δ
M
1
′
∂
Ω
/
∂
M
1
′ along the equal-opening curve.Table 1
Signs of the different terms at different segments.
Segment
∂
Ω
∂
n
1
′
Δ
n
1
′
∂
Ω
∂
n
1
′
Δ
n
1
′
∂
Ω
∂
Q
1
′
Δ
Q
1
′
∂
Ω
∂
Q
1
′
Δ
Q
1
′
∂
Ω
∂
M
1
′
Δ
M
1
′
∂
Ω
∂
M
1
′
Δ
M
1
′
Pump region
AB
+
+
+
+
+
+
+
−
−
BC
+
+
+
+
+
+
+
−
−
Pump braking region
CD
+
+
+
+
+
+
+
+
+
Turbine region
DE
+
+
+
+
+
+
+
+
+
EF
+
+
+
−
+
−
+
+
+
FG
+
+
+
−
+
−
+
−
−
GH
+
+
+
−
−
+
+
−
−
HI
+
−
−
−
−
+
+
−
−
Braking region
IJ
+
−
−
−
−
+
+
−
−
Reverse-pump region
JK
+
−
−
−
−
+
+
−
−
KL
+
+
+
−
−
+
+
−
−


### 3.2. Determination of the Monotonicity ofΩ

#### 3.2.1. Comparison betweenΔ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
M
1
′
∂
Ω
/
∂
M
1
′

AssumingK
1 as the ratio between Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
M
1
′
∂
Ω
/
∂
M
1
′, we let(19)
F
1
=
Δ
n
1
′
Δ
M
1
′
F
2
=
E
D
1
4
X
.According to (18)(20)
K
1
=
F
1
F
2
.The magnitude of E
D
1
4 is usually approximately 10−4, and the magnitude of X is usually in the range of 101–102. Hence, the magnitude of E
D
1
4
X is below 10−2. Therefore, along the entire equal-opening curve, the slopes approach infinity only in the BC segment of the pump region and at points H and K in the turbine region, as shown in Figure 5, where maybe F
1
<
F
2, Δ
n
1
′
∂
Ω
/
∂
n
1
′
<
Δ
M
1
′
∂
Ω
/
∂
M
1
′; however, these areas are also the segments where Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
>
Δ
M
1
′
∂
Ω
/
∂
M
1
′. According to the signs of these three terms in Table 1, we know that ∂
Ω
/
∂
u
>
0. For the other intervals of the equal-opening characteristic curve, Δ
n
1
′
∂
Ω
/
∂
n
1
′
≫
Δ
M
1
′
∂
Ω
/
∂
M
1
′.Figure 5
Variation inF
1
(
n
1
′
,
Q
1
′
).Consequently, except for the BC segment and the neighborhoods of points H and K, the sign of∂
Ω
/
∂
u depends mainly on the signs of Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′, whereas Δ
M
1
′
∂
Ω
/
∂
M
1
′ can be ignored.


#### 3.2.2. Comparison betweenΔ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′

AssumingK
2 as the ratio between Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′, we let(21)
F
3
=
Δ
n
1
′
Δ
Q
1
′
,
F
4
=
n
1
′
-
2
E
D
1
4
X
M
1
′
Q
1
′
2
+
4
C
1
/
C
2
2
D
1
4
.According to (18)(22)
K
2
=
F
3
F
4
.When F
3
>
F
4, Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′.Figure6 shows how the signs of F
3 and F
4 change along the equal-opening curve. F
3 is determined by the shapes of the characteristic curve. It reaches a low value in the BC segment and at points H and K in the entire opening curve while reaching a peak at point D. F
4 reaches its maximum value when Q
1
′
=
-
-
4
C
1
/
C
2
2
D
1
4 and reaches its minimum value near the n
1
′ axis. In the pump braking and turbine operating regions where Q
1
′ smoothly changes (the range of CG or CH, including EG), Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′. At the ends of the pump and reverse-pump regions, Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ is also possible; however, the two positions are usually beyond the calculation range of the transient process. At the other positions, Δ
n
1
′
∂
Ω
/
∂
n
1
′
<
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ exists.Figure 6
Variation inF
3
(
n
1
′
,
Q
1
′
) and F
4
(
n
1
′
,
Q
1
′
).WhenX
=
X
1, opposite signs for Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ occur at two positions. The first position is the curve segment EG where Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
0, Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
<
0, and |
Δ
n
1
′
∂
Ω
/
∂
n
1
′
|
>
|
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
|. Then, ∂
Ω
/
∂
u
>
0. The second position is the curve segment HK where Δ
n
1
′
∂
Ω
/
∂
n
1
′
<
0, Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
>
0, and |
Δ
n
1
′
∂
Ω
/
∂
n
1
′
|
<
|
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
|. Then, ∂
Ω
/
∂
u
>
0. At other positions, Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ are both greater than zero. Therefore, ∂
Ω
/
∂
u
>
0 along the characteristic curve, and Ω is a monotonically increasing function.WhenX
=
X
2, this condition can be demonstrated in the same manner as X
=
X
1. Ω is a monotonically decreasing function in the pump region as well as in the reverse-pump region (where Q
1
′
≤
-
-
4
C
1
/
C
2
2
D
1
4).In summary,Ω is a monotonic function in the search segments; thus, Assumption 1 is verified.


### 3.3. Signs ofΩ in the Search Boundary

#### 3.3.1. WhenC
1
>
0

If the current operating point(
n
1
,
u
=
p
′
,
Q
1
,
u
=
p
′
) is on the left side of the n
1
′
=
0 axis, n
1
,
u
=
0
′
<
n
1
,
u
=
p
′
<
0, and Q
1
,
u
=
0
′
<
Q
1
,
u
=
p
′ can be obtained. According to (7) and (14a), Ω
u
=
0
<
0 if X
u
=
0
>
X
u
=
p and n
u
=
0
<
n
u
=
p
<
0. For the end point (
n
1
,
u
=
1
′
,
Q
1
,
u
=
1
′
) in the reverse-pump region, n
u
=
1
>
0 and Ω
u
=
1
>
0 if n
1
,
u
=
1
′
>
0. As a result, when C
1
>
0, Ω is a monotonic function within the range of u
=
[
0,1
] and Ω
u
=
0
·
Ω
u
=
1
<
0. Therefore, only one point exists where Ω
=
0 along the entire characteristic curve. The same conclusion can also be drawn when the operating point is on the right side of the n
1
′
=
0 axis.


#### 3.3.2. WhenC
1
<
0

Similarly, whenX
=
X
1, Ω
u
=
0
<
0 and Ω
u
=
1
>
0. When X
=
X
2, Ω
u
=
0
>
0 and Ω
u
=
1
<
0.If0
<
a
<
b
<
1, for the boundary point (
n
1
,
u
=
a
′
,
Q
1
,
u
=
a
′
) where Q
1
,
u
=
a
′
=
-
-
4
C
1
/
C
2
2
D
1
4 in the pump region, we can obtain the result X
=
X
1
=
X
2
=
-
C
1 by inserting Q
1
,
u
=
a
′ into (14a) and (14b). Then, inserting X into (17) yields the following result:(23)
Ω
=
n
1
′
-
C
1
D
1
︸
1
+
E
M
1
′
D
1
3
C
1
︸
2
-
n
0
︸
3
-
E
M
t
0
︸
4
.At the position where n
1
′ is relatively larger, the second and fourth terms in (23) are smaller than the other two terms; thus, they can be ignored. By comparing the first term with the third term, Ω
u
=
a
<
0 when C
1
<
-
n
0
2
D
1
2
/
n
1
′
2; otherwise, Ω
u
=
a
>
0. Similarly, at the boundary point (
n
1
,
u
=
b
′
,
Q
1
,
u
=
b
′
) where Q
1
,
u
=
b
′
=
-
-
4
C
1
/
C
2
2
D
1
4 in the reverse-pump region, Ω
u
=
b
>
0 when C
1
<
-
n
0
2
D
1
2
/
n
1
′
2; otherwise, Ω
u
=
b
<
0.In summary, whenC
1
<
-
n
0
2
D
1
2
/
n
1
′
2, Ω
u
=
0
·
Ω
u
=
a
<
0 (or Ω
u
=
b
·
Ω
u
=
1
<
0) in the search boundary of the pump region (or reverse-pump region) only if X
=
X
2. When 0
>
C
1
>
-
n
0
2
D
1
2
/
n
1
′
2, Ω
u
=
0
·
Ω
u
=
a
<
0 (or Ω
u
=
b
·
Ω
u
=
1
<
0) in the search boundary of the pump region (or reverse-pump region) only if X
=
X
1.Based on the systematic discussions ofC
1
>
0 and C
1
<
0, two conclusions can be drawn. (
1
) The signs of Ω at two boundary points in the search segment are opposite; thus, Assumption 2 is verified. (
2
) When C
1
<
-
n
0
2
D
1
2
/
n
1
′
2, X
=
X
2; meanwhile, when C
1
>
-
n
0
2
D
1
2
/
n
1
′
2, X
=
X
1. Thus, the second problem is solved.


## 3.1. Directional Derivative ofΩ

According to the directional derivative theorem, the directional derivative of Ω in line segmentP
i
P
i
+
1 shown in Figure 2 can be expressed as [29](15)
∂
Ω
∂
u
=
∂
Ω
∂
n
1
′
cos
⁡
α
+
∂
Ω
∂
Q
1
′
cos
⁡
β
+
∂
Ω
∂
M
1
′
cos
⁡
γ
.We defineΔ
L
=
Δ
n
1
′
2
+
Δ
Q
1
′
2
+
Δ
M
1
′
2 as the length of the segment, cos
⁡
α
=
Δ
n
1
′
/
Δ
L, cos
⁡
β
=
Δ
Q
1
′
/
Δ
L, and cos
⁡
γ
=
Δ
M
1
′
/
Δ
L as the cosines of the angles separately formed by segment P
i
P
i
+
1 and the three coordinate axes in space. Then(16)
∂
Ω
∂
u
=
1
Δ
L
∂
Ω
∂
n
1
′
Δ
n
1
′
+
∂
Ω
∂
Q
1
′
Δ
Q
1
′
+
∂
Ω
∂
M
1
′
Δ
M
1
′
.Therefore, the sign of∂
Ω
/
∂
u depends on the sign of the sum of Δ
n
1
′
∂
Ω
/
∂
n
1
′, Δ
Q
1
′
∂
Ω
/
∂
Q
1
′, and Δ
M
1
′
∂
Ω
/
∂
M
1
′. When ∂
Ω
/
∂
u
>
0, Ω is a monotonically increasing function within P
i
P
i
+
1; when ∂
Ω
/
∂
u
<
0, Ω is monotonically decreasing. From the calculation method presented in Section 2.2 and (7)–(10), Ω can be expressed as(17)
Ω
=
X
n
1
′
D
1
︸
n
-
E
M
1
′
D
1
3
X
2
+
n
0
+
E
M
t
0
︸
n
∗
.The following parts will considerX
=
X
1 as an example for analysis. X
=
X
2 can be analyzed by the same method. The partial derivatives of (17) are as follows:(18)
∂
Ω
∂
n
1
′
=
X
D
1
,
∂
Ω
∂
Q
1
′
=
-
2
C
2
D
1
2
X
C
2
D
1
2
Q
1
′
2
+
4
C
1
n
1
′
2
D
1
-
E
D
1
3
M
1
′
X
,
∂
Ω
∂
M
1
′
=
E
D
1
3
X
2
.Evidently,∂
Ω
/
∂
n
1
′ and ∂
Ω
/
∂
M
1
′ are constants that are greater than zero, and the sign of ∂
Ω
/
∂
Q
1
′ varies under the following scenarios. (
1
) When the operating point is located in the pump region or pump braking regionwhere n
1
′
<
0 and M
1
′
>
0, ∂
Ω
/
∂
Q
1
′ is positive. (
2
) When the operating point is located in the turbine braking region or reverse-pump region where n
1
′
>
0 and M
1
′
<
0, ∂
Ω
/
∂
Q
1
′ is negative. (
3
) When the operating point is located in the turbine region and n
1
′
<
2
E
D
1
4
M
1
′
X (the estimated magnitude of n
1
′ is not greater than 100), ∂
Ω
/
∂
Q
1
′ is positive; otherwise, ∂
Ω
/
∂
Q
1
′ is negative.The signs ofΔ
n
1
′, Δ
Q
1
′, and Δ
M
1
′ depend on the shape of the equal-opening curve. Because the characteristic curve of pump-turbine is the most complex, this study considers this curve as an example. In particular, the equal-opening curve of pump-turbine has the following characteristics [30] (as shown in Figure 4). Point A is the start point of the pump region. Point B is the inflection point of the pump operating region. Point C is the point where the unit discharge of the pump region is zero. Point D is the point where the unit rotating speed is zero. Point E is the point where n
1
′
=
2
E
D
1
4
M
1
′
X in the turbine region. Point F is the point where the unit torque is the largest in the operating region of the turbine. Point G is the point where the unit discharge is the largest in the operating region of the turbine. Point H is the point where the unit rotating speed is the highest in the operating region of the turbine. Point I is the runaway point. Point J is the point where the unit discharge is zero in the turbine region. Point K is the point where the unit rotating speed is the lowest in the braking region of the turbine. Point L is the end point of the reverse-pump region.Figure 4
Key points of the characteristic curves.Table1 lists the changes in the signs of Δ
n
1
′
∂
Ω
/
∂
n
1
′, Δ
Q
1
′
∂
Ω
/
∂
Q
1
′, and Δ
M
1
′
∂
Ω
/
∂
M
1
′ along the equal-opening curve.Table 1
Signs of the different terms at different segments.
Segment
∂
Ω
∂
n
1
′
Δ
n
1
′
∂
Ω
∂
n
1
′
Δ
n
1
′
∂
Ω
∂
Q
1
′
Δ
Q
1
′
∂
Ω
∂
Q
1
′
Δ
Q
1
′
∂
Ω
∂
M
1
′
Δ
M
1
′
∂
Ω
∂
M
1
′
Δ
M
1
′
Pump region
AB
+
+
+
+
+
+
+
−
−
BC
+
+
+
+
+
+
+
−
−
Pump braking region
CD
+
+
+
+
+
+
+
+
+
Turbine region
DE
+
+
+
+
+
+
+
+
+
EF
+
+
+
−
+
−
+
+
+
FG
+
+
+
−
+
−
+
−
−
GH
+
+
+
−
−
+
+
−
−
HI
+
−
−
−
−
+
+
−
−
Braking region
IJ
+
−
−
−
−
+
+
−
−
Reverse-pump region
JK
+
−
−
−
−
+
+
−
−
KL
+
+
+
−
−
+
+
−
−


## 3.2. Determination of the Monotonicity ofΩ

### 3.2.1. Comparison betweenΔ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
M
1
′
∂
Ω
/
∂
M
1
′

AssumingK
1 as the ratio between Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
M
1
′
∂
Ω
/
∂
M
1
′, we let(19)
F
1
=
Δ
n
1
′
Δ
M
1
′
F
2
=
E
D
1
4
X
.According to (18)(20)
K
1
=
F
1
F
2
.The magnitude of E
D
1
4 is usually approximately 10−4, and the magnitude of X is usually in the range of 101–102. Hence, the magnitude of E
D
1
4
X is below 10−2. Therefore, along the entire equal-opening curve, the slopes approach infinity only in the BC segment of the pump region and at points H and K in the turbine region, as shown in Figure 5, where maybe F
1
<
F
2, Δ
n
1
′
∂
Ω
/
∂
n
1
′
<
Δ
M
1
′
∂
Ω
/
∂
M
1
′; however, these areas are also the segments where Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
>
Δ
M
1
′
∂
Ω
/
∂
M
1
′. According to the signs of these three terms in Table 1, we know that ∂
Ω
/
∂
u
>
0. For the other intervals of the equal-opening characteristic curve, Δ
n
1
′
∂
Ω
/
∂
n
1
′
≫
Δ
M
1
′
∂
Ω
/
∂
M
1
′.Figure 5
Variation inF
1
(
n
1
′
,
Q
1
′
).Consequently, except for the BC segment and the neighborhoods of points H and K, the sign of∂
Ω
/
∂
u depends mainly on the signs of Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′, whereas Δ
M
1
′
∂
Ω
/
∂
M
1
′ can be ignored.


### 3.2.2. Comparison betweenΔ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′

AssumingK
2 as the ratio between Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′, we let(21)
F
3
=
Δ
n
1
′
Δ
Q
1
′
,
F
4
=
n
1
′
-
2
E
D
1
4
X
M
1
′
Q
1
′
2
+
4
C
1
/
C
2
2
D
1
4
.According to (18)(22)
K
2
=
F
3
F
4
.When F
3
>
F
4, Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′.Figure6 shows how the signs of F
3 and F
4 change along the equal-opening curve. F
3 is determined by the shapes of the characteristic curve. It reaches a low value in the BC segment and at points H and K in the entire opening curve while reaching a peak at point D. F
4 reaches its maximum value when Q
1
′
=
-
-
4
C
1
/
C
2
2
D
1
4 and reaches its minimum value near the n
1
′ axis. In the pump braking and turbine operating regions where Q
1
′ smoothly changes (the range of CG or CH, including EG), Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′. At the ends of the pump and reverse-pump regions, Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ is also possible; however, the two positions are usually beyond the calculation range of the transient process. At the other positions, Δ
n
1
′
∂
Ω
/
∂
n
1
′
<
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ exists.Figure 6
Variation inF
3
(
n
1
′
,
Q
1
′
) and F
4
(
n
1
′
,
Q
1
′
).WhenX
=
X
1, opposite signs for Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ occur at two positions. The first position is the curve segment EG where Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
0, Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
<
0, and |
Δ
n
1
′
∂
Ω
/
∂
n
1
′
|
>
|
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
|. Then, ∂
Ω
/
∂
u
>
0. The second position is the curve segment HK where Δ
n
1
′
∂
Ω
/
∂
n
1
′
<
0, Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
>
0, and |
Δ
n
1
′
∂
Ω
/
∂
n
1
′
|
<
|
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
|. Then, ∂
Ω
/
∂
u
>
0. At other positions, Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ are both greater than zero. Therefore, ∂
Ω
/
∂
u
>
0 along the characteristic curve, and Ω is a monotonically increasing function.WhenX
=
X
2, this condition can be demonstrated in the same manner as X
=
X
1. Ω is a monotonically decreasing function in the pump region as well as in the reverse-pump region (where Q
1
′
≤
-
-
4
C
1
/
C
2
2
D
1
4).In summary,Ω is a monotonic function in the search segments; thus, Assumption 1 is verified.


## 3.2.1. Comparison betweenΔ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
M
1
′
∂
Ω
/
∂
M
1
′

AssumingK
1 as the ratio between Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
M
1
′
∂
Ω
/
∂
M
1
′, we let(19)
F
1
=
Δ
n
1
′
Δ
M
1
′
F
2
=
E
D
1
4
X
.According to (18)(20)
K
1
=
F
1
F
2
.The magnitude of E
D
1
4 is usually approximately 10−4, and the magnitude of X is usually in the range of 101–102. Hence, the magnitude of E
D
1
4
X is below 10−2. Therefore, along the entire equal-opening curve, the slopes approach infinity only in the BC segment of the pump region and at points H and K in the turbine region, as shown in Figure 5, where maybe F
1
<
F
2, Δ
n
1
′
∂
Ω
/
∂
n
1
′
<
Δ
M
1
′
∂
Ω
/
∂
M
1
′; however, these areas are also the segments where Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
>
Δ
M
1
′
∂
Ω
/
∂
M
1
′. According to the signs of these three terms in Table 1, we know that ∂
Ω
/
∂
u
>
0. For the other intervals of the equal-opening characteristic curve, Δ
n
1
′
∂
Ω
/
∂
n
1
′
≫
Δ
M
1
′
∂
Ω
/
∂
M
1
′.Figure 5
Variation inF
1
(
n
1
′
,
Q
1
′
).Consequently, except for the BC segment and the neighborhoods of points H and K, the sign of∂
Ω
/
∂
u depends mainly on the signs of Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′, whereas Δ
M
1
′
∂
Ω
/
∂
M
1
′ can be ignored.


## 3.2.2. Comparison betweenΔ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′

AssumingK
2 as the ratio between Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′, we let(21)
F
3
=
Δ
n
1
′
Δ
Q
1
′
,
F
4
=
n
1
′
-
2
E
D
1
4
X
M
1
′
Q
1
′
2
+
4
C
1
/
C
2
2
D
1
4
.According to (18)(22)
K
2
=
F
3
F
4
.When F
3
>
F
4, Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′.Figure6 shows how the signs of F
3 and F
4 change along the equal-opening curve. F
3 is determined by the shapes of the characteristic curve. It reaches a low value in the BC segment and at points H and K in the entire opening curve while reaching a peak at point D. F
4 reaches its maximum value when Q
1
′
=
-
-
4
C
1
/
C
2
2
D
1
4 and reaches its minimum value near the n
1
′ axis. In the pump braking and turbine operating regions where Q
1
′ smoothly changes (the range of CG or CH, including EG), Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′. At the ends of the pump and reverse-pump regions, Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ is also possible; however, the two positions are usually beyond the calculation range of the transient process. At the other positions, Δ
n
1
′
∂
Ω
/
∂
n
1
′
<
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ exists.Figure 6
Variation inF
3
(
n
1
′
,
Q
1
′
) and F
4
(
n
1
′
,
Q
1
′
).WhenX
=
X
1, opposite signs for Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ occur at two positions. The first position is the curve segment EG where Δ
n
1
′
∂
Ω
/
∂
n
1
′
>
0, Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
<
0, and |
Δ
n
1
′
∂
Ω
/
∂
n
1
′
|
>
|
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
|. Then, ∂
Ω
/
∂
u
>
0. The second position is the curve segment HK where Δ
n
1
′
∂
Ω
/
∂
n
1
′
<
0, Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
>
0, and |
Δ
n
1
′
∂
Ω
/
∂
n
1
′
|
<
|
Δ
Q
1
′
∂
Ω
/
∂
Q
1
′
|. Then, ∂
Ω
/
∂
u
>
0. At other positions, Δ
n
1
′
∂
Ω
/
∂
n
1
′ and Δ
Q
1
′
∂
Ω
/
∂
Q
1
′ are both greater than zero. Therefore, ∂
Ω
/
∂
u
>
0 along the characteristic curve, and Ω is a monotonically increasing function.WhenX
=
X
2, this condition can be demonstrated in the same manner as X
=
X
1. Ω is a monotonically decreasing function in the pump region as well as in the reverse-pump region (where Q
1
′
≤
-
-
4
C
1
/
C
2
2
D
1
4).In summary,Ω is a monotonic function in the search segments; thus, Assumption 1 is verified.


## 3.3. Signs ofΩ in the Search Boundary

### 3.3.1. WhenC
1
>
0

If the current operating point(
n
1
,
u
=
p
′
,
Q
1
,
u
=
p
′
) is on the left side of the n
1
′
=
0 axis, n
1
,
u
=
0
′
<
n
1
,
u
=
p
′
<
0, and Q
1
,
u
=
0
′
<
Q
1
,
u
=
p
′ can be obtained. According to (7) and (14a), Ω
u
=
0
<
0 if X
u
=
0
>
X
u
=
p and n
u
=
0
<
n
u
=
p
<
0. For the end point (
n
1
,
u
=
1
′
,
Q
1
,
u
=
1
′
) in the reverse-pump region, n
u
=
1
>
0 and Ω
u
=
1
>
0 if n
1
,
u
=
1
′
>
0. As a result, when C
1
>
0, Ω is a monotonic function within the range of u
=
[
0,1
] and Ω
u
=
0
·
Ω
u
=
1
<
0. Therefore, only one point exists where Ω
=
0 along the entire characteristic curve. The same conclusion can also be drawn when the operating point is on the right side of the n
1
′
=
0 axis.


### 3.3.2. WhenC
1
<
0

Similarly, whenX
=
X
1, Ω
u
=
0
<
0 and Ω
u
=
1
>
0. When X
=
X
2, Ω
u
=
0
>
0 and Ω
u
=
1
<
0.If0
<
a
<
b
<
1, for the boundary point (
n
1
,
u
=
a
′
,
Q
1
,
u
=
a
′
) where Q
1
,
u
=
a
′
=
-
-
4
C
1
/
C
2
2
D
1
4 in the pump region, we can obtain the result X
=
X
1
=
X
2
=
-
C
1 by inserting Q
1
,
u
=
a
′ into (14a) and (14b). Then, inserting X into (17) yields the following result:(23)
Ω
=
n
1
′
-
C
1
D
1
︸
1
+
E
M
1
′
D
1
3
C
1
︸
2
-
n
0
︸
3
-
E
M
t
0
︸
4
.At the position where n
1
′ is relatively larger, the second and fourth terms in (23) are smaller than the other two terms; thus, they can be ignored. By comparing the first term with the third term, Ω
u
=
a
<
0 when C
1
<
-
n
0
2
D
1
2
/
n
1
′
2; otherwise, Ω
u
=
a
>
0. Similarly, at the boundary point (
n
1
,
u
=
b
′
,
Q
1
,
u
=
b
′
) where Q
1
,
u
=
b
′
=
-
-
4
C
1
/
C
2
2
D
1
4 in the reverse-pump region, Ω
u
=
b
>
0 when C
1
<
-
n
0
2
D
1
2
/
n
1
′
2; otherwise, Ω
u
=
b
<
0.In summary, whenC
1
<
-
n
0
2
D
1
2
/
n
1
′
2, Ω
u
=
0
·
Ω
u
=
a
<
0 (or Ω
u
=
b
·
Ω
u
=
1
<
0) in the search boundary of the pump region (or reverse-pump region) only if X
=
X
2. When 0
>
C
1
>
-
n
0
2
D
1
2
/
n
1
′
2, Ω
u
=
0
·
Ω
u
=
a
<
0 (or Ω
u
=
b
·
Ω
u
=
1
<
0) in the search boundary of the pump region (or reverse-pump region) only if X
=
X
1.Based on the systematic discussions ofC
1
>
0 and C
1
<
0, two conclusions can be drawn. (
1
) The signs of Ω at two boundary points in the search segment are opposite; thus, Assumption 2 is verified. (
2
) When C
1
<
-
n
0
2
D
1
2
/
n
1
′
2, X
=
X
2; meanwhile, when C
1
>
-
n
0
2
D
1
2
/
n
1
′
2, X
=
X
1. Thus, the second problem is solved.


## 3.3.1. WhenC
1
>
0

If the current operating point(
n
1
,
u
=
p
′
,
Q
1
,
u
=
p
′
) is on the left side of the n
1
′
=
0 axis, n
1
,
u
=
0
′
<
n
1
,
u
=
p
′
<
0, and Q
1
,
u
=
0
′
<
Q
1
,
u
=
p
′ can be obtained. According to (7) and (14a), Ω
u
=
0
<
0 if X
u
=
0
>
X
u
=
p and n
u
=
0
<
n
u
=
p
<
0. For the end point (
n
1
,
u
=
1
′
,
Q
1
,
u
=
1
′
) in the reverse-pump region, n
u
=
1
>
0 and Ω
u
=
1
>
0 if n
1
,
u
=
1
′
>
0. As a result, when C
1
>
0, Ω is a monotonic function within the range of u
=
[
0,1
] and Ω
u
=
0
·
Ω
u
=
1
<
0. Therefore, only one point exists where Ω
=
0 along the entire characteristic curve. The same conclusion can also be drawn when the operating point is on the right side of the n
1
′
=
0 axis.


## 3.3.2. WhenC
1
<
0

Similarly, whenX
=
X
1, Ω
u
=
0
<
0 and Ω
u
=
1
>
0. When X
=
X
2, Ω
u
=
0
>
0 and Ω
u
=
1
<
0.If0
<
a
<
b
<
1, for the boundary point (
n
1
,
u
=
a
′
,
Q
1
,
u
=
a
′
) where Q
1
,
u
=
a
′
=
-
-
4
C
1
/
C
2
2
D
1
4 in the pump region, we can obtain the result X
=
X
1
=
X
2
=
-
C
1 by inserting Q
1
,
u
=
a
′ into (14a) and (14b). Then, inserting X into (17) yields the following result:(23)
Ω
=
n
1
′
-
C
1
D
1
︸
1
+
E
M
1
′
D
1
3
C
1
︸
2
-
n
0
︸
3
-
E
M
t
0
︸
4
.At the position where n
1
′ is relatively larger, the second and fourth terms in (23) are smaller than the other two terms; thus, they can be ignored. By comparing the first term with the third term, Ω
u
=
a
<
0 when C
1
<
-
n
0
2
D
1
2
/
n
1
′
2; otherwise, Ω
u
=
a
>
0. Similarly, at the boundary point (
n
1
,
u
=
b
′
,
Q
1
,
u
=
b
′
) where Q
1
,
u
=
b
′
=
-
-
4
C
1
/
C
2
2
D
1
4 in the reverse-pump region, Ω
u
=
b
>
0 when C
1
<
-
n
0
2
D
1
2
/
n
1
′
2; otherwise, Ω
u
=
b
<
0.In summary, whenC
1
<
-
n
0
2
D
1
2
/
n
1
′
2, Ω
u
=
0
·
Ω
u
=
a
<
0 (or Ω
u
=
b
·
Ω
u
=
1
<
0) in the search boundary of the pump region (or reverse-pump region) only if X
=
X
2. When 0
>
C
1
>
-
n
0
2
D
1
2
/
n
1
′
2, Ω
u
=
0
·
Ω
u
=
a
<
0 (or Ω
u
=
b
·
Ω
u
=
1
<
0) in the search boundary of the pump region (or reverse-pump region) only if X
=
X
1.Based on the systematic discussions ofC
1
>
0 and C
1
<
0, two conclusions can be drawn. (
1
) The signs of Ω at two boundary points in the search segment are opposite; thus, Assumption 2 is verified. (
2
) When C
1
<
-
n
0
2
D
1
2
/
n
1
′
2, X
=
X
2; meanwhile, when C
1
>
-
n
0
2
D
1
2
/
n
1
′
2, X
=
X
1. Thus, the second problem is solved.


## 4. Engineering Cases

### 4.1. Case One (Theoretical Validation)

To validate the proposed solution, simulations of the load rejection set in a pumped-storage power station were carried out and presented in this section. The details of the pumped-storage power station, along with simulation results of both the traditional method and proposed solution, are presented as follows.The simulated pumped-storage power station has four units (300 MW each and 1200 MW in total) supplied by a single tunnel, and two surge chambers were installed on both upstream and downstream. Figure7 shows the schematic layout of the station. Table 2 lists the basic information of the single machines and the parameters of the calculated case.Table 2
Basic pump-turbine information of the engineering cases.
RatedP (MW)
Ratedn (rpm)
RatedH (m)
RatedQ (m3/s)
Rotational inertia (kg⋅m2)
CalculatedH (m)
CalculatedQ (m3/s)
CalculatedP (MW)
306.1
428.6
419
83.64
6000000
457.00
73.18
306.1Figure 7
Schematic diagram of the piping system of a hydropower station.During the simulation, all turbines simultaneously rejected their loads, except for one that experienced a guide-vane mechanical failure and its opening remained unchanged. For the other three turbines with no malfunction, their guide vanes held their position for 10 s and then completely shut in 15 s. Figures8(a)–8(d) show the difference between the transient parameters using the new solution proposed in this paper and the traditional solution.Figure 8
Comparison of the results of the traditional and new solutions.
(a)
Variation in trajectories
(b)
Variation in flows
(c)
Variation in pressure in spiral case
(d)
Variation in pressure in draft tubeUntil 20 s, the operating trajectories and variations in the flow and pressure at the inlet and outlet of the rotating wheel obtained by traditional method are more logical and credible with no obvious fluctuation. After 20 s, the operational trajectories become thicker at the reverse-S-shaped region of the turbine characteristic curve and even deviate from the scheduled equal-opening curve. Obvious vibrations occur in both the flow and pressure values, indicating that the traditional method is unable to determine the correct operating points. In contrast, no obvious vibration is observed in the calculation results obtained by the proposed method in this paper, and the calculation results until 20 s coincide with those obtained by the traditional solution. Therefore, this new method is clearly superior to the traditional solution and can improve the calculation accuracy of the transient process.The traditional methods is to solve the head first and then seek the operating point. There may be more than one solutions to the one-element cubic equation and the roots obtained by the Newton–Simpson iterative method are correlated to the initial values. Because lack of sufficient accordance, we take the root of last instant as the initial value. But when the right root is the further one, we will not obtain the right root based on the current initial value. Reducing the iterative accuracy or specifying one root according to last instant may result in bigger errors, and the unit parameters may deviate from the equal-opening curve and the transient parameters will jump up and down nonphysically.The new method determines the unit parameters and the head at the same time. In the solving process, the equations of transient process and the operating point of the equal-opening curve are mutually restrained. All roots obtained are consistent with the physical conditions in the calculation of the transition process: the operating points are in the equal-opening curve; the rotational speeds calculated from speed equation and head equation are equal to each other. Therefore the present new solution method avoids the appearance of wrong roots. The results obtained are in line with the physical phenomenon without obvious vibrations.


### 4.2. Case Two (Model Validation)

A pumped-storage station model consisting of model units, piping systems, measuring system, and generator system was established in the laboratory. The schematic of the model is shown in Figure9. Two different model pump-turbines, manufactured by Harbin Electric Corporation, were installed in the model station. The basic information of the pipe system and the pump-turbines is listed in Tables 3 and 4, respectively.Table 3
Basic unit information of the pipe system.
Pipe section
Inlet
Upstream main pipe
Upstream branch pipe
Upstream branch pipe
Downstream main pipe
Outlet
Equivalent diameter (m)
0.597
0.357
0.202
0.300
0.426
0.675
Length (m)
2.251
24.154
2.195
4.891
16.285
1.062Table 4
Basic unit information of the pump-turbines.
N
Q
E (m⋅kW)
Inlet diameter (mm)
Outlet diameter (mm)
Guide-vane height (mm)
Number of blades
Number of guide vanes
Ratedn (rpm)
RatedH (m)
RatedQ (L/s)
Rotational inertia (N⋅m2)
37.91
280
146.34
24.44
9
20
1000
10.54
49.1
66.4Figure 9
Pumped-storage station model.The basic information on the different sensors used in the experiments is presented in Zeng et al. [31]. The positions of these sensors are shown in Figure 9. Various analog signals from these sensors were collected using an HBM Gen7i data acquisition system. The pump-turbine characteristic curve at the rated guide-vane opening was provided by Zeng et al. [31].A load rejection experiment was conducted on the model station. The guide-vane opening of the pump-turbine was kept constant, and the experimental data are shown as the solid line in Figure10. The simulation is also shown in Figure 10 as a dotted line. The high consistency between the result of the experiment and the simulation confirms that the model has good accuracy.Figure 10
Comparison of the simulated parameters and experimental data.
(a)
Pressure in the spiral case
(b)
Pressure in the draft tube
(c)
Rotating speed


### 4.3. Case Three (Field Validation)

A load rejection test of a pumped-storage power station in South China was conducted by ALSTOM. This section presents the comparison of the test results with the simulation results obtained through the mathematical model proposed in Sections2 and 3. The power station has a similar layout with the case presented in Section 4.1. The basic information of the station is listed in Table 5.Table 5
Basic unit information of the pump-turbines.
RatedP (MW)
Ratedn (rpm)
RatedH (m)
RatedQ (m3/s)
Rotational inertia (kg⋅m2)
TestedH (m)
TestedQ (m3/s)
TestedP (MW)
306.1
500
517.4
66.72
6000000
535.03
62.78
298Two pressure sensors were installed on the unit. One was located downstream of the ball valve 3 m away from its center line. The other was located at the ell of the draft tube. The altitude of this sensor was 133.12 m. By comparing the results of field test and simulation, we can clearly see and thus could arrive at the conclusion that the simulation has good accuracy and can reflect the hydraulic transients of the load rejection with acceptable errors. These errors may be caused by the inaccurate characteristic curve, imprecise parameters of the piping system, and errors in the installation location of the sensors and cross-sectional area of the tubes. In addition, the one-dimensional characteristic method cannot simulate the pressure fluctuation, as clearly shown in Figures11(b) and 11(c) at the peak of the spiral case pressure and nadir of the draft tube pressure, leading to inconsistency between the field test and the simulation.Figure 11
Comparison of the simulation and field test.
(a)
Closing law of guide vanes and spherical valve
(b)
Variation in pressure in spiral case
(c)
Variation in pressure in draft tube
(d)
Variation in rotational speed


## 4.1. Case One (Theoretical Validation)

To validate the proposed solution, simulations of the load rejection set in a pumped-storage power station were carried out and presented in this section. The details of the pumped-storage power station, along with simulation results of both the traditional method and proposed solution, are presented as follows.The simulated pumped-storage power station has four units (300 MW each and 1200 MW in total) supplied by a single tunnel, and two surge chambers were installed on both upstream and downstream. Figure7 shows the schematic layout of the station. Table 2 lists the basic information of the single machines and the parameters of the calculated case.Table 2
Basic pump-turbine information of the engineering cases.
RatedP (MW)
Ratedn (rpm)
RatedH (m)
RatedQ (m3/s)
Rotational inertia (kg⋅m2)
CalculatedH (m)
CalculatedQ (m3/s)
CalculatedP (MW)
306.1
428.6
419
83.64
6000000
457.00
73.18
306.1Figure 7
Schematic diagram of the piping system of a hydropower station.During the simulation, all turbines simultaneously rejected their loads, except for one that experienced a guide-vane mechanical failure and its opening remained unchanged. For the other three turbines with no malfunction, their guide vanes held their position for 10 s and then completely shut in 15 s. Figures8(a)–8(d) show the difference between the transient parameters using the new solution proposed in this paper and the traditional solution.Figure 8
Comparison of the results of the traditional and new solutions.
(a)
Variation in trajectories
(b)
Variation in flows
(c)
Variation in pressure in spiral case
(d)
Variation in pressure in draft tubeUntil 20 s, the operating trajectories and variations in the flow and pressure at the inlet and outlet of the rotating wheel obtained by traditional method are more logical and credible with no obvious fluctuation. After 20 s, the operational trajectories become thicker at the reverse-S-shaped region of the turbine characteristic curve and even deviate from the scheduled equal-opening curve. Obvious vibrations occur in both the flow and pressure values, indicating that the traditional method is unable to determine the correct operating points. In contrast, no obvious vibration is observed in the calculation results obtained by the proposed method in this paper, and the calculation results until 20 s coincide with those obtained by the traditional solution. Therefore, this new method is clearly superior to the traditional solution and can improve the calculation accuracy of the transient process.The traditional methods is to solve the head first and then seek the operating point. There may be more than one solutions to the one-element cubic equation and the roots obtained by the Newton–Simpson iterative method are correlated to the initial values. Because lack of sufficient accordance, we take the root of last instant as the initial value. But when the right root is the further one, we will not obtain the right root based on the current initial value. Reducing the iterative accuracy or specifying one root according to last instant may result in bigger errors, and the unit parameters may deviate from the equal-opening curve and the transient parameters will jump up and down nonphysically.The new method determines the unit parameters and the head at the same time. In the solving process, the equations of transient process and the operating point of the equal-opening curve are mutually restrained. All roots obtained are consistent with the physical conditions in the calculation of the transition process: the operating points are in the equal-opening curve; the rotational speeds calculated from speed equation and head equation are equal to each other. Therefore the present new solution method avoids the appearance of wrong roots. The results obtained are in line with the physical phenomenon without obvious vibrations.


## 4.2. Case Two (Model Validation)

A pumped-storage station model consisting of model units, piping systems, measuring system, and generator system was established in the laboratory. The schematic of the model is shown in Figure9. Two different model pump-turbines, manufactured by Harbin Electric Corporation, were installed in the model station. The basic information of the pipe system and the pump-turbines is listed in Tables 3 and 4, respectively.Table 3
Basic unit information of the pipe system.
Pipe section
Inlet
Upstream main pipe
Upstream branch pipe
Upstream branch pipe
Downstream main pipe
Outlet
Equivalent diameter (m)
0.597
0.357
0.202
0.300
0.426
0.675
Length (m)
2.251
24.154
2.195
4.891
16.285
1.062Table 4
Basic unit information of the pump-turbines.
N
Q
E (m⋅kW)
Inlet diameter (mm)
Outlet diameter (mm)
Guide-vane height (mm)
Number of blades
Number of guide vanes
Ratedn (rpm)
RatedH (m)
RatedQ (L/s)
Rotational inertia (N⋅m2)
37.91
280
146.34
24.44
9
20
1000
10.54
49.1
66.4Figure 9
Pumped-storage station model.The basic information on the different sensors used in the experiments is presented in Zeng et al. [31]. The positions of these sensors are shown in Figure 9. Various analog signals from these sensors were collected using an HBM Gen7i data acquisition system. The pump-turbine characteristic curve at the rated guide-vane opening was provided by Zeng et al. [31].A load rejection experiment was conducted on the model station. The guide-vane opening of the pump-turbine was kept constant, and the experimental data are shown as the solid line in Figure10. The simulation is also shown in Figure 10 as a dotted line. The high consistency between the result of the experiment and the simulation confirms that the model has good accuracy.Figure 10
Comparison of the simulated parameters and experimental data.
(a)
Pressure in the spiral case
(b)
Pressure in the draft tube
(c)
Rotating speed


## 4.3. Case Three (Field Validation)

A load rejection test of a pumped-storage power station in South China was conducted by ALSTOM. This section presents the comparison of the test results with the simulation results obtained through the mathematical model proposed in Sections2 and 3. The power station has a similar layout with the case presented in Section 4.1. The basic information of the station is listed in Table 5.Table 5
Basic unit information of the pump-turbines.
RatedP (MW)
Ratedn (rpm)
RatedH (m)
RatedQ (m3/s)
Rotational inertia (kg⋅m2)
TestedH (m)
TestedQ (m3/s)
TestedP (MW)
306.1
500
517.4
66.72
6000000
535.03
62.78
298Two pressure sensors were installed on the unit. One was located downstream of the ball valve 3 m away from its center line. The other was located at the ell of the draft tube. The altitude of this sensor was 133.12 m. By comparing the results of field test and simulation, we can clearly see and thus could arrive at the conclusion that the simulation has good accuracy and can reflect the hydraulic transients of the load rejection with acceptable errors. These errors may be caused by the inaccurate characteristic curve, imprecise parameters of the piping system, and errors in the installation location of the sensors and cross-sectional area of the tubes. In addition, the one-dimensional characteristic method cannot simulate the pressure fluctuation, as clearly shown in Figures11(b) and 11(c) at the peak of the spiral case pressure and nadir of the draft tube pressure, leading to inconsistency between the field test and the simulation.Figure 11
Comparison of the simulation and field test.
(a)
Closing law of guide vanes and spherical valve
(b)
Variation in pressure in spiral case
(c)
Variation in pressure in draft tube
(d)
Variation in rotational speed


## 5. Conclusions

This study has presented a mathematical model of the hydraulic transient process and proposed a new method to solve the turbine boundary. The new solution was based on error function Ω of the rotational speed. According to the change in the values and signs of the error function along the equal-opening characteristic curve, we can discriminate and search the root of the transient equations; whenΩ
=
0, the iteration in the numerical solution was completed.This study has shown thatΩ is a monotonic function in the search segment along the equal-opening characteristic curve and that its signs in the search boundary are opposite to each other. The sign of Ω may be used to determine the search direction of the operating point in the equal-opening characteristic curve in the next time step. Within the search range, only one point exists that could fulfill the condition Ω
=
0; hence, the equation set has only one unique solution.The simulation of the engineering case illustrates that the operating trajectories obtained by the method proposed in this paper are more reasonable and credible, and the pressure at the end of the spiral case and at the entrance of the draft tube obtained by this method does not sharply jump. In contrast with the model and the field test results, the simulation results can accurately reflect the change process of the transient parameters in the load rejection. Therefore, this method is clearly superior to the traditional solution and can improve the simulation accuracy of the transient process of a hydropower station.



---

*Source: 1504659-2016-06-30.xml*