yade-users team mailing list archive

[Bettina Suhr <question269739@xxxxxxxxxxxxxxxxxxxxx>]

New question #269739 on Yade:
https://answers.launchpad.net/yade/+question/269739

Dear all,

I found a nice publication for the verification of the Hertz-Mindlin contact model:
Chung, Y. & Ooi, J. Benchmark tests for verifying discrete element modelling codes at particle impact level Granular Matter, Springer-Verlag, 2011, 13, 643-656

There are 8 simple test defined for verification. For the first test  (normal contact of two identical spheres) my result produced with Yade agrees with the analytical solution. 
The next test is about normal contact between a sphere and a plane and here I get differences to the analytical solution. 
I use the following script:

from yade import plot
import math
import sys
import numpy as np

O.materials.append(FrictMat(young=7.0e10, poisson=0.3, frictionAngle=np.arctan(0.0),    density=2699.,label="alu"))
O.bodies.append([
      sphere(center=(0,0,0.11),radius=0.1, material='alu')
])
O.bodies.append(utils.wall(0, axis=2, sense=1))

#asign initial velocity
O.bodies[0].state.vel=(0,0,-0.2)

O.engines=[
  ForceResetter(),
  #use aabb for approx. contact detection
  InsertionSortCollider([Bo1_Sphere_Aabb(), Bo1_Wall_Aabb(), Bo1_Facet_Aabb()]),
  InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(), Ig2_Facet_Sphere_ScGeom(), Ig2_Wall_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_MindlinPhys()],#ELASTIC = No Damping 
  [Law2_ScGeom_MindlinPhys_Mindlin()]
  ),
  NewtonIntegrator(damping=0.0, gravity=(0,0,0)),#NO DAMPING IN INTEGRATOR
  PyRunner(command='myMonitor()',iterPeriod=1)
]
O.dt=0.1*utils.PWaveTimeStep() # extra save value
print 'time step used: ', O.dt, ' critical time step size ', utils.PWaveTimeStep()
O.trackEnergy=True

t_start_contact=0

def myMonitor():
  global t_start_contact
  if len(O.interactions)>0 and O.interactions[0,1].geom:
    if np.abs(t_start_contact)<1e-10:
      t_start_contact=O.time
    for i in O.interactions:
       plot.addData(d=i.geom.penetrationDepth,Fn=np.abs(i.phys.normalForce[2]), t=O.time-t_start_contact)

plot.plots={'d':('Fn'), 't':('Fn')}

# show the plot on the screen, and update while the simulation runs
plot.plot()



Do I do something wrong? Is this an error in Yade or is this a difference by intention?
Thanks for your help,
Bettina 



-- 
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.