Multi-objective optimization of the airfoil shape of Wells turbine used for wave energy conversion

Energy 36 (2011) 438e446

Contents lists available at ScienceDirect

Energy journal homepage: www.elsevier.com/locate/energy

Multi-objective optimization of the airfoil shape of Wells turbine used for wave energy conversion M.H. Mohamed, G. Janiga, E. Pap, D. Thévenin* Lab. of Fluid Dynamics and Technical Flows, University of Magdeburg “Otto von Guericke”, Germany

a r t i c l e i n f o

a b s t r a c t

Article history: Received 7 July 2010 Received in revised form 1 October 2010 Accepted 9 October 2010 Available online 10 November 2010

Wells turbine is one of the technical systems allowing an efficient use of the power contained in oceans’ and seas’ waves with a relatively low investment level. It converts the pneumatic power of the air stream induced by an Oscillating Water Column into mechanical energy. The standard Wells turbines show several well-known disadvantages: low tangential force, leading to low power output from the turbine; high undesired axial force; usually a low aerodynamic efficiency and a limited range of operation due to stall. In the present work an optimization process is employed in order to increase the tangential force induced by a monoplane Wells turbine using symmetric airfoil blades. The automatic optimization procedure is carried out by coupling an in-house optimization library (OPAL (OPtimization ALgorithms)) with an industrial CFD (Computational Fluid Dynamics) code (ANSYS-Fluent). This multi-objective optimization relying on Evolutionary Algorithms takes into account both tangential force coefficient and turbine efficiency. Detailed comparisons are finally presented between the optimal design and the classical Wells turbine using symmetric airfoils, demonstrating the superiority of the proposed solution. The optimization of the airfoil shape leads to a considerably increased power output (average relative gain of þ11.3%) and simultaneously to an increase of efficiency (þ1%) throughout the full operating range. Ó 2010 Elsevier Ltd. All rights reserved.

Keywords: Wells turbine Wave energy Multi-objective optimization Evolutionary algorithms

1. Introduction In its recent report the United Nations IPCC (“International Panel on Climate Change”) emphasized the challenges associated with climate change and recommended a fast migration toward renewable energy sources. In the US, the European Union and China, policies have been formulated with the objective of rapidly decreasing CO2 emissions. In many nations around the world, policies are being initiated to raise the share of renewable energy as fast as possible. In March 2007, the European Union defined a target of 20% renewable energy for year 2020. In Denmark, a target of 30% renewable energy for year 2025 has just been proposed by the Danish Government [1]. There is a growing perception by society of the risks of dramatic global climate changes due to anthropogenic greenhouse gases, in particular energy related emissions of CO2. This has spurred a renewed interest in carbon-free or carbon-neutral technologies for converting sources of renewable primary energy to electricity and to transportation fuels. However, it takes energy to produce * Corresponding author. Tel.: þ49 391 67 18570; fax: þ49 391 67 12840. E-mail addresses: [email protected] (M.H. Mohamed), janiga@ ovgu.de (G. Janiga), [email protected] (E. Pap), [email protected] (D. Thévenin). 0360-5442/$ e see front matter Ó 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.energy.2010.10.021

energy, even when the primary source is energetically cost-free, such as solar or wind [2]. Wave energy is rightly regarded as one of the renewable energy sources with the greatest potential for development over the course of the next few years. Several technological solutions are currently the object of intensive research. Yet the technological challenge, for all its importance, is not the only requirement for the successful exploitation of wave energy. A thorough resource assessment is also necessary in view of the significant spatial variations of exploitable wave energy. The Atlantic coastlines of Western Europe present a large wave energy resource for two main reasons: their location at the eastern boundary of the Atlantic Ocean and the prevalence of western winds at mid-latitudes, which drive a succession of lows across the ocean that in turn generate powerful swells. As a result, Western European countries figure prominently in wave energy atlases and maps, with annual average wave power exceeding 40 kW m1 in certain regions of Norway, Ireland, UK, France, Spain and Portugal [3,4]. A previous study presented a methodology to assess the possible benefits of the combination of wind energy with the still unexploited, but quite significant in Ireland, wave energy. An analysis of the raw wind and wave resource at certain locations around the coasts of Ireland shows low correlations on the South

M.H. Mohamed et al. / Energy 36 (2011) 438e446

439

R

FA

Rotation

w

vA ut

Flow Flow ut Fig. 1. Operation principle of the Wells turbine, the Oscillating Water Column.

vA

