張 健，王凱民，許俊華，杜成忠

(1. 江蘇科技大學(xué) 船舶與海洋工程學(xué)院，江蘇 鎮(zhèn)江 212003；2. 江蘇新?lián)P子造船有限公司，江蘇 泰州 214500)

0 引 言

隨著全球平均氣溫的升高，冰層融化速度加快，北極資源開發(fā)和與之相關(guān)的北極航線建設(shè)成為當(dāng)前的一個(gè)熱點(diǎn)問(wèn)題。對(duì)我國(guó)而言，正在“一帶一路”倡議背景下，開展北極航道的開發(fā)和利用，打造“冰上絲綢之路”。然而，由于極地海域環(huán)境的復(fù)雜性，船舶與海冰發(fā)生碰撞不可避免，因此破冰船在保障極地資源開發(fā)和海洋活動(dòng)的安全性和可持續(xù)性方面起著至關(guān)重要的作用。

由于船舶破冰是一種復(fù)雜的非線性動(dòng)態(tài)響應(yīng)過(guò)程，再加上海冰本構(gòu)模型的復(fù)雜性，國(guó)內(nèi)外還沒(méi)有一種有效而準(zhǔn)確的冰阻力預(yù)報(bào)方法。冰阻力主要受冰材料特性、船體形狀、船舶運(yùn)動(dòng)和航速等多種因素的影響[1]。目前，計(jì)算冰阻力的方法主要有3 種：基于模型試驗(yàn)或全尺度試驗(yàn)數(shù)據(jù)生成阻力經(jīng)驗(yàn)關(guān)系的經(jīng)驗(yàn)公式[2-5]；基于船舶運(yùn)動(dòng)方程描述船舶與冰動(dòng)態(tài)作用過(guò)程的數(shù)值模擬方法[6-7]；模型或全尺度試驗(yàn)研究方法[8]。

1 冰材料及模型驗(yàn)證

準(zhǔn)確的海冰材料模型是保證數(shù)值模擬結(jié)果準(zhǔn)確度的關(guān)鍵，然而由于冰材料的復(fù)雜性，國(guó)內(nèi)外對(duì)冰材料本構(gòu)模型的研究尚未形成共識(shí)。本文選取的冰材料模型為*MAT_124 等向彈塑性材料[9]。該材料模型中需定義3 種不同的材料曲線：定義2 種不同的屈服應(yīng)力與有效塑性應(yīng)變曲線；一種用于拉伸，另一種用于壓縮；第3 條曲線是應(yīng)變速率與比例壓縮屈服應(yīng)力因子（CYSF）曲線，用于說(shuō)明應(yīng)變速率對(duì)壓縮強(qiáng)度的影響，將等效塑性應(yīng)變值作為失效標(biāo)準(zhǔn)。

表 1 冰材料參數(shù)Tab. 1 The material parameters of ice

圖 1 鋼板與冰錐模型Fig. 1 FE model of ice cone and steel plate

通過(guò)模擬鋼板擠壓圓錐冰體，將碰撞力結(jié)果和現(xiàn)有試驗(yàn)結(jié)果進(jìn)行比較，以驗(yàn)證所選冰材料的適用性。為與試驗(yàn)?zāi)Ｐ椭斜F形狀尺寸基本保持一致，并模擬冰錐固定在試驗(yàn)裝置上，建立直徑為10 cm，錐角為30°的冰錐有限元模型，將冰錐底平面節(jié)點(diǎn)進(jìn)行全約束。由于在鋼板與冰錐的擠壓試驗(yàn)中，鋼板的厚度較大，表面無(wú)明顯的變形，故將鋼板設(shè)為剛體，這樣不僅節(jié)約計(jì)算時(shí)間又不會(huì)產(chǎn)生較大誤差。有限元模型如圖1 所示。為驗(yàn)證該冰材料在不同應(yīng)變率下的適用性，使鋼板分別以1 mm/s 和100 mm/s 的速度擠壓冰錐，圖2 為鋼板以不同速度與冰錐相互作用的冰力-位移曲線，并和實(shí)驗(yàn)結(jié)果進(jìn)行比較，相關(guān)實(shí)驗(yàn)數(shù)據(jù)來(lái)自文獻(xiàn)[10]。從圖中可以看出，冰力隨鋼板位移波動(dòng)上升，反映出鋼板與冰錐碰撞過(guò)程是高度非線性運(yùn)動(dòng)。剛板在低速時(shí)冰力振幅較小，而剛板在速度相對(duì)較高時(shí)，冰力振幅較大，這是因?yàn)殇摪逶诟咚俦鹊退贂r(shí)冰體破碎發(fā)生得更為頻繁。不論是低速時(shí)的計(jì)算結(jié)果還是高速時(shí)的計(jì)算結(jié)果，都能反映出計(jì)算值和試驗(yàn)值的增長(zhǎng)趨勢(shì)基本相符，該海冰材料模型適合應(yīng)用于破冰船破冰過(guò)程的數(shù)值模擬。

