多孔介質中單相流地應力耦合問題的二重網格混合有限元法

摘要： 本文提出了多孔介質中單相流耦合地應力問題的二重網格混合有限元法. 本文首先給出了多孔介質中的單相流的控制方程及地應力問題的控制方程，推導了多孔介質中的微可壓縮流體的單相流耦合地應力問題的非線性控制方程. 然后，為數值求解該方程本文提出了二重網格混合有限元法，將問題轉化為在粗網格上求解小規模非線性問題、在細網格上求解大規模線性問題，以提高計算效率. 最后，數值算例驗證了方法的有效性.

關鍵詞： 多孔介質； 單相流耦合地應力問題； 二重網格法； 混合有限元

中圖分類號： O241. 82 文獻標志碼： A DOI： 10. 19907/j. 0490-6756. 2024. 041002

1 引言

油氣藏數值模擬是一個復雜交叉學科問題，需要地質工程、滲流學、計算數學、熱力學和計算機科學等多學科共同研究. 油氣藏數值模擬需要通過計算機軟件來模擬地下油氣在不同條件下的變化，以便為制定開采方案提供參考，達到油氣高產的目標或預測油氣產量.

目前，多數油氣藏數值模擬都只關注多孔介質中的油氣流動. 但是，在實際的頁巖油氣儲層中，流場和地應力場之間存在重要的相互作用，因而在構建模型及數值模擬時考慮地應力的影響很重要. 對于這個問題，不同領域學者均有研究，這些研究共同推動人們對地應力變化規律及流場與地應力場的耦合相互作用的認識［1-5］. Gudala 等［1］使用單相流與地質力學耦合模型研究了縫隙多孔介質中單相流的流動，并對模型進行了數值求解.Li 等［4］把巖體看作多孔介質，計入溫度和孔隙率的影響，建立了高溫高壓及高孔隙壓環境下的巖體多孔介質流變模型. Chen 等［5］基于計算流體力學和巖土力學，考慮流場和應力場的相互作用建立了流-固耦合模型，等.

在數學模型基礎上，數值模擬還需要合適的計算方法與算法設計. 其中，用混合有限元法數值求解耦合地應力問題已有一些研究［6-9］. Zhang 等［7］使用多尺度混合有限元法研究了可壓縮流體在多孔介質中的滲流. Li 等［8］建立了單相滲流問題的有限元方程，并針對復雜邊界情況用有限元法進行了模擬計算，給出了混合邊界油藏壓力分布的規律. 進一步，鑒于單相流耦合地應力模型的非線性可能導致計算效率低，并考慮到地應力模擬時間較長、數值模擬較慢，研究者提出了二重網格混合有限元法［10-13］. Zhang 等［10］提出了一種求解半線性反應擴散方程的二重網格混合有限元法. Rui 等［11］提出了多孔介質Darcy-Forchheimer 流動的二重網格塊中心差分法，等.

本文研究多孔介質中單相流的耦合地應力問題，提出了一個二重網格混合有限元格式，構造了格式對應的算法，并用數值算例驗證了方法的有效性.

2 模型

在推導問題對應的控制模型前，我們假設：（i）在耦合流動過程中溫度恒定；（ii）考慮基巖和流體的壓縮性時，基巖與裂縫中的流體流動均滿足達西定律；（iii）巖石變形為微小線彈性變形；（iv）不考慮慣性力，并忽略毛細管力與重力作用.

準確模擬多孔介質中單相流耦合地應力問題需要研究流體流動與地應力的耦合問題，以及巖石結構因地應力而產生的線彈性變形問題. 式（1）為多孔介質中單相流的控制方程.

得到線性問題（即方程（13））的二重網格混合有限元法的解（ vn + 1h ，pn + 1h ，un + 1h ）∈ Vh × Wh × Uh.

5 數值算例

考慮10 m × 10 m 的二維區域，區域的右邊界和上邊界牽引力為5 N. 對該區域采用64 × 64 的一致三角形單元離散，如圖1 所示，對應的巖石流體參數如表1 所示. 假設該區域有5 口生產井存在，坐標分別為［5，5］，［2. 5，2. 5］，［2. 5，7. 5］，［7. 5，2. 5］及［7. 5，7. 5］，生產速率為 2 × 10-7 s-1，模擬時間為100 d，時間步長為1440 min，即1 d.

