1. Introduction

Small ships take a portion up to 46% of the whole population of ships (Equasis, 2020), which can be concluded to contribute relatively big amounts of CO₂ emission gases. This is against the IMO goals to reduce the emissions of shipping industries by 50% in 2050 (IMO 2018), (MARPOL, 2020). Greenhouse gas emissions are mainly produced by the engine and propulsion system of the vessel. Therefore, many efforts have been made to solve this problem. The idea is to reduce the vessel resistance by having a better hull and superstructure, improving the efficiency of the propulsor by choosing better design parameters, hull and propulsion system interaction optimization, and also optimizing the strategy for operation (Molland et al., 2014). By operating in the optimum conditions, several improvements could be made to the ship’s performance. One of the aspects of the operation is the trim angle. This issue leads to many studies to optimize the ship’s trim angle which could bring a significant impact to reduce the vessel resistance (Reichel et al., 2014). Energy Saving Devices (ESD) are commonly used to reduce lower propulsive requirements (ITTC, 1999), also trim optimizer devices were studied and proven to reduce the ship’s resistance (Seo, 2010), (Seo et al., 2013). The ship’s trim is majorly affected by the weight distribution of the ship along the longitudinal axis. Therefore, in this study, changes in resistance performance due to the Longitudinal Center of Gravity (LCG) position are studied by empirical and numerical methods.

2. Methodology

In this study, two different speeds are analyzed under 8 different Longitudinal Centers of Gravity (LCG). The empirical analysis is performed using the Savitsky method for planing hulls (Savitsky, 1964), with the aid of Maxsurf Resistance software. While the numerical analysis is performed using the commercial Computational Fluid Dynamics (CFD) software of Star-CCM+. The comparison is made for both analyses and discussed further.

2.1 Target Vessel

A hard chine and flat-roof tunnel stern of a small planing hull form is used as the research vessel. A final 1.972 m length of planing ship corresponding to a 1:9 scale is used to minimize the computational time. The target vessel's hull form and main dimensions are presented in Fig.1 and Table 1 below.

2.2 Study Cases

The variations are made to the Longitudinal Center of Gravity (LCG) and the speed range. The LCG is set within the range of 29.2% - 49.1% of the ship’s overall length (LoA), while the speed range is set as the considered design speed, ranging from Froude Number 0.56 - 0.63. Table 2 will summarize the scope of this study. A total of 16 cases are made to evaluate the resistance performance along with the LCG changes.

3. Empirical Analysis

There are many theoretical and empirical studies of the planing hulls. These studies are mainly focused on the hydrodynamic factors such as hydrodynamic lift, drag, wetted surface area, pressure distributions, impact forces, wake form, spray formation, and dynamic stability are considered to develop the planing hulls (Savitsky, 1964). A summary report on the completed studies was published to present the hydrodynamic lift, drag, wetted surface area, and center of pressure (Korvin-Kroukovsky et al., 1949). Where finally, by Savitsky in 1954, a set of empirical planing equations was developed (Savitsky and Neidinger, 1954). Finally, in 1964, a study to formulate a simple computational procedure or empirical analysis to predict the power requirements of planing hulls by estimating the total resistance, followed by the porpoising stability of the vessel was finished. These studies include technical studies previously done, which were also completed with experimental results, resulting in an empirical formula to evaluate the ship’s total resistance which is commonly known as the Savitsky empirical approach (Savitsky, 1964) for planing hulls.

3.1 Savitsky Empirical Formula and Applicability

The total resistance mainly consists of two main components, namely D_p, or the bottom resistance component due to pressure force, and D_f, or the viscous drag acting tangential to the bottom. Finally, the empirical formula of the resistance estimation is shown in the equation below as follows:

Where Δ is the ship displacement, τ is the trim angle, followed by the Schoenherr turbulent coefficient of C_f, V₁ as the average velocity measured in the bottom section, mean wetted length-beam ratio of λ, b for the breadth of the hull, and the deadrise angle of β.

In the planing regime, a hydrodynamic lift is produced to support the hull displacement. Therefore, an empirical formula is also provided to estimate the lift produced by a given trim and draught, also along with the deadrise angle influence. As the function of beam speed coefficient (C_υ), length-beam ratio (λ), and trim angle (τ), an empirical formula of hydrodynamic lift coefficients for planing hull form with some deadrise angle is presented as follows:

