Influence of Propulsion Type on the Stratified Near Wake of an Axisymmetric Self-Propelled Body

To better understand the influence of swirl on the thermally-stratified near wake of a self-propelled axisymmetric vehicle, three propulsor schemes were considered: a single propeller, contra-rotating propellers (CRP), and a zero-swirl, uniform-velocity jet. The propellers were modeled using an Actuator-Line model in an unsteady Reynolds-Averaged Navier–Stokes simulation, where the Reynolds number is ReL = 3.1× 108 using the freestream velocity and body length. The authors previously showed good comparison to experimental data with this approach. Visualization of vortical structures shows the helical paths of blade-tip vortices from the single propeller as well as the complicated vortical interaction between contra-rotating blades. Comparison of instantaneous and time-averaged fields shows that temporally stationary fields emerge by half of a body length downstream. Circumferentially-averaged axial velocity profiles show similarities between the single propeller and CRP in contrast to the jet configuration. Swirl velocity of the CRP, however, was attenuated in comparison to that of the single propeller case. Mixed-patch contour maps illustrate the unique temperature distribution of each configuration as a consequence of their respective swirl profiles. Finally, kinetic and potential energy is integrated along downstream axial planes to reveal key differences between the configurations. The CRP configuration creates less potential energy by reducing swirl that would otherwise persist in the near wake of a single-propeller wake.


Introduction
Experiments show that propeller-driven wakes evolve from a complicated near wake with discernible propeller-blade features, to a far wake, in which these features have mixed together to form a nearly-axisymmetric field [1,2].Sirviente and Patel [3] show that the near-wake region transitions to the far wake in roughly twelve initial wake diameters, but the development of the far wake can be delayed by appendages on the body [4].This transition is influenced by the Reynolds number, body geometry, and operation of the propulsor [5], which itself has a large impact on the ingested stern boundary layer and downstream turbulence [6,7].The swirling propeller induces helical vortices that are shed from the roots and tips of the individual blades.In the near wake, these vortices break down, which is a topic of extensive study [8].Although experiments show the contribution of swirl [9], its role in the evolution from near to far wake is not well-characterized.
In a stratified wake, a mixed patch is formed by swirl from the propeller, turbulent mixing, and potential effects from the upstream body [10].This mixed patch can further modify the far wake in the event of a mixed-patch collapse when buoyancy forces are large [11][12][13][14].Numerous experiments have explored the interaction between stratification and wake evolution with close observation to the generation of internal gravity waves and coherent structures [15,16].Direct Numerical Simulation (DNS) provides further insight into the physics of the flow, particularly with its turbulence properties [17,18].Background turbulence increases the turbulent kinetic energy and energy transfer in the wake, which in turn lowers the mean velocity and increases horizontal spreading [19].Excess momentum leads to changes in increased turbulent kinetic energy and qualitative changes in the wake dynamics, particularly in downstream vortical structures [20].High levels of stratification in the wake create a non-equilibrium region in which the mean velocity decay is reduced [21,22].By reducing the level of potential energy in the near wake thermal-haline distribution, the effects of buoyancy in the far wake can be reduced.
Originally studied as a disc-with-center-jet [23] and later with self-propelled axisymmetric bodies [24,25], the net-zero-momentum wake functions as a theoretical model of a self-propelled marine vehicle.Beyond experiment, the study of self-propelled wakes includes several numerical methods.Ordered by increasing fidelity and computational expense, these methods include [26]: panel/lattice methods, actuator models, and fully resolved rotating geometry.Generalized actuator models include the Actuator Disk (AD), Actuator Line (AL), and Actuator Surface (AS) models.Each of these models imposes a body force over a volume in a Computational Fluid Dynamics (CFD) simulation to simulate the effects of a propeller on the surrounding fluid.Although a fully resolved propeller may offer more fidelity, its computational requirements are often large, so an actuator model provides a cost-effective alternative [27].
In a self-propelled near wake, the mixed-patch structure and overall potential energy depend largely on the propulsor.A single propeller will mix fluid unopposed within the swirling region of the wake.Contra-rotating propellers of equivalent thrust will modify the initial swirl profile thereby changing the shape of the downstream mixed patch and reducing its potential energy.Contra-rotating propeller blades add additional complexity to the interaction between root and tip vortices and reduce the swirling kinetic energy of the wake.These influences on the near wake may be compared to the simplified case of a zero-swirl, jet-propelled configuration with uniform-velocity, which results in the smallest generation of the potential energy in the wake.
The present study is an extension of Jones and Paterson [28].The unsteady Reynolds-Averaged Navier-Stokes (URANS) equations are solved to examine the near-wake evolution of the stratified, turbulent, net-zero-momentum propeller wake of the axisymmetric Iowa Body using three different propulsion schemes: single propeller, dual contra-rotating propellers (CRP), and a zero-swirl, uniform-velocity jet.The propellers are simulated using the AL model.The Iowa Body hull geometry is chosen for comparison to the non-stratified experiment by Hyun and Patel [2], which is the only known experiment to have phase-averaged propeller data for a self-propelled axisymmetric body.The authors have previously shown good agreement to this experiment for towed and self-propelled configurations [27].Flow visualization reveals the interaction between propeller-root and tip vortices and the additional complexity introduced by CRP.Comparison between instantaneous and time-averaged cross-plane profiles demonstrates the transition from near-to far-wake regions.Circumferentially-averaged profiles of velocity reveal the evolution of momentum, with observations drawn in comparison to the theoretical disc-with-center-jet that is often used in far-wake simulations [11].Mixed-patch velocity and temperature-deviation cross-plane profiles show the structure of kinetic and potential energy in the developed wake.Finally, the relative growth, decay, or persistence of integrated kinetic and potential energy of each propulsion scheme is considered.Compared to the single-propeller configuration, the CRP configuration is more effective at reducing potential and swirling kinetic energy in the wake, with potential energy reductions similar to that of the zero-swirl jet.

