解剖黑油模型源程序-final

1、黑油模型解剖匯報(bào)小組成員：周詩雨、徐晨、李碩軒、徐錫勇、倪威。BLACK OILBLACK OILBLACK OILBLACK OILBLACK OIL黑油模型基本結(jié)構(gòu)(1)GRID1(2)PARM1(3)TABLE(4)UNIT1(5)CODES(6)NOTES(7)TRAN1(8)ORATE(9)SOLMAT(10&11)CGIN&MF1(12)MATBAL(13)PRTPS(14)INTERP(15)INTPVT主程序主程序:從結(jié)構(gòu)上代表了整個(gè)模型，通過調(diào)用不同的子程序，完成全部運(yùn)算1.定義運(yùn)算所需數(shù)組；2.需要調(diào)用和應(yīng)生成的文件的打開和關(guān)閉；3.通過調(diào)用子程序給模型賦基礎(chǔ)數(shù)據(jù)和初始數(shù)據(jù)

2、；4.通過調(diào)用子程序給模型的運(yùn)行確定各種控制；5.在運(yùn)行過程中反復(fù)讀入井點(diǎn)數(shù)據(jù)（包括產(chǎn)量和注入量、井底流壓、流動(dòng)指數(shù)等）以及打印輸出控制碼；6.分層油氣水地質(zhì)儲(chǔ)量計(jì)算；7.進(jìn)行井點(diǎn)產(chǎn)量項(xiàng)處理和形成壓力矩陣；8.通過調(diào)用子程序求解壓力方程；9.進(jìn)行顯式飽和度計(jì)算；10.進(jìn)行過泡點(diǎn)處理；11.根據(jù)壓力和飽和度增量控制，自動(dòng)調(diào)節(jié)計(jì)算時(shí)步長(zhǎng)；12.在每一運(yùn)算時(shí)步末進(jìn)行物質(zhì)平衡檢驗(yàn)；13.打印油井、水井、氣井的分井分層報(bào)告和分井報(bào)告；14.通過調(diào)用子程序進(jìn)行計(jì)算結(jié)果打印輸出和形成文件；15.再啟動(dòng)運(yùn)行方式的選擇和再啟動(dòng)文件的生成和調(diào)用；16.運(yùn)算終止的控制；17.運(yùn)行中途錯(cuò)誤信息的提示。主要功能 （1）

3、、調(diào)用子程序，完成整個(gè)黑油模型計(jì)算； （2）、儲(chǔ)量計(jì)算和運(yùn)行控制； （3）、打印輸出控制。總結(jié):由調(diào)用子程序和自帶部分功能共同完成主程序程序框圖基礎(chǔ)數(shù)據(jù)準(zhǔn)備計(jì)算儲(chǔ)量及相關(guān)參數(shù)壓力求解時(shí)步調(diào)整主程序程序框圖計(jì)算飽和度及相關(guān)參數(shù)計(jì)算分層和總產(chǎn)量指標(biāo)過泡點(diǎn)處理，變泡點(diǎn)計(jì)算 輸出；循環(huán)主程序自帶功能（五部分）一、運(yùn)行參數(shù)的輸入 模型在運(yùn)行過程中，對(duì)每一時(shí)步都要往復(fù)地進(jìn)行運(yùn)算，因此也就要往復(fù)地輸入運(yùn)行參數(shù)，這些數(shù)據(jù)（包括其中的提示符）是：IWLCHG改變工作制度控制碼（1或0）；ICHANG本運(yùn)算階段的時(shí)步數(shù)（0的整數(shù)）；IWLREP油井報(bào)告打印控制碼（1或0）；ISUMRY總報(bào)告打印控制碼（1或0）；

