Wei-wen Zhao, Jian-hua Wang,De-cheng Wan

Computational Marine Hydrodynamics Lab(CMHL),School of Naval Architecture,Ocean and Civil Engineering,State Key Laboratory of Ocean Engineering,Shanghai Jiao Tong University,Shanghai 200240,China

Abstract:In this paper,several commonly used vortex identification methods for marine hydrodynamics are revisited.In order to extract and analyse the vortical structures in marine hydrodynamics,the Q ,2-λ criterion and modified normalized Liutex/Rortex RΩ～method are utilized for vortex identification for propeller open water test,ship drag test,ship propeller-rudder interaction, VIV of a marine riser and VIM of a Spar platform.The limitation of Q and 2-λ criterion is discussed.The Liutex/Rortex RΩ～method is promising for convenient and accurate vortex identification and visualization. However,care should be taken when choosing the small parameter 0b for RΩ～ .We proposed recommended values of 0b for marine hydrodynamic problems.

Key words:Marinehydrodynamics,vortex identification,Omega,Liutex/Rortex,naoe-FOAM-SJTU solver,OpenFOAM

Introduction

Marine hydrodynamic flows are turbulent and disordered,not only because the flow is usually at high Reynolds number,but also due to the complexity of wall-bounded flow especially for those involving complex geometries.These two main factors lead to the ubiquitous vortex structures.For example,the tip and hub vortex of a rotating propeller[1-2],the vortices shed from ship hull when performing large drift angle maneuvering[3-4],the vortices shed from appendages(such asbilge, fin)of a fully appended ship[5-6],vortex shedding of deep-draft column stabilized floaters[7-9],etc..In the past,the research of marine hydrodynamics was focused on macroscopic quantities such as forces,moments and motions due to the limitation of experimental facilities and potential theory methods.However,the advancement of particle imagine velocimetry (PIV)and computational fluid dynamic(CFD)for viscous flow makes it possible to provide detailed flow field information.Researchers are paying more and more attentions on the analysis of vortex structures for hydrodynamics.By extracting and studying the vortex dynamics,researchers can have a deep insight of the detailed flow field such as pressure fluctuation,loads,vibrations and fatigue on structures.

Typically,vortex can be visualized by regions or lines[10].The region method takes vortex as a connected region and uses isosurfaces to identify vortex.The isosurfaces could describe the extended distance from the vortex core.The line method tracks the trajectory of fluid particle through extracting the coresof swirl motion in given area[11].Thesetwo class of visualization methods are complementary to each other for different purposes.Vortex identification methods can be categorized into either Eulerian or Lagrangian.The first type uses Eulerian quantities such as velocity,pressure to help identifying vortex.The second type is usually based on the trajectory of fluid particle motions.

In the current study,we first review several Eulerian vortex identification methods in detail and discuss the relative merits of these methods.Then selective numerical examples of marinehydrodynamic problems are processed by the Q,2-λ criterion and modified normalized Liutex/RortexRΩ～method.The strengths and weaknesses of different methods are discussed based on the extracted and visualized vortex structures.

1.Vortex identification methods

To study the behaviour of vortices,the vortex structures should first be identified and extracted from the flow field.This section outlines the most commonly used vortex identification methods in marine hydrodynamics.

1.1 Vorticity

Vorticity is defined as the curl of velocity.It can be written in thefollowing formula

It is undoubtfully that the magnitude of vorticity is the most widely used quantity to represent vortex cores for free shear flows,especially for 2-D case.However,vorticity magnitude is not necessarily an appropriate vortex identifier for wall-bounded flows because vorticity includes both the near-wall shear motion and swirling motion.The swirling part cannot be extracted from vorticity.For example,in 2-D wall-bounded flow,the maxima and minima of vorticity magnitude occur at the wall,where a vortex obviously does not exist.The large-scale vorticity is due to the strong shear effect of boundary layer in near wall region.In addition,vorticity cannot identify vortex cores in shear flows if the background shear is comparable to the vorticity magnitude[12].Vorticity is a necessary,but not a sufficient condition for vortical motion. Thus vorticity isnot suitable nor accurate to identify vortex.

1.2 Methodsbased on velocity gradient

