aboutsummaryrefslogtreecommitdiffhomepage
path: root/Eigen/src/Core/Fuzzy.h
diff options
context:
space:
mode:
authorGravatar Benoit Jacob <jacob.benoit.1@gmail.com>2010-10-25 10:15:22 -0400
committerGravatar Benoit Jacob <jacob.benoit.1@gmail.com>2010-10-25 10:15:22 -0400
commit4716040703be1ee906439385d20475dcddad5ce3 (patch)
tree8efd3cf3007d8360e66f38e2d280127cbb70daa6 /Eigen/src/Core/Fuzzy.h
parentca85a1f6c5fc33ac382aa2d7ba2da63d55d3223e (diff)
bug #86 : use internal:: namespace instead of ei_ prefix
Diffstat (limited to 'Eigen/src/Core/Fuzzy.h')
-rw-r--r--Eigen/src/Core/Fuzzy.h147
1 files changed, 6 insertions, 141 deletions
diff --git a/Eigen/src/Core/Fuzzy.h b/Eigen/src/Core/Fuzzy.h
index 4a3b5e4dc..0a3620df9 100644
--- a/Eigen/src/Core/Fuzzy.h
+++ b/Eigen/src/Core/Fuzzy.h
@@ -28,8 +28,6 @@
// TODO support small integer types properly i.e. do exact compare on coeffs --- taking a HS norm is guaranteed to cause integer overflow.
-#ifndef EIGEN_LEGACY_COMPARES
-
/** \returns \c true if \c *this is approximately equal to \a other, within the precision
* determined by \a prec.
*
@@ -42,10 +40,10 @@
* \note Because of the multiplicativeness of this comparison, one can't use this function
* to check whether \c *this is approximately equal to the zero matrix or vector.
* Indeed, \c isApprox(zero) returns false unless \c *this itself is exactly the zero matrix
- * or vector. If you want to test whether \c *this is zero, use ei_isMuchSmallerThan(const
+ * or vector. If you want to test whether \c *this is zero, use internal::isMuchSmallerThan(const
* RealScalar&, RealScalar) instead.
*
- * \sa ei_isMuchSmallerThan(const RealScalar&, RealScalar) const
+ * \sa internal::isMuchSmallerThan(const RealScalar&, RealScalar) const
*/
template<typename Derived>
template<typename OtherDerived>
@@ -54,10 +52,10 @@ bool DenseBase<Derived>::isApprox(
RealScalar prec
) const
{
- const typename ei_nested<Derived,2>::type nested(derived());
- const typename ei_nested<OtherDerived,2>::type otherNested(other.derived());
-// std::cerr << typeid(Derived).name() << " => " << typeid(typename ei_nested<Derived,2>::type).name() << "\n";
-// std::cerr << typeid(OtherDerived).name() << " => " << typeid(typename ei_nested<OtherDerived,2>::type).name() << "\n";
+ const typename internal::nested<Derived,2>::type nested(derived());
+ const typename internal::nested<OtherDerived,2>::type otherNested(other.derived());
+// std::cerr << typeid(Derived).name() << " => " << typeid(typename internal::nested<Derived,2>::type).name() << "\n";
+// std::cerr << typeid(OtherDerived).name() << " => " << typeid(typename internal::nested<OtherDerived,2>::type).name() << "\n";
// return false;
return (nested - otherNested).cwiseAbs2().sum() <= prec * prec * std::min(nested.cwiseAbs2().sum(), otherNested.cwiseAbs2().sum());
}
@@ -104,137 +102,4 @@ bool DenseBase<Derived>::isMuchSmallerThan(
return derived().cwiseAbs2().sum() <= prec * prec * other.derived().cwiseAbs2().sum();
}
-#else
-
-template<typename Derived, typename OtherDerived=Derived, bool IsVector=Derived::IsVectorAtCompileTime>
-struct ei_fuzzy_selector;
-
-/** \returns \c true if \c *this is approximately equal to \a other, within the precision
- * determined by \a prec.
- *
- * \note The fuzzy compares are done multiplicatively. Two vectors \f$ v \f$ and \f$ w \f$
- * are considered to be approximately equal within precision \f$ p \f$ if
- * \f[ \Vert v - w \Vert \leqslant p\,\min(\Vert v\Vert, \Vert w\Vert). \f]
- * For matrices, the comparison is done on all columns.
- *
- * \note Because of the multiplicativeness of this comparison, one can't use this function
- * to check whether \c *this is approximately equal to the zero matrix or vector.
- * Indeed, \c isApprox(zero) returns false unless \c *this itself is exactly the zero matrix
- * or vector. If you want to test whether \c *this is zero, use ei_isMuchSmallerThan(const
- * RealScalar&, RealScalar) instead.
- *
- * \sa ei_isMuchSmallerThan(const RealScalar&, RealScalar) const
- */
-template<typename Derived>
-template<typename OtherDerived>
-bool DenseBase<Derived>::isApprox(
- const DenseBase<OtherDerived>& other,
- RealScalar prec
-) const
-{
- return ei_fuzzy_selector<Derived,OtherDerived>::isApprox(derived(), other.derived(), prec);
-}
-
-/** \returns \c true if the norm of \c *this is much smaller than \a other,
- * within the precision determined by \a prec.
- *
- * \note The fuzzy compares are done multiplicatively. A vector \f$ v \f$ is
- * considered to be much smaller than \f$ x \f$ within precision \f$ p \f$ if
- * \f[ \Vert v \Vert \leqslant p\,\vert x\vert. \f]
- * For matrices, the comparison is done on all columns.
- *
- * \sa isApprox(), isMuchSmallerThan(const DenseBase<OtherDerived>&, RealScalar) const
- */
-template<typename Derived>
-bool DenseBase<Derived>::isMuchSmallerThan(
- const typename NumTraits<Scalar>::Real& other,
- RealScalar prec
-) const
-{
- return ei_fuzzy_selector<Derived>::isMuchSmallerThan(derived(), other, prec);
-}
-
-/** \returns \c true if the norm of \c *this is much smaller than the norm of \a other,
- * within the precision determined by \a prec.
- *
- * \note The fuzzy compares are done multiplicatively. A vector \f$ v \f$ is
- * considered to be much smaller than a vector \f$ w \f$ within precision \f$ p \f$ if
- * \f[ \Vert v \Vert \leqslant p\,\Vert w\Vert. \f]
- * For matrices, the comparison is done on all columns.
- *
- * \sa isApprox(), isMuchSmallerThan(const RealScalar&, RealScalar) const
- */
-template<typename Derived>
-template<typename OtherDerived>
-bool DenseBase<Derived>::isMuchSmallerThan(
- const DenseBase<OtherDerived>& other,
- RealScalar prec
-) const
-{
- return ei_fuzzy_selector<Derived,OtherDerived>::isMuchSmallerThan(derived(), other.derived(), prec);
-}
-
-
-template<typename Derived, typename OtherDerived>
-struct ei_fuzzy_selector<Derived,OtherDerived,true>
-{
- typedef typename Derived::RealScalar RealScalar;
- static bool isApprox(const Derived& self, const OtherDerived& other, RealScalar prec)
- {
- EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(Derived,OtherDerived)
- ei_assert(self.size() == other.size());
- return((self - other).squaredNorm() <= std::min(self.squaredNorm(), other.squaredNorm()) * prec * prec);
- }
- static bool isMuchSmallerThan(const Derived& self, const RealScalar& other, RealScalar prec)
- {
- return(self.squaredNorm() <= ei_abs2(other * prec));
- }
- static bool isMuchSmallerThan(const Derived& self, const OtherDerived& other, RealScalar prec)
- {
- EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(Derived,OtherDerived)
- ei_assert(self.size() == other.size());
- return(self.squaredNorm() <= other.squaredNorm() * prec * prec);
- }
-};
-
-template<typename Derived, typename OtherDerived>
-struct ei_fuzzy_selector<Derived,OtherDerived,false>
-{
- typedef typename Derived::RealScalar RealScalar;
- typedef typename Derived::Index Index;
- static bool isApprox(const Derived& self, const OtherDerived& other, RealScalar prec)
- {
- EIGEN_STATIC_ASSERT_SAME_MATRIX_SIZE(Derived,OtherDerived)
- ei_assert(self.rows() == other.rows() && self.cols() == other.cols());
- typename Derived::Nested nested(self);
- typename OtherDerived::Nested otherNested(other);
- for(Index i = 0; i < self.cols(); ++i)
- if((nested.col(i) - otherNested.col(i)).squaredNorm()
- > std::min(nested.col(i).squaredNorm(), otherNested.col(i).squaredNorm()) * prec * prec)
- return false;
- return true;
- }
- static bool isMuchSmallerThan(const Derived& self, const RealScalar& other, RealScalar prec)
- {
- typename Derived::Nested nested(self);
- for(Index i = 0; i < self.cols(); ++i)
- if(nested.col(i).squaredNorm() > ei_abs2(other * prec))
- return false;
- return true;
- }
- static bool isMuchSmallerThan(const Derived& self, const OtherDerived& other, RealScalar prec)
- {
- EIGEN_STATIC_ASSERT_SAME_MATRIX_SIZE(Derived,OtherDerived)
- ei_assert(self.rows() == other.rows() && self.cols() == other.cols());
- typename Derived::Nested nested(self);
- typename OtherDerived::Nested otherNested(other);
- for(Index i = 0; i < self.cols(); ++i)
- if(nested.col(i).squaredNorm() > otherNested.col(i).squaredNorm() * prec * prec)
- return false;
- return true;
- }
-};
-
-#endif
-
#endif // EIGEN_FUZZY_H