使用cublasGemmEx接口在GPU中进行int8型矩阵乘法运算

最新推荐文章于 2024-05-15 15:33:52 发布

沙拉OK

最新推荐文章于 2024-05-15 15:33:52 发布

阅读量3.5k

点赞数 3

分类专栏： GPU比对文章标签： cublas库 GPU比对矩阵乘法

本文链接：https://blog.csdn.net/Hearz/article/details/80636873

版权

GPU比对专栏收录该内容

1 篇文章 0 订阅

订阅专栏

在一些诸如特征比对等等的应用场景中，特征值其实就是一个个的浮点数(float)组成的矩阵，这些浮点数每一个都是4个字节(32位），然后对这两个矩阵进行乘法计算，A * B = C，得出的C矩阵就是能够代表两个特征矩阵相似性的值。毫无疑问矩阵乘法可以在CPU中进行计算，但是，在一些需要实时性的场合下，CPU的运算速度可能就不能满足要求，而通过GPU去做矩阵乘法反而再合适不过，GPU处理浮点数运算的能力非常强，这里记录两种接口来使用GPU计算矩阵乘法。

首先，两个矩阵相乘有意义必须满足：矩阵A为m行n列，矩阵B为n行k列，则矩阵C=AB即C为m行k列。使用cuda的runtime库和cublas库在GPU中进行计算。现在CPU中使用随机浮点数产生两个矩阵h_A和h_B，再使用cublasSetVector函数或cudaMemcpy函数把h_A和h_B从CPU内存中复制到GPU中相应的d_A和d_B中，如果是float型，可以直接使用cublasSegmm函数进行矩阵乘法运算，结果置于d_C中，若为int8型，则可使用cublasGemmEx函数，再使用cublasGetVector函数或cudaMemcpy函数将d_C从GPU内存中复制到CPU内存中。

cublasSegmm函数做FP32型矩阵乘法实现：

   A:M行N列矩阵  B:N行M列矩阵  C:M行M列矩阵

#include "cuda_runtime.h"
#include "cublas_v2.h"
#include <iostream>

using namespace std;

#define M 5
#define N 3

  int main() 
  {   
      // 定义状态变量，矩阵乘法接口的返回值
      cublasStatus_t status;
  
      // 在 CPU内存 中为将要计算的矩阵开辟空间
      float *h_A = (float*)malloc (N*M*sizeof(float));
      float *h_B = (float*)malloc (N*M*sizeof(float));
      
      // 在 CPU内存 中为将要存放运算结果的矩阵开辟空间
      float *h_C = (float*)malloc (M*M*sizeof(float));
  
      // 为待运算矩阵的元素赋予 0-10 范围内的随机数，实际使用时，A矩阵需要做转置，B矩阵不需要转      置，可以手动转置也可以算法内转置
      for (int i=0; i<N*M; i++) 
      {
          h_A[i] = (float)(rand()%10+1);
          h_B[i] = (float)(rand()%10+1);
      }
      
      // 打印待测试的矩阵
      cout << "矩阵 A :" << endl;
      for (int i = 0; i < N * M; i++){
          cout << h_A[i] << " ";
          if ((i + 1) % N == 0) 
              cout << endl;
      }
      cout << endl;
      cout << "矩阵 B :" << endl;
      for (int i = 0; i < N * M; i++){
         cout << h_B[i] << " ";
          if ((i + 1) % M == 0) 
              cout << endl;
      }
      cout << endl;
      
      /*
      ** GPU 计算矩阵相乘
      */
  
      // 创建并初始化 CUBLAS 库对象
      cublasHandle_t handle;
      status = cublasCreate(&handle);
      
      if (status != CUBLAS_STATUS_SUCCESS)
      {
          if (status == CUBLAS_STATUS_NOT_INITIALIZED) {
              cout << "CUBLAS 对象实例化出错" << endl;
          }
          getchar ();
          return EXIT_FAILURE;
      }
  
      float *d_A, *d_B, *d_C;
      // 在 显存 中为将要计算的矩阵开辟空间
      cudaMalloc (
          (void**)&d_A,    // 指向开辟的空间的指针
          N*M * sizeof(float)    //　需要开辟空间的字节数
      );
      cudaMalloc (
          (void**)&d_B,    
          N*M * sizeof(float)    
      );
  
      // 在 显存 中为将要存放运算结果的矩阵开辟空间
      cudaMalloc (
          (void**)&d_C,
          M*M * sizeof(float)    
      );
  
      // 将矩阵数据传递进 显存 中已经开辟好了的空间
      cublasSetVector (
          N*M,    // 要存入显存的元素个数
          sizeof(float),    // 每个元素大小
          h_A,    // 主机端起始地址
          1,    // 连续元素之间的存储间隔
          d_A,    // GPU 端起始地址
          1    // 连续元素之间的存储间隔
      );
     //注意：当矩阵过大时，使用cudaMemcpy是更好地选择：
     //cudaMemcpy(d_A, h_A, sizeof(float)*N*M, cudaMemcpyHostToDevice);

     cublasSetVector (
          N*M, 
          sizeof(float), 
          h_B, 
          1, 
          d_B, 
          1
      );
     //cudaMemcpy(d_B, h_B, sizeof(float)*N*M, cudaMemcpyHostToDevice);
     // 同步函数
     cudaThreadSynchronize();
 
     // 传递进矩阵相乘函数中的参数，具体含义请参考函数手册。
     float a = 1; 
     float b = 0;
     // 矩阵相乘。该函数必然将数组解析成列优先数组
     cublasSgemm (
         handle,    // blas 库对象 
         CUBLAS_OP_T,    // 矩阵 A 属性参数
         CUBLAS_OP_T,    // 矩阵 B 属性参数
         M,    // A, C 的行数 
         M,    // B, C 的列数
         N,    // A 的列数和 B 的行数
         &a,    // 运算式的 α 值
         d_A,    // A 在显存中的地址
         N,    // lda
         d_B,    // B 在显存中的地址
         M,    // ldb
         &b,    // 运算式的 β 值
         d_C,    // C 在显存中的地址(结果矩阵)
         M    // ldc
     );
     
     // 同步函数
     cudaThreadSynchronize();
 
     // 从 显存 中取出运算结果至 内存中去
     cublasGetVector (
         M*M,    //  要取出元素的个数
         sizeof(float),    // 每个元素大小
         d_C,    // GPU 端起始地址
         1,    // 连续元素之间的存储间隔
         h_C,    // 主机端起始地址
         1    // 连续元素之间的存储间隔
     );

     //或使用cudaMemcpy(h_C, d_C, sizeof(float)*M*M, cudaMemcpyDeviceToHost);
     // 打印运算结果
     cout << "计算结果的转置 ( (A*B)的转置 )：" << endl;
 
     for (int i=0;i<M*M; i++)
     {
          cout << h_C[i] << " ";
          if ((i+1)%M == 0) cout << endl;
     }
     
     // 清理掉使用过的内存
     free (h_A);
     free (h_B);
     free (h_C);

     try
     {
         cudaFree (d_A);
         cudaFree (d_B);
         cudaFree (d_C);
     }
     catch(...)
     {
          cout << "cudaFree Error!" << endl;
          // 释放 CUBLAS 库对象
     }

     cublasDestroy (handle);
     return 0;
}

但是GPU的内存一般比较有限，实际项目中，为了能够做更大的矩阵乘法，使用FP32型数据占用空间会很大，所以可以压缩为int8型数据进行矩阵乘法运算，这时使用的接口模型为：

cublasStatus_t cublasGemmEx(cublasHandle_t handle,      //句柄
                            cublasOperation_t transa,   //a矩阵转置选项
                            cublasOperation_t transb,   //b矩阵转置选项
                            int m,                      //a矩阵行数
                            int n,                      //b矩阵列数
                            int k,                      //a矩阵列数兼b矩阵行数
                            const void *alpha,          //乘法因子alpha
                            const void *A,              //A矩阵
                            cudaDataType_t Atype,       //A矩阵的数据类型
                            int lda,                    //A矩阵的列数
                            const void *B,              //B矩阵
                            cudaDataType_t Btype,       //B矩阵的数据类型
                            int ldb,                    //B矩阵的行数
                            const void *beta,           //乘法因子beta
                            void *C,                    //C结果矩阵
                            cudaDataType_t Ctype,       //C结果矩阵数据类型
                            int ldc,                    //C矩阵的行数
                            cudaDataType_t computeType, //计算模式
                            cublasGemmAlgo_t algo)

该矩阵接口和cublasSgemm几乎相同，但多了几个参数：cudaDataType_t Atype A矩阵数据类型、cudaDataType_t Btype B矩阵数据类型、cudaDataType_t Ctype C矩阵数据类型以及cudaDataType_t computeType 计算模式和cublasGemmAlgo_t参数。首先来说，这个接口可以完成的运算表达式如下：

C=αop(A)op(B)+βC

其中，α和β为标量，A、B、C分别是m×k、k×n、m×n的矩阵，当α = 1, β = 0时，即为：

C=op(A)op(B)

也就是做矩阵乘法。α和β分别对应cublasGemmEx中的const void *alpha、const void *beta，这一点与前面cublasSegmm做FP32类型的接口是一样的。

在接口中，cublasOperation_t transa矩阵转置选项有如下三种：

transa == CUBLAS_OP_N表示A矩阵，无任何转置。

transa == CUBLAS_OP_T表示A^T，A的转置。

transa == CUBLAS_OP_C表示A^H，A的共轭转置。

同理，B矩阵也是一样。

在接口中，新增的cudaDataType_t Atype A矩阵数据类型、cudaDataType_t Btype B矩阵数据类型、cudaDataType_t Ctype C矩阵数据类型以及cudaDataType_t computeType这几个参数的取值情况如下表所示：

computeType	Atype/Btype	Ctype
CUDA_R_16F	CUDA_R_16F	CUDA_R_16F
CUDA_R_32I	CUDA_R_8I	CUDA_R_32I
CUDA_R_32F	CUDA_R_16F	CUDA_R_16F
	CUDA_R_8I	CUDA_R_32F
	CUDA_R_16F	CUDA_R_32F
	CUDA_R_32F	CUDA_R_32F
CUDA_R_64F	CUDA_R_64F	CUDA_R_64F
CUDA_C_32F	CUDA_C_8I	CUDA_C_32F
CUDA_C_32F	CUDA_C_32F	CUDA_C_32F
CUDA_C_64F	CUDA_C_64F	CUDA_C_64F

我们这里选择的int8型矩阵乘法，最好选择Atype/Btype为CUDA_R_8I型，计算模式CUDA_R_32I，乘出的结果Ctype为CUDA_R_32I类型。最后的cublasGemmAlgo_t参数可以选择的值如下表所示：

CublasGemmAlgo_t	含义
CUBLAS_GEMM_DEFAULT	Apply Heuristics to select the GEMM algorithm
CUBLAS_GEMM_ALGO0 to CUBLAS_GEMM_ALGO23	Explicitly choose an algorithm
CUBLAS_GEMM_DEFAULT_TENSOR_OP	Apply Heuristics to select the GEMM algorithm while allowing the use of Tensor Core operations if possible
CUBLAS_GEMM_ALGO0_TENSOR_OP to CUBLAS_GEMM_ALGO15_TENSOR_OP	Explicitly choose a GEMM algorithm allowing it to use Tensor Core operations if possible, otherwise falls back to cublas<t>gemmBatched based on computeType

这里可以选择CUBLAS_GEMM_ALGO0~7。

cublasGemmEX函数做int8型矩阵乘法实现：

#include "cuda_runtime.h"
#include "cublas_v2.h"
#include <iostream>
#include <math.h>

using namespace std;

#define M 5
#define N 3
#define BYTE 128

int main() 
{   
      // 定义状态变量，矩阵乘法接口的返回值
      cublasStatus_t status;
  
      // 在 CPU内存 中为将要计算的矩阵开辟空间
      char *h_A = (char *)malloc (N * M * sizeof(char));
      char *h_B = (char *)malloc (N * M * sizeof(char));
      
      // 在 CPU内存 中为将要存放运算结果的矩阵开辟空间
      int *h_C = (int *)malloc ( M * M * sizeof(int));
  
      float *f_A = (float *)malloc(N * M * sizeof(float));
      float *f_B = (float *)malloc(N * M * sizeof(float));
      // 为待运算矩阵的元素赋予 0-10 范围内的随机数，实际使用时，A矩阵需要做转置，B矩阵不需要转置，可以手动转置也可以算法内转置
      for (int i = 0; i < N * M; i++) 
      {
          f_A[i] = (float)(rand()%10+1);
          f_B[i] = (float)(rand()%10+1);
      }

      //FP32转int8型，有精度损失
      for(int i = 0; i < N * M; i++)
      {
          h_A[i] = (char)(round(f_A[i] * BYTE)); 
          h_B[i] = (char)(round(f_B[i] * BYTE)); 
      }
      
      // 打印待测试的矩阵
      /*cout << "矩阵 A :" << endl;
      for (int i = 0; i < N * M; i++){
          cout << (int)h_A[i] << " ";
          if ((i + 1) % N == 0) 
              cout << endl;
      }
      cout << endl;
      cout << "矩阵 B :" << endl;
      for (int i = 0; i < N * M; i++){
         cout << (int)h_B[i] << " ";
          if ((i + 1) % M == 0) 
              cout << endl;
      }
      cout << endl;
      */
      /*
      ** GPU 计算矩阵相乘
      */
  
      // 创建并初始化 CUBLAS 库对象
      cublasHandle_t handle;
      status = cublasCreate(&handle);
      
      if (status != CUBLAS_STATUS_SUCCESS)
      {
          if (status == CUBLAS_STATUS_NOT_INITIALIZED) {
              cout << "CUBLAS 对象实例化出错" << endl;
          }
          getchar ();
          return EXIT_FAILURE;
      }
  
      char *d_A, *d_B;
      int *d_C;
      // 在 显存 中为将要计算的矩阵开辟空间
      cudaMalloc (
          (void **)&d_A,    // 指向开辟的空间的指针
          N * M * sizeof(char)    //　需要开辟空间的字节数
      );
      cudaMalloc (
          (void **)&d_B,    
          N * M * sizeof(char)    
      );
  
      // 在 显存 中为将要存放运算结果的矩阵开辟空间
      cudaMalloc (
          (void **)&d_C,
          M * M * sizeof(int)    
      );
  
      // 将矩阵数据传递进 显存 中已经开辟好了的空间
      cublasSetVector (
          N * M,    // 要存入显存的元素个数
          sizeof(char),    // 每个元素大小
          h_A,    // 主机端起始地址
          1,    // 连续元素之间的存储间隔
          d_A,    // GPU 端起始地址
          1    // 连续元素之间的存储间隔
      );
     //注意：当矩阵过大时，使用cudaMemcpy是更好地选择：
     //cudaMemcpy(d_A, h_A, sizeof(char)*N*M, cudaMemcpyHostToDevice);

     cublasSetVector (
          N * M, 
          sizeof(char), 
          h_B, 
          1, 
          d_B, 
          1
      );
     //cudaMemcpy(d_B, h_B, sizeof(char)*N*M, cudaMemcpyHostToDevice);
     // 同步函数
     cudaThreadSynchronize();
 
     // 传递进矩阵相乘函数中的参数，具体含义请参考函数手册。
     float a = 1.0; 
     float b = 0;
     // 矩阵相乘。该函数必然将数组解析成列优先数组
     cublasSgemm (
         handle,    // blas 库对象 
         CUBLAS_OP_T,    // 矩阵 A 属性参数
         CUBLAS_OP_T,    // 矩阵 B 属性参数
         M,    // A, C 的行数 
         M,    // B, C 的列数
         N,    // A 的列数和 B 的行数
         &a,    // 运算式的 α 值
         d_A,    // A 在显存中的地址
         N,    // lda
         d_B,    // B 在显存中的地址
         M,    // ldb
         &b,    // 运算式的 β 值
         d_C,    // C 在显存中的地址(结果矩阵)
         M    // ldc
     );
    cublasGemmEx (handle,               //句柄
                   CUBLAS_OP_T,         //矩阵 A 属性参数
                   CUBLAS_OP_T,         //矩阵 B 属性参数
                   M,                   //A, C 的行数 
                   M,                   //B, C 的列数
                   N,                   //A 的列数和 B 的行数
                   &a,                  //运算式的 α 值
                   d_A,                 //A矩阵
                   CUDA_R_8I,           //A矩阵计算模式，int8型
                   N,                   //A矩阵的列数
                   d_B,                 //B矩阵
                   CUDA_R_8I,           //B矩阵计算模式，int8型
                   M,                   //B矩阵的行数
                   &b,                  //乘法因子beta
                   d_C,                 //C结果矩阵
                   CUDA_R_32I,          //C矩阵计算模式，int32型
                   M,                   //C矩阵的行数
                   CUDA_R_32I,          //计算模式，int32模式
                   CUBLAS_GEMM_ALGO0    //算法参数
                   )      
     
     // 同步函数
     cudaThreadSynchronize();
 
     // 从 显存 中取出运算结果至 内存中去
     cublasGetVector (
         M * M,    //  要取出元素的个数
         sizeof(int),    // 每个元素大小
         d_C,    // GPU 端起始地址
         1,    // 连续元素之间的存储间隔
         h_C,    // 主机端起始地址
         1    // 连续元素之间的存储间隔
     );

     //或使用cudaMemcpy(h_C, d_C, sizeof(int)*M*M, cudaMemcpyDeviceToHost);
     // 打印运算结果
     cout << "计算结果的转置 ( (A*B)的转置 )：" << endl;
 
     for (int i = 0;i < M * M; i++)
     {
          cout << h_C[i] << " ";
          if ((i+1)%M == 0) 
            cout << endl;
     }
     
     //注意，这里需要进行归一化操作，乘出来的矩阵需要除以128*128，以还原原来的大小。在此就省略这一步。
     // 清理掉使用过的内存
     free (h_C);
     free (h_B);
     free (h_A);
     free (f_B);
     free (f_A);

     try
     {
         cudaFree (d_A);
         cudaFree (d_B);
         cudaFree (d_C);
     }
     catch(...)
     {
          cout << "cudaFree Error!" << endl;
          // 释放 CUBLAS 库对象
     }

     cublasDestroy (handle);
     return 0;
}

以上代码编译命令可以是：g++ GPURetrievalInt8.cpp -o GPURetrievalDemo -L/usr/local/cuda/lib64 -lcudart -lcuda

当运行有误时，可以打印cublasGemmEx的返回值来确认是什么错误，官方定义的错误码如下表所示：

Error Value	Meaning
CUBLAS_STATUS_SUCCESS	the operation completed successfully
CUBLAS_STATUS_NOT_INITIALIZED	the library was not initialized
CUBLAS_STATUS_ARCH_MISMATCH	cublasGemmEx is only supported for GPU with architecture capabilities equal or greater than 5.0
CUBLAS_STATUS_NOT_SUPPORTED	the combination of the parameters Atype, Btype and Ctype or the algorithm, algo is not supported
CUBLAS_STATUS_INVALID_VALUE	the parameters m,n,k<0
CUBLAS_STATUS_EXECUTION_FAILED	the function failed to launch on the GPU

成功时会返回0，也就是CUBLAS_STATUS_SUCCESS。多数情况会返回CUBLAS_STATUS_INVALID_VALUE，这就是m,n,k参数传入不正确所致。

经测试，GPU矩阵乘法可以达到大约1.6亿次/秒的速度，只要提前将向量存入GPU显存中，待比对向量以1个或n个的形式拷贝到GPU显存中，乘法之后选用合适的排序方法如堆排序等，就可以从大到小将相似度最高到最低排列出来。应用前景广泛。

沙拉OK

关注

3
点赞
踩
4

收藏

觉得还不错? 一键收藏
0
评论
使用cublasGemmEx接口在GPU中进行int8型矩阵乘法运算

在一些诸如特征比对等等的应用场景中，特征值其实就是一个个的浮点数(float)组成的矩阵，这些浮点数每一个都是4个字节(32位），然后对这两个矩阵进行乘法计算，A * B = C，得出的C矩阵就是能够代表两个特征矩阵相似性的值。毫无疑问矩阵乘法可以在CPU中进行计算，但是，在一些需要实时性的场合下，CPU的运算速度可能就不能满足要求，而通过GPU去做矩阵乘法反而再合适不过，GPU处...
复制链接

扫一扫