油藏數(shù)值模擬原理黑油模型

1、- 1 - 2006年年4月月13日日 李建芳李建芳 油藏數(shù)值模擬原理油藏數(shù)值模擬原理黑油模型- 2 -黑油模型簡介黑油模型簡介黑油模型是指油質(zhì)較重性質(zhì)的油藏類型。黑油模型是指油質(zhì)較重性質(zhì)的油藏類型。 黑油模型是最完善、最成熟，也是應(yīng)用最黑油模型是最完善、最成熟，也是應(yīng)用最為廣泛的模型。為廣泛的模型。是油藏數(shù)值模擬的基礎(chǔ)，其它模型大都是黑油模是油藏數(shù)值模擬的基礎(chǔ)，其它模型大都是黑油模型的擴展。型的擴展。目的：軟件這樣一個目的：軟件這樣一個“黑箱黑箱”是如何操作的，以是如何操作的，以便指導(dǎo)我們更好地去應(yīng)用軟件。便指導(dǎo)我們更好地去應(yīng)用軟件。- 3 -黑油模型簡介黑油模型簡介1. 黑油模型的基本假設(shè)

2、黑油模型的基本假設(shè)2. 單相滲流黑油模型基本微分方程單相滲流黑油模型基本微分方程3. 三維三相黑油模型基本微分方程三維三相黑油模型基本微分方程4. 初邊值條件初邊值條件5. 網(wǎng)格系統(tǒng)網(wǎng)格系統(tǒng)6. 黑油模型的差分方程黑油模型的差分方程7. 差分方程的線性化差分方程的線性化8. 線性代數(shù)方程組的求解線性代數(shù)方程組的求解9. 黑油模型的主要數(shù)據(jù)流黑油模型的主要數(shù)據(jù)流- 4 -1 、黑油模型基本假設(shè)、黑油模型基本假設(shè)（1） 油藏中的滲流是等溫滲流。油藏中的滲流是等溫滲流。（2） 油藏中最多只有油、氣、水三相，每一相均遵守達油藏中最多只有油、氣、水三相，每一相均遵守達西定律。西定律。（3）油藏?zé)N類只含有

3、油、氣兩個組分。在油藏狀態(tài)下，）油藏?zé)N類只含有油、氣兩個組分。在油藏狀態(tài)下，油氣兩組分可能形成油氣兩相，油組分完全存在于油油氣兩組分可能形成油氣兩相，油組分完全存在于油相內(nèi)，氣組分則可以以自由氣的方式存在于氣相中，相內(nèi)，氣組分則可以以自由氣的方式存在于氣相中，也可以以溶解氣的方式存在于油相中，所以地層內(nèi)油也可以以溶解氣的方式存在于油相中，所以地層內(nèi)油相為油組分和氣組分的某種組合。在常規(guī)油田中，一相為油組分和氣組分的某種組合。在常規(guī)油田中，一般不考慮油組分向氣組分揮發(fā)的現(xiàn)象。般不考慮油組分向氣組分揮發(fā)的現(xiàn)象。 （4）油藏中氣體的溶解和逸出是瞬間完成的，即認為油）油藏中氣體的溶解和逸出是瞬間完成的

4、，即認為油藏中油氣兩相瞬時達到相平衡狀態(tài)。藏中油氣兩相瞬時達到相平衡狀態(tài)。 （5）油水之間不互溶；天然氣也假定不溶于水。）油水之間不互溶；天然氣也假定不溶于水。- 5 -2.基本微分方程基本微分方程單相滲流黑油模型基本微分方程單相滲流黑油模型基本微分方程黑油模型基本微分方程主要是由黑油模型基本微分方程主要是由達西公式達西公式(達西達西定定律律)和和連續(xù)性方程連續(xù)性方程(物質(zhì)守恒定律物質(zhì)守恒定律)所構(gòu)成所構(gòu)成 一維單相滲流的一維單相滲流的Darcy公式為：公式為： 三維單相滲流的三維單相滲流的Darcy公式為：公式為：QK dPU = -Adx()KPgZ U- 6 -滲流速度是一個空間向量，用

5、分量表示為：滲流速度是一個空間向量，用分量表示為： 滲流速度矢量；滲流速度矢量； 壓力；壓力； 勢函數(shù)勢函數(shù) 絕對滲透率；絕對滲透率； 重力加速度常數(shù)；重力加速度常數(shù)； 粘度；粘度； 垂向坐標；垂向坐標； 密度；密度；2. 基本微分方程基本微分方程UK()()()xyzKPgZKxxKPgZKyyKPgZKzz UUUPgZ- 7 -Gringarten, 2002Data from all sources- 8 -2.基本微分方程基本微分方程根據(jù)物質(zhì)守恒原理，在單位時間內(nèi)根據(jù)物質(zhì)守恒原理，在單位時間內(nèi)流入單元內(nèi)的流體質(zhì)量流入單元內(nèi)的流體質(zhì)量-流出單元的流體質(zhì)量流出單元的流體質(zhì)量=單元內(nèi)流體質(zhì)

