cory-master/BabelStream
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Merge branch 'master' into bare_hc

Browse Source
This commit is contained in:
Peter Steinbach 2017-01-30 16:06:34 +01:00
commit c9a45600c8
23 changed files with 533 additions and 301 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
.gitignore vendored
View File

@ -1,6 +1,4 @@
common.h
gpu-stream-cuda
gpu-stream-ocl
gpu-stream-acc

36
ACCStream.cpp
View File

@ -36,13 +36,19 @@ ACCStream<T>::~ACCStream()
}
template <class T>
void ACCStream<T>::write_arrays(const std::vector<T>& h_a, const std::vector<T>& h_b, const std::vector<T>& h_c)
void ACCStream<T>::init_arrays(T initA, T initB, T initC)
{
T *a = this->a;
T *b = this->b;
T *c = this->c;
#pragma acc update device(a[0:array_size], b[0:array_size], c[0:array_size])
{}
unsigned int array_size = this->array_size;
T * restrict a = this->a;
T * restrict b = this->b;
T * restrict c = this->c;
#pragma acc kernels present(a[0:array_size], b[0:array_size], c[0:array_size]) wait
for (int i = 0; i < array_size; i++)
{
a[i] = initA;
b[i] = initB;
c[i] = initC;
}
}
template <class T>
@ -112,6 +118,24 @@ void ACCStream<T>::triad()
a[i] = b[i] + scalar * c[i];
}
}
template <class T>
T ACCStream<T>::dot()
{
T sum = 0.0;
unsigned int array_size = this->array_size;
T * restrict a = this->a;
T * restrict b = this->b;
#pragma acc kernels present(a[0:array_size], b[0:array_size]) wait
for (int i = 0; i < array_size; i++)
{
sum += a[i] * b[i];
}
return sum;
}
void listDevices(void)
{
// Get number of devices

3
ACCStream.h
View File

@ -35,8 +35,9 @@ class ACCStream : public Stream<T>
virtual void add() override;
virtual void mul() override;
virtual void triad() override;
virtual T dot() override;
virtual void write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c) override;
virtual void init_arrays(T initA, T initB, T initC) override;
virtual void read_arrays(std::vector<T>& a, std::vector<T>& b, std::vector<T>& c) override;

2
CMakeLists.txt
View File

@ -29,7 +29,7 @@ include(CheckIncludeFileCXX)
include(CheckCXXCompilerFlag)
set(gpu-stream_VERSION_MAJOR 2)
set(gpu-stream_VERSION_MINOR 1)
set(gpu-stream_VERSION_MINOR 2)
configure_file(common.h.in common.h)
include_directories(${CMAKE_BINARY_DIR})

72
CUDAStream.cu
View File

@ -8,8 +8,6 @@
#include "CUDAStream.h"
#define TBSIZE 1024
void check_error(void)
{
cudaError_t err = cudaGetLastError();
@ -47,6 +45,9 @@ CUDAStream<T>::CUDAStream(const unsigned int ARRAY_SIZE, const int device_index)
array_size = ARRAY_SIZE;
// Allocate the host array for partial sums for dot kernels
sums = (T*)malloc(sizeof(T) * DOT_NUM_BLOCKS);
// Check buffers fit on the device
cudaDeviceProp props;
cudaGetDeviceProperties(&props, 0);
@ -60,29 +61,42 @@ CUDAStream<T>::CUDAStream(const unsigned int ARRAY_SIZE, const int device_index)
check_error();
cudaMalloc(&d_c, ARRAY_SIZE*sizeof(T));
check_error();
cudaMalloc(&d_sum, DOT_NUM_BLOCKS*sizeof(T));
check_error();
}
template <class T>
CUDAStream<T>::~CUDAStream()
{
free(sums);
cudaFree(d_a);
check_error();
cudaFree(d_b);
check_error();
cudaFree(d_c);
check_error();
cudaFree(d_sum);
check_error();
}
template <typename T>
__global__ void init_kernel(T * a, T * b, T * c, T initA, T initB, T initC)
{
const int i = blockDim.x * blockIdx.x + threadIdx.x;
a[i] = initA;
b[i] = initB;
c[i] = initC;
}
template <class T>
void CUDAStream<T>::write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c)
void CUDAStream<T>::init_arrays(T initA, T initB, T initC)
{
// Copy host memory to device
cudaMemcpy(d_a, a.data(), a.size()*sizeof(T), cudaMemcpyHostToDevice);
init_kernel<<<array_size/TBSIZE, TBSIZE>>>(d_a, d_b, d_c, initA, initB, initC);
check_error();
cudaMemcpy(d_b, b.data(), b.size()*sizeof(T), cudaMemcpyHostToDevice);
check_error();
cudaMemcpy(d_c, c.data(), c.size()*sizeof(T), cudaMemcpyHostToDevice);
cudaDeviceSynchronize();
check_error();
}
@ -165,6 +179,48 @@ void CUDAStream<T>::triad()
check_error();
}
template <class T>
__global__ void dot_kernel(const T * a, const T * b, T * sum, unsigned int array_size)
{
extern __shared__ __align__(sizeof(T)) unsigned char smem[];
T *tb_sum = reinterpret_cast<T*>(smem);
int i = blockDim.x * blockIdx.x + threadIdx.x;
const size_t local_i = threadIdx.x;
tb_sum[local_i] = 0.0;
for (; i < array_size; i += blockDim.x*gridDim.x)
tb_sum[local_i] += a[i] * b[i];
for (int offset = blockDim.x / 2; offset > 0; offset /= 2)
{
__syncthreads();
if (local_i < offset)
{
tb_sum[local_i] += tb_sum[local_i+offset];
}
}
if (local_i == 0)
sum[blockIdx.x] = tb_sum[local_i];
}
template <class T>
T CUDAStream<T>::dot()
{
dot_kernel<<<DOT_NUM_BLOCKS, TBSIZE, sizeof(T)*TBSIZE>>>(d_a, d_b, d_sum, array_size);
check_error();
cudaMemcpy(sums, d_sum, DOT_NUM_BLOCKS*sizeof(T), cudaMemcpyDeviceToHost);
check_error();
T sum = 0.0;
for (int i = 0; i < DOT_NUM_BLOCKS; i++)
sum += sums[i];
return sum;
}
void listDevices(void)
{

11
CUDAStream.h
View File

@ -15,16 +15,24 @@
#define IMPLEMENTATION_STRING "CUDA"
#define TBSIZE 1024
#define DOT_NUM_BLOCKS 256
template <class T>
class CUDAStream : public Stream<T>
{
protected:
// Size of arrays
unsigned int array_size;
// Host array for partial sums for dot kernel
T *sums;
// Device side pointers to arrays
T *d_a;
T *d_b;
T *d_c;
T *d_sum;
public:
@ -36,8 +44,9 @@ class CUDAStream : public Stream<T>
virtual void add() override;
virtual void mul() override;
virtual void triad() override;
virtual T dot() override;
virtual void write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c) override;
virtual void init_arrays(T initA, T initB, T initC) override;
virtual void read_arrays(std::vector<T>& a, std::vector<T>& b, std::vector<T>& c) override;
};

