Closed form solution of Eigen frequency of monopile supported offshore wind turbines in deeper waters incorporating stiffness of substructure and SSI

Offshore wind turbines (OWTs) are dynamically loaded structures and therefore the estimation of the natural frequency is an important design calculation to avoid resonance and resonance related effects (such as fatigue). Monopiles are currently the most used foundation type and are also being considered in deeper waters (430 m) where a stiff transition piece will join the monopile and the tapered tall tower. While rather computationally expensive, high fidelity finite element analysis can be carried to find the Eigen solutions of the whole system considering soil–structure interaction; a quick hand calculation method is often convenient during the design optimisation stage or conceptual design stage. This paper proposes a simplified methodology to obtain the first natural frequency of the whole system using only limited data on the WTG (Wind Turbine Generator), tower dimensions, monopile dimensions and the ground. The most uncertain component is the ground and is characterised by two parameters: type of ground profile (i.e. soil stiffness variation with depth) and the soil stiffness at one monopile depth below mudline. In this framework, the fixed base natural frequency of the wind turbine is first calculated and is then multiplied by two non-dimensional factors to account for the foundation flexibility (i.e. the effect of soil–structure interaction). The theoretical background behind the model is the Euler–Bernoulli and Timoshenko beam theories where the foundation is idealised by three coupled springs (lateral, rocking and cross-coupling). 10 wind turbines founded in different ground conditions from 10 different wind farms in Europe (e.g. Walney, Gunfleet sand, Burbo Bank, Belwind, Barrow, Kentish flat, Blyth, Lely, Thanet Sand, Irene Vorrink) have been analysed and the results compared with the measured natural frequencies. The results show good accuracy (errors below 3.5%). A step by step sample calculation is also shown for practical use of the proposed methodology. & 2015 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/).


Introduction
Offshore wind turbines (OWTs) are currently installed in high numbers in Northern Europe and the industry is rapidly developing worldwide. For example, countries such as China, South Korea, Taiwan and Japan are currently developing wind farms and the US has just installed its first pilot offshore wind farm. The Levelised Cost of Energy (LCoE) and SCoE (Society's Cost of Energy) from offshore wind are expected to decrease significantly in the near future and challenging projects in further offshore and deeper water depths (30-50 m), with heavier/larger wind turbines (upwards of 6 MW, with turbines of 10 MW in the testing phase) are underway. OWTs are expected to have a fatigue life of 25-30 years, while foundations are typically designed to last even longer (50 years).
It has been well established that offshore wind turbines are dynamically sensitive structures [10][11][12]3,35], and dynamics of OWTs must be studied to avoid unplanned resonance which may lead to increased fatigue damage through dynamic amplification of responses. Estimation of the natural frequency of the whole system is an important design calculation [25] so as to avoid the excitation frequencies arising from wind turbulence loading, wave loading, the rotational frequency of the turbine (1P) and the blade passing frequency (2P/3P). The installed number of these relatively new structures is likely to grow significantly in deeper waters and will pose additional dynamic challenges, and is therefore the focus of the paper. Specifically, the paper tackles the issue of the extra length of the tower in water (which is generally stiffer than the tower exposed in air) due to deeper water sites which makes the tower relatively more flexible as compared to shallow water installations thereby lowering the natural frequency of the whole system. Fig. 1 shows typical wind and wave power spectral densities for offshore sites which essentially describe the frequency content of dynamic excitations from wind turbulence and waves. The graphs also show the turbine's rotational frequency range (1P) and the turbine's blade passing frequency range (3P) for a range of commercial wind turbines of different capacities (2)(3)(4)(5)(6)(7)(8). A simplified way of translating the wind turbulence spectrum and the wave spectrum into mudline bending moment spectra through linear transfer functions has been presented in [8]. The peak frequencies of loading have been shown to coincide with the peak frequencies of wind turbulence and waves. However, loads at the 1P and 3P frequencies can also be observed through effects such as aerodynamic and mass imbalance of the rotor and rotational sampling of turbulence by the blades, for 1P and 3P loads, respectively [18,8]. Further details on the estimations of the loads can be found in [13].
For the widely used soft-stiff design, the target Eigen frequency (first natural frequency) is a frequency in the gap between 1P and 3P. Fig. 1 also shows measured natural frequencies of a wide array of OWTs from several different wind farms across Europe (see Table 1 for the nomenclature on the wind turbine structure for which the measured natural frequency is shown and Table 5 for the specification of the turbines given in Fig. 1). The figure clearly shows the trend that the target natural frequencies of heavier turbines are indeed closer to the excitation frequencies of the wave and wind making these structures even more sensitive to dynamics. Fig. 2 shows a schematic diagram of a wind turbine system along with the mechanical model relevant to the study of the overall system dynamics which is important for many design calculations such as Eigen frequencies and fatigue. The definitions of support structure, tower, substructure, foundation (monopile considered here), transition piece, mudline and mean sea level are defined in Fig. 2(a) and are consistently used throughout the paper. For clarity, the support structure is the whole structure that supports the heavy turbine i.e. the components below the rotor nacelle assembly (RNA) which includes the tower, substructure and foundation. The foundation is defined as the part of the support structure that is embedded in the ground below the mudline. The tower is typically a tubular tapered column. As the natural frequencies of these systems are very close to the forcing frequencies, the dynamics pose multiple design challenges and the scale of the challenges will vary depending on turbine types and site characteristics. These issues are briefly mentioned in the next section.

