Chin. Sci. Bull. (2014) 59(26):3276–3282 DOI 10.1007/s11434-014-0485-1

csb.scichina.com www.springer.com/scp

Article

Engineering Thermophysics

Numerical study of unsteady cavitation on 2D NACA0015 hydrofoil using free/open source software Victor Hidalgo • Xianwu Luo • Bin Ji Alvaro Aguinaga

•

Received: 25 February 2014 / Accepted: 8 April 2014 / Published online: 17 June 2014 Ó Science China Press and Springer-Verlag Berlin Heidelberg 2014

Abstract The free/open source software as OpenFOAM, Salome meshing and python language under Debian Linux system are evaluated to model unsteady cavitating turbulent flow around a NACA0015 hydrofoil. Based on the cavitation model proposed by Kunz and large eddy simulation (LES) method, we denote the benefits of free software and open source tools as an alternative to proprietary software of computational fluid mechanics, and provided a modified cavitation model to improve numerical accuracy. The simulation results of typical break-off cycle for cavitation shedding are compared to available experimental data, and validated using image processing to find percentage of similarities. The sheet cavity frequency of 7.752 Hz was obtained based on cavitation evolution and pressure fluctuations. The study gives relevant information for CFD software development in the future. Keywords OpenFOAM  Cavitating flow  LES  Open source

1 Introduction Computational fluid dynamics (CFD) is important to solve mathematical problems and adequate use of sources takes advantage in this field, thus international companies as ANSYS on CFD presents several products for quick V. Hidalgo  X. Luo (&)  B. Ji State Key Laboratory of Hydroscience and Engineering, Tsinghua University, Beijing 100084, China e-mail: [email protected] V. Hidalgo  A. Aguinaga Mechanical Engineering Departments, Escuela Politecnica Nacional University, Quito E11-253, Ecuador

123

numerical solutions. However, this kind of commercial software is a limit for technical investigation, due to the fact that it is impossible to change the main code and sometimes is impossible to implement new solvers [1]. Instead of proprietary software, there are free open source softwares (FOSS) as OpenFoam, Salome, Python programing language which are a small group of CFD tools for numerical simulations and other programs as Libre Office and Latex for word processing. The programs are born free, so it is necessary to open the main codes for researchers to take the advantages and this is the philosophy of FOSS. In fact, there is not any kind of locks. Therefore the future of research is going to be in free/open source software development [2]. Figure 1 is a pie chart from data of CFD Online. ‘‘CFD Online is an online center for Computational Fluid Dynamics’’ [3]. This figure shows that OpenFOAM forums for research topics of CFD simulation are in a very important position of the online center, around 14 %. Because OpenFOAM is a library of C?? toolbox for numerical CFD simulations, applications and the main codes are free to change for getting new solvers which can depend on our necessities [4]. Dular et al. [5] carried out the interesting work of a cavitation erosion modelling to consider the effects of cavitation number based on image processing to illustrate experimental results. They mentioned that pit distribution after cavitating flow over hydrofoil was related with flow velocity and cavitation number, and explained those factors were involved with pressure gap. They also depicted the cavitation shedding frequency and the bubble number increased in the time. Usually, the cavity bubbles take a form of cloud, and this cloud is a challenge for numerical simulation and so difficult to be properly simulated by ANSY CFX or any proprietary CFD software. Homogeneous flow is an interesting idea presented by Kubota et al. [6]. This premise is the basis of several

Chin. Sci. Bull. (2014) 59(26):3276–3282

3277

2 Physical description 2.1 Mixture assumptions Cavitating flow is considered as a type of multiphase flow because vapour/gas bubbles tend to move together with liquid phase as travelling cavity which grows in a low pressure region [22]. However, these two phases have presented inter-phase change so that fractions for different phase are needed to analyze. Equations (1)–(3) are fundamental of multiphase flow when there is a change of phase: c¼

Fig. 1 (Color online) Traffic to different services [3]

