Iterative solving of sparse systems of linear equations with (M, N) right-hand size matrix

Question

I would like to solve a sparse linear equations system: A x = b, where A is a (M x M) array, b is an (M x N) array and x is and (M x N) array.

I solve this in three ways using the:

• scipy.linalg.solve(A.toarray(), b.toarray()),
• scipy.sparse.linalg.spsolve(A, b),
• scipy.sparse.linalg.splu(A).solve(b.toarray()) # returns a dense array

I wish to solve the problem using the iterative scipy.sparse.linalg methods:

• scipy.sparse.linalg.cg,
• scipy.sparse.linalg.bicg,
• ...

However, the metods suport only a right hand side b with a shape (M,) or (M, 1). Any ideas on how to expand these methods to (M x N) array b?

Answer 1

A key difference between iterative solvers and direct solvers is that direct solvers can more efficiently solve for multiple right-hand values by using a factorization (usually either Cholesky or LU), while iterative solvers can't. This means that for direct solvers there is a computational advantage to solving for multiple columns simultaneously.

For iterative solvers, on the other hand, there's no computational gain to be had in simultaneously solving multiple columns, and this is probably why matrix solutions are not supported natively in the API of cg, bicg, etc.

Because of this, a direct solution like scipy.sparse.linalg.spsolve will probably be optimal for your case. If for some reason you still desire an iterative solution, I'd just create a simple convenience function like this:

from scipy.sparse.linalg import bicg

def bicg_solve(M, B):
    X, info = zip(*(bicg(M, b) for b in B.T))
    return np.transpose(X), info

Then you can create some data and call it as follows:

import numpy as np
from scipy.sparse import csc_matrix

# create some matrices
M = csc_matrix(np.random.rand(5, 5))
B = np.random.rand(5, 4)

X, info = bicg_solve(M, B)
print(X.shape)
# (5, 4)

Any iterative solver API which accepts a matrix on the right-hand-side will essentially just be a wrapper for something like this.