Evolution of the dynamic response and its effects on the serviceability of offshore wind turbines with stochastic loads and soil degradation

Accepted Manuscript

Evolution of the dynamic response and its effects on the serviceability of offshore wind turbines with stochastic loads and soil degradation Woochul Nam , Ki-Yong Oh , Bogdan I. Epureanu PII: DOI: Reference:

S0951-8320(17)30991-2 10.1016/j.ress.2018.03.017 RESS 6101

To appear in:

Reliability Engineering and System Safety

Received date: Revised date: Accepted date:

24 August 2017 26 January 2018 7 March 2018

Please cite this article as: Woochul Nam , Ki-Yong Oh , Bogdan I. Epureanu , Evolution of the dynamic response and its effects on the serviceability of offshore wind turbines with stochastic loads and soil degradation, Reliability Engineering and System Safety (2018), doi: 10.1016/j.ress.2018.03.017

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

ACCEPTED MANUSCRIPT

Highlights 

Dynamic response of wind turbines is considerably changed due to soil

CR IP T

degradation. The spectral density of stress of soil is studied for offshore wind turbines.



The stochastic behavior of soil stress for offshore wind turbines is

AN US



investigated.

A novel method is proposed to calculate the long-term soil degradation.



The potential of the proposed method for robust design is shown with case

AC

CE

PT

ED

studies.

M



ACCEPTED MANUSCRIPT

Evolution of the dynamic response and its effects on the serviceability of offshore wind turbines with stochastic loads and soil degradation

a

CR IP T

Woochul Nam†a, Ki-Yong Oh†b, and Bogdan I. Epureanu*c

Mechanical Engineering, Chung-Ang University, 84 Heukseok-ro, Dongjak-gu, Seoul, 06974, Republic of Korea [email protected] b

AN US

School of Energy System Engineering, Chung-Ang University, 84 Heukseok-ro, Dongjak-gu, Seoul, 06974, Republic of Korea [email protected]

c

Department of Mechanical Engineering, University of Michigan, 2350 Hayward Street, Ann Arbor, MI 48109-2125, USA [email protected]

M

†These authors contributed equally to this paper

ED

*Corresponding author: Tel.: +1-734-647-6391; fax: +1-734-615-6647

Abstract

AC

CE

PT

Novel methods combined with an integrated simulation platform are suggested for the design of offshore wind turbines (OWTs) and substructures that ensure a 20-year lifespan. These methods enable one to estimate the long-term evolution of the dynamic responses of OWTs due to the degradation of the soil modulus of the foundation under stochastic loading conditions. The results of this study show that random fluctuations of the soil stress caused by stochastic loads (i.e., aerodynamic and hydrodynamic loads acting on OWTs) can be described by a Rayleigh distribution and a Gaussian distribution. By using these probabilistic characteristics, the stochastic fluctuations in the soil stress can be rapidly calculated without using Monte Carlo simulations. Moreover, a new method based on the derivatives of the degradation functions and on the inverse of these functions is also suggested to calculate the mean degradation index. These methods significantly decrease the computational effort, thus overcoming a critical drawback of existing methods. Case studies demonstrate that the dimensions of the substructures significantly affect the evolution of the dynamic response. This suggests that the evolution of the dynamic response should be considered in the design process to secure the serviceability of OWTs and substructures.

ACCEPTED MANUSCRIPT

Keywords offshore wind turbine; soil modulus degradation; stochastic aerodynamic/hydrodynamic load; suction caisson foundation; wind turbine long-term serviceability

CR IP T

1. Introduction

The various advantages of offshore wind power including its high-energy density, low turbulence, low wind shear, and limited civil complaints have encouraged the renewable offshore wind energy sector to enter the large-scale development phase. In Europe, a new offshore wind power system was installed and connected to the grid in 2015, with a capacity of 3 GW, which is twice the capacity installed there in 2013 and 2014 [1]. Moreover, overall, 13 commercial wind farms are under construction, which will have a total capacity of over 4.2 GW [2].

ED

M

AN US

Overall, advances in offshore wind turbine (OWT) technology have resulted in significant increases in commercial wind farms. Particularly, large-size OWTs are economically desirable for wind farms because the foundation costs make up 20% to 45% of the total cost of OWTs, depending on the water depth [3]. For example, the average size of OWTs installed in 2015 and 2016 was 4.2 MW and 4.8 MW, respectively, whereas that of OWTs installed in 2011 was 2.78 MW on average [1-4]. Nonetheless, monopile substructures remain by far the most popular substructure type, even though jacket or tripod structures provide sufficient bearing capacity with relatively low weight. Monopile substructures were deployed on 97% of all OWTs installed in 2015. Moreover, the piling method of installation is still dominant, although the demonstration project in 2002 at Frederikshavn (Denmark) has shown that a suction caisson reduces the steel weight by half compared with a traditional monopile piling solution [5].

PT

The safety and reliability of OWT are not fully ensured because of the harsh offshore load conditions and the mechanical properties of the seabed, which is less robust compared to onshore ground. The design life of a wind turbine is at least 20 years, meaning that OWTs should generate electricity for more than 20 years. Therefore, uncertainty in the design of OWTs needs to be mitigated when constructing commercial offshore wind farms with new technologies.

AC

CE

The risk in the design of substructures for OWTs, especially the design of suction caissons, has motivated many researchers to study new types of foundations and to create many demonstration projects with a variety of different OWTs worldwide [5-9]. The installation of suction caissons in clays and sands has been phenomenologically investigated, and theoretical studies have been validated with several experiments [10-13]. Methods to calculate the friction on the caisson skirt using the cone resistance profile from cone penetration tests and modeling techniques in sand overlaid by clay have also been suggested to predict the suction installation [14,15]. The responses of OWTs under a variety of loading conditions have been studied [16-19]. The responses of caissons to vertical tensile loads have been investigated to prevent caissons from pulling out [20-23]. Many previous studies provide useful information on how to design suction caissons for OWTs. However, studies related to the long-term evolution of the dynamic response are still few, even though the evolution of dynamic response can significantly affect the serviceability of OWTs. Stochastic dynamic loads from wind and waves create periodic stress and strain in soils around the suction caissons of OWTs. Because

ACCEPTED MANUSCRIPT

fluctuations in stress and strain degrade the soil modulus, soil degradation should be investigated in the design of OWTs. Soil and water conditions are different for different sites. For example, the conditions in the Korean peninsula are very different from those in Europe. Most feasible sites in the Korean peninsula that have high wind and potential for large-scale development have a sea depth deeper than 20 m. Moreover, soft clay deposits (N value of 2-3) exist near the seafloor at candidate sites [24].

ED

M

AN US

CR IP T

These challenges motivate us to develop a novel method to predict the long-term evolution of the dynamic response under normal operating conditions. A high-fidelity three-dimensional (3-D) finite element model (FEM) of an OWT consisting of a tower, the substructures, and the soil is developed in this study. Moreover, stochastic aerodynamic and hydrodynamic loads acting on the OWT are created according to the International Electrotechnical Commission (IEC) standards [25,26]. Stochastic aerodynamic loads are created using the power spectral density (PSD) of the wind speed and the coefficient of thrust of the OWT. Stochastic hydrodynamic loads are created using the potential theory and an empirically obtained PSD of the wave elevation. To predict the stress in the soil, which causes the degradation of the soil modulus, the stress PSDs of the soil are calculated for all locations around the OWT foundation. Moreover, an analytical method is proposed to obtain the long-term probabilities of the soil stress. Finally, a novel method is proposed to effectively estimate the nonlinear soil degradation for the design lifespan of OWTs. The proposed method can replace Monte Carlo simulations, which require very intensive computation. In addition, the proposed method for estimating the degradation of the soil modulus in stochastic load conditions is based on an analytical solution. These features of the proposed approach allow significant reductions in the calculation time to predict the long-term evolution of the dynamic response of OWTs under stochastic load conditions. Hence, this new approach enables one to calculate the degradation of the soil over space and time. Case studies with different average wind speeds suggest that the long-term distribution of the wind speed should be considered to obtain the site-specific design of OWTs. Case studies with different dimensions show that an erroneous selection of substructure dimensions significantly affects the long-term serviceability of OWTs. This suggests that soil degradation is a very important factor in designing OWTs.

PT

2. Integrated simulation platform

CE

This section describes the finite element integrated simulation platform we developed. Section 2.1 describes the geometry, dimensions, and material properties of the model. Section 2.2 explains the methods applied to create aerodynamic loads on wind turbines and hydrodynamic loads on substructures.

AC

2.1. Finite element model 2.1.1. Offshore wind turbine and substructure A 3-D finite element model is created using the COMSOL Multiphysics software to predict the evolution of stress on the OWT and the soil caused by long-term stochastic aerodynamic and hydrodynamic loads. The model consists of a tower, a substructure, triple suction caissons, and the soil (Figure 1 (a)). The assembly of rotor blades, nacelle, and hub weighs 105 kg, and its mass moment inertia is 107 kg m2 [2729]. The mass and the mass moment of inertia are lumped at the top of the tower (Figure 1 (b)) as in previous studies [29,30]. The dimensions of the tower (Table 1) are based on the specifications of a 3-

ACCEPTED MANUSCRIPT

CE

PT

ED

M

AN US

CR IP T

