AXPY와 FMA 연산 | Notes

왜 이걸 배우는가?

ChatGPT가 답변을 생성할 때 무슨 일이 벌어질까? GPU 안에서 수십억 번의 곱셈과 덧셈이 동시에 일어난다. 게임 캐릭터가 부드럽게 움직일 때도, 일기예보를 계산할 때도 마찬가지다.

이번에 알아볼 AXPY와 FMA가 바로 그 핵심이다. PyTorch로 모델을 학습시킬 때 Adam 옵티마이저가 파라미터를 업데이트하는 과정, 게임에서 반투명 효과를 만들 때 색상을 섞는 계산, 기상청이 내일 날씨를 예측할 때 대기 흐름을 시뮬레이션하는 모든 과정이 이 연산들의 반복이다.

벡터 덧셈 장에서 실망스러운 결과를 봤다. GPU 커널 자체는 56배 빠른데 메모리 전송 때문에 전체로는 오히려 느렸다. AXPY는 이 문제를 어떻게 극복하는지 보여준다. 같은 메모리로 계산을 더 많이 하면 GPU가 유리해진다. 이것이 GPU 최적화의 핵심 원리다.

FMA는 한 단계 더 들어간다. 하드웨어가 어떻게 "곱하고 더하기"를 한 명령어로 처리하는지 이해하면, 왜 최신 GPU가 딥러닝에 필수인지, 왜 NVIDIA가 Tensor Core를 만들었는지 보인다. 다음 장부터 배울 행렬 곱셈과 컨볼루션도 결국 이 원리의 확장이다.

AXPY: 메모리 병목 극복하기

벡터 덧셈의 실망스러운 결과

이전 장에서 벡터 덧셈을 GPU로 병렬화했다. GPU 커널 단독으로는 CPU보다 56배 빠르다. 하지만 메모리 전송을 포함하면 오히려 느려진다. 전체 시간의 99%가 메모리 전송에 소비된다.

계산 자체는 매우 빠르다. 하지만 데이터를 옮기는 데 시간이 너무 오래 걸린다. GPU의 수천 개 코어가 대부분 놀고 있다. 메모리를 기다리며 시간을 낭비하고 있는 것이다.

해결책: 계산 강도를 높이자

같은 데이터로 더 많은 계산을 하면 메모리 전송 비용을 상쇄할 수 있다. AXPY가 바로 그 예시다.

AXPY = A times X Plus Y = $a \times x + y$

벡터 덧셈은 덧셈만 하지만, AXPY는 곱셈과 덧셈을 모두 한다. 벡터 덧셈은 메모리 12 bytes를 읽고 덧셈 1번만 수행한다(1 FLOP). 반면 AXPY는 같은 12 bytes를 읽지만 곱셈과 덧셈을 모두 수행한다(2 FLOP).

메모리 전송은 같은데 연산량은 2배다. GPU의 계산 유닛을 더 많이 활용할 수 있다. 하지만 2배 정도로는 여전히 부족하다. 진짜 극적인 차이는 다음 장의 행렬 곱셈에서 나타난다.

계산 강도(Arithmetic Intensity)란?
메모리 접근 대비 연산량의 비율이다. 계산 강도 = (연산 횟수) / (메모리 접근 바이트). 높을수록 GPU에 유리하다. 메모리 대역폭 병목이 줄어들고, 계산 유닛 활용률이 올라간다.

벡터 덧셈의 계산 강도를 계산하면 다음과 같다. 메모리는 $A[i]$ , $B[i]$ 읽기 8 bytes, $C[i]$ 쓰기 4 bytes로 총 12 bytes다. 연산은 덧셈 1회로 1 FLOP다.

\text{Arithmetic Intensity} = \frac{1 \text{ FLOP}}{12 \text{ bytes}} = 0.083 \text{ FLOP/byte}

AXPY의 계산 강도는 다음과 같다. 메모리는 $X[i]$ , $Y[i]$ 읽기 8 bytes (스칼라 $a$ 는 캐시됨), $Z[i]$ 쓰기 4 bytes로 총 12 bytes다. 연산은 곱셈 1회, 덧셈 1회로 2 FLOP다.

\text{Arithmetic Intensity} = \frac{2 \text{ FLOP}}{12 \text{ bytes}} = 0.167 \text{ FLOP/byte}

AXPY가 벡터 덧셈보다 2배 높다. 메모리 전송 비용이 상대적으로 줄어든다.

AXPY 루틴이란