Dynamic issues in support structure design
There are two main categories of modern wind turbines based on the operational range: (1) Variable rotational speed machines have an operational speed range (1P range); (2) Constant rotational speed wind turbines that operate at a single rotational speed (fixed 1P).
Variable speed machines have 1P and 2P/3P frequency bands which must be avoided, as opposed to a single frequency (applicable to constant rotational speed machines). They are more restrictive from the point of view of foundation design as the forcing frequency is a band rather than a unique value. The turbines currently used in practice (under operation) are variable speed machines. Fig. 1 shows the 1P and 3P bands for typical offshore wind turbines ranging from 2 MW (not in production anymore) turbines to commercially not yet available 8 MW turbines. It may be observed from Fig. 1 that for most wind turbines, 1P and 3P excitation cover wide frequency bands. To avoid these ranges, the designer of the support structure has three options from the point of view of the stiffness of the tower and the foundation [32]: (1) Soft-soft structures: whereby the natural frequency lies below the 1P frequency range. This design is typical for floating offshore platforms, but are practically impossible to achieve for fixed (grounded) structures. Also, soft-soft fixed structures would likely be subject to high dynamic amplification of wave load response. (2) Soft-stiff structures: whereby the natural frequency lies above the 1P frequency range. This is the typical design choice for practical bottom fixed structures, such as monopile founded OWTs. All wind turbines in Fig. 1 are designed soft-stiff. (3) Stiff-stiff structures: whereby the natural frequency lies above the 3P range. Such designs are typically considered to have a massive support structure and therefore uneconomic. However, the Hornsea twisted jacket foundation supporting a met mast has recently successfully achieved a natural frequency above the 3P band [36].
It is quite clear from Fig. 1 that designing a soft-stiff system which avoids both the 1P and 3P frequency bands is challenging because of the tight tolerance of the target natural frequency. Indeed, some OWTs with a wide rotational speed range, such as the Areva M5000 5 MW in Fig. 1 with 4.5-14.8 rpm do not even have a gap between 1P and 3P and they are overlapping. Modern wind turbines often feature a sophisticated pitch control system Fig. 1. Typical wind and wave spectra, rotational speed (1P) and blade passing (3P) frequency bands for six commercial turbines and measured natural frequencies of the OWTs given in Table 1. designed to leap frog the rotational frequency range that would cause resonance of the structure due to the corresponding 3P frequency. For example, if the rotational speed range is 5-13 rpm, equivalent to 0.083-0.216 Hz, then the 3P frequency band is 0.25-0.65 Hz. If the structural natural frequency is 0.35 Hz, then the pitch control may regulate the rotational speed such that the corresponding 3P frequency avoids 7 10% of the structural natural frequency i.e. the frequency band 0.315-0.385 Hz to comply with the DNV code [25]. This would result in a jump in the rotational speed from 6.3 rpm to 7.7 rpm, avoiding operation between these values. Obviously, this is also not without cost, both in terms of initial costs, maintenance and power production.
The second issue is the top head mass, i.e. the mass of the rotor nacelle assembly (RNA, see Fig. 2 for clarity). With increasing capacity of the turbine, the mass of the rotor-nacelle assembly increases and the natural frequency will also decrease keeping other factors constant. Furthermore, the hub height is also higher due to longer blades (associated with heavier turbines) which increases the flexibility of the taller towers causing further reduction in the natural frequency. This trend can also be seen in Fig. 1, the target natural frequency of wind turbines (in soft-stiff design) with higher power output is typically lower. Fortunately, the large rotor diameter wind turbines also tend to have lower rotational speeds. This is because of the fact that for optimum power production the tip speed ratios (ratio of the speed of the tip of the blades and the incoming wind speed) are maintained at a favourable value, also to avoid blade damage. Furthermore, the top head mass is dependent on the choice of drive (with or without gearbox). Direct drive wind turbines do not need gearboxes, and thus significant weight reduction can be achieved (theoretically). However, the weight reduction is not obvious because the power of an electric generator is proportional to its rotational speed, see Eq. (1) [18] and thus direct drive wind turbines typically have larger diameter (and thus heavier) generators.
where C g is a constant and D g ; L g ; Ω g are the diameter, length and rotational speed of the generator, respectively.
As the total mass of the rotor-nacelle assembly influences the natural frequency, the design of the support structure for wind turbines may be improved by an integrated design.
A further effect is the increasing importance of the flexibility of the monopile, the grouted connection and the transition piece (that is, the substructure) as a consequence of increasing water depth. The longer the monopile, the more effect the substructure has on the natural frequency of the total structure. This is in contrast with shallow water turbines where the tower is overwhelmingly softer than the substructure and the natural frequency of the whole system is critically governed by the bending stiffness of the tower and the soil stiffness.
For large wind turbines the frequency band of the wave excitation (see JONSWAP spectrum in Fig. 1) becomes important as it is very close to the natural frequency of the whole system. This issue gains even more importance as the substructure diameter increases for large wind turbines. If the structure is also installed in deeper water, then both the expected significant wave height and the length of the substructure increase, resulting in dynamic wave loads comparable to or exceeding the wind load.
The natural frequency of the whole system therefore needs to be assessed in the feasibility study phase and also in the early phases of the foundation design. This calls for a simplified methodology that is quick and uses only limited data about the wind turbine and the site. Table 1 Analysed offshore wind farms with the used wind turbines and soil conditions at the sites.

