Numerical Study on Aerodynamic Performance of Different Forms of Adaptive Blades for Vertical Axis Wind Turbines

: The wind energy exploitation technique has been developed very quickly in recent years. The vertical axis wind turbine is a hot research domain due to several advantages: low noise, ﬂexible for installation, ease of maintenance, great safety and credibility, etc. The aerodynamic performances of different forms of airfoils including an active deformation airfoil and a ﬂuid-solid coupling passive airfoil with two-dimensional (2D) and three-dimensional (3D) cases have been investigated numerically in this paper. Firstly, the aerodynamic performances of the airfoils with the maximum deformation amplitudes of their cambers which are 3%, 5% and 7% of the chord length have been discussed, respectively, with the angles of attack in the range of 0 ◦ and 20 ◦ . Secondly, for the angle of attack set at 18 ◦ , the two-way ﬂuid-solid coupling simulations with the Young’s Modulus of 1 Mpa and 2 Mpa have also been investigated. Results show that: (1) for the pseudo 3D and real 3D single active deformation airfoil cases, the lift coefﬁcients increase as the maximum deformation amplitudes augment from 3% to 7% of the chord length at the same angle of attack. With the same maximum deformation amplitude, when the angles of attack increase from 0 ◦ to 20 ◦ , the lift coefﬁcients which increase ﬁrstly and then decrease are bigger than that of the original NACA0012 airfoil. When the maximum deformation amplitude of the pseudo 3D airfoil reaches 5% of the chord length, a relatively good aerodynamic performance with better inhibition effect of vortex generation can be obtained. The 3D vortex distribution demonstrates that the deformable airfoil has a better vortex generation controlling effect at the middle cross-section along the spanwise direction than the non-deformable airfoil. (2) From the aspect of ﬂuid-solid coupling, the lift increases and the drag decreases so that the lift to drag ratio has a big improvement when the Young’s Modulus is equal to 1 Mpa and 2 Mpa. The deformable airfoil can inhibit the generation and the shedding of the surface vortex when the ﬂuid-solid coupling effect is considered.


