CFD prediction of unsteady forces on marine propellers caused by the wake nonuniformity and nonstationarity

Ocean Engineering 104 (2015) 659–672

Contents lists available at ScienceDirect

Ocean Engineering journal homepage: www.elsevier.com/locate/oceaneng

CFD prediction of unsteady forces on marine propellers caused by the wake nonuniformity and nonstationarity N. Abbas a,1, N. Kornev a, I. Shevchuk a, P. Anschau b,2 a b

Chair of Modeling and Simulation, University of Rostock, 18057 Rostock, Germany Potsdam Model Basin, Marquardter Chaussee 100, 14469 Potsdam, Germany

art ic l e i nf o

a b s t r a c t

Article history: Received 29 March 2014 Accepted 6 June 2015 Available online 24 June 2015

Results of computations of the unsteady loadings on marine propellers behind the KVLCC2 tanker using different engineering and numerical methods are presented and analyzed. The hybrid URANS–LES model presented in Kornev et al. (2011) is applied to the arrangement containing both the ship and the rotating propeller with consideration of all interaction effects. Comparison of numerical results with these obtained from engineering methods allows one to validate the empirical estimations based on strong simplifications. Two of three considered engineering methods give standard deviations of the force and moment fluctuations with the accuracy sufficient for practical purposes. However, the hybrid model predicts the existence of strong peak loading (three times larger than the standard deviation for the thrust) that are not detected by any engineering methods. Hybrid computations also show the existence of strong transversal structures in addition to common bilge vortices in the ship stern area in the case of a bare hull. The presence of the rotating propeller destroys these additional structures due to suction effect on the ship boundary layer. & 2015 Elsevier Ltd. All rights reserved.

Keywords: Ship wake Propeller URANS LES Hybrid methods

1. Introduction Ship structure vibrations are caused by a complex combination of different hydrodynamical and mechanical effects. The main source of the hydrodynamically excited vibrations is propeller. Determination of loads on marine propeller is one of the important and challenging problems for the prediction of hull structure and propulsion shafting vibration. The main source of the propeller excited vibrations is the unsteady loading, which results from the rotation of the propeller through nonstationary and nonuniform wake. The nonuniformity of the wake is the dominant effect for the vast majority of ships. Classical engineering methods of the marine propeller forces calculation are based on this effect, approximately assuming the velocity field to be stationary, the so-called “frozen field” (see, for instance, Carlton (2007), Breslin and Andersen (1994), Alte and Baur (1986) and Voitkunski (1985)). These methods allow one to predict the propeller forces fluctuations with dominating blade frequencies proportional to nZ, where Z is the number of blades and n is the frequency of the propeller.

E-mail address: [email protected] (N. Abbas). Tel.: þ49 3814989550. 2 Tel.: þ49 3315671262. 1

http://dx.doi.org/10.1016/j.oceaneng.2015.06.007 0029-8018/& 2015 Elsevier Ltd. All rights reserved.

If the wake is steady the variations of forces and moments are strictly periodical. In mechanics, such processes are referred to as periodic rather than unsteady. In this paper the unsteady nonperiodic loadings are considered that correspond to non-periodic flow within the wake, which is the usual case for turbulent flows. The nonstationarity of the wake plays more important role for ships with large block coefficients (the full-bottomed ships) and is far less studied. The wake of full-bottomed ships contains complicated vortex structures which amplify the unsteady effects in wakes. The existing engineering methods do not take the nonstationarity of the wake into account due to the difficulties connected with its determination. From the measurement point of view it is a big challenge to measure the unsteady velocities in the wake using traditional techniques like the Pitot tube or the modern non-intrusive methods like LDV or PIV. Two latter measurement techniques have certain advantages such as the non-intrusive character of measurements and high accuracy at high particle seeding densities. Although these methods are already widely used in ship hydromechanics, their application for the study of propellers behind the ship (propulsion measurements) is very complicated. First, the whole laser setup should be towed together with the model in a special well-streamlined body at a big distance to prevent the wake perturbations. Second, it is difficult to attain the high particle seeding density, which is necessary for proper flow resolution, under towing tank

660

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

experimental conditions. High costs of measurements is the next disadvantage of both laser diagnostics techniques. From the computational point of view it is difficult to resolve vortex structures responsible for the velocity field fluctuations. The URANS (Unsteady Reynolds Averaged Navies Stokes) method, which is widely used in shipbuilding community see, for instance, Xiao-fei et al. (2009), Park et al. (2013) and Hayati et al. (2013), is not capable of modeling unsteady vortices arising in the ship stern area flow. The modern numerical methods like LES (Large Eddy Simulations) require large computational resources for accurate prediction of flows with large Reynolds numbers which are typical for shipbuilding applications (Jang and Mahesh, 2010). The grid resolution necessary for a pure LES is so huge that it makes the direct application of LES impossible. A practical solution of this problem is the use of a hybrid URANS–LES approach, in which the near body flow region is treated using URANS and far flow regions are treated with LES (see, for instance, Xing et al. (2012) and Stern et al. (2013)). The paper is organized as follows. In the second section we briefly review some methods used in this paper for computations of unsteady propeller loadings and also outline our hybrid method which has been developed in our previous paper (Kornev et al., 2011). The third section is devoted to the validation of the hybrid method using diffuser benchmark test and two propeller test cases. In the fourth section, the geometry of ship configuration and the numerical environment of computations are described. The results of the paper are analyzed in the fifth and sixth sections starting with the case of bare hull. The main results of the paper are presented in the last section, devoted to the unsteady loadings on propellers determined using different methodologies. Finally, we summarize the obtained results in the conclusion.

2. Description of methods for prediction of unsteady loadings on marine propellers In this section we give a short overview of some methods utilized in this paper for computation of forces and moments on a propeller. 2.1. Engineering methods The engineering methods described in this section use different kinds of simplifications and approximations of the available experimental data to estimate the desired quantities. Scheme B is a very simple and at the same time a very efficient method proposed by experts of the Krylov Research Shipbuilding Institute (Voitkunski, 1985). The method is based on the hypothesis of quasi-stationarity and assumes that the thrust and the moment coefficients are known (KT, KQ). These coefficients are determined under the open water conditions as functions of the propeller advance ratio. The force and the moment are then distributed along the blade span using correlations taken from the lifting surface theory. From these distributions it is possible to estimate the local lift and drag coefficients corresponding to each blade profile at a certain radius along the span. At each time instant, the force and the moment arising on this profile are calculated using these coefficients, the local values of the incident velocity and the angle of attack. The total force and the moment are determined by integration of local loadings along the blade span. The nominal wake velocity field including tangential and axial components is represented via the Fourier series as a function of the rotation angle. This allows one to express forces and moments in terms of the Fourier coefficients. The influence of the propeller on the nominal wake is taken into account by the