20
HIPStream.cu
View File

@ -74,15 +74,21 @@ HIPStream<T>::~HIPStream()
check_error();
}
template <class T>
void HIPStream<T>::write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c)
template <typename T>
__global__ void init_kernel(hipLaunchParm lp, T * a, T * b, T * c, T initA, T initB, T initC)
{
// Copy host memory to device
hipMemcpy(d_a, a.data(), a.size()*sizeof(T), hipMemcpyHostToDevice);
const int i = hipBlockDim_x * hipBlockIdx_x + hipThreadIdx_x;
a[i] = initA;
b[i] = initB;
c[i] = initC;
}
template <class T>
void HIPStream<T>::init_arrays(T initA, T initB, T initC)
{
hipLaunchKernel(HIP_KERNEL_NAME(init_kernel), dim3(array_size/TBSIZE), dim3(TBSIZE), 0, 0, d_a, d_b, d_c, initA, initB, initC);
check_error();
hipMemcpy(d_b, b.data(), b.size()*sizeof(T), hipMemcpyHostToDevice);
check_error();
hipMemcpy(d_c, c.data(), c.size()*sizeof(T), hipMemcpyHostToDevice);
hipDeviceSynchronize();
check_error();
}

2
HIPStream.h
View File

@ -37,7 +37,7 @@ class HIPStream : public Stream<T>
virtual void mul() override;
virtual void triad() override;
virtual void write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c) override;
virtual void init_arrays(T initA, T initB, T initC) override;
virtual void read_arrays(std::vector<T>& a, std::vector<T>& b, std::vector<T>& c) override;
};

37
KOKKOSStream.cpp
View File

@ -34,18 +34,18 @@ KOKKOSStream<T>::~KOKKOSStream()
}
template <class T>
void KOKKOSStream<T>::write_arrays(
const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c)
void KOKKOSStream<T>::init_arrays(T initA, T initB, T initC)
{
for(int ii = 0; ii < array_size; ++ii)
View<double*, DEVICE> a(*d_a);
View<double*, DEVICE> b(*d_b);
View<double*, DEVICE> c(*d_c);
parallel_for(array_size, KOKKOS_LAMBDA (const int index)
{
(*hm_a)(ii) = a[ii];
(*hm_b)(ii) = b[ii];
(*hm_c)(ii) = c[ii];
}
deep_copy(*d_a, *hm_a);
deep_copy(*d_b, *hm_b);
deep_copy(*d_c, *hm_c);
a[index] = initA;
b[index] - initB;
c[index] = initC;
});
Kokkos::fence();
}
template <class T>
@ -121,6 +121,23 @@ void KOKKOSStream<T>::triad()
Kokkos::fence();
}
template <class T>
T KOKKOSStream<T>::dot()
{
View<double *, DEVICE> a(*d_a);
View<double *, DEVICE> b(*d_b);
T sum = 0.0;
parallel_reduce(array_size, KOKKOS_LAMBDA (const int index, double &tmp)
{
tmp += a[index] * b[index];
}, sum);
return sum;
}
void listDevices(void)
{
std::cout << "This is not the device you are looking for.";

4
KOKKOSStream.hpp
View File

@ -47,9 +47,9 @@ class KOKKOSStream : public Stream<T>
virtual void add() override;
virtual void mul() override;
virtual void triad() override;
virtual T dot() override;
virtual void write_arrays(
const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c) override;
virtual void init_arrays(T initA, T initB, T initC) override;
virtual void read_arrays(
std::vector<T>& a, std::vector<T>& b, std::vector<T>& c) override;
};

82
OCLStream.cpp
View File

