考虑桩-弹性地基相互作用的单桩风机自振频率

牛文杰

(辽宁工程技术大学力学与工程学院，辽宁 阜新 123000)



考虑桩-弹性地基相互作用的单桩风机自振频率

牛文杰

(辽宁工程技术大学力学与工程学院，辽宁 阜新 123000)

考虑桩-弹性地基相互作用，采用集中质量法和柔度法对泥面线以上的单桩风机结构进行多自由度动力分析，确定单桩风机的自振频率。通过算例给出单桩风机的自振频率值，并与不考虑桩-弹性地基相互作用的单桩风机的自振频率值进行比较。

单桩风机； 自振频率； 桩基础和弹性地基的相互作用； 柔度法； 集中质量法

0 引言

对于水平轴风力发电机[1]，大型风力机塔筒通常都采用变截面的锥形筒体形式，可将其视为顶端受横向力作用的变截面悬臂梁结构[2]。横向力由风力机叶片、轮毂和机舱的迎风阻力组成[3]。风力发电机组运行时，其叶轮上的风荷载和风机偏航引起的荷载通过结构和传动机构作用在塔架顶端[4-5]。已建成的大部分海上风电场都采用单立柱基础[4-6]。风机塔上作用有自重荷载、波浪荷载、流体流过立柱的流体压力荷载[5]以及可能的近海地震作用与浮冰冲击作用。此外，风机叶片的转速一般是常数或在某一指定范围，可通过闭环控制系统实现[7]。因此，对于叶片转速一定的风机，当以塔筒为研究对象时，其上的动力荷载还主要有：轮毂等转子系统施加的动荷载，其频率为一定值，记为f；当风机叶片在转动过程中遮挡住塔筒时，改变了环境流场，会施加给塔筒荷载，如果是三叶片的风机，则叶片转动引起的作用在塔筒上的振动荷载频率是3f[8]。

在上述动力荷载作用下，需要通过工程设计，使现场近海单桩风机的自振频率避开上述动力荷载的激励频率，以避免近海单桩风机在动力荷载作用下发生共振。另外一个工程问题就是：对于近海风机这种高而柔的结构物需要承受各种动力和循环荷载，其在长期的水平荷载作用下，由于桩周土在动应力作用下的材料弱化，风机的自振频率会发生变化。因此需要研究确定近海单桩自振频率的计算方法或实验方法。

确定近海单桩自振频率的实验方法通常采用Bhattacharya等(2011)建议的模型试验，以保证模型与现场风机的相似准数一致[9]。

确定近海单桩自振频率的理论解有以下两方面的研究：

(1) 计算近海单桩风机的自振频率需要考虑桩土相互作用。而关于桩土相互作用方面的研究有：Zaaijer[10]研究了单桩-土体的动力特性，土体被简化为一系列弹簧；Padrn等[11]通过边界元-有限元耦合模型对埋在半空间弹性体中的桩和群桩进行了动力分析，其中桩按照伯努利假设的梁单元进行模拟，而土按照连续的、半无限体、各向同性的线弹性介质通过边界元来进行模拟；Dash等[12]使用数值方法研究了土体中桩发生压-弯耦合时的性状。

(2) 牛文杰等[13]使用集中质量法，采用刚度法对近海风力发电机的塔和泥面线上基础按照2个自由度的动力系统进行了自由振动分析,并得到自振频率。但考虑地基-风机相互作用时的计算模型更加符合实际。Bhattacharya[14]按照无限自由度体系，把地基对塔筒的作用简化为转动弹簧和平动弹簧，弹簧刚度可使用试验测得，并使用解析解计算了风机的自振频率，模型试验和有限元解的结果表明计算结果是合理的。另外，可使用p-y曲线法对海洋平台的桩土结构系统进行模态分析[15]。

上述计算自振频率的方法，或过于简略，或过于复杂，或需要大量的现场测试数据，或需要借助于Ansys软件进行分析计算，因此需要一种快速估算近海单桩风机自振频率的计算方法。

