Adding boundary pores

In this example we will illustrate how boundary pores are added during the network extraction phase. We will use the individual steps of the SNOW algorithm (i.e. porespy.filters.snow_partioning), rather than use the pre-packaged porespy.networks.snow function.

Import packages

import porespy as ps
import matplotlib.pyplot as plt
import numpy as np
import openpnm as op
from edt import edt
from copy import copy

ws = op.Workspace()
ws.settings['loglevel'] = 50
ps.visualization.set_mpl_style()

[01:06:16] ERROR    PARDISO solver not installed, run `pip install pypardiso`. Otherwise,          _workspace.py:56
                    simulations will be slow. Apple M chips not supported.                                         

Create image

im = ps.generators.lattice_spheres(shape=[52, 52], r=10, offset=13, spacing=25)
ps.imshow(im);

snow = ps.filters.snow_partitioning(im)
fig, ax = plt.subplots(1, 3, figsize=[10, 4])
ax[0].imshow(snow.dt, origin='lower', interpolation='none')
ax[1].imshow(snow.peaks > 0, origin='lower', interpolation='none')
ax[2].imshow(snow.regions, origin='lower', interpolation='none');

bds = ps.networks.add_boundary_regions(regions=snow.regions, pad_width=[5, 5])
cm = copy(plt.cm.nipy_spectral)
plt.figure(figsize=[4, 4])
plt.imshow(bds, origin='lower', interpolation='none', cmap=cm);

The image above shows how boundary regions appear. They are 3-voxel thick slices that are padded onto each face and given a unique label to differentiate them from the internal pores to which they connect. There are a few key things to note:

• Boundary regions that touch each other are separated by a small solid obstacle indicated by black. This is to prevent boundary pores from becoming connected to each other when the network extraction function is applied (regions_to_network).

• Corners are also converted to solid obstacles. This is for the same reason as above, but also because it’s not clear which internal pores these corners should be connected to.

• A thickness of 3 voxels was chosen since this was necessary for the marching_cubes surface area calculation. It is also a good idea for the pore centroid to be padded by 1-voxel on each side. This 3-voxel dimension can pose problems during the network extraction stage since most of the physical dimensions will be meaningless. Consequently care must be taken to ensure that these boundary pores are treated correctly. In many cases it may be preferrable to forego the them entirely and just apply boundary conditions on real pores lying on the domain surface.

Passing this image into the regions_to_network function will treat each region as an actual pore:

net = ps.networks.regions_to_network(bds)

Next we’ll open the net dictionary in OpenPNM using the PoreSpy IO class:

pn = op.io.network_from_porespy(net)
fig = plt.figure()
fig = op.visualization.plot_connections(pn, c='w', linewidth=2, ax=fig)
fig = op.visualization.plot_coordinates(pn, c='w', s=200, ax=fig)
plt.imshow(bds.T, cmap=cm, origin='lower')
plt.axis('off');

The approach outlined above is missing one important feature. The boundary pores are not labelled. When the process is performed internally by the regions_to_network function this is taken care of.