极地船舶冰载荷研究方法综述

2023-10-18许育文顾学康赵南王艺陶

装备环境工程 2023年9期

许育文，顾学康，赵南，王艺陶

（中国船舶科学研究中心，江苏无锡 214082）

北极地区拥有丰富的油气资源，是未来各国贸易往来的重要航线。根据美国地质调查局的评估，全球约有13%的石油和30%的天然气储存在北极圈[1]，而近年来随着全球气候变暖，北极地区的冰川大量融化，因而北极航道通行能力也在逐年增加，对于北极地区的探索则越加迫切。在北极地区，船舶及海洋结构物和海冰的相互作用，将会很大程度地影响极地装备的工作性能和安全性能，甚至在严重情况下会导致结构物沉没[2]，如图1 所示。随着极地资源的开发，船舶在极地地区的航行也日益频繁，从而各国也开始逐渐重视和发展船冰相互作用的研究。我国也加大了对于极地研究的投入，涌现了一些杰出的学者，形成了大量有意义的成果。

图1 极地地区航行的冰损事故[2]Fig.1 Ice damage accidents in polar regions[2]

在极地区域航行过程中，船舶与海冰的碰撞无法避免，因此船冰的相互作用对于冰区船舶安全性设计有着重大影响。如何快速准确地分析预测船舶与海洋结构物在极地运营过程中所受的冰载荷，成为了国内外学者研究的重点方向之一。

由于海冰自身物理力学性质非常复杂，与船舶的相互作用过程更加复杂，因此在研究船冰相互作用时的结构冰载荷需要解决诸多问题。首先是对海冰材料特性的研究，海冰主要分为浮冰（漂浮在海面的平整海冰）以及冰山冰（从冰川或冰架上脱落后漂浮的巨大冰块）2 类[2]，如图2 所示。由于海冰在实际形成过程中，其晶粒尺寸、晶体取向、孔隙率、含盐量、温度以及冰的尺度存在着差异，导致海冰的压缩、拉伸以及剪切等力学特性也存在着不同程度的变化，而力学特性对海冰的失效模式有较大影响，决定着冰载荷的大小和形式[3]。因此，研究冰的物理性质、力学特性以及本构关系是冰载荷研究的基础。

图2 海冰类型[2]Fig.2 Type of sea ice[2]:a) floating ice (crushed ice zone);b) iceberg

由于最初研究冰载荷时缺乏实地的观测，以及对海冰物理性质认识的不足，学者们通常采用理论分析与经验相结合的办法，即半经验公式方法对冰载荷进行预测，但其适用范围有限，并且公式的可靠性依赖于对海冰物理力学性质的了解。概率法主要针对船-冰相互作用过程中各种条件的不确定性进行建模，从而对自然环境中冰载荷的随机性进行合理的描述。试验方法主要有模型的冰水池试验以及实船现场测试，其不仅可以直观地观察船冰相互作用的效果、冰层的失效过程，还能掌握实时数据进行载荷分析，但是海冰等模型的制备以及试验所需的环境条件较高，花费的成本很大。随着近些年计算机技术的迅速发展，利用数值仿真的方法模拟船舶在冰区的运动已经成为国内外学者研究的热点方向之一，其主要的研究方法包括有限元法、离散元法等。有限元法对于船冰作用的非线性过程具有良好的模拟效果，但是对于海冰的细观结构以及海冰破坏过程中破碎离散状态的描述具有一定的困难。离散元模型可以采用精细离散化方法对大体积冰进行建模，较好地模拟海冰与结构作用时的破碎规律，但是其开发状态还不够完善，需要合理的参数确定方法，以及高水平的计算效率。

本文主要对船舶及海洋结构物的冰载荷研究方法进行介绍，从理论分析、试验方法以及仿真模拟3个方面入手，总结并展望今后极地船舶冰载荷研究的方向和难点问题。

1 理论分析方法

理论分析方法包括直接计算法以及概率法。在极地船舶冰载荷研究的初期，由于缺少实测数据和对船冰相互作用过程认识不足，研究通常是将理性的理论分析和由经验获得的信息相结合，构建船冰相互作用的理论模型，从而进行冰载荷的直接计算。在后期冰载荷研究的发展过程中，学者们考虑在极地的复杂环境下，海冰的材料性质、几何形状、失效模式以及船冰的作用过程都具有很强的随机性，因此在冰载荷的理论研究中引入概率方法是很合适的。

