39 #include <visp3/core/vpConfig.h>
41 #include <visp3/core/vpColVector.h>
42 #include <visp3/core/vpMath.h>
43 #include <visp3/core/vpMatrix.h>
45 #ifdef VISP_HAVE_EIGEN3
49 #ifdef VISP_HAVE_LAPACK
51 #include <gsl/gsl_linalg.h>
52 #include <gsl/gsl_permutation.h>
56 typedef MKL_INT integer;
58 #ifdef VISP_HAVE_LAPACK_BUILT_IN
59 typedef long int integer;
63 extern "C" int dgetrf_(integer *m, integer *n,
double *a, integer *lda, integer *ipiv, integer *info);
64 extern "C" void dgetri_(integer *n,
double *a, integer *lda, integer *ipiv,
double *work, integer *lwork,
69 #if (VISP_HAVE_OPENCV_VERSION >= 0x020101)
70 #include <opencv2/core/core.hpp>
74 #include <visp3/core/vpException.h>
75 #include <visp3/core/vpMatrixException.h>
136 if (std::fabs(d) < std::numeric_limits<double>::epsilon()) {
146 if (std::fabs(d) < std::numeric_limits<double>::epsilon()) {
151 inv[1][1] = (*this)[0][0] * d;
152 inv[0][0] = (*this)[1][1] * d;
153 inv[0][1] = -(*this)[0][1] * d;
154 inv[1][0] = -(*this)[1][0] * d;
160 if (std::fabs(d) < std::numeric_limits<double>::epsilon()) {
165 inv[0][0] = ((*this)[1][1] * (*this)[2][2] - (*this)[1][2] * (*this)[2][1]) * d;
166 inv[0][1] = ((*this)[0][2] * (*this)[2][1] - (*this)[0][1] * (*this)[2][2]) * d;
167 inv[0][2] = ((*this)[0][1] * (*this)[1][2] - (*this)[0][2] * (*this)[1][1]) * d;
169 inv[1][0] = ((*this)[1][2] * (*this)[2][0] - (*this)[1][0] * (*this)[2][2]) * d;
170 inv[1][1] = ((*this)[0][0] * (*this)[2][2] - (*this)[0][2] * (*this)[2][0]) * d;
171 inv[1][2] = ((*this)[0][2] * (*this)[1][0] - (*this)[0][0] * (*this)[1][2]) * d;
173 inv[2][0] = ((*this)[1][0] * (*this)[2][1] - (*this)[1][1] * (*this)[2][0]) * d;
174 inv[2][1] = ((*this)[0][1] * (*this)[2][0] - (*this)[0][0] * (*this)[2][1]) * d;
175 inv[2][2] = ((*this)[0][0] * (*this)[1][1] - (*this)[0][1] * (*this)[1][0]) * d;
178 #if defined(VISP_HAVE_LAPACK)
180 #elif defined(VISP_HAVE_EIGEN3)
182 #elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
186 "Install Lapack, Eigen3 or OpenCV 3rd party"));
227 return (*
this)[0][0];
229 return ((*
this)[0][0] * (*
this)[1][1] - (*
this)[0][1] * (*
this)[1][0]);
231 return ((*
this)[0][0] * ((*
this)[1][1] * (*
this)[2][2] - (*
this)[1][2] * (*
this)[2][1]) -
232 (*
this)[0][1] * ((*
this)[1][0] * (*
this)[2][2] - (*
this)[1][2] * (*
this)[2][0]) +
233 (*
this)[0][2] * ((*
this)[1][0] * (*
this)[2][1] - (*
this)[1][1] * (*
this)[2][0]));
235 #if defined(VISP_HAVE_LAPACK)
237 #elif defined(VISP_HAVE_EIGEN3)
239 #elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
243 "Install Lapack, Eigen3 or OpenCV 3rd party"));
248 #if defined(VISP_HAVE_LAPACK)
281 #if defined(VISP_HAVE_GSL)
290 unsigned int tda = (
unsigned int)A->tda;
291 for (
unsigned int i = 0; i <
rowNum; i++) {
292 unsigned int k = i * tda;
293 for (
unsigned int j = 0; j <
colNum; j++)
294 A->data[k + j] = (*
this)[i][j];
302 inverse.tda = inverse.size2;
303 inverse.data = Ainv.
data;
307 gsl_permutation *p = gsl_permutation_alloc(
rowNum);
311 gsl_linalg_LU_decomp(A, p, &s);
312 gsl_linalg_LU_invert(A, p, &inverse);
314 gsl_permutation_free(p);
325 integer dim = (integer)
rowNum;
328 integer lwork = dim * dim;
329 integer *ipiv =
new integer[dim + 1];
330 double *work =
new double[lwork];
334 dgetrf_(&dim, &dim, A.
data, &lda, &ipiv[1], &info);
341 dgetri_(&dim, A.
data, &dim, &ipiv[1], work, &lwork, &info);
378 #if defined(VISP_HAVE_GSL)
390 unsigned int tda = (
unsigned int)A->tda;
391 for (
unsigned int i = 0; i <
rowNum; i++) {
392 unsigned int k = i * tda;
393 for (
unsigned int j = 0; j <
colNum; j++)
394 A->data[k + j] = (*
this)[i][j];
397 gsl_permutation *p = gsl_permutation_alloc(
rowNum);
401 gsl_linalg_LU_decomp(A, p, &s);
402 det = gsl_linalg_LU_det(A, s);
404 gsl_permutation_free(p);
416 integer dim = (integer)
rowNum;
419 integer *ipiv =
new integer[dim + 1];
423 dgetrf_(&dim, &dim, A.
data, &lda, &ipiv[1], &info);
429 double det = A[0][0];
430 for (
unsigned int i = 1; i <
rowNum; i++) {
435 for (
int i = 1; i <= dim; i++) {
450 #if (VISP_HAVE_OPENCV_VERSION >= 0x020101)
489 cv::Mat Minv = M.inv(cv::DECOMP_LU);
492 memcpy(A.
data, Minv.data, (
size_t)(8 * Minv.rows * Minv.cols));
532 det = cv::determinant(M);
538 #if defined(VISP_HAVE_EIGEN3)
577 Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> > M(this->
data, this->
getRows(),
579 Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> > A_(A.
data, this->getRows(),
619 Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> > M(this->
data, this->
getRows(),
622 return M.determinant();
unsigned int getCols() const
Type * data
Address of the first element of the data array.
void resize(unsigned int nrows, unsigned int ncols, bool flagNullify=true, bool recopy_=true)
unsigned int rowNum
Number of rows in the array.
unsigned int getRows() const
unsigned int colNum
Number of columns in the array.
error that can be emited by ViSP classes.
Implementation of a matrix and operations on matrices.
vpMatrix inverseByLU() const
vpMatrix inverseByLUEigen3() const
vpMatrix inverseByLUOpenCV() const
double detByLUEigen3() const
double detByLUOpenCV() const
vpMatrix inverseByLULapack() const
double det(vpDetMethod method=LU_DECOMPOSITION) const
double detByLULapack() const