Modulation effect on rotor-stator interaction subjected to fluctuating rotation speed in a centrifugal pump

In the paper, the modulation effect of rotor-stator interaction subjected to fluctuating rotation speed is investigated numerically. A quasi-bivariate variable mode decomposition method for analyzing non-stationary flow fields is proposed for the first time. Different from the previous studies, this paper emphasizes the local features of the transient flow field at fluctuating rotation speed. The results show that the impact of fluctuating rotation speed on the impeller force is limited, but it significantly impacts the pressure pulsation at fluctuating frequency. The correlation between amplitude and fluctuation ratio is positive. Meanwhile, the fluctuating rotation speed also induces sidebands centered on the blade passing frequency and its harmonics. The periodic changes in rotational speed lead to frequency and amplitude modulation of pressure pulsations. The magnitude of the instantaneous frequency and amplitude increases as the fluctuation ratio grows. The change in pressure field is well synchronized with the rotational speed, while the change in velocity field lags behind the rotational speed. The variation of pressure field with rotational speed is mainly attributed to two aspects. One is the flow structure at the fluctuating frequency, and the other is the instantaneous component of the flow structure at blade passing frequency.


Introduction
Centrifugal pumps are widely applied in marine and energy engineering.Along with the larger power and higher speed, people pay more and more attention to the operation stability and vibration of the pump (Zhang et al., 2021a).The rotor-stator interaction (RSI) is the primary source of vibration in the pump (Lin et al., 2021).As early as 1982, (Dring et al., 1982) experimentally investigated the RSI phenomenon in the axial turbine stage, which is divided into two areas.These are potential flow disturbances that extend both upstream and downstream and wake disturbances that propagate only downstream.Subsequently, huge numbers of studies have been performed to investigate the RSI effect of centrifugal pumps over past decades.In centrifugal pumps, the blade passing frequency (BPF) and its harmonics are the central frequencies of the pressure fluctuations (Yu et al., 2022b).The amplitude of pressure pulsation under design conditions is minimized (Zhang et al., 2017b).Compared to the impeller, the RSI effect is more pronounced in the volute (Chalghoum et al., 2018).The pressure distribution around the impeller outlet at CONTACT Chenxing Jiang jiangchx@xmu.edu.cn;Zhijun Shuai shuaizhijun@hrbeu.edu.cnBPF shows a modulation pattern (Barrio et al., 2011).
The pressure fluctuations at the volute tongue are the strongest (Pei et al., 2014).The spatial fluctuation pattern is closely related to the flow rate (Gonzalez et al., 2002).Xu et al (Xu et al., 2021) developed a visual method to investigate the wake and potential interactions.Zhang et al. (Zhang et al., 2016) analyzed the flow structure in the pump using large eddy simulation and captured the impingement of the vortex from the impeller outlet on the volute tongue.They found that the wake flow at the impeller outlet is the main factor inducing pressure pulsation (Zhang et al., 2019).It is crucial to repress the secondary flow inside the impeller and improve the velocity homogeneity at the impeller outlet to reduce hydro-induced vibration (Zhang et al., 2021b;Qian et al., 2020).Huang et al. (Huang et al., 2021), Gangipamula et al. (Gangipamula et al., 2022) and Ye et al. (Ye et al., 2022) reduced pressure fluctuations by modifying the profile of the blade trailing edge.The common ways to optimize pulsation also include increasing the radial gap and adding splitter blades to the impeller (Solis et al., 2011).Meanwhile, the volute tongue plays an essential role in the RSI (Parrondo-Gayo et al., 2002).Zhang et al. (Zhang et al., 2020) tried to optimize the pressure pulsation in the pump by cutting the volute tongue.The bionics is also a useful approach to reduce vibration and noise (Li et al., 2019a;Dai et al., 2021).
Accompanying the increasing complexity of the system, the rotor instability problem is becoming prominent.The non-stationary rotating speed induces excessive vibration, which produces accelerated fatigue and even catastrophic failure.In practice, torsional vibration of the shaft and voltage instability will induce nonstationary rotating speeds (Lv et al., 2017;Zhang et al., 2014).Tsukamoto et al. (Tsukamoto et al., 1995) studied the dynamic characteristics of the pump at fluctuating rotational speeds theoretically and experimentally.The results showed that the dynamic characteristics deviate significantly from the quasi-steady ones as the fluctuation frequency increases.This is consistent with the conclusion obtained by Zhang et al. (Zhang et al., 2014).Meanwhile, they found that the head and rotational speed change synchronously, while the flow rate lags behind the rotational speed.Sun et al. (Sun et al., 2016) found that the torque fluctuations at non-stationary rotating speeds are more significant than those at stationary rotating speeds.Lv et al. (Lv et al., 2017) studied the effect of fluctuating rotation speed on pressure fluctuation characteristics.They concluded that the shaft torsional vibrations induced sidelobe peaks near the blade passing frequency (BPF).In the past few years, scholars have gradually focused on the effect of fluctuating rotation speed on the hydraulic characteristics of centrifugal pumps.However, the RSI mechanism under non-stationary rotating speeds is still unclear.Few studies have emphasized the local features of pressure pulsations.
The pressure pulsations in the pump at fluctuating rotational speeds show complex modulation.Variational mode decomposition (VMD) is suitable for decomposing signals characterized by non-stationary and multicomponent (Sharma & Parey, 2020).The VMD method is proposed by Dragomiretskiy (Dragomiretskiy & Zosso, 2014).The optimal solution of the variational modes is searched by iteration to obtain modal functions with a certain bandwidth.The VMD algorithm is applied widely for fault diagnosis (Li et al., 2019b;Huang et al., 2019;Liu et al., 2016b;Tang et al., 2021;Dibaj et al., 2021;Zhang et al., 2017a).However, it only works with onedimensional signals, and there is little literature applying it to flow field analysis.This paper proposes a quasibivariate VMD method to extract the flow structures with characteristic frequencies in the temporal-spatial flow field.
This paper introduces CFD tools to study centrifugal pumps subjected to fluctuating rotation speeds.Different from the previous studies (Lv et al., 2017;Sun et al., 2016;Zhang et al., 2014), this paper focuses on the transient flow induced by fluctuating rotation speed, with particular attention to the time-frequency characteristics of the pressure pulsation.A quasi-bivariate VMD method is innovatively proposed, which can be applied to obtain reduced-order models of non-stationary flow fields.The relationship between the flow pattern and the fluctuating rotation speed is elaborated.

