International Journal of Rotating Machinery

Volume 2012 (2012), Article ID 474785, 11 pages

doi:10.1155/2012/474785

Research Article

Open-Water Thrust and Torque Predictions of a Ducted Propeller System with a Panel Method

J. Baltazar,¹ J. A. C. Falcão de Campos,¹ and J. Bosschers²

¹Department of Mechanical Engineering, Marine Environment and Technology Center (MARETEC), Instituto Superior Técnico (IST), Technical University of Lisbon, 1049-001 Lisbon, Portugal

²Maritime Research Institute Netherlands (MARIN), Wageningen, The Netherlands

Received 2 December 2011; Accepted 8 March 2012

Academic Editor: Moustafa Abdel-Maksoud

Copyright © 2012 J. Baltazar et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

This paper discusses several modelling aspects that are important for the performance predictions of a ducted propulsor with a low-order Panel Method. The aspects discussed are the alignment of the wake geometry, the influence of the duct boundary layer on the wake pitch, and the influence of a transpiration velocity through the gap. The analysis is carried out for propeller Ka4-70 operating without and inside a modified duct 19A, in which the rounded trailing edge is replaced by a sharp trailing edge. Experimental data for the thrust and torque are used to validate the numerical results. The pitch of the tip vortex is found to have a strong influence on the propeller and duct loads. A good agreement with the measurements is achieved when the wake alignment is corrected for the presence of the duct boundary layer.

1. Introduction

The ability to accurately predict the thrust and torque of a ducted propeller in open-water conditions is very important for a calculation method used in the design stage. RANSE methods have been progressively introduced for the calculation of ducted propeller systems, meeting considerable success in predicting open-water characteristics for the well-known Ka-series [1–3]. However, due to their relative complexity and time requirements, they are not yet routinely used in the design process, which is often still based on the use of inviscid flow methods.

Various numerical methods based on inviscid (potential) flow theory have been proposed for the analysis of ducted propellers. Examples are the combination of a Panel Method, also known as Boundary Element Method, to model the duct with a vortex lattice method for the propeller [4], and a Panel Method for the complete ducted propeller system operating in unsteady flow conditions including blade sheet cavitation [5]. Both methods applied a transpiration velocity model for the gap flow between propeller blade tip and duct inner surface and analysed a duct with a sharp trailing edge.

These references show that the application of inviscid flow models to ducted propellers, albeit of great usefulness, may meet some serious limitations related to the occurrence of flow regions where viscous effects cannot be ignored and have to be modelled in some way for the correct prediction of the ducted propeller thrust and torque. One of such region concerns the gap flow, which has a strong influence on the propeller and duct circulation distribution, and therefore, on the distribution of loading between propeller and duct, as studied in detail by Baltazar and Falcão de Campos [6]. In addition, there may be a considerable interaction between the vorticity shed from the propeller blade tips and the boundary layer developing on the duct inner side, as found in the works of Krasilnikov et al. [3] and Rijpkema and Vaz [7]. This effect has not been studied before with potential flow methods, and its importance is therefore unknown.

The purpose of this paper is to show the importance of an efficient and robust method for the vortex pitch alignment, including the effect of gap modelling and the effect of duct boundary layer. The computational results are compared with open-water data measured at MARIN for the ducted propeller Ka4-70 with operating without and inside a modified duct 19A [8]. The section geometry of this duct, denoted as duct 19Am was obtained by replacing the round trailing edge of the duct 19A [9] by a sharp trailing edge. This sharp trailing edge allows the application of a classical Kutta condition for the prediction of the duct circulation.

Details of the mathematical formulation of the Panel Method are shown in Section 2. The numerical models for the gap region and the interaction with the duct boundary layer are presented in Section 3. This section also discusses the two wake alignment methods that have been investigated, which are a rigid wake model with prescribed geometry and an iterative wake alignment model. The influence of these models on the results and the comparison with the experimental data is shown in Section 4. In Section 5, the main conclusions are drawn.

2. Mathematical Formulation

