problems in making a 2D packing using TriaxialStressController

Asked by Chen Lifan on 2020-10-28

Hello, everyone,
I was making a 2D packing using TriaxialStressController() when I met the problem, i.e., the walls didn't move at all. And I wonder what caused it.

Here is my scripts:
from yade import pack

O.materials.append(CohFrictMat(young=88.e9,alphaKr=2.0,alphaKtw=0.,poisson=.7,frictionAngle=.0,isCohesive=True,normalCohesion=3.e9,shearCohesion=3.e9,momentRotationLaw=True,label='rock'))
O.materials.append(FrictMat(young=2e11,poisson=.4,frictionAngle=0,density=0,label='frictionlessWalls'))
sp = pack.SpherePack()
size = .3
sp.makeCloud(minCorner=(0,0,.05),maxCorner=(size,size,.05),rMean=.005,rRelFuzz=0.2,num=400,periodic=False,seed=1)
O.bodies.append([sphere(center,rad,material='rock') for center,rad in sp])
walls=aabbWalls(thickness=0,material='frictionlessWalls')
wallIds=O.bodies.append(walls)

print len(O.bodies)
for i in range(400):
   p=O.bodies[i]
   p.state.blockedDOFs = 'zXY'
   p.state.mass = 2650 * 0.1 * pi * p.shape.radius**2 # 0.1 = thickness of cylindrical particle
   inertia = 0.5 * p.state.mass * p.shape.radius**2
   p.state.inertia = (.5*inertia,.5*inertia,inertia)

O.dt = .5*utils.PWaveTimeStep()
print O.dt

triax=TriaxialStressController(
   internalCompaction=False,
   stressMask = 3,
   goal1 = -1.e7,
   goal2 = -1.e7,
   goal3 = 0.,
   label = 'triax'
)

O.engines = [
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom6D(),Ig2_Box_Sphere_ScGeom()],
      [Ip2_FrictMat_FrictMat_FrictPhys(),
       Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(
          label='cohesiveIp',
          setCohesionNow=False)],
      [Law2_ScGeom_FrictPhys_CundallStrack(),
       Law2_ScGeom6D_CohFrictPhys_CohesionMoment( # contact law
          useIncrementalForm=False,
          always_use_moment_law=False,
          label='cohesiveLaw')]
   ),
   triax,
   NewtonIntegrator(damping=.2)
]

O.run(1000,True)

Any kind of help is appreciated. Thank you very much!

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Chen Lifan
Solved:
2020-10-28
Last query:
2020-10-28
Last reply:

This question was reopened

Chen Lifan (chenlifan) said : #1

Sorry, I have solved the problem.

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

Hello,

if you ask a question and then come to solution yourself, please next time explain the solution (for further reference and for people coming here with the same problem)

> the walls didn't move at all

depends on definition of walls :-)
the "walls=aabbWalls(thickness=0,material='frictionlessWalls')" indeed does not move
BUT, the triax's "walls" DOES!
try:
b = O.bodies[triax.wall_left_id]
print(b.state.displ()) # the wall does move!
print(O.bodies[triax.wall_left_id].shape) # Sphere (!!)

Either specify wall ids explicitly, or add walls before spheres [1]

cheers
Jan

[1] https://yade-dem.org/doc/yade.wrapper.html#yade.wrapper.TriaxialStressController.wall_back_id

Chen Lifan (chenlifan) said : #3

Hello Jan,

thank you for your detailed instructions. In fact, I solved the problem by using exactly one of the ways you propose above, i.e. specifying the wall ids explicitly. For further reference, here are my scripts:

from yade import pack

O.materials.append(CohFrictMat(young=88.e9,alphaKr=2.0,alphaKtw=0.,poisson=.7,frictionAngle=.0,isCohesive=True,normalCohesion=3.e9,shearCohesion=3.e9,momentRotationLaw=True,label='rock'))
O.materials.append(FrictMat(young=2e11,poisson=.4,frictionAngle=0,density=0,label='frictionlessWalls'))
sp = pack.SpherePack()
size = .3
sp.makeCloud(minCorner=(0,0,.05),maxCorner=(size,size,.05),rMean=.005,rRelFuzz=0.2,num=400,periodic=False,seed=1)
O.bodies.append([sphere(center,rad,material='rock') for center,rad in sp])
walls=aabbWalls(thickness=0,material='frictionlessWalls')
wallIds=O.bodies.append(walls)
#O.bodies.erase(wallIds[4])
#O.bodies.erase(wallIds[5])

#O.cell.hSize = Matrix3(size,0,0, 0,size,0, 0,0,.1)
print len(O.bodies)
for i in range(400):
   p=O.bodies[i]
   p.state.blockedDOFs = 'zXY'
   p.state.mass = 2650 * 0.1 * pi * p.shape.radius**2 # 0.1 = thickness of cylindrical particle
   inertia = 0.5 * p.state.mass * p.shape.radius**2
   p.state.inertia = (.5*inertia,.5*inertia,inertia)

O.dt = .5*utils.PWaveTimeStep()
print O.dt

triax=TriaxialStressController(
   wall_bottom_id=wallIds[2],
   wall_top_id=wallIds[3],
   wall_left_id=wallIds[0],
   wall_right_id=wallIds[1],
   wall_back_id=wallIds[4],
   wall_front_id=wallIds[5],
   internalCompaction=False,
   stressMask = 3, #1 is x, 2 is y, 4 is z, Bitmask, if imposed goal`s are stresses (0 for none, 7 for all, 1 for direction 1, 5 for directions 1 and 3, etc. :ydefault:`7)
   goal1 = -1.e7,
   goal2 = -1.e7, # positive is tension, negative is compression
   goal3 = 0.,
   max_vel=0.01,
   label = 'triax'
)

O.engines = [
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom6D(),Ig2_Box_Sphere_ScGeom()],
      [Ip2_FrictMat_FrictMat_FrictPhys(),
       Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(
          label='cohesiveIp',
          setCohesionNow=False)],
      [Law2_ScGeom_FrictPhys_CundallStrack(),
       Law2_ScGeom6D_CohFrictPhys_CohesionMoment( # contact law
          useIncrementalForm=False,
          always_use_moment_law=False,
          label='cohesiveLaw')]
   ),
   triax,
   NewtonIntegrator(damping=.2)
]
O.run(1000,True)