Turbine-99 III: Proceedings of the third IAHR/ERCOFTAC workshop on draft tube flow

Full text

(1)2005:20. R E S E A R C H R E P O RT Proceedings of the third IAHR/ERCOFTAC Workshop on draft tube flows. Turbine-99 III. M.J. Cervantes T.F. Engström L.H. Gustavsson. Luleå University of Technology Department of Applied Physics and Mechanical Engineering Division of Fluid Mechanics 2005:20 • ISSN:1402-1528 • ISRN: LTU - FR -- 05⁄20 -- SE.

(2) To the memory of Prof. Rolf Karlsson.

(3) Local organizing Committee Dr. M.J. Cervantes, Chairman, Luleå University of Technology, Sweden Dr. T.F. Engström, Scientific-secretary, Luleå University of Technology, Sweden Prof. L.H. Gustavsson, Scientific-secretary, Luleå University of Technology, Sweden T. Gustafsson, Secretary, Luleå University of Technology, Sweden Prof. N. Dahlbäck, Vattenfall Utveckling AB, Sweden Prof. S Lundström, Luleå University of Technology, Sweden Dr. M. Page, Research Institute, Hydro Quebec, Canada. International Scientific Committee Dr. G. Constantinescu, University of Iowa, USA Prof. L. Davidsson, Chalmers University of Technology, Sweden Prof. L.H. Gustavsson, Scientific-secretary, Luleå University of Technology, Sweden.

(4) PREFACE The third IAHR/ERCOFTAC workshop on draft tube flows, Turbine-99 III, is based on the experience gained during the first two workshops and the development of the computational capacities. It is expected to be a step towards better understanding of draft tube flow simulation capabilities. Three cases based on the Turbine-99 benchmark were proposed to the participants: -. steady calculation, unsteady calculation, optimization of the draft tube performance.. More than 30 simulations have been performed by the participants with several turbulence models, near wall treatment, grids and boundary conditions. The complexity of the turbulence models ranges from zero equation model to large eddy simulations. The contribution of the different participants, the protocol of their simulations and a comprehensive comparison of experimental data with the simulations are included in the present document. The two days meeting, December 8-9, 2005, in Porjus will allow a detailed presentation of the participants work as well as discussions.. Dr. Michel Cervantes Project manger Hydro Power University.

(5)

(6) CONTENTS List of participants ................................................................................................................... 3 Workshop program .................................................................................................................. 5 Case description ...................................................................................................................... 7. PAPERS Buntić Ogor I., Dietze S. and Ruprecht A., Numerical simulation of the flow in Turbine-99 draft tube ................................................................................................................................ 23 Cervantes M.J. and Engström F., Eddy viscosity turbulence models and steady draft tube simulations ............................................................................................................................. 37 Tokyay T. and Constantinescu G., CFD simulations of flow in a hydraulic turbine draft tube using near wall RANS and LES models ................................................................................. 45 Galvan S, Page M, Guibault F and Reggio M, Numerical validation of different CFD k-ε turbulence models using FLUENT code ................................................................................57 Galvan S., Guibault F. and Reggio M., Optimization of the inlet velocity profile of Turbine-99 draft tube ................................................................................................................................67 Kurosawa S. and Nakamura K., Unsteady turbulent flow simulation in Turbine-99 draft tube .................................................................................................................................................73 Marjavaara B.D., Kamakoti R., Lundström T.S., Thakur S., Wright J. and W. Shyy, Steady and unsteady CFD calculations of the Turbine-99 draft tube using CFX-5 and Stream .......83 Nilsson H. and Page M., OpenFOAM simulation of the flow in the Hölleforsen draft tube model .....................................................................................................................................101 Page M., Giroux A.-M. and Nicolle, J., Steady and unsteady computations of Turbine-99 draft tube ........................................................................................................................................109 Raisee M., Alemi H. and Iacovides H.,Performance of linear and non-linear low-Re k-ε in prediction of developing turbulent flow through 90° curved ducts ......................................125. 1.