MW class wind turbine [28]. The length of the substructure is 30 m, considering the candidate sites of the Korea peninsula [31], with a mean sea level of 20 m, and other important factors including tide differences, extreme wave height, and water splash [32]. Other dimensions of the substructure and suction caissons are selected according to previous studies [33,34]. The properties of the steel (Young's modulus = 210 GPa, Poisson's ratio = 0.27, and density = 7,850 kg/m3) are used for the tower, substructure, and suction caissons. A damping ratio of 1.5% is used based on estimated values from ‘rotor-stop’ tests with an OWT located in the North Sea [35]. This ratio is reasonable because the structural damping is typically of the order of 0.5% to 1.5%. The lower value of 0.5% accounts for pure material damping, whereas higher values include both structural damping in the material and energy dissipation at the joints [36]. Triangular shell elements are used for all the OWT components, with the thickness listed in Table 1.

AC

Figure 1. Overall OWT model: (a) components of the OWT; (b) 3-D model of the OWT with applied aerodynamic and hydrodynamic loads.

ACCEPTED MANUSCRIPT

4m ~60 m 4 cm 4m 30 m 12 cm 5-7 m 3-5 m 5 cm 1.5 m 5 cm 10 m

CR IP T

𝐷tw 𝐿tw 𝑡tw 𝐷𝑚 𝐿𝑚 𝑡𝑚 𝐷sc 𝐿sc 𝑡sc 𝐷𝑠𝑠 𝑡𝑠𝑠 𝐷tw

AN US

Table 1. Dimensions of the OWT model. Diameter of the tower Length of the tower Thickness of the tower walls Diameter of the main substructure Length of the main substructure Thickness of the main substructure Diameter of the suction caissons Length of the suction caissons Thickness of the suction caisson walls Diameter of the supporting structure Thickness of the supporting structure Distance from the center of the substructure to the suction caissons

2.1.2. Soil behavior

M

The width and the depth of the soil are designed to be large enough to fully account for the stress distribution in the soil. Prior investigations showed that a width 6 times the diameter of the suction caisson and a depth 3.5 times the length of the suction caisson are suitable for monopile foundations to avoid significant effects of the model boundaries on the results [37,38]. Thus, we used a soil computational domain following the same properties.

PT

ED

Tetrahedral solid elements are used to model the soil. The mechanical properties of the soil can be determined by several independent parameters. Table 2 shows the values of these parameters used in this study. The damping ratio of the soil varies over the strain amplitude [39]. The damping ratio 𝜁s of the soil is determined as 15% because this value corresponds to the range of dynamic amplitudes of the strain calculated from the model during the operating conditions. The Poisson’s ratio υs of the soil is estimated to be 0.45.

AC

CE

The mechanical properties of offshore soil can change significantly because of static loads (e.g., the weight of the structure). Because offshore soil is saturated with water, static loads can drain the water, which affects the soil properties. Therefore, changes in the strength and the modulus of the soil due to the weight of the structure should be accounted for to accurately obtain the dynamic responses of OWTs and the soil. It is worth noting that the effect of the weight is considerable, especially for OWTs with suction caissons. Because of the weight of the structure, the soil strength and modulus significantly increase under the suction foot where very large stresses are exerted when the structure vibrates. In the laboratory tests with soil samples, the shear strength 𝑆 of the soil linearly increased with the vertical confining stress σ′ vc [40-43]. Thus, the effective strength 𝑆 of the soil is calculated as 𝑐0 + Δ𝑆1 + Δ𝑆2, where 𝑐0 is the cohesion of the soil, which is determined to be 2 kPa. Δ𝑆1 is the increase in strength caused by the soil weight, and Δ𝑆2 is the increase due to the structure weight. Δ𝑆1 is equal to 𝛼1 σ′ vc , where σ′ vc = 𝛾 ′ 𝑧𝑠 , 𝛾 ′ is the effective unit weight of the soil (chosen as 𝛾′ = 12.6 kPa) and 𝑧𝑠 is the depth of the soil from the seabed. The coefficient 𝛼1 was determined as 0.2 according to experimental observations [40]. Next, a vertical static load corresponding to the weight of the structure is applied to the

ACCEPTED MANUSCRIPT

suction caissons to obtain the stress in the soil caused by the weight. Then, the increase 𝛥𝑆2 is calculated as 𝛥𝑆2 = α2

𝜎𝑥′ +𝜎𝑦′ +𝜎𝑧′ 3

3𝛼1 0 +1

, where 𝛼2 = 2𝐾

and 𝐾0 is the lateral earth pressure coefficient. 𝜎𝑥′ , 𝜎𝑦′ , and 𝜎𝑧′

CR IP T

are the effective normal stresses acting along the x, y, and z directions, respectively. These stresses are calculated by applying static vertical loads corresponding to the weights of the OWT on suction caissons. The lateral earth pressure coefficient can be calculated using 𝐾0 = 1 − sin𝜙, where 𝜙 is the internal friction angle. The internal friction angle for a candidate site varies from 26° to 35° [44]. The value of 35° is used for the internal friction angle to calculate the lateral earth pressure coefficient in this study. The shear modulus 𝐺 of the soil is linearly proportional to its strength [45] as 𝐺 = AGS 𝑆, and thereby AGS is determined as 200 based on the data obtained from the cone penetration test for a candidate site [44]. Figure 2 shows the spatial variation of the soil modulus on the cross section. The modulus of the soil increases over the depth due to the weight of the soil. Significantly, the soil under the suction foot is stiffer than the soil at other locations due to the weight of the structure. Because of this change in the soil properties, the soil under the suction foot considerably contributes to the total bearing capacity.

AC

CE

PT

ED

M

AN US

Bonding conditions are used at the contact surfaces between the soil and the suction caissons, assuming that deflections are small, and hence there is a negligible slip or gaps under normal operating conditions. Note that the IEC suggests verifying the design of OWTs in two cases: extreme load cases and normal operating cases. In extreme load cases, the maximum stress on an OWT should not exceed the design ultimate loads. In normal operating cases, the strain on an OWT should be small enough not to create a slip or gaps. Moreover, the fatigue life and the natural frequency of the OWT installed in the soil should be studied in normal operating cases for the design lifespan [26]. This study focuses on the evolution of the dynamic response over the lifespan (e.g., 20 years) under normal operating conditions. Analysis for the extreme load cases is beyond the scope of this paper. Future studies will include predicting the dynamic responses under the extreme load cases with a modified model that accounts for large deflections and nonlinear behaviors in contact surfaces.

Figure 2. Spatial distribution of the shear modulus in the soils.

ACCEPTED MANUSCRIPT

𝛾′ c0 𝜙 𝜁s υs 𝛼1 AGS

2.3. Stochastic dynamic loads acting on the structure 2.3.1. Aerodynamic loads

12.6 kPa 2 kPa 35° 15% 0.45 0.2 200

CR IP T

Table 2. Properties of the soil. Effective unit weight Cohesion Internal friction angle Damping ratio Poisson’s ratio Ratio of soil strength and vertical compressive stress Ratio of shear stress and strength

𝐹𝑎 (𝑡) =

AN US

When the blades of OWTs rotate because of the wind, the rotation creates a thrust force at the top of the tower. The direction of the aerodynamic load caused by the rotation of the blades is parallel to the rotation axis of the blades, and its magnitude 𝐹𝑎 can be estimated [29] as 1 2 𝑎

2 2

,

(1)

ED

M

where 𝑎 is the air density, is the radial length of the blades, is the wind speed at the hub, and is the coefficient of the thrust. Thus, thrust forces exerted on the structure vary with the site and the design of the turbine. The air density 𝑎 = 1.252 kg/m3 at the height of the hub for a candidate site was measured by the Herald of the Meteorological and Oceanographic Special Research Unit (HeMOSU) [46]. The radial length = 45.6 m of the blades and the value of were provided by Doosan Heavy Industries & Construction for the WinDS3000 site.

AC

CE

PT

The stochastic wind speed at the hub is obtained using the Kaimal spectrum model recommended by the IEC [25]. Details about are provided in the supplementary material. Figure 3 (a1) shows a stochastic wind speed created by using Eq. (S1) when avg is 11 m/s, where avg is the average wind speed at the hub for 10-minute time window measurements. Figure 3 (a2) shows the corresponding aerodynamic loads calculated by using Eq. (1). Note that the thrust force does not exceed a specific value, which is approximately 380 kN in Figure 3 (a2). This bound force corresponds to the rated wind speed for the WinDS3000 because the pitch angle of OWTs is controlled over the rated wind speed to reduce the excessive loads acting on the structure when the wind speed is too high.

M

AN US

CR IP T

ACCEPTED MANUSCRIPT

PT

ED

Figure 3. Dynamic loads acting on an OWT when avg = 11 m/s. (a1) and (a2) show the stochastic wind speed and the aerodynamic loads acting on the hub of an OWT, respectively; (b1) and (b2) show the wave elevation and the hydrodynamic loads acting on the submerged structure; (c1-c5) show the dynamic pressure acting on the main substructure at time t = 3, 25, 38, 55, and 70 seconds, respectively.

2.3.2. Hydrodynamic loads

AC

CE

The stochastic wave elevation can be created by using the empirical equation on the PSD of the wave elevation proposed by Pierson and Moskowitz [47], as shown in Figure 3 (b1). Then, the hydrodynamic dynamic pressure acting on the submerged structure of OWTs due to the wave with frequency f can be calculated as ( 𝑧 𝑡) = − where ̂ =

(

𝑧 )

|

𝐷𝑚

/2 0

= ̂ sin(

𝑐

𝑡) + ̃ c s(

and ̃ =

𝑡),

(2)

0,

