穩態段塞流對鋼懸鏈立管動態響應的影響分析

劉德鵬， 艾尚茂， 孫麗萍， 賈魯生

(1.哈爾濱工程大學 船舶工程學院，黑龍江 哈爾濱 150001； 2.中海油研究總院有限責任公司 工程研究設計院， 北京 100028)

海洋立管的作用是將油氣資源從海底輸送到海上平臺，是海洋油氣行業中的關鍵水下設備。立管內部通常為混合了多種氣體與液體的多相流體，該種流動情況被稱為段塞流。由于立管結構的幾何非線性[1-2]以及段塞流不均勻的質量分布，使得輸送段塞流的立管動態響應預報十分復雜。

考慮到流固耦合方法計算復雜且耗時較長，一般在立管設計中需要建立合適的段塞流數學模型來開展海洋立管總體響應預報。段塞流具有速度和密度隨時間變化的特征[3]，當內部流動速度一致且每部分段塞長度相等時，即為穩態段塞流。早期研究往往將立管內部流動簡化為均勻流考慮[4]。Patel 等[5]提出了一種密度隨位置和時間的變化而呈正弦變化的穩態段塞流模型。當前國內對于段塞流的研究主要集中在嚴重段塞流的成型機理以及影響[6]，穩態段塞流模型相關研究較少。

Shabana[7]提出了一種有限元建模方法，即絕對節點坐標法(absolute nodal coordinate formula, ANCF)。該方法中節點坐標在全局坐標系中定義，有效地避免了復雜的坐標變換。此外，ANCF采用斜率向量代替角坐標的方式來描述結構的局部旋轉[8]，可以準確地描述剛體位移與柔性變形之間的耦合關系。通過將ANCF與任意拉格朗日-歐拉(arbitrary Lagrangian-Eulerian, ALE)描述結合，可以精確高效地描述具有大變形特征的介質[9]。

為預測內部段塞流引起的立管結構振動，本研究在ALE-ANCF框架下建立了海洋立管與內部段塞流的有限元分析模型。以鋼懸鏈立管為例，分析了內部穩態段塞流對立管動態響應的影響。

1 ALE-ANCF有限元理論

1.1 控制體ALE-ANCF描述

ALE采用物質坐標來描述控制體的邊界和沿介質的物質流動，介質上任一點的物質坐標為該點與起始點之間的介質長度。由于控制邊界s1和s2是隨時間變化的，控制區間的位置及長度隨之變化。為便于后續控制區間的描述，定義局部坐標ξ為[10]：

(1)

在任意時刻t時控制體內部任意物質點的空間坐標為：

r=r(ξ,t)

(2)

r=Nw,r′=N′w

(3)

(4)

1.2 穩態段塞流單元

采用Patel等[5]提出的穩態段塞流模型來模擬內部流體密度的時空變化，如圖1所示，其中t、n分別為管道切向與法向，Di為立管內徑，vi為段塞速度，ls為段塞長度。該數學模型中段塞流密度分為恒定分量與時變分量：

ρf=ρ0(1+ηei(ks-ωt+ψ))

(5)

式中：ρf和ρ0分別為為內部流體密度和恒定密度分量；k=2π/ls為段塞流波數；η為波動幅值系數；ω為段塞圓頻率；ψ為初始相位，用于指定初始時刻的段塞流質量分布。

圖1 段塞流模型Fig.1 Configuration of slug flow

段塞流單元可以在管道結構的切向上任意移動，形函數僅約束其在管道結構法向的位移。假定內部流動為穩態段塞流，段塞速度vi=ω/k，此時對于任意段塞流單元，物質的流入速度與流出速度相等。段塞流單元邊界的物質坐標為段塞速度的函數：

s1=-vit,s2=Le-vit

(6)

式中Le為單元長度。將式(6)代入式(4)，即可得到段塞流單元內部任意點的速度與加速度。

段塞流單元的平衡方程為：

(7)

式中：Qi為廣義慣性力；Mf為單元質量陣。根據達朗貝爾原理，Qi與Mf可由虛功方程得到：

(8)

(9)

1.3 立管單元

(10)

式中：Mr、Qf與Qe分別為管道單元的質量矩陣、廣義外力以及廣義粘彈性力。根據達朗貝爾原理，Mr、Qf與Qe由虛功方程得：

