NVIDIA GPU와 CUDA 완전 해부: 왜 GPU가 AI를 지배하는가

CPU 아키텍처 (레이턴시 최적화):
┌────────────────────────────────────────────┐
│  코어 0  │  코어 1  │  코어 2  │  코어 3  │
│ (강력)   │ (강력)   │ (강력)   │ (강력)   │
│ OOO실행  │ OOO실행  │ OOO실행  │ OOO실행  │
├──────────┴──────────┴──────────┴───────────┤
│           L3 캐시 (32MB+)                  │
├────────────────────────────────────────────┤
│           메인 메모리 (DDR5, ~50 GB/s)     │
└────────────────────────────────────────────┘

GPU 아키텍처 (처리량 최적화) — H100 기준:
┌──┬──┬──┬──┬──┬──┬──┬──┬──┬──┬──┐
│SM│SM│SM│SM│SM│SM│SM│SM│SM│SM│...│  132개 SM
│  │  │  │  │  │  │  │  │  │  │   │  각 SM = 128 CUDA 코어
└──┴──┴──┴──┴──┴──┴──┴──┴──┴──┴───┘
       L2 캐시 (50MB)
       HBM3 (80GB, 3.35 TB/s)

GPU Device (H100 SXM5)
│
├── SM 0 (Streaming Multiprocessor)
│   ├── 4개 Warp Scheduler
│   ├── 4개 Dispatch Unit
│   ├── 128개 CUDA Core (FP32)
│   ├── 64개 FP64 Core
│   ├── 4개 Tensor Core (4세대, FP8/FP16/BF16/TF32 지원)
│   ├── 1개 특수함수 유닛 (SFU: sin, cos, sqrt 등)
│   ├── LD/ST 유닛 (Load/Store)
│   ├── L1 캐시 + 공유 메모리 = 228KB (비율 설정 가능)
│   └── 레지스터 파일 (65536개 x 32비트 레지스터)
│
├── SM 1 ...
├── SM 2 ...
│   ...
├── SM 131 ...
│
├── L2 캐시 (50MB)
└── HBM3 메모리 (80GB, 3.35 TB/s 대역폭)

Thread Block (예: 256개 스레드)
├── Warp 0: Thread 0~31   → 동일 명령어 실행 (동기화됨)
├── Warp 1: Thread 32~63  → 동일 명령어 실행 (동기화됨)
├── Warp 2: Thread 64~95  → 동일 명령어 실행 (동기화됨)
├── ...
└── Warp 7: Thread 224~255

// 이 코드는 Warp Divergence를 유발한다
__global__ void divergent_kernel(float* data, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx % 2 == 0) {
        // Warp 내 짝수 스레드만 실행
        data[idx] = data[idx] * 2.0f;  // Step 1
    } else {
        // Warp 내 홀수 스레드만 실행
        data[idx] = data[idx] + 1.0f;  // Step 2
    }
    // 결과: Step 1 실행 시 홀수 스레드는 idle
    //       Step 2 실행 시 짝수 스레드는 idle
    // 실제 처리량: 이론값의 50%
}

// Divergence 없음: 모든 스레드가 같은 경로
__global__ void coherent_kernel(float* data, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    // blockIdx.x 전체가 짝수이거나 홀수이면 divergence 없음
    if (blockIdx.x % 2 == 0) {
        data[idx] = data[idx] * 2.0f;
    } else {
        data[idx] = data[idx] + 1.0f;
    }
}

메모리 계층 (빠름 → 느림):