1.1 直接计算法

学者们简化了海冰的本构方程以及船体板壳等结构的强度理论，结合船冰作用过程中海冰失效模式、结构几何形状和相互作用环境等因素，并依据船冰作用时的试验数据，采用数学推导的方法求得冰载荷。考虑到计算角度和应用重点的不同，可以将极地船舶冰载荷分为全局（Global）角度下船舶在冰区航行时的总冰阻力与局部（Local）角度下不同类型的船体结构受到冰体作用下的冰载荷。

1.1.1 船舶总体冰阻力

早在1956 年，Johnson[4]便给出了船舶在冰区航行时计算船体总冰阻力的简化公式，但是由于当时条件所限，难以得到可靠的海冰物理数据，无法获得精确的冰阻力。Milano[5]假定船体破冰的所有能量都是用来抵消冰阻力的消耗，从而获得了忽略压沉阻力和滑动阻力影响的冰阻力计算公式，但也因此导致计算结果相较于实际情况偏低。Lindqvist[6]则在船舶航行的冰载荷公式中考虑了船型系数、海冰厚度与强度、船冰摩擦系数、船舶航速等因素的影响。Lewis 等[7]对以前研究者们的成果进行了总结，结合实船试验结果推导出了破冰船破冰时的阻力计算公式：

式中：R为冰阻力；h为冰厚；ρw为海水密度；g为重力加速度；B为船宽；σ为冰的弯曲强度；Fn为冰厚的傅氏数，Fn=v/(gh)0.5，其中v为船舶航行速度。

刘源等[8]结合上述的理论公式计算比较了挪威海岸警卫队破冰船破冰过程中的冰阻力，并据此分析了各公式的计算结果、各种方法的优缺点及适用范围。徐义刚等[9]基于冰载荷的解析算法分析了航速、冰厚在不同海冰失效模式下对集装箱船受到总体冰载荷的影响。

1.1.2 局部结构冰载荷

船冰相互作用的形式与海冰特性、结构形状等因素关系紧密，因而研究平整冰层以及碎冰层与不同结构类型相互作用时的冰载荷是近些年来的热点问题之一。其中主要的结构类型有2 种，分别是直面结构和斜面结构。海冰和直面结构作用时，海冰的破坏形式以挤压破碎和屈曲失效为主；海冰和斜面结构作用时，海冰的破坏形式以弯曲断裂为主，进而使得碎冰产生滑移、堆积。

Korzhavin[10]给出了平整冰层与直面结构之间接触时直面结构冰载荷大小的计算公式：

式中：F为结构总体冰载荷；I为局部系数，依据结构宽窄取值；K为接触条件系数；D为结构迎冰面宽度；H为冰层厚度；m为结构形状系数，不同形状的直面结构取值不同（如圆形柱和多边形柱）；σc为海冰单轴抗压强度；v和v0分别为冰速和参考冰速。

随后Croasdale 等[11]和Shi[12]利用一系列海冰和孤立桩作用的试验数据进行了验证，同时也给出了Korzhavin[10]的经验公式中系数的推荐值。ISO 19906[13]中针对直立结构所计算的平均冰载荷为：

式中：PG是结构平均冰载荷；ω为冰与结构的接触宽度；H为冰厚；h1为参考冰厚；CR为海冰强度系数，北极地区建议取值2.8 MPa；β为接触宽度效应的经验系数，取值为–0.16；n为冰厚效应的经验系数，当冰厚小于1 m 时，取值为–0.5+h/5。

除此以外，美国石油工业协会给出的API 规范[14]以及法国船级社给出的冰载荷计算指南[15]进一步总结了冰层和直面结构作用时的冰载荷计算流程。