Introduction
With the rapid development of world economy, the exploitation of fossil energy resources and the corresponding problem of environmental pollution make the search of sustainable energy as an urgent issue for all countries [1]. Wind energy is one such kind of sustainable energy whose conversion technology has been developed very quickly [2][3][4]. The wind turbine is the most important equipment which converts kinetic energy of the coming flow into electrical energy. Two types of wind turbine exist in the industry: the horizontal axis wind turbine (HAWT) and the vertical axis wind turbine (VAWT). The HAWT is widely used in the commercial electrical power generation [5]. In comparison with the HAWT, the VAWT does not need any yaw mechanism and takes the wind from all the directions [6]. Among all the research aspects of VAWT, the study about the deformable blades has attracted many attention of researchers. Daynes [9] proposed a morphing flap equipped with a highly anisotropic cellular structure. This morphing flap enables undergo large deflections and high strains without a large actuation penalty. When compared to a conventional hinged flap, the proposed morphing flap is experimentally validated with a manufactured demonstrator and shown to have reduced actuation requirements. Kerho [10] proposed an adaptive airfoil which helps reduce the negative effects of dynamic stall on rotorcraft blades. In the design of large-scale wind turbine blades, Capuzzi et al. [11] has proposed a novel aero-elastic approach. The annual energy production of turbine is increased by suitably tailoring the elastic response of blade to aerodynamic pressure. Hoogedoorn et al. [12] has demonstrated a computational study about the static aero-elastic response of a two-dimensional (2D) wind turbine airfoil with varying wind conditions. The static aero-elastic effects are found to have the potential to improve the lift and the lift to drag ratio at off-design wind speed conditions. Meanwhile, the flexibility delays stall to a large pitch angle, increasing the operating range of a flexible blade airfoil. Jafaryar et al. [13] proposed a response surface methodology (RSM). This methodology is based on central composite design and can obtain an optimization design for the asymmetric blades geometry of a VAWT. Menter [14] proposed a new kind of design method for the blade of the wind energy power generation based on a topological structure about the central axis of leaf veins. The performance of this kind of bionic design for the flexible blades is studied by using the similarity of structure and working environment between wind turbine blades and plant blades. This kind of flexible blade not only widens the range of wind speed, but also improves the wind power coefficient. Mclaren [15] firstly applies the flying posture and feather structure for the design of the power generation with wind turbines. Due to the streamlined structure on the gull wing surface and the unique feather structure on the wing, there is almost no flow separation on the gull wing surface. Hence, Mclaren extracted the convex configuration data and the curved shape of the gull wing in combination with the requirements of wind turbine blades to design two types of bionic blades Among all the research aspects of VAWT, the study about the deformable blades has attracted many attention of researchers. Daynes [9] proposed a morphing flap equipped with a highly anisotropic cellular structure. This morphing flap enables undergo large deflections and high strains without a large actuation penalty. When compared to a conventional hinged flap, the proposed morphing flap is experimentally validated with a manufactured demonstrator and shown to have reduced actuation requirements. Kerho [10] proposed an adaptive airfoil which helps reduce the negative effects of dynamic stall on rotorcraft blades. In the design of large-scale wind turbine blades, Capuzzi et al. [11] has proposed a novel aero-elastic approach. The annual energy production of turbine is increased by suitably tailoring the elastic response of blade to aerodynamic pressure. Hoogedoorn et al. [12] has demonstrated a computational study about the static aero-elastic response of a two-dimensional (2D) wind turbine airfoil with varying wind conditions. The static aero-elastic effects are found to have the potential to improve the lift and the lift to drag ratio at off-design wind speed conditions. Meanwhile, the flexibility delays stall to a large pitch angle, increasing the operating range of a flexible blade airfoil. Jafaryar et al. [13] proposed a response surface methodology (RSM). This methodology is based on central composite design and can obtain an optimization design for the asymmetric blades geometry of a VAWT. Menter [14] proposed a new kind of design method for the blade of the wind energy power generation based on a topological structure about the central axis of leaf veins. The performance of this kind of bionic design for the flexible blades is studied by using the similarity of structure and working environment between wind turbine blades and plant blades. This kind of flexible blade not only widens the range of wind speed, but also improves the wind power coefficient. Mclaren [15] firstly applies the flying posture and feather structure for the design of the power generation with wind turbines. Due to the streamlined structure on the gull wing surface and the unique feather structure on the wing, there is almost no flow separation on the gull wing surface. Hence, Mclaren extracted the convex configuration data and the curved shape of the gull wing in combination with the requirements of wind turbine blades to design two types of bionic blades with convex and front curved shapes. According to the numerical simulation results, the bionic blades have better aerodynamic performance when compared with the traditional turbine blades.
In the research about the flexible blade of the VAWT, Ying et al. [16] proposed a new kind of VAWT whose blade airfoil geometry deforms automatically according to the change of surface pressure. That is to say: (1) when the surface pressure is relatively high compared with a symmetrical airfoil, the blade surface is pressed inward (such as the pressure surface of NACA2412), and the shape becomes relatively gentle; (2) while the surface pressure of the blade is relatively small compared with a symmetrical airfoil, the blade surface is convex outward (such as the suction surface of NACA2412) and the shape becomes relatively convex. Blades that undergo adaptive changes in surface shape according to the above rules will have a higher capacity in absorbing wind energy, and effectively improve the wind energy extraction efficiency for this type of VAWT. Based on Ying's research, Tong et al. [17] designed and made a blade model suitable for VAWT whose leading edge and trailing edge can have active changes simultaneously with the single input motivation. Firstly, based on the continuous deflection and deformation method of the leading and trailing edges of the airfoil with constant thickness (according to the wing rib concept and NACA0012 airfoil as an example), the blade deformation with different curvatures were achieved by applying a pre-tightening force. Secondly, the relevant wind tunnel force measurement and flow field display platform were built up, and the quasi-steady force measurement and flow field display were performed. Experiments under different deformations were carried out using the smoke line meter, seven-hole probe, and component force balance instrument. Finally, the comparison revealed that the aerodynamic performance characteristics of the blade under different deformation factors, and the results could provide reference for subsequent unsteady experiments. Tong et al. [17] have also conducted numerical simulation concerning the influence of the active deformation movement of the uniform-thickness mid-curve on the aerodynamic performance of the airfoil with NACA0012 as the reference airfoil at different angles of attack before and after the stall. Under the condition that the maximum deformation position of the airfoil was at 0.4 c, the incoming wind speed was of 9 m/s, the airfoil chord length was of 0.4 m, and the Reynolds number was of 2.5 × 10 5 . The study found that the active deformation performances of the airfoil were different at different angles of attack. The deformation amplitude and deformation frequency would bring about different changes in aerodynamic performance.
All the airfoil design mentioned as above are for the numerical simulation or experimental study of the VAWT in 2D configuration. The 3D structural effect should be also taken into account when cases are applied in 3D situations. Hence, this paper conducts the study about the aerodynamic performance of different forms of airfoils including the pseudo 3D adaptive airfoil, the real 3D adaptive airfoil and the fish skeleton airfoil for comparison. The two-way fluid-solid coupling technique is used for the research about the stress characteristic in the fish skeleton deformation process so as to give some advice for the design and development of the adaptive airfoil which are suitable for the VAWT. Section 2 presents the numerical simulation method for the cases involved in this paper, mainly about the governing equations and the fluid-solid coupling equations. These equations are the basic physical principles for the numerical simulation. Section 3 is about the presentation of different cases studied, which includes the geometry, the generation of mesh, the choice of different turbulent models and the fluid-solid coupling setting for the fluid-solid coupling simulation case. Different sub-sections correspond to different cases. Section 4 mainly focuses on the results of the numerical simulations and the corresponding discussions for different cases, where different sub-sections correspond to different cases.

