一、存储器简介
        对于CUDA线程，主要设置的存储器有：局部存储器、共享存储器、全局存储器、只读存储器。         局部存储器：作用域在每个线程的内部，每个线程都有一个自己私有的存储器；         共享存储器：作用域在线程块内，相同块内的线程可访问该存储中的数据，已达到数据共享；         全局存储器：比共享存储器的作用域更大，所有块中的线程都可用访问该存储器中的数据         只读存储器：包括常量存储器和纹理存储器 其中，全局存储器、只读存储器（常量存储器、文理存储器）的生命周期同当前应用程序的当前内核发射器一致！以下是相关存储器的层级结构：

二、共享存储器

        我们可在共享存储器中声明需要共享的变量，用__shared__限定符，可以和__device__限定符一起使用（可选），一个声明为__shared__的变量有以下特点：         一、所属的作用域为当前线程块，也即块内的所有线程都都可访问这一限定符所声明的变量，从而起到了块内线程间的数据共享和通信作用；         二、生命周期和线程保持一致；         用法：                  1、__shared__ float array[10];                        声明了一个大小为10的float数组，块内的线程都共享数据，另外，该数据的大小是固定的，为10
                 2、 extern __shared__ float shared[];                       声明了一个动态的数据，其大小并没有在此指定，通过执行配置的第三个参数Ns指定（表明共享存储器的空间大小，单位字节）      三、利用共享存储器进行矩阵乘法         除了上面提到的，用__shared__声明的变量可在块内线程间起到通信作用和数据共享外，另外一个作用就是可降低对全局存储器的访问频率，从而节省带宽，当然得存在这样的情况：对全局存储器访问频繁，且能够用共享存储器去代替全局存储器。         这里以共享存储器为基础，进行了矩阵的乘法运算，这里进行的是矩阵的分块计算，其算法为：对矩阵进行分块为n*n的大小，每个线程块负责计算一个块大小的数据，而具体到每个线程就只计算对应的数据即可。如下图计算C = A *B
            如果不对矩阵进行分块，就是简单的计算，对于计算矩阵C中的每个元素，都要向矩阵A中读取A.width个元素、向矩阵B中读取B.height个元素，如果每个线程负责计算一个对应的数据，由于并发，在计算期间对全局存储器的访问是多么的繁忙，每个线程都要读取A.width + B.height个元素！如果分块计算，那么对于计算一个块内数据，需要向矩阵A读取B.width / BLOCK_SIZE次，向矩阵B读取A.height / BLOCK_SIZE次。通过增加次数来降低对全局存储器的访问频率，节省带宽。     

四、代码实现

