Development of Combined Load Spectra for Offshore Structures Subjected to Wind, Wave, and Ice Loading

: Fixed offshore wind turbines continue to be developed for high latitude areas where not only wind and wave loads need to be considered, but also moving sea ice. Current rules and regulations for the design of fixed offshore structures in ice-covered waters do not adequately consider effects of ice loading and its stochastic nature on fatigue life of the structure. Ice crushing on such structures results in ice-induced vibrations, which can be represented by loading the structure using a variable-amplitude loading (VAL) sequence. Typical offshore load spectra are developed for wave and wind loading. Thus, a combined VAL spectrum is developed for wind, wave, and ice action. To this goal, numerical models are used to simulate the dynamic ice-, wind-, and wave-structure interaction. The stress time-history at an exemplarily selected critical point in an offshore wind energy monopile support structure is extracted from the model and translated into a VAL sequence, which can then be used as a loading sequence for the fatigue assessment or fatigue testing of welded joints of offshore wind turbine support structures. This study presents the approach to determine combined load spectra and standardized time series for wind, wave, and ice action.


Introduction
As the demand for sustainable development and renewable energy sources surge, so does the need to continually improve technologies, design, and safety in the renewable energy industry [1]. Fixed offshore wind turbines (OWTs) continue to be developed for high latitude areas where not only wind and wave loads need to be considered, but also moving sea ice [2].
Based on current projections of climate change, the decline in icing frequency and sea ice extend will likely continue and tend to benefit the wind energy industry [3,4]; however, while sea ice extents are decreasing, wind parks are progressively moving into the northern hemisphere. In this harsh environment, OWTs are subjected to severe loading. The effects of ice-induced variable loads due to drifting ice sheets and the consequences for fatigue life of offshore structures are not yet fully understood [5][6][7]; nevertheless, it is generally known that loads from ice action on offshore structures can be much higher than those from wind and wave loading [8][9][10][11], and fatigue is the main failure cause for such structures [12].
Fatigue design of OWTs under variable-amplitude loading (VAL) is typically performed by means of stress-based fatigue assessment methods, linear damage accumulation hypothesis, and standardized load spectra or time series [13][14][15][16][17]. In this regard, standardized refers to the process of creating time series that represent a repeated sequence of loads, e.g., one flight journey, a year of loading of an offshore structure etc.
Alternatively, VAL testing can be based on real-time load sequences including load frequency information; however, such sequences are typically too long and complex to be suitable for structures exposed to high numbers of load cycles like ships and offshore structures [18]. Nevertheless, VAL fatigue assessment based on linear damage accumulation hypothesis also has a number of restrictions, see Heim [17] for an extensive summary. Most notably, for scenarios that differ from the usual loading scenarios, calibration using VAL testing is required to verify applicable damage sums. So far, the contribution of ice loads on fatigue damage is based on the knowledge on design for regions without ice loading. Thus, only simple damage superposition for the different loading conditions (i.e., wave and ice) was performed using fatigue design methods and curves not developed with ice loading in mind [19,20].
As the magnitude of loads due to ice action can be much higher than those related to wind and wave loading, differences in fatigue life under combined wind, wave, and ice action are expected. This is exacerbated by the fact that similar effects are well-known from overloading events in tension and compression, e.g., during storms [21][22][23].
To this day, there has never been any fatigue testing for combined wind, wave and ice load spectra. Thus, this study presents a design approach for OWTs that considers VAL due to wind, wave, and ice action. This is achieved by determination of a standardized load sequence from a numerical model that represents a case study for an OWT subjected to a number of drifting ice events during its lifetime. To this goal, the occurrence of different load cases on an OWT in the southern Baltic Sea was simulated to obtain the structural response of an exemplary monopile support structure as load measurements on OWTs under ice action are limited, see Nord et al. [24]. This data was then analyzed using the rainflow counting method to derive load spectra at an exemplary welded joint of the OWT. The presented method and results of this study are expected to serve as a foundation for future investigations on fatigue strength of OWTs under combined wind, wave, and ice action.