@ -16,6 +16,18 @@ std::string kernels{R"CLC(
constant TYPE scalar = startScalar;
kernel void init(
global TYPE * restrict a,
global TYPE * restrict b,
global TYPE * restrict c,
TYPE initA, TYPE initB, TYPE initC)
{
const size_t i = get_global_id(0);
a[i] = initA;
b[i] = initB;
c[i] = initC;
}
kernel void copy(
global const TYPE * restrict a,
global TYPE * restrict c)
@ -50,6 +62,32 @@ std::string kernels{R"CLC(
a[i] = b[i] + scalar * c[i];
}
kernel void stream_dot(
global const TYPE * restrict a,
global const TYPE * restrict b,
global TYPE * restrict sum,
local TYPE * restrict wg_sum,
int array_size)
{
size_t i = get_global_id(0);
const size_t local_i = get_local_id(0);
wg_sum[local_i] = 0.0;
for (; i < array_size; i += get_global_size(0))
wg_sum[local_i] += a[i] * b[i];
for (int offset = get_local_size(0) / 2; offset > 0; offset /= 2)
{
barrier(CLK_LOCAL_MEM_FENCE);
if (local_i < offset)
{
wg_sum[local_i] += wg_sum[local_i+offset];
}
}
if (local_i == 0)
sum[get_group_id(0)] = wg_sum[local_i];
}
)CLC"};
@ -64,9 +102,22 @@ OCLStream<T>::OCLStream(const unsigned int ARRAY_SIZE, const int device_index)
throw std::runtime_error("Invalid device index");
device = devices[device_index];
// Determine sensible dot kernel NDRange configuration
if (device.getInfo<CL_DEVICE_TYPE>() & CL_DEVICE_TYPE_CPU)
{
dot_num_groups = device.getInfo<CL_DEVICE_MAX_COMPUTE_UNITS>();
dot_wgsize = device.getInfo<CL_DEVICE_NATIVE_VECTOR_WIDTH_DOUBLE>() * 2;
}
else
{
dot_num_groups = device.getInfo<CL_DEVICE_MAX_COMPUTE_UNITS>() * 4;
dot_wgsize = device.getInfo<CL_DEVICE_MAX_WORK_GROUP_SIZE>();
}
// Print out device information
std::cout << "Using OpenCL device " << getDeviceName(device_index) << std::endl;
std::cout << "Driver: " << getDeviceDriver(device_index) << std::endl;
std::cout << "Reduction kernel config: " << dot_num_groups << " groups of size " << dot_wgsize << std::endl;
context = cl::Context(device);
queue = cl::CommandQueue(context);
@ -101,10 +152,12 @@ OCLStream<T>::OCLStream(const unsigned int ARRAY_SIZE, const int device_index)
}
// Create kernels
init_kernel = new cl::KernelFunctor<cl::Buffer, cl::Buffer, cl::Buffer, T, T, T>(program, "init");
copy_kernel = new cl::KernelFunctor<cl::Buffer, cl::Buffer>(program, "copy");
mul_kernel = new cl::KernelFunctor<cl::Buffer, cl::Buffer>(program, "mul");
add_kernel = new cl::KernelFunctor<cl::Buffer, cl::Buffer, cl::Buffer>(program, "add");
triad_kernel = new cl::KernelFunctor<cl::Buffer, cl::Buffer, cl::Buffer>(program, "triad");
dot_kernel = new cl::KernelFunctor<cl::Buffer, cl::Buffer, cl::Buffer, cl::LocalSpaceArg, cl_int>(program, "stream_dot");
array_size = ARRAY_SIZE;
@ -120,12 +173,15 @@ OCLStream<T>::OCLStream(const unsigned int ARRAY_SIZE, const int device_index)
d_a = cl::Buffer(context, CL_MEM_READ_WRITE, sizeof(T) * ARRAY_SIZE);
d_b = cl::Buffer(context, CL_MEM_READ_WRITE, sizeof(T) * ARRAY_SIZE);
d_c = cl::Buffer(context, CL_MEM_READ_WRITE, sizeof(T) * ARRAY_SIZE);
d_sum = cl::Buffer(context, CL_MEM_WRITE_ONLY, sizeof(T) * dot_num_groups);
sums = std::vector<T>(dot_num_groups);
}
template <class T>
OCLStream<T>::~OCLStream()
{
delete init_kernel;
delete copy_kernel;
delete mul_kernel;
delete add_kernel;
@ -173,11 +229,29 @@ void OCLStream<T>::triad()
}
template <class T>
void OCLStream<T>::write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c)
T OCLStream<T>::dot()
{
cl::copy(queue, a.begin(), a.end(), d_a);
cl::copy(queue, b.begin(), b.end(), d_b);
cl::copy(queue, c.begin(), c.end(), d_c);
(*dot_kernel)(
cl::EnqueueArgs(queue, cl::NDRange(dot_num_groups*dot_wgsize), cl::NDRange(dot_wgsize)),
d_a, d_b, d_sum, cl::Local(sizeof(T) * dot_wgsize), array_size
);
cl::copy(queue, d_sum, sums.begin(), sums.end());
T sum = 0.0;
for (T val : sums)
sum += val;
return sum;
}
template <class T>
void OCLStream<T>::init_arrays(T initA, T initB, T initC)
{
(*init_kernel)(
cl::EnqueueArgs(queue, cl::NDRange(array_size)),
d_a, d_b, d_c, initA, initB, initC
);
queue.finish();
}
template <class T>

13
OCLStream.h
View File

@ -28,20 +28,30 @@ class OCLStream : public Stream<T>
// Size of arrays
unsigned int array_size;
// Host array for partial sums for dot kernel
std::vector<T> sums;
// Device side pointers to arrays
cl::Buffer d_a;
cl::Buffer d_b;
cl::Buffer d_c;
cl::Buffer d_sum;
// OpenCL objects
cl::Device device;
cl::Context context;
cl::CommandQueue queue;
cl::KernelFunctor<cl::Buffer, cl::Buffer, cl::Buffer, T, T, T> *init_kernel;
cl::KernelFunctor<cl::Buffer, cl::Buffer> *copy_kernel;
cl::KernelFunctor<cl::Buffer, cl::Buffer> * mul_kernel;
cl::KernelFunctor<cl::Buffer, cl::Buffer, cl::Buffer> *add_kernel;
cl::KernelFunctor<cl::Buffer, cl::Buffer, cl::Buffer> *triad_kernel;
cl::KernelFunctor<cl::Buffer, cl::Buffer, cl::Buffer, cl::LocalSpaceArg, cl_int> *dot_kernel;
// NDRange configuration for the dot kernel
size_t dot_num_groups;
size_t dot_wgsize;
public:
@ -52,8 +62,9 @@ class OCLStream : public Stream<T>
virtual void add() override;
virtual void mul() override;
virtual void triad() override;
virtual T dot() override;
virtual void write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c) override;
virtual void init_arrays(T initA, T initB, T initC) override;
virtual void read_arrays(std::vector<T>& a, std::vector<T>& b, std::vector<T>& c) override;
};

111
OMP3Stream.cpp
View File

