極地某構型船首與浮冰碰撞的數值模擬研究

劉 昕 傅昱曉 蔣武杰 倪寶玉＊

（1. 中國人民解放軍92942 部隊 北京 100161；2. 哈爾濱工程大學 船舶工程學院 哈爾濱 150001；3. 上海交通大學 船舶海洋與建筑工程學院 上海 200030）

0 引 言

在地球上，冰川的覆蓋面積約為1 600 萬km2,約占地球表面積的3.1%[1]。隨著全球氣候變暖以及冰川逐漸消融，北極的戰略意義和商業價值正被逐步發掘。與其他海域相比，極地地區由冰川融化而生成的浮冰由于面積更大、冰層更厚，因此與船舶發生碰撞時的危險性也大大提高。這些大尺寸浮冰融化、相互碰撞，又可能會生成更多數量的小塊浮冰。這些浮冰主要集中在冰蓋區與開闊水域中間的邊緣冰區內，這就使船舶在極地航行時，與浮冰碰撞的概率大大增加。因此，模擬冰材料、構建船體模型以及分析冰船碰撞的動態過程，逐漸成為極地領域熱門的研究方向。

基于冰材料物理與力學特性，學者們提出了各種可用于數值模擬的冰材料模型。SCHNEIDER 等[2]依據海冰的失效模式建立了1 種幾何網格方法，用于計算海冰與船體的相互作用。SAWAMURA 等[3]基于瞬態載荷與浮冰作用的過程，得到了浮冰的彎曲特性。EHLERS 等[4]通過把數值模擬結果與早期進行的一系列海冰四點彎曲實驗進行對比，發現數值模擬結果在受力、破壞時間和位移方面與試驗結果吻合較好，驗證了塑性材料可用于模擬冰船碰撞時的海冰材料。馬國寶[5]通過討論海冰性能的影響因素，確定海冰材料在進行數值模擬時應采用彈塑性材料模型，同時模擬圓臺冰撞擊剛性墻的過程，將得到的壓力-面積曲線與ISO 給定曲線進行對比，驗證了模擬冰材料的可行性。

針對船舶與浮冰的作用，LUBBAD 等[6]通過數值模擬方法與實驗方法，研究了冰船碰撞問題，分析船舷承載能力、結構變形、結構穩定性以及板架的承載能力。何菲菲[7]采用理論推導和數值仿真方法，基于非線性有限元法，對破冰船船體運動響應以及破冰能力進行了探索。韓月[8]通過模擬破冰船與海冰的相互作用，討論了不同形狀船首的破冰性能，并基于數值模擬方法預報了船體的運動響應和破冰載荷。張健等[9]提出了連續破冰模式下的船首部冰載荷分析方法，并與數值仿真結果進行對比。

本文首先通過整理海冰的物理力學特性，確定模擬海冰材料的參數，運用LS-DYNA 軟件中的各向同性材料來模擬冰材料；并設置了冰錐與鋼板撞擊過程，通過與ISO 規范提出的P-A 曲線[10]對比，驗證模擬冰材料的可行性；然后基于某極地船型原型構建船首模型，通過模擬船首以不同速度與浮冰碰撞的工況，分析船首速度對于冰船碰撞的影響。

1 模擬海冰材料及其有效性驗證

在研究船舶與海冰碰撞的過程中，LS-DYNA是較常用的軟件。本章將根據海冰材料的物理力學特性，利用LS-DYNA 設置海冰材料的相關參數，并進行冰材料模型的驗證。

1.1 模擬海冰材料

通過對海冰材料的整理，得到其物理力學特性如表1[1]所示。

表1 海冰的物理力學特性統計值

海冰是具有各項異性的彈塑性材料，但考慮到模擬各向異性材料耗費的計算資源較大，且對于研究船舶與浮冰碰撞過程來說，各向異性與各向同性材料的差異性可忽略。因此，在研究浮冰與海上大型船舶相互作用時，可采用LS-DYNA 中的13 號材料（*MAT_ISOTROPIC_ELASTIC_FAILURE）來模擬海冰材料，并輸入參數值：密度為910 kg/m3、剪切模量為2.2 GPa、屈服應力為2.12 MPa、塑性硬化模量為4.26 GPa、體積模量為5.26 GPa、塑性失效應變值為0.35。

該材料是帶有塑性應變失效準則的各項同性彈性斷裂失效模型，即當材料單元的應力或應變達到最大許用值時，單元會失效且無法再承載應力。

1.2 海冰材料的有效性驗證

為驗證模擬冰材料的合理性，本文通過對剛性板撞擊錐形冰體進行模擬，將模擬結果與ISO 收錄的P-A 理論標準曲線[9]進行對比驗證。其中，錐形冰體和剛性板的模型如圖1 所示。

圖1 錐形冰體與剛性板模型

