Numerical Modeling of a Marine Propeller Undergoing Surge and Heave Motion

Abstract

A boundary element method (BEM) and a vortex-lattice method (VLM) are extended in order to predict the unsteady performance of propellers subject to rigid body motions. The methods are applied in the case of prescribed surge and heave motions, and the results are compared with those from other methods.

1. Introduction

Due to the unsteady motions of a ship, the operating propeller is usually subjected to a time-dependent rigid body motion which can drastically change the hydrodynamic performance of the propulsion system and reduce the maneuverability of the ship. Incident waves tend to move the hull in 6 DOF over a range of frequencies which in turn oscillate the propeller in all possible degrees of freedom. In general, a propeller operates in time-varying wake field and thus generates vibratory forces causing severe wear and tear of the propeller system. For all such cases, it becomes important to consider the effect of unsteady contributions to the wake field caused by waves and wave-induced motions. The combined effect of these two contributions was investigated in the past by Sluijs [1] and by Jessup and Boswell [2]. Later, efforts were made by Parsons and Vorus [3] to assess added mass and damping coefficients for a vibrating propeller. Sasajima [4] adapted a quasi-steady approach to estimate propeller bearing forces while Breslin and Andersen [5] in favored fully unsteady approach as opposed to quasisteady approach for calculating time-varying forces on propellers.

In the recent past, Lee and Kinnas [6] developed unsteady wake alignment model to simulate propeller in nonaxisymmetric flows. Politis [7] had modeled the unsteady motion of propeller using a free wake model. Very recently, He [8] developed a full unsteady wake alignment model for hydrofoils and marine propellers which was incorporated by Sharma et al. [9, 10] in modeling the unsteady motion of hydrofoils and marine propellers.

In this paper, a computational approach to analyze the hydrodynamic performance of a marine propeller undergoing unsteady motion is presented. Firstly, an unsteady boundary element method (BEM) model for a general rigid body motion problem is proposed. The kinematic and dynamic boundary conditions of the BEM model are also incorporated into a Vortex Lattice Method (VLM) code. The BEM model is first applied to a pitching hydrofoil followed by validation of numerical results of the current model against those obtained from a commercial RANS solver. The problem of flow around a marine propeller performing surge and heave motion in unbounded fluid is then solved using both VLM and BEM. A fully unsteady wake alignment algorithm is implemented into the VLM code in order to satisfy the force-free condition on the propeller wake surface. Finally, a comparative study of unsteady propeller forces obtained from BEM/VLM, RANS, and other available numerical results is carried out and results are presented. The ultimate goal of the present research work is to develop an integrated 6 DOF propeller-hull motion model by including complete propeller-hull interaction in time. A comprehensive analysis of the flow around such system under most severe motions can be performed by the proposed model.

2. Unsteady BEM Model for an Obstacle with Rigid Body Motion

As shown in Figure 1, an obstacle is subject to rigid body motion under an inertial frame of reference. Consider a fluid particle which is at location 𝐱_𝟏 when 𝑡=𝑡₁. After a small interval, when 𝑡=𝑡₂=𝑡₁+Δ𝑡, the particle is moved to a new location 𝐱_𝟐. At the same time, the point on the body which was at 𝐱_𝟏, when 𝑡=𝑡₁ is moved to a new location 𝐱^′_𝟏, due to rigid body motion. Two velocities, 𝐔_𝐛 and 𝐪, are associated with the movement of the point on the body and the fluid particle, where 𝐔_𝐛 is the velocity due to the rigid body motion, which translated the point on the body from 𝐱_𝟏 to 𝐱^′_𝟏, and 𝐪 is the total velocity which moved the fluid particle from 𝐱_𝟏 to 𝐱_𝟐.

Figure 1 

Schematic plot of an obstacle with rigid body motion.

The velocity flow field can be decomposed into two components: inflow velocity and perturbation velocity, written as 𝐪=𝐔_∞+𝐮, (1) where 𝐪 is the total velocity, 𝐔_∞ is the inflow velocity, and 𝐮 is the perturbation velocity due to the presence of obstacles.

2.1. Governing Equations

