OpenPNM forces manually-passed conns to be in ascending order

[GeorgyBorisochev]

Just discovered a weird issue/feature:

I am generating Conns and Coords separatey as Numpy arrays and then passing them to GenericNetwork to work in OpenPNM.
GenericNetwork seems to swap positions for Conns values so that first entry is always smaller than second entry.

For example, if my generator has a certain throat connecting pores [385 291],
throat.conns will say [291 385].

As far as I can see that does not affect any results downstream, but I am just curious as to why it works that way?

Thanks,
Georgy

[ma-sadeghi]

Hi Georgy,

As you have guessed, that's a feature :) Since in openpnm (at least as of now), we don't support bidirectional throats, sorting network conns in an ascending order ensures that the adjacency matrix is always upper triangular, which can later be leveraged to speed up certain topological queries/manipulations.

Finally, like you said, that doesn't (shouldn't) affect any results downstream, it's just more optimized this way.

Best,
Amin

[GeorgyBorisochev]

Hi Amin,

Thanks for your reply - that makes sense!
I have a follow-up question: been digging in the guts of the ReactiveTransport solver and the way it creates A and b and solves it, and as far as I can see it uses full matrix to pass to solver (pypardiso by default).
You mentioned using upper triangular adjacency matrix, am I right to assume that it is only applicable to topological manipulations and I have to use full A matrix for actual solution?

Georgy

[ma-sadeghi]

Just a quick heads up: we've extensively cleaned up/refactored algorithm objects for v3, to be released soon. So, it might be good idea to also take a look at what we have in the dev branch (if you haven't already).

Regarding your question: the upper triangular trick only works for network topology. It doesn't necessarily work for algorithms. Specifically, the A matrix is not necessarily symmetric for all problems. For instance, for pure diffusion problems (may it be heat transfer, mass transfer, charge transfer, etc.), the A matrix is indeed symmetric, but for problems with advection, your A matrix is asymmetric.

Furthermore, for those problems that "do" produce a symmetric A matrix, we won't lose any performance as most linear solvers (including those that are wrapped in openpnm such as pardiso, umfpack, etc.) already apply different routines for different matrix structures.