Inverse Aerodynamic Optimization Considering Impacts of Design Tip Speed Ratio for Variable-Speed Wind Turbines

Zhiqiang Yang 1, Minghui Yin 1,*, Yan Xu 2, Yun Zou 1, Zhao Yang Dong 3,4 and Qian Zhou 5 1 School of Automation, Nanjing University of Science and Technology, Nanjing 210094, China; yang-zhiqiang@foxmail.com (Z.Y.); zouyun@vip.163.com (Y.Z.) 2 School of Electrical and Electronic Engineering, Nanyang Technological University, Singapore 639798, Singapore; eeyanxu@gmail.com 3 China Southern Power Grid Electric Power Research Institute, Guangzhou 510000, China; zydong@ieee.org 4 School of Electrical and Information Engineering, University of Sydney, Sydney, NSW 2006, Australia 5 Jiangsu Electric Power Company Research Institute, Nanjing 211103, China; xjtu@163.com * Correspondence: ymhui@vip.163.com; Tel.: +86-151-0516-8909


Introduction
Wind energy has been receiving increasing attention as one of the most exploitable renewable energy sources.The efficiency of wind energy conversion is mainly dependent on the aerodynamic shape of the wind turbine rotor [1], an essential component of a wind turbine for harvesting wind energy.Therefore, the aerodynamic optimization of wind turbine rotors plays a crucial role in the design of wind turbines.Generally, the current aerodynamic design for wind turbines is divided into two major categories: the direct method [2][3][4][5][6][7][8][9][10] and the inverse method [11][12][13][14][15][16][17].As compared with the former, the latter is distinguished by its clear principle, analytical process and fast convergence.In the inverse method, some basic design parameters (e.g., design tip speed ratio, blade number, blade radius, hub radius and airfoil distribution) and the desired aerodynamic characteristic (e.g., peak power, rated wind speed and lift coefficient distribution) are firstly specified.Then, the aerodynamic shape parameters (e.g., the blade chord and twist distributions) are directly obtained via analytical calculation.
Variable-speed wind turbines (VSWTs) have become the mainstream of large-scale wind power generation systems in recent years [18,19].They mostly operate in the variable-speed region [20], which corresponds to the wind speed range from the cut-in to the rated speed.In this region, the maximum power point tracking (MPPT) control is employed to regulate rotor speed according to wind speed variation so that the VSWT can be maintained at the design tip speed ratio (TSR) λ dgn (also called the optimal TSR λ opt ) [21].Therefore, in the conventional inverse methods for the aerodynamic optimization of VSWTs, such as the Glauert method [11], the Wilson method [12] and others [13][14][15][16][17] based on the blade element momentum (BEM) theory, the primary objective is to maximize the power coefficient at a specified λ dgn .
As one of the basic parameters in the inverse methods for the aerodynamic optimization of VSWTs, λ dgn exerts a considerable influence upon the aerodynamic shape of rotors [22,23], including the radial chord and twist distributions.To obtain an appropriate λ dgn , the functional relation between λ dgn and the maximum power coefficient C p,max , which is regarded as a type of the turbines' static aerodynamic performance, has been discussed in References [24][25][26][27].According to this relation, the λ dgn corresponding to the maximum value of C p,max is selected and then the chord and twist distributions are optimized based on this λ dgn by an inverse method in Reference [13].
Obviously, in the aforementioned inverse methods, λ dgn is determined only based on the static aerodynamic performance, and the impact of λ dgn on the dynamic process of MPPT is usually neglected.This may result in a considerable reduction in VSWTs' efficiency, as interpreted below: (1) The wind energy capture efficiency, as a closed-loop performance of VSWTs, is actually dependent on not only the static aerodynamic performance but also the MPPT dynamics.It has been pointed out in References [10,[28][29][30][31] that turbines with large inertia usually track λ dgn rather than maintaining at λ dgn .Moreover, the longer the dynamic process lasts, the lower the wind energy capture efficiency is.Especially for the VSWTs designed for low wind speed regions, larger inertia of the wind turbine rotor and weaker aerodynamic torque result in a slower dynamic behavior and a longer dynamic process of MPPT.Therefore, the dynamic process of MPPT and its effect on the wind energy production should be carefully modeled, especially for the low wind speed VSWTs.(2) Besides C p,max , λ dgn also has an effect on the dynamic process of MPPT.With the decrease (increase) of λ dgn , the MPPT process is shortened (lengthened) and correspondingly the wind energy capture efficiency is increased (decreased) [32].(3) To sum up, λ dgn exerts a comprehensive influence upon the wind energy capture efficiency through the static aerodynamic performance and the MPPT dynamics.If a larger λ dgn is conventionally applied just for improving the static aerodynamic performance (i.e., C p,max ), the turbine has to spend more time operating at the off-design TSRs due to the prolonged dynamic process of MPPT and eventually the overall efficiency of the closed-loop VSWT system is reduced.
Therefore, in the inverse aerodynamic design of VSWT rotors, the static aerodynamic performance and the dynamic process of MPPT should be compromised for determining an appropriate λ dgn .In this paper, based on the inverse design program PROPID [33,34] for the aerodynamic optimization of VSWTs, an inverse method fully considering the impact of λ dgn on the static and dynamic behavior of turbines is proposed.Compared to the existing separated design procedure, λ dgn , chord length and twist angle are jointly optimized to maximize the closed-loop performance of VSWTs (i.e., wind energy capture efficiency) in this method.Finally, the proposed method is verified by the simulation results based on the commercial software Bladed [35,36].