(7) COMPARISON BETWEEN EXPERIMENTAL DATA AND SIMULATION RESULTS Summary of the simulations ..................................................................................................137 Engineering quantities (evaluated by the participants) ..........................................................138 Engineering quantities (evaluated by the organizers) ............................................................139 Pressure along the upper centreline (experimental and computational results) ....................144 Pressure along the upper centreline (difference plot experiment and simulation) ................146 Pressure along the lower centreline (experimental and computational results) ....................148 Pressure along the lower centreline (difference plot experiment and simulation) ................150 Cross-section Ia - Radial pressure (experimental and computational results) .......................................152 Cross-section Ib -. Axial velocity (experimental and computational results) .........................................154 Axial velocity (difference plot experiment and simulation) ....................................156 Tangential velocity (experimental and computational results ...................................158 Tangential velocity (difference plot experiment and simulation) .............................160 Secondary flow ..........................................................................................................162. Cross-section II -. Axial velocity (experimental and computational results) .........................................164 Axial velocity (difference plot experiment and simulation) ....................................166 Horizontal velocity (experimental and computational results) .................................168 Horizontal velocity (difference plot experiment and simulation) .............................170 Secondary flow .........................................................................................................172 Clauser plot (lower wall, y = 250 mm) .....................................................................174 Clauser plot (lower wall, y = - 250 mm) ...................................................................176 Velocity profile tangential to the lower wall (y = 250 mm) ......................................178 Velocity profile tangential to the lower wall (y = - 250 mm) ...................................180. Cross section III -. Axial velocity (experimental and computational results) .........................................182 Axial velocity (difference plot experiment and simulation) ....................................184 Horizontal velocity (experimental and computational results) .................................186 Horizontal velocity (difference plot experiment and simulation) ............................188 Secondary flow ..........................................................................................................190. Cross section IVb - Secondary flow ..........................................................................................................192. 2.

(8) List of Participants Authors. Affiliation of the first author. Buntić I Dietze S. Ruprecht A.. Institute of Fluid Mechanics and Hydraulic Machinery University of Stuttgart Pfaffenwaldring 10 D-70550 Stuttgart, Germany Email: buntic@ihs.uni-stuttgart.de. 1. Home code. Cervantes M. Engström F.. Division of Fluid Mechanics Luleå University of Technology 971 87 Luleå, Sweden Michel.cervantes@ltu.se. 2. CFX. Tokyay T. Constantinescu G. .. Civil & Environmental Engineering IIHR-Hydroscience and Engineering The University of Iowa Stanley Hydraulics Laboratory Iowa City, IA 52242 ttokyay@engineering.uiowa.edu. 3. Fluent Home code. Galvan S. Page M. Guibault F. Reggio M.. École Polytechnique de Montréal Département de Génie Mécanique Section aérothermique E-mail sergio@galvan.polymtl.ca. 4,5. Fluent. Kurosawa S. Nakamura K.. Hydraulic Research Laboratory , Toshiba Corporation 20-1 , KANSEI-CHO , TSURUMI-KU,YOKOHAMA 230-0034,JAPAN E-mail : sadao.kurosawa@toshiba.co.jp TEL+81-45-510-6802 FAX+81-45-500-1430. 6. Home code. Marjavaara D. Kamakoti R. Lundström S. Thakur S., Shyy W.. Division of Fluid Mechanics Luleå University of Technology 971 87 Luleå, Sweden daniel.marjavaara@ltu.se. 7. CFX Loci-STREAM. Nilsson H. Page M.. Dept. of Thermo and Fluid Dynamics CHALMERS SE - 412 96 Gothenburg hani@tfd.chalmers.se. 8. OpenFOAM. Page M. Giroux A-M. Nicolle J.. Hydro-Québec, Institut de recherche 1800, boul. Lionel-Boulet, Varennes (Québec), Canada J3X 1S1 email: page.maryse@ireq.ca. 9. CFX Numeca. Raisee M. Alemi H. Iacovides H.. Department of Mechanical Engineering, Faculty of Engineering University of Tehran, Tehran, Iran. 10. Home code. 3. Paper CFD CODE.

(9) 4.

(10) Workshop Program Wednesday 07 December 2005 1500: Bus leaving Luleå Airport for Jokkmokk and Porjus 1900: Dinner at Hotel Porjus 1200-2100: Registration at Hotel Porjus. Thursday 08 December 2005 0845-0900: Welcome address 0900-0930: Lessons learned from Turbine-99 II, Professor H. Gustavsson 0930-1000: Coffee Session 1(Chairman: Dr. Constantinescu) 1000-1030: Numerical simulation of the flow in Turbine-99 draft tube, Buntić Ogor I., Dietze S. and Ruprecht A. 1030-1100: Eddy viscosity turbulence models and steady draft tube simulations, Cervantes M.J. and Engström F. 1100-1130: Numerical validation of different CFD k-ε turbulence models using Fluent code, Galvan S. et al. 1130-1200: OpenFOAM simulation of the flow in the Hölleforsen draft tube model, Nilsson H. and Page M. 1200-1330: Lunch Session 2 (Chairman: Prof. H. Gustavsson) 1330-1400: Performance of linear and non-linear low-Re k-ε in prediction of developing turbulent flow through 90° curved ducts, Raisee M. Alemi H. and Iacovides H. 1400-1430: Optimization of the inlet velocity profile of Turbine-99 draft tube, Galvan S. et al. 1430-1530: Large eddy simulation and unsteady RANS, Professor L. Davidson 1530-1600: Coffee Session 3 (Chairman: Prof. L. Davidsson) 1600-1630: Steady and unsteady CFD calculations of the Turbine-99 draft tube using CFX-5 and Stream, Marjavaara D. et al., 1630-1700: Steady and unsteady computations of Turbine-99 draft tube, Page M. et al. 1700-1730: CFD simulations of flow in a hydraulic turbine draft tube using near wall RANS and LES model, Tokyay T. and Constantinescu G. 1730-1800: Unsteady turbulent flow simulation in Turbine-99 draft tube, Kurosawa S. and Nakamura K. 1800-1900: Visit of the Porjus Hydropower Centre 2000: Dinner at Hotel Porjus. 5.

(11) Friday 09 December 2005 0830-0930: Comparisons of experimental data and simulations results; discussions (the organizers) 0930-1100: Group discussions (coffee) -. Numerics Boundary conditions Turbulence models Flow field characteristics Post processing Optimization Turbine-99 IV. 1100-1230: Presentation of the group discussions 1230-1400: Lunch 1400-1445: Plenary discussion (Prof. L. Davidsson, Prof. H. Gustavsson, Dr. Constantinescu): -. Conclusions Recommendations Further work. 1500: Bus leaving Porjus for Luleå Airport Alternative post-workshop activities. 6.

(12) Case Description The participants were asked to perform the simulation of a draft tube model; cf. figures 5-7 for detailed drawings and the download page at www.turbine-99.org for the CAD files in different format. The model is the same as used in the previous Turbine-99 workshops. Three different types of simulations were proposed.. Case 1: Steady calculation This case had the objective to assess the performance of different turbulence models. To achieve the goal, grid and boundary conditions were specified in detail. A detailed protocol had to be filled out. • • • • •. Grid 1 available in CGNS, Patran Neutral and Plot3D Boundary conditions specified in a section below Turbulence model – free Detailed protocol to be filled out for each simulation (available at the end of March 2005) Data and engineering quantities to be sent to the organization committee, see below. Case 2: Time resolved calculation This case had the objective to determine the differences between steady and unsteady calculations in terms of engineering quantities and flow phenomena. The ambition was also to assess the capability of different turbulence models for time resolved calculations. Grid and boundary conditions were all specified. In order to resolve the blade wakes, a refined grid was proposed. A detailed protocol had to be filled out as well. • • • • • •. Grid 2 available in CGNS, Patran Neutral and Plot3D at www.turbine-99.org Boundary conditions specified in a section below and in AngularResolved.zip available at www.turbine-99.org Turbulence model – free Detailed protocol to be filled out for each simulation (available at the end of March 2005) Data and engineering quantities to be sent to the organization committee, see below Maximum time step and minimum total simulation time specified below. Case 3: Draft tube optimization This case represents the final objective for the industry. Due to assumed economic considerations, the optimization had to fit inside the actual geometry. The objective function was based on the pressure recovery, i.e. the pressure difference between the inlet and the outlet of the draft tube. Examples of possible modifications were the cone and the elbow. • • • •. Geometry – available at www.turbine-99.org Grid – free Steady or unsteady – free Boundary conditions specified below and available at www.turbine-99.org. 7.

(13) • • • •. Turbulence model – free Objective function defined below Detailed protocol to be filled out for each simulation and found in appendix 1 Data and engineering quantities to be sent to the organization committee, see below. Boundary conditions Data were supplied for the operational mode T conducted at 60 % load of a Kaplan draft tube model, which is close to the best efficiency for the system and at the test head (H = 4.5 m). The mode was on-cam, i.e. the top-point (T) on the propeller curve (single runner blade angle curve). The exact settings are: • • • • •. runner speed: N=595 rpm (rotation per minute), flow rate: Q=0.522 m3/s , unit runner speed: DN/√H=140, where D is the runner diameter in meter, unit flow: Q/D2√H=1.00, water temperature: T=15 °C.. Correction for volume flow rate Rescale of the velocity data at the inlet to achieve the correct flow rate. All velocity components including rms data should be scaled with the same factor in order to maintain flow angles etc.. Case 1 Inlet boundary conditions For this case, velocity measurements were performed along a radial line at CS_Ia, cf. fig. 1 and 6 (CS stands for cross section). It was assumed that the velocity components and turbulent quantities were not functions of the azimuthal position, i.e. they were axisymmetric. The mean axial velocity (U), the mean tangential velocity (W) and the turbulent quantities ( u 2 , w 2 , uw ) are presented below in table 3. Since the experimental set of data is not complete, assumptions were made for the radial velocity (V) and the turbulent quantities ( v 2 , uv, vw ) to perform the simulation. The radial velocity (V) at the inlet was defined assuming the flow attached to the runner cone and the draft tube wall. The following relation was used:. V (r ) = U (r ) tan (θ ). θ = θ cone +. θ wall − θ cone Rwall − Rcone. (r − Rcone ). where Rcone ≤ r ≤ Rwall , θ cone = −12.8o and θ wall = +2.8o .. The unknown turbulent quantities ( v 2 , uv, vw ) were related to the measured turbulent quantities as follows: v 2 = w2 uv = vw = uw. 8.

(14) Wall boundary conditions Surface roughness: 10 µm on all walls. Other boundary conditions The other boundary conditions necessary to perform the simulation and not mentioned in the above lines, such as dissipation length scale, were left to the participants to define.. Case 2 Inlet boundary conditions For this case, angular resolved velocity measurements were performed along a radial line at CS_Ia, cf. fig. 1 and 6. The velocity data were divided into 120 bins (∆φ = 3°) for which averages and rms-values were calculated. These should be interpreted into a periodic time series of inlet boundary conditions. The mean axial velocity (U), the mean tangential velocity (W) and the turbulent quantities ( u 2 , w 2 , uw ) are available at www.turbine-99.org in a separate file named AngularResolved.zip. Since the experimental set of data was not complete, assumptions for the radial velocity (V) and the turbulent quantities ( v 2 , uv, vw ) were necessary to perform the simulation. The radial velocity (V) at the inlet was defined assuming the flow attached. The following relation was used: V (r ) = U (r ) tan (θ ) θ −θ θ = θ cone + wall cone (r − Rcone ) Rwall − Rcone where Rcone ≤ r ≤ Rwall , θ cone = −12.8o and θ wall = +2.8o . The unmeasured turbulent quantities ( v 2 , uv, vw ) were related to the measured turbulent quantities as follow: v 2 = w2 uv = vw = uw. A maximum time step (∆t) and a minimum simulation (t) time are given: •. time step:. •. simulation time: the outlet). ∆t ≤ 0.002 s (1/10 of a blade passage) (recommended: ∆t=0.0006 s) t ≥ 1.5 s (average time for a particle from the inlet to (recommended t ≥3 s). Wall boundary conditions Surface roughness: 10 µm on all walls.. 9.

(15) Other boundary conditions The other boundary conditions necessary to perform the simulation and not mentioned above, such as dissipation length scale, were left to the participants to define.. Case 3 The constraint for this case was the flow rate (Q=0.522 m3/s), i.e. any profile may be used. The participants may use the boundary condition of case 1 or 2 for the velocities. The objective function was the pressure recovery:. C prm =. 1 Aout. 1. ∫∫ P dA - A ∫∫ P dA. Aout. in Ain 2. 1 ⎛ Q ⎞ ρ⎜ ⎟ 2 ⎜⎝ Ain ⎟⎠. ,. where P is the pressure, ρ is the density (kg/m3), Q is the flow rate (m3/s), Ain is the area of cross section Ia (fig. 1 and 6) and Aout refers to cross section IVb (fig. 1 and 6). The pressure recovery factor indicates the degree of conversion of kinetic energy into static pressure where a higher value means higher efficiency for the draft tube. The exact value of the pressure recovery factor depends on the whole field solution and can be seen as an integral property of the solution.. Figure 1- Schematic of the draft tube with the measurement cross sections and visualization regions.. Requested data A certain amount of data was requested from the participants of the Turbine-99 III: 1. engineering quantities summarise the numerical results, 2. field data for comparison with experimental data at different cross sections (the cross sections are defined in appendix, figure 6).. 10.

(16) They are described below.. Engineering quantities The engineering quantities to be calculated were:. • Pressure recovery factor Cpr • Mean pressure recovery factor Cprm • Energy loss coefficient ζ Details of how to calculate each one of these quantities are described below. The pressure recovery factor is defined as: C pr =. Pout :wall − Pin:wall 1 ⎛ Q⎞ ρ⎜ ⎟ 2 ⎝ Ain ⎠. 2. ,. where Pout:wall is the outlet averaged static wall pressure at cross section IVb (fig. 1) Pin:wall is the inlet averaged static wall pressure at cross section Ia (fig. 1), ρ is the density (kg/m3), Q is the flow rate (m3/s) and Ain is the area at cross section Ia. Ain is an annular surface with the cone in the center. The inlet averaged static wall pressure, Pin:wall, is a direct result of the whole field solution for all variables (if the pressure is set at the outlet). The mean pressure recovery factor is defined as:. C prm =. 1 Aout. 1. ∫∫ P dA - A ∫∫ P dA in Ain 2. Aout. 1 ⎛ Q ⎞ ρ⎜ ⎟ 2 ⎜⎝ Ain ⎟⎠. ,. where Ain refers to cross section Ia (fig. 1) and Aout refers to cross section IVb (fig. 1). The energy loss coefficient is defined as: 2 2 ⎛ ⎞ ⎛ ⎞ ⎜ P + ρ U ⎟U ⋅ nˆ dA + ⎜ p + ρ U ⎟U ⋅ nˆ dA ∫∫ ∫∫ ⎜ ⎜ 2 ⎟⎠ 2 ⎟⎠ Ain ⎝ Aout ⎝ , ζ = 2 U ρ U ⋅ nˆ dA ∫∫ 2 Ain. where Ain refers to cross section Ia (fig. 1) and Aout refers to cross section IVb (fig. 1). Note that U ⋅ nˆ is mainly negative at Ain. This coefficient is directly coupled to the variation of the mechanical energy in the system. It is not used in experimental work since it requires knowledge about the variables over the whole inlet and outlet cross section. However, in CFD it is easily calculated. A lower value means lower loss for the draft tube. Since, the function of the draft tube is to convert kinetic energy into pressure, the loss factor is not an efficiency indicator (mechanical energy is the sum of the kinetic energy and the pressure).. 11.

(17) The requested data for comparison with experimental data are presented below in Table 1. Table 1- Requested data Position. Section Ia (cf. fig 1). Section Ib (cf. fig 1). Section II (cf. fig 1). Section III (cf. fig 1). Section IVb (cf. fig 1). Centerlines (cf. fig 3a and 3b). Data for Case 1 • x-velocity • y-velocity • z-velocity • pressure • x-velocity • y-velocity • z-velocity • pressure • x-velocity • y-velocity • z-velocity • pressure • x-velocity • y-velocity • z-velocity • pressure • x-velocity • y-velocity • z-velocity • pressure • upper centerline pressure • lower centerline pressure. Data for Case 2 • x-velocity (time averaged) • y-velocity (time averaged) • z-velocity (time averaged) • pressure (time averaged) • x-velocity (time averaged) • y-velocity (time averaged) • z-velocity (time averaged) • pressure (time averaged) • x-velocity (time averaged) • y-velocity (time averaged) • z-velocity (time averaged) • pressure (time averaged) • x-velocity (time averaged) • y-velocity (time averaged) • z-velocity (time averaged) • pressure (time averaged) • x-velocity (time averaged) • y-velocity(time averaged) • z-velocity (time averaged) • pressure (time averaged) • upper centerline pressure (time averaged) • lower centerline pressure, (all time-averaged). Coordinate system: The origin is in the centre of cross section Ia (fig. 1 and 2). The x-axis is pointing in the downstream direction towards the outlet of the draft tube. The z-axis is pointing upwards towards the runner and the y-axis is pointing to the right when watching the draft tube from the outlet (according to the right hand rule).. 12.

(18) Figure 2 - Coordinate system.. Upper wall centreline. Figure 3a - Upper wall centreline of the draft tube.. 13.

(19) Lower wall centreline. Figure 3b - Lower wall centreline of the draft tube.. Measurement details LDV-measurements for case 1 The measurements are done at z=0, (i.e. 127 mm below the runner hub centre) and at α=80° for Ia, cf. figure 4. The coordinate system is shown in fig. 2 and the exact geometry is defined by the CAD-file, which may be found on the download webpage at www.turbine99.org. The location of a single measurement is specified by the radius r, such as ƒ1 r 2 = x 2 + y 2 . The two measured mean velocity components are: • •. U: the axial velocity (in the negative z-direction; cf. Figure 2, side view), W: the tangential velocity (in the negative α-direction; cf. Figure 4, top view). A positive tangential velocity co-rotates with the runner.. Corresponding to these mean values, RMS-valuesƒ2 and the cross correlation between the two components have been determined. u: the RMS(root mean square) valueƒ3 of U w: the RMS valueƒ3 of W uw: the cross correlationƒ3 of U and W. The data are presented below in table 3.. 14.

(20) ƒ1. The mean velocity (X) is defined as. X =. ∑ xτ ∑τ. i i. , where xi is an individual sample and τi is the corresponding weight constant.. i. For these data τi=1. ƒ2. The RMS-values (x’) is defined as. x' =. ∑x τ ∑τ. 2 i i. − X2. , see above.. i. ƒ3. Note that these values include the periodic behaviour of the components and do not only relate to the turbulence intensity.. Figure 4 - Side and top view of the draft tube inlet cone (Read W and not V!). Table 3 - Velocity data for case 1 at section Ia(1). Radius [m]. U [m/s]. W [m/s]. u [m/s]. w [m/s]. 0.0981* 0 6.11 0 0 0.1014 2.44 0.37 0.63 0.49 0.1041 2.55 0.58 0.70 0.51 0.1068 2.72 0.72 0.72 0.61 0.1094 2.86 0.69 0.68 0.48 0.1161 3.07 0.65 0.51 0.39 0.1227 3.06 0.61 0.50 0.36 0.1294 3.14 0.60 0.32 0.37 0.1427 3.27 0.65 0.30 0.37 0.1560 3.38 0.79 0.28 0.38 0.1693 3.41 0.91 0.30 0.37 0.1826 3.57 1.07 0.38 0.47 0.1959 3.69 1.31 0.42 0.70 0.2092 3.71 1.28 0.34 0.52 0.2168 3.65 1.29 0.47 0.49 0.2221 3.55 1.43 0.54 0.63 0.2261 3.55 1.57 0.51 0.80 0.2301 3.54 1.59 0.57 0.83 0.2314 3.48 1.56 0.73 0.79 0.2328 3.40 1.52 0.79 0.73 0.2341 3.16 1.45 0.98 0.65 0.2365* 0 0 0 0 *: wall **: extrapolated values from previous measurements. 15. uw** [m2/s2] 0 -0.011 -0.020 -0.030 -0.038 -0.060 -0.050 -0.030 -0.020 -0.030 -0.040 -0.048 -0.121 -0.047 -0.021 -0.002 0.015 0.033 0.038 0.036 0.031 0.

(21) LDV-measurements for case 2 The LDV measurements for case 2 are available on a separate file named: AngularResolved.zip. The file contains a description of the angular resolved measurements and several files for the data, cf. www.turbine-99.org.. Pressure measurements Together with the geometry downstream the draft tube (cf. figure 3 in appendix) the wall pressure at the outlet is supplied to help the setting of a suitable outlet condition in the calculations. The wall pressure measurements are located at IVa and IVb (cf. figure 1 and 6). The approximate locations of the pressure taps along the circumference of the cross sections are shown in Figure 5. The exact coordinates (in mm) in the data tables are given in the local coordinate system, [yl zl], defined in the same figure. ⎧ y1 = y ⎨ ⎩ z1 = − z + 421.2. The mean value and the standard deviation of the measurements are given:. ƒ. . Pw: mean value of the wall pressure pw: standard deviationƒ of the wall pressure. Note that this value includes the periodic behaviour of the component and the cross correlation with the reference point.. Figure 5 - Upstream view of the pressure tap locations at the outlet of the test draft tube.. 16.

(22) Table 4 - Pressure data at Section IV Section Position yl zl Pw pw (Fig. 5) [mm] [mm] [mBar] [mBar] b1 -250 1063 -1.604 b 1.63 b2 0 1063 0.000 b 1.56 b3 250 1063 -0.039 b 1.62 l1 500 925 -0.711 b 1.92 l2 500 660 -0.693 b 1.85 l3 500 395 -1.210 b l4 500 130 -0.516 b 1.74 r1 -500 130 -2.309 b 1.77 r2 -500 395 -1.807 b 1.66 r3 -500 660 -1.208 b 1.68 r4 -500 925 -1.708 b 1.71. a a a a a a a a. l1 l2 l3 l4 r1 r2 r3 r4. 500 500 500 500 -500 -500 -500 -500. 925 660 395 130 130 395 660 925. -0.330 0.085 -0.800 -2.489 -2.986 -2.174 -1.829 -1.726. 17. 1.76 1.78 1.68 1.62 1.89 1.83 1.66 1.60.

(23) Drawings of the geometry. Figure 5 - Drawing of the draft tube.. 18.

(24) Figure 6 - Location of the different cross sections.. 19.

(25) 31. 2600. 1063 1050. Equalization. 104. 500. 1100. Figure 7 - The downstream tank at the draft tube outlet. The draft tube ends in a large cylindrical tank, having a diameter of 2600 mm, and the outlet centre of the draft tube is located 282 mm above the cylinder centre.. References Daugherty, R.L, Franzini, J.B. and Finnemore, E.J., (1989), Fluid Mechanics with engineering applications, McGraw-Hill Book Company, Singapore. Seeno, Y., Kawaguchi, N. and Tetsuzou, T., (1978), Swirl flow in conical diffusers, Bulletin of the JSME, Vol. 21., No. 151., January, pp. 112-119. Staubli, T. and Deniz, S., (1994), Analysing draft tube kinetics, Int. Water Power & Dam Construction, December, pp. 38-42.. 20.

(26) PAPERS. 21.

(27) 22.

(28) Turbine99 III Proceedings of the third IAHR/ERCOFTAC workshop on draft tube flow 89 December 2005, Porjus, Sweden Paper No. xxx1. NUMERICAL SIMULATION OF THE FLOW IN TURBINE99 DRAFT TUBE Buntić Ogor I., Dietze S., Ruprecht A. Institute of Fluid Mechanics and Hydraulic Machinery University of Stuttgart Pfaffenwaldring 10, D-70550 Stuttgart, Germany buntic@ihs.uni-stuttgart.de. In this paper numerical simulation of the flow in the turbine-99 draft tube using the in-house finite element code FENFLOSS (Finite Element based Numerical FLOw Simulation System) is presented. Besides modeling with standard k-ε model, main focus is on the simulations with the extended k-ε model of Chen-Kim and Very Large Eddy Simulation (VLES) known as promising tool for prediction of unsteady phenomena. Both turbulence modeling approaches are explained in details as well as the applied adaptive filtering technique which can distinguish between numerically resolved and unresolved parts of the flow. Numerical results are analyzed and presented.. problems. Widely spread in industry are Reynoldsaveraged Navier-Stokes (RANS) based turbulence models which are usually not capable to reproduce complex 3D turbulent flows, especially unsteady flows. The highest accuracy for resolving complete turbulence is offered by a Direct Numerical Simulation (DNS). Unfortunately its industrial application is not possible in the foreseeable future. Lately Large Eddy Simulation (LES) starts to be a mature technique despite its high computational cost. Very Large Eddy Simulation (VLES), or also known as Detached Eddy Simulation (DES), starts to expand as a promising compromise for simulation of industrial flow problems with reasonable computational time and costs.. INTRODUCTION. SIMULATION METHOD. ABSTRACT. The increasing consumption and demand on energy again takes into account hydropower as important energy source. Therefore the increase of old hydropower plants efficiency through refurbishment and upgrading is a challenging task. The quality and the design of the draft tube are very important. The power output of a hydraulic turbine is specially affected by the performance of its draft tube whose role is to convert the kinetic energy behind the runner into static pressure. There is a potential of improving the pressure recovery in the draft tube by modifying and optimizing its geometry. Flow in a draft tube is characterized as very intricate turbulent flow followed with appearance of different flow phenomena, e.g. unsteadiness, flow separation, swirling flow etc. Thus its simulation is complicated and time consuming requiring high computational power. Additionally, adequate turbulence modeling is needed which is able to predict such flows satisfactorily. An example of complexity of numerical simulations which are connected to the real industrial applications is Turbine-99 (Engström [1]). Nowadays Computational Fluid Dynamics (CFD) is a powerful tool for simulation and analyzing complex flows such as flows in the elbow draft tubes. It reduces the time and the costs which are necessary for design and optimization of the turbines i.e. draft tubes. The most important part of the modeling is turbulence modeling. It is still one of the fundamental CFD. Governing equations and numerical method The governing equations describing incompressible, viscous and time dependant flow are the NavierStokes equations. They express the conservation of mass and momentum. In the RANS approach, these equations are time or ensemble averaged leading to the well known RANS equations: ∂τ ij ∂U i ∂U i ∂P +Uj =− + ν∇ 2U i − (1) ∂t ∂x j ∂x i ∂x j. ∂U i =0 (2) ∂x i In RANS τij expresses the Reynolds stress tensor which is unknown and has to be modelled. The task of turbulence modelling is the formulation and determination of suitable relations for Reynolds stresses. Details of the new VLES approach are described in next section. The calculations are performed using the program FENFLOSS which is developed at the Institute of Fluid Mechanics and Hydraulic Machinery, University of Stuttgart. It is based on the Finite Element Method. For spatial domain discretisation 8-node hexahedral elements are used. Time discretisation involves a three-level nd fully implicit finite difference approximation of 2 order. For the velocity components and the turbulence quantities a trilinear approximation is. 23.

(29) applied. The pressure is assumed to be constant within element. For advection dominated flow a nd Petrov-Galerkin formulation of 2 order with skewed upwind orientated weighting function is used. For the solution of the momentum and continuity equations a segregated algorithm is used. It means that each momentum equation is handled independently. They are linearised and the linear equation system is solved with the conjugated gradient method BICGSTAB2 with an incomplete LU decomposition (ILU) for preconditioning. The pressure is treated with the modified Uzawa pressure correction scheme (Ruprecht [2]). The pressure correction is performed in a local iteration loop without reassembling the system matrices until the continuity error is reduced to a given order. After solving the momentum and continuity equations, the turbulence quantities are calculated and a new turbulence viscosity is gained. The k and ε equations are also linearised and solved with BICGSTAB2 algorithm with ILU preconditioning. The whole procedure is carried out in a global iteration until convergence is obtained. For unsteady simulation the global iteration has to be performed for each time step. The code is parallelised and computational domain is decomposed using overlapping grids. In that case the linear solver BICGSTAB2 has a parallel performance and the data exchange between the domains is organised on the level of the matrixvector multiplication. The preconditioning is then local on each domain. The data exchange uses MPI (Message Passing Interface) on computers with distributed memory. On the shared memory computers the code applies OpenMP. For more details on the numerical procedure and parallelisation the reader is referred to Maihöfer [3] and Maihöfer et al.[4].. Figure 1 - Degree of turbulence modelling and computational effort for the different approaches.. LES starts to get more practical importance. With LES all anisotropic turbulent structures are resolved in the computation and only the smallest isotropic scales are modelled (schematically shown in Fig. 2). The models used for LES are simple compared to those used for RANS because they only have to describe the influence of the isotropic scales on the resolved anisotropic scales. With increasing Reynolds number the small anisotropic scales strongly decrease becoming isotropic and therefore not resolvable. If there is a gap in the turbulence spectrum between the unsteady mean flow and the turbulent flow, ”classical” RANS i.e. URANS models can be applied, as they are developed for modelling the whole range of turbulent scales (Fig. 2). They show excessive viscous behaviour and very often damp down unsteady motion quite early. It also means that they are not suitable for prediction and analysis of many unsteady vortex phenomena.. Turbulence modeling RANS equations are established as a standard tool for industrial simulations, although it means that the complete turbulence behaviour has to be enclosed within appropriate turbulence model which takes into account all turbulence scales (from the largest eddies to the Kolmogorov scale). DNS is capable to resolve all turbulence scales but it requires a very fine grid resolution. Hence carrying out 3D simulations of the flow with high Reynolds number is very time consuming even for high performance computers (Fig. 1).. Figure 2 - Modelling approaches for RANS and LES.. Figure 3 - Modelling approach in VLES.. 24. Contrary, if there is no spectral gap and even one part of the turbulence can be numerically resolved, VLES can be used. It is very similar to the LES, with.

(30) the difference that a smaller part of the turbulence spectrum is resolved and the influence of a larger part of the spectrum has to be expressed with the model (Fig. 3). It requires the use of an appropriate filtering technique which distinguishes between resolved and modelled part of the turbulence spectrum. Because of its adaptive characteristic it can be applied for the whole range of turbulence modelling approaches from the RANS to the DNS (Fig. 4).. with following coefficients: σk = 0.75, σε = 1.15, c1ε = 1.15, c2ε = 1.15 and c3ε = 0.25. Additionally, these extended k-ε equations need to be filtered. The applied filtering technique is similar to Willems [6]. The smallest resolved length scale Δ used in filter is according to Magnato et al. [7] dependant on the local grid size or the computational time step and local velocity. According to the Kolmogorov theory it can be assumed that the dissipation rate is equal for all scaled. This leads to ε = εˆ (5) This is not acceptable for turbulent kinetic energy. It is filtered according to   Δ  (6) kˆ = k ⋅ 1 − f    L   As a suitable filter 0 for  2/3  f =  Δ 1− for   L  . Figure 4 - Adjustment for adaptive model.. Lately several hybrid methods have been proposed in the literature. All of them are based on the same idea to represent a link between RANS and LES. They try to keep computational efficiency of RANS and the potential of LES to resolve large turbulent structures, even on coarser grids and with high Reynolds number. Various VLES methods slightly differ in filtering techniques, applied model and interpretation of the resolved motion, but broadly speaking they all have a tendency to solve relevant part of the flow and model the rest (Fig. 5)..  u ⋅ Δt Δ = α ⋅ max   hmax. c1ε. ε k. Pk − c 2ε.  ∂ε   ∂x  j. ε2. P  + c 3ε  k  ⋅ Pk k  k2  43 144 4 additional. term. wit h.  ΔV hmax = 3  ΔV. for 2D for 3D. (8). contains model constant α in a range from 1 to 5. Then the Kolmogorov scale L for the whole spectrum is given as k 3/2 L= . (9). ε. Modelled length scales and turbulent viscosity are kˆ 3 / 2 Lˆ = (10) εˆ kˆ 2 νˆt = c µ ⋅ (11) εˆ with cµ = 0.09. The filtering procedure leads to the final equations ∂k ∂k ∂  νˆ  ∂k  ˆ (12) + Uj = ν + t   + Pk − ε ∂t σ k  ∂x j  ∂x j ∂x j  νˆt  ∂ε  ∂ε ∂ε ∂    +Uj = ν +  + ∂t ∂x j ∂x j  σ ε  ∂x j  (13)  Pˆk . ε ˆ ε2 c1ε Pk − c 2ε + c 3ε   ⋅ Pˆk k k  k  with the production term ) ∂U j  ∂U i )  ∂U  Pk = ν t  i + . (14)  ∂x j ∂x i  ∂x j . For more details of the model and its characteristics the reader is referred to Ruprecht [8]. Furthermore, VLES, of course depending on the grid size, is not valid all the way to the wall. Modelling procedure close to the wall is based on standard wall functions, which are most commonly used in industrial practice. Standard wall functions are as well used in case of extended k-ε model of Chen and Kim.. The basis of the adaptive model is the extended k-ε model of Chen and Kim (Chen et al. [5]). It is chosen due to its simplicity and capacity to better handle unsteady flows. Its transport equations for k and ε are given as ∂k ∂k ∂  ν  ∂k  (3) + Uj = ν + t   + Pk − ε ∂x j ∂x j  ∂t σ k  ∂x j .  νt ν + σε . (7). Δ< L. is applied where. Figure 5 - Distinguishing of turbulence spectrum by VLES.. ∂ε ∂ε ∂ +Uj = ∂t ∂x j ∂x j. Δ≥ L.  +  (4). 25.

