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()) {
150 inv[1][1] = (*this)[0][0]*d;
151 inv[0][0] = (*this)[1][1]*d;
152 inv[0][1] = -(*this)[0][1]*d;
153 inv[1][0] = -(*this)[1][0]*d;
160 if (std::fabs(d) < std::numeric_limits<double>::epsilon()) {
164 inv[0][0] = ((*this)[1][1] * (*this)[2][2] - (*this)[1][2] * (*this)[2][1]) * d;
165 inv[0][1] = ((*this)[0][2] * (*this)[2][1] - (*this)[0][1] * (*this)[2][2]) * d;
166 inv[0][2] = ((*this)[0][1] * (*this)[1][2] - (*this)[0][2] * (*this)[1][1]) * d;
168 inv[1][0] = ((*this)[1][2] * (*this)[2][0] - (*this)[1][0] * (*this)[2][2]) * d;
169 inv[1][1] = ((*this)[0][0] * (*this)[2][2] - (*this)[0][2] * (*this)[2][0]) * d;
170 inv[1][2] = ((*this)[0][2] * (*this)[1][0] - (*this)[0][0] * (*this)[1][2]) * d;
172 inv[2][0] = ((*this)[1][0] * (*this)[2][1] - (*this)[1][1] * (*this)[2][0]) * d;
173 inv[2][1] = ((*this)[0][1] * (*this)[2][0] - (*this)[0][0] * (*this)[2][1]) * d;
174 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];
230 return ((*
this)[0][0] * (*
this)[1][1] - (*
this)[0][1] * (*
this)[1][0]);
233 return ((*
this)[0][0] * ((*
this)[1][1] * (*
this)[2][2] - (*
this)[1][2] * (*
this)[2][1]) -
234 (*
this)[0][1] * ((*
this)[1][0] * (*
this)[2][2] - (*
this)[1][2] * (*
this)[2][0]) +
235 (*
this)[0][2] * ((*
this)[1][0] * (*
this)[2][1] - (*
this)[1][1] * (*
this)[2][0]));
237 #if defined(VISP_HAVE_LAPACK)
239 #elif defined(VISP_HAVE_EIGEN3)
241 #elif (VISP_HAVE_OPENCV_VERSION >= 0x020101)
245 "Install Lapack, Eigen3 or OpenCV 3rd party"));
250 #if defined(VISP_HAVE_LAPACK)
283 #if defined(VISP_HAVE_GSL)
292 unsigned int tda = (
unsigned int)A->tda;
293 for (
unsigned int i = 0; i <
rowNum; i++) {
294 unsigned int k = i * tda;
295 for (
unsigned int j = 0; j <
colNum; j++)
296 A->data[k + j] = (*
this)[i][j];
304 inverse.tda = inverse.size2;
305 inverse.data = Ainv.
data;
309 gsl_permutation *p = gsl_permutation_alloc(
rowNum);
313 gsl_linalg_LU_decomp(A, p, &s);
314 gsl_linalg_LU_invert(A, p, &inverse);
316 gsl_permutation_free(p);
327 integer dim = (integer)
rowNum;
330 integer lwork = dim * dim;
331 integer *ipiv =
new integer[dim + 1];
332 double *work =
new double[lwork];
336 dgetrf_(&dim, &dim, A.
data, &lda, &ipiv[1], &info);
343 dgetri_(&dim, A.
data, &dim, &ipiv[1], work, &lwork, &info);
380 #if defined(VISP_HAVE_GSL)
392 unsigned int tda = (
unsigned int)A->tda;
393 for (
unsigned int i = 0; i <
rowNum; i++) {
394 unsigned int k = i * tda;
395 for (
unsigned int j = 0; j <
colNum; j++)
396 A->data[k + j] = (*
this)[i][j];
399 gsl_permutation *p = gsl_permutation_alloc(
rowNum);
403 gsl_linalg_LU_decomp(A, p, &s);
404 det = gsl_linalg_LU_det(A, s);
406 gsl_permutation_free(p);
418 integer dim = (integer)
rowNum;
421 integer *ipiv =
new integer[dim + 1];
425 dgetrf_(&dim, &dim, A.
data, &lda, &ipiv[1], &info);
431 double det = A[0][0];
432 for (
unsigned int i = 1; i <
rowNum; i++) {
437 for (
int i = 1; i <= dim; i++) {
452 #if (VISP_HAVE_OPENCV_VERSION >= 0x020101)
491 cv::Mat Minv = M.inv(cv::DECOMP_LU);
494 memcpy(A.
data, Minv.data, (
size_t)(8 * Minv.rows * Minv.cols));
534 det = cv::determinant(M);
540 #if defined(VISP_HAVE_EIGEN3)
579 Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> > M(this->
data, this->
getRows(),
581 Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> > A_(A.
data, this->getRows(),
621 Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> > M(this->
data, this->
getRows(),
624 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