yade-users team mailing list archive

[Fedor <question702972@xxxxxxxxxxxxxxxxxxxxx>]

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

Fedor gave more information on the question:
Hello,
thank you very much for clarification. But the problem with zero strains and stress has not solved. I put preliminar iteration but I could not reach a confined state (-1e4). 


from __future__ import print_function

import time
import datetime
import os

from yade import pack, plot, export
import numpy as np


#FIXED PARAMETERS

poisson=0.2
R=1e-3 
rate=1e-4
dimcell = 0.03
density= 1e9
young=1e9
frictionAngle=radians(30)
targetPorosity=0.43

#SETTINGS
O.periodic = True
O.cell.hSize = Matrix3(dimcell , 0, 0, 0, dimcell , 0, 0, 0, dimcell )


sp = pack.SpherePack()
sp.makeCloud((0, 0, 0), (dimcell, dimcell, dimcell), rMean=R, rRelFuzz=.1, periodic=True)


pp = O.materials.append(CohFrictMat(
	young=young,
	poisson=poisson,
	frictionAngle=frictionAngle,
	density=density,
	isCohesive=False,
	momentRotationLaw=True,
	etaRoll=.1,label='spheres'
	))


sp.toSimulation(material='spheres')

O.engines = [
        ForceResetter(
               ),
        InsertionSortCollider([Bo1_Sphere_Aabb()]),
        InteractionLoop([Ig2_Sphere_Sphere_ScGeom6D()], [Ip2_CohFrictMat_CohFrictMat_CohFrictPhys()],[Law2_ScGeom6D_CohFrictPhys_CohesionMoment(useIncrementalForm=True,always_use_moment_law=True)]),
	TriaxialStressController(
		stressMask = 7,label='triax'),
	NewtonIntegrator(damping=.2),
        PyRunner(command='addPlotData()', iterPeriod=500),
	]

O.dt = .5 * PWaveTimeStep()
print('time step',O.dt)

triax.goal1=triax.goal2=triax.goal3=-10000


def addPlotData():
	plot.addData(
	        i=O.iter,
	        Ezz=log(O.cell.trsf[2,2]),
		Eyy=log(O.cell.trsf[1,1]),
		Exx=log(O.cell.trsf[0,0]),
		szz=utils.getStress()[2,2],
		syy=utils.getStress()[1,1],
		sxx=utils.getStress()[0,0],
	        u=utils.porosity()
	)


while 1:
	O.run(100, 1)
	print(getStress().trace() / 3.)		  
	if getStress().trace() / 3.<-10000:
  		break

triax.goal1=triax.goal2=triax.goal3=-40000
while triax.porosity>targetPorosity:
	frictionAngle = 0.95*frictionAngle
	setContactFriction(radians(30))
	print (frictionAngle," porosity:",triax.porosity)
	O.run(500,0)
	


# define what to plot
plot.plots = {
        'i ': ('sxx', 'syy', 'szz'),
        ' i': ('Exx', 'Eyy', 'Ezz'),
	' i ':('u')
        # energy plot
        
}
# show the plot
plot.plot()

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