Pump model and parameters
In this paper, the research object is a centrifugal pump with low specific speeds.Figure 1 shows the pump model, which consists of the impeller and the volute.In Table 1, the primary geometry and performance parameters are listed.The impeller inlet and volute outlet are appropriately extended to obtain stable flow.The dimensionless type number Kr is defined as follows (Yu et al., 2022a).
where, g is the gravitational acceleration.The computational domain is decomposed using tetrahedral meshes (as shown in Figure 2).The velocity along the normal wall in the boundary layer varies drastically.Prismatic meshes are used at the walls to reduce   the number of grids.Simultaneously, the impeller blades and the volute tongue are locally refined to capture flow details.

Numerical method
ANSYS CFX 17.2 is used to solve the N-S equations.The SST k-ω model is adopted to simulate turbulence.The boundary condition at the inlet is total pressure.The turbulence intensity is set to 5%.The boundary condition at the outlet is opening.Simultaneously, an additional loss proportional to the velocity is used to simulate the effect of the valve.The advection scheme is set to high resolution.The impeller is set to the rotating domain.
The fluctuating rotation speed is simplified as the rotating speed that varies sinusoidally (Guo et al., 2018).The rotating speed of the impeller is defined in Eq. ( 2).The fluctuation ratio γ is the ratio of ω 1 to ω 0 .
where, ω 1 and f 1 are the fluctuation amplitude and frequency of the rotating speed.In this paper, the fluctuation frequency (FF) is 50 Hz.
In the steady simulation, the frozen rotor is adopted to account for centrifugal effects.In the transient simulation, the transient frozen rotor is adopted to simulate RSI effects.The results obtained from the steady simulation are used as the initial flow field.The convergence criterion is that the monitored variables are periodically varying.
The non-dimensional pressure coefficient (C P ) is defined: . The dependence of the numerical prediction concerning the grid and timestep is investigated, as shown in Figure 3.For this purpose, six grids with different numbers are generated, and four time-steps are selected.The basic unit of the timestep is 5.556×10 −5 s, which corresponds to the time for the impeller to rotate 1°.When the grid number exceeds 232×10 4 and the time step is below one basic unit, the results do not converge further.Therefore, the computation model with 232×10 4 grids is adopted, and the time-step is set to 1 basic unit.
Y plus is the dimensionless wall distance, which is often used as the selection criterion for the solution method in the near-wall region (Chen et al., 2022;Shi et al., 2022).The definition of Y plus is as follows: where, y is the vertical height of the first boundary layer to the wall, and ν, τ w , and ρ are the kinematic viscosity, shear stress, and fluid density, respectively.Figure 4 shows the values of Y plus on stationary and rotating walls.The value on all walls is less than 3.The boundary layer complies with the requirements of the turbulence model.

