Traceback when getting flow rate of a network with a single throat

[Arenhart]

Running the following line:
perm.rate(throats=perm.network.throats("all"), mode="individual")
(where perm is a openpnm.algorithms.StokesFlow of a network with a single throat)
throws the following error:

  File "C:\Work\lib\Python\Lib\site-packages\openpnm\algorithms\_transport.py", line 318, in rate
    R = np.absolute(Qt[throats])
IndexError: too many indices for array: array is 0-dimensional, but 1 were indexed

This seems to be caused by the following line squeezing a ndarray with more than one dimension, but a single value, into a zero-dimension array when it has a single value, while the rest of the code still treats the variable as a one dimension array:
https://github.com/PMEAL/OpenPNM/blob/0efb343f2f564b4cf0189828d49096da1addddda/openpnm/algorithms/_transport.py#L315C46-L315C46
Qt = np.diff(g*X12, axis=1).squeeze()

I fixed the problem locally with the following code below the previous line:

        if Qt.ndim == 0:
            Qt = Qt[np.newaxis]

Is this solution acceptable? Can I make a pull request to fix it in the repository?

[ma-sadeghi]

Hi. Can you share an MWE? (minimum working example that reproduces the bug)

[Arenhart]

Of course.

import openpnm as op
import numpy as np

Nx, Ny, Nz = 1, 1, 2
Lc = 1e-4
pn = op.network.Cubic([Nx, Ny, Nz], spacing=Lc)
pn.add_model_collection(op.models.collections.geometry.spheres_and_cylinders)
pn.regenerate_models()

water = op.phase.Water(network=pn)
water.add_model_collection(op.models.collections.physics.standard)
water["throat.hydraulic_conductance"]=[1,]
perm = op.algorithms.StokesFlow(network=pn, phase=water)

perm.set_value_BC(
        pores=pn.pores("zmin"), 
        values=101325, 
        mode='overwrite'
    )  # pressure in pa: 101325 pa = 1 atm
perm.set_value_BC(
        pores=pn.pores("zmax"), 
        values=0, 
        mode='overwrite'
    )
perm.run()
perm.rate(throats=perm.network.throats("all"), mode="individual")

Running this will throw the following error:

---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
Cell In[47], line 12
      6 perm.set_value_BC(
      7         pores=pn.pores("zmax"), 
      8         values=0, 
      9         mode='overwrite'
     10     )
     11 perm.run()
---> 12 perm.rate(throats=perm.network.throats("all"), mode="individual")

File D:\anaconda3\envs\notebook\lib\site-packages\openpnm\algorithms\_transport.py:318, in Transport.rate(self, pores, throats, mode)
    315 Qt = np.diff(g*X12, axis=1).squeeze()
    317 if throats.size:
--> 318     R = np.absolute(Qt[throats])
    319     if mode == 'group':
    320         R = np.sum(R)

IndexError: too many indices for array: array is 0-dimensional, but 1 were indexed

The error will not occur with larger networks, (increasing Nz to 3 already avoids it). It is definitely a consequence of np.squeeze eliminating all axes with size 1, which is the case for all axes in an array with a single element.

[Arenhart]

Addendum: The bug can alternativelly be solved by replacing line 315 Qt = np.diff(g*X12, axis=1).squeeze() with Qt = np.diff(g*X12, axis=1)[:, 0], which is more compact than the previous suggestion.

[ma-sadeghi]

@Arenhart Thanks for reporting this. It's indeed a bug. Your solution works just fine, but let's use .ravel() instead of [:, 0] (just for the sake of being more Pythonic). We'd really appreciate a PR! Thanks once again!

[Arenhart]

@ma-sadeghi Done
#2798
Thank you for accepting my contribution, and suggesting the more appropriate correction.

[ma-sadeghi]

Thank "you" for contributing to OpenPNM!

If you were interested in further developing OpenPNM, we have a long list of issues and we welcome/appreciate PRs 😊