Influences of wind and rotating speed on the fluid-structure interaction vibration for the offshore wind turbine blade

2.1. Fluid-structure interaction physical model

In order to avoid the deficiency of the wind tunnel test which can only be based on the scaling experiment and truly reflect the effect of the structure of the wind turbine blade and the surrounding flow field, the fluid-structure interaction physical model of the wind turbine blade with the same proportion was established in this work. For the sake of the rotating operation condition of the blade, a rotating inner flow field was arranged at the periphery of the blade to simulate the rotation of the wind turbine blade. At the same time, in order to avoid the influence of the upper wind turbine on the downstream wind turbine, the length of the flow field in the external flow field is 2 times as long as that of the blade and the length of the wake region is 6 times as long as that of the blade. Based on it, the fluid-structure interaction physical model is shown in Fig. 1.

Fig. 1Fluid-structure interaction physical model

2.2. Theoretical model of fluid domain

Lagrange coordinate system is accustomed to being used in solid mechanics and only focuses on the movement of material point. In fluid mechanics, the Euler coordinate system is used and focuses only on the state of the space point. In the numerical calculation of the fluid-structure interaction, in order to solve the problem of moving interface coordination caused by the incongruous coordinates of fluid and solid domains, the ALE method is established to describe the fluid control equations. The letter $\varphi$ is the general variable, and the governing equations of the fluid domain can be expressed as [10]:

where $\rho$ is the air density; for the continuity equation, $\varphi =$ 1; for the momentum equation, $\varphi ={u^{\mathrm{*}}}_j$, ${u^{\mathrm{*}}}_i$ and ${u^{\mathrm{*}}}_j$ are speed components in each direction; for the two-equation $k$-$\epsilon$ turbulence model, $\varphi =\left\{k,\epsilon \right\}$; $\mathrm{\Gamma }$, $S$ are the generalized diffusion coefficient and the generalized source term, which are respectively determined according to $\varphi$.

The finite volume method is used to discretize the governing equations, and the integral operation is applied to the control volume $P$ and time step $\mathrm{\Delta }t$. By combining the Gauss divergence theorem, the general form of the discrete equations can be given as [10]:

here, $a$ and $b$ are coefficients; $f$ is the function to be solved; $d$ is the grid points for the value of the function; $nb$ is the adjacent node.

2.3. Calculation of wind load

Considering the fluid-structure interaction, the relation between wind speed and wind pressure from Bernoulli equation can be converted to the calculation formula of wind-load as [11]:

where, $F$ is the wind load; $C_P$ is the wind pressure distribution coefficient, and is defined as $C_P=(P_{in}^{\mathrm{*}}-P)/(P_{in}^{\mathrm{*}}-P_{out})$, where $P_{in}^{\mathrm{*}}$ represents the total pressure at the mainstream inlet, $P$ is the local static pressure, and $P_{out}$ denotes the static outlet pressure [12]; $A$ is the action area; $\dot{x}$ is the vibration velocity in the blade flapping direction.

2.4. Theoretical model in the structural domain

For the blade, the geometric nonlinear motion equation under the fluid-structure interaction is expressed as [13]:

where, ${\sigma }_{i\mathrm{}j}$ is the stress tensor; $U_i$ is the displacement vector of the blade; $\stackrel{-}{\rho }$ is the material density of the blade; $f_i$ is the vector for wind pressure distribution.

The constitutive equations and the relations between displacement and strain which the geometric nonlinear elastic body meets are respectively described as [13]:

where, $G$, ${\delta }_{ij}$ and $v_a$ is shear modulus, Lame coefficient and unit tensor respectively.

In the rotating state, the governing equations of the wind turbine blade at a time $t+∆t$ is the geometrically nonlinear discrete differential equations of motion, which can be written as:

where $\left[M\right]$ and $\left[C\right]$ are respectively the mass matrix and the damping matrix; ${\left[\dot{x}\right]}_{t+\mathrm{\Delta }t}$ and ${\left[\ddot{x}\right]}_{t+\mathrm{\Delta }t}$ are respectively the speed and the acceleration of the blade at a time $t+\mathrm{\Delta }t$, $\left[x\right]={\left[u,v,w\right]}^T$; $\left[F\left(t\right)\right]$ is the load column-vector related with the blade wind pressure distribution; $\left[K\right]$ is the stiffness matrix of the rotating blade under the fluid-structure interaction, which can be expressed as:

