PFC-FLAC耦合技术在解决这种大尺度变形上就有很高的价值了,将非大变形区域的颗粒用zone来模拟,而我们将颗粒数花费在大变形破坏区域,这样就可以提高模拟的价值了。
一、二维的耦合测试
下方是耦合的测试代码,FLAC2D的关键词并没有在PFC的手册中,各位可以参考ITASCA官网的关键词手册(Itasca Software Guide — Itasca Software 9.0 documentation (itascacg.com))。这里我们用quad关键词生成一个矩形的地基,地基首先进行地应力沉降。注意边界设置,底部固定y向,左右固定x向即可。
沉降完成后,生成一个密度较高的颗粒落在地基上,非常简单的案例,测试足够。
model newdomain extent -2 2zone create quad size 50 5 point 0 0 0 ...point 1 1 0 ...point 2 0 0.2zone gridpoint fix velocity-y range pos-y -0.001 0.001zone gridpoint fix velocity-x range pos-x -0.001 0.001zone gridpoint fix velocity-x range pos-x [1-0.001] [1+0.001]zone initialize density 3e3zone cmodel assign elasticzone property young 10e6 poisson 0.3model gravity 9.8model cycle 1model solvemodel save "diji"model restore "diji"wall-zone createball create position 0.5 0.4 radius 0.02contact cmat default model linear property kn 1e9 ks 1e9 fric 0.5ball attribute density 3e6 damp 0.7model cycle 2000
可以看到耦合效果还是可以的,在落点位置出现较大的位移。

这个测试完成后,就尝试进行一个略复杂的砂土边坡模拟了。
二、颗粒部分生成
我们先生成边坡的颗粒部分,这也是我一般的做法,离散-有限耦合的时候,先完成离散部分的成样,然后再扩展有限元部分的元素。我们将坡脚定义在原点位置,然后定义好边坡的宽和高。里面还有个数值控制着离散元左边界与边坡左边界的距离 width*0.1 ,这个数值越大,离散的部分就越小,这个数值在设置的时候应当充分考虑边坡的大变形范围。
还有一点是distribute成样是先在domain范围内生成颗粒,然后删除range以外的部分,所以domain应当足够小。
model newdef parwidth=2.0height=width*0.5rdmin=0.001rdmax=0.002poro=0.08end@parmodel domain extent [-width*0.5] [width*0.1] [-width*0.1] [height*0.6]wall generate box [-width*0.5+width*0.1] [0] [0] [height*0.5]ball distribute porosity @poro radius @rdmin @rdmax box ...[0] [0] [height*0.5]ball attribute density 3e3 damp 0.7cmat default model linear method deform emod 10e6 kratio 1.5 property fric 0.5model cycle 2000 calm 20model solvemodel save "sample_keli"
三、有限元部分生成
离散有限部分耦合的关键点是要保证材料的一致性,这个其实难度很大,对参数需要足够多的时间去调整。对于这种考虑重力的工况而言,除了常规的力学参数,主要是变形参数保持基本一致外,密度是需要调整的。PFC的密度是给到的ball,同样的式样体积下,zone部分的固体体积和式样体积一致,但是ball部分由于孔隙率和重叠量的存在,总的质量是有浮动的,所以需要调整密度,使得zone部分的地应力或者变形场基本水平即可。
model restore "sample_keli"model mechanical time-total 0ball attribute displacement multiply 0geometry set "zone_range"ball attribute density 2.8e3 damp 0.7zone create quad size 2 10 point 0 [-width*0.5] [-height*0.5] ...point 1 [-width*0.5+width*0.1] [-height*0.5] ...point 2 [-width*0.5] [height*0.5]zone create quad size 10 5 point 0 [-width*0.5+width*0.1] [-height*0.5] ...point 1 [0] [-height*0.5] ...point 2 [-width*0.5+width*0.1] [0]zone create quad size 10 5 point 0 [0] [-height*0.5] ...point 1 [width*0.5] [-height*0.5] ...point 2 [0] [0]zone create quad size 10 5 point 0 [0] [0] ...point 1 [width*0.5] [0] ...point 2 [0] [height*0.5]zone initialize density 3e3zone cmodel assign elasticzone property young 10e6 poisson 0.3cmat default type ball-facet model linear method deform emod 10e6 kratio 1.5 property fric 0 lin_mode 1cmat default type ball-ball model linear method deform emod 10e6 kratio 1.5 property fric 0 lin_mode 1cmat applymodel cleancontact property lin_force 0 0ball attribute force-contact 0 0wall deletewall create id 100 vertices [0] [0] [0] [height*0.6]wall create id 101 vertices [-width*0.5+width*0.1] [0] [-width*0.5+width*0.1] [height*0.6]wall-zone create range pos-x [-width*0.5+width*0.1-0.0001] [0.0001] ...pos-y [-0.00001] [height*0.5+0.0001]zone gridpoint fix velocity-y range pos-y [-height*0.5-0.0001] [-height*0.5+0.0001]zone gridpoint fix velocity-x range pos-x [-width*0.5-0.0001] [-width*0.5+0.0001]zone gridpoint fix velocity-x range pos-x [width*0.5-0.0001] [width*0.5+0.0001]model gravity [9.8*40]model cycle 1000model solve ratio-average 1e-3wall delete walls range id 100 101model solvemodel save "zizhong"

