謝 鵬，常江濤，肖卓臨

(1. 中山大學(xué) 海洋工程與技術(shù)學(xué)院，廣東 珠海 519082; 2. 南方海洋科學(xué)與工程廣東省實(shí)驗(yàn)室(珠海)，廣東 珠海 519082)

天然氣水合物廣泛存在于海底的沉積物中，分布范圍廣、資源量巨大、能量密度高，是一種高效的清潔能源?，F(xiàn)行的水合物開采方法主要有加熱法、降壓法、注入外加劑法等[1]，其基本原理是通過改變水合物所處環(huán)境的溫度和壓力，打破水合物穩(wěn)定存在的平衡條件，使其分解為天然氣。國(guó)際上已經(jīng)有少數(shù)國(guó)家開展了水合物試采，其中產(chǎn)量最高的是我國(guó)于2017年5月10日實(shí)施的在我國(guó)南海神狐海域水合物試采，累計(jì)試采穩(wěn)定產(chǎn)氣60 d，總計(jì)產(chǎn)氣量為30.9×104m3，其中最高日產(chǎn)為3.5×104m3/d，完成了驗(yàn)證我國(guó)水合物開采關(guān)鍵技術(shù)的目標(biāo)[2]。根據(jù)國(guó)家規(guī)劃，預(yù)計(jì)2030年實(shí)現(xiàn)天然氣水合物的產(chǎn)業(yè)化開采，建成產(chǎn)能達(dá)到10×108m3的大型氣田。根據(jù)最高產(chǎn)量3.5×104m3/d的產(chǎn)量進(jìn)行保守估計(jì)，建造一個(gè)10×108m3的大型氣田，日均產(chǎn)量約為273×104m3。為實(shí)現(xiàn)這一目標(biāo)，預(yù)計(jì)需要在氣田內(nèi)布設(shè)78口同等流量的產(chǎn)氣井進(jìn)行同時(shí)開采和聯(lián)合集輸。因此，在大規(guī)模開發(fā)的天然氣水合物開采區(qū)，必須要使用大量海底管道和管匯將各井口串聯(lián)，并進(jìn)行天然氣集輸。

天然氣水合物的開采會(huì)誘發(fā)海床發(fā)生不均勻沉降，沉降會(huì)改變海底管道的受力狀態(tài)，甚至導(dǎo)致管道發(fā)生大變形，造成管道發(fā)生屈曲、斷裂、懸跨等失效[3-4]。大量研究表明，當(dāng)天然氣水合物在儲(chǔ)層內(nèi)分解后，水合物沉積層的結(jié)構(gòu)會(huì)發(fā)生很大程度的擾動(dòng)，其強(qiáng)度和力學(xué)特性也會(huì)發(fā)生改變，當(dāng)水合物的開采量達(dá)到某一程度時(shí)，則會(huì)誘發(fā)海床發(fā)生不均勻沉降、海底滑坡、土層失穩(wěn)等地質(zhì)災(zāi)害[5]，歷史上著名的挪威Storegga第二次海底滑坡，美國(guó)東海岸Carolina岸外的Cape Fear滑坡，西地中海的Balearic巨型濁流層和西非大陸架的海底滑坡體等事故都是由水合物分解造成的[6-7]。萬義釗等[8]基于“熱-流-固-化”四場(chǎng)耦合的水合物開采分析模型，研究了水合物開采早期的海床沉降規(guī)律和影響參數(shù)。研究發(fā)現(xiàn)，在開采第30天時(shí)，開采井位置處的最大沉降量為0.032 m，海底面沉降范圍的半徑約為166 m。開采第60天時(shí)開采井位置的最大沉降量為0.035 m，海底面沉降范圍的半徑約為232 m。金光榮等[9]通過建立南海神狐海域水合物降壓開采模型，研究了更長(zhǎng)時(shí)間開采范圍內(nèi)的海床沉降規(guī)律，獲得了海床沉降量與開采時(shí)間的關(guān)系曲線，并針對(duì)開采井500 m范圍內(nèi)的海床沉降情況進(jìn)行了詳細(xì)分析。結(jié)果顯示，開采第2年時(shí)，在開采井附近位置海床沉降量達(dá)到2.6 m，且隨著開采時(shí)間和開采量的增加，海床沉降量增大。

然而，現(xiàn)有的研究主要針對(duì)天然氣水合物開采過程中水合物的分解過程進(jìn)行研究，并預(yù)測(cè)水合物開采過程引起的海床沉降規(guī)律，并未考慮布設(shè)在開采區(qū)范圍內(nèi)的海底管道等工程設(shè)施，也未研究水合物開采引起的海底管道受力狀態(tài)和失效模式。為解決這一問題，基于有限元分析方法，建立天然氣水合物開采區(qū)內(nèi)“海床-管道”耦合作用模型，以“生死單元”模擬了天然氣水合物開采過程中海床沉降過程，獲得了海床變形形態(tài)及其對(duì)管道應(yīng)力、應(yīng)變、彎矩、懸跨等力學(xué)行為的影響。研究結(jié)果可以為建立天然氣水合物開采區(qū)內(nèi)海底管道的設(shè)計(jì)方法提供依據(jù)，為未來天然氣水合物的大規(guī)模產(chǎn)業(yè)化開采做技術(shù)儲(chǔ)備。

1 “管-土”相互作用的有限元分析

水合物開采過程中，位于井口附近的水合物先分解，這將會(huì)擾動(dòng)開采井口附近區(qū)域的沉積層結(jié)構(gòu)，改變沉積層的強(qiáng)度和力學(xué)特性，誘發(fā)海床發(fā)生不均勻沉降，而沉降最嚴(yán)重的區(qū)域就是開采井口，海床的沉降量以開采井口為中心向外輻射遞減。研究位于水合物開采區(qū)內(nèi)場(chǎng)(in-field)的海底管道，該管道需要與開采井口進(jìn)行連接，在大規(guī)模的水合物開發(fā)過程中預(yù)計(jì)需要在氣田內(nèi)布設(shè)78口同等流量的開采井進(jìn)行同時(shí)開采和聯(lián)合集輸，必須要使用大量的海底管道和管匯將各井口串聯(lián)起來，因此必然導(dǎo)致海底管道經(jīng)過井口附近。當(dāng)管道從開采井口附近穿過時(shí)，由于海床的沉降，管道在自身重力下也會(huì)發(fā)生豎向位移，位移最大處位于開采井處，由于管道發(fā)生變形，管道的應(yīng)力、應(yīng)變、彎矩等力學(xué)特性也會(huì)發(fā)生相應(yīng)的變化。因此，有限元模型的建立以海底管道跨越單井口為研究對(duì)象，模擬在進(jìn)行單井口水合物開采時(shí)引起的海床沉降規(guī)律，及其對(duì)附近管道力學(xué)性能的影響。

1.1 海底管道模型

基于ABAQUS有限元軟件進(jìn)行數(shù)值模擬，管道模型全長(zhǎng)600 m，管道半徑0.31 m, 厚度0.02 m，如圖1所示。管道模型采用X65鋼，主要材料參數(shù)如表1所示。ABAQUS中適用于接觸問題的單元須為減縮積分單元[10]，因此管道采用減縮積分的4節(jié)點(diǎn)曲面殼單元(S4R)，該單元是一種通用的殼單元類型，具有很好的適應(yīng)性，既可以用于厚殼問題的模擬，也可以用于薄殼問題的模擬。減縮積分是相對(duì)于完全積分來說的，減縮積分單元比完全積分單元在每個(gè)方向上少使用一個(gè)積分點(diǎn)，在模型受彎曲時(shí)可以消除計(jì)算中的剪切鎖閉問題且對(duì)計(jì)算結(jié)果的精度影響不大。

圖1 管道模型Fig. 1 Pipeline model

表1 X65 鋼材主要參數(shù)Tab. 1 Main parameters of X65

采用Romberg-Osgood方程建立管道材料的本構(gòu)關(guān)系，其表達(dá)式為：

(1)

式中：ε是材料應(yīng)變；σ是材料應(yīng)力；E為彈性模量；σs是屈服強(qiáng)度；A和B是描述材料塑性應(yīng)變的硬化參數(shù)，A=1.29，B=25.58 。

圖2 海床模型Fig. 2 Seabed model

1.2 海床模型

