yade-users team mailing list archive

[Luc Scholtès <question689058@xxxxxxxxxxxxxxxxxxxxx>]

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

Luc Scholtès posted a new comment:
Hi Jasper,

I think you can achieve what you want with the JCFPM. In fact, you can
reproduced the effect of bending/twisting resistance by increasing the
interaction radius when using the JCFPM. The effect is "similar" in
terms of material behavior/strength. Also, as you know, in such context
(breakage of samples due to impact), fractures/fragments will develop as
a result of the presence of flaws in the material (pre-existing
fractures and/or weakness planes), and the JCFPM was made for that
purpose and I think it might be relevant for your study.

Also, be aware that the way the numerical sample breaks is directly
related to its mechanical properties (and the loading conditions as
well). Maybe you could not achieve a "representative" simulation because
you chose a non adequate set of parameters?

ANW, please find below a "cartoon" simulation that might be useful for
you.

Luc

###

from yade import pack, plot

################# SIMULATIONS DEFINED HERE

#### Material microproperties
intR=1.2 # allows near neighbour interaction (can be adjusted for every packing)

#### material definition
def sphereMat(): return JCFpmMat(type=1,density=3000,young=5e9,poisson=0.2,frictionAngle=radians(15),tensileStrength=5e6,cohesion=5e7)
def wallMat(): return JCFpmMat(type=0,density=3000,young=1e9,poisson=0.2,frictionAngle=radians(15),tensileStrength=0,cohesion=0)

#### create the specimen
pred=pack.inCylinder((0,-0.5,0.5),(0.,0.5,0.5),0.25)
#O.bodies.append(pack.regularHexa(pred,radius=0.025,gap=0.,material=sphereMat))
O.bodies.append(pack.randomDensePack(pred,radius=0.025,rRelFuzz=0.2,useOBB=True,spheresInCell=3000,memoizeDb='/tmp/gts-triax-packings.sqlite',returnSpherePack=False,color=(0.9,0.8,0.6),material=sphereMat))

R=0
Rmax=0
Rmin=1e6
nbSpheres=0.
for o in O.bodies:
 if isinstance(o.shape,Sphere):
   nbSpheres+=1
   R+=o.shape.radius
   if o.shape.radius>Rmax:
     Rmax=o.shape.radius
   if o.shape.radius<Rmin:
     Rmin=o.shape.radius
Rmean=R/nbSpheres

print('nbSpheres=',nbSpheres,' | Rmean=',Rmean, ' | Rmax/Rmin=',
Rmax/Rmin)

O.bodies.append(geom.facetBox((0,0,0),(2,2,1),Quaternion((1,1,0),pi/20.),wallMask=16,color=(1,1,1),wire=False,fixed=1,material=wallMat))

################# ENGINES DEFINED HERE

O.engines=[
    ForceResetter(),
    InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=intR,label='Saabb'),Bo1_Facet_Aabb()]),
    InteractionLoop(
        [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=intR,label='SSgeom'),Ig2_Facet_Sphere_ScGeom()],
        [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1,label='interactionPhys')],
        [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=True,Key='test',label='interactionLaw')]
    ),
    GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=10,timestepSafetyCoefficient=0.9,defaultDt=utils.PWaveTimeStep()),
    NewtonIntegrator(damping=0.4,gravity=(0,0,-100)),
    VTKRecorder(iterPeriod=2000,initRun=True,fileName='fallingBlock-',recorders=['facets','spheres','jcfpm','cracks'],Key='test',label='vtk')
]

################# PREPROCESSING

#### manage interaction detection factor during the first timestep and then set default interaction range
O.step();
### initializes the interaction detection factor
SSgeom.interactionDetectionFactor=-1.
Saabb.aabbEnlargeFactor=-1.

#### coordination number verification
numSSlinks=0
numCohesivelinks=0
for i in O.interactions:
    if isinstance(O.bodies[i.id1].shape,Sphere) and isinstance(O.bodies[i.id2].shape,Sphere):
      numSSlinks+=1
    if i.phys.isCohesive :
      numCohesivelinks+=1

print ("nblinks=", numSSlinks, " | nbCohesivelinks=", numCohesivelinks,
"|| Kcohesive=", 2.0*numCohesivelinks/nbSpheres)

#### delete floating particles
Kmin=3
erased=0
for o in O.bodies:
    if isinstance(o.shape,Sphere):
        nbCont=0
        for i in O.interactions.withBody(o.id) :
            if not i.isReal : continue
            if i.phys.isCohesive:
                nbCont+=1
        if nbCont<Kmin:
            erased+=1
            O.bodies.erase(o.id)
print 'nb of floating particles (erased) = ', erased

################# SIMULATION REALLY STARTS HERE
O.run(14000)

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