contraction of the propeller jet determined using the formula of the actuator disc theory. The Veritec approximation is briefly discussed in the book of Carlton (2007). It is based on the approximation of results of theoretical investigation of the dynamic forces at blade and twice blade frequencies performed for 20 typical ships (Noise and vibration Group, 1985). The results are summarized in the form of dependence of force and moment Fourier coefficients on the blade number (see Table 11.13 in Carlton, 2007). The Wereldsma (1964) approximation is based on the statistical regression analysis of measurements done for 40 ship models with four and five bladed propellers. To get detailed information the reader is referred to the original publication of Wereldsma (1964). Results obtained using engineering methods and presented in the paper are taken from Batrak et al. (2012). 2.2. Hybrid method The hybrid CFD model developed in our previous work (Kornev et al., 2011) is based on the fact that the governing transport equations have the same form in LES and URANS ∂u i ∂ðu i u j Þ ∂p n ∂ðτij þ τij Þ þ ¼ þ ; ∂t ∂xi ∂xj ∂xj l

t

ð1Þ

but the interpretation of the overline symbol is different. In LES it means filtering, and in URANS it stands for the Reynolds averaging (the term “ensemble averaging” is also used in this context). Here we use the standard notation of pn for the pseudo-pressure, and τlij and τtij for the laminar and turbulent stresses respectively. Note that the turbulent stresses are calculated in different ways in LES and URANS regions. The computational domain in our model is dynamically (i.e. at each time step) subdivided into the LES and URANS regions. The key quantities of this decomposition are the integral length scale L and the extended LES filter Δ which are computed for each cell of the mesh. L is determined from the formula of Kolmogorov and Prandtl: L¼C k

3=2

=ε

ð2Þ

where k is the turbulent kinetic energy and ε is the dissipation rate and C is the certain empiric constant. The accuracy of the prediction (2) which, strictly speaking, is valid for high local Reynolds numbers far from the wall, will be proven below in qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 Section 2.2.1. The filter Δ is computed as Δ ¼ 0:5ðdmax þ δ2 Þ, where dmax is the maximal length of the cell edges pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi dmax ¼ maxðdx ; dy ; dz Þ and δ ¼ 3 ðthe cell volumeÞ is the common filter width used in LES. A cell of the mesh belongs to one area or the other depending on the value of L relative to Δ, if L 4 Δ then the cell is in LES area, in other case it is in URANS region. The extended LES filter Δ depends only on the geometry of the mesh and is computed only once, whereas the integral length scale L is changed in space and varies from one time step to another, which results in dynamic decomposition of the computational domain into the LES and URANS regions. The turbulent stresses τtij are calculated from the Boussinesq approximation using the concept of the turbulent viscosity which is considered as the subgrid viscosity. These stresses are computed according to the localized dynamic model of Smagorinsky in the LES region and according to the k   ω SST turbulence model of Menter (1993) in the URANS region. The turbulent kinematic viscosity is smoothed between the LES and URANS regions: νðxÞ ¼ ανt þ ð1  αÞνSGS

ð3Þ

  1  40x x2 þ x1 1 þ þ 10 αðxÞ ¼ arctan π x2 x1 2 x2  x1

ð4Þ

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

where νt is the RANS turbulent kinematic viscosity and νSGS is the LES subgrid viscosity. The factor 40 in the arc tangent function is chosen such that ν  νSGS , when L=Δ 4 1:05 (LES region) and ν  νt , when L=Δ o 0:95 (RANS region) and for 0:95 o L=Δ o1:05 this expression gives a smooth transition of ν between the two regions. The wall functions are used in the near wall URANS region. Comparing this method with Detached Eddy Simulation (DES) one can state that the main difference is way of determination of the interface between URANS and LES. In both methods the interface is not prescribed and determined dynamically. An important parameter in DES is the distance from the wall d. It is shown that if d gets large the URANS models implemented in DES are transformed smoothly into the Smagorinsky model. Between URANS and LES, there is a grey zone in which LES and URANS models are mixed. Within our method the transition between LES and URANS is performed by comparison of integral length L (2) with the grid resolution. There is also a grey zone between URANS and LES determined by the smoothing formula (3). The determination of the interface is the most critical point of our technique. It is based on the integral length L which is calculated using k and ε from URANS simulation assuming that URANS predicts these quantities properly in the whole domain using instantaneous velocities. Strictly speaking, the governing URANS equations are derived for time or ensemble averaged velocities, but not filtered ones like in LES. However, many of the existing LES models are based on the URANS ones with formal replacement of averaged quantities by filtered ones. In our method we suppose that the URANS models written for filtered velocities are valid for calculation of k and ε also in the LES model at least in the area of strong vortex formation and flow fluctuations in the stern area. A possible mistake outside of this area is not critical for the prediction of phenomena which we are interested in. An alternative solution could be a pure preliminary RANS simulation in the whole domain which would then be used for the LES–URANS interface determination and further hybrid simulation with the fixed interface. Unfortunately, this approach did not yield satisfactory results in the area of flow separation and vortex formation, because it could not account for the large-scale flow unsteadiness. The CFD calculations using both URANS and hybrid models were carried out with solvers from OpenFOAM toolkit (Weller et al., 1998).

approximation obtained by Tomas et al. (2009): y  0:3  0:03 L δ ¼ y 2 δ  0:03 0:31 þ δ

661

ð6Þ

where δ is the boundary layer thickness and y is the distance from the wall. The integral length in formula (6) was calculated as the average of four integral lengths obtained from autocorrelation functions for the axial and vertical velocities u02 and w02 as well for the cross correlations of these velocities u0 w0 and w0 u0 . According to Schlichting and Gersten (2000) C is equal to 0.168. This estimation is valid not within the whole boundary layer rather than in close proximity to the wall at 0 oy=δ o 0:2. Analysis of experimental data from Tomas et al. (2009) shows that the constant C in (2) is about 0.35 in the largest part of the boundary layer at 0:2 o y=δ o 0:8. Fig. 2 shows the comparison of the numerical results obtained using the constant C¼0.35 in (5) with experimental data (6). This confirms that within the boundary layer at 0:2 o y=δ o 0:8 the formula (2) can be used and the value of constant C is around 0.35. Within the region of 0:8 o y=δ o1 C is not fixed and changes between 0.35 and 1.0. A proper agreement between measurement and estimation (2) is attained with C¼ 0.5. At the outer border of the boundary layer at y=δ 4 1:3 the constant C can be taken approximately as C  1:0. Thus, the formula of Prandtl Kolmogorov (2) can be used within the boundary layer if the coefficient C is properly chosen depending on y=δ. The integral length scale prediction (2) is in a good agreement with measurements in the outer part of the TBL if the constant C is chosen from the conditions C  1:0 at y=δ 4 1:3 (see. Fig. 3). Since the switch between LES and URANS happens in the outer region of the TBL, it is assumed that the formation of the turbulent boundary layer along the ship and generation of unsteady bilge vortices can be predicted accurately by hybrid method with the constant C¼1.0 in formula (2). The primary task of the present work is the calculation of a propeller located in a ship wake. Thomas et al. (1991) have shown that the value C  1 is also valid for wake flows. Thus, the classic formula (2) with the coefficient C ¼1.0 can be accepted as quite robust one for estimation of the integral length both along the ship boundary layer and in the wake. This value is then used in all calculations presented in this paper.