The methodology applied in this study is based on the inviscid potential flow theory: 𝐮=∇𝜙, (2) where 𝜙 is the perturbation potential. The perturbation velocity field, therefore, is governed by the Laplace Equation: ∇²𝜙=0. (3) If the inflow is also irrotational, we can also define the total potential Φ as 𝐪=∇𝜙+𝑈_∞=∇Φ. (4) Appling Green’s second identity, the Laplace equation ‎(3) can be written as following boundary integrated form: 2𝜋𝜙_𝑝=_𝑆𝐻𝜙_𝑝𝜕𝐺𝑝,𝑝^′𝜕𝐧_𝑝′−𝜕𝜙_𝑝𝜕𝐧_𝑝′𝐺𝑝,𝑝^′+𝑑𝑆_𝑆𝑊𝜙_𝑤𝑦_𝑝′𝜕𝐺𝑝,𝑝^′𝜕𝐧_𝑝′𝑑𝑆, (5) where 𝑆_𝐻 represents the surface of hydrofoil or propeller blades and 𝑆_𝑊 represents the surface of trailing edge wake. 𝐺(𝑝,𝑝^′) is the Green’s function, which is defined as 1/𝑅(𝑝;𝑝^′) in 3D, and 2ln𝑅(𝑝;𝑝^′) in 2D.

2.2. Boundary Conditions

2.2.1. Kinematic Boundary Condition

Without viscosity, the solid surface has to be treated as slip wall, which implies that the fluid particles on the surface of the body can only be sliding on the surface. From Figure 1, we can easily obtain the expression of the kinematic boundary condition on the solid surface: 𝐧⋅𝐪−𝐔_𝐛=𝟎, (6) which can be rewritten as 𝜕𝜙𝐔𝜕𝑛=−_∞−𝐔_𝐛⋅𝐧. (7) It is convenient to define 𝐔_in=𝐔_∞−𝐔_𝐛. (8) Therefore, the kinematic boundary condition can be written as 𝜕𝜙𝜕𝑛=−𝐔_in⋅𝐧. (9)

2.2.2. Dynamic Boundary Condition

The momentum equation for potential flow can be integrated into the form of unsteady Bernoulli’s Equation: 𝑃_∞𝜌+12𝐔_∞⋅𝐔_∞=𝑃𝜌+𝐪⋅𝐪2+𝜕𝜙,𝜕𝑡 (10) where 𝑃_∞ is the reference pressure infinitely far upstream, and 𝑃 is the pressure in the field. It is more convenient to treat the unsteady term in (10) with the material derivative by using following identity: 𝐃𝜙=𝐃𝑡𝜕𝜙𝜕𝑡+𝐪⋅∇𝜙. (11) Equation (10) becomes 𝑃_∞−𝑃𝜌+12𝐔_∞⋅𝐔_∞=𝐪⋅𝐪2+𝐃𝜙𝐃𝑡−𝐪⋅∇𝜙. (12) Following the definition in Figure 1, as Δ𝑡→0, we have 𝐃𝜙=𝜙𝐱𝐃𝑡^𝟐,𝑡₂−𝜙𝐱_𝟏,𝑡₁=𝜙𝐱Δ𝑡_𝟐,𝑡₂−𝜙𝐱^′_𝟏,𝑡₂+𝜙𝐱Δ𝑡^′𝟏,𝑡₂−𝜙𝐱_𝟏,𝑡₁𝑥Δ𝑡=∇𝜙⋅²−𝑥^′₁+Δ𝑡𝛿𝜙𝛿𝑡=∇𝜙⋅𝐪−𝐔_𝐛+𝛿𝜙,𝛿𝑡 (13) where 𝛿𝜙/𝛿𝑡 denotes the rate of change of the perturbation potential associated with a solid particle on the surface of the body.

Substituting (13) into (12), we obtain the dynamic boundary condition on the solid surface: 𝑃_∞−𝑃𝜌+12𝐔_∞⋅𝐔_∞=𝐪⋅𝐪2−∇𝜙⋅𝐔_𝐛+𝛿𝜙.𝛿𝑡 (14) Adding 0.5𝐔_𝐛⋅𝐔_𝐛−𝐔_∞⋅𝐔_𝐛 to both side of the equation above will lead to 𝑃_∞−𝑃𝜌+12𝐔_in⋅𝐔_in=𝐪⋅𝐪2+𝐔_𝐛⋅𝐔_𝐛2−𝐔_∞+∇𝜙⋅𝐔_𝐛+𝛿𝜙,𝛿𝑡 (15) which is 𝑃_∞−𝑃𝜌=12𝐔_in+𝐮²−12𝐔²_in+𝛿𝜙.𝛿𝑡 (16) Equation (16) is a general expression of the unsteady Bernoulli’s equation for an obstacle with any rigid body motion. For a 2D hydrofoil with pitch motion, 𝐔_𝐛 is Ω×𝐫, where Ω angular velocity of pitching, and 𝐫 is the vector from the center of pitching to the points on the surface of the hydrofoil, as shown in Figure 2.

