徐 瑩 ， 胡志強 ， 陳 剛 ， 徐業(yè)峻

（1.上海交通大學海洋工程國家重點實驗室，上海200240；2.高新船舶與深海開發(fā)裝備協(xié)同創(chuàng)新中心，上海200240；3.紐卡斯爾大學 船舶科學與技術(shù)學院，英國 NE1 7RU；4.中國船舶及海洋工程設(shè)計研究院，上海 200011；5.中海油能源發(fā)展采油服務(wù)公司，天津300452）

0 引 言

極地資源的開發(fā)是國際熱點之一。有研究表明北極區(qū)油氣儲量占世界新探明油氣資源的25%[1]，科技的發(fā)展也使得夏季的北極航線成為可能。如圖1所示，北極航線相對于目前的國際貿(mào)易航線縮短許多，能夠有效提高經(jīng)濟效益，并壓縮成本。但是，船舶及海洋結(jié)構(gòu)物與冰的相互作用，嚴重情況下（例如碰撞）可能導致結(jié)構(gòu)物破損甚至沉沒，以及原油天然氣泄露污染等嚴重問題，如圖2。隨著各種海洋平臺和船只在北極地區(qū)的不斷出現(xiàn)，船舶及海洋結(jié)構(gòu)物與冰的相互作用問題成為研究熱點。北歐、俄羅斯、美國和加拿大等環(huán)北極國家對海冰的研究較早，在船冰相互作用問題的研究中處于國際領(lǐng)先地位。近年來，隨著北極開發(fā)的不斷重視，我國也涌現(xiàn)出一些杰出的研究學者和有意義的科研成果。

海冰從產(chǎn)生到發(fā)展經(jīng)歷了復雜的過程，冰的生長環(huán)境、壽命和形態(tài)，以及結(jié)構(gòu)物與其作用時的速率、結(jié)構(gòu)形式和強度等因素都會影響到冰的力學性能。目前各國學者對冰的材料特性以及船冰相互作用問題的了解有限，因此，研究冰的本構(gòu)關(guān)系，以及開發(fā)有效的數(shù)值模型是研究船冰相互作用問題的關(guān)鍵。

圖1 北極航線與目前航線 Fig.1 The arctic routes and current routes

圖2 Reduta Ordona號與冰山相撞[2]Fig.2 Reduta Ordona collided with an iceberg

最早對冰載荷的預測采用經(jīng)驗或半經(jīng)驗公式，但是經(jīng)驗公式法受限于適用范圍。試驗法是研究海冰特性以及船冰作用的最好方法之一，通過試驗可以觀察并研究真實冰的失效過程，獲得直接的數(shù)據(jù)資料，但對試驗條件要求比較嚴格，成本較高。隨著計算機技術(shù)的不斷提高，數(shù)值模擬已經(jīng)成為一種經(jīng)濟有效的研究方法。船冰模擬的數(shù)值方法主要包括有限元法和離散元法。有限元法具有成熟的理論基礎(chǔ)，能夠處理比較復雜的力學模型、邊界條件以及材料和幾何非線性問題，而且商業(yè)軟件較多，方便實現(xiàn)。但由于冰的破碎具有離散的特性，有限元法在模擬介質(zhì)的間斷和離散性方面能力有限，如何完善有限元方法使之能夠更好地模擬冰的離散特性是一個重要課題。離散元方法是一種用途廣泛的海冰模擬方法，在表現(xiàn)海冰的離散特性方面具有一定優(yōu)勢，能夠較好地模擬船冰作用中海冰從完整狀態(tài)到破碎狀態(tài)的過程。目前離散元法開發(fā)得并不完善，如何合理地確定模型相關(guān)參數(shù)（如離散單元間的剛度）以及提高計算效率都是離散元法需要解決的重要問題。由于自然水域中冰的材料，幾何和分布具有很大的隨機性，因此概率法也普遍被認為適于研究冰的某些特征。本文對上述船冰相互作用的主要研究方法進行介紹，包括試驗法、有限元法、離散元法和概率法，并對應(yīng)用上述方法進行的船冰作用研究進行概述。

1 試驗法

船冰相互作用試驗研究方法包括壓痕試驗、壓縮與拉伸試驗、彎曲試驗和大尺度碰撞試驗等。

1.1 壓痕試驗

船冰碰撞時接觸面上冰破壞的物理過程對冰載荷的測量至關(guān)重要，這一過程可通過壓痕試驗來觀察。船冰碰撞時，兩者一旦接觸會在冰面形成若干高壓區(qū)[3]，高壓區(qū)的冰處于三向應(yīng)力狀態(tài)，隨著壓力提升，高壓區(qū)的冰首先出現(xiàn)融化、重結(jié)晶等現(xiàn)象，低壓區(qū)的冰出現(xiàn)微裂紋、彈性變形，之后高壓區(qū)的冰被破壞，并將碎冰擠出，生成新的高壓區(qū)。圖3為冰的碰撞破壞示意圖。

