Tutorial 2 of 3: Digging Deeper into OpenPNM

This tutorial will follow the same outline as Getting Started, but will dig a little bit deeper at each step to reveal the important features of OpenPNM that were glossed over previously.

Learning Objectives

• Explore different network topologies, and learn some handy topological query methods
• Create a heterogeneous domain with different geometrical properties in different regions
• Learn about data exchange between objects
• Utilize pore-scale models for calculating properties of all types
• Propagate changing geometrical and thermo-physical properties to all dependent properties
• Calculate the permeability tensor for the stratified media
• Use the Workspace Manager to save and load a simulation

Building a Cubic Network

As usual, start by importing the OpenPNM and Scipy packages:

Let's generate a cubic network again, but with a different connectivity:

• This Network has pores distributed in a cubic lattice, but connected to diagonal neighbors due to the connectivity being set to 8 (the default is 6 which is orthogonal neighbors). The various options are outlined in the Cubic class's documentation which can be viewed with the Object Inspector in Spyder.

• OpenPNM includes several other classes for generating networks including random topology based on Delaunay tessellations (Delaunay).

• It is also possible to import networks <data_io>_ from external sources, such as networks extracted from tomographic images, or that networks generated by external code.

Initialize and Build Multiple Geometry Objects

One of the main functionalities of OpenPNM is the ability to assign drastically different geometrical properties to different regions of the domain to create heterogeneous materials, such as layered structures. To demonstrate the motivation behind this feature, this tutorial will make a material that has different geometrical properties on the top and bottom surfaces compared to the internal pores. We need to create one Geometry object to manage the top and bottom pores, and a second to manage the remaining internal pores:

• The above statements result in two distinct Geometry objects, each applying to different regions of the domain. geom1 applies to only the pores on the top and bottom surfaces (automatically labeled 'top' and 'bottom' during the network generation step), while geom2 applies to the pores 'not' on the top and bottom surfaces.

• The assignment of throats is more complicated and illustrates the find_neighbor_throats method, which is one of the more useful topological query methods <topology>_ on the Network class. In both of these calls, all throats connected to the given set of pores (Ps1 or Ps2) are found; however, the mode argument alters which throats are returned. The terms 'union' and 'intersection' are used in the "set theory" sense, such that 'union' returns all throats connected to the pores in the supplied list, while 'intersection' returns the throats that are only connected to the supplied pores. More specifically, if pores 1 and 2 have throats [1, 2] and [2, 3] as neighbors, respectively, then the 'union' mode returns [1, 2, 3] and the 'intersection' mode returns [2]. A detailed description of this behavior is given in :ref:topology.

Assign Static Seed values to Each Geometry

In :ref:getting_started we only assigned 'static' values to the Geometry object, which we calculated explicitly. In this tutorial we will use the pore-scale models that are provided with OpenPNM. To get started, however, we'll assign static random seed values between 0 and 1 to each pore on both Geometry objects, by assigning random numbers to each Geometry's 'pore.seed' property:

• Each of the above lines produced an array of different length, corresponding to the number of pores assigned to each Geometry object. This is accomplished by the calls to geom1.Np and geom2.Np, which return the number of pores on each object.
• Every Core object in OpenPNM possesses the same set of methods for managing their data, such as counting the number of pore and throat values they represent; thus, pn.Np returns 1000 while geom1.Np and geom2.Np return 200 and 800 respectively.

Accessing Data Distributed Between Geometries

The segmentation of the data between separate Geometry objects is essential to the management of pore-scale models, although it does create a complication: it's not easy to obtain a single array containing all the values of a given property for the whole network. It is technically possible to piece this data together manually since we know the locations where each Geometry object applies, but this is tedious so OpenPNM provides a shortcut. First, let's illustrate the manual approach using the 'pore.seed' values we have defined:

The following code illustrates the shortcut approach, which accomplishes the same result as above in a single line:

• This shortcut works because the pn dictionary does not contain an array called 'pore.seed', so all associated Geometry objects are then checked for the requested array(s). If it is found, then OpenPNM essentially performs the interleaving of the data as demonstrated by the manual approach and returns all the values together in a single full-size array. If it is not found, then a standard KeyError message is received.
• This exchange of data between Network and Geometry makes sense if you consider that Network objects act as a sort of master object relative Geometry objects. Networks apply to all pores and throats in the domain, while Geometries apply to subsets of the domain, so if the Network needs some values from all pores it has direct access.

Add Pore Size Distribution Models to Each Geometry

Pore-scale models are mathematical functions that are applied to each pore (or throat) in the network to produce some local property value. Each of the modules in OpenPNM (Network, Geometry, Phase and Physics) have a "library" of pre-written models located under "models" (i.e. Geometry.models). Below this level, the models are further categorized according to what property they calculate, and there are typical 2-3 models for each. For instance, under Geometry.models.pore_diameter you will see random, normal and weibull among others.

Pore size distribution models are assigned to each Geometry object as follows:

Pore-scale models tend to be the most complex (i.e. confusing) aspects of OpenPNM, so it's worth dwelling on the important points of the above two commands:

• Both geom1 and geom2 have a models attribute where the parameters specified in the add command are stored for future use if/when needed. The models attribute actually contains a ModelsDict object which is a customized dictionary for storing and managing this type of information.

• The propname argument specifies which property the model calculates. This means that the numerical results of the model calculation will be saved in their respective Geometry objects as geom1['pore.diameter'] and geom2['pore.diameter'].

• Each model stores it's result under the same propname but these values do not conflict since each Geometry object presides over a unique subset of pores and throats.

