Dynamic stress assessment in high head Francis runners

(1)

MASTER’S THESIS

2008:016 CIV

Yasir Dawood Siwani

Dynamic stress assessment in high head Francis runners

MASTER OF SCIENCE PROGRAMME Engineering Physics

Luleå University of Technology

Department of Applied Physics and Mechanical Engineering Division of Fluid mechanics

(2)

Dynamic stress assessment in high head Francis runners

Yasir Dawood Siwani

Master thesis

Department of Applied Physics and Mechanical Engineering Luleå University of Technology

SE-971 87 Luleå Sweden

January 25, 2008

(3)

Preface

This work was conducted under the supervision of Associate Professor Michel Cervantes, at the department of fluid mechanics at luleå university of technology, Sweden and Trond Moltubakk technical leader for HHT mechanical Norway at Rainpower Norway AS.

The work was preformed during 2007/08-12 at HHT mechanical Noway, Rainpower Norway AS in Kjeller, Norway.

There are some people I would like to acknowledge for their help and support during the thesis period. First I would like to thank Trond moltubakk and Michel Cervantes for the time given me both within and outside the thesis, their guidance, patient, criticism and suggestions throughout the work.

I would also like to thank Ståle Risberg, hydraulic design engineer at Rainpower Norway AS, for his help and guidance with the CFD calculations preformed and for the many discussions not only on the hydraulic part of the thesis but also on the mechanical calculations. Ståle Risberg was responsible for obtaining the CFD results to be used for further

(4)

ii

analysis in the thesis.

Special thank to Øystein Gjerde, Mechanical design engineer at Rain- power Norway AS, for his patient and discussions within the mechanical work in the thesis.

Yasir Dawood Siwani January 25, 2008

(5)

Abstract

The deregulation of the energy market in the recent years has lead to an increased willingness to operate hydro turbines over a larger operational range and with increased number of start and stops. This means that the turbine will be subject to a more frequent and larger load variation.

Especially for high head turbines, the impact of dynamic forces has been of considerable interest. High head turbine runners are subjected to dynamic forces originating from various sources. One significant source is the interaction between the runner and guide vanes, also called rotor - stator interaction (RSI). RSI creates periodic pressure fluctuations in the runner. Combined with stress concentration and material defects, the in- duced pressure fluctuations can give an increased risk of fatigue cracks.

For Francis runners, the areas near the trailing edge of the runner blade towards the band and/or crown has been identified as critical areas for fatigue, and for high head runners there has been several examples where fatigue failure have been initiated near the trailing edge toward the band.

The thesis presents a validation of dynamic stress obtained via computational fluid dynamics and the finite element method. A new runner was installed at Tonstad hydropower plant in Norway. Site measurements

(6)

iv

were carried out on the new replacement runner and one of the original design runners. The tests were conducted by Halvard Bjørndal in 2003 [4]. The measurements were partially done to compare the mechanical properties of the new design to the original one, during operation [4], and also for validation of FEA stress prediction.

Measurement results on the new runner were used for validation of the finite element method. A computational fluid dynamic geometry model of the new runner was modeled. Geometrical symmetry allowed for a third of the real runner including guide vanes and an extension of the outlet to be modeled.

Unsteady CFD simulations providing the load for FEA calculations were preformed by Ståle Risberg a hydraulic design engineer at Rainpower Norway AS. The transient CFD solution showed as predicted a fluctuat- ing torque on the blades due to the runner blade - guide vane interaction (or Rotor Stator Interaction (RSI)). The calculations were carried out on 81% load, which gives a mean torque of about 3.3 M N m which is an output of 129.6 M P a.

The CFD analysis resulted in a mean torque of 3.2 M N m and at the guide vane frequency an amplitude of 76.3 kN m. The result yielded a mean torque deviation of 1.8% between the CFD torque and the mean torque calculated based on a nominal output value of 160 M P a.

The dynamic pressure for one cycle, one blade passing a hole guide vane passage, is the dynamic pressure load used in the FEA. The FEA model consisted of a portion of the band and crown, a runner blade and a splitter blade. The mean torque obtained by the FEA was found to be 3.21 M P a. That amounted to a deviation of only 2.71% between the FEA torque and the mean torque calculated for 81% load.

The maximum principle quasi-static mean stress calculated was found to be located at the trailing edge of the blade near the band. The high- est mean stress, at the predefined validation points, was located 15 mm from the trailing edge of the blade on the pressure side. Compared to the measured value, the calculated stress was 32.7% higher. The deviation was 33.8%, 9% and 1% for 25 mm, 200 mm and 300 mm respectively on the pressure side. Looking at the suction side, the mean stress calculated

(7)

v

deviated by 121% at 15 mm and 7.5% 25 mm from the trailing edge.

The stress amplitudes at the guide vane passing frequency (150 Hz) were much lower than the measured values, for four out of the validation points. On the pressure side, the deviation was 67.9%, 63.2%, 32.1%

and 13.3% at 15 mm, 25 mm, 200 mm and 300 mm respectively. The amplitudes on the suction side deviated by 94.3% 15 mm and 71.1%

25 mm from the trailing edge of the blade.

Although the deviation was higher for the amplitudes than for the mean stress, the tendencies were similar. The deviation seamed to decrease for distance further from the trailing edge of the blade. It was concluded that one cause could be, thickness difference between the design model and the real one.

(8)

(9)

List of Figures

1.1 Simplified schematic of a Hydropower plant. . . 3 1.2 Turbine types and there operational head span. The im-

age is gathered from Krivchenko 1994 [16]. . . 4 1.3 Figure showing illustration of a Kaplan turbine and a

schematic of an adjustable blade Kaplan turbine. (The images are printed with permission from Rainpower Nor- way AS) . . . 5 1.4 Illustration of a Pelton turbine.(The images are printed

with permission from Rainpower Norway AS) . . . 6 1.5 Simplified schematic of a Pelton turbine.(The image is an

illustration printed with permission from Rainpower Nor- way AS) . . . 7 1.6 Illustration of a Francis turbine. (The image is an illus-

tration printed with permission from Rainpower Norway AS) . . . 8

(12)

x LIST OF FIGURES 1.7 Simplified schematic of a low head Francis turbine with

terminology. 1. Stayring vanes, 2. Guide vanes, 3. Wicket gate lower ring, 4. Turbine cover, 5. Support bearings of guide vane upper pivot, 6. Attachment of cover of stayring upper band, 7. Shaft flange, 8. Shaft, 9. Runner crown, 10. Runner band, 11. Runner blades, 12, 13, 14. Oper- ating gear, 15- Servomotor, 16. Bearing, 17. Generator thrust bearing support, 18. Runner cone, 19. Runner band seal, 20. Blanching hole. The image is gathered from Krivchenko 1994 [16]. . . 9 1.8 Different designs for Francis runners [2]. a.High head, b.

