variable. La pale de HAWT a un rapport 'aspect trés important, par conséquent, l'écoulement transverse n'a pas assez de distance pour se développer du bord d'attaque au bord de fuite. Au contraire, l'écoulement transverse est capable de croître assez pour conduire à la transition sur l'hélice, où le rapport d'aspect est petit.

Nous soutenons qu'il est nécessaire de travailler dans le référentiel rotatif pour évaluer l'effet de la rotation. Dans le cas du référentiel rotatif, la force centrifuge a une position donnée est une constante, tandis que la force de Coriolis dépend de la vitesse locale. L'amplitude de la vitesse en envergure (cross-flow) est étroitement liée à la force centrifuge et à la force de Coriolis. Par exemple, l'écoulement transverse est généralement le plus forte autour de la séparation de couche limite, où la force de Coriolis change de direction et agit dans le même sens que la force centrifuge.

Three years ago, when I landed in France, the excitement of living and studying in a totally different country makes it less unbearable to leave hometown. Now, with the PhD almost finished, the idea of leaving this beautiful country makes me more and more heavy-hearted as the remaining time elapses.

First of all, I would like to express my thanks to the supervisor of this thesis, HDR Antoine DUCOIN. I still remember vividly that he stood in front of the airport to pick me up when I went out from it. During the past three years, he provides me with excellent cluster computers. He is also quite tolerant to my explorations and mistakes. Almost all my choices were respected. On the other hand, he is quite critical to the works. It is always beneficial to discuss with him. This work won't be possible without his supervision. I am also very much indebted to HDR Caroline BRAUD, who is the co-supervisor of this PhD. Thanks to her, the experimental data on the wind turbine is available. She gives me a lot of advices on the wind turbine simulation. She also generously provides me with many other helps regarding administrative stuff and paper writing etc. I would also like to express my gratitude to the CSI members of the thesis, PU Jacques-André ASTOLFI and MCF Mostafa Safdari SHADLOO. They have been following this works since the very beginning and have given me quite many valuable advices.

The LHEEA lab has a really agreeable working environment. Sijo George, who is my office mate, helps me a lot regarding to the code we use. Gaël Clodic, my another office mate, has enriched my knowledge on the France. The IT department solves all my problems professionally and promptly. The daily talks with MCF Zhe LI and other PhD students are quite enjoyable.

The financial support of this PhD from CSC (China Scholarship Council) is greatly appreciated. I am also very grateful to IDRIS (Institut du Développement et des Ressources en Informatique Scientifique) for the computational hours.

Last but not least, I want to express my thanks to my families. They have been really i List of Figures 1.8 3D boundary layer flow structure. From [START_REF] William S Saric | Stability and transition of three-dimensional boundary layers[END_REF]. y t is the local wall-normal direction. x t is the local streamwise direction and z t is perpendicular to x t and y t . Outside boundary layer, the velocity is totally in streamwise direction. However, there is a velocity component in z t direction inside the boundary layer. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.9 Cross-flow vortices and flow transition (iso-surface of λ2). The inviscid flow is roughly in x direction. From [START_REF] Pc Dörr | Transition control in a three-dimensional boundary layer by direct attenuation of nonlinear crossflow vortices using plasma actuators[END_REF] . . . . . . . . . . . . . . . . . 1.10 (a) The mushroom-like structure in Görtler instability. From [START_REF] Ducoin | Direct numerical simulation of flow instabilities over savonius style wind turbine blades[END_REF]. (b)

The half mushroom structure in cross-flow instability. From [START_REF] Serpieri | Three-dimensional organisation of primary and secondary crossflow instability[END_REF]. . . . . iii LIST OF FIGURES 1.11 Illustration of streamlines in inviscid region and in boundary layer on a rotating disk. Because the wall-normal velocity, which is parallel to the angular velocity vector, does not affect the direction and magnitude of the forces, therefore it is neglected in plotting these streamlines . . . . 12 1.12 The arise of the cross-flow velocity in the rotating disk boundary layer. . [START_REF] Reynolds | An experimental investigation of the circumstances which determine whether the motion of water shall be direct or sinuous, and of the law of resistance in parallel channels[END_REF] 1.13 Self-similar solution of the von Kármán boundary layer. From [START_REF] Lingwood | Absolute instability of the boundary layer on a rotating disk[END_REF] . . . [START_REF] Reynolds | An experimental investigation of the circumstances which determine whether the motion of water shall be direct or sinuous, and of the law of resistance in parallel channels[END_REF] 1. [START_REF] Reynolds | An experimental investigation of the circumstances which determine whether the motion of water shall be direct or sinuous, and of the law of resistance in parallel channels[END_REF] The flow inside a teacup, which is an example of Bödewadt layer. Here the cross-flow velocity is towards the rotating axis. From [START_REF] Peter Alan | Turbulence in rotating, stratified and electrically conducting fluids[END_REF] . . . . . Time is discretized using Forward Euler method. Spatial derivative is approximated using first order upwind finite difference. . . . . . . . . . 4.7 Experiment measurements of K t and K q for different advance ratios J (from Ref. [START_REF] Kuiper | Scale effects on propeller cavitation inception[END_REF]). The symbols are from current simulations. The order 6 (red) and order 8 (black) Sur la pale d'hélice marine, l'écoulement de la couche limite est entièrement tridimensionnel (3D) à cause des effets combinés de la rotation et de la géométrie de la pale. L'instabilité et la transition des 'cross-flow' sont clairement observées. La forme des tourbillons est en bon accord avec la prédiction de la théorie de la stabilité linéaire (LST). Bien qu'il ait été longtemps supposé que la "cross-flow" transition devrait être importante sur les hélices, il s'agit de première observation directe de tels phénomènes à notre connaissance. Parceque l'hélice n'a pas de symétrie de rotation infinie, notre résultat suggère que la couche limite sur les hélices marines est instable par convection.

Cet aspect est différent de l'écoulement de von Kármán, qui est absolument instable. The difference in boundary layer flows and therefore transitions between the HAWT and the marine propeller is likely caused by their shapes. The HAWT blade has a very large aspect ratio, consequently, the cross-flow does not have enough distance to develop from the leading edge to the trailing edge. On the contrary, the cross-flow xiii ABSTRACT velocity is able to grow large enough to lead flow transition on the propeller, where the aspect ratio is small.

We also argue that it is necessary to work in the rotating reference frame in order to evaluate the effect of rotations. In the rotating reference frame, the centrifugal force at one position is a constant, while the Coriolis force depends on the local velocity. The magnitude of the spanwise (cross-flow) velocity is closely related to the relative strength of centrifugal and Coriolis forces. For example, cross-flow is usually the largest around separation bubbles, where Coriolis force changes its direction and acts in the same direction as centrifugal force.

Contents

Acknowledgements i

List of Figures iii

Abstract xi Another broad topic in the LHEEA lab is the naval engineering. The research areas include cavitation of marine propellers [START_REF] Huang | Physical and numerical investigation of cavitating flows around a pitching hydrofoil[END_REF], development of composite hydrofoils [START_REF] Pernod | Experimental and numerical investigation of the fluid-structure interaction on a flexible composite hydrofoil under viscous flows[END_REF],

and numerical methods for hydrodynamic [START_REF] Marrone | δ-sph model for simulating violent impact flows[END_REF][START_REF] Guilmineau | Numerical simulations of the cavitating and non-cavitating flow around the postdam propeller test case[END_REF] etc.

HAWTs harvest energy from the wind. The torque that drives motors is from the lift force on airfoil sections of the turbine blades [START_REF] Martin | Aerodynamics of wind turbines[END_REF]. Whereas marine propellers exert energy to the surrounding fluid. The propulsion is also from the lift force on sections [START_REF] Carlton | Marine propellers and propulsion[END_REF]. Although the HAWTs and marine propellers have quite different geometries and there is another apparent similarity, i.e. they both rotates.

Motivation of the present work

The fact that both HAWTs and marine propellers rotate during operations makes them bear a resemblance to one of the most classic problems in fluid dynamic, which is the laminar-turbulent transition of the rotating disk boundary layer.

The most notable feature of the rotating disk boundary layer is the 3D effect on the boundary layer due to the rotation. This gives rise to a totally different transition mechanism from 2D boundary layers, namely the cross-flow transition [START_REF] Lingwood | Absolute instability of the boundary layer on a rotating disk[END_REF]. However, most studies on HAWTs and marine propellers' boundary layer transition are based on 2D airfoils. To get better understanding of the process on HAWTs and marine propellers, it is necessary to quantify the effect of rotations.

In the rest of this chapter, we will firstly give a brief description of the aerodynamics of 2D airfoils. Then the laminar-turbulent transition of boundary layer is introduced, followed by the description of possible transition scenarios related to the present study.

We then revisit the rotating disk flow and its transition. The rise of cross-flow inside boundary layer is explained. We re-state the similarity between the rotating disk flow and HWATs and marine propellers. At last, the state-of-the-art about the boundary layer transition on rotating blades is reviewed.

Aerodynamics of 2-Dimensional Airfoil

The airfoil is the basis of many practical fluid machinery including HAWTs, propellers, and airplane wings et al. It is a lift generating structure. Figure 1.2 gives an illustration of the flow around a typical 2D airfoil. When the flow passes around the airfoil, there usually exists a pressure difference between the two sides of the airfoil (the side with large pressure is usually referred to as the pressure side, and the one with smaller pressure is referred to as the suction side). This pressure difference is the main source of the force on the airfoil. AoA measures the angle between the incoming flow direction and the direction of airfoil's chord, which is the line connecting the leading edge and the trailing edge of the airfoil. When AoA is small, lift coefficient C l roughly increases linearly with

AoA as predicted by lifting-line theory. However, if AoA further increases, the C l would decrease rapidly after C l having reached maximum. At the same time, the drag increases as well. This phenomenon is called 'stall'.

Stall is generally related to massive flow separations on the airfoil [START_REF] Tani | Low-speed flows involving bubble separations[END_REF][START_REF] Rodríguez | Direct numerical simulation of a naca0012 in full stall[END_REF]. Pressure gradients always present along the airfoil chord. The adversary pressure gradient (where the downstream pressure is larger than upstream pressure) is inevitable due to the nature of airfoils. The adversary pressure gradient would decelerate the flow. When the inertial force is not able to overcome the adversary pressure gradient, the attached boundary layer flow could separate for the wall surface.

When AoA is small, the flow separation does not exist or is very weak. The flow usually reattaches after separations, which is usually close to the trailing edge. As AoA and the adversary pressure gradient increase, the separation points moves upstream on the airfoil. Under stall conditions, massive flow separations appear close to leading edge [START_REF] Rodríguez | Direct numerical simulation of a naca0012 in full stall[END_REF] and the flow does not reattach. The effective shape of the airfoil is profoundly changed and its aerodynamic performance deteriorates.

The boundary layer and its transition 1.3.1 The boundary layer

On the airfoil surface, the flow velocity decelerates to zero due to the non-slip boundary condition. Very close the airfoil wall, there exists a thin layer of fluid, through which the outside inviscid flow velocity decreases to zero on the wall. This fluid layer is referred to as 'boundary layer' and is firstly discovered by Prandtl [START_REF] Schlichting | Boundary-layer theory[END_REF]. The importance of the boundary layer could not be emphasized too much.

δ x ∝ √ x.
When there is a pressure gradient along x, the velocity profile is still self- similar and can be described by the Falkner-Skan solution [START_REF] Schlichting | Boundary-layer theory[END_REF]. Inside the boundary layer, the large velocity gradient in wall-normal direction makes the viscous force nonnegligible.

The boundary layer around airfoil is usually not self-similar. However, many concepts from Blasius solution can be applied to airfoil boundary layer. For example, they both develop from the leading edge toward downstream. The same technique is often used in the derivation and solving of the boundary layer equation for airfoil. What's more, the transition on airfoils is very much similar to flat plate when the boundary layer is attached.

The laminar-turbulent transition

The famous Reynolds' experiment [START_REF] Reynolds | An experimental investigation of the circumstances which determine whether the motion of water shall be direct or sinuous, and of the law of resistance in parallel channels[END_REF] on the pipe flow shows that there are two different flow states depends on a parameter which is now called Reynolds number (Re).

When Re is small, the dye streaks are aligned with each other and the flow is laminar.

When Re is large, the flow become chaotic and turbulent.

The laminar-turbulent transition has many implications on the aerodynamic performances of the bodies. It affects flow separation, increases drags, and introduces fluctuating loads etc. The transition process of boundary layer flow is quite complex and depends on many factors. Figure 1.4 summaries the main routes through which a laminar boundary layer could become turbulent. In short, if the flow is unstable to external disturbance, the disturbance would increase its magnitude and eventually leads to flow transition.

In this thesis, only the natural transition is considered (Path A in figure 1.4). Although natural transitions are more common when the environmental disturbances are small, a good understanding of them is important to predict the boundary layer's responses to large disturbances. In the following subsections, several natural transition modes pertain to airfoil boundary layers are described briefly.

T-S wave instability

The Blasius boundary layer is subjected to the Tollmien-Schlichting (T-S) instability, which is firstly studied theoretically by Tollmien and Schlichting around the late 1920s.

It was confirmed by the experiments of Schubauer and Skramstad [START_REF] Galen | Laminar boundary-layer oscillations and transition on[END_REF][START_REF] Philip | Hydrodynamic stability[END_REF]. When the 3D disturbance becomes large, the flow usually exhibits Λ structures. The 3D flow structure itself is quite unstable and would quickly breakdown to smaller structures. The transition process finishes and the flow becomes fully turbulent.

Flow separation induced transition

The Laminar Separation Bubble (LSB) is commonly seen on the suction side of airfoils. A recirculation zone appears near the wall with part of the fluid moves reversely. After the boundary layer transition to turbulent, the streamlines reattaches.

Depending on many factors such as the geometry of airfoil, Re, AoA, and inflow conditions [START_REF] Muti | Low-reynolds-number separation on an airfoil[END_REF][START_REF] Mm O'meara | Laminar separation bubble characteristics on an airfoil at low reynolds numbers[END_REF], LSBs have different shapes. However, they are generally quite unstable when Re is large. The wall-normal velocity profile of LSB usually has an inflectional point, making it unstable according to inviscid stability theory.

When the separation bubble is small, the unstable T-S wave, which appears in the upstream attached boundary layer, can go through the LSB. In Ref. [START_REF] Rist | On instabilities and transition in laminar separation bubbles[END_REF], the growth rate of the T-S wave in separation region obtained from DNS (Direct Numerical Simulation) agrees quite well that from LST. Inflectional instability can also happen due to the upstream disturbance [START_REF] Sourabh | On the origin of the inflectional instability of a laminar separation bubble[END_REF]. However, this instability is convectively unstable in nature.

When the separation is strong and the reverse flow is large, inviscid instability mech-anism would dominate the transition process. The boundary layer velocity profile with large reverse flow is similar to a mixing layer, where rapid flow transition can be observed due to Kelvin-Helmholtz (K-H) instability. Similar transition processes to K-H instability are often observed around the LSB. The primary stability results in 2D shedding vortices. As the shedding vortices being convected to downstream, secondary instability occurs [START_REF] Marxen | Vortex formation and vortex breakup in a laminar separation bubble[END_REF]. The flow becomes rapidly turbulent. Figure 1.7 shows the structures in a LSB induced transition, where the rollers are clearly visible.

Similar to the K-H instability, LSB is absolute unstable when the reverse flow is large. As a result, the LBS always exhibit unsteadiness. Pauley et al. [START_REF] Pauley | The structure of twodimensional separation[END_REF] use the Stroulhal number based on the boundary-layer momentum thickness at separation and the local free-stream velocity to quantify the shedding frequency, which is in agreement with the prediction of linear inviscid stability theory.

It should be noted that there is no universal criterion to determine whether the LSB is large or small. Generally speaking the flow is absolute unstable only when the reverse flow is as large as 12% -20% of the boundary layer edge velocity [START_REF] Alam | Direct numerical simulation of 'short'laminar separation bubbles with turbulent reattachment[END_REF][START_REF] Rodríguez | The two classes of primary modal instability in laminar separation bubbles[END_REF]. The interactions between the viscous instability and the inflectional instability are also reported in the literatures (for example [START_REF] Joshua | Interaction of viscous and inviscid instability modes in separation-bubble transition[END_REF]).

Cross-flow transition

So far, we have been talking about 2D boundary layer. i.e. the flow is homogeneous in one direction (usually spanwise) and the fluid does not move along that direction. However, 3D boundary layer, where there is no homogeneous direction, is quite common in practical engineering flows. Because cross-flow is zero on the wall and outside the boundary layer, its profile usually has an inflectional point, which makes it unstable according to the Rayleigh theorem. The cross-flow transition can be observed in many typically 3D flows, e.g. the Figure 1.8: 3D boundary layer flow structure. From [START_REF] William S Saric | Stability and transition of three-dimensional boundary layers[END_REF]. y t is the local wall-normal direction. x t is the local streamwise direction and z t is perpendicular to x t and y t . Outside boundary layer, the velocity is totally in streamwise direction. However, there is a velocity component in z t direction inside the boundary layer.

Figure 1.9: Cross-flow vortices and flow transition (iso-surface of λ2). The inviscid flow is roughly in x direction. From [START_REF] Pc Dörr | Transition control in a three-dimensional boundary layer by direct attenuation of nonlinear crossflow vortices using plasma actuators[END_REF] .

rotating disk [START_REF] Rj Lingwood | On the impulse response for swept boundary-layer flows[END_REF], swept wings [START_REF] William S Saric | Stability and transition of three-dimensional boundary layers[END_REF], and the yawed cone [START_REF] Li | High-frequency instabilities of stationary crossflow vortices in a hypersonic boundary layer[END_REF] etc. The saturated cross-flow vortices are subjected to the secondary instability. In this stage, the vortices are unstable with respect to high-frequency oscillations. With the growth of secondary instability, the flow quickly becomes turbulent.

Others

Many other unstable modes also exist depending on the flow configurations. For example, along concave walls, Görtler instability often appears in the boundary layer. The flow along the leading edge of a swept wing is subjected to so-called 'leading edge contamination'. Roughness elements on the wall can also lead to boundary layer transition.

Transition on the 2D airfoil

On 2D airfoil without swept or where the blade spanwise variation is small, the natural transition of the boundary layer can be either triggered by T-S wave or flow separation.