@ -1,111 +0,0 @@
// Copyright (c) 2015-16 Tom Deakin, Simon McIntosh-Smith,
// University of Bristol HPC
//
// For full license terms please see the LICENSE file distributed with this
// source code
#include "OMP3Stream.h"
template <class T>
OMP3Stream<T>::OMP3Stream(const unsigned int ARRAY_SIZE, T *a, T *b, T *c)
{
array_size = ARRAY_SIZE;
this->a = (T*)malloc(sizeof(T)*array_size);
this->b = (T*)malloc(sizeof(T)*array_size);
this->c = (T*)malloc(sizeof(T)*array_size);
}
template <class T>
OMP3Stream<T>::~OMP3Stream()
{
free(a);
free(b);
free(c);
}
template <class T>
void OMP3Stream<T>::write_arrays(const std::vector<T>& h_a, const std::vector<T>& h_b, const std::vector<T>& h_c)
{
#pragma omp parallel for
for (int i = 0; i < array_size; i++)
{
a[i] = h_a[i];
b[i] = h_b[i];
c[i] = h_c[i];
}
}
template <class T>
void OMP3Stream<T>::read_arrays(std::vector<T>& h_a, std::vector<T>& h_b, std::vector<T>& h_c)
{
#pragma omp parallel for
for (int i = 0; i < array_size; i++)
{
h_a[i] = a[i];
h_b[i] = b[i];
h_c[i] = c[i];
}
}
template <class T>
void OMP3Stream<T>::copy()
{
#pragma omp parallel for
for (int i = 0; i < array_size; i++)
{
c[i] = a[i];
}
}
template <class T>
void OMP3Stream<T>::mul()
{
const T scalar = startScalar;
#pragma omp parallel for
for (int i = 0; i < array_size; i++)
{
b[i] = scalar * c[i];
}
}
template <class T>
void OMP3Stream<T>::add()
{
#pragma omp parallel for
for (int i = 0; i < array_size; i++)
{
c[i] = a[i] + b[i];
}
}
template <class T>
void OMP3Stream<T>::triad()
{
const T scalar = startScalar;
#pragma omp parallel for
for (int i = 0; i < array_size; i++)
{
a[i] = b[i] + scalar * c[i];
}
}
void listDevices(void)
{
std::cout << "0: CPU" << std::endl;
}
std::string getDeviceName(const int)
{
return std::string("Device name unavailable");
}
std::string getDeviceDriver(const int)
{
return std::string("Device driver unavailable");
}
template class OMP3Stream<float>;
template class OMP3Stream<double>;

40
OMP3Stream.h
View File

@ -1,40 +0,0 @@
// Copyright (c) 2015-16 Tom Deakin, Simon McIntosh-Smith,
// University of Bristol HPC
//
// For full license terms please see the LICENSE file distributed with this
// source code
#pragma once
#include <iostream>
#include <stdexcept>
#include "Stream.h"
#define IMPLEMENTATION_STRING "Reference OpenMP"
template <class T>
class OMP3Stream : public Stream<T>
{
protected:
// Size of arrays
unsigned int array_size;
// Device side pointers
T *a;
T *b;
T *c;
public:
OMP3Stream(const unsigned int, T*, T*, T*);
~OMP3Stream();
virtual void copy() override;
virtual void add() override;
virtual void mul() override;
virtual void triad() override;
virtual void write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c) override;
virtual void read_arrays(std::vector<T>& a, std::vector<T>& b, std::vector<T>& c) override;
};

108
OMP45Stream.cpp → OMPStream.cpp
View File

