I recently implemented Successive Over Relaxation using Cuda as a part of my course project and was curious to know how I can make the code more efficient.
I'm using Red/Black SOR scheme which is popular version of parallel SOR. To update the value at a grid point, I need to access 4 adjacent values which are generally not contiguous in the memory and hence the GPU is not able to use memory coalescing feature of cuda which makes the algorithm memory bound. To overcome this I have reordered the storage of odd and even grid points in a separate vectors as mentioned in the paper here: http://link.springer.com/chapter/10.1007%2F978-3-642-31464-3_60 which takes advantage of memory coalescing.
Implementing the technique mentioned in the paper I see a speedup of around ~15 (compared to ~10 with traditional Red/Black SOR) when compared with 1 CPU. How much performance speed up I should be expecting? I feel speedup of ~15 is quite low?
I'm calling the kernel functions to update odd values and then even values sequentially from host CPU in following way :
for(size_t it = 0; it < itmax; it++){
cuda_odd_update <<<dimGrid_odd,dimBlock>>> (odd,even);
cuda_even_update <<<dimGrid_even,dimBlock>>> (odd,even);
}
where odd
and even
are the pointers to reordered vectors in device memory. I'm using this approach as I'm not aware of any way to synchronize all the threads in the device. Is there any better way to implement the execution loop?
The code below implements the technique mentioned in the paper and not a traditional Red/Black SOR:
__global__ void cuda_odd_update(double* odd,double* even){
size_t tx = blockIdx.x*blockDim.x + threadIdx.x;
size_t ty = blockIdx.y*blockDim.y + threadIdx.y;
size_t odd_index = tx*height_odd+ty;
size_t even_index = tx*height_even+ty;
if (( (ty == 0 && tx%2 == 0) || (ty == height_odd-1 && tx%2 == 1) || (ty > 0 && ty < height_odd-1) ) && (tx > 0 && tx < width-1 && ty < height_odd) ){
odd[odd_index] = (1.0 - omega)*odd[odd_index] + omega/(2*(1+beta))
* ( (tx%2)*even[even_index]+(1-tx%2)*even[even_index+1]
+ (tx%2)*even[even_index-1]+(1-tx%2)*even[even_index]
+ beta * ( even[even_index-height_even]
+ even[even_index+height_even] )) ;
}
}
__global__ void cuda_even_update(double* odd,double* even){
size_t tx = blockIdx.x*blockDim.x + threadIdx.x;
size_t ty = blockIdx.y*blockDim.y + threadIdx.y;
size_t odd_index = tx*height_odd+ty;
size_t even_index = tx*height_even+ty;
if(( (ty == 0 & tx%2 == 1) || (ty == height_even-1 && tx%2 == 0) || (ty > 0 && ty < height_even-1) ) && ( tx > 0 && tx < width-1 && ty < height_even)){
even[even_index] = (1.0 - omega)*even[even_index] + omega/(2*(1+beta))
* ( (1-tx%2)*odd[odd_index] + (tx%2)*odd[odd_index+1]
+ (1-tx%2)*odd[odd_index-1] + (tx%2)*odd[odd_index]
+ beta * ( odd[odd_index-height_odd]
+ odd[odd_index+height_odd] )) ;
}
}
I'm currently using Intel Xeon CPU E5-2660-2.20GHz (8 cores) and NVIDIA TESLA M2090.
The nvprof
gives the following results (for 1000 iterations):
CPU time : 62.69
==27665== NVPROF is profiling process 27665, command: ./sor
GPU time : 3.73
End!
==27665== Profiling application: ./sor
==27665== Profiling result:
Time(%) Time Calls Avg Min Max Name
49.81% 1.83752s 1000 1.8375ms 1.8333ms 1.8422ms cuda_even_update(double*, double*)
49.67% 1.83226s 1000 1.8323ms 1.8290ms 1.8364ms cuda_odd_update(double*, double*)
0.26% 9.6281ms 2 4.8141ms 4.6306ms 4.9975ms [CUDA memcpy DtoH]
0.26% 9.4569ms 2 4.7284ms 4.7065ms 4.7504ms [CUDA memcpy HtoD]