[Multi-Core & GPU Programming] HW4 - CUDA Conv2d

"CUDA Conv2d"

[Objective]

CUDA를 활용한 직접적인 2D Convolution 구현과 이미지를 행렬로 변환하는 im2col 기법 및 행렬 곱셈(GEMM)을 결합한 방식을 각각 구현하여, GPU 병렬 처리 환경에서 효율적인 메모리 접근 패턴과 고성능 연산 최적화 이해하기

[Theory]

▣ Convolution
· Input 데이터의 특정 영역과 Filter(Kernel)간의 요소별 곱셈 후 합계를 구하는 연산

▣ Parallel Conv2d
· GPU에서 각 출력 Pixel 값을 독립적으로 계산 가능 → 병렬성 극대화
· 효율적인 성능을 내기 위해서는 Shared Memory 활용 및 Coalesced Memory Access 활용 필수

▣ Parallel im2col combined with Matmul
· im2col(Image-To-Column) : Input 이미지에서 Convolution 연산이 수행될 모든 Patch를 추출하여 2D 행렬의 열로 펼치는 과정
    → 다차원 데이터인 이미지가 2차원 행렬 구조로 재구성

· GEMM : 변환된 이미지 행렬과 Filter를 펼친 행렬을 행렬 곱셈 연산으로 계산
    → GPU는 행렬 곱셈에 최적화

[Implementation]

▣ im2col
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
#include "cuda_runtime.h"
#define TILE_WIDTH 8
 
/**
 * @brief Launches the im2col operation, which rearranges image blocks into columns.
 *
 * This function is typically used in convolutional neural networks (CNNs) to
 * transform input image data into a format suitable for matrix multiplication.
 *
 * @param x        Pointer to the input tensor of shape (N, C, H, W), where:
 *                 - N: Batch size
 *                 - C: Number of channels
 *                 - H: Height of the input
 *                 - W: Width of the input
 * @param N        Number of images in the batch.
 * @param C        Number of channels in the input tensor.
 * @param H        Height of the input tensor.
 * @param W        Width of the input tensor.
 * @param kH       Height of the convolution kernel (filter).
 * @param kW       Width of the convolution kernel (filter).
 * @param padH     Padding applied to the height dimension.
 * @param padW     Padding applied to the width dimension.
 * @param strideH  Stride along the height dimension.
 * @param strideW  Stride along the width dimension.
 * @param dilH     Dilation factor for the height dimension.
 * @param dilW     Dilation factor for the width dimension.
 * @param out      Pointer to the output tensor, which stores the rearranged
 *                 image blocks in column format.
*/
__global__ void im2col_kernel(const float* x, int N, int C, int H, int W,
                              int kH, int kW, int padH, int padW, int strideH, int strideW,
                              int dilH, int dilW, int outH, int outW, float* out) {
        int idx = blockIdx.x * blockDim.x + threadIdx.x;
 
        int rows = C * kH * kW;
        int cols = outH * outW;
        int total = rows * cols;
 
        if (idx < total) {
                int w_out = idx % outW;
                int h_out = (idx / outW) % outH;
                int q = (idx / (outH * outW)) % kW;
                int p = (idx / (outH * outW * kW)) % kH;
                int c = idx / (outH * outW * kH * kW);
 
                int h_in = h_out * strideH - padH + p * dilH;
                int w_in = w_out * strideW - padW + q * dilW;
 
                if (h_in >= 0 && h_in < H && w_in >= 0 && w_in < W) {
                        int in_idx = (c * H + h_in) * W + w_in;
                        out[idx] = x[in_idx];
                } else {
                        out[idx] = 0;
                }
        }
}
 
void launch_im2col(const float* x, int N, int C, int H, int W,
                   int kH, int kW, int padH, int padW, int strideH, int strideW,
                   int dilH, int dilW, float* out)
{
        int outH = (H + 2 * padH - dilH * (kH - 1) - 1) / strideH + 1;
        int outW = (W + 2 * padW - dilW * (kW - 1) - 1) / strideW + 1;
 
        int rows = C * kH * kW;
        int cols = outH * outW;
        int total = rows * cols;
        int block = 256;
        int grid = (total + block - 1) / block;
 
        im2col_kernel<<<grid, block>>>(x, N, C, H, W, kH, kW, padH, padW, strideH, strideW, dilH, dilW, outH, outW, out);
}
 