3. Validation 2.2.1. Validation of the Prandtl Kolmogorov formula (2) The expression (2) is the key formula to determine the interface between URANS and LES regions in the present method. To check the accuracy of the prediction (2), we performed a series of calculations of zero gradient turbulent boundary layer (TBL) on a plate using 1420  40  40 cells in a box with sizes of 4  0:1  0:1 m, whereas the mean value of y þ at the center of the wall-adjacent cells was equal to 3.18. Calculations were conducted using simpleFoam solver with k   ω-SST-model. The maximum Re number at the end of the computational domain was Rex ¼ 15  105 for the velocity of incoming flow 3:75 m s  1 . The integral length L is calculated from k   ω-SST-model using following estimations: 3=2

L¼C

k

ε

1=2

;

ε ¼ βn ωk ) L ¼ C

k βn ω

ð5Þ

Accuracy of k and ω prediction in OpenFOAM was proven by comparison with NASA benchmark (NASA, 2003). The results agree well with distributions obtained by NASA, (see Fig. 1). Further methodical investigations were focused on the integral length computations at different distances from the plate leading edge. Numerical results for L were compared with the experimental

3.1. Asymmetric plane diffuser flow Asymmetric plane diffuser is a very popular benchmark test in computational fluid dynamics. There are a lot of measurement data obtained using single component laser Doppler velocimeter (Obi et al., 1993), hot wire and pulsed wire anemometer (Buice and Eaton, 1996 and Buice, 1997), which are available to validate the numerical results. In this study we compare the results obtained with the hybrid method at different grid resolutions with experiments of Buice and Eaton (1996). The geometry of diffuser is presented in Fig. 4. The x-axis is in the streamwise direction, y-axis is in the spanwise direction and the z-axis is the vertical one. The origin of the coordinate system is located at the intersection of the straight and inclined walls as seen in Fig. 5. The computational domain extends downstream to 75H, where H is the width of the inlet channel. In spanwise direction the computational domain had the size of 6.67H. On the left and right boundaries cyclic boundary condition was imposed. The expansion ratio of the diffuser is H out =H ¼ 4:7 with the opening angle of 101. The inlet plane is located at x=H ¼  6 to exclude the upstream influence of the expansion. The inlet channel flow is turbulent and fully

662

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

Fig. 1. Validation of k–ω-SST model implemented in OpenFOAM using NASA benchmark test (NASA, 2003).

Fig. 2. Comparison of turbulent length scales L obtained from the Kolmogorov Prandtl formula (2) and Tomas's empirical estimation (6) for the inner part of the boundary layer.

Fig. 3. Comparison of turbulent length scales L obtained from the Kolmogorov Prandtl formula (2) and Tomas's empirical estimation (6) for the outer part of the boundary layer 1:0 o y=δ o 2:0

Fig. 4. Geometry of diffuser (Buice and Eaton, 1996).

developed one with Reynolds number of Re¼20,000 based on the centreline velocity and the channel height. Three different grids with 5:0  105 , 1:2  106 and 6:5  106 cells are investigated. The coarsest mesh with N x  N y  N z ¼ ð394  30  40Þ cells is shown

in Fig. 6. Other two grids have the following cell distribution in x; y and z directions: ð505  40  60Þ and ð590  100  110Þ. Since the mean value of y þ at the first grid point from the wall is smaller than 1.0 for all three meshes no wall functions were applied. Fig. 7 shows snapshots of a flow visualization obtained on the finest grid using DMM (Dynamic Mixed Model, pure LES) method and the present hybrid model. Mean axial velocity U, root mean squares of the axial u02 and vertical w02 velocities at different sections along the diffuser are presented in Figs. 5 and 8. LES DMM results were obtained only for the finest grid. The initial turbulence is generated using the precursor method, i.e. the flow at x=H ¼ 0 is mapped back to x=H ¼ 6 at each time instant. The inlet profile of the mean axial velocity at x=H ¼  6 is well reproduced for all three grids. At x=H ¼  6 there is a certain discrepancy for fluctuations caused by insufficient resolution of the narrow channel flow. This discrepancy is then traced in the first cross

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

663

Fig. 5. Axial mean velocity profiles x=H þ 10  ux =U b in comparison with the experimental data, Ub stands for bulk velocity (Buice and Eaton, 1996).

Fig. 6. Computational domain of the plane diffuser (x–z view of the coarse mesh).

The open-water diagram for the test propeller VP1356 is presented in Fig. 9 and shows a very good agreement with the experimental data provided by SVA Potsdam. The accuracy of the calculations, conducted using pimpleDyMFoam solver is illustrated in Table 1. A very good accuracy is documented for heavily loaded operating conditions at small advance ratios J. The relative error calculated as ðΦexp  ΦCFD Þ=Φexp (where Φ is the analyzed quantity: K T ; K Q or η ) is getting larger when the advance ratio J increases because the reference values Φexp decrease. At J¼0.8 the maximum relative error is about 8% for KT and 10% for KQ whereas it does not exceed three percent for the efficiency η.

3.3. Validation of velocity distribution for PPTC

Fig. 7. Isosurfaces of instantaneous axial velocity component obtained using hybrid method and LES DMM.

section at x=H ¼ 5:0. Further downstream the agreement between measurements and simulations becomes satisfactory. From x=H ¼ 30 experimental and numerical results agree well. Overall, one can conclude that the hybrid model results show clear trend towards experimental data when the resolution increases. The best agreement is attained with the finest grid containing 6.5M of cells. At moderate resolution used in this work, results of the hybrid model are slightly better than these obtained from LES DMM model because LES usually requires much higher resolution especially in near wall region.

3.2. Validation of the hybrid model for a open water propeller test case The thorough validation of the model for the bare hull has been done in Kornev et al. (2011) for the resistance, the mean velocity and the turbulent kinetic energy fields in the propeller plane. In this paper we present new validation results with rotating propeller geometries. For the interpolation of the solution at the interface between static (hull) and rotating (propeller) meshes Arbitrary Mesh Interface (AMI) was used. This is a relatively new option introduced in OpenFOAM v.2.1.0 and improved in 2.3.0 version, which utilizes the concept of conservative mesh to mesh interpolation by local Galerkin projection. It enables simulation across adjacent, separate static and rotating mesh domains through a special boundary condition. For further details the reader is referred to Farrell and Maddison (2011).