风力发电机的上部塔筒和桩基础都假定为弹性模量相同的钢材制成的套筒。通过张氏法[16-17]得到半无限长桩基础的控制方程；然后通过欧拉-伯努利经典梁理论方程得到泥面线上以上风机上部塔筒的变形控制方程；再通过计算变形的叠加原理与套筒在泥面线处的变形连续性条件，得到泥面线上风机立柱(即套筒)的挠度函数；最后使用柔度法，给出按照多自由度分析的考虑桩-弹性地基相互作用时风机的自振频率值。

1 考虑桩-弹性地基相互作用的自由振动理论

1.1 考虑桩土相互作用时单桩结构物的挠度

图1所示为单桩风力发电机。以下分析中，把风力机简化为截面和弹性模量都相同的单桩结构物(图2)。风机高度为h。

当水平荷载作用于杆件顶端时，对于图2中的单桩结构物即为等截面杆件在泥面线以上的水平挠度y(x)F。这可按照两部分来计算：一部分是桩底端固定时的位移曲线函数y(x)F0；另一部分是假定无外荷载作用，泥面线位置处单桩结构物转角θ0和位移y0引起的桩的水平挠度y(x)F1。最后通过叠加原理可确定y(x)F。

(1)

图1 单桩式风力发电机Fig.1 A mono-pile wind turbine

图2 单桩结构物顶部(杆件)作用水平荷载FFig.2 Mono-pile structure top under horizontal load F

在图2中，以泥面线与杆件的交点为坐标原点，记x是泥面线上杆件某点与原点的距离。根据材料力学，可确定y(x)F0：

(2)

在图2中，泥面线位置处单桩结构物的转角和位移引起的单桩结构物在泥面线以上的水平挠度

(3)

式中：θ0和y0分别是泥面线位置处单桩结构物立柱的转角和挠度(图2)。

综合式(1)、(2)和(3)，可得图2中单桩结构物在顶部作用有水平集中荷载时的水平挠度：

(4)

(1) 泥面线位置处单桩结构物的转角和位移

图2中以泥面线与单桩结构物的交点为坐标原点，记z是泥面线下桩基础某点与原点的距离， 即桩基础的水平挠度。当单桩结构物顶部作用有水平集中力F时，泥面线以下半无限长桩在坐标z处的水平挠度y的控制方程为[16-17]：

(5)

式(5)中物理量的含义可查找文献[16-17]。

泥面线位置处：

(6)

将式(6)代入式(5)，得出泥面线位置处单桩结构物的挠度：

(7)

在单桩结构物发生小挠度的假定下，单桩结构物轴线的转角θ在任一点满足以下关系：

(8)

观察图2可发现：

(9)

由式(5)可得单桩结构物轴线在泥面线以下任一点的转角的正切值为：

(10)

泥面线位置处单桩结构物轴线的纵坐标满足式(6)，将式(6)代入式(10)，得出泥面线位置处单桩结构物轴线转角的正切值：

(11)

由式(9)和(11)得出在图2所示的水平力F作用下，泥面线位置处单桩结构物轴线转角的正切值为：

(12)

根据式(8)和(12)得泥面线位置处单桩结构物轴线转角:

(13)

(2) 单桩结构物的水平挠度y(x)F

结合式(4)、(7)和(13)，得出在图2所示的单桩结构物在顶部作用有水平集中荷载F时，泥面线以上的水平挠度：

(14)

单桩结构物顶端的位移为：

(15)

因此，式(14)代入式(15),得出单桩结构物顶端的水平挠度为：

(16)

1.2 多自由度分析

按照柔度法对图1的结构进行自振频率计算。记机舱、轮毂和叶片的总重量为m0，泥面线以上套筒质量为M0。按照集中质量法[20]，在图3中m1=0.25*M0+m0,m2=0.5*M0。

泥面线以上塔筒高度为h。根据式(16)和(14),图3中在水平集中力F作用下，集中质量位置处的水平挠度分别为：

图3 计算m1位置处的柔度系数Fig.3 Computation of flexibility coefficient at the position of m1

因此在图3中,柔度系数[21]

根据式(16),图4中在水平集中力F的作用下，集中质量位置处的水平挠度为：

(17)

根据式(14), 图4中m2以下位置处的挠度函数为

(18)

