基于高性能離散元方法的極區浮式平臺冰載荷數值分析

孔帥，季順迎，季少鵬，王迎暉，剛旭皓

1 中國船舶科學研究中心, 江蘇 無錫 214082 2 大連理工大學 工業裝備結構分析國家重點實驗室, 遼寧 大連 116024

0 引 言

近 40 年來，受全球氣候急劇變暖的影響，極地海冰覆蓋面積持續減少，引起了各國開發利用北極自然資源的熱情[1-2]。目前，我國在固定式海油鉆井平臺的運營方面已具有較為完備的設計建造基礎[3]。移動式結構中的浮式結構適用于水深較深的海域，更加便于在極區的移動作業[4-6]，但我國在冰區浮式結構運營方面尚處起步階段。因此，建立一套科學、有效的冰載荷預報技術是冰區浮式結構優化設計和冰激疲勞損傷特性分析的關鍵。

近年來，基于單元法的冰載荷預報分析技術已受到各國學者的關注，例如有限單元法[7]、近場動力學[8]和離散元法（discrete element method，DEM）[5,9-11]。離散元法由Cundall 等[12]提出，適于模擬散體材料在準靜態或動態載荷激勵下的運動、變形及破壞過程。寒區海域中的海冰呈現出很強的離散分布特性，離散元模型已被成功應用于碎冰動力過程[5,9-10]及海冰流變學的研究[11]。

為進一步精細化模擬浮式結構?冰相互作用的過程，離散元法需要在海冰的破碎分析、海冰拓撲形態表征及計算效率等方面有所突破。在海冰與浮式結構的相互作用過程中，海冰的破碎屬于典型的脆性破壞過程[13]。Potyondy 等[14-22]基于球形單元提出了平行黏結模型，該模型通過在接觸單元之間建立梁單元模型來傳遞單元之間的力和力矩，并根據梁單元上的最大應力是否達到預設強度來判斷材料是否破碎。平行黏結模型的工程應用性較強，適用于海冰等材料的脆性破壞分析[17-22]。Voronoi 圖是幾何學中重要的圖形概念，該圖可以有效利用限定空間內特征點的分布形式來表征自然狀態下物質的離散特性。該方法已被朱紅日等[22]應用于極區浮碎冰的模擬生成中，具有生成效率高、復現能力強的優點。在離散元數值模型中，接觸判斷及單元搜索占據了較多的計算資源，若采用串行處理技術的傳統離散元數值分析方法，浮式結構冰載荷的分析效率將會極低。龍雪等[17,20-22]運用并行處理化技術，建立了高性能離散元算法并用于分析導管架平臺及船體結構的冰載荷。

為了分析浮式結構冰區作業時的冰載荷，本文將首先基于高性能離散元分析模型及采用平行黏結模型對海冰的破碎過程進行模擬，然后，再模擬構造冰情，分析平臺與海冰的作用過程，最后，利用實測數據對離散元數值結果進行驗證。

1 模擬海冰動力過程的離散元法

基于單元模擬技術的離散元法可以較為便捷地從細觀角度描述海冰單元間的動力行為，且可有效構造復雜的冰型，例如碎冰[22]和冰脊[23]等。對于海冰的破碎過程，可通過引入平行黏結模型進行模擬[14,17-18]。圖形處理器(graphics processing unit，GPU)擁有數量可觀的核且每個核又包含多個線程，近年來，英偉達公司的CUDA-C 產品利用GPU集群處理技術為計算模型的并行化處理提供了新的研究思路。離散元算法中，接觸判斷和單元搜索計算資源占據了絕大部分的計算資源，為此，可通過結合單元鄰居鏈表算法及背景網格搜索法進行并行化加速計算處理，計算效率是串行程序的數十倍[24]，在復雜作業環境模擬、多介質場耦合分析及超大計算域分析等方面具有較強的計算優勢。

1.1 海冰單元接觸模型