(31) COMPUTATIONAL DETAILS . 4,5. 4. Geometry, boundary conditions and simulation setup Turbine-99 draft tube is an elbow sharp-heel draft tube developed for use in Kaplan turbines. Draft tube geometry given by workshop organisers is fixed, as well as defined cross sections (for more details see [9]). Original computational domain starts at cross section CS-Ia and ends at cross section CS-IVb. It also includes a part of the runner hub at the inlet. Computational domain used for simulations presented in this paper is slightly modified. The 2.1m straight extension was added due to the outlet boundary condition. Used computational mesh was also provided by workshop organisers. Computations were performed at the coarsest available mesh (1 million elements, y+ ≈ 50, see Döbbener [10]), due to the interest to investigate behaviour and characteristics of VLES at coarser grids which are computationally “cheaper”. Due to geometry extension at the outlet, real mesh consisted of ca. 1.2 million elements. Inlet boundary conditions are as well predefined by workshop organisers (for details see [9]). They are represented with velocity profiles which are available by means of extensive LDV laboratory measurements. They provide the mean axial (U) and tangential (W) velocity components (Fig. 6 and Fig.. Mean axial velocity [m/s]. 3,5. 3. 2,5. 2. 1,5. 1. 0,5. 0 0,09. 0,11. 0,13. 0,15. 0,17. 0,19. 0,21. 0,23. 0,25. Radius r [m]. . Figure 6 – Profile of measured mean axial velocity. 6. Mean tangential velocity [m/s]. 5. 4. 3. 2. 1. 0 0,09. 0,11. 0,13. 0,15. 0,17. 0,19. 0,21. 0,23. 0,25. Radius r [m]. . Figure 7 - Profile of measured mean tangential velocity.. 7) and three turbulent quantities ( u ' , w ' and u ' w ' ) (Fig. 9) in radial direction at CS-Ia. Radial velocity component (V) is not recorded and therefore no measured data are available. The same is valid for. 0,2. 0,1. 0 0,09. turbulent quantities v ' , u ' v ' and v ' w ' . Therefore assumption for the radial velocity is made using following (Fig. 8): V (r ) = U (r ) ⋅ tan( Θ). 0,11. 0,13. 0,15. 0,17. 0,19. 0,21. 0,23. 0,25. Radial velocity [m/s]. -0,1. -0,2. -0,3. -0,4. Θ − Θ cone (r − R cone ) Θ = Θ cone + wall Rwall − R cone with Rcone≤ r ≤ Rwall, Θcone = -12.8° and Θwall = +2.8°. The unknown turbulent quantities are assumed as follows:. -0,5. -0,6. -0,7 Radius r [m]. . Figure 8 - Profile of assumed mean radial velocity.. v' = w' u ' v ' = v ' w ' = u' w ' As turbulence kinetic energy is needed for the calculations, it is estimated using (Fig. 10) 1 k =  u ' 2 + v ' 2 + w ' 2  .  2 Estimation of dissipation rate (Fig. 11) is performed using k 1.5 ε = cµ l with constant cµ = 0.09 and l = 0.1.. 1,20. 1,00. Turbulent quantities . 0,80. u' w' u'w' u' scale w' scale u'w' scale. 0,60. 0,40. 0,20. 0,00. -0,20 0,09. 0,11. 0,13. 0,15. 0,17. 0,19. 0,21. 0,23. 0,25. Radius r [m]. . Figure 9 - Profile of measure u' , w ' and u' w ' .. . 26. . .

