Simulation of air admission in a propeller hydroturbine during transient events

In this study, multiphysic simulations are carried out in order to model fluid loading and structural stresses on propeller blades during startup and runaway. It is found that air admission plays an important role during these transient events and that biphasic simulations are therefore required. At the speed no load regime, a large air pocket with vertical free surface forms in the centre of the runner displacing the water flow near the shroud. This significantly affects the torque developed on the blades and thus structural loading. The resulting pressures are applied to a quasi-static structural model and good agreement is obtained with experimental strain gauge data.


Introduction
To better characterize turbine residual life, a better understanding of the hydraulic load on the runner blades during damaging events is required. Studies reveal that turbine fatigue damage arises mostly from low occurrence, high amplitude cycles related to transient events such as startups and runaways [1].
In 2014, Hydro-Québec performed an extensive measurement campaign involving strain gauges and accelerometers on the runner blades along with many parameters on the fixed part of one of its machines [2]. The purpose was primarily to assess structural behavior under different operating conditions and to optimize the startup sequence to minimize transient stresses. Part of the instrumentation on the runner blades with protective coating can be seen in figure 1. Despite being very valuable, the costs of on-site measurements along with their punctual nature strongly suggest that these data be supplemented with numerical results. Therefore, validation of numerical tools was also considered as an important part of the process.
Building on previous transient Computational Fluid Dynamic (CFD) work on a Francis turbine, a setup allowing guide vane motion and variable runner speed was readily available [3]. It was decided to numerically reproduce one startup and one runaway scenario for a propeller turbine to further validate the methodology. As it will be shown later in this paper, early tests involving monophasic water simulation (figure 5) failed short to reproduce the expected behavior. After verifying the numerical setup, differences with measured values were deemed important enough to revisit some simplifying assumptions. In this process, aeration was soon identified as one of the possible sources of error since turbine operation was very noisy at speed no load regime due to air admission. In hydraulic turbines, air admission through the upper cover or the turbine shaft can play many roles. It can be used to increase the dissolved oxygen at the tail water, prevent water column separation and dampens/delays dynamic phenomena and vibrations occurring during speed no load, deep part load and sometimes even full load regimes [4]. Although turbine behavior is clearly affected, not much is really understood about what is going on inside the flow. To say that aeration is not well understood is probably an understatement. Part of this knowledge gap is explained by the difficulty to reproduce air admission with scaled models, although the same also applies to numerical models. Since air injection is not usually operational at optimal efficiency, these systems are often overlooked by designers and maintenance staff who mostly see them as an insurance policy for mitigating unwanted dynamic problems. In fact, literature on turbine aeration is scaringly thin for such a common system. Most papers focus on the impact of aeration for dissolved oxygen content for environmental purposes [5][6][7] and with the model homology [8]. Only one CFD study that includes multiphase water/air flow within the turbine was found in the literature [9]. However, in the latter paper, air injection only had a limited influence on flow since normal operation was simulated.

Turbine description
The turbine under study operates under a head of approximately 25m in a run-of-river power plant in Québec. After a short intake, the flow is oriented towards the 24 guide vanes through a semi-spiral casing. Energy is extracted by a six-blade propeller, after which the flow returns to the river through an elbow draft tube with two outlets.

Startup and runaway conditions
Guide vane motion is one of the primary inputs in transient simulations. Since an angular encoder was installed on one of the guide vane, there was no need to convert the servomotor stroke into angular motion. The chosen startup scenario was the default one for this powerhouse. As can be seen on figure 2, after a fast opening, the guide vane reaches a plateau at 30% of the maximal opening. This position is maintained until the runner velocity reaches about 90% of the nominal speed. The guide vanes then fold back to an opening of 15% and the speed controller brings the turbine to synchronous speed.
The runaway represents a load rejection while the turbine operates at full load. This can happen whenever connection to the grid is lost or because the turbine safety system detected a condition that requires immediate shutdown. After detecting the overspeed, guide vanes rapidly close in order to limit the maximal velocity, then a slower slope is used to limit the overpressure in the spiral casing. In this case, after a full closing of the guide vanes, the speed controller regains control over the turbine and brings it back to synchronous speed. This is more representative of a grid connectivity problem than a turbine mechanical issue since the latter commands a complete shutdown.