將海床的不均勻沉降作為海底管道力學(xué)分析的邊界條件，基于有限單元法的實(shí)體單元模擬水合物儲(chǔ)層被開挖后海床的沉降變形過程，構(gòu)造出“等效的”不均勻海床表面變形形態(tài)，不需要對(duì)水合物開采過程中的熱量傳遞、固氣轉(zhuǎn)變等過程進(jìn)行分析。分析時(shí)，海床選擇減縮積分的8節(jié)點(diǎn)六面體單元(C3D8R)，該單元對(duì)模擬非常大的網(wǎng)格扭曲問題有很好的適應(yīng)性，由于模型中海床會(huì)發(fā)生較大的變形，因此選用該單元進(jìn)行網(wǎng)格劃分。海床模型全長(zhǎng)600 m，寬300 m，高200 m，其中天然氣水合物層深度100 m，厚度50 m，海床上表面設(shè)置寬度為20倍管道直徑的條帶作為和管道接觸的區(qū)域，如圖2所示。海床模型材料簡(jiǎn)化為彈塑性材料，模型主要參數(shù)如表2所示。該模型的分析參數(shù)基于我國(guó)南海北部陸坡神狐海域天然氣水合物的成藏規(guī)律確定，研究結(jié)果可以代表未來在該海域內(nèi)海底管道沉降及受力分析的一般情況。

表2 海床主要參數(shù)Tab. 2 Main parameters of seabed

1.3 接觸設(shè)置

管道外壁與海床上表面之間設(shè)置接觸面，采用單純的“主控-從屬”接觸算法。在計(jì)算過程中，為了獲得較好的模擬結(jié)果，選擇管道的外表面作為主接觸面，海床的上表面取20倍管徑寬度的條帶區(qū)域作為從接觸面。接觸作用屬性設(shè)置為法向“硬”接觸，切向采用罰函數(shù)摩擦公式定義切向摩擦力。

2 基于“生死單元”模擬海床不均勻沉降過程

“生死單元”功能是為了克服有限元本身缺陷提出的一種分析方法，被廣泛應(yīng)用于礦產(chǎn)開采、隧道鉆進(jìn)、基坑開挖等巖土工程領(lǐng)域，用來模擬被開挖部分“消失”的情況[11-12]。該單元彌補(bǔ)了有限單元法常規(guī)單元的不足。由于有限元本身是基于連續(xù)介質(zhì)力學(xué)的，所研究的物體需要是連續(xù)的，在這樣的理論框架下，單元本身是不會(huì)消失的。然而在水合物開采過程中，由于水合物分解并經(jīng)開采井運(yùn)輸后，勢(shì)必在原區(qū)域產(chǎn)生短暫的環(huán)空，并造成海床沉降。因此，在進(jìn)行數(shù)值分析時(shí)，需要模擬水合物儲(chǔ)層單元消失或者完全失效的情況。為了模擬這種工況，采用ABAQUS“生死單元”功能進(jìn)行研究。在移除單元之前，ABAQUS自動(dòng)計(jì)算和儲(chǔ)存待移除單元和余下單元交界面上的節(jié)點(diǎn)力，這些節(jié)點(diǎn)力在移除步中逐漸減少為零。在有限元分析過程中，通過分步移除水合物層一定范圍內(nèi)的單元，模擬水合物逐步開采的過程，在此基礎(chǔ)上得到海床的不均勻沉降模型，如圖3所示。具體分析過程可分為13個(gè)分析步。

圖3 海床的不均勻沉降模擬Fig. 3 Simulation of uneven subsidence of seabed

第1分析步：管道左端施加x、y、z三個(gè)方向的約束；右端施加x方向的約束，約束管道的橫向位移。海床底面施加固定約束，海床左右兩端施加x、y、z三個(gè)方向的約束。

第2分析步：為管道施加重力載荷，重力加速度取一較小值0.1 m/s2,使得管道與海床之間的接觸作用平穩(wěn)的建立起來。同時(shí)在管道的右端施加一水平張力T0來模擬管道內(nèi)部的軸向拉力，其值為200 kN，此值主要根據(jù)工程經(jīng)驗(yàn)得到。

第3分析步：為海床施加體力載荷，海床的重度為20 kN/m3,將管道受到的重力加速度改為10 m/s2。

