A Calibration Study with CFD Methodology for Self-Propulsion Simulations at Ship Scale

Abstract

This paper summarises the main findings from the full-scale Computational Fluid Dynamics (CFD) analyses conducted at SINTEF Ocean on the case of MV REGAL, which is one of the benchmark vessels studied in the ongoing joint industry project JoRes. The numerical approach is described in detail, and comparative results are presented regarding the propeller open water characteristics, ship towing resistance, and ship self-propulsion performance. The focus of numerical investigations is on the assessment of the existing simulation best practises applied to a ship-scale case in a blind simulation exercise and the performance thereof with different turbulence modelling methods. The results are compared directly with full-scale performance predictions based on model tests conducted at SINTEF Ocean and sea trials data obtained in the JoRes project.

1. Introduction

Model testing has been, and still is, the dominant method to predict the propulsion performance of ships. However, the use of Computational Fluid Dynamics (CFD) methods for the same purpose is becoming increasingly common and is also becoming an integral part of the ship design process. The greatest advantage of CFD is the possibility to directly perform predictions at full scale, thereby reducing uncertainties related to the influence of Reynolds number (known as scale effect). Numerical simulations are conducted in a strictly controlled environment and are therefore characterized by high repeatability. They can be conducted in a time- and cost-effective manner, especially when best practices are established and automated templates are developed.

Considering the increasing importance of numerical simulations in both research and industrial settings, since as early as 2011, ITTC has been developing and continuously updating the recommended procedures and guidelines for CFD analysis of ship performance [1], with specific sections addressing resistance and flow, propulsion, and manoeuvrability, as well as uncertainty analysis and quality assurance. The recent IMO Resolutions MEPC.350(78) and MEPC.351(78) consider numerical calculations as an acceptable way to derive a ship’s performance in the EEXI regulatory framework [2]. To this end, the International Association of Classification Societies (IACS) has proposed dedicated guidelines concerning the definition of the reference speed used for forecasting the impact of energy efficiency technologies applied to a ship through CFD simulations [3]. The said guidelines include, among other things, three optional approaches to perform calibration of a CFD method used for speed and power prediction: 1-Calibrated CFD with sea trials or model tests of parent hull; 2-Calibrated CFD with model tests of similar ships; 3-Calibrated CFD with sea trials of a set of comparable ships. It is generally assumed that the calibration factor, defined as the ratio between the shaft delivered power found from sea trails (or model test prognosis) and the same quantity found from a CFD calculation, should range between 0.95 and 1.05 to be accepted by a verifier without further technical justification. At the same time, it is understood that successful calibration alone does not necessarily guarantee that the CFD method in question correctly solves the governing equations of the problem or that it accurately represents physical phenomena significant to the problem solution. The classical example is given by the RANS method, which normally shows sufficient accuracy in the prediction of ship resistance and propulsion characteristics but reveals limitations in the analysis of propeller cavitation in the wake behind the ship hull and associated pressure pulses on the hull and radiated noise [4]. To answer the latter questions, verification and validation procedures need to be employed systematically, following formal procedures for error and uncertainty analysis [5,6,7,8].

In the public domain, the experimental material on ship propulsion available for verification, validation, and calibration exercises has largely been limited to the results of model tests, with such benchmark ship cases as KCS, KVLCC1, KVLCC2, JBC, DTMB5415, and DTC. An example of these can be found in the proceedings of the CFD Workshops in Ship Hydrodynamics held in Gothenburg (2010) [9] and Tokyo (2015) [10]. More recently, new benchmark datasets have emerged, dealing with contemporary designs of next generation ships such as, for example, SOBC-1 by SINTEF Ocean [11]. A “blind” comparison of CFD predictions with sea trials data in full-scale has been addressed in the 2016 Workshop on Ship Scale Hydrodynamic Computer Simulation for the case of the general cargo vessel MV REGAL [12]. Presently, the same ship is also used in the ongoing JoRes project [13], where SINTEF Ocean contributes with both the model tests and CFD simulations. The findings from the referred works highlight that numerical results often appear to be sensitive to how the CFD model is constructed, especially the computation grid, near-wall treatment, surface roughness model, boundary conditions, and turbulence modelling approach. Representative examples of the scatter found in practical full-scale predictions are documented in [12]. To improve the predictive capabilities of CFD methods and increase confidence in numerical prognoses, calibration and validation of CFD methods against model tests and sea trials data need to be continued. Such efforts are also seen as important steps towards better understanding the physics behind the effect of Reynolds number and elaboration of scaling procedures employed by testing facilities.

In the application to ship-scale CFD simulations, a significant amount of research has recently been performed to elaborate the guidelines and recommendations for this type of analysis. For instance, grid sensitivity studies with the RANS method are addressed in [14], while grid requirements for LES simulations of full-scale ships are discussed in [15]. The choice and implications of the turbulence modelling approach are considered in works such as [16,17,18]. The revision of general guidelines towards full-scale ship performance prediction is suggested in [19], with specific emphasis on the choice of Y+ and Courant number levels. The latter work makes use of the MV REGAL dataset from the 2016 Workshop [12].

