Numerical investigations of the transient cavitating vortical flow structures over a flexible NACA66 hydrofoil ＊

Ren-fang Huang, Te-zhuan Du, Yi-wei Wang, Cheng-guang Huang

Key Laboratory for Mechanics in Fluid Solid Coupling Systems, Institute of Mechanics, Chinese Academy of Sciences, Beijing 100190, China

Abstract: In this paper, the cavitating flow over a flexible NACA66 hydrofoil is studied numerically by a modified fluid-structure interaction strategy with particular emphasis on understanding the flow-induced vibration and the cavitating vortical flow structures.The modified coupling approaches include (1) the hydrodynamic solution obtained by the large eddy simulation (LES) together with a homogenous cavitation model, (2) the structural deformation solved with a cantilever beam equation, (3) fluid-structural interpolation and volume mesh motion based on the radial basis functions and greedy algorithm. For the flexible hydrofoil, the dominant flow-induced vibration frequency is twice of the cavity shedding frequency. The cavity shedding frequency is same for the rigid and flexible hydrofoils, demonstrating that the structure vibration is not large enough to affect the cavitation evolution. The predicted cavitating behaviors are strongly three-dimensional, that is, the cavity is (a) of a triangular shape near the hydrofoil tip, (b)of a rectangular shape near the hydrofoil root, and (c) with a strong unsteadiness in the middle of the span, including the attached cavity growth, oscillation and shrinkage, break-off and collapse downstream. The unsteady hydroelastic response would strongly affect the cavitation shedding process with small-scale fragments at the cavity rear part. Furthermore, three vortex identification methods (i.e., the vorticity, the Q-criteria and the Ω method) are adopted to investigate the cavitating vortex structures around the flexible hydrofoil. It is indicated that the cavity variation trend is consistent with the vortex evolution. The vortex structures are distributed near the foil trailing edge and in the cavitation region, especially at the cavity-liquid interface. With the transporting downstream the shedding cavities, the vortices gradually increase in the wake flows.

Key words: Cloud cavitation, flexible hydrofoil, Ω method, flow-induced vibration

Introduction

The cavitation occurs when the local pressure drops below the liquid saturation vapor pressure. It is a common phenomenon in water turbines, marine vehicles, propellers and valves, etc. The cavitation oscillation, break-off and collapse process may cause many problems such as the pressure pulsation, the structural vibration, the noise and the surface erosion.Due to its importance in a wide range of fundamental studies and engineering applications, the studies of the cavitation dynamics were comprehensively reviewed in literature[1-4]. Due to the cavitation unsteadiness including the cavity breakdown and collapse, strong instantaneous loads will be produced, causing further hydrodynamic instabilities, even structural failures.Meanwhile, the flow-induced vibration will in turn affect the transient cavitating flows. Therefore, it is of great significance to investigate the flow-induced vibrations in unsteady cavitating flows.

Many experiments[5-7]and numerical simulations[8-10]were conducted to study the cavitation dynamics and the related structure vibrations. Amromin and Kovinskaya[11]analyzed the vibrations of an elastic wing in an attached cavitating flows, with the wing vibration solved based on the beam equation,. It is indicated that the high-frequency band of a significant vibration is related with the elastic resonance and the low-frequency band corresponds to the cavity-volume oscillations governed by the cavity length-based Strouhal number.