where $\left[K_d\right]$ is the stress stiffness matrix generated by the mutual coupling between the blade rigid motion and its elastic deformation; $\left[K_g\right]$ is the geometric stiffness matrix resulted from the in-plane deformation which is caused by the centrifugal inertia force.

For Eq. (7), the displacement, the velocity and the acceleration at $t+∆t$ corresponding to the dynamic response of the blade were calculated by Newmark and Newton-Raphson iterative methods.

2.5. Structural damping calculation

In this work, the influence of structural damping is considered when calculating the dynamic response of the blade. Rayleigh damping model is widely used to describe the orthogonal damping model of the structure damping, whose expression is as follows:

where $\alpha$ is the coefficient proportion to mass and can be ignored in actual engineering; $\beta$ is the coefficient proportional to the stiffness matrix, which can be expressed as:

Here, ${\xi }_i$ is the structural damping ratio, and its value is 0.477 % for this blade [14]; ${\omega }_i$ is the primary response frequency of the structure and the first order natural frequency of the blade ${\omega }_1=$ 0.37559 Hz is chosen.

In engineering, the damping coefficient is obtained by the method of exponential decay, which is generally based on the displacement of any point of the structure [15]. Set the displacement response of the maximum displacement point $\left(x^{\mathrm{*}},y^{\mathrm{*}},\mathrm{}{z}^{\mathrm{*}}\right)$ for the blade as $\stackrel{-}{w}\left(x^{\mathrm{*}},\mathrm{}{y}^{\mathrm{*}},\mathrm{}{z}^{\mathrm{*}},t\right)$, the dynamic equation of the equivalent dynamic system can be expressed as:

where ${\xi }^{\mathrm{*}}$ is the damping coefficient with respect to the displacement response; ${\omega }_i$ is the natural frequency of the blade; ${\stackrel{-}{w}}_k\left(x^{\mathrm{*}},\mathrm{}{y}^{\mathrm{*}},\mathrm{}{z}^{\mathrm{*}},t\right)$ and ${\stackrel{-}{w}}_{k+1}\left(x^{\mathrm{*}},\mathrm{}{y}^{\mathrm{*}},\mathrm{}{z}^{\mathrm{*}},t\right)$ represent the maximum displacements of two adjacent periods $T^{\mathrm{*}}$. Accordingly, the exponential attenuation coefficient can be expressed as:

5.1. Numerical reliability verification

In order to verify the reliability of the theory and numerical model, the wind tunnel test was carried out on the scaling model of 5 MW blade. The scaling model was introduced into the BFSI numerical calculation program in this work, and the first principal stress distribution of the blade surface along the span direction was compared with the wind tunnel test results, which was plotted in Fig. 6. It can be seen from the figure that the numerical results agree well with the experimental ones, which verifies the accuracy of the numerical methods used in this work.

5.2. Effect of rotation and fluid-structure interaction on the displacement and stress of the blade

In view of the shutdown condition ($\stackrel{-}{V}=$ 11.4 m/s and $n=$ 0 rpm) and the rated condition ($\stackrel{-}{V}=$ 11.4 m/s and $n=$ 12.1 rpm), the effects of fluid-structure interaction and rotation on the vibration characteristics of the wind turbine blade were investigated.

The comparison curves of the maximum displacement and Mises stress response under the conditions of shutdown, rated condition, unidirectional fluid-structure interaction (UFSI) and BFSI were given in Fig. 7. It can be found that the change trends of the maximum displacement and Mises stress are basically the same in four cases. Under the shutdown condition, the maximum displacements of the blade when considering BFSI and UFSI are 0.114 m and 0.139 m respectively, and the maximum Mises stresses are 1.02 MPa and 1.25 MPa. Thus, it can be concluded that the effect of BFSI on the vibration response of the blade makes the maximum values of the BFSI vibration displacement and stress response curves are slightly smaller than that of the UFSI vibration curve. Under the rotating condition, the maximum values of the BFSI vibration displacement and stress of the blade are obviously smaller than that of UFSI. It shows that BFSI has a significant influence on the vibration response of the blade under rotating condition. In addition, the attenuation velocity of BFSI vibration response curves is faster, and the occurrence time of the maximum value is later than UFSI in the way that the first peak is 0.7 s late, which is caused by the additional damping generated by the action of BFSI. Whether in UFSI or BFSI, the peak values of displacement and stress under the rated condition are much larger than the corresponding peak values under static condition, which indicates that the influence of rotation on the vibration response of the blade cannot be ignored.

