0 引　言

全球超过80%的海上风能资源都蕴藏在40 m水深以上海域^[1]，我国的深远海风能资源十分丰富，发展海上浮式风机技术是我国风电发展的要求。近年来，随着浮式风机技术快速发展，单机容量增大，考虑浮式风机气动−水动−伺服-弹性全耦合的计算方法受到人们的关注。

国内外学者在浮式风机的时域耦合模拟分析方面展开了大量的研究，着重关注浮式风机的动态响应特征。Marco等^[2]将FAST与OrcaFlex耦合，利用OrcaFlex对系泊缆进行有限元分析，进一步提高了浮式风机耦合系统的计算精度。Wang等^[3]利用气动-水动-伺服-弹性全耦合仿真工具Sima建立了5 MW浮式风机计算模型并进行了时域分析。Rong等^[4]基于全耦合方法对3种不同的单柱式风机进行研究，发现随着外界环境条件恶化，基础响应结果呈上升态势，其中纵荡、纵摇为主要的响应形式。Deng等^[5]利用自编程序研究了浮式基础运动对风电机组动态特性的影响，结果表明基础运动通过改变转子在风场中的姿态和位置对空气动力载荷和叶片变形产生影响。李焱等^[6]考虑系泊系统的拉伸-弯曲-扭转变形研究了不同系缆计算方法对浮式风机动力学响应的影响。张成祥^[7]将FAST软件与自编程序相结合，采用不同的耦合方法对支撑5 MW风机的OC3基础进行了时域计算，结果表明，不同耦合方法的纵荡、纵摇响应结果有明显差异。徐应瑜等^[8]研究了10 MW单柱式浮式风机的运动特性，结果表明风机单机容量增大后，风载荷引起的运动比例增大。Zhang等^[9]对极端工况下系缆的断裂机理进行研究。目前，浮式风机设计中的上部风机和下部基础通常由不同单位设计，基础设计单位由于缺乏上部风机参数而不能进行全耦合计算，戴琼霖等^[10]对此提出了限制性耦合方法，对比了不同耦合方法的差异。

本文进一步提出了限制性强耦合和限制弱耦合方法，并与全耦合计算结果进行了对比，以工作于60 m过渡水深的10 MW浮式风机为例，分析3种耦合方式之间的差异，对浮式基础设计提供理论指导。

1 浮式风机计算理论 1.1 浮式风机系统结构动力学方程

进行耦合计算时，将浮式基础、轮毂、机舱处理为刚体，塔柱、叶片、系泊缆处理为弹性体，采用非线性方法进行同步动态分析，实现风机各子部分间的耦合。考虑风、浪、流载荷，浮式基础的六自由度运动方程如下：

$ \begin{split} & (\boldsymbol{M}+\boldsymbol{m})\ddot{X}\left(t\right)+\boldsymbol{C}\dot{X}\left(t\right)+\boldsymbol{K}X\left(t\right)+\int_{0 }^{ t}\boldsymbol{R}(t-\tau)\dot{X}\left(\tau\right)\mathrm{d}\tau= \\ & \quad\boldsymbol{F}_{\text{wave}}+\boldsymbol{F}_{\text{current}}+\boldsymbol{F}_{\text{moor}}+\boldsymbol{F}_{\text{up}}。\end{split} $

(1)

其中：M为质量矩阵；m为附加质量矩阵；C为阻尼矩阵；K为静水恢复刚度矩阵；F_wave为浮体所受波浪力矩阵；F_current为浮体所受流力矩阵；F_moor为浮体所受系缆力矩阵；F_up为风机上部向下传递的力矩阵；R为速度脉冲函数矩阵。

系泊系统、塔柱、叶片等柔性结构采用梁理论建模并基于有限元方法离散^[11]。

1.2 风载荷

位于水面以上的浮体结构受到风载荷作用。对于塔柱及停机状态的叶片，采用绕流理论计算风载荷：

$ {f}_{L,D}=\frac{1}{2}\rho {{V}_{w}}^{2}{C}_{L,D}。$

(2)

式中：$ {f}_{L,D} $为升力或阻力风压；$ {C}_{L,D} $为截面升力与阻力系数；$ \rho $为空气密度；$ {V}_{w} $为风速。整个塔柱或停机状态的叶片所受风力可对式(2)在长度L上积分获得：