and West Coast, where the waves are dominated by the presence of high energy swells generated by remote westerly wind systems. As a consequence, the integration of wind and waves in combined farms, at these locations, allows the achievement of a more reliable, less variable and more predictable electrical power production [5]. >Dr. Wells, a former professor of civil engineering at Queen’s university of Belfast, proposed in 1976 a self-rectifying axial flow air turbine as a device suitable for wave energy conversion using the OWC (Oscillating Water Column, Fig. 1). In its standard design the air turbine rotor consists of several symmetrical airfoil blades positioned around a hub. Because of its simple and efficient operation, the Wells turbine has already been applied in practice to gain energy from oceans’ waves. For example, one such wave power station with an installed peak capacity of 70 kW has been built in the west of England under the support of the Department of Energy as early as 1991. It employs a twin-rotor symmetrical wing turbine generator. The observed, average electric power generation is 7.5 kW. The Limpet shoreline wave energy concept has been commissioned in December 2000 on the Island of Islay, off the west coast of Scotland. It is intended to enable Islay to replace fossil fuels and become self-sufficient through renewable energy. The OWC feeds a pair of counter-rotating Wells turbines, each of which drives a 250 kW generator, giving a theoretical peak power of 500 kW [6]. The Mighty Whale, a floating OWC, has been installed 1.5 km from the shore in 40-m depth water, just outside the mouth of Gokasho Bay (off the Mie Prefecture) in Japan. The overall peak power capacity is expected at 110 kW [7].

FA

  i h CT ¼ FT = ð1=2Þr v2A þ u2t zbc

(1)

R

  i h CA ¼ Dp0 prt2 = ð1=2Þr v2A þ u2t zbc

(2)

where ut ¼ urt is the peripheral velocity, u is the rotor angular velocity and rt is the tip radius. Furthermore, vA is the axial velocity normal to the plane of rotation, z is the number of blades, b is the blade span, c is the blade chord (see Fig. 3) and Dp0 is the total pressure difference across the rotor. The main advantage of this turbine is that it always rotates in the same direction, normal to the free stream, although the air flow is oscillating. The manufacturing of Wells turbines is easier than for impulse turbines and they are able to work at low flow coefficients. Furthermore, OWC-based energy stations are considerably more compact than other known converters [7]. On the other hand, Wells turbines show some drawbacks compared to alternative systems: low efficiency, poor self-starting characteristics, short operating range, high axial force coefficient and low tangential force coefficient. Therefore, many authors have already tried to identify the best principle of operation and to improve the characteristics of Wells turbines (see in particular [9e11]). Several ideas have been introduced for this purpose. It has been first suggested to place guide vanes before and after the rotor [12e14]. As an alternative, self-pitch controlled blades [15e17] and end plates [18] have been used to improve the performance of the Wells turbine. Two-stage Wells turbines involving symmetric

Rotor hub on shaft

)

rh Sym. Airfoil Blade

