yade-users team mailing list archive

[Leonard <question684106@xxxxxxxxxxxxxxxxxxxxx>]

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

Hi,
I am using example from [1] to achieve target porosity (for short, the MWE), the example works pretty well when ''targetPorosity=0.43". However, when I modify the 'targetPorosity' from 0.43 to 0.7 (because I'd like to generate a loose sand sample), and after the "Compacted state saved", I use "triax.porosity" to check the porosity of the sample, it is only 0.446.
I tried some ways such as changing the number of spheres, but it seems like the largest porosity they can achieve is around 0.45.
Here is the MWE, it may take 30 seconds to run it on Yade 2018.02b, ubuntu 18.04.
#####
from yade import pack, plot
nRead=readParamsFromTable(
	num_spheres=3000,
	compFricDegree = 30,
	key='_triax_base_',
	unknownOk=True
)
from yade.params import table

num_spheres=table.num_spheres
key=table.key
targetPorosity = 0.7 ### 0.43 works well
compFricDegree = table.compFricDegree
finalFricDegree = 30
rate=-0.001
damp=0.2
stabilityThreshold=0.001
young=5e7

mn,mx=Vector3(0,0,0),Vector3(0.7,1.4,0.7)

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])
Gl1_Sphere.quality=3

triax=TriaxialStressController(

	maxMultiplier=1.004,
	finalMaxMultiplier=1.0004,
	thickness = 0,
	stressMask = 7,
	internalCompaction=True,
)

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),
	newton
]

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

triax.goal1=triax.goal2=triax.goal3=-10000
while 1:
  O.run(1000, True)
  unb=unbalancedForce()
  print 'unbF:',unb,' meanStress: ',-triax.meanStress,'top:',-triax.stress(triax.wall_top_id)[1]
  if unb<stabilityThreshold and abs(-10000-triax.meanStress)/10000<0.001:
	break

print "###      Isotropic state saved      ###"

import sys
while triax.porosity>targetPorosity:
	compFricDegree = 0.95*compFricDegree
	setContactFriction(radians(compFricDegree))
	print "\r Friction: ",compFricDegree," porosity:",triax.porosity,
	sys.stdout.flush()
	O.run(500,1)

print "###    Compacted state saved      ###"
###############
I think [1] is a very good example to decrease porosity, and I tried a stupid way to inverse it like:
####
while triax.porosity<targetPorosity:
	compFricDegree = 1.05*compFricDegree
	setContactFriction(radians(compFricDegree))
	print "\r Friction: ",compFricDegree," porosity:",triax.porosity,
	sys.stdout.flush()
	O.run(500,1)
####
And it failed, the porosity hardly changes. Is there any way to reach target porosity if  current porosity < targetPorosity?

Thanks in advance,
Leonard

[1]https://github.com/yade/trunk/blob/master/examples/triax-tutorial/script-session1.py

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