桂俊川，陳 平，馬天壽

“油氣藏地質(zhì)及開發(fā)工程”國家重點(diǎn)實(shí)驗(yàn)室·西南石油大學(xué)，四川 成都 610500

引 言

頁巖氣是一種重要的非常規(guī)能源，受北美頁巖氣革命的影響，中國也加入到了頁巖氣的勘探與開發(fā)熱潮中[1-2]。頁巖氣具有低孔低滲的特點(diǎn)，水平井鉆井技術(shù)和體積壓裂改造是頁巖氣成功開發(fā)的兩大關(guān)鍵技術(shù)[3]。巖石的力學(xué)特性是影響井壁穩(wěn)定，造成井眼垮塌、漏失的重要因素[4]。在水力壓裂過程中，不同空間位置的巖石力學(xué)特性展布會影響裂縫的形成與擴(kuò)展，例如，在膠結(jié)程度較弱的層理面往往先于頁巖本體開裂，使得水力裂縫在延伸過程中沿層理面優(yōu)先擴(kuò)展，從而抑制縫網(wǎng)的形成，降低水力壓裂效率[5]。因此，研究巖石的力學(xué)特性對油氣勘探開發(fā)有著重要的意義。

相比橫向各向同性（VTI）介質(zhì)，在地下巖層中正交各向異性（ORT）介質(zhì)是更為廣泛的存在，尤其在沉積盆地中是最普通的情況之一。目前與沉積巖，如頁巖相關(guān)的巖石力學(xué)研究中，主流觀點(diǎn)是將地層巖石視為宏觀同性介質(zhì)[6-8]。海、河水的方向性流動使得巖石顆粒的定向排列，水平最大最小主應(yīng)力的差異造成兩個(gè)方向上的裂縫、孔隙具有不同程度的擠壓或新裂縫的產(chǎn)生，以及鉆井過程中鉆遇背（向）斜地層等是造成地層表現(xiàn)出正交各向異性的主要因素[9-10]。忽視地層的正交各向異性，利用VTI模型獲取的泊松比、楊氏模量等巖石力學(xué)參數(shù)計(jì)算地應(yīng)力，并用于進(jìn)一步的井壁穩(wěn)定分析和壓裂施工指導(dǎo)等，有時(shí)會導(dǎo)致錯(cuò)誤的結(jié)果。

與VTI介質(zhì)不同，ORT介質(zhì)研究相對較少。Hudson給出裂隙介質(zhì)的本構(gòu)關(guān)系[11]，希望通過研究地震波在裂縫介質(zhì)中的響應(yīng)特征預(yù)測儲層目標(biāo)區(qū)域的裂縫方位、填充物等微觀裂縫參數(shù)。Crampin通過將定向排列的充滿液體的垂直裂隙加入到周期性薄層的橫向各向同性介質(zhì)中，研究了正交各向異性的橫波偏振特征[12]。Bush在巴黎盆地獲取的資料中觀察到了橫波分裂現(xiàn)象，使人們認(rèn)識到沉積盆地中廣泛存在著正交各向異性[13]。Tsvankin將VTI介質(zhì)的Thomsen參數(shù)推廣到了ORT介質(zhì)，便于研究裂縫引起的各向異性[14]。Sehoenberg和Helbig從垂直裂縫介質(zhì)的彈性方程出發(fā)，模擬了ORT介質(zhì)中3個(gè)正交對稱面上縱橫波的波速面[15]。Bakulin等對裂縫誘導(dǎo)的正交各向異性進(jìn)行了研究，利用地震波反射數(shù)據(jù)估計(jì)裂縫參數(shù)[16]。Cheng改進(jìn)了Hudson的一階擾動模型，在Eshelby的基礎(chǔ)上建立了Eshelby-Cheng模型[17]。劉恩儒和曾新吾將多種裂縫模型總結(jié)為3類，推導(dǎo)了3種模型對應(yīng)的柔度解析表達(dá)式，給出統(tǒng)一的解析形式[18]。Franquet和Rodriguez利用測井資料，基于剛度矩陣參數(shù)，建立了正交各向異性地層的地應(yīng)力計(jì)算方法[19]。Alkhalifah利用在聲學(xué)介質(zhì)假設(shè)下導(dǎo)出的色散關(guān)系，獲得了正交介質(zhì)的聲波方程[20]。馬妮等利用地震資料根據(jù)水平應(yīng)力差異比計(jì)算了正交各向異性地層的地應(yīng)力[21]。物理模型實(shí)驗(yàn)方面，Chaedle等最先對正交各向異性樣品進(jìn)行了速度測試，測得了3個(gè)主軸和6個(gè)45°對角方向上的縱橫波速度[22]。李躍等用鋁板和有機(jī)玻璃板構(gòu)成一個(gè)薄互層等效各向異性模型，討論了薄層厚度對各向異性的影響[23]。魏建新通過實(shí)驗(yàn)觀測了正交各向異性介質(zhì)的聲波特征[10，24]。李軍等利用現(xiàn)場實(shí)鉆巖芯進(jìn)行了巖石力學(xué)性質(zhì)正交各向異性實(shí)驗(yàn)，并建立了正交各向異性地層井壁圍巖應(yīng)力模型[25-26]。水平井鉆井和體積壓裂是目前頁巖氣勘探開發(fā)的主要方式，井筒與地層的夾角（井斜角和方位角）隨井眼軌跡的變化而變化。而以上研究成果都是基于巖石本構(gòu)坐標(biāo)系而得到的，沒有考慮巖石力學(xué)特性隨觀測坐標(biāo)系的變化，得到的結(jié)果不便于水平井鉆井和壓裂增產(chǎn)改造的直接應(yīng)用。對于VTI介質(zhì)，由于觀測坐標(biāo)系的不同，造成的巖石力學(xué)特性差異已有一些相關(guān)研究[27-28]，而ORT介質(zhì)則無相應(yīng)報(bào)道。

