强震作用下面板堆石坝跨尺度面板开裂演化分析

孔宪京; 屈永倩; 邹德高; 陈楷; 刘京茂

doi:10.11779/CJGE202006001

强震作用下面板堆石坝跨尺度面板开裂演化分析 English Version

孔宪京^{1, 2,},
屈永倩^{1, 2, ,},
邹德高^{1, 2},
陈楷^{1, 2},
刘京茂^{1, 2}

1.
大连理工大学海岸与近海工程国家重点实验室,辽宁　大连　116024
2.
大连理工大学水利工程学院,辽宁　大连　116024

基金项目:

国家自然科学基金项目 U1965206

国家自然科学基金项目 51779034

中央高校基本科研业务费资助项目 DUT19ZD216

华能集团科技项目 HNKJ18-H25

详细信息

作者简介:
孔宪京(1952—),男,博士,教授,主要从事高土石坝抗震和岩土地震工程方面的研究。E-mail: kongxj@dlut.edu.cn

通讯作者:
屈永倩, E-mail: quyongqian@mail.dlut.edu.cn

中图分类号: TV64
计量
- 文章访问数: 437
- HTML全文浏览量: 23
- PDF下载量: 292
出版历程
- 收稿日期: 2019-09-18
- 网络出版日期: 2022-12-07
- 刊出日期: 2020-05-31

Cross-scale crack evolution analysis for face slab in concrete faced rockfill dams under strong earthquake

KONG Xian-jing^{1, 2,},
QU Yong-qian^{1, 2, ,},
ZOU De-gao^{1, 2},
CHEN Kai^{1, 2},
LIU Jing-mao^{1, 2}

1.
State Key Laboratory of Coastal and Offshore Engineering, Dalian University of Technology, Dalian 116024, China
2.
School of Hydraulic Engineering, Dalian University of Technology, Dalian 116024, China

摘要

摘要: 准确定位面板薄弱区域、定量评价面板破坏程度对面板坝抗震安全评价至关重要。联合非完全点对点界面和四分树网格生成技术,实现了面板坝面板与垫层接触作用的两级跨尺度精细化建模；引入混凝土黏聚力模型并联合筑坝材料广义塑性模型、状态相关的弹塑性接触面模型考虑材料强非线性和破坏过程,建立了强震作用下面板动力弹塑性跨尺度开裂演化分析方法,并研发了显式地震波动分析框架下的SBFEM-FEM耦合计算软件。以200 m级面板坝为例进行了面板动力破坏数值分析,并分别考虑了面板配筋率、竖向地震和坝前水深的影响。结果表明：十分直观地再现面板坝防渗面板的动力开裂演化过程,有助于定位面板局部薄弱区域、定量评价面板破坏程度以及评估抗震措施加固效果,为混凝土面板抗震优化设计和极限抗震能力评估提供技术支持。可拓展用于其它混凝土防渗结构破坏计算,且可容易地扩展至三维分析应用。
- 强震作用 /
- 面板堆石坝 /
- 开裂演化 /
- 黏聚力模型 /
- 弹塑性分析 /
- 跨尺度
Abstract: It is crucially important for seismic safety evaluation of high concrete faced rockfill dams (CFRDs) to accurately locate the weak area of the panel and to quantitatively assess the damage of the face slab. In this study, the cross-scale model for CFRD is established using the interface element with asymmetric nodes and Quadtree for refined simulation of slab and cushion interaction. The cohesive zone model for concrete, the generalized plastic model for rockfill and the state-dependent elasto-plastic interface model are combined and used to describe the strong nonlinearity and failure process. On the above basis, the cross-scale crack evolution analysis method under strong earthquake is established and the coupled SBFEM-FEM analysis software is developed in the explicit earthquake wave motion input method frame. The dynamic failure analyses of slabs are performed for a 200-m-high CFRD considering the reinforcement ratio, vertical earthquake and water level of reservoir. The results indicate that the developed method can visually represent the seismic cracking evolution, conveniently locate weak areas of face slab, quantitatively determine the damage severity, and evaluate the aseismic measures. The research results may provide an effective method for the aseismic design and assessment of ultimate aseismic capacity of concrete slab. The proposed method can be extended to the failure analyses of other concrete structures and three-dimensional investigation and application easily.
- strong earthquake /
- CFRD /
- cracking evolution /
- cohesive zone model /
- elasto-plastic analysis /
- cross-scale

HTML全文

0. 引言

天然土体颗粒间通常具有胶结物质,比如碳酸钙、微生物、黏粒、甲烷水合物等,称为结构性土。结构性土可看成是一种特殊的胶结颗粒材料,颗粒间胶结物质的存在使结构性土的宏观力学性质不同于重塑土^[1-2]。为了分析结构性土工程问题,需要建立合理的结构性土本构模型。在临界状态土力学基础上,相关学者引入状态变量,分别建立了边界面模型^[3]、下负荷面模型^[4]、统一硬化模型^[5]等重塑土本构模型。在加荷过程中,结构性土颗粒间的胶结逐渐破坏并退出承担应力,土体向重塑土转化。因此,在重塑土本构模型基础上,通过引入表征结构性的参量及其演化规律建立结构性土本构模型是一条可行的研究途径。不少学者提出了结构性土本构模型,在模型中假定了不同的结构性表征参量。比如：①通过比较结构性土屈服面大小 ${\bar{p}}_{c}$ 和重塑土屈服面大小 $p_{c}$ 的不同,定义 $Δ p_{c} = {\bar{p}}_{c} - p_{c}$ ^[6]或 $R^{*} = p_{c} /$ ${\bar{p}}_{c}$ ^[7]为结构性表征参量；②通过引入胶结组分（相对完整状态）和摩擦组分（完全调整状态）应力分担参量（扰动参量、结构破损参量）D^[8]或B₀^[9-10]表达结构性衰退程度。由于胶结破损演化难以通过试验测得,结构性表征参量的演化通常需要假设,缺乏微观试验结果验证。

