岩石化学风化力学特性的宏微观定量关联

摘要

为了分析化学风化对岩石力学性能的影响机制，采用密集排列的颗粒集合体模拟岩石材料，将与化学风化时间相关的质量损失率引入颗粒刚度和胶结模型。基于均匀化方法和格构模型，推导获得了风化岩石弹性阶段应力‒应变关系和强度准则以及宏观弹性和强度参数与微观参数（颗粒接触响应特性、风化损伤参数等）的定量关联，并探究了化学风化及其时效性对岩石力学性能的影响规律。结果表明，随着质量损失率的增大以及反应时间的增加，风化岩石的刚度和强度均呈负指数减小；当环境pH值较低时，岩石的峰值黏聚力和峰值内摩擦角随着反应时间的增大均显著降低；当环境pH值较高时，开始阶段岩石的峰值黏聚力显著降低，而峰值内摩擦角变化不大，峰值内摩擦角的显著降低需要更长的风化时间。

关键词

岩石; 化学风化; 宏微观; 时效效应

岩石作为水利水电、交通线路工程中最基础、最重要的岩土介质，其力学性质及演化直接影响工程安全和稳定性^［1］。处于自然环境中的岩石常受到化学风化作用，导致以风化岩石为基础的工程项目容易发生地陷、大沉降等地质灾害，给项目安全带来巨大挑战^［2-3］。微观层面上，化学风化将削弱晶粒之间的胶结或腐蚀矿物颗粒^［4］，岩石力学特性受风化因素和时间的控制。因此，从岩石微观结构演变入手对宏观特性进行研究，可更好地揭示风化岩石宏观力学特性的劣化机制。

借助扫描电子显微镜（SEM）和X射线CT图像扫描（X-ray CT）等先进微观观测技术，研究人员针对各地的风化岩石开展了大量的微观观察试验，探究了岩石风化过程中矿物颗粒间胶结在水作用下的破坏情况^［5］以及微观裂隙的产生和扩展规律^［6］。Li等^［7］、Lin等^［8］、Zhang等^［9］从宏微观角度分析了化学腐蚀对灰岩和砂岩试件孔隙率、孔隙结构、破裂过程、变形及强度的影响，表明水化学溶液的侵蚀对岩石的损伤破坏具有重要影响。室内微观观测试验一定程度上揭示了岩石宏观力学性质劣化的微观机理。然而，室内试验存在制样复杂、耗时长等缺点，更重要的是不能对微观结构变化进行连续观测和统计分析。

由Cundall等^［10］提出的离散元法（DEM）将岩土体视为一定排列方式下的颗粒集合体，通过尽可能真实地描述颗粒间的力学响应来再现宏观力学行为，使材料微观裂纹、应力分布和其他变量的研究和可视化变得更加容易且详细^［11］。Jiang等^［12］将考虑化学风化的岩石粒间接触模型引入离散元软件，探究不同风化程度下岩石的微观胶结破坏过程，指出风化程度与微观胶结的破坏形式密切相关并导致岩石宏观上表现出不同的破坏类型。Li等^［13］将灰岩浸泡在不同pH值的硫酸溶液中以模拟化学腐蚀试验，通过核磁共振技术总结了岩石化学腐蚀的微观损伤法则，并通过改变PFC2D内置平行黏结模型中所有接触的半径乘子模拟了岩石微细观结构的损伤演化过程。

理论模型方面，目前化学风化对于岩石试样力学特性的弱化作用大多通过基于唯象的损伤力学模型来描述^［14-15］，这些模型一定程度上考虑了岩石细观结构变化对宏观力学特性的影响规律，然而损伤变量与风化因素间的关系存在人为假定、缺乏物理意义的情况^［16］，导致考虑化学风化的岩石宏微观关联多限于定性分析阶段。另一方面，已有的研究大多只关注风化因素对某种具体岩石宏观力学行为的定性影响规律，而未对比不同孔隙率下不同岩石类型风化因素影响规律间的差异，岩石风化过程中抗剪强度参数的演化规律也未获得一致性结论^［17-20］。

在颗粒材料细观力学模型的框架下，通过引入能够反映岩石化学风化的微观胶结力学模型，建立考虑岩石化学风化的细观力学模型，推导适用于岩石的宏微观定量关联解答，用以直接探究不同孔隙率下化学风化对岩石宏观力学行为的影响规律，进一步对岩石化学风化的长期时效效应进行分析。

