← Back to team overview

yade-users team mailing list archive

Re: [Question #680783]: 2D clump in the compression test which the mean stress is constant

 

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

Xu Mengqian posted a new comment:
Hi,
The first script is clump consolidation;

from yade import pack,plot,qt,export,utils
O.periodic=True
O.materials.append(FrictMat(young=6.e8,poisson=.8,frictionAngle=radians(26),density=2650))

sp = pack.SpherePack()
size = .01
sp.makeCloud(minCorner=(0,0,.0020),maxCorner=(size,size,.0020),rMean=.0002,rRelFuzz=.4,num=400,periodic=True,seed=1)
sp.toSimulation()
O.cell.hSize = Matrix3(size,0,0, 0,size,0, 0,0,.004)
print len(O.bodies)
for p in O.bodies:
   p.state.blockedDOFs = 'zXY'
   p.state.mass = 2650 * 0.004 * pi * p.shape.radius**2 # 0.1 = thickness of cylindrical particle
   inertia = 0.5 * p.state.mass * p.shape.radius**2
   p.state.inertia = (.5*inertia,.5*inertia,inertia)

O.engines = [
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom()],
      [Ip2_FrictMat_FrictMat_FrictPhys()],
      [Law2_ScGeom_FrictPhys_CundallStrack()]
   ),
   GlobalStiffnessTimeStepper(),
   NewtonIntegrator(damping=0.7),
   PeriTriaxController(
      dynCell=True,
      goal=(-1.e2,-1.e2,0),
      stressMask=3,
      relStressTol=.001,
      maxUnbalanced=.001,
      maxStrainRate=(.3,.3,.0),
      doneHook='Finished()',
      label='biax'
   ),
     PyRunner(command='plotAddData()',iterPeriod=100),
]

def plotAddData():
	plot.addData(
		iter=O.iter,iter_=O.iter,
		sxx=biax.stress[0],syy=biax.stress[1],szz=biax.stress[2],
		exx=O.cell.size[0],eyy=O.cell.size[1],ezz=O.cell.size[2],
		Z=avgNumInteractions(),
		Zm=avgNumInteractions(skipFree=True),
                voidratio=utils.voidratio2D(0.004),
		unbalanced=utils.unbalancedForce(),
		t=O.time,
		gWork=O.energy['gravWork'],
		Ep=O.energy['elastPotential'],
		Edamp=O.energy['nonviscDamp'],
		Ediss=O.energy['plastDissip'],
		Ekin=utils.kineticEnergy(),
      		Etot=O.energy.total(),**O.energy
		
	)
        plot.saveDataTxt('macroFile',vars=('t','exx','eyy','ezz','sxx','syy','szz','Z','Zm','voidratio'))
        plot.saveDataTxt('energyFile',vars=('t','Etot','unbalanced','gWork','Edamp','Ekin'))

O.trackEnergy=True

def Finished():
	O.save('isotropicState.xml')
	print 'Finished'
	print utils.voidratio2D(0.004)
	O.pause()

yade.qt.Controller(), yade.qt.View()  
O.run()
plot.plot(subPlots=True)

The second script is the compression test:

from yade import pack,plot,qt,export
from math import fabs
import numpy as np
O.load('isotropicState.xml')
O.trackEnergy=True
setContactFriction(radians(26.6))

for p in O.bodies:
   p.state.blockedDOFs='zXY'

def servo():
	tol=0.01
	b=1
	sr_x=-.1
	sr_y=0.
	sr_z=0.
        epsilonmax=10
	p=-100
	sx_curr=utils.getStress()[0,0]
	sy_curr=utils.getStress()[1,1]
        sy_req=2*p-sx_curr
        gy=2*epsilonmax/sy_req
	sry_curr=O.cell.velGrad[1,1]
	srz_curr=O.cell.velGrad[2,2]
	if fabs(sy_curr-sy_req)>=tol:
                sr_y=-gy*(sy_req-sy_curr)
	else:
		sr_y=sry_curr
	O.cell.velGrad=Matrix3(sr_x,0,0,0,sr_y,0,0,0,sr_z)

O.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb()]),
	InteractionLoop(
		[Ig2_Sphere_Sphere_ScGeom()],
   		[Ip2_FrictMat_FrictMat_FrictPhys()],
   		[Law2_ScGeom_FrictPhys_CundallStrack()]),
	PyRunner(command='fabric()',iterPeriod=10000),
	NewtonIntegrator(damping=0.2),
        GlobalStiffnessTimeStepper(),
	PyRunner(command='servo()', iterPeriod=1),
	PyRunner(command='plotAddData()',iterPeriod=100),
]

plot.live=True
plot.plots={'iter':('sxx','syy','szz'),'iter_':('exx','eyy','ezz'), ' iter':('voidratio')
}
def plotAddData():
	plot.addData(
		iter=O.iter,iter_=O.iter,
		sxx=utils.getStress()[0,0],
		syy=utils.getStress()[1,1],
		szz=utils.getStress()[2,2],
		exx=O.cell.size[0],
		eyy=O.cell.size[1],
		ezz=O.cell.size[2],
		Z=avgNumInteractions(),
		Zm=avgNumInteractions(skipFree=True),
                voidratio=utils.voidratio2D(0.004),
		unbalanced=utils.unbalancedForce(),
		t=O.time,
		gWork=O.energy['gravWork'],
		Edamp=O.energy['nonviscDamp'],
		Ekin=utils.kineticEnergy(),
      		Etot=O.energy.total(),**O.energy
		
	)
	plot.saveDataTxt('macroFile',vars=('t','exx','eyy','ezz','sxx','syy','szz','Z','Zm','voidratio'))
	plot.saveDataTxt('energyFile',vars=('t','Etot','unbalanced','gWork','Edamp','Ekin'))

O.run()

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