2.1. Potential Flow Problem

Consider a propeller of radius rotating with constant angular velocity inside a duct both advancing with constant axial speed along its axis in an incompressible ideal fluid of constant density and kinematic viscosity otherwise at rest in a domain extending to infinity in all directions. The propeller is made of blades symmetrically distributed around an axisymmetric hub. The duct is also considered to be axisymmetric of inner radius at the propeller plane , which defines a gap height .

In a moving reference frame advancing and rotating with the propeller, we introduce a Cartesian coordinate system , with the positive -axis direction opposite to the propeller axial motion, the -axis coincident with the propeller reference line, passing through the reference point of the root section of the blade , and the -axis completing the right-hand system. We use a cylindrical coordinate system related to the Cartesian system by the transformation:

Figure 1 shows the coordinate system used to describe the propeller geometry and the fluid domain around the ducted propeller.

Figure 1: Propeller coordinate system.

In the reference frame moving with the propeller blades, the flow is steady, and assuming that the perturbation velocity due to the ducted propeller is irrotational then, the velocity field may be described by a perturbation potential in the form:

where

is the undisturbed onset velocity in the moving reference frame.

The perturbation potential satisfies the Laplace equation:

The boundary of the domain consists of the blade surfaces , the duct surface , and the hub surface . The kinematic boundary condition,

is satisfied on the blade, duct, and hub surfaces, where denotes differentiation along the normal and is the unit vector normal to the surface directed outward from the body. At infinity, the flow disturbance due to the ducted propeller vanishes

To allow for the existence of circulation around the propeller blades and the duct, vortex sheets are shed from the trailing edge of the blades and the duct. The boundary conditions on the vortex sheet surfaces are the tangency of the fluid velocity on each side of the sheet and the continuity of the pressure across the sheet. In steady flow, these conditions are

where is the pressure and the indices + and − denote the two sides of the vortex sheets, taken, respectively, on the side of the back (normally the suction side) and face (normally the pressure side) of the blade at the trailing edge, and on the inner side and outer side of the duct at the trailing edge. In the case of the propeller blade, the unit normal to the vortex sheet is defined pointing from the face to the back side of the blade. In the case of the duct, the unit normal to the vortex sheet is defined pointing from the outer to the inner side of the duct.

In order to specify uniquely the circulation around the blades and duct, it is necessary to impose the Kutta condition at the blade trailing edge and at the duct trailing edge. The Kutta condition states that the velocity must remain finite

at a sharp trailing edge. In the theoretical formulation, the Kutta condition, (8), is needed to guarantee that the vortex wakes are shed from the sharp trailing edge, and (7) should be applied at the vortex wake including the trailing edge.

The excess pressure due to the fluid motion (with respect to pure hydrostatic conditions) may be determined from the Bernoulli’s equation for the steady flow in the moving reference frame in the form

where is the pressure of the undisturbed inflow.

Applying Green’s second identity, assuming for the interior region to , , we obtain the integral representation of the perturbation potential at a point on the body surface,

where , is the distance between the field point and the point on the boundary . With the on the surfaces , , and known from the Neumann boundary condition on the body surface, (5), (10) is a Fredholm integral equation of the second kind in the dipole distribution on the surfaces , , and . The Kutta condition, (8), yields the additional relationship between the dipole strength in the wake surfaces and the surface dipole strength at the blade and duct trailing edges. In potential flow theory, the dipole strength at a given point on the wake equals the circulation for a closed contour around the lifting surface (blade or duct) through that point.

2.2. Velocity, Pressure and Forces

The velocity on the surface is obtained by differentiation of the surface potential distribution. From Bernoulli’s equation, (9), the pressure coefficient can be determined from:

where .

The inviscid thrust and torque on the ducted propeller are obtained by integration of the pressure distribution on the blade and duct surfaces. The propeller operation conditions are defined by a single nondimensional parameter: the advance coefficient , where is the rate of revolution and the propeller diameter. The nondimensional thrust and torque of the ducted propeller system are given by the propeller thrust coefficient , the duct thrust coefficient , and the torque coefficient :

