Structural Analysis of a Small H-Darrieus Wind Turbine Using Beam Models: Development and Assessment

The increasing threats of climate change and depletion of fossil fuel energy sources underscore the increasing value of renewable energy. Whilst a substantial proportion of the targets in terms of emission reduction and energy saving may be met through control in demand and large-scale energy generation, increasing interest is paid to distributed generation, especially if installed in new contexts.

Focusing on wind energy, great interest is paid to small and medium size rotors installed in unconventional contexts, e.g., the built environment [1–5].

In particular, Darrieus-type VAWTs are increasingly appreciated in small installations due to several advantages [6–9]:

• They can work effectively with the wind coming from all directions without any orientation system and they are almost noiseless [8,10].

• The power coefficient variation in turbulent and unstructured flows is reduced. Moreover, the efficiency can be even improved in a skewed flow [6,11].

• The electrical apparatus can be positioned on the ground, increasing the reliability in comparison with more conventional small horizontal-axis turbines [10].

The real use of these machines is, however, still limited and most manufacturers are small companies which often produce nearly craftmade products. The road toward a real industrialization of this technology, on the one hand, is prompting the manufacturers to develop new and more appealing design solutions, e.g., by changing the shape of the blades and/or of the supporting struts; on the other hand, efforts are being devoted to increasing the efficiency, in order to make them competitive with respect to more conventional horizontal axis wind turbines (HAWTs), and also to reducing the initial cost of the machines, primarily by means of new materials, e.g., light alloys and plastic [12].

Within this context, the structural analysis of the rotors becomes more and more important: a continuous check between the design solutions needed to maximize the aerodynamic performance of the blades and the structural constraints must be indeed provided. A continuous verification of the static behavior of the structure, therefore, becomes not only a final test but a part of the design process itself: to this purpose, a reliable and fast tool for the structural analysis can provide great advantages to an effective design.

Moreover, while Darrieus rotors are becoming an established industrial technology, proper certifications must be guaranteed. For example, the IEC 61400-2 Standard [13] has been recently suited to small VAWTs (e.g., Ref. [14]); the structural verifications indicated by the standard are indeed numerous and a fast numerical tool is again thought to represent a breakthrough for an effective development of these rotors.

A complete multiphysics analysis is thought to provide the most accurate results and a reliable estimation of the stress conditions in each part of the machine. Its computational time, however, is very high and some experimental validations are required. Even in case of 2D models (e.g., with shell elements) the time required to build and set up the model of a rotor can be extremely long, making this type of analysis often not compatible with an industrial design process, which requires several iterations on the geometry before the final features are assessed [15]. In this study, the application of a beam element approach to simulate the structural behavior of a Darrieus VAWT has been investigated. This modeling scheme is in fact supposed to perform the analysis with very short computing times; conversely, some errors could be introduced, especially in very complex geometries [16,17]. Moreover, if properly calibrated, this type of modeling can be successfully applied also to the analysis of the dynamic response of the wind turbine in real functioning, i.e., with mechanical loads and interactions with the wind flow, needed to accomplish the standards' requirements.

The turbine used for this study is the PRAMAC Revolutionair WT1KW (Fig. 1), a real industrial H-Darrieus rotor [18,19], whose main features are reported in Table 1.

The rated power of 1 kW is declared by the manufacturer at a wind speed of 15 m/s and a TSR of 2.1 corresponding to a revolution speed of 415 rpm.

From a mechanical point of view, the turbine is composed by a metal structure (central tower and six small appendixes to hold the struts) to which plastic components are connected [20]; the plastic material is a custom mix of high-resistance polypropylene and reinforcing fibers, whose characteristics have been calculated with Refs. [21] and [22]. The blades are made of two halves linked at middle span by metal plates sustained by the tie-rods, realized with metal profiles. The blades' core is composed by several stiffening ribs (Fig. 2), closed by thin plastic covers which can be assumed without structural functions (see Ref. [20]).

In order to evaluate the accuracy of the beam model developed in this study, a proper term of comparison was needed. More specifically, the assessment of the proposed modeling was obtained by a benchmark analysis with a 2D model based on shell elements, which had been already developed by the same authors for the industrial customer during the design process of the rotor [23].