numerical cavitation models, considering cavitating flow as multiphase flow, where Rayleigh’s eqaution is used by neglecting the slip between two phases, so the mixture of liquid-vapor is treated as a single fluid for Navier-Stokes equations. Coutier-Delgosha et al. [7] also based on this premise and proposed that condensation and vaporization are controlled by barotropic state law to understand unsteady cavitating flow and cavity shedding in a Venturitype duct. Other methods have been presented basing the cavitation model on mass transport equations, which implemented some vapour or liquid volume fraction and adhering resource term for evaporation-condensation process [8–13]. Reynolds Average Navier Stokes (RANS) equations are used to describe turbulent flows. For most industrial problems, the cavitation phenomena are usually high Reynolds number flow and the turbulence is an important issue. For cavitating flow simulations, RANS is not enough and some modifications are necessary to obtain accurate results [7, 14– 17]. In fact, Li et al. [18] used a modified k-x model combined with mass transport model of Schnerr-Sauer in a multiphase flow RANS solver. The investigation was looking for re-entrant jets, periodic shedding and collapse of cavities, formation, roll-up and transport of cavitating vortex. Nowadays, the computational calculation power has improved rapidly, so large eddy simulations have been conducted to properly model the flow unsteadiness during cavitation [19–21]. Studies of unsteady cavitating flow around a twisted hydrofoil made by Ji et al. [21], show that LES is better option than RANS to find the behaviour of cavity clouds. Based on those studies, the present research tries to validate our proposed methods by using OpenFoam and FOSS to study cavitating flow around a hydrofoil.

81 ; 8

ð1Þ

q ¼ cql þ ð1  cÞqv ;

ð2Þ

l ¼ cll þ ð1  cÞlv ;

ð3Þ

where V and V1, total mixture volume and the liquid volume respectively; c, liquid volume fraction; q, density of mixture fluid; l, molecular viscosity of mixture fluid; l and v, subscripts for liquid and vapour phase. 2.2 Governing equations In the present numerical simulation, continuity and momentum equations are governing equations. Favre filtering operation is applied to the equations to obtain Eqs. (4) and (5). Note that the over-bar denotes filtered dependent variables.   uj oq o q þ ¼ 0; ð4Þ ot oxj   ui uj oðq ui Þ o q o p o þ ¼ þ ½qðR  GÞ; ð5Þ ot oxi oxj oxj where R ¼ 2mSij , filtered viscous stress tensor, Sij ¼   ouj 1 oui i uj , þ 2 oxj oxi , rate of strain tensor, m and G ¼ ui uj  u kinematic viscosity and subgrid stress tensor respectively. Equations (4) to (5) are used in large eddy simulation (LES). 2.3 Cavitation model A mass transport cavitation model proposed by Kunz is implemented in OpenFOAM to simulate cavitating flows, which is based on the work by Merkle et al. [23]. This model considers m? and m- as creation and destruction of liquid [24] that is observed in   oðql cÞ o ql cuj þ ¼ mþ þ m : ð6Þ ot oxj Equation (7) describes the transportation rates in the cavitation model:

123

3278

Chin. Sci. Bull. (2014) 59(26):3276–3282

Fig. 2 (Color online) Illustration of Hydrofoil NACA0015

mþ ¼ m ¼

Cprod qv c2 ð1cÞ ; t1 Cdest qv cmin½0;ppv  ð12ql U12 Þt1

;

ð7Þ

where Cprod = 455 and Cdest = 4100, empirical constants based on the mean flow [25]; U?, undisturbed flow velocity [24, 26]; t? = L/U?, mean flow time scale, in which L is the characteristic length. 2.4 Flow conditions Cavitation number is calculated using pr  pv r¼ 1 2 ; 2 qU1

Fig. 3 Structural meshing around hydrofoil obtained by SALOME

ð8Þ

where pv, vapour pressure; pr, reference pressure.

