Problem with physics

[DMaxJ]

I am having trouble understanding some results of simulations performed on networks extracted using SNOW algorithm. I have used the following approach

1. GEO1: (in the code below): standard pore and throat diameter with throat conduit length
2. GEO2: geometry models using equivalent diameter
3. PHYSICS1: generic physics and then adding Hagen-Poisseuille model
4. PHYSICS2: Standard physics model

import openpnm as op
import os
import numpy as np

ws = op.Workspace()
ws.clear()
fname='IEKXX'
path = 'M:\\...\\TEST_SNOW'
os.chdir(path)

ws.load_project(fname)
proj=ws[fname]
pn=proj.network

Lx = np.amax(pn['pore.coords'][:, 0]) - np.amin(pn['pore.coords'][:, 0])
Ly = np.amax(pn['pore.coords'][:, 1]) - np.amin(pn['pore.coords'][:, 1])
Lz = np.amax(pn['pore.coords'][:, 2]) - np.amin(pn['pore.coords'][:, 2])

# GEO 1
mod = op.models.geometry.throat_endpoints.spherical_pores
pn.add_model(propname='throat.endpoints', model=mod)
mod = op.models.geometry.throat_length.conduit_lengths
pn.add_model(propname='throat.conduit_lengths', model=mod)
mod = op.models.geometry.pore_area.sphere
pn.add_model(propname='pore.area', model=mod)
geo = op.geometry.GenericGeometry(network=pn, pores=pn.Ps, throats=pn.Ts)
geo.regenerate_models()

"""
# GEO 2
mod = op.models.geometry.throat_endpoints.spherical_pores
pn.add_model(propname='throat.endpoints', model= mod, 
             pore_diameter='pore.equivalent_diameter', 
             throat_diameter='throat.equivalent_diameter')
mod = op.models.geometry.throat_length.conduit_lengths
pn.add_model(propname='throat.conduit_lengths', model=mod)
mod = op.models.geometry.pore_area.sphere
pn.add_model(propname='pore.area', model=mod, 
             pore_diameter='pore.equivalent_diameter')
geo = op.geometry.GenericGeometry(network=pn, pores=pn.Ps, throats=pn.Ts)
"""
"""
# PHYSICS 1
air = op.phases.Air(network=pn, name='ph_air')
phys_air = op.physics.GenericPhysics(network=pn, phase=air, geometry=geo, 
                              name='phy_air')
air.add_model(propname='throat.hydraulic_conductance',
               model=op.models.physics.hydraulic_conductance.hagen_poiseuille)
air.regenerate_models()
"""
# PHYSICS 2
air = op.phases.Air(network=pn, name='ph_air')
phys_air = op.physics.Standard(network=pn, phase=air, geometry=geo, 
                              name='phy_air')

Pin=101325
Pout=0
delta_P=Pin-Pout
md=1/0.98e-12*1000

perm = op.algorithms.StokesFlow(network=pn, project=proj, name='stokes')
perm.setup(phase=air)

perm.set_value_BC(pores=pn.pores('left'), values=Pout)
perm.set_value_BC(pores=pn.pores('right'), values=Pin)
perm.run()
air.update(perm.results())

Q = perm.rate(pores=pn.pores('right'), mode='group')
mu = air['pore.viscosity'].max()
K_lr = Q*Lx*mu/(Ly*Lz*delta_P)
K_x = K_lr*md
print('The value of K_lr is:', K_x, 'mD')

Now the results look like this: 1) GEO 1 + PHYSICS 1 61.76 mD 2) GEO 1 + PHYSICS 2 94.05 mD 3) GEO 2 + PHYSICS 1 160.31 mD 4) GEO 2 + PHYSICS 2 88.98 mD

Why am I seeing this variation? I see that the standard physics model also uses Hagen-Poisseulle model. Could anyone please enlighten me? Or is there any mistake?

[ma-sadeghi]

Standard physics also includes a "flow shape factor", which corrects for the fact that pores are not cylinders. To back it up, all our basic conductance models assume that each conduit, i.e. pore + throat + pore, is a series of cylinders. To take into account, for instance the fact that pores are spherical, we use shape factors. Currently, we support two distinct shape factors: ball_and_stick, which assumes pores are spherical and throats are still cylindrical, and conical_frustum_and_stick, which assumes pores are conical frustums (sliced cones) and throats are cylinders. The default shape factor is ball_and_stick, which is by default added in Standard physics, hence the difference in permeability values.

[DMaxJ]

Thanks. It seems I overlooked that aspect. Thanks for the clarification.

[DMaxJ]

When we use generic physics with Hagen-Poisseuille model, then which shape factor is used in OpenPNM? I see that when I use standard physics model with equivalent diameter, the calculated permeability decreases (as compared to generic physics+HP), but the opposite happens when I just use pore diameter and not equivalent diameter.

[ma-sadeghi]

When we use generic physics with Hagen-Poisseuille model, then which shape factor is used in OpenPNM?

No shape factor will be used. Basically, all our classes starting with Generic are templates, meaning that they're empty and only serve as placeholders. You need to manually add models.

I see that when I use standard physics model with equivalent diameter, the calculated permeability decreases (as compared to generic physics+HP), but the opposite happens when I just use pore diameter and not equivalent diameter.

When using shape factors, it's really hard to tell which direction the permeability changes (i.e. increase or decrease). If you exclude shape factors, however, it's easier to predict: if the equivalent diameter is (on average) smaller than pore diameter, you should expect a decrease in permeability.

[DMaxJ]

Thank you. I have the following output

pn['pore.diameter'][10:20]
array([8.24485840e-06, 5.60060999e-06, 8.85239043e-06, 4.69439382e-06,
       3.23719016e-06, 4.36263791e-06, 5.87886151e-06, 9.23880206e-06,
       9.42700852e-06, 5.60060999e-06])

pn['pore.equivalent_diameter'][10:20]
array([1.98855916e-05, 1.50374982e-05, 2.53402125e-05, 1.15277275e-05,
       8.18680566e-06, 9.82535020e-06, 1.35105332e-05, 1.90852088e-05,
       2.56704970e-05, 1.31705242e-05])

It becomes clear from the results that without shape factor the permeability increases by a factor of 3 when equivalent diameter is used. But the same cannot be observed with shape factor. There the calculated permeability value is smaller with equivalent diameter. Your explanation makes it clear.

Any suggestions regarding what could be a/the realistic approach?

[ma-sadeghi]

Can you print the average diameters? I mean print average pore.diameter and average pore.equivalent_diameter?

Any suggestions regarding what could be a/the realistic approach?

@borahdj2010 In my experience (which is limited), you're always better off using "a" shape factor, rather than not. The reason is that pores are almost certainly not a bundle of tubes, but rather much closer to sliced cones or deformed spheres.

@Zohaib-Atiq is much more experienced, he probably has more to offer.

[DMaxJ]

Yes average pore.diameter = 6.175725394373993e-06 average pore.equivalent_diameter = 1.4394915727386392e-05

median pore.diameter = 5.600609985337265e-06 average pore.equivalent_diameter = 1.2571392133939607e-05

[ma-sadeghi]

I'm closing this for now, if you have more questions, feel free to follow up here or open up another issue.