基于數值仿真與經驗公式的船舶破冰阻力預報

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

(1. 江蘇科技大學 船舶與海洋工程學院，江蘇 鎮江 212003；2. 江蘇新揚子造船有限公司，江蘇 泰州 214500)

0 引 言

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

由于船舶破冰是一種復雜的非線性動態響應過程，再加上海冰本構模型的復雜性，國內外還沒有一種有效而準確的冰阻力預報方法。冰阻力主要受冰材料特性、船體形狀、船舶運動和航速等多種因素的影響[1]。目前，計算冰阻力的方法主要有3 種：基于模型試驗或全尺度試驗數據生成阻力經驗關系的經驗公式[2-5]；基于船舶運動方程描述船舶與冰動態作用過程的數值模擬方法[6-7]；模型或全尺度試驗研究方法[8]。

1 冰材料及模型驗證

準確的海冰材料模型是保證數值模擬結果準確度的關鍵，然而由于冰材料的復雜性，國內外對冰材料本構模型的研究尚未形成共識。本文選取的冰材料模型為*MAT_124 等向彈塑性材料[9]。該材料模型中需定義3 種不同的材料曲線：定義2 種不同的屈服應力與有效塑性應變曲線；一種用于拉伸，另一種用于壓縮；第3 條曲線是應變速率與比例壓縮屈服應力因子（CYSF）曲線，用于說明應變速率對壓縮強度的影響，將等效塑性應變值作為失效標準。

表 1 冰材料參數Tab. 1 The material parameters of ice

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

通過模擬鋼板擠壓圓錐冰體，將碰撞力結果和現有試驗結果進行比較，以驗證所選冰材料的適用性。為與試驗模型中冰錐形狀尺寸基本保持一致，并模擬冰錐固定在試驗裝置上，建立直徑為10 cm，錐角為30°的冰錐有限元模型，將冰錐底平面節點進行全約束。由于在鋼板與冰錐的擠壓試驗中，鋼板的厚度較大，表面無明顯的變形，故將鋼板設為剛體，這樣不僅節約計算時間又不會產生較大誤差。有限元模型如圖1 所示。為驗證該冰材料在不同應變率下的適用性，使鋼板分別以1 mm/s 和100 mm/s 的速度擠壓冰錐，圖2 為鋼板以不同速度與冰錐相互作用的冰力-位移曲線，并和實驗結果進行比較，相關實驗數據來自文獻[10]。從圖中可以看出，冰力隨鋼板位移波動上升，反映出鋼板與冰錐碰撞過程是高度非線性運動。剛板在低速時冰力振幅較小，而剛板在速度相對較高時，冰力振幅較大，這是因為鋼板在高速比低速時冰體破碎發生得更為頻繁。不論是低速時的計算結果還是高速時的計算結果，都能反映出計算值和試驗值的增長趨勢基本相符，該海冰材料模型適合應用于破冰船破冰過程的數值模擬。

2 破冰船連續式破冰數值模擬

2.1 破冰船與冰排模型

圖 2 鋼板不同運動速度下的冰力-位移曲線Fig. 2 Curve of ice force-displacement at different steel plate speeds

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

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

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

2.2 船速和冰厚對冰阻力影響

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

為了驗證破冰船船速和冰厚對冰阻力的影響，本文就不同冰厚（0.5 m，1 m，1.5 m），不同船速（1 m/s，1.5 m/s，2 m/s，2.5 m/s）模擬破冰船在連續式破冰場景中的多個工況，得到的冰力曲線如圖4 所示。

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

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

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

通過比較不同工況下冰力時程曲線和冰力表，發現冰力的波動振幅、冰力峰值和平均冰力都隨著船速和冰厚的增大而增大。這是因為較大的速度和較大冰厚都會對冰排產生較大的碰撞力，冰力因冰單元失效而卸載，從而產生較大的振幅，而相比于船速，冰厚對冰阻力影響更大。

3 冰阻力經驗模型

3.1 Lindqvist 公式

該公式是由Lindqvist[2]在波羅的海根據模型試驗和實尺度觀測結果，由冰的失效模式和運動狀態，將冰阻力分成擠壓冰阻力、彎曲破壞形成的冰阻力和冰浸沒阻力3 個部分。Lindqvist 模型的冰阻力公式如下：

3.2 Riska 公式

該公式是由Riska[3]基于航行于波羅的海不同船舶的實船觀測數據，對Lindqvist 公式進行改進，從而提出有關平整冰的冰阻力經驗公式。由于在平整冰環境下敞水阻力很小，通常忽略水動力與冰的相互作用。Riska 模型的冰阻力公式如下：

表 4 Riska 模型經驗系數Tab. 4 Empirical coefficient of Riska model