針對上述問題，考慮實(shí)際頁巖地層中的各向異性特點(diǎn)，從巖石本構(gòu)方程出發(fā)，推導(dǎo)出ORT介質(zhì)巖石力學(xué)參數(shù)的彈性表達(dá)式。然后考慮井眼觀測坐標(biāo)系與巖石本構(gòu)坐標(biāo)系的不同交角，根據(jù)巖石的力學(xué)參數(shù)定義，建立了觀測坐標(biāo)系下ORT介質(zhì)的巖石力學(xué)參數(shù)計(jì)算方法，并進(jìn)行了分析。最后，考慮頁巖與砂巖的互相沉積成巖，研究了不同砂巖含量下頁巖-砂巖互層在觀測坐標(biāo)系下的巖石力學(xué)特性。

1 巖石本構(gòu)方程的各向異性表征

對于一個(gè)一般各向異性的線性彈性固體，應(yīng)力σ與應(yīng)變ε之間存在如下線性關(guān)系[29]

式中：

σ—應(yīng)力張量，GPa；

ε—應(yīng)變張量，無因次；

C—巖石剛度矩陣，GPa。

由于對稱性，C最多具有21個(gè)獨(dú)立常數(shù)。反過來，應(yīng)變也可以表示為應(yīng)力的線性組合

式中：S—巖石的柔度矩陣，GPa-1。

柔度矩陣S和剛度矩陣C互為逆矩陣。對于如圖1所示的正交各向異性巖石，其具有C11、C12、C13、C22、C23、C33、C44、C55、C66共 9 個(gè)獨(dú)立的彈性剛度參數(shù)。

圖1 正交各向異性巖石示意圖（本構(gòu)坐標(biāo)系）Fig.1 Schematic diagram of orthotropic rock（Constitutive coordinate system）

在正交各向異性條件下，剛度矩陣C0和柔度矩陣S0可以分別表示如下

式中：

C0—正交各向異性條件下巖石剛度矩陣，GPa；

Cij—巖石剛度系數(shù)（i=1,2,3,4,5,6，j=1,2,3,4,5,6），GPa。

式中：

S0—正交各向異性條件下巖石柔度矩陣，GPa-1；

E1—沿x軸方向的楊氏模量，GPa；

E2—沿y軸方向的楊氏模量，GPa；

