半潛鉆井平臺輻射問題的高精度算法研究

王科，賀大川，陳彧超，施鵬飛（大連理工大學工業裝備與結構分析國家重點實驗室，工程力學系，遼寧大連116024）

半潛鉆井平臺輻射問題的高精度算法研究

王科，賀大川，陳彧超，施鵬飛
（大連理工大學工業裝備與結構分析國家重點實驗室，工程力學系，遼寧大連116024）

附加質量和阻尼力是決定半潛式鉆井平臺運動響應的重要參數，該研究采用直交邊界邊配置多重節點的新方法去除積分邊界的奇異性。研究中同時應用了高階邊界單元方法近似真實的空間三維復雜結構，依據波浪與結構物相互作用的輻射和繞射理論建立邊界元積分方程，對ISSC型半潛式平臺（Semi-Submersible）的附加質量和阻尼力進行了分析研究。研究表明多重節點與高階邊界元相結合的方法對復雜結構的水動力計算具有很高的精度。

半潛式鉆井平臺；附加質量；阻尼力；雙重和多重節點配置方法；高階邊界元

1 引言

半潛式鉆井平臺（Semi-Submersible）可以用于3 000m水深石油資源的勘探和開發。這種型式的鉆井平臺甲板面積大，性價比高，轉移安裝方便，能適應極端海況條件。隨著石油開采越來越走向深海，研究和開發適應深海石油勘探的新型半潛式鉆井平臺是海洋工程界關注的熱點。半潛式鉆井平臺是由立柱和浮筒組成的大型海工結構物，長度、寬度和深度三個方向的比尺適中，無法象FPSO系統以迎浪狀態工作。半潛式鉆井平臺必須直接承受各方向波浪的入射作用，這種工作方式一大缺點是垂蕩運動過大，為了減小垂蕩運動，增加平臺的附加質量和阻尼力，是方法之一。這也使得對三維實體鉆井平臺的附加質量和阻尼力的研究提出了很高的計算精度要求。

在半潛式鉆井平臺早期的研究過程中，限于當時的計算條件，通常把半潛平臺簡化為圓柱和浮筒的簡單組合形式，同時浮筒也簡化為一種水平橫置的水下圓柱。具有三維復雜的空間結構形式的半潛平臺的附加質量和阻尼力可以用多個獨立圓柱結構的附加質量和阻尼力來近似。由此針對這種相對簡單的結構形式發展了一系列的理論。代表性的工作有：當圓柱的直徑相對于波長較小時，可以忽略波浪繞射的影響，圓柱的波浪力可以用Morison等人[1]的公式計算。對于相對尺度較大的圓柱，Yeung[2]采用特征值展開法，研究了有限水深下直立圓柱的輻射問題。Garrison[3]用三維邊界元法研究了水面直立圓柱的輻射問題及運動響應。

隨著計算方法和計算理論的提高，許多波浪與結構物相互作用的研究工作都是基于求解格林函數的邊界元方法，而且這種趨勢在求解大型復雜結構方面的優勢越來越明顯。例如Chakrabarti等[4]同時應用波浪繞射理論和Morison方程的方法對桁架浮筒型半潛平臺的平板結構對附加質量的影響進行了研究；肖麗娜和索雙武[5]應用挪威船級社（DNV）的SESAM軟件計算作用在深水半潛平臺上的波浪載荷；楊立軍，肖龍飛，楊建民[6]應用HydroD軟件對半潛式平臺垂向運動的低頻響應性能進行了分析；，Matos等人[7]應用WAMIT軟件對深吃水的半潛平臺結構的二階奇模態共振運動響應進行了數值分析和物理模型試驗驗證。上述理論算法或者商用軟件的核心是邊界元或者面元方法和有限水深格林函數的結合，其中典型的面元方法有常數元方法（Hess和Smith[8]），線性元方法（Webster[9]）,以及高階邊界元方法（Lee，Newman[10]；Teng等[11]）。

