Yade

Cylinder powder compaction

Asked by Mithushan Soundaranathan on 2021-01-21

Hi,

I am working on compaction of powder into cylinder shape. I have created a cylinder cloud with spherical particles, which I am trying to compress.
I tried to use PeriTriaxialController, I am getting error message saying that this will not work for aperiodic simulation.
So I tried to use TriaxialCompressionEngine, I am not getting any error message. But problem is that the particle will not compact together, the particles are staying at the same place in the cloud.

My question is, how can I compact particle in a cylinder cloud? Or is it possible to create a cylinder walls around the sphere and in the bottom and apply pressure from the top only (-z direction)?

My code is:

from __future__ import print_function
sigmaIso=-1e5

#import matplotlib
#matplotlib.use('Agg')

# generate loose packing
from yade import pack, qt, plot
from yade import utils

sample_material=CohFrictMat(
    young=4e9
   ,poisson=1
   ,density=1400
   ,frictionAngle=radians(30)
   ,isCohesive=True
   ,normalCohesion=1e8*1.2
   ,shearCohesion=.4e8*1.2
   ,momentRotationLaw=True
   ,label='sample_mat')

O.materials.append(sample_material)

#DATOS PROBETA
radius_tablet = 0.015 #1.5 cm
height_tablet = 0.02 #2 cm
radius_particle = 0.001 #1 mm
radius_std = 0.0005 #0.5 mm

# randomDensePack packing
#pred = pack.inCylinder((0,0,0),(0,0,height_tablet),radius=radius_tablet)
#assembly = pack.randomDensePack(pred,radius=radius_particle,rRelFuzz=0,returnSpherePack=True)
#assembly.toSimulation(material='sample_mat')

sp = pack.SpherePack()
sp.makeCloud((0,0,0),(2,2,2),rMean=.1,num=400)
cyl = pack.inCylinder((1,1,0),(1,1,2),1)
sp = pack.filterSpherePack(cyl,sp,True)
sp.toSimulation()

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_FrictPhys()],
  [Law2_ScGeom_FrictPhys_CundallStrack()]
 ),
 TriaxialCompressionEngine(label='triax',
  # specify target values and whether they are strains or stresses
  goal=(sigmaIso,sigmaIso,sigmaIso),stressMask=7,
  # type of servo-control
  dynCell=True,maxStrainRate=(10,10,10),
  # wait until the unbalanced force goes below this value
  maxUnbalanced=.1,relStressTol=1e-3,
  # call this function when goal is reached and the packing is stable
  doneHook='compactionFinished()'
 ),
 NewtonIntegrator(damping=.2),
 PyRunner(command='addPlotData()',iterPeriod=100),
]
O.dt=.5*PWaveTimeStep()

def addPlotData():
 plot.addData(unbalanced=unbalancedForce(),i=O.iter,
  sxx=triax.stress[0],syy=triax.stress[1],szz=triax.stress[2],
  exx=triax.strain[0],eyy=triax.strain[1],ezz=triax.strain[2],
  # save all available energy data
  Etot=O.energy.total(),**O.energy
 )

# enable energy tracking in the code
O.trackEnergy=True

# define what to plot
plot.plots={'i':('unbalanced',),'i ':('sxx','syy','szz'),' i':('exx','eyy','ezz'),
 # energy plot
 ' i ':(O.energy.keys,None,'Etot'),
}
# show the plot
plot.plot()

def compactionFinished():
 # set the current cell configuration to be the reference one
 O.cell.trsf=Matrix3.Identity
 # change control type: keep constant confinement in x,y, 20% compression in z
 triax.goal=(sigmaIso,sigmaIso,-.2)
 triax.stressMask=3
 # allow faster deformation along x,y to better maintain stresses
 triax.maxStrainRate=(1.,1.,.1)
 # next time, call triaxFinished instead of compactionFinished
 triax.doneHook='triaxFinished()'
 # do not wait for stabilization before calling triaxFinished
 triax.maxUnbalanced=10

def triaxFinished():
 print('Finished')
 O.pause()

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Mithushan Soundaranathan
Solved:
2021-02-04
Last query:
2021-02-04
Last reply:
2021-02-01

This question was reopened

Karol Brzezinski (kbrzezinski) said : #1

If you just need to prepare a cloud prior to your simulation you should check the manual [1] and this example [2].
If the compaction process is within the scope of your analysis, you can create a cylinder from facets, using this function [3].

