KOUKOUVINIS P.， BERGELES G.， GAVAISES M.

Department of Engineering and Mathematical Science， City University London， Northampton Square， London，EC1V0HB， UK， E-mail： Foivos.koukouvinis.1@city.ac.uk

A cavitation aggressiveness index within the Reynolds averaged Navier Stokes methodology for cavitating flows*

KOUKOUVINIS P.， BERGELES G.， GAVAISES M.

Department of Engineering and Mathematical Science， City University London， Northampton Square， London，EC1V0HB， UK， E-mail： Foivos.koukouvinis.1@city.ac.uk

The paper proposes a methodology within the Reynolds averaged Navier Stokes （RANS） solvers for cavitating flows capable of predicting the flow regions of bubble collapse and the potential aggressiveness to material damage. An aggressiveness index is introduced， called cavitation aggressiveness index （CAI） based on the total derivative of pressure which identifies surface areas exposed to bubble collapses， the index is tested in two known cases documented in the open literature and seems to identify regions of potential cavitation damage.

cavitation aggressiveness index （CAI）， cavitation erosion， multiphase cavitating flows

Introduction

Cavitation is a phenomenon commonly associated with the appearance of vapour cavities within the bulk of a liquid flow when pressure drops below the vapour pressure at the local liquid temperature， in such areas due to the existence of impurities in the flow field or surface irregularities， phase-change takes place， initiating from either the material surface or inside the bulk of the flow. The bubbles produced may contain both vapour and non-condensable gas （i.e. air）. Theses bubbles when advected to lower pressure regions， increase in size and if subsequently are brought abruptly into higher pressure regions in the flow， they collapse， during the collapsing process， their energy is transformed into pressure pulses which are emitted into the surrounding， causing material damage and erosion， mainly on non-deformable （metallic） surfaces. There is a substantial effort within the fluid mechanics community in developing a computational tool within the numerical algorithms currently used in cavitating flows， for predicting not only the location of the vapour bubble creation but mainly the region of its destruction， thus aiming into predicting the surface areas most prone to cavitation damage. Cavitation damage is considered to be caused by two fluid mechanics mechanisms： （1） either by the emitted pressure waves impinging on material surfaces during the bubble implosion and its subsequent rebound， according to Brennen and Leighton［1，2］， （2） by the impulse momentum of a penetrating liquid micro jet impinging on the surface and which is created due to asymmetric bubble collapse near the surface and related to a hydraulic jump pressure Hammitt［3，4］. Kedrinskii［5］， questions the potential of a single bubble collapse near a surface to create material damage as the amount of force that can be exerted on a surface by a single collapsing bubble is order of magnitudes smaller than the material hardness while in the case of the liquid micro jet impinging on the surface its kinetic energy is order of magnitudes smaller than the energy stored in the collapsing single bubble. According to Kedrinskii， it is the synchronous collapse of many bubbles which creates a cumulative pressure effect， above a threshold，which can initiate material erosion. Patella et al.［6-9］proposed an erosion model based on the pressure wave power density emitted by the imploding bubbles near a surface. The pressure wave energy emitted is obtained fromthe numerical solution of the Rayleigh-Plesset equation， with compressibility effe-cts taken into account， the time dependent wave energy results were then used for calculating the pit volume due to surface deformation of the material exposed to the pressure pulse， a linear relation was found with the coefficient of linearity depending on the material properties， the reliability of the model was also tested against experimental results on eroded samples of various metals with a good correlation. The wave power density， based on the maximum pressure emitted upon bubble collapse was found to be a good indicator of the erosion aggressiveness of the flow.

Franc and Mitchel［10］estimate the various pressures created upon bubble implosion and they conclude that the impacting pressure due to liquid micro jet may be of the same order of magnitude with the impacting pressure wave due to bubble collapse， however， its duration is orders of magnitude smaller， moreover， the mechanism of collapse of cavitating vortices， found in unsteady hydrofoil flows， is more important creating high impact loads， with order of magnitude higher duration， even higher than the collective synchronous bubble collapse. Dular et al.［11］conducted experiments in a water tunnel and presented erosion results on a cavitating hydrofoil at various flow conditions， the results obtained， mainly after visual inspection with the time of the eroded surface， were used to develop an erosion model based on the concept of the impulse pressure acting on the surface due to the liquid micro jet impingement upon asymmetric bubble collapse near the wall. The measurements indicate that erosion damage is related to cavity unsteadiness， while erosion aggressiveness is correlated to a power law of the nominal flow velocity， besides， increase of the noncondensable gas content in the water of the water tunnel reduces the erosion damage. The proposed erosion model correlates reasonably well， the impulse pressure created by the impinging liquid micro jet with the pit area created on the surface.