Governing Equations
This fluid-flow problem is defined by the unsteady Reynolds-averaged Navier-Stokes (RANS) equations in Boussinesq form with an additional body force term f p to account for the propeller model.
The equations are written in terms of the non-inertial velocity U i .In the equations, t is time, ν is the kinematic viscosity, and ρ is density.The density is expressed as ρ = ρ 0 + ∆ρ, where ρ 0 is a reference value and ∆ρ is the deviation from that value.The gravitational vector g j points downward in the negative z direction, where z is the upward-positive, vertical position.This formulation includes the piezometric pressure, p = p − ρ 0 gz where g is the magnitude of the gravitational vector.
The governing equations are solved using a custom solver written with the CFD framework OpenFOAM.This custom solver takes into account salinity and temperature transport and the corresponding turbulent fluctuations.The transport of temperature T and salinity S in the stratified environment are determined through the following equations with the diffusion coefficients κ T and κ S .

∂T ∂t
∂S ∂t The Reynolds stresses u i u j and turbulent fluxes u j t and u j s are determined using a linear eddy-viscosity closure model.

−u
In these equations, ν t is the eddy viscosity, S ij is the mean rate of strain, and k is the turbulent kinetic energy.For this study, the k − ω SST turbulence model is chosen to compute ν t due to its ease of implementation and relative advantage in computing the attached flow over a body [29].Production terms in the k − ω equations are modified to include buoyancy effects, but in the near wake they are small in comparison to the production due to shear.Wall functions are used in the computation of k and specific turbulence dissipation ω at wall boundaries to relax mesh requirements near the hull in the high Reynolds-number flow.
Density is computed by solving the UNESCO seawater equation of state [30].For the given problem, it is appropriate to approximate the secant bulk modulus as constant at sea-level conditions, even though it is a function of salinity, temperature, and pressure.Thus, the secant bulk modulus is K(S, T, p) = K(0, 20, p atm ) where p atm is atmospheric pressure.Additionally, substituting the hydrostatic pressure for the total pressure, the equation of state becomes, where a n , b n , c n and d 0 terms are empirical coefficients given in Table 1.Because the environment in the present study is isohaline, only changes in temperature from the thermally-stratified background affect changes in density.

