蘇義腦 陳 燁 孫曉峰 閆 鐵 曲晶瑀 段瑞溪

1.東北石油大學(xué)石油工程學(xué)院 2.中國石油海洋工程有限公司

0 引言

我國海域天然氣水合物（以下簡稱水合物）資源豐富，其儲(chǔ)層特點(diǎn)與日本類似，儲(chǔ)層埋深淺，膠結(jié)較弱，未固結(jié)成巖，水合物在砂土中以填充、骨架和膠結(jié)幾種土體賦存方式存在[1-5]，一旦試采破壞了水合物的賦存條件[6-7]，起膠結(jié)和骨架支撐作用的水合物土體就會(huì)失去支撐，儲(chǔ)層砂粒在地層流體的攜帶下流向井筒。水合物儲(chǔ)層巖性以微米級(jí)黏土、粉砂、細(xì)砂為主，南海多個(gè)站位巖心樣品粒度分布測試結(jié)果顯示粒徑小于63 μm的粉砂超過65%[8-10]，有的站位超過70%[11]，甚至有測試結(jié)果顯示小于63 μm的黏土和粉砂占比接近94%[12]。對(duì)于微米級(jí)的泥質(zhì)粉砂，防砂技術(shù)難度較大，目前水合物試采所用的防砂技術(shù)主要包括機(jī)械防砂、加裝防砂網(wǎng)、射孔砂篩防砂、Geoform防砂系統(tǒng)等，這些防砂技術(shù)效果并不理想，實(shí)踐和研究表明，防砂系統(tǒng)很難對(duì)粒徑44 μm以下的砂粒發(fā)揮作用[13]。近年來，陸地（加拿大，2007年）和海域（日本，2013年、2017年）水合物開采，都發(fā)生了因微米級(jí)砂粒突破防砂設(shè)施進(jìn)入井筒，引起井筒砂埋、堵塞，迫使水合物試采作業(yè)提前終止的案例[14-16]。日本共進(jìn)行了兩次大規(guī)模海域水合物開采試驗(yàn)，其中2013年為世界首次海域水合物開采試驗(yàn)，共產(chǎn)氣6天，日產(chǎn)氣2×104m3，后由于大量儲(chǔ)層砂粒突破防砂設(shè)施進(jìn)入井筒，導(dǎo)致試采提前結(jié)束[15,17-19]。2017年，日本總結(jié)第一次水合物試采經(jīng)驗(yàn)，進(jìn)行第二次海域水合物試采作業(yè)，同樣由于井筒沉砂堵塞而終止產(chǎn)氣[20]。因而，水合物開采過程中不僅需要考慮儲(chǔ)層微米級(jí)砂粒突破井筒防砂設(shè)施后井筒設(shè)備的磨損問題，更要考慮微米級(jí)砂粒在井筒內(nèi)的運(yùn)移、沉積、堵塞問題，需要對(duì)微米級(jí)砂粒進(jìn)入井筒后的運(yùn)移沉積規(guī)律特征進(jìn)行研究。

目前，水合物出砂研究主要集中在地質(zhì)儲(chǔ)層中，包括儲(chǔ)層巖石力學(xué)特性（產(chǎn)砂機(jī)理）[21-27]、驅(qū)動(dòng)力[28]、出砂風(fēng)險(xiǎn)[29-33]、評(píng)價(jià)防砂方法有效性[34-37]幾個(gè)方面，而微米級(jí)砂粒隨地層流體進(jìn)入井筒后的運(yùn)移、沉積、堵塞研究尚未見文獻(xiàn)報(bào)道。本文研究的目的是：①通過數(shù)值模擬，獲得微米級(jí)砂粒在節(jié)流管線螺旋段的運(yùn)移沉積特征；②獲得粒徑、水速等參數(shù)對(duì)微米級(jí)砂粒運(yùn)移沉降的影響規(guī)律；③獲得不同情況下的砂粒臨界不沉積水速，建立局部復(fù)雜井段微米級(jí)砂粒運(yùn)移沉積預(yù)測模型，為水合物開采井筒砂埋提供預(yù)測方法；④為減小砂堵風(fēng)險(xiǎn)、確定合理開采工作制度、保證井筒安全提供依據(jù)，為后續(xù)水合物儲(chǔ)層適度防砂工藝的發(fā)展提供技術(shù)支撐。

