Numerical Fluid Simulations of a Multi- Rotor Propeller A validation study for experimental propeller data and identi- fication of vortex sources Bachelor Thesis in Mechanical Engineering Emil Andersson Oliver Sandberg Alrik Sjökvist Examiner: Niklas Andersson Supervisor: Debarshee Gosh Department of Mechanics and Maritime Sciences CHALMERS UNIVERSITY OF TECHNOLOGY Gothenburg, Sweden 2024 www.chalmers.se www.chalmers.se Abstract For centuries, propellers have been used to propel ships, in the beginning through water, and later they began being used to propel aircraft through the air. Pro- pellers increase the pressure in the fluid and are beneficial for low-speed aerospace applications. Today, a common use for propellers is in drones, where the concern for noise emissions arises due to their popularity and in-city use. Therefore, this thesis studied the aerodynamic performance of a 14x8-inch propeller by Advanced Precision Composites and the vortices generated by it using Computational Fluid Dynamics simulations solving steady-state Reynolds Averaged Navier-Stokes equa- tions. The results from these simulations are then compared to experimental data from the University of Illinois Urbana-Champaign. The results of the simulations showed the model most closely following the experi- mental results for low advance ratios with an error for efficiency around 10%. The simulations were performed with a rotational rate of 4038 rpm and inlet velocities ranging from 10.4 m/s to 17.7 m/s. The simulation model also showed the most common occurrence of vortices being from the hub and the tip of the blades. A reduction in tip vortices would be beneficial to reduce the aero-acoustic noise gen- erated by the propeller. i Acknowledgements We want to thank Chalmers University of Technology and the Department of Mech- anics and Maritime Sciences for allowing us to undertake this project, along with the Chalmers Center for Computational Science and Engineering, for providing us with the computational resources. We also want to thank our supervisor, PhD student Debarshee Gosh, for his contin- ued patience and guidance during this project and our examiner, Professor Niklas Andersson, for helping put this project together. Lastly, we thank Professor Thomas Sebastian from MIT for providing information regarding their toroidal propeller. Emil Andersson Oliver Sandberg Alrik Sjökvist Gothenburg, May 27th 2024 ii Contents 1 Introduction 1 1.1 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Purpose and Objectives . . . . . . . . . . . . . . . . . . . . . . . . . . 4 1.3 Limitations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 2 Theory 5 2.1 Propeller physics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.1.1 Coordinate system and Velocity triangles . . . . . . . . . . . . 5 2.1.2 Blade Element Theory . . . . . . . . . . . . . . . . . . . . . . 7 2.1.3 Performance parameters . . . . . . . . . . . . . . . . . . . . . 9 2.1.4 Noise generation . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.2 Airfoil geometry . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.3 Computational Fluid Dynamics . . . . . . . . . . . . . . . . . . . . . 11 2.4 Governing equations . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.5 Moving Reference Frame . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.6 Turbulence Models . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.7 Wall-Y+ Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 2.8 Boundary conditions . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 2.9 Meshing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 2.9.1 Mesh Refinement . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.9.2 Prism Layers . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.9.3 Mesh Quality . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.9.4 Mesh Independence Study . . . . . . . . . . . . . . . . . . . . 17 3 APC 14x8 Propeller 19 3.1 Experimental Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 3.2 Geometry . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 3.3 Experimental Validation . . . . . . . . . . . . . . . . . . . . . . . . . 21 4 Method 22 4.1 Simulation Setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 4.1.1 CFD Domain . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 4.1.2 Mesh Generation . . . . . . . . . . . . . . . . . . . . . . . . . 23 4.1.3 Boundary Conditions . . . . . . . . . . . . . . . . . . . . . . . 25 4.1.4 Physics Models . . . . . . . . . . . . . . . . . . . . . . . . . . 25 4.1.5 Stopping criteria . . . . . . . . . . . . . . . . . . . . . . . . . 25 iii Contents 4.2 Refining Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 4.2.1 Turbulence Models . . . . . . . . . . . . . . . . . . . . . . . . 26 4.2.2 Mesh type . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 4.2.3 Wall-Y+ testing . . . . . . . . . . . . . . . . . . . . . . . . . . 26 4.2.4 Mesh Independence Study . . . . . . . . . . . . . . . . . . . . 27 4.3 Post-processing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 5 Results 29 5.1 Mesh Type Testing . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 5.2 Turbulence Model Testing . . . . . . . . . . . . . . . . . . . . . . . . 29 5.3 Wall-Y+ Testing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 5.4 Mesh Independence Study . . . . . . . . . . . . . . . . . . . . . . . . 30 5.5 Presentation of Results . . . . . . . . . . . . . . . . . . . . . . . . . . 31 5.5.1 Experimental vs Simulated Results . . . . . . . . . . . . . . . 31 5.5.2 Residuals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 5.5.3 Wall-Y+ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 5.5.4 Velocities . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 5.5.5 Pressure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 5.5.6 Swirl Angle . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 5.5.7 Vorticies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 6 Analysis and Discussion 38 6.1 Result Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38 6.2 Methodology Discussion . . . . . . . . . . . . . . . . . . . . . . . . . 40 6.3 Future Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 A Geometry Data I B Experimental propeller graphs II C Experimental propeller data III iv 1 Introduction This chapter introduces propellers, their uses, and the current research. It also includes the purpose of this thesis and possible limitations in its completion. 1.1 Background A propeller is, as the name suggests, a propulsion device. It is found in several domestic, industrial, and military applications. They are used for both marine and aircraft propulsion. In the aerospace industry, propellers are used primarily to propel small-to-medium-sized aircraft. More recently, propellers have also been used to propel drones and other small Unmanned Aerial Vehicles (UAVs)[1]. A propeller is a turbomachine. A turbomachine is a device that transfers energy to or from the flow by the dynamic action of the rotor blade on a continuously moving fluid. Depending on the direction of energy transfer, turbomachines are classified into turbines or pumping devices. A turbine is a turbomachine that produces power by expanding the fluid flowing through it to a lower pressure. A pumping device, instead, absorbs power and increases the pressure or head of the flow. A propeller is a power-consuming turbomachine, i.e. it requires energy to perform work on the flow. The flow may be either a liquid or a gas, like water or air. A propeller consists of two or more rotor blades attached to a common hub, which is connected to a powered shaft. The rotor blades rotate around this hub and displace air by performing work on it, producing thrust. Thrust is a force directed in the direction parallel to the propeller’s rotational axis. The propeller is an old machine part used as a propulsion device for several cen- turies. The regular fan is a commonly used propeller-like structure, often found in everything from household appliances and electronics to a standalone human cooling device. As early as the 18th century, ideas about using four-bladed screws to propel ships were proposed [2]. In the 1830s, the first propeller ships were patented. In Fig. 1.1a an example of a propeller used for marine applications can be seen. In the early days, propellers were primarily used for ships and lighter-than-air aircraft. It was not until the early 20th century and the success of the Wright brothers that a rapid progression was observed in the use of propellers for heavier-than-air aircraft 1 1. Introduction [2], such a propeller can be seen in Fig 1.1b . Several commercial aircraft today still use the Pratt & Whitney Canada PT6 turbo-prop engine, among others, for short-haul flights [3]. (a) Propeller used for marine applic- ations, found outside the SSPA Mari- time Center located in Gothenburg. (b) Propeller used in aerospace applic- ations [4]. Figure 1.1: Different types of propellers. When propelling an aircraft, a choice has to be made between propeller-driven or jet-driven systems. A propeller serves the same purpose as a compressor rotor, i.e. increases the fluid’s pressure flowing through it. A jet engine ingests air and compresses it to a higher pressure and temperature. The high-pressure and high- temperature air is then forced into the combustor, where fuel is combusted. The hot air then passes through a turbine where the flow is expanded to a lower pressure and consequently a higher velocity. The high-velocity gas exits the engine and propels the aircraft forward [5]. Propellers are the better choice for low-speed flights in terms of efficiency [6]. Prob- lems arise when the aircraft’s and, subsequently, the propeller’s speed increases. The rotor tip has a higher velocity than the section close to the hub. The blade velocity (U) is defined as the product of the rotational speed (ω) and blade section radius (r), according to Eq. 1.1. The blade section radius at the hub is smaller than at the blade tip. Consequently, higher blade velocity is achieved at the tip for the same rotational speed compared to the hub, see Fig. 1.2. 2 1. Introduction Figure 1.2: Velocity at different points of the blade. For blade tip speeds with a Mach number greater than 0.3, the flow is considered to be compressible. The Mach number (Ma) is a non-dimensional number, defined as the ratio of the flow velocity (u) to the speed of sound (a), see Eq. 1.2. For tip speeds with Ma > 1, supersonic flow is achieved at the tip and may lead to the formation of shock waves, which may lead to severe losses [6]. U = ω ∗ r (1.1) Ma = u a (1.2) With the current trend of electrification across industries, the efficiency of propellers for aircraft and drones is of the utmost importance. Propellers have around 90% effi- ciency [7]. However, researchers are trying to push that number higher by increasing thrust while minimizing drag [1]. Drones bring another problem to light. More precisely, their use in cities creates more noise pollution. This is due to the noise from the electric motors and the noise generated by the propellers [8]. This has caused research to focus on creating quieter drone propellers that still maintain a high level of aerodynamic efficiency [8]. A rotor’s aerodynamic or aeroacoustic performance may be evaluated either com- putationally or through experiments. Experiments require significant resources and financial investments and are not always feasible. Computational investigations are relatively cheaper and can serve as a valuable tool to predict the aerodynamic and aeroacoustic performance of most fluid systems. However, obtaining a reliable com- putational setup that yields physically accurate solutions can be challenging in its own right. 3 1. Introduction 1.2 Purpose and Objectives As indicated in Sec 1.1, developing low-noise drones is a key area of investigation associated with propeller rotor blade design. Computationally, performing aero- acoustic analysis requires significantly high computational resources and time. Given the resource dependency in this project, a complete aero-acoustic analysis is outside the scope of this work. However, an indication of the main noise sources may be obtained through the aerodynamic analysis of the rotor by identifying noise- producing vortices. The main purpose of this work is to investigate the aerodynamic performance of a multi-rotor propeller computationally. A validation study is performed to compare computational results against numerical data for a multi-rotor propeller. Specific- ally, this work investigates a two-bladed propeller designed by Advanced Precision Composites (APC), which was experimentally tested at the University of Illinois Urbana-Champaign (UIUC). Computations have been performed to evaluate the aerodynamic performance at different operating points of the propeller. 1.3 Limitations Computational simulations for turbomachines take considerable time and have many potential problems. Therefore, this project will limit itself to propeller simulations, and no physical experiments will be conducted. A complete aeroacoustic, Computational Fluid Dynamics (CFD) simulation requires significant computational resources. The high computational demand for aeroacous- tic simulations is a result of the requirement for high spatial resolution, i.e. highly refined grids with the ability to capture sound waves. Additionally, unsteady sim- ulations with significantly small time steps in the order of microseconds need to be performed to capture the sound waves. These resources are unavailable; therefore, this work is limited to aerodynamic simulations attempting to accurately predict the flow field. A reliable aerodynamic analysis is the first step to performing an aeroacoustic analysis. 4 2 Theory The following sections present relevant theories regarding propeller physics, appro- priate fluid dynamics, and CFD. 2.1 Propeller physics A propeller blade has an airfoil cross-section [9]. When the propeller rotates about its axis, it performs work on the surrounding fluid, changing its enthalpy. A propeller absorbs power and increases the pressure of the fluid flowing through it [5]. 2.1.1 Coordinate system and Velocity triangles Since propellers rotate around an axis and can move up and down, a cylinder can be identified, and therefore, cylindrical polar coordinates are appropriate to use [5]. The cylindrical coordinates are defined as radial, axial, and circumferential. Likewise, the velocity components are named in the same manner (cr, cax and cθ), which are illustrated in Fig. 2.1. In reality, a turbomachine has varying velocities in all three directions, but to ease calculations, constant circumferential velocity can be assumed [5], as seen in Fig. 2.1b. As demonstrated in Fig. 2.1a, the flow component along an axisymmetric stream surface is called the meridional velocity, defined in Eq. 2.1. cm = √ c2 ax + c2 r (2.1) The total flow velocity (c) comprises all three velocity components; see Eq. 2.2. However, an axial flow machine has flow primarily parallel to its rotational axis, while a radial flow machine has flow perpendicular to its rotational axis. This means that Eq. 2.1 changes if any of the components are of several orders more significant than the other; the meridional velocity is often defined as just the axial or the radial velocity [5]. 5 2. Theory (a) Meridional view. (b) Axial view. Figure 2.1: Velocity components in a turbomachine. Adapted from [5]. c = √ c2 ax + c2 r + c2 θ = √ c2 m + c2 θ (2.2) In turbomachinery, velocity triangles are often used to evaluate the flow in a two- dimensional plane. Velocity triangles define the velocity components and their angles with respect to the flow direction at the rotor inlet and outlet of a blade section. As mentioned in the previous section, a rotating airfoil section performs work on the fluid [5]. It does this by changing the flow direction, i.e. α2 > α1. An example in Fig. 2.2 shows the velocity diagram of an axial flow rotor for an airfoil cross- section. The blades rotate around the rotational axis with a velocity U (Eq. 1.1). The velocity in the absolute frame of reference at the rotor inlet is defined as c1 and has an angle α1 with respect to the rotational axis. Similarly, the velocity at the rotor outlet is defined as c2 and has an angle α2 with respect to the rotational axis. In the relative frame of reference, that is, a frame of reference that rotates with the same angular velocity as the rotor blades, the velocities at the rotor inlet and rotor outlet are defined as v1 and v2, respectively. The flow angles in the relative frame of reference at the rotor inlet and outlet are β1 and β2, respectively. 6 2. Theory Figure 2.2: Velocity triangle of an axial flow rotor. Adapted from [5]. Both in the absolute and relative frame of reference, the velocities consist of meridi- onal and circumferential velocity components. For the absolute velocity, the angle is called the swirl angle and is defined by the inverse tangent of the components as in Eq. 2.3. αi = tan−1 cθ,i cm,i (2.3) The work rate per second (Ẇc) the blade performs on the fluid is calculated with the Euler work equation (Eq. 2.4) [5]. The mass flow (ṁ) is multiplied by the blade velocity and the difference in the circumferential fluid velocity. Since this is for one blade section, the work rate can be summed or integrated over the entire blade. Ẇc = ṁU∆cθ = ṁU(cθ,2 − cθ,1) (2.4) 2.1.2 Blade Element Theory A propeller has the same characteristics as the rotating blades of an axial flow compressor, meaning the same theory applies to both cases. However, because of various rotational speeds, the blades are slightly twisted, giving, for instance, a lower angle of attack at the tip than at the hub, complicating calculations. A blade can be divided into elements using the cross-sectional fact, and each is individually studied. This is known as Blade Element Theory (BET) [10]. In Fig. 2.3 below, the velocity triangle and the forces acting upon the blade are demonstrated. 7 2. Theory Figure 2.3: Airfoil cross-section of a propeller blade with two-dimensional forces acting upon it. Adapted from [10]. The propeller rotates around the hub with a rotational velocity Ω, giving the blade element a circumferential velocity of vθ = Ωr, r being an arbitrary radius from the hub. The effective air velocity (veff ) is easily calculated with the Pythagorean theorem in Eq. 2.5 from the circumferential, axial (vax) and induced velocity (vi). veff = √ (vax + vi)2 + v2 θ (2.5) The effective velocity is used when calculating the differential lift and drag forces, dFL and dFD, per Eq. 2.6 and 2.7 where c(r) is the element chord, dr the width and CL and CD are the lift and drag coefficients. dFL = 1 2ρv2 effc(r)CLdr (2.6) dFD = 1 2ρv2 effc(r)CDdr (2.7) The differential thrust force and torque can be calculated from dFL and dFD, ac- cording to Eq. 2.8 and Eq. 2.9, thanks to the lift force being perpendicular to the effective velocity and the drag force being normal to it [10]. Notice that the force is multiplied by the radius to get the torque. dT = dFLcos(β − α) − dFDsin(β − α) (2.8) 8 2. Theory dQ = r ∗ dFQ = r ∗ [dFLsin(β + α) + dFDcos(β + α)] (2.9) The whole propeller’s thrust and torque can be obtained by integrating the above equations according to Eq. 2.10 and 2.11, where N is the number of blades. T = N ∫ R Rhub dT (2.10) Q = N ∫ R Rhub dQ (2.11) 2.1.3 Performance parameters Interesting aspects of a propeller are the advance ratio, (J) (Eq. 2.12), which is defined by the velocity (V ), the rotational speed in revolutions per second (n), and the propeller diameter (D). It can be described as the ratio between the upward motion and the rotational motion of the propeller. In the absolute majority of small drone propellers, the advance ratio is in the range of zero and one [11]. It can be bigger than one; however, in order to increase the velocity, more thrust needs to be created, and thus, a higher rotational speed is needed, which in turn lowers the advance ratio. Thrust and power coefficients (CT and CP ) (Eq. 2.13, 2.14) are defined by the thrust and power, respectively, the fluid’s density (ρ), rotational speed in revolutions per second, and the propeller diameter. The efficiency (η) (Eq. 2.15) is the ratio between CT and CP [11]. For propellers, the efficiency usually maxes out at around 0.8 [7]. All four are dimensionless parameters of a propeller’s performance. The power is calculated by multiplying the torque from Eq. 2.10 with the rotational speed (Ω) in radians per second, as per Eq. 2.16 [11]. J = V nD (2.12) CT = T ρn2D4 (2.13) CP = P ρn3D5 (2.14) η = CT CP J (2.15) 9 2. Theory P = Q ∗ Ω (2.16) 2.1.4 Noise generation One of the major contributing factors to propeller noise is the vortices generated by the propeller tips. So if one could reduce the strength of these vortices, one should ideally be able to reduce the amount of sound generated [8]. The tip and hub are where the major vortices for propellers appear [12]. The tip has the biggest impact on noise generation between the two [13]. 2.2 Airfoil geometry An airfoil could be summed up as a stretched oval. This gives a Leading Edge (LE), the most forward part, and a Trailing Edge (TE), the part furthest back of the blade. Between these points exists a straight line called the chord line, with a length called the chord. The airfoil is described as having two-dimensional x- and y-coordinates, with the origin in the LE. Various parameters are calculated at and between different points on the airfoil [10]. They will not be explained in detail, but some of them are: · Thickness: the maximum between two vertical points at the same x- coordinate, · Twist or geometrical pitch: the angle between the chord line and the hori- zontal plane, · Thickness-to-Chord ratio: the thickness divided by the chord. Some of the most used standardized airfoil geometries are the NACA and Clark-Y profiles. The NACA geometry was developed in the 1930s and was completed with a total of 78 different airfoils [14]. Model names range from a four to an eight-digit description. The Clark-Y profile is characterized by a flat lower surface with the chord line and the lower surface merging just before the TE [10]. Figure 2.4 shows the difference between a NACA and a Clark Y airfoil. Figure 2.4: NACA0012 (Green) and Clark-Y (Blue) airfoils. 10 2. Theory 2.3 Computational Fluid Dynamics To obtain the flow field, the governing flow equations, which will be outlined in Sec. 2.4, are discretized in space to obtain a set of linearized equations, which can then be solved using an iterative solver such as Star-CCM+ [15]. As with numerical solving, certain assumptions are made, and boundary conditions must be set. The assumptions can be about the type of fluid flow. CFD programs use discretization methods such as finite-volume and finite-element, allowing calcu- lations on complex geometries and surrounding fluids. The discretization is referred to as meshing. The grid that arises is linear approximations, which are much more easily solved in the form of systems of linear equations [16], [17]. 2.4 Governing equations Figure 2.5: Infinitesimal, cubic Control Volume with equilibrating mass flow. Adapted from [18]. This section describes the fun- damental equations, known as the governing equations, which are the conservation of mass, linear momentum, and energy. Since these equations are in an introductory fluid mechan- ics course at university, they will only be covered superfi- cially, and their derivations will only be discussed. For all the equations, they are defined over an infinitesimal, cubic Control Volume (CV) as seen in Fig. 2.5 [18]. The conservation of mass equa- tion, also known as the continuity equation, is defined in Eq. 2.17. It describes the change in density (ρ) over time and space and the velocities (V) change in each direction [18]. ∂ρ ∂t + ∇ · (ρV) = 0 (2.17) With a steady-state and an incompressible flow assumption [18], the fluid’s density does not change, and the equation reduces to Eq. 2.18. ∇ · V = 0 (2.18) 11 2. Theory The linear momentum equation is defined in Eq. 2.19. The introduction of body forces due to gravity (g) and the surface forces (∇p) due to stresses on the CV surfaces give rise to a viscous stress tensor (τij). ρg − ∇p + ∇ · τij = ρ dV dt (2.19) The conservation of energy equation is derived in the same manner as the continuity and momentum equations and is defined in Eq. 2.20. ρ dû dt + p(∇ · V) = ∇ · (k∇T ) + Φ (2.20) Where dû dt is the change of internal energy, k is the thermal conductivity, ∇T is the change in temperature, and Φ is the viscous-dissipation function. It is important to note that Eq. 2.20 is only valid for a Newtonian fluid under general conditions and neglects radiation heat transfer along with internal heat sources [18]. 2.5 Moving Reference Frame A Moving Reference Frame (MRF) approach can be used to simulate rotation in CFD. The MRF does not move the mesh directly, but instead, a relative velocity term is inserted into the governing equations to simulate the motion [19]. This approach reduces the complexity and computational power required for rotating simulations since the mesh can be kept stationary, making it possible to solve the problem as steady-state [20]. 2.6 Turbulence Models Laminar fluid flow at low Reynolds numbers is quite easily modeled, even analyt- ically. However, many flows around an immersed body are not of this character but instead of a turbulent one with a higher Reynolds number, as is the case with the flow around an object or a propeller pulling and pushing on the fluid medium. The need for a numerical approach to solve these problems is considerable for an engineer. This is done using turbulence models [17], [21]. At the occurrence of turbulence, fluctuations in the fluid’s velocity are significantly noticed. Instead of defining complex functions, the velocity is decomposed in a mean, ū, and the fluctuations as a function of time, u′(t). This is known as Reynolds decomposition, as per Eq. 2.21. 12 2. Theory u(t) = ū + u′(t) (2.21) By applying time averaging on the Reynolds decomposition and, in turn, apply- ing this to Eq. 2.19, time-averaged Navier-Stokes are given. Known as Reynolds Averaged Navier-Stokes (RANS), it is a powerful tool in turbulence modeling. In the RANS equations, extra terms occur due to turbulent fluctuations interacting. These are modeled with known models like the k − ϵ, k − ω, Shear Stress Transport (SST), and Spalart-Allmaras (SA) methods. These turbulence models have different strengths and weaknesses, making them applicable in different situations [17], [18], [21]. In general, the k − ϵ method performs well in computations in the free stream region but poorly in the near-wall region. Conversely, the k − ω does the opposite, performing well in the near-wall region but poorer in the free stream region. The SST method was suggested by Florian Menter and is a combination of both the k − ϵ and Wilcox k − ω methods. The SA method differs from the others in that it involves only one transport equation while the others solve for two unkowns. How these extra terms and equations are derived is on a much higher mathematical level than required to understand this thesis’s purpose and will therefore not be presented here [17], [18], [21]. 2.7 Wall-Y+ Functions A fluid flowing past bodies is exposed to viscous shear forces. This creates an increasing velocity profile in the near-wall region, with the relative velocity of the fluid at the wall being zero. This is called the no-slip condition, and the layer that emerges is called a boundary layer. In the aerodynamics of aviation, the build-up of boundary layers on a surface is of great importance because of the effect these have on flow separation, wake propagation, and ultimately drag and lift performance [18]. A turbulent boundary layer is divided into three regions: the viscous sublayer, the overlap layer, and the outer layer from the nearest to the furthest away from the wall. The viscous sublayer is dominated by viscous shear and is considered laminar; in the outer layer, turbulent shear dominates and is considered turbulent. In the overlap layer, both viscous and turbulent shear are of great notice. The profile of the boundary layer depends on whether the flow is laminar or turbu- lent. In CFD, the boundary layer that is formed as a result of the no-slip condition may be either resolved or modeled. To resolve a boundary layer, a significantly smal- ler grid size has to be used at the walls. The other option is to use mathematical equations to model the velocity profile at the wall. The velocity of the boundary layer is expressed as a dimensionless velocity (u+) defined by the real velocity (u) over the friction velocity (u∗), which in turn is defined by the wall shear stress (τw), and the density (ρ). Equation 2.22 and 2.23 are called the law of the wall and are 13 2. Theory shown in Fig. 2.6. u+ = u u∗ (2.22) u∗ = √ τw ρ (2.23) Figure 2.6: Log graph of the Law of the Wall. Adapted from [22]. The dimensionless wall distance (y+) closest to the wall is then defined by Eq. 2.24, with the real wall distance (y) and kinematic viscosity (ν), up to y+ ≈ 5. It then merges with the overlap layer and the logarithmic Eq. 2.25 takes over after y+ ≈ 30, where κ ≈ 0.41 is the von Kárman constant and B ≈ 5. y+ = yu∗ ν (2.24) u u∗ = 1 κ ln y+ + B (2.25) Depending on the choice of wall-y+, the boundary layer may be either resolved or modeled using mathematical functions. Typically, for wall-y+ values greater than 11, the boundary layer is modeled, and for wall-y+ values below one, it is resolved using wall functions. The details of how wall functions are implemented are outside the scope of this thesis. 14 2. Theory 2.8 Boundary conditions There are seven different kinds of boundary conditions in Star-CCM+. The most common of these are velocity inlet, free stream, wall and pressure outlet. Velocity inlet specifies the fluid velocity when entering the simulation domain. Free stream needs a large enough space between an object and the boundary and is specified by Mach number, static pressure and static temperature. Wall is typically for the surfaces of an object, this could be set to a slip or no-slip condition. However, it could also be used instead of a free stream condition on the outer boundary. A pressure outlet specifies the pressure at the end of the simulation domain. If a reference pressure has been set previously to, for example, 1 atm, then the pressure at the outlet should be set to 0 [20]. 2.9 Meshing The mesh is one of the most important parts of a CFD simulation. The fluid region surrounding the geometry of interest is discretized into a large number of small cells. The governing equations are solved for each of these small cells iteratively [17]. In Star-CCM+, there are mainly three different types of meshes: tetrahedral, poly- hedral, and trimmed cell meshes. The tetrahedral mesh is created from tetrahedral shapes forming a triangular surface boundary like the one in Fig. 2.7a. The poly- hedral mesh, like the one in Fig. 2.7b, uses arbitrary polyhedral shapes created from an underlying tetrahedral mesh with an average of 14 cell faces per cell. The final mesh called trimmed cell mesh that can be seen in Fig. 2.7c, is a square mesh that is based on a hexahedral template mesh from where it trims, hence the name, the cells using the boundary surface [20]. (a) Tetrahedral mesh. (b) Polyhedral mesh. (c) Trimmed cell mesh. Figure 2.7: Different mesh types. 15 2. Theory 2.9.1 Mesh Refinement Certain regions of the domain are of more interest, such as the fluid closely sur- rounding the propeller. In that case, this particular region can be refined to better capture the fluid flow in said region. By using this refinement, the mesh can be kept coarse in the regions where nothing of interest is expected to happen. Doing this reduces the number of cells and, consequently, the computational power and time. [23]. 2.9.2 Prism Layers Figure 2.8: Example of a prism layer. To better capture the boundary layers in the near-wall regions, a prism layer mesh is used. To achieve this, the mesher creates orthogonal prismatic cells like the ones in Fig. 2.8 along the surfaces [24]. There are a few different settings that one needs to consider in order to get good prism layers. The most import- ant ones are gap-fill percentage, minimum thickness percentage, and layer reduction percentage. Gap-fill percentage describes how much of a gap between two surfaces may be filled by prism layer cells. The gap fill percentage is in the range 0-50, where 50 means that the prism layer from each surface may cover 50% of the gap [20]. The minimum thickness percentage determines whether the prism layers should collapse completely at a certain point. This is decided by checking if the prism layer height in an area is less than a specified percentage of the total prism layer height; if so, the prism layer will collapse completely. The layer reduction percentage removes individual layers, as many as required, to fulfill requirements regarding total prism layer height. 2.9.3 Mesh Quality To ensure good simulation results, the mesh needs to be of high quality. This is ensured by inspecting the mesh to make sure the cells are up to the right standards according to certain metrics. These metrics are cell quality, face validity, skewness angle, and volume change [20]. 16 2. Theory (a) Cell quality. (b) Face validity. (c) Skewness angle. (d) Volume change. Figure 2.9: Visualisation of mesh quality metrics. Adapted from [20] Cell quality describes the geometric "fit" to its neighboring cells, see Fig. 2.9a. Per- fect orthogonality would give a metric of 1.0. For StarCCM+, cells with a cell quality of less than 1.0E-5 are considered bad. Face validity is an area-weighted metric that describes if the cell face normal is pointing out from the cell center, see Fig. 2.9b. Values below 1.0 are considered bad, meaning no faces should point slightly inward at all. The Skewness angle is the angle (θ) at the crossing point of a face normal vector between two cells and the vector connecting its two centers, see Fig. 2.9c. It should not be greater than 85 degrees because then it would negatively impact the diffusion of quantities passing these cells. Volume change is the ratio between a cell’s volume and the volume of its largest neighbor, see Fig. 2.9d. This may be the most understandable and logical metric. If this ratio becomes too small, inaccurate jumps in the values of the quantities will occur, giving an inaccurate result. Ratios below 0.01 are considered bad cells. 2.9.4 Mesh Independence Study Only the boundary conditions and physics model should affect the results when do- ing CFD simulations. Due to this, a mesh independence study always needs to be conducted to ensure that the mesh does not affect the results. For a mesh independ- ence study, three meshes of substantially different number of cells are generated. The three meshes are labeled Coarse, Medium, and Fine. The coarsest mesh has the least 17 2. Theory number of cells, and the finest mesh has the most. The number of cells is increased by at least 2.2 times from a coarser to a finer grid. Finally, the mesh independence is established using the Grid Convergence Index (GCI) method [25]. This method compares the results for the different refinements and gives an indication of if the results are within an asymptotic range of convergence. The GCI result describes how far the simulated results are from the asymptotic result when the grid resolu- tion approaches zero. The thrust coefficient (CT ), the pressure coefficient (CP ), and the efficiency (η) are the parameters that are evaluated for the mesh independence. Ideally, mesh independence can be achieved if minimal differences are observed for CT , CP , and η for the Medium and Fine mesh. 18 3 APC 14x8 Propeller The APC 14x8-inch propeller is a low-speed propeller from APC [26]. It has a 14- inch (35.56 cm) diameter, and the 8-inch (20.32 cm) pitch describes the theoretical distance traveled in the axial direction in one revolution. UIUC has performed experimental wind tunnel testing on many of APC’s wide range of propellers [11]. Previous work has recreated the experimental data of some of the propellers, for instance, Ahmed and Rajendran [27]. However, prior validation of the 14x8 propeller has not been found. This chapter presents the experimental data from UIUC on the 14x8 propeller this thesis is recreating. 3.1 Experimental Data In App. B, the experimental results from UIUC are shown. From Fig. B.1, it is gathered that the highest efficiency is between advance ratios 0.55 − 0.60 with a steep drop after 0.65. Figure B.2 shows the continued decline of CT and CP , with the former being more aggressive. The data is shown in tabular format in App. C. 3.2 Geometry A two-bladed multirotor propeller with a Clark-Y cross-section was created using geometry data available in App. A from APC [26]. The different sections were imported into Solidworks, and a solid blade was extruded. Figure 3.1 shows a selection of the different sections. 19 3. APC 14x8 Propeller Figure 3.1: Line plot of the propeller blade. The propeller is 355.6 mm in diameter with a hub diameter of 34 mm. The difference in twist along the propeller is 28.7 degrees. To connect the blade to the hub, an airfoil section of the same size as the innermost section from the data with zero degrees of twist was placed in the center of the hub, and the two sections were twisted together. The distance of 38.3 mm in Fig. 3.1 above is where the hub has completely transitioned into the blade. An ISO view of the propeller can be seen in Fig. 3.2. Figure 3.2: ISO view of the propeller generated from APC data, see App. A. 20 3. APC 14x8 Propeller 3.3 Experimental Validation This section briefly introduces the computation for comparing experimental results and simulations. The error for η, CT , and CP is calculated using Eq. 3.1, 3.2 and 3.3 to the difference in percent. eη = |ηsim − ηexp| ηexp ∗ 100 (3.1) eCT = |CT,sim − CT,exp| CT,exp ∗ 100 (3.2) eCP = |CP,sim − CP,exp| CP,exp ∗ 100 (3.3) 21 4 Method This chapter outlines the steps to meet the objectives outlined in Sec. 1.2. The first step in any computational investigation is to determine the validity and accuracy of the numerical setup. Thus, a numerical validation study is performed. This study compares the experimental data from UIUC against numerical data. The simulations performed within the scope of this work are three-dimensional, steady-state, RANS analyzes. The investigation of the numerical accuracy of a computational setup involves several steps. First, a mesh-type test, testing a poly- hedral and a trimmed cell, is conducted. Then, a turbulence model study is con- ducted. Four different turbulence models, namely the k − ϵ, Wilcox k − ω, SST, and SA turbulence, have been investigated. After, an investigation is made into the effect of wall-y+ on the numerical solution. Lastly, a mesh independence study was performed to ensure the negligible effect of the mesh on the numerical solution. Numerical results for comparison have been obtained for the experimental Best Effi- ciency Point (BEP). Later, with a CFD model reflecting the experimental situation, the rest of the operating points could be simulated. All simulations have been performed using the commercially available software Star- CCM+ - 18.07 [15]. 4.1 Simulation Setup This section describes how the simulation is set up. In CFD, a proper physics model must simulate the right conditions. The different models require different amounts of computational power. For instance, looking at a symmetrical airfoil in nice and calm conditions, the data could be well supplied with two-dimensional, steady-state simulations. Conversely, a turbomachine moves the fluid at high velocities; with changing pressure and complex geometries, a more complex model may be better suited. These considerations are up to the engineer to decide on the task they have been provided with [17]. 22 4. Method 4.1.1 CFD Domain The CFD domain consists of two regions that fulfill the requirements of a MRF. In Fig. 4.1a, a top-side view of the two regions is shown. Both regions are cylindrical. The static domain has a diameter of approximately 11.25Dprop, and the rotating do- main has a diameter of 1.14Dprop. Figure 4.1b displays a side view of the simulation domain. The height of the static domain is 14.06Dprop, and the rotating domain is 0.84Dprop. A boolean subtract was made of the rotating domain from the static domain. Likewise, it was used on the rotating domain and the propeller geometry. (a) Dimensions of the CFD domain in the x − z plane. (b) Dimensions of the CFD domain in the x − y plane. Figure 4.1: Dimensions of the CFD domain. The static domain is set to be stationary, and the rotating domain is set to rotate. An additional frame of reference is created beside the so-called lab reference frame. The rotating domain is set to have an angular velocity of 4038 rpm. Internal interfaces are created between the two domains and between the rotating domain and the propeller. A boolean imprint is also implemented between the two domains. Furthermore, a cylinder is created with a diameter of 1.4Dprop that stretched 2.81Dprop in each direction from the propeller center along the propeller’s central axis. 4.1.2 Mesh Generation A trimmed cell mesh with a base size of 11 mm is used. In the stationary region, where nothing of interest is expected, the mesh is made coarser by setting the target size to 110 mm. Likewise, the mesh is refined around the propeller by setting the target cell size to 2.2 mm. A volumetric control is implemented for the cylinder with a custom isentropic size of 3.3 mm. A prism layer mesh is implemented along the propeller surface. The total thickness 23 4. Method of the prism layers is set to 5.9 mm. The other prism layer settings are presented in Tab. 4.1. Along the boundaries and the interfaces, the prism layers are disabled. Table 4.1: Prism layer controls Control Value Number of prism layers 4 Gap-fill 45 % Minimum thickness percentage 0.01 Layer reduction percentage 0.0 Near wall thickness 1 mm Boundary march angle 50 deg The mesh is checked against the criteria in Sec. 2.9.3 to ensure good mesh quality. In Fig. 4.2, a zoomed-in view of the mesh close to the hub shows the prism layers growing from the wall. Figure 4.2: Mesh in the near wall region. 24 4. Method 4.1.3 Boundary Conditions The inlet of the computational domain is defined as a velocity inlet. The exit of the computational domain is defined as a pressure outlet. The pressure at the outlet is set to atmospheric pressure. All other walls in the domain are set as no-slip walls. A summary of the boundary conditions can be viewed in Tab. 4.2 Table 4.2: Boundary conditions Boundary Condition Inlet Velocity inlet Outlet Pressure outlet Outer Walls No-Slip Wall Propeller No-Slip Wall 4.1.4 Physics Models The physics models chosen for all simulations are three-dimensional, steady-state simulations, based on their suitability for the study’s objectives. The fluid used in the simulations is air, an ideal gas. The gas density is measured with a point probe just below the domain’s inlet. A segregated solver is used to solve the linearized equations obtained from the discretization of the governing equations. Table 4.3 gives an overview of the physics model. Table 4.3: Basic physics model Type Model used Space Three Dimensional Time Steady Material Gas Flow Segregated Equation of State Ideal Gas Energy Segregated Fluid Temperature Viscous Regime Turbulent 4.1.5 Stopping criteria To determine if the simulation converged, the criteria in Tab. 4.4 is used. The continuity (Eq. 2.17) and the three momentum equations (Eq. 2.19) are required to drop by at least three orders of magnitude. Additionally, the standard deviation of the efficiency is required to be less than 0.1 percentage units over at least 1000 iterations. 25 4. Method Table 4.4: Stopping Criteria Measure Type Value Efficiency Standard deviation 0.1 x-momentum Standard deviation 1e-3 y-momentum Standard deviation 1e-3 z-momentum Standard deviation 1e-3 Continuity Standard deviation 1e-3 4.2 Refining Results The following subsections underline the methods used to refine the results and ensure they portray the natural world as closely as possible. 4.2.1 Turbulence Models Four different turbulence models are tested and compared to the experimental case with the highest efficiency from UIUC, which can be studied closer in App. C [11]. The turbulence models compared are k − ϵ, k − ω, SST, and SA. The turbulence model that best matched the experimental data is then selected to be run for further simulations. 4.2.2 Mesh type The two more complex mesh types, the polyhedral mesher and the trimmed cell mesher, are tested to find which yielded the most accurate simulation results. 4.2.3 Wall-Y+ testing Three test simulations are run to find a suitable wall-y+ value with an average wall-y+ of below 1, between 11 and 30, and the last one over 300. Achieving these values required adjustments in the prism layer control settings. On the internet, many wall-y+ calculators are available to perform these calculations. Figure 4.3 demonstrates how the prism layers grow from the wall. 26 4. Method Figure 4.3: Prism layer growth. A cell growth ratio of 1.15 is used for all cases. For the first case, the height of the first prism layer (hwall) is 0.266 mm with four prism layers, giving a total height (htot) of 1.33 mm. In the second case, hwall is 5.92 mm, with three prism layers making htot 20.6 mm. In the last case, hwall is 104 mm, with two prism layers making htot 223 mm. The different settings for the three cases are summarized in Tab. 4.5. Table 4.5: Prism Layer Settings. Value Number of Layers First Cell Height [m] Total Layer Height [m] <1 4 2.66E-4 1.33E-3 [11, 30] 3 5.92E-3 2.06E-2 >300 2 1.04E-1 2.23E-1 4.2.4 Mesh Independence Study The mesh independence study is performed, with each refinement increasing the cell count by at least a factor of 2.2. This refinement is made in steps until the GCI index is deemed satisfactory for each mesh refinement. This is done to ensure that the simulation results are not dependent on any specific mesh. 4.3 Post-processing Different monitors are set up to retrieve results from the simulations and to be able to perform Eq. 2.12-2.16. Plots of the thrust- and torque coefficients and advance ratio are added to visibly verify convergence alongside the stopping criterion listed in Tab. 4.4. To ensure nothing unordinary happened with the simulation, the velocity 27 4. Method magnitude, velocity vectors, and pressure distribution in the surrounding fluid are studied with scalar and vector scenes. After the solution is gathered from the simulation, thrust- and torque coefficients and efficiency are compared to the experimental data. The percentage error is calculated with Eq. 3.1-3.3. 28 5 Results The following sections will cover the results of all the different simulations. 5.1 Mesh Type Testing The results from the tests of the polyhedral and trimmed cell meshes are presented in Tab. 5.1. The table also shows the error compared to the experimental data calculated with Eq. 3.1-3.3. The trimmed cell mesher performs better than the polyhedral mesher, and is deemed the best choice. Table 5.1: Results for the different mesh types. Mesh Type CT [−] eCT [%] CP [−] eCP [%] η[−] eη [%] Polyhedral 0.0373 20.81 0.0305 8.41 0.684 13.53 Trimmed Cell 0.0397 15.71 0.0315 5.41 0.705 10.87 Experimental (UIUC) 0.0471 - 0.0333 - 0.791 - 5.2 Turbulence Model Testing In Tab. 5.2, the results and the errors calculated with Eq. 3.1-3.3 for the different turbulence models are presented. The SA model had the lowest percentage error for the thrust, power, and efficiency coefficients. The SA model is selected for all further investigations. A pattern is identified that CT has the highest error and CP has the lowest error for all models. 29 5. Results Table 5.2: Turbulence model testing with J = 0.558893 and 4038 rpm. Model CT [−] eCT [%] CP [−] eCP [%] η[−] eη [%] k − ϵ 0.0374 20.56 0.0310 6.91 0.675 14.66 k − ω 0.0371 21.23 0.0310 6.91 0.670 15.30 SST 0.0372 21.02 0.0310 6.91 0.671 14.17 SA 0.0397 15.71 0.0315 5.41 0.705 10.87 Experimental (UIUC) 0.0471 - 0.0333 - 0.791 - 5.3 Wall-Y+ Testing The results for the simulations using different wall-y+ values are presented in Tab. 5.3. The simulations were run with the same advance ratio and angular velocity as in the previous tests, and the Spalart-Allmaras turbulence model and the trimmed cell mesher was used. The same table also shows the error in reference to the experi- mental data calculated with Eq. 3.1-3.3. The error has increased from the previous simulations. The numbers still indicate that a moderate wall-y+ gives a better result. Table 5.3: Results for different wall-y+ values. Value CT [−] eCT [%] CP [−] eCP [%] η[−] eη [%] < 1 0.0290 38.43 0.0242 27.33 0.669 15.42 [11, 30] 0.0336 28.66 0.0279 16.22 0.674 14.79 > 300 0.0297 36.94 0.0258 22.52 0.645 18.56 Experimental (UIUC) 0.0471 - 0.0333 - 0.791 - 5.4 Mesh Independence Study The number of cells and subsequent results from the three different meshes in the mesh independence study are shown below in Tab. 5.4. The critical factor used for the independence study was efficiency. Table 5.4: Results from the mesh independence study. Parameter Coarse Medium Fine GCIfine CT 0.0315 0.0396 0.0414 2.31% CP 0.0279 0.0317 0.0319 0.07% η 0.632 0.698 0.725 4.87% 30 5. Results 5.5 Presentation of Results The following sections will outline the results from the simulations. Contours and single-operating point graphs show the simulation done at the BEP for the medium grid size of 11 mm. Graphs comparing experimental and simulation results cover different operating points. 5.5.1 Experimental vs Simulated Results Figure. 5.1 the experimental efficiencies as well as the simulated ones can be seen. The simulated results follow the same shape as the experimental results, although they are a bit lower for low to mid advance ratios. The simulations edge closer to the experimental results and become higher than the experimental results for higher advance ratios. It should be noted that it is for low to mid advance ratios where the simulated results follow the shape of the experimental data. At higher advance ratios, it starts to deviate. Figure 5.1: Simulated efficiency compared to the experimental. Similar to the efficiency, the results for the thrust coefficient in Fig. 5.2a and the power coefficient in Fig. 5.2b, follow the general shape of the experimental results. Both coefficients are lower than the experimental results for lower advanced ratios and higher for the higher advanced ratios. 31 5. Results (a) Simulated CT compared to the exper- imental. (b) Simulated CP compared to the exper- imental. Figure 5.2: Simulated coefficients compared to the experimental. 5.5.2 Residuals Figure 5.3 shows the residuals for continuity, xyz-momentum, and energy for the simulation with a grid size base of 11 mm. All start to asymptote after 1000 itera- tions and deviate less. Note that the scale of the y-axis is in the order of 1e−3, which means the fluctuations are minimal. This convergence is similar for all simulations. Figure 5.3: Residuals for 11 mm base size, first 4000 iterations. 5.5.3 Wall-Y+ The wall-y+ distribution along the propeller in Fig. 5.4 shows that it is higher along the LE, reduces, and is lower at the TE when flow separation has occurred. On the top of the hub, the wall-y+ is lower. 32 5. Results Figure 5.4: Wall-y+ distribution along the propeller. 5.5.4 Velocities In Fig. 5.5 the velocity contour can be seen. From the contours, the speed increases further out on the blades, which follows from a higher velocity further out on the blade, as per Eq. 1.1. Around and below the hub is a low-speed wake area. The hub is more of a blunt object decelerating the flow. Figure 5.5: Velocity magnitude around the propeller. 33 5. Results The radial velocities, as depicted in Fig. 5.6a and 5.6b, exhibit a notable similarity upstream and downstream of the blade with the difference that downstream velocity shifts and accelerates inwards closer at the hub. It also shows the same tendencies for the lowest inlet velocity, highest efficiency, and highest inlet velocity. The radial velocity tends to have a zero velocity center close to the middle of the blade. Shifting closer to the hub for lower inlet velocities and further for greater inlet velocities. (a) Radial velocity upstream. (b) Radial velocity downstream. Figure 5.6: Radial velocity over the propeller. In Fig. 5.7a, the upstream circumferential velocity is close to zero, with the down- stream velocity in Fig. 5.7b peaking at a radius of around 40 mm. The circumfer- ential velocity is also lower for higher inlet velocity. (a) Circumferential velocity upstream. (b) Circumferential velocity downstream. Figure 5.7: Circumferential velocity over the propeller. Figure 5.8a and 5.8b shows the axial velocity downstream is generally greater than upstream, with a dip at the blade’s tip, indicating that the flow has been accelerated by the moving blades. For higher inlet velocities the change of axial velocity over the propeller is less. 34 5. Results (a) Axial velocity upstream. (b) Axial velocity downstream. Figure 5.8: Axial velocity over the propeller. 5.5.5 Pressure In the pressure contour in Fig. 5.9, areas of lower pressure can be seen above the propeller blades as the air gets pulled through the blades and pushed into higher- pressure areas below the blades. It can also be seen that a higher pressure area is present where the incoming air collides with the propeller hub, and as a result of this, there is also a lower pressure area present in the wake of the hub. Figure 5.9: Pressure distribution surrounding the propeller. In Fig. 5.10, the pressure distribution across the top and bottom of the blades can be seen. Per the pressure contour, the pressure is lower on the blade’s top surface 35 5. Results than on the bottom. It can also be seen that the pressure is higher at the LE than at the TE. Figure 5.10: Pressure distribution on the top and bottom of the propeller. The graphs in Fig. 5.11 plot the pressure distribution along the rotor blade against the blade radius. Figure 5.11a shows a mainly negative pressure upstream of the propeller, whereas Fig 5.11b shows it is positive downstream of the blade. It also shows that a higher inlet velocity results in a lower change in pressure between the upstream and downstream of the propeller. (a) Pressure upstream. (b) Pressure downstream. Figure 5.11: Pressure distribution over the propeller. 36 5. Results 5.5.6 Swirl Angle The swirl angle created by the propeller along the blades can be seen in Fig. 5.12. Figure 5.12a and 5.12b show that the flow is angled more by the blade for lower inlet velocities, due to the lower axial velocity seen in Fig. 5.8 (a) Swirl angle upstream. (b) Swirl angle downstream. Figure 5.12: Swirl angle over the propeller. 5.5.7 Vorticies Figure 5.13 shows the primary vortex sources. These are the tips of the propeller blades, along with the hub and the connection between the hub and blades. Small vortices can also be noticed following the TE, but they are insignificant in relation. Figure 5.13: Vortices generated by the propeller. 37 6 Analysis and Discussion In the final chapter of this thesis, the results from the simulations are analyzed. The method is discussed, dealing with issues of the chosen method and possible corrections. Lastly, a proposal for future work related to the topic is conferred. 6.1 Result Analysis The simulation model effectively mirrors the behavior of a real propeller. Compared to the experimental data, the error is at its lowest, around 10 percent for the effi- ciency. This is considered a good result, considering the time and resources available for the project. As anticipated, the discrepancy between numerical and experimental results dimin- ished with a more refined grid around the propeller. Notably, the wall-y+ functions demonstrated superior performance at a moderate value, suggesting that the prism layer settings should be tailored to this condition. Further strengthening is Fig. 5.1, which shows the simulation model following the experimental data for lower advance ratios with an error of around 10% for the lower advance ratios up to the BEP. This is likely due to the increase in inlet velocity leading to more unsteady features in the flow which is hard for the turbulence model to accurately simulate. Unsteady simulations would most likely improve these results. This also means the simulation model should be used for lower advance ratios in general. From this it can also be concluded, in accordance with Eq. 2.15, that the relation between CT and CP , shown in Fig. 5.2a and 5.2b respectively, follow each other best for lower advance ratios. It should also be noted that they both have a lesser decline compared to experimental values. All simulated results have a less irregular curve compared to experimental data due to the nature of the simulations solving equations, while experimental data measure real-world performance. The pressure distribution and velocity magnitude in the surrounding fluid are other expected results that strengthen the model’s integrity. The pressure increase is clearly shown in the contour in Fig. 5.9 and in the graphs found in Fig. 5.11. There is a low-speed area around and below the hub, which should be because the hub can 38 6. Analysis and Discussion not accelerate the fluid. This effect from the hub is also present in the radial velocity graphs in Fig. 5.6. An interesting notation is how the radial velocity centers zero around the middle of the blade. The velocity magnitude upstream and downstream seems to follow the pressure distribution. Upstream, as pressure decreases, the velocity also decreases until the turning point. Downstream, it is the opposite, as the pressure increases, the velocity decreases until the turning point. This could be an explanation for this observed zero-centering effect. However, the rapid pressure decrease close to the hub accelerates the flow inwards into the low-pressure area in the wake of the hub, explaining the negative radial velocity downstream. Common between all downstream velocities is that a higher inlet velocity leads to the distribution being more irregular along the blade than that of the lower inlet velocities, see Fig. 5.6b-5.8b and 5.12b. This is due to the propeller operating at the limit of its capabilities at that point. The circumferential velocity in Fig. 5.7 is almost zero before the fluid passes through the propeller. The blades then push the fluid, and the magnitude is highest around 40 mm from the hub. This is likely due to a combination of the blade’s absolute velocity in this section and the blade’s frontal area pushing on the fluid. The cir- cumferential velocity can be seen to be greater for lower inlet velocities. This could be a result of the fact that the fluid would be in the propeller’s path for a longer time and therefore, it would be pushed more in the rotational direction. The circumferential velocity can then be connected to the different swirl angles along the blade in Fig. 5.12 that have a similar profile to the velocity, see Eq. 2.3, i.e. the blade has changed the flow angle. This means the propeller has performed work on the fluid, which is also strengthened by the fact that the fluid has been accelerated. In summary, the higher the circumferential velocity at a particular section, the more work is being performed and the more the swirl angle changes; note that too large of a change can indicate stall. The largest difference in swirl angle over the propeller is where the geometrical pitch is greatest. Thus, this is where stall will occur first. For the axial velocity in Fig. 5.8, the downstream velocity in Fig. 5.8b is greater than the upstream velocity in Fig. 5.8a, meaning the propeller accelerates the fluid it pulls through. However, it accelerates the fluid less when the inlet velocity is increased due to the propeller not being able to push the fluid anymore. The axial fluid velocity is also higher around the propeller in general than that of the surrounding freestream. It also shows the velocity being highest at 70% of the propeller radius, a common point for manufacturers to display non-dimensional data like Reynold’s number. Figure 5.8b also gives indications on some of the vortices generated by the model, more precisely, the tip vortices. The dip in the downstream velocity at the tip of the blade is an indication of the phenomenon of fluid flowing from the bottom to the top of the propeller past the tip, causing vortices to be created. It also shows how this effect is lessened with a higher inlet velocity. The tip vortices are confirmed by Fig. 5.13 showing their existence, it also showed the prevalence of vortices at 39 6. Analysis and Discussion the hub, caused by abrupt flow separation because of the hub’s non-aerodynamical shape. Thus, mitigating these vortices, especially the tip vortices, improves aero- acoustic performance. However, this thesis only showed the existence of vortices and not their noise generation for this specific propeller. 6.2 Methodology Discussion The methodology of this thesis was first to find suitable settings that best reflected the experimental results. This was done by testing turbulence models, different meshes, and different wall-y+ functions and ensuring the mesh size was independent by doing a mesh independence study. These steps are essential in CFD modeling to provide satisfactory results. More computational resources would allow for further improvements of the simulation model, further ensuring it replicates the natural world as closely as possible. Other points regarding this work can be improved or examined further with more time. Modeling the propeller from APC data proved difficult, and some parameters they provided were ignored to make the creation of a propeller model feasible with the tools available. Certain information, like the transition from blade to hub, was unavailable, and these parts were therefore estimated from product pictures. It would be beneficial in the future to use a geometry directly from APC to remove any uncertainty regarding that part. It should also be noted that this thesis attempts to replicate experimental data without knowing the exact conditions under which these experiments occurred. For example, the temperature and viscosity of the surrounding air during the experi- ments were unknown, so generalized values were used. The CFD software also does not consider the propeller’s possible flex, which is most likely prevalent during high load. 6.3 Future Work This thesis has set up a simulation model for a dual-rotor propeller. For future studies, continuing upon this thesis, transient and unsteady simulations can be simulated to better capture the propeller’s aero-acoustic performance. If time serves, more simulations comparing different wall-y+ functions and more fine-tuned grid settings can be done to improve the simulation result further by reducing the error compared to the experimental result. The same simulation model could also be used in future work to compare the per- formance of the dual-rotor propeller to other types of propellers, like MIT’s toroidal propeller, to see if a toroidal propeller of similar diameter has comparable perform- ance. In addition, aeroacoustic simulations could be performed to compare the sound 40 6. Analysis and Discussion produced by the different types of propellers in relation to the thrust generated. 41 Bibliography [1] V. Patel, K. Nikam, S. Dikshit, M. Agarwalla and C. Zagade, ‘Propeller design and optimization for drones,’ Lecture Notes in Mechanical Engineer- ing, 2022. doi: 10.1007/978-981-16-7059-6_31. [Online]. Available: https: //www.scopus.com/record/display.uri?eid=2- s2.0- 85126230882& doi=10.1007%2f978-981-16-7059-6_31&origin=inward&tx (visited on 29/04/2024). [2] Hartzell Propeller. ‘A short history of the aircraft propeller.’ (2018), [On- line]. Available: https://hartzellprop.com/a-short-history-of-the- aircraft-propeller/ (visited on 27/04/2024). [3] N. Abdullah, P. Brapakaran, M. Sani, M. Marcel and E. Agung, ‘Investiga- tion on the dynamic properties of propeller structure with different number of blades,’ 2020. doi: 10.1088/1757-899X/807/1/012035. [Online]. Avail- able: https: // www.scopus .com /record/ display.uri ?eid =2- s2.0 - 85084073326&doi=10.1088%2f1757-899X%2f807%2f1%2f012035&origin= inward&tx (visited on 27/04/2024). [4] Science Photo Library, ’Light Aircraft Propeller’, [Photograph], n.d. [Online]. Available: https://quest-eb-com.eu1.proxy.openathens.net/images/ search/propeller/detail/132_3125815 (visited on 03/05/2024). [5] S. L. Dixon and C. Hall, Fluid mechanics and thermodynamics of turboma- chinery, 7th ed. Oxford, UK: Elsevier, 2014. [6] N. Hall. ‘Propeller propulsion,’ National Aeronautics and Space Administra- tion. (May 2021), [Online]. Available: https://www.grc.nasa.gov/www/k- 12/airplane/propeller.html (visited on 27/04/2024). [7] Massachusetts Institute of Technology, 11.7 performance of propellers, n.d. [Online]. Available: https : / / web . mit . edu / 16 . unified / www / SPRING / propulsion/notes/node86.html (visited on 06/05/2024). 42 https://doi.org/10.1007/978-981-16-7059-6_31 https://www.scopus.com/record/display.uri?eid=2-s2.0-85126230882&doi=10.1007%2f978-981-16-7059-6_31&origin=inward&tx https://www.scopus.com/record/display.uri?eid=2-s2.0-85126230882&doi=10.1007%2f978-981-16-7059-6_31&origin=inward&tx https://www.scopus.com/record/display.uri?eid=2-s2.0-85126230882&doi=10.1007%2f978-981-16-7059-6_31&origin=inward&tx https://hartzellprop.com/a-short-history-of-the-aircraft-propeller/ https://hartzellprop.com/a-short-history-of-the-aircraft-propeller/ https://doi.org/10.1088/1757-899X/807/1/012035 https://www.scopus.com/record/display.uri?eid=2-s2.0-85084073326&doi=10.1088%2f1757-899X%2f807%2f1%2f012035&origin=inward&tx https://www.scopus.com/record/display.uri?eid=2-s2.0-85084073326&doi=10.1088%2f1757-899X%2f807%2f1%2f012035&origin=inward&tx https://www.scopus.com/record/display.uri?eid=2-s2.0-85084073326&doi=10.1088%2f1757-899X%2f807%2f1%2f012035&origin=inward&tx https://quest-eb-com.eu1.proxy.openathens.net/images/search/propeller/detail/132_3125815 https://quest-eb-com.eu1.proxy.openathens.net/images/search/propeller/detail/132_3125815 https://www.grc.nasa.gov/www/k-12/airplane/propeller.html https://www.grc.nasa.gov/www/k-12/airplane/propeller.html https://web.mit.edu/16.unified/www/SPRING/propulsion/notes/node86.html https://web.mit.edu/16.unified/www/SPRING/propulsion/notes/node86.html Bibliography [8] T. Sebastian, ’Toroidal Propeller’, Powerpoint, Massachusetts Insitute of Tech- nology, 2023. [9] N. Hall. ‘Propeller thrust,’ National Aeronautics and Space Administration. (May 2021), [Online]. Available: https : / / www . grc . nasa . gov / www / k - 12/airplane/propth.html (visited on 07/02/2024). [10] S. Gudmundsson, ‘General aviation aircraft design,’ in 2nd ed., Butterworth- Heinemann, 2022, isbn: 978-0-12-818465-3. [Online]. Available: https://www. sciencedirect.com/book/9780128184653/general-aviation-aircraft- design (visited on 14/02/2024). [11] University of Illinois at Urbana-Champaign, ’UIUC Propeller Data Site’, 2022. [Online]. Available: https://m-selig.ae.illinois.edu/props/propDB. html (visited on 29/02/2024). [12] A. Posa, ‘Dependence of tip and hub vortices shed by a propeller with wing- lets on its load conditions,’ Physics of Fluids, vol. 34, no. 10, 2022. [Online]. Available: https://www.scopus.com/record/display.uri?eid=2-s2.0- 85139944891&doi=10.1063%2f5.0113480&origin=inward&tx (visited on 06/05/2024). [13] H.-D. Yao, Z. Huang, L. Davidson, J. Niu and Z.-W. Chen, Aerospace, vol. 9, no. 12, 2022. [Online]. Available: https://www.scopus.com/record/display. uri ? eid = 2 - s2 . 0 - 85144897085 & doi = 10 . 3390 % 2faerospace9120825 & origin=inward&tx (visited on 06/05/2024). [14] National Aeronautics and Space Administration. ‘NACA Airfoils.’ (2017), [On- line]. Available: https://www.nasa.gov/image-article/naca-airfoils/ (visited on 15/04/2024). [15] SiemensDigital Industries Software, Simcenter StarCCM+, version 2310, Siem- ens, Munich, Germany, 2023. [16] H. H. Hu, ‘Chapter 10 - computational fluid dynamics,’ in Fluid Mechan- ics, P. K. Kundu, I. M. Cohen and D. R. Dowling, Eds., 5th ed., Boston (MA), USA: Academic Press, 2012, pp. 421–472, isbn: 978-0-12-382100-3. [On- line]. Available: https://www.sciencedirect.com/science/article/pii/ B9780123821003100101 (visited on 21/03/2024). [17] H. K. Versteeg and W. Malalasekera, An introduction to computational fluid dynamics: the finite volume method, 2nd ed. Essex, England: Pearson educa- 43 https://www.grc.nasa.gov/www/k-12/airplane/propth.html https://www.grc.nasa.gov/www/k-12/airplane/propth.html https://www.sciencedirect.com/book/9780128184653/general-aviation-aircraft-design https://www.sciencedirect.com/book/9780128184653/general-aviation-aircraft-design https://www.sciencedirect.com/book/9780128184653/general-aviation-aircraft-design https://m-selig.ae.illinois.edu/props/propDB.html https://m-selig.ae.illinois.edu/props/propDB.html https://www.scopus.com/record/display.uri?eid=2-s2.0-85139944891&doi=10.1063%2f5.0113480&origin=inward&tx https://www.scopus.com/record/display.uri?eid=2-s2.0-85139944891&doi=10.1063%2f5.0113480&origin=inward&tx https://www.scopus.com/record/display.uri?eid=2-s2.0-85144897085&doi=10.3390%2faerospace9120825&origin=inward&tx https://www.scopus.com/record/display.uri?eid=2-s2.0-85144897085&doi=10.3390%2faerospace9120825&origin=inward&tx https://www.scopus.com/record/display.uri?eid=2-s2.0-85144897085&doi=10.3390%2faerospace9120825&origin=inward&tx https://www.nasa.gov/image-article/naca-airfoils/ https://www.sciencedirect.com/science/article/pii/B9780123821003100101 https://www.sciencedirect.com/science/article/pii/B9780123821003100101 Bibliography tion, 2007. [Online]. Available: http://200.17.228.88/disciplinas/TM701/ Versteeg_Malalasekera_2ed.pdf (visited on 21/03/2024). [18] F. M White and H. Xue, Fluid Mechanics Ninth Edition, 9th ed. New York (NY), USA: McGraw Hill LLC, 2021, isbn: 9781260575545. [19] A. G. Gardi, ’Moving reference frame and Arbitrary Lagrangian Eulerian ap- proaches for the study of moving domains in Typhon’, 2010. [Online]. Avail- able: https://www.politesi.polimi.it/bitstream/10589/26182/3/ tesiGARDI_corretta.pdf (visited on 29/02/2024). [20] SiemensDigital Industries Software, Simcenter StarCCM+ User Guide, ver- sion 2310, Siemens, Munich, Germany, 2023. [21] L. Davidson. ‘An introduction to turbulence models,’ Chalmers University of Technology. (2022), [Online]. Available: https://www.tfd.chalmers.se/ ~lada/postscript_files/kompendium_turb.pdf (visited on 04/04/2024). [22] CFD Online, Law of the wall, n.d. [Online]. Available: https://www.cfd- online.com/Wiki/Law_of_the_wall (visited on 27/05/2024). [23] J. Tu, G.-H. Yeoh and C. Liu, ‘Chapter 4 - CFD Mesh Generation: A Prac- tical Guideline,’ in Computational Fluid Dynamics (Third Edition), 3rd ed., Oxford, UK: Elsevier, 2018, pp. 125–154, isbn: 978-0-08-101127-0. [Online]. Available: https : / / www . sciencedirect . com / science / article / pii / B9780081011270000040 (visited on 16/04/2024). [24] Siemens. ‘Understanding the subtle difference between the prism layer mesher and the advancing layer mesher.’ (2022), [Online]. Available: https://community. sw.siemens.com/s/article/Understanding-the-subtle-difference- between-the-Prism-Layer-mesher-and-the-Advancing-Layer-Mesher (visited on 15/04/2024). [25] J. W. Slater. ‘Examining Spatial (Grid) Convergence,’ NPARC Alliance CFD Verification and Validation Web Site. (Feb. 2021), [Online]. Available: https: //www.grc.nasa.gov/www/wind/valid/tutorial/spatconv.html (visited on 23/04/2024). [26] APC Propellers, ’14x8’, n.d. [Online]. Available: https://www.apcprop.com/ product/14x8/ (visited on 29/02/2024). 44 http://200.17.228.88/disciplinas/TM701/Versteeg_Malalasekera_2ed.pdf http://200.17.228.88/disciplinas/TM701/Versteeg_Malalasekera_2ed.pdf https://www.politesi.polimi.it/bitstream/10589/26182/3/tesiGARDI_corretta.pdf https://www.politesi.polimi.it/bitstream/10589/26182/3/tesiGARDI_corretta.pdf https://www.tfd.chalmers.se/~lada/postscript_files/kompendium_turb.pdf https://www.tfd.chalmers.se/~lada/postscript_files/kompendium_turb.pdf https://www.cfd-online.com/Wiki/Law_of_the_wall https://www.cfd-online.com/Wiki/Law_of_the_wall https://www.sciencedirect.com/science/article/pii/B9780081011270000040 https://www.sciencedirect.com/science/article/pii/B9780081011270000040 https://community.sw.siemens.com/s/article/Understanding-the-subtle-difference-between-the-Prism-Layer-mesher-and-the-Advancing-Layer-Mesher https://community.sw.siemens.com/s/article/Understanding-the-subtle-difference-between-the-Prism-Layer-mesher-and-the-Advancing-Layer-Mesher https://community.sw.siemens.com/s/article/Understanding-the-subtle-difference-between-the-Prism-Layer-mesher-and-the-Advancing-Layer-Mesher https://www.grc.nasa.gov/www/wind/valid/tutorial/spatconv.html https://www.grc.nasa.gov/www/wind/valid/tutorial/spatconv.html https://www.apcprop.com/product/14x8/ https://www.apcprop.com/product/14x8/ Bibliography [27] H. Ahmed and P. Rajendran, ‘3D CFD Simulation and Experimental Valida- tion of Small APC Slow Flyer Propeller Blade,’ Aerospace, vol. 4, no. 1, Feb. 2017. doi: 10.3390/aerospace4010010. [28] 14x8, Woodland (CA), USA: APC Propellers, 2022. [Online]. Available: https: //www.apcprop.com/technical-information/file-downloads/ (visited on 29/02/2024). 45 https://doi.org/10.3390/aerospace4010010 https://www.apcprop.com/technical-information/file-downloads/ https://www.apcprop.com/technical-information/file-downloads/ A Geometry Data Geometry data for the 14x8 inch propeller. Adapted from [28]. Station [In] Chord [In] Sweep [In] Thickness ratio Twist [Deg] 1.5085 0.9912 0.6106 0.2990 37.5348 1.5988 0.9867 0.6156 0.2892 36.7624 1.6890 0.9835 0.6211 0.2796 35.9843 1.7793 0.9817 0.6269 0.2704 35.1993 1.8695 0.9812 0.6331 0.2615 34.4068 1.9598 0.9819 0.6395 0.2530 33.6060 2.0551 0.9838 0.6466 0.2444 32.7502 2.2290 0.9901 0.6600 0.2296 31.1645 2.4081 0.9998 0.6742 0.2156 29.4965 2.5871 1.0122 0.6885 0.2030 27.8182 2.7662 1.0265 0.7026 0.1917 26.2824 2.9452 1.0421 0.7160 0.1818 24.8974 3.1243 1.0582 0.7286 0.1731 23.6456 3.3033 1.0742 0.7398 0.1658 22.5120 3.4823 1.0894 0.7495 0.1598 21.4838 3.6614 1.1030 0.7572 0.1552 20.5502 3.8404 1.1144 0.7627 0.1518 19.7019 4.0195 1.1228 0.7655 0.1498 18.9307 4.1985 1.1276 0.7654 0.1483 18.2223 4.3776 1.1280 0.7620 0.1468 17.5654 4.5566 1.1232 0.7550 0.1453 16.9544 4.7356 1.1127 0.7440 0.1438 16.3845 4.9147 1.0956 0.7288 0.1423 15.8513 5.0937 1.0713 0.7090 0.1409 15.3510 5.2728 1.0391 0.6842 0.1394 14.8803 5.4518 0.9981 0.6541 0.1379 14.4361 5.6309 0.9478 0.6185 0.1364 14.0153 5.8099 0.8874 0.5769 0.1349 13.6152 5.9889 0.8162 0.5291 0.1334 13.2327 6.1680 0.7334 0.4747 0.1319 12.8165 6.3470 0.6384 0.4133 0.1304 12.3400 6.5261 0.5310 0.3448 0.1289 12.0134 6.7050 0.4101 0.2687 0.1274 11.9298 6.8801 0.2764 0.1853 0.1260 12.0715 7.0000 0.0001 -0.0482 0.1250 8.7923 I B Experimental propeller graphs Figure B.1: Experimental efficiency against advance ratio. Adapted from UIUC [11]. Figure B.2: Experimental thrust and power coefficient against advance ratio. Ad- apted from UIUC [11]. II C Experimental propeller data UIUC wind tunnel data for 4038 rpm. Adapted from [11]. J CT CP η 0.434111 0.069554 0.040846 0.739233 0.459433 0.065386 0.039752 0.755707 0.484630 0.060503 0.038182 0.767935 0.509538 0.055928 0.036611 0.778391 0.535009 0.050954 0.034757 0.784329 0.558893 0.047136 0.033285 0.791466 0.581138 0.042215 0.031044 0.790246 0.605672 0.037089 0.028658 0.783863 0.630828 0.032068 0.025913 0.780670 0.657955 0.025259 0.022517 0.738080 0.682181 0.018815 0.019177 0.669309 0.709614 0.012196 0.015644 0.553236 0.737836 0.005350 0.011812 0.334205 0.758895 -0.001539 0.007785 -0.150048 0.758866 -0.001555 0.007770 -0.151884 III DEPARTMENT OF MECHANICS AND MARITIME SCIENCES CHALMERS UNIVERSITY OF TECHNOLOGY Gothenburg, Sweden www.chalmers.se www.chalmers.se Introduction Background Purpose and Objectives Limitations Theory Propeller physics Coordinate system and Velocity triangles Blade Element Theory Performance parameters Noise generation Airfoil geometry Computational Fluid Dynamics Governing equations Moving Reference Frame Turbulence Models Wall-Y+ Functions Boundary conditions Meshing Mesh Refinement Prism Layers Mesh Quality Mesh Independence Study APC 14x8 Propeller Experimental Data Geometry Experimental Validation Method Simulation Setup CFD Domain Mesh Generation Boundary Conditions Physics Models Stopping criteria Refining Results Turbulence Models Mesh type Wall-Y+ testing Mesh Independence Study Post-processing Results Mesh Type Testing Turbulence Model Testing Wall-Y+ Testing Mesh Independence Study Presentation of Results Experimental vs Simulated Results Residuals Wall-Y+ Velocities Pressure Swirl Angle Vorticies Analysis and Discussion Result Analysis Methodology Discussion Future Work Geometry Data Experimental propeller graphs Experimental propeller data