Medium head and c. Low head. (The image is an illustration printed with permission from Rainpower Norway AS) . . . 10 1.9 The figure shows the location of the strain gauges on the

new runner, blade 7. Picture (a) shows the locations of the gauges on the pressure side. Picture (b) shows the locations on the suction side. The images are prints from the site measurement report [4]. . . 11 1.10 The formation of the rosette on blade 14 of the original

runner. The figure also illustrates the error estimation shown to the right. The image is a print from the site measurement report [4]. . . 12 2.1 Figure showing a bar with a diameter D, length L, cross

section area after deformation is A, the material has the Youngs modulus E. The force F applied on the bar has resulted in an elongation ². b) Shows the nominal stress σ at the cross section. . . . 18 2.2 The figure shows a standard stress strain curve obtained

by tensile testing. The figure to the right shows the elastic contra the plastic region in the curve. On the figure to the left, σ_t is the tensile strength, σ_y is the yield strength, ²_t is the tensile strain and ²_f is the fracture strain. . . 20 2.3 The figure shows the same cross section as in figure 2.1b

but here a defect in the form of a vacant sphere. . . 20 2.4 The figure shows the stress concentration condition for a

small hole with diameter d in an infinitely large and thin plate. . . 21

(13)

LIST OF FIGURES xi 2.5 Figure showing simplistic view of the methodology of the

discretization process. The image is a modified one that originally was found in Felippa 2002 [10]. . . 24 2.6 Figure illustrating the shapes, degrees and node positions

for some elements. . . 28 3.1 Figure of the designed runner model. (a) Showing a trig

view, (b) a side view,(c) a top view and (d) a bottom view. 32 3.2 Figure showing the CFD model used in the analysis. It

represents a third of the total geometry shown in figure 3.1. 33 3.3 Figure of the CFD Mesh. (a) Showing a top view, (b) a

zoom in, for a closer look at the mesh. . . 34 3.4 Figure showing one of the simplifications on the design

model in order to save CPU time. (b) The zoom in design model contour on the band, (c) a zoom in, for a look at simplification done at the area. . . 35 3.5 Figure showing the finite element geometry of model one

used in the analysis. Containing a splitter blade, blade and part of the crown and band. . . 36 3.6 Figure showing the finite element geometry of model two

used in the analysis. a) showing a back view and b) a front. 37 3.7 Figure showing the type of solid mesh element used in

the analysis. It is a ten node tetrahedra element with quadratic displacement behavior. The image is one gathered from the Ansys 11 manual pages [1]. . . 38 3.8 Figure showing the type of surface mesh element used in

the analysis. It is a Second order elements with four (on the right) to eight (on the left) nodes. The image is one gathered from the Ansys 11 manual pages [1]. . . 38 3.9 Figure showing the finite element mesh of model two used

in the analysis. The figure is showing a pressure side view. 38 3.10 Figure showing the boundary conditions for the finite el-

ement analysis. . . 40

(14)

xii LIST OF FIGURES 4.1 Figure showing the results obtained by the CFD calcula-

tion. The figure on the top, figure a, shows the pressure results normalized with respect to the head of the power plant and figure b shows a contour plot of the velocity at the same positions. The velocity is normalized with respect to the the maximum velocity v_max. . . 42 4.2 The figure is showing the pressure results on the surfaces

of the blades. It should be noted the the contour plot presented here is taken at on certain time step out of the transient CFD calculation. . . 43 4.3 Figure showing the results obtained by the FEA calcula-

tion. The figure on the top, figure a, shows the maximum principle stress for a given time step. The plotted contour is of σ/σ_max. For a closer look of the legend values se figures given in the appendix, figure A.3-A.5 . . . 45 4.4 The figure shows the quasi static stress variation on the

pressure side at different point, at various locations from the trailing edge of the blade. . . 46 4.5 The figure shows the quasi static stress variation on the

suction side at two different locations from the trailing edge of the blade. . . 47 4.6 A plot of the mean stress at different locations from the

trailing edge on the blade. The mean stress calculated from the site measurements is here compared to the FEA calculated results. . . 48 4.7 A figure showing plots of the mean stress deviation as a

function of the distance from the trailing edge. . . 49 4.8 A plot of the mean stress at different locations from the

trailing edge on the blade. The mean stress calculated from the site measurements is here compared to the FEA calculated results. . . 50 4.9 A figure showing plots of the stress amplitude deviation

as a function of the distance from the blade trailing edge. 51

(15)

5.1 The figure illustrates the geometry difference at the trailing edge of the blade between the FEA model and the real runner. The figure to the left represents the shape profile of the FEA model and the one on the right is that of the real runner. It should be noted that the figure is only an illustration of the difference and is in no way an exact representation of the geometry. . . 55 A.1 The figure shows the location of the strain gauges on blade

7. The Picture shows the locations of the gauges on the pressure side. The image is a print from the site measurement report[4]. . . 62 A.2 The figure shows the location of the strain gauges on blade

7. The Picture shows the locations on the suction side.

The image is a print from the site measurement report[4]. 63 A.3 Figure showing the results obtained by the FEA calcula-

tion. The figure shows the maximum principle stress for a given time step. The plotted contour is of σ/σ_max. . . . 65 A.4 Figure showing the results obtained by the FEA calcula-

tion. The figure shows a zoom in for a better picture of the trailing edge at the pressure side. The contour plot gives the distribution of the maximum principle stresses for a given time step. The plotted contour is of σ/σ_max. . 66 A.5 Figure showing the results obtained by the FEA calcula-

tion. The figure shows a zoom in for a better picture of the trailing edge at the suction side. The contour plot gives the distribution of the maximum principle stresses for a given time step. The plotted contour is of σ/σ_max. . 67 A.6 Figure showing the results obtained by the FEA calcula-

tion. The contour plot gives the distribution of the maximum principle stresses for a given time step. And the streamlines color code responds to the local velocity of the fluid. The figure was only meant to serv as an illustrative picture of the interaction between fluid and structure and nothing else. . . 68

(16)

xiv LIST OF FIGURES

(17)

List of Tables

B.1 The table presents the normalized mean stress measured.

The stresses are normalized with respect to the maximum measured mean stress. . . 70 B.2 The table presents the normalized stress amplitudes at the

guide vane passing frequency 150Hz. The amplitudes are normalized with respect to the maximum measured stress amplitude. . . 71 B.3 The table presents the normalized FEA principle stresses.

The stresses are normalized with respect to the maximum measured mean stress on site. σ denotes the mean stress and ^∆σ₂ stands for the stress amplitude. It shoould be noted that the ^∆σ₂ was normalized with respect to the maximum measured stress amplitude. . . 73

(18)

CHAPTER 1. INTRODUCTION

Chapter

1 Introduction

The mechanical energy in flowing water has been exploited for centuries.

