Hi There all, I'm working on my first ever Graphics Engine and I'm stuck on the pretty hairy problem of Matrix inversion. I'm using the Building Blocks engine from the book Ultimate 3d Game Engine Design and Architecture as a base and then working off of that. The problem however, is that I have no Idea whatsoever what the code writer was doing when he wrote his matrix inversion code I really don't understand it. I've tried looking up the different methods but I still can't make any sense of it... Here's the code hopefully somebody will be able to tell me how it works or at least a basic idea of it. Thank you!

bool Matrix4x4::inverseMatrix(const Matrix4x4 &m)
{
   float tempMatrix[16] = {0};
   float d12, d13, d23, d24, d34, d41;

	d12 = matrix[2]  * matrix[7]  - matrix[3]  * matrix[6];
	d13 = matrix[2]  * matrix[11] - matrix[3]  * matrix[10];
	d23 = matrix[6]  * matrix[11] - matrix[7]  * matrix[10];
	d24 = matrix[6]  * matrix[15] - matrix[7]  * matrix[14];
	d34 = matrix[10] * matrix[15] - matrix[11] * matrix[14];
	d41 = matrix[14] * matrix[3]  - matrix[15] * matrix[2];

	tempMatrix[0] =   matrix[5] * d34 - matrix[9] * d24 + matrix[13] * d23;
	tempMatrix[1] = -(matrix[1] * d34 + matrix[9] * d41 + matrix[13] * d13);
	tempMatrix[2] =   matrix[1] * d24 + matrix[5] * d41 + matrix[13] * d12;
	tempMatrix[3] = -(matrix[1] * d23 - matrix[5] * d13 + matrix[9]  * d12);

	float determinant = matrix[0] * tempMatrix[0] + matrix[4] * tempMatrix[1] +
                       matrix[8] * tempMatrix[2] + matrix[12] * tempMatrix[3];

	if(determinant == 0.0)
	   {
		   Identity();
         return false;
	   }

	float invDeterminant = 1.0f / determinant;
	
	tempMatrix[0] *= invDeterminant;
	tempMatrix[1] *= invDeterminant;
	tempMatrix[2] *= invDeterminant;
	tempMatrix[3] *= invDeterminant;

	tempMatrix[4] = -(matrix[4] * d34 - matrix[8] * d24 + matrix[12] * d23) * invDeterminant;
	tempMatrix[5] =   matrix[0] * d34 + matrix[8] * d41 + matrix[12] * d13  * invDeterminant;
	tempMatrix[6] = -(matrix[0] * d24 + matrix[4] * d41 + matrix[12] * d12) * invDeterminant;
	tempMatrix[7] =   matrix[0] * d23 - matrix[4] * d13 + matrix[8]  * d12  * invDeterminant;

	d12 = matrix[0]  * matrix[5]  - matrix[1]  * matrix[12];
	d13 = matrix[0]  * matrix[9]  - matrix[1]  * matrix[8];
	d23 = matrix[4]  * matrix[9]  - matrix[5]  * matrix[8];
	d24 = matrix[4]  * matrix[13] - matrix[5]  * matrix[12];
	d34 = matrix[8]  * matrix[13] - matrix[9]  * matrix[12];
	d41 = matrix[12] * matrix[1]  - matrix[13] * matrix[0];

	tempMatrix[8]  =   matrix[7] * d34 - matrix[11] * d24 + matrix[15] * d23 * invDeterminant;
	tempMatrix[9]  = -(matrix[3] * d34 + matrix[11] * d41 + matrix[15] * d13) * invDeterminant;
	tempMatrix[10] =   matrix[3] * d24 + matrix[7]  * d41 + matrix[15] * d12 * invDeterminant;
	tempMatrix[11] = -(matrix[3] * d23 - matrix[7]  * d13 + matrix[11] * d12) * invDeterminant;
	tempMatrix[12] = -(matrix[6] * d34 - matrix[10] * d24 + matrix[14] * d23) * invDeterminant;
	tempMatrix[13] =   matrix[2] * d34 + matrix[10] * d41 + matrix[14] * d13 * invDeterminant;
	tempMatrix[14] = -(matrix[2] * d24 + matrix[6]  * d41 + matrix[14] * d12) * invDeterminant;
	tempMatrix[15] =   matrix[2] * d23 - matrix[6]  * d13 + matrix[10] * d12 * invDeterminant;

   matrix[0]  = tempMatrix[0];
   matrix[1]  = tempMatrix[1];
   matrix[2]  = tempMatrix[2];
   matrix[3]  = tempMatrix[3];

	matrix[4]  = tempMatrix[4];
	matrix[5]  = tempMatrix[5];
   matrix[6]  = tempMatrix[6];
   matrix[7]  = tempMatrix[7];

	matrix[8]  = tempMatrix[8];
	matrix[9]  = tempMatrix[9];
   matrix[10] = tempMatrix[10];
   matrix[11] = tempMatrix[11];

	matrix[12] = tempMatrix[12];
	matrix[13] = tempMatrix[13];
   matrix[14] = tempMatrix[14];
   matrix[15] = tempMatrix[15];

   return true;
}

