相關(guān)系數(shù)控制的致密油藏裂縫隨機(jī)表征和生成方法*

趙玉云 姜漢橋 李俊鍵 成寶洋 范 楨

(1. 浙江理工大學(xué)理學(xué)院數(shù)學(xué)科學(xué)系 浙江杭州 310018； 2. 中國石油大學(xué)(北京)油氣資源與探測國家重點(diǎn)實(shí)驗(yàn)室 北京 102249)

致密油是近年來發(fā)展起來的一種非常規(guī)油氣資源，在我國具有較豐富的儲量，合理有效地開采致密油可以帶來巨大的經(jīng)濟(jì)效益[1-6]。致密油藏中天然裂縫較發(fā)育[7]，從天然露頭和巖心樣品的分析表明裂縫各參數(shù)服從一定的統(tǒng)計學(xué)規(guī)律[8]。同時天然露頭和地質(zhì)統(tǒng)計學(xué)的分析表明，裂縫開度和長度等參數(shù)之間存在一定的相關(guān)性。目前裂縫的隨機(jī)表征方法僅僅考慮了裂縫各參數(shù)的隨機(jī)性，而忽略了這些參數(shù)之間的相關(guān)性[8]。這樣在做致密油藏數(shù)值模擬或者復(fù)雜縫網(wǎng)反演的時候，隨機(jī)生成的裂縫參數(shù)僅僅服從一定的統(tǒng)計分布，而不完全符合實(shí)際情況。

因此，本文引入相關(guān)系數(shù)對各參數(shù)進(jìn)行控制，建立了基于相關(guān)系數(shù)控制的裂縫隨機(jī)表征方法，將裂縫參數(shù)的每一種隨機(jī)分布都轉(zhuǎn)化為均勻分布來生成，用半解析的方法得到了目標(biāo)相關(guān)系數(shù)和初始相關(guān)系數(shù)之間的關(guān)系，從而給出了相關(guān)系數(shù)控制下裂縫各參數(shù)的隨機(jī)生成方法。

1 常用隨機(jī)表征方法及存在問題

致密油藏天然露頭和巖心樣品的分析表明其天然裂縫較發(fā)育，而且裂縫各參數(shù)服從一定的統(tǒng)計學(xué)規(guī)律。

從地球物理的角度來看，裂縫是地層內(nèi)聚力喪失的表現(xiàn)，可以看作是地層破裂的結(jié)果[9]，換言之，是地層不連續(xù)、某些物性參數(shù)發(fā)生突變的地方。從幾何角度來看，二維空間中的裂縫可以簡化為具有有限長度的線段。三維空間中的裂縫可以簡化為具有有限面積的多邊形或橢圓形。

對于二維空間里面的單條裂縫，可以用中心點(diǎn)坐標(biāo)(xc,yc)、長度(L)、方位角(θ)和開度(w)來表征(圖 1)。

圖1 單條裂縫的幾何參數(shù)表征Fig .1 Characterization of geometric parameters of a single fracture

即一條裂縫可以用(xc,yc,θ,L,w)5個參數(shù)唯一確定。同時也可以用裂縫的2個端點(diǎn)坐標(biāo)和開度來表征一條裂縫，即(x1,y1,x2,y2,w)，可見無論何種表征方法，參數(shù)的數(shù)量都是5。

根據(jù)地質(zhì)統(tǒng)計學(xué)，絕大多數(shù)儲層的裂縫各參數(shù)服從的統(tǒng)計學(xué)規(guī)律如表1[8]所示。

表1 裂縫各參數(shù)及其服從的統(tǒng)計分布[8]Table 1 Various parameters of fractures and their statistical distribution[8]

對于一組天然裂縫，傳統(tǒng)的裂縫隨機(jī)表征方法可以表示為如下5組隨機(jī)參數(shù)(axc,bxc)、(ayc,byc)、(μθ,σθ)、(μL,σL)、(μw,σw)。其中，(a,b)為均勻分布的最小值和最大值，(μ,σ)為正態(tài)分布、指數(shù)分布、對數(shù)正態(tài)分布的均值和方差，對于長度的(μL,σL)，當(dāng)分布為均勻分布時，也分別指最小值和最大值。對于具體的一組天然裂縫，其參數(shù)為(xc,yc,θ,L，w)。