Waterwheels driven by the flow of water in rivers were used in ancient India mainly for irrigation of crops. They were used for grinding grain by the ancient Greeks and the Romans were known for using waterwheels for mining. Although useful for many purposes, waterwheels are not very effective in harnessing the full potential of the kinematic energy of flowing water.

The first development of a more modern turbine was done in the 8th century by professor Ján Andrej Segner and is today called the Segner turbine. The unit was an axial machine, much smaller than a water- wheel, that raised the limit for operating heads and flow rate [21]. The development of water turbines took off during the industrial revolution and in 2005 the existing capacity of hydroelectric power was 816 GW which was about 21% of the worlds electricity supply at that time [18].

(19)

1.1. BACKGROUND CHAPTER 1. INTRODUCTION

1.1 Background

The mechanical integrity of turbine runners is put to its very limits with increasing demand for higher efficiency to a lower cost. Designers today have access to design tool that enable them to push the mechanical limits further than ever before. Computational methods have given designers a deeper understanding of the complex turbulent flow and a greater knowledge of the loads they impose on mechanical structures.

Turbines are pushed to be more available and made to operate at a wider rang of operation loads, which means frequent start and stop cycles and a drift at of best efficiency operational conditions. Since turbines are manufactured to work at best efficiency for a longer period of time, they experience loads they are not dimensioned for and thus the risk of fatigue failure is increased. Designers are facet with fatigue, fatigue loading, corrosion, cavitation, corrosion fatigue, stress fatigue and many other challenges when designing a new runner.

Turbine runners are subjected to high dynamic loading that could eventually lead to fatigue failure. Dynamic fatigue loading cycles are of major concern. They are manly categorized in two sections, High Fatigue Cy- cles (HCF) and Low Fatigue Cycles (LCF). LCF originates from Start stop cycles and HCF from fatigue loading during operation. The main source for HCF is thought to be the rotor stator interaction [12].

Fatigue failure generally occurs in areas of high stress concentration where dynamic stress amplitudes are the largest. Structural loading fatigue failure in runners has been observed both in early life and after several years of operation. Fatigue cracking propagate from areas of high stress concentration [12]. For high head Francis runners, the transitional welded area between the blade and band or crown have been identified as critical areas [8]. The risk for fatigue failure is specially high near the trailing edge of the blade. Structural failure originate from various different factors but the major causes are thought to be dynamic hydraulic load fluctuations combined with material defects, especially defects within the weld [8].

With computational methods such as Computational Fluid Dynamics (CFD) and Finite Element Analysis (FEA) we are able to better un-

(20)

CHAPTER 1. INTRODUCTION 1.2. HYDRAULIC MACHINES derstand the cause and that way improve the mechanical properties of future runners.

1.2 Hydraulic machines

There are many different types of hydraulic machines in use at different power and pump stations worldwide. The choice of a machine depends on the application and the installation site. They are chosen mainly based on the head, headwater height over the tailwater, and flow rate at the specific reservoir. Figure 1.1 shows a simplified schematic of a Hy- dropower plant. The most commonly used turbines are Kaplan, Francis and Pelton turbines. The reason why these types are the most predominate today [12], is that they supplement each other. They cover most operation heads and discharge water flows, figure 1.2 shows an approximate head span for different turbine types [16].

Figure 1.1: Simplified schematic of a Hydropower plant.

(21)

1.2. HYDRAULIC MACHINES CHAPTER 1. INTRODUCTION

Figure 1.2: Turbine types and there operational head span. The image is gathered from Krivchenko 1994 [16].

1.2.1 Kaplan turbines

Another classification for Kaplan turbines is axial flow turbine, this is because the fluid flows axially through the runner. It is a low-head machine used for heads ranging approximately from 10 to 70 meters and high flow rates [16]. Figure 1.3(b) shows a schematic of a Kaplan turbine.

As the fluid flows through the runner the static pressure decreases and mechanical energy in transferred to the runner driving the shaft. Kaplan runners with adjustable blades have a larger range of discharges and are used at sites where there is a high head variation [16].

(22)

CHAPTER 1. INTRODUCTION 1.2. HYDRAULIC MACHINES

(a) Illustration of a Kaplan turbine (b) 1. Runner, 2. Guide vane, 3. Outer cover, 4. Lower cover, 5. Runner chamber, 6.

Stay ring, 7. Spiral casing, 8. Vane arm, 9.

Link, 10. Regulating ring, 11. Inner cover, 12. Turbine shaft, 13. Turnťbine bearing, 14. Shaft seal box, 15. Thrust bearing, 16.

Draft tube cone, 17. Drain pump, 18. Inter- mediate shaft, 19. Servo motor.

Figure 1.3: Figure showing illustration of a Kaplan turbine and a schematic of an adjustable blade Kaplan turbine. (The images are printed with permission from Rainpower Norway AS)

(23)

1.2. HYDRAULIC MACHINES CHAPTER 1. INTRODUCTION 1.2.2 Pelton turbines

Figure 1.4: Illustration of a Pelton turbine.(The images are printed with permission from Rainpower Norway AS)

Pelton turbines are used manly for extremely high heads with low discharge [16]. A Pelton turbine operates by one or more nozzles im- pinging jets onto buckets placed circumferentially on the runner disk, se figure 1.5. The buckets change the direction of the flow and the resulting change in momentum exerts a force on the buckets, thus driving the runner [12]. The resulting change in momentum is why Pelton turbines are also classified as impulse machines. The runner is mounted vertically or horizontally and both runner and nozzle are located over the tailwater [16].

(24)

Figure 1.5: Simplified schematic of a Pelton turbine.(The image is an illustration printed with permission from Rainpower Norway AS)

For impulse turbines the static pressure is constant over the runner.

In the case of a reaction turbines the pressure decreases as the fluid flows through the runner. Francis and Kaplan are examples of reaction turbines [12].

(25)

1.2. HYDRAULIC MACHINES CHAPTER 1. INTRODUCTION 1.2.3 Francis turbines

Figure 1.6: Illustration of a Francis turbine. (The image is an illustration printed with permission from Rainpower Norway AS)

Francis turbines are also classified as radial-axial flow turbines. The fluid flows radially, with respect to the shaft, in at the inlet of the runner and exits axially at the outlet. Francis turbines are use for heads typically ranging from 20 to 900m [16]. A simplified scheme of a Francis turbine can be seen in figure 1.7. Francis runners are as the case for Ka- plan runners submerged in water. Newtons third law of motion describes the energy transfer from fluid to runner for these types of turbines, se below [5].

Newtonťs third law:

When two particles interact, the force on one particle is equal and oppo- site to the force on the other.

The water is lead from the reservoir via a pressure conduit, the pen- stock, to the spiral casing. The spiral casing is designed to distribute a uniform flow in to the runner. The fluid passes the stay vanes and goes through the guide vane channels onto the runner blades. Kinematic

(26)

