Repository: PMEAL/OpenPNM
Pull request #1330 (2019-09-26 01:29:13)

Preventing object storage in settings dict, and added test
This PR is on hold until some changes to the PNP algorithm are merged...


Patch:
diff --git a/example_PNP.py b/example_PNP.py
deleted file mode 100644
index 1efb4ba3b..000000000
--- a/example_PNP.py
+++ /dev/null
@@ -1,106 +0,0 @@
-import openpnm as op
-from openpnm.phases import mixtures
-import numpy as np
-ws = op.Workspace()
-proj = ws.new_project()
-# ws.settings['loglevel'] = 20
-
-
-# network, geometry, phase
-np.random.seed(0)
-
-net = op.network.Cubic(shape=[8, 8, 1], spacing=9e-4)
-prs = (net['pore.back'] * net['pore.right'] + net['pore.back'] *
-       net['pore.left'] + net['pore.front'] * net['pore.right'] +
-       net['pore.front'] * net['pore.left'])
-thrts = net['throat.surface']
-op.topotools.trim(network=net, pores=net.Ps[prs], throats=net.Ts[thrts])
-
-
-geo = op.geometry.StickAndBall(network=net, pores=net.Ps, throats=net.Ts)
-pore_d = op.models.misc.constant
-throat_d = op.models.misc.constant
-geo.add_model(propname='pore.diameter', model=pore_d, value=1.5e-4)
-geo.add_model(propname='throat.diameter', model=throat_d, value=1e-4)
-geo.regenerate_models()
-
-sw = mixtures.SalineWater(network=net)
-# Retrieve handles to each species for use below
-Cl, Na, H2O = sw.components.values()
-
-# physics
-phys = op.physics.GenericPhysics(network=net, phase=sw, geometry=geo)
-
-flow = op.models.physics.hydraulic_conductance.hagen_poiseuille
-phys.add_model(propname='throat.hydraulic_conductance',
-               pore_viscosity='pore.viscosity',
-               throat_viscosity='throat.viscosity',
-               model=flow, regen_mode='normal')
-
-current = op.models.physics.ionic_conductance.laplace
-phys.add_model(propname='throat.ionic_conductance',
-               model=current, regen_mode='normal')
-
-eA_dif = op.models.physics.diffusive_conductance.ordinary_diffusion
-phys.add_model(propname='throat.diffusive_conductance.' + Na.name,
-               pore_diffusivity='pore.diffusivity.' + Na.name,
-               throat_diffusivity='throat.diffusivity.' + Na.name,
-               model=eA_dif, regen_mode='normal')
-
-eB_dif = op.models.physics.diffusive_conductance.ordinary_diffusion
-phys.add_model(propname='throat.diffusive_conductance.' + Cl.name,
-               pore_diffusivity='pore.diffusivity.' + Cl.name,
-               throat_diffusivity='throat.diffusivity.' + Cl.name,
-               model=eB_dif, regen_mode='normal')
-
-# algorithms
-sf = op.algorithms.StokesFlow(network=net, phase=sw)
-sf.set_value_BC(pores=net.pores('back'), values=0.01)
-sf.set_value_BC(pores=net.pores('front'), values=0.00)
-sf.settings['rxn_tolerance'] = 1e-12
-sf.run()
-sw.update(sf.results())
-
-p = op.algorithms.OhmicConduction(network=net, phase=sw)
-p.settings['conductance'] = 'throat.ionic_conductance'
-p.settings['quantity'] = 'pore.potential'
-p.set_value_BC(pores=net.pores('left'), values=0.01)
-p.set_value_BC(pores=net.pores('right'), values=0.00)
-p.settings['rxn_tolerance'] = 1e-12
-p.run()
-sw.update(p.results())
-
-eA = op.algorithms.NernstPlanck(network=net, phase=sw, ion=Na.name)
-eA.set_value_BC(pores=net.pores('back'), values=100)
-eA.set_value_BC(pores=net.pores('front'), values=90)
-eA.settings['rxn_tolerance'] = 1e-12
-
-eB = op.algorithms.NernstPlanck(network=net, phase=sw, ion=Cl.name)
-eB.set_value_BC(pores=net.pores('back'), values=100)
-eB.set_value_BC(pores=net.pores('front'), values=90)
-eB.settings['rxn_tolerance'] = 1e-12
-
-ad_dif_mig_Na = op.models.physics.ad_dif_mig_conductance.ad_dif_mig
-phys.add_model(propname='throat.ad_dif_mig_conductance.' + Na.name,
-               model=ad_dif_mig_Na, ion=Na.name,
-               s_scheme='exponential')
-
-ad_dif_mig_Cl = op.models.physics.ad_dif_mig_conductance.ad_dif_mig
-phys.add_model(propname='throat.ad_dif_mig_conductance.' + Cl.name,
-               pore_pressure=sf.settings['quantity'],
-               model=ad_dif_mig_Cl, ion=Cl.name,
-               s_scheme='exponential')
-
-pnp = op.algorithms.ChargeConservationNernstPlanck(network=net, phase=sw)
-pnp.setup(potential_field=p, ions=[eA, eB])
-pnp.settings['max_iter'] = 10
-pnp.settings['tolerance'] = 1e-04
-pnp.settings['charge_conservation'] = 'laplace'
-# Electroneutrality condition does not work with new Mixtures
-# pnp.settings['charge_conservation'] = 'electroneutrality'
-pnp.run()
-
-sw.update(sf.results())
-sw.update(p.results())
-sw.update(eA.results())
-sw.update(eB.results())
diff --git a/openpnm/algorithms/ChargeConservationNernstPlanck.py b/openpnm/algorithms/ChargeConservationNernstPlanck.py
deleted file mode 100644
index faade5077..000000000
--- a/openpnm/algorithms/ChargeConservationNernstPlanck.py
+++ /dev/null
@@ -1,116 +0,0 @@
-import numpy as np
-from openpnm.algorithms import ReactiveTransport
-from openpnm.utils import logging
-from openpnm.models.physics import generic_source_term as gst
-logger = logging.getLogger(__name__)
-
-
-class ChargeConservationNernstPlanck(ReactiveTransport):
-    r"""
-    A subclass of GenericTransport to solve the charge conservation and
-    Nernst-Planck equations.
-    """
-    def __init__(self, settings={}, phase=None, **kwargs):
-        def_set = {'phase': None,
-                   'potential_field': None,
-                   'ions': [],
-                   'charge_conservation': 'electroneutrality',
-                   'tolerance': 1e-4,
-                   'max_iter': 10}
-        super().__init__(**kwargs)
-        self.settings.update(def_set)
-        self.settings.update(settings)
-        if phase is not None:
-            self.setup(phase=phase)
-
-    def setup(self, phase=None, potential_field=None, ions=[],
-              charge_conservation=None, tolerance=None,
-              max_iter=None, **kwargs):
-        r"""
-        """
-        if phase:
-            self.settings['phase'] = phase.name
-        if potential_field:
-            self.settings['potential_field'] = potential_field
-        if ions:
-            self.settings['ions'] = ions
-        if charge_conservation:
-            self.settings['charge_conservation'] = charge_conservation
-        if tolerance:
-            self.settings['tolerance'] = tolerance
-        if max_iter:
-            self.settings['max_iter'] = max_iter
-        super().setup(**kwargs)
-
-    def run(self):
-        # Phase, potential and ions algorithms
-        phase = self.project.phases()[self.settings['phase']]
-        p_alg = self.settings['potential_field']
-        e_alg = self.settings['ions']
-
-        # Define initial conditions (if not defined by the user)
-        try:
-            p_alg[p_alg.settings['quantity']]
-        except KeyError:
-            try:
-                p_alg[p_alg.settings['quantity']] = (
-                    phase[p_alg.settings['quantity']])
-            except KeyError:
-                p_alg[p_alg.settings['quantity']] = np.zeros(
-                    shape=[p_alg.Np, ], dtype=float)
-        for e in e_alg:
-            try:
-                e[e.settings['quantity']]
-            except KeyError:
-                try:
-                    e[e.settings['quantity']] = phase[e.settings['quantity']]
-                except KeyError:
-                    e[e.settings['quantity']] = np.zeros(shape=[e.Np, ],
-                                                         dtype=float)
-
-        # Source term for Poisson or charge conservation (electroneutrality) eq
-        Ps = (p_alg['pore.all'] * np.isnan(p_alg['pore.bc_value']) *
-              np.isnan(p_alg['pore.bc_rate']))
-        mod = gst.charge_conservation
-        phys = p_alg.project.find_physics(phase=phase)
-        phys[0].add_model(propname='pore.charge_conservation', model=mod,
-                          phase=phase, p_alg=p_alg, e_alg=e_alg,
-                          assumption=self.settings['charge_conservation'])
-        p_alg.set_source(propname='pore.charge_conservation', pores=Ps)
-
-        # Define tolerance and initialize residuals
-        tol = self.settings['tolerance']
-        res = {}
-        res['potential'] = 1e+06
-        for e in e_alg:
-            res[e.name] = 1e+06
-
-        # Iterate until solutions converge
-        for itr in range(int(self.settings['max_iter'])):
-            r = str([float(format(i, '.3g')) for i in res.values()])[1:-1]
-            logger.info('Iter: ' + str(itr+1) + ', Residuals: ' + r)
-            convergence = max(i for i in res.values()) < tol
-            if not convergence:
-                # Poisson eq
-                phys[0].regenerate_models()
-                phi_old = p_alg[p_alg.settings['quantity']].copy()
-                p_alg._run_reactive(x=phi_old)
-                phi_new = p_alg[p_alg.settings['quantity']].copy()
-                # Residual
-                res['potential'] = np.sum(np.absolute(phi_old**2 - phi_new**2))
-                # Update phase and physics
-                phase.update(p_alg.results())
-                phys[0].regenerate_models()
-
-                # Ions
-                for e in e_alg:
-                    c_old = e[e.settings['quantity']].copy()
-                    e._run_reactive(x=c_old)
-                    c_new = e[e.settings['quantity']].copy()
-                    # Residual
-                    res[e.name] = np.sum(np.absolute(c_old**2 - c_new**2))
-                    phase.update(e.results())
-
-            if convergence:
-                logger.info('Solution converged: ' + str(res))
-                break
diff --git a/openpnm/algorithms/__init__.py b/openpnm/algorithms/__init__.py
index 95f4a081d..c1ea83e34 100644
--- a/openpnm/algorithms/__init__.py
+++ b/openpnm/algorithms/__init__.py
@@ -26,5 +26,4 @@ from .InvasionPercolation import InvasionPercolation
 from .MixedInvasionPercolation import MixedInvasionPercolation
 from .Porosimetry import Porosimetry
 from .NernstPlanck import NernstPlanck