AXPY는 BLAS(Basic Linear Algebra Subprograms) 라이브러리의 기본 연산이다. 이름은 "A times X Plus Y"를 줄인 것이다.

벡터 형태로 표현하면 다음과 같다.

\begin{aligned} Z &= a \cdot X + Y \\ \\ z_1 &= a \cdot x_1 + y_1 \\ z_2 &= a \cdot x_2 + y_2 \\ z_3 &= a \cdot x_3 + y_3 \\ &\vdots \\ z_n &= a \cdot x_n + y_n \end{aligned}

$X$ , $Y$ , $Z$ 는 $n$ 차원 벡터이고 $a$ 는 스칼라 값이다. 모든 원소에 동일한 $a$ 를 곱한 후 $Y$ 의 해당 원소를 더한다.

CPU 버전 구현

AXPY의 정의와 역사를 확인했다. 실제 코드로 구현하여 성능을 측정한다.

순차 처리 코드

먼저 CPU 기준 구현을 작성한다. 비교 기준점을 만드는 것이다.

for (unsigned i = 0; i < vecSize; ++i) {
    vecZ[i] = saxpy_a * vecX[i] + vecY[i];
}

for (unsigned i = 0; i < vecSize; ++i) {
    vecZ[i] = saxpy_a * vecX[i] + vecY[i];
}

각 원소마다 곱셈 1회, 덧셈 1회를 수행한다.

성능 측정

Intel Core i5-3570에서 527,661 usec가 소요된다. 약 0.53초다.

벡터 덧셈은 459,271 usec였다. AXPY가 약 15% 느리다. 벡터 덧셈은 원소당 1 FLOP, AXPY는 2 FLOP로 연산량이 2배다.

\frac{527{,}661}{459{,}271} = 1.15

연산량이 2배인데 시간은 1.15배만 늘었다. 왜 그럴까?

CPU의 파이프라인과 슈퍼스칼라 아키텍처가 여러 명령어를 동시에 처리한다. 곱셈과 덧셈 유닛이 별도로 존재하므로, 병렬로 실행될 수 있다. 메모리 대역폭이 병목이므로 계산 추가가 전체 시간에 미치는 영향이 작다.

처리량 계산

연산량은 $268{,}435{,}456 \times 2 = 536{,}870{,}912$ FLOP이고, 시간은 $527{,}661$ usec = $0.5277$ sec다.

\text{Throughput} = \frac{536{,}870{,}912}{0.5277} = 1{,}017 \text{ MFLOPS} \approx 1 \text{ GFLOPS}

Intel Core i5-3570의 이론적 최대 성능은 약 100 GFLOPS다. 실제 1%만 활용하고 있다. 메모리 대역폭에 제한되기 때문이다.

CUDA 병렬 구현

CPU 기준 성능을 확인했다. 이제 GPU로 병렬화하여 가속한다.

커널 함수 설계

GPU 커널은 벡터 덧셈과 거의 동일하다. 덧셈 대신 AXPY 연산만 다르다.

__global__ void kernelSAXPY(float* z, const float a,
                             const float* x, const float* y, unsigned n) {
    unsigned i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n) {
        z[i] = a * x[i] + y[i];
    }
}

__global__ void kernelSAXPY(float* z, const float a,
                             const float* x, const float* y, unsigned n) {
    unsigned i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n) {
        z[i] = a * x[i] + y[i];
    }
}

각 스레드가 자신의 인덱스를 계산하고 해당 원소를 처리한다.

성능 측정 결과

GeForce RTX 2070에서 8,130 usec가 소요된다.

CPU 527,661 usec와 비교하면 다음과 같다.

\text{Speedup} = \frac{527{,}661}{8{,}130} = 64.9\times

벡터 덧셈은 56배 가속이었다. AXPY는 65배 가속이다. 계산 강도가 높아지면서 가속 비율도 증가했다.

처리량 분석

GPU 처리량을 계산하면 다음과 같다. 시간은 8,130 usec = 0.00813 sec다.

\text{Throughput} = \frac{536{,}870{,}912}{0.00813} = 66{,}019 \text{ MFLOPS} \approx 66 \text{ GFLOPS}

GeForce RTX 2070의 이론적 최대 성능은 약 7,500 GFLOPS다. 실제 0.9%만 활용하고 있다. 여전히 메모리 대역폭에 제한된다.

하지만 CPU 1 GFLOPS에 비하면 66배 빠르다. 계산 유닛 개수의 차이가 드러난다.

FMA: 하드웨어 가속의 비밀