6、量的變化單元內(nèi)流體質(zhì)量的變化取滲流場中一個微小六面體體積單元來研究取滲流場中一個微小六面體體積單元來研究- 9 -212121212121jjjiiiyyyxxxjjjjjjiiiiiijjjiiiyyyyyyxxxxxxyyyxxx12112112112121212121- 10 -2.基本微分方程基本微分方程該六面體單元的中心點坐標為（該六面體單元的中心點坐標為（x，y，z），），其其每一側(cè)面的質(zhì)量流速均以其側(cè)面的中心點的質(zhì)量每一側(cè)面的質(zhì)量流速均以其側(cè)面的中心點的質(zhì)量來代替。來代替。 方向流入質(zhì)量方向流入質(zhì)量= 方向流出質(zhì)量方向流出質(zhì)量= 方向流入質(zhì)量方向流入質(zhì)量= 方向流出質(zhì)量方向流出質(zhì)

7、量= 方向流入質(zhì)量方向流入質(zhì)量= 方向流出質(zhì)量方向流出質(zhì)量= 單元體內(nèi)流體質(zhì)量的變化單元體內(nèi)流體質(zhì)量的變化= x(, , )2()xxxy zy zu x(,)2()xxxy zyzu yy(, )2()yyx yzxzu (, )2()yyx yzxzu zz(,)2()zzx y zxyu ( , ,)2()zzx y zx yu ()xyzt - 11 -根據(jù)根據(jù)物質(zhì)守恒定律物質(zhì)守恒定律：將上式兩邊同除以將上式兩邊同除以 ，2.基本微分方程基本微分方程(, , )(, , )22( , )( , )22( , ,)( , ,)22()()()()()()()xxxxxy zxy zyyy

8、yx yzx yzzzzzx y zx y zy zuy zux zux zux yux yux y zt x y z (,)(,)22(,)(,)22(,)(,)22()()()()()()()xxxxxy zxy zyyyyxyzxyzzzzzxy zxy zuuxuuyuuzt- 12 -并令并令 ， ，把源匯項（井）加入到方程中，得到把源匯項（井）加入到方程中，得到單相滲流單相滲流連續(xù)性方程連續(xù)性方程 ：2.基本微分方程基本微分方程0 x 0y0z ()()()()yxzuuuxyzt()()()()yxzuuuqxyzt- 13 -用用Hamilton算子算子 表示，并將達西公式代入

9、，得表示，并將達西公式代入，得到，到，單相滲流微分方程（流動方程）：單相滲流微分方程（流動方程）： 生產(chǎn)井生產(chǎn)井 注入井注入井 關(guān)關(guān) 井井2.基本微分方程基本微分方程()()Kqt 0q 0q 0q 流動項源匯項累積項- 14 -這就是地下流體所遵循的滲流規(guī)律。這就是地下流體所遵循的滲流規(guī)律?？傊?，滲流基本微分方程包含三部分內(nèi)容總之，滲流基本微分方程包含三部分內(nèi)容: 流動項流動項、源匯項源匯項以及以及累積項累積項 滲流基本微分方程由兩個定律控制：滲流基本微分方程由兩個定律控制： 達西定律與物質(zhì)守恒定律達西定律與物質(zhì)守恒定律壓力場和飽和度場的分布是流體流動的結(jié)果，流壓力場和飽和度場的分布是流體流

10、動的結(jié)果，流動規(guī)律由達西定律來確定，流量由物質(zhì)守恒來確動規(guī)律由達西定律來確定，流量由物質(zhì)守恒來確定，這兩個規(guī)律是數(shù)值模擬的定，這兩個規(guī)律是數(shù)值模擬的關(guān)鍵關(guān)鍵所在。所在。2.基本微分方程基本微分方程- 15 -將單相滲流達西定律推廣到多相滲流中，得到將單相滲流達西定律推廣到多相滲流中，得到三相流動達西公式：三相流動達西公式： 式中，式中，w，o，g水、油、氣諸相下標。水、油、氣諸相下標。 相對滲透率；相對滲透率；3. 三維三相滲流基本微分方程三維三相滲流基本微分方程()()()rwrwwwwwwwrorooooooorgrgggggggK KK KPgZK KK KPgZK KK KPgZ UU

