Multi-objective optimization of blade profiles for a horizontal axis tidal turbine

A multi-objective optimization method for blade profiles of horizontal axis tidal turbines based on ISIGHT and BEM-CFD model was proposed in this paper. Firstly, the Bezier curve was applied to get the desirable pitch angle and the chord length distribution respectively which were determined by 8 control points. The sample points were extracted according to the Latin hypercube test design, and the objective functions of each sample point were obtained through BEM-CFD model. Then the fitting relationships between the eight variables and the two objective functions were established by applying the response surface method. At last the multi-objective optimization was executed through NSGA-II, and a final optimal solution was chosen from the Pareto frontier solution obtained from the optimization. It was demonstrated that the power coefficient of the turbine was improved by 0.96%, while the trust coefficient was reduced by 2.36%. The automatic multi-objective optimization for blade profiles was realized by introducing the Bezier curve and the BEM-CFD model to the optimization platform of ISIGHT and it can be applicable to engineering practice.


Introduction
Marine currents, as a form of marine renewable energy, have attracted much attention around the world in recent years [1]. And tidal turbine is gradually becoming an economically reliable power generation device due to the cyclical and predictable characteristics of marine currents [2]. The most widely used tidal turbine is horizontal axis tidal turbine (HATT) and vertical axis tidal turbine (VATT). Compared with other kind of tidal turbines, HATTs are more efficient.
As the primary energy conversion mechanism, rotor is an important part of HATT which converts fluid kinetic energy into mechanical energy [3]. Designation of HATTs is considerably similar with wind turbines [4]. Current-state-of-the-art, there are three methods to study the hydrodynamic characteristics of HATTs, including the non-stick model based on two-dimensional airfoil uplift coefficient table such as blade element momentum (BEM) theory, the three-dimensional CFD model based on Navier-Stokes equation and the BEM-CFD model which coupled the two methods mentioned above. The BEM-CFD model can take the non-uniform flow field and free surface into account and is low in computational cost [5] [6]. Considering the computational efficiency and accuracy, the BEM-CFD model was chosen in this study.
The aim of rotor optimization is to get the best hydrodynamic performance and lowest cost of installation and maintenance at the same time. It's a multi-objective optimization. In recent years, genetic algorithm has been widely used in multi-objective optimization of wind turbines, so as to tidal turbines.  Benini E et al. [7] selected three turbines whose power are 600KW, 800KW and 1MW as the objects of this study. Annual energy production and the cost of energy were taken as two objective functions. Evolutionary algorithm (EA) was chosen as the optimization algorithm. The population size was set to 100 and the number of generations was set to 150. This paper made it possible to establish the margins between the improvement of annual energy production and the corresponding cost of energy. Abdulrahman M et al. [8] took a multi-objective optimization for wind farm layout. The MATLAB ga solver was used to obtain the optimum layout. The population size was set as large as possible to widen the random search in order to increase the probability of finding the global optimum solution instead of local optima. At last they got a three-objective optimal solution including output power, the power plant capacity coefficient and the unit output power cost considering commercial turbine model and tower height distribution. Guo et al. [9] applied the multi-objective optimization to a HATT. The quadratic polynomial relationships between the design variables and the objective functions were established according to the experiment design method of BBD and response surface technique. Then the multiobjective genetic algorithm was used to optimize the pitch angle distribution. As a result, the energy capture efficiency was improved while the axial thrust coefficient was reduced.
ISIGHT [10] was used as an optimization platform in this study. It is an excellent multi-disciplinary optimization design platform which is widely used in aviation, power, electronics and other industries. Coupled with BEM-CFD model and Bezier curve, the optimization of a HATT can be executed automatically and rapidly on ISIGHT.
In the first part of this paper, the BEM-CFD model was introduced briefly. The second part gave the multi-objective optimization model including the parametric technology and the multi-objective optimization algorithm of NSGA-Ⅱ. The results and conclusions were shown in part three and four.