加拿大Michel等人[4]采用矩形壓頭對垂直冰板進行壓痕試驗。Devinder等人[5]對淡水冰層進行小尺度壓痕試驗，改變壓頭寬度、冰厚、壓入速度和結(jié)構(gòu)強度。冰——結(jié)構(gòu)接觸面安裝接觸傳感器，測量壓力。結(jié)果表明，低壓痕速度下冰發(fā)生延性變形，高速下發(fā)生連續(xù)脆性壓碎。低應(yīng)變率或壓入速率下，冰的變形主要是蠕變，蠕變應(yīng)變通常比彈性應(yīng)變大；高應(yīng)變率或壓入速度下，斷裂前的冰變形大多是彈性的，因為蠕變需要時間。Li Chuanke等人[6]采用四種尺寸的壓頭在三種速率下進行小尺度壓痕試驗，試驗發(fā)現(xiàn)尺度效應(yīng)與位移率有關(guān)：在高速率下，尺度效應(yīng)比較明顯。2013年Browne[7]采用直徑20mm的壓頭沖擊多晶冰冰塊，發(fā)現(xiàn)低溫下冰易出現(xiàn)微裂紋并導致剝落（spalling）破壞，高溫下冰易出現(xiàn)重結(jié)晶并以破碎（crushing）形式破壞。還有一些采用球形或半球形壓頭進行的多晶冰壓痕試驗見參考文獻[8-10]。2017年，Scopper等人[11-12]采用不同形狀的壓頭，對不同覆蓋層的冰以不同壓入速度進行試驗，以研究這些因素對接觸力以及冰破碎情況的影響。試驗發(fā)現(xiàn)水的浸沒會抑制冰的破碎，冰的表面覆蓋雪或冰顆粒都會提高接觸力。O’Rourke[13]對過往的壓痕試驗進行了回顧，分析和研究了試驗中發(fā)生的冰載荷規(guī)律性振動（冰激振動）。壓頭附近的冰有明顯的損傷層，接觸壓力的規(guī)律性振動主要來源于該處的損傷，中間夾雜的壓力突變一般來源于大塊冰的斷裂和脫落，并伴隨著接觸面積的顯著減小。

圖3 冰的碰撞破壞[3]Fig.3 Ice collision damage

1.2 壓縮與拉伸試驗

壓縮與拉伸試驗一般可分為單軸、雙軸和三軸加載試驗。通過這些試驗可得到冰的抗壓與抗拉強度，復雜應(yīng)力和限制條件下的冰強度測試可得到冰的完整失效面（failure envelope），失效面決定了冰在各種壓力或張力狀態(tài)下的失效狀態(tài)。船冰碰撞中，壓縮破壞是冰破壞的一種主要形式，相關(guān)試驗較多。由于對試驗設(shè)備要求較高，真正的三軸試驗（σ1＞σ2＞σ3）和拉伸試驗較少[14]，許多三軸試驗為圍壓三軸試驗（σ1＞ σ2= σ3）。

1986年，Timco和Frederking[15]進行了真實的三軸壓縮試驗，用不同方式對柱狀冰、多晶冰和柱狀冰混合物進行加載，溫度為-2℃和-12℃，得到了能夠描述各向異性冰的n型多項式屈服函數(shù)。Risk和Frederking[16]對多年冰進行單軸、雙軸和三軸壓縮試驗，并將強度結(jié)果分別應(yīng)用‘Tsai-Wu’和Nadreau破壞準則在p-J2平面上擬合，發(fā)現(xiàn)‘Tsai-Wu’準則能更好地擬合試驗結(jié)果。天津大學馬鍵鈞[17]在冰室內(nèi)進行三軸壓縮試驗，根據(jù)試驗結(jié)果提出了修正的Derradji-Aout破壞準則。

應(yīng)變率對冰的破壞過程和冰載荷有較大影響，許多學者在試驗中就這一問題展開研究。1997年Jones[18]對淡水冰和波羅的海冰進行高應(yīng)變率單軸壓縮試驗，應(yīng)變率范圍10-1-10/s，溫度-11℃，發(fā)現(xiàn)海冰強度隨著應(yīng)變率提升而增大，并給出了冪次關(guān)系式。在應(yīng)變率為10/s時，淡水冰強度是海冰的1.3倍。2003年，Stephen和Gagnon等人[19]對冰山冰進行單軸抗壓試驗，應(yīng)變率范圍擴展到4*10-8/s到6/s，溫度-10℃。圖4為應(yīng)力峰值隨應(yīng)變率的變化曲線。隨著應(yīng)變率的提升，抗壓強度先升高后下降，當應(yīng)變率最高時，強度會回升，Stephen等人推測這與粒狀冰內(nèi)的裂紋擴展速度有關(guān)。

