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.

Advertisements

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[0][0]*matrix[1][1]-matrix[0][1]*matrix[1][0])
        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[0][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….

    Inverse of a matrix: http://www.geeksforgeeks.org/adjoint-inverse-matrix/
    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 [[0]*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
            
    
    

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: