CUDA 二维卷积
由于在二维卷积中卷积核多为横列数为奇数的矩阵，例如：3X3，5X5，本次代码演示只适用横列数为奇数的卷积核。
1、扩边和翻转
【资料图】
在进行二维卷积之前，我们要对源数据进行扩边，即在源数据外围添加 "0"，保证计算结果与原始数据大小一致。具体添加 "0" 的数目与卷积核大小有关，在 host 通过循环实现，代码如下：(bdy_x 和 bdy_y 代表所扩边的大小，例如3X3的卷积内核扩边 "1"，newm 和 newn 代表扩边后矩阵的大小，a是用来指向新数组的一维指针)
for (i = 0; i < newm; i++)
{
for (j = 0; j < newn; j++)
{
if (i < bdy_x || j < bdy_y || i >= m + bdy_x || j >= n + bdy_y)
{ a[i * newn + j] = 0; }
else {a[i * newn + j] = A[i - bdy_x][j - bdy_y];}
}
}
除此之外，还需要对内核进行180° 翻转，在 host 通过循环实现，代码如下：(将180° 翻转分解为：先左右对称，再上下对称，p 和 q 为卷积核大小，B1 和 B2 分别为翻转前、后的卷积核)
for (i = 0; i < p; i++)
{
for (j = 0; j < q; j++)
{ B2[i][j] = B1[p-i-1][q-j-1]; }
}
2、内核计算
首先将扩边后的矩阵根据内核的大小和移动切分成和卷积核大小相同计算单元，并保证要计算的点位于计算单元的中心位置，代码如下：( i 和 j 对应扩边前的下标 = 扩边后计算单元的中心位置，bdy_x 和 bdy_y 代表所扩边的大小，a 为扩边后的源数据，newm 和 newn 代表扩边后矩阵的大小，i 和 j 也是线程下标，共m*n 个线程)
for (i2 = 0; i2 < p; i2++)
{
for (j2 = 0; j2 < q; j2++)
{ newi2 = i2 - bdy_x; newj2 = j2 - bdy_y;
per[i2 * q + j2] = a[(i + bdy_x + newi2) * newn + (j + bdy_y + newj2)]; }
}
切分完成后，将每一个计算单元与卷积核点乘并通过上述循环累加。
mids += per[i2 * q + j2] * N[i2 * q + j2];
3、完整代码如下：
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
const int m = 32, n = 32, p = 3, q = 3;
__constant__ int N[p * q];
__global__ void conv2(int *s, int *a, int bdy_x, int bdy_y, int newn)
{
int i = blockIdx.x * blockDim.x + threadIdx.x;
int j = blockDim.y * blockIdx.y + threadIdx.y;
int i2, j2, newi2, newj2;
int mids = 0;
int per[p * q];
if (i < m && j < n)
{
for (i2 = 0; i2 < p; i2++)
{
for (j2 = 0; j2 < q; j2++)
{ newi2 = i2 - bdy_x; newj2 = j2 - bdy_y;
per[i2 * q + j2] = a[(i + bdy_x + newi2) * newn + (j + bdy_y + newj2)];
mids += per[i2 * q + j2] * N[i2 * q + j2]; }
}
s[i * n + j] = mids; }
}
int main(int argc, char* argv[])
{
int i, j; int * a, * s; int * d_a, * d_s;
int A[m][n], B1[p][q], B2[p][q];
int bdy_x = floor(p / 2), bdy_y = floor(q / 2);
int newm = m + 2 * bdy_x, newn = n + 2 * bdy_y;
//定义一个源数据和卷积核
for (i = 0; i < m; i++)
{
for (j = 0; j < n; j++)
{ A[i][j] = i + j; }
}
for (i = 0; i < p; i++)
{
for (j = 0; j < q; j++)
{ B1[i][j] = i + j; }
}
cudaMallocHost(&a, newm * newn * sizeof(int)); cudaMallocHost(&s, m * n * sizeof(int));
cudaMalloc(&d_a, newm * newn * sizeof(int)); cudaMalloc(&d_s, m * n * sizeof(int));
//扩边
for (i = 0; i < newm; i++)
{
for (j = 0; j < newn; j++)
{ if (i < bdy_x || j < bdy_y || i >= m + bdy_x || j >= n + bdy_y)
{ a[i * newn + j] = 0; }
else
{ a[i * newn + j] = A[i - bdy_x][j - bdy_y]; }
}
}
//内核翻转
for (i = 0; i < p; i++)
{
for (j = 0; j < q; j++)
{ B2[i][j] = B1[p-i-1][q-j-1]; }
}
cudaMemcpyToSymbol(N, B2, p * q * sizeof(int));
cudaMemcpy(d_a, a, newm * newn * sizeof(int), cudaMemcpyHostToDevice);
cudaMemcpy(d_s, s, m * n * sizeof(int), cudaMemcpyHostToDevice);
dim3 threads(8, 8);
dim3 blocks((m + threads.x - 1) / threads.x, (n + threads.y - 1) / threads.y);
conv2 << <blocks, threads >> > (d_s, d_a, bdy_x, bdy_y, newn);
cudaMemcpy(s, d_s, m * n * sizeof(int), cudaMemcpyDeviceToHost);
cudaDeviceSynchronize();
cudaFreeHost(a); cudaFreeHost(s);
cudaFree(d_a); cudaFreeHost(d_s);
}