除了應(yīng)變率，圍壓大小和溫度也是關(guān)鍵因素。Barrette和Jordaan[20]進行三軸壓縮恒定應(yīng)力蠕變試驗，發(fā)現(xiàn)圍壓大于某一數(shù)值時冰活化能顯著增加，冰的強度降低。Sinha[21]對橫向各向同性冰模型進行單軸加壓/卸載蠕變試驗，發(fā)現(xiàn)冰模型在加載之后很短時間內(nèi)只表現(xiàn)出線彈性及滯后彈性行為。更多多軸加載試驗還可參考文獻[22-28]。

圖4 試驗測得應(yīng)力峰值隨應(yīng)變率的變化[19]Fig.4 Peak stress plotted against strain rate from tests

1.3 彎曲試驗

彎曲測試一般會導致非均勻應(yīng)力區(qū)，因此彎曲強度常作為指標測試[29]，試驗方法包括懸臂梁試驗和簡支梁試驗。懸臂梁試驗一般在室外進行，直接將冰梁從自然冰中切出，其中三邊自由一邊與自然冰相連。在冰梁的自由端施加一個漸變的載荷，當載荷達到某一值時連接端斷裂，從而得到相應(yīng)的彎曲強度。簡支梁試驗可在室內(nèi)進行，一般分為三點彎曲試驗和四點彎曲試驗，將載荷平行且反向地作用在冰梁上[14]。

在von Bock等人[30]對模型冰進行的懸臂梁試驗中，斷裂面應(yīng)變不均勻，上部和下部應(yīng)變較大而中性軸應(yīng)變較小，使得冰梁在彎曲時彈性響應(yīng)持續(xù)時間較長。韓紅衛(wèi)和李志軍[31]對淡水冰進行懸臂梁試驗，測量其彎曲強度和有效彈性模量。Frederking等人[32-34]也進行了冰的懸臂梁試驗。

Gagnon和Gammon[35]對冰山冰進行四點梁彎曲試驗，發(fā)現(xiàn)冰梁的彎曲強度隨應(yīng)變率提高而提高，隨溫度升高而降低。Barrette等人[36]對冰山冰進行四點梁彎曲試驗，測量梁的失效強度，建立內(nèi)部結(jié)構(gòu)特性（晶粒、缺陷）間的關(guān)系，發(fā)現(xiàn)平均晶粒尺寸越大，氣泡數(shù)量越多，強度越小。Barrette認為可將缺陷分為三種級別：（1）內(nèi)部空氣（氣泡或其他形狀）；（2）平面缺陷，表現(xiàn)為氣泡構(gòu)成的二維網(wǎng)絡(luò)；（3）形成脈絡(luò)。

大連理工大學的季順迎等人[37]在2008到2010年兩個冬季，對環(huán)渤海沿岸9個測點海冰進行三點彎曲試驗，發(fā)現(xiàn)海冰彎曲強度與鹵水體積的平方根呈負指數(shù)關(guān)系，與加載速率呈線性關(guān)系，并得到了海冰彎曲強度的雙參數(shù)推算公式。在此之后，季順迎等人[38]對2011至2012年現(xiàn)場采集的渤海萊州灣海冰試樣進行室內(nèi)低溫三點彎曲試驗，測量海冰的斷裂韌度，確定了鹵水體積和加載速率共同影響下海冰斷裂韌度關(guān)系式。

1.4 大尺度碰撞試驗

大尺度碰撞試驗以及實船測量能夠更準確地體現(xiàn)真實冰與結(jié)構(gòu)物碰撞，為此學者們進行了一些大尺度碰撞試驗以及總結(jié)工作。1995年，加拿大C-CORE公司進行了真實尺寸冰山與結(jié)構(gòu)物碰撞試驗，通過拖拽小型冰山與垂直固定的平板進行碰撞。之后C-CORE公司[39-43]進行了CCGS Terry Fox破冰船撞擊冰山試驗，冰山尺寸為4-73 m，冰山質(zhì)量為21-22 000 t，冰山表面溫度為0℃，0.6 m深處-15.4℃。Timco[44]總結(jié)了2011年以前公開發(fā)表的大塊浮冰與固定結(jié)構(gòu)物的碰撞載荷，包括大尺度試驗、實地試驗和監(jiān)測數(shù)據(jù)，以及小尺度冰碰撞試驗，綜合以上結(jié)果，Timco認為浮冰碰撞力與動能聯(lián)系緊密且存在冪次關(guān)系，并提出了一個粗略估計浮冰碰撞力的經(jīng)驗公式，即，其中F為碰撞力（kN），M為質(zhì)量（kg），v為浮冰速度（m/s）。M??tt?nen 等人[45]采用 1：3 比例對冰排與不同強度平板進行了碰撞試驗。

