Most efficient matrix inversion in MATLAB

Question

When computing the inverse for some square matrix A in MATLAB, using

Ai = inv(A)
% should be the same as:
Ai = A^-1

MATLAB usually notifies me that this is not the most efficient way of inverting. So what's more efficient? If I have an equation system, using the /,\ operators probably is. But sometimes I need the inverse for other computations.

What's the most efficient way to invert?

Answer 1

If there isn't a clever way to do all your calculations without explicitly forming the inverse then you have to use the "inv" function. You could of course solve a linear system with your matrix and the identity matrix, but there is nothing to be gained by doing that.

Answer 2

I think LU decomposition is more efficient than than inversion (and potentially more stable if you use pivoting). It works especially well if you need to solve for more than one right hand side vector, because once you have the LU decomposition you can do the forward back substitutions for each one as you need it.

I would recommend LU decomposition over a full inverse. I agree if that's what MATLAB is saying.

UPDATE: 3x3 matrix? You can invert that by hand in closed form if you need it. Just check the determinant first to make sure that it's not singular or nearly singular.

Answer 3

I would recommend to use svd (unless you are really absolute sure that your matrix is not ill-conditioned). Then, based on singular values you make your decisions on further actions to take. This may sound like a 'overkill' approach, but in long run it will pay back.

Now if your matrix A is actually invertible, then the pseudo inverse of A coincidence with inv(A), however if you are close to 'singularity' you'll easily make appropriate decision how to proceed to actually make the pseudo inverse. Naturally these decisions will depend on your application.

Added a straightforward example:

> A= randn(3, 2); A= [A A(:, 1)+ A(:, 2)]
A =
  -1.520342  -0.239380  -1.759722
   0.022604   0.381374   0.403978
   0.852420   1.521925   2.374346

> inv(A)
warning: inverse: matrix singular to machine precision, rcond = 0
ans =
   Inf   Inf   Inf
   Inf   Inf   Inf
   Inf   Inf   Inf

> [U, S, V]= svd(A)
U =
  -0.59828  -0.79038   0.13178
   0.13271  -0.25993  -0.95646
   0.79022  -0.55474   0.26040

S =
Diagonal Matrix
  3.6555e+000            0            0
            0  1.0452e+000            0
            0            0  1.4645e-016

V =
   0.433921   0.691650   0.577350
   0.382026  -0.721611   0.577350
   0.815947  -0.029962  -0.577350

> s= diag(S); k= sum(s> 1e-9) % simple thresholding based decision
k =  2

> Ainv= (U(:, 1: k)* diag(1./ s(1: k))* V(:, 1: k)')'
Ainv =
  -0.594055  -0.156258  -0.273302
   0.483170   0.193333   0.465592
  -0.110885   0.037074   0.192290

> A* Ainv
ans =
   0.982633   0.126045  -0.034317
   0.126045   0.085177   0.249068
  -0.034317   0.249068   0.932189

> A* pinv(A)
ans =
   0.982633   0.126045  -0.034317
   0.126045   0.085177   0.249068
  -0.034317   0.249068   0.932189

Answer 4

If you only need the inverse then just do, it will be numerically more stable than inv(A):

inv_A = 1\A;