首页/文章/ 详情

即将直播：岩石的PB模型研究和双轴标定思路（7月14日）

详细信息

作者优秀: 优秀教师/意见领袖/博士学历/特邀专家/独家讲师
平台推荐: 内容稀缺

仿真圈

2年前浏览3003

作者 | lobby 仿真秀优秀讲师

首发 | 仿真秀APP

导读：参数标定是离散元中联系宏观和微观的桥梁，我们模拟的材料是表现出岩石的特征，抑或是金属的特征，就是取决于参数标定的过程。本文基于PFC的平行胶结模型，讲解岩石参数的标定过程。为什么讲解PB模型，而不是从简单的线性模型开始。原因在于，虽然线性模型的参数比较少，但是随之而来的孔隙率、级配等非接触模型参数的影响会放大，导致多因素影响下的标定困难。为了帮助大家更好理解，7月14日20时，笔者将在2022仿真知识周（仿真人自己的节日）第四期讲座带来《岩石的PB模型研究和双轴标定思路》公开直播，详情见后文。

一、构件器

在做参数标定前，我们必须理解我们所使用的模型，只有理解模型的作用机理以及微观参数意义，才能根据每一次结果的反馈去合理的调整参数。

PB模型是一个比较标准的岩石模型，其不光继承了传统线性模型的部分，而且基于微观的物质成分出发，增加了胶结模型部分。为什么其叫平行胶结模型，就是因为这个模型在构件层次有两个法向和切向弹簧。说到构件层次，就得顺便提一下PFC手册中，每个模型都会提供的构件图。

图1.1 构件图

图1.1 展示了PFC中前三个模型的构件图，线性模型、cb模型和pb模型。这里分别介绍一下图中所使用的构件器。

（1）黏壶：用于提供接触间的粘滞性，就好比在两个颗粒间加了棉花一样，给其相对运动提供一定的缓冲，但是不参与静力学意义上的平衡，也就是说速度为0或者几乎为0的情况下其是不起作用的。

（2）

弹簧：用于提**生位移后的力，或者一定力下的位移。这个概念比较简单，我相信各位小时候没有没玩过弹簧的。

（3）分离器：用于仅产生压力，而不产生拉力。上下的横杠杠就像是两个手掌，手掌贴紧的时候，上下手臂才会慢慢产生力，分开后手臂就没有力了。

（4）摩擦器：用于提供切向的摩擦，阻止无限的相对移动。这里的“无限”两个字需要理解一下。因为横向弹簧的存在，所以这个摩擦器并不能阻碍切向的相对运动，但是当摩擦器在一定条件下失效后，这个相对运动就不会和切向弹簧相关了，后面会跟静力学中动摩擦的概念类似。但是比较明显的是，摩擦器生效的时候，这个和静力学中静摩擦的概念还是有区别的。（希望讲这个部分知识的时候，各位的高中物理知识还没还给老师。）

（5）

胶结器：用于提供胶结力，比较明显的是这个力除了会是压力，也有可能是拉力。

在cb模型中，是有一个弹簧，相对于线性模型，胶结器取代分离器，和弹簧是串联关系。这个构件器特征意味着这个模型适用于单一物质的材料。

相对于cb模型，pb模型提供了一个额外的弹簧，并且胶结器和这个额外的弹簧串联。这种特征也是代表这个模型可以代表两个材料。一个材料是被胶结物，一个材料是胶结物。

在这里也回答一下广受争议的一个问题。

我们该如何选择接触模型？

首先肯定是应当从构件器出发的，cb模型胶结器和弹簧是在一起的，也就代表被胶结物和胶结物是一种材料。从此来看，cb模型更加适用于连续性材料。而pb模型区分了被胶结物和胶结物，所以其更加适用于胶结性材料。比如本文所介绍的岩石类材料。

二、变形力学特性

PB模型的变形是线性的，和线性的模型的线性可能还是有一点区别。线性模型中的线性体现在颗粒层次，但是颗粒的相互嵌入，导致在宏观上，线性模型产生的变形和所受到的力在大部分时候并不是线性的。

而PB模型由于胶结的存在，至少在胶结破坏前，颗粒不会产生过大的错动，更不用说嵌入了。这种特性导致PB模型，在宏观上的力位移曲线是线性的。这个其实也是对我们的参数标定带来了便捷。

图2.1 pb模型的法向力位移关系

