岩石冻融损伤细观特征数值模拟研究*

王中文,徐 颖,,谢守冬,焦杨浩楠,于美鲁,谢昊天

(1.安徽理工大学 深部煤矿采动响应与灾害防控国家重点实验室,安徽 淮南 232001;2.安徽理工大学 土木建筑学院,安徽 淮南 232001;3.宏大爆破工程集团有限责任公司,广东 广州 510623)

0 引言

近年来,随着西部寒区露天煤矿资源开发,所面临的冻融灾害问题日益突出,长期冻融循环作用下会引起岩体结构劣化对矿产资源的安全开采造成严重威胁,如边坡和围岩稳定性失效等[1-2]。究其本质是细观冻融损伤引起的力学性能的退化[3-4],最终导致岩体丧失承载能力,引发工程灾害。因此,以细观损伤力学为基础,建立冻融损伤模型,通过数值试验研究岩石内部颗粒细观结构损伤的程度,从细观损伤层面对冻融损伤机理进行可视化分析,对研究岩石冻融损伤特性有重要意义。

针对岩石细观冻融损伤对力学性能的影响,国内外学者们做了大量研究。在理论试验方面,Hsieh等[5]建立细观颗粒组含量与宏观力学特性之间的关系。肖鹏等[6]利用SEM和三轴试验基于细观损伤理论建立本构模型。Park等[7]指出孔隙水频繁冻融是岩石受损重要原因。Khanlari等[8]分析冻融循环后砂岩的强度变化特征,得到岩石孔隙率与其冻融劣化特征之间的关系。汪鑫等[9]基于岩石细观损伤演化特征,揭示冻融损伤劣化的力学机理。

在数值模拟方面,岩石材料在冻融过程中力学变化一直是研究重点。如Zhou等[10]利用PFC研究冻融循环对土石混合体颗粒强度损伤影响。Lin等[11]建立冻融过程中水分迁移的三维模型,以研究低温环境对岩石的影响。Neaupane等[12]将线性应力-应变理论应用在有限元模型,采用热-力-流耦合来模拟冻融循环对岩石力学性能的劣化和岩石的破坏。Feng等[13]通过PFC软件建立宏观参数与细观参数之间的关系。

前人研究为本文通过0,30,80,110次冻融循环室内试验,结合颗粒流数值模拟程序,建立冻融损伤细观力学模型提供理论支持。从颗粒间的接触力、接触量、位移和微裂纹发展过程等细观层面对砂岩冻融损伤机理进行分析,可为分析冻融损伤机理提供一些参考。

1 岩石冻融模拟方法

PFC热分析模块主要针对热传导问题,通过改变颗粒尺寸引入热应变,进而在力学接触模型中引入热应力。

1.1 模型建立

PFC能够模拟由热效应导致变形和力的发展过程。PFC内置接触模型中,只有线性平行黏结模型的黏结键能够考虑热膨胀,同时传递力矩,故本文以线性平行黏结模型作为接触模型进行数值模拟试验。

将饱水岩石简化为岩石颗粒和孔隙水颗粒,3种接触类型:岩石颗粒间的接触、岩石颗粒与水颗粒接触、水颗粒与水颗粒接触。为避免岩样的尺寸效应造成影响,模型尺寸为Φ50 mm×100 mm与砂岩试样尺寸保持一致,然后对文献[14]通过CT所确定的颗粒密度与粒径范围进行参考,并计算模型的分辨率(RES),使其对数值模型的宏观力学参数不产生影响,PFC程序根据颗粒密度与粒径范围生成具体的颗粒数目。分辨率(RES)如式(1)所示:

(1)

式中:R为模型计算直径,mm;Rmax为最大颗粒半径,mm;Rmin为最小颗粒半径,mm。

最终相关参数,如表1所示。

表1 冻融砂岩颗粒参数Table 1 Parameters of freeze-thaw sandstone particles

综上所述,岩石数值模型及细观结构,如图1所示。

图1 数值模型Fig.1 Numerical model

岩石内水冰相变产生的冻胀力是引起损伤的起因[15],则建立砂岩冻融循环模型的核心内容是通过水颗粒的膨胀力等效为冻胀力。颗粒膨胀过程的实现具体如式(2)所示:

ΔR=αRΔT

(2)

式中:α为线性热膨胀系数,(1/℃)10-4;R为颗粒半径,mm;ΔT为温度增量,℃;ΔR为颗粒半径增量,mm。