Problem Description
Actually, the wind energy capture efficiency, as a closed-loop performance of VSWTs, is collectively determined by the static aerodynamic performance as well as the MPPT dynamics.Although a larger design TSR can contribute to better static aerodynamic performance (i.e., higher C p,max ), it also greatly increases the duration of the dynamic process of MPPT so that the overall efficiency of the closed-loop VSWT system can be reduced.However, this comprehensive impact of λ dgn on wind energy capture efficiency through static aerodynamic performance and MPPT dynamics is usually ignored by the conventional aerodynamic design of VSWTs in the literature.

Relationship between λ dgn and C p,max
For single-airfoil wind turbine rotors, the relationship between λ dgn and C p,max characterized by a parabolic function has been derived in References [24][25][26].Furthermore, for multi-airfoil wind turbine rotors, which are generally implemented in practical VSWTs, the 1.5 MW wind partnership for advanced component technologies (WindPACT) turbine [37] developed by the national renewable energy laboratory (NREL), is chosen for analyzing the relation of C p,max with λ dgn .Specifically, keeping the original blade radius, hub radius and airfoils unchanged, the PROPID code [33,34] is applied to design the aerodynamic shape of the WindPACT turbine according to the different λ dgn , which ranges from 5.0 to 10.0 with a step size of 0.5.Then, C p,max of the redesigned wind turbine rotors corresponding to an array of λ dgn are respectively calculated by the software Bladed.
As shown in Figure 1, the variation of C p,max versus λ dgn also exhibits a parabolic shape and the maximum value of C p,max occurs at λ dgn = 8.5.Therefore, λ dgn corresponding to the maximum value of C p,max is selected in the existing inverse methods [11][12][13][14][15].Although a larger design TSR can contribute to better static aerodynamic performance (i.e., higher Cp,max), it also greatly increases the duration of the dynamic process of MPPT so that the overall efficiency of the closed-loop VSWT system can be reduced.However, this comprehensive impact of λdgn on wind energy capture efficiency through static aerodynamic performance and MPPT dynamics is usually ignored by the conventional aerodynamic design of VSWTs in the literature.

Relationship between λdgn and Cp,max
For single-airfoil wind turbine rotors, the relationship between λdgn and Cp,max characterized by a parabolic function has been derived in References [24][25][26].Furthermore, for multi-airfoil wind turbine rotors, which are generally implemented in practical VSWTs, the 1.5 MW wind partnership for advanced component technologies (WindPACT) turbine [37] developed by the national renewable energy laboratory (NREL), is chosen for analyzing the relation of Cp,max with λdgn.Specifically, keeping the original blade radius, hub radius and airfoils unchanged, the PROPID code [33,34] is applied to design the aerodynamic shape of the WindPACT turbine according to the different λdgn, which ranges from 5.0 to 10.0 with a step size of 0.5.Then, Cp,max of the redesigned wind turbine rotors corresponding to an array of λdgn are respectively calculated by the software Bladed.
As shown in Figure 1, the variation of Cp,max versus λdgn also exhibits a parabolic shape and the maximum value of Cp,max occurs at λdgn = 8.5.Therefore, λdgn corresponding to the maximum value of Cp,max is selected in the existing inverse methods [11][12][13][14][15].

The Impact of λdgn on the Maximum Power Point Tracking (MPPT) Process
Besides the aforementioned Cp,max, it has been shown from the viewpoint of a closed-loop model of VSWTs that λdgn has a significant impact on the dynamic process of MPPT [32].

A Closed-Loop Model of Variable-Speed Wind Turbines (VSWTs)
As shown in Figure 2, a closed-loop model of VSWTs consists of turbulent wind, a wind turbine rotor, a drive train, a generator and a MPPT controller.Note that because the electromagnetic response is much faster than the mechanical one, the converter can be neglected and the generator is assumed to instantly follow the torque reference requested by the MPPT controller [30].

The Impact of λ dgn on the Maximum Power Point Tracking (MPPT) Process
Besides the aforementioned C p,max , it has been shown from the viewpoint of a closed-loop model of VSWTs that λ dgn has a significant impact on the dynamic process of MPPT [32].

A Closed-Loop Model of Variable-Speed Wind Turbines (VSWTs)
As shown in Figure 2, a closed-loop model of VSWTs consists of turbulent wind, a wind turbine rotor, a drive train, a generator and a MPPT controller.Note that because the electromagnetic response is much faster than the mechanical one, the converter can be neglected and the generator is assumed to instantly follow the torque reference requested by the MPPT controller [30].In this paper, the commercial software Bladed [35,36] is employed to conduct the time-domain dynamic simulation of the closed-loop model.This software has passed the Germanischer Lloyd's certification and has been widely applied for wind turbine design, analysis and validation.
The aerodynamic model used within Bladed is based on the BEM theory.To sufficiently reproduce the dynamics of rotor wake and induced velocity flow field (i.e., the vorticity trailed into the rotor wake is influenced by the change in blade loading and the time for the change of the induced flow field is finite), the dynamic inflow model [35] is selected.
The commonly-used control strategy for MPPT, known as the optimal torque (OT) control [38], is employed.To track maximum power point, the OT control regulates the generator torque according to the measured rotor speed and a predefined torque versus rotor speed curve that is expressed as: where ω is the rotor speed, ρ is the air density, R is the blade length and λopt is equal to λdgn.