Figure 1.7: Simplified schematic of a low head Francis turbine with terminology. 1. Stayring vanes, 2. Guide vanes, 3. Wicket gate lower ring, 4. Turbine cover, 5. Support bearings of guide vane upper pivot, 6.

Attachment of cover of stayring upper band, 7. Shaft flange, 8. Shaft, 9. Runner crown, 10. Runner band, 11. Runner blades, 12, 13, 14. Op- erating gear, 15- Servomotor, 16. Bearing, 17. Generator thrust bearing support, 18. Runner cone, 19. Runner band seal, 20. Blanching hole.

The image is gathered from Krivchenko 1994 [16].

energy is then absorbed from the fluid by the blades driving the shaft connected to the generator. The water is then lead through the draft tube to the tailwater.

Francis runner

Francis runner designs vary depending on the head. The designs are distinguished between low, medium and high head turbines [2], se firure 1.8. A high head Francis runner may also consist of splitter blades as well as blades to obtain the maximum efficiency possible for higher heads.

(27)

1.3. SITE MEASUREMENTS CHAPTER 1. INTRODUCTION

Figure 1.8: Different designs for Francis runners [2]. a.High head, b.

Medium head and c. Low head. (The image is an illustration printed with permission from Rainpower Norway AS)

1.3 Site measurements

A new unit was installed at Tonstad hydropower plat in Norway. The unit has a nominal head of 430 m, a nominal power of 160 M P a and a rotational speed of 375 rpm. The runner replaced one of the original runners and had a new design.

Measurements where carried out to investigate the difference in mechanical properties between the new and the original design. Strain gauges were mounted near the trailing edge of the blade close to the band, on unit I (the original runner) and unit II (the new runner). The strain gauges were mounted on blade 7 and 14 both unites. The strain gauges mounted at blade seven were installed as a precaution, in case the gauges on blade 14 would destroy during the test [4].

Six strain gauges were installed on each blade, four on the pressure side and two on suction side. Gauges were installed at 15 mm (P S15), 25 mm (P S25), 200 mm (P S200) and 300 mm (P S300) on the pressure side, for the new runner, and 15 mm and 25 mm on the suction side, se figure 1.9. All the strain gauges were mounted at a distance of 25 mm from the band. The strain gauges installed on the original runner had the same positions except for P S200 which was mounted 175 mm from the

(28)

CHAPTER 1. INTRODUCTION 1.3. SITE MEASUREMENTS trailing edge to make room for a rosette strain gauge [4].

(a) (b)

Figure 1.9: The figure shows the location of the strain gauges on the new runner, blade 7. Picture (a) shows the locations of the gauges on the pressure side. Picture (b) shows the locations on the suction side.

The images are prints from the site measurement report [4].

The major principle stress, σ₁, was assumed to have a direction perpendicular to the band. The transverse stress was assumed to be insignifi- cant. That introduced an error which needed to be investigated [4].

The error of measuring the strain in a single direction was estimated by installing a rosette on the original runner. The rosette was installed according to figure 1.10 and was placed on the pressure side 200 mm from the trailing edge on blade 14 on the original runner. As seen to the right in figure 1.10, an error of 0.14% was found after the test [4].

Site measurements on the new runner were used for validating the dynamic stress investigation in the thesis, therefore site measurement results presented in the thesis are from measurements on the new unit.

The measurement for that test was not fulfilled according to the planed test. During the first try, the data logger used, only manege to register results for the first 30 minutes of a 60 minute test session. Other test tries were unsuccessful in improving the first try. The first test mea- surements covered start up and uploading to 290 mm servomotor stroke (81% load), full load is at 343 mm servo stroke [4].

(29)

1.3. SITE MEASUREMENTS CHAPTER 1. INTRODUCTION

Figure 1.10: The formation of the rosette on blade 14 of the original runner. The figure also illustrates the error estimation shown to the right. The image is a print from the site measurement report [4].

The maximum mean stress at 290 mm servomotor stroke was registered to be located at P S15. The mean stresses gradually decreased as the distance from the trailing edge increased. The largest fluctuations oc- curred during startup, at speeds near the natural frequency of one of the blades. That seemed to have triggered interference between mechanical and hydraulic resonance [4].

(30)

CHAPTER 2. THEORY

Chapter

2 Theory

2.1 Fluid mechanics

Fluid mechanics describes the behavior of fluids. It can be classified in three major subbranches. The science describing fluids at rest called fluid statics, in motion (fluid dynamics) and the interaction between solid and fluid [25].

2.1.1 Fluid flows

Laminar flow is typically apparent for high viscosity fluids with low velocity. A laminar flow is smooth and ordered. If the flow is disordered it is classified as turbulent. Turbulent flow usually appear at high veloci- ties and one of the characteristics is velocity fluctuations throughout the domain [25].

A flow is defined as steady if the flow is time independent and unsteady otherwise. As steady implies no change at a point in time, uniform implies no change in location within the domain. Unsteady flow which fluctuate in the same manner about a steady state are defined as periodic [25].

(31)

2.1. FLUID MECHANICS CHAPTER 2. THEORY

Viscosity

When a fluid in contact with another fluid is forced to move, it exerts a shear force through the contact surface in the flow direction on the stationary fluid. The same is apparent if the moving fluid is replace for- instance with a solid plate in motion. This is due to the exitance of a property quantified by what is called viscosity. It is the frictional force developed as the two fluids are forced to move relative to each other.

Viscosity can be seen as a measure of the resistance of deformation in the fluid [25].

If a fluid rate of deformation is proportional to the shear stress on the contact surface, the fluid is classified as Newtonian fluid other wise the fluid is non-Newtonian. The relationship between rate of deformation and the shear stress for a Newtonian fluid in one-dimension is described by equation 2.2 [25]. The relationship is derived and expressed by an experiment with laminar flow of a fluid between two large plates held at a certain distance from each other. One plate is forced to move at a certain velocity while the other is kept stationary. The velocity field between is linear and can be expressed by equation 2.1 [25].

u(y) = y

lV (2.1)

Where u is the the velocity at a certain point, y is the vertical distance from the stand-still plate, l is the distance between the plates and V stands for the velocity of the moving plate.

σ_s= µdu

dy (2.2)

Where σ_sstands for the shear stress acting on the fluid layer and µ stands for the dynamic viscosity. The relation 2.2 can be used to calculate the viscosity of a fluid. Therefore the above described experiment can give a measure of the viscosity [25].

The ratio between the dynamic viscosity and the density of a fluid is a recurring term in fluid mechanics. The term is therefore given a name, kinematic viscosity ν, se equation 2.3.

ν = µ

ρ (2.3)

(32)

