1.Introduction

Wind,solar power and geothermal heat are representative of clean and renewable energy which have the potential to become alternatives to current supplement of fossil fuel sources of energy in the future.While these alternative energy sources have their advantages and drawbacks,wind energy is widely accepted as the cheapest and most economically available one based on current technology.Today,wind energy has proven to be a valuable feature for large-scale future investment in the energy industries worldwide,and many countries install their proper wind turbine generators(WTGs),mainly on land.As offshore wind turbines(OWTs)have gained their popularity,many WTG manufacturers believe that offshore wind energy will play an increasingly important role in the future development.This is supported by the fact that all principal wind turbine manufacturers currently are spending huge amount of money and effort on developing larger offshore WTGs for deeper waters where wind speed is generally higher and steadier,resulting in an increase in energy output.

Although there are many OWTs support options which may range from gravity foundations(for shallow depths of 0-15 m)to floating foundations(for very deep waters of 60-200 m)(Achmus et al.,2009;Lombardi et al.,2013;Damgaard et al.,2015;Abed et al.,2016),most OWTs are supported on monopile foundations,as they are simple structures which are easy and convenient to construct.The accumulated experience from limited monitored data from OWTs over the last 15 years showed that the available design procedures(mostly contained in the API(API and ISO,2011)and DNV(DNV-OS-J101,2004)regulation codes suffer limitations.

The existing methods were established/calibrated by testing small-diameter piles used for supporting offshore platforms in gas and oil industry,often with design criteria and loading conditions which are different from those encountered in an OWT.The inappropriateness of these methods comes from the fact that:

(1)The continuum(soil)is replaced by a series of uncoupled springs.However,reliable results necessitate a rigorous method which can properly account for the true deformation mechanism of soil-monopile interaction.

(2)As they rotate freely,monopiles supporting OWT energy converters undergo severe degradation in the upper soil layer resulting from cyclic loading,whereas offshore jacket piles are significantly restrained against rotation at their heads.

where νsis the Poisson’s ratio,and E is the Young’s modulus.

(4)The API model is calibrated in response to a small number of cycles for offshore fixed platform applications.However,an OWT over its lifetime of 20-25 years may undergo 107-108 cycles of loading.

Due to these complex issues,appropriate determination of the dynamic characteristics of these extremely complex structures through their monopiles head stiffnesses is continuing to challenge designers,as the foundation of an OWT behavior is still not well understood,and also not introduced in the current design guidelines.

The ultimate stress difference(σ1- σ3)ultis defined in terms of the actual failure stress difference by another parameter called‘failure ratio’Rfwhich is given by

Otsmane and Amar Bouzid(2018)formulated a nonlinear pseudo three-dimensional(3D)computation method,combining the FEM and the finite difference method(FDM).They wrote a Fortran computer code called NAMPULAL(nonlinear analysis of monopiles under lateral and axial loadings)to study monopiles under axial,lateral and moment loadings in a medium characterized by the hyperbolic model for representing the stressstrain relationships.In this paper,we attempt to apply NAMPULAL to examining the lateral behavior of monopiles supporting OWTs chosen from five different offshore wind farms in Europe.These offshore wind farms include Lely A2(UK),Irene Vorrink(Netherlands),Kentish Flats(UK),Walney 1(UK)and Noth Hoyle(UK).

To accurately estimate the natural frequency of the OWT structure(tower+substructure)which is a function of monopilesub soil interaction,the monopiles head movements(displacements and rotations)and consequently,the lateral stiffness KL,the rotational stiffness KRand the cross-coupling stiffness KLRare obtained and substituted in the analytical expression of natural frequency for comparison.In general,the results of comparison between the computed and measured natural frequencies showed a good agreement.

2.Natural frequency and modal analysis

OWTs are dynamically sensitive structures,in which the dynamic soil-structure interaction is a pivotal aspect of their design process and consequently,they require accurate soil stiffness estimation in order to ensure that the design frequency matches the actual operational frequency when the wind turbines are constructed.

The natural frequency of the hub-tower-foundation system is the key feature on which the response of an OWT to wind and wave loads depends.This is due to the dynamic nature of the loads on the wind turbine structure and the slenderness of the system.Through determination of the natural frequency,designer can assess the strains produced by loading cycles,through which the fatigue failure of the structure can be ascertained.Therefore,an accurate estimation of this parameter is essential to assess the working life of a wind turbine.

Unlike most large-scale civil Engineering structures,wind turbines are subjected to millions of periodic excitation cycles during their operating life.The rotor spinning at a given velocity induces mass imbalances(gyroscopic effect),causing a frequency known as 1P.In addition to this,the effect of a standard turbine having n blades induces a further excitation due to the blades passing the tower.The frequency of this shadowing effect is nP,where n=3 in most cases.

The modern installed wind turbines are characterized by a range of different velocities in which their rotors are operating.This results in two ranges of operating frequencies around 1P and 3P.In order to avoid resonance,the natural frequency of the tower cannot be in any of these two ranges and must be far from 1P and 3P.

The OWT design can be performed in such a way that the first eigenfrequency lies within three possible ranges:soft-soft,soft stiff and stiff-stiff as shown in Fig. 1.

(1)Soft-soft range:the natural frequency is less than the lower bound of 1P.This implies that the structure is too flexible,and moreover,this is a range where the frequency of waves may lie,therefore leading to resonance.

(2)Stiff-stiff range:this is a range where the tower frequency is higher than the upper bound of blade passing frequency(3P).This range is economically unfeasible as it leads to a too rigid(heavy and expensive)structure,making it inappropriate for design.

(3)Soft-stiff range:in this interval,the natural frequency lies between 1P and 3P.This range is the optimum range for the best possible design.

对9度地震烈度下钢管塔沿塔身高度的位移响应进行研究,分别在4类场地下沿0°方向输入地震波,分别计算5.17,15.5,20.5,26.6,39.0,50.5,71.5,92.0,95.5,101 m塔身高度下的位移响应,得到的钢管塔沿塔身高度的位移响应曲线如图4所示.

where a and b are the empirical constants,which have been obtained by fitting closed-form curves;EIηis the equivalent tower bending stiffness;andηis the soil-foundation interaction coefficient depending on tower’s bending stiffness.The applicability of Eqs.(4)and(5)is conditioned by

In order to satisfy these requirements and to keep the natural frequency of the whole structure in the adequate margin of the soft-stiff range,thus avoiding resonance,a joint effort between foundation designers and turbine manufacturers is performed.Foundation designers need careful site investigations to obtain reliable soil data in order to correctly assess the foundation stiffness.

2.1.Appropriate OWT modeling for dynamic analysis

The natural frequency of a wind turbine is highly dependent on the material properties used in its construction,and is significantly affected by the stiffness of the soil surrounding the monopile.Assessment of foundation stiffness is the key to obtain reliable estimate of system frequency.

In the computation of eigenfrequency f1,most researchers in the past tried to model this complex system according principally to two concepts(Prendergast et al.,2015;Yi et al.,2015).In the first one,Yi et al.(2015)simply considered the soil as a medium having an infinite stiffness.In this regard,Vught(2000)used a model in

Fig. 1.Forcing frequencies against power spectral density for a three-bladed wind turbine(Hz).

which the wind turbine is considered as an inverted pendulum having a flexural rigidity EI,a tower mass per meter mTand a top mass mt.Expression for the first eigenfrequency is given by

where LTis the tower height.The natural frequency expressed by Eq.(1)is based on a uniform tower cross-section.A slight different expression has been proposed by Blevins(2001):

It seems from the first sight that these equations are based on a simple model which ignores the fact that the tower is tubular and conical in shape and generally does not have a constant wall thickness.Additionally,Eqs.(1)and(2)depend only on the tower geometrical and mechanical properties without taking into account the OWT foundation characteristics.This physically does not make sense,as the monopile head movements that occur as a result of the applied loading on the tower can lead to a finite stiffness of the monopile-subsoil system.This obviously has an influence on the value of the first natural frequency and shows that the service limit computations based on Eqs.(1)and(2)are inaccurate.

where γsteelis the steel density usually taken as 7860 kg/m3.

Bearing in mind that the dynamic analysis of the whole system composed of tower-monopile-soil is hard to perform,Prendergast et al.(2015)tried to find a natural frequency expression considering soil stiffness(Zaaijer,2006;Yu et al.,2014;Prendergast et al.,2015).Alternatively,they replaced the subsoil-monopile system by a set of springs through which the tower is connected to the subsoil,as shown in Fig. 2.This figure illustrates a mechanical model,in which the subsoil-monopile interaction is represented by four springs,i.e.a lateral,a rocking,a cross-coupling and a vertical spring whose stiffnesses are KL,KR,KLRand KV,respectively.Most researchers disregarded the axial vibrations since the wind turbines are very stiff vertically.

兔抗人Livin多克隆抗体（武汉博士德生物工程公司,货号：BA1743）、兔抗人STAT-5多克隆抗体（武汉博士德生物工程公司，货号：BA1411）、SP免疫组化染色即用型快速免疫组化MaxVisionTM试剂盒（福州迈新生物技术有限公司，货号：KIT-5001）；DAB显色剂（武汉博士德生物工程公司,货号：AR1021）。

The stiffness of these springs which represent the subsoilmonopile interaction may be estimated from the monopile head load-deformation curves,provided that these curves are obtained by means of a rigorous modeling method,such as FEM.

On the basis of a numerical solution of transcendental frequency equation,Adhikari and Bhattacharya(2011,2012)proposed an exact approach where only lateral and rotational stiffnesses have been included.Furthermore,in order to improve the first natural frequency equation,Arany et al.(2014,2015)derived expressions of natural frequency of OWTs on three-spring flexible foundations by means of two beam models:Bernoulli-Euler and Timoshenko.The natural frequencies in both cases have been obtained numerically from the resulting transcendental equations.They proposed a closed-form expression containing,in addition to KLand KR,the cross-coupling stiffness KLRof the monopile.Their equation for the natural frequency is

炉衬采用砖砌炉衬和耐火纤维炉衬组合形式，砖砌炉衬由耐火层和隔热层组成，根据不同的使用温度选用刚玉、莫来石、氧化铝砖等，炉底灰缝不大于2 mm。由于耐火纤维制品的热导率低、密度小，比热容与普通耐火砖接近，因此与砖砌炉衬相比，其散热损失和蓄热损失大大减小，因此炉体侧墙和炉顶可采用耐火纤维炉衬[4]。

计量方程通过所有制（Own）、行业（Industry）、省份（Province）和时间（Year）等虚拟变量，处理其它固定效应因素。为了分析中间品进口的国际技术溢出对我国整体技术创新的影响，设定四个模型：

where fFBis called fixed base frequency which can be either Eq.(1)or(2).The factorsCRand CLaccount for the stiffness provided by the monopile,and are functions of tower’s geometrical properties.Their analytical expressions are given by

The system stiffness must be such that the natural frequency of the wind turbine does not lie within the rotor frequency excitation bands,as this may induce resonance which could lower the design life significantly.

Fig. 2.The OWT model used:(a)Principal components and(b)Model considering soil stiffness through monopile head springs.

Although Eq.(3)is mathematically attractive as it contains three simple factors,finding its constituting parameters is not an easy task.This equation is very important in the sense that the two first parameters account for the interaction between the soil and the monopile through the spring stiffnesses.Moreover,they incorporate terms related to the equivalent tower stiffness which should be evaluated properly.

The OWTs,especially those installed in deep waters,differ from onshore wind turbines.This difference comes from the fact that these OWTs are considered as slender and heavy structures which require in general three elements in their construction to bear the heavy masses.These include the tower,the monopile overhang and the transition piece which assembles the first two elements(Fig. 2a).

In Fig. 3d,the length of the tower LTaccounts for the distance from the rotor nacelle assembly to the top of the transition piece.The tower has a varying bending stiffness EIT,a thickness tTand top and bottom diameters which are respectively Dtand Db.Neglecting the f l exibilities of the grouted connection,the monopile overhang and transition piece welded together are assumed to constitute one element,called the substructure.The latter has a length of Lswhich is defined as the distance from the mudline(seabed)to the bottom of the tower.The diameter Dsand the thickness tsof the substructure are assumed to have the same values as those of the monopile on which the substructure is founded.Consequently,the bending stiffness of the substructure is the same as that of the monopile EIp.

由表2可以看出，对于G分量图像：对面积在190～250之间的鸡蛋，蛋黄指数大于0.26则判定为双黄鸡蛋；对面积小于190的鸡蛋，蛋黄指数大于0.38则判定为双黄鸡蛋。

As most towers are tapered and tubular,the value of the bending stiffness is used to evaluate the fixed base frequency and consequently the natural frequency is difficult to obtain,although the variation law of EITwith the increasing tower cross-section is easy to establish.In this context,Bhattacharya(2011)studied a tower as a tapered cantilever beam subjected to a concentrated force p applied at its free end(Fig. 3a).Then,by means of beam theory,he computed a parameter fp(m)(termed here as‘tower stiffness coefficient’)as the ratio of the top displacement of a tower having a constant cross-sectionto that of a tower having a linearly varying cross-section This parameter has been determined as

where m is the ratio of bottom diameter to top diameter(Fig. 3d):

However,it is more likely to consider the tower as a tapered beam subjected to an upward tapered load along the whole length(Fig. 3b).Integrating the beam lateral deflection equation twice and setting the suitable boundary conditions,the tower stiffness coefficient corresponding to this load has been obtained as

Fig. 3.OWT model used to evaluate the tower mass and bending stiffness:(a)Tower subjected to p;(b)Tower subjected to downward tapered load;(c)Tower different masses;and(d)Geometrical properties of the tower.msis the mass of the substructure(transition piece+monopile)and Dpis the monopile diameter.

Fig. 4 illustrates the evolution of fp(m)and fqtriang(m)with the increasing values of m.It is clearly seen that for the interval where the ratio m varies from 1 to 2.5,both curves yield the same values.This corresponds to the majority of practical applications.Outside this range,a neat discrepancy is observed and we suggest using the average value if it occurs to find a ratio m greater than 2.5.

Fig. 4.Evolution of tower stiffness coefficient with tower diameters ratio.

The bending stiffness EITof the tower may be evaluated as

读后感：1 8岁的伊福琳·安娜·克露顿说:“信中谦和甚至带着抱歉的口吻，不仅没有刺激我，反而让我觉得敢于向哈佛提出申请，本身就是一件非常值得自豪的事。”

The equivalent bending stiffness of the whole structure(support structure)is given by

1) 在生物量方面,在黏土中,从移栽到DAP30,3个处理间地上部生物量没有显著差异,DAP30到DAP60,甘薯地上部茎叶生物量随着压力的增大而增加,而到了生育后期(DAP118),T2处理高于T3,说明适宜的外源压力有助于甘薯茎叶在中后期的生长．在黏土和砂土中，各生育期甘薯地下部生物量均呈现出一致的趋势．