Motivation behind the study
This study aims to provide a simple, quick and reliable tool for the foundation designer to assess the fundamental natural frequency of the turbine-substructure-foundation-soil system. The simple expressions obtained in this work build on previous research in the dynamics of offshore wind turbines. Analysis of OWTs modelled as Euler-Bernoulli beams and founded on two independent springs (lateral and rocking stiffness) is carried out in [4,3] while soil-structure interaction was experimentally verified in [10]. A major discrepancy between measured and predicted natural frequencies of OWT structures has been reported in a multitude of studies [21,28,32,36,54,56]. The sources of uncertainty related to the discrepancy can be due to theory or physical idealisation of the system or parameters or a combination of the above and are summarised below: (1) Beam model uncertainty i.e. uncertainty in the mathematical idealisation of the physical system: The current methodology is based on the Euler-Bernoulli beam theory, however, the authors have found in [6] that more sophisticated beam models including rotary inertia (Rayleigh beam) and shear deformation (Timoshenko beam) do not provide notable improvements in the results (less than 0.1%). Therefore, the simple Euler-Bernoulli model is satisfactory for wind turbine towers. Further details are given in Appendix A. (2) Uncertainty in the modelling of the substructure: The substructure is modelled by a single equivalent cross section (which is the cross section of the monopile). This approximation may lead to over or underestimation of substructure stiffness. (3) Uncertainty in the modelling of the foundation. The analysis focuses on Eigen solutions of a linear system, and therefore nonlinear soil behaviour is not modelled. The foundation system is not expected to go into the nonlinear regime and therefore the linear approximation is considered acceptable. In the linear regime, the foundation is modelled with three coupled springs (K L , K R and K LR ) resulting in an exact model of the linearised system, see Fig. 2(c). It was found in [6] that the cross coupling spring (K LR ) cannot be excluded from the analysis. (4) The major uncertainty is the well-known issue of measurement of soil stiffness. It is challenging to estimate foundation stiffness with complex soil layers at offshore sites. A sensitivity analysis was carried out [7] and it was found that underestimation of soil stiffness by 30% would typically result in a frequency increase of 0.2-5% while overestimation by 30% causes a natural frequency change of 2-5.6%. The uncertainty analysis was carried out by considering the effect of stiffness changes for the ten wind turbines presented in Table 1, as well as 6 other industrial wind turbines not discussed in this paper. The transcendental equations resulting from the analytical formulation have been solved, using the three springs approximation for the foundation and Timoshenko beam theory for the tower, as given in Appendix B, and also in [6] (a simplified version of the problem has also been given in [4,3]). The equations have been numerically solved for all wind turbines in different soil conditions to analyse the sensitivity. (5) Another aspect is the idealisation of the tapered tower of varying wall thickness with an equivalent diameter and wall thickness as a constant bending stiffness column, as has been used in many of the simplified methods [3,4,50,55,56,6]. (6) A further source of uncertainty is the damping of the real structure. If the damping ratio is ζ, the natural frequency of underdamped vibrationsf 1 can be expressed with the undamped frequency f 1 as The main sources of damping of an offshore wind turbine structure are the aerodynamic damping, soil damping, structural damping (damping of steel), hydrodynamic damping and additional damping from tower motion dampers (e.g. slosh dampers).

Damping of structural vibrations of offshore wind turbines
As the natural frequencies of offshore wind turbines are close to forcing frequencies, damping is critical to restrict damage accumulation and avoid premature maintenance. Therefore discussion is warranted for issues related to damping and for practical design purposes can be idealised in fore-aft and side-to-side vibrations of offshore wind turbines. The main difference between the sway-bending (or rocking) vibrations about two axes (X and Y in Fig. 2) is that in the along-wind direction higher damping is expected due to the high aerodynamic damping caused by the rotating blades interacting with the airflow. On the other hand, for the side-to-side direction, the aerodynamic damping is orders of magnitudes lower.
A non-operational (parked or idling) OWT has similar aerodynamic damping in the fore-aft as in the side-to-side direction. Since the wind load is acting in the along-wind direction, the highest load amplitudes are expected in fore-aft motion. This is because for most wind turbines in water depths less than 30 m, wind loading is the dominant load, while for very large diameter monopiles in medium to deep water, wave loading is expected to have equal or higher magnitude. It is worth noting that due to the yaw mechanism of the wind turbine, the along-wind and crosswind directions are dynamically moving and are not fixed, therefore the foundation is subjected to both cross-wind and along-wind loading in all directions during the lifetime of the OWT. One can conclude that analysing vibrations in both directions is important.
Studies considering the damping of the first bending mode either empirically or theoretically include Camp et al. [19], Tarpjohansen et al. [48], Versteijlen et al. [51], Damgaard and Andersen [21], Damgaard et al. [22] and Shirzadeh et al. [46]. Based on these studies and other estimates and in the absence of other data, the following assessment of damping ratio contributions is recommended: Structural damping: 0.15-1.5%. The value of structural damping depends on the connections in the structure (such as welded connections, grouted connections, etc.) in addition to material damping (usually steel) through energy dissipation in the form of heat (hysteretic damping).
Soil damping: 0.444-1%. The sources of damping resulting from soil-structure interaction (SSI) include hysteretic (material) damping of the soil, wave radiation damping (geometric dissipation) and, to a much lesser extent, pore fluid induced damping. Wave radiation damping and pore fluid induced damping are negligible for excitations below 1 Hz, and therefore hysteretic damping is dominant for the purposes of this study. The soil damping depends on the type of soil and the strain level.
Hydrodynamic damping: 0.07-0.23%. Results from wave radiation and viscous damping due to hydrodynamic drag. In the low frequency vibration of wind turbines the relative velocity of the substructure is low and therefore viscous damping, which is proportional to the square of the velocity is typically very low. The larger contribution results from wave radiation damping, which is proportional to the relative velocity.
Aerodynamic damping: in the fore-aft direction for an operational turbine 1-6%, for a parking turbine or in the crosswind direction 0.06-0.23%. Aerodynamic damping is the result of the relative velocity between the wind turbine structure and the surrounding air. Aerodynamic damping depends on the particular wind turbine, and is inherent in the popular Blade Element Momentum (BEM) theory for aeroelastic analysis of wind turbine rotors. The magnitude for a particular wind turbine also depends on the rotational speed of the turbine.
The total damping of the first mode of vibration is typically between 1% and 4% in side-to-side vibration, or for a parked, stopped or idling turbine. On the other hand, the total damping is between 2% and 8% for an operational wind turbine in the fore-aft direction.

