On Offshore wind turbine fatigue caused by wave influenced wind

A wave influenced wind turbine simulator (WIWiTS) is developed and results from these simulations are used for fatigue analyses. WIWiTS is based on the Simulator fOr Wind Farm Application (SOWFA) developed at NREL. The simulations are transient with an unsteady Reynolds-Averaged Navier-Stokes (URANS) approach. An actuator line representation of the turbine is placed in a domain where the wave and wind are either aligned with each other or opposed to each other. Simulations with four different wave states are compared to a reference case with no waves, but the inlet wind is the same for all cases. The wave will influence the wind field, which in turn affects the damage equivalent load both at the blade root and at the tower base. In a relatively low wind regime (8 m/s in a height of 400 m) our simplified simulations show that the wave influenced wind increase the fatigue damage compared to a situation with no waves, especially for the cases where the wave opposes the wind field.


Introduction
All though it is known that fast moving waves in a low wind regime will influence the whole depth of the marine atmospheric boundary layer (MABL) [1], it is still uncertain to what extent the wave influenced wind field will affect an offshore wind turbine or a wind farm [2].The aim of this study is to investigate if a wave influenced wind field will affect the wind turbine loads and fatigue.This will be done by the use of computational fluid dynamic (CFD) and by introducing a moving wave surface in the actuator line simulations of wind turbine performance and couple this with a structural response code.In the following, this combined setup with integrated wave simulations will be referred to as Wave Influenced Wind Turbine Simulations (WIWiTS).
In this study, we only look at the possible effect the waves will have on the wind.In the MABL there is of cause a close interaction between the wind field and the sea surface where heat and momentum is exchanged.A real representation of the varying wind and wave field is a huge task and very computational requiring.The problem thus had to be simplified.First of all, we only look at how the waves influence the wind field and not the other way around.Hence, the waves are seen as a solid moving wave surface.A neutral atmosphere is anticipated, so no buoyancy effects or heat exchange a Corresponding author: siri.kalvig@stormgeo.comare considered.In addition, the Coriolis force is neglected, and since we are studying processes on a relatively small scale in a short period, this is believed to be a valid assumption.The wave field can be divided into a wind-generated wave field (characterized by short periods and relatively low phase speed) and waves that have propagated away from the source origin or swell (characterized by longer periods and faster phase speed).While wind-waves are normally aligned with the local wind, the local wind direction is not always correlated with the swell direction.Sometimes winds can completely oppose the swell field and this is known to give rise to interesting situations with increased turbulence levels over the sea surface [3].We choose to study cases where the wind is aligned to the wind field and opposing the wind field, and compare these to a control run over a flat surface.The waves chosen will then need to be swell-like waves.These waves can also be fairly well represented by a sinusoidal shape.This wave influenced wind field will then interact with a wind turbine.For the wind turbine modelling part, the Simulator fOr Wind Farm Application (SOWFA) [4], developed at the National Renewable Energy Laboratory (NREL), was used.
In section 2, we give a brief description of the different models that are used and a description of the simulation characteristics.The results are presented and discussed in section 3 and we conclude and give suggestion for further usage in section 4.

Wave generation, actuator line model and FAST
All CFD simulations are done with the open source CFD toolbox OpenFOAM [5].To model the effect of waves on the wind flow there was a need for a method that could resolve the individual waves, including the direction of the waves relative to the wind.For this purpose, CFD with a moving grid approach was chosen.Several sinusoidal waves can be superposed on top of each other and implemented as a boundary condition on a patch in the CFD domain.A transient turbulent solver, able to handle deforming mesh, was used as a starting point.With the help of both OpenCFD b and Acona Flow Technology, a new solver that models flow above a moving wave surface was developed, this is named 'pimpleDyWFoam'.Each grid cell moves now up and down, reaching its maximum elevation at different time increments according to a sinusoidal function, and the movement looks like a wave propagating -much like the surface particles in a real ocean wave.In [6][7][8] a more detailed description of the model is available.
Churchfield [4] implemented Sørensen and Shen's actuator line method [9] in SOWFA.Here, the turbine rotor blades are represented as span-wise sections with airfoil characteristics.The model requires various input parameters.We have used the airfoil characteristics of NREL's 5 MW turbine.This turbine has a hub height of 90 m, rotor radius of 61 m and is assumed to produce power close to 2 MW at a wind speed of 8 m/s and a rotation speed of approximately 9 RPM [10].Using this as guidance other input parameters such as number of span-wise airfoil sections (segments) and the Gaussian width factor (related to the impact area for the force calculated in every blade segments) are determined after best practice, see e.g.[11][12][13].
In SOWFA it is possible to activate a coupling to FAST which is an aeroelastic computer-aided engineering tool for horizontal axis wind turbines [14].The actuator line CFD simulations will then replace the blade element momentum (BEM) part that is usually used with FAST.FAST will calculate the structural response and feed this back into the CFD simulations.Now various output parameters regarding blade loads can be analyzed and used for fatigue calculations.
b OpenCFD Ltd is owned by ESI-OpenCFD and they produce the OpenFOAM® open source CFD toolbox and distribute it through the OpenFOAM Foundation.