Kinetic and Potential Energy
The evolution and transfer of energy in the wake is examined in the form of kinetic and potential energy defined as, The per-unit-volume energy ke, and pe may be integrated over an axial slice of area A in the wake to find energy per-unit-length, KE and PE, as functions of downstream distance.Kinetic energy is computed for the magnitude of velocity and also individually for each component of velocity in cylindrical coordinates.The potential energy per-unit-volume pe follows the formulation of Holliday and McIntyre [31].

Actuator-Line Model
The unsteady propeller for each non-BOR hull form is simulated using an AL model from the Simulator fOr Wind Farm Applications (SOWFA) library [32].The AL model projects a distributed line of force f p in the place of each propeller blade, where F p is the actuator element force composed of contributions from lift F L and drag F D .The distance between CFD cell center and actuator point is r, and ε controls the Gaussian width.This function decays to 1% of its maximum value when ε = 2.15r.If ε is too small, numerical oscillations arise, and if ε is too large, the applied body forces will be smoothed considerably.Troldborg [33] recommends ε ≡ 2∆x where ∆x is the grid spacing at the actuator position.Martınez et al. [34] developed best practices for AL modeling and suggested ε > 2∆x.For the present study, ε ≡ 4∆x was selected because it eliminated the numerical instabilities that arose when ε ≡ 2∆x was assigned.The Cartesian-mesh region of the propeller was refined to a resolution such that, ∆x/∆b = 0.74, where ∆b is the width of each hydrofoil section.Martınez et al. [34] suggests a value smaller than 0.75.Lift and drag at each section are computed from a lookup table of lift and drag coefficients C and C d as functions of α, where ρ is the density, U rel is the local flow speed, c is the chord and w is the width of the actuator section.The relationship between C and C d with α must be predetermined from experiment, simulation, or theory for each hydrofoil section.Figure 1 shows the magnitude of the projected propeller body force | f p |/(ρ 0 R p rps 2 ) on the AL propeller plane of the single-propeller case, where R p is the propeller radius and rps is the propeller rotations per second.

Iowa Body
The axisymmetric Iowa Body, described in the experiment by Hyun and Patel [2], is shown in Figure 2 for the standard, single-propeller case.This geometry is representative of a typical marine vehicle without appendages.Features of this geometry are listed in Table 2 where L is the body length, D is the body diameter, D p is the propeller diameter, and D h is the hub diameter.
Minor modifications are made to the Iowa Body hull for the CRP and jet configurations.For the CRP configuration, the hub is extended by the length of the rotating portion so that a second propeller may be placed directly downstream of the first.For the jet configuration, the hub is truncated at the propeller plane to function as an exhaust port.In effect, the zero-swirl, uniform-velocity jet exhausts with an initial diameter of D h .

Iowa Body Propeller
The Iowa Body propeller is defined by 36 discrete sections to account for variations in radial propeller-blade geometry.Sectional lift C is computed using the analytic expression of Brockett [35] for the NACA 66-Modified foil, where α is the local flow angle of attack, τ is the maximum thickness ratio, and f is the maximum camber ratio.Sectional drag C d is imposed by combining viscous and induced drag at each section, where C d0 is the viscous drag, e is the efficiency factor, and AR is the aspect ratio.
For the present unsteady simulations, α at each section of the propeller blades remains below 3 • at every instant in time.Because α remains small, these analytic expressions do not require additional conditions for stall.Pitch, chord, thickness, and camber distributions for the Iowa Body propeller blade are tabulated in Hyun [1].The Iowa Body propeller has zero rake and zero skew.