Experiment validation
The layout of the closed test system is shown in Figure 5.The hydraulic characteristics of the centrifugal pump are tested.Two fast-response piezoresistive pressure sensors (CYG-401, SQ) with less than ±0.2% uncertainty record static pressure signals at the pump inlet and outlet.A piezoelectric pressure sensor (113B28-PCB) with less than ±0.1% uncertainty records pressure pulsation signals at the outlet.An electromagnetic flowmeter (LDTH-50B) with less than ±0.5% uncertainty is used to measure the flow rate.
The comparison between the experimental data and the simulation results is shown in Figure 6.Both are basically similar, which also demonstrates the correctness of the computation procedure.The simplification of the computational model and boundary conditions is primarily to blame for the imperfect agreement between experimental and simulation results.

Variational mode decomposition
The details of the VMD algorithm are available in Dragomiretskiy's research (Dragomiretskiy & Zosso, 2014), which is briefly described below.The assumption is that the original signal is composed of a certain number of Intrinsic Mode Functions (IMF).Each IMF is defined as an AM-FM signal in Eq. ( 7).
where, A k (t) is the instantaneous amplitude and φ k (t) is the instantaneous phase of the signal u k (t).Subject to the minimum sum of estimated bandwidths for each IMF, the constrained variational model is constructed in Eq. ( 8).where, } are the center frequencies corresponding to the IMFs.The quadratic penalty term and Lagrangian multipliers are introduced to convert the constrained variational problem into an unconstrained one.Constraints are strictly enforced using Lagrangian multipliers.The fidelity of the signal reconstruction can be guaranteed by the quadratic penalty term.The augmented Lagrangian function is as follows.
The above variational problem is solved by the alternate direction method of multipliers.The saddle point of the function is sought by alternately updating u n+1 k , ω n+1 k and λ n+1 .An overview of the complete VMD algorithm is obtained in Algorithm 1.

end for
Dual ascent for all

Hilbert transform
The original signal is decomposed into several IMFs by the VMD.The Hilbert transform is used to demodulate the IMF to determine its frequency and amplitude over time.The Hilbert transform of u k (t) is defined as Eq. ( 13).
Accordingly, the analytic signal is constructed as follows: where, a k (t) is the instantaneous amplitude and θ k (t) is the instantaneous phase, which are written as below: The instantaneous frequency can be further obtained by Eq. ( 17).
The Hilbert spectrum of u k (t) is defined in Eq. ( 18).
The marginal spectrum can be obtained by integrating H(ω, t) over time, as follows.
where, T is the signal length.

Quasi-bivariate variational mode decomposition
The paper proposes a quasi-bivariate VMD method to extract the flow structure at characteristic frequencies.In Figure 7, the quasi-bivariate VMD method's flowchart is depicted.The method includes two decompositions and one correction, which can effectively avoid the mode mixing problem.The particular steps are listed below.First, the flow field is divided into a series of slices along the x-direction, and the first VMD is performed on each slice to obtain IMFs.The mode number K is chosen with reference to the number of prominent peaks of pressure pulsations in the frequency domain.Secondly, the center frequencies corresponding to the obtained IMFs are checked.The mode number K and quadratic penalty term α are adjusted if the decomposition is insufficient or excessive.Subsequently, the VMD is executed again.Until neither exists, the IMF at the specified frequency is extracted.Then, the second VMD is performed on the extracted IMF to filter out the noise part.The obtained signal can be used for time-frequency analysis.Finally, the IMFs extracted from each slice are assembled along the x-direction.

