← Back to team overview

yade-users team mailing list archive

Re: [Question #688685]: permeability - PFV

 

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

Othman Sh posted a new comment:
Thanks Robert for your help.

Bruno, please see the code below.

----------------
# -*- coding: utf-8 -*-


from __future__ import print_function
from builtins import range
from yade import pack, ymport
radiuscyl=.08
heightcyl=.25
dP=4e3 #Pa
visc=10 #Pa.sec
density=1000 #kg/m3
g=9.81 #m/s2
young=1e6


O.materials.append(FrictMat(young=young,poisson=0.2,density=1900,label='spheres'))
O.materials.append(FrictMat(young=young,poisson=0.5,frictionAngle=0,density=0,label='walls'))

yade.qt.View()

############################# spheres #############################

sp=pack.SpherePack()
sp.makeCloud((0,0,0),(.4,.4,2),rMean=0.0083,rRelFuzz=0.1,seed=1)	
### cylinder extraction 										
pred=pack.inCylinder((.2,.2,0),(.2,.2,heightcyl),radiuscyl) 						
spFilter=filterSpherePack(pred,sp,Material=Material, returnSpherePack=True)				
spFilter.toSimulation()

Height=max(utils.aabbDim())


#Fix all particles in their positions. No deformation
for i in O.bodies:
	i.state.blockedDOFs='xyzXYZ'


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()],label="iloop"
	),
	FlowEngine(dead=1,label="flow"),#introduced as a dead engine for the moment, see 2nd section
	GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8),
	NewtonIntegrator(damping=0.2)
]

#B. Activate flow engine and set boundary conditions in order to get permeability
flow.dead=0
flow.useSolver=3
flow.permeabilityFactor=1
flow.viscosity=visc
flow.bndCondIsPressure=[0,0,0,0,1,1]
flow.bndCondValue=[0,0,0,0,0,dP]
flow.boundaryUseMaxMin=[1,1,1,1,1,1]
O.dt=0.1e-8
O.dynDt=False

O.run(1,1)

Q=max(flow.averageVelocity())
Permeability=abs((density*g*Q*Height)/dP)
print ("average velocity: ",Q, "Permeability m/s: ",Permeability,"Permeability in/hr: ", Permeability*141732)

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