E3—沿z軸方向的楊氏模量，GPa；

μ12—垂直于x軸，沿y軸的泊松比，無因次；

μ13—垂直于x軸，沿z軸的泊松比，無因次；

μ21—垂直于y軸，沿x軸的泊松比，無因次；

μ23—垂直于y軸，沿z軸的泊松比，無因次；

μ31—垂直于z軸，沿x軸的泊松比，無因次；

μ32—垂直于z軸，沿y軸的泊松比，無因次；

G12—xy平面的剪切模量，GPa；

G13—xz平面的剪切模量，GPa；

G23—yz平面的剪切模量，GPa。

利用Sonic Scanner[30]測井平臺能夠獲取包括：C33、C44、C55、C66在內(nèi)的 4個(gè)參數(shù)，其余的參數(shù)需要通過實(shí)驗(yàn)獲取。目前，有兩類方法可以獲取這些參數(shù)，巖芯的超聲波實(shí)驗(yàn)測試和室內(nèi)三軸巖石力學(xué)壓縮實(shí)驗(yàn)。

用超聲波波速測試確定這9個(gè)參數(shù)需要測量6個(gè)方向的縱波波速和3個(gè)方向的橫波波速[22，31]，分別包括：（1）沿x、y、z軸方向傳播縱波波速vp-x、vp-y、vp-z，如圖 2a；（2）沿z軸傳播，偏振方向?yàn)閤（y）軸的橫波波速vsv-x、vsv-y和沿y軸方向傳播、偏振方向?yàn)閤軸的橫波波速vsh，如圖2b；（3）在xy（zx、yz）平面與x（z、y）軸成 45°角方向傳播的縱波波速 vp-xy45°、vp-zx45°、vp-yz45°，如圖 2c。

圖2 P-波和S-波波速測量示意Fig.2 The P-wave and S-wave velocities measurement

圖2為縱橫波波速測量示意圖，圖中長虛線單向箭頭表示波的傳播方向，短虛線雙向箭頭表示波的振動方向。結(jié)合巖石的密度ρ和測量的波速，9個(gè)獨(dú)立的彈性剛度參數(shù)可表示如式（5）、式（6）、式（7）所示

式中：

C11，C22，C33剛度系數(shù)，GPa；

ρ—巖石密度，g/cm3；

vp—x沿x軸傳播縱波波速，km/s；

vp—y沿y軸傳播縱波波速，km/s；

vp—z沿z軸傳播縱波波速，km/s。

式中：

C44，C55，C66—剛度系數(shù)，GPa；

vsv-x—沿z軸傳播，偏振方向?yàn)閤軸的橫波波速，km/s；

vsv-y—沿z軸方向傳播，偏振方向?yàn)閥軸的橫波波速，km/s；

vsh—沿y軸方向傳播，偏振方向?yàn)閤軸的橫波波速，km/s。

式中：

C12，C13，C23—剛度系數(shù)，GPa；

vp-xy45°—在xy平面，與x軸夾角成45°的縱波波速，km/s；

vp-zx45°—在zx平面，與z軸夾角成45°的縱波波速，km/s；

vp-yz45°—在yz平面，與y軸夾角成45°的縱波波速，km/s。

如果采用三軸室內(nèi)巖石力學(xué)實(shí)驗(yàn)確定這些參數(shù)，則需要測量沿上述6個(gè)方向取芯巖樣的應(yīng)力應(yīng)變曲線就能完全確定巖石中的各向異性[32]。

2 巖石力學(xué)參數(shù)的彈性表征

由上述可知，由室內(nèi)實(shí)驗(yàn)或現(xiàn)場數(shù)據(jù)可以得到相應(yīng)的彈性剛度系數(shù)，而要獲取對應(yīng)的巖石力學(xué)參數(shù)（泊松比、楊氏模量），則需要建立彈性剛度參數(shù)與巖石力學(xué)參數(shù)的關(guān)系。將式（3）取逆，并結(jié)合式（4），可以推出巖石力學(xué)參數(shù)關(guān)于巖石剛度系數(shù)的表達(dá)式，包括3個(gè)楊氏模量和6個(gè)泊松比。3個(gè)楊氏模量的表達(dá)式為