3.3 Jeong 公式

該公式是由Jeong[4]在Spencer 模型的基礎上，根據模型實驗提出的一種針對標準破冰船模型的冰阻力經驗估算方法，冰阻力由破冰阻力、浮冰阻力和清除冰力3 部分組成。Jeong 模型的冰阻力公式如下：

式中：?分別為浮冰阻力系數，清除冰力系數和破冰阻力系數；和分別為傅汝德數和強度因子；和分別為傅汝德數和強度因子的冪指數。各經驗系數的取值如表5所示。

表 5 各經驗系數的取值Tab. 5 values of each empirical coefficient

3.4 Edwards 公式

Edwards 等[5]通過模型和實船試驗結果推導，并對之前的經驗公式重新進行分析和改進，得到下式：

4 計算結果分析比較

4.1 數值仿真與經驗公式對比分析

將經驗公式計算得到的破冰船在連續破冰場景中不同工況下的平均冰力值，與數值模擬的平均冰力結果進行對比，得到不同工況下經驗公式與數值模擬的冰力對比曲線，如圖5 所示。

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

4.2 參數敏感性分析

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

圖 5 不同工況下經驗公式與數值模擬的冰力比較Fig. 5 Comparison of empirical formula and numerical simulation results under different conditions

表 6 相關主要參數Tab. 6 Related main parameters

通過對冰阻力經驗公式中的參數研究發現，計算結果由經驗公式中各個參數值決定，而每個參數值對結果的影響程度也各有不同。針對在船速為2 m/s，冰厚為1 m，船寬為24 m，冰彎曲強度為1 MPa 的初始條件下，相關參數敏感性分析如圖6 所示。

圖 6 參數敏感性分析Fig. 6 Parameter sensitivity analysis

圖6 反映了冰阻力隨各參數變化率增大而增大的變化趨勢，但每個參數值對冰阻力的影響程度也各有不同。在初始條件下，船寬的變化對Edwards 模型影響最大，若船寬增加20%，Lindqvist，Riska，Jeong，Edwards 的冰阻力計算結果值分別增加了13.1%，5.3%，6.5%和20%；船速的變化對Lindqvist 模型的冰阻力變化率影響最大，對Riska 模型影響最小。若船速增加25%，Lindqvist，Riska，Jeong 和Edwards 的冰阻力計算結果值分別增加了21.3%，9.1%，11.2%和14.2%。冰厚的變化對各經驗模型均有較大影響，其中對Edwards 模型的冰阻力變化率影響最大，對Lindqvist 模型影響最小，若冰厚增加20%，Lindqvist，Riska，Jeong 和Edwards 的冰阻力計算結果值分別增加了30.2%，70.1%，100.2%和108.1%。冰彎曲強度的變化對Jeong 模型冰阻力變化率影響最大，對Edwards 模型影響最小，若冰彎曲強度增加20%，Lindqvist，Jeong 和Edwards 的冰阻力計算結果值分別增加了22.1%，27.4%和16.2%。

5 結 語

本文通過建立破冰船與冰排相關有限元模型，選取驗證合適的冰材料模型，對破冰船在不同冰況下的破冰過程進行數值模擬分析，并對數值模擬和經驗公式計算得到的冰阻力值作了比較分析，同時對各經驗公式中的主要參數作了敏感性分析，得出以下結論：

1）在不同應變率下，由鋼板與冰錐碰撞仿真計算所得的冰力-位移曲線與文獻中的實驗數據結果基本吻合，說明該海冰材料模型適用于破冰過程的數值模擬。

2）由數值仿真和各種經驗公式計算結果，可以發現，在不同船速和冰厚條件下，冰阻力均隨船速和冰厚的增大而增大，而在船速和冰厚2 個參數中，冰厚對冰阻力影響更大。

3）經與數值仿真結果對比分析，不同的經驗公式適用于不同冰厚條件下的冰阻力預報：冰厚為0.5 m時，數值仿真結果與Jeong 公式計算值最為接近；冰厚為1 m 時，數值仿真結果與Edwards 公式計算值最為接近；冰厚為1.5 m 時，數值仿真結果與Lindqvist 公式的計算值最為接近，可以為進行快速冰阻力預報時選取適合的經驗公式提供依據。

4）對不同經驗公式中的不同參數進行敏感性分析，反映出了冰阻力隨各參數變化率增大而增大的變化趨勢，但船寬、船速、冰厚以及冰彎曲強度等參數對各經驗公式冰阻力計算結果的敏感性也各有不同。船寬和冰厚的變化對Edwards 公式影響最大，船速對Lindqvist 公式影響最大，可以為船體結構設計提供參考。