11、UrK- 16 - 對油、水、氣三相，每一相運用物質(zhì)守恒定律對油、水、氣三相，每一相運用物質(zhì)守恒定律，可以建立以下，可以建立以下連續(xù)性方程連續(xù)性方程 3. 三維三相滲流基本微分方程三維三相滲流基本微分方程()()()()()()()()wwwwwoooooooogggggooooggogSqtSqtSSqqttUUUU- 17 - 把達西公式代入連續(xù)性方程，得到把達西公式代入連續(xù)性方程，得到三維三相滲流基本微分方程三維三相滲流基本微分方程 3. 三維三相滲流基本微分方程三維三相滲流基本微分方程()()()()()()()()rwwwwwwwooorooooooorgggroooggogoggg

12、gooK KSqtK KSqtK KK KqqSStt相對滲透率、毛管力是飽和度函數(shù)相對滲透率、毛管力是飽和度函數(shù),而粘度、體積系數(shù)而粘度、體積系數(shù)孔隙度與密度均是壓力的函數(shù)孔隙度與密度均是壓力的函數(shù)- 18 -3. 三維三相滲流基本微分方程三維三相滲流基本微分方程令令 rwwwwoorooooggroooorggggK KTK KTK KTK KT- 19 -上述上述T的物理意義是反映相鄰網(wǎng)格間的流動性，的物理意義是反映相鄰網(wǎng)格間的流動性，又叫傳導(dǎo)率（又叫傳導(dǎo)率（transmissibility) 3. 三維三相滲流基本微分方程三維三相滲流基本微分方程- 20 -3. 三維三相滲流基本微分方

13、程三維三相滲流基本微分方程 oP()()()()()()()()wwwwwoooooooogggggooooggogSTqtSTqtSSTTqqttoSwSoP求解變量求解變量wPgPgS- 21 -3. 三維三相滲流基本微分方程三維三相滲流基本微分方程還要補充一些方程還要補充一些方程飽和度飽和度 ：油水毛管力：油水毛管力：油氣毛管力：油氣毛管力：相對滲透率和毛管力是飽相對滲透率和毛管力是飽和度函數(shù)：和度函數(shù)：cowowPPPcoggoPPP1wogSSS()rwrwwKKS()rgrggKKS(,)rorowgKKSS()cowcowwPPS()cogcoggPPS- 22 -Stone公式

14、：公式：式中，式中， 束縛水飽和度束縛水飽和度 油水兩相流動時油相相對滲透率油水兩相流動時油相相對滲透率 油氣兩相流動時油相相對滲透率油氣兩相流動時油相相對滲透率 3. 三維三相滲流基本微分方程三維三相滲流基本微分方程()() ()()()()()rowwrwwroggrggrorwwrggrowwirKSKSKSKSKKSKSKSwirSrowKrogK- 23 -3. 三維三相滲流基本微分方程三維三相滲流基本微分方程流體流體PVT性質(zhì)性質(zhì)水相：水相： 氣相：氣相： （Z) 油相：隨壓力變化，相態(tài)會發(fā)生變化。油相：隨壓力變化，相態(tài)會發(fā)生變化。 三相狀態(tài)：三相狀態(tài)： 兩相狀態(tài)：兩相狀態(tài)：巖石的

15、高壓物性巖石的高壓物性相滲曲線以及流體相滲曲線以及流體PVT性質(zhì)，以表格形式給出，性質(zhì)，以表格形式給出，稱為靜態(tài)數(shù)據(jù)。稱為靜態(tài)數(shù)據(jù)。wC( )ooooP( )ggooP( )ooP( )wPw( ,)oooosP P( ,)ggoosP P( ,)oosP P( )gPgwBgBgC( )ooCCP( ,)oosCCP P( )rrCC P( )P- 24 -3. 三維三相滲流基本微分方程三維三相滲流基本微分方程實際求解變量實際求解變量初始化輸入數(shù)據(jù)初始化輸入數(shù)據(jù)ljfinit.txtoPwSgS油相流體性質(zhì)與壓力和飽和壓力的關(guān)系油相流體性質(zhì)與壓力和飽和壓力的關(guān)系- 25 -產(chǎn)量項（作源匯項處

16、理）（動態(tài)數(shù)據(jù)）產(chǎn)量項（作源匯項處理）（動態(tài)數(shù)據(jù)）3. 三維三相滲流基本微分方程三維三相滲流基本微分方程()()()()rwwwBHwBHwooroooBHoBHoggroooBHoBHorgggBHgBHgKqPIPPgZZKqPIPPgZZKqPIPPgZZKqPIPPgZZ- 26 -3. 三維三相滲流基本微分方程三維三相滲流基本微分方程其中，其中， 在井基準面深度在井基準面深度 上的井底壓力。上的井底壓力。 節(jié)點等效泄油半徑。節(jié)點等效泄油半徑。 井筒半徑。井筒半徑。X X，Y Y，Z Z 網(wǎng)格節(jié)點網(wǎng)格節(jié)點x x，y y和和z z方向步長方向步長 PI叫采油指數(shù)，叫采油指數(shù)， 通常寫成通