离散单元法^[11]模拟结构性土弥补了室内微观测试在获取胶结破损演化方面的不足,已广泛应用于结构性土体宏微观力学性质研究。Jiang等^[10]通过开展结构性土的二维离散元模拟,对3种结构性土本构模型结构性参量演化假设的合理性进行了验证,基于离散元结果提出了新的胶结破损演化式,并应用到本构模型研究^[12],不足之处在于胶结破损演化式中的参数受加载路径影响。通过PFC软件植入完整胶结接触模型（考虑了颗粒间抗转动和抗扭转能力及胶结尺寸对刚度和强度的影响）建立了结构性土三维离散元试样,开展了侧限压缩、等向压缩、等应力比压缩以及常规三轴和真三轴试验的离散元数值模拟,分析了结构破损参量B₀随等效塑性应变的演化规律,并提出了应力路径相关性低的结构破损参量B_σ。

1. 胶结接触模型

离散元模拟采用软硬复合胶结接触模型框架^[13],考虑了颗粒间抗转动和抗扭转能力^[14]及胶结尺寸对刚度和强度的影响^[15]。对于通过胶结间接传力/力矩的接触（“厚胶结”）直接采用胶结接触公式计算接触力/力矩；对于颗粒直接接触并有胶结传力/力矩的接触（“薄胶结”）,按并联模式通过颗粒接触力/力矩和胶结接触力/力矩之和计算接触力/力矩。

1.1 颗粒接触

假设颗粒间通过一圆形截面接触,半径为

$R_{c} = β R$ ,

(1)

式中, $β$ 为抗转动系数, $R = 2 R_{1} R_{2} / (R_{1} + R_{2})$ , $R_{1}$ , $R_{2}$ 为两颗粒半径。

接触间力/力矩的传递通过下列公式计算：