where is the propeller thrust, the duct thrust, and the propeller torque. The ducted propeller efficiency is given by

3. Numerical Method

3.1. Surface Discretisation

For the numerical solution of the integral equation, (10), we discretise the blade surfaces , the duct surface , the hub surface , and the wake surfaces in bilinear quadrilateral panels. The blade is discretised in the spanwise radial direction by a number of strips, extending chordwise from the blade leading edge to the trailing edge. A cosine formula for the panel spacing in the radial and chordwise directions is used to increase the panel density near the blade extremities. The duct surface is divided into three regions: between, upstream, and downstream of the blades. The regions between and downstream of the blades are discretised along the axial direction with the same panel spacing of the blade tip section and blade wake tip section, respectively. The region upstream of the blades is discretised along the axial direction with a Vinokur stretching function [10]. Equidistant angular spacing is used in circumferential direction. For the discretisation of the hub surface, an elliptical grid generator is used [11]. The blade and duct wake surfaces are discretised in the spanwise and circumferential directions, extending downstream from the trailing edge the corresponding strips on the propeller blade and duct surfaces.

3.2. Solution of the Integral Equation

The integral equation, (10), is solved by the collocation method with the element centre point as collocation point. We assume a constant strength of the dipole and source distributions on each element. The influence coefficients are determined analytically using the formulations of Morino et al. [12].

3.3. Wake Models

Two wake models are considered: a rigid wake model and a wake alignment model for the blade wake where the pitch of the vortex lines are aligned with the local fluid velocity.

In the rigid wake model, the geometry of the wake surfaces is specified empirically. For the blade wake, the pitch of the vortex lines is assumed constant along the axial direction and equal to the blade pitch. However, contraction of the blade wake is prescribed by an empirical formulation, following Hoshino [13] for the open propeller case. The radial variation of the helicoidal lines is occurring in a transition wake region, with a length of , and is kept constant further downstream. For the duct, the wake leaves the trailing edge at the bisector. The value of the dipole strength at the blade and duct trailing edges is determined by the application of an iterative pressure Kutta condition, requiring that the pressure is equal on both sides of the blade and duct panels adjacent to the trailing edge.

In the wake alignment model, the corner points of the blade wake grid panels are displaced with the mean fluid velocity. At the iteration, the geometry in cylindrical coordinates of the wake strip can be determined by using an Euler scheme:

where , , and are the components of the mean vortex sheet velocity along the axial, radial, and circumferential directions, respectively, and is the time step for the Euler vortex convection scheme. The velocity components are calculated from the integral equation of the velocity, obtained by taking the gradient of (10). The velocities are evaluated at a staggered grid of control points, which are located equidistantly from consecutive dipole singularities of the wake. The first wake strip () corresponds to the blade trailing edge. To control the wake alignment stability, the radial coordinates of the blade wake grid are kept constant during the iterative process. Hence,

where

The new axial and circumferential coordinates and are calculated using an Euler scheme:

The nondimensional time step is introduced, which can also be expressed in terms of the number of time steps per propeller revolution . Note that in the wake alignment model, the downstream part of the duct is remeshed at each iteration of the propeller and duct vortex wakes alignment.

3.4. Gap Flow Models

Two different models for the potential flow in the gap region are considered: a closed gap with zero gap width and a gap flow model with transpiration velocity.

In the closed gap model, the blade and duct grids match in the tip region. In this case, the flow is not allowed to pass between the blade tip and the duct inner surface. The blade wake connects to the duct inner surface down to the duct trailing edge.

In the gap model with transpiration velocity, a partial flow between the blade tip and the duct inner surface is allowed to pass in the gap region. The procedure proposed by Hughes [14] is implemented in the current method, where the gap flow is treated as a two-dimensional orifice. The velocity in the gap region can be related to the difference in pressure across the blade tip using Bernoulli’s equation. The reduction in the flow from losses in the orifice is expressed in terms of an empirical discharge coefficient from the total volume flow rate through the gap and the pressure difference across the gap as

