MathFunctions.h

00001 // This file is part of Eigen, a lightweight C++ template library
00002 // for linear algebra.
00003 //
00004 // Copyright (C) 2006-2010 Benoit Jacob <jacob.benoit.1@gmail.com>
00005 //
00006 // This Source Code Form is subject to the terms of the Mozilla
00007 // Public License v. 2.0. If a copy of the MPL was not distributed
00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
00009 
00010 #ifndef EIGEN_MATHFUNCTIONS_H
00011 #define EIGEN_MATHFUNCTIONS_H
00012 
00013 namespace Eigen {
00014 
00015 namespace internal {
00016 
00017 /** \internal \struct global_math_functions_filtering_base
00018   *
00019   * What it does:
00020   * Defines a typedef 'type' as follows:
00021   * - if type T has a member typedef Eigen_BaseClassForSpecializationOfGlobalMathFuncImpl, then
00022   *   global_math_functions_filtering_base<T>::type is a typedef for it.
00023   * - otherwise, global_math_functions_filtering_base<T>::type is a typedef for T.
00024   *
00025   * How it's used:
00026   * To allow to defined the global math functions (like sin...) in certain cases, like the Array expressions.
00027   * When you do sin(array1+array2), the object array1+array2 has a complicated expression type, all what you want to know
00028   * is that it inherits ArrayBase. So we implement a partial specialization of sin_impl for ArrayBase<Derived>.
00029   * So we must make sure to use sin_impl<ArrayBase<Derived> > and not sin_impl<Derived>, otherwise our partial specialization
00030   * won't be used. How does sin know that? That's exactly what global_math_functions_filtering_base tells it.
00031   *
00032   * How it's implemented:
00033   * SFINAE in the style of enable_if. Highly susceptible of breaking compilers. With GCC, it sure does work, but if you replace
00034   * the typename dummy by an integer template parameter, it doesn't work anymore!
00035   */
00036 
00037 template<typename T, typename dummy = void>
00038 struct global_math_functions_filtering_base
00039 {
00040   typedef T type;
00041 };
00042 
00043 template<typename T> struct always_void { typedef void type; };
00044 
00045 template<typename T>
00046 struct global_math_functions_filtering_base
00047   <T,
00048    typename always_void<typename T::Eigen_BaseClassForSpecializationOfGlobalMathFuncImpl>::type
00049   >
00050 {
00051   typedef typename T::Eigen_BaseClassForSpecializationOfGlobalMathFuncImpl type;
00052 };
00053 
00054 #define EIGEN_MATHFUNC_IMPL(func, scalar) Eigen::internal::func##_impl<typename Eigen::internal::global_math_functions_filtering_base<scalar>::type>
00055 #define EIGEN_MATHFUNC_RETVAL(func, scalar) typename Eigen::internal::func##_retval<typename Eigen::internal::global_math_functions_filtering_base<scalar>::type>::type
00056 
00057 /****************************************************************************
00058 * Implementation of real                                                 *
00059 ****************************************************************************/
00060 
00061 template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
00062 struct real_default_impl
00063 {
00064   typedef typename NumTraits<Scalar>::Real RealScalar;
00065   static inline RealScalar run(const Scalar& x)
00066   {
00067     return x;
00068   }
00069 };
00070 
00071 template<typename Scalar>
00072 struct real_default_impl<Scalar,true>
00073 {
00074   typedef typename NumTraits<Scalar>::Real RealScalar;
00075   static inline RealScalar run(const Scalar& x)
00076   {
00077     using std::real;
00078     return real(x);
00079   }
00080 };
00081 
00082 template<typename Scalar> struct real_impl : real_default_impl<Scalar> {};
00083 
00084 template<typename Scalar>
00085 struct real_retval
00086 {
00087   typedef typename NumTraits<Scalar>::Real type;
00088 };
00089 
00090 
00091 /****************************************************************************
00092 * Implementation of imag                                                 *
00093 ****************************************************************************/
00094 
00095 template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
00096 struct imag_default_impl
00097 {
00098   typedef typename NumTraits<Scalar>::Real RealScalar;
00099   static inline RealScalar run(const Scalar&)
00100   {
00101     return RealScalar(0);
00102   }
00103 };
00104 
00105 template<typename Scalar>
00106 struct imag_default_impl<Scalar,true>
00107 {
00108   typedef typename NumTraits<Scalar>::Real RealScalar;
00109   static inline RealScalar run(const Scalar& x)
00110   {
00111     using std::imag;
00112     return imag(x);
00113   }
00114 };
00115 
00116 template<typename Scalar> struct imag_impl : imag_default_impl<Scalar> {};
00117 
00118 template<typename Scalar>
00119 struct imag_retval
00120 {
00121   typedef typename NumTraits<Scalar>::Real type;
00122 };
00123 
00124 /****************************************************************************
00125 * Implementation of real_ref                                             *
00126 ****************************************************************************/
00127 
00128 template<typename Scalar>
00129 struct real_ref_impl
00130 {
00131   typedef typename NumTraits<Scalar>::Real RealScalar;
00132   static inline RealScalar& run(Scalar& x)
00133   {
00134     return reinterpret_cast<RealScalar*>(&x)[0];
00135   }
00136   static inline const RealScalar& run(const Scalar& x)
00137   {
00138     return reinterpret_cast<const RealScalar*>(&x)[0];
00139   }
00140 };
00141 
00142 template<typename Scalar>
00143 struct real_ref_retval
00144 {
00145   typedef typename NumTraits<Scalar>::Real & type;
00146 };
00147 
00148 /****************************************************************************
00149 * Implementation of imag_ref                                             *
00150 ****************************************************************************/
00151 
00152 template<typename Scalar, bool IsComplex>
00153 struct imag_ref_default_impl
00154 {
00155   typedef typename NumTraits<Scalar>::Real RealScalar;
00156   static inline RealScalar& run(Scalar& x)
00157   {
00158     return reinterpret_cast<RealScalar*>(&x)[1];
00159   }
00160   static inline const RealScalar& run(const Scalar& x)
00161   {
00162     return reinterpret_cast<RealScalar*>(&x)[1];
00163   }
00164 };
00165 
00166 template<typename Scalar>
00167 struct imag_ref_default_impl<Scalar, false>
00168 {
00169   static inline Scalar run(Scalar&)
00170   {
00171     return Scalar(0);
00172   }
00173   static inline const Scalar run(const Scalar&)
00174   {
00175     return Scalar(0);
00176   }
00177 };
00178 
00179 template<typename Scalar>
00180 struct imag_ref_impl : imag_ref_default_impl<Scalar, NumTraits<Scalar>::IsComplex> {};
00181 
00182 template<typename Scalar>
00183 struct imag_ref_retval
00184 {
00185   typedef typename NumTraits<Scalar>::Real & type;
00186 };
00187 
00188 /****************************************************************************
00189 * Implementation of conj                                                 *
00190 ****************************************************************************/
00191 
00192 template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
00193 struct conj_impl
00194 {
00195   static inline Scalar run(const Scalar& x)
00196   {
00197     return x;
00198   }
00199 };
00200 
00201 template<typename Scalar>
00202 struct conj_impl<Scalar,true>
00203 {
00204   static inline Scalar run(const Scalar& x)
00205   {
00206     using std::conj;
00207     return conj(x);
00208   }
00209 };
00210 
00211 template<typename Scalar>
00212 struct conj_retval
00213 {
00214   typedef Scalar type;
00215 };
00216 
00217 /****************************************************************************
00218 * Implementation of abs2                                                 *
00219 ****************************************************************************/
00220 
00221 template<typename Scalar>
00222 struct abs2_impl
00223 {
00224   typedef typename NumTraits<Scalar>::Real RealScalar;
00225   static inline RealScalar run(const Scalar& x)
00226   {
00227     return x*x;
00228   }
00229 };
00230 
00231 template<typename RealScalar>
00232 struct abs2_impl<std::complex<RealScalar> >
00233 {
00234   static inline RealScalar run(const std::complex<RealScalar>& x)
00235   {
00236     return real(x)*real(x) + imag(x)*imag(x);
00237   }
00238 };
00239 
00240 template<typename Scalar>
00241 struct abs2_retval
00242 {
00243   typedef typename NumTraits<Scalar>::Real type;
00244 };
00245 
00246 /****************************************************************************
00247 * Implementation of norm1                                                *
00248 ****************************************************************************/
00249 
00250 template<typename Scalar, bool IsComplex>
00251 struct norm1_default_impl
00252 {
00253   typedef typename NumTraits<Scalar>::Real RealScalar;
00254   static inline RealScalar run(const Scalar& x)
00255   {
00256     using std::abs;
00257     return abs(real(x)) + abs(imag(x));
00258   }
00259 };
00260 
00261 template<typename Scalar>
00262 struct norm1_default_impl<Scalar, false>
00263 {
00264   static inline Scalar run(const Scalar& x)
00265   {
00266     using std::abs;
00267     return abs(x);
00268   }
00269 };
00270 
00271 template<typename Scalar>
00272 struct norm1_impl : norm1_default_impl<Scalar, NumTraits<Scalar>::IsComplex> {};
00273 
00274 template<typename Scalar>
00275 struct norm1_retval
00276 {
00277   typedef typename NumTraits<Scalar>::Real type;
00278 };
00279 
00280 /****************************************************************************
00281 * Implementation of hypot                                                *
00282 ****************************************************************************/
00283 
00284 template<typename Scalar>
00285 struct hypot_impl
00286 {
00287   typedef typename NumTraits<Scalar>::Real RealScalar;
00288   static inline RealScalar run(const Scalar& x, const Scalar& y)
00289   {
00290     using std::max;
00291     using std::min;
00292     using std::abs;
00293     using std::sqrt;
00294     RealScalar _x = abs(x);
00295     RealScalar _y = abs(y);
00296     RealScalar p = (max)(_x, _y);
00297     if(p==RealScalar(0)) return RealScalar(0);
00298     RealScalar q = (min)(_x, _y);
00299     RealScalar qp = q/p;
00300     return p * sqrt(RealScalar(1) + qp*qp);
00301   }
00302 };
00303 
00304 template<typename Scalar>
00305 struct hypot_retval
00306 {
00307   typedef typename NumTraits<Scalar>::Real type;
00308 };
00309 
00310 /****************************************************************************
00311 * Implementation of cast                                                 *
00312 ****************************************************************************/
00313 
00314 template<typename OldType, typename NewType>
00315 struct cast_impl
00316 {
00317   static inline NewType run(const OldType& x)
00318   {
00319     return static_cast<NewType>(x);
00320   }
00321 };
00322 
00323 // here, for once, we're plainly returning NewType: we don't want cast to do weird things.
00324 
00325 template<typename OldType, typename NewType>
00326 inline NewType cast(const OldType& x)
00327 {
00328   return cast_impl<OldType, NewType>::run(x);
00329 }
00330 
00331 /****************************************************************************
00332 * Implementation of atanh2                                                *
00333 ****************************************************************************/
00334 
00335 template<typename Scalar, bool IsInteger>
00336 struct atanh2_default_impl
00337 {
00338   typedef Scalar retval;
00339   typedef typename NumTraits<Scalar>::Real RealScalar;
00340   static inline Scalar run(const Scalar& x, const Scalar& y)
00341   {
00342     using std::abs;
00343     using std::log;
00344     using std::sqrt;
00345     Scalar z = x / y;
00346     if (y == Scalar(0) || abs(z) > sqrt(NumTraits<RealScalar>::epsilon()))
00347       return RealScalar(0.5) * log((y + x) / (y - x));
00348     else
00349       return z + z*z*z / RealScalar(3);
00350   }
00351 };
00352 
00353 template<typename Scalar>
00354 struct atanh2_default_impl<Scalar, true>
00355 {
00356   static inline Scalar run(const Scalar&, const Scalar&)
00357   {
00358     EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar)
00359     return Scalar(0);
00360   }
00361 };
00362 
00363 template<typename Scalar>
00364 struct atanh2_impl : atanh2_default_impl<Scalar, NumTraits<Scalar>::IsInteger> {};
00365 
00366 template<typename Scalar>
00367 struct atanh2_retval
00368 {
00369   typedef Scalar type;
00370 };
00371 
00372 /****************************************************************************
00373 * Implementation of pow                                                  *
00374 ****************************************************************************/
00375 
00376 template<typename Scalar, bool IsInteger>
00377 struct pow_default_impl
00378 {
00379   typedef Scalar retval;
00380   static inline Scalar run(const Scalar& x, const Scalar& y)
00381   {
00382     using std::pow;
00383     return pow(x, y);
00384   }
00385 };
00386 
00387 template<typename Scalar>
00388 struct pow_default_impl<Scalar, true>
00389 {
00390   static inline Scalar run(Scalar x, Scalar y)
00391   {
00392     Scalar res(1);
00393     eigen_assert(!NumTraits<Scalar>::IsSigned || y >= 0);
00394     if(y & 1) res *= x;
00395     y >>= 1;
00396     while(y)
00397     {
00398       x *= x;
00399       if(y&1) res *= x;
00400       y >>= 1;
00401     }
00402     return res;
00403   }
00404 };
00405 
00406 template<typename Scalar>
00407 struct pow_impl : pow_default_impl<Scalar, NumTraits<Scalar>::IsInteger> {};
00408 
00409 template<typename Scalar>
00410 struct pow_retval
00411 {
00412   typedef Scalar type;
00413 };
00414 
00415 /****************************************************************************
00416 * Implementation of random                                               *
00417 ****************************************************************************/
00418 
00419 template<typename Scalar,
00420          bool IsComplex,
00421          bool IsInteger>
00422 struct random_default_impl {};
00423 
00424 template<typename Scalar>
00425 struct random_impl : random_default_impl<Scalar, NumTraits<Scalar>::IsComplex, NumTraits<Scalar>::IsInteger> {};
00426 
00427 template<typename Scalar>
00428 struct random_retval
00429 {
00430   typedef Scalar type;
00431 };
00432 
00433 template<typename Scalar> inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random(const Scalar& x, const Scalar& y);
00434 template<typename Scalar> inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random();
00435 
00436 template<typename Scalar>
00437 struct random_default_impl<Scalar, false, false>
00438 {
00439   static inline Scalar run(const Scalar& x, const Scalar& y)
00440   {
00441     return x + (y-x) * Scalar(std::rand()) / Scalar(RAND_MAX);
00442   }
00443   static inline Scalar run()
00444   {
00445     return run(Scalar(NumTraits<Scalar>::IsSigned ? -1 : 0), Scalar(1));
00446   }
00447 };
00448 
00449 enum {
00450   floor_log2_terminate,
00451   floor_log2_move_up,
00452   floor_log2_move_down,
00453   floor_log2_bogus
00454 };
00455 
00456 template<unsigned int n, int lower, int upper> struct floor_log2_selector
00457 {
00458   enum { middle = (lower + upper) / 2,
00459          value = (upper <= lower + 1) ? int(floor_log2_terminate)
00460                : (n < (1 << middle)) ? int(floor_log2_move_down)
00461                : (n==0) ? int(floor_log2_bogus)
00462                : int(floor_log2_move_up)
00463   };
00464 };
00465 
00466 template<unsigned int n,
00467          int lower = 0,
00468          int upper = sizeof(unsigned int) * CHAR_BIT - 1,
00469          int selector = floor_log2_selector<n, lower, upper>::value>
00470 struct floor_log2 {};
00471 
00472 template<unsigned int n, int lower, int upper>
00473 struct floor_log2<n, lower, upper, floor_log2_move_down>
00474 {
00475   enum { value = floor_log2<n, lower, floor_log2_selector<n, lower, upper>::middle>::value };
00476 };
00477 
00478 template<unsigned int n, int lower, int upper>
00479 struct floor_log2<n, lower, upper, floor_log2_move_up>
00480 {
00481   enum { value = floor_log2<n, floor_log2_selector<n, lower, upper>::middle, upper>::value };
00482 };
00483 
00484 template<unsigned int n, int lower, int upper>
00485 struct floor_log2<n, lower, upper, floor_log2_terminate>
00486 {
00487   enum { value = (n >= ((unsigned int)(1) << (lower+1))) ? lower+1 : lower };
00488 };
00489 
00490 template<unsigned int n, int lower, int upper>
00491 struct floor_log2<n, lower, upper, floor_log2_bogus>
00492 {
00493   // no value, error at compile time
00494 };
00495 
00496 template<typename Scalar>
00497 struct random_default_impl<Scalar, false, true>
00498 {
00499   typedef typename NumTraits<Scalar>::NonInteger NonInteger;
00500 
00501   static inline Scalar run(const Scalar& x, const Scalar& y)
00502   {
00503     return x + Scalar((NonInteger(y)-x+1) * std::rand() / (RAND_MAX + NonInteger(1)));
00504   }
00505 
00506   static inline Scalar run()
00507   {
00508 #ifdef EIGEN_MAKING_DOCS
00509     return run(Scalar(NumTraits<Scalar>::IsSigned ? -10 : 0), Scalar(10));
00510 #else
00511     enum { rand_bits = floor_log2<(unsigned int)(RAND_MAX)+1>::value,
00512            scalar_bits = sizeof(Scalar) * CHAR_BIT,
00513            shift = EIGEN_PLAIN_ENUM_MAX(0, int(rand_bits) - int(scalar_bits)),
00514            offset = NumTraits<Scalar>::IsSigned ? (1 << (EIGEN_PLAIN_ENUM_MIN(rand_bits,scalar_bits)-1)) : 0
00515     };
00516     return Scalar((std::rand() >> shift) - offset);
00517 #endif
00518   }
00519 };
00520 
00521 template<typename Scalar>
00522 struct random_default_impl<Scalar, true, false>
00523 {
00524   static inline Scalar run(const Scalar& x, const Scalar& y)
00525   {
00526     return Scalar(random(real(x), real(y)),
00527                   random(imag(x), imag(y)));
00528   }
00529   static inline Scalar run()
00530   {
00531     typedef typename NumTraits<Scalar>::Real RealScalar;
00532     return Scalar(random<RealScalar>(), random<RealScalar>());
00533   }
00534 };
00535 
00536 template<typename Scalar>
00537 inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random(const Scalar& x, const Scalar& y)
00538 {
00539   return EIGEN_MATHFUNC_IMPL(random, Scalar)::run(x, y);
00540 }
00541 
00542 template<typename Scalar>
00543 inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random()
00544 {
00545   return EIGEN_MATHFUNC_IMPL(random, Scalar)::run();
00546 }
00547 
00548 } // end namespace internal
00549 
00550 /****************************************************************************
00551 * Generic math function                                                    *
00552 ****************************************************************************/
00553 
00554 namespace numext {
00555 
00556 template<typename Scalar>
00557 inline EIGEN_MATHFUNC_RETVAL(real, Scalar) real(const Scalar& x)
00558 {
00559   return EIGEN_MATHFUNC_IMPL(real, Scalar)::run(x);
00560 }  
00561 
00562 template<typename Scalar>
00563 inline typename internal::add_const_on_value_type< EIGEN_MATHFUNC_RETVAL(real_ref, Scalar) >::type real_ref(const Scalar& x)
00564 {
00565   return internal::real_ref_impl<Scalar>::run(x);
00566 }
00567 
00568 template<typename Scalar>
00569 inline EIGEN_MATHFUNC_RETVAL(real_ref, Scalar) real_ref(Scalar& x)
00570 {
00571   return EIGEN_MATHFUNC_IMPL(real_ref, Scalar)::run(x);
00572 }
00573 
00574 template<typename Scalar>
00575 inline EIGEN_MATHFUNC_RETVAL(imag, Scalar) imag(const Scalar& x)
00576 {
00577   return EIGEN_MATHFUNC_IMPL(imag, Scalar)::run(x);
00578 }
00579 
00580 template<typename Scalar>
00581 inline typename internal::add_const_on_value_type< EIGEN_MATHFUNC_RETVAL(imag_ref, Scalar) >::type imag_ref(const Scalar& x)
00582 {
00583   return internal::imag_ref_impl<Scalar>::run(x);
00584 }
00585 
00586 template<typename Scalar>
00587 inline EIGEN_MATHFUNC_RETVAL(imag_ref, Scalar) imag_ref(Scalar& x)
00588 {
00589   return EIGEN_MATHFUNC_IMPL(imag_ref, Scalar)::run(x);
00590 }
00591 
00592 template<typename Scalar>
00593 inline EIGEN_MATHFUNC_RETVAL(conj, Scalar) conj(const Scalar& x)
00594 {
00595   return EIGEN_MATHFUNC_IMPL(conj, Scalar)::run(x);
00596 }
00597 
00598 template<typename Scalar>
00599 inline EIGEN_MATHFUNC_RETVAL(abs2, Scalar) abs2(const Scalar& x)
00600 {
00601   return EIGEN_MATHFUNC_IMPL(abs2, Scalar)::run(x);
00602 }
00603 
00604 template<typename Scalar>
00605 inline EIGEN_MATHFUNC_RETVAL(norm1, Scalar) norm1(const Scalar& x)
00606 {
00607   return EIGEN_MATHFUNC_IMPL(norm1, Scalar)::run(x);
00608 }
00609 
00610 template<typename Scalar>
00611 inline EIGEN_MATHFUNC_RETVAL(hypot, Scalar) hypot(const Scalar& x, const Scalar& y)
00612 {
00613   return EIGEN_MATHFUNC_IMPL(hypot, Scalar)::run(x, y);
00614 }
00615 
00616 template<typename Scalar>
00617 inline EIGEN_MATHFUNC_RETVAL(atanh2, Scalar) atanh2(const Scalar& x, const Scalar& y)
00618 {
00619   return EIGEN_MATHFUNC_IMPL(atanh2, Scalar)::run(x, y);
00620 }
00621 
00622 template<typename Scalar>
00623 inline EIGEN_MATHFUNC_RETVAL(pow, Scalar) pow(const Scalar& x, const Scalar& y)
00624 {
00625   return EIGEN_MATHFUNC_IMPL(pow, Scalar)::run(x, y);
00626 }
00627 
00628 // std::isfinite is non standard, so let's define our own version,
00629 // even though it is not very efficient.
00630 template<typename T> bool (isfinite)(const T& x)
00631 {
00632   return x<NumTraits<T>::highest() && x>NumTraits<T>::lowest();
00633 }
00634 
00635 } // end namespace numext
00636 
00637 namespace internal {
00638 
00639 /****************************************************************************
00640 * Implementation of fuzzy comparisons                                       *
00641 ****************************************************************************/
00642 
00643 template<typename Scalar,
00644          bool IsComplex,
00645          bool IsInteger>
00646 struct scalar_fuzzy_default_impl {};
00647 
00648 template<typename Scalar>
00649 struct scalar_fuzzy_default_impl<Scalar, false, false>
00650 {
00651   typedef typename NumTraits<Scalar>::Real RealScalar;
00652   template<typename OtherScalar>
00653   static inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y, const RealScalar& prec)
00654   {
00655     using std::abs;
00656     return abs(x) <= abs(y) * prec;
00657   }
00658   static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar& prec)
00659   {
00660     using std::min;
00661     using std::abs;
00662     return abs(x - y) <= (min)(abs(x), abs(y)) * prec;
00663   }
00664   static inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y, const RealScalar& prec)
00665   {
00666     return x <= y || isApprox(x, y, prec);
00667   }
00668 };
00669 
00670 template<typename Scalar>
00671 struct scalar_fuzzy_default_impl<Scalar, false, true>
00672 {
00673   typedef typename NumTraits<Scalar>::Real RealScalar;
00674   template<typename OtherScalar>
00675   static inline bool isMuchSmallerThan(const Scalar& x, const Scalar&, const RealScalar&)
00676   {
00677     return x == Scalar(0);
00678   }
00679   static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar&)
00680   {
00681     return x == y;
00682   }
00683   static inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y, const RealScalar&)
00684   {
00685     return x <= y;
00686   }
00687 };
00688 
00689 template<typename Scalar>
00690 struct scalar_fuzzy_default_impl<Scalar, true, false>
00691 {
00692   typedef typename NumTraits<Scalar>::Real RealScalar;
00693   template<typename OtherScalar>
00694   static inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y, const RealScalar& prec)
00695   {
00696     return numext::abs2(x) <= numext::abs2(y) * prec * prec;
00697   }
00698   static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar& prec)
00699   {
00700     using std::min;
00701     return numext::abs2(x - y) <= (min)(numext::abs2(x), numext::abs2(y)) * prec * prec;
00702   }
00703 };
00704 
00705 template<typename Scalar>
00706 struct scalar_fuzzy_impl : scalar_fuzzy_default_impl<Scalar, NumTraits<Scalar>::IsComplex, NumTraits<Scalar>::IsInteger> {};
00707 
00708 template<typename Scalar, typename OtherScalar>
00709 inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y,
00710                               const typename NumTraits<Scalar>::Real &precision = NumTraits<Scalar>::dummy_precision())
00711 {
00712   return scalar_fuzzy_impl<Scalar>::template isMuchSmallerThan<OtherScalar>(x, y, precision);
00713 }
00714 
00715 template<typename Scalar>
00716 inline bool isApprox(const Scalar& x, const Scalar& y,
00717                      const typename NumTraits<Scalar>::Real &precision = NumTraits<Scalar>::dummy_precision())
00718 {
00719   return scalar_fuzzy_impl<Scalar>::isApprox(x, y, precision);
00720 }
00721 
00722 template<typename Scalar>
00723 inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y,
00724                                const typename NumTraits<Scalar>::Real &precision = NumTraits<Scalar>::dummy_precision())
00725 {
00726   return scalar_fuzzy_impl<Scalar>::isApproxOrLessThan(x, y, precision);
00727 }
00728 
00729 /******************************************
00730 ***  The special case of the  bool type ***
00731 ******************************************/
00732 
00733 template<> struct random_impl<bool>
00734 {
00735   static inline bool run()
00736   {
00737     return random<int>(0,1)==0 ? false : true;
00738   }
00739 };
00740 
00741 template<> struct scalar_fuzzy_impl<bool>
00742 {
00743   typedef bool RealScalar;
00744   
00745   template<typename OtherScalar>
00746   static inline bool isMuchSmallerThan(const bool& x, const bool&, const bool&)
00747   {
00748     return !x;
00749   }
00750   
00751   static inline bool isApprox(bool x, bool y, bool)
00752   {
00753     return x == y;
00754   }
00755 
00756   static inline bool isApproxOrLessThan(const bool& x, const bool& y, const bool&)
00757   {
00758     return (!x) || y;
00759   }
00760   
00761 };
00762 
00763   
00764 } // end namespace internal
00765 
00766 } // end namespace Eigen
00767 
00768 #endif // EIGEN_MATHFUNCTIONS_H