廖發(fā)林，郭海燕，牛建杰，張 莉

(1. 中國海洋大學 工程學院，山東 青島 266100; 2. 江蘇科技大學 船舶與海洋工程學院，江蘇 鎮(zhèn)江 212003; 3. 山東科技大學 土木工程與建筑學院，山東 青島 266590)

內波作用下海洋立管動力響應數值模擬系統(tǒng)

廖發(fā)林1，郭海燕1，牛建杰2，張 莉3

(1. 中國海洋大學 工程學院，山東 青島 266100; 2. 江蘇科技大學 船舶與海洋工程學院，江蘇 鎮(zhèn)江 212003; 3. 山東科技大學 土木工程與建筑學院，山東 青島 266590)

針對目前內波對海洋立管作用研究缺乏的現狀，基于MATLAB開發(fā)出內波作用下海洋立管動力響應數值模擬系統(tǒng)NSDRI 1.0。該系統(tǒng)主要考慮內孤立波作用，基于mKdV理論求解內孤立波，使用改進的Morison方程計算作用在海洋立管上內孤立波的荷載，進而求解出立管的水平位移以及應力響應。該系統(tǒng)界面清晰友好，操作簡單，適用于頂張力立管。選取一定參數，使用該系統(tǒng)進行計算，并將計算結果與已發(fā)表的成果進行對比，發(fā)現結果較為吻合，從一定程度上驗證了系統(tǒng)的可行性，可對立管的設計提供借鑒。

MATLAB；用戶圖形界面；內孤立波；頂張力立管；動力響應

Abstract: Due to the lack of the study on marine riser under internal wave, Numerical Simulation System of Dynamic Response Under Internal Wave (NSDRI 1.0) based on MATLAB is developed. The effect of internal solitary wave is mainly considered and the theory of mKdV is used to solve internal solitary wave. It uses the modified Morison equation to calculate the load created by internal wave on marine riser. And then, we can get the displacement and stress of marine riser under the internal wave. The system is suitable for top tensioned riser and the interface of the system is so concise and easy to use. The results calculated by this system are consistent with other achievements which have been published. The feasibility of this system is verified to some degree. And the study results may provide some reference to the design of marine riser.

Keywords: MATLAB; GUI; internal solitary wave; top tensioned riser; dynamic response

我國南海內波活動頻繁，對南海油氣資源的開發(fā)有重要影響。1990年，在流花油田的延長測試期間，內波導致纜繩拉斷、船體碰撞[1]。1992年，我國南海石油公司東沙群島附近的半潛式鉆井平臺受到內孤立波的影響，在不到5分鐘的時間內擺動了110°[2]。1998年，中國科學院南海海洋研究所的“實驗三號”科考船在東沙群島南部海域工作時，受到內波的作用，科考船船體產生持續(xù)晃動[3]。眾多事實表明，內波已經成為制約海洋油氣開采和船體航行安全的重要因素，因此，我國南海海域的海洋工程結構物的設計和安裝須要考慮內波的影響。在內波與海洋工程結構物相互作用的研究方面，何景異等[4]研究了兩層流體中Spar平臺在內波作用下的運動響應問題，提出內波作用下Spar平臺運動響應的計算方法。尤云祥等[5]基于mKdV理論，結合Morison公式，給出了內孤立波與TLP平臺相互作用的理論模型，并研究了兩層流體中內孤立波與TLP平臺的相互作用。

立管是連接海底井口和海面平臺的主要連接構件。隨著南海油氣田開發(fā)進程的加快，作為海洋油氣開采中使用最多的頂張力立管，研究其在海洋環(huán)境荷載作用下的響應顯得尤為重要。然而，目前關于內波作用下頂張力立管動力響應的研究較少，如劉碧濤等[6]采用速度入口數值造波及數值耗散消波技術，利用FLUENT建立了內孤立波數值水槽，提出了頂張立管與內孤立波的相互作用問題的數值模擬方法；林忠義等[7]基于mKdV理論，結合改進的Morison公式及考慮內流的頂張力立管動力學方程，給出了內孤立波與頂張立管相互作用的理論模型；張莉等[8]在KdV-mKdV方程的基礎上，對內孤立波作用下的頂張立管進行數值模擬，對內孤立波作用下頂張立管極值響應進行了研究。并且現有的設計規(guī)范和商業(yè)軟件也都未考慮內波的作用。因此，著眼于工程實際應用，建立和發(fā)展工程使用的內波作用下海洋立管動力響應數值模擬系統(tǒng)是現階段急需解決的問題。