by using the potential theory. w is the water density; is the wave number, which depends on the frequency of the wave; and 0 is a parameter of the wave, which is determined by the wave amplitude, wave number, and wave frequency. Details about Eq. (2) and the parameters in this equation are provided in the supplementary material. By using this equation, the dynamic pressure, which varies over the angle

ACCEPTED MANUSCRIPT

, time 𝑡, and depth 𝑧, can be calculated as shown in Figure 3 (c). Figure 3 (b2) shows the total hydrodynamic loads 𝐹 acting on the main substructure, which are obtained by integrating over its surface. Note that the dynamic pressure ̂ is used to predict the responses of the OWT and the soil. ̃ is constant over , and, hence, its integrated effect over is negligible because it is a small pre-stress applied in the radial direction to the substructure.

CR IP T

3. Analysis methods

AC

CE

PT

ED

M

AN US

This section describes the novel methods for predicting the evolution of the dynamic responses of the OWT. These methods provide significant gains in computational efficiency. Figure 4 shows the overall procedure for predicting the degradation of the soil modulus and to consider its effect on the response of the OWT and the soil.

Figure 4. Procedure for predicting the evolution of the dynamic response of an OWT and the soil. By calculating two frequency response functions and by using aerodynamic loads and hydrodynamic pressures, the PSD of the stress tensor is obtained (Section 3.1). Shear stresses vary over the plane and direction, which can be defined by two unit vectors denoted by 𝑛̅ (perpendicular to the plane) and 𝑣̅ (unit vector inside the plane, i.e., perpendicular to 𝑛̅). For every orientation, the fluctuations in the shear stress are characterized and used to calculate the degradation of the soil modulus (Section 3.2). The degradation

ACCEPTED MANUSCRIPT

CR IP T

can be calculated with the degradation coefficient defined in this work. Candidate degradation coefficients corresponding to each pair of vectors (𝑛̅, 𝑣̅ ) are obtained by using a new method to predict the degradation due to stochastic cyclic stress (Section 3.3). Next, the actual degradation coefficient is determined as the maximum value among the candidate degradation coefficients (Section 3.4). The local modulus of the soil is updated after each time interval ∆𝑡 to reflect the soil degradation. The process is continued/repeated sequentially at time intervals of duration ∆𝑡 . The calculations shown in the larger/outer box in Figure 4 are performed at every location in the soil. The calculations shown in the smaller/inner box are performed for each pair of vectors (𝑛̅, 𝑣̅ ). Note that the degradation of the soil modulus is assumed to be much slower than the dynamics of the OWT and the period of the fluctuations in the stress in the soil. Thus, the soil modulus is assumed not to vary over the time duration of ∆𝑡 when the frequency response functions and the PSD of the stress are calculated.

3.1. Dynamic response of the soil and structures

AN US

A time domain analysis requires a heavy computational effort to obtain the long-term evolution of the dynamic response. Thus, a frequency domain approach is used instead of a time domain analysis to predict the dynamic response at all locations in the soil. The stress tensor 𝜎𝑚𝑛 corresponding to frequency f acting on the soil at a location characterized by coordinates 𝑥, 𝑦, and 𝑧 can be obtained as 𝜎𝑚𝑛 (𝑥 𝑦 𝑧 ) = 𝜎𝑚𝑛 𝑎 (𝑥 𝑦 𝑧 ) + 𝜎𝑚𝑛 (𝑥 𝑦 𝑧 ) = 𝐹𝑎 ( )𝐻𝜎

𝑛𝑎

(𝑥 𝑦 𝑧 ) + ̂ ( )𝐻𝜎

𝑛ℎ

(𝑥 𝑦 𝑧 ),

(3)

M

where 𝜎𝑚𝑛 𝑎 denotes the stress tensor caused by aerodynamic loads and 𝜎𝑚𝑛 denotes the stress tensor caused by hydrodynamic pressures. 𝐹𝑎 denotes the aerodynamic force on the structure. 𝐻𝜎𝑖𝑗 𝑎 is the

ED

frequency response functions of the stress for a unit aerodynamic load acting on the structure. All functions in Eq. (3) are complex valued to represent the amplitude and the phase.

CE

PT

The phases of the various frequency components of the aerodynamic loads are equally distributed random variables in [0 π]. Likewise, the phase differences between the various frequency components of the aerodynamic loads and the various frequency components of the hydrodynamic pressures can also be considered as random variables in [0 π]. As a consequence, the amplitude of the response of the OWT and the soil at frequency f at a given location can stochastically vary over time, with values between |𝐹𝑎 𝐻𝑎 | − | ̂ 0 𝐻𝜎 𝑛 ℎ | and |𝐹𝑎 𝐻𝑎 | + | ̂ 0 𝐻𝜎 𝑛 ℎ |.

AC

To characterize the stress acting on the soil over the frequency, the PSD of the stress tensor is calculated. As an example, Figure 5 shows the PSD of 𝜎̃𝑥𝑥 (= 𝜎𝑥𝑥 /𝑆, stress normalized by the strength 𝑆) at a location in the soil near the suction caisson. The aerodynamic loads are dependent on the wind speed and the control strategy of OWTs. When the wind speed is lower than the rated speed, OWTs are controlled to generate the maximum energy. Because of this, the stress increases over the average wind speed when the average speed is small, as shown in Figure 5 (b1) and (b2). The OWT controls the pitch angle of the blades to prevent excessive aerodynamic torque and mechanical loads when the wind speed is faster than the rated wind speed [48]. This prevention attenuates the fluctuations in the stress for a high average speed, as shown in Figure 5 (b3)-(b5). In contrast to the aerodynamic loads, hydrodynamic loads are not actively controlled by the OWT. Thus, the PSD of 𝜎̃𝑥𝑥 caused by hydrodynamic loads monotonically increases over the wind speed because the wave elevation increases with the wind speed. Note that the

ACCEPTED MANUSCRIPT

ED

M

AN US

CR IP T

PSD of the wave elevation has a peak at 0.34, 0.21, 0.16, 0.12, and 0.10 Hz for average wind speeds of 5, 8, 11, 14, 17 m/s, respectively. Thus, the peak of the PSD of 𝜎̃𝑥𝑥 shifts over the average wind speeds (Figure 5 (b)).

CE

PT

Figure 5. PSD of 𝜎̃𝑥𝑥 acting on the soil at a location near the suction caissons with 𝐷sc = 6 m and 𝐷sc = 4 m. (a) depicts the location where the dynamic response is calculated; (b1)-(b5) show the PSD of 𝜎̃𝑥𝑥 𝑎 caused by aerodynamic loads (denoted by the dotted lines) and the PSD of 𝜎̃𝑥𝑥 caused by hydrodynamic loads (denoted by the solid lines) when the average wind speed is 5, 8, 11, 14, and 17 m/s, respectively.

AC

In summary, the PSDs of the stress tensor corresponding to aerodynamic and hydrodynamic loads were calculated for all locations in the soil. Therefore, the stress of all locations in the soil can be predicted fast with the dynamic characteristics of the two loads and the computed frequency response functions.

3.2. Characterization of random stress caused by stochastic loads Aerodynamic and hydrodynamic loads cause stochastically fluctuating stress in the soil. This stress changes the modulus of the soil, which can affect the dynamic response of the OWT. Thus, this stochastic fluctuation of the stress needs to be characterized to ensure the serviceability of the OWT over its lifetime. Particularly, the amplitude of the fluctuating stress should be considered because it was observed

ACCEPTED MANUSCRIPT

CR IP T

to significantly affect the degradation of the soil modulus in previous tests [49-53]. These observations suggest that the dependency of the soil modulus on the mean value of the cyclic stress is also considerable. The amplitude and the mean value are graphically explained in Figure 6. Therefore, both the amplitude and the mean value of the stress acting on the soil around the OWT are characterized using our model. However, our model assumes that the degradation over a cycle is not affected by the frequency of the cyclic stress because these effects are negligible when the frequency is between 0.1 and 1 Hz [54,55], which is the frequency range where aerodynamic and hydrodynamic loads acting on the OWT are considerable. Note that the stress is normalized by the strength 𝑆 of the soil in this study because the modulus of the soil is affected by this ratio of stress and strength. This ratio is referred to as the normalized stress and is denoted as σ ̃ (= σ/𝑆) for normalized stress and τ̃ (= τ/𝑆) for shear stress.

M

AN US

Observations [53-55] of soil samples under cyclic stress cannot be directly used to predict the soil degradation for OWTs because the orientation of the maximum shear stress remains the same over time in the laboratory experiments. Unlike in these experiments, the orientation of the maximum shear stress can stochastically vary over time in an actual OWT. Thus, it is difficult to capture the degradation of the soil around an OWT directly from the experimental data with soil samples. To address this challenge, the following procedure is applied to the soil around the OWT. Stress components (𝜎𝑥𝑥 𝜎𝑦𝑦 𝜎𝑧𝑧 𝜎𝑥𝑦 𝜎𝑦𝑧 , and 𝜎𝑧𝑥 ) at a location in the soil are calculated. Note that these stresses fluctuate over time because of the stochastic loads. Then, shear stresses at this location over various orientations are calculated. Next, the candidate degradation corresponding to each orientation is calculated. Here, an empirically obtained relation of cyclic shear stress and degradation is used. Then, the maximum degradation over the orientations is determined as the actual degradation at that location.

PT

ED

