油藏?cái)?shù)值模擬基礎(chǔ)

油藏?cái)?shù)值模擬根底中國(guó)石油大學(xué)石油工程學(xué)院2021年9月華北油田培訓(xùn)班課程第一章 油藏?cái)?shù)值模擬進(jìn)展油藏?cái)?shù)值模擬的根本概念80年代的油藏?cái)?shù)值模擬進(jìn)展90年代的油藏?cái)?shù)值模擬進(jìn)展

第一節(jié)油藏?cái)?shù)值模擬的根本概念具有同一水動(dòng)力學(xué)系統(tǒng)的油氣聚集體。限于油藏的特殊性：地下，面積大，看不到－－不能直接研究模擬分物理模擬和數(shù)學(xué)模擬。模擬就是用模型來(lái)研究物理過(guò)程，油藏模擬就是用模型來(lái)研究油藏中的流體運(yùn)動(dòng)規(guī)律〔直觀再現(xiàn)〕。物理模擬：根據(jù)相似原理，把自然界中的原型按比例縮小，制成物理模型。然后使原型中的物理過(guò)程按一定的相似關(guān)系在模型中展現(xiàn)。－平板模型、巖心模型不能反映實(shí)際油藏的非均質(zhì)特征，只用來(lái)研究驅(qū)油機(jī)理數(shù)學(xué)模擬

通過(guò)求解某一物理過(guò)程的數(shù)學(xué)方程組來(lái)研究這個(gè)物理過(guò)程變化規(guī)律的方法。一、油藏?cái)?shù)值模擬的定義油藏?cái)?shù)值模擬就是用數(shù)值的方法來(lái)解油藏中流體〔相或者組分〕滲流的偏微分方程組。相：指體系中具有相同成分，相同物理、化學(xué)性質(zhì)的均勻物質(zhì)局部，相與相之間有明顯的界面組分：指混合物體系中的各個(gè)成分。為什么采用數(shù)值求解？

所謂數(shù)值方法是一種近似的解法，即用離散化的方法把連續(xù)函數(shù)轉(zhuǎn)變成離散函數(shù)，用計(jì)算機(jī)來(lái)求解。離散的方法：有限差分法、有限元法、譜分析方法運(yùn)動(dòng)方程、狀態(tài)方程、守恒方程〔質(zhì)量、體積、能量〕優(yōu)點(diǎn)1)能重復(fù)進(jìn)行，可以進(jìn)行所謂的“屢次開(kāi)發(fā)〞。2)可以在短時(shí)間內(nèi)進(jìn)行開(kāi)發(fā)，本錢較低。3)可以模擬各種非均質(zhì)情況及復(fù)雜流體流動(dòng)。油藏?cái)?shù)值模擬的特點(diǎn)缺點(diǎn)1)模擬依賴于油藏描述和生產(chǎn)動(dòng)態(tài)。2)模型本身均有一定的假設(shè)條件，與實(shí)際油藏有一定誤差，簡(jiǎn)化掉某些因素和原理。抓住主要問(wèn)題，滿足現(xiàn)場(chǎng)的需要二、數(shù)值模擬流程

建立計(jì)算機(jī)模型

建立數(shù)值模型

建立數(shù)學(xué)模型不是進(jìn)行具體的數(shù)值模擬工作內(nèi)容1.建立數(shù)學(xué)模型

2.建立數(shù)值模型

通過(guò)離散化，將連續(xù)的偏微分方程組轉(zhuǎn)換成離散的有限差分方程組，再用多種方法將非線性系數(shù)線性化，成為線性代數(shù)方程組，然后求解線性代數(shù)方程組

即建立一套描述油藏中流體滲流的偏微分方程組，包括初、邊值問(wèn)題。偏微分方程組線性代數(shù)方程組得到壓力、飽和度等有限差分方程組離散化線性化解方程組3.建立計(jì)算機(jī)模型將資料〔靜、動(dòng)態(tài)〕輸入，系數(shù)矩陣和常數(shù)項(xiàng)的形成，多種解法和結(jié)果的輸出等，編制成計(jì)算機(jī)程序。數(shù)值模擬的關(guān)鍵是計(jì)算的精度和速度。精度：離散程度、數(shù)值計(jì)算誤差、離散方程穩(wěn)定性速度：計(jì)算機(jī)速度、解法速度、模型準(zhǔn)確性從離散的程度看，精度和速度是矛盾的。三、用途油藏描述油藏動(dòng)態(tài)預(yù)測(cè)驅(qū)油機(jī)理研究1.油藏描述油藏描述是油田開(kāi)發(fā)的根底，是一項(xiàng)系統(tǒng)工程，由多學(xué)科各種方法聯(lián)合研究的結(jié)果。油藏?cái)?shù)值模擬作為一種方法，在油藏描述中起了一定的作用。－不同的方法研究的尺度不同1)孔隙結(jié)構(gòu)研究~10μm級(jí)CT、核磁共振、圖象分析儀、微觀驅(qū)油機(jī)理、毛管壓力實(shí)驗(yàn)2)巖石物理性質(zhì)研究……cm級(jí)利用油層物理實(shí)驗(yàn)室進(jìn)行Ф、K、S、kr等研究。3)測(cè)井研究~m級(jí)利用測(cè)井方法進(jìn)行Ф、K、S等研究。4)井間地震研究1~100m級(jí)利用井間地震進(jìn)行油層連續(xù)性、油層厚度及剩余油分布等研究。5)井間監(jiān)測(cè)m~1000m

利用試井和示蹤劑等方法來(lái)研究主力油層的、k、s等。6)數(shù)值模擬網(wǎng)格

10m~100m