Fatigue Design of Offshore Wind Turbine Support Structures including Ice Loading and Sub-zero Temperatures
Structural design of OWTs is based on international standards like IEC 61400-3-1 [25] and design rules, e.g., DNVGL-ST-0126 [26]. For the estimation of ice loads (for ultimate strength and fatigue assessment), IEC 61400-3 [25] refers to the international standard for Arctic Offshore Structures of the Petroleum and natural gas industry ISO 19906:2019 [27]. This standard, however, focusses on static ice loads for ultimate strength assessment. Consequently, there is a need to test the applicability of state-of-the-art design methods for fatigue assessment of OWTs subjected to combined wind, wave, and ice action.
For design, the most important ice-structure interaction modes are intermittent crushing and the frequency lock-in [24,28]. Their interaction with the structure induces large global loads and displacements; however, their occurrence is typical limited in locations with mild ice conditions such as the Baltic Sea. Hence, designers of OWTs try to avoid such scenarios by preventive measures.
Rare large magnitudes loads, e.g., due to ice-structure interaction, are known to alter fatigue crack propagation by crack retardation or acceleration [29]. This effect is extensively studied for single or periodically load cycles; however, it is also known that events with a series of large cycles may affect crack propagation. A typical example of such over-load effects for ships and offshore structures are storm events, see [21][22][23]. Fatigue design methods for OWTs are supposed to design against crack initiation, i.e., a through thickness crack; however, for large plate thickness, this also includes a period of macro-crack growth, see [13].
Besides the effect of ice loading on fatigue strength, another important aspect for the design of OWTs in cold regions is the effect of sub-zero temperatures on fatigue strength [27,[30][31][32]. In a number of recent studies, an increase in fatigue strength with decreasing temperatures was observed [33][34][35][36][37]. This effect is currently not considered in any design standard. Under variable amplitude loading, small changes of the fatigue limit can have significant effects on fatigue life of OWTs. The current study, however, focusses on the determination of a novel combined VAL spectrum for an offshore wind energy support structure that accounts for wind, wave and ice action. The presented methodology is, nonetheless, applicable for different types of offshore structures and not limited to OWTs. The methodology will subsequently be presented based on a case study of an OWT support structure. Fatigue design of offshore structures is prone to large uncertainties. For OWTs, Velarde et al. [38] found that fatigue damage is more sensitive to uncertainties in wind, wave and soil parameters as compared to uncertainties in structural inputs. This study is not aiming at accounting for all uncertainties in fatigue design, but rather on presenting a general approach how to perform fatigue assessment for combined wind, wave, and ice action. For this purpose, a case study of an OWT located in the southern Baltic Sea was set up and loading conditions were chosen to resemble typical weather and ice conditions for this region.
For the case study, a circumferential butt-weld of an offshore wind turbine support structure at a height of about 75 m above the mean sea level was chosen. At this position, the tower section has a plate thickness of 17 mm and a diameter of 3.5 m. In order to obtain the structural responses at this location, numerical simulations were performed under the influence of the design load cases (DLC). This has been achieved using the structural finite element and aeroelastic model for OWTs (BHawC) and coupling it with the ice crushing model Variation of contact Area model for Numerical Ice Load Level Analyses (VANILLA) [39]. The brittle crushing load estimation of this model is based on the continuous brittle crushing equation of ISO 19906:2019 [27].
The ice model parameters in VANILLA are subsequently scaled to a pre-defined probability of exceeding the peak load FG in a simulated standard crushing event. This is defined by a drifting ice speed vice and an associated crushing length Lcrush to create a continuous brittle crushing on the structure. In the VANILLA project, the defined parameters were based on the Nordströmsgrund light-house (vice = 0.15 m/s and Lcrush = 90 m [39]). The ice-structure interaction in the VANILLA model has been validated with model-scale ice tests at Hamburg Ship Model Basin (HSVA) and full-scale data from Validation of Low Level Ice Forces (LOLEIF) and Measurements on Structures in Ice (STRICE) [40].

Design Load Cases
In the first pre-processing step, a critical load case assumption for the fatigue limit state must be defined. This has been done according to IEC61400-3-1 [25] for an exemplary monopile support structure in the southern Baltic Sea with site-specific wind and wave data, and estimated ice event conditions. These loads have been applied in two different modes of the turbine (idling and in operation). The defined scenarios are originally categorized in 11 DLCs, which are regrouped in four main DLCs. These four load conditions are illustrated in Figure 1. For simplification, waves are considered to be absent during ice events. Due to the differences in loading on the OWT and in aeroelastic damping, large differences in load levels are observed. Exemplarily stress-time responses at the investigated circumferential butt-welded joint are presented in Figure 2. DLC12 (a) results in much higher stresses than during turbine idling (b). This can be traced back to the boundary conditions of the load cases. The idling state (DLC64) represents wind speeds too low for turbine operation. On the other hand, the operational mode simulated high wind speeds. Additionally, a misalignment between wind and wave loading is considered. This leads to a fluctuation in the load history and a higher irregularity. To show a typical VAL sequence for each DLC, the stress-time series of DLC12 and DLC64 that occurs with the highest probability has been selected and presented in Figure  2 together with examples for DLCD4 (operation at an ice speed of 220 mm/s) and DLCD7 (idling at an ice speed of 90 mm/s). All four stress-time series are shortened to present 60s around the samples max stress. Due to the development of frequency lock-in in the second global bending mode of the structure, recognizable by the harmonic oscillation of the stress pattern, the average stress amplitudes are higher in DLCD7 than in DLCD4. It is noted that this is not an effect of just the difference between the aerodynamic damping between the two cases, as the velocity for the example from DLCD4 is so high that a different interaction regime develops. The two peaks in stress in Figure 2(c) are related to periods of time where the wind load causes the relative speed between ice and structure to become close to zero for a while, thereby resulting in a sudden short load build-up.
Furthermore, the substantial higher stress ranges in the idling mode (d) are observed. These differences can be explained by the low aerodynamic damping of the turbine while idling. This leads to higher displacement amplitudes and dynamic amplification in the structure caused by ice induced vibrations.
In total, 15 000 load scenarios with individual duration of 600s were simulated. From each simulation, bending moment histories about x-and y-axis of a defined position of the tower were extracted and used to determine the stress response at the chosen location. Each of the simulations is associated with a weight factor W that defines the occurrence probability of each scenario. These weight factors are based on Weibull distributions of environmental data for wind, wave and ice conditions and their misalignment. The total duration of the simulations sums up to 23.27 years.

Data Processing
Due to different probabilities of load scenarios (for instance different direction of wind and waves), the structural response around the circumference of the butt joint varies significantly. Thus, the methodology proposed by Milaković et al. [5], and presented in Figure 3, was adopted to determine combined load spectra. To this goal, the stress-time series obtained from the numerical simulations are assessed using rainflow counting according to the ASTM E1049-85 [41]. Next, the obtained range pairs are grouped into cumulative distribution functions (CDF) and superimposed. Finally, to shorten fatigue testing based on the combined load spectra, an omission level can be introduced to remove low load cycles that do not contribute significantly to fatigue damage, see [42,43].

Stress-Time Series Analysis
First, the data was reduced from the simulation states which do not occur. Second, the simulations are filtered to reduce the data points to relevant turning points and converting the moments into tensile and compression stresses in the structure. Third, the direction of each stress cycle is determined and weighted by the occurrence probability of the corresponding simulation, see Figure 4.
From Figure 4, the prevailing loading directions can be seen in addition to the number of load cycles for each DLC. Due to the predominate wind and wave direction, a large number of cycles are created in northwest direction within DLC12. Compared to the large number of load cycles for DLC12, the number of cycles associated with the ice load case in DLCD4 (operation) and DLCD7 (idling) are predominantly acting in one direction. This is related to the assumed ice load drift direction for the design location.  The edited load history for the twelve sectors and four load conditions are then assessed using rainflow counting according to the ASTM E1049-85(2017) [41]. In order to find the position around the circumference of the tower with the highest estimated fatigue damage, the loads are subsequently assessed in 30° sectors. To account for load cycles in neighboring sectors, trigonometrical relations are used to calculate the load contribution of each cycle within ±60°, see Figure 5.

Combined load spectra for wind, wave and ice action
Once the stresses are distributed into 12 different positions around the circumference, CDFs of stresses can be calculated for each DLC and 30° sector. Next, the stress spectra obtained for each DLC is added up to determine combined load spectra for wind, wave and ice action for each 30° sector. This is presented in Figure 6 as staircase function with 1 MPa steps.
Similar to Figure 4, the sectors between 150° and 210° are dominated by ice loading. This can be seen from the staircase function being separated in a region with high loads and low number of cycles and one with lower loads but higher number of cycles. These two regions are associated with ice-loading (DLCD4 and DLCD7), and wind and wave loading (DLC12).
For visual purposes and in order to permit fatigue damage calculation based on closed-form solutions, the individual stress spectra for each DLC in Figure 6 are fitted by Weibull functions. Only the parts of the stress spectra above the omission level are used for fitting to maximize accuracy in the relevant parts of the stress spectra.
It is well known that long-term geophysical processes such as stress spectra of ships and offshore structures can be approximated by Weibull functions [43]. To enhance comparability, all stress spectra are normalized by the maximum load cycles observed during the simulated lifetime of the structure. The parameters of the Weibull functions are summarized in Table 1. This simplification leads to some deviation from the original data; however, the fitted Weibull functions match the combined load spectra well in the regions of high loads and high number of cycles. For fatigue testing, the superimposed combined loading can be used to create a random stress-time series; nevertheless, for quick damage calculations, the fitted curves are suitable.  To reduce the number of cycles for the S-N testing an omission level is applied. The omission level was chosen to be 10% of the maximum stress. This significantly reduces the test duration if the stress spectra are used to determine suitable damage sums for this application. There are estimations that omission levels up to 50% of the fatigue limit do not significantly affect fatigue life under VAL [42]; however, for welded joints typically lower omission levels of up to 20% are applied [16,45,46]. For ships and offshore structures, where a large fraction of fatigue damage is related to cycles with low amplitudes, typically lower omission levels are recommended, see Fricke et al. [47]. They applied a level of 11% for the WASH spectrum [48,49]. Thus, an omission of 10% was assumed to be suitable in this study.

Derivation of Stress-Time Series for Fatigue Testing
In order to transfer the spectra into stress-time series that can be used for fatigue testing, random times series are required that reproduce the counted stress cycle. For this purpose, the Markov chain method is applied. This method is a stochastic method based on Markov transition matrices and is well suited to sample stress spectra from cumulative distribution function of stationary processes like wind and wave loading [43,50]. It thus captures the probability of different stress ranges and the number of occurred transitions, see also [49]. Based on this information, a random stochastic VAL series can be constructed. The reproduced load sequence represents upper and lower extrema which can be used as turning points for servo-hydraulic testing machines [51]. The transition matrices obtained for each 30° sector are presented in Figure 7 using 64 classes. Such matrices present distributions of stress ranges based on start class (lower stress) and end class (upper stress).  Again, the ice-related cycles lead to high stresses in the structure with low occurrence probability compared to the frequently acting wind and wave loads. Thus, cycles related to ice loading are found along the diagonal from bottom left to top right. Interestingly, differences in transition matrices can be traced back to the two different ice DLCs. Cycles from DLCD7 (see directions 150° to 210°) show the highest ranges, e.g., from class 10 to class 50, and are located closer to the diagonal from bottom left to top right. Cycles related to DLCD4 lead to the highest stresses in the 30° sector and have higher starting classes (higher mean stress level). This is related to different ice-structure interaction processes developing for idling (DLCD7) and operation (DLCD4) of the turbine. Cycles from wind and wave loading are found close to the diagonal from top left to bottom right.
The calculated matrices include the transition counts from the simulated data, with an applied omission of 10%. Due to the applied omission, the Markov matrices do not contain values along the diagonal from top left to bottom right. The reason for applying the omission before reproducing the random series is the vast number of cycles. This also mitigates the need to remove small cycles once time series were created.
The Markov chain method then creates a discrete sequence after a probability of transferring from one state into another state. A stochastically distributed probability of an occurring state can be selected by normalizing the matrix row by row and creating a single uniformly distributed random number in the interval (0,1) that is used to determine the starting state. After selecting one state, the associated column chooses the next row to find the subsequent stress state. The Markov method is a commonly used technique [51][52][53][54][55] and is also included in the international standard for VAL testing ISO 12110-1:2013 [56]. Figure 8 presents the stress-time series for one year. Therein, large ice loads are clearly visible. Furthermore, there are clear differences in length of the time series which are related to the prevailing wind and wave loading directions in western to northern direction (270° to 0°). In addition, differences in mean stress levels for DLCD7 and DLCD4 lead to predominantly tensile or compressive stress peaks. For example, high tensile stress peaks are clearly visible in the stress-time series for 30°, which are related to DLCD4, i.e., ice loading during operation. In contrast, during idling the global structural response results in a number of distinct compressive load peaks in 180° direction and neighboring sectors.
In order to determine the critical loading direction, each time series has to be evaluated separately. Due to the superimposed stress spectra, it is not straightforward to determine the most critical configuration based on length and height of a time series or spectrum.
These stress-time series could now be used for fatigue testing in order to determine suitable damage sums for combined wind, wave and ice action. This is, however, outside of the scope of this study. The main goal of this study was to present a methodology to determine combined stress spectra, and to highlight differences to typical fatigue loading of OWTs which is only related to wind and wave loading. A number of important aspects for design of OWTs for regions with seasonal ice cover and low temperatures are thus subsequently discussed.  6. Discussion

Sequence Effects of Variable Amplitude Load Sequences
In this study, ice, wind and wave amplitudes are randomly superimposed. This leads to an unnatural sequence of the loads, see Figure 9. This figure presents 50 turning points before and after the maximum load simulated for each direction. Thus, pronounced ice loads are in particular visible between wind and wave-related load cycles for the direction 150° to 210°. Typically, ice loads are only expected during a few winter days and not randomly distributed between wind and wave cycles; nevertheless, randomly distributing loads is the typical approach to achieve conservative test results, see [16,43]. The reason for this conservative approach is that the linear damage accumulation hypothesis does not account for sequence effects [17]. Furthermore, any sequence information is lost during rainflow counting; however, this is the state-of-the-art approach to handle such large datasets of loading for ships and offshore structures. In other engineering industries, VAL fatigue testing is sometimes also performed using combinations of harmonic and stochastic signals, see Decker [18]. This approach is, nevertheless, currently unfeasible due to several reasons. First, data of only a small number of measurement champaigns is available on real ice loading, and ice-structure interaction is known to vary significantly with partially high differences in load frequency, see Nord et al. [24]. This requires high sampling rates to sufficiently determine all load contributions, which is difficult in instrumented real structures. Second, the numerical and test basin representation of ice-structure interaction is an ongoing field of study, see [40,57,58], and third, offshore structures are exposed to load cycles in the order to 10 8 cycles during their lifetime. This can realistically not be applied to fatigue test specimens. Another approach to apply realistic seasonal loading to fatigue test specimens could be exposing the specimens to the four load series successively. This could for example be achieved by a method proposed by Li et al. [21]. Their idea could be adopted to ice load conditions; however, this can only be achieved without superimposing the different load conditions before generating the standardized stress-time series.
The idea of deriving a standardized load sequence with random order of loads was to derive a suitable basis for future fatigue testing of welded joints in OWTs to determine applicable damage sums, and to create a basis for comparison of combined load spectra for OWTs subjected to wind, wave, and ice action. Future tests could investigate the difference between random and successive ice load events for fatigue life of welded joints in OWTs.

Consideration of Temperature Effects on Fatigue Strength
Another important aspect of ships and offshore structures-operating in regions with seasonal ice cover-are temperature effects on material properties and in particular on fatigue strength. There are several recent studies that observed increases in fatigue strength with decreasing temperatures, as long as the temperature remains above the ductile-brittle transition regime, see [30,[32][33][34][35][36]. For example, at -20°C changes of fatigue strength of around 7% were determined from fatigue tests on different welded joints [30]. To this day, this effect is not considered in international standards and guidelines for fatigue design; however, there are studies that present approaches to consider temperature effects during fatigue assessment by means of nominal stress or local fatigue assessment methods, see [59][60][61].
This effect is of particular importance under VAL, as small changes of the fatigue limit can have significant effects on fatigue life of OWTs. The reason for this behavior is that most of the fatigue damage for an OWT is accumulated below the S-N curve kneepoint of its welded joints.

Summary and Conclusions
This study presented the development of a combined load spectra and standardized stress-time sequences for OWTs that account for wind, wave, and ice action. To this goal, the load interaction of an exemplary monopile support structure in the Baltic Sea was numerically assesses. This was achieved by simulating over 3000 different load conditions of 600s duration for four DLCs, and accounting for the occurrence probability of each scenario during the service life of the structure. From the simulations, bending moments at an exemplary location in the support structure of the OWT were determined and used to calculate stress responses.
Next, load spectra were obtained at different locations (30° sectors) around the circumference of the tower. Due to differences in prevailing loading directions, the number and magnitude of load cycles was found to vary significantly around the circumference of the tower. Based on the combined load spectra, it was possible to derive standardized stress-time sequences of one year duration, which can be used in future studies to determine suitable damage sums for OWTs subjected to wind, wave, and ice action. For this purpose, the Markov chain method and an omission level of 10% (based on the maximum stress response) was applied. This level was found to be a good compromise between reducing the cycles and existing knowledge on suitable omission levels. From the current study, the following conclusions are obtained:  Comparing the stress-time series, it becomes apparent that the wind and wave loads dominate the length of the sequences. This phenomenon can be explained by the significantly lower occurrence probability of ice events in the southern Baltic Sea which was used for the case-study. On the contrary, the highest stresses are caused by icestructure interaction. For other regions like the northern Baltic Sea, this behavior might be different.  Thus, combined stress spectra for wind, wave, and ice action on an OWT are very different from standard load spectra without ice contribution. The reason for this are very high ice-related load cycles that can cause overload effects similar to what is known from storm events. Applicable damage sums for such spectra should consequently be determined experimentally.  The sequence of load events is lost by rainflow counting. Thus, it is not possible to exactly reconstruct real-loading sequence; however, the Markov chain method is a suitable tool to create standardized random stress-time series. As high ice-related loads are randomly distributed in between wind-and wave-related loads, conservative results are expected if the stress-time series are used for VAL fatigue testing.  The reported Weibull distributions can be used for comparison fur future projects and to obtain fatigue damage based on closed-form solutions once applicable damage sums are experimentally determined.  Finally, the presented methodology is suitable for different types of offshore structures and not limited to OWTs. The main requirement is a verified numerical model of the structures that accurately depicts wave, wind, and ice action accurately.
The methods and conclusions of this work can be continued by forthcoming investigation of interactions between different load conditions on OWTs in regions with seasonal ice cover. Some influences are pointed out and a method to generate combined VAL series is elaborated. To extend the conclusions, a number of aspects worth investigating are mentioned. Further work will focus on the determination of applicable damage sums for combined wind, wave, and ice action.