$F_{n}^{l} = {\begin{cases} k_{n} δ_{n} - F_{a} (δ_{n} \geq 0) \\ 0 (δ_{n} < 0) \end{cases}$ ,

(2)

$F_{s}^{l} = {(F_{s}^{l})}_{o} - k_{s} Δ δ_{s}$ ,

(3)

$M_{r}^{l} = {(M_{r}^{l})}_{o} - k_{r} Δ θ_{r}$ ,

(4)

$M_{t}^{l} = {(M_{t}^{l})}_{o} - k_{t} Δ θ_{t}$ 。

(5)

式中 $F_{a}$ 为粒间引力,比如毛细力、范德华力、电子引力等； $F_{n}^{l}$ , $F_{s}^{l}$ , $M_{r}^{l}$ , $M_{t}^{l}$ 分别为颗粒法向力、切向力、弯矩和扭矩； $k_{n}$ , $k_{s}$ ,k_r,k_t分别为颗粒法向、切向、抗弯、抗扭刚度, $k_{r} = 0.25 k_{n} R_{c}^{2}$ , $k_{t} = 0.5 k_{s} R_{c}^{2}$ ； $δ_{n}$ 为颗粒重叠量, $Δ δ_{s}$ , $Δ θ_{r}$ 和 $Δ θ_{t}$ 分别为剪切位移、转动角和扭转角增量。

颗粒接触的法切向刚度可按PFC软件手册^[16]按下式计算：

$k_{n} = A E^{*} / L$ ,

(6)

$k_{s} = k_{n} / κ^{*}$ ,

(7)

式中, $A = π R^{2}$ , $E^{*}$ 是颗粒接触的等效模量, $κ^{*}$ 是颗粒接触的法切向刚度比, $L$ 为两颗粒半径之和。

接触抗剪强度、抗弯强度和抗扭强度分别为^[14-15]

$‖ F_{s}^{l} ‖ \leq μ (F_{n} + F_{a})$ ,

(8)

$‖ M_{r}^{l} ‖ \leq 0.25 ζ_{c} R_{c} (F_{n} + F_{a})$ ,

(9)

$M_{t}^{l} \leq 0.65 μ R_{c} (F_{n} + F_{a})$ ,

(10)

式中, $μ$ 是摩擦系数, $ζ_{c}$ 取2.1。颗粒接触力和力矩达到以上强度后颗粒将进入滑移、转动和扭转状态。

1.2 胶结计算

假定胶结物质呈圆形截面连接两颗粒,则胶结法向力 $F_{n}^{b}$ ,切向力 $F_{s}^{b}$ ,弯矩 $M_{r}^{b}$ 和扭矩 $M_{t}^{b}$ 分别为^[16]

$F_{n}^{b} = {(F_{n}^{b})}_{o} - k_{n}^{b} Δ δ_{n}^{}$ ,

(11)

$F_{s}^{b} = {(F_{s}^{b})}_{o} - k_{s}^{b} Δ δ_{s}^{}$ ,

(12)

$M_{r}^{b} = {(M_{r}^{b})}_{o} - k_{r}^{b} Δ θ_{r}^{}$ ,

(13)

$M_{t}^{b} = {(M_{t}^{b})}_{o} - k_{t}^{b} Δ θ_{t}$ 。

(14)

式中 $k_{n}^{b}$ , $k_{s}^{b}$ , $k_{r}^{b}$ , $k_{t}^{b}$ 分别为胶结法向、切向、抗弯和抗扭刚度, $k_{t}^{b} = 0.5 k_{s}^{b} R_{b}^{2}$ , $k_{r}^{b} = 0.25 k_{n}^{b} R_{b}^{2}$ ； $Δ δ_{n}^{}$ 是法向位移增量； $R_{b} = λ_{b} R$ 为胶结半径, $λ_{b}$ 为胶结半径系数。

为了考虑胶结厚度和宽度对胶结刚度的影响,用下式计算胶结法向刚度（图1）^[15]：

图 1 胶结刚度计算示意图

Figure 1. Schematic diagram for calculation of bond stiffness

下载: 全尺寸图片幻灯片

$k_{n}^{b} = 1 / (δ_{A 1} + δ_{A 2} + δ_{B 1} + δ_{B 2} + δ_{C 1} + δ_{C 2})$ ,

(15)

式中, $δ_{A 1}$ , $δ_{B 1}$ , $δ_{C 1}$ 分别为颗粒1的A1_,B1和C1三段的线柔度,

$δ_{A 1} = \frac{R_{1} ξ_{E}}{π E_{1} R_{b}^{2}}$ ,

(16)

$δ_{B 1} = {\begin{matrix} \frac{1}{2 π (E_{1} - E_{2}) R_{1} Λ_{E}} (\ln (\frac{Λ_{E} + 1}{Λ_{E} - 1}) - \ln (\frac{Λ_{E} + ξ_{E}}{Λ_{E} - ξ_{E}})) (E_{1} > E_{2}) \\ \frac{R_{1} - R_{1} ξ_{E}}{π E_{1} R_{b}^{2}} (E_{1} = E_{2}) \end{matrix}$ ,

(17)

$δ_{C 1} = δ_{C 2} = \frac{h_{0}^{b} / 2}{π E_{2} R_{b}^{2}}$ ,

(18)

式中, $E_{1}$ , $E_{2}$ 分别为骨架材料和胶结材料的等效变形模量, $R_{1}$ 为颗粒1的半径, $h_{0}^{b}$ 为胶结厚度,同理可计算 $δ_{A 2}$ , $δ_{B 2}$ 。在模型使用中, $E_{1} = {\bar{E}}^{*}$ 为胶结等效模量, $E_{2} = η_{E} {\bar{E}}^{*}_{}$ , $η_{E}$ 为模量折减系数。 $ξ_{E}$ , $Λ_{E}$ 计算如下：

$ξ_{E} = \sqrt{1 - {(R_{b} / R_{1})}^{2}}$ ,

(19)

$Λ_{E} = \sqrt{1 + \frac{E_{2}}{E_{1} - E_{2}} {(R_{b} / R_{1})}^{2}}$ 。

(20)

胶结切向刚度根据法切向刚度比 $κ_{b}$ 计算：

$k_{s}^{b} = k_{n}^{b} / κ_{b}$ 。

(21)

采用考虑胶结厚度和宽度的胶结强度计算方法^[15,17]。胶结抗压和抗拉强度计算如下：

$R_{nc}^{b} = π R_{b}^{2} σ_{c}^{b} \frac{c_{1}}{λ_{b}^{c_{3}}} \exp (\frac{c_{2} h_{0}^{b}}{2 R_{b}})$ ,

(22)

$R_{nt}^{b} = π R_{b}^{2} σ_{t}^{b}$ ,

(23)

式中, $σ_{c}^{b}$ , $σ_{t}^{b}$ 分别为胶结材料抗压和抗拉强度； $c_{1} = 1.45$ , $c_{2} = - 0.7$ 和 $c_{3} = 1 / 6$ 是反映胶结半径和厚度对抗压强度影响的参数（定值）。

在一定的法向力下（未达到胶结抗压和抗拉强度）,胶结可能承受剪切力、弯矩和扭矩,复合荷载作用下胶结强度符合椭球包面：

${(‖ F_{s}^{b} ‖ / R_{s}^{b})}^{2} + {(‖ M_{r}^{b} ‖ / R_{r}^{b})}^{2} + {(M_{t}^{b} / R_{t}^{b})}^{2} = 1$ ,

(24)

式中, $R_{s}^{b}$ , $R_{r}^{b}$ , $R_{t}^{b}$ 分别为胶结抗剪、抗弯和抗扭强度,

$R_{s}^{b} = S_{s}^{b} (R_{nc}^{b} + R_{nt}^{b})$ ,

(25)

$R_{r}^{b} = S_{r}^{b} (R_{nc}^{b} + R_{nt}^{b}) R_{b}$ ,

(26)

$R_{t}^{b} = S_{t}^{b} (R_{nc}^{b} + R_{nt}^{b}) R_{b}$ ,

(27)

式中, $S_{s}^{b}$ , $S_{r}^{b}$ , $S_{t}^{b}$ 为胶结抗剪、抗弯和抗扭包面的形状,

$S_{i}^{b} = m_{i} f_{n}^{b} {[\ln (1 / f_{n}^{b})]}^{3 / 5}$ ,

(28)

式中,下标i是s,r和t, $m_{s} = 0.5$ , $m_{r} = 0.3$ 和 $m_{t} = 0.36$ 为参数（定值）, $f_{n}^{b} = (F_{n}^{b} + R_{nt}^{b}) / (R_{nc}^{b} + R_{nt}^{b})$ 为正则化胶结法向力。

胶结的存在对结构性土体的宏微观力学性质有很大影响。胶结达到强度包面后会发生胶结破坏,胶结破坏后胶结部分不再承担力和力矩,接触力和力矩退化到残余值（颗粒接触部分）。

2. 离散元模拟步骤

2.1 离散元试样制备

以饱和结构性黄土^[15]为基准制备离散元试样和选取接触参数,为提高研究的适用性开展了不同胶结材料强度离散元试样的数值分析,下文仅以一种胶结材料强度为例进行分析。采用分层欠压法^[18]分5层制备离散元试样,颗粒数目42180个,如图2（a）所示。结构性土离散元试样制备采用的颗粒级配曲线如图2（b）所示。试样生成后在12.5 kPa竖向压力下压缩稳定以反映原位K₀应力状态。试样稳定后若两个颗粒之间的间距小于 $g_{c} R_{b}$ 则生成胶结, $g_{c}$ 为胶结临界厚度系数。

图 2 离散元试样及其粒径级配

Figure 2. DEM sample and grain-size distribution curve

下载: 全尺寸图片幻灯片

颗粒间引力作用能使试样以更松散状态稳定赋存^[19],试样制备过程施加弱吸引力以制备松散试样^[15]：

$F_{a} = σ_{a} d_{50}^{2}$ ,

(29)

式中,d₅₀为中值粒径, $σ_{a}$ 为弱吸引力系数,根据典型黄土的最大稳定孔隙比和饱和重塑土黏聚力,通过试算取4 kPa。颗粒与墙体间仅承担法向力作用,忽略摩擦。颗粒密度=2710 kg/m³,颗粒局部阻尼系数=0.7。颗粒之间的接触参数：颗粒等效模量=800 MPa,法切向刚度比=1.5,摩擦系数=0.5,抗转动系数=0.21,胶结等效模量=200 MPa,胶结模量折减系数=0.2,胶结法切向刚度比=2,胶结半径系数=0.35,胶结材料抗压强度=25MPa,胶结材料抗拉强度=2.5MPa,胶结临界厚度系数=0.1。

2.2 不同应力路径离散元模拟步骤

为了研究不同应力路径结构性土的胶结破损演化,对上述结构性土离散元试样进行了分级加载的侧限压缩试验、等向压缩试验和等应力比压缩试验,以及应变加载的常规三轴试验和等p真三轴试验。

对分级加载压缩试验,上下墙体从12.5 kPa开始分级施加荷载,每次施加的荷载约为上级荷载的 $\sqrt{2}$ 倍；侧墙伺服应力目标通过应力比要求计算（侧限压缩试验使4个侧墙固定）。

对三轴试验,首先将离散元试样等向压缩到不同的围压,然后对侧墙进行伺服控制,使上下两墙相向运动进行剪切试验。在常规三轴试验中,侧墙的伺服目标为围压。在真三轴试验中通过平均应力和中主应力系数b不变的条件计算侧墙伺服目标应力。

3. 离散元模拟结果

为验证离散元数值分析结果的有效性,在研究结构性土的胶结破损演化规律之前,需确保离散元试样能够反映结构性土的主要压缩和剪切特性。考虑到本文针对一般结构性土体（不同胶结强度的胶结颗粒材料）进行破损规律数学建模,因此,下面主要进行定性对比,未对某一特定的结构性土进行标定。

图3给出了结构性土离散元试样在不同压缩路径下的压缩曲线。图3中 $η = q / p$ 为应力比, $CSR = σ_{3} / σ_{1}$ 为主应力比。随着平均应力的增加,离散元试样压缩线先缓慢下降而后逐渐转变为快速下降,存在结构屈服应力,离散元试样能够定性再现结构性土体的主要压缩特性^[20]。离散元试样体积变形特征与真实结构性土定量上差别较大,主要是因为颗粒形状、颗粒破碎机制等微观结构的差别。

图 3 结构性土压缩性质

Figure 3. Compression properties of structured soils

下载: 全尺寸图片幻灯片

图4给出了结构性土离散元试样常规三轴试验应力–应变关系和孔隙比–应变关系。离散元分析能够定性再现结构性土体三轴剪切性质^[2,21]。试样在剪切初始段模量很大,近似弹性段。随着围压的增加,试样抗剪强度增加,孔隙比减小,试样从倾向于软化和剪胀向倾向于硬化和剪缩发展。离散元反映结构性土剪切破坏特性的效果优于反映压缩变形特性。

图 4 结构性土常规三轴剪切性质

Figure 4. Conventional triaxial test results of structured soils

下载: 全尺寸图片幻灯片

图5给出了结构性土离散元试样真三轴试验应力-应变关系和孔隙比–应变关系（以p=300 kPa为例）。可见离散元模拟结果能够反映室内试验中中主应力系数b（ $b = (σ_{2} - σ_{3}) / (σ_{1} - σ_{3})$ ）的影响^[22-23]。同样平均应力下随着b的增加,破坏偏应力下降。b=0时试样倾向于剪缩,随着b的增加试样剪缩程度减弱。天然土试样的变形模量相对离散元试样和人工制备结构性土试样较低,这可能与试样扰动有关。

图 5 结构性土真三轴剪切性质

Figure 5. True triaxial test results of structured soils

下载: 全尺寸图片幻灯片

4. 胶结破损演化规律

4.1 胶结颗粒材料微观力学

土体可以看成一种颗粒集合体,重塑土（无黏结）可以看成是一种散粒体材料而结构性土可以看成是一种胶结颗粒材料。对于理想的散粒体材料,土中的应力通过接触力进行传递,而土中的应变由颗粒间的相对变形产生。基于此,学者通过研究土颗粒之间的接触关系和微观变形机理,建立了理想散粒体材料宏观应力和应变与微观接触力和变形之间的基本方程^[24-25]。Jiang等^[10]在此基础上建立了胶结颗粒材料的宏微观理论框架。

对于颗粒材料,荷载通过颗粒之间的接触力进行传递,因此颗粒材料的代表性单元的平均应力张量可根据下式求得

${\bar{σ}}_{i j}^{} = \frac{1}{V^{}} \sum_{k = 1}^{N^{}} F_{i}^{k} l_{j}^{k}$ ,

(30)

式中,V为代表性单元的体积,N为代表性单元内颗粒间接触点的总数, $F_{i}^{k}$ 为第k个接触点处接触力的分量, $l_{j}^{k}$ 为第k个接触点处方向向量的分量。

对于结构性土体,根据颗粒间接触点是否存在胶结,土体代表性单元平均应力张量 ${\bar{σ}}_{i j}^{}$ 可分解为胶结接触分担应力 ${\bar{σ}}_{i j}^{in}$ 与无胶结接触分担应力 ${\bar{σ}}_{i j}^{d}$ 之和,

$\begin{array}{l} {\bar{σ}}_{i j} = \frac{1}{V^{}} (\sum_{k = 1}^{N^{in}} F_{i}^{k} l_{j}^{k} + \sum_{k = 1}^{N^{d}} F_{i}^{k} l_{j}^{k}) \\ = \frac{1}{V^{}} （ V^{in} {\bar{σ}}_{i j}^{in} + V^{d} {\bar{σ}}_{i j}^{d}) = (1 - λ) {\bar{σ}}_{i j}^{in} + λ {\bar{σ}}_{i j}^{d} \end{array}$ ,

(31)

式中, $N^{in}$ , $N^{d}$ 分别为代表性单元内胶结接触点和无胶结接触点的数目,且 $N = N^{in} + N^{d}$ 。 $V^{in}$ , $V^{d}$ 分别为代表性单元内胶结颗粒所占体积与无胶结颗粒所占体积（含相应的孔隙体积）,且 $V = V^{in} + V^{d}$ ； $λ = V^{d} / V$ 为体积破损率,即代表性单元内无胶结组分所占体积与代表性单元总体积的比值； ${\bar{σ}}_{i j}^{in}$ 和 ${\bar{σ}}_{i j}^{d}$ 可由下式计算：

$\begin{array}{l} {\bar{σ}}_{i j}^{in} = \frac{1}{V^{in}} \sum_{k = 1}^{N^{in}} F_{i}^{k} l_{j}^{k}, \\ {\bar{σ}}_{i j}^{d} = \frac{1}{V^{d}} \sum_{k = 1}^{N^{d}} F_{i}^{k} l_{j}^{k} 。 \end{array}}$

