Stochastic load effect characterization of floating wind turbine support structures

Achieving substantial reductions in Levelized Cost of Energy (LCOE) of floating wind turbines (FWTs) requires robust reliability assessment that accounts for inherent design uncertainties. A key aspect of such reliability assessment is the definition of limit states. In this regard, load effects need to be evaluated accurately. This paper presents a computational framework for evaluating load effects on FWT support structures. The computed load effect is subsequently characterized. A high fidelity finite element model of the National Renewable Energy Laboratory (NREL) 5MW reference turbine mounted on the OC3-Hywind spar buoy was developed and validated for this purpose. The loads from fully coupled time domain aero-hydro-servo-elastic simulations are transferred for detailed finite element (FE) load effect computation in Abaqus. Matlab® and Python are used as the computational tools for automating the whole analysis from start to finish. The initial part of this study addresses the amount of run-in-time to be excluded from response statistics. Based on convergence studies carried out, recommendations are made for run-in-time to be excluded from response statistics. The maximum von Mises stress in the tower as a measure of yielding is the load effect investigated in this study.


Introduction
Floating support structures are indeed the most viable option for deep water deployment of wind turbines. These structures are subjected to random excitation from both wind and waves as well as other loading sources. Designing floating wind turbine systems to withstand these loads throughout their service life at minimal cost requires robust engineering design that ensures the system is neither over-designed nor under-designed but just designed at the optimal level leading to substantial reductions in their LCOE. Support structure design drivers for an ultimate limit state (ULS) design are usually the absolute maximum of key parameters such as tower top/base shear forces and moments and tension at the fairleads and anchors obtained from dynamic simulations. From a reliability standpoint, the load effects in terms of excessive yielding and buckling, deformations, loss of stability etc., which result from the interaction of such loads as computed from coupled dynamic simulations are essential for defining realistic limit state functions. High fidelity finite element analysis (FEA) is very much suited for such robust load effect computation and is used in this study. The maximum observed von Mises stress throughout the tower is taken as the measure of yielding. The computation and subsequent characterization of this load effect is presented in section 4. High fidelity 3-D FEA models have also been employed by various researchers [1][2][3], in obtaining stresses on parts of wind turbines.
Part of this study throws light on the time to be excluded from response statistics in an ULS design where short 10 minutes simulations are carried out. This is crucial as the loads transferred from the aero-hydro-servo-elastic tool to the FEA model should not be characterized by start-up transients. Research findings on fully-coupled stochastic load effect computation and subsequent characterization are then presented. This work can be extended for developing uncertainty models towards more refined reliability assessment and partial safety factor calibration, leading to significant cost reductions in floating support structure design.
The wind turbine used for this study is the NREL offshore 5MW reference wind turbine [4] mounted on OC3-Hywind spar buoy [5]. The time series of von Mises stress is obtained using an interface between NREL aero-hydro-servo-elastic tool, OpenFAST (formerly known as FAST [6]) and FEA solver Abaqus with Matlab R and Python as the computational environment. Given the high computational cost of numerous FE simulations, High Performance Computing (HPC) cluster is used for this work.

Numerical model and tools
2.1. Coupled aero-hydro-servo-elastic time domain solver Dynamic simulations of the floating turbine are performed using the NREL aero-hydro-servoelastic code OpenFAST, which incorporates various modules like AeroDyn [7] for calculation of aerodynamic loads, a hydrodynamics module called HydroDyn [8], mooring dynamics modules (MoorDyn [9] is used in this study), structural dynamics modules SubDyn/ElastoDyn, and control/electrical system dynamics module ServoDyn. The wind loads acting over the rotor, gravitational loading as well as acceleration forces resulting from the contribution of the Rotor and Nacelle Assembly (RNA) in the full turbine system dynamics are reduced to 3 forces and moments at the tower top. The coupled turbine system structural dynamics are solved in OpenFAST which employs Kanes dynamics for solving the turbine system equations of motion. Details of the structural dynamics model employed by OpenFAST can be found in [10][11][12] and are not repeated here for sake of brevity.
2.2. FEA model 2.2.1. Creation of geometry and parts. Based on material properties and geometric dimensions in [5], a 3-D geometry of the OC3-Hywind is generated using Abaqus Python scripting. The model however is a reduced model without mooring system, transition piece and turbine blades. The RNA and platform are represented as rigid bodies with appropriate masses and rotary inertias assigned to their reference points. Only the tower (the part spanning 77.6 m i.e. from the tower top to 10m above the mean water level) is modelled as a deformable body.