Croasdale 等[16]分别分析了完整冰层和破碎冰层与斜面结构的相互作用过程，构建了预测冰层和斜面结构相互作用时冰载荷大小的理论模型。Danys 等[17]提出将冰层理想化为弹性板来处理，由最大拉伸应力来确定冰层是否失效，得到一个用于预测圆锥结构冰载荷大小的理论公式。Ralston[18]将冰层考虑成弹塑性板，利用塑性上限理论研究了冰层和圆锥结构作用时的冰载荷大小。唐百强等[19]基于Ralston 提出的海冰本构模型，构造了平压头接触问题的应力场分布，并依据冰的塑性应力场得出了海冰的塑性极限压力。ISO 19906[13]同样给出了2 种方法计算斜面结构的局部冰载荷，分别是基于Croasdale 的弹性梁分析模型与基于Ralston 的三维塑性分析模型。王国军等[20]依据实测的冰厚、冰速以及海冰破碎尺度等数据，得出了修正的椎体冰载荷的确定性冰力函数。隋吉学[21]采用理想的各向同性弹性材料模拟海冰，假定冰体的失效由弯曲强度主导，提出了冰体与斜平面结构物作用的三维计算模型。

随着冰载荷计算方法的发展，研究者在前人的基础上，结合新的理论方法对冰载荷的计算方法进行了拓展。胡嘉骏等[22]应用内薄层法的基础理论，提出了船体局部结构冰载荷的计算公式，并描述了碰撞结束后冰力下降的卸载过程。王文瑜等[23]利用解析法估算出船冰碰撞中舷侧3 个典型碰撞位置结构变形能，从而得到计及结构变形能的冰载荷计算优化方法。龚榆峰[24]基于损伤力学，建立了计算船体结构受到局部冰载荷的方法，并与基于断裂力学的船体总体冰阻力计算模型进行了对比，提出了船体结构安全的简化评估方法。

1.2 概率法

由于极地区域环境复杂，船体与冰的相互作用具有很强的随机性，因此相比于确定性方法，从统计学的角度出发，利用概率方法去描述冰载荷的随机性，并构建概率模型对冰载荷进行评估分析是更加合适的。

针对船舶总体冰载荷的计算，Kheysin[25]首次使用了概率方法，通过Likhomanov[26]的测量结果，使用泊松分布来表示任意时间间隔内的撞击数量。Kujala 等[27]建立了冰载荷的单参数指数模型，但是其适用范围有限，因此Suyuthi 等[28]将2 种指数分布按照一定比例因子结合构建了三参数指数的模型：

其中0

Taylor 等[29]基于缺陷结构和接触条件的随机性对断裂过程的影响，建立了概率断裂力学模型，模拟冰样的局部断裂。通过应力分析确定冰样的应力状态，然后考虑裂纹尺寸、方向和密度等因素的影响，进行失效概率计算。概率断裂力学方法为改进冰破坏建模提供了一种可靠的选择，其模型主要考虑的是导致脆性断裂的载荷加载速率，未来可以考虑在慢加载速率或者持续加载下随时间变化的裂纹扩展。Sinsabvarodom 等[30]利用蒙特卡罗模拟技术，根据ISO 标准[13]中的公式，对冰-结构相互作用引起的荷载进行了概率评估，针对垂直结构和倾斜结构，估算了冰强度系数以及冰抗弯强度与冰厚度的相关系数，研究了相关系数的变化对冰荷载不确定性的影响。

1.3 小结

综上所述，很多机构和学者都给出了对于结构冰载荷的计算方法，直接计算法能够经济迅速地给出船舶冰区载荷的结果，还能对未来船舶的载荷分析进行一定预测与评估，适用于船舶初步设计阶段。海冰的物理性质与力学性质非常复杂，虽然学者们建立了大量的海冰本构模型，可理论推导过程仍然难以全部涵盖海冰物理力学性质以及结构破坏模式，并且在自然水域，船冰的相互作用被很强的随机性影响，因而概率法可以对分析结果的稳定性和准确性进行进一步优化。但由于概率法依赖于大量的实测数据，对于数据的准确性和实效性要求很高，并且对于特殊情况下的结构冰载荷计算，概率法就不完全合适了。因此，依托于实测数据，构建稳定、准确、应用场景广泛的冰载荷计算模型是未来的重点研究方向之一。

2 试验方法

船-冰相互作用是非线性动态响应过程，仅采用理论推导难以准确描述复杂作用过程中的冰载荷，上述冰载荷理论分析多用来验证其他方法的有效性。在对结构冰载荷的研究中，试验方法是非常有效的手段，并且试验方法还可以提供大量的数据以供理论分析和数值方法参考。试验方法的主要研究方向分为2种：一种是针对海冰的材料性质，分析抗压强度、弯曲强度等强度指标以及冰的变形破坏特点，主要涉及单轴压缩试验、三轴压缩试验、弯曲试验和压痕试验等。另一种是直接进行实体的船冰试验测量冰载荷，包括室内模型试验和实尺度测试。由于本文主要研究的是对冰载荷的测量方法，所以主要介绍的是直接测量冰载荷的试验。