Computational Mesh
The three computational meshes were generated using the software cfMesh [36].Cells are focused near the body, the propulsor region, and in the wake.The hull is located at a depth of one body length.The inlet, outlet, and far-field boundaries are located two body lengths away from the hull.Comparison to simulations from spatially-larger meshes showed that the boundaries of the computational domain did not affect the solution.Mesh design and quality features are listed in Table 3.Because wall functions are used in the computation of turbulence variables, the dimensionless wall distance requirement of y + < 100 can capture the boundary layer effects and the viscous drag of the hull even for boundary cells where y + ≈ 100.A grid-refinement study of the propeller-and wake-region cells showed that 100 cells/D p adequately resolved the AL model and downstream wake cross-plane profiles.These meshes are also visualized in Figure 3. Cutting planes reveal the distribution of cells surrounding the hull and in the wake region.Views of the propulsor region show how the mesh is modified for each configuration.A single AL-modeled propeller is implemented within the highlighted region for the single-propeller case.For the CRP configuration, the hub is extended with one AL-modeled propeller placed behind the first.There is no AL model for the jet configuration since there is no propeller.Instead, fluid is exhausted from the truncated hub.

Numerical Methods
The Navier-Stokes unsteady mass and momentum equations are solved using the Pressure-Implicit with Splitting of Operators (PISO) method [37].This segregated approach decouples operations on pressure and velocity variables.At each time step, the following procedure is followed in the customized OpenFOAM solver.First, the momentum equations are solved to provide velocity by using pressure from the previous time step.Next, the pressure-Poisson equation is solved iteratively with corrections to velocity to conserve mass.Three inner iterations are are used in the present study, each with an additional mesh non-orthogonality correction step.After completion of these inner iterations, turbulence quantities are solved for, followed by salinity and temperature.The time step is then advanced.
Implicit, second-order, backward differencing is used in temporal discretization, while the cell-centered finite volume method is used in spatial discretization.A second-order, linear-upwind scheme is applied to the advective term of the momentum equations.A first-order, upwind scheme is applied to turbulence quantities, and a second-order, linear scheme is applied to all other divergence terms.Laplacian terms are discretized using a second-order, linear scheme that is partially-limited to correct for mesh non-orthogonality.
Two iterative methods are employed to solve the resulting systems of algebraic equations.The pressure equation is solved using the Preconditioned Conjugate Gradient (PCG) method with a residual tolerance of 10 −6 .The momentum, scalar transport, and turbulence equations are solved using the Pre-Conditioned Bi-Conjugate Gradient (PBiCG) scheme with a residual tolerance of 10 −8 .

Initial and Boundary Conditions
Several boundary conditions are employed.Velocity at the inlet is set to the freestream velocity U 0 through a Dirichlet boundary condition.The no-slip condition is set on the hull boundary, and the slip condition is set in the far field.Zero-gradient conditions are specified for velocity and pressure in the outlet.Background turbulence values of k and ω are computed assuming a turbulence intensity of 1% and eddy viscosity ratio ν t /ν of 100.Turbulence variables on the hull boundary are computed with wall functions.Other variables satisfy the zero-gradient Neumann boundary condition.
Initial conditions for pressure and velocity are computed by solving the potential flow equations.The PISO algorithm is then used in the transient simulation.Distributed body forces from the propeller-blades rotate at each time step at the propeller's rotation rate.The simulation is then run until initial-transient flow features advect far downstream and a periodic wake flow field is found.

Flow Field Analysis
This study examines primary flow variables including: deviation of temperature from the background ∆T and the axial, radial, and azimuthal velocities U x , U r , and U θ , respectively.The second invariant of the velocity gradient tensor Q is computed to visualize vortical structures.The cross-plane-integrated kinetic and potential energies are examined, where the kinetic energy is considered exclusively for each of the three components of velocity KE x , KE r , and KE θ .Data are extracted in axial planes within 0.9755 ≤ x/L ≤ 1.5 where x is the downstream distance from the bow of the body and L is the body length.