下面我們先使用牛頓迭代法模擬運行100 d，模擬運行100 d 后區域中的壓力及孔隙率如圖2 所示，速度及剪切應力如圖3 所示，最大主應力及最小主應力如圖4 所示，x 方向位移及y 方向位移如圖5 所示.

下面我們使用二重網格混合有限元方法再次進行計算. 首先，在粗網格上進行8 × 8 的一致三角形單元離散，再劃分同圖1 的64 × 64 細網格.我們先在粗網格上求解非線性問題，再在細網格上求解線性問題，得到數值解. 同樣，我們模擬100 d，時間步長取為1 d.

模擬運行100 d 后，區域中的壓力及孔隙率如圖6 所示，速度及剪切應力如圖7 所示，最大主應力及最小主應力如圖8 所示，x 方向位移及y 方向位移如圖9 所示.

對比2 種方法的結果，可以看到二者幾乎相同. 另外，圖10 展示了2 種方法的計算時間，可見2種方法所需計算時間與模擬天數線性相關，但二重網格混合有限元法的效率明顯高于牛頓迭代法.

6 結論

本文提出用二重網格混合有限元法來求解多孔介質單相流耦合地應力問題，此類問題涉及的模型為單相流耦合地應力模型，模擬的時間尺度大，計算量也大，計算耗時非常長. 我們分別用牛頓迭代法和二重網格混合有限元法對同一問題進行求解，并對比2 種方法的結果和計算時間，發現二重網格混合有限元方法的效率明顯高于牛頓迭代法.

參考文獻：

［1］ Gudala M， Govindarajan S K. Numerical modellingof coupled single-phase fluid flow and geomechanicsin a fractured porous media ［J］. J Petrol Sci Eng，2020， 191： 107215.

［2］ Gao C， Gray K E. A coupled geomechanics and res?ervoir simulator with a staggered grid finite differencemethod［ J］. J Petrol Sci Eng， 2022， 209： 109818.

［3］ Shojaei M J. Pore-scale dynamics of foam flow in porousmedia［ D］. Manchester： University of Manchester，2020.

［4］ Li J， Wang Y， Wang H. A study on the rheologicalmodel of porous media in deep rock mass ［J］. RockSoil Mech， 2008， 29： 2355.

［5］ Chen W， Pei W， Li S， et al. Numerical model andengineering application of fluid structure coupling insaturated porous media ［J］. Chin J Rock Mech Eng，2013（z2）： 3346.

［6］ Zidane A， Firoozabadi A. An implicit numericalmodel for multicomponent compressible two-phaseflow in porous media ［J］. Adv Water Res， 2015，85： 64.

［7］ Zhang N， Yao J. Research on multiscale mixed finiteelement numerical method for compressible fluid flow［J］. Chin J Comput Mech， 2017， 34： 226.

［8］ Li Y， Qu X， Mao Y. Finite element simulation ofseepage pressure distribution in complex boundary oilreservoirs ［J］. J Shandong Inst Petrol Chem Technol，2007， 21： 1.

［9］ Yu J， Ren Y， Cao W， et al. Extended mixed elementmethod for compressible porous media flow［J］. J Shandong Univ（ Nat Sci Ed）， 2017， 52： 25.

［10］ Zhang J， Han H， Yu Y， et al. A new two-gridmixed finite element analysis of semi-linear reaction–diffusion equation ［J］. Comput Math Appl， 2021，92： 172.

［11］ Rui H， Liu W. A two-grid block-centered finite differencemethod for Darcy-Forchheimer flow in porousmedia ［J］. SIAM J Numer Anal， 2015， 53：1941.

［12］ Hu H. Two-grid method for compressible miscibledisplacement problem by mixed finite element methodsand expanded mixed finite element method ofcharacteristics ［J］. Numer Algorithms， 2022，89： 611.

［13］ Gao M. A double grid finite element method for solvingnonlinear porous elastic problems ［D］. Jinan：Shandong University， 2022.

［14］ Rutqvist J， Tsang C F. A study of caprock hydromechanicalchanges associated with CO2-injection into abrine formation［ J］. Environ Geol， 2002， 42： 296.

（責任編輯： 周興旺）

基金項目： 國家自然科學基金（11971337）