Eqs.(1)and(2)for an onshore or offshore wind turbine,where the tower is resting directly on the soil,should be altered to properly find a fixed base natural frequency for an OWT having different masses and different bending stiffnesses,as shown in Fig. 3c.If weadopt Eq.(2),this would have the following form for an OWT composed of both tower and substructure parts:

2)冷水年,菲律宾反气旋的存在使得副高西伸加强,显著加强了其西侧暖湿气流向江南地区输送,高层辐散抽吸作用导致江南地区对流上升运动增强,暖水年相反,表明冷(暖)水年对应着江南雨季降水偏多(少)。

Since the support structure composed of tower and substructure is in contact with soil,the empirical soil-foundation interaction coefficients given in Eq.(6)should be corrected in order to properly take into account the true support structure length and its effective bending stiffness.Thus these expressions are given by

The soil-foundation interaction coefficients given by Eqs.(4)and(5)are evaluated on the basis of Eq.(19).

图3表示测试模式下对低电压SRAM单元进行写0操作来检测PMOS管上的稳定性故障。在测试模式下将低电压SRAM的位线WBLB设置为浮动0,并执行写操作,由于浮动0与QB处的目标值相反,所以成功写入数据变得更加困难,首先QB处受到浮动“0”攻击,若稳定性故障落在PMOS的源极或漏极上,削弱了交叉耦合反相器的上拉能力,则低电压SRAM将无法正确地写入数据,从而检测到PMOS上的稳定性缺陷。

