cuda学习之使用常量内存实现一维模板

文章目录

前言

一、参考链接

二、常量内存

1.概念相关

2.一维模板

3.数据填充

总结

---

前言

学习cuda的教程的时候在常量内存突然看到了一个一维卷积，里面的数据填充的理解消耗了一部分，在这里做一个学习记录。

---

一、参考链接

书：cuda c编程权威指南

bgm^-^：别恐惧 - 黄英华 - 单曲 - 网易云音乐

二、常量内存

1.概念相关

常量内存是一种专用的内存，它用于只读数据和统一访问线程束中线程的数据。常量内存对内核代码而言是只读的，但它对主机而言既是可读又是可写的。

在常量内存中，如果线程束中的所有线程都访问相同的位置，那么这个访问模式就是最优的。如果线程束中的线程访问不同的地址，则访问就需要串行。因此，一个常量内存读取的成本与线程束中线程读取唯一地址的数量呈线性关系。

常量内存的声明方式：

__constant__

2.一维模板

在一维中，在位置x 周围的九点模板会给这些位置上的值应用一些函数：

是实变量函数f在点x 上一阶导数的第八阶中心差分公式

代码如下（示例）：

__global__ void stencil_1d(float * in,float * out)
{
    //1.
    __shared__ float smem[BDIM+2*TEMP_RADIO_SIZE];
    //2.
    int idx=threadIdx.x+blockDim.x*blockIdx.x;
    //3.
    int sidx=threadIdx.x+TEMP_RADIO_SIZE;
    smem[sidx]=in[idx];
    //4.
    if (threadIdx.x<TEMP_RADIO_SIZE)
    {
        if(idx>TEMP_RADIO_SIZE)
            smem[sidx-TEMP_RADIO_SIZE]=in[idx-TEMP_RADIO_SIZE];
        if(idx<gridDim.x*blockDim.x-BDIM)
            smem[sidx+BDIM]=in[idx+BDIM];

    }
    __syncthreads();
    //5.
    if (idx<TEMP_RADIO_SIZE||idx>=gridDim.x*blockDim.x-TEMP_RADIO_SIZE)
        return;
    float temp=.0f;
    //6.
    #pragma unroll
    for(int i=1;i<=TEMP_RADIO_SIZE;i++)
    {
        temp+=coef[i-1]*(smem[sidx+i]-smem[sidx-i]);
    }
    out[idx]=temp;
    //printf("%d:GPU :%lf,\n",idx,temp);
}

main函数：

int main(int argc,char** argv)
{
    printf("strating...\n");
    initDevice(0);
    int dimx=BDIM;
    unsigned int nxy=1<<16;
    int nBytes=nxy*sizeof(float);

    //Malloc
    float* in_host=(float*)malloc(nBytes);
    float* out_gpu=(float*)malloc(nBytes);
    float* out_cpu=(float*)malloc(nBytes);
    memset(out_cpu,0,nBytes);
    initialData(in_host,nxy);

    //cudaMalloc
    float *in_dev=NULL;
    float *out_dev=NULL;

    initialData(in_host,nxy);
    float templ_[]={-1.0,-2.0,2.0,1.0};
    CHECK(cudaMemcpyToSymbol(coef,templ_,TEMP_RADIO_SIZE*sizeof(float)));

    CHECK(cudaMalloc((void**)&in_dev,nBytes));
    CHECK(cudaMalloc((void**)&out_dev,nBytes));
    CHECK(cudaMemcpy(in_dev,in_host,nBytes,cudaMemcpyHostToDevice));
    CHECK(cudaMemset(out_dev,0,nBytes));

    // cpu compute
    double iStart=cpuSecond();
    convolution(in_host,out_cpu,templ_,nxy);
    double iElaps=cpuSecond()-iStart;
    //printf("CPU Execution Time elapsed %f sec\n",iElaps);

    // stencil 1d
    dim3 block(dimx);
    dim3 grid((nxy-1)/block.x+1);
    stencil_1d<<<grid,block>>>(in_dev,out_dev);
    CHECK(cudaDeviceSynchronize());
    iElaps=cpuSecond()-iStart;
    printf("stencil_1d Time elapsed %f sec\n",iElaps);
    CHECK(cudaMemcpy(out_gpu,out_dev,nBytes,cudaMemcpyDeviceToHost));
    checkResult(out_cpu,out_gpu,nxy);
    CHECK(cudaMemset(out_dev,0,nBytes));

    cudaFree(dcoef_ro);
    cudaFree(in_dev);
    cudaFree(out_dev);
    free(out_gpu);
    free(out_cpu);
    free(in_host);
    cudaDeviceReset();
    return 0;
}

可以看到这里定义的grid是一维2048个block；

block定义的是一维32个thread；

一共2^16个数据；

依旧是每个线程负责一个数据；

seme定义了40位[0,39]的共享内存；

sidx范围为[4,35]；

idx指的是全局索引。

3.数据填充

如果模板大小是9，那么我们输出的前4个数据是没办法计算的因为要使用第-1，-2，-3，-4位置的数据，最后4个数据也是不能计算的，因为他要使用 n+1,n+2,n+3,n+4 的数据，这些数据也是没有的，为了保证计算过程中访问不会越界，我们把输入数据两端进行扩充，也就是把上面虽然没有，但是要用的数据填充到输入数据中，当我们我们要处理的是中间的某段数据的时候，那么填充位的数据就来自前面的块对应的输入数据，或者后面线程块对应的输入数据。

关键代码如下：

if (threadIdx.x<TEMP_RADIO_SIZE)
    {
        if(idx>TEMP_RADIO_SIZE)
            smem[sidx-TEMP_RADIO_SIZE]=in[idx-TEMP_RADIO_SIZE];
        if(idx<gridDim.x*blockDim.x-BDIM)
            smem[sidx+BDIM]=in[idx+BDIM];
    }

1. 这里的条件语句指的是“从全局内存中读取数据到共享内存中时，前四个线程负责从左侧和右侧的光环中读取数据到共享内存中”，也就是只让前四个线程负责读取数据这件事

2. 里面的两个if语句的判断条件只是为了保证他们每个下面要执行的命令不会出现错误，因为索引一般都是非负整数，不能出现负值

3. 第一个命令行：

smem[sidx-TEMP_RADIO_SIZE]=in[idx-TEMP_RADIO_SIZE];

是为了从前面的数据块中取4个数据，考虑的是邻居元素的数量是4个，除了这个数据块的32个数据以外，取到左相邻的4个数据。

4. 第二个命令行：

smem[sidx+BDIM]=in[idx+BDIM];

因为条件语句是说的只让前四个thread负责填充数据，所以需要加上一整个block的长度才可以取到后面四个的值。除了这个数据块的32个数据以外，取到右相邻的4个数据。

注意这里数据填充其实是取前后数据块本身就有的数据，而不是给空位赋值。

5. 整个数据最前和最尾的数据，无法做计算的会被return

if (idx<TEMP_RADIO_SIZE||idx>=gridDim.x*blockDim.x-TEMP_RADIO_SIZE)
        return;

---

总结

重点在于理解数据填充方法