yade-users team mailing list archive

[Cloud <question679566@xxxxxxxxxxxxxxxxxxxxx>]

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

    Status: Needs information => Open

Cloud gave more information on the question:
Hi Robert,

Thanks for your reply!

>>Are you simulating fluid?

No. Instead, to carry out the undrained test during triaxial shear, the
volume of the sample was maintained constant. this is a method to
simulate undrained condition in triaxial test without fluid.

We can check confining stress in file 'results_triax_base_' after
running the MWE, then you will find that the confining stress is
changing.

the MWE:

from yade import pack
nRead=readParamsFromTable(
 num_spheres=10000,
 compFricDegree = 30 ,
 key='_triax_base_',
 unknownOk=True
)
from yade.params import table

num_spheres=table.num_spheres
key=table.key
targetPorosity = 0.43
compFricDegree = table.compFricDegree
finalFricDegree = 22
rate=-0.02
damp=0.2
stabilityThreshold=0.01
young=5e6
mn,mx=Vector3(0,0,0),Vector3(10,10,10)
confiningS=400e3

O.materials.append(FrictMat(young=young,poisson=0.5,frictionAngle=radians(compFricDegree),density=2600,label='spheres'))
O.materials.append(FrictMat(young=young,poisson=0.5,frictionAngle=0,density=0,label='walls'))

walls=aabbWalls([mn,mx],thickness=0,material='walls')
wallIds=O.bodies.append(walls)
sp=pack.SpherePack()

sp.makeCloud(mn,mx,-1,0.3333,num_spheres,False, 0.95,seed=1)
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,truncate=1),
 newton
]

Gl1_Sphere.stripes=0
if nRead==0: yade.qt.Controller(), yade.qt.View()

triax.goal1=triax.goal2=triax.goal3=-confiningS

while 1:
 O.run(1000, True)
 unb=unbalancedForce()
 print 'unbalanced force:',unb,' mean stress: ',triax.meanStress
 if unb<stabilityThreshold and abs(-confiningS-triax.meanStress)/confiningS<0.001:
    break

O.save('confinedState'+key+'.yade.gz')
print "### Isotropic state saved ###"

triax.internalCompaction=False
triax.width0=triax.width
triax.depth0=triax.depth
triax.height0=triax.height

setContactFriction(radians(finalFricDegree))

triax.stressMask=0
triax.goal2= rate
triax.goal1=triax.goal3=-rate/2
newton.damping=0.1

O.saveTmp()

from yade import plot

def history():
  plot.addData(
            e22=-triax.strain[1], # axial strain
            e11=-triax.strain[0],
            e33=-triax.strain[2],
       ev=triax.volumetricStrain,
      s11=-triax.stress(triax.wall_right_id)[0],
      s22=-triax.stress(triax.wall_top_id)[1],
      s33=-triax.stress(triax.wall_front_id)[2],
            args=(-triax.stress(triax.wall_right_id)[0]-triax.stress(triax.wall_top_id)[1]-triax.stress(triax.wall_front_id)[2])/3,
            eta = 3*(-triax.stress(triax.wall_top_id)[1]+triax.stress(triax.wall_right_id)[0])/(-triax.stress(triax.wall_right_id)[0]-triax.stress(triax.wall_top_id)[1]-triax.stress(triax.wall_front_id)[2]),
            devi = -triax.stress(triax.wall_top_id)[1]+triax.stress(triax.wall_right_id)[0],
            i=O.iter,
            )

def stoploading():
    axial_strain = abs(triax.strain[1])
    if axial_strain > 0.3:
        print 'axial strain reach Limitation, stop loading!'
        O.pause()

if 1:
 O.engines=O.engines[0:5]+[PyRunner(iterPeriod=20,command='history()',label='recorder')]+O.engines[5:7]
 O.engines=O.engines+[PyRunner(command='stoploading()',iterPeriod=100)]
else:
 O.engines[4]=PyRunner(iterPeriod=20,command='history()',label='recorder')

while 1:
 O.run(1000,True)
 if O.iter%100 == 0:
    print "iter = ", O.iter
         print "Syy = ",abs(triax.strain[1])
        if abs(triax.strain[1]) > 0.3:
  break

plot.plots={'e22':('devi',),' args ':('devi')}
plot.plot()
plot.saveDataTxt('results'+key)
plot.saveGnuplot('plotScript'+key)

Cloud

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