• The model argument contains a handle to the desired function, which is extracted from the models library of the relevant Module (Geometry in this case). Each Geometry object has been assigned a different statistical model, normal and weibull. This ability to apply different models to different regions of the domain is reason multiple Geometry objects are permitted. The added complexity is well worth the added flexibility.

• The remaining arguments are those required by the chosen model. In the above cases, these are the parameters that define the statistical distribution. Note that the mean pore size for geom1 will be 20 um (set by scale) while for geom2 it will be 50 um, thus creating the smaller surface pores as intended. The pore-scale models are well documented regarding what arguments are required and their meaning; as usual these can be viewed with Object Inspector in Spyder.

Now that we've added pore diameter models the each Geometry we can visualize the network in Paraview to confirm that distinctly different pore sizes on the surface regions:

Add Additional Pore-Scale Models to Each Geometry

In addition to pore diameter, there are several other geometrical properties needed to perform a permeability simulation. Let's start with throat diameter:

Instead of using statistical distribution functions, the above lines use the neighbor model which determines each throat value based on the values found 'pore_prop' from it's neighboring pores. In this case, each throat is assigned the minimum pore diameter of it's two neighboring pores. Other options for mode include 'max' and 'mean'.

We'll also need throat length as well as the cross-sectional area of pores and throats, for calculating the hydraulic conductance model later.

Create a Phase Object and Assign Thermophysical Property Models

For this tutorial, we will create a generic Phase object for water, then assign some pore-scale models for calculating their properties. Alternatively, we could use the prewritten Water class included in OpenPNM, which comes complete with the necessary pore-scale models, but this would defeat the purpose of the tutorial.

Note that all Phase objects are automatically assigned standard temperature and pressure conditions when created. This can be adjusted:

A variety of pore-scale models are available for calculating Phase properties, generally taken from correlations in the literature. An empirical correlation specifically for the viscosity of water is available:

Create Physics Objects for Each Geometry

Physics objects are where geometric information and thermophysical properties are combined to produce the pore and throat scale transport parameters. Thus we need to create one Physics object for EACH Phase and EACH Geometry:

Next add the Hagan-Poiseuille model to both:

• The same function (mod) was passed as the model argument to both Physics objects. This means that both objects will calculate the hydraulic conductance using the same function. A model must be assigned to both objects in order for the 'throat.hydraulic_conductance' property be defined everywhere in the domain since each Physics applies to a unique selection of pores and throats.
• The "pore-scale model" mechanism was specifically designed to allow for users to easily create their own custom models. Creating custom models is outlined in :ref:advanced_usage.

Accessing Data Distributed Between Multiple Physics Objects

Just as Network objects can retrieve data from separate Geometries as a single array with values in the correct locations, Phase objects can retrieve data from Physics objects as follows:

• Each Physics applies to the same subset for pores and throats as the Geometries so its values are distributed spatially, but each Physics is also associated with a single Phase object. Consequently, only a Phase object can to request all of the values within the domain pertaining to itself.
• In other words, a Network object cannot aggregate the Physics data because it doesn't know which Phase is referred to. For instance, when asking for 'throat.hydraulic_conductance' it could refer to water or air conductivity, so it can only be requested by water or air.

Pore-Scale Models: The Big Picture

Having created all the necessary objects with pore-scale models, it is now time to demonstrate why the OpenPNM pore-scale model approach is so powerful. First, let's inspect the current value of hydraulic conductance in throat 1 on phys1 and phys2:

Now, let's alter the Geometry objects by assigning new random seeds, and adjust the temperature of water.

So far we have not run the regenerate command on any of these objects, which means that the above changes have not yet been applied to all the dependent properties. Let's do this and examine what occurs at each step:

These two lines trigger the re-calculation of all the size related models on each Geometry object.

This line causes the viscosity to be recalculated at the new temperature. Let's confirm that the hydraulic conductance has NOT yet changed since we have not yet regenerated the Physics objects' models:

Finally, if we regenerate phys1 and phys2 we can see that the hydraulic conductance will be updated to reflect the new sizes on the Geometries and the new temperature on the Phase:

Determine Permeability Tensor by Changing Inlet and Outlet Boundary Conditions

The :ref:getting started tutorial <getting_started> already demonstrated the process of performing a basic permeability simulation. In this tutorial, we'll perform the simulation in all three perpendicular dimensions to obtain the permeability tensor of our heterogeneous anisotropic material.

Set boundary conditions for flow in the X-direction:

The resulting pressure field can be seen using Paraview:

To determine the permeability coefficient we must find the flow rate through the network to use in Darcy's law. The StokesFlow class (and all analogous transport algorithms) possess a rate method that calculates the net transport through a given set of pores:

To find K, we need to solve Darcy's law: Q = KA/(mu*L)(P_in - P_out). This requires knowing the viscosity and macroscopic network dimensions:

The dimensions of the network can be determined manually from the shape and spacing specified during its generation:

The pressure drop was specified as 1 atm when setting boundary conditions, so Kxx can be found as:

We can either create 2 new Algorithm objects to perform the simulations in the other two directions, or reuse alg by adjusting the boundary conditions and re-running it.

The first call to set_boundary_conditions used the overwrite mode, which replaces all existing boundary conditions on the alg object with the specified values. The second call uses the merge mode which adds new boundary conditions to any already present, which is the default behavior.

A new value for the flow rate must be recalculated, but all other parameters are equal to the X-direction:

The values of Kxx and Kyy should be nearly identical since both these two directions are parallel to the small surface pores. For the Z-direction:

The permeability in the Z-direction is about half that in the other two directions due to the constrictions caused by the small surface pores.