17、常寫成KH,可以直接輸入可以直接輸入2lnewK ZPIrrerwr0.2erxy BHPBHzK Z井?dāng)?shù)據(jù)又叫動態(tài)數(shù)據(jù)井?dāng)?shù)據(jù)又叫動態(tài)數(shù)據(jù)ljfmodel.txt- 27 -4 .初邊值條件初邊值條件油氣藏外邊界條件油氣藏外邊界條件: :第一類邊界第一類邊界 條件條件 P=f(x,y,z,t)第二類邊界第二類邊界 條件條件 第三類邊界條件第三類邊界條件: - 28 -4 .初邊值條件初邊值條件4 外邊界外邊界 封閉邊界封閉邊界:尖滅、斷層、圈閉尖滅、斷層、圈閉 流動邊界：邊、底水流動邊界：邊、底水4 內(nèi)邊界，即內(nèi)邊界，即井的生產(chǎn)條件井的生產(chǎn)條件 定產(chǎn)邊界條件（定產(chǎn)邊界條件（rate const

18、raint):給定產(chǎn)液量、產(chǎn)油量給定產(chǎn)液量、產(chǎn)油量、產(chǎn)氣量、注水量、注氣量、產(chǎn)氣量、注水量、注氣量 定壓邊界條件定壓邊界條件(pressure constraint):給定井底流壓給定井底流壓0),(ntzyxpCtzyxp),(- 29 -4 .初邊值條件初邊值條件油藏初始條件油藏初始條件: :( , , ,0)( , , )wwSx y zSx y z( , , ,0)( , , )ggSx y zSx y z- 30 - 5.網(wǎng)格系統(tǒng)網(wǎng)格系統(tǒng)- 31 -在一般情況下在一般情況下流動方程流動方程求不出解析解，只能用離散求不出解析解，只能用離散化方法求數(shù)值解?；椒ㄇ髷?shù)值解。 離散化方法，首

19、先要把求解區(qū)域按一定的網(wǎng)離散化方法，首先要把求解區(qū)域按一定的網(wǎng)格系統(tǒng)進行剖分。網(wǎng)格剖分方法有很多，像直角格系統(tǒng)進行剖分。網(wǎng)格剖分方法有很多，像直角坐標網(wǎng)格、徑向網(wǎng)格、角點網(wǎng)格、坐標網(wǎng)格、徑向網(wǎng)格、角點網(wǎng)格、PEBI網(wǎng)格等。網(wǎng)格等。目前，商用軟件中大都采用目前，商用軟件中大都采用矩形網(wǎng)格系統(tǒng)矩形網(wǎng)格系統(tǒng)。離散。離散化方法主要有：有限差分化方法主要有：有限差分（Finite Difference)、有有限元限元(Finite Element)以及流線法以及流線法(Stream Line)。比較成熟的是比較成熟的是有限差分方法。有限差分方法。有限差分方法大多有限差分方法大多空間變量采用空間變量采用網(wǎng)

20、格塊中心（網(wǎng)格塊中心（Block Center)差分差分, 時時間變量按間變量按向前差分方法向前差分方法。 5.網(wǎng)格系統(tǒng)網(wǎng)格系統(tǒng)- 32 -5.網(wǎng)格系統(tǒng)網(wǎng)格系統(tǒng):矩形網(wǎng)格矩形網(wǎng)格Producer 2Water InjectorProducer 1- 33 -212121212121jjjiiiyyyxxxjjjjjjiiiiiijjjiiiyyyyyyxxxxxxyyyxxx12112112112121212121- 34 -5.網(wǎng)格系統(tǒng)網(wǎng)格系統(tǒng):徑向網(wǎng)格系統(tǒng)徑向網(wǎng)格系統(tǒng)讀取網(wǎng)格數(shù)據(jù)時，R方向循環(huán)最快，然后是 方向，最后是Z方向徑向網(wǎng)格數(shù)據(jù)讀取規(guī)則- 35 -5.網(wǎng)格系統(tǒng)網(wǎng)格系統(tǒng):角點網(wǎng)格角點

21、網(wǎng)格Corner Point- 36 -塊中心網(wǎng)格與角點網(wǎng)格的比較塊中心網(wǎng)格與角點網(wǎng)格的比較塊中心網(wǎng)格表示的含有斷層的傾斜構(gòu)造塊中心網(wǎng)格表示的含有斷層的傾斜構(gòu)造 角點網(wǎng)格表示的含有斷層的傾斜構(gòu)造角點網(wǎng)格表示的含有斷層的傾斜構(gòu)造- 37 -5.網(wǎng)格系統(tǒng)網(wǎng)格系統(tǒng):塊中心網(wǎng)格的數(shù)據(jù)結(jié)構(gòu)塊中心網(wǎng)格的數(shù)據(jù)結(jié)構(gòu)l模型被分為模型被分為20510個網(wǎng)格區(qū)域。個網(wǎng)格區(qū)域。l網(wǎng)格從（網(wǎng)格從（1，1，1）開始向下傾斜。）開始向下傾斜。l各網(wǎng)格塊在各網(wǎng)格塊在X方向上是方向上是300ft，在，在Y方向上是方向上是1000ft。l層厚度從上至下分別是層厚度從上至下分別是32，22，20，4，32，4，26，26，4，28