將上述研究的Ф、K、S等參數(shù)填入數(shù)值模擬網(wǎng)格。7)油藏?cái)?shù)值模擬1000m~10000m

利用油藏?cái)?shù)值方法研究油藏中的P、S分布。8)小井距生產(chǎn)試驗(yàn)區(qū)如大慶、孤東等小井距研究彈性驅(qū)、水驅(qū)和各種提高采收率方法從以上不同的研究尺度可以看出，油藏?cái)?shù)值模擬是適合于油藏大尺度范圍內(nèi)研究的方法，對(duì)生產(chǎn)有直接的應(yīng)用。微觀研究：驅(qū)替機(jī)理和巖石孔隙結(jié)構(gòu)的研究－小尺度：油水運(yùn)動(dòng)規(guī)律和油藏物性參數(shù)－均質(zhì)－不反映油藏非均質(zhì)大尺度：整個(gè)油藏滲流規(guī)律和剩余油分布認(rèn)識(shí)2.油藏動(dòng)態(tài)預(yù)測(cè)在油藏開(kāi)發(fā)的不同階段，利用數(shù)值模擬進(jìn)行動(dòng)態(tài)預(yù)測(cè)，有不同的用途規(guī)劃方案開(kāi)發(fā)方案提高開(kāi)發(fā)效果方案調(diào)整方案開(kāi)發(fā)前期開(kāi)發(fā)初期開(kāi)發(fā)中期開(kāi)發(fā)后期3.油藏驅(qū)油機(jī)理研究1)層內(nèi)油水運(yùn)動(dòng)機(jī)理研究油水運(yùn)動(dòng)的作用力：驅(qū)動(dòng)力：重力：毛管力：正韻律反韻律正韻律0L凹0反韻律反韻律L凸正韻律水驅(qū)特征3)化學(xué)驅(qū)對(duì)油水運(yùn)動(dòng)規(guī)律的影響活性劑驅(qū)主要降低油水界面張力，從而提高油相的相對(duì)滲透率。聚合物驅(qū)主要是增加水的粘度，從而降低油水粘度比。堿驅(qū)可與原油中的酸性物質(zhì)反響生成活性物質(zhì)，可使油水界面張力進(jìn)一步降低，與活性劑驅(qū)有同樣效果。混相驅(qū)可使混相劑與原油作用降低界面張力，減少原油粘度。上述機(jī)理均可在數(shù)值模擬中表達(dá)，從而得到油水運(yùn)動(dòng)規(guī)律?；钚詣?qū)聚合物驅(qū)堿驅(qū)混相驅(qū)思考題什么是數(shù)值模擬，有什么特點(diǎn)和用途？目前國(guó)內(nèi)外比較成熟的軟件有哪些？數(shù)學(xué)模型的分類和推導(dǎo)原那么組分模型〔黑油模型和凝析氣藏〕雙重介質(zhì)模型〔黑油〕注蒸汽熱采模型聚合物驅(qū)模型三元復(fù)合驅(qū)模型水平井模型第二章數(shù)學(xué)模型第一節(jié)數(shù)學(xué)模型的分類和推導(dǎo)原那么一、數(shù)學(xué)模型的分類1.按空間維數(shù)來(lái)分零維---物質(zhì)平衡方程一維－巖心水驅(qū)油，注采井間動(dòng)態(tài)二維－三維2.按流體相數(shù)來(lái)分單相－氣藏、油藏彈性開(kāi)發(fā)兩相－氣藏水驅(qū)，油藏水驅(qū)三相3.按流體組分來(lái)分單組分兩組分…N組分4.按巖石類型來(lái)分單重介質(zhì)〔砂巖〕雙重介質(zhì)〔碳酸鹽巖、低滲透油田〕5.按模型功能來(lái)分黑油模型√√凝析氣藏模型√雙重介質(zhì)模型√熱采模型√根本模型聚合物驅(qū)模型三元復(fù)合驅(qū)模型水平井模型……二、數(shù)學(xué)模型的推導(dǎo)原那么1.質(zhì)量守恒方程單位時(shí)間內(nèi)流入單元體的流體質(zhì)量單位時(shí)間內(nèi)流出單元體的流體質(zhì)量單位時(shí)間內(nèi)從單元體注入或采出的流體質(zhì)量單位時(shí)間內(nèi)單元體中的質(zhì)量增加量－+=推導(dǎo)方法：1)微分方法①直角坐標(biāo)單元體法2)積分法2、運(yùn)動(dòng)方程1)牛頓流體單相或多相或τo為流體的屈服應(yīng)力

4)含有啟動(dòng)壓力－低滲透3、狀態(tài)方程1)液體2)巖石巖石壓縮系數(shù)孔隙壓縮系數(shù)3)氣體理想氣體實(shí)際氣體即dpdClrr1=或流體微可壓縮兩端積分4.能量守恒方程單位時(shí)間內(nèi)流入單元體的能量單位時(shí)間內(nèi)流出單元體的能量單位時(shí)間內(nèi)從單元體注入或采出的能量單位時(shí)間內(nèi)單元體中能量的增加－＋＝式中Hj—流體j相的焓，kJ/kgλ—導(dǎo)熱系數(shù)，kJ/(h.

C.m)Uj—流體的內(nèi)能，kJ/kg(ρC)r—油藏的熱容量，kJ/(m3.

C)qh—源、匯項(xiàng)，kJ/(m3.h)qhl—頂、底層的熱損失，kJ/(m3.h)單元體的能量變化有三局部：熱對(duì)流、熱傳導(dǎo)、內(nèi)部熱源4.能量守恒方程對(duì)流項(xiàng)傳導(dǎo)項(xiàng)式中Hj—流體j相的焓，kJ/kgλ—導(dǎo)熱系數(shù)，kJ/(h.

C.m)Uj—流體的內(nèi)能，kJ/kg