(b

rt

Many research programs attempting to gain energy from waves depend on the OWC [8] as converter mechanism (Fig. 1). The water wave energy is first converted to pneumatic energy in the air, which then passes periodically across a self-rectifying, axial air flow turbine, called Wells turbine. The turbine comprises a number of symmetric airfoils set around the hub radially at 90 stagger angle, with the chord plane normal to the axis of rotation. According to the standard airfoil concept, if the airfoil is set at an angle of incidence a in a fluid flow, it will generate a lift force FL normal to the free stream and a drag force FD in the direction of the free stream. These lift and drag forces can then be combined to get the tangential force FT and the axial force FA as shown in Fig. 2. The tangential force coefficient CT and the axial force coefficient CA are finally calculated as:

Rotation

Fig. 2. Axial and tangential forces acting on a Wells turbine.

sp an

2. Principle of operation

w

Rotation Chord (c)

FA Oscillating flow

FA

FT Forces resolved in direction of rotation

Fig. 3. Schematic description and main characteristic dimensions of a Wells turbine blade.

440

M.H. Mohamed et al. / Energy 36 (2011) 438e446

airfoils have been investigated experimentally and theoretically [19e21] and show an improved performance compared to monoplane turbines. A two-stage Wells turbine constructed from nonsymmetric airfoils based on National Advisory Committee for Aeronautics (NACA) 2421 has also been optimized, mirroring the two stages in order to keep the global symmetry of the turbine [22,23]. Due to the complexity and cost of suitable experiments, the trend is first to predict numerically the performance of Wells turbines [24,25] in order to support new developments. Some of the most important results concerning possible improvements are summarized in Table 1.

3. Purpose of the present work Apart from our own contribution [22,23,26], all the theoretical and experimental investigations listed in the previous section only consider the performance of Wells turbines using standard symmetric airfoils of type NACA 00XX. As an illustration, Fig. 4 shows NACA 0015 and NACA 0021. Most investigations pertaining to Wells turbines have considered NACA 0012, NACA 0015, NACA 0018 and NACA 0021 (e.g., [10,28,29]). The formula for the shape of a NACA 00XX foil, with “XX” being replaced by the percentage of maximum thickness to chord length c, is

y ¼

rffiffiffi  x ct x 0:2969  0:126 0:2 c c x3 x4  x2 þ0:2843 0:1015  0:3517 c c c

(3)

where x is the position along the chord from 0 to c, y is the halfthickness at a given value of x (centerline to external surface), and t is the maximum half-thickness as a fraction of the chord (so that 100t gives a half of the last two digits in the NACA 4-digit denomination). Reference investigations indicated that NACA 0021 airfoil profiles (21% thickness) lead to the best performance for conventional monoplane Wells turbines [10]. There is nevertheless no proof that NACA profiles, as defined by Eq. (3), automatically lead to the best possible performance. An alternative geometry might be much better, in particular for such very specific applications. As a consequence, the present work now concentrates on the optimization of a symmetric airfoil shape, leading to the best possible performance of a Wells turbine (i.e., maximal tangential force coefficient and efficiency). Due to the complexity of the underlying optimization procedure, this first study considers only monoplane Wells turbines (the original design) and a constant turbine solidity. We are indeed well aware of the interest of two-stage Wells turbines [22,26], but Table 1 Main modifications allowing to improve the performance of Wells turbines. Design modification

Gain

Description and comments

Contra-rotating [27]

Improve efficiency by x5% Improve efficiency by x7% Improve efficiency by x5% Improve efficiency by x5% Wider operating range Improve efficiency by x2% Improve efficiency by x0.6%

Double shaft, complex

Pitch setting [15] Guide vanes [14] End plate [18] a

Multi-stage [28] Multi-stageb [22, 23] Multi-stagea [26] a b

Symmetric airfoils. Non-Symmetric airfoils.

such an optimization would require considerably higher computational resources and will be the subject of future publications. 4. Numerical flow simulations It is generally accepted that the air flow frequency in a wave energy device is so small (f < 1 Hz) that dynamic effects are negligible [30]. Steady-state computations are appropriate for such cases. All flow simulations (later denoted also CFD for Computational Fluid Dynamics) presented in this work rely on the industrial software ANSYS-Fluent 6.3. The steady-state ReynoldsAveraged NaviereStokes equations are solved using the SIMPLEC (Semi-Implicit Method for Pressure-linked Equations: Consistent) algorithm for pressure-velocity coupling. The flow variables, the turbulent kinetic energy and the turbulent dissipation rate are discretized in a Finite-Volume formulation using a second-order upwind scheme, while the pressure field is solved using the Presto scheme. The validation and accuracy of these simulations is discussed in a later section. 5. Optimization The mathematical optimization of turbomachines is still a relatively new field of research. Optimization is used to identify geometries and configurations that best satisfy some specified requirements, pertaining here to the performance of the monoplane Wells turbine. Optimization criteria can be considered both for mechanical and for flow properties, as done in the present work. The central goal when designing an improved Wells turbine is to achieve high efficiency and high power output (i.e., high tangential force coefficient). Note that many turbomachine bladings are still designed using direct optimization methods, as in the present project. This situation results mainly from the fact that, 1) the inverse problem is much more difficult to properly formulate as well as to solve than the direct one due to the presence of unknown boundaries, which is the origin of strong nonlinearity and possible ill-posedness, and 2) the inverse design method often leads to airfoil/cascade shapes that are unfeasible from practical design considerations [31,32]. Furthermore, inverse methods are difficult to implement for multiobjective, concurrent optimization problems. Though highly interesting in principle and increasingly promising for future applications, inverse methods are therefore not discussed further in the present article. Readers specifically interesting in inverse methods for turbomachines are referred for instance to [33e36]. It must be kept in mind that turbomachines often operate outside the nominal (or design) conditions. Therefore, the present work considers a typical operating range: an angle of incidence a varying between 5 and 14 ; or equivalently a flow coefficient, i.e., the ratio between axial flow velocity and peripheral velocity of the blade, f ¼ tanðaÞvarying between 0.08 and 0.25. Due to an

For positive small angles, complex Smaller operating range Only for 0:2  f  0:25 Reduce efficiency by x10% Small parameter space Only symmetric airfoils

Fig. 4. Standard airfoils NACA 0015 and NACA 0021.

M.H. Mohamed et al. / Energy 36 (2011) 438e446

441

Table 2 Parameters of the evolutionary algorithm. Parameter

Value

Population size of the first generation, N Number of generations Survival probability Average probability Crossover probability Mutation probability Mutation magnitude

30 40 50% 33.3% 16.7% 100% 30%a (i.e. 15%)

a This value is multiplied by 0.8 at each generation. For example the mutation magnitude is only 4% (2%) after 10 generations. Mutation magnitude must be decreased during the optimization process to stabilize the population.

extremely small power output, there is no practical interest for the operating range below 5 . Attempting operation above a ¼ 14 means working beyond stall conditions, where efficiency and power rapidly decrease. Many publications hence recommended the range between 5 and 14 (e.g., [10,15,21]) and these values have been retained as well in the present study. Optimization methods attempt to determine the n design variables Xi ði ¼ 1.nÞ that maximize a user-defined objective function, denoted OFðUðXi Þ; Xi Þ, where UðXi Þ is here the CFD solution of the flow equations [37]. In the present paper, the free design variables considered for the optimization will be the shape of the blade using a constant solidity ðs ¼ zc=½prt ð1 þ hÞÞ where h is the ratio between hub radius rh and tip radius rt. The objective function contains simultaneously two outputs of the simulation, that should both be maximized as far as possible: the tangential force coefficient CT; and the turbine efficiency h, which is inversely proportional to the axial force coefficient, and is defined for a negligible density change as:

h¼

F T ut

Dp0 Q

(4)

