B9.Unsteady Flow Analysis in Hydraulic Turbomachinery.pdf

Unsteady Flow Analysis in Hydraulic Turbomachinery Albert Ruprecht Institute of Fluid Mechanics and Hydraulic Machinery University of Stuttgart Germany ABSTRACT In the field of hydraulic machinery Computational Fluid Dynamics CFD is routinely used today in re search and development as well as in design At that nearly always steady state simulations are ap plied In this paper however unsteady simulations are shown for different examples The presented examples contain applications with self excited unsteadiness e g vortex shedding or vortex rope in the draft tube as well as applications with externally forced unsteadiness by changing or moving geometries e g rotor stator interactions For these examples the requirements potential and limita tions of unsteady flow analysis assessed Particularly the demands on the turbulence models and the necessary computational efforts are discussed INTRODUCTION For more than a decade Computational Fluid Dynamics CFD is used in the field of hydraulic machinery in research and development as well as in the daily design busi ness Early successful demonstrations are given e g in the GAMM workshop 1 The applications are steadily increasing This is expressed in fig 1 where the per centage of papers dealing with CFD is shown which were presented at the IAHR Symposium on Hydraulic Machinery and Cavitation Starting with Q3D Euler and 3D Euler today usually the Rey nolds averaged Navier Stokes equations together with a robust model of turbulence usually the k model is used It is common practice to apply steady state simulations the unsteadiness in consequence of the rotor stator interactions is addressed by averaging proce dures By this method accurate results are obtained for many questions in the design of com ponents However different problems in turbomachinery arise from un steady flow phenomena In or der to get information on this phenomena or solutions to the problems an unsteady flow analysis is necessary This requires a much higher computational effort roughly a factor 5 10 compared to steady state depending of the problem and of the degree of modeling assumptions With today s computers and software however unsteady problems can be solved Fig 1 Percentage of papers at the IAHR Symposium dealing with CFD Two major groups of unsteady problems can be distinguished The first group are flows with an externally forced unsteadiness This can be caused by unsteady boundary conditions or by changing of the geometry with time Examples are the clo sure of a valve the change of the flow domain in a piston pump or the rotor stator interactions The second group are flows with self excited unsteadiness which are e g turbulent motion vortex shedding Karman vortex street or unsteady vortex be havior e g vortex rope in a draft tube Here the unsteadiness is obtained without any change of the boundary conditions or of the geometry There can also occur a combination of both groups e g flow induced vibrations change of geometry caused by vortex shedding All these phenomena can take place in a turbine or pump and require different solution procedures BASIC EQUATIONS AND NUMERICAL PROCEDURES In hydraulic turbomachinery today usually the Reynolds averaged Navier Stokes equations for an incompressible flow are applied Compared to the steady state the momentum equations contain an additional term prescribing the unsteady change 0 x U x U xx P1 x U U t U ij i j j i jij i j i 1 ij are the Reynolds stresses which are calculated from the turbulence model The continuity equation for incompressible flow reads 0 x U i i 2 and does not contain a time depending term It has to be emphasized that the equations 1 and 2 behaves different in time and in space In space they show elliptic behavior therefore they require boundary conditions on all surfaces In time how ever they are of parabolic nature which mean that there is no feed back from the future to the pres ent or past Because of that no boundary condi tions are required in the future This is schemati cally shown in fig 2 This is the reason why the time discretization is generally carried out in a dif ferent way than the spatial discretization For spa tial discretization usually a Finite Volume or a Finite Element approximation is applied For time discre tization however mostly the Finite Difference method is used A few of the most popular finite difference approximations are shown in fig 3 In addition explicit multi point schemes of Runge Kutta type or predictor corrector schemes are often applied Fig 2 Boundary and initial conditions Fig 3 Time discretization schemes It has to be mentioned that the explicit methods require a restriction of the time step according to stability criteria CFL criteria which depend on the local velocities and the local grid size The implicit methods in contrary are always stable there is no restriction of the time step It can be chosen only according to the physical require ments In order to obtain accurate solutions the time discretization should be at least of 2nd order similar to the spatial discretization Otherwise extremely small time steps would be required The above description of the flow in the Eulerian coordinates can be applied for un steady boundary condition problems as well as for self excited unsteadiness How ever to express problems with moving geometries in Eulerian coordinates is more difficult At the moving boundary a Lagrangian description can be applied very easily since the fluid particles can be traced by this method Combining these two methods an Arbitrary Lagrangian Eulerian ALE method can be utilized This method is suit able for the solution of problems with moving boundaries In the ALE method the ref erence coordinates can be chosen arbitrary In this referential coordinate system the material derivative can be described as j E i jj R i L i x t xf wu t t xf t t xf 3 with the coordinates scooddinateEulerian x scooddinatelreferentia x scooddinateLagrangian x E i R i L i and Wi reference velocity The momentum equations in the ALE formulation can be written as follows 0 x U x U xx P1 x U WU t U ij i j j i jij i jj i 4 The moving of the reference system Wi can be chosen arbitrary If Wi is equal to zero one gets the Eulerian description on the other hand if wi is equal to the velocity of the fluid particle the Lagrangian formulation is obtained The convective term in the transport equations for scalar quantities changes in the same way than in the mo mentum equations This applies also to the k and equations The numerical realization of moving or changing grids can either be obtained by de formation of an existing mesh in each time step For large deformations this requires an automatic grid smoothing algorithm or even an automatic remeshing after a few time steps An other method is the use of different embedded grids which can move against each other In this case a sliding interface between the non matching grids is required This procedure is schematically shown in fig 4 for two different problems namely rotor stator interaction and vibration of a cylinder in a fluid In FENFLOSS the computer code developed at our institute at University of Stutt gart the second approach is applied The interface between the grids is realized by means of dynamic boundary conditions where downstream the node values veloci ties and turbulence quantities are prescribed and upstream pressure and fluxes are introduced as surface conditions A brief over view on the numerical procedures is given in 2 for more details the reader is referred to 3 4 One point has to be em phasized Since the un steady simulations re quire a severe increase of computational effort compared to steady state solutions parallel proce dures are necessary In this case the ALE formulation with moving grids leads to a dynamic change of com munication because the location of exchange boundaries varies with time and can therefore change the computational domain of the processors see 2 In FENFLOSS an implicit solution algorithm is applied As already mentioned this has the advantage that there is no stability limitation for the time step The overall solution procedure including the fluid structure interaction is shown in fig 5 If the movement of the grid does not depend on the flow situation the fluid structure loop vanishes Fig 5 Flow chart of FENFLOSS including fluid structure interaction Fig 4 Moving grid examples APPLICATIONS In the following selected applications are shown and the specific problems for this examples are discussed Firstly some cases with self excited unsteadiness are pre sented Vortex shedding at the inlet of a power plant Problem description The first example shows the flow behavior at the inlet of a low head power plant It is an existing plant with two identical bulb turbines During op eration the inner turbine showed severe bearing problems whereas the outer turbine operates smoothly The reason was expected to be vortex shedding at the inlet By numerical analysis the problem was investigated and it was tried to find a solution to the problem In fig 6 the geometry is shown The calculation has been carried out in 2D as well as in 3D Firstly it was tried to carry out a steady state simulation how ever no converged solution could be obtained Therefore an unsteady simula tion was undertaken The results indicate a strong unsteady motion In fig 7 the velocity distribution at a certain time step is presented Clearly visible are the vor tices shedding from the inlet and moving downstream into the inner turbine This is the reason of the destruction of the bearings In or der to improve the flow behavior a modified ge ometry was sug gested This ge ometry shown in fig 8 has been built in the meantime There are no longer problems with vortex shedding Further details about this application can be found in 5 6 Discussion The physical unsteadi ness of the flow has been indicated by the inability to achieve a con verged steady state solution This is very often the case with flows showing vortex shedding in reality Fig 6 Geometry of power plant inlet Fig 7 Instantaneous velocity vectors vortex shedding at the inlet pier Fig 8 Modified geometry A necessary condition for that is that the numerical scheme does not contain seri ous artificial diffusion which would suppress the unsteady motion The same applies to the used turbulence model The standard k model usually produces a too high eddy viscosity especially in swirling flows and therefore it very often suppresses the unsteady motion This will be discussed again in other applications For many cases at least a streamline curvature correction or even a non linear eddy viscosity formu lation is necessary in order to avoid a too high turbulence production Another point in turbulence modeling is the treatment of the near wall flow It is well known that the use of wall functions usually tends to predict a flow separation too late In case of vortex shedding this can cause a severe reduction of the vortex sizes or even a complete suppression of the vortices More accurate results can be ob tained by solving the flow up to the wall if possible by a low Reynolds or a two layer model The results shown above are achieved by an algebraic turbulence model Baldwin Lomax type where the flow is resolved up to the wall Vortex rope in a draft tube Problem description As an other self excited unsteady flow example the simulation of a vortex rope in a draft tube is shown Here a straight axisymmetrical diffuser is considered The inflow conditions to the diffuser are chosen according to the part load operation of a Francis turbine This means that the flow shows a strong swirl component The inlet velocity distribution and the geometry are presented in fig 9 The instantaneous flow for a certain time step is given in fig 10 where an iso pressure surface as well as the secondary velocity vectors in three cross sections are plotted Clearly the cork screw type flow with an unsymmetrical form is visible although the geometry and the boundary conditions are completely axisymmetrical Fig 9 Geometry and inlet conditions Fig 10 Iso pressure and secondary flow of a vortex rope In fig 11 the secondary velocity and the low pressure region which represents the vortex center is shown in the cross section S indicated in fig 9 for certain time steps Clearly the revolution of the vortex center can be observed This of course causes pressure fluctuations and therefore dynamical forces on the draft tube sur face Fig 11 Secondary motion and low pressure region for different time steps Discussion Concerning the numerical scheme and the turbulence models the dis cussion above also applies here e g application of the standard k model leads to a steady state symmetrical solution This is also reported in 7 The results shown above are achieved by applying the multi scale k model of Kim 8 together with a streamline curvature correction This model shows a much lower eddy viscosity than the standard model especially in swirling flows The application of wall functions does not give any problems here since the flow instability has its origin in the center and is not affected by the prediction of the near wall region Vortex instability in a pipe trifurcation Problem description In the following another problem caused by a vortex instability is shown It is a pipe trifurcation which is es tablished in a power plant in Nepal The tri furcation distributes the water from the pen stock to the three turbine units The prob lem in this plant arises from severe fluctua tions of the power output of the both outer turbines By field measurements the trifur cation was discovered as the reason for the fluctuations By means of CFD and by model tests carried out at ASTROE in Graz the flow behavior should be analyzed and a cure of the problem should be found The geometry of the trifurcation is shown in fig 12 It has a spherical shape The fluctuation in the trifurcation is caused by a strong vortex which tends to be unstable It skips between the two situations sketched in fig 13 In the model tests the secondary velocity of the vortex could be found to be 30 times higher than the transport velocity The reason is that at the top of the sphere there is enough space for a huge vortex to form This vortex concentrates in the side branches and there fore increases the swirl intensity Because of this strong secondary motion there are strong losses at the inlet of the branch which reduces the head of the turbine and therefore causes the reduction of power output During the project it was tried to obtain the unsteady behavior by a k simulation on relatively coarse grids 200 300 000 nodes However these calculations did not show the vortex instability Merely a vortex forms which extends from one side branch to the other The swirl intensity was underpredicted by more than a factor five Because of the low swirl rate the vortex is completely stable and has no tendency of skipping between different stations Even by a dynamical excitation caused by changes of the outlet boundary condition of one branch the predicted vortex did not change its position Only when applying finer grids and another turbulence model the predicted swirl in tensity could be increased Here an algebraic turbulence model with a limitation of the eddy viscosity is applied The used grids consists of about 500 000 nodes As a consequence this leads to an instability of the vortex In the prediction the vortex skips between the two structures shown in fig 14 One of these structures corre sponds quite well with the structure observed in the model tests In the second situa tion the vortex expends from one side branch to the other This complies with the above mentioned stable results The calculated swirl intensity is still more than two times lower compared to the results of the model tests Therefore further investiga Fig 12 Geometry of the trifurcation Fig 13 Vortex structure tions with other turbulence models and with finer grids are necessary and will be car ried out in future Fig 14 Predicted vortex structures For completeness the solution to the problem is shown It consists of the installation of two plates in the upper and lower part of the sphere This is shown in fig 15 Hence no free space is available where the vortex can form Consequently the in tensity of the vortex is dramatically reduced and the vortex is completely stable In the meantime the reconstruction was carried out and the fluctua tion of the power output vanished As a by product the losses in the trifurcation are severely reduced which results in an increase of power output of ap proximately 5 Further details of this problem can be found in 9 10 Discussion As already mentioned the calculations using the k model were not successful It is well known that this model is not able to predict highly swirling flows accurately The unsteady motion of the vortices especially of very slim vortices however very much depends on the swirl intensity In order to prescribe such types of flow with sufficient accuracy it is necessary to have highly sophisticated turbulence models and very fine grids maybe the only way to achieve it is the application of large eddy simulation Rotor stator interaction in an axial tubine The following ex ample belongs to the second group the unsteadiness is forced by mov ing geometries The problem in question is the Fig 15 Modified geometry Fig 16 Geometry of the investigated axial turbine flow in an axial turbine The speciality of this turbine is its relatively low specific speed It has been designed for pressure recuperation in piping systems The ad vantage is that the discharge is nearly independent of the speed because of that the turbine cannot introduce waterhammers in the system The geometry of the turbine is shown in fig 16 It consists of the inlet confuser 12 fixed guide vanes 15 runner blades and the draft tube The stator and rotor part is shown in more detail in fig 17 For the simulation the complete t