edge; h, the maximum thickness as a fraction of chord; y, position along a vertical axis perpendicular to AB line. 3.2 Mesh generation

3 Geometry and solver 3.1 Hydrofoil model The model for CFD simulation is a hydrofoil/airfoil NACA0015. Figure 2 shows this hydrofoil, developed by the National Advisory Committee for Aeronautics (NACA) [27, 28]. The foil section is frequently used for designing hydraulic machine blades and aircraft wings. It should be noted that the hydrofoil geometry in Fig. 2 is plotted using python programming language. For the present study, the hydrofoil is set with attack angle of 8°. Equation (9) was developed by Bertin in the book of Aerodynamics for engineers [29] and used by University of Illinois to make a web software design to draw standard hydrofoils/airfoils [30]. The data from the equation is employed to form the hydrofoil geometry in python and Salome: rffiffiffi   x x2 h x y¼ c 0:296  0:126  0:3516 0:2 c c c x3 x4  þ0:2843 0:1015 ; ð9Þ c c where c, chord length, 200 mm in this paper; x, downstream position along the chord from the hydrofoil leading

123

OpenFOAM use finite volume method to discretize NavierStokes equations. Due to the complexity of NACA series foil shape, mesh generation along the hydrofoil surface is crucial for cavitating flow simulation so as to properly present the cavity motion. Salome version 7.2.0 is used to make structural mesh as shown in Fig. 3. The meshing is an anisotropic adaptation of hexahedrons around the NACA0015 foil. This type of mesh has been studied to especially apply for cavitating flow and to characterize multiphase flow [31, 32]. The computation domain is shown in Fig. 4a. The distances of upstream and downstream of the hydrofoil are selected according to the previous researches [18, 33]. The adoption of hexahedral meshing permits to capture the flow around a hydrofoil perfectly [34]. Principal characteristics of the structural mesh are presented at Table 1. Since OpenFOAM is for three-dimension computational simulations, a span width of 0.002 mm is applied for the present study. Temporarily, the cavitating flow is simplified under the consideration of two dimensional flows. Figure 4b shows final grid result and taper of elements for front face. The homogeneous colour in mesh generation gives information that the size of element changes is proportional and is a guaranty for proper CFD simulation.

Chin. Sci. Bull. (2014) 59(26):3276–3282

3279

Fig. 4 (Color online) Mesh generation: a Computation domain and b Taper analysis of Mesh Table 1 Anisotropic structural meshing

Table 2 Boundary conditions

Object

Mesh

Nodes

50,600

Quadrangles Hexahedrons

Inlet

Outlet

Top and bottom

Front and back

Hydrofoil surface

51,370

Velocity

Pressure

Wall

Symmetry

Wall

24,920

U? = 5.477 m/s

pr = 21.17 kPa

YPlus is a non-dimensionless number that is shown in [35]: us y ; ð10Þ yþ ¼ m

coefficient before ‘‘mag(pSat())’’ is revised from its initial value of 0.01 to 0.001 in order to better visualize the sheet cavity.

where us is the friction velocity at the nearest hydrofoil wall, y is the distance to wall, and m is kinematic viscosity. In OpenFOAM, this number is calculated by writing the yPlusLES option in the main code. The calculated YPlus around the hydrofoil surface is between 0.1 and 6.7, and the distribution is appropriate for LES.

3.4 Computation setup

3.3 Solver OpenFOAM is a group of compiled solvers and temples for CFD simulation. This group have been written in C?? and Fortran language. In this study, the initial solver interPhaseChangeFoam is modified to hyInterPhaseChangeFoam by changing the following code line:  maxðp  pSatðÞ; p0 Þ=maxðp  pSatðÞ; 0:001  magðpSatðÞÞÞ: For the above expression, ‘‘pSat()’’ means the saturated pressure, and ‘‘mag’’ represents the magnitude. Note that the