The shell-elements-based model was satisfactorily validated (both in terms of vibration modes and displacements) by means of experimental tests on the turbine model of Fig. 1, whose results cannot be shown here as they are still covered by industrial secret. In particular, specific tests using an high-precision taster were carried out to evaluate the maximum displacement of the plastic blades (Fig. 2–right) at the middle section of a half-blade, i.e., where the maximum displacement was expected in the plastic: the highest error (∼20%, but low as absolute value) was measured at a rotational speed of 450 rpm (over speed test), whereas the error was constantly lower than 10% for other rotational speeds.

The first step in the FEM modeling of the rotor was a preliminary simplification of the structure with respect to the original CAD model of the real rotor.

Details such as bolts, bearings, and bushings were not modeled (as the connections were considered rigid), while other components with negligible structural tasks (but contributing to the centrifugal load) have been included as concentrated or distributed mass, without contributing to the global model stiffness. The attention was therefore focused on

• the six struts and the three central tie-rods

• the aerodynamic covers of the struts

• the six half-blades

• the central tower

Beam elements were considered the most appropriate solution for the present modeling, as they are the one-dimensional elements having 6DOF for each node. The shaft and the three blades were modeled with elements characterized by their own stiffness dependent on the material and geometric section, their own mass and moment of inertia. In details, the physical properties were evaluated every 5 mm along the vertical span of the blades. In proximity of blades' corners, the properties were instead defined using a bundle of planes having an angular spacing of 3 deg and centered in the curvature center of the ends. Moreover, since the geometry of each section and its orientation notably change along the blades span, due to the stiffening ribs pattern and helix twist (see Fig. 2), an auxiliary axis had to be defined for each blade as the fitting curve containing the positions of the centers of gravity of the sections. The properties were then defined with respect to the intersection of this new virtual axis with the 5-mm-spaced lying planes. Figures 3 and 4 report some details of the model setup for an entire blade and the whole turbine, respectively.

For each section of the blades, the beam elements had two sections at its ends in order to reduce the discontinuities on the nodes: a linear interpolation of the sections' properties was performed. Due to the stiffening ribs, the shear stress center of the sections not always laid within the section itself; as a result, the use of stress data recovery points was introduced. Moreover, although blade covers do not contribute to the total stiffening of the rotor, they actually increase the centrifugal forces, especially as they are positioned quite far from the revolution axis. These components were modeled as nonstructural masses. Finally, the plastic material was considered isotropic [21,22] and a linear behavior was assumed between stress and deformation for the range of deformations attended in the analysis.

The first validation of the proposed model was based on a static analysis under a centrifugal force field. By means of NX NASTRAN^® solver, the rotor was put into revolution with a speed of 450 RPM and the results were compared to those obtained with the shell-based model under the same conditions. To preserve industrial secret, all displacement were normalized by the maximum value experimentally measured (d_M). Figures 5 and 6 report the comparison between the two models in terms of displacements and stresses on the blades, respectively.

Notable agreement can be appreciated between the models in terms of stresses distribution. The beam model indicated a maximum stress of 56.7 MPa, whereas the shell model indicated 59.6 MPa. The maximum displacements were definitely comparable and both slightly overestimated the experimental measurements on the model. As one may notice, the beam model is deemed to be more rigid than the shell one, mainly due to the rigid modeling of the connecting elements in the 1D configuration. In order to quantitatively evaluate the accuracy of the results, a comparison was performed between the stress values (σ = E·ε) at some homologous nodes in the two models (Table 2): in details, the table reports the error on the results of the beam model, calculated as in Eq. (2).

Upon examination of Table 2, it is readily noticeable that a good matching was found between the two approaches, with a maximum error less than 7% and an average error in the order of 3%. As a second step, the analysis was repeated at the same revolution speed but with an unbalancing force of 1 kN applied on the central tower (e.g., reproducing the unbalance due to aerodynamic forces). Stress results are compared in Fig. 7.

Once again, notable agreement was found between the two models, with maximum calculated stresses of a maximum stress of 56.7 MPa for the beam model and 59.4 MPa for the shell model. The same analysis presented in Table 2 was repeated (and not reported for conciseness reasons), obtaining a maximum error of 6.9% and a mean error about 3%.

