yade-users team mailing list archive

[Jibril Coulibaly <question685707@xxxxxxxxxxxxxxxxxxxxx>]

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

    Status: Answered => Open

Jibril Coulibaly is still having a problem:
Hi,

Thank you very much rcaulk and bruno-chareyre for your answers. Using
flow.imposePressure() and flow.averageVelocity() provide sensible
results.

I have a couple remaining questions:

- (1) : If I try to run the script in a loop to obtain permeability in
all 3 directions, I get sensible results on the first run, then I get
zero or nan on subsequent runs. I did not get these issues when using
walls (see script below)

- (2) : Running a similar script, but reading the particles coordinates from an external file instead, I obtain the error message: "GS did not converge in 20k iterations (maybe because the reference pressure is zero?)".
I tried to figure it out from the code but I am not sure what this message is warning me about.

I appreciate your help,
Jibril B. Coulibaly

###########

from yade import pack,Vector3

O.periodic=True
O.cell.hSize=Matrix3(0.1,0,0, 0,0.1,0, 0,0,0.1)
pos_center = Vector3(0.05,0.05,0.05)
sp=pack.SpherePack()
radius=5e-3
num=sp.makeCloud((0,0,0),(.1,.1,.1),radius,.2,1000,periodic=True)
O.bodies.append([sphere(s[0],s[1]) for s in sp])

O.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb()],verletDist=.05*radius),
	InteractionLoop(
		[Ig2_Sphere_Sphere_ScGeom()],
		[Ip2_FrictMat_FrictMat_FrictPhys()],
		[Law2_ScGeom_FrictPhys_CundallStrack()]
	),	GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8),
	NewtonIntegrator(damping=0.),
	PeriodicFlowEngine(dead=0,label="flow")
]

O.dt=0.1e-8
O.dynDt=False

flow.defTolerance=0.3
flow.meshUpdateInterval=300
flow.useSolver=0 # Solver >0 not supported for PBC
flow.permeabilityFactor=1
flow.viscosity=1
flow.imposePressure(pos_center,1.0)

for i in [0,1,2]:
	if i==0:
		gradPress = Vector3(1,0,0)
	elif i==1:
		gradPress = Vector3(0,1,0)
	elif i==2:
		gradPress = Vector3(0,0,1)

	flow.updateTriangulation=True
	flow.gradP=gradPress
	flow.updateTriangulation=True
	O.run(1,1)
	aveVel = flow.averageVelocity()
	print("average velocity=(",aveVel[0],aveVel[1],aveVel[2],")")

#######################

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