其中：

3個(gè)方向6個(gè)泊松比為

楊氏模量和泊松比的關(guān)系滿足

E2μ12=E1μ21、E3μ13=E1μ31、E2μ32=E3μ23

由式（9）可知

μ12≠μ21、μ13≠μ31、μ23≠μ32。

當(dāng)獲得了9個(gè)獨(dú)立的彈性剛度參數(shù)C11、C12、C13、C22、C23、C33、C44、C55、C66之后，將其代入式（8）、式（9），便可得到對應(yīng)的巖石力學(xué)參數(shù)，這樣便建立起了巖石力學(xué)參數(shù)的彈性剛度計(jì)算方法。

對于VTI介質(zhì)，彈性剛度參數(shù)滿足

則式（8）、式（9）退化為

式中：Ev—垂直于層理面的楊氏模量，GPa；

Eh—平行于層理面的楊氏模量，GPa；

μv—垂直于層理面的泊松比，無因次；

μh—平行于層理面的泊松比，無因次。

3 巖石力學(xué)參數(shù)的空間變換

在鉆井過程中，井眼的方位角和井斜角隨時(shí)都在變化，這導(dǎo)致了觀測坐標(biāo)系和巖石本構(gòu)坐標(biāo)系不一致。盡管室內(nèi)實(shí)驗(yàn)可以獲取任意對應(yīng)方位角和井斜角條件下的巖石力學(xué)參數(shù)，但巖芯的獲取較為困難，測試費(fèi)用高，并且鉆取的角度有限。因此，如何利用測量到的巖石本構(gòu)坐標(biāo)系下彈性剛度矩陣，去計(jì)算任意觀測坐標(biāo)系下的巖石力學(xué)參數(shù)，具有很強(qiáng)現(xiàn)實(shí)意義。本文利用Bond變換，首先，將巖石本構(gòu)坐標(biāo)系下的彈性剛度矩陣轉(zhuǎn)換到觀測坐標(biāo)系，然后，根據(jù)巖石力學(xué)參數(shù)的定義，計(jì)算任意方位井斜角條件下的巖石力學(xué)參數(shù)，具體流程如下。

假定地層巖石的本構(gòu)坐標(biāo)系為O-x0y0z0，且Ox0軸，Oy0軸和Oz0軸分別與主應(yīng)力σH、σh和σv方向一致。本構(gòu)坐標(biāo)系O-x0y0z0先以O(shè)z0為軸旋轉(zhuǎn)方位角α，變?yōu)镺-x1y1z1坐標(biāo)系，再以O(shè)y1為軸旋轉(zhuǎn)井斜角β，則得到觀測坐標(biāo)系（井筒坐標(biāo)系）O-xyz，如圖3所示。

圖3 觀測坐標(biāo)系的空間轉(zhuǎn)換Fig.3 Spatial transformation of the observation coordinate system

巖石本構(gòu)坐標(biāo)系（O-x0y0z0）和觀測坐標(biāo)系（O-xyz）之間的坐標(biāo)轉(zhuǎn)換關(guān)系可由坐標(biāo)之間的方向余弦aij表示

式中：

aij—方向余弦（i,j=1,2,3），無因次；

α—方位角，弧度；

β—井斜角，弧度。

假設(shè) σ、ε、C、S和 σ0、ε0、C0、S0分別為觀測坐標(biāo)系和本構(gòu)坐標(biāo)系下的應(yīng)力、應(yīng)變、剛度矩陣和柔度矩陣，則觀測坐標(biāo)系和本構(gòu)坐標(biāo)系下應(yīng)力、應(yīng)變滿足[33-34]

式中：M，N過渡矩陣。

過渡矩陣M、N滿足N-1=M，M矩陣可由aij的表達(dá)式計(jì)算。

將式（11）代入到上述表達(dá)式，可得

觀測坐標(biāo)系及本構(gòu)坐標(biāo)系下的剛度矩陣C和C0，柔度矩陣S和S0分別滿足以下關(guān)系

式中：

MT__矩陣M的轉(zhuǎn)置矩陣，無因次；

NT__矩陣N的轉(zhuǎn)置矩陣，無因次。