图4 计算m2 位置处的柔度系数Fig.4 Computation of flexibility coefficient at the position of m2

(19)

将式(17)和(18)代入式(19)，得：

(20)

因此在图4中，柔度系数[21]

(21)

(22)

上述柔度系数确定后，即可建立集中质点在惯性力作用下的位移方程，进而根据柔度法计算出图3和图4中2个自由度动力系统的自振频率。

2 算例

图1中，风力发电机的塔和基础都假定为弹性模量[18]E=2.1×1011N/m2的钢材制成的套筒。

套筒直径d=4 000 mm,泥面线以上高度h=90 m,桩入土110 m。参考文献[19]中水下桩的桩径和壁厚的比值，图1中套筒的壁厚采用52 mm。根据环形截面惯性矩计算公式得I=1.256 m4。因此钢套筒的环形截面抗弯刚度EI=2.1×1011N/m2×1.256 m4=2.6×1011N·m2。

通过计算，按照钢材密度7.85×103kg/m3,计算出泥面线以上套筒质量M0=4.5×105kg。机舱重量1.0×105kg,轮毂重量0.2×105kg，叶片重量3×0.1×105kg[19]，则机舱、轮毂和叶片的总重量m0=1.5×105kg。

上述设计参数都是根据相关参考文献的设计值组合得到的，与文献[13]的数据相同，但其自振频率分析没有考虑桩基础以及桩土相互作用。这些设计参数是否满足初步设计要求，需根据现场记录的荷载统计数据进行套筒的弯曲正应力和剪应力计算，然后根据材料力学的强度理论进行判断，最后逐步进行精确有限元分析和优化设计。

假设地基土的水平反力系数[16]k=200 kPa，即地基土的水平反力系数比极软淤泥的水平反力系数还小，根据集中质量法及算例参数可知：

m1=0.25*M0+m0=2.625*105kg，

m2=0.5*M0=2.25*105kg,

h=90m，EI=2.6×1011N·m2，β=0.030m-1

根据柔度法[21]，即可计算出图1泥面线以上结构的自振频率(计算过程见附录的fortran源程序)：

第一频率为1.141 Hz,第二频率为8.421 Hz。

3 讨论

3.1 与其他算法计算的自振频率的比较

对于参数完全相同的结构(图1), 当不考虑桩土相互作用时，牛文杰等[13]使用集中质量法，采用刚度法对近海风力发电机的塔和泥面线上基础按照2个自由度的动力系统进行自由振动分析, 得到自振频率，经过验算，文中M应为(0.25*M0+m0),m应为0.5*M0。重新根据刚度法计算，第一频率是2.050 Hz,第二频率是13.559 Hz；根据柔度法计算，第一频率是1.928 Hz,第二频率是13.597 Hz。在地面按照9.8sin(4.626t)的规律运动时(单位m)，质点M的位移幅值AM=-15.084 m，质点m的位移幅值Am=-4.505 m。

而本文考虑桩土相互作用，计算出理论上更加精确的自振频率：第一频率1.141 Hz,第二频率8.421 Hz。

3.2 影响自振频率的因素

在自振频率的计算过程中，通过观察知泥面线以上风机结构的自振频率取决于桩的惯性矩EI、桩径d、地基土的水平反力系数k、桩基础入土深度、风机在泥面线以上的高度h及泥面线以上风机的质量分布。

3.3 改进方向

在自振频率的计算中，文献[9,22]考虑了叶片、机头等部件的重力作用，这使得计算结果更加符合实际。在考虑重力作用的情况下，根据柔度法，单桩结构物的水平挠度会增大。实际上也应该考虑塔筒、桩基础的重力作用，但其计算过程会更加复杂。

另外，对风机整体，当叶片的回转直径达到91.6m时，若风力发电机的叶片、机头等部件简化为一集中质量块进行自振频率的计算，其结果是否合理还需实验检验。风机的支撑结构截面一般都是从下到上逐渐变化的(变小)，还需在今后计算时进行改进。