In the present study, the open benchmark dataset of the vessel MV REGAL, further detailed in Section 2, is chosen by the authors to assess the accuracy of their current CFD simulation practises and study the influence of different turbulence modelling approaches for ship resistance, propeller open water characteristics, and ship self-propulsion calculations.

The range of turbulence modelling approaches includes the conventional Reynolds-averaged Navier–Stokes (RANS) method and Scale-Resolving Simulation (SRS) techniques such as Detached Eddy Simulation (DES) and Large Eddy Simulation (LES). The LES results are only provided for the case of propellers in open water. All simulations are performed in full-scale conditions. The results of resistance and open water calculations are compared with the predictions using SINTEF Ocean’s scaling procedure, which is based on the results of model tests. The results of self-propulsion simulations are compared with both the performance predictions and full-scale trials, considering propeller RPM and shaft delivered power at ship self-propulsion points for several speeds. Moreover, the cavitation patterns on propeller blades are visually compared with full-scale observations on cavitation for selected conditions. Due to the restrictions regarding the data distribution that currently apply in the JoRes Consortium, most of the comparisons are presented in the form of relative differences and plotted without scale.

All CFD results were produced in the spring of 2021 through “blind” simulations, adhering to SINTEF Ocean’s modelling best practises, before either sea trials data or model test results became available. This places the present study in the category of CFD calibration exercises, whereas the employed simulation templates have been subject to earlier verification and validation studies using both the open access data and SINTEF Ocean internal datasets. The said studies included solution sensitivity as regards the influence of mesh resolution (global and local), near-wall treatment, time step, and numerical discretization schemes. The findings from this research suggest the following quantification of uncertainties in CFD predictions when using the task-specific templates employed in this work: towing resistance, −2 to 3%; propeller open water characteristics, −3 to 4%; ship self-propulsion characteristics, −6 to 8%.

2. Benchmark Vessel MV REGAL

MV REGAL, depicted in Figure 1, is a single-screw general cargo vessel that is equipped with a 4-bladed fixed pitch propeller. The ship features a semi-balanced rudder with a fixed horn and does not possess tunnel thrusters. The key specifications of the MV REGAL vessel are listed in

Table 1. Main particulars of MV REGAL.

Hull Particulars			Propeller Particulars
Length between perpendiculars LPP	138	m	Propeller diameter D	5.2	m
Breadth B	23	m	Pitch ratio P(0.7)/D	0.6781
Gross tonnage	11,542	t	Blade area ratio Ae/Ao	0.57

. Additional details concerning its propeller and rudder arrangements can be found in [12,20].

It needs to be remarked that the geometries provided for the JoRes CFD Workshop differ from those utilised in the earlier Lloyd’s Register Workshop organised in 2016 on the same ship case. While the 2016 Workshop used the STL models of the ship hull, propeller, and rudder derived from the 3D laser scan data, the present JoRes geometries were provided as clean solids that were prepared using both the reconstructed 3D laser scan data and original documentation and drawings. Further, such constructive features as bilge keels and protective anodes installed on the ship hull and rudder were removed from the CFD model, and no respective corrections were applied in the performance prediction based on model tests. The approximate geometry models of the superstructure, deck hatches and coamings, and the two cranes on the deck were included in the resistance and self-propulsion CFD simulations to account for the aerodynamic resistance of the ship in a more accurate manner. The initial hydrostatic position of the ship was set according to the draught measurements before the JoRes sea trials, resulting in draught values at the fore and aft marks of 2.97 m and 5.865 m, respectively, which correspond to ballast condition. The presence of hull sagging, hogging, and list is disregarded in the draught measurements. The general view of the MV REGAL self-propulsion simulation setup is depicted in Figure 2.

3. Research Methodology

For the numerical simulations conducted in the present study, the commercial CFD software STAR-CCM+ (version 15.06.007-R8) was employed. A trimmed hexahedral mesher with prism layers on the wall boundaries was employed in all the main types of simulations (towing resistance, propeller in open water, and ship self-propulsion calculations), as well as additional simulations with a flat plate, which were used to assess the influence of surface roughness height. All simulations, except for the flat plate calculations, were carried out in a time-dependent manner using the implicit unsteady segregated flow solver. Simulations involving rotating propeller (open water and self-propulsion) were performed using the sliding mesh technique to fully account for the interaction between the rotating and stationary components in the setup. The properties of water and air used in the simulations were derived from recorded values during the sea trials of MV REGAL. The simulations were performed for full-scale conditions corresponding to sea trials, as specified in the case description for the JoRes CFD Workshop regarding the MV REGAL vessel [20]. Three ship speeds (9, 10.5, and 12 knots) were investigated in the resistance and self-propulsion scenarios, while in the open water case the advance coefficient (J) values of 0.2, 0.3, 0.4, 0.5, and 0.6 were used. The case-specific details of numerical setups are addressed for each individual simulation scenario in their respective subsections below.