The orientation of a shear stress can be defined by a plane where the shear stress is acting and by a unit vector 𝑣̅ , which represents the direction of the shear stress in that plane. A unit vector 𝑛̅ perpendicular to the plane can be used to define the plane. Note that only three independent variables are required to determine both 𝑛̅ and 𝑣̅ , considering that they are orthonormal. The shear stress 𝜏𝑛̅ 𝑣̅ for a given pair of vectors 𝑛̅ and 𝑣̅ can be obtained as 𝜏𝑛̅ 𝑣̅ = 𝑣̅ 𝑇 [𝜎] 𝑛̅, where [𝜎] is the stress tensor. Therefore, the PSD of the shear stress 𝜏𝑛̅ 𝑣̅ can be obtained using the PSDs of the stress components and the unit vectors 𝑛̅ and 𝑣̅ .

AC

CE

To obtain the degradation by fluctuating shear stress, its amplitude and mean value needs to be characterized. Unlike in deterministic systems where the amplitude and mean value of the stress are constant over time, stress acting on the soil around an OWT has a time-varying amplitude 𝜏̃ 𝑛𝑎̅ 𝑣̅ and a timevarying mean value 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ because the OWT is excited by the stochastic loads, as shown in Figure 6. One 𝑎 can capture 𝜏̃ 𝑛̅ 𝑣̅ and 𝜏̃ 𝑛𝑚 ̃ 𝑝1 and τ̃𝑝2 in the time history as 𝜏̃ 𝑛𝑎̅ 𝑣̅ = ̅ 𝑣̅ by using two neighboring peaks τ ̃ 𝑝1 + τ̃𝑝2 )⁄ , as shown in the inset of Figure 6. |τ̃𝑝1 − τ̃𝑝2 |⁄ and 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ = (τ

CR IP T

ACCEPTED MANUSCRIPT

Figure 6. Example of the time history for normalized shear stress acting on the soil.

AN US

The normalized shear stress fluctuates with frequencies of the order of the first natural frequency of the OWT (which is approximately 0.35 Hz), as shown in Figure 6, because it has a narrow bandwidth peak at its first natural frequency. Thus, a large number of peak values need to be calculated to account for the evolution of the dynamics of OWTs in the long-term (approximately 20 years). When using Monte Carlo simulations, intensive computational efforts are indispensable. This intensive computational effort is the main reason why studies on the evolution of the soil modulus for OWTs are still few.

M

By using the characteristics of soil stress around OWTs in the frequency domain, we propose that this intensive computation for Monte-Carlo methods can be avoided by using probability distributions, namely, the Rayleigh distribution and the Gaussian distribution.

PT

ED

The PSD of the stress in soils can be characterized with two peaks, as shown in Figure 5 (b): one peak at the first natural frequency and the other peak at zero frequency. The frequency components around this first natural frequency determine the amplitude of the fluctuating stress because this peak has a narrow bandwidth. Moreover, the PSDs of wind speed and waves dramatically reduce over the frequency for high frequency. Thus, it is assumed that the effects of the high-frequency components (i.e., faster than 1 Hz) of stresses on the degradation of the soil modulus are negligible. Therefore, a PSD corresponding to the frequency components near the first natural frequency can be used to predict the probability distribution function (PDF) for the amplitude of the shear stress.

AC

CE

Next, the peak at near-zero frequency arising from the PSD of the wind speed has very large power at low frequencies. The low-frequency components are attributed to the mean value of the shear stress because they create slowly varying fluctuations in the time domain. Thus, a PSD corresponding to these frequency components can be used to predict the PDF for the mean value of the shear stress. The Rayleigh distribution can capture the PDF for the amplitude of a random variable when it has the PSD of a peak with a narrow bandwidth [56]. Thus, the PDF of the amplitude 𝜏̃ 𝑛𝑎̅ 𝑣̅ of the normalized shear stress can be obtained by using the Rayleigh distribution as PDF( 𝜏̃ 𝑛𝑎̅ 𝑣̅ ) =

𝑎 𝜏̃𝑛 ̅𝑣 ̅

𝑐𝑎

exp (−

𝑎 2 (𝜏̃𝑛 ̅𝑣 ̅)

2𝑐𝑎

),

(4)

where c𝑎 is the variance of 𝜏̃ 𝑛𝑎̅ 𝑣̅ . Generally, the standard deviation of a random signal can be obtained by integrating its PSD over the frequency from 0 to infinity. However, low-frequency components in the PSD of the shear stress create slowly varying fluctuations, which are negligibly attributed to the time-

ACCEPTED MANUSCRIPT

varying amplitude. Instead, they are correlated to the time-varying mean value of the stress. Thus, the variance c𝑎 can be obtained by integrating the PSD of 𝜏̃ 𝑛̅ 𝑣̅ from = 𝑐 to infinity. The frequency 𝑐 is determined as the lower bound of the half-power bandwidth of the peak at the first natural frequency.

PDF(𝜏̃ 𝑛𝑚 ̅ 𝑣̅ ) ≅

1 √2𝜋𝑐

exp (−

2 ( 𝜏̃𝑛 ̅𝑣 ̅ −𝐶0 )

2𝑐

),

CR IP T

The time-varying mean value 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ (𝑡) is determined by the sum of the low-frequency components of 𝜏̃ 𝑛̅ 𝑣̅ . Moreover, these low-frequency components can be considered as independent random variables because of their independent random phases in [0, 2π]. Thus, 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ can be regarded as the sum of several independent random variables. The general central limit theorem states that a PDF for the sum of random variables can be approximated to the Gaussian distribution if the random variables are independent of each other. Note that this theorem holds even when the random variables have different probability distributions. Hence, the PDF for the mean value of 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ can be determined as (5)

AN US

where 0 is the average value of 𝜏̃ 𝑛̅ 𝑣̅ over time and the variance 𝑐𝑚 can be obtained by integrating the PSD of 𝜏̃ 𝑛̅ 𝑣̅ from = 0 to = 𝑐 .

ED

M

To verify that these two distributions are applicable to the stress in the soil, the probabilities of 𝜏̃ 𝑛𝑎̅ 𝑣̅ and 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ obtained from Eqs. (4) and (5), respectively, are compared to the values obtained from the Monte Carlo method. In the Monte Carlo method, the time series are obtained by using the PSD of the shear stress, as shown in Figure 6. Then, the amplitude and the mean values are calculated with the neighboring peak values. Next, their probabilities can be obtained by counting the amplitude and the mean values. Figure 7 (b1) shows the probability of the amplitude of 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ , and Figure 7 (b2) shows the probability of 𝑎 𝜏̃ 𝑛̅ 𝑣̅ . This comparison demonstrates that the probabilities obtained from Eqs. (4) and (5) can successfully estimate the probability calculated using the Monte Carlo method with much faster computation (i.e., approximately 500 times faster than the Monte Carlo method).

PT

Finally, the probability of the amplitude and that of the mean value are integrated. Moreover, the effects of average wind speed avg are incorporated into the probability. Because avg is updated every 10 minutes [46], its probability distribution is required to study the long-term dynamics of OWTs and soil. The Weibull distribution with a shape parameter of 1.83 and a scale parameter of 7.73 is constructed based on the observation at the west coast of Korea [46]. The probability 𝑃𝜏̃𝑛̅ 𝑣̅ 𝑖𝑗 that the shear stress

AC

CE

fluctuates with the amplitude of (𝑖 − 1)∆𝜏̃ < 𝜏̃ 𝑛𝑎̅ 𝑣̅ ≤ 𝑖∆𝜏̃ and with the mean value of (𝑗 − 1)∆𝜏̃ + 𝜏̃ 0 < 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ ≤ 𝑗∆𝜏̃ + 𝜏̃ 0 can be calculated as

𝑃𝜏̃𝑛̅ 𝑣̅ 𝑖𝑗 = ∑𝑙𝑘

𝑉avg 𝑎 𝑖 𝑚 𝑗 𝑃𝑘 𝑃𝑘 , 1 𝑃𝑘

(6)

AN US

CR IP T

ACCEPTED MANUSCRIPT

M

Figure 7. Normalized shear stress acting on the soil near the suction caissons (shown in Figure 5 (a)) with 𝐷sc= 6 m and 𝐷sc = 4 m. (a) represents the PSD of the normalized shear stress 𝜏̃ 𝑛̅ 𝑣̅ ; (b1) and (b2) are the probabilities of the mean value of the normalized shear stress 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ and the amplitude of the normalized shear stress 𝜏̃ 𝑛𝑎̅ 𝑣̅ when the average wind speed avg is 11 m/s. Bars denote the probability calculated by collecting the peak values denoted as circles in (a). Dots denote the probability calculated using Eqs. (4) and (5); (c) shows the probability 𝑃𝜏̃𝑛̅ 𝑣̅ 𝑖𝑗 in Eq. (6). 𝑉

where 𝑃𝑘 avg is the probability of the average wind speed that is the probability that the amplitude is (𝑖 − 1)∆𝜏̃

ED

∆

𝑎𝑖 avg . 𝑃𝑘

avg 0 + ( < 𝜏̃ 𝑛𝑎̅ 𝑣̅ ≤

− 1)∆

avg

𝑖∆𝜏̃ when

< avg

avg

≤

avg 0

+

is larger than 𝑚𝑗

+ ( − 1)∆ avg and less than avg 0 + ∆ avg. This probability can be obtained using Eq. (4). 𝑃𝑘 is the probability that the mean value is (𝑗 − 1)∆𝜏̃ + 𝜏̃ 0 < 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ ≤ 𝑗∆𝜏̃ + 𝜏̃ 0 when avg is larger than avg 0 + ( − 1)∆ avg and less than avg 0 + ∆ avg . This probability can be obtained using Eq. (5). In this study, avg 0 , ∆ avg, ∆𝜏̃ , and 𝜏̃ 0 are determined as 1.5 m/s, 3 m/s, 0.02, and -1, respectively. Figure 7 (c) shows the probability 𝑃𝜏̃𝑛̅ 𝑣̅ 𝑖𝑗 at the location in soils depicted in Figure 5 (a). This calculation is

CE

PT

avg 0

AC

performed for every point of the soil to predict the degradation of the soil, which can be spatially varying.

3.3. Degradation of the soil modulus due to variable-amplitude and variable-mean cyclic stress The degradation index of the soil modulus 𝐺 (= 𝐺 /𝐺1 ) caused by cyclic stress with constant amplitude was measured in several experiments using soil samples [49,53,57-62], where 𝐺1 is the initial shear modulus of the soil before cyclic loads are applied and is the modulus after N-periods of cyclic loading are applied. In the laboratory experiments, various phenomenological degradation functions = ( ) were constructed to consider the effects of the amplitude (which was constant over time) and/or the effects of the mean value (which was constant over time) of the cyclic stress. However, experiments on the degradation of the soil modulus with variable-amplitude (or variable-mean value)

ACCEPTED MANUSCRIPT

CR IP T

cyclic stress are still few, making it difficult to develop an estimation method when variable-amplitude cyclic stresses are applied to structures. Lin and Liao [58] suggested a superposition method and demonstrated its effectiveness by comparing its predictions with the experimental data obtained by Little and Briaud [63]. However, this method requires a heavy computational effort, especially for stochastic loads, suggesting that a new method is needed to reduce the computational time for an OWT where the soil stress varies over time due to the stochastic nature of aerodynamic and hydrodynamic loads (Figure 6). In this section, a simple yet efficient new method is proposed to calculate the degradation of the soil modulus with stochastic loads. Then, the method is compared with the previous method (i.e., superposition method [58]).

3.3.1. Superposition method

M

AN US

The degradation index induced by time-varying cyclic stresses can be calculated by using a superposition method [42,58]. Figure 8 shows an example for the superposition method when two different cyclic stresses are applied. Only the stress amplitude changes in this example, whereas the mean value is assumed to be zero. The degradation function can be expressed as 𝐺 = ( ) = −𝜎 , where 𝜎 is the amplitude of the stress. For example, a degradation function 1 can describe the degradation ratio when the stress fluctuates with a constant amplitude 𝜎1 over the whole cycles. A degradation function 2 can describe the degradation ratio when the stress fluctuates with a constant amplitude 𝜎2 over the whole cycles. Note that this degradation function can represent several experimental data on the cyclic load [53,57,59-61].

PT

ED

First, consider a simple case wherein the amplitude of the stress is changed only once. Specifically, the degradation index 𝐺 is decreased by 1 cycles of cyclic stress with the amplitude 𝜎1 . Subsequently, the amplitude of cyclic stress is changed to 𝜎2 and the stress is applied to the next 2 cycles. The degradation also depends on the current modulus or current degradation index. Thus, when the amplitude of stress changes from 𝜎1 to 𝜎2 , one needs to calculate an equivalent number of cycle eq = 2−1 (𝐺 ts ), where 𝐺 ts = 𝐺 ( 1 ) is the degradation index at the transition instance of the cyclic stress, as shown in Figure 8 (a). Then, the soil is degraded with the degradation function 2 from eq to eq + 2 .

AC

CE

Next, consider the situation where the amplitude of the stress keeps changing and its sequences are random. The probability that amplitude 𝜎 = 𝜎1 is 0.5, and the probability that amplitude 𝜎 = 𝜎2 is also 0.5. Because degradation functions are nonlinear with respect to N, the trend of the degradation index 𝐺 depends on the sequence of the stress (Figure 8 (b1-b3)). The mean degradation index ⟨𝐺 ⟩ can be calculated by averaging the obtained evolution histories (Figure 8 (c)). One hundred sets of random sequences are used to obtain ⟨𝐺 ⟩ and their standard deviation, as shown in Figure 8 (c).

CR IP T

ACCEPTED MANUSCRIPT

AN US

Figure 8. Superposition method for calculating the mean degradation index ⟨𝐺 ⟩ . (a) shows the degradation when a single transition of stress from 𝜎1 to σ2 takes place; (b) shows the evolution of the degradation when the stress keeps changing and its sequence is random. The circles in (b1)-(b3) are examples of the degradation index . The bold line in (c) is the curve obtained by averaging 100 sets of the degradation index 𝐺 . Error bars represent standard deviations.

PT

ED

M

The mean degradation index ⟨𝐺 ⟩ can be a useful metric for predicting the degradation of the soil modulus for OWTs because of the characteristics of the sequence for the random stresses acting on the OWT and the soil. If most of the large-amplitude (i.e., 𝜎2 ) stresses are applied in the first half stage and most of the small-amplitude (i.e., 𝜎1 ) stresses are applied in the latter half stage, the degradation index 𝐺 rapidly decreases over N. In the opposite sequence of stresses, the degradation index slowly decreases over N. However, the stochastic load acting on an OWT does not exhibit these biased sequences but is applied with relatively fair sequences, especially when long-term responses are of interest. Therefore, the soil degradation associated with an OWT can be approximated to the mean degradation index ⟨𝐺 ⟩.

AC

CE

However, calculating ⟨𝐺 ⟩ with the superposition method requires a very large number of calculations because an equivalent number of cycle 𝑒𝑞 should be calculated whenever the amplitude of the stress changes and because several sets of evolution histories should be calculated to obtain the average ⟨𝐺 ⟩. Particularly, a large number (~ × 108 ) of cyclic loads are exerted on an OWT during its design lifetime considering that its first natural frequency is about 0.35 Hz. Moreover, the mean degradation index ⟨ ⟩ has to be calculated for all locations in soils and for all orientations of the shear stress. Therefore, intensive computation is required to account for the degradation of the soil modulus over spaces, orientations, and years, if the superposition method is used for OWTs.

ACCEPTED MANUSCRIPT

3.3.2. Newly proposed method To mitigate this heavy computational effort, a novel method is developed and used to predict the longterm degradation. First, the mean degradation index ⟨𝐺 ⟩ can be obtained by integrating its rate of change per cycle. The probability of each stress can be used as a weighting factor. For example, the amplitude of cyclic stresses is stochastically varying; the amplitude is either σ1 or σ2 , The sequence of σ1 and σ2 is random, but their probabilities are, respectively, 1 and 2 (= 1 − 1 ). The rate of change per cycle can 1

1

+

2

2

. Note that the current degradation index needs to be considered to estimate 1

the degradation. Thus, the rate ( )⟩). Likewise, the rate

