PFC及几种方法简介201130410019一、PFC简介1、pfc背景PFC系列软件是由ITASCA咨询集团(设有ITASCA中国公司)开发的颗粒流分析程序(ParticalFlowCode),分为PFC2D,PFC3D两种特别用于模拟任意性状、大小的二维圆盘或三维球体集合体的运行及其相互作用的强大颗粒分析程序。除了模拟大体积流动和混合材料力学研究,程序更适合于描述在固体材料中细观/宏观裂纹扩展、破坏累积并断裂、破坏冲击和微震响应等高水平课题的深化研究。PFC中的颗粒为刚性体,但在力学关系上允许重叠,以模拟颗粒之间的接触力。颗粒之间的力学关系非常简单,即牛顿第二定律。颗粒之间的接触破坏可以为剪切和张开两种形式,当介质中颗粒间的接触关系(如断开)发生变化时,介质的宏观力学特性受到影响,随着发生破坏的接触数量增多,介质宏观力学特性可以经历从峰前线性到峰后非线性的转化,即介质内颗粒接触状态的变化决定了介质的本构关系。因此,在PFC计算中不需要给材料定义宏观本构关系和对应的参数,这些传统的力学特性和参数通过程序自动获得,而定义它们的是颗粒和水泥的几何和力学参数,如颗粒级配、刚度、摩擦力、粘结介质强度等微力学参数。2、基本假设1)颗粒单元为刚性体;2)接触发生在很小的范围内,即点接触;3)接触特性为柔性接触,接触处允许有一定的“重叠”量;4)“重叠”量的大小与接触力有关,与颗粒大小相比,“重叠”量很小;5)接触处有特殊的连接强度;6)颗粒单元为圆盘形(或球形)。3、优缺点优点:第一、它有潜在的高效率。因为圆形物体间的接触探测比角状物体间的更简单。第二、对可以模拟的位移大小实质上没有限制。第三、由于它们是由粘结的粒子组成,块体可以破裂,不象UDEC和3DEC模拟的块体不能破裂。缺点是:块体的边界不是平的,用户必须接受不平的边界以换取PFC2D提供的优点。4、求解步骤1)定义模拟对象根据模拟意图定义模型的详细程序,假如只对某一力学机制的不同解释作出判断时,可以建立一个比较粗略的模型,只要在模型中能体现要解释的机制即可,对所模拟问题影响不大的特性可以忽略。2)建立力学模型的基本概念首先对分析对象在一定初始特性形成初步概念。为此,应先提出一些问题,如系统是否将变为不稳定系统、问题变形的大小、主要力学特性是否非线性、是否需要定义介质的不连续性、系统边界是实际边界还是无限边界、系统结构有无对称性等3)构造并运行简化模型在建立实际工程模型之前,先构造并运行一系列简化的测试模型,可以提高解题效率。通过这种前期简化模型的运行,可对力学系统的概念有更深入的了解,有时在分析简化模型的结果后(例如所选的接触类型是否有代表性、边界条件对模型结果的影响程度等),还需将第二步加以修改4)补充模拟问题的数据资料模拟实际工程问题需要大量简化模型运行的结果,对于地质力学来说包括:a)几何特性,如地下开挖酮室的形状、地形地貌、坝体形状、岩土结构等;b)地质构造位置,如断层、节理、层面等;c)材料特性,如弹/塑性、后破坏特性等;d)初始条件,如原位应力状态、孔隙压力、饱和度等;e)外荷载,如冲击荷载、开挖应力等。5)模拟运行的进一步准备a)合理确定每一时步所需时间,若运行时间过长,很难得到有意义的结论,所以应该考虑在多台计算机上同时运行。b)模型的运行状态应及时保存,以便在后续运行中调用其结果。例如如果分析中有多次加卸荷过程,要能方便地退回到每一过程,并改变参数后可以继续运行。c)在程序中应设有足够的监控点(如参数变化处、不平衡等),对中间模拟结果随时作出比较分析,并分析颗粒流动状态。6)运行计算模型在模型正式运行之前先运行一些检验模型,然后暂停,根据一些特性参数的试验或理论计算结果来检查模拟结果是否合理,当确定模型运行正确无误时,连接所有的数据文件进行计算。7)解释结果计算结果与实测结果进行分析比较。图形应集中反应要分析的区域如应力集中区,各种计算结果应能方便地输出分析。二、PFC2D计算模型的几种生成方法1、有两个命令可用于生成颗粒流模型:BALL和GENER-ATE,其中,BALL命令是生成单个的颗粒,该命令生成的颗粒可与已存在的颗粒重叠,而GENERATE可生成一系列指定数目的颗粒流,该命令生成的颗粒是不允许重叠的。PFC2D里主要有两种类型的颗粒流:规则排列的和无规则排列的。尽管颗粒的排列是随机的,但在颗粒模型生成后,整个模型的结构特性还是可能会受影响的,比如弱的结构面或各向异性。对于无规律排列的颗粒流模型,一般不可能去描述它的初始接触力的量级大小,这必须在后期要经过一个压缩的过程才可能给予较好的评价。1.1规律排列颗粒流Newdefhexxc=x0yc=y0rc=radiusidc=id_startr2=2.0*radiusyinc=radius*sqrt(3.0)looprow(1,n_row)loopcol(1,n_col)commandballid=idcx=xcy=ycrad=rcend_commandidc=idc+1xc=xc+r2end_loopyc=yc+yincxc=x0+radius*(row-(row/2)*2)end_loopendsetechooffsetx0=0.2y0=0.4radius=0.1setid_start=100n_col=7n_row=8hexsetechoonplotsetcapsize20plotaddaxesblackplotaddballyellowplotshow1.2不规则排列无规则排列,即:对一个给定空隙率的区域,采用颗粒来充填其中需要进行填充的空隙,并确保整个模型保持平衡。对于所能被填充的模型的初始空隙率,是有一个限制值,不能任意小。对于某些空隙率的模型,颗粒的填充可以无接触地排列,对于其它情况的空隙率,颗粒又可以重叠排列。第一种方法,首先建立封闭区域的边界(简称墙体),然后在封闭区域内任意生成一系列无接触的颗粒,最后移动区域的限制墙体,至所需要的空隙率。这种方法有三个缺点:1.区域的几何形状改变;2.收敛速度慢;3.最终的分布趋势是不均匀的第二种方法:运用GENERATE命令生成颗粒体,同时配合关键词高斯分配,即指定颗粒体半径的上下限,然后相应分配一个标准差,同时配合FISH函数来选择颗粒半径,最终生成我们所需要的模型。半径扩展法newdefexpand;---输入数据---n_stiff=1e8;法向连接刚度s_stiff=1e8;剪切连接刚度width=10.0;区域宽height=5.0;区域高tot_vol=width*heightporos=0.12;最终目标空隙率num=300;颗粒体数目rat=1.5;最大最小半径比;---导出所需数据---mult=1.6;初始半径放大系数n0=1.0-(1.0-poros)/mult^2r0=sqrt(height*width*(1.0-n0)/(pi*num))rlo=2.0*r0/(1.0+rat)rhi=rat*rlo_x1=width*(1.0+extend)_y1=0.0commandwallid=1ks=s_stiffkn=n_stiffnodes(_x0,_y0)(_x1,_y1)end_command_x0=width_y0=-extend*height_x1=width_y1=height*(1.0+extend)commandwallid=2ks=s_stiffkn=n_stiffnodes(_x0,_y0)(_x1,_y1)end_command_x0=width*(1.0+extend)_y0=height_x1=-extend*width_y1=heightcommandwallid=3ks=s_stiffkn=n_stiffnodes(_x0,_y0)(_x1,_y1)end_command_x0=0.0_y0=height*(1.0+extend)_x1=0.0_y1=-extend*heightcommandwallid=4ks=s_stiffkn=n_stiffnodes(_x0,_y0)(_x1,_y1)end_command;---generatetheballsandgivethemtheirpropertiescommandgenid=1,numrad=rlo,rhix=0,widthy=0,heightpropdens=1000ks=s_stiffkn=n_stiffend_commandget_porosmult=sqrt((1.0-poros)/(1.0-pmeas))commandiniradmulmultcycle1000propfric0.2cycle2000end_commandenddefget_porossum=0.0bp=ball_headloopwhilebp#nullsum=sum+pi*b_rad(bp)^2bp=b_next(bp)end_looppmeas=1.0-sum/(width*height)endexpandget_porosplotwallballplotshowprintpmeassaveexpand.SAV挤压排斥法指定颗粒体的半径,不限制颗粒的数目,使足够多的颗粒产生来达到所需要的空隙率。但这种方法所带来的缺点是可能在局部区域造成大面积的颗粒重叠,这将会产生很大的挤压力,从而给予颗粒较大的初始速度,这就可能使得颗粒体脱离墙体的限制。为避免此情况的发生,可通过初始的有限步循环计算将颗粒的动能减至零,然后再计算至平衡态。newsetrandom;设定颗粒体数目无限制模式defexploden_stiff=1e8s_stiff=1e8width=10.0height=5.0poros=0.12n_max=1000rlo=0.174rhi=0.261commandwallid1ks=s_stiffkn=n_stiffnodes(0,0)(width,0)wallid2ks=s_stiffkn=n_stiffnodes(width,0)(width,height)wallid3ks=s_stiffkn=n_stiffnodes(width,height)(0,height)wallid4ks=s_stiffkn=n_stiffnodes(0,height)(0,0)end_commandpvol_sum=0.0tot_vol=width*heightcount=0sectionloopn(1,n_max)r_ball=rlo+urand*(rhi-rlo)pvol_new=pvol_sum+pi*r_ball^2if(1.0-pvol_new/tot_vol)porosthenexitsection;(newporositywillbethatspecified)end_ifrb2=r_ball*2.0x_ball=r_ball+urand*(width-rb2)y_ball=r_ball+urand*(height-rb2)commandballx=x_bally=y_ballrad=r_ballend_commandpvol_sum=pvol_newcount=count+1end_loopend_sectiontot_poros=1.0-pvol_sum/tot_volcommandpropdens=1000.0ks=s_stiffkn=n_stiffend_commandii=out(''+string(count)+'ballswerecreated.')enddefquietloopn(1,4)commandcycle5inixvel0yvel0end_commandend_loopendexplodequietcyc980propfric=0.2cyc250plotwallballprinttot_porossaveexplode.SAV三、边界条件PFC2D中有三种边界条件,分别是:墙体边界、颗粒体边界和混合边界。其中颗粒体边界又分为