第4分析步：從第4分析步(即第1步開采)開始模擬水合物的開采過程，開采井口設(shè)置在海床模型上表面的中心點(diǎn)處，開采位置設(shè)置在開采井口下100 m深的水合物層處，開采厚度50 m，開采范圍20 m。開采過程通過ABAQUS中的單元生死功能實(shí)現(xiàn)，通過移除水合物層一定范圍內(nèi)的單元，來模擬逐步開采的過程。

第5～13分析步：重復(fù)第4分析步，開采范圍以20 m為梯度逐步增大，共開采10步，最終開采范圍為200 m。通過這種逐步開采的方式獲得海床的沉降以及管道的力學(xué)響應(yīng)。

3 天然氣水合物開采區(qū)內(nèi)海床沉降規(guī)律及驗(yàn)證

3.1 天然氣水合物開采過程中海床不均勻沉降量

通過ABAQUS進(jìn)行數(shù)值模擬，對(duì)水合物進(jìn)行分步開采，得到開采過程中的每一步海床沉降量及形態(tài)分布結(jié)果，如圖4、圖5所示。提取第5、第7、第9和第10步開采的位移云圖，開采影響范圍分別為100 m、140 m、180 m和200 m，分別如圖5所示。從圖中可以明顯的發(fā)現(xiàn)，隨著開采范圍的增大，海床的位移逐漸增大，第10步開采時(shí)，海床沉降量達(dá)到最大值3.543 m，位移云圖顯示的較大位移區(qū)域逐漸向開采井處收縮。

圖4 海床位移云圖Fig. 4 Displacement nephogram of seabed

圖5 不同開采步下模型的位移云圖Fig. 5 Displacement nephogram of models under different mining steps

3.2 海床沉降量規(guī)律的驗(yàn)證

以“生死單元”法模擬實(shí)體單元消失引起海床沉降變形，等效模擬真實(shí)開采過程中的海床沉降量和海床表面形態(tài)分布規(guī)律，結(jié)果如圖6所示。x=0處為開采井，圖6(a)為海床沉降量在模型長(zhǎng)度方向上隨開采井距離的變化曲線，圖6(b)為海床沉降量在模型寬度方向上隨開采井距離的變化曲線；曲線之間的變化趨勢(shì)相同，海床沉降量隨著與開采井距離的增大而減小。開采過程中開采井處的海床沉降量如表3所示。

圖6 海床沉降量在模型長(zhǎng)度和寬度方向上的變化Fig. 6 Variation of seabed subsidence with the length and width of the model

表3 開采過程中開采井處的海床沉降量Tab. 3 Subsidence of seabed at the mining well in mining process

為驗(yàn)證該模型模擬結(jié)果的準(zhǔn)確性和適用性，與基于“熱-水動(dòng)力-機(jī)械耦合”模型模擬的水合物開采引起的海床沉降情況進(jìn)行對(duì)比和驗(yàn)證。金光榮等[9]通過使用TOUGH2Biot模擬器，擴(kuò)展了Biot合并模型，從而得到了用于水合物開采的“熱-水動(dòng)力-機(jī)械耦合”模型，通過改變開采井底部壓力，獲得了神狐海域天然氣水合物開采過程中海床不均勻沉降的規(guī)律。

將“生死單元”法獲得的海床沉降量與文獻(xiàn)進(jìn)行對(duì)比，結(jié)果如圖7所示，圖7(a)為“生死單元”法海床沉降計(jì)算結(jié)果，圖7(b)為金光榮等[9]通過降壓法開采水合物得到的海床沉降曲線。圖7(a)分別為第5、第7和第9步開采時(shí)海床沉降的結(jié)果，海床最大沉降量分別為1.621 m、2.068 m和2.767 m，沉降發(fā)生的主要區(qū)域位于開采井400 m范圍內(nèi)；圖7(b)分別為開采時(shí)間為第0.5、第1和第2年時(shí)海床沉降的結(jié)果，海床最大沉降量分別為1.52 m、1.80 m和2.50 m，沉降發(fā)生的主要區(qū)域位于開采井500 m范圍內(nèi)。兩種結(jié)果較為接近，海床的最大沉降都位于開采井處，當(dāng)海床與開采井的距離增大時(shí)，沉降量都會(huì)逐漸減小。將圖7(a)中的曲線擬合到圖7(b)中，可以看出文中“生死單元”法模擬的水合物開采引起的海床表面沉降形態(tài)與“熱-水動(dòng)力-機(jī)械耦合”法模擬的海床表面變形形態(tài)能夠很好的吻合，海床沉降的數(shù)值也很接近前人的研究成果。驗(yàn)證了通過使用ABAQUS中單元生死功能模擬水合物開采進(jìn)而引發(fā)海床沉降這一方法是可行的，以該模型獲得的海床沉降規(guī)律作為管道受力分析的邊界條件是合理的。