Flow Coefficients and Case Studies
Several of the important flow coefficients for this propeller-driven flow are the Reynolds number Re L , advance ratio J, thrust coefficient C T , and torque coefficient C Q .An alternate expression for the thrust coefficient C T * is computed for comparison to the jet configuration.
For these expressions, U 0 is the freestream velocity, ν is kinematic viscosity, D p is the diameter of the propeller, R is the radius of the Iowa Body, n is the propeller speed in revolutions per second, F T is the thrust, and F Q is the torque.The thrust-to-drag ratio is F T /F D .The Reynolds number for this study is Re L = 3.1 × 10 8 , a typical operating condition in the ocean.Other coefficients are listed in Table 4.The fore and aft propellers are listed individually for the CRP case, and total thrust is equivalent for all cases.The Froude number Fr provides a measure of the density stratification, where an infinite Fr means zero stratification, and a small Fr means high levels of stratification.The present study considers a linearly varying temperature stratification, typical of an ocean environment, with a Froude number of Fr = 350, where, In these expressions N is the Brunt Väisälä frequency, g is acceleration due to gravity, and z is the vertical coordinate.The influence of buoyancy on the near-wake fluid dynamics is often small and can be quantified by the Richardson number Ri, which is the ratio of buoyancy to flow gradient terms [15].
where T 0 is a reference temperature, T is the mean temperature, and U is the mean velocity.For the single-propeller case, Ri ≈ 2.54 × 10 −3 which indicates that the near-wake inertial forces of the propeller dominate the buoyancy forces.As the local velocity U x decays, Ri increases and buoyancy forces become important further downstream in the far wake, beyond the geometric bounds of these simulations.

Near-Wake Transition
Individual vortices are visualized using the second invariant of the velocity gradient tensor Q. Figure 4 shows contour surfaces of the non-dimensionalized (L/U 0 ) 2 Q = 16.9 with a vertical cutting plane colored by axial velocity defect U x /U 0 − 1 that extends to half of a body length downstream where x/L = 1.5.For the single-propeller case, root and tip vortices induced by the propeller are apparent.These vortices follow a helical path and disappear by x/L ≈ 1.25.For the CRP case, additional complexity is introduced by the interaction between the two propellers.Complicated vortical structures are visible until x/L ≈ 1.35.Negligible vortical structures are found in the case of the zero-swirl jet.This figure illustrates the complexity introduced by contra-rotating propellers in comparison to the other two cases.
To better understand the transition from near to far wake regimes in the propeller-driven cases, axial planes behind in the propulsor are examined.Figure 5 compares the instantaneous and time-averaged axial velocity defect field U x /U 0 − 1 for the single propeller case.Instantaneous fields are taken at an arbitrary time long after initial-transient features have disappeared from the simulation and the flow field has become periodic.Time-averaging occurs over temporal interval of two periods of the propeller.Near the propellers, at x/L = 1.01 and x/L = 1.10, individual propeller-blade wakes can be seen in the instantaneous field.These blade wakes follow the azimuthal motion of the propellers directly upstream.The time-averaged field, by contrast, is axisymmetric.Further downstream at x/L = 1.3,only a small variation is seen between the instantaneous and time-averaged fields, and by x/L = 1.45, the two contour maps are nearly identical showing that the wake is steady and axisymmetric.By half of a body length downstream, the flow is stationary in time and space when viewed from a body-fixed reference frame.
For the CRP case shown in Figure 6, similar observations are drawn.Unsteadiness is apparent for x/L ≤ 1.10, however by x/L = 1.45 the cross-plane profile is temporally stationary.In this case, a unique hexagonal shape is formed due to the interaction between the two opposing three-bladed propellers.This shape is still present at half of a body length downstream of the propulsor.

