NCL Home > Documentation > Functions > General applied math

# inverse_matrix

Computes the inverse of a square matrix using LU factorization.

## Prototype

```	function inverse_matrix (
A [*][*] : numeric
)

return_val [dimsizes(A)] :  float or double
```

## Arguments

A

A two-dimensional (N, N) array.

## Return value

An array of the same size, shape and type A.

## Description

NOTE: NCL does not have a numerical analysis focus. There are better tools for inverting matrices: Matlab, Octave, IDL and Python.

NOTE: Generally, it is best not to invert a matrix to solve an equation. See: Do not invert that matrix.

Perhaps, NCL's solve_linsys would be appropriate.

Using the LAPACK routines DGETRF and DGETRI, this function computes the inverse of a general matrix A(N, M) [here: N=M]using the LU factorization method. This method inverts U and then computes inv(A) by solving the system inv(A) * L = inv(U) for inv(A).

The returned inverse provides one way of solving the set of equations A * x = b (i.e. x = inv(A) * b). Thus, if the inverse is known, the vector x can be solved for any b. However, if efficiency is a concern, this is not the preferred method (see solve_linsys).

## Examples

The identity matrix (I) is found via I = aInv # a. The "#" operator computes the dot product of two arrays.

Example 1: Compute the determinant and inverse of a 3 x 3 matrix:

```   a    = (/(/1.0, -1.0, 2.0/) \
,(/3.0,  0.0, 1.0/) \
,(/1.0,  0.0, 2.0/) /)
write_matrix(a, 3f4.1, False)

aInv = inverse_matrix(a)           ; [3] x [3]
aId  = aInv#a

write_matrix(aInv, "3f4.1", False)
write_matrix(aId , "3f4.1", False)

----------------
qDet=5

1.0 -1.0  2.0      <=== original matrix
3.0  0.0  1.0
1.0  0.0  2.0
===

0.0  0.4 -0.2      <=== inverse
-1.0 -0.0  1.0
0.0 -0.2  0.6
===

1.0  0.0  0.0      <=== identity
0.0  1.0  0.0
0.0  0.0  1.0

```
Example 2: Compute the determinant and inverse of a 4 x 4 matrix:

```   x = (/  (/3.0, 2.0, -1.0,  2.0/), \
(/1.0, 4.0,  0.0,  2.0/), \
(/2.0, 1.0,  2.0, -1.0/), \
(/1.0, 1.0, -1.0,  3.0/) /)
write_matrix(x, 4f4.1, False)

xDet = determinant(x)
xInv = inverse_matrix(x)           ; [4] x [4]
xId  = xInv#x

print("xDet="+xDet)
write_matrix(xInv, "4f4.1", False)
write_matrix(toint(xId) , "4i4" , False)

----------------
xDet=44

3.0  2.0 -1.0  2.0
1.0  4.0  0.0  2.0
2.0  1.0  2.0 -1.0
1.0  1.0 -1.0  3.0
===

0.3 -0.2  0.1 -0.0
0.1  0.3 -0.1 -0.3
-0.5  0.0  0.5  0.5
-0.3 -0.0  0.2  0.6
===

1.0  0.0  0.0  0.0
0.0  1.0  0.0  0.0
0.0  0.0  1.0  0.0
0.0  0.0  0.0  1.0

```
Example 3: Compute the inverse of a NxN matrix (N=20). Due to computer numerics, values such as -0.0 and 0.99999 are returned. These are manipulated using round and toint to return the expected value [0.0 and 1.0].

