PFC模拟随机块体

2年前浏览4077

这里随便写了一串fish，应该可以实现一般意义上的随机块体了。我认为的随机块体，一方面是块体控制点的随机性，一方面是控制点处角度的随机性。控制点的随机性有两个，一个是对于整个块体而言，控制点数目的随机性，一个是对于单个控制点而言，控制点对于原点的距离也有随机性。

先看一下控制点数目的随机函数，这里使用的是高斯随机数math.random.gauss，这个随机数反应的为期望为0，标准差为1的高斯随机分布函数。根据概率论的知识，也比较容易在此基础上进行修正，得到期望为node_num_mean，标准差为node_num_sqrt的分布。为了避免误差，设置上下的阈值，分别为5和30。

控制点处的半径也是一样的道理

每个点处的角度也差不多，不过不同边数块体形状角度期望应该是不一样的，这里将360除以点数作为角度期望。

来看一下全部的代码：

newdomain extent -10 10set random 10001call geo_toolsdef rock_par    node_num_mean=10    node_num_sqrt=2     radius_sqrt=0.2        theta_sqrt=0.2end@rock_par
def get_node_num    temp_gauss=math.random.gauss     temp_num=math.floor(temp_gauss*node_num_sqrt node_num_mean)    if temp_num<5 then        temp_num=5    endif    if temp_num>30 then        temp_num=30    endif    get_node_num=temp_numend
def get_radius    temp_gauss=math.random.gauss     temp_radius=temp_gauss*radius_sqrt 1    if temp_radius<0.1 then        temp_radius=0.1    endif    get_radius=temp_radiusend
def get_theta(theta_local_mean,theta_local_sqrt)    temp_gauss=math.random.gauss     temp_theta=temp_gauss*theta_local_sqrt theta_local_mean    if temp_theta<0.1 then        temp_theta=0.1    endif    get_theta=temp_thetaend
def get_array_sum(array_input)    array_sum=0    array_num=array.size(array_input,1)    loop local array_count(1,array_num)        array_sum =array_input(array_count)    endloop    get_array_sum=array_sumenddef create_rock_geo(num)    loop rock_count(1,num)         geo_name=string.build("rock_%1",rock_count)        rock_geo_collect = geom.set.create(geo_name)         rock_local_poly = geom.poly.create(rock_geo_collect)         node_num=get_node_num        theta_mean=360/float(node_num)        theta_local=0        theta_info=array.create(node_num)        loop node_count(1,node_num-1)            radius_local=get_radius            node_pos_x=math.cos(theta_local*math.degrad)*radius_local            node_pos_y=math.sin(theta_local*math.degrad)*radius_local             theta_incri=get_theta(theta_mean,theta_mean*theta_sqrt)            theta_local =theta_incri            node_point=geom.node.create(rock_geo_collect,vector(node_pos_x,node_pos_y))             e_temp= geom.poly.add.node(rock_geo_collect,rock_local_poly,node_point)        endloop        eer=geom.poly.add.node(rock_geo_collect,rock_local_poly,geom.node.find(rock_geo_collect,1))        pos_count=rock_count-1        x_pos=0.15*pos_count/math.ceiling(math.sqrt(num))        y_pos=0.15*(pos_count-math.sqrt(num)*x_pos/0.15)        if geom.poly.check(rock_local_poly) then            command                clump template create bubblepack distance 80 ratio 0.1 name @geo_name ...                    geometry @geo_name surfcalculate                clump replicate x @x_pos y @y_pos  name @geo_name diameter   0.1            endcommand        endif    endloopend@create_rock_geo(64)
save clump_template

这里的控制点数目的期望为10，看一下效果：

总的来说还是满足我们需求的。

当然有时候也会出现误差，这种误差在生成clump template时候避免了：

比如设置为控制点数为20的时候，错误的形状不会生成template：

这些模板生成好了后，便可以使用generate或者distribute在模板中随机生成。

restore clump_templateclump delete
clump generate size 0.1 0.3 number 100 box -2 2 diameter

微信截图_20220720145846.png

代码&命令科普 PFC

著作权归作者所有，欢迎分享，未经许可，不得转载

首次发布时间：2022-07-20

PFC 入门到精通 81 讲：从 0 掌握离散元颗粒流软件应用与实践

¥89

1条评论

＠_＠

签名征集中

1年前

请问学长，调用的geo_tools文件是软件自带的吗，我应该怎样获取呢

回复 2条回复