CHAPTER 2. THEORY 2.1. FLUID MECHANICS Viscous effects are involved in all fluid flows. The effects are though more significant in some regions of the flow domain than others. In these regions the flow is classified as viscous. Regions where the internal frictional forces are negligibly small compared to for-instance pressure forces, the flow region is called inviscid flow region [25]. This enables approximations simplifying the complexity of the flow, while still keeping the accuracy level of the analysis high.

Compressibility

Compressibility defines the degree in which a fluid changes its density throughout a flow domain. If a volume portion of a fluid varies depending on the pressure at different points in the flow field, the fluid is defined as compressible and incompressible if otherwise. All fluids change in density to some degrees when pressurized, so incompressibility is actually an approximation [25]. A fluid is said to be incompressible if the density remain reasonably constant. Liquids keep there densities nearly constant and therefore are usually called incompressible substances [25]. Gases compressibility is partially model after the mach number, the speed of the gas relative to the speed of sound, se equation 2.4. If the change in density of a gas is under 5%, the flow is approximated as incompressible, usually when M a < 0.3 [25].

M a = u

c (2.4)

u is the velocity of the fluid and c stands for the speed of sound.

2.1.2 Boundary layer equations

Professor Claude-Louis Navier and George Gabriel Stokes formulated the equations of fluid motion, further discussed in section 2.1.3, during the beginning of the 1800s. Even though the Navier-Stokes equation was well known in the mid seventeens century, it is a complex differential equation that could only be solved for simple flow geometries. Scientists were able to solve approximations of the equation (the Euler equation, equation 2.5, and the Bernoulli equation, equation 2.6) [25], but these simplification by them self were not enough to resolve the complete flow, therefore these solutions could not be validated by experiments.

The Euler equation approximation is derived from the Navier-Stokes

(33)

2.1. FLUID MECHANICS CHAPTER 2. THEORY

equation by the following assumption [25].

Euler equation approximation:

Viscous forces << pressure and/or inertial forces

The bernoulli equation eq. 2.6 is derived for stationary and incompressible flows[25].

∂~u

∂t + ~u · ∇~u = −1

ρ∇P + ~F (2.5)

~u · ∇(P ρ +U²

2 + ~F ) = 0 (2.6)

u is the velocity of the fluid, U is the velocity average, P is the pressure and F stands for the external force acting on the fluid.

The assumption means that the viscous term in the Navier-Stokes equation can be neglected. This assumption is not valid close to solid walls where it would mean free slip condition, which does not yield a physical solution. The viscous forces near a solid wall cause the no slip condition, they are therefore significant and can not be neglected.

The breakthrough came in 1904 by the German physicist Ludwig Prandtl.

He suggested that the flow was to be divided into two regions, an inner and an outer flow region. The outer flow regions velocity field would be solved by the Euler and continuity equation. The Bernoulli equation would be used for solving the pressure field for the region, where viscous forces can be neglected. The inner region, called boundary layer, is a very thin layer where viscous forces are not negligible and the flow is rotational. Professor Prandtls approximation is what is known as the boundary layer approximation [25].

The boundary layer is assumed to be very thin and the pressure is shown to be constant across a boundary layer [25]. The boundary layer for incompressible flows is described by equation 2.7.

~u · ∇~u = −1

ρ∇P + ν∇²~u (2.7)

(34)

CHAPTER 2. THEORY 2.2. SOLID MECHANICS

2.1.3 Navier-Stoke and continuity equation

The motion of fluid is mainly described by two differential equations. The conservation of mass also known as the continuity equation, equation 2.8, and the Navier-stokes equation, equation 2.9. These equations combined are able to fully solve the details of the flow for every point in the domain [25]. These differential equations are extremely complex and difficult to solve. The equations are usually combined with other equations and necessary approximation in order to be solved analytically. However with the help of computational fluid dynamics a solution of the equations of motion can be obtain, the method is described further in section 2.3.1.

Dρ

Dt + ρ(∇ · ~u) = 0 (2.8)

D~u Dt = −1

ρ∇P + ν∇²~u + ~F (2.9)

2.2 Solid Mechanics

Material strength testing and theories thrived during the 20th century.

The engineer Alan Arnold Griffith made a profound contribution the the study of materials when he published his first work on the subject in 1920, "The phenomenon of rupture and flow in solids". He made a series of tensile tests with glass bars. The test showed that the tensile strength σ_t, the stress needed to pull a material to the breaking point, seemed to decrees with increasing specimens thickness. The discovery lead him to the following conclusion, there are defects within the material leading to it breaking at much lower stresses than theoretically anticipated. The number of defects is increased with increased bar thickness and that is the cause for lower tensile strength. The conclusion of the study was only valid for extremely brittle materials, since it only took the elastic work into consideration [14].

The theory was further modified in the 1940s by the scientist George Rankine Irwin and the team under him during his time at the US Naval Research Laboratory. He included the plastic work done during tensile testing and developed the concept of stress intensity factors [14].

(35)

2.2. SOLID MECHANICS CHAPTER 2. THEORY

2.2.1 Mechanical stress and strain

The cohesive forces between atoms within a material defines the ideal strength of the material. The word ideal is used here because such strength refers to a material with no discontinuities, defects, and as is known such a material does not exist in reality. The theory is through still applicable in certain loading cases discussed further in the following sections.

The strength is quantified by the material property Youngs modulus E. E is given by Hooks law and it is obtained by calculating the ratio between the normal stress applied and the strain, elongation or elastic deformation, it forces on the material, se equation 2.12 [9]. Figure 2.1 illustrates the standard visualization of the stress-strain relation, the stress is calculated by equation 2.10 and the strain by equation 2.11 [5].

(a)

(b)

Figure 2.1: Figure showing a bar with a diameter D, length L, cross section area after deformation is A, the material has the Youngs modulus E. The force F applied on the bar has resulted in an elongation ². b) Shows the nominal stress σ at the cross section.

σ = F

A (2.10)

(36)

² = ∆L

L (2.11)

E = σ

² (2.12)

Where σ is the stress and ² is the elongation. As seen in equation 2.12, the relation between stress and strain is assumed to be linear. This is approximately true for a material, except for polymers, in the elastic domain before plastic deformation occurs.

In the case of a nonuniform cross section stress, due to discontinuities or geometrical changes in the structure, the following derivation results in an expression for the stress distribution.

Equation 2.10 gives the following for a one dimensional stress distribution:

σ_c= ∆F

∆A (2.13)

σ_c stands for the cross section stress. To look at the stress in one point in the cross section, δA → 0, equation 2.13 becomes:

σ_c= lim

∆A→0

∆F

∆A = dF

dA (2.14)

So the stress can vary at different points on the cross section area. The net force can then be obtained as follows [9].

dF = σ_c(A)dA ⇒ (2.15)

⇒ F_net= Z

σ_c(A)dA (2.16)

Deformation

If a test specimen is loaded with a load that does not exceed the materials yield strength σ_y, the deformation forced on the specimen will disappear when the load is removed, the test subject returns to its original shape.