$ F_{L,D}=\int_{0 }^{L }f_{L,D}\mathrm{d}L。$

(3)

对于作业工况下的叶片，风载荷主要为气动载荷，使用叶素动量理论来求解，具体求解过程与理论可参考文献[7]。

1.3 波浪载荷

对于波浪载荷，基于三维势流理论计算，在Sesam/Hydrod软件计算获得一阶波浪力、平均漂移力、附加质量、静水恢复刚度、势流阻尼的水动力传递函数后，导入Sima进行时域计算。波浪载荷考虑了随机特性且采用JONSWAP波浪谱，其适用于有限风区海域成长阶段的风浪，谱形式为：

$ {S}_{J}\left(w\right)={A}_{\gamma }{S}_{PM}\left(w\right){\gamma }^{\mathrm{exp}\left(-0.5{\left(\frac{w-{w}_{p}}{\sigma {w}_{p}}\right)}^{2}\right)}。$

(4)

式中：γ为无量纲峰形参数对JONSWAP谱的影响见图1；σ为谱宽参数；Aγ = 1−0.287ln(γ)。

1.4 流载荷

考虑流载荷的作用时，浮式基础水下部分受到的流载荷按照如下经验公式计算：

$ F=\frac{1}{2}{C}_{D}{\rho }_{W}{{v}_{流}}^{2}A 。$

(5)

式中：C_D为曳力系数；ρ_w为海水密度；v _流为设计海流流速；A 为投影面积。

1.5 系泊系统

系泊缆受力及位置分析的方法有准静态法和动态法，其中准静态法以悬链线理论为代表，动态求解方法以有限元理论为代表。有限元方法可更加真实地模拟系泊缆截面，且在系泊缆发生大偏移时仍可保持高准确度，有限元方法离散系泊缆的动力方程为：

$ {R}^{I}(r,\ddot{r},t)+{R}^{D}(r,\dot{r},t)+{R}^{S}(r,t)={R}^{E}(r,\dot{r},t)。$

(6)

式中：R^I为惯性力，R^I($ \mathit{r} $,$ \ddot{\mathit{r}} $,t)=$ \ddot{\mathit{r}} $[ M^S+ M^F(r)+ M^H(r)]，M^S为结构质量，M^F(r)为内部流体质量矩阵，M^H(r)为水动力质量矩阵，主考虑莫里森方程的结构加速度项；R^D为阻尼力，R^D=$ \dot{\mathit{r}} $[C^S(r)+ C^H(r)+ C^D($ \mathit{r} $,$ \dot{\mathit{r}} $)]，C^S(r)为内部结构阻尼，C^H(r)为势流阻尼矩阵，C^D($ \mathit{r} $,$ \dot{\mathit{r}} $)为指定的离散阻尼矩阵；R^S为内部结构作用力；R^E为外激力。

2 浮式风机系统参数及模型 2.1 浮式风机系统参数

考虑60 m水深确定10 MW级漂浮式风机系统参数，上部风机采用DTU 10 MW风机^[13]，浮式基础以O-Ostar基础^[14]为母型进行设计，系泊系统为3×1悬链线式系泊系统。

基础几何尺寸如图2所示。为了改进浮式风机基础水动力性能，减小漂浮式风机的横摇、纵摇响应幅值，在浮式风机基础内设置减摇舱室，共3个水舱。减摇舱室的底部距风机基础底部12.5 m，减摇舱室高为5 m，每个减摇舱室内填充水量占舱室体积50%。浮式基础水动力模型及减摇舱室布置示意如图3所示。

由图4可知，布置减摇水舱可有效减小浮式基础的摇摆响应幅值，其纵摇、横摇响应RAO最大分别下降了43%和61%。

浮式风机系统总体参数如表1所示。其中重心坐标所在坐标系均以水线面中心为原点。后续计算中，依据表中参数，对浮式风机进行水动力计算。

本文利用配重块改进浮式风机的系泊系统，在保障系泊性能的前提下，减少了系泊系统的占地面积，2种系泊方案具体参数如表2所示。