With the measurements of the displacements on the free foil tip section using a high speed video camera and the surface velocity vibrations using a Doppler vibrometer, Ducoin et al.[7]studied the fluid-structure interaction on a rectangular cantilevered flexible hydrofoil in the cavitating flows and found that the cavitation greatly increases the vibration level due to the hydrodynamic loading unsteadiness and the change of the modal responses of specific frequencies. Wu et al.[6,8,12]investigated the cavitating flow-induced vibrations of a flexible NACA66 hydrofoil through both experiments and numerical simulations, and it was found that the maximum vibration amplitude keeps relatively small for the inception and sheet cavitation, it increases dramatically for the cloud cavitation and declines for the supercavitation. De La Torre et al.[13]conducted a large range of experiments for a two-dimensional NACA0009 truncated hydrofoil submerged in the air,the still water and in various cavitation regimes, and a linear correlation was found between the added mass and the entrained mass related to the cavity, including the density, the dimension and the location relative to the specific modal shape deformation. Smith et al.[5]experimentally measured the cloud cavitation behaviors around the NACA0009 hydrofoil, and it is demonstrated that the hydrofoil compliance would damp the high frequency force fluctuations and is closely associated with the normal force and the tip deflection. With focus on the effects of the cavitation and the fluid-structure interaction on the mechanism of the vortex generation, Ausoni et al.[14]found that the vortex shedding frequency increases up to 15% in a fully developed cavitation and the cavitation onset will significantly increase the vortex-induced vibration level due to the increase of the vorticity caused by the cavitation. For the cross flow around a two-dimensional elastic cylinder, So et al.[15]investigated the free vibrations and their effects on the near wake flows by using a laser vibrometer to measure the bending displacements and a laser Doppler anemometer to measure the velocities. It was shown that the cylinder vibrations have little or no effect on the mean drag and the normalized mean field,but they enhance the turbulent mixing, with a substantial increase of the turbulent intensities. The studies of Ausoni et al.[14]and So et al.[15], indicated that the structure vibration is closely associated with the large-scale vortical motion. Moreover, the cavitation would promote the vortex production with the contributions from the vortex stretching, the vortex dilatation and the baroclinic torque terms[16-17].

The accurate vortex identification is essential for better understanding the structural vibration characteristics and the transient cavitation patterns, as well as the vortex dynamics. There are various traditional vortex visualization methods including the2λmethod[18], theQ-criterion[19], and the Lagrange method[20-21]. Generally speaking, the2λmethod and theQ-criterion involve the subjective selection of a threshold for the vortex visualization, varying from case to case and with no universally acknowledged guidelines. In view of these problems, Liu and his collaborators recently proposed new vortex identification methods, that is, theΩmethod[22]and the Liutex/Rortex[22]. Zhang et al.[23]reviewed various vortex identification methods with applications of the wake flow around the moving bodies, the atmosphere boundary layer and the reversible pump turbine. It is found that theΩmethod is superior to other methods whenΩ= 0.52, as is recommended to define the vortex boundary. Different vortex identification methods were used by Wang et al.[24]to study the cavitation vortex dynamics in the unsteady sheet/cloud cavitating flows with shock waves,indicating that different methods can identify vortices in different cavitation regions. However, there are few studies about the application of theΩmethod and theQ-criterion to the cavitating flow with consideration of the fluid-structure interaction.

1. Numerical approach

1.1 Governing equations

The unsteady Navier-Stokes equations are solved by using the large eddy simulation (LES) method and a mass transfer cavitation model. In the homogenous framework, the vapor/liquid fluid components are assumed to share the same velocity and pressure fields.The basic governing equations consist of the mass and momentum conservation equations.

whereuiis the velocity component in theidirection,pis the pressure and the laminar viscosityμand the mixture densityρa(bǔ)re defined as

whereαis the volume fraction of the different phases and the subscriptslandvrepresent the liquid water and the water vapor, respectively.

The vapor volume fractionvαis governed by the cavitation model developed by Schnerr and Sauer[25].

The bubble radius is related to the vapor volume fractionvαand the bubble number densityNbas follows

Herein,Nbthe only parameter which needs to be specified, and we letNb=103. The saturated vapor pressure ispv=3540 Pa . The model constants are based on the work of Schnerr and Sauer[25]. This cavitation model has been validated by many cases,such as the cavitating flows around a threedimensional hydrofoil[16,26-27].

Applying the Favre-filtering operation to Eqs. (1)and (2), the LES equations are obtained as:

where the over-bars denote filtered quantities. The non-linear termijτin Eq. (10) is called the sub-grid scale (SGS) stress, which is defined as

The SGS stress is modeled by the eddy viscosity model, where it is assumed that the SGS stresses are proportional to the modulus of the strain rate tensor,, of the filtered large-scale flow.

whereτkkis the isotropic part,is the rate-ofstrain tensor in the resolved scale,tμis the sub-grid scale turbulent viscosity and with the LES wall adapting local eddy-viscosity (WALE) model[28], we have:

where Δsis the sub-grid scale mixing length,kis von Karman’s constant,dis the distance to the closest wall,Vis the volume of the grid cell andCs=0.325 is the WALE constant based on calibrations using freely decaying isotropic homogeneous turbulence[28].

1.2 Fluid-structure interaction (FSI) technique