-from .ChargeConservationNernstPlanck import ChargeConservationNernstPlanck
 from . import metrics
diff --git a/openpnm/utils/misc.py b/openpnm/utils/misc.py
index df79f7edb..94d97e989 100644
--- a/openpnm/utils/misc.py
+++ b/openpnm/utils/misc.py
@@ -102,6 +102,12 @@ class SettingsDict(PrintableDict):
     None
 
     """
+    def __setitem__(self, key, value):
+        if hasattr(value, 'Np'):
+            raise Exception('Cannot store OpenPNM objects in settings, ' +
+                            'store object\'s name instead')
+        super().__setitem__(key, value)
+
     def __missing__(self, key):
         self[key] = None
         return self[key]


Test patch:
diff --git a/tests/unit/utils/UtilsTest.py b/tests/unit/utils/UtilsTest.py
index db6d279c8..13486333c 100644
--- a/tests/unit/utils/UtilsTest.py
+++ b/tests/unit/utils/UtilsTest.py
@@ -46,6 +46,10 @@ class UtilsTest:
         r = D.__repr__()
         assert r == "{'item1': 1, 'item2': 2, 'item3': array([1, 2])}"
 
+    def test_settings_dict(self):
+        with pytest.raises(Exception):
+            self.geo.settings['test'] = self.net
+
 
 if __name__ == '__main__':