The Impact of λdgn on the MPPT Process
The optimal rotor speed opt ω for MPPT is defined as where v is the wind speed.The tracking range of rotor speed, defined as the difference of opt ω corresponding to wind speed variation, is proportionally related to λdgn.Obviously, λdgn exerts a direct influence upon the dynamic process of MPPT.
To illustrate the impact of λdgn on the MPPT process, simulations on the redesigned WindPACT turbines with two different λdgn (i.e., 5.5 and 8.5) excited by a step variation of wind speed are conducted by the Bladed software.The simulation results are compared in Figure 3 and Table 1, respectively.It can be seen that as λdgn increases, the tracking range of rotor speed is widened, which results in a longer dynamic process of MPPT.Correspondingly, the rotor spends more time operating at the off-design TSRs.In this paper, the commercial software Bladed [35,36] is employed to conduct the time-domain dynamic simulation of the closed-loop model.This software has passed the Germanischer Lloyd's certification and has been widely applied for wind turbine design, analysis and validation.
The aerodynamic model used within Bladed is based on the BEM theory.To sufficiently reproduce the dynamics of rotor wake and induced velocity flow field (i.e., the vorticity trailed into the rotor wake is influenced by the change in blade loading and the time for the change of the induced flow field is finite), the dynamic inflow model [35] is selected.
The commonly-used control strategy for MPPT, known as the optimal torque (OT) control [38], is employed.To track maximum power point, the OT control regulates the generator torque according to the measured rotor speed and a predefined torque versus rotor speed curve that is expressed as: where ω is the rotor speed, ρ is the air density, R is the blade length and λ opt is equal to λ dgn .

The Impact of λ dgn on the MPPT Process
The optimal rotor speed ω opt for MPPT is defined as where v is the wind speed.The tracking range of rotor speed, defined as the difference of ω opt corresponding to wind speed variation, is proportionally related to λ dgn .Obviously, λ dgn exerts a direct influence upon the dynamic process of MPPT.
To illustrate the impact of λ dgn on the MPPT process, simulations on the redesigned WindPACT turbines with two different λ dgn (i.e., 5.5 and 8.5) excited by a step variation of wind speed are conducted by the Bladed software.The simulation results are compared in Figure 3 and Table 1, respectively.It can be seen that as λ dgn increases, the tracking range of rotor speed is widened, which results in a longer dynamic process of MPPT.Correspondingly, the rotor spends more time operating at the off-design TSRs.Furthermore, for a given turbulent wind (mean wind speed is 6.5 m/s and the turbulence class is A [39]), the probability distributions of operational TSR corresponding to an array of λdgn (from 5.0 to 10.0) are obtained by the dynamic simulations and plotted in Figure 4.It also can be observed that the larger the λdgn, the smaller the probability of operational TSR around λdgn, namely, the more time the turbine will take to operate at the off-design TSRs.

The Comprehensive Impact of λdgn on Wind Energy Capture Efficiency via Cp,max and the MPPT Process
Because off-design TSRs correspond to the power coefficients less than Cp,max, running at such TSRs results in tracking losses [28][29][30][31].It is inferred that if a turbine operates at off-design TSRs for a long duration due to the increasing of λdgn, its MPPT efficiency is inevitably reduced.Therefore, by combining the discussion mentioned in Sections 2.1 and 2.2, λdgn actually exerts a comprehensive impact on the wind energy production of VSWTs through static aerodynamic performance and MPPT dynamics.That is to say, although a relatively large λdgn yielding the maximum value of Cp,max can contribute to improving the wind energy capture efficiency, the overall efficiency is probably reduced due to the long-duration dynamic process of MPPT.Furthermore, for a given turbulent wind (mean wind speed is 6.5 m/s and the turbulence class is A [39]), the probability distributions of operational TSR corresponding to an array of λ dgn (from 5.0 to 10.0) are obtained by the dynamic simulations and plotted in Figure 4.It also can be observed that the larger the λ dgn , the smaller the probability of operational TSR around λ dgn , namely, the more time the turbine will take to operate at the off-design TSRs.Furthermore, for a given turbulent wind (mean wind speed is 6.5 m/s and the turbulence class is A [39]), the probability distributions of operational TSR corresponding to an array of λdgn (from 5.0 to 10.0) are obtained by the dynamic simulations and plotted in Figure 4.It also can be observed that the larger the λdgn, the smaller the probability of operational TSR around λdgn, namely, the more time the turbine will take to operate at the off-design TSRs.

The Comprehensive Impact of λdgn on Wind Energy Capture Efficiency via Cp,max and the MPPT Process
Because off-design TSRs correspond to the power coefficients less than Cp,max, running at such TSRs results in tracking losses [28][29][30][31].It is inferred that if a turbine operates at off-design TSRs for a long duration due to the increasing of λdgn, its MPPT efficiency is inevitably reduced.Therefore, by combining the discussion mentioned in Sections 2.1 and 2.2, λdgn actually exerts a comprehensive impact on the wind energy production of VSWTs through static aerodynamic performance and MPPT dynamics.That is to say, although a relatively large λdgn yielding the maximum value of Cp,max can contribute to improving the wind energy capture efficiency, the overall efficiency is probably reduced due to the long-duration dynamic process of MPPT.