22、ft- 38 -5.網(wǎng)格系統(tǒng)網(wǎng)格系統(tǒng):塊中心網(wǎng)格的數(shù)據(jù)結(jié)構(gòu)塊中心網(wǎng)格的數(shù)據(jù)結(jié)構(gòu)-前前20個個TOPS定義網(wǎng)格定義網(wǎng)格 (1, 1, 1) 到到 (20, 1, 1)6855.000 6865.0006875.0006885.0006895.0006905.0006915.0006925.0006935.0006945.0007005.0007015.0007025.0007035.0007045.0007055.0007065.0007075.0007085.0007095.000-下面下面20個個TOPS定義網(wǎng)格定義網(wǎng)格 (1, 2, 1) 到到 (20, 2, 1)6930.0006940

23、.0006950.0006960.0006970.0006980.0006990.0007000.0007010.0007020.0007080.0007090.0007100.0007110.0007120.0007130.0007140.0007150.0007160.0007170.000-下面下面20個個TOPS定義網(wǎng)格定義網(wǎng)格 (1, 3, 1) 到到 (20, 3, 1)7030.0007040.0007050.0007060.0007070.0007080.0007090.0007100.0007110.0007120.0007180.0007190.0007200.000721

24、0.0007220.0007230.0007240.0007250.0007260.0007270.000-下面下面20個個TOPS定義網(wǎng)格定義網(wǎng)格 (1, 4, 1) 到到 (20, 4, 1)7130.0007140.0007150.0007160.0007170.0007180.0007190.0007200.0007210.0007220.0007280.0007290.0007300.0007310.0007320.0007330.0007340.0007350.0007360.0007370.000-下面下面20個個TOPS定義網(wǎng)格定義網(wǎng)格 (1, 5, 1) 到到 (20, 5,

25、 1)7205.0007215.0007225.0007235.0007245.0007255.0007265.0007275.0007285.0007295.0007355.0007365.0007375.0007385.0007395.0007405.0007415.0007425.0007435.0007445.000/ 這完成了第一層這完成了第一層TOPS的定義的定義- 39 -角點網(wǎng)格的數(shù)據(jù)結(jié)構(gòu)角點網(wǎng)格的數(shù)據(jù)結(jié)構(gòu)-COORD-此關(guān)鍵字定義坐標線此關(guān)鍵字定義坐標線-X1Y1Z1X2Y2Z2 0.0.6825.0000.0.7023.000 300.0000 0.6835.000300.

26、00000.7033.000 600.0000 0.6845.000600.00000.7043.000 900.0000 0.6855.000900.00000.7053.000 1200.000 0.6865.0001200.0000.7063.000./對于對于20 * 10 的網(wǎng)格模型的網(wǎng)格模型, 需要需要21*11 條坐標線，即條坐標線，即 231條。條。-每一條坐標線都用每一條坐標線都用6個數(shù)來定義，所以關(guān)鍵字個數(shù)來定義，所以關(guān)鍵字COORD包含有包含有1386個數(shù)。個數(shù)。ZCORN-這個關(guān)鍵字定義網(wǎng)格拐點的深度，數(shù)據(jù)按照這個關(guān)鍵字定義網(wǎng)格拐點的深度，數(shù)據(jù)按照X (或或 R) 方向

27、循環(huán)最快，然后是方向循環(huán)最快，然后是Y (或或 THETA) 方向，最后是方向，最后是Z方向方向的順序排列。的順序排列。6825.0006835.0006835.0006845.0006845.0006855.0006855.0006865.0006865.0006875.0006875.0006885.0006885.0006895.0006895.0006905.0006905.0006915.0006915.0006925.0006975.0006985.0006985.0006995.0006995.0007005.0007005.0007015.0007015.0007025.000