1 岩石材料的等效格构模型

从离散元法的基本原理出发，考虑到本文研究内容为风化对宏观特性的影响规律而非颗粒排布方式，因此将岩石概化为由同一粒径圆盘形颗粒以紧密、规则的排列方式形成的胶结颗粒集合体（见图1a）。模型中颗粒代表岩石的微观晶粒，而晶粒间的矿物黏结作用通过能够反映岩石化学风化的微观胶结力学模型描述。不同于离散元法中将连续体材料直接离散化为颗粒体间的相互作用，格构模型将组成材料的基本单元（如岩石中的矿物颗粒）视为质点，将相邻点连接形成格构网络，并基于铁木辛柯梁定义基本单元间的力学特性，该模型被广泛应用于岩石的裂隙扩展、断裂等破坏问题研究^［21］。由于晶粒间的矿物黏结特性与铁木辛柯梁力学响应具有等效性^［22］，因此如图1a所示的岩石概化颗粒模型可以等效为如图1b所示的格构模型，从而将颗粒集合体力学问题简化为杆系力学问题。格构模型与真实岩石的等效性证明见文献［23-24］。

图1  颗粒模型和等效格构模型示意图

Fig.1  Schematic diagram of particle assembly and equivalent lattice model

如图1b所示，不同于离散元法中直接针对大量的离散颗粒体进行显式迭代求解，受连续介质力学中的应力单元启发，本文通过对代表性单元（RVE）进行力学分析以表征材料的宏观力学特性，从而可以方便地获得材料的宏微观定量关联解答。显然，该力学模型的宏观特性将与RVE的结构形式和梁的力学特性（即微观特性）相关。本文的RVE由如图1b所示的4个颗粒所对应的5根格构梁组成，包括上三角和下三角两部分，称之为胞元1（UC-1）和胞元2（UC-2）。

在颗粒排列方式确定后，RVE力学响应主要受控于梁单元的力学特性。考虑到本文的研究对象为风化岩石，因此梁单元力学特性采用考虑化学风化的二维粒间接触模型^［25］，如图2所示。图2中，W、h_min和h_max分别为胶结宽度、最小和最大胶结厚度，R为颗粒半径，D为颗粒直径。根据如图2所示的几何关系，岩石的孔隙率n与W、h_min和R具有如下关系：

$$n=\mathrm{1}-\frac{\mathrm{\pi }{R}^{\mathrm{2}}+\mathrm{3}\left(W(\mathrm{2}R+h_{min}-\sqrt[]{\mathrm{4}{R}^{\mathrm{2}}-W^{\mathrm{2}}})-\mathrm{2}\left(R^{\mathrm{2}}\mathrm{a}\mathrm{r}\mathrm{c}\mathrm{s}\mathrm{i}\mathrm{n}\frac{W}{\mathrm{2}R}-\frac{W}{\mathrm{2}}\sqrt[]{R^{\mathrm{2}}-\frac{W^{\mathrm{2}}}{\mathrm{4}}}\right)\right)}{\mathrm{2}\sqrt[]{\mathrm{3}}(R+h_{\mathrm{m}\mathrm{i}\mathrm{n}}{/\mathrm{2})}^{\mathrm{2}}}$$

（1）

其中，$W=\sqrt[]{\mathrm{4}{R}^{\mathrm{2}}-(\mathrm{2}R+h_{\mathrm{m}\mathrm{i}\mathrm{n}}-h_{\mathrm{m}\mathrm{a}\mathrm{x}}{)}^{\mathrm{2}}}$。通过式（1）可探究孔隙率对风化岩石力学特性的影响规律。接触模型在弹性阶段满足以下线性力‒位移关系：

$$F_{\mathrm{n}}=K_{\mathrm{n}}{x}_{\mathrm{n}}$$， $$F_{\mathrm{v}}=K_{\mathrm{v}}{x}_{\mathrm{v}}$$， $$F_{\mathrm{m}}=K_{\mathrm{m}}{x}_{\mathrm{m}}$$

（2）

式中：F_i、K_i、x_i（i=n，v，m）分别为沿i方向的接触力/力矩、接触刚度以及相对位移/转角；n、v、m分别表示法向、切向、弯转向。接触刚度K_i包括颗粒法向刚度K_np和胶结法向刚度K_nb两部分的贡献。考虑到胶结具有两端向内凹陷的不规则形状特点，胶结法向刚度K_nb通过在整个胶结宽度W上积分得到：

