Aerodynamic Design of Horizontal Axis Wind Turbine Blades

Designing horizontal axis wind turbine (HAWT) blades to achieve satisfactory levels of performance starts with the knowledge of aerodynamic forces acting on the blades. In this paper, a design method based on blade element momentum (BEM) theory is explained for HAWT blades. The method is used to optimize the chord and twist distributions of the blades. Applying this method a 1000W HAWT rotor is designed. A user-interface computer program is written on VISUALBASIC to estimate the aerodynamic performance of the existing HAWT blades and used for the performance analysis of the designed 1000W HAWT rotor. The program gives blade geometry parameters (chord and twist distributions), coefficients of performance and aerodynamic forces (trust and torque) for the following inputs; power required from the turbine, number of blades, design wind velocity and blade profile type (airfoil). The program shows the results with figures. It also gives the three dimensional views of the designed blade elements for visualization after exported to AutoCAD.


MOTIVATION AND BACKGROUND
Energy prices, supply uncertainties, and environmental concerns are driving the developed nations to rethink their energy mix and develop diverse sources of clean, renewable energy.Wind power, as an alternative to fossil fuels, is growing at the rate of 19.2% annually, with a worldwide installed capacity increased from 432,680MW in 2015 to 486,790MW at the end of 2016 [1][2][3].All wind turbines installed by the end of 2016 worldwide can generate around 5% of the world's electricity demands.Wind power is widely used in Europe, Asia, and the United States.China, United States, Germany, India and Spain respectively are top five countries in terms of installed cumulative capacities at the end of 2016.
Wind power is also capable of becoming a major contributor to the developing nations of Africa.
According to 2016 Global Wind Energy Council (GWEC) statistics, all wind turbines installed in Africa had a capacity of 3,906MW [1][2][3].Currently, South Africa is leading the continent with installed capacity of 1,471MW followed by Egypt with 810MW, Morocco with 787MW and Ethiopia with 324MW.Other promising countries in Africa include Tunisia, Kenya, Tanzania, where wind project development is firmly underway [1][2][3].
Wind energy application in Ethiopia has been limited to water pumping in the past.There is now, however, definite plan to exploit wind for power production.According to the Ethiopian Electric Power Corporation (EEPCo), Ethiopia has a capacity of generating more than 10,000MW from wind.Considering this substantial wind resource in the country, the government has committed itself to generate power from wind plants by constructing eight wind farms with total capacities of 1,116MW [3].So far, two wind farms such as Ashegoda, Adama I and Adama II wind farms started producing 324MW.Some of the erected wind turbines in Adama I wind farm are shown in Fig. 1 [3].Wind turbine power production depends on the interaction between the rotor and the wind.The blades of wind turbine receive kinetic energy of the wind, which is then transformed to mechanical energy, electrical forms, depending on our end uses.The objective of this study is to develop a user-interface computer program for Horizontal-Axis Wind turbine (HAWT) blade design and power performance prediction using the Blade Element Momentum (BEM) theory.

Figure 2. Main elements of HAWTs
The primary objective in wind turbine design is to maximize the aerodynamic efficiency, or power extracted from the wind.However, the mechanical strength criteria, and economical aspects are also equally important and therefore must be considered during the design of wind turbines.The design of HAWT blades, to achieve satisfactory level of performance; starts with knowledge of the aerodynamic forces acting on the blades.In this article, HAWT blade design from the aspect of aerodynamics is investigated.The structural design of HAWT blades is also as important as their aerodynamic design.The dynamic structural loads which a rotor will experience play the major role in determining the lifetime of the rotor.Obviously, aerodynamic loads are a major source of dynamic structural behavior, and ought to be accurately determined.In addition, the blade geometry parameters are required for dynamic load analysis of wind turbine rotors [4][5][6][7].

AERODYNAMICS OF HAWTS
Many factors play a role in the design of a wind turbine rotor, including aerodynamics, generator characteristics, blade strength and rigidity, noise levels, etc.But since wind energy conversion system is largely dependent on maximizing its energy extraction, rotor aerodynamics play important role in the minimization of the cost of energy.The rotor consists of the hub and blades of the wind turbine.These are often considered to be its most important components from both performance and overall cost standpoint.