Analysis of flow-induced force in the pump
The non-dimensional force coefficient (C F ) is defined: . Figure 8(a) depicts the orbital diagram of the radial force on the impeller.As the impeller rotates, the force vector shifts in the opposite direction.In one cycle, there are six peaks equal to the blade number.The radial force is phase-shifted at the fluctuating rotation speed in comparison to the stationary rotation speed.Accordingly, the trajectory of the radial force can be divided into lag and lead regions.However, both are similar in magnitude.In Figure 8(b), the radial force in the frequency domain is further analyzed.The primary frequencies of radial forces are the BPF and its harmonics.Fluctuating rotation speed hardly affects the amplitude at BPF but induces FF.The amplitude of the FF is proportional to the fluctuation ratio.
Figure 9(a) shows the pressure pulsation at the pump outlet under different fluctuation ratios, where the instantaneous rotational speed at γ = 0.013 is depicted specifically below for comparative analysis.The pressure   pulsation changes periodically at stationary rotation speeds.As the fluctuation ratio rises, the magnitude of the pressure pulsation gradually increases.In a cycle, pressure pulsation can be divided into two areas: the increased area of pressure pulsation induced by high rotational speed and the decreased area of pressure pulsation corresponding to low rotational speed.Both form a complete pulsation in one cycle.Figure 9(b) shows the pressure pulsation in the frequency domain.Pressure pulsation is primarily characterized by BPF and its harmonics.Pressure pulsations at FF begin to predominate when the fluctuation ratio rises.Similarly, the fluctuating rotational speed has no significant effect on the pressure pulsation at the BPF.
In this paper, standard deviations are adopted to define the pressure pulsation intensity C * p and the force pulsation intensity C * F in Eqs. ( 20) and ( 21), where Cp is the average of C p in one cycle and CF in the same way.C * P0 and C * F0 are the pressure and force pulsation intensities at stationary rotation speeds.
The pressure and force pulsation intensities at different fluctuation ratios are shown in Figure 10.Positive correlation exists between the fluctuation ratio and pressure pulsation intensity.In contrast, the force pulsation intensity is insensitive to the fluctuation ratio.
The force pulsation intensity at fluctuating rotation speed is slightly lower than that at stationary rotation speed.

Analysis of the temporal-spatial features of the flow field
The spatial and temporal distribution of the instantaneous pressure (excluding the average value) around the impeller outlet is shown in Figure 11.If the angular position is fixed, the ordinate shows the pressure changes over time; if the angle is fixed, the abscissa represents the distribution of the pressure around the impeller outlet.At stationary rotation speeds, the pressure distribution exhibits good periodicity in time and space.There are six high and low-pressure areas that advance with the impeller rotation in the circumferential direction.The pressure near the volute tongue varies dramatically.In contrast, the temporal-spatial distribution of pressure at fluctuating rotation speeds is more complex.A whole   high-low pressure area can be captured along the time direction and the corresponding high pressure area is marked.The fluctuating rotation speed is responsible for this situation.High rotational speed induces highpressure pulsations, while low rotational speed does the opposite.
The non-dimensional velocity coefficient (C v ) is defined: 12 shows the spatial and temporal distribution of the instantaneous velocity around the impeller outlet.Overall, the temporal-spatial velocity field is similar to the pressure distribution.However, the two differ in some ways.The effect of the volute tongue on the velocity field is mainly downstream.The change in velocity is significant at angular positions from 60°to 180°.Similarly, fluctuating rotation speeds induce corresponding changes in the velocity field.High rotational speed also leads to high velocity, but with phase lag.This is mainly attributed to the flow inertial.
Figure 13 depicts the distribution of absolute velocities in the mid-plane.The velocity contours for angle = 90°a nd angle = 270°are extracted at fluctuating rotation speeds, which correspond to the moments of maximum and minimum rotation speeds, respectively.For comparative analysis, the velocity contours for angle = 90°at stationary rotation speed are also extracted.The velocity at the trailing edge of the suction side of the blade is relatively high, which is typical of the jet-wake structure.The jet-wake structure is an important source of pressure pulsation (Wu et al., 2021).The strength of the jet-wake structure at high rotation speeds is significantly higher than that at low rotation speeds.The variation of the flow velocity is an intuitive representation of the fluctuating rotation speed, which affects the coupling between the exit angle of the flow and the tongue geometry.This impact can be indirectly evaluated by the relative height of the low-velocity zone at the volute outlet.It can be observed that the velocity distribution at the volute outlet in Case 3 is exactly between Case 1 and Case 2.
Further, the vortex structure in the pump is identified by the Omega method (Liu et al., 2016a).Figure 14 depicts the vortex structure at different moments under the fluctuating rotation speed.The threshold with Omega = 0.52 is adopted and the vortex structure is colored by the turbulent kinetic energy.The vortex shedding in the wake of the blade trailing edge impinges on the volute tongue and subsequently ruptures.A large number of vortex structures exist at the volute outlet.Another contributing factor is the flow separation downstream of the volute tongue.Depending on the rotating speed, the location of the vortex structure downstream of the volute tongue shifts.This is attributed to the change in the inlet angle of the flow at the volute tongue.

