Matrix Determinant And Inverse

May 23, 2017

Today’s exercise is preliminary to the exercise we will have later this week. You are to write programs that calculate the determinant and inverse of a matrix. I won’t go into the math involved behind the matrix arithmetic, as there are many sources on the internet that know far more about the process than I. Google for “matrix determinant” or “matrix inverse”; I used YouTube for my instruction.

Your task is to write programs that calculate the determinant and inverse of a matrix. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

5 Responses to “Matrix Determinant And Inverse”

1. John Cowan said

No worse than addend, summand, subtrahend, minuend, multiplicand.

2. Rutger said

In Python, the determinant part.

```def determinant(matrix):
if len(matrix) == 2:
return (matrix*matrix-matrix*matrix)
else:
negator = 1
result = 0
for col in range(len(matrix)):
sub_matrix = [[] for i in range(len(matrix)-1)]
for sub_row in range(1,len(matrix)):
for sub_col in range(len(matrix)):
if sub_col != col:
sub_matrix[sub_row-1].append(matrix[sub_row][sub_col])
result += negator * matrix[col] * determinant(sub_matrix)
negator *= -1
return result

print(-240 == determinant([[1,-4,9],[-6,7,3],[1,2,3]]))
print(441 == determinant([[6,24,1],[13,16,10],[20,17,15]]))
```
3. Milbrae said

I’ve found these two links with explanations and code. It’s quite a hassle for larger dimensions….

Determinant of a matrix: http://www.geeksforgeeks.org/determinant-of-a-matrix/

Though, I haven’t tried their code yet.

4. Paul said

Here in Python an implementation of an LU decomposition. Once the matrix is decomposed, the determinant and inverse are straightforward.

```from copy import deepcopy
from functools import reduce
from operator import mul

def trans(A): return [list(a) for a in zip(*A)]
def dot(a, b): return sum(ai*bi for ai, bi in zip(a, b))
def matvec(A, b): return [dot(a, b) for a in A]
def matmul(A, B): return trans(matvec(A, b) for b in zip(*B))
def zeros(m, n): return [*n for i in range(m)]

def lu_decompose(A, tol):
""" returns matrix LU with upper and lower matrix
vector P: P[:n] contains the ordering of rows
P[n] = n + numbeer of swaps in P[:n]
"""
LU = deepcopy(A)  # the original matrix is not changed
n = len(LU)
P = list(range(n+1))
for i in range(n):
imax = max(range(i, n), key=lambda m: abs(LU[m][i]))
if abs(LU[i][imax]) < tol:
raise ValueError("Matrix is degenerate")
if imax != i:
P[i], P[imax] = P[imax], P[i]
LU[i], LU[imax] = LU[imax], LU[i]
P[n] += 1
for j in range(i+1, n):
LU[j][i] /= LU[i][i]
for k in range(i+1, n):
LU[j][k] -= LU[j][i] * LU[i][k];
return LU, P

def lu_solve(LU, P, b):
"""solve Ax=b A is te original matrix"""
n = len(LU)
x = [b[p] for p in P[:-1]]
for i in range(n):
x[i] = b[P[i]]
for k in range(i):
x[i] -= LU[i][k] * x[k]
for i in range(n-1, -1, -1):
for  k in range(i + 1, n):
x[i] -= LU[i][k] * x[k]
x[i] /= LU[i][i]
return x

def lu_invert(LU, P):
"""returns the inverse of the original matrix"""
n = len(LU)
IA = zeros(n, n)
for j in range(n):
for i in range(n):
IA[i][j] = 1.0 if P[i] == j else  0.0
for k in range(i):
IA[i][j] -= LU[i][k] * IA[k][j]
for i in range(n-1, -1, -1):
for k in range(i+1, n):
IA[i][j] -= LU[i][k] * IA[k][j]
IA[i][j] /= LU[i][i]
return IA

def lu_determinant(LU, P):
"""returns the determinant of the original matrix"""
n = len(LU)
det = reduce(mul, (LU[i][i] for i in range(n)))
return det if (P[n]-n) % 2 == 0 else -det

def print_matrix(A):
print()
for row in A:
for ai in row:
print("{:8.4f}".format(ai), end=" ")
print()

b = [[6, 24, 1], [13, 16, 10], [20, 17, 15]]
A, P = lu_decompose(b, 1e-8)
ib = lu_invert(A, P)
print_matrix(ib)
print_matrix(matmul(b, ib))
print(lu_determinant(A, P))
```
```inverse
0.1587  -0.7778   0.5079
0.0113   0.1587  -0.1066
-0.2245   0.8571  -0.4898
b*invb
1.0000  -0.0000   0.0000
0.0000   1.0000   0.0000
0.0000   0.0000   1.0000
determinant 440.99999999999994
```
5. Steve said

Klong (for determinant)

```steve@steve-Satellite-L555D:~\$ rlwrap kg
Klong 20161212
:" Beginning of code"
det::{:[2=#x; cal1(x); cal2(x)]}
cal1::{((x:@[0 0])*(x:@[1 1]))-((x:@[0 1])*(x:@[1 0]))}
mults::{[a b]; a::1; b::[]; {b::b,(x*a); a::-a}'*x; b}
makerow::{[a pos rn row t]; a::[]; t::x; rn::y; pos::z; row::t@rn; ((1+&#row):=0,pos){:[x=1; a::a,y; a]}'row; a}
maketbl::{[pos tbl]; tbl::x; pos::y; {makerow(tbl; x; pos)}'(1+!(#tbl)-1)}
maketbls::{[t]; t::x; {maketbl(t; x)}'!#t}
cal2::{[a t]; a::0; t::x; mults(t){a::a+(x*det(y))}'maketbls(t); a}
:" End of code"

t2::[[1 2] [3 4]]
t3::[[1 2 3] [2 3 4] [3 4 5]]
t4::[[1 2 3 4] [2 3 4 5] [3 4 55 666] [7 88 999 1234]]

:" Should be -2"
det(t2)
-2
:" Should be 0"
det(t3)
0
:" Should be 498600"
det(t4)
498600
:" Should be 18 per Wikipedia"
det([[-2 2 -3] [-1 1 3] [2 0 -1]])
18

:" Per @Rutger should be -240"
det([[1 -4 9] [-6 7 3] [1 2 3]])
-240
:" Per @Rutger shuld be 441"
det([[6 24 1] [13 16 10] [20 17 15]])
441

```