MathFunctions.h
Go to the documentation of this file.
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 // Eigen is free software; you can redistribute it and/or
00007 // modify it under the terms of the GNU Lesser General Public
00008 // License as published by the Free Software Foundation; either
00009 // version 3 of the License, or (at your option) any later version.
00010 //
00011 // Alternatively, you can redistribute it and/or
00012 // modify it under the terms of the GNU General Public License as
00013 // published by the Free Software Foundation; either version 2 of
00014 // the License, or (at your option) any later version.
00015 //
00016 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
00017 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
00018 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
00019 // GNU General Public License for more details.
00020 //
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License and a copy of the GNU General Public License along with
00023 // Eigen. If not, see <http://www.gnu.org/licenses/>.
00024 
00025 #ifndef EIGEN_MATHFUNCTIONS_H
00026 #define EIGEN_MATHFUNCTIONS_H
00027 
00028 namespace Eigen {
00029 
00030 namespace internal {
00031 
00052 template<typename T, typename dummy = void>
00053 struct global_math_functions_filtering_base
00054 {
00055   typedef T type;
00056 };
00057 
00058 template<typename T> struct always_void { typedef void type; };
00059 
00060 template<typename T>
00061 struct global_math_functions_filtering_base
00062   <T,
00063    typename always_void<typename T::Eigen_BaseClassForSpecializationOfGlobalMathFuncImpl>::type
00064   >
00065 {
00066   typedef typename T::Eigen_BaseClassForSpecializationOfGlobalMathFuncImpl type;
00067 };
00068 
00069 #define EIGEN_MATHFUNC_IMPL(func, scalar) func##_impl<typename global_math_functions_filtering_base<scalar>::type>
00070 #define EIGEN_MATHFUNC_RETVAL(func, scalar) typename func##_retval<typename global_math_functions_filtering_base<scalar>::type>::type
00071 
00072 
00073 /****************************************************************************
00074 * Implementation of real                                                 *
00075 ****************************************************************************/
00076 
00077 template<typename Scalar>
00078 struct real_impl
00079 {
00080   typedef typename NumTraits<Scalar>::Real RealScalar;
00081   static inline RealScalar run(const Scalar& x)
00082   {
00083     return x;
00084   }
00085 };
00086 
00087 template<typename RealScalar>
00088 struct real_impl<std::complex<RealScalar> >
00089 {
00090   static inline RealScalar run(const std::complex<RealScalar>& x)
00091   {
00092     using std::real;
00093     return real(x);
00094   }
00095 };
00096 
00097 template<typename Scalar>
00098 struct real_retval
00099 {
00100   typedef typename NumTraits<Scalar>::Real type;
00101 };
00102 
00103 template<typename Scalar>
00104 inline EIGEN_MATHFUNC_RETVAL(real, Scalar) real(const Scalar& x)
00105 {
00106   return EIGEN_MATHFUNC_IMPL(real, Scalar)::run(x);
00107 }
00108 
00109 /****************************************************************************
00110 * Implementation of imag                                                 *
00111 ****************************************************************************/
00112 
00113 template<typename Scalar>
00114 struct imag_impl
00115 {
00116   typedef typename NumTraits<Scalar>::Real RealScalar;
00117   static inline RealScalar run(const Scalar&)
00118   {
00119     return RealScalar(0);
00120   }
00121 };
00122 
00123 template<typename RealScalar>
00124 struct imag_impl<std::complex<RealScalar> >
00125 {
00126   static inline RealScalar run(const std::complex<RealScalar>& x)
00127   {
00128     using std::imag;
00129     return imag(x);
00130   }
00131 };
00132 
00133 template<typename Scalar>
00134 struct imag_retval
00135 {
00136   typedef typename NumTraits<Scalar>::Real type;
00137 };
00138 
00139 template<typename Scalar>
00140 inline EIGEN_MATHFUNC_RETVAL(imag, Scalar) imag(const Scalar& x)
00141 {
00142   return EIGEN_MATHFUNC_IMPL(imag, Scalar)::run(x);
00143 }
00144 
00145 /****************************************************************************
00146 * Implementation of real_ref                                             *
00147 ****************************************************************************/
00148 
00149 template<typename Scalar>
00150 struct real_ref_impl
00151 {
00152   typedef typename NumTraits<Scalar>::Real RealScalar;
00153   static inline RealScalar& run(Scalar& x)
00154   {
00155     return reinterpret_cast<RealScalar*>(&x)[0];
00156   }
00157   static inline const RealScalar& run(const Scalar& x)
00158   {
00159     return reinterpret_cast<const RealScalar*>(&x)[0];
00160   }
00161 };
00162 
00163 template<typename Scalar>
00164 struct real_ref_retval
00165 {
00166   typedef typename NumTraits<Scalar>::Real & type;
00167 };
00168 
00169 template<typename Scalar>
00170 inline typename add_const_on_value_type< EIGEN_MATHFUNC_RETVAL(real_ref, Scalar) >::type real_ref(const Scalar& x)
00171 {
00172   return real_ref_impl<Scalar>::run(x);
00173 }
00174 
00175 template<typename Scalar>
00176 inline EIGEN_MATHFUNC_RETVAL(real_ref, Scalar) real_ref(Scalar& x)
00177 {
00178   return EIGEN_MATHFUNC_IMPL(real_ref, Scalar)::run(x);
00179 }
00180 
00181 /****************************************************************************
00182 * Implementation of imag_ref                                             *
00183 ****************************************************************************/
00184 
00185 template<typename Scalar, bool IsComplex>
00186 struct imag_ref_default_impl
00187 {
00188   typedef typename NumTraits<Scalar>::Real RealScalar;
00189   static inline RealScalar& run(Scalar& x)
00190   {
00191     return reinterpret_cast<RealScalar*>(&x)[1];
00192   }
00193   static inline const RealScalar& run(const Scalar& x)
00194   {
00195     return reinterpret_cast<RealScalar*>(&x)[1];
00196   }
00197 };
00198 
00199 template<typename Scalar>
00200 struct imag_ref_default_impl<Scalar, false>
00201 {
00202   static inline Scalar run(Scalar&)
00203   {
00204     return Scalar(0);
00205   }
00206   static inline const Scalar run(const Scalar&)
00207   {
00208     return Scalar(0);
00209   }
00210 };
00211 
00212 template<typename Scalar>
00213 struct imag_ref_impl : imag_ref_default_impl<Scalar, NumTraits<Scalar>::IsComplex> {};
00214 
00215 template<typename Scalar>
00216 struct imag_ref_retval
00217 {
00218   typedef typename NumTraits<Scalar>::Real & type;
00219 };
00220 
00221 template<typename Scalar>
00222 inline typename add_const_on_value_type< EIGEN_MATHFUNC_RETVAL(imag_ref, Scalar) >::type imag_ref(const Scalar& x)
00223 {
00224   return imag_ref_impl<Scalar>::run(x);
00225 }
00226 
00227 template<typename Scalar>
00228 inline EIGEN_MATHFUNC_RETVAL(imag_ref, Scalar) imag_ref(Scalar& x)
00229 {
00230   return EIGEN_MATHFUNC_IMPL(imag_ref, Scalar)::run(x);
00231 }
00232 
00233 /****************************************************************************
00234 * Implementation of conj                                                 *
00235 ****************************************************************************/
00236 
00237 template<typename Scalar>
00238 struct conj_impl
00239 {
00240   static inline Scalar run(const Scalar& x)
00241   {
00242     return x;
00243   }
00244 };
00245 
00246 template<typename RealScalar>
00247 struct conj_impl<std::complex<RealScalar> >
00248 {
00249   static inline std::complex<RealScalar> run(const std::complex<RealScalar>& x)
00250   {
00251     using std::conj;
00252     return conj(x);
00253   }
00254 };
00255 
00256 template<typename Scalar>
00257 struct conj_retval
00258 {
00259   typedef Scalar type;
00260 };
00261 
00262 template<typename Scalar>
00263 inline EIGEN_MATHFUNC_RETVAL(conj, Scalar) conj(const Scalar& x)
00264 {
00265   return EIGEN_MATHFUNC_IMPL(conj, Scalar)::run(x);
00266 }
00267 
00268 /****************************************************************************
00269 * Implementation of abs                                                  *
00270 ****************************************************************************/
00271 
00272 template<typename Scalar>
00273 struct abs_impl
00274 {
00275   typedef typename NumTraits<Scalar>::Real RealScalar;
00276   static inline RealScalar run(const Scalar& x)
00277   {
00278     using std::abs;
00279     return abs(x);
00280   }
00281 };
00282 
00283 template<typename Scalar>
00284 struct abs_retval
00285 {
00286   typedef typename NumTraits<Scalar>::Real type;
00287 };
00288 
00289 template<typename Scalar>
00290 inline EIGEN_MATHFUNC_RETVAL(abs, Scalar) abs(const Scalar& x)
00291 {
00292   return EIGEN_MATHFUNC_IMPL(abs, Scalar)::run(x);
00293 }
00294 
00295 /****************************************************************************
00296 * Implementation of abs2                                                 *
00297 ****************************************************************************/
00298 
00299 template<typename Scalar>
00300 struct abs2_impl
00301 {
00302   typedef typename NumTraits<Scalar>::Real RealScalar;
00303   static inline RealScalar run(const Scalar& x)
00304   {
00305     return x*x;
00306   }
00307 };
00308 
00309 template<typename RealScalar>
00310 struct abs2_impl<std::complex<RealScalar> >
00311 {
00312   static inline RealScalar run(const std::complex<RealScalar>& x)
00313   {
00314     return real(x)*real(x) + imag(x)*imag(x);
00315   }
00316 };
00317 
00318 template<typename Scalar>
00319 struct abs2_retval
00320 {
00321   typedef typename NumTraits<Scalar>::Real type;
00322 };
00323 
00324 template<typename Scalar>
00325 inline EIGEN_MATHFUNC_RETVAL(abs2, Scalar) abs2(const Scalar& x)
00326 {
00327   return EIGEN_MATHFUNC_IMPL(abs2, Scalar)::run(x);
00328 }
00329 
00330 /****************************************************************************
00331 * Implementation of norm1                                                *
00332 ****************************************************************************/
00333 
00334 template<typename Scalar, bool IsComplex>
00335 struct norm1_default_impl
00336 {
00337   typedef typename NumTraits<Scalar>::Real RealScalar;
00338   static inline RealScalar run(const Scalar& x)
00339   {
00340     return abs(real(x)) + abs(imag(x));
00341   }
00342 };
00343 
00344 template<typename Scalar>
00345 struct norm1_default_impl<Scalar, false>
00346 {
00347   static inline Scalar run(const Scalar& x)
00348   {
00349     return abs(x);
00350   }
00351 };
00352 
00353 template<typename Scalar>
00354 struct norm1_impl : norm1_default_impl<Scalar, NumTraits<Scalar>::IsComplex> {};
00355 
00356 template<typename Scalar>
00357 struct norm1_retval
00358 {
00359   typedef typename NumTraits<Scalar>::Real type;
00360 };
00361 
00362 template<typename Scalar>
00363 inline EIGEN_MATHFUNC_RETVAL(norm1, Scalar) norm1(const Scalar& x)
00364 {
00365   return EIGEN_MATHFUNC_IMPL(norm1, Scalar)::run(x);
00366 }
00367 
00368 /****************************************************************************
00369 * Implementation of hypot                                                *
00370 ****************************************************************************/
00371 
00372 template<typename Scalar>
00373 struct hypot_impl
00374 {
00375   typedef typename NumTraits<Scalar>::Real RealScalar;
00376   static inline RealScalar run(const Scalar& x, const Scalar& y)
00377   {
00378     using std::max;
00379     using std::min;
00380     RealScalar _x = abs(x);
00381     RealScalar _y = abs(y);
00382     RealScalar p = (max)(_x, _y);
00383     RealScalar q = (min)(_x, _y);
00384     RealScalar qp = q/p;
00385     return p * sqrt(RealScalar(1) + qp*qp);
00386   }
00387 };
00388 
00389 template<typename Scalar>
00390 struct hypot_retval
00391 {
00392   typedef typename NumTraits<Scalar>::Real type;
00393 };
00394 
00395 template<typename Scalar>
00396 inline EIGEN_MATHFUNC_RETVAL(hypot, Scalar) hypot(const Scalar& x, const Scalar& y)
00397 {
00398   return EIGEN_MATHFUNC_IMPL(hypot, Scalar)::run(x, y);
00399 }
00400 
00401 /****************************************************************************
00402 * Implementation of cast                                                 *
00403 ****************************************************************************/
00404 
00405 template<typename OldType, typename NewType>
00406 struct cast_impl
00407 {
00408   static inline NewType run(const OldType& x)
00409   {
00410     return static_cast<NewType>(x);
00411   }
00412 };
00413 
00414 // here, for once, we're plainly returning NewType: we don't want cast to do weird things.
00415 
00416 template<typename OldType, typename NewType>
00417 inline NewType cast(const OldType& x)
00418 {
00419   return cast_impl<OldType, NewType>::run(x);
00420 }
00421 
00422 /****************************************************************************
00423 * Implementation of sqrt                                                 *
00424 ****************************************************************************/
00425 
00426 template<typename Scalar, bool IsInteger>
00427 struct sqrt_default_impl
00428 {
00429   static inline Scalar run(const Scalar& x)
00430   {
00431     using std::sqrt;
00432     return sqrt(x);
00433   }
00434 };
00435 
00436 template<typename Scalar>
00437 struct sqrt_default_impl<Scalar, true>
00438 {
00439   static inline Scalar run(const Scalar&)
00440   {
00441 #ifdef EIGEN2_SUPPORT
00442     eigen_assert(!NumTraits<Scalar>::IsInteger);
00443 #else
00444     EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar)
00445 #endif
00446     return Scalar(0);
00447   }
00448 };
00449 
00450 template<typename Scalar>
00451 struct sqrt_impl : sqrt_default_impl<Scalar, NumTraits<Scalar>::IsInteger> {};
00452 
00453 template<typename Scalar>
00454 struct sqrt_retval
00455 {
00456   typedef Scalar type;
00457 };
00458 
00459 template<typename Scalar>
00460 inline EIGEN_MATHFUNC_RETVAL(sqrt, Scalar) sqrt(const Scalar& x)
00461 {
00462   return EIGEN_MATHFUNC_IMPL(sqrt, Scalar)::run(x);
00463 }
00464 
00465 /****************************************************************************
00466 * Implementation of standard unary real functions (exp, log, sin, cos, ...  *
00467 ****************************************************************************/
00468 
00469 // This macro instanciate all the necessary template mechanism which is common to all unary real functions.
00470 #define EIGEN_MATHFUNC_STANDARD_REAL_UNARY(NAME) \
00471   template<typename Scalar, bool IsInteger> struct NAME##_default_impl {            \
00472     static inline Scalar run(const Scalar& x) { using std::NAME; return NAME(x); }  \
00473   };                                                                                \
00474   template<typename Scalar> struct NAME##_default_impl<Scalar, true> {              \
00475     static inline Scalar run(const Scalar&) {                                       \
00476       EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar)                                       \
00477       return Scalar(0);                                                             \
00478     }                                                                               \
00479   };                                                                                \
00480   template<typename Scalar> struct NAME##_impl                                      \
00481     : NAME##_default_impl<Scalar, NumTraits<Scalar>::IsInteger>                     \
00482   {};                                                                               \
00483   template<typename Scalar> struct NAME##_retval { typedef Scalar type; };          \
00484   template<typename Scalar>                                                         \
00485   inline EIGEN_MATHFUNC_RETVAL(NAME, Scalar) NAME(const Scalar& x) {                \
00486     return EIGEN_MATHFUNC_IMPL(NAME, Scalar)::run(x);                               \
00487   }
00488 
00489 EIGEN_MATHFUNC_STANDARD_REAL_UNARY(exp)
00490 EIGEN_MATHFUNC_STANDARD_REAL_UNARY(log)
00491 EIGEN_MATHFUNC_STANDARD_REAL_UNARY(sin)
00492 EIGEN_MATHFUNC_STANDARD_REAL_UNARY(cos)
00493 EIGEN_MATHFUNC_STANDARD_REAL_UNARY(tan)
00494 EIGEN_MATHFUNC_STANDARD_REAL_UNARY(asin)
00495 EIGEN_MATHFUNC_STANDARD_REAL_UNARY(acos)
00496 
00497 /****************************************************************************
00498 * Implementation of atan2                                                *
00499 ****************************************************************************/
00500 
00501 template<typename Scalar, bool IsInteger>
00502 struct atan2_default_impl
00503 {
00504   typedef Scalar retval;
00505   static inline Scalar run(const Scalar& x, const Scalar& y)
00506   {
00507     using std::atan2;
00508     return atan2(x, y);
00509   }
00510 };
00511 
00512 template<typename Scalar>
00513 struct atan2_default_impl<Scalar, true>
00514 {
00515   static inline Scalar run(const Scalar&, const Scalar&)
00516   {
00517     EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar)
00518     return Scalar(0);
00519   }
00520 };
00521 
00522 template<typename Scalar>
00523 struct atan2_impl : atan2_default_impl<Scalar, NumTraits<Scalar>::IsInteger> {};
00524 
00525 template<typename Scalar>
00526 struct atan2_retval
00527 {
00528   typedef Scalar type;
00529 };
00530 
00531 template<typename Scalar>
00532 inline EIGEN_MATHFUNC_RETVAL(atan2, Scalar) atan2(const Scalar& x, const Scalar& y)
00533 {
00534   return EIGEN_MATHFUNC_IMPL(atan2, Scalar)::run(x, y);
00535 }
00536 
00537 /****************************************************************************
00538 * Implementation of pow                                                  *
00539 ****************************************************************************/
00540 
00541 template<typename Scalar, bool IsInteger>
00542 struct pow_default_impl
00543 {
00544   typedef Scalar retval;
00545   static inline Scalar run(const Scalar& x, const Scalar& y)
00546   {
00547     using std::pow;
00548     return pow(x, y);
00549   }
00550 };
00551 
00552 template<typename Scalar>
00553 struct pow_default_impl<Scalar, true>
00554 {
00555   static inline Scalar run(Scalar x, Scalar y)
00556   {
00557     Scalar res(1);
00558     eigen_assert(!NumTraits<Scalar>::IsSigned || y >= 0);
00559     if(y & 1) res *= x;
00560     y >>= 1;
00561     while(y)
00562     {
00563       x *= x;
00564       if(y&1) res *= x;
00565       y >>= 1;
00566     }
00567     return res;
00568   }
00569 };
00570 
00571 template<typename Scalar>
00572 struct pow_impl : pow_default_impl<Scalar, NumTraits<Scalar>::IsInteger> {};
00573 
00574 template<typename Scalar>
00575 struct pow_retval
00576 {
00577   typedef Scalar type;
00578 };
00579 
00580 template<typename Scalar>
00581 inline EIGEN_MATHFUNC_RETVAL(pow, Scalar) pow(const Scalar& x, const Scalar& y)
00582 {
00583   return EIGEN_MATHFUNC_IMPL(pow, Scalar)::run(x, y);
00584 }
00585 
00586 /****************************************************************************
00587 * Implementation of random                                               *
00588 ****************************************************************************/
00589 
00590 template<typename Scalar,
00591          bool IsComplex,
00592          bool IsInteger>
00593 struct random_default_impl {};
00594 
00595 template<typename Scalar>
00596 struct random_impl : random_default_impl<Scalar, NumTraits<Scalar>::IsComplex, NumTraits<Scalar>::IsInteger> {};
00597 
00598 template<typename Scalar>
00599 struct random_retval
00600 {
00601   typedef Scalar type;
00602 };
00603 
00604 template<typename Scalar> inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random(const Scalar& x, const Scalar& y);
00605 template<typename Scalar> inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random();
00606 
00607 template<typename Scalar>
00608 struct random_default_impl<Scalar, false, false>
00609 {
00610   static inline Scalar run(const Scalar& x, const Scalar& y)
00611   {
00612     return x + (y-x) * Scalar(std::rand()) / Scalar(RAND_MAX);
00613   }
00614   static inline Scalar run()
00615   {
00616     return run(Scalar(NumTraits<Scalar>::IsSigned ? -1 : 0), Scalar(1));
00617   }
00618 };
00619 
00620 enum {
00621   floor_log2_terminate,
00622   floor_log2_move_up,
00623   floor_log2_move_down,
00624   floor_log2_bogus
00625 };
00626 
00627 template<unsigned int n, int lower, int upper> struct floor_log2_selector
00628 {
00629   enum { middle = (lower + upper) / 2,
00630          value = (upper <= lower + 1) ? int(floor_log2_terminate)
00631                : (n < (1 << middle)) ? int(floor_log2_move_down)
00632                : (n==0) ? int(floor_log2_bogus)
00633                : int(floor_log2_move_up)
00634   };
00635 };
00636 
00637 template<unsigned int n,
00638          int lower = 0,
00639          int upper = sizeof(unsigned int) * CHAR_BIT - 1,
00640          int selector = floor_log2_selector<n, lower, upper>::value>
00641 struct floor_log2 {};
00642 
00643 template<unsigned int n, int lower, int upper>
00644 struct floor_log2<n, lower, upper, floor_log2_move_down>
00645 {
00646   enum { value = floor_log2<n, lower, floor_log2_selector<n, lower, upper>::middle>::value };
00647 };
00648 
00649 template<unsigned int n, int lower, int upper>
00650 struct floor_log2<n, lower, upper, floor_log2_move_up>
00651 {
00652   enum { value = floor_log2<n, floor_log2_selector<n, lower, upper>::middle, upper>::value };
00653 };
00654 
00655 template<unsigned int n, int lower, int upper>
00656 struct floor_log2<n, lower, upper, floor_log2_terminate>
00657 {
00658   enum { value = (n >= ((unsigned int)(1) << (lower+1))) ? lower+1 : lower };
00659 };
00660 
00661 template<unsigned int n, int lower, int upper>
00662 struct floor_log2<n, lower, upper, floor_log2_bogus>
00663 {
00664   // no value, error at compile time
00665 };
00666 
00667 template<typename Scalar>
00668 struct random_default_impl<Scalar, false, true>
00669 {
00670   typedef typename NumTraits<Scalar>::NonInteger NonInteger;
00671 
00672   static inline Scalar run(const Scalar& x, const Scalar& y)
00673   {
00674     return x + Scalar((NonInteger(y)-x+1) * std::rand() / (RAND_MAX + NonInteger(1)));
00675   }
00676 
00677   static inline Scalar run()
00678   {
00679 #ifdef EIGEN_MAKING_DOCS
00680     return run(Scalar(NumTraits<Scalar>::IsSigned ? -10 : 0), Scalar(10));
00681 #else
00682     enum { rand_bits = floor_log2<(unsigned int)(RAND_MAX)+1>::value,
00683            scalar_bits = sizeof(Scalar) * CHAR_BIT,
00684            shift = EIGEN_PLAIN_ENUM_MAX(0, int(rand_bits) - int(scalar_bits))
00685     };
00686     Scalar x = Scalar(std::rand() >> shift);
00687     Scalar offset = NumTraits<Scalar>::IsSigned ? Scalar(1 << (rand_bits-1)) : Scalar(0);
00688     return x - offset;
00689 #endif
00690   }
00691 };
00692 
00693 template<typename Scalar>
00694 struct random_default_impl<Scalar, true, false>
00695 {
00696   static inline Scalar run(const Scalar& x, const Scalar& y)
00697   {
00698     return Scalar(random(real(x), real(y)),
00699                   random(imag(x), imag(y)));
00700   }
00701   static inline Scalar run()
00702   {
00703     typedef typename NumTraits<Scalar>::Real RealScalar;
00704     return Scalar(random<RealScalar>(), random<RealScalar>());
00705   }
00706 };
00707 
00708 template<typename Scalar>
00709 inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random(const Scalar& x, const Scalar& y)
00710 {
00711   return EIGEN_MATHFUNC_IMPL(random, Scalar)::run(x, y);
00712 }
00713 
00714 template<typename Scalar>
00715 inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random()
00716 {
00717   return EIGEN_MATHFUNC_IMPL(random, Scalar)::run();
00718 }
00719 
00720 /****************************************************************************
00721 * Implementation of fuzzy comparisons                                       *
00722 ****************************************************************************/
00723 
00724 template<typename Scalar,
00725          bool IsComplex,
00726          bool IsInteger>
00727 struct scalar_fuzzy_default_impl {};
00728 
00729 template<typename Scalar>
00730 struct scalar_fuzzy_default_impl<Scalar, false, false>
00731 {
00732   typedef typename NumTraits<Scalar>::Real RealScalar;
00733   template<typename OtherScalar>
00734   static inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y, const RealScalar& prec)
00735   {
00736     return abs(x) <= abs(y) * prec;
00737   }
00738   static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar& prec)
00739   {
00740     using std::min;
00741     return abs(x - y) <= (min)(abs(x), abs(y)) * prec;
00742   }
00743   static inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y, const RealScalar& prec)
00744   {
00745     return x <= y || isApprox(x, y, prec);
00746   }
00747 };
00748 
00749 template<typename Scalar>
00750 struct scalar_fuzzy_default_impl<Scalar, false, true>
00751 {
00752   typedef typename NumTraits<Scalar>::Real RealScalar;
00753   template<typename OtherScalar>
00754   static inline bool isMuchSmallerThan(const Scalar& x, const Scalar&, const RealScalar&)
00755   {
00756     return x == Scalar(0);
00757   }
00758   static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar&)
00759   {
00760     return x == y;
00761   }
00762   static inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y, const RealScalar&)
00763   {
00764     return x <= y;
00765   }
00766 };
00767 
00768 template<typename Scalar>
00769 struct scalar_fuzzy_default_impl<Scalar, true, false>
00770 {
00771   typedef typename NumTraits<Scalar>::Real RealScalar;
00772   template<typename OtherScalar>
00773   static inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y, const RealScalar& prec)
00774   {
00775     return abs2(x) <= abs2(y) * prec * prec;
00776   }
00777   static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar& prec)
00778   {
00779     using std::min;
00780     return abs2(x - y) <= (min)(abs2(x), abs2(y)) * prec * prec;
00781   }
00782 };
00783 
00784 template<typename Scalar>
00785 struct scalar_fuzzy_impl : scalar_fuzzy_default_impl<Scalar, NumTraits<Scalar>::IsComplex, NumTraits<Scalar>::IsInteger> {};
00786 
00787 template<typename Scalar, typename OtherScalar>
00788 inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y,
00789                                    typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision())
00790 {
00791   return scalar_fuzzy_impl<Scalar>::template isMuchSmallerThan<OtherScalar>(x, y, precision);
00792 }
00793 
00794 template<typename Scalar>
00795 inline bool isApprox(const Scalar& x, const Scalar& y,
00796                           typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision())
00797 {
00798   return scalar_fuzzy_impl<Scalar>::isApprox(x, y, precision);
00799 }
00800 
00801 template<typename Scalar>
00802 inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y,
00803                                     typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision())
00804 {
00805   return scalar_fuzzy_impl<Scalar>::isApproxOrLessThan(x, y, precision);
00806 }
00807 
00808 /******************************************
00809 ***  The special case of the  bool type ***
00810 ******************************************/
00811 
00812 template<> struct random_impl<bool>
00813 {
00814   static inline bool run()
00815   {
00816     return random<int>(0,1)==0 ? false : true;
00817   }
00818 };
00819 
00820 template<> struct scalar_fuzzy_impl<bool>
00821 {
00822   typedef bool RealScalar;
00823   
00824   template<typename OtherScalar>
00825   static inline bool isMuchSmallerThan(const bool& x, const bool&, const bool&)
00826   {
00827     return !x;
00828   }
00829   
00830   static inline bool isApprox(bool x, bool y, bool)
00831   {
00832     return x == y;
00833   }
00834 
00835   static inline bool isApproxOrLessThan(const bool& x, const bool& y, const bool&)
00836   {
00837     return (!x) || y;
00838   }
00839   
00840 };
00841 
00842 /****************************************************************************
00843 * Special functions                                                          *
00844 ****************************************************************************/
00845 
00846 // std::isfinite is non standard, so let's define our own version,
00847 // even though it is not very efficient.
00848 template<typename T> bool isfinite(const T& x)
00849 {
00850   return x<NumTraits<T>::highest() && x>NumTraits<T>::lowest();
00851 }
00852 
00853 } // end namespace internal
00854 
00855 } // end namespace Eigen
00856 
00857 #endif // EIGEN_MATHFUNCTIONS_H