0 引言
PFC的结果往往都是离散化的结果,我们对于结果大部分都是采用ball的一些显示,比如速度场、位移场等,或者是接触的一些东西。当然我是支持利用离散化的结果去进行描述的,这点也符合离散元的理论特性。
但是当某些情况下,颗粒的离散化比较大的时候,我们需要将一部分颗粒集 合体作为单元来显示,作为某一部分的结果。这时候就是我们测量圆的概念,这点在很多地方也会遇到,比如式样内部的应力分布、空隙分布等。这个是单个或者少部分颗粒无法描述的性状。
这时候或许用有限元那种单元的可视化效果会更加好一点。
1 案例
以一个土体自重沉降的案例为例:
成样:
model new
def par
width = 0.2
hight = width*2
rdmin = 0.009
rdmax = 0.006
poro = 0.4
end
@par
domain extent [-width*2.0] [width*2.0] ...
[-hight*2.0] [hight*2.0]
model random 10001
wall generate box [-width*0.5] [width*0.5] ...
[-hight*0.5] [hight*0.5] expand 1.5
cmat default model linear method deform emod 100e6 kratio 1.5
ball distribute porosity @poro ...
radius [rdmin] [rdmax] ...
box [-width*0.5] [width*0.5] ...
[-hight*0.5] [hight*0.5] ...
[-hight*0.5] [hight*0.5]
ball attribute density 2200 damp 0.7
model cycle 2000 calm 10
model cycle 1
model solve
ball delete range pos-x [-width] [-width*0.5]
ball delete range pos-x [width*0.5] [width]
ball delete range pos-y [-hight] [-hight*0.5]
ball delete range pos-y [hight*0.5] [hight]
ball delete range pos-z [-hight] [-hight*0.5]
ball delete range pos-z [hight*0.5] [hight]
model save "ball_sample"
沉降:
model restore "ball_sample"
model gravity 9.8
wall delete walls range id 2
model cycle 1
model solve
model save "zizhong"
沉降后的结果如图;
可以从力链上分析这个应力是上小下大的,是一个比较典型的力的分布。
3 绘制
这里使用的概念是,先把计算范围网格化,然后以网格中心为测中心,最小的长度作为半径生成测量圆。后面把测量圆的相关数据存到rblock的extra中,ball、clump和rblock这些基础元素都有这种extra属性用来存储用户自定义的数据。
model restore "zizhong"
def create_measure(x_min,x_max,y_min,y_max,z_min,z_max,n_x,n_y,n_z)
x_length=(x_max-x_min)/float(n_x)
y_length=(y_max-y_min)/float(n_y)
z_length=(z_max-z_min)/float(n_z)
rad=math.min(x_length,y_length)
rad=math.min(rad,z_length)
loop local n(1,n_x)
x_pos=x_min x_length*(n-1) x_length*0.5
loop local m(1,n_y)
y_pos=y_min y_length*(m-1) y_length*0.5
loop local k(1,n_z)
z_pos=z_min z_length*(k-1) z_length*0.5
command
measure create position [x_pos] [y_pos] [z_pos] radius [rad*0.5]
rblock create box [x_pos-x_length*0.5] [x_pos x_length*0.5] ...
[y_pos-y_length*0.5] [y_pos y_length*0.5] ...
[z_pos-z_length*0.5] [z_pos z_length*0.5]
endcommand
endloop
endloop
endloop
end
@create_measure([-width*0.5],[width*0.5],[-hight*0.5],[hight*0.5],[-hight*0.5],[hight*0.5],2,4,4)
def get_data
loop foreach mp measure.list
id=measure.id(mp)
rb=rblock.find(id)
rblock.extra(rb,1)=measure.stress.zz(mp)
endloop
end
@get_data
这里存储了z向的应力,测量圆分布如图:
把rblock的相关显示设置为:
与其相对应的rblock的效果为:
隐藏掉measure:
换种显示方式可以得到如下的结果:
配位数的分布如图:
本案例以方形的简单算例进行了计算,针对于复杂的模型,也可以先进行三角形或者多边形的rblock划分,然后再在rblock位置处生成测量圆的方式去做。理论应该是一致的,这里留给大家自己去拓展。