根據(jù)式（1）、式（2），可以得到觀測坐標(biāo)系下應(yīng)力和應(yīng)變的關(guān)系式

通過變換，能夠求得任意方位角和井斜角條件下的柔度矩陣S。利用式（17），根據(jù)楊氏模量和泊松比的定義，即可計(jì)算任意方位角度下的楊氏模量和泊松比。

以觀測坐標(biāo)系下z軸方向的巖石力學(xué)參數(shù)為例。根據(jù)定義[35]，楊氏模量為巖石受到單軸應(yīng)力時(shí)，應(yīng)力相對于應(yīng)變的變化率。而泊松比的定義則是巖石在單軸受壓條件下橫向應(yīng)變與縱向應(yīng)變之比。在觀測坐標(biāo)系下，假定巖石只受到單軸應(yīng)力σz的作用，即

式中：σz—沿z軸的單軸應(yīng)力，GPa。

將式（18）代入式（17），可以得到觀測坐標(biāo)系下與應(yīng)力σ相關(guān)的應(yīng)變ε表達(dá)式，假設(shè)為

式中：εxx—沿x軸的應(yīng)變，無因次；

εyy—沿y軸的應(yīng)變，無因次；

εzz—沿z軸的應(yīng)變，無因次。

則根據(jù)定義，可以求取該方向的楊氏模量和泊松比

若式（17）中的S0用C0表示，則式（20）中解出的楊氏模量和泊松比只與C0、方位角α和井斜角β有關(guān)。類似地，可以求取觀測坐標(biāo)系下x軸和y軸方向的楊氏模量和泊松比。在正交各向異性條件下，觀測坐標(biāo)系中關(guān)于C0、α和β的楊氏模量和泊松比的解析表達(dá)式非常復(fù)雜，這里不給出解析表達(dá)式，下面進(jìn)行數(shù)值分析。

在本構(gòu)坐標(biāo)系下，本文測量了頁巖的剛度系數(shù)：C11=31.89 GPa、C12=10.29 GPa、C13=9.94 GPa、C22=28.18 GPa、C23=10.11 GPa、C33=24.48 GPa、C44=10.64 GPa、C55=10.72 GPa、C66=10.80 GPa。由于測量過程中巖芯露頭有限，沒有測量到C22和C23的數(shù)據(jù)，這里C22取C11和C33的平均值，C23取C12和C13的平均值。利用式（8）、式（9）計(jì)算出本構(gòu)坐標(biāo)系下頁巖的巖石力學(xué)參數(shù)：E1=26.26 GPa，E2=22.63 GPa，E3=19.47 GPa，μ12=0.258，μ13=0.300，μ21=0.222，μ23=0.323，μ31=0.222，μ32=0.278。為方便表述，取Emax=26.26 GPa，Emin=19.47 GPa，μmax=0.323，μmin=0.222。

方位角和井斜角對巖石力學(xué)參數(shù)的影響均是關(guān)于Π的周期函數(shù)。因此，取方位角α和井斜角β的變化范圍均為180°，間隔為15°。這里解出了一個(gè)周期內(nèi)所有方位角、井斜角所對應(yīng)的巖石力學(xué)參數(shù)，為節(jié)省篇幅，分析了由VTI模型和ORT模型計(jì)算的楊氏模量所產(chǎn)生的相對誤差。

圖4a為不同觀測坐標(biāo)系下沿x軸方向的楊氏模量E1，圖4b為采用ORT模型和VTI模型計(jì)算的E1所產(chǎn)生的相對誤差。E1在井斜角30°（150°）和方位角 30°（150°）時(shí)取最大值 26.67 GPa，超過Emax，井斜角為90°時(shí)取最小值19.47 GPa。用VTI模型計(jì)算的楊氏模量在方位角30°～150°和井斜角0°～60°、120°～180°所產(chǎn)生的誤差較大，最大值超過了15%。