The Actuator Disc Theory and The Betz Limit
A simple model, generally attributed to Betz (1926) can be used to determine the power from an ideal turbine rotor, the thrust of the wind on the ideal rotor and the effect of the rotor operation on the local wind field.The simplest aerodynamic model of a wind turbine is known as 'actuator disk model' in which the rotor is assumed like a homogenous disk that removes energy from the wind.
The theory of the ideal actuator disk is based on the following assumptions [8]: Homogenous, incompressible, steady state fluid flow, no frictional drag; The pressure increment or thrust per unit area is constant over the disk; The rotational component of the velocity in the slipstream is zero; There is continuity of velocity through the disk; An infinite number of blades.The analysis of the actuator disk theory assumes a control volume as shown in Fig. 3.
In the control volume below, the only flow is across the ends of the stream tube.The turbine is represented by a uniform actuator disk which creates a discontinuity of pressure in the stream tube of air flowing through it.Note also that this analysis is not limited to any particular type of wind turbine.From the assumption that the continuity of velocity through the disk exists, the velocities at section 2 and 3 are equal to the velocity at the rotor: Applying the conservation low of linear momentum to the control volume enclosing the whole system, the net force can be found on the contents of the control volume.That force is equal and opposite to the thrust, T which is the force of the wind on the wind turbine.Hence from the conservation of linear momentum for a one-dimensional, incompressible, time-invariant flow the thrust is equal and opposite to the change in momentum of air stream: No work is done on either side of the turbine rotor.Thus the Bernoulli function can be used in the two control volumes on either side of the actuator disk.Between the free-stream and upwind side of the rotor (from section 1 to 2 in Fig. 3) and between the downwind side of the rotor and far wake (from section 3 to 4 in Fig. 2) respectively: The thrust can also be expressed as the net sum forces on each side of the actuator disk: By using Equations ( 4) and ( 5), the pressure decrease, p′ can be found as: And by substituting Equation (8) in to Equation ( 6): By equating the thrust values from Equation (3) in which substituting Equation ( 2) in place of m ɺ and Equation (9), the velocity at the rotor plane can be found as: If an axial induction factor, a is defined as the fractional decrease in the wind velocity between the free stream and the rotor plane, then Using Equation (10) and ( 12) one can get: The velocity and pressure distribution are illustrated in Fig. 4. Because of continuity, the diameter of flow field must increase as its velocity decreases and note that there occurs sudden pressure drop at rotor plane which contributes the torque rotating turbine blades.And substituting for R u and w u from Equations ( 12) and (13) to Equations (15): The power performance parameters of a wind turbine can be expressed in dimensionless form, in which the power coefficient, C P is given in the following equation: Using Equation ( 16) and ( 17), the power coefficient C P becomes: The maximum C P is determined by taking the derivative of Equation ( 18) with respect to a and setting it equal to zero yields: when, a = 1/3 This result indicates that if an ideal rotor were designed and operated so that the wind speed at the rotor were 2/3 of the free stream wind speed, then it would be operating at the point of maximum power production.This is known as the Betz limit.
From Equations ( 9) and ( 13) the axial thrust on the disk can be written in the following form: Similar to the power coefficient, thrust coefficient can be defined by the ratio of trust force to dynamic force as shown in the following equation: Using Equation ( 20) and ( 21), the trust coefficient C T becomes: Note that C T has a maximum of 1.0 when a = 0.5 and the downstream velocity (wake velocity) is zero.At maximum power output (a =1/3), C T has a value of 8/9.A graph of the power and trust coefficients for an ideal Betz turbine and the non-dimensionalized downstream wind speed are illustrated in the Fig. 5 [9].
As mentioned above, this idealized model is not valid for axial induction factors greater than 0.5.In practice, as the axial induction factor approaches and exceeds 0.5, complicated flow patterns that are not represented in this simple model result in thrust coefficients that can go as high as 2.0 [10].
In conclusion, the actuator disk theory provides a rational basis for illustrating that the flow velocity at the rotor is different from the free-stream velocity.The Betz limit, C P , max =16/27, shows the maximum theoretically possible rotor power coefficient that can be attained from a wind turbine.In practice, three effects lead to a decrease in the maximum achievable power coefficient: rotation of the wake behind the rotor; finite number of blades and associated tip losses; non-zero aerodynamic drag.Note that the overall turbine efficiency is a function of both the rotor power coefficient and the mechanical (including electrical) efficiency of the wind turbine.Therefore:

The General Momentum Theory
In the previous analysis using linear momentum theory, it was assumed that no rotation was imparted to the flow.The previous analysis can be extended to the case where the rotating rotor generates angular momentum, which can be related to rotor torque.In the case of a rotating wind turbine rotor, the flow behind the rotor rotates in the opposite direction to the rotor, in reaction to the torque exerted by the flow on the rotor.An annular stream tube model of this flow, illustrating the rotation of the wake, is shown in Fig. 6 [9].The generation of rotational kinetic energy in the wake results in less energy extraction by the rotor than would be expected without wake rotation.In general, the extra kinetic energy in the wind turbine wake will be higher if the generated torque is higher.Thus, as will be shown here, slow-running wind turbines (with a low rotational speed and a high torque) experience more wake rotation losses than high-speed wind machines with low torque.Fig. 7 gives a schematic of the parameters involved in this analysis.Subscripts denote values at the crosssections identified by numbers.If it is assumed that the angular velocity imparted to the flow stream, ω, is small compared to the angular velocity, Ω, of the wind turbine rotor, then it can also be assumed that the pressure in the far wake is equal to the pressure in the free stream [10].The analysis that follows is based on the use of an annular stream tube with a radius r and a thickness dr, resulting in a cross-sectional area equal to 2πrdr (see Fig. 7).The pressure, wake rotation, and induction factors are all assumed to be functions of radius.If one uses a control volume that moves with the angular velocity of the blades, the energy equation can be applied in the sections before and after the blades to derive an expression for the pressure difference across the blades [11][12].Note that across the flow disc, the angular velocity of the air relative to the blade increases from Ω to Ω+ ω, while the axial component of the velocity remains constant.The results are: The resulting thrust on an annular element, dT, is: An angular induction factor, a', is then defined as: Note that when wake rotation is included in the analysis, the induced velocity at the rotor consists of not only the axial component, u ∞ a, but also a component in the rotor plane, rΩa'.
The expression for the thrust becomes: Following the previous linear momentum analysis, the thrust on an annular cross-section can also be determined by the following expression that uses the axial induction factor, a, (note that u 1 , the free stream velocity, is designated by u ∞ in this analysis): Equating the two expressions for thrust gives: where λ r is the local tip speed ratio (see below).This result will be used later in the analysis.
The tip speed ratio, λ, defined as the ratio of the blade tip speed to the free stream wind speed, is given by: The tip speed ratio often occurs in the aerodynamic equations for the rotor.The local speed ratio is the ratio of the rotor speed at some intermediate radius to the wind speed: Next, one can derive an expression for the torque on the rotor by applying the conservation of angular momentum.For this situation, the torque exerted on the rotor, Q, must equal the change in angular momentum of the wake.On an incremental annular area element this gives: Since u 2 =u ∞ (1-a) and a'=ω/2Ω, this expression reduces to: The power generated at each element, dP, is given by: dQ dP Ω = (34) Substituting for dQ in this expression and using the definition of the local speed ratio, λ r , (Equation 31), the expression for the power generated at each element becomes: It can be seen that the power from any annular ring is a function of the axial and angular induction factors and the tip speed ratio.The axial and angular induction factors determine the magnitude and direction of the air flow at the rotor plane.The local speed ratio is a function of the tip speed ratio and radius.
The incremental contribution to the power coefficient, dCp, from each annular ring is given by: In order to integrate this expression, one needs to relate the variables a, a', and λ r , [13][14].Solving Equation (29) to express a' in terms of a, one gets: The aerodynamic conditions for the maximum possible power production occur when the term a'(1-a) in Equation ( 37) is at its greatest value.Substituting the value for a' from Equation (38) into a'(1-a) and setting the derivative with respect to a equal to zero yields: This equation defines the axial induction factor for maximum power as a function of the local tip speed ratio in each annular ring.Substituting into Equation (29), one finds that, for maximum power in each annular ring: If Equation ( 39) is differentiated with respect to a, one obtains a relationship between dλ r and da at those conditions that result in maximum power production: Now, substituting the Equations ( 39)-(40) into the expression for the power coefficient (Equation 37) gives: Here the lower limit of integration, a 1 , corresponds to the axial induction factor for λ r = λ h =0 and the upper limit, a 2 , corresponds to the axial induction factor at λ r =λ.Also, from Equation (39): Note that from Equation (39), a 1 =0.25 gives λ r a value of zero.Equation (43) can be solved for the values of a 2 that correspond to operation at tip speed ratios of interest.Note also from Equation (43), a 2 = 1/3 is the upper limit of the axial induction factor, a, giving an infinitely large tip speed ratio.
The definite integral in Equation 42 can be evaluated by changing variables: substituting x for (1-3a) in Equation ( 42) and evaluating the integral, the result is: The result of general momentum theory are graphically represented in Fig. 8 which also shows the Betz limit of an ideal turbine based on the linear momentum analysis performed in the previous section.The result shows that the higher the tip-speed ratio, the greater the maximum power coefficient.Note that Fig. 8, represents the numerical values for C P,max as a function of λ, with corresponding values for the axial induction factor at the tip, a 2 .

