Yade

Size problem in yade simulation

Asked by Huang peilun

I did a gravity deposition test in a cylinder. The following is the original code and it works.

************************************************
from yade import pack

#material
SoilMat=FrictMat(young=25e6,poisson=0.2,density=2650,frictionAngle=0.28)
O.materials.append((SoilMat))

#define the chamber
O.bodies.append(geom.facetCylinder(center=(.5,.5,.5),radius=.5,height=1,wallMask=6,color=(1,0,0)))

#create foundation by making spheres in the cylinder.
sp=pack.SpherePack()
sp.makeCloud((0,0,0),(1,1,1),rMean=.05,rRelFuzz=.3)
cyl=pack.inCylinder((.5,.5,0),(.5,.5,1),.45)
sp=pack.filterSpherePack(cyl,sp,True)
sp.toSimulation(material=SoilMat)

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_FrictPhys()],
  [Law2_ScGeom_FrictPhys_CundallStrack()]
 ),
 # update position using Newton's equations
 NewtonIntegrator(gravity=(0,0,-9.81),damping=0.05),
 # call the checkUnbalanced function (defined below) every 2 seconds
 PyRunner(command='checkUnbalanced()',realPeriod=2,label='checker'),
]
O.dt=.2*PWaveTimeStep()

def checkUnbalanced():
 print('%s: Unbalanced=%s' % (O.iter,unbalancedForce()))
 if O.iter<5000: return
 if unbalancedForce()<.05:
  O.pause()

O.saveTmp()
************************************************

However, when I magnify all units 1000 times, like the following code, some sphere just went through the bottom of the cylinder and I don't know why.

************************************************
from yade import pack

#material
SoilMat=FrictMat(young=25e6,poisson=0.2,density=2650,frictionAngle=0.28)
O.materials.append((SoilMat))

#define the chamber
O.bodies.append(geom.facetCylinder(center=(500,500,500),radius=500,height=1000,wallMask=6,color=(1,0,0)))

#create foundation by making spheres in the cylinder.
sp=pack.SpherePack()
sp.makeCloud((0,0,0),(1000,1000,1000),rMean=50,rRelFuzz=.3)
cyl=pack.inCylinder((500,500,0),(500,500,1000),450)
sp=pack.filterSpherePack(cyl,sp,True)
sp.toSimulation(material=SoilMat)

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_FrictPhys()],
  [Law2_ScGeom_FrictPhys_CundallStrack()]
 ),
 # update position using Newton's equations
 NewtonIntegrator(gravity=(0,0,-9.81),damping=0.05),
 # call the checkUnbalanced function (defined below) every 2 seconds
 PyRunner(command='checkUnbalanced()',realPeriod=2,label='checker'),
]
O.dt=.2*PWaveTimeStep()

def checkUnbalanced():
 print('%s: Unbalanced=%s' % (O.iter,unbalancedForce()))
 if O.iter<5000: return
 if unbalancedForce()<.05:
  O.pause()

O.saveTmp()
************************************************

Thanks in advance!

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Jan Stránský
Solved:
Last query:
Last reply:
Revision history for this message
Robert Caulk (rcaulk) said : 		#1

Have you tried increasing the stiffness of the facets relative to the particles?

Revision history for this message
Huang peilun (hpl16) said : 		#2

Thanks Robert.
I haven't tried increasing the stiffness of the facets because I don't know how. Can you show me how to do that?

Revision history for this message
Robert Caulk (rcaulk) said : 		#3

For example [1].

I will also give you the answer:

O.materials.append(FrictMat(young=desiredStiffness,label='walls'))
O.bodies.append(geom.facetCylinder(center=(.5,.5,.5),radius=.5,height=1,wallMask=6,color=(1,0,0)),material='walls')

However, I still recommend you read [2] carefully to ensure functionality.

Cheers.

Robert

[1]https://gitlab.com/yade-dev/trunk/blob/master/examples/triax-tutorial/script-session1.py#L54
[2]https://yade-dev.gitlab.io/trunk/tutorial-hands-on.html#singles

Revision history for this message
Huang peilun (hpl16) said : 		#4

Thanks!
I increased the stiffness but there still are spheres that go through the bottom of the cylinder. Here's my code.

********************************
from yade import pack

Mat=FrictMat(young=2.06e20)
SoilMat=FrictMat(young=25e6,poisson=0.2,density=2650,frictionAngle=0.28)
O.materials.append((SoilMat,Mat))

O.bodies.append(geom.facetCylinder(center=(500,500,500),radius=500,height=1000,wallMask=6,color=(1,0,0),material=Mat))

sp=pack.SpherePack()
sp.makeCloud((0,0,0),(1000,1000,1000),rMean=50,rRelFuzz=.3)
cyl=pack.inCylinder((500,500,0),(500,500,1000),450)
sp=pack.filterSpherePack(cyl,sp,True)
sp.toSimulation(material=SoilMat)

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_FrictPhys()],
  [Law2_ScGeom_FrictPhys_CundallStrack()]
 ),
 # update position using Newton's equations
 NewtonIntegrator(gravity=(0,0,-9.81),damping=0.05),
 # call the checkUnbalanced function (defined below) every 2 seconds
 PyRunner(command='checkUnbalanced()',realPeriod=2,label='checker'),
]
O.dt=.2*PWaveTimeStep()

def checkUnbalanced():
 print('%s: Unbalanced=%s' % (O.iter,unbalancedForce()))
 if O.iter<5000: return
 if unbalancedForce()<.05:
  O.pause()

O.saveTmp()

Revision history for this message
Best Jan Stránský (honzik) said : 		#5

Hello,

welcome to the playing with physical dimensions and units :-)
I assume you want in both cases "same" results:
- you increase size by a factor f (here f=1000)
- cross sections A increases by factor f^2
- volume and mass and gravity force increases with factor f^3
- penetration depth (pd) increases with factor f
- you left stiffness E unchanged, increasing by factor 1
- repulsive force F=pd*E*A/L increases with f*1*f^2/f = f^2
now it becomes 1000x "smaller" then in the original case! and is not large enough to prevent gravity force pushing particles through the cylinder..
So to "scale" repulsive force the same way as gravity force, you have to scale the stiffness by the factor f, too. Then F=pd*E*A/L=f*f*f^2/f=f^3.

Using
SoilMat=FrictMat(young=25e9,...)
in the second script, I got very similar results as in the first one

cheers
Jan

Revision history for this message
Huang peilun (hpl16) said : 		#6

Thanks Jan!
I tried SoilMat=FrictMat(young=25e9,...) and it works!
Now I understand why these spheres go through the cylinder. It's very interesting!

Revision history for this message
Huang peilun (hpl16) said : 		#7

Thanks Jan Stránský, that solved my question.