圖5a為不同觀測坐標(biāo)系下沿y軸方向的楊氏模量E2，圖5b為采用ORT模型和VTI模型計(jì)算的E2所產(chǎn)生的相對誤差。E2隨著方位角角度的增加呈現(xiàn)增、降、增、降的趨勢，并在方位角為60°（120°）時(shí)取得最大值 26.41 GPa，超過Emax，井斜角的變化對E2沒有影響。用VTI模型計(jì)算的楊氏模量在方位角 0°～60°和 120°～180°所產(chǎn)生的誤差較大，最大值超過了15%。

圖4 沿x軸方向楊氏模量E1及其相對誤差Fig.4 Young′s modulus E1along the x-axis and its relative error

圖5 沿y軸方向楊氏模量E2及其相對誤差Fig.5 Young′s modulus E2along the y-axis and its relative error

圖6a為不同觀測坐標(biāo)系下沿z軸方向的楊氏模量E3，圖6b為采用ORT模型和VTI模型計(jì)算E3所產(chǎn)生的相對誤差。E3隨著方位角和井斜角角度的變化趨勢復(fù)雜，在井斜角為0°（180°）時(shí)取最小值19.47 GPa，井斜角為60°（120°）和方位角為30°（150°）時(shí)取最大值 26.67 GPa，超過Emax。VTI模型計(jì)算的楊氏模量在方位角45°～135°和井斜角45°～135°所產(chǎn)生的誤差較大，最大值超過了15%。

圖7為不同坐標(biāo)系下沿x軸方向的兩個(gè)泊松比，可以看出同一方向的兩個(gè)泊松比具有大致類似的變化趨勢。μ21隨井斜角角度的增加呈現(xiàn)增、降的趨勢，在井斜角 0°（180°）和方位角 45°（135°）時(shí)取最小值 0.193，小于 μmin，在井斜角 30°（150°）和方位角90°時(shí)取最大值0.331，大于μmax。μ31隨井斜角角度的增加呈現(xiàn)降、增、降、增的趨勢，在井斜角45°（135°）和方位角0°（180°）時(shí)取最小值0.167，小于μmin，在井斜角90°和方位角60°（120°）時(shí)取最大值 0.331，大于 μmax。

圖8為不同觀測坐標(biāo)系下沿y軸方向的兩個(gè)泊松比，這兩個(gè)泊松比隨方位角和井斜角的變化并無明顯的對應(yīng)關(guān)系，但具有類似的變化趨勢。μ12隨井斜角和方位角角度的變化關(guān)系復(fù)雜，在井斜角0°（180°）和方位角45°（135°）時(shí)取最小值0.193，小于 μmin，在井斜角 0°（180°）和方位角 45°（135°）時(shí)取最大值0.301。μ32隨井斜角和方位角角度的變化復(fù)雜，在井斜角 90°和方位角 45°（135°）時(shí)取最小值 0.193，在井斜角 45°（135°）和方位角0°（180°）時(shí)取最大值 0.301。

圖6 沿z軸方向楊氏模量E3及其相對誤差Fig.6 Young′s modulus E3along the z-axis and its relative error

圖7 沿x軸方向泊松比μ21和μ31Fig.7 Poisson′s ratiosμ21and μ31along the x-axis

圖8 沿y軸方向泊松比μ12和μ32Fig.8 Poisson′s ratiosμ12and μ32along the y-axis

不同觀測坐標(biāo)系下沿z軸方向的兩個(gè)泊松比隨井斜角和方位角的變化具有一定的相似趨勢，但又有所不同（圖9）。μ13隨井斜角的增加呈現(xiàn)降、增、降、增的趨勢，在井斜角 45°（135°）和方位角0°（180°）時(shí)取最小值 0.167，小于 μmin，在井斜角0°（180°）和方位角 60°（120°）時(shí)取最大值 0.331，大于μmax。μ23隨井斜角的增加呈現(xiàn)降、增的趨勢，在井斜角 90°和方位角 45°（135°）時(shí)取最小值 0.193，小于 μmin，在井斜角 0°（180°）和方位角30°（150°）時(shí)取最大值 0.331，大于 μmax。

圖9 沿z軸方向泊松比μ13和μ23Fig.9 Poisson′s ratiosμ13and μ23along the z-axis

4 層狀巖石的等效彈性性質(zhì)