Methodology and results
In this study, the OWT mechanical model used in [4,3,6] is extended to consider the effect of the stiffness of the substructure which becomes critical for deeper water, see Fig. 2. The rotornacelle assembly (RNA) is modelled as a lumped top head mass m RNA and the tower and the substructure are modelled as beams. The foundation is idealised as three springs (lateral K L , rocking or rotational K R , and cross coupling K LR stiffness). This methodology requires 14 input parameters which are given in Table 1 and the definitions are shown in Fig. 2. Arany et al. [6] demonstrated that the more sophisticated Timoshenko beam model, which takes into account mass moment of inertia and shear deformation in the equations of motion of the tower, provides no notable improvement over the simpler Euler-Bernoulli beam model in natural frequency estimation. It was concluded that the tower can be considered a slender beam and the Euler-Bernoulli beam is sufficiently accurate, therefore it was used to obtain the results presented in this paper.
Appendix B shows the derivation of the frequency equations of the two beam models. The natural frequency in both cases can be obtained numerically from the resulting transcendental equations. Approximate closed form expressions have been fitted to the results to fit the expression given by Eq. (2) which states that the first natural frequency (f 1 ) can be obtained by multiplying the fixed base frequency (cantilever beam frequency f FB ) by two factors C R and C L which account for the flexibility provided by the foundation. In other words, the foundation flexibility coefficients C R and C L are applied to the fixed base (cantilever beam) natural frequency to obtain the natural frequency taking into account the foundation compliance and soil-structure interaction. The formula for the first mode natural frequency is then given as:

Parameters and idealisation
The methodology uses 14 basic parameters for the wind turbine support structure and the site (foundation stiffness) and are summarised in Table 2. It is considered useful to explain the idealisation of the structure in a bit more detail. Fig. 2(b) shows the dimensions of the tapered and equivalent towers. The tapering is assumed to be linear with the tower diameter linearly decreasing from D b at the bottom to D t at the top of the tower. The average diameter is used as the equivalent constant diameter for the idealised tower D T , that is The wall thickness t T is assumed to be constant along the tower. The equivalent mass of the tower is calculated from the geometry and density of the tower material ρ T as m 0 Tower mass m T ton 8 Platform height above mudline L S m 9 Monopile diameter D P m 10 Monopile wall thickness t P m 11 Monopile Young's modulus E P GPa 12 Lateral stiffness of foundation K L GN/m 13 Cross stiffness of foundation K LR GN 14 Rocking stiffness of foundation K R GN m/rad Note, in some cases this mass may not necessarily be exactly the same as the actual mass of the tower m T . However, in cases when the average tower thickness or a range of tower wall thicknesses is not available, the thickness can be back-calculated from the tower mass The model is based on a beam with two different cross section. It is assumed that the cross section of the monopile is the first section of the beam from the mudline to the bottom of the tower, and the second section of the beam is the tower from the tower bottom to the rotor nacelle assembly. In this formulation the flexibilities of the grouted connection and the transition piece are neglected. The first section is of length L S , diameter D P , wall thickness t P and bending stiffness E P I P , while the second cross section is of length L T , diameter D T , wall thickness t T and bending stiffness E T I T . The platform height above mudline L S is defined as the distance from the mudline (seabed) to the bottom of the tower (towertransition piece connection).