Since the data for velocity field of the propeller VP1356 studied above are not available this validation is done for the Potsdam Propeller Test Case (PPTC) proposed by Potsdam Ship Model Basin (SVA) for validating the CFD codes applied to ship propulsion. It was originally prepared for SMP'11 workshop and contains a wide range of experimental data for model propeller VP1304: open water tests, cavitation tests and LDV measurements. Propeller VP1304 is a 5-bladed controllable pitch propeller. Its basic characteristics are listed in Table 2. Since the prediction of integral parameters had already been conducted for VP1356, the aim of this test was to analyze the prediction accuracy of hybrid URANS–LES methods for distributed parameters by comparing with LDV data. For the detailed description of the experimental setup reader is referred to the report 3754 of SVA Mach (2011). The flow around the propeller was studied at J¼1.25, with n¼23 s  1 and U 1 ¼7.20 m s 1. Mesh for the simulation consisted of rotating and static regions. Static mesh region was generated using Ansys ICEM generator, whereas the rotating one was produced using snappyHexMesh generator. The total number of computational cells was 2M. Numerical simulations were performed with the hybrid method and IDDES approach implemented in OpenFOAM. The former yielded KT of 0.242, the latter – 0.233, whereas the experimental value is 0.25. Fig. 10 shows the circumferential distributions of axial, tangential and radial wake fractions at the cross section x=D ¼ 0:094 and radius r=R ¼ 0:7. As one can conclude from the presented plots, the results agree well with experimental data both qualitatively and quantitatively for the axial and radial wake fractions, where the discrepancy is no higher than 10%. The exceptions are the peaks, which are reproduced qualitatively, but the discrepancy reaches considerably higher values. On the contrary, the peaks of wake fraction in tangential direction are captured much better, whereas between them relative discrepancy in some points reaches 50%, which is, however, less than 0.05 in absolute values and can be considered as acceptable. Keeping in mind, that the mesh was quite coarse, the results were considered as satisfactory.

664

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

Fig. 8. Mean velocity fluctuation profiles: (a) x=H þ 20  u0x =U b and (b) x=H þ 50  u0z =U b in comparison with the experimental data (Buice and Eaton, 1996).

Fig. 9. Open-water diagram of the test propeller calculated using different numerical methods.

Table 1 Deviations in percent between measurement and OpenFOAM calculations for the propeller VP1356. J

KT

KQ

η

0.3 0.4 0.5 0.6 0.7 0.8

2.79 2.29 3.54 4.23 5.89 8.36

3.56 4.22 5.19 6.83 7.22 10.24

 1.27  1.72  2.50  2.34  2.25  1.08

Table 2 VP1304 parameters. Diameter Pitch ratio, r/R¼ 0.7 Area ratio Chord length, r/R ¼0.7 Skew Hub ratio Number of blades Sense of rotation

(m) (–) (–) (m) (deg) (–) (–) (–)

0.25 1.64 0.78 0.10 18.837 0.30 5 Right

Fig. 10. Circumferential distributions of mean axial, tangential and radial wake fractions at the cross section x=D ¼ 0:094 and radius r=R ¼ 0:7. SVA propeller test case VP1304.

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

4. Geometry and the numerical environment The doubled model of the KRISO tanker KVLCC2 with the scale 1/58 has been chosen as a test object, since it is widely used in the shipbuilding community for CFD validations (Larsson et al., 2000). The model has length of 5.517 m, breadth of 1 m, draught of 0.359 m and block coefficient of 0.8098. The study of the model has been performed for the constant velocity of u0 ¼ 1:047 m s  1 corresponding to the Reynolds number of Re ¼ 5  106 . The Froude number Fn ¼ 0:142 is small which makes it possible to neglect the water surface deformation effects. The propeller VP1356 used in our test has five blades (Z¼5), the diameter D of 0.17 m, pitch ratio of 0.996, EAR (expanded area ratio) of 0.8, Dhub =D ¼ 0:18 and maximal skewness of 321. The propeller frequency n is 8 rotations per second. The bare hull was studied with four gradually refined grids generated by Ansys ICEM and containing 7M, 13M, 25M and 45M of cells. The first computational point of the grid with 7M of cells was located at y þ ¼ 6 from the wall with a sufficient refinement in the propeller disc (see. Fig. 11). All other grids were obtained from this grid by gradual refinement. For that the cell numbers in all three directions ðx; y; zÞ in each block in the boundary layer and the propeller disc (see Fig. 11(b)) have been increased by the factor of 1.3. The first cell center for 13M, 25M and 45M was, respectively, at y þ ¼ 2:6; 1:6 and 0:5 . The system ship with rotating propeller was studied with coarse (13M) and fine (22.5M) grids (see Fig. 12). The static grids have been generated using the Ansys ICEM software. The coarse static grid contains 11M whereas the fine one 19.5M. The rotating grids have been generated using the snappyHexMesh-mesher provided by OpenFOAM and contain another 2M cells for the coarse grid and 3M for the fine one. The coarse grid of 13M has y þ  0:1–4 in the wall region of aftership and y þ  16 in the foreship area. The computations have been carried out with the fixed maximal Courant number of 3.5, which corresponds to the time step  1:0  10  4 s. The fine grid of 22.5M has y þ  0:1–3:5 in the wall region of aftership and y þ  3 in the foreship area. The computations have been carried out with the fixed maximal Courant number of 4.0, which corresponds to the time step  1:0  10  4 s. For both grids the value of the rotation angle of the propeller for each time step is equal to 0:2881. The spatial discretization of the convective term is performed using the filtered Linear scheme implemented in OpenFOAM. This scheme calculates the face values using blending of linear interpolation with a particular amount of upwind, depending on the ratio of the background (in-cell) gradient and face gradient. The amount of upwind is limited to 20%. This way the high-frequency oscillation modes are filtered and thus the stable solution is obtained without a considerable increase of numerical dissipation, which is undesirable for hybrid and LES simulations. Laplacian term was discretized using the linear scheme with explicit nonorthogonal correction. Pressure gradient was reconstructed using linear scheme based on the Green–Gauss theorem. The equations for k and ω were discretized in the same manner except the convective term, for which a TVD scheme with Sweby flux limiter was applied.

665

The time discretization has been done using the Crank–Nicolson scheme. For the initialization of the flow in the computational domain the steady RANS solutions have been used. Matrix of simulations discussed in the next two sections is presented in Table 3.