cs

> im2col_kernel
① 각 Thread가 im2col 행렬(Unrolled Input Features)의 하나의 원소(c, p, q, h_out, w_out)를 담당
    · idx : 1차원 Block과 Thread Index를 사용하여 현재 Thread의 고유 ID를 계산

② im2col 행렬(Unrolled Input Features) 행렬 크기 계산
    · rows : im2col 결과 행렬(Unrolled Input Features)의 행 개수 → 12(= 3 x 2 x 2)
    · cols : im2col 결과 행렬(Unrolled Input Features)의 열 개수 → 4(= 2 x 2) (Batch Size = 1 → N = 1)
    · total : im2col 결과 행렬(Unrolled Input Features)의 전체 원소 개수 → 48(= 12 x 4)

③ 범위 밖 Thread를 제외한 후, idx를 이용해 c, p, q, h_out, w_out을 계산
    · w_out : Output Feature Map에서의 열 Index
    · h_out : Output Feature Map에서의 행 Index
    · q : Filter 내부의 열 Index
    · p : Filter 내부의 행 Index
    · c : Input Channel 번호

④ 실제 Input Image에서 읽을 위치 계산
    · h_in : Input Image에서 실제로 읽을 행 Index
    · w_in : Input Image에서 실제로 읽을 열 Index

⑤ Input Image 안쪽이면 해당 Index의 값을 복사하고, Input Image 바깥쪽이면 Padding 영역이므로 0을 할당

> launch_im2col
① Output Feature Map 크기 계산
    · outH : Output Feature Map의 Height → 2
    · outW : Output Feature Map의 Width → 2

② im2col 결과 행렬(Unrolled Input Features) 크기 계산
    · rows : im2col 결과 행렬(Unrolled Input Features)의 행 개수 → 12(= 3 x 2 x 2)
    · cols : im2col 결과 행렬(Unrolled Input Features)의 열 개수 → 4(= 2 x 2)
    · total : im2col 결과 행렬(Unrolled Input Features)의 전체 원소 개수 → 48(= 12 x 4)

③ CUDA Block 크기 설정 → 한 Block 안에서 실행할 Thread 개수 설정
④ CUDA Grid 크기 계산 → im2col 결과 행렬의 전체 원소 개수만큼 Thread 필요
⑤ GPU에서 im2col_kernel 실행

▣ matmul
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
/**
 * @brief Launches a matrix multiplication operation on the provided matrices.
 *
 * This function performs the matrix multiplication operation C = A * B, where:
 * - A is an MxK matrix.
 * - B is a KxN matrix.
 * - C is the resulting MxN matrix.
 *
 * @param A Pointer to the first input matrix (MxK).
 * @param B Pointer to the second input matrix (KxN).
 * @param C Pointer to the output matrix (MxN) where the result will be stored.
 * @param M Number of rows in matrix A and matrix C.
 * @param N Number of columns in matrix B and matrix C.
 * @param K Number of columns in matrix A and rows in matrix B.
 */