AXPY를 GPU로 병렬화했다. 계산 강도를 높이는 방법을 배웠다. 이제 하드웨어가 어떻게 계산 자체를 가속하는지 알아본다.

왜 FMA가 필요한가?

AXPY의 핵심 연산은 $a \times x + y$ 다. 이 단순해 보이는 계산에도 최적화 여지가 있다.

전통적인 방식은 두 단계로 나뉜다.

곱셈: $a \times x$ 계산 → 결과를 레지스터에 저장
덧셈: 저장된 값 + $y$ 계산

문제는 중간 결과를 저장하고 다시 읽어야 한다는 것이다. 명령어 2개, 반올림 2번이 필요하다.

FMA(Fused Multiply-Add)는 이 두 단계를 하나의 명령어로 합친다. 중간 저장 없이 곱셈과 덧셈을 한 번에 처리한다.

왜 명령어가 2개에서 1개로 줄어드나?

전통적인 방식은 두 개의 독립적인 명령어로 처리된다. 아래 내용은 컴퓨터 구조 개념을 다루므로 필요한 용어는 blockquote로 설명했다.

명령어 1: 곱셈

plain

mul r1, a, x  // a × x 계산하고 r1에 저장

mul r1, a, x  // a × x 계산하고 r1에 저장

곱셈기가 계산한다. 결과를 FP32로 반올림해서 레지스터 파일에 쓴다.

레지스터 파일이란?
CPU/GPU 내부의 초고속 저장 공간이다. 연산에 필요한 값을 임시로 보관한다. 메모리보다 훨씬 빠르지만 용량이 매우 작다. 각 레지스터는 r0, r1 같은 이름을 가진다.

명령어 2: 덧셈

plain

add r2, r1, y  // r1 읽어서 y와 더함

add r2, r1, y  // r1 읽어서 y와 더함

레지스터 파일에서 r1을 읽는다. 덧셈기가 계산한다. 결과를 FP32로 반올림해서 레지스터 파일에 쓴다.

문제는 중간에 레지스터 파일 접근이 끼어든다는 것이다. 레지스터 읽기/쓰기에 시간이 걸리고, 파이프라인에 stall이 발생할 수 있다.

파이프라인 stall이란?
CPU/GPU는 여러 명령어를 동시에 처리한다(파이프라인). 한 명령어가 이전 명령어 결과를 기다리면 파이프라인이 멈춘다(stall). 성능이 저하된다.

FMA 하드웨어 회로

FMA는 전용 회로를 가진다.

plain

[곱셈기] → (내부 버퍼) → [덧셈기] → 반올림 → 레지스터
   ↓                        ↑
확장 정밀도 유지            y 입력

[곱셈기] → (내부 버퍼) → [덧셈기] → 반올림 → 레지스터
   ↓                        ↑
확장 정밀도 유지            y 입력

곱셈 결과를 레지스터 파일에 쓰지 않는다. 내부 버퍼에 확장 정밀도(FP64 이상)로 유지한다. 덧셈기가 이 확장 정밀도 결과에 y를 더한다. 최종 결과만 FP32로 반올림한다.

하나의 명령어로 처리되므로 명령어 fetch/decode가 1회만 필요하다.

명령어 Fetch/Decode란?
CPU/GPU가 명령어를 실행하는 과정이다. Fetch는 메모리에서 명령어를 가져오는 단계다. Decode는 명령어를 해석해서 어떤 연산인지 파악하는 단계다. 전통 방식은 곱셈 명령어와 덧셈 명령어를 각각 fetch/decode해야 한다. FMA는 1회만 필요하다.

레지스터 파일 접근도 감소한다(중간 쓰기/읽기 없음). 반올림도 1회만 수행하므로 정밀도가 향상된다.

Fused Multiply-Add의 장점

AXPY 연산 $a \times x + y$ 는 곱셈과 덧셈 두 단계로 이루어진다. 전통적인 방식은 곱셈 결과를 레지스터에 저장했다가 다시 읽어야 한다.

\text{temp} = a \times x, \quad \text{result} = \text{temp} + y

중간 결과 temp를 레지스터에 저장했다가 다시 읽어야 한다. 레지스터 압력이 증가하고, 파이프라인에 데이터 의존성이 생긴다.

장점이 세 가지 있다.

속도가 향상된다. 명령어가 2개에서 1개로 줄어들면서 명령어 디코딩 오버헤드가 감소한다. 파이프라인 효율도 증가한다.

정밀도가 향상된다. 반올림이 2번에서 1번으로 줄어들면서 중간 결과의 반올림 오차가 제거된다. 수치 안정성이 높아진다.

