賈英,湯昌福,黃磊
(1.中國石化石油勘探開發(fā)研究院,北京 100083;2.安徽省煤田地質局勘查研究院,安徽 合肥 230088;3.中國石油勘探開發(fā)研究院,北京 100083)
基于迎風差分格式的TVD聚合物驅數(shù)學模擬
賈英1,湯昌福2,黃磊3
(1.中國石化石油勘探開發(fā)研究院,北京 100083;2.安徽省煤田地質局勘查研究院,安徽 合肥 230088;3.中國石油勘探開發(fā)研究院,北京 100083)
聚合物驅替過程中存在飽和度和聚合物質量濃度2個間斷,在2個間斷之間有“油墻”形成。文中分別運用特征線法(MOC)、一階迎風差分(SPU)方法和顯式全變差遞減差分(TVD)方法對聚合物驅替過程進行求解,通過與特征線解進行對比,考察了SPU方法和TVD方法的數(shù)值模擬效果。研究表明:SPU方法的精度低,數(shù)值耗散嚴重,對2個間斷特別是聚合物質量濃度間斷的捕捉能力較差;高階精度TVD方法分辨力高,但在“油墻”處會產生非物理“振蕩”。利用迎風的思想,對傳統(tǒng)TVD方法進行了改進。數(shù)值實驗表明,改進的TVD方法在保持高精確度的同時,還能有效地防止非物理“振蕩”的產生,具有較好的穩(wěn)定性。
聚合物驅;特征線法;一階迎風格式;TVD;非物理“振蕩”
聚合物驅能夠有效改變驅替流度比,擴大掃油面積,提高波及系數(shù)和采出程度。近年來,聚合物驅提高采收率技術得到了大力推廣和應用[1-9]。
聚合物驅是一個復雜的物理過程,在特定的情形下,可以運用特征線法MOC(Method of Characteristics)進行解析求解[2-3],然而,在大部分情況下,聚合物驅都需要借助于數(shù)值方法進行研究[3,7]。在油藏數(shù)值模擬中,對流占優(yōu)滲流方程主要采用單點上游迎風 SPU(Single Point Upstream)方法求解[10]。然而,SPU方法只是一階精度,具有較強數(shù)值耗散作用,驅替過程出現(xiàn)的激波(間斷)可能會被抹平[11-16]。高階精度差分方法雖可提高數(shù)值模擬的精度,但直接使用容易在間斷附近產生非物理“振蕩”[14]??傋儾钸f減(TVD)方法通過通量限制因子控制模擬過程中總變差的變化,在獲取高階精度的同時,還能有效防止非物理“振蕩”的出現(xiàn)[13-14];因而,TVD方法在油藏數(shù)值模擬中已得到推廣和應用[11-15]。
本文首先運用MOC方法求解一維聚合物驅過程,通過與MOC解進行對比,進一步分析了聚合物驅SPU和TVD數(shù)值解的效果,提出了改進的TVD數(shù)值通量的重構方法,在保證了傳統(tǒng)TVD方法的高階精度的同時,還防止了非物理“振蕩”的出現(xiàn)。
聚合物驅主要通過增加水相黏度減小驅替流度比來提高采收率。在聚合物驅替過程中,還伴有吸附、水相滲透率下降等現(xiàn)象[1-3]。本文主要考慮了水相黏度的變化及巖石吸附對聚合物驅替的影響。并且假設:均質不可壓等溫流動;油水兩相流動,聚合物溶于水中;驅替過程對流占優(yōu),不考慮擴散作用,且忽略重力和毛細管力對驅替過程的影響。聚合物溶液的黏度是聚合物質量濃度的函數(shù),可用修正的水相黏度方程表示:
聚合物驅的吸附為等溫吸附過程,可以用線性吸附或Langmuir吸附表示。為研究方便,假設線性吸附過程,即聚合物吸附是聚合物質量濃度的線性函數(shù):
一維聚合物驅數(shù)學模型的水相連續(xù)方程和聚合物質量濃度方程分別為
將式(3)、式(4)無量綱化,得:
將式(6)、式(7)二元雙曲型方程組寫成矩陣形式:
式(8)的邊界條件和初始條件分別為
聚合物溶液的黏度隨聚合物質量濃度呈線性變化。關于聚合物質量濃度的多項式H(C)為
研究中,取a1=1 cm3/mg,a2=9 cm3/mg,μw=1 mPa·s,μo=10 mPa·s。
采用Corey型相滲曲線,即:
取Swi=0.2,Sor=0.3,Kw=0.6,Ko=1.0,no=2,nw=2。
設注入聚合物溶液質量濃度Cin=1 mg/cm3,吸附相關常數(shù)D=0.1。油、水相對滲透率曲線及水相分流量曲線見圖1。
圖1 油、水相對滲透率及水相分流量曲線
式(8)和式(9)組成了一維聚合物驅的Riemann問題,可以采用MOC求解。
Riemann問題的解可以由簡單波描述,簡單波由稀疏波和激波(或接觸間斷)組成。式(8)存在2個特征值,分別為
圖2(τ=0.1)為一維聚合物驅Riemann問題MOC模擬結果??梢钥闯?,聚合物驅過程中存在兩族波,分別為S族和C族。其中:S族波為飽和度間斷,聚合物質量濃度C為常數(shù);C族波由接觸間斷和稀疏波組成,在接觸間斷處,聚合物質量濃度C和飽和度Sw都是間斷的,但其特征速度λc為常數(shù)。飽和度間斷和聚合物質量濃度間斷分別決定了聚合物驅見水和見聚的時間。在兩族波之間存在“油墻”,飽和度保持為常數(shù)。
圖2 MOC模擬飽和度和質量濃度空間分布情況
一維聚合物驅Riemann問題的MOC解具體可參考文獻[2-3]。該計算結果與巖心實驗CT掃描結果相一致(見圖3)。
圖3 水飽和度曲線及CT分析數(shù)據(jù)[16]
當初始條件和注入條件保持穩(wěn)定時,一維聚合物驅可以運用特征線法(MOC)快速準確地求解。然而,實際應用過程中,飽和度和聚合物質量濃度的初始分布一般都不均勻,且注入條件也經常發(fā)生變化,此時常需要采用數(shù)值方法求解,所采用的數(shù)值方法應能夠準確地捕捉聚合物驅中出現(xiàn)的物理現(xiàn)象。
式(8)的顯式差分方法守恒差分格式可以表示為
不同的通量限制因子(Φ(r))構成不同的差分方法。
3.1 一階迎風(SPU)差分方法
在油藏數(shù)值模擬中,由于具有明確的物理特點,經常采用單點上游加權迎風格式SPU[10]求解,即數(shù)值流通量采用其上游的網格值計算:
圖4(τ=0.1)為不同網格數(shù)時SPU解與MOC解析解之間的對比??梢钥闯觯壕W格數(shù)較少時,聚合物驅過程中的2個間斷——特別是質量濃度間斷——都被不同程度抹平了,SPU的數(shù)值模擬結果嚴重偏離精確解;隨著網格數(shù)的增加,SPU數(shù)值解趨于MOC解,此時為獲得一個較為精確的差分解,需要更多網格。
圖4 一維聚合物驅不同網格數(shù)SPU差分模擬結果
上述結果表明,SPU方法對聚合物驅過程中的2個間斷特別是質量濃度間斷的捕捉能力較弱。這是因為:接觸間斷與激波不同,它不具有“自銳利”(Selfsharpening)行為,精確解的一些特征被差分方法具有的數(shù)值耗散掩蓋,需要采用其他方法(例如粒子追蹤方法、高階精度差分方法等[11,17])才能改進模擬的結果。
3.2 高分辨率顯式全變差遞減(TVD)差分方法
TVD差分方法能提高數(shù)值模擬的的精度,然而直接采用該方法容易在間斷附近產生非物理“振蕩”。為了定量研究“振蕩”情形,定義時間步n的總變差[11]為
TVD差分方法通過限制總變差的增加,防止非物理“振蕩”的產生,即:
此時式(14)中的通量限制因子Φ為梯度r的函數(shù),使得在模擬區(qū)域的不同部分需采用不同的差分格式。目前,學者們提出了各種限制因子,如Min mod限制因子、Fromm限制因子、Van Leer限制因子以及Leonard限制因子等[8]。本文采用Van Leer限制因子的TVD差分方法研究一維聚合物驅過程:
圖5(τ=0.1)為一維聚合物驅TVD差分方法模擬結果。通過與SPU解以及MOC解進行對比可以發(fā)現(xiàn),聚合物驅數(shù)值模擬的效果得到顯著改善,在較少的網格數(shù)時就能獲得較為精確的數(shù)值解 (N=400)。因此,TVD差分方法捕捉間斷的能力要強于SPU方法,聚合物驅差分方法應優(yōu)先選擇TVD差分方法,保證在較少的網格時得到更加精確的數(shù)值解。
圖5 一維聚合物驅TVD差分模擬結果
然而,從圖5還可以發(fā)現(xiàn),TVD差分方法在提高精度的同時,卻在 “油墻”處產生了非物理 “振蕩”。Thiele等研究發(fā)現(xiàn),增加時間差分精度、控制CFL的大小或者采用耗散性較強的限制因子(如Min mod)可以減緩“振蕩”,但并不能達到消除“振蕩”的目的。
3.3 迎風TVD差分方法
傳統(tǒng)TVD差分方法在“油墻”處會出現(xiàn)非物理“振蕩”?;谟筒財?shù)值模擬中常用的迎風思想,本文提出了新的數(shù)值流通量的計算方法,以抑制“振蕩”的發(fā)生。在新的TVD差分方法中,飽和度方程中的數(shù)值流通量采用迎風格式,而質量濃度方程中的總變差通過質量濃度的TVD來限制,即:
聚合物質量濃度方程中數(shù)值通量依式(21)計算:
圖6(τ=0.1)給出了改進的TVD解(U-TVD)與MOC解及傳統(tǒng)TVD解的對比??梢园l(fā)現(xiàn),迎風TVD差分方法在保持TVD差分方法高階精度的同時,還有效地抑制了“油墻”處非物理“振蕩”的出現(xiàn)。
1)一階迎風(SPU)方法模擬聚合物驅能力較弱。網格數(shù)較少時,數(shù)值解嚴重偏離精確解;雖然可通過增加網格數(shù)改善數(shù)值模擬結果,但增加的計算代價巨大。
2)TVD差分方法能以較少網格獲得較精確的數(shù)值解,聚合物驅數(shù)值模擬應盡量采用高分辨率差分方法。
3)傳統(tǒng)TVD差分方法會在飽和度和質量濃度間斷之間存在的“油墻”處產生不合理的非物理“振蕩”。基于迎風思想,改進了TVD數(shù)值通量計算方法,在保證TVD方法的高精度同時,還能防止“振蕩”產生,使得TVD方法更為穩(wěn)定。該方法為下一步開展三維聚合物驅數(shù)值模擬奠定了基礎。
圖6 一維聚合物驅迎風TVD差分模擬結果
μp,μw,μo分別為聚合物溶液、水相和油相黏度,mPa·s;C為聚合物在水相中的質量濃度,mg/cm3;Ca為聚合物吸附質量濃度,mg/cm3;Cin為聚合物溶液質量濃度,mg/cm3;a為吸附常數(shù);x為空間坐標,m;t為時間坐標,s;ρr,ρw分別為巖石和水相密度,kg/m3;q為注入流量,m3/s;Sw為水飽和度;Swi,Soi分別為束縛水飽和度和殘余油飽和度;φ為孔隙度;A為截面積,m2;f為水相分流量函數(shù) (即水飽和度和聚合物質量濃度的函數(shù));Krw,Kro分別為水相和油相的相對滲透率;Kw,Ko分別為殘余油飽和度端點水相滲透率和束縛水飽和度端點油相滲透率;τ為注入孔隙體積倍數(shù);Vp為孔隙體積,m3;L為模型長度,m;xD為無量綱長度;D是吸附相關常數(shù);a1,a2為常數(shù) (通過實驗取值);nw,no分別為水相、油相相對滲透率曲線Coery系數(shù);n和i分別表示時間和空間節(jié)點(或網格位置);N為空間網格數(shù);和分別為網格i的左、右數(shù)值流通量;Φ()r為通量限制因子;r為梯度;下標分別表示上一網格和下一網格;表示在網格位置的質量濃度差比值。
[1]NEEDHAM R B,DOE P H.Polymer flooding review[J].Journal of Petroleum Technology,1987,39(12):1503-1507.
[2]POPE GARY A.The application of fractional flow theory to enhanced oil recovery[J].SPE Journal 1980,20(20):191-205.
[3]PATTON J T,COATS K H,COLEGROVE G T.Prediction of polymer flood performance[R].SPE 2546,1971.
[4]蔣明.齊40塊稠油油藏聚合物驅可行性室內實驗研究[J].斷塊油氣田,1998,5(1):45-49.
[5]宋考平,楊二龍,王錦梅,等.聚合物驅提高驅油效率機理及驅油效果分析[J].石油學報,2004,25(3):71-74.
[6]方文超.陸上稠油油田多輪次聚合物驅提高采收率方法[J].斷塊油氣田,2015,22(5):619-622.
[7]蘆文生.綏中36-1油田聚合物驅數(shù)值模擬研究[J].中國海上油氣,2002,16(5):333-340.
[8]周叢叢.大慶油田一、二類油層聚合物驅注采指標變化規(guī)律[J].斷塊油氣田,2015,22(5):610-613.
[9]侯健,李振泉,王玉斗,等.考慮擴散和吸附作用的聚合物驅替過程滲流數(shù)值模擬[J].計算物理,2003,20(3):239-244.
[10]MALLISON B T,GERRITSEN M G,JESSEN K,et al.Higher-order upwind schemes for two-phase,multi-component Flow[R].SPE 79691,2005.
[11]THIELE M R,EDWARDS M.Physically based higher order godunov schemes for compositional simulation[R].SPE 66403,2001.
[12]LEER V B.Towards the Ultimate conservative difference Scheme,V-A second-order sequel to Godunov′s method[J].Journal of Computational Physics,1979,32(1):101-136.
[13]HARTEN A.High resolution schemes for hyperbolic conservation laws[J].Journal of Computational Physics,1983,49(3):357-393.
[14]RUBIN B,EDWARDS M G.Extension of the TVD mid-point scheme to higher order accuracy in time[R].SPE 25265,1993.
[15]BRENIER Y,JAFFRé J.Upstream differencing for multiphase flow in reservoir simulation[J].Siam Journal on Numerical Analysis,1991,28(3):685-696.
[16]SU S.Laboratory investigation of polymer-flooding,ASP-flooding,and foam-flooding in heavy oil system[D].Calgary:University of Calgary,2012.
[17]JUANES R,LIE K A.Numerical modeling of multiphase first-contact miscible flows.Ⅱ.Front-tracking/streamline simulation[J].Transport Porous Media,2008,72(1):97-120.
(編輯 李宗華)
Mathematical simulation of polymer flooding with upstream-TVD solution
JIA Ying1,TANG Changfu2,HUANG Lei3
(1.Petroleum Exploration and Production Research Institute,SINOPEC,Beijing 100083,China;2.Exploration Research Institute, Anhui Provincial Bureau of Coal Geology,Anhui 230088,China;3.Research Institute of Petroleum Exploration&Development, PetroChina,Beijing 100083,China)
For polymer flooding,there are two discontinuities of saturation and concentration,between which"oil bank"exists.The polymer flooding was investigated with method of characteristics(MOC),single point upstream(SPU)finite difference method and explicit total variation diminishing(TVD)finite difference method respectively.Based on the solution of MOC,performance of simulation results of both SPU and TVD were evaluated.It is suggested TVD exhibits a better numerical performance while the SPU scheme suffers numerical diffusion badly.The conventional TVD method is modified on the calculation of numerical flux,so that non-physical oscillation can be prevented on the"oil bank"between the two discontinuities.Numerical tests demonstrate that the upstream-TVD can achieve better stability than the conventional TVD method.
polymer flooding;MOC;SPU;TVD;non-physical oscillation
TE357.46+1
A
10.6056/dkyqt201703017
2016-12-12;改回日期:2017-03-09。
賈英,女,1980年生,高級工程師,博士,油氣田開發(fā)工程專業(yè)。E-mail:jiaying.syky@sinopec.com。
賈英,湯昌福,黃磊.基于迎風差分格式的TVD聚合物驅數(shù)學模擬[J].斷塊油氣田,2017,24(3):373-377.
JIA Ying,TANG Changfu,HUANG Lei.Mathematical simulation of polymer flooding with upstream-TVD solution[J].Fault-Block Oil&Gas Field,2017,24(3):373-377.