A high-Reynolds near-wall treatment method was employed in the analyses. While using fine near-wall meshes with Y+ < 5 offers advantages in accurately predicting the frictional component of forces and moments, as well as modelling boundary layer separation/detachment, especially in Scale-Resolving Simulations, it becomes impractical in a full-scale case due to the excessively large mesh size and the small time step required to maintain the desired Courant number level in the areas of mesh refinement. In this regard, one needs to remember that it is not sufficient to increase mesh density only in the direction normal to the wall. Eddies developing in the near-wall region also require fine mesh resolution in the spanwise and streamwise directions. Failing to capture those eddies may compromise the overall simulation quality. Therefore, avoiding the resolution of viscous sub-layers may be a more reliable approach for solving high-Reynolds flows, even with such techniques as DES and LES. Another reason for choosing the high Y+ near-wall treatment is the inclusion of surface roughness, which relies on the use of roughness-modified wall functions. The present setup employs blended wall functions, supporting the so-called “All Y+ Treatment” method.

Separately, as a part of the JoRes project, SINTEF Ocean has conducted a model test campaign with the MV REGAL vessel. The campaign included towing resistance, open water, and propulsion tests in calm water at the scale of 1:23.111. The results of model tests were used in full-scale performance prediction for the sea trial conditions, following the standard procedure applied at SINTEF Ocean for single-screw ships [21]. The results of performance prediction were compared against the sea trials data, which was post-processed according to the ISO15016 standard [22], as well as the CFD calculations conducted at full-scale. The comparisons focused on propeller RPM and propeller shaft power, P_D. Additional comparisons were made between the CFD calculations and model test results extrapolated to full-scale conditions, examining the towing resistance of the ship, its dynamic position, and propeller characteristics in open water.

3.1. Open Water Propeller Simulations

The open water simulations were performed in full-scale using the same propeller setup as in self-propulsion conditions. This means that, unlike a conventional open water model test setup, the propeller was not driven from downstream but from upstream. As a result, it operated in a pushing mode behind the ship hull, similar to the setup in propulsion conditions, and had the same hub cap. This setup is illustrated in Figure 3, which also provides an overview of the overall mesh refinement pattern around the propeller.

The cylindrical rotating propeller region used in the sliding mesh calculation has the dimensions 0.125D (upstream), 0.288D (downstream), and 0.5385D (radius), measured from the propeller plane, where D represents the propeller diameter. These dimensions are smaller than those typically applied in the SINTEF Ocean standard open water CFD setup, because, in the self-propulsion case studied in this work, the propeller region had to be accommodated within a tight space between the ship hull and the rudder, as shown in Figure 2. The rationale was then to use the same region dimensions in open water calculations to avoid the influence of region size when deriving propulsion factors. The respective dimensions of the main fluid region (also cylindrical in this case) were 5D (inlet), 15D (outlet), and 5D (radius), measured from the propeller centre. Additional simulations were performed to assess the influence of the downstream extension of the sliding mesh region on propeller open water characteristics and the quality of propeller slipstream resolution.

The cell size in the propeller region and the first (finest) volumetric control around the propeller and slipstream (as depicted in Figure 3) were set to 1.3% of the base size, which corresponds to the propeller diameter, D. The volumetric control extends to the distance of 3D downstream. On the propeller blades, the maximum target cell size is 0.65% of base, and it is reduced to the minimum size of 0.0203125% of base at the blade edges and tip. The edge and tip refinements are achieved using the surface patches extracted from the initial CAD geometry of the propeller, as shown in Figure 4. In the absence of such patches, similar refinement can be achieved by means of volumetric controls in the shape of a tube following the leading edge feature curve, typically obtained from blade surface wireframe. The propeller model included a gap between the rotating propeller hub and the stationary shaft, measuring 17 mm in size, which is consistent with the self-propulsion setup. The inclusion of the hub gap in the numerical model allows for more accurate computation of forces acting on the propeller by avoiding uncertainty related to the integration of pressure on the side of the propeller hub facing the shaft. The total number of cells in the open water setup was 22.1 million, of which 18.6 million were accommodated in the propeller region.

Figure 5 illustrates the arrangement of prism layers on propeller blades. The prism layer mesh consists of 10 layers. The height of the first near-wall cell is selected to target an average Wall Y+ approximately 60 in the middle section of the blade. This choice implies the use of wall functions, which offer computational efficiency and inclusion of surface roughness. The total thickness of prism layers is 0.235% of the propeller diameter. The present settings result in a layer stretch factor of approximately 1.35 and Wall Y+ distribution shown in Figure 6. With a coarse near-wall mesh, one cannot achieve a uniform distribution of Y+ over the whole blade. However, as depicted in the figure, the applied settings provide a favourable Y+ range between 30 and 100, avoiding buffer zones everywhere except the gap and hub vortex separation area. The same Y+ range on the propeller is also met at other advance coefficients and in self-propulsion calculations.

In the open water scenario, several turbulence modelling approaches were investigated. These included (i) the traditional unsteady RANS method using the k-$\mathsf{\omega }$ SST model with linear constitutive relation [23]; (ii) the Improved Delayed Detached Eddy Simulation (IDDES) method [24], which incorporates a subgrid length-scale dependence on the wall distance. This allows the RANS part of the solution to be utilised in the thin near-wall region where the wall distance is smaller than the boundary layer thickness. The DES formulation with the k-$\mathsf{\omega }$ SST turbulence model in the RANS zones was employed [25]; (iii) Large Eddy Simulation (LES) method with the Smagorinsky Subgrid Scale model [26] and Modified Van Driest damping function [27]; (iv) Scale-Resolving Hybrid (SRH) turbulence model [28], which is a continuous hybrid RANS-LES technique. It switches continuously (unlike the DES method) from the RANS model to the LES model when the mesh resolution is fine and the time step is small. Similar to the DES solution, the SRH method was applied with the k-$\mathsf{\omega }$ SST model in the RANS part of the solution.

