Study on vibration of the offshore large-scale wind turbine blade under fluid structure interaction at different wind speeds

2.1. Physical model

Fig. 1Coupling schematic diagram between three dimensional flow field and structure

FSI problems are mainly divided into UFSI and BFSI. The mechanism of BFSI is that the fluid generates pressure on solid surface, and the pressure results in deformation, which conversely affects the fluid flow at the same time. In other words, the interaction includes both the force effect and the deformation occurrence. The mechanism of UFSI is that the fluid action force on structure at the fluid-solid interface can be obtained by the calculation of fluid domain, and then the deformation of solid is calculated by applying the fluid action force to structure. The FSI model of three dimensional flow field is shown in Fig. 1.

2.2. Control and discrete equations of fluid domain

The basic conservation laws of fluid domain include the mass conservation law, the momentum conservation law and the energy conservation law. The general form of governing equations is established to analyze the above equations conveniently and solve them by the same calculation program. Denote $\mathbf{\varphi }$ as the general variable which represents the directions of velocity components $u^{*}$, $v^{*}$, $w^{*}$, and then one can get [12]:

where $\rho$, $\mathrm{\Gamma }$ and $S$ are respectively the air density, the generalized diffusion coefficient and the generalized source term. Eq. (1) is solved for $\varphi =1$ (continuity equation), $\varphi =u_j^{*}$ (momentum equation, $u_i^{*}$ and $u_j^{*}$ stand for velocity components in each direction, and the indexes of $i$ and $j$ range from 1 to 3), and $\mathbf{\varphi }=\left\{k\text{,}\epsilon \right\}$ ($k\text{-}\epsilon$ two-equation turbulence model). $\mathrm{\Gamma }$ and $S$ are determined according to $\mathbf{\varphi }$.

After the discretization of the governing Eq. (1) is conducted by using finite volume method, integral operation is done in the control volume $P$ and at time foot $\mathrm{\Delta }t$ according to the computational grid. By further applying Gauss divergence theorem, the general form of the discrete equations can be obtained as [12]:

where $nb$ is adjacent nodes, $nb=\{w,e,s,n,b,t\}$, $a_p$ is a coefficient, and $b$ can be got from the source term.

Iterative solution of Eq. (2) is carried out to obtain the velocity field which satisfies the continuity equation, so that one can calculate the pressure distribution on the surface of blade.

2.3. Nonlinear dynamic equations

The motion equations of the blade with geometric nonlinear taken into account can be written as [13]:

where ${\mathbf{\sigma }}_{ij}$ is stress tensor, ${\mathbf{U}}_i$ is displacement vector, $\stackrel{-}{\rho }$ is material density, and ${\mathbf{f}}_i$ is body force vector. The constitutive equation and the relation between displacement and strain for the geometrical nonlinear elastic body are respectively expressed as:

where $G$ and $\lambda$ respectively stand for shear modulus and Lame coefficients, and ${\mathbf{\sigma }}_{ij}$ is unit tensor.

Generally, the discrete differential motion equation of structure domain is written as [14]:

where $\left[\mathbf{M}\right]$ and $\left[\mathbf{C}\right]$ represent mass matrix and damping matrix respectively, $\left[\mathbf{x}\right]={[u,v,w]}^T$, ${\left[\mathbf{F}\right]}_{t+\mathrm{\Delta }t}$ calculated by data transfer of dynamic mesh is column vector of the wind load related to the wind pressure distribution of the blade at time $t+\mathrm{\Delta }t$. In addition, the nonlinear term ${\left[\stackrel{-}{\mathbf{P}}\right]}_{t+\mathrm{\Delta }t}$ of the equation can be expressed as:

in which, ${\left[{\mathbf{K}}_T\right]}_t$ is the tangential stiffness matrix computed by nodal displacements at time $t$.

For Eq. (6), the displacement, the velocity and the acceleration at $t+\mathrm{\Delta }t$ corresponding to the blade dynamic response are solved by Newmark method [15] and Newton-Raphson iterative method.

3.1. Solid modeling and parameter settings

In the process of finite element analysis by using ANSYS, the establishment of solid model is the key link which directly affects the accuracy and the reliability of the numerical results. In this study, the design parameters for 1.5 MW wind turbine are listed in Table 1. The simplified solid model of the blade was set up by Pro/E, as shown in Fig. 2.

Fig. 2The solid model of 1.5 MW wind turbine blade

Currently, the glass-fiber reinforce plastic (GRP) is mainly applied to the commercial wind turbine blade, which has the characteristics of light weight, high strength and good rigidity etc. Specific material parameters are shown in Table 2 [16].

3.2. Finite element mesh generation

The basis of differential equation discretization and numerical solution is model meshing, which directly affects the accuracy of the subsequent numerical calculation. The overall mesh generation including fluid domain and structure domain can be seen in Fig. 3.

Fig. 3The meshing generation of fluid field and structure field

5.1. Model verification

To demonstrate the reliability of the theory, solid modeling and finite element simulation method, the experimental model in reference [17] is chosen to carry out the verification calculation. 3D model and numerical results are respectively shown in Fig. 5 and Fig. 6.

Fig. 5The verification model of wind blade