2.1 室内模型试验

室内模型试验通过建立模型冰来模拟海冰工况，在试验过程中，需要模拟实际情况中不同的极地海洋环境以及海冰几何构造，测量结构的冰载荷响应数据，并通过采集的测试数据推导出实尺度下的结构冰载荷。因此，拥有合适的试验环境和试验设备便显得十分重要，我国也正在加大对于室内冰水池建设的投入。

对于室内试验来说，首先要考虑模型冰的制备，Zong 等[31]利用聚丙烯蜡模拟浮冰，分析了海冰的几何尺度和物理性质对结构冰载荷的影响。Barker 等[32]为了研究海冰与圆锥结构相互作用时的结构冰载荷，采用了7 种构型的模型结构进行试验，并通过改变冰速，冰厚、结构形状、结构上的水位、结构刚度和结构固有频率等参数，来分析其对模型所受冰载荷的影响。Cui 等[33]在辽东湾滨海咸水湖取冰，开展了冰块冲击载荷的试验研究，试验得到了不同冲击头、不同冲击速度以及不同尺寸冰样对冰冲击荷载的影响，如图3 所示。Guo 等[34]利用石蜡为原材料的非冻结模型冰，开展了船舶在不同密集度工况下碎冰区航行时的载荷试验。然而，他们都未考虑在冲击过程中“尺寸效应”的影响，而模型试验的相似是要重点分析的问题，将来也可以进一步讨论冲击角度、曲面曲率等关键船-冰接触参数对冰载荷的影响。

图3 冲击载荷对比[33]Fig.3 Comparison of impact load[33]:a) different speed;b) impact heads in different shapes

船舶在浮冰区和平整冰区的行驶过程中，对船体所受的总体冰载荷进行测量是模型试验的重点。黄焱等[35-37]利用天津大学低温水池开展了大型运输船在浮冰区和平整冰区的冰阻力试验研究，获得了船体在碎冰航道中航行时碎冰块的运动过程，以及船体在平整冰区航行时船体结构周围的海冰失效形式。此外，他们还对船首与浮冰的碰撞过程进行了模拟试验，利用船首表面的触觉传感器获得船冰碰撞的时变载荷，如图4 所示。在这些试验中得到的参数数据可以作为数值模拟的输入，从而验证模拟过程和理论分析的准确性，但是得考虑好试验时的尺寸和试验条件。

图4 整体碰撞冰载荷的时空演变历程[37]Fig.4 Temporal and spatial variation of the global ice impact load[37]

同样，对船体局部结构所受冰载荷的模型分析也是必不可少的，Riska 等[38]和Aksne[39]开展了室内模型试验，测量了船体不同结构与平整冰相互作用的局部冰载荷，并对比验证了数值分析的结果。Li 等[40]利用模型船试验观察了冰破碎过程和碎冰堆积，测量了载荷大小和碎冰尺寸，分析了船速、航向角和冰厚对冰屑形成和堆积的影响，并研究了船舶的局部结构在平整冰区航行时不同航速、航向角和冰厚条件下的冰载荷。Polach 等[41]研究了船舶垂荡和纵摇运动对船体冰载荷的影响。

2.2 实尺度实测

由于室内环境与真实环境差距较大，室内模型试验的推导结果具有局限性，因此在实际环境中的实尺度实测便在结构冰载荷的研究中必不可少。实船现场测试通过在关键区域布置传感器，从而测量出冰激响应，之后通过载荷识别技术对数据进行统计分析，进一步获得冰载荷时程，得到真实环境下冰载荷时历。在冰载荷的现场测量领域，测力传感器可以直接测量结构冰载荷。