(32)

为了研究结构性土体的结构性演化,需要在本构模型中引入（或蕴含）胶结破损参量。体积破损率就是一个简单直观的胶结破损参量,考虑到胶结组分和无胶结组分的体积难以划分,可用无胶结接触数量与总接触数量比值近似表示,

$λ = \frac{V^{d}}{V} \approx \frac{N^{d}}{N}$ 。

(33)

无胶结组分分担的应力与总应力的比值也可以作为胶结破损参量,

$B_{0} = \frac{λ {\bar{σ}}_{i j}^{d}}{{\bar{σ}}_{i j}} = \frac{\sum_{k = 1}^{N^{d}} F_{i}^{k} l_{j}^{k}}{\sum_{k = 1}^{N} F_{i}^{k} l_{j}^{k}}$ 。

(34)

胶结破损参量表征胶结破损程度或结构性衰退程度。当胶结破损参量为0时,土体结构性完好,土体应力全部由胶结组分承担；随着胶结破损参量的增加,土体应力由胶结组分和无胶结组分共同分担（胶结组分承担应力逐渐变小）；当胶结破损参量为1时,土体结构性完全破坏,土体应力全部由无胶结组分承担。

胶结破损参量 $λ$ , $B_{0}$ 难以通过室内试验获得,离散单元法可以用来研究 $λ$ , $B_{0}$ 的演化。尽管 $λ$ , $B_{0}$ 参量的大小有一定不同,但二者演化规律是相似的,因此本文以破损参量 $B_{0}$ 为例进行研究。