Governing Equations
In this part, the governing equations for the simulation cases involved in this paper are presented. For cases which involved pseudo 3D and real 3D simulation, the form of the governing equations keeps the same except that the z direction momentum equation needs to be extended.
The numerical simulation uses the 2D non-steady turbulence model, whose governing equations can be denoted as follows: Continuous equation: Momentum equation: x direction: y direction: where Γ k and Γ ω indicate the effective diffusion coefficients k and ω, respectively. S k and S ω are the source terms defined by the users (N/m 3 ). Besides, G k and G ω indicate that the generation of the turbulent kinetic is caused by the gradient of average speed (Pa.s), which also means that the diffusion of k is caused by the turbulence: The turbulent viscosity is calculated as follows: The strain rate can be formulated as: Based on these basic equations, different solution methods for the simulation of fluid flow are used. Several different commercial software (FLUENT, CFX and Comsol, etc.) and open-source algorithm (OpenFoam) are available for the engineers to conduct the simulation, in which Fluent is adopted in this paper.

The Two-Way Fluid-Solid Coupling Theory and Equations
In this paper, a fluid-solid problem is studied. Hence the two-way fluid-solid coupling theory and Equations are described in this section. The two-way fluid-solid coupling theory and Equations take both the fluid and solid theory into consideration which helps the simulation results to approach the real physical phenomenon. The fluid-solid coupling equations are slightly different when compared with the traditional fluid dynamic equa- tions. A brief introduction from several different aspects is made: the geometric nonlinear equations, the equilibrium equations, the stress-strain constitutive relation, the boundary conditions, and the fluid-solid coupling equations.
(1) Geometric nonlinear equations The study about the strain is the key objective of solid mechanics. The formulation of strain components is highly correlated with the physical deformation of the object studied. When studying a problem of relative large deformation, there exists six basic strain components for each point in the nonlinear elastic body, denoted as following:  (11) where ε xx , ε yy , ε zz , γ xy , γ yz and γ zx are 6 basic strain components for each point in the nonlinear elastic body.
(2) Equilibrium equations Inside the elastic body, if the deformation is relatively large, then the relationship between the external force acting on it and the internal force of the elastic body will be very complicated. When the deformation amplitude is fixed, the equilibrium equations can be formulated as: where σ i (i = x, y, z) is the normal stress, τ i (i = xy, yz, xz) is the shear stress, and F i (i = x, y, z) is the body force.
(3) The stress-strain constitutive relation The stress-strain constitutive relationship can be denoted as the following formulation: where θ is the first strain invariant and E is the modulus of elasticity.