with Q the volumetric flow-rate through the turbine. 6. Optimization methodology Until recently, the denomination “optimization” was mostly used in the engineering literature to describe a trial-and-error, manual procedure, at the difference of a real, mathematical optimization. This is now changing rapidly. In the present project, direct

Fig. 5. Grid-independence study for the tangential force coefficient.

Fig. 6. Influence of the turbulence model on the tangential force coefficient, compared to experimental results [43].

optimization will be considered to obtain the optimal shape geometry. In this case, an appropriate algorithm must be chosen. A considerable experience is available in our group concerning the mathematical optimization relying on CFD evaluations [37]. We therefore employ our own optimization library, OPAL (OPtimization ALgorithms), containing many different optimization techniques. OPAL has already been coupled in the past with different CFD solvers (in-house codes, ANSYS-Fluent, ANSYS-CFX) and has been employed successfully to improve a variety of applications, for example heat exchangers [38], burners [39] or turbomachines [22,40e42]. For configurations involving concurrent objectives, Evolutionary Algorithms are particularly robust and have therefore been used in the present study. The employed optimization parameters are listed in Table 2. A fully automatic optimization takes place, using OPAL (decisionmaker for the configurations to investigate), the commercial tool Gambit for geometry and grid generation (including quality check) and the industrial CFD code ANSYS-Fluent to compute the flow field around the modified Wells turbine. As a result of the CFD computation the tangential force coefficient and the efficiency are determined, and both are stored in a result file. The procedure is automated using journal scripts (Gambit, Fluent) and a master program written in C, calling all codes in the right sequence. By checking the values stored in the result file, OPAL is able to decide how to modify the rotor geometry, before starting a new iteration. The fully coupled optimization procedure is a complex task, which has been described in detail in previous publications. We thus refer the interested reader to [37] for a complete description of the procedure.

Fig. 7. Allowed parameter space for the moving points P2 to P12.

442

M.H. Mohamed et al. / Energy 36 (2011) 438e446

Table 3 Parameter space for the moving points P2 to P12. Point

Parameter

Minimum

Maximum

P2

XP2 =c YP2 =c XP3 =c YP3 =c XP4 =c YP4 =c XP5 =c YP5 =c XP6 =c YP6 =c XP7 =c YP7 =c XP8 =c YP8 =c XP9 =c YP9 =c XP10 =c YP10 =c XP11 =c YP11 =c XP12 =c YP12 =c

0.01 0.023 0.015 0.025 0.05 0.034 0.1 0.043 0.2 0.054 0.35 0.052 0.45 0.043 0.55 0.039 0.65 0.034 0.75 0.026 0.85 0.015

0.015 0.043 0.035 0.066 0.1 0.11 0.2 0.143 0.3 0.154 0.45 0.152 0.55 0.143 0.65 0.11 0.75 0.094 0.85 0.066 0.95 0.035

P3 P4 P5 P6 P7 P8 P9 P10 P11 P12

a

0.6

Efficiency

0.5

0.4

0.3 Optimization output Monoplane Wells turbine NACA 0021

0.2

0.1 0

0.04

0.08

0.12

0.16

Tangential force coefficient (CT)

b

(B)

(A)

0.52

Efficiency

Fig. 9. Input parameters of the optimization and objectives represented using parallel coordinates. The parameters of the optimal shape are connected with an orange line going through gray circles. The standard design (NACA 0021) is shown with a blue dashed line going through gray squares; a) X-coordinates of the variable points (P2.P12); b) Ycoordinates of the variable points (P2.P12) (for interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article).

0.525

0.515

7. Results and discussions 0.51

CFD-based optimization can only deliver quantitative, meaningful results, when the evaluations of the different configurations carried out with CFD are of high quality. Therefore, we first discuss the quality of the CFD model, starting with a grid-independence study.

0.505 Optimization output Optimum tang. force coeff. Optimum efficiency

0.5 0.495 0.1275

0.129

0.1305

0.132

0.1335

Tangential force coefficient (CT) Fig. 8. Objectives of the optimization; a) for all computed configurations b) for the best configurations (i.e., zoom on the upper-right part (dashed square) of a).

7.1. Grid dependency Several different two-dimensional grids of increasing density and quality, composed of 2120 up to 119,000 cells, have been tested for the baseline, symmetric blade configuration NACA 0021. All

M.H. Mohamed et al. / Energy 36 (2011) 438e446

other parameters of the CFD are unchanged. Corresponding results are shown in Fig. 5. It is easy to notice that the nine coarsest grids are associated with a large variation of the objective functions (here, the tangential force coefficient is represented). On the other hand, all remaining grids employing more than 60,000 cells lead to a variation of the target variables smaller than 1.08%. Since the cost of a CFD evaluation obviously increases rapidly with the number of grid cells, the intermediate grid range between 60,000 and 80,000 cells is retained for all further results shown in the present paper.