There are some limitations of the empirical formula, especially the geometry of the vessel and speed range. The formula is applied by assuming that the hull form is having a constant deadrise, trim angle, and forward speed. In addition, the empirical formula is also developed based on hard-chine planing hull-form vessels. With this formula, the non-dimensional beam speed coefficient (C_υ) is used to evaluate the limitations. The reason for this approach can be justified since the planing ships’ length (L_WL) varies with trim angle (τ), and speed (V), while the beam (b) is mostly constant. This empirical formula mostly covers planing hull form which has a beam speed coefficient in the range of 0.6 - 1.3. (0.6 ≤ C_υ ≤ 13), deadrise angle of up to 30 degrees (β ≤ 30^o), trim angle within the range of 2 - 15 degrees (2^o ≤ τ ≤ 15^o), and a wetted length-beam ratio of less than 4 (λ ≤ 4).

3.2 Empirical Analysis Results

The total resistance estimated by the Savitsky approach is calculated using the aid of Maxsurf Resistance software (Bentley, 2022). These calculations are done by following the manual of the software to ensure its reliability. The calculations are started by importing the 3D geometry of the small planing ship, which previously has been set into a given draught. After all surfaces are measured, the next step is the speed definition, which in this study is set within the range of 0.56 - 0.63. Savitsky planing method is used since it covers the goal of this study by providing the LCG as one of the parameters to be input. These LCG changes were input to the hydrostatic parameters windows. The LCG input in the windows is used by Maxsurf within the algorithm.

The initial trim angle is assumed, and the equilibrium trim angle is interpolated between the two initial trim that have been chosen. Where the displacement acted as a negative force, and the D_f or viscous component of drag, parallel to the keel line acted as a positive value. The equilibrium condition is defined as the zero-sum of these two forces. The resistance is calculated on the equilibrium trim condition, and the result is presented in the table and plotted in the graph.

By the Savitsky approach, it is agreed for both cases that lower resistance can be achieved by shifting LCG towards the bow section, compared to the stern section. The results give a reasonable result with a low computational cost and time, which takes around 30 minutes - 1 hour.

Resulting in a maximum of total resistance reduction of up to 53.29% for Froude Number 0.56 and up to 51.29% for Froude Number 0.63. The total resistance estimated by the Savitsky approach tends to keep going lower within the considered cases. The optimum LCG has achieved on LCG around 49.1% LoA from the transom for Froude Number 0.56 and Froude Number 0.63. The resistance is within a huge range of 23.750 - 50.850 N for the Froude Number 0.56, and 27.160 - 55.760 N for the Froude Number 0.63. The results are summarized in Table 3 and presented in Fig. 2.

4. Numerical Analysis

