知识点
上节课问题部分答疑
- 一个block只有1024个线程。
- 最佳参数需要自己尝试
- 每个并行程序输出乱序,是因为输出buff导致。
- 注意判断矩阵数据越界。
存储单元
与CPU(host)交互的内存Global、Constant、Texture Memory
此外线程可与其他的shared memory、 register、local memory等内存单元读写。
内存申请操作
GPU与CPU数据异步运行,且需要提前在GPU上申请数据空间,主要功能函数为cudaMemcpy()
、cudaFree()
。
gpu和cpu异步运行,需要同步函数
矩阵相乘
cpu:
CPU上计算矩阵乘法需要重复计算乘法并累加,计算繁杂,花费较大
void cpu_matrix_mult(int *h_a, int *h_b, int *h_result, int m, int n, int k)
{
for (int i = 0; i < m; ++i)
{
for (int j = 0; j < k; ++j)
{
int tmp = 0;
for (int h = 0; h < n; ++h)
{
tmp += h_a[i * n + h] * h_b[h * k + j];
}
h_result[i * k + j] = tmp;
}
}
}
如果矩阵a,b维度1000,就要作1000^3次tmp += h_a[i * n + h] * h_b[h * k + j];
GPU
利用矩阵索引替代for循环中的i,j
矩阵乘法=每个元素乘积累加的运算。每个无关元素运算互不影响
GPU提高速度:利用空间换取时间。多个核心替代运算步骤
实际上,二维矩阵数据储存在一维的连续空间上,所以可以根据x,y坐标找到二维矩阵中的元素。
index索引计算示例如下图。threadidx(2,0)
为网格单元中所示单元数字, blockidx(1,1)
为所示单元所在的block坐标, blockDim(4,3)
为每个block的大小。
实验实例
矩阵乘法示例
![
对于 sum+=M[row*n_size+i]*N[i*k_size+col];
M是行乘,所有元素内存位置相邻,故是
行数*宽度(起始位置)+i(当行迭代索引)
N是列乘,每个元素内存位置都相差一个行元素的位置,故是
i*列数(起始位置)+col(高度),内存中离散位置
#include <stdio.h>
#include <cuda_runtime.h>
#include "device_launch_parameters.h"
#include <math.h>
#define BLOCK_SIZE 16
__global__ void gpu_matrix_mult(int *M,int *N, int *P,int m_size,int n_size, int k_size){
//计算当前线程坐标
int row=blockIdx.x*blockDim.x+threadIdx.x;
int col=blockIdx.y*blockDim.y+threadIdx.y;
//乘积运算
int sum=0;
if(row<m_size&&col<k_size){
for(int i=0;i<n_size;i++){
sum+=M[row*n_size+i]*N[i*k_size+col];
}
P[row*k_size+col]=sum;
}
}
int main(int argc, char const *argv[])
{
int m=100;
int n=100;
int k=100;
// 申请计算内存
int *h_a, *h_b, *h_c, *h_cc;
cudaMallocHost((void **) &h_a, sizeof(int)*m*n);
cudaMallocHost((void **) &h_b, sizeof(int)*n*k);
cudaMallocHost((void **) &h_c, sizeof(int)*m*k);
cudaMallocHost((void **) &h_cc, sizeof(int)*m*k);
for (int i = 0; i < m; ++i) {
for (int j = 0; j < n; ++j) {
h_a[i * n + j] = rand() % 1024;
}
}
for (int i = 0; i < n; ++i) {
for (int j = 0; j < k; ++j) {
h_b[i * k + j] = rand() % 1024;
}
}
//申请结果内存
int *d_a, *d_b, *d_c;
cudaMalloc((void**)&d_a,sizeof(int)*m*n);
cudaMalloc((void**)&d_b,sizeof(int)*k*n);
cudaMalloc((void**)&d_c,sizeof(int)*m*k);
//转移数据
cudaMemcpy(d_a,h_a,sizeof(int)*m*n,cudaMemcpyHostToDevice);
cudaMemcpy(d_b,h_b,sizeof(int)*k*n,cudaMemcpyHostToDevice);
//计算
unsigned int grid_rows = (m + BLOCK_SIZE - 1) / BLOCK_SIZE;
unsigned int grid_cols = (k + BLOCK_SIZE - 1) / BLOCK_SIZE;
dim3 dimGrid(grid_cols, grid_rows);
dim3 dimBlock(BLOCK_SIZE, BLOCK_SIZE);
gpu_matrix_mult<<<dimGrid,dimBlock>>>(d_a,d_b,d_c,m,n,k);
cudaMemcpy(h_c,d_c,sizeof(int)*m*k,cudaMemcpyDeviceToHost);
// free memory
cudaFree(d_a);
cudaFree(d_b);
cudaFree(d_c);
cudaFreeHost(h_a);
cudaFreeHost(h_b);
cudaFreeHost(h_c);
cudaFreeHost(h_cc);
return 0;
}
利用nvprof查看
size=100
size=1000
总结
- 两个指针是指针的指针,cuda内存用指针指向数据,如要修改指针故要用两个星号。
- gpu和cpu异步运行,需要同步函数
- blockid要根据grid的分块判断,上图就把grid分为一个块,所以整个idx=(0,0),Dim=(12,6),threadidx
=(6,3),所以整个结果为(0,0)*(12,6)+(6,3)=(6,3)
- 通过工具分析,其实增加size=1000后计算所用时间变化并不大,说明了GPU并行计算的高校性。
- GPU大部分时间花在计算上,Host到Device的时间相比device到Host花销更大。