考慮內孤立波對頂張力立管的作用，基于MATLAB平臺開發(fā)出內波作用下海洋立管動力響應數值模擬系統(tǒng)NSDRI 1.0，用于分析計算頂張力立管的內波作用響應。隨著研究深入，將進一步升級開發(fā)出計算鋼懸線立管等柔性立管在內孤立波作用下的動力響應計算軟件，以期能在工程實際中應用。

1 計算原理

圖1 內孤立波與頂張力立管示意Fig. 1 Internal solitary wave and top tensioned riser

1.1頂張力立管

假定頂張力立管在初始時刻處于垂直狀態(tài)，頂張力立管的上端安裝有張緊器提供張力，并與水面平臺鉸接，下端與海底井口的萬向節(jié)鉸接，因此將其簡化為兩端受張力的簡支梁。建立坐標系如圖1所示。直角坐標系OXYZ用來描述內波的運動，OXY平面位于未擾動界面上。直角坐標系oxyz用來描述立管的動力響應，坐標原點位于立管下端的萬向節(jié)上，oz軸過立管的中心軸。假定內波初始位置到立管中心軸距離為d，內波傳播方向水平向右。兩個坐標系之間的轉換關系如下：

X=x+d，Y=y，Z=z-h2

1.1.1 控制方程的建立

根據達朗貝爾原理可以得到頂張力立管順流向的運動控制方程[9]：

式中：mr為立管單位長度質量，mf為管內流動的不可壓縮流體單位質量，ma為外部流體附加質量，C為結構阻尼，C′為水動力阻尼，V為管內流體流速，Te為有效張力，f′為流體作用力。

1.1.2 有限元離散

將立管從底端到頂端平均分成n個單元，每個單元長度為l，采用經典Hermit梁單元并利用Galerkin方法對方程(1)進行離散，可以得到整體坐標系下立管矩陣形式的控制方程：

采用Newmark-β法，在時域內對立管的控制方程進行逐步積分，得到立管各個單元每一時刻的位移、速度和加速度，用于進一步計算立管的位移和應力。

1.2內波

1.2.1 內波波面方程

設內孤立波振幅為η0，在自由面為剛蓋的情況下，使用mKdV理論對波面方程進行求解[10]：

其中，

η0為振幅，當h1h2時，取正值，此時兩層流體的界面上凸，稱為上凸型內孤立波，cm為內孤立波mKdV理論解的相速度。

波面方程給出后，由水動力學和流體連續(xù)方程可以得到內孤立波引起的水平方向速度和加速度[9]。

1)內孤立波引起的水平速度U：

其中，Δ=kX-cmt，以下公式中相同。

1.2.2 內波作用力

采用改進的Morison公式來模擬內孤立波對海洋立管的作用力[11]：

2 用戶界面設計

以MATLAB7.6.0為平臺，使用圖形用戶界面開發(fā)環(huán)境(GUIDE)，開發(fā)出內波作用下海洋立管動力響應數值模擬系統(tǒng)——NSDRI 1.0，系統(tǒng)簡潔清晰，使用簡單，便于數據修改與保存，實現了對內波作用下海洋立管動力響應計算結果的圖形、數值以及動畫的輸出。

在本系統(tǒng)的設計中，界面的布局設計采用自頂向下的設計方法，即先設計出主界面，再設計每個子界面。根據運算需要，將系統(tǒng)主要設置為三大模塊：參數輸入模塊；分析計算模塊以及結果輸出模塊，每個模塊分工明確，條理清晰，易于操作。

2.1參數輸入模塊