沉積巖在成巖過程中，通常都是兩種或兩種以上的巖石互相沉積而成。以頁巖為例，淺色的粉砂巖和暗色的泥頁巖互層（如圖1所示），韻律交替發(fā)育是造成頁巖成層分布的主要原因[36-37]。不同的巖性具有不同的彈性力學(xué)性質(zhì)，研究砂-頁巖互層組成的巖石的彈性性質(zhì)，及其對應(yīng)的巖石力學(xué)特性，為鉆井和壓裂提供基礎(chǔ)數(shù)據(jù)和指導(dǎo)依據(jù)，具有重要的現(xiàn)實(shí)意義和工程意義。VTI介質(zhì)的層狀巖石的等效彈性性質(zhì)可以由Backus理論進(jìn)行模擬[38]，對于ORT介質(zhì)，則可以用Gerrard提出的模型來分析多層正交各向異性巖石的彈性性質(zhì)[39]。假定巖石是由頁巖和砂巖兩種巖體互相沉積而成，兩種巖體的體積分別為V1和V2（V1+V2=1），并且兩種巖石本構(gòu)坐標(biāo)系的對稱軸互相重合。

該模型需要滿足兩點(diǎn)假設(shè)：（1）與任何單個(gè)層的厚度相比，層系統(tǒng)的規(guī)模較大；（2）層的厚度和材料特性相對于它們在系統(tǒng)內(nèi)的各自位置而隨機(jī)變化。模型的具體表達(dá)式如下[39-40]

其中：

5 算例分析

假定砂巖的剛度系數(shù)為C11=49.27 GPa、C12=18.00 GPa、C13=16.24 GPa、C22=45.38 GPa、C23=17.12 GPa、C33=41.49 GPa、C44=14.94 GPa、C55=15.29 GPa、C66=15.63 GPa。由式（8）、式（9）計(jì)算出本構(gòu)坐標(biāo)系下砂巖的巖石彈性參數(shù)：E1=39.58 GPa，E2=35.34 GPa，E3=32.91 GPa，μ12=0.295，μ13=0.270，μ21=0.263，μ23=0.310，μ31=0.224，μ32=0.288。取砂巖體積V2分別為 0，0.1，…，1.0，計(jì)算出混合巖石的等效巖石力學(xué)參數(shù)。

圖10為楊氏模量隨砂巖體積含量的變化關(guān)系，隨著砂巖體積含量的增加，楊氏模量呈現(xiàn)增加的趨勢，其中，E1和E2呈線性增加，而E3以一種近線性的方式增加。

圖10 楊氏模量隨砂巖體積含量的變化Fig.10 Variation of Young′s modulus with volume content of sandstone

圖11為泊松比隨砂巖體積含量的變化關(guān)系，與楊氏模量不同，泊松比隨砂巖體積含量的變化關(guān)系更為復(fù)雜。μ12和μ21隨砂巖體積含量的增加以近似線性的方式增加，μ13和μ23隨砂巖體積含量的增加以近似線性的方式降低，μ31和μ32隨砂巖體積含量的增加先降低后增加。

為研究不同砂巖含量條件下頁巖和砂巖互層的巖石力學(xué)參數(shù)在觀測坐標(biāo)系下的變化特征，取砂巖體積V2為0.2、0.4、0.6、0.8，利用上述流程，模擬了不同方位角和井斜角條件下的巖石力學(xué)參數(shù)的分布特征。考慮到E2μ12=E1μ21、E3μ13=E1μ31、E2μ32=E3μ23。本文只給出了3個(gè)楊氏模量E1、E2、E3和3個(gè)泊松比μ12、μ23、μ31在不同砂巖體積含量下隨方位角和井斜角的變化關(guān)系。

圖12為不同觀測坐標(biāo)系下楊氏模量E1在不同砂巖體積含量條件下的變化圖。可以看出，隨著砂巖體積含量的增加，楊氏模量隨方位角和井斜角變化的曲面形態(tài)基本不變，數(shù)值上逐漸增加。圖12a中紅色圓圈部分相比于圖12d中的對應(yīng)部分的楊氏模量更為平緩。