在圖1 中，上部圓錐半徑為0.5 m、高為0.289 m，下部圓柱半徑為0.5 m、高為0.1 m。剛性板為1.2 m×1.2 m 的正方形板，厚度為0.1 m。錐形冰的網格尺寸為0.04 ~ 0.05 m，剛性板網格尺寸為0.03 m。錐形冰前端距剛性板表面距離為0.02 m，剛性板以0.1 m/s 的速度向錐形冰撞去。剛性板材料選用了LS-DYNA 材料庫中第20 號材料“020-RIGID”。該材料不會變形且沒有應力、應變等，相當于固定的剛體，因此在計算中，僅需采集冰錐的受力情況。為保證碰撞的連續發生，錐形冰底部所有節點的6個自由度運動均被約束。

在剛性板向錐形冰撞擊的過程中，定義剛性板與錐形冰的接觸為侵蝕接觸，即冰單元體失效時,程序會自動搜索下一個接觸節點。其中，從面定義為浮冰模型，主面定義為剛性板模型。計算過程中,僅包含從面上的作用力。將靜摩擦系數設置為1.5,動摩擦系數設置為1.0。在對稱平面的選擇上，不包括有正常邊界約束的面；在內部節點的選擇上,儲存器的分配使腐蝕接觸可以發生，否則，單元格被侵蝕后將沒有接觸。在對固體元素進行相鄰材料的處理上，設定實體元素面位于材料子集的邊界上。圖2 為剛性板撞擊錐形冰體的應力云圖。

圖2 剛性板撞擊錐形冰體的應力云圖

由圖2 可見，在剛性板與錐形冰碰撞的過程中，隨著剛性板的運動，錐形冰逐漸斷裂，冰單元失效且失效面不平整。從時間分別為0.125 s、0.15 s、0.17 s 和0.19 s 的錐形冰體應力云圖可以看出，接觸面的失效應力分布不均勻，存在明顯的高應力區和低應力區。其中，最大應力分別為12.950 MPa、38.785 MPa、36.400 MPa 和39.701 MPa。

為驗證冰材料的合理性，需要將剛性板-錐形冰體碰撞的壓強-面積曲線與ISO 標準理論曲線進行對比。國際標準化委員會在ISO/CD 19906（2010）中，推薦了冰山碰撞載荷下造成嚴重事故時的冰體材料P-A 曲線[9]，這是1 條反比例曲線。

本文所模擬的剛性板與錐形冰碰撞中，接觸面的壓強與面積關系式為：

式中：P為接觸面的壓強，MPa；A為碰撞接觸面積，m2。

剛性板與錐形冰碰撞的壓強-面積曲線和ISO標準理論曲線的對比情況如圖3 所示。

圖3 壓強-面積曲線對比情況

由圖3 可見，初始階段壓強值都比較大，隨著碰撞接觸面積的增加，壓強值迅速減小并趨于平穩。當面積達到0.4 m2時，本文模擬冰材料對應的曲線產生小幅波動。但整體來看，兩者吻合度較高，所以本文可以選取該材料模型進行極地海冰與船舶相互作用的數值模擬。

2 船冰碰撞動態響應的工況分析

2.1 典型工況分析與主要船首部件受力分析

本文所用的船體模型是基于某極地船型建立的，且該船型非破冰船型，在極區航道航行時，會盡量避免與大塊浮冰、冰脊碰撞，浮冰厚度主要參考較薄的一年冰或多年冰來選取。本文設置的冰厚為0.6 m，浮冰形狀為正方形，邊長為10 m。遇到較大尺寸的浮冰時，建議繞行。綜合考慮計算精度與計算性能，冰單元的大小取200 mm×200 mm×100 mm。為模擬浮冰在海上漂浮的狀態，設置浮冰模型初始狀態為靜止，邊界條件為自由邊界，可在6 個自由度上進行自由運動。

考慮到整船結構建模計算量巨大，且船冰碰撞主要發生在船首區域，因此只選擇船首進行建模，如圖4 所示，網格大小為300 mm×400 mm。船首模型除外部殼體，還包括甲板、縱骨、縱桁、橫艙壁、肋骨、肋板和舭部等24 個殼部件和47 個梁部件，船首與浮冰距離為0.5 m。

圖4 船首有限元模型

模型船首材料采用高強度鋼，當應力超出材料的彈性極限時，它便由彈性狀態變為塑性狀態。本文用LS-DYNA 材料庫中的“001_ELASTIC”來模擬鋼材料。該材料是一種各向同性的彈性材料，當所受應變超出其承受范圍時會導致結構變形失效，對應參數選擇為楊氏模量2 060 GPa、密度為7 850 kg/m3、泊松比為0.3。

選擇船首速度6 kn 為典型工況，整個碰撞過程持續約5 s，各時間點船首與浮冰碰撞的應力云圖及浮冰模型的應力云圖如下頁圖5 所示。