$$K_{\mathrm{n}\mathrm{b}}=\mathrm{2}{\int }_{\mathrm{0}}^{W/\mathrm{2}}\frac{E_{\mathrm{0}}}{h_i}\mathrm{d}{W}_i=E_{\mathrm{0}}{\int }_{\mathrm{0}}^{W/\mathrm{2}}\frac{\mathrm{1}}{(-\sqrt[]{R^{\mathrm{2}}-W_i^{\mathrm{2}}}+R+h_{\mathrm{m}\mathrm{i}\mathrm{n}}/\mathrm{2})}\mathrm{d}{W}_i$$

（3）

式中，E₀为胶结弹性模量。粒间接触的切向刚度K_v和弯转向刚度K_m通过切向‒法向刚度比λ=K_v/K_n和弯转向‒法向刚度比ξ=K_m/K_n定义。

图2  胶结接触模型示意图

Fig.2  Schematic diagram of bonded contact model

根据室内胶结铝棒单轴拉伸、压缩以及压剪扭复杂加载路径试验，结合前人提出的胶结强度包线表达式^［25］，定义如下简化的峰值拉伸、压缩、剪切和扭转荷载以及复合压‒剪‒扭破坏包线表达式：

$$F_{\mathrm{t}}={\sigma }_{\mathrm{t}}WL{\zeta }^{\mathrm{0.474}}{\mathrm{e}}^{\mathrm{0.794}\beta }$$

（4）

$$F_{\mathrm{c}}={\sigma }_{\mathrm{c}}WL{\zeta }^{-\mathrm{0.58}}{\mathrm{e}}^{-\mathrm{0.611}\beta }$$

（5）

$$F_{\mathrm{s}}=\mathrm{0.34}{\zeta }^{\mathrm{0.42}}{\mathrm{e}}^{-\mathrm{0.611}\beta }{F}_{\mathrm{c}}(f+\mathrm{2}{\mathrm{e}}^{-\mathrm{1}}\sqrt[]{\mathrm{1}-\mathrm{4}{(f-\mathrm{0.5})}^{\mathrm{2}}})$$

（6）

$$F_{\mathrm{r}}=\mathrm{0.1}W{F}_{\mathrm{c}}(f+\mathrm{2}{\mathrm{e}}^{-\mathrm{1}}\sqrt[]{\mathrm{1}-\mathrm{4}{(f-\mathrm{0.5})}^{\mathrm{2}}})$$

（7）

$$\frac{F_{\mathrm{v}}^{\mathrm{2}}}{F_{\mathrm{s}}^{\mathrm{2}}}+\frac{F_{\mathrm{m}}^{\mathrm{2}}}{F_{\mathrm{r}}^{\mathrm{2}}}=\mathrm{1}$$

（8）

式中：F_t、F_c、F_s、F_r分别为胶结的峰值拉伸、压缩、剪切、扭转荷载；σ_c、σ_t分别为胶结单轴抗压和抗拉强度；ζ、β分别为胶结延伸率和长细比，ζ=W/D，β=h_min/W；f为表征法向力水平的法向应力比，f=（F_n+F_t）/（F_c+F_t）；L为z方向单位长度。

根据室内微观观测结果，化学风化将削弱晶粒之间的胶结或引起矿物颗粒腐蚀^［3］。反映到如图2所示的胶结接触模型中表现为胶结宽度W以及颗粒半径R的减小。因此根据Ciantia等^［16］的研究，将质量损失率α作为定量反映化学风化损伤程度的物理指标。α由试样在化学风化过程中损失的质量M_loss与风化前质量M₀之比表示，即α=M_loss/M₀。根据如图2所示的胶结模型示意图，胶结物质的质量M与面积S有关。在h_min一定的情况下，胶结宽度W与胶结面积S存在非线性关系，因此假定胶结宽度和颗粒刚度与质量损失率α存在以下关系：

$$W=W_{\mathrm{0}}{\left(\mathrm{1}-\frac{\alpha }{{\alpha }_{\mathrm{c}\mathrm{r}\mathrm{i}}}\right)}^{\mathrm{2}}$$， $$K_{\mathrm{n}\mathrm{p}}=K_{\mathrm{n}\mathrm{p},\mathrm{0}}{\left(\mathrm{1}-\frac{\alpha }{{\alpha }_{\mathrm{c}\mathrm{r}\mathrm{i}}}\right)}^{\mathrm{2}}$$