When the Reynolds number is large, the boundary layer close to the leading edge would be unstable to T-S wave [START_REF] Würz | Detuned resonances of tollmien-schlichting waves in an airfoil boundary layer: Experiment, theory, and direct numerical simulation[END_REF][START_REF] Wu Mengmeng | A transition prediction method for flow over airfoils based on high-order dynamic mode decomposition[END_REF]. The T-S wave on the airfoil is very much similar to the flat plate case describe above. However, there is usually a pressure gradient along the streamwise direction on the airfoil. Similar to the Falkner-Skan boundary layer [START_REF] Ar Wazzan | Spatial and temporal stability charts for the falkner-skan boundary-layer profiles[END_REF], the adversary pressure gradient has a dis-stabilizing effect whereas the favorable pressure gradient makes the boundary layer be more stable.

When the chord Reynolds number is small (for example, the unmanned aerial vehicles), the laminar boundary layer can maintain a long distance on the airfoil. In that case, the boundary layer transition is usually triggered by the LSB, which appears when 1.5 The rotating disk flow

Cross-flow velocity

The laminar rotating disk flow has an exact self-similar solution of the Navier-Stokes (NS) equations, which is first discovered by von Kármán [START_REF] Th Von Karman | Über laminare und turbulente reibung[END_REF]. It can be classified into a family of NS solutions, where rotations play an important role [START_REF] Go | Note on a class of solutions of the navier-stokes equations representing steady rotationally-symmetric flow[END_REF][START_REF] Rj Lingwood | The effects of surface mass flux on the instability of the bek system of rotating boundary-layer flows[END_REF]. The essential feature of von Kármán flow is the appearance of cross-flow velocity inside boundary layer, which can be explained by the relationship between centrifugal and Coriolis forces in the rotating reference frame.

Let's consider the flow on a rotating disk in figure 1.11, which rotates with an angular velocity ω = (0, 0, ω z ). It is natural to work under the rotating reference frame in which the disk is stationary. Because rotating reference frames are non-inertial frames, two additional fictitious forces arise. The centrifugal force -ω × (r × ω) on a fluid parcel at r = (r, θ, z) (r, θ, z are radius, angle, and depth coordinates respectively in the cylindrical system, which is in the rotating reference frame) is (rω 2 z , 0, 0). And the Coriolis force -2ω × U on it would be (-2u θ ω z , 2u r ω z , 0), where U = (u r , u θ , u z ) is the velocity of the fluid parcel with respect to the rotating reference frame.

In the rotating reference frame, we can treat the fictitious forces like real forces, and pretend we are in an inertial frame [START_REF] Louis | Analytical mechanics[END_REF]. For a fluid parcel far away from the disk, it is not affected by the presence of the disk (the wall-normal velocity induced by the disk is neglected, which would not change the analysis) and does circular motion with radius r and u θ = rω z . The centripetal force needed for this circular motion is (-rω 2 z , 0, 0), which is exactly the combination of the corresponding centrifugal force and Coriolis force.

For a fluid parcel at the same r but inside the boundary layer, its circumferential velocity u θ is smaller than rω z because of non-slip boundary condition. Therefore, the Coriolis force, which depends on the velocity, is smaller than the inviscid case.

However, the centrifugal velocity force, which depends only on r, is the same as the inviscid case. As a result, the combination of the centrifugal force and Coriolis force is no longer large enough to keep the fluid parcel doing circular motion with radius r, but with a larger radius. This means that a radial velocity appears inside the boundary layer. The cross-flow velocity profile introduces an inward viscous stress in radial direction, which partially counteract the outward motion to the extent that an equilibrium is reached.

Figure 1.12 summaries the above process. In short, the radial (cross-flow) velocity can be seen as a secondary flow induced by the primary flow in circumferential direction. Later we will see that the process illustrated in the figure can explain the more complex phenomena on rotating blade boundary layer. For example, when flow separates, the Coriolis force changes its direction and acts in the same direction as centrifugal force. On the other hand, because the separation bubble is usually thicker than the attached boundary layer, the viscous force is small. These two factors would result a large radial flow around separation regions.

Figure 1.13 plots the velocity profiles of von Kármán solution from [START_REF] Lingwood | Absolute instability of the boundary layer on a rotating disk[END_REF]. The circumferential velocity is similar to the Blasius boundary layer. However, the radial velocity inside boundary layer is not zero as explained above. Because the fluid is pumped outwards near the wall, there is a wall-normal velocity so the continuity equation is sat-

1.
In the rotating reference frame, the velocity inside boundary layers decreases due to the none-slip boundary condition.

2. The radial Coriolis force inside boundary layers decreases because it linearly depends on the circumferential velocity.

3. The fluid particle can no longer do the perfect circular motion. i.e. a radial (cross-flow) velocity appears inside boundary layers.

The cross-flow velocity profile induces an inwards radial viscous force, which partially counteracts the Coriolis force

to the extent that an equilibrium is reached. Here the cross-flow velocity is towards the rotating axis. From [START_REF] Peter Alan | Turbulence in rotating, stratified and electrically conducting fluids[END_REF] isfied. It should also be noted that the Re in the rotating disk flow is usually defined as r ω/ν, where ν is the kinematic viscosity.

Another example of the rotating flow system is the Bödewadt layer, where the disk is at rest while the inviscid flow rotates. It can be analyzed in the rotating reference frame similarly. Like the von Kármán flow, there is a cross-flow velocity close to the disk. However, the cross-flow velocity in Bödewadt layer is towards the rotating axis.

A daily example is the flow inside a teacup. After stirring, the tea-leaf at the cup bottom would alway come to the center. The cross-flow instability has received much attention because of its importance in modern airliner swept wings [START_REF] William S Saric | Stability and transition of three-dimensional boundary layers[END_REF]. As a prototype for cross-flow transition, the rotating disk boundary layer has been under extensive theoretical, experimental, and numerical studies for its geometry simplicity. In the experiment of Gregory et al. [START_REF] Gregory | On the stability of threedimensional boundary layers with application to the flow due to a rotating disk[END_REF] around 30 stationary spiral vortices were observed on the disk. The vortices propagate in the direction predicted by the inviscid stability theory. Malik et al. [START_REF] Mujeeb R Malik | Instability and transition in rotating disk flow[END_REF] showed that the and the experiment (Re c = 297) of Kobayashi [START_REF] Kobayashi | Spiral vortices in boundary layer transition regime on a rotating disk[END_REF]. In the two important works of Lingwood [START_REF] Lingwood | Absolute instability of the boundary layer on a rotating disk[END_REF][START_REF] Rj Lingwood | An experimental study of absolute instability of the rotating-disk boundary-layer flow[END_REF], the author found that the boundary layer on the rotating disk is absolutely unstable when Re is larger than 510. This value corresponds well with the turbulence onset Re observed in several experiments. However, later it is shown that the cross-flow instability are only convective unstable in Falkner-Skan-Cooke (FSC) flows and swept wing boundary layer [START_REF] Rj Lingwood | On the impulse response for swept boundary-layer flows[END_REF][START_REF] Taylor | The long-time behaviour of incompressible sweptwing boundary layers subject to impulsive forcing[END_REF]. The global instability of rotating disks has also been discussed in many literatures [START_REF] Yim | Transition to turbulence in the rotating disk boundary layer of a rotor-stator cavity[END_REF].

Implication

As mentioned earlier, both HAWTs and marine propellers rotate during operation. As a result, the same reasoning in figure 1.12 can be readily applied. There should be a cross-flow velocity in their boundary layers. Therefore, a natural question to ask is whether cross-flow instability and transition can be observed on HAWTs and marine propeller blades.

HAWTs and marine propellers

The boundary layer's state could have a big implication on the aerodynamic of the bodies. A well-know example is the dimples on golf balls, which prompt the boundary layer transition. The turbulent boundary layer is usually harder to separate from the wall. The drag is therefore smaller and the golf ball can fly farther.

A good knowledge of the transition process is needed to model it in numerical simulations, which are widely used to predict and optimize the performance of turbines or propellers. Much effort has been done to understand the transition on rotating blades in the past decades.

HAWTs blades boundary layer

Effect of the rotation HAWT blades consist of a series of 2D airfoil sections with different thicknesses, chord lengths, and twist angles. In the design and optimization of HAWTs, Blade Element Theory (BET) is widely used. After dividing the blade into sections (or elements) in the spanwise direction, BET assumes that the flow around each section is locally 2dimensional. The lift and drag forces for each section can be obtained from momentum theory or calculated from 2D potential flow code such as XFOIL [START_REF] Drela | Xfoil: An analysis and design system for low reynolds number airfoils[END_REF]. The total torque of the blade can be then obtained by summing local sections' performance together [START_REF] Martin | Aerodynamics of wind turbines[END_REF][START_REF] Polat | Aerodynamic shape optimization of wind turbine blades using a parallel genetic algorithm[END_REF].

Nonetheless, the flow around the HAWT blade is not strictly 2D. As early as 1945, Himmelskamp [START_REF] Himmelskamp | Profile investigations on a rotating airscrew[END_REF] observed that the lift coefficient for the rotating blade is increased compared with the non-rotating case, and the stall is postponed. This phenomenon is referred to as rotational augmentation in the later literature. Fogarty [START_REF] Le Fogarty | The laminar boundary layer on a rotating blade[END_REF] solved the boundary-layer equations for rotating blade and found that the effect of rotation is small, which contradicts with engineers' impression. He suggests one of possible reasons might be that his analysis is only valid for attached flow region. McCroskey's flow visualization of rotating blade shows that the boundary layer transition on rotating blade resembles non-rotating blade, but there is a significant outward radial velocity in laminar separation bubble and in trailing-edge separated flow [START_REF] Mccroskey | Measurements of boundary layer transition, separation and streamline direction on rotating blades[END_REF]. Snel [START_REF] Snel | Sectional prediction of 3d effects for separated flow on rotating blade[END_REF] performed a scale magnitude analysis for rotating blade boundary layer equations, where he considered flow separation. His results show an improved C l prediction compared to 2D calculation. Corten [START_REF] Gp Corten | Flow separation on wind turbines blades[END_REF] pointed out that the boundary layer assumption in Snel's analysis is not necessary for separation regions. He also shows that the separated flow will move in outward radial direction due to centrifugal force.

Nowadays, it is well recognized that rotational augmentation is closely related to the large spanwise velocities in separation regions. On the one hand, the radial velocity pump fluid from separation region to spanwise direction, leading to volume reduction of separation bubble [START_REF] Pk Chaviaropoulos | Investigating three-dimensional and rotational effects on wind turbine blades by means of a quasi-3d navier-stokes solver[END_REF]. On the other hand, the radial velocity induces a streamwise Coriolis force, which partially counteracts the adverse pressure gradient [START_REF] Zhu | Combined effect of rotational augmentation and dynamic stall on a horizontal axis wind turbine[END_REF]. In section 1.5, we pointed out the similarity between the boundary layer flows on HAWT blades and on the rotating disk. Here we argue that the radial velocities in separation regions of HAWT blades can be analyzed similarly. Here we repeat the argument.

In figure 1.16 a blade is rotating in the plane (x, y), i.e., the angular velocity ω of the rotating reference frame fixed to the blade is (0, 0, ω z ). Considering a fluid parcel at r = (r, θ, z) (r, θ, z are radius, angle, and depth in the cylindrical coordinate system respectively), its velocity relative to the rotating reference frame U is (u r , u θ , u z ). The centrifugal force -ω × (ω × r) on this fluid parcel would be (rω 2 z , 0, 0) (in local radial, circumferential, and depth directions respectively). And the Coriolis force -2ω × U on it would be (-2u θ ω z , 2u r ω z , 0).

Far away from the blade, the fluid parcel experience circular motion with u θ = rω z .

The centripetal force needed for this circular motion is (-rω 2 z , 0, 0), which is exactly the combination of centrifugal and Coriolis forces. When flow detaches, the circumferential velocity decreases and reverses its direction. As a result, Coriolis force changes its direction and acts in the same direction as centrifugal force, which results in large radial velocity u r .

It is obvious to see that there should be a radial flow velocity inside the attached boundary layer as well. The argument on the rotating disk flow applies to HAWT blades without any changes (please refer to section 1.5). Indeed, this radial velocity was observed in [START_REF] Bosschers | Influence of blade rotation on the sectional aerodynamics of rotational blades[END_REF], where the boundary layer equations on a rotating blade was solved.

Boundary layer transition on rotating blades Due to the rotation and their large scale, it is very challenging to detect and measure the unstable waves HAWT blades in similarly to the flat plate. However, the transition locations on full-scale HAWTs are obtained using different experimental techniques.

In Reichstein et al. [START_REF] T Reichstein | Investigation of laminarturbulent transition on a rotating wind-turbine blade of multimegawatt class with thermography and microphone array[END_REF], the thermo-graphic imaging was used as a non-intrusive method for detecting flow transition on a MW wind turbine. The transition location from thermo-graphic imaging is in good agreement with results from microphones and transitional CFD simulations. Figure 1.17 is taken from [START_REF] T Reichstein | Investigation of laminarturbulent transition on a rotating wind-turbine blade of multimegawatt class with thermography and microphone array[END_REF]. It can be seen that when the wind turbine is in regular operation mode, the transition appears quite near the leading edge (x/c = 0.05).

In the DanAero MW projects [START_REF] Troldborg | Danaero mw: final report[END_REF][START_REF] Ha Madsen | Transition characteristics measured on a 2mw 80m diameter wind turbine rotor in comparison with transition data from wind tunnel measurements[END_REF], the transition location on a MW turbine is derived from surface pressure fluctuations measured by microphones. The transition location on the full scale turbine blade is advanced as compared to wind tunnel experiments. This difference of transition locations is caused by different characteristics of incoming turbulence in wind tunnel and field measurements.

Cavitation and marine propeller boundary layer

The marine propeller is another important type of flow machinery. Similar to HAWT blades, propeller blades are based on the hydrodynamic of 2D blade. However, propeller sections are defined on a series of concentric cylindrical surfaces as shown in Similar to the HAWTs, transition locations on propeller blades can be obtained experimentally. In the classic work of Kuiper [START_REF] Kuiper | Cavitation inception on ship propeller models[END_REF], he shows that although the exact transition location on propeller blades depends on a number of factors, the flow regimes can be generally divided into laminar and turbulent regions according to the direction of paint streaks. In the laminar flow region, the paint streaks point outwards towards the tip, whereas they are aligned to the circumferential direction in the turbulent region.

This streamline pattern difference is still the most commonly-used criteria to distinguish the flow status in marine propeller experiments [START_REF] Sánchez-Caja | Scale effects on tip loaded propeller performance using a ranse solver[END_REF][START_REF] Bhattacharyya | Scale effects on open water characteristics of a controllable pitch propeller working within different duct designs[END_REF].

Using oil-film interferometry, Schülein et al. [START_REF] Schuelein | Transition detection and skin friction measurements on rotating propeller blades[END_REF] obtained the friction coefficient on a whole blade of a high-speed aeronautic propeller. Some fine flow structures, not only in the attached flow but also in the separation bubble, are also clearly visualized in their paper. However, these works are more focused on transition locations on propeller blades instead of the transition process itself.

There are also numerous of CFD studies about the transition on propeller blade.

At present, almost all those studies rely on the RANS, where transition models are usually used to capture the boundary layer transition. Transition models in the RANS usually include empirical correlations and their effectiveness depend on the transition scenarios [START_REF] Robin | Correlation-based transition modeling for unstructured parallelized computational fluid dynamics codes[END_REF][START_REF] Medida | A new crossflow transition onset criterion for rans turbulence models[END_REF]. Therefore, the understanding of possible transition mechanisms is vital for modeling the transition problem. One of the most popular transition models is γ -Re θ model developed by Menter et al. [START_REF] Florian R Menter | A correlation-based transition model using local variables-part i: model formulation[END_REF]. Bhattacharyya et al. [START_REF] Bhattacharyya | Scale effects on open water characteristics of a controllable pitch propeller working within different duct designs[END_REF] shows that the inclusion of γ -Re θ model results a better agreement with experimental paint test in terms of flow patterns on the propeller. Pawar and Brizzolara [START_REF] Pawar | Relevance of transition turbulent model for hydrodynamic characteristics of low reynolds number propeller[END_REF] show that γ -Re θ is able to correctly predict the complex flow phenomena on propellers like leading edge separation. However, the inlet condition need to be tuned RANS in order to get a better agreement with experiment. In the Moran-Guerrero et al. [START_REF] Moran-Guerrero | On the influence of transition modeling and crossflow effects on open water propeller simulations[END_REF], they investigated the boundary layer transition on a propeller using γ -Re θ that is able to take cross-flow into account. They show that the addition of cross-flow terms in the γ -Re θ model can promote the transition onset and result a larger turbulent region.

Interaction with cavitation

The cavitation is perhaps the most important and challenging problem on ship propellers. As mentioned above, propellers work based on the pressure difference between the two sides of the blades. When the pressure on the suction side decreases too much locally, the fluid could become vapor rapidly. When the bubble ( or cavity ) filled with vapor is convected to downstream where the ambient pressure is large, it collapses. During the implosion, there is usually a high speed water jet which can reach supersonic speed. If the water jet hits the propeller, it exerts a considerable force on the body. Because the cavitation can generate numerous bubble, it is a main cause of propeller wear and failure. In the transition flow region, the pressure fluctuates at a significant magnitude. Previous studies show that the cavitation might occur either at the location of minimum C p or at the transition region [START_REF] Young | The influence of hydrofoil oscillation on boundary layer transition and cavitation noise[END_REF]. The fluctuation in the laminar separation bubble might affect the inception of sheet cavitation [START_REF] Vh Arakeri | Viscous effects in the inception of cavitation on axisymmetric bodies[END_REF][START_REF] Croci | Attached cavitation in laminar separations within a transition to unsteadiness[END_REF][START_REF] Che | Control effect of micro vortex generators on leading edge of attached cavitation[END_REF] whereas the turbulent spot during transition can affect the bubble cavitation [START_REF] Van Terwisga | Tu delft: Opencourseware: Cavitation on ship propellers[END_REF].

Summary

At the beginning of this chapter, we introduced some basic fluid dynamic concepts such as airfoil aerodynamics, boundary layers and the laminar-turbulent transition. Then we explained the origin of cross-flow on the rotating disk and its transition.

After having pointed out the similarity between the rotating disk and HAWTs and marine propellers, we reviewed the related works on rotating blades' boundary layer.