(ρC)r—油藏的熱容量，kJ/(m3.

C)qh—源、匯項(xiàng)，kJ/(m3.h)qhl—頂、底層的熱損失，kJ/(m3.h)熱源項(xiàng)流體項(xiàng)巖石項(xiàng)同理，可用微分法或積分法導(dǎo)出5.數(shù)學(xué)模型定解條件完整的數(shù)學(xué)模型：泛定方程(組)：描述油層流體運(yùn)動(dòng)規(guī)律的滲流方程(組)，可以求解，但是解有假設(shè)干個(gè)。定解條件：邊界條件，和初始條件。1.邊界條件

邊界條件：外邊界條件和內(nèi)邊界條件。外邊界條件：指油藏外邊界所處的狀態(tài)；內(nèi)邊界條件：指油水井所處的狀態(tài)。外邊界條件三種形式：定壓外邊界、定流量邊界和混合邊界。

定壓外邊界

油藏的邊界上某一時(shí)刻的壓力是已知的，這種邊界條件又稱為第一類邊界條件或狄利克里條件。例如具有較大天然供水區(qū)的油藏，油水邊界的壓力是一定的.定流量邊界油藏的邊界上有流量流過(guò)，且流量是的，這種邊界條件又稱為第二類邊界條件或紐曼條件。最簡(jiǎn)單且最常見(jiàn)的定流量邊界條件是封閉邊界，油藏邊界上無(wú)流量通過(guò)，例如油藏尖滅邊界、封閉斷層邊界和規(guī)那么注采井網(wǎng)的對(duì)稱線等均可看成封閉邊界條件.混合邊界

混合邊界是指油藏邊界條件是壓力和壓力導(dǎo)數(shù)的線性組合形式，這種邊界條件又稱為第三類邊界條件

對(duì)于具有邊、底水的油藏，油水邊界處的條件在整個(gè)開(kāi)發(fā)過(guò)程中是不斷變化的，既不定壓，也不定流量，稱為混合邊界條件。但的具體形式不易直接表達(dá)出來(lái)，對(duì)于實(shí)際問(wèn)題，這種邊界條件在某一特定的時(shí)間內(nèi)可以近似為第一類或第二類邊界條件，所以混合邊界條件在油藏?cái)?shù)值模擬中比較少見(jiàn)。內(nèi)邊界條件又稱為井點(diǎn)條件，通常有定產(chǎn)量〔或注入量〕條件和定井底流壓條件。⑴定產(chǎn)量(或注入量)條件當(dāng)井的產(chǎn)量(注入量)給定時(shí)，滲流方程中需要加上源匯質(zhì)量項(xiàng)。注水井：定注水量生產(chǎn)井：油相產(chǎn)量、水相產(chǎn)量、總產(chǎn)液量層位：給定的是總的井口產(chǎn)量，各相產(chǎn)量需要內(nèi)部劈產(chǎn)，考慮滲透率和相對(duì)滲透率。⑵定流壓條件由于井眼幾何尺寸遠(yuǎn)遠(yuǎn)小于油藏，因此油井或注水井通常作為點(diǎn)匯或點(diǎn)源，所以油藏?cái)?shù)值模擬中通常把定井底流壓條件轉(zhuǎn)化為源匯項(xiàng)處理.2.初始條件初始條件是指在初始時(shí)刻，油藏內(nèi)壓力分布、飽和度分布或溫度分布是的。(1)初始?jí)毫?/p>

給定油水界面或者油氣界面的深度，根據(jù)流體的密度和高差計(jì)算非平衡初始化平衡初始化假設(shè)深度，位于氣頂區(qū)內(nèi)，那么：假設(shè)深度，位于底水區(qū)內(nèi)，那么假設(shè)深度，位于含油區(qū)內(nèi)，那么：GOCWOCThesecondclass〔2)初始飽和度非平衡初始化－適用于沒(méi)有邊底水的油藏

平衡初始化〔3)初始溫度分布對(duì)于等溫滲流開(kāi)發(fā)過(guò)程，一般不需要給出初始溫度分布，但對(duì)于熱力采油等非等溫滲流過(guò)程，那么需要指定油藏的初始溫度分布。確定初始溫度分布時(shí)需要給出某一基準(zhǔn)深度的溫度，然后根據(jù)溫度梯度計(jì)算出初始溫度分布，即：第二節(jié)組分模型一、組分模型的一般式假設(shè)條件2.組分質(zhì)量守恒方程)符合達(dá)西滲流定律2)等溫滲流3)油、氣、水三相，N個(gè)組分Cig表示氣相中i組分的質(zhì)量分?jǐn)?shù)〔氣相合計(jì)1〕Cio表示油相中i組分的質(zhì)量分?jǐn)?shù)Ciw表示水相中i組分的質(zhì)量分?jǐn)?shù)4)油相和氣相隨著壓力變化而發(fā)生相態(tài)變化5)流體和巖石均微可壓縮6)油藏非均質(zhì)和各向異性7)考慮毛管力和重力3.輔助方程

1)2)3)4)5)6)7)任何兩相之間都有比例常數(shù)4.未知數(shù)和方程數(shù)再加上邊界條件和初始條件組分模型的通用模型，可以根據(jù)油藏條件簡(jiǎn)化假設(shè)條件

1)符合達(dá)西滲流定律

2)等溫滲流

3)油、氣、水三相,甲烷、重?zé)N和水三個(gè)組分。重?zé)N組分完全溶解在油相中，而甲烷組分可以自由氣形式存在于氣相中，也可以溶解氣形式存在于油相中。氣體的溶解和逸出在泡點(diǎn)瞬間完成。水組分存在于水相中,也可以考慮氣組分溶解在水相中。－黑油模型