During the fluid-structure interaction (FSI), the flexible structure will deform under the hydrodynamic forces, and the structure deformations will in turn alter the surrounding flows and change the hydrodynamic forces, until a static equilibrium is reached. As the FSI technique is the key issue in the CFD-CSD coupling approach[29], Huang et al.[30]developed a coupling strategy including the Fluid-structure interpolation and the volume mesh motion schemes based on radial basis functions. Their coupling strategy enjoys a good performance in studying the aeroelastic effects on the grid fin aerodynamics in transonic and supersonic regimes. Inspired by their work, a modified fluid-structure technique is proposed in this paper to solve this static hydro-elastic problem by incorporating a user defined function (UDF) code in the commercial software FLUENT, where the hydrodynamic forces and the structural motions are solved separately. The iterative scheme is as follows:

(1) Compute the hydrodynamic forces using the governing equations given in Section 1.1.

(2) Interpolate the forces from the hydrodynamic nodes onto the structural nodes. It is noted that the hydrodynamic nodes are not the same as the structural nodes at the fluid-structure interface. Therefore,according to the physical laws including the conserva-tion of the total force, torque and energy, it is necessary to realize the transfer of the force and the displacement between the hydrodynamic nodes and the structural nodes by using the radial basis function(RBF) interpolation[30].

(3) Calculate the structural deformation using the flexibility method. The elastic hydrofoil here is represented by a cantilever beam. The natural frequencyiωand its corresponding normalized shape functionare as follows:

(4) Interpolate the deformation from the structural nodes onto the hydrodynamic nodes at the interface and deform the volume grids in the fluid computational domain by using the RBF interpolation method. In order to reduce the computational complexity of the mesh motion and improve the efficiency, the greedy algorithm proposed by Rendall and Allen[31]is adopted to reduce the number of control points, with an adequate accuracy of the deformation at the fluid-structure interface, the maximum interpolation error being in the order of 10-3.The assessment of the RBF interpolation accuracy was made in literature[30].

(5) If the static equilibrium is reached, move to the next time-step, else repeat (1)-(4).

1.3 Physical and simulation setup

A NACA66 hydrofoil is used in the present study,with its computational domain as shown in Fig. 1. The hydrofoil has a maximum thickness-to-chord ratio of 12% at the position 45% chord length from the leading edge. The hydrofoil chord length isC= 0.15m , the hydrofoil span length isl= 0.191m,and the angle of attack is 8 degrees. The hydrofoil is fixed in the water tunnel and the test section is 14C(length)× 10C(height)× 2C(width), with the free stream velocityU∞=5 m/s yielding a Reynolds number ofRe=ρU∞C/μ=5.59×107. The outlet pressure is set according to the cavitation number as. Figure 2 shows the three-dimensional fluid mesh. An O-Y type mesh is generated in the computational domain with a sufficient refinement near the hydrofoil surface in order to satisfy the condition that y+=/1yτμν≈,whereyis the distance from the first cell to the hydrofoil surface,τμis the wall frictional velocity.Based on the mesh independence study, the number of the final grid nodes is 3 937 911.

Fig. 1 Computational domain for NACA66

Fig. 2 Fluid mesh grids

Table 1 The parameters of the hydrofoils

The unsteady governing equations are discretized in both space and time. The second order implicit scheme is used for transient formulations. For the spatial discretization of the convection terms, the PRESTO! scheme is used for the pressure equation,and the second order upwind scheme is used for the other convection-diffusion equations. The SIMPLEC algorithm is used for the pressure-velocity coupling solution, to obtain a converged solution quickly with the under-relaxation factor set to 1.0. The numerical simulation of the unsteady cavitating induced vibrations are started from an unsteady cavitating flow field. Subsequently, the in-house UDF code is compiled and turned on in the dynamic mesh panel.The time step is set to 1×10-4s so that the courant number is about 1.

In order to estimate the accuracy of the present fluid-structure interaction technique, we simulate the fully-wetted hydrofoil in the still water with an initial force and analyze the vibration response by monitoring the vibrations at the hydrofoil tip. On the other hand, the ANSYS workbench commercial software is used to perform the modal analysis incorporated with the acoustic method in order to obtain the first-order mode of the hydrofoil in the water. Figure 3 shows the time dependent displacement at the foil tip. The hydrofoil oscillates periodically and its displacement gradually decreases due to the energy dissipation. In Fig. 4, the calculated natural frequency in the still water isby using the present modified FSI technique, and this is the same as predicted by the ANSYS workbench,demonstrating that the modified fluid-structure interaction technique in the present study is reliable.