Mesh generation and model validation.
Few elements are generated for the rigid parts (mesh sizes greater than 2 m) while the flexible tower has a mesh size of 0.5 m. The tower mesh size of 0.5 m is chosen based on good agreement between the tower natural frequencies extracted in Abaqus modal analysis and results obtained using BModes [13]. BModes is chosen for validation of the FEA model because data from BModes modal analysis is used in calculating polynomial coefficients needed in OpenFAST input files. The tower can be classed as a thinwalled structural member and therefore shell elements are used as the tower mesh element type. The shell thickness is represented by an analytical field defined by a simple mathematical expression of thickness as a function of the distance from the tower top.
It is imperative that the modal dynamics of the finite element (FE) model is similar to the model used in OpenFAST simulations from which loads are to be transferred. The loads or displacements outputted by OpenFAST already account for the effect of hydrodynamic the rigid RNA base plate and the tower top nodes. The tower base is assigned a fixed boundary condition. To account for the weight of the tower along with the tower inertial forces, the tower translational and rotational accelerations from OpenFAST together with gravitational acceleration of 9.80665 m/s 2 is applied to the tower as body forces in Abaqus. Similar analysis methodology was employed by Andrew et al [2] for evaluating stress utilization and buckling constraints using FEA in the optimization of composite towers for floating wind turbine concepts. Since the dynamics of the FWT system is already captured by OpenFAST, the problem is reduced to solving the equilibrium of all forces acting on the tower as represented by equation 1.
where K is the tower stiffness matrix, {D} is the vector of tower displacements and {F} is the vector of the forces acting on the tower. The displacements can then be calculated followed by strains and stresses. If tower displacements are known from aero-hydro-servo-elastic simulations, equation 1 can also be solved by specifying displacement boundary conditions.

Reference site and metocean conditions
The site chosen for this study is the East Coast of the United States of America. The realistic metocean data was created by Stewart et al [14]. Water depth at this site is assumed to be 320 m. The expected values of significant wave height (H s ) and median peak-spectral period (T p ) for each average wind speed is used as published in [15].

