CFD configurations for hydraulic turbine startup

This paper presents various numerical setups for modelling Francis turbine startups involving moving meshes and variable runner speed in order to help define best practices. During the accelerating phase of the startup, the flow is self-similar between channels, thus making single sector configuration appropriate. Adding the draft tube improves the results by allowing pressure recovery midway during in the startup. At the speed no-load regime, a rotating stall phenomenon occurs and can only be capted with the full runner included in the simulation. Comparison with experimental data, such as runner speed and strain gauge measurements, generally shows good agreement.

Historically, hydraulic turbine structural analyses were focused on the static and dynamic loading on the blade near best efficiency operating point. Current-day researches are pointing out that low cycle fatigue transient events such as startups can greatly influence their life expectancy [1,2,3]. In order to optimize maintenance and operation of units while limiting the cumulative damage, electrical utilities need to understand the full impact of start/stop cycles. However, reproducing these transient events numerically remains a challenge due to the lack of general guidelines and the complexity of the physics involved.
A first attempt at modeling the startup phase of a Francis turbine from rest to speed no-load raised a lot of questions on the influence of the numerical setup and the required components to be included in the simulation [4]. In this paper, the results obtained from different configurations and parameters are analyzed and compared to the available experimental data in order to help define the best practice for hydraulic turbine startup simulations including the speed no-load regime. The Francis turbine under investigation is shown in Figure 1. It has a semi-spiral casing, 24 guide vanes, 13 runner blades and an elbow draft tube. It operates under a nominal head of 23 m, has a diameter of 5,4 m and a nominal output of 48 MW.
Two startup scenarios, presented in the right part of figure 1, are studied here. The first one involves an opening of the guide vanes up to 10 degrees where they stay until the group reaches a speed of about 70 rpm. The speed governor then takes control of the group and brings the machine to synchronization speed of 75 rpm. The startup 2 follows the same logic but with an opening up to 16 degrees. It is considered more aggressive since the opening is higher and consequently the plateau is shorter. While startup2 exhibits an overshoot in velocity, both scenarios have the same opening slope and require the same amount of time to reach synchronization. Total pressure corresponding to the nominal head was specified at the inlet along with a pressure boundary condition at the outlet so that the flowrate was calculated at each timestep. Transient rotor stator interfaces were used between the rotating and stationary parts, involving some scaling between the components when different pitches were involved. The numerical setup used two junction box routines with CFX 14.5 to perform simulations as illustrated in figure 2. The first one updated the mesh around the guide vane at each timestep according to a prescribed motion and the second updated the runner velocity based on the calculated torque. A companion script, running alongside the CFD solver, was sometimes used to pre-generate and clean up the meshes. More details about boundary conditions and the numerical setup can be found in [4].  The first simulations were performed to check the influence of mesh density and numerical schemes on the results. These tests included 2 nd order scheme for spatial discretisation, high resolution turbulence, 2 nd order backward scheme for updating runner velocity and doubling of the mesh size presented in table 1. All these simulations were performed with a simplified setup including only a sector of the distributor and runner (DR). As can be seen in After determining that numerical sensitivity was low with respect to the quantities of interest, the next step was to investigate which components were to be included in the simulation in order to capture the right physics of the startup. Questions about whether or not to include the draft tube, the complete runner and the spiral casing were raised since each additional component increases the numerical cost of the simulation. The tested configurations can be seen in figure 3 from the simplest to the more complete.   Table 3 summarizes the results and provides more information on the computational cost (run time, number of cores and disk space). The draft tube can be seen to have a significant impact on both the maximum and final velocities. Use of the complete geometry raises the numerical cost by two orders of magnitude when compared to the simplest configuration.   (1), becomes significant around 10 seconds. At this time, the combination of the flowrate and runner velocity provides an inlet profile to the draft tube which has limited swirl (figure 4 right), allowing kinetic energy to be converted into pressure by the draft tube. This in turn increases the apparent head seen by the turbine and thus affects the runner torque and acceleration.  figure 5 for the DRA configuration. The agreement is very good except at the beginning of the startup. This could be explained by the biphasic nature of the flow that was neglected in the calculation. Initially, the runner is not filled with water as it is above the tail water level, so calculated torque is probably overestimated. Possible air ingestion during the startup could also influence the results. As the curves for the DRA360 and DRAS 360 are similar, they are not shown here.
As stated previously, the aim of this study was to evaluate the impact of startup on turbine life expectancy. Before performing cumulative damage analyses, the pressure loading has to be converted into mechanical stress. There are multiple models that could be used for fluid-structure interaction, including fully coupled two-way interaction. Here, a simple quasi-static model was used. This assumes that the structure is not deformed significantly by the flow and that the timescale of the deformation is less than the fluid timescale. For the duration of the startup, a pressure load was extracted every second and used as a load for static structural analysis in Ansys Mechanical. Gravity, centrifugal forces and rigid coupling with the turbine shaft completed the boundary conditions for this model. The structural mesh was made out of 133K tetrahedral elements for a periodic sector. Since a strain gauge was installed on one blade, which provided experimental data for comparison purposes, it was reproduced with a surface bounded element in the numerical model ( figure 6). This region was previously identified as the most constrained one. A comparison of the deformation along the principal axis of the strain gauge showed reasonable agreement for the two startup scenarios. High frequency fluctuations could not be reproduced by the structural model (figure 7) since it was quasi-static. The rapid increase of strain at the beginning is consistent with the CFD acceleration that is overestimated for this period and abstraction of the biphasic nature of the flow in the simulation is suspected to be responsible for this. Some of the discrepancies between experimental and numerical results could also be related to strain gauge position that is not precisely known.
During most of the startup, the simulations including the full runner and distributor showed that the individual dynamic torque on the blades remarquably collapse on the average curve (figure 8), meaning that the flow is self-similar between channels. Small oscillations appear around 15s, correlated with the maximum pressure recorded on the leading edge, but are rapidly damped. However, large torque oscillations appear quickly after reaching the speed no load regime. This was rather unexpected since this turbine did not report any synchronization problem. Initially, these fluctuations were thought to come from the draft tube but the flow in this region was very stable.  Inspection of the runner confirms that the flow was self-similar between the different channels until 25s (figure 9, left). At this time, two counter rotating vortices were present at the leading edge of the turbine. However they rapidly started to oscillate, leading to an unstable state where 6 or 7 secondary flow cells appeared in the runner (figure 9, right).
For the 25 to 40s period, the numerical signal was compared with the one from an experimental pressure probe located at the leading edge (figure 10). The simulated frequency of 1.12 Hz in the rotating frame was also found on the experimental signal. This corresponds to 70% of the runner frequency and seems to be consistent with the rotating stall reported by Hasmatuchi [5]. No clear explanation has yet been found to explain why the numerical simulation overestimates the amplitude of this phenomenon. It is possible that the k-epsilon turbulence model that was used, which is known to be diffusive, could only allow for major structures of the flow to exist. Also, as one of the countermeasures mentioned in the literature for such phenomena [6,7] is to offset some guide vane positions, it is possible that small geometrical discrepancies between the blades and the guide vanes in  the real turbine could be assimilated to mistuning and could therefore reduce the fluctuation amplitude. This is obviously not the case in the numerical model where all the sectors are identical. Another possibility is that air pockets located in the runner hub are brought back to the leading edge area by the large pumping vortex at the trailing edge and add some damping to the system. Figure 11 provides an overall view of the phenomenon. On the left side, the radial velocity contours and vectors clearly show the flow ejection originating from stall cells. The right part shows that the stall cells are not fully contained within the runner but extend to the vaneless space downstream of the guide vanes. Here the high total pressure zones are regions for which the turbine had transferred energy through the pumping motion occurring on the suction side of the runner channel.  Most analyses of rotating stall phenomenon are presented using a constant guide vane opening for determining whether or not the machine will be in an unstable regime. For the 27 to 40s period, the opening angle is almost constant since the speed governor has taken control of the turbine prior to grid synchronization. A similar analysis can then be performed. The characteristic curves for this period ( figure 12) show an inflection point where a positive slope can clearly be seen. This criterion is often mentioned as the onset of instability for pump turbine [6]. While the region that showed limited fluctuations, around 15s, is also highlighted in this figure, nothing of particular interest seems to be happening based on the graph.
The rotating stall occurring at the leading edge of the runner involves an interaction between two major counter-rotating vortices. Figure 13 shows the Z component of vorticity in the runner reference frame. Based on this illustration, the pumping flow from the suction side is seen to be the driving mechanism that generates vortices and partially blocks channel or lets the flow through. At first there S-Shape Characteristic curve. The rotating stall phenomenon quickly appears after the flowrate and torque curves show a positive slope at the end of the startup.
Rotating stall at the runner leading edge colored by the vorticity Z-component (red is conterclockwise and blue is clockwise). This figure shows how stall cells evolve from one channel to the other by alternatively blocking/allowing the flow through. is an ejection of the flow originating from the suction side of the runner into the vaneless space (A). This blocks the passage and pushes the stall cell to the adjacent channel to the left (B and C). The flow returns to the channel (D and E) until a new stall cell arrives and is strengthened by the pumping on the suction side (F), after which the whole process starts again.
Having a system that can reproduce all the machine's transients can be very helpful for predicting a turbine's life expectancy. To model the startup, a configuration that includes one distributor sector, one runner blade and the draft tube is enough to capture the proper behaviour since the flow is selfsimilar between the hydraulic channels. This is therefore the recommended configuration for capturing the main physics of the startup. The inclusion of the draft tube improves the results by allowing pressure recovery midway in the startup when a proper combination of flowrate and runner velocity provides a swirl-free profile at the runner outlet. A quasi-static structural analysis showed acceptable prediction of the average quantities. Although these simulations were not aimed at studying the speed no-load regime, a rotating stall was captured at the end of the startup. To reproduce the phenomenon, the complete runner must be included which substantially increases the numerical cost of the simulation. During this regime, about six or seven rotating recirculation cells were identified in the vaneless space. An analysis of the characteristic curve for this turbine correlates the onset of fluctuations with an inflection point. This can be associated with reported unstable S-shape characteristic curves. The simulated frequency was found on the experimental signal, although with less amplitude.
More data with different types of turbines will be necessary to further validate the methodology used to simulate transient processes. Experimental campaigns on prototype runners, designed with transient operation in mind, are actually underway and should provide adequate data for this purpose. Investigations to improve the structural coupling and analyses are also currently ongoing.