線性元和常數元由于對物體的形狀采用面元來近似，對于復雜的結構需要的節點和單元數量都很多，限于邊界元法形成的矩陣為滿陣的特點，對于復雜結構計算速度慢，精度差，實際上已不適用。高階邊界元由于其形函數能精確的近似真實的物體幾何形狀，因此其所用的節點和單元數量都少，同時物理變量也可以采用與幾何近似一致的形狀函數。目前典型的高階邊界元方法主要有兩類，一類是采用樣條單元（Lee，Newman[10]），一個單元的幾何形狀和物理量由臨近單元的樣條函數組合逼近而成。另一類高階單元是采用在三角形單元或者四邊形單元上布置6個節點或者8個節點，形成只跟本單元節點有關的多節點高階等參單元（Teng等[11]）。

由于復雜的空間結構存在大量的由直交邊界面形成的直角邊，這類邊界在數值計算中存在著幾何奇異性，如何正確地處理這類奇異性，是高階邊界元方法能否準確預測海工結構物水動力特性的關鍵。這類邊界一個共同的特征是邊界物理量的法向導數不連續，樣條高階邊界元方法（Lee，Newman[10]）處理直角邊的措施是在極小的范圍內應用特殊樣條函數對直角邊進行光滑處理，實際上是用數學意義上逼近的虛擬邊界近似真實的幾何邊界。這種方法的主要優點是通過模型邊界平滑方法去除了幾何奇異性，缺點是模型中的幾何形狀不是真實的幾何形狀，而是由樣條函數逼近而成的。高階等參單元（Teng等[11]）方法由于單元各節點的物理量僅跟本單元有關，因此不存在把直角邊進行光滑處理的問題，其處理直角邊物理量的方法是將位于同一位置但屬于不同單元直角邊節點上的物理量取算術平均值。高階等參單元的一個主要缺點是對于復雜結構由于本身計算法向導數的誤差會在直角邊節點物理量取平均值的過程中不斷累積，使得幾何奇異性的特征在局部沒有得到有效改善，造成即便增加網格單元密度，也對計算精度提高不大的情況出現。

在本研究中，我們在高階等參單元的基礎上引進了一個新的改進方法-多重節點匹配方法，很好地緩解了復雜結構水動力計算的直角邊幾何奇異性問題，而且由于該方法是直接針對海工結構物的有限元模型進行處理，具有很好的工程應用前景。通過對ISSC型半潛平臺的有限元模型進行附加質量和阻尼力的計算分析，表明與常數單元方法相比，該方法計算精度更高，完全可以解決復雜的半潛平臺的水動力問題。

本文第二章闡述了應用三維有限水深格林函數求解波浪與結構物相互作用的輻射問題的基本理論和數值離散方法，闡述了多重節點方法的具體應用。第三章計算了ISSC型半潛鉆井平臺的附加質量和阻尼力,并對相關的影響因素進行了定量分析和討論，得到了具有參考價值的結果，結論及建議在第四章給出。

圖1 計算示意圖Fig.1 Calculation sketch

2 基本方程及公式

2.1 控制方程與邊界條件

一個任意形狀的三維物體位于水深為h的流體域內，定義笛卡爾坐標系如圖1所示，該坐標系相對于無擾動的自由表面和物體處于靜止狀態，原點取在自由表面上，x軸和y軸沿水平面方向，z軸垂直向上為正。在研究物體的輻射問題時，假定物體在靜平衡位置附近作頻率為ω的微幅振動并且產生向外傳播的波浪，物體的這種運動也稱為輻射運動。此時物體周圍的流體運動可由線性波浪理論描述。這樣所有的需要考慮的動力學變量均可表示為A e-iωt的實部，這里A（x,y,z）為復數幅值。物面上任一點,y, z）處的法向向量幅值Vn可寫成：

假定流體為無粘性、不可壓縮并且運動是無旋的理想流體，則物體作輻射運動時，周圍流體的運動可表示為輻射速度勢Φ=Re｝的梯度，即▽Φ。復數輻射勢φ滿足的流體運動控制方程可寫成如下形式：