(32) boundary conditions. Steady cases were performed for calculations with k-ε model and extended k-ε model of Chen and Kim. In both cases previously mentioned wall function was used. Unsteady calculations were performed for further investigation of VLES. Used time steps were 0.1s and 0.01s. Important was to investigate the influence of time step on the filtering approach used in VLES.. 1,2. Turbulence kinetic energy [m2/s2]. 1. 0,8. 0,6. 0,4. 0,25. Steady calculation of k-ε model showed good and quick convergence. Velocity distributions in symmetry and cut planes are shown in Fig. 12 and Fig. 13, respectively. Streamlines of the flow behind the runner and at the outlet of the draft tube are presented in Fig.14 and Fig.15.. 0,25. Figure 12 - Velocity distribution in symmetry plane, steady k-ε. 0,2. 0 0,09. 0,11. 0,13. 0,15. 0,17. 0,19. 0,21. 0,23. Radius r [m]. Figure 10 - Turbulence kinetic energy profile. 1. 0,9. 0,8. Dissipation rate [m2/s3]. 0,7. 0,6. 0,5. 0,4. 0,3. 0,2. 0,1. 0 0,09. 0,11. 0,13. 0,15. 0,17. 0,19. 0,21. 0,23. Radius r [m]. Figure 11 - Dissipation rate profile.. Correction of values was not performed due to the fact of significant discrepancy of ca. 6.5% between stated flow rate and the one integrated from measured velocity. Consequently, stated Reynolds number is 1.7 million and Reynolds number in simulation is 1.6 million. Difference can be seen in Fig 6. to Fig. 11. (lines with filled square are original not scaled data and lines with not filled square are scaled data). Slight modification of the profile at the runner point in axial direction was necessary due to implementation of rotating wall boundary condition and stability of convergence. Boundary condition set at the runner hub is rotating wall with runner speed of N = 595 rpm. For all other cone and side walls no-slip wall boundary condition in conjunction with standard wall function was set. At the outlet section constant pressure was imposed. Simulated operational point was operational mode T (defined by workshop organisers) corresponding to 60% load which is close to the best efficiency of the given system. Defined flow rate in mode T is 3 Q = 0.522 m /s.. Figure 13 - Velocity distribution in planes Ib, II, III and IVb, steady k-ε. RESULTS AND DISCUSSION One of the main aims of the investigation is assessment of FENFLOSS code and implemented turbulence models to predict the flow phenomena and pressure recovery in Turbine-99 draft tube. Valuable knowledge is also general understanding of the flow which appears in draft tubes. Calculations of Turbine-99 were performed as steady and unsteady calculations with already mentioned. Figure 14 – Streamlines behind the runner, steady k-ε. 27.

