0.   引言

堆石料作为一种广泛应用的工程材料,其遇水后产生湿化变形一直是工程界关注的焦点问题。工程实践表明,堆石坝的上游堆石体在初次蓄水及库水位变动时会产生较大的沉降变形；或是降雨入渗、下游尾水位的抬高,也会导致下游堆石体产生明显的沉降变形^[1]。国内外堆石坝工程中,存在不少因湿化变形导致的安全事故,例如瀑布沟心墙堆石坝因上游堆石体的湿化变形导致坝顶出现纵向裂缝^[2]；降雨入渗使天生桥面板堆石坝上游堆石体变形加剧,面板脱空愈加严重^[3]。因此,湿化变形对堆石坝的安全运行影响重大,研究堆石体的湿化变形特性和细观机理十分必要。

近年来,国内外学者为研究堆石体湿化变形特性进行了大量的室内物理试验。Baud等^[4]对不同湿化程度的砂岩进行常规三轴压缩试验,认为颗粒强度降低的原因为孔隙中的承压水会使岩石变弱变脆,且孔隙水会使破碎能量降低。Ham等^[5]进行了花岗岩在干湿两种状态下的单颗粒破碎及一维压缩试验,结果表明由于水的存在,单颗粒初始破碎强度降低,强度变异性增大。魏松等^[6]对粗粒料进行了三轴湿化变形的试验研究,发现湿化变形与颗粒破碎密切相关,可以从颗粒破碎角度解释湿化变形产生的机理。丁艳辉等^[7-8]进行了堆石体在侧限条件下的非饱和湿化试验,结果表明颗粒小孔隙中的水对湿化变形起主要作用,而大孔隙中的水的影响不大,此后进一步进行了弱风化花岗岩的三轴湿化试验,认为堆石体的湿化变形是浸水导致的物态弱化变形,并将湿化变形分为湿化瞬态变形和湿态流变变形。贾宇峰等^[9-10]进行了等应力比路径下的湿化变形试验和复杂应力路径下的排水三轴剪切试验,并统计了颗粒破碎情况,结果表明应力路径、应力水平、应力强度和颗粒破碎率决定了等应力比路径下的湿化变形。遇水湿化不仅会导致颗粒强度减小和颗粒破碎。Ulusay等^[11]通过倾斜试验分别测得不同岩石在干湿状态下的摩擦系数,结果表明湿化使得岩石的摩擦系数降低。在堆石体湿化试验结果的基础上,学者们提出了不同的湿化变形模拟方法,如单线法,双线法,模量软化法等,且不断改进堆石体湿化变形的计算方法和湿化本构模型,并用于堆石坝的湿化变形计算^[12-16]。

试验研究揭示了堆石体的湿化变形特性和影响因素,与此同时,细观数值模拟也被用来验证堆石体湿化变形的细观机理。黄绪武等^[17]采用连续-离散耦合分析方法（combined finite and discrete element method,FDEM）,进行堆石体的单线法双轴压缩数值试验,结果表明湿化变形主要是颗粒间摩擦系数和颗粒强度的降低引起的。杨贵等^[18]采用三维离散元方法模拟了粗粒料在干湿条件下的三轴剪切试验,得到不同围压下试样强度的变化规律。王蕴嘉等^[19]运用二维离散元模拟粗粒料的单、双线法侧限湿化试验,研究了湿化路径、颗粒形状及软化系数对湿化变形的影响,结果表明单线法湿化应变大于双线法,低应力湿化变形主要原因为颗粒滑移,而高应力湿化变形则主要是受颗粒破碎的影响。然而该研究仅考虑了颗粒的二维形态,未能分析三维空间的颗粒形状对湿化变形的影响。Shao等^[20]运用离散元法模拟恒定围压和应力比条件下不同粒径试样的三轴湿化试验,模拟的应力–应变曲线与干湿颗粒的单轴压缩试验结果相吻合,尤其是在高应力比和围压条件下,大尺寸试样湿化变形明显大于小尺寸试样。