圖7 海床位移曲線Fig. 7 Displacement curve of seabed

4 海床不均勻沉降引起的管道力學(xué)性能變化

基于獲得的海床沉降規(guī)律，研究海床沉降對(duì)管道造成的力學(xué)行為變化，分析管道的豎向位移、曲率和彎矩、應(yīng)力應(yīng)變分布、管道懸跨等問題。

4.1 管道豎向位移

在水合物開采過程中，管道豎向位移隨著海床的沉降發(fā)生變化。圖8(a)為管道位置變化示意，圖8(b)為每一開采步的管道豎向位移量，表4為開采過程中開采井處管道的沉降量。結(jié)果顯示，管道在開采井附近的豎向位移最大，為3.407 m，當(dāng)管道位置遠(yuǎn)離開采井時(shí)，由于海床的沉降量逐漸減小，管道的豎向位移也逐漸減小。管道的位移曲線關(guān)于開采井位置對(duì)稱且沉降主要發(fā)生在開采井周圍400 m范圍內(nèi)，沉降的最大值總是出現(xiàn)在開采井處。這一形態(tài)與海床的表面沉降形態(tài)是一致的。

圖8 管道的豎向位移Fig. 8 Vertical displacement of pipeline

表4 開采過程中開采井處的管道沉降量Tab. 4 Subsidence of pipeline at the mining well in mining process

4.2 管道的曲率及彎矩分布

圖9 管道曲率Fig. 9 Curvature of pipeline

圖10 管道彎矩Fig. 10 Bending moment of pipeline

4.3 管道的應(yīng)力應(yīng)變分布

由于天然氣水合物的開采導(dǎo)致管道發(fā)生較大的沉降，盡管沉降的跨度較大，但是依然會(huì)對(duì)管道的應(yīng)力應(yīng)變?cè)斐捎绊憽Ｓ捎诤４渤两抵饕l(fā)生在開采井附近區(qū)域，距離開采井附近150 m區(qū)域內(nèi)的沉降量較大，對(duì)管道應(yīng)力應(yīng)變的影響更為顯著，因此只需選擇該區(qū)域進(jìn)行分析即可。如圖11(a)所示，管道應(yīng)力和應(yīng)變的最大值均出現(xiàn)在管道中點(diǎn)的底部。每一步開采后，管道的最大應(yīng)力和應(yīng)變值如圖11(b)所示，隨著水合物開采范圍的擴(kuò)大，管道的最大應(yīng)力和應(yīng)變都隨之增加，應(yīng)力值的范圍為9.97～24.37 MPa，應(yīng)變值的范圍為4.38×10-5～1.17×10-4。

圖11 管道的最大應(yīng)力、應(yīng)變Fig. 11 Maximum stress and strain of pipeline

為了探究該150 m范圍內(nèi)管道上應(yīng)力應(yīng)變隨水平位置的變化趨勢(shì)，須對(duì)每一步開采所對(duì)應(yīng)的管道的應(yīng)力應(yīng)變進(jìn)行數(shù)據(jù)提取，得到的結(jié)果如圖12所示。

圖12 開采過程中管道的應(yīng)力、應(yīng)變分布Fig. 12 Stress and strain distribution of pipeline in mining process

