壓裂井組非線性滲流模型求解

黃迎松

（中國石化勝利油田分公司勘探開發研究院，山東東營 257015）

低滲透油藏中流體的流動不再符合線性達西滲流規律，只有當壓力梯度大于啟動壓力梯度后流體才會流動。隨著驅替壓力梯度增大，流體在低滲透地層中的滲流過程經歷非線性和擬線性2個滲流階段。目前，在研究低滲透油藏滲流時，大都將非線性和擬線性滲流階段簡化為超過擬啟動壓力梯度后的擬線性滲流段，忽略了非線性滲流階段對流體流動的影響[1-3]。尹芝林等人[4]基于動態滲透率的概念，采用統一的運動形式描述非線性和擬線性滲流階段的滲流規律，認為與采用擬線性滲流規律相比，采用非線性滲流規律計算出的壓力變化更為平緩。綜合考慮前人的研究成果，筆者提出分段描述非線性、擬線性滲流規律，同時考慮壓裂裂縫中高速流動特征的求解方法，以精細描述低滲透油藏的流動規律，精確刻畫地層中的壓力分布，提高低滲透壓裂注采井組生產指標的計算精度。

開發低滲透、特低滲透油藏時一般采用水力壓裂，學者們對此進行了很多研究[5-7]，有部分學者將裂縫和地層看成統一系統，但這種處理會造成模擬不準確；壓裂裂縫中流體的流速快，可能會出現高速非達西滲流，需要根據雷諾數判斷流動形態[8]。此外，裂縫的導流能力沿縫長及隨裂縫中壓力的變化而變化[9-10]。基于以上特征，筆者考慮非線性和擬線性滲流段的滲流特征，建立了低滲透三維油水兩相達西滲流和高速非達西滲流耦合的數學模型，采用有限差分法進行了求解，并分析了計算結果。

1 壓裂注采井組耦合數學模型的建立

1.1 基本假設

建立模型前，進行以下假設：1）油藏內流體的流動為等溫流動，油藏外邊界封閉；2）地層巖石和流體微可壓縮；3）三維地層中有油水兩相參與滲流，分油藏和壓裂裂縫2個區域分別建立滲流方程；4）油藏區域考慮非線性和擬線性2種滲流的特征，認為油相和水相的啟動壓力梯度為常數；5）裂縫系統考慮達西和高速非達西流動，且考慮裂縫導流能力的時變特征；6）考慮重力和毛細管力的影響。

1.2 油藏系統的滲流方程

圖1為低滲透地層中流體流速與壓力梯度的關系曲線（圖中，GA為流體的啟動壓力梯度，GB為線性流動階段的擬啟動壓力梯度，GC為開始呈現線性流動時的壓力梯度，單位均為10-1MPa/cm）。

當驅替壓力梯度不小于GA且不大于GC時，流體流速-壓力梯度曲線呈現下凹的非線性段，可用二次函數近似描述非線性段：

圖1 低滲透地層中流體流速與壓力梯度的關系Fig. 1 Relationship between fluid flow velocity and pressure gradient in low permeability formation

當驅替壓力梯度大于GC時，流體流速-壓力梯度關系曲線為一條直線，可描述為：

式中：p為壓力，10-1MPa；v為流體流速，cm/s；K為滲透率，D；μ為流體黏度，mPa·s；a，b和c為二次函數方程的系數。

考慮三維油水兩相流動，且油藏區域與裂縫之間存在交互流動項，根據式（1）和式（2），可得到油藏系統的滲流方程。

式中：D為油藏深度，cm；Kr為相對滲透率；ρ為流體密度，g/cm3；φ 為油藏孔隙度；S 為飽和度；τlmf為油藏和壓裂裂縫系統交互流動項，g/（cm3·s）；qlm為單位時間單位體積的產量項，g/（cm3·s）；下標 l=o，w；o，w分別代表油相和水相；下標m代表油藏；下標f代表裂縫；為壓力梯度，10-1MPa/cm；t為時間，s；g 為重力加速度，m/s2。

1.3 裂縫系統的滲流方程

因為裂縫寬度較小，可以忽略流體在寬度y方向的流動，建立坐標系 Ox′z′，x′軸沿裂縫延伸方向，z′軸與油藏坐標系的z軸相同，當其與油藏系統坐標系不產生混淆的情況下也可將裂縫坐標系記為（x，z）。裂縫中流體的流動形態用卡迪雷夫雷諾數來判斷。