레지스터를 절약한다. 중간 결과를 저장할 필요가 없으므로 레지스터 압력이 감소한다. 더 많은 병렬 작업이 가능해진다.

반올림 오차 비교

전통적인 MAC/MAD 명령어는 두 번 반올림한다.

\text{temp} = \text{round}(a \times x), \quad \text{result} = \text{round}(\text{temp} + y)

FMA는 한 번만 반올림한다.

\text{result} = \text{round}(a \times x + y)

정밀도 향상 예시를 보자.

plain

a = 1.0
x = 1.0e-20
y = 1.0

MAC/MAD:
  temp = round(1.0 * 1.0e-20) = 1.0e-20
  result = round(1.0e-20 + 1.0) = 1.0  (작은 값이 사라짐)

FMA:
  result = round(1.0 * 1.0e-20 + 1.0) = 1.0 + 1.0e-20  (정확히 표현 가능하면 보존)

a = 1.0
x = 1.0e-20
y = 1.0

MAC/MAD:
  temp = round(1.0 * 1.0e-20) = 1.0e-20
  result = round(1.0e-20 + 1.0) = 1.0  (작은 값이 사라짐)

FMA:
  result = round(1.0 * 1.0e-20 + 1.0) = 1.0 + 1.0e-20  (정확히 표현 가능하면 보존)

큰 수와 작은 수를 더할 때 FMA가 더 정확하다.

CUDA FMA 함수

CUDA math 라이브러리는 FMA 함수를 제공한다.

float fmaf(float a, float x, float y);
double fma(double a, double x, double y);

float fmaf(float a, float x, float y);
double fma(double a, double x, double y);

반환값은 (a * x + y)를 FMA 명령어로 계산한 결과다. CPU와 GPU 모두 하드웨어 명령어로 구현되어 있다.

PTX(Parallel Thread Execution) 어셈블리에서는 fma.rn.f32 명령어로 나타난다. rn은 round to nearest를 의미한다.

FMA를 사용한 SAXPY

FMA 함수를 배웠다. AXPY 커널에 적용하여 성능을 비교한다.

커널 수정

커널 함수의 계산 부분을 FMA로 교체한다.

__global__ void kernelSAXPY(float* z, const float a,
                             const float* x, const float* y, unsigned n) {
    unsigned i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n) {
        z[i] = fmaf(a, x[i], y[i]);
    }
}

__global__ void kernelSAXPY(float* z, const float a,
                             const float* x, const float* y, unsigned n) {
    unsigned i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n) {
        z[i] = fmaf(a, x[i], y[i]);
    }
}

a * x[i] + y[i] 대신 fmaf(a, x[i], y[i])를 사용한다. 의미는 동일하지만 하드웨어 구현이 다르다.

성능 측정

256M 원소 처리 시 8,121 usec가 소요된다.

일반 연산 8,130 usec와 거의 차이가 없다. 왜 그럴까?

컴파일러가 이미 최적화를 수행했기 때문이다. NVCC는 -fmad=true 플래그가 기본값이다. 일반 곱셈-덧셈 표현식을 자동으로 FMA로 변환한다.

하지만 명시적으로 FMA를 사용하면 장점이 있다.

컴파일러 의존성을 제거한다. 최적화 설정과 무관하게 FMA 사용이 보장된다. -fmad=false로 컴파일해도 FMA가 적용된다.

정밀도를 보장한다. 정확히 한 번만 반올림함을 명시할 수 있다. 수치 알고리즘에서 중요한 부분이다.

코드 의도를 명확히 한다. FMA 사용 의도를 명시적으로 표현한다. 코드 리뷰 시 이해하기 쉬워진다.

FMA 활용 예제: 선형 보간

FMA의 기본 사용법을 배웠다. 실제 응용 사례로 선형 보간을 구현한다.

선형 보간이란

선형 보간(Linear Interpolation)은 두 값 사이의 중간값을 계산하는 연산이다. 애니메이션, 물리 시뮬레이션, 그래픽에서 필수적이다.

f(t) = (1 - t) \cdot v_0 + t \cdot v_1

$t=0$ 이면 $v_0$ 을, $t=1$ 이면 $v_1$ 을 반환한다. $t=0.5$ 면 중간값이 나온다.

예를 들어 v₀=10, v₁=20이면 다음과 같다.

plain

f(0.0) = 1.0 * 10 + 0.0 * 20 = 10
f(0.5) = 0.5 * 10 + 0.5 * 20 = 15
f(1.0) = 0.0 * 10 + 1.0 * 20 = 20