The Comprehensive Impact of λ dgn on Wind Energy Capture Efficiency via C p,max and the MPPT Process
Because off-design TSRs correspond to the power coefficients less than C p,max , running at such TSRs results in tracking losses [28][29][30][31].It is inferred that if a turbine operates at off-design TSRs for a long duration due to the increasing of λ dgn , its MPPT efficiency is inevitably reduced.Therefore, by combining the discussion mentioned in Sections 2.1 and 2.2, λ dgn actually exerts a comprehensive impact on the wind energy production of VSWTs through static aerodynamic performance and MPPT dynamics.That is to say, although a relatively large λ dgn yielding the maximum value of C p,max can contribute to improving the wind energy capture efficiency, the overall efficiency is probably reduced due to the long-duration dynamic process of MPPT.
To demonstrate such a comprehensive impact of λ dgn on wind energy capture efficiency, the relationship between the efficiency and λ dgn is obtained via the Bladed-based simulation on the WindPACT turbine.Note that the wind energy capture efficiency P favg [20] is in fact a closed-loop performance indicator of VSWTs operating under turbulent wind.Its expression and calculation are provided in Section 3.1.As depicted in Figure 5, P favg first increases and then decreases with λ dgn .This indicates that with the increase of λ dgn , the prolonged MPPT process rather than the raised C p,max gradually plays a dominant role in the MPPT performance.To demonstrate such a comprehensive impact of λdgn on wind energy capture efficiency, the relationship between the efficiency and λdgn is obtained via the Bladed-based simulation on the WindPACT turbine.Note that the wind energy capture efficiency Pfavg [20] is in fact a closed-loop performance indicator of VSWTs operating under turbulent wind.Its expression and calculation are provided in Section 3.1.As depicted in Figure 5, Pfavg first increases and then decreases with λdgn.This indicates that with the increase of λdgn, the prolonged MPPT process rather than the raised Cp,max gradually plays a dominant role in the MPPT performance.Hence, there is a tradeoff, which can be optimized by λdgn, between the static aerodynamic performance (i.e., increasing Cp,max) and the dynamic process of MPPT (i.e., reducing the tracking range of rotor speed).They should be simultaneously considered in determining λdgn so that the wind energy capture efficiency can be eventually improved.

Aerodynamic Optimization Considering the Comprehensive Impact of λdgn
Based on the inverse design method, this paper proposes a new method for the aerodynamic optimization of VSWTs.Different from the separated procedure applied in the conventional methods, the design TSR and the distributions of chord and twist are jointly optimized for maximizing the wind energy capture efficiency in the proposed method so that the effect of λdgn on both Cp,max and the MPPT process can be integrally considered.Essentially, the new method aims to improve the closed-loop performance of VSWTs through the coordination between the static aerodynamic performance and the MPPT dynamics.

Wind Energy Capture Efficiency
The wind energy capture efficiency is approximated as the ratio of the captured power to the available wind power Pfavg [20], i.e., where N is the total of the simulation steps, J is the rotor inertia, i v , cap i P , , ω ei are the wind speed, captured wind power, inflow wind power, rotor speed, generator torque and generator speed at the i-th step, respectively.
Before calculating Pfavg according to the Equation ( 3), the dynamic simulation on the closed-loop model of VSWTs is conducted for a given turbulent wind.The detailed procedure for calculating Pfavg is described as follows: Hence, there is a tradeoff, which can be optimized by λ dgn , between the static aerodynamic performance (i.e., increasing C p,max ) and the dynamic process of MPPT (i.e., reducing the tracking range of rotor speed).They should be simultaneously considered in determining λ dgn so that the wind energy capture efficiency can be eventually improved.

Aerodynamic Optimization Considering the Comprehensive Impact of λ dgn
Based on the inverse design method, this paper proposes a new method for the aerodynamic optimization of VSWTs.Different from the separated procedure applied in the conventional methods, the design TSR and the distributions of chord and twist are jointly optimized for maximizing the wind energy capture efficiency in the proposed method so that the effect of λ dgn on both C p,max and the MPPT process can be integrally considered.Essentially, the new method aims to improve the closed-loop performance of VSWTs through the coordination between the static aerodynamic performance and the MPPT dynamics.

Wind Energy Capture Efficiency
The wind energy capture efficiency is approximated as the ratio of the captured power to the available wind power P favg [20], i.e., where N is the total of the simulation steps, J is the rotor inertia, v i , P i cap , P i inflow , ω i , T e,i and ω e,i are the wind speed, captured wind power, inflow wind power, rotor speed, generator torque and generator speed at the i-th step, respectively.
Before calculating P favg according to the Equation ( 3), the dynamic simulation on the closed-loop model of VSWTs is conducted for a given turbulent wind.The detailed procedure for calculating P favg is described as follows: Step 1 Model the concerning wind turbine in the Bladed software, including wind turbine rotor, drive train, tower, and etc.
Step 2 Call the aerodynamic analysis module provided by the Bladed software to calculate the rotor's optimal TSR λ opt and the corresponding maximum power coefficient C p,max .
Step 3 Configure the parameters of the MPPT controller by the Equation (1).
Step 4 Perform the dynamic simulation module in the Bladed software with the given turbulent wind and obtain the simulated trajectories of the wind turbine.
Step 5 Calculate P favg according to the Equation (3) based on the simulated trajectories.