4)油相和氣相隨壓力而發(fā)生相態(tài)變化

5)流體和巖石均微可壓縮；

6)油藏非均質(zhì)和各向異性；

7)考慮毛管力和重力。2.組分質(zhì)量守恒方程水相僅存在水組分氣相僅存在氣組分地面單位體積油相在地下中可同時(shí)存在氣組分和油組分代入組分的質(zhì)量守恒方程水組分油組分氣組分3.輔助方程

1)2)

3)

4)

5)水的物性參數(shù)可以視為不變必要的方程4.未知數(shù)和方程數(shù)再加上邊界條件和初始條件

由假設(shè)條件知道了每個(gè)組分在各相中的存在分?jǐn)?shù)，所以在這里不需要再進(jìn)行求解，未知數(shù)減少第三節(jié)雙重介質(zhì)的黑油模型一、雙重孔隙介質(zhì)油藏的有關(guān)概念

雙重介質(zhì)油藏是存在天然裂縫的油藏，這種油藏在實(shí)際分析中，常視為由兩種孔隙介質(zhì)組成，即基質(zhì)巖塊介質(zhì)和裂縫介質(zhì)，且兩種介質(zhì)均勻分布，油藏中任何一個(gè)體積單元都存在著這兩個(gè)系統(tǒng)。由于兩種孔隙介質(zhì)具有不同的儲(chǔ)油性和滲透性，因此當(dāng)油井生產(chǎn)時(shí)壓力波的擴(kuò)散和地下流體的滲流規(guī)律將與均質(zhì)油藏完全不同。在雙重介質(zhì)中的任何一點(diǎn)應(yīng)同時(shí)引進(jìn)兩個(gè)壓力(即裂縫中的壓力和基質(zhì)巖塊中的壓力)參數(shù)，同時(shí)也將存在兩個(gè)滲流場(chǎng)。

由于裂縫系統(tǒng)的滲透率比基巖系統(tǒng)的滲透率大得多，認(rèn)為原地下流體由基質(zhì)巖塊到裂縫系統(tǒng)，然后由裂縫系統(tǒng)流到井筒，這種模型稱為雙孔單滲模型，此外還有雙孔雙滲模型。

由于裂縫系統(tǒng)的滲透率比基巖系統(tǒng)的滲透率大得多，認(rèn)為地下流體由基質(zhì)巖塊到裂縫系統(tǒng)，然后由裂縫系統(tǒng)流到井筒，忽略由基質(zhì)巖塊系統(tǒng)直接流入井筒的過(guò)程，即：基巖系統(tǒng)裂縫系統(tǒng)井筒?；|(zhì)巖塊系統(tǒng)裂縫系統(tǒng)基質(zhì)巖塊系統(tǒng)km裂縫系統(tǒng)kf雙孔單滲模型雙孔雙滲模型竄流：擬穩(wěn)態(tài)竄流

假設(shè)基巖系統(tǒng)壓力處處相等，即不考慮基巖內(nèi)部的流體流動(dòng)，竄流只與兩個(gè)系統(tǒng)的壓力有關(guān)。大多數(shù)竄流屬于擬穩(wěn)態(tài)竄流。不穩(wěn)態(tài)竄流假設(shè)基巖系統(tǒng)壓力處處不相等，即考慮基巖內(nèi)部的流體流動(dòng)。一般情況下竄流在早期屬于不穩(wěn)態(tài)竄流。一、假設(shè)條件除巖石性質(zhì)外，其余條件與組分模型的黑油模型相同。巖石性質(zhì)如下：1〕巖石孔隙結(jié)構(gòu)由基巖和裂縫組成。2〕基巖是儲(chǔ)油的主要空間，裂縫是流油的主要通道。假設(shè)基巖不發(fā)生流動(dòng)，稱為雙孔單滲介質(zhì)；假設(shè)基巖發(fā)生流動(dòng)，稱為雙孔雙滲介質(zhì)；3〕基巖和裂縫之間產(chǎn)生質(zhì)量交換，為擬穩(wěn)態(tài)竄流。二、雙孔雙滲黑油模型1.組分質(zhì)量守恒方程-分裂縫和基質(zhì)兩個(gè)系統(tǒng)列方程1)裂縫系統(tǒng)水組分油組分氣組分2)基質(zhì)系統(tǒng)水組分油組分氣組分2.輔助方程3.未知數(shù)和方程數(shù)1.組分質(zhì)量守恒方程1)裂縫系統(tǒng)同雙孔雙滲黑油模型。2)基質(zhì)系統(tǒng)根據(jù)定義，基質(zhì)系統(tǒng)中沒(méi)有滲流，只有向裂縫的竄流。油組分水組分氣組分2.輔助方程同雙孔雙滲黑油模型。3.未知數(shù)和方程數(shù)同雙孔雙滲黑油模型。作業(yè)完整的數(shù)值模擬方程組包括哪些局部？其中都涉及到哪些方程？黑油模型的假設(shè)條件是什么？考慮的是什么情況的流動(dòng)？雙重介質(zhì)模型與單重介質(zhì)有哪些不同？

第三章油藏?cái)?shù)值模擬計(jì)算方法

主要內(nèi)容2.1差分格式2.2網(wǎng)格排列方法2.3差分方程組的求解方法2.4預(yù)處理共扼梯度求解方法第1節(jié)差分格式1.1離散網(wǎng)格系統(tǒng)

數(shù)值模擬中離散包括：