2.2.Procedures to estimate the monopile head stiffnesses

In dealing with monopiles supporting OWTs,design engineers need to compute KL,KRand KLR.Two ways are often considered to compute these stiffnesses.The fi rst way is to model the monopile using the Winkler concept.This procedure which is alternatively called p-y approach is assumed to be sufficiently accurate for monopile diameter Dp?2 m,as p-y curves have been established for small-diameter and slender piles in offshore gas and oil industry.However,several investigations indicated that the pile deflections of large-diameter monopiles are underestimated for service loads and overestimated for small operational loads,which has been confirmed in a separate work(Otsmane and Amar Bouzid,2018).The second way is to directly employ values of KL,KRand KLR given in the existing standards(EC8 for example,where pile-head stiffness of flexible pile is provided).Although the expressions containing these coefficients have been determined for various soil profiles(three profiles in most cases:constant soil stiffness,linear variation of soil stiffness with depth,and variation of soil stiffness with square root of depth),they encompass monopile-soil Young’s modulus ratio Ep/Es.However,recent research(Higgins et al.,2013;Abed et al.,2016;Aissa et al.,2017)confirmed that the lateral behavior of large-diameter monopile fundamentally depends on monopile slenderness Lp/Dprather than monopile-soil relative stiffness Ep/Es.

A different class of researchers suggested that the stiffness coefficients could be obtained from the elastic behavior of soilmonopile system under lateral loading where both soil and monopile are elastic.Indeed,some researchers(e.g.Carter and Kulhawy,1992;Higgins et al.,2013;Abed et al.,2016;Aissa et al.,2017)performed parametric studies using FEM,in which the ratios Ep/Eswere varied and load-deflection curves were drawn.These studies confirmed that the short monopile head stiffness for monopiles embedded in elastic media depends only on the monopile slenderness rather than the stiffness ratio.Results from the aforementioned references are given in Table 1 for homogeneous soils.The corresponding values for Gibson soil are given in Table 2,and Table 3 provides values for soils where stiffness varies with square root of depth.For simplicity,results presented in these tables are restricted to Poisson’s ratio equal to 0.4.

