-1
$\begingroup$

I'm going to build a C-library so I can do linear algebra at embedded systems. It's most for machine learning. https://github.com/DanielMartensson/EmbeddedAlgebra/

Anyway, I need to compute inverse of matrix $A$. What is the absolute best algorithm to use for that? Remember it must be a numerical algorithm.

Which one should I use?

Analytic solution

$\mathbf{A}^{-1}={1 \over \begin{vmatrix}\mathbf{A}\end{vmatrix}}\mathbf{C}^{\mathrm{T}}={1 \over \begin{vmatrix}\mathbf{A}\end{vmatrix}} \begin{pmatrix} \mathbf{C}_{11} & \mathbf{C}_{21} & \cdots & \mathbf{C}_{n1} \\ \mathbf{C}_{12} & \mathbf{C}_{22} & \cdots & \mathbf{C}_{n2} \\ \vdots & \vdots & \ddots & \vdots \\ \mathbf{C}_{1n} & \mathbf{C}_{2n} & \cdots & \mathbf{C}_{nn} \\ \end{pmatrix}$

Cholesky decomposition

$\mathbf{A}^{-1} = (\mathbf{L}^{*})^{-1} \mathbf{L}^{-1}$

Eigendecomposition

$\mathbf{A}^{-1}=\mathbf{Q}\mathbf{\Lambda}^{-1}\mathbf{Q}^{-1}$

Cayley–Hamilton method

$\mathbf{A}^{-1} = \frac{1}{\det (\mathbf{A})}\sum_{s=0}^{n-1}A^s$

Newton's method

$X_{k+1} = 2X_k - X_k A X_k$

Neuman series $\lim_{n \to \infty} (\mathbf I - \mathbf A)^n = 0$

Gaussian elimination

$A = \begin{bmatrix} 2 & -1 & 0 \\ -1 & 2 & -1 \\ 0 & -1 & 2 \end{bmatrix}$

To find the inverse of this matrix, one takes the following matrix augmented by the identity, and row reduces it as a 3 × 6 matrix: $ [ A | I ] = \left[ \begin{array}{rrr|rrr} 2 & -1 & 0 & 1 & 0 & 0\\ -1 & 2 & -1 & 0 & 1 & 0\\ 0 & -1 & 2 & 0 & 0 & 1 \end{array} \right]$

By performing row operations, one can check that the reduced row echelon form of this augmented matrix is:

$[ I | B ] = \left[ \begin{array}{rrr|rrr} 1 & 0 & 0 & \frac34 & \frac12 & \frac14\\[3pt] 0 & 1 & 0 & \frac12 & 1 & \frac12\\[3pt] 0 & 0 & 1 & \frac14 & \frac12 & \frac34 \end{array} \right]$

$\endgroup$
  • 5
    $\begingroup$ As a non-mathematical comment: for the love of god, don't. Firstly, you probably don't want to invert a matrix at all (see, for example, here). Secondly, use one of the many impementations that already exist. $\endgroup$ – user3482749 Jan 11 at 20:16
  • $\begingroup$ So what should I do? $\endgroup$ – Daniel Mårtensson Jan 11 at 21:24
  • $\begingroup$ In almost all cases, you should just solve $Ax = b$, since that's probably what you want to do anyway, and it's vastly quicker than computing $A^{-1}$ then calculating $A^{-1}b$ (indeed, it's quicker than either of those steps alone). If you really must invert $A$, just use one of the implementations that already exist in other people's algorithms. $\endgroup$ – user3482749 Jan 11 at 21:33
  • 1
    $\begingroup$ See here. The correct, answer, however, is almost certainly "just use somebody else's solution". $\endgroup$ – user3482749 Jan 11 at 22:35
  • 1
    $\begingroup$ All six of them are in C. $\endgroup$ – user3482749 Jan 12 at 11:27
3
$\begingroup$

It truly depends on the type of matrix you're going to compute the inverse from. Some methods are better for some classes of matrices than other.

But more importantly, why do you want to invert matrices? In many problems, you don't need to invert matrices, but only need to apply the inverse to some vectors. The latter problem is much easier to tackle, especially from a computational complexity standpoint (e.g. if your matrix is very large) and stability point of view (a matrix could have bad conditioning, be numerically close to non-invertible etc...). Does your application require that compute the matrix inverse?

$\endgroup$
  • $\begingroup$ No. But I need to find the inverse. Compute it or not. $\endgroup$ – Daniel Mårtensson Jan 11 at 21:26

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.