E3S Web of Conferences
The bases for these CFD simulations are the Navier-Stokes equations and the continuity equation for an incompressible Newtonian fluid [15].Our wave generation method was built from the OpenFOAM specific solver pimpleDyMFoam, and SOWFA is built from the OpenFOAM specific solver pisoFoam.Originally, SOWFA was set up for Large Eddy Simulations (LES), but we use an Unsteady Reynold Average Navier-Stokes (URANS) approach.The original solver for SOWFA was modified to include our wave generating method and in order to run using the URANS approach.As turbulent closure, we have used the standard k-epsilon model [16].These modifications of existing solvers and implementation of a new wave generation method have resulted in the combined WIWiTS setup.

Simulation characteristics
Several waves can be superimposed on top of each other to create different sea states, but we used only one wave, a swell like wave.We varied the period (T p ) and length (L).Following the dispersion relation in deep water, the wave speed c then depends only on the wavelength.We studied four different waves in combination with a wind field with wind speed of 8 m/s in 400 m height.The amplitude of the wave was held constant at 2 m.This is 30.5 times less than the rotor radius.The simulation cases are listed in Table 1.
Table 1.WIWiTS cases.The wave speed c will be aligned with the wind (+) or opposed with the wind (-).
Case (named after the period, T p ) Wave parameters The computational domain was 700 m x 260 m x 400 m with the turbine placed 550 m from the inlet (Figure 1a and 1b).Grid dependency studies performed on a two-dimensional setup (without turbine representation) showed that in order to have grid independent solutions, a very fine mesh is preferable.Also, the domain must be large to minimize boundary effects.This resulted in hundreds of millions of cells when used on an equivalent three-dimensional case.Simulations on such a mesh were not feasible.The mesh was instead constructed with a background-graded mesh having a refinement towards the wave surface and a refined area around the turbine rotor.This resulted in a cell size of approximately one meter near the wave surface inside the refined area.The domain length and possibly also the height are believed to be too short to avoid all boundary effects, but the simulations will nevertheless give indications on the relative differences in fatigue for different wave conditions compared to a no-wave situation.For pressure and velocities calculations, WIWiTS uses PIMPLE, which is a hybrid version of the PISO-SIMPLE algorithm [15].For the calculation of the gradient term and the Laplacian term in the Navier-Stoks equations, Gauss linear discretization schemes were used and for the convective term, the Gauss upwind scheme was used.The boundary condition used for the different paches are listed in Table 2.

Fatigue calculations
Fatigue damage is based on a combination of material properties and load history.The material property gives information regarding the maximum number of cycles that the material can withstand with a given stress amplitude.This is normally given as an S-N curve [14]; (1) where N is the predicted number of cycles to failure for the stress range ∆ߪ, ܽ ത is the intercept of the design S-N curve with the log N axis and m is the negative inverse slope of the S-N curve.A glass epoxy material for the blade is assumed, where the slope m = 9, and intercepts the curve, ܽ ഥ at 70 MPa.
The load history from the simulation will provide the bending moment, from which one can estimate the stress.From the stress history, the number of amplitudes at a given stress range can be estimated.For large simulations, the stress ranges are divided into blocks.The number of cycles with a stress range within the block will be n.The fatigue damage, D, is calculated using the Palmer-Miner rule: (2) where k is the number of stress blocks, ݊ is the number of stress cycles in block i and ܰ is the number of cycles to failure at a constant stress range ∆ߪ .A combination of the probability of the wind conditions over a year combined with the fatigue damage from a simulation will give the expected fatigue damage over a year, and the inverse of the damage would be the expected lifetime.In this study, the aim is not to investigate the lifetime of the structure, but rather the effect of the change in fatigue damage due to the wave influencing the wind profile.The damage equivalent load is a measurement of the fatigue damage load relative to a reference damage load.We will use the no 2nd Symposium on OpenFOAM in Wind Energy 02001-p. 5 wave simulation to estimate the reference damage load (see Table 1.).The definition of the damage equivalent load is: (3) where ܰ ாொ is the equivalent number of cycles in the simulation.Here, ܰ ாொ is the number of cycles needed for the damage equivalent load in the reference simulation (no waves) to be one.The stress at the blade root, σ, is estimated using the following relation: where M is the bending moment of the blade root, y is the outer radius of the blade root, and I is the moment of inertia of the circular cylinder at the blade root.Only the bending moments in flapwise direction are evaluated, as the gravitational force is dominating the bending moment in edgewise direction.The first 200 seconds of the simulation is discarded due to start-up transients.