2. Results and discussions

2.1 Flow-induced vibration characteristics

As the vibration characteristics of the flexible hydrofoil, the dry frequency (in the air) isand the fully-wetted frequency (in the still water) is, as shown in Table 2. It is shown that a 20.4% reduction is resulted in the natural frequency due to the added mass in the still water. In the air, the added mass of the air can be ignored, while in the still water, the added mass due to the water must be considered.

Fig. 3 Predicted vibration displacement in still water by using the modified FSI technique

Fig. 4 Predicted natural frequency in still water by using the modified FSI technique ( fwC/U∞=1.17)

Table 2 Predicted frequencies of flexible hydrofoil in air,still water and cavitating flows

Fig. 5 (Color online) Frequency spectrum of the vapor volume evolution around the flexible hydrofoil and rigid hydrofoil

In the present work, we mainly focus on the flow-induced vibration characteristics of the flexible hydrofoil, so it is necessary to carry out a preliminary analysis for the rigid hydrofoil around the cavitating flows. Figure 5 shows the frequency spectrum of the vapor volume evolution around the flexible and rigid hydrofoils. The cavity shedding frequency for the rigid hydrofoil isfc1C/U∞=0.51, which is the same as that for the flexible hydrofoil. This demonstrates that the structure vibration is not large enough to affect the cavitation evolution, so the cavitation around the flexible hydrofoil keeps the same pace as that around the rigid hydrofoil. Since 4 864 samples are taken every 1×10-4s during the unsteady calculations, the frequency resolution of the vapor volume evolution isf*C/U∞=0.06, which could be improved by increasing either the sampling time or the sampling frequency.

For the flexible hydrofoil, the cavitation unsteadiness would cause the vibrations of the flexible structure, as illustrated in Fig. 6. It is shown that the hydrofoil tip vibrates evidently with small amplitudes in the time history. As the first derivative of the vibration displacement, the vibration velocity represents the vibration intensity and the displacement can be visually observed at the hydrofoil tip. The vibration displacement ranges from 0.74 mm to 0.82 mm, and the average displacement is 0.77 mm with the standard deviation of 1.5×10-5. This is in the same order of magnitude as predicted in literature[6,8,32-33].The vibration velocity fluctuates ranging from-9.0×10-3m/s to 9.1×10-3m/s, and the average value is -1.3×10-4m/s with the standard deviation of 3.4×10-3.

Fig. 6 Time-dependent vibrations of displacement, velocity and acceleration at the flexible hydrofoil tip

To investigate the effect of the cavitation behavior on the flow-induced vibration, the Fast Fourier Transform is carried out for the vibration displacement and the results are analyzed together with the vapor volume oscillations, as shown in Fig. 7.It is found that the dominant flow-induced vibration frequencyfcsC/U∞=0.99 is approximately twice the cavity shedding frequency. It is observed that the other harmonics due to the cavitation oscillations also make contributions to the responses of the flexible hydrofoil.

Fig. 7 (Color online)The displacement and vapor volume oscillations of the flexible hydrofoil in cavitating flows

The cavitation evolution is closely associated with changes in the lift coefficient, which is defined as. Therefore, it is reasonable to analyze the amplitude responses of the lift coefficient along with the cavitation oscillations for the flexible hydrofoil as shown in Fig. 8. The lift coefficientfluctuates periodically ranging from 0.43 to 0.53, and the average value is 0.50 with the standard deviation of 0.010. TheCLamplitude is significantly amplified at the cavity shedding frequency.

2.2 The transient cavitation patterns associated with the flow-induced vibrations

Figure 9 shows the time-history of the vapor volume and the vibration displacement of the flexible hydrofoil within several cycles. Eight instants in one typical cavitation cycle are selected from Fig. 9 to illustrate the transient cavitating flow patterns around the flexible hydrofoil in Fig. 10 by using the iso-surface ofαv=0.1 overlaid with the streamwise velocity. The two-dimensional contours represent the distributions of the vapor volume fraction at plane I,60 mm away from the solid wall as shown in Fig. 11.

Fig. 8 (Color online) The lift (CL ) and vapor volume oscillations of the flexible hydrofoil in cavitating flows

