This commit is contained in:
Thomas Capricelli 2009-12-01 21:07:20 +01:00
commit 291fd4f8da
40 changed files with 292 additions and 202 deletions

View File

@ -152,6 +152,9 @@ add_subdirectory(doc EXCLUDE_FROM_ALL)
include(CTest) include(CTest)
enable_testing() # must be called from the root CMakeLists, see man page enable_testing() # must be called from the root CMakeLists, see man page
include(EigenTesting)
ei_init_testing()
if(EIGEN_LEAVE_TEST_IN_ALL_TARGET) if(EIGEN_LEAVE_TEST_IN_ALL_TARGET)
add_subdirectory(test) # can't do EXCLUDE_FROM_ALL here, breaks CTest add_subdirectory(test) # can't do EXCLUDE_FROM_ALL here, breaks CTest
else() else()
@ -164,6 +167,7 @@ add_subdirectory(demos EXCLUDE_FROM_ALL)
add_subdirectory(blas EXCLUDE_FROM_ALL) add_subdirectory(blas EXCLUDE_FROM_ALL)
# must be after test and unsupported, for configuring buildtests.in
add_subdirectory(scripts EXCLUDE_FROM_ALL) add_subdirectory(scripts EXCLUDE_FROM_ALL)
# TODO: consider also replacing EIGEN_BUILD_BTL by a custom target "make btl"? # TODO: consider also replacing EIGEN_BUILD_BTL by a custom target "make btl"?

View File

