// (Taken from http://savingyoutime.wordpress.com/2009/09/21/c-matrix-inversion-boostublas/) #ifndef INVERSION_H #define INVERSION_H #include #include #include /* Matrix inversion routine. Uses lu_factorize and lu_substitute in uBLAS to invert a matrix */ template bool InvertMatrix(const boost::numeric::ublas::matrix& input, boost::numeric::ublas::matrix& inverse) { using namespace boost::numeric::ublas; typedef permutation_matrix pmatrix; // create a working copy of the input matrix A(input); // create a permutation matrix for the LU-factorization pmatrix pm(A.size1()); // perform LU-factorization int res = lu_factorize(A, pm); if (res != 0) return false; // create identity matrix of "inverse" inverse.assign(identity_matrix (A.size1())); // backsubstitute to get the inverse lu_substitute(A, pm, inverse); return true; } #endif // INVERSION_H