天津大學是目前國內(nèi)少數(shù)擁有海冰試驗池的院校之一。該試驗冰池進行了大量海洋結(jié)構(gòu)物與海冰相互作用的試驗研究，史慶增等人[46]在冰池中模擬了遼東灣抗冰平臺與冰排碰撞，黃焱[47]對冰激柔性錐體及直立樁結(jié)構(gòu)振動問題進行了冰池模型試驗研究。

2 有限元法

有限元法是常用的數(shù)值方法，具有成熟的理論基礎(chǔ)，商業(yè)軟件較多，便于實現(xiàn)。開發(fā)廣泛有效的冰模型對于數(shù)值模擬是非常重要的，但需要遵循冰的物理性質(zhì)和失效方式。

2.1 本構(gòu)模型與有限元

復雜的環(huán)境和形成過程使冰具有多種形態(tài)和內(nèi)部結(jié)構(gòu)，采用一種本構(gòu)模型模擬所有冰的特點是比較困難的，目前的本構(gòu)模型往往針對某種冰或在某種情形下有效。本部分將冰的本構(gòu)模型和有限元模擬結(jié)合起來進行討論。

1978年，Sinha等人[48]開發(fā)了一種針對單軸受壓的粘彈性本構(gòu)模型，應(yīng)變包括彈性、延遲彈性和粘性應(yīng)變。Karr和Choi等人[49]考慮了微裂紋對延遲彈性應(yīng)變和粘性應(yīng)變的影響。Ralston[50]將一種應(yīng)用于巖土力學的塑性屈服函數(shù)引入柱狀晶體冰。Jorddan等人[51]提出了兩種方法完善應(yīng)用于連續(xù)性分析的本構(gòu)模型，一個是基于流變模型的方法，另一個是基于指數(shù)函數(shù)的方法，兩種方法都能較好地模擬延遲彈性應(yīng)變。Lade[52-53]將冰看做一種摩擦材料，通過硬化/軟化塑性模型和非關(guān)聯(lián)流動法則模擬不可恢復的變形。

多晶冰的自身損傷與孔隙率是兩個重要變量。已有的損傷會導致微裂紋、重結(jié)晶與晶粒邊界的產(chǎn)生，從而提高冰的延性；孔隙也會提高冰的延性。Singh和Jorddan[54]討論了自身損傷和孔隙率的影響，用損傷函數(shù)模擬裂紋和晶粒邊界的變形，并模擬了孔隙率的變化，冰的變形包括彈性、延遲彈性和粘性蠕變變形。Keller等人[55]建立了一個針對多晶冰的粘彈性本構(gòu)模型，使用分布式纖維束模型構(gòu)建損傷函數(shù)，將損傷變量耦合入延遲彈性變形中。

前文提到，應(yīng)變率對冰的破壞過程有很大影響，本構(gòu)模型的構(gòu)造也應(yīng)考慮這一點。Jorddan等人[56]認為，在低應(yīng)變率下，蠕變與微裂紋對冰影響較大，此時冰體現(xiàn)粘彈特性，應(yīng)變率較高時容易發(fā)生脆性失效。Derradji-Aouat[57]針對淡水橫向各向同性柱狀晶S-2冰，建立三維彈塑性本構(gòu)模型。該模型由非裂紋子模型和裂紋子模型組成，適用于模擬冰在低應(yīng)變率下的延性行為，不能模擬高應(yīng)變率下的脆性失效。

試驗?zāi)Ｐ捅从沉艘环N冰的理想狀態(tài)（未變形時各向同性），其生長和冷卻過程是已知的，可通過測量得到數(shù)值模型需要的重要參數(shù)。2013年，Polach[58]針對芬蘭阿爾托大學冰池的模型冰——細粒冰，建立了模型冰的數(shù)值模型。發(fā)現(xiàn)該模型冰屈服強度相對較低，不可恢復應(yīng)變發(fā)生較早，模型中未考慮粘性和溫度的影響。