但是，對于真實(shí)的裂縫，各參數(shù)之間存在一定的相關(guān)性，比如長度和開度，一般來說具有正相關(guān)性，即裂縫越長，其開度也越大。傳統(tǒng)的裂縫隨機(jī)表征方法忽略了這種參數(shù)之間的相關(guān)性。如表 2顯示的一組裂縫參數(shù)，雖然長度和開度均服從對數(shù)正態(tài)分布，但是對于編號為1的裂縫，其長度最長，而開度卻最小，不太符合實(shí)際情況。

因此傳統(tǒng)的裂縫隨機(jī)表征方法存在一定的缺陷，需要進(jìn)一步完善。

表2 無相關(guān)系數(shù)控制下隨機(jī)生成的裂縫參數(shù)Table 2 Parameters of the fractures randomly generated without the control of correlation coefficient

2 相關(guān)系數(shù)控制的裂縫隨機(jī)表征和生成方法

2.1 相關(guān)系數(shù)控制的裂縫隨機(jī)表征方法

在統(tǒng)計學(xué)中，相關(guān)系數(shù)也叫皮爾遜積-矩相關(guān)系數(shù)，用來度量2個隨機(jī)變量或2組統(tǒng)計數(shù)據(jù)之間的相關(guān)性，其值在-1和1之間。其定義如下:

(1)

或

(2)

若ρXY=1則隨機(jī)變量X和Y呈完全正相關(guān)，具體表現(xiàn)為X越大則Y越大，X越小則Y越小。若ρXY=-1則隨機(jī)變量X和Y呈完全負(fù)相關(guān)，表現(xiàn)為X越大則Y越小，X越小則Y越大。類似的，若r=1則統(tǒng)計數(shù)據(jù)X和Y呈完全正相關(guān)，若r=-1則統(tǒng)計數(shù)據(jù)X和Y呈完全負(fù)相關(guān)。

對于裂縫的5個參數(shù)，其兩兩之間的相關(guān)系數(shù)分別如表 3所示。

表3 裂縫5個參數(shù)之間的相關(guān)系數(shù)Table 3 Correlation coefficients between 5 parameters of fractures

基于相關(guān)系數(shù)的控制，一組天然裂縫可以用如下參數(shù)表征:axc、bxc、ayc、byc、μθ、σθ、μL、σL、μw、σw、ρ(xc，yc)、ρ(xc，θ)、ρ(xc，L)、ρ(xc，w)、ρ(yc，θ)、ρ(yc，L)、ρ(yc，w)、ρ(θ，L)、ρ(θ，w)、ρ(L,w)。

這種表征方法比傳統(tǒng)的表征方法多了相關(guān)系數(shù)的參數(shù)，可以反映前面各參數(shù)之間的相關(guān)性。

2.2 基于統(tǒng)計分布變換的裂縫隨機(jī)生成方法

由于目前統(tǒng)計得到的裂縫各參數(shù)服從的分布只有均勻分布、正態(tài)分布、指數(shù)分布和對數(shù)正態(tài)分布4種分布，因此只需研究這4種分布的隨機(jī)數(shù)的生成即可。

在統(tǒng)計學(xué)中，關(guān)聯(lián)結(jié)構(gòu)是處理隨機(jī)變量相關(guān)性問題的一種方法。通過關(guān)聯(lián)結(jié)構(gòu)來確定2個隨機(jī)變量聯(lián)合分布方法的思想為:將每個邊緣分布都轉(zhuǎn)換為平均分布，即由均勻分布的隨機(jī)數(shù)可以通過一定的變換法則得到任意分布的隨機(jī)數(shù)，具體的變換法則如下：

x=F-1(u)～F

(3)

式(3)中：u是服從均勻分布的隨機(jī)數(shù)，F(xiàn)是想要得到的隨機(jī)數(shù)服從的分布函數(shù)，則由變換法則x=F-1(u)得到的隨機(jī)數(shù)x服從分布F。因此為了得到服從分布F和G的隨機(jī)數(shù)，只需要生成2組服從均勻分布的隨機(jī)數(shù)u1和u2，再做變換x=F-1(u1)和y=G-1(u2)即可。但是，這樣得到的隨機(jī)數(shù)x和y是沒有相關(guān)系數(shù)控制的。為了得到相關(guān)系數(shù)為ρobj的分別服從分布F和G的隨機(jī)數(shù)，初始生成的服從均勻分布的隨機(jī)數(shù)u1和u2應(yīng)該滿足一定的相關(guān)系數(shù)ρu。

事實(shí)上，要得到相關(guān)系數(shù)為ρu的隨機(jī)數(shù)并不是一件容易的事情，更多的時候需要借助服從正態(tài)分布的隨機(jī)數(shù)來實(shí)現(xiàn)。而由正態(tài)分布的隨機(jī)數(shù)得到均勻分布隨機(jī)數(shù)的變換法則恰好是變換法則(式(3))的逆變換，即