Due to the shortcoming of vorticity for identifying vortex,several methods based on the local characteristics of velocity gradient are proposed.As the velocity gradient tensor is independent of different Galilean frames of reference,these methods are Galilean invariant.

The characteristic polynomial equation of velocity gradient tensor ?u can be given by the following formula

1.2.4 Liutex/Rortex method

As the skew symmetric tensor Ωcannot represent the rigid rotational motion,the Liutex/ Rortex method[16-20]is proposed to address the issues of how to extract rigid rotational part from the fluid motion.The strict definition of Liutex vector can be summarized as follows:

(1)The direction of Liutex vector r

If the velocity gradient tensor ?u has three real eigenvalues,that means the fluid motion only contains stretch and deforming.In such condition,the magnitude of Liutex vector is zero.On the contrary,if ?u has one real and two conjugate complex eigenvalues,there exists vortex in this local region.The Liutex vector is aligned with the direction of eigenvector corresponding to the only real eigenvalue.The normalized eigenvector ismarked as r.

(2)Themagnitude of Liutex vector R

By transforming the original coordinate system xyz to a new coordinate system XYZ with a transform matrix Q,the original -z axisis rotated to the direction of r .In the new coordinate system,fluid will only rotate in the XY plane,i.e.,the-Z direction derivatives ofXU andYU are zero.After rotation, the velocity gradient tensor becomes

whereciλ is the imaginary part of eigenvectors corresponding to the conjugate complex eigenvalues of ?u.

The author claims that by using the explicit formula,the calculation efficiency for Liutex can be dramatically improved by 37%[21],which is rather impressive.From the formula,it is easy to judge that Liutex vector involves the eigenvalues,eigenvectors of velocity gradient tensor,as well as local vorticity.Furthermore,Liutex vector implies not only the rotation strength but also the rotation axis of the local fluid motion.The physical meaning is rather clear,it decomposes vorticity into the vortical part R and non-vortical part -Rω and only uses the vortical part to represent vortex.

1.2.5-Ω Liutex method

By combining the advantages of Ωand Liutex methods,Dong et al.[19]proposed the normalized Liutex method which defines a scalarRΩwith the following formula

The combined method shows promising results.It can not only measure the relative rotation strength on the plane perpendicular to the local rotation axis,but also separate rotational vortices from shear layers,discontinuity structures and other non-physical structures.However,Liu and Liu[20]found that the extracted isosurfaces formed byRΩis not smooth enough and contain many bulges.They improved it and proposed the modified normalizedRΩ～method,whoseformula isgiven by

1.3 Other methods

1.3.1 Helicity

The helicity is introduced by Levy et al.[22]to extract vortex core lines.The normalized helicity is given by a scalar quantity defined by

1.3.2 Pressure method

For 2-D flow,pressure minimum always occurs at the rotation centre of the fluid.It is therefore reasonable to extract rotation of fluid by locating the minima of pressure.By tracing the sectionally local minimum of pressure, the disconnected vortex filament can be extracted[23]. However,this method is not suitable for unsteady,viscous and 3-D flow.

2.Numerical approaches

In this section,we choose several typical flows in marine hydrodynamics to demonstrate the characteristics of several vortex identification methods.All the selected examples are numerically computed using the naoe-FOAM-SJTU solver[24-26],which is developed on top of the OpenFOAM framework.It utilizes the data structures and low-level infrastructures of OpenFOAM and consists of several specialized features for marine hydrodynamics. Figure 1 shows the functional modulesof the solver,in which theblue where the blue frame and red frame represent developed and developing modules,respectively.

The calculation of Q and2λare provided in OpenFOAM’s function object framework as postprocessing utilities.Similarly,we implemented Ω ,Liutex and -Ω Liutex methods under the frame work of function object as postprocessing utilities.In such a way these vortex identification methods can be applied to arbitrary unstructured polyhedral computational mesh with complex geometries.

Fig.1 (Color online) Main framework of naoe-FOAM-SJTU solver[27]

The comparisons between Ω ,Liutex and -Ω Liutex have been discussed extensively.In this paper,we compare Q and2λmethods with the most recently proposedRΩ～method.Recalling that the small parameter0b inRΩ～method is defined to prohibit the computational noise and avoid dividing by zero[20].The default value of0b is recommended to set as a value between 0.001 and 0.002.However,this value is not suitable for typical marine hydrodynamic flows.We will discuss the applicable values of0b for marine hydrodynamic flows in Section 4.

