yade-users team mailing list archive

[Lukas Wasner <question703081@xxxxxxxxxxxxxxxxxxxxx>]

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

    Status: Answered => Open

Lukas Wasner is still having a problem:
Hello Jan,
thank you for the answer.
Regarding to point [1], here is the "get to the right porosity" Script:

readParamsFromTable(rMean=.0012, sig3=-100)
from yade import pack, plot
from yade.params.table import *

mn,mx=Vector3(0,0,0),Vector3(.1,.1,.1)
compFricDegree=30
damp=0.2
young=5e6
stabilityThreshold=0.01
targetPorosity=.4

O.materials.append(FrictMat(young=young,poisson=0.5,frictionAngle=radians(compFricDegree),density=2500,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,rMean=rMean,seed=1)
sp.makeCloud(mn,mx,rMean=.5*rMean,seed=1)
sp.toSimulation(material='spheres')

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

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
]

triax.goal1=triax.goal2=triax.goal3=sig3

while 1:
	O.run(1000, True)
	unb=unbalancedForce()
	print('unbalanced force:',unb,' mean stress: ',triax.meanStress)
	if unb<stabilityThreshold and abs(sig3-triax.meanStress)/sig3<0.001:
		print('Compression equilibrium reached')
		break

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

if triax.porosity<targetPorosity:
		sp.fromSimulation()
		sp.save(f'packtest_rMean{rMean}_sig{sig3}_n_{triax.porosity}.txt')
		print('Pack ist exportiert')
		O.pause()

waitIfBatch()

---------------------------------------------------------------------
Regarding point [2]:

I have tried the volume argument and unfortunately get different results:
-------------------------------------------------------------------------------
For script 3 (Predicate is a box).
I get a similar porosity when I insert the volume of the predicate into the porosity function.
.5672 without argument
.56898 with Predicate Volume as argument

from yade import pack

sp=SpherePack()
sp.load('packtest_rMean0.0012_sig-100.txt')

print('target porosity is 0.38542')

pred=pack.inAlignedBox(minAABB=(-.01,-.01,-.01), maxAABB=(.01, .01, .01))
O.bodies.append(geom.facetBox((0, 0, 0), (.01, .01, .01), wallMask=31))


spf=pack.filterSpherePack(predicate=pred,spherePack=sp,returnSpherePack=True)
spf.toSimulation()
print('porosity without argument',porosity())
print('porosity with argument volume',porosity(volume=(.02*.02*.02)))

---------------------------------------------------------------------------------------------------------------

For script 2 (Predicate is a cylinder).
I get a lower porosity if I insert the volume of the predicate into the porosity function. Nevertheless, this is much too high.

----------------------------------------------------------------------------------------------------------------

from yade import pack

sp=SpherePack()
sp.load('packtest_rMean0.0012_sig-100.txt')
print('target porosity is 0.38542')

O.bodies.append(geom.facetCylinder(center=(0, 0, .01), radius=.01, height=.02,segmentsNumber=50, wallMask=6))
pred=pack.inCylinder(centerBottom=(0,0,0), centerTop=(0,0,.02), radius=.01)

spf2=pack.filterSpherePack(predicate=pred,spherePack=sp,returnSpherePack=True)
spf2.toSimulation()
print('porosity without argument',porosity())
print('porosity with argument', porosity(volume=((pi*.01**2)*.02)))

------------------------------------------------------------------------------------------------------------

Unfortunately, at this point I do not know how to proceed, since now only the proof has been given that the porosity function uses the volume of the predicate, which it is supposed to do.
Why there is a larger deviation with the cylinder predicate, I cannot explain myself.

I would still like to find a way to work with the porosity function,
since I would like to map the pore number over the stress later in the
oedometer experiment.

I conclude that the filtered section of the imported sphere packing has
a smaller total particle volume in relation to the new (predicate)
volume than the total particle volume in relation to the final
compressed particle cloud (old volume) from [1].

Since in [1] the particle cloud was compressed isotropically, I assumed
that the porosity is equally distributed over the entire volume and the
cropped part accordingly also has the same porosity as the entire sphere
packing volume.

How do I manage to get the same porosity in my predicate as it is in the
old volume (imported pack)?

Thanks in advance,
Lukas

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