Postal adress: Visit adress Telephone Fax Bank
Høgskolen i Ålesund Larsgårdsvegen 2 70 16 12 00 70 16 13 00 7694 05 00636
N-6025 Ålesund Internett E-mail Enterprise no.
Norway www.hials.no [email protected] NO 971 572 140
TITLE:
CFD Investigation in Scale Effects on Propellers with Different Blade Area Ratio
CANDIDATE NAME:
Zicheng Chen
DATE: COURSE CODE: COURSE TITLE: RESTRICTION: 29/5/2015 IP501909 MSc Thesis
STUDY PROGRAM: PAGES/APPENDIX: LIBRARY NO.:
MSc-Ship Design 63/3
SUPERVISOR(S):
Karl Henning Halse
ABSTRACT:
The performance prediction of marine propulsor has long been assessed by means of model test. Due to the larger difference in Reynolds number between model and full scale, the traditional extrapolation of measured force to full scale may be questionable, which implies knowledge investment into the scale effect on propellers. In the report, a number of propellers with blade area ratio variation were analyzed systematically by a RANSE method, the obtained numerical results agreed with the towing tank test fairly well. Parts Based Meshing technique is utilized for the highly automated mesh generation in model and full scale. The concerning scale effects due to blade area ratio variation are demonstrated and explained by comparing open water characteristics including thrust, torque, and pressure distribution, etc. The study indicates that trailing edge radial flow behavior contributes a lot to the propeller scale effect.
1
PREFACE
This report is the result of my Master’s thesis project for the two-year MSc program Ship Design at Ålesund University College (HIÅ).
The motivation for me to cover this topic is to fully develop the capability of dealing with marine engineering problems with the state-of-the-art CFD (Computational Fluid Dynamics) tools. In this report, a number of systematically varied propellers were analyzed numerically and validated against experimental results, with primary focus on the problem of scale effect.
This MSc thesis is originated from the research project PROPSCALE (Full Scale Performance Prediction for Energy Efficient Ship Design), 2013-2016. Led by MARINTEK, the project involves research partners including CSSRC, NTNU, HIÅ and TUHH, in collaboration with the leading Norwegian industrial partners like Rolls-Royce Marine AS, etc. Its primary objective is to acquire theoretical knowledge and develop practical methods for the accurate prediction of vessels equipped with novel energy efficient propulsion systems in full scale, and to quantify scale effect on propulsive characteristics of these vessels through the development and application of CFD methods.
I wish to express my sincere gratitude to my supervisor, Professor Karl Henning Halse, for providing suggestion and help throughout the project. I am also very grateful to Dr. Vladimir Krasilnikov, manager of the project PROPSCALE. Maine application of CFD method requires large amount of best practices, I have greatly benefited from his generous support and valuable comments during the research.
Ålesund, May 25, 2015.
Zicheng Chen
2
ABSTRACT
The performance prediction of marine propulsor has long been assessed by means of model test.
Due to the larger difference in Reynolds number between model and full scale, the traditional extrapolation of measured force to full scale may be questionable, which implies knowledge investment into the scale effect on propellers. In the report, a number of propellers with blade area ratio variation were analyzed systematically by a RANSE method, the obtained numerical results agreed with the towing tank test fairly well. Parts Based Meshing technique is utilized for the highly automated mesh generation in model and full scale. The concerning scale effects due to blade area ratio variation are demonstrated and explained by comparing open water characteristics including thrust, torque, and pressure distribution, etc. The study indicates that trailing edge radial flow behavior contributes a lot to the propeller scale effect.
Keywords
Propeller, Scale effect, Blade area ratio, RANSE
3
Tables of contents
TERMINOLOGY ... 4
1. INTRODUCTION ... 5
2. BACKGROUND AND THEORETICAL BASIS ... 6
2.1 Propeller theory ... 6
2.2 CFD introduction ... 10
2.3 STAR-CCM+ ... 18
3. METHODOLOGY ... 19
3.1 Propeller geometry ... 19
3.2 Simulation setup ... 23
4. RESULTS ... 36
4.1 Plot and Scene ... 36
4.2 Propeller open water characteristics ... 47
4.3 Experimental results VS numerical results ... 52
4.4 Pressure and friction components ... 53
5. DISCUSSION ... 55
5.1 Comparison between numerical and experimental results ... 55
5.2 General propeller open water characteristics ... 57
5.3 Scale effects ... 58
6. CONCLUSION ... 65
7. FUTURE WORK ... 66
REFERENCES ... 67
APPENDIX ... 68
4
TERMINOLOGY
Symbols
Γ Circulation
𝛽𝑖 Section hydrodynamic pitch angel [deg]
𝑇 Propeller thrust [N]
Q Propeller torque [N·m]
𝜀 Lift to Drag ratio Z Number of blades 𝜂0 Propeller efficiency D Propeller diameter [m]
n Propeller rotational speed [Hz]
𝑉𝐴 Advance speed [m/s]
𝜌 Water density [kg/m3]
µ Water dynamic viscosity [Pa·s]
𝜐 Water kinematic viscosity [m2∙ s−1] 𝐾𝑇 Propeller thrust coefficient
𝐾𝑇
Propeller torque coefficient J
Advance Coefficient Re Reynolds number div Vector operator 𝑈⃗⃗
Velocity vector
𝜙 General scalar conserved property 𝑓𝑜 Maximum section camber 𝑥𝑟
Propeller section rake 𝑏𝑟
Chord length 𝑒𝑜 Max thickness 𝑇𝑠 Propeller section skew 𝐶𝑝 Pressure Coefficient
𝜁𝑥, 𝜁𝑦, 𝜁𝑧 Vorticity component in three direction 𝛿 Boundary layer thickness [m]
Abbreviations
BAR Blade area ratio
DNS
Direct Numerical Simulation
LES
Large eddy simulation
RANS
Reynolds Averaged Navier-Stokes CPP Controllable pitch propeller
FPP Fixed pitch propeller
TKE Turbulent kinetic energy CFL Courant number
5
1. INTRODUCTION
The global shipping transportation has been increased significantly in the past decades. This raises significant concerns in the society regarding the environmental impact of sea operations. As a result, the shipbuilding industry is facing an ever more challenging regulatory leading to drastic reductions of emissions of pollutants. Therefore, it requires further development of new energy efficient designs of vessels and propulsion systems. This, in its turn, requires accurate prediction of vessels equipped with novel energy efficient propulsion systems in full scale.
Nowadays, the performance prediction of propulsion system is still mainly based on the results of towing tests, open water tests and self-propulsion tests in model scale. Due to the larger difference in Reynolds number between model and full scale, the traditional extrapolation procedure may be inadequate to capture the scale effect on characteristics of propellers, which can be extremely complicated by propeller geometry and the interaction among blades, duct and hull structure in some cases.
In recent years, significant progress has been made in the scale effect on marine propulsor by means of CFD methods. Krasilnikov et al. [1], has studied the influence of blade skew on scale effect using FLUENT solver with the SST k-ω model. Investigations on the scale effect of ducted propellers was carried out in Abdel-Maksoud et al. [2] and Krasilnikov et al. [3].
In the present study, numerical simulation with systematically varied propellers are performed in the commercial code STAR-CCM+, which solves RANS equations in their integral form, by means of Finite Volumes methods. Highly automatically mesh generation is also realized, due to the Parts Based Operation functionality in STAR-CCM+. The numerical simulation is performed as implicit unsteady, velocities and pressures are solved in a segregated flow, and then coupled by means of the SIMPLE algorithm. SST k-ω turbulence model is selected for its better performance compared to other two- equation models when boundary layer flow is under adverse pressure gradients (e.g. separation), which is expected under heavy loading condition.
Scale effect on propeller open water characteristics were studied in this project, nevertheless open water characteristics are not always sufficient for an adequate propeller design which should account for the interaction between propeller and ship hull, however this is not addressed in this paper. When applying CFD methods with assumption of fully turbulent flow to predict propeller performance in model scale, one should be aware that laminar flow zone can exist during model test, this will contribute to the difference between numerical and experimental results. It requires more advanced transition model to capture this effect. When applying the numerical results in full scale to practice use, one has to take consideration of roughness effects. Besides, other effects like propeller cavitation has a considerable impact on the propeller performance under certain conditions, this requires additional investigation in future.
6
2. BACKGROUND AND THEORETICAL BASIS
The propeller, whose name comes from the Latin “propellare” (to drive forward) is a very old idea.
It was however at the beginning of the 19th century, when a suitable power source-the steam engine-was developed that an efficient screw propeller was produced, AB Volvo Penta [4]. Most marine vessels are fitted with propulsor to overcome resistance as it moves forward. Meanwhile, operations like maneuvering and station-keeping also require this kind of robust driving source.
2.1 Propeller theory
2.1.1 Working Principles
A propeller creates pulling power by moving water and creating a column of water behind it. Hence, the accelerated water will generate a counterforce on the blade, which forms the so-called blade force. As shown in Fig.2.1, the propeller’s blade force which is caused by the pressure difference between the blade’s pressure and suction side, can be split into two parts. Axial component makes up the thrust, and the other tangential component constitutes the torque. The thrust and torque in their non-dimensional form are among the most important values describing the propeller’s performance.
Figure 2.1 Blade force diagram Figure 2.2 Water column in propeller slipstream
Behind the propeller, the increased water cylinder moves both axially (equivalent to the propeller’s thrust) and in a rotational direction (equivalent to the propeller’s torque). Hence, the water cylinder or column has a kinetic energy. See Fig.2.2. The kinetic energy in the outgoing column is an important part of the propeller’s energy diss. Other losses are friction when the blades cut through the water, flow separation due to vortex and cavitation bubbles at the blade.
There are several terminologies to define the propeller's characteristics such as: diameter, pitch, blade area ratio, rake, skew etc. All these characteristics are calculated and optimized to design the propeller accordingly to specific needs of the customer and the ship characteristics. The propeller’s blade area ratio BAR is a parameter used to relate the size of a blade to its diameter, see Figure 2.3. It is a term used to denote the ratio of expanded area 𝐴𝐸 to disc area 𝐴𝑜. If one knows the blade width at 70% of the propeller radius then the blade area ratio is expressed as:
7 𝐵𝐴𝑅 =𝐴𝐸
𝐴𝑜= 0.43 × 𝑁𝑜. 𝑏𝑙𝑎𝑑𝑒𝑠 ×𝐵𝑙𝑎𝑑𝑒 𝑤𝑖𝑑𝑡ℎ
𝐷𝑖𝑎𝑚𝑒𝑡𝑒𝑟 (2.1.1)
BAR is critical to the control of cavitation and changes of BAR affect propeller’s efficiency a lot. In the present paper, three types of propellers with systematically varied BAR will be investigated.
Figure 2.3 Blade area and disc area
2.1.2 Blade element theory
Blade element theory is a mathematical process originally designed to determine the behavior of propellers. It involves breaking a blade down into several small sections then determining the forces on each of these small blade elements.
The lift force generated on each blade section due to the pressure differential on the suction and pressure side of propeller, can be derived from the famous Kutta–Joukowski theorem:
𝑑𝐿 = 𝜌𝑉𝑅𝛤(𝑟)𝑑𝑟 (2.1.2)
where Γ is the circulation, 𝑉𝑅 is the flow velocity. At the same time, there will be a drag force with the direction normal to the lift force and coincident with incoming flow. The Drag component is caused by viscosity and will be affected greatly by Reynolds number. According to blade element theory, a typical propeller blade section at the radius r and the corresponding velocity diagram are shown blow.
Figure 2.4 Blade element theory
From the plot, one can easily find the lift force 𝑑𝐿 and drag force 𝑑𝐷 with their components in axial and tangential direction. The resulting thrust and torque can be expressed as:
𝑑𝑇 = 𝑑𝐿𝑎− 𝑑𝐷𝑎= 𝑑𝐿 𝑐𝑜𝑠 𝛽𝑖− 𝑑𝐷 𝑠𝑖𝑛 𝛽𝑖 (2.1.3) 𝑑𝐹 = 𝑑𝐿𝑡+ 𝑑𝐷𝑡= 𝑑𝐿 𝑠𝑖𝑛 𝛽𝑖+ 𝑑𝐷 𝑐𝑜𝑠 𝛽𝑖 (2.1.4)
8
where 𝑑𝐿 = 𝜌𝑉𝑅𝛤(𝑟)𝑑𝑟, 𝛽𝑖 is the section hydrodynamic pitch angel. The related torque for the blade section 𝑑𝑄 = 𝑟𝑑𝐹. Finally, the total thrust of the propeller will be the integration of axial lift vectors for the sections from root to tip.
𝑇 = 𝑍 ∙ ∫ 𝐹𝑅 𝐿∙ 𝑐𝑜𝑠 𝛽𝑖∙ (1 − 𝜀 𝑡𝑎𝑛 𝛽𝑖) ∙ 𝑑𝑟
𝑟𝐻 (2.1.5)
𝑄 = 𝑍 ∙ ∫ 𝐹𝑅 𝐿∙ 𝑠𝑖𝑛 𝛽𝑖∙ (1 + 𝜀 𝑐𝑜𝑡 𝛽𝑖) ∙ 𝑟 ∙ 𝑑𝑟
𝑟𝐻 (2.1.6)
where 𝜀 is the lift to drag ratio, 𝜀 = 𝐹𝐷/𝐹𝐿 and Z is the number of blades.
When propeller is operating at the advance speed of 𝑉𝐴 with rotational speed of n, the effective power is 𝑇 ∙ 𝑉𝐴, and the received power from engine is defined as 2πnQ. Therefore propeller efficiency is :
𝜂0= 𝑇𝑉𝐴
2πnQ (2.1.7)
2.1.3 General open water characteristics
As we may know, the thrust and torques produced by the propeller are expressed in their most fundamental form in terms of a series of non-dimensional characteristics. To establish these expressions, the principle of dimensional similarity can be applied to geometrically similar propellers. One can imagine, the hydrodynamic performance (thrust and torque) of a propeller is depend upon the following parameters:
(a) The propeller diameter (D).
(b) The rotational speed (n).
(c) The speed of advance (𝑉𝐴).
(d) The density of the fluid (𝜌).
(e) The dynamic viscosity of the fluid (µ).
(f) The gravity (g).
Therefore, the thrust (T) can be assumed to be proportional to D, n, 𝑉𝐴, 𝜌, µ and g.
𝑇 = 𝑘𝐷𝑎 𝑛𝑏 𝑉𝐴𝑐 𝜌𝑑 𝜇𝑒 𝑔𝑓 (2.1.8) where k is proportional coefficient, and a, b, c, d, e, f are unknown index.
Since the above equation must be dimensionally correct, and the parameters involved can be replaced by equating indices for M (mass), L (length) and T (time).
𝑀𝐿
𝑇2 = 𝑘𝐿𝑎 (1 𝑇)𝑏 (𝐿
𝑇)𝑐 (𝑀 𝐿3)𝑑 (𝐿2
𝑇)𝑒 (𝐿
𝑇2)𝑓 (2.1.9)
From the above
𝑀: 1 = 𝑑
𝐿: 1 = 𝑎 + 𝑐 − 3𝑑 + 2𝑒 + 𝑓 𝑇: − 2 = −𝑏 − 𝑐 − 𝑒 − 2𝑓} →
𝑑 = 1 𝑎 = 4 − 𝑐 − 2𝑒 − 𝑓
𝑏 = 2 − 𝑐 − 𝑒 − 2𝑓} (2.1.10) Further
𝑇 = 𝑘𝐷4−𝑐−2𝑒−𝑓 𝑛2−𝑐−𝑒−2𝑓 𝑉𝐴𝑐 𝜌1 𝜇𝑒 𝑔𝑓= 𝑘𝜌𝑛2𝐷4(𝑉𝐴
𝑛𝐷)𝑐( 𝜐
𝑛𝐷2)𝑒( 𝑔𝐷
𝑛2𝐷2)𝑓 (2.1.11) One more universal format of the above equation is
𝑇 = 𝜌𝑛2𝐷4∙ 𝑓1(𝑉𝐴 𝑛𝐷,𝑛𝐷2
𝜐 ,𝑛2𝐷2
𝑔𝐷 ) (2.1.12)
9
These non-dimensional groups are known by the following:
Thrust coefficient 𝐾𝑇
𝐾𝑇= 𝑇
𝜌𝑛2𝐷4= 𝑓1(𝑉𝐴
𝑛𝐷,𝑛𝐷2 υ ,𝑛2𝐷2
𝑔𝐷 ) (2.1.13)
where 𝑉𝐴
𝑛𝐷 is the advance coefficient J, 𝑛𝐷
2 𝜐 and 𝑛
2𝐷2
𝑔𝐷 are Reynolds number Re and Froude number Fr respectively.
From experience, when a marine propeller is working sufficiently far away from the free surface so as not to cause surface waves, the influence of Froude number on propeller thrust can be ignored. In conclusion, if the propeller is submerged deeply enough, the open water characteristics are depend only upon advance coefficient J and Reynolds number Re.
𝐾𝑇= 𝑓1(𝐽, 𝑅𝑒) (2.1.14)
2.1.4 Propeller Scale effects
To be able to accurately predict a ship’s the propulsion efficiency, a good understanding of scale effects is a prerequisite. Traditionally, prediction of open water characteristics of a propeller design has been assessed by model scale tests. The experimental methodology is based on the similarity theory, which postulates that two propellers with geometric and dynamic similarity will have the identical hydrodynamic characteristics. The premise of dynamic similarity requires the equivalence of advance coefficient J, and Reynolds number Re of the target two propellers.
Assuming 𝑉𝐴𝑠, 𝑛𝑠, 𝐷𝑠, 𝜐𝑠 an 𝑉𝐴𝑚, 𝑛𝑚, 𝐷𝑚, 𝜐𝑚 represents the advance speed, rotational speed, diameter and dynamic viscosity of the fluid of the full-scale and mode-scale propeller respectively. 𝜆 is the scale factor and 𝜆 = 𝐷𝑠/𝐷𝑚.
To keep J and Re to be identical for model and full scale propeller, two equations must be fulfilled at the same time.
𝑉𝐴𝑚
𝑛𝑚𝐷𝑚= 𝑉𝐴𝑠 𝑛𝑠𝐷𝑠 𝑛𝑚𝐷𝑚2
𝜐𝑚 =𝑛𝑠𝐷𝑠2 𝜐𝑠 }
(2.1.15)
Obviously, the requirements above can’t be realized in reality due to the unrealistic rotational speed and advance speed for the model propeller. Hence, in model test only the advance coefficient J is kept identical, whilst 𝑅𝑒>𝑅𝑒𝑐.
The different propeller performance characteristics between model and full scale resulted from different Reynolds number, is so-called scale effects.
2.1.5 Scaling method
The scale effects affecting performance characteristics are essentially viscous in nature, and as such are primarily due to boundary layer phenomena dependent on Reynolds number. Due to the methods of testing model propellers and the consequent changes in Reynolds number between model and full scale, there can arise a different boundary layer structure to the flow over the blades. Whilst it is generally recognized that most full-scale propellers will have a primarily turbulent flow over the blade
10
surface this need not be the case for the model where characteristics related to laminar flow can prevail over significant parts of the blade. In order to quantify the effect of scale on the performance characteristics of a propeller, the 1978 ITTC recommended scaling procedure is presented.
𝐾𝑇𝑠= 𝐾𝑇𝑚− ∆𝐾𝑇
𝐾𝑄𝑠 = 𝐾𝑄𝑚− ∆𝐾𝑄
} (2.1.16)
where the scale corrections ∆𝐾𝑇 and ∆𝐾𝑄 are given by
∆𝐾𝑇= −0.3∆𝐶𝐷(𝑃 𝐷)(𝑐𝑍
𝐷) (2.1.17)
∆𝐾𝑄 = 0.25∆𝐶𝐷(𝑐𝑍
𝐷) (2.1.18)
The term ∆𝐶𝐷 relates to the change in drag coefficient introduced by the differing flow regimes at model and full scale, and is formally written as
∆𝐶𝐷= 𝐶𝐷𝑚− 𝐶𝐷𝑠 (2.1.19)
where
𝐶𝐷𝑚= 2 (1 +2𝑡
𝐶) [0.044 𝑅𝑛𝑥1/6− 5
𝑅𝑛𝑥2/3] (2.1.20)
𝐶𝐷𝑠= 2 (1 +2𝑡
𝐶) [1.89 + 1.62 log10(𝑐 𝐾𝑝)]
−2.5
(2.1.21) In these relationships t/c is the section thickness to chord ratio; P/D is the pitch ratio; c is the section chord length and 𝑅𝑛𝑥 is the local Reynolds number, all relating to the section located 0.75R. The blade roughness 𝐾𝑝 is taken as 30×10−6m.
It is well acknowledged that the ITTC-procedure method can provide a reliable and practical solution to propeller testing. However, the ITTC-procedure is unable to consider the local flow condition. In the meantime, it do not reflect correctly the effect of Reynolds number on propeller characteristics, in particular, as far as high skew and balanced skew propeller designs are concerned.
2.2 CFD introduction
Computational fluid dynamics, or CFD, is a computational technology that enables scientists and engineers to perform numerical experiment in a virtual laboratory. It employs numerical methods and algorithms to solve the equations that describe fluid flow. CFD method has being widely applied in marine industry for the prediction of ship performance in waves, design and analysis of propulsion system, investigation on behaviors of floating structure, as well as for many other applications. For the study of scale effects on marine propellers, CFD method can provide an insight into detailed flow patterns that are difficult, expensive or impossible to study using traditional experimental techniques.
Before the investigation, it is necessary to introduce the underlying physical problems and mathematical principles employed by the CFD method.
11
2.2.1 Governing equation
The equations that describe the dynamics of fluid represent fundamental laws of physics stating conservation of mass, momentum end energy. It should be noted that energy equation is normally not included in the marine CFD as water temperature is assumed to be fixed.
Mass conservation
For any studied elemental fluid volume, the mass conservation law states that the mass of fluid is conserved, i.e. rate of increase of mass in fluid element equals to the net rate of mass flow into the fluid element. Based on this law, the continuity equation can be derived as:
𝜕𝜌
𝜕𝑡+𝜕(𝜌𝑢)
𝜕𝑥 +𝜕(𝜌𝑣)
𝜕𝑦 +𝜕(𝜌𝑤)
𝜕𝑧 = 0 (2.2.1)
or in the mathematical definition for 𝑑𝑖𝑣 of a vector property:
𝜕𝜌
𝜕𝑡+ 𝑑𝑖𝑣(𝜌𝑈⃗⃗ ) = 0 (2.2.2)
The incompressible isothermal flow is the model we adopt in the marine CFD simulation for the majority of applications. Under this definition, the continuity equations becomes:
𝜕𝑢
𝜕𝑥+𝜕𝑣
𝜕𝑦+𝜕𝑤
𝜕𝑧 = 𝑑𝑖𝑣𝑈⃗⃗ = 0 (2.2.3)
Momentum conservation
According to the Newton’s second law, the net force applied on the fluid element equals its mass time the acceleration of the element. For a moving fluid element, the momentum equation in x, y and z- direction can be written as:
𝜌𝐷𝑢
𝐷𝑡 = −𝜕𝑝
𝜕𝑥+𝜕𝜏𝑥𝑥
𝜕𝑥 +𝜕𝜏𝑦𝑥
𝜕𝑦 +𝜕𝜏𝑧𝑥
𝜕𝑧 + 𝜌𝑓𝑥
𝜌𝐷𝑣
𝐷𝑡 = −𝜕𝑝
𝜕𝑦+𝜕𝜏𝑥𝑦
𝜕𝑥 +𝜕𝜏𝑦𝑦
𝜕𝑦 +𝜕𝜏𝑧𝑦
𝜕𝑧 + 𝜌𝑓𝑦
𝜌𝐷𝑤
𝐷𝑡 = −𝜕𝑝
𝜕𝑧+𝜕𝜏𝑥𝑧
𝜕𝑥 +𝜕𝜏𝑦𝑧
𝜕𝑦 +𝜕𝜏𝑧𝑧
𝜕𝑧 + 𝜌𝑓𝑧
(2.2.4)
For the isotropic Newtonian fluids, the viscous stress in the momentum equation can be related to the rates of linear deformations of the fluid element, and the latter are expressed through the velocity components. The viscous stress are given as:
𝜏𝑥𝑥 = 2𝜇𝜕𝑢𝜕𝑥; 𝜏𝑦𝑦= 2𝜇𝜕𝑣𝜕𝑦; 𝜏𝑧𝑧= 2𝜇𝜕𝑤𝜕𝑧; 𝜏𝑥𝑦= 𝜏𝑦𝑥= 𝜇 (𝜕𝑣𝜕𝑥+𝜕𝑢𝜕𝑦);
𝜏𝑥𝑧= 𝜏𝑧𝑥= 𝜇 (𝜕𝑢𝜕𝑧+𝜕𝑤𝜕𝑥); 𝜏𝑦𝑧= 𝜏𝑧𝑦 = 𝜇 (𝜕𝑤𝜕𝑦+𝜕𝑣𝜕𝑦)
(2.2.5)
Then substitute the above equations for viscous stress into momentum equation. For the x-, y- and z- momentum, one can write down:
𝜌𝐷𝑢
𝐷𝑡 = −𝜕𝑝
𝜕𝑥+ 𝜕
𝜕𝑥(2𝜇𝜕𝑢
𝜕𝑥) + 𝜕
𝜕𝑦[𝜇 (𝜕𝑣
𝜕𝑥+𝜕𝑢
𝜕𝑦)] + 𝜕
𝜕𝑧[𝜇 (𝜕𝑢
𝜕𝑧+𝜕𝑤
𝜕𝑥)] + 𝑆𝑀𝑥 𝜌𝐷𝑣
𝐷𝑡 = −𝜕𝑝
𝜕𝑦+ 𝜕
𝜕𝑥[𝜇 (𝜕𝑢
𝜕𝑦+𝜕𝑣
𝜕𝑥)] + 𝜕
𝜕𝑦(2𝜇𝜕𝑣
𝜕𝑦) + 𝜕
𝜕𝑧[𝜇 (𝜕𝑢
𝜕𝑧+𝜕𝑤
𝜕𝑥)] + 𝑆𝑀𝑦 𝜌𝐷𝑤
𝐷𝑡 = −𝜕𝑝
𝜕𝑧+ 𝜕
𝜕𝑥[𝜇 (𝜕𝑢
𝜕𝑧+𝜕𝑤
𝜕𝑥)] + 𝜕
𝜕𝑦[𝜇 (𝜕𝑣
𝜕𝑧+𝜕𝑤
𝜕𝑦)] + 𝜕
𝜕𝑧(2𝜇𝜕𝑤
𝜕𝑧) + 𝑆𝑀𝑧
(2.2.6)
The above three equations are the x, y and Z components, respectively, of the momentum equation.
12
They are scalar equations and are called Navier-Stokes Equations in honor of the two men-the Frenchman M.Navier and the Englishman G.Stokes-who independently obtained the equations in the first half of the nineteenth century, J.D Anderson [5].
2.2.2 Governing equations in conservative and integral forms
Using the unfolded expressions for the substantive derivative, mathematical definition of div and grad, the governing transport equations for the mass and momentum can be written as follows:
𝑑𝑖𝑣𝑈⃗⃗ = 0 𝜌𝜕𝑢
𝜕𝑡 + 𝜌𝑑𝑖𝑣(𝑢𝑈⃗⃗ ) = −𝜕𝑝
𝜕𝑥+ 𝜇𝑑𝑖𝑣(𝑔𝑟𝑎𝑑𝑢) + 𝑆𝑀𝑥
𝜌𝜕𝑣
𝜕𝑡 + 𝜌𝑑𝑖𝑣(𝑣𝑈⃗⃗ ) = −𝜕𝑝
𝜕𝑦+ 𝜇𝑑𝑖𝑣(𝑔𝑟𝑎𝑑𝑣) + 𝑆𝑀𝑦
𝜌𝜕𝑤
𝜕𝑡 + 𝜌𝑑𝑖𝑣(𝑤𝑈⃗⃗ ) = −𝜕𝑝
𝜕𝑧+ 𝜇𝑑𝑖𝑣(𝑔𝑟𝑎𝑑𝑤) + 𝑆𝑀𝑧
(2.2.7)
The equations are in the so-called conservative or divergence form. In view of obvious commonalities of these equations, one can write a general conserved form for fluid property 𝜙:
𝜌𝜕𝜙
𝜕𝑡 + 𝜌𝑑𝑖𝑣(𝜙𝑈⃗⃗ ) = 𝑑𝑖𝑣(Γ𝜙𝑔𝑟𝑎𝑑𝜙) + 𝑆𝜙 (2.2.8) The conservation forms the basis of computational procedure in the finite volume method. It clearly reflects different contributions into the transport of a fluid property. In the left-hand side, the first term expresses the rate id change of 𝜙 the fluid element, and the second term expresses the net of flow 𝜙 due to convection. In the right-hand side, the first term expresses the rate of change of 𝜙 due to diffusion with the corresponding diffusion coefficient Γ𝜙, and the second term expresses the rate of increase of 𝜙 due to sources. The key procedure in development of the finite volume method is the integration of Eq.2.2.8 over a 3D control volume CV, this can be written as follows:
𝜌 ∫𝜕𝜙
𝜕𝑡𝑑𝑉 + 𝜌
𝐶𝑉
∫ 𝑑𝑖𝑣(𝜙𝑈⃗⃗ )𝑑𝑉 = ∫ 𝑑𝑖𝑣(Γ𝜙𝑔𝑟𝑎𝑑𝜙)𝑑𝑉
𝐶𝑉
+ ∫ 𝑆𝜙𝑑𝑉
𝐶𝑉 𝐶𝑉
(2.2.9) The volume integrals of the convective and diffusive terms can be transformed into the surface integrals over the surface A bounding the control volume using the divergence theorem by Gauss-Ostrogradsky which states
∫ 𝑑𝑖𝑣𝑎 𝑑𝑉 = ∫ 𝑛⃗ 𝑎
𝐴
𝑑𝐴
𝐶𝑉
(2.2.10) Where 𝑎 is an arbitrary vector property and 𝑛⃗ is the normal to the surface element 𝑑𝐴. In worlds, the divergence theorem states that the outward flux of a vector field through a closed surface is equal to the volume integral of the divergence of this vector field on the volume boundary by the surface. After the application of the divergence theorem and changing the order of integration and differentiation in the rate of the change term, the integral form of the transport equations can be derived:
𝜌 𝜕
𝜕𝑡( ∫ 𝜙𝑑𝑉
𝐶𝑉
) + 𝜌 ∫ 𝑛⃗ ∙ (𝜙𝑈⃗⃗ )𝑑𝐴 = ∫ 𝑛⃗ ∙ (Γ𝜙𝑔𝑟𝑎𝑑𝜙)𝑑𝐴
𝐴
+ ∫ 𝑆𝜙𝑑𝑉
𝐶𝑉 𝐴
(2.2.11)
The integral form of the transport equations represents the statement of conservation of a fluid property for a finite size control volume. This is the principle difference of the integral form from the conservative from which expresses the same conservation principle for an infinitely small fluid element.
13
2.2.3 Discretization and Solution
Discretization method
To solve the transport equations in the computer, one has to transfer them into discretized form. This process is so-called discretization. The typical discretization methods are finite difference, finite element and finite volume methods. Due to limited space, only finite volume method is discussed here, because of its widely application in the CFD commercial solver.
In finite volume method, discretization schemes are used for the approximation of surface and volume integrals that represent different terms of the equations governing transport of solution variables. Such an approximation allows one to convert a general scalar transport equation to an algebraic equation that can be solved numerically. For a general scalar conserved property 𝜙, the integral form of the transport equation has been written in Eq.2.2.9. After the computational mesh is built and each control volume (CV) is defined, this governing equation can be integrated about each CV, resulting in discrete equations that conserve property 𝜙 on a control-volume base. The integral transport equation applies to each CV, as well as to the whole solution domain. Therefore, if one sums the discrete equations for all CVs, one obtains the global transport equations, because the surface integrals over inner CV faces cancel out.
Equation Solution method
Pressure-based segregated algorithm
The pressure-based solver uses a solution algorithm where the governing equations are solved sequentially (i.e., segregated from one another). Because the governing equations are non-linear and coupled, the solution loop must be carried out iteratively in order to obtain a converged numerical solution. In the segregated algorithm, the individual governing equations for the solution variables (e.g.
u, v, w, p, k, e etc.) are solved one after another. Each governing equation, while being solved, is
"decoupled" or "segregated" from other equations, hence its name. The segregated algorithm is memory-efficient, since the discretized equations need only be stored in the memory one at a time.
However, the solution convergence is relatively slow, inasmuch as the equations are solved in a decoupled manner.
Pressure-based coupled algorithm
In a coupled algorithm, the momentum and continuity equations are solved in a closely coupled manner.
Hence, the rate of solution convergence significantly improves when compared to the segregated algorithm. However, the memory requirement increases by 1.5 - 2 times that of the segregated algorithm since the discrete system of all momentum and pressure-based continuity equations needs to be stored in the memory when solving for the velocity and pressure fields.
2.2.4 Computation mesh and boundary
The discretization scheme of finite volume method implies that the control volume (CV) is the prerequisite for the numerical computation. These control volume or subdomains of the computation domain are often called elements or cells, and the collection of all elements or cells is called a mesh or grid. Modern mesh generation technique have evolved to the point where very complex computation domains can be meshed efficiently with a variety of mesh types to ensure high quality simulation. Based upon the connectivity of the mesh, the mesh can be classified into three types.
14 Structured Meshes
A structured mesh is characterized by regular connectivity that can be expressed as a two or three dimensional array. This restricts the element choices to quadrilaterals in 2D or hexahedra in 3D. It is robust in calculation, however obviously not effective for very complicated geometry.
Hexahedron Tetrahedron Pyramid Prism Polyhedron Figure 2.5 Cell types used in modern CFD codes
Unstructured Meshes
An unstructured mesh is characterized by irregular connectivity is not readily expressed as a two or three dimensional array in computer memory. Compared to structured meshes, the storage requirements for an unstructured mesh can be substantially larger since the neighborhood connectivity must be explicitly stored.
Hybrid Meshes
A hybrid mesh is a mesh that contains structured portions and unstructured portions. Those parts of the geometry that are regular can have structured grids and those that are complex can have unstructured grids. It combines both of merits of structured mesh and unstructured mesh.
Boundary condition
The equations of fluid motion are solved for a certain computation domain. In view of the nature of the problems solved and equations that describe them, as discussed in Section 2.2.1 and 2.2.2, this domain has to be restricted by boundaries. On the boundaries one has to specify the boundary condition and initial conditions under which the equations are solved. These conditions will tell the solver what physical processes take place and should be accounted for at the boundaries of computation domain. The typical boundary conditions in CFD are Wall condition, Symmetrical boundary condition, Inlet, outlet boundary condition and Periodic boundary condition.
Wall represent the impenetrable surfaces that bound the fluid. The most common type is No-slip boundary condition in viscous flow simulation which ensure that fluid sticks to the wall and moves with the same velocity as the wall.
The inlet and outlet boundary are the surface through which fluid enters and leaves the computation domain. The most common inlet condition include velocity inlet, pressure inlet and mass flow inlet.
The condition imposed on the outlet boundaries can be of pressure outlet, pressure far-field and outflow.
Symmetry boundaries allow one to take benefit of physical flow symmetry, in order to reduce the size of computation domain thus save memory and time.
Periodic boundaries allow for the account of periodically repeating nature of the flow in the simulation. This boundary condition is often introduced in the case of propeller.
15
Figure 2.6 Typical boundaries in the CFD simulation of vessel
2.2.5 Turbulence modelling approach
Turbulent flow is the most likely flow condition that can be met in the simulation of marine tasks. In CFD, the random nature of turbulent flow complicates their numerical simulation greatly. Extensive theoretical and experimental research in the mechanisms of turbulence has been made in recent decades. Due to the complexity of these work and related results, only key features of turbulent flow and its numerical models are included in this part.
Phenomenon of turbulence
The phenomenon of turbulence was discovered by Osborne Reynolds. In 1883 he performed a very illustrative experiment with the water flow in a pipe where he injected a painted jet in the middle of the pipe. The observed views show that there is a laminar-turbulent transition between laminar and turbulent flow, and this is affected by the Reynolds number. The majority of industrial flows are high-Reynolds- number flows and are, consequently, turbulent. Turbulent flow has the following features, according to Krasilnikov [6].
Turbulent flow are unsteady, their properties are random function of time.
Turbulent flows are three dimensional, even if the averaged velocity field has a dominant 2D direction, the instantaneous field fluctuates in 3D space.
Turbulent flows are chaotic vertical flows. In which, vortices (eddies) move in a chaotic fashion.
Turbulent flows fluctuate on a broad range of length and time scales.
The large eddies in the turbulent flow carry the main portion of turbulent kinetic energy.
Approaches to turbulence
Direct Numerical Simulation (DNS)
DNS methods solve the Navier-Stokes equations directly for all scales of turbulent motions without any turbulence modeling. At the initial step of the solution, the solver generates small disturbances which begin to grow in strength and amplitude as solution proceeds. So, in fact, in the numerical solution, just like in the experiment, one can obtain transition from laminar to turbulent flow regime. Thus, the chaotic motion and interaction of turbulent eddies of different scales is simulated by mathematical principles.
However turbulent eddies can vary significantly in scale, and a great part of them, especially in the
16
beginning are of microscopic length scale. In order to capture such vortices in the numerical solution, the grid resolution must be so fine that only supercomputers can handle it. Thus, it can hardly be expected that DNS will become an engineering tool.
Large eddy simulation (LES)
LES methods provide solution for the largest scale motions of turbulent flow. In other words, the largest and energetically most important turbulence eddies are computed directly, while the effect of smaller eddies, which are not resolved is accounted for through additional stresses obtained from the turbulence theory. The expressions for these stress are much simpler compared to turbulence models in RANS, and for them it is easier to build a consistent theory, because small-scale turbulence is isotropic. The large eddies carry the main portion of turbulent kinetic energy, and in LES they are computed directly, which is a big advantage. Nowadays, LES methods did not find wide use in engineering application, basically the same reasons as DNS. While only the large-scale eddies are calculated, LES methods still require great computational resources.
Figure 2.7 Illustration of turbulence scales resolved by DNS, LES and RANS method
Reynolds Averaged Navier-Stokes (RANS)
RANS methods are the most commonly used at present for practical calculations of viscous turbulent flows in marine CFD. The RANS approach is based on time averaging of general transport equations and representation of total flow characteristics (velocity and pressure) as a sum of averaged and fluctuating values. The turbulent stresses are modeled by turbulence model, either empirical or semi- empirical using experimental and statistical data. The RANS method resolves only the vortices of largest scale comparable with the size of flow domain, while the rest of turbulence is accounted for through a turbulence model. The RANS method provides the main tool for the engineering simulation of turbulent flows. In Reynolds averaging, the solution variables in the instantaneous Navier-Stokes equations are decomposed into the mean (ensemble-averaged or time-averaged) and fluctuating components. For the velocity components:
𝑢𝑖= 𝑢̅𝑖+ 𝑢𝑖′ (2.2.12)
where 𝑢̅𝑖 and 𝑢𝑖′ are the mean and fluctuating velocity components(𝑖 = 1,2,3) Likewise, for pressure and other scalar quantities:
𝜙 = 𝜙̅ + 𝜙′ (2.2.13)
where 𝜙 denotes a scalar such as pressure, energy, or species concentration. Substituting expressions of this form for the flow variables into the instantaneous continuity and momentum equations and taking a time average yields the ensemble-averaged momentum equations. They can be written in Cartesian tensor form as:
17
𝜕𝜌
𝜕𝑡+𝜕(𝜌𝑢𝑖)
𝜕𝑥𝑖 = 0 (2.2.14)
𝜕(𝜌𝑢𝑖)
𝜕𝑡 +𝜕(𝜌𝑢𝑖𝑢𝑗)
𝜕𝑥𝑗 = −𝜕𝜌
𝜕𝑥𝑖+ 𝜕
𝜕𝑥𝑗[𝜇 (𝜕𝑢𝑖
𝜕𝑥𝑗+𝜕𝑢𝑗
𝜕𝑥𝑖−2 3𝛿𝑖𝑗𝜕𝑢𝑙
𝜕𝑥𝑙)] + 𝜕
𝜕𝑥𝑗(−𝜌𝑢̅̅̅̅̅̅̅) 𝑖̇′𝑢𝑗̇′ (2.2.15) These equations are called Reynolds-averaged Navier-Stokes (RANS) equations. They have the same general form as the instantaneous Navier-Stokes equations, with the velocities and other solution variables now representing ensemble-averaged (or time-averaged) values. Additional terms now appear that represent the effects of turbulence. These turbulent stresses, −𝜌𝑢̅̅̅̅̅̅̅𝑖̇′𝑢𝑗̇′ must be modeled in order to close equation.
Turbulence model for RANS method
The Reynolds-averaged approach to turbulence modeling requires that the turbulent stresses be appropriately modeled. How to relate the turbulent stress to the dynamic of the flow is one of the central problems of the RANS method and it constitutes the subject of turbulence model development.
A common method employed by many turbulence models is the Boussinesq hypothesis, which implies the analogy between the process of turbulent mixing and molecular diffusion, which is the cause of viscous stresses. According to the hypothesis, a parameter named turbulent or eddy viscosity 𝜇𝑇 is introduced by direct analogy with conventional dynamic viscosity 𝜇, and the components of turbulent stress are expressed as a product of 𝜇𝑇 and rates of linear deformations.
−𝜌𝑢̅̅̅̅̅̅̅ = 𝜇𝑖̇′𝑢𝑗̇′ 𝑇(𝜕𝑢𝑖
𝜕𝑥𝑗+𝜕𝑢𝑗
𝜕𝑥𝑖) −2
3(𝜌𝑘 + 𝜇𝑡𝜕𝑢𝑘
𝜕𝑥𝑘) 𝛿𝑖𝑗 (2.2.16)
The Boussinesq hypothesis is used in the Spalart-Allmaras model, the k-𝜀 turbulence models, and the k-ω models. In the case of the Spalart-Allmaras model, only one additional transport equation (representing turbulent viscosity) is solved. In the case of k-𝜀 and k-ω models (two-equation model), two additional transport equations are solved for turbulent kinetic energy k[𝑚2/𝑠2] and its dissipation rate 𝜀 and 𝜔 [𝑚2/𝑠2]. The first variable, k, determines the energy in the turbulence, whereas the second variable can be thought of as the variable that determines the scale of the turbulence (length-scale or time-scale). Advantage of this approach is the relatively low computational cost associated with the computation of the turbulent viscosity. Disadvantage of this hypothesis is that it assumes an isotropic scalar quantity, which is not strictly true. When it comes to complex flows, like flows with strong curvature, or strongly accelerated or decelerated flows the Boussinesq assumption is simply not valid. This give two equation models inherent problems to predict the flow separation and even attached flow with adverse pressure gradients.
SST k-ω turbulence model
The SST k-ω model was introduced in 1994 by F.R. Menter to deal with the strong free stream sensitivity of the k-ω turbulence model and improve the predictions of adverse pressure gradients, Reference [7].
The model combines both of the merits of k-ω and k-𝜀 turbulence model. The use of a k-ω formulation in the inner parts of the boundary layer makes the model directly usable all the way down to the wall through the viscous sub-layer, hence the SST k-ω model can be used as a low-Re turbulence model without any extra damping functions. The SST formulation switches to a k-ε behavior in the free- stream and thereby avoids the common k-ω problem that the model is too sensitive to the inlet free- stream turbulence properties.
18
2.2.6 Near-Wall Treatment
Turbulent flows are significantly affected by the presence of walls. Numerous experiments have shown that the viscous-affected region can be largely made up of three layers with their corresponding wall 𝑦+, namely the:
Viscous sublayer (𝑦+< 5)
Buffer layer or blending region (5 ≤ 𝑦+≤ 30)
Fully turbulent or log-law region (𝑦+> 30 to 60)
The wall 𝑦+ is the distance to the way, made dimensionless with the friction velocity 𝑈𝑡 and kinematic viscosity 𝜐, 𝑦+ is similar to local Reynolds number.
𝑦+=𝑦 ∙ 𝑈𝑡
𝜐 (2.2.17)
Very close to the wall, viscous damping reduces the tangential velocity fluctuations, while kinematic blocking reduces the normal fluctuations. Toward the outer part of the near-wall region, however, the turbulence is rapidly augmented by the production of turbulence kinetic energy due to the large gradients in mean velocity. Finally, the interim buffer region between the viscous sublayer and the fully turbulent layer where the effects of molecular viscosity and turbulence are equally important. Detailed flow measurements and DNS calculation indicate that the values of 𝑦+ parameter of the buffer layer lie in the range 5 ≤ 𝑦+≤ 30.
Traditionally, there are two approaches to modeling the near-wall region. In one approach, the viscosity- affected inner region (viscous sublayer and buffer layer) is not resolved. Instead, semi-empirical formulas called "wall functions'' are used to bridge this region between the wall and the fully-turbulent region. In another "near-wall modeling'' approach, the turbulence models are modified to enable the viscosity-affected region to be resolved with a mesh all the way to the wall, including the viscous sublayer.
In most high-Re-number flows, the wall function approach substantially saves computational resources.
In the viscosity-affected near-wall region, in which the solution variables change most rapidly, does not need to be resolved. The wall function approach is popular because it is economical, robust, and reasonably accurate. It’s a practical option for the near-wall treatments for industrial flow simulations.
2.3 STAR-CCM+
STAR-CCM+ is one of the largest commercial CFD solver in the world, with its unrivalled ability to tackle problems involving multi-physics and complex geometries, it provides an entire engineering process for solving problems involving flow (of fluids or solids), heat transfer and stress. The solver come packaged with pre- and post-processors which provide everything a user needs to go from the raw CAD geometry to a final flow analysis and visualization. The result is a code that offers outstanding ease-of-use delivered through an object-based tree-structured GUI, which guides even the most novice user through the set-up and analysis of a CFD problem. Utilizing the latest numerical algorithms, physical models and state-of-the-art software coding, it provides the user with a toolset capable of tackling the most complex multi-disciplinary engineering problems without compromising ease-of-use over capability or accuracy.
To ensure that users are constantly updated with the very latest advances throughout the product, there are three major releases of STAR-CCM+ every year, Reference [8].
19
3. METHODOLOGY
Commercial code STAR-CCM+ has been used for the investigation in the scale effects on propellers with systematically varied blade area ratio. Highly automated mesh generation for model and full scale is realized, due to Parts Based Operation functionality in STAR-CCM+. The CFD code solves RANS equations in their integral form, by means of Finite Volumes methods. The spatial discretization of the convective terms is done with a second order upwind based scheme, whereas the diffusive terms are discretized with second order central scheme. The numerical simulation is performed as implicit unsteady, velocities and pressures are solved in a segregated flow, and then coupled by means of the SIMPLE algorithm. The SST k-𝜔 turbulence model is chosen for turbulence closure. The rotation of the propeller is modelled used a moving reference frame system. The numerical results will be compared to the experimental data.
3.1 Propeller geometry
To be able to carry out the numerical simulations in STAR-CCM+, one has to obtain the systematically varied propeller geometry first. This section will introduce the basic modeling procedures for parent propeller P1374 and its variants.
3.1.1 Systematically varied propellers
In the PROPSCALE project, the original parent propeller P1374 in model tests is a controllable pitch propeller (CPP), with its parameters can be found in Table 3.1. In the systematic CFD analyses, propellers derived from P1374 are divided into the four series: "CPP Pitch Series", "Blade Area Ratio Series"," Skew Series", and "Blade Number Series". In all series, except the first one ("CPP Pitch Series"), propellers are considered as fixed-pitch propellers (FPP). "Blade Area Ratio Series" with varied blade area ratio including 0.4, 0.6 and 0.8 are investigated in this paper.
Table 3.1 Parameters of parent propeller P1374
Blades Diameter P(0.7)/D BAR Skew Hub ratio
4 0.25 [m] 1.1 0.6 23 [deg] 0.24
Preliminary calculations done with the parent propeller P1374 have shown quite heavy loading on the outer blade sections, resulting in strong tip vortex. This result is thought to be related to the radial distributions of chord length and pitch at the outer blade sections, which may not be typical for conventional open propeller designs (it should be remembered that propeller P1374 was conceived as a compromise design to be used in the tests with both open and ducted propulsors). Obviously, the aforementioned phenomena may have considerable influence on scale effects. Therefore, it is planned to include in the investigations one alternative distribution of chord length.
20
3.1.2 Geometry modeling
The GAMBIT software is used to produce the initial blade model as a solid (Parasolid format). Then, this initial model can be re-meshed in STAR-CCM+ by a special Java macro, in order to fix minor surface flaws and ensure a better quality of subsequent Boolean operations with geometry parts. Table 3.2 lists all the geometrical elements used to model the parent propeller P1374, using cylindrical section definition.
When generating the geometry of fixed pitch series propellers the following conditions apply:
Distributions of max section camber 𝑓𝑜/R and rake 𝑥𝑟/R are kept unchanged through series;
Chord length distribution 𝑏𝑟/b(0.7) is kept unchanged through series;
Distribution of max thickness to chord length ratio 𝑒𝑜/𝑏𝑟 is kept unchanged through series (whereas 𝑒𝑜//R varies according to given 𝐴𝐸/𝐴𝑜 and Z);
Distributions of pitch 𝑃𝑟/P(0.7) and skew 𝑇𝑠/𝑇𝑠(max) are based on the parent propeller.
Table 3.2 Geometrical elements of parent propeller P1374
r/R 𝒃𝒓/R 𝒆𝒐/R cs/R 𝑷𝒓/R 𝒇𝒐/R
0.240 0.2668 0.0761 0.0000 2.1560 0.0024
0.250 0.2827 0.0750 0.0036 2.1579 0.0074
0.300 0.3599 0.0694 0.0217 2.1667 0.0141
0.350 0.4327 0.0640 0.0393 2.1745 0.0175
0.400 0.5008 0.0588 0.0557 2.1813 0.0196
0.500 0.6208 0.0490 0.0814 2.1917 0.0216
0.600 0.7151 0.0400 0.0899 2.1979 0.0214
0.700 0.7758 0.0318 0.0707 2.2000 0.0198
0.800 0.7879 0.0244 0.0115 2.1764 0.0169
0.900 0.7128 0.0178 -0.1014 2.0753 0.0125
0.975 0.4903 0.0134 -0.2302 1.9322 0.0065
0.990 0.3674 0.0125 -0.2612 1.8958 0.0041
1.000 0.1500 0.0120 -0.2829 1.8700 0.0000
r/R 𝒃𝒓/b(0.7) 𝒆𝒐/𝒃𝒓 𝑻𝒔 [deg] 𝑷𝒓/P(0.7) 𝑻𝒔/𝑻𝒔(max)
0.240 0.3439 0.2852 0.0000 0.9800 0.00000
0.250 0.3644 0.2653 0.4856 0.9809 0.02111
0.300 0.4639 0.1928 2.7202 0.9849 0.11827
0.350 0.5577 0.1479 4.5747 0.9884 0.19890
0.400 0.6455 0.1174 6.0255 0.9915 0.26198
0.500 0.8002 0.0789 7.6501 0.9962 0.33261
0.600 0.9217 0.0559 7.4164 0.9991 0.32245
0.700 1.0000 0.0410 5.1755 1.0000 0.22502
0.800 1.0156 0.0310 0.7558 0.9893 0.03286
0.900 0.9188 0.0250 -6.0601 0.9433 -0.26348
0.975 0.6319 0.0273 -12.9012 0.8783 -0.56092
0.990 0.4736 0.0340 -14.4602 0.8617 -0.62870
1.000 0.1933 0.08000 -15.5355 0.8500 -0.67546
21
Fig.3.1-3.4 present the characteristics of the propeller P1374, including maximum camber distribution, pitch distribution, etc. For propellers having alternative chord length distribution, maximum thickness to chord length ratio 𝑒𝑜/b is modified compared to the parent propeller, in order to ensure a more realistic blade geometry at the outer blade sections; Fig.3.3 and 3.4 plots the difference. Detailed geometrical data of this variant propeller can be found in Appendix A.
Figure 3.1 Max camber distribution Figure 3.2 Pitch distribution
Figure 3.3 Chord length distribution Figure 3.4 Max thickness distribution Propellers with varied BAR are derived by changing the chord length to radius ratio 𝑏𝑟/R. Max thickness is altered, however max thickness to chord length ratio 𝑒𝑜/b is same through the series. Apart from these, all other geometrical elements were kept identical. The obtained geometry of the series propellers in a general view can be found in Table 3.3.
Table 3.3 General view of propeller with varied blade area ratio
BAR=0.4 BAR=0.6 BAR=0.8
D=0.25m; Z=4; Sew 23°; P(0.7)/D=1.1
0.000 0.005 0.010 0.015 0.020 0.025
0.2 0.4 0.6 0.8 1.0
parent r/R
0.0 0.5 1.0 1.5 2.0 2.5
0.2 0.4 0.6 0.8 1.0
parent r/R
0 0.2 0.4 0.6 0.8 1
0.2 0.4 0.6 0.8 1.0
parent alternative r/R
0.00 0.02 0.04 0.06 0.08
0.2 0.4 0.6 0.8 1.0
parent alternative r/R
22
Figure 3.5 presents the geometry of parent propeller and its variant with alternative chord length distribution, one can easily tell the difference between the two propellers.
Figure 3.5 Parent propeller with BAR 0.6 and Propeller with alternative chord length distribution
Blade Tip
As one can find in Fig.3.6, there is a sharp edge in the blade tip area marked by feature curves. This is because the blade tip used in the simulations is just the cylindrical blade section at the radius 1.0R. For realistic propellers, the blade tip is normally faired to round off sharp edges. However, since the definition of such rounding is quite uncertain (and will surely be different for propellers having different blade area ratio and skew) the initial tip geometry is used in these systematic calculations. The feature curves describing the sharp edge is critical to control the mesh quality in the blade tip area. It will be discussed later in Part 3.2.2 meshing->feature curve.
Figure 3.6 Blade tip Figure 3.7 Propeller hub
Propeller Hub
The propeller hub is defined according to the hub ratio with diameter D=0.06m in model scale. In order to ensure correct intersection between blade and hub, the blade was extended towards inner radii on purpose. One can find the hub geometry in Figure 3.7.
23
3.2 Simulation setup
This section covers the procedure for carrying out simulations of the series propellers. Great thanks to Vladimir, from MARINTEK, who provided many necessary guidelines and advices in this part.
3.2.1 Domain modeling
When modelling propeller in straight flow in open water conditions, one can take advantage of flow axial symmetry, and use only one blade passage domain with appropriate periodic boundaries. The most straightforward setup for one blade passage flow simulation implies the use of a sector, having angular dimension 360/Z deg. (where Z is the number of propeller blades), which is cut from a cylinder and includes one whole blade, as shown in Fig.3.8. Such a setup is also simplest for post-processing.
However, if propeller blades are wide, they may not be entirely accommodated in the domain as described above.
Figure 3.8 One-blade-passage domain Figure 3.9 Two-blades-split domain
The simplest way to solve the problem is to use an alternative one blade passage setup that includes the same cylindrical sector, but instead splits two neighbouring blades, as shown in Fig.3.9. Such a setup will ensure that complete blade geometry will be accommodated in the one blade passage domain, and flow periodicity will be observed. The present work utilized One-blade-passage domain to simulate all the cases. At the same time, for comparison, propeller with BAR 0.6 in model scale will be simulated once again in the Two-blades-split domain.
A gap between propeller hub and shaft was created to replicate the physical setup of the open water model test. Practice shows that only by including the gap between the rotating propeller hub and stationary shaft housing one can achieve adequate prediction of the axial force acting on the hub. In the present simulations, propeller thrust includes blade thrust and hub thrust, exactly as measured in the tests. Blade thrust is less dependent on whether one includes the gap in the simulation or not. However, hub force may not be predicted correctly in the setup without gap (blades on shaft, part of which is separated to form the hub surface) because of the uncertainty in the definition of gauge pressure on open surface.
24
3.2.2 Meshing
A mesh is the discretized representation of the computational domain. STAR-CCM+ provides meshers and tools that can be used to generate a quality mesh for various geometries and applications, Reference [9]. In the present simulation, the Surface Remesher, Polyhedral Mesher, Prism Layer Mesher and Extruder are applied to generate the volume mesh. Base size in model scale is set to 0.25m equals to propeller diameter. While in full scale, it has to be augmented by the scale factor 20.
The idea behind the base size is that one can set other meshing values like surface minimum / target size, Prism layer thickness e.g. proportionally to the base value. All the mesh values which are defined by the ratio to the base size will be changed automatically in full scale.
Surface Remesher
In order to improve the overall quality of an existing surface and optimize it for the volume mesh models, the surface remesher was used to re-triangulate the surface. Parameters including relative target size and relative minimum size are used to control the surface mesh size in particular region.
Table 3.4 Parameters of surface remesher (% to base size) Relative target size Relative minimum size
Blade surface 0.5% 0.125%
Tip 0.125% 0.125%
Trailing edge 0.125% 0.125%
Hub and Shaft 3% 1%
Polyhedral Mesher
It generates a volume mesh that is composed of polyhedral-shaped cells. The polyhedral meshing model utilizes an arbitrary polyhedral cell shape in order to build the core mesh. The polyhedral core mesh density can be increased or decreased by using the volume mesh blending factor, which is set to be 0.5 at the present case.
Prism Layer Mesher
The prism layer mesh model is used with a core volume mesh to generate orthogonal prismatic cells next to wall surfaces or boundaries. Prism layers are mainly used to resolve flow boundary layers and they are critical to improve the accuracy of the flow solution.
Table 3.5 Parameters of Prism layer mesh
Base size Prism layer thickness Number of prism layers Prism layer stretching
0.25m 0.25% base size 10 1.4
The near wall 𝑦+ is very sensitive to the prism layer mesh setting. The automatic mesh setup used in the simulations should allow for adequate near-wall (boundary layer) flow treatment in both the model scale and full scale conditions without modification of Prism Layer Mesher setting. In model scale, it will result in wall 𝑦+<5, while in full scale, it will result in wall 30<𝑦+300. For this reason, the present mesh