在離散元法中，通常把散體看作具有一定形狀和質量顆粒單元的集合，單元之間的運動方程相對獨立，但可以發生接觸并產生相互作用力。在顆粒的相互作用過程中，考慮因單元之間相對速度和彈性變形而引起的作用力。顆粒單元之間的碰撞過程可采用彈簧?阻尼器?滑塊的唯象模型進行模擬[19]，并且單元之間的接觸力可解耦為法向分量和切向分量，解耦后的法向與切向接觸模型如圖1 所示。因此，第i個接觸對中單元之間因接觸而產生的力Fi可 由法向接觸力Fni和切向接觸力Fsi疊加計算，即

圖1 單元間法向與切向接觸模型Fig. 1 Normal and tangential contact model between elements

1.2 模擬海冰破碎的平行黏結模型

海冰的破碎過程直接影響船體結構冰載荷的頻率和數值，海冰在船體的持續擠壓過程中將變為尺寸更小的冰塊。通過分析觀測得到的北極海冰資料，可知其破碎冰多呈多邊形的幾何形狀。圖2所示為采用Voronoi 算法中的種子點分布特性控制海冰幾何規則度而生成的碎冰域，該算法可有效控制碎冰場中的海冰面積及其分布特征。考慮破碎冰的再破碎，本文采用平行黏結模型模擬冰層的凍結作用。海冰內部不僅存在力的傳遞，還存在力矩的傳遞，因此，對海冰進行離散元分析時常采用平行黏結單元[17-22]。平行黏結是將2 個球體黏結在一起，其黏結單元不僅可以傳遞力，還可以傳遞力矩，如圖2 左下角所示。

圖2 模擬海冰凍結效應的平行黏結模型Fig. 2 Parallel bond model for simulating the freezing effect of sea ice

2 錨泊系統作用下浮式結構的冰載荷分析

冰區采油平臺或者采油船等浮式結構通常采用錨泊系統在冰區作業，同時，周邊會有破冰船協助破冰。Kulluk 石油鉆井浮式平臺隸屬于殼牌公司，其曾長期服役于阿拉斯加等北極地區，在冰區積累了大量的冰情及冰載荷實測資料[5]。

2.1 浮式結構運動求解

Kulluk 浮式平臺的整體結構為雙圓錐體，向下的傾斜角為31.4°，該傾斜角的設計使得海冰在此易發生彎曲破碎而降低冰載荷；其上甲板處直徑為81 m，吃水為11.5 m，水線處直徑70 m，排水量28 000 t。平臺現場工作圖及對應的有限元模型如圖3 所示，平臺結構外表面由三角形單元構成。冰區作業中，浮式結構會受到海水浮力、海流拖曳力及海浪的作用，本文在采用離散元法模擬浮式結構與海冰相互作用時，將只考慮浮力與拖曳力對浮式結構的影響，而忽略海浪的作用。浮式結構在整體坐標軸方向平動和沿局部坐標軸轉動的動力方程可以表示為：

圖3 Kulluk 浮式平臺Fig. 3 Kulluk floating platform

浮式結構的錨泊系統可以采用張緊式錨鏈將錨泊點與浮式平臺連接起來[25]，具體形式如圖4(a)所示。錨鏈系統采用12 根錨鏈將平臺進行約束，錨鏈的具體分布形式如圖4(b) 所示。圖中：D，d分別為Kulluk 平臺雙圓錐結構上、下表面的外直徑；L為錨鏈長度；h為水深；xG，yG分別為橫、縱坐標軸；θ 為錨鏈間夾角。

圖4 錨鏈布置圖Fig. 4 Layout of the mooring lines

張緊式錨泊系統采用線性彈簧模擬單一錨鏈的恢復力T：

2.2 對比分析實測數據

根據在Kulluk 浮式平臺現場實測的數據，可知平臺周圍主要以破碎冰為主。利用Voronoi算法生成浮冰，浮冰區域面積為 300 m×200 m，破碎冰平均面積為 200 m2，破碎冰的運動由速度控制在0.5 m/s 的固定邊界推進。離散單元其他模擬參數如表1 所示。

表1 Kulluk 浮式平臺與碎冰作用時的離散單元計算參數Table 1 Computational parameters of DEM when Kulluk floating platform interacts with pack ice