图2.1为PB模型的胶结法向力位移曲线，其中gs为颗粒间的重叠量，当这个值为正时，颗粒是分开的，颗粒接触上是拉力。这个值为负值时，接触上是压力。图上可以看出，胶结的变形是线性的，当正向拉力大于某个值的时候，其会发生破坏，导致胶结消失。而压力方向上，这个线是无穷的，代表着胶结并不会因为压力而发生破坏。

图2.2 pb模型的切向力位移关系

图2.2为PB模型的胶结切向力位移曲线，其实ds为切向位移。有一些同学会问：这切向位移怎么只有正的，没有负的。关于正负的问题是很多同学去纠结的，可能也不会是在这个地方。

这里也要借鲁迅先生的嘴说句话：正负本身并没有意义，人给他定义多了后，才有了意义。

对于切向来说，正向和反向对于接触来说是一样的，所以正负是没有必要去定义，这里直接给个正值代表切向位移就够了。
对于切向来说，其变形也是线性的，当切向力达到某个值后会发生胶结的破坏。
对于二维来说，弯转也会有破坏，三维的时候，扭转也会有破坏。这里不去过多的介绍，机理和切向类似。
需要介绍的是，正向的力是反应在正应力上，切向力弯转力都反应在切向力上。

PB模型的破坏机理是正向拉应力和切向应力达到某个值是会发生破坏的。

三、强度力学特性

下面就介绍一下这个值是怎么确定的。图3.1(a)为pb模型胶结的强度包线，x向为正应力，正向为压应力，y向为切应力。这个曲线，学土木的学生应该会觉得眼熟，我们学习的第一个关于破坏的本构模型就是摩尔库伦模型。图3.1(b)为我们学习的摩尔库伦模型，可以发现的是，这两个在切应力的角度都是一条斜直线，所以这两个的公式也是非常相近的。

不过相对于摩尔库伦模型，PB胶结的强度包线多了一个抗拉强度，反应到图上就是左边的竖直线。做成这样的原因是，斜直线的强度包线是一个比较理想的情况。很多实验表明，岩石等胶结材料的强度包线是一个类似于子弹头形状的，这里我在做直播的时候看看能不能找到那个图，在很多弹性力学的教材中是有介绍到的。

相对于完全的斜直线，用抗拉强度封闭左边的强度，更加吻合那种子弹头的强度包线，并且相对于曲线，直线判断起来会更加容易一点。

图3.1 pb模型的强度包线

四、微观参数

介绍完了PB模型基本的力学特性，这里简单展示一下上面力学特性所对应的微观参数了。

（1）颗粒部分

变形相关：kn\ks 或者emod\kratio 力位移曲线的斜率
强度相关：fric 决定切向力的峰值，相当于静力学的动摩擦系数

（2）胶结部分

变形相关：pb_kn\pb_ks 或者pb_emod\pb_kratio 力位移曲线的斜率
强度相关：pb_ten 胶结抗拉强度 pb_coh 切应力曲线的截距 pb_fa 切应力曲线的斜角度

五、参数敏感性分析

这里以双轴为参数标定使用的单元实验，先进行参数的敏感性分析。这里主要针对于上述使用的参数。

（1）Emod

相对于kn，我更喜欢用emod来决定颗粒的变形特性，原因在之前的课程中介绍过，直播中也会再次阐述，欢迎大家来battle。

图5.1为不同emod对应的应力应变关系，可以发现的是emod对曲线斜率的影响是比较显著的，图5.2也给出了20MPa对应的弹性模量的关系，可以看到几乎是线性相关的。图5.3反映了emod不光对变形有影响，也会引起强度的下降。

图5.1 不同emod对应的应力应变曲线

图5.2 emod和弹性模量E的关系

图5.3 emod和强度的关系

（2）PB_Emod

图5.4为不同pb_emod对应的应力应变曲线，可以看到其效果和emod类似，对刚度的影响也是正相关的，且几乎为直线。而其对强度的影响却几乎相反，越大的pb_emod式样容易产生更大的强度。

图5.4 不同pb_emod对应的应力应变曲线

图5.5 pb_emod和弹性模量E的关系

图5.6 pb_emod和强度的关系

（3）PB_Coh

图5.7为不同pb_coh对应的应力应变曲线，可以看到其对曲线斜率几乎没有影响，对强度的影响是显著的，几乎是线性的。

图5.7 不同pb_coh对应的应力应变曲线