Control equations
The CFD model is based on the solution of Navier-Stokes equations. For steady-state fluids, the mass and momentum conservation equations are as follows: where is the fluid density, is the i-th component of velocity vector, and represent the dynamic laminar and turbulent viscosities respectively, and represents other additional source(e.g., source due to rotor rotation).
Compared with other turbulence models, − ε model is the most widely used turbulence model because it is simple, stable and cost lowly in computation. It solves the effect of turbulence on the average flow at sub-grid scale. The main limitation of this model is using only the length scale to calculate the turbulence viscosity and assuming that the turbulence is isotropic. While in nature, the turbulence diffusion occurs over a wide range of length and time scales and the vortices may become anisotropic in the rotational flow. High-order turbulence model can predict turbulence structures better, such as the Reynolds stress model, large eddy simulation et al, but they require higher computational costs. The main object of this study was turbine rotor rather than the fluid flow structure at the high turbulence region, so − ε model was chosen. Two equations are included in this model-equation and ε equation-which represents the turbulence energy and the dissipation of this energy respectively: .
And turbulent viscosity is calculated: In equations (3)(4)(5), , , 1 , 2 and are constants, and represents the turbulence generation rate. The effective viscosity used in the diffusion equation in the momentum equation, equation and ε equation is sum of the laminar and turbulent viscosities.

Blade element momentum theory
The blade element momentum theory combines the blade element theory and the momentum theory to predict the hydrodynamic performance of rotors. In this method, the rotor plane is simplified by actuator disk theory and the influence on a turbine is time-averaged for a long period of time. The force exerted by the rotor is in same at a specific radius on a given axial plane and its magnitude is determined by the blade geometry and hydrodynamic characteristics.
The BEM theory and CFD simulation are coupled by source terms produced by turbine's rotation. The source terms can be obtained through BEM theory based on lift line theory rather than the geometric parameters of blade. Thus the characteristics of blade are described by source terms which can improve the quality of mesh in use. However, the transient flow characteristics such as tip vortices, fluid separation and turbulence formation at downstream can't be calculated because of Reynolds timeaveraged equation. Figure 1(a) shows the discretization of a rotor. Usually a rotor is discretized into 10~20 elements. For each blade element, the geometric parameters include the chord length , the thickness and radial width , as shown in figure 1 Here, and are the lift and drag coefficients respectively, and: where and ′ are the axial and tangential induction factor, 1 is the inflow velocity and Ω is the angular velocity.
The axial and tangential forces on the blade can be determined by projecting the lift and drag forces onto the axial and tangential directions as shown in figure 1(c): where is the flow inclination angle defined by: Substituting equation (6) (7) into (9) (10): These two components are source terms when converted into force per volume.

Process of BEM-CFD model
The process of the BEM-CFD method is shown in figure 2. First, the CFD model is set with the appropriate initial and boundary conditions and then the velocity and pressure field of the fluid in domain can be calculated by the completed model. In the meantime, the source terms gotten through BEM theory are introduced to the CFD model by momentum equation. These two models have a close coupling relationship in each iteration throughout the simulation.

Definitions
The tip speed ratio (TSR) is the ratio of the tip velocity to the inflow velocity: The rotor performance can be expressed by power coefficient and thrust coefficient , defined as:

Cases
The simulations presented in this paper were based on the tow-tank experiments of a three-bladed tidal turbine taken by Bahaj et al [11] [12]. The diameter of the rotor is 0.8m and the domain of CFD model is a rectangle, as shown in figure 3(a). The width, depth and the length are 2.4m(3D), 1.4m(1.75D) and 8m(10D) including 1.6m(2D) of upstream and 6.4m(8D) of downstream. The mesh was generated by ICEM, as shown in figure 3(b). There were 5.7 × 10 5 nodes, 5.5 × 10 5 hexahedral meshes and 4.3 × 10 4 tetrahedral meshes.

Parametric modeling of blades
There are five airfoils chosen in this case, NACA63-812 to NACA63-824. The parameters of the rotor were shown in table 1. The pitch angle and chord length distribution were parameterized by Bezier curve. The expression of n-order Bezier curve is: where ( ) represents one point on the curve, t∈(0,1), , is basis function of Bernstein and are control points. The fitting curve of pitch angle and chord length distribution were calculated by fourorder Bezier curve in this paper. Results were shown in figure 4. It showed that four-order Bezier curve fitted well.

Genetic algorithm
Genetic algorithm (GA) is a metaheuristic inspired by the process of natural selection which belongs to the larger class of evolutionary algorithms (EA). Genetic algorithm is commonly used to generate highquality solutions to optimization and search problems by relying on bio-inspired operators such as mutation, crossover and selection [13]. Among the genetic algorithms, the second generation of non-dominated sorting algorithm(NSGA-II) has good performance in exploration. Deb, K et al. [14] find that NSGA-II is able to maintain a better spread of solutions and converged better in the obtained non-dominated front compared with other two elitist MOEAs-PAES and SPEA. And the diversity preserving mechanism used in NSGA-II is found to be the best among the three approaches. In particular, on a problem which has strong parameter interactions, NSGA-II has been able to come closer to the true front than the other two approaches. So the NSGA-II was chosen as the multi-objective genetic optimization algorithm in this study.