Results and discussion
The relative difference in direction between the wind and the wave will affect the wind profile over the rotor swept area.In Figure 2 only horizontal wind profiles for one of the four wave cases are shown (T p =10 sec) as well as wind profiles for the no wave case.Profiles are sampled 300 meters from the inlet and at the rotor plane.The profiles in the no wave case are as expected, but the profiles in the situation with waves presents illustrates that the waves have an effect on the wind field and this will be dependent on the wind direction relative to the wave direction (blue and red profiles in figure 2).Since the aim of this study is to examine the possible effect of this fluctuating wind field on the wind turbine, we will in the following investigate stress on the blades and the tower.Profile of horizontal wind (m/s) 300 m from the inlet (upper) and in the rotor plane (lower) from case with wave period of 10 sec and from the no wave case (see Table 1).Wind and wave in the same direction (blue), wind and wave in the opposite direction (red) and no wave case (black).Profiles are sampled over approximately one wavelength.

E3S Web of Conferences
2nd Symposium on OpenFOAM in Wind Energy 02001-p.7 Figure 3 shows the stress at the blade root due to the flapwise moment after removing the start-up transients.The top graph is for the opposed wave condition, while the bottom graph is the aligned wave condition.Both graphs only consider the last 100 seconds in the fatigue calculations.It is clear that the stress variations are larger in the opposed wave condition relative to the aligned wave condition.This is true for all wave conditions with opposed wind investigated.However, regarding the aligned wind cases the stress is not notably affected.The reference case with no waves is included in both graphs, and it coincides very well with the simulations considering aligned wind.Therefore, we expect that the difference in fatigue will not be large for the aligned case.
The opposed wave condition shows larger variations in the stress range than the aligned condition.Especially the waves with long periods show larger variations in stress than shorter periods.Our simulations show a periodicity in the lower range of the stress for some of the wave conditions, as illustrated in Figure 5.
Figure 4 shows the stress at the base of the tower.Similar to the results for the blade root moment, the aligned case has smaller stress ranges compared to the opposed wave condition.The mean value correspond relatively well with the stress in the reference case with no waves.The mean stress from opposed wave condition is higher than the aligned wave condition and the no wave condition for the stress at the tower base.The stress range has a large variance for some wave periods, e.g.T p = 8 sec.

E3S Web of Conferences
02001-p.8The stress at the blade root and at the tower root will have different periods.This is illustrated in Figure 5, where the stress at the root of each of the three blades is shown in the upper graph.In the lower graph in Figure 5, it is shown that the local minima of the stress will oscillate with a period of 56 sec.This period is the multiple of the period of the rotor rotation, 7 seconds, and the wave period, T p , 8 sec.If the local stress minima for all three blades are considered, the oscillating period will be 8 sec, which is equal to the waveperiod.
2nd Symposium on OpenFOAM in Wind Energy 02001-p.9The frequency spectrum for the blade root stress in the flapwise direction is shown in Figure 6 for all wave frequencies.The frequency of one rotation, 1P, is 0.14 Hz. 2P and 3P are integers of this frequency.It is evident that the rotational frequency, 1P, contains the most energy, while the higher multiples of the rotational frequency, i.e. 2P and 3P, have less energy.1P has the highest energy as this is where the mean profile is sampled, while the higher multiples of the rotational sampling are due by the rotational sampling of the ambient turbulence.The wave frequencies applied vary from 0.10 Hz to 0.17 Hz.The frequency plot for the stress in the tower base is shown in Figure 7.The wave frequencies are visible in both the aligned and the opposed cases (wave period of 6-10 seconds corresponds to frequencies of 0.1-0.17Hz).The intensity of the response spectrum varies with the period of the wave applied, and most energy is found for wave periods close to the rotational frequency 1P (here: 0.14 Hz).The only peak that is found in the reference case with no waves is at the 3P (0.42 Hz).This peak is due to the vertical shear of the mean wind profile, and correspond to the peak at 1P (0.14 Hz) for the bending moment in flapwise direction in Figure 6.The intensity of this peak is relatively similar for all simulations.
The small peak at around 0.3 Hz in Figure 7 may be the eigenfrequency of the tower.For a fully dynamic aero-hydro-elastic-servo response analysis, the response at this frequency is normally larger and will sometimes dominate the spectra.The analysis presented here is limited to an analysis with little ambient turbulence and no hydrodynamic load, and the excitation of this mode is therefore limited.The equivalent fatigue damages, shown in Table 3 and Table 4, are estimated using Equation 3. The length of the simulation is very limited, only 100 seconds.The equivalent damage is estimated using the no wave case as a reference load case.
2nd Symposium on OpenFOAM in Wind Energy