時(shí)間區(qū)域的離散，空間區(qū)域的離散產(chǎn)量時(shí)間物性參數(shù)非線性→線性第一類邊界條件〔定壓〕－－點(diǎn)中心網(wǎng)格第二、三類邊界條件－－－塊中心網(wǎng)格實(shí)際油藏的邊界比較復(fù)雜，一般都將邊界處理成封閉的邊界，所以采用塊中心網(wǎng)格?！虨榱搜芯繂?wèn)題的方便，特做如下約定：離散空間坐標(biāo)用作標(biāo)記，空間步長(zhǎng)為，網(wǎng)格交界處用作標(biāo)記；離散時(shí)間坐標(biāo)用n作標(biāo)記，時(shí)間步長(zhǎng)為；空間坐標(biāo)標(biāo)記為下標(biāo)，時(shí)間坐標(biāo)標(biāo)記為上標(biāo)，例如表示在一維空間上i離散點(diǎn)n時(shí)間步的參數(shù)值，表示在一維空間上i與i+1離散網(wǎng)格交界處n+1時(shí)間步的參數(shù)值。1.2差分格式

1-i

i

1+i

xxD-

x

xxD+

根據(jù)Taylor展開(kāi)公式為一階后差分為一階前差分一階中心差分兩式相加為二階差分兩式相減常系數(shù)偏微分方程的差分格式－古典顯式差分格式

離散網(wǎng)格點(diǎn)〔i，n)建立差分，且右端對(duì)時(shí)間的差分取一階前差分令，整理上式可得