Johnston 等[42]利用MOTAN 惯性测量系统，得到了船体总体冰荷载的观测曲线。Timco[43]通过实地试验对固定结构物与浮冰的碰撞载荷进行了总结，得到了一个浮冰与结构物碰撞的载荷计算公式。然而应用更广泛的是通过结构应变来对船体局部冰载荷进行测量，即利用载荷识别算法对结构冰载荷进行反演，得到所需要的载荷数据。Leira 等[44]将KV Svalbard号上的载荷与实测应变联系起来，通过实测应变对作用在KV Svalbard 船体上的冰载荷进行了估计，并且通过时间序列预测方法来预报结构应力和载荷。但该预测方法需要冰荷载有关的参数是恒定的，并且对于预测载荷需要全面的校准。Suominen 等[45]利用极地船舶Agulhas II 在波罗的海进行海冰试验，通过剪切应变仪对船尾肩部开展了全尺寸测量，结合统计学模型获得了不同长度冰荷载的发生频率及大小的发生概率。

不同载荷长度下冰载荷分布的拟合如图5 所示。可以看出，Weibull 概率分布模型对实测载荷的拟合效果最好，但是指数分布和对数正态分布在某些情况下可以接受，即短加载情况下载荷分布更符合指数分布，较长加载情况下分布形状更符合对数正态分布。Kim 等[46-47]对破冰船“ARAON”进行了海冰现场试验，将破冰船在浮冰区和层冰区的试验结果与模型试验结果进行了比较分析。Jeon 等[48]用安装在船体内板和船艏横架上的应变片测量了作用在“ARAON”船艏的局部冰载荷，并对由剪切应变数据估计的局部冰载荷和由船体板应力所估计的局部冰载荷进行了对比，认为通过剪切应变估计冰载荷不仅同样精确，仪器成本和计算程序也更加简单。

图5 不同载荷长度下冰载荷分布拟合[45]Fig.5 Fitting of ice load distribution under different load length[45]

在实船试验的结构冰载荷监测系统中，载荷识别算法作为监测系统的核心，近年来被广泛关注。现在主流的冰载荷识别方法为影响系数矩阵法、时域反卷积算法和支持向量机法，其中影响系数矩阵法的优势在于能够反映出船冰碰撞的时程特征，但是难以反映冰载荷的动力效应；时域反卷积算法则能够解决这个问题，只是结果稳定性和快速性不足；支持向量机算法对于非线性问题有很好的应用效果，并且机器学习算法的发展具有更大的前景，能够为适应更加复杂环境的载荷识别提供更加有效的选择。大连理工大学项目组[49-51]分别基于以上3 种算法建立了冰载荷反演识别模型，在我国“雪龙号”破冰船的船艏区域安装应变片（如图6 所示），并在“雪龙号”执行极地航行任务时完成了试验测量，通过反演识别模型研究了船体破冰时载荷的时程特征以及分布规律。何帅康等[52]对内河轻型破冰船的冰载荷进行了监测，利用影响系数矩阵法对监测区域的冰压力进行了识别，并将河冰冰载荷特性与海冰进行对比。孔帅等[53-54]基于形函数的载荷识别算法建立了船体结构动冰载荷快速识别模型，该模型不仅考虑了冰载荷的动力效应，并且保障了求解快速性。

图6 船体结构海冰载荷测量的应变位置及现场测量示意图[49]Fig.6 Location of strain sensors for ice load on ship hull and schematic diagram of field measurement[49]

2.3 小结

我国在冰区船舶的模型试验方面目前还处于起步阶段，包括室内冰水池的建设、冰区船舶的模型试验技术和模型冰的精确制备技术也需要进一步的发展，并且国内尚未形成系统的制备模型冰的流程，将模型冰的制备标准化，并测定在标准化流程下的模型冰的材料参数，是未来研究的重点方向之一。相比于室内模型试验，实尺度实测更能直接准确地提供冰载荷数据，但是现场环境极为复杂，很容易出现难以控制的天气与海冰冰情，并且试验成本较高，可重复性较低，这些都对获得普适性的冰载荷模型造成很大的困难。

3 仿真模拟

随着计算机技术的快速发展，仿真模拟方法在船冰相互作用中的应用也愈发完善。仿真模拟方法通过简化船舶与海冰之间的相互作用过程，建立船冰的相互作用模型，并运用数值方法计算出船舶所受到的冰载荷。现阶段主流的仿真模拟方法为有限元法、离散元法以及各类仿真方法的耦合应用。

3.1 有限元方法