參數輸入模塊相當于前處理階段，主要作用是將用戶輸入的數值轉化成字符串并存儲起來，供分析計算模塊的運算使用。參數輸入模塊由四部分組成，包括立管參數、流體參數、內波參數以及波浪及浮體參數，界面如圖2～5所示。

圖2 立管參數輸入界面Fig. 2 Riser parameter

圖3 流體參數輸入界面Fig. 3 Flow parameter

圖4 內波參數輸入界面Fig. 4 Internal wave parameter

圖5 波浪及浮體參數輸入界面Fig. 5 Wave and Vessel motion parameter

2.2分析計算模塊

分析計算模塊是系統(tǒng)至關重要的組成部分，這個模塊決定著計算能否順利進行、計算速度的快慢以及計算結果的正確與否。在計算內波時，考慮的是非線性較為顯著(即振幅與水深比大于0.1)的內波，因此使用mKdV理論進行波面方程的計算。在時域內使用Newmark-β法對立管控制方程進行逐步積分，計算立管每一時刻的位移、速度和加速度，再利用改進的Morison公式計算內孤立波作用力，最后求解立管在內孤立波作用下的應力及位移。

2.3結果輸出模塊

2.3.1 圖形結果

圖形結果界面可以查看五個不同的圖像，分別是任意一點位移時程曲線、任意一點應力時程曲線、任意時刻位移曲線、最大位移包絡圖以及最大應力包絡圖。這五個圖為互斥關系，一次只能顯示一個圖像，可以根據需要逐一查看。點擊界面上的編輯圖形按鈕，可以進入MATLAB圖形編輯界面，在這里可以對圖形進行各種后處理。

2.3.2 數值結果

為了使計算結果更加準確地呈現出來，使用VideoSoft FlexArray Control控件，以表格的形式顯示數值結果[14]，數值結果包括了各個節(jié)點的最大位移與最大應力，并將所有節(jié)點中的極值加粗表示。數值結果只能查看不能編輯。

3 實例分析

3.1驗證計算

文獻[15]研究了頂張力立管在內孤立波與非均勻海流共同作用下的動力響應問題。利用NSDRI 1.0的系統(tǒng)，對文獻[15]的算例參數進行計算，得到立管順流向動力特性并與文獻對比，結果如圖6所示。

圖6 計算結果比較Fig. 6 Comparison of consequences

對比圖6(a)和圖6(b)可知，立管各時刻位移圖相對文獻[15]有一些差別，經分析，主要原因：一是使用莫里森方程計算單位長度立管上的順流向作用力時，文獻考慮了內波在順流向、橫向以及垂向產生的速度分量，而本文作為簡化，只考慮了內波在順流向產生的速度分量；二是對立管動力方程求解時，文獻使用的是顯式積分法中的四階Runge-Kutta法，而本文使用的是隱式積分法中的Newmark-β法。這兩個原因綜合導致了位移圖的差異。但同時可以看出，立管位移最大值以及到達最大值的時刻與文獻[15]吻合較好，從一定程度上說明了計算軟件的可行性。

3.2實例計算

1998年，中國科學院南海海洋研究所科考船“實驗三號”在進行錨定觀測時，獲取了內孤立波經過時的資料。蔡樹群等[16]通過對資料進行分析，得到內孤立波的相關特征，參照這些特征，選取的計算實例參數如下：水深 472 m，上下層水深分別為60 m和412 m，密度分別為1 025 kg/m3和1 028 kg/m3，內孤立波振幅為-75 m，立管長度為 472 m，頂部張力為748 kN，外徑 0.26 m，內徑 0.2 m，管內流體密度998 kg/m3，立管密度7 850 kg/m3，彈性模量E=2.07×1011N/m2，拖曳力系數CD=0.7，附加質量系數Ca=1，內孤立波距離立管中心軸的初始距離為1 250 m，利用NSDRI 1.0進行計算，結果如圖7～12所示。

圖7 任意一點位移時程曲線Fig. 7 Displacement-time curve

圖8 任意一點應力時程曲線Fig. 8 Stress-time curve

圖9 任意時刻位移曲線Fig. 9 Displacement curve

圖10 最大位移包絡圖Fig.10 Envelope of displacement