┌─────────────────────────────────────────────────────┐
│  레지스터 (32-bit, 65536개/SM)                      │
│  레이턴시: ~1 사이클    대역폭: 극대 (로컬 접근)     │
│  범위: 단일 스레드만 접근 가능                       │
├─────────────────────────────────────────────────────┤
│  L1 캐시 / 공유 메모리 (228KB/SM on H100)           │
│  레이턴시: ~32 사이클   대역폭: ~19 TB/s             │
│  범위: Thread Block 내 모든 스레드                   │
├─────────────────────────────────────────────────────┤
│  L2 캐시 (50MB on H100)                             │
│  레이턴시: ~200 사이클  대역폭: ~12 TB/s             │
│  범위: 모든 SM 공유                                  │
├─────────────────────────────────────────────────────┤
│  글로벌 메모리 HBM3 (80GB on H100)                  │
│  레이턴시: ~700 사이클  대역폭: 3.35 TB/s            │
│  범위: 모든 스레드 접근 가능                         │
└─────────────────────────────────────────────────────┘

// Naive 행렬 곱셈 커널 — 느린 이유를 이해하기 위한 버전
__global__ void matmul_naive(float* A, float* B, float* C, int N) {
    // 각 스레드가 출력 행렬의 한 원소를 담당
    int row = blockIdx.y * blockDim.y + threadIdx.y;
    int col = blockIdx.x * blockDim.x + threadIdx.x;

    if (row < N && col < N) {
        float sum = 0.0f;
        for (int k = 0; k < N; k++) {
            // 매 반복마다 글로벌 메모리 접근 2번 → 레이턴시 ~700 사이클
            sum += A[row * N + k] * B[k * N + col];
        }
        C[row * N + col] = sum;
    }
}

// 호출 예시
dim3 blockDim(16, 16);  // 256 스레드/블록
dim3 gridDim((N + 15) / 16, (N + 15) / 16);
matmul_naive<<<gridDim, blockDim>>>(A, B, C, N);

#define TILE_SIZE 16

// 타일링된 행렬 곱셈 — 공유 메모리 활용
__global__ void matmul_tiled(float* A, float* B, float* C, int N) {
    // 공유 메모리에 타일 선언 (L1 수준 레이턴시)
    __shared__ float tile_A[TILE_SIZE][TILE_SIZE];
    __shared__ float tile_B[TILE_SIZE][TILE_SIZE];

    int tx = threadIdx.x, ty = threadIdx.y;
    int row = blockIdx.y * TILE_SIZE + ty;
    int col = blockIdx.x * TILE_SIZE + tx;
    float sum = 0.0f;

    // K 차원을 타일 단위로 순회
    for (int t = 0; t < N / TILE_SIZE; t++) {
        // 단계 1: 각 스레드가 타일 하나씩 공유 메모리에 로드
        tile_A[ty][tx] = A[row * N + t * TILE_SIZE + tx];
        tile_B[ty][tx] = B[(t * TILE_SIZE + ty) * N + col];

        // 단계 2: 모든 스레드가 로드 완료할 때까지 동기화
        __syncthreads();

        // 단계 3: 공유 메모리에서 타일 내 내적 계산
        // 이 루프는 전부 공유 메모리(~32 사이클) 접근
        for (int k = 0; k < TILE_SIZE; k++) {
            sum += tile_A[ty][k] * tile_B[k][tx];
        }

        // 단계 4: 다음 타일 로드 전에 동기화
        __syncthreads();
    }

    C[row * N + col] = sum;
}

타일링 원리:

A 행렬:         B 행렬:
┌──┬──┬──┬──┐  ┌──┬──┬──┬──┐
│T0│T1│T2│T3│  │T0│T1│T2│T3│
├──┼──┼──┼──┤  ├──┼──┼──┼──┤
│  │  │  │  │  │  │  │  │  │
└──┴──┴──┴──┘  └──┴──┴──┴──┘

각 Thread Block은 C의 16×16 타일 하나를 담당.
A의 해당 행 타일들과 B의 해당 열 타일들을 순서대로
공유 메모리에 올려서 계산 → 글로벌 메모리 접근 최소화.

// WMMA API를 이용한 Tensor Core 직접 호출
#include <mma.h>
using namespace nvcuda::wmma;

