海底滑坡對置于海床表面管線作用力的CFD模擬

王忠濤，王寒陽，張宇

(1.大連理工大學 海岸和近海工程國家重點實驗室，遼寧 大連 116024；2.大連理工大學 巖土工程研究所，遼寧 大連 116024)

?

海底滑坡對置于海床表面管線作用力的CFD模擬

王忠濤1，2，王寒陽1，2，張宇1，2

(1.大連理工大學 海岸和近海工程國家重點實驗室，遼寧 大連 116024；2.大連理工大學 巖土工程研究所，遼寧 大連 116024)

海底管線是海洋工程中用于傳輸原油和天然氣等的重要通道，通常放置在海床表面或處于懸跨狀態。本文采用計算流體動力學CFD法模擬了不同沖擊角度下海底滑坡對置于海床表面的海底管線的作用，得到了管線所受的軸向荷載和法向荷載與滑坡沖擊角度之間的關系。同時分析了沿沖擊方向管線截面形狀與管線所受阻力之間的關系。對已有研究進行拓展延伸，豐富了不同工況下軸向阻力系數和法向阻力系數的計算成果，得出了海底滑坡對置于海床表面管線沖擊力的計算公式。

海底滑坡；海底管線；CFD模擬；非牛頓流體；沖擊截面形狀

1　引言

海洋中蘊藏著豐富的礦產資源，已成為資源開發較為活躍的區域。由于承受著各種地質環境條件的制約，因而海底通常處于不穩定狀態[1]。而海底滑坡則是海底中最常見、出現頻率最高的地質災害[2]。海底管線作為海上運輸資源的主要手段，其數量也隨著開發活動的愈加頻繁而日益增加[3]。海底滑坡能破壞海底石油平臺與海洋工程設施，甚至能夠引起諸如海嘯等自然災害，給經濟社會造成重大損失[4]。因此，評價海底的穩定性并進行滑坡-管線相互作用分析，不僅對資源開發防避相關地質災害具有重要的科學意義和應用價值，而且可為深水區油氣勘探平臺的穩定性分析和天然氣水合物勘探開發提供重要的科學依據[5]。

Zakeri等[6]在2007年通過常規模型試驗模擬了管線受垂直向水下滑坡的沖擊作用，研究了不同泥漿含量、泥漿流速對沖擊荷載的影響，建立了雷諾數與法向阻力系數CD-90間的數學表達關系式。隨后Zakeri[7]于2009年采用計算流體動力學CFD法模擬了不同沖擊角度下水下滑坡對懸跨管線的作用力，計算中采用Herschel-Bulkley模型，擴展了之前模型試驗的結果，提出了軸向阻力系數CD-0與雷諾數的計算公式，并修正CD-90的計算公式。Liu等[8]采用Power-Law模型，對不同沖擊角度下懸跨管線的受力進行計算，并通過對比土力學方法與流體力學方法的結果，提出法向阻力系數與沖擊角度相關，并進一步修正了CD-90和CD-0的計算公式。

Liu等、Zakeri等研究均模擬了懸跨狀態下的管線受力，尚未考慮不同沖擊角度滑坡對置于海床表面管線的受力。又由于底床粗糙度是研究河口海岸沉積物運移和水流結構的重要因素[9]，故在模擬時將采用粗糙底面以區別于Liu等的研究。另外李少華認為管線所受到的阻力與管線表面粗糙程度、管線沿沖擊方向的形狀等因素有關[10]，因此將進一步考慮沖擊截面形狀對管線阻力的影響。

本研究采用CFD法，模擬不同沖擊角度下海底滑坡對置于海床表面管線的作用力，討論沖擊角度、泥漿含量、滑移速度的影響，提出置于海床表面工況下管線法向阻力系數CD-90和軸向阻力系數CD-0的計算公式，并整理了CD-0與沖擊角度的相關性。同時討論了管線所受阻力與沖擊截面形狀之間的關系。

2　CFD數值模擬

2.1流變模型

海底滑坡從初始狀態失穩滑移，逐漸轉變成泥漿狀態、到形成混濁流，最終穩定沉積下來。剪切速率與切應力不滿足線性關系，故不能看作牛頓流體。在流體力學研究中，將海底滑坡的流動物質看成僅由一種物質組成的黏滯力很高的流體，可以用Bingham、Herschel-Bulkley或其他Power-Law等非牛頓流體的流變模型來描述。下面僅介紹本文模擬采用的Herschel-Bulkley流變模型。

