Periodic flow and permeability

Hello Jibril,

I'm not sure you'll find solace in my answer, but a quick review of your post lead me to these two points:

1st point: In the MWE you present, you are indexing cells before a triangulation is created. FlowEngine builds the triangulation during the first time step that it is active. So it would be necessary to make a step (O.run(1,1)) before assigning pressures. Alternatively, you could use flow.emulateAction()[1] to create the triangulation without running a mechanical timestep.

2nd point:
># Solver >0 not supported for PBC

Should not be true for yadedaily versions and compiled trunk versions.

Cheers,

Robert

[1]https://yade-dem.org/doc/yade.wrapper.html#yade.wrapper.PeriodicFlowEngine.emulateAction

Hello Robert,

Thank you very much for your help and prompt answer.
I am running Yade 2018.02b that I installed running "sudo apt-get install yade", that must be why I do not have access to "solver>0".

Regarding the building of the triangulation, moving the "O.run(1,1)" before assigning pressures indeed stops the error. However, it still yields zero boundary flux, so I must be missing a piece of the setup still.

I tried running a second time using "O.run(1,1)" after assigning the pressures, hoping it would run the system with the properly assigned pressures but I get a segmentation fault instead.

I appreciate your guidance,
Jibril B. Coulibaly

This question was expired because it remained in the 'Open' state without activity for the last 15 days.

Setting to answered to avoid launchpad janitor from deleting thread.

Hi,
A periodic problem has no boundaries and therefore "flow.getBoundaryFlux(1)" is meaningless (Body #1 is not a boundary anyway).
What you want is "flow.averageVelocity()".

In order to impose pressure at (x,y,z) coordinates you need flow.imposePressure(). Better avoid playing with mesh elements directly since they change all the time.
Cheers
Bruno

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

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

Hi,
Your script reveals that writing "flow.updateTriangulation=True" repeatedly, at each iteration, leads to problems with periodic solver. It needs to be investigated. Most likely some internals are not updated correctly after assigning gradP multiple times.

In fact changing the boundary conditions repeatedly at one point in time (i.e. in 0 iteration) is not something supposed to happen in a concrete simulation and probably it was not anticipated.

Nevertheless, you can always use three engines instead of one, so you are sure there is no such issue. Below is a shorter, working, version of your script. Results with useSolver=3 (default):

average velocity= Vector3(-4.537650575699425e-6,-4.7346547328556296e-8,-5.300733352662439e-8)
average velocity= Vector3(-4.734654732786302e-8,-4.190262426966041e-6,5.93967279201477e-8)
average velocity= Vector3(-5.300733352643585e-8,5.9396727919603324e-8,-4.479154269170494e-6)

Bruno

from yade import pack

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

for p in [Vector3(1,0,0),Vector3(0,1,0),Vector3(0,0,1)]:
 flow=PeriodicFlowEngine()
 flow.imposePressure(pos_center,0)
 flow.gradP=p
 flow.emulateAction()
 flow.saveVtk()
 print("average velocity=",flow.averageVelocity())

The problem is reported here: https://gitlab.com/yade-dev/trunk/issues/141