A typical value as used by Hughes [14] is . From (18), the mean velocity through the gap at a given chordwise location can be expressed in terms of the difference in pressure between the pressure and suction sides of the blade at that location:

where is the mean relative flow velocity through the gap and is the pressure difference at the blade tip in terms of the pressure coefficient defined by (11).

To incorporate this expression for the mean velocity into the Panel Method, an additional strip of panels along the blade tip is introduced which closes the gap. This wake strip continues past the duct trailing edge, connecting the duct wake with the blade wake. The dipole strength on the last wake strip is also determined from the Kutta condition. In the case of a closed gap, the boundary condition on the gap panels sets the source strength to cancel the normal component of the inflow velocity, (5). To allow for the existence of a transpiration velocity through the gap, the boundary condition on the gap strip becomes

where is the unit normal vector to the mean camber line at the gap strip on the same chordwise position of the panel.

The value of is obtained by an iterative procedure. First, the potential flow problem is solved as if the gap were completely closed and (5) is used to impose the kinematic boundary condition on the gap strip. The pressure distribution is then computed by differentiating the potential using a finite difference scheme. The potential flow problem is solved again with the kinematic boundary condition on the gap strip specified by (20). The pressures on the gap strip are then recomputed and used to update the boundary condition on the gap panels. The solution process will be repeated until a specified convergence criterion is met.

3.5. A Simple Model for the Interaction of the Blade Wake with the Duct Boundary Layer

A simple model for the interaction between the blade wake and the boundary layer on the duct inner side is implemented as well. Due to the duct boundary layer, a reduction in the axial velocity may be taken into account in the convection of blade vorticity in the wake alignment model. Considering as the duct boundary layer thickness and assuming a power law distribution for the velocity profile, we have

To correct the potential flow velocities in the blade wake, a linear variation of the axial velocity is considered in the gap region, . The corrected axial velocity at the duct surface, , is calculated by extrapolation from its value at the gap, , and at the edge of the duct boundary layer, . The velocity at the edge of the duct boundary layer is identical to the velocity on the duct surface in the potential flow method. A correction factor for the velocity at the duct surface can be derived from (21) and is given by

4. Results

4.1. General

Results are presented for the propeller Ka4-70 without and inside the duct 19Am. The propeller Ka4-70 is a four-bladed propeller of the Kaplan type with a large chord at the blade tip. A pitch-diameter ratio of is considered. The duct 19Am is based on duct 19A, modified such that the trailing edge is almost sharp. The particulars of the Ka-series and the duct 19A section geometry can be found in Kuiper [15]. The section geometry of the duct 19Am is given in Bosschers and van der Veeken [8], which also contains the results of the experiments performed at MARIN in the Netherlands. The gap between the duct inner surface and the blade tip is uniform and equal to of the propeller radius. The profiles of the duct 19A and 19Am are compared in Figure 2. The comparison of the experimental open-water characteristics of the ducted propeller system with duct 19A [15] and with duct 19Am [8] is presented in Figure 3.

Figure 2: Profiles of duct 19A and 19Am.

Figure 3: Experimental open-water characteristics of propeller Ka4-70 in duct 19A and duct 19Am.

Figure 4 shows a typical panel arrangement of propeller Ka4-70 inside the duct 19Am with an open-water hub. Calculations were carried out with the rigid wake model and wake alignment model. The wake geometries were obtained after 5 iterations of the wake alignment model using time steps per revolution, which leads to an angular step of 4 degrees. For the ducted propeller case, a transpiration velocity gap model with a discharge coefficient equal to [14] and a closed gap model were considered. A specified tolerance of on the gap strip control points was applied as convergence criterion, corresponding to approximately 20 iterations for the transpiration velocity. In general, the iterative pressure Kutta condition converged after 5 iterations to a precision of . In the iterative procedure to deal with nonlinearities of the numerical method, the outer loop is the wake alignment model, the middle loop is the transpiration velocity gap model, and the inner loop is the iterative pressure Kutta condition.