（9）

式中：W₀、k_np，0分别为化学风化前胶结宽度和颗粒法向刚度；α_cri为临界质量损失率。质量损失率达到α_cri时试样强度很小，可认为胶结物已经全部溶解。将式（9）代入2节的弹性和强度参数宏微观定量关联即可得到考虑化学风化影响的质量损失率对岩石弹性和强度参数的定量影响规律。

2 岩石模型弹性和强度参数的宏微观定量关联

2.1　弹性参数的宏微观关联

小变形假定下，均匀化理论认为离散颗粒体内部的应力场或者应变场均匀，即著名的Reuss’s假设（静态假设）^［26］和Voigt’s假设（动态假设）^［27］。动态假设从宏观平均应变出发，将微观的粒间位移表示为应变的函数，基于微观接触模型求得微观粒间接触力，再由粒间接触力均匀化求得宏观平均应力，完成从宏观应变到宏观应力的关系建立，而静态假设采用相反的思路建立从宏观应变到宏观应力的关系。以胞元UC-1为例，根据动态假设法，将3个颗粒间的相对位移代入式（2）中的粒间本构关系可得到3个颗粒间的接触力，并根据能量平衡原理^［21］将粒间接触力均匀化为宏观应力，最终得到UC-1的宏观等效应力‒应变关系，即：

$$\left\{\begin{array}{c}{\sigma }_{\mathrm{11}}\\ {\sigma }_{\mathrm{22}}\end{array}\right\}=\frac{\sqrt[]{\mathrm{3}}}{\mathrm{4}}\left[\begin{array}{cc}\mathrm{3}{K}_{\mathrm{n}}+K_{\mathrm{v}}& K_{\mathrm{n}}-K_{\mathrm{v}}\\ K_{\mathrm{n}}-K_{\mathrm{v}}& \mathrm{3}{K}_{\mathrm{n}}+K_{\mathrm{v}}\end{array}\right]\left\{\begin{array}{c}{\epsilon }_{\mathrm{11}}\\ {\epsilon }_{\mathrm{22}}\end{array}\right\}$$

（10）

式中：${\sigma }_{\mathrm{11}}$、${\sigma }_{\mathrm{22}}$分别为x、y方向的应力；${\epsilon }_{\mathrm{11}}$、${\epsilon }_{\mathrm{22}}$分别为x、y方向的应变。以胞元UC-2为例可以得到同样的结果，具体推导过程见Chang等^［22］的研究。将式（10）与平面应力下各向同性柯西连续体模型的本构方程进行对比，可得弹性模量E和泊松比ν与粒间接触刚度之间的宏微观定量关联，即：

$$E=\frac{\mathrm{2}\sqrt[]{\mathrm{3}}{K}_{\mathrm{n}}(\mathrm{1}+\lambda )}{\mathrm{3}+\lambda }$$， $$\nu =\frac{\mathrm{1}-\lambda }{\mathrm{3}+\lambda }$$

（11）

式（11）以显式解的形式直接展示微观参数对宏观力学行为的影响规律。可见，弹性模量E为粒间法向刚度K_n的线性函数并随着刚度比λ的增大而增大。由式（3）中胶结法向刚度K_nb与胶结宽度W和胶结弹性模量E₀的关系可以推断，弹性模量E随胶结宽度W的减小而负指数减小，随着胶结弹性模量E₀的增大而线性增大。泊松比ν为刚度比λ的单值函数并随着λ的增大而非线性减小。

2.2　强度准则及强度参数的宏微观关联

为获得RVE的强度（即宏观强度），在图3中节点B/D和A/C处施加荷载F_x和F_y，通过结构力学分析得到格构梁“c”（c=1~5）跨中内力f_ic（i=n，v，m）。将内力代入式（4）~（8）定义的强度准则，即可得到梁破坏时对应的荷载F_x和F_y，并换算为应力得到RVE的强度，推导中考虑了格构梁的渐进破坏过程。如图3a中虚线所示，格构梁1首先被拉坏，将此时对应的外部荷载定义为材料的初始屈服强度，该强度类似于岩石典型应力‒应变曲线^［28］中的起裂应力，此时岩石晶粒间的胶结以拉坏为主；随后，如图3b中点划线所示，格构梁2~5被剪坏，对应的外部荷载定义为材料的峰值强度，此时RVE彻底失去承载能力。当RVE受y方向单轴拉伸应力时，如图3c所示，只存在一种破坏情况，即斜杆直接受拉剪破坏，此时对应的外部荷载定义为材料的单轴抗拉强度。