__global__ void matmul_kernel(const float* A, const float* B, float* C, int M, int N, int K) {
        __shared__ float tileA[TILE_WIDTH][TILE_WIDTH];
        __shared__ float tileB[TILE_WIDTH][TILE_WIDTH];
 
        int col = blockIdx.x * blockDim.x + threadIdx.x;
        int row = blockIdx.y * blockDim.y + threadIdx.y;
        float sum = 0;
 
        for (int t = 0; t < (K + TILE_WIDTH - 1) / TILE_WIDTH; t++) {
                if (row < M && t * TILE_WIDTH + threadIdx.x < K) {
                        tileA[threadIdx.y][threadIdx.x] = A[row * K + t * TILE_WIDTH + threadIdx.x];
                } else {
                        tileA[threadIdx.y][threadIdx.x] = 0;
                }
 
                if (t * TILE_WIDTH + threadIdx.y < K && col < N) {
                        tileB[threadIdx.y][threadIdx.x] = B[(t * TILE_WIDTH + threadIdx.y) * N + col];
                } else {
                        tileB[threadIdx.y][threadIdx.x] = 0;
                }
                __syncthreads();
 
                for (int k = 0; k < TILE_WIDTH; k++) {
                        sum += tileA[threadIdx.y][k] * tileB[k][threadIdx.x];
                }
                __syncthreads();
        }
        if (row < M && col < N) {
                C[row * N + col] = sum;
        }
}
// Tiling X Version
/*__global__ void matmul_kernel(const float* A, const float* B, float*C, int M, int N, int K) {
        int col = blockIdx.x * blockDim.x + threadIdx.x;
        int row = blockIdx.y * blockDim.y + threadIdx.y;
        float sum = 0;
        if (row < M && col < N) {
                for (int k = 0; k < K; k++) {
                        sum += A[row * K + k] * B[k * N + col];
                }
                C[row * N + col] = sum;
        }
}*/
 
 
void launch_matmul(const float* A, const float* B, float* C, int M, int N, int K)
{
        dim3 block(TILE_WIDTH, TILE_WIDTH);
        dim3 grid((N + TILE_WIDTH - 1) / TILE_WIDTH, (M + TILE_WIDTH - 1) / TILE_WIDTH);
 
        matmul_kernel<<<grid, block>>>(A, B, C, M, N, K);
}
 
cs

> matmul_kernel
① Shared Memory 선언
② 각 Thread가 담당할 결과 행렬 C의 Column과 Row Index를 계산하고, 곱셈 결과를 누적할 변수 선언
③ A, B 행렬의 일부 Tile을 Shared Memory에 Load한 후, 동기화 수행
④ Shared Memory에 저장된 Tile 값을 이용해 부분 내적을 수행하고 sum에 누적한 후, 다음 Tile 계산 전에 동기화 수행
⑤ 누적된 sum을 마지막에 결과 행렬 C에 저장

> launch_matmul
① Block 하나에 Thread를 (TILE_WIDTH x TILE_WIDTH)개 배치 → Block 하나가 결과 행렬 C의 Tile 영역을 담당
② 결과 행렬 C 전체를 처리할 수 있는 Block 개수를 계산
③ GPU에서 matmul_kernel 실행

▣ direct_conv2d
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
/**
 * @brief Launches a 2D convolution operation using the direct convolution method.
 *
 * @param x Pointer to the input tensor of shape (N, C, H, W), where:
 *          - N: Batch size
 *          - C: Number of input channels
 *          - H: Height of the input
 *          - W: Width of the input
 * @param w Pointer to the weight tensor of shape (K, C, kH, kW), where:
 *          - K: Number of output channels
 *          - C: Number of input channels
 *          - kH: Height of the kernel
 *          - kW: Width of the kernel
 * @param y Pointer to the output tensor of shape (N, K, outH, outW), where:
 *          - outH: Computed output height
 *          - outW: Computed output width
 * @param N Number of input batches.
 * @param C Number of input channels.
 * @param H Height of the input tensor.
 * @param W Width of the input tensor.
 * @param K Number of output channels.
 * @param kH Height of the convolution kernel.
 * @param kW Width of the convolution kernel.
 * @param padH Padding applied to the height dimension.
 * @param padW Padding applied to the width dimension.
 * @param strideH Stride applied to the height dimension.
 * @param strideW Stride applied to the width dimension.
 * @param dilH Dilation applied to the height dimension of the kernel.
 * @param dilW Dilation applied to the width dimension of the kernel.
 */
__global__ void conv2d_direct_kernel(const float* x, const float* w, float* y,
                                     int N, int C, int H, int W,
                                     int K, int kH, int kW,
                                     int padH, int padW, int strideH, int strideW,
                                     int dilH, int dilW, int outH, int outW) {
        int idx = blockIdx.x * blockDim.x + threadIdx.x;
        int total = K * outH * outW;
 
        if (idx < total) {
                int w_out = idx % outW;
                int h_out = (idx / outW) % outH;
                int k = idx / (outH * outW);
 
                float sum = 0;
 
                for (int c = 0; c < C; c++) {
                        for (int p = 0; p < kH; p++) {
                                for (int q = 0; q < kW; q++) {
                                        int h_in = h_out * strideH - padH + p * dilH;
                                        int w_in = w_out * strideW - padW + q * dilW;
 
                                        if (h_in >= 0 && h_in < H && w_in >= 0 && w_in < W) {
                                                int x_idx = (c * H + h_in) * W + w_in;
                                                int weight_idx = ((k * C + c) * kH + p) * kW + q;
 
                                                sum += x[x_idx] * w[weight_idx];
                                        }
                                }
                        }
                }
                y[idx] = sum;
        }
}
 