1 計(jì)算方法

1.1 計(jì)算模型

降壓開采法是目前水合物開采的主流方法之一。前蘇聯(lián)麥索亞哈、加拿大麥肯齊、美國阿拉斯加北坡、日本南海海槽東部、日本南海海槽、中國南海神狐都采用過降壓法開采水合物[20]，日本2013年在南海海槽東部采用降壓開采法進(jìn)行了第一次海域水合物試采作業(yè)，其開采系統(tǒng)如圖1-a所示。

降壓破壞了水合物穩(wěn)定賦存的壓力條件，使水合物發(fā)生分解，產(chǎn)生甲烷氣和水，分解后產(chǎn)生的氣液兩相在井下分離并通過獨(dú)立的管路產(chǎn)出，其中的甲烷等烴類氣體通過油管通道舉升到地面，分離水通過節(jié)流管線單獨(dú)輸運(yùn)到井口。分離水在井下分離后經(jīng)過C型通道進(jìn)入節(jié)流管線，由于節(jié)流管線下部為剛性管路，且上下接箍位置不在同一側(cè)，因而存在一段螺旋段，如圖1-b所示。可見，水流通道既存在C型段，也存在螺旋段，是流道復(fù)雜砂粒易沉降的位置。為此，對(duì)水合物開采局部復(fù)雜管路微米級(jí)砂粒運(yùn)移規(guī)律進(jìn)行數(shù)值模擬研究，建立計(jì)算的幾何模型如圖1-c所示。由于節(jié)流管線為水流通道，為此選定水為連續(xù)相，微米級(jí)砂粒為分散相。

1.2 控制方程

連續(xù)性方程為：

式中下標(biāo)l和s分別表示水和微米級(jí)砂粒；α表示濃度；表示速度向量；ρ表示密度。

每個(gè)計(jì)算單元相互滲透，因此他們濃度的和是一致的。水和砂粒的動(dòng)量方程可以表示為：

圖1 水合物開采系統(tǒng)、節(jié)流管線及計(jì)算模型示意圖

當(dāng)αl≤0.8時(shí)，水—砂粒的交換系數(shù)可以表示為：

當(dāng)αl＞0.8時(shí)，水—砂粒的交換系數(shù)可以表示為：

拖拽系數(shù)（CD）可表示為：

砂粒的雷諾數(shù)定義為：

式中ds表示砂粒粒徑。

各相的應(yīng)力—應(yīng)變張量可表示為：

Realizablek—ε模型中的湍動(dòng)能（k），耗散率（ε）為：

式中σk和σ?分別表示湍動(dòng)能方程k和耗散率方程ε的紊流普朗特?cái)?shù)；C2表示常數(shù)。

Gk表示由于層流速度梯度而產(chǎn)生的湍流動(dòng)能，可以表示如下：

1.3 邊界條件及網(wǎng)格劃分

水和均勻分散的微米級(jí)砂粒從節(jié)流管線底部進(jìn)入，經(jīng)C型段、螺旋段后進(jìn)入垂直段而后到達(dá)地面，如圖2所示。在入口處選擇速度入口邊界條件，出口選擇壓力出口邊界條件，壁面采用無滑移壁面條件。

適當(dāng)?shù)木W(wǎng)格劃分軟件和高精度網(wǎng)格質(zhì)量是CFD計(jì)算的必要先決條件。為提高計(jì)算準(zhǔn)確性，加快計(jì)算和收斂速度，使用ICEM軟件采用“自下而上”方式對(duì)計(jì)算流場劃分六面體結(jié)構(gòu)網(wǎng)格，圖2為整體計(jì)算流場及局部位置三維網(wǎng)格劃分示意圖。設(shè)置邊界層網(wǎng)格，管路徑向第一層網(wǎng)格高度的選取符合y＋≈30，網(wǎng)格增長率取1.2，沿徑向第一層網(wǎng)格高度采用式（16）計(jì)算

式中Δy表示第一層網(wǎng)格高度，mm；Re表示雷諾數(shù)；D表示水力特征長度，mm。