Calculation procedure for simple natural frequency estimation
The natural frequency estimation process can be summarised in three steps: Step 1: Calculate fixed base natural frequency f FB of the variable cross section beam using the monopile flexibility coefficient C MP .
Step 3: Calculate and apply foundation flexibility coefficients C L and C R and apply them to obtain the first natural frequency on the flexible foundation, i.e. f 1 ¼ C L C R f FB The steps are described in more detail along with the methodology and calculation steps together with references for further reading.
Calculation procedure for Step 1: The fixed base natural frequency can be calculated by the simple natural frequency formula of the cantilever beam. The derivation can be found in Appendix A.
where for the tower and the RNA only, the equivalent stiffness of the beam k based on the tower bending stiffness E T I T and length L T , and the equivalent mass (considering the mass of the tower m 0 T and the mass of the rotor-nacelle assembly m RNA ) are given by: Substituting back into the formula, the following is given for the fixed base natural frequency.
In these expressions the second moment of area I T and the mass m 0 T of a constant diameter tower are considered. That is, if the average diameter of the tower is D T and the average wall thickness of the tower is t T then the equivalent second moment of area I T and the equivalent tower mass m 0 T are A similar expression is given specifically for wind turbine towers by van der Tempel and Molenaar [50].
In this paper the simple expression in Eq. (8) is used. The flexibility of the monopile is taken into account by a monopile flexibility coefficient C MP expressed in terms of two nondimensional numbers, bending stif fness ratio : where E T I T is the bending stiffness of the constant diameter tower, E P I P is the bending stiffness of the monopile, L S is the platform height above the seabed (distance from the mudline to the top of the transition piece) and L T is the height of the tower.
The fixed base natural frequency of the wind turbinesubstructure system is Fig . 3 shows the effect of the parameters χ and ψ on the fixed base natural frequency. To make these parameters more physically meaningful, one can say that increasing χ represents a stiffer monopile, while increasing ψ represents deeper water.
Calculation procedure for Step 2: In the second step, the foundation stiffness is non-dimensionalised using the tower's length L T and the equivalent bending stiffness of the tapered tower EI η . The non-dimensional lateral, cross-coupling and rotational stiffness values, η L ; η LR ; η R respectively, are necessary for this calculation These formulae use the equivalent bending stiffness for a tapered tower as given in [7][8][9]. The derivation is shown in Appendix C. Only the final formulae are presented here for brevity: where EI t is the bending stiffness at the top of the tower.
Calculation procedure for Step 3: The empirical foundation flexibility factors are applied to the natural frequency where a; b are empirical constants. The values of these constants were found to be a % 0:6 and b % 0:5 and were obtained by fitting closed form curves to the solutions obtained by solving the transcendental equations resulting from the analytical formulation given in Appendix B. Further details can be found in [6] and a simplified version of the problem has also been solved in [4,3]. A total of 16 wind turbines have been used to find best-fit values of a and b minimising the mean squared error between the numerical and closed form results. The curve fit reproduces the natural frequency for all 16 wind turbines with excellent accuracy ( o1%) for practical parameter combinations. The applicability of the formulae in Eq. (15) is limited to With these the natural frequency of the OWT on a flexible foundation can be calculated as given in Eq. (2).

Example step by step calculation: Blyth Offshore Wind Farm Vestas V66 2 MW
The input data for the calculations are shown in Table 3 (number III). Foundation stiffness is one of the most uncertain quantities and Appendix D discusses a method to obtain foundation stiffness in the absence of detailed site data.
Step 1: The parameters in Table 2 are used to calculate the fixed base natural frequency of the tower In this formulation D T is the average tower diameter The tower wall thickness is determined from the mass of the tower as and I T is the second moment of area of the equivalent constant diameter tower cross section The bending stiffness ratio and length ratio are calculated as where the monopile's second moment of area is The natural frequency of the monopile supported wind turbine on a fixed base is given as Step 2: The equivalent bending stiffness needed for the nondimensional stiffness parameters is calculated.
See details on the derivation of these expressions in Appendix B. Then the non-dimensional stiffness parameters can be obtained from Step 3: The foundation flexibility factors are calculated the natural frequency is then obtained from The measured frequency at Blyth wind farm is 0.488 Hz, and thus there is an excellent match with an error of 0.2%. The natural frequency on the monopile foundation is about 4.9% lower than it would be on a completely stiff fixed base. This difference is small because at the Blyth site the monopile is grouted into very stiff bedrock. The influence of the foundation flexibility is higher for most sites with values ranging from 4% to 15% for wind turbines analysed in this paper.

Application of the methodology to installed wind turbines and discussion
To demonstrate the general applicability of this method, the natural frequency of several wind turbines will be determined. The analysed offshore wind farms and turbines are presented in Table 1, as well as references to the sources of data. The input data of the turbines is summarised in  Table 4 and the results are given in Table 5 and compared with the measured data. The measured natural frequency values have been reported by various researchers and have been obtained by different methodologies of signal processing. Further details on the case studies can be found in the references given in Table 1. The results show that the empirical formula presented in this paper approximates the natural frequency within À3.2% to 3.1% of the measured frequency. This can be considered an excellent match given such a simplified methodology. The difference between the fixed base natural frequency and that for the flexible foundation using the present method is typically within the range 3-15%, which is considered typical for most offshore wind turbines. The flexibility (percentage reduction of the first natural frequency) introduced by the compliance of the foundation is reported for each wind turbine in Table 5.
As the foundation stiffness parameter ðK L ; K R Þ increases, the foundation flexibility coefficients ðC L ; C R Þ also increase and the natural frequency approaches the fixed base natural frequency f FB . Most of the contribution to the frequency change (4-15%) results from the rotational foundation flexibility coefficient C R and the lateral foundation stiffness coefficient C L has very limited influence (of less than 1% on the natural frequency) of all analysed turbines. Consequently, the most important foundation stiffness parameter is the rotational stiffness K R .
In Fig. 4 three turbines are used to illustrate the effect of the non-dimensional rotational stiffness η R on the natural frequency in terms of the rotational foundation flexibility coefficient C R . It should be noted that the curves flatten out for increasing rotational stiffness. The foundation designer should choose the foundation stiffness such that the OWT structure is in the flat part of this curve. This ensures that even if the foundation stiffness was estimated with significant errors or if the stiffness changes during the lifetime of the turbine, the natural frequency change is limited and does not affect the foundation's ability to meet the Fatigue Limit State (FLS) and Serviceability Limit State (SLS) criteria. On the dropping part of the curve, however, a slight increase or decrease of the stiffness might cause substantial change in the natural frequency, which may lead to increased fatigue damage and reduced service life.
It is important to emphasise that the foundation stiffnesses K L ; K LR ; K R and thus the non-dimensional foundation stiffnesses η L ; η LR ; η R are not independent parameters. Appendix D provides a brief discussion about the available methods for estimating the foundation stiffness, and simple expressions are provided for the simple cases of cohesive soils, cohesionless soils and rock. Following these formulae, the independent parameters of foundation stiffness are the geometry of the pile (diameter D P , wall thickness t P and embedded length L P as defined in Fig. 2(b)) and soil parameters (expressed using either the modulus of subgrade reaction k h , the coefficient of subgrade reaction n h or the effective shear modulus G Ã ). Furthermore, as discussed in Section 1, the damping in the first sway-bending mode also influences the measured natural frequency. The damping is lower for side-to-side direction and for a parked/idling turbine than in the fore-aft direction for an operational turbine due to the significant contribution of aerodynamic damping. The total damping of the system was estimated between 2% and 8% in Section 1, which amounts to a frequency change of about 0.02-0.32%.