Dular and Coutier-Delgosha［12］used a CFD model based on RANS finite volume methodology for time dependent flows coupled with a barotropic model to account for cavitating flows， a modification was also employed to the k-εRNG model for computing turbulent viscosity which takes in an heuristic way into account mixture compressibility， the CFD results were post processed and information about pressure distribution was coupled to the erosion model， proposed in Ref.［11］. The CFD methodology predicted the various phases of cavity cloud collapse， the associated pressure signatures in time and the resulting erosion as manifested by the formation of the pits in time. Even though the predictions are not very well correlated to measurements， the coupling of CFD post processing and the erosion model gives a direction for developing such an approach for engineering design. Terwisga et al.［13］reviewed the various erosion models reported in the literature and drew the conclusion that an initial bubble implosion synchronizes the implosion of a bubble cloud and that the synchronous implosion of many bubbles either per se or after breaking up of the traveling cavitation vortices may be the physical mechanism leading to cavitation erosion. Li［14］tested various erosion functions like pressure， partial derivative of pressure with time， the rate of change of volume vapour fraction and their time integrals based on CFD results after post processing. She concluded that the maximum rate of change of pressure in time is a better criterion for cavitation damage and it seems that high values of the rate of change of vapour volume fraction do not correlate with erosion damage. Li and Terwisga［15］investigated the time derivative of pressure as an erosion risk index in a flow around a hydrofoil at 8oangle of attack and concluded that an erosion intensity function which is the time mean value of the time derivatives of pressure above a threshold correlates well with the experimental evidence of erosion risk regions on the hydrofoil surface.

1. The cavitation aggressiveness index （CAI）

The paper proposes a methodology for predicting the region of bubble collapse and its cavitation aggressiveness. The idea for the proposed methodology is based on the concept that for a vapour bubble to collapse two conditions must be met： （1） the total derivative of the vapour volume fractiona must be negative（as bubbles must have decreasing volume） or equally the total derivative of the mixture density must be positive and （2） the total derivative of the pressurep must be positive （bubbles collapse at regions of increasing pressure），

The idea of using total derivatives （and not just partial derivatives） in estimating cavitation damage originates from the fact that bubbles at the final stages of collapse follow the streamlines due to their small size and collapse along them， the use of the total derivative makes the procedure applicable not only to unsteady but also to quasi steady-state flow calculations， flow regions of high positive total pressure， calculated with a steady-state RANS methodology， might be potential areas for cavitation damage in flows with intermittent cavitation regions， thus indicating secondary areas of cavitation damage. The previous two conditions，Eqs.（1） and （2）， define a topology in the flow field where cavitation damage might appear.

1.1 The concept of the CAI

In a cavitating flow there is a region where an unsteady cavity is formed in the flow， inside of which and at some region， the total derivative of the vapour volume fraction is negative， periodically also due to the interaction of the re-entrant jet and the oscillating cavity， the cavity breaks up and a vapour cloud is advected by the flow which collapses further downstream at pressure recovery regions， i.e.， in regions characterized also by negative total derivative of vapour volume fraction. These two regions are regions vulnerable to erosion damage. The total derivative of pressure indicates the rate at which the pressure changes as a fluid particle moves along the streamline， a bubble following the streamlines having a maximum radius of Rmaxobtained at a flow region of minimum pressure Pmin（which inside the vapor cavity is equal to the vapor pressure） reduces in size and finally collapses under the influence of the pressure difference between the external pressure （Pout）， i.e.， the recovery pressure far outside the cavity， and its vapour pressure，the strength of the emitted pressure pulses depends on the imposed pressure difference and a characteristic time scale defined by the time required for bubble collapse when the environmental pressure is abruptly imposed. The time scale characteristic of the bubble collapse， is given by［1］（Rayleigh bubble collapse time），

where a=0.917. Therefore a time derivative of pressure difference characteristic of bubble collapse is of the order of

The maximum bubble radius at equilibrium can be calculated on the assumption of isothermal bubble growth and applying Eq.（2.3） of Franc and Michel［10］，is estimated to be equal to

where the initial pressure Pg0of the non-condensable gas inside the bubble upon inception（R0）is assumed to be

P0being the total pressure of the flow. From the above it can be estimated that the characteristic rate of pressure for bubble collapse is

（a is taken equal to unity）.

Similarly， within the cavity， a threshold hydrodynamic rate pressure difference due to cavity oscillations could be defined as

