CN111444641A - 一种考虑冻融环境下的岩体工程稳定性分析方法 - Google Patents
一种考虑冻融环境下的岩体工程稳定性分析方法 Download PDFInfo
- Publication number
- CN111444641A CN111444641A CN202010084632.6A CN202010084632A CN111444641A CN 111444641 A CN111444641 A CN 111444641A CN 202010084632 A CN202010084632 A CN 202010084632A CN 111444641 A CN111444641 A CN 111444641A
- Authority
- CN
- China
- Prior art keywords
- damage
- freeze
- stress
- thaw
- rock mass
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 239000011435 rock Substances 0.000 title claims abstract description 103
- 238000004458 analytical method Methods 0.000 title claims abstract description 29
- 238000010257 thawing Methods 0.000 title claims description 31
- 230000006378 damage Effects 0.000 claims abstract description 122
- 238000000034 method Methods 0.000 claims abstract description 44
- 230000000694 effects Effects 0.000 claims abstract description 14
- 230000008878 coupling Effects 0.000 claims abstract description 13
- 238000010168 coupling process Methods 0.000 claims abstract description 13
- 238000005859 coupling reaction Methods 0.000 claims abstract description 13
- 230000014509 gene expression Effects 0.000 claims description 32
- 238000012937 correction Methods 0.000 claims description 25
- 230000008569 process Effects 0.000 claims description 24
- 238000011068 loading method Methods 0.000 claims description 22
- 238000004364 calculation method Methods 0.000 claims description 17
- 230000005484 gravity Effects 0.000 claims description 13
- 239000000463 material Substances 0.000 claims description 13
- 230000011514 reflex Effects 0.000 claims description 11
- IAZDPXIOMUYVGZ-UHFFFAOYSA-N Dimethylsulphoxide Chemical compound CS(C)=O IAZDPXIOMUYVGZ-UHFFFAOYSA-N 0.000 claims description 3
- 230000002706 hydrostatic effect Effects 0.000 claims description 3
- 238000005516 engineering process Methods 0.000 abstract 1
- 238000007710 freezing Methods 0.000 description 18
- 230000008014 freezing Effects 0.000 description 17
- 238000010586 diagram Methods 0.000 description 9
- 230000001133 acceleration Effects 0.000 description 5
- 230000015556 catabolic process Effects 0.000 description 5
- 238000006731 degradation reaction Methods 0.000 description 5
- 239000011159 matrix material Substances 0.000 description 5
- 238000004422 calculation algorithm Methods 0.000 description 4
- 230000008859 change Effects 0.000 description 4
- 238000011160 research Methods 0.000 description 4
- 230000003313 weakening effect Effects 0.000 description 4
- 230000009471 action Effects 0.000 description 3
- 230000003902 lesion Effects 0.000 description 3
- 238000012360 testing method Methods 0.000 description 3
- 238000009825 accumulation Methods 0.000 description 2
- 238000005056 compaction Methods 0.000 description 2
- 238000012669 compression test Methods 0.000 description 2
- 230000001808 coupling effect Effects 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 230000007774 longterm Effects 0.000 description 2
- 238000013507 mapping Methods 0.000 description 2
- 230000009467 reduction Effects 0.000 description 2
- 239000002689 soil Substances 0.000 description 2
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 2
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 1
- 208000027418 Wounds and injury Diseases 0.000 description 1
- 230000002411 adverse Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000000556 factor analysis Methods 0.000 description 1
- 208000014674 injury Diseases 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 229920006395 saturated elastomer Polymers 0.000 description 1
- 230000001932 seasonal effect Effects 0.000 description 1
- 238000010187 selection method Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000000638 solvent extraction Methods 0.000 description 1
Images
Landscapes
- Investigating Strength Of Materials By Application Of Mechanical Stress (AREA)
Abstract
本发明实施例公开了一种考虑冻融环境影响下的岩体工程稳定性分析方法，其包括：S1、创建表达冻融次数与岩体损伤之间关系的冻融损伤模型；S2、创建表达岩体塑性引起的损伤效应的岩体荷载损伤模型；S3、创建表达冻融荷载综合作用的冻融荷载耦合损伤演化模型；S4、基于Hoek‑Brown岩体非线性强度准则，创建弹塑性损伤模型；S5、求解所述弹塑性损伤模型对应的应力；S6、获取达到极限平衡状态时岩体工程对应的安全系数。本发明解决了传统技术中冻融岩体稳定性分析存在误差，准确率低的问题。
Description
技术领域
本发明实施例涉及数字处理分析领域，尤其涉及一种考虑冻融环境影响下 的岩体工程稳定性分析方法。
背景技术
季冻区岩体工程稳定性受到冻融和荷载的综合影响。一方面，岩体中的水 冰相变导致岩体冻胀力，导致岩体微裂纹扩展；另一方面，外部荷载作用改变 了岩体的应力，也会引起岩体形变和破坏的产生。因此，研究冻融和荷载综合 作用下岩体损伤规律和计算方法，对保证岩体工程安全具有重要的意义。但是 目前冻融损伤模型的研究仅限于线性强度准则，不能很好地适应岩体节理、强 度非线性的特点。
具体的：已有的研究有的采用摩尔库仑体作为塑性准则，建立硬岩冻融损 伤的长期变形模型构建方法。还有的研究建立了一种等围岩下的岩石冻融损伤 本构，也是基于DP的塑性准则。但是摩尔库仑体和DP准则都是线性强度准则， 不太适合节理岩体的非线性强度特征；
同时由于传统的岩体稳定性分析方法大都采用理想弹塑性模型和线性强度 准则，因此带来很大的误差。如，已有的硬岩冻融损伤的长期变形模型和等围 岩下的岩石冻融损伤本构目前大多数是一维的本构，其并不能用于三维岩体的 数值计算。现有冻融环境下H-B参数选取方法不能反映冻融-荷载耦合作用；因 此可以说稳定性分析是岩土工程中一个重要难题，而考虑冻融影响的三维岩体 安全系数计算更是一个有待解决的技术问题。
发明内容
基于此，为解决现有技术所存在的不足，特提出了一种考虑冻融环境影响 下的岩体工程稳定性分析方法。
一种考虑冻融环境影响下的岩体工程稳定性分析方法，其特征在于，包括：
S1、创建表达冻融次数与岩体损伤之间关系的冻融损伤模型；
S2、创建表达岩体塑性引起的损伤效应的岩体荷载损伤模型；
S3、创建表达冻融荷载综合作用的冻融荷载耦合损伤演化模型；
S4、基于Hoek-Brown岩体非线性强度准则，创建弹塑性损伤模型；
S5、求解所述弹塑性损伤模型对应的应力；
S6、获取达到极限平衡状态时岩体工程对应的安全系数。
可选的，在其中一个实施例中，所述冻融损伤模型对应的表达式：
其中，Dt表示冻融损伤变量，N表示冻融循环次数，Nf表示实验最大冻融次 数，Df表示与Nf对应的损伤值；q、r均为常数；
其中，n0和n分别表示岩石冻融前后孔隙率，VP和V′P分别表示岩石冻融前 后纵波波速。
可选的，在其中一个实施例中，所述岩体荷载损伤模型对应的表达式为：
其中，
表示等效塑性应变，Dm表示受荷损伤变量，α与β均为经验系数；
等效塑性应变
对应的表达式为：
式中：εp1、εp2、εp3分别为x轴、y轴、z轴三个方向上的主塑性应变；
可选的，在其中一个实施例中，所述冻融荷载耦合损伤演化模型对应的表达 式为：
Dc＝Dt+Dm-DtDm (5)
其中，Dc表示岩石的冻融-荷载总损伤变量；DtDm表示耦合项。
可选的，在其中一个实施例中，所述弹塑性损伤模型对应的表达式为：
式中：p表示静水压力，J2表示第二偏应力不变量，θ表示罗德角，σci表示完 整岩石单轴抗压强度，mb，s和a均表示经验常数；
其中，
式(7)-(9)中：GSI表示地质强度指标；mi表示岩石的软硬程度；
对应的剪切模量G(Dc)、体积模量K(Dc)表达式分别为
式中，G0表示初始剪切模量，K0表示初始体积模量，E0表示未损 伤材料的弹性模量，ν表示泊松比。
可选的，在其中一个实施例中，求解所述弹塑性损伤模型对应的应 力的过程包括：
S51、在主应力空间中，对所述弹塑性损伤模型进行弹性预测分析以获得所 述模型自弹性状态进入塑性状态对应的极限参数；
S52、进行塑性修正分析，即对应力空间进行划分，判断应力回映位置，并 根据回映位置，建立更新应力及多个或单一塑性因子对应的求解式以对应力进 行更新；所述应力回映位置包括单一屈服面、双屈服面相交的棱线或者多屈服 面相交的尖点处中的一种或者多种组合；
S53、进行塑性修正分析，即完成应力更新后，计算出等效塑性应变后对损 伤变量进行更新并获得所对应的应力张量；
S54、获取一致切线模量以给定迭代求解过程所需的模量。
可选的，在其中一个实施例中，所述进行塑性修正分析的过程包括：
S521、对应力空间进行划分即在主应力空间中给定多个屈服函数
f1＝σ1-σ3-σci[(1-Dc)(mbσ3/σci+s)]a
f2＝σ2-σ3-σci[(1-Dc)(mbσ3/σci+s)]a
f3＝σ2-σ1-σci[(1-Dc)(mbσ1/σci+s)]a
f4＝σ3-σ1-σci[(1-Dc)(mbσ1/σci+s)]a
f5＝σ3-σ2-σci[(1-Dc)(mbσ2/σci+s)]a
f6＝σ1-σ2-σci[(1-Dc)(mbσ2/σci+s)]a
其中，塑性势函数gi与对应的屈服函数形式相同，由参数
sg、ag分别 替换mb、s、a即可；则多屈服面本构模型的流动法则表达式为：
式中：dεp为塑性应变增量；dλi为塑性因子增量；m为大于零的屈服函数的 个数，dλi≥0，fi≤0，dλifi(σ，γ，ω)≤0；
S522、塑性修正过程中，保持应变增量Δεn+1及损伤变量(Dc)n不变，则在塑性 状态下，对应的应变增量中包含弹性应变增量Δεe和塑性应变增量Δεp，其中， 塑性应变增量表达式为：
式中，m为所述屈服函数的函数值大于零的方程个数，应力修正量Δσp的表 达式为：
式中：
为屈服面fi＞0上的塑性应力修正方向；
S523、进行应力返回位置的判断并求得更新应力值σn+1：则更新后的应力和 内变量表达式为：
可选的，在其中一个实施例中，tn+1时刻的所对应的应力张量σ′n+1对应的计 算公式为：
可选的，在其中一个实施例中，通过重力加载安全系数分析方法获取达到 极限平衡状态时岩体工程对应的安全系数，重力加载法下岩体安全系数K的表 达式为：
式中：Sstep为计算所述模型破坏时最大加载步数；Δg离心加载系数，Δg＝0.01；γ为所述模型的重度。
实施本发明实施例，将具有如下有益效果：
本发明提出了一种考虑强度非线性准则的、考虑强度退化和刚度退化的冻 融岩体工程稳定性的分析方法，其能够有效提高冻融环境的岩体工程安全性分 析的准确率。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施 例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述 中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付 出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
其中：
图1为一个实施例中实施温冻融区域岩体稳定性分析方法技术流程图；
图2为一个实施例中损伤引起的刚度退化示意图；
图3a为一个实施例中荷载作用下损伤演化规律示意图(β对损伤演化规律 的影响(α＝0.01))；
图3b为一个实施例中荷载作用下损伤演化规律示意图(α对损伤演化规律 的影响(β＝0.3))；
图4为一个实施例中H-B准则在主应力空间的区域划分示意图；
图5为一个实施例中本文的弹塑性损伤模型应力回映算法求解过程示意图；
图6为一个实施例中某片麻岩损伤变量随冻融次数变化规律示意图；
图7为一个实施例中不同冻融次数下岩石应力-应变曲线示意图；
图8为一个实施例中边坡计算模型示意图；
图9a为一个实施例中不同加载步下边坡损伤演化规律(冻融次数100)示意 图(K＝1.30)；
图9b为一个实施例中不同加载步下边坡损伤演化规律(冻融次数100)示意 图(K＝1.60)。
图9c为一个实施例中不同加载步下边坡损伤演化规律(冻融次数100)示意 图(K＝1.70)。
图9d为一个实施例中不同加载步下边坡损伤演化规律(冻融次数100)示意 图(K＝1.94)。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实 施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅 仅用以解释本发明，并不用于限定本发明。
除非另有定义，本文所使用的所有的技术和科学术语与属于本发明的技术 领域的技术人员通常理解的含义相同。本文中在本发明的说明书中所使用的术 语只是为了描述具体的实施例的目的，不是旨在限制本发明。可以理解，本发 明所使用的术语“第一”、“第二”等可在本文中用于描述各种元件，但这些元 件不受这些术语限制。这些术语仅用于将第一个元件与另一个元件区分。举例 来说，在不脱离本申请的范围的情况下，可以将第一元件称为第二元件，且类 似地，可将第二元件为第一元件。第一元件和第二元件两者都是元件，但其不 是同一元件。
在本实施例中，特提出了一种考虑冻融环境影响下的岩体工程稳定性分析 方法，如图1所示，该方法包括：S1、创建表达冻融次数与岩体损伤之间关系 的冻融损伤模型；S2、创建表达岩体塑性引起的损伤效应的岩体荷载损伤模型； S3、创建表达冻融荷载综合作用的冻融荷载耦合损伤演化模型；S4、基于 Hoek-Brown岩体非线性强度准则，创建弹塑性损伤模型；S5、求解所述弹塑性 损伤模型对应的应力；S6、获取达到极限平衡状态时岩体工程对应的安全系数。 其中，Hoek-Brown岩体非线性强度准则简称H-B准则，是非线性强度准则，其 适合节理岩体的强度特征，可以应用到冻融损伤分析中，即可以基于所述原理 建立冻融损伤本构模型并求解，但是其在求解过程中在应力空间中还存在奇异 点的等问题，为此，本发明首先考虑冻融-荷载耦合作用对岩体强度和刚度的弱 化作用，建立基于H-B准则的岩体弹塑性耦合损伤模型并推导了其本构方程； 同时为解决模型数值求解过程中的奇异点问题，给出了其在主应力空间中的完 全隐式求解算法。另在此基础上采用离心加载法，以在计算过程中保持c、
不 变，逐步增加重力加速度G的值直到边坡发生失稳破坏，将此时的重力加速度 值与实际重力加速度值之比定义为安全系数的方法实现分析过程。
在一些具体的实施例中，所述创建表达冻融次数与岩体损伤之间关系的冻 融损伤模型；根据冻融次数与损伤的关系即以纵波波速和孔隙率双物理参数为 变量表示动弹性模量，并在此基础上，所述冻融损伤模型对应的表达式：
其中，Dt表示冻融损伤变量，N表示冻融循环次数，Nf表示实验最大冻融次 数，Df表示与Nf对应的损伤值；q、r均为材料常数，通过实验获取即可；
其中，n0和n分别表示岩石冻融前后孔隙率，VP和V′P分别表示岩石冻融前 后纵波波速。
在一些具体的实施例中，所述S2在考虑荷载作用下，创建表达岩体塑性引 起的损伤效应的岩体荷载损伤模型；在外力荷载的作用下，岩石内部同样会产 生损伤，试验研究表明：损伤随着等效塑性应变
的累积而加剧，受荷损伤变 量Dm可以表示为
的幂指数形式，即所述岩体荷载损伤模型对应的表达式为：
其中，
表示等效塑性应变，Dm表示受荷损伤变量，α与β均为经验系数，通 过实验获取即可；通常α取值范围为[0，+∞]，决定了损伤后岩石材料软化曲线的 初始斜率；β取值范围为[0，1]，决定了岩石最大损伤值，不同α，β值下的损伤 变量变化曲线如图3所示。从图3中可以看出，α越大，损伤演化速率越慢；β 越大，岩石的最终损伤值越小。
其中，等效塑性应变
对应的表达式为：
式中：εp1、εp2、εp3分别为x轴、y轴、z轴三个方向上的主塑性应变。
在一些具体的实施例中，所述S3创建表达冻融荷载综合作用的冻融荷载耦合 损伤演化模型；所述冻融荷载耦合损伤演化模型对应的表达式为：
Dc＝Dt+Dm-DtDm (5)
式中：Dc为岩石的冻融-荷载总损伤变量；DtDm为耦合项。由式(6)可知，耦 合损伤使冻融-荷载总损伤有所弱化；冻融-荷载耦合损伤变量表示的岩石应力- 应变关系：σ＝(1-Dc)Eε。
在一些具体的实施例中，所述S4基于Hoek-Brown岩体非线性强度准则， 创建弹塑性损伤模型；
所述弹塑性损伤模型对应的表达式为：
式中：p表示静水压力，J2表示第二偏应力不变量，θ表示罗德角，σci表示完 整岩石单轴抗压强度，mb，s和a均表示经验常数，均用于反映岩体质量和岩体 类型；
其中，
式(7)-(9)中：GSI表示地质强度指标；mi表示岩石的软硬程度；DC为 冻融-荷载共同作用的损伤变量，反映了冻融与荷载综合作用对强度的影响；
对应的剪切模量G(Dc)、体积模量K(Dc)表达式可以用初始剪切模量G0和 初始体积模量K0表示，即对应分别为
式中，G0表示初始剪切模量，K0表示初始体积模量，E0表示未损伤材料的 弹性模量，ν表示泊松比。同时其体现了损伤对材料刚度的影响，当应变达到一 定的阈值时，弹性模量值随着损伤的累积开始减小，材料应力-应变关系显现出 非线性特征，此为损伤引起刚度退化，见图2。
在一些具体的实施例中，所述S5求解所述弹塑性损伤模型对应的应力；由 于H-B准则在应力空间中存在棱线和尖点等奇异点，在这些奇异点处屈服函数 外法线方向不唯一，导数不连续，使得数值实施存在一定的困难；此外在主应 力空间中考虑问题能够减少应力维数(从六维减少到三维)，在几何上也更加直 观。因此本文推导了主应力空间中耦合损伤模型的完全隐式求解算法，求解过 程分为弹性预测，塑性修正和损伤修正。在塑性修正过程中，首先对应力空间 进行划分，判断应力回映点的位置(单一屈服面、双屈服面相交的棱线或者多 屈服面相交的尖点处)，根据回映位置的不同，建立更新应力及多个或单一塑性 因子的Newton-Raphson求解式，很好地解决了空间角点问题。
基于上述原理可知，弹性预测与塑性修正均在有效应力空间进行，最后通 过损伤修正得到名义应力张量，则求解所述弹塑性损伤模型对应的应力的过程 包括：
S51、在主应力空间中，对所述弹塑性损伤模型进行弹性预测分析以获得所 述模型自弹性状态进入塑性状态对应的极限参数；即对所述弹塑性损伤模型进 行弹性预测分析过程：已知tn时刻应变增量Δεn+1、应力σn和损伤变量(Dc)n，则tn+1时刻的弹性预测应力为：
式中：
为损伤弹性矩阵；计算
并代入f1中对应力状态进行判断：即若f1＜0，材料仍处于弹性阶段，对变量进行更新：
若f1＞0，则进入塑性修正阶段。
S52、进行塑性修正分析，即对应力空间进行划分，判断应力回映位置，并 根据回映位置，建立更新应力及多个或单一塑性因子对应的求解式以对应力进 行更新；所述应力回映位置包括单一屈服面、双屈服面相交的棱线或者多屈服 面相交的尖点处中的一种或者多种组合；所述进行塑性修正分析的过程包括： S521、基于3个主应力之间的大小关系，对应力空间进行划分即在主应力空间 中给定多个屈服函数，对应的公式为
f1＝σ1-σ3-σci[(1-Dc)(mbσ3/σci+s)]a
f2＝σ2-σ3-σci[(1-Dc)(mbσ3/σci+s)]a
f3＝σ2-σ1-σci[(1-Dc)(mbσ1/σci+s)]a
f4＝σ3-σ1-σci[(1-Dc)(mbσ1/σci+s)]a
f5＝σ3-σ2-σci[(1-Dc)(mbσ2/σci+s)]a
f6＝σ1-σ2-σci[(1-Dc)(mbσ2/σci+s)]a
其中，塑性势函数gi与对应的屈服函数形式相同，由参数
sg、ag分别 替换mb、s、a即可；则多屈服面本构模型的流动法则表达式为：
式中：dεp为塑性应变增量；dλi为塑性因子增量；m为式mb公式中大于零 的屈服函数的方程个数，由Kuhn-Tucker加卸载准则可得 dλi≥0，fi≤0，dλifi(σ，γ，ω)≤0；
S522、塑性修正过程中，保持应变增量Δεn+1及损伤变量(Dc)n不变，则在塑性 状态下，对应的应变增量中包含弹性应变增量Δεe和塑性应变增量Δεp，其中， 塑性应变增量表达式为：
式中，m为所述屈服函数的函数值大于零的方程个数，此时的更新应力位 于屈服面以外，需要对其进行修正，应力修正量Δσp的表达式为：
式中：
为屈服面fi＞0上的塑性应力修正方向；
S523、进行应力返回位置的判断并求得更新应力值σn+1：则更新后的应力和 内变量表达式为：
进行应力返回位置的判断过程包括：在f1与f2的交线l1上有σ1＝σ2(为方便 描述，对下标n进行省略)，因此直线方程为：
同理，f1和f6的交线l6上有σ2＝σ3，直线方程为：
当更新应力只位于屈服面f1上时，由于塑性修正应力增量
的方 向h1：
式中：v为泊松比，k的表达式如下：
如，图4对应力空间进行了区域划分。由图4可以看出，h1与f1的两条边 界线组成的边界面确定了应返回至f1面的应力区域；h1、h2和f1、f2的交线确定 了返回至l1线的区域，h1、h6和f1、f2的交线确定了返回值l6线的区域(h2，h6表 达式见下文)。
在l1上取一点σl1＝(σ1，σ1，σ3)，l6上取一点σl6＝(σ1，σ3，σ3)，若同时满足：
则更新应力位于f1面上。式中：
rl1、rl6为棱线l1和l6的方向向量。
若
则更新应力位于l1线上，同理若
则更新应 力位于l6线上。由图3可以看到，由h1、h2组成的平面将尖点位置与l1分割开来， 该平面的法向量为：
na1＝h1×h2
由h1、h6组成的平面将尖点位置与l6区域分割开来，平面法向量为：
na6＝h6×h1
尖点处应力值为σapex＝(σa，σa，σa)，σa＝sσci/mb。若有：
和
同时成立，则更新应力位于尖点处。
应力求解过程的关键是对下述方程组进行求解，即：
式中：R为残差值。
S53、进行塑性修正分析，即完成应力更新后，通过式(2)计算出等效塑性应 变后代入式(1)对损伤变量进行更新，即
并将上式 代入式所述演化公式求得当前迭代步下的总损伤变量Dc，继而对mb，s进行修 正，获得所对应的应力张量；tn+1时刻的所对应的应力张量
对应的计算公式 为：
S54、获取一致切线模量以给定迭代求解过程所需的模量。为了保证有限元 方程组整体迭代求解过程中具有二阶收敛速度，需要给出一致切线模量的表达 式：
式中：n为屈服函数值大于零的屈服函数个数；A为n阶方阵：
δij为Kronecker符号，αi的表达式为：
式中：I为单位矩阵；至此完成对上述弹塑性损伤模型应力回映求解过程，具体 见图5。
在一些具体的实施例中，所述S6获取达到极限平衡状态时岩体工程对应的 安全系数，即通过重力加载的方法对岩体工程进行计算，直到达到极限平衡状 态，重力加载系数就是该岩体工程对应得安全系数。
首先建立岩体工程的模型，基于S2-S4的模型和算法，通过对逐步增加重力 加速度的过程基于进行数值模拟，直至岩体工程达到临界破坏的状态，此时加 重力加速度增加的倍数即为岩体工程的安全系数；其中，重力加载法下岩体安 全系数K的表达式为：
式中：Sstep为计算模型破坏时最大加载步数；Δg离心加载系数，本文选取Δg＝0.01； γ为模型重度。
下面以实际案例对所述方法进行进一步说明以及验证；为了验证所建冻融- 荷载耦合损伤模型，进行室内冻融三轴压缩试验：岩样取自吉林省某边坡工程， 主要组成为混合片麻岩。首先，根据工程所在地区气候情况，确定冻融循环温 度为-18℃-22℃，冻融循环周期定为12h(-18℃保持6h，22℃保持6h)，分别 对试样进行0、10、20、40和80次冻融。冻融前将试样放在饱和器中饱和6h， 然后取出试样放入盛水容器，再将其放入冻融箱中进行循环，冻融过程中保持 试样的饱和状态；其次，分别对冻融前后的岩石试样进行声波波速和孔隙率测 量，将所得结果代入式(10)，计算冻融损伤变量Dt，并按照式(9)进行拟合，所 得结果如图6所示。并对未冻融的试样进行不同围压下的压缩试验，获取力学 参数如下：重度γ＝20kN/m3，弹性模量E＝15GPa，泊松比v＝0.23；确定H-B强 度参数GSI＝90，mi＝16.4，σci＝130MPa，s＝0.03后，通过对峰后曲线拟合的方 式得到损伤参数a＝0.02，β＝0.8，(Dc)ini＝Dt。再次分别对比无冻融，冻融次数为 10、40和80次下的岩石应力应变曲线的试验值与计算值，结果如图7所示。 由图7可知，在非线性压密段二者存在一定的偏差，原因在于本文模型没有考 虑材料的压密特性，但却很好的体现了其弹性模量随冻融次数的增大而发生弱 化的现象。经历0、10、40和80次冻融的岩石试样，由试验所得峰值强度分别 为123.58MPa、107.79MPa、94.71MPa和81.49MPa，数值计算所得峰值强度为 125.40MPa、112.75MPa、98.96MPa和83.90MPa，与试验值吻合较好。当试样 处于弹性阶段时，耦合损伤变量Dc＝Dt，随着冻融次数的增加，Dc值不断增大， 通过式(2)、(3)和(5)对岩石强度参数mb、s和刚度参数E起到弱化作用，使得试 样峰值强度和弹性模量不断降低；进入峰后段时，试样开始产生塑性应变并不 断累积，由式(8)可知荷载损伤逐步发育，导致试样强度不断降低，最终发生破 坏(计算不收敛)。综上，本发明所建模型能较好的反映冻融-荷载耦合损伤对 岩体强度和刚度的弱化作用。
以下再用一个例子来说明本发明方法的应用效果。吉林省某边坡工程，四季 变化明显，冬季昼夜温差较大。为研究冻融次数对工程稳定性的影响规律，本 文做出如下假定：入口段围岩计算参数与前节试验室相同；从最不利的方面考 虑问题，假定整个研究模型都处于冻融范围内。根据设计图纸建立边坡计算模 型如图8所示。模型共划分为11511个节点和11275个单元，左右两侧施加水平 边界约束，底部施加固定边界约束。
由于所建模型参数较多，利用强度折减法计算安全系数较为困难，因此本 文选用离心加载法对其进行求解。求解过程中，保持其余各参数不变，逐步增 加模型的容重，直至计算出现不收敛现象即判定边坡发生破坏。计算结果如表1 所示。
表1边坡安全系数
由表1可知，离心加载法计算所得安全系数与传统Spencer法求得的安全系 数相比略大，这是因为边坡破坏形式主要为剪切破坏，而岩土材料的抗剪强度 与法线应力成正比例关系，不断增加岩体容重的同时，其抗剪强度也在增长， 因此得到的安全系数较大。考虑荷载损伤时，边坡安全系数大幅降低，这是由 于在容重加载过程中，局部材料发生软化现象，强度降低，使得计算过早出现 不收敛。上述计算结果一方面证明了所述有限元分析方法的正确性，同时也说 明了荷载损伤对岩土工程的稳定性具有较大的影响。
另本例记录冻融循环次数为100(Dt＝0.172)时，离心加载过程中，边坡损伤场 演化规律如图9a-d所示。由图9各图可知，由于冻融循环的劣化作用，边坡存 在初始损伤场；随着重力荷载的增加，坡体内部产生荷载损伤，与冻融损伤共 同作用，造成岩体材料刚度和强度的退化，最终导致坡体发生破坏。
以上所述实施例仅表达了本申请的几种实施方式，其描述较为具体和详细， 但并不能因此而理解为对本申请专利范围的限制。应当指出的是，对于本领域 的普通技术人员来说，在不脱离本申请构思的前提下，还可以做出若干变形和 改进，这些都属于本申请的保护范围。因此，本申请专利的保护范围应以所附 权利要求为准。
Claims (9)
1.一种考虑冻融环境影响下的岩体工程稳定性分析方法，其特征在于，包括：
S1、创建表达冻融次数与岩体损伤之间关系的冻融损伤模型；
S2、创建表达岩体塑性引起的损伤效应的岩体荷载损伤模型；
S3、创建表达冻融荷载综合作用的冻融荷载耦合损伤演化模型；
S4、基于Hoek-Brown岩体非线性强度准则，创建弹塑性损伤模型；
S5、求解所述弹塑性损伤模型对应的应力；
S6、获取达到极限平衡状态时岩体工程对应的安全系数。
2.根据权利要求1所述的方法，其特征在于，所述冻融损伤模型对应的表达式：
其中，Dt表示冻融损伤变量，N表示冻融循环次数，Nf表示实验最大冻融次数，Df表示与Nf对应的损伤值；q、r均为常数；
其中，n0和n分别表示岩石冻融前后孔隙率，VP和V′P分别表示岩石冻融前后纵波波速。
3.根据权利要求1所述的方法，其特征在于，所述岩体荷载损伤模型对应的表达式为：
其中，
表示等效塑性应变，Dm表示受荷损伤变量，α与β均为经验系数；
等效塑性应变
对应的表达式为：
式中：εp1、εp2、εp3分别为x轴、y轴、z轴三个方向上的主塑性应变。
4.根据权利要求1所述的方法，其特征在于，所述冻融荷载耦合损伤演化模型对应的表达式为：
Dc＝Dt+Dm-DtDm (5)
其中，Dc表示岩石的冻融-荷载总损伤变量；DtDm表示耦合项。
5.根据权利要求1所述的方法，其特征在于，所述弹塑性损伤模型对应的表达式为：
式中：p表示静水压力，J2表示第二偏应力不变量，θ表示罗德角，σci表示完整岩石单轴抗压强度，mb，s和a均表示经验常数；
其中，
式(7)-(9)中：GSI表示地质强度指标；mi表示岩石的软硬程度；
对应的剪切模量G(Dc)、体积模量K(Dc)表达式分别为
式中，G0表示初始剪切模量，K0表示初始体积模量，E0表示未损伤材料的弹性模量，ν表示泊松比。
6.根据权利要求1所述的方法，其特征在于，求解所述弹塑性损伤模型对应的应力的过程包括：
S51、在主应力空间中，对所述弹塑性损伤模型进行弹性预测分析以获得所述模型自弹性状态进入塑性状态对应的极限参数；
S52、进行塑性修正分析，即对应力空间进行划分，判断应力回映位置，并根据回映位置，建立更新应力及多个或单一塑性因子对应的求解式以对应力进行更新；所述应力回映位置包括单一屈服面、双屈服面相交的棱线或者多屈服面相交的尖点处中的一种或者多种组合；
S53、进行塑性修正分析，即完成应力更新后，计算出等效塑性应变后对损伤变量进行更新并获得所对应的应力张量；
S54、获取一致切线模量以给定迭代求解过程所需的模量。
7.根据权利要求6所述的方法，其特征在于，所述进行塑性修正分析的过程包括：
S521、对应力空间进行划分即在主应力空间中给定多个屈服函数
f1＝σ1-σ3-σci[(1-Dc)(mbσ3/σci+s)]a
f2＝σ2-σ3-σci[(1-Dc)(mbσ3/σci+s)]a
f3＝σ2-σ1-σci[(1-Dc)(mbσ1/σci+s)]a
f4＝σ3-σ1-σci[(1-Dc)(mbσ1/σci+s)]a
f5＝σ3-σ2-σci[(1-Dc)(mbσ2/σci+s)]a
f6＝σ1-σ2-σci[(1-Dc)(mbσ2/σci+s)]a
其中，塑性势函数gi与对应的屈服函数形式相同，由参数mbg、sg、ag分别替换mb、s、a即可；则多屈服面本构模型的流动法则表达式为：
式中：dεp为塑性应变增量；dλi为塑性因子增量；m为大于零的屈服函数的个数，dλi≥0，fi≤0，dλifi(σ，γ，ω)≤0；
S522、塑性修正过程中，保持应变增量Δεn+1及损伤变量(Dc)n不变，则在塑性状态下，对应的应变增量中包含弹性应变增量Δεe和塑性应变增量Δεp，其中，塑性应变增量表达式为：
式中，m为所述屈服函数的函数值大于零的方程个数，应力修正量Δσp的表达式为：
式中：
为屈服面fi＞0上的塑性应力修正方向；
S523、进行应力返回位置的判断并求得更新应力值σn+1：则更新后的应力和内变量表达式为：
8.根据权利要求6所述的方法，其特征在于，tn+1时刻的所对应的应力张量σ′n+1对应的计算公式为：
9.根据权利要求1所述的方法，其特征在于，通过重力加载安全系数分析方法获取达到极限平衡状态时岩体工程对应的安全系数，重力加载法下岩体安全系数K的表达式为：
式中：Sstep为计算所述模型破坏时最大加载步数；Δg离心加载系数，Δg＝0.01；γ为所述模型的重度。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010084632.6A CN111444641A (zh) | 2020-02-10 | 2020-02-10 | 一种考虑冻融环境下的岩体工程稳定性分析方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010084632.6A CN111444641A (zh) | 2020-02-10 | 2020-02-10 | 一种考虑冻融环境下的岩体工程稳定性分析方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN111444641A true CN111444641A (zh) | 2020-07-24 |
Family
ID=71652477
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010084632.6A Pending CN111444641A (zh) | 2020-02-10 | 2020-02-10 | 一种考虑冻融环境下的岩体工程稳定性分析方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN111444641A (zh) |
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112945750A (zh) * | 2020-12-15 | 2021-06-11 | 大连海事大学 | 一种低温岩石蠕变试验装置及系统 |
| CN113032955A (zh) * | 2021-02-05 | 2021-06-25 | 中国科学院武汉岩土力学研究所 | 一种适用于地震荷载下岩石动态本构模型的构建方法 |
| CN113505514A (zh) * | 2021-08-04 | 2021-10-15 | 大连海事大学 | 一种复杂扰动条件下岩体弹塑性损伤-渗流耦合计算方法 |
| CN115235884A (zh) * | 2022-06-22 | 2022-10-25 | 重庆交通大学 | 一种饱和砂岩的损伤状况判别及剩余强度预估方法 |
| CN115389363A (zh) * | 2022-07-28 | 2022-11-25 | 四川路桥桥梁工程有限责任公司 | 一种测量混凝土盐冻后强度的预测模型 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107515291A (zh) * | 2017-08-23 | 2017-12-26 | 西安科技大学 | 一种等围压作用下岩石冻融损伤本构模型的构建方法 |
| CN108829916A (zh) * | 2018-04-25 | 2018-11-16 | 中铁二院工程集团有限责任公司 | 硬岩冻融损伤长期变形模型的构建方法 |
| CN109522611A (zh) * | 2018-10-25 | 2019-03-26 | 长江大学 | 新型岩体损伤本构模型构建方法及装置 |
| CN110705165A (zh) * | 2019-10-08 | 2020-01-17 | 中国石油大学(华东) | 一种构建岩石材料弹塑性-损伤耦合力学本构模型的方法 |
-
2020
- 2020-02-10 CN CN202010084632.6A patent/CN111444641A/zh active Pending
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107515291A (zh) * | 2017-08-23 | 2017-12-26 | 西安科技大学 | 一种等围压作用下岩石冻融损伤本构模型的构建方法 |
| CN108829916A (zh) * | 2018-04-25 | 2018-11-16 | 中铁二院工程集团有限责任公司 | 硬岩冻融损伤长期变形模型的构建方法 |
| CN109522611A (zh) * | 2018-10-25 | 2019-03-26 | 长江大学 | 新型岩体损伤本构模型构建方法及装置 |
| CN110705165A (zh) * | 2019-10-08 | 2020-01-17 | 中国石油大学(华东) | 一种构建岩石材料弹塑性-损伤耦合力学本构模型的方法 |
Non-Patent Citations (4)
| Title |
|---|
| QUANSHENG LIU 等: "a prediction model for uniaxial compressive strength of deteriorated rocks due to freeze-thaw", 《COLD REGIONS SCIENCE AND TECHNOLOGY》 * |
| 唐春安 等: "岩土破裂过程分析RFPA 离心加载法", 《岩土工程学报》 * |
| 张慧梅 等: "冻融与荷载耦合作用下岩石损伤模型的研究", 《岩石力学与工程学报》 * |
| 许梦飞 等: "基于Hoek-Brown准则的岩体弹塑性损伤模型及其应力回映算法研究", 《工程力学》 * |
Cited By (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112945750A (zh) * | 2020-12-15 | 2021-06-11 | 大连海事大学 | 一种低温岩石蠕变试验装置及系统 |
| CN113032955A (zh) * | 2021-02-05 | 2021-06-25 | 中国科学院武汉岩土力学研究所 | 一种适用于地震荷载下岩石动态本构模型的构建方法 |
| CN113505514A (zh) * | 2021-08-04 | 2021-10-15 | 大连海事大学 | 一种复杂扰动条件下岩体弹塑性损伤-渗流耦合计算方法 |
| CN113505514B (zh) * | 2021-08-04 | 2024-01-05 | 大连海事大学 | 一种复杂扰动条件下岩体弹塑性损伤-渗流耦合计算方法 |
| CN115235884A (zh) * | 2022-06-22 | 2022-10-25 | 重庆交通大学 | 一种饱和砂岩的损伤状况判别及剩余强度预估方法 |
| CN115389363A (zh) * | 2022-07-28 | 2022-11-25 | 四川路桥桥梁工程有限责任公司 | 一种测量混凝土盐冻后强度的预测模型 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN111444641A (zh) | 一种考虑冻融环境下的岩体工程稳定性分析方法 | |
| CN111695285B (zh) | 一种各向异性岩体应力-损伤-渗流耦合数值模拟方法 | |
| JP5911077B2 (ja) | 空気と水と土骨格の連成計算装置および連成計算方法並びに連成計算プログラム | |
| Valliappan et al. | Seismic analysis of arch dams—a continuum damage mechanics approach | |
| Dujc et al. | Quadrilateral finite element with embedded strong discontinuity for failure analysis of solids | |
| Orduña | Non-linear static analysis of rigid block models for structural assessment of ancient masonry constructions | |
| Réthoré et al. | A discrete model for the dynamic propagation of shear bands in a fluid‐saturated medium | |
| Liao et al. | Crack propagation modelling using the weak form quadrature element method with minimal remeshing | |
| Desai et al. | Numerical algorithms and mesh dependence in the disturbed state concept | |
| Lilja et al. | Finite-discrete element modelling of sea ice sheet fracture | |
| CN113505514B (zh) | 一种复杂扰动条件下岩体弹塑性损伤-渗流耦合计算方法 | |
| Luo et al. | Viscoelastic analysis of the creep characteristics of interlayered rock specimens under uniaxial compression | |
| Lee et al. | Structural analysis of pipe-framed greenhouses using interface elements for cross-over connections | |
| Rafi et al. | Analytical modeling of concrete beams reinforced with carbon FRP bars | |
| CN114662341A (zh) | 岩体临界滑动面极限分析方法 | |
| Chen et al. | Accelerated and stabilized meshfree method for impact-blast modeling | |
| Yao et al. | Material point method for numerical simulation of failure phenomena in multiphase porous media | |
| Al-Shayea et al. | A plastic-damage model for stress-strain behavior of soils | |
| Li et al. | Explicit integration algorithm of the bounding surface model based on swell–shrink rules for cyclic behaviors of clay | |
| Liu et al. | A structured grid based B-Spline finite elements method combining local isogeometry analysis technique for crack problems | |
| Lundqvist et al. | A plane-stress plasticity model for masonry for the explicit finite element time integration scheme | |
| Wackerfuß | Efficient finite element formulation for the analysis of localized failure in beam structures | |
| CN114595495B (zh) | 一种沥青混凝土心墙-土石料接触面数值模拟方法 | |
| Cattaneo et al. | A comparison between numerical integration algorithms for unsaturated soils constitutive models | |
| Peng et al. | FEM simulation of deformation and strength characteristics in geogrid-reinforced sand retaining wall under the change of loading rate |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| RJ01 | Rejection of invention patent application after publication |
Application publication date: 20200724 |
|
| RJ01 | Rejection of invention patent application after publication |