Figure 4: Panel arrangement for propeller Ka4-70 inside duct 19Am. 50 × 25 panels on each blade; 142 × 160 panels on the duct; 55 × 80 panels on the hub.

4.2. Ka4-70 without Duct

First, the numerical results are presented for the Ka4-70 propeller without duct in uniform flow. In this case, the rigid wake model and the same wake vortex pitch alignment model as used for the ducted propeller is applied in combination with, or without, an empirical contraction. Figure 5 shows the comparison of the thrust and torque coefficients with experimental data from the open-water tests. A section viscous drag coefficient of 0.007 and a lower limiting value for the chordwise component of the blade section force between sections and were used for all computations. The limiting value models the effect of flow separation which eliminates the unphysical suction peaks at the leading edge in the potential flow theory. The limiting value has been set to zero. For a flat plate, this would imply the suppression of the leading edge suction force. By applying a leading edge suction force correction, a significant increase of the torque coefficient at low advance ratio and a decrease in the thrust coefficient are obtained. Small differences between the wake models are seen with the exception of the low advance ratios. In this case, a reduction of the force coefficients is observed with wake alignment models for low advance ratios. These differences are related to the difficulty in obtaining converged and smooth aligned wake surfaces when large induced velocities are present. In general, a fair-to-good agreement is seen between the numerical predictions and the experimental data. The reasons for the underprediction of the force coefficients at low advance ratios are probably the limitations of the wake alignment model and the flow separation from the propeller blades, which is modelled by a simple suppression of the suction force.

Figure 5: Thrust and torque coefficients for propeller Ka4-70 showing the influence of the wake alignment model and the leading edge suction force (SF) correction.

Figure 6 presents the panel arrangement including the aligned wake for the advance ratio . To control the convergence and smoothness of the wake grid, the tip vortex roll-up phenomenon is not modelled with the present alignment method. Figure 7 shows the blade circulation obtained from the wake dipole strength and the chordwise pressure distribution at the radial section . Similar circulation and pressure distributions are obtained between all wake models. Please note that large suction peaks near the blade leading edge are obtained in the inviscid flow solution, which are physically not realistic.

Figure 6: Panel arrangement for propeller Ka4-70 at for the wake alignment model with empirical contraction. Note that the only one blade wake grid (400 × 25) is shown.

(a)

(b)

Figure 7: Radial circulation (a) and chordwise pressure (b) distributions. Propeller Ka4-70 at . Influence of the wake model.

4.3. Ka4-70 inside Duct 19Am

Secondly, the effect of the gap and the wake alignment models are investigated for propeller Ka4-70 inside duct 19Am. Calculations are presented for the rigid wake model and for the wake alignment model without and with correction of the axial velocity due to the duct boundary layer. The thickness of the duct boundary layer is roughly assumed to be , corresponding to . A power law velocity profile with is also assumed. The boundary layer thickness is a critical parameter and, since there is yet no information on suitable values to be used in this model, a sensitivity study is presented at the end of this section.

Figure 8 shows the comparison of the predicted thrust and torque coefficients with experimental data from the open-water tests. A section viscous drag coefficient of 0.007 was used for all computations for the propeller. Lacking any adequate method to estimate the duct viscous drag under the influence of a propeller, no viscous drag correction to the duct thrust has been applied in the present comparison. The differences between the closed and transpiration gap models are very small in this case. The comparison with the experiments shows that the numerical results obtained with the rigid wake model and wake alignment model without duct boundary layer correction over-predict the propeller thrust and torque coefficients. By applying the axial velocity correction at the tip in the wake alignment model, a significative reduction of the propeller thrust and torque coefficients is observed. A good agreement with the measurements may be achieved. The duct thrust coefficient agrees well with the measurements for low advance coefficients. The correction for the duct boundary layer results in an increase of the duct thrust. However, for high advance ratios, all models over-predict the duct thrust, which is probably due to the occurrence of flow separation on the outer side of the duct, and it is not modelled in the present method. The differences between the measured and the predicted ducted propeller efficiencies are around for the advance ratios lower than .