4.2 胶结破损演化规律

随着塑性应变的发展,胶结逐渐破损。由于弹性应变较小,塑性应变采用如下简单方法计算。离散元试样的弹性能 $E_{k}$ 为颗粒弹性能 $E_{k}^{p}$ 和胶结弹性能 $E_{k}^{b}$ 之和：

$E_{k}^{p} = \frac{1}{2} (\frac{{(F_{n}^{l})}^{2}}{k_{n}} + \frac{{‖ F_{s}^{l} ‖}^{2}}{k_{s}} + \frac{{‖ M_{r}^{l} ‖}^{2}}{k_{r}} + \frac{{(M_{t}^{l})}^{2}}{k_{t}})$ ,

(35)

$E_{k}^{b} = \frac{1}{2} (\frac{{(F_{n}^{b})}^{2}}{k_{n}^{b}} + \frac{{‖ F_{s}^{b} ‖}^{2}}{k_{s}^{b}} + \frac{{‖ M_{r}^{b} ‖}^{2}}{k_{r}^{b}} + \frac{{(M_{t}^{b})}^{2}}{k_{t}^{b}})$ ,

(36)

$E_{k} = E_{k}^{p} + E_{k}^{b} = \frac{1}{2} (p ε_{v}^{e} + q ε_{s}^{e}) V_{DEM}$ ,

(37)

$K = p / ε_{v}^{e}$ ,

(38)

$G = q / ε_{s}^{e}$ ,

(39)

式中, $V_{DEM}$ 为离散元试样的体积, $ε_{v}^{e}$ 为弹性体应变, $ε_{s}^{e}$ 为弹性偏应变,假设体积模量K和剪切模量G相等,联立以上各式可求得 $ε_{v}^{e}$ , $ε_{s}^{e}$ ,从而求出塑性体应变和塑性偏应变。

图6给出了破损参量随塑性体应变和塑性偏应变的变化规律。对于三轴剪切试验,塑性体应变较小且变化复杂,而对于等应力比压缩试验（应力比较小时）,塑性偏应变较小,二者均不适合单独用于描述试样破损参量演化。因此,胶结的破损既与塑性体应变有关,又与塑性偏应变有关。

图 6 不同应力路径下破损参量B₀随塑性体应变和塑性偏应变的演化规律

Figure 6. Evolution of degradation variable B₀ versus plastic volumetric and deviator strain for DEM sample under different stress paths

下载: 全尺寸图片幻灯片

为了考虑两种塑性应变对胶结破损参量的影响,采用下式简单等效塑性应变研究胶结破损演化。

$ε_{d}^{p} = \sqrt{{(ε_{v}^{p})}^{2} + {(ε_{s}^{p})}^{2}}$ ,

(40)

式中, $ε_{v}^{p}$ 为塑性体应变, $ε_{s}^{p}$ 为塑性偏应变。

图7给出结构性土离散元试样胶结破损参量 $B_{0}$ 随等效塑性应变的变化规律。从图7中可见,对不同应力路径,破损参量 $B_{0}$ 随等效塑性应变的增加先较快增加,而后缓慢增加。但不同应力路径的破损参量 $B_{0}$ 演化在数值上有较大不同。对不同应力比压缩试验,随着应力比的增加, $B_{0}$ 表征的胶结破损速率加快；对不同中主应力系数b的真三轴试验,随着b的增加, $B_{0}$ 表征的胶结破损速率减小；对不同围压常规三轴试验,围压也对 $B_{0}$ 表征的胶结破损速率影响显著。

图 7 不同应力路径下破损参量B₀随等效塑性应变的演化规律

Figure 7. Evolution of degradation variable B₀ versus equivalent plastic strain for DEM sample under different stress paths

下载: 全尺寸图片幻灯片

考虑到胶结破损跟胶结强度相关,在宏观上,胶结破损速率跟目前应力状态和强度包线之间的相对关系相关,结构参数 $B_{0}$ 破损的速率随应力比的增加而增加。因此,将应力比和强度包线引入胶结破损参量 $B_{0}$ ,通过应力比与峰值应力比的比值 $η / M_{f}$ 调整胶结破损速率,构建一个新的胶结破损参量如下式,从而蕴含应力路径对胶结破损的影响,

$B_{σ} = \frac{2 B_{0}}{(1 + η^{2} / M_{f}^{2})}$ ,

(41)

式中, $M_{f} = q_{f} / p$ 为峰值应力比,本文通过常规三轴试验结果采用莫尔–库仑准则得到峰值偏应力 $q_{f} =$ $0.9275 p + 123.8 kPa$ 。