α是1个微观属性,但可以由连续固体材料的宏观线性热膨胀系数αt来设定,如式(3)所示:

(3)

(4)

1.2 细观参数标定

采用颗粒离散元方法进行相关研究时,通过确定与宏观物理参数相对应的细观参数,建立两者联系,本文将冻融模型参数分为2个步骤:

1)对未冻融砂岩模型进行标定。包括:线性部分有效模量(Em)、平行黏结部分的有效模量(Pbe)、法向与切向刚度比(Pbk)、法向强度(Pbt)、切向强度(Pbc)。最终获得的参数如表2所示。

表2 未冻融砂岩细观参数Table 2 Mesoscopic parameters of sandstone without freeze-thaw

2)确定冻融循环所需的参数。颗粒线性热膨胀系数(α)、热接触单位长度热阻(η)、颗粒恒定体积比热容(Cv)。热膨胀系数参考不可恢复变形占总变形比例(0.115～0.24)[16],如表3所示。

表3 砂岩冻融参数Table 3 Freeze-thaw parameters of sandstone

1.3 数值模型验证

将进行不同冻融循环次数的试件在MTS-816型电液伺服试验机试验,控制速率为0.001 mm/s,应力-应变曲线如图2所示。

图2 不同冻融循环次数下应力-应变曲线Fig.2 Stress-strain curves under different numbers of freeze-thaw cycles

为探究冻融循环下,冻融损伤引起的力学性能的劣化原因,进行离散元模拟时,将峰值应力和应变作为标定吻合的主要因素,由表2～3所示的细观参数,得到模拟结果,与室内试验对比,如图3所示。在细观参数不变的情况下,进行不同冻融循环次数的模拟试验,得到的峰值应力与峰值应变均与室内试验大致吻合,则离散元模拟结果可用来研究砂岩细观冻融损伤所引起的力学性能的劣化。

图3 不同冻融次数试验与模拟应力-应变曲线对比Fig.3 Comparison of experimental and simulated stress-strain curves under different numbers of freeze-thaw cycles

数值模拟在不同冻融循环次数下,得到的峰值强度分别为60.17,58.14,48.57,42.17 MPa;峰值应变分别为0.763%,0.768%,0.826%,0.838%;与表4所示的室内试验结果相比,数值模拟试验结果拟合较好,如图4所示。

图4 峰值应力和峰值应变模拟与试验结果对比Fig.4 Comparison on simulated and experimental results of peak stress and peak strain

表4 不同冻融循环次数下力学参数Table 4 Mechanical parameters under different numbers of freeze-thaw cycles

2 冻融损伤岩石细观特征模拟分析

基于冻融模型,从颗粒间的接触力、接触量、位移场、微裂纹等细观特征来分析砂岩在冻融循环的作用下损伤原因及程度。

2.1 冻融砂岩颗粒接触力分析

在PFC模型中,颗粒之间由于相互挤压和剪切而产生接触力[17]。当接触力大于黏结键的拉伸或剪切强度时,黏结键就会断裂,导致接触量的减少,记录颗粒间的接触力变化可以从细观特征方面定量地分析砂岩在冻融循环过程损伤程度,如图5所示。

图5 不同冻融循环次数颗粒间的接触力Fig.5 Contact force between particles under different numbers of freeze-thaw cycles

为定量分析砂岩冻融劣化程度,在不同冻融循环次数下岩石颗粒接触力如表5所示。冻融循环次数与接触力之间的关系如图6所示。

表5 不同冻融循环下颗粒之间的接触力Table 5 Contact force between particles under different numbers of freeze-thaw cycles

由表5可知,最大接触力由5.795 N减少为5.740 N,减小比例为0.94%,接触量由140 524个减少为140 322个,减小比例为0.144%。这是因为水颗粒低温膨胀,岩石颗粒受到挤压和剪切,岩石颗粒之间黏结逐渐被破坏,而接触面积的变化会继续影响岩石颗粒之间的接触力和强度的变化,从而导致岩石内部细观结构的变化和力学性质的劣化。如表4所示,岩石的峰值强度减小,峰值应变的增大。

2.2 冻融砂岩颗粒位移场分布特征分析

颗粒的位移和位置变化可以反映冻融循环作用下岩石的变形[18]。不同冻融循环次数下的颗粒位移最大值如表6所示。拟合得到二者的函数关系如图7所示。

