4 by 4 matrix inverse, need help.

[serengeor]

I was writing my own 4 by 4 matrix class and started writing inversion code. Though it seams that I can't get it coded right.
It seams that I only calculate the determinant right and the 1 value of inverse matrix (Also 14th value seams right, only that it has the wrong sign).
Here is the code for inverse function:

Code: Select all

 
 
    bool inverse()
    {
        /**
        indexes for array
            1   2   3   4
         _____________________
       1 |   0   1   2   3   |
       2 |   4   5   6   7   |
       3 |   8   9   10  11  |
       4 |   12  13  14  15  |
        **/
        T a00=          M[0]* ( M[5]*( (M[10]*M[15]) - (M[11]*M[14])) - M[6] * ( M[9]*M[15] - M[11]*M[13] ) + M[7] *(M[9]*M[14] - M[10]*M[13]) );
        T a01=-1*       M[1]* ( M[4]*( (M[10]*M[15]) - (M[11]*M[14])) - M[6] * ( M[8]*M[15] - M[11]*M[12] ) + M[7] *(M[8]*M[14] - M[10]*M[12]) );
        T a02=          M[2]* ( M[4]*( (M[9] *M[15]) - (M[11]*M[13])) - M[5] * ( M[8]*M[15] - M[11]*M[12] ) + M[7] *(M[8]*M[13] - M[9] *M[12]) );
        T a03=-1*       M[3]* ( M[4]*( (M[9] *M[14]) - (M[10]*M[13])) - M[5] * ( M[8]*M[14] - M[10]*M[12] ) + M[6] *(M[8]*M[13] - M[9] *M[12]) );
        T det = a00 + a01 + a02 + a03;
 
        printf("Determinant=%f\n",(float)det);
 
        if(det!=(T)0)
        {
            Matrix4<T> m;
            m.M[0]=a00;
            m.M[1]=a01;
            m.M[2]=a02;
            m.M[3]=a03;
 
            m.M[4]=-1*    M[4]* ( M[1] * ((M[10]*M[15]) - (M[11]*M[14])) - M[2] * ( (M[9]*M[15]) - (M[11]*M[13]) ) + M[3] *( (M[9]*M[14]) - (M[10]*M[13]) ) );
            m.M[5]=       M[5]* ( M[0] * ((M[10]*M[15]) - (M[11]*M[14])) - M[2] * ( (M[8]*M[15]) - (M[11]*M[12]) ) + M[3] *( (M[8]*M[14]) - (M[10]*M[12]) ) );
            m.M[6]=-1*    M[6]* ( M[0] * ((M[9] *M[15]) - (M[11]*M[13])) - M[1] * ( (M[8]*M[15]) - (M[11]*M[12]) ) + M[3] *( (M[8]*M[13]) - (M[9] *M[12]) ) );
            m.M[7]=       M[7]* ( M[0] * ((M[9] *M[14]) - (M[10]*M[13])) - M[1] * ( (M[8]*M[14]) - (M[10]*M[12]) ) + M[2] *( (M[8]*M[13]) - (M[9] *M[12]) ) );
 
            m.M[8]=       M[8]* ( M[1] * ((M[6] *M[15]) - (M[7] *M[14])) - M[2] * ( (M[5]*M[15]) - (M[7] *M[13]) ) + M[3] *( (M[5]*M[14]) - (M[6] *M[13]) ) );
            m.M[9]=-1*    M[9]* ( M[0] * ((M[6] *M[15]) - (M[7] *M[14])) - M[2] * ( (M[4]*M[15]) - (M[7] *M[12]) ) + M[3] *( (M[4]*M[14]) - (M[6] *M[12]) ) );
            m.M[10]=      M[10]*( M[0] * ((M[5] *M[15]) - (M[7] *M[13])) - M[1] * ( (M[4]*M[15]) - (M[7] *M[12]) ) + M[3] *( (M[4]*M[13]) - (M[5] *M[12]) ) );
            m.M[11]=-1*   M[11]*( M[0] * ((M[5] *M[14]) - (M[6] *M[13])) - M[1] * ( (M[4]*M[14]) - (M[6] *M[12]) ) + M[2] *( (M[4]*M[13]) - (M[5] *M[12]) ) );
 
 
            m.M[12]=-1*   M[12]*( M[1] * ((M[6] *M[11]) - (M[7] *M[10])) - M[2] * ( (M[5]*M[11]) - (M[7] *M[9])  ) + M[3] *( (M[5]*M[10]) - (M[6] *M[9]) ) );
            m.M[13]=      M[13]*( M[0] * ((M[6] *M[11]) - (M[7] *M[10])) - M[2] * ( (M[4]*M[11]) - (M[7] *M[8])  ) + M[3] *( (M[4]*M[10]) - (M[6] *M[8]) ) );
            m.M[14]=-1*   M[14]*( M[0] * ((M[5] *M[11]) - (M[7] *M[9] )) - M[1] * ( (M[4]*M[11]) - (M[7] *M[8])  ) + M[3] *( (M[4]*M[9] ) - (M[5] *M[8]) ) );
            m.M[15]=      M[15]*( M[0] * ((M[5] *M[10]) - (M[6] *M[9] )) - M[1] * ( (M[4]*M[10]) - (M[6] *M[8])  ) + M[2] *( (M[4]*M[9] ) - (M[5] *M[8]) ) );
            m.transpose();
            det=(T)1.0/det;
            *this=(m*det);
            return true;
        }
        return false;
    }
 
 

Here is the code for the whole matrix class: http://pastebin.com/vYU2cKzi
Here is the test code: http://pastebin.com/tStuDzmt