2种系泊方案均采用三缆分布形式，其分布形式如图5(a)所示，系泊方案2所占面积为系泊方案1所占面积的4倍。此外，方案2中的系泊缆为均匀缆，未设置配重块；方案1中的系泊缆设置有配重块^[15]，布置于系缆中部，如图5(b)所示。

针对此2种不同的系泊方案，进行额定、极限工况的时域计算，对比其系泊缆受力，额定、极限工况环境参数如表3所示。

额定、极限工况下，风浪方向均为0°方向，故对比系泊性能时仅考虑与风浪同向的1号系泊缆受力。

由图6可知，系泊方案1和系泊方案2的时域曲线规律一致，其均值、最值无明显差异。针对本文中设计的系泊缆参数，利用CCS《材料与焊接规范》^[16]中的估算公式对系泊缆的破断张力进行估算，本文系泊缆的等级对应规范中的R3：

$ F=0.0223k，$

(7)

$ k={d}^{2}\cdot(44-0.08d)。$

(8)

式中：F为系泊缆的破断张力，kN；d为系泊缆直径，mm。

计算可得系泊缆破断张力为13828 kN，额定工况下的1号系缆最大系泊张力为5201.02 kN，仅达到破断张力的36.9%，结合CCS《海上浮动设施入籍规范》^[17]中系泊索设计中对张力安全系数F的规定：

$ F=\frac{{P}_{B}}{{T}_{{m}}}。$

(9)

式中：P_B为系泊缆的最小额定拉断强度，kN；T_m为系泊缆最大张力，kN，完整自存工况下动力安全系数规定为1.67，本文设计系泊缆的安全系数为2.66，符合规范设计要求。

故采用了配重块的系泊系统不但可以满足浮体系泊要求，而且有效减小了系泊系统的占地面积。

2.2 浮式风机耦合分析方法及模型设置

设定3种耦合计算方法，即全耦合、限制性强耦合、限制性弱耦合方法。其中，全耦合方法即气动-水动-弹性-伺服全耦合计算^[3]。相对于全耦合，限制性耦合方法不建立上部风机模型，用塔顶的输入载荷替代上部风机作用并计算浮式基础运动。根据载荷的来源又可分为限制性强耦合和限制性弱耦合。其中，限制性强耦合的上部输入载荷来自于全耦合计算结果；限制性弱耦合的上部载荷来自于相同单机容量的固定式风机的计算结果。

为实现不同的耦合功能，建立了3个计算模型，包括全耦合模型、固定式风机模型和带塔筒的浮式基础，如图7所示，限制性耦合模型的排水量、重心高度等浮体参数与全耦合模型的设置完全相同。

限制性耦合计算过程如图8所示。对于限制性强耦合模型，先由全耦合模型，进行时域计算并提取塔顶六自由度受力，然后将塔顶六自由度受力施加到塔顶位置。对于限制性弱耦合模型，先由固定式风机模型，进行时域计算并提取塔顶六自由度受力，然后将固定式风机模型塔顶六自由度受力施加到塔顶位置。

3 计算结果及分析 3.1 自由衰减结果

首先对浮式风机系统进行静水自由衰减计算以分析风机的固有特性。设定浮式风机在六自由度下有一定的初始位移，释放浮体自由衰减，记录时域曲线，对衰减数据进行傅里叶变化即可得到浮式风机各自由度运动的固有频率。

风机六自由度自由衰减曲线见图9，六自由度固有周期结果见表4。可知，六自由度下固有周期最小的运动为垂荡，其固有周期为20.7 s，而波浪能量集中的周期范围多小于20 s，据文献[18]，周期小于20 s的波浪占比大于98%。故本风机六自由度的固有周期均明显远离此波浪范围，风机具有较好稳性，不易发生共振。

3.2 计算工况

浮式风机受到的环境载荷主要包括风、浪、流载荷3种，对于风载荷考虑湍流风且风谱采用海洋工程领域传统的NPD谱（又称ISO 19901-1风谱）；波浪载荷考虑了随机特性且采用JONSWAP波浪谱；流模型为剪切流，分为3层，包括海面海流、中层海流、近泥面海流，依次对应水深0.5、36、60 m，表5中每3个海况采用同样的环境载荷参数，分别对应全耦合，限制性强耦合和限制性弱耦合，工况参数如表5所示，考虑了3种典作业工况，即低风速工况、额定风速工况、高风速工况。此外，风浪流载荷方向一致，均为0°方向。