(11)

(12)

(13)

式中：ρr為立管材料密度；Ar為立管截面積；f為單元外力矢量；ε和κ分別為單元的軸向應變和曲率，定義為：

(14)

單元的剛度矩陣可由Qe對w直接求導得到[11-12]:

(15)

式(15)右端兩部分分別表示剛度矩陣的拉伸分量和彎曲分量。

施加在立管單元上的外力包括重力、浮力、拖曳力和慣性力，后兩項水動力通過莫里森公式計算：

(16)

(17)

1.4 數值求解方法

靜態平衡方程的求解采用Newton-Raphson方法，迭代形式為[13]：

w(n+1)=w(n)+Δw(n)

Δw(n)=[K(w(n))](Qe-Qf+Qr)

(18)

式中：n為迭代次數，當Δw滿足：

‖Δw(n)‖≤w0

(19)

時迭代終止。

在t+Δt時刻，輸送段塞流的立管的平衡方程為：

(20)

式(20)可采用Newmark法進行求解，在該方法中，t+Δt時刻的位移、速度與加速度為：

(21)

將式(21)代入t+Δt時刻的平衡方程，可以得到遞推公式：

(22)

其中：

2 穩態段塞流對立管的影響分析

2.1 立管與段塞流參數

以鋼懸鏈立管(SCR)為例，分析穩態段塞流對立管的影響。圖2給出了鋼懸鏈立管模型的示意圖，該立管總長1 260 m，工作水深為900 m，水平跨距700 m，頂部懸掛點位于水面下10 m。立管外徑D=0.4 m與內徑Di=0.36 m，空氣中單位長度的質量為189.1 kg/m，法向拖曳力系數與慣性力系數分別為1與2。

圖2 鋼懸鏈立管Fig.2 A schematic presentation of the SCR model

鋼懸鏈立管內部輸送流體為API 48原油和甲烷氣體，內部流體平均密度選為632.135 kg/m3，液相密度與氣相密度分別為790 kg/m3與0.675 kg/m3，波動幅值系數為0.25，段塞速度為2.5 m/s，段塞長度Ls=80 m。

為了探討穩態段塞流對鋼懸鏈立管的影響，本研究還對輸送內部均勻流的立管進行了模擬，該內部均勻流模型具有與穩態段塞流相同的平均密度與流動速度，2種內部流動情況下的內部流體總質量相同。

2.2 靜態分析

在分析穩態段塞流對立管的動態響應影響之前，先簡要分析內部均勻流與段塞流情況下立管的靜態特性。依據單個段塞的長度，將鋼懸鏈立管劃分為315個等長單元。采用土壤剛度來描述立管與海床之間的相互作用[14]。圖3給出了2種內部流動情況下的立管張力與彎矩特性曲線。

通過觀察圖3可以發現2種內部流動情況下的立管張立曲線幾乎完全重合，這說明段塞流的不均勻密度分布對立管的張力影響很小。2種內部流動情況下的彎矩曲線具有相似的趨勢，但是內部穩態段塞流情況下的彎矩曲線呈現出明顯的局部波動。其原因是鋼懸鏈立管的截面抗彎剛度極大，導致很小的曲率變形就會引起很大的彎矩變化。

圖3 立管靜態張力與彎矩Fig.3 Static tension and bending moment along the SCR

2.3 動態分析

在動態仿真中，通過對立管的頂部節點施加X方向的強迫位移來模擬頂部平臺的波頻運動。激勵位移rx可以表示為：

(23)

式中：頂部激勵位移幅值λ和周期Τm分別取為0.5 m和9 s。動態仿真時長設置為2 000 s, 時間步長為0.05 s。

圖4給出了鋼懸鏈立管頂部張力的時歷曲線。可以看出在內部穩態段塞流的情況下，立管頂部張力的變化呈現出雙頻疊加的特性。與內部均勻流情況相比，穩態段塞流情況下頂部張力的波動幅值有所增大，2種情況下頂部張力幅值之間的差值為3.78 kN，僅為內部均勻流情況下頂部靜態張力的0.27%，這說明穩態段塞流對鋼懸鏈立管的動態張力影響很小。

圖4 頂部張力時歷曲線Fig.4 Time histories of the top tension