桩实际上是有一定长度的，桩端的支撑条件对风机体系的自振频率存在一定影响。本文将桩作为半无限长考虑是一种工程假定，实际近海风机桩基础的入土深度需要根据近海地质资料、海洋环境荷载和正常使用极限等因素确定。实际工程中地基大多数是由不同岩土层组成的(例如海砂-软土互层)，因此地基对桩的横向支撑刚度也必然是随深度变化的，本文没有考虑这一变化，与实际情况不符，需在今后加以考虑，例如使用p-y曲线而非张氏法考虑桩土相互作用。桩周土在复杂动应力[23]作用下会发生应力-应变关系的变化，进而影响单桩风机的自振频率，也需在今后计算中加以考虑[24]。

4 结论

以泥面线以上的风机上部结构为研究对象，给出快速估算近海单桩风机自振频率的计算方法。按照集中质量法，把泥面线以上的风机结构按照多自由度动力系统进行分析，计算出更加精确的自振频率值。

观察自振频率的计算过程，可知泥面线以上风机结构的自振频率取决于桩的惯性矩EI、桩径d、地基土的水平反力系数k、桩基础入土深度、风机在泥面线以上的高度h及泥面线以上风机结构的质量分布。

)

[1] 赵先民.水平轴风力机动力特性和流场的数值模拟[D].上海:同济大学,2007.

ZHAOXian-min.DynamicCharacteristicsofaHorizontalAxisWindTurbineandNumericalSimiluationofFlowField[D].Shanghai:TongjiUniversity,2007.(inChinese)

[2] Wang F,Chen Q,Yu G-C.Research on Tower Stiffness of Large Scale Wind Turbine Unit[J].New Energy and Technology,1997,20(6):38-39.

[3] 汤炜梁，袁奇，韩中合.风力机塔筒抗台风设计[J].太阳能学报,2008,29(4):422-427.

TANG Wei-liang,YUAN Qi,HAN Zhong-he.Withstanding Typhoon Design of Wind Turbine Tower[J].Acta Energiae Solaris Sinica,2008,29(4):422-427.(in Chinese)

[4] 王懿，段梦兰，尚景宏，等.海上风机基础结构力学分析[J].中国海洋平台，2009,24(4):14-20.

WANG Yi，DUAN Meng-lan，SHANG Jing-hong，et al. Analysis of Mechanics for Offshore Wind Turbines Founadtion[J].China Offshore P1atform,2009,24(4):14-20.(in Chinese)

[5] Merle C Potter,David C Wiggert.Mechanics of Fluids[M].Beijing:China Machine Press,2003.

[6] 尤汉强,杨敏.循环荷载下海洋风机单桩基础简化分析模型[J].岩土工程学报,2010,32(增刊2):13-16.

YOU Han-qiang,YANG Min.Simplified Model for Single-pile of Offshore Wind Turbine Under Cyclic Loading[J].Chinese Journal of Geotechnical Engineering,2010,32(Supp.2):13-16.(in Chinese)

[7] Jan van der Tempel,David-pieter Molenaar.Wind Turbine Structural Dynamics——A Review of the Principles for Modern Power Generation,Onshore and Offshore[J].Wind Engineering,2002,26(4):211-220.

[8] Subhamoy Bhattacharya,James A Cox,Domenico Lombardi,et al.Dynamics of Offshore Wind Turbines Supported on Two Foundations[J].Geotechnical Engineering,2013,166(2):159-169.

[9] Subhamoy Bhattacharya，Domenico Lombardi，David Muir Wood.Similitude Relationships for Physical Modelling of Monopile-supported Offshore Wind Turbines[J].International Journal of Physical Modelling in Geotechnics,2011,11(2)：58-68.

[10] Zaaijerm M B,Foundation Modeling to Assess Dynamic Behaviour of Offshore Wind Turbines[J].Applied Ocean Research,2006,28(1):45-57.

[12] Suresh R Dash,Subhamoy Bhattacharya,Anthony Blakeborough.Bending-buckling Interaction as a Failure Mechanism of Piles in Liquefiable Soils[J].Soil Dynamic and Earthquake Engineering,2010,30(1/2):32-39.

[13] 牛文杰, 王振宇, 李洪然.海上风机-基础水平振动的多自由度分析[J].土木建筑与环境工程, 2011, 33 (增刊1):121-124.