(33) Figure 18 - Streamlines after the runner, steady Chen and Kim model. Figure 15 – Streamlines at cross section IVb, steady k-ε. Due to the additional terms, extended k-ε model of Chen and Kim is expected to show in study case slightly different result i.e less damping of swirl behind the runner. As it can be seen in Fig. 16 to Fig. 18 swirling flow behind the runner is a bit brighter and the region of higher velocity on the right side (in the flow direction) is a bit longer. It can be also seen that the structure of the vortex developing behind runner is different (Fig. 14 and Fig. 18). Velocity distributions in cross sections II, III and IVb do not show significant difference for investigated models. In both cases vortices are suppressed after the elbow. At the place where flow enters the horizontal part of the draft tube, Chen and Kim model shows still existing small vortex on the right side. On the other hand, both cases show the appearance of the recirculation zone on the right side (in the flow direction) at the cross section IVb (Fig. 15 and Fig. 19).. Figure 19 - Streamlines at cross section IVb, steady Chen and Kim model. Further investigations were concentrated on unsteady calculations using VLES. In presented VLES approach unsteady calculations were necessary due to the filtering technique. Two time steps were used 0.1s and 0.01s to investigate the influence of time step on implemented filtering. According to Helmrich [11] model constant α was set to 2. Due to the previously presented definition of the filtering approach it is expected that for coarse mesh filter function becomes 0 and thus whole turbulence is expected to be modelled with Chen and Kim turbulence model. On the other hand for fine mesh filtering function approaches 1 resulting in complete resolving of the turbulence. Used mesh is relatively fine, times steps are coarse and relatively fine. Therefore the question is which part influences the filtering - mesh size or u ⋅Δt .. Figure 16 - Velocity distribution in symmetry plane, steady Chen and Kim model. It was expected that time step 0.1s would lead to the more modeling of the whole turbulence spectrum and that result would develop toward the results calculated with Chen and Kim model. Unfortunately, for this calculation convergence problems occurred. As initialization for VLES with time step 0.01s a short unsteady calculation with Chen and Kim with the same time step was used. It showed the convergence to the same steady state which is previously presented. VLES calculation using time step 0.01s was time averaged and presented in Fig. 20 to Fig. 23.. Figure 17 - Velocity distribution in planes Ib, II, III and IVb, steady Chen and Kim model. 28.

