1. Introduction
Transient analyses of hydropower systems often are carried out via numerical simulation schemes based on time and space discretization, such as the method of characteristics (MOC) [
1], the finite difference/volume method [
2,
3], and the transfer function/differential equations [
4,
5]. The basic idea of these time domain numerical schemes is to transform continuous systems into discrete systems [
6], thus reducing the degree of freedom in the system model and calculating the time evolutions of system states with satisfactory computational cost and numerical precision. However, most time domain numerical schemes are under the constraint of the Courant condition [
7], where the Courant number is required to be close to one to restrain the numerical dissipation. The strong constraint will bring a substantial computational burden to the complex high-order systems. In addition, due to the nonlinearities, such as the dissipative terms and the uncertain parameters in these models [
8], it is difficult to be implemented in the system stability analysis. Therefore, frequency domain methods, as another alternative, were extensively explored for stability analysis [
2,
9].
The conventional frequency domain schemes, e.g., the hydraulic impedance method [
10,
11,
12,
13,
14] and the transfer matrix method [
2], are widely used in applied hydraulic transients. The hydraulic impedance, defined as the complex head deviation divided by the complex flow rate deviation, is always formulated as a function of the complex frequency
s. The relationship between the piezometric head and the discharge at an arbitrary longitudinal position along a pipe or a channel can be explicitly described with the model. This method shows potential for the leakage calibration [
15], blockage detection [
16,
17], and the resonance analysis [
18]. Suo and Wylie [
8] obtained the response of a steady oscillatory flow with the linearized Saint-Venant equations of a pressurized pipe. Kim [
11] formulated an impedance matrix model for complicated pipe networks and discussed the superiority of the impulse response method over other time domain schemes. In addition, he also discussed the dynamic matrix computation method [
12] and incorporated the genetic algorithm into the address-oriented impedance matrix for the generic calibration of the system parameters of the pipe networks with complex topology [
13]. Zhou et al. [
14] investigated the self-excited hydraulic and mechanical vibration characteristics of the pumped storage power station based on the hydraulic impedance method. The free and forced oscillation responses were obtained and analyzed through this approach. Yang [
9] systematically revealed the mathematical expressions of the hydraulic impedances of different components in the hydropower station, including the hydraulic turbine, the pressurized pipe, and the reservoir. Moreover, the stability of the hydropower station under different operating conditions is investigated with the hydraulic impedance method.
The traditional continuous hydraulic impedance can handle complex frequency domain modeling of simple pipe networks, such as pipes in series, pipe networks with single branches, and single loops. However, it cannot address the modeling of high-dimensional hydraulic networks with multi-loops very well. To fully take advantage of the well-developed circuit theory, the methodology of the complex frequency domain equivalent circuit-based hydraulic impedance model is inspired to calculate the specific impedance of any given location along the complicated pipeline by considering the pipe network as a discrete system. The basic idea of this resolution, inspired by an analogy between the hydraulic circuit and the electrical conductor, was first mentioned by Paynter [
19] and Jaeger [
20]. Souza Jr. [
21] performed the numerical analysis of hydraulic transients in a hydropower plant using this methodology. Dr. Nicolet [
22] extended the equivalent circuit model (ECM) to the mathematical modeling of different hydraulic facilities and validated its effectiveness in the simulation of hydraulic transients in engineering practice. Over the past few years, Zhao et al. [
23] applied ECM to the dynamic analysis of operating condition conversion processes in a pumped storage plant. Zheng et al. [
24] expanded the ECM to describe hydraulic transients of a hydraulic system with pipe and open channel flows and achieved satisfactory simulation precision. To the best of the authors’ knowledge, most of the existing literature about this approach only focused on its privilege in the time domain numerical simulation [
19,
20,
21,
22,
23,
24]. Seldom was the equivalent circuit-based complex frequency-domain analysis reported. Although Nicolet [
25] discussed the application of the frequency domain ECM to obtain the frequency responses in the reservoir–pipe–valve system, the system discussed is too simple to fully reflect the complex frequency domain description capability of the ECM. Therefore, it is necessary to further explore the features of the complex frequency domain ECM in more complicated hydraulic systems.
This paper proposes an equivalent circuit-based discrete hydraulic impedance model for the power station system. The equivalent hydraulic circuit is established by an analogy between the water hammer wave propagation in pressurized pipes and the electromagnetic wave propagation in conductors. The discrete hydraulic impedance at any discrete node is obtained by calculating the total specific impedance from one boundary side to the other according to the circuit theory. Subsequently, the detailed procedure of the free oscillation response analysis of the system via the proposed discrete impedance model is introduced. According to the oscillation analysis process, the variation trends of the decay coefficients of different orders of oscillation for the diversion and tailrace pipelines are discussed. Furthermore, the influence of turbine impedance on the decay coefficients and the frequencies corresponding to different orders of oscillation are investigated, respectively.
The rest of this paper is organized as follows: First, the methodology of the equivalent circuit-based hydraulic impedance for different hydraulic components, such as the pressurized pipe and the hydraulic turbine is introduced in
Section 2. In
Section 3, the overall equivalent circuit topology for the hydraulic system of the studied hydropower station is presented, and the detailed numerical solution process of the free oscillation response is illustrated. Then, numerical simulation results of the free oscillation responses for the hydropower station system are presented, and the variation trends of the system oscillation responses are systematically analyzed in
Section 4. Eventually, the conclusions are given in
Section 5.
4. Simulation Result Analyses
Simulation experiments were conducted to reveal the features of different orders of oscillation in the hydropower station system. Since there are two pipelines in the system plant, i.e., the diversion pipeline and the tailrace pipeline, the oscillation characteristics of either pipe may be influenced by the other. In addition, the influence of the hydraulic turbine impedance on the system’s oscillation patterns needs further investigation. Note that the friction losses along the pipelines are neglected to explore the damping effect of the hydraulic turbine.
4.1. Oscillation Mode Characteristics of the Overall Hydropower Station System
It is known that the continuous model of the hydropower station system has infinite degrees of freedom. Therefore, the oscillation order for this system can be as high as infinite. However, as the proposed ECM is a discrete approximation to the system, only the characteristics of the lower orders of oscillation are precise in the modeling, as discussed in
Section 2.1. Firstly, we investigate the decay coefficients of different orders of oscillation. It is found that the damping characteristics of the system oscillations when
and
are quite different from each other. The variation curves of the decay coefficients of the oscillations under different hydraulic turbine operating conditions with relatively small values of
ZT (
) are illustrated in
Figure 7.
Figure 7 shows that the all the decay coefficients gradually move towards zero when the angular frequency is less than 15 rad/s. At the 5th order of oscillation (about 15.58 rad/s), the decay coefficient reaches the crest and is quite close to zero (but is still negative). When the angular frequency is within 15 rad/s to 32 rad/s, the decay coefficient first decreases gradually as the angular frequency increases. At the 10th order of oscillation (about 31.16 rad/s), the decay coefficient reaches the nadir and its value is much smaller than those of the neighboring orders. Then, the decay coefficient gradually increases and reaches another crest at the 15th order of oscillation (about 46.73 rad/s), where the value at this crest is also negative and close to zero.
The decay coefficients of the oscillations of the 1st, 2nd, 3rd, 4th, 6th, 7th, 8th, 9th, 11th, 12th, 13th, 14th, and 15th orders are within the range of [−2, 0], but the value of the decay coefficient of the 10th-order oscillation varies from −22 to 0, which is quite different from the oscillations of other orders. The corresponding angular frequency ranges from 31.1 rad/s to 31.4 rad/s, which is very close to twice the base natural frequency of pipe #2 given in
Table 1.
Different from that of the hydraulic turbine operating conditions with smaller values of
ZT (
ZT ≤ 3.5), the variation curve of the decay coefficient of the system frequency responses under the hydraulic turbine operating conditions with larger values of
ZT (
ZT ≥ 3.6) is illustrated in
Figure 8.
Figure 8 shows that all the decay coefficients gradually increase to zero when the angular frequency is less than 15 rad/s, which is the same trend compared to that in
Figure 7. At the 6th-order oscillation (about 15.70 rad/s), the decay coefficient suddenly decreases sharply to the nadir (which is much smaller than the value of neighboring orders) and then rapidly surges to a normal value (about σ = −0.28) at the 7th-order oscillation. The decay coefficient at the 6th-order oscillation ranges from −18.11 to −2.32 when
ZT decreases from 8 to 3.6. While the angular frequency is within 16 rad/s to 33 rad/s, the decay coefficient at first decreases gradually as the angular frequency increases. The decay coefficient reaches the nadir at the 11th-order oscillation (about 32.88 rad/s). Then, the decay coefficient gradually increases and reaches another crest at the 15th order of oscillation (about 46.72 rad/s), where its value is negative and close to zero. Afterward, the decay coefficient suddenly decreases sharply to the unusual nadir again at the 16th-order oscillation (about 47.11 rad/s) and rapidly rises to the normal crest (within the range of [−0.26, −0.21]) at the 17th-order oscillation (about 49.85 rad/s). From
Figure 8, we can also observe the angular frequencies of the two unusual oscillations (the 6th-order and the 16th-order) are 15.70 rad/s and 47.11 rad/s, respectively. Their values are close to the 1st-order (15.71 rad/s) and 3rd-order (47.12) natural frequencies of pipe #2, respectively.
It is proven that in a simple reservoir–pipe–valve hydraulic system [
25], the frequencies of the free oscillation responses at the valve inlet are even harmonics when the hydraulic impedance of the valve
Zv is smaller than the characteristic impedance of the pipe
Zc (i.e.,
Zv <
Zc). In contrast, the frequencies of the free oscillation responses at the valve inlet are odd harmonics when
Zv >
Zc. This result is quite related to our aforementioned findings. In this study, the oscillation frequency distribution obtained in the hydropower station system also presents obvious odd harmonics performance when
ZT is relatively large and even harmonics performance when
ZT is relatively small. The hydraulic turbine acts as the valve in the reservoir–pipe–valve system. However, the system plant is more complex and consists of multiple pipelines, the coupling effect of the pipes greatly influences the frequency responses of the overall system. The critical hydraulic turbine impedance
ZT to determine the two different oscillation modes is located within the range of [3.5, 3.6]. According to the analysis above, the 10th-order oscillation when
ZT ≤ 3.5 and the 6th-order and 16th-order oscillations when
ZT ≥ 3.6 are considered to correspond to the natural frequencies of pipe #2.
4.2. Variation Trends of the Decay Coefficients of Oscillations of the Two Pipes
Apart from the three oscillations related to pipe #2, all other oscillations are considered related to pipe #1. The decay coefficients corresponding to different orders of oscillation related to pipe #1 under the different conditions of
ZT are depicted in
Figure 9. As seen in
Figure 9, the decay coefficients of the 5th-order oscillation and the 14th-order oscillation are kept close to zero for all values of
ZT. The angular frequencies of the two orders are 15.58 rad/s and 46.73 rad/s, respectively. These two frequencies are close to the 1st-order and 3rd-order oscillations of pipe #2. The 4th-order, 6th-order, 13th-order, and 15th-order oscillations of pipe #1 have relatively weak damping capability because the absolute values of their decay coefficients for all
ZT are relatively small. The reason why they have slow convergence characteristics is that the 4th-order and the 6th-order oscillations are next to the 5th-order oscillation, and the 13th-order and 15th-order oscillations are next to the 14th-order oscillation. On the contrary, the 1st-order, 9th-order, and 10th-order oscillations have a relatively strong damping capability because the absolute values of their decay coefficients for all
ZT are relatively large. The reason why they have fast convergence characteristics is that the angular frequencies of the 9th-order oscillation and the 10th-order oscillation are next to the 2nd-order oscillation of pipe #2, and the 1st-order oscillation is next to the 0th-order oscillation (i.e., the DC component).
In addition, the decay coefficients of the oscillation modes related to pipe #1 are within the range of [−1.6, 0]. For each order of oscillation, the decay coefficient value decreases to the corresponding nadir as ZT increases when ZT is small, then the decay coefficient increases as ZT continues to increase. However, the critical hydraulic turbine impedance of the decay coefficient nadir for each oscillation order is different.
The decay coefficient variation curves of the first three orders of oscillation related to pipe #2 with different
ZT are shown in
Figure 10. It is seen that the 1st-order oscillation and the 3rd-order oscillation can be observed when
and their decay coefficients σ nearly have the same variation trend with each other. σ increases gradually as the value of
ZT increases. On the contrary, the 2nd-order oscillation is observed when
, and σ gradually decreases as the value of
ZT increases.
4.3. Influence of the Turbine Impedance on the Oscillations Frequencies Distribution
In
Section 4.1 and
Section 4.2, the decay coefficients corresponding to different orders of oscillation are investigated. In this subsection, we reveal the variation trend of the oscillation frequencies of the system. From
Figure 7 and
Figure 8, the oscillation modes of the hydropower system are separated into two parts, i.e., the oscillations related to pipe #1 and the oscillations related to pipe #2. As discussed in
Section 4.1, the oscillation frequencies’ distribution of the hydropower station system tends to present odd harmonics performance when
ZT is relatively large and even harmonics performance when
ZT is relatively small. The critical turbine impedance of the two oscillation patterns is within the range of [3.5, 3.6]. Therefore, the variation trends of the frequencies of different orders of oscillation of the system with different
ZT are illustrated in
Figure 11.
Both frequency responses of pipe #1 and pipe #2 when ZT ≤ 2 can be considered even harmonics, although the frequencies of some oscillation orders of pipe #1 in this region are greatly affected by the coupling effect of the two pipes. On the contrary, the frequency responses of pipe #1 and pipe #2 when ZT ≥ 5 are odd harmonics for pipe #1 and pipe #2. Furthermore, an obvious frequency transition region can be observed between ZT = 2 and ZT = 5. The oscillation modes of pipe #1 gradually shift to the odd harmonics from the even harmonics as ZT increases. During this transition, the oscillation modes of pipe #1 are between the two oscillation modes. However, the oscillation modes of pipe #2 suddenly change at a critical turbine impedance (within the region of ZT ∈ [3.5, 3.6]).
5. Conclusions
An equivalent circuit-based hydraulic impedance model for a hydropower station is proposed in this paper. The free oscillation analysis of the hydraulic system of the hydropower station with this method is carried out. Several conclusions are condensed below:
(1) Compared with the analytical frequency responses of the traditional continuous hydraulic impedance method, the frequency responses of the proposed equivalent circuit-based discrete impedance method for higher frequency impulses are more prone to being inaccurate. The overall modeling accuracy is sensitive to the precision of the spatial meshing.
(2) There exists a critical hydraulic turbine impedance ZTC to determine the oscillation modes of the system. Performances of the system’s oscillation pattern for a relatively small ZT and a relatively large ZT are different. When ZT < ZTC, the oscillation modes are even harmonics, and when ZT > ZTC, the oscillation modes are odd harmonics.
(3) An obvious transition region of the frequency response for pipe #1 can be observed around the critical turbine impedance. The oscillation frequency modes of pipe #1 gradually shift to the odd harmonics from the even harmonics as ZT increases. However, the oscillation modes of pipe #2 suddenly change at a critical turbine impedance, and no transition process is observed.