Joint Aerodynamic Optimization Based on the Inverse Design
The schematic diagrams of the conventional inverse design and the joint aerodynamic optimization proposed in this paper are illustrated in Figure 6.In the former method, as shown in Figure 6a, λ dgn is first determined according to the relationship between λ dgn and C p,max , then the aerodynamic shape is derived through analytical calculation for the specified λ dgn and other basic parameters (e.g., hub radius, rotor radius and airfoil distribution).Obviously, λ dgn and the aerodynamic shape are separately obtained and the comprehensive effect of λ dgn on the wind energy capture efficiency cannot be considered.To overcome this disadvantage, this paper proposes a new aerodynamic optimization method based on the inverse design.In this method, as represented in Figure 6b, by introducing P favg as a performance indicator, λ dgn and the aerodynamic shape are jointly optimized so that the static aerodynamic performance and the MPPT dynamics, which are dependent on λ dgn , are coordinated.
Energies 2016, 9, 1023 7 of 15 Step 1 Model the concerning wind turbine in the Bladed software, including wind turbine rotor, drive train, tower, and etc.
Step 2 Call the aerodynamic analysis module provided by the Bladed software to calculate the rotor's optimal TSR λopt and the corresponding maximum power coefficient Cp,max.
Step 3 Configure the parameters of the MPPT controller by the Equation (1).
Step 4 Perform the dynamic simulation module in the Bladed software with the given turbulent wind and obtain the simulated trajectories of the wind turbine.
Step 5 Calculate Pfavg according to the Equation (3) based on the simulated trajectories.

Joint Aerodynamic Optimization Based on the Inverse Design
The schematic diagrams of the conventional inverse design and the joint aerodynamic optimization proposed in this paper are illustrated in Figure 6.In the former method, as shown in Figure 6a, λdgn is first determined according to the relationship between λdgn and Cp,max, then the aerodynamic shape is derived through analytical calculation for the specified λdgn and other basic parameters (e.g., hub radius, rotor radius and airfoil distribution).Obviously, λdgn and the aerodynamic shape are separately obtained and the comprehensive effect of λdgn on the wind energy capture efficiency cannot be considered.To overcome this disadvantage, this paper proposes a new aerodynamic optimization method based on the inverse design.In this method, as represented in Figure 6b, by introducing Pfavg as a performance indicator, λdgn and the aerodynamic shape are jointly optimized so that the static aerodynamic performance and the MPPT dynamics, which are dependent on λdgn, are coordinated.Specifically, in the proposed method, the PROPID code [32,33] is employed for inverse aerodynamic design; Pfavg is calculated based on the Bladed-based dynamic simulation; the interval quartering algorithm [40] is adopted to search the aerodynamic design with maximum Pfavg from a number of trial λdgn.The flowchart of the joint aerodynamic optimization for VSWTs is shown in Figure 7 and the detailed procedure is described as follows: Step 1 Initialization.
Step 1.1 Initialize the aerodynamic parameters according to the baseline wind turbine, including blade number B, blade length R, rotor hub radius Rhub, airfoils and original distributions of chord Specifically, in the proposed method, the PROPID code [32,33] is employed for inverse aerodynamic design; P favg is calculated based on the Bladed-based dynamic simulation; the interval quartering algorithm [40] is adopted to search the aerodynamic design with maximum P favg from a number of trial λ dgn .The flowchart of the joint aerodynamic optimization for VSWTs is shown in Figure 7 and the detailed procedure is described as follows: Step 1 Initialization.
Step 1.1 Initialize the aerodynamic parameters according to the baseline wind turbine, including blade number B, blade length R, rotor hub radius R hub , airfoils and original distributions of chord and twist.B, R, R hub and airfoils (including the relative thickness) remain unchanged during the optimization procedure.
Step 1.2 Initialize the PROPID code with B, R, R hub , airfoils and original distributions of chord and twist determined in Step 1.1.
Step 1.3 Considering the constraints of material cost and noise, λ dgn usually ranges from 5.0 to 9.0 [22,23].Thus, set the search interval of λ dgn as [5.0, 9.0] in this paper.
Step 2 Divide the search interval equally into four connected subintervals and generate three independent boundaries of all the subintervals.
Step 3 Choose one of the boundaries as a trial λ dgn .
Step 4 Call the PROPID code to perform the inverse design for the chosen λ dgn and obtain the aerodynamic shape parameters.
Step 4.1 Specify the distributions of lift coefficient and axial induction factor.It is assumed that the lift coefficient distributes around the maximum lift-to-drag ratio of airfoils along the blade span, as listed in Appendix A, and the axial induction factor of each blade is 0.333, except the sections near the blade root [33].
Step 4.2 Run the PROPID code to derive the distributions of chord and twist corresponding to the chosen λ dgn .Thus, a newly designed wind turbine with the λ dgn chosen at Step 3 is obtained.
Step 5 Calculate P favg for the newly designed turbine according to the procedure presented in Section 3.1.
Step 6 Check whether all the three boundaries are chosen.If yes, go to Step 7; otherwise, return to Step 3.
Step 7 Check whether the length of the subinterval is less than a termination threshold (which is determined as 0.1 in this paper by comprehensively considering the optimization result and the manufacturing precision of wind turbine blades).If yes, go to Step 8; otherwise, define the two adjacent subintervals sharing the boundary corresponding to the maximum P favg as the new search interval, and then go to Step 2.
Step 8 Output the jointly optimized wind turbine rotor corresponding to the maximum P favg .
Energies 2016, 9, 1023 8 of 15 and twist.B, R, Rhub and airfoils (including the relative thickness) remain unchanged during the optimization procedure.
Step 1.2 Initialize the PROPID code with B, R, Rhub, airfoils and original distributions of chord and twist determined in Step 1.1.
Step 1.3 Considering the constraints of material cost and noise, λdgn usually ranges from 5.0 to 9.0 [22,23].Thus, set the search interval of λdgn as [5.0, 9.0] in this paper.
Step 2 Divide the search interval equally into four connected subintervals and generate three independent boundaries of all the subintervals.
Step 3 Choose one of the boundaries as a trial λdgn.
Step 4 Call the PROPID code to perform the inverse design for the chosen λdgn and obtain the aerodynamic shape parameters.
Step 4.1 Specify the distributions of lift coefficient and axial induction factor.It is assumed that the lift coefficient distributes around the maximum lift-to-drag ratio of airfoils along the blade span, as listed in Appendix A, and the axial induction factor of each blade section is 0.333, except the sections near the blade root [33].
Step 4.2 Run the PROPID code to derive the distributions of chord and twist corresponding to the chosen λdgn.Thus, a newly designed wind turbine with the λdgn chosen at Step 3 is obtained.
Step 5 Calculate Pfavg for the newly designed turbine according to the procedure presented in Section 3.1.
Step 6 Check whether all the three boundaries are chosen.If yes, go to Step 7; otherwise, return to Step 3.
Step 7 Check whether the length of the subinterval is less than a termination threshold (which is determined as 0.1 in this paper by comprehensively considering the optimization result and the manufacturing precision of wind turbine blades).If yes, go to Step 8; otherwise, define the two adjacent subintervals sharing the boundary corresponding to the maximum Pfavg as the new search interval, and then go to Step 2.
Step 8 Output the jointly optimized wind turbine rotor corresponding to the maximum Pfavg.
7. The length of the subinterval is less than 0.1?Note that performing the inverse design by PROPID code, as a step of the whole optimization procedure, just aims to obtain the blade shape corresponding to a candidate λdgn.Then, Pfavg for each candidate rotor obtained by the Bladed-based dynamic simulation can reflect the flow dynamics Note that performing the inverse design by PROPID code, as a step of the whole optimization procedure, just aims to obtain the blade shape corresponding to a candidate λ dgn .Then, P favg for each candidate rotor obtained by the Bladed-based dynamic simulation can reflect the flow dynamics predicted by the aeroelastic code.Moreover, the maximum P favg corresponds to an appropriate coordination between the static aerodynamic performance and MPPT dynamics.