計算結果表明，在t=752 s，距水面56 m時，內孤立波作用使立管產生的水平位移達到最大，約為18.6倍的立管外徑；在距水面28 m處，立管產生的應力達到最大，約為127.4 MPa。

使用控制變量法，保持其余參數不變，計算頂張力為708 kN以及788 kN時立管的應力及位移，并將計算結果進行對比，結果如表1所示。

表1 不同頂張力下的計算結果Tab. 1 Results under different top tensions

由表1可知，隨著頂張力的增大，立管的位移減小，且變化明顯，最大位移的節(jié)點位置不變。立管的應力值由拉伸應力與彎曲應力組成。頂張力增大時，拉伸應力增大，由于位移減小，導致彎曲應力減小，在本文設定的頂張力范圍內，彎曲應力減小的程度大于拉伸應力增大的程度，因此應力值總體上呈現減小的趨勢，最大應力的節(jié)點位置不變。由此可見，頂張力對立管的作用相當顯著。因此，在設計立管時，應在保證結構安全的前提下結合成本因素，合理設置頂張力大小。

4 結 語

考慮內孤立波作用，基于MATLAB開發(fā)出《內波作用下海洋立管動力響應數值模擬系統(tǒng)——NSDRI》，該系統(tǒng)主要用于分析內孤立波作用下頂張力立管的動力響應，得到立管在內孤立波經過時的變形曲線及受力曲線。利用NSDRI對文獻[15]的數據進行計算并對比，結果基本吻合，說明了本系統(tǒng)的可行性。參照蔡樹群等人的數據進行計算，所得結果說明內孤立波對海洋立管有重大威脅；通過變量控制法進行計算，發(fā)現頂張力越大，立管在內孤立波作用下的響應越小。該系統(tǒng)界面清晰友好，使用簡單，可為實際工程提供參考。

[1] 陳景輝.南海流花11-1油田開發(fā)工程[J].中國海洋平臺，1996，11(1):43-45. (CHEN Jinghui. The development of liuhua 11-1 oilfield in South China Sea[J].China Offshore Platform, 1996,11(1):43-45. (in Chinese))

[2] BOLE J B, EBBESMEYER C C, ROMEA R D. Soliton currents in the South China Sea: measurements and theoretical modeling[C]//Proceedings of the Offshore Technology Conference. 1994: 367-376.

[3] 蔡樹群，甘子均. 南海北部孤立子內波的研究進展[J]. 地球科學進展，2001，16(2):215-219. (CAI Shuqun, GAN Zijun. Progress in the study of the internal soliton in the northern South China Sea[J].Advances in Earth Sciences, 2001,16(2):215-219. (in Chinese))

[4] 何景異，尤云祥，石強，等.兩層流體中內波作用下Spar平臺運動響應[J].海洋工程，2008，26(1):1-11. (HE Jingyi, YOU Yunxiang, SHI Qiang, et al. The motion responses of the spar platform induced by internal waves in a two-layer fluid[J]. The Ocean Engineering, 2008,26(1)1-11. (in Chinese))

[5] 尤云祥，李巍，時忠民，等.海洋內孤立波中張力腿平臺的水動力特性[J].上海交通大學學報，2010，44(1):56-61. (YOU Yunxiang, LI Wei, SHI Zhongmin, et al.Hydrodynamic characteristics of tension leg platforms in ocean internal solitary waves[J]. Journal of Shanghai Jiaotong University, 2010,44(1):56-61. (in Chinese))

[6] 劉碧濤，李巍，尤云祥，等.內孤立波與深海立管相互作用數值模擬[J].海洋工程，2011，29(4):1-7. (LIU Bitao, LI Wei, YOU Yunxiang, et al.Numerical simulation of interaction of internal solitary waves with deep-sea risers[J]. The Ocean Engineering, 2011, 29(4):1-7. (in Chinese))

[7] 林忠義，尤云祥，曲衍，等.內孤立波作用下頂部張緊式立管動力特性[J].海洋工程，2012，30(2):20-25. (LIN Zhongyi, YOU Yunxiang, QU Yan, et al. Dynamic response of the top tension riser in internal solitary waves[J]. The Ocean Engineering, 2012, 30(2): 20-25. (in Chinese))

