Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
cuda_atomic.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2016 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_cuda_atomic_h
17 #define dealii_cuda_atomic_h
18 
19 #include <deal.II/base/config.h>
20 
21 #ifdef DEAL_II_COMPILER_CUDA_AWARE
22 
23 DEAL_II_NAMESPACE_OPEN
24 
25 namespace LinearAlgebra
26 {
27  namespace CUDAWrappers
28  {
34  inline __device__ float
35  atomicAdd_wrapper(float *address, float val)
36  {
37  return atomicAdd(address, val);
38  }
39 
40 
41 
47  inline __device__ double
48  atomicAdd_wrapper(double *address, double val)
49  {
50  // Use native instruction for CUDA 8 on Pascal or newer architecture
51 # if __CUDACC_VER_MAJOR__ >= 8 && \
52  (!defined(__CUDA_ARCH__) || __CUDA_ARCH__ >= 600)
53  return atomicAdd(address, val);
54 # else
55 
56  unsigned long long int *address_as_ull =
57  reinterpret_cast<unsigned long long int *>(address);
58  unsigned long long int old = *address_as_ull, assumed;
59  do
60  {
61  assumed = old;
62  old = atomicCAS(address_as_ull,
63  assumed,
64  __double_as_longlong(val +
65  __longlong_as_double(assumed)));
66  }
67  while (assumed != old);
68 
69  return __longlong_as_double(old);
70 # endif
71  }
72 
73 
74 
80  inline __device__ float
81  atomicMax_wrapper(float *address, float val)
82  {
83  int *address_as_int = reinterpret_cast<int *>(address);
84  int old = *address_as_int, assumed;
85  do
86  {
87  assumed = old;
88  old = atomicCAS(address_as_int,
89  assumed,
90  atomicMax(address_as_int, __float_as_int(val)));
91  }
92  while (assumed != old);
93 
94  return __longlong_as_double(old);
95  }
96 
97 
98 
104  inline __device__ double
105  atomicMax_wrapper(double *address, double val)
106  {
107  unsigned long long int *address_as_ull =
108  reinterpret_cast<unsigned long long int *>(address);
109  unsigned long long int old = *address_as_ull, assumed;
110  do
111  {
112  assumed = old;
113  old = atomicCAS(address_as_ull,
114  assumed,
115  atomicMax(address_as_ull,
116  static_cast<unsigned long long int>(
117  __double_as_longlong(val))));
118  }
119  while (assumed != old);
120 
121  return __longlong_as_double(old);
122  }
123  } // namespace CUDAWrappers
124 } // namespace LinearAlgebra
125 
126 DEAL_II_NAMESPACE_CLOSE
127 
128 #endif
129 
130 #endif
float atomicMax_wrapper(float *address, float val)
Definition: cuda_atomic.h:81
float atomicAdd_wrapper(float *address, float val)
Definition: cuda_atomic.h:35