图3  RVE荷载条件以及格构梁破坏模式

Fig.3  Loading conditions of RVE and failure modes of beams

限于篇幅，本文不再展示具体的推导过程，详细推导过程和验证可参考已有成果^［24］。材料的初始屈服强度和峰值强度分别为：

$${\sigma }_{\mathrm{1}}=\frac{\sqrt[]{\mathrm{3}}}{\mathrm{2}}\frac{(\mathrm{1}-\mathrm{2}{A}_{\mathrm{c}})}{B_{\mathrm{c}}}{\sigma }_{\mathrm{3}}+\frac{\mathrm{1}}{\mathrm{2}{B}_{\mathrm{c}}}\frac{F_{\mathrm{t}}}{D}$$

（12）

$${\sigma }_{\mathrm{1},\mathrm{p}\mathrm{e}\mathrm{a}\mathrm{k}}=-\frac{\mathrm{1}}{\mathrm{2}{D}_{\mathrm{c}}}\left(\left(P_{\mathrm{c}}{\sigma }_{\mathrm{3}}+\frac{F_{\mathrm{t}}}{D}{Q}_{\mathrm{c}}\right)+\sqrt[]{{\left(P_{\mathrm{c}}{\sigma }_{\mathrm{3}}+\frac{F_{\mathrm{t}}}{D}{Q}_{\mathrm{c}}\right)}^{\mathrm{2}}-\mathrm{4}{D}_{\mathrm{c}}\left(L_{\mathrm{c}}{\sigma }_{\mathrm{3}}^{\mathrm{2}}+\frac{F_{\mathrm{t}}}{D}{M}_{\mathrm{c}}{\sigma }_{\mathrm{3}}+\frac{F_{\mathrm{t}}^{\mathrm{2}}}{D^{\mathrm{2}}}{N}_{\mathrm{c}}\right)}\right)$$

（13）

式中：σ₁、σ_1，peak分别表示材料在初始围压σ₃下的初始屈服强度和峰值强度；A_c、B_c为量纲一的初始屈服强度参数，是刚度比λ和ξ的函数；D_c、P_c、Q_c、L_c、M_c、N_c为量纲一的峰值强度参数，是胶结峰值拉压荷载比$\eta (=F_{\mathrm{t}}/F_{\mathrm{c}})$以及胶结内摩擦系数μ_b的函数。初始屈服强度包线满足Mohr-Coulomb（MC）强度准则，受刚度比λ和ξ以及微观粒间抗拉强度F_t/D的控制。峰值强度包线在主应力空间上与岩石常用的Hoek-Brown（HB）强度准则类似，由线性部分和非线性部分组成。根据MC强度准则以及HB强度准则的定义，初始屈服强度参数为：

$${\sigma }_{\mathrm{u},\mathrm{i}\mathrm{n}\mathrm{i}}=\frac{\mathrm{1}}{\mathrm{2}{B}_{\mathrm{c}}}\frac{F_{\mathrm{t}}}{D}$$

（14）

$$c_{\mathrm{i}\mathrm{n}\mathrm{i}}=\frac{(\sqrt[]{{(p-\mathrm{1})}^{\mathrm{2}}+\mathrm{4}p}-p+\mathrm{1})}{\mathrm{8}{B}_{\mathrm{c}}\sqrt[]{p}}\frac{F_{\mathrm{t}}}{D}$$

（15）

$${\phi }_{\mathrm{i}\mathrm{n}\mathrm{i}}=\mathrm{a}\mathrm{r}\mathrm{c}\mathrm{t}\mathrm{a}\mathrm{n}((p-\mathrm{1})/\mathrm{2}\sqrt[]{p})$$

（16）

式中：σ_u，ini、c_ini、φ_ini分别为初始屈服强度对应的初始单轴抗压强度、初始黏聚力和初始内摩擦角；$p=\sqrt[]{\mathrm{3}}(\mathrm{1}-\mathrm{2}{A}_{\mathrm{c}})/\mathrm{2}{B}_{\mathrm{c}}$。峰值强度参数为：

$${\sigma }_{\mathrm{T}}=\frac{Y_{\mathrm{t}}-\sqrt[]{Y_{\mathrm{t}}^{\mathrm{2}}-\mathrm{4}{X}_{\mathrm{t}}{Z}_{\mathrm{t}}}}{\mathrm{2}{X}_{\mathrm{t}}}\frac{F_{\mathrm{t}}}{D}$$

