하나를 들이고 여섯을 재다 — Go의 SIMD 실험에서 Mere가 가져온 것

9월에 Go 팀이 SIMD 실험에 관한 글을 썼다. Go 1.26은 amd64용 SIMD를 GOEXPERIMENT=simd 아래에서 내놓았고, 그 API는 CPU의 모양을 따른 것이었다. Go 1.27은 이를 arm64 NEON과 WebAssembly까지 아우르는, CPU에 의존하지 않는 simd 패키지로 만들었다. 아직 실험 기능이고, 어떻게 생각해서 그렇게 했는지가 꼼꼼하게 적혀 있다.

Mere는 내가 나를 위해 만드는 작은 언어다. 인터프리터와, 컴파일하는 백엔드 네 개 —— C, LLVM IR, WebAssembly, RISC-V 기계어 —— 를 갖고 있고, SIMD도 이미 있었다. 평범한 루프를 clang이 벡터화할 수 있게 해 주는 컴파일러 패스와, 128비트 타입 세 개(f64x2, f32x4, u8x16)다 (수치는 이전 글에 있다). 그래서 Go의 글이 던진 질문은 “Mere에 SIMD가 필요한가”가 아니라 “Go의 어떤 판단을 들이면 Mere의 SIMD가 나아지는가”였다.

Go의 설계를 항목 하나씩 보고, 결정하기 전에 각각 측정했다. 들어온 것은 하나뿐이다. 이 글은 그 하나와, 들어오지 않은 여섯과, 확인하는 도중에 나온 것들의 이야기다. 수치는 모두 arm64(Apple clang 21, node 24)에서 측정했다.

표부터

Go의 기능	Mere	결정한 근거
MulAdd(융합 곱셈-덧셈)	들였다: fma, f64x2_fma	레인으로 쓴 정점 변환이 140 ms → 100 ms
루프 끝자락을 다루는 LoadVPart / StorePart	들이지 않음	이를 쓰는 유일한 프로그램에서의 상한이 0.16 s → 0.15 s
레인 폭마다의 마스크 타입과 IfElse	아직 들이지 않음	select가 필요한 곳이 없고, 여분의 명령은 타입을 늘리지 않고 지울 수 있다
“반드시 돌아가게” 하는 에뮬레이션	원칙으로 삼지 않음	정확히 재현할 수 있는 것만 소프트웨어로 계산한다. 지금은 Wasm의 fma뿐
타입에서 폭을 없애기	들이지 않음	Mere의 타깃은 모두 128비트 폭
feature variants(1.28 예정)	들이지 않음	실행 시에 CPU를 가려내는 문제로, 미리 빌드하는 Mere에는 생기지 않는다
합산 계열(Go 1.27에는 없음)	이미 있음	reduce_add, any_true, first_true

처음 추천은 틀렸다

글을 읽은 시점에는 끝자락 읽기를 1순위로 추천할 생각이었다. Go의 LoadFloat32sPart(x[i:])는 루프의 마지막 짧은 벡터를 읽어 들이고, 모자란 부분을 0으로 채우고, 몇 요소를 읽었는지 돌려준다. 그래서 루프의 끝까지 레인 그대로 처리할 수 있다. Mere의 바이트 단위 프로그램은 마지막 몇 바이트를 위해 루프를 스칼라로 이중으로 쓰고 있을 테니, 그것이 사라지리라 —— 고 생각했다.

만들기 전에 확인해 보니, 그 전제는 파일 하나를 연 것만으로 무너졌다. UTF-8 검증기는 끝자락을 이미 사용자 코드에서 0으로 채워 한 블록으로 처리하고 있었다. 0은 ASCII라서 채워도 답이 바뀌지 않는다. 끝자락과 16바이트 미만의 줄을 스칼라로 처리하는 것은 mgrep(Mere로 쓴 grep)의 literal prefilter 하나뿐이었다.

그래서 그것으로 상한을 쟀다. 81 MB의 소스 코드 줄(170만 줄, 평균 49바이트)에서 스칼라 경로를 지나는 바이트는 전체의 14%다.

mgrep 변형	첫 바이트가 드묾	첫 바이트가 흔함
지금 그대로(끝자락은 스칼라)	0.16 s	0.31 s
끝자락을 사용자 코드에서 0으로 채움	0.26 s	0.38 s
끝자락을 통째로 생략 —— 답은 틀림, 상한	0.15 s	0.27 s

내장 끝자락 읽기로 도달할 수 있는 것은 잘해야 마지막 행까지다. 게다가 그 행은 끝자락에서의 대조까지 생략하고 있어서 실제보다 좋아 보인다. 1순위로 추천하려던 기능은, 그것을 쓸 수 있는 유일한 프로그램에게 6% 정도의 가치였다.