图8给出了结构性土离散元试样胶结破损参量 $B_{σ}$ 随等效塑性应变的破损规律。对不同应力路径,胶结破损参量 $B_{σ}$ 随等效塑性应变逐渐增加,但破损速率逐渐降低。除等向压缩试验（ $η = 0$ ）外,试验应力比、围压和中主应力系数对 $B_{σ}$ 表征的胶结破损演化影响很小。由此可见,胶结破损参量 $B_{σ}$ 应力路径相关性明显降低,可作为状态变量引入结构性土体二元介质类本构模型（比如扰动状态概念模型、损伤模型等）表征胶结破损情况,使模型的结构破损规律更为合理。通过数据分析,建议采用指数函数拟合破损参量 $B_{σ}$ 与等效塑性应变的关系,

图 8 不同应力路径下破损参量B_σ随等效塑性应变的演化规律

Figure 8. Evolution of degradation variable B_σversus equivalent plastic strain for DEM sample under different stress paths

下载: 全尺寸图片幻灯片

$B_{σ} = 1 - \exp [- c_{a} {(ε_{d}^{p})}^{c_{b}}]$ ,

(42)

式中, $c_{a}$ , $c_{b}$ 为胶结破损相关的参数,可通过侧限压缩或常规三轴试验确定。

5. 结论

通过将相对完备的三维胶结接触模型植入粒间接触,反映结构性土粒间胶结的关键特征,建立了结构性土离散元模型；开展了侧限压缩、等向压缩、等应力比压缩以及常规三轴、真三轴等不同加载路径下的离散元数值分析；再现了结构性土体的主要宏观力学性质（尤其是结构性土的破坏特性）；分析了胶结破损参量 $B_{0}$ , $B_{σ}$ 随等效塑性应变的演化规律。

（1）胶结破损参量（ $λ$ , $B_{0}$ , $B_{σ}$ ）表征胶结破损程度或结构性衰退程度。当胶结破损参量为0时,土体结构性完好,土体应力全部由胶结组分承担；当胶结破损参量为1时,土体结构性完全破坏,土体应力全部由无胶结组分承担。

（2）胶结破损参量 $B_{0}$ 的演化受应力比、围压和中主应力系数等的影响,具有应力路径相关性。将应力比和强度包线引入胶结破损参量 $B_{0}$ ,构建了一个新的胶结破损参量 $B_{σ}$ ,应力路径相关性明显降低。在本构模型应用中,胶结破损参量 $B_{σ}$ 随等效塑性应变的演化规律可用指数函数表示,并通过侧限压缩或常规三轴试验确定相关参数。

不同应力路径下 $B_{σ}$ 的演化仍然具有一定的离散性,今后将考虑改进等效塑性应变的定义和强度包面公式。因为胶结破损规律 $B_{σ}$ 基于离散元结果提出,适用于理想胶结颗粒材料（理想结构性土体）,对真实结构性土体的适用性需要进一步验证。

图 1 黏聚力模型和界面单元

Figure 1. Cohesive zone model and interface elements

下载: 全尺寸图片幻灯片

图 2 牵引力和位移关系

Figure 2. Relationship between traction and displacement

下载: 全尺寸图片幻灯片

图 3 混凝土梁几何模型

Figure 3. Geometry of concrete beams

下载: 全尺寸图片幻灯片

图 4 测得的冲击荷载时间曲线

Figure 4. Time histories of measured loads

下载: 全尺寸图片幻灯片

图 5 加载点的荷载位移曲线

Figure 5. Relationship between load and displacement at loading point

下载: 全尺寸图片幻灯片

图 6 裂纹尖端扩展

Figure 6. Time histories of crack extension

下载: 全尺寸图片幻灯片

图 7 测点应变时程

Figure 7. Time histories of strain at observing points

下载: 全尺寸图片幻灯片

图 8 面板堆石坝跨尺度分析模型

Figure 8. Cross-scale model for CFRDs

下载: 全尺寸图片幻灯片

图 9 地震动加速度时程

Figure 9. Acceleration histories of input ground motion

下载: 全尺寸图片幻灯片

图 10 地震过程中裂缝位置和最大宽度

Figure 10. Location and maximum width of crack during earthquake

下载: 全尺寸图片幻灯片

图 11 面板迎水面开裂单元的裂缝宽度、刚度退化因子、混凝土顺坡向应力和钢筋应力时程

Figure 11. Time histories of axial stress of steel, concrete stress along slope, degeneration factor, and crack width of at upstream side during earthquake

下载: 全尺寸图片幻灯片

图 12 面板震后残余裂缝

Figure 12. Residual crack after earthquake

下载: 全尺寸图片幻灯片

图 13 面板顺坡向静动叠加应力

Figure 13. Stresses along slope of face slab at certain time

下载: 全尺寸图片幻灯片

图 14 考虑竖向地震后,面板的开裂位置,损伤时刻和最大裂缝宽度

Figure 14. Crack location, initial cracking moment and maximum crack width of concrete face slab considering vertical earthquake

下载: 全尺寸图片幻灯片

图 15 考虑竖向地震后,面板残余裂缝分布

Figure 15. Distribution of residual crack considering vertical earthquake

下载: 全尺寸图片幻灯片

图 16 坝体震后变形轮廓（放大50倍）

Figure 16. Residual deformations of dam (magnified 50 times)

下载: 全尺寸图片幻灯片

表 1 面板和基岩线弹性模型参数

Table 1 Parameters of slab and bedrock

名称 E/GPa ρ/(kg·m^-3) $ν$

面板 31 2500 0.17
基岩 10 2400 0.25

下载: 导出CSV

表 2 配筋率对面板动力开裂的影响

Table 2 Influences of reinforcement ratio on crack width

配筋率/% 地震中最大裂缝宽度/mm 震后残余裂缝宽度/mm

0.3 15.0 1.7
0.6 11.3 1.4
0.9 0.1 0.0

下载: 导出CSV

表 3 坝前水深对面板动力开裂的影响

Table 3 Influences of water level on crack width

坝前水深/m 地震中最大裂缝宽度/mm 震后残余裂缝宽度/mm

190 20.0 1.2
170 0.1 0.0
140 0.0 0.0

下载: 导出CSV

参考文献(29)

