yade-users team mailing list archive

[Xu Mengqian <question680783@xxxxxxxxxxxxxxxxxxxxx>]

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

Hi all:
I tried to simulate the compression test which the mean stress is constant,
The mean stress remains constant at the beginning,
but as  the particles continue to compress,
the mean stress can't remain constant later...

This is my script:
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))

filename1='Packing2D_area.txt'
a=np.loadtxt(filename1)

for p in O.bodies:
#if isinstance(b.shape,Clump):
   p.state.blockedDOFs='zXY'

def servo():
	tol=0.01
	b=1
	sr_x=-.1
	sr_y=0.
	sr_z=0.
        epsilonmax=5
	p=-100.

	sx_curr=utils.getStress()[0,0]
	sy_curr=utils.getStress()[1,1]
#	sz_curr=utils.getStress()[2,2]

#	sy_req=b*sx_curr+(1-b)*(2*p-sx_curr*(1+b))/(2-b)
        sy_req=2*p-sx_curr
#      sz_req=(2*p-sx_curr*(1+b))/(2-b)

        gy=2*epsilonmax/sy_req
#     gz=2*epsilonmax/sz_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),
]

# plotting 
plot.live=True
plot.plots={'iter':('sxx','syy','szz'),'iter_':('exx','eyy','ezz'), ' iter':('voidratio')
#			, ' iter':('unbalanced'), 
#            ' iter ':(O.energy.keys,None,'Etot')
}

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=(O.cell.size[0]*O.cell.size[1]-a)/a,
		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.run()
#O.cell.trsf=Matrix3.Identity
#O.cell.velGrad=Matrix3(sr_x,0,0,0,sr_y,0,0,0,sr_z)


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