裂縫中流體流動的卡迪雷夫雷諾數計算式為：

裂縫中流體的運動方程為：

β為非達西因子，由介質參數孔隙度和滲透率決定，可表示為：

裂縫系統的連續性方程為：

τlmf為油藏和裂縫系統的交互流動項，可表示為：

σ是基質塊形狀因子，由基質巖塊形狀的維數及其特征長度決定：

式中：d為裂縫面的維數；L為基質部分的特征長度，m。

根據雙重介質關于形狀因子的計算方法，利用達西公式推導油藏和裂縫系統的交互流動項的計算公式：

式中：Re為流體雷諾數；Kf為裂縫滲透率，D；φf為裂縫孔隙度；wf為裂縫的寬度，m；Df為裂縫穿過油藏網格的長度，m；Dx為油藏網格x方向的步長，m；Dy為油藏網格y方向的步長，m。

M. Y. Soliman[9]研究發現，裂縫的導流能力隨著縫長增長呈線性降低或指數降低。水力裂縫隨裂縫中壓力的變化張開或閉合，該過程中裂縫的導流能力也會發生改變。因此，裂縫的導流能力是縫長和壓力的函數，函數形式與裂縫性質有關。

式中：i為生產過程中的時刻序號，i0為初始生產時刻。

1.4 輔助方程

考慮油水兩相流動，必須滿足以下輔助方程。

油藏系統：

1.5 定解條件

油藏系統的初始條件為：

裂縫系統的初始條件為：

油藏外邊界為封閉條件：

井底邊界采用定壓條件：

式中：Q為油井地面產量，cm3/s；B為流體體積系數。

2 數學模型的求解

采用有限差分法求解壓裂注采井組耦合數學模型，采用順序解法求解壓力和飽和度。為提高計算結果的精度，計算時要減小時間步長，以達到工程計算精度的要求。非線性方程組先進行線性化處理，得到線性差分方程組，避免使用迭代法求解非線性差分方程組，以減少計算量和提高計算速度，并給出了油藏系統和裂縫系統的具體差分格式。

2.1 油藏系統的差分格式

式中：Cp和Cl為巖石和流體的壓縮系數，10-1MPa-1。

以第一項為例說明差分處理方法：

滲流方程后2項的差分格式為：

y和z方向的差分處理方式同上，不再贅述。

2.2 裂縫系統的差分格式

裂縫系統在x方向上采用非等距網格，根據裂縫的長度和方位及壓裂井的坐標，可求得裂縫系統通過的油藏網格和所穿過的長度，每個油藏網格對應裂縫系統的一個網格，油藏和裂縫之間的流動交互項只存在于相應的網格，通過交互項建立耦合方程組。

從而，裂縫系統滲流方程x方向的差分格式為：

同理可以寫出z方向的差分格式，結合竄流、產量項以及右端項的差分可以得到裂縫高速非達西滲流的差分格式。油藏系統和裂縫系統相交的網格上才會有竄流交互項，竄流項的差分按照竄流公式離散即可。

由于主要考慮的是油水兩相流動，在求解飽和度時可以只求解油相或者水相的飽和度，然后根據輔助方程求出另外一相的飽和度。水相的飽和度可以采用下式求解：

式中：ΔAwn為采用顯示處理方法得到的n+1時刻壓力項前的系數；GWWT為式（27）中的所有n時刻項的組合；Vp為網格塊的孔隙體積，Vp= ΔxΔyΔzφ。

2.3 求解框圖

根據上述差分模型，可以將之轉化為計算機模型進行求解分析，求解框圖見圖2。

3 算例分析

3.1 算例計算結果

圖2 求解過程框圖Fig. 2 Solution process block diagram

選擇五點法井網中的1個壓裂注采井組，只壓裂中心注水井，不壓裂采油井，假設裂縫為雙翼對稱裂縫，裂縫延伸方向與地層最大水平主應力方向一致。井組的基本參數為：儲層厚度8.50 m，孔隙度0.21，地層深度1 000.00 m，原始地層壓力10.0 MPa，儲層平均滲透率4.0 mD，原油黏度3. 5 mPa·s，地層水黏度0.45 mPa·s，井距200.00 m，裂縫半長 80.00 m。根據上文建立的壓裂注采井組耦合數學模型編制求解程序，計算該注采井組的動態。該注采井組的網格系統平面如圖3所示。

