Parametric Study of the Hydrodynamic Characteristics of the Pumpjet Propulsor for the SUBOFF Submarine

: Submarines with pumpjet propulsors have recently been used in many countries to improve their propulsion and noise performance. The pumpjet has the advantage of improving noise and cavitation performance by increasing the pressure inside the duct, and it has good efﬁciency at a low advance ratio. This study analyzed the propulsive performance of a pumpjet propulsor afﬁxed to a SUBOFF submarine as a function of variations in the design parameters using computational ﬂuid dynamics analysis. The incidence angle and camber of the duct and the pitch angle of the stator were selected as the design parameters. To investigate the propulsion performance as a function of the design parameter variations, each case’s thrust, torque, and efﬁciency were compared, and the velocity and pressure ﬁelds of each case were analyzed. In addition, the efﬁciency was analyzed using the non-dimensional mass ﬂow rate and area ratio difference for each case. The duct incidence angle contributed most dominantly to the dimensionless ﬂow rate and the difference in area ratio, and these two factors resulted in high efﬁciency at certain values. It is expected that further research will be conducted in the near future that takes into account cavitation inception speed and cavitation.


Introduction
Recently, pumpjet propulsors (PJPs) have been used in submarines for military purposes.The PJP is renowned for its outstanding cavitation performance.The PJP has improved cavitation performance, with propulsion efficiency comparable to that of a conventional propeller.The PJP consists of a rotor, a stator, and a duct surrounding the two aforementioned components.By reducing the flow velocity in the duct, the pressure of the flow field in the duct is increased.This improves cavitation performance and increases resistance to noise caused by cavitation.In a PJP, the function of the stator is to recover the rotational energy by controlling the tangential velocity of the flow entering the rotor, which is the same function as in an energy-saving device.The stator also improves cavitation performance by correcting the downstream swirl of the rotor.There is a slight difference in the role of the stator depending on whether it is located upstream or downstream of the rotor.The PJP with the stator in front of the rotor is called a stator-rotor PJP (S-R PJP), whereas the PJP with the stator behind the rotor is referred to as the rotor-stator PJP (R-S PJP) [1].The S-R and R-S PJP are shown schematically in Figure 1.Because the stator of the R-S PJP is responsible for approximately 25% of the total thrust, it serves to reduce the load on the rotor.This reduces the noise and cavitation generated by the rotor.However, the stator of the S-R PJP effectively improves the noise performance of the rotor by improving Typically, an optimized PJP has a smaller diameter than a conventional propell particular, the ratio of the propeller diameter to the ship's hull is smaller for a subm than for a torpedo.This means that the PJP of the submarine is highly dependent o boundary layer flow.Thurston and Evanbar [2] derived the energy relationship defi the efficiency of a propeller-inducing boundary layer flow in a body of revolution published the quantitative interrelationships of the corresponding parameters.demonstrated that a substantial gain in propulsion efficiency can be achieved by indu a boundary layer flow as opposed to a free stream flow and claimed that signific higher propulsion efficiency was possible with the PJP with submerged bodies of re tion.However, the limitations and compromises imposed by the advance ratio, cavita and structure were not analyzed in detail in this study.The design of a PJP is challen due to the difficulty of analyzing the turbulent boundary layer flow and the compl of the interactions between the components of the PJP.
The PJP was designed using methods based on the potential flow and simplifie draulic analogies and charts.McCormick and Eisenhuth [3] described a design proce for PJPs that takes cavitation and efficiency into account.They mentioned that the should be designed first because it regulates the velocity and pressure of the flow ent the rotor and has a large effect on the load at the rotor tip.In this study, the duct replaced by a ring vortex, and momentum theory was used to determine the velocity pressure upstream of the rotor section as well as the required pressure increase.The c line of the duct should be aligned with the velocity determined by the hull and rot avoid negative pressure peaks at the leading edge of the duct.The duct section shou determined after all design procedures, taking into account the possibility of flow se tion on the outer surface of the duct.Henderson et al. [4] investigated the design of d for PJPs.To design the duct, the flow between the duct and the central body of the pulsor was simplified as a one-dimensional flow.Bruce et al. [5] published a design cedure for a wake-adapted pumpjet attached to the stern of an axisymmetric hull in t of energy and resistance and used this procedure to design the PJP of the Akron air In this study, the mass flow of the propulsor was chosen to minimize kinetic energy through the duct.Furuya and Chiang [6] developed a quasi-three-dimensional PJP d method combining blade-to-blade flow theory and blade-through flow theory to Typically, an optimized PJP has a smaller diameter than a conventional propeller.In particular, the ratio of the propeller diameter to the ship's hull is smaller for a submarine than for a torpedo.This means that the PJP of the submarine is highly dependent on the boundary layer flow.Thurston and Evanbar [2] derived the energy relationship defining the efficiency of a propeller-inducing boundary layer flow in a body of revolution and published the quantitative interrelationships of the corresponding parameters.They demonstrated that a substantial gain in propulsion efficiency can be achieved by inducing a boundary layer flow as opposed to a free stream flow and claimed that significantly higher propulsion efficiency was possible with the PJP with submerged bodies of revolution.However, the limitations and compromises imposed by the advance ratio, cavitation, and structure were not analyzed in detail in this study.The design of a PJP is challenging due to the difficulty of analyzing the turbulent boundary layer flow and the complexity of the interactions between the components of the PJP.
The PJP was designed using methods based on the potential flow and simplified hydraulic analogies and charts.McCormick and Eisenhuth [3] described a design procedure for PJPs that takes cavitation and efficiency into account.They mentioned that the duct should be designed first because it regulates the velocity and pressure of the flow entering the rotor and has a large effect on the load at the rotor tip.In this study, the duct was replaced by a ring vortex, and momentum theory was used to determine the velocity and pressure upstream of the rotor section as well as the required pressure increase.The chord line of the duct should be aligned with the velocity determined by the hull and rotor to avoid negative pressure peaks at the leading edge of the duct.The duct section should be determined after all design procedures, taking into account the possibility of flow separation on the outer surface of the duct.Henderson et al. [4] investigated the design of ducts for PJPs.To design the duct, the flow between the duct and the central body of the propulsor was simplified as a one-dimensional flow.Bruce et al. [5] published a design procedure for a wake-adapted pumpjet attached to the stern of an axisymmetric hull in terms of energy and resistance and used this procedure to design the PJP of the Akron airship.In this study, the mass flow of the propulsor was chosen to minimize kinetic energy loss through the duct.Furuya and Chiang [6] developed a quasi-three-dimensional PJP design method combining blade-to-blade flow theory and blade-through flow theory to overcome the limitations of previous approaches.Whang et al. [7] presented a panel method using a tip-leakage vortex model for analyzing the steady-state hydrodynamic performance of the PJP.The PJP was separated into two independent systems, the rotor-hub and stator-duct, and the interaction between the systems was evaluated based on the induced velocity.
Computational fluid dynamics (CFD) analysis has been used in multiple studies to simulate the flow around the PJP and estimate its performance because the potential flow theory has theoretical limitations for propeller problems inside the duct [1,[8][9][10][11].Suryanarayana et al. [9][10][11] designed a PJP for a torpedo using CFD, and its performance was validated by experiments.Qin et al. [12] performed CFD analysis based on an improved delayed detached eddy (IDDE) simulation to study the flow around a PJP.CFD analysis was performed for the R model (rotor only), R-D model (rotor within the duct), and PJP model (rotor and stator within the duct) to investigate the role of each component.Li et al. [13] analyzed the ambient flow and propulsion performance of a PJP attached to a submarine using the same CFD analysis method as the study by Qin; the results showed that the performance of the PJP was significantly different from the performance in open-water conditions due to the hull-retarded flow.
As mentioned earlier, the PJP consists of a duct, stator, and rotor; consequently, there are many design parameters.Therefore, it is time-consuming and expensive to analyze the effects of the parameters through experiments.In the majority of studies, the effects of the parameters on the PJP have been analyzed numerically.Due to the narrow distance between the inner surface of the duct and the tip of the rotor, the performance of the PJP, including cavitation and efficiency, is sensitive to variations in the distance.Many studies have examined the hydrodynamic characteristics as a function of the gap [8,14,15].As the tip clearance increases, the propulsion and cavitation performances of the PJP decreases.Lu et al. [16] compared the results for cavitating and non-cavitating conditions and discovered that at high advance ratios, the difference in efficiency was smaller at small clearances for the non-cavitating condition and was bigger at small clearances for the cavitating condition.In addition to tip clearance, numerous parametric studies have been conducted on other parameters.Yu et al. [17] investigated the variation of the propulsion performance of a PJP as a function of stator parameters, such as stator pitch angle, chord length, and rotorstator spacing, using CFD analysis.As the stator pitch angle increased, the circumferential velocity of the rotor inflow and the overall thrust and propulsion efficiency of the PJP increased.However, the noise performance may degrade due to the sharp increase in unsteady pressure fluctuations.The rotor-stator spacing had no significant effect on the propulsion performance.However, if the rotor-stator spacing is appropriately selected, the fluctuation amplitudes of the unsteady force are significantly reduced.The chord length of the stator had little effect on the propulsion performance and maximum fluctuation amplitude of the unsteady force.Whang et al. [18] and Huang et al. [19] investigated the effects of duct parameters on the PJP.Whang et al. [18] used the same analytical model as the previous study to analyze the propulsion performance of a PJP according to the gap clearance, camber, and angle of attack of the duct.Huang et al. [19] systematically analyzed the effects of duct parameters (length-diameter ratio, incidence angle, shrinkage ratio at the duct inlet, expansion ratio at the duct outlet, and tip clearance) in the PJP using CFD analysis.In this study, the duct parameters were analyzed in terms of openwater performance, thrust fluctuation, pressure field and its fluctuation, velocity field, and vortices.
The effects of duct parameters have also been investigated in studies of typical ducted propellers.Theoretical calculations performed by Oosterveld [20] demonstrated the effects of a decelerating duct on the length and thickness of the duct.The risk of flow separation at the outer surface of the duct was reduced by increasing the length or decreasing the thickness of the decelerating duct.Huyer and Dropkin [21] performed a CFD analysis of a two-dimensional axisymmetric duct behind a ship's hull to determine the effects of duct parameters.This study did not parameterize the geometry of the duct because it focused on the effect of the flow inside the duct.The study presented the velocity distribution of the flow within the duct as well as the pressure distribution on the duct based on the difference between the mean flow velocities within the duct.The results showed that the drag on the hull decreased, while the drag on the duct increased with an increase in the difference between the mean flow velocities.As the difference in mean flow velocities increased, the mass flow rate inside the duct decreased.Bontempo et al. [22] applied axial momentum theory and a nonlinear semi-analytical actuator disk model to investigate the flow around a decelerating duct with varying length, camber, and thicknesses of the duct.The study showed that the length had the smallest effect on the drag of the duct compared to its camber and thickness.Although the cavitation behavior of the rotor can be improved by increasing the camber and length of duct or decreasing the thickness of duct, increasing the camber causes cavitation in the mid-chord of the duct, and decreasing the thickness causes cavitation at the leading edge of the duct at very high advance ratios.Gaggero et al. [23] presented an optimization approach for the design of accelerating and decelerating ducts.Multiple control points were used to modify the geometry of the duct, which was represented by a B-spline curve.This method can be performed with an optimization process using diverse data.
Previous studies have shown that ducts have a significant impact on PJP performance.This is because the duct has a significant effect on flow entering the stator and duct.As mentioned earlier, potential analysis is a limitation in analyzing the closed-tip clearance in the duct-propeller problem.In addition, the majority of previous studies have analyzed the characteristics of the PJP under open-water conditions.However, the flow behind bodies of revolution, such as submarines, is quite different from the flow under open water conditions due to the large inclination angle of the astern body hull form.
In this study, a parametric analysis of the PJP of the SUBOFF submarine was performed.The incidence angle and camber were selected as the duct parameters and the pitch angle was selected as the stator parameter.The hydrodynamic characteristics of the PJP were investigated for each parameter, and the results were summarized and analyzed according to the variations of the area ratio difference between the inlet and outlet stations of the duct and the operating station of the rotor.CFD was used for most of the computations and the initial design of the rotor and stator was performed using the potential analysis.
The in-house developed potential code (PASTA) was used for the design [24], and Star CCM+ v15.04 was used for CFD computations.The CFD analysis was verified by a grid dependence test, and a comparison with experimental results from the David Taylor Research Center (DTRC) and the Korea Research Institute of Ship & Ocean Engineering (KRISO) was performed for validation.
For the proper design of a submarine-PJP, the current study focuses on a parametric study that takes favorable flow conditions and hydrodynamic efficiency into account.Although not specified, favorable flow conditions may also contribute to the improvement of cavitation performance.It is expected that cavitation performance and noise will be studied when cavitation inception speed (CIS) and cavitation conditions are explicitly specified, and future research will investigate optimization with two objective functions of efficiency and cavitation, including CIS.The current study may be the first step toward achieving this objective.
The first section of the current study presents the specifications of the target submarine.It describes the types of design variables of the PJP used in the analysis and their definitions, and the dimensionless coefficients used for the performance analysis.The second section describes the numerical method used in the current study.The numerical methods include governing equations, turbulence model, grids, and boundary conditions.The third section presents the results of the grid dependence tests and comparison with experimental results for validation and describes the results of the parametric study.

Submarine
The submerged target vehicle was the SUBOFF submarine (Model No. 5470) of the Defense Advanced Research Projects Agency (DARPA) conducted by the DTRC, and the arrangement of appendages corresponded to AFF-8.It consists of a fairwater and four stern appendages attached to an axisymmetric body [25].The shape and dimensions of this submarine are shown in Figure 2 and Table 1, respectively [26].The total length and maximum diameter of the submarine model were 4.356 and 0.508 m, respectively.r.Sci.Eng.2023, 11, x FOR PEER REVIEW stern appendages attached to an axisymmetric body [25].The shap this submarine are shown in Figure 2 and Table 1, respectively [26].maximum diameter of the submarine model were 4.356 and 0.508 m,  Because the hub diameter of the rotor was smaller than that of rine in the rotor plane, separate hubs and caps were attached to moun the total body length of the SUBOFF submarine was slightly reduced

Pumpjet Propulsor
In the current study, the S-R PJP was chosen for the paramet higher efficiency than the R-S PJP at a low advanced coefficient.This noise performance due to the low speed of revolution, despite the ciency compared to the R-S PJP [12,27].
Figure 3 shows the geometry and cross-sectional profile of the P the duct and stator parameters.The incidence angle and camber of t by the definition of the cross-section of the annular airfoil.The flow divided into the inlet, outlet, and working stations.The reference po coordinates are the leading and trailing edges of the duct and rotor The inlet, outlet, and working areas were defined as the areas betwe hull.Stator parameters are defined similarly to those of typical prop Table 1.The specifications of the DARPA SUBOFF submarine.The data were from [26].

Item Specification
Fore-body Length Because the hub diameter of the rotor was smaller than that of the SUBOFF submarine in the rotor plane, separate hubs and caps were attached to mount the PJP.As a result, the total body length of the SUBOFF submarine was slightly reduced to 4.287 m.

Pumpjet Propulsor
In the current study, the S-R PJP was chosen for the parametric analysis.It has a higher efficiency than the R-S PJP at a low advanced coefficient.This leads to a beneficial noise performance due to the low speed of revolution, despite the low maximum efficiency compared to the R-S PJP [12,27].
Figure 3 shows the geometry and cross-sectional profile of the PJP used to represent the duct and stator parameters.The incidence angle and camber of the duct are reflected by the definition of the cross-section of the annular airfoil.The flow through the PJP was divided into the inlet, outlet, and working stations.The reference points along the x-axis coordinates are the leading and trailing edges of the duct and rotor plane, respectively.The inlet, outlet, and working areas were defined as the areas between the duct and the hull.Stator parameters are defined similarly to those of typical propeller parameters, except for the stator's local coordinate system.Thus, the pitch angle is defined as the angle formed by the coordinates of the stator axis and the chord line of the stator.Table 3 lists the design parameters of the PJP used in this study.For simplicity, the duct incidence angle, duct camber, and stator pitch angle are abbreviated as DIA, DC, and SPA, respectively.The DIA was chosen over a wide range to observe the effects of changes at the inlet and outlet of the duct, whereas the SPA was chosen over a wide range to observe the effects of tangential energy recovery and decreasing thrust.The DC was chosen over a range from 0% to 5.6% to analyze the effects of the deceleration duct in the case where the DIA is aligned with the streamline of the flow around the PJP.In the current study, the term "negative thrust" was used instead of "resistance".This reduces terminology confusion, as ducts can either generate thrust or resistance.It is also because the thrust of the rotor, stator, and ducts are used to define the performance of the PJP.The non-dimensional coefficients used to evaluate the propulsion performance of the PJP are defined as follows: A parametric study was performed by varying the incidence angle and camber of the duct and the pitch angle of the stator with respect to the reference PJP.The parameters such as the rotor diameter, rotor and stator blade numbers, and duct length were decided based on the authors' experience.Those parameters are not thought to be the optimum values, but they do not appear to be too far from optimal values because the performance of the reference PJP is not bad.The specifications of the reference PJP are listed in Table 2.The rotor diameter is 0.21 m, and the tip clearance represents the distance between the rotor tip and the duct expressed as a percentage of the rotor diameter.The rotor and stator each had seven and eleven blades, respectively.Table 3 lists the design parameters of the PJP used in this study.For simplicity, the duct incidence angle, duct camber, and stator pitch angle are abbreviated as DIA, DC, and SPA, respectively.The DIA was chosen over a wide range to observe the effects of changes at the inlet and outlet of the duct, whereas the SPA was chosen over a wide range to observe the effects of tangential energy recovery and decreasing thrust.The DC was chosen over a range from 0% to 5.6% to analyze the effects of the deceleration duct in the case where the DIA is aligned with the streamline of the flow around the PJP.In the current study, the term "negative thrust" was used instead of "resistance".This reduces terminology confusion, as ducts can either generate thrust or resistance.It is also because the thrust of the rotor, stator, and ducts are used to define the performance of the PJP.The non-dimensional coefficients used to evaluate the propulsion performance of the PJP are defined as follows: where, m are the speed of the submarine, the working area of the rotor, and the mass flow rate, respectively.V A is the advance velocity of the thruster, which is calculated by averaging the velocity at the rotor plane through CFD analysis.The subscripts r, s, d denote the rotor, stator, and duct, respectively.

Numerical Approach
To compute the flow around the submarine and PJP, the incompressible fluid field was modeled using the Reynolds Averaged Navier-Stokes (RANS) equations.The continuity and momentum equations are expressed as follows: where, U is velocity, x is direction, t is time, ρ is density, µ is viscosity, and −ρu i u j is Reynolds stress.The Reynolds stress generated during the averaging process of the momentum equation consists of six unknown components.As there were unknown velocity components and pressures, the number of unknowns was ten.Therefore, six additional equations are required to close this equation.An Elliptic Blending Reynolds Stress Model (EB-RSM) was used to perform the calculations.The RSM predicts complex flows more accurately than eddy viscosity models because the transport equations for Reynolds stresses naturally account for the effects of turbulence anisotropy, streamline curvature, swirl rotation, and high strain rates [28].The EB-RSM is based on an inhomogeneous near-wall formulation of a quasilinear quadratic pressure-strain term.The blending function was used to blend the viscous sub-layer with the log-layer formulation of the pressure-strain term.This approach requires the solution of an elliptic equation for the blending function [29].It is known that the EB-RSM has a relatively low sensitivity to the grid, which is a limitation of the RSM [30].
The numerical analysis was performed using the popular software Star-CCM+.The finite volume method (FVM) was used to compute the governing equation, and the semiimplicit method for pressure-linked equations (SIMPLE) was used as the pressure-velocity coupling algorithm.

Grid and Numerical Setup
The CFD analysis was performed for the resistance and self-propulsion tests.The properties of the fluid were determined by administering freshwater with a temperature of 23.4 • C as recommended by the International Towing Tank Conference (ITTC).The CFD analysis for the resistance test was performed for the steady-state fluid and the velocity condition range was between 1.58 and 9.14 m/s.The corresponding Reynolds number ranged from 7.43 × 10 6 to 4.31 × 10 7 .The CFD analysis for the self-propulsion test was performed for the unsteady-state fluid and the velocity conditions of 2.63, 3.15, and 3.68 m/s.The towing forces used for the CFD parametric study were 26.39, 36.28, and 47.50 N for each speed.However, the towing forces in the experiments conducted by KRISO were 31.43,43.42, and 57.04 N. Therefore, the results of the self-propulsion test were compared between the CFD analysis and the model test using the towing force used in the experiment.The reason is that the PJP was designed using parametric research.The time step for the simulation was the time it takes for the rotor to rotate 2 • .
Figure 4 shows the flow domain and boundary conditions for the CFD analysis.The global coordinate system was defined in the flow direction along the x-axis, starboard direction along the y-axis, and top direction along the z-axis.The boundaries of the flow domain satisfied the velocity inlet and pressure outlet conditions at the inlet and outlet, while the remaining boundaries satisfied the symmetry conditions.The shape of the flow domain is hexahedral and the distance to the inlet is 2.5 LOA, the distance to the outlet is 7.5 LOA, and the distance in the up, down, left, and right directions is 2.5 LOA, with respect to the rotor center.
between the CFD analysis and the model test using the towing force used in the experiment.The reason is that the PJP was designed using parametric research.The time step for the simulation was the time it takes for the rotor to rotate 2°.
Figure 4 shows the flow domain and boundary conditions for the CFD analysis.The global coordinate system was defined in the flow direction along the x-axis, starboard direction along the y-axis, and top direction along the z-axis.The boundaries of the flow domain satisfied the velocity inlet and pressure outlet conditions at the inlet and outlet, while the remaining boundaries satisfied the symmetry conditions.The shape of the flow domain is hexahedral and the distance to the inlet is 2.5 LOA, the distance to the outlet is 7.5 LOA, and the distance in the up, down, left, and right directions is 2.5 LOA, with respect to the rotor center.Figure 5 shows the coarse, medium, and fine grid systems of the submarine with the PJP used to verify the grid dependence and grid configuration inside and outside the PJP.Table 4 lists the number of grids in each system.The length of the grids in the coarse, medium, and fine grid systems was reduced by a factor of √2, and the total number of grid elements is 12.7 M, 19.0 M, and 34.0 M, respectively.The height of the first cell layer close to the wall was set so that y+ equals one.Although the length of the grid element increased due to the grid dependence test, y+ was maintained at less than two.The y+ value for the CFD analysis is shown in Figure 6. Figure 5 shows the coarse, medium, and fine grid systems of the submarine with the PJP used to verify the grid dependence and grid configuration inside and outside the PJP.Table 4 lists the number of grids in each system.The length of the grids in the coarse, medium, and fine grid systems was reduced by a factor of √ 2, and the total number of grid elements is 12.7 M, 19.0 M, and 34.0 M, respectively.The height of the first cell layer close to the wall was set so that y+ equals one.Although the length of the grid element increased due to the grid dependence test, y+ was maintained at less than two.The y+ value for the CFD analysis is shown in Figure 6.Table 4.The number of grid elements required to validate the current computations for each grid system.Table 4.The number of grid elements required to validate the current computations for each grid system.

Grid Dependence and Validation Results
Grid dependence tests were performed with three different grid systems to verify the CFD analysis.The results of the resistance acting on the hull, as well as the thrust and

Grid Dependence and Validation Results
Grid dependence tests were performed with three different grid systems to verify the CFD analysis.The results of the resistance acting on the hull, as well as the thrust and torque coefficients of the PJP at the three velocity conditions are listed in Table 5.All the results obtained from the unsteady analysis were averaged over one rotation of the rotor after the flow had sufficiently converged.The three cases with increasing grid density converged to the similar values of residuals.The medium grid system was chosen for the current CFD computation because the difference between the results of the fine and medium grid systems was less than 0.2% and the number of grid elements in the fine grid system was considerably larger.The CFD computation was validated by comparing experimental results conducted by the DTRC and KRISO for the SUBOFF submarine, as shown in Figure 7 [25,31].The trends of the three cases are nearly identical, and the resistance differences between the experiment and CFD are 4.11% and 7.32% at the target speed of 3.15 m/s for the DTRC and KRISO, respectively.The results of the self-propulsion test showed a difference of 2.38% for K T and 6.22% for K Q at the target speed of 3.15 m/s.Due to grid dependence and comparison with experimental results from reputable institutions, the current simulation system can be deemed trustworthy.The CFD computation was validated by comparing experimental results conducted by the DTRC and KRISO for the SUBOFF submarine, as shown in Figure 7 [25,31].The trends of the three cases are nearly identical, and the resistance differences between the experiment and CFD are 4.11% and 7.32% at the target speed of 3.15 m/s for the DTRC and KRISO, respectively.The results of the self-propulsion test showed a difference of 2.38% for   and 6.22% for   at the target speed of 3.15 m/s.Due to grid dependence and comparison with experimental results from reputable institutions, the current simulation system can be deemed trustworthy.

Propulsion Performance Results
The verified CFD analysis method was used to analyze the propulsion performance as a function of the variation in the selected parameters of the PJP at a design speed of 3.15 m/s.The selected parameters of the PJP were the duct incidence angle (DIA), duct camber

Propulsion Performance Results
The verified CFD analysis method was used to analyze the propulsion performance as a function of the variation in the selected parameters of the PJP at a design speed of 3.15 m/s.The selected parameters of the PJP were the duct incidence angle (DIA), duct camber (DC), and stator pitch angle (SPA).Figure 8 shows the thrust coefficients of the rotor, stator, duct, and overall propulsor, in addition to the torque coefficient of the rotor and the efficiency of the PJP.The results for the DIA and DC show that the highest efficiency is obtained when the duct produces a positive thrust.In addition, the thrusts of the rotor and duct exhibit opposite trends, whereas the thrust of the stator increased marginally but not considerably.In the case of DC, the total thrust tends to resemble the thrust of the stator because the thrusts of the rotor and duct cancel one another out.However, in the case of DIA, the total thrust decreases substantially.The reason for this is that the change in the DIA significantly affects the drag on the hull, so the required thrust also changes significantly.The SPA results also show an increase in total thrust, but the thrust coefficient increases as the rotor speed decreases, and the efficiency is highest at approximately 15 • .At this the angle, the recovery of tangential energy and the reduction of stator thrust are properly balanced.
case of DIA, the total thrust decreases substantially.The reason for this is that the change in the DIA significantly affects the drag on the hull, so the required thrust also changes significantly.The SPA results also show an increase in total thrust, but the thrust coefficient increases as the rotor speed decreases, and the efficiency is highest at approximately 15 °.At this the angle, the recovery of tangential energy and the reduction of stator thrust are properly balanced.Figures 9 and 12 show the results of variations in DIA.A small DIA causes flow separation near the leading edge inside the duct by creating a stagnation point on the outer surface of the leading edge of the duct.The radial velocity distribution inside the duct changes sharply near a 0.7 of radius due to flow separation.The separated flow passes through the stator and negatively affects the stator downstream.In addition, these downstream flows have a negative effect on the rotor and strengthen the rotor tip vortex.In addition, the change in radial velocity distribution within the duct affects the angle of attack of the rotor and stator, causing the pressure in the leading edge of the rotor and stator to change rapidly in the span direction.This tendency gradually disappears according to the increase in DIA.As the DIA increases, the axial velocity within the duct decreases, which causes the rotor with an increased angle of attack to increase the pressure on the pressure side and decrease the pressure on the suction side.In contrast, the trailing edge of the rotor has increased the pressure on both the pressure and suction sides.The low-pressure area at the rotor tip becomes larger according to the increase in the angle of attack at the rotor tip.The pressure on the duct is divided into three parts according to its characteristics: the duct`s inner surface and the leading and trailing edge of the duct`s outer surface.It decreases at the leading edge of the outer surface and increases overall at the inner surface and remains unchanged at the trailing edge of the outer surface accord- Figures 9-11 show the pressure distributions on the surface of the rotor, stator, and duct corresponding to the variations of DIA, DC, and SPA, respectively.The pressure distribution is represented by C P , which is defined as C P = (P − P ∞ )/ 0.5ρV 2 ∞ .Figures 12-14 show the streamline at the leading edge of the duct and the axial velocity and vorticity magnitude distributions around the PJP on the y = 0 plane corresponding to the variations of DIA, DC, and SPA, respectively.
Figures 9 and 12 show the results of variations in DIA.A small DIA causes flow separation near the leading edge inside the duct by creating a stagnation point on the outer surface of the leading edge of the duct.The radial velocity distribution inside the duct changes sharply near a 0.7 of radius due to flow separation.The separated flow passes through the stator and negatively affects the stator downstream.In addition, these downstream flows have a negative effect on the rotor and strengthen the rotor tip vortex.In addition, the change in radial velocity distribution within the duct affects the angle of attack of the rotor and stator, causing the pressure in the leading edge of the rotor and stator to change rapidly in the span direction.This tendency gradually disappears according to the increase in DIA.As the DIA increases, the axial velocity within the duct decreases, which causes the rotor with an increased angle of attack to increase the pressure on the pressure side and decrease the pressure on the suction side.In contrast, the trailing edge of the rotor has increased the pressure on both the pressure and suction sides.The low-pressure area at the rotor tip becomes larger according to the increase in the angle of attack at the rotor tip.The pressure on the duct is divided into three parts according to its characteristics: the duct's inner surface and the leading and trailing edge of the duct's outer surface.It decreases at the leading edge of the outer surface and increases overall at the inner surface and remains unchanged at the trailing edge of the outer surface according to the increase of DIA.When the DIA is between 0 • and 6 • , the duct acts as an acceleration duct because the internal pressure was lower than the external pressure.A low DIA is not suitable for a PJP because the duct acts as an acceleration duct, and the radial pressure distributions of the rotor and stator change rapidly.However, if the DIA increases too much, the risk of flow separation on the outer surface of the duct at the trailing edge increases and a vortex is generated at the duct trailing edge due to the velocity difference between flow inside and outside the duct.Figures 10 and 13 show the results of variations in DC.High DC forms a stagnation point on the outer surface of the duct leading edge and causes flow separation, although it is smaller than low DIA.However, unlike low DIA, high DC reduces the axial flow velocity within the duct, thereby increasing the angle of attack of the rotor and increasing the CP difference between the pressure and suction sides near the rotor's leading edge.Also, the pressure within the duct and the pressure at the leading edge of the duct's outer surface increase according to the increase in DC.The increased pressure inside the duct somewhat suppresses the increase in the low-pressure area at the rotor tip caused by the increased angle of attack of the rotor.High DC increases the risk of flow separation at the trailing edge of the duct and strengthens the vortex at the trailing edge of the duct.Low DC and low DIA reduce the duct inlet, while high DC and high DIA reduce the duct outlet.Therefore, the hydrodynamic features at the leading and trailing edges of the duct are comparable.
Because low DC causes the duct to behave as an acceleration duct, as does low DIA, it is not suitable for PJP.
(c) Figures 11 and 14 show the results of variations in SPA.High SPA increases the angle of attack of the stator and rotor, which increases the pressure difference between the pressure and suction sides of the rotor and stator.As the angle of attack increases, so does the low-pressure area at the rotor tip.High SPA causes flow separation at the high radius of the stator and increases the resistance of the stator.However, proper SPA reduces the tangential velocity of the rotor downstream, thereby weakening the rotor's downstream vortex and increasing the rotor's efficiency.

Parametric Study Discussion
The design parameters were summarized in terms of the non-dimensional mass flow rate and the area ratio difference to provide an intuitive comparison between the design parameters of the PJP.The area ratio difference is the difference in the non-dimensional area ratios, which is the working area of the rotor divided by the inlet and outlet areas of the duct.These were denoted as        Figures 10 and 13 show the results of variations in DC.High DC forms a stagnation point on the outer surface of the duct leading edge and causes flow separation, although it is smaller than low DIA.However, unlike low DIA, high DC reduces the axial flow velocity within the duct, thereby increasing the angle of attack of the rotor and increasing the CP difference between the pressure and suction sides near the rotor's leading edge.Also, the pressure within the duct and the pressure at the leading edge of the duct`s outer surface increase according to the increase in DC.The increased pressure inside the duct somewhat suppresses the increase in the low-pressure area at the rotor tip caused by the increased angle of attack of the rotor.High DC increases the risk of flow separation at the trailing edge of the duct and strengthens the vortex at the trailing edge of the duct.Low Velocity and vorticity distribution around the PJP according to 0, 2.4, 4.0, 5.6% of DC at plane x = 0 (top: streamlines and axial velocity (V X /V) around the duct leading edge; middle: axial velocity (V X /V) distribution around the PJP; bottom: vorticity (ν M ) distribution around the PJP).Figures 10 and 13 show the results of variations in DC.High DC forms a stagnation point on the outer surface of the duct leading edge and causes flow separation, although it is smaller than low DIA.However, unlike low DIA, high DC reduces the axial flow velocity within the duct, thereby increasing the angle of attack of the rotor and increasing the CP difference between the pressure and suction sides near the rotor's leading edge.Also, the pressure within the duct and the pressure at the leading edge of the duct`s outer surface increase according to the increase in DC.The increased pressure inside the duct somewhat suppresses the increase in the low-pressure area at the rotor tip caused by the increased angle of attack of the rotor.High DC increases the risk of flow separation at the trailing edge of the duct and strengthens the vortex at the trailing edge of the duct.Low  of SPA at plane x = 0 (top: streamlines and axial velocity (V X /V) around the duct leading edge; middle: axial velocity (V X /V) distribution around the PJP; bottom: vorticity (ν M ) distribution around the PJP).
Figure 15 shows the efficiency of each case based on the non-dimensional mass flow rate and the difference in area ratio.As the flow rate within the duct is constant, the difference in the area ratio becomes the difference in axial velocity ratio by applying the continuity equation (mass conservation).The change in axial velocity is the component that produces the overall thrust; therefore, the difference in area-ratio correlates to efficiency.Because efficiency is the ratio between the energy transferred by the rotor and the axial energy gain, the results between the design parameters do not match perfectly, but the correlation between the difference in area ratio and efficiency is notable.The non-dimensional mass flow rate according to variation of the DIA ranges from 0.6 to 0.85, the DC case ranges from 0.63 to 0.78, and the SPA case ranges from 0.67 to 0.7.According to the variations in DIA and DC, the highest efficiency was obtained when the mass flow rate was approximately 0.77.However, in the case of SPA, the highest efficiency was achieved at the lowest flow rate.This is because the duct has a significant effect on the axial and radial velocities, whereas the stator primarily affects the circumferential and axial velocities.As the efficiency depends on the variation of the area ratio difference, the variation of the area ratio difference due to the variation of DC is small compared to the variation of DIA, and SPA has no influence on the area ratio difference because it cannot change the configuration of the duct.The relationship between the area ratio difference and efficiency shows that the efficiency is the highest at an area ratio difference of approximately 0.1 according to the variation of DIA.It also shows that if the area ratio difference is in the positive range, a lower area ratio difference leads to higher efficiency.However, when the area ratio difference is reduced to a certain point, the duct acts as an accelerating duct and the pressure within the duct decreases accordingly.Therefore, it is recommended to maintain an area ratio difference of 0.3 or more to consider CIS or cavitation.
The area ratio difference can be used to determine the inner profile of the duct in terms of efficiency.However, the difference in the outer profile of the duct changes the As the efficiency depends on the variation of the area ratio difference, the variation of the area ratio difference due to the variation of DC is small compared to the variation of DIA, and SPA has no influence on the area ratio difference because it cannot change the configuration of the duct.The relationship between the area ratio difference and efficiency shows that the efficiency is the highest at an area ratio difference of approximately 0.1 according to the variation of DIA.It also shows that if the area ratio difference is in the positive range, a lower area ratio difference leads to higher efficiency.However, when the area ratio difference is reduced to a certain point, the duct acts as an accelerating duct and the pressure within the duct decreases accordingly.Therefore, it is recommended to maintain an area ratio difference of 0.3 or more to consider CIS or cavitation.
The area ratio difference can be used to determine the inner profile of the duct in terms of efficiency.However, the difference in the outer profile of the duct changes the flow inside and outside the duct, which has a substantial effect on the thrust of the rotor and duct.For ducts, the pressure distributions on the inner and outer surfaces changed significantly.The change of SPA has little effect on the thrust of the duct but has a significant effect on the performance of the rotor itself, affecting the pressure distribution on the rotor, stator, inner surface of the duct, and downstream vortex.Therefore, the interior profile of the duct must be designed based on the difference between the flow rate and area ratio, and the outer profile of the duct, as well as the rotor and stator, must be designed based on the pressure distribution of the PJP.

Conclusions
In this study, a parametric study of a PJP was conducted from a hydrodynamic point of view.The PJP was designed for the SUBOFF submarine, and its propulsion performance was evaluated for the design parameters (duct incidence angle (DIA), duct camber (DC), and stator pitch angle (SPA)).The propulsion performance for each design parameter was also analyzed using the non-dimensional mass flow rate, C m, and area ratio difference ( A W A O − A W A I ).Star-CCM+ was used as a tool for the parametric study.The CFD analysis method was validated by a grid dependence test and the experimental results from DTRC and KRISO were compared.For a parametric study of the PJP, each case's propulsion performance at the self-propulsion point was compared.The propulsion performance was also analyzed using a non-dimensional mass flow rate, C m, and area ratio difference ( A W A O − A W A I ).DIA is the most dominant influence on the hydrodynamic characteristics and efficiency of the PJP.This is because of the large variation in the non-dimensional mass flow rate and area ratio difference and the dominant effect on flow separation at the leading and trailing edge of the duct.DC has the same effect as DIA, but it is expected to change the rotor advance ratio mainly because it changes the non-dimensional mass flow rate more than the area ratio difference.SPA changes the performance of the rotor itself because it changes the tangential velocity, which causes a change in the angle of attack of the rotor.The highest efficiency was attained when the non-dimensional mass flow rate was approximately 0.77 and the area ratio difference was approximately 0.1.However, a PJP with high non-dimensional flow rate and low area ratio difference are not suitable for a PJP because the duct functions as an accelerating duct.Therefore, to design a PJP that is optimal in terms of efficiency and cavitation, the flow rate and area ratio of the PJP must be properly determined.In addition, the outer contour of the duct must be determined by taking into account the flow separation that occurs at the leading and trailing edge of the duct.
This study is expected to be the initial step toward the optimal design of PJPs for submarines.It is also anticipated that additional parameters, including the diameter, number of stator and rotor blades, and length of the duct, will be incorporated into a broader spectrum of applications in the near future.

Figure 2 .
Figure 2. The geometry of the DARPA SUBOFF submarine.

Figure 2 .
Figure 2. The geometry of the DARPA SUBOFF submarine.

Figure 3 .
Figure 3.The schematic configuration used to define the parameters: (a) the geometry and profile of the PJP; (b) the profile of the stator section.

Figure 3 .
Figure 3.The schematic configuration used to define the parameters: (a) the geometry and profile of the PJP; (b) the profile of the stator section.

Figure 4 .
Figure 4. Domain and boundary conditions for the CFD analysis.

Figure 4 .
Figure 4. Domain and boundary conditions for the CFD analysis.

Figure 5 .
Figure 5.The grid system of the SUBOFF submarine and PJP: (a) the coarse, medium, and fine grids; (b) the grid within the PJP in the case of the medium grid; (c) the detailed submarine and PJP grid in the case of the medium grid.

Figure 6 .
Figure 6.The y+ value for the entire submarine and within the PJP: (a) the entire submarine; (b) the PJP.

Figure 5 .
Figure 5.The grid system of the SUBOFF submarine and PJP: (a) the coarse, medium, and fine grids; (b) the grid within the PJP in the case of the medium grid; (c) the detailed submarine and PJP grid in the case of the medium grid.

Figure 5 .
Figure 5.The grid system of the SUBOFF submarine and PJP: (a) the coarse, medium, and fine grids; (b) the grid within the PJP in the case of the medium grid; (c) the detailed submarine and PJP grid in the case of the medium grid.

Figure 6 .
Figure 6.The y+ value for the entire submarine and within the PJP: (a) the entire submarine; (b) the PJP.

Figure 6 .
Figure 6.The y+ value for the entire submarine and within the PJP: (a) the entire submarine; (b) the PJP.

Figure 7 .
Figure 7.Comparison of resistance and self-propulsion test results between EFD and CFD: (a) resistance test results between the DTRC and KRISO and CFD analysis; (b) self-propulsion test resultsbetween the KRISO and CFD analysis.The data were form[31,32].

Figure 7 .
Figure 7.Comparison of resistance and self-propulsion test results between EFD and CFD: (a) resistance test results between the DTRC and KRISO and CFD analysis; (b) self-propulsion test resultsbetween the KRISO and CFD analysis.The data were form[31,32].

Figure 8 .
Figure 8. Self-propulsion performance as a function of changes in design parameters: (a) the DIA; (b) the DC; (c) the SPA.Figures 9-11 show the pressure distributions on the surface of the rotor, stator, and duct corresponding to the variations of DIA, DC, and SPA, respectively.The pressure distribution is represented by   , which is defined as   = ( −  ∞ )/(0.5∞ 2 ).Figures 12-14 show the streamline at the leading edge of the duct and the axial velocity and vorticity magnitude distributions around the PJP on the  = 0 plane corresponding to the variations of DIA, DC, and SPA, respectively.Figures9 and 12show the results of variations in DIA.A small DIA causes flow separation near the leading edge inside the duct by creating a stagnation point on the outer surface of the leading edge of the duct.The radial velocity distribution inside the duct changes sharply near a 0.7 of radius due to flow separation.The separated flow passes through the stator and negatively affects the stator downstream.In addition, these downstream flows have a negative effect on the rotor and strengthen the rotor tip vortex.In addition, the change in radial velocity distribution within the duct affects the angle of attack of the rotor and stator, causing the pressure in the leading edge of the rotor and stator to change rapidly in the span direction.This tendency gradually disappears according to the increase in DIA.As the DIA increases, the axial velocity within the duct decreases, which causes the rotor with an increased angle of attack to increase the pressure on the pressure side and decrease the pressure on the suction side.In contrast, the trailing edge of the rotor has increased the pressure on both the pressure and suction sides.The low-pressure area at the rotor tip becomes larger according to the increase in the angle of attack at the rotor tip.The pressure on the duct is divided into three parts according to its characteristics: the duct`s inner surface and the leading and trailing edge of the duct`s outer surface.It decreases at the leading edge of the outer surface and increases overall at the inner surface and remains unchanged at the trailing edge of the outer surface accord-

Figure 8 .
Figure 8. Self-propulsion performance as a function of changes in design parameters: (a) the DIA; (b) the DC; (c) the SPA.

Figure 9 .
Figure 9. Pressure coefficient (C P ) distribution according to 0, 6, 12, 18 • of DIA; the black dashed line is a reference line for geometrical comparison: (a) the rotor at the top position (top: pressure side, bottom: outside); (b) the duct (top: inside, bottom: outside); (c) the stator.

Figure 12 .
Figure 12.Axial velocity and vorticity magnitude distribution around the PJP according to 0, 6, 18 ° of DIA at x = 0 plane (top: streamlines and axial velocity (  /) around the duct leading ed middle: axial velocity (  /) distribution around the PJP; bottom: vorticity magnitude distribut around the PJP).

Figure 12 .
Figure12.Axial velocity and vorticity magnitude distribution around the PJP according to 0, 6, 12, 18 • of DIA at x = 0 plane (top: streamlines and axial velocity (V X /V) around the duct leading edge; middle: axial velocity (V X /V) distribution around the PJP; bottom: vorticity magnitude distribution around the PJP).

Figure 13 .
Figure13.Velocity and vorticity distribution around the PJP according to 0, 2.4, 4.0, 5.6% of DC at plane x = 0 (top: streamlines and axial velocity (V X /V) around the duct leading edge; middle: axial velocity (V X /V) distribution around the PJP; bottom: vorticity (ν M ) distribution around the PJP).

Figure 14 .
Figure 14.Velocity and vorticity distribution around the PJP according to 0, 10, 20, 30• of SPA at plane x = 0 (top: streamlines and axial velocity (V X /V) around the duct leading edge; middle: axial velocity (V X /V) distribution around the PJP; bottom: vorticity (ν M ) distribution around the PJP).

Figure 15 .
Figure 15.The efficiency of each case based on the non-dimensional mass flow rate and the area ratio difference: (a) the non-dimensional mass flow rate; (b) the area ratio difference.

Figure 15 .
Figure 15.The efficiency of each case based on the non-dimensional mass flow rate and the area ratio difference: (a) the non-dimensional mass flow rate; (b) the area ratio difference.

Table 1 .
The specifications of the DARPA SUBOFF submarine.The data wer

Table 2 .
The specifications of the reference PJP.

Table 3 .
Design parameters of the PJP.

Table 2 .
The specifications of the reference PJP.

Table 3 .
Design parameters of the PJP.
the advance ratio, thrust coefficient, torque coefficient, thruster efficiency, and non-dimensional mass flow rate and T, Q, n, D, ρ represent thrust, torque, number of revolutions of the rotor, diameter of the rotor, and water density, respectively.V ∞ , A W , .

Table 4 .
The number of grid elements required to validate the current computations for each grid system.

Table 5 .
The comparison of resistance, thrust, and torque depending on variation of the grid systems.

Table 5 .
The comparison of resistance, thrust, and torque depending on variation of the grid systems.
Table 6 lists the non-dimensional mass flow rate and area ratio for each design parameter.

Table 6 .
The non-dimensional mass flow rate and area ratio of duct inlet and outlet to rotor working surface as design parameters.