基于多功能CPTU测试的基坑开挖扰动深度确定方法

李赞; 刘松玉; 吴恺; 蔡国军; 童立元; 刘文亮

doi:10.11779/CJGE202101021

基于多功能CPTU测试的基坑开挖扰动深度确定方法 English Version

李赞^{1, 2,},
刘松玉^{1, 2, ,},
吴恺^{1, 2},
蔡国军^{1, 2},
童立元^{1, 2},
刘文亮^{1, 2}

1.
东南大学岩土工程研究所,江苏　南京　211189
2.
江苏省城市地下工程与环境安全重点实验室,江苏　南京　211189

基金项目:

国家重点研发计划项目 2016YFC0800201

国家自然科学基金项目 41972269

国家自然科学基金项目 51878157

江苏省交通工程建设局科技项目 CX-2019GC02

详细信息

作者简介:
李赞(1989— ),男,博士研究生,主要从事原位测试技术、软土地基处理等方面的研究工作。E-mail: lizan@seu.edu.cn

通讯作者:
刘松玉, E-mail: liusy@seu.edu.cn

中图分类号: TU43
计量
- 文章访问数: 291
- HTML全文浏览量: 36
- PDF下载量: 152
出版历程
- 收稿日期: 2020-01-07
- 网络出版日期: 2022-12-04
- 刊出日期: 2020-12-31

Determination of disturbance depth due to excavation using multifunctional CPTU tests

LI Zan^{1, 2,},
LIU Song-yu^{1, 2, ,},
WU Kai^{1, 2},
CAI Guo-jun^{1, 2},
TONG Li-yuan^{1, 2},
LIU Wen-liang^{1, 2}

1.
Institute of Geotechnical Engineering, Southeast University, Nanjing210096, China
2.
Jiangsu Key Laboratory of Urban Underground Engineering and Environmental Safety, Southeast University, Nanjing210096, China

摘要

摘要: 基坑开挖卸荷扰动区深度的确定是基坑设计的重要指标。现有的关于扰动区深度的确定大都采用室内试验方法,无法反映土体的原位应力状态。基于多功能孔压静力触探(CPTU)原位测试,对南京、无锡、常州3个基坑工程进行开挖前、后的原位测试,分析了开挖卸荷对锥尖阻力q_t、侧壁摩阻力f_s、摩阻比R_f、电阻率ρ的影响；进一步结合室内卸荷回弹试验,提出了基于原位测试的基坑开挖扰动深度确定方法。结果表明,对于粉土和粉砂地基,开挖卸荷会使锥尖阻力和侧壁摩阻力降低,且随着土体埋深增加,其影响程度逐渐减小；对摩阻比和电阻率影响则较小；基坑开挖扰动区深度为CPTU锥尖阻力q_t衰减20%时对应的深度距坑底开挖面的距离。
- 多功能孔压静力触探 /
- 开挖卸荷 /
- 锥尖阻力 /
- 侧壁摩阻力 /
- 扰动深度
Abstract: The disturbance depth of excavations is an important parameter for their design. Most of the existing methods for determining the disturbed depth are laboratory tests, which cannot reflect the in-situ stress state of soils. The in-situ multifunctional piezocone penetration tests (CPTU) are carried out before and after the excavation of three projects in Nanjing, Wuxi and Changzhou. The influences of unloading induced by excavation on the cone resistance (q_t), sleeve resistance (f_s), friction ratio (R_f) and resistivity (ρ) are also analyzed. Further combining with the laboratory unloading resilience tests, the disturbance depth of excavation based on in-situ tests is proposed. The results show that the unloading effect may reduce the values of cone resistance and sleeve resistance for silt and silty sand foundations, and the corresponding influences on the amplitudes decrease with the increase of depth. However, the excavation has little effect on the friction ratio and resistivity. The disturbance depth of excavation can be determined by the depth beneath the excavation bottom where the cone resistance (q_t) is reduced by 20% as compared to that prior to excavation.
- multifunctional CPTU /
- excavation unloading /
- cone resistance /
- sleeve resistance /
- disturbance depth

HTML全文

0. 引言

柱状节理常见于玄武岩岩体,节理面独特的空间结构分布使得岩体强度表现出强各向异性。已有学者针对柱状节理岩体开展了现场原位试验和大量的相似材料室内试验^[1-5],验证了柱状节理岩体的强度各向异性并总结了材料强度参数的各向异性分布规律。但目前对于柱状节理岩体各向异性本构研究的报道还不多,而已有的一些各向同性模型则不能很好的描述柱状节理岩体的强度各向异性特性。

柱状节理岩体为岩块和结构面的综合体,其表征单元体较室内试验可研究的尺度要大许多,很难通过室内试验来研究柱状节理岩体的变形和强度特征。而原位试验虽然可以在一定程度上反映节理岩体的工程力学特性,但受限于场地、设备、成本等因素,一些物理力学参数难以获取。由于柱状节理岩体具有结构周期性,可以应用数值均匀化方法开展多尺度分析计算,如多重有限元（FE²）或跨尺度有限元离散元耦合（FEM/DEM）。该类方法避免了复杂本构方程的建立,适用于满足均匀化条件的复杂细观结构材料。但多尺度方法计算量巨大,难以直接应用于大型水电工程的计算。因此建立可以准确描述柱状节理岩体力学特性的本构模型具有重要意义。柱状节理岩体显著的各向异性是本构的建模的难点。微结构张量作为材料微观结构的一种度量,可用来描述岩体介质中孔隙或裂纹的空间分布。本文在传统力学模型的基础上引入微结构张量,建立柱状节理岩体各向异性本构模型。

Kanatani^[6-7]最早系统研究了材料空间分布与微结构张量的关系,初步建立了空间分布函数的表达式。Pietruszczak等^[8]通过组合强度参数的空间分布,把经典各向同性准则拓展为各向异性,建立了基于微结构张量的各向异性屈服准则。Gao等^[9]将微结构张量参数引入屈服函数,提出了一种具有交叉各向异性的岩土材料的广义破坏准则,该准则在描述大范围材料的强度随加载方向的变化时具有通用性和鲁棒性。钟世英等^[10]采用多组节理面局部坐标系与整体坐标系的投影关系定义了材料空间的微结构张量表达,将其引入莫尔库仑强度准则,得到了多组节理岩体的各向异性强度准则。

然而现阶段基于微结构张量的各向异性本构模型尚存在不足,或以页岩等层状岩石为研究对象^[8-9],模型不能较好反映柱状节理岩体物理力学试验特性,或模型较为复杂,参数不易确定^[10],实际应用较为困难。