@ -5,26 +5,33 @@
// For full license terms please see the LICENSE file distributed with this
// source code
#include "OMP45Stream.h"
#include "OMPStream.h"
template <class T>
OMP45Stream<T>::OMP45Stream(const unsigned int ARRAY_SIZE, T *a, T *b, T *c, int device)
OMPStream<T>::OMPStream(const unsigned int ARRAY_SIZE, T *a, T *b, T *c, int device)
{
omp_set_default_device(device);
array_size = ARRAY_SIZE;
#ifdef OMP_TARGET_GPU
omp_set_default_device(device);
// Set up data region on device
this->a = a;
this->b = b;
this->c = c;
#pragma omp target enter data map(to: a[0:array_size], b[0:array_size], c[0:array_size])
#pragma omp target enter data map(alloc: a[0:array_size], b[0:array_size], c[0:array_size])
{}
#else
// Allocate on the host
this->a = (T*)malloc(sizeof(T)*array_size);
this->b = (T*)malloc(sizeof(T)*array_size);
this->c = (T*)malloc(sizeof(T)*array_size);
#endif
}
template <class T>
OMP45Stream<T>::~OMP45Stream()
OMPStream<T>::~OMPStream()
{
#ifdef OMP_TARGET_GPU
// End data region on device
unsigned int array_size = this->array_size;
T *a = this->a;
@ -32,35 +39,64 @@ OMP45Stream<T>::~OMP45Stream()
T *c = this->c;
#pragma omp target exit data map(release: a[0:array_size], b[0:array_size], c[0:array_size])
{}
#else
free(a);
free(b);
free(c);
#endif
}
template <class T>
void OMP45Stream<T>::write_arrays(const std::vector<T>& h_a, const std::vector<T>& h_b, const std::vector<T>& h_c)
void OMPStream<T>::init_arrays(T initA, T initB, T initC)
{
unsigned int array_size = this->array_size;
#ifdef OMP_TARGET_GPU
T *a = this->a;
T *b = this->b;
T *c = this->c;
#pragma omp target update to(a[0:array_size], b[0:array_size], c[0:array_size])
{}
#pragma omp target teams distribute parallel for simd map(to: a[0:array_size], b[0:array_size], c[0:array_size])
#else
#pragma omp parallel for
#endif
for (int i = 0; i < array_size; i++)
{
a[i] = initA;
b[i] = initB;
c[i] = initC;
}
}
template <class T>
void OMP45Stream<T>::read_arrays(std::vector<T>& h_a, std::vector<T>& h_b, std::vector<T>& h_c)
void OMPStream<T>::read_arrays(std::vector<T>& h_a, std::vector<T>& h_b, std::vector<T>& h_c)
{
#ifdef OMP_TARGET_GPU
T *a = this->a;
T *b = this->b;
T *c = this->c;
#pragma omp target update from(a[0:array_size], b[0:array_size], c[0:array_size])
{}
#else
#pragma omp parallel for
for (int i = 0; i < array_size; i++)
{
h_a[i] = a[i];
h_b[i] = b[i];
h_c[i] = c[i];
}
#endif
}
template <class T>
void OMP45Stream<T>::copy()
void OMPStream<T>::copy()
{
#ifdef OMP_TARGET_GPU
unsigned int array_size = this->array_size;
T *a = this->a;
T *c = this->c;
#pragma omp target teams distribute parallel for simd map(to: a[0:array_size], c[0:array_size])
#else
#pragma omp parallel for
#endif
for (int i = 0; i < array_size; i++)
{
c[i] = a[i];
@ -68,14 +104,18 @@ void OMP45Stream<T>::copy()
}
template <class T>
void OMP45Stream<T>::mul()
void OMPStream<T>::mul()
{
const T scalar = startScalar;
#ifdef OMP_TARGET_GPU
unsigned int array_size = this->array_size;
T *b = this->b;
T *c = this->c;
#pragma omp target teams distribute parallel for simd map(to: b[0:array_size], c[0:array_size])
#else
#pragma omp parallel for
#endif
for (int i = 0; i < array_size; i++)
{
b[i] = scalar * c[i];
@ -83,13 +123,17 @@ void OMP45Stream<T>::mul()
}
template <class T>
void OMP45Stream<T>::add()
void OMPStream<T>::add()
{
#ifdef OMP_TARGET_GPU
unsigned int array_size = this->array_size;
T *a = this->a;
T *b = this->b;
T *c = this->c;
#pragma omp target teams distribute parallel for simd map(to: a[0:array_size], b[0:array_size], c[0:array_size])
#else
#pragma omp parallel for
#endif
for (int i = 0; i < array_size; i++)
{
c[i] = a[i] + b[i];
@ -97,22 +141,51 @@ void OMP45Stream<T>::add()
}
template <class T>
void OMP45Stream<T>::triad()
void OMPStream<T>::triad()
{
const T scalar = startScalar;
#ifdef OMP_TARGET_GPU
unsigned int array_size = this->array_size;
T *a = this->a;
T *b = this->b;
T *c = this->c;
#pragma omp target teams distribute parallel for simd map(to: a[0:array_size], b[0:array_size], c[0:array_size])
#else
#pragma omp parallel for
#endif
for (int i = 0; i < array_size; i++)
{
a[i] = b[i] + scalar * c[i];
}
}
template <class T>
T OMPStream<T>::dot()
{
T sum = 0.0;
#ifdef OMP_TARGET_GPU
unsigned int array_size = this->array_size;
T *a = this->a;
T *b = this->b;
#pragma omp target teams distribute parallel for simd reduction(+:sum) map(tofrom: sum)
#else
#pragma omp parallel for reduction(+:sum)
#endif
for (int i = 0; i < array_size; i++)
{
sum += a[i] * b[i];
}
return sum;
}
void listDevices(void)
{
#ifdef OMP_TARGET_GPU
// Get number of devices
int count = omp_get_num_devices();
@ -125,6 +198,9 @@ void listDevices(void)
{
std::cout << "There are " << count << " devices." << std::endl;
}
#else
std::cout << "0: CPU" << std::endl;
#endif
}
std::string getDeviceName(const int)
@ -136,5 +212,5 @@ std::string getDeviceDriver(const int)
{
return std::string("Device driver unavailable");
}
template class OMP45Stream<float>;
template class OMP45Stream<double>;
template class OMPStream<float>;
template class OMPStream<double>;

11
OMP45Stream.h → OMPStream.h
View File

@ -14,10 +14,10 @@
#include <omp.h>
#define IMPLEMENTATION_STRING "OpenMP 4.5"
#define IMPLEMENTATION_STRING "OpenMP"
template <class T>
class OMP45Stream : public Stream<T>
class OMPStream : public Stream<T>
{
protected:
// Size of arrays
@ -29,15 +29,16 @@ class OMP45Stream : public Stream<T>
T *c;
public:
OMP45Stream(const unsigned int, T*, T*, T*, int);
~OMP45Stream();
OMPStream(const unsigned int, T*, T*, T*, int);
~OMPStream();
virtual void copy() override;
virtual void add() override;
virtual void mul() override;
virtual void triad() override;
virtual T dot() override;
virtual void write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c) override;
virtual void init_arrays(T initA, T initB, T initC) override;
virtual void read_arrays(std::vector<T>& a, std::vector<T>& b, std::vector<T>& c) override;

38
RAJAStream.cpp
View File

@ -21,12 +21,6 @@ RAJAStream<T>::RAJAStream(const unsigned int ARRAY_SIZE, const int device_index)
d_a = new T[ARRAY_SIZE];
d_b = new T[ARRAY_SIZE];
d_c = new T[ARRAY_SIZE];
forall<policy>(index_set, [=] RAJA_DEVICE (int index)
{
d_a[index] = 0.0;
d_b[index] = 0.0;
d_c[index] = 0.0;
});
#else
cudaMallocManaged((void**)&d_a, sizeof(T)*ARRAY_SIZE, cudaMemAttachGlobal);
cudaMallocManaged((void**)&d_b, sizeof(T)*ARRAY_SIZE, cudaMemAttachGlobal);
@ -50,12 +44,17 @@ RAJAStream<T>::~RAJAStream()
}
template <class T>
void RAJAStream<T>::write_arrays(
const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c)
void RAJAStream<T>::init_arrays(T initA, T initB, T initC)
{
std::copy(a.begin(), a.end(), d_a);
std::copy(b.begin(), b.end(), d_b);
std::copy(c.begin(), c.end(), d_c);
T* a = d_a;
T* b = d_b;
T* c = d_c;
forall<policy>(index_set, [=] RAJA_DEVICE (int index)
{
a[index] = initA;
b[index] = initB;
c[index] = initC;
});
}
template <class T>
@ -115,6 +114,23 @@ void RAJAStream<T>::triad()
});
}
template <class T>
T RAJAStream<T>::dot()
{
T* a = d_a;
T* b = d_b;
RAJA::ReduceSum<reduce_policy, T> sum(0.0);
forall<policy>(index_set, [=] RAJA_DEVICE (int index)
{
sum += a[index] * b[index];
});
return T(sum);
}
void listDevices(void)
{
std::cout << "This is not the device you are looking for.";

6
RAJAStream.hpp
View File

@ -18,11 +18,13 @@
typedef RAJA::IndexSet::ExecPolicy<
RAJA::seq_segit,
RAJA::omp_parallel_for_exec> policy;
typedef RAJA::omp_reduce reduce_policy;
#else
const size_t block_size = 128;
typedef RAJA::IndexSet::ExecPolicy<
RAJA::seq_segit,
RAJA::cuda_exec<block_size>> policy;
typedef RAJA::cuda_reduce<block_size> reduce_policy;
#endif
template <class T>
@ -49,9 +51,9 @@ class RAJAStream : public Stream<T>
virtual void add() override;
virtual void mul() override;
virtual void triad() override;
virtual T dot() override;
virtual void write_arrays(
const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c) override;
virtual void init_arrays(T initA, T initB, T initC) override;
virtual void read_arrays(
std::vector<T>& a, std::vector<T>& b, std::vector<T>& c) override;
};