#include #include #include #include using namespace std; //下面定义矩阵分块的信息，用于计算矩阵C，C = A * B //注意为了处理简单设置A和B的行数和列数相等 const int g_blockSize = 16; //定义每块的大小 const int g_numBlock = 100; //块数 const int g_matrixSize = g_blockSize * g_numBlock; //矩阵的大小 struct Matrix { int width; //矩阵宽度 int height; //矩阵高度 int stride; //每行的步长，方便操作子矩阵,非子矩阵：width = stride float * pElements; //数据元素的首地址 }; //创建一个矩阵matrix，其大小为size*size的方阵 //随机数初始化 void InitHostMatrix(Matrix & matrix,int size) { matrix.width = size; matrix.height = size; matrix.stride = size; matrix.pElements = new float[matrix.stride * matrix.height]; int n = matrix.stride * matrix.height; for (int i = 0; i < n; ++i) { matrix.pElements[i]= 1.f * rand() / RAND_MAX; } } //销毁一个host矩阵 void DestroyHostMatrix(Matrix & matrix) { delete[] matrix.pElements; matrix.pElements = nullptr; } //用一个host矩阵创初始化device矩阵 //参数copyContent决定了是否将host矩阵的内容拷贝到device矩阵中,defalt:true void InitDeviceMatrix(Matrix& deviceMatrix, const Matrix& hostMatrix,bool copyContent = true) { assert(hostMatrix.width > 0 && hostMatrix.height > 0 && hostMatrix.stride && hostMatrix.pElements); deviceMatrix.width = hostMatrix.width; deviceMatrix.height = hostMatrix.height; deviceMatrix.stride = hostMatrix.stride; int totalSize = sizeof(float) * hostMatrix.stride * hostMatrix.height; cudaError err = cudaMalloc(&deviceMatrix.pElements, totalSize); if (err != cudaSuccess) { std::cerr << "call cudaMalloc fail!" << std::endl; } if (copyContent) { err = cudaMemcpy(deviceMatrix.pElements, hostMatrix.pElements, totalSize, cudaMemcpyHostToDevice); if (err != cudaSuccess) { std::cerr << "call cudaMemcpy fail!" << std::endl; } } } //销毁一个Device矩阵 void DestroyDeviceMatrix(Matrix & deviceMatrix) { if (deviceMatrix.pElements) { cudaFree(deviceMatrix.pElements); deviceMatrix.pElements = nullptr; } } //获取矩阵的子矩阵,根据blockSize来划分块 //声明为__device,表明该函数由设备负责调用 //这里的row、col单位为块，并不是以具体的元素为单位 __device__ Matrix GetSubMatrix(const Matrix matrix, int row, int col) { Matrix resultMatrix; resultMatrix.width = g_blockSize; resultMatrix.height = g_blockSize; resultMatrix.stride = matrix.stride; //！！！这里设置为父矩阵的步长，以便用GetElement函数能够获取到正确的元素 resultMatrix.pElements = matrix.pElements + row * g_blockSize * matrix.stride + col * g_blockSize; return resultMatrix; } //获取矩阵中指定行、列的元素 //这里row、col单位为元素个数，不为块 __device__ float GetElement(const Matrix matrix,int row,int col) { return *(matrix.pElements + matrix.stride * row + col); } //设置矩阵指定行、列的元素 __device__ void SetElment(Matrix matrix, int row, int col,float value) { *(matrix.pElements + matrix.stride * row + col) = value; } //计算matC,matC = matA * matB //每个线程负责计算matC中的一个元素，唯一 __global__ void MatMulKernel(const Matrix matA, const Matrix matB, Matrix matC) { //当前计算矩阵块编号 const int blockRow = blockIdx.y; //注意x、y是要颠倒的 const int blockCol = blockIdx.x; const int threadRow = threadIdx.y; const int threadCol = threadIdx.x; //定义共享子矩阵,每个线程负责从matC中复制数据 __shared__ float sharedSubA[g_blockSize][g_blockSize]; __shared__ float sharedSubB[g_blockSize][g_blockSize]; Matrix subC = GetSubMatrix(matC, blockRow, blockCol); float result = 0.f; //矩阵matA、matB都被分成了g_numBlock*g_numBlock块 for (int b = 0; b < g_numBlock; ++b) { //1、获取当前块的数据 //根据当前块索引，获取matA、matB对应的子矩阵,并将其内容复制到共享存储器中（每个线程只复制一个数据元素） Matrix subA = GetSubMatrix(matA, blockRow, blockCol); Matrix subB = GetSubMatrix(matB, blockRow, blockCol); //将数据复制到共享存储器中 sharedSubA[threadRow][threadCol] = GetElement(subA, threadRow, threadCol); sharedSubB[threadRow][threadCol] = GetElement(subB, threadRow, threadCol); //同步线程 __syncthreads(); //2、计算 for (int e = 0; e < g_blockSize; ++e) { result += sharedSubA[threadRow][e] * sharedSubB[e][threadCol]; //累计结果 } __syncthreads(); } //将结果写入矩阵C对应的子矩阵中 SetElment(subC, threadRow, threadCol, result); } int main() { //定义主机上的矩阵 Matrix A, B, C; InitHostMatrix(A,g_matrixSize); //初始化host矩阵，行列分别为g_matrixSize InitHostMatrix(B,g_matrixSize); InitHostMatrix(C,g_matrixSize); Matrix dA, dB, dC; InitDeviceMatrix(dA, A); InitDeviceMatrix(dB, B); InitDeviceMatrix(dC, C, false); //不拷贝C的内容到dC中 dim3 gridSize(g_numBlock, g_numBlock); //这里预先指定了g_matrixSize = g_numBlock * g_blockSize dim3 blockSize(g_blockSize, g_blockSize); MatMulKernel << > >(dA, dB, dC); //MatMulKernel << > >(dA, dB, dC); cudaError err = cudaGetLastError(); if (err != cudaSuccess) { std::cerr << "call MatMulKernel!" << std::endl; } //将计算结构从device拷贝到host err = cudaMemcpy(C.pElements, dC.pElements, sizeof(float)*dC.stride * dC.height, cudaMemcpyDeviceToHost); if (err != cudaSuccess) { std::cerr << "call cudaMemcpy fail!" << std::endl; } //销毁矩阵 DestroyHostMatrix(A); DestroyHostMatrix(B); DestroyHostMatrix(C); DestroyDeviceMatrix(dA); DestroyDeviceMatrix(dB); DestroyDeviceMatrix(dC); cin.get(); return 0; }