28、7025.0007035.0007035.0007045.0007045.0007055.0007055.0007065.0007065.0007075.000./ 對于對于1000個網(wǎng)格，需要個網(wǎng)格，需要 8000 個個ZCORN值值-這些只是角點網(wǎng)格模型的數(shù)據(jù)體中很少的一部分這些只是角點網(wǎng)格模型的數(shù)據(jù)體中很少的一部分- 40 -塊中心網(wǎng)格與角點網(wǎng)格的比較塊中心網(wǎng)格與角點網(wǎng)格的比較BC（塊中心網(wǎng)格）CP（角點網(wǎng)格）網(wǎng)格描述簡單網(wǎng)格描述復(fù)雜前處理軟件不是必須的必須使用前處理軟件與其它許多模擬器兼容與很少的其它模擬器兼容很難為非常規(guī)構(gòu)造建立模型能夠準確的為非常規(guī)構(gòu)造建立模型網(wǎng)格數(shù)據(jù)量小網(wǎng)格數(shù)據(jù)量

29、大無法區(qū)分傾斜面和斷層區(qū)分傾斜面和斷層很難如實的建立尖滅和剝蝕面的模型如實的建立尖滅和剝蝕面的模型沿著斷層面建立的網(wǎng)格連接不正確，需要手工進行修改。能夠正確的模擬沿著斷層面層之間的接觸關(guān)系容易建立徑向網(wǎng)格模型若沒有前處理軟件，則很難建立徑向網(wǎng)格模型。- 41 -x1ii正常連接的網(wǎng)格關(guān)系正常連接的網(wǎng)格關(guān)系粗細網(wǎng)格交界處的粗細網(wǎng)格交界處的流量計算流量計算非正常連接技術(shù)1, 11/2riii iiik kyzQBx 兩點流量公式：兩點流量公式：- 42 -4 傳統(tǒng)的結(jié)構(gòu)網(wǎng)格傳統(tǒng)的結(jié)構(gòu)網(wǎng)格-即規(guī)則矩形剖即規(guī)則矩形剖分的網(wǎng)格分的網(wǎng)格,具有如下不足具有如下不足:4 存在較嚴重的網(wǎng)格取向效應(yīng)存在較嚴重的網(wǎng)

30、格取向效應(yīng)4 不能保證每口井都在網(wǎng)格中心不能保證每口井都在網(wǎng)格中心4 不能很好地描述油藏形態(tài)不能很好地描述油藏形態(tài)4 粗細網(wǎng)格間的過渡銜接性差粗細網(wǎng)格間的過渡銜接性差4 研究新的網(wǎng)格技術(shù)勢在必行研究新的網(wǎng)格技術(shù)勢在必行.結(jié)構(gòu)網(wǎng)格的不足xy好好差差井井好好開發(fā)軟件- 43 -Gurpinar, 2001Castellini, 2001Wolfsteiner et al., 2002Prevost 2003- 44 -4 常用的非結(jié)構(gòu)網(wǎng)格常用的非結(jié)構(gòu)網(wǎng)格4 各網(wǎng)格中心之間連線構(gòu)成三角各網(wǎng)格中心之間連線構(gòu)成三角網(wǎng)網(wǎng)4 PEBI網(wǎng)格（又叫垂直平分網(wǎng)格網(wǎng)格（又叫垂直平分網(wǎng)格Perpendicular B

31、isectors):各三各三角形邊的垂直平分線構(gòu)成網(wǎng)格角形邊的垂直平分線構(gòu)成網(wǎng)格邊界邊界4 CVFE網(wǎng)格網(wǎng)格-控制體有限元網(wǎng)格控制體有限元網(wǎng)格(Control volume finite element):各三角形邊的中點與重心的連各三角形邊的中點與重心的連線構(gòu)成網(wǎng)格邊界線構(gòu)成網(wǎng)格邊界4 待解決的問題待解決的問題4 非結(jié)構(gòu)網(wǎng)格條件下代數(shù)方程組非結(jié)構(gòu)網(wǎng)格條件下代數(shù)方程組的快速求解的快速求解CVFE網(wǎng)格PEBI網(wǎng)格開發(fā)軟件- 45 -5.網(wǎng)格系統(tǒng)網(wǎng)格系統(tǒng): PEBI網(wǎng)格網(wǎng)格PEBI- 46 -Locally different grids- 47 -New Faulted GridHeavily

32、Faulted ReservoirNew Fault Traces- 48 - 49 -Geological Model of a Heavily Faulted ReservoirSimulation Model- 50 -1. 1. 概述概述 局部網(wǎng)格加密局部網(wǎng)格加密(Local Grid Refinement(Local Grid Refinement，簡稱簡稱LGR)LGR)技術(shù)主技術(shù)主要分：要分：靜態(tài)局部網(wǎng)格加密靜態(tài)局部網(wǎng)格加密(Static Local Grid Refinement)(Static Local Grid Refinement)局部加局部加密區(qū)域在整個計算過程中固定