129
SYCLStream.cpp
View File

@ -18,14 +18,6 @@ std::vector<device> devices;
void getDeviceList(void);
program * p;
/* Forward declaration of SYCL kernels */
namespace kernels {
class copy;
class mul;
class add;
class triad;
}
template <class T>
SYCLStream<T>::SYCLStream(const unsigned int ARRAY_SIZE, const int device_index)
{
@ -38,24 +30,39 @@ SYCLStream<T>::SYCLStream(const unsigned int ARRAY_SIZE, const int device_index)
throw std::runtime_error("Invalid device index");
device dev = devices[device_index];
// Determine sensible dot kernel NDRange configuration
if (dev.is_cpu())
{
dot_num_groups = dev.get_info<info::device::max_compute_units>();
dot_wgsize = dev.get_info<info::device::native_vector_width_double>() * 2;
}
else
{
dot_num_groups = dev.get_info<info::device::max_compute_units>() * 4;
dot_wgsize = dev.get_info<info::device::max_work_group_size>();
}
// Print out device information
std::cout << "Using SYCL device " << getDeviceName(device_index) << std::endl;
std::cout << "Driver: " << getDeviceDriver(device_index) << std::endl;
std::cout << "Reduction kernel config: " << dot_num_groups << " groups of size " << dot_wgsize << std::endl;
queue = new cl::sycl::queue(dev);
/* Pre-build the kernels */
p = new program(queue->get_context());
p->build_from_kernel_name<kernels::copy>();
p->build_from_kernel_name<kernels::mul>();
p->build_from_kernel_name<kernels::add>();
p->build_from_kernel_name<kernels::triad>();
p->build_from_kernel_name<init_kernel>();
p->build_from_kernel_name<copy_kernel>();
p->build_from_kernel_name<mul_kernel>();
p->build_from_kernel_name<add_kernel>();
p->build_from_kernel_name<triad_kernel>();
p->build_from_kernel_name<dot_kernel>();
// Create buffers
d_a = new buffer<T>(array_size);
d_b = new buffer<T>(array_size);
d_c = new buffer<T>(array_size);
d_sum = new buffer<T>(dot_num_groups);
}
template <class T>
@ -64,6 +71,7 @@ SYCLStream<T>::~SYCLStream()
delete d_a;
delete d_b;
delete d_c;
delete d_sum;
delete p;
delete queue;
@ -76,11 +84,11 @@ void SYCLStream<T>::copy()
{
auto ka = d_a->template get_access<access::mode::read>(cgh);
auto kc = d_c->template get_access<access::mode::write>(cgh);
cgh.parallel_for<kernels::copy>(p->get_kernel<kernels::copy>(),
cgh.parallel_for<copy_kernel>(p->get_kernel<copy_kernel>(),
range<1>{array_size}, [=](item<1> item)
{
auto id = item.get();
kc[id[0]] = ka[id[0]];
auto id = item.get()[0];
kc[id] = ka[id];
});
});
queue->wait();
@ -94,11 +102,11 @@ void SYCLStream<T>::mul()
{
auto kb = d_b->template get_access<access::mode::write>(cgh);
auto kc = d_c->template get_access<access::mode::read>(cgh);
cgh.parallel_for<kernels::mul>(p->get_kernel<kernels::mul>(),
cgh.parallel_for<mul_kernel>(p->get_kernel<mul_kernel>(),
range<1>{array_size}, [=](item<1> item)
{
auto id = item.get();
kb[id[0]] = scalar * kc[id[0]];
auto id = item.get()[0];
kb[id] = scalar * kc[id];
});
});
queue->wait();
@ -112,11 +120,11 @@ void SYCLStream<T>::add()
auto ka = d_a->template get_access<access::mode::read>(cgh);
auto kb = d_b->template get_access<access::mode::read>(cgh);
auto kc = d_c->template get_access<access::mode::write>(cgh);
cgh.parallel_for<kernels::add>(p->get_kernel<kernels::add>(),
cgh.parallel_for<add_kernel>(p->get_kernel<add_kernel>(),
range<1>{array_size}, [=](item<1> item)
{
auto id = item.get();
kc[id[0]] = ka[id[0]] + kb[id[0]];
auto id = item.get()[0];
kc[id] = ka[id] + kb[id];
});
});
queue->wait();
@ -131,28 +139,81 @@ void SYCLStream<T>::triad()
auto ka = d_a->template get_access<access::mode::write>(cgh);
auto kb = d_b->template get_access<access::mode::read>(cgh);
auto kc = d_c->template get_access<access::mode::read>(cgh);
cgh.parallel_for<kernels::triad>(p->get_kernel<kernels::triad>(),
cgh.parallel_for<triad_kernel>(p->get_kernel<triad_kernel>(),
range<1>{array_size}, [=](item<1> item)
{
auto id = item.get();
ka[id] = kb[id[0]] + scalar * kc[id[0]];
auto id = item.get()[0];
ka[id] = kb[id] + scalar * kc[id];
});
});
queue->wait();
}
template <class T>
void SYCLStream<T>::write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c)
T SYCLStream<T>::dot()
{
auto _a = d_a->template get_access<access::mode::write, access::target::host_buffer>();
auto _b = d_b->template get_access<access::mode::write, access::target::host_buffer>();
auto _c = d_c->template get_access<access::mode::write, access::target::host_buffer>();
for (int i = 0; i < array_size; i++)
queue->submit([&](handler &cgh)
{
_a[i] = a[i];
_b[i] = b[i];
_c[i] = c[i];
auto ka = d_a->template get_access<access::mode::read>(cgh);
auto kb = d_b->template get_access<access::mode::read>(cgh);
auto ksum = d_sum->template get_access<access::mode::write>(cgh);
auto wg_sum = accessor<T, 1, access::mode::read_write, access::target::local>(range<1>(dot_wgsize), cgh);
size_t N = array_size;
cgh.parallel_for<dot_kernel>(p->get_kernel<dot_kernel>(),
nd_range<1>(dot_num_groups*dot_wgsize, dot_wgsize), [=](nd_item<1> item)
{
size_t i = item.get_global(0);
size_t li = item.get_local(0);
size_t global_size = item.get_global_range()[0];
wg_sum[li] = 0.0;
for (; i < N; i += global_size)
wg_sum[li] += ka[i] * kb[i];
size_t local_size = item.get_local_range()[0];
for (int offset = local_size / 2; offset > 0; offset /= 2)
{
item.barrier(cl::sycl::access::fence_space::local_space);
if (li < offset)
wg_sum[li] += wg_sum[li + offset];
}
if (li == 0)
ksum[item.get_group(0)] = wg_sum[0];
});
});
T sum = 0.0;
auto h_sum = d_sum->template get_access<access::mode::read, access::target::host_buffer>();
for (int i = 0; i < dot_num_groups; i++)
{
sum += h_sum[i];
}
return sum;
}
template <class T>
void SYCLStream<T>::init_arrays(T initA, T initB, T initC)
{
queue->submit([&](handler &cgh)
{
auto ka = d_a->template get_access<access::mode::write>(cgh);
auto kb = d_b->template get_access<access::mode::write>(cgh);
auto kc = d_c->template get_access<access::mode::write>(cgh);
cgh.parallel_for<init_kernel>(p->get_kernel<init_kernel>(),
range<1>{array_size}, [=](item<1> item)
{
auto id = item.get()[0];
ka[id] = initA;
kb[id] = initB;
kc[id] = initC;
});
});
queue->wait();
}
template <class T>
@ -244,5 +305,5 @@ std::string getDeviceDriver(const int device)
// TODO: Fix kernel names to allow multiple template specializations
//template class SYCLStream<float>;
template class SYCLStream<float>;
template class SYCLStream<double>;