where lrefand urefare characteristic length and velocity （time） scales of the vapor cavity （i.e.， length of cavity and vapor velocity at the middle of the cavity）andPout，cavis the pressure outside the cavity. The two time scales introduced above， one due to bubble collapse and the other due to cavity oscillation are maybe orders of magnitude different （i.e.， μs and ms respectively）， therefore in order to bridge these two regions a logarithmic scale is proposed and a cavitation aggressiveness index （CAI） on a ten-scale is thus defined as

Values of the CAI above zero （related to cavity unsteadiness to 10 （related to bubble collapse） indicate areas of hydrodynamic cavitation aggressiveness from cavity oscillations and break up to vapor cloud collapse.

1.2Preliminary investigation of the CAI behavior in a flow

The turbulent flow around a planar-convex hydrofoil at zero and four degrees angle of attack has been chosen as a preliminary case for testing the proposed CAI， the flow has been documented experimentally by Le et al.［16］and studied numerically by Franket al.［17］using the ANSYS CFX commercial code with a homogeneous multiphase model. In both cases the flow velocity is small and cavitation damage is not expected， thus these two cases were chosen to test the behavior of the proposed Cavitation Aggressiveness Index.

Fig.1（a） Streamlines around the hydrofoil at 0oand pressure distribution

Fig.1（b） Distribution of vapour volume fraction （VOF）-two vapour cavities are formed

1.2.1Turbulent flow around the planar-convex hydrofoil at zero degree angle of attack， Le et al.［11］

The flow around the hydrofoil has been numerically investigated using the commercial software Fluent with a fine grid around the hydrofoil supplemented with local grid refinement， the k-εturbulence model was used， curvature correction and pressure gradients effects in the wall functions were also included and theyplusvalues were around 50. The inlet flow velocity was set to 5 m/s and the cavitation number equal to 0.4， the Reynolds number of the flow was 106. Details of the experimental set up can be found in Ref.［16］. For the particular flow conditions studied，according to Le et al.［16］a steady attached very thin cavity develops at the leading edge of the hydrofoil with a length around 0.4 the hydrofoil chord. Figure 1（a） shows the steady state streamlines around the hydrofoil area indicating two regions of low pressure， a small one near the leading edge of the suction side and a larger one at the pressure side near the trailing edge of the hydrofoil， the low static pressure inside these regions is associated with cavitation and Fig.1（b） indicates the corresponding fraction of the vapour volume in the two cavity regions. The vapour cavity is shorter than the flow separation region， at the suction side of the hydrofoil the vapour cavity length is in fair agreement compared to the experiments （experimental uncertainty from 2% to 10%） and to the corresponding predictions of Ref.［17］.

Fig.1（c） Distribution of positive values of total derivative of mixture density （Da/Dt＜0）

Fig.1（d） Distribution of positive values of total derivative of pressure （Dp/Dt＞0）

Fig.1（e） Regions of potential cavitation erosion-intersection of Da/Dt＜0and Dp/Dt＞0

Fig.1（f） CAI values on and around the hydrofoil

Figure 1（c） indicates the regions of negative values of the total derivative of the vapour volume fraction. Figure 1（d） indicates regions of positive values of total derivative. The Boolean intersection of the re-gions in Fig.1（c） and 1（d）， shown in Fig.1（e）， indicate the areas of bubble collapse under pressure. The corresponding CAI values， Fig.1（f）， on the hydrofoil surface are very small， apart from the tiny region at the leading edge， indicating that the hydrofoil might be free of cavitation erosion and the cloud collapse takes place at regions far from the hydrofoil surface. Time dependent calculations indicate that the flow exhibits some limited time unsteadiness， mainly at the cavity region of the pressure side with the cavity at the suction side remaining steady.

Fig.2（a） Streamlines around the hydrofoil at 4oangle of attack and pressure distribution

Fig.2（b） Distribution of VOF

Fig.2（c） Distribution of positive values of total derivative of mixture density （Da/Dt＜0）

1.2.2Turbulent flow around the plano-convex hydrofoil at 4oangle of attack， Le et al.［16］

For the experimental conditions described in Section 1.2.1 and for a cavitation number of 0.5 and angle of attack 4o， Fig.2（a） shows the time mean streamlines around the hydrofoil obtained after a time dependent solution of the flow field over many cycles，in contrast to the previous case， the suction side separation region has increased inheight whilst the pressure side separation region has disappeared in the time mean picture， the low static pressure in the region of separated flow is associated with cavitation and Fig.2（b） indicates the corresponding fraction of the vapour volume in the cavity. The time mean length of the cavity at the suction side of the hydrofoil is in qualitatively agreement with the measurements and the prediction of Frank et al.［17］but slightly underestimated.