5. Results for the bare hull In our previous paper (Kornev et al., 2012) we have already studied the bare hull flow and the unsteady loadings on a propeller. The weakness of this work was a relatively low resolution. The bare hull was calculated with  2M and the hull with rotating propeller with  4M cells. In this paper the influence of grid resolution on fluctuations and vortex structures in wake as well as the grid dependency of unsteady propeller loadings will be presented. 5.1. Influence of grid resolution on velocity fluctuations and the resolution of the kinetic energy in wake Investigations were performed in the propeller plane at different radii from the propeller axis (see Fig. 13). Circumferential distributions of the r.m.s. of axial velocity fluctuations referred to the mean axial velocity in percent are shown at different radii and grid resolutions. The statistical data were gathered within 10 s of model time. The r.m.s. distributions are strongly fluctuating and the convergence in a classical sense is difficult to recognize. However, one can conclude that 7M cells grid is too coarse to resolve the fluctuations because of smoothing effect of big cells. When the cell number increases from 13M to 45M the r.m.s. level is stabilized at r=R r 0:6. The r.m.s. distributions fluctuate around approximately the same level. At r=R ¼ 1 there is the discrepancy between medium grid 13M and fine grids close to the peak region at Φ  501 and Φ  3101. Since the peak region is very narrow, one can expect that this discrepancy is negligible for propeller loadings prediction. The only unsatisfactory results are documented at 501 o Φ o 1251 and 2251 oΦ o 3101 at r=R ¼ 0:8. The level of fluctuations is gradually reduced when the resolution grows. Taking these results into account further investigations with rotating propellers are performed for 11M and 19.5M grid cells around the hull excluding the propeller area. With the resources available to the authors the use of finer grids is unfeasible for the hull with rotating propeller because of huge necessary computational time. However, in the next subsection we will show that despite unsatisfactory convergence of fluctuations at certain areas in the propeller plane, the convergence for integral forces and moments on the propeller can be attained at moderate resolutions. Assessment of LES implies the check of resolution of the turbulent kinetic energy. Quality of LES is considered as satisfactory if the ratio between the resolved and the total kinetic energy (TKE) K res =K tot is larger than 80 percent. This ratio was studied in the region of our primary interest, i.e. in the propeller plane. The modeled SGS kinetic energy was estimated from the approximation

Fig. 11. Computational grid of KVLCC2 with 7 million of cells. (a) Computational grid of KVLCC2. (b) Blocks of the boundary layer and the propeller disc.

666

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

Fig. 12. Computational grid containing static and rotating (yellow) grids. (a) 13 million cells. (b) 22.5 million cells. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.) Table 3 Simulation matrix. Configuration

Grid

URANS

Hybrid

Grid study

Study of TKE resolution

Study of vortex structures

Unsteady loadings

Hull

7M 13M 25M 45M 13M 22.5M

þ þ   þ 

þ þ þ þ þ þ

þ þ þ þ þ þ

þ þ þ   

þ þ   þ þ

    þ þ

Hull þpropeller

of Mason and Callen (1986): K SGS mod

2

2

¼ ðcD ΔÞ jSj =0:3

ð7Þ

where cD is the dynamic constant of the Smagorinsky model and Sij is the strain rate tensor. In the overlapping areas of URANS and LES the modeled kinetic energy is calculated as the weighted sum of the SGS kinetic energy and the URANS part URANS K mod ¼ ð1  αÞK SGS mod þ αK mod

ð8Þ

where α is the blending factor (see (4)). As seen in Fig. 14 the ratio K res =K tot is more than 90 percent already for moderate resolution of 13M grid. Resolution of TKE is generally not satisfactory for 7M grid although at radius r=R r 0:8 there are regions around the angles Φ  1801 and Φ  01 where this ratio is quite acceptable. The reason can be found from the analysis of λ2 field which shows that these regions are populated by strong large scale vortices which are properly resolved on the coarse 7M grid. 5.2. Vortex structures of the bare ship hull Ship hulls are designed as well-streamlined bodies to prevent the separations and to reduce the resistance as much as possible. That is why the vortices created on ship hulls are not so pronounced as those behind the bluff bodies like, for instance, ship superstructures or offshore constructions. The study of ship vortex structures has attracted the attention of shipbuilders mainly due to inhomogeneity of propulsors' inflow caused by vortices shed from the hull. The vortices entering the propeller disc have strong impact on the (Stern et al., 2013) vibration. Recently, the interest to the ship vortical flows has increased due to new trends in the propulsive efficiency improvements and the development of stealth technologies by Navies (Gorski, 2001). The strongest vortex structures are created in the form of the longitudinal vortices when the ship performs intensive transversal motions like maneuvering at large drift angles or rolling. However, during the longitudinal motion the ship hull is also rich of vortices created at the bow, on bilges and on the appendages. They are influenced by wave and free surface effects. The most significant are the bilge vortices entering the propeller. On tankers and bulk carriers they are formed behind the tubular hull part due to the same flow mechanisms like those in the case of the leading edge vortices of a delta wing. Typical ship vortices are thoroughly discussed and schematically shown in Gorski (2001).

Within the present work the vortex structures have been identified using three different criteria: Q, λ2 (Jeong and Hussain, 1995) and λci (Zhou et al., 1999). We have not observed any substantial difference in vortex structures detected with these criteria. Therefore, the results below are presented only for the λ2 criterion. The λ2 criterion applied to the instantaneous and averaged velocity field reveals common bilge vortices shown in Fig. 15. These structures can be clearly seen in both numerical simulations using URANS (Fig. 15(a) and (c)), and hybrid URANS– LES methods (Fig. 15(b) and (d)). The most important result of the present work and our previous paper (Kornev et al., 2012) is the observation of strong transversal structures in the ship stern area looking like vortex rings shed from the deadwood. In Fig. 15(b) and (d) these structures are visualized by applying the λ2 criterion to the instantaneous velocity field. They appear due to the instability of the shear layer arising at the boundary between the slow wake flow jet and the outer flow. This mechanism is similar to that of creation of transversal vortices at the free jet boundary which is well described in the literature. Being unsteady, these vortices are smoothed out and disappear from the averaged flow pattern (see Fig. 15(b) and (d)). The transversal vortex structures are detected only from the hybrid URANS–LES calculations and cannot be observed in pure URANS calculations due to strong diffusion typical for URANS technique, which is illustrated in Fig. 15 (a) and (c). It is also necessary to point out, that the transversal vortex structures presented in Fig. 15(b) and (d) are observed only behind the bare hull. Presence of working propeller leads to the suction effect and acceleration of the flow inside of the wake jet. The vortices in the shear layer become weaker which prevents the creation of ring-like concentrated structures, see Fig. 16. An important question to be addressed here is the identification of structures when the LES resolution gets higher. A vorticity field consists of structures of different scales. Even the tip vortex of a wing, which is one of the most stable vortex structures, is not a solid structure rather than is formed by a set of vortices with different sizes. If the LES resolution is coarse, the coherent structures can be easily identified if the physical mechanisms responsible for their creation are not damped by large artificial viscosity typical for coarse grids. When the resolution grows the vorticity field shows many fine details which make the structure identification more complicated. For that, it is necessary to apply either the Proper Orthogonal Decomposition (POD) analysis or a

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