2 破冰船連續(xù)式破冰數(shù)值模擬

2.1 破冰船與冰排模型

圖 2 鋼板不同運(yùn)動(dòng)速度下的冰力-位移曲線Fig. 2 Curve of ice force-displacement at different steel plate speeds

在冰厚不超過(guò)1.5 m 時(shí)，破冰船通過(guò)船體首部特殊形狀將冰層局部破碎，進(jìn)行“連續(xù)式”破冰。本文以某極地破冰船在連續(xù)式破冰場(chǎng)景為對(duì)象，研究不同船速，不同冰厚對(duì)破冰阻力的影響。該船的主要參數(shù)包括：船長(zhǎng)約123 m，型寬約24 m，吃水7.8 m，船首角23°以及外飄角58°。相關(guān)材料參數(shù)如表2 所示。

表 2 極地破冰船材料參數(shù)Tab. 2 Polar icebreaker material parameters

船體首部區(qū)域采用精細(xì)網(wǎng)格，在首部以后采用簡(jiǎn)化處理方法，對(duì)碰撞影響不大及遠(yuǎn)離碰撞區(qū)域的船體結(jié)構(gòu)設(shè)為剛體。采用體單元建立120 m×60 m 的冰層有限元模型，冰層前端面處于自由狀態(tài)，左、右端面和后端面三面全約束?？紤]到冰層浮于水面，根據(jù)冰層厚度的不同，使冰層上表面的略高于水線面。對(duì)冰層施加重力的同時(shí)在冰層底部施加相應(yīng)的向上的浮力。其接觸控制使用侵蝕接觸，其關(guān)鍵字*CONTACT_ERODING_SURFACE_TO_SURFACE。為了不考慮破冰過(guò)程對(duì)船體運(yùn)動(dòng)的影響，對(duì)船體除x（定速運(yùn)動(dòng)方向）方向外的其他5 個(gè)自由度進(jìn)行剛性約束。破冰船與冰排模型如圖3 所示。

2.2 船速和冰厚對(duì)冰阻力影響

圖 3 破冰船和冰排模型Fig. 3 The model of icebreaker and ice row

為了驗(yàn)證破冰船船速和冰厚對(duì)冰阻力的影響，本文就不同冰厚（0.5 m，1 m，1.5 m），不同船速（1 m/s，1.5 m/s，2 m/s，2.5 m/s）模擬破冰船在連續(xù)式破冰場(chǎng)景中的多個(gè)工況，得到的冰力曲線如圖4 所示。

圖 4 不同船速冰厚下的冰力曲線Fig. 4 Ice force curve under different ship speeds and ice thicknesses

從圖中整體冰力時(shí)程曲線可以看出，冰力曲線整體呈鋸齒狀波動(dòng)上升或下降。在破冰船不斷向前破冰過(guò)程中，起初由于船冰接觸面積較小，因而冰力較小，隨著船冰接觸面積的逐漸增大，導(dǎo)致船冰摩擦力也增大，冰力也逐漸上升。但同時(shí)船舶破冰過(guò)程，冰排發(fā)生彎曲破壞，并伴隨冰體不斷失效，因而會(huì)產(chǎn)生不同程度的冰力卸載，冰力表現(xiàn)出鋸齒狀波動(dòng)。不同工況下冰力值如表3 所示。

表 3 不同工況下冰力值Tab. 3 Ice force values under different working conditions

通過(guò)比較不同工況下冰力時(shí)程曲線和冰力表，發(fā)現(xiàn)冰力的波動(dòng)振幅、冰力峰值和平均冰力都隨著船速和冰厚的增大而增大。這是因?yàn)檩^大的速度和較大冰厚都會(huì)對(duì)冰排產(chǎn)生較大的碰撞力，冰力因冰單元失效而卸載，從而產(chǎn)生較大的振幅，而相比于船速，冰厚對(duì)冰阻力影響更大。

3 冰阻力經(jīng)驗(yàn)?zāi)Ｐ?/h2>

3.1 Lindqvist 公式

該公式是由Lindqvist[2]在波羅的海根據(jù)模型試驗(yàn)和實(shí)尺度觀測(cè)結(jié)果，由冰的失效模式和運(yùn)動(dòng)狀態(tài)，將冰阻力分成擠壓冰阻力、彎曲破壞形成的冰阻力和冰浸沒(méi)阻力3 個(gè)部分。Lindqvist 模型的冰阻力公式如下：

