Phase-field modeling of frost propagation of cracks for rock mass under frost action
-
摘要: 寒区岩体冻结时，裂隙水相变成冰体积膨胀，冰-岩相互作用常引起裂隙冻胀扩展。为探索冻结岩体裂隙冻胀扩展的模拟和预测方法并深入认识不同条件下裂隙冻胀扩展规律，引入以弥散化标量场表征裂隙的相场法描述裂隙冻胀扩展，利用等效热膨胀系数法施加裂隙冰相变膨胀，在COMSOL Multiphysics软件平台实现对岩体裂隙冰-岩相互作用力学场控制方程和裂隙冻胀扩展相场演化控制方程的求解。基于冻结岩体冻胀扩展试验，对多种工况下裂隙冻胀扩展进行相场法模拟。模拟所得不同倾角、外荷载作用下单裂隙冻胀扩展，及不同岩桥倾角、裂隙夹角工况时双裂隙冻胀扩展均与试验获得裂隙冻胀扩展形态接近，说明文中建立的相场法模型能够较准确地模拟预测冻结岩体的裂隙冻胀扩展。受外荷载作用时，单裂隙冻胀扩展向荷载作用方向偏转。对不同岩桥倾角的双裂隙，因裂隙之间相互作用，裂隙内尖端冻胀扩展总是向另一条裂隙方向偏转，裂隙外尖端则产生共面冻胀扩展。对不同夹角的双裂隙，平行双裂隙会形成两条独立呈“蝶状”的冻胀扩展裂隙，而具有倾斜交角的双裂隙会形成交于一点的辐射状冻胀扩展裂隙。Abstract: As the rock masses in cold regions freeze, water in cracks turns into ice and expands in volume, and the mechanical interaction between ice and rock may lead to the frost propagation of cracks. To study the prediction method for the frost propagation of cracks and further cognize the laws of the frost propagation under different conditions, the phase-field model which represents cracks in a diffusive way with a scalar field is introduced to simulate the frost propagation, and the method for the equivalent thermal expansion coefficient is utilized to simulate the volume expansion of ice in cracks due to phase transition. Moreover, the governing equations for stress field for ice-rock interaction and the governing equations of phase-field evolution for the frost propagation of cracks are solved through the COMSOL Multiphysics software. Numerical simulations with the phase-field model on the frost propagation are conducted based on a series of frost propagation experiments on the rock masses. The numerical results are similar to the experimental ones for both the frost propagation of a single crack under different dip angles and external loads, and the frost propagation of double cracks under different dip angles of rock bridge and different intersection angles. It is indicated that the phase-field model established can accurately simulate the frost propagation of cracks in the rock masses. Furthermore, when subjected to an external load, the frost propagation of a single crack deflects towards the direction of the load. For the double cracks with different dip angles of rock bridge, the frost propagation of inner tips always deflects towards the adjacent crack due to the interaction of two cracks, while the outer tips propagate approximately along the coplanar direction. For double cracks with different intersection angles, two independent new frost cracks in a butterfly shape will form when the two cracks are parallel to each other, while new frost cracks in a radiation shape will form for the double cracks with an inclined intersection angle.
-
Keywords:
- rock mass /
- frost propagation /
- phase-field method /
- single crack /
- double cracks
-
0. 引言
中国寒区分布广泛，寒区含水岩体冻结时，裂隙水相变成冰体积膨胀，产生作用于裂隙壁的冻胀力，冻胀力增大常引起裂隙冻胀扩展，甚至岩体断裂，诱发寒区岩体工程冻害问题[1]。因此，寒区岩体冻融劣化与裂隙冻胀扩展演化问题引起了广泛关注。
冻胀力是诱发岩体冻融劣化和驱动裂隙冻胀扩展的力学基础。Davidson等[2]基于光弹性效应测量裂隙水冻结成冰时产生的冻胀力，黄诗冰等[3]和乔趁等[4]均采用薄膜压力传感器测量了岩体裂隙水冰相变过程中的冻胀力变化。尽管学者们进行了一系列测量裂隙内冻胀力的试验[5]，确定冻胀力的量值仍存在困难，因为它受冻结速率、未冻水含量、冻结方向[6]、裂隙几何尺寸、岩石力学性质[7]等诸多因素的影响。在岩体裂隙冻胀力的理论计算方面，Tan等[8]基于复变函数方法建立了椭圆形裂隙冻胀力解析解，刘泉声等[9]采用等效热膨胀系数法热力耦合模拟了含水裂隙冻结时的冻胀力，并建立柱形裂隙中冻胀力演化模型，研究了考虑和忽略水分迁移时冻胀力量值变化[10]。
当冻胀力超过岩石抗拉强度，则会导致岩体裂隙冻胀扩展。在裂隙冻胀扩展的试验研究方面，Huang等[11]研究了不同倾角、不同长度单裂隙的冻胀扩展模式以及因裂隙扩展导致的岩体强度损失。此外，Wang等[12]利用声发射和CT扫描技术，研究了含不同夹角双裂隙岩体岩桥的冻胀断裂特征。李平等[13]进行了含双裂隙岩体冻胀试验和单轴压缩试验，研究了不同岩桥倾角双裂隙的冻胀扩展规律和岩体强度劣化规律。
在冻结岩体裂隙冻胀扩展的理论计算方面，早期主要结合断裂力学分析不同裂隙形态时冻胀扩展的产生条件[14]。Huang等[15]基于复变函数方法探讨了椭圆裂隙时最大拉应力准则和应力强度因子准则对岩体裂隙冻胀扩展的适用性。近年来，因可同时将冻结岩体的热传递、水冰相变、冰-岩相互作用等过程进行耦合，THM多场耦合数值方法在冻结裂隙岩体的研究中应用较广[16]，但基于常规有限元方法的THM耦合在模拟裂隙冻胀扩展时遇到困难。于是，学者们开始利用扩展有限元方法处理岩体裂隙冻胀扩展问题[17]，但扩展有限元往往需要特有的算法和计算软件平台[18]。
相场法引入弥散式连续标量场表征裂隙[19]，避免了裂隙形成和扩展引起的不连续问题给数值计算带来的极大困难，裂隙扩展后无需重新划分网格。且相场法基于能量变分原理，由能量变化驱动相场演化，即裂隙扩展，无需预设裂隙扩展方向。因此，相场法在岩石等脆性材料断裂研究中[20-21]，尤其是岩石的水力劈裂方面[22-23]，得到了深入的应用。冻结岩体裂隙冻胀扩展与岩石水力劈裂均为岩石在内部压力，即冰压或水压作用下，产生断裂的现象，具有相似性。因此，可将相场法应用于冻结岩体裂隙冻胀扩展的研究，以发挥相场法在描述裂隙扩展方面的独特优势。而目前在冻结多孔介质力学分析中，相场法理念多被用来描述水冰相变过程[24]，尚未见用于描述裂隙冻胀扩展。
因此，本文利用相场法描述裂隙冻胀扩展，采用等效热膨胀系数法施加裂隙冰相变膨胀，在COMSOL Multiphysics软件平台实现冻结岩体裂隙冰-岩相互作用力学场控制方程和裂隙冻胀扩展相场演化控制方程的求解。基于含单裂隙、双裂隙岩体冻胀扩展试验，对不同倾角、外荷载作用下岩体单裂隙冻胀扩展以及不同岩桥倾角、裂隙夹角条件下岩体双裂隙冻胀扩展进行相场法模拟，以探索冻结岩体裂隙冻胀扩展模拟和预测方法，深入认识冻结岩体裂隙冻胀扩展规律。
1. 岩体裂隙冻胀扩展的控制方程
1.1 裂隙冻胀扩展的相场法描述
研究对象岩体Ω以∂Ω为外边界，冻结过程中，裂隙水相变成裂隙冰Ω i，体积膨胀，导致岩体内初始裂隙发生冻胀扩展，产生扩展裂隙Γ，见图 1（a）。在相场法中，为描述裂隙扩展引入标量场ϕ，场量ϕ∈[0,1]，ϕ=0表示未发生裂隙扩展，而ϕ=1表示裂隙扩展，新裂隙产生，见图 1（b）。标量l0为长度参数，控制着相场过渡区域范围，反映了裂隙扩展的宽度。l0越大，相场中裂隙扩展弥散区域越宽；而l0越小，相场中裂隙扩展区域越尖锐。根据相场ϕ的定义，单位体积内冻胀扩展裂隙表面密度可表示为[19]
γ(ϕ,∇ϕ)=ϕ22l0+l02|∇ϕ|2。 (1) 设岩石材料的临界断裂能密度为Gc，裂隙冻胀扩展所需的断裂能Ψf可近似为[20]
Ψf=∫ΓGcdS=∫ΩGc(ϕ22l0+l02|∇ϕ|2)dΩ。 (2) 1.2 相场演化的控制方程
岩体冻结过程中，冻胀扩展裂隙Γ的弹性表面能由岩体的弹性应变能转化而来，即岩体的弹性应变能驱动着相场演化[21]。冻胀扩展裂隙为Ⅰ型裂纹，仅由张拉荷载驱动，即弹性应变能的张拉部分驱动着相场演化，压缩部分不会导致裂隙扩展，故需将岩体的弹性应变能分解为张拉和压缩两部分。本文采用Miehe等[19]提出的分解方法，首先将应变张量分解如下：
ε±=d∑a=1⟨εa⟩±na⊗na (3) 式中：ε+，ε−分别为拉应变张量和压应变张量；d为几何维度，平面问题中d=2；εa为主应变，na为主应变的方向向量；⟨εa⟩±=(εa±|εa|)/2。然后，利用分解后的应变张量将弹性应变能密度分解为[19]
ψ±ε(ε)=λ2⟨tr(ε)⟩2±+μtr(ε2±)。 (4) 式中：λ，μ为拉梅常数。
假定相场演化只影响张拉部分的弹性应变能密度，可引入式（5）表征裂隙冻胀扩展时岩体刚度衰减和弹性应变能密度变化[25]：
ψε(ε)=[(1−k)(1−ϕ)2+k]ψ+ε(ε)+ψ−ε(ε)。 (5) 式中：k为为防止ϕ→1时ψ+ε(ε)项消失产生数值计算奇异而引入的模型参数，文中取k=1.0×10-9[20]。
根据Griffith能量平衡理论，裂隙冻胀扩展过程中，岩体Ω的总势能由弹性应变能、断裂能和外力做功构成[26]：
L=∫Ω{[(1−k)(1−ϕ)2+k]ψ+ε(ε)+ψ−ε(ε)}dΩ+ ∫ΩGc(ϕ22l0+l02|∇ϕ|2)dΩ−∫∂Ωf⋅udS。 (6) 式中：f为面力；u为位移向量。
本文仅考虑面力做功，忽略体力做功。根据弹性力学的变分原理，令能量函数L的一阶变分δL=0[27]，可得应力场、相场的控制方程和边界条件分别为[25]
Div(σ)=0 ，[2l0(1−k)ψ+ϵGc+1]ϕ−l20∇2ϕ=2l0(1−k)ψ+ϵGc ，} (7) σ⋅m=f,∇ϕ⋅m=0。} (8) 式中：m为边界外法向量；柯西应力张量σ可表示为
σ=[(1−k)(1−ϕ)2+k][λ⟨tr(ε)⟩+I+2με+]+λ⟨tr(ε)⟩−I+2με−。 (9) 式中：I为单位阵。
1.3 裂隙冰相变膨胀的控制方程
冻结过程中，岩体裂隙冻胀扩展是裂隙水相变成冰体积膨胀所致。将裂隙冰Ω i（图 1）视为线弹性材料，考虑相变体积膨胀，其物理方程可表示为
εij=1Ei[(1+νi)σij−νiσkkδij]+εfδij。 (10) 式中：Ei，νi分别为裂隙冰的弹性模量和泊松比；εij，σij为分别为裂隙冰Ω i域的应变张量和应力张量，εf为裂隙水冰相变膨胀线应变，δij为Kronecker符号。
不考虑体力作用，裂隙冰Ω i应力场控制方程为
∂σij∂xi=0。 (11) 冻结过程中，裂隙冰与岩体相互作用，在裂隙冰Ω i与岩体Ω的交界面上，满足位移和应力连续条件。
1.4 控制方程的求解
冻结岩体裂隙冻胀扩展的相场法模拟中，待求解控制方程包括岩体Ω域内的应力场和相场控制方程，即式（7），以及裂隙冰Ω i域内的应力场控制方程，即式（11）。利用COMSOL软件实现控制方程的求解。
由于裂隙冻胀扩展是不可逆的，用来描述裂隙扩展的相场ϕ应是单调递增函数。因此，引入应变历史函数H(x,t)，并将其定义为[28]
H(x,t)=max (12) 将式（7）中 \psi _\varepsilon ^ + 替换为 H(\boldsymbol{x},t) ，可保证相场的演化是不可逆的。可利用COMSOL中Domain ODEs and DAEs模块求解式（12），以获得应变历史函数 H(\boldsymbol{x},t) 。
应力场控制方程的求解利用COMSOL中的固体力学（Solid Mechanics）模块，在此模块下建立两个线弹性材料子块，分别对应岩体 \mathit{\Omega} 域和裂隙冰 {\mathit{\Omega} _{{\text{ }}i}} 域内的应力场。在岩体 \mathit{\Omega} 域内，随着裂隙冻胀扩展，材料的弹性矩阵需不断更新，弹性矩阵D可表示为[20]
\begin{array}{r} \boldsymbol{D}=\frac{\partial \boldsymbol{\sigma}}{\partial \boldsymbol{\varepsilon}}=2 \mu\left\{\left[(1-k)(1-\phi)^2+k\right] P_{i j k l}^{+}+P_{i j k l}^{-}\right\}+ \\ \lambda\left\{\left[(1-k)(1-\phi)^2+k\right] H_{\varepsilon}(\operatorname{tr}(\varepsilon))+H_{\varepsilon}(-\operatorname{tr}(\varepsilon))\right\} J_{i j k l} 。 \end{array} (13) 式中： {J_{ijkl}} = {\delta _{ij}}{\delta _{kl}} ； {H_\varepsilon }(x) 为Heaviside函数； P_{ijkl}^ \pm 可表示为
\begin{aligned} P_{i j k l}^{ \pm}= & \sum\limits_{a=1}^3 \sum\limits_{b=1}^3 H_{\varepsilon}\left( \pm \varepsilon_a\right) \delta_{a b} n_{a i} n_{a j} n_{b k} n_{b l}+ \\ & \sum\limits_{a=1}^3 \sum\limits_{b \neq a}^3 \frac{1}{2} \frac{\left\langle\varepsilon_a\right\rangle_{ \pm}-\left\langle\varepsilon_b\right\rangle_{ \pm}}{\varepsilon_a-\varepsilon_b} n_{a i} n_{b j}\left(n_{a k} n_{b l}+n_{b k} n_{a l}\right), \end{aligned} (14) 其中， {n_{ai}} 为方向向量 \boldsymbol{n}_a 的第i个分量。
在裂隙冰 {\mathit{\Omega} _{{\text{ }}i}} 域内，将冰体视为各向同性线弹性材料，水冰相变膨胀线应变 {\varepsilon _{\text{f}}} 采用等效热膨胀系数的方式施加，将其等效为
{\varepsilon _{\text{f}}} = {\alpha _{\text{T}}}(T - {T_0}) \text{，} (15) 式中： {\alpha _{\text{T}}} 为等效热膨胀系数， T ， {T_0} 分别为当前温度与初始温度。将 T 设置为时间t的函数，即可获得冻结过程中随时间逐渐增大的相变膨胀线应变。
岩体 \mathit{\Omega} 域内相场控制方程利用Helmholtz方程模块，因为相场控制方程与Helmholtz方程具有相同的方程形式，可通过调整Helmholtz方程的系数，将其转化为相场控制方程，从而输入到COMSOL软件中。
计算过程中，首先固体力学模块中裂隙冰 {\mathit{\Omega} _{{\text{ }}i}} 域内体积膨胀对岩体 \mathit{\Omega} 域产生力的作用，求得 \mathit{\Omega} 域的主应变和应变能，更新应变历史函数 H(\boldsymbol{x},t) ；然后，相场 \phi 在与 H(\boldsymbol{x},t) 相应能量的作用下发生演化；最后，基于演化后的相场更新弹性矩阵D，进入下一个计算循环。计算时采用广义 \alpha 法进行非线性求解，并设置各个场之间采用分离式计算模式，以实现稳定求解。
2. 冻结岩体单裂隙冻胀扩展
2.1 不同倾角单裂隙冻胀扩展
刘艳章等[29]进行了单裂隙倾角分别为0°，30°，45°，60°，90°5种工况岩样的冻胀扩展试验，裂隙长2 cm，宽2 mm，高3.5 cm，岩样长12 cm，宽6 cm，高4 cm。基于此试验进行冻结岩体裂隙冻胀扩展的相场法模拟。因裂隙宽度相对高度较小，可视为平面应变问题。试算发现将试件长度方向两端各截去2 cm，建立8 cm×6 cm的计算模型，与建立12 cm×6 cm的模型，结果基本相同。为减少模型单元数，提高计算效率，采用8 cm×6 cm计算模型，如图 2，后续章节计算也进行类似简化。计算模型底边界为固定约束，其他边界无约束。岩石密度为2.17 g/cm3，弹性模量为6.1 GPa，泊松比为0.12，临界断裂能密度 {G_{\text{c}}} 为10 J/m2。裂隙冰的弹性模量为580 MPa，泊松比为0.3。 {l_0} 取0.25 mm，划分网格时，最大单元格尺寸取等于 {l_0} ，采用自由三角形单元（Free Triangular）进行网格划分。计算初始阶段，相场ϕ值较小，计算收敛性强，设置每个时步裂隙冰的膨胀线应变增量为57.78 με；计算至100时步时，相场ϕ值接近1，计算收敛性下降，将每个时步裂隙冰的膨胀线应变增量降低至28.89 με。
图 3为计算所得倾角0°裂隙岩体拉（压）应力分布、冻胀扩展过程与试验裂隙扩展照片对比。因文中相场云图图例相同，仅在图 3中提供图例，后续图中不再重复。图 3中，当裂隙冰冻胀体应变 \varepsilon _{\text{v}}^f 达到2.253%，在裂隙冰压力作用下，裂隙角部应力集中区域拉应力超过2 MPa，应力发展驱动相场演化，裂隙开始发生冻胀扩展，上侧两个角部获得扩展优势。随 \varepsilon _{\text{v}}^{\text{f}} 增至2.773%，裂隙角部拉应力集中区域向远处发展，冻胀扩展裂隙继续发展至约长10 mm。当 \varepsilon _{\text{v}}^{\text{f}} 达到3.293%，裂隙角部拉应力集中区域发展至接近边界，冻胀扩展裂隙接近贯通岩体，扩展裂隙与拉应力集中区域形态相似、位置对应重合。对比图 3，计算所得裂隙冻胀扩展形态与试验裂隙扩展照片接近，说明文中建立的相场法模型能够较准确地模拟冻结岩体裂隙冻胀扩展。
![]()
图 3 倾角0°应力分布与裂隙冻胀扩展Figure 3. Stresses and frost propagation of crack with dip angle of 0°图 4~7分别为裂隙倾角30°，45°，60°，90°时，裂隙冻胀扩展相场法模拟结果与试验照片对比。与倾角0°裂隙相比， \varepsilon _{\text{v}}^{\text{f}} 达到3.293%时，倾斜裂隙冻胀扩展长度明显变短，说明初始裂隙端部距离边界越远，裂隙所受约束作用越强，冻胀扩展越困难。倾角30°，45°裂隙工况下，岩体裂隙冻胀扩展规律相似，表现为裂隙沿距离边界较近的矩形初始裂隙对角处冻胀扩展，并逐渐向约束较弱的边界方向发展，扩展裂隙接近贯通岩体，与冻胀扩展试验结果相似。倾角60°，90°裂隙工况下，裂隙不再沿对角处发生冻胀扩展，而是沿矩形初始裂隙单侧角部发生冻胀扩展，与冻胀扩展试验结果相似。由图 4~7可见，不同倾角裂隙条件下，文中建立的相场法模型计算所得裂隙冻胀扩展规律与试验裂隙扩展形态照片接近。
2.2 不同外荷载作用下单裂隙冻胀扩展
Tang等[30]进行了无外荷载、竖向荷载3，6 MPa 3种工况下单裂隙岩体的冻胀扩展试验，裂隙长1 cm，宽1 mm，倾角为60°。基于此试验进行外荷载与冻结共同作用下裂隙冻胀扩展的相场法模拟。计算模型的几何尺寸和边界条件如图 8。当外荷载为竖向荷载时，计算模型底边界为固定约束，顶边界施加竖向荷载p，其他边界无约束；当外荷载为水平荷载时，计算模型底边界和一侧边界法向位移为零，另一侧边界施加水平荷载p，顶边界无约束。Tang等[30]采用的试验材料与刘艳章等[29]相同，因此，材料参数、 {l_0} 取值、网格划分、计算过程设置与2.1节相同。
图 9为无外荷载工况下裂隙冰冻胀体应变 \varepsilon _{\text{v}}^{\text{f}} =6.933%时，冻结岩体裂隙冻胀扩展计算结果与试验照片对比。计算结果与试验照片的裂隙冻胀扩展规律一致，在无外荷载作用下，冻胀扩展裂隙沿初始裂隙长度方向发展，发生共面裂隙扩展。
图 10为1 MPa竖向与水平荷载分别作用下 \varepsilon _{\text{v}}^{\text{f}} = 6.933%时冻结岩体裂隙冻胀扩展计算结果。在1 MPa竖向荷载作用下，裂隙不再发生共面扩展，由于竖向荷载作用一定程度地抑制了冻胀力引起的拉应力在竖直方向的分量，冻胀扩展向竖直方向偏转，扩展裂隙倾角为74°，初始裂隙倾角为60°，冻胀扩展向竖直方向偏转了14°。在1 MPa水平荷载作用下，水平荷载部分地抑制了冻胀力引起的拉应力在水平方向的分量，冻胀扩展向水平方向偏转，扩展裂隙倾角为40°，冻胀扩展向水平方向偏转了20°。尽管水平荷载与竖向荷载大小相等，由于水平荷载边界距离初始裂隙比竖向荷载边界近，水平荷载对裂隙冻胀扩展的约束作用更强，导致扩展裂隙向水平方向偏转的角度比向竖直方向大，且相同裂隙冰冻胀体应变 \varepsilon _{\text{v}}^{\text{f}} 时，水平荷载作用下裂隙冻胀扩展距离明显更短。
![]()
图 10 外荷载1 MPa工况裂隙冻胀扩展Figure 10. Frost propagation of crack with external load of 1 MPa图 11，12分别为竖向荷载3，6 MPa工况 \varepsilon _{\text{v}}^{\text{f}} = 6.933%时裂隙冻胀扩展计算结果与试验照片对比。图 11中，竖向荷载为3 MPa时，竖向荷载更大程度地抑制了冻胀力引起的拉应力在竖直方向的分量，裂隙冻胀扩展方向相对于1 MPa时进一步向竖直方向偏转，冻胀扩展裂隙倾角约为85°，计算得裂隙冻胀扩展形态与试验照片接近。由图 12可见，竖向荷载为6 MPa时，竖向荷载作用完全抑制了冻胀力引起的拉应力在竖直方向的分量，裂隙冻胀扩展笔直地沿竖直方向发展，计算结果与试验照片规律一致。
![]()
图 11 竖向外荷载3 MPa工况裂隙冻胀扩展Figure 11. Frost propagation of crack with vertical load of 3 MPa![]()
图 12 竖向外荷载6 MPa工况裂隙冻胀扩展Figure 12. Frost propagation of crack with vertical load of 6 MPa3. 冻结岩体双裂隙冻胀扩展
3.1 不同岩桥倾角的双裂隙冻胀扩展
李平等[13]进行了岩桥倾角α分别为45°，90°，135°，180°4种工况下双裂隙岩体的冻胀扩展试验，两条裂隙长均为1 cm，宽均为1 mm，倾角均为45°，岩桥中心点与岩体形心重合，岩桥长度为1 cm，如图 13。基于此试验进行冻结岩体双裂隙冻胀扩展的相场法模拟。计算模型尺寸与边界条件如图 13，底边界为固定约束，其他边界无约束。李平等[13]采用的试验材料与刘艳章等[29]相同，因此，材料参数、 {l_0} 取值、网格划分、计算过程设置与2.1节相同。
图 14，15分别为岩桥倾角45°，90°工况双裂隙冻胀扩展相场法模拟结果与试验照片对比。岩桥倾角45°时，岩桥与双裂隙共面，岩桥受两个裂隙冻胀扩展的共同作用，当裂隙冰冻胀体应变 \varepsilon _{\text{v}}^{\text{f}} 达到3.207%时，裂隙内尖端之间1 cm长的岩桥即将冻胀扩展贯通，而裂隙外尖端处的冻胀扩展距离仅约2.3 mm。当 \varepsilon _{\text{v}}^{\text{f}} 达到6.760%时，冻胀扩展沿双裂隙共面方向延伸几乎贯通岩体，扩展裂隙在接近边界处产生一定程度的向边界的弯曲，计算所得裂隙冻胀扩展形态与试验裂隙扩展情况相似。岩桥倾角90°时，岩桥与初始裂隙成135°夹角，裂隙外尖端处沿裂隙共面方向发生冻胀扩展，裂隙内尖端处冻胀扩展受两裂隙之间相互作用影响，扩展方向均向另一裂隙中心位置方向偏转，呈弧形发展。当 \varepsilon _{\text{v}}^{\text{f}} 达到4.507%时，裂隙内外尖端的冻胀扩展距离均约为6 mm；当 \varepsilon _{\text{v}}^{\text{f}} 达到7.973%时，裂隙外尖端的共面冻胀扩展接近贯通岩体，裂隙内尖端的弧形冻胀扩展裂隙接近贯通岩桥。计算所得裂隙冻胀扩展形态与试验照片相似。
图 16，17分别为岩桥倾角135°，180°时双裂隙冻胀扩展模拟结果与试验照片对比，模拟与试验裂隙冻胀扩展形态相似。岩桥倾角135°时，两条裂隙在同一竖直线上，裂隙外尖端仍沿裂隙共面方向冻胀扩展，裂隙内尖端冻胀扩展均向另一裂隙外尖端位置偏转，沿竖直方向发展。 \varepsilon _{\text{v}}^{\text{f}} 达到4.507%时，裂隙外尖端的冻胀扩展距离约为8mm，而内尖端的冻胀扩展长度仅约为4.5 mm； \varepsilon _{\text{v}}^{\text{f}} 达到7.973%时，裂隙外尖端的共面冻胀扩展接近边界位置，裂隙内尖端的冻胀扩展继续向着另一裂隙的外尖端沿竖直方向发展。另外，裂隙内尖端位置出现了裂隙冻胀分支扩展的现象。岩桥倾角180°时，裂隙外尖端仍沿裂隙共面方向冻胀扩展，裂隙内尖端出现了比135°工况更明显的冻胀分支扩展现象，其中一条分支扩展向另一初始裂隙外尖端位置偏转，另一条分支则向边界方向发展。
从上述4种岩桥倾角工况下双裂隙冻胀扩展规律可见，尽管裂隙倾角45°保持不变，不同岩桥倾角的双裂隙冻胀扩展差异显著。在无外荷载作用下，裂隙外尖端均沿裂隙共面方向发生冻胀扩展，基本不受两裂隙之间的相互作用影响；而裂隙内尖端的冻胀扩展受两裂隙间的相互作用影响显著，冻胀扩展方向总是向着另一条裂隙方向发展。岩桥倾角从90°~180°，随着裂隙中心距离逐渐变小，裂隙之间的相互作用变强，裂隙内尖端的冻胀扩展方向偏转程度也逐渐增大。
3.2 外荷载作用下的双裂隙冻胀扩展
在3.1节不同岩桥倾角双裂隙冻胀扩展基础上，在以图 13为例的计算模型的顶边界施加竖向荷载p，p取1，3 MPa两种工况，研究外荷载作用下的双裂隙冻胀扩展。参数取值、计算过程等与3.1节相同。
图 18为外荷载作用下岩桥倾角45°工况 \varepsilon _{\text{v}}^{\text{f}} = 6.760%时双裂隙冻胀扩展相场模拟结果。与图 14（b）对比可见，在竖向荷载作用下裂隙内尖端不再发生共面冻胀扩展，而是先向竖向扩展，然后再向另一条裂隙偏转，形成弧形冻胀扩展裂隙。1 MPa荷载时，双裂隙之间相互作用和外荷载对裂隙冻胀扩展影响均较显著，弧形扩展裂隙最终都弯向另一条初始裂隙内尖端，两条弧形扩展裂隙形成一个椭圆。3 MPa荷载时，外荷载对裂隙冻胀扩展的影响比双裂隙之间相互作用更显著，裂隙内尖端冻胀扩展向竖向偏转明显，仅稍向另一条裂隙弯曲成弧状。裂隙外尖端冻胀扩展基本不受双裂隙之间相互作用影响，在竖向荷载作用下，外尖端冻胀扩展均向竖直方向偏转，荷载越大偏转角度越大，这与荷载作用下单裂隙冻胀扩展规律一致。
图 19~21分别为外荷载作用下岩桥倾角90°，135°，180°工况 \varepsilon _{\text{v}}^{\text{f}} = 7.973%时双裂隙冻胀扩展相场法模拟结果。对比图 19，15（b）可见，此时竖向荷载与双裂隙之间的相互作用形成合效应，均引起冻胀扩展向竖直方向偏转，在竖向荷载作用下裂隙内尖端冻胀扩展进一步向竖直方向偏转，扩展裂隙弯曲程度变小。1 MPa荷载时，裂隙内尖端冻胀扩展裂隙间距比无外荷载时更近，岩桥接近贯通；3 MPa荷载时，裂隙内尖端沿竖直方向冻胀扩展，岩桥直接贯通。对比图 16（b），20可见，此时两条裂隙在同一竖直线上，竖向荷载与双裂隙之间的相互作用均引起裂隙内尖端沿竖直方向冻胀扩展，竖向荷载没有显著改变裂隙内尖端的冻胀扩展形态，但竖向荷载抑制了裂隙内尖端位置出现的裂隙冻胀分支扩展现象，尤其是3 MPa荷载时，冻胀分支扩展现象完全消失。对比图 17（b），21可见：①竖向荷载减弱了裂隙内尖端向另一条裂隙外尖端冻胀扩展的趋势，使裂隙内尖端冻胀扩展偏转角度变小；②竖向荷载也抑制了裂隙内尖端出现的裂隙冻胀分支扩展现象，1 MPa荷载时，裂隙内尖端冻胀分支扩展现象明显减弱，3 MPa荷载时，冻胀分支扩展现象基本消失。
3.3 不同夹角的双裂隙冻胀扩展
Tang等[30]进行了不同夹角双裂隙岩体冻胀扩展试验，两裂隙长均为1 cm，宽均为1 mm，左侧裂隙倾角90°，右侧裂隙保持形心位置不变，与左侧裂隙形成夹角β分别为0°，30°，45°3种工况，如图 22。双裂隙中心点距离为1 cm，中心点连线的中点与岩体形心重合。基于此试验进行冻结岩体双裂隙冻胀扩展的相场法模拟。计算模型尺寸与边界条件如图 22，计算模型底边界为固定约束，其他边界无约束。Tang等[30]采用的试验材料与刘艳章等[29]相同，因此，材料参数、 {l_0} 取值、网格划分、计算过程设置与2.1节相同。
图 23~25分别为裂隙夹角0°，30°，45°3种工况下双裂隙冻胀扩展相场法模拟结果与试验照片对比。整体上，计算所得裂隙冻胀扩展形态与试验裂隙扩展情况相似。夹角0°时，两裂隙相互平行，由于双裂隙之间区域受到两侧裂隙冻胀的挤压作用，冻胀扩展从裂隙尖端向背离另一条裂隙的方向偏转，呈“蝶状”扩展形态。
当夹角β增大至30°时，右侧倾斜裂隙靠近左侧竖直裂隙的尖端向左侧裂隙冻胀扩展，偏转角度达到32°，偏转的冻胀扩展裂隙与左侧裂隙的冻胀扩展裂隙交叉，形成一条主扩展裂隙，而远离左侧裂隙的尖端接近共面冻胀扩展，稍向边界偏转。β=45°时，右侧裂隙靠近左侧裂隙的尖端向左侧裂隙冻胀扩展，很快与左侧裂隙的扩展裂隙交叉，形成一条主扩展裂隙。对比3种工况的双裂隙冻胀扩展规律可见，平行双裂隙会以初始裂隙为桥梁形成两条独立的扩展裂隙，而具有倾斜交角的双裂隙冻胀扩展会交与一点，形成以扩展裂隙交点为中心的3条辐射状冻胀扩展裂隙。
4. 结论
（1）相场法模拟计算所得不同倾角、外荷载作用下冻结岩体单裂隙冻胀扩展，以及不同岩桥倾角、裂隙夹角工况时冻结岩体双裂隙冻胀扩展规律均与试验裂隙冻胀扩展形态接近，说明文中建立的相场法模型能够较准确地模拟冻结岩体的裂隙冻胀扩展。
（2）不同倾角的单裂隙表现为沿距离边界较近的初始裂隙角部冻胀扩展，并逐渐向约束较弱的边界方向发展。受外荷载作用时，单裂隙冻胀扩展向荷载作用方向偏转，荷载越大，偏转角度越大。
（3）对不同岩桥倾角的双裂隙，因裂隙之间的相互作用，裂隙内尖端冻胀扩展总是向另一条裂隙方向偏转，且裂隙距离越近，偏转角度越大，裂隙外尖端产生共面冻胀扩展。受外荷载作用时，双裂隙外尖端冻胀扩展向荷载方向偏转，内尖端冻胀扩展的偏转方向和程度受荷载大小和双裂隙间相互作用共同控制。
（4）对不同夹角双裂隙，平行双裂隙会以初始裂隙为桥梁形成两条独立呈“蝶状”的冻胀扩展裂隙，而具有倾斜交角的双裂隙冻胀扩展会交于一点，形成以扩展裂隙交点为中心的三条辐射状冻胀扩展裂隙。
-
![]()
图 3 倾角0°应力分布与裂隙冻胀扩展
Figure 3. Stresses and frost propagation of crack with dip angle of 0°
![]()
图 10 外荷载1 MPa工况裂隙冻胀扩展
Figure 10. Frost propagation of crack with external load of 1 MPa
![]()
图 11 竖向外荷载3 MPa工况裂隙冻胀扩展
Figure 11. Frost propagation of crack with vertical load of 3 MPa
![]()
图 12 竖向外荷载6 MPa工况裂隙冻胀扩展
Figure 12. Frost propagation of crack with vertical load of 6 MPa
-
[1] 刘泉声, 黄诗冰, 康永水, 等. 裂隙岩体冻融损伤研究进展与思考[J]. 岩石力学与工程学报, 2015, 34(3): 452-471. doi: 10.13722/j.cnki.jrme.2015.03.003 LIU Quansheng, HUANG Shibing, KANG Yongshui, et al. Advance and review on freezing-thawing damage of fractured rock[J]. Chinese Journal of Rock Mechanics and Engineering, 2015, 34(3): 452-471. (in Chinese) doi: 10.13722/j.cnki.jrme.2015.03.003
[2] DAVIDSON G, NYE J. A photoelastic study of ice pressure in rock cracks[J]. Cold Regions Science and Technology, 1985, 11(2): 141-153. doi: 10.1016/0165-232X(85)90013-8
[3] 黄诗冰, 刘泉声, 程爱平, 等. 低温岩体裂隙冻胀力与冻胀扩展试验初探[J]. 岩土力学, 2018, 39(1): 78-84. https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201801011.htm HUANG Shibing, LIU Quansheng, CHENG Aiping, et al. Preliminary experimental study of frost heaving pressure in crack and frost heaving propagation in rock mass under low temperature[J]. Rock and Soil Mechanics, 2018, 39(1): 78-84. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201801011.htm
[4] 乔趁, 王宇, 宋正阳, 等. 饱水裂隙花岗岩周期冻胀力演化特性试验研究[J]. 岩土力学, 2021, 42(8): 2141-2150. https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX202108010.htm QIAO Chen, WANG Yu, SONG Zhengyang, et al. Experimental study on the evolution characteristics of cyclic frost heaving pressure of saturated fractured granite[J]. Rock and Soil Mechanics, 2021, 42(8): 2141-2150. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX202108010.htm
[5] BOST M, POUYA A. Stress generated by the freeze–thaw process in open cracks of rock walls: empirical model for tight limestone[J]. Bulletin of Engineering Geology and the Environment, 2017, 76(4): 1491-1505. doi: 10.1007/s10064-016-0955-6
[6] LÜ Zhitao, XIA Caichu, WANG Yuesong, 等. Frost heave and freezing processes of saturated rock with an open crack under different freezing conditions[J]. 结构与土木工程前沿, 2020, 14(4): 947-960. https://www.cnki.com.cn/Article/CJFDTOTAL-PGJX202310015.htm LÜ Zhitao, XIA Caichu, WANG Yuesong, et al. Frost heave and freezing processes of saturated rock with an open crack under different freezing conditions[J]. Frontiers of Structural and Civil Engineering, 2020, 14(4): 947-960. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-PGJX202310015.htm
[7] JIA H, LEITH K, KRAUTBLATTER M. Path-dependent frost-wedging experiments in fractured, low-permeability granite[J]. Permafrost and Periglacial Processes, 2017, 28(4): 698-709. doi: 10.1002/ppp.1950
[8] TAN X, CHEN W, LIU H, et al. A unified model for frost heave pressure in the rock with a penny-shaped fracture during freezing[J]. Cold Regions Science and Technology, 2018, 153: 1-9. doi: 10.1016/j.coldregions.2018.04.016
[9] 刘泉声, 黄诗冰, 康永水, 等. 低温冻结岩体单裂隙冻胀力与数值计算研究[J]. 岩土工程学报, 2015, 37(9): 1572-1580. doi: 10.11779/CJGE201509003 LIU Quansheng, HUANG Shibing, KANG Yongshui, et al. Numerical and theoretical studies on frost heaving pressure in a single fracture of frozen rock mass under low temperature[J]. Chinese Journal of Geotechnical Engineering, 2015, 37(9): 1572-1580. (in Chinese) doi: 10.11779/CJGE201509003
[10] 刘泉声, 黄诗冰, 康永水, 等. 裂隙冻胀压力及对岩体造成的劣化机制初步研究[J]. 岩土力学, 2016, 37(6): 1530-1542. doi: 10.16285/j.rsm.2016.06.002 LIU Quansheng, HUANG Shibing, KANG Yongshui, et al. Preliminary study of frost heave pressure and its influence on crack and deterioration mechanisms of rock mass[J]. Rock and Soil Mechanics, 2016, 37(6): 1530-1542. (in Chinese) doi: 10.16285/j.rsm.2016.06.002
[11] HUANG S, LIU Y, GUO Y, et al. Strength and failure characteristics of rock-like material containing single crack under freeze-thaw and uniaxial compression[J]. Cold Regions Science and Technology, 2019, 162: 1-10. doi: 10.1016/j.coldregions.2019.03.013
[12] WANG Y, FENG W, WANG H, et al. Rock bridge fracturing characteristics in granite induced by freeze-thaw and uniaxial deformation revealed by AE monitoring and post-test CT scanning[J]. Cold Regions Science and Technology, 2020, 177: 103115. doi: 10.1016/j.coldregions.2020.103115
[13] 李平, 唐旭海, 刘泉声, 等. 双裂隙类砂岩冻胀断裂特征与强度损失研究[J]. 岩石力学与工程学报, 2020, 39(1): 115-125. https://www.cnki.com.cn/Article/CJFDTOTAL-YSLX202001012.htm LI Ping, TANG Xuhai, LIU Quansheng, et al. Experimental study on fracture characteristics and strength loss of intermittent fractured quasi-sandstone under freezing and thawing[J]. Chinese Journal of Rock Mechanics and Engineering, 2020, 39(1): 115-125. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YSLX202001012.htm
[14] THARP T. Conditions for crack propagation by frost wedging[J]. Bulletin of the Geological Society of America, 1987, 99(1): 94-102. doi: 10.1130/0016-7606(1987)99<94:CFCPBF>2.0.CO;2
[15] HUANG S, LIU Q, LIU Y, et al. Frost heaving and frost cracking of elliptical cavities (fractures) in low-permeability rock[J]. Engineering Geology, 2018, 234: 1-10. doi: 10.1016/j.enggeo.2017.12.024
[16] 黄诗冰, 刘泉声, 刘艳章, 等. 低温热力耦合下岩体椭圆孔(裂)隙中冻胀力与冻胀开裂特征研究[J]. 岩土工程学报, 2018, 40(3): 459-467. doi: 10.11779/CJGE201803009 HUANG Shibing, LIU Quansheng, LIU Yanzhang, et al. Frost heaving pressure and characteristics of frost cracking in elliptical cavity (crack) of rock mass under coupled thermal-mechanical condition at low temperature[J]. Chinese Journal of Geotechnical Engineering, 2018, 40(3): 459-467. (in Chinese) doi: 10.11779/CJGE201803009
[17] ZHOU Y. Numerical simulation of fracture propagation in freezing rocks using the extended finite element method (XFEM)[J]. International Journal of Rock Mechanics and Mining Sciences, 2021, 148: 104963. doi: 10.1016/j.ijrmms.2021.104963
[18] TAO S, TANG X, RUTQVIST J, et al. The influence of stress anisotropy and stress shadow on frost cracking in rock[J]. Computers and Geotechnics, 2021, 133: 103967. doi: 10.1016/j.compgeo.2020.103967
[19] MIEHE C, HOFACKER M, WELSCHINGER F. A phase field model for rate-independent crack propagation: Robust algorithmic implementation based on operator splits[J]. Computer Methods in Applied Mechanics and Engineering, 2010, 199(45-48): 2765-2778.
[20] ZHOU S, ZHUANG X, ZHU H, et al. Phase field modelling of crack propagation, branching and coalescence in rocks[J]. Theoretical and Applied Fracture Mechanics, 2018, 96: 174-192.
[21] ZHOU S, ZHUANG X, RABCZUK T. A phase-field modeling approach of fracture propagation in poroelastic media[J]. Engineering Geology, 2018, 240: 189-203.
[22] 刘嘉, 薛熠, 高峰, 等. 层理页岩水力裂缝扩展规律的相场法研究[J]. 岩土工程学报, 2022, 44(3): 464-473. doi: 10.11779/CJGE202203008 LIU Jia, XUE Yi, GAO Feng, et al. Propagation of hydraulic fractures in bedded shale based on phase-field method[J]. Chinese Journal of Geotechnical Engineering, 2022, 44(3): 464-473. (in Chinese) doi: 10.11779/CJGE202203008
[23] ZHOU S, ZHUANG X, RABCZUK T. Phase field method for quasi-static hydro-fracture in porous media under stress boundary condition considering the effect of initial stress field[J]. Theoretical and Applied Fracture Mechanics, 2020, 107: 102523.
[24] SWEIDAN A, HEIDER Y, MARKERT B. A unified water/ice kinematics approach for phase-field thermo-hydro- mechanical modeling of frost action in porous media[J]. Computer Methods in Applied Mechanics and Engineering, 2020, 372: 113358.
[25] BORDEN M, VERHOOSEL C, SCOTT M, et al. A phase-field description of dynamic brittle fracture[J]. Computer Methods in Applied Mechanics and Engineering, 2012, 217: 77-95.
[26] 丁遂栋. 断裂力学[M]. 北京: 机械工业出版社, 1997. DING Suidong. Fracture Mechanics[M]. Beijing: China Machine Press, 1997. (in Chinese)
[27] 吴家龙. 弹性力学[M]. 2版. 北京: 高等教育出版社, 2011. WU Jialong. Elasticity[M]. 2nd ed. Beijing: Higher Education Press, 2011. (in Chinese)
[28] MIEHE C, WELSCHINGER F, HOFACKER M. Thermodynamically consistent phase-field models of fracture: Variational principles and multi-field FE implementations[J]. International Journal for Numerical Methods in Engineering, 2010, 83(10): 1273-1311.
[29] 刘艳章, 郭赟林, 黄诗冰, 等. 冻融作用下裂隙类砂岩断裂特征与强度损失研究[J]. 岩土力学, 2018, 39(增刊2): 62-71. https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX2018S2010.htm LIU Yanzhang, GUO Yunlin, HUANG Shibing, et al. Study of fracture characteristics and strength loss of crack quasi-sandstone under freeze-thaw cycles[J]. Rock and Soil Mechanics, 2018, 39(S2): 62-71. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX2018S2010.htm
[30] TANG X, TAO S, LI P, et al. The propagation and interaction of cracks under freeze-thaw cycling in rock-like material[J]. International Journal of Rock Mechanics and Mining Sciences, 2022, 154: 105112.
-
期刊类型引用(4)
1. 刘乃飞，李宁，汪双杰，宋战平，王莉平，徐拴海. 含相变低温岩体水–热–力特性研究进展与展望. 岩石力学与工程学报. 2025(03): 521-542 . 
百度学术
2. 乔雄，杨小龙，冯勇. 高海拔严寒地区隧道冻害研究现状及展望. 隧道建设(中英文). 2024(02): 233-256 . 百度学术
3. 李聪然，马明星，潘懿，宋军涛，杜晓飞，刘文进，岳祖润. 季冻区露天矿山边坡稳定性分析及加固措施研究. 河北地质大学学报. 2024(03): 73-79 . 百度学术
4. CHEN Wenhua，XIANG Tian. Frost heaving and crack initiation characteristics of tunnel rock mass in cold regions under low-temperature degradation. Journal of Mountain Science. 2024(08): 2844-2854 . 必应学术
其他类型引用(2)
-
其他相关附件
下载:
计量
- 文章访问数: 343
- HTML全文浏览量: 54
- PDF下载量: 113
- 被引次数: 6