The boundary conditions are given as Table 2 in accordance to the computation domain in Fig. 4a. The reference pressure designated at the domain outlet corresponds to cavitation number of r = 1.2. Because the unsteady simulation is conducted, the time step is set as 2910-5 s, the start time is 0 s and end time is 6.5910-1 s.

3.5 Python image processing Comparison of images between pictures in the literature and graphs from the present research is used to validate CFD simulation by OpenFOAM. Therefore, a program is developed in Python language to process images. Image processing is based on ITU-R 601-2 Luma Transform. The first step is to convert an image to monochrome picture [36]. Then, the picture is filtered and

123

3280

Chin. Sci. Bull. (2014) 59(26):3276–3282

4 Results and discussion The simulation is conducted with the cavitation number of r =1.2, which is the same as many previous studies [39]. In Fig. 5, the simulation results listed at the left column, are compared to experimental results from Yakushiji et al. [40] listed at the right column. Figures 5 and 6 show a regular break-off cycle of cavitation from 0.141 to 0.270 s. The following milestones describe this process: (1)

(2)

(3) Fig. 5 (Color online) Typical cycle of sheet cavity

Fig. 6 (Color online) The cloud of bubbles

Once cavity shedding occurs, the pressure close to hydrofoil wall presents peaks and fluctuates because of cavitation dynamics. Thus, a time period known as breakoff cycle can be defined. Figure 7 shows the typical breakoff cycle for the predicted pressure fluctuations at the point close to leading edge along the hydrofoil wall, i.e. x/c=0.2. For better understanding, the reference pressure level of pr and inferior limit pv are shown in Fig. 7, and six instants are highlighted in Figs. 5, 6 and 7. Based on these outcomes, the follows can be seen: (1)

divided into many areas with different gray scale, which corresponds to the movement of the cavities. The last step is to compare these areas. The algorithm is generated using the package/libraries as Python Imaging Library (PIL), Numerical Python (numpy) and Matplotlib [37, 38].

Cavity expansion at rear part induces a re-entrant jet, which is observed as small cavity separation from the hydrofoil suction wall at instant t =0.141 s. When the re-entrant jet approaches the hydrofoil leading edge, it breaks off the sheet cavity and makes the cavity shedding. Then, cloud cavity appears at instant t =0.228 s. At 0.270 s, a new sheet cavity starts to grow, while the cloud cavity moves downstream with the main flow.

(2)

There is a pressure pulse near the instants of 0.155 s, which represents a small cavity shedding near hydrofoil leading edge as shown in Fig. 5. There are several pressure pulses in the period of 0.24–0.26 s, and the highest peak occurs before 0.26 s with the value of 63.0 kPa. These peaks probably

Fig. 7 (Color online) Predicted pressure fluctuations in a typical break-off cycle

123

Chin. Sci. Bull. (2014) 59(26):3276–3282

3281

Fig. 8 Cavity comparison between present result and the experiment [40], at an instant of t = 0.228 s

mean the pressure pulses induced by the collapses of cavity cloud downstream the hydrofoil. (3) The break-off cycle is defined to explain the evolution of a sheet cavity and cloud cavity by f ¼ 1=Dtcycle ; ð11Þ where f, frequency. Dtcycle, the estimated time that is close to 0.155 s. Thus, the frequency of cavitation evolution is 7.752 Hz. Figure 8 shows the cavities at the instant of t =0.228 s. By using the program described in subsection 3.5, the image comparison between the numerical prediction and the experiment [40] is conducted to check the similarity. It is noted that the differences between the simulation and the experiment are mainly found at the sheet and cloud cavity areas. The total difference in area is 14.07 %. Thus, the cavitation evolution over the hydrofoil is reasonably simulated. Nevertheless, the numerical accuracy of the present method should be improved, because both pressure fluctuation and cavitation erosion are very sensitive to the cavity evolution. In the future, the authors would like to modify cavitation model by changing the constants Cprod, Cdest, and to refine the mesh around the hydrofoil surface. Further, the mathematic treatments applied by commercial CFD codes should be implemented to enhance the robustness and to reduce the numerical uncertainty for the present method.