(4) Boundary conditions
Boundary conditions are of great importance in the numerical simulation of the fluidsolid coupling cases. Different boundary conditions lead to totally different simulation results. In this paper, it is assumed that there exists force at the surface of the object, then the stress boundary conditions at the interface can be denoted as follows:    f 1 = σ 11 n 1 + σ 21 n 2 + σ 31 n 3 f 2 = σ 12 n 1 + σ 22 n 2 + σ 32 n 3 f 3 = σ 13 n 1 + σ 23 n 2 + σ 33 n 3 (14) where f 1 , f 2 and f 3 are the surface force components, and n 1 , n 2 and n 3 are the cosine of the out direction normal values.

(5) The fluid-solid coupling equations
In the fluid-solid coupling calculation, it is needed to obey the corresponding conservation principles. The following equations need to be obeyed at the interface: where τ is the stress, d is the displacement, q is the heat flux, and T is the temperature. The above are the conservation equations when conducting the fluid-solid coupling analysis. These equations can be solved by setting the required parameters, the corresponding initial conditions, the boundary conditions. There exists two different kinds of solution methods: the direction method and the separation method. The direct method involves the solution of the following equations: where k is the iteration time step, A f f , ∆X k f and B f are respectively the system matrix, to be solved values and the outsider forces in the flow field. A ss , ∆X k s and B s are respectively the system matrix, to be solved values and the outsider forces in the solid field. A s f and A f s are the fluid-solid coupling matrix.
As the direct solution method involves the synchronous solution of the fluid and solid, there will not be any time lag. However, in real solution process, it is relatively hard for realizing the synchronous solution of the CFD and CSM technique. Meanwhile, the synchronous solution does not show any effect in the convergence of solving the equation. The difference between the separation method and the direct method is that the above mentioned controlling equations won't be needed. The solution for the controlling equation of fluid and solid are conducted separately. The whole solution process involves only the exchange of the solution results of the fluid field and the solid field. Meanwhile, the separation method involves less computer resources and is applied widely in the fluid-solid coupling calculations.

Validation of the Numerical Method
In this section, details about the numerical method are presented, such as the geometry, the mesh generation and the turbulence model. The validation of the numerical method chosen for the cases is also studied in this paper.

Camber Varying Model with Fixed Chord Length
In this part, how to generate the camber varying model with fixed chord length is introduced. From the self-defined functions, the coordinates of the mesh nodes of the airfoil surface can be exported. The node coordinates have been updated in real time according to the given motion law, so as to realize the airfoil deformation motion. The deformation of the airfoil can be achieved by keeping the chord length of the airfoil unchanged and changing the length of the mean camber line. The corresponding motion equation can be denoted as: In this equation, y 1 (x, t) = Acx(x − c) means the deformation amplitude of the mean camber line, where A is the maximum deformation amplitude in the thickness direction of the airfoil, c is the chord length of the airfoil, x is the coordinate in the chord direction, and g(t) = sin(2π f t) is the time-varying function. The airfoil deformation curve can be denoted as in Figure 2:

Cl/Cd
Formatted Table   Figure 2. The airfoil deformation profile.
With this model, an airfoil with the deformation obeying certain function can be obtained so as to conduct the following simulation.