Velocity Profiles
The evolution of the near wake may also be described by circumferentially-averaged velocity profiles.Figure 7 shows the circumferentially-averaged velocity defect profiles U x /U 0 for the three self-propelled configurations at the downstream positions x/L = 1.01, x/L = 1.3, and x/L = 1.5.Just behind the propulsor at x/L = 1.01, the jet is shown to have uniform, positive velocity leaving the exhaust port, while negative velocity due to drag appears for r/R p > 0.2.By x/L = 1.3 and further at x/L = 1.5, the jet profile appears as a classical net-zero-momentum wake and may be described using the analytical formulation of the disc-with-center-jet.
Circumferentially averaged profiles of the swirl component of velocity U θ are shown in Figure 8.The jet configuration profile exhibits zero swirl because there are no sources of swirl for this case.The jet exhaust contains uniform axial velocity, the body is axisymmetric, and buoyancy forces are relatively small.The single-propeller configuration at x/L = 1.01, shows swirl imparted by the propeller and rotating hub.By x/L = 1.3 and x/L = 1.5, most of the momentum due to swirl exists in a region centered at r/R p ≈ 0.4.In the case of the CRP, regions of positive and negative swirl develop due to interference between the opposing, contra-rotating blades.Throughout the near wake, the CRP swirl magnitude is attenuated, remaining less than half of that of the single propeller, which is explained by the interactions of opposing azimuthal forces of the two propellers.x/L =1.30 x/L =1.50 single CRP jet The propeller-driven cases show their own unique profiles.Positive momentum from the propellers exists in a region near r/R p ≈ 0.6, while negative momentum due to drag from the body exists near the center and further outward.These two circumferentially-averaged profiles are nearly equivalent and decay at similar rates, which is explained by the similar distributions of axial momentum.They also may be defined analytically using a process described in Jones and Paterson [27].For all configurations, the positive momentum decays more quickly than the negative momentum, which is a feature of idealized wakes described by Tennekes and Lumley [38].The theoretical, axisymmetric drag wake decays with the power of −2/3, while the theoretical, axisymmetric jet decays with the power of −1.

Velocity and Temperature Fields in the Mixed-Patch
The wake is mixed by half of a body length downstream of the stern.Unsteadiness from the propulsor has disappeared and axial gradients are small in comparison to the transverse.This location is significant because a cross-plane profile may be considered as the initial data plane (IDP) for further far-wake simulations, which is beyond the scope of this paper.Unique features of the IDP cross-plane profiles of and temperature deviation are presented.

Velocity Field
Cross-plane contour maps of the axial velocity defect U x /U 0 − 1 at x/L = 1.5 are shown for the single propeller and CRP cases in Figure 9.For the single propeller, axial velocity is axisymmetric and has previously been fit to an analytical curve as a function of radial distance [27].The CRP, however, is not axisymmetric, and a steady, hexagonal profile is formed.The geometric shape is attributed to the two three-bladed propellers interacting with one another.Unlike the single-propeller profile, the CRP profile is a function of both radial and azimuthal positions.To create an analytical expression as a function of radial position alone, the profile must first be circumferentially averaged.
The swirl component of velocity U θ /U 0 is shown in Figure 10.Again, the single-propeller velocity is axisymmetric, while the CRP velocity has discernible geometry.The magnitude of swirl velocity is much higher for the single propeller case, since the CRP propulsor imparts opposing azimuthal forces.Swirl velocity from the CRP is less than half of that of the single propeller and varies in azimuthal sign.