近幾年內(nèi)，采用有限元方法的船冰研究還有很多。2017年，Kolari[59]基于wing-crack模型開發(fā)了一個三維各向異性損傷模型，用以描述冰的脆性失效。將該模型嵌入ABAQUS中對冰的單軸壓縮和拉伸進行數(shù)值模擬并得到了相應(yīng)的失效模式和強度。Wong[60]針對層冰的彎曲斷裂建立了三維本構(gòu)模型，并采用有限元對圓錐形橋墩與層冰的相互作用進行了分析。Liu Zhenhui[61-62]構(gòu)建了一種塑性模型模擬冰山與結(jié)構(gòu)物的碰撞，并進行了船側(cè)和船首碰撞分析。Taylor[8]采用單元移除方法，模擬壓力軟化和破裂過程。Gagnon等人[63]采用一種硬的可壓碎泡沫建立了冰破碎的數(shù)值模型，該模型能夠模擬冰壓碎過程的兩個重要特性：快速融化和剝離。2012年，Metrikin[64]對相對速度較小的碎冰——結(jié)構(gòu)物相互作用進行模擬，結(jié)合了剛體運動與有限元方法。紐芬蘭紀念大學[65]和C-CORE[66]采用有限元法對冰山與海底結(jié)構(gòu)的相互作用進行了分析，以保證水下結(jié)構(gòu)的安全性。

哈爾濱工程大學的任慧龍等人[67]探討了破冰引起的船舯彎矩計算，翟帥帥等人[68]探討了不同本構(gòu)模型下的船冰相互作用。大連理工大學武文華等人[69]采用LS-DYNA分析了JZ20-2NW平臺的抗冰性能，韓雷等人[70]采用ANASYS分析了冰與錐形結(jié)構(gòu)物的相互作用，張大勇等人[71]對海洋結(jié)構(gòu)物的冰阻力以及冰致疲勞問題進行研究。上海交通大學石礎(chǔ)[72]建立了一個考慮溫度影響的本構(gòu)模型，并嵌入LS-DYNA模擬船冰碰撞，分析船體結(jié)構(gòu)的損傷和變形。天津大學的盧安[73]、大連海事大學的宋艷平等人[74]也采用有限元方法開展了冰與結(jié)構(gòu)物相互作用的研究。

2.2 冰的失效方式與準則

失效準則是船冰作用問題的一個重要準則，冰的失效模式包括蠕變失效，輻射裂紋，軸向斷裂，屈曲，剝落和壓碎，Timco建立了冰的失效圖，如圖5所示。

中低限制情況下（如雪），摩擦力對斷層的形成影響很大，稱為庫倫斷層；高限制情況下，摩擦滑動被抑制，易發(fā)生脆性失效，由于塑性變形作用比較明顯，稱之為塑性斷層。冰出現(xiàn)庫倫斷層或塑性斷層也與限制程度有關(guān)，Golding[75]和 Schulson[76]對這一問題進行了探討。Renshaw等人[77]開發(fā)了一種符合試驗結(jié)果的塑性斷層失效應(yīng)力定量模型，并采用該模型對淡水冰的壓碎失效圖進行完善。

近幾年來，加拿大的Derradji-Aouat等人開發(fā)的針對冰的多失效面方法得到較為廣泛的應(yīng)用。2000年，Derradji-Aouat[78-80]針對淡水各向同性冰和冰山，根據(jù)Mroz’s的多失效面理論（multi-surface failure theory），建立了針對高應(yīng)變率下的脆性失效準則，之后結(jié)合試驗數(shù)據(jù)將失效模型擴展至鹽水冰。應(yīng)用該方法對冰山冰進行模擬和驗證，能夠捕捉到高應(yīng)變率（大于10-3/s）下冰的彈性響應(yīng)和脆性失效。Martonen等人[81]根據(jù)Derradji-Aouat的失效準則模擬了冰與結(jié)構(gòu)物的相互作用問題，并與全尺度測量結(jié)果進行對比，在某些情況下吻合良好。

Schulson[82]總結(jié)了多晶冰的失效機制，尤其是脆性失效。對于冰的拉伸性能采用裂紋成核和擴展描述，對多軸加載的情況采用摩擦裂紋滑動和剪切斷裂機制來描述。對于冰在壓力作用下由延性到脆性的過渡，可以用裂紋末端蠕變、裂紋擴展進行解釋。Kolari等人[83]引入了一種網(wǎng)格更新方法來預測冰的失效，由于該方法不斷對網(wǎng)格進行更新，且不預定裂紋擴展路徑，容易導致模擬時間過長。Jordaan等人[84]認為損傷力學可以較好地描述冰的軟化等特性。Goldstein[85]開發(fā)了一種針對大接觸區(qū)的冰失效模型。

圖5 冰的失效模式與失效圖Fig.5 Ice failure mode and failure map

2.3 裂紋的模擬——粘合單元法