33、不變的，常用于斷層、尖滅密區(qū)域在整個計算過程中固定不變的，常用于斷層、尖滅、裂縫以及井附近的加密。該技術(shù)已經(jīng)成熟，在多加軟件、裂縫以及井附近的加密。該技術(shù)已經(jīng)成熟，在多加軟件公司的軟件產(chǎn)品中都有這一功能。公司的軟件產(chǎn)品中都有這一功能。動態(tài)局部網(wǎng)格加密動態(tài)局部網(wǎng)格加密(Dynamic Local Grid Refinement)(Dynamic Local Grid Refinement)局部加局部加密區(qū)域在整個計算過程中隨流體的運動位置而變化，主要密區(qū)域在整個計算過程中隨流體的運動位置而變化，主要用于追蹤驅(qū)替前緣。該方法在研究之中，最早是用于追蹤驅(qū)替前緣。該方法在研究之中，最早是Heinema

34、nnHeinemann將動態(tài)局部網(wǎng)格用于水驅(qū)油藏數(shù)值模擬將動態(tài)局部網(wǎng)格用于水驅(qū)油藏數(shù)值模擬 。局部網(wǎng)格加密局部網(wǎng)格加密(Local Grid Refinement)- 51 -Productive Area 1Productive Area 2Sealing InterbeddingAquiferAquifer- 52 -AquiferProductive AreaVertical Refinement Zone- 53 -隨驅(qū)替前緣移動的動態(tài)局部網(wǎng)格加密系統(tǒng)隨驅(qū)替前緣移動的動態(tài)局部網(wǎng)格加密系統(tǒng)- 54 -開發(fā)軟件- 55 -Single Well BlockHybrid Grid Refin

35、ement)Horizontal Well GridComputed Results- 56 -Independent time stepping At any locationSeveral different gridding options are available within one model in the same position- 57 -Window GridBasic GridWindow polygonW indow G rid B asic G rid- 58 -Automatically and implicitly calculated in/out-flux

36、will be used as a von Neumann boundary condition- 59 -Good match required windowing- 60 - 差分方程差分方程 6.差分方程差分方程(Difference Equations)11111111111()()() ()()() ()()()() nnnnwwwwwwwnnooononooooooonnnnggooggoggngnooggooggVTqSStVTqSStTTqqVSSSStxxxyyyzzzTTTT - 61 - 6. 差分方程差分方程111122()()xxxiiiiiiTTT 111122()

37、()yyyjjjjjjTTT 111122()()zzzkkkkkkTTT - 62 - 該差分方程組無論是左端還是右端都是非線該差分方程組無論是左端還是右端都是非線性的，是一組性的，是一組非線性代數(shù)方程組非線性代數(shù)方程組。非線性最強的。非線性最強的參數(shù)是相對滲透率和毛管力。這兩個參數(shù)在時間參數(shù)是相對滲透率和毛管力。這兩個參數(shù)在時間方面的處理不同，對計算的穩(wěn)定性以及計算結(jié)果方面的處理不同，對計算的穩(wěn)定性以及計算結(jié)果影響最大，也是形成形形色色解法的關(guān)鍵所在。影響最大，也是形成形形色色解法的關(guān)鍵所在。 根據(jù)隱式程度的不同，數(shù)值模擬的算法有根據(jù)隱式程度的不同，數(shù)值模擬的算法有IMPES (Impli

38、cit Pressure Explicit Saturation) 方方法、法、半隱式半隱式(Semi-Implicit Method)方法方法和和全隱式全隱式方法方法(Fully Implicit Method )。IMPES：隱式求解壓力顯式求解飽和度。對流度隱式求解壓力顯式求解飽和度。對流度中與時間有關(guān)的量均取上一時間步的值。中與時間有關(guān)的量均取上一時間步的值。6.差分方程差分方程- 63 -IMPES方法的優(yōu)點是速度快，缺點是穩(wěn)定性差方法的優(yōu)點是速度快，缺點是穩(wěn)定性差,尤其是錐進問題以及過泡點問題不能很好解決。尤其是錐進問題以及過泡點問題不能很好解決。因此時間步長不能太大。目前，因此時

39、間步長不能太大。目前，IMPES方法在化方法在化學(xué)驅(qū)模型，尤其是流線法中發(fā)揮著作用。學(xué)驅(qū)模型，尤其是流線法中發(fā)揮著作用。全隱式方法全隱式方法主要采用了主要采用了Newton-Raphson迭代技迭代技術(shù)術(shù), 對方程中各項系數(shù)，包括井系數(shù)全部隱式處對方程中各項系數(shù)，包括井系數(shù)全部隱式處理。全隱式方法較之半隱式方法，具有更好的穩(wěn)理。全隱式方法較之半隱式方法，具有更好的穩(wěn)定性，時間步長更長，徹底解決了過泡點問題，定性，時間步長更長，徹底解決了過泡點問題，是八十年代數(shù)值模擬技術(shù)走向成熟的關(guān)鍵技術(shù)，是八十年代數(shù)值模擬技術(shù)走向成熟的關(guān)鍵技術(shù)，也是大型油藏數(shù)值模擬的主流方法。也是大型油藏數(shù)值模擬的主流方法。

