State of the art in numerical simulation of high head Francis turbines

This is an O Abstract. The Francis-99 test case consists in a high head Francis turbine model, which geometry together with meshes and detailed experimental measurements is freely available at www.francis-99.org. Three workshops were initially planned to exchange experience on numerical investigations of the test case concerning steady state operating conditions, transient operating conditions and fluid structure analysis. The first workshop was held in Trondheim, Norway in December 2014. Some results of the 14 contributions are presented. They are concerned with the influence of the near wall space discretization and turbulence modelling in order to capture hydraulic efficiency, torque, pressure and velocity with a good uncertainty at three operating conditions.


Introduction
Computational fluid dynamic (CFD) is used to study the flow inside hydraulic turbines since some decades, Keck [1].Numerical modelling of hydraulic turbines is challenging because a specific modelling approach applied to investigate an operating condition does not obligatory work for another and several flow phenomena may appear simultaneously.An open test case allowing the researcher to interact on such issues is necessary to develop CFD technique applied to hydraulic turbines.A test case named Francis-99 was created.The Francis-99 test case consists of a high head Francis turbine model, which geometry together with detailed experimental measurements are freely available at www.francis-99.org.Three workshops were initially planned.
The first workshop aimed to determine the state of the art in simulation of high head Francis turbines under steady operations: part load (PL), best efficiency point (BEP) and high load (HL).It was held in Trondheim, Norway, in December 2014.These paper summaries some numerical results presented in the first workshop.All the papers presented at the workshop are freely available in the Journal of Physics, Volume 579.

Model test rig
Experimental measurements were conducted on a model Francis turbine test rig.The hydraulic system is capable to develop a head around 12 m for the open loop and 100 m for the close loop.The open loop and close loop systems were used for the pressure and velocity measurements, respectively.Head across the turbine was maintained constant for both types of measurements.The operating head and discharge were 12 m and 0.2 m 3 s À1 at the BEP, respectively.The turbine includes, spiral casing, distributor with 14 stay vanes integrated inside the spiral casing and 28 guide vanes, runner with 15 blades and 15 splitters, and an elbow type draft tube.The guide vanes can be operated from 0°to 14°.Minimum angle of 4°and the maximum of 14°correspond to the part load and full load operating conditions, respectively.

Instrumentation
The model Francis turbine was instrumented with sensors to perform model test according to IEC 60193 [2].Additionally, six pressure sensors were mounted at different locations inside the turbine.Figure 1 shows the locations of the pressure sensors.A miniature pressure sensor VL01 was flush mounted in the vaneless space to acquire pressure pulsations developed by the rotor stator interactions.Two miniature pressure sensors P42 and P71 were mounted on the blade pressure side and trailing edge, respectively.The other miniature pressure sensor S51 was mounted on the blade suction side.The remaining two pressure sensors DT11 and DT21 were mounted on the wall of the draft tube cone.Laser Doppler velocimetry (LDV) and two-dimensional particle image velocimetry (2D PIV) were used for the velocity measurement in the draft tube cone, Bergan, 2015 [3].

Turbine performance
Steady state measurements on the model Francis turbine were carried out to determine the turbine performances.The maximum hydraulic efficiency was around 93.4% at the BEP, n ED = 0.18, and Q ED = 0.15 [4].Three operating points were investigated for the Francis-99 workshop, part load (PL), BEP, and high load (HL).The guide vanes angular positions for the PL, BEP and HL correspond to 3.9°, 9.9°, and 12.4°, respectively.The discharge values were 0.071, 0.203, and 0.221 m 3 s À1 for the PL, BEP and HL, respectively.

Geometry and mesh
The three-dimensional geometry was prepared using available AutoCAD drawing of the complete turbine.The geometry includes spiral casing, distributor, runner, and draft tube.The geometry was prepared using the CAD software ANSYS ICEM CFD.A hexahedral mesh of the complete turbine was prepared using ICEM CFD.Threedimensional blocking technique was applied to create hexahedral elements in the computational domains.The turbine was divided into three domains, spiral casing with distributor, runner, and draft tube.Total nodes were around 12 million in the complete turbine.Over 5 million nodes were created in the runner and around 4 million nodes in the distributor.The mesh in the draft tube was relatively coarse, around three million nodes.Fine mesh close to the boundaries and in the complex passages was created.The pre-mesh blocks were provided to the researchers, which allowed creating finer mesh.
3 Results from the first workshop