4、IPMAP網(wǎng)格壓力打印控制碼（1或0）；ISOMAP網(wǎng)格含油飽和度打印控制碼（1或0）；ISWMAP網(wǎng)格含水飽和度打印控制碼（1或0）；ISGMAP網(wǎng)格含氣飽和度打印控制碼（1或0）；IPBMAP網(wǎng)格飽和壓力打印控制碼（1或0）。(2)時(shí)步長(zhǎng)DAY給定時(shí)步長(zhǎng)（天）DTMIN最小時(shí)步長(zhǎng)（天）DTMAX最大時(shí)步長(zhǎng)（天）工作制度和打印控制代碼、與時(shí)間有關(guān)參數(shù)主程序自帶功能二、原始地質(zhì)儲(chǔ)量計(jì)算 在主程序中的第一時(shí)步內(nèi)要進(jìn)行原始的油水氣（包括游離氣和溶解氣）的儲(chǔ)量計(jì)算，程序自動(dòng)計(jì)算出分層的油水氣儲(chǔ)量，也計(jì)算總的油水氣儲(chǔ)量。首先是計(jì)算出各網(wǎng)格點(diǎn)的孔隙體積VP，然后根據(jù)各網(wǎng)格點(diǎn)壓力值內(nèi)插求其體積系數(shù)BO、B

5、W和BG，然后計(jì)算各層的油水氣（分別按游離氣和溶解氣計(jì)算）：OOIP原油原始地質(zhì)儲(chǔ)量（地面方）；OWIP地層水原始地質(zhì)儲(chǔ)量（地面方）；ODGIP溶解氣原始地質(zhì)儲(chǔ)量（地面方）；OFGIP游離氣原始地質(zhì)儲(chǔ)量（地面方）。 根據(jù)分層的油氣水地質(zhì)儲(chǔ)量再計(jì)算全油藏地質(zhì)儲(chǔ)量，即TOOIP全油藏原油原始地質(zhì)儲(chǔ)量（地面方）；TOWIP全油藏地層水原始地質(zhì)儲(chǔ)量（地面方）；TDGIP全油藏溶解氣原始地質(zhì)儲(chǔ)量（地面方）；TFGIP全油藏游離氣原始地質(zhì)儲(chǔ)量（地面方）。主程序自帶功能三、計(jì)算時(shí)步長(zhǎng)的自動(dòng)調(diào)節(jié) 時(shí)間步長(zhǎng)的自動(dòng)調(diào)節(jié)是為了按照用戶要求和題目本身的運(yùn)行情況，避免由于時(shí)間步長(zhǎng)過大、引起壓力值增量和飽和度增量過大，從

6、而造成誤差過大，也避免由于步長(zhǎng)取得過小，從而造成不必要的計(jì)算時(shí)間過長(zhǎng)的情況?？傊菫榱吮ＷC在必要計(jì)算精度的同時(shí)，以最大的時(shí)步長(zhǎng)進(jìn)行模型的運(yùn)行。PPM網(wǎng)格壓力變化最大值；SOM網(wǎng)格含油飽和度變化最大值；SWM網(wǎng)格含水飽和度變化最大值；SGM網(wǎng)格含氣飽和度變化最大值；DSMC網(wǎng)格飽和度變化最大值；DPMC網(wǎng)格壓力變化最大絕對(duì)值；FACT1時(shí)間步長(zhǎng)增長(zhǎng)系數(shù)；FACT2時(shí)間步長(zhǎng)縮短系數(shù)。根據(jù)需要調(diào)整合適的時(shí)間步長(zhǎng)，避免過大或過小主程序自帶功能四、單井報(bào)告的輸出和數(shù)據(jù)存儲(chǔ) 每一計(jì)算階段末的分井分層的生產(chǎn)數(shù)據(jù)是通過主程序進(jìn)行輸出或存盤。通過油井報(bào)告在屏幕或打印機(jī)輸出的數(shù)據(jù)有：FT運(yùn)算階段末的模擬生產(chǎn)時(shí)間（

