衡 帅,杨春和,,张保平,郭印同,王 磊,魏元龙
(1.中国科学院武汉岩土力学研究所岩土力学与工程国家重点实验室,湖北 武汉 430071;2.重庆大学 煤矿灾害动力学与控制国家重点实验室,重庆 400044;3.中国石油化工股份有限公司石油工程技术研究院,北京 100101)
页岩气是指主要赋存于低孔超低渗、富含有机质的暗色泥页岩或高碳泥页岩层系中,在海相、海陆过渡相和陆相都有分布的一种非常规能源,其分布广泛,开采潜力巨大,是常规能源的理想接替[1-2]。页岩气井钻进完成后,由于地层的低渗透率,只有少数天然裂缝特别发育的井可直接投入生产,而大多数井都要经过酸化、压裂等储层改造后才能获得比较理想的产量[3-4]。但近年来,随着钻井工艺技术的不断提高,水平井由于其泄油面积大、单井产量高和穿透度大等优点,在油气开采中越来越多地被采用。对页岩气井,由于水平段一般较长,且地层层理发育,力学性质硬脆,部分还具有高陡倾角,这使井壁容易发生坍塌、漏失等严重的失稳问题[5]。井壁失稳不仅会影响钻进时效,更重要的是会降低固井质量,影响压裂作业的安全进行,由于目前对层理性页岩地层的基础理论和试验研究还相对较薄弱,井壁失稳很难得到解决。而国内外水平井压裂现场施工中经常出现的破裂压力高、裂缝压不开等现象[6]表明,目前对水平井水力裂缝的起裂机制还没有较深刻的认识。因此,根据我国页岩气藏储层特点,深入开展与水力压裂开采相关的岩石力学特性研究,探索出我国页岩气开发的合适途径具有重要的现实意义。
受层理面影响,页岩地层的力学性质、强度特征和破裂模式均表现出明显的各向异性,这在分析水平井井壁稳定性和水力裂缝的扩展问题时应予以足够重视。目前,国内外对各向异性岩石的强度特征和破裂模式进行了大量研究。Jaeger 等[7-11]对层状岩体弹性参数的各向异性进行了试验和理论研究,并建立了横观各向同性岩体的破坏准则。Jung-Woo 等[12]通过不同角度下片麻岩、页岩和片岩的单轴压缩和巴西劈裂试验,研究了这3 种岩石弹性参数及强度的各向异性。Niandou 等[13]研究了Tournemire 页岩的力学各向异性,并将其破裂模式分为剪切破坏和张拉破坏。赵文瑞[14]研究了砂岩强度随层面和主应力轴夹角变化时的各向异性,并得出了当主应力轴与弱面呈30°角时强度最低的结论。高春玉等[15]研究了砂板岩力学特性的各向异性,并分析了层理倾角对岩石变形特性、强度特征及其参数的影响。刘运思等[16]通过不同层理方向板岩的巴西劈裂试验,分析了破坏模式、抗拉强度、弹性模量和泊松比的各向异性。刘胜利等[17]通过试验研究了绿泥石片岩的各向异性特征。由此可知,虽然国内外对岩石强度和破裂模式的各向异性研究较多,但对页岩各向异性的研究成果还较少,而页岩的力学参数和强度特征不仅直接关系到页岩气压裂开采的参数设计,还对水平井井壁的稳定性分析和水力裂缝的扩展规律影响显著。因此,对页岩的强度和破裂模式的各向异性进行研究显得尤为必要。
本文基于不同层理角度页岩的单轴和三轴压缩试验,研究了层理面影响下页岩的力学特性、强度特征和破裂模式的各向异性,分析了其破坏机制的各向异性,并给出了将页岩视为横观各向同性体的5个独立材料参数,为页岩气水力压裂开采中水平钻井和压裂施工参数的设计等提供了必要的技术基础。
试验所用岩芯取自重庆彭水页岩气区块储层自然延伸的石柱县漆辽海相志留统龙马溪组露头页岩。该地层为黑色-深黑色碳质页岩,薄层-中厚层平行交互,层理面发育,层间黏结力小且易风化开裂,地层倾角为70°。清除表面风化层后,选取底部扰动相对较小的页岩,采集原始尺寸均远大于300 mm×300 mm×300 mm 的岩块,从而保证试样加工的完整性。
为了解层理面影响下页岩的力学特性、强度特征和破裂模式的各向异性,在取芯时钻取方向与层理面的夹角依次为0°、30°、60°和90°,如图1 所示。加工好的圆柱体试样直径为50 mm,长度为100 mm,误差±0.5 mm,端面平行度±0.02 mm。为避免试验结果的离散性,每组试验至少做3个试样,并取平均值。
图1 页岩定向取芯示意图Fig.1 Directional coring diagram of shale
本次试验是在美国产MTS815 Flex Test GT岩石力学测试系统上进行,该设备轴向最大加载可达2 800 kN,测试精度高,性能稳定。试验时采用轴向位移控制,其加载速率为0.18 mm/min,加载过程中采集试样的轴向力、轴向位移和环向位移。
不同层理角度页岩单轴和三轴压缩试验强度和变形参数结果的平均值如表1 所示。
表1 不同层理角度页岩单轴和三轴压缩试验结果Table 1 Results from uniaxial and triaxial compression tests for shale with different bedding orientations
由表1 得不同层理角度页岩在不同围压下弹性模量和泊松比的变化规律如图2、3 所示。
图2 弹性模量随层理角度的变化图Fig.2 Variation of elastic modulus with orientation of bedding plane
由图2 可知,单轴压缩时,随着层理角度的增大,弹性模量逐渐减小,即平行层理方向最大,垂直层理方向最小。三轴压缩时,随着围压的增大,各角度页岩的弹性模量均逐渐增大,但增加幅度逐渐减小;同一围压下,弹性模量的增加速率随层理角度的增大逐渐减小;对60°和90°页岩,围压的增大对弹性模量的变化几乎没有影响。为分析围压对弹性模量各向异性的影响,定义弹性模量各向异性度为
式中: RE为弹性模量的各向异性度;E0°为0°页岩的弹性模量;E90°为90°的弹性模量。不同围压下弹性模量的各向异性度如表2 所示。
表2 不同围压下弹性模量的各向异性度Table 2 Degrees of anisotropy of elastic modulus at different confining pressures
由表可知,该页岩地层的弹性模量各向异性度为1.768。随着围压的增加,弹性模量的各向异性度不断增加,但总体上增加速率不断减小。由0°和90°弹性模量的变化规律和表2 知,平行层理方向的弹性模量随围压的增加变化明显,这可能是由围压对层理面间的孔隙和微裂隙的压密作用而抑制侧向变形引起的;而在垂直层理方向,层理间孔隙和微裂隙的压密作用受围压影响较小,因而对轴向变形影响较小,使弹性模量随围压几乎没有变化。
图3 泊松比随层理角度的变化图Fig.3 Variation of Poisson’s ratio with orientation of bedding plane
由图3 可知,总体上,0°、30°和60°、90°页岩的泊松比随围压的增加呈现出了相反的变化规律。0°和30°页岩的泊松比随围压的增加不断增大,而60°和90°页岩的泊松比随围压的增加不断减小。这可能是层理间孔隙和微裂隙较发育,压力对水平层理方向变形影响较小,对垂直层理方向变形影响较大引起的。而低围压时,泊松比表现出的不规则变化特征是层理间孔隙和微裂隙瓦解的结果[18],高围压时,该瓦解作用受到抑制,泊松比随层理角度的变化表现出与弹性模量相似的变化规律。泊松比的该变化特征进一步说明了层理面黏结力相对较弱,为页岩地层的薄弱面,在钻井过程中,钻井液滤液易沿层理间的微裂隙浸入地层,使页岩强度降低,在井下钻具的扰动下,极易出现掉块、坍塌等井壁失稳现象,因此,层理间孔隙和微裂隙的良好发育也是页岩地层井壁易失稳的重要原因之一。
由图3 可知,在层理角度为50°~60°间大致存在一个角度,该方向页岩的泊松比将几乎始终保持为固定值,并不随围压的增加而变化,而页岩泊松比该特殊性质的可靠性,还需进一步深入的试验验证。
不同层理角度页岩在不同围压下压缩强度的变化规律如图4 所示。
图4 压缩强度随层理角度的变化图Fig.4 Variation of compressive strength with orientation of bedding plane for shale
由图4 可知,单轴压缩时,页岩的抗压强度在0°和90°时最高,30°时最低,呈现出两边高、中间低的U 型变化规律。三轴压缩时,随着围压的不断升高,0°页岩的压缩强度增加速率较快,为最高值,90°次之,30°为最小值,总体上仍呈现出两边高、中间低的U 型变化规律。因此,受层理面影响,页岩的压缩强度呈现出显著的各向异性特征。定义强度各向异性度为
式中:Rc为强度各向异性度;σcmax和 σcmin分别为压缩强度的最大值和最小值。不同围压下压缩强度的各向异性度如表3 所示。
表3 不同围压下压缩强度的各向异性度Table 3 Degrees of anisotropy of compressive strength at different confining pressures for shale
由表可知,随着围压的增加,压缩强度的各向异性度不断降低。而对该页岩地层,强度各向异性度cR=2.185,为中等各向异性[9]。
弹性模量和压缩强度的各向异性度随围压的变化规律如图5 所示。
图5 页岩各向异性度变化图Fig.5 Variation of degree of anisotropy with orientation of bedding plane for shale
由图可知,总体上,弹性模量和压缩强度的各向异性度随围压的增加呈现出了相反的变化规律。弹性模量各向异性度的增加是由压力对层理间孔隙和微裂隙的压密作用较显著引起的,而压缩强度各向异性度的减小是由围压抑制层理面开裂后破裂机制的改变引起的。因此,对页岩地层高地应力时,如果忽略其弹性模量的各向异性,将对工程实际问题的分析和设计带来较大误差;而层理面的易开裂性质受到较大抑制后,强度各向异性度将较小,较利于井壁的稳定性等。
针对Mohr-Coulomb 准则仅适用于岩体剪切破坏的缺陷,Hoek-Brown 在大量试验数据的基础上,提出了Hoek-Brown 经验强度准则[19]。该准则综合考虑了岩体强度、结构面强度和岩体结构等因素的影响,能对岩体的拉伸和剪切破坏机制进行描述,较适用于各向异性较明显的页岩。
Hoek-Brown 经验强度准则为
式中:σc为岩石的单轴抗压强度;m、s为经验参数。m 反映岩石的软硬程度,其取值范围为0.001~25,对完整坚硬岩体取25;s 与岩石内部颗粒间抗拉强度和啮合程度有关,其取值范围为0~1,对破碎岩石取0,完整岩体取1;α为与岩体特征有关的常数,对完整岩体取0.5。
对页岩,通过考虑不同层理角度的Hoek-Brown强度准则来描述其强度的各向异性。为此,引入能反映岩体强度各向异性的新参数kβ来考虑Hoek-Brown 强度准则的各向异性特征[20]。定义kβ为
式中:mβ为层理角度β 时参数m 值;mi为层理角度为90°时参数m 的值。
由式(2)、(4)知,对层理角度为90°的页岩,kβ=1。由式(3)知,不同角度的Hoek-Brown强度准则可表示为
式中:σcβ为层理角度β 时页岩的单轴抗压强度。
根据式(5),拟合不同层理角度页岩的压缩强度,得到的强度包络线如图6 所示。
图6 不同层理角度页岩的强度包络线Fig.6 Failure envelopes of shale with different bedding orientations
根据曲线拟合,不同层理角度页岩的强度参数如表4 所示。
表4 不同层理角度页岩的Hoek-Brown 强度准则参数表Table 4 Parameters in Hoek-Brown failure criterion of shale with different bedding orientations
由表4 知,不同层理角度页岩的Hoek-Brown强度准则参数差别较大,但总体上仍大致呈现了类似强度特征的U 型变化规律,各向异性特征明显。
页岩力学特性和强度的各向异性与其破裂模式密切相关。图7 展示了不同层理角度页岩在单轴和三轴压缩时的典型破裂形态。
单轴压缩时,当轴向应力达到峰值强度后,伴随着能量的突然释放,多个宏观裂缝迅速贯穿试样,使岩样失去继续承载能力,形成多个拉伸、剪切破裂面,其破裂形态具有明显的层状硬脆性岩石破坏特征,各向异性较明显。各角度岩样破坏特征如下:
图7 典型的页岩试样破裂样式图Fig.7 Typical failure modes of shale specimens
(1)β =0°:试样沿层理面的张拉劈裂破坏。破坏的岩样存在多个平行于层理面且贯通岩样两端面的张拉破裂面,这些破裂面将岩样分成多个薄板状岩块,由于裂后的岩块还能继续承载,在继续加载的过程中,岩板受压而弯曲,直至部分发生屈曲失稳而折断。
(2)β =30°:试样沿层理面形成单一剪切破坏面。破坏的岩样沿30°层理面形成贯穿整个试样的平整破裂面,试样发生明显的沿层理面的剪切滑移。
(3)β =60°:试样形成贯穿层理的多个剪切破坏面。破坏的岩样自两端形成大角度的剪切破裂面,破裂面贯穿多个层理面向岩样中部扩展,最终通过60°层理面连接,形成近似Z 型的多剪切破坏面。
(4)β =90°:试样形成贯穿层理的张拉破裂面。由于岩样端面与试验机压头的摩擦作用抑制了两端部的侧向变形,试样中部在较大的侧向张力作用下形成了贯穿层理面的张拉破坏,而该破坏又使岩样沿层理面开裂为平行的几部分。
三轴压缩时,页岩的破裂模式主要为剪切破坏,且随着围压的升高,破裂面的数量逐渐减少,脆性破裂特征逐渐减弱。其破裂模式可分为:
(1)单剪切面破坏。破坏的岩样均有一宏观主剪切面,且该剪切面基本都贯穿试样两端面,但随层理角度的变化,剪切面的平整度也有所变化。不同围压下,30°岩样剪切破裂面均为层理面,较平整;60°岩样剪切破裂面有明显的弯曲,且有与主剪切面不相交的层理开裂现象;90°岩样剪切破裂面有一定程度的弯曲,但部分开裂的层理面与主剪切面相交。
(2)共轭剪切破坏。岩样破坏后,有两个以上的多剪切破裂面,且大致形成两组相互平行的剪切面,该两组破裂面的交叉贯穿将岩样分为较多的块体,从而形成共轭剪切破裂面。0°试样在围压10~30 MPa 时均为此破裂模式。
对比不同层理角度页岩的强度特征和破裂模式,可知:层理角度30°左右时,页岩发生沿层理的剪切滑移破坏,强度较低,这表明页岩的层理面为地层中的薄弱面,是页岩力学特性、强度特征和破裂模式各向异性的根源。沿层理面的剪切滑移破坏是页岩地层井壁易失稳的主要原因之一。而在水力压裂过程中,层理面过弱时,压裂液易沿层理进入储层,而首先压开地层中的层理面,难以形成裂缝网络,达不到良好的压裂效果。
通过对不同层理角度页岩的破裂面与层理面及加载方向的关系(见表5)进行分析,可得出不同层理角度页岩破裂的主控因素(见表6)。
单轴压缩时,其破坏机制分为4 种类型,表现出了明显的各向异性。0°页岩为层理弱面主控的沿层理的张拉劈裂破坏;30°为层理弱面主控的沿层理的剪切滑移破坏;60°为基质体和层理弱面共同控制的贯穿层理和沿层理的剪切破坏;90°为基质体和层理弱面共同控制的贯穿层理的张拉破坏。无论哪种破坏机制,层理弱面均起到了重要作用。因此,总体上看,层理弱面的存在是引起页岩单轴压缩破坏机制各向异性的主要原因。
三轴压缩时,其破坏机制分为3 种类型,也表现出了较强的各向异性。0°页岩为基质体主控的共轭剪切破坏;30°为层理弱面主控的沿层理的剪切滑移破坏;60°和90°为基质体主控的贯穿层理的剪切破坏。对比单轴压缩时的破坏机制可知,三轴压缩时,层理弱面对破坏机制的影响已显著减小,而层理面倾角和围压效应对破坏机制的影响较大。
表5 页岩破裂模式与层理倾角和围压的关系Table 5 Relationships of failure mode,bedding orientation with confining pressure
表6 不同层理角度页岩破裂的主控因素Table 6 Dominate factors controlling the failure of shale specimens with different orientations
总之,产生页岩单轴和三轴压缩破坏机制各向异性的根源为页岩的层状沉积结构和层理面间的弱胶结作用,且其压缩强度的各向异性特征是由破坏机制的各向异性控制的。
页岩这种在平行层理和垂直层理表现出明显不同力学效应的岩体,在工程上通常被视为横观各向同性体。如图8 所示,在整体坐标系xyz 中,其应力-应变关系为[21]:
式中:ε为应变张量;S为柔度矩阵;σ为应力张量。
矩阵形式表示为
式(7)中共有5个独立的弹性参数,E1为横观各向同性面内的弹性模量;E2为垂直于横观各向同性面的弹性模量;ν1、ν2为泊松比,ν1为y 方向单独作用正应力而无其他应力分量时,x 方向应变与y 方向应变之比的负值,ν2为z 方向单独作用正应力而无其他应力分量时,x 方向应变与z 方向应变之比的负值;G12为垂直于横观各向同性面内的剪切模量。
图8 整体坐标系内的横观各向同性材料示意图Fig.8 Schematic diagram of transversely isotropic material in global coordinate system
要确定可视为横观各向同性材料的页岩的5个独立材料参数,除平行和垂直层理两方向的试样外,至少还需一个任意方向,而本文选择加载方向与层理呈30°和60°角的试样。加载方向与层理面呈β 角的弹性模量[21]为
式中:Eβ为加载方向与层理呈β 角时的弹性模量。
根据试验结果及式(8)的曲线拟合知,页岩的5个独立材料参数如表7 所示。
表7 页岩各向异性材料参数汇总表Table 7 Summary of anisotropic parameters of shale
在确定了页岩地层的各向异性材料参数后,就可以进一步分析层理性页岩地层的地应力、水平井井壁围岩应力状态、井壁稳定性、地层破裂压力及水力裂缝起裂及扩展规律等。虽然页岩各向异性问题较复杂繁琐,但国内外学者对该问题已有了足够重视,且已有了一定的研究成果,这可为下一步分析页岩水力压裂的相关问题提供参考。
(1)龙马溪组页岩地层的力学参数表现出了明显的各向异性特征。平行层理方向弹性模量最大,垂直层理方向最小;随着围压的增加,同一角度页岩弹性模量的增加速率逐渐减小。0°、30°和60°、90°页岩的泊松比随围压的增加呈现出了相反的变化规律,这可能是由层理间孔隙和微裂缝的良好发育引起的。
(2)不同围压下,0°页岩的强度最高,90°次之,30°最低,总体上呈现出两边高、中间低的U 型变化规律。不同角度的Hoek-Brown 强度准则参数也大致呈现了U 型变化规律,能较好地反映页岩强度的各向异性特征。
(3)页岩破裂模式的各向异性与层理倾角和围压的大小密切相关。破裂模式的各向异性是由破裂机制的各向异性引起的,而强度的各向异性是由破裂机制的各向异性控制的。单轴压缩时,0°页岩为沿层理的张拉劈裂破坏;30°为沿层理的剪切滑移破坏;60°为贯穿层理和沿层理的剪切破坏;90°为贯穿层理的张拉劈裂破坏。三轴压缩时,0°页岩为贯穿层理的共轭剪切破坏;30°为沿层理的剪切滑移破坏;60°和90°为贯穿层理的剪切破坏。
(4)层理面为页岩地层的薄弱面,其层状沉积结构和层间的弱胶结作用是造成力学特性、强度特征和破裂模式各向异性的主要原因。沿层理面的剪切滑移是页岩地层井壁易失稳的重要原因之一。在水力压裂过程中,层理面过弱时,压裂液易沿层理进入储层,而首先压开地层中的层理面,难以形成裂缝网络,达不到良好的压裂效果。
[1]肖钢,唐颖.页岩气及其勘探开发[M].北京:高等教育出版社,2012.
[2]张金川,金之钧,袁明生.页岩气成藏机制及分布[J].天然气工业,2004,24(7):15-18.ZHANG Jin-chuan,JIN Zhi-jun,YUAN Ming-sheng.Shale gas forming mechanism and its distribution[J].Natural Gas Industry,2004,24(7):15-18.
[3]唐颖,唐玄,王广源,等.页岩气开发水力压裂技术综述[J].地质通报,2011,30(2-3):393-399.TANG Ying,TANG Xuan,WANG Guang-yuan,et al.Summary of hydraulic fracturing technology in shale gas development[J].Geological Bulletin of China,2011,30(2-3):393-399.
[4]邹才能,董大忠,王社较,等.中国页岩气形成机制、地质特征及资源潜力[J].石油勘探与开发,2010,37(6):641-653.ZOU Cai-neng,DONG Da-zhong,WANG She-jiao,et al.Geological characteristics,formation mechanism and resource potential of shale gas in China[J].Petroleum Exploration and Development,2010,37(6):641-653.
[5]赵凯,袁俊亮,邓金根,等.层理产状对页岩气水平井井壁稳定性的影响[J].科学技术与工程,2013,13(3):580-583.ZHAO Kai,YUAN Jun-liang,DENG Jin-gen,et al.Effect of bedding plane occurrence on horizontal shale gas wellbore stability[J].Science Technology and Engineering,2013,13(3):580-583.
[6]SOLIMAN M Y,HALLIBURTON S.Interpretation of pressure behavior of fractured,deviated,and horizontal wells[J].SPE,1990,10:7-15.
[7]JAEGER J C.Shear failure of anisotropic rock[J].Geological Magazine,1960,97(1):65-72.
[8]TALIERCIO A,LANDRIANI G S.A failure condition for layered rock[J].International Journal of Rock Mechanics and Mining Sciences &Geomechanics Abstracts,1988,25(5):299-305.
[9]RAMAMURTHY T.Strength,modulus responses of anisotropic rocks[J].Compressive Rock Engineering,1993,(1):313-329.
[10]TIEN Y M,KUO M C.A failure criterion for transversely isotropic rocks[J].International Journal of Rock Mechanics and Mining Sciences,2001,38(3):399-412.
[11]NASSERI M H B,RAO K S,RAMAMURTHY T.Anisotropic strength and deformational behavior of Himalayan schists[J].International Journal of Rock Mechanics and Mining Sciences,2003,40(1):3-23.
[12]JUNG-WOO CHO,HANNA KIM,SEOKWON JEON,et al.Deformation and strength anisotropy of Asan gneiss,Boryeong shale,and Yeoncheon schist[J].International Journal of Rock Mechanics and Mining Sciences,2012,50:158-169.
[13]NIANDOU H,SHAO J F,HENRY J P,et al.Laboratory investigation of the mechanical behaviour of Tournemire shale[J].International Journal of Rock Mechanics and Mining Sciences,1997,34(1):3-16.
[14]赵文瑞.泥质粉砂岩各向异性强度特征[J].岩土工程学报,1984,6(1):32-36.ZHAO Wen-rui.Strength properties of anisotropic rock of an argillaceous siltstone[J]Chinese Journal of Geotechnical Engineering,1984,6(1):32-36.
[15]高春玉,徐进,李忠洪,等.雪峰山隧道砂板岩各向异性力学特性的试验研究[J].岩土力学,2011,32(5):1360-1364.GAO Chun-yu,XU Jin,LI Zhong-hong,et al.Experimental study of anisotropically mechanical characteristics of sandy slate in Xuefeng mountain tunnel[J].Rock and Soil Mechanics,2011,32(5):1360-1364.
[16]刘运思,傅鹤林,饶军应,等.不同层理方位影响下板岩各向异性巴西圆盘劈裂试验研究[J].岩石力学与工程学报,2012,31(4):785-791.LIU Yun-si,FU He-lin,RAO Jun-ying,et al.Research on Brazilian disc splitting tests for anisotropy of slate under influence of different bedding orientations[J].Chinese Journal of Rock Mechanics and Engineering,2012,31(4):785-791.
[17]刘胜利,陈善雄,余飞,等.绿泥石片岩各向异性特性研究[J].岩土力学,2012,33(12):3616-3623.LIU Sheng-li,CHEN Shan-xiong,YU Fei,et al.Anisotropic properties study of chlorite schist[J].Rock and Soil Mechanics,2012,33(12):3616-3623.
[18]刘斌,席道瑛,葛宁洁,等.不同围压下岩石中泊松比的各向异性[J].地球物理学报,2002,45(6):880-890.LIU Bin,XI Dao-ying,GE Ning-jie,et al.Anisotropy of Poisson’s ratio in rock samples under confining pressures[J].Chinese Journal of Geophysics,2002,45(6):880-890.
[19]宋建波,张倬元,于远忠,等.岩体经验强度准则及其在地质工程中的应用[M].北京:地质出版社,2002.
[20]SAROGLOU H,TSIAMBAOS G.A modified Hoek-Brown failure criterion for anisotropic intact rock[J].International Journal of Rock Mechanics and Mining Sciences,2008,45(2):223-234.
[21]沈观林,胡更开.复合材料力学[M].北京:清华大学出版社,2006.