yade-users team mailing list archive

[VG <question290950@xxxxxxxxxxxxxxxxxxxxx>]

Question #290950 on Yade changed:
https://answers.launchpad.net/yade/+question/290950

VG posted a new comment:
Thanks Jan!
I modified the generation of top and bottom boxes as you suggested. The boxes seem to be generated the way they should be. However, as soon as the simulation starts, those boxes just disappear. I am not sure what is causing this behavior. Here is the updated MWE script:


# Generate a periodic cell with few super-particles.
# Each superparticle consists of smaller sub-particles of 75 microns
# bonded with cohesive bonds.

from yade import pack,qt,plot,utils,export,ymport
from math import *
from pylab import rand
import datetime
import os,shutil

#############################################################################
# Set up run
#############################################################################
run_name="PeriodicCohesive_MWE"
data_root_dir="."
if not os.path.exists(data_root_dir):
  print "The data root directory you specified:"
  print data_root_dir
  print "does not exist. Exiting..."
  exit()
timestamp=datetime.datetime.now().strftime("%Y.%m.%d-%H%M%S")
run_dir_name=timestamp+'-'+run_name
data_dir_path=os.path.join(data_root_dir,run_dir_name)
if not os.path.exists(data_dir_path):
  os.makedirs(data_dir_path)
else:
  print("Something is really wrong if you get this message. Exiting...\n")
  exit()

run_name_base_path=os.path.join(data_dir_path,run_name+'-')
# opts.script is defined in yade, which inserts the script file (e.g.
# this file) into itself.
this_script=os.path.abspath(opts.script)
shutil.copy(this_script,data_dir_path)


#############################################################################
# Materials
#############################################################################
plate_material=CohFrictMat(
    young=200e9
   ,poisson=0.3
   ,density=8000
   ,frictionAngle=radians(30)
   ,normalCohesion=1e10
   ,shearCohesion=1e10
   ,momentRotationLaw=True
   ,label='plate_mat')
O.materials.append(plate_material)


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


#############################################################################
# Component dimensions and operating condition
#############################################################################
# Granular material dimension
sample_diameter=2e-4
sample_radius=sample_diameter/2.0
# Sub-particle dimension
particle_diameter=74.e-6
particle_radius=particle_diameter/2.

############################################################################# 
# grinding plate dimension
#############################################################################

rotvel=2./3.*pi*(1.5+0.5)*.254

#############################################################################
# Periodic Geometry
#############################################################################


# Set up periodic boundary conditions
O.periodic=True
xExt=4*sample_diameter
yExt=3.*sample_diameter*2 #to block the periodicity in y direction
zExt=xExt
xLim=xExt
yLim=yExt/4
zLim=zExt
O.cell.hSize=Matrix3(
  xExt, 0, 0,
  0, yExt, 0,
  0, 0, zExt)


length=xExt
height=yExt
width=zExt

# Top and bottom plate thickness
thickness=0.1*height


#########
# Bottom
#########
#bottom=box(
#    center=(length/2.0,yLim - thickness/2.0,width/2.0)
#   ,extents=(4*length,thickness/2.0,width)
#   ,wire=False
#   ,material='plate_mat'
#   )

bottomBoxes = []
for ix in (0,1,2):
 for iz in (0,1,2):
   bottomBoxes.append(box( # create 3x3 boxes with 1/3 cell size
   center=(xExt/6.*(1+2*ix),yLim - thickness/2.0,zExt/6.*(1+2*iz))
   ,extents=(xExt/6.,thickness/2.0,zExt/6.)
   ,wire=False
   ,material='plate_mat'
   ))

bottom_id,bottom_ids = O.bodies.appendClumped(bottomBoxes) # bottom_id
is the clump id,


O.bodies[bottom_id].state.blockedDOFs='xyzXYZ'


#############################################################################
# Particle Packing
#############################################################################

min_corner= (0,yLim,0)
max_corner= (xLim, yExt-yLim, zLim)

sp=pack.SpherePack()
sp.makeCloud( min_corner,max_corner, rMean=sample_radius, periodic=False)

print "Generated ",len(sp)," aggregates"

###########################################
# Sample
###########################################
for s in sp:
  sphere=pack.inSphere((s[0][0],s[0][1],s[0][2]),s[1])
  sp1=pack.randomDensePack(
      sphere
     ,spheresInCell=2000
     ,radius=particle_radius
     ,memoizeDb='/tmp/triaxPackCache.sqlite'
     ,returnSpherePack=True
     )

  sp1.toSimulation(material='sample_mat',color=(0.9,0.8,0.6))
  print 'Generated ',len(sp1),' particles'

Gl1_Sphere(stripes=True)

#########
# Top
#########
#top=box(
#    center=(length/2.0,yExt - yLim + thickness*2,width/2.0)
#   ,extents=(0.16*length,thickness/2.0,0.16*width)
#   ,wire=False
#   ,material='plate_mat'
#   )

#top_id=O.bodies.append(top)

#########
 # Top
#########
topBoxes = []
for ix in (0,1,2):
 for iz in (0,1,2):
   topBoxes.append(box( # create 3x3 boxes with 1/3 cell size
   center=(xExt/6.*(1+2*ix),yExt - yLim + thickness/2.0,zExt/6.*(1+2*iz))
   ,extents=(xExt/6.,thickness/2.0,zExt/6.)
   ,wire=False
   ,material='plate_mat'
   ))

top_id,top_ids = O.bodies.appendClumped(topBoxes) # top_id is the clump
id,


O.bodies[top_id].state.blockedDOFs='xzXYZ'


plate_downforce=-0.0036

O.forces.addF(top_id,(0,plate_downforce,0),permanent=True)

O.bodies[top_id].state.vel[0]= rotvel


#############################################################################
# Run the simulation
#############################################################################


O.dt=0.5*PWaveTimeStep()

O.engines=[
  ForceResetter(),
  InsertionSortCollider([
    Bo1_Sphere_Aabb()
   ,Bo1_Box_Aabb()
   ], allowBiggerThanPeriod=True
   ),
  InteractionLoop(
    [Ig2_Sphere_Sphere_ScGeom6D()
    ,Ig2_Box_Sphere_ScGeom6D() 
    ],
    [
     Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(setCohesionNow=True,label="cohesiveIp")
    ],
    [
     Law2_ScGeom6D_CohFrictPhys_CohesionMoment()
    ]
  ), # End InteractionLoop 
  NewtonIntegrator(damping=0.8,gravity=(0.,0.,0.)),

]

O.saveTmp()

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.