Sacado Package Browser (Single Doxygen Collection)  Version of the Day
Sacado_DynamicArrayTraits.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Sacado Package
5 // Copyright (2006) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // This library is free software; you can redistribute it and/or modify
11 // it under the terms of the GNU Lesser General Public License as
12 // published by the Free Software Foundation; either version 2.1 of the
13 // License, or (at your option) any later version.
14 //
15 // This library is distributed in the hope that it will be useful, but
16 // WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 // Lesser General Public License for more details.
19 //
20 // You should have received a copy of the GNU Lesser General Public
21 // License along with this library; if not, write to the Free Software
22 // Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23 // USA
24 // Questions? Contact David M. Gay (dmgay@sandia.gov) or Eric T. Phipps
25 // (etphipp@sandia.gov).
26 //
27 // ***********************************************************************
28 // @HEADER
29 
30 #ifndef SACADO_DYNAMICARRAYTRAITS_HPP
31 #define SACADO_DYNAMICARRAYTRAITS_HPP
32 
33 #include <new>
34 #include <cstring>
35 
36 #include "Sacado_Traits.hpp"
37 #if defined(HAVE_SACADO_KOKKOSCORE)
38 #include "Kokkos_Core.hpp"
39 #endif
40 
41 namespace Sacado {
42 
46  template <typename T, bool isScalar = IsScalarType<T>::value>
47  struct ds_array {
48 
50  static T* my_alloc(const int sz) {
51 #if defined(__CUDACC__) && defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_USE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
52  T* m;
53  CUDA_SAFE_CALL( cudaMallocManaged( (void**) &m, sz*sizeof(T), cudaMemAttachGlobal ) );
54 #else
55  T* m = static_cast<T* >(operator new(sz*sizeof(T)));
56 #endif
57  return m;
58  }
59 
61  static void my_free(T* m, int sz) {
62  if (sz > 0) {
63 #if defined(__CUDACC__) && defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_USE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
64  CUDA_SAFE_CALL( cudaFree(m) );
65 #else
66  operator delete((void*) m);
67 #endif
68  }
69  }
70 
73  static T* get(int sz) {
74  if (sz > 0) {
75  T* m = my_alloc(sz);
76  T* p = m;
77  for (int i=0; i<sz; ++i)
78  new (p++) T();
79  return m;
80  }
81  return NULL;
82  }
83 
86  static T* get_and_fill(int sz) {
87  if (sz > 0) {
88  T* m = my_alloc(sz);
89  T* p = m;
90  for (int i=0; i<sz; ++i)
91  new (p++) T(0.0);
92  return m;
93  }
94  return NULL;
95  }
96 
102  static T* get_and_fill(const T* src, int sz) {
103  if (sz > 0) {
104  T* m = my_alloc(sz);
105  T* p = m;
106  for (int i=0; i<sz; ++i)
107  new (p++) T(*(src++));
108  return m;
109  }
110  return NULL;
111  }
112 
118  static T* strided_get_and_fill(const T* src, int stride, int sz) {
119  if (sz > 0) {
120  T* m = my_alloc(sz);
121  T* p = m;
122  for (int i=0; i<sz; ++i) {
123  new (p++) T(*(src));
124  src += stride;
125  }
126  return m;
127  }
128  return NULL;
129  }
130 
133  static void copy(const T* src, T* dest, int sz) {
134  for (int i=0; i<sz; ++i)
135  *(dest++) = *(src++);
136  }
137 
140  static void strided_copy(const T* src, int src_stride,
141  T* dest, int dest_stride, int sz) {
142  for (int i=0; i<sz; ++i) {
143  *(dest) = *(src);
144  dest += dest_stride;
145  src += src_stride;
146  }
147  }
148 
151  static void zero(T* dest, int sz) {
152  for (int i=0; i<sz; ++i)
153  *(dest++) = T(0.);
154  }
155 
158  static void strided_zero(T* dest, int stride, int sz) {
159  for (int i=0; i<sz; ++i) {
160  *(dest) = T(0.);
161  dest += stride;
162  }
163  }
164 
167  static void destroy_and_release(T* m, int sz) {
168  T* e = m+sz;
169  for (T* b = m; b!=e; b++)
170  b->~T();
171  my_free(m, sz);
172  }
173  };
174 
179  template <typename T>
180  struct ds_array<T,true> {
181 
183  static T* my_alloc(const int sz) {
184 #if defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_USE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
185  T* m;
186  CUDA_SAFE_CALL( cudaMallocManaged( (void**) &m, sz*sizeof(T), cudaMemAttachGlobal ) );
187 #else
188  T* m = static_cast<T* >(operator new(sz*sizeof(T)));
189 #endif
190  return m;
191  }
192 
194  static void my_free(T* m, int sz) {
195  if (sz > 0) {
196 #if defined( CUDA_VERSION ) && ( 6000 <= CUDA_VERSION ) && defined(KOKKOS_USE_CUDA_UVM) && !defined( __CUDA_ARCH__ )
197  CUDA_SAFE_CALL( cudaFree(m) );
198 #else
199  operator delete((void*) m);
200 #endif
201  }
202  }
203 
206  static T* get(int sz) {
207  if (sz > 0) {
208  T* m = my_alloc(sz);
209  return m;
210  }
211  return NULL;
212  }
213 
216  static T* get_and_fill(int sz) {
217  if (sz > 0) {
218  T* m = my_alloc(sz);
219 #ifdef __CUDACC__
220  for (int i=0; i<sz; ++i)
221  m[i] = 0.0;
222 #else
223  std::memset(m,0,sz*sizeof(T));
224 #endif
225  return m;
226  }
227  return NULL;
228  }
229 
235  static T* get_and_fill(const T* src, int sz) {
236  if (sz > 0) {
237  T* m = my_alloc(sz);
238  for (int i=0; i<sz; ++i)
239  m[i] = src[i];
240  return m;
241  }
242  return NULL;
243  }
244 
250  static T* strided_get_and_fill(const T* src, int stride, int sz) {
251  if (sz > 0) {
252  T* m = my_alloc(sz);
253  for (int i=0; i<sz; ++i)
254  m[i] = src[i*stride];
255  return m;
256  }
257  return NULL;
258  }
259 
262  static void copy(const T* src, T* dest, int sz) {
263  if (sz > 0)
264 #ifdef __CUDACC__
265  for (int i=0; i<sz; ++i)
266  dest[i] = src[i];
267 #else
268  std::memcpy(dest,src,sz*sizeof(T));
269 #endif
270  }
271 
274  static void strided_copy(const T* src, int src_stride,
275  T* dest, int dest_stride, int sz) {
276  for (int i=0; i<sz; ++i) {
277  *(dest) = *(src);
278  dest += dest_stride;
279  src += src_stride;
280  }
281  }
282 
285  static void zero(T* dest, int sz) {
286  if (sz > 0)
287 #ifdef __CUDACC__
288  for (int i=0; i<sz; ++i)
289  dest[i] = T(0.);
290 #else
291  std::memset(dest,0,sz*sizeof(T));
292 #endif
293  }
294 
297  static void strided_zero(T* dest, int stride, int sz) {
298  for (int i=0; i<sz; ++i) {
299  *(dest) = T(0.);
300  dest += stride;
301  }
302  }
303 
306  static void destroy_and_release(T* m, int sz) {
307  my_free(m, sz);
308  }
309  };
310 
311 } // namespace Sacado
312 
313 #endif // SACADO_DYNAMICARRAY_HPP
static KOKKOS_INLINE_FUNCTION T * strided_get_and_fill(const T *src, int stride, int 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, int sz)
Destroy array elements and release memory.
static KOKKOS_INLINE_FUNCTION void copy(const T *src, T *dest, int sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION T * my_alloc(const int sz)
static KOKKOS_INLINE_FUNCTION void my_free(T *m, int sz)
static KOKKOS_INLINE_FUNCTION void strided_zero(T *dest, int stride, int sz)
Zero out array dest of length sz.
expr true
static KOKKOS_INLINE_FUNCTION void destroy_and_release(T *m, int sz)
Destroy array elements and release memory.
static KOKKOS_INLINE_FUNCTION void strided_copy(const T *src, int src_stride, T *dest, int dest_stride, int sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION T * my_alloc(const int sz)
#define KOKKOS_INLINE_FUNCTION
#define T
Definition: Sacado_rad.hpp:573
static KOKKOS_INLINE_FUNCTION void my_free(T *m, int sz)
static KOKKOS_INLINE_FUNCTION T * get_and_fill(int sz)
Get memory for new array of length sz and fill with zeros.
static KOKKOS_INLINE_FUNCTION void zero(T *dest, int sz)
Zero out array dest of length sz.
static KOKKOS_INLINE_FUNCTION void strided_copy(const T *src, int src_stride, T *dest, int dest_stride, int sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION void copy(const T *src, T *dest, int sz)
Copy array from src to dest of length sz.
static KOKKOS_INLINE_FUNCTION T * strided_get_and_fill(const T *src, int stride, int sz)
Get memory for new array of length sz and fill with entries from src.
static KOKKOS_INLINE_FUNCTION void strided_zero(T *dest, int stride, int sz)
Zero out array dest of length sz.
static KOKKOS_INLINE_FUNCTION T * get_and_fill(const T *src, int sz)
Get memory for new array of length sz and fill with entries from src.
static KOKKOS_INLINE_FUNCTION T * get_and_fill(int sz)
Get memory for new array of length sz and fill with zeros.
static KOKKOS_INLINE_FUNCTION T * get_and_fill(const T *src, int sz)
Get memory for new array of length sz and fill with entries from src.
Dynamic array allocation class that works for any type.
static KOKKOS_INLINE_FUNCTION void zero(T *dest, int sz)
Zero out array dest of length sz.