In a second step, the full numerical model and in particular the employed turbulence model are validated by comparison with published experimental results for a standard, monoplane Wells turbine at a flow Reynolds number Re ¼ 2.4 105 using the chord as characteristic length. Published studies usually consider a range Re ¼ 1 105 to Re ¼ 5 105, since this corresponds to realistic conditions for employing a Wells turbine. As already stated before, solidity is constant and equal to s ¼ 0.67, which corresponds to a number of blades equal to 8. This solidity value has been retained since it has been proved in a number of studies to be highly suitable for Wells turbines (for instance [15,16,21]). Furthermore, it has been shown in a previous study of our group [23] that the optimal solidity is found for 0:64  s  0:70 in this configuration, taking into account the mutual interaction effect between the blades. Considering this small span, it was decided to keep the usual value of 0.67 as fixed. The influence of the turbulence model is shown in Fig. 6. These results demonstrate the excellent agreement obtained between CFD and experiments for this standard configuration, in particular when using the Realizable k-3 turbulence model. This model thus appears suitable to predict the performance of the turbine in the later investigated range of operation (flow coefficient varying from f ¼ 0.08 to 0.25) and is now kept for all further simulations. In this particular case the Reynolds-Stress Model (RSM) leads to a considerably longer computing time but surprisingly to a worse agreement than the k-3 models. This is probably due to the low turbulence level and to a larger influence of the inflow turbulence boundary conditions, which are not characterized in the experiment. Inlet boundary conditions for the RSM model have been implemented using different possibilities, prescribing either k and 3 or the Table 4 Optimum shape parameters. Point

Parameter

Value

P2

XP2 =c YP2 =c XP3 =c YP3 =c XP4 =c YP4 =c XP5 =c YP5 =c XP6 =c YP6 =c XP7 =c YP7 =c XP8 =c YP8 =c XP9 =c YP9 =c XP10 =c YP10 =c XP11 =c YP11 =c XP12 =c YP12 =c

0.0124947 0.033155 0.0269924 0.045006 0.0678583 0.0608499 0.155326 0.078402 0.248998 0.0877814 0.405838 0.0751153 0.469922 0.0666932 0.56625 0.0493063 0.723109 0.037401 0.798623 0.0315091 0.8891079 0.0214618

P4 P5 P6 P7 P8 P9 P10 P11 P12

Optimum shape Polynomial fit of optimum shape NACA 0021

0.12c 0.08c 0.04c 0 0

0.4c

0.2c

0.8c

0.6c

c

Fig. 10. Comparison between the original profile NACA 0021 (solid line), the optimal airfoil shape described by splines (black squares showing the position of the control points) and the corresponding polynomial fit (Eq. (5), dashed line).

7.2. CFD validation

P3

443

turbulence intensity together with a length scale. Nevertheless, it has been impossible to obtain a better agreement [22]. Therefore, the RSM model appears to be inappropriate for this configuration, associated with a low but unknown inflow turbulence level. 7.3. Optimization of airfoil shape After having checked the accuracy of an individual evaluation relying on CFD, it is now possible to start the optimization procedure. As explained previously, only symmetric blades are considered in what follows, based on the profile NACA 0021 for a first guess. To illustrate the optimization, a fixed angle of incidence a ¼ 8 (flow coefficient f ¼ 0.14) is considered. Twenty two free parameters are considered simultaneously by the OPAL optimizer, explaining the difficulty of the process. Increasing the number of free parameters leads obviously to a corresponding increase in algorithm complexity and computational time. In the present case the outer boundary of the airfoil (or airfoil shape) is constructed with thirteen points; two fixed points (P1 and P13) and eleven variable points (P2, P3, P4, P5, P6, P7, P8, P9, P10, P11 and P12) as shown in Fig. 7. Knowing the exact position of these 13 points, the full profile is finally reconstructed for one face of the airfoil using standard splines (Nonuniform rational B-splines, NURBS). In the present version of the software, only such splines might be used to parametrize the geometry. The order of a NURBS curve defines the number of nearby control points that influence any given point on the curve. The curve is represented mathematically by a polynomial of degree one less than the order of the curve; this means that the spline order is 13 in our case. Then, the obtained face is mirrored to obtain the full symmetric airfoil. Every point P2 to P12 has two coordinates (Xpi, Ypi), where i ¼ 2.12. The parameter space considered in the optimization has been defined as documented in Table 3 and illustrated in Fig. 7. The corresponding parameter spaces have been selected to cover all usual NACA 00XX, while avoiding collisions between reference points and keeping acceptable geometries. The reference point is point P1 (0,0), origin of the cartesian coordinate system. As a whole, the optimization process thus involves twenty two parameters (or degrees of freedom) Xpi and Ypi with i ¼ 2.12 and two objectives (efficiency and tangential force coefficient) that should be simultaneously maximized in a concurrent manner. The results presented in Fig. 8(a) and (b) indicate that the two considered objectives are indeed considerably influenced by the airfoil shape. Fig. 8(a) shows all evaluation results. As a whole, 615

Table 5 Polynomial coefficients of optimal airfoil shape (best fit). A

B

D

E

H

K

1.6958588

5.4277515

6.560073

3.737973

0.898334

0.01409362

444

M.H. Mohamed et al. / Energy 36 (2011) 438e446