Herschel-Bulkley模型被應用于Zakeri分析泥漿的組成特性中，本文模擬保持一致，流變方程為：

(1)

(2)

非牛頓流體雷諾數Renon-Newtonian計算公式[5]為：

(3)

式中，τ值可根據公式(1)計算獲得，不同含量泥漿的流變模型參數列于表1。本文泥漿體流動過程中計算得到的雷諾數Re均小于350，可視為層流[12]，與工程實際中雷諾數的范圍也較為接近[5]。

2.2CFD數值模型

CFD法是計算技術與數值計算技術的結合體，是流體試驗用數值模擬方法求解的過程[13]。通過ANSYS CFX 13.0軟件進行模擬計算，對模型進行前、后處理，得到所需要的數據。

本文采用水和泥漿二相自由表面流，且泥漿組成與Zakeri模型試驗試樣保持一致，具體見表1。共模擬7個不同的沖擊角度，包括0°、15°、30°、45°、60°、75°、90°。其中，90°是指垂直于管線軸線方向的沖擊，0°是指平行于管線軸向方向的沖擊。由Liu等[8]研究結果得到，管線直徑的大小對管線受力并無顯著影響，本研究中為與Zakeri試驗數據對比，管線外徑取為25 mm，管線具體的平面布置情況如圖1所示。模型的邊界條件與Zakeri[7]、Liu等[8]相同，管線采用銅管，銅管表面應用粗糙無滑移邊界條件，粗糙率ks取0.001 5 mm。四周邊壁采用光滑壁面條件，底面采用粗糙無滑移邊界條件，ks取0.5 mm。銅管與底面的粗糙率均與Zakeri的模型試驗保持一致。網格由ANSYS ICEM CFD 13.0繪制，其中一個算例網格見圖2所示。圖2a左側藍色的矩形為泥漿進口，在模擬開始后泥漿以預先設定的恒定速度流入，泥漿出口在右端，邊界條件設為開放式邊界條件。在管線周圍進行網格加密，加密區網格最大尺寸為7.5 mm。由于圓管表面處的流速變化較大，故在表面設置邊界層，邊界層的厚度取為3 mm，共5層。其他區域的網格最大尺寸為50 mm。

表1　泥漿的組成和流變特性[6]

模擬中泥漿流動采用連續流體、層流、自由表面流，而水流動產生的湍流采用雷諾數平均法RANS中的雙方程k-ε模型，即速度與長度分開求解的傳輸模型的一種，適合絕大多數的工程湍流模型，其中k為湍動能，定義為速度波動變化量，單位為m2/s2；ε為湍動能耗散，即速度波動耗散的速率[13]。Liu等[8]等利用低流速產生較小雷諾數，并考慮流速低時泥漿對管線的覆蓋問題，限定最小流速為0.3 m/s。本文模擬中，流速低到0.2 m/s甚至更小時，泥漿仍能完全覆蓋管線，實現了更低流速的工況。為了消除或減小邊壁對管線受力的影響，管線所計算的陰影部分至少距離模型邊壁4倍管徑長度，且管線并沒有完全延伸至另一側壁面[6—7]，如圖1所示。假定在沖擊過程中，管線不發生變形。

圖1　各角度模型尺寸平面圖Fig.1　Plan view of domain dimensions and criteria of model

圖2　沖擊角度為60°的網格及局部放大圖Fig.2　The local enlarged grid of impact angle of 60°

3　CFD數值結果及分析

3.1流體力學方法基礎

根據CFD模擬結果，Zakeri[7]提出管線處于懸跨狀態工況下計算管線軸向阻力系數和法向阻力系數的計算公式，如下：

(4)

(5)

Kundu和Cohen[14]提出軸向阻力與法向阻力可分別表示為：

(6)

(7)

式中，FD-90為與管線軸向方向垂直的阻力(法向力)，FD-0為與管線軸線方向平行的阻力(軸向力)，A為計算長度管線的投影面面積，θ為流動方向與管線軸線方向的夾角。

3.2流體力學方法分析

計算得到作用于管線上的力Fx、Fz分解可得任意角度管線軸向和法向的沖擊力計算公式 ：

(8)

(9)

其中，對90°管線，FD-90=Fx，FD-0=-Fz。在計算中，垂直沖擊時實際測得z方向上所受的阻力Fz的值非常小，即Fz≈0。所以，當沖擊方向與管線軸線方向垂直時，只需計算法向阻力。同理，0°管線只需計算軸向阻力。