[8] 張莉，郭海燕，李效民.南海內孤立波作用下頂張力立管極值響應研究[J].振動與沖擊，2013，32(10):100-104. (ZHANG Li, GUO Haiyan, LI Xiaomin. Extreme response of top tensioned riser under internal solitary wave in South China Sea[J]. Journal of Vibration an Shock,2013,32(10):100-104. (in Chinese))

[9] 張莉.深海立管內孤立波作用的動力特性及動力響應研究[D]. 青島：中國海洋大學，2013. (ZHANG Li. Study on dynamic characteristics and dynamic responses of deepwater riser under internal solitary wave[D]. Qingdao: Ocean University of China,2013. (in Chinese))

[10] MICHALLET H, BARTHELEMY E. Experimental study of interfacial solitary waves[J]. Journal of Fluid Mech., 1998,336:159-177．

[11] HAN S M. Non-linear coupled transverse and axialvibration of a compliant structure. Part2: Force vibration[J]. Journal of Sound and Vibration, 2000, 237(5)：875-900.

[12] 王樹青，梁丙臣.海洋工程波浪力學[M].青島：中國海洋大學出版社，2013:156-157. (WANG Shuqing, LIANG Bingchen. Wave mechanics for ocean engineering[M]. Qingdao: China Ocean University Press,2013:156-157. (in Chinese))

[13] 黃文昊，尤云祥，王旭，等.圓柱型結構內孤立波載荷實驗及其理論模型[J].力學學報, 2013, 45(5):716-727. (HUANG Wenhao, YOU Yunxiang, WANG Xu, et al. Internal solitary wave loads experiments and its theoretical model for a cylindrical structure[J]. Chinese Journal of Theoretical and Applied Mechanics, 2013, 45(5):716-727. (in Chinese))

[14] 李洪春,郭海燕,李效民.基于MATLAB的海洋立管渦激振動數值模擬系統(tǒng)研究[J].中國海洋大學學報，2010，40(S): 207-212. (LI Hongchun, GUO Haiyan, LI Xiaomin. Research on numerical simulation system of vortex-induced vibration of marine riser based on MATLAB[J].Periodical of Ocean University of China,2010,40(S): 207-212. (in Chinese))

[15] 蔣武杰, 林忠義, 尤云祥,等. 內孤立波與非均勻海流共同作用下頂張立管動力特性[J].水動力學研究與進展，2012，27(4)： 424-435. (JIANG Wujie, LIN Zhongyi, YOU Yunxiang, et al. Dynamic characteristics of a top tension riser under combined internal solitary wave and non-uniform current[J].Chinese Journal of Hydrodynamics,2012,27(4):424-435. (in Chinese))

[16] 蔡樹群，甘子鈞，龍小敏. 南海北部孤立子內波的一些特征和演變[J]. 科學通報，2001，46(15)：1245-1250. (CAI Shuqun, GAN Zijun, LONG Xiaomin. Characteristics and evolution of the internal solitary wave in the northern South China Sea [J].Chinese Science Bulletin,2001,46(15):1245-1250. (in Chinese))

Numerical simulation system of dynamic response of marine riser under internal wave

LIAO Falin1, GUO Haiyan1, NIU Jianjie2, ZHANG Li3

(1. College of Engineering , Ocean University of China , Qingdao 266100, China; 2. School of Naval Architecture & Ocean Engineering, Jiangsu University of Science and Technology, Zhenjiang 212003, China; 3. College of Civil Engineering and Architecture, Shandong University of Science and Technology, Qingdao 266590, China)

P751

A

10.16483/j.issn.1005-9865.2017.03.010

1005-9865(2017)03-0076-07

2016-03-28

國家自然科學基金(51279187)；山東省科技重大專項(2015ZDZX04003)

廖發(fā)林(1992-)，男，福建人，碩士研究生，主要從事海洋立管結構動力分析研究。E-mail：851268483@qq.com

郭海燕，女，教授，博士生導師。E-mail：hyguo@ouc.edu.cn