(a)

(b)

Figure 8: Comparison between predicted thrust and torque coefficients for the Ka4-70 propeller inside duct 19Am. Predictions are presented for the transpiration velocity gap model (a) and the closed gap model (b).

Figure 9 presents the panel arrangement including the aligned wake with the duct boundary layer correction. A strong reduction of the blade wake pitch is seen at the tip. The pitch angle of the blade wake at the axial positions and is illustrated in Figure 10. The correction in the axial velocity due to the duct boundary layer introduces a large reduction of the vortex pitch at the blade wake tip, which is responsible for the significant decrease of propeller thrust and torque and increase of duct thrust.

(a)

(b)

Figure 9: Panel arrangement for propeller Ka4-70 inside duct 19Am at . Gap model with transpiration velocity. Wake alignment model with duct boundary layer correction. The duct grid is not shown in the bottom figure.

(a)

(b)

Figure 10: Vortex pitch distributions at and . Influence of gap and wake models. Ka4-70 inside duct 19Am at .

Figure 11 presents the circulation distributions on the propeller blade and duct. The duct circulation is shown as function of the angular position between two consecutive blade wakes at the duct trailing edge. The pressure distributions are shown in Figures 12 and 13 for the gap model with transpiration velocity and closed gap model, respectively. The pressure distributions are shown as function of the chordwise position , where is defined along the cylindrical section chord for the propeller blade and as the axial distance along the meridional duct section, from the leading edge to the trailing edge. In the closed gap model, the blade circulation at the tip is equal to the duct circulation discontinuity. Note that a finite circulation at the blade tip is also obtained in the gap model with transpiration velocity. Similar results are obtained with the two gap models. Small differences are seen in the blade tip and duct circulations. In the closed gap model, the flow is not allowed to pass through the gap and increases the blade tip and duct circulations. For the pressure distributions and propeller forces, similar results are obtained between the two gap models. From the comparison between the wake models, a decrease in the blade circulation is obtained with the wake alignment model, which has the opposite effect on the duct circulation. Although the correction due to the duct boundary layer affects only the blade wake geometry near the tip, a very strong influence on the propeller and duct circulation distributions is obtained. For the blade pressure distribution, no suction peak is seen in the wake alignment model with duct boundary layer correction, indicating a different incidence angle to the blade section. The differences in the duct pressure distributions on the duct inner side in front of the propeller explain the differences in duct thrust.

(a)

(b)

Figure 11: Blade (a) and duct (b) circulation distributions. Influence of gap and wake models. Ka4-70 inside duct 19Am at .

(a)

(b)

Figure 12: Blade chordwise pressure distribution at (a). Duct chordwise pressure distribution at degrees (b). Influence of wake model. Ka4-70 inside duct 19Am at . Gap model with transpiration velocity.

(a)

(b)

Figure 13: Blade chordwise pressure distribution at (a). Duct chordwise pressure distribution at degrees (b). Influence of wake model. Ka4-70 inside duct 19Am at . Closed gap model.

A sensitivity study for the influence of the duct boundary layer thickness on the predicted thrust and torque coefficients has been performed using the transpiration gap model. The relative differences between the numerical results and the experimental open-water thrust and torque values are presented in Table 1 for various values of the duct boundary layer thickness. It is observed that the influence is very significant and that an increase in boundary layer thickness always leads to a reduction of propeller thrust and torque.

Table 1: Variation with the duct boundary layer thickness of the relative differences between the predicted and experimental (exp) thrust and torque coefficients. Ka4-70 inside duct 19Am at . Transpiration velocity gap model.

0.0%
19.8%
17.2%
11.9%

1.0%
8.4%
9.0%
3.2%

2.0%
0.4%
4.0%
−3.0%

3.0%
−1.5%
2.2%
−4.3%

4.0%
−3.3%
0.2%
−5.9%

5.0%
−6.1%
−0.8%
−8.1%

5. Conclusions