Case-1: The Pseudo 3D NACA0012 Airfoil
The numerical simulation method of the adaptive blades with active deformation as depicted in Section 3.1 shares the same settings as that of the traditional fixed blades. In this paper, the NACA0012 airfoil is chosen for the following calculation as it is one of the most used in the literature.
In the calculation of the pseudo 3D NACA0012 airfoil (chord length is 0.601 m and spanwise length is 0.914 m), the k−ω−SST model with the incoming flow at Re = 1 × 10 6 is chosen. Daróczy et al. [18] systematically analyzed and compared the characteristic curves of Darrieus wind turbines calculated by adopting different computational fluid dynamics software such as Star-CCM+ and Fluent and by using different turbulence models. The 2D airfoils suitable for Darrieus wind turbines with different turbulence models are simulated and calculated. The simulation results are compared with four experimental results. It is concluded that the k-ε realizable and the k-ω-SST models are the optimum turbulence models for 2D simulation. By comparing Spalart-Allmaras model with k-ω-SST model, Gosselin et al. [19] concluded that k-ω-SST model is the most suitable model to simulate wind turbine characteristic curve. The k-ε and Spalart-Allmaras models can not capture and predict the development of the flow, especially in the laminar separation bubble [20,21]. Therefore, this paper adopts the k-ω-SST model for the corresponding research. The total mesh size is 457,600, just as shown in Figure 3. In order to validate the numerical method, different mesh sizes to calculate the lift coefficient of the airfoil at the angle of attack set as 10 • have been conducted, just as shown in Table 1. Through the analysis of Table 1, it is found that the lift coefficients for these two meshes are nearly the same, which means that the mesh size of 457,600 is enough for the following calculation. By comparing the result based on the mesh size of 457,600 and that of Ladson [22], as shown in Figure 4, it is found that these two results coincide well with each other, thus further proves that this mesh can be applied to the numerical simulation of subsequent deformed airfoil.  [18] systematically analyzed and compared the characteristic curves of Darrieus wind turbines calculated by adopting different computational fluid dynamics software such as Star-CCM+ and Fluent and by using different turbulence models. The 2D airfoils suitable for Darrieus wind turbines with different turbulence models are simulated and calculated. The simulation results are compared with four experimental results. It is concluded that the k-ε realizable and the k-ω-SST models are the optimum turbulence models for 2D simulation. By comparing Spalart-Allmaras model with k-ω-SST model, Gosselin et al. [19] concluded that k-ω-SST model is the most suitable model to simulate wind turbine characteristic curve. The k-ε and Spalart-Allmaras models can not capture and predict the development of the flow, especially in the laminar separation bubble [20,21]. Therefore, this paper adopts the k-ω-SST model for the corresponding research. The total mesh size is 457,600, just as shown in Figure 3. In order to validate the numerical method, different mesh sizes to calculate the lift coefficient of the airfoil at the angle of attack set as 10° have been conducted, just as shown in Table 1. Through the analysis of Table 1, it is found that the lift coefficients for these two meshes are nearly the same, which means that the mesh size of 457,600 is enough for the following calculation. By comparing the result based on the mesh size of 457,600 and that of Ladson [22], as shown in Figure 4, it is found that these two results coincide well with each other, thus further proves that this mesh can be applied to the numerical simulation of subsequent deformed airfoil.   Cl-exp Cl-sim-pseudo Cl-sim-2d Cd-exp Cd-sim-pseudo Cd-sim-2d

Case-2: The Real 3D NACA0012 Airfoil
The real 3D NACA0012 airfoil simulation share nearly the same numerical method as that of the pseudo 3D NACA0012 airfoil. The geometry structure and the mesh structure are showed in Figure 5. Figure 5a is the geometry structure of the airfoil, Figure 5b is the mesh at the end part of the airfoil, Figure 5c is the view of the whole grid zone, and Figure 5d is the airfoil circumferential grid distribution. As shown in Table 1, it is found that the mesh sizes of 3,140,000 and 3,500,000 give nearly the same results with a relative error of less than 1%, thus the mesh size of 3,140,000 is used for the following calculation.

