# Create a pack in which the cylindrical particles and spherical particles are randomly distributed.

Asked by Pengfei Tan on 2019-12-05

Hello,

I would like to simulate a process in which a roller is used to spread cylindrical particles and spherical particles. So, first I need to create a pack which includes the random distribution of the cylindrical particles and spherical particles. I know the cylindrical particles and spherical particles can be generated using commands "cylinder()" and "pack.SpherePack.makeCloud()", respectively.

Could you please give me some suggestions about how to create a pack in which the cylindrical particles and spherical particles are randomly distributed?

## Question information

Language:
English Edit question
Status:
Solved
For:
Assignee:
No assignee Edit question
Solved by:
Pengfei Tan
Solved:
2019-12-06
Last query:
2019-12-06
2019-12-05
 Jan Stránský (honzik) said on 2019-12-05: #1

Hello,

a quick solution would be:
1) create ordinary sphere packing using makeCloud (possible called multiple times with different parameters for spheres and spheres to be replaced by cylinders)
2) replace respective spheres by cylinders
3) (if needed) compress the packing

For more detailed answer, please provide more detailed description of the problem:
- if/what you have already tried
- definition of "random distribution" (positions? orientation? volume fraction? ...)
- ...

cheers
Jan

 Pengfei Tan (tpf516) said on 2019-12-05: #2

(1) I already created the cylindrical particles. My code is as follows:

##############################################

# encoding: utf-8
"An example showing how to create cylinders with random length."

import gts, os.path, locale

locale.setlocale(locale.LC_ALL, 'en_US.UTF-8') # Note: gts is locale-dependent. If, for example, german locale is used, gts.read()-function does not import floats normally

#### Parameters ####
Lmin=3. # minimum length
Lmax=6. # maximum length
n=10 # number of cylinders in one direction
r=0.5 # radius of the cylinder element
phi=30. # friction angle
E=1e6 # Young's modulus
wire=False
fixed=False
z=-1.2
color=[0,0,1]
Ns=200
L=5 #length of the cylinder
dx=0.5
dy=0
nx=20
ny=1
z=2

#### Engines ####
O.engines=[
ForceResetter(),
InsertionSortCollider([
Bo1_GridConnection_Aabb(),
Bo1_PFacet_Aabb(),
Bo1_Sphere_Aabb(),
]),
InteractionLoop([
Ig2_GridNode_GridNode_GridNodeGeom6D(),
Ig2_Sphere_Sphere_ScGeom(),
Ig2_Sphere_GridConnection_ScGridCoGeom(),
Ig2_GridConnection_GridConnection_GridCoGridCoGeom(), # cylinder-cylinder interaction
Ig2_Sphere_PFacet_ScGridCoGeom(), # needed for GridNode-pFacet interaction (why is this not included in Ig2_GridConnection_PFacet_ScGeom???)
Ig2_GridConnection_PFacet_ScGeom(), # Cylinder-pFcet interaction
],
[
Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(setCohesionNow=True,setCohesionOnNewContacts=False),
Ip2_FrictMat_FrictMat_FrictPhys()
],
[
Law2_ScGeom_FrictPhys_CundallStrack(), # contact law for sphere-sphere
Law2_ScGeom6D_CohFrictPhys_CohesionMoment(), # contact law for "internal" cylider forces
Law2_ScGridCoGeom_FrictPhys_CundallStrack(), # contact law for Cylinder-pFacet interaction
Law2_GridCoGridCoGeom_FrictPhys_CundallStrack() # contact law for cylinder-cylinder interaction
]
),
GlobalStiffnessTimeStepper(timestepSafetyCoefficient=0.1,label='ts'),
NewtonIntegrator(gravity=(0.,0,-10),damping=0.5,label='newton'),
RotationEngine(angularVelocity=1,rotationAxis=(0,1,0),ids=[0,1,2])+TranslationEngine(translationAxis=[1,0,0],velocity=0.5,ids=[0,1,2]),
#VTKRecorder(fileName='3d-vtk-',recorders=['all'],iterPeriod=5000),

]
O.dt=0.1*PWaveTimeStep()