Figure 2 

A hydrofoil with pitch motion.

As shown in Figure 3, in the case of a propeller under surge and heave motion, 𝐔_𝐛 and 𝐔_in are defined as 𝐔_𝑏=𝝎_prop×𝐫+𝐔_𝑏,surge+𝐔_𝑏,heave,𝐔_in=𝐔_∞−𝝎_prop×𝐫−𝐔_𝑏,surge−𝐔_𝑏,heave, (17) where 𝝎_prop is the angular velocity of the rotation of the propeller; 𝐫 is the radial vector from the axis of the propeller to the points on the blade surface; 𝐔_𝑏,surge and 𝐔_𝑏,heave are the surge and heave velocity of the propeller, respectively. The surge and heave motion of a propeller are incorporated into a BEM code, PROPCAV. It has to be pointed out that the formulation presented above is based on the assumption of inviscid and unseparated flow. If the amplitude of motion of the moving hydrofoils or propellers becomes too large, the boundary layer may separate even at the leading edge of the hydrofoil or the blade of propeller and then form leading edge vortex, which is beyond the scope of application of current models.

Figure 3 

A propeller with surge and heave motion.

3. Unsteady VLM Model for Propeller with Rigid Body Motion

The kinematic and dynamic boundary conditions described for the BEM model can be easily applied to vortex lattice method (VLM). In this study, a VLM code, MPUF3A, is modified to incorporate the correct kinematic and dynamic boundary conditions for a propeller with either surge or heave motions. He [8] developed a fully unsteady wake alignment model based on MPUF3A. The effects of the unsteady wake model on propellers undergoing heave motions are also investigated. The details on VLM and the unsteady wake alignment can be found in He [8].

4. Application to 2D Pitching Hydrofoils

The BEM model presented in this study is first applied to a pitching hydrofoil, which has an NACA66 thickness form with 5% thickness and 0% camber over chord length. The angular frequency, 𝑘 of the pitching motion is 𝜋, and the rotation-axis origin is set as the leading edge of the hydrofoil. The angle between the hydrofoil and the inflow, or, the instantaneous angle of attack, 𝛼, is defined as follows: 𝛼=𝛼₀sin𝑘𝑡, (18) where 𝛼₀ denotes the amplitude of the angle between the hydrofoil and the inflow. RANS and also inviscid finite volume method (FVM) simulations are performed with commercial software FLUENT for the same configuration. The inviscid FVM simulation is carried out in order to assess viscous effects. The standard 𝑘-𝜀 turbulence model with standard wall function is used to solve RANS equations. Two-dimensional unstructured grids were generated with structured prism boundary layer cells attached on the hydrofoil. In total, 407K number of cells is distributed over the entire domain for pitch motion while 390K cells are provided for heave motion. Obtained range of 𝑦⁺ for both the motions is between 40 and 120. Dynamic mesh is used to simulate the pitch motion. Computational (CPU) time of 17 hours is taken for the RANS simulation on 4 Intel Xeon 2.54 GHz CPUs whereas 2D BEM model takes only 3 minutes on a single CPU for simulating the complete flow.

Figures 4 and 5 compare the unsteady lift and drag coefficients from 2D BEM model and inviscid FVM/RANS where 𝐶_𝐷=𝐷/(0.5𝜌𝑈²𝑐). Results from both of these methods show close level of agreement. However, a slight shift between inviscid and viscous drag coefficient can be noticed. This shift is attributed to the viscous effects which are not included in the inviscid FVM simulation and 2D BEM model. The shift of the drag coefficient appears to be a constant value, which indicates that it can be modeled with a time-independent constant skin friction coefficient. In Figures 6 and 7, pressure distribution on the foil surface is compared at different flow time. The results suggest that except at the trailing edge of the hydrofoil, both 2D BEM model and RANS agree well with each other. In the unsteady BEM model, an iterative pressure Kutta (IPK) condition, proposed by Kinnas and Hsin [11], has been implemented to obtain zero loading at the trailing edge of the hydrofoil.