26
SYCLStream.h
View File

@ -15,6 +15,16 @@
#define IMPLEMENTATION_STRING "SYCL"
namespace sycl_kernels
{
template <class T> class init;
template <class T> class copy;
template <class T> class mul;
template <class T> class add;
template <class T> class triad;
template <class T> class dot;
}
template <class T>
class SYCLStream : public Stream<T>
{
@ -27,6 +37,19 @@ class SYCLStream : public Stream<T>
cl::sycl::buffer<T> *d_a;
cl::sycl::buffer<T> *d_b;
cl::sycl::buffer<T> *d_c;
cl::sycl::buffer<T> *d_sum;
// SYCL kernel names
typedef sycl_kernels::init<T> init_kernel;
typedef sycl_kernels::copy<T> copy_kernel;
typedef sycl_kernels::mul<T> mul_kernel;
typedef sycl_kernels::add<T> add_kernel;
typedef sycl_kernels::triad<T> triad_kernel;
typedef sycl_kernels::dot<T> dot_kernel;
// NDRange configuration for the dot kernel
size_t dot_num_groups;
size_t dot_wgsize;
public:
@ -37,8 +60,9 @@ class SYCLStream : public Stream<T>
virtual void add() override;
virtual void mul() override;
virtual void triad() override;
virtual T dot() override;
virtual void write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c) override;
virtual void init_arrays(T initA, T initB, T initC) override;
virtual void read_arrays(std::vector<T>& a, std::vector<T>& b, std::vector<T>& c) override;
};

3
Stream.h
View File

@ -29,9 +29,10 @@ class Stream
virtual void mul() = 0;
virtual void add() = 0;
virtual void triad() = 0;
virtual T dot() = 0;
// Copy memory between host and device
virtual void write_arrays(const std::vector<T>& a, const std::vector<T>& b, const std::vector<T>& c) = 0;
virtual void init_arrays(T initA, T initB, T initC) = 0;
virtual void read_arrays(std::vector<T>& a, std::vector<T>& b, std::vector<T>& c) = 0;
};

9
common.h.in
View File

@ -1,9 +0,0 @@
// Copyright (c) 2015-16 Tom Deakin, Simon McIntosh-Smith,
// University of Bristol HPC
//
// For full license terms please see the LICENSE file distributed with this
// source code
#define VERSION_STRING "@gpu-stream_VERSION_MAJOR@.@gpu-stream_VERSION_MINOR@"

69
main.cpp
View File