如圖12所示，每一步開采后，管道的應(yīng)力應(yīng)變隨水平位置的變化趨勢(shì)是一致的，當(dāng)管道位置靠近開采井時(shí)，管道的位移和變形逐漸增大，導(dǎo)致管道的應(yīng)力應(yīng)變逐漸增大；當(dāng)管道位置遠(yuǎn)離開采井時(shí)，管道的位移和變形逐漸減小，導(dǎo)致管道的應(yīng)力應(yīng)變逐漸減小。管道的應(yīng)力應(yīng)變曲線關(guān)于開采井位置對(duì)稱且最大值總是出現(xiàn)在開采井處。從圖11和圖12可以看出，天然氣水合物開采導(dǎo)致海床的不均勻沉降會(huì)改變管道的應(yīng)力應(yīng)變狀態(tài)，開采范圍的擴(kuò)大會(huì)導(dǎo)致管道內(nèi)部的應(yīng)力和應(yīng)變逐漸增大，但是其數(shù)值都還在安全范圍以內(nèi)，因此水合物的開采對(duì)管道的應(yīng)力應(yīng)變影響相對(duì)較小，在應(yīng)力應(yīng)變的角度上，管道是安全可靠的。

4.4 管道的懸跨

在水合物開采過程中，管道雖然隨著海床的沉降發(fā)生豎向位移，但由于管道本身具有彎曲剛度且受到軸向拉力的約束，在豎向位移發(fā)生到一定程度后，則不再發(fā)生變化。因此，當(dāng)海床的沉降量達(dá)到一定值時(shí)，管道不再繼續(xù)沉降，而是與海床發(fā)生分離，在海床表面形成懸跨，且逐步增大，如圖13所示。在開采過程中海床和管道的相對(duì)位移變化量如圖14所示。如圖14(a)所示，天然氣水合物開采在前4步即開采范圍在80 m以內(nèi)時(shí)，海床與管道的最大沉降量相等，這表示海床與管道是始終接觸的，當(dāng)天然氣水合物開采到第5步即開采范圍達(dá)到100 m時(shí)，海床與管道的最大沉降量分別為1.621 m和1.576 m，出現(xiàn)0.045 m的懸空，此時(shí)海床和管道已發(fā)生分離。隨著開采范圍的繼續(xù)增大，海床與管道之間的懸空距離整體呈上升趨勢(shì)，如圖14(b)所示。當(dāng)海水流經(jīng)懸跨管道時(shí)，在一定流速條件下，懸跨管道兩側(cè)會(huì)出現(xiàn)旋渦，并以一定的頻率交替脫落，從而在管道表面形成周期載荷，將會(huì)使懸跨管道在順流方向及橫流方向上發(fā)生振動(dòng)，即渦激振動(dòng)，很可能對(duì)海底管道造成疲勞損傷乃至破壞[13-15]。

圖13 管道與海床分離產(chǎn)生懸跨(開采步10)Fig. 13 Suspension span caused by separation of pipeline and seabed (Mining step 10)

圖14 開采過程中管道和海床的沉降量及懸跨高度變化Fig. 14 Change of seabed subsidence and height of suspension span of pipeline during mining process

5 結(jié) 語

經(jīng)過研究和計(jì)算分析，可以得到以下結(jié)論：

1) 有限單元中“生死單元”法可以等效模擬天然氣水合物開采過程中海床的變形規(guī)律，為海底管道的沉降分析提供邊界條件，對(duì)于模擬海底管道在海床沉降邊界下的安全校核有重要的指導(dǎo)意義。

2) 天然氣水合物開采引起的海床沉降區(qū)域主要集中在開采井附近400 m的范圍內(nèi)，其對(duì)管道的應(yīng)力、應(yīng)變、位移、曲率和彎矩均產(chǎn)生影響。由海床沉降引起的管道應(yīng)力應(yīng)變基本處于彈性范圍內(nèi)，不會(huì)造成管道本身的強(qiáng)度破壞；但是當(dāng)海床沉降達(dá)到一定程度時(shí)，鋪設(shè)在開采區(qū)內(nèi)的海底管道將與海床發(fā)生分離，形成懸跨，并引發(fā)渦激振動(dòng)風(fēng)險(xiǎn)，對(duì)管道造成疲勞破壞。

3) 基于天然氣水合物“單井開采”模型進(jìn)行海床沉降研究，與實(shí)際工程中的“多井聯(lián)合開采”模式仍有一定差距。后續(xù)研究中可進(jìn)一步考慮多井聯(lián)合開采工況下海床的變形規(guī)律及其影響，為天然氣水合物的大規(guī)模開發(fā)提供技術(shù)支撐。