diff --git a/src/Core/Dot.h b/src/Core/Dot.h index 88b180c4c..158b852eb 100644 --- a/src/Core/Dot.h +++ b/src/Core/Dot.h @@ -84,7 +84,7 @@ typename NumTraits::Real Object::norm2() const template typename NumTraits::Real Object::norm() const { - return Sqrt(norm2()); + return NumTraits::Real>::sqrt(norm2()); } template diff --git a/src/Core/Fuzzy.h b/src/Core/Fuzzy.h index eadffcff7..0334e6d96 100644 --- a/src/Core/Fuzzy.h +++ b/src/Core/Fuzzy.h @@ -34,7 +34,7 @@ bool Object::isApprox(const OtherDerived& other) const { return((*this - other).norm2() <= std::min(norm2(), other.norm2()) - * Abs2(NumTraits::epsilon())); + * NumTraits::epsilon2()); } else { @@ -50,7 +50,7 @@ bool Object::isNegligble(const Scalar& other) const { if(IsVector) { - return(norm2() <= Abs2(other) * Abs2(NumTraits::epsilon())); + return(norm2() <= NumTraits::abs2(other) * NumTraits::epsilon2()); } else { @@ -67,7 +67,7 @@ bool Object::isNegligble(const Object& ot { if(IsVector) { - return(norm2() <= other.norm2() * Abs2(NumTraits::epsilon())); + return(norm2() <= other.norm2() * NumTraits::epsilon2()); } else { diff --git a/src/Core/Numeric.h b/src/Core/Numeric.h index 1bd48044b..8437dfb2e 100644 --- a/src/Core/Numeric.h +++ b/src/Core/Numeric.h @@ -39,6 +39,7 @@ template<> struct NumTraits static const bool HasFloatingPoint = false; static int epsilon() { return 0; } + static int epsilon2() { return 0; } static int real(const int& x) { return x; } static int imag(const int& x) { EI_UNUSED(x); return 0; } static int conj(const int& x) { return x; } @@ -47,10 +48,23 @@ template<> struct NumTraits static int abs2(const int& x) { return x*x; } static int rand() { - // "rand()%21" would be bad. always use the high-order bits, not the low-order bits. - // note: here (gcc 4.1) static_cast seems to round the nearest int. - // I don't know if that's part of the standard. - return -10 + static_cast(std::rand() / ((RAND_MAX + 1.0)/20.0)); + // "rand() % n" is bad, they say, because the low-order bits are not random enough. + // However here, 21 is odd, so rand() % 21 uses the high-order bits + // as well, so there's no problem. + return (std::rand() % 21) - 10; + } + static bool negligible(const int& a, const int& b) + { + EI_UNUSED(b); + return(a == 0); + } + static bool approx(const int& a, const int& b) + { + return(a == b); + } + static bool lessThanOrApprox(const int& a, const int& b) + { + return(a <= b); } }; @@ -64,6 +78,7 @@ template<> struct NumTraits static const bool HasFloatingPoint = true; static float epsilon() { return 1e-5f; } + static float epsilon2() { return epsilon() * epsilon(); } static float real(const float& x) { return x; } static float imag(const float& x) { EI_UNUSED(x); return 0; } static float conj(const float& x) { return x; } @@ -74,6 +89,18 @@ template<> struct NumTraits { return std::rand() / (RAND_MAX/20.0f) - 10.0f; } + static bool negligible(const float& a, const float& b) + { + return(abs(a) <= abs(b) * epsilon()); + } + static bool approx(const float& a, const float& b) + { + return(abs(a - b) <= std::min(abs(a), abs(b)) * epsilon()); + } + static bool lessThanOrApprox(const float& a, const float& b) + { + return(a <= b || approx(a, b)); + } }; template<> struct NumTraits @@ -86,6 +113,7 @@ template<> struct NumTraits static const bool HasFloatingPoint = true; static double epsilon() { return 1e-11; } + static double epsilon2() { return epsilon() * epsilon(); } static double real(const double& x) { return x; } static double imag(const double& x) { EI_UNUSED(x); return 0; } static double conj(const double& x) { return x; } @@ -96,6 +124,18 @@ template<> struct NumTraits { return std::rand() / (RAND_MAX/20.0) - 10.0; } + static bool negligible(const double& a, const double& b) + { + return(abs(a) <= abs(b) * epsilon()); + } + static bool approx(const double& a, const double& b) + { + return(abs(a - b) <= std::min(abs(a), abs(b)) * epsilon()); + } + static bool lessThanOrApprox(const double& a, const double& b) + { + return(a <= b || approx(a, b)); + } }; template struct NumTraits > @@ -109,6 +149,7 @@ template struct NumTraits > static const bool HasFloatingPoint = NumTraits::HasFloatingPoint; static Real epsilon() { return NumTraits::epsilon(); } + static Real epsilon2() { return epsilon() * epsilon(); } static Real real(const Complex& x) { return std::real(x); } static Real imag(const Complex& x) { return std::imag(x); } static Complex conj(const Complex& x) { return std::conj(x); } @@ -122,48 +163,16 @@ template struct NumTraits > { return Complex(NumTraits::rand(), NumTraits::rand()); } + static bool negligible(const Complex& a, const Complex& b) + { + return(abs2(a) <= abs2(b) * epsilon2()); + } + static bool approx(const Complex& a, const Complex& b) + { + return(NumTraits::approx(std::real(a), std::real(b)) + && NumTraits::approx(std::imag(a), std::imag(b))); + } + // lessThanOrApprox wouldn't make sense for complex numbers }; -template typename NumTraits::Real Real(const T& x) -{ return NumTraits::real(x); } - -template typename NumTraits::Real Imag(const T& x) -{ return NumTraits::imag(x); } - -template T Conj(const T& x) -{ return NumTraits::conj(x); } - -template typename NumTraits::FloatingPoint Sqrt(const T& x) -{ return NumTraits::sqrt(x); } - -template typename NumTraits::RealFloatingPoint Abs(const T& x) -{ return NumTraits::abs(x); } - -template typename NumTraits::Real Abs2(const T& x) -{ return NumTraits::abs2(x); } - -template T Rand() -{ return NumTraits::rand(); } - -template bool Negligible(const T& a, const T& b) -{ - return(Abs(a) <= Abs(b) * NumTraits::epsilon()); -} - -template bool Approx(const T& a, const T& b) -{ - if(NumTraits::IsFloat) - return(Abs(a - b) <= std::min(Abs(a), Abs(b)) * NumTraits::epsilon()); - else - return(a == b); -} - -template bool LessThanOrApprox(const T& a, const T& b) -{ - if(NumTraits::IsFloat) - return(a < b || Approx(a, b)); - else - return(a <= b); -} - #endif // EI_NUMERIC_H diff --git a/src/Core/Random.h b/src/Core/Random.h index 9827f5585..a90d69d5b 100644 --- a/src/Core/Random.h +++ b/src/Core/Random.h @@ -51,7 +51,7 @@ template class Random { EI_UNUSED(row); EI_UNUSED(col); - return Rand(); + return NumTraits::rand(); } protected: diff --git a/test/main.cpp b/test/main.cpp index 8610cf363..f5c758961 100644 --- a/test/main.cpp +++ b/test/main.cpp @@ -25,7 +25,7 @@ #include "main.h" -genTest::genTest() +EigenTest::EigenTest() { unsigned int t = (unsigned int) time( NULL ); qDebug() << "Initializing random number generator with seed" diff --git a/test/main.h b/test/main.h index 47ea1b928..85968134c 100644 --- a/test/main.h +++ b/test/main.h @@ -29,7 +29,7 @@ #include #include "../src/Core.h" -USING_EIGEN_DATA_TYPES +using namespace Eigen; #include #include @@ -62,6 +62,9 @@ template bool TestNegligible(const T& a, const T& b) return(Abs(a) <= Abs(b) * TestEpsilon()); } +//template +//bool TestNegligible + template bool TestApprox(const T& a, const T& b) { if(Eigen::NumTraits::IsFloat)