Bugzilla – Attachment 556 Details for
Bug 949
Add static assertion for incompatible scalar types in decompositions
Home
|
New
|
Browse
|
Search
|
[?]
|
Reports
|
Requests
|
Help
|
Log In
[x]
|
Forgot Password
Login:
[x]
This bugzilla service is closed. All entries have been migrated to
https://gitlab.com/libeigen/eigen
[patch]
Static assertions on wrong scalar types
bug_949_static_assert_on_scalar_type_for_solvers.diff (text/plain), 22.05 KB, created by
Gael Guennebaud
on 2015-03-09 21:36:06 UTC
(
hide
)
Description:
Static assertions on wrong scalar types
Filename:
MIME Type:
Creator:
Gael Guennebaud
Created:
2015-03-09 21:36:06 UTC
Size:
22.05 KB
patch
obsolete
>diff --git a/Eigen/src/Cholesky/LDLT.h b/Eigen/src/Cholesky/LDLT.h >--- a/Eigen/src/Cholesky/LDLT.h >+++ b/Eigen/src/Cholesky/LDLT.h >@@ -221,16 +221,21 @@ template<typename _MatrixType, int _UpLo > > #ifndef EIGEN_PARSED_BY_DOXYGEN > template<typename RhsType, typename DstType> > EIGEN_DEVICE_FUNC > void _solve_impl(const RhsType &rhs, DstType &dst) const; > #endif > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ } > > /** \internal > * Used to compute and store the Cholesky decomposition A = L D L^* = U^* D U. > * The strict upper part is used during the decomposition, the strict lower > * part correspond to the coefficients of L (its diagonal is equal to 1 and > * is not stored), and the diagonal entries correspond to D. > */ > MatrixType m_matrix; >@@ -419,16 +424,18 @@ template<typename MatrixType> struct LDL > > } // end namespace internal > > /** Compute / recompute the LDLT decomposition A = L D L^* = U^* D U of \a matrix > */ > template<typename MatrixType, int _UpLo> > LDLT<MatrixType,_UpLo>& LDLT<MatrixType,_UpLo>::compute(const MatrixType& a) > { >+ check_template_parameters(); >+ > eigen_assert(a.rows()==a.cols()); > const Index size = a.rows(); > > m_matrix = a; > > m_transpositions.resize(size); > m_isInitialized = false; > m_temporary.resize(size); >diff --git a/Eigen/src/Cholesky/LLT.h b/Eigen/src/Cholesky/LLT.h >--- a/Eigen/src/Cholesky/LLT.h >+++ b/Eigen/src/Cholesky/LLT.h >@@ -165,16 +165,22 @@ template<typename _MatrixType, int _UpLo > > #ifndef EIGEN_PARSED_BY_DOXYGEN > template<typename RhsType, typename DstType> > EIGEN_DEVICE_FUNC > void _solve_impl(const RhsType &rhs, DstType &dst) const; > #endif > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ } >+ > /** \internal > * Used to compute and store L > * The strict upper part is not used and even not initialized. > */ > MatrixType m_matrix; > bool m_isInitialized; > ComputationInfo m_info; > }; >@@ -372,16 +378,18 @@ template<typename MatrixType> struct LLT > * \returns a reference to *this > * > * Example: \include TutorialLinAlgComputeTwice.cpp > * Output: \verbinclude TutorialLinAlgComputeTwice.out > */ > template<typename MatrixType, int _UpLo> > LLT<MatrixType,_UpLo>& LLT<MatrixType,_UpLo>::compute(const MatrixType& a) > { >+ check_template_parameters(); >+ > eigen_assert(a.rows()==a.cols()); > const Index size = a.rows(); > m_matrix.resize(size, size); > m_matrix = a; > > m_isInitialized = true; > bool ok = Traits::inplace_decomposition(m_matrix); > m_info = ok ? Success : NumericalIssue; >diff --git a/Eigen/src/Eigenvalues/ComplexEigenSolver.h b/Eigen/src/Eigenvalues/ComplexEigenSolver.h >--- a/Eigen/src/Eigenvalues/ComplexEigenSolver.h >+++ b/Eigen/src/Eigenvalues/ComplexEigenSolver.h >@@ -229,16 +229,22 @@ template<typename _MatrixType> class Com > > /** \brief Returns the maximum number of iterations. */ > Index getMaxIterations() > { > return m_schur.getMaxIterations(); > } > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ } >+ > EigenvectorType m_eivec; > EigenvalueType m_eivalues; > ComplexSchur<MatrixType> m_schur; > bool m_isInitialized; > bool m_eigenvectorsOk; > EigenvectorType m_matX; > > private: >@@ -246,16 +252,18 @@ template<typename _MatrixType> class Com > void sortEigenvalues(bool computeEigenvectors); > }; > > > template<typename MatrixType> > ComplexEigenSolver<MatrixType>& > ComplexEigenSolver<MatrixType>::compute(const MatrixType& matrix, bool computeEigenvectors) > { >+ check_template_parameters(); >+ > // this code is inspired from Jampack > eigen_assert(matrix.cols() == matrix.rows()); > > // Do a complex Schur decomposition, A = U T U^* > // The eigenvalues are on the diagonal of T. > m_schur.compute(matrix, computeEigenvectors); > > if(m_schur.info() == Success) >diff --git a/Eigen/src/Eigenvalues/EigenSolver.h b/Eigen/src/Eigenvalues/EigenSolver.h >--- a/Eigen/src/Eigenvalues/EigenSolver.h >+++ b/Eigen/src/Eigenvalues/EigenSolver.h >@@ -294,16 +294,23 @@ template<typename _MatrixType> class Eig > { > return m_realSchur.getMaxIterations(); > } > > private: > void doComputeEigenvectors(); > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ EIGEN_STATIC_ASSERT(!NumTraits<Scalar>::IsComplex, NUMERIC_TYPE_MUST_BE_REAL); >+ } >+ > MatrixType m_eivec; > EigenvalueType m_eivalues; > bool m_isInitialized; > bool m_eigenvectorsOk; > ComputationInfo m_info; > RealSchur<MatrixType> m_realSchur; > MatrixType m_matT; > >@@ -361,16 +368,18 @@ typename EigenSolver<MatrixType>::Eigenv > } > return matV; > } > > template<typename MatrixType> > EigenSolver<MatrixType>& > EigenSolver<MatrixType>::compute(const MatrixType& matrix, bool computeEigenvectors) > { >+ check_template_parameters(); >+ > using std::sqrt; > using std::abs; > using numext::isfinite; > eigen_assert(matrix.cols() == matrix.rows()); > > // Reduce to real Schur form. > m_realSchur.compute(matrix, computeEigenvectors); > >diff --git a/Eigen/src/Eigenvalues/GeneralizedEigenSolver.h b/Eigen/src/Eigenvalues/GeneralizedEigenSolver.h >--- a/Eigen/src/Eigenvalues/GeneralizedEigenSolver.h >+++ b/Eigen/src/Eigenvalues/GeneralizedEigenSolver.h >@@ -258,16 +258,23 @@ template<typename _MatrixType> class Gen > */ > GeneralizedEigenSolver& setMaxIterations(Index maxIters) > { > m_realQZ.setMaxIterations(maxIters); > return *this; > } > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ EIGEN_STATIC_ASSERT(!NumTraits<Scalar>::IsComplex, NUMERIC_TYPE_MUST_BE_REAL); >+ } >+ > MatrixType m_eivec; > ComplexVectorType m_alphas; > VectorType m_betas; > bool m_isInitialized; > bool m_eigenvectorsOk; > RealQZ<MatrixType> m_realQZ; > MatrixType m_matS; > >@@ -285,16 +292,18 @@ template<typename _MatrixType> class Gen > // // TODO > // return matV; > //} > > template<typename MatrixType> > GeneralizedEigenSolver<MatrixType>& > GeneralizedEigenSolver<MatrixType>::compute(const MatrixType& A, const MatrixType& B, bool computeEigenvectors) > { >+ check_template_parameters(); >+ > using std::sqrt; > using std::abs; > eigen_assert(A.cols() == A.rows() && B.cols() == A.rows() && B.cols() == B.rows()); > > // Reduce to generalized real Schur form: > // A = Q S Z and B = Q T Z > m_realQZ.compute(A, B, computeEigenvectors); > >diff --git a/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h b/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h >--- a/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h >+++ b/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h >@@ -342,16 +342,21 @@ template<typename _MatrixType> class Sel > /** \brief Maximum number of iterations. > * > * The algorithm terminates if it does not converge within m_maxIterations * n iterations, where n > * denotes the size of the matrix. This value is currently set to 30 (copied from LAPACK). > */ > static const int m_maxIterations = 30; > > protected: >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ } >+ > MatrixType m_eivec; > RealVectorType m_eivalues; > typename TridiagonalizationType::SubDiagonalType m_subdiag; > ComputationInfo m_info; > bool m_isInitialized; > bool m_eigenvectorsOk; > }; > >@@ -377,16 +382,18 @@ EIGEN_DEVICE_FUNC > static void tridiagonal_qr_step(RealScalar* diag, RealScalar* subdiag, Index start, Index end, Scalar* matrixQ, Index n); > } > > template<typename MatrixType> > EIGEN_DEVICE_FUNC > SelfAdjointEigenSolver<MatrixType>& SelfAdjointEigenSolver<MatrixType> > ::compute(const MatrixType& matrix, int options) > { >+ check_template_parameters(); >+ > using std::abs; > eigen_assert(matrix.cols() == matrix.rows()); > eigen_assert((options&~(EigVecMask|GenEigMask))==0 > && (options&EigVecMask)!=EigVecMask > && "invalid option parameter"); > bool computeEigenvectors = (options&ComputeEigenvectors)==ComputeEigenvectors; > Index n = matrix.cols(); > m_eivalues.resize(n,1); >diff --git a/Eigen/src/LU/FullPivLU.h b/Eigen/src/LU/FullPivLU.h >--- a/Eigen/src/LU/FullPivLU.h >+++ b/Eigen/src/LU/FullPivLU.h >@@ -385,16 +385,22 @@ template<typename _MatrixType> class Ful > > #ifndef EIGEN_PARSED_BY_DOXYGEN > template<typename RhsType, typename DstType> > EIGEN_DEVICE_FUNC > void _solve_impl(const RhsType &rhs, DstType &dst) const; > #endif > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ } >+ > MatrixType m_lu; > PermutationPType m_p; > PermutationQType m_q; > IntColVectorType m_rowsTranspositions; > IntRowVectorType m_colsTranspositions; > Index m_det_pq, m_nonzero_pivots; > RealScalar m_maxpivot, m_prescribedThreshold; > bool m_isInitialized, m_usePrescribedThreshold; >@@ -429,16 +435,18 @@ FullPivLU<MatrixType>::FullPivLU(const M > m_usePrescribedThreshold(false) > { > compute(matrix); > } > > template<typename MatrixType> > FullPivLU<MatrixType>& FullPivLU<MatrixType>::compute(const MatrixType& matrix) > { >+ check_template_parameters(); >+ > // the permutations are stored as int indices, so just to be sure: > eigen_assert(matrix.rows()<=NumTraits<int>::highest() && matrix.cols()<=NumTraits<int>::highest()); > > m_isInitialized = true; > m_lu = matrix; > > const Index size = matrix.diagonalSize(); > const Index rows = matrix.rows(); >diff --git a/Eigen/src/LU/PartialPivLU.h b/Eigen/src/LU/PartialPivLU.h >--- a/Eigen/src/LU/PartialPivLU.h >+++ b/Eigen/src/LU/PartialPivLU.h >@@ -204,16 +204,22 @@ template<typename _MatrixType> class Par > m_lu.template triangularView<UnitLower>().solveInPlace(dst); > > // Step 3 > m_lu.template triangularView<Upper>().solveInPlace(dst); > } > #endif > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ } >+ > MatrixType m_lu; > PermutationType m_p; > TranspositionType m_rowsTranspositions; > Index m_det_p; > bool m_isInitialized; > }; > > template<typename MatrixType> >@@ -420,16 +426,18 @@ void partial_lu_inplace(MatrixType& lu, > ::blocked_lu(lu.rows(), lu.cols(), &lu.coeffRef(0,0), lu.outerStride(), &row_transpositions.coeffRef(0), nb_transpositions); > } > > } // end namespace internal > > template<typename MatrixType> > PartialPivLU<MatrixType>& PartialPivLU<MatrixType>::compute(const MatrixType& matrix) > { >+ check_template_parameters(); >+ > // the row permutation is stored as int indices, so just to be sure: > eigen_assert(matrix.rows()<NumTraits<int>::highest()); > > m_lu = matrix; > > eigen_assert(matrix.rows() == matrix.cols() && "PartialPivLU is only for square (and moreover invertible) matrices"); > const Index size = matrix.rows(); > >diff --git a/Eigen/src/QR/ColPivHouseholderQR.h b/Eigen/src/QR/ColPivHouseholderQR.h >--- a/Eigen/src/QR/ColPivHouseholderQR.h >+++ b/Eigen/src/QR/ColPivHouseholderQR.h >@@ -393,16 +393,22 @@ template<typename _MatrixType> class Col > > #ifndef EIGEN_PARSED_BY_DOXYGEN > template<typename RhsType, typename DstType> > EIGEN_DEVICE_FUNC > void _solve_impl(const RhsType &rhs, DstType &dst) const; > #endif > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ } >+ > MatrixType m_qr; > HCoeffsType m_hCoeffs; > PermutationType m_colsPermutation; > IntRowVectorType m_colsTranspositions; > RowVectorType m_temp; > RealRowVectorType m_colSqNorms; > bool m_isInitialized, m_usePrescribedThreshold; > RealScalar m_prescribedThreshold, m_maxpivot; >@@ -431,16 +437,18 @@ typename MatrixType::RealScalar ColPivHo > * the factorization is stored into \c *this, and a reference to \c *this > * is returned. > * > * \sa class ColPivHouseholderQR, ColPivHouseholderQR(const MatrixType&) > */ > template<typename MatrixType> > ColPivHouseholderQR<MatrixType>& ColPivHouseholderQR<MatrixType>::compute(const MatrixType& matrix) > { >+ check_template_parameters(); >+ > using std::abs; > Index rows = matrix.rows(); > Index cols = matrix.cols(); > Index size = matrix.diagonalSize(); > > // the column permutation is stored as int indices, so just to be sure: > eigen_assert(cols<=NumTraits<int>::highest()); > >diff --git a/Eigen/src/QR/FullPivHouseholderQR.h b/Eigen/src/QR/FullPivHouseholderQR.h >--- a/Eigen/src/QR/FullPivHouseholderQR.h >+++ b/Eigen/src/QR/FullPivHouseholderQR.h >@@ -375,16 +375,22 @@ template<typename _MatrixType> class Ful > > #ifndef EIGEN_PARSED_BY_DOXYGEN > template<typename RhsType, typename DstType> > EIGEN_DEVICE_FUNC > void _solve_impl(const RhsType &rhs, DstType &dst) const; > #endif > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ } >+ > MatrixType m_qr; > HCoeffsType m_hCoeffs; > IntDiagSizeVectorType m_rows_transpositions; > IntDiagSizeVectorType m_cols_transpositions; > PermutationType m_cols_permutation; > RowVectorType m_temp; > bool m_isInitialized, m_usePrescribedThreshold; > RealScalar m_prescribedThreshold, m_maxpivot; >@@ -414,16 +420,18 @@ typename MatrixType::RealScalar FullPivH > * the factorization is stored into \c *this, and a reference to \c *this > * is returned. > * > * \sa class FullPivHouseholderQR, FullPivHouseholderQR(const MatrixType&) > */ > template<typename MatrixType> > FullPivHouseholderQR<MatrixType>& FullPivHouseholderQR<MatrixType>::compute(const MatrixType& matrix) > { >+ check_template_parameters(); >+ > using std::abs; > Index rows = matrix.rows(); > Index cols = matrix.cols(); > Index size = (std::min)(rows,cols); > > m_qr = matrix; > m_hCoeffs.resize(size); > >diff --git a/Eigen/src/QR/HouseholderQR.h b/Eigen/src/QR/HouseholderQR.h >--- a/Eigen/src/QR/HouseholderQR.h >+++ b/Eigen/src/QR/HouseholderQR.h >@@ -191,16 +191,22 @@ template<typename _MatrixType> class Hou > > #ifndef EIGEN_PARSED_BY_DOXYGEN > template<typename RhsType, typename DstType> > EIGEN_DEVICE_FUNC > void _solve_impl(const RhsType &rhs, DstType &dst) const; > #endif > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ } >+ > MatrixType m_qr; > HCoeffsType m_hCoeffs; > RowVectorType m_temp; > bool m_isInitialized; > }; > > template<typename MatrixType> > typename MatrixType::RealScalar HouseholderQR<MatrixType>::absDeterminant() const >@@ -343,16 +349,18 @@ void HouseholderQR<_MatrixType>::_solve_ > * the factorization is stored into \c *this, and a reference to \c *this > * is returned. > * > * \sa class HouseholderQR, HouseholderQR(const MatrixType&) > */ > template<typename MatrixType> > HouseholderQR<MatrixType>& HouseholderQR<MatrixType>::compute(const MatrixType& matrix) > { >+ check_template_parameters(); >+ > Index rows = matrix.rows(); > Index cols = matrix.cols(); > Index size = (std::min)(rows,cols); > > m_qr = matrix; > m_hCoeffs.resize(size); > > m_temp.resize(cols); >diff --git a/Eigen/src/SVD/SVDBase.h b/Eigen/src/SVD/SVDBase.h >--- a/Eigen/src/SVD/SVDBase.h >+++ b/Eigen/src/SVD/SVDBase.h >@@ -212,16 +212,22 @@ public: > > #ifndef EIGEN_PARSED_BY_DOXYGEN > template<typename RhsType, typename DstType> > EIGEN_DEVICE_FUNC > void _solve_impl(const RhsType &rhs, DstType &dst) const; > #endif > > protected: >+ >+ static void check_template_parameters() >+ { >+ EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar); >+ } >+ > // return true if already allocated > bool allocate(Index rows, Index cols, unsigned int computationOptions) ; > > MatrixUType m_matrixU; > MatrixVType m_matrixV; > SingularValuesType m_singularValues; > bool m_isInitialized, m_isAllocated, m_usePrescribedThreshold; > bool m_computeFullU, m_computeThinU; >@@ -235,17 +241,19 @@ protected: > * Default constructor of SVDBase > */ > SVDBase() > : m_isInitialized(false), > m_isAllocated(false), > m_usePrescribedThreshold(false), > m_computationOptions(0), > m_rows(-1), m_cols(-1), m_diagSize(0) >- {} >+ { >+ check_template_parameters(); >+ } > > > }; > > #ifndef EIGEN_PARSED_BY_DOXYGEN > template<typename Derived> > template<typename RhsType, typename DstType> > void SVDBase<Derived>::_solve_impl(const RhsType &rhs, DstType &dst) const >diff --git a/failtest/CMakeLists.txt b/failtest/CMakeLists.txt >--- a/failtest/CMakeLists.txt >+++ b/failtest/CMakeLists.txt >@@ -42,16 +42,28 @@ ei_add_failtest("swap_1") > ei_add_failtest("swap_2") > > ei_add_failtest("sparse_ref_1") > ei_add_failtest("sparse_ref_2") > ei_add_failtest("sparse_ref_3") > ei_add_failtest("sparse_ref_4") > ei_add_failtest("sparse_ref_5") > >+ei_add_failtest("partialpivlu_int") >+ei_add_failtest("fullpivlu_int") >+ei_add_failtest("llt_int") >+ei_add_failtest("ldlt_int") >+ei_add_failtest("qr_int") >+ei_add_failtest("colpivqr_int") >+ei_add_failtest("fullpivqr_int") >+ei_add_failtest("jacobisvd_int") >+ei_add_failtest("bdcsvd_int") >+ei_add_failtest("eigensolver_int") >+ei_add_failtest("eigensolver_cplx") >+ > if (EIGEN_FAILTEST_FAILURE_COUNT) > message(FATAL_ERROR > "${EIGEN_FAILTEST_FAILURE_COUNT} out of ${EIGEN_FAILTEST_COUNT} failtests FAILED. " > "To debug these failures, manually compile these programs in ${CMAKE_CURRENT_SOURCE_DIR}, " > "with and without #define EIGEN_SHOULD_FAIL_TO_BUILD.") > else() > message(STATUS "Failtest SUCCESS: all ${EIGEN_FAILTEST_COUNT} failtests passed.") > message(STATUS "") >diff --git a/failtest/bdcsvd_int.cpp b/failtest/bdcsvd_int.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/bdcsvd_int.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/SVD" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR int >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ BDCSVD<Matrix<SCALAR,Dynamic,Dynamic> > qr(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+} >diff --git a/failtest/colpivqr_int.cpp b/failtest/colpivqr_int.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/colpivqr_int.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/QR" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR int >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ ColPivHouseholderQR<Matrix<SCALAR,Dynamic,Dynamic> > qr(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+} >diff --git a/failtest/eigensolver_cplx.cpp b/failtest/eigensolver_cplx.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/eigensolver_cplx.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/Eigenvalues" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR std::complex<double> >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ EigenSolver<Matrix<SCALAR,Dynamic,Dynamic> > eig(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+} >diff --git a/failtest/eigensolver_int.cpp b/failtest/eigensolver_int.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/eigensolver_int.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/Eigenvalues" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR int >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ EigenSolver<Matrix<SCALAR,Dynamic,Dynamic> > eig(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+} >diff --git a/failtest/fullpivlu_int.cpp b/failtest/fullpivlu_int.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/fullpivlu_int.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/LU" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR int >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ FullPivLU<Matrix<SCALAR,Dynamic,Dynamic> > lu(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+} >diff --git a/failtest/fullpivqr_int.cpp b/failtest/fullpivqr_int.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/fullpivqr_int.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/QR" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR int >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ FullPivHouseholderQR<Matrix<SCALAR,Dynamic,Dynamic> > qr(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+} >diff --git a/failtest/jacobisvd_int.cpp b/failtest/jacobisvd_int.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/jacobisvd_int.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/SVD" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR int >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ JacobiSVD<Matrix<SCALAR,Dynamic,Dynamic> > qr(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+} >diff --git a/failtest/ldlt_int.cpp b/failtest/ldlt_int.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/ldlt_int.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/Cholesky" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR int >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ LDLT<Matrix<SCALAR,Dynamic,Dynamic> > ldlt(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+} >diff --git a/failtest/llt_int.cpp b/failtest/llt_int.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/llt_int.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/Cholesky" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR int >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ LLT<Matrix<SCALAR,Dynamic,Dynamic> > llt(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+} >diff --git a/failtest/partialpivlu_int.cpp b/failtest/partialpivlu_int.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/partialpivlu_int.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/LU" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR int >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ PartialPivLU<Matrix<SCALAR,Dynamic,Dynamic> > lu(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+} >diff --git a/failtest/qr_int.cpp b/failtest/qr_int.cpp >new file mode 100644 >--- /dev/null >+++ b/failtest/qr_int.cpp >@@ -0,0 +1,14 @@ >+#include "../Eigen/QR" >+ >+#ifdef EIGEN_SHOULD_FAIL_TO_BUILD >+#define SCALAR int >+#else >+#define SCALAR float >+#endif >+ >+using namespace Eigen; >+ >+int main() >+{ >+ HouseholderQR<Matrix<SCALAR,Dynamic,Dynamic> > qr(Matrix<SCALAR,Dynamic,Dynamic>::Random(10,10)); >+}
You cannot view the attachment while viewing its details because your browser does not support IFRAMEs.
View the attachment on a separate page
.
View Attachment As Diff
View Attachment As Raw
Actions:
View
|
Diff
Attachments on
bug 949
:
537
|
555
| 556