3.Applicationsof vortex identification methods

3.1 Propeller open water

Fig.2 (Color online)Propeller operating at =0.9J .Isosurfaces of different 2-λ criterion,coloured by vorticity magnitude

Fig.3 (Color online)Propeller operating at =0.9J .Isosurfacesof different RΩ～ ,coloured by vorticity magnitude

Figure 2 shows isosurfaces of different2λ values coloured with vorticity magnitude.The tip vorticesof the propeller are resolved clearly.However,using different Q values as the threshold for contour,various length of tip vortices are extracted.For small 2λvalues shown in Fig.2(a),more detailed wake structures are identified and extracted.As threshold of Q increases,those weak vortex (fluid regions with smaller vorticity magnitude)cannot be identified.The structures of vortex filament are also distinct.The vortex tube extracted by larger Q threshold value is thinner.Figure 3 shows the propeller wake structures which are identified and extracted by different values of modified normalized RortexRΩ～and coloured by vorticity magnitude.The contour threshold value 0.52 is the recommended value in the original paper[20].Besides this,we also choose three other threshold values for comparison.As depicted in the pictures,although the vortex lengths in wake region are different,the vortex tube structures extracted by differentRΩ～values and the corresponding vorticities remain keep almost the same.

3.2 Ship drag test

Figures 4 and 5 illustrate the vortical structures extracted by Q-criterion and modified normalized Rortex method,respectively.The vortex above free surfaces are eliminated for sake of simplicity.At first glance,the two figures shows similar vortex structures.In the vicinity of stern,massively separation flow due to the sharp curvature change of the hull surface is observed.The captured structures is,however,very distinctive.In Fig.4,the vortices identified by Q contain redundant shear motions near the hull surface.These deformations is excluded in the vortices extracted byRΩ～ .

Fig.4 (Color online) Vortical structuresin wakeregion presented by isosurfaceswith =5Q and coloured by vorticity magnitude

Fig.5 (Color online) Vortical structuresin wakeregion presented by isosurfaces with Ω～R =0.52 and coloured by vorticity magnitude

3.3 Ship propeller-rudder interaction

The propeller-rudder interaction of the ONRT surface combatant during a zigzag maneuver is visualized.Detailed results of the maneuvering study will be published on the SIMMAN2020 workshop.The simulation follows the SIMMAN2020 case 5.2.2,corresponding to a 20/20 zigzag test,starting to starboard side.The model is scaled at a ratio of 1:48.935 and is self-propelled with free running,which meansall six degrees of freedomsfor the model are released,at model point in calm water.The corresponding Froude and Reynolds number are Fr =0.200 and Re=3.39 ×106,respectively.In the 20/20 zigzag maneuver test the rudder is turned to the opposite direction at a rate of 35°/s to a rudder angle of 20°every time the ship heading reaches ±20°.The computational grid for moving components such as hull,propellers and rudders are generated separately and are then assembled into a composite overset grid.All component grids can translate and rotate freely with respect to other grids.

Figure 6 gives the instantaneous flow visualizations of hull-propeller-rudder interaction for ONRT during zigzag maneuvering.The vortical structuresare identified by different criterions and coloured by vorticity magnitude.Similar to the results in Section 3.1,the dominant wake structures are the clockwise and anti-clockwise tip vortices generated by the contra-rotating twin-screw propellers.Strong hub vortices are also observed.Owing to the numerical dissipation and coarse grid resolution in the rudder downstream wake regions,the vorticity of hub vortices become smaller.These weak vortices is not identified by the Q-criterion.However,theRΩ～method can capture both the weak and strong tip vorticessimultaneously

Fig.6 (Color online)Hull-propeller-rudder interaction of the fully appended ONRTduring 20/20 zigzag maneuver test.Isosurfaces by (a), (b) and coloured by vorticity magnitude