has to be evaluated at an equivalent number of cycle 2

has to be evaluated at an equivalent number of cycle

( )⟩). Then, the mean degradation index ⟨ r ⟩ at ⟨

(

∗ )⟩

= ∫1

∗

⟨𝐺𝑟 ⟩

∗

𝑑 + 1 = ∫1 [

1

1

∗

= +

| 1

However, Eq. (7) still requires intensive computation because of change have to be obtained at every N.

eq 1

=

eq 2

=

can be calculated as 2

2

] 𝑑 + 1,

|

(7)

2

AN US

−1 1 (⟨𝐺 −1 2 (⟨𝐺

CR IP T

be obtained as

and

eq 1

eq 2

and the corresponding rate

This problem can be addressed by converting the -N graph into an N- graph (Figure 9). Note that Figure 9 (b) is generated by rotating the graph shown in Figure 9 (a). Then, the mean degradation index ⟨𝐺 ⟩ at = ∗ can be related with ∗ as ⟨𝐺𝑟 ( ∗ )⟩

𝑝1 ⁄

⟨𝐺𝑟 ⟩ 1

|

1 ⟨𝐺𝑟 ⟩

𝑑⟨𝐺 ⟩ + 1=∫1

𝑝 + 2 ⁄

1

M

⟨𝐺𝑟 ( ∗ )⟩

= ∫1

⟨𝐺𝑟 ( ∗ )⟩

= ∫1

𝑝1

|

2 ⟨𝐺𝑟 ⟩

ℎ1 |

⟨𝐺𝑟 ( ∗ )⟩

𝑑⟨𝐺 ⟩ + 1 = ∫1

ED

∗

1

+𝑝2

ℎ2 |

𝑝1 −1 ⁄ ℎ1 |

𝑑⟨𝐺 ⟩ + 1

2

1

⟨𝐺𝑟 ⟩ ⟨𝐺 ⟩ 𝑟

𝑝 + 2 −1 ⁄ ℎ2 |

𝑑⟨𝐺 ⟩ + 1 ,

(8)

⟨𝐺𝑟 ⟩ ⟨𝐺 ⟩ 𝑟

CE

PT

where 1−1 and 2−1 are the inverse functions of 1 and 2 , respectively. This method significantly reduces the computational time because it simply uses the degradation index ⟨𝐺 ( )⟩ as an integrating variable instead of calculating eq 1 and eq 2 at every cycle. Moreover, this approach is also beneficial because the integrating variable ⟨𝐺 ( )⟩ is the same for σ1 and σ2 .

AC

Eq. (8) can be generalized to consider the degradation due to 𝑛-amplitudes of stresses 𝜎𝑖 as ∗

⟨𝐺𝑟 ( ∗ )⟩

= ∫1

1 𝑑⟨𝐺 ∑𝑖 𝑝𝑖 𝑖′ ( 𝑖−1 (⟨𝐺𝑟 ⟩))

⟩ + 1,

(9)

where 𝑖 is the degradation function corresponding to an amplitude σ𝑖 , 𝑖 is the probability that the stress amplitude is equal to σ𝑖 , and 𝑖′ = 𝑑 𝑖 ⁄𝑑 . To derive this equation, a formula on the derivative of the inverse function (i.e., 1/( −1 )′ = ′( −1 )) is used.

CR IP T

ACCEPTED MANUSCRIPT

AN US

Figure 9. An example of the proposed method for calculating the rate of change for ⟨𝐺 ⟩. (a) shows the calculation of the rate of change using the rate of change of the degradation functions 1 and 2 ; (b) shows the rate of change calculated using the inverse functions of the degradation functions.

PT

ED

M

To verify the accuracy of the proposed method, the mean degradation index ⟨𝐺 ⟩ obtained from Eq. (9) is compared to the value obtained from the superposition method when 10 different amplitudes of cyclic stress are applied. The values for 𝑖 are also different over the amplitudes. ⟨𝐺 ⟩ obtained from Eq. (9) shows a good agreement with the value obtained from the superposition method, as shown in Figure 10 (a). To demonstrate the effectiveness of the proposed method, the computation time (𝑡prev ) for the superposition method and the time (𝑡new ) for the proposed method are calculated, and then the ratio (𝑡prev /𝑡new) is shown over the cycles in Figure 10 (b). The reduction in the computation time is more noticeable as the number of cycles increases. The number of cycles in the stress applied to OWTs per day is about 3 × 104 when its first natural frequency is 0.35 Hz. For this number of cycles, the proposed method is 5.8 × 104 times faster than the superposition method. Furthermore, the reduction in the computation time grows as the number of cycles, as shown in Figure 10 (b). Therefore, the effectiveness of this new method is more remarkable when it is applied to the problem of OWTs that have a long lifespan.

CE

Generally, the effect of amplitude and that of number of cycles are coupled in the function 𝑖′ ( 𝑖−1 ). However, if these effects can be decoupled as 𝑖′ ( 𝑖−1 ) = 𝐹⟨𝐺𝑟 ⟩ 𝐹𝜎𝑖 , where 𝐹⟨𝐺𝑟 ⟩ is a function of ⟨𝐺 ⟩ only,

AC

and 𝐹𝜎𝑖 is a function of 𝜎𝑖 only, Eq. (9) can be simplified as ⟨𝐺𝑟 ⟩ 1 1 ∫1 𝐹 (⟨𝐺 ⟩) 𝑑⟨𝐺 𝑝 𝐹 ⟨𝐺 ⟩ 𝑖 𝜎 𝑟 𝑖 1 𝑟 𝑖

= ∑𝑛

⟩ + 1.

This equation enables one to further reduce the calculation time if 1⁄𝐹⟨𝐺𝑟 ⟩ is analytically integrable.

(10)

CR IP T

ACCEPTED MANUSCRIPT

Figure 10. Comparison of the mean degradation index ⟨𝐺 ⟩ calculated from the superposition and the proposed method. (a) shows the calculated degradation histories; (b) shows the ratio of their computation times. The dotted line denotes the number of cycles that can be applied during a day.

AN US

The newly proposed method is applied to the OWT and the soil to predict the effects of the degradation in the soil modulus on the dynamic characteristics of OWTs. The degradation function should consider both the amplitude of the normalized shear stress 𝜏̃ 𝑛𝑎̅ 𝑣̅ and the mean value of the normalized shear stress 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ . Stress-controlled cyclic load tests conducted by Andersen et al. [53] are used because they observed the effects of amplitudes and mean values of the cyclic stresses on the degradation of the soil samples. Using their experimental data, we developed a degradation function as (11)

1.78

where 𝑎 = (4.73|𝜏̃ 𝑛𝑚 ̅ 𝑣̅ |

M

( ) = exp(−𝑎√ − 1)

+ 1)(𝜏̃ 𝑛𝑎̅ 𝑣̅ )5. For this degradation function,

′ ( −1 )

=

𝑎 2 ⟨𝐺𝑟 ⟩ , 2 ln⟨𝐺𝑟 ⟩

which can be

ED

split into a function 𝐹𝜎𝑖 (= 𝑎2 / ) and another function 𝐹⟨𝐺𝑟 ⟩ (= ⟨𝐺 ⟩⁄ln⟨𝐺 ⟩) . Moreover, 1⁄𝐹⟨𝐺𝑟 ⟩ is analytically integrable. Therefore, Eq. (11) can be transformed into (12)

PT

𝑛 ⟨𝐺 ⟩ = exp (−√∑𝑚 𝑎2 ( − 1)). 𝑗 1 ∑𝑖 1 𝑃𝜏̃𝑛 ̅𝑣 ̅ 𝑖𝑗 𝑖𝑗

𝑛 Here, 𝐷𝐺 is defined as 𝐷𝐺 = ∑𝑚 𝑎2 , and this coefficient is referred to as the degradation 𝑗 1 ∑𝑖 1 𝑃𝜏̃𝑛 ̅𝑣 ̅ 𝑖𝑗 𝑖𝑗

CE

coefficient because it determines the degradation. Note that this degradation coefficient 𝐷𝐺 is a spatially varying coefficient, and it depends on the stress acting on each location in the soil and on the orientation of the shear stress. Again, 𝑃𝜏̃𝑛̅ 𝑣̅ 𝑖𝑗 is the probability that 𝜏̃ 𝑛𝑎̅ 𝑣̅ is (𝑖 − 1)∆𝜏̃ < 𝜏̃ 𝑛𝑎̅ 𝑣̅ ≤ 𝑖∆𝜏̃ and 𝜏̃ 𝑛𝑚 ̅ 𝑣̅ is

AC

1.78 (𝑗 − 1)∆𝜏̃ + 𝜏̃ 0 < 𝜏̃ 𝑛𝑚 + 1]((𝑖 − 0.5)∆𝜏̃ )5 . ̅ 𝑣̅ ≤ 𝑗∆𝜏̃ + 𝜏̃ 0 , and 𝑎𝑖𝑗 = [4.73|(𝑗 − 0.5)∆𝜏̃ + 𝜏̃ 0 |

4. Case studies and discussion 4.1. Long-term degradation of the soil modulus The long-term degradation of the soil modulus is predicted with the proposed methods through the integrated simulation platform. It is observed that the degradation of soil modulus stops when the modulus decreased to 10-30% of its initial modulus [64]. Thus, it is assumed that the degradation of the

ACCEPTED MANUSCRIPT

modulus stops when the modulus decreases to 20% of its initial modulus. Moreover, the wind direction of a candidate site is assumed not to change significantly. This is a very critical load condition for an OWT because the soil around the caisson located opposite to the wind direction suffers from significant modulus degradation.

CR IP T

The degradation of the modulus is calculated for two different conditions: 1) for wind conditions where its average speed with 10-minute time window measurements is almost constant and 2) for wind conditions where its average speed with 10-minute time window measurements is stochastically varying.