軸向阻力系數CD-0、法向阻力系數CD-90與雷諾數的關系見圖3、圖4。可以看出CD-90與Re以及CD-0與Re的關系曲線的數據點均比較分散。原因在于獲取數據點的管線與泥漿流動方向夾角不同，即沖擊角度對法向阻力系數CD-90和軸向阻力系數CD-0均有一定影響。數據結果列于表2。

圖3　法向阻力系數CD-90與雷諾數Re的關系曲線Fig.3　Relationship curves of Reynolds number and normal drag coefficient

圖4　軸向阻力系數CD-0與雷諾數Re的關系曲線Fig.4　Relationship curves of Reynolds number and axial drag coefficient

表2　各個角度管線的阻力系數與雷諾數的關系式

考慮泥漿流動方向對阻力系數的影響，將每個阻力系數乘(除)以關于沖擊角度θ的函數，可得到：

(10)

(11)

(12)

(13)

將式(12)和式(13)代入式(6)和(7)中，可得管線阻力的表達式為：

(14)

(15)

分別代入Herschel-Bulkley模型公式，則管線阻力表達式可改寫為：

(16)

(17)

式中，τy可由表1得到。可見，作用在管線的力只與泥漿的基本特性、沖擊速度以及管線的幾何尺寸有關。

圖6　總軸向阻力系數CD-0與雷諾數Re的關系曲線Fig.6　Relationship curve of total axial drag coefficient and Reynolds number

4　沖擊截面形狀對阻力影響

4.1沖擊截面的高寬比

在工程中，顆粒的形態通常用顆粒的“長寬比α”來表示[15]。長寬比的定義式為：

(18)

式中，L為顆粒投影的最大弦長，W為與顆粒投影面積相等的與最大弦長相對應的短邊長。

已有研究表明[15]，顆粒的長寬比能夠有效地描述顆粒的形狀。因假設在沖擊過程中管線不發生變形，引入截面形狀作為影響管線所受阻力的因素進行研究。本文借鑒“長寬比”概念，提出管線截面“高寬比”來表征泥漿在不同角度沖擊管線工況下截面形狀的參數。定義高寬比λr為沖擊截面在豎直方向上的長度與水平方向的長度的比值，即：

(19)

式中，V表示管線截面在豎直方向上的長度值，H表示管線截面在水平方向上的長度值，見圖7。

圖7　沖擊角度為45°的截面形狀對比Fig.7　Comparision with the impact angle of 45°

本研究中，沿泥漿流動方向，不同角度管線的截面為橢圓形，沖擊截面的高寬比λr即為橢圓截面長短兩軸之比。當λr=1時，即泥漿沖擊方向與管線軸線垂直，沖擊截面為圓形。當λr<1時，即當管線與流動方向的夾角為θ時，沖擊截面豎向方向上的長度值小于水平方向的長度值。隨著沖擊角度的增加，豎向和水平向兩方向上截面長度的比值越來越小，即高寬比λr隨著沖擊角度增加而減小。當λr>1時，即截面豎向方向上的長度值大于水平方向的長度值，對應實際工程中豎向傾斜的管線。另外，由于0°管線的沖擊截面比較特殊，故沒有計入本次截面形狀的討論中。

4.2阻力系數與截面形狀關系

考慮不同沖擊角度的管線，將截面的高寬比進行分別表示。沖擊截面水平方向的長度即為橢圓的長軸，可以表示為：

H=2a，

(20)

式中，a值可由沖擊角度θ與管線直徑D表示出來。不同角度管線的具體截面的大小可見圖8，利用三角函數關系進行推導可得到關系式：

(21)

圖8　沖擊角度為θ的管線平面圖Fig.8　The pipeline plan of the impact angle θ

結合式(20)、(21)得到截面水平方向的長度值H，截面豎向的長度值即為管徑D，即V=D。因此，整理以上各式并利用定義式(19)得到不同沖擊角度管線截面的高寬比為：

(22)

利用式(22)對本文模擬數據的結果進行分析，包括流體力學方法中的法向阻力系數CD-90與軸向阻力系數CD-0。結合式(12)、(13)，可以得到法向阻力系數與軸向阻力系數與高寬比關系式為：

(23)

(24)

當沖擊截面為純圓形時，對應泥漿垂直管線軸向方向沖擊工況，根據式(24)計算得到CD-0為0，即管線軸向受力為0，這與模擬計算結果及分析結論一致。沖擊角度為0°的管線即為λr=0的情況。因此，軸向阻力計算式(24)中，給定高寬比的適用范圍為，0≤λr≤1。

