I am trying to get the inverse of a 4x4 square matrix for opengl, so it is a column major matrix. Also I would like to avoid discussions about making my code into loops, it's quite challenging enough to follow without loops.
I wrote some code that can do most basic matrix operations but matrix inverse was kinda rough, I found this answer: inverting a 4x4 matrix
I implemented it in my code
void mat4inverse(mat4* out, const mat4* in)
{
double inv[16], det;
int i;
inv[0] = in->m[5] * in->m[10] * in->m[15] -
in->m[5] * in->m[11] * in->m[14] -
in->m[9] * in->m[6] * in->m[15] +
in->m[9] * in->m[7] * in->m[14] +
in->m[13] * in->m[6] * in->m[11] -
in->m[13] * in->m[7] * in->m[10];
inv[4] = -in->m[4] * in->m[10] * in->m[15] +
in->m[4] * in->m[11] * in->m[14] +
in->m[8] * in->m[6] * in->m[15] -
in->m[8] * in->m[7] * in->m[14] -
in->m[12] * in->m[6] * in->m[11] +
in->m[12] * in->m[7] * in->m[10];
inv[8] = in->m[4] * in->m[9] * in->m[15] -
in->m[4] * in->m[11] * in->m[13] -
in->m[8] * in->m[5] * in->m[15] +
in->m[8] * in->m[7] * in->m[13] +
in->m[12] * in->m[5] * in->m[11] -
in->m[12] * in->m[7] * in->m[9];
inv[12] = -in->m[4] * in->m[9] * in->m[14] +
in->m[4] * in->m[10] * in->m[13] +
in->m[8] * in->m[5] * in->m[14] -
in->m[8] * in->m[6] * in->m[13] -
in->m[12] * in->m[5] * in->m[10] +
in->m[12] * in->m[6] * in->m[9];
inv[1] = -in->m[1] * in->m[10] * in->m[15] +
in->m[1] * in->m[11] * in->m[14] +
in->m[9] * in->m[2] * in->m[15] -
in->m[9] * in->m[3] * in->m[14] -
in->m[13] * in->m[2] * in->m[11] +
in->m[13] * in->m[3] * in->m[10];
inv[5] = in->m[0] * in->m[10] * in->m[15] -
in->m[0] * in->m[11] * in->m[14] -
in->m[8] * in->m[2] * in->m[15] +
in->m[8] * in->m[3] * in->m[14] +
in->m[12] * in->m[2] * in->m[11] -
in->m[12] * in->m[3] * in->m[10];
inv[9] = -in->m[0] * in->m[9] * in->m[15] +
in->m[0] * in->m[11] * in->m[13] +
in->m[8] * in->m[1] * in->m[15] -
in->m[8] * in->m[3] * in->m[13] -
in->m[12] * in->m[1] * in->m[11] +
in->m[12] * in->m[3] * in->m[9];
inv[13] = in->m[0] * in->m[9] * in->m[14] -
in->m[0] * in->m[10] * in->m[13] -
in->m[8] * in->m[1] * in->m[14] +
in->m[8] * in->m[2] * in->m[13] +
in->m[12] * in->m[1] * in->m[10] -
in->m[12] * in->m[2] * in->m[9];
inv[2] = in->m[1] * in->m[6] * in->m[15] -
in->m[1] * in->m[7] * in->m[14] -
in->m[5] * in->m[2] * in->m[15] +
in->m[5] * in->m[3] * in->m[14] +
in->m[13] * in->m[2] * in->m[7] -
in->m[13] * in->m[3] * in->m[6];
inv[6] = -in->m[0] * in->m[6] * in->m[15] +
in->m[0] * in->m[7] * in->m[14] +
in->m[4] * in->m[2] * in->m[15] -
in->m[4] * in->m[3] * in->m[14] -
in->m[12] * in->m[2] * in->m[7] +
in->m[12] * in->m[3] * in->m[6];
inv[10] = in->m[0] * in->m[5] * in->m[15] -
in->m[0] * in->m[7] * in->m[13] -
in->m[4] * in->m[1] * in->m[15] +
in->m[4] * in->m[3] * in->m[13] +
in->m[12] * in->m[1] * in->m[7] -
in->m[12] * in->m[3] * in->m[5];
inv[14] = -in->m[0] * in->m[5] * in->m[14] +
in->m[0] * in->m[6] * in->m[13] +
in->m[4] * in->m[1] * in->m[14] -
in->m[4] * in->m[2] * in->m[13] -
in->m[12] * in->m[1] * in->m[6] +
in->m[12] * in->m[2] * in->m[5];
inv[3] = -in->m[1] * in->m[6] * in->m[11] +
in->m[1] * in->m[7] * in->m[10] +
in->m[5] * in->m[2] * in->m[11] -
in->m[5] * in->m[3] * in->m[10] -
in->m[9] * in->m[2] * in->m[7] +
in->m[9] * in->m[3] * in->m[6];
inv[7] = in->m[0] * in->m[6] * in->m[11] -
in->m[0] * in->m[7] * in->m[10] -
in->m[4] * in->m[2] * in->m[11] +
in->m[4] * in->m[3] * in->m[10] +
in->m[8] * in->m[2] * in->m[7] -
in->m[8] * in->m[3] * in->m[6];
inv[11] = -in->m[0] * in->m[5] * in->m[11] +
in->m[0] * in->m[7] * in->m[9] +
in->m[4] * in->m[1] * in->m[11] -
in->m[4] * in->m[3] * in->m[9] -
in->m[8] * in->m[1] * in->m[7] +
in->m[8] * in->m[3] * in->m[5];
inv[15] = in->m[0] * in->m[5] * in->m[10] -
in->m[0] * in->m[6] * in->m[9] -
in->m[4] * in->m[1] * in->m[10] +
in->m[4] * in->m[2] * in->m[9] +
in->m[8] * in->m[1] * in->m[6] -
in->m[8] * in->m[2] * in->m[5];
det = in->m[0] * inv[0] + in->m[1] * inv[4] + in->m[2] * inv[8] + in->m[3] * inv[12];
if (det == 0)
printf("oh oh!");
det = 1.0 / det;
//printf("\n\ndeterminant:%f\n\n",det);
for (i = 0; i < 16; i++)
out->m[i] = inv[i] * det;
}
I create a matrix that looks like this:
[ 1.000000 3.000000 4.000000 9.000000 ]
[ 5.000000 6.000000 7.000000 2.000000 ]
[ 8.000000 8.000000 8.000000 9.000000 ]
[ 0.000000 0.000000 0.000000 1.000000 ]
after passing it to my transpose function that I listed above the results are:
[ -1.000000 1.000000 -0.375000 10.375000 ]
[ 2.000000 -3.000000 1.625000 -26.625000 ]
[ -1.000000 2.000000 -1.125000 15.125000 ]
[ 0.000000 0.000000 0.000000 1.000000 ]
I am not sure but that seems a bit wrong to me. This is for opengl rendering and my matrix is a 1D 16 double array.
Am I getting correct results?
additional info Here is my math over on code review: https://codereview.stackexchange.com/questions/98692/4d-matrix-math-library-for-use-with-opengl