As in the first static analysis, the beam model appears more rigid. Errors in the calculated stress values were indeed attended due to the fact that they are calculated in the 1D model by means of stress data recovery, i.e., the stresses are not calculated exactly on the node but in the space around it.

As a general remark, however, the errors both on displacement and on stresses have to be considered extremely satisfactory based on the enormously different modeling.

In case of rotating structures, a correct prediction of the dynamic response is of capital importance. A modal analysis was therefore carried out on the two models.

The calculated modes are consistent between the two models, with the same predicted behavior of the turbine. For example, Figures 8 and 9 report the comparison between the first and the second vibration modes for the two modes.

Upon examination of the figures, it is worth noticing that the first mode, which is the closest to the operating range of the rotor and hence the most dangerous one, mainly excites the central tower, making it bending laterally, whereas the second mode is a flexional one involving both the blades and the tie-rods. As a general remark, as the analysis referred to the entire structure of the turbine, all the predicted modes were due to a combination of flexional and torsional components. All the modes after the second one were indeed differing from each other mainly in the direction of deformations.

With the same approach used for the static validation, based on the evidence that the same types of modes had been predicted, an estimation of the frequency error (Eq. (3)) between the two models was carried out (Table 3).

In particular, Table 3 reports the comparison between the vibrations modes within 25 Hz, which were of particular interest for the manufacturer. Some duplicated modes coming from the analysis were here not displayed. The percent error was a little bit higher than expectations only for the second mode, although the absolute error is only 1.2 Hz, which is definitely a good result if one considers the approximation introduced by the beam model. In particular, the industrial partner was mainly concerned about achieving a sufficiently accurate estimation of the vibration frequencies in order to prevent the turbine from functioning in that frequency range. In this view, the present frequency resolution was considered acceptable for this analysis.

Once the 1D model had been assessed, the attention was focused on evaluating the effects of wind load on the behavior of the rotor. This type of analysis is in fact not conventional in small Darrieus turbines, especially with helix-shaped blades.

The aerodynamic stresses were calculated by means of the VARDAR code of the University of Florence [15,24,25].

The code makes use of the blade element momentum (BEM) theory, by which the performance of the rotor are calculates coupling the momentum equation in the mainstream direction of the wind and the aerodynamic analysis of the interactions between the airfoils in motion and the oncoming flow on the rotor [10,15].

In particular, the VARDAR code has been specifically developed for H-Darrieus wind turbines using an improved version of a double multiple streamtubes approach with variable interface factors [10,26] embedded with some advanced routines to account for the main secondary and parasitic effects [27].

The capabilities of the VARDAR code have been validated with several test campaigns in the wind tunnel on 1:1 models of H-Darrieus rotors. In particular, the version of the code used in this study was preliminary verified by comparing numerical simulations with the experimental data collected on the study turbine in the wind tunnel (for further details on the experimental campaign, see Refs. [15], [24], and [27]). The comparison, presented in a dimensionless form to preserve the nondisclosure agreement with the industrial partner, highlighted notable agreement between expectations and experiments (Fig. 10), confirming that the numerical code can be reliably used to estimate the aerodynamic loads on the rotor.

The results at the nominal wind speed of 15 m/s are reported in Fig. 11.

Upon examination of the figure, one can notice that the tangential force presents the typical trend of a highly loaded turbine, with the greater energy extraction located in the upwind half of the cycle; moreover, a lack in the torque production can be noticed around ϑ = 270 deg due to the shadowing effect of the supporting tower, which has a quite high impact in such a small-diameter turbine. On the other hand, it is worth noticing that the normal force changes direction every half cycle, thus representing a high-amplitude oscillating load for the supporting tower.

With the same approach used for the static validation of the beam model, displacements and stresses on the rotor due to the predicted aerodynamics loads are presented in Fig. 12.

Since at each section the resultant of the wind force is not applied on the center of gravity, it actually introduces both a force and a torsional moment reported to the center of gravity, as effectively modeled in this approach.

A maximum displacement of approximately one tenth of the global measured one and a maximum stress of 9.23 MPa were observed. These values appear notably lower than those obtained in case of the only centrifugal forces (see Figs. 5 and 6), confirming that the stresses due to the very high rotational speeds of small rotors are by far the predominant component to be considered. However, since the aerodynamic component is a nonaxisymmetric load on the rotor, it has a direct influence on the natural frequencies of the structure and can introduce fatigue problems during the operation: further analyses will be devoted at investigating these aspects.