Although many works has been done, they are more focused on the transition location on HAWT and propeller blades. Little is known about the transition process itself on rotating blades.

Indeed, the complex operating conditions, twisted geometries and unsteadiness make the accurate measurements of the boundary layer profile and transition by experiments to be very difficult. On the other hand, CFD offers a feasible approach, which can get the full information of the flow fields. To the best of our knowledge, there is no numerical simulation which fully resolve the laminar-turbulent transition on rotating HAWT or propeller blades. This work aims partially fill this gap.

Chapter 2

Methodology

In this chapter, we will first present the stability theory. Then we offer a description of Spectral Element Method, which is the main numerical tool used in this thesis. However, both subjects are so developed that it is impossible to give an extensive explanation. Therefore, we will focus on the basic concepts and terminology. References are given as much as possible.

Hydrodynamic instability

Stability theory studies the response of a system when it is disturbed. The toy in figure 2.1 together with a horizontal plane define a stable system. The base state of the toy is simply stationary. Any disturbance added to it would decay due to frictions and the system would finally return to the base state. It likes the Gömböc which would always return to the same position.

The laminar-turbulent transition is usually accompanied by the amplification of small disturbances. The hydrodynamic instability theory attempts studies the evolution of small disturbances in the laminar flow. If the disturbance increases in its amplitude, the laminar base flow is said to be unstable and could transition to turbulence. The development history of the linear instability theory can be found in the monograph of Drazin and Reid [START_REF] Philip | Hydrodynamic stability[END_REF]. It is an example of a stable system. From taobao.com (b) The Gömböc, which is a homogeneous body having just one stable and one unstable point of equilibrium. From Wikipedia.

The Linearized NS equations

The low-speed motion of Newtonian fluid is governed by the incompressible Navier-Stokes equations, which can be written as:

∇ • U b = 0 ∂U b ∂t + U b • ∇U b = -∇p + 1 Re ∆U b (2.1)
Suppose that we have a laminar base flow field (U b , p b ) which satisfies eq. 2.1.

After introducing a perturbation (U , p ) to the base flow, the total flow field (U + U b , p b + p ) also satisfies eq.2.1, which now is written as:

∇ • (U b + U ) = 0 ∂(U b + U ) ∂t + (U b + U ) • ∇(U b + U ) = -∇(p b + p ) + 1 Re ∆(U b + U ) (2.2)
Subtracting eq. 2.1 satisfied by the (U b , p b ) from eq. 2.2, we obtain the governing equations of the disturbance (U , p ). Because the amplitude of the disturbance is generally very small, the nonlinear terms of the disturbance can be neglected. The Linearized NS equation can be written as:

∇ • U = 0 ∂U ∂t + U • ∇U b + U b • ∇U = -∇p + 1 Re ∆U (2.3)
By solving eq. 2.3, we can obtain the evolution of a perturbation in the laminar base flow. If the amplitude (or energy) of the perturbation increases, then we say that the base flow is unstable. However, eq. 2.3 are still a set of Partial Differential Equations and are elliptic. Solving them is as expensive as solving the original NS equation.

Linear Stability Theory

With additional assumptions, eq. 2.3 can be further reduced. Consider the plane Poiseuille flow in a channel formed by two parallel walls, which is shown in figure 2.2 and has a steady base flow solution:

U b = (U b x , U b y , U b z ) = (1 -y 2 , 0, 0) p b = -2x/Re + C. (2.4)
The Re is based on the centerline velocity and the half height of the channel.

Suppose the disturbance had a wave form in the streamwise and the spanwise direc-tions, it can be then written as: U = Û (y)e i(βz+αx-ωt) + c.c. p = p(y)e i(βz+αx-ωt) + c.c.

(2.5) where α and β are complex wave-numbers of the disturbance in streamwise and spanwise (or radial) directions respectively, and ω is a complex frequency. Û and p are the shape of the disturbance along y, or, eigenfunction. If α i < 0, β i < 0, or ω i > 0, the wave-amplitude in eq. 2.5 would grow exponentially either with time or space.

Substitute eq. 2.4 and eq. 2.5 into equations 2.3 and rearrange the result, we can obtain the O-S (Orr-Sommerfield) equation [START_REF] Peter | Stability and Transition in Shear Flows[END_REF]:

{(D 2 -k 2 ) 2 -iαRe[(U b x -ω/α)(D 2 -k 2 ) -D 2 U b x ]}U y = 0 (2.6)
where D is the differential operator in y direction.

k = α 2 + β 2 .
The O-S equation is a fourth order Ordinary Differential Equations (ODEs). It defines an eigenvalue problem in that the solution U y (eigenfunction) only exists with certain values of ω, α, and β. In the spatial mode analysis, ω i is set to be 0. If the imaginary part of the streamwise wave-number α i is negative, the flow is unstable and the disturbance grows along streamwise direction. The 3D boundary layers usually have more component than just the streamwise velocity. In that case, the final governing equation of linear stability analysis is usually a system of ODEs instead of the O-S equation ( [START_REF] Mujeeb | Numerical methods for hypersonic boundary layer stability[END_REF]). However, all the concepts in O-S equation are still valid and the final system is an eigenvalue problem too.

Further notes on LST

The parallel assumption on base flow

The above derivation of O-S equation relies on two assumptions. Firstly, the base flow only has a streamwise component and is uniform along z and x directions. The boundary layer flow, e.g. Blasius solution, is usually not uniform in streamwise direction.

However, the streamwise variation of the flow is generally small and can be neglected (the local parallel assumption). Therefore, the O-S equation is still applicable.

LST is a local analysis method in that it depends only on the velocity profile at one location. In certain situations, the base flow's variation in streamwise or spanwise direction could not be neglected (e.g. flow around roughness elements). For such flows, the Global Stability theory is developed (see e.g. [START_REF] Theofilis | Global linear instability[END_REF]), where the base flow and the eigenfunction are 2D or 3D.

The normal mode assumption on the disturbance

The second assumption in LST is that the disturbance can be written in wave forms. This is referred to as normal mode analysis. It implies that the small disturbances with different wave-numbers grow independently and do not interact with each other.

Whereas this is true for the T-S wave and cross-flow transitions, other growing mechanism can play an important role in other flows (for instance Couette flow [START_REF] Lloyd N Trefethen | Hydrodynamic stability without eigenvalues[END_REF]).

Comparison of different simulation strategies for transitional flows

With the advent of high performance computation, the numerical simulation becomes an indispensable research method in fluid dynamics. The numerical simulation provides many observations and insights on the boundary layer transition. In stability analysis, the numerical simulation is also needed to obtain the base flows around complex geometries.

Transitional flows and turbulent flows share many similarities. Although we focus on the laminar-turbulent transition of the boundary layer, it is necessary and useful to firstly introduce the basic concepts in turbulent simulation. One of the biggest challenges to simulate turbulence is its multi-scales nature. From the perspective of the energy, the kinematic energy in large eddies is successively transferred to smaller ones as they break into smaller eddies. This process continues until the energy being dissipated by the viscous effect at the smallest scales. Figure 2.4 illustrate this process.

It can be shown that the ratio between the large eddy length scales and dissipation length scale is Re 3/4 [START_REF] Stephen B Pope | Turbulent flows[END_REF]. This implicates that if one wants to resolve all the flow details in 3D, the total mesh point N 3 is of the order of Re 9/4 . Considering that the time-step dt usually has to be smaller when the mesh is finer, the total operation is O(Re 3 (log Re 3 )) [START_REF] Gary | A primer on direct numerical simulation of turbulence-methods, procedures and guidelines[END_REF]. If the Re doubles, the computational cost increases by a factor of 11! For a MegaWatts (MW) wind turbine, the Re based on total blade length is O(1 × 10 7 ). The Re of a modern airliner is O(1 × 10 8 ). Even with modern cluster computers, there is no way to directly simulate the turbulent flow with such large Re. In short, the demand for faster calculation is insatiable in the fluid dynamic.

RANS

For engineering applications, we are often only interested in the time-averaged forces and moments induced by the flows, therefore there is no need to resolve all the flow details. The Reynolds Averaged Navier-Stokes(RANS) equations are based on this idea.

RANS equations describes the time-averaged behavior of turbulent flows. However, RANS equations are not closed. The relation between the Reynolds stress R ij and the fluctuation velocity need to be described by empirical models [START_REF] Joel H Ferziger | Computational methods for fluid dynamics[END_REF]. There are many models exist but their effectiveness generally depends on the flows. Most of RANS models are based on the concept of eddy viscosity, which relates the Reynolds stress to the mean flow quantities.

As shown in figure 2.4, RANS only resolves the very large scale motions of turbulent flow. Any unsteady small scale motion is averaged out. As a result, it is not able to capture the laminar-turbulent transition process. To take into account transition, extra transition models are often employed in RANS simulations. Almost all the tran-sition models are based on some empirical relations between local flow quantities and transitions.

LES

Large Eddy Simulation (LES), as its name indicates, aims to resolve the large energycontaining motion of turbulent flows (figure 2.4). It is firstly introduced in the simulation of atmospheric flow [START_REF] Smagorinsky | General circulation experiments with the primitive equations: I. the basic experiment[END_REF]. As the development of more and more powerful computers, it is making its way towards industry applications.

Although the large motions of the turbulence are geometry dependent, according to the Kolmogorov theory, the small scale motion of turbulent is universal when Re is sufficiently large. Therefore we can only capture the geometry-dependent large motions, while leave the small scale motion's effect to be modeled. Since their effect is mainly dissipating the energy, the modeling can be realized by introducing a kind of artificial viscosity.

A finer grid than RANS is needed in LES, but the gird size is usually larger than the dissipation scale. Therefore the LES models are sometimes called sub-grid models because they emulate the effect of the eddies smaller than the grid. Generally speaking, there are two broad classes of LES models.

Explicit LES

The traditional LES (or Explicit LES) solves the filtered NS equations [START_REF] Misra | A vortex-based subgrid stress model for largeeddy simulation[END_REF]. Similar to RANS, filtered NS equations are not closed. Explicit LES try to model the subgrid stress τ ij from the resolved flow field. In the classic Smagorinsky model, an eddy viscosity term ν r is added to the kinematic viscosity:

ν r = (C s ∆) 2 S (2.7)
where ∆ is filter width, S = 2S ij S ij is magnitude of strain-rate tensor S ij of the resolved flow. The coefficient C s is Smagorinsky coefficient and is around 0.1 -0.2.

A constant C s , does not apply for transitional flow because the subgrid stress (unresolved stress) is different for laminar and turbulent flows. In the dynamic Smagorinsky model, the C s is determined by local resolved flow field and adjust accordingly [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF] [90]. Because the result coefficient is usually quite volatile across the field, averaging techniques are usually needed in the dynamic Smagorinsky model [START_REF] Meneveau | A lagrangian dynamic subgrid-scale model of turbulence[END_REF].

Implicit LES

Implicit-LES (ILES) solve the unmodified NS equations directly. However, the mesh resolution in ILES is not fine enough to capture the smallest scale of the turbulence.

The numerical dissipation is used as a subgrid model.

The differential equations after discretization are not the same as the original ones.

For relatively simple schemes like finite difference method, it can be shown that the truncate error of a convective problem includes even-order derivatives, which behaves like diffusion terms [START_REF] Joel H Ferziger | Computational methods for fluid dynamics[END_REF].

Consider the following convection problem, which is adapted from [START_REF] William | Numerical partial differential equations: finite difference methods[END_REF]: As can be seen, the amplitude of the wave is damped as it is being convected.