Fig. 11. Static pressure distribution (in Pa) around a) the standard airfoil (NACA 0021) and b) the optimum configuration, together with Cp along the blade at different flow coefficients f.

different configurations have been finally tested by the optimizer, leading to 15 days of total computing time on a standard PC. In Fig. 8(a) the performance of the standard airfoil NACA 0021 (tangential force coefficient and efficiency (CT,h)¼(0.11632; 0.5109)) is also plotted for comparison. Globally, the two considered objectives are not fully concurrent but increase simultaneously, which is not a complete surprise since the tangential force appears on the numerator in Eq. (4) defining the efficiency. When considering now only the best configurations of Fig. 8(a), located in the upper-right corner (marked by a black square in dashed line), a more complex picture appears, as documented in Fig. 8(b). For the last percent of performance improvement, the two objectives (tangential force coefficient and efficiency) become indeed slightly concurrent and cannot be optimized simultaneously. Similar effects have been observed already in previous studies [22,23,26]. Two optimal conditions are finally found: Point A, ðCTA ; hA Þ ¼ ð0:1325; 0:519Þ (highest tangential force); and Point B, ðCTB ; hB Þ ¼ ð0:1281; 0:520Þ (highest efficiency). By analyzing in detail the resulting geometries and considering daily engineering purposes, the increase in tangential force coefficient (higher power output) appears to be more significant and valuable than the very

slightly increased efficiency. Therefore, the most interesting point is globally Point A with ðCTA ; hA Þzð0:1325; 0:519Þ. The results of the optimization process can be usefully visualized in a different manner using parallel coordinates (Fig. 9). Here again, the performance of the standard airfoil NACA 0021 is also plotted for comparison, close to the middle of the parameter space (thick dashed gray line). Fig. 9(a) and (b) indicate by parallel coordinates the X and Y coordinates of the eleven moving points (P2.P12), together with the two objectives. This figure demonstrates that very different shapes have been evaluated on the way toward the optimal solution. The optimum configuration (Point A) is indicated by a thick gray line. It can be seen that the optimal airfoil shape leads only to a slightly higher efficiency (þ0.78%) compared to the standard airfoil (NACA 0021). However, the tangential force coefficient CT is at the same time increased by 0.0162. This means a relative increase of 12.2% for the present flow coefficient, equal to 0.14. The geometrical parameters corresponding to the optimal shape are listed in Table 4. The resulting shape of the optimal airfoil in comparison with the standard NACA 0021 is shown in Fig. 10.

M.H. Mohamed et al. / Energy 36 (2011) 438e446

0.5

50 % NACA 0021 Spline optimum shape Fitting optimum shape %Relative increase

0.4

40

0.3

30

0.2

20

0.1

10

0 0.05

0.1

0.15 0.2 Flow coefficient

0.25

0.3

0.54 0.52 0.5

Efficiency

tangential force coefficient is observed throughout for all values of f, compared to the conventional turbine based on standard airfoils NACA 0021. The absolute gain for CT increases even slightly with f. The relative increase is higher than 8.8% throughout the useful operating range, with an average gain of 11.3% (Fig. 12 a). At the same time the efficiency of the optimized shape is always higher than for the conventional design, the difference being lower for large flow coefficients. The corresponding gain varies between 0.2% and up to 3.2%, with an average increase of 1% (Fig. 12 b). No significant difference is observed in Fig. 12 between the performance of the exact profile described by splines and the associated polynomial fit, Eq. (5). 8. Conclusions

0

b

% Relative increase

Tangential force coefficient (CT)

a

445

0.48 0.46 0.44

NACA 0021 Spline optimum shape Fitting optimum shape

0.42 0.4 0.05

0.1

0.15 0.2 Flow coefficient

0.25

0.3

Fig. 12. Performance of the spline optimal configuration (black empty square), fitting optimal one (Eq. (5), black cross) compared to the conventional Wells turbine relying on the NACA 0021 profile (filled squares). The corresponding relative increase is shown with stars; a) tangential force coefficient; b) efficiency.

Knowing all points P1 to P13, the full profile is again reconstructed using standard splines of order 13. Nevertheless, a simple polynomial description of this profile would be helpful for practical purposes. Furthermore, such polynomials efficiently remove possible surface oscillations induced by the spline-based description of the geometry. An excellent fit (average residual error of 0.38%) has been obtained with following polynomial description:

 3  2    5  4 Y X X X X X þK þB þD þE þH ¼ A c c c c c c

(5)

with the constants A to H listed in Table 5. The static pressure distribution around the optimal airfoil described by Eq. (5) is shown in Fig. 11 compared with the pressure distribution around the standard NACA airfoil for different values of f, showing large modifications. In particular, the obtained pressure difference is much larger with the optimal geometry. Finally, it is important to check how the gain induced by the new airfoil shape would change as a function of the flow coefficient f, since such a turbine must be able to work also for off-design conditions. Therefore, the performance of the optimal shape has been finally computed for the full range of useful f-values, as shown in Fig. 12. These results demonstrate that the improvement of