Fig. 9 (Color online) Time-history of the vapor volume and hydrofoil displacement

Fromt1tot3, the attached cavity is generated from the leading edge and grows on the suction side,meanwhile the cloud cavity continuously moves downstream with the decrease of the cloud cavity volume until it collapses att4. It is observed that the cavity rear part begins to oscillate and breaks into small-scale cloud cavities att4. The displacement goes up firstly to a maximum and then keeps decreasing to a minimum att3. During this process,as shown in Fig. 12, a re-entrant flow develops at the rear of the attached cavity due to the strong adverse pressure gradient and moves toward the leading edge.The re-entrant flow breaks through the liquid-vapor interface and makes the attached cavity detached from the hydrofoil, causing the attached cavity to gradually shrink and break the rear part into the cloud cavities fromt4tot7. The displacement goes up fromt3tot5, then gradually decreases to a minimum displacement corresponding to the maximum vapor volume and tends to increase untilt7. Att=t7, the attached cavity shrinks to the leading edge with the presence of several medium-scale cloud cavities. Further in this cycle, the cloud cavities continuously move downstream as illustrated att=t8and begins to collapse where the pressure is higher than the saturated vapor pressure. At the same time, the cavity starts to grow,followed by a new cycle.

As shown by the cavitation patterns in Fig. 10,the cavitation structures are strongly three-dimensional due to the existence of the gap flow and the solid wall (referring to the root end of the hydrofoil).The sheet cavity near the hydrofoil tip is very stable with a triangular shape due to the gap liquid flow.Besides, there is a rectangular cavity att=t1andt6-t8due to the existence of the low-velocity boundary layer at the solid wall of the hydrofoil root.The strong unsteadiness of the sheet/cloud cavities is mainly observed in the middle of the span, and the unsteady cavitating patterns are close to the hydrofoil root, which indicates that the gap flow has a greater effect on the cavitation behavior than the solid wall of the hydrofoil root.

To investigate the effect of the hydrofoil flexibility on the cavitating flows, the evolutions of the cavitation patterns for the rigid and flexible hydrofoils are compared in Fig. 13. It is observed the cavitation evolutions for the rigid and flexible hydrofoils are similar, including the cavity growth,development, shedding and collapse processes. For the flexible hydrofoil, although the vibration displacement is very small, the cavitation features are more complex during the cavity developing and shedding processes. From Figs. 13(b) and 13(c), it is seen that the cavity trailing edge (CTE) of the flexible hydrofoil is longer than that of the rigid hydrofoil. In the shedding process, the attached cavity around the flexible hydrofoil starts to oscillate and break into small-scale cloud cavities at the cavity rear part, as shown in Figs. 13(d)-13(f), while the cavity around the rigid hydrofoil is much more stable. Subsequently,instead of shedding of the large-scale cloud cavities as shown in Figs. 13(g) and 13(h), the cavity for the flexible hydrofoil breaks into medium-scale cloud cavities and moves downstream. It is found that the cavitation shedding process is closely associated with the unsteady hydroelastic response, and the cavity rear part appears to be fragmented.

Fig. 10 The variation of the three-dimensional cavitation structures during a typical cycle depicted by the iso-surface of α v=0.1 overlaid with the streamwise velocity. The two-dimensional contours represent the distributions of the vapor volume fraction at plane I, as shown in Fig. 11

Fig. 11 (Color online) The position of Plane I, 60 mm away from the solid wall

Fig. 12 (Color online) Velocity vectors near the leading edge at t = 25%Tref

2.3 Analysis of unsteady vortex structures obtained with different vortex identification methods

From the above results, the transition of the attached cavity to the cloud cavity is strongly unstable.According to the cavitation-vortex interaction mechanism proposed by Ji et al.[16-17], the cavitation affects the vorticity distribution by the vortex stretching, the vortex dilation and the baroclinic torque terms. Inspired by their work, three vortex identification methods, i.e., the vorticity, theQ-criteria and theΩmethod, are adopted to further investigate the unsteady vortex features around the flexible hydrofoil.