Figure 4 

Lift coefficient versus time.

Figure 5 

Drag coefficient versus time.

Figure 6 

Pressure distribution on the pitching hydrofoil at 𝑡/𝑇=0.3845.

Figure 7 

Pressure distribution on the pitching hydrofoil at 𝑡/𝑇=1.482.

Figure 8 shows the comparison between unsteady wake profiles obtained from 2D BEM model and RANS. Both the trajectories and rollups are predicted accurately by 2D BEM model. As shown in Figure 9, range of 𝑦⁺ is between 40 and 100, which is acceptable for standard wall function.

Figure 8 

Comparison of wake pattern for pitch motion where 𝑘=𝜋 and pitch amplitude = 0.1 rad. BEM is shown by solid curve, and RANS by colormap.

Figure 9 

𝑦+ on the hydrofoil surface.

5. Application to Propellers

Marine propellers are subjected to spatially nonuniform unsteady inflow, which causes periodic loads on the propeller blades and shaft system. The unsteadiness in the inflow is primarily responsible for blade pressure fluctuations and cavity inception whose effect on the propulsion system could be significant. In addition to this, the transmission of unsteady forces through shafting into the vehicle can cause severe wear and tear of the system. Therefore, it becomes important to analyze the effects of unsteadiness in the inflow on the propulsive performance of a propeller for an efficient design of the hull propeller system.

5.1. Propeller Undergoing Pure Surge Motion

A model propeller DTNSRDC P4381 undergoing surge motion respectively is simulated. The geometry of this propeller is available to the author through the technical report by Boswell [12]. Fully wetted unsteady runs are performed at design advance coefficients (𝐽_𝑠=0.889).

In this study, the surge motion is defined as harmonic oscillation: 𝐔_𝑏,surge=𝐔_surgesin𝜔_surge𝑡, (19) where 𝐔_surge denotes the amplitude of the surge velocity, and 𝜔_surge is the angular frequency of the surge motion. Table 1 shows the parameters of the surge motion in this case.

Table 1 

Surge motion parameters.

In VLM, the propeller blade is discretized using 20 panels in chordwise and 18 panels in spanwise direction while in BEM, discretization is done using 60 panels in chordwise and 20 panels in spanwise direction. Grids are distributed uniformly along the spanwise direction while in the chordwise direction full cosine distribution is provided. Thickness/loading coupling has been implemented in VLM to incorporate the thickness.

RANS simulation for the same propeller undergoing the same surge motion is also carried out. Modeling of surge motion of propeller in RANS solver is performed by leveraging moving mesh feature of FLUENT and using UDF to give surge motion to the complete domain. Total numbers of 1.75 million structured hexahedral cells are provided in a periodic domain. 192 hours are taken by RANS to simulate 9 seconds of flow on 24 Intel Xeon 2.54 GHz CPUs, on the other hand, VLM takes only 12 minutes and BEM takes 15 minutes of CPU time for the complete simulation.

A detailed convergence study of both the BEM and VLM codes can be found in Sharma et al. [10].

Figures 10 and 11 show the comparison between time variation of the thrust and torque coefficients among different numerical methods. BEM and VLM compare well with each other but show some differences with RANS. It is also noticed that the time average KT over one cycle of surge oscillation is growing for RANS, while it is constant for BEM and VLM. Similarly, mean KQ for RANS is growing, but it remains constant for BEM or VLM.

Figure 10 

Unsteady thrust coefficient.

Figure 11 

Unsteady torque coefficient.

5.2. Propeller Undergoing Pure Heave Motion

A model propeller DTMB P4718 undergoing heave motion is simulated. The performance of this propeller underdoing heave motion was investigated by Politis [7] independently. Fully wetted unsteady runs by VLM are performed at design advance coefficients (𝐽_𝑠=0.751). In this study, the heave motion is in the form of harmonic oscillation: 𝐔_𝑏,heave=𝐔_heavesin𝜔_heave𝑡, (20) where 𝐔_heave denotes the amplitude of the heave velocity, and 𝜔_heave is the angular frequency of the heave motion. Table 2 shows the parameters of the heave motion in this case.

Table 2 

Heave motion parameters.