为此,本文提出了一种既保留了数学的严谨性又具有一定实用性的方法。首先,结合微结构张量理论推导了柱状节理岩体各向异性弹塑性本构积分公式,然后基于FLAC^3D对自定义本构进行了二次开发,并结合柱状节理相似结构材料三轴岩石力学试验对本文提出的自定义本构模型进行了验证。应用所开发各向异性弹塑性本构模型对白鹤滩水电站坝基开挖工程进行了数值模拟分析,数值模拟结果与坝基岩体开挖工程监测结果基本相符。研究结果对白鹤滩水电工程建设具有一定应用价值。

1. 各向异性弹塑性本构

1.1 弹性刚度矩阵及空间旋转变换

理想弹性体的本构关系可表示为

$ε = C σ$ 。

(1)

柱状节理岩体（如图1）可简化为横观各向同性材料,由广义胡克定律,在材料空间的局部坐标系下其柔度矩阵可表示为

图 1 全局坐标系中的柱状节理材料

Figure 1. Cylindrical joint materials in global coordinate system

下载: 全尺寸图片幻灯片

$\begin{matrix} C = [\begin{matrix} \frac{1}{E_{1}} & - \frac{ν_{12}}{E_{1}} & - \frac{ν_{13}}{E_{3}} & 0 & 0 & 0 \\ - \frac{ν_{12}}{E_{1}} & \frac{1}{E_{1}} & - \frac{ν_{13}}{E_{3}} & 0 & 0 & 0 \\ - \frac{ν_{13}}{E_{3}} & - \frac{ν_{13}}{E_{3}} & \frac{1}{E_{3}} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1}{G_{13}} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{G_{13}} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{2 (1 + ν_{12})}{E_{1}} \end{matrix}] \\ = {[D^{e}]}^{- 1} \end{matrix}$ ,

(2)

式中, $E_{1}$ 为平行于各向同性面的弹性模量, $E_{3}$ 为垂直于各向同性面的弹性模量, $G_{13}$ 为垂直平面的剪切模量, $ν_{12}$ 为平行面内的泊松比, $ν_{13}$ 为垂直方向的泊松比, $D^{e}$ 为弹性刚度矩阵。