With the development of computational analysis, a numerical simulation could be performed in order to simulate the towing test experiment conditions. The idea is to reduce the computational cost and time for the experiment test which is limited to small changes in the ship design process. Therefore, numerical analysis to estimate the total resistance of small planing ships is performed to complete this study. In this study, a commercial Computational Fluid Dynamics (CFD）software namely Simcenter Star-CCM+ code is used.

4.1 Methodology

The methodology used in this study is based on recent studies for resistance evaluation. Implicit unsteady time-step is applied along with Reynolds-averaged Navier-Stokes (RANS) equation as the governing equations for the incompressible flow, and viscous fluids in a finite volume within a time control, which represents in the equation (4), (5), and (6).

The turbulence flow is simulated by the SST-κ - ω, which is commonly used for ship hydrodynamics simulations, due to the limitations of the k-ε model which may not be appropriately modeled the flow separation (Pena et al., 2020), which also confirmed by (Pena and Huang, 2021) that the SST-κ - ω turbulence model is a solid turbulence modeling strategy which capable to model pressure gradients and flow separation. Ship resistance by numerical simulations was studied (Elkafas et al., 2019) by using the SST-κ - ω model gives high-accuracy modeling of the boundary layer and gives accurate predictions of the amount of flow separation where turbulence is present, making it the best in marine engineering applications. The use of the SST-κ - ω model for resistance calculation (Pacuraru et al., 2022) (Zha and Ye, 2014) is satisfactory which is indicated by the accuracy of numerical computation.

The following equations (7) and (8) are used for the turbulence model.

To separate the multi-fluid phase of water and fluid at the free-surface area, a volume of fluid (VOF) method is used, along with High-Resolution Interface Capturing (HRIC). The light fluid of air and heavy fluid of fresh water is applied for the upper and below waterline region, respectively. The air-water interface is defined by a volume fraction of F, between the range of 0 to 1 (0 < F < 1), or within the time range that could be represented in the equation (9).

In order to see the effect of LCG changes, the small planing ship motion in the calm water needs to be modeled especially for longitudinal motion, namely heave and pitch. To perform the motions, Dynamic Fluid Body Interaction (DFBI) method is used, which could simulate the interaction between the ship’s hull and the multiphase fluid of water and air. The LCG changes could be set within the simulations under these settings.

The methodology listed in Table 4 needs to be verified before being used in the targeted ship of small planing ship. In order to do the verification, a simple geometry of planing hull form is used. The planing hull used is set into the numerical simulation as close as possible to the main goal of the targeted vessel of the small planing ship. The principal dimension of the simple planing hull (Warped-2) is listed in Table 5.

The Warped-hull 2 is generally a single V-type planing hull that experimentally towing tested in the Tank of the University of Naples Federico II (Begovic and Bertorello, 2012). With a provided towing test result and hull geometry, validation work is performed to see the appropriateness of the methods.

The results are shown in the Fig. 3 and 4 below, where it can be concluded that the methods are applicable to be performed for this study due to the reasonable accuracy with average difference about 5.80% for CFD simulation compared to the towing test results.

4.2 Computational Domain and Boundaries

After a proper methodology is chosen, a virtual towing tank is made to evaluate the resistance performance of the target vessel, considering ITTC recommendations (ITTC, 2014). The virtual towing tank is set with the small planing ship length overall (L) as the parameters. A final computational domain and boundaries are presented and summarized in the Fig. 5 and Table 6.

4.3 Mesh Generation

Mesh configuration is one of the most important factors in numerical analysis, due to its sensitivity to the simulation result while considering the computational time. Therefore, the mesh study is first performed for the y+ value to capture the flow interaction in the boundary layer. All y+ treatment method is used by having 15 layers of prism mesh and a final y+ value of greater than 40 (y⁺ > 40) is set considering the convergence of resistance and other hydrodynamic parameters such as trim and wetted surface area. By maintaining the targeted y+ value, the mesh configurations are now set to various sizes to get a proper size to ensure accuracy as well as to have an optimized computational time. Richardson’s Extrapolation method (Richardson, 1927) is used to get the proper mesh configurations. The time step is chosen at 0.001 s to satisfy the CFL numbers of less than 1.

Three mesh configurations namely coarse, medium, and fine mesh is produced, and the results are summarized in Table 6. By having a reasonable Grid Convergence Index (GCI). Grid uncertainty analysis is a method of estimating the virtual solution when the grid space is zero by extrapolating the results of several grid spacing (Roache, 1994). A medium-mesh configuration is chosen in this study, since it gives a reasonable GCI as summarized in Table 5. Volumetric control regions are added for the free surface area to capture the wake and spray formation better. A final mesh consisting of approximately 2.57 million cells is used in this study as shown in Fig. 6.

4.4 Numerical Analysis Results

The total resistance estimated by numerical analysis tends to keep going lower within the considered cases, where at some point, an optimum LCG is achieved. The optimum LCG is chosen based on the lowest resistance that could be achieved in the considered Froude Number.

As shown in Fig. 7, the optimum LCG has achieved on LCG around 46.2% LoA from the transom for Froude Number 0.56, while around 43.4% for Froude Number 0.63. By numerical analysis of CFD, at some point, the resistance trend starts to increase. Noting that 29.2% LoA as a reference LCG point, a maximum of total resistance reduction up to 45.16% and 41.12% for Froude Number 0.56 and 0.63 is achieved within the LCG changes range as summarized in Table 8. The analysis shows a sensitivity of the simulations with a relatively longer computational time of 5 - 6 days.

With the CFD simulations, the total resistance is deeply analyzed by considering the shear and pressure force components. Within the resistance reduction of up to 41.12% - 45.16%, the main reduction is concluded as the result of pressure force reduction within the LCG changes. These results are well proven by separating the total resistance into two main components, namely shear and pressure force. Despite the fact that the resistance is greatly reduced along the LCG changes, the shear force is slightly increased. However, the increase is greatly countered by the significant reduction of the pressure force.

Table 9 summarizes the analysis, where can be seen that the shear force is increased by about 0.34 - 2.10 N compared to case (1) for Froude Number 0.56, and increased by about 0.65 N - 3.23 N compared to case (10) for Froude Number 0.63. However, the pressure force is greatly reduced by 8.68 N - 23.72 N and 8.81 N - 23.20 N for the equivalent cases, respectively. So it can be concluded that the total resistance reduction is dominated by the reduction of the pressure force along the LCG changes.

5. Results and Discussion

Both empirical and numerical analyses are performed to evaluate the effect of LCG changes on the total resistance performance of the small planing ship. The results are then compared to see the reliability of each method. The figure below represents both results for either CFD and Savitsky empirical formula, for each Froude Number considered.

From the plot presented in Fig.8, and the results summarized in Table 10, it can be concluded that for all the test conditions, the CFD analysis gives a lower total resistance estimation compared to the Savitsky empirical formula. It can be seen that the triangle and circle marked graphs (CFD results) are almost lower in call cases compared to the square and x-marked graphs (Savitsky approach) for both speeds considered. The first 5 cases of Froude Number 0.56, and 0.63 with LCG points ranging from 29.2% - 40.6% LoA shows a similar trend of resistance with a relatively small differences ranging from 0.85% - 7.41%.

These results gap keeps getting wider for the rest test cases. Shifting the LCG further to the bow section gives as high as 16.58% differences occurred in the analysis. For each speed considered, the lowest resistance is achieved in case 7 and case 14 by the CFD, but the resistance keeps going lower even until case 8 and 16 for the Savitsky method. This anomaly is then investigated, and two factors are believed to be the supporting arguments for this phenomenon could happen.

First, one of the causes would be the limitations of the Savitsky empirical formula, especially in the hydrodynamic lift estimations. The hydrodynamic lift of equation (3) is derived from the lift of flat planing surfaces in equation (2) to include the influence of deadrise angle by a regression of the previous studies, where the vessels mostly have a length-beam ratio lower than 4 (λ ≤ 4). LCG located at the stern section, tends to reduce the length of the waterline because of the trim by stern occurred, which lifts the bow section of the vessel, making the target vessel within the range of the length-beam ratio of the vessels used to derive the formula. However, shifting the LCG towards the bow section will somehow increase the length of the waterline because of deeper immersion of the bow section, which in other words, tends to break the limit of the length-beam ratio of the Savitsky empirical formula, making it gives a higher hydrodynamic lift coefficient values which in fact, reduce the total resistance performance.

Second, in this study, the LCG changes will affect the given trim of the vessel, which by equation (1) could give a variation of the total resistance. Using the Maxsurf Resistance to help the calculations, the changes of LCG area input into the software by changing the hydrostatic data of LCG, while maintaining other parameters to be the same within the interface. However, in the real conditions of the vessel’s operation, changing LCG is not only affect the trim conditions, but also the displacement, length of the waterline, or even beam and deadrise angle (if the geometry varies along the hull length). These hydrostatic parameters could be one of the reasons the resistance keeps getting lower since the main parameter that changes from equation (1) is only the trim angle, while in fact, most of the variables are subject to change.

6. Conclusion

This study concludes the Longitudinal Center of Gravity (LCG) changes to the total resistance performance of the small planing ship by empirical and numerical methods. As a result, both methods answer the sensitivity of trim conditions for small planing ship to the resistance. By operating in an optimum trim condition, several improvements to reduce resistance which resulting lower power requirements could be achieved.

Both analyses gives a reasonable sensitivity considering the computational time. A total resistance reduction of 45.16% is achieved by operating in the optimum trim conditions at the Froude Number 0.56 based on the CFD analysis which takes about 5 - 6 days to complete. While with other methods of Savitsky empirical formula which has relatively lower computational time of 10 - 30 minutes, the resistance trend keeps going lower, as low as 53.29% total resistance reduction in the Froude Number 0.56.

With the CFD simulations, no significant differences in shear force are achieved within all the LCG changes. Concluding it the main differences of the resistance are achieved due to the reduction of the pressure resistance. Even though a different trim condition results in different wetted surface areas and lengths of the waterline, in fact, the shear force analysis shows no significant differences in shear force between all the cases, and the pressure resistance is the main factor that is changed in the vessel's resistance performance. This might happened due to the different positions of the center of pressure, and also the wave height produced under different trim conditions. Resulting in various resistance in every LCG position conducted in this study.