Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Stokhos_DynArrayTraits.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Stokhos Package
4 //
5 // Copyright 2009 NTESS and the Stokhos contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef STOKHOS_DYN_ARRAY_TRAITS_HPP
11 #define STOKHOS_DYN_ARRAY_TRAITS_HPP
12 
13 #include <new>
14 #include <cstring>
15 
16 #include "Kokkos_Core_fwd.hpp"
17 
18 namespace Stokhos {
19 
21 
25  template <typename T> struct IsScalarType2 {
26  static const bool value = false;
27  };
28 
30 #define STOKHOS_BUILTIN_SPECIALIZATION(t) \
31  template <> struct IsScalarType2< t > { \
32  static const bool value = true; \
33  };
34 
39 
40 #undef STOKHOS_BUILTIN_SPECIALIZATION
41 
46  template <typename T, typename device_t,
47  bool isScalar = IsScalarType2<T>::value>
48  struct DynArrayTraits {
49 
50  typedef T value_type;
51  typedef device_t execution_space;
52 
54  static
55  KOKKOS_INLINE_FUNCTION
56  void copy(const volatile T* src, volatile T* dest, std::size_t sz) {
57  //if (sz > 0) std::memcpy(dest,src,sz*sizeof(T));
58  for (std::size_t i=0; i<sz; ++i)
59  *(dest++) = *(src++);
60  }
61 
63  static
64  KOKKOS_INLINE_FUNCTION
65  void copy(const volatile T* src, T* dest, std::size_t sz) {
66  //if (sz > 0) std::memcpy(dest,src,sz*sizeof(T));
67  for (std::size_t i=0; i<sz; ++i)
68  *(dest++) = *(src++);
69  }
70 
72  static
73  KOKKOS_INLINE_FUNCTION
74  void copy(const T* src, volatile T* dest, std::size_t sz) {
75  //if (sz > 0) std::memcpy(dest,src,sz*sizeof(T));
76  for (std::size_t i=0; i<sz; ++i)
77  *(dest++) = *(src++);
78  }
79 
81  static
82  KOKKOS_INLINE_FUNCTION
83  void copy(const T* src, T* dest, std::size_t sz) {
84  //if (sz > 0) std::memcpy(dest,src,sz*sizeof(T));
85  for (std::size_t i=0; i<sz; ++i)
86  *(dest++) = *(src++);
87  }
88 
90  static
91  KOKKOS_INLINE_FUNCTION
92  void zero(T* dest, std::size_t sz) {
93  if (sz > 0) std::memset(dest,0,sz*sizeof(T));
94  }
95 
97  static
98  KOKKOS_INLINE_FUNCTION
99  void zero(volatile T* dest, std::size_t sz) {
100  //if (sz > 0) std::memset(dest,0,sz*sizeof(T));
101  for (std::size_t i=0; i<sz; ++i)
102  *(dest++) = T(0.);
103  }
104 
106  static
107  KOKKOS_INLINE_FUNCTION
108  void fill(T* dest, std::size_t sz, const T& v) {
109  //if (sz > 0) std::memset(dest,v,sz*sizeof(T));
110  for (std::size_t i=0; i<sz; ++i)
111  *(dest++) = v;
112  }
113 
115  static
116  KOKKOS_INLINE_FUNCTION
117  void fill(volatile T* dest, std::size_t sz, const T& v) {
118  //if (sz > 0) std::memset(dest,v,sz*sizeof(T));
119  for (std::size_t i=0; i<sz; ++i)
120  *(dest++) = v;
121  }
122 
124  static
125  KOKKOS_INLINE_FUNCTION
126  T* get_and_fill(std::size_t sz, const T& x = T(0.0)) {
127  T* m = 0;
128  if (sz > 0) {
129  m = static_cast<T* >(operator new(sz*sizeof(T)));
130  //std::memset(m,x,sz*sizeof(T));
131  for (std::size_t i=0; i<sz; ++i)
132  m[i] = x;
133  }
134  return m;
135  }
136 
141  static
142  KOKKOS_INLINE_FUNCTION
143  T* get_and_fill(const T* src, std::size_t sz) {
144  T* m = 0;
145  if (sz > 0) {
146  m = static_cast<T* >(operator new(sz*sizeof(T)));
147  for (std::size_t i=0; i<sz; ++i)
148  m[i] = src[i];
149  }
150  return m;
151  }
152 
157  static
158  KOKKOS_INLINE_FUNCTION
159  T* get_and_fill(const volatile T* src, std::size_t sz) {
160  T* m = 0;
161  if (sz > 0) {
162  m = static_cast<T* >(operator new(sz*sizeof(T)));
163  for (std::size_t i=0; i<sz; ++i)
164  m[i] = src[i];
165  }
166  return m;
167  }
168 
170  static
171  KOKKOS_INLINE_FUNCTION
172  void destroy_and_release(T* m, std::size_t sz) {
173  if (sz > 0) operator delete((void*) m);
174  }
175 
177  static
178  KOKKOS_INLINE_FUNCTION
179  void destroy_and_release(volatile T* m, std::size_t sz) {
180  if (sz > 0) operator delete((void*) m);
181  }
182  };
183 
187  template <typename T, typename device_t>
188  struct DynArrayTraits<T, device_t, false> {
189 
190  typedef T value_type;
191  typedef device_t execution_space;
192 
194  static
195  KOKKOS_INLINE_FUNCTION
196  void fill(T* dest, std::size_t sz, const T& v) {
197  for (std::size_t i=0; i<sz; ++i)
198  *(dest++) = v;
199  }
200 
202  static
203  KOKKOS_INLINE_FUNCTION
204  void fill(volatile T* dest, std::size_t sz, const T& v) {
205  for (std::size_t i=0; i<sz; ++i)
206  *(dest++) = v;
207  }
208 
210  static
211  KOKKOS_INLINE_FUNCTION
212  void copy(const volatile T* src, volatile T* dest, std::size_t sz) {
213  for (std::size_t i=0; i<sz; ++i)
214  *(dest++) = *(src++);
215  }
216 
218  static
219  KOKKOS_INLINE_FUNCTION
220  void copy(const volatile T* src, T* dest, std::size_t sz) {
221  for (std::size_t i=0; i<sz; ++i)
222  *(dest++) = *(src++);
223  }
224 
226  static
227  KOKKOS_INLINE_FUNCTION
228  void copy(const T* src, volatile T* dest, std::size_t sz) {
229  for (std::size_t i=0; i<sz; ++i)
230  *(dest++) = *(src++);
231  }
232 
234  static
235  KOKKOS_INLINE_FUNCTION
236  void copy(const T* src, T* dest, std::size_t sz) {
237  for (std::size_t i=0; i<sz; ++i)
238  *(dest++) = *(src++);
239  }
240 
242  static
243  KOKKOS_INLINE_FUNCTION
244  void zero(T* dest, std::size_t sz) {
245  for (std::size_t i=0; i<sz; ++i)
246  *(dest++) = T(0.);
247  }
248 
250  static
251  KOKKOS_INLINE_FUNCTION
252  void zero(volatile T* dest, std::size_t sz) {
253  for (std::size_t i=0; i<sz; ++i)
254  *(dest++) = T(0.);
255  }
256 
258  static
259  KOKKOS_INLINE_FUNCTION
260  T* get_and_fill(std::size_t sz, const T& x = T(0.0)) {
261  T* m = 0;
262  if (sz > 0) {
263  m = static_cast<T* >(operator new(sz*sizeof(T)));
264  T* p = m;
265  for (std::size_t i=0; i<sz; ++i)
266  new (p++) T(x);
267  }
268  return m;
269  }
270 
275  static
276  KOKKOS_INLINE_FUNCTION
277  T* get_and_fill(const T* src, std::size_t sz) {
278  T* m = 0;
279  if (sz > 0) {
280  m = static_cast<T* >(operator new(sz*sizeof(T)));
281  T* p = m;
282  for (std::size_t i=0; i<sz; ++i)
283  new (p++) T(*(src++));
284  }
285  return m;
286  }
287 
292  static
293  KOKKOS_INLINE_FUNCTION
294  T* get_and_fill(const volatile T* src, std::size_t sz) {
295  T* m = 0;
296  if (sz > 0) {
297  m = static_cast<T* >(operator new(sz*sizeof(T)));
298  T* p = m;
299  for (std::size_t i=0; i<sz; ++i)
300  new (p++) T(*(src++));
301  }
302  return m;
303  }
304 
306  static
307  KOKKOS_INLINE_FUNCTION
308  void destroy_and_release(T* m, std::size_t sz) {
309  T* e = m+sz;
310  for (T* b = m; b!=e; b++)
311  b->~T();
312  operator delete((void*) m);
313  }
314 
316  static
317  KOKKOS_INLINE_FUNCTION
318  void destroy_and_release(volatile T* m, std::size_t sz) {
319  T* e = m+sz;
320  for (T* b = m; b!=e; b++)
321  b->~T();
322  operator delete((void*) m);
323  }
324  };
325 
326 #if defined(KOKKOS_ENABLE_CUDA)
327 
334  template <typename T>
335  struct DynArrayTraits<T,Kokkos::Cuda,true> {
336 
337  typedef T value_type;
338  typedef Kokkos::Cuda execution_space;
339 
341  static
342  KOKKOS_INLINE_FUNCTION
343  void copy(const volatile T* src, volatile T* dest, std::size_t sz) {
344  //if (sz > 0) std::memcpy(dest,src,sz*sizeof(T));
345  for (std::size_t i=0; i<sz; ++i)
346  *(dest++) = *(src++);
347  }
348 
350  static
351  KOKKOS_INLINE_FUNCTION
352  void copy(const volatile T* src, T* dest, std::size_t sz) {
353  //if (sz > 0) std::memcpy(dest,src,sz*sizeof(T));
354  for (std::size_t i=0; i<sz; ++i)
355  *(dest++) = *(src++);
356  }
357 
359  static
360  KOKKOS_INLINE_FUNCTION
361  void copy(const T* src, volatile T* dest, std::size_t sz) {
362  //if (sz > 0) std::memcpy(dest,src,sz*sizeof(T));
363  for (std::size_t i=0; i<sz; ++i)
364  *(dest++) = *(src++);
365  }
366 
368  static
369  KOKKOS_INLINE_FUNCTION
370  void copy(const T* src, T* dest, std::size_t sz) {
371  //if (sz > 0) std::memcpy(dest,src,sz*sizeof(T));
372  for (std::size_t i=0; i<sz; ++i)
373  *(dest++) = *(src++);
374  }
375 
377  static
378  KOKKOS_INLINE_FUNCTION
379  void zero(T* dest, std::size_t sz) {
380  if (sz > 0) std::memset(dest,0,sz*sizeof(T));
381  }
382 
384  static
385  KOKKOS_INLINE_FUNCTION
386  void zero(volatile T* dest, std::size_t sz) {
387  //if (sz > 0) std::memset(dest,0,sz*sizeof(T));
388  for (std::size_t i=0; i<sz; ++i)
389  *(dest++) = T(0.);
390  }
391 
393  static
394  KOKKOS_INLINE_FUNCTION
395  void fill(T* dest, std::size_t sz, const T& v) {
396  //if (sz > 0) std::memset(dest,v,sz*sizeof(T));
397  for (std::size_t i=0; i<sz; ++i)
398  *(dest++) = v;
399  }
400 
402  static
403  KOKKOS_INLINE_FUNCTION
404  void fill(volatile T* dest, std::size_t sz, const T& v) {
405  //if (sz > 0) std::memset(dest,v,sz*sizeof(T));
406  for (std::size_t i=0; i<sz; ++i)
407  *(dest++) = v;
408  }
409 
411  static
412  KOKKOS_INLINE_FUNCTION
413  T* get_and_fill(std::size_t sz, const T& x = T(0.0)) {
414  T* m = 0;
415  if (sz > 0) {
416 #if defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_ENABLE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
417  cudaMallocManaged( (void**) &m, sz*sizeof(T), cudaMemAttachGlobal );
418 #else
419  m = static_cast<T* >(operator new(sz*sizeof(T)));
420 #endif
421  //std::memset(m,x,sz*sizeof(T));
422  for (std::size_t i=0; i<sz; ++i)
423  m[i] = x;
424  }
425  return m;
426  }
427 
432  static
433  KOKKOS_INLINE_FUNCTION
434  T* get_and_fill(const T* src, std::size_t sz) {
435  T* m = 0;
436  if (sz > 0) {
437 #if defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_ENABLE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
438  cudaMallocManaged( (void**) &m, sz*sizeof(T), cudaMemAttachGlobal );
439 #else
440  m = static_cast<T* >(operator new(sz*sizeof(T)));
441 #endif
442  for (std::size_t i=0; i<sz; ++i)
443  m[i] = src[i];
444  }
445  return m;
446  }
447 
452  static
453  KOKKOS_INLINE_FUNCTION
454  T* get_and_fill(const volatile T* src, std::size_t sz) {
455  T* m = 0;
456  if (sz > 0) {
457 #if defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_ENABLE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
458  cudaMallocManaged( (void**) &m, sz*sizeof(T), cudaMemAttachGlobal );
459 #else
460  m = static_cast<T* >(operator new(sz*sizeof(T)));
461 #endif
462  for (std::size_t i=0; i<sz; ++i)
463  m[i] = src[i];
464  }
465  return m;
466  }
467 
469  static
470  KOKKOS_INLINE_FUNCTION
471  void destroy_and_release(T* m, std::size_t sz) {
472 #if defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_ENABLE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
473  cudaFree(m);
474 #else
475  if (sz > 0) operator delete((void*) m);
476 #endif
477  }
478 
480  static
481  KOKKOS_INLINE_FUNCTION
482  void destroy_and_release(volatile T* m, std::size_t sz) {
483 #if defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_ENABLE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
484  cudaFree(m);
485 #else
486  if (sz > 0) operator delete((void*) m);
487 #endif
488  }
489  };
490 
496  template <typename T>
497  struct DynArrayTraits<T, Kokkos::Cuda, false> {
498 
499  typedef T value_type;
500  typedef Kokkos::Cuda execution_space;
501 
503  static
504  KOKKOS_INLINE_FUNCTION
505  void fill(T* dest, std::size_t sz, const T& v) {
506  for (std::size_t i=0; i<sz; ++i)
507  *(dest++) = v;
508  }
509 
511  static
512  KOKKOS_INLINE_FUNCTION
513  void fill(volatile T* dest, std::size_t sz, const T& v) {
514  for (std::size_t i=0; i<sz; ++i)
515  *(dest++) = v;
516  }
517 
519  static
520  KOKKOS_INLINE_FUNCTION
521  void copy(const volatile T* src, volatile T* dest, std::size_t sz) {
522  for (std::size_t i=0; i<sz; ++i)
523  *(dest++) = *(src++);
524  }
525 
527  static
528  KOKKOS_INLINE_FUNCTION
529  void copy(const volatile T* src, T* dest, std::size_t sz) {
530  for (std::size_t i=0; i<sz; ++i)
531  *(dest++) = *(src++);
532  }
533 
535  static
536  KOKKOS_INLINE_FUNCTION
537  void copy(const T* src, volatile T* dest, std::size_t sz) {
538  for (std::size_t i=0; i<sz; ++i)
539  *(dest++) = *(src++);
540  }
541 
543  static
544  KOKKOS_INLINE_FUNCTION
545  void copy(const T* src, T* dest, std::size_t sz) {
546  for (std::size_t i=0; i<sz; ++i)
547  *(dest++) = *(src++);
548  }
549 
551  static
552  KOKKOS_INLINE_FUNCTION
553  void zero(T* dest, std::size_t sz) {
554  for (std::size_t i=0; i<sz; ++i)
555  *(dest++) = T(0.);
556  }
557 
559  static
560  KOKKOS_INLINE_FUNCTION
561  void zero(volatile T* dest, std::size_t sz) {
562  for (std::size_t i=0; i<sz; ++i)
563  *(dest++) = T(0.);
564  }
565 
567  static
568  KOKKOS_INLINE_FUNCTION
569  T* get_and_fill(std::size_t sz, const T& x = T(0.0)) {
570  T* m = 0;
571  if (sz > 0) {
572 #if defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_ENABLE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
573  cudaMallocManaged( (void**) &m, sz*sizeof(T), cudaMemAttachGlobal );
574 #else
575  m = static_cast<T* >(operator new(sz*sizeof(T)));
576 #endif
577  T* p = m;
578  for (std::size_t i=0; i<sz; ++i)
579  new (p++) T(x);
580  }
581  return m;
582  }
583 
588  static
589  KOKKOS_INLINE_FUNCTION
590  T* get_and_fill(const T* src, std::size_t sz) {
591  T* m = 0;
592  if (sz > 0) {
593 #if defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_ENABLE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
594  cudaMallocManaged( (void**) &m, sz*sizeof(T), cudaMemAttachGlobal );
595 #else
596  m = static_cast<T* >(operator new(sz*sizeof(T)));
597 #endif
598  T* p = m;
599  for (std::size_t i=0; i<sz; ++i)
600  new (p++) T(*(src++));
601  }
602  return m;
603  }
604 
609  static
610  KOKKOS_INLINE_FUNCTION
611  T* get_and_fill(const volatile T* src, std::size_t sz) {
612  T* m = 0;
613  if (sz > 0) {
614 #if defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_ENABLE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
615  cudaMallocManaged( (void**) &m, sz*sizeof(T), cudaMemAttachGlobal );
616 #else
617  m = static_cast<T* >(operator new(sz*sizeof(T)));
618 #endif
619  T* p = m;
620  for (std::size_t i=0; i<sz; ++i)
621  new (p++) T(*(src++));
622  }
623  return m;
624  }
625 
627  static
628  KOKKOS_INLINE_FUNCTION
629  void destroy_and_release(T* m, std::size_t sz) {
630  T* e = m+sz;
631  for (T* b = m; b!=e; b++)
632  b->~T();
633  operator delete((void*) m);
634  }
635 
637  static
638  KOKKOS_INLINE_FUNCTION
639  void destroy_and_release(volatile T* m, std::size_t sz) {
640  T* e = m+sz;
641  for (T* b = m; b!=e; b++)
642  b->~T();
643  operator delete((void*) m);
644  }
645  };
646 
647 #endif
648 
649 } // namespace Stokhos
650 
651 #endif // STOKHOS_DYN_ARRAY_TRAITS_HPP
static KOKKOS_INLINE_FUNCTION void copy(const T *src, volatile T *dest, std::size_t sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION T * get_and_fill(const T *src, std::size_t sz)
Get memory for new array of length sz and fill with entries from src.
static KOKKOS_INLINE_FUNCTION void destroy_and_release(T *m, std::size_t sz)
Destroy array elements and release memory.
static KOKKOS_INLINE_FUNCTION void fill(volatile T *dest, std::size_t sz, const T &v)
Fill array dest of length sz with value v.
static KOKKOS_INLINE_FUNCTION T * get_and_fill(const volatile T *src, std::size_t sz)
Get memory for new array of length sz and fill with entries from src.
static KOKKOS_INLINE_FUNCTION void destroy_and_release(volatile T *m, std::size_t sz)
Destroy array elements and release memory.
static KOKKOS_INLINE_FUNCTION void copy(const volatile T *src, volatile T *dest, std::size_t sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION void copy(const T *src, volatile T *dest, std::size_t sz)
Copy array from src to dest of length sz.
Base template specification for IsScalarType.
static KOKKOS_INLINE_FUNCTION T * get_and_fill(std::size_t sz, const T &x=T(0.0))
Get memory for new array of length sz and fill with zeros.
static KOKKOS_INLINE_FUNCTION void zero(T *dest, std::size_t sz)
Zero out array dest of length sz.
static KOKKOS_INLINE_FUNCTION void fill(T *dest, std::size_t sz, const T &v)
Fill array dest of length sz with value v.
static KOKKOS_INLINE_FUNCTION void zero(T *dest, std::size_t sz)
Zero out array dest of length sz.
static KOKKOS_INLINE_FUNCTION void fill(T *dest, std::size_t sz, const T &v)
Fill array dest of length sz with value v.
static KOKKOS_INLINE_FUNCTION T * get_and_fill(std::size_t sz, const T &x=T(0.0))
Get memory for new array of length sz and fill with zeros.
#define STOKHOS_BUILTIN_SPECIALIZATION(t)
Specialization of above classes to built-in types.
static KOKKOS_INLINE_FUNCTION void copy(const volatile T *src, T *dest, std::size_t sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION T * get_and_fill(const volatile T *src, std::size_t sz)
Get memory for new array of length sz and fill with entries from src.
static KOKKOS_INLINE_FUNCTION void copy(const volatile T *src, T *dest, std::size_t sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION void copy(const T *src, T *dest, std::size_t sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION void fill(volatile T *dest, std::size_t sz, const T &v)
Fill array dest of length sz with value v.
static KOKKOS_INLINE_FUNCTION void copy(const volatile T *src, volatile T *dest, std::size_t sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION void copy(const T *src, T *dest, std::size_t sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION void zero(volatile T *dest, std::size_t sz)
Zero out array dest of length sz.
Dynamic array allocation class that is specialized for scalar i.e., fundamental or built-in types (fl...
static KOKKOS_INLINE_FUNCTION void destroy_and_release(volatile T *m, std::size_t sz)
Destroy array elements and release memory.
static KOKKOS_INLINE_FUNCTION T * get_and_fill(const T *src, std::size_t sz)
Get memory for new array of length sz and fill with entries from src.
static KOKKOS_INLINE_FUNCTION void zero(volatile T *dest, std::size_t sz)
Zero out array dest of length sz.
static KOKKOS_INLINE_FUNCTION void destroy_and_release(T *m, std::size_t sz)
Destroy array elements and release memory.