u=Φ(z)～U

(4)

其中：z是服從標(biāo)準(zhǔn)正態(tài)分布的隨機(jī)數(shù)；Φ是標(biāo)準(zhǔn)正態(tài)分布的分布函數(shù)

(5)

因此基于關(guān)聯(lián)結(jié)構(gòu)的裂縫隨機(jī)生成方法包含如下3個步驟：①生成相關(guān)系數(shù)為ρ的正態(tài)分布隨機(jī)數(shù)z1和z2。②通過式(4)得到相關(guān)系數(shù)為ρu的服從均勻分布的隨機(jī)數(shù)u1和u2。③通過式(3)得到相關(guān)系數(shù)為ρobj的分別服從分布F和G的隨機(jī)數(shù)x和y。

根據(jù)關(guān)聯(lián)結(jié)構(gòu)理論中的相關(guān)結(jié)論[10-11]，當(dāng)ρobj和ρmin、ρmax之間存在線性關(guān)系時，ρ和ρobj之間存在如下關(guān)系：

(6)

(7)

式(6)～(7)中：ρmin和ρmax分別為ρobj可能取到的最小值和最大值。本文通過半解析的方法確定了ρmin和ρmax，從而得到了ρ和ρobj的關(guān)系。

1) 均勻分布和正態(tài)分布。

設(shè)服從標(biāo)準(zhǔn)正態(tài)分布，記作Z～N(0,1)，則-Z～N(0,1)，且Φ(Z)服從[0,1]之間的均勻分布，記作Φ(Z)～U[0,1]。則有

(8)

(9)

易知

ρ(a+(b-a)Φ(Z),μ-σZ)=

(10)

ρ(a+(b-a)Φ(Z),μ+σZ)=

(11)

所以只需要求ρ(Φ(Z),Z)的值。直接用解析的方法求ρ(Φ(Z),Z)比較困難，分別生成了200組服從標(biāo)準(zhǔn)正態(tài)分布的隨機(jī)數(shù)，每一組有2 000個值，可以求得200個ρ(Φ(Z),Z)，如圖2所示，其均值為0.977。

因此對于任意的均勻分布和正態(tài)分布，有

(12)

2) 均勻分布和均勻分布。

易知對于任意均勻分布有ρmin=-1，ρmax=1，所以對于任意均勻分布

(13)

3) 均勻分布和指數(shù)分布。

設(shè)u服從標(biāo)準(zhǔn)均勻分布，即u～U[0,1]，則1-u～U[0,1]，X=F-1(u)服從指數(shù)分布，其中，F(xiàn)是參數(shù)為λ的指數(shù)分布的分布函數(shù)，其定義為

F(x)=1-e-λx

(14)

圖2 隨機(jī)生成200個ρ(Φ(Z),Z)的統(tǒng)計直方圖Fig .2 Statistical histogram of 200 ρ(Φ(Z),Z) randomly generatedF-1(x)=-ln(1-x)

(15)

對于任意的均勻分布隨機(jī)變量Y=a+(b-a)u～U[a,b]和任意的指數(shù)分布隨機(jī)變量，易知

ρmin=ρ(a+(b-a)u,F-1(1-u))

(16)

ρmax=ρ(a+(b-a)u,F-1(u))

(17)

其中，F(xiàn)0為參數(shù)λ=1的指數(shù)分布函數(shù)，即

F0(x)=1-e-x

(18)

因此對于任意的均勻分布和指數(shù)分布，有

(19)

圖3 隨機(jī)生成200個和的統(tǒng)計直方圖Fig .3 Statistical histogram of 200 randomly generated

4) 均勻分布和對數(shù)正態(tài)分布。

設(shè)服從標(biāo)準(zhǔn)對數(shù)正態(tài)分布，記作Z～LN(0,1)，根據(jù)標(biāo)準(zhǔn)正態(tài)分布的定義有l(wèi)nZ～N(0,1)，所以-lnZ～N(0,1)，1/Z～LN(0,1)，ΦL(Z)和ΦL(1/Z)均服從[0,1]之間的均勻分布，ΦL為標(biāo)準(zhǔn)對數(shù)正態(tài)分布的分布函數(shù)，其定義為

(20)

(21)

(22)

同樣采用數(shù)值模擬的方法，分別生成了200組服從標(biāo)準(zhǔn)對數(shù)正態(tài)分布的隨機(jī)數(shù)，每一組有2 000個值，可以分別求得200個ρ(ΦL(1/Z),Z)和ρ(ΦL(Z),Z)，如圖4所示，其均值分別為-0.70和0.70。

