Abstract
A 1.5 MW horizontal-axis wind turbine blade and fluid field model are established to study the difference in the unsteady flow field and structural vibration of the wind turbine blade under one- and two-way fluid-structure interactions. The governing equations in fluid field and the motion equations in structural were developed, and the corresponding equations were discretized with the Galerkin method. Based on ANSYS CFX fluid dynamics and mechanical structural dynamics calculation software, the effects of couplings on the aerodynamic and vibration characteristics of the blade are compared and analyzed in detail. Results show that pressure distributions at different sections of the blade are concentrated near the leading edge, and the leeward side of two-way coupling is slightly higher than that of one-way coupling. Deformation along the blade span shows a nonlinear change under the coupling effect. The degree of amplitude attenuation in two-way coupling is significantly greater than that in one-way coupling because of the existence of aerodynamic damping. However, the final amplitude is still higher than the one-way coupling. The Mises stress fluctuation in the windward and leeward sides is more obvious than one-way coupling, and the discrepancy must not be ignored.
1. Introduction
As wind power technology continuously grows, the generating power of wind turbines has increased from the initial 500 W to the current ≥2 MW, and the diameter of the wind wheel from the original 3 m has currently extended to ≥100 m [1]. However, due to the long span, short chord, and good flexibility of the blades, it becomes more susceptible to aeroelastic issues caused by inertial forces, elastic forces, and aerodynamic loads during the actual operation. In some severe cases, the blade may be deformed too much and may cause fracture damage, which will bring huge losses to the wind farm [2]. Therefore, an accurate method for fluid-structure interaction analysis must be established for wind turbine blades, which offers important engineering significance for blade design and performance analysis.
According to the data transfer on fluid-structure interfaces, coupling can be categorized into one-way and two-way solutions [3]. The calculation process of one-way coupling is shown in Figure 1. It reveals that the solution of flow field and structure domain is independent, and the data transmission is carried out only when the convergence conditions are reached in each solver. However, the structural model does not feedback the deformation of the mesh after it receives aerodynamic forces generated by the fluid field, which reduces accuracy but saves considerable time and resources. Compared with one-way coupling, two-way coupling is more complex and involves additional computational time and resources. The accuracy of the solution is also improved. The solution is shown in Figure 2. In the first step of transient analysis, the convergence of the fluid field is transmitted as forces to the structure model, which uses the convergence as a boundary condition to calculate the deformation in the structure. The deformation is passed into the fluid field, causing the fluid grid to re-mesh. The process will iterate until the fluid forces and structural displacement reach the preset convergence criterion.
Researchers have done a lot of analysis in numerical simulation and experimental study to reveal the dynamic response and the mechanism of the blade under fluid-structure interaction. Lee et al. [4] developed an analytical process considering the effects of aerodynamic pressure and finite element analysis to determine the structural response of FRP wind turbine blades in fluid-structure interaction. Yu and Kwon [5] used a coupled CFD-CSD method to predict the unsteady aerodynamic load behaviors and blade deformation considering tower interference. Hsu and Bazilevs [6] and Bazilevs et al. [7] employed a nonmatching interface discretization approach with shells NURBS-based isogeometric analysis (IGA). It revealed that coupling FEM and IGA gives a good combination of efficiency, accuracy, and flexibility of the computational procedures for wind turbine coupling. Korobenko et al. [8] and Miao et al. [9] presented a simulation of two back-to-back horizontal-axis wind turbines at full scale to study the wake and its effects on downstream wind turbine. Kumar and Wurm [10] studied the propeller composite blade on the fluid-structure coupling. Dynamic mesh technology was used to analyze the structure and flow field deformation, and the results were in good agreement with the experiment. Lanzafame et al. [11] used a moving reference frame model to simulate NREL PHASE VI and compared the result with NREL’s design reference data which showed good capabilities of the modified transitional model. Rafiee et al. [12] used the modified BEM theory and CFD method to analysis horizontal-axis wind turbine, combined with ANSYS and MATLAB, to construct fluid-structure interaction in validation and blade dynamic stall characteristics analysis. Wilson et al. [13] used and shear stress transport turbulence to investigate the interactions of horizontal-axis wind turbine models in the neutrally stratified atmospheric boundary layer (ABL) and compared flow characteristics of downstream for ADM, ALM, and FRM rotor models. Leble and Barakos [14] took a 10 MW offshore floating wind turbine as the research object. The mooring system dynamics and the smooth particle hydrodynamics (SPH) coupling were used to analyze the force and moment of the mooring system and the pressure distribution on the blade. Zhang et al. [15] developed governing equations in fluid field and the motion equations in structural with geometric nonlinearity based on ALE description and analyzed the blade vibration under FSI and wind shear effect.
Although the dynamic response of the blades and the flow field characteristics of one- or two-way coupling have been studied, few works have considered the difference between the two types of coupling. Therefore, in this work, the two types of coupling for 1.5 MW horizontal-axis wind turbine blades are employed under rated conditions by using the commercial software ANSYS CFX and mechanical structural solver. The blade’s deformation for each time step is considered in one- and two-way FSI analysis, to make the solutions more reliable. The torque and blade surface pressure on different sections are evaluated by the CFD method. The stress and displacement characteristics of the blade, which can provide the reference value for the structural stability of the blades, are obtained.
2. Mathematical and Numerical Methods
2.1. Fluid Mechanics Governing Equations
Assuming fluid is uniform, irrotational, inviscid, and limited to linearly small perturbations, the fluid equation can be written as
The continuity condition of the fluid is as follows:where are the displacement component along the direction. is the density of the fluid; is a fluid dynamic pressure. is the compression modulus of the fluid. is the compression wave velocity of the fluid.
Discretized with the Galerkin method, the flow field at any point of the pressure distribution can be approximated aswhere is the shape function vector and is the pressure vector.
To bring the boundary conditions into the equation, equation (3) can be written aswhere is the normal acceleration and represent the surface area of the fluid-structure coupling interface, free surface, fixed boundary, and infinity border, respectively.
The normal acceleration on the fluid interface can also be discretized in Galerin, which is presented aswhere is the insertion function vector of the structural system, is the normal acceleration vector of the node, is the structural displacement vector, and is the coordinate transformation matrix.
Therefore, the discretized fluid equation can be expressed aswhere is the input excitation vector.
2.2. Structural Dynamics Equations
The structure equation of motion can be obtained by the finite element method, which is shown as follows:where is the displacement vector, is the mass matrix of structure, is the damping matrix of structure, is the stiffness matrix of structure, is the node vector of fluid force on the fluid-structure interface, and is the external excitation vector except .
On the fluid-structure coupling surface, the pressure distribution within the fluid element can be approximately discretized as
If there is a normal node virtual displacement vector at the interface, it can be expressed aswhere is the shape function of structure element.
The virtual work can be described as
Thus, the generalized force vector of the normal direction is presented as
After coordinate transformation, the generalized force on the overall coordinates is given as
Summarize the contribution of each fluid element:
The equation of motion of the structure under the fluid can be expressed as
Equations (6) and (14) form a fluid-structure coupling system between blades and flow field.
3. Solid Modeling and Boundary Conditions
3.1. Solid Modeling
ANSYS Design Modeler is used to build the blade solid model and fluid domain to ensure the compatibility of the model and the accuracy of the simulation results, after the airfoil coordinate is transformed. The main design parameters of the 1.5 MW wind turbine blade studied in this work are listed in Table 1, and the model shown in Figures 3(a) and 3(b) depicts the external stationary and internal rotating fluid domains, which are established by the Boolean operation. The internal rotating domain is slightly larger than the diameter of the rotor and can thus enclose the entire rotor and capture the fluid flowing around the blades. The ratio of the length, width, and height of the external stationary domain is 3 D : 7 D : 2.5 D, where D refers to the rotor diameter. The wind turbine is placed at the central axis of the stationary domain, and the distance from the inlet is 105 m. The following material properties of the blades are defined as GRP.
(a)
(b)
3.2. Computational Fluid Dynamics (CFD) Mesh
The quality of the grid determines the accuracy of the simulation results. In this work, MESHING in ANSYS Workbench is employed. The internal rotating and external stationary domains use a tetrahedral mesh, considering that the generation of tetrahedral meshes is fast and adaptable to complex shapes. Advanced options for the MESHING function are employed to encrypt the meshes around the blades. The mesh sensitivity study is presented in Table 2; three pairs of different mesh sizes at the blade surface and rotated domain are investigated, respectively, to evaluate the grid independently. To balance the computing resources and precision, the mesh size of blade surface and rotated domain are selected as 0.3 m and 1 m, which is only 0.73% relative error to the designed power of 1.5 MW. The fluid domain consists of 1,799,882 elements, whereby the structure comprises 216,802 elements. The overall meshing is shown in Figure 4.
(a)
(b)
3.3. Boundary Conditions
The inlet flow rate is set at 12. The outlet is defined as open, which allows the fluid to freely flow through the boundary. Around the stationary domain is defined as a symmetry, which can reduce the size of the outflow field and ignore the real wall effect. The ground is selected as the nonslip wall, indicating that the fluid sticks to the wall and moves with the same velocity at the wall. One- and two-way coupling solutions are defined at a transient state, and the moving grid of the rotating domain must be turned on in two-way coupling so that it can receive mesh deformation from the blades. The steady-state calculation results are set as the initial value for the transient, which can reduce calculation time. The solution time is set to 10 s, and the time step is 0.05 s. The frozen rotor method is selected to handle the rotating, static interface, which can achieve the same accuracy in an unsteady calculation but consumes less time than slid or dynamic meshes. Interface grid node matching is selected and set at GGI mode, and the turbulence model is selected as the model.
4. Results and Discussions
4.1. Model Reliability Validation
To verify the accuracy of the mathematical model and the calculation method, the Mises stress distribution in one-way coupling is carried out by using the same parameters of the 600 W wind turbine blade entity model in reference [16]. The calculated value is in good agreement with reference [16], and the maximum Mises stress relative error does not exceed 5%. So, the mathematical model and calculation method in this work can be used to predict vibration change of large size wind turbine blades. The corresponding calculation results are shown in Figure 5.
4.2. Torque and Pressure Distribution under FSI
The torque variations in one- and two-way coupling solutions under uniform wind speeds of 4, 8, and 12 m/s are shown in Figure 6. At 4 m/s, the torque variation in one-way coupling is relatively stable. While in two-way coupling, it rises first and then decreases, and the fluctuation amplitude is large. With the increase of wind speed, the torque variation in two-way coupling is generally higher than that in one-way coupling. At the rated speed, the torque fluctuations in both FSIs are similar, but those in two-way coupling are slightly higher than those of the one-way coupling.
(a)
(b)
(c)
The spanwise pressure distribution at different blade cross sections for 12 m/s wind speed is shown in Figure 7. In both FSIs, the pressure value on the leeward side is greater than the windward of the blade, which fits the fundamental characteristic of the wind turbine blade. From the leading edge to the trailing edge, the pressure at each cross-section initially increases and then decreases. When the position is closer to the trailing edge, the difference between the windward and leeward sides is smaller. The pressure changes in both FSIs are similar on the windward side and higher in two-way coupling than those in one-way coupling on the leeward side, though the discrepancy is not very evident.
(a)
(b)
(c)
4.3. Analysis of Structural Modeling
4.3.1. Blade Vibration
Fluid-structure coupling inevitably leads to blade deformation. Figure 8 depicts the distribution of the blade deformation and the vibration change along the blade span at a rated wind speed. A nonlinear change in blade span vibration can be observed. From the root to the tip, the deformation gradually increases and reaches the maximum at the tip position. According to the vibration curve of the blade span, when the position is closer to the root, the difference in amplitude is smaller between one- and two-way coupling solutions. In contrast, when close to the tip, the difference between the two couplings is magnified. At 0.35 s, the coupling is in the initial stage, and two-way coupling exhibits higher coupling effect than one-way coupling. At 2.75 s, the vibration amplitude difference reaches the maximum under one- and two-way coupling solutions within the whole coupling time. At 8.35 s, the coupling tends to be stable, and two-way coupling demonstrates larger amplitude difference than that of one-way coupling.
The three points a, b, and c are, respectively, obtained at 10, 20, and 34 m of the blade section to analyze the vibration change under the couplings, which is shown in Figure 9. Under the coupling effect, the deformation of the three points involves a decaying sine wave, and the vibration amplitude of two-way coupling is greater than that of one-way coupling at the beginning. However, due to the existence of aerodynamic and structural damping, two-way coupling exhibits a continuously attenuating vibration and a lower amplitude than one-way coupling. As time passed, two-way coupling reaches a steady value earlier than one-way coupling. With the height of the points on the blade increasing, the time to reach the stable value is also prolonged. From the trend of the vibration curve, it can be seen that the amplitude of the blade for two-way coupling is eventually higher than that for one-way coupling.
4.3.2. Stress Distribution on Windward and Leeward Sides
Figure 10 depicts the stress contour of the Mises and the distribution on the windward and leeward sides for wind turbine blades in one- and two-way coupling solutions at 0.35 s, 2.75 s, and 8.35 s, respectively.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(k)
(l)
The stress contour reveals that, in both FSIs, the stress is mainly concentrated in the roots and middle of the blade, gradually decreases from the root to the tip, and achieves the maximum at the transition areas of the root to the airfoil. However, distinguishing the discrepancy from the contour map is difficult. Thus, the stress distribution curve along the span length at different times is analyzed, as shown in Figure 10, which can represent the main differences on the blade surface.
The stress distribution curve in Figure 11 depicts that the one- and two-way coupling stress along the span length presents nonlinear variation. The windward side stress in one-way coupling exhibits a slightly greater variation than the leeward side stress during the whole computing time. However, the two-way coupling is significantly different. In the initial period of coupling time at 0.35 s, the stress distribution in the windward side at two-way coupling is close to that of one-way coupling, but the leeward side is much lower than that of one-way coupling. Thus, the vibration of the blade is higher in two-way coupling than in one-way coupling. At the mid-coupling period time of 2.75 s, the maximum discrepancy of stress between the one- and two-way coupling solutions is reached. Therefore, the vibration amplitude difference between the two coupling solutions is at the maximum at this moment. At the end of the coupling period time of 8.35 s, two-way coupling tends to be stable and its windward side stress is greater than the one-way coupling stress. But the stress in the leeward side is similar to that in the one-way coupling, which eventually makes the blade vibration in two-way coupling higher than that in one-way coupling.
(a)
(b)
(c)
5. Conclusions
In this work, the 1.5 MW horizontal-axis wind turbine blade is analyzed under one- and two-way coupling, and the aerodynamic characteristics and vibration changes of the blades are compared in detail. The conclusions are summarized as follows:(1)Under low wind speed, the torque variation is relatively stable in one-way coupling but changes significantly in two-way coupling. With the increase of wind speed, the torque amplitude of the two-way coupling is higher than that of the one-way coupling.(2)In one- and two-way coupling solutions, the pressure of each cross section of the blade initially increases and then decreases from the leading edge to the trailing edge. When the position is closer to the trailing edge of the blade, the difference between the windward and leeward sides is smaller. However, on the leeward side, the pressure distribution of two-way coupling is higher than that of one-way coupling but the difference is not very evident.(3)The deformation along the span length presents nonlinear variation, gradually increasing from the root to the tip under the couplings. Three typical positions at the blade are obtained, and the deformations present a decaying sine wave. Due to the existence of aerodynamic and structural damping, two-way coupling exhibits a continuously attenuating vibration and a lower amplitude than one-way coupling and reaches the steady state earlier.(4)Under the two couplings, the Mises stress gradually decreases from the blade root to the tip and achieves the maximum at the transition areas of the root to the airfoil. The Mises stress distribution curve between the windward and the leeward surface is obviously fluctuant, but stress distribution is close at the one-way coupling. The difference between the two couplings cannot be ignored.
Nomenclature
| : | Compression modulus of fluid |
| : | Compression wave velocity of fluid |
| : | Shape function vector |
| : | Normal acceleration |
| : | Insertion function vector of the structural system |
| : | Normal acceleration vector of the node |
| : | Structural displacement vector |
| : | Coordinate transformation matrix |
| : | Input excitation vector |
| : | Displacement vector |
| : | Mass matrix of structure |
| : | Damping matrix of structure |
| : | Stiffness matrix of structure |
| : | Node vector of fluid force on fluid-structure interface |
| : | External excitation vector except |
| : | Normal node virtual displacement vector |
| : | Shape function of structure element. |
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This work was supported by the National Natural Science Foundation of China (Grant no. 11502140) and the “Innovation Action Program” for the Local Colleges and Universities of Capacity Building Projects in 2015 (Grant no. 15110501000).