Python left multiplication of of a matrix with inverse of a sparse matrix
I'm trying to calculate an expression of the form K = P*C.T*S^-1 (implementation of a Kalman filter)
All the involved matrices are sparse, and I would of course like to avoid calculating the actual inverse.
I tried using
import scipy.sparse.linalg as spln
self.K = self.P.dot(spln.spsolve(S.T, C).T)
The problem is that spsolve expects it's second argument to be a vector and not a matrix.
edit: Clarificat开发者_运维知识库ion, the problem could in Matlab be solved by K = P * (C / S), so what I'm looking for is a method similar to spsolve but which can accept a matrix as its second argument. This could of course be done by splitting C into a number of column vectors c1..cn and solving the problem for each of them and then reassembling them into a matrix, but I suspect doing that will be both cumbersome and inefficient.
edit2&3: The dimensions of the matrices will typically be around P~10⁶x10^6, S~100x100, C=100x10⁶. P diagonal and S symmetric and C will only have one element per row. It will be used for an implementation of a Kalman filter using sparse matrices, see
http://en.wikipedia.org/wiki/Kalman_filter#The_Kalman_filter
As a workaround can do
import numpy as np
from scipy.sparse.linalg import splu
def spsolve2(a, b):
a_lu = splu(a)
out = np.empty((A.shape[1], b.shape[1]))
for j in xrange(b.shape[1]):
out[:,j] = a_lu.solve(b[:,j])
return out
self.K = self.P.dot(spsolve2(S.T, C).T)
But yes, it's a bug that spsolve
does not accept matrices.
However, as your S
is not very large, you can as well use a dense inverse.
精彩评论