This type of deformation is classified as elastic. Plastic deformation occurs if σ_y is exceeded. In this case, when the load is released the material will recover the elastic deformation but keep the plastically deformed length [23, 14]. Figure 2.2 shows a standard stress-strain curve plotted for a material after tensile testing. The figure also gives an idea of the information such a curve reveals about the material properties of the specimen.

(37)

2.2. SOLID MECHANICS CHAPTER 2. THEORY

Figure 2.2: The figure shows a standard stress strain curve obtained by tensile testing. The figure to the right shows the elastic contra the plastic region in the curve. On the figure to the left, σ_t is the tensile strength, σ_y is the yield strength, ²_t is the tensile strain and ²_f is the fracture strain.

Stress concentration factor

The definition of stress concentration is, a local elevation of stress levels.

Consider the setup in figure 2.1. If a deficiency is introduced within the bar, in the form of a vacant sphere of radius R, take out a small sphere of material as sown in figure 2.3. The same loading condition is applied.

Figure 2.3: The figure shows the same cross section as in figure 2.1b but here a defect in the form of a vacant sphere.

The stress sheared by the present material in the previous case is now sheared by the material around the defect. This leads to elevated stress levels around the sphere, an area of stress concentration appears [9]. The cross section stress distribution is no longer homogeneous. The stress distribution seen in figure 2.1b is called the nominal stress and is noted as σ_n. The stress distribution around the sphere is much more complex and

(38)

three dimensional than illustrated in figure 2.3 [9].

Lets look at the case for an infinite and thin plate, loaded with the nominal tensile stress σ_n. The plate has a hole with the diameter d, see figure 2.4.

Figure 2.4: The figure shows the stress concentration condition for a small hole with diameter d in an infinitely large and thin plate.

The stress distribution near the hole is the same as in the previous case, it is complex, but the stress at the circle curve is shown to be described by equation 2.17 [9].

σ_θ = (1 − 2cos(2θ))σ_n (2.17) Equation 2.17 gives σ_θ = 3σ_n for θ = 90^◦ and θ = 270^◦, which is noted as σ_max. The equation also gives σ_θ = −1σ_n for θ = 0^◦ and θ = 180^◦, which is the minimum stress at the circle. The negative value means that the stresses in these areas are contractive. The stress concentration factor is derived from the maximum stress and is formulated by equation 2.18 [12], it should be noted that the equation is only valid for a circular hole. In this case the concentration factor is K_t= 3.

K_t= σ_max

σ_n (2.18)

The expressions of K_t differ for different types of discontinuities. Dia- grams, tables and exact expressions of the solution of K_tare today given in books available for construction engineers and scientists.

2.2.2 Fatigue

Mechanical structures have been known to reach fracture failure at load levels below the tensile strength of the material, fractures have also oc- curred at stress levels even below the yield strength. This phenomenon is

(39)

2.2. SOLID MECHANICS CHAPTER 2. THEORY a progressive structural damage known as fatigue. A structure subjected to load cycles can fail if the load amplitude is high enough to achieve fracture, it is called fatigue failure or fatigue fracture [12].

Fatigue loading

Fatigue life assessment are often characterized by S − N curves, also known as Wöhler curves. The curve plots the stress (S) against the number of cycles (N ) needed to achieve fatigue failure [12]. Fatigue fail- ure often starts with crack nucleation (initiation) and then as the load cycles continue, the crack grows continuously eventually leading to fracture [12].

There are two main types of fatigue loading, Low cycle fatigue (LCF) and high cycle fatigue (HCF). Low cycle fatigue, is defined at approxi- mately N = 10³, and occur where the stress levels are high enough to cause plastic deformation. LCF is often characterized by Coffin-Manson relation given by equation 2.19.

∆²_p

2 = ²⁰_f(2N )^c (2.19)

Where ^∆²₂^p is the plastic strain amplitude, ²⁰_f is the empirical constant fatigue ductility coefficient and c is the empirical constant fatigue ductility exponent. HCF is usually defined for cycles ranging between N = 10⁴− 10⁸. HCF occur for load cycles with stress levels below the yield strength where the deformation is primarily elastic.

Fatigue crack and crack growth

G. R. Irwin extended theory of A. A. Griffiths work lead Irwin and his team to the formulation of the stress intensity factor, K. Griffith formulated the concept, that an initiated crack in a brittle environment will grow during loading if the total energy of the system is lower by the crack propagation. Írwens work based on Griffiths findings developed a theory for crack propagation in ductile material and formulated equation 2.20 [12].

K = Cσ_n√

πa (2.20)

Here a is the crack length and C is a dimensionless constant derived for the crack shape and the geometry around the crack. Exact solutions of

(40)

CHAPTER 2. THEORY 2.3. NUMERICS the crack tip stress field for an infinite plate and relatively simple load cases are available in engineering books.

The solutions for stress intensity factors available are often calculated for simpler load cases. Mechanical structures inhabit more complex stress fields originating from various sources, like residual stress fields and stress field at notches (geometrical discontinuities). Due to these discontinuities stress concentration is evident, which lead to gradient stress fields.

The stress gradient, χ, given by equation 2.21, for propagating cracks are in most cases not constant [12]. So solving the stress intensity factor is a complex and non trivial task.

χ = 1 σ_max

¯¯

¯¯δσ_y δx

¯¯

x=x0

(2.21) The stress gradient in equation 2.21 is given for a spherical cavity in an infinite body in tension in the y axes direction. Here x₀ is where the maximum stress is located. The crack is shown to grow in the direction of the largest stress gradient and perpendicular to the maximum mean stress [12].

To accurately predict the crack propagation of a crack in a mechani- cal structure, the geometry factor, F , and the stress gradient are key components. The stress field at the propagation front is influenced by the stress gradient and visa versa. This means that in order to obtain a solution as accurate as possible, χ is preferably recalculated after each crack growing load cycle [12].

2.3 Numerics

Numerical analysis is a representation of continues mathematics by algorithms. Numerical mathematics has been around for a long time.

The earliest numeric mathematical find dates back to somewhere be- tween 1800 − 1600BC in ancient Babylonia. The find was a clay tablet with an engraved numerical approximation of the square root of 2, √

2.

Great mathematicians and physicists like Isaac Newton, Joseph Louis Lagrange, Carl Friedrich Gauss and Leonhard Euler have formulated numerical algorithms used for a wide range of computational applications [19].

(41)

2.3. NUMERICS CHAPTER 2. THEORY

Computational numerics can be seen as the discretization of the physical setup of the problem at hand. The methodology and process involved in applying such calculations on a physical problem is summarized by the figure 2.5 [10].

The process called idealization describes going form the actual physical

Figure 2.5: Figure showing simplistic view of the methodology of the discretization process. The image is a modified one that originally was found in Felippa 2002 [10].