将色条设为一致,可以发现这里的密度还不是一个非常好的状态,颗粒底部的位移明显小了点,说明密度还需要增大。而颗粒内部的连续并不需要过于在意,由于颗粒洒落的时候并不是连续变形,会导致初始时候中间颗粒的位移较大。
而应力也会有相应的变化,一般而言并不会完完全全的呈现水平带状分布,调整参数到差不多的状态就可以了。

四、开挖
这里开挖主要删除两个部分,第一是右上方的zone,第二十根据坡度产生的斜半平面以上的颗粒。常规有限元的生死单元和设置模型是null的方法,在这里只需要删除颗粒就可以了,这个还是比较好的。
注意一下初始位移量清零,还有一点是ball-ball接触参数应和耦合部分的ball-facet参数保持一致。
model restore "zizhong"model mechanical time-total 0model domain extent [-width*0.5] [width*0.6] [-height*0.6] [height*0.6]ball attribute displacement multiply 0zone gridpoint initialize displacement 0 0zone delete range pos-x 0 [width*0.5] pos-y 0 [height*0.5]wall-zone create range pos-x [0] [width*0.6] ...pos-y [-0.2] [2]wall property "fric" 0.5ball property "fric" 0.5ball delete range plane origin 0 0 normal ...[math.cos(math.pi/6.0)] [math.sin(math.pi/6.0)][final_time=2][baocunpinlv=final_time/20.0][time_record_savefile=-1][count=0]def savefiletime=mech.time.totalif time-time_record_savefile >= baocunpinlv thenfilename=string.build("jieguo%1",count)commandmodel save @filenameendcommandtime_record_savefile=timecount +=1endifendfish callback add @savefile -0.5model cycle 1model solve time @final_time
测试的时候把收敛时间设置为1s,各位可以根据计算情况继续往下面算。下面为1s后的y向位移场,在分析的时候应当将色卡设置一致,可以看到总的来说,变形的连续性还尚可,调整变形参数应当会有更好的结果。

做成动图效果也是非常好看的。

七、Lobby的PFC6.0入门学习讲义
在此向读者朋友推荐我今年新上线的新课程《Lobby的PFC6.0入门学习讲义》适合新手入门。

本课程为PFC6.0软件入门,常用的案例技术解析;全字幕的课程,由简入深的讲解,用户可以获得多达二十种案例,丰富的技术点讲解,完全从零开始的案例实操,覆盖clump/cluster/rblock/离散连续耦合/流固耦合/柔性膜等内容,让用户掌握基本的后处理技能。
以下是我的课程安排




扫码立即试看
(完)