Open water calculations were performed with both smooth and rough propeller surfaces. The influence of surface roughness was accounted for by using roughness-modified wall functions, which shift the log layer of the inner part of the boundary layer closer to the wall. Mathematically, it is achieved by means of a roughness function [29] that modifies the log law offset coefficient depending on the equivalent sand-grain roughness height, viscosity, and velocity scale. The influence of roughness was investigated only with the RANS and DES methods, considering the two values of roughness height: 8.68 $\mathsf{\mu }$m resulting from the JoRes propeller surface roughness measurements and 30.0 $\mathsf{\mu }$m, which is a standard value of sand-grain roughness used at SINTEF Ocean in CFD calculations on older ship propellers in service that have undergone cleaning and polishing before trials.

The time-accurate solution is achieved by rotating the sliding mesh region (propeller region) about the shaft axis with a specific angular step. The time step corresponding to 2 deg of propeller rotation was applied in both the open water and self-propulsion calculations, which in the authors’ experience is sufficient in most practical cases. However, an additional test with a step of 1 deg was conducted to assess the sensitivity of the LES and SRH models to time step size. Representative levels of Courant number obtained with the time step of 2 deg are shown in Figure 7. Reducing the time step to 1 deg lowers the Courant number by a factor of 2.

The open water simulations were performed using the multi-phase flow formulation with the Volume of Fluid (VOF) model, as in the self-propulsion calculations, but with the reference pressure set to atmospheric conditions to prevent occurrence of cavitation.

3.2. Hull Resistance Simulations

The resistance simulations were performed in full scale using a geometry model consisting of the ship’s hull, rudder, and propeller hub. The calculation matrix included both cases with and without superstructure and cranes to evaluate their influence on the ship’s resistance and dynamic position. The initial hydrostatic position for the resistance simulations was the same as in the self-propulsion case, and it corresponded to the draught marks provided by JoRes as specified earlier.

The dimensions of the rectangular computation domain were assigned according to existing practises to minimise the influence of the outer boundaries on the quality of the numerical solution: 4LPP in X (inlet and outlet) and Y (side boundaries) directions from the ship’s aft perpendicular/centre plane, 2LPP in the Z direction to the bottom boundary, and 1LPP in the Z direction to the top boundary, from the base. As an additional measure to mitigate wave reflections, the VOF damping method described in [30] was employed, with the damping zone extending to a distance of 2LPP from the inlet, outlet, and side boundaries. The resistance simulations were performed using the full hull model.

The Trimmer hexahedral mesh employed in the resistance case was designed to ensure refinement near the free surface, with a size equal to 0.1% of LPP in the Z direction. In the bow and stern regions, this size also applies to the isotropic cells around the hull. In the mid-ship area, the cells are anisotropic, with aspect ratios varying from 2 to 4. Conventional Kelvin’s wake refinement controls were used to capture the wave systems generated by the ship. Figure 8 and Figure 9 give illustrations of the mesh around the ship hull with superstructure and cranes, and in the area of free surface.

Additional mesh refinement controls are implemented around the propeller and rudder locations to ensure an adequate level of refinement for capturing the key characteristics of the separated hull wake. In this specific area, the refinement pattern is isotropic, and the cell size is set to 1.35% of the propeller diameter, which matches the size applied in the propeller slipstream in open water and self-propulsion calculations.

The prism layer mesh on the hull consists of 16 layers, with the height of the first near-wall cell chosen to target an average Wall Y+ of approximately 130 on the hull and 50 on the rudder. The total boundary layer thickness is chosen based on the consideration of stretch factor (which varies between 1.30 and 1.35) and reasonably smooth transition to the core mesh, which is particularly relevant in the stern area. A typical distribution of Wall Y+ on the ship is presented in Figure 10a. As in the case of the propeller, it is impossible to provide a uniform distribution of Y+, but its values are kept above the buffer region (i.e., above 30) everywhere except in the flow separation zones.

The number of prism layers on the deck, superstructure, and cranes is reduced to 6. Because these parts are subject to intensive separation of the air flow with multiple stall areas, the Y+ varies significantly. However, the range of Y+ between 30 and 240 is well preserved, as shown in Figure 10b. The total number of cells in the resistance setup is 15.8 million, which is higher than the usual count in routine resistance calculations. This is primarily due to the inclusion of on-deck features and the finer resolution of propulsor area.

The numerical solution for the free surface is obtained using the VOF method with the Flat VOF Wave model and the blended High-Resolution Interface Capturing (HRIC) scheme [31]. The Courant number limits in the blended HRIC scheme are set to high values (Co_l = 200 and Co_u = 250) to ensure that the HRIC is used irrespective of the time step. These settings mitigate solution dependency on the chosen time step, which is essential when processing the results of resistance and self-propulsion calculations, such as deriving the thrust deduction factor. A representative distribution of Courant number around the ship observed in the resistance simulations is shown in Figure 11.