Blade Element Theory
The momentum theories, which have been developed in the previous sections, are based on a consideration of the mean axial and rotational velocity in the slipstream and determine the thrust and torque of a blade from the rate of decrease of momentum of the fluid.The theories determine an upper limit to the power coefficient of any blade, depending on the free-stream wind velocity and on the power extracted, but they restrict the understanding of the effect of rotor geometry (i.e.blade airfoil section, chord and twist).The blade-element theory is an alternative method of analyzing the behavior of blades due to their motion through air.
As shown in Fig. 9, for this analysis, it is assumed that the blade is divided into N sections (or elements) and the aerodynamic force acting on each blade element can be estimated as the force on suitable airfoil characteristics of the same cross-section adopted for the blade elements.Finally assuming that the behavior of each element is not affected by the adjacent elements of the same blade, the force on the whole blade can be derived by adding the contributions of all the elements along the blade [15].A diagram showing the developed blade element at radius r and the velocities and forces acting on this element is given in Fig. 10.The relative wind velocity u rel is the vector sum of the wind velocity at the rotor u ∞ (1-a) (the vector sum of the free-stream wind velocity, u ∞ and the induced axial velocity -au ∞ ) and the wind velocity due to rotation of the blade.And this rotational component is the vector sum of the blade section velocity, Ωr and the induced angular velocity a'Ωr.Hence the relative wind velocity will be as shown on the velocity diagram in Fig. 10.The minus sign in the term u ∞ (1-a) is due to the retardation of flow while the air approaching the rotor and the plus sign in the term Ωr(1+a') as shown in Fig. 10 is due to the flow of air in the reverse direction of blade rotation after air particles hit the blades and so give torque.From the velocity diagram in Fig. 10, the following relationships can be determined: (1 ) ( sin( ) cos( )) sin ( ) Thus, from blade element theory, two Equations (57 and 58) have been obtained.They define the normal force (thrust) and the tangential force (torque) on an annular rotor section as a function of the flow angles at the blades and airfoil characteristics.Here it is convenient to consider in turn the following assumptions which are the bases of developing blade element theory [15]: The assumption that the behavior of an element is not affected by the adjacent elements of the same blade The airfoil characteristics to be adopted for the element