有限元法是应用极广泛的数值方法之一，它可以处理相对复杂的几何构型、力学物理模型和边界条件，同时其理论基础上的可靠性和计算的高效性也是值得信任的。但对于船冰相互作用的数值模拟，建立合适的海冰模型是非常重要的，特别是建立合适的海冰本构模型。同时，明确海冰的失效准则也是有限元方法的关键，包括屈曲、剥落、压碎、蠕变失效、辐射裂纹和轴向断裂等失效过程。本文介绍的有限元法还包括基于传统有限元法优化的内聚力单元法以及扩展有限元法。

3.1.1 传统有限元法

Aksnes[39]基于弹性梁理论建立了海冰的一维有限元模型，其失效是由最大截面弯矩控制，据此模拟了船冰的相互作用。Ehler 等[55]则采用塑性模型和单元删除方法模拟了冰梁的弯曲失效，但是未考虑海冰的破碎和弯折等失效模式。Su 等[56]在考虑了上述失效模式后，提出了一个用于模拟船舶连续破冰情况下船舶整体和局部海冰的数值模型，得到了船舶在平整冰区航行中随时间变化的冰载荷。Sazidy[57]应用塑性板模型和单元删除法模拟了楔形冰板的弯断，同时分析了在楔形冰失效过程中应变率的影响。Lubbad 等[58]基于弹性基础板理论，将海冰当作各向同性的弹脆性材料开发了一套仿真程序，认为船舶受到的冰载荷主要是由于海冰失效以及海冰运动所导致的，并使用Phyx引擎模拟了破冰船破冰过程中的总体响应过程。与Lubbad 等[58]不同的是，Kolari[59-61]在描述冰的脆性失效时建立了一个三维各向异性损伤模型，利用该模型对柱状冰进行了双向加载的数值模拟，得到了相应的失效模式和强度。柱状冰脆性破坏模式如图7 所示，在单轴压缩条件下，沿加载方向发生劈裂破坏；在中低约束条件下，加载面上的破坏模式为剪切型断裂；在更高的约束下，柱形冰的失效模式为剥落失效。

同样，徐莹[62]利用有限元方法模拟了Pond Inlet冰山压痕试验，将海冰材料设定为黏弹塑性材料，在冰山碰撞过程中删除压溃单元，并将数值模拟所得的冰载荷和压力-面积曲线与试验结果进行对比，验证了数值方法对压溃冰载荷预报的有效性和准确性。王刚等[63]通过有限元模型研究了冰层和锥体结构相互作用时，接触宽度对于冰层弯曲破坏模式的影响。Wong[64]通过三维有限元模型模拟了冰层的弯曲断裂，并分析了圆锥形桥墩与冰层的相互作用。Liu 等[65]基于连续介质力学，建立了船-冰山碰撞的综合弹塑性模型，来模拟冰山与结构物的相互作用，并研究了接触面压力的演变规律、船艏损伤分布和冰山破碎情况。王陈阳等[66]以极地物探船为对象，在LS-DYNA软件中建立了物探船-水-冰耦合有限元模型，获得了冰载荷曲线，并据此对该极地物探船的局部结构进行了冰致疲劳损伤评估。常顺[67]建立了船-冰碰撞的2种有限元模型，分别是基于流固耦合方法和基于浮冰载荷二次开发的简化方法的模型，并对比分析了2 种方法模拟的结果。岳前进等[68-70]采用有限元方法分析了冰与锥形结构物相互作用的破坏过程，JZ20-2NW平台的抗冰性能以及海洋结构物的冰致疲劳问题。Kim 等[71-72]考虑了冰-水、冰-结构和冰-冰相互作用对于船舶冰区航行的影响，采用有限元模拟计算出船-冰接触压力和水动力，并据此研究了由冰载荷引起的船体结构疲劳损伤。

3.1.2 内聚力有限元法

内聚力单元法（Cohesive Element Method）近年来在海洋工程方面得到了广泛应用，该方法基于弹塑性断裂力学及传统有限元法，采用厚度接近0 的单元表示为内聚力单元，将其插入传统单元之间，并将材料损伤全部聚集在内聚力单元中，当内聚力单元失效时，使得传统单元分离，从而模拟裂纹的产生和扩展。Gribanov 等[73]利用具体试验数据校准材料特性和模型参数，采用内聚力模型来模拟冰圆柱试样在单轴加载条件下的断裂行为（如图8 所示），获得了应力-应变和损伤累积曲线。Gürtner 等[74-75]利用内聚力单元法来处理冰的特征断裂，其中内聚单元放置在冰的有限元网格之间，建立了模拟冰-碎石相互作用的计算模型。之后基于内聚力单元法建立了灯塔有限元模型，模拟灯塔与海冰相互作用，并将实测的灯塔-海冰作用过程与基于内聚力单元法的模拟过程进行了比较，验证了内聚力单元法的适用性。