Mesh
Some of the research groups used their own software and techniques to create meshes.Meshes with different densities were created ranging from 3 to 48 million nodes.Three modelling approaches have been applied: (1) complete turbine modelling, (2) component modelling, and (3) passage modelling.In the complete turbine modelling, all the components, i.e., spiral casing, stay vanes, guide vanes, runner with blades and splitters, and draft tube, of the Francis turbine have been modeled.In the component modelling approach, runner and draft tube have been modeled.In the passage modelling approach, a passage of distributor and runner blade has been modeled.Majority of the simulations have been performed using the complete turbine modelling approach.The maximum mesh density in the complete turbine was around 48 million nodes, Čelič [5].It was challenging to create a good quality mesh in the turbine because the guide vanes and blades have complex geometrical profiles.
Nicolle and Cupillard [6] investigated the influence of the y+ values and mesh resolution near the boundaries.The y+ value and mesh resolution have an effect on flow field.Using shear stress transport (SST) model, the runner output power decreased by nearly 10% as the y+ value decreased from ∼60 to a value below 1.The result obtained without wall function was closer to the experimental one pointing out the uncertainty in using wall function to determine losses in hydraulic turbines.
Meshes were created in the labyrinth seals to investigate the flow leakage losses in the turbine by some of the researchers: Amstutz [7], Jošt [8] and Čelič [5].Close to the labyrinth seals wall, fine mesh was created for improving the velocity profile from the solid boundary by Čelič [5], who used 63 million nodes.The volumetric loss due to the seals needs to be quantified especially at PL as its relative amount becomes significant.

Hydraulic efficiency
Most turbulence models overestimated the hydraulic efficiency for all three operating points; PL, BEP, and HL.The prediction was best at the BEP and HL conditions; within ±3%.The challenge was to model the flow field at the PL condition.Reynolds stress model (RSM), SST, and zonal large eddy simulation (ZLES) models show small difference at the PL condition.The efficiency was obtained with variation up to 12%.The SST and ZLES require fine mesh, which increases computational cost.The RSM shows reasonable result between computational power and required numerical accuracy.For the computation of the overall turbine performance, application of RANS models may be reasonable but these models are incapable to resolve the flow unsteadiness in the draft tube, Jošt [8].The discrepancy at PL was mainly attributed to the labyrinth seals, which contribution to the efficiency was not accounted for.Čelič [3] simulated the labyrinth seals and got 2% difference for the efficiency at PL.
Figure 2 shows the comparison of the experimental and numerical values of the hydraulic efficiencies at the PL, BEP, and HL conditions, Stoessel [9].The best performance of the numerical models is observed at the BEP and yield an error of 1.16%, the difference to the experimental data increases to 2.72% at high load and 10.67% at part load.The hydraulic efficiency obtained with steady state simulations show larger difference than the results obtained with unsteady simulations.
The hydraulic efficiency may be sometimes misleading because it is an integral value relating several variables and does not provide any detailed information about individual variables, Trivedi [3].Based on the results presented in the workshop, the difference in the quantities such as head, discharge, and torque was larger than the difference in hydraulic efficiency at the corresponding operating condition compared to the experimental data.The torque and head values were over-predicted by 13.9% and 19.5% at the BEP, respectively Lenarcic [10].The hydraulic efficiency showed a difference of À3.9%.Similarly, at the PL, the torque, head and the hydraulic efficiency were overpredicted by 28%, 14.1%, and 8.7% than the corresponding experimental values, respectively.These results indicate a difficulty for the turbulence models to predict accurately the losses mainly due to the use of wall function.

Pressure loading
The results indicate that the average pressure variation in the turbine might not be considerably influenced by the mesh density, discretization scheme, and turbulence model.One may use RANS models with extended wall function and relatively coarse mesh if the boundary layer and turbulent resolution do not matter significantly.However, the RANS models have some limitations to capture the unsteady flow phenomena in the draft tube including unsteady pressure pulsations, Nicolle [6].The SST model was unable to capture the unsteady pressure pulsations in the draft tube as captured by the SAS-SST model.
Figure 3 shows instantaneous pressure variation in the vaneless space at the BEP, Buron [11].The selected turbulence model, Spalart-Allmaras, was able to capture the flow unsteadiness in the vaneless space and the runner, not presented here.Spectral analysis of the pressure pulsations was also carried out by Stefan [12].Figure 4 shows a fast Fourier transform analysis (FFT) of the signal acquired at the vaneless space location VL01 during PL operating conditions with FLUENT and OpenFOAM   solvers.Both Solvers have underestimated the harmonic frequency of 15 by 50%.At the BEP, both solvers show difference less than 10% but, as compared to the experimental value, the amplitude of blade passing frequency (f n = 30) are over estimated by 50%.The numerical solvers failed to capture harmonic frequency of 15 in the vaneless space at the BEP, observed in the measurement.The pressure amplitudes were found to be sensitive to the mesh density, time step sizing, turbulence modelling, discretization order, and boundary conditions; Čelič [5], Nicolle [6], Amstutz [7], Jošt [8], Buron [11], Stefan [12] and Minakov [13].