I'm sure if you could find my mistakes in a few of those lines I would work my way out

[serengeor]

Sorry for wasting your time , I just found my stupid mistake.
I had made all calculations almost correctly, except I didn't need to multiply some of the values like I did in determinant calculation.
The corrected code now looks like this:

Code: Select all

 
bool inverse()
    {
        /**
        indexes for array
            1   2   3   4
         _____________________
       1 |   0   1   2   3   |
       2 |   4   5   6   7   |
       3 |   8   9   10  11  |
       4 |   12  13  14  15  |
        **/
        T a00=   ( M[5]*( (M[10]*M[15]) - (M[11]*M[14])) - M[6] * ( M[9]*M[15] - M[11]*M[13] ) + M[7] *(M[9]*M[14] - M[10]*M[13]) );
        T a01=-1*( M[4]*( (M[10]*M[15]) - (M[11]*M[14])) - M[6] * ( M[8]*M[15] - M[11]*M[12] ) + M[7] *(M[8]*M[14] - M[10]*M[12]) );
        T a02=   ( M[4]*( (M[9] *M[15]) - (M[11]*M[13])) - M[5] * ( M[8]*M[15] - M[11]*M[12] ) + M[7] *(M[8]*M[13] - M[9] *M[12]) );
        T a03=-1*( M[4]*( (M[9] *M[14]) - (M[10]*M[13])) - M[5] * ( M[8]*M[14] - M[10]*M[12] ) + M[6] *(M[8]*M[13] - M[9] *M[12]) );
        T det = a00*M[0] + a01*M[1] + a02*M[2] + a03*M[3];
 
        printf("Determinant=%f\n",(float)det);
 
        if(det!=(T)0)
        {
            Matrix4<T> m;
            m.M[0]=a00;
            m.M[1]=a01;
            m.M[2]=a02;
            m.M[3]=a03;
 
            m.M[4]=-1*    ( M[1] * ((M[10]*M[15]) - (M[11]*M[14])) - M[2] * ( (M[9]*M[15]) - (M[11]*M[13]) ) + M[3] *( (M[9]*M[14]) - (M[10]*M[13]) ) );
            //m.M[4]=-1*    M[4]* ( (M[1] * M[10]*M[15]) + (M[9]*M[14]*M[3]) + (M[2]*M[11]*M[13]) - (M[3]*M[10]*M[13]) - (M[9]*M[2]*M[15]) - (M[14]*M[11]*M[1]));
            m.M[5]=       ( M[0] * ((M[10]*M[15]) - (M[11]*M[14])) - M[2] * ( (M[8]*M[15]) - (M[11]*M[12]) ) + M[3] *( (M[8]*M[14]) - (M[10]*M[12]) ) );
            m.M[6]=-1*    ( M[0] * ((M[9] *M[15]) - (M[11]*M[13])) - M[1] * ( (M[8]*M[15]) - (M[11]*M[12]) ) + M[3] *( (M[8]*M[13]) - (M[9] *M[12]) ) );
            m.M[7]=       ( M[0] * ((M[9] *M[14]) - (M[10]*M[13])) - M[1] * ( (M[8]*M[14]) - (M[10]*M[12]) ) + M[2] *( (M[8]*M[13]) - (M[9] *M[12]) ) );
 
            m.M[8]=       ( M[1] * ((M[6] *M[15]) - (M[7] *M[14])) - M[2] * ( (M[5]*M[15]) - (M[7] *M[13]) ) + M[3] *( (M[5]*M[14]) - (M[6] *M[13]) ) );
            m.M[9]=-1*    ( M[0] * ((M[6] *M[15]) - (M[7] *M[14])) - M[2] * ( (M[4]*M[15]) - (M[7] *M[12]) ) + M[3] *( (M[4]*M[14]) - (M[6] *M[12]) ) );
            m.M[10]=      ( M[0] * ((M[5] *M[15]) - (M[7] *M[13])) - M[1] * ( (M[4]*M[15]) - (M[7] *M[12]) ) + M[3] *( (M[4]*M[13]) - (M[5] *M[12]) ) );
            m.M[11]=-1*   ( M[0] * ((M[5] *M[14]) - (M[6] *M[13])) - M[1] * ( (M[4]*M[14]) - (M[6] *M[12]) ) + M[2] *( (M[4]*M[13]) - (M[5] *M[12]) ) );
 
 
            m.M[12]=-1*   ( M[1] * ((M[6] *M[11]) - (M[7] *M[10])) - M[2] * ( (M[5]*M[11]) - (M[7] *M[9])  ) + M[3] *( (M[5]*M[10]) - (M[6] *M[9]) ) );
            m.M[13]=      ( M[0] * ((M[6] *M[11]) - (M[7] *M[10])) - M[2] * ( (M[4]*M[11]) - (M[7] *M[8])  ) + M[3] *( (M[4]*M[10]) - (M[6] *M[8]) ) );
            m.M[14]=-1*   ( M[0] * ((M[5] *M[11]) - (M[7] *M[9] )) - M[1] * ( (M[4]*M[11]) - (M[7] *M[8])  ) + M[3] *( (M[4]*M[9] ) - (M[5] *M[8]) ) );
            m.M[15]=      ( M[0] * ((M[5] *M[10]) - (M[6] *M[9] )) - M[1] * ( (M[4]*M[10]) - (M[6] *M[8])  ) + M[2] *( (M[4]*M[9] ) - (M[5] *M[8]) ) );
            m.transpose();
            det=(T)1.0/det;
            *this=(m*det);
            return true;
        }
        return false;
    }
 

Maybe it will be useful to somebody.