667

Fig. 14. Ratio between the resolved and the total kinetic energy at different radii in the propeller plane.

Fig. 13. Axial velocity fluctuations along circles at different radii in the propeller plane.

local time or spatial averaging. For instance, in our previous investigations of vortex structures on dimpled heat exchanger (Voskoboinick et al., 2013; Turnow et al., 2012), we observed the presence of long period oscillations in velocity and pressure measurements. The reason for this phenomenon was found from URANS simulations which revealed the big vortex structure

changing its orientation. LES confirmed the long period oscillations. However, the vortex field was filled with small vortices chaotically distributed within the dimple. Only local time averaging allowed one to recognize, that small vortices acts in the same manner as the big solid structure. With the other words, their statistically averaged collective effect can be treated as the effect of a single coherent vortex structure. We give this sample to show that simulations with different level of resolution are often necessary in order not to miss important flow features. In our previous work the transversal structures behind the bare hull were obtained on coarse grid with 1.8M. As seen from Fig. 15

668

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

Fig. 15. Visualization of vortex structures by λ2 ¼  1 criterion applied to the velocity field in the ship stern area. (a) k  ω-SST: Vortex structures in ship stern area (7 M grid). (b) Hybrid method: Vortex structures in ship stern area (7 M grid). (c) k  ω-SST: Vortex structures in ship stern area (13 M grid). (d) Hybrid method: Vortex structures in ship stern area (13 M grid).

Fig. 17. Coordinate system and forces directions convention (Batrak et al., 2012).

6. Results for the hull with rotating propeller The hybrid method has been applied to the complex system with both the ship hull and the propeller VP1356 taking all interaction effects into account. The coordinate system and directions of forces are shown in Fig. 17.

6.1. Grid influence

Fig. 16. Vortex structures around the hull with rotating propeller visualized by λ2 ¼  50 criterion. (a) k   ω SST : 13 M grid. (b) Hybrid: 13M grid. (c) Hybrid: 22.5M grid.

(b) and (d) the transversal structures remain in fine grid simulations without any additional filtering although they become more irregular due to fine details appearing with increase of grid resolution. Results of these calculations confirmed the presence of transversal vortex structures obtained first in our coarse LES simulations (Kornev et al., 2012).

Fig. 18 shows standard deviations of thrust Px and moment Mx versus time using 13M and 22.5M grids. Both the thrust and moment standard deviations converge in time. Despite grid dependency problems for velocity fluctuations at certain places within the propeller plane mentioned in Section 5.1, the thrust and the moment are proved to be independent on the grid resolution. This can be due to two effects. First, the fluctuations presented in Fig. 13 are computed in the nominal wake. Due to the suction effect of propeller, separation effects are diminished, vortices shed from the hull becomes weaker in the effective wake and the velocity fluctuations decrease. Second, the propeller has an integrating effects levelling the differences in velocity fluctuations.

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

669

Fig. 18. Convergence of the thrust and torque ðM x Þ solution.

6.2. Unsteady loadings on propeller Fig. 19 presents the spectra of the thrust fluctuations E(f). As seen five dominating blade passing frequencies (BPF) up to 5 nZ  200 Hz are resolved in our simulations. Comparison between the hybrid and the URANS k   ω SST results shows that the magnitude of the thrust fluctuations in the hybrid simulation is much higher than in the k   ω SST one. Both spectra are continuous although at frequencies larger than the first BPF both spectra can be approximated as the discrete ones. The peaks corresponding to BPF in the k   ω SST spectrum are narrower than in the hybrid one. The most remarkable difference is the presence of low frequency modes of the thrust fluctuations at f o nZ ¼ 40 in hybrid simulations. As shown in the subfigure the most contribution to the full r.m.s. of the thrust fluctuations 1 Rf R defined as the ratio Iðf Þ ¼ EðχÞ dχ= EðχÞ dχ is made by the low 0

0

frequency modes. It is almost 70% in the range f onZ ¼ 40 whereas the contribution of the first blade passing frequency f ¼nZ is about 25%. The contributions of the next BPFs is about 5%. Certainly, this result raises questions about its plausibility because of the restricted time of statistical processing and the transitional process which numerical solution needs to develop from an initial state. The latter effect is easy to prove by comparison of results obtained for different periods of statistical processing. From our estimations the processing time period of 10 s used further does not contain the transitional process. The smallest captured frequency is 2π=10 s  1 whereas the Nyquist frequency is 104 s  1 . The aliasing effects from high frequencies f 4 104 s  1 for fast decaying discrete spectra should be negligible. Taking these facts into account, we believe that the long period oscillations are realistic and caused by the vortex structures shed from the hull which are not resolved by URANS. The presence of low frequency modes for resistance and transversal force was also revealed in DES simulations by Xing et al. (2012) who found the relation between these modes and typical vortex structures shed into the wake. In the ship hydromechanics the common procedure to quantify the unsteady propeller loadings is based on the force representation in form of the series on BPF mZn: P x ¼ P x þ N P ðam cos ðmZntÞ þ bm sin ðmZntÞÞ, where P x is the mean thrust.

Fig. 19. Spectra of the thrust E(f) obtained from the hybrid (13M grid, grey lines) Rf and k  ω SST (13M grid, black lines) models. Iðf Þ ¼ 1 EðχÞ dχ, r.m.s stands for the full R 0 root mean square of the thrust fluctuations r:m:s ¼ Eðf Þ df . (For interpretation of 0 the references to color in this figure caption, the reader is referred to the web version of this paper.)

BPF analysis can result in a substantial underestimation of r.m.s. especially for continuous spectra like for this from the hybrid simulations case. To the knowledge of the authors this problem was not addressed in the literature so far, perhaps, because the long period fluctuations of the propeller thrust are not critical from the strength point of view. However, this problem requires additional investigation in the future. The results for standard deviations of forces and moments are summarized in Fig. 20. Standard deviations based on engineering methods were calculated using amplitudes of different harmonics published in Batrak et al. (2012): P 0x ¼ ððP x  P x Þ2 Þ1=2 2 !2 31=2 ZT N X 1 4 ¼ lim Ak sin ðknZt þ ϕk Þ dt 5 T-1T k¼1

ð9Þ

0

m¼1

It implies the discrete spectrum of the thrust neglecting the low qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 frequency modes f o nZ. Three first magnitudes cm ¼ a2m þ bm of the BPF series are distributed as ðc1 : c2 : c3 Þ=P x ¼ 1:84 %: 0:47 %: 0:18% for the hybrid simulation and ðc1 : c2 : c3 Þ=P x ¼ 0:6 %: 0:06 %: 0:047% for k   ω SST one. Since the peaks around BPF are smoothed the coefficients cm were calculated from the condition that the energy of single BPF modes are equal to the energy of corresponding smoothed peaks. Our analysis shows clearly that