并且滿足以下邊界條件：

φz-Kφ=0,z=0在自由表面上

（3）式中K=ω2/g為深水波數，g為重力加速度。

總輻射勢φ可分解為由縱蕩、橫蕩、垂蕩、橫搖、縱搖和首搖6個運動模態產生的輻射速度勢之和，即：

其中：ξj，j=1,2,…,6為物體的第j個方向運動的復數振幅分量，φj為單位振幅輻射勢，對應于物體的第j個運動模態。將方程（4）代入到方程（3）中可得輻射勢φj滿足的物面條件：

上述關于物體表面的輻射速度勢的邊界值問題（1-5）可通過格林定理在物面邊界上建立積分方程式的方法求解。其中物面上的輻射勢φj滿足的積分方程為：

式中：J0（）為零階貝塞爾函數，在公式（7）的積分中，積分變量k的積分路徑為k值在實軸上的實根上半部分，這樣的積分路徑使格林函數滿足波浪的遠場輻射條件。積分方程（6）可有效去除浮在水面上的物體波浪力計算中出現的不規則頻率現象（Lee，Newman和Zhu[12]），格林函數的計算方法可參照（Newman[13-14]）。

由物體的輻射運動模態引起的作用在物體上的總的流體力可由物面上各點壓力積分求得，物面上各點的水動壓力P由Bernoulli方程計算：

這里aij定義為附加質量力，bij定義為阻尼力。

2.2 數值求解

數值求解時，首先將物體表面進行離散化處理，積分方程（6）在物面上的積分可近似為在物面上一系列離散單元上的積分，這里的物面包括物體表面和內域水線面。在高階邊界元方法中，離散單元可采用6節點三角形單元（圖2）或者8節點四邊形單元（圖3），其形函數可表示為：

圖2 六節點形函數Fig2 6 node shape function

圖3 八節點形函數Fig.3 8 node shape function

則上述單元中任意一點的物理量及其導數，例如速度勢φ，可由單元節點值表示如下：

這里s表示單元上的節點數目，s=6或者8。則積分方程式（6）離散如下：

2.3 直角邊的多重節點配置方法

很多商用軟件例如ANSYS和GAMBIT等專業建模軟件能很好地構建精細而復雜的三維真實空間結構模型，并成功地應用在結構分析和強度計算中。但這種建模軟件提供的有限元數據模型不能直接為邊界元所用，一個重要原因是有限元模型中提供的直交物面形成的直角邊上的控制節點處，各物理量的法向導數不連續，邊界積分方程在該處存在幾何奇異性。高階等參邊界元方法如果直接應用有限元數據模型，通常方法是在直角邊控制節點處的物理量采用具有同一節點但不同單元數值結果的平均值。這種方法對于簡單形狀例如圓柱等，計算誤差較小，但對于半潛平臺這樣極為復雜的空間正交結構，直角邊的空間幾何奇異性非常明顯，靠取平均值和加密直角邊附近的網格單元，對計算精度提高不大。特別是當存在復雜空間曲面相交時，空間曲面法向導數由于網格質量剖分不佳而引起的計算誤差會導致計算錯誤。為解決上述問題，本研究提出了一種直角邊的雙重或者多重節點配置方法可完全去除積分邊界的幾何奇異性，提高計算結果的準確性。其主要思路是通過對半潛平臺的邊界積分方程進行整體分析，首先對具有直角邊的相交曲面進行合理分類，其次對有限元模型數據進行重新編號使得在原來直角邊上具有同一位置但屬于不同物面的任一個控制點處進行雙重或者多重節點配置，同時不同的重節點屬于唯一的高階單元，其方向導數也由該單元唯一確定。以直交圓柱為例，圓柱側面與底面相交直角邊，其上控制點原來的法向導數為側面與底面法向導數的平均值，大小和方向如圖4a所示，配置雙重或者多重節點后，相交直角邊上重節點的法向導數如圖4b所示。圖4c為配置多重節點的半潛平臺法向導數局部放大圖，其中相交直角邊包括立柱側面與底面相交線和浮箱側面與上下底面相交線。配置雙重節點實際上也增加了控制點數量，但由于重節點數量占整體節點比例不大，計算速度影響不明顯但卻顯著地提高了計算的精度。