Lu 等[76]建立了冰层与斜平面结构的内聚力有限元模型，并分析在了不同网格尺寸以及单元形状情况下海冰裂纹扩展路径的改变。此外，Lu 等[77]分别采用单元删除法、内聚力单元法（CEM）、扩展有限元法（XFEM）模拟了海冰的弯曲变形过程，并讨论了几种方法的优缺点。Lu 等认为，单元删除法效率最高，但模拟结果的准确性取决于模型的完整度和精细度；CEM 可使用相对简洁的冰体单元材料模型，但插入内聚力单元会改变结构的整体刚度；XFEM 可以较好地模拟裂纹初始状态和扩展，但该方法目前仍不成熟，需要进一步的发展。

3.1.3 扩展有限元法

扩展有限元法（Extended Finite Element Method）对有限元法进行了改进，结合断裂力学的相关理论，令裂纹能够穿越单元内部，从而模拟裂纹的发生和扩展过程。Lu 等[78]建立了一种模拟冰-结构相互作用过程中浮冰断裂的模型，针对劈裂破坏模式，提出了一种基于扩展有限元法的模拟劈裂裂纹扩展的数值方法。模拟的裂纹扩展路径与现实比对如图9 所示。可以看出，裂纹扩展轨迹在视觉上与实际情况很相似，但是模拟裂缝路径与现场观测结果的精确匹配仍然是一大挑战。最后Lu 利用实测浮冰的几何数据，可以数值预测非线性劈裂裂纹的扩展路径。

图9 浮冰裂纹路径模拟与现实对比[78]Fig.9 Comparison between simulated and real floating ice crack path[78]

Ji 等[79]通过扩展动态有限元建立了冰-结构相互作用模型，从反馈机制和能量机制2 个方面分析了冰-结构相互作用的物理机理。Xu 等[80]基于扩展有限元法，采用内聚区概念描述了层冰断裂以及裂纹扩展过程，并通过模拟波罗的海登陆艇艏与层冰碰撞的现场试验，验证了数值方案的正确性。

综上所述，有限元能够较好地处理比较复杂的力学模型以及几何非线性问题，但是冰的破碎具有离散特性，很难精确地处理海冰破坏时的离散状态。有限元法在模拟船舶结构的动态响应具有一定的优势，但在模拟介质间断和离散方面有一定的困难，因此如何使用有限元模拟海冰的失效模式有一定的挑战。

3.2 离散元方法

有限元方法虽然能很好地描述结构的动态响应，但是很难准确模拟海冰的破碎失效过程。离散元法具有描述黏结和失效的功能，能对海冰破坏、重叠和堆积等过程进行有效的模拟，且该方法已经在岩石、陶瓷、玻璃等脆性材料的破碎过程分析中得到了大量应用，因此研究者们基于离散元方法也开展了一定的研究。

Løset[81]建立了二维离散圆盘的浮碎冰模型，并基于此对其进行了动力学的数值模拟。Hanse 等[82]则利用二维离散元模型对浮冰区的船舶响应进行了模拟计算。杨斌等[83]通过构建海冰黏弹性模型，计算了船舶航行过程中受到的总体碎冰载荷。

海冰的离散单元模型可分为颗粒离散单元和块体离散单元2 种，这2 种离散单元的构建如图10 所示。在颗粒离散单元模型中，海冰依据冰晶结构特点由多个球形颗粒黏接而成，由颗粒材料性质及单元间黏接强度描述海冰的破坏模式，块体离散单元模型则能更合适地模拟海冰重叠和堆积过程。季顺迎等[84]将海冰离散为若干个规则排列且具有黏接-破碎功能的颗粒单元，并据此进行了冰与海上直立平台的模拟计算，其作用过程如图11 所示。Bateman 等[85]利用海冰的离散模型模拟了在弯曲等载荷作用下海冰的断裂破坏。

