cpu_backend.h - sphere - GPU-based 3D discrete element method algorithm with optional fluid coupling
HTML git clone git://src.adamsgaard.dk/sphere
DIR Log
DIR Files
DIR Refs
DIR LICENSE
---
cpu_backend.h (6507B)
---
1 // cpu_backend.h -- CUDA execution-model emulation for the OpenMP CPU backend.
2 // Supplies the vector types, qualifier macros, thread-index globals and CUDA
3 // runtime API shims needed to compile the kernel sources as plain C++.
4 #ifndef CPU_BACKEND_H_
5 #define CPU_BACKEND_H_
6
7 #ifdef SPHERE_GPU
8 #error "cpu_backend.h must not be included in CUDA builds"
9 #endif
10
11 #include <cstdlib>
12 #include <cstring>
13 #include <cmath>
14 #include <chrono>
15
16 //// Vector types normally provided by the CUDA toolkit headers
17
18 struct double2 { double x, y; };
19 struct double3 { double x, y, z; };
20 struct double4 { double x, y, z, w; };
21 struct float3 { float x, y, z; };
22 struct float4 { float x, y, z, w; };
23 struct int3 { int x, y, z; };
24 struct uint2 { unsigned int x, y; };
25 struct uint3 { unsigned int x, y, z; };
26
27 struct dim3 {
28 unsigned int x, y, z;
29 dim3(unsigned int x_ = 1, unsigned int y_ = 1, unsigned int z_ = 1)
30 : x(x_), y(y_), z(z_) {}
31 };
32
33 inline double2 make_double2(double x, double y)
34 {
35 double2 v; v.x = x; v.y = y; return v;
36 }
37
38 inline double3 make_double3(double x, double y, double z)
39 {
40 double3 v; v.x = x; v.y = y; v.z = z; return v;
41 }
42
43 inline double4 make_double4(double x, double y, double z, double w)
44 {
45 double4 v; v.x = x; v.y = y; v.z = z; v.w = w; return v;
46 }
47
48 inline float3 make_float3(float x, float y, float z)
49 {
50 float3 v; v.x = x; v.y = y; v.z = z; return v;
51 }
52
53 inline float4 make_float4(float x, float y, float z, float w)
54 {
55 float4 v; v.x = x; v.y = y; v.z = z; v.w = w; return v;
56 }
57
58 inline int3 make_int3(int x, int y, int z)
59 {
60 int3 v; v.x = x; v.y = y; v.z = z; return v;
61 }
62
63 inline uint2 make_uint2(unsigned int x, unsigned int y)
64 {
65 uint2 v; v.x = x; v.y = y; return v;
66 }
67
68 inline uint3 make_uint3(unsigned int x, unsigned int y, unsigned int z)
69 {
70 uint3 v; v.x = x; v.y = y; v.z = z; return v;
71 }
72
73 // Subset of helper_math.h operators used by the kernels
74 inline int3 operator+(int3 a, int3 b)
75 {
76 return make_int3(a.x + b.x, a.y + b.y, a.z + b.z);
77 }
78
79 // float3/float4 operators used by the raytracer
80 inline float4 make_float4(float3 v, float w)
81 {
82 return make_float4(v.x, v.y, v.z, w);
83 }
84
85 inline float3 operator-(const float3& a)
86 {
87 return make_float3(-a.x, -a.y, -a.z);
88 }
89
90 inline float3 operator+(float3 a, float3 b)
91 {
92 return make_float3(a.x + b.x, a.y + b.y, a.z + b.z);
93 }
94
95 inline float3 operator-(float3 a, float3 b)
96 {
97 return make_float3(a.x - b.x, a.y - b.y, a.z - b.z);
98 }
99
100 inline float3 operator*(float3 a, float3 b)
101 {
102 return make_float3(a.x*b.x, a.y*b.y, a.z*b.z);
103 }
104
105 inline float3 operator*(float3 a, float b)
106 {
107 return make_float3(a.x*b, a.y*b, a.z*b);
108 }
109
110 inline float3 operator*(float b, float3 a)
111 {
112 return make_float3(b*a.x, b*a.y, b*a.z);
113 }
114
115 inline float3 operator/(float3 a, float b)
116 {
117 return make_float3(a.x/b, a.y/b, a.z/b);
118 }
119
120 inline float dot(float3 a, float3 b)
121 {
122 return a.x*b.x + a.y*b.y + a.z*b.z;
123 }
124
125 inline float3 cross(float3 a, float3 b)
126 {
127 return make_float3(a.y*b.z - a.z*b.y,
128 a.z*b.x - a.x*b.z,
129 a.x*b.y - a.y*b.x);
130 }
131
132 inline float length(float3 v)
133 {
134 return sqrtf(dot(v, v));
135 }
136
137 inline float3 normalize(float3 v)
138 {
139 return v * (1.0f/sqrtf(dot(v, v)));
140 }
141
142 //// Qualifier macros and thread-index emulation
143
144 #define __global__
145 #define __device__
146 #define __host__
147 #define __constant__
148 #define __forceinline__ inline
149 #define __inline__ inline
150 #define __shared__ static
151
152 extern thread_local uint3 threadIdx;
153 extern thread_local uint3 blockIdx;
154 extern dim3 blockDim; // uniform per launch
155 extern dim3 gridDim;
156
157 // Safe as a no-op: the only kernels relying on intra-block cooperation
158 // (reorderArrays, summation) have dedicated CPU implementations.
159 inline void __syncthreads() {}
160
161 //// CUDA runtime API shims
162
163 enum cudaMemcpyKind {
164 cudaMemcpyHostToHost,
165 cudaMemcpyHostToDevice,
166 cudaMemcpyDeviceToHost,
167 cudaMemcpyDeviceToDevice
168 };
169
170 typedef int cudaError_t;
171 const cudaError_t cudaSuccess = 0;
172
173 // calloc, not malloc: some kernels (e.g. the Darcy solvers reading
174 // dev_darcy_div_v_p) consume arrays no kernel ever writes, and rely on
175 // fresh CUDA allocations being zeroed pages.
176 inline cudaError_t cudaMalloc(void** p, size_t n)
177 {
178 *p = calloc(1, n);
179 return cudaSuccess;
180 }
181
182 inline cudaError_t cudaFree(void* p)
183 {
184 free(p);
185 return cudaSuccess;
186 }
187
188 inline cudaError_t cudaMemcpy(void* dst, const void* src, size_t n,
189 cudaMemcpyKind)
190 {
191 memcpy(dst, src, n);
192 return cudaSuccess;
193 }
194
195 inline cudaError_t cudaMemset(void* p, int v, size_t n)
196 {
197 memset(p, v, n);
198 return cudaSuccess;
199 }
200
201 inline cudaError_t cudaDeviceSynchronize() { return cudaSuccess; }
202 inline cudaError_t cudaThreadSynchronize() { return cudaSuccess; }
203 inline cudaError_t cudaGetLastError() { return cudaSuccess; }
204 inline cudaError_t cudaSetDevice(int) { return cudaSuccess; }
205 inline cudaError_t cudaDeviceReset() { return cudaSuccess; }
206
207 #define cudaMemcpyToSymbol(sym, src, ...) memcpy(&(sym), (src), sizeof(sym))
208 #define cudaMemcpyFromSymbol(dst, sym, ...) memcpy((dst), &(sym), sizeof(sym))
209
210 // Event shims for the PROFILING timers
211 typedef std::chrono::high_resolution_clock::time_point cudaEvent_t;
212
213 inline void cudaEventCreate(cudaEvent_t*) {}
214
215 inline void cudaEventRecord(cudaEvent_t& e, int = 0)
216 {
217 e = std::chrono::high_resolution_clock::now();
218 }
219
220 inline void cudaEventSynchronize(cudaEvent_t&) {}
221
222 inline void cudaEventElapsedTime(float* ms, cudaEvent_t& tic, cudaEvent_t& toc)
223 {
224 *ms = std::chrono::duration<float, std::milli>(toc - tic).count();
225 }
226
227 inline void cudaEventDestroy(cudaEvent_t&) {}
228
229 //// CPU launch driver: one OpenMP-parallel loop over emulated blocks*threads
230
231 template<typename F>
232 void cpuKernelLaunch(const dim3 grid, const dim3 block, size_t /*smem*/, F f)
233 {
234 gridDim = grid;
235 blockDim = block;
236 const long nb = (long)grid.x * grid.y * grid.z;
237 const long nt = (long)block.x * block.y * block.z;
238 // The if clause keeps tiny launches serial and cheap
239 #pragma omp parallel for schedule(static) collapse(2) if (nb*nt > 512)
240 for (long b = 0; b < nb; ++b) {
241 for (long t = 0; t < nt; ++t) {
242 blockIdx.x = b % grid.x;
243 blockIdx.y = (b / grid.x) % grid.y;
244 blockIdx.z = b / ((long)grid.x * grid.y);
245 threadIdx.x = t % block.x;
246 threadIdx.y = (t / block.x) % block.y;
247 threadIdx.z = t / ((long)block.x * block.y);
248 f();
249 }
250 }
251 }
252
253 #endif
254 // vim: tabstop=8 expandtab shiftwidth=4 softtabstop=4