[1] https://yade-dem.org/doc/user.html#sphere-packings
[2] https://gitlab.com/yade-dev/trunk/blob/master/examples/concrete/brazilian.py
[3] https://yade-dem.org/doc/yade.geom.html?highlight=facetcylinder#yade.geom.facetCylinder

Mithushan Soundaranathan (mithushan93) said : #2

Hi,

thank you very much, I have managed to fill a cylinder with a cloud. But have problem when I try to compress the spheres inside the cylinder only in the z-direction. I am not sure how to do it?

Best,
Mithu

Mithushan Soundaranathan (mithushan93) said : #3

Here is my code:

from __future__ import print_function
sigmaIso=-1e5

#import matplotlib
#matplotlib.use('Agg')

# generate loose packing
from yade import pack, qt, plot
from yade import utils

sample_material=CohFrictMat(
    young=4e9
   ,poisson=1
   ,density=1400
   ,frictionAngle=radians(30)
   ,isCohesive=True
   ,normalCohesion=1e8*1.2
   ,shearCohesion=.4e8*1.2
   ,momentRotationLaw=True
   ,label='sample_mat')

O.materials.append(sample_material)

#DATOS PROBETA
radius_tablet = 0.015 #1.5 cm
height_tablet = 0.02 #2 cm
radius_particle = 0.001 #1 mm
radius_std = 0.0005 #0.5 mm

# randomDensePack packing
#pred = pack.inCylinder((0,0,0),(0,0,height_tablet),radius=radius_tablet)
#assembly = pack.randomDensePack(pred,radius=radius_particle,rRelFuzz=0,returnSpherePack=True)
#assembly.toSimulation(material='sample_mat')

sp = pack.SpherePack()
sp.makeCloud((0,0,0),(2.1,2.1,2.1),rMean=.1,num=500)
#cyl = pack.inCylinder((1,1,0),(1,1,2),1)
#sp = pack.filterSpherePack(cyl,sp,True,material=sample_material)
sp.toSimulation(material=sample_material)

id_wall = O.materials.append(FrictMat(density=10000,young=1e11,poisson=.5,frictionAngle=atan(0.5),label="wallmat"))
O.materials[id_wall]

### WALLS

O.bodies.append(yade.geom.facetCylinder((1,1,1),radius=1.6,height=3,segmentsNumber=20,wallMask=6,material=id_wall))

O.bodies.append(yade.geom.facetBox((1,1,1),(3,3,3)))

Karol Brzezinski (kbrzezinski) said : #4

Hi,

the first thing is, you need to add "Bo1_Wall_Aabb()" and "Bo1_Facet_Aabb()" in the collider, so the balls could interact with walls and facets. The same story with "Ig2_Wall_Sphere_ScGeom(), Ig2_Facet_Sphere_ScGeom()" in the interaction loop.

I erased the box from your example, and added one wall at the top of your cloud. It moves down and pushes the balls. The rest is up to you. Please find working code below.

Cheers,
Karol

##############################
from __future__ import print_function
sigmaIso=-1e5

#import matplotlib
#matplotlib.use('Agg')

# generate loose packing
from yade import pack, qt, plot
from yade import utils

sample_material=CohFrictMat(
    young=4e9
   ,poisson=1
   ,density=1400
   ,frictionAngle=radians(30)
   ,isCohesive=True
   ,normalCohesion=1e8*1.2
   ,shearCohesion=.4e8*1.2
   ,momentRotationLaw=True
   ,label='sample_mat')

O.materials.append(sample_material)

#DATOS PROBETA
radius_tablet = 0.015 #1.5 cm
height_tablet = 0.02 #2 cm
radius_particle = 0.001 #1 mm
radius_std = 0.0005 #0.5 mm

# randomDensePack packing
#pred = pack.inCylinder((0,0,0),(0,0,height_tablet),radius=radius_tablet)
#assembly = pack.randomDensePack(pred,radius=radius_particle,rRelFuzz=0,returnSpherePack=True)
#assembly.toSimulation(material='sample_mat')

sp = pack.SpherePack()
sp.makeCloud((0,0,0),(2.1,2.1,2.1),rMean=.1,num=500)
#cyl = pack.inCylinder((1,1,0),(1,1,2),1)
#sp = pack.filterSpherePack(cyl,sp,True,material=sample_material)
sp.toSimulation(material=sample_material)

