Sunday, September 30, 2012

Matrix Multiplication: 3x3 and 3x1


Formula:

| a11 a12 a13 |    | b1 |    | a11*b1 + a12*b2 + a13*b3 |
| a21 a22 a23 | x | b2 | = | a21*b1 + a22*b2 + a23*b3 |
| a31 a32 a33 |    | b3 |    | a31*b1 + a32*b2 + a33*b3 |

For arrays:

| a11 a12 a13 |      | a[0][0] a[0][1] a[0][2] |
| a21 a22 a23 |  =  | a[1][0] a[1][1] a[1][2] |
| a31 a32 a33 |      | a[2][0] a[2][1] a[2][2] |

| b1 |    | b[0] |
| b2 | = | b[1] |
| b3 |    | b[2] |

| c1 |    | a11 a12 a13 |    | b1 |    | a[0][0]*b[0] + a[0][1]*b[1] + a[0][2]*b[2] |
| c2 | = | a21 a22 a23 | x | b2 | = | a[1][0]*b[0] + a[1][1]*b[1] + a[1][2]*b[2] |
| c3 |    | a31 a32 a33 |    | b3 |    | a[2][0]*b[0] + a[2][1]*b[1] + a[2][2]*b[2] |

Code:   

c[0] = a[0][0]*b[0] + a[0][1]*b[1] + a[0][2]*b[2]
c[1] = a[1][0]*b[0] + a[1][1]*b[1] + a[1][2]*b[2]
c[2] = a[2][0]*b[0] + a[2][1]*b[1] + a[2][2]*b[2]

Matrix Inverse - 3x3

Operations of matrices of 3x3 are very common in day to day applications. Because the dimensions are already known, its better to avoid using a for-loop because of the extra computation involved in finding the subscript values.
Following is a quick formula for the inverse of a 3x3 matrix (implemented as a 2D array).

Formula

The inverse of a 3x3 matrix:
| a11 a12 a13 |-1                   |   a33a22-a32a23  -(a33a12-a32a13)   a23a12-a22a13   | 
| a21 a22 a23 |    =  1/DET *  | -(a33a21-a31a23)   a33a11-a31a13  -(a23a11-a21a13) |
| a31 a32 a33 |                      |   a32a21-a31a22  -(a32a11-a31a12)   a22a11-a21a12   |

where DET is the determinant of the matrix, i.e.
DET  =  a11(a33a22 - a32a23) - a21(a33a12 - a32a13) + a31(a23a12 - a22a13)

For 2D array

Moving on from matrices to 2D arrays
| a11 a12 a13 |        | mat[0][0] mat[0][1] mat[0][2] |
| a21 a22 a23 |   =   | mat[1][0] mat[1][1] mat[1][2] |
| a31 a32 a33 |        | mat[2][0] mat[2][1] mat[2][2] |
and,
| a11 a12 a13 |-1     | inv[0][0] inv[0][1] inv[0][2] |
| a21 a22 a23 |   =   | inv[1][0] inv[1][1] inv[1][2] |
| a31 a32 a33 |        | inv[2][0] inv[2][1] inv[2][2] |

Hence for code, the assignments are:
inv[0][0] = mat[2][2] * mat[1][1] - mat[2][1] * mat[1][2]
inv[0][1] = mat[2][1] * mat[0][2] - mat[2][2] * mat[0][1]
inv[0][2] = mat[1][2] * mat[0][1] - mat[1][1] * mat[0][2]
inv[1][0] = mat[2][0] * mat[1][2] - mat[2][2] * mat[1][0]
inv[1][1] = mat[2][2] * mat[0][0] - mat[2][0] * mat[0][2]
inv[1][2] = mat[1][0] * mat[0][2] - mat[1][2] * mat[0][0]
inv[2][0] = mat[2][1] * mat[1][0] - mat[2][0] * mat[1][1]
inv[2][1] = mat[2][0] * mat[0][1] - mat[2][1] * mat[0][0]
inv[2][2] = mat[1][1] * mat[0][0] - mat[1][0] * mat[0][1] 

DET =   mat[0][0]*inv[0][0] + mat[1][0]*inv[1][0] + mat[2][0]*inv[0][2]