Blade Element-momentum (BEM)Theory
BEM theory refers to the determination of a wind turbine blade performance by combining the equations of general momentum theory and blade element theory.
In this case by equating the elemental thrust force equations from general momentum theory and blade element theory (Equations 28 and 57 respectively) the following relationship is obtained: 4sin ( ) and equating the elemental torque derived in both general momentum theory and blade element theory (Equations 33 and 58 respectively): [ ] Equation 60 can be rearranged by using equation 47, which relates a, a', φ and λ r based on the geometric considerations: [ ] ) 4 cos( ) In the calculation of induction factors a and a', accepted practice is to set C D zero for the purpose of determining induction factors independently from airfoil characteristics.For airfoils with low drag coefficient, this simplification introduces negligible errors [10].So equation 59, 60 and 61 can be rewritten considering By using these three equations the following useful relationships result after some algebraic manipulation; )) cos( ) (sin(

Solution Methods
Two solution methods will be proposed using these equations to determine the flow conditions and forces at each blade section.The first one uses the measured airfoil characteristics and the BEM equations to solve directly for C L and a.This method can be solved numerically, but it also lends itself to a graphical solution that clearly shows the flow conditions at the blade and the existence of multiple solutions [9].The second solution is an iterative numerical approach that is most easily extended for flow conditions with large axial induction factors.
The first method focuses in solving for C L and α.Since φ=θ+α, for a given blade geometry and operating conditions, there are two unknowns in Equation 65, C L and α at each section.In order to find these values, one can use the empirical C L vs. α curves for the chosen airfoil.One then finds the C L and α from the empirical data that satisfy Equation 65.This can be done either numerically or graphically.Once C L and α have been found, a' and a can be determined from any two of Equations 62 through 68.It should be verified that the axial induction factor at the intersection point of the curves is less than 0.5 to ensure that the result is valid.
The other equivalent solution method starts with guesses for a and a', from which flow conditions and new induction factors are calculated.Specifically: 1. Guess values of a and a'. 2. Calculate the angle of the relative wind (Equation 47). 3. Calculate the angle of attack from Equation 45 and then C L and C D .4. Update a and a' from Equations 62 and 63 or 66 and 67.
The process is then repeated until the newly calculated induction factors are within some acceptable tolerance of the previous ones.This method is especially useful for highly loaded rotor conditions, as described in Section 3.2.

Calculation of Power Coefficient
Once a has been obtained from each section, the overall rotor power coefficient may be calculated using the expression for the elemental power from Equation 34, elemental torque from Equation 58 together with Equation 36 as: Finally by using Equations 62 and 68 in equation 69, the general form of power coefficient expression can be obtained as: where λ h is the local speed ratio at the hub.Note that when, C D ≈0, the equation above for P C is the same as the one derived from the general momentum theory (Equation 37).An alternative expression for the power coefficient can be derived after performing the tedious algebra and by inserting equations 65 and 66 into equation 69; sin ( )(cos( ) sin( ))(sin( ) 8 cos( )) 1 ( / )cot( ) Note that even though the axial and angular induction factors were determined assuming the C D =0, the drag is included in the power coefficient calculation.By the same token, the thrust coefficient C T can be found beginning from the definition of thrust coefficient in the equation 2. sin ( )(cos( ) sin( ))(sin( ) 4 cos( )) 1 ( / )cot( )

HAWT BLADE DESIGN
In this chapter, the application of BEM theory on the HAWT blade design and analyzing the aerodynamic performance of a rotor will be explained.The concept of tip-loss factor, the flow states in which HAWTs are operating and the introduction of airfoil selection criteria in HAWT blade design will be discussed.After giving all these necessary knowledge for a blade design, the blade design procedure for an optimum rotor and power performance prediction procedure is given.It should be noted here that various methods [16][17][18][19][20] for HAWT blade design and predicting performance of a rotor have been studied.

The Tip-loss Factor
There was an assumption; the rotor has an infinite number of blades which was used in all theories discussed in the previous chapter.With the aid of this assumption radial velocity of the flow across the rotor plane and in the wake has been always neglected and by so the derivations of governing equations have been established.But near the boundary of the slipstream the air tends to flow around the edges of the blade tips and acquires an important radial velocity.Because the pressure on the suction side of a blade is lower than on the pressure side which causes the air to flow around the tip from the lower to upper surface, reducing lift and hence power production near the tip.
A number of methods have been suggested for including the effect of tip losses.An approximate method of estimating the effect of this radial flow and hence including the effect of tip losses has been given by L.Prandtl.The expression obtained by Prandtl for tip-loss factor is given as: The application of this equation for the losses at the blade tips is to provide an approximate correction to the system of Equations (27, 33, 62-67 and 71) summarized in the previous sections for predicting rotor performance and blade design.
Thus, to order of approximation of the analysis the correct form of the axial momentum equation (Equation 27) will be taken to be [21]; Similarly, also the angular momentum equation (Equation 33) will be assumed to be [18]; Thus the effect of the tip-loss is to reduce slightly the thrust and torque contributed by the elements near the tips of the blades.
Equations derived in Section 2.3 are all based on the definition of the forces used in blade element theory and remain unchanged.When the forces from the general momentum theory are set equal using the method of BEM theory as performed before, the derivation of the flow condition is changed.

HAWT Flow States
Measured wind turbine performance closely approximates the results of BEM theory at low values of the axial induction factors but general momentum theory is no longer valid at axial induction factors greater that 0.5, because according to the Equation 12, the wind velocity in the far wake would be negative.In practice, as the axial induction factor increases above 0.5, the flow patterns through the wind turbine become much more complex than those predicted by the general momentum theory.The thrust coefficient represented by the Equation 20 can be used to characterize the different flow states of a rotor.Fig. 11 shows flow states and thrust force vectors T associated with a wide range of axial induction factors.According to Fig. 11, for negative induction factors (a < 0) it is simple to continue the analysis to show that the device will act as a propeller producing an upwind force (i.e.C T <0) that adds energy to the wake.This is typical of the propeller state.The operating states relevant to HAWTs are designated by the windmill state and the turbulent wake state.The windmill state is the normal operating state.The windmill state is characterized by the flow conditions described by general momentum theory for axial induction factors less than about 0.5.As illustrated by the data in Fig. 11 obtained on wind turbines, above a=0.5, rotor thrust increases up to 2 with increasing induction factor in the turbulent wake state, instead of decreasing as predicted by the Equation 12.While general momentum theory no longer describes the turbine behavior, Glauert's empirical formula for axial induction factor from 0.4 to 1.0 are often used in HAWT rotor design for predicting wind turbine flow states [22].
When the induction somewhat over unity, the flow regime is called the vortex ring state and when a >2.0 the rotor reverses the direction of flow which is termed propeller brake state with power being added to flow to create downwind thrust on the rotor.
In the turbulent wake state, as stated before, a solution can be found by using the Glauert empirical relationship between the axial induction factor, a and the thrust coefficient, C T in conjunction with the blade element theory.The empirical relationship developed by Glauert, including tip losses is [21]; (1 / ) 0.143 0.0203 0.6427(0.88 ) This equation is valid for a > 0.4 or equivalently for C T <0.96.The Glauert empirical relationship was determined for the overall thrust coefficient for a rotor.It is customary, however, to assume that it applies equally to equivalently local thrust coefficients for each blade section.The local thrust coefficient C Tr can be defined for each annular rotor section as: From the Equation 57 for the elemental thrust force from blade element theory together with Equation (84), the local thrust coefficient becomes: (1 ) ( cos sin ) / sin

Airfoil selection in HAWT Blade Design
Designing HAWT blade depends on knowledge of the properties of airfoils.The most significant flow factor influencing the behavior of airfoils is that of viscosity which is characterized by the Reynolds number of the airfoil/fluid combination.The Reynolds number Re is defined by: Airfoils in use on modern wind turbines range in representative chord size from about 0.3 m on a smallscale turbine over 2 m on a megawatt-scale rotor.Tip speeds typically range from approximately 45 to 90 m/s.Then for HAWT airfoils Reynolds number range from about 0.5 million to 10 million.This implies that turbine airfoils generally operate beyond sensitive.It should be noted that there are significant differences in airfoil behavior at different Reynolds numbers.For that reason it must be made sure that appropriate Reynolds number data are available for the blade design.
The lift and drag characteristics of airfoils show also significant aspect-ratio dependence at angles of attack larger than 30 deg.But in the fully-attached region in which HAWT is operating under windmill state are not greatly affected by aspect-ratio, so that two-dimensional (i.e.infinite aspect-ratio) data can be used in blade design at low angles of attack.However, when twodimensional data are used, tip-loss factor must be added as described in equation 73.
There are evidently many engineering requirements into the selection of a wind turbine airfoil.These include primary requirements related to aerodynamic performance, structural, strength and stiffness, manufacturability and maintainability.
The usual assumption is that high lift and low drag are desirable for an airfoil and that the drag-to-lift ratio γ which is known as glide ratio as given below is a critical consideration [22]; Airfoils for HAWTs are often designed to be used at low angles of attack, where lift coefficients are fairly high and drag coefficients are fairly low(i.e.fairly low glide ratio).In this thesis, NACA4412 airfoil characteristics has been obtained and included in the airfoil database of the BLADE DESIGN program as an example.

Blade Design procedure
Designing a blade shape from a known airfoil type for an optimum rotor means determining the blade shape parameters; chord length distribution and twist distribution along the blade length for a certain tipspeed ratio at which the power coefficient of the rotor is maximum.
As it can be seen from the Equation 71, the overall power coefficient, P C depends on the relative wind angle (φ), local tip-speed ratio λ r , the glide ratio (C D /C L ) and the tip loss factor (F).To get maximum P C value from this equation is only possible make the elemental power coefficient maximum for each blade element.In other words, the term in the integral of the mentioned equation should be maximum for each blade element in order to get maximum overall power coefficient from the summation of each.
As it can be seen from the Equation 82, the overall power coefficient, P C depends on the relative wind angle (φ), local tip-speed ratio λ r , the glide ratio (C D /C L ) and the tip loss factor (F).To get maximum C P value from this equation is only possible to make the elemental power coefficient maximum for each blade element.In other words, the term in the integral of the mentioned equation should be maximum for each blade element in order to get maximum overall power coefficient from the summation of each.
For a minimum glide ratio of selected airfoil type a relationship can be established between the relative wind angle and local tip-speed ratio to determine the optimum relative wind angle, φ opt for a certain local tipspeed ratio.The condition is given as: When the optimum relative wind angle values are plotted with respect to the corresponding local tip-speed ratio values at which the elemental power coefficient is maximum for a wide range of glide ratios, the relationship has be found to be nearly independent of glide ratio and tip-loss factor.Hence the condition given before can be rearranged as: Therefore a general relationship can be found between optimum relative wind angle and local tipspeed ratio which will be applicable for any airfoil type by taking the partial derivative of the term above, i.e.: { } 2 sin (cos sin )(sin cos ) Equation 90 reveals after some algebra [18]; Having found the solution of determining the optimum relative wind angle for a certain local tipspeed ratio, the rest is nothing but to apply the equations from equation 79 to 82, which were derived from the blade-element momentum theory and modified including the tip loss factor, to define the blade shape and to find out the maximum power coefficient for a selected airfoil type.
The procedure of blade design begins with dividing the blade length into N elements.The local tip-speed ratio for each blade element can then be calculated with the use of Equation 31 as given below: Then according to Equation 91 the optimum relative wind angle for each blade element is determined as: From the Equation 73 the tip loss factor for each blade element can be found as: The chord-length distribution can then be calculated for each blade element by using Equations 79 and 56 as: (sin cos ) where, C L,design , is chosen such that the glide ratio is minimum at each blade element.The solidity can be calculated as follows: The twist distribution can be determined from the following equation which can be easily derived from the velocity diagram shown in Figure 9 as: where, α design , is again the design angle of attack at which C L,design is obtained.Finally, the power coefficient is determined using a sum approximating the integral in Equation 82 as: As a result, for a selected airfoil type and for a specified tip-speed ratio and blade length (i.e.rotor radius), the blade shape can be designed for optimum rotor and from the calculation of power coefficient the maximum power that can be extracted from the wind can then be found for any average wind velocity.During the design process, from the required power and wind velocity, the rotor radius, R can be estimated by including the effect of probable C P and efficiencies.Then, choose the design tip speed ratio and the number of blades according to the type of application of the wind turbine and then select the design aerodynamic conditions, C L,design and α design , such that C D,design /C L,design is at a minimum for each blade section.Finaly, divide the blade into N elements and use the optimum rotor theory to estimate the shape of the i th blade with a midpoint radius of r i and Equations (92-98).

SAMPLE BLADE DESIGN ON BLADE DESIGN PROGRAM
As it was stated before, the main objective of this paper is to develop a user interface computer program on HAWT blade design.The main purpose of constructing such a program is to collect and represent all the studies of this paper in a visual and more understandable way and also to provide a design program for the users who are dealing with HAWT blade design.The general view of the program and its working flow schematic are illustrated in Fig. 12 and Fig. 13 respectively.The main features of the program can be summarized as below: The program was written on VISUALBASIC.
Input part requires the following information; turbine power required, design wind velocity, the number of blades and airfoil type to be used for the blade profile.The output parts are composed of design condition output, blade geometry output, 3-D visualization output and figures output.The design condition output part gives the information for designed blade.These information are design power coefficient, corresponding design tip-speed ratio and rotor diameter.Blade geometry output part shows the chord-length in meter and setting angle in degree for each blade element.The 3-D visualization output part gives three dimensional views of the designed blade on after exporting to AutoCAD windows.
A sample blade design is demonstrated on blade design program according to the following inputs; Turbine power required, P=1000W, Design wind velocity, u ∞ =8 m/s, Number of blades, B=3, Airfoil type -NASA63206.The first task is to find the design power coefficient which is dependent on the selected airfoil type.Then the design tip-speed ratio at which the maximum power coefficient (i.e.design power coefficient) can be obtained for the three-bladed rotor.This design tipspeed ratio can be used to determine the rotational velocity of the rotor at design condition.Finally rotor diameter is found out with the known values of design power coefficient and design wind velocity and required power.
Wind turbines are most commonly classified by their rated power at a certain rated wind speed however annual energy output is actually a more important measure for evaluating a wind turbine's value at a given site.Multiplying the rated power output by the rough capacity factor and the number of hours in a year, annual energy production can be estimated.Capacity factor is nothing but the wind turbine's actual energy output for the year divided by the energy output if the machine operated at its rated power output for the entire year.If differences < 0.01 Take updated a and a' as actual values and for considered blade element and calculate C P (-), dT (KNm) and dQ (KNm)  In Fig. 13 the chord-length distribution with respect to radial location of each blade element both of which are normalized with blade radius is shown.Similarly twist distribution with respect to radial location is illustrated in the same figure.The blade chord-length and twist distribution for an optimum three-bladed rotor at the design tip-speed ratio λ d =4 is tabulated in Table 1 for the airfoil NACA4412 whose lift coefficient and drag coefficient values are taken at Re=1x10 6 .As can be seen from Fig. 13, in the output part the blade design program, the radius of the rotor, rotational speed of the rotor, the trust coefficient and power coefficient are given.Finally, Fig. 14 shows the blade elements of NACA 4412 in isometric view for the designed blade.

CONCLUSION
In this study, aerodynamic design of horizontal-axis wind turbine blades was investigated and a user-interface computer program called as BLADE DESIGN PROGRAM was written for the use of its outputs in further studies.All the studies on HAWT blade design were presented on a user interface computer program written on VISUALBASIC.User of the blade design program gives the required power output, number of blades, blade profile and design wind velocity as input and the program gives design power coefficient, design tip-speed ratio and rotor diameter as output.In addition to that, blade geometries will be listed.Three dimensional views of the blade (the blade model) were what the program does lastly after exported to AUTOCAD.
Structural design of HAWT blades is as important as their aerodynamic design.The dynamic structural loads which a rotor will experience play the major role in determining the lifetime of the rotor.Obviously, aerodynamic loads are a major source of loading and must be well understood before the structural response can be accurately determined and also the blade geometry parameters are required for dynamic load analysis of wind turbine rotors.So such a study on the dynamic load analysis of HAWT blades might also use the outputs of BLADE DESIGN PROGRAM.

Figure 1 .
Figure 1.Some of the erected wind turbines in Adama I wind farm in Ethiopia

Figure 3 .
Figure 3. Idealized flow through a wind turbine represented by a non-rotating actuator disc flow, air mass flow rate through the disk can be written as:

Figure 4 .
Figure 4. Velocity and Pressure distribution along the stream tubeThe power output, P is equal to the thrust times the velocity at the rotor plane:

Figure 5 .
Figure 5. Operating parameters for a Betz turbine

Figure 6 .
Figure 6.Stream tube model of flow behind rotating wind turbine blade

Figure 7 .
Figure 7. Geometry for rotor analysis

Figure 8 .
Figure 8. Theoretical maximum power coefficient as a function of tip speed ratio for an ideal HAWT with and without wake rotation

Figure 9 .
Figure 9. Schematic of blade elements

Figure 10 .
Figure 10.Blade geometry for analysis of a HAWT the rotor has B number of blades, the total normal (trust) and tangential force on the element at a distance r by rearranging Equation 50 and 51 with the use of equation 46, 48 and 49: due to the tangential forces, dL operating at a distance r from the center is given by: And by defining solidity ratio, σ as following; 46 and 56 into Equations 52 and 55, the general form of elemental thrust and torque)

Figure 11 .
Figure 11.Relationship between the axial induction factor, flow state and thrust of a rotor [19]

Figure 12 .
Figure 12.Flow chart of the iteration procedure for HAWT blade design

Figure 13 .
Figure 13.Output from the program performed for the sample blade design case Table 1.Variation of blade chord, twist distribution, normal (thrust) force and torque along the blade r(m) c (m) Twist (deg.)Trust (N) Torque (Nm)