Figure 7 depicts another view of the instantaneous vortical structures around the fully appended ONRT model during zigzag maneuver test.Although the simulations are performed using the two-phase VOF model,only vortices in the water are displayed for the sake of simplicity.The rudders have been executed and turned to the port side,which provide lift rudder force.A consequence turning moment turns the ship towards left.The “fake”vortices owing to the boundary layer shear motions near the hull surface such as bow,bilge keel and skeg vortices are removed from the isosurfaces of Ω～R=0.52.Moreover,isosurfaces of Q=100 can only capture propeller and rudder vortices,with some redundant non-physical vortices near hull,since these are strong vortices and have large vorticity.Other weak vortices generated by free surface ship waves,bilge keel and skeg are failed to be captured by isosurfaces of Q=100,but can be successfully observed in the isosurfaces of Ω～R=0.52.

Fig.7 (Color online) Vortical structures around the fully appen ded ONRT during 20/20 zigzag maneuver test.Isosurfaces by (a),(b)and coloured by vorticity magnitude

3.4 VIV of a flexible riser

Vortex induced vibrations (VIV)is the main cause of marine risers’ fatigue damage.It is very hard to predict VIV of marine risers due to the nonlinearity and instability of the flow.By taking the flexible structural deforming into account,the fluid-structure interaction problem become more complex.In this subsection,a riser under stepped flow is simulated.The diameter and length of the riser is 0.028 m and 13.12 m,respectively. Lower part (45%of total length)of the riser is submerged in water with speed of 0.16 m/s.The top is tensioning with a pretension force of 405 N.This case setup is in accordance with the experiment which is performed by Chaplin et al.[31-32].The VIV of marine riser is simulated by 2.5-D strip method,i.e.,the fluid is solved in a finite set of 2-D strip sections and the structure is obtained with 3-D model.In the current solver’s implementation,the flow fields in 2-D domains are solved using finite volume method provided by OpenFOAM,and the dynamic structural response is solved using finite element method with Euler-Bernoulli beam model.The hydrodynamic force is obtained from fluid and passed to structure and the displacements are calculated in structure and passed to moving boundaries of fluid mesh in each section.

Figure 8 depicts the 2-D vortices of VIV for a flexible marine riser in fluid sections,represented and coloured by different methods.In both pictures,clear Kármán vortex streets are observed in each section.As been elaborated in previous section,the maxima of vorticity magnitude occurs at the near wall region (Fig.8(a)).It is convinced that the disadvantages of vorticity magnitude in representing vortices are obvious.Different vortices in the vicinity of wall are identified as a single connected region which is non-physical. Furthermore,the vorticity magnitude is so small in the far wake region that it is hard to distinguish vortical structures.On the other hand,RΩ～ ,unlike vorticity magnitude,can clearly resolve the individual vortex cores in both near wake region and far wake region.

Fig.8 (Color online)Vortices of VIV for a flexible marine riser in 2-D strip sections.The vorticesare represented and coloured by (a),(b)

3.5 VIM of a spar platform

Fig.9 (Color online) Vortical structures of VIM for a Spar with helical strakes. Isosurfaces by (a), (b)and coloured by vorticity magnitude

4.Discussion on the sensitivity of 0b

Fig.10 (Color online)Hull-propeller-rudder interaction of the self-propelled ONRT at model point. Isosurfaces by Ω～R =0.52

5.Conclusions

In this paper,we first revisit the most commonly used vortex identification methods for marine hydrodynamic flow problems.After that,we introduce the in-house solver naoe-FOAM-SJTU and several numerical examples which are computed by the solver.These classical examples are selected to demonstrate the advantages and disadvantages of different vortex identification methods.After the discussions,the following conclusionscan be drawn.

(1)Vortex identification methods such as Q and2-λ criterion is quite sensitive to the threshold value of isosurface.Hence,it is hard to extract the correct vortical structures for further vortex dynamics analysis.

(2)The-Ω Liutex (modified normalized Rortex)method shows promising results in vortex identification and visualization for marine hydrodynamics.It removes vorticity due to shear motion and only keeps the pure rotational motion part.Furthermore,it can identify strong and weak vortices simultaneously and is insensitive to isosurface threshold value.

(3)For marine hydrodynamic flows which involve strong absolute vortices strength,0b should be set to a small value.We recommend that0b is set to 10-6.

Acknowledgements

This work was supported by the Chang Jiang Scholars Program (Grant No.T2014099),the Shanghai Excellent Academic Leaders Program (Grant No.17XD1402300) and the Innovative Special Project of Numerical Tank of Ministry of Industry and Information Technology of China (2016-23/09).