How to set the PeriodicFlowEngine

[1]https://link.springer.com/article/10.1007/s11242-011-9915-6

Hi,

1/ No, if it is periodic there is just no walls. You simply impose a macroscopic pressure gradient. (well... it really depends what you are trying to achieve, not clear from the question).

2/ As for "should I set", it really depends on what you are trying to simulate, so I cannot answer.

Bruno

Dear Bruno,
Thank you for your response. In my mind, the flow problem is a boundary value problem, so the PeriodicFlowEngine should get the position and value of the boundary. I have read the documentation on PeriodicFlowEngine, I don’t know if my understanding is wrong. In my opinion, boundaryUseMaxMin is used to set the boundary position. The default Boolean vector of boundaryUseMaxMin is true on six faces, which means using the max-min coordinates of X, Y, and Z to add spheres as the boundary. Otherwise, the wall position is used. I saw the related problem in Answer #3 in https://answers.launchpad.net/yade/+question/688195.

another question, I have read your paper 'Pore-Scale Modeling of Viscous Flow and Induced Forces in Dense Sphere Packings', in that paper, the total force generated on one particle by the fluid includes the effects of absolute pressure and viscous stress. So I set flow.pressureForce=True, flow.viscousShear=True, flow.viscousNormalBodyStress=True, flow.viscousShearBodyStress=True, flow.normalLubrication=True, flow.shearLubrication=True. Is that correct? The complete code looks like this:

from yade import pack,qt,utils,ymport
from yade import export
from yade import timing
from yade import plot
from builtins import range
import numpy
import time
from math import *

O.load('prova.xml')

num_spheres=5000
young=750.e6
finalFricDegree = 0.35
graindensity=2600.0
poisson=0.17
en=0.9
sigmaIso=-150.e3
radiusR=0.001/2.
D50=radiusR*2.0

#O.materials.append(FrictViscoMat(young=young,poisson=poisson,frictionAngle=compFricDegree,density=graindensity,betan=0.0336,label='spheres'))

triax=PeriTriaxController(
 dynCell=True,
 mass=0.2,
 goal=[sigmaIso,sigmaIso,sigmaIso],
 stressMask=7,
 maxStrainRate=[0.0001,0.0001,0.0001],
 label='triax')

newton=NewtonIntegrator(damping=0.1)

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb()]),
 InteractionLoop(
         [Ig2_Sphere_Sphere_ScGeom()],
         [Ip2_FrictViscoMat_FrictViscoMat_FrictViscoPhys()],
         [Law2_ScGeom_FrictViscoPhys_CundallStrackVisco(traceEnergy=True)]
 ),
 PeriodicFlowEngine(dead=1,label="flow"),
 GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.5),
 triax,
 newton
]

while 1:
 O.run(1000, True)
 unb=unbalancedForce()
 if unb<0.00001:
  break

triax.dead=0
flow.dead=0

flow.pressureForce=True
flow.viscousShear=True
flow.viscousNormalBodyStress=True
flow.viscousShearBodyStress=True
flow.normalLubrication=True
flow.shearLubrication=True

flow.defTolerance=0.1
flow.meshUpdateInterval=200
flow.useSolver=3
flow.permeabilityFactor=1
flow.viscosity=0.001
flow.bndCondIsPressure=[0,0,0,0,1,1]
flow.bndCondValue=[0,0,0,0,0,0]
GlobalStiffnessTimeStepper.dead=True
O.dt=min(0.5*PWaveTimeStep(),0.5*1./1200.*pi/flow.viscosity*graindensity*(D50/2.)**2)
O.dynDt=False

O.run(10,1)

while 1:
 O.run(5000, True)
 unb=unbalancedForce()
 if unb<0.0001:
  if O.cell.velGrad[0,1]<0.01:
   O.cell.velGrad=O.cell.velGrad+Matrix3(0,0.001,0, 0,0,0, 0,0,0)
   O.dt=min(0.5*PWaveTimeStep(),0.5*1./1200.*pi/flow.viscosity*graindensity*(D50/2.)**2)
  if O.cell.velGrad[0,1]>0.01:
   break

Thanks Bruno Chareyre, that solved my question.