Discussion on the parameters for foundation stiffness estimation
There are two groups of parameters necessary to determine foundation stiffness, categorised as pile data and soil data, as tabulated in Tables 6 and 7 for all the turbines listed in Table 1. The pile parameters are the pile diameter D P , the pile wall thickness t P , and the pile embedded length L P , as defined in Fig. 2(b). The diameter of monopiles is typically constant along the length. However, the pile wall thickness may change. In this paper the wall thickness was chosen as the pile thickness in the region below the mudline, because the top layers are considered more    important from the point of view of pile head deflection/rotation and pile stiffness. As can be seen based on Appendix D, the stiffness of the soil is described by different parameters for different soil types. For cohesive soils where the modulus of subgrade reaction is considered constant with depth below mudline, the critical parameter is the modulus of subgrade reaction k h . For cohesionless soils where k h increases linearly with depth, the key parameter is the coefficient of subgrade reaction n h , which describes the rate at which the modulus of subgrade reaction increases with depth. The parameter n h , however, is a simple linear function of the unit weight γ S and thus the density ρ S of the cohesionless soil using the formula of Terzaghi [49] as given in Eq.
(D5) (in Appendix D). In case of rocks the important parameters are the shear modulus G S and the Poisson's ratio ν S , from which the effective shear modulus of Eq. (D7) can be calculated. As mentioned in Section 3, the foundation stiffness parameters K L ; K R and K LR are not independent. It is apparent from Table D1 (in Appendix D) that the independent parameters required to determine the foundation flexibility factor F FF ¼ C R C L are different for each soil type and pile slenderness combination. The simplest cases are for cohesionless soils where only two independent parameters are necessary. For slender piles these are the density of soil ρ S and the bending stiffness of the monopile E P I P , and for rigid piles they are the soil density ρ S and the pile embedded length L P .
These two cases are shown in Fig. 5. There are more independent parameters for clays and rock, and those cases therefore do not allow for easy visualisation.

Conclusions
A simple methodology has been presented to calculate the first natural frequency of an offshore wind turbine founded on a monopile. The methodology is based on multiplying the fixed base cantilever beam natural frequency by two foundation flexibility factors to include the compliance of the foundation as well as the flexibility of the substructure. The monopile foundation is modelled by three coupled springs (lateral, rotational and cross coupling stiffness). The closed form formulae presented in this paper were obtained from fitting curves to the natural frequency results obtained from numerically solving the transcendental equations. Several conclusions can be drawn from this study: It is shown through the study of 10 wind turbines that the natural frequency can be predicted with accurately (the error range is 73.5%).
For the wind turbines considered, foundation flexibility reduced the fixed base natural frequency typically by 4-15%.
The foundation flexibility factor is very sensitive to the rotational stiffness K R of the monopile. On the other hand, the natural frequency change from the fixed base frequency due to lateral stiffness of the monopile ( K L ) is limited to about 1%.
The simple framework for the calculation of the foundation stiffness following the work of Poulos and Davis [42] Randolph [43] and Carter and Kulhawy [20] used in this study was found to give reasonable results for the 10 wind turbines studied where three different soil types (cohesive/clay, cohesionless/ sand and rock) and two limit case approximations of slender/ infinitely long and rigid/infinitely stiff piles were analysed. It was found that most real monopiles fall between the limiting cases of infinitely long and rigid piles.
The effect of water depth and monopile stiffness were also analysed and it is found that the methodology can correctly represent the flexibility of the monopile and therefore the approximation is applicable both in shallow and deep water.
While the final design verification can be carried out using detailed finite element investigations, the method presented in the paper can be used to estimate the natural frequency using limited information about the wind turbine and the site and may be a useful tool for initial analyses and conceptual design.