AN US

First, the degradation is predicted when the average wind speed does not change over time. Note that instantaneous wind speed stochastically varies as Eq. (S1) even though the average wind speed ( avg) is fi xed. This calculation is performed to characterize the effects of average wind speed on the degradation of the modulus. The degradation is predicted using three different average wind speeds: 5, 11, and 17 m/s. Figure 11 shows the evolutions of the soil modulus around the suction caisson on the right when the OWT has three suctions with a diameter of 5 m and a length of 3 m. This suction caisson supports a higher load than the other suction caissons because of the direction of the wind. Thus, soil degradation is faster around this caisson than that around other caissons.

M

The effects of average wind speed on the rate of degradation are not monotonic. When the average wind speed is determined as 5 m/s, the degradation is very slow because dynamic loads are small, as shown in Figure 11. For an average wind speed of 11 m/s, the aerodynamic load is significant. As a result, the degradation rate is the fastest in Figure 11. The dynamic responses for the higher average wind speed than the rated wind speed (e.g., 17 m/s) are smaller than those corresponding to avg = 11 m/s because of the control strategy of OWTs as aforementioned. Thus, the degradation for avg = 17 m/s is slower than that for avg = 11 m/s, but is faster than the degradation for avg = 5 m/s.

AC

CE

PT

ED

The spatial variation of degradation can also be characterized from the proposed method. It is obvious that the stress acting on the soil near the suction caissons is larger than the stress in the soil away from the caissons. Among the soil near the caissons, the soil located farther away from the centerline of the tower is more critical to the degradation, as shown in Figure 11 (b1) and (b5). The same behavior is predicted in Figure 11 (b2) and (b7). Second, the rate of degradation depends on the depth, as shown in Figure 11 (b5)-(b7), because the stress and strength vary over the depth. The degradation is fast near the seabed (Figure 11 (b5)) because the strength of the soil is small. Note that the modulus of the soil in Figure 11 (b5) rapidly converges to the minimum value (i.e., 0.2). Soils around the tip shows relatively slow degradation (Figure 11 (b7)) due to the weight of the soils above them. Next, the degradation occurs faster for soils outside of the suction skirt than for soils inside of the suction skirt (Figure 11 (b2) and (b3)). The same characteristics can be observed in Figure 11 (b4) and (b6). This observation can be explained by the fact that the static compression from the weight of the structures increases the modulus of the soil inside the suction skirt, as shown in Figure 2.

AN US

CR IP T

ACCEPTED MANUSCRIPT

Figure 11. Evolutions 𝐺(𝑡)/𝐺(𝑡 = 0) in the soil modulus over the design life at locations around the suction caisson.

PT

ED

M

