diff --git a/openpnm/algorithms/_transport.py b/openpnm/algorithms/_transport.py index ac004e856b..2150252d13 100644 --- a/openpnm/algorithms/_transport.py +++ b/openpnm/algorithms/_transport.py @@ -311,6 +311,7 @@ def rate(self, pores=[], throats=[], mode='group'): if g.size == self.Nt: g = np.tile(g, (2, 1)).T # Make conductance an Nt by 2 matrix # The next line is critical for rates to be correct + # We could also do "g.T.flatten()" or "g.flatten('F')" g = np.flip(g, axis=1) Qt = np.diff(g*X12, axis=1).ravel() diff --git a/openpnm/models/geometry/pore_surface_area/_funcs.py b/openpnm/models/geometry/pore_surface_area/_funcs.py index 11e9612d47..0fc9846e66 100644 --- a/openpnm/models/geometry/pore_surface_area/_funcs.py +++ b/openpnm/models/geometry/pore_surface_area/_funcs.py @@ -40,7 +40,7 @@ def sphere( R = network[pore_diameter] / 2 value = 4 * _np.pi * R**2 Tca = network[throat_cross_sectional_area] - _np.subtract.at(value, network.conns.flatten(), _np.repeat(Tca, repeats=2)) + _np.subtract.at(value, network.conns.T.flatten(), _np.repeat(Tca, repeats=2)) return value @@ -69,7 +69,7 @@ def circle( """ value = _np.pi * network[pore_diameter] Tca = network[throat_cross_sectional_area] - _np.subtract.at(value, network.conns.flatten(), _np.repeat(Tca, repeats=2)) + _np.subtract.at(value, network.conns.T.flatten(), _np.repeat(Tca, repeats=2)) return value @@ -95,7 +95,7 @@ def cube( D = network[pore_diameter] value = 6.0 * D**2 Tca = network[throat_cross_sectional_area] - _np.subtract.at(value, network.conns.flatten(), _np.repeat(Tca, repeats=2)) + _np.subtract.at(value, network.conns.T.flatten(), _np.repeat(Tca, repeats=2)) return value @@ -121,5 +121,5 @@ def square( D = network[pore_diameter] value = 4.0 * D Tca = network[throat_cross_sectional_area] - _np.subtract.at(value, network.conns.flatten(), _np.repeat(Tca, repeats=2)) + _np.subtract.at(value, network.conns.T.flatten(), _np.repeat(Tca, repeats=2)) return value diff --git a/openpnm/network/_network.py b/openpnm/network/_network.py index 1d78a0b566..6bdf793ce1 100644 --- a/openpnm/network/_network.py +++ b/openpnm/network/_network.py @@ -359,15 +359,16 @@ def create_incidence_matrix(self, weights=None, fmt='coo', # Check if provided data is valid if weights is None: weights = np.ones((self.Nt,), dtype=int) - elif np.shape(weights)[0] != self.Nt: + if np.shape(weights)[0] == self.Nt: + weights = np.append(weights, weights) + elif np.shape(weights)[0] != 2*self.Nt: raise Exception('Received dataset of incorrect length') conn = self['throat.conns'] - row = conn[:, 0] - row = np.append(row, conn[:, 1]) + row = conn[:, 1] + row = np.append(row, conn[:, 0]) col = np.arange(self.Nt) col = np.append(col, col) - weights = np.append(weights, weights) temp = sprs.coo.coo_matrix((weights, (row, col)), (self.Np, self.Nt)) diff --git a/tests/unit/models/geometry/PoreSurfaceAreaTest.py b/tests/unit/models/geometry/PoreSurfaceAreaTest.py index d86c04d07b..fc4c52354d 100644 --- a/tests/unit/models/geometry/PoreSurfaceAreaTest.py +++ b/tests/unit/models/geometry/PoreSurfaceAreaTest.py @@ -6,64 +6,66 @@ class PoreSurfaceAreaTest: def setup_class(self): - self.net = op.network.Cubic(shape=[5, 5, 5], spacing=1.0) - self.net['pore.diameter'] = 1 - self.net['throat.cross_sectional_area'] = 0.1 + self.net3D = op.network.Cubic(shape=[5, 5, 5], spacing=1.0) + self.net3D['pore.diameter'] = 1 + self.net3D['throat.cross_sectional_area'] = 0.1 + self.net2D = op.network.Cubic(shape=[5, 5, 1], spacing=1.0) + self.net2D['pore.diameter'] = 1 + self.net2D['throat.cross_sectional_area'] = 0.1 def test_sphere(self): - self.net.add_model(propname='pore.surface_area', - model=mods.sphere, - regen_mode='normal') + self.net3D.add_model(propname='pore.surface_area', + model=mods.sphere, + regen_mode='normal') a = np.array([2.54159265, 2.64159265, 2.74159265, 2.84159265]) - b = np.unique(self.net['pore.surface_area']) + b = np.unique(self.net3D['pore.surface_area']) assert_allclose(a, b) def test_circle(self): - self.net.add_model(propname='pore.surface_area', + self.net2D.add_model(propname='pore.surface_area', model=mods.circle, regen_mode='normal') - a = np.array([2.54159265, 2.64159265, 2.74159265, 2.84159265]) - b = np.unique(self.net['pore.surface_area']) + a = np.array([2.74159265, 2.84159265, 2.94159265]) + b = np.unique(self.net2D['pore.surface_area']) assert_allclose(a, b) def test_cube(self): - self.net.add_model(propname='pore.surface_area', - model=mods.cube, - regen_mode='normal') + self.net3D.add_model(propname='pore.surface_area', + model=mods.cube, + regen_mode='normal') a = np.array([5.4, 5.5, 5.6, 5.7]) - b = np.unique(self.net['pore.surface_area']) + b = np.unique(self.net3D['pore.surface_area']) assert_allclose(a, b) def test_square(self): - self.net.add_model(propname='pore.surface_area', + self.net2D.add_model(propname='pore.surface_area', model=mods.square, regen_mode='normal') - a = np.array([3.4, 3.5, 3.6, 3.7]) - b = np.unique(self.net['pore.surface_area']) + a = np.array([3.6, 3.7, 3.8]) + b = np.unique(self.net2D['pore.surface_area']) assert_allclose(a, b) def test_circle_multi_geom(self): - self.net = op.network.Cubic(shape=[10, 1, 1], spacing=1.0) - self.net['pore.left'] = False - self.net['pore.left'][[0, 1, 2]] = True - self.net['pore.diameter'] = 1 - self.net['throat.cross_sectional_area'] = 0.1 - self.net['throat.cross_sectional_area'][[0, 1, 2]] = 0.3 - self.net['pore.right'] = ~self.net['pore.left'] - self.net.add_model(propname='pore.surface_area', - model=mods.circle, - domain='left', - regen_mode='normal') - self.net.add_model(propname='pore.surface_area', - model=mods.circle, - domain='right', - regen_mode='normal') - a = np.array([2.84159265, 2.54159265, 2.54159265]) - b = self.net['pore.surface_area@left'] - c = np.array([2.74159265, 2.94159265, 2.94159265, - 2.94159265, 2.94159265, 2.94159265, - 3.04159265]) - d = self.net['pore.surface_area@right'] + net = op.network.Cubic(shape=[10, 1, 1], spacing=1.0) + net['pore.left'] = False + net['pore.left'][[0, 1, 2]] = True + net['pore.diameter'] = 1 + net['throat.cross_sectional_area'] = 0.1 + net['throat.cross_sectional_area'][[0, 1, 2]] = 0.3 + net['pore.right'] = ~net['pore.left'] + net.add_model(propname='pore.surface_area', + model=mods.circle, + domain='left', + regen_mode='normal') + net.add_model(propname='pore.surface_area', + model=mods.circle, + domain='right', + regen_mode='normal') + a = np.array([2.84159265, 2.74159265, 2.74159265]) + b = net['pore.surface_area@left'] + c = np.array([2.74159265, 2.74159265, 2.74159265, 2.94159265, + 2.94159265, 2.94159265, 3.04159265]) + d = net['pore.surface_area@right'] assert_allclose(a, b) assert_allclose(c, d)