id_wall = O.materials.append(FrictMat(density=10000,young=1e11,poisson=.5,frictionAngle=atan(0.5),label="wallmat"))
O.materials[id_wall]

### WALLS

O.bodies.append(yade.geom.facetCylinder((1,1,1),radius=1.6,height=3,segmentsNumber=20,wallMask=6,material=id_wall))

#find max Z coordinate of your cloud
max_z = aabbExtrema()[1][2]
# add top wall
pushWallId = O.bodies.append(wall((0,0,max_z ), 2, material = O.materials[id_wall])) #

#set vertical downward speed to the wall
O.bodies[pushWallId].state.vel = (0,0,-1)

################### ENGINES
O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Wall_Aabb(),Bo1_Facet_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom(), Ig2_Facet_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_FrictPhys()],
  [Law2_ScGeom_FrictPhys_CundallStrack()]
 ),
 NewtonIntegrator(damping=.2),
]
O.dt=.5*PWaveTimeStep()

Mithushan Soundaranathan (mithushan93) said : #5

Hi Karol,

Thank you very much, that was exactly what I was looking for.

Best regards,
Mithu

Mithushan Soundaranathan (mithushan93) said : #6

Hi,

A quick question, the wall moving past the bottom of the cylinder. How do I stop it?

Best,
Mithu

Jan Stránský (honzik) said : #7

a quick answer, see [1] :-)
cheers
Jan

https://yade-dem.org/doc/user.html#stop-conditions

Mithushan Soundaranathan (mithushan93) said : #8

Hi Jan,

Thanks very much, is it possible to give a certain porosity as a stop condition?
How can I save the spheres position and radius after simulation?

Best,
Mithu

Jan Stránský (honzik) said : #9

> is it possible to give a certain porosity as a stop condition?

sure, see [1] :-) the checkUnbalanced() approach, just instead of unbalanced force you check porosity

> How can I save the spheres position and radius after simulation?

the easiest way is to use export.text [2]

cheers
Jan

[1] https://yade-dem.org/doc/user.html#stop-conditions
[2] https://yade-dem.org/doc/yade.export.html#yade.export.text

Mithushan Soundaranathan (mithushan93) said : #10

Hi Jan,

Thank you very much :-)

Best,
Mithu

Mithushan Soundaranathan (mithushan93) said : #11

Hi,

I was wondering if it possible to measure the volume of the particle packing at different time during simulation and also "count"the number of particles initially?

I tried to use the num function in makeCloud, but did get the desired number. F.eks when num=100 I only got around 50 particles, and for num=1000 I only got around 400 particles?

Best regards,
Mithushan

Karol Brzezinski (kbrzezinski) said : #12

If you do not trust makeCloud (I guess the number of the particles could be smaller if there is not enough space in predictor), why don't you count the spheres in the loop:

count = 0
for b in O.bodies:
    if isinstance(b.shape, Sphere):
        count +=1

Approximate size of your cloud in the current shape could be obtained using tesselation wrapper [1]:
TW=TesselationWrapper()
TW.triangulate()
TW.computeVolumes()
totVolume = 0
for i in range(count):
    totVolume += TW.volume(i)

Don't forget to triangulate and computeVolumes again to 'refresh' the state (when the shape of your cloud changes).

Karol

[1] https://yade-dem.org/doc/user.html?highlight=tesselation

Mithushan Soundaranathan (mithushan93) said : #13

Hi Karol,

Thank you very much, but problem is that I am measuring volume during compression and this method does account for particle deformation.

Best,
Mithu

Karol Brzezinski (kbrzezinski) said : #14

Hi Mithushan,

I think that I may have trouble understanding the problem.

As far as I know your sample code, you use 'regular' spheres, thus particles don't deform. The shape of your cloud deforms under compression. I expected that you would like to take this deformation into account, and that is why I proposed the tesselation wrapper. If you want to know the volume of the cloud before simulation, you can compute it from its dimensions e.g.
aabbDim()[0]*aabbDim()[1]*aabbDim()[2]

However, I don't think it is the problem, because the tesselation wrapper would give you the same result (before compression).

Could you try to rephrase the problem?

Karol

Mithushan Soundaranathan (mithushan93) said : #15

Hi Karol,

Thank you for your help.

Best,
Mithushan