图10 海冰的离散单元模型[84,86]Fig.10 Discrete element model of sea ice[84,86]:a) parallel bonding model for 2 units;b) ice blocks generated based on Voronoi cutting algorithm

图11 离散单元模拟的海冰和直立结构相互作用过程[84]Fig.11 Discrete element simulation of the interaction process between sea ice and vertical structures[84]

块状离散单元模型的应用则更加广泛。刘璐等[86]建立了基于闵可夫斯基和方法的扩展多面体单元模型，描述了非规则颗粒单元的几何形态，并采用Voronoi 算法获得了碎冰的初始随机分布状态，分析了冰块尺度和速度对圆桩上冰载荷的影响。Li 等[87]则通过三维圆盘单元模型模拟了碎冰的离散状态，分析了海冰间、海冰与船间作用力，并研究了航速和海冰参数的影响。Long 等[88]通过离散元法模拟了核电站取水口浮冰的堆积以及海冰与核动力浮式平台的相互作用过程，预判了海冰的堆积高度，从而降低了核电站取水口浮冰堆积的可能危险。

离散元方法在海冰与椎体结构相互作用的数值模拟上也取得了非常广泛的应用。刘圆等[89]在季顺迎等[84]的海冰离散模型基础上，计算了海上风机结构与平整冰层相互作用时结构的冰载荷时历。Wang 等[90]模拟分析了平整冰和碎冰对锥体风电结构的冰载荷以及冰激结构的振动响应。同样，针对锥体结构和海冰的作用过程，Long 等[91-92]构建了海冰强度和微观离散元参数的关系，确定了参数的合理取值范围，最后提出了锥体受到静态冰载荷的计算公式，并分析了锥体结构与海冰作用的破坏过程，认为倒锥体结构相比于正锥体结构具有更好的抗冰性能。针对锥体结构-冰的离散元研究，大连理工大学项目组进行了大量详细的工作，其成果不仅有效地证明了离散元方法应用在海冰与锥体结构作用过程中的可靠性，也为冰区海上风电单桩结构的抗冰锥设计和海洋工程结构的抗冰设计提供了参考依据。

3.3 方法耦合

随着数值模拟在冰载荷的研究中的发展，越来越多的研究人员考虑将不同的数值方法进行结合，来达到更加优化的结果。如通过有限法构建船体结构的数值模型，再通过离散元法分析海冰的破碎过程，将2种方法进行耦合，既解决了有限元法无法针对冰体破碎和冰堆积模拟的困境，同样解决了离散元模型在进行结构响应计算时不够成熟的问题。Lilja 等[93-94]利用了三维有限元-离散元模型来模拟层冰，由同向旋转、黏性阻尼的梁有限元组成的平面内梁格与刚性离散单元形成实际冰层质心相连，通过Voronoi 算法生成并网格化薄片，在考虑中心垂直受载和边缘垂直受载情况下进行了冰层挠度计算。Shao 等[95]采用离散元方法模拟了海冰的破碎过程，并通过有限元模型计算了海洋平台结构的动力响应，之后通过对海冰和海洋平台结构的相互作用进行耦合分析，建立了冰激海洋平台结构振动的离散元-有限元模型，并基于王帅霖等[96]提出的海冰与海洋平台结构在接触面上的接触算法，研究了平台的结构响应和海冰荷载特性。

不仅有FEM-DEM 模型，石玉云等[97]建立了开敞水域海洋平台计算流体动力学数值模型以及海冰的离散元模型，研究了流速和碎冰密集度对海洋平台结构冰载荷的影响。闫允鹤等[98]也利用了CFD-DEM 耦合方法研究了碎冰航道内上浮潜艇所承受的冰载荷特性。

3.4 小结

随着船舶结构冰载荷研究的深入，仿真模拟方法占据了愈发重要的地位，结合结构冰载荷的研究需求，有限元方法从传统方法发展了内聚力单元法和扩展有限元法，离散元法也与各个仿真方法进行了耦合分析。有限元法在分析结构失效、断裂等过程中有一定的优势，并且理论基础深厚，使用方便；离散元法在模拟冰的离散和破碎过程中更加合适，但对结构的变形和损伤的模拟则不如有限元法快速简洁。因而，将各个仿真模拟方法进行耦合，构造合适的船-冰相互作用的数值仿真模型具有深厚的潜力。