只要給出物理問(wèn)題的初始條件和邊界條件，直接由上式就可以計(jì)算出各離散點(diǎn)的具體數(shù)值，因此稱為顯式差分格式，其截?cái)嗾`差為。常系數(shù)偏微分方程的差分格式－古典隱式差分格式

在離散網(wǎng)格點(diǎn)〔i,n+1)建立差分，且右端對(duì)時(shí)間的差分取一階后差分假設(shè)物理問(wèn)題的初始條件和邊界條件，實(shí)際上是以離散點(diǎn)參數(shù)值為未知量的封閉代數(shù)方程組，對(duì)于點(diǎn)中心網(wǎng)格，當(dāng)時(shí)n＝0，定壓外邊界的方程組形式為

求解方程組可以得到各離散點(diǎn)的具體數(shù)值，各離散點(diǎn)的數(shù)值是聯(lián)立計(jì)算出來(lái)的，因此稱為隱式差分格式，其截?cái)嗾`差為當(dāng)時(shí)n＝0，封閉外邊界的方程組形式為：求解方程組可以得到各離散點(diǎn)的具體數(shù)值，各離散點(diǎn)的數(shù)值是聯(lián)立計(jì)算出來(lái)的，因此稱為隱式差分格式，其截?cái)嗾`差為未知未知變系數(shù)偏微分方程的差分格式

由于油層巖石的非均質(zhì)性、油藏流體的多相流動(dòng)特征、變滲流截面，油藏?cái)?shù)學(xué)方程的二階偏微分通常是變系數(shù)的，例如非均質(zhì)油層單相微可壓縮流體變截面滲流的數(shù)學(xué)方程為

由于變系數(shù)二階偏微分方程所反映的物理問(wèn)題為非均勻體系，建立差分方程時(shí)的離散空間和離散時(shí)間都可能是非均勻的，通常將二階偏微分轉(zhuǎn)化為一階偏微分的方法進(jìn)行處理。令對(duì)于三維網(wǎng)格系統(tǒng)有每個(gè)網(wǎng)格都可以得到一個(gè)方程1.3差分格式的適定性

由于差分方程是偏微分方程的微分利用差分替代后得到的，這種近似必然帶來(lái)一定的誤差，所以微分方程是適定的，而差分方程未必是適定的。差分方程的適定取決于：誤差的大小，差分格式的相容性、收斂性和穩(wěn)定性第二節(jié)網(wǎng)格排列方法微分方程采用差分近似以后形成了差分方程，只要給出物理問(wèn)題的初始條件和邊界條件，所有離散點(diǎn)(網(wǎng)格)對(duì)應(yīng)的隱式差分格式組合起來(lái)形成一個(gè)代數(shù)方程組，方程組系數(shù)矩陣具有如下特點(diǎn)：⑴空間離散步長(zhǎng)不能太大，因此離散格點(diǎn)較多，而每一個(gè)離散格點(diǎn)對(duì)應(yīng)一個(gè)未知數(shù)，因此系數(shù)矩陣是高階的；⑵某一離散點(diǎn)的差分格式中，只出現(xiàn)本離散格點(diǎn)和周圍離散格點(diǎn)未知數(shù)，與其余大局部離散格點(diǎn)無(wú)關(guān)，矩陣的每一行上都存在大量零元素，所以系數(shù)矩陣是稀疏的；⑶系數(shù)矩陣中非零元素出現(xiàn)的形式與空間網(wǎng)格的離散未知數(shù)排列順序有關(guān)。一、一維離散網(wǎng)格排序012345

014235

(a)

(b)

對(duì)于〔a)排列格式3

-1

0

0

-1

3

-1

0

0

-1

3

-1

0

0

-1

3

三對(duì)角端點(diǎn)0,5對(duì)于〔b)排列格式012345

014235

(a)

(b)

兩種排列格式的計(jì)算結(jié)果是相同的計(jì)算速度和計(jì)算工作量可能不同目的：系數(shù)矩陣的結(jié)構(gòu)便于計(jì)算，且節(jié)省計(jì)算存儲(chǔ)量。

二、二維離散網(wǎng)格排序⒈標(biāo)準(zhǔn)排列格式

標(biāo)準(zhǔn)排列格式是由離散網(wǎng)格的自然下標(biāo)按照遞增順序和循環(huán)形成的一維未知數(shù)順序。⑴行標(biāo)準(zhǔn)排列格式行標(biāo)準(zhǔn)列標(biāo)準(zhǔn)行標(biāo)準(zhǔn)排列行標(biāo)準(zhǔn)排列⑵列標(biāo)準(zhǔn)排列格式兩種排列的方式都形成5對(duì)角矩陣，但是形狀不同，非零系數(shù)的集中程度是不相同，導(dǎo)致需要的存儲(chǔ)量和計(jì)算速度不同。集中程度？按小的方向首先排列，帶寬小。⒉D2對(duì)角排列格式⒊A3點(diǎn)交替排列格式

⒋D4排列格式該排列格式是D2和A3組合而成的排列格式

三、三維離散網(wǎng)格排序⒈標(biāo)準(zhǔn)排列格式對(duì)于三維滲流問(wèn)題，標(biāo)準(zhǔn)排列格式的系數(shù)矩陣非零元素為七對(duì)角，假設(shè)三維滲流區(qū)域網(wǎng)格剖分后三個(gè)方向的網(wǎng)格數(shù)分別為，那么首先以縱向?qū)?yīng)的方向進(jìn)行排列，所得到系數(shù)矩陣的帶寬最小。除剖面地質(zhì)模型之外，油藏?cái)?shù)值模擬中一般總是最小。三、三維離散網(wǎng)格排序⒉D4排列格式對(duì)于三維滲流問(wèn)題，根據(jù)以前離散點(diǎn)標(biāo)號(hào)的約定，排列格式的原那么如下：⑴三維網(wǎng)格下標(biāo)之和為奇數(shù)者，按奇數(shù)遞增順序優(yōu)先排序；然后再按偶數(shù)排列；？？？原因？三、三維離散網(wǎng)格排序⒉D4排列格式對(duì)于三維滲流問(wèn)題，根據(jù)以前離散點(diǎn)標(biāo)號(hào)的約定，排列格式的原那么如下：⑴三維網(wǎng)格下標(biāo)之和為奇數(shù)者，按奇數(shù)遞增順序優(yōu)先排序；然后再按偶數(shù)排列；⑵三維網(wǎng)格下標(biāo)之和相同者，按K值遞減順序排序，即K值大者優(yōu)先排序；⑶三維網(wǎng)格下標(biāo)之和相同者，且K值相同，按J值遞減順序排序，即J值大者優(yōu)先排序。如何排列？1,1,1＝32,1,1＝41,2,1＝42,2,1＝51,1,2＝42,1,2＝51,2,2＝52,2,2＝61,1,3＝52,1,3＝61,2,3＝62,2,3＝7①②③④⑤⑥⑦⑧⑨⑩(11)(12)①②③④⑤⑥⑦⑧⑨⑩(11)(12)四、局部網(wǎng)格加密與死節(jié)點(diǎn)消除排序⒈局部網(wǎng)格加密與死節(jié)點(diǎn)數(shù)值模擬中網(wǎng)格的劃分?jǐn)?shù)量(或節(jié)點(diǎn)數(shù))與數(shù)值解的精度有著直接的關(guān)系，為提高解的精度，希望網(wǎng)格劃分得愈多愈好，但從另一方面來(lái)講，網(wǎng)格數(shù)量又至少以一次冪的關(guān)系影響著計(jì)算機(jī)的內(nèi)存量和計(jì)算時(shí)間。而在注采井較密的部位和物性變化較快的部位需要密網(wǎng)格，其他部位較少的網(wǎng)格即滿足計(jì)算精度。－－網(wǎng)格加密另一方面，在劃分網(wǎng)格的矩形區(qū)域內(nèi)，不可防止地要框入—些無(wú)效格塊，將這類格塊視為死節(jié)點(diǎn)。在網(wǎng)格排序時(shí)避開(kāi)這些死節(jié)點(diǎn)，可減少節(jié)點(diǎn)數(shù)目。⒉排序方法和系數(shù)矩陣的生成局部網(wǎng)格加密是一個(gè)在原劃分網(wǎng)格根底上增加節(jié)點(diǎn)數(shù)目的方法，而死節(jié)點(diǎn)的消除是一個(gè)減少節(jié)點(diǎn)數(shù)目的方法，兩種方法的目的各不相同，但它們的原理相同。隨著加密網(wǎng)格個(gè)數(shù)不同、加密點(diǎn)和死節(jié)點(diǎn)出現(xiàn)的位置不同，其系數(shù)矩陣的階數(shù)和結(jié)構(gòu)也要發(fā)生相應(yīng)的改變。對(duì)于網(wǎng)格局部加密和死節(jié)點(diǎn)消除通常采用自然網(wǎng)格排序方法。⑴局部加密當(dāng)排到需進(jìn)行局部加密的格塊后，仍按著原來(lái)的序號(hào)，用順序排序法對(duì)加密網(wǎng)格進(jìn)行排序，排完后，保持這個(gè)序號(hào)繼續(xù)向下排。⑵死節(jié)點(diǎn)消除當(dāng)排到死節(jié)點(diǎn)處時(shí)，將該格塊跳過(guò)去，直到遇到某一格塊不為死節(jié)點(diǎn)為止，然后按原序號(hào)繼續(xù)排下去。斷層

網(wǎng)格的排列方向和順序決定了矩陣的形式；矩陣的形式?jīng)Q定了求解過(guò)程中存儲(chǔ)量和求解速度不同系數(shù)矩陣形式可以采用不同的求解方式。第三節(jié)差分方程組的求解方法采用隱式差分格式建立的數(shù)值模型為一個(gè)方程組對(duì)于常系數(shù)偏微分方程，差分方程組是線性的對(duì)于變系數(shù)微分方程，差分方程組可能是非線性的，這主要取決于傳導(dǎo)系數(shù)中非線性項(xiàng)的處理方法。一、線性方程組的求解方法⒈三對(duì)角系數(shù)矩陣方程求解方法LU分解法⑴根據(jù)矩陣乘法規(guī)那么計(jì)算元素由的第一行的元素得

根據(jù)矩陣乘法規(guī)那么，計(jì)算其它元素整理以上三式得

⑵計(jì)算未知變量由計(jì)算過(guò)渡矩陣Y

由計(jì)算X適用于一維問(wèn)題的求解作業(yè)1，有限差分的格式：向前、向后、中間差分。2誤差產(chǎn)生的原因，受什么因素影響？3數(shù)值模擬方程系數(shù)矩陣的特點(diǎn)，產(chǎn)生原因？4矩陣的排列格式對(duì)系數(shù)矩陣有什么影響？5線性方程組的解法？6非線性方程組的解法？一維油水兩相水驅(qū)油的數(shù)值模擬方法一維徑向單相流的數(shù)值模擬方法第四章一維油藏的數(shù)值模擬方法用于研究均勻井網(wǎng)和單井的動(dòng)態(tài)第一節(jié)一維油水兩相水驅(qū)油的數(shù)值模擬方法一、數(shù)學(xué)模型1.假設(shè)條件1)符合達(dá)西滲流定律2)等溫滲流3)油、水兩相及油水兩組分4)一維流動(dòng)5)流體和巖石不可壓縮；6)油藏巖石性質(zhì)〔k,〕沿一維非均質(zhì)7)不考慮毛管力和重力2.組分質(zhì)量守恒方程由組分質(zhì)量守恒方程的一般式逐步簡(jiǎn)化到上述假設(shè)條件。1)三維油氣水三相N個(gè)組分2)三維油水兩相兩組分3)一維油水兩相兩組分4)假設(shè)不考慮流體和巖石壓縮性令，將質(zhì)量流量轉(zhuǎn)變?yōu)轶w積流量

不考慮重力假設(shè)不考慮油水兩相之間的毛管力－油相、水相壓力相同水組分油組分3.輔助方程4.未知數(shù)和方程數(shù)(1)(2)(3)0L水油＋水5.初始條件和邊界條件

在巖心中飽和油和束縛水，然后在左端注入水，右端先出油，后出油和水，要求巖心中各點(diǎn)壓力、飽和度隨時(shí)間的變化。I.CB.C二、差分方程組的建立1.預(yù)備知識(shí)1)方程解法問(wèn)題

順序求解(Sequential)，先求P，再求S

隱式壓力顯式飽和度IMPES

即ImplicitPressureExplicitSaturation

隱式壓力隱式飽和度IMPIMS

即ImplicitPressureImplicitSaturation

所謂隱式，即用一個(gè)線性代數(shù)方程組求解一組未知函數(shù)；所謂顯式，即用一個(gè)線性代數(shù)計(jì)算求解一個(gè)未知數(shù)。

聯(lián)立求解(Simultaneous)，P、S同時(shí)求解

半隱式(Semi-Implicit)

全隱式(Fully-Implicit)2)方程非線性系數(shù)的顯式和隱式處理

Krl隨Sl而變化，而Sl又是未知函數(shù)，隨時(shí)間變化，因此Kr有以下幾種處理：顯式，即Sl取n時(shí)刻，為值。半隱式，即krl隨時(shí)間而變化，用Taylor級(jí)數(shù)展開(kāi)，取前二項(xiàng)。隱式，即Sl取n+1

時(shí)刻，因Sln+1為未知，采用迭代法，常用的為Newton－Raphson方法。簡(jiǎn)單計(jì)算采用顯式的處理方法3)方程非線性系數(shù)項(xiàng)取上游權(quán)i-1ii+1x流動(dòng)由i-1到i流動(dòng)由i到i－1誤差比較大上游權(quán)的意義在于把流動(dòng)滯后由時(shí)間上〔采用顯式處理導(dǎo)致n+1－n〕所造成的誤差轉(zhuǎn)移到空間上進(jìn)行補(bǔ)償，因此上游權(quán)處理方法是穩(wěn)定的。隱式壓力顯式飽和度IMPES流動(dòng)由i-1到i流動(dòng)由i到i－1顯式處理2、隱式求壓力的方程－IMPES為消除SL項(xiàng)，〔1〕＋〔2〕可得〔4〕對(duì)〔4〕式進(jìn)行離散化，采用塊中心網(wǎng)格，且網(wǎng)格大小相等，均為x采用二階隱式差分格式后得：i=123n-1n

x〔5〕i=123n-1n

x〔5〕式可以分以下三種情況來(lái)討論①對(duì)于第2至n-1個(gè)網(wǎng)格，因無(wú)注入和采出qv=0,其中系數(shù)項(xiàng)采用顯式處理，并采用上游權(quán)原那么：〔6)〔7〕i=123n-1n

x②對(duì)于第一個(gè)網(wǎng)格i=1,注入為qv,(5)式中第二項(xiàng)，由于沒(méi)有流體從0流到1網(wǎng)格，因此無(wú)此項(xiàng)，第一項(xiàng)系數(shù)采用上游權(quán)顯示處理。變形，并令(8)i=123n-1n

x③對(duì)于第n個(gè)網(wǎng)格，i=n，產(chǎn)出為qv，(5)式中第一項(xiàng)，由于沒(méi)有流體從n流到n+1網(wǎng)格，因此無(wú)此項(xiàng)，第二項(xiàng)系數(shù)采用上游權(quán)顯式處理：〔6〕、〔7〕、〔8〕構(gòu)成了i從1到n的線性代數(shù)方程組。總的方程式00矩陣方程如下：系數(shù)矩陣為三對(duì)角矩陣，可用LU分解法求解i=1,2,,n壓力值?！?…………………..…..3.顯式求飽和度方程方程〔1〕采用二階隱式差分后，系數(shù)項(xiàng)采用顯式處理和上游權(quán)原那么，可得〔9〕式中Pn＋1已從隱式壓力方法求解得到，而求Sw可分以下三種情況：a:2,3,…n-1b:1c:n(9)1〕i=2，3，，n-1,qv=0(9)式可寫(xiě)為：

〔10〕2〕i=1,〔9〕式中只有第一項(xiàng)，可得：兩端乘以，并令Qv=qvA△x那么〔11〕只注水初始含水飽和度增加快，然后上升速度逐漸減慢3)i=n，〔9〕式中只有第二項(xiàng)，可得：兩端乘以Ax，并令n-1n產(chǎn)水量，出口端為？見(jiàn)水后？油水界面沒(méi)有到達(dá)井底前，Qw為0，井底含水飽和度不變，水的流度項(xiàng)為0含水率利用〔10〕、〔11〕、〔12〕式即可求得i=1,2,,n的的值?！?2〕如何求解含油飽和度？?jī)煞N方法：方法1：根據(jù)飽和度方程求解方法2：根據(jù)油相方程，按水相飽和度的方法求解同樣，對(duì)于1,2～n-1，n采用不同的方程第1個(gè)網(wǎng)格：qov＝？？？？為0，因?yàn)槭亲⑺?~n網(wǎng)格：qov＝？沒(méi)有產(chǎn)量項(xiàng)為0第n個(gè)網(wǎng)格：qov＝？產(chǎn)油項(xiàng)三、計(jì)算框圖讀入SwcPi不同Sw下的KrwKro值〔線性差值〕計(jì)算λwλoλ值計(jì)算壓力方程的系數(shù)矩陣和右邊項(xiàng)c(I),a(I),b(I),d(I)用分解法求P(I)

顯式求Sw(I)So(I)T<Tmax打印TP(I)Sw(I)SO(I)

結(jié)束程序開(kāi)始讀入μw,μo,K,φ,L,△x

讀入S(I),Krw(I),Kro(I),I=1……M讀入Tmax,△t,QvT=0YNT+t

程序一形成求解相對(duì)滲透率的子程序或函數(shù)在一維水驅(qū)油巖心中：φ=0.3;Swc=0.2;Pi=0MPa;A=10cm2;L=1m;Qv=0.1cm3/s;Pn=0MPa油水相對(duì)滲透率數(shù)據(jù)如下表：Tmax=500s△t=10s取塊中心網(wǎng)格N=40△x=2.5cm要求：編制一維水驅(qū)油數(shù)值模擬程序。1)打印100、200、300、400、500s時(shí)的壓力和飽和度隨巖心長(zhǎng)度上的分布，并繪圖表示。2)求出水突破時(shí)間。第二節(jié)一維徑向單相流的數(shù)值模擬方法一、數(shù)學(xué)模型1.假設(shè)條件1)符合達(dá)西滲流定律2)等溫滲流3)單相流體4)一維徑向向井流動(dòng)5)巖石不可壓縮，流體微可壓縮6)油藏巖石性質(zhì)〔k,〕沿徑向不發(fā)生變化7)不考慮重力2.質(zhì)量守恒方程由三維柱坐標(biāo)質(zhì)量守恒方程逐步簡(jiǎn)化到上述假設(shè)條件1)三維單相非均質(zhì)油藏可壓縮流體和巖石2)三維單相均質(zhì)油藏可壓縮流體和巖石

k=常數(shù)＝常數(shù)，不考慮粘度變化，單相3)一維單相均質(zhì)油藏可壓縮流體和巖石

〔1〕〔2〕〔3〕2.質(zhì)量守恒方程4)一維單相均質(zhì)油藏可壓縮流體，不可壓縮巖石

〔4〕〔5〕右端相：左端相：按微可壓縮處理3、初始條件和邊界條件假設(shè)園形邊界中心一口井，單相流體向井底流動(dòng)，求在各種內(nèi)外邊界條件下的壓力分布。rerwrerw定壓封閉初始條件1)外邊界2)內(nèi)邊界定產(chǎn)定流壓可求以下問(wèn)題：1.定壓外邊界條件下:1)內(nèi)邊界定產(chǎn)，求不同時(shí)間沿徑向壓力分布和井底流壓。

2)內(nèi)邊界定流壓，求不同時(shí)間沿徑向壓力分布和產(chǎn)量。2.封閉外邊界條件下

1)內(nèi)邊界定產(chǎn)，求不同時(shí)間沿徑向壓力分布和井底流壓。

2)內(nèi)邊界定流壓，求不同時(shí)間沿徑向壓力分布和產(chǎn)量。二、差分方程組的建立1.預(yù)備知識(shí)—不等距徑向網(wǎng)格的建立必要性當(dāng)單相流體向井底流動(dòng)時(shí)，由于外邊界附近，流動(dòng)截面大，流速小，壓力變化也小，因此網(wǎng)格尺寸要大；而在井底附近，流動(dòng)截面小，流速大，壓力變化也大，因此網(wǎng)格尺寸要小。因此沿徑向可以采用不等距網(wǎng)格。等比級(jí)數(shù)變化令…將不等距網(wǎng)格r坐標(biāo)轉(zhuǎn)換成等距網(wǎng)格x坐標(biāo)，即將等比級(jí)數(shù)取對(duì)數(shù)。…xxxrwr1r2r3rnx0x