According to the Helmholtz velocity decomposition, the velocity gradient tensor ?ucan be decomposed into the symmetric partAand the anti-symmetric partBas in Eq. (21). Different vortex identification methods are defined based on a general understanding that the symmetric part represents the deformation and the anti-symmetric part represents the rotation. The vorticity is defined as the anti-symmetric tensor of the velocity gradient tensor, i.e., the partBin Eq. (23). In theQ-criterion[19],Qis defined as in Eq. (24), whereaandbare the squares of the Frobenius norm ofAandBgiven in Eq. (22). The part withQ＞0 defines as the vortex region. Liu et al.[22]proposed a new vortex identification method ofΩgiven in Eq.(24), which is a ratio of the vortical vorticity over the total vorticity.Ω=0 represents the pure deformation andΩ=1 indicates the rigidly rotational flow.ξis a small positive number used to avoid division by zero and is recommended to take a value in the level of 10-3[34-35].

Fig. 13 (Color online) Comparison of cavitation features between the rigid and flexible hydrofoils

To have a further insight into the cavitation-vortex interaction, Fig. 14 shows the relation between the cavitation and the vortex evolution. Sincebin Eq.(22) represents the vortical vorticity, it can be treated as the vortex strength in this paper, which is also thus selected by Zhang et al.[36]. It is observed that the variation trend of the cavity is consistent with that of the vortex.

Figure 15 depicts the vapor volume fraction distribution and the contours of the vorticity, theQ-criteria and theΩmethod on the plane I. These results are used to (1) evaluate the capability of different vortex identification methods, (2) analyze the effect of the cavitation evolution on the vortex structures.

Fig. 14 (Color online) Time-dependent vapor volume and vortex strength represented by b

Fort= 17%Tref, the attached cavity develops along the suction surface with the cloud cavity transporting downstream. The cavity interface is represented by the black line with the contour line ofαv=0.1. The vorticity （ωz) on the plane I concentrates around the foil surface with negative values on the suction surface and positive values on the pressure surface. Negative vorticity means the clockwise rotation and positive vorticity means the counter-clockwise rotation. The clockwise vortex is mainly located in the sheet cavity region, and a strip counter-clockwise vortex is observed on the foil pressure surface. A vortex pair is observed in the foil trailing. Both theQ-criteria and theΩmethod identify strong vortices in the sheet/cloud cavitation region and weak vortices in the foil trailing edge region and the foil wake flow region. The differences are (1) theΩmethod identified vortices mainly locate outside the cavity interface and strong vortices are observed in the foil wake flows, while (2) theQ-criteria identifies vortices inside and outside the sheet cavity with weak vortices in the foil wake flows.The large-scale vortex structures are identified by experiments in the cavity sheet rear part and in the wake region[37-38]. Therefore, it is shown that theQ-criteria and theΩmethod can identify vortices in the cavity region very well but identify vortices in the foil wake with different positions and magnitudes.

In the shrinkage stage (t= 50%Tref), a re-entrant jet flow is observed beneath the attached cavity and it causes the cavity oscillate in the rear part and break into small-scale cloud cavities. With the vorticity identification, the vortex structures are shown to be distributed around the foil surface in the clockwise direction in the re-entrant jet region and in the counter-clockwise direction in the cavity region.There is also a positive strip vortex on the pressure surface and a vortex pair near the foil trailing edge.Both theQ-criteria and theΩmethod can identify the vortex structures in the cavitation region and near the foil trailing edge. TheQ-criteria not only can identify cavity-interface-like vortices which are easily distinguished by theΩmethod but also can identify vortices inside the sheet cavity region. In contrast,strong vortices in the foil wake flow are identified by theΩmethod.

Att=Tref, the cloud cavity is shedding, moving downstream and the sheet cavity starts to develop again. The vorticity method identifies clockwise vortices near the suction surface and counter-clockwise vortices with a strip shape on the pressure surface. A vortex pair is identified by thezωmethod at the foil trailing edge. The vortex structures in the cavitation region and near the foil trailing edge are visualized by both theQ-criteria and theΩmethod, but vortices in the wake flow can only be clearly identified by theΩmethod.

From these results, it is seen that the vortices identified by the vorticity method are varied in the transient cavitation evolution due to the strong vortex-cavitation interaction, which is extensively studied by Ji et al.[17]. It is noted that the vorticity shows similar distributions, that is, the clockwise vorticity near the suction surface and counterclockwise vorticity with a strip shape on the pressure surface. The vortex structures identified by theQ-criteria and theΩmethod are distributed in the cavitation region and near the foil trailing edge,especially at the cavity-liquid interface. Due to the downstream transportation of the shedding cavities,the vortices increase and become intensely unstable which can only be identified by theΩmethod.