粘合單元法（cohesive element method）是一種描述裂紋產(chǎn)生和擴展的方法，可在有限元中實現(xiàn)，近年來得到了廣泛應(yīng)用。粘合單元法的主要假設(shè)是粘合單元符合某種力-位移關(guān)系曲線，可稱為軟化曲線，該曲線的構(gòu)造應(yīng)符合材料的具體特性?；緟?shù)有裂紋初始應(yīng)力和斷裂能。裂紋初始應(yīng)力是裂紋開始發(fā)生的應(yīng)力，斷裂能是完全打開粘結(jié)裂紋所需要的外部能量，即曲線下方的面積。這兩個參數(shù)可以通過試驗獲得，從而構(gòu)建出軟化函數(shù)。該方法模擬裂紋發(fā)展的基本過程是，有限元模型中包含普通單元和粘合單元，粘合單元按照軟化曲線關(guān)系發(fā)展，直至到達最大位移時被刪除，從而形成真正的裂紋。因此材料分為三個區(qū)域：（1）無損傷區(qū)域；（2）裂紋過渡區(qū)，有應(yīng)力在裂紋交界面?zhèn)鬟f；（3）真正的裂紋，材料完全損壞，沒有應(yīng)力在裂紋交界面?zhèn)鬟f。

Gürtner等人[86]在漢堡水池進行了相關(guān)試驗校準的冰模型試驗，這種基于能量的斷裂準則對能量釋放率的確定有較高要求。Hilding等人[87]對應(yīng)用于連續(xù)破碎的粘結(jié)單元法進行研究，并認為采取該方法可以對全尺度船冰作用下冰的連續(xù)破碎進行模擬。Lu Wenjun[88]采用粘合單元法模擬冰與錐形結(jié)構(gòu)的作用，討論了粘合單元法面臨的挑戰(zhàn)，例如網(wǎng)格依賴性，并采用兩種方法減輕其影響。王峰[89]模擬了傾斜海洋結(jié)構(gòu)物和層冰的相互作用。采用彈塑性軟化模型來代表冰的本構(gòu)關(guān)系，通過內(nèi)聚力單元的失效模擬出層冰的彎曲斷裂，并討論了網(wǎng)格和參數(shù)敏感性問題。Gribanov等人[90]建立了一個內(nèi)聚力模型來模擬柱形冰樣本的斷裂，通過有限元分析得到應(yīng)力——應(yīng)變和損傷累積的結(jié)果。

3 離散元法

有限元方法可以處理相對復雜的力學模型和邊界條件，但需要明確海冰的本構(gòu)模型和計算參數(shù)。此外，有限元方法對脆性破裂、壓碎和冰碎片的大尺度移動模擬能力不足。離散元方法是一種用途廣泛的模擬海冰特性的數(shù)值方法。離散元法可以在細觀尺度下描述海冰材料的結(jié)構(gòu)，能夠模擬冰從完整狀態(tài)到破裂狀態(tài)的過程，在宏觀上又能夠合理地模擬海冰與結(jié)構(gòu)作用時的破碎規(guī)律，獲得冰荷載的冰力時程，具有一定優(yōu)勢。

離散元模型分為兩種：硬模型和軟模型。硬模型中，粒子遵循非彈性碰撞，認為粒子間的相互作用是瞬間的，遵守動量定理，需滿足非滲透約束，通過沖力實現(xiàn)動量交換。軟接觸中，粒子可以輕微重疊，通過接觸關(guān)系計算接觸力，如彈簧-阻尼模型，分別代表動量交換和能量耗散。接觸力的計算基于粒子之間的重疊量，與粒子被賦予的材料性質(zhì)有關(guān)，如彈性、粘性、塑性參數(shù)和摩擦系數(shù)，典型的接觸力模型如圖6所示。粒子間的接觸可以持續(xù)，代價是計算時間較長。離散元法主要有兩個方面最為耗時：一是檢測粒子單元的接觸狀態(tài)，在計算接觸力之前，首先要判斷單元間的接觸狀態(tài)，這一過程需要反復的接觸探測，非常耗時；二是動力學方程求解。

多晶巖土材料，如冰和石頭，由晶體和顆粒構(gòu)成，內(nèi)部存在一些缺陷。當樣品尺度較小時，只有較少的缺陷，其失效方式主要是新裂紋的發(fā)生和擴展，而材料尺寸較大時，已有的關(guān)鍵位置的缺陷會產(chǎn)生較大影響，這是一種尺度效應(yīng)。Selvadurai等人[91]采用二維離散元方法對船冰相互作用進行模擬，引入了尺寸相關(guān)破碎應(yīng)力（size-dependent fragmentation stress）的概念，并假設(shè)應(yīng)力狀態(tài)為壓力主導的區(qū)域為粘塑性失效，應(yīng)力狀態(tài)為張力主導的區(qū)域為脆性破裂。該模型體現(xiàn)了冰的彈性、粘塑性、脆性斷裂等多種特性。