Next, the degradation is calculated when the 10-minute average wind speed stochastically changes. The probability of avg is required to calculate the probability of the stress using Eq. (6). As aforementioned, the probability is obtained from the Weibull distribution with a shape parameter of 1.83 and a scale parameter of 7.73. The degradation rate with varying avg is between the rate for a constant average wind speed of 11 m/s and the rate for 17 m/s. This suggests that the effects of average wind speed on the degradation is considerable. Therefore, an accurate statistical study on the wind speed is critical to the site-specific optimization design of OWTs. Examples of the evolutions in the modulus over 3-D spaces at different times are shown in Figure 12 when the OWT has three suctions with a diameter of 6 m and a length of 4 m. The change in the modulus

CE

is normalized with its initial value as

Δ𝐺( ) 𝐺( 0)

=

𝐺(

0)−𝐺( ) . 𝐺( 0)

Note that the changes in the modulus are

AC

calculated from the actual mean degradation index. The maximum change is 80% because the degradation of the soil modulus converges to 20% of the initial modulus. Because of the wind direction, the modulus of soils around the suction caisson located in the opposite direction degraded faster than that in other locations. The soil near the seabed shows a faster degradation because of its low strength and high stress concentration. Note again that the evolution over 20 years can be predicted with our method, which significantly reduces the computation.

CR IP T

ACCEPTED MANUSCRIPT

AN US

Figure 12. Spatial distribution of the degradation of the soil modulus over 20 years with stochastically varying average wind speed.

4.2. Serviceability of OWTs

PT

ED

M

To guarantee the reliability and safety of an OWT, we need to ensure that its natural frequency is not overlapped with the rotational frequency (1P) range of its turbine over its lifetime (i.e., 20 years). The deg radation of the soil modulus needs to be considered in the design of OWTs because the degradation of the soil modulus can shift their natural frequencies. To demonstrate the effectiveness of the proposed approach, the evolution in the natural frequency over time is predicted with different suction caissons (i.e., different lengths and diameters) (Figure 13). The rotating frequency of the blades for a 3-MW wind turbine ranges from 8.6 to 18.4 rpm [28]. Thus, the natural frequency should be higher than 0.31 Hz and be lower than 0.43 Hz. The initial natural frequencies are determined to be the same (i.e., 0.355 Hz) for different suction caissons by varying the length of the tower. Thus, the OWTs investigated in this work are reliable in their dynamic perspective when they are installed in the seabed. Figure 13 shows the evolution of the natural frequency over a year.

AC

CE

The evolutions of the natural frequency are different over the size of the suction caissons and over the average wind speed, although their initial natural frequencies are the same. First, the OWT with smaller suction caissons show a larger shift in the natural frequency because the stress concentration occurrs more severely in soils around the suctions. Second, the decrease in the natural frequency is negligible for a fixed average wind speed of 5 m/s, whereas the change is very significant for an average wind speed of 11 m/s. This effect of average wind speed is remarkable, especially for a small dimension of the suction caissons (i.e., 𝐷sc = 5 m and 𝐿sc = 3 m). The OWT with a suction caisson of 𝐷𝑠𝑐 = 5 m and 𝐿𝑠𝑐 = 3 m is not reliable after 8 years because its natural frequency overlaps with the 1P operating range of the OWT, as shown in Figure 13 (d). This result suggests that the OWT can fail after several years even if it operates properly at the initial and early stages. In contrast, two other dimensions show that the natural frequency of the OWT does not change significantly, meaning that these dimensions are reliable at the dynamic response perspective. Therefore, the proposed methodologies can be effectively used for the design and optimization of the OWT. Note that the design of the OWT needs an intensive verification for a variety of design load cases according to

ACCEPTED MANUSCRIPT

AN US

CR IP T

the IEC [26]. This study is focused on the dynamic load analysis and fatigue life of OWTs under normal operating conditions.

M

Figure 13. Changes in the natural frequency over the 20-year lifespan of an OWT.

5. Conclusion

AC

CE

PT

ED

The evolution of dynamic responses, i.e., the shift of the natural frequency, due to long-term cyclic loads is one of the important design metrics for OWTs. The overlap between the natural frequency of OWTs and the operating frequencies (1P and 3P) results in resonance and magnifies a significant stress and strain in OWTs, thereby reducing the serviceability of OWTs. To investigate the effects of soil degradation on the evolution of dynamic responses and the serviceability of OWTs, an integrated simulation platform with novel methods is developed. The proposed approach and methods significantly reduce the computation time, which is the main obstacle in predicting the long-term reliability of OWTs. First, frequency domain analysis and statistical methods are used to calculate the probability distribution of the soil stress. By investigating the aerodynamic and hydrodynamic loads and the structural characteristics of OWTs, we have revealed that the stochastic fluctuation of the stress can be modeled with a Rayleigh distribution and a normal distribution. Moreover, a novel method is suggested to calculate the nonlinear soil degradation. These methods reduce the calculation time significantly, thereby allowing us to predict the long-term evolution of the dynamic responses of OWTs. In addition, the effectiveness of this model is demonstrated with case studies. These studies suggest that the dimension of the substructures should be determined based on the prediction of the soil degradation. Therefore, this integrated simulation platform can be useful for enhancing the estimation of the lifespan of OWTs.

ACCEPTED MANUSCRIPT

Acknowledgements

AC

CE

PT

ED

M

AN US

CR IP T

This research was supported by the New & Renewable Energy Program of the Korea Institute of Energy Technology Evaluation and Planning (KETEP) grant (no. 20143010024330) funded by the Korean Government Ministry of Knowledge Economy and by the Chung-Ang University Research Grants in 2017.

ACCEPTED MANUSCRIPT

References [1] The European offshore wind industry – key trends and statistics. The European Wind Energy Association. (2016). Available from: https://windeurope.org/about-wind/statistics/offshore/keytrends-2015/.

CR IP T

[2] The European offshore wind industry – key trends and statistics 1st half 2016. Wind Europe. (2016). Available from: https://windeurope.org/aboutwind/statistics/off-shore/european-offshore-windindustry-key-trends-statistics-1st-half-2016/. [3] The European offshore wind industry – key 2011 trends and statistics. The European Wind Energy Association. (2012). Available from: https://windeurope.org/about-wind/statistics/offshore/keytrends-2011/. [4] Europe’s onshore and offshore wind energy potential: an assessment of environmental and economic constraints. European Environment Agency (2009). EEA Technical report No 6.

AN US

[5] L. B. Ibsen, M. Liingaard, and S. A. (2005) Nielsen. Bucket Foundation, a status. Copenhagen offshore wind. Copenhagen, Denmark. [6] B. Skaare, T. D. Hanson, and F.G. Nielsen. (2007). Importance of control strategies on fatigue life of floating wind turbines. 26th International conference on offshore mechanics and arctic engineering, San Diego, USA.

M

[7] A. Myhr and K. J. Maus. The twenty-first international offshore and polar engineering conference. Maui, USA; 2011: pp. 353-360. [8] H. Shin. (2011). Model Test of the OC3-Hywind Floating Offshore Wind Turbine. The twenty-first international offshore and polar engineering conference. Maui, USA, pp. 361-366.

ED

[9] T. D. Hanson, B. Skaare, R. Yttervik, F. G. Nielsen, and O. Havmøller, (2011). Comparison of measured and simulated responses at the first full scale floating wind turbine Hywind. EWEA Offshore. Amsterdam, Netherland.

PT

[10] G. T. Houlsby and B. W. Byrne. (2000). Suction caisson foundations for offshore wind turbines and anemometer mast, Wind engineering, vol. 24 (4), pp. 249-255.

CE

[11] B. W. Byrne, G. T. Houlsby, C. Martin, and P. Fish. (2002). Suction caisson foundations for offshore wind turbines, Wind Engineering, vol. 26 (3), pp. 145-155.

AC

[12] G. T. Houlsby and B. W. Byrne. (2005). Design procedures for installation of suction caissons in sand, Geotechnical Engineering, vol. 158, pp. 135-144. [13] G. T. Houlsby and B. W. Byrne (2005). Design procedures for installation of suction caissons in clay and other materials. Geotechnical Engineering, vol. 158, pp. 75-82. [14] M. Senders and M. F. Randolph. (2009). CPT-Based method for the installation of suction caissons in sand, Journal of Geotechnical and Geoenvironmental Engineering, vol. 135: pp. 14-25. [15] M. Senders, M. Randolph and C. Gaudin. (2007). Theory for the installation of suction caissons in sand overlaid by clay, International Offshore Site Investigation and Geotechnics Conference. London, UK, pp. 429-438.

ACCEPTED MANUSCRIPT

[16] J. P. Doherty, G. T. Houlsby, and A. J. Deeks. (2005). Stiffness of flexible caisson foundations embedded in nonhomogeneous elastic soil. Journal of Geotechnical and Geoenvironmental Engineering, vol. 131 (12): pp. 1498-1508. [17] M. Achmus, C. T. Akdag, K. Thieken. (2013). Load-bearing behavior of suction bucket foundations in sand, Applied Ocean Research. vol. 43, pp. 157-165.

CR IP T

[18] C. P. Aubeny, S. W. Han, and J. D. Murff. (2003). Inclined load capacity of suction caissons, International Journal for Numerical and Analytical Methods in Geomechanics, vol. 27, pp. 1235-1254. [19] J. H. Zhang, Z. Y. Chen, and F. Li. (2010). Three dimensional limit analysis of suction bucket foundations, Ocean Engineering, vol. 37, pp. 790-799. [20] W. Deng and J. P. Carter. (2002). A theoretical study of the vertical uplift capacity of suction caissons, International Journal of Offshore and Polar Engineering, vol. 12 (2), pp. 89-97.