Fig.2（d） Distribution of positive values of total derivative of pressure

Fig.2（e） Regions of potential cavitation damage

Figure 2（c） indicates the regions of positive values of the total derivative of mixture density （negative values of the total derivative of vapour volume fraction）. Figure 2（d） indicates regions of positive values of total derivative of pressure. The Boolean intersection of the regions in Figs.2（c） and 2（d） indicate the areas of bubble collapse under pressure. Figure 2（e）indicates the values of CAI distribution where it is found that cloud collapse takes place far from the hydrofoil surface except at the leading edge where high CAI values are calculated， therefore the hydrofoil under the present operating conditions should not face erosion problems apart at the leading edge. However，the flow field at this flow condition according to Le［16］，i.e.，4oangle of attack， is unsteady and Figs.2（f）， 2（g），2（h）， 2（i） show four indicative time instants of the flow development， the maximum cavity length is greater than half the hydrofoil chord.

For these particular instances the CAI values are shown in Fig.2（j） it is concluded that at the leading edge the CAI index is high indicating cavitation aggressiveness， whilst a smaller value is found during thetime dependent calculations not computed from the time mean values of the flow field as it is associated with the oscillating cavity.

Figs.2（f）-2（i） Flow development and vapor cavity at various time instants

Fig.2（j） Instantaneous distribution of CAI on the hydrofoil suction surface， for the instances above

Fig.3（a） The experimental set up， Franc et al.［18］

Fig.3（b） The erosion zones measured

2. The LEGI erosion investigation， Franc et al.［18，19］

Erosion tests on specimens have been systematically conducted at the specially build hydraulics tunnel in the university of Grenoble （LEGI）， the particular case to be examined in the present investigation is shown in Fig.3（a）. Water with a velocity of 40.5 m/s is flowing into the round pipe， where it is forced to turn abruptly into the small radial passage， upon turn， cavitation takes place and a large cavity attached to the upper part of the disk is formed. Three erosion zones are found， one just at the bent， i.e.， entrance to the narrow passage， a second one attached to the upper disk and a third one on the lower disk opposite the end of the second erosion zone. The measured erosion areas are shown in Fig.3（b）.

The numerical solution of the turbulent flow field was obtained using the Fluent software， it has been assumed that the flow is steady， incompressible and axisymmetric. The standard k-εturbulence model hasbeen utilized with non-equilibrium wall functions to account for pressure gradient effects in the wall functions. A telescopic grid refinement technique was adopted that led to grid independent solution and the ypluswas kept around 50， in the fully turbulent region. Figure 3（c） indicates pressure distribution where the pressure region below the liquid vapour pressure is evident， a vapour cavity is formed attached to the upper disk， as shown in Fig.3（d）， a large recirculation region is also formed attached to the upper disk with a length larger than the vapour cavity extending up to position 0.026 m， ending in the region of Zone 3 of erosion of the target-lower disk， the vapor cavity ends at a distance within the erosion Zone 2.

Fig.3（c） Pressure distribution and details of the streamlines structure

Fig.3（d） Distribution of the vapour volume fraction（）a

Fig.3（e） Distribution of positive values of total derivative of mixture density （Da/Dt＜0）

Fig.3（f） Distribution of positive values of total pressure gradient

Fig.3（g） Regions of potential cavitation damage

Fig.3（h） CAI values in the flow field and on the wall

Figure 3（e） presents the regions of bubble collapse， designated by the positive values of the total derivative of the mixture density and Fig.3（g） indicates the regions of possible cavitation damage on the condition of Da/Dt＜0（Boolean intersection of flow areas of Figs.3（e） and 3（f）. Figure 3（h） indicates the values of CAI. High values of CAI are found at the corner of the upper disk， location A， as also at the lower disk at two locations， locations C and D， atY= 0.021 and Y=0.027respectively， lower values of CAI around zero are found atY=0.016at the upper disk， location B， all four predicted erosion prone locations are within the erosion zones of the experiments，as shown in Fig.3（b）. The CAI value at the upper diskis associated with the vapor cavity “oscillation”， whilst the two locations on the lower disk are associated with cloud collapse. In reality the investigated flow exhibits unsteadiness and the eroded area would spread around these locations； a detailed unsteady flow calculation should indicate the extent of the areas prone to erosion， the interested reader is also addressed to the recent publications of Gavaises et al.［19，20］.

3. Conclusion

The paper has proposed a CAI based on the total derivative of pressure， properly non-dimensionalized which proved to be a reliable tool capable of identifying flow regions of potentialmaterial damage due to cavitation. The index being based on post processing results from CFD calculations using RANS solvers is suitable both for steady and unsteady flows， its possible correlation with energy release concepts upon bubble collapse could create an index of erosion rate of material.