圖4 a圓柱側面與底面相交直角邊控制點原法向導數方向Fig.4a Original normal vector along sharp edge between cylinder body and bottom surface

圖4 b配置雙重節點后相交直角邊上控制點新法向導數方向Fig.4b New normal vector along sharp edge between cylinder body and bottom surface

圖4 c配置雙重節點后半潛平臺法向導數局部分布圖Fig.4c Local view of normal vector distribution on Semi-Submersible based on new nodes relocationmethod

3 附加質量與阻尼力計算結果分析

圖5 ISSC型半潛平臺計算視圖Fig.5 Calculation views of ISSC type Semi-Submersible

應用上述格林函數理論和雙重及多重節點配置技術對半潛平臺的附加質量和阻尼力進行了研究。半潛平臺的附加質量和阻尼力受波浪的頻率和物體的形狀影響很大，本文以ISSC型半潛平臺為例研究該型平臺的附加質量和阻尼力，為波浪荷載及運動響應設計提供依據。該型式半潛平臺的空間模型結構及計算網格如圖5所示，半潛平臺基本參數如表1所示。在網格剖分中，采用6節點和8節點混合單元進行剖分。其中兩端立柱側面有228個節點，104個單元；中間立柱側面有224個節點，102個單元；圓筒頂面有874個節點，364個單元；浮筒側面有1 032個節點，416個單元；浮筒底面有861個節點，380個單元；合計共3 895個節點，1 674個單元。

表1 半潛平臺基本參數Tab.1 Basic parameters of Sem i-Submersible

圖6至圖13為半潛式采油平臺附加質量和阻尼力計算結果值，其中圖6至圖11為附加質量和阻尼系數陣主對角線上的值。圖12、13為考慮耦合運動時，橫蕩—橫搖，縱蕩—縱搖附加質量和阻尼力的值。圖中橫坐標為深水波浪頻率ω，縱坐標為附加質量和阻尼力aij，bij。當半潛平臺作平動時，由圖6、7、8可見，在波浪頻率趨于零時，縱蕩、垂蕩和橫蕩的附加質量都趨于不同的有限值，但此時垂蕩的附加質量最大；隨著波浪頻率的增加，縱蕩、垂蕩和橫蕩的附加質量都出現不同程度的波動，其中以縱蕩波動最大，垂蕩波動最小。由圖6可以看出，當波浪頻率ω<0.6時，縱蕩附加質量幾乎不變，而當波浪頻率增加到ω=1.0附近時，縱蕩附加質量會出現大幅度的波動，說明縱蕩附加質量對ω=1.0附近的波浪頻率比較敏感。橫蕩附加質量先是在ω<0.6范圍內小幅度波動，但是當波浪頻率進一步增加時，橫蕩附加質量迅速增加，達到極大值1.68×105kg，然后又快速減小，變化幅度較大，當ω=0.9時=7.0×104kg。垂蕩方向的附加質量相對變化較小，隨著波浪頻率的增加，呈波動式變化，當ω= 0.9時a33min=1.6×105kg。就阻尼力而言，縱蕩、橫蕩和垂蕩的阻尼力變化趨勢與各自的附加質量類似，但達到各自極值的頻率不同。當ω=0.95時7.0×104kg/s，當ω=0.8時，=7.0×104kg/s，當ω=0.9時3.0×104kg/s。但在ω→0和ω→∞時三者都趨于零。這主要是因為ω→0時物體不可能產生向外傳播的波浪，物體周圍流體的能量無法被波浪帶走，而阻尼力與物體向遠方傳播的波浪能量有關。

圖6 縱蕩運動附加質量與阻尼力Fig.6 Added mass and damping for surge