void launch_conv2d_direct(const float* x, const float* w, float* y,
                          int N, int C, int H, int W,
                          int K, int kH, int kW,
                          int padH, int padW, int strideH, int strideW,
                          int dilH, int dilW)
{
        int outH = (H + 2 * padH - dilH * (kH - 1) - 1) / strideH + 1;
        int outW = (W + 2 * padW - dilW * (kW - 1) - 1) / strideW + 1;
 
        int block = 256;
        int grid = (K * outH * outW + block - 1) / block;
 
        conv2d_direct_kernel<<<grid, block>>>(x, w, y, N, C, H, W, K, kH, kW, padH, padW, strideH, strideW, dilH, dilW, outH, outW);
}
 
cs

> conv2d_direct_kernel
① 각 Thread가 담당할 1D Index 계산
② 전체 출력 원소 개수 계산
③ 현재 Thread가 계산할 출력의 Column과 Row, Output Channel Index를 계산하고, 곱셈 결과를 누적할 변수 선언
④ 모든 Input Channel을 반복하면서 Input Image의 Column과 Row Index를 계산
⑤ Padding 때문에 입력 범위를 벗어난 경우를 제외하고, Input과 Weight의 1차원 Index 계산 후 서로 곱하여 누적
⑥ 최종 결과를 Output에 저장

> launch_conv2d_direct
① 출력의 Row와 Column 개수를 계산
② Block 하나에 Thread 256개를 배치
③ 전체 결과 원소를 처리할 수 있는 Block 개수 계산
④ GPU에서 conv2d_direct_kernel 실행

[Result]

· Direct Kernel은 Shared Memory Reuse가 없어 Global Memory에 반복적으로 접근하지만, im2col 행렬 생성 비용 X
    → im2col + matmul 방식보다 빠름
· Filter Size가 7x7일때도 Direct Conv2d보다 im2col + matmul의 성능이 더 낮음
    ∵ Batch Size가 1로 고정되어 있기 때문에 im2col + matmul의 효율↓

∴ Batch Size가 증가
    → im2col Matrix Column↑
    → GEMM 크기↑
    → GPU 병렬성↑
    → im2col + matmul의 성능↑

[Discussion]

· HW3과 마찬가지로 실행중인 Process를 Kill하여 여유 Process를 보유하고 있어야 함

<im2col_kernel에서 for문 사용>
· 각 Thread가 Output Feature의 한 원소를 담당 = im2col의 한 Column 담당
    → 각 Thread가 Kernel 전체를 다 돌며, 각 Thread당 (kH x kW)번의 연산 반복

<matmul_kernel에서 Tiling X>
· Naïve 방식이 Tiling 방식보다 빠름
    → 연산량에 비해 Shared Memory로 복사하여 Load하고, 모든 Thread가 대기하도록 동기화하는 것에서 Overhead 발생
    → Tiling에서 for 구문 안의 if-else로 인해 Warp Divergence가 반복적으로 발생하여 연산 Unit이 낭비되어 성능 저하


<TILE_WIDTH = 16>                    <TILE_WIDTH = 32>
· im2col은 K차원이 매우 작아서 TILE_WIDTH를 늘리면 Thread 낭비 발생

[Reference]

· HW4 - CUDA Conv2d (MGP) - Yongjun Park
· 7_cuda_dnn_mgp_2026 (MGP) - Yongjun Park
· [Multi-Core & GPU Programming] HW3 - CUDA LoRA

댓글

댓글 쓰기

이 블로그의 인기 게시물

[Mini-NPU RTL] NN Reference Model