이 측정으로 mgrep이 정말로 손해를 보는 곳도 알게 됐다. Mere의 SIMD 읽기는 bytes를 읽고, 줄은 str 이라서 줄마다 한 번씩 복사한다. 보통 길이의 줄에서는 그 복사가 얻을 수 있는 효과의 4할쯤을 먹고 있었다. 이것은 Mere 자신의 문제로, Go의 설계가 답해 주는 이야기가 아니어서 따로 이름을 붙여 남겨 두었다.

마스크와 에뮬레이션, 짧게

Mere에서는 두 u8x16을 비교하면 일치한 레인이 0xFF인 u8x16이 돌아온다. NEON과 WebAssembly의 방식이다. RISC-V에서는 비교가 마스크 레지스터에 비트 마스크를 만들기 때문에, 그것을 바이트열로 되돌리는 데 명령 두 개가 든다. WebAssembly에서는 처음 일치한 레인을 물으면, 진짜 마스크 타입이 있을 때보다 명령이 두 개 더 든다. 한편 UTF-8 검증기는 비교 결과를 데이터로 쓴다. 1과 AND를 취해 연속 바이트를 센다. 마스크를 별도 타입으로 만들면 거기에 변환이 필요하다. 그리고 select가 필요한 프로그램은 아직 없다. 자주 쓰는 모양을 백엔드에서 패턴으로 잡으면 타입을 늘리지 않고 여분의 명령을 지울 수 있다. 이것은 select가 필요한 JSON 스캐너나 base64 코덱이 나타날 때까지 기다린다.

Go는 하드웨어에 없는 SIMD 연산을 모두 에뮬레이션하므로, simd를 쓰는 Go 프로그램은 어디서나 돌아간다. Mere는 내릴 수 없는 연산을 그 백엔드에서 이름을 들어 거부한다. 예를 들어 RISC-V 백엔드는 f64x2를 거부한다. 게다가 Mere에는 컴파일된 모든 백엔드를 대조하는 기준으로 인터프리터가 있다. Go의 GODEBUG=simd=0이 테스트 도구로서 주는 것이다. 원칙으로서의 에뮬레이션은 Mere에는 별로 효과가 없다. 하지만 연산마다의 판단으로서의 에뮬레이션은 중요해졌다. 그 연산이 Mere가 들인 유일한 기능이었다.

융합 곱셈-덧셈

fma a b c는 a * b + c를 반올림 한 번으로 계산한다. 풀어 쓴 a * b + c는 곱 다음과 합 다음에 두 번 반올림한다. 두 답은 마지막 비트에서 다를 수 있고, 곱과 c가 거의 상쇄될 때 가장 다르다. fma 0.1 10.0 (-1.0)은 5.55e-17이고, 0.1 * 10.0 - 1.0은 정확히 0.0이 된다.

Mere는 이것을 일부러 갖지 않고 있었다. C 백엔드는 #pragma STDC FP_CONTRACT OFF를 출력한다. arm64의 clang은 -O2에서 a * b + c를 제멋대로 fma 명령으로 묶기 때문이다. 그런 일이 일어났을 때, 내적이 인터프리터와 최적화한 바이너리에서 다른 비트를 냈다. Mere에는 “최적화로 답을 바꾸면 안 된다”는 규칙이 있다. 묶기를 끄면 그 규칙은 지킬 수 있지만, 곱셈-덧셈의 연쇄를 빠르게 하는 유일한 명령도 내려놓게 된다.

Go의 글은 다른 질문을 세울 계기가 됐다. 융합한 답을 프로그램이 이름으로 부탁할 수 있게 하면 어떨까. 그러면 규칙은 유지된다. a * b + c는 어디서나 두 번 반올림하고, fma는 어디서나 한 번 반올림한다.

한 줄도 쓰기 전에 얼마나 가치가 있는지 쟀다. 정점 변환(정점 20,000개를 4x4 행렬에 통과시키기를 5,000번)에 대해, Mere가 출력한 C의 pragma를 켠 것이다.

커널	pragma OFF(출하 상태)	pragma ON	융합 명령 수
스칼라	0.20 s	0.15 s	0 → 348
f64x2 레인	0.14 s	0.10 s	0 → 72

둘 다 1.4배. 답은 바뀐다 —— 바뀌어야 한다 —— 그런데 융합한 스칼라 판과 융합한 레인 판은 서로 같은 비트를 냈다. fma는 결정적이라서, 융합한 답은 반올림의 우연이 아니라 제대로 정의된 답이 된다. 내장 함수를 쓴 뒤 f64x2_fma로 고쳐 쓴 레인 판은 C 백엔드에서 140 ms 대신 100 ms로 돈다. clang은 레인 fma를 두 개씩 fmla.2d 하나로 묶는다.