圖5給出了2種內部流動情況下沿管長彎矩分布時歷圖。比較圖5(a)、(b)可以發現彎矩在穩態段塞流情況下的變化幅值明顯大于內部均勻流情況。立管觸地段的彎矩在2種流動情況下的變化趨勢相近，而在立管懸掛段，彎矩的變化頻率與段塞流頻率相近。

圖5 沿管長彎矩分布時歷Fig.5 Distribution of bending moments along the arc length of the riser

為更明顯地觀察彎矩變化規律，分別選擇觸地點附近的節點35(距錨點140 m處)與懸掛段彎矩最大處附近的節點75(距錨點300 m處)作為觀測點，圖6給出了上述2節點處彎矩變化(距錨點140 m處)的時歷曲線。可以看出，穩態段塞流情況下上述兩節點的彎矩變化均呈現出雙頻疊加的特性，節點35處的彎矩變化以頂部激勵頻率為主，而節點75處的彎矩變化以穩態段塞流的頻率為主，這說明穩態段塞流對懸掛段彎矩的影響更為顯著。

圖6 節點35與75處彎矩時歷曲線Fig.6 Time histories of bending moment at node 35 and 75

位移計算結果如圖7所示。結果顯示，段塞流對節點35的X向位移幾乎沒有影響，而其Z向位移會隨著段塞流的不均勻質量移動發生波動，2種內部流動情況下的最大Z向位移值相差0.01D。穩態段塞流情況下節點75的X與Z向位移均出現了雙頻疊加的特性，與均勻流動情況下的最大位移值相差分別為0.07D和0.08D。上述分析表明，穩態段塞流對立管觸地段影響極小，但是會改變其Z向的運動特性。而對于立管懸掛段，穩態段塞流的影響顯著，這也是該段彎矩變化特性發生改變的原因。

圖7 節點35與75處位移時歷曲線Fig.7 Time history of displacements at node 35 and 75

3 段塞長度參數分析

為分析段塞長度對立管動態響應的影響，對輸送不同長度段塞流的鋼懸鏈立管進行了時域仿真。段塞長度與管徑之間的比值選取為Ls/D=100，160，220，280。內部流體平均密度與流動速度保持不變圖8和圖9分別給出了不同段塞長度下的節點35與節點75處的位移時歷曲線。

圖8 不同段塞長度下節點35處位移時歷曲線Fig.8 Time histories of displacements at node 35 under different slug length

由圖8(a)可知，段塞長度的改變對觸地點X方向的位移沒有明顯影響。對于圖8(b)，可以發現隨著段塞長度的逐漸增加，觸地點的Z向位移幅值顯著增大。此外，在液塞段經過時，觸地點的Z向位移受到抑制，而氣塞段經過時則與之相反，這一現象是由液塞與氣塞交替流動導致的立管內部流體質量不斷變化而引起的。由圖9(a)、(b)可以看出，段塞長度的改變會導致立管懸掛段位移幅值的變化，隨著段塞長度的逐漸增加，懸掛段X與Z方向位移的雙頻疊加特性愈為明顯。

圖9 不同段塞長度下節點75處位移時歷曲線Fig.9 Time histories of displacements at Node 75 under different slug length

不同段塞長度下的彎矩時歷曲線如圖10所示。由圖10(a)可以看出，隨著段塞長度的增加觸地點彎矩明顯增大，其變化的雙頻疊加特性更加明顯，當Ls/D=100,280時，彎矩平均值分別為101.06 kN與127.53 kN,兩者之間差值為觸地點靜態彎矩的36%。對于圖10(b)，隨著段塞長度的增加，懸掛段彎矩的波動幅值逐漸增大，其變化的主導頻率與段塞流頻率相接近。

圖10 不同段塞長度下的彎矩時歷曲線Fig.10 Time histories of bending moment under different slug lengths

4 結論

1) 穩態段塞流的不均勻密度分布對立管的張力影響很小，但是對立管懸掛段的彎矩影響顯著，彎矩變化主導頻率與段塞流頻率相近。

2) 穩態段塞流會導致鋼懸鏈立管位移幅值的增大，液塞與氣塞交替流動導致的內流質量變化會改變觸地點垂向的運動特性。

3) 隨著段塞長度的增加，立管懸掛段位移的雙頻疊加特性更為明顯，觸地點處的彎矩顯著增大。