首页
材料特性
使用命令 PROPERTY 可以分配材料特性参数信息。难点在于选择能反映真实材料的特
性参数。连续介质程序中的材料参数是宏观参数，可以由实验室得到，而PFC2D 则是采用
颗粒之间的微观参数。如果材料的微观参数已知，则可以直接应用PFC2D 模拟；如果微观
参数未知，则需要采用逆向模拟，即进行许多试验和模拟，建立宏观参数和微观参数的关系。
当不连续体（如裂缝、断层、节理等）要设置在材料中进行模拟时，可以使用节理生成器引
入弱面。
颗粒和接触的特性
除了与粘结相关的特性以外，都要使用PROPERTY 命令给颗粒分配特性参数。有些参
数，如摩擦系数和接触刚度，与接触行为相关，但是也必须分配给颗粒。这些接触行为由组
成接触的两个颗粒的特性得到。
特性在空间的变化
参数沿空间线性变化：prop s_bond 1e6 grad -1e4 range x 0 50。即，接触粘结剪切强度由
x=0处的1e6 线性变化到x=50处的5e5（10e6 − 1e4*50）。
参数按照更复杂的函数变化则需要通过FISH 函数完成。
参数按组分配：使用命令GROUP和RANGE GROUP。GROUP将某空间范围内的颗粒
定义为一组后，即使以后这些颗粒移动到初始定义的RANGE 空间范围外，仍然保持在一个
组内，仍然可以对其同时分配或改变参数。需要注意的是，通过GROUP设置接触特性时，
只有GROUP内部球体间形成的接触可以被改变，而属于不同GROUP的颗粒间形成的接触
不受影响。如果需要改变不同GROUP之间的颗粒形成的接触的特性，则需要使用JSET命
令，这将在后面的节理平面部分给出详细说明。
基于已知微观特性的直接模拟
如果模拟的材料由圆颗粒组成，则颗粒的特性可以直接输入PFC2D。各参数的相对重
要性取决于模拟的类型和目的。如果颗粒材料小应变时的整体模量很重要，则应该使用
HERTZ-MINDLIN 模型，输入剪切模量和孔隙比；如果材料的强度比模量更重要，则使用
计算效率高的线性接触模型，这种情况下，如果颗粒的重叠在合理范围（比如小于平均半径
的5%），则切向和法向接触刚度就不重要了；如果没用颗粒间粘结，则材料的强度决定于
摩擦系数。
对于矩形边界的颗粒组，可以建立模量和强度的解析表达式。在解析表达式中，存在
一些普遍的关系：第一，二维颗粒组的整体模量与颗粒刚度成正比，但与颗粒半径无关；第
二，对于使用接触粘结模型的颗粒组，粘结对整体强度的贡献（区别于摩擦的贡献，即抗拉
强度）与平均颗粒半径成反比。平行粘结对整体行为的影响更加复杂。
未知微观特性时的逆向模拟
先做一系列物理试验，得到材料的各种宏观参数，然后不断调整数值试验的输入参数，
使数值模拟结果与物理试验结果相吻合，然后建立物理试验参数与数值试验参数之间的关
系，用于数值方法模拟更多的大型问题。这个过程有点想是撞运气，但是我们可以根据经验
缩小参数的范围。
双轴试验
任何合成材料的整体特性可以通过对颗粒材料进行一系列试验获得。这些试验以数值形
势进行并模拟物理试验。
双轴试验PFC模型
合成材料试样由一组圆颗粒组成，其局限性主要有两点：首先，PFC2D 模型只发生平
面应力或平面应变，这与真实材料的三轴试样中的响应不同；其次，为了计算应力，假定颗
粒都是单位厚度。
该双轴试验的模拟中，用四个墙限制矩形试样，顶墙和底墙模拟加荷板，左墙和右墙
模拟试样侧限。试样应变控制形式加荷，指定顶板和底板的速度。在试验的所有阶段，左右
墙的速度都由数值伺服机制自动控制，以使侧向应力恒定。应力和应变的计算方法是，将所
有作用在相应墙上的力相加，求得宏观上的平均应力和应变。材料的各应力和应变的响应，
可以使用HISTORY 来进行跟踪观察。
试样制备
方法详见《PFC2D 学习笔记之颗粒生成细节》
计算并控制应力状态
应力和应变状态的有FISH 函数GET_SS确定。应力是作用在每组相对的墙上的应力平
均值，而作用在每个墙上的应力则是用作用在该墙上的总的力除以试样长度。应变的计算公
式为2*（L-L0）/(L0+L)。get_ss的核心程序如下：
wadd1 = find_wall(1)
wadd2 = find_wall(2)
wadd3 = find_wall(3)
wadd4 = find_wall(4)
xdif = w_x(wadd2) - w_x(wadd4)
ydif = w_y(wadd3) - w_y(wadd1)
new_xwidth = width + xdif
new_height = height + ydif
wsxx = 0.5 * (w_xfob(wadd4) - w_xfob(wadd2)) / (new_height * 1.0)
wsyy = 0.5 * (w_yfob(wadd1) - w_yfob(wadd3)) / (new_xwidth * 1.0)
wexx = 2.0 * xdif / (width + new_xwidth)
weyy = 2.0 * ydif / (height + new_height)
wevol = wexx + weyy
加载过程中，要使用数值伺服机制，通过调整侧墙的速度使侧向应力保持不变，该伺服机
制通过FISH函数SERVO和GET_GAIN实现。SERVO每个CYCLE中调用一次，它通过GET_SS
函数得到应力，并使用伺服机制调节墙的速度以减小测量应力与目标应力之差。变量Y_SERVO
作为开关，若为零，则伺服机制不作用在顶板和底板。GET_gain函数用于获得GAIN参数”G”。
GET_GAIN的核心代码如下：
alpha = 0.5 ; relaxation factor
count = 0
avg_stiff = 0
cp = contact_head ; find avg. number of contacts on x-walls
loop while cp # null
if c_ball1(cp) = wadd2
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball1(cp) = wadd4
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball2(cp) = wadd2
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball2(cp) = wadd4
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
cp = c_next(cp)
end_loop
nxcount = count / 2.0
avg_stiff = avg_stiff / count
gx = alpha * (height * 1.0) / (avg_stiff * nxcount * tdel)
count = 0
avg_stiff = 0
cp = contact_head ; find avg. number of contacts on y-walls
loop while cp # null
if c_ball1(cp) = wadd1
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball1(cp) = wadd3
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball2(cp) = wadd1
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball2(cp) = wadd3
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
cp = c_next(cp)
end_loop
nycount = count / 2.0
avg_stiff = avg_stiff / count
gy = alpha * (width * 1.0)/ (avg_stiff * nycount * tdel)
SERVO的核心程序如下
while_stepping
get_ss ; compute stresses & strains
udx = gx * (wsxx - sxxreq)
w_xvel(wadd4) = udx
w_xvel(wadd2) = -udx
if y_servo = 1 ; switch stress servo on or off
udy = gy * (wsyy - syyreq)
w_yvel(wadd1) = udy
w_yvel(wadd3) = -udy
end_if
微观特性选择方法
线性接触模型
1、初始弹摸月接触刚度线性相关。平行粘结能增加接触刚度，但接触粘结不能。
2、泊松比取决于取决于压密状态和剪切接触刚度与法向接触刚度的比值。
3、峰值强度取决于摩擦系数和粘结强度。如果只有摩擦系数，材料表现出塑性或微弱软化
行为；如果粘结强度也加上，材料将表现出更大的峰值强度但是也有更大的应变软化。
4、当围压增加时，摩擦强度部分将提高，但是粘结强度部分不变。因此，高围压时材料将
表现出更大的塑性。
5、如果使用接触粘结模型，峰值后的卸载-加载模量仅比初始模量稍低；如果使用平行粘结
模型，则该模量会随着应变的增加而减小；随着平行粘结破坏时，材料不断积累损伤。
6、如果平行粘结刚度与接触刚度的比值增大，损伤累积的速率随着应变的发展儿增大。
7、如果强度值围绕平均值分布而不是只设为单一值，峰值的形状会更加平滑和宽大。
8、试样孔隙比越低，则在峰值后的体积膨胀越大。
9、粘结可以在任何平均应力值处设置；在施加剪切应力前，该应力将增加或降低。在平均
应力或偏应力处设置粘结的作用是，为了锁闭接触力和内能。这将对材料行为有影响。
节理平面
使用 JSET 命令，用关键字DIP（倾角）和ORIGIN（平面上一点）指定方向和位置，
可以创建单个或一组节理平面。设置节理平面后，会自动识别出落在节理平面上的接触（组
成该接触的两个BALL位于节理平面两侧），赋予节理号，并分配不同于其他接触的接触模
型和特性。节理号可以用PRINT CONTACT命令显示。
JSET 命令只能应用在颗粒组达到初始平衡状态后，且只有法向力为非零或具有粘结的
接触才可以转变成节理平面接触。可以在备样时采用消除漂浮颗粒程序将这些漂浮颗粒消
除。
创建单个连续节理平面的命令：
JSET id=1 dip=-45.0 dd=20.0 origin=(0.0,0.0,2.0)
创建一组连续节理平面的命令：
jset id=2 dip=0 origin=(3.0,2.0) number=2 spacing=5.0
也可以创建由一组线段组成的平面，使用关键字RADIUS 和AREA_RATIO。RADIUS 为线段长
度的一半，AREA_RATIO 为线段总长与区域总长度的比值，但没有处理重叠的线段，所以
AREA_RATIO只是节理平面百分率的近似值。
关键字DIP、SPACING、AREA_RATIO、RADIUS可以使用均值和方差，默认为服从均匀分布，
如果使用关键字GAUSS则服从高斯正态分布。如下：
jset id 1 dip 50,2 number 5 spacing 50,3，其中spacing为沿着平均平面方向的
间距。、
为节理面接触分配接触模型和特性：
property n bond= 0 s bond= 0 range jset 1
property fric 0.25 range jset 1,2
注意：用JSET创建的一个节理平面并分配较弱的特性一般并不能使材料沿该平面滑动，因为该
平面是有粗糙度的。如果要指定一个滑面，需要自定义FISH函数创建一个至少三个颗粒宽度的
弱性质条带。
使用命令 PROPERTY 可以分配材料特性参数信息。难点在于选择能反映真实材料的特
性参数。连续介质程序中的材料参数是宏观参数，可以由实验室得到，而PFC2D 则是采用
颗粒之间的微观参数。如果材料的微观参数已知，则可以直接应用PFC2D 模拟；如果微观
参数未知，则需要采用逆向模拟，即进行许多试验和模拟，建立宏观参数和微观参数的关系。
当不连续体（如裂缝、断层、节理等）要设置在材料中进行模拟时，可以使用节理生成器引
入弱面。
颗粒和接触的特性
除了与粘结相关的特性以外，都要使用PROPERTY 命令给颗粒分配特性参数。有些参
数，如摩擦系数和接触刚度，与接触行为相关，但是也必须分配给颗粒。这些接触行为由组
成接触的两个颗粒的特性得到。
特性在空间的变化
参数沿空间线性变化：prop s_bond 1e6 grad -1e4 range x 0 50。即，接触粘结剪切强度由
x=0处的1e6 线性变化到x=50处的5e5（10e6 − 1e4*50）。
参数按照更复杂的函数变化则需要通过FISH 函数完成。
参数按组分配：使用命令GROUP和RANGE GROUP。GROUP将某空间范围内的颗粒
定义为一组后，即使以后这些颗粒移动到初始定义的RANGE 空间范围外，仍然保持在一个
组内，仍然可以对其同时分配或改变参数。需要注意的是，通过GROUP设置接触特性时，
只有GROUP内部球体间形成的接触可以被改变，而属于不同GROUP的颗粒间形成的接触
不受影响。如果需要改变不同GROUP之间的颗粒形成的接触的特性，则需要使用JSET命
令，这将在后面的节理平面部分给出详细说明。
基于已知微观特性的直接模拟
如果模拟的材料由圆颗粒组成，则颗粒的特性可以直接输入PFC2D。各参数的相对重
要性取决于模拟的类型和目的。如果颗粒材料小应变时的整体模量很重要，则应该使用
HERTZ-MINDLIN 模型，输入剪切模量和孔隙比；如果材料的强度比模量更重要，则使用
计算效率高的线性接触模型，这种情况下，如果颗粒的重叠在合理范围（比如小于平均半径
的5%），则切向和法向接触刚度就不重要了；如果没用颗粒间粘结，则材料的强度决定于
摩擦系数。
对于矩形边界的颗粒组，可以建立模量和强度的解析表达式。在解析表达式中，存在
一些普遍的关系：第一，二维颗粒组的整体模量与颗粒刚度成正比，但与颗粒半径无关；第
二，对于使用接触粘结模型的颗粒组，粘结对整体强度的贡献（区别于摩擦的贡献，即抗拉
强度）与平均颗粒半径成反比。平行粘结对整体行为的影响更加复杂。
未知微观特性时的逆向模拟
先做一系列物理试验，得到材料的各种宏观参数，然后不断调整数值试验的输入参数，
使数值模拟结果与物理试验结果相吻合，然后建立物理试验参数与数值试验参数之间的关
系，用于数值方法模拟更多的大型问题。这个过程有点想是撞运气，但是我们可以根据经验
缩小参数的范围。
双轴试验
任何合成材料的整体特性可以通过对颗粒材料进行一系列试验获得。这些试验以数值形
势进行并模拟物理试验。
双轴试验PFC模型
合成材料试样由一组圆颗粒组成，其局限性主要有两点：首先，PFC2D 模型只发生平
面应力或平面应变，这与真实材料的三轴试样中的响应不同；其次，为了计算应力，假定颗
粒都是单位厚度。
该双轴试验的模拟中，用四个墙限制矩形试样，顶墙和底墙模拟加荷板，左墙和右墙
模拟试样侧限。试样应变控制形式加荷，指定顶板和底板的速度。在试验的所有阶段，左右
墙的速度都由数值伺服机制自动控制，以使侧向应力恒定。应力和应变的计算方法是，将所
有作用在相应墙上的力相加，求得宏观上的平均应力和应变。材料的各应力和应变的响应，
可以使用HISTORY 来进行跟踪观察。
试样制备
方法详见《PFC2D 学习笔记之颗粒生成细节》
计算并控制应力状态
应力和应变状态的有FISH 函数GET_SS确定。应力是作用在每组相对的墙上的应力平
均值，而作用在每个墙上的应力则是用作用在该墙上的总的力除以试样长度。应变的计算公
式为2*（L-L0）/(L0+L)。get_ss的核心程序如下：
wadd1 = find_wall(1)
wadd2 = find_wall(2)
wadd3 = find_wall(3)
wadd4 = find_wall(4)
xdif = w_x(wadd2) - w_x(wadd4)
ydif = w_y(wadd3) - w_y(wadd1)
new_xwidth = width + xdif
new_height = height + ydif
wsxx = 0.5 * (w_xfob(wadd4) - w_xfob(wadd2)) / (new_height * 1.0)
wsyy = 0.5 * (w_yfob(wadd1) - w_yfob(wadd3)) / (new_xwidth * 1.0)
wexx = 2.0 * xdif / (width + new_xwidth)
weyy = 2.0 * ydif / (height + new_height)
wevol = wexx + weyy
加载过程中，要使用数值伺服机制，通过调整侧墙的速度使侧向应力保持不变，该伺服机
制通过FISH函数SERVO和GET_GAIN实现。SERVO每个CYCLE中调用一次，它通过GET_SS
函数得到应力，并使用伺服机制调节墙的速度以减小测量应力与目标应力之差。变量Y_SERVO
作为开关，若为零，则伺服机制不作用在顶板和底板。GET_gain函数用于获得GAIN参数”G”。
GET_GAIN的核心代码如下：
alpha = 0.5 ; relaxation factor
count = 0
avg_stiff = 0
cp = contact_head ; find avg. number of contacts on x-walls
loop while cp # null
if c_ball1(cp) = wadd2
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball1(cp) = wadd4
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball2(cp) = wadd2
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball2(cp) = wadd4
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
cp = c_next(cp)
end_loop
nxcount = count / 2.0
avg_stiff = avg_stiff / count
gx = alpha * (height * 1.0) / (avg_stiff * nxcount * tdel)
count = 0
avg_stiff = 0
cp = contact_head ; find avg. number of contacts on y-walls
loop while cp # null
if c_ball1(cp) = wadd1
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball1(cp) = wadd3
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball2(cp) = wadd1
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
if c_ball2(cp) = wadd3
count = count + 1
avg_stiff = avg_stiff + c_kn(cp)
end_if
cp = c_next(cp)
end_loop
nycount = count / 2.0
avg_stiff = avg_stiff / count
gy = alpha * (width * 1.0)/ (avg_stiff * nycount * tdel)
SERVO的核心程序如下
while_stepping
get_ss ; compute stresses & strains
udx = gx * (wsxx - sxxreq)
w_xvel(wadd4) = udx
w_xvel(wadd2) = -udx
if y_servo = 1 ; switch stress servo on or off
udy = gy * (wsyy - syyreq)
w_yvel(wadd1) = udy
w_yvel(wadd3) = -udy
end_if
微观特性选择方法
线性接触模型
1、初始弹摸月接触刚度线性相关。平行粘结能增加接触刚度，但接触粘结不能。
2、泊松比取决于取决于压密状态和剪切接触刚度与法向接触刚度的比值。
3、峰值强度取决于摩擦系数和粘结强度。如果只有摩擦系数，材料表现出塑性或微弱软化
行为；如果粘结强度也加上，材料将表现出更大的峰值强度但是也有更大的应变软化。
4、当围压增加时，摩擦强度部分将提高，但是粘结强度部分不变。因此，高围压时材料将
表现出更大的塑性。
5、如果使用接触粘结模型，峰值后的卸载-加载模量仅比初始模量稍低；如果使用平行粘结
模型，则该模量会随着应变的增加而减小；随着平行粘结破坏时，材料不断积累损伤。
6、如果平行粘结刚度与接触刚度的比值增大，损伤累积的速率随着应变的发展儿增大。
7、如果强度值围绕平均值分布而不是只设为单一值，峰值的形状会更加平滑和宽大。
8、试样孔隙比越低，则在峰值后的体积膨胀越大。
9、粘结可以在任何平均应力值处设置；在施加剪切应力前，该应力将增加或降低。在平均
应力或偏应力处设置粘结的作用是，为了锁闭接触力和内能。这将对材料行为有影响。
节理平面
使用 JSET 命令，用关键字DIP（倾角）和ORIGIN（平面上一点）指定方向和位置，
可以创建单个或一组节理平面。设置节理平面后，会自动识别出落在节理平面上的接触（组
成该接触的两个BALL位于节理平面两侧），赋予节理号，并分配不同于其他接触的接触模
型和特性。节理号可以用PRINT CONTACT命令显示。
JSET 命令只能应用在颗粒组达到初始平衡状态后，且只有法向力为非零或具有粘结的
接触才可以转变成节理平面接触。可以在备样时采用消除漂浮颗粒程序将这些漂浮颗粒消
除。
创建单个连续节理平面的命令：
JSET id=1 dip=-45.0 dd=20.0 origin=(0.0,0.0,2.0)
创建一组连续节理平面的命令：
jset id=2 dip=0 origin=(3.0,2.0) number=2 spacing=5.0
也可以创建由一组线段组成的平面，使用关键字RADIUS 和AREA_RATIO。RADIUS 为线段长
度的一半，AREA_RATIO 为线段总长与区域总长度的比值，但没有处理重叠的线段，所以
AREA_RATIO只是节理平面百分率的近似值。
关键字DIP、SPACING、AREA_RATIO、RADIUS可以使用均值和方差，默认为服从均匀分布，
如果使用关键字GAUSS则服从高斯正态分布。如下：
jset id 1 dip 50,2 number 5 spacing 50,3，其中spacing为沿着平均平面方向的
间距。、
为节理面接触分配接触模型和特性：
property n bond= 0 s bond= 0 range jset 1
property fric 0.25 range jset 1,2
注意：用JSET创建的一个节理平面并分配较弱的特性一般并不能使材料沿该平面滑动，因为该
平面是有粗糙度的。如果要指定一个滑面，需要自定义FISH函数创建一个至少三个颗粒宽度的
弱性质条带。
» 求助~询问。。PFC2D 求助ing
» 求助~询问。。PFC2D
» PFC是怎么用颗粒模拟岩体的
» 希望大家踊跃发言，积极讨论PFC 技术与问题~
» PFC2D学习笔记之流固耦合
» PFC2D学习笔记之使用细则
» PFC学习初期总结
» PFC2D学习笔记之颗粒生成