4.3豎向傾斜管線截面形狀分析

在實際的海洋工程中，海底管線可能由于海床表面的凹凸不平等原因而處于豎向傾斜的狀態，對應上述高寬比λr>1情況。將兩種取值范圍(λr>1和λr<1)的高寬比截面形狀進行對比(圖7)，圖7a、7b分別表示管線在豎向平面內傾斜角度45°時的截面形狀與沖擊角度為45°的水平管線的截面形狀。

當管線為豎向傾斜時，截面形狀變成豎向橢圓形，同水平管線分析方法給出此時高寬比表達式為：

(25)

流體力學分析方法具體計算結果見圖9。得到法向阻力系數CD-90與Re的數據擬合結果為：

(26)

式(26)的擬合度R2為0.998 9。

圖9　豎向傾斜角度為45°的阻力系數關系線Fig.9　The relationship of 45° in vertical inclination

當管線豎向傾斜的角度為45°時，可由式(25)計算得到高寬比的值為λr=1.414 2。利用高寬比與阻力系數計算式(23)，代入高寬比，得到相應的系數為15.834 8，這一數值與模擬得到的15.34非常接近，誤差僅在3%左右。由于要實現分析軸向阻力系數的豎向橢圓截面比較復雜，并未分析討論，因此本文給出式(24)的適用范圍僅為0≤λr≤1。

5　結論

本文針對置于海床表面工況的管線，利用計算流體動力學CFD法，采用與Zakeri的試驗條件一致的粗糙不滑移的底面邊界條件，進行法向阻力系數與軸向阻力系數的研究，討論了滑坡沖擊角度、沖擊速度以及沖擊截面形狀關于阻力系數的影響規律，得到如下結論：

(1)得到管線放置在海床表面的工況下，法向阻力系數CD-90與軸向阻力系數CD-0的表達式。由表達式可見，作用在管線的力只與泥漿的基本特性、沖擊速度以及管線的幾何尺寸有關。

(2)提出截面“高寬比λr”用于表征沖擊截面形狀，并得到了法向阻力系數CD-90與軸向阻力系數CD-0隨高寬比λr變化的表達式，并分析討論了表達式的的適用范圍。

(3)通過CFD算例驗證，針對豎向傾斜45°工況的管線，數值計算得到的結果與本文法向阻力系數CD-90隨高寬比λr表達式中相應高寬比下的結果僅相差3%，證明了阻力系數隨高寬比變化的表達式是有效的。

[1]杜軍， 李培英， 李萍， 等. 基于海洋災害地質評價基礎上的我國近海海底穩定性區劃[J]. 海洋學報， 2014， 36(5):124-129．

Du Jun， Li Peiying， Li Ping， et al. The seabed stability zonation based on the marine geohazards evaluation in China[J]. Haiyang Xuebao， 2014， 36(5):124-129．

[2]徐元芹， 劉樂軍， 李培英， 等. 我國典型海島地質災害類型特征及成因分析[J]. 海洋學報， 2015， 37(9):71-83．

Xu Yuanqin， Liu Lejun， Li Peiying， et al. Geology disaster feature and genetic analysis of typical islands，China[J]. Haiyang Xuebao， 2015， 37(9):71-83．

[3]陳東景， 李培英， 劉樂軍， 等. 海底地質災害對社會經濟發展影響的特點與趨勢[J]. 海洋開發與管理， 2010， 27(6):80-84．

Chen Dongjing， Li Peiying ， Liu Lejun， et al. On the characteristic and developing trends of submarine geological hazard impact on society and economy of China[J]. Ocean Development and Management， 2010， 27(6):80-84．

[4]李細兵， 李家彪， 吳自銀， 等. 北呂宋海槽深海滑坡沉積及其分布特征[J]. 海洋學報， 2010， 32(5):17-24．

Li Xibing， Li Jiabiao， Wu Ziyin， et al. Deepwater landslide deposition and distribution characteristics in the northern Luzon Trough[J]. Haiyang Xuebao， 2010， 32(5):17-24．

[5]張丙坤， 李三忠， 夏真， 等. 南海北部深水區新生代巖漿巖分布規律及其與海底地質災害的相關性[J]. 海洋學報， 2014， 36(11):90-100．

Zhang Bingkun， Li Sanzhong， Xia Zhen， et al. Distribution of Cenozoic igneous rocks and its relation to submarine geological hazards in the deepwater area of the northern South China Sea[J]. Haiyang Xuebao， 2014， 36(11):90-100．