Quote:

bool Matrix4x4::inverseMatrix(const Matrix4x4 &m)
{
float tempMatrix[16] = {0};
float d12, d13, d23, d24, d34, d41;

d12 = matrix[2] * matrix[7] - matrix[3] * matrix[6];
d13 = matrix[2] * matrix[11] - matrix[3] * matrix[10];
d23 = matrix[6] * matrix[11] - matrix[7] * matrix[10];
d24 = matrix[6] * matrix[15] - matrix[7] * matrix[14];
d34 = matrix[10] * matrix[15] - matrix[11] * matrix[14];
d41 = matrix[14] * matrix[3] - matrix[15] * matrix[2];

tempMatrix[0] = matrix[5] * d34 - matrix[9] * d24 + matrix[13] * d23;
tempMatrix[1] = -(matrix[1] * d34 + matrix[9] * d41 + matrix[13] * d13);
tempMatrix[2] = matrix[1] * d24 + matrix[5] * d41 + matrix[13] * d12;
tempMatrix[3] = -(matrix[1] * d23 - matrix[5] * d13 + matrix[9] * d12);

float determinant = matrix[0] * tempMatrix[0] + matrix[4] * tempMatrix[1] +
matrix[8] * tempMatrix[2] + matrix[12] * tempMatrix[3];

if(determinant == 0.0)
{
Identity();
return false;
}

The above part calculates the determinant of the matrix. If a matrix has a determinant equal to zero the matrix is said to be Singular. Singular matrices are not invertible.

Quote:
float invDeterminant = 1.0f / determinant;

tempMatrix[0] *= invDeterminant;
tempMatrix[1] *= invDeterminant;
tempMatrix[2] *= invDeterminant;
tempMatrix[3] *= invDeterminant;

tempMatrix[4] = -(matrix[4] * d34 - matrix[8] * d24 + matrix[12] * d23) * invDeterminant;
tempMatrix[5] = matrix[0] * d34 + matrix[8] * d41 + matrix[12] * d13 * invDeterminant;
tempMatrix[6] = -(matrix[0] * d24 + matrix[4] * d41 + matrix[12] * d12) * invDeterminant;
tempMatrix[7] = matrix[0] * d23 - matrix[4] * d13 + matrix[8] * d12 * invDeterminant;

d12 = matrix[0] * matrix[5] - matrix[1] * matrix[12];
d13 = matrix[0] * matrix[9] - matrix[1] * matrix[8];
d23 = matrix[4] * matrix[9] - matrix[5] * matrix[8];
d24 = matrix[4] * matrix[13] - matrix[5] * matrix[12];
d34 = matrix[8] * matrix[13] - matrix[9] * matrix[12];
d41 = matrix[12] * matrix[1] - matrix[13] * matrix[0];

tempMatrix[8] = matrix[7] * d34 - matrix[11] * d24 + matrix[15] * d23 * invDeterminant;
tempMatrix[9] = -(matrix[3] * d34 + matrix[11] * d41 + matrix[15] * d13) * invDeterminant;
tempMatrix[10] = matrix[3] * d24 + matrix[7] * d41 + matrix[15] * d12 * invDeterminant;
tempMatrix[11] = -(matrix[3] * d23 - matrix[7] * d13 + matrix[11] * d12) * invDeterminant;
tempMatrix[12] = -(matrix[6] * d34 - matrix[10] * d24 + matrix[14] * d23) * invDeterminant;
tempMatrix[13] = matrix[2] * d34 + matrix[10] * d41 + matrix[14] * d13 * invDeterminant;
tempMatrix[14] = -(matrix[2] * d24 + matrix[6] * d41 + matrix[14] * d12) * invDeterminant;
tempMatrix[15] = matrix[2] * d23 - matrix[6] * d13 + matrix[10] * d12 * invDeterminant;

matrix[0] = tempMatrix[0];
matrix[1] = tempMatrix[1];
matrix[2] = tempMatrix[2];
matrix[3] = tempMatrix[3];

matrix[4] = tempMatrix[4];
matrix[5] = tempMatrix[5];
matrix[6] = tempMatrix[6];
matrix[7] = tempMatrix[7];

matrix[8] = tempMatrix[8];
matrix[9] = tempMatrix[9];
matrix[10] = tempMatrix[10];
matrix[11] = tempMatrix[11];

matrix[12] = tempMatrix[12];
matrix[13] = tempMatrix[13];
matrix[14] = tempMatrix[14];
matrix[15] = tempMatrix[15];

return true;
}