Case-2: The Real 3D NACA0012 Airfoil
The real 3D NACA0012 airfoil simulation share nearly the same numerical method as that of the pseudo 3D NACA0012 airfoil. The geometry structure and the mesh structure are showed in Figure 5. Figure 5a is the geometry structure of the airfoil, Figure 5b is the mesh at the end part of the airfoil, Figure 5c is the view of the whole grid zone, and Figure 5d is the airfoil circumferential grid distribution. As shown in Table 1, it is found that the mesh sizes of 3,140,000 and 3,500,000 give nearly the same results with a relative error of less than 1%, thus the mesh size of 3,140,000 is used for the following calculation.

Case-3: The Two-way Fluid-solid Coupling Case
The two-way fluid-solid coupling case involves the NACA0012 airfoil with chord length of 1 m and spanwise length of 0.2 m. This section focuses on the study about the aerodynamic performance for the airfoil at the angle of attack 18° and the Young's modulus of 1 Mpa and 2 Mpa. In this paper, airfoil skeleton models with different shapes, as shown in Figure 6, have been studied by using the − − model, and it is found that the airfoil shape in Figure 6d which is sampled from the shape of fish seems to be more in line with the bionic principles. The boundary condition of the solid zone, the whole grid structure, the airfoil circumferential grid distribution and the grid distribution of the solid zone are depicted in Figure 7. To ensure the accuracy of calculation, the lift coefficients for different mesh sizes including 400,000, 650,000, 800,000 and 1,500,000 have

Case-3: The Two-Way Fluid-Solid Coupling Case
The two-way fluid-solid coupling case involves the NACA0012 airfoil with chord length of 1 m and spanwise length of 0.2 m. This section focuses on the study about the aerodynamic performance for the airfoil at the angle of attack 18 • and the Young's modulus of 1 Mpa and 2 Mpa. In this paper, airfoil skeleton models with different shapes, as shown in Figure 6, have been studied by using the k−ω−SST model, and it is found that the airfoil shape in Figure 6d which is sampled from the shape of fish seems to be more in line with the bionic principles. The boundary condition of the solid zone, the whole grid structure, the airfoil circumferential grid distribution and the grid distribution of the solid zone are depicted in Figure 7. To ensure the accuracy of calculation, the lift coefficients for different mesh sizes including 400,000, 650,000, 800,000 and 1,500,000 have been calculated, as shown in Table 1. Finally, the mesh with the size of 650,000 is chosen to keep a balance between saving the computational resources and keeping a relatively high level of accuracy.

Case-1: Simulation Results for the Pseudo 3D NACA0012
In this section, the lift-drag curves of the NACA0012 airfoil with the maximum deformation amplitudes of 3%, 5% and 7% of chord length at the angle of attack from 0 • to 20 • have been calculated and plotted in Figure 8. From this figure, it is found that the lift coefficients firstly increase and then decrease as the angles of attack augment from 0 • to 20 • . Compared with the original lift coefficient curve with non-deformation, all the lift coefficient curves with deformation have demonstrated better lift performance. At the same angle of attack, bigger the maximum deformation amplitude, better the lift performance.

Case-1: Simulation Results for the Pseudo 3D NACA0012
In this section, the lift-drag curves of the NACA0012 airfoil with the maximum deformation amplitudes of 3%, 5% and 7% of chord length at the angle of attack from 0° to 20° have been calculated and plotted in Figure 8. From this figure, it is found that the lift coefficients firstly increase and then decrease as the angles of attack augment from 0° to 20°. Compared with the original lift coefficient curve with non-deformation, all the lift coefficient curves with deformation have demonstrated better lift performance. At the same angle of attack, bigger the maximum deformation amplitude, better the lift performance.  Coefficient α (º) Cd-sim-pseudo