f(0.0) = 1.0 * 10 + 0.0 * 20 = 10
f(0.5) = 0.5 * 10 + 0.5 * 20 = 15
f(1.0) = 0.0 * 10 + 1.0 * 20 = 20

0과 1 사이를 선형으로 보간한다.

수식 변형

식을 전개하면 다음과 같다.

\begin{aligned} f(t) &= (1 - t) \cdot v_0 + t \cdot v_1 \\ &= v_0 - t \cdot v_0 + t \cdot v_1 \\ &= v_0 + t \cdot (v_1 - v_0) \end{aligned}

마지막 형태가 FMA에 적합하다. 한 번의 FMA로 계산할 수 있다.

일반 구현

일반 연산으로 구현하면 다음과 같다.

__global__ void kernel_lerp(float* z, const float t,
                             const float* x, const float* y, unsigned n) {
    unsigned i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n) {
        z[i] = (1.0F - t) * x[i] + t * y[i];
    }
}

__global__ void kernel_lerp(float* z, const float t,
                             const float* x, const float* y, unsigned n) {
    unsigned i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n) {
        z[i] = (1.0F - t) * x[i] + t * y[i];
    }
}

곱셈 2회, 덧셈 2회, 뺄셈 1회가 필요하다. 총 5개 연산이다.

FMA 최적화

FMA를 적용하면 다음과 같이 바뀐다.

__global__ void kernel_lerp_fma(float* z, const float t,
                                 const float* x, const float* y, unsigned n) {
    unsigned i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n) {
        float diff = y[i] - x[i];
        z[i] = fmaf(t, diff, x[i]);
    }
}

__global__ void kernel_lerp_fma(float* z, const float t,
                                 const float* x, const float* y, unsigned n) {
    unsigned i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n) {
        float diff = y[i] - x[i];
        z[i] = fmaf(t, diff, x[i]);
    }
}

v₀ + t * (v₁ - v₀) 형태로 계산한다. 뺄셈 1회, FMA 1회로 줄었다.

더 최적화하면 다음과 같다.

z[i] = fmaf(t, y[i], fmaf(-t, x[i], x[i]));

z[i] = fmaf(t, y[i], fmaf(-t, x[i], x[i]));

안쪽 FMA는 x[i] - t * x[i] = x[i] * (1 - t)를 계산한다. 바깥쪽 FMA는 그 결과에 t * y[i]를 더한다.

전체적으로 (1 - t) * x[i] + t * y[i]가 된다. FMA 2개로 계산 완료다.

성능 비교

512M 원소 처리 시 측정 결과는 다음과 같다.

plain

일반 연산: 15,419 usec
FMA 버전: 15,415 usec

일반 연산: 15,419 usec
FMA 버전: 15,415 usec

성능 차이는 0.026%로 미미하다. 컴파일러 최적화 덕분이다.

하지만 정밀도는 향상된다. 반올림 오차가 줄어들어 수치 안정성이 높아진다.

다음 단계로

여기까지 AXPY와 FMA의 핵심을 배웠다. AXPY로 계산 강도 개념을 배웠다. 같은 데이터로 더 많은 계산을 수행하면 GPU의 장점이 커진다. FMA는 현대 GPU의 핵심 명령어로, 딥러닝과 과학 계산에 필수적이다.

다음 장에서는 GPU 내부 동작을 더 깊이 이해한다. Warp 스케줄링과 SM 구조를 배우면 왜 병렬 처리가 효율적인지 명확히 알 수 있다.

그 다음에는 2D 데이터인 행렬 연산으로 확장한다. 행렬 곱셈은 계산 강도가 매우 높아 GPU의 진가를 발휘하기에 이상적이다.

왜 이걸 배우는가?

AXPY: 메모리 병목 극복하기

벡터 덧셈의 실망스러운 결과

해결책: 계산 강도를 높이자

AXPY 루틴이란

CPU 버전 구현

순차 처리 코드

성능 측정

처리량 계산

CUDA 병렬 구현

커널 함수 설계

성능 측정 결과

처리량 분석

FMA: 하드웨어 가속의 비밀

왜 FMA가 필요한가?

왜 명령어가 2개에서 1개로 줄어드나?

Fused Multiply-Add의 장점

반올림 오차 비교

CUDA FMA 함수

FMA를 사용한 SAXPY

커널 수정

성능 측정

FMA 활용 예제: 선형 보간

선형 보간이란

수식 변형

일반 구현

FMA 최적화

성능 비교

다음 단계로

관련 문서