Using the values of CRand CL,the natural frequency was found very close to unity,i.e.close to fixed base.Precisely,it may be argued that they do not reflect the soil actual stiffness and they do not bring an actual alteration to the natural frequency which remains very close to the fixed base frequency.

As an alternative procedure,it is appropriate to evaluate the stiffness coefficients on the basis of initial stiffness(tangential values at the origin)of the monopile head-deformations curves,resulting from the study of the nonlinear behavior of the subsoil in which the monopile is embedded.Assuming that the monopile head movements and applied efforts are expressed in function of fl exibility coefficients,this can be given in matrix form as where H and M are the shear force and overturning moment applied at the monopile head,respectively;uLandθRare the lateral displacement and rotation of the monopile head,respectively;and IL,IR,and ILRare the lateral,rotational and cross-coupling flexibility coefficients,respectively.

铁路建设是集企业利益和社会效益于一体的系统工程，线路方案、车站的选择涉及企业、社会各方的利益。由于很难兼顾各方的全部意见，选择推荐合理、可实施的方案，总体除需具备很高的政策水平和很强的业务能力外，还需具备出众的组织沟通协调能力，特别是要与业主保持良好的沟通，真心实意当好业主的参谋。

As the aim is to find the stiffness coefficients,it is easy to reverse the matrix in Eq.(20)to obtain:

The stiffness coefficients are related to flexibility ones by the following terms:

Table 1 Short monopile stiffness coefficients proposed by different researchers in homogeneous soils.

Note:EsDis the soil Young’s modulus at one monopile diameter depth,and Lpis the monopile length.

Source KL/(EsDDp) KLR/（EsDD2p） KR/（EsDD3p）Carter and Kulhawy(1992)1.884(Lp/Dp)0.627 -1.048(Lp/Dp)1.483 1.91(Lp/Dp)2.049 Higgins et al.(2013)2.426(Lp/Dp)0.71 -1.44(Lp/Dp)1.67 1.789(Lp/Dp)2.459 Aissa et al.(2017) 2.756(Lp/Dp)0.668 -1.595(Lp/Dp)1.6361.731(Lp/Dp)2.495

Table 2 Stiffness coefficients for short monopiles proposed for Gibson soils.

Source KL/(EsDDp) KLR/（EsDD2p） KR/（EsDD3p）Higgins et al.(2013)0.929(Lp/Dp)2.041 -0.633(Lp/Dp)3.0610.672(Lp/Dp)3.941 Abed et al.(2016) 1.708(Lp/Dp)1.661 -1.233(Lp/Dp)2.6551.153(Lp/Dp)3.605

Table 3 Stiffness coefficients for short monopiles proposed for soils whose stiffness increases with square root of depth.

Source KL/(EsDDp) KLR/（EsDD2p） KR/（EsDD3p）Abed et al.(2016) 2.841(Lp/Dp)0.977 -2.933(Lp/Dp)1.767 3.894(Lp/Dp)2.562

Determination of the values KL,KRand KLRis not straightforward in the FEM analyses controlled by forces.The flexibility coefficients should be determined first,and then be reversed to obtain the stiffness coefficients of Eq.(22).To do so,an arbitrary pure horizontal load(H≠0 and M=0)is applied at the monopile head at the mudline level,and a plot depicting the increasing values of H versus the corresponding values of monopile head displacements u is illustrated(Fig. 5a).The parameter 1/ILis then obtained bysimply computing the slope of the resulting curve at the origin.The parameter 1/ILRis computed from the curve giving the variation of H in function of rotationθissued from the same analysis(Fig. 5c).

As the rocking flexibility coefficient needs a pure bending,the monopile-soil system is analyzed under an overturning moment(M≠0 and H=0)applied at the topof the pile at the mudline level.From the curve portraying the M increasing values against the obtained rotationsθ,the reciprocal of flexibility coefficient 1/IRis evaluated by simply computing the slope of curve tangent at the origin(Fig. 5b).

This procedure is followed in this paper,when OWT monopiles of the different wind farms are considered in the next sections.

3.Numerical methodology:the computer program NAMPULAL

A pseudo 3D FEM model has been performed to study soilstructure interaction problems in nonlinear media.This procedure,called nonlinear finite element vertical slices model(NFEVSM),involves the combination of the FEM and the FDM for capturing the behavior of the embedded structure and its surrounding soil being considered to obey the hyperbolic model as proposed by Duncan and Chang(1970).The 3D soil-structure problem plotted in Fig. 6 shows a soil-structure interaction problem example(Fig. 6a)and the vertical slices model where different slices are acted upon by external forces and body forces(Fig. 6b).

The stress and deformation analyses in each slice are conducted by the conventional FEM,using two-dimensional(2D)finite elements.According to the standard formulation in the displacement based FEM,the element stiffness matrix in slice i can be written as

where B and BTare the strain field-nodal displacement matrix and its transpose,respectively;N and NTare the shape function matrix and its transpose,respectively;piis the external force vector to which the slice i is subjected;and biis the body force vector which has the following compact form:

where

Fig. 5.Monopile head load-movement curves permitting to obtain(a)spring lateral stiffness,(b)spring rocking stiffness,and(c)spring coupling stiffness.

Fig. 6.(a)Real-world soil-structure interaction problem,and(b)The vertical slices model showing the interacting slices subjected to external and body forces.

where the superscripts pc,pr and f lst and for proper contribution of the slice itself,contribution of the preceding slice,and contribution of the subsequent slice,respectively;ai-1,aiand ai+1are the element nodal displacement vectors of slices i-1,i and i+1,respectively;I is the identity matrix;and Gi-1,Giand Gi+1are the shear moduli at slices i-1,i and i+1,respectively.

In Eq.(22),the matrix Dpscorresponds to a problem of plane stresses,which may be given as