Acknowledgments

References

［1］BRENNEN C. E. Cavitation and bubble dynamics［M］. New York， USA： Oxford University Press， 1995.

［2］LEIGHTON T. G. The acoustic bubble［M］. London，UK： Academic Press， 1994.

［3］HAMMITT F. G. Cavitation and multiphase flow phenomena［M］. New York， USA： McGraw-Hill， 1980.

［4］KIM K. H.， CHAHINE G. and FRANC J. P. et al. Advanced experimental and numerical techniques for cavitation erosion prediction［M］. Dordrecht， The Netherlands： Springer， 2014， 3-20.

［5］KEDRINSKII V. K. Bubble cluster， cumulative jets，and cavitation erosion［J］. Journal of Applied Mechanics and Technical Physics， 1996， 37（4）： 476-483.

［6］PATELLA R. F.， CHALLIER G. and REBOUD J. L. Cavitation erosion mechanism： numerical simulations of the interaction between pressure waves and solid boundaries［C］. CAV2001： Fourth International Symposium on Cavitation. Pasadena， CA， USA， 2001.

［7］PATELLA F. R.， CHOFFAT T. and REBOUD J.-L. et al. Mass loss simulation in cavitation erosion： Fatigue criterion approach［J］. Wear， 2013， 300（1-2）： 205-215.

［8］PATELLA R. F.， ARCHER A. and FLAGEUL C. Numerical and experimental investigations on cavitation erosion［C］. 26th IAHR Symposium on Hydraulic Machinery and Systems. Beijing， China， 2012.

［9］PATELLA R. F.， REBOUD J. L. and ARCHER A. Cavitation damage measurement by 3D laser profilometry［J］. Wear， 2000， 246（1-2）： 59-67.

［10］FRANC J. P.， MICHEL J. M. Fundamentals of cavitation［M］. New York， USA： Springer-Verlag， 2004， 328.

［11］DULAR M.， STOFFEL B. and ?IROK B. Development of a cavitation erosion model［J］. Wear， 2006， 261（5-6）：642-655.

［12］DULAR M.， COUTIER-DELGOSHA O. Numerical modelling of cavitation erosion［J］. International Journal for Numerical Methods in Fluids， 2009， 61（12）：1388-1410.

［13］TERWISGA T. V.， ZIRU L. and FITZSIMMONS P. et al. Cavitation erosion-A review of physical mechanisms and erosion risk models［C］. Proceedings of the 7th International Symposium on Cavitation. Ann Arbor， Michigan， USA， 2009.

［14］LI Z. R. Assessment of cavitation erosion with a multiphase Reynolds-averaged Navier-Stokes method［D］. Doctoral Thesis， Delft， The Netherlands： Delft University of Technology， 2012.

［15］LI Z.， TERWISGA T. J. C. On the capability of a RANS method to assess the cavitation erosion risk on a hydrofoil［C］. 8th International Symposium on Cavitation （CAV）. Singapore， 2012.

［16］LE Q.， FRANC J. P. and MICHEL J. M. Partial cavities：Global behavior and mean pressure distribution［J］. Journal of Fluids Engineering， 1993， 115（2）： 243-248.

［17］FRANK T.， LIFANTE C. and JEBAUER S. et al. CFD simulation of cloud and tip vortex cavitation on hydrofoils［C］. 6th International Conference on Multiphase Flow （ICMF 2007）. Leipzig， Germany， 2007.

［18］FRANC J. P.， RIONDET M. and KARIMI A. et al. Material and velocity effects on cavitation erosion pitting［J］. Wear， 2012， 274（3）： 248-259.

［19］GAVAISES M.， VILLA F. and KOUKOUVINIS P. et al. Visualisation and LES simulation of cavitation cloud formation and collapse in an axisymmetric geometry［J］. International Journal of Multiphase Flow， 2015， 68：14-26.

［20］KOUKOUVINIS P.， BERGELES G. and GAVAISES M. A new methodology for estimating cavitation erosion： Application on a high speed cavitation rig［C］. 6th European Conference on Computational Fluid Dynamics. Barcelona， Spain， 2014.

（February 20， 2015， Revised May 2， 2015）

* Biography： KOUKOUVINIS P. （1986-）， Male， Ph. D.

GAVAISES M.，

E-mail： M.Gavaises@city.ac.uk

The research leading to these results has funding from the People Programme （Marie Curie Actions） of the European Union's Seventh Framework Programme FP7/2007-2013/ under REA grant agreement No. 324313.