The Courant number remains below 2.0 for the greatest part of the domain, but it increases to 7.0–8.0 in the regions where the ship’s bow and stern wave systems are generated and to 9.0–10.0 in the area where the stern wave interacts with the rudder. These higher values of Courant number are a natural consequence of the high induced velocities and finer mesh size in these specific areas. In propulsion simulations, due to the small time step ($\Delta$t∼2 deg), the Courant number remains well below 1.0 throughout the domain. Special tests were conducted to verify that the mentioned differences in Courant number have little, if any, impact on the computed ship resistance.

The Dynamic Fluid–Body Interaction (DFBI) model with two degrees of freedom in heave and pitch was used to solve the hydrodynamic position of the ship. The DFBI solver was applied to the whole domain. The equilibrium body motion option was chosen to accelerate convergence towards the sought-after steady-state solution. On the cautionary, it needs to be noted that, depending on the case, the equilibrium solution may lead to large variations in the body’s position at intermediate time steps and introduce disturbances that remain in the numerical solution for a long time. This effect is particularly noticeable in the DFBI solutions conducted without mesh morphing.

The turbulence model investigations primarily focused on two solutions: the RANS method with the k-$\omega$ SST model and the IDDES method with the k-$\omega$ SST model in the RANS domain. Preliminary studies on the case without superstructure also included the LES method, which was found to provide a good prognosis of total resistance of the smooth hull but had larger discrepancies with the RANS and DES solutions in terms of pressure and viscous resistance components. Furthermore, the LES model is not applicable with rough surfaces, which limits its use in the present case where the influence of hull roughness is essential. The results obtained with the SRH model revealed dependency on time step, because as the time step decreases, a larger part of the solution is treated by LES, and eventually when the time step is small (as in the self-propulsion case), the SRH solution is found to be equivalent to that of LES.

Assigning an appropriate value of sand-grain roughness to the ship’s hull required separate investigations. Initially, a surface roughness of 3.244 mm was provided in the JoRes case description from measurements on the hull using an underwater roughness scanner. Such a high reported value of average hull roughness may be caused by the following main factors: vessel age, quality of underwater hull cleaning, and method used for conversion of locally measured values. Furthermore, the available photos of the hull surface revealed numerous patches, pits, and bumps, which can be considered macro-scale roughness measured in millimeters rather than micrometers. At a later stage, JoRes provided another estimation of hull roughness using an alternative approach, resulting in a value of 440 $\mathsf{\mu }$m. Converting the measured value of technical averaged hull roughness (AHR) to the equivalent sand-grain roughness height used in the wall functions is not a straightforward process. Therefore, the value of this parameter was derived from additional calculations with a flat plate model.

In these calculations, the computational domain consisted of a rectangular box of length equal to the ship’s LPP. The bottom surface of the domain was set as the no-slip wall boundary, representing the flat plate. The upstream boundary was designated as the velocity inlet, while the downstream and top boundaries were set as pressure outlets. The side boundaries were imposed with a symmetry boundary condition. The computation mesh was constructed to replicate the prism layer mesh settings used in the ship resistance calculations and provided the same target Wall Y+. At the inlet, a Blasius velocity profile was utilized, which represents the theoretical solution for the streamwise and normal velocities in the laminar boundary layer over a smooth plate. This was done following the recommendation from the JoRes case description [20] in an attempt to mitigate a slight acceleration in the boundary layer that may occur locally when using a constant velocity profile. The Blasius profile applied at the inlet was calculated at a distance of x/L = 0.005 from the plate leading edge. A comparison between the two methods to set up the inlet velocity revealed that local flow acceleration downstream of the inlet is reduced and the displacement velocity decreases somewhat faster along the plate when the Blasius profile is applied. However, only very minor differences between the two solutions were found in the computed distributions of the local friction coefficient and in its integral value. The calculations were performed with a smooth surface and an equivalent roughness height (r_g) of 50, 100, 15, and 200 $\mathsf{\mu }$m. Figure 12 illustrates the computation mesh and the field of streamwise (axial) velocity in the near-wall region of the smooth flat plate.

In

Table 2. Computed values of the integral friction coefficients (C_F) of the flat plate compared with the ITTC friction lines.

	CF × 103
Velocity [knots]	9.0	10.5	12.0
ITTC57 (smooth)	1.6164	1.5851	1.5587
ITTC57 + ITTC78 (rough)	1.6854	1.6795	1.6742
Österlund (exp., smooth)	1.5666	1.5368	1.5116
RANS (smooth)	1.5847	1.5542	1.5287
RANS (rg = 50 $\mathsf{\mu }$m)	1.6314	1.6200	1.6128
RANS (rg = 100 $\mathsf{\mu }$m)	1.7644	1.7659	1.7690
RANS (rg = 150 $\mathsf{\mu }$m)	1.8798	1.8863	1.8922
RANS (rg = 200 $\mathsf{\mu }$m)	1.9746	1.9790	1.9681

, the computed values of the integral friction coefficient (C_F) are compared with the ITTC friction lines used in ship performance prediction procedures and the experimental correlations by Österlund [32]. The calculation results for the smooth plate are approximately 2% below the ITTC57 friction line [33] and 1% above the experimental results by Österlund. The computed streamwise distributions of the local friction coefficient (cf) are presented in Figure 13, where they are compared with different theoretical solutions and experimental correlations. The results show good agreement with Österlund’s data.