Analysis of flow-induced force in the pump
The pressure pulsation at the pump outlet is decomposed using the VMD.In Figure 15, the signal (SIG) is decomposed into a number of IMFs and a residual part.T is the impeller rotation period.The center frequencies corresponding to IMF1-IMF3 are the BPF and its harmonics.The central frequency of IMF4 is the FF.It is worth noting that the instantaneous amplitude variations of IMF1-IMF3 at fluctuating rotation speed are more significant than those at stationary rotation speed.
The spectrum of each IMF is shown in Figure 16.There is no mode mixing problem among the IMFs.A pair of sideband frequencies centered on the BPF and its harmonics can be captured at fluctuating rotation speeds.The higher the harmonics, the more pronounced the sidelobe peak.The difference between the upper and lower sideband frequencies and the center frequency is the fluctuation frequency.The sidebands are mainly attributed to the modulation between the non-stationary rotation speed and the rotor-stator interaction.
The instantaneous frequency and amplitude of each IMF are obtained by the Hilbert transform.The frequency variation of each IMF with time at fluctuating rotation speed is shown in Figure 17.The line's color indicates the relative magnitude of the instantaneous amplitude.The instantaneous amplitude from high to low is indicated by red to blue.The instantaneous frequency of IMF1-IMF3 fluctuates periodically.As the center frequency increases, the frequency pulsation becomes greater.The instantaneous frequency and amplitude are associated with the rotation speed.High rotation speed induces high instantaneous frequencies.Nevertheless, the high instantaneous amplitude appears in the region of the increasing rotation speed.In contrast, there is basically no frequency change in IMF4.The maximum instantaneous amplitude corresponds to the low rotational speed.
Figure 18 shows the instantaneous frequency and amplitude of IMF1 at different fluctuation ratios (excluding the mean value).The instantaneous frequency and amplitude do not change basically at stationary rotation speeds.As the fluctuation ratio rises, the pulsation range of the two increases.The instantaneous frequency  changes almost synchronously with the rotation speed.However, the instantaneous frequency lags behind the instantaneous amplitude.The pressure pulsation at the pump outlet stems mainly from the impingement of the wakes on the volute tongue.The high speed of the impeller blade approaching the volute tongue tend to cause high amplitudes of pressure pulsations.Similar conclusions are obtained in our previous study (Yu et al., 2022b).The reason for the phase lag can be attributed to the following, the instantaneous frequency mainly depends on the rotation speed, while the instantaneous amplitude is related to the RSI effect.The marginal spectrum represents the accumulation of signal amplitudes in the time domain.The marginal spectrum at fluctuating rotation speed is normalized, as shown in Figure 19.A single peak at the FF indicates no frequency fluctuation.In contrast, the BPF and its harmonics are continuous spectra over a certain range.As the harmonic increases, the range of the continuous spectrum also enlarges.Compared with the Fourier spectrum, the Hilbert marginal spectrum is able to eliminate the spurious frequency components.It can accurately reflect the local features of non-stationary signals.
The flow structure at the characteristic frequencies is extracted by quasi-bivariate VMD. Figure 20 shows the temporal-spatial pressure field around the impeller outlet under BPF.The impeller blade moves circumferentially as it turns, with the high-low pressure side moving forward.The pressure distribution in space is distorted by the volute tongue.At stationary rotation speed, the pressure shows a periodic pulsation with time.In contrast, fluctuations in rotational speed induce corresponding changes in pressure.The areas of increased and decreased pressure can be captured in the time domain.
The flow structure around the impeller outlet at the FF is shown in Figure 21.The pressure distribution in space is similar.In time, the pressure field varies synchronously with the rotational speed.Given the above, the pressure variation induced by fluctuating rotation speed in Figure 11 includes two forms.One is the flow structure at fluctuation frequency.The other is the change in the instantaneous amplitude of the flow structure under the BPF and its harmonics.The former is the main component.