where Ak is the amplitude of kth mode, N is the number of harmonic modes (usually N ¼ 2) and ϕk is the phase displacement. Amplitudes Ak are published in Batrak et al. (2012) for different blade numbers. Each amplitude is represented as the sum of a mean value which is valid for all possible ships and a deviation, i.e. dev Ak ¼ Amean 7 Adev accounts for the variation of Ak depending k k . Ak on ship types. Particularly, for full bottomed ships this reads mean Amax ¼ Amean þAdev are referred in k k k . Results obtained using Ak Fig. 20 to as “mean”, whereas the designation “max” corresponds

670

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

Fig. 20. Standard deviations of the forces and moments calculated using different approaches: V: Veritec, W: Wereldsma, S: Scheme B and CFD hybrid computations with 13M and 22M grids. All force standard deviations are referred to the mean thrust, whereas the moments deviations to the mean torque. Other designations stand for: Px: thrust, Py: vertical force, Pz: horizontal force, Mx: torque, My: horizontal moment, Mz: vertical moment.

Fig. 21. Variation of the thrust in time (13M grid).

Fig. 22. Cumulative distribution function of the thrust.

to Amax . Since the phase displacement is not known, it is set to zero k to get upper estimation for standard deviations. The same procedure is applied to all forces and moments. All force standard deviations are referred to the mean thrust, whereas the moments deviations to the mean torque. Although the representation (9) is based on the BPF analysis and can be not quite accurate this disadvantage is diminished due to adjustment of the coefficients Ak to match the experimental data. The standard deviation of the thrust varies in the range between 1.86% and 4%, the vertical forces between 0.35% and 3.38%, the horizontal force between 0.83% and 3.74%, all referred to the mean thrust. For the torque, the results of all methods, except the Veritec approximation, are in a good agreement with each other. Taking into account that the safety factor commonly used for the forces and moments variations in practical shaft

calculations is taken a few times larger than that obtained from simple estimations, the discrepancy between the methods can be considered as insignificant. The situation with the horizontal M 0y ¼ ððM y  M y Þ2 Þ1=2 and the vertical M 0z ¼ ððM z  M z Þ2 Þ1=2 moment variations is more complicated. Results by the Veritec method substantially deviates from the results of the other methods. The estimation for the horizontal moment fluctuations around 2–3.2% percent has been obtained from hybrid methods and Scheme B approach, whereas the Wereldsma's estimation is around 7%. The vertical moment fluctuations are around 5% of the mean torque according to Wereldsma (1964), Scheme B and hybrid CFD computations, whereas the Veritec method predicts up to 17 percent, which is more than three times as large as the others. This analysis shows that the accuracy of predictions of standard deviations of forces and moments provided by simple engineering methods like Wereldsma (1964) and Scheme B can be considered as quite acceptable for practical purposes. The best agreement with CFD is attained for the Scheme B method. All engineering methods are not capable of predicting the peak loadings. Fig. 21 presents the time history of the thrust within 10 s. While the standard deviation is moderate (2.9% for 13M and 2.5% for 22.5M), a few strong peak force loadings up to 13% act on the propeller shaft. The peaks of vertical and horizontal forces are approximately twice as small as that of the thrust. The peaks of the torque, horizontal and vertical moment fluctuations are, respectively, 10, 7.83 and 18.25% referred to the mean torque. The occurrence frequency of peak loadings of a certain magnitude can be estimated from the cumulative probability density distribution diagram presented in Fig. 22 for the thrust variations.

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

The horizontal axis shows the thrust related to the mean thrust P x =P x whereas the vertical axis illustrates the probability of the event that the instantaneous thrust is smaller than P x =P x . For instance, the probability of the event that the thrust is in the range between 90% and 103% of the mean thrust is 0.85. In other words, the probability that the thrust fluctuation is larger than its standard deviation is almost 15%. The amplitude of thrust fluctuations can be of order of three standard deviations. URANS computations with k   ω SST model show regular periodic oscillations of small amplitude which is much less than that from hybrid method and engineering estimations. While the averaged propeller loadings and the standard deviations are reliably predicted by engineering methods, hybrid techniques based on the combination of LES and URANS approaches are necessary to detect the peak loadings on marine propellers of full bottomed ships. It can be supposed from the analysis of the literature (Vartdal et al., 2009) that the unsteady loadings are getting much larger during the maneuvering motions resulting in the damage of the shaft bearings. To our opinion, the LES or hybrid approaches are especially important for simulation of propeller hydrodynamics at transient motion conditions.

7. Conclusion The paper presents recent achievements in the development and application of the hybrid LES–URANS technique (Kornev et al., 2011) for the prediction of the flow in the stern area and unsteady loadings on propellers of full bottomed ships. Detailed analysis of the wake flow performed using λ2 criterion revealed the existence of strong transversal vortex structures in the wake of a bare ship hull which have not been documented in the literature. They disappear for the whole arrangement consisting both the propeller and the hull due to suction effect caused by the propeller. While our previous work Kornev et al. (2011) was dedicated to the calculations of the bare hull the focus of the present paper is the calculation of the whole system containing both the propeller and the hull with consideration of all interaction effects. The unsteady loadings on the propeller of the KVLCC2 tanker have been calculated with numerical URANS and hybrid methods as well as with a few engineering approaches. Analysis of spectra shows that the magnitude of the thrust fluctuations in the hybrid simulation is much higher than in the k   ω SST one. Both spectra are continuous although at frequencies larger than the first BPF both spectra can be approximated as the discrete ones. The most remarkable difference is the presence of low frequency modes of the thrust fluctuations at f o nZ ¼ 40 in hybrid simulations which are caused by vortex structures shed from the hull. Results show that the accuracy of predictions for standard deviations of forces and moments on propellers provided by simple engineering methods like Wereldsma (1964) and Scheme B ones can be considered as quite acceptable for practical purposes. However, all these methods are not capable of predicting the peak loadings which occur frequently. For instance, the thrust peak fluctuation can attain the value of up to 13% of the mean thrust. Hybrid techniques based on the combination of LES and URANS approaches are necessary to detect the peak loadings on marine propellers of full bottomed ships. At present, the existence of such peaks has been revealed only in numerical investigations and should be confirmed by experimental data which are still missing.

