Eigen  3.3.0
 
Loading...
Searching...
No Matches
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
46namespace Eigen {
47
48struct half;
49
50namespace half_impl {
51
52#if !defined(EIGEN_HAS_CUDA_FP16)
53
54// Make our own __half definition that is similar to CUDA's.
55struct __half {
56 EIGEN_DEVICE_FUNC __half() {}
57 explicit EIGEN_DEVICE_FUNC __half(unsigned short raw) : x(raw) {}
58 unsigned short x;
59};
60
61#endif
62
63EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC __half raw_uint16_to_half(unsigned short x);
64EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC __half float_to_half_rtne(float ff);
65EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC float half_to_float(__half h);
66
67struct half_base : public __half {
68 EIGEN_DEVICE_FUNC half_base() {}
69 EIGEN_DEVICE_FUNC half_base(const half_base& h) : __half(h) {}
70 EIGEN_DEVICE_FUNC half_base(const __half& h) : __half(h) {}
71};
72
73} // namespace half_impl
74
75// Class definition.
76struct half : public half_impl::half_base {
77 #if !defined(EIGEN_HAS_CUDA_FP16)
78 typedef half_impl::__half __half;
79 #endif
80
81 EIGEN_DEVICE_FUNC half() {}
82
83 EIGEN_DEVICE_FUNC half(const __half& h) : half_impl::half_base(h) {}
84 EIGEN_DEVICE_FUNC half(const half& h) : half_impl::half_base(h) {}
85
86 explicit EIGEN_DEVICE_FUNC half(bool b)
87 : half_impl::half_base(half_impl::raw_uint16_to_half(b ? 0x3c00 : 0)) {}
88 template<class T>
89 explicit EIGEN_DEVICE_FUNC half(const T& val)
90 : half_impl::half_base(half_impl::float_to_half_rtne(static_cast<float>(val))) {}
91 explicit EIGEN_DEVICE_FUNC half(float f)
92 : half_impl::half_base(half_impl::float_to_half_rtne(f)) {}
93
94 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(bool) const {
95 // +0.0 and -0.0 become false, everything else becomes true.
96 return (x & 0x7fff) != 0;
97 }
98 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(signed char) const {
99 return static_cast<signed char>(half_impl::half_to_float(*this));
100 }
101 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(unsigned char) const {
102 return static_cast<unsigned char>(half_impl::half_to_float(*this));
103 }
104 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(short) const {
105 return static_cast<short>(half_impl::half_to_float(*this));
106 }
107 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(unsigned short) const {
108 return static_cast<unsigned short>(half_impl::half_to_float(*this));
109 }
110 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(int) const {
111 return static_cast<int>(half_impl::half_to_float(*this));
112 }
113 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(unsigned int) const {
114 return static_cast<unsigned int>(half_impl::half_to_float(*this));
115 }
116 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(long) const {
117 return static_cast<long>(half_impl::half_to_float(*this));
118 }
119 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(unsigned long) const {
120 return static_cast<unsigned long>(half_impl::half_to_float(*this));
121 }
122 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(long long) const {
123 return static_cast<long long>(half_impl::half_to_float(*this));
124 }
125 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(unsigned long long) const {
126 return static_cast<unsigned long long>(half_to_float(*this));
127 }
128 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(float) const {
129 return half_impl::half_to_float(*this);
130 }
131 EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(double) const {
132 return static_cast<double>(half_impl::half_to_float(*this));
133 }
134
135 EIGEN_DEVICE_FUNC half& operator=(const half& other) {
136 x = other.x;
137 return *this;
138 }
139};
140
141namespace half_impl {
142
143#if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 530
144
145// Intrinsics for native fp16 support. Note that on current hardware,
146// these are no faster than fp32 arithmetic (you need to use the half2
147// versions to get the ALU speed increased), but you do save the
148// conversion steps back and forth.
149
150__device__ half operator + (const half& a, const half& b) {
151 return __hadd(a, b);
152}
153__device__ half operator * (const half& a, const half& b) {
154 return __hmul(a, b);
155}
156__device__ half operator - (const half& a, const half& b) {
157 return __hsub(a, b);
158}
159__device__ half operator / (const half& a, const half& b) {
160 float num = __half2float(a);
161 float denom = __half2float(b);
162 return __float2half(num / denom);
163}
164__device__ half operator - (const half& a) {
165 return __hneg(a);
166}
167__device__ half& operator += (half& a, const half& b) {
168 a = a + b;
169 return a;
170}
171__device__ half& operator *= (half& a, const half& b) {
172 a = a * b;
173 return a;
174}
175__device__ half& operator -= (half& a, const half& b) {
176 a = a - b;
177 return a;
178}
179__device__ half& operator /= (half& a, const half& b) {
180 a = a / b;
181 return a;
182}
183__device__ bool operator == (const half& a, const half& b) {
184 return __heq(a, b);
185}
186__device__ bool operator != (const half& a, const half& b) {
187 return __hne(a, b);
188}
189__device__ bool operator < (const half& a, const half& b) {
190 return __hlt(a, b);
191}
192__device__ bool operator <= (const half& a, const half& b) {
193 return __hle(a, b);
194}
195__device__ bool operator > (const half& a, const half& b) {
196 return __hgt(a, b);
197}
198__device__ bool operator >= (const half& a, const half& b) {
199 return __hge(a, b);
200}
201
202#else // Emulate support for half floats
203
204// Definitions for CPUs and older CUDA, mostly working through conversion
205// to/from fp32.
206
207EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator + (const half& a, const half& b) {
208 return half(float(a) + float(b));
209}
210EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator * (const half& a, const half& b) {
211 return half(float(a) * float(b));
212}
213EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator - (const half& a, const half& b) {
214 return half(float(a) - float(b));
215}
216EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator / (const half& a, const half& b) {
217 return half(float(a) / float(b));
218}
219EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator - (const half& a) {
220 half result;
221 result.x = a.x ^ 0x8000;
222 return result;
223}
224EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half& operator += (half& a, const half& b) {
225 a = half(float(a) + float(b));
226 return a;
227}
228EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half& operator *= (half& a, const half& b) {
229 a = half(float(a) * float(b));
230 return a;
231}
232EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half& operator -= (half& a, const half& b) {
233 a = half(float(a) - float(b));
234 return a;
235}
236EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half& operator /= (half& a, const half& b) {
237 a = half(float(a) / float(b));
238 return a;
239}
240EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator == (const half& a, const half& b) {
241 return float(a) == float(b);
242}
243EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator != (const half& a, const half& b) {
244 return float(a) != float(b);
245}
246EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator < (const half& a, const half& b) {
247 return float(a) < float(b);
248}
249EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator <= (const half& a, const half& b) {
250 return float(a) <= float(b);
251}
252EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator > (const half& a, const half& b) {
253 return float(a) > float(b);
254}
255EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool operator >= (const half& a, const half& b) {
256 return float(a) >= float(b);
257}
258
259#endif // Emulate support for half floats
260
261// Division by an index. Do it in full float precision to avoid accuracy
262// issues in converting the denominator to half.
263EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half operator / (const half& a, Index b) {
264 return half(static_cast<float>(a) / static_cast<float>(b));
265}
266
267// Conversion routines, including fallbacks for the host or older CUDA.
268// Note that newer Intel CPUs (Haswell or newer) have vectorized versions of
269// these in hardware. If we need more performance on older/other CPUs, they are
270// also possible to vectorize directly.
271
272EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC __half raw_uint16_to_half(unsigned short x) {
273 __half h;
274 h.x = x;
275 return h;
276}
277
278union FP32 {
279 unsigned int u;
280 float f;
281};
282
283EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC __half float_to_half_rtne(float ff) {
284#if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 300
285 return __float2half(ff);
286
287#elif defined(EIGEN_HAS_FP16_C)
288 __half h;
289 h.x = _cvtss_sh(ff, 0);
290 return h;
291
292#else
293 FP32 f; f.f = ff;
294
295 const FP32 f32infty = { 255 << 23 };
296 const FP32 f16max = { (127 + 16) << 23 };
297 const FP32 denorm_magic = { ((127 - 15) + (23 - 10) + 1) << 23 };
298 unsigned int sign_mask = 0x80000000u;
299 __half o;
300 o.x = static_cast<unsigned short>(0x0u);
301
302 unsigned int sign = f.u & sign_mask;
303 f.u ^= sign;
304
305 // NOTE all the integer compares in this function can be safely
306 // compiled into signed compares since all operands are below
307 // 0x80000000. Important if you want fast straight SSE2 code
308 // (since there's no unsigned PCMPGTD).
309
310 if (f.u >= f16max.u) { // result is Inf or NaN (all exponent bits set)
311 o.x = (f.u > f32infty.u) ? 0x7e00 : 0x7c00; // NaN->qNaN and Inf->Inf
312 } else { // (De)normalized number or zero
313 if (f.u < (113 << 23)) { // resulting FP16 is subnormal or zero
314 // use a magic value to align our 10 mantissa bits at the bottom of
315 // the float. as long as FP addition is round-to-nearest-even this
316 // just works.
317 f.f += denorm_magic.f;
318
319 // and one integer subtract of the bias later, we have our final float!
320 o.x = static_cast<unsigned short>(f.u - denorm_magic.u);
321 } else {
322 unsigned int mant_odd = (f.u >> 13) & 1; // resulting mantissa is odd
323
324 // update exponent, rounding bias part 1
325 f.u += ((unsigned int)(15 - 127) << 23) + 0xfff;
326 // rounding bias part 2
327 f.u += mant_odd;
328 // take the bits!
329 o.x = static_cast<unsigned short>(f.u >> 13);
330 }
331 }
332
333 o.x |= static_cast<unsigned short>(sign >> 16);
334 return o;
335#endif
336}
337
338EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC float half_to_float(__half h) {
339#if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 300
340 return __half2float(h);
341
342#elif defined(EIGEN_HAS_FP16_C)
343 return _cvtsh_ss(h.x);
344
345#else
346 const FP32 magic = { 113 << 23 };
347 const unsigned int shifted_exp = 0x7c00 << 13; // exponent mask after shift
348 FP32 o;
349
350 o.u = (h.x & 0x7fff) << 13; // exponent/mantissa bits
351 unsigned int exp = shifted_exp & o.u; // just the exponent
352 o.u += (127 - 15) << 23; // exponent adjust
353
354 // handle exponent special cases
355 if (exp == shifted_exp) { // Inf/NaN?
356 o.u += (128 - 16) << 23; // extra exp adjust
357 } else if (exp == 0) { // Zero/Denormal?
358 o.u += 1 << 23; // extra exp adjust
359 o.f -= magic.f; // renormalize
360 }
361
362 o.u |= (h.x & 0x8000) << 16; // sign bit
363 return o.f;
364#endif
365}
366
367// --- standard functions ---
368
369EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool (isinf)(const half& a) {
370 return (a.x & 0x7fff) == 0x7c00;
371}
372EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool (isnan)(const half& a) {
373#if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 530
374 return __hisnan(a);
375#else
376 return (a.x & 0x7fff) > 0x7c00;
377#endif
378}
379EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bool (isfinite)(const half& a) {
380 return !(isinf EIGEN_NOT_A_MACRO (a)) && !(isnan EIGEN_NOT_A_MACRO (a));
381}
382
383EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half abs(const half& a) {
384 half result;
385 result.x = a.x & 0x7FFF;
386 return result;
387}
388EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half exp(const half& a) {
389 return half(::expf(float(a)));
390}
391EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half log(const half& a) {
392#if defined(EIGEN_HAS_CUDA_FP16) && defined __CUDACC_VER__ && __CUDACC_VER__ >= 80000 && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 530
393 return Eigen::half(::hlog(a));
394#else
395 return half(::logf(float(a)));
396#endif
397}
398EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half log1p(const half& a) {
399 return half(numext::log1p(float(a)));
400}
401EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half log10(const half& a) {
402 return half(::log10f(float(a)));
403}
404EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half sqrt(const half& a) {
405 return half(::sqrtf(float(a)));
406}
407EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half pow(const half& a, const half& b) {
408 return half(::powf(float(a), float(b)));
409}
410EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half sin(const half& a) {
411 return half(::sinf(float(a)));
412}
413EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half cos(const half& a) {
414 return half(::cosf(float(a)));
415}
416EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half tan(const half& a) {
417 return half(::tanf(float(a)));
418}
419EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half tanh(const half& a) {
420 return half(::tanhf(float(a)));
421}
422EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half floor(const half& a) {
423 return half(::floorf(float(a)));
424}
425EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half ceil(const half& a) {
426 return half(::ceilf(float(a)));
427}
428
429EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half (min)(const half& a, const half& b) {
430#if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 530
431 return __hlt(b, a) ? b : a;
432#else
433 const float f1 = static_cast<float>(a);
434 const float f2 = static_cast<float>(b);
435 return f2 < f1 ? b : a;
436#endif
437}
438EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC half (max)(const half& a, const half& b) {
439#if defined(EIGEN_HAS_CUDA_FP16) && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 530
440 return __hlt(a, b) ? b : a;
441#else
442 const float f1 = static_cast<float>(a);
443 const float f2 = static_cast<float>(b);
444 return f1 < f2 ? b : a;
445#endif
446}
447
448EIGEN_ALWAYS_INLINE std::ostream& operator << (std::ostream& os, const half& v) {
449 os << static_cast<float>(v);
450 return os;
451}
452
453} // end namespace half_impl
454
455// import Eigen::half_impl::half into Eigen namespace
456// using half_impl::half;
457
458namespace internal {
459
460template<>
461struct random_default_impl<half, false, false>
462{
463 static inline half run(const half& x, const half& y)
464 {
465 return x + (y-x) * half(float(std::rand()) / float(RAND_MAX));
466 }
467 static inline half run()
468 {
469 return run(half(-1.f), half(1.f));
470 }
471};
472
473template<> struct is_arithmetic<half> { enum { value = true }; };
474
475} // end namespace internal
476
477template<> struct NumTraits<Eigen::half>
478 : GenericNumTraits<Eigen::half>
479{
480 EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half epsilon() {
481 return half_impl::raw_uint16_to_half(0x0800);
482 }
483 EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half dummy_precision() { return Eigen::half(1e-2f); }
484 EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half highest() {
485 return half_impl::raw_uint16_to_half(0x7bff);
486 }
487 EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half lowest() {
488 return half_impl::raw_uint16_to_half(0xfbff);
489 }
490 EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half infinity() {
491 return half_impl::raw_uint16_to_half(0x7c00);
492 }
493 EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE Eigen::half quiet_NaN() {
494 return half_impl::raw_uint16_to_half(0x7c01);
495 }
496};
497
498} // end namespace Eigen
499
500// C-like standard mathematical functions and trancendentals.
501EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half fabsh(const Eigen::half& a) {
502 Eigen::half result;
503 result.x = a.x & 0x7FFF;
504 return result;
505}
506EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half exph(const Eigen::half& a) {
507 return Eigen::half(::expf(float(a)));
508}
509EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half logh(const Eigen::half& a) {
510#if defined __CUDACC_VER__ && __CUDACC_VER__ >= 80000 && defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 530
511 return Eigen::half(::hlog(a));
512#else
513 return Eigen::half(::logf(float(a)));
514#endif
515}
516EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half sqrth(const Eigen::half& a) {
517 return Eigen::half(::sqrtf(float(a)));
518}
519EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half powh(const Eigen::half& a, const Eigen::half& b) {
520 return Eigen::half(::powf(float(a), float(b)));
521}
522EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half floorh(const Eigen::half& a) {
523 return Eigen::half(::floorf(float(a)));
524}
525EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half ceilh(const Eigen::half& a) {
526 return Eigen::half(::ceilf(float(a)));
527}
528
529namespace std {
530
531#if __cplusplus > 199711L
532template <>
533struct hash<Eigen::half> {
534 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE std::size_t operator()(const Eigen::half& a) const {
535 return static_cast<std::size_t>(a.x);
536 }
537};
538#endif
539
540} // end namespace std
541
542
543// Add the missing shfl_xor intrinsic
544#if defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 300
545__device__ EIGEN_STRONG_INLINE Eigen::half __shfl_xor(Eigen::half var, int laneMask, int width=warpSize) {
546 return static_cast<Eigen::half>(__shfl_xor(static_cast<float>(var), laneMask, width));
547}
548#endif
549
550// ldg() has an overload for __half, but we also need one for Eigen::half.
551#if defined(__CUDA_ARCH__) && __CUDA_ARCH__ >= 350
552EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Eigen::half __ldg(const Eigen::half* ptr) {
553 return Eigen::half_impl::raw_uint16_to_half(
554 __ldg(reinterpret_cast<const unsigned short*>(ptr)));
555}
556#endif
557
558
559#if defined(__CUDA_ARCH__)
560namespace Eigen {
561namespace numext {
562
563template<>
564EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
565bool (isnan)(const Eigen::half& h) {
566 return (half_impl::isnan)(h);
567}
568
569template<>
570EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
571bool (isinf)(const Eigen::half& h) {
572 return (half_impl::isinf)(h);
573}
574
575template<>
576EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE
577bool (isfinite)(const Eigen::half& h) {
578 return (half_impl::isfinite)(h);
579}
580
581} // namespace Eigen
582} // namespace numext
583#endif
584
585#endif // EIGEN_HALF_CUDA_H
Namespace containing all symbols from the Eigen library.
Definition: Core:287
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_sign_op< typename Derived::Scalar >, const Derived > sign(const Eigen::ArrayBase< Derived > &x)
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:33
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_isfinite_op< typename Derived::Scalar >, const Derived > isfinite(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_isnan_op< typename Derived::Scalar >, const Derived > isnan(const Eigen::ArrayBase< Derived > &x)