A low-order Panel Method has been used to predict the thrust and torque of propeller Ka4-70 without and inside a modified duct 19A in open-water conditions. Calculations were carried out with a rigid wake model and an iterative wake alignment model for the blade wake pitch. For the open propeller, a fair-to-good agreement is seen between the numerical predictions and the experiments for all wake models. For the ducted propeller, the loading predictions of the duct propeller system were found to be critically dependent on the blade wake pitch, especially at the tip. A model to incorporate the influence of the duct boundary layer on the pitch of tip vorticity was found necessary to correctly predict the open-water characteristics of the ducted propeller. This mechanism was not studied before using panel methods but seems plausible from the physical point of view. The quantitative aspects of this mechanism will be investigated in more detail by comparing the results of the present panel method with RANS simulations.

References

A. Sanchez-Caja, P. Rautaheimo, and T. Siikonen, “Simulation of incompressible viscous flow around a ducted propeller using a RANS equation solver,” in Proceedings of the 23rd Symposium on Naval Hydrodynamics, Val de Reuil, France, 2000.

M. Abdel-Maksoud and H.-J. Heinke, “Scale effects on ducted propellers,” in Proceedings of the 24th Symposium on Naval Hydrodynamics, Fukuoka, Japan, 2003.

V. I. Krasilnikov, J. Y. Sun, Z. Zhang, and F. Hong, “Mesh generation technique for the analysis of ducted propellers using a commercial RANSE solver and its application to scale effect study,” in Proceedings of the 10th Numerical Towing Tank Symposium (NuTTS’07), Hamburg, Germany, 2007.

J. E. Kerwin, S. A. Kinnas, J.-T. Lee, and W.-Z. Shih, “A surface panel method for the hydrodynamic analysis of ducted propellers,” Transactions of Society of Naval Architects and Marine Engineers, vol. 95, 1987.

H. Lee and S. A. Kinnas, “Prediction of cavitating performance of ducted propellers,” in Proceedings of the 6th International Symposium on Cavitation, Wageningen, The Netherlands, 2006.

J. Baltazar and J. A. C. Falcão de Campos, “On the modelling of the flow in ducted propellers with a panel method,” in Proceedings of the 1st International Symposium on Marine Propulsors, Trondheim, Norway, 2009.

D. Rijpkema and G. Vaz, “Viscous flow computations on propulsors: verification, validation and scale effects,” in Proceedings of the Developments in Marine CFD, London, UK, 2011.

J. Bosschers and R. van der Veeken, “Open water tests for propeller Ka4-70 and duct 19A with a sharp trailing edge,” MARIN Report 224457-2-VT, 2008.

J. D. vanManen and M. W. C. Oosterveld, “Analysis of ducted-propeller design,” Transactions of Society of Naval Architects and Marine Engineers, vol. 74, pp. 522–562, 1966.

M. Vinokur, “On one-dimensional stretching functions for finite-difference calculations,” Journal of Computational Physics, vol. 50, no. 2, pp. 215–234, 1983.

R.L. Sorensen, “Three-dimensional elliptic grid generation about fighter aircraft for zonal finite-difference computations,” in Proceedings of the AIAA 24th Aerospace Sciences Conference, Reno, Nev, USA, 1986.

L. Morino, L. T. Chen, and E. O. Suciu, “Steady and oscillatory subsonic and supersonic aerodynamics around complex configurations,” AIAA Journal, vol. 13, no. 3, pp. 368–374, 1975.

T. Hoshino, “Hydrodynamic analysis of propellers in steady flow using a surface panel method 2nd report: flow field around propeller,” Journal of The Society of Naval Architects of Japan, vol. 166, pp. 79–92, 1989.

M. J. Hughes, “Implementation of a special procedure for modeling the tip clearance flow in a panel method for ducted propellers,” in Proceedings of the Propellers/Shafting ’97 Symposium, Virginia Beach, Va, USA, 1997.

G. Kuiper, The Wageningen Propeller Series, no. 92-001, MARIN Publication, 1992.