Acknowledgments The authors acknowledge gratefully the support from the Federal Ministry for Economic Affairs and Energy (BMWi) (Grant

671

03SX261A). A part of this work has been performed within the framework of the project ShipLES monitored by Dr. A. Nitz. CFD calculations have been performed on IBM pSeries 690 Supercomputer at the North German Alliance for the Advancement of High-Performance Computing (HLRN). The authors wish to acknowledge the community developing the OpenFOAM. Special thanks are to Andrey Taranov and Evgeny Shchukin who developed the initial version of the hybrid code. We are also very grateful to Dr. Yuri Batrak et al. (2012) from Intellectual Maritime Technologies for explanation of his results (Batrak et al. (2012)) and many helpful discussions.

References Alte, B., Baur, M., 1986. Propulsion, Handbuch der Werften, Hansa Verlag, Hamburg, pp. 183–192. Batrak, Y.A., Shestopal, V.P., Batrak, R.Y., 2012. Propeller hydrodynamic loads in relation to propulsion shaft alignment and vibration calculations. In: Proceedings of the Propellers–Shafting 2012 Symposium, 10.1–10.14. Breslin, J.P., Andersen, P., 1994. Hydrodynamics of Ship Propellers. Cambridge Ocean Technology Series. Buice, C.U., Eaton, J.K., 1996. Experimental investigation of flow through an asymmetric plane diffuser, CTR Annual Research Briefs, Center for Turbulence Research. Stanford University, Stanford, pp. 243–248. Buice, C.U., 1997. Experimental investigation of flow through an asymmetric plane diffuser (Ph.D. thesis). Stanford University, August 1997. Carlton, J.S., 2007. Marine Propellers and Propulsion. Butterworth – Heinemanns, Jordan Hill, Oxford OX2 8DP, USA. Farrell, P.E., Maddison, J.R., 2011. Conservative interpolation between volume meshes by local Galerkin projection. Comput. Methods Appl. Mech. Eng. 200, 89–100. Gorski, J.J., 2001. Marine vortices and their computation. In: RTO Applied Vehicle Technology Panel (AVT) Symposium, Leon, Norway. Hayati, A.N., Hashemi, S.M., Shams, M., 2013. A study on the behind-hull performance of marine propellers astern autonomous underwater vehicles at diverse angles of attack. Ocean Eng. 59, 152–163. Jang, H., Mahesh, K., 2010. Large Eddy Simulation of Marine Propulsors in Crashback. In: 28th Symposium on Naval Hydrodynamics, Pasadena, California, 12– 17 September. Jeong, J., Hussain, F., 1995. On the identification of a vortex. J. Fluid Mech. 285, 69–94. Kornev, N., Taranov, A., Shchukin, E., Kleinsorge, L., 2011. Development of hybrid urans–les methods for flow simulation in the ship stern area. Ocean Eng. 38, 1831–1838. Kornev, N., Taranov, A., Shchukin, E., Springer, J., Palm, M., Batrak Y., Development, application and validation of hybrid URANS–LES methods for flow simulation in the ship stern area, In: Proceedings of the 29th Symposium on Naval Hydrodynamics Gothenburg, Sweden, August 2012, 26–31. Larsson, L., Stern, F., Bertram, V., (Eds.), 2000. Gothenburg 2000: A Workshop on Numerical Ship Hydrodynamics, Chalmers University of Technology, Gothenburg. Mach, K.-P., 2011. Report 3754: LDV velocity measurements with the model propeller VP1304, Potsdam Ship Model Basin. Mason, P.J., Callen, N.S., 1986. On the magnitude of the subgrid-scale Eddy coefficient in large-Eddy simulations of turbulent channel flow. J. Fluid Mech. 162, 439–462. Menter, F.R., 1993. Zonal two equation k – ω turbulence models for aerodynamic flows. AIAA Papers, 2893–2906. Weller, H.G., Tabor, G., Jasak, H., Fureby, C., 1998. A tensorial approach to CFD using object oriented techniques. Comput. Phys. 12 (6), 620–631. NASA, 2003, 〈http://turbmodels.larc.nasa.gov/flatplate.html〉. Noise and vibration Group, 1985. Vibration Control in Ships. VERITEC :Marine Technology Consultants (Norway), Oslo. Obi, S., Aoki, K., Masuda, S., 1993. Experimental and computational study of turbulent separating flow in and asymmetric plane diffuser. In: Ninth Symposium Turbulent Shear Flows, Kyoto, Japan, August 16–19, 1993, pp. 305. Park, S., Park, S.W., Rhee, S.H., Lee, S.B., Choi, J.-E., Kang, S.H., 2013. Investigation on the wall function implementation for the prediction of ship resistance. Int. J. Nav. Archit. Ocean Eng. 5, 33–46. Schlichting, H., Gersten, K., 2000. Boundary-layer theory. Springer Verlag, Berlin, Heidelberg. Stern, M., Yang, J., Wang, Z., Hosseini, H.-S., Mousaviraad, M., Bhushan, S., Xing, T., 2013. Computational ship hydrodynamics: nowadays and way forward. Int Shipbuild. Prog. 60 (1–4), 3–105. Tomas, S., Eiff, O., Masson, V., 2009. Experimental study of the turbulent structure of a boundary layer developing over a rough surface. 19éme Congrés Francais de Mécanique. Thomas, Russell H., Schetz, Joesph A., Pelletier, Dominique H., 1991. Threedimensional finite element method analysis of turbulent flow over selfpropelled slender bodies. J. Propuls. Power 7.2, 281–287.

672

N. Abbas et al. / Ocean Engineering 104 (2015) 659–672

Turnow, J., Kornev, N., Hassel, E., Zhdanov, V., 2012. Flow structures and heat transfer on dimpled surfaces. Int. J. Heat Fluid Flow 35, 168–175. Vartdal, B.J., Gjestland, T., Arvidsen, T.I., 2009. Lateral propeller forces and their effects on shaft bearings. In: First International Symposium on Marine Propulsion SMP09, Trondheim, Norway. Voskoboinick, V., Kornev, N., Turnow, J., 2013. Study of Near Wall Coherent Flow Structures on Dimpled Surfaces Using Unsteady Pressure Measurements. Flow Turbul. Combust. 90 (4), 709–722. Voitkunski (Ed.), 1985. Handbook on Ship Theory, vol. 1. pp. 564–572.

Wereldsma, R., 1964. Propeller excited shaft and hull vibrations of single screw ships. Int. Shipbuild. Prog. 11 (124), 547–553. Xiao-fei, H., Zhen-yu, H., Fang-wen, H., 2009. Unsteady hydrodynamics forces of propeller predicted with viscous CFD. J. Hydrodyn. (Ser. A) 2009-06. Xing, T., Bhushan, S., Stern, F., 2012. Vortical an turbulent structures for KVLCC2 at drift angle 0, 12, and 30 degrees. J. Ocean Eng. 55, 23–43. Zhou, J., Adrian, R.J., Balachandar, S., Kendall, T.M., 1999. Mechanisms for generating coherent packets of hairpin vortices in channel flow. J. Fluid Mech. 387, 353–396.