7、天）；WELLD井名；IQN1井點(diǎn)X坐標(biāo)；IQN2井點(diǎn)Y坐標(biāo)；K層號(hào)；PWEC計(jì)算的分井分層流動(dòng)壓力值（MPa）；PWF給定的分井分層流動(dòng)壓力值（MPa）；PID分層流動(dòng)指數(shù)；QOO分層日生油量（m3/d）；QGG分層日產(chǎn)氣量（ m3/d ）；QWW分層日產(chǎn)水量（ m3/d ）；GOR分層生產(chǎn)油氣比（m3/m3）；WOR分層生產(chǎn)水油比（m3/m3）；CUMO分層累積產(chǎn)油量（m3）；CUMG分層累積產(chǎn)氣量（m3）；CUMW分層累積產(chǎn)水量（m3）；TOR全油藏日產(chǎn)油量（ m3/d ）；TGR全油藏日產(chǎn)氣量（ m3/d ）；TWR全油藏日產(chǎn)水量（ m3/d ）；TOC全油藏累積產(chǎn)油量（m3）；TGC

8、全油藏累積氣量（m3）；TWC全油藏累積水量（m3）；主程序自帶功能四、單井報(bào)告的輸出和數(shù)據(jù)存儲(chǔ) 通過再啟動(dòng)文件存入盤中的分井分層累積產(chǎn)量，包括J井號(hào)K層號(hào)CUMO累積產(chǎn)油量；CUMG累積產(chǎn)氣量；CUMW累積產(chǎn)水量；五、黑油模型運(yùn)行結(jié)束方式 黑油模型運(yùn)行結(jié)束方式有下列幾種：1.最大時(shí)步數(shù)NMAX2.最大模擬時(shí)間TMAX3.井?dāng)?shù)據(jù)卡片讀完4.最大水油比WORMAX5.最高油氣比GORMAX6.最大平均壓力PAMAX7.最低平均壓力PAMIN一、Grid1子程序參數(shù)說明II，JJ，KK：分別為X，Y，Z方向網(wǎng)格數(shù)； KDX，KDY，KDZ：分別為DX，DY，DZ輸入方式控制代碼；DX，DY，DZ：

9、分別為網(wǎng)格點(diǎn)X，Y，Z方向尺寸；NUMDX，NUMDY，NUMDZ：分別為網(wǎng)格尺寸DX，DY， DZ修正系數(shù)；IDCODE：網(wǎng)格尺寸修正值打印控制代碼；KEL：地層深度計(jì)算控制代碼；ELEV：到網(wǎng)格頂部的深度值。功能介紹 根據(jù)油藏均質(zhì)程度，幾何形狀及井網(wǎng)規(guī)則程度和用戶要求：進(jìn)行三維網(wǎng)格劃分輸入網(wǎng)格尺寸數(shù)據(jù)并修正計(jì)算各網(wǎng)格中部深度網(wǎng)格信息Grid1流程圖Y,Z方向尺寸的讀入類似X方向參數(shù)說明KPH:控制孔隙度數(shù)據(jù)輸入的代碼；KKX:控制X方向滲透率數(shù)據(jù)輸入的代碼；KKY :控制Y方向滲透率數(shù)據(jù)輸入的代碼；KKZ :控制Z方向滲透率數(shù)據(jù)輸入的代碼；NUMP :需要修正孔隙度值的網(wǎng)格數(shù)目；NUMKX

10、 :需要修正X方向滲透率值的網(wǎng)格數(shù)目；NUMKY :需要修正Y方向滲透率值的網(wǎng)格數(shù)目；NUMKZ :需要修正Z方向滲透率值的網(wǎng)格數(shù)目。功能介紹 在油氣藏網(wǎng)格大小已經(jīng)確定的情況下，根據(jù)油層均質(zhì)程度，把已確定出的各網(wǎng)格點(diǎn)的孔隙度和滲透率值，賦予模型各網(wǎng)格點(diǎn)。PARM1子程序孔隙度、滲透率PARM1流程圖參數(shù)說明SAT(I):相對(duì)滲透率數(shù)據(jù)表中的飽和度值；KRO(I) :相對(duì)滲透率數(shù)據(jù)表中的油相相對(duì)滲透率值；PBO:原始飽和壓力；VSLOPE:在高于飽和壓力下，地層原油粘度對(duì)壓力的斜率；RSLOPE:在高于飽和壓力下，地層原油溶解油氣比對(duì)壓力的斜率；POT(I) :原油高壓物性數(shù)據(jù)表中的壓力值；MU