(34) well, strong recirculation zone in cross section IVb is observed. Fig. 24. presents the zone of low pressure (light colored left side in flow direction) where unsteadiness seems to appear and continues to spread through draft tube always on the right side.. Figure 20 - Velocity distribution in symmetry plane, VLES. Figure 24 - Vortex streamlines. Calculated pressure recovery factor (based on the wall pressure) and mean pressure recovery factor (based on the entire inlet cross section) for investigate cases are presented in Table 1.. Figure 21 - Velocity distribution in planes Ib, II, III and IVb, VLES. Table 1. Pressure recovery factors Case Cpr 1.047 k-ε steady Chen & Kim steady 1.036 VLES 0.1s VLES 0.01s 1.056. Cprm 0.790 0.731 0.714. As the simulation of the flow in draft tube is time consuming, most of the calculations are performed on CRAY Opteron and NEC Xeon clusters at High Performance Computing Center HLRS, University of Stuttgart. Figure 22 - Streamlines after the runner, VLES. CONCLUSIONS In this paper simulations with standard k-ε model, extended k-ε model of Chen and Kim and VLES are presented. Predefined mesh and boundary conditions were used with the exception of the outlet region. All steady state calculations converged well. Slight difference of results between k-ε and Chen and Kim model can be noticed, especially in vortex behind the runner. Contrary to both k-ε models, VLES, which is based on Chen and Kim k-ε model, shows much more unsteady phenomena throughout the draft tube. It means that one part of the flow is really resolved and other part modelled.. Figure 23 - Streamlines at cross section IVb, VLES. As it can be seen VLES shows different flow picture compared to the two previous cases. The flow shows a much more unsteady and turbulent characteristic. Appeared vortex behind the runner does not show in cross section Ib circular shape. Vortex stretches throughout whole draft tube after the elbow, especially on the right side (in the flow direction). As. ACKNOWLEDGEMENT The authors would like to acknowledge the help of Mr. Volker Brost in data transformation and also thank other institute’s colleagues for their valuable comments and suggestions.. 29.

(35) [10] DöbbenerG., Grotjans H., 2005, Mesh Generation for Turbine 99 Workshop, ANSYS Germany [11] Helmrich T., Simulation Instationärer Wirbelstrukturen in Hydraulischen Maschinen, Ph.D. thesis, University of Stuttgart. NOMENCLATURE f hmax k L Pk u Ui U P τij α ∆ ∆t ε ν ν t ∆V N Q U V W Cpr Cprm ^ ~ i. filter function local grid size, m 2 2 turbulent kinetic energy, m /s Kolmogorov length scale, m production term local velocity, m/s filtered velocity, m/s averaged velocity, m/s averaged pressure, Pa Reynolds stresses, Pa model constant resolved length scale, m time step, s 2 3 dissipation rate, m /s 2 kinematic viscosity, m /s 2 turbulent viscosity, m /s 2 3 size of the local element, m or m runner speed, rpm 3 flow rate, m /s axial velocity, m/s radial velocity, m/s tangential velocity, m/s pressure recovery factor mean pressure recovery factor modelled resolved covariant indices, i = 1,2,3. REFERENCES [1] Engström T.F., Gustavsson H., Karlsson R. I., 2002, Turbine-99 – Workshop 2 on Draft Tube Flow, Proceedings of the XXIst IAHR Symposium on Hydraulic Machinery and Systems, Laussane, Switzerland [2] Ruprecht A., 1989, Finite Elemente zur Berechnung dreidimensionaler turbulenter Strömungen in komplexen Geometrien, Ph.D. thesis, University of Stuttgart [3] Maihöfer M., 2002, Effiziente Verfahren zur Berechnung dreidimensionaler Strömungen mit nichtpassenden Gittern, Ph.D. thesis, University of Stuttgart [4] Maihöfer M., Ruprecht A., 2003, A Local Grid Refinement Algorithm on Modern High-Performance Computers, Proceedings of Parallel CFD 2003, Elsavier, Amsterdam [5] Chen Y.S., Kim S.W., 1987, Computation of turbulent flows using an extended k-ε .turbulence closure model, NASA CR179204 [6] Willems W., 1997, Numerische Simulation turbulenter Scherströmungen mit einem Zwei-Skalen Turbulenzmodell, Ph.D. thesis, Shaker Verlag, Aachen [7] Magnato F., Gabi M., 2000, A new adaptive turbulence model for unsteady flow fields in rotating machinery, Proceedings of the 8th International Symposium on Transport Phenomena and Dynamics of Rotating Machinery (ISROMAC 8) [8] Ruprecht A., 2005 Numerische Strömungssimulation am Beispiel hydraulischer Strömungsmaschinen, Habilitation thesis, University of Stuttgart [9] Turbine-99 Workshop III - Description. 30.

(36) 31.

(37) 32.

(38) 33.

(39) 34.

(40) 35.

(41) 36.