圖4所示為該井組注水開發90 d后的壓力和飽和度計算結果。從圖4可以看出：壓力和飽和度圍繞壓裂裂縫形成等值線，沿著裂縫兩翼呈現對稱分布；裂縫附近的等壓線密集，距離裂縫越遠，等壓線越稀疏；靠近油井井底時，等壓線又逐漸變密。同樣，注入水也是沿著壓裂裂縫逐漸向外擴散。計算結果顯示的壓力和飽和度的分布與地層實際情況相符，說明建立的模型是正確的。

3.2 壓裂裂縫對注采井組動態的影響

圖3 五點法井網網格系統平面示意Fig. 3 Schematic diagram of the five-spot well pattern system

圖4 壓裂注采井組壓力和飽和度的計算結果Fig. 4 Calculation results of the pressure and saturation of fractured injection-production well group

為了研究壓裂裂縫對注采井組生產動態的影響，首先計算了考慮裂縫與不考慮裂縫時的注采動態，結果見圖5。從圖5可以看出：注水井未壓裂時，油井的產油量迅速降低，維持在很低的水平；注水井壓裂后，產油量先降低后升高，與注水井未壓裂時相比產油量明顯增高，表明壓裂注水井具有明顯的增產效果，其主要原因是壓裂裂縫提高了注水井的吸水能力，油藏能量得到補充，保持地層壓力不降低，生產壓差穩定。

圖5 注水井壓裂對注采井組產量的影響Fig. 5 Effect of fracturing in water injection well on the production of injection-production well group

3.3 壓裂裂縫導流能力對注采井組動態的影響

由文獻[10]可知，壓裂裂縫導流能力是地層閉合壓力的函數。目前，壓裂井進行數值模擬時多認為裂縫導流能力隨時間變化[11-16]，但實際上壓裂裂縫導流能力的變化主要還是由閉合壓力變化引起的，因此，筆者將導流能力設計成裂縫中壓力的函數。在此基礎上對比文獻算法和本文算法的計算結果，結果見圖6。

圖6 裂縫導流能力對壓裂注采井組產量的影響Fig. 6 Effect of fracture flow conductivity on the production of fractured injection-production well group

由圖6可以看出，2種算法計算的2條產量-生產時間曲線具有相同的趨勢，但本文算法計算出的產量要略低于文獻算法。文獻算法中，注水井壓裂裂縫的導流能力隨生產時間增長而逐漸降低；而本文算法中，注水井壓裂裂縫導流能力的變化取決于壓裂裂縫中壓力的變化：裂縫中的壓力大于閉合壓力時，裂縫導流能力隨裂縫中壓力降低而降低；裂縫中壓力小于等于閉合壓力時，此時裂縫的導流能力趨近于零。注水井壓裂裂縫的導流能力決定了注采井組補充地層能量的能力，因此也影響了整個注采井組的產量。這也是本文計算結果與文獻計算結果存在一定的差異的原因。

4 結 論

1）根據低滲透油藏非線性滲流特征和壓裂裂縫導流能力的變化規律，建立了低滲透油藏壓裂井組油藏和裂縫耦合的非線性數學模型。該模型綜合考慮了低滲透油藏非線性和擬線性滲流規律，以及壓裂裂縫出現達西和非達西滲流的情況，模擬結果與地層實際情況相符，說明該模型正確有效。

2）提出了利用Taylor展開將油藏系統非線性滲流方程轉化為線性差分方程組的方法，編制了計算機求解程序，計算結果證明該方法有效，能夠用來模擬低滲透油藏壓裂井注采組的生產動態，與其他非線性模型相比，解法簡單。

3）低滲透油藏具有非線性滲流特征，綜合考慮非線性段和擬線性段滲流比僅考慮擬線性滲流更加符合地層實際滲流情況，所以不能忽略非線性段滲流的影響。

4）五點法井網壓裂注采井組的計算結果表明，注水井壓裂可以產生明顯的增產效果，壓裂裂縫的導流能力隨著裂縫中壓力的變化而變化，裂縫導流能力降低，造成注水井注水效果變差，導致油井產量降低。