system to a simplified mathematical model. The model is designed ex- tracting certain aspects of behavior of the physical system, depending on the application. Therefore the idealization process is the most sensitive and important step, since there are many factors that need to be taken into account. Mathematical modeling is basically a filtering step, where physical details which are not of interest for the analysis are filtered out [10]. These models are usually a system of coupled partial differential equation in space and time. They are often complex not easy to solve. The variables of such a system are called degrees freedom (DOF).

A mathematical idealized model has an infinite number of DOFs, it is what is called a continues model [10].

It is necessary to reduce the number of DOFs to a finite number, in order to gain a numerically handleable and solvable system, this step is called discretization (transferring continues equations into discreet counterparts) [10]. A discreet model is obtained by spatial and time discretization. The discreet model is only an approximate model of the mathematical one, so validation of the process is of upmost importance.

(42)

CHAPTER 2. THEORY 2.3. NUMERICS

Every step is an error source introducing errors that need to be un- derstood and accounted for in the final solution. The only way to really do that is by studying the real physical systems behavior, what is meant by that statement is validation against experimental results.

2.3.1 Computational Fluid Dynamics (CFD)

Computational fluid dynamics has made a profound impact on the design process for engineer. This section will briefly go through the solution methodology and theory of CFD calculations. The CFD method is used for numerical analysis of the flow of fluids. The equations of motion are as mentioned earlier a set of partial differential equation (PDEs) that describe the flow, equation 2.8 and 2.9.

Obtaining a solution for a turbulent flow is not an easy task. The finer features the flow are unsteady and three-dimensionally random. There are randomly swirling vortical structures within a turbulent flow called turbulent Eddies. They are of various sizes and time scales and add to the difficulties of the calculation [25]. Even though the complexity of a turbulent flow is high, a meaningful, physical and increasingly accurate solution can still be obtained. Such a solution requirers the total devo- tion to understanding the setup. Knowledge of the mathematical and physical ground of the specific problem is essential. The computational time can be largely reduced by appropriate approximations and assump- tions.

There are many different methods, of different difficulties and accuracy, used to solve the PDEs. The direct numerical simulation (DNS) method is one of those. The method solves the unsteady motions of the flow directly for all scale sizes. The method gives a highly accurate solution.

However the size and time scale difference between large and small Ed- dies make the DNS method computationally heavy. It requires fine grids and a great deal of CPU time [11].

A simpler technique is to only resolve the large Eddies and model the small Eddies, the method is called Large Eddy simulation (LES). The model of the small Eddies basically assumes that they are isotropic (independent of direction) and behave in a statistically similar manner. This

(43)

2.3. NUMERICS CHAPTER 2. THEORY method is thus not as heavy as DNS and require less computational time.

The method is still time consuming [25].

Turbulence models

Modeling all the turbulence features of the flow with what is called turbulence models is the most commercially used. Here mathematical models of all the turbulent unsteady features, such as the mixing and diffusion caused by the turbulent Eddies, are made. The calculation difficulty is considerably reduces by this solving method, while the accuracy is held high depending on the chosen model.

For CFD calculation the Navier-stokes equation is time, velocity and pressure averaged via equation 2.22 and 2.23. The resulting equation is the Reynolds-averaged Navier-stokes equation (RANS), given in tensor form by equation 2.24 [11].

P = ¯P + p⁰ (2.22) u_i= ¯U_i+ u⁰_i (2.23) δU_i

δt + U_jδU_i δx_j = −1

ρ δP

δx_i + νδ²U_i δx²_i +δτ_ij

δx_i (2.24)

Where U_i and P are statistical averages, u⁰_i and p⁰ are the fluctuations for the velocity and the pressure and τ_ij is the Reynolds stress tensor [11].

Equation 2.24 is then modeled to be solved. There are many turbulence models in use today including one-equation, two-equation and Reynolds stress models (RSM). The first mentioned is mostly validated for aerodynamical flows and thus primarily used within that area. The Reynolds stress models are some of the more complex turbulence models in common use. They are based on modeling dynamic equations for the Reynolds stress tensor. The improvement in solution accuracy, said to be gained with RSM models, is still to be proven.

Two-equation models are the predominate turbulence models in use today. This is because they offer a balance between numerical effort and computational accuracy. The mainly used two-equation models are the k −² and shear stress transport (SST) models. whereas the -˛² is the more widely used. The SST model combines the use of two models, the k − ² and the k − ω model. Both models are modeled relating the Reynolds

(44)

CHAPTER 2. THEORY 2.3. NUMERICS stresses to the mean velocity gradients and the turbulent viscosity [11].

The k − ² models is shown to handel the free stream flow outside the boundary layer with greater accuracy, but is less efficient in the boundary layer region, where the Reynolds number is low. The model needs finer grids near the wall to account for the flow formulation of the bound- ary layer. The k − ω on the other hand gives a better resolution of the boundary layer but is said to be less stable within the free-stream flow.

The model is shown to be sensitive to the inlet free-stream turbulence properties [11]. The SST model combines, as mentioned, the two previ- ously discussed. It uses the k − ω formulation in areas of low Reynolds number flows (near walls) and the k − ² model further away from the walls, in the free-stream. Although the method offers a higher accuracy solution, it requires an enormous amount of CPU time and therefore usually chosen after an initial solution by k − ², if required. The k − ² is by far the most used model, it is shown to be very accurate for many diffident applications.

2.3.2 Finite Element Method (FEM)

The finite element method is a discretization method used in many applications. The method accuracy is though most frequently validated for structural mechanic calculations. FEM analysis is carried out by decom- posing the whole physical domain into a finite number of subdomains, finite elements, and calculating an approximate solution for each subdo- main.

The finite elements come in different shapes and they possess various properties. They are of different special dimensions one, two and three dimensions. The elements are also of different complexity. Elements of the same classification can be of different mathematical orders, linear or of higher order. The choice of element directly impacts the accuracy and the convergence of the solution.

For the finite element method, the elements do not overlap in space as for some other discretization methods like the finite difference method.

This property is called disjoint support. Each element consists of one (in the case of a one dimensional elements) or several distinguished points called nodes. Each node serve a purpose. The classifications are geo- metric nodes for nodes defining the element geometry and connection

(45)

2.3. NUMERICS CHAPTER 2. THEORY nodes for nodes home for the defined degrees of freedom [10]. For many studies the two definitions can define one and the same node. Figure 2.6 illustrates the shapes, degrees and node positions for some elements.

Figure 2.6: Figure illustrating the shapes, degrees and node positions for some elements.

An approximation of the mathematical models describing the behavior of the physical problem is represented by suitable functions in each element and solved. The mathematical models of the physical setup are usually a set of partial differential equations or integral equations. The element equations are then assembled and a discrete solution for the system of equations is then obtained, meaning the values for every node is solved. That way approximate solution for the whole domain is obtained [3]. The accuracy of the FEM solution refers to the mathematical models ability to predict the behavior of interest of the real system. The term convergence in FEM analysis sense refers to the accuracy of the discretization of the mathematical model [22].

