yade-users team mailing list archive

[maquankun <question673124@xxxxxxxxxxxxxxxxxxxxx>]

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

I tried to make the displacement field of the second stage of triaxial compression. I have this answer in the following question, but I have not implemented this method. Here is my code, please tell me where the error is.
https://answers.launchpad.net/yade/+question/478212.
from yade import pack,plot,qt,export
nRead=readParamsFromTable(
        num_spheres=1000,
        compFricDegree = 30,
        key='_triax_base_',
        unknownOk=True
)
from yade.params import table
num_spheres=table.num_spheres# number of spheres
key=table.key
targetPorosity = 0.42
compFricDegree = table.compFricDegree
finalFricDegree = 30
rate=-0.02
damp=0.2
stabilityThreshold=0.01
young=5e8
mn,mx=Vector3(0,0,0),Vector3(1,2,1) # corners of the initial packing
weiya=100000
poisson=0.1
O.materials.append(FrictMat(young=young,poisson=poisson,frictionAngle=radians(compFricDegree),density=2650,label='spheres'))
O.materials.append(FrictMat(young=young,poisson=poisson,frictionAngle=0,density=0,label='walls'))
walls=aabbWalls([mn,mx],thickness=0,material='walls')
wallIds=O.bodies.append(walls)
psdSizes,psdCumm =[0.19,0.191,0.195,0.196,0.2],[0.,0.1,0.5,0.6,1.]
                    #[0.01,0.03,0.11,0.13,0.2],[0.,0.1,0.5,0.6,1.]
                    # [0.05,0.07,0.14,0.155,0.2],[0.,0.1,0.5,0.6,1.]
                    #[0.1,0.11,0.16,0.17,0.2],[0.,0.1,0.5,0.6,1.]
                    #[0.15,0.16,0.184,0.19,0.2],[0.,0.1,0.5,0.6,1.]
                    #[0.19,0.191,0.195,0.196,0.2],[0.,0.1,0.5,0.6,1.]
sp=pack.SpherePack();
sp.makeCloud(mn,mx,num=num_spheres, psdSizes=psdSizes, psdCumm=psdCumm,distributeMass=True,porosity=0.42)
O.bodies.append([sphere(center, rad, material='spheres') for center, rad in sp])
triax=TriaxialStressController(
        maxMultiplier=1.+2e4/young,
        finalMaxMultiplier=1.+2e3/young,
        thickness = 0,
        stressMask = 7,
        internalCompaction=False,
)
newton=NewtonIntegrator(damping=damp)
O.engines=[
        ForceResetter(),
        InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
        InteractionLoop(
                [Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
                [Ip2_FrictMat_FrictMat_FrictPhys()],
                [Law2_ScGeom_FrictPhys_CundallStrack()]
        ),
        GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8),
        triax,
        TriaxialStateRecorder(iterPeriod=100,file='WallStresses'+table.key),
        newton,
]
yade.qt.Controller(), yade.qt.View()
triax.goal1=triax.goal2=triax.goal3=-weiya
while 1:
  O.run(1000, True)
  unb=unbalancedForce()
  print 'unbalanced force:',unb,' mean stress: ',triax.meanStress
  if unb<stabilityThreshold and abs(-weiya-triax.meanStress)/weiya<0.001:
    break
print "###      Isotropic state saved      ###"
import sys
while triax.porosity>targetPorosity:
	compFricDegree = 0.95*compFricDegree
	setContactFriction(radians(compFricDegree))
	print "\r Friction: ",compFricDegree," porosity:",triax.porosity,
	sys.stdout.flush()
	O.run(500,1)
print 'Box Volume: ', triax.boxVolume
print "###    Compacted state saved      ###"
triax.internalCompaction=False
setContactFriction(radians(finalFricDegree))
triax.stressMask = 5
triax.goal1=-weiya
triax.goal2=rate
triax.goal3=-weiya
newton.damping=0.1
O.saveTmp()
vtk = export.VTKExporter('/home/ma/dspl/jan-')
vtk.exportSpheres(useRef=True,what=[('dspl','b.state.displ()')])

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