fma 명령이 없는 백엔드

다섯 중 넷은 순조로웠다. 인터프리터는 OCaml의 Float.fma, C는 fma(3), LLVM은 llvm.fma를 부른다. RISC-V 백엔드는 배정밀도 연산기가 없어서 이름으로 거부한다. 진짜 판단은 WebAssembly였다. fma 명령이 없고, JavaScript 호스트에도 Math.fma가 없다.

선택지는 둘이었다. Mere가 늘 하듯이 거기서 거부하거나, Go가 하듯이 소프트웨어로 계산하거나. 소프트웨어를 고른 이유는 하나다. fma는 IEEE-754가 올바른 반올림을 요구하는 연산 중 하나로, +, -, *, /, sqrt와 마찬가지로 올바른 답이 딱 하나 있고 그것은 재현할 수 있다. 거부해도 아무도 문제에서 벗어나지 못한다. fma를 쓰고 싶은 프로그램은 WebAssembly 판에서만 a * b + c라고 쓰고, 거기서 다른 답을 얻게 된다. 그것이야말로 pragma가 막으려던 것이다.

처음 만든 소프트 fma는 Mere로 15줄쯤이었다. Boldo와 Melquiond의 기법을 쓴다. 곱을 double 두 개로 정확히 구하고, 더하고, 중간에 한 번만 “홀수로의 반올림”을 한다. 네 백엔드 모두에서 3만 건에 대해 하드웨어와 일치했다. 하지만 그것이 성립하는 것은 중간에 아무것도 오버플로나 언더플로를 일으키지 않는 동안뿐이고, 확인한 것은 2^-100에서 2^100 정도의 범위뿐이었다. 게다가 WebAssembly에서는 부를 때마다 할당을 해서 200만 번쯤에서 메모리가 바닥났다. 대신 출하한 것은 WebAssembly 런타임 안에 정수로 쓴 것으로, musl의 fma.c와 같은 구성이다.

1. 가수 세 개를 54비트로 정규화한다.
2. 108비트 곱을 32비트씩 나눠 만든다. WebAssembly에는 64비트끼리 곱해 128비트를 돌려주는 명령이 없기 때문이다.
3. c를 곱에 맞춘다. 밀려난 비트는 sticky 비트로 모은다.
4. 더하거나 뺀다.
5. 반올림은 정확히 한 번. 상위 63비트를 double로 변환하는 곳이 반올림이 된다. 결과가 비정규화 수가 될 때는 먼저 비정규화 수의 마지막 자리에서 반올림해 두어, 마지막 2의 거듭제곱 배에서 다시 반올림하지 않게 한다.

확인하는 쪽이 쓰는 것보다 오래 걸렸다.

• 먼저 C로, 하드웨어와 2억 5,000만 건. 특수값, 비정규화 수, 지수의 위아래 끝, 상쇄, 타이에서 만든 입력으로, 불일치는 0.
• 분기를 하나씩 일부러 망가뜨렸다. 열 가지 중 여덟 가지에서 불일치가 나왔다. 남은 두 가지는 어떻게 해도 답을 바꿀 수 없었다. 이유도 알아냈다. 53비트끼리의 곱 한가운데에 64비트 이상의 0이 늘어서고 그 양쪽에 1이 서 있어야 하는데, 그런 곱은 존재하지 않는다. 그래서 이 둘은 “확인했다”고 하지 않고, “관측할 수 없는” 이유를 글로 남겼다.
• 다음으로 WAT 자체. 출력된 모듈에서 잘라내 node 위에서 하드웨어의 답과 5,000만 건을 더 대조했다. 불일치는 0. node에서의 대조 자체도, WAT를 망가뜨려 떨어지는 것을 보고 확인했다.

이 중 하나에는 하마터면 속을 뻔했다. 밀려난 c의 sticky 비트를 떨어뜨리는 망가뜨리기는, 무작위 입력 2,000만 건에서는 통과해 버렸다. 무작위 입력은, 그 비트만이 반올림을 결정하는 유일한 경우인 타이에 거의 닿지 않는다. 일부러 만든 타이 계열을 더하자 —— 정확히 절반에서 반올림되는 홀수 54비트의 정확한 곱과, 그 훨씬 아래에서 방향을 정하는 작은 c —— 같은 망가뜨리기가 2,000만 건 중 10,272건에서 떨어졌다. 타이는 상설 게이트에 넣어 두었다. 컴파일러의 WebAssembly 런타임에서 그 비트를 빼면, 게이트는 계열의 이름을 들며 떨어진다. Wasm differs from fma(3) in: ties.