(3)Monopiles are relatively shorter and rigid piles with a length to diameter ratio(Lp/Dp)in the range of 2-6 and a diameter(Dp)of up to 8 m envisaged for the next generation of turbines,whereas offshore piled foundations in the offshore oil and gas industry have a length to diameter ratio(Lp/Dp)of over 30 and relevant recommendations have been set on the basis of full-scale loading tests on long,slender and flexible piles with a diameter of 0.61 m(Reese et al.,1974).

From this fact,and unlike a fully 3D or plane strain problem,the value of Poisson’s ratio equal to 0.5 is no longera singular value.It is clear from Eq.(24)that the fictitious body forces applied to a slice i depend essentially on its own nodal displacements and on those of slices sandwiching it.The numerical analysis of the vertical slicing model has led to the familiar equations of a pseudo plane stress problem with body forces representing the interaction between the slices,forming the structure and its surrounding medium.Substituting Eqs.(24)and(25)into Eq.(23)gives a more detailed governing equation:

In a more compact form,Eq.(29)becomes

This equation cannot be solved straight-fully,since the right hand terms are not available explicitly at the same time.Consequently,an updating iterative process is needed:

where j denotes the iteration number and jmaxis the maximum number of iterations allowed in the numerical process.

为了进一步验证交互双模自适应无迹卡尔曼滤波算法的性能，测速电机主轴的运动状态采用式(19)表示的机动性更强的变速模型M2和恒速模型M1交替的形式，采用蒙特卡洛方法仿真200 ms(其中，41-90 ms以及111-160 ms采用变速模型，变速因子ζ分别为1和-1，其余步采用恒速模型。仿真结果如图8至图10所示。

不同课程体系下临床医学本科专业学生心理压力对比分析………………………迟天缘，何思林，范琳琳，缪 倩，陈红霞，黄玉巧，王丽霞，胡吉富（94）

The nonlinearity in vertical slices model stems from the implementation of the hyperbolic model proposed by Duncan and Chang(1970)for modeling the soil.In fact,they found out that both tangential modulus Eiand ultimate stress difference(σ1- σ3)ultare dependent on the minor principal stressσ3.More precisely,they suggested for the initial tangent modulus the following formula:

A Fortran computer program called NAMPULAL for the analysis of axially and laterally loaded single monopiles has been written.Although approximate,the computer code NAMPULALisa coherent tool which exhibits many advantages over other numerical codes in dealing with nonlinear soil-structure interaction problems.For further details,the reader can refer to Otsmane and Amar Bouzid(2018)and only the features of this computer program are given here.

Concerning accurate prediction of the monopile head stiffnesses,numerical analysis using the finite element method(FEM)constitutes an excellent alternative to capture the real behavior of this type of foundations and hence to accurately estimate the dynamic characteristics of an OWT.

Using Mohr-Coulomb failure criterion where the envelope is considered as a straight line,the principal stress difference at failure is related to the confining pressureσ3as

where c is the cohesion andφis the internal Friction angle.The tangent modulus Etis given by

For unloading and reloading cycles,Duncan and Chang(1970)proposed the following expression:

where Euris the unloading-reloading modulus and Kuris the corresponding modulus number.

A thorough literature investigation has been performed by Otsmane and Amar Bouzid(2018)to keep the hyperbolic modeling parameters sufficiently accurate and to make their use practical for solving the soil-structure interaction problems.Theauthorsexamined a large number of well-established correlations between soil physical parameters especially those of sandy deposits whose behaviors are mainly governed by their internal friction angles,and proposed relationships between the sand relative density and the confining pressure.This has been achieved using mainly the recommendations made by well-known researchers who carried out a great number of careful experiments.These parameters are listed in Tables 4 and 5 along with references of their origin.

Table 4 Soil stiffness parameters in terms of soil friction angle and confining pressure.

Note:Esis the modulus of elasticity of sand,k0is the earth pressure coefficient at rest,Dris the relative density andσv0is the overburden pressure.

Equation Source?Es=1025e2.93Dr 1+2k0 3 σv0?0.51 Otsmane and Amar Bouzid(2018)k0=1-sinφ Jâky(1944)

Table 5 Parameters governing the hyperbolic model according to correlations and recommendations.

Equation Source n=0.51 Otsmane and Amar Bouzid(2018)Rf=0.7 Wong and Duncan(1974)Kur=1025e2.93Dr??0.51 1+2k03k0p-0.49a Otsmane and Amar Bouzid(2018)K=0.667KurDuncan and Wong(1999)

Equations in both Tables 4 and 5 have been implemented in the FEM computer code NAMPULAL which will be described in the next paragraph for evaluating soil model parameters related to the five wind farm sites considered in Section 4.In the Duncan-Chang’s basic model,the Poisson’s ratio νswas assumed to be constant throughout the whole process.

where K is the dimensionless factor termed as ‘modulus number’,n is a dimensionless parameter called ‘modulus exponent’,and pais the atmospheric pressure used to make K and n dimensionless.

Although the computational process in NAMPULAL is naturally iterative to fulfill slices equilibrium,it does not require a significant number of iterations to reach convergence.For the problems analyzed so far,a number of 20 iterations are generally sufficient to reach accurate solutions within acceptable margins.

A number of 20 slices have been implemented in NAMPULAL.This number,which has been set on the basis of parametric study involving many monopile behavior parameters(Amar Bouzid et al.,2005),has been found sufficient to accurately model many soilstructure interaction problems(Amar Bouzid et al.,2005;Otsman and Amar Bouzid,2018).

（25）带叶苔 Pallavicinia lyellii（Hook.）Gray，Nat.赵文浪等（2002）；刘胜祥等（1999）；李粉霞等（2011）

Unlike most implemented elastoplastic constitutive models,which necessitate a significant number of iterations to subdue the unbalanced forces,the implemented hyperbolic model in NAMPULAL requires only two iterations.This fact alleviates considerably the whole process of solution and makes it easier to find fast solutions even for the most complex soil-structure interaction problems.

In addition to the rectangular cross-sectional monopiles that are automatically considered due to the shape of the vertical slice,the solid circular or tubular cross-sectional monopiles are easily dealt with by prescribing the effective bending stiffness.Hence,an equivalent Young’s modulus is adopted according to the following formula:

This equation has been set on the assumption that the square cross-sectional monopile under consideration in NAMPULAL has the same cross-sectional area as the effective circular crosssectional monopile.

The performances of this computer code have been assessed against analysis of the behavior of a number of OWT monopiles where other commercial packages are used,such as FLAC3D,ABAQUS®and PLAXIS(Otsmane and Amar Bouzid,2018).The results were in excellent agreement with those of the aforementioned powerful numerical tools.This computer code is employed to determine the monopile head stiffnesses for the OWTs examined in this paper.

Table 6 List of the five OWTs with soil conditions at the sites.

No. Wind farm name Country Soil conditions at the site Sources providing data and measured natural frequencies 1 Lely A2 offshore wind farm UK Soft clay in the uppermost layer to dense and very dense sand layers below Zaaijer(2002),Arany et al.(2016),Amar Bouzid(2016)2 Irene Vorrink offshore wind farm Netherlands Soft layers of silt and clay in the upper seabed to dense sand and very dense sand below Zaaijer(2002),Arany et al.(2016),Amar Bouzid(2016)3 Kentish Flats offshore wind farm UK Layers of dense sand and f i rm clay Arany et al.(2016),Amar Bouzid(2016)4 Walney 1 offshore wind farm UK Medium and dense sand layers Abed et al.(2016)5 North Hoyle wind farm UK Sand and clay layers Leblanc(2009),Arany et al.(2016)

Table 7 Input parameters for the five OWTs chosen for this study.

OWT component dimension Symbol(unit) Lely A2 Irene Vorrink Kentish Flats Walney 1 North Hoyle Tower height LT(m) 37.9 44.5 60.06 67.3 67 Substructure height Ls(m) 12.1 5.2/6 16 37.3 7 Structure height L(m) 50 49.7/50.5 76.06 104.6 74 Tower top diameter Dt(m) 1.9 1.7 2.3 3 2.3 Tower bottom diameter Db(m) 3.2 3.5 4.45 5 4 Tower wall thickness tT(mm) 13 13 22 40 35 Substructure diameter Ds(m) 3.2 3.5 4.3 6 4 Substructure wall thickness ts(mm) 35 28 45 80 50 Young’s modulus of tower material ET(GPa) 210 210 210 210 210 Tower mass mT(t) 31.44 37 108 260 130 Top mass mt(t) 32 35.7 130.8 234.5 100 Monopile diameter Dp(m) 3.2 3.5 4.3 6 4 Monopile wall thickness tp(mm) 35 28 45 80 50 Young’s modulus of monopile material Ep(GPa) 210 210 210 210 210 Monopile depth Lp(m) 13.5 19 29.5 23.5 33 Shear modulus of the soil Gs(MPa) 140 55 60 70 230 Poisson’s ratio of the soil νs 0.4 0.5 0.4 0.4 0.4 Young’s modulus of the soil Es(MPa) 392 165 168 196 644 Measured frequency f1(Hz) 0.634 0.546/0.563 0.339 0.35 0.35

Table 8 Masses and bending stiffnesses of the OWTs constitutive elements.

Wind farm name Tower mass,mT(t) Substructure mass,ms(t) Support structure mass,mTs(t) Structure bending stiffness,EITs(GN m2)Monopile bending stiffness,EIp(GN m2)Lely A2 30.866 33.09 63.956 40.072 142.105 Irene Vorrink 36.955 14.403/12.482 51.358/49.437 35.222/34.234 96.65 Kentish Flats 109.399 75.649 185.048 136.833 285.916 Walney 1 263.234 436.207 699.441 664.58 1369.032 North Hoyle 180.374 34.138 214.512 139.191 254.161

Table 9 Adopted soil deformation and strength parameters as well as hyperbolic model parameters for the OWTs chosen.

Note:γis the soil unit weight.

Wind farm name Es(MPa) νs c(kPa) φ(?) Rf K n Kur Dr γ(kN/m3) k0 Lely A2 392 0.4 0 40 0.7 524.905 0.51 787.358 0.6 21 0.3572 Irene Vorrink 165 0.5 0 40 0.7 524.905 0.51 787.358 0.6 21 0.3572 Kentish Flats 168 0.4 0 40 0.7 524.905 0.51 787.358 0.6 21 0.3572 Walney 1 196 0.4 0 40 0.7 524.905 0.51 787.358 0.6 21 0.3572 North Hoyle 644 0.4 0 40 0.7 524.905 0.51 787.358 0.6 21 0.3572

Fig. 7.Monopile head displacement against applied horizontal load for different OWTs.

Fig. 8.Monopile head rotation against applied horizontal load for different OWTs.

Fig. 9.Monopile head rotation against applied overturning moment for different OWTs.

4.Computed and measured first natural frequencies for five different offshore wind turbines

In order to assess the performances of the computer code NAMPULAL for a wide range of geotechnical applications,in terms of the accuracy,utilityand potential of the NFEVSM,five OWTs have been selected from five wind farm sites.These are Lely A2(UK),Irene Vorrink(Netherlands),Kentish Flats(UK),Walney 1(UK)and North Hoyle(UK).These wind turbines have been chosen for the full availability of their data,especially the measured first natural frequency.Soil conditions at the site and sources from which the OWT data are adopted are summarized in Table 6.

OWTs structural data along with some soil deformation characteristics and measured natural frequencies are presented in Table 7.Since these data come directly from the OWT manufacturers,it is not possible to check their accuracy,except data relevant to the tower mass,provided that the tower height is correct.

Slight differences between computed values of mTand those provided in the reference(Arany et al.,2016)are noticed.Thus computed data in Table 8 are used in the coming computations.

Although the OWT structural data were available which enable the users to compute any structural behavior parameter,the parameters relevant to soil behavior were not found.However,only the different strata of each site are given,but nothing about strength and deformation parameters.

The site investigations indicate that almost all the OWTs chosen in this paper are installed through deep layers of dense sand.The pertinent hyperbolic parameters have been computed according to the prescribed values and relationships given in Tables 4 and 5.These are given in Table 9.

A comprehensive mesh study has been performed to find the optimal finite element mesh that captures the behaviors of monopiles under lateral loading in a nonlinear medium characterized by the hyperbolic model as a yield criterion.A mesh of 20 times monopile diameter Dpin both sides of the monopile and one monopile length Lpunder the monopile tip has been adopted for the study of all OWTs considered here.Furthermore,35 finite elements in both sides of the monopile and 36 finite elements in vertical direction as well as 20 slices have been chosen to analyze the pseudo 3D medium under consideration.

As the monopile head stiffness does not depend on the loading level,a horizontal load H of 1000 kN is applied in 10 increments at the top of each monopile in the five wind farms considered,aiming to compute the monopile head flexibility coefficients ILand ILR.Fig. 7 shows the evolution of monopile head displacements as a function of the increasing lateral load H.

The evolution of rotations is a function of applied lateral load and is plotted in Fig. 8.This figure is used to determine the cross coupling flexibility coefficient ILRfor all monopiles considered here.As the flexibility coefficient IRrequires a pure bending,an applied moment M at the top of monopile of 20,000 kN m in magnitude has been considered and the corresponding rotations are plotted in Fig. 9.

Table 10 Flexibility coefficients IL,ILRand IRand their corresponding stiffness coefficients KL,KLRand KRrelevant to monopiles in the OWTs chosen.

Wind farm name IL(m/GN) IR(rad/(GN m)) ILR(GN-1) KL(GN/m) KR(GN m/rad) KLR(GN)Lely A2 5.788 0.115 -0.571 0.339 17.049 -1.682 Irene Vorrink 4.993 0.132 -0.498 0.321 12.169 -1.213 Kentish Flats 3.413 0.056 -0.268 0.472 28.975 -2.278 Walney 1 2.522 0.018 -0.148 0.755 103.625 -6.096 North Hoyle 3.571 0.0603 -0.289 0.459 27.199 -2.208

Table 11 Fixed base natural frequencies for different OWTs.

Wind farm name fFB(Hz)Lely A2 0.719 Irene Vorrink 0.659-0.669 Kentish Flats 0.368 Walney 1 0.333 North Hoyle 0.404

The almost linear relationships between H and u and H andθon one hand and that between M andθon the other hand somewhat make it easy to compute IL,ILRand IRwhich can be performed by simply inverting the slopes of their corresponding load deformation curves.Then using Eq.(22),the stiffness coefficients are obtained.Flexibility and stiffness coefficients are respectively shown in Tables 10 and 11 for all turbines considered in this paper.

Fig. 10.Values of KL,KLRand KRgiven by NAMPULAL against those developed by Arany et al.(2016).

Table 12

CRand CLcomputed values for the OWTs considered in the current study.

Wind farm name CR CL Lely A2 0.867 0.996 Irene Vorrink 0.867 0.997 Kentish Flats 0.857 0.998 Walney 1 0.837 0.997 North Hoyle 0.849 0.998

Table 13

Predicted and measured natural frequencies of all OWTs.Wind farm namePredicted frequency(fη=CRCLfFB)(Hz)

Measured frequency(Hz)Error(%)

Lely A2 0.621 0.634 2.05

Irene Vorrink 0.570/0.579 0.546-0.563 4.39/2.84 Kentish Flats 0.315 0.339 7.08

Walney 1 0.277 0.35 20.86

North Hoyle 0.342 0.35 2.28

As the monopile head stiffness coefficients play an important role in the correct assessment of the natural frequency,which is in turn a significant parameter in the design of any OWT,it is useful to compare the values listed in Table 10 with those provided by other methods.On the basis of the formulas developed by Poulos and Davis(1980),Randolph(1981)and Carter and Kulhawy(1992).Arany et al.(2016)determined the values of the monopile stiffness coefficients which are added to the histograms of Fig. 10 for comparison.

One important point can emerge from the close examination of the histograms shown in Fig. 10.NAMPULAL’s results are approximately half those given by Arany et al.(2016)for the OWTs whose supporting monopiles are driven in dense sand.We believe that the NFEVSM results are more accurate than those of Arany et al.(2016).This is probably due to the fact that these authors used data from works performed on slender piles using the Winkler model for which many questions had been raised about its applicability to large-diameter monopiles.

Eq.(15),whose different constitutive parts are evaluated using Eqs.(16)and(17),is employed here to give the fixed base natural frequency.This expression,which depends only on the OWT structure properties,gives the values of the fixed base natural frequencies for different turbines shown in Table 11.

The factors CRand CLdepend on values of IL,IRand ILR.Table 12 shows the values of CRand CLfor the five OWTs considered in this paper.These values make it quite clear that CRis the dominant factor that can bring the value of the fixed base frequency to the measured one.However,CLis very close to unity,and hence its influence in changing the value of fFBis very small.This has been also noticed by Arany et al.(2016).

The natural frequency which is simply obtained by multiplying the flexibility coefficients by the fixed base frequency for each OWT is given in Table 13.Also shown are errors between the measured and the computed natural frequencies.5.Conclusions

In this article,a nonlinear finite element computer code NAMPULAL developed for soil-pile interaction has been used to analyze five different monopiles from five European wind farms.

As the natural frequency of the whole wind turbine structure is of paramount importance in the design of OWTs,developing reliable methods for its determination is an active area of research.Stiffness of foundation is the key of natural frequency calculation and this work is based on three-spring model where the monopilesoil interaction is modeled by lateral spring KL,cross-coupling KLR and rotational spring KR.NAMPULAL code has been adopted to find the foundation stiffness.The code has the capability to incorporate nonlinear soil model and in this study,Duncan-Chang hyperbolic model has been used.These stiffnesses in turn were used to obtain natural frequency of the whole wind turbine system and the results obtained were compared with the measurements.Good agreement was noted between prediction and observation.Conflict of interest

The authors wish to confirm that there are no known Conflicts of interest associated with this publication and there has been no significant financial support for this work that could have influenced its outcome.

References

Abed Y,Amar Bouzid Dj,Bhattacharya S,Aissa MH.Static impedance functions for monopiles supporting offshore wind turbines in non-homogeneous soils emphasis on soil/monopile interface characteristics.Earthquakes and Structures 2016;10(5):1143-79.

Achmus M,Kuo YS,Abdel-Rahman K.Behavior of monopile foundations under cyclic lateral load.Computers and Geotechnics 2009;36(5):725-35.

Adhikari S,Bhattacharya S.Dynamic analysis of wind turbine towers on flexible foundations.Shock and Vibrations 2012;19(1):37-56.

Adhikari S,Bhattacharya S.Vibrations of wind-turbines considering soil-structure interaction.Wind and Structures 2011;14(2):85-112.

Aissa MH,Amar Bouzid Dj,Bhattacharya S.Monopile head stiffness for serviceability limit state calculations in assessing the natural frequency of offshore wind turbines.International Journal of Geotechnical Engineering 2017;35(7).https://doi.org/10.1080/19386362.2016.1270794.

Amar Bouzid Dj.Monopile head stiffness and natural frequency assessment of some installed OWTs using a pseudo 3D nonlinear FE model.In:World congress on Advances in Civil,Environmental and Materials research(ACEM16).Guildford,Surrey,UK:The University of Surrey;2016.p.1-15.http://epubs.surrey.ac.uk/id/eprint/844988.

Amar Bouzid Dj,Vermeer PA,Tiliouine B.Finite element vertical slices model:validation and application to an embedded square footing under combined loading.Computers and Geotechnics 2005;32(2):72-91.

American Petroleum Institute(API),International Organization for Standardization(ISO).ANSI/API specification RP 2GEO:geotechnical considerations and foundation design for offshore structures.Washington,D.C.,USA:API;2011.

Arany L,Bhattachary S,Hogan SJ,Macdonald J.Dynamic soil-structure interaction issues of offshore wind turbines.In:Proceedings of the 9th international conference on structural dynamics,Porto,Portugal;2014.p.3611-7.

Arany L,Bhattacharya S,Adhikari S,Hogan SJ,Mcdonald JHG.An analytical model to predict the natural frequency of offshore wind turbines on three-spring flexible foundations using two different beam models.Soil Dynamics and Foundation Engineering 2015;74:40-5.

Arany L,Bhattacharya S,Adhikari S,Mcdonald JHG,Hogan SJ.Closed form solution of Eigen frequency of monopile supported offshore wind turbines in deeper waters stiffness of substructure and SSI.Soil Dynamics and Foundation Engineering 2016;83:18-32.

Bhattacharya S.SDOWT:simplif i ed dynamics of wind turbines.User’s Manual.Bristol,UK:Bristol Laboratory for Advanced Dynamics Engineering;2011.

Blevins RD.Formulas for frequencies and mode shapes.Krieger Pub Co;2001.

Carter JP,Kulhawy FH.Analysis of laterally loaded shafts in rock.Journal of Geotechnical Engineering 1992;118(6):839-55.

Damgaard M,Bayat M,Andersen LV,Ibsen LB.Dynamic response sensitivity of an offshore wind turbine for varying subsoil conditions.Ocean Engineering 2015;101:227-34.

DNV-OS-J101.Offshore standard:design of offshore wind turbine structures.Hellerup,Denmark:Det Norske Veritas(DNV);2004.

Duncan JM,Chang CY.Nonlinear analysis of stress and strain in soils.Journal of Soil Mechanics and Foundation Division 1970;96(5):1629-53.

Duncan JM,Wong KS.Soil properties manual.User’s manual for SAGE,vol.II.Blacksburg,Virginia,USA:Centre for Geotechnical Practice and Research,Virginia Polytechnic Institute and State University;1999.

Higgins W,Vasquez C,Basu D,Griffiths DV.Elastic solutions for laterally loaded piles.Journal for Geotechnical and Geoenvironmental Engineering 2013;139(7):1096-103.

Jâky JA.Nyugalmi nyomâs tényezöje(The coefficient of earth pressure at rest).Agyar Mérnok és Epitész Egylet Közlönye(Journal for Society of Hungarian Architects and Engineers)1944:355-8(in Hungarian).

Leblanc C.Design of offshore wind turbine support structures-selected topics in the field of geotechnical engineering.PhD Thesis.Aalborg,Denmark:Aalborg University;2009.

Lombardi D,Bhattacharya S,Wood DM.Dynamic soil-structure interaction of monopile supported wind turbines in cohesive soil.Soil Dynamics and Earthquake Engineering 2013;49:165-80.

Otsmane L,Amar Bouzid D.An efficient FE model for SSI:theoretical background and assessment by predicting the response of large diameter monopiles supporting OWECs.Computers and Geotechnics 2018;97:155-66.

Poulos HG,Davis EH.Pile foundation analysis and design.New York,USA:Wiley;1980.

Prendergast LJ,Gavin K,Doherty P.An investigation into the effect of scour on the natural frequency of an offshore wind turbine.Ocean Engineering 2015;101:1-11.

Randolph MF.The response of flexible piles to lateral loading.Géotechnique 1981;31(2):247-59.

Reese L,Cox WR,Koop FD.Analysis of laterally loaded piles in sand.In:Proceedings of the offshore technology conference,Houston,USA;1974.https://doi.org/10.4043/2080-MS.

Vught JH.Considerations on the dynamics of support structures for an offshore wind energy converter.PhD Thesis.Delft,the Netherlands:Delft University of Technology;2000.

Wong KS,Duncan JM.Hyperbolic stress-strain parameters for non-linear finite element analyses of stresses and movements in soil masses.Report No.TE-74-3.Berkeley,USA:University of California;1974.

Yi JH,Kim SB,Yoon GL,Andersen LV.Natural frequency of bottom-fixed offshore wind turbines considering pile-soil-interaction with material uncertainties and scouring depth.Wind and Structures 2015;21(6):625-39.

Yu L,Bhattacharya S,Li L,Guo Z.Dynamic characteristics of offshore wind turbines on different types of foundations.Electronic Journal of Geotechnical Engineering 2014;19:2917-36.

Zaaijer MB.Design methods for offshore wind turbines at exposed sites(OWTES)-sensitivity analysis for foundations of offshore wind turbines.Technical report.Delft,The Netherlands:Delft University of Technology;2002.

Zaaijer MB.Foundation modeling to assess dynamic behavior of offshore wind turbines.Applied Ocean Research 2006;28(1):45-57.