NIU Wen-jie,WANG Zhen-yu,LI Hong-ran.The Multiple-degree-of-freedom Analysis of Offshore Wind Turbine and Foundation in Horizontal Vibration[J].Journal of Civil, Architectural & Environmental Engineering,2011,33(Supp.1):121-124.(in Chinese))

[14] Bhattacharya S,Adhikari S.Experimental Validation of Soil-structure Interaction of Offshore Wind Turbines[J].Soil Dynamics and Earthquake Engineering,2011,31(5):805-816.

[15] 杨进，刘书杰.ANSYS在海洋石油工程中的应用[M].北京：石油工业出版社,2010.

YANG Jin,LIU Shu-jie.The Application of ANSYS in Offshore Oil Engineering[M].Beijing:Petroleum Industry Press,2010.

[16] 刘金砺.桩基础设计与计算[M].北京：中国建筑工业出版社，1990.

LIU Jin-Li.Design and Calculation of Pile Foundation[M].Beijing:China Architecture & Building Press,1990.(in Chinese)

[17] 卢世深,林亚超.桩基础的计算和分析[M].北京:人民交通出版,1987.

LU Shi-shen,LIN Ya-chao.Calculation and Analysis of Pile Foundation[M].Beijing:China Communications Press,1987.

[18] Abdel-rahman K,Achmus M.Finite Element Modelling of Horizontally Loaded Monopile Foundations for Offshore Wind Energy Converters in Germany[C]// International Symposium on Frontiers in Offshore Geotechnics.Perth,Australia,2005:309-396.

[19] 郇彩云.海上桩式风机基础结构设计与研究[D].大连：大连理工大学,2009.

XUAN Cai-yun.Design,Research of Pile-style Foundation Structures for Offshore Wind Turbine[D].Dalian:Dalian University of Technology,2009.(in Chinese)

[20] 龙驭球,包世华.结构力学[M].北京：高等教育出版社,1994.

LONG Yu-qiu,BAO Shi-hua.Structural Mechanics[M].Beijing:Higher Education Press,1994.(in Chinese)

[21] 周竞欧,朱伯钦,许哲明.结构力学[M].上海:同济大学出版社,2004.

ZHOU Jing-ou,ZHU Bo-qin,XU Zhe-ming.Structural Mechanics[M].Shanghai:Tongji University Press,2004.(in Chinese)

[22] 余璐庆.海上风机桶形基础安装与支撑结构动力特性研究[D].杭州：浙江大学,2014.

YU Lu-qing.Study on the Installation Behavior of Suction Caisson and the Dynamic Properties of Offshore Wind Turbine Structures[D].Hangzhou:Zhejiang University,2014.(in Chinese)

[23] 左照坤,童朝霞.近海风机单桩基础桩周土应力特征分析[J].地震工程学报,2014,36(3):549-554.

ZUO Zhao-kun,TONG Zhao-xia.Stress Characteristics of Soil around the Pile of a Monopile Foundation in Offshore Wind Turbines[J].China Earthquake Engineering Journal,2014,36(3):549-554.(in Chinese)

[24] 余璐庆，王立忠，Bhattacharya.海上风机支撑结构动力特性模型试验研究[J].地震工程学报,2014,36(4):797-803.

YU Lu-qing,WANG Li-zhong,Bhattacharya.Scaled Model Test Study of the Dynamic Behavior of an Offshore Wind Turbine Support Structure[J].China Earthquake Engineering Journal,2014,36(4):797-803.(in Chinese)

附录 计算图1泥面线以上结构自振频率的FORTRAN源程序

REAL M1,M2,H,EI,BETA,F11,F21,F22,F12,LMDA1,LMDA2,OMIGA1,OMIGA2

M1=2.625*100000

M2=2.25*100000

H=90

EI=2.6*100000*100000*10

BETA=0.030

F11=-H*H*H/(6*EI)+H*H*H/(2*EI)+((1+BETA*H)/(2*EI*BETA*BETA)+H/(2*EI*BETA))*H+(1+BETA*H)/(2*EI*BETA*BETA*BETA)

