Eigen  3.3.3
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 // source: http://www.geom.uiuc.edu/~huberty/math5337/groupe/digits.html
00014 // TODO this should better be moved to NumTraits
00015 #define EIGEN_PI 3.141592653589793238462643383279502884197169399375105820974944592307816406L
00016 
00017 
00018 namespace Eigen {
00019 
00020 // On WINCE, std::abs is defined for int only, so let's defined our own overloads:
00021 // This issue has been confirmed with MSVC 2008 only, but the issue might exist for more recent versions too.
00022 #if EIGEN_OS_WINCE && EIGEN_COMP_MSVC && EIGEN_COMP_MSVC<=1500
00023 long        abs(long        x) { return (labs(x));  }
00024 double      abs(double      x) { return (fabs(x));  }
00025 float       abs(float       x) { return (fabsf(x)); }
00026 long double abs(long double x) { return (fabsl(x)); }
00027 #endif
00028 
00029 namespace internal {
00030 
00051 template<typename T, typename dummy = void>
00052 struct global_math_functions_filtering_base
00053 {
00054   typedef T type;
00055 };
00056 
00057 template<typename T> struct always_void { typedef void type; };
00058 
00059 template<typename T>
00060 struct global_math_functions_filtering_base
00061   <T,
00062    typename always_void<typename T::Eigen_BaseClassForSpecializationOfGlobalMathFuncImpl>::type
00063   >
00064 {
00065   typedef typename T::Eigen_BaseClassForSpecializationOfGlobalMathFuncImpl type;
00066 };
00067 
00068 #define EIGEN_MATHFUNC_IMPL(func, scalar) Eigen::internal::func##_impl<typename Eigen::internal::global_math_functions_filtering_base<scalar>::type>
00069 #define EIGEN_MATHFUNC_RETVAL(func, scalar) typename Eigen::internal::func##_retval<typename Eigen::internal::global_math_functions_filtering_base<scalar>::type>::type
00070 
00071 /****************************************************************************
00072 * Implementation of real                                                 *
00073 ****************************************************************************/
00074 
00075 template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
00076 struct real_default_impl
00077 {
00078   typedef typename NumTraits<Scalar>::Real RealScalar;
00079   EIGEN_DEVICE_FUNC
00080   static inline RealScalar run(const Scalar& x)
00081   {
00082     return x;
00083   }
00084 };
00085 
00086 template<typename Scalar>
00087 struct real_default_impl<Scalar,true>
00088 {
00089   typedef typename NumTraits<Scalar>::Real RealScalar;
00090   EIGEN_DEVICE_FUNC
00091   static inline RealScalar run(const Scalar& x)
00092   {
00093     using std::real;
00094     return real(x);
00095   }
00096 };
00097 
00098 template<typename Scalar> struct real_impl : real_default_impl<Scalar> {};
00099 
00100 #ifdef __CUDA_ARCH__
00101 template<typename T>
00102 struct real_impl<std::complex<T> >
00103 {
00104   typedef T RealScalar;
00105   EIGEN_DEVICE_FUNC
00106   static inline T run(const std::complex<T>& x)
00107   {
00108     return x.real();
00109   }
00110 };
00111 #endif
00112 
00113 template<typename Scalar>
00114 struct real_retval
00115 {
00116   typedef typename NumTraits<Scalar>::Real type;
00117 };
00118 
00119 /****************************************************************************
00120 * Implementation of imag                                                 *
00121 ****************************************************************************/
00122 
00123 template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
00124 struct imag_default_impl
00125 {
00126   typedef typename NumTraits<Scalar>::Real RealScalar;
00127   EIGEN_DEVICE_FUNC
00128   static inline RealScalar run(const Scalar&)
00129   {
00130     return RealScalar(0);
00131   }
00132 };
00133 
00134 template<typename Scalar>
00135 struct imag_default_impl<Scalar,true>
00136 {
00137   typedef typename NumTraits<Scalar>::Real RealScalar;
00138   EIGEN_DEVICE_FUNC
00139   static inline RealScalar run(const Scalar& x)
00140   {
00141     using std::imag;
00142     return imag(x);
00143   }
00144 };
00145 
00146 template<typename Scalar> struct imag_impl : imag_default_impl<Scalar> {};
00147 
00148 #ifdef __CUDA_ARCH__
00149 template<typename T>
00150 struct imag_impl<std::complex<T> >
00151 {
00152   typedef T RealScalar;
00153   EIGEN_DEVICE_FUNC
00154   static inline T run(const std::complex<T>& x)
00155   {
00156     return x.imag();
00157   }
00158 };
00159 #endif
00160 
00161 template<typename Scalar>
00162 struct imag_retval
00163 {
00164   typedef typename NumTraits<Scalar>::Real type;
00165 };
00166 
00167 /****************************************************************************
00168 * Implementation of real_ref                                             *
00169 ****************************************************************************/
00170 
00171 template<typename Scalar>
00172 struct real_ref_impl
00173 {
00174   typedef typename NumTraits<Scalar>::Real RealScalar;
00175   EIGEN_DEVICE_FUNC
00176   static inline RealScalar& run(Scalar& x)
00177   {
00178     return reinterpret_cast<RealScalar*>(&x)[0];
00179   }
00180   EIGEN_DEVICE_FUNC
00181   static inline const RealScalar& run(const Scalar& x)
00182   {
00183     return reinterpret_cast<const RealScalar*>(&x)[0];
00184   }
00185 };
00186 
00187 template<typename Scalar>
00188 struct real_ref_retval
00189 {
00190   typedef typename NumTraits<Scalar>::Real & type;
00191 };
00192 
00193 /****************************************************************************
00194 * Implementation of imag_ref                                             *
00195 ****************************************************************************/
00196 
00197 template<typename Scalar, bool IsComplex>
00198 struct imag_ref_default_impl
00199 {
00200   typedef typename NumTraits<Scalar>::Real RealScalar;
00201   EIGEN_DEVICE_FUNC
00202   static inline RealScalar& run(Scalar& x)
00203   {
00204     return reinterpret_cast<RealScalar*>(&x)[1];
00205   }
00206   EIGEN_DEVICE_FUNC
00207   static inline const RealScalar& run(const Scalar& x)
00208   {
00209     return reinterpret_cast<RealScalar*>(&x)[1];
00210   }
00211 };
00212 
00213 template<typename Scalar>
00214 struct imag_ref_default_impl<Scalar, false>
00215 {
00216   EIGEN_DEVICE_FUNC
00217   static inline Scalar run(Scalar&)
00218   {
00219     return Scalar(0);
00220   }
00221   EIGEN_DEVICE_FUNC
00222   static inline const Scalar run(const Scalar&)
00223   {
00224     return Scalar(0);
00225   }
00226 };
00227 
00228 template<typename Scalar>
00229 struct imag_ref_impl : imag_ref_default_impl<Scalar, NumTraits<Scalar>::IsComplex> {};
00230 
00231 template<typename Scalar>
00232 struct imag_ref_retval
00233 {
00234   typedef typename NumTraits<Scalar>::Real & type;
00235 };
00236 
00237 /****************************************************************************
00238 * Implementation of conj                                                 *
00239 ****************************************************************************/
00240 
00241 template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
00242 struct conj_impl
00243 {
00244   EIGEN_DEVICE_FUNC
00245   static inline Scalar run(const Scalar& x)
00246   {
00247     return x;
00248   }
00249 };
00250 
00251 template<typename Scalar>
00252 struct conj_impl<Scalar,true>
00253 {
00254   EIGEN_DEVICE_FUNC
00255   static inline Scalar run(const Scalar& x)
00256   {
00257     using std::conj;
00258     return conj(x);
00259   }
00260 };
00261 
00262 template<typename Scalar>
00263 struct conj_retval
00264 {
00265   typedef Scalar type;
00266 };
00267 
00268 /****************************************************************************
00269 * Implementation of abs2                                                 *
00270 ****************************************************************************/
00271 
00272 template<typename Scalar,bool IsComplex>
00273 struct abs2_impl_default
00274 {
00275   typedef typename NumTraits<Scalar>::Real RealScalar;
00276   EIGEN_DEVICE_FUNC
00277   static inline RealScalar run(const Scalar& x)
00278   {
00279     return x*x;
00280   }
00281 };
00282 
00283 template<typename Scalar>
00284 struct abs2_impl_default<Scalar, true> // IsComplex
00285 {
00286   typedef typename NumTraits<Scalar>::Real RealScalar;
00287   EIGEN_DEVICE_FUNC
00288   static inline RealScalar run(const Scalar& x)
00289   {
00290     return real(x)*real(x) + imag(x)*imag(x);
00291   }
00292 };
00293 
00294 template<typename Scalar>
00295 struct abs2_impl
00296 {
00297   typedef typename NumTraits<Scalar>::Real RealScalar;
00298   EIGEN_DEVICE_FUNC
00299   static inline RealScalar run(const Scalar& x)
00300   {
00301     return abs2_impl_default<Scalar,NumTraits<Scalar>::IsComplex>::run(x);
00302   }
00303 };
00304 
00305 template<typename Scalar>
00306 struct abs2_retval
00307 {
00308   typedef typename NumTraits<Scalar>::Real type;
00309 };
00310 
00311 /****************************************************************************
00312 * Implementation of norm1                                                *
00313 ****************************************************************************/
00314 
00315 template<typename Scalar, bool IsComplex>
00316 struct norm1_default_impl
00317 {
00318   typedef typename NumTraits<Scalar>::Real RealScalar;
00319   EIGEN_DEVICE_FUNC
00320   static inline RealScalar run(const Scalar& x)
00321   {
00322     EIGEN_USING_STD_MATH(abs);
00323     return abs(real(x)) + abs(imag(x));
00324   }
00325 };
00326 
00327 template<typename Scalar>
00328 struct norm1_default_impl<Scalar, false>
00329 {
00330   EIGEN_DEVICE_FUNC
00331   static inline Scalar run(const Scalar& x)
00332   {
00333     EIGEN_USING_STD_MATH(abs);
00334     return abs(x);
00335   }
00336 };
00337 
00338 template<typename Scalar>
00339 struct norm1_impl : norm1_default_impl<Scalar, NumTraits<Scalar>::IsComplex> {};
00340 
00341 template<typename Scalar>
00342 struct norm1_retval
00343 {
00344   typedef typename NumTraits<Scalar>::Real type;
00345 };
00346 
00347 /****************************************************************************
00348 * Implementation of hypot                                                *
00349 ****************************************************************************/
00350 
00351 template<typename Scalar>
00352 struct hypot_impl
00353 {
00354   typedef typename NumTraits<Scalar>::Real RealScalar;
00355   static inline RealScalar run(const Scalar& x, const Scalar& y)
00356   {
00357     EIGEN_USING_STD_MATH(abs);
00358     EIGEN_USING_STD_MATH(sqrt);
00359     RealScalar _x = abs(x);
00360     RealScalar _y = abs(y);
00361     Scalar p, qp;
00362     if(_x>_y)
00363     {
00364       p = _x;
00365       qp = _y / p;
00366     }
00367     else
00368     {
00369       p = _y;
00370       qp = _x / p;
00371     }
00372     if(p==RealScalar(0)) return RealScalar(0);
00373     return p * sqrt(RealScalar(1) + qp*qp);
00374   }
00375 };
00376 
00377 template<typename Scalar>
00378 struct hypot_retval
00379 {
00380   typedef typename NumTraits<Scalar>::Real type;
00381 };
00382 
00383 /****************************************************************************
00384 * Implementation of cast                                                 *
00385 ****************************************************************************/
00386 
00387 template<typename OldType, typename NewType>
00388 struct cast_impl
00389 {
00390   EIGEN_DEVICE_FUNC
00391   static inline NewType run(const OldType& x)
00392   {
00393     return static_cast<NewType>(x);
00394   }
00395 };
00396 
00397 // here, for once, we're plainly returning NewType: we don't want cast to do weird things.
00398 
00399 template<typename OldType, typename NewType>
00400 EIGEN_DEVICE_FUNC
00401 inline NewType cast(const OldType& x)
00402 {
00403   return cast_impl<OldType, NewType>::run(x);
00404 }
00405 
00406 /****************************************************************************
00407 * Implementation of round                                                   *
00408 ****************************************************************************/
00409 
00410 #if EIGEN_HAS_CXX11_MATH
00411   template<typename Scalar>
00412   struct round_impl {
00413     static inline Scalar run(const Scalar& x)
00414     {
00415       EIGEN_STATIC_ASSERT((!NumTraits<Scalar>::IsComplex), NUMERIC_TYPE_MUST_BE_REAL)
00416       using std::round;
00417       return round(x);
00418     }
00419   };
00420 #else
00421   template<typename Scalar>
00422   struct round_impl
00423   {
00424     static inline Scalar run(const Scalar& x)
00425     {
00426       EIGEN_STATIC_ASSERT((!NumTraits<Scalar>::IsComplex), NUMERIC_TYPE_MUST_BE_REAL)
00427       EIGEN_USING_STD_MATH(floor);
00428       EIGEN_USING_STD_MATH(ceil);
00429       return (x > Scalar(0)) ? floor(x + Scalar(0.5)) : ceil(x - Scalar(0.5));
00430     }
00431   };
00432 #endif
00433 
00434 template<typename Scalar>
00435 struct round_retval
00436 {
00437   typedef Scalar type;
00438 };
00439 
00440 /****************************************************************************
00441 * Implementation of arg                                                     *
00442 ****************************************************************************/
00443 
00444 #if EIGEN_HAS_CXX11_MATH
00445   template<typename Scalar>
00446   struct arg_impl {
00447     static inline Scalar run(const Scalar& x)
00448     {
00449       EIGEN_USING_STD_MATH(arg);
00450       return arg(x);
00451     }
00452   };
00453 #else
00454   template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
00455   struct arg_default_impl
00456   {
00457     typedef typename NumTraits<Scalar>::Real RealScalar;
00458     EIGEN_DEVICE_FUNC
00459     static inline RealScalar run(const Scalar& x)
00460     {
00461       return (x < Scalar(0)) ? Scalar(EIGEN_PI) : Scalar(0); }
00462   };
00463 
00464   template<typename Scalar>
00465   struct arg_default_impl<Scalar,true>
00466   {
00467     typedef typename NumTraits<Scalar>::Real RealScalar;
00468     EIGEN_DEVICE_FUNC
00469     static inline RealScalar run(const Scalar& x)
00470     {
00471       EIGEN_USING_STD_MATH(arg);
00472       return arg(x);
00473     }
00474   };
00475 
00476   template<typename Scalar> struct arg_impl : arg_default_impl<Scalar> {};
00477 #endif
00478 
00479 template<typename Scalar>
00480 struct arg_retval
00481 {
00482   typedef typename NumTraits<Scalar>::Real type;
00483 };
00484 
00485 /****************************************************************************
00486 * Implementation of log1p                                                   *
00487 ****************************************************************************/
00488 
00489 namespace std_fallback {
00490   // fallback log1p implementation in case there is no log1p(Scalar) function in namespace of Scalar,
00491   // or that there is no suitable std::log1p function available
00492   template<typename Scalar>
00493   EIGEN_DEVICE_FUNC inline Scalar log1p(const Scalar& x) {
00494     EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar)
00495     typedef typename NumTraits<Scalar>::Real RealScalar;
00496     EIGEN_USING_STD_MATH(log);
00497     Scalar x1p = RealScalar(1) + x;
00498     return ( x1p == Scalar(1) ) ? x : x * ( log(x1p) / (x1p - RealScalar(1)) );
00499   }
00500 }
00501 
00502 template<typename Scalar>
00503 struct log1p_impl {
00504   static inline Scalar run(const Scalar& x)
00505   {
00506     EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar)
00507     #if EIGEN_HAS_CXX11_MATH
00508     using std::log1p;
00509     #endif
00510     using std_fallback::log1p;
00511     return log1p(x);
00512   }
00513 };
00514 
00515 
00516 template<typename Scalar>
00517 struct log1p_retval
00518 {
00519   typedef Scalar type;
00520 };
00521 
00522 /****************************************************************************
00523 * Implementation of pow                                                  *
00524 ****************************************************************************/
00525 
00526 template<typename ScalarX,typename ScalarY, bool IsInteger = NumTraits<ScalarX>::IsInteger&&NumTraits<ScalarY>::IsInteger>
00527 struct pow_impl
00528 {
00529   //typedef Scalar retval;
00530   typedef typename ScalarBinaryOpTraits<ScalarX,ScalarY,internal::scalar_pow_op<ScalarX,ScalarY> >::ReturnType result_type;
00531   static EIGEN_DEVICE_FUNC inline result_type run(const ScalarX& x, const ScalarY& y)
00532   {
00533     EIGEN_USING_STD_MATH(pow);
00534     return pow(x, y);
00535   }
00536 };
00537 
00538 template<typename ScalarX,typename ScalarY>
00539 struct pow_impl<ScalarX,ScalarY, true>
00540 {
00541   typedef ScalarX result_type;
00542   static EIGEN_DEVICE_FUNC inline ScalarX run(ScalarX x, ScalarY y)
00543   {
00544     ScalarX res(1);
00545     eigen_assert(!NumTraits<ScalarY>::IsSigned || y >= 0);
00546     if(y & 1) res *= x;
00547     y >>= 1;
00548     while(y)
00549     {
00550       x *= x;
00551       if(y&1) res *= x;
00552       y >>= 1;
00553     }
00554     return res;
00555   }
00556 };
00557 
00558 /****************************************************************************
00559 * Implementation of random                                               *
00560 ****************************************************************************/
00561 
00562 template<typename Scalar,
00563          bool IsComplex,
00564          bool IsInteger>
00565 struct random_default_impl {};
00566 
00567 template<typename Scalar>
00568 struct random_impl : random_default_impl<Scalar, NumTraits<Scalar>::IsComplex, NumTraits<Scalar>::IsInteger> {};
00569 
00570 template<typename Scalar>
00571 struct random_retval
00572 {
00573   typedef Scalar type;
00574 };
00575 
00576 template<typename Scalar> inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random(const Scalar& x, const Scalar& y);
00577 template<typename Scalar> inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random();
00578 
00579 template<typename Scalar>
00580 struct random_default_impl<Scalar, false, false>
00581 {
00582   static inline Scalar run(const Scalar& x, const Scalar& y)
00583   {
00584     return x + (y-x) * Scalar(std::rand()) / Scalar(RAND_MAX);
00585   }
00586   static inline Scalar run()
00587   {
00588     return run(Scalar(NumTraits<Scalar>::IsSigned ? -1 : 0), Scalar(1));
00589   }
00590 };
00591 
00592 enum {
00593   meta_floor_log2_terminate,
00594   meta_floor_log2_move_up,
00595   meta_floor_log2_move_down,
00596   meta_floor_log2_bogus
00597 };
00598 
00599 template<unsigned int n, int lower, int upper> struct meta_floor_log2_selector
00600 {
00601   enum { middle = (lower + upper) / 2,
00602          value = (upper <= lower + 1) ? int(meta_floor_log2_terminate)
00603                : (n < (1 << middle)) ? int(meta_floor_log2_move_down)
00604                : (n==0) ? int(meta_floor_log2_bogus)
00605                : int(meta_floor_log2_move_up)
00606   };
00607 };
00608 
00609 template<unsigned int n,
00610          int lower = 0,
00611          int upper = sizeof(unsigned int) * CHAR_BIT - 1,
00612          int selector = meta_floor_log2_selector<n, lower, upper>::value>
00613 struct meta_floor_log2 {};
00614 
00615 template<unsigned int n, int lower, int upper>
00616 struct meta_floor_log2<n, lower, upper, meta_floor_log2_move_down>
00617 {
00618   enum { value = meta_floor_log2<n, lower, meta_floor_log2_selector<n, lower, upper>::middle>::value };
00619 };
00620 
00621 template<unsigned int n, int lower, int upper>
00622 struct meta_floor_log2<n, lower, upper, meta_floor_log2_move_up>
00623 {
00624   enum { value = meta_floor_log2<n, meta_floor_log2_selector<n, lower, upper>::middle, upper>::value };
00625 };
00626 
00627 template<unsigned int n, int lower, int upper>
00628 struct meta_floor_log2<n, lower, upper, meta_floor_log2_terminate>
00629 {
00630   enum { value = (n >= ((unsigned int)(1) << (lower+1))) ? lower+1 : lower };
00631 };
00632 
00633 template<unsigned int n, int lower, int upper>
00634 struct meta_floor_log2<n, lower, upper, meta_floor_log2_bogus>
00635 {
00636   // no value, error at compile time
00637 };
00638 
00639 template<typename Scalar>
00640 struct random_default_impl<Scalar, false, true>
00641 {
00642   static inline Scalar run(const Scalar& x, const Scalar& y)
00643   { 
00644     typedef typename conditional<NumTraits<Scalar>::IsSigned,std::ptrdiff_t,std::size_t>::type ScalarX;
00645     if(y<x)
00646       return x;
00647     // the following difference might overflow on a 32 bits system,
00648     // but since y>=x the result converted to an unsigned long is still correct.
00649     std::size_t range = ScalarX(y)-ScalarX(x);
00650     std::size_t offset = 0;
00651     // rejection sampling
00652     std::size_t divisor = 1;
00653     std::size_t multiplier = 1;
00654     if(range<RAND_MAX) divisor = (std::size_t(RAND_MAX)+1)/(range+1);
00655     else               multiplier = 1 + range/(std::size_t(RAND_MAX)+1);
00656     do {
00657       offset = (std::size_t(std::rand()) * multiplier) / divisor;
00658     } while (offset > range);
00659     return Scalar(ScalarX(x) + offset);
00660   }
00661 
00662   static inline Scalar run()
00663   {
00664 #ifdef EIGEN_MAKING_DOCS
00665     return run(Scalar(NumTraits<Scalar>::IsSigned ? -10 : 0), Scalar(10));
00666 #else
00667     enum { rand_bits = meta_floor_log2<(unsigned int)(RAND_MAX)+1>::value,
00668            scalar_bits = sizeof(Scalar) * CHAR_BIT,
00669            shift = EIGEN_PLAIN_ENUM_MAX(0, int(rand_bits) - int(scalar_bits)),
00670            offset = NumTraits<Scalar>::IsSigned ? (1 << (EIGEN_PLAIN_ENUM_MIN(rand_bits,scalar_bits)-1)) : 0
00671     };
00672     return Scalar((std::rand() >> shift) - offset);
00673 #endif
00674   }
00675 };
00676 
00677 template<typename Scalar>
00678 struct random_default_impl<Scalar, true, false>
00679 {
00680   static inline Scalar run(const Scalar& x, const Scalar& y)
00681   {
00682     return Scalar(random(real(x), real(y)),
00683                   random(imag(x), imag(y)));
00684   }
00685   static inline Scalar run()
00686   {
00687     typedef typename NumTraits<Scalar>::Real RealScalar;
00688     return Scalar(random<RealScalar>(), random<RealScalar>());
00689   }
00690 };
00691 
00692 template<typename Scalar>
00693 inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random(const Scalar& x, const Scalar& y)
00694 {
00695   return EIGEN_MATHFUNC_IMPL(random, Scalar)::run(x, y);
00696 }
00697 
00698 template<typename Scalar>
00699 inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random()
00700 {
00701   return EIGEN_MATHFUNC_IMPL(random, Scalar)::run();
00702 }
00703 
00704 // Implementatin of is* functions
00705 
00706 // std::is* do not work with fast-math and gcc, std::is* are available on MSVC 2013 and newer, as well as in clang.
00707 #if (EIGEN_HAS_CXX11_MATH && !(EIGEN_COMP_GNUC_STRICT && __FINITE_MATH_ONLY__)) || (EIGEN_COMP_MSVC>=1800) || (EIGEN_COMP_CLANG)
00708 #define EIGEN_USE_STD_FPCLASSIFY 1
00709 #else
00710 #define EIGEN_USE_STD_FPCLASSIFY 0
00711 #endif
00712 
00713 template<typename T>
00714 EIGEN_DEVICE_FUNC
00715 typename internal::enable_if<internal::is_integral<T>::value,bool>::type
00716 isnan_impl(const T&) { return false; }
00717 
00718 template<typename T>
00719 EIGEN_DEVICE_FUNC
00720 typename internal::enable_if<internal::is_integral<T>::value,bool>::type
00721 isinf_impl(const T&) { return false; }
00722 
00723 template<typename T>
00724 EIGEN_DEVICE_FUNC
00725 typename internal::enable_if<internal::is_integral<T>::value,bool>::type
00726 isfinite_impl(const T&) { return true; }
00727 
00728 template<typename T>
00729 EIGEN_DEVICE_FUNC
00730 typename internal::enable_if<(!internal::is_integral<T>::value)&&(!NumTraits<T>::IsComplex),bool>::type
00731 isfinite_impl(const T& x)
00732 {
00733   #ifdef __CUDA_ARCH__
00734     return (::isfinite)(x);
00735   #elif EIGEN_USE_STD_FPCLASSIFY
00736     using std::isfinite;
00737     return isfinite EIGEN_NOT_A_MACRO (x);
00738   #else
00739     return x<=NumTraits<T>::highest() && x>=NumTraits<T>::lowest();
00740   #endif
00741 }
00742 
00743 template<typename T>
00744 EIGEN_DEVICE_FUNC
00745 typename internal::enable_if<(!internal::is_integral<T>::value)&&(!NumTraits<T>::IsComplex),bool>::type
00746 isinf_impl(const T& x)
00747 {
00748   #ifdef __CUDA_ARCH__
00749     return (::isinf)(x);
00750   #elif EIGEN_USE_STD_FPCLASSIFY
00751     using std::isinf;
00752     return isinf EIGEN_NOT_A_MACRO (x);
00753   #else
00754     return x>NumTraits<T>::highest() || x<NumTraits<T>::lowest();
00755   #endif
00756 }
00757 
00758 template<typename T>
00759 EIGEN_DEVICE_FUNC
00760 typename internal::enable_if<(!internal::is_integral<T>::value)&&(!NumTraits<T>::IsComplex),bool>::type
00761 isnan_impl(const T& x)
00762 {
00763   #ifdef __CUDA_ARCH__
00764     return (::isnan)(x);
00765   #elif EIGEN_USE_STD_FPCLASSIFY
00766     using std::isnan;
00767     return isnan EIGEN_NOT_A_MACRO (x);
00768   #else
00769     return x != x;
00770   #endif
00771 }
00772 
00773 #if (!EIGEN_USE_STD_FPCLASSIFY)
00774 
00775 #if EIGEN_COMP_MSVC
00776 
00777 template<typename T> EIGEN_DEVICE_FUNC bool isinf_msvc_helper(T x)
00778 {
00779   return _fpclass(x)==_FPCLASS_NINF || _fpclass(x)==_FPCLASS_PINF;
00780 }
00781 
00782 //MSVC defines a _isnan builtin function, but for double only
00783 EIGEN_DEVICE_FUNC inline bool isnan_impl(const long double& x) { return _isnan(x)!=0; }
00784 EIGEN_DEVICE_FUNC inline bool isnan_impl(const double& x)      { return _isnan(x)!=0; }
00785 EIGEN_DEVICE_FUNC inline bool isnan_impl(const float& x)       { return _isnan(x)!=0; }
00786 
00787 EIGEN_DEVICE_FUNC inline bool isinf_impl(const long double& x) { return isinf_msvc_helper(x); }
00788 EIGEN_DEVICE_FUNC inline bool isinf_impl(const double& x)      { return isinf_msvc_helper(x); }
00789 EIGEN_DEVICE_FUNC inline bool isinf_impl(const float& x)       { return isinf_msvc_helper(x); }
00790 
00791 #elif (defined __FINITE_MATH_ONLY__ && __FINITE_MATH_ONLY__ && EIGEN_COMP_GNUC)
00792 
00793 #if EIGEN_GNUC_AT_LEAST(5,0)
00794   #define EIGEN_TMP_NOOPT_ATTRIB EIGEN_DEVICE_FUNC inline __attribute__((optimize("no-finite-math-only")))
00795 #else
00796   // NOTE the inline qualifier and noinline attribute are both needed: the former is to avoid linking issue (duplicate symbol),
00797   //      while the second prevent too aggressive optimizations in fast-math mode:
00798   #define EIGEN_TMP_NOOPT_ATTRIB EIGEN_DEVICE_FUNC inline __attribute__((noinline,optimize("no-finite-math-only")))
00799 #endif
00800 
00801 template<> EIGEN_TMP_NOOPT_ATTRIB bool isnan_impl(const long double& x) { return __builtin_isnan(x); }
00802 template<> EIGEN_TMP_NOOPT_ATTRIB bool isnan_impl(const double& x)      { return __builtin_isnan(x); }
00803 template<> EIGEN_TMP_NOOPT_ATTRIB bool isnan_impl(const float& x)       { return __builtin_isnan(x); }
00804 template<> EIGEN_TMP_NOOPT_ATTRIB bool isinf_impl(const double& x)      { return __builtin_isinf(x); }
00805 template<> EIGEN_TMP_NOOPT_ATTRIB bool isinf_impl(const float& x)       { return __builtin_isinf(x); }
00806 template<> EIGEN_TMP_NOOPT_ATTRIB bool isinf_impl(const long double& x) { return __builtin_isinf(x); }
00807 
00808 #undef EIGEN_TMP_NOOPT_ATTRIB
00809 
00810 #endif
00811 
00812 #endif
00813 
00814 // The following overload are defined at the end of this file
00815 template<typename T> EIGEN_DEVICE_FUNC bool isfinite_impl(const std::complex<T>& x);
00816 template<typename T> EIGEN_DEVICE_FUNC bool isnan_impl(const std::complex<T>& x);
00817 template<typename T> EIGEN_DEVICE_FUNC bool isinf_impl(const std::complex<T>& x);
00818 
00819 template<typename T> T generic_fast_tanh_float(const T& a_x);
00820 
00821 } // end namespace internal
00822 
00823 /****************************************************************************
00824 * Generic math functions                                                    *
00825 ****************************************************************************/
00826 
00827 namespace numext {
00828 
00829 #ifndef __CUDA_ARCH__
00830 template<typename T>
00831 EIGEN_DEVICE_FUNC
00832 EIGEN_ALWAYS_INLINE T mini(const T& x, const T& y)
00833 {
00834   EIGEN_USING_STD_MATH(min);
00835   return min EIGEN_NOT_A_MACRO (x,y);
00836 }
00837 
00838 template<typename T>
00839 EIGEN_DEVICE_FUNC
00840 EIGEN_ALWAYS_INLINE T maxi(const T& x, const T& y)
00841 {
00842   EIGEN_USING_STD_MATH(max);
00843   return max EIGEN_NOT_A_MACRO (x,y);
00844 }
00845 #else
00846 template<typename T>
00847 EIGEN_DEVICE_FUNC
00848 EIGEN_ALWAYS_INLINE T mini(const T& x, const T& y)
00849 {
00850   return y < x ? y : x;
00851 }
00852 template<>
00853 EIGEN_DEVICE_FUNC
00854 EIGEN_ALWAYS_INLINE float mini(const float& x, const float& y)
00855 {
00856   return fminf(x, y);
00857 }
00858 template<typename T>
00859 EIGEN_DEVICE_FUNC
00860 EIGEN_ALWAYS_INLINE T maxi(const T& x, const T& y)
00861 {
00862   return x < y ? y : x;
00863 }
00864 template<>
00865 EIGEN_DEVICE_FUNC
00866 EIGEN_ALWAYS_INLINE float maxi(const float& x, const float& y)
00867 {
00868   return fmaxf(x, y);
00869 }
00870 #endif
00871 
00872 
00873 template<typename Scalar>
00874 EIGEN_DEVICE_FUNC
00875 inline EIGEN_MATHFUNC_RETVAL(real, Scalar) real(const Scalar& x)
00876 {
00877   return EIGEN_MATHFUNC_IMPL(real, Scalar)::run(x);
00878 }
00879 
00880 template<typename Scalar>
00881 EIGEN_DEVICE_FUNC
00882 inline typename internal::add_const_on_value_type< EIGEN_MATHFUNC_RETVAL(real_ref, Scalar) >::type real_ref(const Scalar& x)
00883 {
00884   return internal::real_ref_impl<Scalar>::run(x);
00885 }
00886 
00887 template<typename Scalar>
00888 EIGEN_DEVICE_FUNC
00889 inline EIGEN_MATHFUNC_RETVAL(real_ref, Scalar) real_ref(Scalar& x)
00890 {
00891   return EIGEN_MATHFUNC_IMPL(real_ref, Scalar)::run(x);
00892 }
00893 
00894 template<typename Scalar>
00895 EIGEN_DEVICE_FUNC
00896 inline EIGEN_MATHFUNC_RETVAL(imag, Scalar) imag(const Scalar& x)
00897 {
00898   return EIGEN_MATHFUNC_IMPL(imag, Scalar)::run(x);
00899 }
00900 
00901 template<typename Scalar>
00902 EIGEN_DEVICE_FUNC
00903 inline EIGEN_MATHFUNC_RETVAL(arg, Scalar) arg(const Scalar& x)
00904 {
00905   return EIGEN_MATHFUNC_IMPL(arg, Scalar)::run(x);
00906 }
00907 
00908 template<typename Scalar>
00909 EIGEN_DEVICE_FUNC
00910 inline typename internal::add_const_on_value_type< EIGEN_MATHFUNC_RETVAL(imag_ref, Scalar) >::type imag_ref(const Scalar& x)
00911 {
00912   return internal::imag_ref_impl<Scalar>::run(x);
00913 }
00914 
00915 template<typename Scalar>
00916 EIGEN_DEVICE_FUNC
00917 inline EIGEN_MATHFUNC_RETVAL(imag_ref, Scalar) imag_ref(Scalar& x)
00918 {
00919   return EIGEN_MATHFUNC_IMPL(imag_ref, Scalar)::run(x);
00920 }
00921 
00922 template<typename Scalar>
00923 EIGEN_DEVICE_FUNC
00924 inline EIGEN_MATHFUNC_RETVAL(conj, Scalar) conj(const Scalar& x)
00925 {
00926   return EIGEN_MATHFUNC_IMPL(conj, Scalar)::run(x);
00927 }
00928 
00929 template<typename Scalar>
00930 EIGEN_DEVICE_FUNC
00931 inline EIGEN_MATHFUNC_RETVAL(abs2, Scalar) abs2(const Scalar& x)
00932 {
00933   return EIGEN_MATHFUNC_IMPL(abs2, Scalar)::run(x);
00934 }
00935 
00936 template<typename Scalar>
00937 EIGEN_DEVICE_FUNC
00938 inline EIGEN_MATHFUNC_RETVAL(norm1, Scalar) norm1(const Scalar& x)
00939 {
00940   return EIGEN_MATHFUNC_IMPL(norm1, Scalar)::run(x);
00941 }
00942 
00943 template<typename Scalar>
00944 EIGEN_DEVICE_FUNC
00945 inline EIGEN_MATHFUNC_RETVAL(hypot, Scalar) hypot(const Scalar& x, const Scalar& y)
00946 {
00947   return EIGEN_MATHFUNC_IMPL(hypot, Scalar)::run(x, y);
00948 }
00949 
00950 template<typename Scalar>
00951 EIGEN_DEVICE_FUNC
00952 inline EIGEN_MATHFUNC_RETVAL(log1p, Scalar) log1p(const Scalar& x)
00953 {
00954   return EIGEN_MATHFUNC_IMPL(log1p, Scalar)::run(x);
00955 }
00956 
00957 #ifdef __CUDACC__
00958 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
00959 float log1p(const float &x) { return ::log1pf(x); }
00960 
00961 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
00962 double log1p(const double &x) { return ::log1p(x); }
00963 #endif
00964 
00965 template<typename ScalarX,typename ScalarY>
00966 EIGEN_DEVICE_FUNC
00967 inline typename internal::pow_impl<ScalarX,ScalarY>::result_type pow(const ScalarX& x, const ScalarY& y)
00968 {
00969   return internal::pow_impl<ScalarX,ScalarY>::run(x, y);
00970 }
00971 
00972 template<typename T> EIGEN_DEVICE_FUNC bool (isnan)   (const T &x) { return internal::isnan_impl(x); }
00973 template<typename T> EIGEN_DEVICE_FUNC bool (isinf)   (const T &x) { return internal::isinf_impl(x); }
00974 template<typename T> EIGEN_DEVICE_FUNC bool (isfinite)(const T &x) { return internal::isfinite_impl(x); }
00975 
00976 template<typename Scalar>
00977 EIGEN_DEVICE_FUNC
00978 inline EIGEN_MATHFUNC_RETVAL(round, Scalar) round(const Scalar& x)
00979 {
00980   return EIGEN_MATHFUNC_IMPL(round, Scalar)::run(x);
00981 }
00982 
00983 template<typename T>
00984 EIGEN_DEVICE_FUNC
00985 T (floor)(const T& x)
00986 {
00987   EIGEN_USING_STD_MATH(floor);
00988   return floor(x);
00989 }
00990 
00991 #ifdef __CUDACC__
00992 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
00993 float floor(const float &x) { return ::floorf(x); }
00994 
00995 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
00996 double floor(const double &x) { return ::floor(x); }
00997 #endif
00998 
00999 template<typename T>
01000 EIGEN_DEVICE_FUNC
01001 T (ceil)(const T& x)
01002 {
01003   EIGEN_USING_STD_MATH(ceil);
01004   return ceil(x);
01005 }
01006 
01007 #ifdef __CUDACC__
01008 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01009 float ceil(const float &x) { return ::ceilf(x); }
01010 
01011 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01012 double ceil(const double &x) { return ::ceil(x); }
01013 #endif
01014 
01015 
01018 inline int log2(int x)
01019 {
01020   eigen_assert(x>=0);
01021   unsigned int v(x);
01022   static const int table[32] = { 0, 9, 1, 10, 13, 21, 2, 29, 11, 14, 16, 18, 22, 25, 3, 30, 8, 12, 20, 28, 15, 17, 24, 7, 19, 27, 23, 6, 26, 5, 4, 31 };
01023   v |= v >> 1;
01024   v |= v >> 2;
01025   v |= v >> 4;
01026   v |= v >> 8;
01027   v |= v >> 16;
01028   return table[(v * 0x07C4ACDDU) >> 27];
01029 }
01030 
01039 template<typename T>
01040 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01041 T sqrt(const T &x)
01042 {
01043   EIGEN_USING_STD_MATH(sqrt);
01044   return sqrt(x);
01045 }
01046 
01047 template<typename T>
01048 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01049 T log(const T &x) {
01050   EIGEN_USING_STD_MATH(log);
01051   return log(x);
01052 }
01053 
01054 #ifdef __CUDACC__
01055 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01056 float log(const float &x) { return ::logf(x); }
01057 
01058 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01059 double log(const double &x) { return ::log(x); }
01060 #endif
01061 
01062 template<typename T>
01063 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01064 typename NumTraits<T>::Real abs(const T &x) {
01065   EIGEN_USING_STD_MATH(abs);
01066   return abs(x);
01067 }
01068 
01069 #ifdef __CUDACC__
01070 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01071 float abs(const float &x) { return ::fabsf(x); }
01072 
01073 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01074 double abs(const double &x) { return ::fabs(x); }
01075 
01076 template <> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01077 float abs(const std::complex<float>& x) {
01078   return ::hypotf(x.real(), x.imag());
01079 }
01080 
01081 template <> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01082 double abs(const std::complex<double>& x) {
01083   return ::hypot(x.real(), x.imag());
01084 }
01085 #endif
01086 
01087 template<typename T>
01088 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01089 T exp(const T &x) {
01090   EIGEN_USING_STD_MATH(exp);
01091   return exp(x);
01092 }
01093 
01094 #ifdef __CUDACC__
01095 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01096 float exp(const float &x) { return ::expf(x); }
01097 
01098 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01099 double exp(const double &x) { return ::exp(x); }
01100 #endif
01101 
01102 template<typename T>
01103 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01104 T cos(const T &x) {
01105   EIGEN_USING_STD_MATH(cos);
01106   return cos(x);
01107 }
01108 
01109 #ifdef __CUDACC__
01110 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01111 float cos(const float &x) { return ::cosf(x); }
01112 
01113 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01114 double cos(const double &x) { return ::cos(x); }
01115 #endif
01116 
01117 template<typename T>
01118 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01119 T sin(const T &x) {
01120   EIGEN_USING_STD_MATH(sin);
01121   return sin(x);
01122 }
01123 
01124 #ifdef __CUDACC__
01125 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01126 float sin(const float &x) { return ::sinf(x); }
01127 
01128 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01129 double sin(const double &x) { return ::sin(x); }
01130 #endif
01131 
01132 template<typename T>
01133 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01134 T tan(const T &x) {
01135   EIGEN_USING_STD_MATH(tan);
01136   return tan(x);
01137 }
01138 
01139 #ifdef __CUDACC__
01140 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01141 float tan(const float &x) { return ::tanf(x); }
01142 
01143 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01144 double tan(const double &x) { return ::tan(x); }
01145 #endif
01146 
01147 template<typename T>
01148 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01149 T acos(const T &x) {
01150   EIGEN_USING_STD_MATH(acos);
01151   return acos(x);
01152 }
01153 
01154 #ifdef __CUDACC__
01155 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01156 float acos(const float &x) { return ::acosf(x); }
01157 
01158 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01159 double acos(const double &x) { return ::acos(x); }
01160 #endif
01161 
01162 template<typename T>
01163 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01164 T asin(const T &x) {
01165   EIGEN_USING_STD_MATH(asin);
01166   return asin(x);
01167 }
01168 
01169 #ifdef __CUDACC__
01170 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01171 float asin(const float &x) { return ::asinf(x); }
01172 
01173 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01174 double asin(const double &x) { return ::asin(x); }
01175 #endif
01176 
01177 template<typename T>
01178 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01179 T atan(const T &x) {
01180   EIGEN_USING_STD_MATH(atan);
01181   return atan(x);
01182 }
01183 
01184 #ifdef __CUDACC__
01185 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01186 float atan(const float &x) { return ::atanf(x); }
01187 
01188 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01189 double atan(const double &x) { return ::atan(x); }
01190 #endif
01191 
01192 
01193 template<typename T>
01194 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01195 T cosh(const T &x) {
01196   EIGEN_USING_STD_MATH(cosh);
01197   return cosh(x);
01198 }
01199 
01200 #ifdef __CUDACC__
01201 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01202 float cosh(const float &x) { return ::coshf(x); }
01203 
01204 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01205 double cosh(const double &x) { return ::cosh(x); }
01206 #endif
01207 
01208 template<typename T>
01209 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01210 T sinh(const T &x) {
01211   EIGEN_USING_STD_MATH(sinh);
01212   return sinh(x);
01213 }
01214 
01215 #ifdef __CUDACC__
01216 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01217 float sinh(const float &x) { return ::sinhf(x); }
01218 
01219 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01220 double sinh(const double &x) { return ::sinh(x); }
01221 #endif
01222 
01223 template<typename T>
01224 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01225 T tanh(const T &x) {
01226   EIGEN_USING_STD_MATH(tanh);
01227   return tanh(x);
01228 }
01229 
01230 #if (!defined(__CUDACC__)) && EIGEN_FAST_MATH
01231 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01232 float tanh(float x) { return internal::generic_fast_tanh_float(x); }
01233 #endif
01234 
01235 #ifdef __CUDACC__
01236 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01237 float tanh(const float &x) { return ::tanhf(x); }
01238 
01239 template<> EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01240 double tanh(const double &x) { return ::tanh(x); }
01241 #endif
01242 
01243 template <typename T>
01244 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01245 T fmod(const T& a, const T& b) {
01246   EIGEN_USING_STD_MATH(fmod);
01247   return fmod(a, b);
01248 }
01249 
01250 #ifdef __CUDACC__
01251 template <>
01252 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01253 float fmod(const float& a, const float& b) {
01254   return ::fmodf(a, b);
01255 }
01256 
01257 template <>
01258 EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
01259 double fmod(const double& a, const double& b) {
01260   return ::fmod(a, b);
01261 }
01262 #endif
01263 
01264 } // end namespace numext
01265 
01266 namespace internal {
01267 
01268 template<typename T>
01269 EIGEN_DEVICE_FUNC bool isfinite_impl(const std::complex<T>& x)
01270 {
01271   return (numext::isfinite)(numext::real(x)) && (numext::isfinite)(numext::imag(x));
01272 }
01273 
01274 template<typename T>
01275 EIGEN_DEVICE_FUNC bool isnan_impl(const std::complex<T>& x)
01276 {
01277   return (numext::isnan)(numext::real(x)) || (numext::isnan)(numext::imag(x));
01278 }
01279 
01280 template<typename T>
01281 EIGEN_DEVICE_FUNC bool isinf_impl(const std::complex<T>& x)
01282 {
01283   return ((numext::isinf)(numext::real(x)) || (numext::isinf)(numext::imag(x))) && (!(numext::isnan)(x));
01284 }
01285 
01286 /****************************************************************************
01287 * Implementation of fuzzy comparisons                                       *
01288 ****************************************************************************/
01289 
01290 template<typename Scalar,
01291          bool IsComplex,
01292          bool IsInteger>
01293 struct scalar_fuzzy_default_impl {};
01294 
01295 template<typename Scalar>
01296 struct scalar_fuzzy_default_impl<Scalar, false, false>
01297 {
01298   typedef typename NumTraits<Scalar>::Real RealScalar;
01299   template<typename OtherScalar> EIGEN_DEVICE_FUNC
01300   static inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y, const RealScalar& prec)
01301   {
01302     return numext::abs(x) <= numext::abs(y) * prec;
01303   }
01304   EIGEN_DEVICE_FUNC
01305   static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar& prec)
01306   {
01307     return numext::abs(x - y) <= numext::mini(numext::abs(x), numext::abs(y)) * prec;
01308   }
01309   EIGEN_DEVICE_FUNC
01310   static inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y, const RealScalar& prec)
01311   {
01312     return x <= y || isApprox(x, y, prec);
01313   }
01314 };
01315 
01316 template<typename Scalar>
01317 struct scalar_fuzzy_default_impl<Scalar, false, true>
01318 {
01319   typedef typename NumTraits<Scalar>::Real RealScalar;
01320   template<typename OtherScalar> EIGEN_DEVICE_FUNC
01321   static inline bool isMuchSmallerThan(const Scalar& x, const Scalar&, const RealScalar&)
01322   {
01323     return x == Scalar(0);
01324   }
01325   EIGEN_DEVICE_FUNC
01326   static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar&)
01327   {
01328     return x == y;
01329   }
01330   EIGEN_DEVICE_FUNC
01331   static inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y, const RealScalar&)
01332   {
01333     return x <= y;
01334   }
01335 };
01336 
01337 template<typename Scalar>
01338 struct scalar_fuzzy_default_impl<Scalar, true, false>
01339 {
01340   typedef typename NumTraits<Scalar>::Real RealScalar;
01341   template<typename OtherScalar> EIGEN_DEVICE_FUNC
01342   static inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y, const RealScalar& prec)
01343   {
01344     return numext::abs2(x) <= numext::abs2(y) * prec * prec;
01345   }
01346   EIGEN_DEVICE_FUNC
01347   static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar& prec)
01348   {
01349     return numext::abs2(x - y) <= numext::mini(numext::abs2(x), numext::abs2(y)) * prec * prec;
01350   }
01351 };
01352 
01353 template<typename Scalar>
01354 struct scalar_fuzzy_impl : scalar_fuzzy_default_impl<Scalar, NumTraits<Scalar>::IsComplex, NumTraits<Scalar>::IsInteger> {};
01355 
01356 template<typename Scalar, typename OtherScalar> EIGEN_DEVICE_FUNC
01357 inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y,
01358                               const typename NumTraits<Scalar>::Real &precision = NumTraits<Scalar>::dummy_precision())
01359 {
01360   return scalar_fuzzy_impl<Scalar>::template isMuchSmallerThan<OtherScalar>(x, y, precision);
01361 }
01362 
01363 template<typename Scalar> EIGEN_DEVICE_FUNC
01364 inline bool isApprox(const Scalar& x, const Scalar& y,
01365                      const typename NumTraits<Scalar>::Real &precision = NumTraits<Scalar>::dummy_precision())
01366 {
01367   return scalar_fuzzy_impl<Scalar>::isApprox(x, y, precision);
01368 }
01369 
01370 template<typename Scalar> EIGEN_DEVICE_FUNC
01371 inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y,
01372                                const typename NumTraits<Scalar>::Real &precision = NumTraits<Scalar>::dummy_precision())
01373 {
01374   return scalar_fuzzy_impl<Scalar>::isApproxOrLessThan(x, y, precision);
01375 }
01376 
01377 /******************************************
01378 ***  The special case of the  bool type ***
01379 ******************************************/
01380 
01381 template<> struct random_impl<bool>
01382 {
01383   static inline bool run()
01384   {
01385     return random<int>(0,1)==0 ? false : true;
01386   }
01387 };
01388 
01389 template<> struct scalar_fuzzy_impl<bool>
01390 {
01391   typedef bool RealScalar;
01392   
01393   template<typename OtherScalar> EIGEN_DEVICE_FUNC
01394   static inline bool isMuchSmallerThan(const bool& x, const bool&, const bool&)
01395   {
01396     return !x;
01397   }
01398   
01399   EIGEN_DEVICE_FUNC
01400   static inline bool isApprox(bool x, bool y, bool)
01401   {
01402     return x == y;
01403   }
01404 
01405   EIGEN_DEVICE_FUNC
01406   static inline bool isApproxOrLessThan(const bool& x, const bool& y, const bool&)
01407   {
01408     return (!x) || y;
01409   }
01410   
01411 };
01412 
01413   
01414 } // end namespace internal
01415 
01416 } // end namespace Eigen
01417 
01418 #endif // EIGEN_MATHFUNCTIONS_H
 All Classes Functions Variables Typedefs Enumerations Enumerator Friends