為了檢測網(wǎng)格獨(dú)立性，劃分3種數(shù)量網(wǎng)格（158 347，200 816，411 075）用于流場（44 μm，1%濃度）模擬，模擬結(jié)果如表1所示?？梢钥闯銎湎鄬?duì)誤差均小于4.5%。綜合考慮計(jì)算資源消耗、計(jì)算準(zhǔn)確性，選擇網(wǎng)格數(shù)200816進(jìn)行數(shù)值模擬計(jì)算。

通過有限體積法對(duì)控制方程進(jìn)行離散化之后，將相位耦合的SIMPLE算法應(yīng)用于壓力—速度耦合。動(dòng)量方程使用二階迎風(fēng)格式進(jìn)行離散化以提高計(jì)算精度。對(duì)每個(gè)尺度殘差分量采用10－5的收斂準(zhǔn)則來確定兩個(gè)連續(xù)迭代之間的相對(duì)誤差。時(shí)間步長取0.005 s，每個(gè)步長迭代20步。

選取3種粒徑砂粒（44 μm，10 μm，4 μm），3種出砂濃度（1%，5%，10%），多種流速進(jìn)行數(shù)值模擬計(jì)算，獲得不同流速下微米級(jí)砂粒在復(fù)雜井段的運(yùn)移沉積情況，以進(jìn)一步確定不同粒徑砂粒在不同濃度下運(yùn)移的臨界不沉積水速。

圖2 整體計(jì)算區(qū)域及局部位置3D網(wǎng)格劃分圖

表1 網(wǎng)格獨(dú)立性檢驗(yàn)表

2 模擬結(jié)果及分析

2.1 微米級(jí)砂粒沉積分析

圖3顯示了44 μm粒徑、10%出砂濃度在不同水速下微米級(jí)砂粒濃度分布云圖。

由圖3可以看出水速較小時(shí)，微米級(jí)砂粒在管路中的沉積較為明顯，主要集中在C型段和螺旋段的下部（圖3-a）；局部砂粒濃度最高達(dá)0.55。隨著流速的增加，砂粒沉積情況改善較大，C型段、螺旋段下部沉積減少，局部砂粒濃度最大為0.24，集中在螺旋段中上部，管路其余部分多為0.12左右（圖3-b），繼續(xù)增大水速，微米級(jí)砂粒沉積越來越小，只在C型段轉(zhuǎn)角留存有部分沉積（圖3-c）?？梢?，微米級(jí)砂粒對(duì)水速較為敏感，較小的水速變化就能引起較大的砂粒濃度的變化。微米級(jí)砂粒較易沉積的部位主要在C型段的拐角處及螺旋段，這是由于流向變化和重力作用引起的。

圖3 44 μm粒徑、10%出砂濃度在不同水速下微米級(jí)砂粒濃度分布云圖

2.2 微米級(jí)砂粒臨界不沉積流速的確定

對(duì)3種粒徑、3種出砂濃度下微米級(jí)砂粒在復(fù)雜井段內(nèi)的運(yùn)移沉積情況進(jìn)行數(shù)值模擬，共獲得100多組模擬數(shù)據(jù)，模擬結(jié)果如圖4所示。

圖4 3種粒徑、3種濃度下流速與復(fù)雜管段內(nèi)砂粒沉積濃度關(guān)系圖

以粒徑44 μm顆粒1%濃度下，管路內(nèi)微米級(jí)砂粒濃度隨水速變化規(guī)律為例，由圖4-a可以看出，當(dāng)水速較大時(shí)，復(fù)雜管段內(nèi)砂粒濃度無限趨近于地層出砂濃度1%，說明管路內(nèi)幾乎沒有砂粒沉降，曲線變化平緩。隨著水速的逐漸減小，管路內(nèi)微米級(jí)砂粒濃度逐漸增加，當(dāng)水速減小到一定程度后，小幅度的速度減小將引起管路內(nèi)砂粒濃度的急劇升高，濃度呈指數(shù)增長趨勢。因此可將曲線變化大體劃分為3個(gè)部分：平緩區(qū)、指數(shù)增長區(qū)以及兩者交匯的過渡區(qū)。臨界不沉積水速位于兩種變化趨勢交匯的過渡區(qū)內(nèi)。如何確定臨界不沉積水速，筆者提出了一種方法：由于臨界不沉積水速位于過渡區(qū)內(nèi)，可分別對(duì)平緩區(qū)和指數(shù)增長區(qū)進(jìn)行直線擬合，擬合曲線分別表示兩種區(qū)域的變化趨勢，則這兩條直線的交點(diǎn)即為臨界點(diǎn)，其橫坐標(biāo)即為臨界不沉積水速。水速大于這個(gè)值，管路內(nèi)砂粒沉積量較小，水速小于這個(gè)值，管路內(nèi)的砂粒沉積量將急劇增加，且水速越小，沉積量越大。由此可得到44 μm，1%濃度下微米級(jí)砂粒在水合物試采局部復(fù)雜井段的臨界不沉積水速為0.106 5 m/s。同理，可得到其他工況下的臨界不沉積水速，如圖4-b～i所示。