Regarding the case of a plate with roughness, the ITTC57 friction line results for C_F corrected with the ITTC78 roughness allowance according to Townsin [34,35] are found to be between the CFD results for the flat plate with sand-grain roughness heights of 50 $\mathsf{\mu }$m and 100 $\mathsf{\mu }$m, as shown in Table 2. Based on these findings, a sand roughness height of 80 $\mathsf{\mu }$m was chosen to be applied in the main resistance and self-propulsion simulations for comparison with the predictions based on model test data. To assess the influence of higher roughness, and considering the high roughness values measured on MV REGAL, additional simulations were performed with r_g = 150 $\mathsf{\mu }$m. The same values of roughness height were applied to the ship’s hull and rudder surfaces.

3.3. Self-Propulsion Simulations

As for the topology of the computation domain, mesh regions, and mesh settings on the ship and propeller, the self-propulsion simulation setup is largely a combination of the setups used in the resistance and open water simulations. Figure 14 provides an illustration of the mesh in the aftship area, where the sliding mesh propeller region is fitted between the ship hull and the rudder. Similar to the open water case, this region accommodates the entire propeller hub and cap. The total number of cells in the self-propulsion setup is 32.8 million, with 13.3 million cells in the main fluid region and 19.5 million cells in the propeller region.

The self-propulsion calculation is performed in several steps. In the first step, the propeller region is fixed, and the solution is performed using the Moving Reference Frame (MRF) method until convergence is attained for the free surface flow. In the second step, the sliding mesh region is set in motion with the initial value of propeller RPM. The RPM is subsequently adjusted to determine the vessel’s self-propulsion point, where the force balance Equation (1) is satisfied:

$$|R_{SP}+\Delta R-T_{PX}|\le \epsilon ,$$

(1)

where $R_{SP}$ is the ship resistance with an operating propeller, $\Delta R$ is the correction accounting for additional resistance components that are not modelled in the simulation (set to 0 in the present case), $TPX$ is the X-component of propeller thrust, and $\epsilon$ is the solution tolerance (set to 0.5% of propeller thrust).

Finally, in the third step, the phase transfer solver with the cavitation model by Scherr-Sauer [36] is activated. The saturation pressure is gradually ramped through the first propeller revolution to its value specified in water properties to address the possible occurrence of cavitation. During the cavitation calculation, it is assumed that thrust-resistance balance is not violated, and therefore propeller RPM is fixed to its value determined in the second step.

In the present study, the ship’s position was fixed in the self-propulsion analysis using the dynamic sinkage and trim values obtained from the DFBI resistance calculation. Both RANS and DES turbulence modeling approaches were investigated, with an equivalent sand grain roughness height of 80 $\mathsf{\mu }$m applied to the ship hull and rudder, and a roughness height of 30 $\mathsf{\mu }$m to the propeller. An additional run with the RANS method using a roughness height of 150 $\mathsf{\mu }$m on the hull and rudder was conducted after the fashion of resistance analyses.

4. Results and Comparisons

In this section, the results of the resistance, open water, and self-propulsion simulations are presented and compared with the experimental data. In all the scenarios, the comparisons are made in full scale, using the results of model tests at SINTEF Ocean extrapolated to full-scale conditions. For the self-propulsion cases, comparisons are complemented by the results of sea trials conducted in the JoRes project. Respecting the current restrictions regarding the data distribution in the JoRes Consortium, the comparisons are presented in the form of relative differences and as plots without scale.

4.1. Results of Open Water Calculations

Figure 15 presents comparisons between the values of propeller thrust coefficient, propeller torque coefficient, and open water efficiency computed using different turbulence modeling approaches and experimental data. The numerical results presented in these figures correspond to a smooth propeller surface. It is important to mention three specific aspects regarding the experimental results. Firstly, in accordance with SINTEF Ocean’s performance prediction procedures, the propeller open water characteristics obtained in model tests are not scaled. Secondly, the model test characteristics consider only the propeller blades, excluding hub thrust and torque. Thirdly, as described in Section 3.1, while the geometry of the propeller blades and hub was exactly the same in the simulations and model tests, different shaft arrangements and hub caps were used. It is mainly the third aspect that explains the differences of 4–5% in propeller thrust and torque in the range of advance coefficient (J) from 0.2 to 0.5. The presence of the dynamometer shaft downstream of the propeller results in an increase of propeller loading by affecting the contraction of propeller slipstream. At the highest J = 0.6, which is already beyond the point of maximum efficiency, the influence of Reynolds number becomes more prominent, explaining the larger difference in propeller torque and hence efficiency. The different turbulence models used in the CFD calculations show good agreement in terms of predicted propeller characteristics. Only slightly larger differences are observed for the LES and SRH solutions, particularly at the highest, J = 0.6.

