From 217d984abc2dd619d9ba0c585f27065e40380fe3 Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Wed, 11 May 2016 10:22:15 -0700 Subject: Fixed a typo in my previous commit --- test/CMakeLists.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'test/CMakeLists.txt') diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 7bed6a45c..3c0f9b685 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -147,7 +147,7 @@ ei_add_test(nomalloc) ei_add_test(first_aligned) ei_add_test(nullary) ei_add_test(mixingtypes) -ei_add_test(packetmath) +ei_add_test(packetmath "-DEIGEN_FAST_MATH=1") ei_add_test(unalignedassert) ei_add_test(vectorization_logic) ei_add_test(basicstuff) -- cgit v1.2.3 From 5a71eb598585f0413ed630f5007f3c0295c7ba60 Mon Sep 17 00:00:00 2001 From: Gael Guennebaud Date: Wed, 18 May 2016 14:03:03 +0200 Subject: Big 1213: add regression unit test. --- test/CMakeLists.txt | 2 ++ test/bug1213.cpp | 13 +++++++++++++ test/bug1213.h | 8 ++++++++ test/bug1213_main.cpp | 18 ++++++++++++++++++ 4 files changed, 41 insertions(+) create mode 100644 test/bug1213.cpp create mode 100644 test/bug1213.h create mode 100644 test/bug1213_main.cpp (limited to 'test/CMakeLists.txt') diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 3c0f9b685..9d49f1e97 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -259,6 +259,8 @@ ei_add_test(dense_storage) ei_add_test(ctorleak) ei_add_test(mpl2only) +add_executable(bug1213 bug1213.cpp bug1213_main.cpp) + check_cxx_compiler_flag("-ffast-math" COMPILER_SUPPORT_FASTMATH) if(COMPILER_SUPPORT_FASTMATH) set(EIGEN_FASTMATH_FLAGS "-ffast-math") diff --git a/test/bug1213.cpp b/test/bug1213.cpp new file mode 100644 index 000000000..581760c1a --- /dev/null +++ b/test/bug1213.cpp @@ -0,0 +1,13 @@ + +// This anonymous enum is essential to trigger the linking issue +enum { + Foo +}; + +#include "bug1213.h" + +bool bug1213_1(const Eigen::Vector3f& x) +{ + return bug1213_2(x); +} + diff --git a/test/bug1213.h b/test/bug1213.h new file mode 100644 index 000000000..040e5a470 --- /dev/null +++ b/test/bug1213.h @@ -0,0 +1,8 @@ + +#include + +template +bool bug1213_2(const Eigen::Matrix& x); + +bool bug1213_1(const Eigen::Vector3f& x); + diff --git a/test/bug1213_main.cpp b/test/bug1213_main.cpp new file mode 100644 index 000000000..4802c0003 --- /dev/null +++ b/test/bug1213_main.cpp @@ -0,0 +1,18 @@ + +// This is a regression unit regarding a weird linking issue with gcc. + +#include "bug1213.h" + +int main() +{ + return 0; +} + + +template +bool bug1213_2(const Eigen::Matrix& ) +{ + return true; +} + +template bool bug1213_2(const Eigen::Vector3f&); -- cgit v1.2.3 From 32a41ee659686fe1fb76156f7a55287acf14d4bb Mon Sep 17 00:00:00 2001 From: Gael Guennebaud Date: Mon, 4 Jul 2016 15:13:35 +0200 Subject: bug #707: add inplace decomposition through Ref<> for Cholesky, LU and QR decompositions. --- Eigen/src/Cholesky/LDLT.h | 21 ++++- Eigen/src/Cholesky/LLT.h | 16 +++- Eigen/src/LU/FullPivLU.h | 44 +++++++--- Eigen/src/LU/PartialPivLU.h | 62 ++++++++++++-- Eigen/src/QR/ColPivHouseholderQR.h | 37 ++++++--- Eigen/src/QR/CompleteOrthogonalDecomposition.h | 47 +++++++---- Eigen/src/QR/FullPivHouseholderQR.h | 27 ++++-- Eigen/src/QR/HouseholderQR.h | 36 ++++++-- test/CMakeLists.txt | 1 + test/inplace_decomposition.cpp | 110 +++++++++++++++++++++++++ 10 files changed, 337 insertions(+), 64 deletions(-) create mode 100644 test/inplace_decomposition.cpp (limited to 'test/CMakeLists.txt') diff --git a/Eigen/src/Cholesky/LDLT.h b/Eigen/src/Cholesky/LDLT.h index 538aff956..a31b3d6aa 100644 --- a/Eigen/src/Cholesky/LDLT.h +++ b/Eigen/src/Cholesky/LDLT.h @@ -52,7 +52,6 @@ template class LDLT enum { RowsAtCompileTime = MatrixType::RowsAtCompileTime, ColsAtCompileTime = MatrixType::ColsAtCompileTime, - Options = MatrixType::Options & ~RowMajorBit, // these are the options for the TmpMatrixType, we need a ColMajor matrix here! MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime, MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime, UpLo = _UpLo @@ -61,7 +60,7 @@ template class LDLT typedef typename NumTraits::Real RealScalar; typedef Eigen::Index Index; ///< \deprecated since Eigen 3.3 typedef typename MatrixType::StorageIndex StorageIndex; - typedef Matrix TmpMatrixType; + typedef Matrix TmpMatrixType; typedef Transpositions TranspositionType; typedef PermutationMatrix PermutationType; @@ -97,6 +96,7 @@ template class LDLT /** \brief Constructor with decomposition * * This calculates the decomposition for the input \a matrix. + * * \sa LDLT(Index size) */ template @@ -110,6 +110,23 @@ template class LDLT compute(matrix.derived()); } + /** \brief Constructs a LDLT factorization from a given matrix + * + * This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref. + * + * \sa LDLT(const EigenBase&) + */ + template + explicit LDLT(EigenBase& matrix) + : m_matrix(matrix.derived()), + m_transpositions(matrix.rows()), + m_temporary(matrix.rows()), + m_sign(internal::ZeroSign), + m_isInitialized(false) + { + compute(matrix.derived()); + } + /** Clear any existing decomposition * \sa rankUpdate(w,sigma) */ diff --git a/Eigen/src/Cholesky/LLT.h b/Eigen/src/Cholesky/LLT.h index 19578b216..ad163c749 100644 --- a/Eigen/src/Cholesky/LLT.h +++ b/Eigen/src/Cholesky/LLT.h @@ -54,7 +54,6 @@ template class LLT enum { RowsAtCompileTime = MatrixType::RowsAtCompileTime, ColsAtCompileTime = MatrixType::ColsAtCompileTime, - Options = MatrixType::Options, MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime }; typedef typename MatrixType::Scalar Scalar; @@ -95,6 +94,21 @@ template class LLT compute(matrix.derived()); } + /** \brief Constructs a LDLT factorization from a given matrix + * + * This overloaded constructor is provided for inplace solving when + * \c MatrixType is a Eigen::Ref. + * + * \sa LLT(const EigenBase&) + */ + template + explicit LLT(EigenBase& matrix) + : m_matrix(matrix.derived()), + m_isInitialized(false) + { + compute(matrix.derived()); + } + /** \returns a view of the upper triangular matrix U */ inline typename Traits::MatrixU matrixU() const { diff --git a/Eigen/src/LU/FullPivLU.h b/Eigen/src/LU/FullPivLU.h index 2d01b18c6..113b8c7b8 100644 --- a/Eigen/src/LU/FullPivLU.h +++ b/Eigen/src/LU/FullPivLU.h @@ -97,6 +97,15 @@ template class FullPivLU template explicit FullPivLU(const EigenBase& matrix); + /** \brief Constructs a LU factorization from a given matrix + * + * This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref. + * + * \sa FullPivLU(const EigenBase&) + */ + template + explicit FullPivLU(EigenBase& matrix); + /** Computes the LU decomposition of the given matrix. * * \param matrix the matrix of which to compute the LU decomposition. @@ -105,7 +114,11 @@ template class FullPivLU * \returns a reference to *this */ template - FullPivLU& compute(const EigenBase& matrix); + FullPivLU& compute(const EigenBase& matrix) { + m_lu = matrix.derived(); + computeInPlace(); + return *this; + } /** \returns the LU decomposition matrix: the upper-triangular part is U, the * unit-lower-triangular part is L (at least for square matrices; in the non-square @@ -459,25 +472,28 @@ FullPivLU::FullPivLU(const EigenBase& matrix) template template -FullPivLU& FullPivLU::compute(const EigenBase& matrix) +FullPivLU::FullPivLU(EigenBase& matrix) + : m_lu(matrix.derived()), + m_p(matrix.rows()), + m_q(matrix.cols()), + m_rowsTranspositions(matrix.rows()), + m_colsTranspositions(matrix.cols()), + m_isInitialized(false), + m_usePrescribedThreshold(false) { - check_template_parameters(); - - // the permutations are stored as int indices, so just to be sure: - eigen_assert(matrix.rows()<=NumTraits::highest() && matrix.cols()<=NumTraits::highest()); - - m_lu = matrix.derived(); - m_l1_norm = m_lu.cwiseAbs().colwise().sum().maxCoeff(); - computeInPlace(); - - m_isInitialized = true; - return *this; } template void FullPivLU::computeInPlace() { + check_template_parameters(); + + // the permutations are stored as int indices, so just to be sure: + eigen_assert(m_lu.rows()<=NumTraits::highest() && m_lu.cols()<=NumTraits::highest()); + + m_l1_norm = m_lu.cwiseAbs().colwise().sum().maxCoeff(); + const Index size = m_lu.diagonalSize(); const Index rows = m_lu.rows(); const Index cols = m_lu.cols(); @@ -557,6 +573,8 @@ void FullPivLU::computeInPlace() m_q.applyTranspositionOnTheRight(k, m_colsTranspositions.coeff(k)); m_det_pq = (number_of_transpositions%2) ? -1 : 1; + + m_isInitialized = true; } template diff --git a/Eigen/src/LU/PartialPivLU.h b/Eigen/src/LU/PartialPivLU.h index ac2902261..c862d9692 100644 --- a/Eigen/src/LU/PartialPivLU.h +++ b/Eigen/src/LU/PartialPivLU.h @@ -26,6 +26,17 @@ template struct traits > }; }; +template +struct enable_if_ref; +// { +// typedef Derived type; +// }; + +template +struct enable_if_ref,Derived> { + typedef Derived type; +}; + } // end namespace internal /** \ingroup LU_Module @@ -102,8 +113,29 @@ template class PartialPivLU template explicit PartialPivLU(const EigenBase& matrix); + /** Constructor for inplace decomposition + * + * \param matrix the matrix of which to compute the LU decomposition. + * + * If \c MatrixType is an Eigen::Ref, then the storage of \a matrix will be shared + * between \a matrix and \c *this and the decomposition will take place in-place. + * The memory of \a matrix will be used througrough the lifetime of \c *this. In + * particular, further calls to \c this->compute(A) will still operate on the memory + * of \a matrix meaning. This also implies that the sizes of \c A must match the + * ones of \a matrix. + * + * \warning The matrix should have full rank (e.g. if it's square, it should be invertible). + * If you need to deal with non-full rank, use class FullPivLU instead. + */ + template + explicit PartialPivLU(EigenBase& matrix); + template - PartialPivLU& compute(const EigenBase& matrix); + PartialPivLU& compute(const EigenBase& matrix) { + m_lu = matrix.derived(); + compute(); + return *this; + } /** \returns the LU decomposition matrix: the upper-triangular part is U, the * unit-lower-triangular part is L (at least for square matrices; in the non-square @@ -251,6 +283,8 @@ template class PartialPivLU EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); } + void compute(); + MatrixType m_lu; PermutationType m_p; TranspositionType m_rowsTranspositions; @@ -284,7 +318,7 @@ PartialPivLU::PartialPivLU(Index size) template template PartialPivLU::PartialPivLU(const EigenBase& matrix) - : m_lu(matrix.rows(), matrix.rows()), + : m_lu(matrix.rows(),matrix.cols()), m_p(matrix.rows()), m_rowsTranspositions(matrix.rows()), m_l1_norm(0), @@ -294,6 +328,19 @@ PartialPivLU::PartialPivLU(const EigenBase& matrix) compute(matrix.derived()); } +template +template +PartialPivLU::PartialPivLU(EigenBase& matrix) + : m_lu(matrix.derived()), + m_p(matrix.rows()), + m_rowsTranspositions(matrix.rows()), + m_l1_norm(0), + m_det_p(0), + m_isInitialized(false) +{ + compute(); +} + namespace internal { /** \internal This is the blocked version of fullpivlu_unblocked() */ @@ -470,19 +517,17 @@ void partial_lu_inplace(MatrixType& lu, TranspositionType& row_transpositions, t } // end namespace internal template -template -PartialPivLU& PartialPivLU::compute(const EigenBase& matrix) +void PartialPivLU::compute() { check_template_parameters(); // the row permutation is stored as int indices, so just to be sure: - eigen_assert(matrix.rows()::highest()); + eigen_assert(m_lu.rows()::highest()); - m_lu = matrix.derived(); m_l1_norm = m_lu.cwiseAbs().colwise().sum().maxCoeff(); - eigen_assert(matrix.rows() == matrix.cols() && "PartialPivLU is only for square (and moreover invertible) matrices"); - const Index size = matrix.rows(); + eigen_assert(m_lu.rows() == m_lu.cols() && "PartialPivLU is only for square (and moreover invertible) matrices"); + const Index size = m_lu.rows(); m_rowsTranspositions.resize(size); @@ -493,7 +538,6 @@ PartialPivLU& PartialPivLU::compute(const EigenBase diff --git a/Eigen/src/QR/ColPivHouseholderQR.h b/Eigen/src/QR/ColPivHouseholderQR.h index e847bc434..db50b5675 100644 --- a/Eigen/src/QR/ColPivHouseholderQR.h +++ b/Eigen/src/QR/ColPivHouseholderQR.h @@ -51,7 +51,6 @@ template class ColPivHouseholderQR enum { RowsAtCompileTime = MatrixType::RowsAtCompileTime, ColsAtCompileTime = MatrixType::ColsAtCompileTime, - Options = MatrixType::Options, MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime, MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime }; @@ -59,7 +58,6 @@ template class ColPivHouseholderQR typedef typename MatrixType::RealScalar RealScalar; // FIXME should be int typedef typename MatrixType::StorageIndex StorageIndex; - typedef Matrix MatrixQType; typedef typename internal::plain_diag_type::type HCoeffsType; typedef PermutationMatrix PermutationType; typedef typename internal::plain_row_type::type IntRowVectorType; @@ -135,6 +133,27 @@ template class ColPivHouseholderQR compute(matrix.derived()); } + /** \brief Constructs a QR factorization from a given matrix + * + * This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref. + * + * \sa ColPivHouseholderQR(const EigenBase&) + */ + template + explicit ColPivHouseholderQR(EigenBase& matrix) + : m_qr(matrix.derived()), + m_hCoeffs((std::min)(matrix.rows(),matrix.cols())), + m_colsPermutation(PermIndexType(matrix.cols())), + m_colsTranspositions(matrix.cols()), + m_temp(matrix.cols()), + m_colNormsUpdated(matrix.cols()), + m_colNormsDirect(matrix.cols()), + m_isInitialized(false), + m_usePrescribedThreshold(false) + { + computeInPlace(); + } + /** This method finds a solution x to the equation Ax=b, where A is the matrix of which * *this is the QR decomposition, if any exists. * @@ -453,21 +472,19 @@ template template ColPivHouseholderQR& ColPivHouseholderQR::compute(const EigenBase& matrix) { - check_template_parameters(); - - // the column permutation is stored as int indices, so just to be sure: - eigen_assert(matrix.cols()<=NumTraits::highest()); - - m_qr = matrix; - + m_qr = matrix.derived(); computeInPlace(); - return *this; } template void ColPivHouseholderQR::computeInPlace() { + check_template_parameters(); + + // the column permutation is stored as int indices, so just to be sure: + eigen_assert(m_qr.cols()<=NumTraits::highest()); + using std::abs; Index rows = m_qr.rows(); diff --git a/Eigen/src/QR/CompleteOrthogonalDecomposition.h b/Eigen/src/QR/CompleteOrthogonalDecomposition.h index 398e1aa77..967eb35dd 100644 --- a/Eigen/src/QR/CompleteOrthogonalDecomposition.h +++ b/Eigen/src/QR/CompleteOrthogonalDecomposition.h @@ -48,16 +48,12 @@ class CompleteOrthogonalDecomposition { enum { RowsAtCompileTime = MatrixType::RowsAtCompileTime, ColsAtCompileTime = MatrixType::ColsAtCompileTime, - Options = MatrixType::Options, MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime, MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime }; typedef typename MatrixType::Scalar Scalar; typedef typename MatrixType::RealScalar RealScalar; typedef typename MatrixType::StorageIndex StorageIndex; - typedef Matrix - MatrixQType; typedef typename internal::plain_diag_type::type HCoeffsType; typedef PermutationMatrix PermutationType; @@ -114,10 +110,27 @@ class CompleteOrthogonalDecomposition { explicit CompleteOrthogonalDecomposition(const EigenBase& matrix) : m_cpqr(matrix.rows(), matrix.cols()), m_zCoeffs((std::min)(matrix.rows(), matrix.cols())), - m_temp(matrix.cols()) { + m_temp(matrix.cols()) + { compute(matrix.derived()); } + /** \brief Constructs a complete orthogonal decomposition from a given matrix + * + * This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref. + * + * \sa CompleteOrthogonalDecomposition(const EigenBase&) + */ + template + explicit CompleteOrthogonalDecomposition(EigenBase& matrix) + : m_cpqr(matrix.derived()), + m_zCoeffs((std::min)(matrix.rows(), matrix.cols())), + m_temp(matrix.cols()) + { + computeInPlace(); + } + + /** This method computes the minimum-norm solution X to a least squares * problem \f[\mathrm{minimize} ||A X - B|| \f], where \b A is the matrix of * which \c *this is the complete orthogonal decomposition. @@ -165,7 +178,12 @@ class CompleteOrthogonalDecomposition { const MatrixType& matrixT() const { return m_cpqr.matrixQR(); } template - CompleteOrthogonalDecomposition& compute(const EigenBase& matrix); + CompleteOrthogonalDecomposition& compute(const EigenBase& matrix) { + // Compute the column pivoted QR factorization A P = Q R. + m_cpqr.compute(matrix); + computeInPlace(); + return *this; + } /** \returns a const reference to the column permutation matrix */ const PermutationType& colsPermutation() const { @@ -354,6 +372,8 @@ class CompleteOrthogonalDecomposition { EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); } + void computeInPlace(); + /** Overwrites \b rhs with \f$ \mathbf{Z}^* * \mathbf{rhs} \f$. */ template @@ -384,20 +404,16 @@ CompleteOrthogonalDecomposition::logAbsDeterminant() const { * CompleteOrthogonalDecomposition(const MatrixType&) */ template -template -CompleteOrthogonalDecomposition& CompleteOrthogonalDecomposition< - MatrixType>::compute(const EigenBase& matrix) { +void CompleteOrthogonalDecomposition::computeInPlace() +{ check_template_parameters(); // the column permutation is stored as int indices, so just to be sure: - eigen_assert(matrix.cols() <= NumTraits::highest()); - - // Compute the column pivoted QR factorization A P = Q R. - m_cpqr.compute(matrix); + eigen_assert(m_cpqr.cols() <= NumTraits::highest()); const Index rank = m_cpqr.rank(); - const Index cols = matrix.cols(); - const Index rows = matrix.rows(); + const Index cols = m_cpqr.cols(); + const Index rows = m_cpqr.rows(); m_zCoeffs.resize((std::min)(rows, cols)); m_temp.resize(cols); @@ -443,7 +459,6 @@ CompleteOrthogonalDecomposition& CompleteOrthogonalDecomposition< } } } - return *this; } template diff --git a/Eigen/src/QR/FullPivHouseholderQR.h b/Eigen/src/QR/FullPivHouseholderQR.h index e21966056..5c2f57d04 100644 --- a/Eigen/src/QR/FullPivHouseholderQR.h +++ b/Eigen/src/QR/FullPivHouseholderQR.h @@ -60,7 +60,6 @@ template class FullPivHouseholderQR enum { RowsAtCompileTime = MatrixType::RowsAtCompileTime, ColsAtCompileTime = MatrixType::ColsAtCompileTime, - Options = MatrixType::Options, MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime, MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime }; @@ -135,6 +134,26 @@ template class FullPivHouseholderQR compute(matrix.derived()); } + /** \brief Constructs a QR factorization from a given matrix + * + * This overloaded constructor is provided for inplace solving when \c MatrixType is a Eigen::Ref. + * + * \sa FullPivHouseholderQR(const EigenBase&) + */ + template + explicit FullPivHouseholderQR(EigenBase& matrix) + : m_qr(matrix.derived()), + m_hCoeffs((std::min)(matrix.rows(), matrix.cols())), + m_rows_transpositions((std::min)(matrix.rows(), matrix.cols())), + m_cols_transpositions((std::min)(matrix.rows(), matrix.cols())), + m_cols_permutation(matrix.cols()), + m_temp(matrix.cols()), + m_isInitialized(false), + m_usePrescribedThreshold(false) + { + computeInPlace(); + } + /** This method finds a solution x to the equation Ax=b, where A is the matrix of which * \c *this is the QR decomposition. * @@ -430,18 +449,16 @@ template template FullPivHouseholderQR& FullPivHouseholderQR::compute(const EigenBase& matrix) { - check_template_parameters(); - m_qr = matrix.derived(); - computeInPlace(); - return *this; } template void FullPivHouseholderQR::computeInPlace() { + check_template_parameters(); + using std::abs; Index rows = m_qr.rows(); Index cols = m_qr.cols(); diff --git a/Eigen/src/QR/HouseholderQR.h b/Eigen/src/QR/HouseholderQR.h index 03bc8e6cd..f2a9cc080 100644 --- a/Eigen/src/QR/HouseholderQR.h +++ b/Eigen/src/QR/HouseholderQR.h @@ -47,7 +47,6 @@ template class HouseholderQR enum { RowsAtCompileTime = MatrixType::RowsAtCompileTime, ColsAtCompileTime = MatrixType::ColsAtCompileTime, - Options = MatrixType::Options, MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime, MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime }; @@ -102,6 +101,24 @@ template class HouseholderQR compute(matrix.derived()); } + + /** \brief Constructs a QR factorization from a given matrix + * + * This overloaded constructor is provided for inplace solving when + * \c MatrixType is a Eigen::Ref. + * + * \sa HouseholderQR(const EigenBase&) + */ + template + explicit HouseholderQR(EigenBase& matrix) + : m_qr(matrix.derived()), + m_hCoeffs((std::min)(matrix.rows(),matrix.cols())), + m_temp(matrix.cols()), + m_isInitialized(false) + { + computeInPlace(); + } + /** This method finds a solution x to the equation Ax=b, where A is the matrix of which * *this is the QR decomposition, if any exists. * @@ -151,7 +168,11 @@ template class HouseholderQR } template - HouseholderQR& compute(const EigenBase& matrix); + HouseholderQR& compute(const EigenBase& matrix) { + m_qr = matrix.derived(); + computeInPlace(); + return *this; + } /** \returns the absolute value of the determinant of the matrix of which * *this is the QR decomposition. It has only linear complexity @@ -203,6 +224,8 @@ template class HouseholderQR { EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); } + + void computeInPlace(); MatrixType m_qr; HCoeffsType m_hCoeffs; @@ -354,16 +377,14 @@ void HouseholderQR<_MatrixType>::_solve_impl(const RhsType &rhs, DstType &dst) c * \sa class HouseholderQR, HouseholderQR(const MatrixType&) */ template -template -HouseholderQR& HouseholderQR::compute(const EigenBase& matrix) +void HouseholderQR::computeInPlace() { check_template_parameters(); - Index rows = matrix.rows(); - Index cols = matrix.cols(); + Index rows = m_qr.rows(); + Index cols = m_qr.cols(); Index size = (std::min)(rows,cols); - m_qr = matrix.derived(); m_hCoeffs.resize(size); m_temp.resize(cols); @@ -371,7 +392,6 @@ HouseholderQR& HouseholderQR::compute(const EigenBase::run(m_qr, m_hCoeffs, 48, m_temp.data()); m_isInitialized = true; - return *this; } #ifndef __CUDACC__ diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 9d49f1e97..27a9ec4ac 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -258,6 +258,7 @@ ei_add_test(rvalue_types) ei_add_test(dense_storage) ei_add_test(ctorleak) ei_add_test(mpl2only) +ei_add_test(inplace_decomposition) add_executable(bug1213 bug1213.cpp bug1213_main.cpp) diff --git a/test/inplace_decomposition.cpp b/test/inplace_decomposition.cpp new file mode 100644 index 000000000..4900312f1 --- /dev/null +++ b/test/inplace_decomposition.cpp @@ -0,0 +1,110 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2016 Gael Guennebaud +// +// This Source Code Form is subject to the terms of the Mozilla +// Public License v. 2.0. If a copy of the MPL was not distributed +// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. + +#include "main.h" +#include +#include +#include + +// This file test inplace decomposition through Ref<>, as supported by Cholesky, LU, and QR decompositions. + +template void inplace(bool square = false, bool SPD = false) +{ + typedef typename MatrixType::Scalar Scalar; + typedef Matrix RhsType; + typedef Matrix ResType; + + Index rows = MatrixType::RowsAtCompileTime==Dynamic ? internal::random(2,EIGEN_TEST_MAX_SIZE/2) : MatrixType::RowsAtCompileTime; + Index cols = MatrixType::ColsAtCompileTime==Dynamic ? (square?rows:internal::random(2,rows)) : MatrixType::ColsAtCompileTime; + + MatrixType A = MatrixType::Random(rows,cols); + RhsType b = RhsType::Random(rows); + ResType x(cols); + + if(SPD) + { + assert(square); + A.topRows(cols) = A.topRows(cols).adjoint() * A.topRows(cols); + A.diagonal().array() += 1e-3; + } + + MatrixType A0 = A; + MatrixType A1 = A; + + DecType dec(A); + + // Check that the content of A has been modified + VERIFY_IS_NOT_APPROX( A, A0 ); + + // Check that the decomposition is correct: + if(rows==cols) + { + VERIFY_IS_APPROX( A0 * (x = dec.solve(b)), b ); + } + else + { + VERIFY_IS_APPROX( A0.transpose() * A0 * (x = dec.solve(b)), A0.transpose() * b ); + } + + // Check that modifying A breaks the current dec: + A.setRandom(); + if(rows==cols) + { + VERIFY_IS_NOT_APPROX( A0 * (x = dec.solve(b)), b ); + } + else + { + VERIFY_IS_NOT_APPROX( A0.transpose() * A0 * (x = dec.solve(b)), A0.transpose() * b ); + } + + // Check that calling compute(A1) does not modify A1: + A = A0; + dec.compute(A1); + VERIFY_IS_EQUAL(A0,A1); + VERIFY_IS_NOT_APPROX( A, A0 ); + if(rows==cols) + { + VERIFY_IS_APPROX( A0 * (x = dec.solve(b)), b ); + } + else + { + VERIFY_IS_APPROX( A0.transpose() * A0 * (x = dec.solve(b)), A0.transpose() * b ); + } +} + + +void test_inplace_decomposition() +{ + EIGEN_UNUSED typedef Matrix Matrix43d; + for(int i = 0; i < g_repeat; i++) { + CALL_SUBTEST_1(( inplace >, MatrixXd>(true,true) )); + CALL_SUBTEST_1(( inplace >, Matrix4d>(true,true) )); + + CALL_SUBTEST_2(( inplace >, MatrixXd>(true,true) )); + CALL_SUBTEST_2(( inplace >, Matrix4d>(true,true) )); + + CALL_SUBTEST_3(( inplace >, MatrixXd>(true,false) )); + CALL_SUBTEST_3(( inplace >, Matrix4d>(true,false) )); + + CALL_SUBTEST_4(( inplace >, MatrixXd>(true,false) )); + CALL_SUBTEST_4(( inplace >, Matrix4d>(true,false) )); + + CALL_SUBTEST_5(( inplace >, MatrixXd>(false,false) )); + CALL_SUBTEST_5(( inplace >, Matrix43d>(false,false) )); + + CALL_SUBTEST_6(( inplace >, MatrixXd>(false,false) )); + CALL_SUBTEST_6(( inplace >, Matrix43d>(false,false) )); + + CALL_SUBTEST_7(( inplace >, MatrixXd>(false,false) )); + CALL_SUBTEST_7(( inplace >, Matrix43d>(false,false) )); + + CALL_SUBTEST_8(( inplace >, MatrixXd>(false,false) )); + CALL_SUBTEST_8(( inplace >, Matrix43d>(false,false) )); + } +} -- cgit v1.2.3 From b64b9d017267b7a0e057a41171574e7dc9cd428b Mon Sep 17 00:00:00 2001 From: Gael Guennebaud Date: Wed, 20 Jul 2016 15:20:14 +0200 Subject: Add a unit test to stress our solvers with Boost.Multiprecision --- test/CMakeLists.txt | 10 ++++ test/boostmultiprec.cpp | 137 ++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 147 insertions(+) create mode 100644 test/boostmultiprec.cpp (limited to 'test/CMakeLists.txt') diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 27a9ec4ac..abc847794 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -327,6 +327,16 @@ if(EIGEN_TEST_EIGEN2) message(WARNING "The Eigen2 test suite has been removed") endif() +# boost MP unit test +find_package(Boost) +if(Boost_FOUND) + include_directories(${Boost_INCLUDE_DIRS}) + ei_add_test(boostmultiprec "" "${Boost_LIBRARIES}") + ei_add_property(EIGEN_TESTED_BACKENDS "Boost.Multiprecision, ") +else() + ei_add_property(EIGEN_MISSING_BACKENDS "Boost.Multiprecision, ") +endif() + # CUDA unit tests option(EIGEN_TEST_CUDA "Enable CUDA support in unit tests" OFF) diff --git a/test/boostmultiprec.cpp b/test/boostmultiprec.cpp new file mode 100644 index 000000000..a75711ff4 --- /dev/null +++ b/test/boostmultiprec.cpp @@ -0,0 +1,137 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2016 Gael Guennebaud +// +// This Source Code Form is subject to the terms of the Mozilla +// Public License v. 2.0. If a copy of the MPL was not distributed +// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. + +#ifdef EIGEN_TEST_MAX_SIZE +#undef EIGEN_TEST_MAX_SIZE +#endif + +#define EIGEN_TEST_MAX_SIZE 50 + +#ifdef EIGEN_TEST_PART_1 +#include "cholesky.cpp" +#endif + +#ifdef EIGEN_TEST_PART_2 +#include "lu.cpp" +#endif + +#ifdef EIGEN_TEST_PART_3 +#include "qr.cpp" +#endif + +#ifdef EIGEN_TEST_PART_4 +#include "qr_colpivoting.cpp" +#endif + +#ifdef EIGEN_TEST_PART_5 +#include "qr_fullpivoting.cpp" +#endif + +#ifdef EIGEN_TEST_PART_6 +#include "eigensolver_selfadjoint.cpp" +#endif + +#ifdef EIGEN_TEST_PART_7 +#include "jacobisvd.cpp" +#endif + +#ifdef EIGEN_TEST_PART_8 +#include "bdcsvd.cpp" +#endif + +#include + +#undef min +#undef max +#undef isnan +#undef isinf +#undef isfinite + +#include +#include + +namespace mp = boost::multiprecision; +typedef mp::number, mp::et_off> Real; + +namespace Eigen { + template<> struct NumTraits : GenericNumTraits { + static inline Real dummy_precision() { return 1e-50; } + }; + + template<> + Real test_precision() { return 1e-50; } +} + +namespace boost { +namespace multiprecision { + // to make ADL works as expected: + using boost::math::isfinite; + + // some specialization for the unit tests: + inline bool test_isMuchSmallerThan(const Real& a, const Real& b) { + return internal::isMuchSmallerThan(a, b, test_precision()); + } + + inline bool test_isApprox(const Real& a, const Real& b) { + return internal::isApprox(a, b, test_precision()); + } + + inline bool test_isApproxOrLessThan(const Real& a, const Real& b) + { return internal::isApproxOrLessThan(a, b, test_precision()); } + + Real get_test_precision(const Real*) { + return test_precision(); + } + + Real test_relative_error(const Real &a, const Real &b) { + return Eigen::numext::sqrt(Real(Eigen::numext::abs2(a-b))/Real((Eigen::numext::mini)(Eigen::numext::abs2(a),Eigen::numext::abs2(b)))); + } +} +} + +namespace Eigen { + +} + +void test_boostmultiprec() +{ + typedef Matrix Mat; + + std::cout << "NumTraits::epsilon() = " << NumTraits::epsilon() << std::endl; + std::cout << "NumTraits::dummy_precision() = " << NumTraits::dummy_precision() << std::endl; + std::cout << "NumTraits::lowest() = " << NumTraits::lowest() << std::endl; + std::cout << "NumTraits::highest() = " << NumTraits::highest() << std::endl; + + for(int i = 0; i < g_repeat; i++) { + int s = internal::random(1,EIGEN_TEST_MAX_SIZE); + + CALL_SUBTEST_1( cholesky(Mat(s,s)) ); + + CALL_SUBTEST_2( lu_non_invertible() ); + CALL_SUBTEST_2( lu_invertible() ); + + CALL_SUBTEST_3( qr(Mat(internal::random(1,EIGEN_TEST_MAX_SIZE),internal::random(1,EIGEN_TEST_MAX_SIZE))) ); + CALL_SUBTEST_3( qr_invertible() ); + + CALL_SUBTEST_4( qr() ); + CALL_SUBTEST_4( cod() ); + CALL_SUBTEST_4( qr_invertible() ); + + CALL_SUBTEST_5( qr() ); + CALL_SUBTEST_5( qr_invertible() ); + + CALL_SUBTEST_6( selfadjointeigensolver(Mat(s,s)) ); + + TEST_SET_BUT_UNUSED_VARIABLE(s) + } + + CALL_SUBTEST_7(( jacobisvd(Mat(internal::random(EIGEN_TEST_MAX_SIZE/4, EIGEN_TEST_MAX_SIZE), internal::random(EIGEN_TEST_MAX_SIZE/4, EIGEN_TEST_MAX_SIZE/2))) )); + CALL_SUBTEST_8(( bdcsvd(Mat(internal::random(EIGEN_TEST_MAX_SIZE/4, EIGEN_TEST_MAX_SIZE), internal::random(EIGEN_TEST_MAX_SIZE/4, EIGEN_TEST_MAX_SIZE/2))) )); +} + -- cgit v1.2.3 From d075d122ea142c9c4ee2e7837493396f54916bc0 Mon Sep 17 00:00:00 2001 From: Gael Guennebaud Date: Fri, 22 Jul 2016 14:34:19 +0200 Subject: Move half unit test from unsupported to main tests --- test/CMakeLists.txt | 1 + test/half_float.cpp | 242 ++++++++++++++++++++++++++++++++++++ unsupported/test/CMakeLists.txt | 1 - unsupported/test/cxx11_float16.cpp | 244 ------------------------------------- 4 files changed, 243 insertions(+), 245 deletions(-) create mode 100644 test/half_float.cpp delete mode 100644 unsupported/test/cxx11_float16.cpp (limited to 'test/CMakeLists.txt') diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index abc847794..e42e850ca 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -259,6 +259,7 @@ ei_add_test(dense_storage) ei_add_test(ctorleak) ei_add_test(mpl2only) ei_add_test(inplace_decomposition) +ei_add_test(half_float) add_executable(bug1213 bug1213.cpp bug1213_main.cpp) diff --git a/test/half_float.cpp b/test/half_float.cpp new file mode 100644 index 000000000..289e4f662 --- /dev/null +++ b/test/half_float.cpp @@ -0,0 +1,242 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// This Source Code Form is subject to the terms of the Mozilla +// Public License v. 2.0. If a copy of the MPL was not distributed +// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. + +#include + +#include "main.h" + +#include + +using Eigen::half; + +void test_conversion() +{ + using Eigen::half_impl::__half; + + // Conversion from float. + VERIFY_IS_EQUAL(half(1.0f).x, 0x3c00); + VERIFY_IS_EQUAL(half(0.5f).x, 0x3800); + VERIFY_IS_EQUAL(half(0.33333f).x, 0x3555); + VERIFY_IS_EQUAL(half(0.0f).x, 0x0000); + VERIFY_IS_EQUAL(half(-0.0f).x, 0x8000); + VERIFY_IS_EQUAL(half(65504.0f).x, 0x7bff); + VERIFY_IS_EQUAL(half(65536.0f).x, 0x7c00); // Becomes infinity. + + // Denormals. + VERIFY_IS_EQUAL(half(-5.96046e-08f).x, 0x8001); + VERIFY_IS_EQUAL(half(5.96046e-08f).x, 0x0001); + VERIFY_IS_EQUAL(half(1.19209e-07f).x, 0x0002); + + // Verify round-to-nearest-even behavior. + float val1 = float(half(__half(0x3c00))); + float val2 = float(half(__half(0x3c01))); + float val3 = float(half(__half(0x3c02))); + VERIFY_IS_EQUAL(half(0.5f * (val1 + val2)).x, 0x3c00); + VERIFY_IS_EQUAL(half(0.5f * (val2 + val3)).x, 0x3c02); + + // Conversion from int. + VERIFY_IS_EQUAL(half(-1).x, 0xbc00); + VERIFY_IS_EQUAL(half(0).x, 0x0000); + VERIFY_IS_EQUAL(half(1).x, 0x3c00); + VERIFY_IS_EQUAL(half(2).x, 0x4000); + VERIFY_IS_EQUAL(half(3).x, 0x4200); + + // Conversion from bool. + VERIFY_IS_EQUAL(half(false).x, 0x0000); + VERIFY_IS_EQUAL(half(true).x, 0x3c00); + + // Conversion to float. + VERIFY_IS_EQUAL(float(half(__half(0x0000))), 0.0f); + VERIFY_IS_EQUAL(float(half(__half(0x3c00))), 1.0f); + + // Denormals. + VERIFY_IS_APPROX(float(half(__half(0x8001))), -5.96046e-08f); + VERIFY_IS_APPROX(float(half(__half(0x0001))), 5.96046e-08f); + VERIFY_IS_APPROX(float(half(__half(0x0002))), 1.19209e-07f); + + // NaNs and infinities. + VERIFY(!(numext::isinf)(float(half(65504.0f)))); // Largest finite number. + VERIFY(!(numext::isnan)(float(half(0.0f)))); + VERIFY((numext::isinf)(float(half(__half(0xfc00))))); + VERIFY((numext::isnan)(float(half(__half(0xfc01))))); + VERIFY((numext::isinf)(float(half(__half(0x7c00))))); + VERIFY((numext::isnan)(float(half(__half(0x7c01))))); + +#if !EIGEN_COMP_MSVC + // Visual Studio errors out on divisions by 0 + VERIFY((numext::isnan)(float(half(0.0 / 0.0)))); + VERIFY((numext::isinf)(float(half(1.0 / 0.0)))); + VERIFY((numext::isinf)(float(half(-1.0 / 0.0)))); +#endif + + // Exactly same checks as above, just directly on the half representation. + VERIFY(!(numext::isinf)(half(__half(0x7bff)))); + VERIFY(!(numext::isnan)(half(__half(0x0000)))); + VERIFY((numext::isinf)(half(__half(0xfc00)))); + VERIFY((numext::isnan)(half(__half(0xfc01)))); + VERIFY((numext::isinf)(half(__half(0x7c00)))); + VERIFY((numext::isnan)(half(__half(0x7c01)))); + +#if !EIGEN_COMP_MSVC + // Visual Studio errors out on divisions by 0 + VERIFY((numext::isnan)(half(0.0 / 0.0))); + VERIFY((numext::isinf)(half(1.0 / 0.0))); + VERIFY((numext::isinf)(half(-1.0 / 0.0))); +#endif +} + +void test_numtraits() +{ + std::cout << "epsilon = " << NumTraits::epsilon() << std::endl; + std::cout << "highest = " << NumTraits::highest() << std::endl; + std::cout << "lowest = " << NumTraits::lowest() << std::endl; + std::cout << "inifinty = " << NumTraits::infinity() << std::endl; + std::cout << "nan = " << NumTraits::quiet_NaN() << std::endl; + +} + +void test_arithmetic() +{ + VERIFY_IS_EQUAL(float(half(2) + half(2)), 4); + VERIFY_IS_EQUAL(float(half(2) + half(-2)), 0); + VERIFY_IS_APPROX(float(half(0.33333f) + half(0.66667f)), 1.0f); + VERIFY_IS_EQUAL(float(half(2.0f) * half(-5.5f)), -11.0f); + VERIFY_IS_APPROX(float(half(1.0f) / half(3.0f)), 0.33333f); + VERIFY_IS_EQUAL(float(-half(4096.0f)), -4096.0f); + VERIFY_IS_EQUAL(float(-half(-4096.0f)), 4096.0f); +} + +void test_comparison() +{ + VERIFY(half(1.0f) > half(0.5f)); + VERIFY(half(0.5f) < half(1.0f)); + VERIFY(!(half(1.0f) < half(0.5f))); + VERIFY(!(half(0.5f) > half(1.0f))); + + VERIFY(!(half(4.0f) > half(4.0f))); + VERIFY(!(half(4.0f) < half(4.0f))); + + VERIFY(!(half(0.0f) < half(-0.0f))); + VERIFY(!(half(-0.0f) < half(0.0f))); + VERIFY(!(half(0.0f) > half(-0.0f))); + VERIFY(!(half(-0.0f) > half(0.0f))); + + VERIFY(half(0.2f) > half(-1.0f)); + VERIFY(half(-1.0f) < half(0.2f)); + VERIFY(half(-16.0f) < half(-15.0f)); + + VERIFY(half(1.0f) == half(1.0f)); + VERIFY(half(1.0f) != half(2.0f)); + + // Comparisons with NaNs and infinities. +#if !EIGEN_COMP_MSVC + // Visual Studio errors out on divisions by 0 + VERIFY(!(half(0.0 / 0.0) == half(0.0 / 0.0))); + VERIFY(half(0.0 / 0.0) != half(0.0 / 0.0)); + + VERIFY(!(half(1.0) == half(0.0 / 0.0))); + VERIFY(!(half(1.0) < half(0.0 / 0.0))); + VERIFY(!(half(1.0) > half(0.0 / 0.0))); + VERIFY(half(1.0) != half(0.0 / 0.0)); + + VERIFY(half(1.0) < half(1.0 / 0.0)); + VERIFY(half(1.0) > half(-1.0 / 0.0)); +#endif +} + +void test_basic_functions() +{ + VERIFY_IS_EQUAL(float(numext::abs(half(3.5f))), 3.5f); + VERIFY_IS_EQUAL(float(abs(half(3.5f))), 3.5f); + VERIFY_IS_EQUAL(float(numext::abs(half(-3.5f))), 3.5f); + VERIFY_IS_EQUAL(float(abs(half(-3.5f))), 3.5f); + + VERIFY_IS_EQUAL(float(numext::floor(half(3.5f))), 3.0f); + VERIFY_IS_EQUAL(float(floor(half(3.5f))), 3.0f); + VERIFY_IS_EQUAL(float(numext::floor(half(-3.5f))), -4.0f); + VERIFY_IS_EQUAL(float(floor(half(-3.5f))), -4.0f); + + VERIFY_IS_EQUAL(float(numext::ceil(half(3.5f))), 4.0f); + VERIFY_IS_EQUAL(float(ceil(half(3.5f))), 4.0f); + VERIFY_IS_EQUAL(float(numext::ceil(half(-3.5f))), -3.0f); + VERIFY_IS_EQUAL(float(ceil(half(-3.5f))), -3.0f); + + VERIFY_IS_APPROX(float(numext::sqrt(half(0.0f))), 0.0f); + VERIFY_IS_APPROX(float(sqrt(half(0.0f))), 0.0f); + VERIFY_IS_APPROX(float(numext::sqrt(half(4.0f))), 2.0f); + VERIFY_IS_APPROX(float(sqrt(half(4.0f))), 2.0f); + + VERIFY_IS_APPROX(float(numext::pow(half(0.0f), half(1.0f))), 0.0f); + VERIFY_IS_APPROX(float(pow(half(0.0f), half(1.0f))), 0.0f); + VERIFY_IS_APPROX(float(numext::pow(half(2.0f), half(2.0f))), 4.0f); + VERIFY_IS_APPROX(float(pow(half(2.0f), half(2.0f))), 4.0f); + + VERIFY_IS_EQUAL(float(numext::exp(half(0.0f))), 1.0f); + VERIFY_IS_EQUAL(float(exp(half(0.0f))), 1.0f); + VERIFY_IS_APPROX(float(numext::exp(half(EIGEN_PI))), 20.f + float(EIGEN_PI)); + VERIFY_IS_APPROX(float(exp(half(EIGEN_PI))), 20.f + float(EIGEN_PI)); + + VERIFY_IS_EQUAL(float(numext::log(half(1.0f))), 0.0f); + VERIFY_IS_EQUAL(float(log(half(1.0f))), 0.0f); + VERIFY_IS_APPROX(float(numext::log(half(10.0f))), 2.30273f); + VERIFY_IS_APPROX(float(log(half(10.0f))), 2.30273f); +} + +void test_trigonometric_functions() +{ + VERIFY_IS_APPROX(numext::cos(half(0.0f)), half(cosf(0.0f))); + VERIFY_IS_APPROX(cos(half(0.0f)), half(cosf(0.0f))); + VERIFY_IS_APPROX(numext::cos(half(EIGEN_PI)), half(cosf(EIGEN_PI))); + //VERIFY_IS_APPROX(numext::cos(half(EIGEN_PI/2)), half(cosf(EIGEN_PI/2))); + //VERIFY_IS_APPROX(numext::cos(half(3*EIGEN_PI/2)), half(cosf(3*EIGEN_PI/2))); + VERIFY_IS_APPROX(numext::cos(half(3.5f)), half(cosf(3.5f))); + + VERIFY_IS_APPROX(numext::sin(half(0.0f)), half(sinf(0.0f))); + VERIFY_IS_APPROX(sin(half(0.0f)), half(sinf(0.0f))); + // VERIFY_IS_APPROX(numext::sin(half(EIGEN_PI)), half(sinf(EIGEN_PI))); + VERIFY_IS_APPROX(numext::sin(half(EIGEN_PI/2)), half(sinf(EIGEN_PI/2))); + VERIFY_IS_APPROX(numext::sin(half(3*EIGEN_PI/2)), half(sinf(3*EIGEN_PI/2))); + VERIFY_IS_APPROX(numext::sin(half(3.5f)), half(sinf(3.5f))); + + VERIFY_IS_APPROX(numext::tan(half(0.0f)), half(tanf(0.0f))); + VERIFY_IS_APPROX(tan(half(0.0f)), half(tanf(0.0f))); + // VERIFY_IS_APPROX(numext::tan(half(EIGEN_PI)), half(tanf(EIGEN_PI))); + // VERIFY_IS_APPROX(numext::tan(half(EIGEN_PI/2)), half(tanf(EIGEN_PI/2))); + //VERIFY_IS_APPROX(numext::tan(half(3*EIGEN_PI/2)), half(tanf(3*EIGEN_PI/2))); + VERIFY_IS_APPROX(numext::tan(half(3.5f)), half(tanf(3.5f))); +} + +void test_array() +{ + typedef Array ArrayXh; + Index size = internal::random(1,10); + Index i = internal::random(0,size-1); + ArrayXh a1 = ArrayXh::Random(size), a2 = ArrayXh::Random(size); + VERIFY_IS_APPROX( a1+a1, half(2)*a1 ); + VERIFY( (a1.abs() >= half(0)).all() ); + VERIFY_IS_APPROX( (a1*a1).sqrt(), a1.abs() ); + + VERIFY( ((a1.min)(a2) <= (a1.max)(a2)).all() ); + a1(i) = half(-10.); + VERIFY_IS_EQUAL( a1.minCoeff(), half(-10.) ); + a1(i) = half(10.); + VERIFY_IS_EQUAL( a1.maxCoeff(), half(10.) ); + + std::stringstream ss; + ss << a1; +} + +void test_half_float() +{ + CALL_SUBTEST(test_conversion()); + CALL_SUBTEST(test_numtraits()); + CALL_SUBTEST(test_arithmetic()); + CALL_SUBTEST(test_comparison()); + CALL_SUBTEST(test_basic_functions()); + CALL_SUBTEST(test_trigonometric_functions()); + CALL_SUBTEST(test_array()); +} diff --git a/unsupported/test/CMakeLists.txt b/unsupported/test/CMakeLists.txt index c50ea53c5..6188b421a 100644 --- a/unsupported/test/CMakeLists.txt +++ b/unsupported/test/CMakeLists.txt @@ -117,7 +117,6 @@ ei_add_test(special_functions) # when using visual studio. We should make the check more strict to enable the tests for # newer versions of MSVC. if (NOT CMAKE_CXX_COMPILER_ID STREQUAL "MSVC") -ei_add_test(cxx11_float16) ei_add_test(cxx11_tensor_dimension) ei_add_test(cxx11_tensor_map) ei_add_test(cxx11_tensor_assign) diff --git a/unsupported/test/cxx11_float16.cpp b/unsupported/test/cxx11_float16.cpp deleted file mode 100644 index d2b48772e..000000000 --- a/unsupported/test/cxx11_float16.cpp +++ /dev/null @@ -1,244 +0,0 @@ -// This file is part of Eigen, a lightweight C++ template library -// for linear algebra. -// -// This Source Code Form is subject to the terms of the Mozilla -// Public License v. 2.0. If a copy of the MPL was not distributed -// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. - -#define EIGEN_TEST_NO_LONGDOUBLE -#define EIGEN_TEST_NO_COMPLEX -#define EIGEN_TEST_FUNC cxx11_float16 - -#include - -#include "main.h" - -#include - -using Eigen::half; - -void test_conversion() -{ - // Conversion from float. - VERIFY_IS_EQUAL(half(1.0f).x, 0x3c00); - VERIFY_IS_EQUAL(half(0.5f).x, 0x3800); - VERIFY_IS_EQUAL(half(0.33333f).x, 0x3555); - VERIFY_IS_EQUAL(half(0.0f).x, 0x0000); - VERIFY_IS_EQUAL(half(-0.0f).x, 0x8000); - VERIFY_IS_EQUAL(half(65504.0f).x, 0x7bff); - VERIFY_IS_EQUAL(half(65536.0f).x, 0x7c00); // Becomes infinity. - - // Denormals. - VERIFY_IS_EQUAL(half(-5.96046e-08f).x, 0x8001); - VERIFY_IS_EQUAL(half(5.96046e-08f).x, 0x0001); - VERIFY_IS_EQUAL(half(1.19209e-07f).x, 0x0002); - - // Verify round-to-nearest-even behavior. - float val1 = float(half(__half(0x3c00))); - float val2 = float(half(__half(0x3c01))); - float val3 = float(half(__half(0x3c02))); - VERIFY_IS_EQUAL(half(0.5f * (val1 + val2)).x, 0x3c00); - VERIFY_IS_EQUAL(half(0.5f * (val2 + val3)).x, 0x3c02); - - // Conversion from int. - VERIFY_IS_EQUAL(half(-1).x, 0xbc00); - VERIFY_IS_EQUAL(half(0).x, 0x0000); - VERIFY_IS_EQUAL(half(1).x, 0x3c00); - VERIFY_IS_EQUAL(half(2).x, 0x4000); - VERIFY_IS_EQUAL(half(3).x, 0x4200); - - // Conversion from bool. - VERIFY_IS_EQUAL(half(false).x, 0x0000); - VERIFY_IS_EQUAL(half(true).x, 0x3c00); - - // Conversion to float. - VERIFY_IS_EQUAL(float(half(__half(0x0000))), 0.0f); - VERIFY_IS_EQUAL(float(half(__half(0x3c00))), 1.0f); - - // Denormals. - VERIFY_IS_APPROX(float(half(__half(0x8001))), -5.96046e-08f); - VERIFY_IS_APPROX(float(half(__half(0x0001))), 5.96046e-08f); - VERIFY_IS_APPROX(float(half(__half(0x0002))), 1.19209e-07f); - - // NaNs and infinities. - VERIFY(!(numext::isinf)(float(half(65504.0f)))); // Largest finite number. - VERIFY(!(numext::isnan)(float(half(0.0f)))); - VERIFY((numext::isinf)(float(half(__half(0xfc00))))); - VERIFY((numext::isnan)(float(half(__half(0xfc01))))); - VERIFY((numext::isinf)(float(half(__half(0x7c00))))); - VERIFY((numext::isnan)(float(half(__half(0x7c01))))); - -#if !EIGEN_COMP_MSVC - // Visual Studio errors out on divisions by 0 - VERIFY((numext::isnan)(float(half(0.0 / 0.0)))); - VERIFY((numext::isinf)(float(half(1.0 / 0.0)))); - VERIFY((numext::isinf)(float(half(-1.0 / 0.0)))); -#endif - - // Exactly same checks as above, just directly on the half representation. - VERIFY(!(numext::isinf)(half(__half(0x7bff)))); - VERIFY(!(numext::isnan)(half(__half(0x0000)))); - VERIFY((numext::isinf)(half(__half(0xfc00)))); - VERIFY((numext::isnan)(half(__half(0xfc01)))); - VERIFY((numext::isinf)(half(__half(0x7c00)))); - VERIFY((numext::isnan)(half(__half(0x7c01)))); - -#if !EIGEN_COMP_MSVC - // Visual Studio errors out on divisions by 0 - VERIFY((numext::isnan)(half(0.0 / 0.0))); - VERIFY((numext::isinf)(half(1.0 / 0.0))); - VERIFY((numext::isinf)(half(-1.0 / 0.0))); -#endif -} - -void test_numtraits() -{ - std::cout << "epsilon = " << NumTraits::epsilon() << std::endl; - std::cout << "highest = " << NumTraits::highest() << std::endl; - std::cout << "lowest = " << NumTraits::lowest() << std::endl; - std::cout << "inifinty = " << NumTraits::infinity() << std::endl; - std::cout << "nan = " << NumTraits::quiet_NaN() << std::endl; - -} - -void test_arithmetic() -{ - VERIFY_IS_EQUAL(float(half(2) + half(2)), 4); - VERIFY_IS_EQUAL(float(half(2) + half(-2)), 0); - VERIFY_IS_APPROX(float(half(0.33333f) + half(0.66667f)), 1.0f); - VERIFY_IS_EQUAL(float(half(2.0f) * half(-5.5f)), -11.0f); - VERIFY_IS_APPROX(float(half(1.0f) / half(3.0f)), 0.33333f); - VERIFY_IS_EQUAL(float(-half(4096.0f)), -4096.0f); - VERIFY_IS_EQUAL(float(-half(-4096.0f)), 4096.0f); -} - -void test_comparison() -{ - VERIFY(half(1.0f) > half(0.5f)); - VERIFY(half(0.5f) < half(1.0f)); - VERIFY(!(half(1.0f) < half(0.5f))); - VERIFY(!(half(0.5f) > half(1.0f))); - - VERIFY(!(half(4.0f) > half(4.0f))); - VERIFY(!(half(4.0f) < half(4.0f))); - - VERIFY(!(half(0.0f) < half(-0.0f))); - VERIFY(!(half(-0.0f) < half(0.0f))); - VERIFY(!(half(0.0f) > half(-0.0f))); - VERIFY(!(half(-0.0f) > half(0.0f))); - - VERIFY(half(0.2f) > half(-1.0f)); - VERIFY(half(-1.0f) < half(0.2f)); - VERIFY(half(-16.0f) < half(-15.0f)); - - VERIFY(half(1.0f) == half(1.0f)); - VERIFY(half(1.0f) != half(2.0f)); - - // Comparisons with NaNs and infinities. -#if !EIGEN_COMP_MSVC - // Visual Studio errors out on divisions by 0 - VERIFY(!(half(0.0 / 0.0) == half(0.0 / 0.0))); - VERIFY(half(0.0 / 0.0) != half(0.0 / 0.0)); - - VERIFY(!(half(1.0) == half(0.0 / 0.0))); - VERIFY(!(half(1.0) < half(0.0 / 0.0))); - VERIFY(!(half(1.0) > half(0.0 / 0.0))); - VERIFY(half(1.0) != half(0.0 / 0.0)); - - VERIFY(half(1.0) < half(1.0 / 0.0)); - VERIFY(half(1.0) > half(-1.0 / 0.0)); -#endif -} - -void test_basic_functions() -{ - VERIFY_IS_EQUAL(float(numext::abs(half(3.5f))), 3.5f); - VERIFY_IS_EQUAL(float(abs(half(3.5f))), 3.5f); - VERIFY_IS_EQUAL(float(numext::abs(half(-3.5f))), 3.5f); - VERIFY_IS_EQUAL(float(abs(half(-3.5f))), 3.5f); - - VERIFY_IS_EQUAL(float(numext::floor(half(3.5f))), 3.0f); - VERIFY_IS_EQUAL(float(floor(half(3.5f))), 3.0f); - VERIFY_IS_EQUAL(float(numext::floor(half(-3.5f))), -4.0f); - VERIFY_IS_EQUAL(float(floor(half(-3.5f))), -4.0f); - - VERIFY_IS_EQUAL(float(numext::ceil(half(3.5f))), 4.0f); - VERIFY_IS_EQUAL(float(ceil(half(3.5f))), 4.0f); - VERIFY_IS_EQUAL(float(numext::ceil(half(-3.5f))), -3.0f); - VERIFY_IS_EQUAL(float(ceil(half(-3.5f))), -3.0f); - - VERIFY_IS_APPROX(float(numext::sqrt(half(0.0f))), 0.0f); - VERIFY_IS_APPROX(float(sqrt(half(0.0f))), 0.0f); - VERIFY_IS_APPROX(float(numext::sqrt(half(4.0f))), 2.0f); - VERIFY_IS_APPROX(float(sqrt(half(4.0f))), 2.0f); - - VERIFY_IS_APPROX(float(numext::pow(half(0.0f), half(1.0f))), 0.0f); - VERIFY_IS_APPROX(float(pow(half(0.0f), half(1.0f))), 0.0f); - VERIFY_IS_APPROX(float(numext::pow(half(2.0f), half(2.0f))), 4.0f); - VERIFY_IS_APPROX(float(pow(half(2.0f), half(2.0f))), 4.0f); - - VERIFY_IS_EQUAL(float(numext::exp(half(0.0f))), 1.0f); - VERIFY_IS_EQUAL(float(exp(half(0.0f))), 1.0f); - VERIFY_IS_APPROX(float(numext::exp(half(EIGEN_PI))), 20.f + float(EIGEN_PI)); - VERIFY_IS_APPROX(float(exp(half(EIGEN_PI))), 20.f + float(EIGEN_PI)); - - VERIFY_IS_EQUAL(float(numext::log(half(1.0f))), 0.0f); - VERIFY_IS_EQUAL(float(log(half(1.0f))), 0.0f); - VERIFY_IS_APPROX(float(numext::log(half(10.0f))), 2.30273f); - VERIFY_IS_APPROX(float(log(half(10.0f))), 2.30273f); -} - -void test_trigonometric_functions() -{ - VERIFY_IS_APPROX(numext::cos(half(0.0f)), half(cosf(0.0f))); - VERIFY_IS_APPROX(cos(half(0.0f)), half(cosf(0.0f))); - VERIFY_IS_APPROX(numext::cos(half(EIGEN_PI)), half(cosf(EIGEN_PI))); - //VERIFY_IS_APPROX(numext::cos(half(EIGEN_PI/2)), half(cosf(EIGEN_PI/2))); - //VERIFY_IS_APPROX(numext::cos(half(3*EIGEN_PI/2)), half(cosf(3*EIGEN_PI/2))); - VERIFY_IS_APPROX(numext::cos(half(3.5f)), half(cosf(3.5f))); - - VERIFY_IS_APPROX(numext::sin(half(0.0f)), half(sinf(0.0f))); - VERIFY_IS_APPROX(sin(half(0.0f)), half(sinf(0.0f))); - // VERIFY_IS_APPROX(numext::sin(half(EIGEN_PI)), half(sinf(EIGEN_PI))); - VERIFY_IS_APPROX(numext::sin(half(EIGEN_PI/2)), half(sinf(EIGEN_PI/2))); - VERIFY_IS_APPROX(numext::sin(half(3*EIGEN_PI/2)), half(sinf(3*EIGEN_PI/2))); - VERIFY_IS_APPROX(numext::sin(half(3.5f)), half(sinf(3.5f))); - - VERIFY_IS_APPROX(numext::tan(half(0.0f)), half(tanf(0.0f))); - VERIFY_IS_APPROX(tan(half(0.0f)), half(tanf(0.0f))); - // VERIFY_IS_APPROX(numext::tan(half(EIGEN_PI)), half(tanf(EIGEN_PI))); - // VERIFY_IS_APPROX(numext::tan(half(EIGEN_PI/2)), half(tanf(EIGEN_PI/2))); - //VERIFY_IS_APPROX(numext::tan(half(3*EIGEN_PI/2)), half(tanf(3*EIGEN_PI/2))); - VERIFY_IS_APPROX(numext::tan(half(3.5f)), half(tanf(3.5f))); -} - -void test_array() -{ - typedef Array ArrayXh; - Index size = internal::random(1,10); - Index i = internal::random(0,size-1); - ArrayXh a1 = ArrayXh::Random(size), a2 = ArrayXh::Random(size); - VERIFY_IS_APPROX( a1+a1, half(2)*a1 ); - VERIFY( (a1.abs() >= half(0)).all() ); - VERIFY_IS_APPROX( (a1*a1).sqrt(), a1.abs() ); - - VERIFY( ((a1.min)(a2) <= (a1.max)(a2)).all() ); - a1(i) = half(-10.); - VERIFY_IS_EQUAL( a1.minCoeff(), half(-10.) ); - a1(i) = half(10.); - VERIFY_IS_EQUAL( a1.maxCoeff(), half(10.) ); - - std::stringstream ss; - ss << a1; -} - -void test_cxx11_float16() -{ - CALL_SUBTEST(test_conversion()); - CALL_SUBTEST(test_numtraits()); - CALL_SUBTEST(test_arithmetic()); - CALL_SUBTEST(test_comparison()); - CALL_SUBTEST(test_basic_functions()); - CALL_SUBTEST(test_trigonometric_functions()); - CALL_SUBTEST(test_array()); -} -- cgit v1.2.3 From 9908020d36952b8d82d561bd2b018dd69ef5042e Mon Sep 17 00:00:00 2001 From: Gael Guennebaud Date: Mon, 25 Jul 2016 14:25:56 +0200 Subject: Add minimal support for Array, and fix Tensor --- Eigen/src/Core/NumTraits.h | 21 ++++++++++++++++++ test/CMakeLists.txt | 1 + test/array_of_string.cpp | 32 +++++++++++++++++++++++++++ unsupported/Eigen/CXX11/src/Tensor/TensorIO.h | 5 ----- 4 files changed, 54 insertions(+), 5 deletions(-) create mode 100644 test/array_of_string.cpp (limited to 'test/CMakeLists.txt') diff --git a/Eigen/src/Core/NumTraits.h b/Eigen/src/Core/NumTraits.h index 53821e8dc..42cffbd3b 100644 --- a/Eigen/src/Core/NumTraits.h +++ b/Eigen/src/Core/NumTraits.h @@ -234,6 +234,27 @@ struct NumTraits > static inline RealScalar dummy_precision() { return NumTraits::dummy_precision(); } }; +template<> struct NumTraits + : GenericNumTraits +{ + enum { + RequireInitialization = 1, + ReadCost = HugeCost, + AddCost = HugeCost, + MulCost = HugeCost + }; + + static inline int digits10() { return 0; } + +private: + static inline std::string epsilon(); + static inline std::string dummy_precision(); + static inline std::string lowest(); + static inline std::string highest(); + static inline std::string infinity(); + static inline std::string quiet_NaN(); +}; + } // end namespace Eigen #endif // EIGEN_NUMTRAITS_H diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index e42e850ca..926b284e6 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -260,6 +260,7 @@ ei_add_test(ctorleak) ei_add_test(mpl2only) ei_add_test(inplace_decomposition) ei_add_test(half_float) +ei_add_test(array_of_string) add_executable(bug1213 bug1213.cpp bug1213_main.cpp) diff --git a/test/array_of_string.cpp b/test/array_of_string.cpp new file mode 100644 index 000000000..e23b7c59e --- /dev/null +++ b/test/array_of_string.cpp @@ -0,0 +1,32 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2016 Gael Guennebaud +// +// This Source Code Form is subject to the terms of the Mozilla +// Public License v. 2.0. If a copy of the MPL was not distributed +// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. + +#include "main.h" + +void test_array_of_string() +{ + typedef Array ArrayXs; + ArrayXs a1(3), a2(3), a3(3), a3ref(3); + a1 << "one", "two", "three"; + a2 << "1", "2", "3"; + a3ref << "one (1)", "two (2)", "three (3)"; + std::stringstream s1; + s1 << a1; + VERIFY_IS_EQUAL(s1.str(), std::string(" one two three")); + a3 = a1 + std::string(" (") + a2 + std::string(")"); + VERIFY((a3==a3ref).all()); + + a3 = a1; + a3 += std::string(" (") + a2 + std::string(")"); + VERIFY((a3==a3ref).all()); + + a1.swap(a3); + VERIFY((a1==a3ref).all()); + VERIFY((a3!=a3ref).all()); +} diff --git a/unsupported/Eigen/CXX11/src/Tensor/TensorIO.h b/unsupported/Eigen/CXX11/src/Tensor/TensorIO.h index f3a3a1b88..a901c5dd4 100644 --- a/unsupported/Eigen/CXX11/src/Tensor/TensorIO.h +++ b/unsupported/Eigen/CXX11/src/Tensor/TensorIO.h @@ -13,11 +13,6 @@ namespace Eigen { namespace internal { -template<> -struct significant_decimals_impl - : significant_decimals_default_impl -{}; - // Print the tensor as a 2d matrix template -- cgit v1.2.3 From 616a7a19129dfe1b64a4e578cd35e84861d8046b Mon Sep 17 00:00:00 2001 From: Benoit Steiner Date: Mon, 3 Oct 2016 17:09:33 -0700 Subject: Improved support for compiling CUDA code with clang as the host compiler --- test/CMakeLists.txt | 2 +- unsupported/test/CMakeLists.txt | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) (limited to 'test/CMakeLists.txt') diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 926b284e6..e17985107 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -355,7 +355,7 @@ if(CUDA_FOUND) set(CUDA_PROPAGATE_HOST_FLAGS OFF) if("${CMAKE_CXX_COMPILER_ID}" STREQUAL "Clang") - set(CUDA_NVCC_FLAGS "-ccbin /usr/bin/clang" CACHE STRING "nvcc flags" FORCE) + set(CUDA_NVCC_FLAGS "-ccbin ${CMAKE_C_COMPILER}" CACHE STRING "nvcc flags" FORCE) endif() if(EIGEN_TEST_CUDA_CLANG) set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11 --cuda-gpu-arch=sm_30") diff --git a/unsupported/test/CMakeLists.txt b/unsupported/test/CMakeLists.txt index 9eac6ec73..a1823beaa 100644 --- a/unsupported/test/CMakeLists.txt +++ b/unsupported/test/CMakeLists.txt @@ -203,7 +203,7 @@ if(CUDA_FOUND AND EIGEN_TEST_CUDA) message(STATUS "Flags used to compile cuda code: " ${CMAKE_CXX_FLAGS}) if("${CMAKE_CXX_COMPILER_ID}" STREQUAL "Clang") - set(CUDA_NVCC_FLAGS "-ccbin /usr/bin/clang" CACHE STRING "nvcc flags" FORCE) + set(CUDA_NVCC_FLAGS "-ccbin ${CMAKE_C_COMPILER}" CACHE STRING "nvcc flags" FORCE) endif() if(EIGEN_TEST_CUDA_CLANG) set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11 --cuda-gpu-arch=sm_${EIGEN_CUDA_COMPUTE_ARCH}") -- cgit v1.2.3