堆石体湿化变形问题十分复杂,目前从细观层面对湿化变形机理的研究仍有待进一步的探索。因此,本文采用连续–离散耦合分析方法,参照室内试验成果率定了干燥与饱和试样的细观参数,进而开展了单颗粒破碎和颗粒集合体单轴压缩的双线法湿化数值试验,通过比较干燥与饱和试样的加载过程与破碎规律,探讨了堆石料的湿化变形特性和机理。

1.   堆石体湿化对堆石坝变形的影响

以瀑布沟、毛尔盖两个典型的堆石坝工程为例,在堆石坝施工、蓄水和运行过程中,由于库水位的变化（初次蓄水、水位波动）,上游较高高程的部分堆石体会经历多次干湿交替。如图1所示,堆石坝测点的沉降变形增量与初次蓄水及库水位的周期变化息息相关,在首次达到正常蓄水位及运行期间水位波动上升阶段,变形增量显著增加^[21]。由此可见,堆石坝的变形与堆石体的湿化变形特性密切相关,并会在一定程度上影响堆石坝工程的稳定性与安全运行。

图  1  大坝沉降与水位变化关系图

Figure  1.  Relationship between dam subsidence and change of water level

下载: 全尺寸图片 幻灯片

2.   湿化对颗粒破碎的影响

2.1   干湿状态的单颗粒破碎试验

首先进行颗粒在干燥和饱和状态下的单颗粒压缩试验。试验选用板岩颗粒,其密度为2240 kg/m³,弹性模量为13.2 GPa,泊松比为0.2。颗粒形状以扁平状为主,局部棱角相对分明。将60个板岩颗粒放入真空干燥箱中烘干24 h。对干燥后的30个单颗粒进行单轴压缩试验,将另外30个颗粒放入水中浸泡24 h,保证颗粒充分饱和。为了防止水分蒸发,将浸泡24 h后的30个颗粒直接进行单轴压缩试验。使用MTS810试验系统进行单颗粒压缩试验（图2）,试验加载速度为0.01 mm/s,保证试验在准静态条件下进行。在加载初期,变形主要发生在颗粒和加载板接触处。由于接触处的变形,荷载–位移曲线表现为Hertz形式的接触行为。在达到峰值后,出现贯穿颗粒的大裂纹,随后荷载急剧减小直至试验结束。干燥与饱和颗粒的典型荷载–位移曲线如图3所示,饱和颗粒到达峰值荷载时的位移提前,且最大荷载值明显小于干燥颗粒。

图  2  单颗粒压缩试验

Figure  2.  Single particle compression test

下载: 全尺寸图片 幻灯片

图  3  干燥与饱和颗粒的典型荷载–位移曲线

Figure  3.  Typical load-displacement curves for dry and saturated particles

下载: 全尺寸图片 幻灯片

对干燥与饱和颗粒的表面进行SEM电镜扫描,结果表明遇水饱和会导致颗粒的岩石矿物成分和微观结构发生变化,如图4所示,其孔隙、微裂纹等缺陷增多,因此颗粒破碎强度明显降低。

图  4  干燥与饱和试样的SEM电镜扫描图

Figure  4.  SEM images of dried and saturated samples

下载: 全尺寸图片 幻灯片

2.2   基于FDEM的单颗粒破碎模拟

在单颗粒压缩试验的基础上,开展基于FDEM的单颗粒破碎数值试验,从细观层面进一步研究遇水湿化对颗粒破碎的影响。在连续–离散耦合分析方法中,颗粒被划分为二阶四面体实体单元,实体单元之间插入无厚度的界面单元。实体单元只产生弹性变形,而裂纹萌生、扩展直至失效的全过程,则通过基于内聚力模型的界面单元来模拟^[22-24]。

界面单元在不同应力状态下可能发生拉伸、剪切和拉剪破坏,如图5所示。图中D为损伤因子,当D由0增大至1,界面单元完全失效,实体单元直接接触,裂纹向前扩展。当法向应力$t_{\text{n}}$大于界面抗拉强度$f_{\text{n}}$时,发生拉伸破坏；当切向应力$t_{\text{s}}$大于界面抗剪强度$f_{\text{s}}$时,发生剪切破坏。界面抗剪强度$f_{\text{s}}$由带拉断的Mohr– Coulomb准则确定$f_{\text{s}}=c-t_{\text{n}}\text{tan}\phi$,c为界面凝聚力,φ为材料内摩擦角,法向应力$t_{\text{n}}$以拉为正。考虑拉伸和剪切两种断裂机制以及两者的相互作用,采用二次应力准则作为起裂准则：

