Eigen  3.2.92
MathFunctions.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2006-2010 Benoit Jacob <jacob.benoit.1@gmail.com>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #ifndef EIGEN_MATHFUNCTIONS_H
11 #define EIGEN_MATHFUNCTIONS_H
12 
13 // source: http://www.geom.uiuc.edu/~huberty/math5337/groupe/digits.html
14 #define EIGEN_PI 3.141592653589793238462643383279502884197169399375105820974944592307816406
15 
16 namespace Eigen {
17 
18 // On WINCE, std::abs is defined for int only, so let's defined our own overloads:
19 // This issue has been confirmed with MSVC 2008 only, but the issue might exist for more recent versions too.
20 #if EIGEN_OS_WINCE && EIGEN_COMP_MSVC && EIGEN_COMP_MSVC<=1500
21 long abs(long x) { return (labs(x)); }
22 double abs(double x) { return (fabs(x)); }
23 float abs(float x) { return (fabsf(x)); }
24 long double abs(long double x) { return (fabsl(x)); }
25 #endif
26 
27 namespace internal {
28 
49 template<typename T, typename dummy = void>
50 struct global_math_functions_filtering_base
51 {
52  typedef T type;
53 };
54 
55 template<typename T> struct always_void { typedef void type; };
56 
57 template<typename T>
58 struct global_math_functions_filtering_base
59  <T,
60  typename always_void<typename T::Eigen_BaseClassForSpecializationOfGlobalMathFuncImpl>::type
61  >
62 {
63  typedef typename T::Eigen_BaseClassForSpecializationOfGlobalMathFuncImpl type;
64 };
65 
66 #define EIGEN_MATHFUNC_IMPL(func, scalar) Eigen::internal::func##_impl<typename Eigen::internal::global_math_functions_filtering_base<scalar>::type>
67 #define EIGEN_MATHFUNC_RETVAL(func, scalar) typename Eigen::internal::func##_retval<typename Eigen::internal::global_math_functions_filtering_base<scalar>::type>::type
68 
69 /****************************************************************************
70 * Implementation of real *
71 ****************************************************************************/
72 
73 template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
74 struct real_default_impl
75 {
76  typedef typename NumTraits<Scalar>::Real RealScalar;
77  EIGEN_DEVICE_FUNC
78  static inline RealScalar run(const Scalar& x)
79  {
80  return x;
81  }
82 };
83 
84 template<typename Scalar>
85 struct real_default_impl<Scalar,true>
86 {
87  typedef typename NumTraits<Scalar>::Real RealScalar;
88  EIGEN_DEVICE_FUNC
89  static inline RealScalar run(const Scalar& x)
90  {
91  using std::real;
92  return real(x);
93  }
94 };
95 
96 template<typename Scalar> struct real_impl : real_default_impl<Scalar> {};
97 
98 template<typename Scalar>
99 struct real_retval
100 {
101  typedef typename NumTraits<Scalar>::Real type;
102 };
103 
104 /****************************************************************************
105 * Implementation of imag *
106 ****************************************************************************/
107 
108 template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
109 struct imag_default_impl
110 {
111  typedef typename NumTraits<Scalar>::Real RealScalar;
112  EIGEN_DEVICE_FUNC
113  static inline RealScalar run(const Scalar&)
114  {
115  return RealScalar(0);
116  }
117 };
118 
119 template<typename Scalar>
120 struct imag_default_impl<Scalar,true>
121 {
122  typedef typename NumTraits<Scalar>::Real RealScalar;
123  EIGEN_DEVICE_FUNC
124  static inline RealScalar run(const Scalar& x)
125  {
126  using std::imag;
127  return imag(x);
128  }
129 };
130 
131 template<typename Scalar> struct imag_impl : imag_default_impl<Scalar> {};
132 
133 template<typename Scalar>
134 struct imag_retval
135 {
136  typedef typename NumTraits<Scalar>::Real type;
137 };
138 
139 /****************************************************************************
140 * Implementation of real_ref *
141 ****************************************************************************/
142 
143 template<typename Scalar>
144 struct real_ref_impl
145 {
146  typedef typename NumTraits<Scalar>::Real RealScalar;
147  EIGEN_DEVICE_FUNC
148  static inline RealScalar& run(Scalar& x)
149  {
150  return reinterpret_cast<RealScalar*>(&x)[0];
151  }
152  EIGEN_DEVICE_FUNC
153  static inline const RealScalar& run(const Scalar& x)
154  {
155  return reinterpret_cast<const RealScalar*>(&x)[0];
156  }
157 };
158 
159 template<typename Scalar>
160 struct real_ref_retval
161 {
162  typedef typename NumTraits<Scalar>::Real & type;
163 };
164 
165 /****************************************************************************
166 * Implementation of imag_ref *
167 ****************************************************************************/
168 
169 template<typename Scalar, bool IsComplex>
170 struct imag_ref_default_impl
171 {
172  typedef typename NumTraits<Scalar>::Real RealScalar;
173  EIGEN_DEVICE_FUNC
174  static inline RealScalar& run(Scalar& x)
175  {
176  return reinterpret_cast<RealScalar*>(&x)[1];
177  }
178  EIGEN_DEVICE_FUNC
179  static inline const RealScalar& run(const Scalar& x)
180  {
181  return reinterpret_cast<RealScalar*>(&x)[1];
182  }
183 };
184 
185 template<typename Scalar>
186 struct imag_ref_default_impl<Scalar, false>
187 {
188  EIGEN_DEVICE_FUNC
189  static inline Scalar run(Scalar&)
190  {
191  return Scalar(0);
192  }
193  EIGEN_DEVICE_FUNC
194  static inline const Scalar run(const Scalar&)
195  {
196  return Scalar(0);
197  }
198 };
199 
200 template<typename Scalar>
201 struct imag_ref_impl : imag_ref_default_impl<Scalar, NumTraits<Scalar>::IsComplex> {};
202 
203 template<typename Scalar>
204 struct imag_ref_retval
205 {
206  typedef typename NumTraits<Scalar>::Real & type;
207 };
208 
209 /****************************************************************************
210 * Implementation of conj *
211 ****************************************************************************/
212 
213 template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
214 struct conj_impl
215 {
216  EIGEN_DEVICE_FUNC
217  static inline Scalar run(const Scalar& x)
218  {
219  return x;
220  }
221 };
222 
223 template<typename Scalar>
224 struct conj_impl<Scalar,true>
225 {
226  EIGEN_DEVICE_FUNC
227  static inline Scalar run(const Scalar& x)
228  {
229  using std::conj;
230  return conj(x);
231  }
232 };
233 
234 template<typename Scalar>
235 struct conj_retval
236 {
237  typedef Scalar type;
238 };
239 
240 /****************************************************************************
241 * Implementation of abs2 *
242 ****************************************************************************/
243 
244 template<typename Scalar,bool IsComplex>
245 struct abs2_impl_default
246 {
247  typedef typename NumTraits<Scalar>::Real RealScalar;
248  EIGEN_DEVICE_FUNC
249  static inline RealScalar run(const Scalar& x)
250  {
251  return x*x;
252  }
253 };
254 
255 template<typename Scalar>
256 struct abs2_impl_default<Scalar, true> // IsComplex
257 {
258  typedef typename NumTraits<Scalar>::Real RealScalar;
259  EIGEN_DEVICE_FUNC
260  static inline RealScalar run(const Scalar& x)
261  {
262  return real(x)*real(x) + imag(x)*imag(x);
263  }
264 };
265 
266 template<typename Scalar>
267 struct abs2_impl
268 {
269  typedef typename NumTraits<Scalar>::Real RealScalar;
270  EIGEN_DEVICE_FUNC
271  static inline RealScalar run(const Scalar& x)
272  {
273  return abs2_impl_default<Scalar,NumTraits<Scalar>::IsComplex>::run(x);
274  }
275 };
276 
277 template<typename Scalar>
278 struct abs2_retval
279 {
280  typedef typename NumTraits<Scalar>::Real type;
281 };
282 
283 /****************************************************************************
284 * Implementation of norm1 *
285 ****************************************************************************/
286 
287 template<typename Scalar, bool IsComplex>
288 struct norm1_default_impl
289 {
290  typedef typename NumTraits<Scalar>::Real RealScalar;
291  EIGEN_DEVICE_FUNC
292  static inline RealScalar run(const Scalar& x)
293  {
294  EIGEN_USING_STD_MATH(abs);
295  return abs(real(x)) + abs(imag(x));
296  }
297 };
298 
299 template<typename Scalar>
300 struct norm1_default_impl<Scalar, false>
301 {
302  EIGEN_DEVICE_FUNC
303  static inline Scalar run(const Scalar& x)
304  {
305  EIGEN_USING_STD_MATH(abs);
306  return abs(x);
307  }
308 };
309 
310 template<typename Scalar>
311 struct norm1_impl : norm1_default_impl<Scalar, NumTraits<Scalar>::IsComplex> {};
312 
313 template<typename Scalar>
314 struct norm1_retval
315 {
316  typedef typename NumTraits<Scalar>::Real type;
317 };
318 
319 /****************************************************************************
320 * Implementation of hypot *
321 ****************************************************************************/
322 
323 template<typename Scalar>
324 struct hypot_impl
325 {
326  typedef typename NumTraits<Scalar>::Real RealScalar;
327  static inline RealScalar run(const Scalar& x, const Scalar& y)
328  {
329  EIGEN_USING_STD_MATH(abs);
330  EIGEN_USING_STD_MATH(sqrt);
331  RealScalar _x = abs(x);
332  RealScalar _y = abs(y);
333  Scalar p, qp;
334  if(_x>_y)
335  {
336  p = _x;
337  qp = _y / p;
338  }
339  else
340  {
341  p = _y;
342  qp = _x / p;
343  }
344  if(p==RealScalar(0)) return RealScalar(0);
345  return p * sqrt(RealScalar(1) + qp*qp);
346  }
347 };
348 
349 template<typename Scalar>
350 struct hypot_retval
351 {
352  typedef typename NumTraits<Scalar>::Real type;
353 };
354 
355 /****************************************************************************
356 * Implementation of cast *
357 ****************************************************************************/
358 
359 template<typename OldType, typename NewType>
360 struct cast_impl
361 {
362  EIGEN_DEVICE_FUNC
363  static inline NewType run(const OldType& x)
364  {
365  return static_cast<NewType>(x);
366  }
367 };
368 
369 // here, for once, we're plainly returning NewType: we don't want cast to do weird things.
370 
371 template<typename OldType, typename NewType>
372 EIGEN_DEVICE_FUNC
373 inline NewType cast(const OldType& x)
374 {
375  return cast_impl<OldType, NewType>::run(x);
376 }
377 
378 /****************************************************************************
379 * Implementation of round *
380 ****************************************************************************/
381 
382 #if EIGEN_HAS_CXX11_MATH
383  template<typename Scalar>
384  struct round_impl {
385  static inline Scalar run(const Scalar& x)
386  {
387  EIGEN_STATIC_ASSERT((!NumTraits<Scalar>::IsComplex), NUMERIC_TYPE_MUST_BE_REAL)
388  using std::round;
389  return round(x);
390  }
391  };
392 #else
393  template<typename Scalar>
394  struct round_impl
395  {
396  static inline Scalar run(const Scalar& x)
397  {
398  EIGEN_STATIC_ASSERT((!NumTraits<Scalar>::IsComplex), NUMERIC_TYPE_MUST_BE_REAL)
399  EIGEN_USING_STD_MATH(floor);
400  EIGEN_USING_STD_MATH(ceil);
401  return (x > Scalar(0)) ? floor(x + Scalar(0.5)) : ceil(x - Scalar(0.5));
402  }
403  };
404 #endif
405 
406 template<typename Scalar>
407 struct round_retval
408 {
409  typedef Scalar type;
410 };
411 
412 /****************************************************************************
413 * Implementation of arg *
414 ****************************************************************************/
415 
416 #if EIGEN_HAS_CXX11_MATH
417  template<typename Scalar>
418  struct arg_impl {
419  static inline Scalar run(const Scalar& x)
420  {
421  EIGEN_USING_STD_MATH(arg);
422  return arg(x);
423  }
424  };
425 #else
426  template<typename Scalar, bool IsComplex = NumTraits<Scalar>::IsComplex>
427  struct arg_default_impl
428  {
429  typedef typename NumTraits<Scalar>::Real RealScalar;
430  EIGEN_DEVICE_FUNC
431  static inline RealScalar run(const Scalar& x)
432  {
433  return (x < Scalar(0)) ? Scalar(EIGEN_PI) : Scalar(0); }
434  };
435 
436  template<typename Scalar>
437  struct arg_default_impl<Scalar,true>
438  {
439  typedef typename NumTraits<Scalar>::Real RealScalar;
440  EIGEN_DEVICE_FUNC
441  static inline RealScalar run(const Scalar& x)
442  {
443  EIGEN_USING_STD_MATH(arg);
444  return arg(x);
445  }
446  };
447 
448  template<typename Scalar> struct arg_impl : arg_default_impl<Scalar> {};
449 #endif
450 
451 template<typename Scalar>
452 struct arg_retval
453 {
454  typedef typename NumTraits<Scalar>::Real type;
455 };
456 
457 /****************************************************************************
458 * Implementation of log1p *
459 ****************************************************************************/
460 template<typename Scalar, bool isComplex = NumTraits<Scalar>::IsComplex >
461 struct log1p_impl
462 {
463  static inline Scalar run(const Scalar& x)
464  {
465  EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar)
466  typedef typename NumTraits<Scalar>::Real RealScalar;
467  EIGEN_USING_STD_MATH(log);
468  Scalar x1p = RealScalar(1) + x;
469  return ( x1p == Scalar(1) ) ? x : x * ( log(x1p) / (x1p - RealScalar(1)) );
470  }
471 };
472 
473 #if EIGEN_HAS_CXX11_MATH
474 template<typename Scalar>
475 struct log1p_impl<Scalar, false> {
476  static inline Scalar run(const Scalar& x)
477  {
478  EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar)
479  using std::log1p;
480  return log1p(x);
481  }
482 };
483 #endif
484 
485 template<typename Scalar>
486 struct log1p_retval
487 {
488  typedef Scalar type;
489 };
490 
491 /****************************************************************************
492 * Implementation of pow *
493 ****************************************************************************/
494 
495 template<typename Scalar, bool IsInteger>
496 struct pow_default_impl
497 {
498  typedef Scalar retval;
499  static inline Scalar run(const Scalar& x, const Scalar& y)
500  {
501  EIGEN_USING_STD_MATH(pow);
502  return pow(x, y);
503  }
504 };
505 
506 template<typename Scalar>
507 struct pow_default_impl<Scalar, true>
508 {
509  static inline Scalar run(Scalar x, Scalar y)
510  {
511  Scalar res(1);
512  eigen_assert(!NumTraits<Scalar>::IsSigned || y >= 0);
513  if(y & 1) res *= x;
514  y >>= 1;
515  while(y)
516  {
517  x *= x;
518  if(y&1) res *= x;
519  y >>= 1;
520  }
521  return res;
522  }
523 };
524 
525 template<typename Scalar>
526 struct pow_impl : pow_default_impl<Scalar, NumTraits<Scalar>::IsInteger> {};
527 
528 template<typename Scalar>
529 struct pow_retval
530 {
531  typedef Scalar type;
532 };
533 
534 /****************************************************************************
535 * Implementation of random *
536 ****************************************************************************/
537 
538 template<typename Scalar,
539  bool IsComplex,
540  bool IsInteger>
541 struct random_default_impl {};
542 
543 template<typename Scalar>
544 struct random_impl : random_default_impl<Scalar, NumTraits<Scalar>::IsComplex, NumTraits<Scalar>::IsInteger> {};
545 
546 template<typename Scalar>
547 struct random_retval
548 {
549  typedef Scalar type;
550 };
551 
552 template<typename Scalar> inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random(const Scalar& x, const Scalar& y);
553 template<typename Scalar> inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random();
554 
555 template<typename Scalar>
556 struct random_default_impl<Scalar, false, false>
557 {
558  static inline Scalar run(const Scalar& x, const Scalar& y)
559  {
560  return x + (y-x) * Scalar(std::rand()) / Scalar(RAND_MAX);
561  }
562  static inline Scalar run()
563  {
564  return run(Scalar(NumTraits<Scalar>::IsSigned ? -1 : 0), Scalar(1));
565  }
566 };
567 
568 enum {
569  meta_floor_log2_terminate,
570  meta_floor_log2_move_up,
571  meta_floor_log2_move_down,
572  meta_floor_log2_bogus
573 };
574 
575 template<unsigned int n, int lower, int upper> struct meta_floor_log2_selector
576 {
577  enum { middle = (lower + upper) / 2,
578  value = (upper <= lower + 1) ? int(meta_floor_log2_terminate)
579  : (n < (1 << middle)) ? int(meta_floor_log2_move_down)
580  : (n==0) ? int(meta_floor_log2_bogus)
581  : int(meta_floor_log2_move_up)
582  };
583 };
584 
585 template<unsigned int n,
586  int lower = 0,
587  int upper = sizeof(unsigned int) * CHAR_BIT - 1,
588  int selector = meta_floor_log2_selector<n, lower, upper>::value>
589 struct meta_floor_log2 {};
590 
591 template<unsigned int n, int lower, int upper>
592 struct meta_floor_log2<n, lower, upper, meta_floor_log2_move_down>
593 {
594  enum { value = meta_floor_log2<n, lower, meta_floor_log2_selector<n, lower, upper>::middle>::value };
595 };
596 
597 template<unsigned int n, int lower, int upper>
598 struct meta_floor_log2<n, lower, upper, meta_floor_log2_move_up>
599 {
600  enum { value = meta_floor_log2<n, meta_floor_log2_selector<n, lower, upper>::middle, upper>::value };
601 };
602 
603 template<unsigned int n, int lower, int upper>
604 struct meta_floor_log2<n, lower, upper, meta_floor_log2_terminate>
605 {
606  enum { value = (n >= ((unsigned int)(1) << (lower+1))) ? lower+1 : lower };
607 };
608 
609 template<unsigned int n, int lower, int upper>
610 struct meta_floor_log2<n, lower, upper, meta_floor_log2_bogus>
611 {
612  // no value, error at compile time
613 };
614 
615 template<typename Scalar>
616 struct random_default_impl<Scalar, false, true>
617 {
618  static inline Scalar run(const Scalar& x, const Scalar& y)
619  {
620  typedef typename conditional<NumTraits<Scalar>::IsSigned,std::ptrdiff_t,std::size_t>::type ScalarX;
621  if(y<x)
622  return x;
623  std::size_t range = ScalarX(y)-ScalarX(x);
624  std::size_t offset = 0;
625  // rejection sampling
626  std::size_t divisor = (range+RAND_MAX-1)/(range+1);
627  std::size_t multiplier = (range+RAND_MAX-1)/std::size_t(RAND_MAX);
628 
629  do {
630  offset = ( (std::size_t(std::rand()) * multiplier) / divisor );
631  } while (offset > range);
632 
633  return Scalar(ScalarX(x) + offset);
634  }
635 
636  static inline Scalar run()
637  {
638 #ifdef EIGEN_MAKING_DOCS
639  return run(Scalar(NumTraits<Scalar>::IsSigned ? -10 : 0), Scalar(10));
640 #else
641  enum { rand_bits = meta_floor_log2<(unsigned int)(RAND_MAX)+1>::value,
642  scalar_bits = sizeof(Scalar) * CHAR_BIT,
643  shift = EIGEN_PLAIN_ENUM_MAX(0, int(rand_bits) - int(scalar_bits)),
644  offset = NumTraits<Scalar>::IsSigned ? (1 << (EIGEN_PLAIN_ENUM_MIN(rand_bits,scalar_bits)-1)) : 0
645  };
646  return Scalar((std::rand() >> shift) - offset);
647 #endif
648  }
649 };
650 
651 template<typename Scalar>
652 struct random_default_impl<Scalar, true, false>
653 {
654  static inline Scalar run(const Scalar& x, const Scalar& y)
655  {
656  return Scalar(random(real(x), real(y)),
657  random(imag(x), imag(y)));
658  }
659  static inline Scalar run()
660  {
661  typedef typename NumTraits<Scalar>::Real RealScalar;
662  return Scalar(random<RealScalar>(), random<RealScalar>());
663  }
664 };
665 
666 template<typename Scalar>
667 inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random(const Scalar& x, const Scalar& y)
668 {
669  return EIGEN_MATHFUNC_IMPL(random, Scalar)::run(x, y);
670 }
671 
672 template<typename Scalar>
673 inline EIGEN_MATHFUNC_RETVAL(random, Scalar) random()
674 {
675  return EIGEN_MATHFUNC_IMPL(random, Scalar)::run();
676 }
677 
678 // Implementatin of is* functions
679 
680 // std::is* do not work with fast-math and gcc, std::is* are available on MSVC 2013 and newer, as well as in clang.
681 #if (EIGEN_HAS_CXX11_MATH && !(EIGEN_COMP_GNUC_STRICT && __FINITE_MATH_ONLY__)) || (EIGEN_COMP_MSVC>=1800) || (EIGEN_COMP_CLANG)
682 #define EIGEN_USE_STD_FPCLASSIFY 1
683 #else
684 #define EIGEN_USE_STD_FPCLASSIFY 0
685 #endif
686 
687 template<typename T>
688 EIGEN_DEVICE_FUNC
689 typename internal::enable_if<internal::is_integral<T>::value,bool>::type
690 isnan_impl(const T&) { return false; }
691 
692 template<typename T>
693 EIGEN_DEVICE_FUNC
694 typename internal::enable_if<internal::is_integral<T>::value,bool>::type
695 isinf_impl(const T&) { return false; }
696 
697 template<typename T>
698 EIGEN_DEVICE_FUNC
699 typename internal::enable_if<internal::is_integral<T>::value,bool>::type
700 isfinite_impl(const T&) { return true; }
701 
702 template<typename T>
703 EIGEN_DEVICE_FUNC
704 typename internal::enable_if<(!internal::is_integral<T>::value)&&(!NumTraits<T>::IsComplex),bool>::type
705 isfinite_impl(const T& x)
706 {
707  #if EIGEN_USE_STD_FPCLASSIFY
708  using std::isfinite;
709  return isfinite EIGEN_NOT_A_MACRO (x);
710  #else
711  return x<NumTraits<T>::highest() && x>NumTraits<T>::lowest();
712  #endif
713 }
714 
715 template<typename T>
716 EIGEN_DEVICE_FUNC
717 typename internal::enable_if<(!internal::is_integral<T>::value)&&(!NumTraits<T>::IsComplex),bool>::type
718 isinf_impl(const T& x)
719 {
720  #if EIGEN_USE_STD_FPCLASSIFY
721  using std::isinf;
722  return isinf EIGEN_NOT_A_MACRO (x);
723  #else
724  return x>NumTraits<T>::highest() || x<NumTraits<T>::lowest();
725  #endif
726 }
727 
728 template<typename T>
729 EIGEN_DEVICE_FUNC
730 typename internal::enable_if<(!internal::is_integral<T>::value)&&(!NumTraits<T>::IsComplex),bool>::type
731 isnan_impl(const T& x)
732 {
733  #if EIGEN_USE_STD_FPCLASSIFY
734  using std::isnan;
735  return isnan EIGEN_NOT_A_MACRO (x);
736  #else
737  return x != x;
738  #endif
739 }
740 
741 #if (!EIGEN_USE_STD_FPCLASSIFY)
742 
743 #if EIGEN_COMP_MSVC
744 
745 template<typename T> EIGEN_DEVICE_FUNC bool isinf_msvc_helper(T x)
746 {
747  return _fpclass(x)==_FPCLASS_NINF || _fpclass(x)==_FPCLASS_PINF;
748 }
749 
750 //MSVC defines a _isnan builtin function, but for double only
751 EIGEN_DEVICE_FUNC inline bool isnan_impl(const long double& x) { return _isnan(x); }
752 EIGEN_DEVICE_FUNC inline bool isnan_impl(const double& x) { return _isnan(x); }
753 EIGEN_DEVICE_FUNC inline bool isnan_impl(const float& x) { return _isnan(x); }
754 
755 EIGEN_DEVICE_FUNC inline bool isinf_impl(const long double& x) { return isinf_msvc_helper(x); }
756 EIGEN_DEVICE_FUNC inline bool isinf_impl(const double& x) { return isinf_msvc_helper(x); }
757 EIGEN_DEVICE_FUNC inline bool isinf_impl(const float& x) { return isinf_msvc_helper(x); }
758 
759 #elif (defined __FINITE_MATH_ONLY__ && __FINITE_MATH_ONLY__ && EIGEN_COMP_GNUC)
760 
761 #if EIGEN_GNUC_AT_LEAST(5,0)
762  #define EIGEN_TMP_NOOPT_ATTRIB EIGEN_DEVICE_FUNC inline __attribute__((optimize("no-finite-math-only")))
763 #else
764  // NOTE the inline qualifier and noinline attribute are both needed: the former is to avoid linking issue (duplicate symbol),
765  // while the second prevent too aggressive optimizations in fast-math mode:
766  #define EIGEN_TMP_NOOPT_ATTRIB EIGEN_DEVICE_FUNC inline __attribute__((noinline,optimize("no-finite-math-only")))
767 #endif
768 
769 template<> EIGEN_TMP_NOOPT_ATTRIB bool isnan_impl(const long double& x) { return __builtin_isnan(x); }
770 template<> EIGEN_TMP_NOOPT_ATTRIB bool isnan_impl(const double& x) { return __builtin_isnan(x); }
771 template<> EIGEN_TMP_NOOPT_ATTRIB bool isnan_impl(const float& x) { return __builtin_isnan(x); }
772 template<> EIGEN_TMP_NOOPT_ATTRIB bool isinf_impl(const double& x) { return __builtin_isinf(x); }
773 template<> EIGEN_TMP_NOOPT_ATTRIB bool isinf_impl(const float& x) { return __builtin_isinf(x); }
774 template<> EIGEN_TMP_NOOPT_ATTRIB bool isinf_impl(const long double& x) { return __builtin_isinf(x); }
775 
776 #undef EIGEN_TMP_NOOPT_ATTRIB
777 
778 #endif
779 
780 #endif
781 
782 // The following overload are defined at the end of this file
783 template<typename T> bool isfinite_impl(const std::complex<T>& x);
784 template<typename T> bool isnan_impl(const std::complex<T>& x);
785 template<typename T> bool isinf_impl(const std::complex<T>& x);
786 
787 } // end namespace internal
788 
789 /****************************************************************************
790 * Generic math functions *
791 ****************************************************************************/
792 
793 namespace numext {
794 
795 #ifndef __CUDA_ARCH__
796 template<typename T>
797 EIGEN_DEVICE_FUNC
798 EIGEN_ALWAYS_INLINE T mini(const T& x, const T& y)
799 {
800  EIGEN_USING_STD_MATH(min);
801  return min EIGEN_NOT_A_MACRO (x,y);
802 }
803 
804 template<typename T>
805 EIGEN_DEVICE_FUNC
806 EIGEN_ALWAYS_INLINE T maxi(const T& x, const T& y)
807 {
808  EIGEN_USING_STD_MATH(max);
809  return max EIGEN_NOT_A_MACRO (x,y);
810 }
811 #else
812 template<typename T>
813 EIGEN_DEVICE_FUNC
814 EIGEN_ALWAYS_INLINE T mini(const T& x, const T& y)
815 {
816  return y < x ? y : x;
817 }
818 template<>
819 EIGEN_DEVICE_FUNC
820 EIGEN_ALWAYS_INLINE float mini(const float& x, const float& y)
821 {
822  return fmin(x, y);
823 }
824 template<typename T>
825 EIGEN_DEVICE_FUNC
826 EIGEN_ALWAYS_INLINE T maxi(const T& x, const T& y)
827 {
828  return x < y ? y : x;
829 }
830 template<>
831 EIGEN_DEVICE_FUNC
832 EIGEN_ALWAYS_INLINE float maxi(const float& x, const float& y)
833 {
834  return fmax(x, y);
835 }
836 #endif
837 
838 
839 template<typename Scalar>
840 EIGEN_DEVICE_FUNC
841 inline EIGEN_MATHFUNC_RETVAL(real, Scalar) real(const Scalar& x)
842 {
843  return EIGEN_MATHFUNC_IMPL(real, Scalar)::run(x);
844 }
845 
846 template<typename Scalar>
847 EIGEN_DEVICE_FUNC
848 inline typename internal::add_const_on_value_type< EIGEN_MATHFUNC_RETVAL(real_ref, Scalar) >::type real_ref(const Scalar& x)
849 {
850  return internal::real_ref_impl<Scalar>::run(x);
851 }
852 
853 template<typename Scalar>
854 EIGEN_DEVICE_FUNC
855 inline EIGEN_MATHFUNC_RETVAL(real_ref, Scalar) real_ref(Scalar& x)
856 {
857  return EIGEN_MATHFUNC_IMPL(real_ref, Scalar)::run(x);
858 }
859 
860 template<typename Scalar>
861 EIGEN_DEVICE_FUNC
862 inline EIGEN_MATHFUNC_RETVAL(imag, Scalar) imag(const Scalar& x)
863 {
864  return EIGEN_MATHFUNC_IMPL(imag, Scalar)::run(x);
865 }
866 
867 template<typename Scalar>
868 EIGEN_DEVICE_FUNC
869 inline EIGEN_MATHFUNC_RETVAL(arg, Scalar) arg(const Scalar& x)
870 {
871  return EIGEN_MATHFUNC_IMPL(arg, Scalar)::run(x);
872 }
873 
874 template<typename Scalar>
875 EIGEN_DEVICE_FUNC
876 inline typename internal::add_const_on_value_type< EIGEN_MATHFUNC_RETVAL(imag_ref, Scalar) >::type imag_ref(const Scalar& x)
877 {
878  return internal::imag_ref_impl<Scalar>::run(x);
879 }
880 
881 template<typename Scalar>
882 EIGEN_DEVICE_FUNC
883 inline EIGEN_MATHFUNC_RETVAL(imag_ref, Scalar) imag_ref(Scalar& x)
884 {
885  return EIGEN_MATHFUNC_IMPL(imag_ref, Scalar)::run(x);
886 }
887 
888 template<typename Scalar>
889 EIGEN_DEVICE_FUNC
890 inline EIGEN_MATHFUNC_RETVAL(conj, Scalar) conj(const Scalar& x)
891 {
892  return EIGEN_MATHFUNC_IMPL(conj, Scalar)::run(x);
893 }
894 
895 template<typename Scalar>
896 EIGEN_DEVICE_FUNC
897 inline EIGEN_MATHFUNC_RETVAL(abs2, Scalar) abs2(const Scalar& x)
898 {
899  return EIGEN_MATHFUNC_IMPL(abs2, Scalar)::run(x);
900 }
901 
902 template<typename Scalar>
903 EIGEN_DEVICE_FUNC
904 inline EIGEN_MATHFUNC_RETVAL(norm1, Scalar) norm1(const Scalar& x)
905 {
906  return EIGEN_MATHFUNC_IMPL(norm1, Scalar)::run(x);
907 }
908 
909 template<typename Scalar>
910 EIGEN_DEVICE_FUNC
911 inline EIGEN_MATHFUNC_RETVAL(hypot, Scalar) hypot(const Scalar& x, const Scalar& y)
912 {
913  return EIGEN_MATHFUNC_IMPL(hypot, Scalar)::run(x, y);
914 }
915 
916 template<typename Scalar>
917 EIGEN_DEVICE_FUNC
918 inline EIGEN_MATHFUNC_RETVAL(log1p, Scalar) log1p(const Scalar& x)
919 {
920  return EIGEN_MATHFUNC_IMPL(log1p, Scalar)::run(x);
921 }
922 
923 template<typename Scalar>
924 EIGEN_DEVICE_FUNC
925 inline EIGEN_MATHFUNC_RETVAL(pow, Scalar) pow(const Scalar& x, const Scalar& y)
926 {
927  return EIGEN_MATHFUNC_IMPL(pow, Scalar)::run(x, y);
928 }
929 
930 template<typename T> EIGEN_DEVICE_FUNC bool (isnan) (const T &x) { return internal::isnan_impl(x); }
931 template<typename T> EIGEN_DEVICE_FUNC bool (isinf) (const T &x) { return internal::isinf_impl(x); }
932 template<typename T> EIGEN_DEVICE_FUNC bool (isfinite)(const T &x) { return internal::isfinite_impl(x); }
933 
934 template<typename Scalar>
935 EIGEN_DEVICE_FUNC
936 inline EIGEN_MATHFUNC_RETVAL(round, Scalar) round(const Scalar& x)
937 {
938  return EIGEN_MATHFUNC_IMPL(round, Scalar)::run(x);
939 }
940 
941 template<typename T>
942 EIGEN_DEVICE_FUNC
943 T (floor)(const T& x)
944 {
945  EIGEN_USING_STD_MATH(floor);
946  return floor(x);
947 }
948 
949 template<typename T>
950 EIGEN_DEVICE_FUNC
951 T (ceil)(const T& x)
952 {
953  EIGEN_USING_STD_MATH(ceil);
954  return ceil(x);
955 }
956 
957 // Log base 2 for 32 bits positive integers.
958 // Conveniently returns 0 for x==0.
959 inline int log2(int x)
960 {
961  eigen_assert(x>=0);
962  unsigned int v(x);
963  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 };
964  v |= v >> 1;
965  v |= v >> 2;
966  v |= v >> 4;
967  v |= v >> 8;
968  v |= v >> 16;
969  return table[(v * 0x07C4ACDDU) >> 27];
970 }
971 
972 } // end namespace numext
973 
974 namespace internal {
975 
976 template<typename T>
977 bool isfinite_impl(const std::complex<T>& x)
978 {
979  return (numext::isfinite)(numext::real(x)) && (numext::isfinite)(numext::imag(x));
980 }
981 
982 template<typename T>
983 bool isnan_impl(const std::complex<T>& x)
984 {
985  return (numext::isnan)(numext::real(x)) || (numext::isnan)(numext::imag(x));
986 }
987 
988 template<typename T>
989 bool isinf_impl(const std::complex<T>& x)
990 {
991  return ((numext::isinf)(numext::real(x)) || (numext::isinf)(numext::imag(x))) && (!(numext::isnan)(x));
992 }
993 
994 /****************************************************************************
995 * Implementation of fuzzy comparisons *
996 ****************************************************************************/
997 
998 template<typename Scalar,
999  bool IsComplex,
1000  bool IsInteger>
1001 struct scalar_fuzzy_default_impl {};
1002 
1003 template<typename Scalar>
1004 struct scalar_fuzzy_default_impl<Scalar, false, false>
1005 {
1006  typedef typename NumTraits<Scalar>::Real RealScalar;
1007  template<typename OtherScalar> EIGEN_DEVICE_FUNC
1008  static inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y, const RealScalar& prec)
1009  {
1010  EIGEN_USING_STD_MATH(abs);
1011  return abs(x) <= abs(y) * prec;
1012  }
1013  EIGEN_DEVICE_FUNC
1014  static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar& prec)
1015  {
1016  EIGEN_USING_STD_MATH(abs);
1017  return abs(x - y) <= numext::mini(abs(x), abs(y)) * prec;
1018  }
1019  EIGEN_DEVICE_FUNC
1020  static inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y, const RealScalar& prec)
1021  {
1022  return x <= y || isApprox(x, y, prec);
1023  }
1024 };
1025 
1026 template<typename Scalar>
1027 struct scalar_fuzzy_default_impl<Scalar, false, true>
1028 {
1029  typedef typename NumTraits<Scalar>::Real RealScalar;
1030  template<typename OtherScalar> EIGEN_DEVICE_FUNC
1031  static inline bool isMuchSmallerThan(const Scalar& x, const Scalar&, const RealScalar&)
1032  {
1033  return x == Scalar(0);
1034  }
1035  EIGEN_DEVICE_FUNC
1036  static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar&)
1037  {
1038  return x == y;
1039  }
1040  EIGEN_DEVICE_FUNC
1041  static inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y, const RealScalar&)
1042  {
1043  return x <= y;
1044  }
1045 };
1046 
1047 template<typename Scalar>
1048 struct scalar_fuzzy_default_impl<Scalar, true, false>
1049 {
1050  typedef typename NumTraits<Scalar>::Real RealScalar;
1051  template<typename OtherScalar>
1052  static inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y, const RealScalar& prec)
1053  {
1054  return numext::abs2(x) <= numext::abs2(y) * prec * prec;
1055  }
1056  static inline bool isApprox(const Scalar& x, const Scalar& y, const RealScalar& prec)
1057  {
1058  return numext::abs2(x - y) <= numext::mini(numext::abs2(x), numext::abs2(y)) * prec * prec;
1059  }
1060 };
1061 
1062 template<typename Scalar>
1063 struct scalar_fuzzy_impl : scalar_fuzzy_default_impl<Scalar, NumTraits<Scalar>::IsComplex, NumTraits<Scalar>::IsInteger> {};
1064 
1065 template<typename Scalar, typename OtherScalar> EIGEN_DEVICE_FUNC
1066 inline bool isMuchSmallerThan(const Scalar& x, const OtherScalar& y,
1067  typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision())
1068 {
1069  return scalar_fuzzy_impl<Scalar>::template isMuchSmallerThan<OtherScalar>(x, y, precision);
1070 }
1071 
1072 template<typename Scalar> EIGEN_DEVICE_FUNC
1073 inline bool isApprox(const Scalar& x, const Scalar& y,
1074  typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision())
1075 {
1076  return scalar_fuzzy_impl<Scalar>::isApprox(x, y, precision);
1077 }
1078 
1079 template<typename Scalar> EIGEN_DEVICE_FUNC
1080 inline bool isApproxOrLessThan(const Scalar& x, const Scalar& y,
1081  typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision())
1082 {
1083  return scalar_fuzzy_impl<Scalar>::isApproxOrLessThan(x, y, precision);
1084 }
1085 
1086 /******************************************
1087 *** The special case of the bool type ***
1088 ******************************************/
1089 
1090 template<> struct random_impl<bool>
1091 {
1092  static inline bool run()
1093  {
1094  return random<int>(0,1)==0 ? false : true;
1095  }
1096 };
1097 
1098 template<> struct scalar_fuzzy_impl<bool>
1099 {
1100  typedef bool RealScalar;
1101 
1102  template<typename OtherScalar> EIGEN_DEVICE_FUNC
1103  static inline bool isMuchSmallerThan(const bool& x, const bool&, const bool&)
1104  {
1105  return !x;
1106  }
1107 
1108  EIGEN_DEVICE_FUNC
1109  static inline bool isApprox(bool x, bool y, bool)
1110  {
1111  return x == y;
1112  }
1113 
1114  EIGEN_DEVICE_FUNC
1115  static inline bool isApproxOrLessThan(const bool& x, const bool& y, const bool&)
1116  {
1117  return (!x) || y;
1118  }
1119 
1120 };
1121 
1122 
1123 } // end namespace internal
1124 
1125 } // end namespace Eigen
1126 
1127 #endif // EIGEN_MATHFUNCTIONS_H
Definition: LDLT.h:16
Definition: StdDeque.h:58
Definition: Eigen_Colamd.h:54