Kulluk 浮式平臺在冰厚為0.8 m 的碎冰區作業時，縱向、橫向和垂向的冰載荷F時程曲線如圖5 所示。從該冰力時程曲線上可以看出，平臺的縱向、垂向冰力更連續且數值較大。其中，縱向冰力對平臺漂移產生作用，垂向冰載荷是影響結構升沉運動的主要因素，而在橫向，由于計算冰域與結構的對稱性，導致其冰載荷均值接近0，大體成對稱性分布。

圖5 3 個方向的船體冰載荷時程曲線Fig. 5 Time histories of ice loads in three directions

在Kulluk 浮式平臺報告中，冰載荷數據為作用于結構上的總體冰載荷，即圖5 中3 個方向冰載荷的合力。另外，報告中未有效分清海冰冰況的密集度、冰型和冰塊尺寸等信息，僅將其復雜的工況簡化成了冰厚?冰載荷限值形式（式（22）和式（23）)[5]，其對應冰載荷的2 個限值方程為：

限值Fα對應于平臺周圍海冰清理不良時的冰載荷上限值，此時海冰冰況嚴重，包括冰脊、重疊冰等復雜冰型；Fβ對應于平臺周圍海冰清理較好時的冰載荷上限值，此時海冰冰況較輕，平臺附近會有破冰船協助破冰，有很多開闊的水域。圖6示出了采用離散元法模擬的1.1 m 厚碎冰與Kulluk浮式平臺相互作用的過程，其中圖6(a)顯示的是碎冰在平臺前發生堆積，而平臺兩側的碎冰則繼續向前運動并在平臺后方形成一段開闊的水域(圖6(b))。從中可以看出，采用離散元法模擬的海冰冰情介于實測報告中所提的清理良好與清理不良這2 種工況之間。

圖6 離散元法模擬中破碎冰與Kulluk 浮式平臺間的相互作用Fig. 6 Interaction between pack ice and Kulluk floating platform simulated with DEM

圖7 所示為不同冰厚(hi=0.5～2.0 m) 下平臺總體冰載荷時程均值（Fi）與實測報告中數據（以噸計位）的對比。從中可以看出，由離散元法模擬得到的結果位于實測冰載荷范圍內，說明離散元冰載荷模型可以較好地適用于浮式結構的冰載荷分析。離散元數值結果顯示，在冰厚較小時黏結單元較易失效，大塊海冰碎裂成小尺寸冰塊。此時，海冰不宜在平臺迎冰面發生堆積轉而向平臺兩側滑移，進而導致冰力降低，可見，平臺冰載荷受冰厚的影響較大。海冰管理作業可以及時、有效地減弱平臺周邊的冰情，因此在實際平臺作業過程中，可以根據觀測的冰厚信息及時調配破冰船進行海冰管理作業，破碎其周邊海冰以提升結構安全性[5,26]。

圖7 離散元法模擬結果與實測數據的對比[26]Fig. 7 Comparison between the DEM results and field test data[26]

3 結 論

本文利用離散元法建立了用于浮式結構冰載荷分析的數值預報模型，并通過對比Kulluk 浮式平臺的實測數據，對該模型進行了驗證性分析，得到以下主要結論及未來需要研究的方向：

1） 基于CUDA-C 并行算法和Voronoi 分割算法圖形處理技術，可以較好地模擬冰區浮式結構?冰的相互作用過程，且該數值預報模型的有效性已通過了Kulluk 浮式平臺多年實測數據的檢驗，因此本文提出的數值模型可以為冰區浮式平臺冰載荷快速預報提供可行的研究手段。

2） 數值分析結果及實測數據均表明，冰厚對浮式結構冰載荷的影響較大，將直接影響平臺的作業安全性，因此，可通過海冰管理等方式對浮式結構的冰載荷進行控制。然而，國內有關冰區工程管理的經驗多集中于固定式平臺結構的風險分析與生命周期管理，如何建立科學調配破冰船隊、氣象預報及作業單元的管理體系成為浮式平臺安全運營的關鍵。

3） 冰區浮式結構的錨泊系統模擬技術大多采用線性彈簧來模擬張緊式錨泊系統，然而錨鏈自身的大變形和非線性特征也會對浮式結構運動響應特性帶來較大的影響，因此，下一步將采用更加符合真實錨鏈力學特性的力學模型。