@ -30,7 +30,7 @@ template<typename T> inline typename NumTraits<T>::Real epsilon()
return std::numeric_limits<typename NumTraits<T>::Real>::epsilon(); return std::numeric_limits<typename NumTraits<T>::Real>::epsilon();
} }
template<typename T> inline typename NumTraits<T>::Real precision(); template<typename T> inline typename NumTraits<T>::Real dummy_precision();
template<typename T> inline T ei_random(T a, T b); template<typename T> inline T ei_random(T a, T b);
template<typename T> inline T ei_random(); template<typename T> inline T ei_random();
@ -55,7 +55,7 @@ template<typename T> inline typename NumTraits<T>::Real ei_hypot(T x, T y)
*** int *** *** int ***
**************/ **************/
template<> inline int precision<int>() { return 0; } template<> inline int dummy_precision<int>() { return 0; }
inline int ei_real(int x) { return x; } inline int ei_real(int x) { return x; }
inline int& ei_real_ref(int& x) { return x; } inline int& ei_real_ref(int& x) { return x; }
inline int ei_imag(int) { return 0; } inline int ei_imag(int) { return 0; }
@ -92,15 +92,15 @@ template<> inline int ei_random()
{ {
return ei_random<int>(-ei_random_amplitude<int>(), ei_random_amplitude<int>()); return ei_random<int>(-ei_random_amplitude<int>(), ei_random_amplitude<int>());
} }
inline bool ei_isMuchSmallerThan(int a, int, int = precision<int>()) inline bool ei_isMuchSmallerThan(int a, int, int = dummy_precision<int>())
{ {
return a == 0; return a == 0;
} }
inline bool ei_isApprox(int a, int b, int = precision<int>()) inline bool ei_isApprox(int a, int b, int = dummy_precision<int>())
{ {
return a == b; return a == b;
} }
inline bool ei_isApproxOrLessThan(int a, int b, int = precision<int>()) inline bool ei_isApproxOrLessThan(int a, int b, int = dummy_precision<int>())
{ {
return a <= b; return a <= b;
} }
@ -109,7 +109,7 @@ inline bool ei_isApproxOrLessThan(int a, int b, int = precision<int>())
*** float *** *** float ***
**************/ **************/
template<> inline float precision<float>() { return 1e-5f; } template<> inline float dummy_precision<float>() { return 1e-5f; }
inline float ei_real(float x) { return x; } inline float ei_real(float x) { return x; }
inline float& ei_real_ref(float& x) { return x; } inline float& ei_real_ref(float& x) { return x; }
inline float ei_imag(float) { return 0.f; } inline float ei_imag(float) { return 0.f; }
@ -140,15 +140,15 @@ template<> inline float ei_random()
{ {
return ei_random<float>(-ei_random_amplitude<float>(), ei_random_amplitude<float>()); return ei_random<float>(-ei_random_amplitude<float>(), ei_random_amplitude<float>());
} }
inline bool ei_isMuchSmallerThan(float a, float b, float prec = precision<float>()) inline bool ei_isMuchSmallerThan(float a, float b, float prec = dummy_precision<float>())
{ {
return ei_abs(a) <= ei_abs(b) * prec; return ei_abs(a) <= ei_abs(b) * prec;
} }
inline bool ei_isApprox(float a, float b, float prec = precision<float>()) inline bool ei_isApprox(float a, float b, float prec = dummy_precision<float>())
{ {
return ei_abs(a - b) <= std::min(ei_abs(a), ei_abs(b)) * prec; return ei_abs(a - b) <= std::min(ei_abs(a), ei_abs(b)) * prec;
} }
inline bool ei_isApproxOrLessThan(float a, float b, float prec = precision<float>()) inline bool ei_isApproxOrLessThan(float a, float b, float prec = dummy_precision<float>())
{ {
return a <= b || ei_isApprox(a, b, prec); return a <= b || ei_isApprox(a, b, prec);
} }
@ -157,7 +157,7 @@ inline bool ei_isApproxOrLessThan(float a, float b, float prec = precision<float
*** double *** *** double ***
**************/ **************/
template<> inline double precision<double>() { return 1e-12; } template<> inline double dummy_precision<double>() { return 1e-12; }
inline double ei_real(double x) { return x; } inline double ei_real(double x) { return x; }
inline double& ei_real_ref(double& x) { return x; } inline double& ei_real_ref(double& x) { return x; }
@ -189,15 +189,15 @@ template<> inline double ei_random()
{ {
return ei_random<double>(-ei_random_amplitude<double>(), ei_random_amplitude<double>()); return ei_random<double>(-ei_random_amplitude<double>(), ei_random_amplitude<double>());
} }
inline bool ei_isMuchSmallerThan(double a, double b, double prec = precision<double>()) inline bool ei_isMuchSmallerThan(double a, double b, double prec = dummy_precision<double>())
{ {
return ei_abs(a) <= ei_abs(b) * prec; return ei_abs(a) <= ei_abs(b) * prec;
} }
inline bool ei_isApprox(double a, double b, double prec = precision<double>()) inline bool ei_isApprox(double a, double b, double prec = dummy_precision<double>())
{ {
return ei_abs(a - b) <= std::min(ei_abs(a), ei_abs(b)) * prec; return ei_abs(a - b) <= std::min(ei_abs(a), ei_abs(b)) * prec;
} }
inline bool ei_isApproxOrLessThan(double a, double b, double prec = precision<double>()) inline bool ei_isApproxOrLessThan(double a, double b, double prec = dummy_precision<double>())
{ {
return a <= b || ei_isApprox(a, b, prec); return a <= b || ei_isApprox(a, b, prec);
} }
@ -206,7 +206,7 @@ inline bool ei_isApproxOrLessThan(double a, double b, double prec = precision<do
*** complex<float> *** *** complex<float> ***
*********************/ *********************/
template<> inline float precision<std::complex<float> >() { return precision<float>(); } template<> inline float dummy_precision<std::complex<float> >() { return dummy_precision<float>(); }
inline float ei_real(const std::complex<float>& x) { return std::real(x); } inline float ei_real(const std::complex<float>& x) { return std::real(x); }
inline float ei_imag(const std::complex<float>& x) { return std::imag(x); } inline float ei_imag(const std::complex<float>& x) { return std::imag(x); }
inline float& ei_real_ref(std::complex<float>& x) { return reinterpret_cast<float*>(&x)[0]; } inline float& ei_real_ref(std::complex<float>& x) { return reinterpret_cast<float*>(&x)[0]; }
@ -224,15 +224,15 @@ template<> inline std::complex<float> ei_random()
{ {
return std::complex<float>(ei_random<float>(), ei_random<float>()); return std::complex<float>(ei_random<float>(), ei_random<float>());
} }
inline bool ei_isMuchSmallerThan(const std::complex<float>& a, const std::complex<float>& b, float prec = precision<float>()) inline bool ei_isMuchSmallerThan(const std::complex<float>& a, const std::complex<float>& b, float prec = dummy_precision<float>())
{ {
return ei_abs2(a) <= ei_abs2(b) * prec * prec; return ei_abs2(a) <= ei_abs2(b) * prec * prec;
} }
inline bool ei_isMuchSmallerThan(const std::complex<float>& a, float b, float prec = precision<float>()) inline bool ei_isMuchSmallerThan(const std::complex<float>& a, float b, float prec = dummy_precision<float>())
{ {
return ei_abs2(a) <= ei_abs2(b) * prec * prec; return ei_abs2(a) <= ei_abs2(b) * prec * prec;
} }
inline bool ei_isApprox(const std::complex<float>& a, const std::complex<float>& b, float prec = precision<float>()) inline bool ei_isApprox(const std::complex<float>& a, const std::complex<float>& b, float prec = dummy_precision<float>())
{ {
return ei_isApprox(ei_real(a), ei_real(b), prec) return ei_isApprox(ei_real(a), ei_real(b), prec)
&& ei_isApprox(ei_imag(a), ei_imag(b), prec); && ei_isApprox(ei_imag(a), ei_imag(b), prec);
@ -243,7 +243,7 @@ inline bool ei_isApprox(const std::complex<float>& a, const std::complex<float>&
*** complex<double> *** *** complex<double> ***
**********************/ **********************/
template<> inline double precision<std::complex<double> >() { return precision<double>(); } template<> inline double dummy_precision<std::complex<double> >() { return dummy_precision<double>(); }
inline double ei_real(const std::complex<double>& x) { return std::real(x); } inline double ei_real(const std::complex<double>& x) { return std::real(x); }
inline double ei_imag(const std::complex<double>& x) { return std::imag(x); } inline double ei_imag(const std::complex<double>& x) { return std::imag(x); }
inline double& ei_real_ref(std::complex<double>& x) { return reinterpret_cast<double*>(&x)[0]; } inline double& ei_real_ref(std::complex<double>& x) { return reinterpret_cast<double*>(&x)[0]; }
@ -261,15 +261,15 @@ template<> inline std::complex<double> ei_random()
{ {
return std::complex<double>(ei_random<double>(), ei_random<double>()); return std::complex<double>(ei_random<double>(), ei_random<double>());
} }
inline bool ei_isMuchSmallerThan(const std::complex<double>& a, const std::complex<double>& b, double prec = precision<double>()) inline bool ei_isMuchSmallerThan(const std::complex<double>& a, const std::complex<double>& b, double prec = dummy_precision<double>())
{ {
return ei_abs2(a) <= ei_abs2(b) * prec * prec; return ei_abs2(a) <= ei_abs2(b) * prec * prec;
} }
inline bool ei_isMuchSmallerThan(const std::complex<double>& a, double b, double prec = precision<double>()) inline bool ei_isMuchSmallerThan(const std::complex<double>& a, double b, double prec = dummy_precision<double>())
{ {
return ei_abs2(a) <= ei_abs2(b) * prec * prec; return ei_abs2(a) <= ei_abs2(b) * prec * prec;
} }
inline bool ei_isApprox(const std::complex<double>& a, const std::complex<double>& b, double prec = precision<double>()) inline bool ei_isApprox(const std::complex<double>& a, const std::complex<double>& b, double prec = dummy_precision<double>())
{ {
return ei_isApprox(ei_real(a), ei_real(b), prec) return ei_isApprox(ei_real(a), ei_real(b), prec)
&& ei_isApprox(ei_imag(a), ei_imag(b), prec); && ei_isApprox(ei_imag(a), ei_imag(b), prec);
@ -281,7 +281,7 @@ inline bool ei_isApprox(const std::complex<double>& a, const std::complex<double
*** long double *** *** long double ***
******************/ ******************/
template<> inline long double precision<long double>() { return precision<double>(); } template<> inline long double dummy_precision<long double>() { return dummy_precision<double>(); }
inline long double ei_real(long double x) { return x; } inline long double ei_real(long double x) { return x; }
inline long double& ei_real_ref(long double& x) { return x; } inline long double& ei_real_ref(long double& x) { return x; }
inline long double ei_imag(long double) { return 0.; } inline long double ei_imag(long double) { return 0.; }
@ -304,15 +304,15 @@ template<> inline long double ei_random()
{ {
return ei_random<double>(-ei_random_amplitude<double>(), ei_random_amplitude<double>()); return ei_random<double>(-ei_random_amplitude<double>(), ei_random_amplitude<double>());
} }
inline bool ei_isMuchSmallerThan(long double a, long double b, long double prec = precision<long double>()) inline bool ei_isMuchSmallerThan(long double a, long double b, long double prec = dummy_precision<long double>())
{ {
return ei_abs(a) <= ei_abs(b) * prec; return ei_abs(a) <= ei_abs(b) * prec;
} }
inline bool ei_isApprox(long double a, long double b, long double prec = precision<long double>()) inline bool ei_isApprox(long double a, long double b, long double prec = dummy_precision<long double>())
{ {
return ei_abs(a - b) <= std::min(ei_abs(a), ei_abs(b)) * prec; return ei_abs(a - b) <= std::min(ei_abs(a), ei_abs(b)) * prec;
} }
inline bool ei_isApproxOrLessThan(long double a, long double b, long double prec = precision<long double>()) inline bool ei_isApproxOrLessThan(long double a, long double b, long double prec = dummy_precision<long double>())
{ {
return a <= b || ei_isApprox(a, b, prec); return a <= b || ei_isApprox(a, b, prec);
} }
@ -321,7 +321,7 @@ inline bool ei_isApproxOrLessThan(long double a, long double b, long double prec
*** bool *** *** bool ***
**************/ **************/
template<> inline bool precision<bool>() { return 0; } template<> inline bool dummy_precision<bool>() { return 0; }
inline bool ei_real(bool x) { return x; } inline bool ei_real(bool x) { return x; }
inline bool& ei_real_ref(bool& x) { return x; } inline bool& ei_real_ref(bool& x) { return x; }
inline bool ei_imag(bool) { return 0; } inline bool ei_imag(bool) { return 0; }
@ -334,15 +334,15 @@ template<> inline bool ei_random()
{ {
return (ei_random<int>(0,1) == 1); return (ei_random<int>(0,1) == 1);
} }
inline bool ei_isMuchSmallerThan(bool a, bool, bool = precision<bool>()) inline bool ei_isMuchSmallerThan(bool a, bool, bool = dummy_precision<bool>())
{ {
return !a; return !a;
} }
inline bool ei_isApprox(bool a, bool b, bool = precision<bool>()) inline bool ei_isApprox(bool a, bool b, bool = dummy_precision<bool>())
{ {
return a == b; return a == b;
} }
inline bool ei_isApproxOrLessThan(bool a, bool b, bool = precision<bool>()) inline bool ei_isApproxOrLessThan(bool a, bool b, bool = dummy_precision<bool>())
{ {
return int(a) <= int(b); return int(a) <= int(b);
} }

View File

@ -558,27 +558,27 @@ template<typename Derived> class MatrixBase
template<typename OtherDerived> template<typename OtherDerived>
bool isApprox(const MatrixBase<OtherDerived>& other, bool isApprox(const MatrixBase<OtherDerived>& other,
RealScalar prec = precision<Scalar>()) const; RealScalar prec = dummy_precision<Scalar>()) const;
bool isMuchSmallerThan(const RealScalar& other, bool isMuchSmallerThan(const RealScalar& other,
RealScalar prec = precision<Scalar>()) const; RealScalar prec = dummy_precision<Scalar>()) const;
template<typename OtherDerived> template<typename OtherDerived>
bool isMuchSmallerThan(const MatrixBase<OtherDerived>& other, bool isMuchSmallerThan(const MatrixBase<OtherDerived>& other,
RealScalar prec = precision<Scalar>()) const; RealScalar prec = dummy_precision<Scalar>()) const;
bool isApproxToConstant(const Scalar& value, RealScalar prec = precision<Scalar>()) const; bool isApproxToConstant(const Scalar& value, RealScalar prec = dummy_precision<Scalar>()) const;
bool isConstant(const Scalar& value, RealScalar prec = precision<Scalar>()) const; bool isConstant(const Scalar& value, RealScalar prec = dummy_precision<Scalar>()) const;
bool isZero(RealScalar prec = precision<Scalar>()) const; bool isZero(RealScalar prec = dummy_precision<Scalar>()) const;
bool isOnes(RealScalar prec = precision<Scalar>()) const; bool isOnes(RealScalar prec = dummy_precision<Scalar>()) const;
bool isIdentity(RealScalar prec = precision<Scalar>()) const; bool isIdentity(RealScalar prec = dummy_precision<Scalar>()) const;
bool isDiagonal(RealScalar prec = precision<Scalar>()) const; bool isDiagonal(RealScalar prec = dummy_precision<Scalar>()) const;
bool isUpperTriangular(RealScalar prec = precision<Scalar>()) const; bool isUpperTriangular(RealScalar prec = dummy_precision<Scalar>()) const;
bool isLowerTriangular(RealScalar prec = precision<Scalar>()) const; bool isLowerTriangular(RealScalar prec = dummy_precision<Scalar>()) const;
template<typename OtherDerived> template<typename OtherDerived>
bool isOrthogonal(const MatrixBase<OtherDerived>& other, bool isOrthogonal(const MatrixBase<OtherDerived>& other,
RealScalar prec = precision<Scalar>()) const; RealScalar prec = dummy_precision<Scalar>()) const;
bool isUnitary(RealScalar prec = precision<Scalar>()) const; bool isUnitary(RealScalar prec = dummy_precision<Scalar>()) const;
template<typename OtherDerived> template<typename OtherDerived>
inline bool operator==(const MatrixBase<OtherDerived>& other) const inline bool operator==(const MatrixBase<OtherDerived>& other) const
@ -722,13 +722,13 @@ template<typename Derived> class MatrixBase
ResultType& inverse, ResultType& inverse,
typename ResultType::Scalar& determinant, typename ResultType::Scalar& determinant,
bool& invertible, bool& invertible,
const RealScalar& absDeterminantThreshold = precision<Scalar>() const RealScalar& absDeterminantThreshold = dummy_precision<Scalar>()
) const; ) const;
template<typename ResultType> template<typename ResultType>
void computeInverseWithCheck( void computeInverseWithCheck(
ResultType& inverse, ResultType& inverse,
bool& invertible, bool& invertible,
const RealScalar& absDeterminantThreshold = precision<Scalar>() const RealScalar& absDeterminantThreshold = dummy_precision<Scalar>()
) const; ) const;
Scalar determinant() const; Scalar determinant() const;

View File

@ -171,7 +171,7 @@ EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(_Scalar,_AmbientDim)
* determined by \a prec. * determined by \a prec.
* *
* \sa MatrixBase::isApprox() */ * \sa MatrixBase::isApprox() */
bool isApprox(const AlignedBox& other, typename NumTraits<Scalar>::Real prec = precision<Scalar>()) const bool isApprox(const AlignedBox& other, typename NumTraits<Scalar>::Real prec = dummy_precision<Scalar>()) const
{ return m_min.isApprox(other.m_min, prec) && m_max.isApprox(other.m_max, prec); } { return m_min.isApprox(other.m_min, prec) && m_max.isApprox(other.m_max, prec); }
protected: protected:

View File

@ -146,7 +146,7 @@ public:
* determined by \a prec. * determined by \a prec.
* *
* \sa MatrixBase::isApprox() */ * \sa MatrixBase::isApprox() */
bool isApprox(const AngleAxis& other, typename NumTraits<Scalar>::Real prec = precision<Scalar>()) const bool isApprox(const AngleAxis& other, typename NumTraits<Scalar>::Real prec = dummy_precision<Scalar>()) const
{ return m_axis.isApprox(other.m_axis, prec) && ei_isApprox(m_angle,other.m_angle, prec); } { return m_axis.isApprox(other.m_axis, prec) && ei_isApprox(m_angle,other.m_angle, prec); }
}; };
@ -165,7 +165,7 @@ template<typename QuatDerived>
AngleAxis<Scalar>& AngleAxis<Scalar>::operator=(const QuaternionBase<QuatDerived>& q) AngleAxis<Scalar>& AngleAxis<Scalar>::operator=(const QuaternionBase<QuatDerived>& q)
{ {
Scalar n2 = q.vec().squaredNorm(); Scalar n2 = q.vec().squaredNorm();
if (n2 < precision<Scalar>()*precision<Scalar>()) if (n2 < dummy_precision<Scalar>()*dummy_precision<Scalar>())
{ {
m_angle = 0; m_angle = 0;
m_axis << 1, 0, 0; m_axis << 1, 0, 0;

View File

@ -50,7 +50,7 @@ MatrixBase<Derived>::eulerAngles(int a0, int a1, int a2) const
Matrix<Scalar,3,1> res; Matrix<Scalar,3,1> res;
typedef Matrix<typename Derived::Scalar,2,1> Vector2; typedef Matrix<typename Derived::Scalar,2,1> Vector2;
const Scalar epsilon = precision<Scalar>(); const Scalar epsilon = dummy_precision<Scalar>();
const int odd = ((a0+1)%3 == a1) ? 0 : 1; const int odd = ((a0+1)%3 == a1) ? 0 : 1;
const int i = a0; const int i = a0;

View File

@ -257,7 +257,7 @@ public:
* determined by \a prec. * determined by \a prec.
* *
* \sa MatrixBase::isApprox() */ * \sa MatrixBase::isApprox() */
bool isApprox(const Hyperplane& other, typename NumTraits<Scalar>::Real prec = precision<Scalar>()) const bool isApprox(const Hyperplane& other, typename NumTraits<Scalar>::Real prec = dummy_precision<Scalar>()) const
{ return m_coeffs.isApprox(other.m_coeffs, prec); } { return m_coeffs.isApprox(other.m_coeffs, prec); }
protected: protected:

View File

@ -123,7 +123,7 @@ public:
* determined by \a prec. * determined by \a prec.
* *
* \sa MatrixBase::isApprox() */ * \sa MatrixBase::isApprox() */
bool isApprox(const ParametrizedLine& other, typename NumTraits<Scalar>::Real prec = precision<Scalar>()) const bool isApprox(const ParametrizedLine& other, typename NumTraits<Scalar>::Real prec = dummy_precision<Scalar>()) const
{ return m_origin.isApprox(other.m_origin, prec) && m_direction.isApprox(other.m_direction, prec); } { return m_origin.isApprox(other.m_origin, prec) && m_direction.isApprox(other.m_direction, prec); }
protected: protected:

View File

@ -163,7 +163,7 @@ public:
* *
* \sa MatrixBase::isApprox() */ * \sa MatrixBase::isApprox() */
template<class OtherDerived> template<class OtherDerived>
bool isApprox(const QuaternionBase<OtherDerived>& other, RealScalar prec = precision<Scalar>()) const bool isApprox(const QuaternionBase<OtherDerived>& other, RealScalar prec = dummy_precision<Scalar>()) const
{ return coeffs().isApprox(other.coeffs(), prec); } { return coeffs().isApprox(other.coeffs(), prec); }
/** return the result vector of \a v through the rotation*/ /** return the result vector of \a v through the rotation*/
@ -514,7 +514,7 @@ inline Derived& QuaternionBase<Derived>::setFromTwoVectors(const MatrixBase<Deri
// under the constraint: // under the constraint:
// ||x|| = 1 // ||x|| = 1
// which yields a singular value problem // which yields a singular value problem
if (c < Scalar(-1)+precision<Scalar>()) if (c < Scalar(-1)+dummy_precision<Scalar>())
{ {
c = std::max<Scalar>(c,-1); c = std::max<Scalar>(c,-1);
Matrix<Scalar,2,3> m; m << v0.transpose(), v1.transpose(); Matrix<Scalar,2,3> m; m << v0.transpose(), v1.transpose();
@ -590,7 +590,7 @@ template <class OtherDerived>
Quaternion<typename ei_traits<Derived>::Scalar> Quaternion<typename ei_traits<Derived>::Scalar>
QuaternionBase<Derived>::slerp(Scalar t, const QuaternionBase<OtherDerived>& other) const QuaternionBase<Derived>::slerp(Scalar t, const QuaternionBase<OtherDerived>& other) const
{ {
static const Scalar one = Scalar(1) - precision<Scalar>(); static const Scalar one = Scalar(1) - dummy_precision<Scalar>();
Scalar d = this->dot(other); Scalar d = this->dot(other);
Scalar absD = ei_abs(d); Scalar absD = ei_abs(d);
if (absD>=one) if (absD>=one)

View File

@ -121,7 +121,7 @@ public:
* determined by \a prec. * determined by \a prec.
* *
* \sa MatrixBase::isApprox() */ * \sa MatrixBase::isApprox() */
bool isApprox(const Rotation2D& other, typename NumTraits<Scalar>::Real prec = precision<Scalar>()) const bool isApprox(const Rotation2D& other, typename NumTraits<Scalar>::Real prec = dummy_precision<Scalar>()) const
{ return ei_isApprox(m_angle,other.m_angle, prec); } { return ei_isApprox(m_angle,other.m_angle, prec); }
}; };

View File

@ -107,7 +107,7 @@ public:
* determined by \a prec. * determined by \a prec.
* *
* \sa MatrixBase::isApprox() */ * \sa MatrixBase::isApprox() */
bool isApprox(const UniformScaling& other, typename NumTraits<Scalar>::Real prec = precision<Scalar>()) const bool isApprox(const UniformScaling& other, typename NumTraits<Scalar>::Real prec = dummy_precision<Scalar>()) const
{ return ei_isApprox(m_factor, other.factor(), prec); } { return ei_isApprox(m_factor, other.factor(), prec); }
}; };

View File

@ -424,7 +424,7 @@ public:
* determined by \a prec. * determined by \a prec.
* *
* \sa MatrixBase::isApprox() */ * \sa MatrixBase::isApprox() */
bool isApprox(const Transform& other, typename NumTraits<Scalar>::Real prec = precision<Scalar>()) const bool isApprox(const Transform& other, typename NumTraits<Scalar>::Real prec = dummy_precision<Scalar>()) const
{ return m_matrix.isApprox(other.m_matrix, prec); } { return m_matrix.isApprox(other.m_matrix, prec); }
/** Sets the last row to [0 ... 0 1] /** Sets the last row to [0 ... 0 1]

View File

@ -154,7 +154,7 @@ public:
* determined by \a prec. * determined by \a prec.
* *
* \sa MatrixBase::isApprox() */ * \sa MatrixBase::isApprox() */
bool isApprox(const Translation& other, typename NumTraits<Scalar>::Real prec = precision<Scalar>()) const bool isApprox(const Translation& other, typename NumTraits<Scalar>::Real prec = dummy_precision<Scalar>()) const
{ return m_coeffs.isApprox(other.m_coeffs, prec); } { return m_coeffs.isApprox(other.m_coeffs, prec); }
}; };

View File

@ -234,8 +234,9 @@ template<typename _MatrixType> class FullPivLU
* who need to determine when pivots are to be considered nonzero. This is not used for the * who need to determine when pivots are to be considered nonzero. This is not used for the
* LU decomposition itself. * LU decomposition itself.
* *
* When it needs to get the threshold value, Eigen calls threshold(). By default, this calls * When it needs to get the threshold value, Eigen calls threshold(). By default, this
* defaultThreshold(). Once you have called the present method setThreshold(const RealScalar&), * uses a formula to automatically determine a reasonable threshold.
* Once you have called the present method setThreshold(const RealScalar&),
* your value is used instead. * your value is used instead.
* *
* \param threshold The new value to use as the threshold. * \param threshold The new value to use as the threshold.
@ -303,7 +304,7 @@ template<typename _MatrixType> class FullPivLU
inline int dimensionOfKernel() const inline int dimensionOfKernel() const
{ {
ei_assert(m_isInitialized && "LU is not initialized."); ei_assert(m_isInitialized && "LU is not initialized.");
return m_lu.cols() - rank(); return cols() - rank();
} }
/** \returns true if the matrix of which *this is the LU decomposition represents an injective /** \returns true if the matrix of which *this is the LU decomposition represents an injective
@ -316,7 +317,7 @@ template<typename _MatrixType> class FullPivLU
inline bool isInjective() const inline bool isInjective() const
{ {
ei_assert(m_isInitialized && "LU is not initialized."); ei_assert(m_isInitialized && "LU is not initialized.");
return rank() == m_lu.cols(); return rank() == cols();
} }
/** \returns true if the matrix of which *this is the LU decomposition represents a surjective /** \returns true if the matrix of which *this is the LU decomposition represents a surjective
@ -329,7 +330,7 @@ template<typename _MatrixType> class FullPivLU
inline bool isSurjective() const inline bool isSurjective() const
{ {
ei_assert(m_isInitialized && "LU is not initialized."); ei_assert(m_isInitialized && "LU is not initialized.");
return rank() == m_lu.rows(); return rank() == rows();
} }
/** \returns true if the matrix of which *this is the LU decomposition is invertible. /** \returns true if the matrix of which *this is the LU decomposition is invertible.
@ -402,6 +403,7 @@ FullPivLU<MatrixType>& FullPivLU<MatrixType>::compute(const MatrixType& matrix)
m_nonzero_pivots = size; // the generic case is that in which all pivots are nonzero (invertible case) m_nonzero_pivots = size; // the generic case is that in which all pivots are nonzero (invertible case)
m_maxpivot = RealScalar(0); m_maxpivot = RealScalar(0);
for(int k = 0; k < size; ++k) for(int k = 0; k < size; ++k)
{ {
// First, we need to find the pivot. // First, we need to find the pivot.
@ -418,10 +420,10 @@ FullPivLU<MatrixType>& FullPivLU<MatrixType>::compute(const MatrixType& matrix)
// if the pivot (hence the corner) is exactly zero, terminate to avoid generating nan/inf values // if the pivot (hence the corner) is exactly zero, terminate to avoid generating nan/inf values
if(biggest_in_corner == RealScalar(0)) if(biggest_in_corner == RealScalar(0))
{ {
// before exiting, make sure to initialize the still uninitialized row_transpositions // before exiting, make sure to initialize the still uninitialized transpositions
// in a sane state without destroying what we already have. // in a sane state without destroying what we already have.
m_nonzero_pivots = k; m_nonzero_pivots = k;
for(int i = k; i < size; i++) for(int i = k; i < size; ++i)
{ {
rows_transpositions.coeffRef(i) = i; rows_transpositions.coeffRef(i) = i;
cols_transpositions.coeffRef(i) = i; cols_transpositions.coeffRef(i) = i;
@ -617,7 +619,7 @@ struct ei_solve_retval<FullPivLU<_MatrixType>, Rhs>
* Step 4: result = Q * c; * Step 4: result = Q * c;
*/ */
const int rows = dec().matrixLU().rows(), cols = dec().matrixLU().cols(), const int rows = dec().rows(), cols = dec().cols(),
nonzero_pivots = dec().nonzeroPivots(); nonzero_pivots = dec().nonzeroPivots();
ei_assert(rhs().rows() == rows); ei_assert(rhs().rows() == rows);
const int smalldim = std::min(rows, cols); const int smalldim = std::min(rows, cols);

View File

@ -74,7 +74,8 @@ template<typename _MatrixType> class ColPivHouseholderQR
ColPivHouseholderQR(const MatrixType& matrix) ColPivHouseholderQR(const MatrixType& matrix)
: m_qr(matrix.rows(), matrix.cols()), : m_qr(matrix.rows(), matrix.cols()),
m_hCoeffs(std::min(matrix.rows(),matrix.cols())), m_hCoeffs(std::min(matrix.rows(),matrix.cols())),
m_isInitialized(false) m_isInitialized(false),
m_usePrescribedThreshold(false)
{ {
compute(matrix); compute(matrix);
} }
@ -153,54 +154,63 @@ template<typename _MatrixType> class ColPivHouseholderQR
/** \returns the rank of the matrix of which *this is the QR decomposition. /** \returns the rank of the matrix of which *this is the QR decomposition.
* *
* \note This is computed at the time of the construction of the QR decomposition. This * \note This method has to determine which pivots should be considered nonzero.
* method does not perform any further computation. * For that, it uses the threshold value that you can control by calling
* setThreshold(const RealScalar&).
*/ */
inline int rank() const inline int rank() const
{ {
ei_assert(m_isInitialized && "ColPivHouseholderQR is not initialized."); ei_assert(m_isInitialized && "ColPivHouseholderQR is not initialized.");
return m_rank; RealScalar premultiplied_threshold = ei_abs(m_maxpivot) * threshold();
int result = 0;
for(int i = 0; i < m_nonzero_pivots; ++i)
result += (ei_abs(m_qr.coeff(i,i)) > premultiplied_threshold);
return result;
} }
/** \returns the dimension of the kernel of the matrix of which *this is the QR decomposition. /** \returns the dimension of the kernel of the matrix of which *this is the QR decomposition.
* *
* \note Since the rank is computed at the time of the construction of the QR decomposition, this * \note This method has to determine which pivots should be considered nonzero.
* method almost does not perform any further computation. * For that, it uses the threshold value that you can control by calling
* setThreshold(const RealScalar&).
*/ */
inline int dimensionOfKernel() const inline int dimensionOfKernel() const
{ {
ei_assert(m_isInitialized && "ColPivHouseholderQR is not initialized."); ei_assert(m_isInitialized && "ColPivHouseholderQR is not initialized.");
return m_qr.cols() - m_rank; return cols() - rank();
} }
/** \returns true if the matrix of which *this is the QR decomposition represents an injective /** \returns true if the matrix of which *this is the QR decomposition represents an injective
* linear map, i.e. has trivial kernel; false otherwise. * linear map, i.e. has trivial kernel; false otherwise.
* *
* \note Since the rank is computed at the time of the construction of the QR decomposition, this * \note This method has to determine which pivots should be considered nonzero.
* method almost does not perform any further computation. * For that, it uses the threshold value that you can control by calling
* setThreshold(const RealScalar&).
*/ */
inline bool isInjective() const inline bool isInjective() const
{ {
ei_assert(m_isInitialized && "ColPivHouseholderQR is not initialized."); ei_assert(m_isInitialized && "ColPivHouseholderQR is not initialized.");
return m_rank == m_qr.cols(); return rank() == cols();
} }
/** \returns true if the matrix of which *this is the QR decomposition represents a surjective /** \returns true if the matrix of which *this is the QR decomposition represents a surjective
* linear map; false otherwise. * linear map; false otherwise.
* *
* \note Since the rank is computed at the time of the construction of the QR decomposition, this * \note This method has to determine which pivots should be considered nonzero.
* method almost does not perform any further computation. * For that, it uses the threshold value that you can control by calling
* setThreshold(const RealScalar&).
*/ */
inline bool isSurjective() const inline bool isSurjective() const
{ {
ei_assert(m_isInitialized && "ColPivHouseholderQR is not initialized."); ei_assert(m_isInitialized && "ColPivHouseholderQR is not initialized.");
return m_rank == m_qr.rows(); return rank() == rows();
} }
/** \returns true if the matrix of which *this is the QR decomposition is invertible. /** \returns true if the matrix of which *this is the QR decomposition is invertible.
* *
* \note Since the rank is computed at the time of the construction of the QR decomposition, this * \note This method has to determine which pivots should be considered nonzero.
* method almost does not perform any further computation. * For that, it uses the threshold value that you can control by calling
* setThreshold(const RealScalar&).
*/ */
inline bool isInvertible() const inline bool isInvertible() const
{ {
@ -226,13 +236,80 @@ template<typename _MatrixType> class ColPivHouseholderQR
inline int cols() const { return m_qr.cols(); } inline int cols() const { return m_qr.cols(); }
const HCoeffsType& hCoeffs() const { return m_hCoeffs; } const HCoeffsType& hCoeffs() const { return m_hCoeffs; }
/** Allows to prescribe a threshold to be used by certain methods, such as rank(),
* who need to determine when pivots are to be considered nonzero. This is not used for the
* QR decomposition itself.
*
* When it needs to get the threshold value, Eigen calls threshold(). By default, this
* uses a formula to automatically determine a reasonable threshold.
* Once you have called the present method setThreshold(const RealScalar&),
* your value is used instead.
*
* \param threshold The new value to use as the threshold.
*
* A pivot will be considered nonzero if its absolute value is strictly greater than
* \f$ \vert pivot \vert \leqslant threshold \times \vert maxpivot \vert \f$
* where maxpivot is the biggest pivot.
*
* If you want to come back to the default behavior, call setThreshold(Default_t)
*/
ColPivHouseholderQR& setThreshold(const RealScalar& threshold)
{
m_usePrescribedThreshold = true;
m_prescribedThreshold = threshold;
}
/** Allows to come back to the default behavior, letting Eigen use its default formula for
* determining the threshold.
*
* You should pass the special object Eigen::Default as parameter here.
* \code qr.setThreshold(Eigen::Default); \endcode
*
* See the documentation of setThreshold(const RealScalar&).
*/
ColPivHouseholderQR& setThreshold(Default_t)
{
m_usePrescribedThreshold = false;
}
/** Returns the threshold that will be used by certain methods such as rank().
*
* See the documentation of setThreshold(const RealScalar&).
*/
RealScalar threshold() const
{
ei_assert(m_isInitialized || m_usePrescribedThreshold);
return m_usePrescribedThreshold ? m_prescribedThreshold
// this formula comes from experimenting (see "LU precision tuning" thread on the list)
// and turns out to be identical to Higham's formula used already in LDLt.
: epsilon<Scalar>() * m_qr.diagonalSize();
}
/** \returns the number of nonzero pivots in the QR decomposition.
* Here nonzero is meant in the exact sense, not in a fuzzy sense.
* So that notion isn't really intrinsically interesting, but it is
* still useful when implementing algorithms.
*
* \sa rank()
*/
inline int nonzeroPivots() const
{
ei_assert(m_isInitialized && "LU is not initialized.");
return m_nonzero_pivots;
}
/** \returns the absolute value of the biggest pivot, i.e. the biggest
* diagonal coefficient of U.
*/
RealScalar maxPivot() const { return m_maxpivot; }
protected: protected:
MatrixType m_qr; MatrixType m_qr;
HCoeffsType m_hCoeffs; HCoeffsType m_hCoeffs;
PermutationType m_cols_permutation; PermutationType m_cols_permutation;
bool m_isInitialized; bool m_isInitialized, m_usePrescribedThreshold;
RealScalar m_precision; RealScalar m_prescribedThreshold, m_maxpivot;
int m_rank; int m_nonzero_pivots;
int m_det_pq; int m_det_pq;
}; };
@ -259,61 +336,84 @@ ColPivHouseholderQR<MatrixType>& ColPivHouseholderQR<MatrixType>::compute(const
{ {
int rows = matrix.rows(); int rows = matrix.rows();
int cols = matrix.cols(); int cols = matrix.cols();
int size = std::min(rows,cols); int size = matrix.diagonalSize();
m_rank = size;
m_qr = matrix; m_qr = matrix;
m_hCoeffs.resize(size); m_hCoeffs.resize(size);
RowVectorType temp(cols); RowVectorType temp(cols);
m_precision = epsilon<Scalar>() * size;
IntRowVectorType cols_transpositions(matrix.cols()); IntRowVectorType cols_transpositions(matrix.cols());
int number_of_transpositions = 0; int number_of_transpositions = 0;
RealRowVectorType colSqNorms(cols); RealRowVectorType colSqNorms(cols);
for(int k = 0; k < cols; ++k) for(int k = 0; k < cols; ++k)
colSqNorms.coeffRef(k) = m_qr.col(k).squaredNorm(); colSqNorms.coeffRef(k) = m_qr.col(k).squaredNorm();
RealScalar biggestColSqNorm = colSqNorms.maxCoeff();
for (int k = 0; k < size; ++k) RealScalar threshold_helper = colSqNorms.maxCoeff() * ei_abs2(epsilon<Scalar>()) / rows;
{
int biggest_col_in_corner;
RealScalar biggestColSqNormInCorner = colSqNorms.end(cols-k).maxCoeff(&biggest_col_in_corner);
biggest_col_in_corner += k;
// if the corner is negligible, then we have less than full rank, and we can finish early m_nonzero_pivots = size; // the generic case is that in which all pivots are nonzero (invertible case)
if(ei_isMuchSmallerThan(biggestColSqNormInCorner, biggestColSqNorm, m_precision)) m_maxpivot = RealScalar(0);
for(int k = 0; k < size; ++k)
{ {
m_rank = k; // first, we look up in our table colSqNorms which column has the biggest squared norm
for(int i = k; i < size; i++) int biggest_col_index;
RealScalar biggest_col_sq_norm = colSqNorms.end(cols-k).maxCoeff(&biggest_col_index);
biggest_col_index += k;
// since our table colSqNorms accumulates imprecision at every step, we must now recompute
// the actual squared norm of the selected column.
// Note that not doing so does result in solve() sometimes returning inf/nan values
// when running the unit test with 1000 repetitions.
biggest_col_sq_norm = m_qr.col(biggest_col_index).end(rows-k).squaredNorm();
// we store that back into our table: it can't hurt to correct our table.
colSqNorms.coeffRef(biggest_col_index) = biggest_col_sq_norm;
// if the current biggest column is smaller than epsilon times the initial biggest column,
// terminate to avoid generating nan/inf values.
// Note that here, if we test instead for "biggest == 0", we get a failure every 1000 (or so)
// repetitions of the unit test, with the result of solve() filled with large values of the order
// of 1/(size*epsilon).
if(biggest_col_sq_norm < threshold_helper * (rows-k))
{ {
cols_transpositions.coeffRef(i) = i; m_nonzero_pivots = k;
m_hCoeffs.coeffRef(i) = Scalar(0); m_hCoeffs.end(size-k).setZero();
} m_qr.corner(BottomRight,rows-k,cols-k)
.template triangularView<StrictlyLowerTriangular>()
.setZero();
break; break;
} }
cols_transpositions.coeffRef(k) = biggest_col_in_corner; // apply the transposition to the columns
if(k != biggest_col_in_corner) { cols_transpositions.coeffRef(k) = biggest_col_index;
m_qr.col(k).swap(m_qr.col(biggest_col_in_corner)); if(k != biggest_col_index) {
std::swap(colSqNorms.coeffRef(k), colSqNorms.coeffRef(biggest_col_in_corner)); m_qr.col(k).swap(m_qr.col(biggest_col_index));
std::swap(colSqNorms.coeffRef(k), colSqNorms.coeffRef(biggest_col_index));
++number_of_transpositions; ++number_of_transpositions;
} }
// generate the householder vector, store it below the diagonal
RealScalar beta; RealScalar beta;
m_qr.col(k).end(rows-k).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta); m_qr.col(k).end(rows-k).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta);
// apply the householder transformation to the diagonal coefficient
m_qr.coeffRef(k,k) = beta; m_qr.coeffRef(k,k) = beta;
// remember the maximum absolute value of diagonal coefficients
if(ei_abs(beta) > m_maxpivot) m_maxpivot = ei_abs(beta);
// apply the householder transformation
m_qr.corner(BottomRight, rows-k, cols-k-1) m_qr.corner(BottomRight, rows-k, cols-k-1)
.applyHouseholderOnTheLeft(m_qr.col(k).end(rows-k-1), m_hCoeffs.coeffRef(k), &temp.coeffRef(k+1)); .applyHouseholderOnTheLeft(m_qr.col(k).end(rows-k-1), m_hCoeffs.coeffRef(k), &temp.coeffRef(k+1));
// update our table of squared norms of the columns
colSqNorms.end(cols-k-1) -= m_qr.row(k).end(cols-k-1).cwise().abs2(); colSqNorms.end(cols-k-1) -= m_qr.row(k).end(cols-k-1).cwise().abs2();
} }
m_cols_permutation.setIdentity(cols); m_cols_permutation.setIdentity(cols);
for(int k = 0; k < size; ++k) for(int k = 0; k < m_nonzero_pivots; ++k)
m_cols_permutation.applyTranspositionOnTheRight(k, cols_transpositions.coeff(k)); m_cols_permutation.applyTranspositionOnTheRight(k, cols_transpositions.coeff(k));
m_det_pq = (number_of_transpositions%2) ? -1 : 1; m_det_pq = (number_of_transpositions%2) ? -1 : 1;
@ -330,13 +430,12 @@ struct ei_solve_retval<ColPivHouseholderQR<_MatrixType>, Rhs>
template<typename Dest> void evalTo(Dest& dst) const template<typename Dest> void evalTo(Dest& dst) const
{ {
const int rows = dec().rows(), cols = dec().cols(); const int rows = dec().rows(), cols = dec().cols(),
nonzero_pivots = dec().nonzeroPivots();
dst.resize(cols, rhs().cols()); dst.resize(cols, rhs().cols());
ei_assert(rhs().rows() == rows); ei_assert(rhs().rows() == rows);
// FIXME introduce nonzeroPivots() and use it here. and more generally, if(nonzero_pivots == 0)
// make the same improvements in this dec as in FullPivLU.
if(dec().rank()==0)
{ {
dst.setZero(); dst.setZero();
return; return;
@ -346,28 +445,26 @@ struct ei_solve_retval<ColPivHouseholderQR<_MatrixType>, Rhs>
// Note that the matrix Q = H_0^* H_1^*... so its inverse is Q^* = (H_0 H_1 ...)^T // Note that the matrix Q = H_0^* H_1^*... so its inverse is Q^* = (H_0 H_1 ...)^T
c.applyOnTheLeft(householderSequence( c.applyOnTheLeft(householderSequence(
dec().matrixQR().corner(TopLeft,rows,dec().rank()), dec().matrixQR(),
dec().hCoeffs().start(dec().rank())).transpose() dec().hCoeffs(),
); true
));
if(!dec().isSurjective())
{
// is c is in the image of R ?
RealScalar biggest_in_upper_part_of_c = c.corner(TopLeft, dec().rank(), c.cols()).cwise().abs().maxCoeff();
RealScalar biggest_in_lower_part_of_c = c.corner(BottomLeft, rows-dec().rank(), c.cols()).cwise().abs().maxCoeff();
// FIXME brain dead
const RealScalar m_precision = epsilon<Scalar>() * std::min(rows,cols);
if(!ei_isMuchSmallerThan(biggest_in_lower_part_of_c, biggest_in_upper_part_of_c, m_precision*4))
return;
}
dec().matrixQR() dec().matrixQR()
.corner(TopLeft, dec().rank(), dec().rank()) .corner(TopLeft, nonzero_pivots, nonzero_pivots)
.template triangularView<UpperTriangular>() .template triangularView<UpperTriangular>()
.solveInPlace(c.corner(TopLeft, dec().rank(), c.cols())); .solveInPlace(c.corner(TopLeft, nonzero_pivots, c.cols()));
for(int i = 0; i < dec().rank(); ++i) dst.row(dec().colsPermutation().indices().coeff(i)) = c.row(i);
for(int i = dec().rank(); i < cols; ++i) dst.row(dec().colsPermutation().indices().coeff(i)).setZero(); typename Rhs::PlainMatrixType d(c);
d.corner(TopLeft, nonzero_pivots, c.cols())
= dec().matrixQR()
.corner(TopLeft, nonzero_pivots, nonzero_pivots)
.template triangularView<UpperTriangular>()
* c.corner(TopLeft, nonzero_pivots, c.cols());
for(int i = 0; i < nonzero_pivots; ++i) dst.row(dec().colsPermutation().indices().coeff(i)) = c.row(i);
for(int i = nonzero_pivots; i < cols; ++i) dst.row(dec().colsPermutation().indices().coeff(i)).setZero();
} }
}; };
@ -376,7 +473,7 @@ template<typename MatrixType>
typename ColPivHouseholderQR<MatrixType>::HouseholderSequenceType ColPivHouseholderQR<MatrixType>::matrixQ() const typename ColPivHouseholderQR<MatrixType>::HouseholderSequenceType ColPivHouseholderQR<MatrixType>::matrixQ() const
{ {
ei_assert(m_isInitialized && "ColPivHouseholderQR is not initialized."); ei_assert(m_isInitialized && "ColPivHouseholderQR is not initialized.");
return HouseholderSequenceType(m_qr, m_hCoeffs.conjugate()); return HouseholderSequenceType(m_qr, m_hCoeffs.conjugate(), false);
} }
#endif // EIGEN_HIDE_HEAVY_CODE #endif // EIGEN_HIDE_HEAVY_CODE

View File

@ -193,7 +193,7 @@ SVD<MatrixType>& SVD<MatrixType>::compute(const MatrixType& matrix)
int i,its,j,k,l,nm; int i,its,j,k,l,nm;
Scalar anorm, c, f, g, h, s, scale, x, y, z; Scalar anorm, c, f, g, h, s, scale, x, y, z;
bool convergence = true; bool convergence = true;
Scalar eps = precision<Scalar>(); Scalar eps = dummy_precision<Scalar>();
Matrix<Scalar,Dynamic,1> rv1(n); Matrix<Scalar,Dynamic,1> rv1(n);
g = scale = anorm = 0; g = scale = anorm = 0;

View File

@ -300,7 +300,7 @@ class AmbiVector<_Scalar>::Iterator
* In practice, all coefficients having a magnitude smaller than \a epsilon * In practice, all coefficients having a magnitude smaller than \a epsilon
* are skipped. * are skipped.
*/ */
Iterator(const AmbiVector& vec, RealScalar epsilon = RealScalar(0.1)*precision<RealScalar>()) Iterator(const AmbiVector& vec, RealScalar epsilon = RealScalar(0.1)*dummy_precision<RealScalar>())
: m_vector(vec) : m_vector(vec)
{ {
m_epsilon = epsilon; m_epsilon = epsilon;

View File

@ -185,7 +185,7 @@ class CompressedStorage
return m_values[id]; return m_values[id];
} }
void prune(Scalar reference, RealScalar epsilon = precision<RealScalar>()) void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar>())
{ {
size_t k = 0; size_t k = 0;
size_t n = size(); size_t n = size();

View File

@ -208,7 +208,7 @@ class DynamicSparseMatrix
inline void finalize() {} inline void finalize() {}
void prune(Scalar reference, RealScalar epsilon = precision<RealScalar>()) void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar>())
{ {
for (int j=0; j<outerSize(); ++j) for (int j=0; j<outerSize(); ++j)
m_data[j].prune(reference,epsilon); m_data[j].prune(reference,epsilon);

View File

@ -94,7 +94,7 @@ class SparseLDLT
: m_flags(flags), m_status(0) : m_flags(flags), m_status(0)
{ {
ei_assert((MatrixType::Flags&RowMajorBit)==0); ei_assert((MatrixType::Flags&RowMajorBit)==0);
m_precision = RealScalar(0.1) * Eigen::precision<RealScalar>(); m_precision = RealScalar(0.1) * Eigen::dummy_precision<RealScalar>();
} }
/** Creates a LDLT object and compute the respective factorization of \a matrix using /** Creates a LDLT object and compute the respective factorization of \a matrix using
@ -103,7 +103,7 @@ class SparseLDLT
: m_matrix(matrix.rows(), matrix.cols()), m_flags(flags), m_status(0) : m_matrix(matrix.rows(), matrix.cols()), m_flags(flags), m_status(0)
{ {
ei_assert((MatrixType::Flags&RowMajorBit)==0); ei_assert((MatrixType::Flags&RowMajorBit)==0);
m_precision = RealScalar(0.1) * Eigen::precision<RealScalar>(); m_precision = RealScalar(0.1) * Eigen::dummy_precision<RealScalar>();
compute(matrix); compute(matrix);
} }

View File

@ -54,7 +54,7 @@ class SparseLLT
SparseLLT(int flags = 0) SparseLLT(int flags = 0)
: m_flags(flags), m_status(0) : m_flags(flags), m_status(0)
{ {
m_precision = RealScalar(0.1) * Eigen::precision<RealScalar>(); m_precision = RealScalar(0.1) * Eigen::dummy_precision<RealScalar>();
} }
/** Creates a LLT object and compute the respective factorization of \a matrix using /** Creates a LLT object and compute the respective factorization of \a matrix using
@ -62,7 +62,7 @@ class SparseLLT
SparseLLT(const MatrixType& matrix, int flags = 0) SparseLLT(const MatrixType& matrix, int flags = 0)
: m_matrix(matrix.rows(), matrix.cols()), m_flags(flags), m_status(0) : m_matrix(matrix.rows(), matrix.cols()), m_flags(flags), m_status(0)
{ {
m_precision = RealScalar(0.1) * Eigen::precision<RealScalar>(); m_precision = RealScalar(0.1) * Eigen::dummy_precision<RealScalar>();
compute(matrix); compute(matrix);
} }

View File

@ -59,7 +59,7 @@ class SparseLU
SparseLU(int flags = 0) SparseLU(int flags = 0)
: m_flags(flags), m_status(0) : m_flags(flags), m_status(0)
{ {
m_precision = RealScalar(0.1) * Eigen::precision<RealScalar>(); m_precision = RealScalar(0.1) * Eigen::dummy_precision<RealScalar>();
} }
/** Creates a LU object and compute the respective factorization of \a matrix using /** Creates a LU object and compute the respective factorization of \a matrix using
@ -67,7 +67,7 @@ class SparseLU
SparseLU(const MatrixType& matrix, int flags = 0) SparseLU(const MatrixType& matrix, int flags = 0)
: /*m_matrix(matrix.rows(), matrix.cols()),*/ m_flags(flags), m_status(0) : /*m_matrix(matrix.rows(), matrix.cols()),*/ m_flags(flags), m_status(0)
{ {
m_precision = RealScalar(0.1) * Eigen::precision<RealScalar>(); m_precision = RealScalar(0.1) * Eigen::dummy_precision<RealScalar>();
compute(matrix); compute(matrix);
} }

View File

@ -341,7 +341,7 @@ class SparseMatrix
} }
} }
void prune(Scalar reference, RealScalar epsilon = precision<RealScalar>()) void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar>())
{ {
int k = 0; int k = 0;
for (int j=0; j<m_outerSize; ++j) for (int j=0; j<m_outerSize; ++j)

View File

@ -467,32 +467,32 @@ template<typename Derived> class SparseMatrixBase : public AnyMatrixBase<Derived
template<typename OtherDerived> template<typename OtherDerived>
bool isApprox(const SparseMatrixBase<OtherDerived>& other, bool isApprox(const SparseMatrixBase<OtherDerived>& other,
RealScalar prec = precision<Scalar>()) const RealScalar prec = dummy_precision<Scalar>()) const
{ return toDense().isApprox(other.toDense(),prec); } { return toDense().isApprox(other.toDense(),prec); }
template<typename OtherDerived> template<typename OtherDerived>
bool isApprox(const MatrixBase<OtherDerived>& other, bool isApprox(const MatrixBase<OtherDerived>& other,
RealScalar prec = precision<Scalar>()) const RealScalar prec = dummy_precision<Scalar>()) const
{ return toDense().isApprox(other,prec); } { return toDense().isApprox(other,prec); }
// bool isMuchSmallerThan(const RealScalar& other, // bool isMuchSmallerThan(const RealScalar& other,
// RealScalar prec = precision<Scalar>()) const; // RealScalar prec = dummy_precision<Scalar>()) const;
// template<typename OtherDerived> // template<typename OtherDerived>
// bool isMuchSmallerThan(const MatrixBase<OtherDerived>& other, // bool isMuchSmallerThan(const MatrixBase<OtherDerived>& other,
// RealScalar prec = precision<Scalar>()) const; // RealScalar prec = dummy_precision<Scalar>()) const;
// bool isApproxToConstant(const Scalar& value, RealScalar prec = precision<Scalar>()) const; // bool isApproxToConstant(const Scalar& value, RealScalar prec = dummy_precision<Scalar>()) const;
// bool isZero(RealScalar prec = precision<Scalar>()) const; // bool isZero(RealScalar prec = dummy_precision<Scalar>()) const;
// bool isOnes(RealScalar prec = precision<Scalar>()) const; // bool isOnes(RealScalar prec = dummy_precision<Scalar>()) const;
// bool isIdentity(RealScalar prec = precision<Scalar>()) const; // bool isIdentity(RealScalar prec = dummy_precision<Scalar>()) const;
// bool isDiagonal(RealScalar prec = precision<Scalar>()) const; // bool isDiagonal(RealScalar prec = dummy_precision<Scalar>()) const;
// bool isUpperTriangular(RealScalar prec = precision<Scalar>()) const; // bool isUpperTriangular(RealScalar prec = dummy_precision<Scalar>()) const;
// bool isLowerTriangular(RealScalar prec = precision<Scalar>()) const; // bool isLowerTriangular(RealScalar prec = dummy_precision<Scalar>()) const;
// template<typename OtherDerived> // template<typename OtherDerived>
// bool isOrthogonal(const MatrixBase<OtherDerived>& other, // bool isOrthogonal(const MatrixBase<OtherDerived>& other,
// RealScalar prec = precision<Scalar>()) const; // RealScalar prec = dummy_precision<Scalar>()) const;
// bool isUnitary(RealScalar prec = precision<Scalar>()) const; // bool isUnitary(RealScalar prec = dummy_precision<Scalar>()) const;
// template<typename OtherDerived> // template<typename OtherDerived>
// inline bool operator==(const MatrixBase<OtherDerived>& other) const // inline bool operator==(const MatrixBase<OtherDerived>& other) const

View File

@ -201,7 +201,7 @@ class SparseVector
EIGEN_DEPRECATED void endFill() {} EIGEN_DEPRECATED void endFill() {}
inline void finalize() {} inline void finalize() {}
void prune(Scalar reference, RealScalar epsilon = precision<RealScalar>()) void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar>())
{ {
m_data.prune(reference,epsilon); m_data.prune(reference,epsilon);
} }

View File

@ -107,7 +107,10 @@ endmacro(ei_add_test_internal)
# #
# Again, ctest -R allows to run all matching tests. # Again, ctest -R allows to run all matching tests.
macro(ei_add_test testname) macro(ei_add_test testname)
set(cmake_tests_list "${cmake_tests_list}${testname}\n") get_property(EIGEN_TESTS_LIST GLOBAL PROPERTY EIGEN_TESTS_LIST)
set(EIGEN_TESTS_LIST "${EIGEN_TESTS_LIST}${testname}\n")
set_property(GLOBAL PROPERTY EIGEN_TESTS_LIST "${EIGEN_TESTS_LIST}")
file(READ "${testname}.cpp" test_source) file(READ "${testname}.cpp" test_source)
set(parts 0) set(parts 0)
string(REGEX MATCHALL "CALL_SUBTEST_[0-9]+|EIGEN_TEST_PART_[0-9]+" string(REGEX MATCHALL "CALL_SUBTEST_[0-9]+|EIGEN_TEST_PART_[0-9]+"
@ -204,10 +207,12 @@ macro(ei_init_testing)
define_property(GLOBAL PROPERTY EIGEN_TESTED_BACKENDS BRIEF_DOCS " " FULL_DOCS " ") define_property(GLOBAL PROPERTY EIGEN_TESTED_BACKENDS BRIEF_DOCS " " FULL_DOCS " ")
define_property(GLOBAL PROPERTY EIGEN_MISSING_BACKENDS BRIEF_DOCS " " FULL_DOCS " ") define_property(GLOBAL PROPERTY EIGEN_MISSING_BACKENDS BRIEF_DOCS " " FULL_DOCS " ")
define_property(GLOBAL PROPERTY EIGEN_TESTING_SUMMARY BRIEF_DOCS " " FULL_DOCS " ") define_property(GLOBAL PROPERTY EIGEN_TESTING_SUMMARY BRIEF_DOCS " " FULL_DOCS " ")
define_property(GLOBAL PROPERTY EIGEN_TESTS_LIST BRIEF_DOCS " " FULL_DOCS " ")
set_property(GLOBAL PROPERTY EIGEN_TESTED_BACKENDS "") set_property(GLOBAL PROPERTY EIGEN_TESTED_BACKENDS "")
set_property(GLOBAL PROPERTY EIGEN_MISSING_BACKENDS "") set_property(GLOBAL PROPERTY EIGEN_MISSING_BACKENDS "")
set_property(GLOBAL PROPERTY EIGEN_TESTING_SUMMARY "") set_property(GLOBAL PROPERTY EIGEN_TESTING_SUMMARY "")
set_property(GLOBAL PROPERTY EIGEN_TESTS_LIST "")
endmacro(ei_init_testing) endmacro(ei_init_testing)
if(CMAKE_COMPILER_IS_GNUCXX) if(CMAKE_COMPILER_IS_GNUCXX)
@ -218,7 +223,7 @@ if(CMAKE_COMPILER_IS_GNUCXX)
set(COVERAGE_FLAGS "") set(COVERAGE_FLAGS "")
endif(EIGEN_COVERAGE_TESTING) endif(EIGEN_COVERAGE_TESTING)
if(EIGEN_TEST_RVALUE_REF_SUPPORT OR EIGEN_TEST_C++0x) if(EIGEN_TEST_RVALUE_REF_SUPPORT OR EIGEN_TEST_C++0x)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++0x") set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=gnu++0x")
endif(EIGEN_TEST_RVALUE_REF_SUPPORT OR EIGEN_TEST_C++0x) endif(EIGEN_TEST_RVALUE_REF_SUPPORT OR EIGEN_TEST_C++0x)
if(CMAKE_SYSTEM_NAME MATCHES Linux) if(CMAKE_SYSTEM_NAME MATCHES Linux)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${COVERAGE_FLAGS} -g2") set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${COVERAGE_FLAGS} -g2")

View File

@ -1,3 +1,6 @@
get_property(EIGEN_TESTS_LIST GLOBAL PROPERTY EIGEN_TESTS_LIST)
configure_file(buildtests.in ${CMAKE_BINARY_DIR}/buildtests)
configure_file(check.in ${CMAKE_BINARY_DIR}/check) configure_file(check.in ${CMAKE_BINARY_DIR}/check)
configure_file(debug.in ${CMAKE_BINARY_DIR}/debug) configure_file(debug.in ${CMAKE_BINARY_DIR}/debug)
configure_file(release.in ${CMAKE_BINARY_DIR}/release) configure_file(release.in ${CMAKE_BINARY_DIR}/release)

View File

@ -3,12 +3,13 @@
if [ $# == 0 -o $# -ge 3 ] if [ $# == 0 -o $# -ge 3 ]
then then
echo "usage: ./buildtests regexp [jobs]" echo "usage: ./buildtests regexp [jobs]"
echo " makes tests matching the regexp, with <jobs> concurrent make jobs" echo " makes tests matching the regexp, with [jobs] concurrent make jobs"
exit 0 exit 0
fi fi
TESTSLIST="${cmake_tests_list}" TESTSLIST="${EIGEN_TESTS_LIST}"
targets_to_make=`echo "$TESTSLIST" | grep "$1" | sed s/^/test_/g | xargs echo`
targets_to_make=`echo "$TESTSLIST" | egrep "$1" | sed s/^/test_/g | xargs echo`
if [ $# == 1 ] if [ $# == 1 ]
then then

View File

@ -4,9 +4,9 @@
if [ $# == 0 -o $# -ge 3 ] if [ $# == 0 -o $# -ge 3 ]
then then
echo "usage: ./check regexp [jobs]" echo "usage: ./check regexp [jobs]"
echo " makes and runs tests matching the regexp, with <jobs> concurrent make jobs" echo " makes and runs tests matching the regexp, with [jobs] concurrent make jobs"
exit 0 exit 0
fi fi
# TODO when ctest 2.8 comes out, honor the jobs parameter # TODO when ctest 2.8 comes out, honor the jobs parameter
./buildtests $* && ctest -R $1 ./buildtests "$1" "$2" && ctest -R "$1"

View File

@ -1 +1 @@
This file is just there as a signature to help identify directories containing Eigen3. When writing for a script looking for Eigen3, just look for this file. This is especially useful to help disambiguate with Eigen2... This file is just there as a signature to help identify directories containing Eigen3. When writing a script looking for Eigen3, just look for this file. This is especially useful to help disambiguate with Eigen2...

View File

@ -3,10 +3,6 @@ add_custom_target(buildtests)
add_custom_target(check COMMAND "ctest") add_custom_target(check COMMAND "ctest")
add_dependencies(check buildtests) add_dependencies(check buildtests)
include(EigenTesting)
ei_init_testing()
option(EIGEN_SPLIT_LARGE_TESTS "Split large tests into smaller executables" ON) option(EIGEN_SPLIT_LARGE_TESTS "Split large tests into smaller executables" ON)
find_package(GSL) find_package(GSL)
@ -169,5 +165,3 @@ if(CMAKE_COMPILER_IS_GNUCXX)
endif(CMAKE_COMPILER_IS_GNUCXX) endif(CMAKE_COMPILER_IS_GNUCXX)
ei_add_property(EIGEN_TESTING_SUMMARY "CXX_FLAGS: ${CMAKE_CXX_FLAGS}\n") ei_add_property(EIGEN_TESTING_SUMMARY "CXX_FLAGS: ${CMAKE_CXX_FLAGS}\n")
ei_add_property(EIGEN_TESTING_SUMMARY "Sparse lib flags: ${SPARSE_LIBS}\n") ei_add_property(EIGEN_TESTING_SUMMARY "Sparse lib flags: ${SPARSE_LIBS}\n")
configure_file(buildtests.in ${CMAKE_BINARY_DIR}/buildtests)

View File

@ -27,7 +27,7 @@
#include <Eigen/QR> #include <Eigen/QR>
template<typename Derived1, typename Derived2> template<typename Derived1, typename Derived2>
bool areNotApprox(const MatrixBase<Derived1>& m1, const MatrixBase<Derived2>& m2, typename Derived1::RealScalar epsilon = precision<typename Derived1::RealScalar>()) bool areNotApprox(const MatrixBase<Derived1>& m1, const MatrixBase<Derived2>& m2, typename Derived1::RealScalar epsilon = dummy_precision<typename Derived1::RealScalar>())
{ {
return !((m1-m2).cwise().abs2().maxCoeff() < epsilon * epsilon return !((m1-m2).cwise().abs2().maxCoeff() < epsilon * epsilon
* std::max(m1.cwise().abs2().maxCoeff(), m2.cwise().abs2().maxCoeff())); * std::max(m1.cwise().abs2().maxCoeff(), m2.cwise().abs2().maxCoeff()));

View File

@ -28,11 +28,11 @@
template<typename MatrixType> void qr() template<typename MatrixType> void qr()
{ {
// int rows = ei_random<int>(20,200), cols = ei_random<int>(20,200), cols2 = ei_random<int>(20,200); int rows = ei_random<int>(2,200), cols = ei_random<int>(2,200), cols2 = ei_random<int>(2,200);
int rows=3, cols=3, cols2=3;
int rank = ei_random<int>(1, std::min(rows, cols)-1); int rank = ei_random<int>(1, std::min(rows, cols)-1);
typedef typename MatrixType::Scalar Scalar; typedef typename MatrixType::Scalar Scalar;
typedef typename MatrixType::RealScalar RealScalar;
typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, MatrixType::RowsAtCompileTime> MatrixQType; typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, MatrixType::RowsAtCompileTime> MatrixQType;
typedef Matrix<Scalar, MatrixType::ColsAtCompileTime, 1> VectorType; typedef Matrix<Scalar, MatrixType::ColsAtCompileTime, 1> VectorType;
MatrixType m1; MatrixType m1;
@ -44,19 +44,11 @@ template<typename MatrixType> void qr()
VERIFY(!qr.isInvertible()); VERIFY(!qr.isInvertible());
VERIFY(!qr.isSurjective()); VERIFY(!qr.isSurjective());
MatrixType r = qr.matrixQR();
MatrixQType q = qr.matrixQ(); MatrixQType q = qr.matrixQ();
VERIFY_IS_UNITARY(q); VERIFY_IS_UNITARY(q);
// FIXME need better way to construct trapezoid MatrixType r = qr.matrixQR().template triangularView<UpperTriangular>();
for(int i = 0; i < rows; i++) for(int j = 0; j < cols; j++) if(i>j) r(i,j) = Scalar(0); MatrixType c = q * r * qr.colsPermutation().inverse();
MatrixType b = qr.matrixQ() * r;
MatrixType c = MatrixType::Zero(rows,cols);
for(int i = 0; i < cols; ++i) c.col(qr.colsPermutation().indices().coeff(i)) = b.col(i);
VERIFY_IS_APPROX(m1, c); VERIFY_IS_APPROX(m1, c);
MatrixType m2 = MatrixType::Random(cols,cols2); MatrixType m2 = MatrixType::Random(cols,cols2);
@ -80,15 +72,8 @@ template<typename MatrixType, int Cols2> void qr_fixedsize()
VERIFY(!qr.isInvertible()); VERIFY(!qr.isInvertible());
VERIFY(!qr.isSurjective()); VERIFY(!qr.isSurjective());
Matrix<Scalar,Rows,Cols> r = qr.matrixQR(); Matrix<Scalar,Rows,Cols> r = qr.matrixQR().template triangularView<UpperTriangular>();
// FIXME need better way to construct trapezoid Matrix<Scalar,Rows,Cols> c = qr.matrixQ() * r * qr.colsPermutation().inverse();
for(int i = 0; i < Rows; i++) for(int j = 0; j < Cols; j++) if(i>j) r(i,j) = Scalar(0);
Matrix<Scalar,Rows,Cols> b = qr.matrixQ() * r;
Matrix<Scalar,Rows,Cols> c = MatrixType::Zero(Rows,Cols);
for(int i = 0; i < Cols; ++i) c.col(qr.colsPermutation().indices().coeff(i)) = b.col(i);
VERIFY_IS_APPROX(m1, c); VERIFY_IS_APPROX(m1, c);
Matrix<Scalar,Cols,Cols2> m2 = Matrix<Scalar,Cols,Cols2>::Random(Cols,Cols2); Matrix<Scalar,Cols,Cols2> m2 = Matrix<Scalar,Cols,Cols2>::Random(Cols,Cols2);
@ -118,7 +103,7 @@ template<typename MatrixType> void qr_invertible()
ColPivHouseholderQR<MatrixType> qr(m1); ColPivHouseholderQR<MatrixType> qr(m1);
m3 = MatrixType::Random(size,size); m3 = MatrixType::Random(size,size);
m2 = qr.solve(m3); m2 = qr.solve(m3);
VERIFY_IS_APPROX(m3, m1*m2); //VERIFY_IS_APPROX(m3, m1*m2);
// now construct a matrix with prescribed determinant // now construct a matrix with prescribed determinant
m1.setZero(); m1.setZero();
@ -150,7 +135,7 @@ template<typename MatrixType> void qr_verify_assert()
void test_qr_colpivoting() void test_qr_colpivoting()
{ {
for(int i = 0; i < 1; i++) { for(int i = 0; i < g_repeat; i++) {
CALL_SUBTEST_1( qr<MatrixXf>() ); CALL_SUBTEST_1( qr<MatrixXf>() );
CALL_SUBTEST_2( qr<MatrixXd>() ); CALL_SUBTEST_2( qr<MatrixXd>() );
CALL_SUBTEST_3( qr<MatrixXcd>() ); CALL_SUBTEST_3( qr<MatrixXcd>() );

View File

@ -188,7 +188,7 @@ template<typename _Scalar> class AlignedVector3
} }
template<typename Derived> template<typename Derived>
inline bool isApprox(const MatrixBase<Derived>& other, RealScalar eps=precision<Scalar>()) const inline bool isApprox(const MatrixBase<Derived>& other, RealScalar eps=dummy_precision<Scalar>()) const
{ {
return m_coeffs.template start<3>().isApprox(other,eps); return m_coeffs.template start<3>().isApprox(other,eps);
} }

View File

@ -272,7 +272,7 @@ void MatrixExponential<MatrixType>::computeUV(float)
} else { } else {
const float maxnorm = 3.925724783138660f; const float maxnorm = 3.925724783138660f;
m_squarings = std::max(0, (int)ceil(log2(m_l1norm / maxnorm))); m_squarings = std::max(0, (int)ceil(log2(m_l1norm / maxnorm)));
MatrixType A = *m_M / std::pow(Scalar(2), m_squarings); MatrixType A = *m_M / std::pow(Scalar(2), Scalar(m_squarings));
pade7(A); pade7(A);
} }
} }
@ -291,7 +291,7 @@ void MatrixExponential<MatrixType>::computeUV(double)
} else { } else {
const double maxnorm = 5.371920351148152; const double maxnorm = 5.371920351148152;
m_squarings = std::max(0, (int)ceil(log2(m_l1norm / maxnorm))); m_squarings = std::max(0, (int)ceil(log2(m_l1norm / maxnorm)));
MatrixType A = *m_M / std::pow(Scalar(2), m_squarings); MatrixType A = *m_M / std::pow(Scalar(2), Scalar(m_squarings));
pade13(A); pade13(A);
} }
} }

View File

@ -46,7 +46,7 @@ public:
* flags \a flags. */ * flags \a flags. */
SkylineInplaceLU(MatrixType& matrix, int flags = 0) SkylineInplaceLU(MatrixType& matrix, int flags = 0)
: /*m_matrix(matrix.rows(), matrix.cols()),*/ m_flags(flags), m_status(0), m_lu(matrix) { : /*m_matrix(matrix.rows(), matrix.cols()),*/ m_flags(flags), m_status(0), m_lu(matrix) {
m_precision = RealScalar(0.1) * Eigen::precision<RealScalar > (); m_precision = RealScalar(0.1) * Eigen::dummy_precision<RealScalar > ();
m_lu.IsRowMajor ? computeRowMajor() : compute(); m_lu.IsRowMajor ? computeRowMajor() : compute();
} }

View File

@ -589,7 +589,7 @@ public:
m_data.squeeze(); m_data.squeeze();
} }
void prune(Scalar reference, RealScalar epsilon = precision<RealScalar > ()) { void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar > ()) {
//TODO //TODO
} }

View File

@ -206,7 +206,7 @@ public:
memset(m_lowerProfile, 0, m_diagSize * sizeof (int)); memset(m_lowerProfile, 0, m_diagSize * sizeof (int));
} }
void prune(Scalar reference, RealScalar epsilon = precision<RealScalar>()) { void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar>()) {
//TODO //TODO
} }

View File

@ -1,8 +1,3 @@
include(EigenTesting)
enable_testing()
find_package(Adolc) find_package(Adolc)
include_directories(../../test) include_directories(../../test)

View File

@ -47,7 +47,11 @@ complex<long double> promote(long double x) { return complex<long double>( x);
cerr <<"idx\ttruth\t\tvalue\t|dif|=\n"; cerr <<"idx\ttruth\t\tvalue\t|dif|=\n";
for (size_t k0=0;k0<size_t(fftbuf.size());++k0) { for (size_t k0=0;k0<size_t(fftbuf.size());++k0) {
complex<long double> acc = 0; complex<long double> acc = 0;
#ifdef _GNU_SOURCE
long double phinc = -2.*k0* M_PIl / timebuf.size(); long double phinc = -2.*k0* M_PIl / timebuf.size();
#else
long double phinc = -2.*k0* M_PI / timebuf.size();
#endif
for (size_t k1=0;k1<size_t(timebuf.size());++k1) { for (size_t k1=0;k1<size_t(timebuf.size());++k1) {
acc += promote( timebuf[k1] ) * exp( complex<long double>(0,k1*phinc) ); acc += promote( timebuf[k1] ) * exp( complex<long double>(0,k1*phinc) );
} }