The general equation of the finite element method is given by equation

(46)

CHAPTER 2. THEORY 2.3. NUMERICS

2.25 [22], where an approximate solution for u is of interest.

u ≈ Xn j=1

u_jψ_j+ Xm j=1

c_jφ_j (2.25)

n is the number of nodes, u_j is the values of u at each element node, ψ_j stand for the interpolation functions, c_j are the nodeless coefficients and φ_j are the approximation derived by using concepts of interpolation theory functions [22].

The solution given in equation 2.25 leads to the formation of the known matrix expression 2.26 [22].

[K]{u} = {f } + {Q} (2.26)

Where [K] is the coefficient matrix and {f } is the source vector, they are called the stiffness matrix and the force vector respectively in solid and structural mechanics. {u} and {Q} are the primary and the secondary element nodal degrees of freedom [22]. The expression describes the FEM problem to be solved.

(47)

(48)

CHAPTER 3. MODELING

Chapter

3 Modeling

When creating a computer model, approximations are made. As mentioned in the theory chapter each simplification introduces an error that need to be marginalized as much as possible.

The geometry used for the calculations was, as mentioned in previous chapters, the new design runner. A design geometry (Pre production model, the model prior to manufacture) was used for creating a CFD and FEA model, figure 3.1 shows the design model. It should be noted that the design model does not include the fillet at the blade-band and crown T-joint section.

(49)

CHAPTER 3. MODELING

(a) (b)

(c) (d)

Figure 3.1: Figure of the designed runner model. (a) Showing a trig view, (b) a side view,(c) a top view and (d) a bottom view.

(50)

CHAPTER 3. MODELING 3.1. CFD MODEL

Figure 3.2: Figure showing the CFD model used in the analysis. It represents a third of the total geometry shown in figure 3.1.

3.1 CFD model

3.1.1 Geometry

The fluid domain is extracted from the design model. Geometrical symmetry allowed for a reduction of the model by two thirds, and thereby gaining a reduction in computational time. The CFD model is shown in figure 3.2. The model contains five blades, five splitter blades, eight guide vanes and an extension of the outlet. The geometry was extracted using the contour of the band and crown. The extension to the outlet was added to avoid any effect it would have had on the flow near the areas of interest. It should be noted that the extension is not a part of the draft tube.

3.1.2 Preprocessing

The mesh used in the CFD calculation is a tetrahedra mesh with a size of approximately 710, 000 nodes, se figure 3.3. The simulation setup was a transient rotor stator simulation. The turbulence model chosen for the simulation was the two equation k − ² model. The advection scheme was high resolution with a second order backward Euler transient scheme. The inlet boundary conditions were chosen to be of the type total pressure, set to the pressure equal to the head hight. The flow direction at the inlet was specified to flow in at an appropriate angle

(51)

3.2. FEM MODEL CHAPTER 3. MODELING

(a) (b)

Figure 3.3: Figure of the CFD Mesh. (a) Showing a top view, (b) a zoom in, for a closer look at the mesh.

with respect to the angle of the stay vanes. The outlet conditions were set to average static pressure 0. A time step of 0.5^◦ was decidedly used and a maximum of 10 iterations per loop was set. The torque of the blade and splitter was monitored during the simulation. The monitoring of the blades was done to enable a manual stop of the simulation once the torque reached a steady periodicity.

3.2 FEM model

As mentioned earlier some approximation was made in order to keep the computational time limited. The original sketch had complex band and crown contours which needed to be smoothed or in some cases remade.

When changing the geometry, there is a constant struggle between simplifying to gain a manageable simulation time and keeping the changes to areas where they do not affect the results. Figure 3.4 shows an example of a simplification made on the band geometry.

(52)

CHAPTER 3. MODELING 3.2. FEM MODEL

(a)

(b) (c)

Figure 3.4: Figure showing one of the simplifications on the design model in order to save CPU time. (b) The zoom in design model contour on the band, (c) a zoom in, for a look at simplification done at the area.

(53)

Figure 3.5: Figure showing the finite element geometry of model one used in the analysis. Containing a splitter blade, blade and part of the crown and band.

3.2.1 Geometry

Two different FEA models were made. Model one is the simulation model used in the analysis. Model two was used to investigate wether the boundary conditions being close to the areas of interest had any effect on the results, especially near the trailing edge of the blade where the boundary conditions are closest to the blade.

Model one

The FEA model created had a portion of the band and crown, a blade and a splitter blade, see figure 3.5. Since the runner consists of 15 blades and splitters, model one was cut out of the total runner with an angle of 360^◦/15 = 24^◦. The cut to the crown and band, as seen in figure 3.5, follows the blades shape.

Model two

The model was made so that the cyclic symmetry boundary condition, explained in section 3.2.3, was put further away from the areas of interest.

(54)

CHAPTER 3. MODELING 3.2. FEM MODEL

(a) (b)

Figure 3.6: Figure showing the finite element geometry of model two used in the analysis. a) showing a back view and b) a front.

The geometry is shown in Figure 3.6. As can be seen the model has the same shape as FEA model one, but here a splitter blade is added on one side and a blade on the other.

3.2.2 Mesh

Model one and two were meshed with second order tetrahedra mesh elements with quadratic displacement behavior [1], se figure 3.7. They are ten node elements with three DOF at each node, translation in x, y and z.

Second order surface elements with four to eight nodes were created, figure 3.8. The elements were used to apply the pressure load. They were created on the surfaces of the blade and splitter blade.

The mesh can be seen in figure 3.9. The mesh was made uniform with one average element side length obtained from a mesh study. The mesh study was preformed to minimize the numerical error in the simulation model.

(55)

Figure 3.7: Figure showing the type of solid mesh element used in the analysis. It is a ten node tetrahedra element with quadratic displacement behavior. The image is one gathered from the Ansys 11 manual pages [1].

Figure 3.8: Figure showing the type of surface mesh element used in the analysis. It is a Second order elements with four (on the right) to eight (on the left) nodes. The image is one gathered from the Ansys 11 manual pages [1].

Figure 3.9: Figure showing the finite element mesh of model two used in the analysis. The figure is showing a pressure side view.

Dynamic stress assessment in high head Francis runners

MASTER’S THESIS

Yasir Dawood Siwani

Dynamic stress assessment in high head Francis runners

Dynamic stress assessment in high head Francis runners

Yasir Dawood Siwani

Preface

Abstract

Contents

List of Figures

List of Tables

1

Introduction

1.1 Background

1.2 Hydraulic machines

1.3 Site measurements

2

Theory

2.1 Fluid mechanics

2.2 Solid Mechanics

2.3 Numerics

3

Modeling

3.1 CFD model

3.2 FEM model