纯砂土的研究目前基本上趋于饱和了，目前离散元在散体材料部分的研究大都是还原砾砂石材料的形状。在此之前，离散元开发了很多考虑抗转动的模型，来表征材料的球度对其转动特性的影响。随着计算能力的提高以及研究的深化，直接生成现实形状的材料是当前主要的工作。于是PFC6.0出现了块体计算元素，本文主要讲解clump来模拟块体，而ball模拟砂土的分层压缩法。只将体积大的用clump也可以有效率的进行计算。

    一般来说我们生成的clump的数目不会特别多，这样就对其均匀性要求比较高了。和砂土不一样，砂土的均匀性体现在孔隙率的均匀上，而块石的均匀性体现在分布上。如果不分层的话，很容易出现某些地方clump很多，而有些地方clump很少，这样必然对破坏模式产生影响。
 

    这部分实现原理也比较简单，只要在生成颗粒的时候进行分流就可以了。先看一下本文的四种粒径：
 

def par    width=0.15    height= width*2         y_vel=10.0        cengshu=5.0    pengzhangxishu=5.0    tbjipei=table.create("jipei")    table(tbjipei,0.003)=0.3    table(tbjipei,0.004)=0.5    table(tbjipei,0.01)=0.7    table(tbjipei,0.015)=1        poro=0.12    ;dengchaxian=0.02end

    其中0.003和0.004粒径的我们选择用Ball去模拟，和0.01和0.015粒径的我们使用clump去模拟。
 

首先我们需要准备好形状文件并导入成clump模板。
 

geometry import suishi.dxfclump template create name suishi ...                         geometry suishi ...                            bubblepack ratio 0.3 distance 130 ...                            surfcalculate 

之后在之前分层压缩代码的基础上进行分流：
 

if keligid=3 then                 command                    clump generate azimuth 0 360 templates 1 suishi size [table.x(tbjipei,keligid)] diameter number ...                    [math.ceiling(kelinum)] box [-width*0.5] [width*0.5] ...                    [dibuzuobiao] [pengzhangGaoDu] group @cenggName slot 1 group @ballgroup slot 2 tries 20000                                    clump attribute density 2.7e3 damp 0.7                endcommand            else if keligid=4 then                 command                    clump generate azimuth 0 360 templates 1 suishi size [table.x(tbjipei,keligid)] diameter number ...                    [math.ceiling(kelinum)] box [-width*0.5] [width*0.5] ...                    [dibuzuobiao] [pengzhangGaoDu] group @cenggName slot 1 group @ballgroup slot 2 tries 20000                                   clump attribute density 2.7e3 damp 0.7                endcommand            else                command                 ball generate radius [table.x(tbjipei,keligid)*0.5] number [math.ceiling(kelinum)] box [-width*0.5] [width*0.5] ...                    [dibuzuobiao] [pengzhangGaoDu] group @cenggName slot 1 group @ballgroup slot 2 tries 20000                ball attribute density 2.7e3 damp 0.7                endcommand            endif

    这里使用clump generate进行clump的生成，注意使用azimuth指定clump的旋转角度，这样生成的clump就不会都是一样的方向。

    我们看一下最后的结果：
 

其实用surface来显示clump会好些：
 

当然四种粒径的颗粒也可以看出来：
 

完整代码如下，用户需自己提供suishi.dxf的形状文件。

new def par    width=0.15    height= width*2         y_vel=10.0        cengshu=5.0    pengzhangxishu=5.0    tbjipei=table.create("jipei")    table(tbjipei,0.003)=0.3    table(tbjipei,0.004)=0.5    table(tbjipei,0.01)=0.7    table(tbjipei,0.015)=1        poro=0.12    ;dengchaxian=0.02end@pardomain extent [-height*10] [height*10]set random 10001wall generate box [-width*0.5] [width*0.5] [-height*0.5] [height*(pengzhangxishu-1)-height*0.5] expand 1.5cmat default model linear method deformability emod 100e6 kratio 1.5 geometry import suishi.dxfclump template create name suishi ...                         geometry suishi ...                            bubblepack ratio 0.3 distance 130 ...                            surfcalculate def create_sample    tb_ceng=table.create("ceng_poro")    loop n(1,cengshu)        dibuzuobiao=-height*0.5+height*(n-1)/cengshu        dingbuzuobiao=-height*0.5+height*(n)/cengshu                pengzhangGaoDu=dingbuzuobiao+height*(pengzhangxishu-1)/cengshu           command            wall delete walls range id 3            wall create id 3 vertices [-width*0.5*1.5] [pengzhangGaoDu] [width*0.5*1.5] [pengzhangGaoDu]         endcommand        poroceng=poro;-(n-(cengshu+1)*0.5)*dengchaxian/(cengshu-1)        table(tb_ceng,n)=poroceng                 keligeoup=table.size(tbjipei)        Vs=height*width*(1-poroceng)/cengshu                cenggName=string.build("ceng%1",n)        loop m(1,keligeoup)            keligid=keligeoup+1-m            vslocal=0            if keligid=1 then                vslocal=Vs*table.y(tbjipei,keligid)            else                            vslocal=Vs*(table.y(tbjipei,keligid)-table.y(tbjipei,keligid-1))            endif            kelinum=vslocal/(math.pi*(table.x(tbjipei,keligid)*0.5)^2)            ballgroup=string.build("ball%1",keligid)                        if keligid=3 then                 command                    clump generate azimuth 0 360 templates 1 suishi size [table.x(tbjipei,keligid)] diameter number ...                    [math.ceiling(kelinum)] box [-width*0.5] [width*0.5] ...                    [dibuzuobiao] [pengzhangGaoDu] group @cenggName slot 1 group @ballgroup slot 2 tries 20000                                    clump attribute density 2.7e3 damp 0.7                endcommand            else if keligid=4 then                 command                    clump generate azimuth 0 360 templates 1 suishi size [table.x(tbjipei,keligid)] diameter number ...                    [math.ceiling(kelinum)] box [-width*0.5] [width*0.5] ...                    [dibuzuobiao] [pengzhangGaoDu] group @cenggName slot 1 group @ballgroup slot 2 tries 20000                                   clump attribute density 2.7e3 damp 0.7                endcommand            else                command                 ball generate radius [table.x(tbjipei,keligid)*0.5] number [math.ceiling(kelinum)] box [-width*0.5] [width*0.5] ...                    [dibuzuobiao] [pengzhangGaoDu] group @cenggName slot 1 group @ballgroup slot 2 tries 20000                ball attribute density 2.7e3 damp 0.7                endcommand            endif        endloop                command            wall attribute yvel [-y_vel] range id 3            solve time [(pengzhangGaoDu-dingbuzuobiao)/y_vel]                        wall attribute yvel 0 range id 3            cycle 2000        endcommand            endloopend@create_samplesolvesave sample