The influence of surface roughness on propeller performance in open water is illustrated in Figure 16. Only the RANS and DES simulations included the roughness model. Both numerical solutions predict similar magnitudes of the roughness effect, which is mainly seen in the increase of propeller torque and, consequently, the decrease of propeller efficiency. These influences become more significant at higher J values, where the relative contribution of the frictional component is greater. The increase in surface friction is primarily responsible for the observed roughness effect, as is evident from the plots of skin friction distribution on the propeller blade shown in Figure 17. Surface roughness also leads to a reduction in blade section lift, but it is counteracted by local changes in the flow pattern near the trailing edge, where roughness delays flow separation. This explains a minor increase in propeller thrust shown by the calculations including roughness at J = 0.2 to 0.4. Further studies are needed to understand whether this result is entirely physical or partly caused by the modified wall functions in use.

Scale-resolving turbulence methods provide a better basis for the detailed resolution of vortical structures generated by the propeller, even when using a fairly modest mesh, as in the present study. This is illustrated by Figure 18, which displays the field of vorticity magnitude in the propeller slipstream. The DES method is found to improve the resolution of the tip vortex in the near-field and the hub vortex in the large part of the domain. The LES and SRH methods offer superior resolution of tip vortices compared to both the RANS and DES methods in the whole domain. It is important to emphasize that, with the present solution settings, the accuracy of prediction of the integral propeller characteristics is comparable for all methods. Another aspect revealed by the present study is related to the influence of the sliding mesh interfaces. It was found that the presence of a sliding interface downstream of the propeller facilitates premature diffusion of the tip vortex and introduces changes in the dynamic behavior of the hub vortex. While this aspect may be of lesser importance for the calculation of ship propulsion performance, it is certainly relevant to the prediction of cavitation and noise characteristics, as well as unsteady loads on the rudder. The vorticity fields presented in Figure 18 were obtained with the long propeller region, extending 2.5D downstream of the propeller. However, it is impossible to use such a region in self-propulsion simulations with the sliding mesh method. The influence of the sliding mesh region extension on the integral characteristics of the propeller is negligible.

4.2. Results of Resistance Calculations

For the reasons mentioned in Section 3.2, the ship towing resistance simulations were performed only with the RANS and DES methods, including both the cases of smooth hull and hull with an equivalent sand-grain roughness of 80 $\mathsf{\mu }$m. An additional calculation with the RANS method was performed using a higher value of roughness height 150 $\mathsf{\mu }$m. In Figure 19, the total resistance of ship hull, superstructure, and cranes predicted by these simulations is compared to the predictions based on model test data. Figure 20 shows the same comparison as relative differences, in percent, for selected ship speeds.

The calculations conducted with smooth surfaces of hull and rudder underpredict the resistance by approximately 9%. A good agreement (within 1–2% depending on speed) is obtained with both the RANS and DES methods when using the sand-grain roughness height of 80 $\mathsf{\mu }$m, as recommended from the flat plate studies described in Section 3.2. The DES method predicts higher values of resistance, which is due to a larger pressure component of hull resistance. This is explained by the DES capturing more accurately the flow separation at the aftship. The frictional components of hull resistance and rudder resistance predicted by DES are slightly lower compared to the RANS predictions.

Figure 21 presents the results for the aerodynamic resistance of superstructure and cranes. These values were obtained through separate calculations without and with the aforementioned on-deck features. According to both the RANS and DES methods, the superstructure and cranes add approximately 2.5% to the total resistance of the ship. They have a relatively small effect on the dynamic position of the vessel, reducing both sinkage and positive trim (by bow up). The computed values of sinkage and trim are compared with experimental predictions in Figure 22.

While the RANS and DES methods generally provide similar predictions for ship resistance and dynamic position, the choice of turbulence modeling method has a more noticeable impact on the depiction of the nominal wake on the propeller plane, as shown in Figure 23.

As expected, the DES solution results in a heavier wake field (nominal wake fraction 0.518 compared to 0.495 predicted by RANS) with better-resolved separation areas on the portside, starboard, and below the stern tube. The wake images shown in Figure 23 were obtained by time averaging of the wake field data recorded over the last 50 s of simulation time. However, it is in the time-varying wake field where the differences are most prominent. While the wake field in the RANS method hardly shows any changes with time, the wake field resulting from the DES simulation is highly unsteady. The aforementioned separation zones downstream of the shaft tube contain a collection of eddies with varying sizes and intensities. These eddies continuously interact, causing the shape and extent of the separation zones to change over time. To illustrate this, Figure 24 presents instantaneous snapshots of the wake field from the DES simulation taken with a time interval of 5 s.

4.3. Results of Self-Propulsion Calculations

The RPM and shaft delivered power of propeller found, obtained from the self-propulsion calculations using the RANS and DES methods, are compared with the model test predictions and sea trials data provided by the JoRes project in Figure 25. The relative differences for these values, compared to the model test predictions and sea trials data, are shown in Figure 26 and Figure 27, respectively. Similar to the resistance case, both the RANS and DES solutions were computed with equivalent sand-grain roughness on ship hull and rudder surfaces equal to 80 $\mathsf{\mu }$m. An additional calculation with the RANS method was performed using a higher value of roughness height of 150 $\mathsf{\mu }$m. The propeller surface roughness value of 30 $\mathsf{\mu }$m was adopted in these analyses. The presented values are obtained by time averaging over the last five propeller revolutions. The total number of propeller revolutions varied in each case, depending on how fast the self-propulsion condition was achieved. The accuracy of satisfaction of the self-propulsion condition given by the Equation (1) was within 0.2% in all cases.