Eliminating Start-up transients
The initial part of dynamic simulations is often characterized by start-up transients especially for floating platforms. These transients exist due to the influence of gravity and rotor rotation on the structural displacements (assumed to be undisplaced at start-up) [16]. However, with structural damping, these start-up transients die out over time. It is pertinent to ensure that the response statistics are truly representative of the structural behaviour of the FWT before they are transferred to Abaqus for load effect computation or as such any other post-processing. IEC 61400-3 [17] recommends the removal of the first 5 s or longer from the simulation statistics to reduce the impact of start-up transients. This 5 s run-in-time clearly would not suffice considering the long natural periods of floating support structures. Although 60 s was suggested by Haid et al [15] if proper initial conditions (ICs) are set, no explanation or backing was given. Guzmn et al [18] investigated the use of moving average and backwards analysis of standard deviation with the later recommended for checking converged statistics. The simulations carried out in [18] were long simulations of more than 3600 s with a recommendation of 600 s as the run-in-time to be expunged. For ULS simulations that require short simulations of 600 s and where the averages of the time series are of less importance than the absolute maximum values, it becomes quite challenging to ascertain the amount of run-in-time to be discarded. To investigate this, simulations based on Design Load Case 1.2 (DLC 1.2) from IEC 61400-3 are performed. The turbine control systems are active for all simulations. The wind turbine is assumed to be class II as described in IEC 61400-1 [19]. A total of 11 wind bins are simulated (4 m/s to 24 m/s wind and associated sea states conditioned on mean wind speed). IEC Kaimal spectral model and a turbulence intensity of 0.14 is used for generation of turbulent wind realizations with TurbSim [20] and JONSWAP spectrum for generation of random seas. The TurbSim input files were modified to output 700 s turbulent wind time series comprised of repeated 100 s realization. The wave elevation time series also comprised repeated 100 s wave elevation realizations. This divides a 700 s simulation into 7 windows of which corresponding wind and wave inputs are equal. This enables convergence towards stationarity of responses within short 700 s simulations. The last window is used as the reference window for paired comparison as described by equation 2.
where d is the percentage difference between response x of window n, (n=1,2,..   figure 2b, the surge response converges to less than ±15% in less than 60 s. Similarly, for pitch response shown in figures 3a & 3b, more than 150 s is required for convergence within ±20% of the reference window for zero ICs simulations as opposed to proper ICs case of around 54 s (if the spikes caused by controller action in simulations of wind speed bins greater than the rated wind speed are ignored). Although the tensions in the fairleads converge quicker, the effect of proper ICs is very visible as seen in figures 4a & 4b.
Using ±20% as the convergence criterion, 50-60 s would suffice as the amount of run-in-time to be excluded from response statistics if proper ICs are set. The effect of second-order waves was neglected in this study considering second-order forces do not greatly influence the response of the OC3-Hywind spar platform [21]. For other floating support structures where second-order hydrodynamics cannot be neglected, longer window length should be used for setting ICs and more than 60 s run-in-time is suggested.

Statistics of load effect and characterization
For every time step in the OpenFAST simulations, loads and accelerations are transferred to Abaqus for von Mises stress computation based on the analysis model described in subsection 2.3. The maximum von Mises stress in the deformable tower for every time step is extracted using a Python subroutine. Results presented here are for DLC 1.2, mean wind speed, U = 12m/s, T p = 10s, and Hs = 6m. In accordance with IEC 61400-3 guidelines, 6 unique 600 s realizations of wind and wave are used, making up 1 hr long simulation time. The statistics of the von Mises stress for the 6 realizations are presented in table 2. Figure 5 shows the time series of maximum von Mises stress on the tower at each time step for realization 1.   The von Mises stress distribution on the tower at one time step can be seen in figure 7 with an indication on the location where this maximum value occurred. The aft region of the tower base experiences the most extreme von Mises stresses. It must be noted that the simulations in this work was for 0 • wind and wave direction.
The load effect is treated as a stochastic variable characterized by a Weibull distribution as shown in figure 6. For 6 realizations of 48001 time steps, a total of 288,006 Abaqus jobs were performed in this study. The highest observed von Mises stress in 1 hr simulation was 179.51 M P a corresponding to a stress utilization of 0.5716 if the tower allowable stress is taken as 314 M P a.

Conclusions
In this paper, fully coupled dynamic load effect computation was carried out using a validated high fidelity FE representation of the OC3-Hywind FWT. The initial part of this work focused on developing a methodology for automatically setting ICs for aero-hydro-servo-elastic time domain simulations for floating support structure concepts. Paired comparison with a reference window was used to check for convergence towards stationarity. It was found that with proper ICs, the chosen periodic output channels converged to less than ±20% of the reference values in less than 60 s. Consequently, about 50-60s can be considered sufficient run-in-time to be excluded from response statistics if proper ICs are set. Setting ICs for tower-top FA and SS initial displacements and the effect of second-order hydrodynamics was not part of this study and will be investigated in future.
With the computational interface presented in this paper, load effect comprising 288,006 data points was computed and fitted with a Weibull distribution. From a reliability standpoint, this study presents an approach that treats load effects as stochastic variables and could be  8 used in establishing uncertainty models for robust reliability assessment leading to calibration of currently used partial safety factors and cost reduction.