Description
I am trying to solve several independent systems of equations at the same time using numpy.linalg.solve, e.g., a has shape (N, M, M). The challenge I'm running into is how to deal with the case when, for some values of N, the last two dimensions comprise a singular matrix.
Consider this example
import numpy as np
arr = [[[1, 0], [0, 1]], [[1, 0], [0, 1]], [[0, 0], [0, 0]]]
np.linalg.solve(arr, 1)On numpy 1.9.2, this will raise a LinAlgError, but I still would like to know the solutions for all the full-rank matrices. (I would be okay with getting back NaNs for the singular cases.) I'm aware I can compute the singular values and use fancy indexing to slice the array, but in my algorithm I have to do several slicing/filtering steps and I would prefer not to lose the alignment of a with b -- my real-world a typically has 6-7 dimensions and all the indexing arrays make the code hard to follow.
My ideal solution seems to involve masked arrays since I end up doing several other filtering steps, but it appears that most (all?) functions in the linalg family ignore array masks.
Is there an alternative solution already available, or perhaps a suggestion on how I could contribute a solution that I could package as a PR?