Conclusions
In this paper, the hydraulic characteristics of centrifugal pumps at fluctuating rotation speeds are studied by numerical methods.The time-frequency features of the pressure pulsation are analyzed by VMD and Hilbert transform.A quasi-bivariate VMD method is proposed to analyze the non-stationary flow field.The mechanism of rotor-stator interaction at non-stationary rotating speed is revealed.The experiment verifies the correctness of the calculation procedure.The main conclusions are as follows: (1) In comparison to the stationary rotation speed, the fluctuating rotation speed causes flow-induced force under FF.The amplitude and fluctuation ratio are positively correlated in a certain range.The effect  Briefly, the study makes clear the dynamic characteristics when the rotation speed changes periodically with small amplitude around its mean value.And special interest lies on the time-frequency features of pressure pulsations.The relevant conclusions could provide an essential reference for vibration control and fault diagnosis of centrifugal pumps.

Improvements and future direction
However, there are several imperfections in the paper.First, only studies under rotation speeds subjected to sinusoidal variations are conducted.Second, the object of this paper is a centrifugal pump with a low specific speed.If possible, extending conclusions obtained in the paper to centrifugal pumps with medium and high specific speed.The next phase of work will focus on studying more complex fluctuation patterns and perfecting the differences between different pump types.

Disclosure statement
No potential conflict of interest was reported by the author(s).

Figure 6 .
Figure 6.Comparison of experimental and simulated results: (a) pump performance; (b) pressure fluctuation.

Figure 9 .
Figure 9.Effect of fluctuating rotation speed on pressure pulsation: (a) time domain; (b) frequency domain.

Figure 10 .
Figure 10.Effect of fluctuating rotation speed on pressure pulsation and radial intensity.

Figure 11 .
Figure 11.Time history of the pressure pulsation around the impeller outlet: (a) stationary rotation speed; (b) fluctuating rotation speed.

Figure 12 .
Figure 12.Time history of the velocity pulsation around the impeller outlet: (a) stationary rotation speed; (b) fluctuating rotation speed.

Figure 14 .
Figure 14.Vortex structure comparisons in different moments.

Figure 16 .
Figure 16.Decomposition results of pressure pulsation in the frequency domain: (a) stationary rotation speed; (b) fluctuating rotation speed.

Figure 18 .
Figure 18.Effect of fluctuating rotation speed on instantaneous amplitude and frequency.

Figure 21 .
Figure 21.Time history of pressure pulsation under FF.

b 1 =
outlet width of impeller, mm C p = pressure coefficient C F = force coefficient C v = velocity coefficient C * p = pressure pulsation intensity C * F = force pulsation intensity Cp = average pressure coefficient CF = average force coefficient C * P0 = pressure pulsation intensity at stationary rotation speed C * F0 = force pulsation intensity at stationary rotation speed D 1 = diameter of inlet, mm D 2 = diameter of outlet, mm D 3 = outside diameter of impeller, mm D 4 = base circle diameter of volute, mm F = force, N f n = rotational frequency, Hz f 1 = fluctuation frequency of the rotating speed, Hz g = gravity, m 2 /s H = head, m H r = pump design head, m K r = type number K = mode number N = signal length P = power, W p = pressure, Pa Q = flow rate, m 3 /s Q r = design flow rate, m 3 /h T = impeller rotation period t = time, s v = velocity, m/s y = vertical height of the first boundary layer to the wall, m

Table 1 .
Geometry and performance parameters.