3.2 Riska 公式

該公式是由Riska[3]基于航行于波羅的海不同船舶的實(shí)船觀測(cè)數(shù)據(jù)，對(duì)Lindqvist 公式進(jìn)行改進(jìn)，從而提出有關(guān)平整冰的冰阻力經(jīng)驗(yàn)公式。由于在平整冰環(huán)境下敞水阻力很小，通常忽略水動(dòng)力與冰的相互作用。Riska 模型的冰阻力公式如下：

表 4 Riska 模型經(jīng)驗(yàn)系數(shù)Tab. 4 Empirical coefficient of Riska model

3.3 Jeong 公式

該公式是由Jeong[4]在Spencer 模型的基礎(chǔ)上，根據(jù)模型實(shí)驗(yàn)提出的一種針對(duì)標(biāo)準(zhǔn)破冰船模型的冰阻力經(jīng)驗(yàn)估算方法，冰阻力由破冰阻力、浮冰阻力和清除冰力3 部分組成。Jeong 模型的冰阻力公式如下：

式中：?分別為浮冰阻力系數(shù)，清除冰力系數(shù)和破冰阻力系數(shù)；和分別為傅汝德數(shù)和強(qiáng)度因子；和分別為傅汝德數(shù)和強(qiáng)度因子的冪指數(shù)。各經(jīng)驗(yàn)系數(shù)的取值如表5所示。

表 5 各經(jīng)驗(yàn)系數(shù)的取值Tab. 5 values of each empirical coefficient

3.4 Edwards 公式

Edwards 等[5]通過(guò)模型和實(shí)船試驗(yàn)結(jié)果推導(dǎo)，并對(duì)之前的經(jīng)驗(yàn)公式重新進(jìn)行分析和改進(jìn)，得到下式：

4 計(jì)算結(jié)果分析比較

4.1 數(shù)值仿真與經(jīng)驗(yàn)公式對(duì)比分析

將經(jīng)驗(yàn)公式計(jì)算得到的破冰船在連續(xù)破冰場(chǎng)景中不同工況下的平均冰力值，與數(shù)值模擬的平均冰力結(jié)果進(jìn)行對(duì)比，得到不同工況下經(jīng)驗(yàn)公式與數(shù)值模擬的冰力對(duì)比曲線，如圖5 所示。

可以看出，在不同船速和不同冰厚條件下，由各種經(jīng)驗(yàn)公式和數(shù)值仿真得到的計(jì)算結(jié)果，都反映出了冰阻力隨船速和冰厚的增大而增大的變化趨勢(shì)。但各個(gè)經(jīng)驗(yàn)解析公式和數(shù)值模擬對(duì)該破冰船在連續(xù)式破冰場(chǎng)景中的不同工況下的冰阻力預(yù)報(bào)，都存在不同程度的差異。Riska 公式在各種工況下所計(jì)算出的冰阻力都最小。其中，在冰厚為0.5 m 時(shí)，數(shù)值仿真結(jié)果與Lindqvist 的冰阻力的計(jì)算結(jié)果偏差最大，而與Jeong 公式計(jì)算值最為接近；在冰厚為1 m 時(shí)，數(shù)值仿真結(jié)果與Lindqvist 公式和Jeong 公式的冰阻力計(jì)算結(jié)果偏差最大，而與Edwards 公式計(jì)算值最為接近；在冰厚為1.5 m 時(shí)，數(shù)值仿真結(jié)果與Jeong 公式的冰阻力的計(jì)算結(jié)果偏差最大，而與Lindqvist 公式計(jì)算值最為接近。

4.2 參數(shù)敏感性分析

敏感性分析是分別將每個(gè)自變量在一定范圍內(nèi)進(jìn)行變化，得到相應(yīng)變化趨勢(shì)和變化程度，最后得出自變量對(duì)結(jié)果的影響程度，而敏感性系數(shù)便是衡量影響程度的大小，敏感性系數(shù)越大，則該自變量對(duì)結(jié)果的影響便越大[11]。表6 列出了不同冰阻力經(jīng)驗(yàn)公式中涉及的主要參數(shù)，符號(hào)“√”表示冰阻力經(jīng)驗(yàn)公式中所包含的相應(yīng)參數(shù)。

圖 5 不同工況下經(jīng)驗(yàn)公式與數(shù)值模擬的冰力比較Fig. 5 Comparison of empirical formula and numerical simulation results under different conditions

表 6 相關(guān)主要參數(shù)Tab. 6 Related main parameters