It can be concluded that both the RANS and DES methods provide comparable prognoses of the ship’s propulsion performance when using a hull roughness of 80 $\mathsf{\mu }$m, which is close to the predictions derived from model test data. The DES results are somewhat closer to model test predictions in terms of shaft power, with relative differences of approximately 2.5% at all ship speeds. The differences in terms of propeller RPM are also within 2.5%.

The RANS results are closer to model test data in terms of RPM (approximately 1%), but they reveal larger deviations in shaft power (2–6% depending on speed). Both the CFD calculations and model test predictions underestimate shaft power and RPM compared to the sea trials data post-processed according to the ISO15016 standard. For the CFD results using the DES method and a hull roughness of 80 $\mathsf{\mu }$m, the underpredictions amount to 3–4% in terms of RPM and 17% in terms of shaft power. The RANS calculations performed with the increased value of hull surface roughness of 150 $\mathsf{\mu }$m reduce these differences to approximately 2% and 12%, respectively.

The high values of hull surface roughness and numerous surface imperfections documented on MV REGAL in the JoRes project indicate that the hull surface conditions may indeed be one of the factors explaining the deviations between the predictions and sea trials. However, in the authors’ opinion, an even more likely reason is related to the presence of bilge keels on the ship hull and sacrificial anodes on both the hull and rudder, which were not considered in the CFD simulations and model test predictions. According to the authors’ experience, the combined contribution of these constructive features can easily amount to 11–13% of increased hull resistance in propulsion conditions. This increase may be even higher in the case of bilge keel misalignment. When combined with the surface conditions of a 25-year-old ship hull, these influences may result in an increased power demand of 15%. Therefore, the differences between the performance predictions and sea trials data found in the present case are not at all surprising and even quite expected.

When comparing these findings with the results of simulations conducted on MV REGAL during the 2016 Workshop [12], it is worth noting that a much closer agreement with sea trials data in terms of shaft delivered power (within 1.5%) was obtained by the authors while applying very similar settings in the CFD model. The 2016 exercises dealt with a different loading condition of the ship corresponding to the draught values of 4.899 m and 5.597 m at the fore and aft marks, instead of 2.97 m and 5.865 m, respectively, applied in the present case. Further, as already mentioned earlier in Section 2, the 2016 simulations used the geometries of ship hull, propeller, and rudder derived from the 3D laser scan data (with some of the surface imperfection and constructive features naturally present), whereas the present simulations used clean CAD geometries reconstructed from drawings. This latter aspect may be critical for direct comparisons with sea trial data.

For the speed of 10.5 knots, the RANS self-propulsion simulation conducted with a hull roughness of 150 $\mathsf{\mu }$m was extended with the cavitation simulation using a phase change model. Fifteen additional propeller revolutions were performed at the same propeller RPM determined from self-propulsion analysis. The cavitation images obtained from the simulation were compared with borescope photography taken during the sea trials. The comparison, depicting selected positions of the blade on the starboard and portside of the ship, is presented in Figure 28.

The experimental images reveal a considerable presence of bubbles swept by the flow from the ship’s bow, under the hull, and onto the propeller. Further, the propeller loading at the observation conditions was higher than that achieved in the CFD simulation. The difference in propeller loading and use of the RANS method presumably explain the smaller extent of tip vortex cavitation shown by the simulations, while an overall cavitation pattern is captured realistically. A better resolution of the tip vortex flow would be achieved with the DES method, but the DES simulations were not performed for the case of increased hull roughness at the time of preparation of the present manuscript.

5. Conclusions

To calibrate the existing modeling practices for ship self-propulsion analyses, CFD calculations at full scale have been conducted on the benchmark ship MV REGAL, after the fashion of a “blind” simulation exercise. Comparisons have been made with the experimental predictions based on model tests and the results of sea trials performed within the joint industry project JoRes.

The analysis of comparative results demonstrates that the applied CFD modeling practices are mature and capable of predicting a ship’s performance characteristics with the accuracy required for practical applications. In particular, the results of the DES method applied in this paper are found to be in good agreement with the prognosis based on model tests. In this case, the differences in terms of shaft delivered power and propeller RPM do not exceed 2.5%, which is well within the accepted range of CFD/EFD calibration factor (0.95 to 1.05), according to [3]. The DES method provides a good compromise between the computational demand, accuracy of prediction of propeller and hull forces, and resolution of flow details. It also supports the inclusion of surface roughness and shows little sensitivity to the simulation time step.

Both the CFD calculations and model test prognosis underpredict shaft delivered power by 12 to 17% and propeller RPM by 2 to 4%, depending on the applied value of equivalent hull roughness, when compared to the sea trials data. These differences are thought to be caused by the additional resistance of bilge keels and sacrificial anodes, which were not included in the numerical model, as well as the surface conditions of the ship’s hull. The use of accurately scanned hull, propeller, and rudder geometries is therefore deemed highly important for direct comparisons between full-scale CFD predictions and sea trial data on old ships in service. The accuracy of roughness measurements, conversion of the measured values to equivalent sand-grain roughness height applied in CFD, and roughness distribution pattern on the hull require closer investigations concerning both their impacts on hull resistance and wake field on the propeller.