[6]Zakeri A， H?eg K， Nadim F. Submarine debris flow impact on pipelines——Part Ⅰ: Experimental investigation [J]. Coastal Engineering， 2008， 55(12):1209-1218．

[7]Zakeri A. Submarine debris flow impact on suspended (free-span) pipelines: Normal and longitudinal drag forces[J]. Ocean Engineering， 2009， 36(6): 489-499．

[8]Liu Jun， Tian Jianlong， Yi Ping. Impact forces of steady submarine landslides on suspended pipelines[J]. Ocean Engineering， 2015， 95(2):116-127．

[9]柏秀芳， 龔德俊， 李思忍， 等. 根據流速剖面估計海底粗糙長度的研究[J]. 海洋學報， 2008， 30(4):176-180．

Bai Xiufang， Gong Dejun， Li Siren， et al. Estimating a seabed roughness length from current profiles[J]. Haiyang Xuebao， 2008， 30(4):176-180．

[10]李少華. 海洋平臺水動力系數反演方法及試驗研究[D]. 長沙：湖南大學， 2010．

Li Shaohua. Identification method and experimental analysis of hydrodynamic coefficients for ocean platform[D]. Changsha：Hunan University， 2010．

[11]Locat J. Normalized rheological behaviour of fine muds and their flow properties in a pseudoplastic regime[C]//Proceedings of the First International Conference. ASCE， New York， 1997: 260-269．

[12]陳懋章. 黏性流體動力學基礎[M]. 北京:高等教育出版社， 2004．

Chen Maozhang. Fundamentals of Viscous Fluid Dynamics[M]. Beijing：Higher Education Press， 2004．

[13]謝龍漢， 趙新宇， 張炯明. ANSYS CFX 流體分析及仿真[M]. 北京:電子工業出版社， 2012．

Xie Longhan， Zhao Xinyu， Zhang Jiongming. ANSYS CFX fluid Analysis and Simulation[M]. Beijing：Publishing House of Electronics Industry， 2012.

[14]Kundu P K， Cohen I M. Fluid Mechanics[M].Amsterdam：Elsevier Academic Press， 2004.

[15]高國瑞. 近代土質學[M]. 南京:東南大學出版社， 1990.

Gao Guorui. Modern Soil Science[M]. Nanjing: Southeast University Press， 1990．

陸雪駿， 程和琴， 周權平，等. 強潮流作用下橋墩不對稱“雙腎”型沖刷地貌特征與機理[J]. 海洋學報， 2016， 38(9):118-125， doi:10.3969/j.issn.0253-4193.2016.09.012

Lu Xuejun， Cheng Heqin， Zhou Quanping， et al. Features and mechanism of asymmetric double-kidneys scoured geomorphology of pier in tidal estuary[J]. Haiyang Xuebao， 2016， 38(9): 118-125， doi：10.3969/j.issn.0253-4193.2016.09.012

CFD numerical analysis of submarine landslides impact on laid-on-seafloor pipeline

Wang Zhongtao1，2， Wang Hanyang1，2， Zhang Yu1，2

(1.StateKeyLaboratoryofCoastalandOffshoreEngineering，DalianUniversityofTechnology，Dalian116024，China; 2.InstituteofGeotechnicalandGeomechanicsEngineering，DalianUniversityofTechnology，Dalian116024，China)

Submarine pipelines， in state of suspended or laid on the seafloor， were important transmissible channels of crude oil and natural gas in marine engineering. For the laid-on-seafloor pipelines， in order to get the axial and normal drag force of the pipeline under different impact angles， the CFD method was used to simulate the action of submarine landslide affecting on the submarine pipeline. Further， it was analysed that the relationship between the cross section shape of the pipeline in the impact direction and the drag force of the pipeline.The research extended the predecessors’ work， enriched the calculated results of the longitudinal and the normal drag coefficient， and obtained the formula for estimating the impact force of the debris flow on the seabed pipelines．

submarine landslide; laid-on-seafloor pipeline; CFD method; non-newtonian fluid; cross section

2016-04-14；

2016-05-09。

國家自然科學基金項目(41572252)；高等學校博士學科點專項科研基金資助課題(20120041130002)。

王忠濤(1974—)，男，遼寧省大連市人，副教授，主要從事海洋巖土工程等方面的教學和科研。E-mail：zhongtao@dlut.edu.cn

P642.22

A

0253-4193(2016)09-0110-08