5 Conclusions Based on the present study, the following conclusions can be drawn: (1)

(2)

The cavitating turbulent flow simulation in this investigation presented fairly good results for hydrofoil cavitation. The simulation predicts a typical break-off cycle with a frequency of 7.752 Hz. In this cycle, there are several peaks and valleys of pressure fluctuations, which would relate with the change of phases.

(3) (4)

It is observed that re-entrant jet has a great influence on the cavity shedding from the hydrofoil wall. The reasonable numerical simulation of cavitating turbulent flow around the hydrofoil verifies that free open source software, such as OpenFOAM, Salome, etc. are good computational fluid dynamics tools, which can reduce the research budget and allow any customized changes for cavitation studies.

Acknowledgments This work was financially supported by the National Natural Science Foundation of China (51206087, 51179091 and 51376100) and the Major National Scientific Instrument and Equipment Development Project (2011YQ07004901). Conflict of interest Victor Hidalgo, Xianwu Luo, Bin Ji and Alvaro Aguinaga declare that they have no conflict of interest.

References 1. Fitzgerald B (2006) The transformation of open source software. MIS Quarterly 30:587–598 2. Scacchi W (2010) The future of research in free/open source software development. In: Proceedings of the FSE/SDP workshop on Future of software engineering research, pp 315–320 3. Jones P (2013) Web server statistics of CFD-Online. cfdonline.com 4. Sasongko N, Arif M (2009) Open source computational fluid dynamic: challenges and its future. Conference on Open Source, Jakarta, Indonesia 5. Dular M, Stoffel B, Sirok B (2006) Development of a cavitation erosion model. Wear 261:642–655 6. Kubota A, Kato H, Yamaguchi H et al (1992) A new modelling of cavitating flows: a numerical study of unsteady cavitation on a hydrofoil section. J Fluid Mech 240:59–96 7. Coutier-Delgosha O, Reboud J, Delannoy Y (2003) Numerical simulation of the unsteady behaviour of cavitating flows. Int J Numer Methods Fluids 42:527–548 8. Ji B, Luo X, Arndt R et al (2014) Numerical simulation of three dimensional cavitation shedding dynamics with special emphasis on cavitation-vortex interaction. Ocean Eng. http://dx.doi.org/10. 1016/j.oceaneng.2014.05.005 9. Singhal A, Athavale M, Li H et al (2002) Mathematical basis and validation of the full cavitation model. J Fluids Eng 124:617–624 10. Kunz R, Boger D, Stinebring D et al (2000) A preconditioned NavierStokes method for two-phase flows with application to cavitation prediction. Comput Fluids 29:849–875

123