[1]	孔宪京, 邹德高, 刘京茂. 高土石坝抗震安全评价与抗震措施研究进展[J]. 水力发电学报, 2016, 35(7): 1-14. https://www.cnki.com.cn/Article/CJFDTOTAL-SFXB201607001.htm KONG Xian-jing, ZOU De-gao, LIU Jing-mao. Developments in seismic safety evaluation methods and aseismic measures for high rockfill dams[J]. Journal of Hydroelectric Engineering, 2016, 35(7): 1-14. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-SFXB201607001.htm
[2]	陈生水, 阎志坤, 傅中志, 等. 特高面板砂砾石坝结构安全性论证[J]. 岩土工程学报, 2017, 39(11): 1949-1958. doi: 10.11779/CJGE201711001 CHEN Sheng-shui, YAN Zhi-kun, FU Zhong-zhi, et al. Evaluation of safety performance of extremely high slab-faced gravel dams[J]. Chinese Journal of Geotechnical Engineering, 2017, 39(11): 1949-1958. (in Chinese) doi: 10.11779/CJGE201711001
[3]	刘京茂, 孔宪京, 邹德高. 接触面模型对面板与垫层间接触变形及面板应力的影响[J]. 岩土工程学报, 2015, 37(4): 700-710. https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201504019.htm LIU Jing-mao, KONG Xian-jing, ZOU De-gao. Effects of interface models on deformation of interface between slab and cushion layer and slab stress of concrete faced rock fill dam[J]. Chinese Journal of Geotechnical Engineering, 2015, 37(4): 700-710. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201504019.htm
[4]	QU Yong-qian, ZOU De-gao, KONG Xian-jing, et al. Seismic damage performance of the steel fiber reinforced face slab in the concrete-faced rockfill dam[J]. Soil Dynamics and Earthquake Engineering, 2019, 119: 320-330. doi: 10.1016/j.soildyn.2019.01.018
[5]	XU Bin, ZOU De-gao, KONG Xian-jing, et al. Dynamic damage evaluation on the slabs of the concrete faced rockfill dam with the plastic-damage model[J]. Computers and Geotechnics, 2015, 65: 258-265. doi: 10.1016/j.compgeo.2015.01.003
[6]	DAKOULAS P. Longitudinal vibrations of tall concrete faced rockfill dams in narrow canyons[J]. Soil Dynamics and Earthquake Engineering, 2012, 41: 44-58. doi: 10.1016/j.soildyn.2012.05.010
[7]	HILLERBORG A, MODEER M, PETERSON E. Analysis of crack formation and crack growth in concrete by means of fracture mechanics and finite elements[J]. Cement and Concrete Research, 1976, 6: 773-82. doi: 10.1016/0008-8846(76)90007-7
[8]	ARICI Y. Investigation of the cracking of CFRD face plates[J]. Computers and Geotechnics, 2011, 38(7): 905-916. doi: 10.1016/j.compgeo.2011.06.004
[9]	CEN Wei-jun, WEN Lang-sheng, ZHANG Zi-qi, et al. Numerical simulation of seismic damage and cracking of concrete slabs of high concrete face rockfill dams[J]. Water Science and Engineering, 2016, 9(3): 205-211. doi: 10.1016/j.wse.2016.09.001
[10]	BARENBLATT G I. The mathematical theory of equilibrium cracks in brittle fracture[J]. Advances in Applied Mechanics, 1962(7): 55-129.
[11]	DUGDALE D S. Yielding of steel sheets containing slits[J]. Journal of the Mechanics and Physics of Solids, 1960, 8(2): 100-104. doi: 10.1016/0022-5096(60)90013-2
[12]	DAI Qing-li, NG K. 2D cohesive zone modeling of crack development in cementitious digital samples with microstructure characterization[J]. Construction and Building Materials, 2014, 54: 584-595. doi: 10.1016/j.conbuildmat.2013.12.095
[13]	PAN Jian-wen, ZHANG Chu-han, XU Yan-jie, et al. A comparative study of the different procedures for seismic cracking analysis of concrete dams[J]. Soil Dynamics and Earthquake Engineering, 2011, 31(11): 1594-1606. doi: 10.1016/j.soildyn.2011.06.011
[14]	SU Xiang-ting, YANG Zhen-jun, LIU Guo-hua. Finite element modelling of complex 3D static and dynamic crack propagation by embedding cohesive elements in Abaqus[J]. Acta Mechanica Solida Sinica, 2010, 23(3): 271-282. doi: 10.1016/S0894-9166(10)60030-4
[15]	QU Yong-qian, ZOU De-gao, KONG Xian-jing, et al. A novel interface element with asymmetric nodes and its application on concrete-faced rockfill dam[J]. Computers and Geotechnics, 2017, 85: 103-116. doi: 10.1016/j.compgeo.2016.12.013
[16]	CHEN Kai, ZOU De-gao, KONG Xian-jing, et al. Global concurrent cross-scale nonlinear analysis approach of complex CFRD systems considering dynamic impervious panel-rockfill material-foundation interactions[J]. Soil Dynamics and Earthquake Engineering, 2018, 114: 51-68. doi: 10.1016/j.soildyn.2018.06.027
[17]	周墨臻, 张丙印, 张宗亮, 等. 超高面板堆石坝面板挤压破坏机理及数值模拟方法研究[J]. 岩土工程学报, 2015, 37(8): 1426-1432. https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201508014.htm ZHOU Mo-zhen, ZHANG Bing-yin, ZHANG Zong-liang, et al. Mechanisms and simulation methods for extrusion damage of concrete faces of high concrete-faced rockfill dams[J]. Chinese Journal of Geotechnical Engineering, 2015, 37(8): 1426-1432. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201508014.htm
[18]	CHEN Kai, ZOU De-gao, KONG Xian-jing, et al. A novel nonlinear solution for the polygon scaled boundary finite element method and its application to geotechnical structures[J]. Computers and Geotechnics, 2017, 82: 201-210. doi: 10.1016/j.compgeo.2016.09.013
[19]	CHEN Kai, ZOU De-gao, KONG Xian-jing, et al. An efficient nonlinear octree SBFEM and its application to complicated geotechnical structures[J]. Computers and Geotechnics, 2018, 96: 226-245. doi: 10.1016/j.compgeo.2017.10.021
[20]	CHEN Kai, ZOU De-gao, KONG Xian-jing, et al. Elasto- plastic fine-scale damage failure analysis of metro structures based on coupled SBFEM-FEM[J]. Computers and Geotechnics, 2019, 108: 280-294. doi: 10.1016/j.compgeo.2018.12.030
[21]	QU Yong-qian, ZOU De-gao, KONG Xian-jing, et al. A flexible various-scale approach for soil-structure interaction and its application in seismic damage analysis of the underground structure of nuclear power plants[J]. Science China Technological Sciences, 2018, 61(7): 1092-1106. doi: 10.1007/s11431-017-9269-7
[22]	GONG Jin, ZOU De-gao, KONG Xian-jing, et al. An extended meshless method for 3D interface simulating soil-structure interaction with flexibly distributed nodes[J]. Soil Dynamics and Earthquake Engineering, 2019, 125: 1-15.
[23]	DU J, YON J H, HAWKINS N M, et al. Fracture process zone for concrete for dynamic loading[J]. ACI Material Journal 1992, 89(3): 252-258.
[24]	PASTOR M, ZIENKIEWICZ O C. A generalized plasticity, hierarchical model for sand under monotonic and cyclic loading[C]//Proceedings of the 2nd International Symposium on Numerical Models in Geomechanic, 1986, Ghent: 131-149.
[25]	XU Bin, ZOU De-gao, LIU Hua-bei. Three-dimensional simulation of the construction process of the Zipingpu concrete face rockfill dam based on a generalized plasticity model[J]. Computers and Geotechnics, 2012, 43: 143-154. doi: 10.1016/j.compgeo.2012.03.002
[26]	LIU Jing-mao, ZOU De-gao, KONG Xian-jing. A three-dimensional state-dependent model of soil–structure interface for monotonic and cyclic loadings[J]. Computers and Geotechnics, 2014, 61: 166-177. doi: 10.1016/j.compgeo.2014.05.012
[27]	刘晶波, 吕彦东. 结构–地基动力相互作用问题分析的一种直接方法[J]. 土木工程学报, 1998, 31(3): 55-64. https://www.cnki.com.cn/Article/CJFDTOTAL-TMGC199803008.htm LIU Jing-bo, LÜ Yan-dong. A direct method for analysis of dynamic soil-structure interaction[J]. China Civil Engineering Journal, 1998, 31(3): 55-64. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-TMGC199803008.htm
[28]	余翔, 孔宪京, 邹德高, 等. 覆盖层上土石坝非线性动力响应分析的地震波动输入方法[J]. 岩土力学, 2018, 39(5): 1858-1866, 1876. https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201805039.htm YU Xiang, KONG Xian-jing, ZOU De-gao, et al. Seismic wave input method for nonlinear dynamic analysis of earth dam built on overburden[J]. Rock and Soil Mechanics, 2018, 39(5): 1858-1866, 1876. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201805039.htm
[29]	水电工程水工建筑物抗震设计规范：NB 35047—2015[S]. 2015. Code for Seismic Design of Hydraulic Structures of Hydropower Project: NB 35047—2015[S]. 2015. (in Chinese)