#### Creat materials ####
O.materials.append( CohFrictMat( young=E,poisson=0.3,density=1000,frictionAngle=radians(phi),normalCohesion=1e10,shearCohesion=1e10,momentRotationLaw=True,label='cMat' ) ) #material to create the gridConnections
O.materials.append( FrictMat( young=E,poisson=0.3,density=1000,frictionAngle=radians(phi),label='fMat' ) ) # material for general interactions

#cylindrical roller: the IDs of cylinder must be [0, 1, 2] which are used for translation and rotation engines.
nodesIds=[]
cylIds=[]

#### Spheres ####
O.materials.append(FrictMat(young=4.0e6,poisson=0.5,frictionAngle=frictionAngle1,density=1600,label='spheremat'))
sp=pack.SpherePack()
Ns=sp.makeCloud(Vector3(-5,-5,3.0),Vector3(+5,+5,+6.0),0.3,0.2,Ns,False,0.8)

nodesIds1=[]
cylIds1=[]
#### Create cylinders ####
for i in range(0,nx):
x=-5+i*dx
for j in range(0,ny):
y=-2.5+j*dy
color=[1,0,0],fixed=False,intMaterial='cMat',extMaterial='fMat')

# collect history of data which will be plotted
# each item is given a names, by which it can be the unsed in plot.plots
# the **O.energy converts dictionary-like O.energy to plot.addData arguments
plot.saveDataTxt('bbb.txt.bz2')

####create Pfacet from gts file####

#### Creat ground with pFacets ####
color=[255./255.,102./255.,0./255.]
n0=O.bodies.append( gridNode([-10,-10,0],r,wire=False,fixed=True,material='cMat',color=color) )
n1=O.bodies.append( gridNode([10,-10,0],r,wire=False,fixed=True,material='cMat',color=color) )
n2=O.bodies.append( gridNode([10,10,0],r,wire=False,fixed=True,material='cMat',color=color) )
n3=O.bodies.append( gridNode([-10,10,0],r,wire=False,fixed=True,material='cMat',color=color) )
O.bodies.append( gridConnection(n0,n1,r,color=color,material='fMat') )
O.bodies.append( gridConnection(n1,n2,r,color=color,material='fMat') )
O.bodies.append( gridConnection(n2,n0,r,color=color,material='fMat') )
O.bodies.append( gridConnection(n2,n3,r,color=color,material='fMat') )
O.bodies.append( gridConnection(n3,n0,r,color=color,material='fMat') )
O.bodies.append( pfacet(n0,n1,n2,wire=False,material='fMat',color=color) )
O.bodies.append( pfacet(n0,n2,n3,wire=False,material='fMat',color=color) )

#### For viewing ####
qt.View()
Gl1_Sphere.stripes=True

#### Allows to reload the simulation ####
O.saveTmp()

#####################################################

(2). "random distribution" means the positions and orientations.

Because my spherical particles and cylindrical particles have total different shapes, how to replace the spheres with cylinders?

Thank you a lot!

 Jan Stránský (honzik) said on 2019-12-05: #3

a MWE:
#########
def randomOrientation(): # most likely not very well random, but for this example it is ok
from random import random as r
return Quaternion((r(),r(),r()),2*pi*r())

# create loose spherical packing using makeCloud
mi,ma = (0,0,0),(20,20,20)
nCyls,nSphs = 40,30
sp.makeCloud(mi,ma,rMean=2,rRelFuzz=.5,num=nCyls) # makeCloud for cylinders
sp.makeCloud(mi,ma,rMean=1,rRelFuzz=.5,num=nSphs) # makeCloud for spheres
if i < nCyls: # add cylinder
ori = randomOrientation()
O.bodies.append(box(pos,extents,ori))