（17）

$${\sigma }_{\mathrm{u},\mathrm{p}\mathrm{e}\mathrm{a}\mathrm{k}}=\frac{-Q_{\mathrm{c}}-\sqrt[]{Q_{\mathrm{c}}^{\mathrm{2}}-\mathrm{4}{D}_{\mathrm{c}}{N}_{\mathrm{c}}}}{\mathrm{2}{D}_{\mathrm{c}}}\frac{F_{\mathrm{t}}}{D}$$

（18）

$$c_{\mathrm{p}\mathrm{e}\mathrm{a}\mathrm{k}}({\sigma }_{\mathrm{3}})=\frac{\mathrm{1}-\mathrm{s}\mathrm{i}\mathrm{n}{\phi }_{\mathrm{p}\mathrm{e}\mathrm{a}\mathrm{k}}({\sigma }_{\mathrm{3}})}{\mathrm{2}\mathrm{c}\mathrm{o}\mathrm{s}{\phi }_{\mathrm{p}\mathrm{e}\mathrm{a}\mathrm{k}}({\sigma }_{\mathrm{3}})}{\sigma }_{\mathrm{u},\mathrm{p}\mathrm{e}\mathrm{a}\mathrm{k}}$$

（19）

$${\phi }_{\mathrm{p}\mathrm{e}\mathrm{a}\mathrm{k}}({\sigma }_{\mathrm{3}})=\mathrm{a}\mathrm{r}\mathrm{c}\mathrm{s}\mathrm{i}\mathrm{n}((H_c^{}-\mathrm{1})/(H_c^{}+\mathrm{1}))$$

（20）

$$H_{\mathrm{c}}=-\frac{\mathrm{1}}{\mathrm{2}{D}_{\mathrm{c}}}\left(P_{\mathrm{c}}+\frac{P_{\mathrm{c}}(P_{\mathrm{c}}{\sigma }_{\mathrm{3}}+F_{\mathrm{t}}{Q}_{\mathrm{c}}/D)-\mathrm{2}{D}_{\mathrm{c}}(\mathrm{2}{L}_{\mathrm{c}}{\sigma }_{\mathrm{3}}+F_{\mathrm{t}}{M}_{\mathrm{c}}/D)}{\sqrt[]{(P_{\mathrm{c}}{\sigma }_{\mathrm{3}}+F_{\mathrm{t}}{Q}_{\mathrm{c}}{/D)}^{\mathrm{2}}-\mathrm{4}{D}_{\mathrm{c}}(L_{\mathrm{c}}{\sigma }_{\mathrm{3}}^{\mathrm{2}}+F_{\mathrm{t}}{M}_{\mathrm{c}}{\sigma }_{\mathrm{3}}/D+F_{\mathrm{t}}^{\mathrm{2}}{N}_{\mathrm{c}}/D^{\mathrm{2}})}}\right)$$

式中：σ_T为材料的单轴抗拉强度；σ_u，peak、c_peak、φ_peak分别为峰值强度对应的峰值单轴抗压强度、峰值黏聚力和峰值内摩擦角；X_t、Y_t、Z_t为量纲一的峰值强度参数，是胶结峰值拉压荷载比η（=F_t/F_c）以及胶结内摩擦系数μ_b的函数。从式（14）~（20）定义的宏微观定量关联可以看到：单轴抗压强度、黏聚力以及单轴抗拉强度均为胶结抗拉强度F_t/D的线性函数，而由式（4）可知，宏观单轴抗压强度、黏聚力以及单轴抗拉强度均随着胶结宽度W的增大而负指数减小，随着胶结单轴抗拉强度σ_t的增大而线性增大。

3 化学风化对砂岩宏观力学特性的影响规律

基于宏微观定量关联式，探究风化对某一典型砂岩宏观力学参数及时效性的影响规律。表1为本节所采用的微观胶结参数。根据式（1）中孔隙率与最小胶结厚度的关系，h_min取如表1所示的不同值以探究孔隙率对宏观力学特性的影响规律。表1同时展示了当h_min=1 μm时，基于本文宏微观定量解得到的风化前宏观力学参数，结果满足典型砂岩的宏观物理力学指标（表1中室内试验）。