施引文献(29)

期刊类型引用(19)

1.	沈书拯，江巍，欧阳晔，谭亿虹，吴剑，陈勇. 基于三维连续-离散耦合模型的门架式双排抗滑桩受力分析. 防灾减灾工程学报. 2025(01): 84-94 . 百度学术
2.	王洪涛，李建华，殷允腾，张红军，宣兆腾，赵明珠. 承压水作用下岩质地层基坑抗突涌安全厚度研究. 建筑结构学报. 2024(06): 213-223 . 百度学术
3.	马勤. 超挖对软土基坑工程稳定性的影响分析. 现代矿业. 2024(04): 121-123 . 百度学术
4.	张瑞程. 新型斜撑-排桩基坑支护体系工程应用分析. 施工技术(中英文). 2024(12): 162-166 . 百度学术
5.	王强，宋歌，廖利，刘国伟，蔡崇德. 全地下空间枢纽深基坑支护结构变形特征控制. 铁路技术创新. 2024(03): 178-185 . 百度学术
6.	赵升峰，汪敏营，邬喜春，李明东. 不同深度基坑中桩加支撑支护特性研究. 江苏建筑. 2024(06): 112-115+129 . 百度学术
7.	赵玉成，代兆宇，夏瑞萌，娄海成. 雄安新区深基坑双排桩支护的实践与探讨. 粉煤灰综合利用. 2024(06): 79-83+95 . 百度学术
8.	蔡红岩. 供水泵站基坑排桩支护局部失效引发连续破坏机理研究. 陕西水利. 2023(01): 100-102 . 百度学术
9.	曾诚. 基坑局部设置双排桩桩长桩径对连续变形控制效果影响研究. 中国新技术新产品. 2023(07): 81-85 . 百度学术
10.	梁文辉. 深基坑双排桩式围护结构支护施工技术研究. 江西建材. 2023(09): 244-245+248 . 百度学术
11.	黄礼明. 某邻近既有建筑的异形软土基坑变形控制. 福建建筑. 2023(11): 84-88 . 百度学术
12.	徐琳，黄程翔，邵根才，杨璞，李宝山. 软土地区基坑超挖对围护结构变形的影响. 建筑安全. 2023(12): 33-36 . 百度学术
13.	魏焕卫，李传斌，种记鑫，郑晓. 局部超挖对内支撑结构的施工效应. 计算机辅助工程. 2023(04): 56-65 . 百度学术
14.	常强. 双排桩治理边坡的数值模拟研究. 黑龙江水利科技. 2022(04): 37-39+179 . 百度学术
15.	郑刚，程雪松，周海祚，张天奇，于晓旋，刁钰，王若展，衣凡，张文彬，郭伟. 岩土与地下工程结构韧性评价与控制. 土木工程学报. 2022(07): 1-38 . 百度学术
16.	唐小军，袁仕贵. 铁路桥台基坑开挖支护方案及稳定性计算. 中国水运. 2022(08): 146-148 . 百度学术
17.	唐小军，袁仕贵. 铁路桥台基坑开挖支护方案及稳定性计算. 中国水运. 2022(15): 146-148 . 百度学术
18.	孙永梅，程海涛，米春荣，孙贝. 土体超挖对深基坑变形影响分析. 山西建筑. 2021(19): 61-63 . 百度学术
19.	曾锦秀. 基于上限法的边坡双排抗滑桩变形计算方法. 福建工程学院学报. 2021(06): 538-544 . 百度学术