about cpm material

Asked by xjin

I write codes like this:
------------------------------------------------------------------------------------
from yade import pack,geom,qt
from yade.gridpfacet import *
from pylab import *

O.engines=[
ForceResetter(),
InsertionSortCollider([
Bo1_Sphere_Aabb(),Bo1_GridConnection_Aabb(),Bo1_Facet_Aabb(),
]),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(),
Ig2_Facet_Sphere_ScGeom(),
Ig2_Sphere_GridConnection_ScGridCoGeom(),
Ig2_GridNode_GridNode_GridNodeGeom6D(),
Ig2_GridConnection_GridConnection_GridCoGridCoGeom(),],
[Ip2_CpmMat_CpmMat_CpmPhys(),
Ip2_FrictMat_CpmMat_FrictPhys(),
Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(setCohesionNow=True,setCohesionOnNewContacts=False)],
[Law2_ScGeom_CpmPhys_Cpm(),
Law2_ScGeom_FrictPhys_CundallStrack(),
Law2_ScGeom6D_CohFrictPhys_CohesionMoment(),
Law2_ScGridCoGeom_FrictPhys_CundallStrack()]
),
NewtonIntegrator(gravity=(0,0,-10),damping=0.1,label='newton')]

O.materials.append(CohFrictMat(young=3e2,poisson=0.3,density=1e1,frictionAngle=10,normalCohesion=1e7,shearCohesion=1e7,momentRotationLaw=True,label='spheremat'))

O.materials.append(FrictMat(young=3e2,poisson=.25,frictionAngle=0.5,density=1e1,label='face'))

O.materials.append(CpmMat(density=1e1,young=3e2,poisson=0.25,frictionAngle=0.5,sigmaT=2e1,epsCrackOnset=1e-4,relDuctility=30,label='spheresh'))
### Parameters of a rectangular grid ###
L=0.1 #length [m]
l=0.05 #width [m]
nbL=10 #number of nodes for the length [#]
nbl=5 #number of nodes for the width [#]
r=L/100.0 #radius
color=[255./255.,102./255.,0./255.]
nodesIds=[]
#Create all nodes first :
for i in range(0,nbL):
for j in range(0,nbl):
nodesIds.append( O.bodies.append(gridNode([i*L/nbL,j*l/nbl,0],r,wire=False,fixed=False,material='spheremat',color=color)) )

#Create connection between the nodes
for i in range(0,len(nodesIds)):
for j in range(i+1,len(nodesIds)):
dist=(O.bodies[i].state.pos - O.bodies[j].state.pos).norm()
if(dist<=L/nbL*1.01):
O.bodies.append( gridConnection(i,j,r,color=color) )

#Set a fixed node
O.bodies[0].dynamic=False
O.bodies[4].dynamic=False
O.bodies[45].dynamic=False
O.bodies[49].dynamic=False

O.bodies.append(sphere([0.05,0.025,0.01],0.0061,material='spheresh'))
O.bodies.append(sphere([0.04,0.023,0.01],0.0021,material='spheresh'))
O.bodies.append(sphere([0.03,0.028,0.008],0.0061,material='spheresh'))
S1=len(O.bodies)

O.bodies.append(geom.facetBox((L/2.0,l/2.0,-0.01),(L*1.2/2,l*1.2/2,L),wallMask=63,material='face'))

S2=len(O.bodies)

for i in range(S1,S2):
O.bodies[i].dynamic=False

O.dt=1e-06

-----------------------------------------------------------------------------
After some timr it runs, I get the message:
--------------------------------------------------------------------------
FATAL /home/dj/yade/trunk/pkg/dem/ConcretePM.cpp:380 go: Verification `!std::isnan(epsN)' failed!
FATAL /home/dj/yade/trunk/pkg/dem/ConcretePM.cpp:380 go: in interaction #137+#83
terminate called without an active exception
--------------------------------------------------------------------------------
can you tell me what's wrong with it?Thanks a lot!

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
Best Jan Stránský (honzik) said :
#1

Hello,

> After some timr it runs

next time please try to be more specific.. a few time steps? hundreds of thousand time steps? ... just to know what to expect during testing

as the error says, epsN (normal strain) of an interaction is NaN. This is maybe because the in CpmMat the implementation of Sphere-GridConnection is missing [1] and refLength, needed for correct computations, is not set.

Sepcifying spheremat (otherwise their material is the last O.materials.appended, i.e. CpmMat) of gridConnections (which sounds logical to me) solves the problem.
O.bodies.append( gridConnection(i,j,r,color=color,material='spheremat') )

cheers
Jan

[1] https://github.com/yade/trunk/blob/master/pkg/dem/ConcretePM.cpp#L300

Revision history for this message
xjin (jpeng22) said :
#2

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