表1  砂岩化学风化分析微观参数取值及其对应的宏观力学参数

Tab.1  Values of micro parameters and corresponding macro properties used in analysis of sandstone chemical weathering

微观参数	数值	宏观参数	数值（室内试验［29］）
颗粒直径，D/mm	1	弹性模量，E/GPa	3.64（3.60）
初始胶结宽度，W0/mm	0.3	峰值单轴抗压强度，σu，peak/MPa	68.1（63.0）
胶结弹性模量，E0/MPa	45	单轴抗拉强度，σT/MPa	4.2
颗粒法向刚度，Knp/（MN∙m-1）	100	峰值黏聚力，cpeak/MPa	17.2
切向‒法向刚度比，λ	0.1	内摩擦角，φ/（°）	45.7
弯转向‒法向刚度比，ξ/mm2	0.15	孔隙率，n/%	7.81
胶结单轴抗拉强度，σt/MPa	12
胶结单轴抗压强度，σc/MPa	150
最小胶结厚度，hmin/μm	1、5、10、20、50、100

3.1　风化对力学参数的影响规律

图4为不同孔隙率下弹性模量比β_E=E_b（α）/E_b随质量损失率α的变化，同时与室内试验^［16］结果进行对比。E_b（α）和E_b分别为质量损失率α时和未受到风化时的弹性模量。由图4可见，弹性模量比β_E随质量损失率α的增加而非线性减小。不同孔隙率下的变化规律存在明显差异：当孔隙率较小时，弹性模量比变化率先增大后减小，而当孔隙率较大时，弹性模量比变化率一直减小，并且该情况下理论解与室内试验结果在定性上吻合更好。相同质量损失率α下，孔隙率越大β_E越小。

图4  不同孔隙率下弹性模量比随质量损失率的变化

Fig.4  Variation of elastic modulus ratio with mass loss ratio under different porosities

图5为不同孔隙率下抗压强度随质量损失率α的变化。初始屈服强度比β_c，ini=σ_u，ini（α）/σ_u，ini（0）和初始黏聚力c_ini随质量损失率α的增加而非线性减小，减小的速率随质量损失率α的增加而减小；孔隙率越大，β_c，ini和c_ini越大，质量损失率的影响越明显。峰值强度比β_c，peak=σ_u，peak（α）/σ_u，peak（0）随质量损失率α增加而减小，减小的速率随α的增加而减小。相反地，相同质量损失率α下，孔隙率越大β_c，peak越小。β_c，peak理论解与室内试验结果^［16］吻合较好。

图5  不同孔隙率下抗压强度随质量损失率的变化

Fig.5  Variation of compressive strength with mass loss ratio under different porosities

图6a为不同孔隙率下围压σ₃=0时对应的峰值内摩擦角φ_peak和峰值黏聚力c_peak随质量损失率α的变化。φ_peak随α先增大后减小并收敛于30°附近，增大阶段增长率逐渐减小；孔隙率越小，质量损失率α的影响越明显；相同质量损失率α下，孔隙率越大φ_peak越小。图6b中，峰值黏聚力c_peak随质量损失率α增大而非线性减小，变化率也逐渐减小；孔隙率越小，质量损失率α的影响越明显；相同质量损失率α下，孔隙率越大c_peak越小。

图6  不同孔隙率下强度参数随质量损失率的变化

Fig.6  Variation of strength parameters with mass loss ratio under different porosities

图7为不同孔隙率下单轴抗拉强度比β_T= σ_T（α）/σ_T（0）随质量损失率α的变化以及理论解与室内试验结果对比^［19，30］。β_T随α增大而非线性减小，变化率也逐渐减小，孔隙率对β_T的影响较小。理论解与室内试验结果定性上吻合较好，但室内试验结果中质量损失率对单轴抗拉强度的影响更为显著，这是由于室内试验中试样在化学侵蚀下会产生不均匀问题，理论解以均匀性简化，因此有所偏差。

图7  不同孔隙率下单轴抗拉强度比随质量损失率的变化

Fig.7  Variation of uniaxial tensile strength ratio with mass loss ratio under different porosities

3.2　风化的时间效应分析

酸性环境下，岩石中的钙质主要与氢离子H⁺发生反应，因此氢离子浓度是反应速率的决定因素。根据岩石风化的室内试验结果^［29］，岩石中由酸性环境所导致的钙质溶解量与风化时间具有如下规律：

$$m_t=m_{\mathrm{0}}-M\gamma tV$$

