1 //============================================================================
2 // Copyright (c) Kitware, Inc.
3 // All rights reserved.
4 // See LICENSE.txt for details.
5 //
6 // This software is distributed WITHOUT ANY WARRANTY; without even
7 // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
8 // PURPOSE. See the above copyright notice for more information.
9 //============================================================================
10 // **** DO NOT EDIT THIS FILE!!! ****
11 // This file is automatically generated by Math.h.in
12
13 #ifndef vtk_m_Math_h
14 #define vtk_m_Math_h
15
16 #include <vtkm/TypeTraits.h>
17 #include <vtkm/Types.h>
18 #include <vtkm/VecTraits.h>
19
20 #include <limits> // must be found with or without CUDA.
21 #ifndef VTKM_CUDA
22 #include <cmath>
23 #include <cstring>
24 #include <limits.h>
25 #include <math.h>
26 #include <stdlib.h>
27 #endif // !VTKM_CUDA
28
29 #if !defined(VTKM_CUDA_DEVICE_PASS)
30 #define VTKM_USE_STL
31 #include <algorithm>
32 #endif
33
34 #ifdef VTKM_MSVC
35 #include <intrin.h> // For bitwise intrinsics (__popcnt, etc)
36 #include <vtkm/internal/Windows.h> // for types used by MSVC intrinsics.
37 #ifndef VTKM_CUDA
38 #include <math.h>
39 #endif // VTKM_CUDA
40 #endif // VTKM_MSVC
41
42 #define VTKM_CUDA_MATH_FUNCTION_32(func) func##f
43 #define VTKM_CUDA_MATH_FUNCTION_64(func) func
44
45 // clang-format off
46 namespace vtkm
47 {
48
49 //-----------------------------------------------------------------------------
50 namespace detail
51 {
52 template <typename T>
53 struct FloatingPointReturnType
54 {
55 using ctype = typename vtkm::VecTraits<T>::ComponentType;
56 using representable_as_float_type = std::integral_constant<bool,
57 ((sizeof(ctype) < sizeof(float)) || std::is_same<ctype, vtkm::Float32>::value)>;
58 using Type = typename std::conditional<representable_as_float_type::value,
59 vtkm::Float32,
60 vtkm::Float64>::type;
61 };
62 } // namespace detail
63
64 /// Returns the constant 2 times Pi in float32.
65 ///
TwoPif()66 static constexpr inline VTKM_EXEC_CONT vtkm::Float32 TwoPif()
67 {
68 return 6.28318530717958647692528676655900576f;
69 }
70
71 /// Returns the constant Pi in float32.
72 ///
Pif()73 static constexpr inline VTKM_EXEC_CONT vtkm::Float32 Pif()
74 {
75 return 3.14159265358979323846264338327950288f;
76 }
77
78 /// Returns the constant Pi halves in float32.
79 ///
Pi_2f()80 static constexpr inline VTKM_EXEC_CONT vtkm::Float32 Pi_2f()
81 {
82 return 1.57079632679489661923132169163975144f;
83 }
84
85 /// Returns the constant Pi thirds in float32.
86 ///
Pi_3f()87 static constexpr inline VTKM_EXEC_CONT vtkm::Float32 Pi_3f()
88 {
89 return 1.04719755119659774615421446109316762f;
90 }
91
92 /// Returns the constant Pi fourths in float32.
93 ///
Pi_4f()94 static constexpr inline VTKM_EXEC_CONT vtkm::Float32 Pi_4f()
95 {
96 return 0.78539816339744830961566084581987572f;
97 }
98
99 /// Returns the constant Pi one hundred and eightieth in float32.
100 ///
Pi_180f()101 static constexpr inline VTKM_EXEC_CONT vtkm::Float32 Pi_180f()
102 {
103 return 0.01745329251994329547437168059786927f;
104 }
105
106 /// Returns the constant 2 times Pi in float64.
107 ///
TwoPi()108 static constexpr inline VTKM_EXEC_CONT vtkm::Float64 TwoPi()
109 {
110 return 6.28318530717958647692528676655900576;
111 }
112
113 /// Returns the constant Pi in float64.
114 ///
Pi()115 static constexpr inline VTKM_EXEC_CONT vtkm::Float64 Pi()
116 {
117 return 3.14159265358979323846264338327950288;
118 }
119
120 /// Returns the constant Pi halves in float64.
121 ///
Pi_2()122 static constexpr inline VTKM_EXEC_CONT vtkm::Float64 Pi_2()
123 {
124 return 1.57079632679489661923132169163975144;
125 }
126
127 /// Returns the constant Pi thirds in float64.
128 ///
Pi_3()129 static constexpr inline VTKM_EXEC_CONT vtkm::Float64 Pi_3()
130 {
131 return 1.04719755119659774615421446109316762;
132 }
133
134 /// Returns the constant Pi fourths in float64.
135 ///
Pi_4()136 static constexpr inline VTKM_EXEC_CONT vtkm::Float64 Pi_4()
137 {
138 return 0.78539816339744830961566084581987572;
139 }
140
141 /// Returns the constant Pi one hundred and eightieth in float64.
142 ///
Pi_180()143 static constexpr inline VTKM_EXEC_CONT vtkm::Float64 Pi_180()
144 {
145 return 0.01745329251994329547437168059786927;
146 }
147
148 /// Returns the constant 2 times Pi.
149 ///
150 template <typename T>
151 static constexpr inline VTKM_EXEC_CONT auto TwoPi() -> typename detail::FloatingPointReturnType<T>::Type
152 {
153 using FT = typename detail::FloatingPointReturnType<T>::Type;
154 using RAFT = typename detail::FloatingPointReturnType<T>::representable_as_float_type;
155 return static_cast<FT>(RAFT::value ? TwoPif() : TwoPi());
156 }
157
158 /// Returns the constant Pi.
159 ///
160 template <typename T>
161 static constexpr inline VTKM_EXEC_CONT auto Pi() -> typename detail::FloatingPointReturnType<T>::Type
162 {
163 using FT = typename detail::FloatingPointReturnType<T>::Type;
164 using RAFT = typename detail::FloatingPointReturnType<T>::representable_as_float_type;
165 return static_cast<FT>(RAFT::value ? Pif() : Pi());
166 }
167
168 /// Returns the constant Pi halves.
169 ///
170 template <typename T>
171 static constexpr inline VTKM_EXEC_CONT auto Pi_2() -> typename detail::FloatingPointReturnType<T>::Type
172 {
173 using FT = typename detail::FloatingPointReturnType<T>::Type;
174 using RAFT = typename detail::FloatingPointReturnType<T>::representable_as_float_type;
175 return static_cast<FT>(RAFT::value ? Pi_2f() : Pi_2());
176 }
177
178 /// Returns the constant Pi thirds.
179 ///
180 template <typename T>
181 static constexpr inline VTKM_EXEC_CONT auto Pi_3() -> typename detail::FloatingPointReturnType<T>::Type
182 {
183 using FT = typename detail::FloatingPointReturnType<T>::Type;
184 using RAFT = typename detail::FloatingPointReturnType<T>::representable_as_float_type;
185 return static_cast<FT>(RAFT::value ? Pi_3f() : Pi_3());
186 }
187
188 /// Returns the constant Pi fourths.
189 ///
190 template <typename T>
191 static constexpr inline VTKM_EXEC_CONT auto Pi_4() -> typename detail::FloatingPointReturnType<T>::Type
192 {
193 using FT = typename detail::FloatingPointReturnType<T>::Type;
194 using RAFT = typename detail::FloatingPointReturnType<T>::representable_as_float_type;
195 return static_cast<FT>(RAFT::value ? Pi_4f() : Pi_4());
196 }
197 /// Returns the constant Pi one hundred and eightieth.
198 ///
199 template <typename T>
200 static constexpr inline VTKM_EXEC_CONT auto Pi_180() -> typename detail::FloatingPointReturnType<T>::Type
201 {
202 using FT = typename detail::FloatingPointReturnType<T>::Type;
203 using RAFT = typename detail::FloatingPointReturnType<T>::representable_as_float_type;
204 return static_cast<FT>(RAFT::value ? Pi_180f() : Pi_180());
205 }
206
207 /// Compute the sine of \p x.
208 ///
209
Sin(vtkm::Float32 x)210 inline VTKM_EXEC_CONT vtkm::Float32 Sin(vtkm::Float32 x)
211 {
212 #ifdef VTKM_CUDA
213 return VTKM_CUDA_MATH_FUNCTION_32(sin)(x);
214 #else
215 return std::sin(x);
216 #endif
217 }
218
Sin(vtkm::Float64 x)219 inline VTKM_EXEC_CONT vtkm::Float64 Sin(vtkm::Float64 x)
220 {
221 #ifdef VTKM_CUDA
222 return VTKM_CUDA_MATH_FUNCTION_64(sin)(x);
223 #else
224 return std::sin(x);
225 #endif
226 }
227 template <typename T>
Sin(const T & x)228 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Sin(const T& x)
229 {
230 using RT = typename detail::FloatingPointReturnType<T>::Type;
231 return vtkm::Sin(static_cast<RT>(x));
232 }
233 template <typename T, vtkm::IdComponent N>
Sin(const vtkm::Vec<T,N> & x)234 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Sin(
235 const vtkm::Vec<T, N>& x)
236 {
237 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
238 for (vtkm::IdComponent index = 0; index < N; index++)
239 {
240 result[index] = vtkm::Sin(x[index]);
241 }
242 return result;
243 }
244 template <typename T>
Sin(const vtkm::Vec<T,4> & x)245 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Sin(
246 const vtkm::Vec<T, 4>& x)
247 {
248 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
249 vtkm::Sin(x[0]), vtkm::Sin(x[1]), vtkm::Sin(x[2]), vtkm::Sin(x[3]));
250 }
251 template <typename T>
Sin(const vtkm::Vec<T,3> & x)252 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Sin(
253 const vtkm::Vec<T, 3>& x)
254 {
255 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
256 vtkm::Sin(x[0]), vtkm::Sin(x[1]), vtkm::Sin(x[2]));
257 }
258 template <typename T>
Sin(const vtkm::Vec<T,2> & x)259 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Sin(
260 const vtkm::Vec<T, 2>& x)
261 {
262 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Sin(x[0]),
263 vtkm::Sin(x[1]));
264 }
265
266 /// Compute the cosine of \p x.
267 ///
268
Cos(vtkm::Float32 x)269 inline VTKM_EXEC_CONT vtkm::Float32 Cos(vtkm::Float32 x)
270 {
271 #ifdef VTKM_CUDA
272 return VTKM_CUDA_MATH_FUNCTION_32(cos)(x);
273 #else
274 return std::cos(x);
275 #endif
276 }
277
Cos(vtkm::Float64 x)278 inline VTKM_EXEC_CONT vtkm::Float64 Cos(vtkm::Float64 x)
279 {
280 #ifdef VTKM_CUDA
281 return VTKM_CUDA_MATH_FUNCTION_64(cos)(x);
282 #else
283 return std::cos(x);
284 #endif
285 }
286 template <typename T>
Cos(const T & x)287 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Cos(const T& x)
288 {
289 using RT = typename detail::FloatingPointReturnType<T>::Type;
290 return vtkm::Cos(static_cast<RT>(x));
291 }
292 template <typename T, vtkm::IdComponent N>
Cos(const vtkm::Vec<T,N> & x)293 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Cos(
294 const vtkm::Vec<T, N>& x)
295 {
296 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
297 for (vtkm::IdComponent index = 0; index < N; index++)
298 {
299 result[index] = vtkm::Cos(x[index]);
300 }
301 return result;
302 }
303 template <typename T>
Cos(const vtkm::Vec<T,4> & x)304 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Cos(
305 const vtkm::Vec<T, 4>& x)
306 {
307 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
308 vtkm::Cos(x[0]), vtkm::Cos(x[1]), vtkm::Cos(x[2]), vtkm::Cos(x[3]));
309 }
310 template <typename T>
Cos(const vtkm::Vec<T,3> & x)311 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Cos(
312 const vtkm::Vec<T, 3>& x)
313 {
314 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
315 vtkm::Cos(x[0]), vtkm::Cos(x[1]), vtkm::Cos(x[2]));
316 }
317 template <typename T>
Cos(const vtkm::Vec<T,2> & x)318 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Cos(
319 const vtkm::Vec<T, 2>& x)
320 {
321 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Cos(x[0]),
322 vtkm::Cos(x[1]));
323 }
324
325 /// Compute the tangent of \p x.
326 ///
327
Tan(vtkm::Float32 x)328 inline VTKM_EXEC_CONT vtkm::Float32 Tan(vtkm::Float32 x)
329 {
330 #ifdef VTKM_CUDA
331 return VTKM_CUDA_MATH_FUNCTION_32(tan)(x);
332 #else
333 return std::tan(x);
334 #endif
335 }
336
Tan(vtkm::Float64 x)337 inline VTKM_EXEC_CONT vtkm::Float64 Tan(vtkm::Float64 x)
338 {
339 #ifdef VTKM_CUDA
340 return VTKM_CUDA_MATH_FUNCTION_64(tan)(x);
341 #else
342 return std::tan(x);
343 #endif
344 }
345 template <typename T>
Tan(const T & x)346 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Tan(const T& x)
347 {
348 using RT = typename detail::FloatingPointReturnType<T>::Type;
349 return vtkm::Tan(static_cast<RT>(x));
350 }
351 template <typename T, vtkm::IdComponent N>
Tan(const vtkm::Vec<T,N> & x)352 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Tan(
353 const vtkm::Vec<T, N>& x)
354 {
355 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
356 for (vtkm::IdComponent index = 0; index < N; index++)
357 {
358 result[index] = vtkm::Tan(x[index]);
359 }
360 return result;
361 }
362 template <typename T>
Tan(const vtkm::Vec<T,4> & x)363 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Tan(
364 const vtkm::Vec<T, 4>& x)
365 {
366 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
367 vtkm::Tan(x[0]), vtkm::Tan(x[1]), vtkm::Tan(x[2]), vtkm::Tan(x[3]));
368 }
369 template <typename T>
Tan(const vtkm::Vec<T,3> & x)370 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Tan(
371 const vtkm::Vec<T, 3>& x)
372 {
373 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
374 vtkm::Tan(x[0]), vtkm::Tan(x[1]), vtkm::Tan(x[2]));
375 }
376 template <typename T>
Tan(const vtkm::Vec<T,2> & x)377 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Tan(
378 const vtkm::Vec<T, 2>& x)
379 {
380 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Tan(x[0]),
381 vtkm::Tan(x[1]));
382 }
383
384 /// Compute the arc sine of \p x.
385 ///
386
ASin(vtkm::Float32 x)387 inline VTKM_EXEC_CONT vtkm::Float32 ASin(vtkm::Float32 x)
388 {
389 #ifdef VTKM_CUDA
390 return VTKM_CUDA_MATH_FUNCTION_32(asin)(x);
391 #else
392 return std::asin(x);
393 #endif
394 }
395
ASin(vtkm::Float64 x)396 inline VTKM_EXEC_CONT vtkm::Float64 ASin(vtkm::Float64 x)
397 {
398 #ifdef VTKM_CUDA
399 return VTKM_CUDA_MATH_FUNCTION_64(asin)(x);
400 #else
401 return std::asin(x);
402 #endif
403 }
404 template <typename T>
ASin(const T & x)405 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type ASin(const T& x)
406 {
407 using RT = typename detail::FloatingPointReturnType<T>::Type;
408 return vtkm::ASin(static_cast<RT>(x));
409 }
410 template <typename T, vtkm::IdComponent N>
ASin(const vtkm::Vec<T,N> & x)411 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> ASin(
412 const vtkm::Vec<T, N>& x)
413 {
414 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
415 for (vtkm::IdComponent index = 0; index < N; index++)
416 {
417 result[index] = vtkm::ASin(x[index]);
418 }
419 return result;
420 }
421 template <typename T>
ASin(const vtkm::Vec<T,4> & x)422 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> ASin(
423 const vtkm::Vec<T, 4>& x)
424 {
425 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
426 vtkm::ASin(x[0]), vtkm::ASin(x[1]), vtkm::ASin(x[2]), vtkm::ASin(x[3]));
427 }
428 template <typename T>
ASin(const vtkm::Vec<T,3> & x)429 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> ASin(
430 const vtkm::Vec<T, 3>& x)
431 {
432 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
433 vtkm::ASin(x[0]), vtkm::ASin(x[1]), vtkm::ASin(x[2]));
434 }
435 template <typename T>
ASin(const vtkm::Vec<T,2> & x)436 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> ASin(
437 const vtkm::Vec<T, 2>& x)
438 {
439 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::ASin(x[0]),
440 vtkm::ASin(x[1]));
441 }
442
443 /// Compute the arc cosine of \p x.
444 ///
445
ACos(vtkm::Float32 x)446 inline VTKM_EXEC_CONT vtkm::Float32 ACos(vtkm::Float32 x)
447 {
448 #ifdef VTKM_CUDA
449 return VTKM_CUDA_MATH_FUNCTION_32(acos)(x);
450 #else
451 return std::acos(x);
452 #endif
453 }
454
ACos(vtkm::Float64 x)455 inline VTKM_EXEC_CONT vtkm::Float64 ACos(vtkm::Float64 x)
456 {
457 #ifdef VTKM_CUDA
458 return VTKM_CUDA_MATH_FUNCTION_64(acos)(x);
459 #else
460 return std::acos(x);
461 #endif
462 }
463 template <typename T>
ACos(const T & x)464 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type ACos(const T& x)
465 {
466 using RT = typename detail::FloatingPointReturnType<T>::Type;
467 return vtkm::ACos(static_cast<RT>(x));
468 }
469 template <typename T, vtkm::IdComponent N>
ACos(const vtkm::Vec<T,N> & x)470 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> ACos(
471 const vtkm::Vec<T, N>& x)
472 {
473 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
474 for (vtkm::IdComponent index = 0; index < N; index++)
475 {
476 result[index] = vtkm::ACos(x[index]);
477 }
478 return result;
479 }
480 template <typename T>
ACos(const vtkm::Vec<T,4> & x)481 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> ACos(
482 const vtkm::Vec<T, 4>& x)
483 {
484 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
485 vtkm::ACos(x[0]), vtkm::ACos(x[1]), vtkm::ACos(x[2]), vtkm::ACos(x[3]));
486 }
487 template <typename T>
ACos(const vtkm::Vec<T,3> & x)488 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> ACos(
489 const vtkm::Vec<T, 3>& x)
490 {
491 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
492 vtkm::ACos(x[0]), vtkm::ACos(x[1]), vtkm::ACos(x[2]));
493 }
494 template <typename T>
ACos(const vtkm::Vec<T,2> & x)495 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> ACos(
496 const vtkm::Vec<T, 2>& x)
497 {
498 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::ACos(x[0]),
499 vtkm::ACos(x[1]));
500 }
501
502 /// Compute the arc tangent of \p x.
503 ///
504
ATan(vtkm::Float32 x)505 inline VTKM_EXEC_CONT vtkm::Float32 ATan(vtkm::Float32 x)
506 {
507 #ifdef VTKM_CUDA
508 return VTKM_CUDA_MATH_FUNCTION_32(atan)(x);
509 #else
510 return std::atan(x);
511 #endif
512 }
513
ATan(vtkm::Float64 x)514 inline VTKM_EXEC_CONT vtkm::Float64 ATan(vtkm::Float64 x)
515 {
516 #ifdef VTKM_CUDA
517 return VTKM_CUDA_MATH_FUNCTION_64(atan)(x);
518 #else
519 return std::atan(x);
520 #endif
521 }
522 template <typename T>
ATan(const T & x)523 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type ATan(const T& x)
524 {
525 using RT = typename detail::FloatingPointReturnType<T>::Type;
526 return vtkm::ATan(static_cast<RT>(x));
527 }
528 template <typename T, vtkm::IdComponent N>
ATan(const vtkm::Vec<T,N> & x)529 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> ATan(
530 const vtkm::Vec<T, N>& x)
531 {
532 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
533 for (vtkm::IdComponent index = 0; index < N; index++)
534 {
535 result[index] = vtkm::ATan(x[index]);
536 }
537 return result;
538 }
539 template <typename T>
ATan(const vtkm::Vec<T,4> & x)540 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> ATan(
541 const vtkm::Vec<T, 4>& x)
542 {
543 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
544 vtkm::ATan(x[0]), vtkm::ATan(x[1]), vtkm::ATan(x[2]), vtkm::ATan(x[3]));
545 }
546 template <typename T>
ATan(const vtkm::Vec<T,3> & x)547 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> ATan(
548 const vtkm::Vec<T, 3>& x)
549 {
550 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
551 vtkm::ATan(x[0]), vtkm::ATan(x[1]), vtkm::ATan(x[2]));
552 }
553 template <typename T>
ATan(const vtkm::Vec<T,2> & x)554 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> ATan(
555 const vtkm::Vec<T, 2>& x)
556 {
557 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::ATan(x[0]),
558 vtkm::ATan(x[1]));
559 }
560
561 /// Compute the arc tangent of \p x / \p y using the signs of both arguments
562 /// to determine the quadrant of the return value.
563 ///
ATan2(vtkm::Float32 x,vtkm::Float32 y)564 static inline VTKM_EXEC_CONT vtkm::Float32 ATan2(vtkm::Float32 x, vtkm::Float32 y)
565 {
566 #ifdef VTKM_CUDA
567 return VTKM_CUDA_MATH_FUNCTION_32(atan2)(x, y);
568 #else
569 return std::atan2(x, y);
570 #endif
571 }
ATan2(vtkm::Float64 x,vtkm::Float64 y)572 static inline VTKM_EXEC_CONT vtkm::Float64 ATan2(vtkm::Float64 x, vtkm::Float64 y)
573 {
574 #ifdef VTKM_CUDA
575 return VTKM_CUDA_MATH_FUNCTION_64(atan2)(x, y);
576 #else
577 return std::atan2(x, y);
578 #endif
579 }
580
581 /// Compute the hyperbolic sine of \p x.
582 ///
583
SinH(vtkm::Float32 x)584 inline VTKM_EXEC_CONT vtkm::Float32 SinH(vtkm::Float32 x)
585 {
586 #ifdef VTKM_CUDA
587 return VTKM_CUDA_MATH_FUNCTION_32(sinh)(x);
588 #else
589 return std::sinh(x);
590 #endif
591 }
592
SinH(vtkm::Float64 x)593 inline VTKM_EXEC_CONT vtkm::Float64 SinH(vtkm::Float64 x)
594 {
595 #ifdef VTKM_CUDA
596 return VTKM_CUDA_MATH_FUNCTION_64(sinh)(x);
597 #else
598 return std::sinh(x);
599 #endif
600 }
601 template <typename T>
SinH(const T & x)602 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type SinH(const T& x)
603 {
604 using RT = typename detail::FloatingPointReturnType<T>::Type;
605 return vtkm::SinH(static_cast<RT>(x));
606 }
607 template <typename T, vtkm::IdComponent N>
SinH(const vtkm::Vec<T,N> & x)608 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> SinH(
609 const vtkm::Vec<T, N>& x)
610 {
611 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
612 for (vtkm::IdComponent index = 0; index < N; index++)
613 {
614 result[index] = vtkm::SinH(x[index]);
615 }
616 return result;
617 }
618 template <typename T>
SinH(const vtkm::Vec<T,4> & x)619 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> SinH(
620 const vtkm::Vec<T, 4>& x)
621 {
622 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
623 vtkm::SinH(x[0]), vtkm::SinH(x[1]), vtkm::SinH(x[2]), vtkm::SinH(x[3]));
624 }
625 template <typename T>
SinH(const vtkm::Vec<T,3> & x)626 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> SinH(
627 const vtkm::Vec<T, 3>& x)
628 {
629 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
630 vtkm::SinH(x[0]), vtkm::SinH(x[1]), vtkm::SinH(x[2]));
631 }
632 template <typename T>
SinH(const vtkm::Vec<T,2> & x)633 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> SinH(
634 const vtkm::Vec<T, 2>& x)
635 {
636 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::SinH(x[0]),
637 vtkm::SinH(x[1]));
638 }
639
640 /// Compute the hyperbolic cosine of \p x.
641 ///
642
CosH(vtkm::Float32 x)643 inline VTKM_EXEC_CONT vtkm::Float32 CosH(vtkm::Float32 x)
644 {
645 #ifdef VTKM_CUDA
646 return VTKM_CUDA_MATH_FUNCTION_32(cosh)(x);
647 #else
648 return std::cosh(x);
649 #endif
650 }
651
CosH(vtkm::Float64 x)652 inline VTKM_EXEC_CONT vtkm::Float64 CosH(vtkm::Float64 x)
653 {
654 #ifdef VTKM_CUDA
655 return VTKM_CUDA_MATH_FUNCTION_64(cosh)(x);
656 #else
657 return std::cosh(x);
658 #endif
659 }
660 template <typename T>
CosH(const T & x)661 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type CosH(const T& x)
662 {
663 using RT = typename detail::FloatingPointReturnType<T>::Type;
664 return vtkm::CosH(static_cast<RT>(x));
665 }
666 template <typename T, vtkm::IdComponent N>
CosH(const vtkm::Vec<T,N> & x)667 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> CosH(
668 const vtkm::Vec<T, N>& x)
669 {
670 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
671 for (vtkm::IdComponent index = 0; index < N; index++)
672 {
673 result[index] = vtkm::CosH(x[index]);
674 }
675 return result;
676 }
677 template <typename T>
CosH(const vtkm::Vec<T,4> & x)678 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> CosH(
679 const vtkm::Vec<T, 4>& x)
680 {
681 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
682 vtkm::CosH(x[0]), vtkm::CosH(x[1]), vtkm::CosH(x[2]), vtkm::CosH(x[3]));
683 }
684 template <typename T>
CosH(const vtkm::Vec<T,3> & x)685 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> CosH(
686 const vtkm::Vec<T, 3>& x)
687 {
688 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
689 vtkm::CosH(x[0]), vtkm::CosH(x[1]), vtkm::CosH(x[2]));
690 }
691 template <typename T>
CosH(const vtkm::Vec<T,2> & x)692 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> CosH(
693 const vtkm::Vec<T, 2>& x)
694 {
695 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::CosH(x[0]),
696 vtkm::CosH(x[1]));
697 }
698
699 /// Compute the hyperbolic tangent of \p x.
700 ///
701
TanH(vtkm::Float32 x)702 inline VTKM_EXEC_CONT vtkm::Float32 TanH(vtkm::Float32 x)
703 {
704 #ifdef VTKM_CUDA
705 return VTKM_CUDA_MATH_FUNCTION_32(tanh)(x);
706 #else
707 return std::tanh(x);
708 #endif
709 }
710
TanH(vtkm::Float64 x)711 inline VTKM_EXEC_CONT vtkm::Float64 TanH(vtkm::Float64 x)
712 {
713 #ifdef VTKM_CUDA
714 return VTKM_CUDA_MATH_FUNCTION_64(tanh)(x);
715 #else
716 return std::tanh(x);
717 #endif
718 }
719 template <typename T>
TanH(const T & x)720 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type TanH(const T& x)
721 {
722 using RT = typename detail::FloatingPointReturnType<T>::Type;
723 return vtkm::TanH(static_cast<RT>(x));
724 }
725 template <typename T, vtkm::IdComponent N>
TanH(const vtkm::Vec<T,N> & x)726 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> TanH(
727 const vtkm::Vec<T, N>& x)
728 {
729 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
730 for (vtkm::IdComponent index = 0; index < N; index++)
731 {
732 result[index] = vtkm::TanH(x[index]);
733 }
734 return result;
735 }
736 template <typename T>
TanH(const vtkm::Vec<T,4> & x)737 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> TanH(
738 const vtkm::Vec<T, 4>& x)
739 {
740 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
741 vtkm::TanH(x[0]), vtkm::TanH(x[1]), vtkm::TanH(x[2]), vtkm::TanH(x[3]));
742 }
743 template <typename T>
TanH(const vtkm::Vec<T,3> & x)744 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> TanH(
745 const vtkm::Vec<T, 3>& x)
746 {
747 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
748 vtkm::TanH(x[0]), vtkm::TanH(x[1]), vtkm::TanH(x[2]));
749 }
750 template <typename T>
TanH(const vtkm::Vec<T,2> & x)751 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> TanH(
752 const vtkm::Vec<T, 2>& x)
753 {
754 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::TanH(x[0]),
755 vtkm::TanH(x[1]));
756 }
757
758 /// Compute the hyperbolic arc sine of \p x.
759 ///
760
ASinH(vtkm::Float32 x)761 inline VTKM_EXEC_CONT vtkm::Float32 ASinH(vtkm::Float32 x)
762 {
763 #ifdef VTKM_CUDA
764 return VTKM_CUDA_MATH_FUNCTION_32(asinh)(x);
765 #else
766 return std::asinh(x);
767 #endif
768 }
769
ASinH(vtkm::Float64 x)770 inline VTKM_EXEC_CONT vtkm::Float64 ASinH(vtkm::Float64 x)
771 {
772 #ifdef VTKM_CUDA
773 return VTKM_CUDA_MATH_FUNCTION_64(asinh)(x);
774 #else
775 return std::asinh(x);
776 #endif
777 }
778 template <typename T>
ASinH(const T & x)779 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type ASinH(const T& x)
780 {
781 using RT = typename detail::FloatingPointReturnType<T>::Type;
782 return vtkm::ASinH(static_cast<RT>(x));
783 }
784 template <typename T, vtkm::IdComponent N>
ASinH(const vtkm::Vec<T,N> & x)785 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> ASinH(
786 const vtkm::Vec<T, N>& x)
787 {
788 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
789 for (vtkm::IdComponent index = 0; index < N; index++)
790 {
791 result[index] = vtkm::ASinH(x[index]);
792 }
793 return result;
794 }
795 template <typename T>
ASinH(const vtkm::Vec<T,4> & x)796 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> ASinH(
797 const vtkm::Vec<T, 4>& x)
798 {
799 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
800 vtkm::ASinH(x[0]), vtkm::ASinH(x[1]), vtkm::ASinH(x[2]), vtkm::ASinH(x[3]));
801 }
802 template <typename T>
ASinH(const vtkm::Vec<T,3> & x)803 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> ASinH(
804 const vtkm::Vec<T, 3>& x)
805 {
806 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
807 vtkm::ASinH(x[0]), vtkm::ASinH(x[1]), vtkm::ASinH(x[2]));
808 }
809 template <typename T>
ASinH(const vtkm::Vec<T,2> & x)810 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> ASinH(
811 const vtkm::Vec<T, 2>& x)
812 {
813 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::ASinH(x[0]),
814 vtkm::ASinH(x[1]));
815 }
816
817 /// Compute the hyperbolic arc cosine of \p x.
818 ///
819
ACosH(vtkm::Float32 x)820 inline VTKM_EXEC_CONT vtkm::Float32 ACosH(vtkm::Float32 x)
821 {
822 #ifdef VTKM_CUDA
823 return VTKM_CUDA_MATH_FUNCTION_32(acosh)(x);
824 #else
825 return std::acosh(x);
826 #endif
827 }
828
ACosH(vtkm::Float64 x)829 inline VTKM_EXEC_CONT vtkm::Float64 ACosH(vtkm::Float64 x)
830 {
831 #ifdef VTKM_CUDA
832 return VTKM_CUDA_MATH_FUNCTION_64(acosh)(x);
833 #else
834 return std::acosh(x);
835 #endif
836 }
837 template <typename T>
ACosH(const T & x)838 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type ACosH(const T& x)
839 {
840 using RT = typename detail::FloatingPointReturnType<T>::Type;
841 return vtkm::ACosH(static_cast<RT>(x));
842 }
843 template <typename T, vtkm::IdComponent N>
ACosH(const vtkm::Vec<T,N> & x)844 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> ACosH(
845 const vtkm::Vec<T, N>& x)
846 {
847 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
848 for (vtkm::IdComponent index = 0; index < N; index++)
849 {
850 result[index] = vtkm::ACosH(x[index]);
851 }
852 return result;
853 }
854 template <typename T>
ACosH(const vtkm::Vec<T,4> & x)855 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> ACosH(
856 const vtkm::Vec<T, 4>& x)
857 {
858 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
859 vtkm::ACosH(x[0]), vtkm::ACosH(x[1]), vtkm::ACosH(x[2]), vtkm::ACosH(x[3]));
860 }
861 template <typename T>
ACosH(const vtkm::Vec<T,3> & x)862 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> ACosH(
863 const vtkm::Vec<T, 3>& x)
864 {
865 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
866 vtkm::ACosH(x[0]), vtkm::ACosH(x[1]), vtkm::ACosH(x[2]));
867 }
868 template <typename T>
ACosH(const vtkm::Vec<T,2> & x)869 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> ACosH(
870 const vtkm::Vec<T, 2>& x)
871 {
872 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::ACosH(x[0]),
873 vtkm::ACosH(x[1]));
874 }
875
876 /// Compute the hyperbolic arc tangent of \p x.
877 ///
878
ATanH(vtkm::Float32 x)879 inline VTKM_EXEC_CONT vtkm::Float32 ATanH(vtkm::Float32 x)
880 {
881 #ifdef VTKM_CUDA
882 return VTKM_CUDA_MATH_FUNCTION_32(atanh)(x);
883 #else
884 return std::atanh(x);
885 #endif
886 }
887
ATanH(vtkm::Float64 x)888 inline VTKM_EXEC_CONT vtkm::Float64 ATanH(vtkm::Float64 x)
889 {
890 #ifdef VTKM_CUDA
891 return VTKM_CUDA_MATH_FUNCTION_64(atanh)(x);
892 #else
893 return std::atanh(x);
894 #endif
895 }
896 template <typename T>
ATanH(const T & x)897 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type ATanH(const T& x)
898 {
899 using RT = typename detail::FloatingPointReturnType<T>::Type;
900 return vtkm::ATanH(static_cast<RT>(x));
901 }
902 template <typename T, vtkm::IdComponent N>
ATanH(const vtkm::Vec<T,N> & x)903 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> ATanH(
904 const vtkm::Vec<T, N>& x)
905 {
906 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
907 for (vtkm::IdComponent index = 0; index < N; index++)
908 {
909 result[index] = vtkm::ATanH(x[index]);
910 }
911 return result;
912 }
913 template <typename T>
ATanH(const vtkm::Vec<T,4> & x)914 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> ATanH(
915 const vtkm::Vec<T, 4>& x)
916 {
917 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
918 vtkm::ATanH(x[0]), vtkm::ATanH(x[1]), vtkm::ATanH(x[2]), vtkm::ATanH(x[3]));
919 }
920 template <typename T>
ATanH(const vtkm::Vec<T,3> & x)921 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> ATanH(
922 const vtkm::Vec<T, 3>& x)
923 {
924 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
925 vtkm::ATanH(x[0]), vtkm::ATanH(x[1]), vtkm::ATanH(x[2]));
926 }
927 template <typename T>
ATanH(const vtkm::Vec<T,2> & x)928 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> ATanH(
929 const vtkm::Vec<T, 2>& x)
930 {
931 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::ATanH(x[0]),
932 vtkm::ATanH(x[1]));
933 }
934
935 //-----------------------------------------------------------------------------
936 /// Computes \p x raised to the power of \p y.
937 ///
Pow(vtkm::Float32 x,vtkm::Float32 y)938 static inline VTKM_EXEC_CONT vtkm::Float32 Pow(vtkm::Float32 x, vtkm::Float32 y)
939 {
940 #ifdef VTKM_CUDA
941 return VTKM_CUDA_MATH_FUNCTION_32(pow)(x, y);
942 #else
943 return std::pow(x, y);
944 #endif
945 }
Pow(vtkm::Float64 x,vtkm::Float64 y)946 static inline VTKM_EXEC_CONT vtkm::Float64 Pow(vtkm::Float64 x, vtkm::Float64 y)
947 {
948 #ifdef VTKM_CUDA
949 return VTKM_CUDA_MATH_FUNCTION_64(pow)(x, y);
950 #else
951 return std::pow(x, y);
952 #endif
953 }
954
955 /// Compute the square root of \p x.
956 ///
957
Sqrt(vtkm::Float32 x)958 inline VTKM_EXEC_CONT vtkm::Float32 Sqrt(vtkm::Float32 x)
959 {
960 #ifdef VTKM_CUDA
961 return VTKM_CUDA_MATH_FUNCTION_32(sqrt)(x);
962 #else
963 return std::sqrt(x);
964 #endif
965 }
966
Sqrt(vtkm::Float64 x)967 inline VTKM_EXEC_CONT vtkm::Float64 Sqrt(vtkm::Float64 x)
968 {
969 #ifdef VTKM_CUDA
970 return VTKM_CUDA_MATH_FUNCTION_64(sqrt)(x);
971 #else
972 return std::sqrt(x);
973 #endif
974 }
975 template <typename T>
Sqrt(const T & x)976 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Sqrt(const T& x)
977 {
978 using RT = typename detail::FloatingPointReturnType<T>::Type;
979 return vtkm::Sqrt(static_cast<RT>(x));
980 }
981 template <typename T, vtkm::IdComponent N>
Sqrt(const vtkm::Vec<T,N> & x)982 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Sqrt(
983 const vtkm::Vec<T, N>& x)
984 {
985 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
986 for (vtkm::IdComponent index = 0; index < N; index++)
987 {
988 result[index] = vtkm::Sqrt(x[index]);
989 }
990 return result;
991 }
992 template <typename T>
Sqrt(const vtkm::Vec<T,4> & x)993 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Sqrt(
994 const vtkm::Vec<T, 4>& x)
995 {
996 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
997 vtkm::Sqrt(x[0]), vtkm::Sqrt(x[1]), vtkm::Sqrt(x[2]), vtkm::Sqrt(x[3]));
998 }
999 template <typename T>
Sqrt(const vtkm::Vec<T,3> & x)1000 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Sqrt(
1001 const vtkm::Vec<T, 3>& x)
1002 {
1003 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1004 vtkm::Sqrt(x[0]), vtkm::Sqrt(x[1]), vtkm::Sqrt(x[2]));
1005 }
1006 template <typename T>
Sqrt(const vtkm::Vec<T,2> & x)1007 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Sqrt(
1008 const vtkm::Vec<T, 2>& x)
1009 {
1010 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Sqrt(x[0]),
1011 vtkm::Sqrt(x[1]));
1012 }
1013
1014 /// Compute the reciprocal square root of \p x. The result of this function is
1015 /// equivalent to <tt>1/Sqrt(x)</tt>. However, on some devices it is faster to
1016 /// compute the reciprocal square root than the regular square root. Thus, you
1017 /// should use this function whenever dividing by the square root.
1018 ///
1019 #ifdef VTKM_CUDA
RSqrt(vtkm::Float32 x)1020 static inline VTKM_EXEC_CONT vtkm::Float32 RSqrt(vtkm::Float32 x)
1021 {
1022 return rsqrtf(x);
1023 }
RSqrt(vtkm::Float64 x)1024 static inline VTKM_EXEC_CONT vtkm::Float64 RSqrt(vtkm::Float64 x)
1025 {
1026 return rsqrt(x);
1027 }
1028 template <typename T>
RSqrt(T x)1029 static inline VTKM_EXEC_CONT vtkm::Float64 RSqrt(T x)
1030 {
1031 return rsqrt(static_cast<vtkm::Float64>(x));
1032 }
1033 #else // !VTKM_CUDA
RSqrt(vtkm::Float32 x)1034 static inline VTKM_EXEC_CONT vtkm::Float32 RSqrt(vtkm::Float32 x)
1035 {
1036 return 1 / vtkm::Sqrt(x);
1037 }
RSqrt(vtkm::Float64 x)1038 static inline VTKM_EXEC_CONT vtkm::Float64 RSqrt(vtkm::Float64 x)
1039 {
1040 return 1 / vtkm::Sqrt(x);
1041 }
1042 template <typename T>
RSqrt(T x)1043 static inline VTKM_EXEC_CONT vtkm::Float64 RSqrt(T x)
1044 {
1045 return 1 / static_cast<vtkm::Float64>(x);
1046 }
1047 #endif // !VTKM_CUDA
1048
1049 template <typename T, vtkm::IdComponent N>
RSqrt(const vtkm::Vec<T,N> & x)1050 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> RSqrt(
1051 const vtkm::Vec<T, N>& x)
1052 {
1053 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1054 for (vtkm::IdComponent index = 0; index < N; index++)
1055 {
1056 result[index] = vtkm::RSqrt(x[index]);
1057 }
1058 return result;
1059 }
1060 template <typename T>
RSqrt(const vtkm::Vec<T,4> & x)1061 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> RSqrt(
1062 const vtkm::Vec<T, 4>& x)
1063 {
1064 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1065 vtkm::RSqrt(x[0]), vtkm::RSqrt(x[1]), vtkm::RSqrt(x[2]), vtkm::RSqrt(x[3]));
1066 }
1067 template <typename T>
RSqrt(const vtkm::Vec<T,3> & x)1068 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> RSqrt(
1069 const vtkm::Vec<T, 3>& x)
1070 {
1071 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1072 vtkm::RSqrt(x[0]), vtkm::RSqrt(x[1]), vtkm::RSqrt(x[2]));
1073 }
1074 template <typename T>
RSqrt(const vtkm::Vec<T,2> & x)1075 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> RSqrt(
1076 const vtkm::Vec<T, 2>& x)
1077 {
1078 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::RSqrt(x[0]),
1079 vtkm::RSqrt(x[1]));
1080 }
1081
1082 /// Compute the cube root of \p x.
1083 ///
1084
Cbrt(vtkm::Float32 x)1085 inline VTKM_EXEC_CONT vtkm::Float32 Cbrt(vtkm::Float32 x)
1086 {
1087 #ifdef VTKM_CUDA
1088 return VTKM_CUDA_MATH_FUNCTION_32(cbrt)(x);
1089 #else
1090 return std::cbrt(x);
1091 #endif
1092 }
1093
Cbrt(vtkm::Float64 x)1094 inline VTKM_EXEC_CONT vtkm::Float64 Cbrt(vtkm::Float64 x)
1095 {
1096 #ifdef VTKM_CUDA
1097 return VTKM_CUDA_MATH_FUNCTION_64(cbrt)(x);
1098 #else
1099 return std::cbrt(x);
1100 #endif
1101 }
1102 template <typename T>
Cbrt(const T & x)1103 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Cbrt(const T& x)
1104 {
1105 using RT = typename detail::FloatingPointReturnType<T>::Type;
1106 return vtkm::Cbrt(static_cast<RT>(x));
1107 }
1108 template <typename T, vtkm::IdComponent N>
Cbrt(const vtkm::Vec<T,N> & x)1109 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Cbrt(
1110 const vtkm::Vec<T, N>& x)
1111 {
1112 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1113 for (vtkm::IdComponent index = 0; index < N; index++)
1114 {
1115 result[index] = vtkm::Cbrt(x[index]);
1116 }
1117 return result;
1118 }
1119 template <typename T>
Cbrt(const vtkm::Vec<T,4> & x)1120 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Cbrt(
1121 const vtkm::Vec<T, 4>& x)
1122 {
1123 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1124 vtkm::Cbrt(x[0]), vtkm::Cbrt(x[1]), vtkm::Cbrt(x[2]), vtkm::Cbrt(x[3]));
1125 }
1126 template <typename T>
Cbrt(const vtkm::Vec<T,3> & x)1127 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Cbrt(
1128 const vtkm::Vec<T, 3>& x)
1129 {
1130 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1131 vtkm::Cbrt(x[0]), vtkm::Cbrt(x[1]), vtkm::Cbrt(x[2]));
1132 }
1133 template <typename T>
Cbrt(const vtkm::Vec<T,2> & x)1134 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Cbrt(
1135 const vtkm::Vec<T, 2>& x)
1136 {
1137 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Cbrt(x[0]),
1138 vtkm::Cbrt(x[1]));
1139 }
1140
1141 /// Compute the reciprocal cube root of \p x. The result of this function is
1142 /// equivalent to <tt>1/Cbrt(x)</tt>. However, on some devices it is faster to
1143 /// compute the reciprocal cube root than the regular cube root. Thus, you
1144 /// should use this function whenever dividing by the cube root.
1145 ///
1146 #ifdef VTKM_CUDA
RCbrt(vtkm::Float32 x)1147 static inline VTKM_EXEC_CONT vtkm::Float32 RCbrt(vtkm::Float32 x)
1148 {
1149 return rcbrtf(x);
1150 }
RCbrt(vtkm::Float64 x)1151 static inline VTKM_EXEC_CONT vtkm::Float64 RCbrt(vtkm::Float64 x)
1152 {
1153 return rcbrt(x);
1154 }
1155 template <typename T>
RCbrt(T x)1156 static inline VTKM_EXEC_CONT vtkm::Float64 RCbrt(T x)
1157 {
1158 return rcbrt(static_cast<vtkm::Float64>(x));
1159 }
1160 #else // !VTKM_CUDA
RCbrt(vtkm::Float32 x)1161 static inline VTKM_EXEC_CONT vtkm::Float32 RCbrt(vtkm::Float32 x)
1162 {
1163 return 1 / vtkm::Cbrt(x);
1164 }
RCbrt(vtkm::Float64 x)1165 static inline VTKM_EXEC_CONT vtkm::Float64 RCbrt(vtkm::Float64 x)
1166 {
1167 return 1 / vtkm::Cbrt(x);
1168 }
1169 template <typename T>
RCbrt(T x)1170 static inline VTKM_EXEC_CONT vtkm::Float64 RCbrt(T x)
1171 {
1172 return 1 / vtkm::Cbrt(static_cast<vtkm::Float64>(x));
1173 }
1174 #endif // !VTKM_CUDA
1175
1176 template <typename T, vtkm::IdComponent N>
RCbrt(const vtkm::Vec<T,N> & x)1177 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> RCbrt(
1178 const vtkm::Vec<T, N>& x)
1179 {
1180 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1181 for (vtkm::IdComponent index = 0; index < N; index++)
1182 {
1183 result[index] = vtkm::RCbrt(x[index]);
1184 }
1185 return result;
1186 }
1187 template <typename T>
RCbrt(const vtkm::Vec<T,4> & x)1188 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> RCbrt(
1189 const vtkm::Vec<T, 4>& x)
1190 {
1191 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1192 vtkm::RCbrt(x[0]), vtkm::RCbrt(x[1]), vtkm::RCbrt(x[2]), vtkm::RCbrt(x[3]));
1193 }
1194 template <typename T>
RCbrt(const vtkm::Vec<T,3> & x)1195 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> RCbrt(
1196 const vtkm::Vec<T, 3>& x)
1197 {
1198 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1199 vtkm::RCbrt(x[0]), vtkm::RCbrt(x[1]), vtkm::RCbrt(x[2]));
1200 }
1201 template <typename T>
RCbrt(const vtkm::Vec<T,2> & x)1202 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> RCbrt(
1203 const vtkm::Vec<T, 2>& x)
1204 {
1205 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::RCbrt(x[0]),
1206 vtkm::RCbrt(x[1]));
1207 }
1208
1209 /// Computes e**\p x, the base-e exponential of \p x.
1210 ///
1211
Exp(vtkm::Float32 x)1212 inline VTKM_EXEC_CONT vtkm::Float32 Exp(vtkm::Float32 x)
1213 {
1214 #ifdef VTKM_CUDA
1215 return VTKM_CUDA_MATH_FUNCTION_32(exp)(x);
1216 #else
1217 return std::exp(x);
1218 #endif
1219 }
1220
Exp(vtkm::Float64 x)1221 inline VTKM_EXEC_CONT vtkm::Float64 Exp(vtkm::Float64 x)
1222 {
1223 #ifdef VTKM_CUDA
1224 return VTKM_CUDA_MATH_FUNCTION_64(exp)(x);
1225 #else
1226 return std::exp(x);
1227 #endif
1228 }
1229 template <typename T>
Exp(const T & x)1230 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Exp(const T& x)
1231 {
1232 using RT = typename detail::FloatingPointReturnType<T>::Type;
1233 return vtkm::Exp(static_cast<RT>(x));
1234 }
1235 template <typename T, vtkm::IdComponent N>
Exp(const vtkm::Vec<T,N> & x)1236 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Exp(
1237 const vtkm::Vec<T, N>& x)
1238 {
1239 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1240 for (vtkm::IdComponent index = 0; index < N; index++)
1241 {
1242 result[index] = vtkm::Exp(x[index]);
1243 }
1244 return result;
1245 }
1246 template <typename T>
Exp(const vtkm::Vec<T,4> & x)1247 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Exp(
1248 const vtkm::Vec<T, 4>& x)
1249 {
1250 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1251 vtkm::Exp(x[0]), vtkm::Exp(x[1]), vtkm::Exp(x[2]), vtkm::Exp(x[3]));
1252 }
1253 template <typename T>
Exp(const vtkm::Vec<T,3> & x)1254 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Exp(
1255 const vtkm::Vec<T, 3>& x)
1256 {
1257 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1258 vtkm::Exp(x[0]), vtkm::Exp(x[1]), vtkm::Exp(x[2]));
1259 }
1260 template <typename T>
Exp(const vtkm::Vec<T,2> & x)1261 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Exp(
1262 const vtkm::Vec<T, 2>& x)
1263 {
1264 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Exp(x[0]),
1265 vtkm::Exp(x[1]));
1266 }
1267
1268 /// Computes 2**\p x, the base-2 exponential of \p x.
1269 ///
1270
Exp2(vtkm::Float32 x)1271 inline VTKM_EXEC_CONT vtkm::Float32 Exp2(vtkm::Float32 x)
1272 {
1273 #ifdef VTKM_CUDA
1274 return VTKM_CUDA_MATH_FUNCTION_32(exp2)(x);
1275 #else
1276 return std::exp2(x);
1277 #endif
1278 }
1279
Exp2(vtkm::Float64 x)1280 inline VTKM_EXEC_CONT vtkm::Float64 Exp2(vtkm::Float64 x)
1281 {
1282 #ifdef VTKM_CUDA
1283 return VTKM_CUDA_MATH_FUNCTION_64(exp2)(x);
1284 #else
1285 return std::exp2(x);
1286 #endif
1287 }
1288 template <typename T>
Exp2(const T & x)1289 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Exp2(const T& x)
1290 {
1291 using RT = typename detail::FloatingPointReturnType<T>::Type;
1292 return vtkm::Exp2(static_cast<RT>(x));
1293 }
1294 template <typename T, vtkm::IdComponent N>
Exp2(const vtkm::Vec<T,N> & x)1295 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Exp2(
1296 const vtkm::Vec<T, N>& x)
1297 {
1298 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1299 for (vtkm::IdComponent index = 0; index < N; index++)
1300 {
1301 result[index] = vtkm::Exp2(x[index]);
1302 }
1303 return result;
1304 }
1305 template <typename T>
Exp2(const vtkm::Vec<T,4> & x)1306 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Exp2(
1307 const vtkm::Vec<T, 4>& x)
1308 {
1309 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1310 vtkm::Exp2(x[0]), vtkm::Exp2(x[1]), vtkm::Exp2(x[2]), vtkm::Exp2(x[3]));
1311 }
1312 template <typename T>
Exp2(const vtkm::Vec<T,3> & x)1313 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Exp2(
1314 const vtkm::Vec<T, 3>& x)
1315 {
1316 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1317 vtkm::Exp2(x[0]), vtkm::Exp2(x[1]), vtkm::Exp2(x[2]));
1318 }
1319 template <typename T>
Exp2(const vtkm::Vec<T,2> & x)1320 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Exp2(
1321 const vtkm::Vec<T, 2>& x)
1322 {
1323 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Exp2(x[0]),
1324 vtkm::Exp2(x[1]));
1325 }
1326
1327 /// Computes (e**\p x) - 1, the of base-e exponental of \p x then minus 1. The
1328 /// accuracy of this function is good even for very small values of x.
1329 ///
1330
ExpM1(vtkm::Float32 x)1331 inline VTKM_EXEC_CONT vtkm::Float32 ExpM1(vtkm::Float32 x)
1332 {
1333 #ifdef VTKM_CUDA
1334 return VTKM_CUDA_MATH_FUNCTION_32(expm1)(x);
1335 #else
1336 return std::expm1(x);
1337 #endif
1338 }
1339
ExpM1(vtkm::Float64 x)1340 inline VTKM_EXEC_CONT vtkm::Float64 ExpM1(vtkm::Float64 x)
1341 {
1342 #ifdef VTKM_CUDA
1343 return VTKM_CUDA_MATH_FUNCTION_64(expm1)(x);
1344 #else
1345 return std::expm1(x);
1346 #endif
1347 }
1348 template <typename T>
ExpM1(const T & x)1349 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type ExpM1(const T& x)
1350 {
1351 using RT = typename detail::FloatingPointReturnType<T>::Type;
1352 return vtkm::ExpM1(static_cast<RT>(x));
1353 }
1354 template <typename T, vtkm::IdComponent N>
ExpM1(const vtkm::Vec<T,N> & x)1355 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> ExpM1(
1356 const vtkm::Vec<T, N>& x)
1357 {
1358 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1359 for (vtkm::IdComponent index = 0; index < N; index++)
1360 {
1361 result[index] = vtkm::ExpM1(x[index]);
1362 }
1363 return result;
1364 }
1365 template <typename T>
ExpM1(const vtkm::Vec<T,4> & x)1366 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> ExpM1(
1367 const vtkm::Vec<T, 4>& x)
1368 {
1369 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1370 vtkm::ExpM1(x[0]), vtkm::ExpM1(x[1]), vtkm::ExpM1(x[2]), vtkm::ExpM1(x[3]));
1371 }
1372 template <typename T>
ExpM1(const vtkm::Vec<T,3> & x)1373 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> ExpM1(
1374 const vtkm::Vec<T, 3>& x)
1375 {
1376 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1377 vtkm::ExpM1(x[0]), vtkm::ExpM1(x[1]), vtkm::ExpM1(x[2]));
1378 }
1379 template <typename T>
ExpM1(const vtkm::Vec<T,2> & x)1380 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> ExpM1(
1381 const vtkm::Vec<T, 2>& x)
1382 {
1383 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::ExpM1(x[0]),
1384 vtkm::ExpM1(x[1]));
1385 }
1386
1387 /// Computes 10**\p x, the base-10 exponential of \p x.
1388 ///
1389 #ifdef VTKM_CUDA
Exp10(vtkm::Float32 x)1390 static inline VTKM_EXEC_CONT vtkm::Float32 Exp10(vtkm::Float32 x)
1391 {
1392 return exp10f(x);
1393 }
Exp10(vtkm::Float64 x)1394 static inline VTKM_EXEC_CONT vtkm::Float64 Exp10(vtkm::Float64 x)
1395 {
1396 return exp10(x);
1397 }
1398 template <typename T>
Exp10(T x)1399 static inline VTKM_EXEC_CONT vtkm::Float64 Exp10(T x)
1400 {
1401 return exp10(static_cast<vtkm::Float64>(x));
1402 }
1403 #else // !VTKM_CUDA
Exp10(vtkm::Float32 x)1404 static inline VTKM_EXEC_CONT vtkm::Float32 Exp10(vtkm::Float32 x)
1405 {
1406 return vtkm::Pow(10, x);
1407 }
Exp10(vtkm::Float64 x)1408 static inline VTKM_EXEC_CONT vtkm::Float64 Exp10(vtkm::Float64 x)
1409 {
1410 return vtkm::Pow(10, x);
1411 }
1412 template <typename T>
Exp10(T x)1413 static inline VTKM_EXEC_CONT vtkm::Float64 Exp10(T x)
1414 {
1415 return vtkm::Pow(10, static_cast<vtkm::Float64>(x));
1416 }
1417 #endif // !VTKM_CUDA
1418
1419 template <typename T, vtkm::IdComponent N>
Exp10(const vtkm::Vec<T,N> & x)1420 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Exp10(
1421 const vtkm::Vec<T, N>& x)
1422 {
1423 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1424 for (vtkm::IdComponent index = 0; index < N; index++)
1425 {
1426 result[index] = vtkm::Exp10(x[index]);
1427 }
1428 return result;
1429 }
1430 template <typename T>
Exp10(const vtkm::Vec<T,4> & x)1431 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Exp10(
1432 const vtkm::Vec<T, 4>& x)
1433 {
1434 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1435 vtkm::Exp10(x[0]), vtkm::Exp10(x[1]), vtkm::Exp10(x[2]), vtkm::Exp10(x[3]));
1436 }
1437 template <typename T>
Exp10(const vtkm::Vec<T,3> & x)1438 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Exp10(
1439 const vtkm::Vec<T, 3>& x)
1440 {
1441 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1442 vtkm::Exp10(x[0]), vtkm::Exp10(x[1]), vtkm::Exp10(x[2]));
1443 }
1444 template <typename T>
Exp10(const vtkm::Vec<T,2> & x)1445 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Exp10(
1446 const vtkm::Vec<T, 2>& x)
1447 {
1448 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Exp10(x[0]),
1449 vtkm::Exp10(x[1]));
1450 }
1451
1452 /// Computes the natural logarithm of \p x.
1453 ///
1454
Log(vtkm::Float32 x)1455 inline VTKM_EXEC_CONT vtkm::Float32 Log(vtkm::Float32 x)
1456 {
1457 #ifdef VTKM_CUDA
1458 return VTKM_CUDA_MATH_FUNCTION_32(log)(x);
1459 #else
1460 return std::log(x);
1461 #endif
1462 }
1463
Log(vtkm::Float64 x)1464 inline VTKM_EXEC_CONT vtkm::Float64 Log(vtkm::Float64 x)
1465 {
1466 #ifdef VTKM_CUDA
1467 return VTKM_CUDA_MATH_FUNCTION_64(log)(x);
1468 #else
1469 return std::log(x);
1470 #endif
1471 }
1472 template <typename T>
Log(const T & x)1473 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Log(const T& x)
1474 {
1475 using RT = typename detail::FloatingPointReturnType<T>::Type;
1476 return vtkm::Log(static_cast<RT>(x));
1477 }
1478 template <typename T, vtkm::IdComponent N>
Log(const vtkm::Vec<T,N> & x)1479 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Log(
1480 const vtkm::Vec<T, N>& x)
1481 {
1482 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1483 for (vtkm::IdComponent index = 0; index < N; index++)
1484 {
1485 result[index] = vtkm::Log(x[index]);
1486 }
1487 return result;
1488 }
1489 template <typename T>
Log(const vtkm::Vec<T,4> & x)1490 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Log(
1491 const vtkm::Vec<T, 4>& x)
1492 {
1493 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1494 vtkm::Log(x[0]), vtkm::Log(x[1]), vtkm::Log(x[2]), vtkm::Log(x[3]));
1495 }
1496 template <typename T>
Log(const vtkm::Vec<T,3> & x)1497 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Log(
1498 const vtkm::Vec<T, 3>& x)
1499 {
1500 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1501 vtkm::Log(x[0]), vtkm::Log(x[1]), vtkm::Log(x[2]));
1502 }
1503 template <typename T>
Log(const vtkm::Vec<T,2> & x)1504 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Log(
1505 const vtkm::Vec<T, 2>& x)
1506 {
1507 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Log(x[0]),
1508 vtkm::Log(x[1]));
1509 }
1510
1511 /// Computes the logarithm base 2 of \p x.
1512 ///
1513
Log2(vtkm::Float32 x)1514 inline VTKM_EXEC_CONT vtkm::Float32 Log2(vtkm::Float32 x)
1515 {
1516 #ifdef VTKM_CUDA
1517 return VTKM_CUDA_MATH_FUNCTION_32(log2)(x);
1518 #else
1519 return std::log2(x);
1520 #endif
1521 }
1522
Log2(vtkm::Float64 x)1523 inline VTKM_EXEC_CONT vtkm::Float64 Log2(vtkm::Float64 x)
1524 {
1525 #ifdef VTKM_CUDA
1526 return VTKM_CUDA_MATH_FUNCTION_64(log2)(x);
1527 #else
1528 return std::log2(x);
1529 #endif
1530 }
1531 template <typename T>
Log2(const T & x)1532 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Log2(const T& x)
1533 {
1534 using RT = typename detail::FloatingPointReturnType<T>::Type;
1535 return vtkm::Log2(static_cast<RT>(x));
1536 }
1537 template <typename T, vtkm::IdComponent N>
Log2(const vtkm::Vec<T,N> & x)1538 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Log2(
1539 const vtkm::Vec<T, N>& x)
1540 {
1541 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1542 for (vtkm::IdComponent index = 0; index < N; index++)
1543 {
1544 result[index] = vtkm::Log2(x[index]);
1545 }
1546 return result;
1547 }
1548 template <typename T>
Log2(const vtkm::Vec<T,4> & x)1549 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Log2(
1550 const vtkm::Vec<T, 4>& x)
1551 {
1552 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1553 vtkm::Log2(x[0]), vtkm::Log2(x[1]), vtkm::Log2(x[2]), vtkm::Log2(x[3]));
1554 }
1555 template <typename T>
Log2(const vtkm::Vec<T,3> & x)1556 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Log2(
1557 const vtkm::Vec<T, 3>& x)
1558 {
1559 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1560 vtkm::Log2(x[0]), vtkm::Log2(x[1]), vtkm::Log2(x[2]));
1561 }
1562 template <typename T>
Log2(const vtkm::Vec<T,2> & x)1563 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Log2(
1564 const vtkm::Vec<T, 2>& x)
1565 {
1566 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Log2(x[0]),
1567 vtkm::Log2(x[1]));
1568 }
1569
1570 /// Computes the logarithm base 10 of \p x.
1571 ///
1572
Log10(vtkm::Float32 x)1573 inline VTKM_EXEC_CONT vtkm::Float32 Log10(vtkm::Float32 x)
1574 {
1575 #ifdef VTKM_CUDA
1576 return VTKM_CUDA_MATH_FUNCTION_32(log10)(x);
1577 #else
1578 return std::log10(x);
1579 #endif
1580 }
1581
Log10(vtkm::Float64 x)1582 inline VTKM_EXEC_CONT vtkm::Float64 Log10(vtkm::Float64 x)
1583 {
1584 #ifdef VTKM_CUDA
1585 return VTKM_CUDA_MATH_FUNCTION_64(log10)(x);
1586 #else
1587 return std::log10(x);
1588 #endif
1589 }
1590 template <typename T>
Log10(const T & x)1591 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Log10(const T& x)
1592 {
1593 using RT = typename detail::FloatingPointReturnType<T>::Type;
1594 return vtkm::Log10(static_cast<RT>(x));
1595 }
1596 template <typename T, vtkm::IdComponent N>
Log10(const vtkm::Vec<T,N> & x)1597 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Log10(
1598 const vtkm::Vec<T, N>& x)
1599 {
1600 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1601 for (vtkm::IdComponent index = 0; index < N; index++)
1602 {
1603 result[index] = vtkm::Log10(x[index]);
1604 }
1605 return result;
1606 }
1607 template <typename T>
Log10(const vtkm::Vec<T,4> & x)1608 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Log10(
1609 const vtkm::Vec<T, 4>& x)
1610 {
1611 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1612 vtkm::Log10(x[0]), vtkm::Log10(x[1]), vtkm::Log10(x[2]), vtkm::Log10(x[3]));
1613 }
1614 template <typename T>
Log10(const vtkm::Vec<T,3> & x)1615 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Log10(
1616 const vtkm::Vec<T, 3>& x)
1617 {
1618 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1619 vtkm::Log10(x[0]), vtkm::Log10(x[1]), vtkm::Log10(x[2]));
1620 }
1621 template <typename T>
Log10(const vtkm::Vec<T,2> & x)1622 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Log10(
1623 const vtkm::Vec<T, 2>& x)
1624 {
1625 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Log10(x[0]),
1626 vtkm::Log10(x[1]));
1627 }
1628
1629 /// Computes the value of log(1+x) accurately for very small values of x.
1630 ///
1631
Log1P(vtkm::Float32 x)1632 inline VTKM_EXEC_CONT vtkm::Float32 Log1P(vtkm::Float32 x)
1633 {
1634 #ifdef VTKM_CUDA
1635 return VTKM_CUDA_MATH_FUNCTION_32(log1p)(x);
1636 #else
1637 return std::log1p(x);
1638 #endif
1639 }
1640
Log1P(vtkm::Float64 x)1641 inline VTKM_EXEC_CONT vtkm::Float64 Log1P(vtkm::Float64 x)
1642 {
1643 #ifdef VTKM_CUDA
1644 return VTKM_CUDA_MATH_FUNCTION_64(log1p)(x);
1645 #else
1646 return std::log1p(x);
1647 #endif
1648 }
1649 template <typename T>
Log1P(const T & x)1650 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Log1P(const T& x)
1651 {
1652 using RT = typename detail::FloatingPointReturnType<T>::Type;
1653 return vtkm::Log1P(static_cast<RT>(x));
1654 }
1655 template <typename T, vtkm::IdComponent N>
Log1P(const vtkm::Vec<T,N> & x)1656 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Log1P(
1657 const vtkm::Vec<T, N>& x)
1658 {
1659 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
1660 for (vtkm::IdComponent index = 0; index < N; index++)
1661 {
1662 result[index] = vtkm::Log1P(x[index]);
1663 }
1664 return result;
1665 }
1666 template <typename T>
Log1P(const vtkm::Vec<T,4> & x)1667 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Log1P(
1668 const vtkm::Vec<T, 4>& x)
1669 {
1670 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
1671 vtkm::Log1P(x[0]), vtkm::Log1P(x[1]), vtkm::Log1P(x[2]), vtkm::Log1P(x[3]));
1672 }
1673 template <typename T>
Log1P(const vtkm::Vec<T,3> & x)1674 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Log1P(
1675 const vtkm::Vec<T, 3>& x)
1676 {
1677 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
1678 vtkm::Log1P(x[0]), vtkm::Log1P(x[1]), vtkm::Log1P(x[2]));
1679 }
1680 template <typename T>
Log1P(const vtkm::Vec<T,2> & x)1681 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Log1P(
1682 const vtkm::Vec<T, 2>& x)
1683 {
1684 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Log1P(x[0]),
1685 vtkm::Log1P(x[1]));
1686 }
1687
1688 //-----------------------------------------------------------------------------
1689 /// Returns \p x or \p y, whichever is larger.
1690 ///
1691 template <typename T>
1692 static inline VTKM_EXEC_CONT T Max(const T& x, const T& y);
1693 #ifdef VTKM_USE_STL
Max(vtkm::Float32 x,vtkm::Float32 y)1694 static inline VTKM_EXEC_CONT vtkm::Float32 Max(vtkm::Float32 x, vtkm::Float32 y)
1695 {
1696 return (std::max)(x, y);
1697 }
Max(vtkm::Float64 x,vtkm::Float64 y)1698 static inline VTKM_EXEC_CONT vtkm::Float64 Max(vtkm::Float64 x, vtkm::Float64 y)
1699 {
1700 return (std::max)(x, y);
1701 }
1702 #else // !VTKM_USE_STL
Max(vtkm::Float32 x,vtkm::Float32 y)1703 static inline VTKM_EXEC_CONT vtkm::Float32 Max(vtkm::Float32 x, vtkm::Float32 y)
1704 {
1705 #ifdef VTKM_CUDA
1706 return VTKM_CUDA_MATH_FUNCTION_32(fmax)(x, y);
1707 #else
1708 return std::fmax(x, y);
1709 #endif
1710 }
Max(vtkm::Float64 x,vtkm::Float64 y)1711 static inline VTKM_EXEC_CONT vtkm::Float64 Max(vtkm::Float64 x, vtkm::Float64 y)
1712 {
1713 #ifdef VTKM_CUDA
1714 return VTKM_CUDA_MATH_FUNCTION_64(fmax)(x, y);
1715 #else
1716 return std::fmax(x, y);
1717 #endif
1718 }
1719 #endif // !VTKM_USE_STL
1720
1721 /// Returns \p x or \p y, whichever is smaller.
1722 ///
1723 template <typename T>
1724 static inline VTKM_EXEC_CONT T Min(const T& x, const T& y);
1725 #if defined(VTKM_USE_STL) && !defined(VTKM_HIP)
Min(vtkm::Float32 x,vtkm::Float32 y)1726 static inline VTKM_EXEC_CONT vtkm::Float32 Min(vtkm::Float32 x, vtkm::Float32 y)
1727 {
1728 return (std::min)(x, y);
1729 }
Min(vtkm::Float64 x,vtkm::Float64 y)1730 static inline VTKM_EXEC_CONT vtkm::Float64 Min(vtkm::Float64 x, vtkm::Float64 y)
1731 {
1732 return (std::min)(x, y);
1733 }
1734 #else // !VTKM_USE_STL OR HIP
Min(vtkm::Float32 x,vtkm::Float32 y)1735 static inline VTKM_EXEC_CONT vtkm::Float32 Min(vtkm::Float32 x, vtkm::Float32 y)
1736 {
1737 #ifdef VTKM_CUDA
1738 return VTKM_CUDA_MATH_FUNCTION_32(fmin)(x, y);
1739 #else
1740 return std::fmin(x, y);
1741 #endif
1742 }
Min(vtkm::Float64 x,vtkm::Float64 y)1743 static inline VTKM_EXEC_CONT vtkm::Float64 Min(vtkm::Float64 x, vtkm::Float64 y)
1744 {
1745 #ifdef VTKM_CUDA
1746 return VTKM_CUDA_MATH_FUNCTION_64(fmin)(x, y);
1747 #else
1748 return std::fmin(x, y);
1749 #endif
1750 }
1751 #endif // !VTKM_USE_STL
1752
1753 namespace detail
1754 {
1755
1756 template <typename T>
Max(T x,T y,vtkm::TypeTraitsScalarTag)1757 static inline VTKM_EXEC_CONT T Max(T x, T y, vtkm::TypeTraitsScalarTag)
1758 {
1759 return (x < y) ? y : x;
1760 }
1761
1762 template <typename T>
Max(const T & x,const T & y,vtkm::TypeTraitsVectorTag)1763 static inline VTKM_EXEC_CONT T Max(const T& x, const T& y, vtkm::TypeTraitsVectorTag)
1764 {
1765 using Traits = vtkm::VecTraits<T>;
1766 T result;
1767 for (vtkm::IdComponent index = 0; index < Traits::NUM_COMPONENTS; index++)
1768 {
1769 Traits::SetComponent(
1770 result, index, vtkm::Max(Traits::GetComponent(x, index), Traits::GetComponent(y, index)));
1771 }
1772 return result;
1773 }
1774
1775 template <typename T>
Min(T x,T y,vtkm::TypeTraitsScalarTag)1776 static inline VTKM_EXEC_CONT T Min(T x, T y, vtkm::TypeTraitsScalarTag)
1777 {
1778 return (x < y) ? x : y;
1779 }
1780
1781 template <typename T>
Min(const T & x,const T & y,vtkm::TypeTraitsVectorTag)1782 static inline VTKM_EXEC_CONT T Min(const T& x, const T& y, vtkm::TypeTraitsVectorTag)
1783 {
1784 using Traits = vtkm::VecTraits<T>;
1785 T result;
1786 for (vtkm::IdComponent index = 0; index < Traits::NUM_COMPONENTS; index++)
1787 {
1788 Traits::SetComponent(
1789 result, index, vtkm::Min(Traits::GetComponent(x, index), Traits::GetComponent(y, index)));
1790 }
1791 return result;
1792 }
1793
1794 } // namespace detail
1795
1796 /// Returns \p x or \p y, whichever is larger.
1797 ///
1798 template <typename T>
Max(const T & x,const T & y)1799 static inline VTKM_EXEC_CONT T Max(const T& x, const T& y)
1800 {
1801 return detail::Max(x, y, typename vtkm::TypeTraits<T>::DimensionalityTag());
1802 }
1803
1804 /// Returns \p x or \p y, whichever is smaller.
1805 ///
1806 template <typename T>
Min(const T & x,const T & y)1807 static inline VTKM_EXEC_CONT T Min(const T& x, const T& y)
1808 {
1809 return detail::Min(x, y, typename vtkm::TypeTraits<T>::DimensionalityTag());
1810 }
1811
1812 /// Clamp \p x to the given range.
1813 ///
1814
Clamp(vtkm::Float32 x,vtkm::Float32 lo,vtkm::Float32 hi)1815 inline VTKM_EXEC_CONT vtkm::Float32 Clamp(vtkm::Float32 x, vtkm::Float32 lo, vtkm::Float32 hi)
1816 {
1817 return x > lo ? (x < hi ? x : hi) : lo;
1818 }
1819
Clamp(vtkm::Float64 x,vtkm::Float64 lo,vtkm::Float64 hi)1820 inline VTKM_EXEC_CONT vtkm::Float64 Clamp(vtkm::Float64 x, vtkm::Float64 lo, vtkm::Float64 hi)
1821 {
1822 return x > lo ? (x < hi ? x : hi) : lo;
1823 }
1824 //-----------------------------------------------------------------------------
1825
1826 //#ifdef VTKM_CUDA
1827 #define VTKM_USE_IEEE_NONFINITE
1828 //#endif
1829
1830 #ifdef VTKM_USE_IEEE_NONFINITE
1831
1832 namespace detail
1833 {
1834
1835 union IEEE754Bits32 {
1836 vtkm::UInt32 bits;
1837 vtkm::Float32 scalar;
1838 };
1839 #define VTKM_NAN_BITS_32 0x7FC00000U
1840 #define VTKM_INF_BITS_32 0x7F800000U
1841 #define VTKM_NEG_INF_BITS_32 0xFF800000U
1842 #define VTKM_EPSILON_32 1e-5f
1843
1844 union IEEE754Bits64 {
1845 vtkm::UInt64 bits;
1846 vtkm::Float64 scalar;
1847 };
1848 #define VTKM_NAN_BITS_64 0x7FF8000000000000ULL
1849 #define VTKM_INF_BITS_64 0x7FF0000000000000ULL
1850 #define VTKM_NEG_INF_BITS_64 0xFFF0000000000000ULL
1851 #define VTKM_EPSILON_64 1e-9
1852
1853 template <typename T>
1854 struct FloatLimits;
1855
1856 template <>
1857 struct FloatLimits<vtkm::Float32>
1858 {
1859 using BitsType = vtkm::detail::IEEE754Bits32;
1860
1861 VTKM_EXEC_CONT
1862 static vtkm::Float32 Nan()
1863 {
1864 BitsType nan = { VTKM_NAN_BITS_32 };
1865 return nan.scalar;
1866 }
1867
1868 VTKM_EXEC_CONT
1869 static vtkm::Float32 Infinity()
1870 {
1871 BitsType inf = { VTKM_INF_BITS_32 };
1872 return inf.scalar;
1873 }
1874
1875 VTKM_EXEC_CONT
1876 static vtkm::Float32 NegativeInfinity()
1877 {
1878 BitsType neginf = { VTKM_NEG_INF_BITS_32 };
1879 return neginf.scalar;
1880 }
1881
1882 VTKM_EXEC_CONT
1883 static vtkm::Float32 Epsilon() { return VTKM_EPSILON_32; }
1884 };
1885
1886 template <int N>
1887 struct FloatLimits<vtkm::Vec<vtkm::Float32, N>>
1888 {
1889 using BitsType = vtkm::detail::IEEE754Bits32;
1890
1891 VTKM_EXEC_CONT
1892 static vtkm::Vec<vtkm::Float32, N> Nan()
1893 {
1894 BitsType nan = { VTKM_NAN_BITS_32 };
1895 return vtkm::Vec<vtkm::Float32, N>(nan.scalar);
1896 }
1897
1898 VTKM_EXEC_CONT
1899 static vtkm::Vec<vtkm::Float32, N> Infinity()
1900 {
1901 BitsType inf = { VTKM_INF_BITS_32 };
1902 return vtkm::Vec<vtkm::Float32, N>(inf.scalar);
1903 }
1904
1905 VTKM_EXEC_CONT
1906 static vtkm::Vec<vtkm::Float32, N> NegativeInfinity()
1907 {
1908 BitsType neginf = { VTKM_NEG_INF_BITS_32 };
1909 return vtkm::Vec<vtkm::Float32, N>(neginf.scalar);
1910 }
1911
1912 VTKM_EXEC_CONT
1913 static vtkm::Vec<vtkm::Float32, N> Epsilon()
1914 {
1915 return vtkm::Vec<vtkm::Float32, N>(VTKM_EPSILON_32);
1916 }
1917 };
1918
1919 template <>
1920 struct FloatLimits<vtkm::Float64>
1921 {
1922 using BitsType = vtkm::detail::IEEE754Bits64;
1923
1924 VTKM_EXEC_CONT
1925 static vtkm::Float64 Nan()
1926 {
1927 BitsType nan = { VTKM_NAN_BITS_64 };
1928 return nan.scalar;
1929 }
1930
1931 VTKM_EXEC_CONT
1932 static vtkm::Float64 Infinity()
1933 {
1934 BitsType inf = { VTKM_INF_BITS_64 };
1935 return inf.scalar;
1936 }
1937
1938 VTKM_EXEC_CONT
1939 static vtkm::Float64 NegativeInfinity()
1940 {
1941 BitsType neginf = { VTKM_NEG_INF_BITS_64 };
1942 return neginf.scalar;
1943 }
1944
1945 VTKM_EXEC_CONT
1946 static vtkm::Float64 Epsilon() { return VTKM_EPSILON_64; }
1947 };
1948
1949 template <int N>
1950 struct FloatLimits<vtkm::Vec<vtkm::Float64, N>>
1951 {
1952 using BitsType = vtkm::detail::IEEE754Bits64;
1953
1954 VTKM_EXEC_CONT
1955 static vtkm::Vec<vtkm::Float64, N> Nan()
1956 {
1957 BitsType nan = { VTKM_NAN_BITS_64 };
1958 return vtkm::Vec<vtkm::Float64, N>(nan.scalar);
1959 }
1960
1961 VTKM_EXEC_CONT
1962 static vtkm::Vec<vtkm::Float64, N> Infinity()
1963 {
1964 BitsType inf = { VTKM_INF_BITS_64 };
1965 return vtkm::Vec<vtkm::Float64, N>(inf.scalar);
1966 }
1967
1968 VTKM_EXEC_CONT
1969 static vtkm::Vec<vtkm::Float64, N> NegativeInfinity()
1970 {
1971 BitsType neginf = { VTKM_NEG_INF_BITS_64 };
1972 return vtkm::Vec<vtkm::Float64, N>(neginf.scalar);
1973 }
1974
1975 VTKM_EXEC_CONT
1976 static vtkm::Vec<vtkm::Float64, N> Epsilon()
1977 {
1978 return vtkm::Vec<vtkm::Float64, N>(VTKM_EPSILON_64);
1979 }
1980 };
1981
1982 #undef VTKM_NAN_BITS_32
1983 #undef VTKM_INF_BITS_32
1984 #undef VTKM_NEG_INF_BITS_32
1985 #undef VTKM_EPSILON_32
1986 #undef VTKM_NAN_BITS_64
1987 #undef VTKM_INF_BITS_64
1988 #undef VTKM_NEG_INF_BITS_64
1989 #undef VTKM_EPSILON_64
1990
1991 } // namespace detail
1992
1993 /// Returns the representation for not-a-number (NaN).
1994 ///
1995 template <typename T>
1996 static inline VTKM_EXEC_CONT T Nan()
1997 {
1998 return detail::FloatLimits<T>::Nan();
1999 }
2000
2001 /// Returns the representation for infinity.
2002 ///
2003 template <typename T>
2004 static inline VTKM_EXEC_CONT T Infinity()
2005 {
2006 return detail::FloatLimits<T>::Infinity();
2007 }
2008
2009 /// Returns the representation for negative infinity.
2010 ///
2011 template <typename T>
2012 static inline VTKM_EXEC_CONT T NegativeInfinity()
2013 {
2014 return detail::FloatLimits<T>::NegativeInfinity();
2015 }
2016
2017 /// Returns the difference between 1 and the least value greater than 1
2018 /// that is representable.
2019 ///
2020 template <typename T>
2021 static inline VTKM_EXEC_CONT T Epsilon()
2022 {
2023 return detail::FloatLimits<T>::Epsilon();
2024 }
2025
2026 #else // !VTKM_USE_IEEE_NONFINITE
2027
2028 /// Returns the representation for not-a-number (NaN).
2029 ///
2030 template <typename T>
2031 static inline VTKM_EXEC_CONT T Nan()
2032 {
2033 return std::numeric_limits<T>::quiet_NaN();
2034 }
2035
2036 /// Returns the representation for infinity.
2037 ///
2038 template <typename T>
2039 static inline VTKM_EXEC_CONT T Infinity()
2040 {
2041 return std::numeric_limits<T>::infinity();
2042 }
2043
2044 /// Returns the representation for negative infinity.
2045 ///
2046 template <typename T>
2047 static inline VTKM_EXEC_CONT T NegativeInfinity()
2048 {
2049 return -std::numeric_limits<T>::infinity();
2050 }
2051
2052 /// Returns the difference between 1 and the least value greater than 1
2053 /// that is representable.
2054 ///
2055 template <typename T>
2056 static inline VTKM_EXEC_CONT T Epsilon()
2057 {
2058 return std::numeric_limits<T>::epsilon();
2059 }
2060 #endif // !VTKM_USE_IEEE_NONFINITE
2061
2062 /// Returns the representation for not-a-number (NaN).
2063 ///
2064 static inline VTKM_EXEC_CONT vtkm::Float32 Nan32()
2065 {
2066 return vtkm::Nan<vtkm::Float32>();
2067 }
2068 static inline VTKM_EXEC_CONT vtkm::Float64 Nan64()
2069 {
2070 return vtkm::Nan<vtkm::Float64>();
2071 }
2072
2073 /// Returns the representation for infinity.
2074 ///
2075 static inline VTKM_EXEC_CONT vtkm::Float32 Infinity32()
2076 {
2077 return vtkm::Infinity<vtkm::Float32>();
2078 }
2079 static inline VTKM_EXEC_CONT vtkm::Float64 Infinity64()
2080 {
2081 return vtkm::Infinity<vtkm::Float64>();
2082 }
2083
2084 /// Returns the representation for negative infinity.
2085 ///
2086 static inline VTKM_EXEC_CONT vtkm::Float32 NegativeInfinity32()
2087 {
2088 return vtkm::NegativeInfinity<vtkm::Float32>();
2089 }
2090 static inline VTKM_EXEC_CONT vtkm::Float64 NegativeInfinity64()
2091 {
2092 return vtkm::NegativeInfinity<vtkm::Float64>();
2093 }
2094
2095 /// Returns the difference between 1 and the least value greater than 1
2096 /// that is representable.
2097 ///
2098 static inline VTKM_EXEC_CONT vtkm::Float32 Epsilon32()
2099 {
2100 return vtkm::Epsilon<vtkm::Float32>();
2101 }
2102 static inline VTKM_EXEC_CONT vtkm::Float64 Epsilon64()
2103 {
2104 return vtkm::Epsilon<vtkm::Float64>();
2105 }
2106
2107 //-----------------------------------------------------------------------------
2108 /// Returns true if \p x is not a number.
2109 ///
2110 template <typename T>
2111 static inline VTKM_EXEC_CONT bool IsNan(T x)
2112 {
2113 #ifndef VTKM_CUDA
2114 using std::isnan;
2115 #endif
2116 return (isnan(x) != 0);
2117 }
2118
2119 /// Returns true if \p x is positive or negative infinity.
2120 ///
2121 template <typename T>
2122 static inline VTKM_EXEC_CONT bool IsInf(T x)
2123 {
2124 #ifndef VTKM_CUDA
2125 using std::isinf;
2126 #endif
2127 return (isinf(x) != 0);
2128 }
2129
2130 /// Returns true if \p x is a normal number (not NaN or infinite).
2131 ///
2132 template <typename T>
2133 static inline VTKM_EXEC_CONT bool IsFinite(T x)
2134 {
2135 #ifndef VTKM_CUDA
2136 using std::isfinite;
2137 #endif
2138 return (isfinite(x) != 0);
2139 }
2140
2141 //-----------------------------------------------------------------------------
2142 /// Round \p x to the smallest integer value not less than x.
2143 ///
2144
2145 inline VTKM_EXEC_CONT vtkm::Float32 Ceil(vtkm::Float32 x)
2146 {
2147 #ifdef VTKM_CUDA
2148 return VTKM_CUDA_MATH_FUNCTION_32(ceil)(x);
2149 #else
2150 return std::ceil(x);
2151 #endif
2152 }
2153
2154 inline VTKM_EXEC_CONT vtkm::Float64 Ceil(vtkm::Float64 x)
2155 {
2156 #ifdef VTKM_CUDA
2157 return VTKM_CUDA_MATH_FUNCTION_64(ceil)(x);
2158 #else
2159 return std::ceil(x);
2160 #endif
2161 }
2162 template <typename T>
2163 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Ceil(const T& x)
2164 {
2165 using RT = typename detail::FloatingPointReturnType<T>::Type;
2166 return vtkm::Ceil(static_cast<RT>(x));
2167 }
2168 template <typename T, vtkm::IdComponent N>
2169 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Ceil(
2170 const vtkm::Vec<T, N>& x)
2171 {
2172 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
2173 for (vtkm::IdComponent index = 0; index < N; index++)
2174 {
2175 result[index] = vtkm::Ceil(x[index]);
2176 }
2177 return result;
2178 }
2179 template <typename T>
2180 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Ceil(
2181 const vtkm::Vec<T, 4>& x)
2182 {
2183 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
2184 vtkm::Ceil(x[0]), vtkm::Ceil(x[1]), vtkm::Ceil(x[2]), vtkm::Ceil(x[3]));
2185 }
2186 template <typename T>
2187 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Ceil(
2188 const vtkm::Vec<T, 3>& x)
2189 {
2190 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
2191 vtkm::Ceil(x[0]), vtkm::Ceil(x[1]), vtkm::Ceil(x[2]));
2192 }
2193 template <typename T>
2194 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Ceil(
2195 const vtkm::Vec<T, 2>& x)
2196 {
2197 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Ceil(x[0]),
2198 vtkm::Ceil(x[1]));
2199 }
2200
2201 /// Round \p x to the largest integer value not greater than x.
2202 ///
2203
2204 inline VTKM_EXEC_CONT vtkm::Float32 Floor(vtkm::Float32 x)
2205 {
2206 #ifdef VTKM_CUDA
2207 return VTKM_CUDA_MATH_FUNCTION_32(floor)(x);
2208 #else
2209 return std::floor(x);
2210 #endif
2211 }
2212
2213 inline VTKM_EXEC_CONT vtkm::Float64 Floor(vtkm::Float64 x)
2214 {
2215 #ifdef VTKM_CUDA
2216 return VTKM_CUDA_MATH_FUNCTION_64(floor)(x);
2217 #else
2218 return std::floor(x);
2219 #endif
2220 }
2221 template <typename T>
2222 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Floor(const T& x)
2223 {
2224 using RT = typename detail::FloatingPointReturnType<T>::Type;
2225 return vtkm::Floor(static_cast<RT>(x));
2226 }
2227 template <typename T, vtkm::IdComponent N>
2228 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Floor(
2229 const vtkm::Vec<T, N>& x)
2230 {
2231 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
2232 for (vtkm::IdComponent index = 0; index < N; index++)
2233 {
2234 result[index] = vtkm::Floor(x[index]);
2235 }
2236 return result;
2237 }
2238 template <typename T>
2239 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Floor(
2240 const vtkm::Vec<T, 4>& x)
2241 {
2242 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
2243 vtkm::Floor(x[0]), vtkm::Floor(x[1]), vtkm::Floor(x[2]), vtkm::Floor(x[3]));
2244 }
2245 template <typename T>
2246 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Floor(
2247 const vtkm::Vec<T, 3>& x)
2248 {
2249 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
2250 vtkm::Floor(x[0]), vtkm::Floor(x[1]), vtkm::Floor(x[2]));
2251 }
2252 template <typename T>
2253 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Floor(
2254 const vtkm::Vec<T, 2>& x)
2255 {
2256 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Floor(x[0]),
2257 vtkm::Floor(x[1]));
2258 }
2259
2260 /// Round \p x to the nearest integral value.
2261 ///
2262
2263 inline VTKM_EXEC_CONT vtkm::Float32 Round(vtkm::Float32 x)
2264 {
2265 #ifdef VTKM_CUDA
2266 return VTKM_CUDA_MATH_FUNCTION_32(round)(x);
2267 #else
2268 return std::round(x);
2269 #endif
2270 }
2271
2272 inline VTKM_EXEC_CONT vtkm::Float64 Round(vtkm::Float64 x)
2273 {
2274 #ifdef VTKM_CUDA
2275 return VTKM_CUDA_MATH_FUNCTION_64(round)(x);
2276 #else
2277 return std::round(x);
2278 #endif
2279 }
2280 template <typename T>
2281 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Round(const T& x)
2282 {
2283 using RT = typename detail::FloatingPointReturnType<T>::Type;
2284 return vtkm::Round(static_cast<RT>(x));
2285 }
2286 template <typename T, vtkm::IdComponent N>
2287 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> Round(
2288 const vtkm::Vec<T, N>& x)
2289 {
2290 vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, N> result;
2291 for (vtkm::IdComponent index = 0; index < N; index++)
2292 {
2293 result[index] = vtkm::Round(x[index]);
2294 }
2295 return result;
2296 }
2297 template <typename T>
2298 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4> Round(
2299 const vtkm::Vec<T, 4>& x)
2300 {
2301 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 4>(
2302 vtkm::Round(x[0]), vtkm::Round(x[1]), vtkm::Round(x[2]), vtkm::Round(x[3]));
2303 }
2304 template <typename T>
2305 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3> Round(
2306 const vtkm::Vec<T, 3>& x)
2307 {
2308 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 3>(
2309 vtkm::Round(x[0]), vtkm::Round(x[1]), vtkm::Round(x[2]));
2310 }
2311 template <typename T>
2312 static inline VTKM_EXEC_CONT vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2> Round(
2313 const vtkm::Vec<T, 2>& x)
2314 {
2315 return vtkm::Vec<typename detail::FloatingPointReturnType<T>::Type, 2>(vtkm::Round(x[0]),
2316 vtkm::Round(x[1]));
2317 }
2318
2319 //-----------------------------------------------------------------------------
2320 /// Computes the remainder on division of 2 floating point numbers. The return
2321 /// value is \p numerator - n \p denominator, where n is the quotient of \p
2322 /// numerator divided by \p denominator rounded towards zero to an integer. For
2323 /// example, <tt>FMod(6.5, 2.3)</tt> returns 1.9, which is 6.5 - 2*2.3.
2324 ///
2325 static inline VTKM_EXEC_CONT vtkm::Float32 FMod(vtkm::Float32 x, vtkm::Float32 y)
2326 {
2327 #ifdef VTKM_CUDA
2328 return VTKM_CUDA_MATH_FUNCTION_32(fmod)(x, y);
2329 #else
2330 return std::fmod(x, y);
2331 #endif
2332 }
2333 static inline VTKM_EXEC_CONT vtkm::Float64 FMod(vtkm::Float64 x, vtkm::Float64 y)
2334 {
2335 #ifdef VTKM_CUDA
2336 return VTKM_CUDA_MATH_FUNCTION_64(fmod)(x, y);
2337 #else
2338 return std::fmod(x, y);
2339 #endif
2340 }
2341
2342 /// Computes the remainder on division of 2 floating point numbers. The return
2343 /// value is \p numerator - n \p denominator, where n is the quotient of \p
2344 /// numerator divided by \p denominator rounded towards the nearest integer
2345 /// (instead of toward zero like FMod). For example, <tt>FMod(6.5, 2.3)</tt>
2346 /// returns -0.4, which is 6.5 - 3*2.3.
2347 ///
2348 #ifdef VTKM_MSVC
2349 template <typename T>
2350 static inline VTKM_EXEC_CONT T Remainder(T numerator, T denominator)
2351 {
2352 T quotient = vtkm::Round(numerator / denominator);
2353 return numerator - quotient * denominator;
2354 }
2355 #else // !VTKM_MSVC
2356 static inline VTKM_EXEC_CONT vtkm::Float32 Remainder(vtkm::Float32 x, vtkm::Float32 y)
2357 {
2358 #ifdef VTKM_CUDA
2359 return VTKM_CUDA_MATH_FUNCTION_32(remainder)(x, y);
2360 #else
2361 return std::remainder(x, y);
2362 #endif
2363 }
2364 static inline VTKM_EXEC_CONT vtkm::Float64 Remainder(vtkm::Float64 x, vtkm::Float64 y)
2365 {
2366 #ifdef VTKM_CUDA
2367 return VTKM_CUDA_MATH_FUNCTION_64(remainder)(x, y);
2368 #else
2369 return std::remainder(x, y);
2370 #endif
2371 }
2372 #endif // !VTKM_MSVC
2373
2374 /// Returns the remainder on division of 2 floating point numbers just like
2375 /// Remainder. In addition, this function also returns the \c quotient used to
2376 /// get that remainder.
2377 ///
2378 template <typename QType>
2379 static inline VTKM_EXEC_CONT vtkm::Float32 RemainderQuotient(vtkm::Float32 numerator,
2380 vtkm::Float32 denominator,
2381 QType& quotient)
2382 {
2383 int iQuotient;
2384 // See: https://github.com/ROCm-Developer-Tools/HIP/issues/2169
2385 #if defined(VTKM_CUDA) || defined(VTKM_HIP)
2386 const vtkm::Float32 result =
2387 VTKM_CUDA_MATH_FUNCTION_32(remquo)(numerator, denominator, &iQuotient);
2388 #else
2389 const vtkm::Float32 result = std::remquo(numerator, denominator, &iQuotient);
2390 #endif
2391 quotient = static_cast<QType>(iQuotient);
2392 return result;
2393 }
2394 template <typename QType>
2395 static inline VTKM_EXEC_CONT vtkm::Float64 RemainderQuotient(vtkm::Float64 numerator,
2396 vtkm::Float64 denominator,
2397 QType& quotient)
2398 {
2399 int iQuotient;
2400 #ifdef VTKM_CUDA
2401 const vtkm::Float64 result =
2402 VTKM_CUDA_MATH_FUNCTION_64(remquo)(numerator, denominator, &iQuotient);
2403 #else
2404 const vtkm::Float64 result = std::remquo(numerator, denominator, &iQuotient);
2405 #endif
2406 quotient = static_cast<QType>(iQuotient);
2407 return result;
2408 }
2409
2410 /// Gets the integral and fractional parts of \c x. The return value is the
2411 /// fractional part and \c integral is set to the integral part.
2412 ///
2413 static inline VTKM_EXEC_CONT vtkm::Float32 ModF(vtkm::Float32 x, vtkm::Float32& integral)
2414 {
2415 // See: https://github.com/ROCm-Developer-Tools/HIP/issues/2169
2416 #if defined(VTKM_CUDA) || defined(VTKM_HIP)
2417 return VTKM_CUDA_MATH_FUNCTION_32(modf)(x, &integral);
2418 #else
2419 return std::modf(x, &integral);
2420 #endif
2421 }
2422 static inline VTKM_EXEC_CONT vtkm::Float64 ModF(vtkm::Float64 x, vtkm::Float64& integral)
2423 {
2424 #if defined(VTKM_CUDA)
2425 return VTKM_CUDA_MATH_FUNCTION_64(modf)(x, &integral);
2426 #else
2427 return std::modf(x, &integral);
2428 #endif
2429 }
2430
2431 //-----------------------------------------------------------------------------
2432 /// Return the absolute value of \p x. That is, return \p x if it is positive or
2433 /// \p -x if it is negative.
2434 ///
2435 static inline VTKM_EXEC_CONT vtkm::Int32 Abs(vtkm::Int32 x)
2436 {
2437 return abs(x);
2438 }
2439 static inline VTKM_EXEC_CONT vtkm::Int64 Abs(vtkm::Int64 x)
2440 {
2441 #if VTKM_SIZE_LONG == 8
2442 return labs(x);
2443 #elif VTKM_SIZE_LONG_LONG == 8
2444 return llabs(x);
2445 #else
2446 #error Unknown size of Int64.
2447 #endif
2448 }
2449 static inline VTKM_EXEC_CONT vtkm::Float32 Abs(vtkm::Float32 x)
2450 {
2451 #ifdef VTKM_CUDA
2452 return VTKM_CUDA_MATH_FUNCTION_32(fabs)(x);
2453 #else
2454 return std::fabs(x);
2455 #endif
2456 }
2457 static inline VTKM_EXEC_CONT vtkm::Float64 Abs(vtkm::Float64 x)
2458 {
2459 #ifdef VTKM_CUDA
2460 return VTKM_CUDA_MATH_FUNCTION_64(fabs)(x);
2461 #else
2462 return std::fabs(x);
2463 #endif
2464 }
2465 template <typename T>
2466 static inline VTKM_EXEC_CONT typename detail::FloatingPointReturnType<T>::Type Abs(T x)
2467 {
2468 #ifdef VTKM_CUDA
2469 return VTKM_CUDA_MATH_FUNCTION_64(fabs)(static_cast<vtkm::Float64>(x));
2470 #else
2471 return std::fabs(static_cast<vtkm::Float64>(x));
2472 #endif
2473 }
2474 template <typename T, vtkm::IdComponent N>
2475 static inline VTKM_EXEC_CONT vtkm::Vec<T, N> Abs(const vtkm::Vec<T, N>& x)
2476 {
2477 vtkm::Vec<T, N> result;
2478 for (vtkm::IdComponent index = 0; index < N; index++)
2479 {
2480 result[index] = vtkm::Abs(x[index]);
2481 }
2482 return result;
2483 }
2484 template <typename T>
2485 static inline VTKM_EXEC_CONT vtkm::Vec<T, 4> Abs(const vtkm::Vec<T, 4>& x)
2486 {
2487 return vtkm::Vec<T, 4>(vtkm::Abs(x[0]), vtkm::Abs(x[1]), vtkm::Abs(x[2]), vtkm::Abs(x[3]));
2488 }
2489 template <typename T>
2490 static inline VTKM_EXEC_CONT vtkm::Vec<T, 3> Abs(const vtkm::Vec<T, 3>& x)
2491 {
2492 return vtkm::Vec<T, 3>(vtkm::Abs(x[0]), vtkm::Abs(x[1]), vtkm::Abs(x[2]));
2493 }
2494 template <typename T>
2495 static inline VTKM_EXEC_CONT vtkm::Vec<T, 2> Abs(const vtkm::Vec<T, 2>& x)
2496 {
2497 return vtkm::Vec<T, 2>(vtkm::Abs(x[0]), vtkm::Abs(x[1]));
2498 }
2499
2500 /// Returns a nonzero value if \p x is negative.
2501 ///
2502 static inline VTKM_EXEC_CONT vtkm::Int32 SignBit(vtkm::Float32 x)
2503 {
2504 #ifndef VTKM_CUDA
2505 using std::signbit;
2506 #endif
2507 return static_cast<vtkm::Int32>(signbit(x));
2508 }
2509 static inline VTKM_EXEC_CONT vtkm::Int32 SignBit(vtkm::Float64 x)
2510 {
2511 #ifndef VTKM_CUDA
2512 using std::signbit;
2513 #endif
2514 return static_cast<vtkm::Int32>(signbit(x));
2515 }
2516
2517 /// Returns true if \p x is less than zero, false otherwise.
2518 ///
2519 static inline VTKM_EXEC_CONT bool IsNegative(vtkm::Float32 x)
2520 {
2521 return (vtkm::SignBit(x) != 0);
2522 }
2523 static inline VTKM_EXEC_CONT bool IsNegative(vtkm::Float64 x)
2524 {
2525 return (vtkm::SignBit(x) != 0);
2526 }
2527
2528 /// Copies the sign of \p y onto \p x. If \p y is positive, returns Abs(\p x).
2529 /// If \p y is negative, returns -Abs(\p x).
2530 ///
2531 static inline VTKM_EXEC_CONT vtkm::Float32 CopySign(vtkm::Float32 x, vtkm::Float32 y)
2532 {
2533 #ifdef VTKM_CUDA
2534 return VTKM_CUDA_MATH_FUNCTION_32(copysign)(x, y);
2535 #else
2536 return std::copysign(x, y);
2537 #endif
2538 }
2539 static inline VTKM_EXEC_CONT vtkm::Float64 CopySign(vtkm::Float64 x, vtkm::Float64 y)
2540 {
2541 #ifdef VTKM_CUDA
2542 return VTKM_CUDA_MATH_FUNCTION_64(copysign)(x, y);
2543 #else
2544 return std::copysign(x, y);
2545 #endif
2546 }
2547
2548 template <typename T, vtkm::IdComponent N>
2549 static inline VTKM_EXEC_CONT vtkm::Vec<T, N> CopySign(const vtkm::Vec<T, N>& x,
2550 const vtkm::Vec<T, N>& y)
2551 {
2552 vtkm::Vec<T, N> result;
2553 for (vtkm::IdComponent index = 0; index < N; index++)
2554 {
2555 result[index] = vtkm::CopySign(x[index], y[index]);
2556 }
2557 return result;
2558 }
2559
2560 /// Decompose floating poing value
2561 ///
2562
2563 inline VTKM_EXEC_CONT vtkm::Float32 Frexp(vtkm::Float32 x, vtkm::Int32 *exponent)
2564 {
2565 // See: https://github.com/ROCm-Developer-Tools/HIP/issues/2169
2566 #if defined(VTKM_CUDA) || defined(VTKM_HIP)
2567 return VTKM_CUDA_MATH_FUNCTION_32(frexp)(x, exponent);
2568 #else
2569 return std::frexp(x, exponent);
2570 #endif
2571 }
2572
2573 inline VTKM_EXEC_CONT vtkm::Float64 Frexp(vtkm::Float64 x, vtkm::Int32 *exponent)
2574 {
2575 #ifdef VTKM_CUDA
2576 return VTKM_CUDA_MATH_FUNCTION_64(frexp)(x, exponent);
2577 #else
2578 return std::frexp(x, exponent);
2579 #endif
2580 }
2581
2582 inline VTKM_EXEC_CONT vtkm::Float32 Ldexp(vtkm::Float32 x, vtkm::Int32 exponent)
2583 {
2584 #ifdef VTKM_CUDA
2585 return VTKM_CUDA_MATH_FUNCTION_32(ldexp)(x, exponent);
2586 #else
2587 return std::ldexp(x, exponent);
2588 #endif
2589 }
2590
2591 inline VTKM_EXEC_CONT vtkm::Float64 Ldexp(vtkm::Float64 x, vtkm::Int32 exponent)
2592 {
2593 #ifdef VTKM_CUDA
2594 return VTKM_CUDA_MATH_FUNCTION_64(ldexp)(x, exponent);
2595 #else
2596 return std::ldexp(x, exponent);
2597 #endif
2598 }
2599
2600 // See: https://randomascii.wordpress.com/2012/01/23/stupid-float-tricks-2/ for why this works.
2601 inline VTKM_EXEC_CONT vtkm::UInt64 FloatDistance(vtkm::Float64 x, vtkm::Float64 y)
2602 {
2603 static_assert(sizeof(vtkm::Float64) == sizeof(vtkm::UInt64), "vtkm::Float64 is incorrect size.");
2604 static_assert(std::numeric_limits<vtkm::Float64>::has_denorm == std::denorm_present, "FloatDistance presumes the floating-point type has subnormal numbers.");
2605
2606 if (!vtkm::IsFinite(x) || !vtkm::IsFinite(y)) {
2607 return 0xFFFFFFFFFFFFFFFFL;
2608 }
2609
2610 // Signed zero is the sworn enemy of this process.
2611 if (y == 0) {
2612 y = vtkm::Abs(y);
2613 }
2614 if (x == 0) {
2615 x = vtkm::Abs(x);
2616 }
2617
2618 if ( (x < 0 && y >= 0) || (x >= 0 && y < 0) )
2619 {
2620 vtkm::UInt64 dx, dy;
2621 if (x < 0) {
2622 dy = FloatDistance(0.0, y);
2623 dx = FloatDistance(0.0, -x);
2624 }
2625 else {
2626 dy = FloatDistance(0.0, -y);
2627 dx = FloatDistance(0.0, x);
2628 }
2629
2630 return dx + dy;
2631 }
2632
2633 if (x < 0 && y < 0) {
2634 return FloatDistance(-x, -y);
2635 }
2636
2637 // Note that:
2638 // int64_t xi = *reinterpret_cast<int64_t*>(&x);
2639 // int64_t yi = *reinterpret_cast<int64_t*>(&y);
2640 // also works, but generates warnings.
2641 // Good option to have if we get compile errors off memcpy or don't want to #include <cstring> though.
2642 // At least on gcc, both versions generate the same assembly.
2643 vtkm::UInt64 xi;
2644 vtkm::UInt64 yi;
2645 memcpy(&xi, &x, sizeof(vtkm::UInt64));
2646 memcpy(&yi, &y, sizeof(vtkm::UInt64));
2647 if (yi > xi) {
2648 return yi - xi;
2649 }
2650 return xi - yi;
2651 }
2652
2653 inline VTKM_EXEC_CONT vtkm::UInt64 FloatDistance(vtkm::Float32 x, vtkm::Float32 y)
2654 {
2655 static_assert(sizeof(vtkm::Float32) == sizeof(vtkm::Int32), "vtkm::Float32 is incorrect size.");
2656 static_assert(std::numeric_limits<vtkm::Float32>::has_denorm == std::denorm_present, "FloatDistance presumes the floating-point type has subnormal numbers.");
2657
2658 if (!vtkm::IsFinite(x) || !vtkm::IsFinite(y)) {
2659 return 0xFFFFFFFFFFFFFFFFL;
2660 }
2661
2662 if (y == 0) {
2663 y = vtkm::Abs(y);
2664 }
2665 if (x == 0) {
2666 x = vtkm::Abs(x);
2667 }
2668
2669 if ( (x < 0 && y >= 0) || (x >= 0 && y < 0) )
2670 {
2671 vtkm::UInt64 dx, dy;
2672 if (x < 0) {
2673 dy = FloatDistance(0.0f, y);
2674 dx = FloatDistance(0.0f, -x);
2675 }
2676 else {
2677 dy = FloatDistance(0.0f, -y);
2678 dx = FloatDistance(0.0f, x);
2679 }
2680 return dx + dy;
2681 }
2682
2683 if (x < 0 && y < 0) {
2684 return FloatDistance(-x, -y);
2685 }
2686
2687 vtkm::UInt32 xi_32;
2688 vtkm::UInt32 yi_32;
2689 memcpy(&xi_32, &x, sizeof(vtkm::UInt32));
2690 memcpy(&yi_32, &y, sizeof(vtkm::UInt32));
2691 vtkm::UInt64 xi = xi_32;
2692 vtkm::UInt64 yi = yi_32;
2693 if (yi > xi) {
2694 return yi - xi;
2695 }
2696 return xi - yi;
2697 }
2698
2699 // Computes ab - cd.
2700 // See: https://pharr.org/matt/blog/2019/11/03/difference-of-floats.html
2701 template<typename T>
2702 inline VTKM_EXEC_CONT T DifferenceOfProducts(T a, T b, T c, T d)
2703 {
2704 T cd = c * d;
2705 T err = std::fma(-c, d, cd);
2706 T dop = std::fma(a, b, -cd);
2707 return dop + err;
2708 }
2709
2710 // Solves ax² + bx + c = 0.
2711 // Only returns the real roots.
2712 // If there are real roots, the first element of the pair is <= the second.
2713 // If there are no real roots, both elements are NaNs.
2714 // The error should be at most 3 ulps.
2715 template<typename T>
2716 inline VTKM_EXEC_CONT vtkm::Vec<T, 2> QuadraticRoots(T a, T b, T c)
2717 {
2718 if (a == 0)
2719 {
2720 if (b == 0)
2721 {
2722 if (c == 0)
2723 {
2724 // A degenerate case. All real numbers are roots; hopefully this arbitrary decision interacts gracefully with use.
2725 return vtkm::Vec<T,2>(0,0);
2726 }
2727 else
2728 {
2729 return vtkm::Vec<T,2>(vtkm::Nan<T>(), vtkm::Nan<T>());
2730 }
2731 }
2732 return vtkm::Vec<T,2>(-c/b, -c/b);
2733 }
2734 T delta = DifferenceOfProducts(b, b, 4*a, c);
2735 if (delta < 0)
2736 {
2737 return vtkm::Vec<T,2>(vtkm::Nan<T>(), vtkm::Nan<T>());
2738 }
2739
2740 T q = -(b + vtkm::CopySign(vtkm::Sqrt(delta), b)) / 2;
2741 T r0 = q / a;
2742 T r1 = c / q;
2743 if (r0 < r1)
2744 {
2745 return vtkm::Vec<T,2>(r0, r1);
2746 }
2747 return vtkm::Vec<T,2>(r1, r0);
2748 }
2749
2750 /// Bitwise operations
2751 ///
2752
2753 /// Find the first set bit in @a word, and return its position (1-32). If no
2754 /// bits are set, returns 0.
2755 #ifdef VTKM_CUDA_DEVICE_PASS
2756 // Need to explicitly mark this as __device__ since __ffs is device only.
2757 inline __device__
2758 vtkm::Int32 FindFirstSetBit(vtkm::UInt32 word)
2759 {
2760 // Output is [0,32], with ffs(0) == 0
2761 return __ffs(static_cast<int>(word));
2762 }
2763 #else // CUDA_DEVICE_PASS
2764 inline VTKM_EXEC_CONT
2765 vtkm::Int32 FindFirstSetBit(vtkm::UInt32 word)
2766 {
2767 # if defined(VTKM_GCC) || defined(VTKM_CLANG)
2768
2769 // Output is [0,32], with ffs(0) == 0
2770 return __builtin_ffs(static_cast<int>(word));
2771
2772 # elif defined(VTKM_MSVC)
2773
2774 // Output is [0, 31], check return code to see if bits are set:
2775 vtkm::UInt32 firstSet;
2776 return _BitScanForward(reinterpret_cast<DWORD*>(&firstSet), word) != 0
2777 ? static_cast<vtkm::Int32>(firstSet + 1) : 0;
2778
2779 # elif defined(VTKM_ICC)
2780
2781 // Output is [0, 31], undefined if word is 0.
2782 return word != 0 ? _bit_scan_forward(word) + 1 : 0;
2783
2784 # else
2785
2786 // Naive implementation:
2787 if (word == 0)
2788 {
2789 return 0;
2790 }
2791
2792 vtkm::Int32 bit = 1;
2793 while ((word & 0x1) == 0)
2794 {
2795 word >>= 1;
2796 ++bit;
2797 }
2798 return bit;
2799
2800 # endif
2801 }
2802 #endif // CUDA_DEVICE_PASS
2803
2804 /// Find the first set bit in @a word, and return its position (1-64). If no
2805 /// bits are set, returns 0.
2806 #ifdef VTKM_CUDA_DEVICE_PASS
2807 // Need to explicitly mark this as __device__ since __ffsll is device only.
2808 inline __device__
2809 vtkm::Int32 FindFirstSetBit(vtkm::UInt64 word)
2810 {
2811
2812 // Output is [0,64], with ffs(0) == 0
2813 return __ffsll(static_cast<long long int>(word));
2814 }
2815 #else // CUDA_DEVICE_PASS
2816 inline VTKM_EXEC_CONT
2817 vtkm::Int32 FindFirstSetBit(vtkm::UInt64 word)
2818 {
2819 # if defined(VTKM_GCC) || defined(VTKM_CLANG) || defined(VTKM_ICC)
2820
2821 // Output is [0,64], with ffs(0) == 0
2822 return __builtin_ffsll(static_cast<long long int>(word));
2823
2824 # elif defined(VTKM_MSVC)
2825
2826 // Output is [0, 63], check return code to see if bits are set:
2827 vtkm::UInt32 firstSet;
2828 return _BitScanForward64(reinterpret_cast<DWORD*>(&firstSet), word) != 0
2829 ? static_cast<vtkm::Int32>(firstSet + 1) : 0;
2830
2831 # else
2832
2833 // Naive implementation:
2834 if (word == 0)
2835 {
2836 return 0;
2837 }
2838
2839 vtkm::Int32 bit = 1;
2840 while ((word & 0x1) == 0)
2841 {
2842 word >>= 1;
2843 ++bit;
2844 }
2845 return bit;
2846
2847 # endif
2848 }
2849 #endif // CUDA_DEVICE_PASS
2850
2851 /// Count the total number of bits set in @a word.
2852 #ifdef VTKM_CUDA_DEVICE_PASS
2853 // Need to explicitly mark this as __device__ since __popc is device only.
2854 inline __device__
2855 vtkm::Int32 CountSetBits(vtkm::UInt32 word)
2856 {
2857 return __popc(word);
2858 }
2859 #else // CUDA_DEVICE_PASS
2860 inline VTKM_EXEC_CONT
2861 vtkm::Int32 CountSetBits(vtkm::UInt32 word)
2862 {
2863 # if defined(VTKM_GCC) || defined(VTKM_CLANG)
2864
2865 return __builtin_popcount(word);
2866
2867 # elif defined(VTKM_MSVC)
2868
2869 return static_cast<vtkm::Int32>(__popcnt(word));
2870
2871 # elif defined(VTKM_ICC)
2872
2873 return _popcnt32(static_cast<int>(word));
2874
2875 # else
2876
2877 // Naive implementation:
2878 vtkm::Int32 bits = 0;
2879 while (word)
2880 {
2881 if (word & 0x1)
2882 {
2883 ++bits;
2884 }
2885 word >>= 1;
2886 }
2887 return bits;
2888
2889 # endif
2890 }
2891 #endif // CUDA_DEVICE_PASS
2892
2893 /// Count the total number of bits set in @a word.
2894 #ifdef VTKM_CUDA_DEVICE_PASS
2895 // Need to explicitly mark this as __device__ since __popcll is device only.
2896 inline __device__
2897 vtkm::Int32 CountSetBits(vtkm::UInt64 word)
2898 {
2899 return __popcll(word);
2900 }
2901 #else // CUDA_DEVICE_PASS
2902 inline VTKM_EXEC_CONT
2903 vtkm::Int32 CountSetBits(vtkm::UInt64 word)
2904 {
2905 # if defined(VTKM_GCC) || defined(VTKM_CLANG)
2906
2907 return __builtin_popcountll(word);
2908
2909 # elif defined(VTKM_MSVC)
2910
2911 return static_cast<vtkm::Int32>(__popcnt64(word));
2912
2913 # elif defined(VTKM_ICC)
2914
2915 return _popcnt64(static_cast<vtkm::Int64>(word));
2916
2917 # else
2918
2919 // Naive implementation:
2920 vtkm::Int32 bits = 0;
2921 while (word)
2922 {
2923 if (word & 0x1)
2924 {
2925 ++bits;
2926 }
2927 word >>= 1;
2928 }
2929 return bits;
2930
2931 # endif
2932 }
2933 #endif // CUDA_DEVICE_PASS
2934
2935 } // namespace vtkm
2936 // clang-format on
2937
2938 #endif //vtk_m_Math_h
2939