(23)

圖4 隨機(jī)生成200個ρ(ΦL(1/Z),Z)和ρ(ΦL(Z),Z)的統(tǒng)計直方圖Fig .4 Statistical histogram of 200 ρ(ΦL(1/Z),Z) and ρ(ΦL(Z),Z) randomly generated

5) 正態(tài)分布和指數(shù)分布。

設(shè)Z服從標(biāo)準(zhǔn)正態(tài)分布，即Z～N(0,1)，則Φ(Z)～U[0,1]，X=F-1(Φ(Z))服從指數(shù)分布，其中，F(xiàn)是參數(shù)為λ的指數(shù)分布的分布函數(shù)，其定義見式(12)。

對于任意的正態(tài)分布隨機(jī)變量X=μ+σZ～N(μ,σ2)和任意的指數(shù)分布隨機(jī)變量，易知

ρmin=ρ(μ+σZ,F-1(Φ(-Z)))=

(24)

ρmax=ρ(μ+σZ,F-1(Φ(Z)))=

(25)

其中，F(xiàn)0的定義見式(18)。

因此對于任意的正態(tài)分布和指數(shù)分布，有

(26)

6) 正態(tài)分布和對數(shù)正態(tài)分布。

(27)

(28)

圖5 隨機(jī)生成200個和的統(tǒng)計直方圖Fig .5 Statistical histogram of 200 randomly generated

(29)

7) 指數(shù)分布和對數(shù)正態(tài)分布。

公式(6)對于指數(shù)分布和對數(shù)正態(tài)分布的情況不再適用，本文采用多項(xiàng)式直接對ρobj和ρ進(jìn)行擬合，擬合結(jié)果如圖7a所示，擬合公式為

(30)

8) 對數(shù)正態(tài)分布和對數(shù)正態(tài)分布。

公式(7)對于對數(shù)正態(tài)分布和對數(shù)正態(tài)分布的情況也不適用，本文采用多項(xiàng)式直接對ρobj和ρ進(jìn)行

圖6 隨機(jī)生成200個和的統(tǒng)計直方圖Fig .6 Statistical histogram of 200 randomly generated

圖7 ρobj和ρ的多項(xiàng)式擬合結(jié)果Fig .7 Polynomial fitting results for ρobj and ρ

(31)

3 模型驗(yàn)證及算例分析

主要針對前6種情況進(jìn)行驗(yàn)證。將半解析公式得到的結(jié)果和200次隨機(jī)模擬取平均值得到的結(jié)果進(jìn)行對比，并計算R2，計算得到的結(jié)果如圖8所示。從圖8可以看出，R2均大于0.99，擬合效果較好。

在實(shí)際的應(yīng)用中，通常通過對致密油藏的巖心或天然露頭的裂縫進(jìn)行統(tǒng)計分析，得到各參數(shù)的均值方差以及它們之間的相關(guān)系數(shù)，再用本文建立的半解析方法進(jìn)行相關(guān)系數(shù)控制下的裂縫隨機(jī)生成，代入到致密油藏數(shù)值模擬器進(jìn)行數(shù)值模擬計算，進(jìn)行產(chǎn)量的預(yù)測等。

裂縫真實(shí)參數(shù)的獲取是一件較為麻煩的事情，為了對相關(guān)系數(shù)控制下的裂縫生成結(jié)果進(jìn)行分析，研究不同相關(guān)系數(shù)對縫網(wǎng)形態(tài)的影響規(guī)律，本文構(gòu)建了一個理想的二維模型，假設(shè)其裂縫各參數(shù)服從的統(tǒng)計分布及其參數(shù)如表 4所示，相關(guān)系數(shù)如表 5所示。

在400 m×200 m的區(qū)域，在表5所示的相關(guān)系數(shù)控制下生成300條裂縫，其各參數(shù)之間的相關(guān)系數(shù)的相對誤差如表 7所示。

從表7可以看出，相對誤差都控制在了25%以內(nèi)，對于相關(guān)系數(shù)絕對值大于0.5的分量，相對誤差均控制在了6%以內(nèi)，6次隨機(jī)模擬得到的裂縫圖像如圖9所示。

圖8 200次隨機(jī)模擬取平均值結(jié)果和半解析結(jié)果的對比Fig .8 Comparison of semi-analytical results and the average of 200 random simulations