荷蘭的Dorival等人[92]采用二維離散元模型模擬冰蓋和剛性結(jié)構(gòu)的相互作用。模型劃分為六邊形網(wǎng)絡(luò)，如圖7所示，軸向彈簧——阻尼器連接位于單元中心的質(zhì)量點，失效準則為彈簧——阻尼元素的伸長量，斷裂后相應(yīng)元素從模型中移除，最終裂紋表現(xiàn)為沿某曲線的多個元素斷裂。圖8為模擬出的冰碰撞損壞后的離散元模型。為考慮冰的幾何和材料隨機性，假設(shè)彈簧剛度、阻尼系數(shù)、最大伸長量等參數(shù)服從高斯分布。實例分析表明該模型能夠捕捉到裂紋的產(chǎn)生和擴展，但仍需通過試驗得到較準確的相關(guān)參數(shù)。芬蘭的Riikil?等人[93]采用大小隨機的二維圓盤粒子，粒子之間通過無質(zhì)量彈性梁連接，模擬海冰的粘彈性變形和斷裂。梁的斷裂采用能量法，并考慮阻尼。該模型不僅具備顆粒模型的優(yōu)勢，如對接觸和斷裂的處理比較直接，同時也能夠描述冰的連續(xù)粘彈性。Tuhkuri等人[94-98]采用離散元法對層冰的彎曲、拉伸和相應(yīng)的冰載荷和冰阻力進行了很多研究，并于2018年發(fā)表了關(guān)于離散元方法在船冰相互作用中應(yīng)用進展的文章。

圖6 彈簧、阻尼器和庫倫摩擦力構(gòu)成的接觸模型[100]Fig.6 Contact model consisted of spring，dash-pot and Coulomb friction

圖7 六邊形點陣模型[92] Fig.7 Hexagonal lattice model

圖8 離散元法模擬冰撞擊結(jié)構(gòu)物后的破碎[92]Fig.8 Ice sheet damaged after crashing against the structure

2012年，Sun Shanshan等人[99]采用三維盤單元模擬荷葉冰，研究了荷葉冰對于固定圓柱樁的載荷，考慮了冰塊受到的拖曳力、附加質(zhì)量力、浮力等。由于冰載荷的振動特性，參照有義波高的定義，將碰撞力和彎矩的最大三分之一的平均值定義為有義力和有義彎矩。模擬過程中認為荷葉冰尺寸一致，并忽略相鄰浮冰之間的凍結(jié)。

Hansen等人[100]對錨泊定位于碎冰區(qū)的海洋結(jié)構(gòu)物進行響應(yīng)分析，每塊浮冰用圓盤單元表示，單元——單元和單元——結(jié)構(gòu)物的碰撞模擬為線粘——彈性，考慮庫倫摩擦接觸以防止切向相對滑動。Lau等人[101]采用基于離散元的軟件DECICE進行多項分析，并采用新的接觸檢測序列和平行算法提高軟件計算性能，縮短計算時間。Yulmetov等人[102]采用非光滑離散元方法模擬了碎冰區(qū)的冰山自由漂移和拖曳。此外還開發(fā)了一種數(shù)值冰區(qū)生成方法，能夠根據(jù)現(xiàn)場測量的分布規(guī)律模擬出隨機形狀的多邊形冰塊。Hopskins和L?set等人[103-105]也對離散元方法進行了研究。

大連理工大學采用離散元方法對海冰進行了多年研究，并應(yīng)用在船舶與海洋結(jié)構(gòu)物方面。季順迎[106]總結(jié)討論了不同尺度下海冰的離散分布特性以及塊體、圓盤、顆粒等不同形態(tài)的離散單元的應(yīng)用。2013年，季順迎等人[107]采用離散元方法模擬碎冰區(qū)海冰與船舶結(jié)構(gòu)的相互作用，碎冰由三維圓盤單元構(gòu)成，考慮了碎冰受到的浮力、拖曳力和附加質(zhì)量。同年，對海冰與直立海洋結(jié)構(gòu)作用的破碎過程進行了離散元分析[108]，將海冰單元間的粘接強度設(shè)為溫度和鹽度的函數(shù)，確定了不同鹵水體積下海冰發(fā)生脆性破壞的單軸壓縮強度，得到了冰力時程與結(jié)構(gòu)冰振響應(yīng)。此外，季順迎等人[109]還模擬了冰層與錐形結(jié)構(gòu)物之間的相互作用，采用球形離散單元，建立法向力、切向力與粒子位移和速度的關(guān)系，切向力極限為庫倫摩擦力；通過應(yīng)力修正，考慮溫度和鹽度對強度造成的影響。