图  5  界面单元本构模型示意图

Figure  5.  Constitutive relations for interface elements

下载: 全尺寸图片 幻灯片

式中,$\langle \rangle$为Macaulay符号,表示当$t_{\text{n}}<0$时不会导致损伤的发生。

采用光学扫描技术获取真实颗粒的表面点云数据,然后通过表面点云密度粗化、点云三角网格化、有限元网格划分、插入无厚度界面单元等操作,生成FDEM数值试验所需模型。本文共生成了48个不规则岩石颗粒对应的数值模型,试验组数满足进行单颗粒破碎强度统计的要求^[25]。为避免颗粒破碎的尺寸效应,通过缩放将颗粒的等效粒径（定义为与颗粒等体积的圆球的直径）均调整为60 mm。采用5个形状参数来定量表征不规则岩石颗粒的形状特征：伸长率EI、扁平率FI、Domokos因子S_f、凸度C_X和球度Ψ_3D,每个形状参数对应的平均值和标准差如表1所示。颗粒的形状分布箱形图如图6所示。

表  1  颗粒形状参数表

Table  1.  Parameters of particle shape

形状参数	公式	平均值	标准差
伸长率	EI=I/L	0.80	0.115
扁平率	FI=S/I	0.79	0.099
Domokos因子	$\begin{array}{l}{S}_{\text{f}}=(1/S+1/I+1/L)\times \\ \sqrt{S^2+I^2+L^2}/\sqrt{3}\end{array}$	3.21	0.088
凸度	$C_{\text{X}}=\frac{V}{V_{\text{ch}}}$	0.88	0.031
球度	${\psi }_{\text{3D}}=\frac{\sqrt[3]{36\text{π}{V}^2}}{A}$	0.81	0.031
注：S,I,L分别为颗粒的短轴、中轴、长轴,V和A为颗粒的实际体积和表面积,$V_{\text{ch}}$为颗粒外围凸包体的体积。

下载: 导出CSV 

| 显示表格

图  6  颗粒形状分布图

Figure  6.  Distribution of particle shape

下载: 全尺寸图片 幻灯片

在单颗粒压缩数值试验中,上下加载板均为刚性,先让颗粒在重力作用下自由下落到下加载板上,使其达到稳定状态,尽量避免颗粒的尖角破碎对颗粒破碎强度和破碎模式的影响。加载时,将下加载板固定,上加载板以5×10^-5 mm/步的速率进行位移控制式加载,模拟准静态加载过程。部分单颗粒试样的数值模型及加载示意图如图7所示。

图  7  部分颗粒和单颗粒压缩试验加载示意图

Figure  7.  Partial particle samples and diagram of loading for single particle crushing tests

下载: 全尺寸图片 幻灯片

岩石单轴压缩FDEM数值试验所需的细观参数可分为实体单元参数和界面单元参数。其中实体单元参数以室内单轴压缩试验的成果作为参照对象。而数值试验中界面单元等细观参数的率定,则需要通过不断地试错选优,直至数值试验与室内试验的应力–应变曲线的结果基本一致,如图8所示。其中,为反映颗粒材料内部微裂纹和微孔洞所导致的材料强度的非均匀性,给定抗压强度均值$f_{\text{c}}{}^{\text{e}}$,假定界面单元抗压强度$f_{\text{c}}$服从Weibull分布,其抗拉强度$f_{\text{t}}$则通过给定的$f_{\text{c}}/f_{\text{t}}$比值确定。根据参数率定的结果,干燥与饱和颗粒的FDEM细观参数如表2所示。

图  8  岩石单轴压缩室内及数值试验应力–应变曲线

Figure  8.  Stress-strain curves of rock in uniaxial compression physical and numerical tests