Finally, the study was devoted at analyzing the aeroelastic behavior of the rotor, i.e., the dynamic behavior of the structure when interacting with the wind.

The aeroelastic analysis of the rotor was carried out using the GAROS software by AeroFEM GmbH [28,29].

GAROS is a program for the aeroelastic and rotordynamic analysis of wind turbines, both with horizontal and vertical-axis. The program is based on a modal coupling method between rotor and tower [30]. The structures are idealized by finite 1D elements. The subroutine for modal analysis contains the fully coupled rotor dynamic matrices, including the geometric stiffness matrix. In the stability analysis, the quasi-steady linearized aerodynamic stiffness, damping, and mass matrices are considered.

The static deformation of the blades can be accounted for in order to calculate the stability behavior at different operating conditions. The true unsteady aerodynamic forces, including the Theodorsen function, which describes the time lag between motion and force, can be included in the simulation [30]; in a similar case, a nonlinear iteration is applied for the solution of the complex eigenvalues. In addition to the stability analysis, dynamic response analysis can also be investigated. Aerodynamic forces due to wind gradients, tower influence, oblique flows, gusts, and turbulence can be accounted for; centrifugal and gravity excitation forces are also included. For further details on the software architecture, see Ref. [29].

Even though the GAROS software has an internal program to define the 1D model of the turbine, in this study the already verified model was used. The beam model needed, however, to be slightly modified to fit the requirements of GAROS software. More specifically, two concentrated masses, connected by a rigid element, were added in correspondence with each node of the beam model (Fig. 13); this solution allows indeed the software to evaluate the orientation of the section during the aeroelastic interactions.

In order to preserve the structural characteristics of the model, however, very light masses were selected, whose almost null influence on the dynamic behavior of the rotor was verified with a preliminary analysis in NX NASTRAN^® [31].

Finally, the aerodynamic characteristics of the airfoils were introduced in the GAROS database. In particular, wind tunnel data for the custom airfoil were used [15,32], including lift, drag, and moment coefficients over the entire range of incidence angles between 0 deg and 180 deg.

The first analysis carried out with GAROS was focused on the stability characterization, either with or without accounting for the flow-structure interaction. Figure 14 reports the Campbell's diagram of the turbine with mechanical loads only.

The results of Fig. 14 were consistent with the outcomes of the preliminary analysis with NX NASTRAN^® (see Table 3), even though slight modifications can be noticed in the predicted frequencies. The little discrepancy was mainly due to some simplifications that were needed to fit the model to the GAROS requirements.

On the other hand, the modes types remained the same already described for the preliminary 1D analysis. The overall agreement was then considered fully satisfactory. Upon examination of the results, the vibration modes are independent from the rotational speed only up to approximately 160 rpm, whereas the first mode (the flexional one of the tower) has a remarkable dependence from ω. In case of mechanical loads only, no instabilities were noticed in the operating range of the turbine. The excitations introduced are only the multiples of the rotational speed since the only periodic interaction of the rotating parts is with the central tower.

The same analysis was then repeated introducing the aerodynamic contribution coming from the interaction between the wind and the structure (Fig. 15). As readily noticeable from Fig. 15, the aerodynamic component influences the dependency between the vibrational modes and the revolution speed. More specifically, the aerodynamic interactions have a detrimental effect on the structure's stability, as the natural frequencies are shifted toward the excitations of the system. On the other hand, although proper experimental verifications are needed, the absence of predicted instabilities of the system was considered very positive for the structural integrity of the rotor. The second set of tests was focused on evaluating the possibilities offered by this new software, properly set with aerodynamic coefficients and a verified beam model, in terms of fulfilling some of the IEC 61400-2 requirements.

First, the response of the rotor during a revolution speed ramp (Fig. 16) was evaluated, as shown in Fig. 17, where the calculated displacement of a node at middle-span of a blade is reported.

It is readily noticeable that the mean displacement of the blade coherently followed the revolution speed trend; moreover, the predicted maximum displacement in these conditions was consistent with the predicted value on the 2D FEM model (absolute error less than 0.01 mm).

