Scipy Sparse: Singular Matrix Warning after SciPy/NumPy Update

Question

my problem results from nodal analysis for a large system of resistors. I am basically setting up a large sparse matrix A, my solution vector b, and I'm trying to solve the linear equation A * x = b. In order to do so, I'm using the scipy.sparse.linalg.spsolve method.

Until recently, everything worked fine until I upgraded SciPy from v0.13.3 to v0.19.1 (which also included a NumPy Upgrade to v1.13.1). I'm running Python 2.7.6. When using the same code as before the update, I get errors, especially for systems that produce matrices > 10000 x 10000 . The warnings are:

SparseEfficiencyWarning: splu requires CSC matrix format
  warn('splu requires CSC matrix format', SparseEfficiencyWarning)
MatrixRankWarning: Matrix is exactly singular
  warn("Matrix is exactly singular", MatrixRankWarning)

spsolve is then - sometimes - unable to find a solution.

As I am performing nodal analysis, a singular matrix is expected since the position of the ground potential is generally not well-defined. However, before the update, a solution was found in 99% of the cases, maybe more. Now, I'm at 10% for large systems at best. I have not changed the algorithm and for a few tests, I have used identical code as before. Here is how I set up my calculation:

1. I generate a random three-dimensional network of resistors (I realize that I could accidentally create unsolvable networks but the percentages above should not change that drastically). The only SciPy/NumPy functions used here is np.random
2. I create a sparse lil-matrix which I fill with conductance values extracted from my resistor network. I also create a solution vector which is not sparse.
3. I convert the conductance matrix to csr-format and use the spsolve method. This is where my code lately fails.

Could it be the method that has changed?

Is spsolve maybe even inappropriate? The matrices I create are generally symmetric and in a block tridiagonal form. Is there a more efficient way to solve the linear equation than spsolve?

Every sort of help is greatly appreciated! Thanks for reading.

Here is how my matrices look like in 'spy'-representation

Answer 1

Your previous scipy-version is quite old and it used umfpack for this task.

Due to licensing-issues (GPL is not compatible with scipy and i think umfpack switched the license at some point), this lib was removed and now superlu is used. Many people observed slow-downs (and robustness issues), but evaluating performance might not be that easy (superlu can be fast and robust too).

Read also this.

You probably got two options:

• Tune superlu's parameters (read the official superlu documentation and scipy's docs on how to pass these options)
  • Pivoting and Ordering is very important!
  • There is also a symmetric-mode (not really a highly-tuned symmetric-matrix solver, but possibly better pivoting-rules)
  • Maybe iterative-refinement can help too (unsure!)
• If the license-stuff is not a problem for you: use scikit-umfpack to make scipy use umfpack again!

If your matrix is PSD, cholmod, available within scikit-sparse (currently unmaintained!) might be the lib to use (again: license)!