Fig. 6Comparison between the calculated and experimental first principal stress distribution

Fig. 7Maximum displacement and Mises stress response curves in the spanwise direction

a) Displacement response

b) Mises stress response

The maximum value of response curves and the distribution curves of the blade section’s maximum displacement and von-Mises stress are depicted in Fig. 8 respectively in four cases. It can be drawn from the curve graphs that the maximum displacement of the blade section in the span-wise direction occurs at the tip, and the maximum stress occurs at the relative wingspan of 0.55. Under the rated condition, the difference of the maximum displacement and von-Mises stress was more significant than that under the static condition, which means the influence of BFSI on the blade’s vibration characteristic under the rated speed is greater than that in the case of a shutdown. The difference between the maximum displacements and stresses under rotating and static conditions in UFSI is bigger than the one in BFSI, which illustrates that the contribution of the rotation under UFSI is more apparent than that of the rotation under BFSI.

Fig. 8Change of maximum displacement and Mises stress of blade section along the span-wise direction

a) Displacement distribution

b) Mises stress distribution

5.3. Effect of rotating speed on the vibration characteristics of the blade

In order to ascertain the influence of the rotating speed on the dynamic response of the blade, four different rotating speeds are compared and analyzed in the 12.1 rpm range of the rated speed, with 2.6 rpm as the speed increment equal.

Under the rated inflow wind speed and different rotating speeds, the displacement and Mises stress responses of the blade are shown in Fig. 9. It is not difficult to conclude that the occurrence time of the maximum displacement and stress is 2.3 s, 2.2 s, 2.1 s and 2.0 s. That is, the greater the average wind speed is, the earlier the maximum appears. With the equal increase of rotating speed, the increasing amplitude of maximum values of the displacement peaks are 0.38 m, 0.52 m and 0.56 m respectively, and the increased amplitude of the maximum values of the stress peaks are 3.50 MPa, 4.81 MPa, 5.20 MPa respectively, both presenting the nonlinear increasing trend and indicating that the contribution of rotational speed to displacement and stress increases with the increase of the rotating speed. The displacement and stress response curves at different rotational speeds show a nonlinear decreasing trend, which is more evident with the increase of rotating speed.

Fig. 9Maximum displacement and Mises stress response curves under different rotating speed

a) Displacement response

b) Mises stress response

Fig. 10 depicts the distribution of the blade’s displacement at the time when the displacement is maximum in Fig. 9. It can be seen that the maximum displacement region of the blade occurs in the tip position. The displacement of the blade at the first half is relatively small, while the one at the last half is significant. The displacement along the span-wise direction presents the nonlinear growth. The increased amplitude shows the nonlinear increase trend along the span-wise direction under different rotating speed, indicating that the effect of speed on the displacement along the span-wise direction is more and more apparent.

Fig. 10Displacement distribution at the maximum displacement moment under different rotating speed

a) Displacement distribution contour

b) Displacement change along the span-wise direction

At the time when the stress is maximum in Fig. 9(b), the distribution of the stress of the blade is shown in Fig. 11. It is well known that the maximum stress region appears in the vicinity of the center of the blade, the stress is greater in the relative wingspan of 0.4-0.7, the maximum Mises-stress appears in the relative wingspan of 0.55. Corresponding to the stress distribution contour, the stress curve has a number of extreme value. With the equal increase of the rotating speed, the increased amplitude of the maximum stress of the blade sections mainly presents the attenuation trend from the central area to both sides. On the whole, the contribution of the rotational speed to the stress is gradually reduced from the centre to both sides of the blade.

Fig. 11Mises stress distribution at the maximum stress moment under different rotating speed

a) Stress distribution contour

b) Stress change along the span-wise direction

5.4. Effect of inflow wind speed on the vibration characteristics of the blade