This work investigates the aerodynamic performance of a monoplane Wells turbine consisting of symmetric airfoils, using numerical simulations. Due to the potential importance of the airfoil shape, an aggressive mathematical optimization procedure has been carried out considering simultaneously twenty two free parameters (coordinates of eleven geometrical points employed to characterize the shape). Two concurrent objectives (efficiency and tangential force coefficient) have been maximized in a concurrent manner. The optimization relies on Evolutionary Algorithms, all geometrical configurations are evaluated in an automatic manner by CFD, taking into account the influence of the mutual interaction between the blades. This optimization procedure is able to identify a considerably better configuration than the standard design relying on NACA 0021. A relative increase of the tangential force coefficient exceeding 8.8% (as a mean, 11.3%) is obtained for the full operating range. At the same time, the efficiency improves also by at least 0.2% and up to 3.2% (as a mean, 1%). In future, these findings will be checked by experimental measurements in the wind tunnel at our institute. Furthermore, a similar optimization process involving two-stage Wells turbines with non-symmetric airfoils will be attempted. Acknowledgments The Ph.D. work of Mr. Mohamed is supported financially by a bursary of the Egyptian government. References [1] Lund H, Mathiesen BV. Energy system analysis of 100% renewable energy systems, the case of Denmark in years 2030 and 2050. Energy 2009;34 (5):524e31. [2] Gonçalves da Silva C. The fossil energy/climate change crunch: can we pin our hopes on new energy technologies? Energy 2010;35(3):1312e6. [3] Iglesias G, Carballo R. Wave energy potential along the death coast (Spain). Energy 2009;35(5):1963e75. [4] Iglesias G, Carballo R. Offshore and inshore wave energy assessment: asturias (N Spain). Energy 2010;35(11):1964e72. [5] Fusco F, Nolan G, Ringwood JV. Variability reduction through optimal combination of wind/wave resources: an Irish case study. Energy 2010;35 (1):314e25. [6] Folley M, Curran R, Whittaker T. Comparison of LIMPET contra-rotating Wells turbine with theoretical and model test predictions. Ocean Eng 2006;33 (8e9):1056e69. [7] Brooke J. Wave energy conversion. Elsevier Science and Technology; 2003. [8] El Marjani A, Castro Ruiz F, Rodriguez MA, Parra Santos MT. Numerical modelling in wave energy conversion systems. Energy 2008;33(1):1246e53. [9] Takao M, Thakker A, Abdulhadi R, Setoguchi T. Effect of blade profile on the performance of large-scale Wells turbine. Int J Sustainable Energy 2006;25 (1):53e61. [10] Raghunathan S, Setoguchi T, Kaneko K. Aerodynamics of monoplane Wells turbine e a review. Int J Offshore Polar Eng 1994;4(1):68e75. [11] Brito-Melo A, Gato LMC, Sarmento AJNA. Analysis of wells turbine design parameters by numerical simulation of the OWC performance. Ocean Eng 2002;29(1):1463e77.

446

M.H. Mohamed et al. / Energy 36 (2011) 438e446