```  b     = random_normal( 0,10, (/N,N/))
write_matrix(b, "20f6.1", False)

bDet  = determinant(b)
bInv  = inverse_matrix(b)
bId   = bInv#b

print("bDet="+bDet)
bId   = where(abs(bId).lt.1e-5, 0.0, bId)    ; cleaner output
write_matrix(toint(round(bID,0),"20i3",False)             ; printing convenience

-----------------
bDet=-6.023e+28

8.7   8.7  19.5   9.8   3.5  -1.3  -8.6   8.8   9.3 -11.7   9.2   7.4 -28.2   0.1   0.5   7.5   1.8  -6.2  10.2  -5.0
11.7  -6.7  -3.3  16.9  15.1   0.3   3.5   1.7 -13.8   8.1  20.8   9.1 -11.3  -2.9  17.3   3.2  14.2  18.0  -5.2   1.6
-6.2  -0.5  -4.1   3.5  13.9  -0.8 -12.3  11.1   3.1  -5.4 -28.3 -11.8  15.9  -1.9 -11.1  20.9 -16.3 -11.3   5.6  -6.6
4.2  -6.1   6.0  -5.0  11.1  13.9  -1.2  10.6  18.6   2.2   3.1   1.3   3.2  15.5  -9.8 -13.7  -1.6   1.8  10.8 -11.3
9.9  -7.2  -1.1  10.3   8.0   7.1  -2.6  -7.7   5.7   8.7   4.3  16.8   9.5  -6.0  -3.5   1.7  -5.0  11.6  -4.6  -0.0
-9.0  -3.2  -2.0 -16.4   4.8   2.2  -3.0  -6.9   8.5  15.7  -0.9   7.6  29.0   5.1  -1.6 -10.0   6.8 -11.9  -3.7  -4.8
-4.1   6.2 -14.3  11.4 -10.9   0.9  -3.4   7.2  10.2   3.7  -0.1   4.3  12.5 -11.2  11.9   6.7   2.5   4.8  12.5   9.0
-5.6  -7.2   5.4   9.7 -22.7   2.7   5.6   7.6   4.8  -9.7   1.7  -6.2  24.3  -7.8  -4.9  -3.3   1.0 -21.3 -18.4  17.8
-10.5  16.5  19.0  -6.5 -14.1  -2.3  13.7  -0.4   8.5  -3.2   3.0 -23.0   7.0  -5.5  20.9  19.1   0.1 -12.8  -6.7  -5.5
9.8   8.8  17.0  -1.0  -4.1   3.7  17.7   3.6  22.6 -25.5 -13.8   8.1  -2.8  11.1  -3.0   0.1   7.9   5.4   9.3  15.2
-3.9   4.1 -10.6  10.8   3.4   7.3 -13.2   6.7 -11.0  11.6   3.0  14.7   0.2  -5.8 -22.6  -8.8 -18.7  10.9  -6.1   0.5
-7.4  11.4  -3.8 -21.2  -6.0  -9.6  -4.3   5.4   9.4   5.4  -7.1   4.2   6.8  -7.8  12.3  -0.8  -3.5  -1.6 -11.6  -3.8
2.0  16.6  -2.9  -6.6  10.9   5.7   2.8  12.2   4.9   4.6  11.2 -10.3  -9.7 -19.4   1.3 -14.0 -14.1   8.9   8.1   7.1
4.7   1.4  13.6  -8.6 -12.7  -1.1  -1.9  -1.8   7.4  -7.2  -4.0  -7.4   3.8  -6.3  -4.9   6.2  -5.2 -14.1  10.0  -0.5
9.8 -13.3  -3.5   9.7  -6.9 -10.0 -11.4   6.2  -3.9   1.0 -10.4  -7.5  17.6   1.9  -8.9  -0.3  -1.1  13.4  -8.4   6.0
12.5  -1.7  -0.3   2.0  -0.8  -2.2   4.4  -5.1  -2.6 -14.9   0.2  13.9   5.7   2.4 -13.7  -0.2  -9.4   0.0   7.1   6.5
4.4  11.5 -10.8  10.7   5.6  12.0  -9.2   2.7   9.2   8.8   0.7   9.5   7.6   3.2   5.8  -6.1  -7.9   7.1   2.3  -5.0
-1.9  -7.9   4.8   7.9 -19.1   7.7   3.0 -11.7  -7.8  -7.5  -6.1  -7.0 -15.2  24.6  -0.8  -8.1 -17.4 -11.2  -4.6   5.0
4.3  -5.7   7.2   3.0   2.9  -3.0   3.6 -15.1 -14.3   6.6  10.8  -3.5 -12.7  -3.2  -3.7 -14.5   0.3  -8.8  16.2 -19.0
-3.1  -2.4  -1.7  -8.2 -11.0  24.6  -9.3  -3.8   9.3  -5.2  -6.6  -9.9  19.2  -0.2  16.6  15.4  -7.3 -16.3  -3.9  -7.9

-0.0450  0.0339  0.0001  0.0072 -0.0108 -0.0399 -0.0412  0.0143 -0.0226 -0.0149 -0.0232  0.0040 -0.0220  0.1120 -0.0261 -0.0096  0.0795 -0.0107 -0.0256 -0.0314
-0.0055 -0.0049 -0.0075 -0.0195 -0.0247  0.0011 -0.0207 -0.0025  0.0088  0.0031  0.0059 -0.0152 -0.0048  0.0213 -0.0061  0.0000  0.0503 -0.0139 -0.0107 -0.0106
0.0304  0.0221  0.0084 -0.0108  0.0034  0.0467  0.0114 -0.0118  0.0194  0.0372  0.0350  0.0072  0.0187 -0.0222  0.0322 -0.0067 -0.0483  0.0176  0.0275  0.0206
0.0106 -0.0091  0.0070 -0.0059  0.0111 -0.0037  0.0075  0.0114  0.0063  0.0027 -0.0041 -0.0093 -0.0024 -0.0255  0.0058 -0.0058  0.0120  0.0020  0.0143 -0.0059
0.0190 -0.0032  0.0147 -0.0109  0.0089  0.0243 -0.0068 -0.0031  0.0003  0.0018 -0.0129 -0.0097  0.0220 -0.0397  0.0055  0.0088 -0.0116  0.0013  0.0019  0.0049
-0.0165  0.0147 -0.0030  0.0002 -0.0076 -0.0061 -0.0100  0.0015 -0.0095  0.0175  0.0214 -0.0112 -0.0027  0.0241 -0.0141 -0.0219  0.0064 -0.0052 -0.0008  0.0205
-0.0522  0.0244  0.0104  0.0258  0.0046 -0.0398 -0.0092  0.0164  0.0027 -0.0002 -0.0038  0.0132 -0.0223  0.0504 -0.0434 -0.0056  0.0241 -0.0025 -0.0047 -0.0320
-0.0199  0.0482  0.0157  0.0310 -0.0304 -0.0156  0.0033  0.0137 -0.0049  0.0057  0.0215  0.0273 -0.0129  0.0457 -0.0139 -0.0065  0.0043  0.0013 -0.0039 -0.0112
0.0083 -0.0431 -0.0094  0.0127  0.0373 -0.0193  0.0013  0.0064 -0.0020 -0.0225 -0.0369 -0.0064 -0.0001 -0.0144 -0.0077 -0.0051  0.0122 -0.0057 -0.0179 -0.0152
-0.0434  0.0397  0.0089  0.0107  0.0027 -0.0100 -0.0108  0.0024 -0.0094 -0.0096  0.0045  0.0014 -0.0209  0.1093 -0.0348 -0.0420  0.0420  0.0104 -0.0262 -0.0436
0.0223 -0.0198 -0.0207  0.0109  0.0019  0.0086  0.0077 -0.0036  0.0183 -0.0308 -0.0094 -0.0194  0.0100 -0.0328  0.0073  0.0353 -0.0137 -0.0050 -0.0153  0.0037
0.0035  0.0232  0.0067  0.0023  0.0015  0.0096  0.0127  0.0008 -0.0050  0.0158  0.0218  0.0260 -0.0109  0.0031 -0.0113  0.0074 -0.0150  0.0069  0.0103  0.0030
0.0079  0.0046  0.0024 -0.0027 -0.0083  0.0274  0.0100 -0.0030  0.0165  0.0096  0.0106 -0.0024  0.0082 -0.0268  0.0234  0.0222 -0.0120 -0.0018  0.0132  0.0143
-0.0013  0.0079 -0.0025  0.0108 -0.0189  0.0117 -0.0032 -0.0102  0.0090 -0.0097 -0.0034 -0.0106 -0.0095  0.0090  0.0008  0.0128  0.0156  0.0135 -0.0217 -0.0120
0.0181  0.0237  0.0108 -0.0073 -0.0014  0.0245  0.0211 -0.0053  0.0036  0.0172 -0.0001  0.0251  0.0182 -0.0342  0.0178  0.0065 -0.0253  0.0214  0.0232  0.0213
-0.0165  0.0085  0.0019  0.0077  0.0047 -0.0180 -0.0042 -0.0075  0.0061 -0.0196 -0.0018 -0.0127 -0.0209  0.0459 -0.0237  0.0009  0.0101 -0.0060 -0.0339 -0.0148
-0.0081 -0.0041 -0.0122 -0.0093 -0.0178 -0.0066 -0.0139  0.0055 -0.0135  0.0063 -0.0036 -0.0185 -0.0178  0.0189 -0.0060 -0.0247  0.0229 -0.0255 -0.0051 -0.0018
0.0152 -0.0172 -0.0128 -0.0004  0.0101  0.0094  0.0177 -0.0231  0.0188  0.0167  0.0190 -0.0018  0.0131 -0.0566  0.0360  0.0077 -0.0411 -0.0016  0.0133  0.0237
0.0055  0.0116  0.0063  0.0030 -0.0120  0.0212  0.0313 -0.0158  0.0022  0.0109  0.0109 -0.0050  0.0087  0.0026  0.0068  0.0082 -0.0237  0.0093  0.0138  0.0061
0.0083  0.0126  0.0024 -0.0211 -0.0006  0.0338  0.0072 -0.0056 -0.0090  0.0048 -0.0024 -0.0182  0.0231  0.0161  0.0028 -0.0047 -0.0127  0.0197 -0.0243 -0.0050

1  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  1  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  1  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  1  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  0  1  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  1  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  1  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  1  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  1  0  0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  1  0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  1  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  1  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  1  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  1  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  1  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  1  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  1  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  1  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  1  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  1

```

## Errors

The LAPACK subrourines return the following information:

info
= 0: successful exit

< 0: if info = -i, the i-th argument had an illegal value

> 0: if info = i, U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, and division by zero will occur if it is used to solve a system of equations.