图7 冻融循环次数与颗粒位移关系Fig.7 Relationship between numbers of freeze-thaw cycles and particle displacement

表6 不同冻融循环次数下颗粒的最大位移Table 6 Maximum displacement of particles under different numbers of freeze-thaw cycles

由表6和图7可知,随着冻融循环次数的增加,岩石颗粒间的最大位移值增大。孔隙水的冻结膨胀和裂缝的形成使得砂岩试样发生变形,随着冻融循环次数的增加,会造成颗粒位移的增加。不同冻融循环次数颗粒的位移,如图8所示。

图8 不同冻融循环次数颗粒间的位移Fig.8 Displacement between particles under different numbers of freeze-thaw cycles

2.3 冻融砂岩颗粒裂隙演化特征

岩石颗粒间的黏结键断裂,会导致微裂纹的出现[19]。颗粒之间的黏结键断裂之后,模型内部的应力会进行重新分布,进而使微裂纹之间发生连通、扩展,最终导致岩石的宏观破坏,如图9所示。

图9 不同冻融循环次数裂隙数量及类型Fig.9 Numbers and types of cracks under different numbers of freeze-thaw cycles

由图9可知,裂纹主要集中于试样边界区域附近,这是因为靠近岩样中心的岩石受到外层岩石的约束,而靠近样品表面的岩石受到外层岩石的约束很少或没有。因此,裂纹更有可能发生在试样表面而不是试样内部。记录不同冻融循环次数下砂岩的裂纹数量、类型,如表7所示。

表7 不同冻融循环下裂纹数量Table 7 Numbers of cracks under different numbers of freeze-thaw cycles

由表7可知,随着冻融循环次数的增加,拉伸裂纹的数量占比分别为76.67%,71.58%,67.51%,表明冻融循环过程中岩石试样主要是由于拉伸破坏引起的,这是因为水颗粒在低温下变为冰颗粒挤压着周围岩石颗粒,而岩石颗粒对冰颗粒具有约束作用,冰颗粒承受压力,岩石颗粒更容易承受拉力,从而产生拉伸裂纹,但剪切裂纹占比在逐渐增多。这是因为冻融循环过程中颗粒尺寸的变化在力学接触模型中引入热应力,这些应力会引起裂纹的形成、扩展,随着冻融次数的增加,砂岩模型内部的损伤程度逐渐变大,内部应力也随之变化,岩石的表面和内部之间的应力差异会导致剪切裂纹的形成和扩展,同时在水冰膨胀与收缩不断相互作用下拉伸裂纹形成和扩展受到抑制,最终导致模型内剪切裂纹的增加。

岩体在冻融循环作用下失稳破坏实质是微裂纹的形成和扩展导致力学性能劣化。为更深入了解裂纹的形成和扩展趋势,引入裂纹倾向进行研究,裂纹倾向为裂纹在平面上的投影所指的方向,用360°方位角表示。如图10所示,由图10(a)可知,在冻融循环次数较少时,拉伸和剪切裂纹的倾向相对集中一致,皆在0°～180°,随着冻融循环次数增多,不同倾向上的裂纹数量均不同。结合上文接触力和位移的矢量图分析,可知某一倾向内颗粒的接触力越小位移越大,裂纹的数量越多,而在相同倾向上,裂纹数量反映着在冻融过程中应力累积量以及宏观裂隙的扩张方向和类型的可能性大小。

3 结论

1)从细观角度定量分析冻融循环对砂岩试样的损伤机理。通过可视化的位移图和破坏模式图,发现水颗粒的膨胀使岩石颗粒间的接触恶化,导致试样中出现微裂纹,是砂岩力学性能下降的原因。

2)随着冻融循环次数的增加,颗粒间的接触量及最大接触力在逐渐的减小。最大接触力由5.79 N减少为5.74 N,最大接触力降低0.94%;接触量由140 524个减少为140 322个,接触量减小0.14%。同时颗粒的位置随着冻融循环次数在变化,试样尺寸会略有增加。

3)冻融循环过程中,试样上会产生微裂纹。微裂纹的数量会随冻融循环的次数不断增加,相比中心区域,边界区域微裂纹更容易发生且数量在不断增加。并以拉伸裂纹为主,但随着冻融次数的增加,剪切裂纹的占比在逐渐增大,分别达到23.33%,28.42%,32.49%。