2.3 臨界不沉積流速結(jié)果分析

微米級(jí)砂粒臨界不沉積水速隨粒徑及地層出砂濃度的變化規(guī)律如圖5、6所示。

圖5 不同地層出砂濃度下臨界流速隨粒徑變化規(guī)律圖

由圖5可以看出，隨著粒徑的增大，微米級(jí)砂粒的臨界不沉積水速逐漸增大，粒徑較小時(shí)，出砂濃度的變化對(duì)臨界不沉積水速的影響較小，粒徑較大時(shí)，出砂濃度差別造成的臨界不沉積水速的變化變得更加明顯。

圖6 不同粒徑下臨界流速隨地層出砂濃度變化規(guī)律圖

由圖6可以看出，隨著地層出砂濃度的增加，微米級(jí)砂粒的臨界不沉積水速逐漸增加，且濃度較大時(shí)，粒徑變化造成的臨界不沉積水速的變化更大。造成這種現(xiàn)象的原因是因?yàn)殡S著粒徑的增大，自重增加，砂粒懸浮需要更大的流速，另一方面，砂粒較小時(shí)，比表面積更大，砂粒更易懸浮在溶液中，不易沉降。而且地層出砂濃度越大，井筒內(nèi)含砂量越多，沉積更容易發(fā)生。粒徑和地層出砂濃度的雙重作用導(dǎo)致粒徑較大時(shí)，濃度差別造成的臨界不沉積水速變化更加明顯，出砂濃度較大時(shí)，粒徑差別造成的臨界不沉積水速變化更大。

2.4 微米級(jí)砂粒沉積預(yù)測模型建立

由于試采局部井段管路復(fù)雜，室內(nèi)實(shí)驗(yàn)存在較大難度，數(shù)值模擬計(jì)算需要耗費(fèi)較長的時(shí)間和計(jì)算資源，為獲得一種相對(duì)便捷且成本較低的計(jì)算方式，對(duì)上述模擬數(shù)據(jù)進(jìn)行分析，水合物試采局部復(fù)雜管路微米級(jí)砂粒的濃度（C）是地層出砂濃度（C0）、液流速度（v）、液流密度（ρl）、液流黏度（μl）、砂粒粒徑（ds）及螺旋段井筒直徑（Dlx）的函數(shù)，對(duì)這些變量進(jìn)行分析，如表2所示，這些變量之間的數(shù)學(xué)關(guān)系可表示為：

表2 獨(dú)立變量分析表

應(yīng)用Buckingham-Π定理，對(duì)上述變量進(jìn)行無量綱分析，組合成無量綱群，將復(fù)雜問題簡化為較小組合變量之間的關(guān)系。為此定義4個(gè)具有實(shí)際意義的無量綱量，對(duì)所研究變量進(jìn)行無量綱化。第1個(gè)無量綱量為復(fù)雜管路井筒內(nèi)微米級(jí)砂粒濃度（π1，C），是因變量，其他變量為自變量；第2個(gè)無量綱量為地層出砂濃度（π2，C0），表示地層出砂速率對(duì)砂粒沉積運(yùn)移的影響；第3個(gè)無量綱量為（π3，ds/Dlx），代表了粒徑對(duì)微米級(jí)砂粒運(yùn)移沉積的影響；第4個(gè)無量綱量為（π4，ρ1vDlx/μ1）代表了螺旋管雷諾數(shù)，式中的特征長度為螺旋管的水力直徑，具體如表3所示。