When the numerical program in this work is degenerated into the same blade structure as reference [17] in which the maximum Mises stress at the blade root is 42.70 MPa, the simulated result is 43.53 MPa. The relative error is only 1.94 %, which confirms the validities of our model.

Fig. 6The contour of maximum Mises stress

5.2. Grid-independent verification

With the increase of mesh density and grid unit of fluid domain and structural domain, the maximum brandish displacement listed in Table 3 is convergent. In this study, on the premise that the precision and the speed of calculation are ensured, grid unit numbers in fluid domain and structural domain are respectively set as 430,847 and 9656 in the end.

5.3. Wind pressure distribution on the blade

In the process of displacement response calculation under BFSI, it is found that the maximum displacement occurs at the time of 1.3 s. Therefore, the wind pressure distribution, displacement and stress distributions in this study all occurred at this moment.

Under condition that wind turbine stops and that inflow wind speed is 20 m/s with FSI considered, the pressure distributions in cross section of the blade middle part at 1.3 s are plotted in Fig. 7. It can be seen that pressure distribution away from the blade structure is almost not affected by the blade vibration. Moreover, the flow field near the blade is more complicated, high region and low region of wind pressure respectively appears at windward side and leeward side, while high pressure region moves downstream to low pressure region gradually. After comparing Fig. 7(a) and (b), it is found that the influence region of the pressure distribution under UFSI is much larger, and the maximum value of wind pressure under UFSI is higher than the one under BFSI.

Fig. 7Wind pressure distribution in cross section of the blade middle part at 1.3 s

a) BFSI

b) UFSI

5.4. The influence of FSI on displacement and stress

According to reference [12], the vibration of wind turbine blade under stationary state is mainly concentrated in brandish direction. Therefore, the vibration analysis in the direction of brandish is carried out in this work, and we take 20 m/s average wind speed as an example.

5.4.1. Response analysis on displacement and stress

The response curves of maximum displacement and maximum Mises stress under UFSI and BFSI are compared, as shown in Fig. 8(a) and (b). As it can be seen, the variations of displacement and stress over time identically present attenuated trend. In addition, the peak values of displacement and stress under BFSI are both less than the ones obtained from UFSI, which results from the additional damping generated by BFSI.

Fig. 8Compared curves under BFSI and UFSI

a) Displacement response

b) Mises stress response

5.4.2. Distribution of displacement and stress

Fig. 9 displays the displacement distribution of the blade at 1.3 s under UFSI and BFSI. It is easy to find that maximum displacement both occurs at the blade tip. By comparison, the peak value of maximum displacement under BFSI is lower.

Fig. 10 exhibits variation curves of maximum displacement of the blade airfoil along span direction under the same conditions as used in Fig. 9. It can be found from Fig. 10 that two curves both increase non-linearly along span direction, and the displacement difference between BFSI and UFSI is becoming greater from the root to the tip.

Fig. 9Displacement distribution of the blade at 1.3 s

a) BFSI

b) UFSI

Fig. 10Variation curves of displacement along the span direction under FSI at 1.3 s

Fig. 11Mises stress distribution of the blade at 1.3 s

a) BFSI

b) UFSI

Fig. 11 and 12 give the distributions of Mises stress and the variation curves of maximum Mises stress along the span direction under the same conditions as Fig. 9 and 10 respectively. From Fig. 11, the maximum Mises stress under UFSI is relatively bigger, while the maximum Mises stress concentrates in the middle region of the blade rather than the root, and the reason is that section modulus in bending of the blade varies from root to tip. From Fig. 12, the Mises stress along span direction almost shows logarithmic normal distribution trends, namely, the stress is diminishing respectively from the peak point of the blade middle region to both sides.

Fig. 12Variation curves of Mises stress along the span direction under FSI at 1.3 s

6.1. Comparative analysis of displacement and stress response

Fig. 13 plots the response curves of displacement and stress of the blade under FSI at inflow wind speeds 12 m/s, 16 m/s and 20 m/s. As can be seen from Fig. 13(a), the maximum displacement of the blade increases with inflow wind speed, the displacement difference between three different kinds of inflow wind speeds are 0.01 m, 0.05 m and 0.13 m respectively, and the increased amplitude which show rapidly non-linear increasing trend, which implies that the effect of BFSI on the blade displacement increases obviously with inflow wind speed. It is found from Fig. 13(b) that maximum Mises stress also increases with inflow wind speed, and the difference between three conditions are respectively 0.11 MPa, 0.47 MPa and 1.17 MPa, which indicates that BFSI also has a significant impact on the blade stress response with the increase of inflow wind speed.

Fig. 13Compared response curves under FSI at different wind speeds

a) Displacement

b) Mises stress

6.2. Distribution of displacement and stress

Fig. 14 reveals the distribution laws of maximum displacement and maximum Mises stress along the span direction of the blade under FSI at different inflow wind speeds. It can be observed from Fig. 14(a) that the maximum displacement rises non-linearly from the root to the tip, the effect of FSI on the displacement increases with inflow wind speed. Fig. 14(b) shows that Mises stress presents the non-linearly decreasing trend and the FSI effect is gradually reducing from the middle region to both sides of the blade, and that the effect of BFSI on stress response is also significant with increasing inflow wind speed.

Fig. 14Compared curves along the span direction at 1.3 s under FSI at different wind speeds

a) Maximum displacement

b) Maximum Mises stress