Simulation Results
The inverse design method and the proposed one are respectively applied to optimize the aerodynamic shape of the baseline rotor.Comparison and analysis are conducted to verify the effectiveness of the joint optimization.

Baseline Wind Turbine and the Parameters of Turbulent Wind
The 1.5 MW WindPACT turbine [37] is chosen as the baseline and its blade geometry is listed in Table 2.The parameters for the Bladed software to generate a three-dimensional turbulent wind are summarized in Table 3.They are determined in compliance with the IEC61400-1 regulation [39].Note that in order to represent the common turbulence characteristic of a wind farm, this turbulent wind should be generated according to the specified mean wind speed and turbulence intensity, which can be statistically determined by long-term wind data of a wind farm.To highlight the comprehensive impact of λ dgn on wind energy production, the simulation case aims at low wind speed regions (such as the southeast of China [41,42]).Furthermore, by referring to the specifications of the commercial wind turbines manufactured by Envision Energy [43] (see Table 4), the mean wind speed of turbulence is set to 6.5 m/s.Therefore, the wind turbine rotor optimized by the proposed method is suitable for low wind speed VSWTs.

Determination of λ dgn for the Inverse Design Method
As mentioned above, λ dgn and the aerodynamic shape are separately obtained in the conventional inverse design method.Moreover, the primary aim of these methods is usually to improve C p,max at λ dgn [11][12][13][14][15][16] based on the implicit assumption that VSWTs can be maintained at λ dgn by the MPPT controller.Therefore, λ dgn that can contribute to the increasing of C p,max is of high priority.
To the authors' knowledge, λ dgn is commonly determined based on the empirical analysis [14,15].Furthermore, λ dgn corresponding to the maximum C p,max is selected according to the functional relation between λ dgn and C p,max , which, however, is derived based on single-airfoil wind turbine rotors [24][25][26][27].Because it is hard to deduce the analytical relationship between λ dgn and C p,max for multi-airfoil wind turbine rotors, in implementing the inverse design method in this paper, λ dgn is chosen as 8.5 according to the analysis result presented in Section 2.1.

Performance Comparison of the Optimized Wind Turbine Rotors
The optimized wind turbine rotors obtained by the inverse design method (hereinafter referred to as "separately-optimized rotor") and the proposed method (hereinafter referred to as "jointly-optimized rotor") are compared from the following aspects.

