Abstract
The internal flow of reactor coolant pump (RCP) is much more complex than the flow of a general mixed-flow pump due to high temperature, high pressure, and large flow rate. The pressure pulsation that is induced by rotor-stator interaction (RSI) has significant effects on the performance of pump; therefore, it is necessary to figure out the distribution and propagation characteristics of pressure pulsation in the pump. The study uses CFD method to calculate the behavior of the flow. Results show that the amplitudes of pressure pulsation get the maximum between the rotor and stator, and the dissipation rate of pressure pulsation in impellers passage is larger than that in guide vanes passage. The behavior is associated with the frequency of pressure wave in different regions. The flow rate distribution is influenced by the operating conditions. The study finds that, at nominal flow, the flow rate distribution in guide vanes is relatively uniform and the pressure pulsation amplitude is the smallest. Besides, the vortex shedding or backflow from the impeller blade exit has the same frequency as pressure pulsation but there are phase differences, and it has been confirmed that the absolute value of phase differences reflects the vorticity intensity.
1. Introduction
The reactor coolant pump is one of the key equipment of Pressurized Water Reactor (PWR). The pump is characterized by large flow rate, high head, working in high temperature, and high-pressure environment and requires safe and efficient operation for about 40–60 years. The flow in reactor coolant pump is quite complex. The relative movement between rotating impellers and stationary guide vanes [1], as well as the trailing edge vortex shedding [2], results into highly irregular pressure pulsation. The pressure pulsation is known in literature as spreading in a form of waves in the whole hydraulic system [3]. Not only will the pump head, flow rate, and efficiency be influenced, but also potentially hydraulic exciting forces [4] and hydraulic components resonance will be induced [5] or even induce fatigue cracks of the construction [6, 7]. Rzentkowski and Zbroja found that the pump acoustic behavior is exclusively governed by the pressure term, and the excitation frequency agrees well with the pressure pulsation frequency [8]. Dan et al. investigated unsteady flow characteristics in a mixed-flow reactor coolant pump using large-eddy simulation method. Results show that, at the nominal flow rate, in two particular diffuser passages near the spherical casing discharge nozzle, the flow structures are uneven compared with that in the other flow passages. Large-scale flow separation and backflow structures easily occur at the regions near these two passages [9]. However, the pressure pulsation in RCP still needs more investigation.
During the past several decades, researches have been carried out experimentally and numerically to study the effects of impeller or guide vane geometries and operating conditions on pressure pulsations in a mixed-flow pump. Haban et al. explored a precise method for solving the unsteady potential flow; results show that applying the second viscosity into the pressure oscillation analysis, the amplitudes calculated correspond much better to the ones obtained from model measurement, compared with the previous mathematical modeling [10]. According to the flow field distribution in the vanless space, Nicolet et al. set up a simplified acoustic model based on the grid drive excitation source between stationary and moving parts, and the resonance conditions on pressure pulsation were obtained [11]. Zhu et al. calculated pressure pulsation of a reactor coolant pump under low flow conditions. Results showed that the backflow in the flow passage caused the unsteady pressure pulsation in the impellers and guide vanes passage and the circumferential secondary flow. Backflow mainly exists at the inlet and outlet of the flow passages; therefore the pressure in these spaces fluctuates remarkably and behaves as unobvious periodicity [12]. Mixed-flow pump with impellers and guide vanes is adopted in the third-generation reactor coolant pump design. Because of high temperature, high pressure, and large flow rate, the pressure pulsation due to rotor-stator interaction which will apply an alternating load to blades is much stronger than that in general mixed-flow pumps. Some investigations have been conducted to study the optimization design of the mixed-flow pump [13] and flows instabilities [14]. However, most of them focus on the unsteady flow field under different working conditions [15] or the hydraulic performance of the pump [16]. The unsteady flow, especially in the rotor and stator interface regions, is very important to the safety analysis of a nuclear reactor, as it could generate serious axial force and vibration, threatening the integrity of blades and bearings. Therefore, a comprehensive analysis and prediction of the pressure pulsation induced by RSI are essential for the design.
This paper used the commercial CFD tool, Fluent, to carry out unsteady 3D turbulent numerical calculation in a reactor coolant pump. Transient rotor stator interface was used for unsteady flow calculation. Monitoring points were mounted to achieve pressure pulsation characteristics in different direction of the interface. The calculation results are processed by using the Fast Fourier Transform (FFT); thus this paper quite fully revealed the characteristics of pressure pulsation in the region between rotor and stator, and corresponding analysis was given as well.
2. Model and Method
2.1. Model and Mesh
First a mixed-flow reactor coolant pump was designed in laboratory with the main parameters listed in Table 1. The overall flow area consists of inlet pipe, impellers, guide vanes, and spherical casing as shown in Figure 1.
The mesh quality has a critical influence on the numerical calculation accuracy; ANSYS-ICEM was used to construct the grid of the model pump. To guarantee the accuracy of calculation, unstructured mesh was used due to the complexity of computational domain. And grid cells near surface were encrypted as shown in Figure 2, especially on the surface where the pressure gradients are large and flow separation easily occurs [17]. The value of maximum layers of boundary layer was set to 20; the thickness of the first layer of boundary layer was less than 0.01 mm in order to ensure the value + ≤ 1 near the wall. With the grid encrypted, it could provide adequate resolution in the critical regions of the computation domain. What is more, a grid independence study was performed to verify whether the number of nodes of the mesh model used in this study is enough to justify no further improvement on the mesh model [18]. The grid independence check with four different sets of mesh numbers was carried out as listed in Table 2. After the grid independence check, it is noted that the value of the pump head approximately remains a constant when the total numbers of mesh are larger than 9597689; as a result, for the balance of accuracy and calculation speed, case 3 is selected for numerical calculation considering the current computational ability.
2.2. Expected RSI Patterns in Mixed-Flow Pump with Guide Vanes
In the rotor frame, it is generally believed that the nonuniform flow field at the impellers outlet caused by wake effect and blade loading generates a periodic flow pattern. In the stator frame as well, the pressure field attached to the rotating impeller blade induces incoming flow field distortions [19]. The periodic flow of the rotating and stationary frame can be expressed by the following Fourier series: and are pressures of rotating and stationary frames, and are harmonic order, and are amplitudes for the th harmonic in rotating frame and th harmonic in stationary frame, and are phases for the th harmonic in rotating frame and th harmonic in stationary frame, and are angles in the rotating and stationary frames, is the number of impeller blades, and is the number of guide vane blades.
The resulting pressure field, combining the pressure field of guide vanes and impeller blades, is characterized by a strong modulation process as illustrated in Figure 3 [19]. The pressure in the vaneless gap between rotor and stator can be expressed as a superposition of the pressure in rotating frame and stationary frame. And considering that there is a relation between the angle in rotating frame with a rotating speed and that in stationary frame, ; the pressure in vaneless gap can be expressed in stationary frame as follows:
The modulated pressure equation that describes the pressure in vaneless gap is a function of time and space; this pressure field has the following maximum and minimum diametrical pressure modes [20]:
The diametrical mode number k indicates the number of high/low pressure zones of specific frequency component in the circumferential direction. The pressure field of specific frequency component rotates with time due to the rotation of the impellers, and the sign of the diametrical mode number k shows that the rotating direction of pressure field is in the same direction as that of impellers when k is positive or counterrotating when k is negative. Generally, it is believed that the higher the harmonic order, the lower the pressure pulsation energy and the lower amplitudes correspondingly. According to the theory of Kubota et al. [21], the influence of diametrical mode number can be neglected. As results shown in Table 3, the lowest value corresponds to and . This analysis leads us to expect the maximum value of pressure pulsation amplitude resulting from modulations of rotor-stator interactions occurring at 2 times the impeller blade passing frequency in the stationary frame: , where is the pressure pulsation frequency and is the impeller rotating frequency.
2.3. Flow-Controlled Equations
In this paper, the flow field inside the pump is described based on the Reynolds time-averaged three-dimensional Navier-Stokes equations, and turbulent effects are included by means of the Renormalization group (RNG) - turbulence model and standard log-wall functions to solve the flow near solid boundaries [22]. The RNG - turbulence model can predict the behavior of complex flow of separation, secondary flow, and vortex in wall region accurately. Three-dimensional unsteady numerical calculation of the whole flow field of the reactor coolant pump is carried out.
The RNG - turbulence model equations consist of two equations; the first equation is the turbulence kinetic energy equation and the second equation is the dissipation rate equation:
In these equations, represents the generation of turbulence kinetic energy due to the mean velocity gradients. is the generation of turbulence kinetic energy due to buoyancy. represents the contribution of the fluctuating dilatation in compressible turbulence to the overall dissipation rate; , , and are constants. and are the reciprocal of turbulent Prandtl numbers for and , respectively. and are user-defined source terms.
The turbulent (or eddy) viscosity, , is computed by combining and as follows:
In these equations, the constant , , , , and .
2.4. Solution Parameters and Boundary Conditions
The RNG - turbulence model has better ability to predict the behavior of vortex flow field compared with standard - turbulence model with good stability and reasonable precision. First of all, a steady calculation should be achieved, and then the steady results are set as the initial boundary of the transient calculation; the absolute convergence criterion of transient calculation is that the residual is less than . In order to have enough resolution for the transient calculation, time step size is set to 111.1 s. Thus, each impeller revolution is calculated in a time sequence of 360 time steps, which means that every time step the impeller revolves 1 degree [23]. And the convergence of transient calculation has a great influence on the results, so at least thirty periods have been computed for unsteady analysis to guarantee the numerical accuracy.
The boundary conditions are as follows: mass flow inlet at the inlet pipe, average relative pressure of 0 Pa at the spherical casing outlet, frame motion (steady calculation) and mesh motion (transient calculation) for the impellers domain, and standard wall functions for the solid surfaces.
3. Results of the Transient Calculation and Corresponding Analysis
3.1. Validation of the Numerical Method
In order to validate the accuracy of the current numerical calculation, experimental and numerical investigation of the model pump is carried out in a closed-loop test rig. The test used pump components are shown in Figure 4. When the efficiency is calculated in CFD method, the mechanical loss and volume loss are ignored, so the numerical result of efficiency must be modified. Here the mechanical efficiency and volume efficiency are defined as 0.96 and 0.96 [24]. Figure 5 shows the comparisons between experiment results and modified numerical results for - and -. As observed, there is a satisfactory agreement for all the test mass flow. The trend of the curve is basically the same as experiments in the simulations for all 3 turbulence models, so the CFD method can predict the internal flow of pump accurately with admissible error. In the vicinity of the optimal condition , the relative errors of efficiency and head are 5.3% and 6.2% based on the standard - turbulence model, 1.8% and 1.9% based on the RNG k- turbulence mode, and 3.3% and 4.6% based on the SST - turbulence model, respectively. When pump works at partial flow condition , the relative errors of efficiency and head are 4.9% and 5.4% based on the standard - turbulence model, 2% and 1.2% based on the RNG k- turbulence mode, and 2.9% and 1.3% based on the SST - turbulence model, respectively. And when pump works at overload flow condition , the relative errors of efficiency and head are 3.1% and 11.8% based on the standard - turbulence model, 2.4% and 4.2% based on the RNG - turbulence mode, and 2.3% and 8.5% based on the SST - turbulence model, respectively. In general, the external characteristics of the pump calculated by RNG - turbulence model are in good agreement with experimental value. And the error value satisfies the engineering application requirements. Comparing the CFD data with the experimental data, it can be concluded that the current numerical method could capture the mainly unsteady flow structures of the model pump.
(a)
(b)
3.2. Pressure Pulsation Analyses in Time Domain
6 monitoring points, P3, P2, P1, Q3, Q2, and Q1, mounted on the flow passage are used to figure out the characteristics of the pressure pulsation in the streamwise direction, and all the points are locating at the mid-span plane of the passage. 6 monitoring points, M1~M6, are mounted on the interface from shroud to hub, and 6 monitoring points, S1~S6, are mounted on the interface from the pressure side to the suction side. The following results are calculated at nominal flow rate; the monitoring points are shown in Figure 6.
(a)
(b)
(c)
Figure 7(a) is the pressure pulsation graph in time domain of monitoring points P1~P3 and Q1~Q3. It shows that there exists obvious pressure pulsation and that this pulsation gets stronger near the interface. For monitor points of the P series, P1 shows the highest pressure in the average; pressure at monitor P3 is negative because of high-speed rotating impellers forming negative relative pressure at the impellers inlet. Negative pressure can never exist in RCP; here the negative relative pressure exists because the outflow pressure is set to 0 for calculation, while the actual pressure is very high. For monitors of the Q series, the average value of pressure amplitude roughly stays the same. The pulsation amplitude of P1 and Q1 is the maximum with respect to corresponding series.
(a) Monitors P and Q series pressure pulsation in time domain
(b) Monitors P and Q series pressure pulsation amplitude
(c) Monitors M1~M6 pressure pulsation in time domain
(d) Monitors M1~M6 pressure pulsation amplitude
(e) Monitors S1~S6 pressure pulsation in time domain
(f) Monitors S1~S6 pressure pulsation amplitude
The pressure pulsation along the mid-span centerline (P3 to Q3) is shown in Figure 7(b). The high velocity of the fluid at the impellers outlet generates strong impact on the guide vanes; therefore the maximum pulsation amplitude occurs at the guide vanes side near the interface. With the formation of the maximum pressure pulsation, it spreads upstream of the impellers passage and downstream of the guide vanes passage in form of pressure wave, respectively. With the spreading, the energy of the pressure pulsation wave is dissipated by highly turbulent flow and the pressure pulsation amplitude decays rapidly.
From Figure 7(d), it is clearly shown that the amplitude of the pressure pulsation from hub to shroud varies greatly; when in the space close to the hub, the amplitude slightly rises and then decreases. From the graphic, monitor M6 has the maximum pressure pulsation amplitude. This is mainly caused by three factors: (i) Figure 6 shows that M6 is the closest to guide vanes; the pressure wave caused by rotor-stator interaction spreads a short distance, so the energy dissipation is little. (ii) Monitor M6 is located near the impeller shroud and the rotation radius is large, so the circular velocity is large and the high-speed fluid particles strongly impact on the guide vane blades. (iii) The impeller blades are three-dimensional twisted blades, and the guide vanes are straight blades; the guide vanes inlet stagger angle is designed according to the outlet angle of the impellers mid-span streamline, so the fluid particle from point M6 will form a large attack angle at the guide vanes inlet edge. The rotor-stator interaction effect is more significant.
As can be seen from Figure 7(f), from pressure side to suction side, the pressure pulsation amplitude shows small change for most monitors. However, the pressure pulsation amplitude of S6 which locates in the suction side suddenly increases and gets the maximum 9 kPa, the pressure pulsation amplitude in the pressure side is about 52 kPa, and the amplitude in the suction side is 1.8 times larger than that in the pressure side.
3.3. Pressure Pulsation Analysis in Frequency Domain
Fast Fourier Transform (FFT) technology is adopted to acquire the frequency spectrogram from the data shown in Figure 7.
As can be seen from Figure 8, the frequencies of all monitors in impellers passage are 275 Hz and corresponding multiplier. And the frequencies in guide vanes passage are 150 Hz and corresponding multiplier. These frequencies can be predicted from the formula , where n represents the rotation revolution, is the blade number, and is the frequency multiple. For P series, the characteristic frequency can be calculated as Hz. For series, the characteristic frequency can be calculated as Hz. It is fully consistent with the predicted results. The two spectrum graphics clearly show that if the noise signal was filtered, the spectrum of the pressure pulsation induced by rotor-stator interaction would be seen as discrete spectrum. A major feature of discrete spectrum is that its energy is relatively concentrated, so in the process of designing of the pump, the characteristic frequency should be avoided.
With the increase of frequency, the amplitudes of pressure pulsation decay rapidly. It is worth noting that the dissipation rate of the pressure pulsation in the guide vanes passage is lower than that in the impellers passage. It means that the energy of pressure pulsation in impellers passage is more concentrated at the fundamental frequency; on the contrary, the energy in guide vanes passage is more dispersed in the frequency domain. For different frequency components with the same amplitude in pressure field, high-frequency components carry more pulsation energy, and the turbulent dissipation rate caused by fluid viscosity is proportional to the square of frequency [25]. So, the turbulent dissipation of energy in impellers passage is more concentrated. Referring to the Stokes-Kirchhoff equation (5) about classical acoustic absorption coefficient [26], the pressure wave dissipation rate is proportional to the square of frequency.
represents the frequency, is the fluid density, represents the acoustic speed, is the viscous coefficient, is the coefficient of thermal conductivity, and and are the specific heat capacities at constant volume and constant pressure. The frequency in impellers passage is 275 Hz, and the frequency in guide vanes passage is 150 Hz. Therefore, the dissipation rate in impellers passage is times that in guide vanes passage.
From Figure 9 that shows frequency graph of monitor M series, the dominant frequency of monitors M1, M2, M5, and M6 which are close to shroud or hub is 275 Hz. However, the dominant frequency of monitors M3 and M4 which are in the center region is not 275 Hz but the double, that is, 550 Hz. As mentioned previously [21], this phenomenon can be explained by the diameter mode number theory proposed by Kubota et al. and analyzed referring to the analysis method of Franke et al. [20]. In the center region, the fluid flow is quite free from constraint; the flow pattern agrees with the diameter mode theory, so the dominant frequency is 2 times the impeller blade passing frequency. Ming demonstrated that the natural frequency of zero nodal diameter of RCP impellers under wet mode is 552 Hz, which is very close to 550 Hz [27]. It is likely to cause resonance. In the area closing to the boundary, nevertheless, its flow pattern is more complicated and subjected to the boundary condition; the theory is not applicable to this situation any more, and the dominant frequency is the blade frequency.
From the frequency graph of monitors S1~S6 in Figure 10, it is figured out that the dominant frequency of monitors close to blade is 275 Hz and that in the center region is doubling blade frequency. This can be explained by the Kubota theory as analyzed in Section 2.2.
3.4. Pressure Pulsation at Different Operating Conditions
The nuclear coolant pump chamber requires structure intensity, so the nuclear coolant pump uses spherical casing that makes a more uniform force distribution. However, the spherical casing structure forces the redistribution of the pressure flied and the velocity field inside the guide vanes passage; thus the flow structure characteristics in circumferential direction are changed, and the pressure pulsation causes the variations in flow rate. In order to investigate the flow rate variations in different operation conditions, a comparative analysis of the unsteady characteristics of the flow rate distribution in guide vanes passage under the condition of 0.6–1.2 is carried out.
The unsteady results of flow rate in Figure 11 show that, at partial flow 0.6, the flow rate of the guide vanes passage in circumferential direction distributes unevenly, and at the moments , , , and , the flow distribution of the guide vanes shows a periodic change with time. At partial flow 0.8, the nonuniformity of the internal flow distribution in the guide vanes is weakened, and, at design condition, the flow rate distribution is relatively uniform. So, this is the reason why the flow rate pulsations are negligible in previous section. However, at overload flow 1.2, the flow rate distribution trends to be nonuniform again.
(a)
(b)
(c)
(d)
In order to investigate the pressure pulsation characteristics at different operating conditions, pressure coefficient is defined as follows:where is reference pressure, ρ is fluid density, and u2 is tangential velocity at the impeller outlet.
Monitors P3, P2, and P1 and Q3, Q2, and Q1 that are mounted along the streamwise direction are chosen as research objects; Figure 12 shows the comparison of the pressure coefficient at different operating conditions. The results indicate that, for each monitoring point, the pressure pulsation amplitude is the smallest at design condition, and, at partial condition 0.6, the pressure pulsation amplitude varies significantly along the streamwise direction. Especially in the vicinity of the rotor-stator interface, its amplitude is the largest among the four operating conditions; at 0.6 condition, the pressure amplitude of P1 is 2.8 times the design condition. This phenomenon might be caused by large-scale flow separation at partial flow. This large-scale flow separation increased the pressure pulsation in pump system. And Figure 13 indicates that the blade passing frequency is the dominant frequency for all operating conditions.
3.5. Phase Difference Analysis of the Pressure Pulsation Wave
The pressure pulsation can lead to blade vibration; studies show that the unsteady force is influenced by the phase difference of blade vibration [28]. When analyzed in time domain as shown in Figure 7(c), the pressure at monitors M3, M4, and M5 comprises two pressure waves that have the same frequency but there exists phase difference. This phenomenon also exists in monitor S series. The phase difference is associated with complex flow structure at this region, including the vortex shedding from the impeller blade exit and even backflow near the outlet of flow passage. Power spectral density function (PSDF) can represent the distribution of signal energy in different frequency components; the frequency spectrum function represents the amplitude and the phase of each frequency. The frequency of signal can be determined in the amplitude curves of cross-spectral function, and the phase difference can be read in the phase curves of cross-spectral function [29].
From the amplitude curves of PSDF as shown in Figure 14, it can be identified that the extreme frequency of M2 is Hz, Hz, and Hz and that of S2 is Hz, Hz, and Hz. The dominant frequency in different direction is slightly different but almost around 1500/60 = 25 Hz, that is, the impeller rotating frequency. The other frequencies are the impeller blade passing frequency and harmonics. In Table 4, it is shown that, from hub to shroud, the absolute value of phase difference decreases. What is more, in Table 5, from the pressure side to the suction side, the absolute value of phase difference decreases as well. Furthermore, in Figure 7(c), it is shown that, in the center region of the flow passage, the wave superposition phenomenon is obvious; however, in the region close to the boundary, the phenomenon is not so obvious. Substantial interaction between the periodic unsteadiness and the aperiodic turbulence may occur in this case; there are random phase differences between these two kinds of motion at the same frequency. Near the wall in the backflow region, the turbulence structure progressively lags the ensemble-averaged flow oscillation [30]. Figure 15 shows that when the span is 0.05, large scale of backflow or vortex shedding exists, and the phase difference of M1 is the maximum. As a result, comparing Table 4 with Figure 15, the conclusion is that the larger the span is, the less vortex shedding or backflow exists and the smaller the phase difference is.
(a) Phase difference of monitor M2
(b) Phase difference of monitor S2
(a)
(b)
(c)
4. Conclusion
Numerical calculation method is adopted to study the time domain and frequency domain information of the pressure pulsation induced by rotor-stator interaction; this paper draws the following important conclusions.(1)From impellers inlet to guide vanes outlet, the amplitude of pressure pulsation first increases and then decreases and gets the maximum between rotor and stator. The dissipation rate of pressure pulsation wave in stator passage is lower than that in the rotor passage due to low frequency in stator passage.(2)The diametrical modes theory and the calculation results demonstrated that the dominant frequency of the monitors that are close to the boundary is blade passing frequency, that is, 275 Hz, and for those monitors which are in the center region, the dominant frequency is 550 Hz. 550 Hz is close to the natural frequency of the impellers; the designers should keep away from this frequency to avoid resonance.(3)At design condition, the flow distribution in guide vanes is relatively uniform, and the pressure pulsation amplitude is the smallest. When the operation of pump deviates from the design conditions, the flow rate distributes nonuniformly in circumferential direction and shows a periodic change with time.(4)There is wave superposing phenomenon of monitors located between rotor and stator; through the analysis of amplitude curves and phase curves of PSDF, the phase differences of every monitor can be identified. The vorticity structure induced by vortex shedding or backflow progressively lags the ensemble-averaged flow oscillation. There is notable strong positive correlation between phase difference and vorticity. Large amount of vorticity exists in hub and the pressure side of blade.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This project is supported by the National Basic Research Program of China [Grant no. 2015CB057301] and the Science Fund for Creative Research Groups of the National Natural Science Foundation of China [Grant no. 51521064].