AN US

[21] M. Zeinoddini and M. Nabipour. (2009). A parametric study on the pull-out response of suction caissons, Journal of Marine Engineering, vol. 5 (9), pp. 1/E-16/E. [22] K. Thieken, M. Achmus, and C. Schröder. (2014). On the behavior of suction buckets in sand under tensile loads, Computers and Geotechnics, vol. 60, pp. 88-100. [23] M. Achmus and K. Thieken (2014). Numerical Simulation of the Tensile Resistance of suction buckets in sand, Journal of Ocean and Wind Energy, vol. 1 (4), pp. 231-239.

M

[24] M. S. Ryu, J.S. Lee, D. J. Kwag, and S. Bang. (2016). Comparison of two meteorological tower foundations for off-shore wind turbines. Proceedings of the Twenty-sixth International Ocean and Polar Engineering Conference. Rhodes, Greece. [25] IEC 61400-1 Wind turbines-part 1: design requirements. 3rd ed. IEC. (2005).

ED

[26] IEC 61400-3 Wind turbinesd-Part 3: design requirements for offshore wind turbines. IEC. (2009).

PT

[27] S. Malhotra (2011). Selection, design and construction of offshore wind turbine foundations, INTECH Open Access Publisher. [28] Brochure of Vestas V90 3MW. Available from: https://www.vestas.com/en/products/turbines/v903_0_mw#.

CE

[29] S. Bisoi, S. Haldar. (2014). Dynamic analysis of offshore wind turbine in clay considering soilmonopile-tower interaction. Soil dynamics and earthquake engineering, vol. 63, pp. 19-35.

AC

[30] S. Bhattacharya, S. Adhikari. (2011). Experimental validation of soil-structure interaction of offshore wind turbines. Soil Dynamics and Earthquake Engineering, vol. 31 (5), pp. 805-816. [31] J.-Y. Kim, K.-Y. Oh, K.-S. Kang, and J.-S. Lee. (2013). Site selection of offshore wind farms around the Korean Peninsula through economic evaluation. Renewable Energy, vol. 54: pp. 189-195. [32] K.-Y. Oh, J.-Y. Kim, and J.-S. Lee. (2013). Preliminary evaluation of monopile foundation dimensions for an offshore wind turbine by analyzing hydrodynamic load in the frequency domain, Renewable Energy, vol. 54, pp. 211-218. [33] G. T. Houlsby and B. W. Byrne. (2004). Calculation procedures for installation of suction caissons. University of Oxford. Report No. OUEL2268/04.

ACCEPTED MANUSCRIPT

[34] L. B. Ibsen, M. Liingaard, and S. A. Nielsen. (2005). Bucket Foundation, a status,.Copenhagen offshore wind. Copenhagen, Denmark. [35] M. Damgaard, J. K. F. Andersen, L. B. Ibsen, and L. V. Andersen. (2012). Natural frequency and damping estimation of an offshore wind turbine structure, International offshore and polar engineering conference. pp. 300-307.

CR IP T

[36] R. Shirzadeh, C. Devriendt, M. A. Bidakhvidi, and P.Guillaume. (2013). Experimental and computational damping estimation of an offshore wind turbine on a monopile foundation, Journal of Wind Engineering and Industrial Aerodynamics. vol. 120, pp. 96-106. [37] M. Achmus, C. T. Akdag, and K. Thieken. (2013). Load-bearing behavior of suction bucket foundations in sand, Applied Ocean Research. vol. 43, pp. 157-165. [38] K. Thieken, M. Achmus, C. Schröder. (2014). On the behavior of suction buckets in sand under tensile loads, Computers and Geotechnics. vol. 60: pp. 88-100.

AN US

[39] J. Zhang, and R. D. Andrus, and C. H. Juang. (2005). Normalized shear modulus and material damping ratio relationships, Journal of Geotechnical and Geoenvironmental Engineering. vol. 131 (4), pp. 453-464. [40] C. C. Ladd and D. J. DeGroot. (2003). Recommended practice for soft ground site characterization: Arthur Casagrande Lecture. 12th Panamerican Conference on Soil Mechanics and Geotechnical Engineering. vol. 1, pp. 1-57.

M

[41] G. Mesri. (1975). New design procedure for stability of soft clays, Journal of Geotechnical and Geoenvironmental Engineering. vol. 100 (7), pp. 763-786.

ED

[42] Lefebvre, G and Ladd, CC and Mesri, G and Tavenas, F. (1983). Report of the testing committee, Committee of Specialists on Sensitive Clays on the NBR Complex, SEBJ, Montreal, Annexe I. [43] Vaid, YP and Campanella, RG. (1974). Triaxial and plane strain behavior of natural clay, Journal of Geotechnical and Geoenvironmental Engineering. vol. 100 (3), pp. 207-224.

PT

[44] K.-Y. Oh, W. Nam, M. S. Ryu, J.-Y. Kim, and B. I. Epureanu. A review of foundations of offshore wind energy convertors: Current status and future perspectives, Renewable and Sustainable Energy Reviews. under review.

CE

[45] Terzaghi, Karl and Peck, Ralph B and Mesri, Gholamreza. (1996). Soil mechanics in engineering practice, John Wiley & Sons.

AC

[46] K.-Y. Oh, J.-Y. Kim, J.-K Lee, M.-S. Ryu, and J.-S. Lee. (2012). An assessment of wind energy potential at the demonstration offshore wind farm in Korea. Energy. vol. 46 (1), pp. 555-563. [47] W. J. Pierson and L. Moskowitz. (1964). A proposed spectral form for fully developed wind seas based on the similarity theory of SA Kitaigorodskii, Journal of geophysical research. vol. 69 (24), pp. 5181-5190. [48] J. Jonkman, S. Butterfield, W. Musial, and G. Scott. (2009). Definition of a 5 MW reference wind turbine for offshore system development, Technical Report, NREL/TP-500-38060. NREL. [49] G. G. Moses, S. N. Rao, and P. N. Rao. (2003). Undrained strength behaviour of a cemented marine clay under monotonic and cyclic loading, Ocean engineering. vol. 30 (14), pp. 1765-1789.

ACCEPTED MANUSCRIPT

[50] A. S. Azzouz, A. M. Malek, and M. M. Baligh. (1989). Cyclic behavior of clays in undrained simple shear, Journal of Geotechnical Engineering. vol. 115 (5), pp. 637-657. [51] G. Lefebvre and P. Pfendler. (1996). Strain rate and preshear effects in cyclic resistance of soft clay, Journal of Geotechnical Engineering. vol. 122 (1), pp. 21-26. [52] Stewart and E. Harry. (1986). Permanent strains from cyclic variable-amplitude loadings, Journal of Geotechnical Engineering. vol. 112 (6): pp. 646-660.

CR IP T

[53] K. H. Andersen, A. Kleven, and D. Heien. (1988). Cyclic soil data for design of gravity structures, Journal of Geotechnical Engineering. vol. 114 (5), pp. 517-539. [54] A. M. Ansal and A. Erken. (1989). Undrained behavior of clay under cyclic shear stresses, Journal of Geotechnical Engineering. vol. 115 (7), pp. 968-983. [55] A. F. L. Hyde, K. Yasuhara, and K. Hirao. (1993). Stability criteria for marine clay under one-way cyclic loading, Journal of Geotechnical Engineering. ASCE. vol. 119 (11), pp. 1771–1889.

AN US

[56] J. S. Bendat and A. G. Piersol. (2011). Random data: analysis and measurement procedures, John Wiley and Sons. [57] I. M. Idriss, R. Dobry, E. H. Doyle, and R. D. Singh. (1976). Behavior of soft clays under earthquake loading conditions, Proceeding of Eighth Annual Offshore Techenical Conference, Offshore Technology Conference, Dallas.

M

[58] S.-S. Lin and J.-C. Liao. (1999). Permanent strains of piles in sand due to cyclic lateral loads, Journal of Geotechnical and Geoenvironmental Engineering. vol. 125 (9), pp. 798-802. [59] M. Vucetic and R. Dobry. (1988). Degradation of marine clays under cyclic loading, Journal of Geotechnical Engineering. vol. 114 (2), pp. 133-149.

ED

[60] A. R. Mortezaie and M. Vucetic. (2013). Effect of Frequency and Vertical Stress on Cyclic Degradation and Pore Water Pressure in Clay in the NGI Simple Shear Device, Journal of Geotechnical and Geoenvironmental Engineering. vol. 139 (10), pp. 1727-1737.

CE

PT

[61] K. Tabata and M. Vucetic. (2010). Threshold shear strain for cyclic degradation of three clays, Proceedings of 5th International Conference on Recent Advances in Geotechnical Earthquake Engineering and Soil Dynamics, Missouri Univ. of Science and Technology, Rolla, MO. [62] K. H. Andersen, W. F. Rosenbrand, S. F. Brown, and J. H. Pool. (1980). Cyclic and static laboratory tests on Drammen clay, Journal of the Geotechnical Engineering Division. vol. 105 (5), pp. 499-529.

AC

[63] R. L. Little and J. L. Briaud. (1988). Full scale cyclic lateral load tests on six single pile in sand, Miscellaneous paper GL-88-27, Geotechnical Division, Texas A&M University. [64] T. D. Stark and G. Mesri. (1992). Undrained shear strength of liquefied sands for stability analysis, Journal of Geotechnical Engineering. vol. 118 (1), pp. 1727-1747.