基于2.2节中的计算模型，根据表4进行时域计算，对响应结果进行分析，对比3种不同计算方法的结果差异。因风浪流方向一致，因此，浮式基础运动仅讨论纵荡、垂荡、纵摇运动，系缆力对比仅考虑1号缆和2号缆。

3.3 浮式基础响应结果对比

图10所示为额定工况下3种不同方法得到的浮式基础纵荡时域曲线和频谱图。图11对比了9个工况下的纵荡响应统计结果，其中柱形图表示纵荡均值，散点位置代表纵荡极值，误差棒高度表示标准差，后续统计图的符号含义与此相同。

从图10(a)可以看出，限制性强耦合结果略高于全耦合计算结果，而限制性弱耦合计算结果则低于全耦合计算结果。具体而言，限制性强耦合与全耦合计算结果差异较小，其极值、均值和标准差的差异均在10%以内；而限制性弱耦合与全耦合结果的差异均超过10%，尤其在工况1、工况2下，标准差分别相差41%和32%。在3个作业工况，浮式风机在额定工况下的纵荡响应最为显著。此外，图10(b)的频谱图表明，对于工况2不同耦合方法的频率组成成分一致，主要包括低频纵荡响应，同时激发了纵荡方向的固有频率成分和由波浪载荷引起的波浪频率成分，但在低频及纵荡固有频率下，限制性弱耦合响应幅值明显小于全耦合及限制性强耦合的响应幅值。

由图12可知，全耦合和限制性强耦合方法得到的垂荡均值，最值，标准差的差异微乎其微，而限制性弱耦合方法得到的垂荡均值，最值，标准差与全耦合垂荡结果差异较大。尽管不同耦合方法垂荡的标准差差异均在5%以内，但由于各方法垂荡均值结果存在差异，各耦合计算方式得到的垂荡最值会呈现出相应的变化。

由图13可知，全耦合方法与限制性强耦合方法得到的纵摇基本一致，其最值、均值、标准差的差异均在2%以内；但限制性弱耦合方法与全耦合方法的计算结果相差较大，特别是纵摇均值结果差异在20%以上。

3.4 系泊缆张力响应结果

由图14可知，该缆的张力变化规律与纵荡基本一致，限制性强耦合计算结果较全耦合结果稍大但二者差异较小约5%，而限制性弱耦合计算结果较全耦合更小且二者差别较大，最大差别约12%。不同耦合方法的频率成分一致，包括低频响应、纵荡固有频率成分、波浪载荷激励引起的波频运动响应、风力机叶轮的转动频率成分（1P），全耦合、限制性强耦合在低频、纵荡固有频率、1P下幅值明显高于限制性弱耦合。

由图15可知，工况1、工况2下，限制性弱耦合结果的最值、均值最大，全耦合结果次之，限制性强耦合结果最小。浮式风机由1、2、3号系泊缆共同组成系泊系统并提供恢复力，在浮式风机基础不发生横向运动时，其系缆力响应主要受基础的纵荡影响，前文结果显示，限制性强耦合计算得到的纵荡大对应1号缆响应值偏大，则2号缆松弛，系泊缆张力也偏小。在工况3下，全耦合模型的2号缆最值较限制性强、弱耦合更大，差异约10%。

4 结　语

1）设计了适用于中等水深的10 MW浮式风机基础和系泊系统，浮式基础中增加了减摇水舱，有效地减小了纵摇、横摇运动；通过在系泊缆上增加配重块，既保证系泊缆的整体性能，又较好地控制了系泊缆的整体长度，减小系泊半径，避免与其他设备干扰。

2）提出将限制性耦合细分为限制性强耦合与限制性弱耦合，对比3种不同耦合方法之间的差异。结果表明，限制性强耦合得到的纵荡与全耦合计算结果差异在10%以内，垂荡差异在5%以内，纵摇差异在2%以内，1号缆系缆力差异在5%以内。

3）对于关键响应参数，限制性强耦合计算结果稍大于全耦合计算结果，且结果可靠，采用风机厂家提供的全耦合塔顶载荷作为输入数据进行浮式基础设计可以获得较为准确且可靠的结果。