Acknowledgement
The first author is supported by an EPSRC DTP studentship via the University of Bristol.
Appendix A. -Natural frequency of a cantilever beam with variable cross section The motion of the cantilever beam can be described as a single degree of freedom mass-spring system. The free vibration of this system is given by Assuming harmonic vibration, the following equation can be obtained: and from the circular frequency the Hertz frequency is easily obtained The flexibility of the substructure expresses the dependence of the natural frequency on the water depth, that is, the flexibility of the monopile above the mudline and that of the transition piece. For the sake of simplicity, the model used in this paper assumes that the monopile's bending stiffness continues throughout the water depth and up to the root of the tower. For clarity, see Fig. A1.
Bending stif fness ratio : χ ¼ Platf orm=tower length ratio : ψ ¼ Castigliano's second theorem for a linearly elastic 1 DoF structure can be written as where U is the strain energy, q is the generalised displacement, and Q is the generalised force. For the particular problem, the theorem can be used to calculate the top head deflection (the total deflection at the hub wð0Þ) due to a horizontal force F acting at the hub as where the moment distribution along the structure caused by the horizontal force F is given as The stiffness of the 1DoF system is then given as From this the natural frequency is calculated as

Dimensionless group Formula Dimensionless group Formula
Non-dimensional lateral stiffness K L ; K R ; K LR are the lateral, rotational and cross stiffness of the foundation, respectively; EIη is the equivalent bending stiffness of the tapered tower; L T is the hub height above the bottom of the tower; P Ã is the modified axial force, m RNA is the mass of the rotor-nacelle assembly; m T is the mass of the tower; J is the rotary inertia of the top mass; μ is the equivalent mass per unit length of the tower.
n The rotary inertia is taken to be zero for all wind turbines considered as information is not available in the referenced literature.
where m is the generalised mass of the 1DoF system. The stiffness can be written as where k T is the stiffness of the tower without the substructure, given as The fixed-base natural frequency of the tower-substructure system (excluding foundation stiffness) is given using the bending stiffness ratio χ and the platform/tower length ratio ψ as Appendix B. The Euler-Bernoulli beam model equation The equation of motion using the Euler-Bernoulli beam model for a beam with an axial force is where EIðzÞ is the bending stiffness distribution along the axial coordinate z, μðzÞ is the distribution of mass per unit length, P Ã is the axial force acting on the beam due to the top head mass and the self-weight of the tower, p z; t ð Þ is the excitation of the beam, wðz; tÞ is the deflection profile.
Using constant equivalent values for the axial force, bending stiffness and mass per length, and considering free harmonic vibration of the beam with separation of variables w z; t ð Þ¼W z ð Þ•e iωt , the equation can be reduced to the following using the non-dimensional parameters of Table A1 and the dimensionless axial coordinate ξ ¼ z=L: is the non-dimensional circular frequency. Using the non-dimensional numbers as defined above, the boundary conditions can be written for the bottom of the tower ðξ ¼ 0Þ: and the top of the tower (ξ ¼1): The parameters used in the boundary conditions are defined in Table A1. The characteristic equation for the equation of motion can be written as The four solutions are then with which the solution is in the form which can be transformed using Euler's identity to with : (B11) Substituting this form of the solution into the boundary conditions, one obtains four equations, written in matrix form as with and Looking for nontrivial solutions of this equation one obtains from which one can obtain the non-dimensional circular frequency Ω, and from that the natural frequency using The equation that has to be solved is transcendental and therefore solutions can only be obtained numerically.