Aerodynamic Shape
The distributions of chord length and twist angle for the optimized and the original wind turbine rotors are respectively plotted in Figure 8.It can be observed from Figure 8a that compared to the original rotor, the chord length of the jointly-optimized rotor is averagely increased by 16.5%, while the chord length of the separately-optimized rotor is obviously decreased for the most part of the blade except the segments near blade root.As shown in Figure 8b, the twist angle of the jointly-optimized rotor is larger than the original rotor at the segments near the blade root and tip, while the twist angle of the separately-optimized rotor is reduced near the root and increased near the tip.It is important to note that to highlight the effect of λdgn on the closed-loop performance of VSWTs through the dynamic process of MPPT, some constraints on rotor loads and material cost are not included in this paper, which results in a relatively large adjustment of aerodynamic shape.Therefore, to improve the applicability of the proposed work, the joint aerodynamic optimization considering the above constraints is worth developing further in the future.

Power Coefficient versus TSR Curve (Static Aerodynamic Performance)
The power coefficient Cp versus TSR curves of the optimized wind turbine rotors are compared in Figure 9.Although the shapes of the two curves are very similar, the optimal TSR, namely λdgn, of the separately-optimized rotor (which equals to 8.5) is notably larger than that of the jointlyoptimized rotor (which equals to 6.4) and Cp,max of the former rotor is 0.97% greater than for the latter rotor, which is consistent with the analysis in Section 2.1.By comparing the probability distribution of the operational TSR corresponding to the two optimized rotors, as shown in Figure 10, the jointly-optimized rotor operates around the optimal TSR more frequently, which agrees with the analysis in Section 2.2.It is important to note that to highlight the effect of λ dgn on the closed-loop performance of VSWTs through the dynamic process of MPPT, some constraints on rotor loads and material cost are not included in this paper, which results in a relatively large adjustment of aerodynamic shape.Therefore, to improve the applicability of the proposed work, the joint aerodynamic optimization considering the above constraints is worth developing further in the future.

Power Coefficient versus TSR Curve (Static Aerodynamic Performance)
The power coefficient C p versus TSR curves of the optimized wind turbine rotors are compared in Figure 9.Although the shapes of the two curves are very similar, the optimal TSR, namely λ dgn , of the separately-optimized rotor (which equals to 8.5) is notably larger than that of the jointly-optimized rotor (which equals to 6.4) and C p,max of the former rotor is 0.97% greater than for the latter rotor, which is consistent with the analysis in Section 2.1.It is important to note that to highlight the effect of λdgn on the closed-loop performance of VSWTs through the dynamic process of MPPT, some constraints on rotor loads and material cost are not included in this paper, which results in a relatively large adjustment of aerodynamic shape.Therefore, to improve the applicability of the proposed work, the joint aerodynamic optimization considering the above constraints is worth developing further in the future.

Power Coefficient versus TSR Curve (Static Aerodynamic Performance)
The power coefficient Cp versus TSR curves of the optimized wind turbine rotors are compared in Figure 9.Although the shapes of the two curves are very similar, the optimal TSR, namely λdgn, of the separately-optimized rotor (which equals to 8.5) is notably larger than that of the jointlyoptimized rotor (which equals to 6.4) and Cp,max of the former rotor is 0.97% greater than for the latter rotor, which is consistent with the analysis in Section 2.1.

Probability Distribution of Operational TSR (MPPT Dynamics)
By comparing the probability distribution of the operational TSR corresponding to the two optimized rotors, as shown in Figure 10, the jointly-optimized rotor operates around the optimal TSR more frequently, which agrees with the analysis in Section 2.2.By comparing the probability distribution of the operational TSR corresponding to the two optimized rotors, as shown in Figure 10, the jointly-optimized rotor operates around the optimal TSR more frequently, which agrees with the analysis in Section 2.2.

Closed-Loop Performance
The static and closed-loop performance of the optimized wind turbine rotors are summarized in Table 5.To estimate the dynamic annual energy production (AEP) [5], the dynamic simulation on the closed-loop VSWT system is firstly conducted for a different turbulent wind with a different mean wind speed from the cut-in (3 m/s) to cut-out limit (25 m/s), and the wind power generation corresponding to a different mean wind speed is then weighted with the given Weibull probability density function of wind speed.The shape and scale parameters in the probability density function are set to 1.97 and 7.35 [42].It is revealed from Table 5 that although Cp,max of the separately-optimized rotor is 0.97% larger than that of the jointly-optimized one, the latter rotor can capture 2.76% more wind energy from the turbulence (mean wind speed is 6.5 m/s) and correspondingly increase in annual energy yield is 247.76MWh (6.21%).This indicates that shortening the tracking range and increasing the frequency of the wind turbine rotor operating around the optimal TSR by decreasing λdgn can improve wind energy production more effectively than increasing Cp,max.In other words, the MPPT dynamics play a more prominent role in the closed-loop performance than the static aerodynamic performance.Through the coordination between the static aerodynamic performance and the MPPT dynamics performed by the proposed method, higher wind energy capture efficiency of VSWTs under turbulence can be achieved.

Conclusions
This paper firstly reveals that λdgn exerts a comprehensive influence upon the overall efficiency of the closed-loop VSWT system through static aerodynamic performance and the dynamic process of MPPT.With the increase of λdgn, though Cp,max is improved, the tracking range of rotor speed as well as the dynamic process of MPPT become longer, which might eventually reduce the wind energy capture efficiency.Therefore, in the inverse aerodynamic design for VSWTs, the impact of λdgn on the static performance and the MPPT dynamic process should be compromised for determining an appropriate λdgn.
Based on the inverse design program PROPID for the aerodynamic optimization of VSWTs, this paper then proposes a new method, fully considering the impact of λdgn on the static and dynamic behavior of wind turbines.Different from the conventional separated design procedure, the design TSR, chord length and twist angle are jointly optimized to maximize the closed-loop performance of