This is clearly illustrated in Fig. 16, where the three-dimensional vortex structures are visualized by theQ-criteria and theΩmethod with different thresholds att= 50%Tref. It is observed that both methods can identify the vortex ring at the foil trailing edge and the tip-leakage vortex at the foil tip although with different thresholds. However, from Figs. 16(a)and 16(b), the vortex structures in the attached cavity region are more sensitive to the thresholds of theΩmethod, that is, theΩthresholds mainly affect the vortex visualization in the attached cavitation region.In contrast, the thresholds of theQ-criteria have major effects on the vortex identification in the wake flows as shown in Figs. 16(c) and 16(d).

Furthermore, the vortex strength is analyzed quantitatively. Based on the vortex structures identified by the iso-surface ofΩ= 0.52 in Fig. 16(a)and the iso-surface ofQ=200 s-2in Fig. 16(b), we have 121 122 points in theΩvortex and 109 463 points in theQvortex. Figure 17 shows the vortex strength (b) distributions of the two vortex structures att= 50%Tref, where the horizontal axis is scaled by taking twice lg10 of thebvalues in order to clearly display thebdistributions. TheΩmethod can identify vortices in regions I and II while theQ-criteria can only identify vortices in region III.Note that in regions II and III largebvalues mean strong vortices and in region I smallbvalues mean weak vortices. Therefore, it is demonstrated that theΩmethod can identify all scales of vortices including strong and weak vortices and theQ-criteria prefers to identify strong vortices.

Fig. 15 (Color online) Distributions of the vapor volume fraction (αv) and vortex visualization with vorticity, Q-criteria and Ω method at three instants. The black line represents the contour line of αv=0.1

Fig. 16 (Color online) Three-dimensional vortex visualizations using Q-criteria and Ω method at t = 50%Tref

Fig. 17 (Color online) The statistical distribution of the vortex regions identified by the Q = 200 s-2 and Ω =0.52 at t = 50%Tref indicating that Ω method can identify all scales of vortices including strong and weak vortices and Q- criteria can identify strong vortices

3. Conclusions

In this paper, numerical simulations of the unsteady cavitating flows around a flexible NACA66 hydrofoil are carried out by using the modified partitioned FSI approach. The modi- fied coupling strategy includes: (1) the hydrodynamic solution using a LES together with a homogenous cavitation model, (2) the structural deformation solved as for a cantilever beam,(3) fluid-structural interpolation and volume mesh motion based on the radial basis functions and greedy algorithm. The flow-induced vibration characteristics together with the transient cavitation behaviors are investigated. Furthermore, the unsteady vortex features due to the hydroelastic response are analyzed by using three vortex identification methods (vorticity,Q-criteria andΩmethod). Several observations are as follows:

(1) The dominant flow-induced vibration frequency is twice the cavity shedding frequency.Note that the cavity shedding frequency is the same for the rigid and flexible hydrofoils, demonstrating that the structure vibration is not large enough to affect the cavitation evolution.

(2) The predicted cavitating behaviors are strongly three-dimensional, that is, the cavity is (a) in a triangular shape near the hydrofoil tip, (b) in a rectangular shape near the hydrofoil root, and (c) with a strong unsteadiness in the middle of the span,including the attached cavity growth, oscillation and shrinkage, break-off and collapse downstream. The unsteady hydroelastic response will affect the cavitation shedding process with small-scale fragments at the cavity rear part.

(3) For the flexible hydrofoil, the cavity variation trend is consistent with the vortex evolution. The vortex structures are distributed near the foil trailing edge and in the cavitation region, especially at the cavity-liquid interface. With the transporting downstream the shedding cavities, the vortices gradually increase in the wake flows.

(4) Regarding to the vortex identification methods, both theQ-criteria and theΩmethod can identify the vortex structures near the foil trailing edge and in the cavitation region. TheΩmethod can visualize the vortex increase in the wake flows with the transporting downstream the shedding cavities. Besides, theΩmethod can identify all scales of vortices including strong and weak vortices and theQ-criteria prefers to identify strong vortices.

It is noted that the elastic hydrofoil is solved as a cantilever beam in this paper without considering the effect of twist on the cavity behaviors. The twist deformation is also very important for the pressure distribution and the structural displacement, and this can be incorporated into the present FSI technique based on the three-dimensional structural modes.