Case-2: Results for the Real 3D NACA0012 Simulation
In this section, the lift-drag curves of the real 3D NACA0012 airfoil with the maximum deformation amplitudes of 3%, 5% and 7% of chord length at the angle of attack from 0 • to 20 • have been calculated and depicted in Figure 13. From this figure, it can be found that the lift coefficient firstly increases and then decreases as the angles of attack augment from 0 • to 20 • . Compared with the original lift coefficient curve with non-deformation, all the lift coefficient curves with deformation have demonstrated better lift performance. At the same angle of attack, bigger the maximum deformation amplitude, better the lift performance.
Energies 2021, 14, x. https://doi.org/10.3390/xxxxx www.mdpi.com/journal/energies Figure 8. Comparison of the lift and drag coefficients of the NACA0012 airfoil between the maximum deformation amplitudes of 3%, 5% and 7% of chord length at the angles of attack of 0° to 20°. The flow field and vortex are demonstrated in Figures 14 and 15. In Figure 14, 0.1B, 0.5B and 0.9B represent different cross-sections along 10%, 50% and 90% of the airfoil spanwise direction, respectively. From this figure, it can be found out that the restrain of the vortex for 0.1B is not obvious, the restrain of the vortex for 0.5B is very obvious and the vortex has nearly disappeared. No obvious change can be found for the 0.9B case.

Case-2: Results for the Real 3D NACA0012 Simulation
In this section, the lift-drag curves of the real 3D NACA0012 airfoil with the maximum deformation amplitudes of 3%, 5% and 7% of chord length at the angle of attack from 0° to 20° have been calculated and depicted in Figure 13. From this figure, it can be found that the lift coefficient firstly increases and then decreases as the angles of attack augment from 0° to 20°. Compared with the original lift coefficient curve with non-deformation, all the lift coefficient curves with deformation have demonstrated better lift performance. At the same angle of attack, bigger the maximum deformation amplitude, better the lift performance.    of 0.05c. From this figure, it can be found that the surface vortex of the middle cross-section of the airfoil has been well controlled.

Case-3: Results for the Two-Way Fluid-Solid Coupling Case
In this section, the aerodynamic performance of the fish skeleton airfoil at the angle of attack of 18 • with the Young's modulus of 1 Mpa and 2 Mpa is discussed. Table 2 shows the lift, drag, lift to drag ratio of the airfoil with the Young's modulus of 1 Mpa, 2 Mpa and the original NACA0012 airfoil. From this table, it can be seen that the airfoils with the Young's modulus of 1 Mpa and 2 Mpa have a higher lift coefficient and lower drag coefficient when compared with the original NACA0012 airfoil. In Figure 16, the surface vortex of the airfoil with the Young's Modulus of 2 Mpa is greatly reduced in scale when compared with that of the original airfoil. and the original NACA0012 airfoil. From this table, it can be seen that the airfoils with the Young's modulus of 1 Mpa and 2 Mpa have a higher lift coefficient and lower drag coefficient when compared with the original NACA0012 airfoil. In Figure 16, the surface vortex of the airfoil with the Young's Modulus of 2 Mpa is greatly reduced in scale when compared with that of the original airfoil.
From Figures 16 and 17, it can be found that the pressure in the suction side of the airfoil increases and the pressure at the pressure side is reduced. The gap between the suction surface and the pressure surface is greatly increased, thus leading to the improvement of the aerodynamic performance of the airfoil. Figure 18 shows , it can be seen that the stress is mainly located in the first half part of the suction surface. While, the stress in the pressure surface is very small. By comparing Figures 21  and 22, even though the deformation amplitude of the airfoil with the Young's modulus at 1Mpa is larger than that of 2 Mpa, the corresponding stress is smaller, which is more suitable for the structure of the deformable airfoil. Original Airfoil From Figures 16 and 17, it can be found that the pressure in the suction side of the airfoil increases and the pressure at the pressure side is reduced. The gap between the suction surface and the pressure surface is greatly increased, thus leading to the improvement of the aerodynamic performance of the airfoil. Figure 18 shows the deformation shape of the airfoil in one cycle. Figures 19 and 20 show the structural stress of the airfoil with the Young's modulus at 1 Mpa and 2 Mpa. The maximum deformation amplitudes of the airfoil are 104.31 mm and 55.20 mm, which equal to 10% and 5.5% of the whole chord length respectively. Their maximum deformation posi-tions are both located in the middle part of the chord length direction. Figures 21 and 22 demonstrate equal stress of the airfoil with the Young's modulus at 1 Mpa and 2 Mpa. The maximum stress is 0.14302 Mpa and 0.14581 Mpa respectively. From the figures, it can be seen that the stress is mainly located in the first half part of the suction surface. While, the stress in the pressure surface is very small. By comparing Figures 21 and 22, even though the deformation amplitude of the airfoil with the Young's modulus at 1Mpa is larger than that of 2 Mpa, the corresponding stress is smaller, which is more suitable for the structure of the deformable airfoil.

