7.4. 混合颗粒堆积

接着上一个香蕉颗粒堆积算例 颗粒堆积 ， 我们加入梨子颗粒混在一起进行混合堆积 (图 7.5)，即两种不同形状颗粒的堆积。因此，我们需要引入两个形状模板，并分别构建颗粒。以下为具体的建模脚本。

图 7.5 香蕉梨子混合堆积

7.4.1. 基本堆积脚本

以下脚本和上一个例子保持一致。我们首先在格子点上生成颗粒的位置信息。

import numpy as np
import sys, random
randSeed = 12345654 #to keep the same config.
random.seed(randSeed)
gui = True
sim = None
from dem.core import *
from dem.rtdem import *
from pySudoMath import *
if gui:
    app.core.setAnalysisType(AnalysisType="RTDEM", randSeed = randSeed)
    sim=app.core.sim
else:
    sim = DEMSim()

#define a lattice grid
#R: distance between two neighboring nodes
#num_x: number of nodes along x axis
#num_y: number of nodes along y axis
#num_z: number of nodes along z axis
#return a list of the postions of all nodes
scales = list()
def GridInitial(R,num_x=10,num_y=10,num_z=20, shift=[0,0,0]):#this part is slow, but just for test
    pos = list()
    for k in range(num_z):
        for j in range(num_y):
            for i in range(num_x):
                x = i*R*2.0*1.2+R+shift[0]+random.uniform(0, 0.1*R)
                y = j*R*2.0*1.2+R+shift[1]+random.uniform(0, 0.1*R)
                z = k*R*2.0*1.2+R+shift[2]
                pos.append(Vector3f(x,y,z))
                scales.append(random.uniform(0.6,1.2))
    return pos

psize = 0.09
Num = 5
Num_z = 20
pos = GridInitial(psize,num_x=Num,num_y=Num,num_z=Num_z,shift = [0.0*psize,0.0*psize,3.0*psize])
pnum = len(pos)
ori = list()
for i in range(pnum):
    #ang = random.uniform(0,np.pi)
    #vec = Vector3f()
    ori.append(Quaternionr(-3.1415926/4, Vector3f(0, 1.0, 1.0)))

7.4.2. 生成颗粒和容器

首先，我们随机打乱颗粒位置信息（在数组pos中的位置）；然后选择前一部分作为香蕉颗粒的位置信息，后一部分作为梨子颗粒的位置信息。

• 生成颗粒

  # 随机打乱位置信息
  pos = random.sample(pos,len(pos))
  
  # 生成香蕉颗粒
  poly=PolyhedronData()
  poly.setPosOri( pos[:400],  ori[:400])
  poly.loadSTL("./Banana-lowResol.stl", 0.002)
  poly.computGeometry()
  span = poly.bounds.span()
  print(span)
  #poly.setScaleFactors(np.array( scales ))
  sim.pdata.addObject(poly)
  
  # 生成梨子颗粒
  poly2=PolyhedronData()
  poly2.setPosOri( pos[400:],  ori[400:])
  poly2.loadSTL("./Pear.stl", 0.02)
  poly2.computGeometry()
  span = poly2.bounds.span()
  print(span)
  #poly.setScaleFactors(np.array( scales ))
  sim.pdata.addObject(poly2)
  

• 生成容器

  # 生成容器
  box=PolyhedronData()
  box.loadSTL("./box.stl", Vector3f(0.6,0.6,0.6), False)
  box.setPosition( [Vector3f(0.5,0.5,1.7)])
  box.setMovable(False)
  sim.pdata.addObject(box)
  

7.4.3. 其它模型设置

与上一个例子类似。

# 颗粒材料
mat = Material()
mat.density = 2650.0
mat.kn = 1e4*1e2
mat.ks = 1e4
mat.friction = 0.3
sim.pdata.material = mat
# 容器材料
wmat = Material()
wmat.density = 1000.0
wmat.kn = 1e5*1e3
wmat.ks = 1e5
wmat.friction = 0.3
sim.pdata.wallMaterial = wmat

sim.integrator.gravity = Vector3f(0,0,-10.0)
sim.integrator.damping = 0.3
sim.integrator.viscousDamping = 0.01


sim.init()

sim.dt = 1e-4