表3 獨(dú)立變量無量綱化表

4組具有實(shí)際物理意義的無量綱量之間的數(shù)學(xué)關(guān)系可以描述如下：

利用OriginPro 2019軟件非線性擬合工具對(duì)其進(jìn)行分析，確定擬合方程如下：

式中a1、a2、a3、a4表示系數(shù)。對(duì)數(shù)據(jù)進(jìn)行擬合，擬合結(jié)果如表4所示。

則最終的復(fù)雜管路微米級(jí)砂粒沉積預(yù)測方程為：

對(duì)數(shù)據(jù)點(diǎn)模型預(yù)測值進(jìn)行統(tǒng)計(jì)誤差分析來檢驗(yàn)?zāi)Ｐ偷臏?zhǔn)確性。檢驗(yàn)?zāi)Ｐ蜏?zhǔn)確性的統(tǒng)計(jì)學(xué)參數(shù)包括：平均誤差（ME）、平均絕對(duì)誤差（MAE）、標(biāo)準(zhǔn)差（SD）、平均絕對(duì)百分比誤差（MAPE）、均方根誤差（RMSE），這些統(tǒng)計(jì)參數(shù)定義如式（21）～（25）所示[39]。統(tǒng)計(jì)參數(shù)分析結(jié)果如表5所示?？梢钥闯鯩E、MAE、RMSE均在0.7以下，SD和MAPE也在9.1以內(nèi)，結(jié)合表4中的相關(guān)系數(shù)、殘差平方和及擬合優(yōu)度，可以看出模型誤差較小。

表5 統(tǒng)計(jì)參數(shù)分析結(jié)果表

式中n表示模擬的組數(shù)；ecal、esim分別表示計(jì)算、數(shù)值模型的井筒砂粒濃度。

x的表達(dá)式為：

為了更好地描述水合物試采局部復(fù)雜井段沉砂的嚴(yán)重程度，方便計(jì)算和預(yù)測不同條件下的砂粒臨界不沉積水速，定義砂沉積濃度比這一概念，即：試采后井筒內(nèi)的含砂濃度C與地層出砂濃度C0（與降壓方式、地層出砂速率等有關(guān)）之比。即

當(dāng)e＝1時(shí)，說明進(jìn)入井筒內(nèi)的砂量等于運(yùn)移出井筒的砂量，此時(shí)的水速即對(duì)應(yīng)為臨界不沉積水速；當(dāng)e＞1時(shí)，說明井筒內(nèi)有沉積，且e值越大，沉積越嚴(yán)重。

定義砂沉積濃度比后，利用式（20），可計(jì)算得到井筒內(nèi)的沉砂濃度，進(jìn)一步利用式（27）計(jì)算砂沉積濃度比，直觀顯示井筒內(nèi)的沉砂情況。另一方面，令e＝1，可以反推計(jì)算得到不同情況下復(fù)雜管路中的砂粒臨界不沉積水速，從而為安全合理安排水合物生產(chǎn)制度、降壓幅度提供參考。

3 結(jié)論

1）螺旋段內(nèi)微米級(jí)砂粒沉積情況隨著水速的增加而逐漸改善，其中螺旋段上部的砂粒清潔難度要大于螺旋段下部。

2）提出了一種微米級(jí)砂粒臨界不沉積水速的確定方法，得到了3種粒徑3種濃度下水合物試采局部復(fù)雜管路微米級(jí)砂粒的臨界不沉積水速。

3）微米級(jí)砂粒臨界不沉積水速隨著粒徑和地層出砂濃度的增大而增大，且粒徑較大時(shí)濃度差別造成的臨界不沉積水速變化更加明顯，濃度較大時(shí)粒徑差別造成的臨界不沉積水速變化更大。

4） 提出了砂沉積濃度比的概念，當(dāng)濃度比為1時(shí)，井筒內(nèi)沒有砂粒沉積，當(dāng)濃度比大于1時(shí)，井筒內(nèi)有沉積，且濃度比越大沉積越嚴(yán)重。

5）得到了水合物試采局部復(fù)雜井段井筒沉砂濃度預(yù)測模型，計(jì)算不同條件下的砂粒臨界不沉積水速，為安全合理安排水合物生產(chǎn)制度、降壓幅度提供有益參考。