Appendix C. -Tower idealisation
The towers of offshore wind turbines are tapered towers with diameters decreasing from the bottom to the top. Typically, the wall thickness of the tower also decreases with height. However, some small and medium sized turbines have constant wall thickness. The formulation presented in this paper replaces this tower shape with an equivalent constant diameter, constant wall thickness tower. The average tower diameter is used in combination with an equivalent tower wall thickness t T . Note that if the average wall thickness is determined from a range ðB14Þ of wall thicknesses of the tower or from the mean of the top and bottom wall thicknesses, then the tower mass, as calculated from the idealised tower geometry, may not be the same as the actual tower mass m T . In the case that information about tower wall thickness is not available, the equivalent thickness can be chosen such that the actual tower mass is maintained, that is, The non-dimensional stiffness parameters are normalised with the length L T and the bending stiffness E T I T of the tower. When calculating these non-dimensional stiffness parameters, the equivalent bending stiffness is calculated such that the deflection at the tower top due to a force acting perpendicular to the tower at the tower top is the same for the equivalent constant diameter tower as that of the tapered tower. The derivation is given here.
The diameter varies along the structure as with z ¼ 0 at the top of the tower and positive downwards. The second moment of area is then given as Moment curvature relation is written as where F is the horizontal force at the hub. One can write and from that w is obtained via integration The boundary conditions are used to calculate the constants.
The deflection at the end of the column is obtained by substituting z ¼ 0 into the equation of deflection The stiffness is then given as Verification: for a cantilever beam of constant diameter lim q-1 In the paper the following notification is used -Guidance on the calculation of foundation stiffness The stiffness of the foundation is used as input in the calculations shown in this paper. However, determining the stiffness is the most challenging task. In the absence of careful (very expensive) and detailed site measured data, existing formulations available in the literature may be used to estimate the stiffness of the foundation theoretically. For the sake of brevity only a brief summary of some common methods are presented here.
For static analysis, two methods are given in [42], one methodology based on a subgrade analysis approach partially based on Broms [16] and Broms [17] and another based on an elastic continuum approach. Results are provided for both cohesive and cohesionless soils, with simple analytical results available for slender (infinitely  ρ θ ! long) and rigid (infinitely stiff) piles. Another static stiffness approach is recommended in Eurocode 8 Part 5 [26] based on Gazetas [27], developed for seismic analysis of slender piles. Dynamic analysis using a subgrade reaction approach of a beam on an elastic foundation was developed in terms of p-y curves [37,39,41,44,5], which is the currently accepted design procedure suggested by e.g. the DNV code for offshore wind turbines [25]. This methodology is based on Winkler's approach [53] and utilises nonlinear springs to model the elastic foundation. Analytical solutions for dynamic analysis have been published by e.g. Nogami and Novak [38,40]. Green's function based models have also been developed for dynamic stiffness of piles, e.g. [29,31]. Recent work regarding dynamic stiffness and damping of monopiles of offshore wind turbines include Shadlou and Bhattacharya [45], Zania [56] and Damgaard et al. [23]. For seismic analysis, a three-step solution was suggested by Kausel et al. [30], including simplified approximate expressions that make the method more attractive than finite element analysis. Finite element based methods for analysis of piles were developed by e.g. Blaney et al. [15].
In this appendix, however, the simple formulae in [42,43,20] are used. These formulations neglect the frequency dependence of foundation stiffness, which can be justified for dynamic loading of offshore wind turbines because the frequencies of excitation are so low. For seismic analysis, however, this frequency dependence should be taken into account. These approximations provide quick and easy solutions for foundation stiffness, achieving sufficiently accurate results for the natural frequency as shown in the paper.
The two main parameters to decide the calculation methodology for foundation stiffness are: (1) Soil condition and ground profile at the site and (2) Pile slenderness/rigidity. The ground type/ soil determines the main soil parameter for estimating soil stiffness, while determining whether the pile can be considered slender or rigid. This allows for simplification of the foundation stiffness estimation and for the use of closed form solutions instead of graphs or numerical methods.
Ground profile and soil conditions at the site There are three main categories to consider from the point of view of the current analysis.
(1) Cohesive soils: the horizontal modulus of subgrade reaction k h is considered constant with depth below the mudline. This is typically used for over-consolidated clayey soils which is often encountered in offshore conditions, see Bhattacharya et al., [14]) for typical North Sea soils. For normally consolidated cohesive soils, the subgrade reaction increases linearly with depth. (2) Cohesionless soils: the horizontal modulus of subgrade reaction k h is assumed to increase with the square root of the depth below the mudline. This can be typically used for loose to medium dense sand and gravels. (3) Bedrock: The foundation stiffness is determined by the shear modulus of the soil. This is typically used for weathered bedrock and very dense sand. (4) Complex layered soils: For real sites the soil is obviously not always as simple as the above categories and different layers are often observed. One of the above categories may be chosen for such complex sites, bearing in mind that from the point of view of pile head deflection/rotation and pile stiffness the upper layers of the soil are of higher importance, and the upper layers should be weighted accordingly. The soil types at the sites of all the wind turbines considered in this paper are listed in Table 1 in Section 1. Alternatively, the stiffness can be more accurately obtained from p-y curves.
The constant horizontal modulus of subgrade reaction k h for cohesive soils can be determined using the expression of Vesic [52], using the parameters from Tables 6 and 7: where I P is the second moment of area of the pile cross section E S is the elastic modulus of the soil, and ν S is the soil's Poisson's ratio.
Broms [17] provides another expression for clays based on the secant Young's modulus E 50 which can be used in combination with the formula of Skempton [47] or the more conservative formula of Davisson [24] as quoted in [42], in terms of the undrained shear strength c u Skempton : E 50 ¼ ð50 À 200Þ Á c u and k h ¼ ð80 À 320Þ Á c u D P The linearly varying modulus of subgrade reaction of cohesionless soils can be written following [42] as: where n h is the coefficient of subgrade reaction, calculated for sand after [49] as where γ sand is the specific weight of sand and A ¼ 100 À300 for loose sand, A ¼ 300 À 1000 for medium sand, and A ¼ 1000 À 2000 for dense sand. The modulus of subgrade reaction k h for cohesive soils and the coefficient of subgrade reaction n h for cohesionless soils are used to calculate the foundation stiffness, as given in Table D1. In Table D1 ρ and θ are mudline pile head deflection and rotation, respectively, F and M are horizontal force and overturning moment at the mudline.

Slenderness/rigidity of the pile
Simple closed form solutions are readily available for the simplified cases obtained by assuming either a slender pile or a rigid pile.
(1) Slender pile: The monopile is idealised as 'slender' or 'infinitely long' assuming that the pile flexibly deflects and that the pile fails first by yielding through a plastic hinge (as opposed to failure of the soil). (2) Rigid pile: The monopile is idealised as 'rigid' or 'infinitely stiff' assuming that the pile undergoes rigid body rotation (the soil fails first). However significant bending moment may be generated in the pile Two main methods are presented here to determine whether a pile can be slender or rigid. The first one is given in [42] and is based on the modulus of subgrade reaction k h and the bending stiffness of the pile E P I P . The slenderness parameter is calculated as The pile is considered slender or infinitely long if βL P 42:5, and considered rigid if βL P o1:5: The second method is that of Randolph [43] and Carter and Kulhawy [20], which is based on the equivalent Young's modulus E e of the pile, and the effective shear strength of the soil G Ã : where G S is the shear modulus of the soil. The pile is considered rigid; if L P D P r 0:05 This methodology is used for rocks by Carter and Kulhawy [20]. The formulae summarised in Table D1 were used to approximate the foundation stiffness in the present study and were found to provide good approximations in terms of natural frequency. The difference between the rigid and slender pile approximations is typically 1-6%.