Meshes
Considering the simulation setup, one of the main challenges is to cover the full range of the guide vane openings for the runaway scenario by using a single structured topology. Some compromises are obviously required to maintain good overall grid quality while covering a span of almost 50 degrees ( figure 3). This includes stretching of the mesh in the distributor channel at maximum opening. A small extrusion on the distributor ceiling is included to represent the air intake. The exact shape of the air ducts is not reproduced since a porous region is used instead. Four sectors of the distributor are included in the simulations to match the pitch of the runner domain which only contains a single blade. The runner blade tip gap is modeled with 9 nodes. The cavitation lip is omitted. Finally, the draft tube is also included in the simulation. More details about the meshes are given in table 1.

Numerical parameters
The measured head is imposed with a constant total pressure at the distributor inlet and an average static pressure at the draft tube outlet. The air intake is defined as an opening at the atmospheric pressure. Interfaces between the rotating and stationary domains are defined as transient rotor stator. The water/air mixture is modelled with a homogenous model and an interface compression scheme is used to obtain a sharp interface between the two phases. Gravity is included in the model to allow for the air to go up through buoyancy. The reference pressure is set at the draft tube tail water level. Notwithstanding the hydrostatic pressure gradient, air is defined as an incompressible fluid. This is considered realistic enough for the air pocket around the runner. Cavitation is not modelled. Two porous regions are used to replace some complex physics. The first one allows for the 0 to 1 degree opening of the guide vanes and is present only during 0.5s of the startup. Its resistance coefficient was adjusted empirically to allow for a linear flowrate increase. The second porous region is used to control the amount of air entering the turbine and represents the losses in the air duct system. It is the only user-controlled parameter in the simulation since the guide vane motion and the head are imposed. Its value was set using steady simulations in order to obtain the synchronous speed at measured speed no load opening. The calibration of this coefficient is not based on experimental results and directly relates to the size of the air pocket around the runner.
A blended numerical scheme with 75% second-order accuracy is used for pressure and velocity equations. A standard k-epsilon model with scalable wall functions together with a first order scheme is used for turbulence. Time marching is done using a second-order backward scheme. The timestep was found to be one of the most important parameters to control with regards to calculation stability. The first instants (5s for startup and 10s for runway) required a timestep of 5E-4 s. Ather that initial period, it was increased to 2E-3s to accelerate the run. This allowed the rms Courant number to remain below 5 and corresponds to a runner rotation of about 1 degree/timestep at synchronous speed. Although residual convergence remains a constant challenge in biphasic simulations, engineering quantities stabilized within 10 iterations of each inner loop before moving to the next timestep. As mentioned previously, the mesh motion of the distributor and the variable speed of the runner were handled by junction box routines within Ansys CFX 16.0[3].
The initial condition for the startup is the measured leakage flow between the guide vanes. In this condition, the runner is mostly dewatered since it is located above the tailwater level. Figure 4 shows this condition with a water volume fraction of 0.5. For the runaway condition, the initial condition is simply the maximal load operating point.   Figure 5 shows the results of transient simulations both with and without the presence of air. Two problems are found in the monophasic water simulation. First, there is a too rapid acceleration of the turbine at the beginning, but more importantly there is an offset of more than 10% in the final velocity. The inclusion of air in the simulation resolves the two problems together. The slower startup is explained by the time required to fill out the vaneless space between the guide vane and the runner which delays the torque increase. For the final runner velocity, it is not surprising to see that there is a match in the final velocity since the air intake porous losses are calibrated to obtain such a result. However, the fact that the runner velocity evolution is properly captured means that turbine forces/inertia ratio is well reproduced. Since the water flowrate is very similar with or without air, the difference between the two flows must be found elsewhere. The speed no load regime is a stable operating point where the turbine operates globally at no torque. In order for the net torque to be zero, part of the blade must be pumping while the other part produces power. It is important to keep this in mind since when aeration is present, an air pocket naturally forms near the hub where stands a recirculation in the monophasic simulation. This pocket rapidly grows in size to occupy a significant volume in the runner and affects the water flow distribution on the blades by displacing it towards the shroud. By concentrating the flow at a higher radius on the blade, it naturally increases the torque produced by the runner, thus explaining why rotation velocity is higher in presence of air than without it. A view of the pressure loading at the speed no load regime is presented in figure 6. In the right part, one can clearly see that a constant pressure region is present within the grey air pocket. Despite the effect of gravity, the air cavity has an almost vertical surface. Figure 7 shows what happens during startup with the help of the water volume fraction in a meridional view. Initially, since the runner is above the tailwater, a large air region is present behind the guide vanes. At t=2s, while the guide vanes are still moving, a vertical water front advances quickly to fill out the volume. At t=6s, the runner is completely surrounded by water and most of the air has now moved into the draft tube. Between 8 and 10s, a new air pocket, fed by the air intake, starts to form around the hub. At t=12s, some of the air continues on its way to the draft tube exit, while the rest is caught within the recirculation region below the hub and returns to the runner. This largely explains the rapid growth of the air pocket than can be seen at the final speed no load regime (15s).