Taking 4.2 m/s as the average wind speed increment, 5 kinds of different average wind speeds were compared and analyzed in the range of 19.8 m/s. Under the rated rotating speed, the displacement and Mises stress response curves of the blade to different inflow wind speed are plotted in Fig. 12. It can be found that the occurrence time of maximum values of displacement and stress peak values was 2.2 s, 2.1 s, 2.0 s, 1.9 s and 1.8 s respectively, meaning that the maximum value appears earlier as the average wind speed becomes larger. With the equal increase of the average wind speed, the maximum values of the displacement and stress peak values increase, but the increase amplitude is decreased between adjacent wind velocities. The increase amplitudes of the displacement are 0.69 m, 0.39 m, 0.26 m and 0.19 m in turn and that of the stress are 6.39 MPa, 3.59 MPa, 2.44 MPa and 1.77 MPa respectively, which both present the nonlinear decreasing trend, that is, the contribution of the average wind speed to the displacement and stress decreases with the increase of the wind speed. The displacement and stress response curves of different average wind speeds show the nonlinear attenuation trend, and their attenuation rate increases with the increase of the average wind speed.

Fig. 12Maximum displacement and Mises stress response curves under different average wind speed

a) Displacement response

b) Mises stress response

Fig. 13 depicts the displacement distribution of the blade at the time corresponding to the maximum displacement of the blade. It can be observed that the displacement contour showing the displacement from the root to the tip presents a gradual increasing trend. Under the same average wind speed, the displacement along the span-wise direction shows a non-linear increase trend, and this trend is more obvious with the increase of the wind speed. The displacement difference under different wind speed increases along the span-wise direction and reaches the maximum at the tip of the blade, suggesting that the contribution of the average wind speed on the displacement is more and more obvious along the span-wise direction. With the equal increase of the average wind speed, the increased amplitude of the maximum displacement in each section decreases unequally, which means that the influence of the wind speed on the maximum displacement of the blade decreases with the increase of average wind speed.

Fig. 13Displacement distribution at the maximum displacement moment at different average wind speeds

a) Displacement distribution contour at 15.6 m/s

b) Displacement change along the span-wise direction

Corresponding to the maximum stress moment in Fig. 13(b), the stress distribution of blade Mises is shown in Fig. 14. The stress distribution contour indicates that the maximum stress appears in the central region and at the root of the blade. The maximum stress occurred in a relative wing-span of 0.55 under different average wind speed. The difference of the maximum stress of the blade section under different average wind speed descends from the centre of the blade to both sides, showing that the effect of the average wind speed on the stress of the blade is larger in the centre of the blade, and smaller in blade’s both sides. What’s more, similar to the situation in the blade's displacement, the contribution of the average wind speed to the increase of the blade stress decreases with the increase of wind speed.

Fig. 14Mises stress distribution at the maximum stress moment at different average wind speed

a) Stress distribution contour

b) Stress change along the span-wise direction

5.5. Determination of equivalent damping coefficients of rotating blades

The maximum displacement attenuation curve of the blade under the rated condition was calculated under BFSI, as shown in Fig. 15. It is not difficult to draw the conclusion that the nonlinear attenuation trend of the blade displacement is quite obvious, and the basic exponential attenuation law is presented. Moreover, under the fluid-structure interaction effect, the blade vibration is experienced 50 s, that is, the maximum displacement response of the blade is stable after 18 cycles.

Fig. 15Attenuation curve of the blade maximum displacement under BFSI effect

Fig. 16Variation curves of the equivalent damping coefficient with the number of cycles

The displacement response data in Fig. 15 was extracted, and the equivalent damping coefficients of GRP blade were calculated by Eq. (12). The variation curve of the equivalent damping coefficients with the number of cycles was plotted in Fig. 16. As can be seen from the graph, the equivalent damping coefficient of the blade presents the nonlinear decreasing trend with the increase of the number of cycles: the decreasing velocity gradually decreases, and the equivalent damping coefficient is gradually close to zero. Based on it, if the damping coefficient is set as a constant, large errors will appear in the blade dynamics calculation under actual engineering. Only taking the damping coefficient change with the vibration cycle number into consideration, can the wind turbine blade dynamic attenuation characteristics be more truly reflected.