yade-users team mailing list archive

[Luc Scholtès <question664499@xxxxxxxxxxxxxxxxxxxxx>]

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

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

Did you actually try to run the script that I gave you in #12 of the
following thread:

https://answers.launchpad.net/yade/+question/661250

?

If not, I strongly suggest that you do so and I copy paste it here once
again (please read the entire message since I add some comments for you
after the script):

--- BEGINNING OF SCRIPT ---

# -*- coding: utf-8 -*-
from yade import pack, plot

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

#### packing (previously constructed)
OUT='compressionTest_JCFPM'

#### Simulation Control
rate=-0.01 #deformation rate
iterMax=10000 # maximum number of iterations
saveVTK=2000 # saving output files for paraview

#### Material microproperties
intR=1.1 # allows near neighbour interaction (can be adjusted for every packing)
DENS=2500
YOUNG=20e9
FRICT=7
ALPHA=0.1
TENS=1e6
COH=1e6

#### material definition
def sphereMat(): return JCFpmMat(type=1,density=DENS,young=YOUNG,poisson=ALPHA,frictionAngle=radians(FRICT),tensileStrength=TENS,cohesion=COH)

#### create the specimen
pred=pack.inCylinder((0,0,0),(0,1,0),0.25)
O.bodies.append(pack.regularHexa(pred,radius=0.025,gap=0.,material=sphereMat))

R=0
Rmax=0
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
Rmean=R/nbSpheres

print 'nbSpheres=',nbSpheres,' | Rmean=',Rmean

#### boundary condition (see utils.uniaxialTestFeatures
bb=utils.uniaxialTestFeatures()
negIds,posIds,longerAxis,crossSectionArea=bb['negIds'],bb['posIds'],bb['axis'],bb['area']

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

O.engines=[
 ForceResetter(),
        InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=intR,label='Saabb')]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=intR,label='SSgeom')],
  [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1,label='interactionPhys')],
  [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=True,Key=OUT,label='interactionLaw')]
 ),
 UniaxialStrainer(strainRate=rate,axis=longerAxis,asymmetry=0,posIds=posIds,negIds=negIds,crossSectionArea=crossSectionArea,blockDisplacements=1,blockRotations=1,setSpeeds=0,stopStrain=0.1,dead=1,label='strainer'),
 GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=10,timestepSafetyCoefficient=0.5, defaultDt=utils.PWaveTimeStep()),
 NewtonIntegrator(damping=0.4,label='newton'),
 PyRunner(iterPeriod=int(100),initRun=True,command='recorder()',label='data'),
        VTKRecorder(iterPeriod=int(saveVTK),initRun=True,fileName=OUT+'-',recorders=['spheres','jcfpm','cracks'],Key=OUT,label='vtk')
]

################# RECORDER DEFINED HERE

def recorder():
    yade.plot.addData({'i':O.iter,
         'eps':strainer.strain,
         'sigma':strainer.avgStress,
         'tc':interactionLaw.nbTensCracks,
         'sc':interactionLaw.nbShearCracks,
         'te':interactionLaw.totalTensCracksE,
         'se':interactionLaw.totalShearCracksE,
         'unbF':utils.unbalancedForce()})
    plot.saveDataTxt(OUT)

################# 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

################# SIMULATION REALLY STARTS HERE
strainer.dead=0
O.run(iterMax)

--- END OF SCRIPT ---

This script corresponds to a uniaxial test on a cylinder (which is what
you want to achieve if I understood well). At the end of the simulation,
you will obtain a text file "compressionTest_JCFPM" containing the
simulation data (stress, strain, number of cracks, ...) that you can
plot with the second script I gave you in
(https://answers.launchpad.net/yade/+question/661250)  and some vtu/vtk
files which you can then use to visualize your results with Paraview
(you will have access to the particles, the cracks and the interaction
network). The "streamlines" you mention are, I guess, the displacement
of the particles. For that, there is a procedure in Paraview which
consists in comparing 2 vtu files (basically, you have to compute the
displacement by substracting the positions at t2 and t1). As mentioned
by Bruno, this is a different question.

Luc

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