圖7 橫蕩運動附加質量與阻尼力Fig.7 Added mass and damping for sway

圖8 垂蕩運動附加質量與阻尼力Fig.8 Added mass and damping for heave

圖9 橫搖運動附加質量與阻尼力Fig.9 Added mass and damping for roll

圖10 縱搖運動附加質量與阻尼力Fig.10 Added mass and damping for pitc

圖11 艏搖運動附加質量與阻尼力Fig.11 Addedmass and damping for yaw

圖12 橫蕩—橫搖耦合運動附加質量與阻尼力Fig.12 Added mass and damping for coupled sway-roll

圖13 縱蕩—縱搖耦合運動附加質量與阻尼力Fig.13 Added mass and damping for coupled surge-pitch

橫搖、縱搖運動時附加質量與垂蕩運動的變化類似，都呈逐步震蕩下降的趨勢，但縱搖附加質量受垂蕩運動影響更明顯。阻尼系數都是在ω<0.3時很小，近似為0值，隨著ω逐漸增加，阻尼系數出現不同程度的波動，ω=0.6x=5.0×107kg/s，ω=0.85，=1.2×108kg/s。需要特別注意的是艏搖時附加質量和阻尼力的變化。在ω<0.4時，附加質量a66變化不大，大小為8.0×108kg·m2。隨著波浪頻率的進一步增加，附加質量迅速上升，在ω=0.75時達到極大值=1.05×109kg·m2，隨后在0.9<ω<1.05的區間，a66振動下降，在ω=1.0時，4.0×108kg·m2。艏搖的阻尼系數在ω<0.65范圍內變化很小，尤其是在ω<0.4時阻尼系數幾乎為0。當波浪頻率ω>0.7時，阻尼系數迅速增加，當ω=0.85時4.1×108kg/s。當0.85<ω<1.0時，阻尼系數振蕩下降。艏搖時附加質量和阻尼系數在ω較大的時候會出現振蕩變化，說明在艏搖時高頻短波對這種結構形式影響明顯。

圖12為橫蕩—橫搖耦合運動的附加質量a42和阻尼力b42。由圖可以看出，當波浪頻率ω較小的時候，波浪頻率的變化對附加質量和阻尼系數的影響不大，但是當ω>0.3時，附加質量和阻尼系數變化劇烈，這種強烈波動是由于入射波浪與半潛平臺中間水體相互作用的結果。

圖13為縱蕩—縱搖耦合運動的附加質量和阻尼力。由圖可以看出，附加質量和阻尼系數隨波浪頻率變化趨勢基本相同，但是阻尼系數的變化相對于附加質量的變化有些滯后。ω<0.25時，阻尼系數為0值，當ω繼續增加時，阻尼系數呈振蕩式變化，分別在ω=0.45和ω=0.85產生兩個極大值。附加質量在ω<0.4范圍內小幅度變化，隨著ω逐步增加，附加質量也會出現劇烈波動，在ω=0.75時達到最大值4.2×106kg·m。ω=0.6和ω=1.0時?2.0×106kg·m。

4 結論

本文對深海石油開采中的重要裝備:半潛式鉆井平臺的附加質量和阻尼力進行了計算分析，獲得了具有實際參考價值的荷載數據。研究發現：

（1）由于半潛式鉆井平臺空間三維幾何比尺適中，平動和轉動六個方向的附件質量和阻尼力受波浪的影響均很明顯。附件質量比阻尼力對半潛平臺的影響大。

（2）在求解速度勢的過程中，采用在直交的邊界邊布置雙重和多重節點的方法，消除了速度勢法向導數的幾何奇異性，能更準確地描述速度勢場的變化規律，節省了計算時間，提高了計算精度。

（3）采用高階邊界元方法和新的積分方程，消除了波浪不規則頻率的影響，計算結果具有很好的收斂性。

[1]Morison JR,O'Brien M P,Johnson JW,Schaaf SA.The force exerted by surface waves on piles[J].Petroleum Transactions,AIME,1950,18:149-157.