Conclusions
This paper studies the aerodynamic performances of different forms of airfoils including the pseudo 3D adaptive airfoil, the real 3D adaptive airfoil and the fish skeleton airfoil. The following conclusions can be drawn: Firstly, with the increase of the angle of attack, for the pseudo 3D NACA0012 airfoil with the airfoil maximum deformation amplitude of 3%, 5% and 7% of the chord length, the lift coefficients firstly increase and then decrease. The lift coefficients of these deformed airfoils are all bigger than that of the original airfoil. At the same angle of attack, bigger the maximum deformation amplitude, higher the lift coefficient. The pseudo 3D NACA0012 airfoil reduces the scale of the surface vortex in some extend when compared with that of the deformed airfoil. When the maximum deformation amplitude reaches 0.05c, the best vortex controlling effect can be achieved and the airfoil demonstrates good aerodynamic performance .
Secondly, with the increase of the angle of attack, the lift coefficients of the real 3D NACA0012 airfoil at the maximum deformation amplitude of 3%, 5% and 7% of chord

Conclusions
This paper studies the aerodynamic performances of different forms of airfoils including the pseudo 3D adaptive airfoil, the real 3D adaptive airfoil and the fish skeleton airfoil. The following conclusions can be drawn: Firstly, with the increase of the angle of attack, for the pseudo 3D NACA0012 airfoil with the airfoil maximum deformation amplitude of 3%, 5% and 7% of the chord length, the lift coefficients firstly increase and then decrease. The lift coefficients of these deformed airfoils are all bigger than that of the original airfoil. At the same angle of attack, bigger the maximum deformation amplitude, higher the lift coefficient. The pseudo 3D NACA0012 airfoil reduces the scale of the surface vortex in some extend when compared with that of the deformed airfoil. When the maximum deformation amplitude reaches 0.05c, the best vortex controlling effect can be achieved and the airfoil demonstrates good aerodynamic performance.
Secondly, with the increase of the angle of attack, the lift coefficients of the real 3D NACA0012 airfoil at the maximum deformation amplitude of 3%, 5% and 7% of chord length firstly increase and then decrease. All of these airfoils have a higher lift coefficient than the original airfoil. The 3D vortex figure demonstrates that the deformable airfoil has a better vortex generation controlling effect at the middle cross-section along the spanwise direction than the non-deformable airfoil.
Thirdly, a new fish skeleton-like airfoil based on the NACA0012 airfoil has been proposed. At the angle of attack of 18 • , the aerodynamic performance of the fish skeleton like airfoil with the Young's modulus of 1 Mpa and 2 Mpa have been studied by using the fluid-solid coupling method. Results show that the passive deformation airfoil shares nearly the same characteristics as the active deformation airfoil. With a smaller Young's modulus value, a smaller stress for the airfoil surface is tended to achieve, thus leading to a better aerodynamic performance. Therefore, this kind of new fish skeleton-like airfoil is expected to improve the aerodynamic performance of the VAWT and also augment its energy extraction efficiency performance.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest.