Predicting absolute permeability#
The example explains absolute permeabilty calculations on a cubic network. Note that permeability calculation for an extracted network from PoreSpy follows similar steps in assigning phase, algorithm and calculating permeability.
import numpy as np
import openpnm as op
op.visualization.set_mpl_style()
np.random.seed(10)
%matplotlib inline
np.set_printoptions(precision=5)
Create a random cubic network#
pn = op.network.Cubic(shape=[15, 15, 15], spacing=1e-6)
pn.add_model_collection(op.models.collections.geometry.spheres_and_cylinders)
pn.regenerate_models()
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
Cell In[2], line 3
1 pn = op.network.Cubic(shape=[15, 15, 15], spacing=1e-6)
2 pn.add_model_collection(op.models.collections.geometry.spheres_and_cylinders)
----> 3 pn.regenerate_models()
File /opt/hostedtoolcache/Python/3.10.18/x64/lib/python3.10/site-packages/openpnm/core/_models.py:476, in ModelsMixin2.regenerate_models(self, propnames, exclude)
474 for item in propnames:
475 try:
--> 476 self.run_model(item)
477 except KeyError as e:
478 msg = (f"{item} was not run since the following property"
479 f" is missing: {e}")
File /opt/hostedtoolcache/Python/3.10.18/x64/lib/python3.10/site-packages/openpnm/core/_models.py:507, in ModelsMixin2.run_model(self, propname, domain)
505 if item.startswith(propname+"@"):
506 _, _, domain = item.partition("@")
--> 507 self.run_model(propname=propname, domain=domain)
508 else: # domain was given explicitly
509 domain = domain.split('.', 1)[-1]
File /opt/hostedtoolcache/Python/3.10.18/x64/lib/python3.10/site-packages/openpnm/core/_models.py:521, in ModelsMixin2.run_model(self, propname, domain)
519 if 'domain' not in inspect.getfullargspec(mod_dict['model']).args:
520 _ = kwargs.pop('domain', None)
--> 521 vals = mod_dict['model'](self, **kwargs)
522 if isinstance(vals, dict): # Handle models that return a dict
523 for k, v in vals.items():
File /opt/hostedtoolcache/Python/3.10.18/x64/lib/python3.10/site-packages/openpnm/models/misc/_neighbor_lookups.py:123, in from_neighbor_pores(target, prop, mode, ignore_nans)
121 network = target.network
122 throats = target.Ts
--> 123 P12 = network.find_connected_pores(throats)
124 pvalues = target[prop][P12]
125 if ignore_nans:
File /opt/hostedtoolcache/Python/3.10.18/x64/lib/python3.10/site-packages/openpnm/network/_network.py:443, in Network.find_connected_pores(self, throats, flatten, mode)
391 r"""
392 Return a list of pores connected to the given list of throats
393
(...)
440
441 """
442 Ts = self._parse_indices(throats)
--> 443 pores = topotools.find_connected_sites(bonds=Ts, network=self,
444 flatten=flatten, logic=mode)
445 return pores
File /opt/hostedtoolcache/Python/3.10.18/x64/lib/python3.10/site-packages/openpnm/topotools/_graphtools.py:37, in find_connected_sites(bonds, **kwargs)
36 def find_connected_sites(bonds, **kwargs):
---> 37 return queries.find_connected_nodes(inds=bonds, **kwargs)
File /opt/hostedtoolcache/Python/3.10.18/x64/lib/python3.10/site-packages/openpnm/_skgraph/queries/_funcs.py:167, in find_connected_nodes(network, inds, flatten, logic)
165 temp = temp.astype(float)
166 temp[inds] = np.nan
--> 167 temp = np.reshape(a=temp, newshape=[len(edges), 2], order='F')
168 neighbors = temp
169 else:
TypeError: reshape() got some positional-only arguments passed as keyword arguments: 'a'
Create phase object#
It is assumed that a generic phase flows through the porous medium. As absolute permeability is the porous medium property and not the fluid property, any other fluid with an assigned viscosity value can be used as the phase.
Permeability of the network is then found by applying Stokes Flow algorithm. To simulate the fluid flow algorithm, hydraulic conductance of the conduits must be defined. Here, the model collection is used to assign basic pore-scale models including generic_hydraulic to the phase.
phase = op.phase.Phase(network=pn)
phase['pore.viscosity']=1.0
phase.add_model_collection(op.models.collections.physics.basic)
phase.regenerate_models()
[14:35:54] WARNING throat.entry_pressure was not run since the following _models.py:480 property is missing: 'throat.surface_tension'
WARNING throat.diffusive_conductance was not run since the _models.py:480 following property is missing: 'throat.diffusivity'
Apply Stokes flow#
To calculate permeability in x direction, a constant pressure boundary condition is applied on the left and right side of the network. Note that a similar procedure can be followed to find the permeability in y and z directions.
inlet = pn.pores('left')
outlet = pn.pores('right')
flow = op.algorithms.StokesFlow(network=pn, phase=phase)
flow.set_value_BC(pores=inlet, values=1)
flow.set_value_BC(pores=outlet, values=0)
flow.run()
phase.update(flow.soln)
Note
The Solution attribute
`flow.soln` is a `dict` with the `quantity` as the `key` (i.e. `'pore.pressure'`) and the solution as the `value` (i.e an `ndarray`). The last line in the cell above updates the `phase` with the new computed values of `pore.pressure` from solving the Stokes flow transport algorithm.
We can visulalize the pressure in the network:
ax = op.visualization.plot_connections(pn)
ax = op.visualization.plot_coordinates(pn, ax=ax, color_by=phase['pore.pressure'])
/home/amin/Code/OpenPNM/openpnm/visualization/_plottools.py:341: UserWarning: The figure layout has changed to tight
fig.tight_layout()
Calculate permeability#
Calculate the permeability using Darcy’s law:
\[ K_{abs}= \frac{Q}{A} \frac{\mu L} {\Delta P} \]
where \(Q\) is the inlet flow rate, \(A\) is the inlet area, and \(L\) is the distance between inlet and outlet. As pressure difference and viscosity were assumed to be 1, we have a simplified equation.
# NBVAL_IGNORE_OUTPUT
Q = flow.rate(pores=inlet, mode='group')[0]
A = op.topotools.get_domain_area(pn, inlets=inlet, outlets=outlet)
L = op.topotools.get_domain_length(pn, inlets=inlet, outlets=outlet)
# K = Q * L * mu / (A * Delta_P) # mu and Delta_P were assumed to be 1.
K = Q * L / A
print(f'The value of K is: {K/0.98e-12*1000:.2f} mD')
[14:35:56] WARNING Attempting to estimate inlet area...will be low _topotools.py:1042
WARNING Attempting to estimate domain length...could be low if _topotools.py:1086 boundary pores were not added
The value of K is: 0.07 mD
Note
Finding permeability of extracted network
1) The methods of finding domain area in `topotools` is based on Scipy's `ConvexHull`, where a convexhull that includes the inlet pores will be created to approximate the inlet area. Both `get_domain_area` and `get_domain_length` can be a useful approximation in estimating area and length in extracted networks. In extracted networks without boundary pores where inlet/outlet pores are not necessarily located on an almost flat plane, the estimated value could be low.
2) In this example we assumed the network has spherical pores and cylindrical throats. Different geometrical shapes for pores and throats are defined in `geometry` collection. Each geometry collection includes pore-scale size factor models that are necessary for finding hydraulic conductance of the conduits and applying transport algorithm, accordingly. These size factor models can alternatively be assigned using `add_model` method to the network and choosing the model from `op.models.geometry.hydraulic_size_factors`. For more information on these pore-scale models see size factor example notebook.