이미지
"NN Reference Model" 1. Neuron 단 1개의 Neuron만 존재하며 해당 Neuron이 가중합을 잘 구하고 활성화 함수를 제대로 통과시키는지 확인 → 28x28=784개의 입력 데이터와 784개의 가중치를 각각 곱한 후 ∑, 이후 활성화 함수를 통과시킨 뒤, 1개의 값 출력 → "test_data.txt"에서 데이터를 읽어 Neuron에 직접 넣음 → 가중치와 편향이 "w_1_0.mif", "b_1_0.mif"로부터 정적으로 고정되어 Load neuron_tb.sv) include.v) 공통 상수를 미리 정해두는 파일로, 숫자를 고치면 Network 깊이, Neuron 개수 등을 조절 ReLU.v) 공통 활성화 함수인 ReLU를 하드웨어로 구현한 코드 Sig_ROM.v) 공통 Sigmoid 함수 연산을 실시간으로 계산하지 않고  미리 계산된 결과값 을 메모리에서 꺼내 쓰는  LookUp Table  방식의 구현 Weight_Memory.v) 공통 가중치를 저장하는 메모리 장치로 가중치를 읽고 쓴다. Neuron.v) 공통 입력 데이터를 받아 가중치와 곱하고, 누적 한 후, 편향을 더하고, 활성화 함수를 거쳐 출력하는 연산 Unit 2. Network 수많은 Neuron이 합쳐진 전체 시스템의 데이터 흐름과 제어 Logic이 제대로 작동하는지 확인 → 28x28=784개의 입력 데이터가 30개의 Neuron에게 동시에 전달되면 784x30개의 가중치를 곱해 30개의 결과 출력     (첫 번째 입력이 1 Cycle에 서로 다른 30개의 가중치와 병렬 연산으로 곱해진다. ∴ 784 Cycle) → 30개의 입력이 다음 Layer의 20개 Neuron에 전달되어 30x20개의 가중치를 곱해 20개의 결과 출력 → 20개의 입력이 다음 Layer의 10개 Neuron에 전달되어 20x10개의 가중치를 곱해 10개의 결과 출력 → 10개의 출력은 입력 이미지 0~9의 점수를 나...
자세한 내용 보기

[Mini-NPU RTL] TPU (Study Paper)

이미지
"In-Datacenter Performance Analysis of a Tensor Processing Unit" 1. Introduction to Neural Networks   NN은 입력 값들의 가중 합에 비선형 함수(ex. max(0, value))를 적용하는 단순한 인공 뉴런에 기초한다.   ※ 가중합 : 입력 데이터에 중요도인 가중치 를 곱해서 더한다.   ※ 비선형 함수 : 계산된 합을 그대로 쓰지 않고 변형하는 함수(활성화 함수)   ※ max(0, value) : ReLU(Rectified Linear Unit) 함수이며, 음수면 0, 양수면 해당 값을 그대로 내보낸다.   NN은 계층(Layer) 으로 모이며, 한 계층의 출력은 다음 계층의 입력이 된다.   DNN(Deep Neural Networks) 은 NN에서 계층의 수가 많은 것 을 의미한다.   NN은 개발 단계인 Training(학습) 단계와 실전 단계인 Inference(추론) 단계가 있다.   개발자는 계층의 수 와 NN의 유형 을 선택하고, 학습 과정을 통해 가중치 를 결정한다.   학습이 끝난 후, 추론 단계에서는 소수점 아래 긴 자리까지 필요 없는 경우가 많다.   ∴ Quantization(양자화) 단계에서 단순한 정수( 8-bit )로 바꿔 용량↓ & 속도↑         → 8-bit 정수 곱셈 : 16-bit 부동소수점 곱셈보다 에너지 소모 6배↓ & 칩 면적 6배↓         → 8-bit 정수 덧셈 : 16-bit 부동소수점 덧셈보다 에너지 소모 13배↓ & 칩 면적 38배↓   NN의 유형은 3가지로 볼 수 있다.   ·  MLP(Multi-Layer Perceptrons) : 각 계층은 이전 계층의 모든 출력 값(Fully Conn...
자세한 내용 보기