11、OT(I) :原油高壓物性數(shù)據(jù)表中的粘度值；BOT(I) :原油高壓物性數(shù)據(jù)表中的體積系數(shù)值。功能介紹 將數(shù)據(jù)文件中有關(guān)巖石性質(zhì)（相對(duì)滲透率和壓縮吸水），地層流體數(shù)據(jù)（高壓物性、油氣水密度）等以一定規(guī)則輸入。TABLE子程序相滲、PVTTABLE流程圖RSOT(I) :原油高壓物性數(shù)據(jù)表中的溶解氣油比值；MSAT:相對(duì)滲透率表中的飽和度個(gè)數(shù)；MPOT:原油高壓物性表中的壓力值個(gè)數(shù)；RHOSCO:標(biāo)準(zhǔn)狀況下原油密度；RHOSCW:標(biāo)準(zhǔn)狀況下水的密度；RHOSCG:標(biāo)準(zhǔn)狀況下天然氣密度；BOPT(I) :不同壓力下原油體積系數(shù)的斜率。參數(shù)說明（續(xù)）UNIT1：給模型賦初值開始KRESTKPI讀入P

12、I,PGOC,WOC,GOC計(jì)算BTP,BBO,RSO計(jì)算PI,BBW,RSW計(jì)算PGOC,BBG計(jì)算PN(I,J,K)SO=SO1，SW=SW1,SG=SG1讀入SO,SW,SG=1-SO-SW結(jié)束重啟動(dòng)讀入OOIP,OWIP,ODGIP,OFGIP讀入NVQN1活動(dòng)井?dāng)?shù)讀入第I口井在第K層的產(chǎn)油量CUMO，產(chǎn)水量CUMW，產(chǎn)氣量CUMG讀入重啟運(yùn)算起始時(shí)間ETIR,再啟動(dòng)起始步數(shù)NRESTT,再啟動(dòng)初始平均壓PAVGO讀入P(I,J,K), SO(I,J,K),SW(I,J,K),SG(I,J,K),PBOT(I,J,K)初始啟動(dòng)：主要以不同方式給各個(gè)網(wǎng)格點(diǎn)賦初始?jí)毫惋柡投取PI=0，

13、原始?jí)毫Ω鶕?jù)油水界面壓力、高度差和流體密度折算得到； KPI=1，則各網(wǎng)格點(diǎn)原始?jí)毫τ蓴?shù)據(jù)文件直接讀入。KSI=0，則各層油氣水均為定值； KSI=1，則各網(wǎng)格點(diǎn)原始油氣水飽和度不同，由數(shù)據(jù)文件直接讀入。再啟動(dòng)方式：飽和度和壓力值由再啟動(dòng)文件讀入，還應(yīng)讀入再啟動(dòng)起始點(diǎn)的儲(chǔ)量、累計(jì)產(chǎn)量，起始時(shí)刻及步數(shù)等=1=0=0=0=1WOCGOCEL(I,J,K)KSI初始數(shù)據(jù)的輸入分為：初始啟動(dòng)和再啟動(dòng)結(jié)束CODES：輸入運(yùn)行控制變量開始讀入KSN1,KSM1,KCO1,KTR,KCOFF輸出NN,FACT1,FACT2,TMAX,WORMAX,GORMAX,PAMIN,PAMAX,MITER,KNK,O

14、MEGA,TOL,TOL1,DSMAX,DPMAXPAMIN=9.869*PAMINPAMAX=9.869*PAMAXDPMAX=9.869*DPMAXKSN=KSN1KSM=KSM1KCO=KCO1結(jié)束FACT1：在時(shí)間步長(zhǎng)自動(dòng)控制的情況下，增大時(shí)間步長(zhǎng)的參數(shù)；FACT2：在時(shí)間步長(zhǎng)自動(dòng)控制的情況下，縮小時(shí)間步長(zhǎng)的參數(shù)；TMAX：最大模擬時(shí)間；KNK：共軛梯度法中，正交化次數(shù)的限定；MITER：每一時(shí)步內(nèi)共軛梯度法的最大迭代次數(shù)；OMEGA：加速收斂因子；TOL：在迭代過程中，兩次迭代之間最大的壓力變化值；TOL1：求解過程中，改變松弛因子的參數(shù)；DSMAX：最大飽和度變化允許值；DPMAX