[2]Yeung RW.Added mass and damping of a vertical cylinder in finite-depth waters[J].Applied Ocean Research,1981,3 (3):119-133.

[3]Garrison C J.Hydrodynamics of large objects in the sea,Part II.Motion of free-floating bodies[J].Journal of Hydronautics,1975,9(2):58-63.

[4]Chakrabarti S,Barnett J,Kanchi H,Mehta A,Yim J.Design analysis of a truss pontoon semi-submersible concept in deep water[J].Ocean Engineering,2007,34:621-629.

[5]肖麗娜,索雙武.深水半潛平臺的環境載荷分析[J].船海工程,2012,41(1):85-87. Xiao Lina,Suo Shuangwu.Analysis of environmental load of semi-submerged platform in deep water[J].2012,41(1):85-87.

[6]楊立軍,肖龍飛,楊建民.半潛式平臺垂向運動低頻響應特性[J].海洋工程,2010,28(2):1-7. Yang Lijun,Xiao Longfei,Yang Jianmin.Low-frequency response of verticalmotions of a semi-submersible[J].The O-cean Engineering,2010,28(2):1-7.

[7]Matos V L F,Simos A N,Sphaier SH.Second-order resonant heave,roll and pitch motions of a deep-draft semi-submersible:Theoretical and experimental results[J].Ocean Engineering,2011,38:2227-2243.

[8]Hess JL,Smith A M O.Calculation of non-lifting potential flow about arbitrary three-dimensional bodies[J].JShip Res., 1964,8:22-44.

[9]Lee CH,Newman JN.Computation ofwave effects using the panelmethod[M].in:Chakrabarti,S.(Ed.),NumericalModeling in Fluid-Structure Interaction,MIT Press,2004.

[10]Eatock Taylor R,Chau F P.Wave diffraction theory-some developments in linear and non-linear theory[J].Journalof Offshore Mechanics and Arctic Engineering,1992,114:185-194.

[11]Teng B,Eatock Taylor R.New higher-order boundary elementmethods forwave diffraction/radiation[J].Appl.Ocean Res, 1995,17:71-77.

[12]Lee C H,Newman JN,Zhu X.An extended boundary-integral-equation method for the removal of irregular-frequency effects[J].International Journal for NumericalMethods in Fluids,1996,23:637-660.

[13]Newman JN.Algorithms for the free-surface Green function[J].Journal of Engineering Mathematics,1985,19:57-67.

[14]Newman JN.Distribution of sources and dipoles over a quadrilateral[J].Journal of Engineering Mathematics,1986,20: 113-126.

An efficient algorithm to study radiation problem of ISSC type Sem i-Submersible

WANG Ke,HE Da-chuan,CHEN Yu-chao,SHIPeng-fei
(State Key Laboratory of Structural Analysis for Industrial Equipment,Departmentof Engineering Mechanics, Dalian University of Technology,Dalian 116024,China)

Addedmass and damping are critical for hydrodynamic response of Semi-Submersible.This paper employs a new double ormultiple nodes relocationmethod on sharp edges to remove geometrical singularity. High order boundary elementmethod is utilized to approximate real three dimensional complicated structure, the added mass and damping of ISSC Semi-Submersible are investigated based on linear wave radiation and diffraction theory.It is found that the new algorithm combingmultiple node relocationmethod and high order boundary elementmethod is efficient for hydrodynamic problem of complicated structure.

Semi-Submersible;addedmass;damping;double and multiple nodes relocationmethod;higher-order boundary elementmethod

U674.38+1

A

10.3969/j.issn.1007-7294.2014.10.006

1007-7294（2014）10-1204-09

2014-04-21

國家重點基礎研究發展計劃資助（2013CB036101）；大連理工大學基礎科研費資助（DUT10LK43）

王科（1970-），男，大連理工大學工程力學系副教授，E-mail:kwang@dlut.edu.cn；賀大川（1986-），男，博士研究生。