通過(guò)對(duì)冰阻力經(jīng)驗(yàn)公式中的參數(shù)研究發(fā)現(xiàn)，計(jì)算結(jié)果由經(jīng)驗(yàn)公式中各個(gè)參數(shù)值決定，而每個(gè)參數(shù)值對(duì)結(jié)果的影響程度也各有不同。針對(duì)在船速為2 m/s，冰厚為1 m，船寬為24 m，冰彎曲強(qiáng)度為1 MPa 的初始條件下，相關(guān)參數(shù)敏感性分析如圖6 所示。

圖 6 參數(shù)敏感性分析Fig. 6 Parameter sensitivity analysis

圖6 反映了冰阻力隨各參數(shù)變化率增大而增大的變化趨勢(shì)，但每個(gè)參數(shù)值對(duì)冰阻力的影響程度也各有不同。在初始條件下，船寬的變化對(duì)Edwards 模型影響最大，若船寬增加20%，Lindqvist，Riska，Jeong，Edwards 的冰阻力計(jì)算結(jié)果值分別增加了13.1%，5.3%，6.5%和20%；船速的變化對(duì)Lindqvist 模型的冰阻力變化率影響最大，對(duì)Riska 模型影響最小。若船速增加25%，Lindqvist，Riska，Jeong 和Edwards 的冰阻力計(jì)算結(jié)果值分別增加了21.3%，9.1%，11.2%和14.2%。冰厚的變化對(duì)各經(jīng)驗(yàn)?zāi)Ｐ途休^大影響，其中對(duì)Edwards 模型的冰阻力變化率影響最大，對(duì)Lindqvist 模型影響最小，若冰厚增加20%，Lindqvist，Riska，Jeong 和Edwards 的冰阻力計(jì)算結(jié)果值分別增加了30.2%，70.1%，100.2%和108.1%。冰彎曲強(qiáng)度的變化對(duì)Jeong 模型冰阻力變化率影響最大，對(duì)Edwards 模型影響最小，若冰彎曲強(qiáng)度增加20%，Lindqvist，Jeong 和Edwards 的冰阻力計(jì)算結(jié)果值分別增加了22.1%，27.4%和16.2%。

5 結(jié) 語(yǔ)

本文通過(guò)建立破冰船與冰排相關(guān)有限元模型，選取驗(yàn)證合適的冰材料模型，對(duì)破冰船在不同冰況下的破冰過(guò)程進(jìn)行數(shù)值模擬分析，并對(duì)數(shù)值模擬和經(jīng)驗(yàn)公式計(jì)算得到的冰阻力值作了比較分析，同時(shí)對(duì)各經(jīng)驗(yàn)公式中的主要參數(shù)作了敏感性分析，得出以下結(jié)論：

1）在不同應(yīng)變率下，由鋼板與冰錐碰撞仿真計(jì)算所得的冰力-位移曲線與文獻(xiàn)中的實(shí)驗(yàn)數(shù)據(jù)結(jié)果基本吻合，說(shuō)明該海冰材料模型適用于破冰過(guò)程的數(shù)值模擬。

2）由數(shù)值仿真和各種經(jīng)驗(yàn)公式計(jì)算結(jié)果，可以發(fā)現(xiàn)，在不同船速和冰厚條件下，冰阻力均隨船速和冰厚的增大而增大，而在船速和冰厚2 個(gè)參數(shù)中，冰厚對(duì)冰阻力影響更大。

3）經(jīng)與數(shù)值仿真結(jié)果對(duì)比分析，不同的經(jīng)驗(yàn)公式適用于不同冰厚條件下的冰阻力預(yù)報(bào)：冰厚為0.5 m時(shí)，數(shù)值仿真結(jié)果與Jeong 公式計(jì)算值最為接近；冰厚為1 m 時(shí)，數(shù)值仿真結(jié)果與Edwards 公式計(jì)算值最為接近；冰厚為1.5 m 時(shí)，數(shù)值仿真結(jié)果與Lindqvist 公式的計(jì)算值最為接近，可以為進(jìn)行快速冰阻力預(yù)報(bào)時(shí)選取適合的經(jīng)驗(yàn)公式提供依據(jù)。

4）對(duì)不同經(jīng)驗(yàn)公式中的不同參數(shù)進(jìn)行敏感性分析，反映出了冰阻力隨各參數(shù)變化率增大而增大的變化趨勢(shì)，但船寬、船速、冰厚以及冰彎曲強(qiáng)度等參數(shù)對(duì)各經(jīng)驗(yàn)公式冰阻力計(jì)算結(jié)果的敏感性也各有不同。船寬和冰厚的變化對(duì)Edwards 公式影響最大，船速對(duì)Lindqvist 公式影響最大，可以為船體結(jié)構(gòu)設(shè)計(jì)提供參考。