In this case, the spatial discretization on the propeller blade is the same as that in the case of surge motion. The results from VLM are compared with those obtained from Politis [7], as shown in Figure 12. The unsteady single-blade thrust coefficients from the present methods appear to be fairly close to each other and also close to the Politis’ results. In the case of heave motion, the VLM code with fully unsteady wake alignment (FWA) is also tested. Although in many situations, such as inclined shaft propeller and propeller at low advance ratio, the fully aligned wake model could significantly improve the results, in this case, it does not make a big difference. The reason could be that the propeller is still working near the design condition, in which the pitch angle of the wake is relatively small.

Figure 12 

Comparison of KT from MPUF-3A and Politis [7]. Unlike other KT plots, direction of thrust is considered positive, that is, along the direction of ship’s forward motion.

6. Conclusions and Future Work

A time-dependent 2D BEM model has been developed and validated against RANS solver to analyze the nature of unsteady forces, pressure distributions, and wake profiles for an oscillating hydrofoil. Simulations of heaving and pitching motion are performed for a set of parameters to avoid inception of leading edge vortex (LEV) which can drastically change pressure distribution on the surface of hydrofoil. In the future, the LEV model developed by Tian and Kinnas [13] will be incorporated to analyze flapping hydrofoil for a wider range of parameters.

The unsteady hydrodynamic analysis of propeller under surge and heave motion has been performed using vortex lattice method (VLM) and boundary element method (BEM). A fully unsteady wake alignment algorithm is implemented in the VLM code. On comparing the results from the present method to those from Politis [7], fair-to-moderate level of agreement has been found between the two. Surge motion of a propeller has successfully been simulated using RANS solver; however, on comparing time-varying KT and KQ with those obtained from BEM or VLM, differences were found for some cases. Discrepancies observed while comparing the results from the VLM and the BEM codes need to be investigated further.

The numerical methods proposed in this study are not intended as a substitution of RANS simulations. The effects of turbulence and vortical separated flows are difficult to be handled by the methods presented in this paper. In these cases, RANS or LES could give results with more physical representations of the problem, at the same time, with significantly increased computational cost. The uncertainties in the numerical solutions have to be addressed as well. The study by Sharma [9] showed that current methods converged to the same results as refining the spatial and temporal resolution. Direct correlation of the results from current methods with any experimental data would be greatly helpful to validate the current methods. In the current study, only two simple harmonic motions are considered. In the future, it is worthwhile to analyze the propeller undergoing surge and heave motions with more harmonic components and understand the corresponding propeller performance in the frequency domain. The present method will also be extended to analyze pitch, yaw, and sway motion of the propeller.

Having a strong background in simulating unsteady motion of hydrofoils and marine propellers, a more comprehensive approach is required in dealing with the minute details of flow around the oscillating hull propeller system. The long-term objective of the current research work is to integrate the present model with hull motion model followed by coupling with Reynolds Averaged Navier Stokes (RANSs) method to develop a complete 6 DOF hull-propeller system motion model. Coupling with RANS by Kinnas et al. [14] can be improved to provide accurate propeller-hull interaction in time. In this approach, the propeller forces are represented by spatially distributed and time-dependent body forces in the RANS code to account for the effects of propeller-induced flows. In the cases that wake alignment is important, unsteady wake alignment could be implemented using the velocity field from RANS. The treatment of the different time step sizes of the hull motion and the propeller sub-system has to be investigated further. Development of the required model will analyze the complete flow around hull-propeller system and predict forces, moments, and pressure distribution with utmost accuracy.

Acknowledgments

Support for this research was provided by the U.S. Office of Naval Research (Contract nos. N00014-07-1-0616 and N00014-10-1-0931) and Phases V and VI of the Consortium on Cavitation Performance of High Speed Propulsors with the following current members: American Bureau of Shipping, Daewoo Shipbuilding and Marine Engineering Co. Ltd., Kawasaki Heavy Industry Ltd., Rolls-Royce Marine AB, Rolls-Royce Marine AS, Samsung Heavy Industries Co. Ltd., SSPA AB, Sweden, VA Tech Escher Wyss GmbH, Wärtsilä Propulsion Netherlands B.V., Wärtsilä Propulsion Norway AS, Wärtsilä Lips Defense S.A.S., and Wärtsilä CME Zhenjiang Propeller Co. Ltd.