Temperature Field
Because these simulations take place in a thermally-stratified environment, any vertical redistribution of the fluid will generate potential energy.Mixing in the wake plays an important role in the redistribution of temperature T. Given that the background temperature field T b is initially linearly stratified, mixing from the wake develops a temperature deviation ∆T = T − T b .This field is non-dimensionalized by the linear change in temperature over the depth of one propeller-blade length ∆T R p .
Figure 11 shows ∆T/∆T R p for the single-propeller case.Additional radial profiles help to visualize how the field varies in polar coordinates.A unique cross-plane profile shape is formed that is steady in time.Colder fluid has been driven to the top, while warmer fluid has been driven to the bottom of the wake.A maximum is found near the center of the warm region and a minimum is found near the center of the colder region.Two "tails" are shown trailing off of the warm and cold regions as a result of the counter-clockwise swirling motion of the fluid due to the propeller.
The mixed-patch ∆T/∆T R p cross-plane profile of the CRP case is shown in Figure 12.Compared to the single-propeller case, the magnitude of ∆T/∆T R p is less than half.The profile is split between an inner region where clockwise-swirling fluid dominates and an exterior region where counter-clockwise-swirling fluid dominates as shown previously in Figure 8.The complexity and lower ∆T/∆T R p magnitude in the profile arise directly from the initial interactions of the opposing azimuthal forces of the contra-rotating blades that drive the swirling fluid.While a single propeller can transport the temperature field across the swirling wake region unimpeded, the addition of an opposing propeller directly counters this effect.The net-swirl in the wake is reduced and regions of both positive and negative swirl exist.The inner, negative-swirl region forms a ∆T/∆T R p profile that mirrors the single-propeller case because of the sign difference in swirl.The outer, positive swirl region shares the same sign of ∆T/∆T R p as the single-propeller case.In effect, the net loss in swirl reduces the overall potential energy.Finally, the mixed-patch ∆T/∆T R p cross-plane profile of the jet case is shown in Figure 13.For this configuration, the magnitude of ∆T/∆T R p is the smallest due to the absence of swirl.Instead, the shearing of axial momentum and potential effects from the body control the shape of the profile.The positive-momentum "jet" core entrains fluid from the negative-momentum "drag" periphery of the wake.Potential effects from the body further influence the temperature distribution.The resulting transport of ∆T/∆T R p suspends warmer fluid above colder fluid in the center, and the reverse in sign in the periphery.Without swirl, the distribution of ∆T/∆T R p in the central core of the jet wake is opposite in sign to that of the single-propeller case.

Comparison to a Perfectly-Mixed Temperature Field
The simulated cross-plane profiles may be compared to the idealized, perfectly-mixed profile shown in Figure 14.This conceptual profile assumes perfect mixing such that the temperature distribution T is uniform up until the boundary of the wake disc, beyond which T = T b .The temperature deviation ∆T = T − T b , however, varies because of changes from the background stratification within the disc.Given the linear background stratification, ∆T also varies linearly in the vertical direction.The single propeller case relates most closely to this idealized profile, but, because it does not perfectly-mix T, differences can be observed.The upper and lower regions of colder and warmer fluid are shifted from the centerline, and the maximum temperature deviations are not on the wake boundary but instead closer inward.Results from the single propeller case show that T is not mixed uniformly within the disc of the swirling wake.

Potential and Kinetic Energy Evolution
Kinetic and potential energy is integrated along axial planes downstream for the three configurations, as shown in Figure 15.Kinetic energy KE is computed individually for the three components of velocity, namely radial, swirl, and axial as KE r , KE θ , and KE x , respectively, as well as for the velocity magnitude, KE.Downstream distance is described both by x/L measured from the bow and by x /D p measured from the stern.For all three cases, KE r decays more rapidly than KE x and KE θ .The swirl component notably exhibits the slowest decay in the near wake, an observation consistent with Sirviente and Patel [9].For the single-propeller case, this relative persistence leads to a rise in PE due to expansion of the wake and entrainment of the surrounding passive scalar T, indicating a change in density.For the CRP, the PE does not grow due to the opposing regions of positive and negative swirl velocity.Instead, the PE decays with a rate similar to the zero-swirl jet.This result shows that the contra-rotating blades can effectively reduce PE, which will reduce the strength of buoyancy effects in the far wake.
Comparing the swirl component for the three cases, KE θ of the CRP is an order of magnitude lower than that of the single propeller.Counteracting azimuthal momentum leads to a reduction in the swirl kinetic energy.The reduced KE θ of the CRP is consistent with its reduced PE.The jet KE θ is small because of the initial absence of swirl.
Additionally, comparison of the total KE shows that the CRP is the most-effective configuration for the total reduction of energy by x/L = 1.5.This increased decay rate suggests that the CRP far wake will decay more quickly than the single propeller and jet cases.The single propeller is less effective than the jet due to persisting KE θ from its unidirectional swirling velocity.

