# inverting a 4x4 matrix

### Question

i am looking for a sample code implementation on how to invert a 4x4 matrix. i know there is gaussian eleminiation, LU decomposition, etc. but instead of looking at them in detail i am really just looking for the code to do this.

language ideally C++, data is available in array of 16 floats in cloumn-major order.

thank you!

1
74
6/27/2013 12:24:59 PM

here:

``````bool gluInvertMatrix(const double m[16], double invOut[16])
{
double inv[16], det;
int i;

inv[0] = m[5]  * m[10] * m[15] -
m[5]  * m[11] * m[14] -
m[9]  * m[6]  * m[15] +
m[9]  * m[7]  * m[14] +
m[13] * m[6]  * m[11] -
m[13] * m[7]  * m[10];

inv[4] = -m[4]  * m[10] * m[15] +
m[4]  * m[11] * m[14] +
m[8]  * m[6]  * m[15] -
m[8]  * m[7]  * m[14] -
m[12] * m[6]  * m[11] +
m[12] * m[7]  * m[10];

inv[8] = m[4]  * m[9] * m[15] -
m[4]  * m[11] * m[13] -
m[8]  * m[5] * m[15] +
m[8]  * m[7] * m[13] +
m[12] * m[5] * m[11] -
m[12] * m[7] * m[9];

inv[12] = -m[4]  * m[9] * m[14] +
m[4]  * m[10] * m[13] +
m[8]  * m[5] * m[14] -
m[8]  * m[6] * m[13] -
m[12] * m[5] * m[10] +
m[12] * m[6] * m[9];

inv[1] = -m[1]  * m[10] * m[15] +
m[1]  * m[11] * m[14] +
m[9]  * m[2] * m[15] -
m[9]  * m[3] * m[14] -
m[13] * m[2] * m[11] +
m[13] * m[3] * m[10];

inv[5] = m[0]  * m[10] * m[15] -
m[0]  * m[11] * m[14] -
m[8]  * m[2] * m[15] +
m[8]  * m[3] * m[14] +
m[12] * m[2] * m[11] -
m[12] * m[3] * m[10];

inv[9] = -m[0]  * m[9] * m[15] +
m[0]  * m[11] * m[13] +
m[8]  * m[1] * m[15] -
m[8]  * m[3] * m[13] -
m[12] * m[1] * m[11] +
m[12] * m[3] * m[9];

inv[13] = m[0]  * m[9] * m[14] -
m[0]  * m[10] * m[13] -
m[8]  * m[1] * m[14] +
m[8]  * m[2] * m[13] +
m[12] * m[1] * m[10] -
m[12] * m[2] * m[9];

inv[2] = m[1]  * m[6] * m[15] -
m[1]  * m[7] * m[14] -
m[5]  * m[2] * m[15] +
m[5]  * m[3] * m[14] +
m[13] * m[2] * m[7] -
m[13] * m[3] * m[6];

inv[6] = -m[0]  * m[6] * m[15] +
m[0]  * m[7] * m[14] +
m[4]  * m[2] * m[15] -
m[4]  * m[3] * m[14] -
m[12] * m[2] * m[7] +
m[12] * m[3] * m[6];

inv[10] = m[0]  * m[5] * m[15] -
m[0]  * m[7] * m[13] -
m[4]  * m[1] * m[15] +
m[4]  * m[3] * m[13] +
m[12] * m[1] * m[7] -
m[12] * m[3] * m[5];

inv[14] = -m[0]  * m[5] * m[14] +
m[0]  * m[6] * m[13] +
m[4]  * m[1] * m[14] -
m[4]  * m[2] * m[13] -
m[12] * m[1] * m[6] +
m[12] * m[2] * m[5];

inv[3] = -m[1] * m[6] * m[11] +
m[1] * m[7] * m[10] +
m[5] * m[2] * m[11] -
m[5] * m[3] * m[10] -
m[9] * m[2] * m[7] +
m[9] * m[3] * m[6];

inv[7] = m[0] * m[6] * m[11] -
m[0] * m[7] * m[10] -
m[4] * m[2] * m[11] +
m[4] * m[3] * m[10] +
m[8] * m[2] * m[7] -
m[8] * m[3] * m[6];

inv[11] = -m[0] * m[5] * m[11] +
m[0] * m[7] * m[9] +
m[4] * m[1] * m[11] -
m[4] * m[3] * m[9] -
m[8] * m[1] * m[7] +
m[8] * m[3] * m[5];

inv[15] = m[0] * m[5] * m[10] -
m[0] * m[6] * m[9] -
m[4] * m[1] * m[10] +
m[4] * m[2] * m[9] +
m[8] * m[1] * m[6] -
m[8] * m[2] * m[5];

det = m[0] * inv[0] + m[1] * inv[4] + m[2] * inv[8] + m[3] * inv[12];

if (det == 0)
return false;

det = 1.0 / det;

for (i = 0; i < 16; i++)
invOut[i] = inv[i] * det;

return true;
}
``````