表4 理想模型裂縫各參數(shù)服從的統(tǒng)計分布及其參數(shù)Table 4 Statistical distribution of the fracture parameters of the ideal model

表5 裂縫5個參數(shù)之間的相關(guān)系數(shù)Table 5 Correlation coefficients between 5 parameters of fractures

表6 當(dāng)時的對數(shù)正態(tài)分布與其他各分布的ρ和ρobj的關(guān)系Table 6 Relationship of ρ and ρobj between the lognormal distribution and others when

表7 按半解析方法得到的50條裂縫各參數(shù)之間的相關(guān)系數(shù)的相對誤差Table 7 Relative errors of correlation coefficients between parameters of 50 fractures obtained by semi-analyticalmethod %

如果只按照表所示的參數(shù)而忽視了表所示的相關(guān)系數(shù)的信息，隨機(jī)生成的裂縫各參數(shù)之間的相關(guān)系數(shù)誤差會很大。在無相關(guān)系數(shù)控制的情況下隨機(jī)生成了6組天然裂縫，并計算了相關(guān)系數(shù)的相對誤差的平均值，得到的裂縫圖像如圖 10所示，結(jié)果如表 8所示。結(jié)果表明，無相關(guān)系數(shù)控制下得到的相對誤差均大于95%，縫網(wǎng)的幾何形態(tài)也有較大的差異。

如果把裂縫各參數(shù)的相關(guān)系數(shù)改為如表9所示，得到的裂縫圖像如圖11所示。

從圖 9和圖 11可以很直觀地看出，即使是同一組隨機(jī)參數(shù)，對于不同的相關(guān)系數(shù)控制，得到的天然裂縫的縫網(wǎng)結(jié)構(gòu)在幾何形態(tài)上有著很大的差異。裂縫中心點(diǎn)橫坐標(biāo)xc和縱坐標(biāo)yc的相關(guān)性決定了裂縫富集的方向，當(dāng)xc和yc呈正相關(guān)時，裂縫向西南-東北對角線富集(圖 11)，當(dāng)xc和yc呈負(fù)相關(guān)時，裂縫向西北-東南對角線富集(圖 9)。裂縫方位角和中心點(diǎn)橫坐標(biāo)xc的相關(guān)性決定了裂縫沿x軸方向的走向趨勢線形狀，當(dāng)θ和xc呈正相關(guān)時，xc越大，θ越大，裂縫沿x軸方向的走向趨勢線呈上凸形狀(圖 11)，當(dāng)θ和xc呈負(fù)相關(guān)時，xc越大，θ越小，裂縫沿x軸方向的走向趨勢線呈下凹形狀(圖 9)。

圖9 表5所示的相關(guān)系數(shù)控制下6次隨機(jī)模擬生成的300條裂縫Fig .9 Six images of 300 fractures generated with the control of the correlation coefficients shown in Table 5

圖10 無相關(guān)系數(shù)控制下生成的300條裂縫Fig .10 Image of 300 fractures generated without the control of correlation coefficients

表8 無相關(guān)系數(shù)控制下生成的300條裂縫各參數(shù)之間的相關(guān)系數(shù)的相對誤差Table 8 Relative errors of the correlation coefficients between the parameters of 300 fractures generated without the control of correlation coefficients %

表9 裂縫各參數(shù)呈正相關(guān)的相關(guān)系數(shù)矩陣Table 9 Correlation coefficient matrix with positive correlation of fracture parameters

圖11 表9所示的相關(guān)系數(shù)控制下生成的300條裂縫Fig .11 Image of 300 fractures generated with the control of the correlation coefficients shown in Table 9

4 結(jié)論

本文引入相關(guān)系數(shù)來定量表征裂縫各參數(shù)的相關(guān)性，建立了基于相關(guān)系數(shù)控制的致密油藏裂縫隨機(jī)表征方法；并且將裂縫參數(shù)的每一種隨機(jī)分布都轉(zhuǎn)化為均勻分布來生成，從而給出了相關(guān)系數(shù)控制下的裂縫各參數(shù)的隨機(jī)生成方法。本文用半解析的方法給出了目標(biāo)相關(guān)系數(shù)和初始相關(guān)系數(shù)之間的函數(shù)關(guān)系，為相關(guān)系數(shù)控制下的裂縫隨機(jī)生成提供了更加快捷有效的方法。

本文提供的方法可以有效的控制裂縫各參數(shù)的相關(guān)性，得到的裂縫參數(shù)更加符合實(shí)際情況，為致密油藏裂縫系統(tǒng)的數(shù)值模擬和縫網(wǎng)反演提供了基礎(chǔ)。