Gauss

int Matrix::inverse()
{
    int d, i, j, k;
    for(j = 0; j < RS_K; j++)
        for(k = 0; k < RS_K; k++)
            inv[j][k] = j==k ? 1 : 0;
    RSD f;
    for(d = 0; d < RS_K; d++)
    {
        k = d;
        if(!a[d][k])
        {
            for(j = d+1; !a[j][k] && j < RS_K; j++);
            if(j==RS_K)
                return RS_K-d;
            swap_lines(a, d, j);
            swap_lines(inv, d, j);
        }
        f = a[d][k];
        for(i = 0; i < RS_K; i++)
            a[d][i] /= f, inv[d][i] /= f;
        for(j = 0; j < RS_K; j++)
        {
            if(j == d)  continue;
            f = a[j][d];
            for(k = 0; k < RS_K; k++)
                a[j][k] -= f*a[d][k], inv[j][k] -= f*inv[d][k];
        }
    }
    return 0;
}