aboutsummaryrefslogtreecommitdiffhomepage
diff options
context:
space:
mode:
Diffstat
-rw-r--r--Eigen/Core2
-rw-r--r--Eigen/src/Cholesky/Cholesky.h2
-rw-r--r--Eigen/src/Cholesky/CholeskyWithoutSquareRoot.h4
-rw-r--r--Eigen/src/Core/MatrixBase.h4
-rwxr-xr-xEigen/src/Core/SolveTriangular.h (renamed from Eigen/src/Core/InverseProduct.h)10
-rw-r--r--Eigen/src/LU/LU.h6
-rw-r--r--Eigen/src/QR/SelfAdjointEigenSolver.h6
-rw-r--r--Eigen/src/Sparse/SparseMatrixBase.h2
-rw-r--r--Eigen/src/Sparse/TriangularSolver.h2
-rw-r--r--bench/benchCholesky.cpp14
-rw-r--r--bench/benchVecAdd.cpp86
-rw-r--r--bench/btl/libs/blitz/btl_tiny_blitz.cpp1
-rw-r--r--bench/btl/libs/eigen2/eigen2_interface.hh3
-rw-r--r--doc/snippets/MatrixBase_marked.cpp2
-rw-r--r--test/submatrices.cpp10
-rw-r--r--test/triangular.cpp8
16 files changed, 90 insertions, 72 deletions
diff --git a/Eigen/Core b/Eigen/Core
index db2fde30d..0cb89bb60 100644
--- a/Eigen/Core
+++ b/Eigen/Core
@@ -44,7 +44,7 @@ namespace Eigen {
#include "src/Core/Dot.h"
#include "src/Core/Product.h"
#include "src/Core/DiagonalProduct.h"
-#include "src/Core/InverseProduct.h"
+#include "src/Core/SolveTriangular.h"
#include "src/Core/MapBase.h"
#include "src/Core/Map.h"
#include "src/Core/Block.h"
diff --git a/Eigen/src/Cholesky/Cholesky.h b/Eigen/src/Cholesky/Cholesky.h
index 48145fa46..af5dfb430 100644
--- a/Eigen/src/Cholesky/Cholesky.h
+++ b/Eigen/src/Cholesky/Cholesky.h
@@ -141,7 +141,7 @@ typename Derived::Eval Cholesky<MatrixType>::solve(const MatrixBase<Derived> &b)
const int size = m_matrix.rows();
ei_assert(size==b.rows());
- return m_matrix.adjoint().template part<Upper>().inverseProduct(matrixL().inverseProduct(b));
+ return m_matrix.adjoint().template part<Upper>().solveTriangular(matrixL().solveTriangular(b));
}
/** \cholesky_module
diff --git a/Eigen/src/Cholesky/CholeskyWithoutSquareRoot.h b/Eigen/src/Cholesky/CholeskyWithoutSquareRoot.h
index 1eb9904fc..69dcdb8a2 100644
--- a/Eigen/src/Cholesky/CholeskyWithoutSquareRoot.h
+++ b/Eigen/src/Cholesky/CholeskyWithoutSquareRoot.h
@@ -148,9 +148,9 @@ typename Derived::Eval CholeskyWithoutSquareRoot<MatrixType>::solve(const Matrix
ei_assert(size==b.rows());
return m_matrix.adjoint().template part<UnitUpper>()
- .inverseProduct(
+ .solveTriangular(
(matrixL()
- .inverseProduct(b))
+ .solveTriangular(b))
.cwise()/m_matrix.diagonal()
);
}
diff --git a/Eigen/src/Core/MatrixBase.h b/Eigen/src/Core/MatrixBase.h
index 25a4c8b08..8ad54e62a 100644
--- a/Eigen/src/Core/MatrixBase.h
+++ b/Eigen/src/Core/MatrixBase.h
@@ -320,10 +320,10 @@ template<typename Derived> class MatrixBase
Derived& operator*=(const MatrixBase<OtherDerived>& other);
template<typename OtherDerived>
- typename OtherDerived::Eval inverseProduct(const MatrixBase<OtherDerived>& other) const;
+ typename OtherDerived::Eval solveTriangular(const MatrixBase<OtherDerived>& other) const;
template<typename OtherDerived>
- void inverseProductInPlace(MatrixBase<OtherDerived>& other) const;
+ void solveTriangularInPlace(MatrixBase<OtherDerived>& other) const;
template<typename OtherDerived>
diff --git a/Eigen/src/Core/InverseProduct.h b/Eigen/src/Core/SolveTriangular.h
index cd27c4ddf..a9867e63c 100755
--- a/Eigen/src/Core/InverseProduct.h
+++ b/Eigen/src/Core/SolveTriangular.h
@@ -212,13 +212,13 @@ struct ei_trisolve_selector<Lhs,Rhs,Upper,ColMajor>
}
};
-/** "in-place" version of MatrixBase::inverseProduct() where the result is written in \a other
+/** "in-place" version of MatrixBase::solveTriangular() where the result is written in \a other
*
- * \sa inverseProduct()
+ * \sa solveTriangular()
*/
template<typename Derived>
template<typename OtherDerived>
-void MatrixBase<Derived>::inverseProductInPlace(MatrixBase<OtherDerived>& other) const
+void MatrixBase<Derived>::solveTriangularInPlace(MatrixBase<OtherDerived>& other) const
{
ei_assert(derived().cols() == derived().rows());
ei_assert(derived().cols() == other.rows());
@@ -244,10 +244,10 @@ void MatrixBase<Derived>::inverseProductInPlace(MatrixBase<OtherDerived>& other)
*/
template<typename Derived>
template<typename OtherDerived>
-typename OtherDerived::Eval MatrixBase<Derived>::inverseProduct(const MatrixBase<OtherDerived>& other) const
+typename OtherDerived::Eval MatrixBase<Derived>::solveTriangular(const MatrixBase<OtherDerived>& other) const
{
typename OtherDerived::Eval res(other);
- inverseProductInPlace(res);
+ solveTriangularInPlace(res);
return res;
}
diff --git a/Eigen/src/LU/LU.h b/Eigen/src/LU/LU.h
index a7d5cd2e8..59f1051be 100644
--- a/Eigen/src/LU/LU.h
+++ b/Eigen/src/LU/LU.h
@@ -257,7 +257,7 @@ void LU<MatrixType>::computeKernel(Matrix<typename MatrixType::Scalar,
m_lu.corner(TopLeft, m_rank, m_rank)
.template marked<Upper>()
- .inverseProductInPlace(y);
+ .solveTriangularInPlace(y);
for(int i = 0; i < m_rank; i++)
result->row(m_q.coeff(i)) = y.row(i);
@@ -309,7 +309,7 @@ bool LU<MatrixType>::solve(
l.setZero();
l.corner(Eigen::TopLeft,rows,smalldim)
= m_lu.corner(Eigen::TopLeft,rows,smalldim);
- l.template marked<UnitLower>().inverseProductInPlace(c);
+ l.template marked<UnitLower>().solveTriangularInPlace(c);
// Step 3
if(!isSurjective())
@@ -326,7 +326,7 @@ bool LU<MatrixType>::solve(
d(c.corner(TopLeft, m_rank, c.cols()));
m_lu.corner(TopLeft, m_rank, m_rank)
.template marked<Upper>()
- .inverseProductInPlace(d);
+ .solveTriangularInPlace(d);
// Step 4
result->resize(m_lu.cols(), b.cols());
diff --git a/Eigen/src/QR/SelfAdjointEigenSolver.h b/Eigen/src/QR/SelfAdjointEigenSolver.h
index 2d833bf13..d366bafc7 100644
--- a/Eigen/src/QR/SelfAdjointEigenSolver.h
+++ b/Eigen/src/QR/SelfAdjointEigenSolver.h
@@ -230,17 +230,17 @@ compute(const MatrixType& matA, const MatrixType& matB, bool computeEigenvectors
Cholesky<MatrixType> cholB(matB);
// compute C = inv(U') A inv(U)
- MatrixType matC = cholB.matrixL().inverseProduct(matA);
+ MatrixType matC = cholB.matrixL().solveTriangular(matA);
// FIXME since we currently do not support A * inv(U),
// let's do (inv(U') A')' :
- matC = (cholB.matrixL().inverseProduct(matC.adjoint())).adjoint();
+ matC = (cholB.matrixL().solveTriangular(matC.adjoint())).adjoint();
compute(matC, computeEigenvectors);
if (computeEigenvectors)
{
// transform back the eigen vectors: evecs = inv(U) * evecs
- m_eivec = cholB.matrixL().adjoint().template marked<Upper>().inverseProduct(m_eivec);
+ m_eivec = cholB.matrixL().adjoint().template marked<Upper>().solveTriangular(m_eivec);
}
}
diff --git a/Eigen/src/Sparse/SparseMatrixBase.h b/Eigen/src/Sparse/SparseMatrixBase.h
index cc43e7e07..e36e7e760 100644
--- a/Eigen/src/Sparse/SparseMatrixBase.h
+++ b/Eigen/src/Sparse/SparseMatrixBase.h
@@ -142,7 +142,7 @@ class SparseMatrixBase : public MatrixBase<Derived>
}
template<typename OtherDerived>
- OtherDerived inverseProduct(const MatrixBase<OtherDerived>& other) const;
+ OtherDerived solveTriangular(const MatrixBase<OtherDerived>& other) const;
protected:
diff --git a/Eigen/src/Sparse/TriangularSolver.h b/Eigen/src/Sparse/TriangularSolver.h
index 41361a471..45679fe86 100644
--- a/Eigen/src/Sparse/TriangularSolver.h
+++ b/Eigen/src/Sparse/TriangularSolver.h
@@ -155,7 +155,7 @@ struct ei_sparse_trisolve_selector<Lhs,Rhs,Upper,ColMajor>
template<typename Derived>
template<typename OtherDerived>
-OtherDerived SparseMatrixBase<Derived>::inverseProduct(const MatrixBase<OtherDerived>& other) const
+OtherDerived SparseMatrixBase<Derived>::solveTriangular(const MatrixBase<OtherDerived>& other) const
{
ei_assert(derived().cols() == other.rows());
ei_assert(!(Flags & ZeroDiagBit));
diff --git a/bench/benchCholesky.cpp b/bench/benchCholesky.cpp
index dd0c7f83b..f64b61b71 100644
--- a/bench/benchCholesky.cpp
+++ b/bench/benchCholesky.cpp
@@ -18,7 +18,7 @@ using namespace Eigen;
#endif
#ifndef TRIES
-#define TRIES 4
+#define TRIES 10
#endif
typedef float Scalar;
@@ -29,6 +29,13 @@ __attribute__ ((noinline)) void benchCholesky(const MatrixType& m)
int rows = m.rows();
int cols = m.cols();
+ int cost = 0;
+ for (int j=0; j<rows; ++j)
+ {
+ int r = std::max(rows - j -1,0);
+ cost += 2*(r*j+r+j);
+ }
+
int repeats = (REPEAT*1000)/(rows*rows);
typedef typename MatrixType::Scalar Scalar;
@@ -70,7 +77,8 @@ __attribute__ ((noinline)) void benchCholesky(const MatrixType& m)
std::cout << "fixed ";
std::cout << covMat.rows() << " \t"
<< (timerNoSqrt.value() * REPEAT) / repeats << "s \t"
- << (timerSqrt.value() * REPEAT) / repeats << "s";
+ << (timerSqrt.value() * REPEAT) / repeats << "s "
+ << "(" << 1e-6 * cost*repeats/timerSqrt.value() << " MFLOPS)\n";
#ifdef BENCH_GSL
@@ -108,7 +116,7 @@ __attribute__ ((noinline)) void benchCholesky(const MatrixType& m)
int main(int argc, char* argv[])
{
- const int dynsizes[] = {/*4,6,8,12,16,24,32,49,64,67,128,129,130,131,132,*/256,257,258,259,260,512,0};
+ const int dynsizes[] = {/*4,6,8,12,16,24,32,49,64,67,128,129,130,131,132,*/256,257,258,259,260,512,900,0};
std::cout << "size no sqrt standard";
#ifdef BENCH_GSL
std::cout << " GSL (standard + double + ATLAS) ";
diff --git a/bench/benchVecAdd.cpp b/bench/benchVecAdd.cpp
index aa211dce0..ef85b3ced 100644
--- a/bench/benchVecAdd.cpp
+++ b/bench/benchVecAdd.cpp
@@ -22,7 +22,7 @@ int main(int argc, char* argv[])
int size = SIZE * 8;
int size2 = size * size;
Scalar* a = ei_aligned_malloc<Scalar>(size2);
- Scalar* b = ei_aligned_malloc<Scalar>(size2);
+ Scalar* b = ei_aligned_malloc<Scalar>(size2+4)+1;
Scalar* c = ei_aligned_malloc<Scalar>(size2);
for (int i=0; i<size; ++i)
@@ -33,22 +33,22 @@ int main(int argc, char* argv[])
BenchTimer timer;
timer.reset();
- for (int k=0; k<3; ++k)
+ for (int k=0; k<10; ++k)
{
timer.start();
benchVec(a, b, c, size2);
timer.stop();
}
std::cout << timer.value() << "s " << (double(size2*REPEAT)/timer.value())/(1024.*1024.*1024.) << " GFlops\n";
-
+ return 0;
for (int innersize = size; innersize>2 ; --innersize)
{
if (size2%innersize==0)
{
int outersize = size2/innersize;
- MatrixXf ma = MatrixXf::map(a, innersize, outersize );
- MatrixXf mb = MatrixXf::map(b, innersize, outersize );
- MatrixXf mc = MatrixXf::map(c, innersize, outersize );
+ MatrixXf ma = Map<MatrixXf>(a, innersize, outersize );
+ MatrixXf mb = Map<MatrixXf>(b, innersize, outersize );
+ MatrixXf mc = Map<MatrixXf>(c, innersize, outersize );
timer.reset();
for (int k=0; k<3; ++k)
{
@@ -60,9 +60,9 @@ int main(int argc, char* argv[])
}
}
- VectorXf va = VectorXf::map(a, size2);
- VectorXf vb = VectorXf::map(b, size2);
- VectorXf vc = VectorXf::map(c, size2);
+ VectorXf va = Map<VectorXf>(a, size2);
+ VectorXf vb = Map<VectorXf>(b, size2);
+ VectorXf vc = Map<VectorXf>(c, size2);
timer.reset();
for (int k=0; k<3; ++k)
{
@@ -95,40 +95,40 @@ void benchVec(Scalar* a, Scalar* b, Scalar* c, int size)
for (int k=0; k<REPEAT; ++k)
for (int i=0; i<size; i+=PacketSize*8)
{
- a0 = ei_pload(&a[i]);
- b0 = ei_pload(&b[i]);
- a1 = ei_pload(&a[i+1*PacketSize]);
- b1 = ei_pload(&b[i+1*PacketSize]);
- a2 = ei_pload(&a[i+2*PacketSize]);
- b2 = ei_pload(&b[i+2*PacketSize]);
- a3 = ei_pload(&a[i+3*PacketSize]);
- b3 = ei_pload(&b[i+3*PacketSize]);
- ei_pstore(&a[i], ei_padd(a0, b0));
- a0 = ei_pload(&a[i+4*PacketSize]);
- b0 = ei_pload(&b[i+4*PacketSize]);
-
- ei_pstore(&a[i+1*PacketSize], ei_padd(a1, b1));
- a1 = ei_pload(&a[i+5*PacketSize]);
- b1 = ei_pload(&b[i+5*PacketSize]);
-
- ei_pstore(&a[i+2*PacketSize], ei_padd(a2, b2));
- a2 = ei_pload(&a[i+6*PacketSize]);
- b2 = ei_pload(&b[i+6*PacketSize]);
-
- ei_pstore(&a[i+3*PacketSize], ei_padd(a3, b3));
- a3 = ei_pload(&a[i+7*PacketSize]);
- b3 = ei_pload(&b[i+7*PacketSize]);
-
- ei_pstore(&a[i+4*PacketSize], ei_padd(a0, b0));
- ei_pstore(&a[i+5*PacketSize], ei_padd(a1, b1));
- ei_pstore(&a[i+6*PacketSize], ei_padd(a2, b2));
- ei_pstore(&a[i+7*PacketSize], ei_padd(a3, b3));
+// a0 = ei_pload(&a[i]);
+// b0 = ei_pload(&b[i]);
+// a1 = ei_pload(&a[i+1*PacketSize]);
+// b1 = ei_pload(&b[i+1*PacketSize]);
+// a2 = ei_pload(&a[i+2*PacketSize]);
+// b2 = ei_pload(&b[i+2*PacketSize]);
+// a3 = ei_pload(&a[i+3*PacketSize]);
+// b3 = ei_pload(&b[i+3*PacketSize]);
+// ei_pstore(&a[i], ei_padd(a0, b0));
+// a0 = ei_pload(&a[i+4*PacketSize]);
+// b0 = ei_pload(&b[i+4*PacketSize]);
+//
+// ei_pstore(&a[i+1*PacketSize], ei_padd(a1, b1));
+// a1 = ei_pload(&a[i+5*PacketSize]);
+// b1 = ei_pload(&b[i+5*PacketSize]);
+//
+// ei_pstore(&a[i+2*PacketSize], ei_padd(a2, b2));
+// a2 = ei_pload(&a[i+6*PacketSize]);
+// b2 = ei_pload(&b[i+6*PacketSize]);
+//
+// ei_pstore(&a[i+3*PacketSize], ei_padd(a3, b3));
+// a3 = ei_pload(&a[i+7*PacketSize]);
+// b3 = ei_pload(&b[i+7*PacketSize]);
+//
+// ei_pstore(&a[i+4*PacketSize], ei_padd(a0, b0));
+// ei_pstore(&a[i+5*PacketSize], ei_padd(a1, b1));
+// ei_pstore(&a[i+6*PacketSize], ei_padd(a2, b2));
+// ei_pstore(&a[i+7*PacketSize], ei_padd(a3, b3));
-// ei_pstore(&a[i+2*PacketSize], ei_padd(ei_pload(&a[i+2*PacketSize]), ei_pload(&b[i+2*PacketSize])));
-// ei_pstore(&a[i+3*PacketSize], ei_padd(ei_pload(&a[i+3*PacketSize]), ei_pload(&b[i+3*PacketSize])));
-// ei_pstore(&a[i+4*PacketSize], ei_padd(ei_pload(&a[i+4*PacketSize]), ei_pload(&b[i+4*PacketSize])));
-// ei_pstore(&a[i+5*PacketSize], ei_padd(ei_pload(&a[i+5*PacketSize]), ei_pload(&b[i+5*PacketSize])));
-// ei_pstore(&a[i+6*PacketSize], ei_padd(ei_pload(&a[i+6*PacketSize]), ei_pload(&b[i+6*PacketSize])));
-// ei_pstore(&a[i+7*PacketSize], ei_padd(ei_pload(&a[i+7*PacketSize]), ei_pload(&b[i+7*PacketSize])));
+ ei_pstore(&a[i+2*PacketSize], ei_padd(ei_ploadu(&a[i+2*PacketSize]), ei_ploadu(&b[i+2*PacketSize])));
+ ei_pstore(&a[i+3*PacketSize], ei_padd(ei_ploadu(&a[i+3*PacketSize]), ei_ploadu(&b[i+3*PacketSize])));
+ ei_pstore(&a[i+4*PacketSize], ei_padd(ei_ploadu(&a[i+4*PacketSize]), ei_ploadu(&b[i+4*PacketSize])));
+ ei_pstore(&a[i+5*PacketSize], ei_padd(ei_ploadu(&a[i+5*PacketSize]), ei_ploadu(&b[i+5*PacketSize])));
+ ei_pstore(&a[i+6*PacketSize], ei_padd(ei_ploadu(&a[i+6*PacketSize]), ei_ploadu(&b[i+6*PacketSize])));
+ ei_pstore(&a[i+7*PacketSize], ei_padd(ei_ploadu(&a[i+7*PacketSize]), ei_ploadu(&b[i+7*PacketSize])));
}
}
diff --git a/bench/btl/libs/blitz/btl_tiny_blitz.cpp b/bench/btl/libs/blitz/btl_tiny_blitz.cpp
index eff3c1881..9fddde752 100644
--- a/bench/btl/libs/blitz/btl_tiny_blitz.cpp
+++ b/bench/btl/libs/blitz/btl_tiny_blitz.cpp
@@ -23,7 +23,6 @@
#include "action_matrix_vector_product.hh"
#include "action_matrix_matrix_product.hh"
#include "action_axpy.hh"
-#include "timers/x86_perf_analyzer.hh"
BTL_MAIN;
diff --git a/bench/btl/libs/eigen2/eigen2_interface.hh b/bench/btl/libs/eigen2/eigen2_interface.hh
index 64314ab3b..f3ba5e8af 100644
--- a/bench/btl/libs/eigen2/eigen2_interface.hh
+++ b/bench/btl/libs/eigen2/eigen2_interface.hh
@@ -134,7 +134,7 @@ public :
}
static inline void trisolve_lower(const gene_matrix & L, const gene_vector& B, gene_vector& X, int N){
- X = L.template marked<Lower>().inverseProduct(B);
+ X = L.template marked<Lower>().solveTriangular(B);
}
static inline void cholesky(const gene_matrix & X, gene_matrix & C, int N){
@@ -146,6 +146,7 @@ public :
static inline void lu_decomp(const gene_matrix & X, gene_matrix & C, int N){
C = X.lu().matrixLU();
+// C = X.inverse();
}
static inline void tridiagonalization(const gene_matrix & X, gene_matrix & C, int N){
diff --git a/doc/snippets/MatrixBase_marked.cpp b/doc/snippets/MatrixBase_marked.cpp
index f367def43..d4c11c6fc 100644
--- a/doc/snippets/MatrixBase_marked.cpp
+++ b/doc/snippets/MatrixBase_marked.cpp
@@ -6,4 +6,4 @@ n.part<Eigen::Lower>() *= 2;
cout << "Here is the matrix n:" << endl << n << endl;
cout << "And now here is m.inverse()*n, taking advantage of the fact that"
" m is upper-triangular:" << endl
- << m.marked<Eigen::Upper>().inverseProduct(n);
+ << m.marked<Eigen::Upper>().solveTriangular(n);
diff --git a/test/submatrices.cpp b/test/submatrices.cpp
index f0475f99f..cc00ae810 100644
--- a/test/submatrices.cpp
+++ b/test/submatrices.cpp
@@ -66,6 +66,7 @@ template<typename MatrixType> void submatrices(const MatrixType& m)
m2 = MatrixType::Random(rows, cols),
m3(rows, cols),
mzero = MatrixType::Zero(rows, cols),
+ ones = MatrixType::Ones(rows, cols),
identity = Matrix<Scalar, MatrixType::RowsAtCompileTime, MatrixType::RowsAtCompileTime>
::Identity(rows, rows),
square = Matrix<Scalar, MatrixType::RowsAtCompileTime, MatrixType::RowsAtCompileTime>
@@ -121,6 +122,14 @@ template<typename MatrixType> void submatrices(const MatrixType& m)
Matrix<Scalar,Dynamic,Dynamic> b = m1.template block<BlockRows,BlockCols>(3,3);
VERIFY_IS_APPROX(b, m1.block(3,3,BlockRows,BlockCols));
}
+
+ // stress some basic stuffs with block matrices
+ VERIFY(ones.col(c1).sum() == Scalar(rows));
+ VERIFY(ones.row(r1).sum() == Scalar(cols));
+
+ VERIFY(ones.col(c1).dot(ones.col(c2)) == Scalar(rows));
+ std::cerr << ones.row(r1).dot(ones.row(r2)) << " == " << cols << "\n";
+ VERIFY(ones.row(r1).dot(ones.row(r2)) == Scalar(cols));
}
void test_submatrices()
@@ -131,5 +140,6 @@ void test_submatrices()
CALL_SUBTEST( submatrices(MatrixXcf(3, 3)) );
CALL_SUBTEST( submatrices(MatrixXi(8, 12)) );
CALL_SUBTEST( submatrices(MatrixXcd(20, 20)) );
+ CALL_SUBTEST( submatrices(MatrixXf(20, 20)) );
}
}
diff --git a/test/triangular.cpp b/test/triangular.cpp
index a1e5383bc..de3b85537 100644
--- a/test/triangular.cpp
+++ b/test/triangular.cpp
@@ -80,15 +80,15 @@ template<typename MatrixType> void triangular(const MatrixType& m)
// test back and forward subsitution
m3 = m1.template part<Eigen::Lower>();
- VERIFY(m3.template marked<Eigen::Lower>().inverseProduct(m3).cwise().abs().isIdentity(test_precision<RealScalar>()));
+ VERIFY(m3.template marked<Eigen::Lower>().solveTriangular(m3).cwise().abs().isIdentity(test_precision<RealScalar>()));
m3 = m1.template part<Eigen::Upper>();
- VERIFY(m3.template marked<Eigen::Upper>().inverseProduct(m3).cwise().abs().isIdentity(test_precision<RealScalar>()));
+ VERIFY(m3.template marked<Eigen::Upper>().solveTriangular(m3).cwise().abs().isIdentity(test_precision<RealScalar>()));
// FIXME these tests failed due to numerical issues
// m1 = MatrixType::Random(rows, cols);
- // VERIFY_IS_APPROX(m1.template part<Eigen::Upper>().eval() * (m1.template part<Eigen::Upper>().inverseProduct(m2)), m2);
- // VERIFY_IS_APPROX(m1.template part<Eigen::Lower>().eval() * (m1.template part<Eigen::Lower>().inverseProduct(m2)), m2);
+ // VERIFY_IS_APPROX(m1.template part<Eigen::Upper>().eval() * (m1.template part<Eigen::Upper>().solveTriangular(m2)), m2);
+ // VERIFY_IS_APPROX(m1.template part<Eigen::Lower>().eval() * (m1.template part<Eigen::Lower>().solveTriangular(m2)), m2);
VERIFY((m1.template part<Eigen::Upper>() * m2.template part<Eigen::Upper>()).isUpper());
generated by cgit v1.2.3 (git 2.39.1) at 2024-05-23 08:44:17 +0000