52 #include <visp3/core/vpCPUFeatures.h>
53 #include <visp3/core/vpColVector.h>
54 #include <visp3/core/vpConfig.h>
55 #include <visp3/core/vpDebug.h>
56 #include <visp3/core/vpException.h>
57 #include <visp3/core/vpMath.h>
58 #include <visp3/core/vpMatrix.h>
59 #include <visp3/core/vpTranslationVector.h>
61 #if defined(VISP_HAVE_SIMDLIB)
62 #include <Simd/SimdLib.h>
65 #ifdef VISP_HAVE_LAPACK
67 #include <gsl/gsl_eigen.h>
68 #include <gsl/gsl_linalg.h>
69 #include <gsl/gsl_math.h>
70 #elif defined(VISP_HAVE_MKL)
72 typedef MKL_INT integer;
74 void vpMatrix::blas_dsyev(
char jobz,
char uplo,
unsigned int n_,
double *a_data,
unsigned int lda_,
double *w_data,
75 double *work_data,
int lwork_,
int &info_)
77 MKL_INT n =
static_cast<MKL_INT
>(n_);
78 MKL_INT lda =
static_cast<MKL_INT
>(lda_);
79 MKL_INT lwork =
static_cast<MKL_INT
>(lwork_);
80 MKL_INT info =
static_cast<MKL_INT
>(info_);
82 dsyev(&jobz, &uplo, &n, a_data, &lda, w_data, work_data, &lwork, &info);
86 #if defined(VISP_HAVE_LAPACK_BUILT_IN)
87 typedef long int integer;
91 extern "C" integer dsyev_(
char *jobz,
char *uplo, integer *n,
double *a, integer *lda,
double *w,
double *WORK,
92 integer *lwork, integer *info);
94 void vpMatrix::blas_dsyev(
char jobz,
char uplo,
unsigned int n_,
double *a_data,
unsigned int lda_,
double *w_data,
95 double *work_data,
int lwork_,
int &info_)
97 integer n =
static_cast<integer
>(n_);
98 integer lda =
static_cast<integer
>(lda_);
99 integer lwork =
static_cast<integer
>(lwork_);
100 integer info =
static_cast<integer
>(info_);
102 dsyev_(&jobz, &uplo, &n, a_data, &lda, w_data, work_data, &lwork, &info);
104 lwork_ =
static_cast<int>(lwork);
105 info_ =
static_cast<int>(info);
110 #if !defined(VISP_USE_MSVC) || (defined(VISP_USE_MSVC) && !defined(VISP_BUILD_SHARED_LIBS))
111 const unsigned int vpMatrix::m_lapack_min_size_default = 0;
112 unsigned int vpMatrix::m_lapack_min_size = vpMatrix::m_lapack_min_size_default;
119 unsigned int ncols,
double svThreshold,
vpMatrix &Ap,
int &rank_out,
int *rank_in,
122 Ap.
resize(ncols, nrows,
true,
false);
125 double maxsv = sv[0];
129 unsigned int sv_size = sv.
size();
130 for (
unsigned int i = 0; i < sv_size; ++i) {
131 if (sv[i] > (maxsv * svThreshold)) {
136 unsigned int rank =
static_cast<unsigned int>(rank_out);
138 rank =
static_cast<unsigned int>(*rank_in);
141 for (
unsigned int i = 0; i < ncols; ++i) {
142 for (
unsigned int j = 0; j < nrows; ++j) {
143 for (
unsigned int k = 0; k < rank; ++k) {
144 Ap[i][j] += (V[i][k] * U[j][k]) / sv[k];
153 for (
unsigned int i = 0; i < nrows; ++i) {
154 for (
unsigned int j = 0; j < rank; ++j) {
155 (*imA)[i][j] = U[i][j];
162 imAt->
resize(ncols, rank);
163 for (
unsigned int i = 0; i < ncols; ++i) {
164 for (
unsigned int j = 0; j < rank; ++j) {
165 (*imAt)[i][j] = V[i][j];
172 kerAt->
resize(ncols - rank, ncols);
174 unsigned int v_rows = V.
getRows();
175 for (
unsigned int k = 0; k < (ncols - rank); ++k) {
176 unsigned j = k + rank;
177 for (
unsigned int i = 0; i < v_rows; ++i) {
178 (*kerAt)[k][i] = V[i][j];
193 if (((r + nrows) > M.
rowNum) || ((c + ncols) > M.
colNum)) {
195 "Cannot construct a sub matrix (%dx%d) starting at "
196 "position (%d,%d) that is not contained in the "
197 "original matrix (%dx%d)",
201 init(M, r, c, nrows, ncols);
204 #if (VISP_CXX_STANDARD >= VISP_CXX_STANDARD_11)
267 vpMatrix::vpMatrix(
unsigned int nrows,
unsigned int ncols,
const std::initializer_list<double> &list)
342 unsigned int rnrows = r + nrows;
343 unsigned int cncols = c + ncols;
353 resize(nrows, ncols,
false,
false);
355 if (this->
rowPtrs ==
nullptr) {
358 for (
unsigned int i = 0; i < nrows; ++i) {
359 memcpy((*
this)[i], &M[i + r][c], ncols *
sizeof(
double));
406 unsigned int rnrows = r + nrows;
407 unsigned int cncols = c + ncols;
419 M.
resize(nrows, ncols,
false,
false);
420 for (
unsigned int i = 0; i < nrows; ++i) {
421 memcpy(M[i], &(*
this)[i + r][c], ncols *
sizeof(
double));
450 for (
unsigned int i = 0; i <
rowNum; ++i) {
451 for (
unsigned int j = 0; j <
colNum; ++j) {
489 for (
unsigned int i = 0; i <
rowNum; ++i) {
490 for (
unsigned int j = 0; j <
colNum; ++j) {
491 At[j][i] = (*this)[i][j];
496 #if defined(VISP_HAVE_SIMDLIB)
500 const int tileSize = 32;
501 for (
unsigned int i = 0; i <
rowNum; i += tileSize) {
502 for (
unsigned int j = 0; j <
colNum; j++) {
503 for (
unsigned int b = 0; b < static_cast<unsigned int>(tileSize) && i + b <
rowNum; b++) {
504 At[j][i + b] = (*this)[i + b][j];
544 bool useLapack = ((
rowNum > vpMatrix::m_lapack_min_size) || (
colNum > vpMatrix::m_lapack_min_size));
545 #if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
550 #if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
551 const double alpha = 1.0;
552 const double beta = 0.0;
553 const char transa =
't';
554 const char transb =
'n';
556 vpMatrix::blas_dgemm(transa, transb,
rowNum,
rowNum,
colNum, alpha,
data,
colNum,
data,
colNum, beta, B.
data,
562 for (
unsigned int i = 0; i <
rowNum; ++i) {
563 for (
unsigned int j = i; j <
rowNum; ++j) {
569 for (
unsigned int k = 0; k <
colNum; ++k) {
570 ssum += *(pi++) * *(pj++);
600 bool useLapack = ((
rowNum > vpMatrix::m_lapack_min_size) || (
colNum > vpMatrix::m_lapack_min_size));
601 #if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
606 #if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
607 const double alpha = 1.0;
608 const double beta = 0.0;
609 const char transa =
'n';
610 const char transb =
't';
612 vpMatrix::blas_dgemm(transa, transb,
colNum,
colNum,
rowNum, alpha,
data,
colNum,
data,
colNum, beta, B.
data,
617 for (
unsigned int i = 0; i <
colNum; ++i) {
619 for (
unsigned int j = 0; j < i; ++j) {
622 for (
unsigned int k = 0; k <
rowNum; ++k) {
623 s += (*(ptr + i)) * (*(ptr + j));
631 for (
unsigned int k = 0; k <
rowNum; ++k) {
632 s += (*(ptr + i)) * (*(ptr + i));
692 #if (VISP_CXX_STANDARD >= VISP_CXX_STANDARD_11)
695 if (
this != &other) {
711 other.rowPtrs =
nullptr;
713 other.data =
nullptr;
747 if (
dsize !=
static_cast<unsigned int>(list.size())) {
748 resize(1,
static_cast<unsigned int>(list.size()),
false,
false);
751 std::copy(list.begin(), list.end(),
data);
781 unsigned int nrows =
static_cast<unsigned int>(lists.size()), ncols = 0;
782 for (
auto &l : lists) {
783 if (
static_cast<unsigned int>(l.size()) > ncols) {
784 ncols =
static_cast<unsigned int>(l.size());
788 resize(nrows, ncols,
false,
false);
789 auto it = lists.begin();
790 for (
unsigned int i = 0; i <
rowNum; ++i, ++it) {
791 std::copy(it->begin(), it->end(),
rowPtrs[i]);
812 for (
unsigned int i = 0; i <
rowNum; ++i) {
813 for (
unsigned int j = 0; j <
colNum; ++j) {
822 resize(1, 1,
false,
false);
871 unsigned int rows = A.
getRows();
875 for (
unsigned int i = 0; i < rows; ++i) {
876 (*this)[i][i] = A[i];
914 for (
unsigned int i = 0; i < min_; ++i) {
931 unsigned int rows = A.
getRows();
932 DA.
resize(rows, rows,
true);
934 for (
unsigned int i = 0; i < rows; ++i) {
952 for (
unsigned int j = 0; j < 3; ++j) {
956 for (
unsigned int j = 0; j < 3; ++j) {
958 for (
unsigned int i = 0; i < 3; ++i) {
959 t_out[i] +=
rowPtrs[i][j] * tj;
996 bool useLapack = ((A.
rowNum > vpMatrix::m_lapack_min_size) || (A.
colNum > vpMatrix::m_lapack_min_size));
997 #if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
1002 #if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
1008 vpMatrix::blas_dgemv(trans, A.
colNum, A.
rowNum, alpha, A.
data, A.
colNum, v.
data, incr, beta, w.
data, incr);
1013 for (
unsigned int j = 0; j < A.
colNum; ++j) {
1015 for (
unsigned int i = 0; i < A.
rowNum; ++i) {
1047 bool useLapack = ((A.
rowNum > vpMatrix::m_lapack_min_size) || (A.
colNum > vpMatrix::m_lapack_min_size) ||
1048 (B.
colNum > vpMatrix::m_lapack_min_size));
1049 #if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
1054 #if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
1055 const double alpha = 1.0;
1056 const double beta = 0.0;
1057 const char trans =
'n';
1058 vpMatrix::blas_dgemm(trans, trans, B.
colNum, A.
rowNum, A.
colNum, alpha, B.
data, B.
colNum, A.
data, A.
colNum, beta,
1064 const unsigned int BcolNum = B.
colNum;
1065 const unsigned int BrowNum = B.
rowNum;
1066 double **BrowPtrs = B.
rowPtrs;
1067 for (
unsigned int i = 0; i < A.
rowNum; ++i) {
1068 const double *rowptri = A.
rowPtrs[i];
1070 for (
unsigned int j = 0; j < BcolNum; ++j) {
1072 for (
unsigned int k = 0; k < BrowNum; ++k) {
1073 s += rowptri[k] * BrowPtrs[k][j];
1097 "Cannot multiply (%dx%d) matrix by (%dx%d) matrix as a "
1102 const unsigned int BcolNum = B.
colNum;
1103 const unsigned int BrowNum = B.
rowNum;
1104 double **BrowPtrs = B.
rowPtrs;
1105 for (
unsigned int i = 0; i < A.
rowNum; ++i) {
1106 const double *rowptri = A.
rowPtrs[i];
1108 for (
unsigned int j = 0; j < BcolNum; ++j) {
1110 for (
unsigned int k = 0; k < BrowNum; ++k) {
1111 s += rowptri[k] * BrowPtrs[k][j];
1134 "Cannot multiply (%dx%d) matrix by (%dx%d) matrix as a "
1143 bool useLapack = ((A.
rowNum > vpMatrix::m_lapack_min_size) || (A.
colNum > vpMatrix::m_lapack_min_size) ||
1144 (B.
colNum > vpMatrix::m_lapack_min_size));
1145 #if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
1150 #if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
1151 const double alpha = 1.0;
1152 const double beta = 0.0;
1153 const char trans =
'n';
1154 vpMatrix::blas_dgemm(trans, trans, B.
colNum, A.
rowNum, A.
colNum, alpha, B.
data, B.
colNum, A.
data, A.
colNum, beta,
1160 const unsigned int BcolNum = B.
colNum;
1161 const unsigned int BrowNum = B.
rowNum;
1162 double **BrowPtrs = B.
rowPtrs;
1163 for (
unsigned int i = 0; i < A.
rowNum; ++i) {
1164 const double *rowptri = A.
rowPtrs[i];
1166 for (
unsigned int j = 0; j < BcolNum; ++j) {
1168 for (
unsigned int k = 0; k < BrowNum; ++k) {
1169 s += rowptri[k] * BrowPtrs[k][j];
1220 unsigned int RcolNum = R.
getCols();
1221 unsigned int RrowNum = R.
getRows();
1222 for (
unsigned int i = 0; i <
rowNum; ++i) {
1225 for (
unsigned int j = 0; j < RcolNum; ++j) {
1227 for (
unsigned int k = 0; k < RrowNum; ++k) {
1228 s += rowptri[k] * R[k][j];
1250 const unsigned int McolNum = M.
getCols();
1251 const unsigned int MrowNum = M.
getRows();
1252 for (
unsigned int i = 0; i <
rowNum; ++i) {
1253 const double *rowptri =
rowPtrs[i];
1255 for (
unsigned int j = 0; j < McolNum; ++j) {
1257 for (
unsigned int k = 0; k < MrowNum; ++k) {
1258 s += rowptri[k] * M[k][j];
1285 bool useLapack = ((
rowNum > vpMatrix::m_lapack_min_size) || (
colNum > vpMatrix::m_lapack_min_size) ||
1286 (V.
colNum > vpMatrix::m_lapack_min_size));
1287 #if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
1292 #if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
1293 const double alpha = 1.0;
1294 const double beta = 0.0;
1295 const char trans =
'n';
1296 vpMatrix::blas_dgemm(trans, trans, V.
colNum,
rowNum,
colNum, alpha, V.
data, V.
colNum,
data,
colNum, beta, M.
data,
1301 #if defined(VISP_HAVE_SIMDLIB)
1304 unsigned int VcolNum = V.
getCols();
1305 unsigned int VrowNum = V.
getRows();
1306 for (
unsigned int i = 0; i <
rowNum; ++i) {
1309 for (
unsigned int j = 0; j < VcolNum; ++j) {
1311 for (
unsigned int k = 0; k < VrowNum; ++k) {
1312 s += rowptri[k] * V[k][j];
1341 bool useLapack = ((
rowNum > vpMatrix::m_lapack_min_size) || (
colNum > vpMatrix::m_lapack_min_size) ||
1342 (V.
getCols() > vpMatrix::m_lapack_min_size));
1343 #if !(defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL))
1348 #if defined(VISP_HAVE_LAPACK) && !defined(VISP_HAVE_LAPACK_BUILT_IN) && !defined(VISP_HAVE_GSL)
1349 const double alpha = 1.0;
1350 const double beta = 0.0;
1351 const char trans =
'n';
1352 vpMatrix::blas_dgemm(trans, trans, V.
getCols(),
rowNum,
colNum, alpha, V.
data, V.
getCols(),
data,
colNum, beta,
1357 #if defined(VISP_HAVE_SIMDLIB)
1360 unsigned int VcolNum = V.
getCols();
1361 unsigned int VrowNum = V.
getRows();
1362 for (
unsigned int i = 0; i <
rowNum; ++i) {
1365 for (
unsigned int j = 0; j < VcolNum; ++j) {
1367 for (
unsigned int k = 0; k < VrowNum; ++k) {
1368 s += rowptri[k] * V[k][j];
1401 double **ArowPtrs = A.
rowPtrs;
1402 double **BrowPtrs = B.
rowPtrs;
1403 double **CrowPtrs = C.
rowPtrs;
1405 for (
unsigned int i = 0; i < A.
rowNum; ++i) {
1406 for (
unsigned int j = 0; j < A.
colNum; ++j) {
1407 CrowPtrs[i][j] = (wB * BrowPtrs[i][j]) + (wA * ArowPtrs[i][j]);
1432 double **ArowPtrs = A.
rowPtrs;
1433 double **BrowPtrs = B.
rowPtrs;
1434 double **CrowPtrs = C.
rowPtrs;
1436 for (
unsigned int i = 0; i < A.
rowNum; ++i) {
1437 for (
unsigned int j = 0; j < A.
colNum; ++j) {
1438 CrowPtrs[i][j] = BrowPtrs[i][j] + ArowPtrs[i][j];
1466 double **ArowPtrs = A.
rowPtrs;
1467 double **BrowPtrs = B.
rowPtrs;
1468 double **CrowPtrs = C.
rowPtrs;
1470 for (
unsigned int i = 0; i < A.
rowNum; ++i) {
1471 for (
unsigned int j = 0; j < A.
colNum; ++j) {
1472 CrowPtrs[i][j] = BrowPtrs[i][j] + ArowPtrs[i][j];
1514 double **ArowPtrs = A.
rowPtrs;
1515 double **BrowPtrs = B.
rowPtrs;
1516 double **CrowPtrs = C.
rowPtrs;
1518 for (
unsigned int i = 0; i < A.
rowNum; ++i) {
1519 for (
unsigned int j = 0; j < A.
colNum; ++j) {
1520 CrowPtrs[i][j] = ArowPtrs[i][j] - BrowPtrs[i][j];
1548 double **ArowPtrs = A.
rowPtrs;
1549 double **BrowPtrs = B.
rowPtrs;
1550 double **CrowPtrs = C.
rowPtrs;
1552 for (
unsigned int i = 0; i < A.
rowNum; ++i) {
1553 for (
unsigned int j = 0; j < A.
colNum; ++j) {
1554 CrowPtrs[i][j] = ArowPtrs[i][j] - BrowPtrs[i][j];
1579 double **BrowPtrs = B.
rowPtrs;
1581 for (
unsigned int i = 0; i <
rowNum; ++i) {
1582 for (
unsigned int j = 0; j <
colNum; ++j) {
1583 rowPtrs[i][j] += BrowPtrs[i][j];
1598 double **BrowPtrs = B.
rowPtrs;
1599 for (
unsigned int i = 0; i <
rowNum; ++i) {
1600 for (
unsigned int j = 0; j <
colNum; ++j) {
1601 rowPtrs[i][j] -= BrowPtrs[i][j];
1623 double **ArowPtrs = A.
rowPtrs;
1624 double **CrowPtrs = C.
rowPtrs;
1626 for (
unsigned int i = 0; i < A.
rowNum; ++i) {
1627 for (
unsigned int j = 0; j < A.
colNum; ++j) {
1628 CrowPtrs[i][j] = -ArowPtrs[i][j];
1647 for (
unsigned int i = 0; i <
rowNum; ++i) {
1648 for (
unsigned int j = 0; j <
colNum; ++j) {
1670 if (std::fabs(x - 1.) < std::numeric_limits<double>::epsilon()) {
1674 unsigned int Brow = B.
getRows();
1675 unsigned int Bcol = B.
getCols();
1678 C.
resize(Brow, Bcol,
false,
false);
1680 for (
unsigned int i = 0; i < Brow; ++i) {
1681 for (
unsigned int j = 0; j < Bcol; ++j) {
1682 C[i][j] = B[i][j] * x;
1695 if (std::fabs(x - 1.) < std::numeric_limits<double>::epsilon()) {
1702 for (
unsigned int i = 0; i <
rowNum; ++i) {
1703 for (
unsigned int j = 0; j <
colNum; ++j) {
1714 if (std::fabs(x - 1.) < std::numeric_limits<double>::epsilon()) {
1718 if (std::fabs(x) < std::numeric_limits<double>::epsilon()) {
1725 double xinv = 1 / x;
1727 for (
unsigned int i = 0; i <
rowNum; ++i) {
1728 for (
unsigned int j = 0; j <
colNum; ++j) {
1729 C[i][j] =
rowPtrs[i][j] * xinv;
1739 for (
unsigned int i = 0; i <
rowNum; ++i) {
1740 for (
unsigned int j = 0; j <
colNum; ++j) {
1751 for (
unsigned int i = 0; i <
rowNum; ++i) {
1752 for (
unsigned int j = 0; j <
colNum; ++j) {
1766 if (std::fabs(x - 1.) < std::numeric_limits<double>::epsilon()) {
1770 for (
unsigned int i = 0; i <
rowNum; ++i) {
1771 for (
unsigned int j = 0; j <
colNum; ++j) {
1782 if (std::fabs(x - 1.) < std::numeric_limits<double>::epsilon()) {
1786 if (std::fabs(x) < std::numeric_limits<double>::epsilon()) {
1790 double xinv = 1 / x;
1792 for (
unsigned int i = 0; i <
rowNum; ++i) {
1793 for (
unsigned int j = 0; j <
colNum; ++j) {
1815 double *optr = out.
data;
1816 for (
unsigned int j = 0; j <
colNum; ++j) {
1817 for (
unsigned int i = 0; i <
rowNum; ++i) {
1872 #if defined(VISP_HAVE_SIMDLIB)
1875 for (
unsigned int i = 0; i <
dsize; ++i) {
1891 unsigned int r1 = m1.
getRows();
1892 unsigned int c1 = m1.
getCols();
1893 unsigned int r2 = m2.
getRows();
1894 unsigned int c2 = m2.
getCols();
1896 out.
resize(r1 * r2, c1 * c2,
false,
false);
1898 for (
unsigned int r = 0; r < r1; ++r) {
1899 for (
unsigned int c = 0; c < c1; ++c) {
1900 double alpha = m1[r][c];
1901 double *m2ptr = m2[0];
1902 unsigned int roffset = r * r2;
1903 unsigned int coffset = c * c2;
1904 for (
unsigned int rr = 0; rr < r2; ++rr) {
1905 for (
unsigned int cc = 0; cc < c2; ++cc) {
1906 out[roffset + rr][coffset + cc] = alpha * *(m2ptr++);
1929 unsigned int r1 = m1.
getRows();
1930 unsigned int c1 = m1.
getCols();
1931 unsigned int r2 = m2.
getRows();
1932 unsigned int c2 = m2.
getCols();
1935 out.
resize(r1 * r2, c1 * c2,
false,
false);
1937 for (
unsigned int r = 0; r < r1; ++r) {
1938 for (
unsigned int c = 0; c < c1; ++c) {
1939 double alpha = m1[r][c];
1940 double *m2ptr = m2[0];
1941 unsigned int roffset = r * r2;
1942 unsigned int coffset = c * c2;
1943 for (
unsigned int rr = 0; rr < r2; ++rr) {
1944 for (
unsigned int cc = 0; cc < c2; ++cc) {
1945 out[roffset + rr][coffset + cc] = alpha * *(m2ptr++);
2134 #if defined(VISP_HAVE_LAPACK)
2136 #elif defined(VISP_HAVE_EIGEN3)
2138 #elif defined(VISP_HAVE_OPENCV)
2203 #if defined(VISP_HAVE_LAPACK)
2205 #elif defined(VISP_HAVE_EIGEN3)
2207 #elif defined(VISP_HAVE_OPENCV)
2213 "Install Lapack, Eigen3 or OpenCV 3rd party"));
2279 #if defined(VISP_HAVE_LAPACK)
2281 #elif defined(VISP_HAVE_EIGEN3)
2283 #elif defined(VISP_HAVE_OPENCV)
2289 "Install Lapack, Eigen3 or OpenCV 3rd party"));
2345 #if defined(VISP_HAVE_LAPACK)
2347 #elif defined(VISP_HAVE_EIGEN3)
2349 #elif defined(VISP_HAVE_OPENCV)
2354 "Install Lapack, Eigen3 or OpenCV 3rd party"));
2410 #if defined(VISP_HAVE_LAPACK)
2412 #elif defined(VISP_HAVE_EIGEN3)
2414 #elif defined(VISP_HAVE_OPENCV)
2419 "Install Lapack, Eigen3 or OpenCV 3rd party"));
2423 #ifndef DOXYGEN_SHOULD_SKIP_THIS
2424 #if defined(VISP_HAVE_LAPACK)
2463 unsigned int nrows =
getRows();
2464 unsigned int ncols =
getCols();
2473 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr,
nullptr,
nullptr,
nullptr);
2520 unsigned int nrows =
getRows();
2521 unsigned int ncols =
getCols();
2530 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr,
nullptr,
nullptr,
nullptr);
2532 return static_cast<unsigned int>(rank_out);
2583 unsigned int nrows =
getRows();
2584 unsigned int ncols =
getCols();
2589 Ap.
resize(ncols, nrows,
false,
false);
2594 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr,
nullptr,
nullptr,
nullptr);
2596 return static_cast<unsigned int>(rank_out);
2708 unsigned int nrows =
getRows();
2709 unsigned int ncols =
getCols();
2713 if (nrows < ncols) {
2714 U.
resize(ncols, ncols,
true);
2723 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr, &imA, &imAt, &kerAt);
2725 return static_cast<unsigned int>(rank_out);
2766 unsigned int nrows =
getRows();
2767 unsigned int ncols =
getCols();
2769 double svThreshold = 1e-26;
2777 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in,
nullptr,
nullptr,
nullptr);
2830 unsigned int nrows =
getRows();
2831 unsigned int ncols =
getCols();
2833 double svThreshold = 1e-26;
2841 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in,
nullptr,
nullptr,
nullptr);
2901 unsigned int nrows =
getRows();
2902 unsigned int ncols =
getCols();
2904 double svThreshold = 1e-26;
2907 Ap.
resize(ncols, nrows,
false,
false);
2912 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in,
nullptr,
nullptr,
nullptr);
3031 unsigned int nrows =
getRows();
3032 unsigned int ncols =
getCols();
3034 double svThreshold = 1e-26;
3037 if (nrows < ncols) {
3038 U.
resize(ncols, ncols,
true);
3047 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, &imA, &imAt, &kerAt);
3053 #if defined(VISP_HAVE_EIGEN3)
3092 unsigned int nrows =
getRows();
3093 unsigned int ncols =
getCols();
3101 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr,
nullptr,
nullptr,
nullptr);
3148 unsigned int nrows =
getRows();
3149 unsigned int ncols =
getCols();
3157 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr,
nullptr,
nullptr,
nullptr);
3159 return static_cast<unsigned int>(rank_out);
3211 unsigned int nrows =
getRows();
3212 unsigned int ncols =
getCols();
3216 Ap.
resize(ncols, nrows,
false,
false);
3221 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr,
nullptr,
nullptr,
nullptr);
3223 return static_cast<unsigned int>(rank_out);
3335 unsigned int nrows =
getRows();
3336 unsigned int ncols =
getCols();
3340 if (nrows < ncols) {
3341 U.
resize(ncols, ncols,
true);
3350 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr, &imA, &imAt, &kerAt);
3352 return static_cast<unsigned int>(rank_out);
3393 unsigned int nrows =
getRows();
3394 unsigned int ncols =
getCols();
3396 double svThreshold = 1e-26;
3403 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in,
nullptr,
nullptr,
nullptr);
3456 unsigned int nrows =
getRows();
3457 unsigned int ncols =
getCols();
3459 double svThreshold = 1e-26;
3463 Ap.
resize(ncols, nrows,
false,
false);
3468 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in,
nullptr,
nullptr,
nullptr);
3528 unsigned int nrows =
getRows();
3529 unsigned int ncols =
getCols();
3531 double svThreshold = 1e-26;
3534 Ap.
resize(ncols, nrows,
false,
false);
3539 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in,
nullptr,
nullptr,
nullptr);
3658 unsigned int nrows =
getRows();
3659 unsigned int ncols =
getCols();
3661 double svThreshold = 1e-26;
3664 if (nrows < ncols) {
3665 U.
resize(ncols, ncols,
true);
3673 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, &imA, &imAt, &kerAt);
3679 #if defined(VISP_HAVE_OPENCV)
3718 unsigned int nrows =
getRows();
3719 unsigned int ncols =
getCols();
3727 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr,
nullptr,
nullptr,
nullptr);
3774 unsigned int nrows =
getRows();
3775 unsigned int ncols =
getCols();
3780 Ap.
resize(ncols, nrows,
false,
false);
3785 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr,
nullptr,
nullptr,
nullptr);
3787 return static_cast<unsigned int>(rank_out);
3838 unsigned int nrows =
getRows();
3839 unsigned int ncols =
getCols();
3843 Ap.
resize(ncols, nrows,
false,
false);
3848 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr,
nullptr,
nullptr,
nullptr);
3850 return static_cast<unsigned int>(rank_out);
3962 unsigned int nrows =
getRows();
3963 unsigned int ncols =
getCols();
3967 if (nrows < ncols) {
3968 U.
resize(ncols, ncols,
true);
3977 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out,
nullptr, &imA, &imAt, &kerAt);
3979 return static_cast<unsigned int>(rank_out);
4020 unsigned int nrows =
getRows();
4021 unsigned int ncols =
getCols();
4023 double svThreshold = 1e-26;
4030 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in,
nullptr,
nullptr,
nullptr);
4083 unsigned int nrows =
getRows();
4084 unsigned int ncols =
getCols();
4086 double svThreshold = 1e-26;
4090 Ap.
resize(ncols, nrows,
false,
false);
4095 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in,
nullptr,
nullptr,
nullptr);
4155 unsigned int nrows =
getRows();
4156 unsigned int ncols =
getCols();
4158 double svThreshold = 1e-26;
4161 Ap.
resize(ncols, nrows,
false,
false);
4166 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in,
nullptr,
nullptr,
nullptr);
4285 unsigned int nrows =
getRows();
4286 unsigned int ncols =
getCols();
4288 double svThreshold = 1e-26;
4291 if (nrows < ncols) {
4292 U.
resize(ncols, ncols,
true);
4301 compute_pseudo_inverse(U, sv, V, nrows, ncols, svThreshold, Ap, rank_out, &rank_in, &imA, &imAt, &kerAt);
4372 #if defined(VISP_HAVE_LAPACK)
4374 #elif defined(VISP_HAVE_EIGEN3)
4376 #elif defined(VISP_HAVE_OPENCV)
4383 "Install Lapack, Eigen3 or OpenCV 3rd party"));
4455 #if defined(VISP_HAVE_LAPACK)
4457 #elif defined(VISP_HAVE_EIGEN3)
4459 #elif defined(VISP_HAVE_OPENCV)
4466 "Install Lapack, Eigen3 or OpenCV 3rd party"));
4546 return pseudoInverse(Ap, sv, svThreshold, imA, imAt, kerAt);
4737 #if defined(VISP_HAVE_LAPACK)
4739 #elif defined(VISP_HAVE_EIGEN3)
4741 #elif defined(VISP_HAVE_OPENCV)
4751 "Install Lapack, Eigen3 or OpenCV 3rd party"));
4864 #if defined(VISP_HAVE_LAPACK)
4866 #elif defined(VISP_HAVE_EIGEN3)
4868 #elif defined(VISP_HAVE_OPENCV)
4878 "Install Lapack, Eigen3 or OpenCV 3rd party"));
4930 for (
unsigned int i = 0; i < column_size; ++i) {
4931 c[i] = (*this)[i_begin + i][j];
5056 if (((j_begin + row_size) >
colNum) || (i >=
rowNum)) {
5060 if ((r.
data !=
nullptr) && (
data !=
nullptr)) {
5061 memcpy(r.
data, (*
this)[i] + j_begin, row_size *
sizeof(
double));
5104 unsigned int min_size = std::min<unsigned int>(
rowNum,
colNum);
5108 diag.resize(min_size,
false);
5110 for (
unsigned int i = 0; i < min_size; ++i) {
5111 diag[i] = (*this)[i][i];
5150 unsigned int nra = A.
getRows();
5151 unsigned int nrb = B.
getRows();
5161 std::cerr <<
"A and C must be two different objects!" << std::endl;
5166 std::cerr <<
"B and C must be two different objects!" << std::endl;
5172 if ((C.
data !=
nullptr) && (A.
data !=
nullptr) && (A.
size() > 0)) {
5177 if ((C.
data !=
nullptr) && (B.
data !=
nullptr) && (B.
size() > 0)) {
5214 std::cerr <<
"A and C must be two different objects!" << std::endl;
5253 std::cerr <<
"A and C must be two different objects!" << std::endl;
5338 unsigned int nca = A.
getCols();
5339 unsigned int ncb = B.
getCols();
5348 if ((B.
getRows() == 0) || ((nca + ncb) == 0)) {
5349 std::cerr <<
"B.getRows() == 0 || nca+ncb == 0" << std::endl;
5384 typedef std::string::size_type size_type;
5389 std::vector<std::string> values(m * n);
5390 std::ostringstream oss;
5391 std::ostringstream ossFixed;
5392 std::ios_base::fmtflags original_flags = oss.flags();
5394 ossFixed.setf(std::ios::fixed, std::ios::floatfield);
5396 size_type maxBefore = 0;
5397 size_type maxAfter = 0;
5399 for (
unsigned int i = 0; i < m; ++i) {
5400 for (
unsigned int j = 0; j < n; ++j) {
5402 oss << (*this)[i][j];
5403 if (oss.str().find(
"e") != std::string::npos) {
5405 ossFixed << (*this)[i][j];
5406 oss.str(ossFixed.str());
5409 values[(i * n) + j] = oss.str();
5410 size_type thislen = values[(i * n) + j].
size();
5411 size_type p = values[(i * n) + j].find(
'.');
5413 if (p == std::string::npos) {
5424 size_type totalLength = length;
5428 maxAfter = std::min<size_type>(maxAfter, totalLength - maxBefore);
5430 if (!intro.empty()) {
5433 s <<
"[" << m <<
"," << n <<
"]=\n";
5435 for (
unsigned int i = 0; i < m; ++i) {
5437 for (
unsigned int j = 0; j < n; ++j) {
5438 size_type p = values[(i * n) + j].find(
'.');
5439 s.setf(std::ios::right, std::ios::adjustfield);
5440 s.width(
static_cast<std::streamsize
>(maxBefore));
5441 s << values[(i * n) + j].substr(0, p).c_str();
5444 s.setf(std::ios::left, std::ios::adjustfield);
5445 if (p != std::string::npos) {
5446 s.width(
static_cast<std::streamsize
>(maxAfter));
5447 s << values[(i * n) + j].substr(p, maxAfter).c_str();
5450 s.width(
static_cast<std::streamsize
>(maxAfter));
5460 s.flags(original_flags);
5462 return static_cast<int>(maxBefore + maxAfter);
5503 unsigned int this_rows = this->
getRows();
5504 unsigned int this_col = this->
getCols();
5506 for (
unsigned int i = 0; i < this_rows; ++i) {
5507 for (
unsigned int j = 0; j < this_col; ++j) {
5508 os << (*this)[i][j] <<
", ";
5510 if (this->
getRows() != (i + 1)) {
5511 os <<
";" << std::endl;
5514 os <<
"]" << std::endl;
5550 unsigned int this_rows = this->
getRows();
5551 unsigned int this_col = this->
getCols();
5552 os <<
"([ " << std::endl;
5553 for (
unsigned int i = 0; i < this_rows; ++i) {
5555 for (
unsigned int j = 0; j < this_col; ++j) {
5556 os << (*this)[i][j] <<
", ";
5558 os <<
"]," << std::endl;
5560 os <<
"])" << std::endl;
5593 unsigned int this_rows = this->
getRows();
5594 unsigned int this_col = this->
getCols();
5595 for (
unsigned int i = 0; i < this_rows; ++i) {
5596 for (
unsigned int j = 0; j < this_col; ++j) {
5597 os << (*this)[i][j];
5598 if (!(j == (this->
getCols() - 1))) {
5644 os <<
"vpMatrix " << matrixName <<
" (" << this->
getRows() <<
", " << this->
getCols() <<
"); " << std::endl;
5646 unsigned int this_rows = this->
getRows();
5647 unsigned int this_col = this->
getCols();
5648 for (
unsigned int i = 0; i < this_rows; ++i) {
5649 for (
unsigned int j = 0; j < this_col; ++j) {
5651 os << matrixName <<
"[" << i <<
"][" << j <<
"] = " << (*this)[i][j] <<
"; " << std::endl;
5654 for (
unsigned int k = 0; k <
sizeof(double); ++k) {
5655 os <<
"((unsigned char*)&(" << matrixName <<
"[" << i <<
"][" << j <<
"]) )[" << k <<
"] = 0x" << std::hex
5656 <<
static_cast<unsigned int>(((
unsigned char *)&((*
this)[i][j]))[k]) <<
"; " << std::endl;
5680 unsigned int rowNumOld =
rowNum;
5712 if (r.
size() == 0) {
5716 unsigned int oldSize =
size();
5721 memcpy(
data + oldSize, r.
data,
sizeof(
double) * r.
size());
5753 if (c.
size() == 0) {
5758 unsigned int oldColNum =
colNum;
5763 for (
unsigned int i = 0; i <
rowNum; ++i) {
5764 memcpy(
data + (i *
colNum), tmp.
data + (i * oldColNum),
sizeof(
double) * oldColNum);
5788 unsigned int a_rows = A.
getRows();
5789 for (
unsigned int i = r; i < (r + a_rows); ++i) {
5847 for (
unsigned int i = 0; i <
rowNum; ++i) {
5848 for (
unsigned int j = 0; j <
rowNum; ++j) {
5850 if (std::fabs(At_A[i][j]) > std::numeric_limits<double>::epsilon()) {
5856 #if defined(VISP_HAVE_LAPACK)
5857 #if defined(VISP_HAVE_GSL)
5859 gsl_vector *eval = gsl_vector_alloc(
rowNum);
5862 gsl_eigen_symmv_workspace *w = gsl_eigen_symmv_alloc(
rowNum);
5865 unsigned int Atda = (
unsigned int)m->tda;
5866 for (
unsigned int i = 0; i <
rowNum; i++) {
5867 unsigned int k = i * Atda;
5868 for (
unsigned int j = 0; j <
colNum; j++)
5869 m->data[k + j] = (*
this)[i][j];
5871 gsl_eigen_symmv(m, eval, evec, w);
5873 gsl_eigen_symmv_sort(eval, evec, GSL_EIGEN_SORT_ABS_ASC);
5875 for (
unsigned int i = 0; i <
rowNum; i++) {
5876 evalue[i] = gsl_vector_get(eval, i);
5879 gsl_eigen_symmv_free(w);
5880 gsl_vector_free(eval);
5882 gsl_matrix_free(evec);
5886 const char jobz =
'N';
5887 const char uplo =
'U';
5894 lwork =
static_cast<int>(wkopt);
5902 "You should install Lapack 3rd party"));
5966 for (
unsigned int i = 0; i <
rowNum; ++i) {
5967 for (
unsigned int j = 0; j <
rowNum; ++j) {
5969 if (std::fabs(At_A[i][j]) > std::numeric_limits<double>::epsilon()) {
5979 #if defined(VISP_HAVE_LAPACK)
5980 #if defined(VISP_HAVE_GSL)
5982 gsl_vector *eval = gsl_vector_alloc(
rowNum);
5985 gsl_eigen_symmv_workspace *w = gsl_eigen_symmv_alloc(
rowNum);
5988 unsigned int Atda = (
unsigned int)m->tda;
5989 for (
unsigned int i = 0; i <
rowNum; i++) {
5990 unsigned int k = i * Atda;
5991 for (
unsigned int j = 0; j <
colNum; j++)
5992 m->data[k + j] = (*
this)[i][j];
5994 gsl_eigen_symmv(m, eval, evec, w);
5996 gsl_eigen_symmv_sort(eval, evec, GSL_EIGEN_SORT_ABS_ASC);
5998 for (
unsigned int i = 0; i <
rowNum; i++) {
5999 evalue[i] = gsl_vector_get(eval, i);
6001 Atda = (
unsigned int)evec->tda;
6002 for (
unsigned int i = 0; i <
rowNum; i++) {
6003 unsigned int k = i * Atda;
6004 for (
unsigned int j = 0; j <
rowNum; j++) {
6005 evector[i][j] = evec->
data[k + j];
6009 gsl_eigen_symmv_free(w);
6010 gsl_vector_free(eval);
6012 gsl_matrix_free(evec);
6016 const char jobz =
'V';
6017 const char uplo =
'U';
6024 lwork =
static_cast<int>(wkopt);
6033 "You should install Lapack 3rd party"));
6058 unsigned int nbline =
getRows();
6059 unsigned int nbcol =
getCols();
6064 V.
resize(nbcol, nbcol,
false);
6070 if (nbline < nbcol) {
6071 U.
resize(nbcol, nbcol,
true);
6074 U.
resize(nbline, nbcol,
false);
6083 for (
unsigned int i = 0; i < nbcol; ++i) {
6084 if (sv[i] > maxsv) {
6089 unsigned int rank = 0;
6090 for (
unsigned int i = 0; i < nbcol; ++i) {
6091 if (sv[i] > (maxsv * svThreshold)) {
6096 kerAt.
resize(nbcol - rank, nbcol);
6097 if (rank != nbcol) {
6098 for (
unsigned int j = 0, k = 0; j < nbcol; ++j) {
6100 if ((sv[j] <= (maxsv * svThreshold)) &&
6101 (std::fabs(V.
getCol(j).
sumSquare()) > std::numeric_limits<double>::epsilon())) {
6102 unsigned int v_rows = V.
getRows();
6103 for (
unsigned int i = 0; i < v_rows; ++i) {
6104 kerAt[k][i] = V[i][j];
6132 unsigned int nbrow =
getRows();
6133 unsigned int nbcol =
getCols();
6138 V.
resize(nbcol, nbcol,
false);
6144 if (nbrow < nbcol) {
6145 U.
resize(nbcol, nbcol,
true);
6148 U.
resize(nbrow, nbcol,
false);
6156 double maxsv = sv[0];
6158 unsigned int rank = 0;
6159 for (
unsigned int i = 0; i < nbcol; ++i) {
6160 if (sv[i] > (maxsv * svThreshold)) {
6165 kerA.
resize(nbcol, nbcol - rank);
6166 if (rank != nbcol) {
6167 for (
unsigned int j = 0, k = 0; j < nbcol; ++j) {
6169 if (sv[j] <= (maxsv * svThreshold)) {
6170 for (
unsigned int i = 0; i < nbcol; ++i) {
6171 kerA[i][k] = V[i][j];
6178 return (nbcol - rank);
6199 unsigned int nbrow =
getRows();
6200 unsigned int nbcol =
getCols();
6201 unsigned int dim_ =
static_cast<unsigned int>(dim);
6206 V.
resize(nbcol, nbcol,
false);
6212 if (nbrow < nbcol) {
6213 U.
resize(nbcol, nbcol,
true);
6216 U.
resize(nbrow, nbcol,
false);
6223 kerA.
resize(nbcol, dim_);
6225 unsigned int rank = nbcol - dim_;
6226 for (
unsigned int k = 0; k < dim_; ++k) {
6227 unsigned int j = k + rank;
6228 for (
unsigned int i = 0; i < nbcol; ++i) {
6229 kerA[i][k] = V[i][j];
6234 double maxsv = sv[0];
6235 unsigned int rank = 0;
6236 for (
unsigned int i = 0; i < nbcol; ++i) {
6237 if (sv[i] > (maxsv * 1e-6)) {
6241 return (nbcol - rank);
6300 #ifdef VISP_HAVE_GSL
6302 double *b =
new double[size_];
6303 for (
size_t i = 0; i < size_; i++)
6305 gsl_matrix_view m = gsl_matrix_view_array(this->
data,
rowNum, colNum);
6306 gsl_matrix_view em = gsl_matrix_view_array(b,
rowNum,
colNum);
6307 gsl_linalg_exponential_ss(&m.matrix, &em.matrix, 0);
6311 memcpy(expA.
data, b, size_ *
sizeof(
double));
6332 for (
unsigned int i = 0; i <
rowNum; ++i) {
6334 for (
unsigned int j = 0; j <
colNum; ++j) {
6335 sum += fabs((*
this)[i][j]);
6337 if ((
sum > nA) || (i == 0)) {
6346 double sca = 1.0 / pow(2.0, s);
6349 v_expE = (c * exp) + v_eye;
6350 v_expD = (-c * exp) + v_eye;
6351 for (
int k = 2; k <= q; ++k) {
6352 c = (c * (
static_cast<double>((q - k) + 1))) / (
static_cast<double>(k * (((2 * q) - k) + 1)));
6353 v_expcX = exp * v_expX;
6355 v_expcX = c * v_expX;
6356 v_expE = v_expE + v_expcX;
6358 v_expD = v_expD + v_expcX;
6361 v_expD = v_expD - v_expcX;
6366 exp = v_expX * v_expE;
6367 for (
int k = 1; k <= s; ++k) {
6400 unsigned int m_rows = M.
getRows();
6401 unsigned int m_col = M.
getCols();
6402 for (
unsigned int i = 0; i < col; ++i) {
6403 for (
unsigned int j = 0; j < row; ++j) {
6404 M_comp[i][j] = M[i][j];
6406 for (
unsigned int j = row + 1; j < m_rows; ++j) {
6407 M_comp[i][j - 1] = M[i][j];
6410 for (
unsigned int i = col + 1; i < m_col; ++i) {
6411 for (
unsigned int j = 0; j < row; ++j) {
6412 M_comp[i - 1][j] = M[i][j];
6414 for (
unsigned int j = row + 1; j < m_rows; ++j) {
6415 M_comp[i - 1][j - 1] = M[i][j];
6432 unsigned int nbline =
getRows();
6433 unsigned int nbcol =
getCols();
6438 V.
resize(nbcol, nbcol,
false);
6444 if (nbline < nbcol) {
6445 U.
resize(nbcol, nbcol,
true);
6448 U.
resize(nbline, nbcol,
false);
6457 for (
unsigned int i = 0; i < nbcol; ++i) {
6458 if (sv[i] > maxsv) {
6464 unsigned int rank = 0;
6465 for (
unsigned int i = 0; i < nbcol; ++i) {
6466 if (sv[i] > (maxsv * svThreshold)) {
6472 double minsv = maxsv;
6473 for (
unsigned int i = 0; i < rank; ++i) {
6474 if (sv[i] < minsv) {
6479 if (std::fabs(minsv) > std::numeric_limits<double>::epsilon()) {
6480 return maxsv / minsv;
6483 return std::numeric_limits<double>::infinity();
6501 unsigned int h_col = H.
getCols();
6502 for (
unsigned int i = 0; i < h_col; ++i) {
6503 HLM[i][i] += alpha * H[i][i];
6517 for (
unsigned int i = 0; i <
dsize; ++i) {
6518 double x = *(
data + i);
6535 if (this->
dsize != 0) {
6544 unsigned int maxRank = std::min<unsigned int>(this->
getCols(), this->
getRows());
6547 unsigned int boundary = std::min<unsigned int>(maxRank, w.
size());
6552 for (
unsigned int i = 0; i < boundary; ++i) {
6577 for (
unsigned int i = 0; i <
rowNum; ++i) {
6579 for (
unsigned int j = 0; j <
colNum; ++j) {
6580 x += fabs(*(*(
rowPtrs + i) + j));
6597 double sum_square = 0.0;
6600 for (
unsigned int i = 0; i <
rowNum; ++i) {
6601 for (
unsigned int j = 0; j <
colNum; ++j) {
6603 sum_square += x * x;
6609 #if defined(VISP_BUILD_DEPRECATED_FUNCTIONS)
6657 for (
unsigned int j = 0; j <
getCols(); ++j) {
6658 c[j] = (*this)[i - 1][j];
6684 for (
unsigned int i = 0; i <
getRows(); ++i) {
6685 c[i] = (*this)[i][j - 1];
6698 for (
unsigned int i = 0; i <
rowNum; ++i) {
6699 for (
unsigned int j = 0; j <
colNum; ++j) {
6701 (*this)[i][j] = val;
Implementation of a generic 2D array used as base class for matrices and vectors.
unsigned int getCols() const
double * data
Address of the first element of the data array.
double ** rowPtrs
Address of the first element of each rows.
void insert(const vpArray2D< Type > &A, unsigned int r, unsigned int c)
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 dsize
Current array size (rowNum * colNum)
unsigned int size() const
Return the number of elements of the 2D array.
unsigned int getRows() const
unsigned int colNum
Number of columns in the array.
Implementation of column vector and the associated operations.
void resize(unsigned int i, bool flagNullify=true)
error that can be emitted by ViSP classes.
@ functionNotImplementedError
Function not implemented.
@ dimensionError
Bad dimension.
@ divideByZeroError
Division by zero.
Implementation of an homogeneous matrix and operations on such kind of matrices.
static Type maximum(const Type &a, const Type &b)
Implementation of a matrix and operations on matrices.
void svdLapack(vpColVector &w, vpMatrix &V)
vpColVector eigenValues() const
vpMatrix pseudoInverseOpenCV(double svThreshold=1e-6) const
std::ostream & cppPrint(std::ostream &os, const std::string &matrixName="A", bool octet=false) const
vpMatrix & operator<<(double *)
double cond(double svThreshold=1e-6) const
vpMatrix pseudoInverseEigen3(double svThreshold=1e-6) const
vpMatrix hadamard(const vpMatrix &m) const
vpMatrix operator-() const
vp_deprecated void setIdentity(const double &val=1.0)
vpMatrix & operator/=(double x)
Divide all the element of the matrix by x : Aij = Aij / x.
int print(std::ostream &s, unsigned int length, const std::string &intro="") const
std::ostream & maplePrint(std::ostream &os) const
unsigned int kernel(vpMatrix &kerAt, double svThreshold=1e-6) const
vpMatrix inverseByLU() const
vpMatrix operator*(const vpMatrix &B) const
vpMatrix operator/(double x) const
Cij = Aij / x (A is unchanged)
void stack(const vpMatrix &A)
vp_deprecated void stackMatrices(const vpMatrix &A)
vpMatrix & operator-=(const vpMatrix &B)
Operation A = A - B.
vpMatrix & operator*=(double x)
Multiply all the element of the matrix by x : Aij = Aij * x.
vpMatrix operator*(const double &x, const vpMatrix &B)
static vpMatrix juxtaposeMatrices(const vpMatrix &A, const vpMatrix &B)
vp_deprecated double euclideanNorm() const
vpColVector stackColumns()
vp_deprecated vpColVector column(unsigned int j)
static void sub2Matrices(const vpMatrix &A, const vpMatrix &B, vpMatrix &C)
void svd(vpColVector &w, vpMatrix &V)
vpMatrix operator+(const vpMatrix &B) const
void diag(const double &val=1.0)
vpRowVector getRow(unsigned int i) const
void svdOpenCV(vpColVector &w, vpMatrix &V)
vpMatrix & operator+=(const vpMatrix &B)
Operation A = A + B.
static void add2WeightedMatrices(const vpMatrix &A, const double &wA, const vpMatrix &B, const double &wB, vpMatrix &C)
double inducedL2Norm() const
double infinityNorm() const
vpMatrix & operator=(const vpArray2D< double > &A)
double frobeniusNorm() const
vpColVector getDiag() const
static void createDiagonalMatrix(const vpColVector &A, vpMatrix &DA)
void solveBySVD(const vpColVector &B, vpColVector &x) const
static void computeHLM(const vpMatrix &H, const double &alpha, vpMatrix &HLM)
vpMatrix pseudoInverseLapack(double svThreshold=1e-6) const
vpMatrix pseudoInverse(double svThreshold=1e-6) const
vpMatrix & operator,(double val)
static void multMatrixVector(const vpMatrix &A, const vpColVector &v, vpColVector &w)
vp_deprecated vpRowVector row(unsigned int i)
vpColVector getCol(unsigned int j) const
double det(vpDetMethod method=LU_DECOMPOSITION) const
vp_deprecated void init()
static void add2Matrices(const vpMatrix &A, const vpMatrix &B, vpMatrix &C)
void insert(const vpMatrix &A, unsigned int r, unsigned int c)
void svdEigen3(vpColVector &w, vpMatrix &V)
void kron(const vpMatrix &m1, vpMatrix &out) const
vpMatrix transpose() const
static void mult2Matrices(const vpMatrix &A, const vpMatrix &B, vpMatrix &C)
static void negateMatrix(const vpMatrix &A, vpMatrix &C)
std::ostream & csvPrint(std::ostream &os) const
unsigned int nullSpace(vpMatrix &kerA, double svThreshold=1e-6) const
std::ostream & matlabPrint(std::ostream &os) const
vpMatrix extract(unsigned int r, unsigned int c, unsigned int nrows, unsigned int ncols) const
Implementation of a rotation matrix and operations on such kind of matrices.
Implementation of row vector and the associated operations.
void resize(unsigned int i, bool flagNullify=true)
Class that consider the case of a translation vector.