[12] Setoguchi T, Takao M, Kaneko K, Inoue M. Effect of guide vanes on the performance of a Wells turbine for wave energy conversion. Int J Offshore Polar Eng 1998;8(2):155e60. [13] Mohamed MH. Aerodynamic performance of the Wells turbine used in the sea wave conversion. Master’s Thesis, Helwan University, Cairo: Egypt; 2003. [14] Gato LMC, Falcao AFO. Performance of the Wells turbine with a double row of guide vanes. JSME Int J 1990;33(2):265e71. Ser 2. [15] Setoguchi T, Santhakumar S, Takao M, Kim TH, Kaneko K. A modified Wells turbine for wave energy conversion. Renew Energ 2003;28(1):79e91. [16] Setoguchi T, Santhakumar S, Takao M, Kim TH, Kaneko K. Effect of guide vane shape on the performance of a Wells turbine. Renew Energ 2001;23(1):1e15. [17] Webster M, Gato LMC. The effect of rotor blade sweep on the performance of the Wells turbine. Int J Offshore Polar Eng 1999;9(3):233e9. [18] Takao M, Setoguchi T, Kinoue Y, Kaneko K. Wells turbine with end plates for wave energy conversion. Ocean Eng 2007;34(11e12):1790e5. [19] Raghunathan S. Prediction of performance of biplane Wells turbine. In: Proc. of the 3d Int. Offshore and Polar Engineering Conference, Singapore; 1993, 167e171. [20] Kinoue Y, Kim TH, Setoguchi T, Mohammad M, Kaneko K, Inoue M. Hysteretic characteristics of monoplane and biplane Wells turbine for wave power conversion. Energy Convers Manag 2004;45(9e10):1617e29. [21] Setoguchi T, Takao M. Current status of self rectifying air turbines for wave energy conversion. Energy Convers Manag 2006;47(15e16):2382e96. [22] Mohamed MH, Janiga G, and Thévenin D. Performance optimization of a modified Wells turbine using non-symmetric airfoil blades. In: ASME Turbo Expo Conference, (GT2008-50815), Berlin: Germany; 2008. [23] Mohamed MH, Janiga G, Pap E and Thévenin D. Optimal shape of a modified Wells turbine considering mutual interaction between the blades. In: First Int. Conference of Energy Engineering ICEE-1, Aswan: Egypt, 2008. [24] Lee H, Kim J, Lee Y, Setoguchi T and Kang C. Numerical analysis of flow characteristics in a Wells turbine for wave power conversion. In: Proc. of the 10th Int. Offshore and Polar Engineering Conference, Seattle, USA; 2000, 4:427e433. [25] Kim T, Lee Y, Park I, Setoguchi T and Kang C. Numerical analysis for unsteady flow characteristics of the Wells turbine. In: Proceedings of the 12th International Offshore and Polar Engineering Conference, Kitakyushu, Japan; 2002, 12:694e699. [26] Mohamed MH, Janiga G, Pap E, Thévenin D. Airfoil shape optimization of a two-stage Wells turbine. Miskolc, Hungary. In: Bikfalvi P, editor. MicroCAD10 international scientific conference; 2010. p. 51e6. [27] Curran R, Gato LMC. The energy conversion performance of several types of Wells turbine designs. Proc IME J Power Energ 1997;211(2):133e45. [28] Raghunathan S. A methodology for Wells turbine design for waves energy conversion. Proc IME J Power Energ 1995;209(A3):221e32. [29] Raghunathan S, Tan CP. Effect of blade profile on the performance of Wells turbine. Int J Heat Fluid Flow 1985;6:17e22. [30] Raghunathan S, Beattie WC. Aerodynamic performance of contra-rotating Wells turbine for wave energy conversion. Proc IME J Power Energ 1996;210 (6):431e47. [31] Van den Braembussche R. Numerical optimization for advanced turbomachinery design. In: Thévenin D, Janiga G, editors. Optimization and computational fluid dynamics. Berlin, Heidelberg: Springer-Verlag; 2008. p. 147e88. [32] Liu GL. A new generation of inverse shape design problem in aerodynamics and aerothermoelasticity: concepts, theory and methods. Aircraft Eng Aerospace Tech Int J 2000;72(4):334e44. [33] Giannakoglou KC, Papadimitriou DI. Adjoint methods for shape optimization. In: Thévenin D, Janiga G, editors. Optimization and computational fluid dynamics. Berlin, Heidelberg: Springer-Verlag; 2008. p. 79e108.

[34] Mohammadi B, Pironneau O. Applied shape optimization for fluids. Oxford University Press; 2001. [35] Bonaiuti D, Zangeneh M. On the coupling of inverse design and optimization techniques for the multiobjective, multipoint design of turbomachinery blades. J Turbomachinery 2009;131(2):021014. [36] Ramamurthy R and Ghaly W. Dual point redesign of an axial compressor airfoil using a viscous inverse design method. Proc of ASME Turbo Expo 2010, GT2010e23400. Glasgow, UK; 2010. [37] Thévenin D, Janiga G. Optimization and computational fluid dynamics. Berlin, Heidelberg: Germany: Springer-Verlag; 2008. [38] Hilbert R, Janiga G, Baron R, Thévenin D. Multi-objective shape optimization of a heat exchanger using parallel genetic algorithms. Int J Heat Mass Tran 2006;49:2567e77. [39] Janiga G, Thévenin D. Reducing the CO emissions in a laminar burner using different numerical optimization methods. Proc IME J Power Energ 2007;221 (5):647e55. [40] Mohamed MH, Janiga G, Pap E, Thévenin D. Optimal performance of a Savonius turbine using frontal guiding plates. Budapest: Hungary. In: Vad J, editor. In: Proc. 14th Int. conference on modelling fluid flow (CMFF’09); 2009. p. 871e8. [41] Mohamed MH, Janiga G, Pap E, Thévenin D. Optimization of Savonius turbines using an obstacle shielding the returning blade. Renew Energ 2010;35(11): 2618e26. [42] Mohamed MH, Janiga G, Pap E, Thévenin D. Optimal blade shape of a modified Savonius turbine using an obstacle shielding the returning blade. Energy Conversion and Management 2011;52(1):236e42. [43] Setoguchi T, Takao M, Kaneko K. Hysteresis on Wells turbine characteristics in reciprocating flow. Int J Rotating Machinery 1998;1(4):17e24.

Nomenclature b: Blade span (m) c: Blade chord (m) CA: Axial force coefficient CT: Tangential force coefficient Dt: Turbine diameter (m) Dp0 : Total pressure difference (Pa) h: Hub to tip ratio FD: Drag force (N) FL: Lift force (N) FT: Tangential force (N) FA: Axial force (N) N: Rotational speed of rotor (rpm) Q: Volumetric flow-rate (m3/s) R: Resultant force (N) rt: Tip blade radius (m) rh: Hub blade radius (m) s: Blade solidity T: Output torque (N m) ut: Tip blade speed (m/s) vA: Axial air velocity (m/s) w: Relative velocity (m/s) z: Number of blades a: Incidence angle ( ) h: Aerodynamic efficieny f: Flow coefficient r: Density (kg/m3)