Conclusions
The influence of swirl on the evolution of self-propelled, stratified near wakes and the development of the mixed patch has not previously been well-characterized.In this study, the linearly stratified near wake of the Iowa body was investigated with three separate propulsor configurations: single propeller, contra-rotating propellers, and a zero-swirl, uniform-velocity jet.Unsteady, rotating propeller blades were simulated using an AL model in a URANS computation.Comparison between the configurations revealed unique differences in the evolution of the near wake.
While clear root and tip vortices were visible in the single-propeller case, the CRP disrupted these structures, introducing additional complexity in the wake evolution.Nevertheless, by half of a body length downstream, the wake flow fields were steady in time.The single-propeller and CRP cases shared similar circumferentially-averaged axial velocity defect profiles due to similar spanwise loading in the propulsor.Swirl velocity, however, varied between the two propeller-driven cases, with the CRP introducing both positive and negative swirl regions exhibiting half of the magnitude of the single-propeller case.Furthermore, by half of a body length downstream, the magnitude of the temperature deviation ∆T/∆T R p for the CRP was less than half of that of the single propeller.The jet ∆T/∆T R p magnitude was the smallest, due to the absence of swirl.Contour maps of velocity revealed that the single propeller has an axisymmetric profile, whereas the CRP exhibits a unique hexagonal structure as a result of its two three-bladed propellers.The evolution of kinetic and potential energy varied as a direct result of the swirl imparted by each propulsor.Because of the interaction of positive and negative swirl, the CRP configuration showed an order of magnitude lower swirling kinetic energy compared to the single propeller configuration.Additionally, its potential energy was similar in decay and magnitude to that of the swirl-free jet, and the total kinetic energy decayed most rapidly out of the three propulsion schemes.These results indicate that the CRP can effectively reduce potential energy that would otherwise develop from a single-propeller configuration.By removing potential energy, buoyancy effects in the far wake will be weakened.

Figure 1 .
Figure 1.Non-dimensional body force on mesh slice at propeller plane | f p |/(ρ 0 R p rps 2 ) for the single-propeller case.
(a) Vertical cutting plane through standard Iowa Body mesh with nested refinement.(b) Boundary-layer cells.(c) Propulsor region for standard mesh.(d) Propulsor region for CRP mesh.(e)Propulsor region for jet mesh.

Figure 3 .
Figure 3. Computational meshes generated for each configuration.

Figure 7 .
Figure 7. Circumferentially-averaged velocity defect U x /U 0 − 1 profiles for each configuration at various distances downstream.

Figure 8 .
Figure 8. Circumferentially-averaged swirl velocity U θ /U 0 profiles for each configuration at various distances downstream.
Profile in radial directions.
Profile in radial directions.
Profile in radial directions.

Figure 14 .
Figure 14.Idealized mixed-patch profile of ∆T/∆T R p for a wake region of constant T.

Figure 15 .
Figure 15.Integrated energy evolution downstream of the vehicle for each configuration, with x measured from the bow of the hull and x measured from the stern.

Table 1 .
Coefficients in UNESCO equation of state for seawater.

Table 3 .
Mesh design and quality features.