A second test was then carried out by analyzing the response of the turbine to a variation of the wind speed while the turbine itself maintains its revolution speed (Fig. 18).

With similar conditions, the average of the blade's displacements in time was attended to be constant, as it is basically function of the rotor's revolution speed.

On the other hand, the amplitude of the displacements was attended to follow the wind speed due to the different intensity of the aerodynamic forces acting on the structure. The results of the analysis are presented in Fig. 19.

The theoretical predictions were fully accomplished and the effects of the wind speed reduction are clearly noticeable.

Finally, the GAROS software was exploited to reproduce one of the most challenging analyses indicated by the IEC 61400-2 Standard, i.e., the response of the turbine to an oncoming wind with a given turbulence in terms of intensity and power spectrum [13]; more specifically, among the different models proposed by the Standard, a von Karman turbulence was here considered, whose resulting wind trend in time is presented in Fig. 20.

The turbine's response to this oncoming wind, evaluated as the displacement at middle blade, is presented in Fig. 21 and compared to that attended in nominal conditions (Fig. 22), i.e., with constant wind. Please note that the small transitory in Fig. 22 is only due to the rising ramp of the wind that is needed as an input is the GAROS software.

A comparative analysis of the two responses clearly shows that a strongly turbulent wind, as that prescribed by the IEC 61400-2 Standard, has a notable impact on the structural behavior of a rotor like the one investigated in this study, whose light material induces notable deformations. The developed methodology is in fact thought to be able to correctly predict the turbine behavior under the majority of load cases imposed by the standard, accomplishing the industrial need of a fast and effective tool to verify the rotors.

As a final remark, it is worth noticing that the development of a reliable 1D model of the rotor allowed reducing the total degrees of freedom of the structure to 6084 from the original number of 751965 for the 2D shell model.

The relative reduction of the calculation time was impressive, going from a few days horizon to only few hours.

In this study, an extensive structural analysis of a novel Darrieus wind turbine was presented.

A 1D model of the rotor, based on beam elements was first developed and assessed in comparison to a more complex shell model and experimental data of the industrial partner.

Once the reliability of the simplified model had been verified, the effects of the wind flow interaction were evaluated.

First, the net effects of the aerodynamic loads, calculated with a BEM code specifically calibrated on the turbine, were evaluated: the analysis showed that in small wind turbines, rotating very fast, the stresses due to the aerodynamic loads are small in comparison to centrifugal contributions, although the unbalance of the resulting force can have an impact on the stability of the rotor.

The 1D model was then set to fulfill the requirements of the GAROS software, which allows an aeroelastic analysis of the Darrieus turbines. The software was tested in some study cases, including some test indicated by the IEC 61400-2 Standard, obtaining a very good agreement with theoretical expectations.

The proposed modeling ensures a notable reduction of the calculation resources in comparison to more sophisticated techniques and is then thought to represent an interesting tool for industrial manufacturers who want to effectively design and verify their models.

Thanks are due to PRAMAC Spa for the industrial cooperation. The authors would also like to acknowledge Dr. Lorenzo Follo for his contribution to the numerical analysis.

Nomenclature

AR =

aspect ratio

c =

blades' chord (m)

c_L, c_D =

lift and drag coefficients

CAD =

computer aided design

D =

turbine's diameter (m)

d_M =

maximum measured displacement (m)

f =

frequency (Hz)

F =

force per unit length (N)

f_err =

frequency error (%)

FEM =

finite elements method

H =

turbine's height (m)

I_z =

momentum of inertia (kg·m²)

N =

blades' number

P =

turbine's power (W)

R =

turbine's radius (m)

TSR =

tip-speed ratio

U =

wind speed (m/s)

W =

relative speed (m/s)

Greek Symbols

α =

incidence angle (deg)

ε =

stress error (%)

ϑ =

azimuthal position of the blade (deg)

ρ =

air density (kg/Nm³)

σ =

stress (MPa)

Σ =

turbine's solidity

χ =

blade's helix angle (deg)

ω =

revolution speed (rad/s)

Subscripts

n =

normal

REF =

reference value

t =

tangential

∞ =

undisturbed conditions

Superscript

* =

dimensionless value