@ -15,7 +15,8 @@
#include <iomanip>
#include <cstring>
#include "common.h"
#define VERSION_STRING "3.0"
#include "Stream.h"
#if defined(CUDA)
@ -34,10 +35,8 @@
#include "ACCStream.h"
#elif defined(SYCL)
#include "SYCLStream.h"
#elif defined(OMP3)
#include "OMP3Stream.h"
#elif defined(OMP45)
#include "OMP45Stream.h"
#elif defined(OMP)
#include "OMPStream.h"
#endif
// Default size of 2^25
@ -47,7 +46,7 @@ unsigned int deviceIndex = 0;
bool use_float = false;
template <typename T>
void check_solution(const unsigned int ntimes, std::vector<T>& a, std::vector<T>& b, std::vector<T>& c);
void check_solution(const unsigned int ntimes, std::vector<T>& a, std::vector<T>& b, std::vector<T>& c, T& sum);
template <typename T>
void run();
@ -63,13 +62,11 @@ int main(int argc, char *argv[])
parseArguments(argc, argv);
// TODO: Fix SYCL to allow multiple template specializations
#ifndef SYCL
// TODO: Fix Kokkos to allow multiple template specializations
#ifndef KOKKOS
if (use_float)
run<float>();
else
#endif
#endif
run<double>();
@ -86,9 +83,9 @@ void run()
std::cout << "Precision: double" << std::endl;
// Create host vectors
std::vector<T> a(ARRAY_SIZE, startA);
std::vector<T> b(ARRAY_SIZE, startB);
std::vector<T> c(ARRAY_SIZE, startC);
std::vector<T> a(ARRAY_SIZE);
std::vector<T> b(ARRAY_SIZE);
std::vector<T> c(ARRAY_SIZE);
std::streamsize ss = std::cout.precision();
std::cout << std::setprecision(1) << std::fixed
<< "Array size: " << ARRAY_SIZE*sizeof(T)*1.0E-6 << " MB"
@ -97,6 +94,9 @@ void run()
<< " (=" << 3.0*ARRAY_SIZE*sizeof(T)*1.0E-9 << " GB)" << std::endl;
std::cout.precision(ss);
// Result of the Dot kernel
T sum;
Stream<T> *stream;
#if defined(CUDA)
@ -131,20 +131,16 @@ void run()
// Use the SYCL implementation
stream = new SYCLStream<T>(ARRAY_SIZE, deviceIndex);
#elif defined(OMP3)
// Use the "reference" OpenMP 3 implementation
stream = new OMP3Stream<T>(ARRAY_SIZE, a.data(), b.data(), c.data());
#elif defined(OMP45)
// Use the "reference" OpenMP 3 implementation
stream = new OMP45Stream<T>(ARRAY_SIZE, a.data(), b.data(), c.data(), deviceIndex);
#elif defined(OMP)
// Use the OpenMP implementation
stream = new OMPStream<T>(ARRAY_SIZE, a.data(), b.data(), c.data(), deviceIndex);
#endif
stream->write_arrays(a, b, c);
stream->init_arrays(startA, startB, startC);
// List of times
std::vector<std::vector<double>> timings(4);
std::vector<std::vector<double>> timings(5);
// Declare timers
std::chrono::high_resolution_clock::time_point t1, t2;
@ -176,11 +172,17 @@ void run()
t2 = std::chrono::high_resolution_clock::now();
timings[3].push_back(std::chrono::duration_cast<std::chrono::duration<double> >(t2 - t1).count());
// Execute Dot
t1 = std::chrono::high_resolution_clock::now();
sum = stream->dot();
t2 = std::chrono::high_resolution_clock::now();
timings[4].push_back(std::chrono::duration_cast<std::chrono::duration<double> >(t2 - t1).count());
}
// Check solutions
stream->read_arrays(a, b, c);
check_solution<T>(num_times, a, b, c);
check_solution<T>(num_times, a, b, c, sum);
// Display timing results
std::cout
@ -192,15 +194,16 @@ void run()
std::cout << std::fixed;
std::string labels[4] = {"Copy", "Mul", "Add", "Triad"};
size_t sizes[4] = {
std::string labels[5] = {"Copy", "Mul", "Add", "Triad", "Dot"};
size_t sizes[5] = {
2 * sizeof(T) * ARRAY_SIZE,
2 * sizeof(T) * ARRAY_SIZE,
3 * sizeof(T) * ARRAY_SIZE,
3 * sizeof(T) * ARRAY_SIZE
3 * sizeof(T) * ARRAY_SIZE,
2 * sizeof(T) * ARRAY_SIZE
};
for (int i = 0; i < 4; i++)
for (int i = 0; i < 5; i++)
{
// Get min/max; ignore the first result
auto minmax = std::minmax_element(timings[i].begin()+1, timings[i].end());
@ -224,12 +227,13 @@ void run()
}
template <typename T>
void check_solution(const unsigned int ntimes, std::vector<T>& a, std::vector<T>& b, std::vector<T>& c)
void check_solution(const unsigned int ntimes, std::vector<T>& a, std::vector<T>& b, std::vector<T>& c, T& sum)
{
// Generate correct solution
T goldA = startA;
T goldB = startB;
T goldC = startC;
T goldSum = 0.0;
const T scalar = startScalar;
@ -242,6 +246,9 @@ void check_solution(const unsigned int ntimes, std::vector<T>& a, std::vector<T>
goldA = goldB + scalar * goldC;
}
// Do the reduction
goldSum = goldA * goldB * ARRAY_SIZE;
// Calculate the average error
double errA = std::accumulate(a.begin(), a.end(), 0.0, [&](double sum, const T val){ return sum + fabs(val - goldA); });
errA /= a.size();
@ -249,6 +256,7 @@ void check_solution(const unsigned int ntimes, std::vector<T>& a, std::vector<T>
errB /= b.size();
double errC = std::accumulate(c.begin(), c.end(), 0.0, [&](double sum, const T val){ return sum + fabs(val - goldC); });
errC /= c.size();
double errSum = fabs(sum - goldSum);
double epsi = std::numeric_limits<T>::epsilon() * 100.0;
@ -264,6 +272,13 @@ void check_solution(const unsigned int ntimes, std::vector<T>& a, std::vector<T>
std::cerr
<< "Validation failed on c[]. Average error " << errC
<< std::endl;
// Check sum to 8 decimal places
if (errSum > 1.0E-8)
std::cerr
<< "Validation failed on sum. Error " << errSum
<< std::endl << std::setprecision(15)
<< "Sum was " << sum << " but should be " << goldSum
<< std::endl;
}