15、：最大壓力變化允許值；WORMAX：最大水油比；GORMAX：最大油氣比；PAMIN：最低油層平均壓力；PAMAX：最高油層平均壓力；參數(shù)說明KIP=1：定產(chǎn)量生產(chǎn)井KIP=2：定產(chǎn)量注水井KIP=3：定產(chǎn)量注氣井KIP=-1：定壓力生產(chǎn)井KIP=-2：定壓力注水井KIP=-3：定壓力注氣井KIP=-11：隱式定壓力生產(chǎn)井KIP=-12：隱式定壓力注水井KIP=-13：隱式定壓力注氣井NVQN總井?dāng)?shù)NINJ注水井?dāng)?shù)NPROD生產(chǎn)井?dāng)?shù)WELLID井號(hào)IQN1井所在X方向網(wǎng)格號(hào)IQN2井所在Y方向網(wǎng)格號(hào)IQN3射孔井段頂部Z方向網(wǎng)格號(hào)LAYER射開層數(shù)KIP井型QVO井產(chǎn)油量QVW井產(chǎn)水量QVG井

16、產(chǎn)氣量PID分井分層流動(dòng)指數(shù)PWF分井分層流動(dòng)壓力NODES-讀入井點(diǎn)工作制度的運(yùn)行數(shù)據(jù)，包括井型、井坐標(biāo)、生產(chǎn)數(shù)據(jù)等。START讀入井號(hào)、網(wǎng)格號(hào)、井標(biāo)示、生產(chǎn)數(shù)據(jù)、流動(dòng)指數(shù)判別井類型（KIP值），輸出相應(yīng)生產(chǎn)信息ENDTRAN1-計(jì)算有限差分方程中流動(dòng)項(xiàng)各方向的傳導(dǎo)系數(shù)網(wǎng)格I-1和網(wǎng)格I之間在X方向的傳導(dǎo)系數(shù)網(wǎng)格J-1和網(wǎng)格J之間在Y方向的傳導(dǎo)系數(shù)網(wǎng)格K-1和網(wǎng)格K之間在Z方向的傳導(dǎo)系數(shù)START計(jì)算各個(gè)網(wǎng)格的滲透率修正因子FACX，F(xiàn)ACY，F(xiàn)ACZ計(jì)算各個(gè)網(wǎng)格的滲透率KX, KY, KZ計(jì)算各個(gè)網(wǎng)格傳導(dǎo)系數(shù)TX，TY，TZ讀入要修正傳導(dǎo)率的網(wǎng)格參數(shù)（NUMT*）和修正值打印和重啟動(dòng)代碼判

17、斷輸出修正值方式END調(diào)和平均NUMTX ：修正X方向傳導(dǎo)系數(shù)的網(wǎng)格數(shù)NUMTY ：修正Y方向傳導(dǎo)系數(shù)的網(wǎng)格數(shù)NUMTZ：修正Z方向傳導(dǎo)系數(shù)的網(wǎng)格數(shù)ITCODE ：打印控制代碼KREST ：重啟動(dòng)控制代碼流程圖根據(jù)不同的生產(chǎn)條件（定流量、定壓-顯式、定壓-隱式），根據(jù)流動(dòng)系數(shù)確定分井分層油氣水產(chǎn)量及全井流動(dòng)壓力和全井總產(chǎn)量ORATE-計(jì)算各井分層分層及全井生產(chǎn)指標(biāo)STARTNRATE=0NRATE=1讀入壓力數(shù)據(jù)，計(jì)算油、氣、水流度讀入壓力數(shù)據(jù)，計(jì)算CIP值定油生產(chǎn)井定液生產(chǎn)井定注水量井定注氣量井油井定流量定流壓水井氣井定流壓油井水井氣井END計(jì)算公式見P194-196顯式隱式SOLMAT求解

18、（I,J,K），(I+1,J,K) (I-1,J,K)，(I,J+1,K) (I,J-1,K)， (I,J,K+1) (I,J,K-1)七點(diǎn)處的參數(shù)值開始計(jì)算參數(shù)I=I+1,I=IIJ=J+1,J=JJK=K+1,KLz=y-z Ux=z-x被誰調(diào)用：CGIN算法描述：x=x0;dx=100;N=1;while dx=EFC & NKNK r=b-Ax;q=M-1r; do k=0,KNK-1 ak=(r,M-1r)/(q,Aq) x=x+aq; r=r-aAq; dx=ak*q; if(dxKNK) exit. b=(r,M-1r)/(r,M-1r) q=M-1r+bq N=N+1參數(shù)說明P

19、AVGO:時(shí)步初油藏平均壓力 WOR :水油比PAVG :時(shí)步末油藏平均壓力 GOR :氣油比OP :產(chǎn)油量 CWOR :累積水油比WP :產(chǎn)水量 功能介紹 對(duì)每一步末的計(jì)算結(jié)果進(jìn)行物質(zhì)平衡檢驗(yàn)，確定計(jì)算誤差。這一計(jì)算分別對(duì)油、氣、水進(jìn)行。MATBAL子程序CGOR :累積氣油比WI :注水量 MBEO :原油物質(zhì)平衡誤差百分?jǐn)?shù)GP :產(chǎn)氣量 MBEW :水的物質(zhì)平衡誤差百分?jǐn)?shù)GI :注氣量MBEG :氣的物質(zhì)平衡誤差百分?jǐn)?shù)MATBAL流程圖COP-累積產(chǎn)油量OPR:全油藏日產(chǎn)油量CWP-累計(jì)產(chǎn)水量；CGP-累積產(chǎn)氣量；CWI-累積注水量；CGI-累積注氣量；參數(shù)說明（續(xù)）計(jì)算值誤差計(jì)算ETN

20、:總生產(chǎn)時(shí)間（天）； NVQN:總井?dāng)?shù)（口）NPRD:總生產(chǎn)井?dāng)?shù)（口； NINJ:總注入井?dāng)?shù)（口）PAVG:平均油藏壓力（兆帕）； OPR:全油藏日產(chǎn)油量（方/日）GPR:全油藏日產(chǎn)氣量（方/日）； WPR:全油藏日產(chǎn)水量（方/日）WIR:全油藏日注水量（方/日）； GORM:全油藏瞬時(shí)生產(chǎn)油氣比（米3/米3）FW:全油藏綜合含水率（%）； COP:全油藏累積產(chǎn)油量（方）CWP:全油藏累積產(chǎn)水量（方） CGP:全油藏累積產(chǎn)氣量（方）CWI:全油藏累積注水量（方） VPR:折算采油速度（%）COP/TOOIP:采出程度（%）再啟動(dòng)運(yùn)行所需的初始數(shù)據(jù)的生成、包括：ETN:再啟動(dòng)起始時(shí)刻（天）；NLM1:再啟動(dòng)時(shí)步數(shù)（1的整數(shù)）；PAVG:再啟動(dòng)初始時(shí)刻油藏平均壓力（兆帕）；CGP:再啟動(dòng)初始時(shí)刻油藏累積產(chǎn)氣量（米3）；CGI:再啟動(dòng)初始時(shí)刻油藏累積注氣量（米3）；CWI:再啟動(dòng)初始時(shí)刻油藏累積注水量（米3）；P(I,J,K):再啟動(dòng)初始時(shí)刻油藏內(nèi)各點(diǎn)壓力分布；SO(I,J,K):再啟動(dòng)初始時(shí)刻油藏內(nèi)各點(diǎn)含油飽和度分布；SW(I,J,K):再啟動(dòng)初始時(shí)刻油藏內(nèi)各點(diǎn)含水飽和度變化；SG(I,J,K):再啟動(dòng)初始時(shí)刻油藏內(nèi)各點(diǎn)含氣飽和度分布；PBOT(I,J,K):再啟動(dòng)初始時(shí)刻油藏內(nèi)