图5.8 pb_coh和强度关系

（4）PB_Fa

可以看到pb_fa 对强度的影响也是有的，相对于pb_coh，其也会增加强度包线的斜率。

图5.9 pb_fa和强度关系

总结：

(1) 强度参数对变形无影响；
(2) 变形参数会影响刚度，而且也会对强度有一定的影响；
(3) Pb_coh会影响强度的大小，对包线的斜率无影响；
(4) Pb_fa会影响强度的大小，也会影响包线的斜率。

六、岩石参数标定思路

（1）标定思路

我们的参数标定思路也是从上述的结论出发的。

1、Emod\pb_emod根据得到的弹性模量进行调整
2、Pb_coh根据得到的强度均值进行调整
3、Pb_fa根据得到的包线斜率进行调整

这里假设三个围压 w1、w2、w3对应的目标轴压为z1、z2、z3。目标弹性模量为E。提取一下强度均值：QP=(z1 z2 z3)/3 强度包线比率：R=z3/z1

每一次调整的流程为：

1、获得三个围压下的轴压x1、x2、x3。及其对应的应变e1、e2、e3。

2、变形参数调整。

计算当前三个弹性模量：E1=(x1-w1)/e1 E2=(x2-w2)/e2 E3=(x3-w3)/e3
平均弹性模量为：EP=(E1 E2 E3)/3
EmodNext=EmodLast*(EP)/E
PB_EmodNext=PB_EmodLast*(EP)/E
这里的EmodLast和PB_EmodLast是上一次使用的弹性参数，EmodNext和PB_EmodNext是下一次使用的。

3、强度参数调整

计算三个围压的强度均值：QPN=(x1 x2 x3)/3 强度包线比率：RN=x3/x1
调整pb_coh：PB_CohNext=PB_CohLast*QP/QPN
调整pb_fa： PB_FaNext=PB_FaLast*R/RN

反复进行上述过程，直到所模拟的曲线满足要求。

（2）案例

目标弹性模量5GPa、三个围压5MPa、10MPa、20MPa对应的轴压为52.74MPa、69.13MPa、96.62MPa。

默认参数为：

结果反馈

图6.1 标定过程应力应变曲线

可以看到在第二次的时候，弹性模量就有了比较大的修正。后面更多的是对强度的修正。模糊的来说，第5次的结果就足够作为最终的标定结果了。按此流程标定的次数越多，得到的结果越精确，第10次已经有比较好的结果了，标定到第19次的时候，与目标结果基本是吻合的。感兴趣的同学可以在根据曲线计算一下弹性模量，图中可以直接看出强度峰值，和目标值对比一下。

七、我的公开课

参数标定是很费时间的，这一点毋庸置疑。特别是满足的应力路径和破坏模式越多，对参数的要求也会越来越高。我也有师兄花两三年的时候标定一套满足各种应力路径的参数，这个是需要时间和经验积累的。

本文主要针对于双轴实验，基于微观参数特性出发，基本给出了满足双轴实验的标定流程。这个流程重复的次数越多，得到的最后结果也是越精确的，相对于没有经验的凑参数，在本文的框架下至少可以在足够的时间付出后，可以得到满足条件的参数的。

为了帮助大家更好理解，7月14日20时，笔者将在2022仿真知识周（仿真人自己的节日）第四期讲座带来《岩石的PB模型研究和双轴标定思路》公开直播，以下是课程安排

八、关于仿真知识周

迄今为止，仿真秀已连续举办了三届仿真知识周-仿真人自己的节日，今年7月11日，第四届仿真知识周将如期而至。抱着回馈广大仿真学习者与从业者的热情，仿真秀请来30位行业专家为大家奉上知识盛宴，平台的专栏作者们也纷纷慷慨让利，释出大量优惠，让这个专属仿真人的节日比炎炎夏日更加火热！

首先，平台面向所有用户赠出每人200元全场通用优惠券，可与平台的限时秒杀同时叠加使用，多重优惠叠加，点击“阅读原文”或者扫码进入活动主页立刻领券。在本届仿真知识周中，我们汇集了30位专家，推出20场线上直播，其中有知名软件原厂与行业大咖的技术沙龙，也有全网稀缺学习资源的讲师在线分享经验。专题涵盖：软件原厂技术沙龙、汽车行业、电磁仿真、土木建筑、制造工艺、工业设备。以下是直播安排：