Delaunay and Voronoi Tessalation#

Generate Random Networks Based on Delaunay Triangulations, Voronoi Tessellations, or Both#

A random network offers several advantages over the traditional Cubic arrangement: the topology is more ‘natural’ looking, and a wider pore size distribution can be achieved since pores are not constrained by the lattice spacing. Random networks can be tricky to generate however, since the connectivity between pores is difficult to determine or define. One surprisingly simple option is to use Delaunay triangulation to connect base points (which become pore centers) that are randomly distributed in space. The Voronoi tessellation is a complementary graph that arises directly from the Delaunay graph which also connects essentially randomly distributed points in space into a network. OpenPNM offers both of these types of network, plus the ability to create a network containing both including interconnections between Delaunay and Voronoi networks via the DelaunayVoronoiDual class. In fact, creating the dual, interconnected network is the starting point, and any unwanted elements can be easily trimmed.

import openpnm as op
import matplotlib.pyplot as plt
pn =, shape=[1, 1, 1])
―――――――――――――――――――――――――――――――――――――――――――――――――――――――――――――――――――――――――――――― : net_01
#     Properties                                    Valid Values
1     pore.coords                                     786 / 786
2     throat.conns                                   4422 / 4422
#     Labels                                        Assigned Locations
1     pore.all                                      786
2     pore.back                                     62
3     pore.bottom                                   74
4     pore.boundary                                 325
5     pore.delaunay                                 223
6     pore.front                                    68
7     pore.internal                                 461
8     pore.left                                     59
9     pore.right                                    53
10    pore.surface                                  190
11                                      65
12    pore.voronoi                                  563
13    throat.all                                    4422
14    throat.boundary                               972
15    throat.delaunay                               660
16    throat.interconnect                           2640
17    throat.internal                               2629
18    throat.surface                                821
Parameters                          Values

The above line of code is deceptively simple. The returned network (pn) contains a fully connected Delaunay network, its complementary Voronoi network, and interconnecting throats (or bonds) between each Delaunay pore (node) and its neighboring Voronoi pores. Such a highly complex network would be useful for modeling pore phase transport (i.e. diffusion) on one network (i.e. Delaunay), solid phase transport (i.e. heat transfer) on the other network (i.e. Voronoi), and exchange of a species (i.e. heat) between the solid and void phases via the interconnecting bonds. Each pore and throat is labelled accordingly (i.e. ‘pore.delaunay’, ‘throat.voronoi’), and the interconnecting throats are labelled ‘throat.interconnect’. Moreover, pores and throats lying on the surface of the network are labelled ‘surface’.

A quick visualization of this network can be accomplished using OpenPNM’s built-in graphing tool. The following shows only the Voronoi connections that lie on the surface of the cube:

Ts = pn.throats(['voronoi', 'boundary'], mode='and')
op.topotools.plot_connections(network=pn, throats=Ts)
<mpl_toolkits.mplot3d.art3d.Line3DCollection at 0x7f309ada9c40>

One central feature of these networks are the flat boundaries, which are essential when performing transport calculations since they provide well-defined control surfaces for calculating flux. This flat surfaces are accomplished by reflecting the base points across each face prior to performing the tessellations.

Plotting the internal Voronoi throats with a different color gives a good idea of the topology:

fig, ax = plt.subplots()
Ts = pn.throats(['voronoi', 'boundary'], mode='and')
op.topotools.plot_connections(network=pn, throats=Ts, alpha=0.5, ax=ax)
Ts = pn.throats(['voronoi', 'internal'], mode='and')
op.topotools.plot_connections(network=pn, throats=Ts, c='g', ax=ax)
Ts = pn.throats(['voronoi', 'surface'], mode='and')
op.topotools.plot_connections(network=pn, throats=Ts, c='r', ax=ax)
<mpl_toolkits.mplot3d.art3d.Line3DCollection at 0x7f309aed2a60>

The green lines are internal connections, and red lines are connections between internal notes and boundary nodes.

Delaunay Network#

As the name suggests, the VoronoiDelaunayDual contains both the Delaunay triangulation and the Voronoi tessellation within the same topology. It is simple to delete one network (or the other) by trimming all of the other network’s pores, which also removes all connected throats including the interconnections:

Ps = pn.pores(['voronoi'])
op.topotools.trim(network=pn, pores=Ps)
Ts = pn.throats(['surface'])
op.topotools.trim(network=pn, throats=Ts)
<mpl_toolkits.mplot3d.art3d.Line3DCollection at 0x7f309afa9490>

Create Random Networks of Spherical or Cylindrical Shape#

Many porous materials come in spherical or cylindrical shapes, such as catalyst pellets. The DelaunayVoronoiDual Network class can produce these geometries by specifying the domain_size in cylindrical [r, z] or spherical [r] coordinates:

cyl =, shape=[2, 5])
op.topotools.plot_connections(network = cyl, throats=cyl.throats('boundary'))
<mpl_toolkits.mplot3d.art3d.Line3DCollection at 0x7f309ad2d700>
sph =, shape=[2])
op.topotools.plot_connections(network = sph, throats=sph.throats('surface'))
<mpl_toolkits.mplot3d.art3d.Line3DCollection at 0x7f3098353280>

Note that the cylindrical and spherical networks don’t look very nice when too few points are used, so at least about 200 is recommended.

Assign Pore Sizes to the Random Network#

With pore centers randomly distributed in space it becomes challenging to know what pore size to assign to each location. Assigning pores that are too large results in overlaps, which makes it impossible to properly account for porosity and transport lengths. OpenPNM includes a Geometry model called largest_sphere that solves this problem. Let’s assign the largest possible pore size to each Voronoi node in the sph network just created:

Ps = sph.pores('voronoi')
Ts = sph.throats('voronoi')
geom = op.geometry.GenericGeometry(network=sph, pores=Ps, throats=Ts)
mod = op.models.geometry.pore_size.largest_sphere
geom.add_model(propname='pore.diameter', model=mod)
mod = op.models.geometry.throat_length.ctc
geom.add_model(propname='throat.length', model=mod)
mod = op.models.geometry.throat_size.from_neighbor_pores
geom.add_model(propname='throat.diameter', model=mod)
geom.show_hist(['pore.diameter', 'throat.length', 'throat.diameter'])
AttributeError                            Traceback (most recent call last)
Input In [10], in <cell line: 3>()
      1 mod = op.models.geometry.pore_size.largest_sphere
      2 geom.add_model(propname='pore.diameter', model=mod)
----> 3 mod = op.models.geometry.throat_length.ctc
      4 geom.add_model(propname='throat.length', model=mod)
      5 mod = op.models.geometry.throat_size.from_neighbor_pores

AttributeError: module 'openpnm.models.geometry.throat_length' has no attribute 'ctc'

The resulting geometrical properties can be viewed with geom.plot_histograms() (note that each realization will differ slightly):