This was lifted from MESA implementation of the GLU library.

91
2/15/2012 9:35:13 PM

If anyone looking for more costumized code and "easier to read", then I got this

``````var A2323 = m.m22 * m.m33 - m.m23 * m.m32 ;
var A1323 = m.m21 * m.m33 - m.m23 * m.m31 ;
var A1223 = m.m21 * m.m32 - m.m22 * m.m31 ;
var A0323 = m.m20 * m.m33 - m.m23 * m.m30 ;
var A0223 = m.m20 * m.m32 - m.m22 * m.m30 ;
var A0123 = m.m20 * m.m31 - m.m21 * m.m30 ;
var A2313 = m.m12 * m.m33 - m.m13 * m.m32 ;
var A1313 = m.m11 * m.m33 - m.m13 * m.m31 ;
var A1213 = m.m11 * m.m32 - m.m12 * m.m31 ;
var A2312 = m.m12 * m.m23 - m.m13 * m.m22 ;
var A1312 = m.m11 * m.m23 - m.m13 * m.m21 ;
var A1212 = m.m11 * m.m22 - m.m12 * m.m21 ;
var A0313 = m.m10 * m.m33 - m.m13 * m.m30 ;
var A0213 = m.m10 * m.m32 - m.m12 * m.m30 ;
var A0312 = m.m10 * m.m23 - m.m13 * m.m20 ;
var A0212 = m.m10 * m.m22 - m.m12 * m.m20 ;
var A0113 = m.m10 * m.m31 - m.m11 * m.m30 ;
var A0112 = m.m10 * m.m21 - m.m11 * m.m20 ;

var det = m.m00 * ( m.m11 * A2323 - m.m12 * A1323 + m.m13 * A1223 )
- m.m01 * ( m.m10 * A2323 - m.m12 * A0323 + m.m13 * A0223 )
+ m.m02 * ( m.m10 * A1323 - m.m11 * A0323 + m.m13 * A0123 )
- m.m03 * ( m.m10 * A1223 - m.m11 * A0223 + m.m12 * A0123 ) ;
det = 1 / det;

return new Matrix4x4() {
m00 = det *   ( m.m11 * A2323 - m.m12 * A1323 + m.m13 * A1223 ),
m01 = det * - ( m.m01 * A2323 - m.m02 * A1323 + m.m03 * A1223 ),
m02 = det *   ( m.m01 * A2313 - m.m02 * A1313 + m.m03 * A1213 ),
m03 = det * - ( m.m01 * A2312 - m.m02 * A1312 + m.m03 * A1212 ),
m10 = det * - ( m.m10 * A2323 - m.m12 * A0323 + m.m13 * A0223 ),
m11 = det *   ( m.m00 * A2323 - m.m02 * A0323 + m.m03 * A0223 ),
m12 = det * - ( m.m00 * A2313 - m.m02 * A0313 + m.m03 * A0213 ),
m13 = det *   ( m.m00 * A2312 - m.m02 * A0312 + m.m03 * A0212 ),
m20 = det *   ( m.m10 * A1323 - m.m11 * A0323 + m.m13 * A0123 ),
m21 = det * - ( m.m00 * A1323 - m.m01 * A0323 + m.m03 * A0123 ),
m22 = det *   ( m.m00 * A1313 - m.m01 * A0313 + m.m03 * A0113 ),
m23 = det * - ( m.m00 * A1312 - m.m01 * A0312 + m.m03 * A0112 ),
m30 = det * - ( m.m10 * A1223 - m.m11 * A0223 + m.m12 * A0123 ),
m31 = det *   ( m.m00 * A1223 - m.m01 * A0223 + m.m02 * A0123 ),
m32 = det * - ( m.m00 * A1213 - m.m01 * A0213 + m.m02 * A0113 ),
m33 = det *   ( m.m00 * A1212 - m.m01 * A0212 + m.m02 * A0112 ),
};
``````

I don't write the code, but my program did. I made a small program to make a program that calculate the determinant and inverse of any N-matrix.

I do it because once in the past I need a code that inverses 5x5 matrix, but nobody in the earth have done this so I made one.

Take a look about the program here.

EDIT: The matrix layout is row-by-row (meaning `m01` is in the first row and second column). Also the language is C#, but should be easy to convert into C.