yade-users team mailing list archive

[Yuri Bezmenov <question706817@xxxxxxxxxxxxxxxxxxxxx>]

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

Hi, 
I am simulating Triaxial Test using TriaxialStressController. 
I noticed some particles flying out of the boundary in the simulation. i.e. all particles should be within the boundary but some particles are flying beyond the walls. 

here is my code 



from yade import pack, qt, plot, export
compFricDegree=30  # contact friction during the confining phase
finalFricDegree = 30  # contact friction during the deviatoric loading
rate = -0.001  # loading rate (strain rate)
stabilityThreshold = 0.01  # we test unbalancedForce against this value in different loops 
young = 5e6  # contact stiffness
mn, mx = Vector3(0, 0, 0), Vector3(1, 1, 1)
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)
# psdSizes, psdCumm = [75e-6, 150e-6, 212e-6, 300e-6, 425e-6, 600e-6, 1000e-6, 2000e-6, 3350e-6, 4750e-6], [0, 0.14, 0.22, 0.39, 0.48, 0.62, 0.72, 0.85, 0.97, 1]
sp = pack.SpherePack()
sp.makeCloud(mn, mx, rMean=.07, rRelFuzz=.03, periodic=True)
# sp.makeCloud(mn, mx, psdSizes=psdSizes, psdCumm=psdCumm, distributeMass=True, num=20000)
O.bodies.append([sphere(center, rad, material='spheres') for center, rad in sp])

triax = TriaxialStressController(
	thickness=0,
	stressMask=7,
	internalCompaction=False,
)

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,
	NewtonIntegrator(damping=0.2),
	# VTKRecorder(fileName='/home/files/RigidTriax/vtkFiles', recorders=['all'], iterPeriod=500, label="visul", dead=True),
	PyRunner(command='vtk()', iterPeriod=100, label="visul", dead=True), 
	PyRunner(command='addPlotData()', iterPeriod=100, label="graph", dead=True),
]

def addPlotData():
	plot.addData(
		e11=-triax.strain[0], e22=-triax.strain[1], 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],
		i=O.iter,
		n=triax.porosity
		# Etot=O.energy.total(),
		# **O.energy
	)
	
O.run(1000,True)

vtkExporter = export.VTKExporter('/home/files/RigidTriax/test')
plot.plots = {
        'i':('e11','e22','e33'),
        'i ':('s11','s22','s33'),
        'e22':('s11','s22','s33'),
	  ' i ': ('n','ev')
        # energy plot
        # ' i ': (O.energy.keys, None, 'Etot', 'ev'),
}

plot.plot()

def vtk():
	vtkExporter.exportSpheres()
	vtkExporter.exportInteractions()
	vtkExporter.exportFacets()


###   APPLYING CONFINING PRESSURE   ###
triax.goal1 = triax.goal2 = triax.goal3 = -100000
while 1:
	O.run(1000, True)
	unb=unbalancedForce()
	print ('unbalanced force:',unb,' mean stress: ',triax.meanStress)
	if unb<stabilityThreshold and abs(-100000-triax.meanStress)/100000<0.0005:
		break
# O.save('confinedState.yade.gz')
print ("###      Isotropic state reached      ###")

triax.depth0=triax.depth
triax.height0=triax.height
triax.width0=triax.width
O.step()

visul.dead=False
graph.dead=False

###   DEVIATORIC LOADING   ###
# triax.internalCompaction=False
setContactFriction(radians(finalFricDegree))
triax.stressMask = 5
triax.goal2=rate
triax.goal1=-100000
triax.goal3=-100000
while -triax.strain[1] <0.15:
	O.run(10000,1)		######change to 10000
	print ("10000 step foreward	Strain=",-triax.strain[1]*100,"%")
print ("###    deviatoric loading finished      ###")
plot.saveDataTxt('triaxial.txt')

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