下载: 全尺寸图片 幻灯片

表  2  细观参数取值

Table  2.  Values of mesoscopic parameters

细观参数		干燥试样	饱和试样
实体单元	密度$\rho$/(kg·m-3)	2240	2400
	弹性模量$E$/GPa	13.2	11.88
	泊松比$v$	0.2	0.2
界面单元	抗压强度$f_{\text{c}}{}^{\text{e}}$ /MPa	110	96
	压拉强度比$f_{\text{c}}/f_{\text{t}}$	18.96	21.30
	Weibull模数$m$	10	10
	内摩擦角$\phi$/(°)	40	37
	Ⅰ型断裂能$G_{{\rm I}\text{C}}$/(N·m-1)	85.5	72.0
	Ⅱ型断裂能$G_{\mathrm{Ⅱ}\text{C}}$/(N·m-1)	427.5	360.0
接触模型	摩擦系数${\mu }_{\text{f}}$	0.57700	0.50925
	法向刚度$k_{\text{n}}$/(N·m-3)	9.9×1013	8.91×1013
	切向刚度$k_{\text{s}}$/(N·m-3)	4.125×1013	3.7125×1013

下载: 导出CSV 

| 显示表格

干燥与饱和颗粒的典型荷载–位移关系曲线和颗粒渐进破坏过程基本相同。如图9所示,荷载–位移曲线存在多个峰值,其中荷载的最大峰值通常被认为与颗粒的整体破坏有关。荷载–位移曲线可分为3个主要阶段：加载初期,曲线快速震荡上升,该阶段的局部峰值是由局部尖角破碎引起的；达到峰值后,曲线迅速软化,颗粒承载能力降低,局部碎片增多；加载后期,曲线出现不同程度的波动,这与加载板接触部位的局部破碎更加剧烈有关,大块碎片也出现二次破碎,曲线再次出现较大峰值。界面单元破坏数量的变化趋势与荷载–位移曲线演化基本同步。

图  9  典型荷载–位移曲线和颗粒渐进破坏过程

Figure  9.  Typical load-displacement curve and progressive failure process of particles

下载: 全尺寸图片 幻灯片

2.3   颗粒破碎强度统计分析

分析图9可知,颗粒的峰值荷载与整体破碎密切相关,故将峰值荷载F作为颗粒破碎的特征值。将上下加载板之间的距离d作为颗粒的名义粒径,则颗粒的破碎强度可表示为^[26]

试验研究、数值模拟和理论分析表明颗粒的破碎强度服从Weibull分布^[27-29]。基于Weibull分布模型,对于任意名义粒径为d、破碎强度为$\sigma$的颗粒,其累计幸存概率$P_{\text{s}}$可表示为

式中,m为Weibull模数,$d_0$为基准粒径,${\sigma }_0$为37%累计幸存概率的颗粒特征强度,下文中均称为特征破碎强度。

当进行有限次的单颗粒破碎试验时,试验中某个颗粒的累计幸存概率可根据其破碎强度大小的位次进行估算：

式中,N为试验颗粒总数,i为按升序排列的颗粒破碎强度列表中该颗粒的序号。

当名义粒径相同,即$d=d_0$时,则得出

在$\ln [\ln (1/P_{\text{s}})]$-$\ln \sigma$坐标系中,Weibull模数m和特征破碎强度${\sigma }_0$分别为最优拟合线的斜率和横轴上的截距。

对干燥与饱和颗粒的物理试验和FDEM数值试验结果进行破碎强度统计分析后,得到相应的Weibull拟合线,如图10所示,除少数低破碎强度颗粒外,数值试验与物理试验的数据点均能够较好地满足Weibull拟合线。其中物理试验的干燥与饱和颗粒的特征破碎强度分别为5.16,4.04 MPa,Weibull模数分别为3.537,3.123；数值试验的干燥与饱和颗粒的特征破碎强度分别为4.82,3.73 MPa,Weibull模数分别为3.455,3.244。对比可知,数值试验与物理试验结果较为接近,表明连续-离散耦合分析方法和表2所列参数能反映颗粒在干燥和饱和状态下的破碎特性。