Conclusions and further suggestions
It has been challenging to validate the results.To our knowledge there does not exist similar work that we can compare our results with, nor observations.We have reason to believe that the domain size is too small and we cannot guarantee that boundary effects are not present.In addition, grid independency was not completely reached.Nevertheless, we think that these results demonstrate that the wave will influence the wind field, which in turn affect the equivalent fatigue damage at both the blade root and the tower base.In a relatively low wind regime (reference wind of 8 m/s in 400 m height) the wave influenced wind increases the fatigue damage compared to a situation with no waves, especially for the cases where the wave opposes the wind field.Of the four wave cases (eight simulations) studied, the larger wave periods (8 and 10 seconds) give rise to the highest damage equivalent loads, and these cases also result in the highest peaks in the frequency spectrum.It is currently not possible to state whether the effect of wave influenced wind gives a significant response on the wind turbine compared to other effects that causes fluctuations in the wind, i.e. varying stratification and natural turbulence fluctuations.There is a need for longer stimulation times and more simulations in order to be able to conclude more specifically regarding the wave influenced wind impact on turbine performance and fatigue.Several sea states and wind regimes must be investigated, and it is also preferable to further develop the method here presented.Grid independency should be established and one should ensure that the computational domain is large enough.We have used an incompressible solver and one should consider implementing WIWiTS on a compressible solver.More realistic atmospheric conditions need to be implemented with buoyancy effect present.Different discretization schemes should be tested.The upwind schemes that we have used are known to be quite dissipative in URANS simulations, and a higher order scheme (e.g Quick) could generate more turbulence and hence change the results.

Figure 1a .
Figure 1a.Illustration of the WIWiTS domain.The WIWiTS mesh is graded near the surface and refined in a region around the turbine.

Figure 1b .
Figure 1b.An actuator line representation of the NREL 5 MW turbine is introduced in the moving wave domain.Shaded area shows the grid cells around the turbine.The horizontal velocities are shown as colour contours (m/s).

Figure 2 .
Figure 2.Profile of horizontal wind (m/s) 300 m from the inlet (upper) and in the rotor plane (lower) from case with wave period of 10 sec and from the no wave case (see Table1).Wind and wave in the same direction (blue), wind and wave in the opposite direction (red) and no wave case (black).Profiles are sampled over approximately one wavelength.

Figure 3 .
Figure 3. Stress at the blade root due to the flapwise bending moment at the blade root.

Figure 4 :
Figure 4: Stress at the tower base due to the tower bending in fore-aft direction.

Figure 5 .
Figure 5.The stress for the three blades for opposed wave and period 8 seconds.The upper figure shows the stress in solid lines for each blade, with the local minima marked with asterisks.In the lower figure, all the local minima are shown with asterisks, and a sinusoidal curve with period 56 sec is fitted to the minima of blade 1.For all local minima, a sinusoidal curve of 8 seconds is used.

Figure 6 .
Figure 6.Frequency spectrum for the blade root bending moment in the flapwise direction.

Figure 7 .
Figure 7. Frequency spectrum for the tower base bending moment in the fore-aft direction.

Table 2 .
Boundary conditions on the different patches for the WIWiTS, with OpenFoam specific naming.On the ground patch, the boundary conditions for U will be different in the case of a flat surface (fixedValue) compared to a moving wave surface (movingWallVelocity). *Full openFOAM specific name: atmBoundaryLayerInletVelocity (logarithmic wind profile with U 400m = 8 m/s and roughness length, z 0 = 0.0002 m)