F21=-H*H*H*0.5*0.5*0.5/(6*EI)+H*0.5*H*0.5*H/(2*EI)

F21=F21+((1+BETA*H)/(2*EI*BETA*BETA)+H/(2*EI*BETA))*0.5*H+(1+BETA*H)/(2*EI*BETA*BETA*BETA)F22=H*H*H*0.5*0.5*0.5/(3*EI)+((1+BETA*0.5*H)/(2*EI*BETA*BETA)+0.5*H/(2*EI*BETA))*0.5*H+(1+BETA*0.5*H)/(2*EI*BETA*BETA*BETA)

F12=H*H*H*0.5*0.5*0.5/(3*EI)+H*0.5*H*0.5*0.5*H/(2*EI)

F12=F12+((1+BETA*0.5*H)/(2*EI*BETA*BETA)+0.5*H/(2*EI*BETA))*H+(1+BETA*0.5*H)/(2*EI*BETA*BETA*BETA)

LMDA1=F11*M1+F22*M2

LMDA1=LMDA1+SQRT((F11*M1+F22*M2)*(F11*M1+F22*M2)-4*(F11*F22-F12*F12)*M1*M2)

LMDA1=LMDA1/2

LMDA2=F11*M1+F22*M2

LMDA2=LMDA2-SQRT((F11*M1+F22*M2)*(F11*M1+F22*M2)-4*(F11*F22-F12*F12)*M1*M2)

LMDA2=LMDA2/2

OMIGA1=SQRT(1/LMDA1)

OMIGA2=SQRT(1/LMDA2)

WRITE(*,*)F11

WRITE(*,*)F21

WRITE(*,*)F22

WRITE(*,*)F12

WRITE(*,*)LMDA1

WRITE(*,*)LMDA2

WRITE(*,*)OMIGA1

WRITE(*,*)OMIGA2

end program

Natural Frequency of Mono-pile Wind Turbine Considering Interaction between Pile and Elastic Foundation

NIU Wen-jie

(SchoolofMechanicsandEngineering,LiaoningTechnicalUniversity,Fuxin123000,Liaoning,China)

In this paper, we analyze the mono-pile wind turbine, comprising blade bodies, a hub, a cabin, a conical shell tower, a transition piece, and a mono-pile foundation. Our study objective was to determine the natural frequency of the mono-pile wind turbine above the mudline. We assumed the conical shell tower, transition piece, and mono-pile foundation to be a vertical steel tube (VST) with a uniform cross section and with the same elastic modulus. We used the lumped-mass method to simplify the mono-pile wind turbine above the mudline as two concentrated masses. One concentrated mass at the VST top is the summation of the mass of the blade bodies, hub, cabin and one-fourth the mass of the VST above the mudline. The other concentrated mass at the middle of the VST is one-half the mass of the VST above the mudline. Taking the mono-pile foundation below the mudline as a semi-infinite long pile, we used Zhang’s method to obtain the horizontal stiffness and rotational stiffness of the pile top. In this way, we can determine the deflection of the VST under horizontal load at any point of the wind turbine above the mudline. The flexibility coefficient of the unit horizontal load at any position of the VST can be determined accordingly. Then, using both the flexibility and lumped-mass methods, we can determine the natural frequency of the mono-pile wind turbine above the mudline. Lastly, with an example, we retrieved the natural frequency value of the mono-pile wind turbine above the mudline. We then compared the analysis result with the mono-pile wind turbine natural frequency obtained without considering the interaction between the pile and elastic groundwork, and found little difference between them. The new theory we propose in this paper is more accurate theoretically and should be validated experimentally.

mono-pile wind turbine; natural frequency; interaction between pile and elastic foundation; flexibility method; lumped-mass method

2016-04-26

辽宁工程技术大学博士科研启动基金(11-415)

牛文杰(1982-)，男，河南漯河人，副教授，博士，从事近海单桩风机工程中岩土力学、结构力学和流体力学问题的研究。E-mail: nwj1982@sohu.com。

TU352.2; TU411.3; TU473.1

A

1000-0844(2016)05-0713-07

10.3969/j.issn.1000-0844.2016.05.0713