1. 定义
a. 在变量前加上 __shared__, 此变量存储于共享内存中. 1-D或2-D都可以.
b. CUDA在每一个block中都copy了一份此变量.
c. 同一个block中的threads共享此内存, 但是无法读写其他block的拷贝.
d. 共享内存的访问延迟远低于常见的buffer.
e. 需要同步机制(synchronization)来协调threads的读写行为. 防止计算没有完成就执行后续指令.
f. 所有的thread都应该做一样的事. 如果做的事不一样(有if判断), 有些thread可能永远无法执行到同步阶段.
2. 使用共享内存求内积
#include "../common/book.h"#define imin(a,b) (a<b?a:b)const int N = 33 * 1024;
const int threadPerBlock = 256;
const int blockPerGrid = imin( 32, (N+threadPerBlock-1) / threadPerBlock );__global__ void dot( float *a, float *b, float *c)
{//共享内存, 每个block都有一份拷贝__shared__ float cache[threadPerBlock];// thread的索引int tid = threadIdx.x + blockIdx.x * blockDim.x;// 共享内存的索引,每个block都有cache, 故只用threadIdx.x即可int cacheIdx = threadIdx.x;float temp = 0;while(tid<N){//当前tid的thread负责把tid,和tid间隔threadIdx总量整数倍的向量做乘-加操作.temp += a[tid] * b[tid];tid += blockDim.x * gridDim.x;}// 完成求和之后,当前thread把和放在对应的cache中cache[cacheIdx] = temp;// 在当前block内做同步操作, 等所有thread都完成乘-加运算之后才能做reduction.__syncthreads();//reduction, 向量缩减.//缩减后的结果在cache[0]里.int i = blockDim.x/2;while (i!=0){if (cacheIdx<i){cache[cacheIdx] += cache[cacheIdx + i];}//同步, 等所有thread都完成了当次缩减了才能做下一次的缩减.//书上说: 同步不能放在if里面, 否则报错.//经过试验没有报错, 结果正确.__syncthreads();i /= 2;}// 一个block输出一个值,即cache[0]. 所以c的长度和block数量相同.// 限制cacheIdx == 0是为了只做一次赋值操作,节省时间.if (cacheIdx == 0){c[blockIdx.x] = cache[0];}// 没有做剩下的累加操作是因为在CPU上做小批量的累加更加有效.
}int main(void)
{float *a, *b, c, *partial_c;float *dev_a, *dev_b, *dev_partial_c;//分配CPU端的内存a = (float *)malloc( N*sizeof(float) );b = (float *)malloc( N*sizeof(float) );partial_c = (float *)malloc( blockPerGrid*sizeof(float));//分配GPU端的内存HANDLE_ERROR( cudaMalloc( (void**)&dev_a, N*sizeof(float)));HANDLE_ERROR( cudaMalloc( (void**)&dev_b, N*sizeof(float)));HANDLE_ERROR( cudaMalloc( (void**)&dev_partial_c, blockPerGrid*sizeof(float)));//将主机内存填入数据for (int i=0; i<N; i++){a[i] = i;b[i] = i*2;}//将向量a和b拷入GPUHANDLE_ERROR( cudaMemcpy( dev_a, a, N*sizeof(float), cudaMemcpyHostToDevice));HANDLE_ERROR( cudaMemcpy( dev_b, b, N*sizeof(float), cudaMemcpyHostToDevice));//GPU上做点积运算dot<<<blockPerGrid, threadPerBlock>>>(dev_a, dev_b, dev_partial_c);//将向量拷入主机HANDLE_ERROR( cudaMemcpy( partial_c, dev_partial_c, blockPerGrid*sizeof(float), cudaMemcpyDeviceToHost));//剩余CPU运算, 求累加和c = 0;for (int i=0; i<blockPerGrid; i++){c += partial_c[i];}//验证结果是否正确
#define sum_square(x) (x*(x+1)*(2*x+1)/6)printf( "Does GPU value %.6g = %.6g?\n",c,2 * sum_square( (float)(N-1) ) );//释放内存cudaFree( dev_a );cudaFree( dev_b );cudaFree( dev_partial_c);free( a );free( b );free( partial_c);return 0;
}3. 不使用同步指令的问题
不使用同步指令, 可能共享内存的计算还没有完成就去执行下面的步骤了. 会造成错误.
#include "cuda.h"
#include "../common/book.h"
#include "../common/cpu_bitmap.h"#define DIM 1024
#define PI 3.1415926535897932f__global__ void kernel( unsigned char *ptr)
{int x = threadIdx.x + blockIdx.x * blockDim.x;int y = threadIdx.y + blockIdx.y * blockDim.y;int offset = x + y * blockDim.x * gridDim.x;__shared__ float shared[16][16];const float period = 128.0f;shared[threadIdx.x][threadIdx.y] = 255 * (sinf(x*2.0f*PI/period) + 1.0f)*(sinf(y*2.0f*PI/period) + 1.0f)/4.0f;// 必须加上同步指令.否则有的计算没有完成.__syncthreads();ptr[offset*4 + 0] = 0;ptr[offset*4 + 1] = shared[15-threadIdx.x][15-threadIdx.y];ptr[offset*4 + 2] = 0;ptr[offset*4 + 3] =255;
}
int main(void)
{CPUBitmap bitmap(DIM, DIM);unsigned char *dev_bitmap;HANDLE_ERROR( cudaMalloc( (void**)&dev_bitmap,bitmap.image_size() ) );dim3 grids(DIM/16, DIM/16);dim3 threads(16, 16);kernel<<<grids, threads>>>( dev_bitmap );HANDLE_ERROR( cudaMemcpy( bitmap.get_ptr(), dev_bitmap,bitmap.image_size(),cudaMemcpyDeviceToHost ) );bitmap.display_and_exit();cudaFree( dev_bitmap );}有同步指令的输出(部分图):
没有同步指令的输出(部分图, 错误):