（21）

式中：γ为反应速率；m_t、m₀分别为风化后t时刻和初始时刻钙物质的质量；M为钙元素摩尔质量；V为溶液体积。根据式（21）以及质量损失比α=M_loss/M₀的定义，可对岩石化学风化时效性进行分析。

图8为不同pH值下弹性及强度参数随反应时间的变化，并与室内试验结果^［20，31］进行对比。如图8a所示，弹性模量比β_E=E_b（t）/E_b及其变化率随着反应时间增大而减小，相同时间段内减小程度随pH值的减小而增大，即相同反应时间下，pH值越小试样的弹性模量比越小。理论解与室内试验结果^［31］定性上一致，但室内试验中pH值对弹性模量比的影响程度不如理论解明显，可能由于室内试验采用的是结构更为致密、原生裂隙更少的页岩，而本文时效性分析中引入的反应速率式（21）是基于钙质砂岩的室内试验结果^［29］。

图8  不同pH值下弹性及强度参数随反应时间的变化

Fig.8  Variation of elastic and strength parameters with reaction time under different pH values

图8b中，峰值强度比β_c，peak=σ_u，peak（t）/σ_{u， peak}（0）及其变化率随着反应时间的增大而减小，相同时间段内减小程度随pH值的减小而增大，即相同反应时间下pH值越小试样的峰值强度比越小。β_c，peak的理论解与室内试验结果^［31］定性上吻合较好，定量上相同风化时间内理论解中峰值强度比的减小程度更高，并且pH值对理论解的影响程度远大于室内试验结果。图8c中，峰值黏聚力c_peak随反应时间的变化规律与峰值强度比类似，理论解与室内试验结果^［20］吻合较好。

图8d为不同pH值下试样峰值内摩擦角φ_peak随反应时间的变化。可以看到，随着反应时间的增大，峰值内摩擦角先轻微增大，达到峰值不久后迅速减小，当质量损失率达到100%时，峰值内摩擦角跌落到30°附近。pH值越小，峰值内摩擦角增大的时间越短，越早达到峰值并迅速跌落，pH=7时甚至要到100 d左右的风化时间才能完成峰值内摩擦角的轻微增大过程。室内试验中^［20］峰值内摩擦角在60 d的风化时间内轻微减小，同样表现出变化不大的特点。对于pH值较低溶液中的岩石，岩石的峰值黏聚力和峰值内摩擦角随着反应时间的增大均显著降低；对于pH值较高溶液中的岩石，开始阶段岩石的峰值黏聚力显著降低，而峰值内摩擦角变化不大，峰值内摩擦角的显著降低需要更长的反应时间，该结论与针对灰岩的室内试验结果^［18］一致。

4 结论

（1） 随着化学风化过程的进行，试样弹性模量比随着质量损失率的增大呈负指数减小。风化开始时，孔隙率越大弹性模量比变化率越大；当孔隙率较小时，弹性模量比变化率先增大后减小，而当孔隙率较大时，弹性模量比变化率始终减小。相同质量损失率α下，孔隙率越大弹性模量比越小。

（2） 初始屈服强度比、峰值强度比及初始黏聚力随着质量损失率的增大呈负指数减小。相同质量损失率α下，孔隙率越大，初始屈服强度比越大，峰值强度比越小，初始黏聚力越大。峰值内摩擦角随α先增大后减小并迅速收敛于30°附近，相同质量损失率α下孔隙率越大，峰值内摩擦角越小。单轴抗拉强度及其变化率随着质量损失率的增大而非线性减小，孔隙率对单轴抗拉强度的影响不明显。

（3） 弹性模量比、峰值强度比和峰值黏聚力及其三者的变化率随着化学风化反应时间t增大呈负指数减小，相同时间段内减小程度随pH值的减小而增大。峰值内摩擦角随着化学风化反应时间t先轻微增大后迅速减小并跌落到30°附近。对于pH值较低溶液中的岩石，岩石的峰值黏聚力和峰值内摩擦角随着反应时间的增大均显著降低；对于pH值较高溶液中的岩石，开始阶段岩石的峰值黏聚力显著降低，而峰值内摩擦角变化不大，峰值内摩擦角的显著降低需要更长的反应时间。

作者贡献声明

周志豪：公式推导，参数分析，论文撰写。

王华宁：项目负责人，论文修改。

蒋明镜：项目负责人，研究思路指导，论文修改。