4 概率法

冰的失效往往是多種失效現(xiàn)象同時發(fā)生，且由于環(huán)境和生長過程不同，導致冰在材料和幾何方面存在不均勻性，因此還可以從統(tǒng)計學角度進行分析和預測。近年來船冰問題的概率研究也得到了很大發(fā)展。

2006年，Jorddan等人[110-111]的研究表明局部壓力值具有隨機性，但具備空間與時間上的關(guān)聯(lián)性。CCORE的Li Chuanke和Jordaan等人[112]建立了概率模型，模擬局部壓力和整體壓力之間的概率平均效應(yīng)，并采用向量自回歸模型展示局部壓力之間的空間和時間關(guān)系，通過得到的壓力——面積關(guān)系和壓力——方面比關(guān)系估算總的冰載荷。Kujala等人[113]認為冰載荷僅在少數(shù)情況下服從單參數(shù)指數(shù)模型。NTNU的Suyuthi[114]建立三參數(shù)指數(shù)模型描述碎冰區(qū)作用于船體的冰載荷，該模型是兩種指數(shù)分布按照一定比例因子的結(jié)合體，并且更符合真實的數(shù)據(jù)分布規(guī)律。

C-CORE的Stuckey等人[115]搭建了一個基于概率算法的模型對冰山設(shè)計載荷進行評估。該模型采用貝葉斯技術(shù)確定冰山尺寸和碰撞速度的分布，并改良了面積——滲透深度模型和離心模型。Jordann等人[116]建立了一個概率斷裂模型以研究冰的破壞過程。首先采用ABAQUS進行應(yīng)力分析，將結(jié)果導入Matlab的PFM程序中，PFM程序能夠區(qū)分不同的應(yīng)力狀態(tài)，根據(jù)應(yīng)力狀態(tài)選擇恰當?shù)牧鸭y模型，然后對不同的裂紋模型建立失效概率計算方程，并計及裂紋的尺寸、方向和密度等因素的影響。Thijssen等人[117]采用軟件SILS（Sea Ice Loads Software），模擬結(jié)構(gòu)物與亞北極區(qū)的浮冰、冰脊間的相互作用，進行冰載荷估算和敏感性分析。該軟件采用的是ISO19906建議的概率方法，并結(jié)合蒙特卡羅算法。

5 其他方法

上述數(shù)值方法還可以互相結(jié)合以適用于更加廣泛的問題。Paavilainen等人[118-121]采用有限元——離散元結(jié)合的方法（FEM-DEM）模擬較寬傾斜結(jié)構(gòu)物周圍發(fā)生的冰塊堆積過程。邵帥等人[122-123]采用FEM-DEM方法對海洋結(jié)構(gòu)物的冰激振動進行了研究。Timoshenko梁單元兩端連接離散單元，梁間碰撞接觸力也采用離散元的方法。數(shù)值方法分析船冰問題一般工作量較大且比較耗時，Ehlers等人[124]曾采用數(shù)值和解析相結(jié)合的方法降低計算成本。此外實地監(jiān)測也是研究的重要手段，可以獲得最直接的數(shù)據(jù)。

薛彥卓等人[125-127]采用近場動力學（peridynamics）對冰的三點彎曲試驗、冰與加筋板的接觸以及船舶在冰區(qū)的航行進行模擬。葉禮裕等人[128]采用近場動力學對冰——槳相互作用進行了模擬。Robb等人[129]結(jié)合平滑粒子流體動力學模型（SPH）與離散元方法模擬了含有漂浮固體物質(zhì)的自由表面液體流動。

6 結(jié) 語

試驗法是冰研究的重要方法，也是數(shù)值模擬的基礎(chǔ)。有限元方法具有深厚的理論基礎(chǔ)，應(yīng)用廣泛，相應(yīng)的商業(yè)軟件較多，方便實現(xiàn)；離散元方法的優(yōu)點在于更適合模擬冰的破碎和離散特性；概率法則從統(tǒng)計學角度解釋冰的隨機性。船冰相互作用主要有兩個問題：一是本構(gòu)模型的開發(fā)與完善，二是數(shù)值方法的改進。對本構(gòu)模型的開發(fā)一般需要試驗數(shù)據(jù)的支持，也可參考相似材料和領(lǐng)域中的方法，例如混凝土和某些巖石。對數(shù)值方法的改進，例如離散元法，已有很多應(yīng)用但目前開發(fā)得并不完善，此外，與冰相關(guān)的某些子問題可以采用統(tǒng)計學的方法，并與有限元或離散元方法相結(jié)合。