Velocity distribution in the draft tube
Most turbulence models were applied to investigate the velocity field and vortex rope at the off-design conditions, i.e., PL.Some of the simulations were able to capture the vortex rope; Minakov [13] and Jošt [8] whereas others presented steady state behavior; Stoessel [9], Stefan [12] and Yaping [14].Two equation RANS models have shown difficulty to capture the flow unsteadiness, as expected.
Figures 5 and 6 show velocity profiles at the upper (L1) and lower (L2) sections of the draft tube cone during PL operating condition, Mossinger [15].Numerical simulations were carried out with three turbulence models and compared with the provided experimental data.All the turbulence models show almost similar trend at both sections.The standard k-e model shows difficulty to capture the tangential velocity at the upper section.Overall, SST and RSM models overestimated the axial velocity and standard k-e model underestimated the velocity at both the sections.Interestingly, these turbulence models are not able to capture accurately the tangential velocity at the lower section.Further, near the wall, small variation are seen, where the measurement shows a velocity drop at the upper section.
Numerical studies on different modelling approaches indicate that the one passage modelling and inclusion/ exclusion of spiral casing with distributor vanes does not have significant influence on the draft tube flow, Mossinger [15].Furthermore, improving the mesh density has no effect on the axial velocity at the part load operating condition; Čelič [5] and Jošt [8].However, small improvement in the velocity profile at the BEP and high load operating conditions may be observed with mesh refinement, Minakov [13].Flow modelling in a draft tube is a challenging task and expensive.One may simplify the draft tube flow modelling based on flow conditions to be investigated or interested.Thus, time and effort required to model the hydraulic turbine can be saved.Francis-99 is a complete test case available to investigate CFD on Francis turbines.Complete geometry, meshes and experimental measurements are available.The first workshop was held in Trondheim, Norway in December 2014.The workshop focused on steady state operation.The workshop has provided an open platform to the academic and industrial researchers to investigate a high head Francis turbine and explore the knowledge and their capabilities within hydropower.
Globally, the efficiency is well captured as well as the mean pressure and velocity profiles.This result is misleading as the torque and head may present strong deviation from the experimental results.This difficulty is the result of a poor estimation of the losses in the turbine.This limitation is attributed to the extensive use of wall function to resolve the near wall flow.The pressure amplitude due to rotor-stator interaction and vortex breakdown are other parameters difficult to capture accurately.
The next workshop will be held in December 2016.It will focus on transient operation, load variation and start-stop.

Fig. 1 .
Fig. 1.Locations of the mounted pressure sensors in the model Francis turbine during measurements at the BEP, HL, and PL; VL01vanless space, P42-blade pressure side, P71-blade trailing edge, S51-blade suction side, DT11 and DT21-draft tube; L1 and L2 correspond to the measurement sections for LDV.

Fig. 2 .
Fig. 2. Comparison of experimental and numerical values of hydraulic efficiency at the part load, BEP, and high load operating conditions, Stoessel [9].

Fig. 4 .
Fig. 4. Comparison of pressure amplitudes and frequencies in the vaneless space (VL01) during part load operating condition; the dimensionless frequency of 30 correspond to the blade passing frequency, Stefan [12].

Fig. 5 .
Fig.5.Comparison of axial (c ax ) and tangential (c u ) velocity profiles at the upper measurement section (L1) with the experimental velocity profiles at the part load, Mossinger[15]; c 0 is the average velocity at the measurement plane.

Fig. 6 .
Fig.6.Comparison of axial (c ax ) and tangential (c u ) velocity profiles at the lower measurement section (L2) with the experimental velocity profiles at the part load, Mossinger[15]; c 0 is the average velocity at the measurement plane.