Eigen  3.2.93
Half.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // This Source Code Form is subject to the terms of the Mozilla
5 // Public License v. 2.0. If a copy of the MPL was not distributed
6 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
7 //
8 // The conversion routines are Copyright (c) Fabian Giesen, 2016.
9 // The original license follows:
10 //
11 // Copyright (c) Fabian Giesen, 2016
12 // All rights reserved.
13 // Redistribution and use in source and binary forms, with or without
14 // modification, are permitted.
15 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
16 // “AS IS” AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
17 // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
18 // A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
19 // HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
20 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
21 // LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
22 // DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
23 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
24 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
25 // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26 
27 
28 // Standard 16-bit float type, mostly useful for GPUs. Defines a new
29 // type Eigen::half (inheriting from CUDA's __half struct) with
30 // operator overloads such that it behaves basically as an arithmetic
31 // type. It will be quite slow on CPUs (so it is recommended to stay
32 // in fp32 for CPUs, except for simple parameter conversions, I/O
33 // to disk and the likes), but fast on GPUs.
34 
35 
36 #ifndef EIGEN_HALF_CUDA_H
37 #define EIGEN_HALF_CUDA_H
38 
39 #if __cplusplus > 199711L
40 #define EIGEN_EXPLICIT_CAST(tgt_type) explicit operator tgt_type()
41 #else
42 #define EIGEN_EXPLICIT_CAST(tgt_type) operator tgt_type()
43 #endif
44 
45 
46 namespace Eigen {
47 
48 namespace half_impl {
49 
50 #if !defined(EIGEN_HAS_CUDA_FP16)
51 
52 // Make our own __half definition that is similar to CUDA's.
53 struct __half {
54  EIGEN_DEVICE_FUNC __half() {}
55  explicit EIGEN_DEVICE_FUNC __half(unsigned short raw) : x(raw) {}
56  unsigned short x;
57 };
58 
59 #endif
60 
61 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC __half raw_uint16_to_half(unsigned short x);
62 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC __half float_to_half_rtne(float ff);
63 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC float half_to_float(__half h);
64 
65 // Class definition.
66 struct half : public __half {
67  EIGEN_DEVICE_FUNC half() {}
68 
69  EIGEN_DEVICE_FUNC half(const __half& h) : __half(h) {}
70  EIGEN_DEVICE_FUNC half(const half& h) : __half(h) {}
71 
72  explicit EIGEN_DEVICE_FUNC half(bool b)
73  : __half(raw_uint16_to_half(b ? 0x3c00 : 0)) {}
74  template<class T>
75  explicit EIGEN_DEVICE_FUNC half(const T& val)
76  : __half(float_to_half_rtne(static_cast<float>(val))) {}
77  explicit EIGEN_DEVICE_FUNC half(float f)
78  : __half(float_to_half_rtne(f)) {}
79 
80  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(bool) const {
81  // +0.0 and -0.0 become false, everything else becomes true.
82  return (x & 0x7fff) != 0;
83  }
84  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(signed char) const {
85  return static_cast<signed char>(half_to_float(*this));
86  }
87  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(unsigned char) const {
88  return static_cast<unsigned char>(half_to_float(*this));
89  }
90  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(short) const {
91  return static_cast<short>(half_to_float(*this));
92  }
93  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(unsigned short) const {
94  return static_cast<unsigned short>(half_to_float(*this));
95  }
96  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(int) const {
97  return static_cast<int>(half_to_float(*this));
98  }
99  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(unsigned int) const {
100  return static_cast<unsigned int>(half_to_float(*this));
101  }
102  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(long) const {
103  return static_cast<long>(half_to_float(*this));
104  }
105  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(unsigned long) const {
106  return static_cast<unsigned long>(half_to_float(*this));
107  }
108  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(long long) const {
109  return static_cast<long long>(half_to_float(*this));
110  }
111  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(unsigned long long) const {
112  return static_cast<unsigned long long>(half_to_float(*this));
113  }
114  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(float) const {
115  return half_to_float(*this);
116  }
117  EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(double) const {
118  return static_cast<double>(half_to_float(*this));
119  }
120 
121  EIGEN_DEVICE_FUNC half& operator=(const half& other) {
122  x = other.x;
123  return *this;
124  }
125 };
126 
127 #if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 530
128 
129 // Intrinsics for native fp16 support. Note that on current hardware,
130 // these are no faster than fp32 arithmetic (you need to use the half2
131 // versions to get the ALU speed increased), but you do save the
132 // conversion steps back and forth.
133 
134 __device__ half operator + (const half& a, const half& b) {
135  return __hadd(a, b);
136 }
137 __device__ half operator * (const half& a, const half& b) {
138  return __hmul(a, b);
139 }
140 __device__ half operator - (const half& a, const half& b) {
141  return __hsub(a, b);
142 }
143 __device__ half operator / (const half& a, const half& b) {
144  float num = __half2float(a);
145  float denom = __half2float(b);
146  return __float2half(num / denom);
147 }
148 __device__ half operator - (const half& a) {
149  return __hneg(a);
150 }
151 __device__ half& operator += (half& a, const half& b) {
152  a = a + b;
153  return a;
154 }
155 __device__ half& operator *= (half& a, const half& b) {
156  a = a * b;
157  return a;
158 }
159 __device__ half& operator -= (half& a, const half& b) {
160  a = a - b;
161  return a;
162 }
163 __device__ half& operator /= (half& a, const half& b) {
164  a = a / b;
165  return a;
166 }
167 __device__ bool operator == (const half& a, const half& b) {
168  return __heq(a, b);
169 }
170 __device__ bool operator != (const half& a, const half& b) {
171  return __hne(a, b);
172 }
173 __device__ bool operator < (const half& a, const half& b) {
174  return __hlt(a, b);
175 }
176 __device__ bool operator <= (const half& a, const half& b) {
177  return __hle(a, b);
178 }
179 __device__ bool operator > (const half& a, const half& b) {
180  return __hgt(a, b);
181 }
182 __device__ bool operator >= (const half& a, const half& b) {
183  return __hge(a, b);
184 }
185 
186 #else // Emulate support for half floats
187 
188 // Definitions for CPUs and older CUDA, mostly working through conversion
189 // to/from fp32.
190 
191 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator + (const half& a, const half& b) {
192  return half(float(a) + float(b));
193 }
194 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator * (const half& a, const half& b) {
195  return half(float(a) * float(b));
196 }
197 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator - (const half& a, const half& b) {
198  return half(float(a) - float(b));
199 }
200 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator / (const half& a, const half& b) {
201  return half(float(a) / float(b));
202 }
203 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator - (const half& a) {
204  half result;
205  result.x = a.x ^ 0x8000;
206  return result;
207 }
208 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half& operator += (half& a, const half& b) {
209  a = half(float(a) + float(b));
210  return a;
211 }
212 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half& operator *= (half& a, const half& b) {
213  a = half(float(a) * float(b));
214  return a;
215 }
216 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half& operator -= (half& a, const half& b) {
217  a = half(float(a) - float(b));
218  return a;
219 }
220 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half& operator /= (half& a, const half& b) {
221  a = half(float(a) / float(b));
222  return a;
223 }
224 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator == (const half& a, const half& b) {
225  return float(a) == float(b);
226 }
227 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator != (const half& a, const half& b) {
228  return float(a) != float(b);
229 }
230 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator < (const half& a, const half& b) {
231  return float(a) < float(b);
232 }
233 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator <= (const half& a, const half& b) {
234  return float(a) <= float(b);
235 }
236 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator > (const half& a, const half& b) {
237  return float(a) > float(b);
238 }
239 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator >= (const half& a, const half& b) {
240  return float(a) >= float(b);
241 }
242 
243 #endif // Emulate support for half floats
244 
245 // Division by an index. Do it in full float precision to avoid accuracy
246 // issues in converting the denominator to half.
247 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator / (const half& a, Index b) {
248  return half(static_cast<float>(a) / static_cast<float>(b));
249 }
250 
251 // Conversion routines, including fallbacks for the host or older CUDA.
252 // Note that newer Intel CPUs (Haswell or newer) have vectorized versions of
253 // these in hardware. If we need more performance on older/other CPUs, they are
254 // also possible to vectorize directly.
255 
256 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC __half raw_uint16_to_half(unsigned short x) {
257  __half h;
258  h.x = x;
259  return h;
260 }
261 
262 union FP32 {
263  unsigned int u;
264  float f;
265 };
266 
267 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC __half float_to_half_rtne(float ff) {
268 #if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 300
269  return __float2half(ff);
270 
271 #elif defined(EIGEN_HAS_FP16_C)
272  __half h;
273  h.x = _cvtss_sh(ff, 0);
274  return h;
275 
276 #else
277  FP32 f; f.f = ff;
278 
279  const FP32 f32infty = { 255 << 23 };
280  const FP32 f16max = { (127 + 16) << 23 };
281  const FP32 denorm_magic = { ((127 - 15) + (23 - 10) + 1) << 23 };
282  unsigned int sign_mask = 0x80000000u;
283  __half o;
284  o.x = static_cast<unsigned short>(0x0u);
285 
286  unsigned int sign = f.u & sign_mask;
287  f.u ^= sign;
288 
289  // NOTE all the integer compares in this function can be safely
290  // compiled into signed compares since all operands are below
291  // 0x80000000. Important if you want fast straight SSE2 code
292  // (since there's no unsigned PCMPGTD).
293 
294  if (f.u >= f16max.u) { // result is Inf or NaN (all exponent bits set)
295  o.x = (f.u > f32infty.u) ? 0x7e00 : 0x7c00; // NaN->qNaN and Inf->Inf
296  } else { // (De)normalized number or zero
297  if (f.u < (113 << 23)) { // resulting FP16 is subnormal or zero
298  // use a magic value to align our 10 mantissa bits at the bottom of
299  // the float. as long as FP addition is round-to-nearest-even this
300  // just works.
301  f.f += denorm_magic.f;
302 
303  // and one integer subtract of the bias later, we have our final float!
304  o.x = static_cast<unsigned short>(f.u - denorm_magic.u);
305  } else {
306  unsigned int mant_odd = (f.u >> 13) & 1; // resulting mantissa is odd
307 
308  // update exponent, rounding bias part 1
309  f.u += ((unsigned int)(15 - 127) << 23) + 0xfff;
310  // rounding bias part 2
311  f.u += mant_odd;
312  // take the bits!
313  o.x = static_cast<unsigned short>(f.u >> 13);
314  }
315  }
316 
317  o.x |= static_cast<unsigned short>(sign >> 16);
318  return o;
319 #endif
320 }
321 
322 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC float half_to_float(__half h) {
323 #if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 300
324  return __half2float(h);
325 
326 #elif defined(EIGEN_HAS_FP16_C)
327  return _cvtsh_ss(h.x);
328 
329 #else
330  const FP32 magic = { 113 << 23 };
331  const unsigned int shifted_exp = 0x7c00 << 13; // exponent mask after shift
332  FP32 o;
333 
334  o.u = (h.x & 0x7fff) << 13; // exponent/mantissa bits
335  unsigned int exp = shifted_exp & o.u; // just the exponent
336  o.u += (127 - 15) << 23; // exponent adjust
337 
338  // handle exponent special cases
339  if (exp == shifted_exp) { // Inf/NaN?
340  o.u += (128 - 16) << 23; // extra exp adjust
341  } else if (exp == 0) { // Zero/Denormal?
342  o.u += 1 << 23; // extra exp adjust
343  o.f -= magic.f; // renormalize
344  }
345 
346  o.u |= (h.x & 0x8000) << 16; // sign bit
347  return o.f;
348 #endif
349 }
350 
351 // --- standard functions ---
352 
353 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool (isinf)(const half& a) {
354  return (a.x & 0x7fff) == 0x7c00;
355 }
356 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool (isnan)(const half& a) {
357 #if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 530
358  return __hisnan(a);
359 #else
360  return (a.x & 0x7fff) > 0x7c00;
361 #endif
362 }
363 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool (isfinite)(const half& a) {
364  return !(isinf EIGEN_NOT_A_MACRO (a)) && !(isnan EIGEN_NOT_A_MACRO (a));
365 }
366 
367 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half abs(const half& a) {
368  half result;
369  result.x = a.x & 0x7FFF;
370  return result;
371 }
372 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half exp(const half& a) {
373  return half(::expf(float(a)));
374 }
375 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half log(const half& a) {
376  return half(::logf(float(a)));
377 }
378 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half log10(const half& a) {
379  return half(::log10f(float(a)));
380 }
381 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half sqrt(const half& a) {
382  return half(::sqrtf(float(a)));
383 }
384 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half pow(const half& a, const half& b) {
385  return half(::powf(float(a), float(b)));
386 }
387 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half sin(const half& a) {
388  return half(::sinf(float(a)));
389 }
390 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half cos(const half& a) {
391  return half(::cosf(float(a)));
392 }
393 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half tan(const half& a) {
394  return half(::tanf(float(a)));
395 }
396 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half tanh(const half& a) {
397  return half(::tanhf(float(a)));
398 }
399 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half floor(const half& a) {
400  return half(::floorf(float(a)));
401 }
402 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half ceil(const half& a) {
403  return half(::ceilf(float(a)));
404 }
405 
406 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half (min)(const half& a, const half& b) {
407 #if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 530
408  return __hlt(b, a) ? b : a;
409 #else
410  const float f1 = static_cast<float>(a);
411  const float f2 = static_cast<float>(b);
412  return f2 < f1 ? b : a;
413 #endif
414 }
415 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half (max)(const half& a, const half& b) {
416 #if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 530
417  return __hlt(a, b) ? b : a;
418 #else
419  const float f1 = static_cast<float>(a);
420  const float f2 = static_cast<float>(b);
421  return f1 < f2 ? b : a;
422 #endif
423 }
424 
425 EIGEN_ALWAYS_INLINE std::ostream& operator << (std::ostream& os, const half& v) {
426  os << static_cast<float>(v);
427  return os;
428 }
429 
430 } // end namespace half_impl
431 
432 // import Eigen::half_impl::half into Eigen namespace
433 using half_impl::half;
434 
435 namespace internal {
436 
437 template<>
438 struct random_default_impl<half_impl::half, false, false>
439 {
440  static inline half run(const half& x, const half& y)
441  {
442  return x + (y-x) * half(float(std::rand()) / float(RAND_MAX));
443  }
444  static inline half run()
445  {
446  return run(half(-1.f), half(1.f));
447  }
448 };
449 
450 template<> struct is_arithmetic<half_impl::half> { enum { value = true }; };
451 
452 } // end namespace internal
453 
454 template<> struct NumTraits<Eigen::half_impl::half>
455  : GenericNumTraits<Eigen::half_impl::half>
456 {
457  EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half_impl::half epsilon() {
458  return half_impl::raw_uint16_to_half(0x0800);
459  }
460  EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half_impl::half dummy_precision() { return half_impl::half(1e-2f); }
461  EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half_impl::half highest() {
462  return half_impl::raw_uint16_to_half(0x7bff);
463  }
464  EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half_impl::half lowest() {
465  return half_impl::raw_uint16_to_half(0xfbff);
466  }
467  EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half_impl::half infinity() {
468  return half_impl::raw_uint16_to_half(0x7c00);
469  }
470  EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half_impl::half quiet_NaN() {
471  return half_impl::raw_uint16_to_half(0x7c01);
472  }
473 };
474 
475 } // end namespace Eigen
476 
477 // C-like standard mathematical functions and trancendentals.
478 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half fabsh(const Eigen::half& a) {
479  Eigen::half result;
480  result.x = a.x & 0x7FFF;
481  return result;
482 }
483 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half exph(const Eigen::half& a) {
484  return Eigen::half(::expf(float(a)));
485 }
486 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half logh(const Eigen::half& a) {
487  return Eigen::half(::logf(float(a)));
488 }
489 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half sqrth(const Eigen::half& a) {
490  return Eigen::half(::sqrtf(float(a)));
491 }
492 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half powh(const Eigen::half& a, const Eigen::half& b) {
493  return Eigen::half(::powf(float(a), float(b)));
494 }
495 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half floorh(const Eigen::half& a) {
496  return Eigen::half(::floorf(float(a)));
497 }
498 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half ceilh(const Eigen::half& a) {
499  return Eigen::half(::ceilf(float(a)));
500 }
501 
502 namespace std {
503 
504 #if __cplusplus > 199711L
505 template <>
506 struct hash<Eigen::half> {
507  EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE std::size_t operator()(const Eigen::half& a) const {
508  return static_cast<std::size_t>(a.x);
509  }
510 };
511 #endif
512 
513 } // end namespace std
514 
515 
516 // Add the missing shfl_xor intrinsic
517 #if defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 300
518 __device__ EIGEN_STRONG_INLINE Eigen::half __shfl_xor(Eigen::half var, int laneMask, int width=warpSize) {
519  return static_cast<Eigen::half>(__shfl_xor(static_cast<float>(var), laneMask, width));
520 }
521 #endif
522 
523 // ldg() has an overload for __half, but we also need one for Eigen::half.
524 #if defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 350
525 EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half __ldg(const Eigen::half* ptr) {
526  return Eigen::internal::raw_uint16_to_half(
527  __ldg(reinterpret_cast<const unsigned short*>(ptr)));
528 }
529 #endif
530 
531 #endif // EIGEN_HALF_CUDA_H
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_tanh_op< typename Derived::Scalar >, const Derived > tanh(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_isfinite_op< typename Derived::Scalar >, const Derived > isfinite(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_sqrt_op< typename Derived::Scalar >, const Derived > sqrt(const Eigen::ArrayBase< Derived > &x)
Namespace containing all symbols from the Eigen library.
Definition: Core:271
Definition: Half.h:502
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_ceil_op< typename Derived::Scalar >, const Derived > ceil(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_isnan_op< typename Derived::Scalar >, const Derived > isnan(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_cos_op< typename Derived::Scalar >, const Derived > cos(const Eigen::ArrayBase< Derived > &x)
const Product< MatrixDerived, PermutationDerived, AliasFreeProduct > operator*(const MatrixBase< MatrixDerived > &matrix, const PermutationBase< PermutationDerived > &permutation)
Definition: PermutationMatrix.h:543
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: XprHelper.h:35
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_floor_op< typename Derived::Scalar >, const Derived > floor(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_isinf_op< typename Derived::Scalar >, const Derived > isinf(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_abs_op< typename Derived::Scalar >, const Derived > abs(const Eigen::ArrayBase< Derived > &x)
Definition: Eigen_Colamd.h:50
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_log_op< typename Derived::Scalar >, const Derived > log(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_tan_op< typename Derived::Scalar >, const Derived > tan(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_sign_op< typename Derived::Scalar >, const Derived > sign(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_sin_op< typename Derived::Scalar >, const Derived > sin(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_exp_op< typename Derived::Scalar >, const Derived > exp(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_log10_op< typename Derived::Scalar >, const Derived > log10(const Eigen::ArrayBase< Derived > &x)