__global__ void tensor_core_matmul(half* a_ptr, half* b_ptr,
                                    float* c_ptr, int M, int N, int K) {
    // 각 Warp가 16x16 출력 타일 하나를 담당
    // fragment = Tensor Core가 다루는 행렬 조각
    fragment<matrix_a, 16, 16, 16, half, row_major> a_frag;
    fragment<matrix_b, 16, 16, 16, half, col_major> b_frag;
    fragment<accumulator, 16, 16, 16, float> c_frag;

    // 누산기 초기화
    fill_fragment(c_frag, 0.0f);

    // K 차원 순회하며 타일 로드 후 곱셈
    for (int k = 0; k < K; k += 16) {
        // 글로벌 메모리에서 fragment 로드
        load_matrix_sync(a_frag, a_ptr + /* offset */, K);
        load_matrix_sync(b_frag, b_ptr + /* offset */, N);

        // 단 하나의 명령어로 16x16x16 행렬 곱 누산
        // 내부에서 Tensor Core 하드웨어 직접 호출
        mma_sync(c_frag, a_frag, b_frag, c_frag);
    }

    // 결과를 글로벌 메모리에 저장
    store_matrix_sync(c_ptr + /* offset */, c_frag, N, mem_row_major);
}

CUDA 실행 계층:

Grid (전체 커널 실행)
├── Block (0,0): 256 스레드
│   ├── Warp 0: Thread 0~31 ── SM에서 동시 실행
│   ├── Warp 1: Thread 32~63
│   └── ...
├── Block (1,0): 256 스레드
├── Block (2,0): 256 스레드
└── ...

각 Block은 하나의 SM에 할당됨.
SM은 동시에 여러 Block을 처리할 수 있음 (occupancy).

import torch

A = torch.randn(4096, 4096, device='cuda', dtype=torch.float16)
B = torch.randn(4096, 4096, device='cuda', dtype=torch.float16)

# 이 한 줄 내부에서 일어나는 일:
# 1. cuBLAS 또는 CUTLASS 커널 선택
# 2. Grid: (256, 256), Block: (16, 16) 등으로 구성
# 3. Tensor Core 활용한 tiled GEMM 실행
# 4. 총 ~4096개 Thread Block, 각 16*16=256 스레드
result = torch.matmul(A, B)  # → 내부적으로 cublasSgemmEx 호출

SM 리소스 한계 (H100 기준):
- 최대 동시 Thread Block: 32개
- 최대 동시 스레드: 2048개
- 레지스터: 65536개 (스레드 수 × 스레드당 레지스터로 나눔)
- 공유 메모리: 228KB (Block들이 나눠 씀)

Occupancy = 실제 활성 Warp 수 / 최대 가능 Warp 수
목표: 50% 이상 (메모리 레이턴시 숨기기 위해)

import torch
import time

# 1. dtype 선택: FP16 또는 BF16 사용 (Tensor Core 활성화)
A = torch.randn(8192, 8192, device='cuda', dtype=torch.bfloat16)
B = torch.randn(8192, 8192, device='cuda', dtype=torch.bfloat16)

# 2. cuBLAS를 위한 워밍업 (첫 실행은 커널 컴파일/선택에 시간 소요)
for _ in range(5):
    _ = torch.matmul(A, B)
torch.cuda.synchronize()

# 3. 측정
start = time.perf_counter()
for _ in range(100):
    C = torch.matmul(A, B)
torch.cuda.synchronize()
end = time.perf_counter()

# 4. TFLOPS 계산
# 8192x8192 행렬 곱: 2 * N^3 = 2 * 8192^3 = 1.1e12 FLOP
flops = 2 * 8192**3 * 100
elapsed = end - start
tflops = flops / elapsed / 1e12
print(f"달성 성능: {tflops:.1f} TFLOPS")
# 목표: 700+ TFLOPS (이론값 989의 70%+ 달성 가능)