The last part does the actual inversion defined by the following equation:

A^-1 = (1 / det(A)) * adj(A)

Notice how we multiply the inverse determinant times the ajoint of matrix A. This definition of inversion goes along way to explain why if determinate is zero the matrix doesn't have an inverse.

The only other interesting thing about the above definition is the adjoint. You can Google or refer to a math text on how to calculate the ajoint of a matrix.


As a side note... I'm not entirely sure why you would want to invert a 4x4 matrix (presuming it is only to contain a 3D affine transformation). 4x4 transformation matrices are frequently non-invertible. For example consider the matrix with zero translation, that would produce a row with all zeros. It is known that a matrix with any row or column with all zeros has a determinant equal to zero, and is therefore non-invertible:

|1 0 0 0 |
|0 1 0 0 |
|0 0 1 1 |
|0 0 0 0 |

We would rather perform inversion by inverting the 3x3 linear transformation matrix in the upper left, and negating the translation portion in the bottom row of the 4x4 matrix.

Thats just my two-cents; I'm curious what other people have to say on that matter.

Quote:Original post by fpsgamer
As a side note... I'm not entirely sure why you would want to invert a 4x4 matrix (presuming it is only to contain a 3D affine transformation). 4x4 transformation matrices are frequently non-invertible. For example consider the matrix with zero translation, that would produce a row with all zeros. It is known that a matrix with any row or column with all zeros has a determinant equal to zero, and is therefore non-invertible:

|1 0 0 0 |
|0 1 0 0 |
|0 0 1 1 |
|0 0 0 0 |

We would rather perform inversion by inverting the 3x3 linear transformation matrix in the upper left, and negating the translation portion in the bottom row of the 4x4 matrix.

Thats just my two-cents; I'm curious what other people have to say on that matter.

This isn't quite right.

First of all, a 4x4 homogeneous transform matrix with a translation of (0, 0, 0) does not have a row of all zeros, but rather a row of [0 0 0 1] (remember that a 3D point represented in 4D homogeneous coordinates has a w value of 1).

Second of all, there are some cases in basic 3D math where a general 4x4 inversion is useful. Although you can special-case just about any inversion, a general inversion function can be useful for, say, inverting projection matrices for picking purposes (such inversions can also be special-cased, but IMO it's not really worth the trouble).

Also, although inversion of various affine transforms can certainly be special-cased, inverting the linear portion and negating the translation won't always give the correct results (in fact, it usually won't).

Thanks for all your help! I was getting so confused whilst looking at the code and then looking at implementing a gaussian inversion method and then looking at the theory of cramer's rule. Hopefully I'll be able to implement my game engine in such a way that it actually runs. I'll keep the forums posted on my results.

-Jorick