Skip to content

Commit

Permalink
FIxing badges, and updating BC example
Browse files Browse the repository at this point in the history
  • Loading branch information
jgostick committed Mar 6, 2016
1 parent 37aeed5 commit 729e9a8
Show file tree
Hide file tree
Showing 4 changed files with 78 additions and 77 deletions.
2 changes: 1 addition & 1 deletion README.md
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
![https://api.travis-ci.org/PMEAL/OpenPNM-Examples.svg](https://travis-ci.org/PMEAL/OpenPNM-Examples)
![](https://api.travis-ci.org/PMEAL/OpenPNM-Examples.svg)

# OpenPNM-Examples
A collection of scripts illustrating the use of OpenPNM
Expand Down
2 changes: 1 addition & 1 deletion Simulations/bond_and_site_percolation.md
Original file line number Diff line number Diff line change
Expand Up @@ -91,4 +91,4 @@ The stored results can be visualized with Matplotlib on the same axes as the bon

The percolation threshold for bond percolation is lower that for site percolation. The theoretical values are [24.88 and 31.16 respectively](http://en.wikipedia.org/wiki/Percolation_threshold#Thresholds_on_3D_lattices), which agrees with the present results considering the small size of the Network used here.

![http://imgur.com/BZTSD6e.png]()
![]](http://imgur.com/BZTSD6e.png)
145 changes: 73 additions & 72 deletions Simulations/exploring_boundary_conditions.md
Original file line number Diff line number Diff line change
@@ -1,131 +1,132 @@
.. _boundary_conditions_example:
# Diffusion Simulations with Various Boundary Conditions

===============================================================================
Diffusion Simulations with Various Boundary Conditions
===============================================================================
This example outlines how to perform a variety of transport simulations on pore networks. It uses the Fickian diffusion algorithm as the basis of the example and demonstrates a number of different boundary conditions specifications.

+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Generating Network, Geometry, Phases and Physics
+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
## Generating Network, Geometry, Phases and Physics

Start by generating a basic cubic network and the other required components:

.. code-block:: python
``` python
>>> import OpenPNM
>>> pn = OpenPNM.Network.Cubic(name='net',shape=[10,10,10])
>>> pn.add_boundaries()

import OpenPNM
pn = OpenPNM.Network.Cubic(name='net',shape=[10,10,10])
pn.add_boundaries()
```

In the last call, the method ``add_boundaries()`` is called, which means that a layer of boundary pores around the network is generated. These boundary pores will be used in the following calculations. Next we generate a geometry for the network and the phase, in this case air. A geometry can span over a part of the network only, so we need to specify to which pores and throats this geometry object should apply. For this example, we want it to apply to all pores and throats of the network. To do so, we can get all pore and throat indices of the network with the ``pn.pores()`` and ``pn.throats()`` calls, and pass these to the geometry object.

.. code-block:: python
``` python
>>> Ps = pn.pores()
>>> Ts = pn.throats()
>>> geo = OpenPNM.Geometry.Stick_and_Ball(network=pn, pores=Ps, throats=Ts)

Ps = pn.pores()
Ts = pn.throats()
geo = OpenPNM.Geometry.Stick_and_Ball(network=pn,name='basic',pores=Ps,throats=Ts)
```

Next we want to change the pore and throat sizes and throat lengths that were generated by the creation of the ``geo`` object. A change like this requires a regeneration of the geometry object, which is accomplished by calling the ``geo.regenerate()`` method.

.. code-block:: python
``` python
>>> geo['throat.diameter'] = geo['throat.diameter']*1e-5
>>> geo['pore.diameter'] = geo['pore.diameter']*1e-5
>>> geo['throat.length'] = geo['throat.length']*1e-5
>>> geo.regenerate()

geo['throat.diameter'] = geo['throat.diameter']*1e-5
geo['pore.diameter'] = geo['pore.diameter']*1e-5
geo['throat.length'] = geo['throat.length']*1e-5
geo.regenerate()
```

Then the phase object is created and a custom value is set.

.. code-block:: python
``` python
>>> air = OpenPNM.Phases.Air(network=pn)
>>> air['pore.Dac'] = 1e-7 # Add custom properties directly

air = OpenPNM.Phases.Air(network=pn)
air['pore.Dac'] = 1e-7 # Add custom properties directly
```

In the next step, a physics object is instantiated. A physics objects can span over several geometries, so we need to specify again to which pores and throats it should apply.

.. code-block:: python
``` python
>>> phys = OpenPNM.Physics.Standard(network=pn,phase=air,geometry=geo)

```

phys = OpenPNM.Physics.Standard(network=pn,phase=air,geometry=geo)
Now we tell the physics object to use the 'bulk_diffusion' model from the 'diffusive_conductance' library.

Now we tell the physics object to use the 'bulk_diffusion' of the 'diffusive_conductance' model.
``` python
>>> mode = OpenPNM.Physics.models.diffusive_conductance.bulk_diffusion
>>> phys.add_model(model= mod,
... propname='throat.gdiff_ac',
... pore_diffusivity='pore.Dac')

.. code-block:: python
```

phys.add_model(model=OpenPNM.Physics.models.diffusive_conductance.bulk_diffusion,
propname='throat.gdiff_ac',
pore_diffusivity='pore.Dac')
## Generate an Algorithm Object

+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Generate an Algorithm Object
+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
All algorithms in OpenPNM are independent objects. The Fickian diffusion algorithm is instantiated as follows:

.. code-block:: python
``` python
>>> alg = OpenPNM.Algorithms.FickianDiffusion(network=pn, phase=air)

alg = OpenPNM.Algorithms.FickianDiffusion(network=pn,phase=air)
```

-------------------------------------------------------------------------------
Apply Dirichlet Conditions to Two Faces
-------------------------------------------------------------------------------
### Apply Dirichlet Conditions to Two Faces

Now this algorithm needs to know about the boundary conditions which are to be applied. Let's start by defining Dirichlet conditions on two opposite faces. This is done by first finding the pore indices that correspond to the two faces. The generation of cubic networks automatically adds pores to the network with the label of the different faces. Let's use the 'top_boundary' and 'bottom_boundary' for this and apply Dirichlet boundary conditions to both and apply a numerical value for the boundary conditions:

.. code-block:: python
``` python
>>> BC1_pores = pn.pores(labels=['top_boundary'])
>>> alg.set_boundary_conditions(bctype='Dirichlet', bcvalue=0.6, pores=BC1_pores)
>>> BC2_pores = pn.pores(labels=['bottom_boundary'])
>>> alg.set_boundary_conditions(bctype='Dirichlet', bcvalue=0.4, pores=BC2_pores)

BC1_pores = pn.pores(labels=['top_boundary'])
alg.set_boundary_conditions(bctype='Dirichlet', bcvalue=0.6, pores=BC1_pores)
BC2_pores = pn.pores(labels=['bottom_boundary'])
alg.set_boundary_conditions(bctype='Dirichlet', bcvalue=0.4, pores=BC2_pores)
```

The above code adds the Dirichlet boundary conditions to both the pores at the 'top_boundary' and the 'bottom_boundary'. The Fickian algorithm looks for this specific label when analyzing and setting up the problem. Note that the boundary condition labels will only be applied to this specific algorithm. This is designed to allow multiple algorithms to exist simultaneously without interfering with each other.
The above code adds the Dirichlet boundary conditions to both the pores at the ``'top_boundary'`` and the ``'bottom_boundary'``. The **Fickian** algorithm looks for this specific label when analyzing and setting up the problem. Note that the boundary condition labels will only be applied to this specific algorithm. This is designed to allow multiple algorithms to exist simultaneously without interfering with each other.

Once the boundary conditions are specified, the algorithm can be run quite simply as:

.. code-block:: python

alg.run(conductance='throat.diffusive_conductance')
``` python
>>> alg.run(conductance='throat.diffusive_conductance')

```

This runs the algorithm using 'throat.diffusive_conductance' as the conductance. Using this argument will be important if another diffusive conductance is defined for the phase. The results are then stored on the Algorithm object. This is done to prevent simultaneous objects from interfering with each other. If and when the results of an Algorithm are required by the network model they must be explicitly sent *out* using:
This runs the algorithm using ``throat.diffusive_conductance`` as the conductance. Using this argument will be important if another diffusive conductance is defined for the phase. The results are then stored on the Algorithm object. This is done to prevent simultaneous objects from interfering with each other. If and when the results of an Algorithm are required by the network model they must be explicitly sent *out* using:

.. code-block:: python
``` python
>>> alg.return_results()

alg.return_results()
```

Each Algorithm must subclass the `return_results()` method so that it sends the correct information out the network and/or phase. In the case of the Fickian Algorithm, the 'mole_fraction' of the phase is stored on the Phase object in question. Running a different version of the Algorithm and calling `return_results()` will overwrite any previous values. The results of this simulation should produce the following visualization (done in Paraview):
Each Algorithm must subclass the `return_results()` method so that it sends the correct information out the network and/or phase. In the case of the Fickian Algorithm, the 'mole_fraction' are stored on the Phase object in question. Running a different version of the Algorithm and calling `return_results()` will overwrite any previous values. The results of this simulation should produce the following visualization (done in Paraview):

.. image:: http://imgur.com/Ae9cG0D.png
![](http://imgur.com/Ae9cG0D.png)

-------------------------------------------------------------------------------
Apply Neumann Conditions to a Group of Internal Pores
-------------------------------------------------------------------------------
### Apply Neumann Conditions to a Group of Internal Pores

The code below sets the total rate leaving a group of pores cumulatively. Note that the same Algorithm object is used (`alg`), so the Dirichlet boundary conditions applied in the previous step still exist. The lines below define a group of 10 pores which are generating mass at a set rate, which is accomplished by creating a 'Neumann_group' boundary condition and placing the numerical value of the rate in 'bcvalue'.
The code below sets the total rate leaving a group of pores cumulatively. Note that the same Algorithm object is used (`alg`), so the Dirichlet boundary conditions applied in the previous step still exist. The lines below define a group of 10 pores which are generating mass at a set rate, which is accomplished by setting a ``'Neumann_group'`` boundary condition and placing the numerical value of the rate in ``'bcvalue'``.

.. code-block:: python
``` python
>>> BC3_pores = [50,51,52,53,54,40,41,42,43,44]
>>> alg.set_boundary_conditions(bctype='Neumann_group', bcvalue=-5e-3, pores=BC3_pores)
>>> alg.run(conductance='throat.diffusive_conductance')
>>> alg.return_results()

BC3_pores = [50,51,52,53,54,40,41,42,43,44]
alg.set_boundary_conditions(bctype='Neumann_group', bcvalue=-5e-3, pores=BC3_pores)
alg.run(conductance='throat.diffusive_conductance')
alg.return_results()
```

This results in the image below, where a region of high concentration can be seen in the core of the domain due to the mass production:

.. image:: http://imgur.com/px45ANu.png
![](http://imgur.com/px45ANu.png)

-------------------------------------------------------------------------------
Apply Neumann Conditions in Several Pores Individually
-------------------------------------------------------------------------------
### Apply Neumann Conditions in Several Pores Individually

One of the options for specifying Neumann conditions is to apply the same rate to multiple pores. Begin by removing some of the conditions applied above, then set a few pores on the 'bottom' face to each have the same specific rate.

.. code-block:: python
``` python
>>> alg.set_boundary_conditions(bctype='Neumann_group', pores=BC3_pores, mode='remove') # This removes label from pores
>>> alg.set_boundary_conditions(bctype='Dirichlet',pores=BC2_pores, mode='remove')
>>> alg.set_boundary_conditions(bctype='Neumann',pores=BC2_pores, bcvalue=1e-10)
>>> alg.run(conductance='throat.diffusive_conductance')
>>> alg.return_results()

alg.set_boundary_conditions(bctype='Neumann_group', pores=BC3_pores, mode='remove') # This removes label from pores
alg.set_boundary_conditions(bctype='Dirichlet',pores=BC2_pores, mode='remove')
alg.set_boundary_conditions(bctype='Neumann',pores=BC2_pores, bcvalue=1e-10)
alg.run(conductance='throat.diffusive_conductance')
alg.return_results()
```

This results in image below. Notice that the concentration on the inlet face is not uniform, and that the smaller pores have a somewhat higher concentration (darker red), which is necessary if their flux is the be the same as larger, more conductive pores.

.. image:: http://imgur.com/50hGves.png
![](http://imgur.com/50hGves.png)
6 changes: 3 additions & 3 deletions Simulations/simulating_capillary_pressure_curves.md
Original file line number Diff line number Diff line change
Expand Up @@ -54,7 +54,7 @@ Then add the desired models to this object using:

This means that ``phys`` will now have a model called 'capillary_pressure', that when called will calculate throat entry pressures using the 'washburn' model. The Washburn model requires that the **Phase** object it's attached (mercury in this case) has the necessary physical properties of surface tension and contact angle.

| **Note** : Both surface tension and contact angle are actually 'phase system' properties, rather than solely water properties. It is an open problem in OpenPNM to figure out how to treat these sort of properties more rigorously. For the present time, they must be entered a single phase properties, in this case in ``Hg``.
> **Note** : Both surface tension and contact angle are actually 'phase system' properties, rather than solely water properties. It is an open problem in OpenPNM to figure out how to treat these sort of properties more rigorously. For the present time, they must be entered a single phase properties, in this case in ``Hg``.
## Run a Drainage Simulation

Expand Down Expand Up @@ -83,7 +83,7 @@ Next, we specify through which pores mercury enters the network.

In MIP experiment this is all pores on the surface which are actually the "boundary" pores we created above. These are found using the *wildcard* operator with the ``'boundary'`` label.

| **Advanced Features**: It is also possible to call the ``set_outlets`` method to specify through pores which the defending phase exits the network, but this is not relevant to an MIP simulation since the sample is evacuated of air. Moreover, it's possible to simulate secondary drainage by using the ``set_residual`` method, but this requires knowing the locations of the residual non-wetting phase from some other simulation.
> **Advanced Features**: It is also possible to call the ``set_outlets`` method to specify through pores which the defending phase exits the network, but this is not relevant to an MIP simulation since the sample is evacuated of air. Moreover, it's possible to simulate secondary drainage by using the ``set_residual`` method, but this requires knowing the locations of the residual non-wetting phase from some other simulation.
Now the ``MIP`` algorithm is ready to ``run``:

Expand All @@ -106,4 +106,4 @@ This algorithm produces 4 data arrays and stores them on the ``MIP`` object. ``

It is possible using the information stored on the ``MIP`` object to reproduce the capillary pressure curve manually. Since this is such a common operation, however, the **Drainage** class has methods already available for doing this. The raw capillary pressure curve data can be obtained using the ``get_drainage_data`` method, which returns a table of data for plotting in an external program of your choice. Alternatively, a plot can be created directly with ``MIP.plot_drainage_curve()``.

.. image:: http://i.imgur.com/ZxuCict.png
![](http://i.imgur.com/ZxuCict.png)

0 comments on commit 729e9a8

Please sign in to comment.