图  10  数值与物理试验结果及Weibull拟合曲线

Figure  10.  Numerical and physical test results and Weibull fitting curves

下载: 全尺寸图片 幻灯片

相较于干燥颗粒,物理试验中饱和颗粒的特征破碎强度降低了21.7%,数值试验中饱和颗粒的特征破碎强度降低了22.6%。该结果表明,在物理试验和数值试验中,遇水湿化对颗粒破碎强度的降低程度相对一致。Weibull模数m反映颗粒破碎强度的离散程度,m越小,离散性越高。两组试样的m值较为接近,因此,遇水湿化对颗粒破碎强度离散性的影响较小。已有研究表明^[30-31],颗粒破碎强度的Weibull模数与颗粒的形状特征有关,本文的湿化数值试验采用形状相同的颗粒,得到近似相等的Weibull模数与之相符。

3.   单轴压缩湿化数值试验

3.1   数值试样制备与加载

进一步采用FDEM数值试验研究湿化对颗粒集合体压缩特性的影响。湿化试验的研究思路主要有单线法和双线法两种。采用单线法进行湿化变形数值试验时,需要把干燥试样加载至某一应力状态,然后模拟通水饱和过程,并通过改变细观参数反映遇水饱和对颗粒的影响。基于单线法的湿化数值试验,模拟难度较大且计算成本较高,不便进行干燥和饱和状态的对比研究。故采用双线法思路开展湿化数值试验,分别对干燥和饱和试样进行加载试验,湿化变形即为两者在同等应力状态下的变形差值。

在FDEM单轴压缩数值试验中,数值试样包括1150个不规则形状颗粒,其初始孔隙比$e_0$为0.614,试样尺寸为300 mm×310 mm（直径×高度）。为分析试验中颗粒级配的演化过程,试样初始级配为近似均一级配,颗粒的等效粒径均为27～29 mm。将数值试样置于圆桶形刚性容器中,对底部及侧面环向刚性加载板施加全约束,上加载板以1.33×10^-4 mm/步的速率进行位移控制式加载,模拟颗粒集合体的准静态单轴压缩。颗粒集合体数值试样及加载示意图如图11所示。

图  11  颗粒集合体数值试样和单轴压缩试验加载示意图

Figure  11.  Schematic diagram of numerical samples and uniaxial compression test loading of particle assembly

下载: 全尺寸图片 幻灯片

3.2   压缩过程

由图12（a）的$e$–$\lg p$压缩曲线可得,饱和试样的屈服应力略小于干燥试样,这是由于遇水饱和导致颗粒强度降低,宏观上表现为更小的屈服应力。轴向应力相同时,饱和试样的孔隙比更小,而且两组压缩曲线在屈服后的斜率基本一致,即压缩指数近似,表明湿化会增大堆石体屈服前的压缩性,但影响较小。这与McDowell等^[32]、Ma等^[33]进行的具有不同颗粒破碎强度的集合体单轴压缩特性类似。

图  12  湿化变形特征曲线

Figure  12.  Characteristic curves of wetting deformation

下载: 全尺寸图片 幻灯片

干燥与饱和试样对应的界面单元破坏率曲线（界面单元破坏数/总数）如图12（b）所示。由图可知,发生屈服前,界面单元破坏率很小,表明加载初期主要是由颗粒的弹性变形和位置重排引起的变形。饱和试样产生附加的湿化变形说明水起润滑作用,使得颗粒间的摩擦系数降低,体现在通过参数率定得到的干燥和饱和状态下的FDEM计算参数（表2）,但其作用相对有限。发生屈服后,破坏率迅速增加,开始发生大范围的颗粒破碎。等轴向应力条件下,饱和试样具有更高的破坏率,表明饱和试样更易发生颗粒破碎。

如图12（c）所示,湿化变形的进程还可以通过干燥与饱和试样的孔隙比差值和颗粒破碎率$p_{\text{c}}$（母颗粒破碎数量与初始数量的比值）来体现。曲线均可大致分为3个阶段：低应力少量破碎阶段、剧烈破碎阶段及高应力缓慢破碎阶段。低应力阶段大多数颗粒的应力水平不满足破碎准则；剧烈破碎阶段,颗粒破碎率迅速上升,湿化变形明显加剧；高应力阶段,破碎率减缓,湿化变形逐渐放缓。