WebAssembly에서의 비용은 한 번에 약 12 ns. 한 번의 반복에서 fma를 8단 중첩하면 130 ms이고, 융합하지 않는 식은 30 ms였다. 반복마다 float를 할당하는 루프 안에서는 이 차이가 잡음에 묻힌다.

게이트가 찾아낸, fma가 아닌 것

상설 게이트는 여섯 계열의 입력을 인터프리터, 두 최적화 수준의 C, LLVM, WebAssembly에 통과시키고, 각각을 C 라이브러리의 fma(3)와 비트 단위로 비교한다. 이 게이트는 처음 초록이 되기 전에, fma와 무관한 LLVM 백엔드의 버그를 두 개 찾아냈다.

x != x가 NaN을 “NaN이 아니다”라고 했다. IEEE-754는 NaN의 페이로드를 정하지 않으므로, 게이트는 NaN 결과를 하나의 값으로 모아서 비교한다. NaN인지는 r != r로 본다. LLVM 쪽은 하드웨어와 한 계열에서만 어긋났다. NaN이 나오는 유일한 계열이다. 백엔드는 float의 !=를 순서 있는 비교로 하고 있어서, 한쪽이 NaN이면 항상 false가 됐다. 인터프리터와 C와 WebAssembly는 모두 true를 돌려준다. float의 경계 사례를 보는 parity 테스트는 NaN에 대해 ==, <, >는 물었지만, !=는 한 번도 묻지 않았다.

내장 함수의 이름을 지역적으로 묶으면, 그것이 프로그램 전체로 샜다. fma의 parity 테스트는 fma라는 이름의 지역 함수가 그 스코프 안에만 속하는지 확인한다. LLVM에서는 그것이 파일 안의 모든 fma를 가로챘다. 최상위 코드가, 마지막으로 출력한 함수의 “밖으로 끌어올린 안쪽 함수” 표를 남긴 채 출력되고 있었던 것이다. 같은 버그로, 어느 함수가 자기 atan2를 정의하면 atan2 1.0 1.0이 2.0이 됐다. 처음 쓴 회귀 테스트는 버그가 남은 채로 통과했다. 새는 것은 마지막 함수의 표뿐이었기 때문이다. 지금의 테스트는 덮어쓰는 함수를 맨 마지막에 둔다.

!= 버그 옆을 보다가 세 번째도 나왔다. f_min과 f_max에 NaN을 넘기면 인터프리터와 WebAssembly는 NaN을 돌려주고, C와 LLVM은 인자 순서에 따라 다른 쪽 값을 돌려줬다. 이것은 알려진 실패 목록에 넣고 조금 뒤에 고쳤다. 고친 방법은 OCaml 자신의 Float.min을 절차 그대로 모든 백엔드에 옮기는 것이다. 첫 판은 NaN을 만드는 데 a + b를 썼고, 비트를 비교하는 게이트가 잡아냈다. signaling NaN을 넘기면 OCaml은 그대로 돌려주지만, 덧셈은 quiet NaN으로 만든다.

시작하기 전의 나에게 할 말

• 기능을 만들기 전에 상한을 잰다. 끝자락 읽기는 백엔드 다섯 개에 내장 함수를 하나 더해서 6%였다. 끝자락을 통째로 생략한 변형 —— 답이 틀리는 것을 알면서 —— 이 반나절 만에 결론을 냈다.
• 규칙은 “fma를 쓰지 않는다”가 아니라 “최적화로 답을 바꾸지 않는다”였다. pragma가 왜 있는지를 읽은 덕분에 1.4배가 “금지”에서 “쓸 수 있음”으로 바뀌었다.
• 정확히 재현할 수 있는 것은 에뮬레이션하고, 나머지는 거부한다. 이것은 전체 방침이 아니라 연산마다의 판단이고, 처음 그것을 통과한 것이 fma였다.
• 무작위 테스트는 타이에 닿지 않는다. 일부러 만든다.
• 게이트의 수확은 흔히 대상의 옆에 있다. fma는 첫날부터 맞았다. 빨개진 것은 게이트가 결과를 모으는 데 쓴 NaN 판정과, parity 테스트가 스코프를 확인하는 데 쓴 지역 바인딩이었다.

fma에는 테스트와 위의 측정 말고는 아직 실제 사용자가 없다. 가장 자연스러운 것은 작은 3D 엔진의 선형대수로, 그곳은 이미 올바르게 반올림하는 연산에 대해 모든 백엔드가 같은 비트를 내는지 확인하고 있다. 판단과 측정을 모두 적은 자세한 기록은 Mere 저장소의 docs/simd-go-comparison.md에 있다.