3

I'm writing a program with Armadillo C++ (4.400.1)

I have a matrix that has to be sparse and complex, and I want to calculate the inverse of such matrix. Since it is sparse it could be the pseudoinverse, but I can guarantee that the matrix has the full diagonal.

In the API documentation of Armadillo, it mentions the method .i() to calculate the inverse of any matrix, but sp_cx_mat members do not contain such method, and the inv() or pinv() functions cannot handle the sp_cx_mat type apparently.

sp_cx_mat Y;

/*Fill Y ensuring that the diagonal is full*/

sp_cx_mat Z = Y.i();

or

sp_cx_mat Z = inv(Y);

None of them work.

I would like to know how to compute the inverse of matrices of sp_cx_mat type.

Unapiedra
  • 15,037
  • 12
  • 64
  • 93
Santi Peñate-Vera
  • 1,053
  • 4
  • 33
  • 68
  • What are the dimensions? X * X or X * Y? – Surt Nov 01 '14 at 11:38
  • The matrix is square, so X by X – Santi Peñate-Vera Nov 01 '14 at 19:08
  • 2
    The inverse of a sparse matrix is not necessarily sparse. Do you really need the inverse of a sparse matrix? All literature I've seen after a quick search recommends to solve the underlying problem in a different way. I.e. use iterative methods (or others) to solve Ax = b. – Unapiedra Nov 05 '14 at 09:02
  • You could compute the SVD. Let `U S V == Y`, then `Y.pinv() == U S.pinv() V`. Because the Y is of full-rank (because of that diagonal) `Armadillo::svd_econ()` doesn't make sense. So, computing the SVD will be significantly slower than other methods for solving the pseodo-inverse. – Unapiedra Nov 05 '14 at 09:21
  • Yes, I need the inverse, because that inverted matrix will be used for the solution of (A^-1)x=b systems. Not calculating the inverted matrix would make me reformulate all the electrical models that fill the initial non-inverted matrix. – Santi Peñate-Vera Nov 05 '14 at 09:26
  • indeed Y.pinv() is nice, but it doesn't exist in the library for complex saparse matrices. – Santi Peñate-Vera Nov 05 '14 at 09:26
  • If you want to solve (A^-1)x=b, why don't you solve x=Ab instead? – dani Nov 05 '14 at 10:00
  • Because, some items in x and b are missing, therefore the system is solved using newton raphson methods. – Santi Peñate-Vera Nov 05 '14 at 10:51

1 Answers1

5

Sparse matrix support in Armadillo is not complete and many of the factorizations/complex operations that are available for dense matrices are not available for sparse matrices. There are a number of reasons for this, the largest being that efficient complex operations such as factorizations for sparse matrices is still very much an open research field. So, there is no .i() function available for cx_sp_mat or other sp_mat types. Another reason for this is lack of time on the part of the sparse matrix developers (...which includes me).

Given that the inverse of a sparse matrix is generally going to be dense, then you may simply be better off turning your cx_sp_mat into a cx_mat and then using the same inversion techniques that you normally would for dense matrices. Since you are planning to represent this as a dense matrix anyway, then it's a fair assumption that you have enough RAM to do that.

ryan
  • 839
  • 4
  • 13