全局坐标系下的柔度矩阵 $C^{'}$ 可由局部坐标系下的柔度矩阵 $C$ 旋转得到

$C^{'} = Q^{T} C Q$ ,

(3)

式中, $Q$ 为一个包含局部坐标轴与全局坐标轴夹角 $α$ 的二阶张量,可表示为

$Q = [\begin{matrix} \cos^{2} α & \sin^{2} α & 0 & \frac{\sin 2 α}{- 2} & 0 & 0 \\ \sin^{2} α & \cos^{2} α & 0 & \frac{\sin 2 α}{2} & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 \\ \sin 2 α & - \sin 2 α & 0 & \cos 2 α & 0 & 0 \\ 0 & 0 & 0 & 0 & \cos α & - \sin α \\ 0 & 0 & 0 & 0 & \sin α & \cos α \end{matrix}]$ 。

(4)

对于横观各向同性材料,需要给出两个夹角,分别为材料的倾向 $α$ 和倾角 $β$ 。后续推导若不加说明,默认在材料空间中的局部坐标系展开。

1.2 各向异性屈服函数及塑性势函数

定义微结构张量 $a_{i j}$ 主轴与材料空间的局部坐标系同轴,三维空间中微结构张量 $a_{i j}$ 可由单位向量 $e_{i}$ 表示为

$\begin{matrix} a_{i j} = a_{1} e_{i}^{(1)} e_{j}^{(1)} + a_{2} e_{i}^{(2)} e_{j}^{(2)} + e_{i}^{(3)} e_{j}^{(3)} \\ = a_{1} m_{i j}^{(1)} + a_{2} m_{i j}^{(2)} + a_{3} m_{i j}^{(3)} \end{matrix}$ ,.

(5)

式中, $a_{i}$ 为微结构张量的主值, $m_{i j}$ 为结构取向张量。

定义广义加载单位向量 $l_{i}$ 为

$l_{i} = \frac{L_{i}}{\sqrt{L_{k} L_{k}}} = \frac{L_{i}}{\sqrt{σ_{i j} σ_{i j}}} = \frac{L_{i}}{\sqrt{tr σ^{2}}}$ ,

(6)

$L_{i} = \sum_{x = 1}^{3} L_{x} e_{i}^{(x)}$ ,

(7)

$L_{x} = \sqrt{σ_{1 x}^{2} + σ_{2 x}^{2} + σ_{3 x}^{2}}$ 。

(8)

考虑微结构张量 $a_{i j}$ 在 $l_{i}$ 方向上的投影,即

$\begin{matrix} η = a_{i j} l_{i} l_{j} = a_{1} \frac{tr(m^{(1)} σ^{2})}{tr σ^{2}} + a_{2} \frac{tr(m^{(2)} σ^{2})}{tr σ^{2}} + \\ a_{3} \frac{tr(m^{(3)} σ^{2})}{tr σ^{2}} \end{matrix}$ 。

(9)

式（9）所定义参数 $η$ 是一个无量纲的标量,综合考虑了材料结构面取向对加载空间的影响,反映了材料强度的方向依赖性,因此可用来描述材料强度的各向异性空间分布。式（9）也可表示为偏量的形式：

$η = \hat{η} (1 + A_{i j} l_{i} l_{j})$ ,

(10)

$\hat{η} = a_{k k} / 3, A_{i j} = (a_{i j} - \hat{η} δ_{i j}) / \hat{η}$ 。

(11)

式（10）最早由Kanatani^[6-7]提出。更一般的,该式可表示为

$\begin{matrix} η = \hat{η} (1 + A_{i j} l_{i} l_{j} + A_{i j k l} l_{i} l_{j} l_{k} l_{l} + \dots) \\ = \hat{η} (1 + A_{i j} l_{i} l_{j} + b_{1} {(A_{i j} l_{i} l_{j})}^{2} + \dots) \end{matrix}$ 。

(12)

本文采用一种修正的莫尔库仑准则^[11]来描述岩石柱体的塑性行为,并考虑各向异性和塑性硬化：

$\begin{matrix} f = \sqrt{{(q K (θ))}^{2} + {(e_{f} c \cos ϕ)}^{2}} + \\ η_{f} h (ξ) (p \sin ϕ - c \cos ϕ) \end{matrix}$ ,

(13)

$\begin{matrix} p = \frac{σ_{k k}}{3},_{}^{} {_{}^{}}_{}^{}_{}^{} q = \sqrt{\frac{3}{2} J_{2}}, \\ J_{2} = s_{i j} s_{i j},_{}^{} {_{}^{}}_{}^{}_{}^{} s_{i j} = σ_{i j} - \frac{σ_{k k}}{3} δ_{i j} \end{matrix}$ 。

(14)

式中 $c$ 为材料黏聚力, $ϕ$ 为材料内摩擦角,可通过试验确定； $e_{f}$ 为子午面上拟合屈服线的双曲偏心率,可通过调整 $e_{f}$ 来调整抗拉强度大小； $η_{f}$ 为微结构张量参数, $h (ξ)$ 为塑性硬化函数, $ξ$ 为塑性内变量； $p$ 为平均应力, $q$ 为偏应力, $J_{2}$ 为应力偏量第二不变量, $θ$ 为Lode角, $K (θ)$ 为关于Lode角的函数。

Lode角函数采用分段形式对屈服面角点进行平滑修正：

$K (θ) = {\begin{matrix} \cos θ - \frac{\sin θ \cos ϕ}{\sqrt{3}} & | θ | \leq θ_{t} \\ K_{1} - K_{2} \sin 3 θ & | θ | > θ_{t} \end{matrix}$ ,

(15)

$\begin{matrix} K_{1} = \frac{1}{3} \cos θ_{t} (3 + \tan θ_{t} \tan 3 θ_{t} + \\ \frac{1}{\sqrt{3}} (\tan 3 θ_{t} - 3 \tan θ_{t}) sign(θ) \sin ϕ) \end{matrix}$ ,

(16)

$K_{2} = \frac{1}{3 \cos 3 θ_{t}} (sign(θ) \sin θ_{t} + \frac{1}{\sqrt{3}} \sin θ_{t} \sin ϕ)$ ,

(17)

式中, $θ_{t}$ 为容许Lode角,当 $| θ | \leq θ_{t}$ 时,与经典Mohr–Coulomb准则一致,当 $| θ | > θ_{t}$ 时,对屈服函数迹线进行平滑处理。

取等效塑性剪应变 ${\bar{γ}}^{p}$ 作为塑性内变量,定义塑性硬化函数为^[12]

$h (ξ) = h_{0} - (h_{1} - h_{0}) \frac{ξ}{B + ξ}$ ,

(18)

$ξ = {\bar{γ}}^{p} = \int d {\bar{γ}}^{p},_{}^{} {_{}^{}}_{}^{} d {\bar{γ}}^{p} = \sqrt{\frac{2}{3} d e_{i j}^{p} d e_{i j}^{p}}$ ,

(19)

式中, $h_{0}$ 为初始屈服阈值, $h_{1}$ 为峰值屈服阈值, $B$ 为描述塑性硬化率大小的参数。

塑性势函数与屈服函数形式一致,即

$g = \sqrt{q^{2} K^{2} (θ) + {(e_{f} c \cos ψ)}^{2}} + η_{f} p \sin ψ$ ,

(20)

式中, $ψ$ 为膨胀角,当 $ψ \neq ϕ$ 时即为非关联流动。

1.3 屈服函数及塑性势函数求导

由于本构数值积分中将涉及到屈服函数及塑性势函数的偏导数,故先行给出偏导数的求解。首先计算屈服函数对应力张量的偏导数,根据式（13）有

$\frac{\partial f}{\partial σ_{i j}} = (\frac{\partial f}{\partial I_{1}} \frac{\partial I_{1}}{\partial σ_{i j}} + \frac{\partial f}{\partial J_{2}} \frac{\partial J_{2}}{\partial σ_{i j}} + \frac{\partial f}{\partial θ} \frac{\partial θ}{\partial σ_{i j}}) + (\frac{\partial f}{\partial η_{f}} \frac{\partial η_{f}}{\partial σ_{i j}})$ ,

(21)

式中,

$\begin{matrix} \frac{\partial I_{1}}{\partial σ_{i j}} = δ_{i j},_{}^{} {_{}^{}}_{}^{} \frac{\partial J_{2}}{\partial σ_{i j}} = s_{i j}, \\ \frac{\partial θ}{\partial σ_{i j}} = \frac{\sqrt{3}}{2 \sqrt{J_{2}^{3}} \cos 3 θ} (\frac{3 J_{3}}{2 J_{2}} s_{i j} - s_{i k} s_{k j} + \frac{2}{3} J_{2} s_{i j}) \end{matrix}$ ,

(22)

$\frac{\partial η_{f}}{\partial σ_{i j}} = \frac{2 η_{0} (A_{k i} σ_{k j} σ_{p q} σ_{p q} - A_{p k} σ_{p q} σ_{k q} σ_{i j})}{{(σ_{m n} σ_{m n})}^{2}}$ ,

(23)

$\frac{\partial f}{\partial I_{1}} = \frac{1}{3} η_{f} h \sin ϕ, \frac{\partial f}{\partial J_{2}} = \frac{1}{2} α_{f} K, \frac{\partial f}{\partial θ} = α_{f} J_{2} \frac{d K}{d θ}$ ,

(24)

$\frac{\partial f}{\partial η_{f}} = h (ξ) (p \sin ϕ - c \cos ϕ)$ ,

(25)

$α_{f} = 3 K / \sqrt{3 J_{2} K^{2} + {(e_{f} c \cos ϕ)}^{2}}$ ,

(26)

$\frac{d K}{d θ} = {\begin{matrix} - \sin θ - \frac{\cos θ \sin ϕ}{\sqrt{3}} & | θ | \leq θ_{t} \\ - 3 K_{2} \cos 3 θ & | θ | > θ_{t} \end{matrix}$ 。

(27)

然后计算塑性势函数对应力张量的偏导数。由于式（20）采取了与屈服函数相一致的形式,因此塑性势函数对应力张量的一阶偏导数只需将内摩擦角 $ϕ$ 替换为膨胀角 $ψ$ ,即

$\frac{\partial g}{\partial σ_{i j}} = (\frac{\partial g}{\partial I_{1}} \frac{\partial I_{1}}{\partial σ_{i j}} + \frac{\partial g}{\partial J_{2}} \frac{\partial J_{2}}{\partial σ_{i j}} + \frac{\partial g}{\partial θ} \frac{\partial θ}{\partial σ_{i j}}) + (\frac{\partial g}{\partial η_{f}} \frac{\partial η_{f}}{\partial σ_{i j}})$ 。

(28)

式中,

$\frac{\partial g}{\partial I_{1}} = \frac{1}{3} η_{f} h \sin ψ, \frac{\partial g}{\partial J_{2}} = \frac{1}{2} α_{g} K, \frac{\partial g}{\partial θ} = α_{g} J_{2} \frac{d K}{d θ}$ ,

(29)

$\frac{\partial g}{\partial η_{f}} = h (ξ) p \sin ψ$ ,

(30)

$α_{g} = 3 K / \sqrt{3 J_{2} K^{2} + {(e_{f} c \cos ψ)}^{2}}$ 。

(31)

1.4 弹塑性刚度矩阵

根据塑性一致性条件,对于式（13）中的屈服函数,其增量形式满足：

$d f = \frac{\partial f}{\partial σ_{i j}} d σ_{i j} + \frac{\partial f}{\partial h} \frac{\partial h}{\partial ξ} d ξ = 0$ 。

(32)

此外,考虑到应力增量和应变增量满足

$d σ_{i j} = D_{i j k l}^{e} d ε_{k l}^{e} = D_{i j k l}^{e} (d ε_{k l} - d ε_{k l}^{p})$ 。

(33)

根据塑性力学,塑性应变增量可由塑性因子 $λ$ 以流动法则表示为

$d ε_{k l}^{p} = d λ \frac{\partial g}{\partial σ_{i j}}$ 。

(34)

等效塑性剪应变 ${\bar{γ}}^{p}$ 也可以流动法则表示为

$d ξ = d {\bar{γ}}^{p} = d λ \sqrt{\frac{2}{3} dev (\frac{\partial g}{\partial σ_{i j}}) dev (\frac{\partial g}{\partial σ_{i j}})}$ 。

(35)

将式（32）～（35）联立,可解得

$d λ = H^{- 1} \frac{\partial f}{\partial σ_{i j}} D_{i j k l}^{e} d ε_{k l}^{p}$ ,

(36)

$\begin{matrix} H = \frac{\partial f}{\partial σ_{i j}} D_{i j k l}^{e} \frac{\partial g}{\partial σ_{k l}} - \\ \frac{\partial f}{\partial h} \frac{\partial h}{\partial ξ} \sqrt{\frac{2}{3} dev (\frac{\partial g}{\partial σ_{i j}}) dev (\frac{\partial g}{\partial σ_{i j}})} \end{matrix}$ ,

(37)

式中,H为塑性硬化模量。

将式（36）回代入式（33）,可得应变增量与应力增量的一致性表达：

$d σ_{i j} = D_{i j k l}^{e} (d ε_{k l} - d λ \frac{\partial g}{\partial σ_{i j}}) = D_{i j k l}^{ep} d ε_{k l}$ ,

(38)

$D_{i j k l}^{ep} = D_{i j k l}^{e} - H^{- 1} D_{i j p q}^{e} \frac{\partial f}{\partial σ_{p q}} \frac{\partial g}{\partial σ_{m n}} D_{m n k l}^{e}$ ,

(39)

式中, $D_{i j k l}^{ep}$ 即为一致性弹塑性刚度矩阵。

2. 数值算例

基于建立的柱状节理岩体的各向异性弹塑性本构模型,结合FLAC^3D软件对该模型进行了二次开发,并对柱状节理相似结构材料三轴岩石力学试验成果进行了模拟。根据国际岩石力学学会（ISRM）推荐标准,圆柱样直径为50 mm,柱体高度为100 mm,网格数量10000个,底面固定,顶面控制加载速度为5×10^-7 m/s,柱体四周施加固定围压4 MPa。

可由不同空间方向的柱状节理岩体的峰值强度进行拟合确定微结构张量系数。图2给出了强度参数关于方向角的微结构张量参数拟合曲线。

图 2 不同拟合阶数下单轴抗压强度随柱体倾角变化规律

Figure 2. Variation of uniaxial compressive strength under different fitting orders with column inclination

下载: 全尺寸图片幻灯片

为了解释固有的各向异性,假设单轴抗压强度 $f_{c}$ 受样本方向的影响,并结合与式（12）相似的分布函数来描述其变化。为了能反映材料各向异性屈服强度,并防止欠拟合和过拟合,本文采用二阶微结构张量来描述柱状节理岩体强度参数的空间分布屈服函数：

$f_{c} = \hat{f} {1 + A (1 - 3 \cos^{2} β) + b {[A (1 - 3 \cos^{2} β)]}^{2}}$ ,

(40)

式中,

$\hat{f} = 12.76 MPa, A = 1.894, b = 0.3526$ 。

(41)

其余计算参数可根据室内试验曲线^[1]数据,应用单纯形最优化方法在MATLAB中构造寻优函数编程计算得到。模型计算参数如表1所示。

表 1 柱状节理岩体基本力学参数表

Table 1. Mechanical parameters of columnar joints

编号	^#1	^#2	^#3	^#4	^#5	^#6	^#7	^#8	^#9	^#10	^#11	^#12	^#13	^#14	^#15	^#16
参数	E₁/MPa	E₂/MPa	G₁₃/MPa	$v_{12}$	$v_{13}$	c/MPa	$φ$ /(°)	$e_{f}$	$\overset{⌢}{η}$	A	b	$h_{0}$	$h_{1}$	$e_{f}$	B	$θ_{t}$
数值	9.81	22.1	4.69	0.22	0.06	0.72	35.2	0.75	0.23	1.89	0.35	0.77	1.21	0.53	0.0018	25°

下载: 导出CSV

| 显示表格

应用表1给出的力学参数,在FLAC^3D中可以调用自定义本构模型来模拟柱状节理岩体的弹塑性变形行为。数值模拟结果如图3和图4所示。

图 3 三轴数值模拟位移云图（倾角75°,围压4 MPa）

Figure 3. Displacement nephogram by triaxial numerical simulation (angle of 75° and confining pressure of 4 MPa)

下载: 全尺寸图片幻灯片

图 4 三轴试验与模拟结果对比（倾角75°,围压4 MPa）^[13]

Figure 4. Comparison between triaxial test and simulated results (angle of 75°and confining pressure of 4 MPa)

下载: 全尺寸图片幻灯片

本文将柱状节理岩体概化为横观各向同性连续体,建立描述岩体力学特性的宏观本构模型。该模型采用了一种考虑抗拉强度和屈服面角点平滑的修正Mohr–Coulomb准则,将其扩展到了空间各向异性,包含一个各向异性标量参数,该参数用应力和结构取向张量的混合不变量表示。

结果表明,柱体倾角为75°时节理岩体的变形具有明显的方向性,模拟所得应力–应变曲线与试验数据^[13]基本吻合。

应该强调的是,在对问题进行分析时,式（18）中硬化特性仅通过调用一个标量参数B来描述,由于式（18）仅限于硬化状态,不包括应变软化分支。当然也可以采用一种更复杂的方法,并假设这些特征受到围压的影响,然而这样将不可避免地生成一组额外的参数及公式。此外,该方法目的在于研究工程尺度柱状节理岩体的物理力学特性,应用了一种折中的处理策略,当研究尺度较小时,则需要进一步考虑对岩体内部的结构面开展精细模拟。

3. 白鹤滩高坝坝基开挖工程数值分析

3.1 白鹤滩坝区工程地质特征及网格建模

白鹤滩水电站坐落于金沙江下游云南巧家县镜内,是金沙江下游干流河段梯级开发的第二个梯级电站。主坝为混凝土双曲拱坝,最大坝高289 m,计划2022年建设完工。白鹤滩水电站坝区属中山峡谷地貌,坝址为单斜地层,主要为二叠系上统峨眉山组玄武岩,如图5,可分为11个岩流层,岩层走向N30～50°E,倾向SE,倾角15°～25°。岩层岩性主要分为斜斑玄武岩、隐晶玄武岩、柱状节理玄武岩、杏仁状玄武岩、角砾熔岩和凝灰岩等。右岸坡顶顶部出露三叠系下统飞仙关组泥质粉砂岩、粉砂质泥岩及少量砂岩^[2]。其中出露 $P_{2} β_{3}$ 层为柱状节理玄武岩,柱体长度为2～3 m,柱体倾角70°～85°,柱体截面主要为不规则五边形和四边形,截面边长为10～20 cm。柱体内部隐节理发育,开挖前节理闭合紧密^[14]。

图 5 拱坝轴线地质剖面示意图

Figure 5. Illustration of engineering geological section of arch axis

下载: 全尺寸图片幻灯片

图6（a）为开挖模拟所用网格模型,坐标原点位于拱坝轴线,建模范围为1300 m×1500 m×1296 m。在模型局部区域,对单元进行了细化剖分,网格平均边长约为1.5 m。采用Hypermesh剖分三维实体,并形成网格文件,然后将数据转换为相应的软件格式,建立了数值计算网格模型,单元总数145万个,节点26万个。坝基开挖边坡及柱状节理玄武岩体为计算分析重点关注的部位。除了层间错动带,该网格模型还考虑了断层F16、F17,以及柱状节理玄武岩层 $P_{2} β_{3}$ 中的层内错动带LS3319、LS3319-1、RS3315和RS336,如图6（b）所示。

图 6 白鹤滩坝基开挖计算网格模型

Figure 6. Calculation grid model for excavation of foundation of …Baihetan Dam

下载: 全尺寸图片幻灯片

3.2 模型参数及工况设置

对模型前后左右4个侧面施加法向约束,对底面施加x,y,z三向约束作为边界条件。模型自下而上共划分有8个地层,模型计算参数根据工程地质条件和岩石力学参数试验等确定,主要岩体力学参数如表2所示。

表 2 主要岩体基本物理力学参数

Table 2. Basic mechanical parameters of main rock mass

岩层分类	$ρ$ /(kg·m^-3)	E/GPa	$ν$	c/MPa	$φ$ /(°)
C₄层面以上	2700	11.0	0.23	1.10	47.7
C₄~C₃层面	2700	13.2	0.23	1.37	51.8
C_3-1~柱状节理顶面	2600	10.2	0.24	1.07	48.0
柱状节理层	2700	10.0	0.25	0.90	47.7
角砾熔岩层	2600	9.0	0.24	1.08	47.2
微新岩岩层	2850	22.0	0.22	2.10	57.2
C₄	2130	0.4	0.35	0.10	28.0
C₃	2130	0.4	0.34	0.10	28.0
C_3-1	2050	0.7	0.35	0.04	28.0
F₁₆	2130	1.5	0.35	0.20	29.2
F₁₇风化上段	2000	0.7	0.35	0.04	15.6
F₁₇风化下段	2130	1.2	0.35	0.05	20.8
F₁₇微新岩体	2130	2.0	0.35	0.20	29.2
LS₃₃₁	2150	1.5	0.35	0.30	30.0
LS₃₃₁₉	2110	0.3	0.35	0.10	26.0
LS_3319-1	2000	0.3	0.35	0.15	27.0
RS₃₃₆	2000	0.3	0.35	0.15	27.0
RS₃₃₁₅	2000	0.3	0.35	0.15	26.0

下载: 导出CSV

| 显示表格

考虑左右岸坝基岩石工程现场实际施工工序,依照表3对边坡开挖过程设置了10个开挖步。

表 3 模型开挖工况设置

Table 3. Setting of model excavation conditions

开挖步编号	开挖台阶高程	现场开挖完成时间
^#1	左岸800 m	2013.10.24
^#2	左岸700 m	2014.08.19
^#3	右岸800 m	2014.11.10
^#4	左岸650 m	2014.12.03
^#5	右岸700 m	2015.03.30
^#6	右岸650 m	2015.06.24
^#7	右岸600 m	2016.01.27
^#8	左岸600 m	2016.06.03
^#9	600 m高程以下	2016.07.26
^#10	保护层	2016.08.14

下载: 导出CSV

| 显示表格

3.3 计算结果分析

图7为白鹤滩开挖模型左岸坝基不同高程处测点位移变化曲线图。左岸坝基建基面岩体位移随开挖进程逐渐增大；开挖完成时,岩体位移随岩体所在高程减小而先增大后减小,坝基底部柱状节理岩体位移大于上部非柱状节理岩体。

图 7 左岸坝基不同高程处测点位移变化曲线图

Figure 7. Curves of displacement change at measurement points at different elevations of dam foundation at left bank

下载: 全尺寸图片幻灯片

图8为考虑柱状节理岩体各向异性时白鹤滩坝基边坡开挖过程位移变化云图。

图 8 考虑柱状节理岩体各向异性模型开挖过程位移变化云图

Figure 8. Cloud charts of displacement change in excavation process considering anisotropic model of columnar jointed rock mass

下载: 全尺寸图片幻灯片

开挖完成时,建基面最大位移超过50 mm,平均位移为25～40 mm；坝基开挖全过程中,最大位移出现在各开挖步开挖梯段的边坡坡脚及开挖台阶面上；岩体变形范围和变形程度随开挖进程逐渐增大；表层岩体的变形较大,变形程度随岩体深度增加而逐渐减小；岩体变形围绕断层和错动带发展,断层和错动带周围的岩体变形显著大于不含断层和错动带的岩体（实际现场开挖中对坝基出露断层采取了置换处理）。

图9为左岸坝基658 m高程处多点变位计Mzjc-2监测所得位移时程曲线,以及分别采用各向异性模型和各向同性模型计算所得的变形曲线。从图中可以看到,相比于只采用各向同性模型,采用各向异性模型的位移计算结果与现场监测结果更为接近。

图 9 左岸坝基多点变位计Mzjc-2监测结果与计算结果对比

Figure 9. Comparison between monitoring and calculated results of Mzjc-2 multi-point displacement meter

下载: 全尺寸图片幻灯片

表4为坝基开挖完成时刻左岸坝基柱状节理试验区及左岸坝肩下游边坡处多点变位计变形监测结果与数值模型计算结果的对比。

表 4 坝基开挖结束时监测变形与计算位移对比

Table 4. Comparison of monitored deformations and calculated displacements at completion end of excavation of dam foundation

测点编号	监测结果/mm	各向异性模型		各向同性模型
测点编号	监测结果/mm	计算值/mm	误差/%	计算值/mm	误差/%
Mzjc-2	25.4	23.7	6.7	19.2	24.4
Mzjc-3	24.9	21.8	12.4	18.0	27.7
M_ZBJ-3	6.6	5.6	15.2	4.9	25.8

下载: 导出CSV

| 显示表格

结果表明：各向异性模型位移计算误差明显小于各向同性模型,这在柱状节理区域表现更为显著。

为生成质量较高的三维网格模型,本文所用开挖模型中开挖梯段高差（50～100 m）与实际工程施工台阶（5～10 m）存在较大差距,这在模拟计算结果中主要体现在变形的时间节点与实际监测结果存在一定出入。另一方面,本文计算所采用的各向异性弹塑性本构模型不能完整体现节理岩体卸荷松弛的时效性。在后续的工作中,考虑引入各向异性时效损伤以本文的相关工作进一步深化完善。

4. 结论

（1）结合微结构张量理论扩展了经典莫尔库仑模型,完整推导了柱状节理岩体的各向异性弹塑性本构模型。考虑各向异性的本构模型可以直观反映岩体强度的各向异性空间分布,计算所需参数均可由试验结果计算得到。

（2）结合FLAC^3D二次开发接口对自定义本构模型进行了开发。对所开发模型进行了数值验证,数值计算结果与柱状节理相似材料室内三轴压缩试验结果基本一致,表明该模型可以体现柱状节理材料的空间各向异性特性。

（3）应用所开发本构模型对白鹤滩水电站坝基边坡开挖过程进行了模拟。计算结果表明,与未考虑各向异性相比,考虑柱状节理岩体强度各向异性时坝基边坡位移更接近于现场多点位移计监测结果,本文研究工作对白鹤滩水电工程建设有一定参考价值。

图 1 卸荷影响深度

Figure 1. Depths affected by unloading

下载: 全尺寸图片幻灯片

图 2 粉砂（南京场地）卸荷回弹试验结果

Figure 2. Results of unloading rebound tests on silty sand at Nanjing site

下载: 全尺寸图片幻灯片

图 3 粉土（无锡场地）卸荷回弹试验结果

Figure 3. Results of unloading rebound tests on silty sand at Wuxi site

下载: 全尺寸图片幻灯片

图 4 CPTU现场测试

Figure 4. In-situ piezocone penetration tests

下载: 全尺寸图片幻灯片

图 5 南京场地电阻率CPTU测试结果

Figure 5. Results of piezocone penetration tests on resistivity at Nanjing site

下载: 全尺寸图片幻灯片

图 6 无锡场地电阻率CPTU测试结果

Figure 6. Results of piezocone penetration tests on resistivity at Wuxi site

下载: 全尺寸图片幻灯片

图 7 常州场地CPT测试结果

Figure 7. Results of cone penetration tests at Changzhou site

下载: 全尺寸图片幻灯片

图 8 南京场地锥尖阻力衰减

Figure 8. Decrease of cone resistance at Nanjing site

下载: 全尺寸图片幻灯片

图 9 无锡场地锥尖阻力衰减

Figure 9. Decrease of cone resistance at Wuxi site

下载: 全尺寸图片幻灯片

图 10 常州场地锥尖阻力衰减

Figure 10. Decrease of cone resistance at Changzhou site

下载: 全尺寸图片幻灯片

表 1 现有开挖扰动区深度确定的室内试验方法

Table 1 Existing laboratory test methods for determination of disturbance depth

文献	试验类型	极限卸荷比R_u	扰动区深度h_u
秦爱芳等^[11]	直剪试验	0.67	0.50D
潘林有等^[12]		0.64	0.56D
程玉梅^[15]		0.66	0.52D
秦爱芳等^[11]	卸荷回弹试验	0.80	0.25D
张淑朝等^[13]		0.90	0.11D
邓指军等^[14]		0.81	0.23D
潘林有等^[16]		0.67	0.25D

下载: 导出CSV

表 2 场地主要土层的物理力学性质

Table 2 Physical and mechanical properties of soils

场地名称	土层	层厚/m	$γ$ /(kN·m^-3)	$γ^{'}$ /(kN·m^-3)	G_s	w/%	e₀
南京	淤泥质粉质黏土	2.4	17.8	8.07	2.68	38.3	1.081
南京	粉砂	4.6	19.1	9.58	2.65	25.0	0.723
无锡	粉质黏土	5.0	19.5	9.45	2.71	30.3	0.810
无锡	粉土	6.4	19.4	9.64	2.71	27.1	0.774
常州	素填土	2.3	19.1	9.26	2.75	31.4	0.890
	淤泥质粉质黏土	3.6	18.6	8.51	2.72	38.1	1.020
	⑤₁粉土	4.7	19.4	9.48	2.70	29.0	0.794
	粉砂	4.4	19.7	9.80	2.69	26.0	0.724
	⑤₃粉土	10.0	19.0	9.14	2.71	31.0	0.870
注：表中, $γ$ 为天然土层的重度, $γ^{'}$ 为浮重度,G_s为土粒相对质量密度,w为含水率,e₀为原状土孔隙比。

下载: 导出CSV

表 3 室内试验与原位测试的结果对比

Table 3 Comparison of results between laboratory and in-situ tests

场地名称	土层	D/m	h_u/D		R_u
场地名称	土层	D/m	原位	室内	原位	室内
南京	粉砂	5.8	0.21	0.25	0.83	0.80
无锡	粉土	8.0	0.34	0.25	0.76	0.80
常州	粉砂	11.5	0.26	—	0.79	—
常州	⑤₃粉土	18.0	0.22	—	0.82	—

下载: 导出CSV

参考文献(23)

[1]	孙钧, 周健, 龚晓南, 等. 受施工扰动影响土体环境稳定理论与变形控制[J]. 同济大学学报(自然科学版), 2004, 32(10): 1261-1269. doi: 10.3321/j.issn:0253-374X.2004.10.001 SUN Jun, ZHOU Jian, GONG Xiao-nan, et al. Theory of soil environment stability and deformation control under influence of construction disturbance[J]. Journal of Tongji University (Natural Science Edition), 2004, 32(10): 1261-1269. (in Chinese) doi: 10.3321/j.issn:0253-374X.2004.10.001
[2]	孙钧. 城市环境土工学[M]. 上海: 上海科学技术出版社, 2005. SUN Jun. Geoenvironmental Engineering[M]. Shanghai: Shanghai Science & Technology Press, 2005. (in Chinese)
[3]	ROCHELLE P L, SARRAILH J, TAVENAS F, et al. Causes of sampling disturbance and design of a new sampler for sensitive soils[J]. Canadian Geotechnical Journal, 1981, 18(1): 52-66. doi: 10.1139/t81-006
[4]	SANTAGATA M C, GERMAINE J T. Sampling disturbance effects in normally consolidated clays[J]. Journal of Geotechnical and Geo-environmental Engineering, 2002, 128(12): 997-1006. doi: 10.1061/(ASCE)1090-0241(2002)128:12(997)
[5]	HONG Z S, HAN J. Evaluation of sample quality of sensitive clay using intrinsic compression concept[J]. Journal of Geotechnical and Geo-environmental Engineering, 2007, 133(1): 83-90. doi: 10.1061/(ASCE)1090-0241(2007)133:1(83)
[6]	CARROLL R, LONG M. Sample disturbance effects in silt[J]. Journal of Geotechnical and Geo-environmental Engineering, 2017, 143(9): 04017061. doi: 10.1061/(ASCE)GT.1943-5606.0001749
[7]	李忠超. 软黏土中轨道交通地下深开挖工程变形及稳定性研究[D]. 杭州: 浙江大学, 2015. LI Zhong-chao. Deformation and Stability Investigation of Underground Excavation of Subway Transit System in Soft Clay[D]. Hangzhou: Zhejiang University, 2015. (in Chinese)
[8]	徐永福, 陈建山, 傅德明. 盾构掘进对周围土体力学性质的影响[J]. 岩石力学与工程学报, 2003, 22(7): 1174-1179. doi: 10.3321/j.issn:1000-6915.2003.07.024 XU Yong-fu, CHEN Jian-shan, FU De-ming. Effect of shield tunneling on mechanical properties of soils[J]. Chinese Journal of Rock Mechanics and Engineering, 2003, 22(7): 1174-1179. (in Chinese) doi: 10.3321/j.issn:1000-6915.2003.07.024
[9]	邓永锋, 刘松玉, 季署月. 取样扰动对固结系数的影响研究[J]. 岩土力学, 2007, 28(12): 2687-2690. doi: 10.3969/j.issn.1000-7598.2007.12.040 DENG Yong-feng, LIU Song-yu, JI Shu-yue. Effect of sample disturbance on consolidation coefficient[J]. Rock and Soil Mechanics, 2007, 28(12): 2687-2690. (in Chinese) doi: 10.3969/j.issn.1000-7598.2007.12.040
[10]	TERZAGHI, K. Theoretical Soil Mechanics[M]. New York: John Wiley & Sons, Inc, 1943.
[11]	秦爱芳, 胡中雄, 彭世娟. 上海软土地区受卸荷影响的基坑工程被动区土体加固深度研究[J]. 岩土工程学报, 2008, 30(6): 935-940. doi: 10.3321/j.issn:1000-4548.2008.06.025 QIN Ai-fang, HU Zhong-xiong, PENG Shi-juan. Depth of soil stabilization in passive area of foundation pits for Shanghai soft clay[J]. Chinese Journal of Geotechnical Engineering, 2008, 30(6): 935-940. (in Chinese) doi: 10.3321/j.issn:1000-4548.2008.06.025
[12]	潘林有, 程玉梅, 胡中雄. 卸荷状态下黏性土强度特性试验研究[J]. 岩土力学, 2001, 22(4): 490-493. doi: 10.3969/j.issn.1000-7598.2001.04.030 PAN Lin-you, CHENG Yu-mei, HU Zhong-xiong. Experimental study on the shear strength of clay under the unloading state[J]. Rock and Soil Mechanics, 2001, 22(4): 490-493. (in Chinese) doi: 10.3969/j.issn.1000-7598.2001.04.030
[13]	张淑朝, 张建新, 张阳, 等. 基坑开挖卸荷土体回弹实验研究[J]. 岩土工程学报, 2008, 30(增刊1): 426-429. https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC2008S1092.htm ZHANG Shu-chao, ZHANG Jian-xin, ZHANG Yang, et al. Experimental study on rebound of clay due to foundation excavatiom[J]. Chinese Journal of Geotechnical Engineering, 2008, 30(S1): 426-429. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC2008S1092.htm
[14]	邓指军, 贾坚. 地铁车站深基坑卸荷回弹影响深度的试验[J]. 城市轨道交通研究, 2008, 11(3): 52-55. doi: 10.3969/j.issn.1007-869X.2008.03.014 DENG Zhi-jun, JIA Jian. On unloading resilience depth in deep station foundation pit[J]. Urban Mass Transit, 2008, 11(3): 52-55. (in Chinese) doi: 10.3969/j.issn.1007-869X.2008.03.014
[15]	程玉梅. 基坑开挖卸荷影响深度的确定[J]. 大坝观测与土工测试, 2001, 25(5): 36-38. doi: 10.3969/j.issn.1671-3893.2001.05.012 Determination of the influence depth due to excavation unloading[J]. Dam Observation and Geotechnical Tests, 2001, 25(5): 36-38. (in Chinese) doi: 10.3969/j.issn.1671-3893.2001.05.012
[16]	潘林有, 胡中雄. 深基坑卸荷回弹问题的研究[J]. 岩土工程学报, 2002, 24(1): 101-104. doi: 10.3321/j.issn:1000-4548.2002.01.023 PAN Lin-you, HU Zhong-xiong. Experimental study on the resilience of pit under unloading[J]. Chinese Journal of Geotechnical Engineering, 2002, 24(1): 101-104. (in Chinese) doi: 10.3321/j.issn:1000-4548.2002.01.023
[17]	沈珠江. 原状取土还是原位测试——土质参数测试技术发展方向刍议[J]. 岩土工程学报, 1996, 18(5): 90-91. doi: 10.3321/j.issn:1000-4548.1996.05.017 SHEN Zhu-jiang. Soil sampling on site or in-situ testing—A preliminary research on development of soil parameters testing technology[J]. Chinese Journal of Geotechnical Engineering, 1996, 18(5): 90-91. (in Chinese) doi: 10.3321/j.issn:1000-4548.1996.05.017
[18]	BUDHU M, WU C S. Numerical analysis of sampling disturbances in clay soils[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 1992, 16(7): 467-492. doi: 10.1002/nag.1610160702
[19]	刘松玉, 蔡国军, 童立元. 现代多功能CPTU技术理论与工程应用[M]. 北京: 科学出版社, 2013. LIU Song-yu, CAI Guo-jun, TONG Li-yuan. Research and engineering application of the digital multifunction piezocone penetration test system (CPTU)[M]. Beijing: Science Press, 2013. (in Chinese)
[20]	陈云敏, 胡琦, 陈仁朋. 杭州地铁湘湖车站基坑坍塌引起的基底土深层扰动与沉降分析[J]. 土木工程学报, 2014, 47(7): 110-117. https://www.cnki.com.cn/Article/CJFDTOTAL-TMGC201407016.htm CHEN Yun-min, HU Qi, CHEN Ren-peng. Soil disturbance by the collapse of retaining wall for a pit excavation and the induced additional settlement: a case study of Hangzhou Metro Xianghu Station[J]. China Civil Engineering Journal, 2014, 47(7): 110-117. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-TMGC201407016.htm
[21]	李建民. 超深超大基坑回弹变形计算方法的试验研究[D]. 北京: 中国建筑科学研究院, 2010. LI Jian-min. Research on the Calculation Method of Resilient Deformation in Ultra-Deep Foundation Pit[D]. Beijing: China Academy of Building Research, 2010. (in Chinese)
[22]	刘国华, 王振宇, 黄建平. 土的电阻率特性及其工程应用研究[J]. 岩土工程学报, 2004, 26(1): 83-87. https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC200401016.htm LIU Guo-hua, WANG Zhen-yu, HUANG Jian-ping. Research on electrical resistivity feature of soil and it's application[J]. Rock and Soil Mechanics, 2004, 26(1): 83-87. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC200401016.htm
[23]	韩立华, 刘松玉, 杜延军. 温度对污染土电阻率影响的试验研究[J]. 岩土力学, 2007, 28(6): 1151-1155. doi: 10.3969/j.issn.1000-7598.2007.06.016 HAN Li-hua, LIU Song-yu, DU Yan-jun. Experiment study on effect of temperature on electrical resistivity of contaminated soils[J]. Rock and Soil Mechanics, 2007, 28(6): 1151-1155. (in Chinese) doi: 10.3969/j.issn.1000-7598.2007.06.016

施引文献(35)

期刊类型引用(24)

1.	夏晶晶，贺姣姣. 基于直剪试验土石混合体路基的稳定性问题研究. 山西建筑. 2025(01): 125-128+172 . 百度学术
2.	郜力君，郝捷. 含石率及坡度对土石混合路堑边坡稳定性的影响分析. 交通世界. 2025(Z2): 154-156 . 百度学术
3.	邹弈，朱碧堂，吴颖彪，周宇航. 断层破碎带黏性夹泥岩体与基岩接触界面大型直剪试验研究. 华东交通大学学报. 2025(01): 45-51 . 百度学术
4.	李利萍，余泓浩，李秋雨，潘一山. 砂岩不同含水特性对超低摩擦效应影响试验研究. 力学学报. 2025(03): 687-700 . 百度学术
5.	程虎，李重情，穆朝民. 冻结温度对不同粒径冻土石混合体劈裂特性的影响. 煤矿安全. 2024(01): 160-166 . 百度学术
6.	代兴先，宋杨，张志彬，齐子怡，刘棋瑞. 土石混合料力学特性和颗粒破碎研究. 河北水利电力学院学报. 2024(01): 1-8 . 百度学术
7.	李慎刚，石云方，刘晋宁，蒋琛. 碎石土路基填料压实及渗透特性. 工程科学学报. 2024(05): 918-926 . 百度学术
8.	SHI Yunfang，LI Shengang，JIANG Chen，LIU Jinning. Gravel hardness effect on compaction characteristics of gravelly soil. Journal of Mountain Science. 2024(04): 1432-1443 . 必应学术
9.	谢周州，赵炼恒，李亮，黄栋梁，张子健，周靖. 基于振动台试验的不同含石率土-石混合体边坡地震动响应差异性研究. 岩土力学. 2024(08): 2324-2337 . 百度学术
10.	王楠楠，高霞，张吉哲，张保勇，吴强. 基于平行黏结模型的含瓦斯水合物煤体宏细观力学性质研究. 煤炭学报. 2024(S1): 314-326 . 百度学术
11.	魏东旭，彭雄志，张佳，范帅，冯瑨. 石灰岩碎石土大型直剪试验研究. 勘察科学技术. 2024(04): 1-5+43 . 百度学术
12.	李刚，尹小涛. 基于数字化施工的山区公路土石混合弃渣工程尺度参数测定方法研究. 水利与建筑工程学报. 2024(06): 76-82 . 百度学术
13.	陈小翔. 碎石土路基填料压实及渗透特性研究. 江西建材. 2024(10): 52-54 . 百度学术
14.	杨忠平，李进，刘浩宇，张益铭，刘新荣. 土石混合体-基岩界面剪切力学特性块石尺寸效应. 岩土力学. 2023(04): 965-974 . 百度学术
15.	刘建平，周花玉，何天奎，余镜南，张坤，潘玉丛，刘泉声. 含根量对根–土石复合体的抗剪强度影响试验初探——以垂丝海棠为例. 岩石力学与工程学报. 2023(S1): 3618-3628 . 百度学术
16.	刘飞禹，孔剑捷，姚嘉敏. 含石量和压实度对格栅-土石混合体界面剪切特性的影响. 岩土工程学报. 2023(05): 903-911 . 本站查看
17.	韩志洋，曹志翔，黄开放. 基于离散元模拟的土石混合体剪切与变形特性研究. 中国农村水利水电. 2023(05): 238-244 . 百度学术
18.	唐丽云，黄涛，汪卫兵，金龙，孙强，李国玉，罗滔. 冻融循环下土石混合体-混凝土界面剪切特性及孔隙结构演化特征试验研究. 中南大学学报(自然科学版). 2023(05): 1954-1969 . 百度学术
19.	张俊云，张乐，高福洲，唐永吉，何卓岭，王鹰. 干湿循环下红层土石混合料强度及变形特性的试验研究. 西南交通大学学报. 2023(06): 1394-1404 . 百度学术
20.	程晓颖，乔婷，秦建敏，季顺迎. 复杂形态碎石颗粒的三维离散元模拟及试验验证. 计算机辅助工程. 2023(04): 40-47 . 百度学术
21.	YANG Zhong-ping，LI Shi-qi，TIAN Xin，HU Yuan-xin，LI Wan-kun. Cumulative damage effect on debris slopes under frequent microseisms. Journal of Mountain Science. 2022(03): 781-797 . 必应学术
22.	崔熙灿，张凌凯，王建祥. 高堆石坝砂砾石料的细观参数反演及三轴试验模拟. 农业工程学报. 2022(04): 113-122 . 百度学术
23.	宋颖能. 朱家涧水库大坝土石方及渗控工程实施探讨. 内蒙古煤炭经济. 2022(08): 160-162 . 百度学术
24.	覃国强. 不同含石率及坡度条件下土石混合体路堑高边坡失稳机制分析. 福建交通科技. 2022(06): 21-25 . 百度学术