v t + v x = 0 x ∈ (0, 1) v(x, 0) = sin 40 (πx) x ∈ (0, 1) v(0, t) = v(1, t) t ≥ 0 (2.
At the same time, the wave is diffused to a lager range as well. However, as mentioned above, eq. 2.8 is a pure convection problem, the wave shape should stay the same as being convected. The dumping and diffusion are due to numerical effect. If we define the energy of the wave as E = 1 0 v 2 /2, its value decreases as the wave being dumped. Figure 2.6 shows the energy variation as the waves being convected. For smaller mesh size, the numerical dissipation is smaller.

For more sophisticated numerical schemes, the numerical dissipation also exits. After tuning, the numerical dissipation is able to serve as a sub-grid model in LES ( [START_REF] Adams | Implicit large-eddy simulation: theory and application[END_REF][START_REF] Moura | Linear dispersiondiffusion analysis and its application to under-resolved turbulence simulations using discontinuous galerkin spectral/hp methods[END_REF]). This approach is refereed to as Implicit-LES because there is no explicit LES models in the solver. It is also called as 'Under-resolved DNS (UDNS)' because the only difference between ILES and DNS is the mesh resolution. If the ILES mesh is successively refined, the solutions should be closer and closer to DNS result.

Because there is no explicit filtering operation, ILES is faster than classic LES.

However, the effectiveness of ILES often depends on the numerical method and meshes [START_REF] Domaradzki | Effective eddy viscosities in implicit modeling of decaying high reynolds number turbulence with and without rotation[END_REF], [START_REF] Zhou | Estimating the effective reynolds number in implicit large-eddy simulation[END_REF]. Good familiarity and understanding about the numerical method are often needed for ILES.

DNS

DNS (Direct Numerical Simulation) solves the NS equations without any turbulence modeling. Since DNS has to resolve the smallest dissipative scale, it requires very fine mesh. Although DNS is able to tackle flow problems with higher and higher Re since its birth (see [START_REF] Ishihara | Study of high-reynolds number isotropic turbulence by direct numerical simulation[END_REF][START_REF] Lee | Direct numerical simulation of turbulent channel flow up to τ ≈ 5200[END_REF]), its applications are still restricted to relatively simple geometry and relatively low Re. There are very few DNSs which deal with practical industrial flow problems.

However, DNS provides a non-intrusive way to obtain all the flow information. It is an important tool to valid theories and develop new turbulent models for RANS and LES methods ( [START_REF] Moin | Direct numerical simulation: a tool in turbulence research[END_REF][START_REF] Wu | Reynolds-averaged navierstokes equations with explicit data-driven reynolds stress closure can be illconditioned[END_REF]).

The nature of DNS requires low-dissipation schemes such as compact schemes or spectral method, which will be the subject of the following section.

Choice of the numerical method

In this thesis, we are more interested in the transition process instead of the turbulence.

The DNS up to Kolmogorov scale is not necessary as long as the transition is well resolved. Therefore, we chose to use LES as the main tool.

Another issue is the choice between the Explicit LES and the ILES. We extended the original explicit LES model in Nek5000 to general 3D problems. The details are given in the appendix. However, our tests show that the Explicit LES doesn't possess significant advantages over ILES in terms of transition studies. One of the reasons is that Explicit LES would result some eddy viscosity as long as there are small structures in flow fields. However, during the transition, the flow before fully breakdowns has fine patterns but it is still laminar technically. We don't want any artificial viscosity in this situation.

On the other hand, the Explicit LES involves many operations such as filtering, which can be costly. As a result, the following simulations are based on ILES (or UDNS). That said, the meshes used are finer than the commonly-seen LES because the unstable waves during transition are quite small.

Spectral Element Method

Spectral Method

In spectral method, unknown solutions of differential equations are approximated by the sum of 'basis functions'. The commonly-used basis functions include Fourier series, Legendre series, and Chebyshev series. Spectral method has two important advantages.

It has low numerical dissipation and converges fast.

Perhaps the most commonly-seen spectral method is Fourier Spectral Method, where the trigonometric functions are used as basis functions. It is usually used in periodic problems. For instance, DNS of 3D homogeneous isotropic turbulence is usually based on Fourier Spectral Method. Fourier series can be also used only in certain directions instead of in all of 3 directions. For example, in turbulent channel simulation, one can use Fourier series in the two directions parallel to the wall and use other methods in the wall-normal direction [START_REF] Lee | Direct numerical simulation of turbulent channel flow up to τ ≈ 5200[END_REF].

For non-periodic problem, the Fourier Spectral Method can be used with modifications (e.g. the Fourier expansion is used in stream-wise direction of flat plate boundary layer in [START_REF] Chevalier | Simson-a pseudo-spectral solver for incompressible boundary layer flow[END_REF]). But basis functions like Chebyshev polynomials and Legendre polynomials are more suitable choices in these cases. One of the earliest application of Chebyshev method is in hydrodynamic instability [START_REF] Steven | Accurate solution of the orr-sommerfeld stability equation[END_REF], where the Orr-Sommerfield equation is solved numerically. However, both the Fourier Spectral Method and the Chybeshev Method are usually used for single-domain problems with simple geometry. Their lack of geometry flexibility severely restricts their application in more complex flow problem. This shortcoming is overcome by Spectral Element Method (SEM).

Spectral Element Method

Spectral Element Method (SEM) is firstly introduced by Patera [START_REF] Anthony T Patera | A spectral element method for fluid dynamics: laminar flow in a channel expansion[END_REF]. It is a hybrid method combining the geometry flexibility of Finite Element Method and the high accuracy of Spectral Method. The model problem 2.8 is used again to illustrate the basic idea of FEM. Recall the model problem 2.8:

v t + v x = 0 x ∈ (0, 1) v(x, 0) = sin 40 (πx) x ∈ (0, 1) v(0, t) = v(1, t) t ≥ 0 (2.9)
In SEM, the exact solution of 2.9 is approximated by a series expansion v(t, x) ≈ N i=1 c i (t)φ i (x), where N is the number of basis functions used in the expansion; φ i (x) is the basis function (or trial function) on [0, 1], c i (t) is the time-dependent coefficient for the ith basis function. Since this expansion is only an approximation of the exact solution, there is an error when one plug the expansion into the differential equation:

N i=1 c i (t + dt) -c i (t) dt φ i (x) + N i=1 c i (t)φ i (x) = R(x) (2.10)
Where the time derivative is approximated using Eular scheme.

We want to minimize the residual R(x) on the domain [0, 1]. This is done with the help of another function series called test (weight) function ξ(x) j . We require that the R(x) be orthogonal to ξ(x) j for every j.

1 0 ξ j (x)R(x) = 1 0 ξ j (x) N i=1 c i (t + dt) - 1 0 c i (t) dt φ i (x)+ 1 0 ξ j (x) N i=1 c i (t)φ i (x) = 0 (2.11)
Different choices of the test function lead to different formulations. In collocation method, the test function ξ j is the Dirac delta function δ(x j ), where x j is a set of point.

In other word, collocations method requires R(x) to be zero at the points x j .

A more common choice is to let the test function ξ(x) be the same as the trail function φ(x). This formulation is refereed to as Galerkin method. Rearrange eq. 2.11 and we have:

N i=1 c i (t + dt) -c i (t) dt 1 0 φ j (x)φ i (x) + N i=1 c i (t) 1 0 φ j (x)φ i (x) = 0 (2.12)
For each test function, we have an equation like 2.12. There are totally N equations.

So we can write them in matrix form:

A c(t + dt) = A c(t) -dtB c(t) (2.13)
where In the Spectral Element Method, the Lagrange polynomials l i used in elements is based on Gauss-Lobatto-Legendre (GLL) points. This is the most important difference with Finite Element Method, which uses even-spaced point. Figure 2.8 shows the five 4th order l i s based on GLL points. GLL points are in the domain [-1, 1], but they can be readily transformed into the domain of any local element. Compared with evenly distributed points, the space between points near the boundaries is smaller than the center. This distribution can eliminate the Runge's phenomenon, which appears on even-spaced grids. Therefore high order schemes are possible in SEM.

A ij = 1 0 φ i (x)φ j (x) is the mass matrix, B ij = 1 0 φ i (x) φ j (x
Using GLL Lagrange polynomials is refereed to as the nodal approach because the coefficients c i equals to the value of v at corresponding nodes. SEM can be also based on the modal basis, where the functions in each element are the combinations of Legendre polynomials. Although these two approaches are equivalent, modal approach could be more convenient when explicit filtering is needed. dissipation of in the SEM method is quite small. On the other hand, similar to FEM, SEM has much larger geometry flexibility than Spectral Method.

For 2D and 3D problem, the basic idea of SEM is the same. However, the formulas in 2D and 3D are much more complex than 1D. Many special numerical treatments are also necessary. For more information about the formulation and application of SEM in fluid dynamic, please refer to [START_REF] Karniadakis | Spectral/hp element methods for computational fluid dynamics[END_REF], [START_REF] Mo Deville | High-order methods for incompressible fluid flow[END_REF].

The SEM belongs to hp methods. The mesh convergence study and refinement can be done either by increasing the order of spectrum in the elements (p refinement) or generating new mesh (h refinement). In this thesis, we rely on p refinement after obtaining the meshes for elements.

Nek5000

The NS solver used in the present study is Nek5000 [START_REF] Fischer | Open source spectral element cfd solver[END_REF], which is a highly scalable solver based on SEM. It is widely used in the studies of hydrodynamic instability and receptivity problems because of its high accuracy, geometry flexibility, as well as extensibility.

Nek5000 was developed by Argonne National Laboratory in U.S. several decades ago. It has been actively maintained and now supports many useful features like deforming mesh, conjugate heat transfer, and Magnetohydrodynamics (MHD) equations etc. Here we offer a description of the two features that is relevant to the present study. elements has to be hexahedron. As a result, the meshes are usually generated as structured ones in meshing softwares. Then they are transformed into the format needed by Nek5000.

For certain complex geometries, it could be very difficult to generated one single structured mesh even with block-structured mesh. Nek5000 provides an elegant solution to this kind of problem, which is the overlapping overset meshes techniques. Walsh's solution [START_REF] Walsh | Eddy solutions of the navier-stokes equations[END_REF] is calculated using the overlapping meshes in figure 2.10a.

The result x velocity at t = 0.04 was given in the contour 2.10b. It can be seen that the flow field is continuous across the boundaries. It looks like that there is only one single domain.

In chapter 5, we will show why this technique is necessary for the marine blade simulation.

NS equations in the rotating reference frame

For the HWAT and marine propeller blades simulations, it is a natural choice to let the meshes (therefore the geometry) to be stationary during the calculations. So the moving mesh is avoided. The effect of rotation is taken into account by solving the NS equations in the rotating reference frame.

Suppose the velocity of the fluid is W for a stationary observer in inertial frame.

Then W satisfies the NS equations 2.1.

Suppose there is another observer who is stationary to the HWAT or the propeller (i.e. he or she rotates with the HWATs or propellers). The velocity of the fluid with respect to this observer is U . W and U are related by:

W = ω × r + U (2.14)
Where r is the position vector of the fluid with respect to the second observer. Substitute 2.14 into 2.1 satisfied by W , we will get:

∇ • U = 0 ∂U ∂t + U • ∇U = -∇p + 1 Re ∆U -ω × (ω × r) Centrifugal term -2ω × U Coriolis term (2.15)
The second observer defines a rotating reference frame, in which the NS equations In the following sections, all the velocities refer to the velocity in rotating reference frames. As mentioned before, Nek5000 has a good extensibility. The centrifugal and Coriolis terms is modeled by the volume force in the simulations. It should be also noted that the rotation axis of the HAWT blade and the propeller blade are aligned to z direction. I.e. ω = (0, 0, ω z ). As a result, The velocity component in z direction U z , which is parallel to the angular velocity vector, does not affect the Coriolis force.

Similarly, z coordinate does not affect the centrifugal force.

IDRIS

All the simulations in this thesis are realized on the cluster computers in the IDRIS

Summary

This chapter described the main tools at our disposal for the boundary layer transition study. We first introduced the linear stability theory and shortly discussed its limitation.

Then we described the different strategies in the turbulent modeling with emphasis on the correlation between numerical dissipation and ILES.

We demonstrated the basic idea of SEM using a simple model problem. It is a

Galerkin method so has a good geometry adaptivity. On the other hand, it is a spectral method and therefore has low numerical dissipation. The code Nek5000 used in the present study is introduced. At last, we briefly introduced the super computer we used.

Chapter 3

Rotating HAWT blade

This chapter presents the results on a HAWT blades. We will first describe the set-ups of the simulations such as the blade, the parameters as well as the meshes. Then the main results from LES and stability analysis are given. Followed by discussion and implication of the present results.

Set-ups of the simulations

The HAWT blade

The model under investigation is the LM38.8 blade from the DanAero MW project, whose section is a series of NACA63-4XX profiles (XX is the maximum thicknesses).

It has a 1.24-meter extension in the root, so the total span length R is 40.04 meter. The simulation parameters are taken from the first round of calculations of DanAero MW experiments: the incoming wind velocity U ∞ is uniform and equals to 6.1m/s, the blade rotates with an angular velocity ω z = 12.3rpm. However, to make the simulation feasible, we reduced the Reynolds number by one order. For that purpose, the kinematic viscosity ν is set to 2 * 10 -4 m 2 /s, which is about one order larger than that of air. The induction effect of the wind turbine is not taken into consideration in the present study.

However, it can be included by changing the normal incoming velocity U ∞ to 9.15m/s if the induction factor a = 0.33 for the ideal wind turbine is used.

The meshes and boundary conditions

Even with the reduced Reynolds number, it is still impossible to resolve the whole blade with our computational capacity (please see figure 3.2). So we restrict the simulation domain from y = 15m to y = 35m, where the chord Reynolds' number is roughly constant (figure 3.1). For Section5 and Section8 marked in figure 3.1 (the names follow the convention in DanAero MW project), 2D airfoil (uniform span and without rotation)

wind tunnel measurements of pressure coefficient C p are also available for each section.

In the simulation, 42 sensors are put along Section1 to record the variation of flow quantities with time. The total element number is 0.67 million (wake mesh included), which means the total degree of freedom is 0.49 billion if order 10 is used in each element.

On the blade surface, non-slip wall BC is applied. Outflow BC is used on the vertical surface behind the trailing edge. Dirichlet BC is used for the surface around the blade. The BC on the two ends of spanwise direction is hard to define because the computational domain is truncated from a whole blade and the flow quantities on those boundaries are not known. We simply impose Dirichlet BCs on them. On all the Dirichlet BCs, the velocities at position (x, y, z) are specified as U x = ω z y, U y = -ω z x, and The mesh for the airfoil simulation is uniform in the spanwise direction and has a roughly same resolution as the 3D blade simulation. The blade rotation is not considered. When airfoil sections are tested in wind tunnels, it usually has a uniform span and does not rotate. The airfoil simulation serves the purpose of replicating wind tunnel experiments.

U z = U ∞ ,

Results

As mentioned in the section 2.3.3, the NS equations are solved in the open source code Nek5000. For the 3D blade simulations, the Coriolis and centrifugal terms are included.

The solver NEK5000 is based on the Cartesian coordinate. However, when circum-ferential and radial velocities are refereed to, the velocities in x -y plane are transformed into the cylindrical coordinate system that shares the same z as the Cartesian frame. Obviously, without the blade, radial velocity u r would be zero.

A proper local coordinate system is also needed when velocity profiles are concerned. In figure 1.16, n is the wall-normal direction of a point on the blade wall.

Streamwise velocity is obtained by projecting circumferential velocity and U z to wall tangential direction (perpendicular to n). The local radial velocity is the same as that in cylindrical coordinate. Another note is that all the velocities are normalized by the incoming flow velocity U ∞ .

The influence of BCs at two ends

To evaluate the influence of BCs at span end on the flow field in the center blade, simulation of a more truncated case is performed. Its span range is from y = 18m to y = 32m. All the set-ups for this smaller case is the same as the standard case except for the domain size.

Figure 3.4 compares the velocity profiles at center blade. Both simulations are UDNS with spectrum order 8. The streamwise velocities at pressure side agree so well that they are indistinguishable. On suction side, a small difference of the two results can be observed near the wall. There is a weak separation for standard domain results, whereas the smaller domain profile is barely attached.

The radial velocity profiles in figure 3.4b are typical cross-flow profile, i.e., the maximum velocity appears inside boundary layer. So an inflectional point exists and the flow is inherently unstable. There are some small differences between the two simulations in the figure 3.4b, but their shapes and magnitudes agree quite well. It can be concluded that the domain size is large enough to avoid that the BCs at the span end affect the flow field in center of the blade.

Mesh Convergence

For the standard domain 3D rotating blade, simulations have been performed with three different orders 6, 8, 10. The order 6 results are omitted for the sake of clarity. Figure 3.5a plots the pressure coefficients along Section8 from the simulations and the 2D wind tunnel experiment. It can be seen that the numerical simulation results overlap with each other in most region. The only discrepancy appears when there is sudden pressure increase (x/c l = 0.4 on suction side and x/c l = 0.8 on pressure side). The order 8 result fails to capture the small pressure plateau on the mid-chord of suction side, which is caused by separation bubble. On the pressure side, the separation happens quite late and is near the trailing edge.

The chord Reynolds' number of wind tunnel experiments is 1.5 million, which is five times larger than the Re c of numerical simulation. Despite the Reynolds numbers difference, the C p agrees well with each other. In fact, when the Re c is large enough, the boundary layer is thin enough to not affect C p very much. lution does not fully capture the turbulent profile. Overall, the mesh resolution is fine enough to give a mesh-independent laminar solution.

Figure 3.6 compares the spectrum of the velocity signals from the order 8 and order 10 results. At x/c l = 0.36, although the flow is still laminar, the TS wave has developed a significant amplitude (which will be analyzed further in the next section). Therefore there is a frequency peak in both simulations (figure 3.6a). However, the frequency peaks differ slightly in their values, which means that the unsteady transitional flow is not fully mesh independent. At x/c l = 0.63, the flow is fully turbulent. Both simula-tions show similar spectrum distributions. As expected, the order 10 resolves more high frequency components than order 8 in the figure 3.6b.

2D simulations with different orders are also performed. In the follows sections, only the highest resolution (order 10) results will be shown. The radial velocities inside the boundary layer are typical cross-flow profiles, i.e., the maximum velocity appears inside the boundary layer. So an inflectional point exists and the flow is inherently unstable. On the suction side, the radial velocity is small near the leading edge, where the flow is fully attached. It progressively increases until reaching its maximum value (U r = 0.2 ) just before the separation bubble (x/c = 0.4).

Flow fields

However, the largest radial velocity is lower than 5% of the corresponding streamwise velocity. Although the cross-flow profile is inherently unstable, its magnitude must be large enough to overcome the viscous effect. In rotating disk and swept wing boundary layers, the cross-flow velocity has to reach 10% -20% of the streamwise velocity to intrigue cross-flow transition. Furthermore, the cross-flow component has to maintain a rather long distance to let the disturbance grow. The current result shows that cross-flow transition is unlikely to happen in rotating HAWT blade since cross-flow's magnitude is too small and it is restricted to the separation region which is quite short. The first noticeable feature is the spanwise streaks around mid-chord of the blade, which are formed by TS waves during laminar-turbulent transition. After the transition, small structures appear in the turbulent region. An interesting feature is that the TS waves are not perfectly aligned to the spanwise direction, i.e., the spanwise wave-number of the From the coherent structures, the flow field can be divided into three regions. In the two regions near the tip and root, the mesh resolution in y direction is not enough. As a result, the coherent structures show some patterns related to the mesh (please refer to figure 3.3). For the center portion from around y = 23 to y = 27, the laminar to turbulent transition and the small turbulent structures are captured.

Boundary layer instability and transition

The skin friction on the suction side of the blade is given in figure 3.11. The c f near the tip and the root is effected by the non-physical boundary condition and is unusually large. But their regions are restricted to the first two layers of the mesh in spanwise direction, which is another evidence that the non-physical BC has a neglected effect on the blade center.

Along the spanwise direction, The c f roughly increases at large y because that the local chord Reynolds' number increases. The boundary layer is relatively thiner at larger y. Along one section, the c f is large near the leading edge because that the the TS waves incept around x/c l = 0.4. They grow quickly and break down after three waves length. After the breakdown, small structures appear and the flow becomes turbulent. For the airfoil simulation, the vortex shedding is also visible. The TS waves in the 2D simulation are 2D (the streak are perfectly aligned to spanwise direction and therefore the spanwise wave number is zero), whereas the TS waves in 3D rotating blade simulation are 3D waves. As a consequence, one end of TS waves can be seen.

However, the small difference between 2D and 3D flow is hardly of any importance from an engineering perspective since the transitional flow region is rather small. When the base flow is 2D (no cross-flow component), the most unstable wave is 2D

(β r = 0), and 3D unstable wave's growth rate is symmetric along β r = 0. However, as shown in figure 3.14b, when cross-flow exist, the most unstable waves β r = 0. This might be the reason that the TS wave in figure 3.10 is not aligned to spanwise.

The flow structures on the pressure side of the blade is given in figure 3.15. As indicated in the C p distribution, the flow transition happens later than the suction side. This is a result of the more favorable pressure distribution along the chord. Similar to the suction side, the weak separation bubble is basically 2D.

Summary

In this chapter, numerical simulations resolving the boundary layer laminar-turbulent transition are performed for the flow around the rotating LM38.8 blade. A corresponding airfoil simulation with a uniform span is also performed to assess the effect of rotation on this laminar to turbulent transition. It is shown that for the current configuration, the local flow field and pressure distribution around the 3D rotating blade can be well approximated by the corresponding 2D airfoil flow except for the small radial In the rotating 3D blade boundary layer, a radial velocity arises because the Coriolis and centrifugal forces cannot keep fluid doing a perfect circular motion. This is very similar to the rotating disk flow introduced in section 1.5. The radial flow is larger in separation bubbles than in attached flow regions. This is in agreement with conclusions from rotating helicopter blades experiments [START_REF] Mccroskey | Measurements of boundary layer transition, separation and streamline direction on rotating blades[END_REF]. The radial velocity in separation can be analyzed similarly to the rotating disk flow. There are two reasons responsible for the larger velocity: firstly, the Coriolis force acts in the same direction as the centrifugal force when flow reverses; Secondly, the viscous effect is small since the separation bubble is thick.

However, the radial velocity on the rotating HAWT blade is not large enough to affect the flow transition significantly. It is shown that the transition on the rotating blade is caused by TS wave, which resembles the TS wave on the 2D airfoil.

There are several limitations in the current study as well. First, the region investigated does not include tip and root, where the 3D effect is stronger. The parameter chosen in this study (incoming flow and angular velocity) is near-optimal, so there is no massive separation. However, the rotation effect would be more significant when large flow separations exist. In a wind turbine farm, most wind turbines operate in the wakes of others. At the same time, wind turbines usually operate in atmospheric boundary layers. The incoming flow is rather complex and depends on the upstream turbines, topology of the terrain, and thermal effects. As a result, the transition observed in experiments [START_REF] Ha Madsen | Transition characteristics measured on a 2mw 80m diameter wind turbine rotor in comparison with transition data from wind tunnel measurements[END_REF] [58] happens quite near the leading edge and is likely triggered by the bypass transition instead of natural transitions. So another limitation of this study is the uniform inflow and therefore the lack of incoming free stream turbulence. Nevertheless, as it is shown that the effect of rotation on the boundary layer flow is quite small in the attached boundary layer, it is reasonable to believe that the bypass transition on 3D blade is similar to airfoil flow. This study might offer some justifications to study the transition on a HAWT blade using airfoil wind tunnel experiments, as long as the incoming flow condition is well controlled.

turbulent transition was measured in the context of investigating the cavitation inception in the original experiments. This propeller is chosen because of its large blade area ratio, and relatively small pressure gradient on the bulk of the blade. Therefore, cross-flow has a sufficient distance to develop from the leading edge to the trailing edge. The blade geometry data can be found in Kuiper's paper.

In the simulations, solely one phase flow is considered, i.e. there is no free surface or cavitation. The blade is in 3D Cartesian coordinates system (x, y, z) and does not move during the simulation. Its rotation is taken into account by solving the NS equations in the rotating reference frame. To take advantage of the circumferential periodicity of the propeller, the flow around a single blade is simulated. The shaft of the blade, therefore the incoming flow direction (or the ship advancing direction), is in z axis. The directrix of the blade is aligned with the y axis.

In the simulations, length scales are normalized by the radius of the blade R, and velocities are normalized by the incoming flow U z . The angular velocity of the blade is ω z = 4.33 radian per unit time. This gives an advance ratio J = U z /nD = 0.724, where the D = 2 is the propeller's diameter, and n = 0.69 is the number of revolutions per unit time. This J is the same as the Kuiper's original experiment. With the computational resources available, the Reynolds number is reduced to make a DNS on the whole blade to be possible. The kinematic viscosity ν of the fluid is set to be 1/80000. The diameter Reynolds' number Re N = nD 2 /ν based on this ν is 0.22 × 10 6 , which is one third of the experiment Re N in Ref [START_REF] Kuiper | Scale effects on propeller cavitation inception[END_REF]. Another widely-used Reynolds' number in the rotating disk studies is Re r = r ω z /ν. Re r linearly depends on the radius r and is usually interpreted as the position. In this study, Re r = 588 at the blade tip (r = R).

The Navier-Stokes equations in the rotating reference frame are solved using the open source code Nek5000 [START_REF] Fischer | Open source spectral element cfd solver[END_REF]. Nek5000 is a highly scalable CFD solver based on the Spectral Element Method (SEM), which is a hybrid method combining the geometry flexibility of Finite Element Method and the high accuracy of Spectral Method. It has been used in numbers of hydrodynamic stability studies [START_REF] Loiseau | Investigation of the roughness-induced transition: global stability analyses and direct numerical simulations[END_REF][START_REF] Appelquist | Global linear instability of the rotating-disk flow investigated through simulations[END_REF]. More information on the applications of SEM in fluid dynamic can be found in Refs [START_REF] Mo Deville | High-order methods for incompressible fluid flow[END_REF][START_REF] Karniadakis | Spectral/hp element methods for computational fluid dynamics[END_REF] et al.

Overlapping overset meshes

Propeller blade sections are usually defined on a series of concentric cylindrical surfaces. In the present case, we expand the cylindrical surfaces along the corresponding straight geodesics that are on the upper side of the plane x = 0. After expanding the cylindrical surfaces into z -x planes, the blade sections on these planes, which are simply airfoil profiles with a given pitch, form a new expanded blade. The mesh is generated on this expanded blade first, and then folded to obtain the final mesh. Although Nek5000 uses unstructured meshes, the elements have to be hexahedron.

Therefore the mesh for Nek5000 is usually generated as structured ones in the softwares. As shown in the figure 4.1, the pitch of the sections and periodic boundaries in the circumferential direction make it quite challenging to generate a mesh around the propeller. We decide to use the overlapping overset meshes technique in Nek5000 to However, as we will see, this mesh is sufficient to resolve the whole transition process.

On the blade surface, the non-slip BC is applied. The two boundaries in the cir- which is the summation of the unperturbed rotating velocity and the incoming velocity.

The uniform incoming flow U z = 1 means there is no inlet turbulence. This is in agreement of the original experiments of Kuiper, which is performed in a towing tank with the fluid at rest.

Definition of local flow directions

As mentioned above, velocities in Nek5000 are based on the Cartesian coordinate system and are represented by (U x , U y , U z ) in the solver. However, the cylindrical coordinate system is a better choice for the results analysis because of the nature of propellers.

After the velocity solutions are obtained, they are transformed into the cylindrical coordinate system which shares the same z axis with the Cartesian system used in the Because the LST is a local analysis method, proper local flow directions are also needed so that the wave-numbers and growth rates make sense. For general 3D flow, the local streamwise direction is usually defined as the direction of streamline in the inviscid region [START_REF] William S Saric | Stability and transition of three-dimensional boundary layers[END_REF]. However, this means that the streamwise and cross-flow directions experience spatial variation, which makes measuring wave-number and growth rate quite difficult along one direction.

For the present propeller case, it is natural to work on the cylindrical surfaces on which the blade sections are defined. Figure 4.4 shows the expanded propeller section on y = 0.6. If we neglect the geometry variation in the radial direction, the local flow directions definition is straightforward. Similar to 2D airfoil case: the local wall-normal direction is perpendicular to the wall and the streamwise is in wall tangential direction.

To obtain the streamwise velocity U s and the wall-normal velocity U n , U θ and U z need to be projected into the corresponding directions on the expanded plane. The local radial direction and therefore U r are the same as in the cylindrical coordinates system.

LNS simulation mesh

A body-fitted mesh is also generated for LNS simulation. Because we focus on the boundary layer instability, computational domains for LNS simulation needs not be so large as to enclose the whole blade. The structured LNS mesh consists of only one mesh block, where the grid points can be indexed by i, j, k, corresponding to the streamwise, wall-normal, and radial directions respectively.

This mesh is also generated firstly on the expanded blade and then folded. Figure 4.5 shows the first k -i slice in j direction of the final mesh, which is on the pressure side of the blade. The grid lines in j direction are aligned to the wall-normal direction defined in the section 4.1.3.

On each side of the blade, any point on the blade can be located by two cylindrical coordinates (r, θ). In this thesis, θ is measured with respect to y axis and is positive toward the trailing edge. On the blade, the gird points share the same k index have the same θ coordinate, the gird points share the same i index have the same r coordinate.

This facilitates the measurement of growth rates as well as wave-numbers.

The LNS mesh covers the region -7.6 • < θ < 32 • and 0.6 < r < 0.9 on the pressure side. The domain's height in the wall-normal direction is 0.06R. Totally 80 × 80×20 (radial, streamwise, and wall-normal directions respectively) elements are used.

The flow field on the LNS mesh is interpolated from the full mesh solution with spectral accuracy. The interpolated solution can be then used as the base flow U b in LNS. In the simulations, the homogeneous BCs are used on the four boundaries in radial and wall-normal directions, which means that the velocity disturbances are zero on them.

Outflow BC is applied to the streamwise boundary near the trailing edge. Velocity inlet BC is used on the streamwise boundary near the leading edge. The disturbance, which is introduced by velocity inlet BC, will be discussed later.

A similar mesh is also generated for the suction side of the blade for the analysis of DNS result. The mesh covers a slightly different region 0 • < θ < 28 • and 0.65 < r < 0.95.

Results

Overlapping overset meshes

To save computational resources, the solution on the outside mesh is obtained firstly with a coarse inside mesh (0.5 million elements), which does not resolve boundary layer flow. This permits using a large time step to get a solution for the inviscid flow region. After the outside solution is obtained, the inside mesh simulations are performed without marching the outside mesh solution.

Mesh convergence and velocity profiles

SEM is a type of h/p method. Mesh convergence study can be done either by using different meshes (h refinement), or the same mesh with different spectral orders (p refinement). For the inside mesh solution, two simulations with order 6 and order 8 are performed to make sure the solution is mesh independent. Totally around 3 million cpuhours on the Jean Zay supercomputer of IDRIS (France) are spent to get the converged solutions.

Thrust coefficient K t and torque coefficient K q are usually used to measure the propeller's performance. They are defined as:

K t = T /ρn 2 D 4 K q = Q/ρn 2 D 5 (4.1)
where T is the total thrust on the propeller and Q is the total torque.

The variation K t and K q with advance ratios J is given in figure 4.7. The values of the coefficients calculated from the present simulations are marked out by symbols. Despite the Re difference between the experiments and simulations, K t and K q from simulations are quite close to the experiment values with the same J. The relative differences between experiments and simulations are less than 10% for both K t and K q . The difference between the order 6 and order 8 results are so small that they are barely separable in the figure. This is smaller than that on the rotating disk (around 18% of inviscid circumferential velocity), but it is large enough to trigger cross-flow transition on the swept wing [START_REF] Wang | Control of crossflow instability over a swept wing using dielectric-barrier-discharge plasma actuators[END_REF].

At most locations, the cross-flow velocities are slightly negative in the inviscid region, which means that the flow is towards the blade root outside the boundary layer.

On the rotating disk, the boundary layer is homogeneous in the circumferential direction. However, there is no such infinite rotational symmetry for propellers. At the same radial position r, the boundary layer develops from the leading edge to the trailing edge. On the suction side, the boundary layer profiles exhibit similar characteristics to the pressure side. At θ = 0, r = 0.95, there is a significant disagreement between order 6 and order 8 results. The reason is that this location is close to the tip and the flow is close to transition. In the following section, the order 8 DNS result will be used exclusively.

Laminar-turbulent transition on the blade

The However, the turbulent tip vortex does not fully detach from the blade and it interacts with the boundary layer close to the tip.

On the blade surface, two different laminar to turbulent processes are observed:

1. At about 0.75 < r < 0.98, a typical natural transition can be observed.

2. At about 0.4 < r < 0.7, transition caused by the laminar separation bubble near the trailing edge.

The separation-induced transition near the trailing edge will be discussed in the section 4.2.4. In this section, we focus on the cross-flow transition. Between 0.7 < r < 0.95, there are numbers of clearly observable streamwise streaks on the blade. These Further downstream at θ = 26 • , the flow has a similar pattern as the previous position. However, as shown in figure 4.9, the cross-flow velocities are larger downstream.

As a result, the radial wave appears earlier than the previous θ position. The transition direction. As a result, Coriolis force changes its direction, which means that it now acts in the same direction as centrifugal force. On the other hand, the large thickness of the separation region means that the viscous force is relatively small. These two effects would cause a large radial velocity in the separation region on the propeller.

On other rotating blades like helicopter rotor and Horizontal Axis Wind Turbines, the large radial (cross-flow) velocities around flow separations are well-documented [START_REF] Mccroskey | Measurements of boundary layer transition, separation and streamline direction on rotating blades[END_REF][START_REF] Pk Chaviaropoulos | Investigating three-dimensional and rotational effects on wind turbine blades by means of a quasi-3d navier-stokes solver[END_REF]. This radial velocity would induce a Coriolis force component in the circumferential direction, which partially suppresses flow separations and leads to so-called rotational augmentation. However, few authors pointed out the relationship between the circumferential separations and the large radial velocities. propellers [START_REF] Sánchez-Caja | Scale effects on tip loaded propeller performance using a ranse solver[END_REF][START_REF] Bhattacharyya | Scale effects on open water characteristics of a controllable pitch propeller working within different duct designs[END_REF], where the high-quality images are available.

Flow transition's influence on the surface streamlines

For comparison, figure 4.17 gives the RANS results of skin friction and streamline.

The streamline in different flow regime regions have similar patterns as our result and the experiment. However, even though that the Re in the RANS simulation is the same as the experiment (i.e. it is three times larger than the present simulation), the turbulent region is restricted in the small tip region. This shows the necessity of including the proper transition model in RANS.

The surface streamlines show that the near-wall cross-flow velocity is much smaller in the laminar region than in the turbulent region. Kuiper [START_REF] Kuiper | Cavitation inception on ship propeller models[END_REF] attributed this to the wall stress difference in the two regions. For the laminar region, where the wall shear stress are weak, the centrifugal force dominates and tend to push to fluid outward. In contrast, the large wall stress in turbulent boundary layer would pull the streamline back to the circumferential direction. However, we argue that this explanation is only partially true.

It is necessary to see the phenomena in rotating reference frame. Close to the wall, the circumferential velocity of a turbulent boundary layer increases faster than a laminar one, so the Coriolis force is less decreased. As a result, the volume force which pushes the fluid outward is smaller in turbulent boundary layer. On the other hand, the large velocity gradient means that the viscous effect is strong. A small cross-flow velocity would produce a large viscous force in radial direction, making the fluid difficult to be pushed outward. In this sense, turbulent boundary layer is exactly opposite to the separated flow, where the cross-flow is unusually large.

In rotating disk boundary layer, the surface streamline direction experiences a similar change. In the laminar region, the angle between the surface streamlines and circumferential direction is around 40 • . Whereas this angle is around 12 • ∼ 17 • in turbulent region [START_REF] Appelquist | Turbulence in the rotating-disk boundary layer investigated through direct numerical simulations[END_REF][START_REF] Itoh | Turbulent boundary layer on a rotating disk in infinite quiescent fluid[END_REF].

Stability Analysis

Linear Stability Analysis

To further verify that the streamwise streaks observed in the DNS are indeed cross-flows vortices, linear stability analysis is performed. The prediction of LST is compared with the observed cross-flow vortices in terms of wave-numbers and wave-angle.

Both LST and LNS require a steady, unperturbed base flow. However, it is rather difficult to separate the total DNS field into a base flow field and a perturbation field because both of them are unknown and the perturbation are usually quite small. We chose to use part of the total DNS solution, where there is no significant presence of perturbations, as the base flow for LST and LNS. Most of the boundary layer flow on the suction side is affected by the existence of cross-flow waves. Whereas the flow on the pressure side is more 'laminar', which makes it more suitable as the base flow. In the following section, the base flow always refers to the flow on the pressure side of the blade. The radial wave-number of the most unstable stationary wave predicted by LST is around β r ≈ 3.1. Although a slightly different value of β r = 3.7 is observed in DNS, this wave is still quite unstable according to LST. Its growth rate predicted by LST is quite close to the most unstable wave. For the streamwise wavenumber, the difference between DNS result and LST prediction is quite small. From the wave-numbers in radial and streamwise directions, the wave angle φ = arctan(α i /β r ) can be calculated.

The local φ is around 8 • , which is much smaller than the rotating disk value of 14 • .

The reason is that the strength of cross-flow is smaller on the propeller than on rotating disks. As a result, the unstable wave is less convected in the radial direction.

Linearized Navier-Stokes simulation

Because the growth rate of the cross-flow vortices from DNS is hard to measure, we perform a LNS simulation to quantify the growth of the cross-flow wave on the blade.

The mesh and BCs used for LNS have been described in section 4. When the disturbance enters the domain through the inlet, an adjustment stage is usually needed. The large growth rate of LNS results at s < 10 is a result of this adjustment.

Further downstream, the LST and LNS results show a similar trend, but LST overestimates the growth rate because neither the streamwise curvature nor Coriolis force is included.

Convective/absolute nature of the instability

The above results demonstrate that flow transition on propeller blades share many similarities with rotating disk flow transition. Since it is well known that the rotating disk boundary layer flow is absolutely unstable when Re is large enough, a natural question to ask is that can the instability of the propeller blade boundary layer be absolutely unstable.

However, if the flow were locally absolutely unstable, it is would be hard to obtain a laminar solution by DNS. Any disturbance during the simulation would persist and lead to the transition. As Landau et al. [START_REF] Landau | Fluid mechanics[END_REF] put it, those (solutions) which do (exist in nature) must not only obey the equation of fluid dynamics but also be stable. At least on Nonetheless, without any known perturbation source, the cross-flow vortices appear on the pressure side and are not convected downstream. It is worth to make a simple investigation. A LNS simulation with the same settings as the previous one is performed.

However, the inlet disturbance lasts only a limited time t ≤ 0.064. When t > 0.064, the boundary condition is set to be homogeneous, so the source of disturbance disappears.

The U fields at different times are given in figure 4.22. The convective nature of the instability is clearly shown in the figure. After the initial pulse disturbance is removed, the wave-package keeps growing in its amplitude. At the same time, it is convected downstream until it goes out of the simulation domain.

Summary

The laminar-turbulent transition on marine propellers is a complex problem. The LNS simulation also shows that the boundary layer is absolutely unstable in the radial direction. However, since the periodicity is broken on the propeller, the unstable wave exhibits convective natures in the streamwise direction.

Chapter 5

Conclusion and Perspectives

Conclusion

In this thesis, we investigated the boundary layer transition of two typical rotating blade machineries: a HAWT and a marine propeller. Numerical simulations which fully resolve the transition process are performed with high order SEM method. The results provide some direct observations on the transition process on them.

The most prominent feature of HAWTs and propellers is the rotation. The effect of rotation on the boundary layer flow can be analyzed similarly to von Kármán swirling flow, where there is an exact solution. We argue that it is necessary to see the problems This is reflected in the surface streamline directions and should be explained through the relationship between the Coriolis force and the circumferential velocity.

Our results imply the transition on marine propeller can not be based on 2D airfoil tests or simulations, where the 3D effects are neglected. For RANS simulation of marine propeller, it is important to use a transition model that is capable to take cross-flow into account.

The difference between HAWT and marine propeller blades is likely caused by their shapes. On the rotating disk, the only relevant quantity is the Reynolds number Re = r ω/ν. However, on rotating blade, the distance to the leading edge is also important because the 3D effect need time to develop. For HAWT blades, which is very thin, the cross-flow velocity does not have an enough long distance to grow before leaving the blade from the trailing edge. As a result, the 3D effect due to the rotation is rather small.

Perspectives

Although we have identified the different transition scenarios on HAWT and propeller blades, more works could be done regarding the boundary layer transition on rotating blades. Following are several points that we believe is worth pursuing:

The influence of free stream turbulence

In this thesis, the inflow for both cases is uniform, which means that there is no free stream turbulence. This is not a realistic assumption for HAWT and propellers, which usually operate under complex inflow conditions.

The free stream turbulence provides a source of disturbance. When its intensity increase, the primary instability appears earlier and the transition finishes earlier. When the inflow turbulence is large enough, a different transition scenario (bypass transition) could happen.

For the HAWT blade, because the flow is very close to 2D airfoil flow, it is reasonable to believe that the bypass transition would be similar to 2D flow as well. However, for marine propellers, the 3D effect due to the rotation could play a role in the bypass transition.

The 3D separation bubble

We have shown that when flow separates on rotating blades, a much larger cross-flow velocity appears around the separation bubbles than in the attached boundary layer. In that case, the separation bubbles have essentially 3D structures.

As we mentioned, this 3D separation should be quite unstable. At present, almost all the studies on separation-induce transition are based on 2D LSBs (i.e. the streamlines of the flow is in 2D planes). There is seldom any study on the instability of 3D LSBs.

However, this kind of flow can be quite important from engineering perspectives (for instance, leading edge separation of marine propellers, which plays an important role in the cavitation inceptions) and therefore needs to be further studied.

The aeronautic propellers

There are other rotating blades machineries like aeronautic propellers which are not studied in the present thesis. The aspect ratios of aeronautic propeller are much smaller than HAWT blades. On the other hands, they are usually larger than marine propeller blades. Therefore it is rather difficult to reach any conclusion about its boundary layer transition based on the present results. Further works are necessary to get a better understanding.

model, the coefficient C s in equation 2.7 can be obtained using the Lilly contraction:

c s = M ij L ij M kl M kl (A.2)
where the tensor L ij and M ij can be obtained from the resolved flow. Their definition can be found in [START_REF] Stephen B Pope | Turbulent flows[END_REF][START_REF] Lilly | A proposed modification of the germano subgrid-scale closure method[END_REF].

A.1 Filtering Commonly-seen filtering include the box filter, the cut-off filter, and the Gaussian filter. For spectral method, the cut-off filter is a natural choice [START_REF] John | Two comments on filtering (artificial viscosity) for chebyshev and legendre spectral and spectral element methods[END_REF], where the physical quantity is transformed to the frequency domain. The high frequency components are cut off before being transformed back.

As mentioned before, the SEM can be switched back and forth from the nodal and the modal formulation. Supposing we have an one-dimension velocity u i on the GLL point x i (i = 0, ..., N ). Since the Legendre polynomials P j (x) form a complete orthog- [START_REF] Hugh | Spectral element filtering techniques for large eddy simulation with dynamic estimation[END_REF] the flow. Therefore the Lagrangian averaging idea does not apply.

Inspired by [START_REF] Chapelier | A spectral-element dynamic model for the large-eddy simulation of turbulent flows[END_REF], we spatially average the viscosity inside each element. The final viscosity field is obtained by interpolations. The details can be found in [START_REF] Chapelier | A spectral-element dynamic model for the large-eddy simulation of turbulent flows[END_REF]. However, we use a quadratic interpolation instead of linear one.

A.3 Test cases

Tests on flows such as the flat plate boundary layers show that the above LES model results some improvements compared with ILES (UDNS). However, for flows with long transitional regions, although the model correctly distinguishes the laminar and turbulent regions, the eddy viscosity is not zero in the transitional region, which is technically laminar. Due to the artificial viscosity, the equivalent Reynolds number is decreased.

Because this thesis is focused on the flow transition, we do not want any artificial viscosity in transitional regions, which can be quite long in certain cases. Therefore we chose to used UDNS (ILES) as the main tool. Nonetheless, we give part of the Exlicit LES results on two test cases.

1. 1

 1 (a) NM80 wind turbine in the DanAero MW project. (b) A marine propeller, from Wikipedia. . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Flow around a 2D airfoil. The incoming flow is in x direction. . . . . . 1.3 Blasius boundary layer on a flat plate. . . . . . . . . . . . . . . . . . . 1.4 Paths to transition with increasing environmental disturbances. In this thesis, we only consider the path A in the figure. From [15] . . . . . . . 1.5 Top view and side view of the T-S wave transition on the flat plate. From [16] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.6 Flow structures around a typical 2D separation bubble. From [19] . . . 1.7 LSB induced transition on a flat plate. The flow separation is caused by adversary pressure gradient. The contour is based on the normalized streamwise velocity. . . . . . . . . . . . . . . . . . . . . . . . . . . . .

15 1.

 15 15 (a) Boundary layer transition on the rotating disk. The spirals are crossflow vortices. (b) Close-up view of the cross-flow vortices. The highfrequency oscillations due to secondary instability can be clearly seen. From [42]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 1.16 Inside the separation bubble, streamwise velocity changes its direction (red arrow), which causes the Coriolis force (green arrow) acts in the same direction as the centrifugal force (black arrow). . . . . . . . . . . 19 1.17 Transition location on the center of a HAWT blade (where there is a glove). Visualized using thermal images. The bright color is laminar flow whereas dark color is turbulent. (a) Start-up stage. (b) Operating. From [58] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 1.18 An example of transition location (red line) derived from pressure fluctuations, which varies with time. From the DanAero MW project. For the meaning of the symbols in the figure, please refer to [62] . . . . . . 21 1.19 The Propeller C from Kuiper [65]. The green surface is one of the cylindrical surfaces on which the blade sections are defined. . . . . . . 22 1.20 Illustration of the suction side C p , flow structure and the cavitation on hydrofoil sections. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 1.21 (a)Sheet cavitation (b) Bubble cavitation. From [75] . . . . . . . . . . 24 2.1 (a) A traditional Chinese roly-poly toy (A-person-won't-fall). It is an example of a stable system. From taobao.com (b) The Gömböc, which is a homogeneous body having just one stable and one unstable point of equilibrium. From Wikipedia. . . . . . . . . . . . . . . . . . . . . . . 2.2 Plane Poiseuille flow. The height is normalized by the half wall distance. The velocities is normalized by the centerline velocity. . . . . . . 2.3 Example waves with different wave-vector. The inviscid flow is in x direction. (a) k = (1, 0) This wave pattern is usually seen in the T-S wave transition. (b) k = (1, 1) This pattern can be observed in the oblique transition. (c) k = (0, 1) This pattern is usually seen in the cross-flow transition. . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.4 Energy cascade in turbulent flow. k is the wave-number. The larger k is, the smaller the length scale. The green bars represent the resolved scale ranges and red bars modeled scale range in different simulation strategies. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.5 Solution of the model problem 2.8 at different times. ∆x = 0.005.

2. 6

 6 Variation of wave energy with time for different mesh sizes. . . . . . . 2.7 The global trail function φ i (x) on [0, 1]. The domain is divided into three elements. Inside each element, 4th order Lagrange polynomial is used to construct the piecewise global trail function. . . . . . . . . . . . 2.8 4th order Lagrange polynomial on GLL points, l i is shifted upwards by 2i. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.9 Spectral Element Method solution of the problem 2.9. Six elements, 15th order in each elements. . . . . . . . . . . . . . . . . . . . . . . . v LIST OF FIGURES 2.10 (a) An example of overlapping meshes. (This is just a illustration because there is no need to use two meshes for this simple geometry.) (b) The x velocity of Walsh's eddy solution at t = 0.04 obtained using the meshes in figure (a) . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 3.1 Variation of Re c , U r , U t along the span (Velocities are normalized by U ∞ ). 48 3.2 The first version of the mesh. The domain encloses the whole blade. However, this mesh is too large and lots of elements are wasted because we are only interested in the boundary layer. . . . . . . . . . . . . . . . 48 3.3 The mesh for the turbine blade . . . . . . . . . . . . . . . . . . . . . . 50 3.4 Boundary layer velocity profile at x/c l = 0.4 y = 25.24m(y/R = 0.63), (a) streamwise velocity (b) radial velocity . . . . . . . . . . . . . 52 3.5 Comparison of results obtained by different orders (a) pressure coefficient along Section8 (y/R = 0.75), normalized by local incoming flow, x coordinate is normalized by local chord length (b) velocity profile at x/c l = 0.5, y = 25.24m(y/R = 0.63) . . . . . . . . . . . . . . . . . . 53 3.6 The spectrum of the circumferential velocity signal (a) x/c l = 0.36, Section1, Suction side. (b) x/c l = 0.63, Section1, Suction side. . . . . . 53 3.7 Averaged velocity contour (a) 3D blade simulation, circumferential velocity on r = 25.24m (b) The same as (a), radial velocity (c) 2D simulation U x (c) 2D simulation, spanwise velocity U y . . . . . . . . . . . . 54 3.8 Comparison of pressure coefficients of 3D simulation, Airfoil simulation and wind tunnel measurement . . . . . . . . . . . . . . . . . . . . 55 3.9 Velocity profile. red-dashed line: 3D rotating blade; black-solid line: Airfoil (a) Streamwise velocity, Suction side (b) Streamwise velocity, Pressure side (c) Cross-flow velocity, Suction side (d) Cross-flow velocity, Pressure side . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 3.10 Coherent structure on suction side. Iso-contour of λ2, colored by circumferential velocity . (The blade is not to scale) . . . . . . . . . . . . 58 3.11 Skin friction contour. Normalized using U ∞ . . . . . . . . . . . . . . . 59 vi LIST OF FIGURES 3.12 Close up view of coherent structure on suction side (a) 3D rotating blade simulation (b) Airfoil simulation . . . . . . . . . . . . . . . . . . . . . 3.13 Spatio-temporal diagram of velocity signal (a) circumferential velocity from 3D simulation at y = 25.24m (the sensors are along Section1) (b) Airfoil simulation (the sensors are located along the mid-span section) .

3. 14

 14 Variation of growth rate of unstable wave with (a) ω r ( Normalized by 1ms. spanwise wave β r number is zero. The arrows mark out the observed ones in the simulations) (b) spanwise wave number β r (ω r = 0.4) 3.15 Coherent structure on the pressure side. . . . . . . . . . . . . . . . . . 4.1 The expanded blade section (in red) on the z -x plane at y = 0.8. Parts of the adjacent blades are also plotted. The blue lines are domain boundaries in the circumferential direction. The periodic points pairs have the same z coordinate. . . . . . . . . . . . . . . . . . . . . . . . . 4.2 The layout of overset overlapping mesh. The outside mesh is in black and the inside mesh is red.(a)z -x slice of the expanded meshes. (b) x -y slice of the final meshes. . . . . . . . . . . . . . . . . . . . . . . 4.3 Topology of the structured mesh blocks on the blade (the pressure side). The numbers are the element numbers along the corresponding edges. . 4.4 Definition of wall-normal and streamwise directions on the expanded z -x plane. The geometry non-uniformity in spanwise direction is neglected. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.5 First LNS mesh layer on the pressure side of the blade. . . . . . . . . . 4.6 U z contour on the expanded z -x slices at y = 0.5 and y = 0.8. The 2D streamlines are based on U θ and U z , i.e. the radial motion is neglected.

  results are quite close that they overlap with each other. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.8 Estimation of y+ on the suction side. . . . . . . . . . . . . . . . . . . . vii LIST OF FIGURES 4.9 Velocity profiles, the arrows indicate either increasing r or increasing θ. (a) Pressure side, r = 0.8, θ = -7.6 • , 7 • , and 21.2 • (b) Pressure side, θ = 0 • , r = 0.6, 0.75, and 0.9 (c) Suction side, r = 0.80, θ = 0.0 • , 10.8 • , and 21.5 • (d) Suction side, θ = 0 • , r = 0.65, 0.80, and 0.95 77 4.10 Iso-surface of lambda2 criterion (λ2 = -400) on the suction side, colored by velocity in z direction. . . . . . . . . . . . . . . . . . . . . . . 78 4.11 U z contour on the expanded section r = 0.95, which cut through the center of the small radial aligned structures. The boundary layer is fully attached. The only plausible explanation of the oscillation marked out is T-S waves. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78 4.12 U z contour at different θ locations (the suction side). d is distance-towall. These planes are j -k slices of LNS mesh. (a) θ = 0 • (b) θ = 13 • (c) θ = 26 • . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79 4.13 Iso-surface of lambda2 criterion on pressure side. . . . . . . . . . . . . 81 4.14 U z contour for different θ locations. (a) θ = -7 • (b) θ = 11 • (c) θ = 30 • 82 4.15 2D streamlines and the contour of radial velocity U r near the trailing edge on the expanded z -x plane y = 0.45. The radial movement of the streamlines are neglected. Because the blade is not uniform in spanwise, the wall-normal velocity is not exactly zero. As a result, some of the streamlines end on the wall. . . . . . . . . . . . . . . . . . 83 4.16 (a)Streamlines based on the current DNS. They are generated as close to the wall as possible. (b) Paint streaks from the experiment of [65] . . 84 4.17 Skin friction and streamlines obtained by RANS simulations. From [114]. 85 4.18 α i and α r of stationary cross-flow waves with different β r . The horizontal line marks the α r measured from the DNS. The vertical line marks the β r measured from the DNS. All the wave-numbers are normalized using the boundary layer scale 0.01R . . . . . . . . . . . . . . . . . . . 88 4.19 Inlet disturbance function (a) a(r) (b) b(d) . . . . . . . . . . . . . . . . 89 viii LIST OF FIGURES 4.20 LNS simulation (a) Iso-surfaces of U z . red: 2 × 10 -5 , blue: -2 × 10 -5 (b) U s along streamwise direction on the section r = 0.75. s is the local curve length to the inlet on the airfoil. Each red line represents one different wall-normal distance, the envelope of the wave is interpolated from the crests and troughs of the wave. . . . . . . . . . . . . . . . . . 4.21 Growth rates along streamwise direction at the section r = 0.75. For LST, the radial wave-number β r = 3.5, which is the same as the inlet disturbance of LNS simulation . . . . . . . . . . . . . . . . . . . . . . 4.22 Iso-contour of perturbed velocity at different time (a) t = 0.128 (b) t = 0.256(c) t = 0.384(d) t = 0.512 . . . . . . . . . . . . . . . . . . . A.1 An example of the Gaussian filter . . . . . . . . . . . . . . . . . . . . A.2 An example filter kernel function for N = 7 . . . . . . . . . . . . . . . A.3 Smagorinsky length along the channel height, the influence of elements boundaries can be clearly observed. From [120] . . . . . . . . . . . . . A.4 Illustration of flat plate flow set-up (not to scale) . . . . . . . . . . . . . A.5 Iso-surface of λ2, colored by distance to wall. (a)DNS (b)LES . . . . . A.6 Comparison of mean profile (solid lines) and RMS (dashed line) with that from [122](a) Re θ = 670, (b) Re θ = 1000 . . . . . . . . . . . . . . A.7 Contour of instantaneous eddy-viscosity from the LES model on a horizontal plane. The viscosity is 1/Re = 0.0016. . . . . . . . . . . . . . . A.8 Same as figure 3.5. . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.9 Same as figure 3.13. . . . . . . . . . . . . . . . . . . . . . . . . . . . A.10 (a) Flow structure on the suction side of the center HAWT blade. The dynamic Smagorinsky model is turned on. (b) The iso-surface of eddyviscosity obtained by the model. . . . . . . . . . . . . . . . . . . . ¸ais Les éoliennes à axe horizontal (HAWT) et les hélices marines sont deux pplications très étudiés dans le domaine de la recherche en mécanique des fluides. Cependant, la transition de couche limite sur celles-ci n'est pas encore entièrement comprise. Dans cette thèse, des Simulations Numériques Direct sont effectuées sur ces deux cas, à l'aide de supercalculateurs. Le principal l'objectif est d'étudier l'effet de la rotation sur la transition de couche limite. Les couches limites des HAWT et des hélices marines partagent un point commun avec l'écoulement de von Kármán, qui est créé par un disque rotatif immergé dans un fluide. L'écoulement de von Kármán peut induire ce que l'on appelle une 'cross flow' transition. L'objet de la présente étude est d'étudier la possibilité d'une transition d'écoulement transversal sur les HAWT et les hélices marines. Cette étude montre que la transition naturelle de couche limite sur les HAWTs et les hélices marine est induite par des mécanismes distinctement différents. Le résultat de l'écoulement autour d'un pale de HAWT montre que le profil de la couche limite est très proche d'un profil bidimensionel. Sur la pale, la vitesse dans le sens de l'envergure est faible lorsque la couche limite est attachée. En conséquence, la transition naturelle est très similaire au profil 2D et est due aux ondes de Tollmien-Schlichting (T-S).

  Les différences observée sur les écoulements de couche limite sur les HAWTs et les hélices marines sont probablement causées par leur degré de complexité géométrique xi The main objective is to study the effect of rotations on the boundary layer transition. The boundary layers on HAWTs and marine propellers share an apparent common point with von Kármán swirling flow, which is created by a rotating disk in the otherwise still fluid. The von Kármán swirling flow is the prototype of cross-flow transition. Therefore one focus of the present study is the possibility of the cross-flow transition on HAWTs and marine propellers. The present study shows that the natural boundary layer transitions on the HAWT and the marine propeller are induced by distinctively different mechanisms. The numerical result of a HAWT blade shows that the boundary layer profile on it is very close to the 2-Dimensional (2D) airfoil flow. On the blade, the velocity in spanwise direction is small in the attached boundary layer. As a result, the natural transition on HAWT blade is very similar to the 2D airfoil flow and is due to Tollmien-Schlichting (T-S) waves.On the marine propeller blade, the boundary layer flow is fully 3-Dimensional (3D) due to the rotation. The cross-flow instability and transition are clearly observed. The shapes of the cross-flow vortices are in good agreement with the prediction of Linear Stability Theory (LST). Although it has been long assumed that cross-flow transition should be important for propellers, this is the first direct observation of such phenomena as far as we know. Because the propeller does not have infinite rotational symmetry, our result suggests the boundary layer on the marine propeller is convectively unstable. This is different with the von Kármán boundary layer flow, which is absolutely unstable.

1 Background and Introduction 1 . 1

 11 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.1 Motivation of the present work . . . . . . . . . . . . . . . . . . 1.2 Aerodynamics of 2-Dimensional Airfoil . . . . . . . . . . . . . . . . . 1.3 The boundary layer and its transition . . . . . . . . . . . . . . . . . . . 1.3.1 The boundary layer . . . . . . . . . . . . . . . . . . . . . . . . 1.3.2 The laminar-turbulent transition . . . . . . . . . . . . . . . . . 1.3.3 T-S wave instability . . . . . . . . . . . . . . . . . . . . . . . . 1.3.4 Flow separation induced transition . . . . . . . . . . . . . . . . 1.3.5 Cross-flow transition . . . . . . . . . . . . . . . . . . . . . . . 1.3.6 Others . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4 Transition on the 2D airfoil . . . . . . . . . . . . . . . . . . . . . . . . 1.5 The rotating disk flow . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.5.1 Cross-flow velocity . . . . . . . . . . . . . . . . . . . . . . . . 1.5.2 Its transition . . . . . . . . . . . . . . . . . . . . . . . . . . .

Figure 1 . 1 :

 11 Figure 1.1: (a) NM80 wind turbine in the DanAero MW project. (b) A marine propeller, from Wikipedia.

Figure 1 . 2 :

 12 Figure 1.2: Flow around a 2D airfoil. The incoming flow is in x direction.

Figure 1 .

 1 Figure 1.3 shows the Blasius profile, which describes the boundary layer formed by flow over a flat plate. It is a self-similar solution in that the velocity profiles at different streamwise locations x is the same after being normalized. The boundary layer thickness of Blasius solution increases along streamwise direction according to

Figure 1 . 3 :

 13 Figure 1.3: Blasius boundary layer on a flat plate.

Figure 1 . 4 :

 14 Figure 1.4: Paths to transition with increasing environmental disturbances. In this thesis, we only consider the path A in the figure. From[START_REF] William S Saric | Stability and transition of three-dimensional boundary layers[END_REF]

Figure 1 . 5 showsFigure 1 . 5 :

 1515 Figure 1.5: Top view and side view of the T-S wave transition on the flat plate. From [16]

Figure 1 .

 1 Figure 1.6 depicts the flow structures around a 2D LSB. The adversary pressure gradient in stream-wise makes the initially wall-attached streamlines to separate from the wall.

Figure 1 . 6 :Figure 1 . 7 :

 1617 Figure 1.6: Flow structures around a typical 2D separation bubble. From[START_REF] Philip | Laminar separation bubbles in two and three dimensional incompressible flow[END_REF]

Figure 1 .

 1 Figure 1.8 shows the typical flow structure inside a 3D boundary layer. Inside the boundary layer, there is a velocity component which is perpendicular to the inviscid flow direction. This velocity profile is refereed to as 'cross-flow'.

Figure 1 .

 1 Figure1.9 shows the process of cross-flow transition. Similar to the T-S wave transition, the whole process can be roughly divided into linear, nonlinear (saturation), and breakdown stages. The small disturbance, which usually form a wave pattern in the cross-flow direction, firstly grows exponentially in the linear stage. As goes to downstream, the wave develops into streamwise vortices with larger amplitude. Both stationary and traveling cross-flow vortices can exist. When their amplitudes are around 1% ∼ 10% of the streamwise inviscid velocities, nonlinearity comes into effect. The vortices grow with smaller growth rates and eventually saturate. The saturated crossflow vortices can maintain a rather long distance. Their flow structures are usually referred to as "half mushroom" in contrast to the "mushroom" structures observed in Görtler instability.

Figure 1 .

 1 Figure 1.10: (a) The mushroom-like structure in Görtler instability. From [32]. (b) The half mushroom structure in cross-flow instability. From [33].

Figure 1 . 11 :

 111 Figure 1.11: Illustration of streamlines in inviscid region and in boundary layer on a rotating disk. Because the wall-normal velocity, which is parallel to the angular velocity vector, does not affect the direction and magnitude of the forces, therefore it is neglected in plotting these streamlines

Figure 1 . 12 :

 112 Figure 1.12: The arise of the cross-flow velocity in the rotating disk boundary layer.

Figure 1 . 13 :Figure 1 . 14 :

 113114 Figure1.13: Self-similar solution of the von Kármán boundary layer. From[START_REF] Lingwood | Absolute instability of the boundary layer on a rotating disk[END_REF]

Figure 1 .

 1 Figure 1.15: (a) Boundary layer transition on the rotating disk. The spirals are crossflow vortices. (b) Close-up view of the cross-flow vortices. The high-frequency oscillations due to secondary instability can be clearly seen. From [42].

CHAPTER 1 .

 1 BACKGROUND AND INTRODUCTIONCoriolis and curvature terms, which have a stabilizing effect, need to be included in linear stability analysis. The critical Reynolds' Re c number obtained by Linear Stability Theory (LST) is around 287, which is in agreement with their experiment (Re c = 294)

Figure 1 . 16 :

 116 Figure 1.16: Inside the separation bubble, streamwise velocity changes its direction (red arrow), which causes the Coriolis force (green arrow) acts in the same direction as the centrifugal force (black arrow).

Figure 1 . 17 :

 117 Figure 1.17: Transition location on the center of a HAWT blade (where there is a glove). Visualized using thermal images. The bright color is laminar flow whereas dark color is turbulent. (a) Start-up stage. (b) Operating. From[START_REF] T Reichstein | Investigation of laminarturbulent transition on a rotating wind-turbine blade of multimegawatt class with thermography and microphone array[END_REF]

Figure 1 . 18 :

 118 Figure 1.18: An example of transition location (red line) derived from pressure fluctuations, which varies with time. From the DanAero MW project. For the meaning of the symbols in the figure, please refer to[START_REF] Bak C Døssing | High frequency microphone measurements for detection of transition on airfoils[END_REF]

Figure 1 . 19 :

 119 Figure 1.19: The Propeller C from Kuiper [65]. The green surface is one of the cylindrical surfaces on which the blade sections are defined.

Figure 1 .

 1 Figure 1.20 shows the typical C p distribution on the suction side of a blade section.

Figure 1 . 20 :

 120 Figure 1.20: Illustration of the suction side C p , flow structure and the cavitation on hydrofoil sections.

Figure 1 .

 1 Figure 1.21: (a)Sheet cavitation (b) Bubble cavitation. From[START_REF] Van Terwisga | Tu delft: Opencourseware: Cavitation on ship propellers[END_REF]

Figure 2 . 1 :

 21 Figure 2.1: (a) A traditional Chinese roly-poly toy (A-person-won't-fall). It is an example of a stable system. From taobao.com (b) The Gömböc, which is a homogeneous body having just one stable and one unstable point of equilibrium. From Wikipedia.

Figure 2 . 2 :

 22 Figure 2.2: Plane Poiseuille flow. The height is normalized by the half wall distance. The velocities is normalized by the centerline velocity.

Figure 2 . 3 :

 23 Figure 2.3: Example waves with different wave-vector. The inviscid flow is in x direction. (a) k = (1, 0) This wave pattern is usually seen in the T-S wave transition. (b) k = (1, 1) This pattern can be observed in the oblique transition. (c) k = (0, 1) This pattern is usually seen in the cross-flow transition.

Figure 2 . 4 :

 24 Figure 2.4: Energy cascade in turbulent flow. k is the wave-number. The larger k is, the smaller the length scale. The green bars represent the resolved scale ranges and red bars modeled scale range in different simulation strategies.

8 )

 8 Equation 2.8 is a pure convection problem. The initial condition v 0 = sin 40 (πx) creates a localized wave-shape. The exact solution of this equation is v = v 0 (x -t), which always has the same wave shape as the initial condition. Because the domain is periodic at the boundaries, any wave exits the domain from the right boundary would re-enter the domain from the left.

Figure 2 .

 2 Figure 2.5 plots the numerical solution of eq. 2.8 obtained by finite difference method.As can be seen, the amplitude of the wave is damped as it is being convected.

Figure 2 . 5 :

 25 Figure 2.5: Solution of the model problem 2.8 at different times. ∆x = 0.005. Time is discretized using Forward Euler method. Spatial derivative is approximated using first order upwind finite difference.

Figure 2 . 6 :

 26 Figure 2.6: Variation of wave energy with time for different mesh sizes.

  ) is the advection matrix, c(t + dt) is the solution vector.Although the trial function and the test function are defined on the global domain, they are constructed locally.

Figure 2 .

 2 [START_REF] Guilmineau | Numerical simulations of the cavitating and non-cavitating flow around the postdam propeller test case[END_REF] shows the trial functions on the global domain [0, 1] when 3 elements are used. It can be seen that the trial functions φ i (x)s are piecewise functions and each φ i (x) is nonzero in only one or two elements. Inside each element, Lagrange polynomials l n are used to construct the global trial functions. Once the trial functions are available, the mass matrix A and advection matrix B can be evaluated by Gauss quadrature rule of integration. In fact, A and B are evaluated locally in each element and then assembled to get the final ones.

Figure 2 . 7 :

 27 Figure 2.7: The global trail function φ i (x) on [0, 1]. The domain is divided into three elements. Inside each element, 4th order Lagrange polynomial is used to construct the piecewise global trail function.

Figure 2 .

 2 Figure 2.9 plots the SEM solution of the advection problem 2.9. Compared with

Figure 2 . 8 :

 28 Figure 2.8: 4th order Lagrange polynomial on GLL points, l i is shifted upwards by 2i.

Figure 2 . 9 :

 29 Figure 2.9: Spectral Element Method solution of the problem 2.9. Six elements, 15th order in each elements.

1 .Figure 2 .

 12 Figure 2.10: (a) An example of overlapping meshes. (This is just a illustration because there is no need to use two meshes for this simple geometry.) (b) The x velocity of Walsh's eddy solution at t = 0.04 obtained using the meshes in figure (a) .

Figure 2 .

 2 Figure 2.10a shows an example overlapping meshes. The total domain is a square with periodic Boundary Conditions (BCs) in both directions. It is filled by two partially overlapping meshes. The outside mesh is in black and the inside in red. During the simulation, there are two sessions running simultaneously. Onefor each mesh. For the boundaries located in the another mesh, velocity BC is used and the value specified is interpolated from the another session. For the implementation details, please refer to[START_REF] Mittal | Nonconforming schwarz-spectral element methods for incompressible flow[END_REF].

become 2 . 15 .

 215 Compared with the NS equation in the inertial reference frame, two additional terms arise because of rotation. The centrifugal force ω × (ω × r) depends on the ω and r. The Coriolis force 2ω × U depends on ω as well as the local flow velocity U .

Figure 3 . 1 :

 31 Figure 3.1: Variation of Re c , U r , U t along the span (Velocities are normalized by U ∞ ).

Figure 3 . 2 :

 32 Figure 3.2: The first version of the mesh. The domain encloses the whole blade. However, this mesh is too large and lots of elements are wasted because we are only interested in the boundary layer.

Figure 3 .

 3 Figure 3.3 shows the computational domain and mesh. As can be seen in figure 3.3a, the mesh in span direction is concentrated on the center of the blade. The coarse mesh near the root and tip should dissipate any large nonphysical flow structures induced by the inaccurate BCs. Local mesh refinement in wall-normal direction for the boundary layer is shown in figure 3.3b. There are 230 elements in spanwise direction, 25 elements in wall-normal direction, and 78 elements around the airfoil in the streamwise direction.

  which are the sum of unperturbed rotation and wind velocities. A smaller domain simulation shows that the boundary conditions at two span ends have negligible influence on the flow in the center region.

Figure 3 . 3 :

 33 Figure 3.3: The mesh for the turbine blade

Figure 3 . 4 :

 34 Figure 3.4: Boundary layer velocity profile at x/c l = 0.4 y = 25.24m(y/R = 0.63), (a) streamwise velocity (b) radial velocity

Figure 3 .Figure 3 . 5 :Figure 3 . 6 :

 33536 Figure 3.5b compares the velocity profile at the blade center. The pressure side velocity profiles are laminar flow and both simulations give almost the same result. On suction side, where the flow is turbulent, the results do not match each other. Although the boundary layer's thickness and edge velocity are similar, it seems that order 8 so-

Figure 3 . 7 :Figure 3 . 7

 3737 Figure 3.7: Averaged velocity contour (a) 3D blade simulation, circumferential velocity on r = 25.24m (b) The same as (a), radial velocity (c) 2D simulation U x (c) 2D simulation, spanwise velocity U y

Figure 3 . 8 :

 38 Figure 3.8: Comparison of pressure coefficients of 3D simulation, Airfoil simulation and wind tunnel measurement

Figure 3 .

 3 Figure 3.10 gives an overall view of coherent structures at the suction side of the blade.

Figure 3 . 10 :

 310 Figure 3.10: Coherent structure on suction side. Iso-contour of λ2, colored by circumferential velocity . (The blade is not to scale)

Figure 3 . 11 :

 311 Figure 3.11: Skin friction contour. Normalized using U ∞

Figure 3 .

 3 Figure 3.12 gives a close-up view of the coherent structures at the suction side of the blade mid-span as well as the airfoil case. The spanwise vortices are formed by TS waves. 3D and 2D simulations show similar characteristics. In both simulations,

Figure 3 .

 3 Figure 3.13 gives the spatio-temporal distribution of streamwise velocity signals recorded during the simulations. Again, the 2D and 3D simulations are similar. The sensors are inside boundary layers and have a roughly constant distance to the wall. As a result, the velocity signal is large near the leading edge since the boundary layer is thin. It decreases until x/c l = 0.3 as the boundary layer grows thicker and thicker. The

Figure 3 . 12 :

 312 Figure 3.12: Close up view of coherent structure on suction side (a) 3D rotating blade simulation (b) Airfoil simulation

Figure 3 .Figure 3 . 13 :

 3313 Figure 3.14 shows linear stability analysis results. The base flow used is at x/c = 0.3 on suction side of the 2D and 3D simulations. The Coriolis and centrifugal terms are not included in the linear stability equation. Figure 3.14a shows the variation of growth rate of a 2D wave (spanwise wave number β r = 0) with frequency. The unstable wave found is obviously TS wave instead of cross-flow mode since β r = 0. The frequency observed in simulations are marked with arrows. Although there is a difference between the growth rates between 2D and 3D simulations, they follow the same

Figure 3 . 14 :

 314 Figure 3.14: Variation of growth rate of unstable wave with (a) ω r ( Normalized by 1ms. spanwise wave β r number is zero. The arrows mark out the observed ones in the simulations) (b) spanwise wave number β r (ω r = 0.4)

Figure 3 . 15 :

 315 Figure 3.15: Coherent structure on the pressure side.

Figure 4 . 1 :

 41 Figure 4.1: The expanded blade section (in red) on the z -x plane at y = 0.8. Parts of the adjacent blades are also plotted. The blue lines are domain boundaries in the circumferential direction. The periodic points pairs have the same z coordinate.

Figure 4 .

 4 Figure 4.1 shows the expanded blade section on the z -x plane at y = 0.8. The local incoming flow U t is the sum of the incoming (or ship advancing) velocity U z and the velocity U θ arising from the propeller rotation. Because marine propellers rotate fast, a large pitch of the blade section is necessary to ensure the angle of attacks is appropriate in operation conditions.

Figure 4 . 2 :

 42 Figure 4.2: The layout of overset overlapping mesh. The outside mesh is in black and the inside mesh is red.(a)z -x slice of the expanded meshes. (b) x -y slice of the final meshes.

Figure 4 .

 4 Figure 4.2 shows two slices of the layout of the inside and outside meshes. The blade is mounted on an infinitely long shaft, whose radius is 0.313R. The radial range of the computational domain is from the shaft to 3R, but the inside mesh only ranges from 0.313R to 1.1R. The full domain in z (the shaft axis) direction extends to 1.8R upstream and 3.2R downstream. The outside mesh has totally 0.22 million elements. 0.95 million elements are used in the inside mesh.Figure 4.3 shows the layout of the

Figure 4 .

 4 [START_REF] Torres | Statistical characteristics of interacting wind turbine wakes from a 7-month lidar measurement campaign[END_REF] shows the layout of the mesh blocks on the pressure side of the blade. There are 25 elements in the wall-normal direction, among which at least 10 ∼ 12 are inside the boundary layer in average. We do not intend to perform a DNS which resolves the Kolmogorov scale of the turbulence.

Figure 4 . 3 :

 43 Figure 4.3: Topology of the structured mesh blocks on the blade (the pressure side). The numbers are the element numbers along the corresponding edges.

Figure 4 . 4 :

 44 Figure 4.4: Definition of wall-normal and streamwise directions on the expanded z -x plane. The geometry non-uniformity in spanwise direction is neglected.

Figure 4 . 5 :

 45 Figure 4.5: First LNS mesh layer on the pressure side of the blade.

Figure 4 .

 4 Figure 4.6 shows contour of U z on two different z -x slices. Before expanding, half of the domain is rotated 90 degrees towards the another half. As a result, the pressure and suction sides of the blade form a passage similar to a compressor stator passage.The cyclic periodic boundaries meet each other in the center of the passage. In figure4.6, the velocity contours are continuous across the meshes and the periodic boundary, which is in the center of the passage. This means that not only the cyclic periodic BC works properly, but also the overlapping overset meshes technique.

Figure 4 . 6 :

 46 Figure 4.6: U z contour on the expanded z -x slices at y = 0.5 and y = 0.8. The 2D streamlines are based on U θ and U z , i.e. the radial motion is neglected.

Figure 4 . 7 :

 47 Figure 4.7:Experiment measurements of K t and K q for different advance ratios J (from Ref.[START_REF] Kuiper | Scale effects on propeller cavitation inception[END_REF]). The symbols are from current simulations. The order 6 (red) and order 8 (black) results are quite close that they overlap with each other.

Figure 4 . 8 :

 48 Figure 4.8: Estimation of y+ on the suction side.

Figure 4 .

 4 Figure 4.8 gives the wall unit y+ contour based on the order 8 mesh. Because the grid spacings are not uniform across one element, the y+ is underestimated slightly in the figure. As we will see, the flow is mostly laminar on the blade. The friction velocity is therefore small. As a result, the y+ is less than 0.05 in most regions. Even for the turbulent region, the maximum y+ is still less than 0.2.

Figure 4 .

 4 Figure 4.9 shows the velocity profiles from the two simulations at different locations. Both streamwise and cross-flow velocities are plotted against distance-to-wall d.The order 6 and order 8 results agree quite well with each other. The difference between them is barely distinguishable at most positions.

Figure 4 .

 4 9a shows the velocity profiles at the same r but with different θ on the pressure side. The larger θ is, the closer to trailing edge is the position. The thickness of the boundary layer increases as the position goes downstream. The strength of cross-flow velocity increases as well towards the trailing edge.The pressure side boundary layer velocities at the same circumferential location θ but different radial locations are plotted in figure4.9b. As r increases, the inviscid streamwise velocity increases, and the boundary thickness decreases. This is also different from the rotating disk flow, where the boundary layer thickness is a constant along the radial direction. The strength of cross-flow decreases because the location with a large r is closer to the leading edge.

  instantaneous coherent structures on the suction sides of the blade are visualized by the lambda2 criterion in figure 4.10. The flow regimes on the blade are different depending on the radius r. Firstly, there is a clearly visible tip vortex on the top edge of the blade. It is laminar close to the leading edge, where it incepts at first. As it goes to downstream, the laminar tip vortex rapidly becomes turbulent. The turbulent tip vortex then slightly detaches from the blade, making the shade difference in the figure.

Figure 4 . 9 :

 49 Figure 4.9: Velocity profiles, the arrows indicate either increasing r or increasing θ. (a) Pressure side, r = 0.8, θ = -7.6 • , 7 • , and 21.2 • (b) Pressure side, θ = 0 • , r = 0.6, 0.75, and 0.9 (c) Suction side, r = 0.80, θ = 0.0 • , 10.8 • , and 21.5 • (d) Suction side, θ = 0 • , r = 0.65, 0.80, and 0.95

Figure 4 . 10 :

 410 Figure 4.10: Iso-surface of lambda2 criterion (λ2 = -400) on the suction side, colored by velocity in z direction.

Figure 4 .

 4 Figure 4.11: U z contour on the expanded section r = 0.95, which cut through the center of the small radial aligned structures. The boundary layer is fully attached. The only plausible explanation of the oscillation marked out is T-S waves.

Figure 4 . 12 :

 412 Figure 4.12: U z contour at different θ locations (the suction side). d is distance-to-wall. These planes are j -k slices of LNS mesh. (a) θ = 0 • (b) θ = 13 • (c) θ = 26 •

Figure 4 .

 4 Figure 4.12 plot the U z contour on several wall-normal slices along the radial direction at different θ locations. At θ = 0 • , the flow is all laminar in the shown range. There is no visible waves in radial direction when r < 0.9. Close to the tip 0.90 < r < 0.94, a short wave in radial direction can be observed, which increases the boundary layer thickness.Downstream at θ = 13 • , there is a clear distinction between the laminar and turbulent regions. The flow is fully laminar for r < 0.85. When r increases, the boundary layer experience periodic oscillations with larger and larger amplitudes. After about 3 wavelengths, the boundary layer becomes turbulent flow. The transition process finishes completely at r = 0.92.

Figure 4 . 13 :

 413 Figure 4.13: Iso-surface of lambda2 criterion on pressure side.

Figure 4 .- 7 •Figure 4 . 14 :

 47414 Figure 4.[START_REF] Schlichting | Boundary-layer theory[END_REF] shows the flow structures on the pressure side of the blade. Although several cross-flow vortices can be observed near the trailing edge at r > 0.8R, they are restricted to a smaller region than suction side. The flow transition does not fully complete until the vortices leave the blade. Closer to the root, the flow separation is very weak compared the suction side. This is more or less expected because the pressure gradient on the pressure side is usually more favorable than the suction side.

Figure 4 .

 4 Figure 4.15: 2D streamlines and the contour of radial velocity U r near the trailing edge on the expanded z -x plane y = 0.45. The radial movement of the streamlines are neglected. Because the blade is not uniform in spanwise, the wall-normal velocity is not exactly zero. As a result, some of the streamlines end on the wall.

Figure 4 .Figure 4 .

 44 Figure 4.15 shows the 2D streamlines and the cross-flow velocity contour around the trailing edge of the expanded section r = 0.45. On the suction side, the directions of streamlines indicate that the flow separates and form a Laminar Separation Bubble(LSB). When reverse flows inside the bubble is large enough, the boundary layer can be absolutely unstable[START_REF] Rodríguez | The two classes of primary modal instability in laminar separation bubbles[END_REF][START_REF] Embacher | Direct numerical simulations of laminar separation bubbles: investigation of absolute instability and active flow control of transition to turbulence[END_REF] and the flow would quickly become turbulent. Because the instability of LSB is similar to the Kelvin-Helmholtz instability, large vortex shedding structures are also observed in figure4.10.

Figure 4 .

 4 Figure 4.16a shows the surface streamlines based on the current DNS. It can be seen that they point to different directions depending on the flow regimes. Around the separation region, the flow is almost totally in the radial direction due to the large radial velocity discussed in the section 4.2.4. In the laminar region, the streamlines are deviated significantly outwards. After the flow become turbulent, they are pulled back to the

Figure 4 . 17 :

 417 Figure 4.17: Skin friction and streamlines obtained by RANS simulations. From [114].

Figure 4 .

 4 Figure 4.18 shows the linear stability analysis result at r = 0.80, θ = 11 • . At this location, the disturbance in the DNS result is quite small compared to downstream, but its amplitude is large enough for measuring the wave-numbers. In the LST, only the stationary cross-flow waves are considered. Variation of the growth rate α i and the streamwise wave-number α r are plotted as a function radial wave-number β r . The wave-numbers of the vortices observed in DNS are also given in the figure. Since the waves are almost aligned to streamwise direction in DNS, the α r is small and it is hard to measure a local value. The α r marked out in the figure is calculated from two streamwise wavelengths.

1 . 4 .

 14 Figure 4.18: α i and α r of stationary cross-flow waves with different β r . The horizontal line marks the α r measured from the DNS. The vertical line marks the β r measured from the DNS. All the wave-numbers are normalized using the boundary layer scale 0.01R

(4. 2 )Figure 4 .

 24 Figure 4.20a shows the U z iso-contours of the converged LNS field. The inlet disturbance develops into cross-flow streaks downstream. The shape of the streaks matches the DNS results. They orient roughly in circumferential direction but deviate slightly outward in the radial direction. As moving downstream, the disturbance propagates in both outward and inward radial directions.

Figure 4 . 19 :

 419 Figure 4.19: Inlet disturbance function (a) a(r) (b) b(d)

Figure 4 . 20 :

 420 Figure 4.20: LNS simulation (a) Iso-surfaces of U z . red: 2 × 10 -5 , blue: -2 × 10 -5 (b) U s along streamwise direction on the section r = 0.75. s is the local curve length to the inlet on the airfoil. Each red line represents one different wall-normal distance, the envelope of the wave is interpolated from the crests and troughs of the wave.

Figure 4 . 21 :Figure 4 .

 4214 Figure 4.21: Growth rates along streamwise direction at the section r = 0.75. For LST, the radial wave-number β r = 3.5, which is the same as the inlet disturbance of LNS simulation

Figure 4 . 22 :

 422 Figure 4.22: Iso-contour of perturbed velocity at different time (a) t = 0.128 (b) t = 0.256(c) t = 0.384(d) t = 0.512

  The boundary layer flow often interacts with cavitation and different transition scenarios could coexist. In this thesis, DNSs are performed to identify potential mechanisms that are responsible for the boundary layer transition on propellers. The computational domain used follows the pitch direction of the blade. Overset overlapping meshes are used to ensure cyclic boundary conditions. The transitional boundary layer is full-resolved in this study.As far as we know, this is a first direct observation of cross-flow vortices on propellers. It is show that the radial velocity in the attached boundary layer is large enough to trigger cross-flow instability. The cross-flow vortices observed in the present results are similar to that on the rotating disk but with depths in z axis. LST can predict the shape of the unstable cross-flow waves but overestimates the growth rate of the cross-flow wave without Coriolis force and streamwise curvature included. Our results show that cross-flow transition is an important transition route on the rotating propeller blades, as has been expected previously.It is also argued that it is necessary to see the boundary layer flows under the rotating reference frame. Close to the blade root and trailing edge, the transition location is very close to the flow separation region. An exceptionally large cross-flow presents around the separation points, which is essentially caused by the same reason as the cross-flow in the attached flow region. However, when flow separates, the Coriolis force in the rotating reference frame acts outwards in the radial direction, causing the large radial velocity. Opposite to the separation region, the cross-flow velocity in the turbulent region is small, which causing the direction deviation when streamlines go into the turbulent region from the laminar region.

  in rotating reference frame. The cross-flow (spanwise or radial velocity) can be seen as a secondary flow induced by the primary flow in the circumferential direction. The Coriolis force relates the circumferential velocity with the cross-flow velocity. The numerical results show this argument is not only valid for attached laminar boundary layer, but also in separation bubbles and turbulent boundary layers.The HAWT results show that the boundary layer around a rotating HAWT blade is very close to the corresponding 2D airfoil case. Their boundary layer velocity profiles are almost identical in the attached laminar flow region. As a result, the natural transition on the 3D blade is induced by T-S wave. The transition location is almost the same for the blade and the 2D airfoil. This suggests that the 2D wind tunnel test can be served as an approximation of the flow around a blade section.On the marine propeller blade, the boundary layer transition is induced by the crossflow. As far as we know, this is the first direct observation of cross-flow vortices on marine propeller. The shape of the cross-flow vortices agrees well with the prediction of linear stability theory. For different flow regimes, the flow directions are different.

Figure A. 1 :

 1 Figure A.1: An example of the Gaussian filter

Figure A. 3 :

 3 Figure A.3: Smagorinsky length along the channel height, the influence of elements boundaries can be clearly observed. From[START_REF] Hugh | Spectral element filtering techniques for large eddy simulation with dynamic estimation[END_REF]

  HAWTs blades boundary layer . . . . . . . . . . . . . . . . . . 18 1.6.2 Cavitation and marine propeller boundary layer . . . . . . . . . 21 1.7 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 Hydrodynamic instability . . . . . . . . . . . . . . . . . . . . . . . . . 27 2.1.1 The Linearized NS equations . . . . . . . . . . . . . . . . . . . 28 2.1.2 Linear Stability Theory . . . . . . . . . . . . . . . . . . . . . . 29 2.1.3 Further notes on LST . . . . . . . . . . . . . . . . . . . . . . . 31 Flow fields . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2.4 Boundary layer instability and transition . . . . . . . . . . . . . 3.3 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Methodology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1.1 Blade geometry and parameters . . . . . . . . . . . . . . . . . 4.1.2 Overlapping overset meshes . . . . . . . . . . . . . . . . . . . 4.1.3 Definition of local flow directions . . . . . . . . . . . . . . . . The wind turbine blade . . . . . . . . . . . . . . . . . . . . . . 106

	2 Methodology 2.1 4 Marine propeller 4.1 Bibliography Chapter 1	27 109
	Background and Introduction
	1.1 Background
	de Nantes) is the wind energy. The specific topics under active studies in the lab include
	the aerodynamic of blade sections [1, 2], turbine wakes [3], as well as floating wind
	turbines etc.
	xviii

1.5.3 Implication . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.6 HAWTs and marine propellers . . . . . . . . . . . . . . . . . . . . . . xv CONTENTS 1.6.1 2.2 Comparison of different simulation strategies for transitional flows . . . 32 2.2.1 RANS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 2.2.2 LES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 2.2.3 DNS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 2.2.4 Choice of the numerical method . . . . . . . . . . . . . . . . . 37 2.3 Spectral Element Method . . . . . . . . . . . . . . . . . . . . . . . . . 38 2.3.1 Spectral Method . . . . . . . . . . . . . . . . . . . . . . . . . 38 2.3.2 Spectral Element Method . . . . . . . . . . . . . . . . . . . . . 39 2.3.3 Nek5000 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 2.3.4 NS equations in the rotating reference frame . . . . . . . . . . 45 2.3.5 IDRIS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 2.4 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 3 Rotating HAWT blade 47 3.1 Set-ups of the simulations . . . . . . . . . . . . . . . . . . . . . . . . . 47 3.1.1 The HAWT blade . . . . . . . . . . . . . . . . . . . . . . . . . 47 3.1.2 The meshes and boundary conditions . . . . . . . . . . . . . . 49 3.1.3 The 2D airfoil simulation . . . . . . . . . . . . . . . . . . . . . 50 3.2 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50 3.2.1 The influence of BCs at two ends . . . . . . . . . . . . . . . . 51 3.2.2 Mesh Convergence . . . . . . . . . . . . . . . . . . . . . . . . 52 xvi CONTENTS 3.2.3 4.1.4 LNS simulation mesh . . . . . . . . . . . . . . . . . . . . . . . 4.2 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2.1 Overlapping overset meshes . . . . . . . . . . . . . . . . . . . 4.2.2 Mesh convergence and velocity profiles . . . . . . . . . . . . . 4.2.3 Laminar-turbulent transition on the blade . . . . . . . . . . . . 4.2.4 Separation-induced transition and the influence of rotation . . . 4.2.5 Flow transition's influence on the surface streamlines . . . . . . 4.3 Stability Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.1 Linear Stability Analysis . . . . . . . . . . . . . . . . . . . . . 4.3.2 Linearized Navier-Stokes simulation . . . . . . . . . . . . . . . 4.3.3 Convective/absolute nature of the instability . . . . . . . . . . . 4.4 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 Conclusion and Perspectives 5.1 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.2 Perspectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A The Explicit LES model in the Nek5000 A.1 Filtering . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.2 Averaging of the eddy viscosity . . . . . . . . . . . . . . . . . . . . . . A.3 Test cases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.3.1 The flat plate . . . . . . . . . . . . . . . . . . . . . . . . . . . xvii CONTENTS A.3.2 As a source of renewable energy, the installation capacity of wind turbines has been increasing steadily around the world. It is projected that this trend would last at least for the foreseeable future. What's more, the recent development of offshore wind turbines offers a new promising possibility. As one of the fundamental problems, the aerodynamic of wind turbines and wind farms has attracted many researchers. As indicates by its name, one of the many research fields of the LHEEA (Laboratory in Hydrodynamics, Energetics and Atmospheric Environment) in ECN ( École Centrale

(

  Institut du Développement et des Ressources en Informatique Scientifique) located in Paris, France. It has 1528 HPE XA730i nodes and the peak power can reach 4,89 PFlop/s. It also has 1292 NVIDIA GPUs. At the time of writing, IDRIS ranks 54th on the top500 list of the fastest computer in the world.

Together with free-stream turbulence, the boundary layer transition happens earlier in the experiment. After the transition, the flow re-attaches and there is no separation.

Marine propeller

Despite the apparent similarities between the rotating disk and propellers, there are very few investigations on propeller boundary layers instability and transition. The works reviewed in Chapter 1 are more concerned with the transition location. Although it has been pointed out that cross-flow instability is relevant for the rotating blade boundary layer transition [START_REF] Moran-Guerrero | On the influence of transition modeling and crossflow effects on open water propeller simulations[END_REF][START_REF] Schuelein | Transition detection and skin friction measurements on rotating propeller blades[END_REF], to the present authors best knowledge, there is no concrete numerical or experimental evidence that shows cross-flow instability exists on marine propellers. This chapter aims to provide a numerical simulation that captures the whole laminar-turbulent transition process on a marine propeller and to find the similarities and differences between the boundary layers on the rotating disk and propellers.

The structure of this chapter is similar to the last one. We will first give the descriptions of the numerical set-ups. Then the results from DNS and stability analysis are given. We also address the problem of convective/absolute instability. At last, it ends with a summary.

Methodology

Blade geometry and parameters

The propeller investigated in the present study is the Propeller C from Kuiper [START_REF] Kuiper | Scale effects on propeller cavitation inception[END_REF],

which was tested in the towing tank at MARIN institute in the late 1970s. The laminar to However, an interesting feature is that the radial range of the wave-package keeps increasing as being convected downstream. The disturbance propagates both outward and inward radial directions, which means that the flow in absolutely unstable in the radial direction. Similar results are obtained on the FSC flows and swept wing boundary layer flows [START_REF] Taylor | The long-time behaviour of incompressible sweptwing boundary layers subject to impulsive forcing[END_REF][START_REF] Rj Lingwood | On the impulse response for swept boundary-layer flows[END_REF], where the flows are absolutely unstable in one direction while are convective unstable in another direction.

As pointed out in Ref [START_REF] Bertolotti | Absolute/convective instability investigation of primary and secondary crossflow vortices[END_REF], because of the periodicity in the circumferential direction of rotating disk flow, only the radial group velocity has to vanish for an absolute instability to exist. In another word, on the rotating disk, although an unstable wavepackage is convected in the circumferential direction, it will return to the location where it begins with a larger amplitude as long as it is not convected in the radial direction.

Obviously propellers do not have this periodicity, as a result, unstable waves are convected downstream and eventually detach from the blade.

If the same definition of Re = r ω z /ν as the rotating disk flow is used, the Re range of the base flow in LNS is between 353 and 529. However, in rotating disk flow, the absolute instability does not appear until Re > 510, which is much higher than the current propeller blade case. This difference is likely caused by the different cross-flow velocities in the inviscid region. In rotating disk flow, the cross-flow in the inviscid region is zero. However, as have been seen in the figure 4.9, the cross-flow in the inviscid region is negative in the blade case, which means that the fluid is flowing inward in the inviscid region. This provides a mechanism for the unstable wave to propagate toward the blade root.

The above study on the convective/absolute nature of the instability is only a very first attempt and not conclusive. An interesting observation in figure 4.13 is that the cross-flow vortices are very close to the trailing edge. Similar observation have been made near the edge of the rotating disk in Yim et. al. [START_REF] Yim | Transition to turbulence in the rotating disk boundary layer of a rotor-stator cavity[END_REF] The superposition of so called 'elephant global modes', which exists near the boundary of the local convective instability and the local absolute instability, forms a front.

Appendix A The Explicit LES model in the Nek5000

This appendix addresses some specific problems when implementing the classic Smagorinsky LES model in the Spectral Element Method. Therefore the theoretical basis of the LES and the derivation of governing equations are only introduced very briefly. For the details, please refer to the textbook of Pope [START_REF] Stephen B Pope | Turbulent flows[END_REF] and the references given below.

As mentioned in the section 2.2.2, the traditional LES solves the filtered NS equations:

where over-bar represents the resolved flow quantities using the given meshes (which can be viewed as a filter). The subgrid stress τ ij need to be modeled using the resolved quantities.

To account for the dissipation due to the unresolved small scale motions, artificial viscosity is needed in modeling the subgrid stress τ ij . In the dynamic Smagorinsky onal basis, u i can be expanded as a series:

where u * j is the coefficient corresponding to jth Legendre polynomial P j (x).

Before transforming back, the coefficient u * j is multiplied by a kernel function l j , which is 1 for small i and 0 for large i. The high-frequency oscillation is removed in this way. Figure A.2 shows an example l i for N = 7.

In the implementation of the filtering operation, the matrix-vector multiplication is used to simplify the above calculation. The equation A.3 can be written in the following matrix form [START_REF] Hugh | Spectral element filtering techniques for large eddy simulation with dynamic estimation[END_REF]:

Where B ij = P j (x i ). So:

Define matrix L = diag(l 0 , ..., l N ). Equation A.4 can be written as:

In order to conserve the continuity across the elements boundary in the filtered flow field, the basis function in the above expansion is based on following polynomials instead of P j (x).

φ i = P j -P j-2 (f or i greater than 1)

φ i equals 0 at two boundaries when i > 1. Therefore the filtering using φ i does not affect the value at boundaries. If the original flow field is continuous, the filtered one is continuous as well.

A.2 Averaging of the eddy viscosity

The above model is available in the turbChannel example case in Nek5000. In the application of dynamic Smagorinsky model, it is found that C s varies volatilely across the domain. Therefore some averaging techniques are need to keep simulations stable.

For flow with homogeneous direction, the averaging can be done across that direction.

This method is used in the turbChannel case.

To extend the LES model in the turbChannel case to general 3D flows, another averaging technique other than averaging along homogeneous directions is need. One brilliant idea is from [START_REF] Meneveau | A lagrangian dynamic subgrid-scale model of turbulence[END_REF], where the dynamic eddy viscosity is time-averaged in a Lagrangian frame. In another word, the averaging is along flow tracelines. A weighting function is used to give more weights to the near position.

However, this averaging method is not suitable for the SEM. The filtering operation introduced above is not unif orm or equal across the elements. As a result, the final C s exhibits some mesh-dependent patterns as shown in A.3, which are not convected with

A.3.1 The flat plate

The first case is a transitional flat plate boundary layer flow with zero pressure gradient.

The computational domain is 1000δ 0 ×120δ 0 ×62.24δ 0 in x, y and z directions, where the δ 0 is the displacement thickness at the inlet x 0 . The Reynolds number Re δ 0 based on δ 0 is 650.

A.3.2 The wind turbine blade

The second test case of the Explicit LES is the wind turbine blade. The mesh and parameters are the same as in the chapter 3. The order is 8. Although its been long assumed that cross-flow transition should be important for propellers, this is the first direct observation of such phenomena as far as we know. Because the propeller does not have infinite rotational symmetry, our result suggests the boundary layer on the marine propeller is convectively unstable. This is different with von Kármán boundary layer flow, which is absolutely unstable.