在PFC里跑单轴压缩这件事,入门很简单,跑出“岩石味”却要命地难。我早期用平行黏结模型做出来的试件,裂纹齐刷刷沿着加载方向一路劈到底,和文献里那种矿物尺度上的沿晶破裂、穿晶破裂混合形态完全对不上。后来我才意识到,问题不在于接触本构参数没调好,而在于模型里根本没有矿物晶粒和晶界这个结构层级。这也是为什么我后来把重心全转到GBM(Grain-Based Model,基于晶粒的模型)上。
这个算例的完整目标,是在PFC5.0里搭一个由三种矿物晶体组成的类岩石材料试件,跑一遍单轴压缩2d全流程,得到应力应变曲线、裂纹演化过程和最终破坏形态。三种矿物怎么布、晶界怎么设、参数怎么标定、加载时哪些坑一定会踩,这篇文章把整个从零搭建到出结果的思路完整拆开来讲。对于刚从“会跑PB”上升到“想搞懂GBM”的朋友,这应该能帮你省下不少瞎试的时间。
1. 平行黏结模型跑不出“岩石味”,GBM到底改了什么
1.1 为什么均匀颗粒模型让人越跑越怀疑
很多初学PFC的人,第一次单轴压缩都是拿平行黏结模型(Parallel Bond Model)来做的。操作确实简单:生成一堆圆球,给接触赋上平行黏结参数,然后压到破坏。但跑过几次之后你会发现,这种模型有几个很难受的特征。
第一,破坏路径太规整。颗粒是随机排列的,但裂纹一旦形成,往往直接连成一条宏观裂缝,沿着加载方向拉通。这种形态在真实岩石试件里很少见,除非是特别致密均质的细粒岩石。第二,宏观力学响应缺乏“过程感”。真实岩石在峰前有明显的非线性变形、微破裂萌生、峰后软化,但简单的平行黏结模型经常在峰值前一直线性,到了峰值瞬间崩掉,应力应变曲线形态和实验对不上。第三,也是最根本的问题:模型里没有矿物矿物学结构,你无法回答“裂纹是从石英内部穿过还是沿着长石晶界扩展”这种问题。而岩石破坏,恰恰是从这些细观结构边界开始的。
我当时反复调平行黏结的抗拉强度和内聚力,把峰值强度调得和实验一致了,结果破裂模式还是一根劈裂裂缝;去对比文献里的声发射定位图,完全没有那种“早期离散微破裂逐步连成宏观破裂面”的演化规律。这才意识到这是模型结构的天花板,不是参数问题。
1.2 GBM把岩石当成“多晶体复合体”
GBM的核心思想,是把岩石试件还原成一个由若干矿物晶粒拼接而成的多晶体。每一个晶粒内部是一种材料属性,晶粒与晶粒之间是强度偏弱的晶界。在受压过程中,应力在矿物之间的传递和重新分配,破裂同时存在两种模式:裂纹穿过晶粒内部,叫穿晶破裂;裂纹沿着晶界走,叫沿晶破裂。两种模式的比例,直接决定了宏观破坏是脆性劈裂、剪切带还是更复杂的复合破坏。
相比平行黏结模型,GBM等于给模型增加了两个结构层级:一是矿物晶粒这个“块体”单元,二是晶界这个“弱面”系统。晶粒内部用强度较高的接触模型模拟矿物晶体,晶界用强度更低的接触模拟,这样应力会优先在晶界处集中和释放。宏观上看到的现象就是:荷载还很低时,晶界处已经出现微破裂;随着荷载增加,微破裂逐渐聚集成宏观裂纹;最终破坏面往往是绕过强硬矿物颗粒、穿过软弱矿物,沿着晶界走出一条不规则的路径。这条路径看着“扭曲”,恰恰就是真实岩石最典型的破坏形态。
对于只熟悉平行黏结模型的人,你可以把GBM理解成“给颗粒堆加了配筋”:每个“配筋区域”代表一个矿物晶体,区域之间的接缝是天然的裂缝引导线。这样一来,微观裂纹的起裂位置、扩展方向、贯通路径都变得有据可依,而不是靠颗粒随机排列碰运气。
1.3 三种矿物组合对应什么岩石材料
标题里特别强调“三种矿物组成”,这不是随便凑数。天然岩石绝大多数是多矿物的,其中最常见的组合有三种:
- 石英 + 钾长石 + 黑云母:典型的花岗岩类组合,也是GBM论文里最常见的对象。
- 石英 + 长石 + 白云母/黏土矿物:更接近砂岩、片麻岩类。
- 方解石 + 白云石 + 黏土矿物:碳酸盐岩类。
这个算例选的是石英、长石和云母,对应的是中酸性岩浆岩或变质岩的矿物组成。三种矿物在弹性模量、强度、颗粒形状上差异明显:石英最硬最脆,长石居中,云母是典型的片状软弱矿物。这种差异越大,GBM的价值越能体现出来。如果三种矿物力学性质差不多,那和均质模型没什么区别,GBM的优势就没了。
对类岩石材料来说也一样,比如水泥基材料里的骨料、硬化水泥浆体、界面过渡区,本质上也是“三种组分”。所以这个模型结构完全可以迁移去模拟水泥砂浆、石膏、混凝土类材料,只要把“石英、长石、云母”替换成“骨料、浆体、界面”即可。
2. PFC5.0里搭GBM的两条路线,我选了哪条
2.1 路线一:Voronoi几何+DFN+接触模型替换
PFC5.0已经支持DFN(离散裂隙网络),这是实现GBM最天然的一条路。思路如下。
第一步,先生成颗粒集合体并压实到目标孔隙率,这一步和普通模型没有区别。第二步,借用Voronoi多边形把试件区域切分成若干胞元区域,每个胞元代表一个矿物晶粒。因为PFC本身没有一键生成Voronoi的命令,实际工程中通常有三种做法:用外部软件生成Voronoi图形后导入DXF/DWG作为Geometry;用Fish自己写凸多边形生成算法;或者用DFN的节理段组合去逼近晶粒边界。第三步,遍历所有接触,判断接触两侧的两个颗粒是否属于同一个Voronoi胞元:如果是,这个接触保留为“晶粒内部”接触,赋高强度的接触模型;如果不是,说明这个接触跨越了晶界,把它替换为SmoothJoint(光滑节理)接触模型,赋低强度的晶界参数。
第四步才是矿物物性的差异化。将Voronoi胞元按照面积占比随机分配为石英、长石、云母三类,并打上不同的group标签。这样遍历接触时,除了判断是否跨晶界,还要判断晶粒内部属于哪一组矿物,再赋对应的矿物力学参数。
这条路的好处在于,晶粒的几何形态接近真实岩石的多边形晶粒,晶界几何坐标明确,后续统计沿晶破裂和穿晶破裂的时候能很清晰地区分。另外PFC5.0对DFN有比较成熟的显示和统计工具,晶界上的裂纹可以直接用DFN节点关联起来。
2.2 路线二:Clump/Cluster聚合矿物晶粒
另一类常见做法是用Clump或者Cluster把一组颗粒“粘”成一个不可破坏的矿物晶体,再用另一个接触模型连接不同的Clump。比如先用ball distribute生成随机颗粒,再按矿物区域给ball分组,然后用cluster或clump命令让同一组内的颗粒形成一个刚体或强粘结体。
这条路实现起来在代码上更直接,因为PFC的clump命令是内置的,不需要外部几何文件。但实际用起来有几个麻烦。Clump本身是刚体,内部不能破裂,所以“穿晶破裂”在Clump内部永远不可能出现;Cluster虽然保留了颗粒间接触,但晶粒的边界形态取决于最初如何给颗粒分组,往往不是标准的多边形,而是锯齿状颗粒轮廓。晶粒边界不规则程度太高,结果是沿晶破裂的路径会被颗粒排列的几何形态干扰,削弱GBM模型本应展示的矿物学特征。再加上矿物区域划分如果太细碎,Cluster的生成和计算都会变得非常不稳定。
2.3 我的选型结论和理由
这两条路线我实际都试过。我的结论是,如果是做研究、发论文,或者想真正还原矿物级别的破坏机理,第一路线优先。原因很简单:Voronoi+DFN的晶粒几何可控、物理意义清晰,而且穿晶、沿晶两种破裂可以在后处理中通过“裂纹是否落在晶界上”直接判别。这个优势在分析破坏模式时几乎是决定性的。
Clump/Cluster路线更适合快速验证,尤其是当你的目标只是粗略看不同矿物占比对强度的影响,并不关心裂纹具体从哪个晶粒穿过。它胜在代码量小、入门快,但它的“天花板”很低——穿晶破裂机制的缺失会让结果在细观尺度上看不太可信。当然,如果模型只是“类岩石”而不追求矿物学精度,比如模拟均匀砂浆材料,那Cluster路线完全够用。
还有一条折中路线:用Voronoi几何做晶粒划分,但不用DFN,直接通过判断接触两侧颗粒的group是否相同来决定赋什么接触模型。这个方法不需要额外设置DFN几何,实现难度大幅降低。唯一要小心的是,没有DFN实体的话,晶界上裂纹的可视化需要额外做映射处理。我这次给的代码框架用的是Voronoi+Voronoi几何判断的思路,属于在原生能力和代码复杂度之间取了一个平衡点。
3. 代码核心:从空模型到可加载的三矿物试件
3.1 模型尺寸、颗粒粒径和孔隙率的设定逻辑
模型尺寸的选择会直接影响计算时间和尺寸效应。单轴压缩2d模型,建议宽高比取1:2,比如宽50mm、高100mm,这个比例能保证试件中部破裂形态少受端部约束影响。颗粒半径不建议太小,否则颗粒数量爆炸;也不建议太大,否则破裂路径太粗糙。我常用的是半径0.25mm到0.4mm的均匀随机分布,配合12%到15%的孔隙率。颗粒数大致在8000到15000之间,PFC5.0单核跑起来压力不大。如果你用更高端的CPU,还能开多线程。
粒径比例max/min不能太大,控制在1.6:1以内比较合适。比例太大,大颗粒会主导破裂路径,小颗粒只能填充孔隙,矿物晶粒的随机性会被颗粒尺寸的不均匀性淹没。配位数在压实后稳定在3.5到4.2左右就行,配位数过低说明颗粒太稀疏,应力传递网络不完整;过高则模型偏刚性,强度偏高。
3.2 三种矿物区域的几何划分与物性赋值
矿物区域划分的做法是:先将模型区域用Voronoi多边形分成若干胞元,然后按面积比例把这些胞元随机指派给石英、长石、云母三种矿物。这里要用到ball的group属性,遍历每个ball,判断其中心落在哪个矿物胞元内,就给它打到哪个group里。
关键点在于,接触模型的赋值要同时区分“矿物内部接触”和“晶界接触”。内部接触的参数根据颗粒所在group决定,晶界接触则统一用弱化的SmoothJoint参数。用一个Fish循环遍历所有contact,读取contact两侧ball的group,然后对接触模型和参数进行差分赋值。
推荐的比例可以这样设:石英40%、长石35%、云母25%,这是花岗岩里比较典型的三矿物体积分数。如果你的研究需要不同矿物配比,直接改这个概率比例重新打标签即可。
下面是一个简化版的Fish逻辑框架,实际运行时要根据你的模型做变量替换。
; 假设已经生成颗粒并压实到目标孔隙率 ; 第一步:按矿物胞元给颗粒分组 loop foreach b ball.list if geometry.in.polygon(b, polyQuartz) then ball.group(b, 'quartz') else if geometry.in.polygon(b, polyFeldspar) then ball.group(b, 'feldspar') else if geometry.in.polygon(b, polyBiotite) then ball.group(b, 'biotite') endif endloop ; 第二步:遍历接触并赋值 loop foreach ct contact.list b1 = contact.end1(ct) b2 = contact.end2(ct) g1 = ball.group(b1, 1) g2 = ball.group(b2, 1) if g1 == g2 then ; 晶粒内部接触,使用平行黏结 contact.group(ct, 'intra') if g1 == 'quartz' then ; 赋石英参数 else if g1 == 'feldspar' then ; 赋长石参数 else ; 赋云母参数 endif else ; 跨晶界接触,使用SmoothJoint弱面参数 contact.group(ct, 'inter') contact.model(ct, 'smoothjoint') ; 赋晶界参数 endif endloop这一段是整套代码的核心。很多人在这一步容易顺手用“按接触的均值属性赋值”来简化,结果就是晶粒边界体现不出来,GBM直接变回了准均质模型。正确的做法必须显式遍历接触、逐条判断group。代码看起来略啰嗦,但这是GBM模型信息量的根本来源。
3.3 伺服加载与单轴压缩的FISH实现
单轴压缩采用位移控制方式。顶部墙和底部墙作为加载墙,两侧墙在加载开始前删除或设为无摩擦约束。PFC里最常见的稳定方案是用伺服机制动态调整墙速——每次时间步根据当前墙体接触力去更新墙速,实现恒应变率加载。
伺服控制的核心逻辑是:设定一个目标加载速度;每一步根据墙体的实时接触力或当前应变率,用比例控制更新墙体的速度。简单点也可以直接给顶部墙一个恒定速度,比如0.05m/s,但注意这个速度是墙的速度,实际施加在模型上的应变率要除以试件高度。100mm高的试件,墙速0.05m/s对应应变率0.5/s,对准静态加载来说已经明显偏快了。建议准静态模拟时墙速控制在0.005m/s量级,也就是应变率0.05/s以下,模拟出来的曲线更接近实验室准静态结果。
一个简易的伺服FISH代码如下:
def servo ; 计算当前墙体接触力 f_curr = wall.force.contact.y(wall.id('top')) ; 计算目标接触力 f_target = stress_target * width ; 比例控制 v_gain = 0.3 delta_v = v_gain * (f_target - f_curr) / f_target wall.vel.y(wall.id('top')) = wall.vel.y(wall.id('top')) + delta_v end伺服代码要反复调试,v_gain太大容易振荡,太小则加载速率上不去、时间步数暴涨。实际中我一般从0.2到0.5起步,根据应力曲线的平滑程度微调。
3.4 应力应变、裂纹和破坏形态的记录方案
单轴压缩要记录的核心数据是轴向应力、轴向应变、裂纹数量、裂纹类型以及最终破坏形态。轴向应力通过墙体的接触力除以截面积得到,二维情况下厚度取单位1,截面积就是试件宽度。应变用上下墙的位移差除以试件初始高度。
PFC5.0里stress和strain的history需要自己写Fish函数,没有内置的“单轴压缩”一键输出。建议用history命令把应力、应变、裂纹数都记录下来。裂纹统计要分两类:穿晶破裂(intra-granular crack)和沿晶破裂(inter-granular crack)。这两个数据的区分方式,是在裂纹产生时将裂纹所属的接触类型(内部接触还是晶界接触)写入crack的group标签。这样在后面处理时,可以一键按group统计两种裂纹的占比。
破坏形态的查看,最直观的是在PFC5.0的View中用crack显示裂纹,并按裂纹的剪切/拉伸属性着色。拉伸裂纹通常用白色短线表示,剪切裂纹用黄色短线表示。三矿物模型的试件破裂后,典型的输出是:边缘云母区域先出现分散的微裂纹,随后沿晶界扩展成弧形裂缝,最终在长石和石英晶粒附近形成穿晶锐利裂纹,两条系统汇合成主破坏面。看到这种形态时,你的GBM才算真正搭建成功了。
4. 矿物参数怎么标,晶界参数怎么调
4.1 参考参数区间与物理含义
GBM最耗时的环节是标定,而标定的第一步是先知道各矿物参数的合理区间。这里给出一组我常用的参考范围,单位统一用Pa和N/m,弹模、刚度和强度的换算要注意PFC里面的力和位移单位一致。
| 参数项 | 石英 | 长石 | 云母 | 晶界 |
|---|---|---|---|---|
| 密度 (kg/m³) | 2650 | 2600 | 2800 | - |
| 接触模量 (GPa) | 60~70 | 45~55 | 15~25 | 5~10 |
| 刚度比 kn/ks | 1.5~2.0 | 1.5~2.0 | 1.0~1.5 | 1.0左右 |
| 抗拉强度 (MPa) | 80~130 | 50~90 | 15~30 | 5~15 |
| 内聚力 (MPa) | 120~180 | 80~120 | 25~50 | 10~20 |
石英的接触模量和抗拉强度最高,这个不用怀疑,矿物学硬度摆在那里。长石略低,云母则明显偏软偏弱。晶界参数整体比矿物强度低一个量级,这样才能保证裂纹优先在晶界产生。
但要特别提醒,PFC里的微观接触参数和宏观弹性模量不是一一对应的,接触模量设了60GPa不代表模型宏观模量就是60GPa。宏观模量还受孔隙率、配位数、颗粒粒径影响。所以标定时不能只看接触参数数值大小,要以宏观应力应变曲线的斜率为准反推。
4.2 标定的三层流程
三层标定是我用下来最稳定的流程,每一层都在上一层的参数基础上微调。
第一层,标单矿物。先跑一个只有石英颗粒的纯矿物试件,通过单轴压缩反演出石英的接触模量、刚度比、抗拉强度、内聚力。目标是用这些微观参数跑出来的宏观弹性模量和强度符合石英岩的参考值。长石、云母各自重复这个过程。
第二层,标晶界。把三种矿物拼起来,固定矿物内部参数不变,只调晶界参数。晶界强度调高,模型表现为脆性劈裂,破坏形态像均质材料;晶界强度调低,模型表现为大量沿晶微破裂,峰前非线性明显,甚至还没到峰值就已经产生大量离散裂纹。调整的目标是让模型的破坏模式呈现出“硬矿物穿晶破裂和软矿物沿晶破裂混合”的状态。
第三层,标组合。将完整三矿物试件的宏观强度和破坏形态与实验对比。如果峰值强度偏高,优先降低云母强度和晶界强度,而不要动石英参数;如果破坏形态太碎,通常是晶界强度太低,要适当上调。依次反馈迭代,通常三轮之后能匹配到可接受的误差范围。
4.3 这组参数对结果的敏感性表现
参数敏感性分析很多人会忽略,但它在调参时非常管用。我做过一组对比试验:固定其他参数,只把晶界抗拉强度从8MPa提到16MPa,宏观峰值强度会明显提升,而且破坏形态从“弥散型微破裂”变成“集中型剪切破坏”,这说明晶界强度对宏观破坏模式的控制力极强。
云母的接触模量对宏观弹性模量影响不大,但对峰值强度影响不小。这是因为云母作为软弱矿物,在应力重新分配过程中承担了更多塑性变形和局部破坏,相当于“应力缓冲带”。把云母强度调高一点,应力缓冲效果变差,穿晶破裂比例上升,宏观脆性增强。反之,云母强度调低,模型会表现出更显著的延性特征,峰前非线性段变长,峰后跌落变缓。
这些敏感性关系是GBM模型的核心语言。掌握之后,你不需要再盲目地六七个参数一起调,而是知道“破坏形态不对应该动哪个参数”“弹性模量不对应该动哪个参数”“峰后跌落太陡应该往哪个方向调”。
5. 踩坑实录:五处最容易翻车的地方
5.1 加载墙速度过大,曲线全在抖
伺服加载刚开始跑时,应力应变曲线如果出现高频锯齿状波动,先别急着怀疑接触参数。十有八九是加载速度太大,模型内部惯性力不可忽略,系统始终处于非平衡状态。
判断方法很简单:在计算过程中实时监测不平衡力比(unbalanced force ratio),如果这个比值长期高于1e-3,说明模型内部力传递远没有达到准静态平衡,应力波的扰动正在主导破坏过程。解决办法是把墙速降下来,例如从0.05m/s降到0.005m/s,然后重新看曲线。准静态单轴模拟的要求是,加载过程中的动能与应变能之比控制在1%以内,否则裂纹起裂的位置和扩展路径都会失真。
5.2 晶界“拉链式”开裂:几何规则化惹的祸
第一次看到自己的GBM试件压坏时,我满心期待看到不规则的破裂路径,结果看到的却是沿晶界整圈“拉链式”开裂,裂纹沿着每一条晶界扫过,试件被分成很多个小多面体碎块。这种破坏形态显然不合理。
问题的根源在晶粒几何的规则化。如果用Voronoi生成晶粒的时候没有加随机扰动,晶粒尺寸接近,晶界近乎等长且规则,那么每一条晶界的受力状态也接近,结果就出现同步破坏。真实岩石的晶粒尺寸、形状都有显著随机性,局部晶界受力不均,有些位置的晶界先裂,另一部分则作为应力桥梁继续承担荷载。
解决办法是给Voronoi种子点加随机扰动,或者用多个随机数种子生成多个晶粒方案做对比。晶粒尺寸系数(也就是平均晶粒面积)也要和颗粒粒径匹配好,一个晶粒至少要包含几十个颗粒,否则颗粒离散误差会把晶粒边界抹平。
5.3 矿物刚度差异导致的应力集中失真
三种矿物的弹性模量差异特别大时,比如石英70GPa、云母20GPa,交界处的刚度突变会引发严重的局部应力集中。这个现象在物理上是合理的,但数值上容易被过度放大,导致云母区域还没到正常强度就提前崩解掉了。
我踩过这个坑之后,建议在赋参数时不要直接照搬矿物宏观看实测值,而是做一次“过渡处理”。可以把云母的接触模量从20GPa适当提到25GPa,让模量对比度稍微缓和;或者让云母的强度不要低到离谱的程度,保持一定的残余承载力。毕竟在PFC里颗粒是圆形的,和真实云母片状晶体的嵌锁效应不同,数值模型需要一定的参数补偿才不会让软弱矿物变成“真空区”。
5.4 裂纹统计口径混乱,导出数据对不上
GBM模型通常同时存在本体裂纹和晶界处裂纹,PFC的crack记录默认只输出裂纹位置和断键类型。如果你没有提前给不同group的接触做标记,后处理时你就分不清哪些裂纹是穿晶的、哪些是沿晶的。
正确做法是,在遍历接触赋值的时候,同时给接触设置一个自定义属性,比如用contact.group(或额外定义contact的额外属性)记录intra_qtz、intra_fsp、inter_biot等标签。裂纹生成时会继承它断裂母体接触的属性,这样导出的crack数据就能按标签直接分类统计。另一个容易忽略的细节是,PFC的DFN节点不一定和裂纹一一对应,靠节点去数裂纹会漏算或重复计算。统计时以crack的接触属性为准,不要用DFN节点数来替代。
5.5 算力不够,调参只能干等
单轴压缩GBM模型的计算量比普通平行黏结模型高出不少,尤其是在遍历接触赋参数、裂纹大量产生的阶段。如果计算机配置一般,一个工况可能跑大半天,调参效率极低。
我试过比较有效的几个提速手段。颗粒粒径适当放大,比如最小半径从0.2mm提到0.3mm,颗粒数量能减少将近一半,对宏观强度的尺寸效应影响在几个百分点内,可以通过后标定补偿回来。另外把模型的domain范围设置紧凑,别留太多空白区,减少颗粒搜索的计算量。还有就是PFC5.0支持多线程,把计算线程数调到与CPU核心数匹配,效果立竿见影。最后建议把不需要的测量圆和history监控全部关掉,它们对计算性能的拖累往往被低估。
6. 这套代码拿到手之后,还能往哪些方向扩展
6.1 换加载路径:巴西劈裂、三轴、循环荷载
单轴压缩只是GBM模型的基础测试。一旦你有了这个三矿物试件框架,换加载方式只是改边界条件的事。巴西劈裂只需把加载墙左右放置,并沿径向加载;三轴则需要给侧向边界加围压伺服控制,用PFC5.0的wall或ball-based伺服都能实现。很多做地下工程的同行还会把单轴改成循环加卸载,用GBM观察循环荷载下的疲劳裂纹扩展规律,重点研究裂纹何时在晶界处萌生、如何随循环次数累积。
注意:换加载路径后,原有的矿物参数不一定能直接复用。巴西劈裂对抗拉强度参数敏感,三轴对抗剪强度参数敏感,你可能需要在原参数基础上各做一个局部重新标定。但三矿物试件本身的几何和参数框架不用动,改起来成本很低。
6.2 提取AE事件和能量耗散
GBM模型与实验室的声发射(AE)监测有天然对应关系。每一条裂纹产生,就对应一次微破裂事件,而经典声发射定位实验记录的正是岩石内部微破裂的位置和时序。
把这个思路落到代码里,最简单的方式是在循环中实时监测crack的产生,将每一轮的裂纹坐标按时间窗口聚类,作为一个“AE事件”,统计事件的空间位置和能量。能量可以用裂纹断裂时的接触应变能释放量的近似值来估算。有了AE事件的时间和空间信息,你就可以对照实验室中试样破坏前的声发射定位图,看AE事件是否先在软矿物区域聚集,再向硬矿物区域迁移,最后在破坏前密集连成核。这种验证方式很有说服力。
6.3 微观参数的尺寸效应与工程参数换算
PFC里有一个老问题:微观参数标定结果依赖颗粒尺寸。同一套接触参数,颗粒放大之后宏观强度会变。GBM模型同样存在这个问题,甚至因为晶粒尺寸和颗粒尺寸的双重影响,尺寸效应更明显。
扩展方向之一是系统研究晶粒尺寸与颗粒粒径的比值对宏观强度的影响规律。一般建议晶粒内部至少包含30到60个颗粒,小于这个阈值时破裂路径会出现明显的“网格依赖”。另一个方向是把微观参数标定结果和工程岩体参数建立关联,比如通过Take a series of不同尺度模型的模拟结果,外推出REV(代表性体积单元)尺寸和对应的宏观强度参数,这对工程尺度的数值计算非常有帮助。
还有一个实用技巧:在做参数敏感性分析的时候,用GBM模型做“虚拟实验”,代替一部分昂贵的室内实验。先在一个合理的参数空间里生成几十组模型,跑出宏观强度范围,再把实验值投影进去,反找最匹配的参数组合。这样比反复试错标定高效得多。
这套代码我前后调了两周,其中一半时间都花在参数标定和排查晶界赋值bug上。现在回头看,GBM真正的难点不在PFC命令本身,而在于对“矿物晶粒—晶界”这个结构层级的理解。刚开始如果你拿着代码跑出来的结果和预期不符,不要怀疑模型框架错了,先按第五部分的踩坑清单逐一排查。把三矿物试件的单轴压缩跑通之后,你会发现自己对PFC、对岩石破裂过程的理解都会上一个台阶。