Runaway
For the runaway condition, overall results are also qualitatively good as shown in figure 8. The maximum velocity of 150% is well reproduced both with and without air. However, while it took about 5 seconds to reach that maximal velocity in the measurements and it was maintained for a few seconds, it is reached earlier in the simulation and no such plateau is observed. Here, it is suspected that cavitation could explain the different behaviour. Figure 9 shows the regions with a negative absolute pressure near the hub and on the pressure side leading edge highlighted in red. These regions indicate where cavitation should be observed and where blade loading should be affected accordingly.
In figure 8, differences between the setups with and without water become more visible in the slowdown process when air starts to enter massively into the machine. Once again, the monophasic water simulation proves to be too reactive. Although air entrainment helps to approach the experiments, some mismatch is still visible in the slowdown process. Looking at the water volume fraction evolution during runaway (figure 10), it is clear that the air pocket around the runner hub comes solely from the upstream air intake. Initially, there is no air in the machine. The turbine starts breathing prior to reaching the maximum speed between 2 and 4s (2.9s), probably helped by the low pressure region highlighted earlier on the hub. Growth of this region continues until the turbine reaches the speed no load regime near 15s. In the last frame at 30s, some air bubbles are visible which are being shed from the main air cavity and have drifted into the draft tube. This explains why an air supply is continuously required for this regime to be stable. Figure 10. Contours of water volume fraction during runaway (air = red, water = blue)

Structural analysis
As mentioned previously, the measurements performed on the runner blades only included strain gauges and accelerometers. While global results such as rotating speed were in good agreement, fluid numerical results had to be transposed to a structural model in order to obtain a more direct comparison. A numerical model of the runner containing the six blades with the runner cone was therefore built. Gauge locations were individually measured since they differ slightly from blade to blade and this partly explains why a complete structural runner was used. Fillets on the hub were included but the cavitation lip was omitted since sensitivity tests showed it has a very limited impact on the blade strains. The finite element model mesh was made of 850k elements and 1.4m nodes ( figure 11, left). Additional refinement near the gauge sites is also visible. The transient calculated pressure presented previously was duplicated and imposed on the blades as a quasi-static loading with a period of 1s (figure 11, right). Rotation speed was adjusted accordingly at every timestep. Gravity and coupling with the turbine shaft completed the setup. Structural analyses were performed with ANSYS Mechanical using default parameters. Once again, the results were very good, as can be shown in figure 12, which displays the measured deformation on a site located near the hub (site B, figure 1) against the numerical results. This site was located on one of the high-stress regions identified from the static analysis during the design phase. The amplitude of the deformation, significant for the low-cycle fatigue, was accurately captured. This simplified setup thus allows the startup to be optimized using numerical simulations. High-frequency fluctuations are not expected since a quasi-static structural analysis was performed.

Conclusion
The initial goal of this project was to validate the transient numerical setup developed previously on a Francis turbine with a propeller and to obtain stresses for comparison with on-site measurements. Startup and runaway scenarios were chosen since they represent some of the most challenging and damaging conditions, thus allowing nearly all possible operating regimes to be simulated. The first simulations with water as the only fluid led to unsatisfactory results. Those were significantly improved when air was included. While aeration is often present in hydraulic turbines to mitigate dynamic phenomena, its impact on transient behaviour with numerical CFD tools was not studied until now. A major influence of aeration can be expected whenever air is injected in a recirculating region since this leads to the rapid growth of an air cavity thus affecting the water flow topology. Such regions are common in off-design and transient operations. Performing transient biphasic simulations in a hydraulic turbine is by no means trivial, but with a careful selection of numerical parameters, very useful information can be obtained on turbine behaviour. A better understanding of the general flow evolution during transients is essential to accurately assess the residual life of the machine. The high-amplitude stress cycles on the machine can be obtained numerically. This opens the door to numerically develop scenarios for limiting transient stresses on the runner by modifying guide vane motion together with turbine aeration.
Despite the good results obtained here, a lot remains to be done to further improve the model. Cavitation in particular would be worth including since it seems to play an important role during the runaway process. A better assessment of the air flowrate is also required to accurately characterize its influence during transient events. This might require modelling of the air ducts and a better understanding of the interphase mass transfer between air and water to prevent the use of a homogenous model and the free surface interface compression scheme. Finally, while much can be done with numerical tools, some improvements must also come from experiments or operationally. This involves at least knowing when aeration systems are operational, but also how much air is admitted.