Optimization model
There are two objectives demonstrating the performance of the rotor -power coefficient and thrust coefficient which have already been illustrated. The former represents the ability of a rotor to convert energy contained in flow field into mechanical energy and it is an important parameter to judge the hydrodynamic performance of a rotor. The latter represents the axial thrust of the flow applied to the rotor, which affects the structure stability. This study optimized the rotor performance to simultaneously maximize while minimizing . Pitch angle and chord length distribution were two geometry parameters of rotor which were chosen as the variables. And these two distribution curves were controlled by four control points respectively, the range of which were shown in Table 2. The process of the optimization model was shown in figure 5. Firstly, 120 sample points were extracted in the range of variables according to Latin hypercube test design. For each sample point, new pitch angle and chord length distribution were fitted by Bezier curve. And the objective functions with new blade profiles of each sample point were calculated through BEM-CFD. After that the fitting relationships between the eight variables and the two objective functions were established by applying the response surface method. The multi-objective optimization of the blade profiles was executed using the algorithm NSGA-II, and a final optimal solution was chosen from the Pareto frontier solution obtained from the optimization.

Results and discussion
In order to verify the accuracy of BEM-CFD model, the predicting performances in different operating conditions were compared with the experimental results, as shown in figure 6. It can be seen that the predictions from BEM-CFD were in good agreement with experiment.
The objective function values of the 120 sample points calculated by BEM-CFD were shown in figure 7, where the abscissa was the serial number of sample N and the ordinates were the objective function and . Then the quadratic polynomial relations between the objective functions and the design variables were established according to the 120 results by applying response surface method.  The main effects of pitch angle variables in the polynomial were shown in figure 8. As for power coefficient , it indicated that the control points 2 had a relatively greater impact while 1 was on the contrast. And as for thrust coefficient , 2 also had a relatively greater impact while 4 affect less. The final optimal population distribution was shown in figure 9 after several generations through NSGA-II and a final optimal solution was chosen. Table 3 showed the geometry and performance comparison between prototype and the optimum. It can be seen that the power coefficient was increased by 0.96%, while the axial thrust coefficient was reduced by 2.36%. Figure 10 showed the prototype and the optimal distribution of pitch angle and chord length and the three-dimensional geometry. The optimal pitch angle increased at the hub while decrease at the tip. The overall distribution span increased to 21 degree from 15 degree resulting in a higher blade twist. The chord length had a similar increase along the blade.   Optimized tidal turbine was analyzed again by both BEM-CFD and CFD methods over a series of TSRs. From figure 11 and 12 it can be seen that the power coefficients have been improved at low TSR 9 1234567890 ''"" region and the thrust coefficients have a slightly improvement at the same time. And in the high TSR region, the power and the thrust coefficients have an obviously improvement. Figure 11. Comparison of the power and thrust coefficients between the original and the optimal blades by BEM-CFD method. Figure 12. Comparison of the power and thrust coefficients between the original and the optimal blades by CFD method.
In summary, the hydrodynamic performance and the stability of the rotor were enhanced in the same time by optimizing the pitch angle and chord length distribution. This method had a value in engineering.

Conclusion
In this paper, the BEM-CFD model was used to predict the performance of a rotor, and the multiobjective optimization algorithm of NSGA-II was adopted to optimize the pitch angle and chord length. The conclusions are as follows: The Bezier curve is used to parameterize the geometry of a rotor whose precision is high in fitting. It reduces design variables obviously which can not only meet the requirements of the design optimization but also decrease the time of optimization.
The BEM-CFD model is used to predict the performance of rotor. It takes the effect of nacelle, upright, and free-surface into account and the characteristics of blades are described by source terms. The BEM-CFD model improves the accuracy of prediction and is a good foundation to automatic optimization. ISIGHT provides a platform which couples all the method mentioned above. The optimal rotor has a better hydrodynamic and stability performance. The power coefficient is increased by 0.96%, while the thrust coefficient is decreased by 2.36%. The main effect graph shows that the control points 2 located in a quarter from the hub is most influential to the performance of the rotor.