1. Introduction
The significance of underwater radiated noise (URN) is increasingly emphasized within the submarine industry. The URN generated by propulsors is recognized as a critical factor affecting the overall marine environment besides being an engineering performance issue. Previous studies have reported that the dominant frequency components radiated from ships may overlap with the communication and detection frequency bands of marine mammals and fish. Such interference can potentially degrade communication capability, cause changes in migration routes, and lead to disturbances in spawning behavior. Furthermore, if these effects accumulate over a long period, they can disrupt the equilibrium of marine ecosystems. Consequently, studies characterizing noise under cavitating and noncavitating conditions and those developing prediction methods are being actively conducted. For reliable noise prediction, prior validation is required to assess the accuracy of computational fluid dynamics (CFD) in reproducing the complex flow around a submarine.
This study aims to validate the propeller open-water (POW) characteristics and self-propulsion performance of an existing submarine geometry to establish a CFD method that can be applied for future analyses of submarine noise and cavitation. The DARPA SUBOFF submarine, introduced by Grove et al. (1989), was selected as the target for this analysis. SUBOFF is widely utilized as a benchmark geometry for comparing and validating the prediction performance of resistance and flow fields across turbulence models and it is supported by the relatively well established model test and local flow measurement data. Representative SUBOFF model test data include the resistance test results of Crook (1990) and the maneuvering model test results of Roddy (1990). In addition, Huang & Liu (1994) provided local flow measurement data around SUBOFF, establishing a basis for validating numerical results in terms of the resistance and flow field.
Thus far, several numerical studies have reported on self-propulsion under SUBOFF and propulsor conditions. For example, Byeon et al. (2018) performed resistance and self-propulsion simulations for the SUBOFF submarine with the INSEAN E1619 propeller, comparing and analyzing the effects of turbulence model selection on the prediction of resistance and self-propulsion performance. Meanwhile, Chase (2012) conducted high-resolution numerical simulations for the SUBOFF submarine and its propulsor to analyze the self-propulsion performance and validate the applicability of the self-propulsion analysis procedure through comparison with towing tank experiments. Sezen et al. (2018) predicted the resistance and self-propulsion performance of SUBOFF using Reynolds-averaged Navier–Stokes (RANS)-based simulations and examined the applicability and limitations of self-propulsion analysis methods by comparing the virtual disk approach with actual propeller modeling results. In addition, Alin et al. (2010) compared and discussed the RANS, detached eddy simulation (DES), and large eddy simulation (LES) results for SUBOFF flow and evaluated the applicability of DES and LES under limited grid conditions.
Meanwhile, Park et al. (2019) numerically analyzed the resistance and self-propulsion performance of the SUBOFF submarine model under large cavitation tunnel (LCT) conditions and compared the results with the KRISO LCT model test data. The study established a cavitation analysis procedure for the submarine–propulsor system in the LCT environment. In addition, studies addressing the hull–appendage–propulsor interaction and noise prediction have been reported. For example, Choi et al. (2019) performed numerical flow simulations considering hull–appendage–propulsor interaction, and they predicted propulsor radiated noise under cavitating and noncavitating conditions based on the simulation results. Seo & Park (2021) conducted DES-based simulations to compare and analyze flow fields under varying LCT model test conditions.
In this study, based on the POW characteristics of INSEAN E1619 and the self-propulsion performance of SUBOFF, the results obtained using RANS and LES are compared to quantitatively assess the difference in prediction accuracy depending on the turbulence model. In the RANS simulations, the representative eddy viscosity-based shear stress transport (SST) k–ω model and the Reynolds stress model (RSM), which account for anisotropic turbulence effects, are applied to examine their influences on POW and self-propulsion predictions. Furthermore, through comparison with LES results, the impact of resolving unsteady turbulent structures on the prediction of resistance and self-propulsion performance is discussed. The findings of this study are expected to serve as fundamental data for selecting appropriate turbulence models and analysis strategies in future studies on flow-induced phenomena such as submarine noise and cavitation.
2. Numerical Method
Numerical simulations in this study were performed using STAR-CCM+ version 18. The governing equations are the incompressible Navier–Stokes equations discretized using the finite volume method. Pressure–velocity coupling was addressed using the semi-implicit method for pressure-linked equations algorithm for both steady and unsteady calculations. Time integration was performed using a second-order accurate scheme, and the convection terms were computed using a second-order spatial discretization scheme.
Under model test conditions, free-surface effects were excluded in both the POW and self-propulsion analyses. The POW performance analysis was conducted under steady and unsteady flow conditions, while the submarine self-propulsion analysis was performed only under unsteady flow conditions. The moving reference frame (MRF) approach was applied for the steady flow calculations, and the sliding mesh technique was used for the unsteady flow calculations. RANS-based steady (MRF) simulations were executed for 2000 iterations. For the sliding mesh simulations, unsteady calculations were performed after the flow field converged under steady conditions using a time step corresponding to a propeller rotation angle of 2°. LES calculations were initially conducted with a time step corresponding to a rotation angle of 20°, which was gradually reduced to a final time step corresponding to 0.5° considering numerical stability and convergence characteristics. This time step setting was adopted to maintain the CFL number (Courant-Friedrichs-Lewy number) below 1 based on the grid resolution.
3. Propeller Open-Water Characteristics
3.1 Target Propeller
The seven-bladed E1619 propeller designed by INSEAN in 2009 is used in this study. The principal particulars of the E1619 propeller are summarized in Table 1. The diameters (D) of the propellers used in the computations are 0.485 m and 0.262 m, and the hub ratio is 0.226. In addition, the pitch ratio (P/D) at 0.7R is 1.150, and the chord lengths at 0.75R are 0.0068 m (D = 0.485 m) and 0.0037 m (D = 0.262 m). The extended area ratio is 0.608. The numerical results were validated by comparison with the thrust coefficient (KT), torque coefficient (KQ), and open water efficiency (η0) measured in the model tests. Fig. 1 shows the geometry of the E1619 propeller.
3.2 Boundary Conditions
The computational domain and boundary conditions for the POW performance analysis are shown in Fig. 2. The domain was configured as a cylinder, and a velocity inlet condition was applied at the inflow boundary with a diameter of 16D. A pressure outlet condition was applied at the outflow boundary, and its diameter was maintained at 16D, identical to that of the inflow boundary. The axial length between the inlet and outlet boundaries was set to 20D. The cylindrical side surface was treated as a symmetry boundary. The computational domain was divided into a rotating region containing the propeller and a stationary region outside it.
3.3 Turbulence Model and Wall Treatment
For turbulence modeling, SST k–ω and RSM were applied in the RANS simulations, while a subgrid-scale model was used for the LES simulations. In all computations, the near-wall grid resolution was set such that y+ remained below 5 to directly resolve the wall-bounded flow. Steady RANS simulations were performed using the MRF approach; the INSEAN model test conditions, wherein the inflow velocity is fixed and the rotational speed is varied, were applied directly. In this case, the Reynolds number (Re) varies significantly with changes in rotational speed, and therefore, y+ ranged from 1 to 4. In contrast, unsteady RANS (URANS) simulations were conducted using the sliding mesh technique with the same grid system, and y+ was maintained below 1.
3.4 Grid System and GCI
The grid system comprised three levels: coarse, medium, and fine. Grid convergence was verified using the SST k–ω turbulence model. A polyhedral mesh generation method was applied in the rotating region to improve the geometric fidelity of the propeller, while a Trimmer mesh was used in the stationary region. The rotating region is indicated in red in Fig. 3.
The grid convergence index (GCI) was evaluated through RANS simulations with the MRF approach to assess grid convergence at J = 0.2081, 0.4625, 0.6458, 0.7400, 0.8005, and 0.9661. The GCI calculation was performed following the discretization error estimation procedure proposed by Celik et al. (2008), and the INSEAN model test results were referenced for comparison and validation. The GCI for the grid is formulated as
where Q1, Q2, and Q3 represent the values calculated on the fine, medium, and coarse grids, respectively. The grid refinement ratio r is set to
2
The convergence index p is calculated as
R is defined as
The GCI results, number of grid cells, J, KT, and KQ are summarized in Table 2. The R(G) values for both KT and KQ corresponded to monotonic convergence. In addition, the numerical uncertainties of KT and KQ were evaluated to be within 2% for all J conditions, confirming that the grid system used in this study achieved sufficient grid convergence.
Fig. 4 compares the GCI results with the INSEAN model test results and the results of Byeon et al. (2018)
3.5 Results According to the Calculation Method
The GCI results confirmed that sufficiently accurate results can be obtained even with the medium grid system. The POW performance analysis was conducted using the MRF approach, and therefore, additional calculations are required to examine the difference from the sliding mesh method. In the URANS simulations, the time step was set such that a propeller rotation angle of 2° corresponds to one step for each rotational speed condition. Tables 3 and 4 compare of the results from the two calculation methods. The comparison shows within no significant difference in performance prediction was observed between the MRF approach and the sliding mesh method in the range of analysis conditions considered in this study. This indicates that the calculation method has a negligible effect on the results. Consequently, subsequent simulations with different turbulence models were also performed using the MRF approach.
These results are compared with the INSEAN model test results and those of Byeon et al. (2018) in Fig. 5.
3.6 Results According to the Turbulence Model
The numerical results are summarized in Tables 5 and 6; the CFD results are presented graphically in Fig. 6 for comparison with the INSEAN model test results and those of Byeon et al. (2018). Byeon et al. (2018) used approximately 15.0M grid cells in their study. The comparison results revealed that the SST k–ω model exhibited errors within 7% for KT and within 5% for KQ and η0, while the RSM showed differences within 5% for KT, KQ, and η0 compared with the experimental results. Therefore, under the conditions of this study, RSM produces smaller discrepancies with the INSEAN experimental results. In addition, compared with Byeon et al.’s (2018) results, which were obtained with the same turbulence models and numerical methods, RSM showed a tendency to better follow the experimental results, particularly in the high advance ratio region.
For the experimental results, all turbulence models predicted relatively lower thrust compared to the torque over the entire range of advance ratios, and this causes differences in η0. Both SST k–ω and RSM showed a trend in which KT and KQ decrease at a nearly constant rate with an increasing advance ratio, similar to the experimental results. Byeon et al. (2018) showed a similar decreasing trend; however, in the RSM results of this study, the reduction in KT and KQ over the range from J = 0.8005 to 0.9661 was larger than that reported by Byeon et al. (2018).
The same boundary conditions as in the RANS simulations were applied to the LES, and approximately 80.0M grid cells were used. This increase in the total number of grid cells was intended to resolve large-scale vortex structures inherent to LES. In addition, the near-wall grid resolution was configured so that y+ remained less than or equal to 1. Unsteady simulations were performed because steady-state analysis is not suitable for LES; the operating conditions were set with reference to the INSEAN model test conditions. The numerical results are summarized in Table 7 and presented graphically in Fig. 7. For comparison, the results were compared not only with the experimental data but also with those of Chase (2012) based on CFDShip-Iowa. Chase (2012) used 73.4M grid cells with y+ = 1 and applied the delayed DES (DDES) turbulence model.
The LES results showed an overall difference of ~ 5% compared to the experimental data. Except for the condition J = 0.2081, KT and KQ were predicted to be lower than the experimental values for most advance ratios. Compared to the results of Chase (2012), based on J ≈ 0.8, the LES results showed relatively lower and higher values in the low and high advance ratio regions, respectively. The LES results exhibited a thrust-to-torque ratio similar to the experimental data over the entire range of advance ratios, resulting in a negligible difference in η0.
The MRF analysis performed in RANS for the propeller with a diameter of 0.485 m applied the same conditions as the INSEAN model test; however, the Reynolds number (Re) varied within the range of 3.19E+4 to 1.36E+5. Therefore, even with the same grid system, the y+ value can vary significantly depending on the operating conditions. In the self-propulsion test, a propeller with a diameter of 0.262 m was used and the sliding mesh method was applied. The propeller with a diameter of 0.485 m propeller was scaled down to the 0.262 m condition and the grid system was constructed to have the same number of grid cells as in the original case to examine the prediction characteristics attributed to changes in scale and numerical method. In addition, to minimize variations in Re, the rotational speed was fixed at 50 RPS, and numerical simulations were conducted within the range of Re = 1.19E+5 to 1.29E+5.
The sliding mesh approach was used as the numerical method, consistent with the self-propulsion simulation, and the time step was set to correspond to a propeller rotation angle of 2° for obtaining converged results. The results are summarized in Tables 8 and 9 and compared with those of Byeon et al. (2018) in Fig. 8. The comparison results show that, in the range of J = 0.2081–0.8005, all turbulence models exhibit differences of ~5% compared to the experimental data, whereas at J = 0.9661, the difference increases to ~10%. This can be attributed to the rate of thrust reduction decreasing while the rate of torque reduction increasing in the MRF analysis with an increase in the advance ratio, which in turn results in a relatively larger discrepancy. Consequently, the error tends to increase as the reduction in torque becomes more pronounced compared to the experimental data. However, in the sliding mesh results, the ratio of thrust to torque gradually increased with the increasing advance ratio, and unlike the MRF results, the error in η0 became closer to the experimental values in the high advance ratio region.
3.7 Flow Field Comparison
A comparison of flow field differences according to turbulence models and calculation methods is presented below. Fig. 9 shows a graph comparing the overall POW results.
Fig. 10 compares the axial velocity distribution measured at x/R = 0.17 downstream of the propeller plane under the condition J = 0.7400 using laser doppler velocimeter (LDV) in model tests conducted by Liefvendahl et al. (2010) along with the numerical results of Chase (2012), Byeon et al. (2018), and the present study. The results of Chase (2012), which applied DDES, confirm that the velocity distribution induced by the blades and the contrast between high- and low-velocity regions formed by the tip vortex are reproduced similar to that in the model test measurements. Further, a similar tendency is observed in the LES results of the present study.
In contrast, the RANS-based results show more pronounced low-velocity regions. Among the RANS models, RSM exhibits a relatively stronger contrast between the low- and high-velocity regions compared to the SST k–ω model. This results in a more pronounced velocity distribution. Further, compared to the case computed using the sliding mesh method at D = 0.262 m, the case computed using the MRF approach at 0.485 m shows more distinct low-velocity regions consistently in the RANS models (SST k–ω and RSM). These characteristics of the RANS distributions have been reported in the RSM-based RANS results of Byeon et al. (2018).
Fig. 11 compares vortex structures according to turbulence models and calculation methods at the advance ratio J = 0.7400. A Q-criterion value of 10,000 was applied to visualize the vortex structures. In the RANS simulations using the MRF approach, the grid resolution in the rotating region is relatively lower than that in the sliding mesh simulations, and thus, both SST k–ω and RSM show smaller vortex sizes and shorter vortex lengths compared to those of the sliding mesh results. In contrast, when the same grid system is applied, similar vortex structures are observed regardless of the calculation method. Further, although the grid system was nondimensionalized using the propeller diameter D and set to be identical, the absolute grid size is smaller in the case of D = 0.262 m than in the case of D = 0.485 m, which enables relatively finer vortex structures to be resolved.
In the LES results, the vortices are formed with longer structures, while the thickness of the vortex core is relatively thinner compared to that in the RANS results. Further, the strongest vortices are observed near the blade tip region where the low-velocity regions are formed, and this is consistent with all RANS models. However, in the LES results, continuous vortex structures are additionally observed not only near the tip region but also in the vicinity of the trailing edge.
4. Self-Propulsion Performance
4.1 Target Submarine Model
The DARPA SUBOFF 5470 hull equipped with the INSEAN E1619 propeller is used as the target submarine model to perform self-propulsion simulations. The model is shown in Fig. 12, and the principal particulars are summarized in Table 10. The main particulars of the DARPA SUBOFF 5470 model are as follows: length overall (LOA) of 4.356 m and length between perpendiculars (LPP) of 4.261 m. The maximum hull radius (Rmax) is 0.254 m, displacement (∇) is 0.718 m3, and wetted surface area (WSA) is 6.338 m2. These are identical to the specifications used in the resistance analysis. In addition, the propeller radius was set to 0.262 m to maintain scale consistency.
4.2 Boundary Conditions and Grid System
The boundary conditions and computational domain for the self-propulsion analysis are as follows. A velocity inlet condition was applied 2.0LPP upstream of the hull, and a pressure outlet condition was applied 3.0LPP downstream of the hull. The remaining boundaries were treated as symmetry conditions, and the lateral, upper, and lower boundaries were placed at a distance of 1.5LPP from the hull.
For the grid system, a polyhedral mesh was applied in the rotating region, while the remaining regions were constructed using the Trimmer mesh generation method. In addition, the grid system validated through GCI verification in the resistance and POW analyses was applied identically to both the rotating and stationary regions. A total of 10.2M cells (5.2M and 5.0M cells in the stationary and rotation regions, respectively) were used for the RANS models (SST k–ω and RSM), while a total of 139.2M cells were used (86.1M and 53.1M cells in the stationary and rotation regions, respectively) for the LES model.
4.3 Self-Propulsion Results
For the RANS models, the grid was configured such that y+ remained below 1 over the entire domain, and the simulations were performed accordingly. For the LES model, a wall-modeled LES condition was applied to the submarine hull with y+ ≥ 50, while the propeller was resolved directly with y+ ≤ 1. This setup was applied with reference to the self-propulsion simulation conditions used in An et al. (2025) (hull y+ ≥ 30 and propeller y+ ≤ 1).
The self-propulsion simulations for RANS were initially performed under steady-state conditions using the MRF approach to sufficiently converge the flow field. Subsequently, they were switched to an unsteady analysis using the sliding mesh method in which the propeller generated by the grid was directly rotated. For LES, the time step based on the propeller rotation angle was gradually reduced in the sequence of 20°, 10°, 5°, 2°, 1°, and 0.5° to ensure convergence stability in the initial transient region. The final time step was set based on the same criterion as the POW analysis, which corresponded to a propeller rotation angle of 2° for RANS and 0.5° for LES. The simulations were performed until sufficient convergence was achieved for each condition. Further, the self-propulsion performance was evaluated using the thrust identity method, wherein the propeller rotational speed was varied under a given inflow velocity condition until the hull resistance and propeller thrust were balanced, following the same procedure as in the model tests. Subsequently, the advance ratio J and self-propulsion factors were determined based on the rotational speed at the self-propulsion point. This procedure can be understood as a model-scale self-propulsion factor evaluation based on the ITTC propulsion test and 1978 ITTC performance prediction method; however, full-scale extrapolation was not performed in this study, and only model-scale self-propulsion factors were obtained.
The self-propulsion condition was set at a ship speed of 2.75 m/s, and the hull efficiency, relative rotative efficiency, and advance ratio were compared at the obtained self-propulsion point. The POW efficiency and advance ratio were calculated using the numerical POW results from this study and model test-based POW results. The results of this study are summarized in Table 11 with the previous numerical results of Chase (2012), Sezen et al. (2018), and Byeon et al. (2018). The RSM results showed an overall trend similar to those of Byeon et al. (2018) and Chase (2012). KT showed the closest agreement with previous studies in the case of RSM, whereas SST k–ω and LES predicted relatively lower KT values. The hull resistance calculated under the present self-propulsion condition (2.75 m/s) was 73.60, 88.0, and 80.0 N for SST k–ω, RSM, and LES, respectively. The differences in thrust and torque in the POW analysis were not as large as the differences in resistance, which are reflected in the prediction trends of KT and KQ. In addition, the relative rotative efficiency (ηR) based on INSEAN EFD shows a similar trend when compared with the results of Sezen et al. (2018) using the same turbulence models, even when differences exist in KT prediction. This can be interpreted as the self-propulsion factors being determined by the hull–propeller interaction than by individual coefficients (KT, KQ). Further validation of the present results through comparison with model test data is required in future work.
5. Conclusions
This study was conducted to provide fundamental data for improving future submarine self-propulsion test procedures and for cavitation and noise research. The self-propulsion analysis of the DARPA SUBOFF submarine and the POW analysis of the INSEAN E1619 propeller were performed and validated through comparisons with experimental data and previous numerical results. The variations in the POW performance, wake characteristics, and self-propulsion performance based on the turbulence model were examined using RANS-based turbulence models (SST k–ω and RSM) and the higher-fidelity turbulence model LES.
In the POW analysis, steady simulations based on the MRF approach were performed under the INSEAN model test conditions, and an improved agreement compared to that reported in previous studies was confirmed. Further, based on the grid system validated through GCI verification, the differences in POW results according to turbulence models, calculation methods (MRF/sliding mesh), and Reynolds number variations were compared quantitatively. The LES results showed a different trend in the thrust–torque relationship compared to the RANS results, indicating the possibility of variations in self-propulsion factors (η0, ηR, etc.) in the self-propulsion analysis.
Unlike in the previous studies, the self-propulsion analysis was performed considering a test procedure in which a propeller scaled down from 0.485 m to 0.262 m in diameter was applied. Self-propulsion factors were evaluated using the POW results obtained from model tests and those from the present study to analyze the resulting changes in prediction characteristics. Although a direct comparison with the model test results was not conducted, the validity of the computational results was assessed through comparison with previous numerical studies. Consequently, the RSM showed an overall good agreement with the RSM and DDES results of the previous studies, while SST k–ω showed some differences in KT but exhibited a similar trend in η0 compared to the previous SST k–ω results. Based on this, the characteristics of variations in self-propulsion factors according to turbulence models were summarized through comparison with LES results.
Based on model test results already obtained in the large cavitation tunnel (LCT) of the Korea Research Institute of Ships and Ocean Engineering (KRISO), the findings of this study will be used in this future to perform CFD simulations under identical test conditions and compare and validate the results. The resistance, self-propulsion, and cavitation test results obtained under LCT conditions will be quantitatively compared with CFD results to evaluate prediction accuracy, and changes in resistance and self-propulsion characteristics caused by wall effects between the LCT and free-surface conditions will be analyzed. These validation results are expected to be used as fundamental data for cavitation and noise research.
PDF Links
PubReader
ePub Link
Full text via DOI
Download Citation
Print