(42) Turbine-99 III Proceedings of the third IAHR/ERCOFTAC workshop on draft tube flow 8-9 December 2005, Porjus, Sweden Paper No. 2. EDDY VISCOSITY TURBULENCE MODELS AND STEADY DRAFT TUBE SIMULATIONS Cervantes M.J. and Engström F. Division of Fluids Mechanics Luleå University of Technology 971 87 Luleå, Sweden michel.cervantes@ltu.se. equation (ZE), standard k-ε and shear stress transport (SST) models. The use of the zero equation model may be surprising at first. It is motivated by the small sensitivity of the engineering quantities to the dissipation length scale with the standard k-ε model, cf. Cervantes and Engström [1]. Detailed analysis of the resulting flows was performed to compare the different turbulence models and possible improvements have been identified.. ABSTRACT Computations of the Turbine-99 benchmark have been performed for two dimensional steady inlet boundary conditions. Three different turbulence models were used: zero equation model, k-ε and shear stress model (SST). The results from the engineering quantities indicate small differences on the mean pressure recovery and the loss factor, while larger differences appear for the wall pressure recovery.. NUMERIC. INTRODUCTION. The commercial code ANSYS CFX 10.0 was used to perform the simulations. The code used the finite volume method and has a coupled unstructured solver.. Draft tube flow simulations are challenging for the numerical community due to the different flow phenomena appearing simultaneously such e.g. turbulence, separation, swirl and unsteadiness. Turbine-99 aims to determine the state of the art within the area by proposing a test case, which consists of a model draft tube of the Hölleforsen hydropower station with detailed pressure and velocity measurements at different cross-sections.. Turbulence models Three different turbulence models have been used: zero equation model (ZE), standard k-ε and shear stress transport (SST). The eddy viscosity assumption is the base of the three different models used in the work. The eddy viscosity assumption assumes that. The first workshop was held in Porjus in 1999. A large scatter in the pressure recovery (±50%) was obtained. The scatter was attributed to the large number of free parameters such as grid topology and inlet radial velocity. The second workshop was held in Älvkarleby in 2001. There, grid, turbulence model and boundary conditions were set by the organizers. The scatter encountered during the first workshop disappeared. However, some discrepancies still appeared between the participants results. The reason was attributed to the implementation of the boundary conditions and the determination of the engineering quantities, i.e. pre and post processing. After discussion, the elaboration of a protocol was proposed to make simulations details such as max residual or scheme order available to be able to explain discrepancies between different results.. the Reynolds stress. ( − ρ uiu j ) is. proportional to the. strain rate tensor ( Sij ) according to 2 kδ ij (1), 3 where k represents the turbulent kinetic energy. − ρ ui u j = µT Sij −. The zero equation model models the turbulence viscosity µT as the product of a turbulent velocity scale UT and a turbulence length scale lT such as. µT = ρ f µ U T lT. (2),. 1 3 D. where l T = V / 7 and VD is the fluid volume domain, in the present case lT = 0.19 m . UT is the maximum velocity in the fluid domain. This model was proposed by Prandtl and Kolmogorov and is based on an analogy with gas viscosity. Such a formulation may be correct if the flow had one length scale and time scale, which is not the case for internal turbulent flows. Nevertheless, the results may be of interest for comparison with more complex turbulence models.. The third IAHR/ERCOFTAC Turbine-99 workshop is a natural continuation of the two first workshops. For this workshop, three cases are proposed: steady calculation, unsteady calculation and optimization. The present work presents steady flow simulations performed with 3 different turbulence models: zero. 37.

(43) The standard k-ε where a transport energy k and one ε are introduced. turbulent velocity turbulent viscosity:. The reason for the lack of convergence with second order accurate scheme for the momentum and the continuity equation with the SST model is unclear. It may be attributed to the unsteady behavior of the flow in the draft tube.. model is a two-equation model, equation for the turbulent kinetic for the turbulent eddy dissipation These two variables yield the and length scales to form the. Table 1 – Scheme order used. k2. µT = Cµ ρ. ε. ZE. (3),. Continuity Momentum Ke Ed Tef. SBF 1 SBF 1 SBF: specified blend factor HR: high resolution. where Cµ is a constant. The transport equation for the turbulent kinetic energy fails to describe the flow near the wall. The use of wall functions, damped viscosity or two-layer models is necessary to overcome the failure. In the present work, scalable wall function was used [2].. k. ω. SST SBF 0.5 SBF 0.5 HR HR. RESULTS. The shear stress transport (k-ω) based model is a two-equation model similar to the k-ε model. The transport equation for the turbulent dissipation is replaced by an equation for the turbulent frequency ω. Turbulent kinetic energy and frequency are related through the turbulent viscosity such as. µT = ρ. k-ε SBF 1 SBF 1 HR HR -. Convergence The residuals obtained for the simulations are presented in table 2. The minimum, area average and maximum values of y+ on the wall and the runner hub are presented in table 3. Higher values are obtained on the runner hub, especially for the SST model. The grid is not enough fine to allow y+ around 1. Nevertheless, the values obtained are satisfying compared to those obtained during the second workshop; see Cervantes and Engström [3].. (4).. The k-ω formulation is advantageous for near wall treatment compared to the standard k-ε. The SST model is a development of the Wilcox and baseline (BSL) models. It is known to give accurate predictions of the onset and the amount of flow separation under adverse pressure gradients [2].. Table 2 - Residual reached. ZE Boundary conditions. Umom.rms Umom.max Vmom.rms Vmom.max Wmom.rms Wmom.max P-Mass. rms P-Mass max Ke rms Ke max Ed rms Ed max Tef rms Tef max. The steady boundary conditions proposed for case 1 were used for the calculations. The geometry proposed by the organizers ends immediately after the straight diffuser. Thus, recirculation is expected. The outlet boundary condition opening was used to allow flow in both directions, see [2] for more details. Grid The grids furnished by the organizers for case 1 were used. For the ZE and k-ε models, the grid for y+=50 was used, while the grid for y+=1 was used for the SST model. The characteristics of the grids are: 1002360 nodes, 981424 hexahedral elements and a minimum face Angle of 20.8°. The first grid has a maximum edge length ratio of 4585 and a maximum element volume ratio of 10.5, while the second has a maximum edge length ratio of 132 and a maximum volume ratio of 8.9.. k-ε. SST. 9.8E-10. 3.6E-8. 1.1E-5. 2.8E-7. 1.8E-5. 2.8E-4. 1.3E-9. 6.7E-9. 1.2E-5. 4.0E-7. 2.4E-7. 3.7E-4. 1.0E-9. 1.1E-6. 2.0E-5. 3.0E-7. 2.4E-7. 3.4E-4. 5.0E-8. 1.2E-8. 9.0E-7. 1.1E-5. 3.5E-7. 8.3E-5. -. 3.4E-8. 3.9E-5. -. 2.0E-5. 5.7E-4. -. 1.8E-8. -. -. 9.7E-6. -. -. -. 1.1E-6. -. -. 3.6E-5. ZE. k-ε. SST. 9 202 258 3 50 133. 19 152 217 1 32 126. 30 148 225 0 1 6. Table 3 – y+ values y+ min hub y+ mean hub y+ max hub y+ min wall y+ mean wall y+ max wall. Discretization The schemes used for the discretization of the different equations are presented in table 1 for the three models. Convergence difficulties were encountered with the SST model. A blend factor of 0.5 was necessary to obtain convergence, i.e. hybrid scheme.. Engineering quantities The engineering quantities are presented in table 4. They are defined by:. 38.