3282 11. Huang B, Wang G (2011) A modified density based cavitation model for time dependent turbulent cavitating flow computations. Chin Sci Bull 56:1985–1992 12. Zhang X, Wu Z, Xiang S et al (2013) Modeling cavitation flow of cryogenic fluids with thermodynamic phase-change theory. Chin Sci Bull 58:567–574 13. Cao X, Zhang X, Qiu L et al (2009) Validation of full cavitation model in cryogenic fluids. Chin Sci Bull 54:1633–1640 14. Hunag B, Wang G (2011) Partially averaged navier-stokes method for time-dependent turbulent cavitating flows. J Hydrodyn 23:26–33 15. Ji B, Luo X, Peng X et al (2012) Numerical analysis of cavitation evolution and excited pressure fluctuation around a propeller in non-uniform wake. Int J Multiphase Flow 43:13–21 16. Huang B, Wang G, Zhao Y et al (2013) Physical and numerical investigation on transient cavitating flows. Sci China Tech Sci 56:2207–2218 17. Ji B, Luo X, Wu Y et al (2013) Numerical analysis of unsteady cavitating turbulent flow and shedding horse-shoe vortex structure around a twisted hydrofoil. Int J Multiphase Flow 51:33–43 18. Li D, Grekula M, Lindell P (2010) Towards numerical prediction of unsteady sheet cavitation on hydrofoils. J Hydrodyn 22: 741–746 19. Luo X, Ji B, Peng X et al (2012) Numerical simulation of cavity shedding from a three-dimensional twisted hydrofoil and induced pressure fluctuation by large-eddy simulation. J Fluids Eng 134:041202 20. Lu N, Bensow R, Bark G (2010) LES of unsteady cavitation on the delft twisted foil. J Hydrodyn 22:784–791 21. Ji B, Luo X, Peng X et al (2013) Three-dimensional large eddy simulation and vorticity analysis of unsteady cavitating flow around a twisted hydrofoil. J Hydrodyn 25:510–519 22. Bensow R, Bark G (2010) Simulating cavitating flows with LES in OpenFoam. V European conference on computational fluid dynamics, pp 14–17 23. Merkle C, Feng J, Buelow P (1998) Computational modeling of the dynamics of sheet cavitation. 3rd International symposium on cavitation, Grenoble, France, vol 2, pp 47–54 24. Morgut M, Nobile E (2011) Influence of the mass transfer model on the numerical prediction of the cavitating flow around a marine propeller. Proceedings of 2nd international symposium on Marine Propulsors, pp 1–8

123

Chin. Sci. Bull. (2014) 59(26):3276–3282 25. Lu N (2008) Solve cavitating flow around a 2D hydrofoil using a user modified version of interPhaseChangeFoam. Tutorial, SE 412:96 26. Vallier A, Nilsson H, Revstedt J (2010) Mass transfer cavitation model with variable density of nuclei. 7th International conference on multiphase flow, ICMF, Tampa, USA 27. Cengel Y, Boles M, Kanonuglu M (2011) Thermodynamics: an engineering approach, vol 5. McGraw-Hill, New York 28. Habis T (1933) The 7 by 10 foot wind tunnel of the national advisory committee for aeronautics 29. Bertin J, Smith M (1989) Aerodynamics for engineers, vol. 1. Prentice Hall, Englewood 30. Islam M, Aim R, Carriveau R et al (2009) Investigation of low Reynolds number airfoils for fixed-pitch straight-bladed. VAWR, AIAA aerospace sciences meeting including the new horizons forums and aerospace exposition, pp 5–8 31. Frey P, Alauzet F (2005) Anisotropic mesh adaptation for CFD computations. Comput Methods Appl Mech Eng 194:5068–5082 32. Muller J (2002) Anisotropic adaptation and multigrid for hybrid grids. Int J Numer Meth Fluids 40:445–455 33. Ji B, Luo X, Wu Y et al (2012) Partially-averaged Navier-Stokes method with modified k-epsilon model for cavitating flow around a marine propeller in a non-uniform wake. Int J Heat Mass Tran 55:6582–6588 34. Shengwei M (2011) Top 5 misunderstandings on (good) mesh. caewatch.com 35. Ari M, Salim S, Cheah S (2009) Wall y? approach for dealing with turbulent flow over a surface mounted cube: part 1—low Reynolds number. Proceedings of 7th international conference on CFD in the minerals and process industries, Melbourne, Australia 36. Jacquin A (1992) Image coding based on a fractal theory of iterated contractive image transformations. IEEE Trans Image Process 1:18–30 37. Bird S, Klein E, Loper E (2009) Natural language processing with Python. O’Reilly Media, Inc., Sebastopol 38. Hunter J (2007) Matplotlib: a 2D graphics environment. Comput Sci Eng 9:90–95 39. Dular M, Bachert R (2009) The issue of strouhal number definition in cavitating flow. J Mech Eng Slovak Acad Sci 55: 666–674 40. Yakushiji R, Yamaguchi H, Maeda S et al (2001) Investigation for unsteady cavitation and re-entrant jet on a foil section. Jpn J Soc Naval Archit 190:61–74