圖13為不同觀測坐標(biāo)系下楊氏模量E2在不同砂巖體積含量條件下的變化圖。隨著砂巖體積含量的增加，E2隨方位角和井斜角變化的曲面形態(tài)完全一致，砂巖體積含量越高，E2越大。

圖14為不同觀測坐標(biāo)系下楊氏模量E3在不同砂巖體積含量下的變化圖。楊氏模量隨方位角和井斜角變化的曲面形態(tài)隨著砂巖體積含量的增加基本保持不變，數(shù)值上逐漸增加。相比于圖14a中紅色圓圈部分，圖14d中的對應(yīng)部分的楊氏模量更為陡峭。

圖13 不同砂巖體積含量條件下的E2Fig.13 Young′s modulus E2under different sandstone volume contents

圖14 不同砂巖體積含量條件下的E3Fig.14 Poisson′s ratio E3under different sandstone volume contents

圖15為在不同砂巖體積含量條件下泊松比μ12隨方位角和井斜角的變化，不同砂巖體積含量下μ12曲面形態(tài)基本一致，數(shù)值上略有增加。砂巖體積含量越高，μ12曲面的變化越趨于平緩。圖14a中紅色圓圈部分，4個(gè)曲面角向下折，而對應(yīng)的圖15d中4個(gè)曲面角向上翹。

圖16為在不同砂巖體積含量條件下泊松比μ23隨方位角和井斜角的變化，不同砂巖體積含量下μ23曲面形態(tài)基本一致。砂巖體積含量越高，μ23曲面的變化更趨于平緩，并且在方位角90°（圖16a中紅色橢圓部分）附近出現(xiàn)較為明顯的增加，其余部分則變化不大。

圖15 不同砂巖體積含量條件下的μ12Fig.15 Poisson′s ratio μ12under different sandstone volume contents

圖16 不同砂巖體積含量條件下的μ23Fig.16 Poisson′s ratio μ23under different sandstone volume contents

圖17為泊松比μ31在不同砂巖體積含量條件下隨方位角和井斜角的變化，不同砂巖體積含量下μ31曲面形態(tài)基本一致。砂巖體積含量的增加使得μ31曲面的變化更趨于平緩，這種變化主要體現(xiàn)在相同方位角不同井斜角的剖面上。而相同井斜角下，不同方位角對曲面形態(tài)的影響較小。

圖17 不同砂巖體積含量條件下的μ31Fig.17 Poisson′s ratio μ31under different sandstone volume contents

6 結(jié) 論

（1）頁巖的楊氏模量隨方位角和井斜角的變化趨勢復(fù)雜，3個(gè)楊氏模量E1、E2、E3隨方位角和井斜角變化的最大值均超過了本構(gòu)坐標(biāo)系下楊氏模量的最大值26.26 GPa，通過比較由VTI模型和ORT模型計(jì)算的楊氏模量所產(chǎn)生的相對誤差，發(fā)現(xiàn)誤差的最大值均超過了15%，正交各向異性條件下的巖石力學(xué)參數(shù)不能用VTI模型近似。

（2）頁巖的泊松比隨方位角和井斜角的變化并無明顯的對應(yīng)關(guān)系，但同一方向的泊松比隨方位角和井斜角的變化具有大致相似的變化趨勢。在觀測坐標(biāo)系下，所有泊松比的最小值均小于本構(gòu)坐標(biāo)系下的最小值0.222，除了μ12、μ32以外，其余的4個(gè)泊松比隨方位角和井斜角變化的最大值均大于本構(gòu)坐標(biāo)系下的最大值0.323。

（3）在本構(gòu)坐標(biāo)系下，隨砂巖和頁巖互層砂巖體積含量的增加，楊氏模量增加，泊松比μ12和μ21增加，μ13和μ23降低，μ31和μ32則先降低后增加。

（4）相同井斜角和方位角條件下，砂巖和頁巖互層楊氏模量隨砂巖含量的增加而有所增加，而泊松比則變化不一。巖石力學(xué)參數(shù)的曲面形態(tài)除了在某些方位角和井斜角有所變化外，基本保持一致。砂巖體積含量越高，則巖石力學(xué)參數(shù)曲面變化越趨于平緩。