(44) 1 Aout Cp mean =. ∫∫. Aout. 1 P dAAin. 1 ⎛ Q ⎞ ρ⎜ ⎟ 2 ⎜⎝ Ain ⎟⎠. 1. (5),. 0.8. CP. P − Pin:wall Cp wall = out :wall 2 ⎛ Q ⎞ 1 ρ ⎜⎜ ⎟⎟ 2 ⎝ Ain ⎠. ∫∫ P dA. 0.6. 0.4. Ain. (6),. 2. 0.2. 0. 2 ⎛ ⎜P + ρ U 2 ⎜⎜ Ain ⎝. ∫∫ ζ =. 2 ⎞ ⎛ ⎟ U ⋅ ndA ⎜P + ρ U ˆ − 2 ⎟⎟ ⎜⎜ ⎠ Aout ⎝. ∫∫. ∫∫ ρ. U. ⎞ ⎟ U ⋅ ndA ˆ ⎟⎟ ⎠. ZE k−ε SST 0. 0.2. 0.4 0.6 Centerline L. 0.8. 1. Figure 1 – Pressure recovery along the upper wall (plane XZ, y=0).. (7),. 2. 2. ˆ U ⋅ ndA. Ain. 1.2. 1. The ZE model predicts a relatively high pressure recovery at the wall. The value of 1.397 for the k-ε model is identical to the value obtained for the second workshop with CFX-4 and a different grid; see Cervantes and Engström [3]. The SST model has a value similar to the experimental one: Cp-exp=1.12. The mean pressure recovery is similar for the turbulence models investigated. The loss factor is as well of the same order for the three models.. CP. 0.8. 0.6. 0.4. 0.2. 0. The pressure recovery along the upper and lower center lines is presented in figures 1 and 2. The three models give similar results, indicating the small importance of turbulence. The pressure recovery is dominated by area changes.. ZE k−ε SST 0. 0.2. 0.4 0.6 Centerline L. 0.8. 1. Figure 2 – Pressure recovery along the lower wall (plane XZ, y=0).. Flow visualization Figure 3 represents the velocity contours below the cone for the turbulence models investigated. The separation on the cone is not captured by the zero equation model, while the results are qualitatively similar for the other two models. Assuming a more accurate result for the SST model, the results of this simulation shows an early separation on the cone, which is negative to the pressure recovery. It should occur further downstream.. Table 4 – Engineering quantities. ZE. k-ε. SST. Cp wall. 1.487. 1.397. 1.241. Cp mean. 0.884. 0.896. 0.906. ζ. 0.281. 0.290. 0.272. The streamlines issued from the hub are presented in figure 4 as well as the velocity contours at sections Ib, II, III and IVb. Similar pattern is obtained for the vortex rope with ZE and k-ε after the elbow. For the SST model, the rope does not deviate to the left; it is fairly straight. The difference is also present on the velocity contours at the outlet plane of the draft tube. The flow near the vortex rope has a relatively low kinetic energy. On the left side looking downstream, the fluid has a higher kinetic energy. This indicates a poor use of the straight diffuser due to the vortex rope.. 39.

(45) Figure 3 a – Velocity contours near the cone, zero equation model (plane XZ, y=0).. Figure 4 a – Streamlines issued from the hub and velocity contours at Ib, II, III, IVb (zero equation model).. Figure 3 b – Velocity contours near the cone, k-ε model (plane XZ, y=0).. Figure 4 b – Streamlines issued from the hub and velocity contours at Ib, II, III, IVb (k-ε model).. Figure 4 c – Streamlines issued from the hub and velocity contours at Ib, II, III, IVb (SST model).. Figure 3 c – Velocity contours near the cone, SST mode (plane XZ, y=0)l.. 40.

(46) The mechanical energy (E) defined as the sum of the static pressure, the mean kinetic energy and the turbulent kinetic energy: E =P+. 1 1 ρU 2 + ρ u 2 2 2. CONCLUSION Three eddy viscosity turbulence models were used to simulate the Turbine-99 test case. The results indicate qualitatively similar values for the mean engineering quantities such as the mean pressure recovery and the lost factor. This indicate the weak influence of the turbulence on the main flow, which is dominated by inertial effects.. (8),. at the plane XZ for y=0 are presented in figure 5. A stronger decrease of the mechanical energy is obtained at the inlet of the draft tube with ZE. The lowest mechanical energy is observed below the runner hub, indicating low mean velocity and pressure. The region is due to the separation on the cone, which creates a recirculation bubble below the runner hub, see figure 3. An extension of the runner hub with a lower angle may be beneficial for the present regime to decrease the separation zone. Experiments performed by Vekve [3] on a Francis model at part load indicate that an extension of the cone length may be positive to the efficiency.. Flow visualization has identified an early separation on the runner hub as well as a vortex rope issued from the runner hub, which negatively affects the function of the draft tube by blockage effects. ACKNOWLEDGEMENTS The authors are thankful to D. Marjavaara for constructive discussions.. REFERENCES [1] Cervantes M.J. and Engström F.E., Factorial Design Applied to CFD. Journal of Fluid Engineering, Vol. 126, Issue 5, pp 791-798. [2] ANSYS CFX 10.0, help manual. [3] Engström T.F., Gustavsson L.H. and Karlsson R.I., Proceedings of Turbine-99 - Workshop 2, The second ERCOFTAC Workshop on Draft Tube Flow, Älvkarleby, Sweden, June 18-20 2001. Figure 5 a – Mechanical energy contours, zero equation model (plane XZ, y=0).. [4] Vekve T., An Experimental Investigation of Draft Tube Flow, Norvegian University of Science and Technology, Doctoral theses at NTNU 2004:36.. Figure 5 b – Mechanical energy contours, k-ε model (plane XZ, y=0).. Figure 5 c – Mechanical energy contours, SST model (plane XZ, y=0).. 41.

(47) 42.

(48) 43.

(49) 44.

(50) Turbine-99 III Proceedings of the third IAHR/ERCOFTAC workshop on draft tube flow 8-9 December 2005, Porjus, Sweden Paper No. 3. CFD SIMULATIONS OF FLOW IN A HYDRAULIC TURBINE DRAFT TUBE USING NEAR WALL RANS AND LES MODELS Tokyay T. and Constantinescu G. IIHR Hydroscience and Engineering, Department of Civil and Environmental Engineering The University of Iowa, Iowa City, IA 52242 sconstan@engineering.uiowa.edu. strong induced adverse pressure gradients in the diffuser, formation and breakdown of vortices, swirl and cavitation. All these phenomena are challenging to turbulence modeling. Thus, accurate prediction of losses associated with the undesirable flow phenomena or design of draft tubes with minimized flow losses needs tools and models that can correctly and realistically account for these phenomena. Both geometric conditions (e.g., the continuous cross sectional changes form circular at the inlet to rectangular at the exit causes complexities in grid generation) and the previously discussed flow phenomena make the numerical simulation of a draft tube very challenging.. ABSTRACT This paper describes the use of near-wall RANS and LES models that do not employ wall functions to predict the flow in the draft tube of a Kaplan hydraulic turbine which was selected as a benchmark test (case 1 corresponding to mean steady inflow conditions) for the 2005 IAHR/ERCOFTAC Third Workshop on Draft Tube Flow. Fluent is used to perform the RANS simulations while a finite-volume non-dissipative LES solver that can use hybrid unstructured meshes is used to perform the LES simulation at a lower Reynolds number. The numerical results from RANS are presented in terms of relevant integral quantities of engineering interest and distribution of mean velocity, vorticity and turbulent kinetic energy (TKE) fields at relevant sections. The flow structures observed in the RANS solution and in the instantaneous LES flow fields are also discussed and compared.. The main objective of this work is to put in perspective the capabilities of two different widely used RANS models (k-ε and SST) in versions that do not employ the popular wall function approach in simulating the flow in a draft tube of a hydraulic turbine by fine enough meshes to resolve the near wall flow. The second objective of this ongoing study is try to obtain an LES solution for the same draft tube geometry at a relatively lower Reynolds number to investigate the role of coherent structures. As this flow is dominated by high swirling motions, adverse pressure gradients, potential flow separation and unsteadiness, the use of eddy resolving techniques like LES would be ideal. However due to the very high Reynolds numbers at which the draft tubes are operating, use of well resolved (no wall functions) LES at these Reynolds numbers is very expensive. So one can either use LES with wall functions or do a well resolved LES simulation at a lower Reynolds number. In the case of RANS use of near wall models without wall functions is not anymore a problem especially if parallel flow solvers are used. Interestingly, there are very few attempts to use near-wall RANS models to simulate the flow in draft tubes, though avoiding the use of wall functions was shown to improve flow predictions significantly in many complex turbulent flows.. INTRODUCTION The flow in draft tubes is of interest because of its highly complex nature and due to the significant function that draft tubes fulfill in hydro-power systems as energy extractors. A draft tube consists of a short conical diffuser followed by a strongly curved 90º elbow of a varying cross section (circular to elliptical to rectangular) and then by a rectangular diffuser section. The flow entering the draft tube is characterized as three dimensional, turbulent and highly swirling in the azimuthal direction relative to the diffuser due to the wake of the turbine runner. Moreover, the very sharp bend can cause complex 3D flow separation on the inner wall and in the corner on the outer wall side. Since energy is the main concern in draft tube flows, minimizing the losses is the big challenge in the design of draft tubes. The relative importance of the losses varies depending on the head. Particularly in the low head plants, the losses in the draft tubes become relatively large and the design and redesign of existing draft tubes become a critical issue.. In the present paper we will report results of the RANS simulations and some preliminary results from the LES simulation using a dynamic Smagorinsky model (Lilly, [4]) on a mesh that allows avoiding the. Poor inflow conditions in a draft tube can be associated with some undesirable flow phenomena such as flow reversal, unsteadiness, presence of. 45.

No results found