Closed-Loop Performance
The static and closed-loop performance of the optimized wind turbine rotors are summarized in Table 5.To estimate the dynamic annual energy production (AEP) [5], the dynamic simulation on the closed-loop VSWT system is firstly conducted for a different turbulent wind with a different mean wind speed from the cut-in (3 m/s) to cut-out limit (25 m/s), and the wind power generation corresponding to a different mean wind speed is then weighted with the given Weibull probability density function of wind speed.The shape and scale parameters in the probability density function are set to 1.97 and 7.35 [42].It is revealed from Table 5 that although C p,max of the separately-optimized rotor is 0.97% larger than that of the jointly-optimized one, the latter rotor can capture 2.76% more wind energy from the turbulence (mean wind speed is 6.5 m/s) and correspondingly the increase in annual energy yield is 247.76MWh (6.21%).This indicates that shortening the tracking range and increasing the frequency of the wind turbine rotor operating around the optimal TSR by decreasing λ dgn can improve wind energy production more effectively than increasing C p,max .In other words, the MPPT dynamics play a more prominent role in the closed-loop performance than the static aerodynamic performance.Through the coordination between the static aerodynamic performance and the MPPT dynamics performed by the proposed method, higher wind energy capture efficiency of VSWTs under turbulence can be achieved.

Conclusions
This paper firstly reveals that λ dgn exerts a comprehensive influence upon the overall efficiency of the closed-loop VSWT system through static aerodynamic performance and the dynamic process of MPPT.With the increase of λ dgn , though C p,max is improved, the tracking range of rotor speed as well as the dynamic process of MPPT become longer, which might eventually reduce the wind energy capture efficiency.Therefore, in the inverse aerodynamic design for VSWTs, the impact of λ dgn on the static performance and the MPPT dynamic process should be compromised for determining an appropriate λ dgn .
Based on the inverse design program PROPID for the aerodynamic optimization of VSWTs, this paper then proposes a new method, fully considering the impact of λ dgn on the static and dynamic behavior of wind turbines.Different from the conventional separated design procedure, the design TSR, chord length and twist angle are jointly optimized to maximize the closed-loop performance of VSWTs, namely wind energy capture efficiency, in this method.Through the coordination between the static aerodynamic performance and the MPPT dynamics performed by the proposed method, higher wind energy capture efficiency of VSWTs under turbulence can be achieved, which is validated by the simulation results based on the Bladed software.

Figure 1 .
Figure 1.The variation of Cp,max versus λdgn obtained from simulations on the WindPACT turbine.

Figure 1 .
Figure 1.The variation of C p,max versus λ dgn obtained from simulations on the WindPACT turbine.

Figure 3 .
Figure 3. Trajectories of the wind turbines with different λdgn excited by a same step wind speed.

Figure 4 .
Figure 4.The probability distribution of operational tip speed ratio (TSR) corresponding to an array of λdgn.

Figure 3 .
Figure 3. Trajectories of the wind turbines with different λ dgn excited by a same step wind speed.

Figure 3 .
Figure 3. Trajectories of the wind turbines with different λdgn excited by a same step wind speed.

Figure 4 .
Figure 4.The probability distribution of operational tip speed ratio (TSR) corresponding to an array of λdgn.

Figure 4 .
Figure 4.The probability distribution of operational tip speed ratio (TSR) corresponding to an array of λ dgn .

Figure 6 .
Figure 6.The schematic diagrams of the conventional inverse design and the proposed aerodynamic design.(a) Separated optimization in the inverse design method; (b) Joint optimization in the proposed method.

Figure 6 .
Figure 6.The schematic diagrams of the conventional inverse design and the proposed aerodynamic design.(a) Separated optimization in the inverse design method; (b) Joint optimization in the proposed method.

Figure 7 .
Figure 7. Flowchart of the aerodynamic optimization of the VSWT rotor.

Figure 8 .
Figure 8. Distributions of chord length and twist angle for the original and optimized blades.(a) Chord length; (b) Twist angle.

Figure 9 .
Figure 9. Power coefficient versus TSR curves of the optimized rotors.

Figure 8 .
Figure 8. Distributions of chord length and twist angle for the original and optimized blades.(a) Chord length; (b) Twist angle.

Figure 8 .
Figure 8. Distributions of chord length and twist angle for the original and optimized blades.(a) Chord length; (b) Twist angle.

Figure 9 .
Figure 9. Power coefficient versus TSR curves of the optimized rotors.

Figure 9 .
Figure 9. Power coefficient versus TSR curves of the optimized rotors.

Figure 10 .
Figure 10.Probability distribution of operating TSR of the optimized rotors.

Figure 10 .
Figure 10.Probability distribution of operating TSR of the optimized rotors.

Table 1 .
Response of the wind turbines with different λdgn excited by a same step wind speed.

Table 1 .
Response of the wind turbines with different λ dgn excited by a same step wind speed.

Table 1 .
Response of the wind turbines with different λdgn excited by a same step wind speed.

Table 2 .
Geometry of the NREL 1.5 MW wind turbine blade.

Table 4 .
Specifications of the low wind speed wind turbines developed by Envision Energy.

Table 5 .
Comparison results of the static and closed-loop performance between different rotors.

Table 5 .
Comparison results of the static and closed-loop performance between different rotors.