可以看出，在船首與浮冰發生第1 次碰撞時，船首的應力開始由首端向后傳遞，但船首端部與浮冰接觸的區域仍然是應力最大區域。浮冰也因與船首碰撞而發生破壞，但此時碰撞速度較慢（即加載速率較小），浮冰模型表面雖有明顯裂紋，但裂紋并未擴張，且碰撞過程中，船首的一部分動能轉化為浮冰的動能和內能。所以在第1 次碰撞后，浮冰有了速度，并發生位移和翻轉。

在t=3 s 附近，船首與浮冰發生第2 次碰撞。浮冰由于第1 次碰撞時已具備一定速度，所以在這次碰撞中，船首和浮冰的碰撞應力相對較小，浮冰模型只發生了局部破壞，并沒有明顯的裂紋擴張，船首所受應力也有所減少。

對于殼部件而言，船首底部甲板是碰撞區域主要的支撐部件，船首端部結構是與浮冰直接接觸的部件。在船首區域，船底板距離碰撞區域很近，而且也是承受總縱彎曲的主要部件之一。所以，殼部件主要研究這3 個部件。目標船型為縱骨架式船，總縱彎曲主要由縱向梁分擔，所以選取的梁部件均為縱向梁，包括上層建筑主要縱桁、主甲板縱骨和船底甲板縱骨，如下頁圖6 所示。

圖6 船首主要部件選取

圖7 給出了船首上述各主要受力部件所包含單元，在每一時間步所受應力的最大值組成的曲線，表2 則給出了此過程中上述各部件的最大應力及對應時間。

圖7 船首主要部件上單元的最大應力

表2 各部件最大應力及對應時間

可以看出：縱向梁部件中承受應力較大的是上層建筑主要縱桁和船底甲板縱骨，其應力最大值分別達到了26.95 MPa 和32.82 MPa；對于殼部件而言，應力在達到最大值之后波動減小，其中承受應力最大的是直接與浮冰接觸的船首端部構件，在0.4 s 達到最大值105.75 MPa。

由于各部件相對碰撞區域的位置不同，因此應力的峰值不同，并且距碰撞區域越遠的構件，應力峰值出現的時間也相對較晚。各部件結構不同，縱向梁部件碰撞峰值后的應力波動情況相對較平緩，而殼部件則呈現間接脈沖式分布。

2.2 不同速度下船首受力情況對比

破冰船在浮冰水域航行時，既要完成破冰，還會受到碎冰和小型浮冰的阻力，因此航行速度較慢。為研究不同速度對船首與浮冰碰撞模擬結果的影響，本文選取船首航行方向（x軸方向）速度分別為5 kn（工況一）、6 kn（工況二）和7 kn（工況三），其余方向的運動都被約束。浮冰的速度設置為0 m/s，6 個自由度的運動均為自由運動，尺寸保持10 m×10 m×0.6 m 不變。

圖8 所示為工況一、工況二和工況三對應的船首所受合力-時間圖像。可以看出，上述各工況下均出現2 個較明顯的波峰，說明船首與浮冰發生了2 次較劇烈的碰撞（參見圖5 工況）。

圖8 所選工況下的船首所受合力-時間圖

隨著速度增大，波峰的出現逐漸提前，第2 次碰撞的波峰情況尤其明顯。通過比較2 個峰值，可以看出船首速度對于第1 次碰撞的冰力影響更為明顯。

由于冰船相距較近，浮冰初始狀態是靜止的，因此第1 次波峰出現的時間比較接近，且隨著船速增大，各工況下對應的船首所受合力的最大值隨之增大，如表3 所示。對于冰船第2 次碰撞，隨著船速增大，第2 次的峰值也隨之增大，且出現時間隨之提前。

表3 所選工況下的合力極大值

圖9 給出了各工況下的船首與浮冰碰撞的結果。

圖9 工況一、工況二和工況三的碰撞結果

由圖9 可見，速度對于浮冰失效方式的影響，在工況一（5 kn）和工況二（6 kn）時，海冰的加載速率較小，浮冰模型表面的裂紋沒有明顯擴張；在工況三（7 kn）時，隨著船速增大，海冰的加載速率變大，浮冰模型的裂紋明顯擴張，最終變為兩半，出現了劈裂破壞。

3 結 語

本文根據冰材料常見的物理及力學特性，基于LS-DYNA 中的各向同性材料設置了冰材料模型，模擬錐形冰撞擊鋼板的過程，通過與ISO 規范提出的壓力-面積曲線對比，驗證了模擬冰材料的可行性。

文中設置3 組速度分別為5 kn、6 kn 和7 kn時的工況，分析船首與浮冰碰撞過程。船首速度決定了外力對浮冰的加載速度。船速越大、冰載荷也越大，速度為7 kn 時的合力達到2.1 MN，且碰撞后浮冰會發生翻轉。隨著船首前進，冰與船會發生二次碰撞。由于本船型船首瘦削，所以浮冰容易產生劈裂破壞。在碰撞過程中，船首端部構件承受的應力最大，建議根據實際航行工況決定是否對該部分構件進行加固。