40、 6. 差分方程差分方程- 64 - 全隱式方法和半隱式方法的根本差別在于：全隱式方法和半隱式方法的根本差別在于： 半隱式方法半隱式方法是將非線性方程直接線性化形成是將非線性方程直接線性化形成線性代數(shù)方程組，來求解線性代數(shù)方程組，來求解n+1時刻的值。時刻的值。 全隱式方法全隱式方法采用采用Newton-Raphson方法，迭方法，迭代地求解非線性方程組，得到代地求解非線性方程組，得到n+1時刻的變量。時刻的變量。 用一維非線性方程表示：用一維非線性方程表示： 給定初始近似給定初始近似 對函數(shù)對函數(shù) 在點處在點處 進行一階泰勒展開進行一階泰勒展開, 可近可近似地表示為似地表示為6.差分方程差分

41、方程000f( ) f( ) f( )() 0 xxx x xf( ) 0 x 0 xf( ) 0 x 0 x- 65 -半隱式：半隱式：全隱式方法：全隱式方法：6.差分方程差分方程10()()nnnnf xxxf xx為初值000( )( )f xxxf x- 66 -下面介紹全隱式方法。為了減小計算機舍入誤下面介紹全隱式方法。為了減小計算機舍入誤差的影響，差的影響， 模型中不直接求壓力和飽和度，而模型中不直接求壓力和飽和度，而是求其增量。是求其增量。其中：其中：x代表求解變量代表求解變量 n時間迭代步時間迭代步 牛頓迭代步牛頓迭代步6.差分方程差分方程1nnxxx1llxxx11nllxx

42、xxl- 67 -求解變量：求解變量： ， ，對于二相狀態(tài)：對于二相狀態(tài)： 對于三相狀態(tài)：對于三相狀態(tài)： 6. 差分方程差分方程1111111111111()()()()()()()()()()()()(llllnwwwwwwwwwlolololonooooooooollllglglooggogglgngooggooggoVTqSSStVTqSSStTTqqVSSSSt )oggSSPSwSxSxSoSxPs- 68 -對方程組進行線性化處理。對所有變量進行對方程組進行線性化處理。對所有變量進行Taylor展開，保留線性部分。變量進行運算后展開，保留線性部分。變量進行運算后也保留線性部分，舍去

43、高階項。形成線性代數(shù)也保留線性部分，舍去高階項。形成線性代數(shù)方程組。方程組。7. 差分方程的線性化差分方程的線性化1llTTTTTPSwSxPSwSx1llwPSwSxPSwSx - 69 -全隱式方法形成的線性代數(shù)方程組的結(jié)構(gòu)全隱式方法形成的線性代數(shù)方程組的結(jié)構(gòu) 7. 差分方程的線性化差分方程的線性化xbA矩陣結(jié)構(gòu)示意圖矩陣結(jié)構(gòu)示意圖- 70 - 油藏數(shù)值模擬的線性代數(shù)方程組是大型稀疏矩陣油藏數(shù)值模擬的線性代數(shù)方程組是大型稀疏矩陣，其求解效率在油藏數(shù)值模擬計算中舉足輕重。，其求解效率在油藏數(shù)值模擬計算中舉足輕重。因為，其計算量占整個計算量約因為，其計算量占整個計算量約70%80%。因。因此，

44、研制快速求解大型稀疏矩陣的新方法成為油此，研制快速求解大型稀疏矩陣的新方法成為油藏數(shù)值模擬發(fā)展的一個重要方向。藏數(shù)值模擬發(fā)展的一個重要方向。 線性代數(shù)方程組的求解方法主要有三類：線性代數(shù)方程組的求解方法主要有三類： 1）直接法：）直接法：Gauss消元、消元、D4排序的排序的Gauss消元消元 2）迭代法：）迭代法：SOR(松弛迭代法）松弛迭代法） 3）共軛梯度類法、預(yù)處理共軛梯度類法（）共軛梯度類法、預(yù)處理共軛梯度類法（PCG方法）方法） 8. 線性代數(shù)方程組的求解線性代數(shù)方程組的求解- 71 - 目前最先進的方法是預(yù)處理共軛梯度類方目前最先進的方法是預(yù)處理共軛梯度類方法。法。BO軟件中采用軟件中采用RSVP方法。該方法計算快方法。該方法計算快，精度高。，精度高。8. 線性代數(shù)方程組的求解線性代數(shù)方程組的求解- 72 -黑油模型油藏數(shù)值模擬器流程黑油模型油藏數(shù)值模擬器流程二、黑油模型簡介二、