为进一步分析颗粒破碎情况,图13统计了两组试样在压缩过程中的新增颗粒数量。发生屈服前,基本没有新增颗粒,再次证明该阶段变形主要为颗粒弹性变形与位置重排,基本不发生破碎；发生屈服后,颗粒开始发生破碎,新增颗粒数量与轴向应力成正相关。在加载结束时,干燥与饱和试样的最终颗粒数量分别为67537,64489。干燥试样在等应变条件下的轴向应力更大,故最终颗粒数量更多。在轴向应力相等时,饱和试样的颗粒数量明显更多,遇水湿化会显著增加颗粒破碎量。

图  13  单轴压缩过程中的新增颗粒数量

Figure  13.  Number of new particles during uniaxial compression

下载: 全尺寸图片 幻灯片

为比较轴向应力相等时两组试样的颗粒尺寸分布情况,现提取轴向应力${\sigma }_{\text{v}}$=11.25 MPa时的试样级配曲线,如图14所示,插图为颗粒尺寸频数分布。相同应力水平下,饱和试样的级配曲线上移,该试样中小尺寸颗粒多,大尺寸颗粒少,表明湿化会降低颗粒的破碎强度,破碎更易发生,故级配曲线上移。

图  14  试样级配曲线(${\sigma }_{\text{v}}$=11.25 MPa)

Figure  14.  Gradation curves of samples (${\sigma }_{\text{v}}$=11.25 MPa)

下载: 全尺寸图片 幻灯片

4.   湿化变形机理分析

水对堆石颗粒具有润滑作用,使得颗粒间摩擦系数降低,促进颗粒滑移与位置重排,是低应力水平下产生湿化变形的主要原因,但因摩擦系数降低而产生的湿化变形量较小,FDEM模拟表明湿化使颗粒间的摩擦系数仅降低了11.7%。水对堆石颗粒具有软化作用,颗粒的矿物成分和微观结构发生明显的变化,使其变形模量和破碎强度降低,FDEM模拟表明湿化使颗粒破碎强度降低了22.6%,湿化使堆石体的压缩模量降低了13.1%。遇水湿化之后的颗粒集合体在颗粒接触处产生更大的接触变形和更多颗粒破碎,宏观上表现出湿化变形,并且颗粒破碎导致的变形占堆石体湿化变形的主要部分。

5.   结论

本文采用室内试验和连续–离散耦合分析方法进行了湿化对堆石颗粒破碎特性的影响研究,通过SEM电镜扫描和数值模拟,揭示了湿化导致颗粒破碎强度降低的机理,并进行了集合体在干燥和饱和状态下的压缩特性研究,主要结论如下：

（1）干燥与饱和颗粒的SEM电镜扫描结果表明,遇水饱和导致颗粒的岩石矿物成分和微观结构发生变化,其孔隙、微裂纹等缺陷增多,颗粒破碎强度和表面摩擦系数降低。单颗粒破碎物理试验和FDEM数值试验均表明,干燥与饱和试样的破碎强度服从Weibull分布,且Weibull模数相近,两组饱和颗粒的特征破碎强度的降低程度基本一致,其中FDEM模拟的饱和颗粒的特征破碎强度降低了22.6%,表明湿化对颗粒强度的劣化作用明显。

（2）在颗粒集合体的单轴压缩数值试验中,湿化降低了颗粒材料的屈服应力和压缩模量,增大了其压缩性。颗粒破碎主要发生在屈服后,与干燥试样相比,饱和试样在等应力水平下的颗粒破碎程度更大。

（3）堆石体湿化变形的细观机理可总结为,遇水湿化对堆石体有润滑、软化和劣化作用,使堆石颗粒的摩擦系数、变形模量和破碎强度均出现了不同程度的降低。其中,破碎强度降低导致的颗粒破碎是堆石体产生湿化变形的最主要原因。