GHCデバッグ日誌
概要
GHC 9.0.1-alpha1で自作ライブラリーをテストしていたらGHC付属ライブラリーのバグを見つけた。再現性がないためにデバッグに苦労したが、数日かけて原因を突き止めた。
発端
GHC 9.0.1-alpha1が出たので自作ライブラリーをテストした。しかし、fusedMultiplyAddに関するテストが一部失敗する。QuickCheckで使う乱数のシードを固定しても、毎回異なる値で失敗する。
Haskellは純粋な言語で、同じ計算は同じ結果を返すはずなので、毎回異なる値で失敗する再現性のないバグというのはかなり嫌なこと（メモリ安全性の破壊）が起きていることを意味する。
問題のコードはざっくり抜き出すとこんな感じである：
-- FMAのリファレンス実装
fusedMultiplyAdd_viaRational :: (RealFloat a) => a -> a -> a -> a
fusedMultiplyAdd_viaRational x y z
| isFinite x && isFinite y && isFinite z =
case toRational x * toRational y + toRational z of
0 -> x * y + z
r -> fromRational r
| isFinite x && isFinite y = z + z -- x * is finite, but z is Infinity or NaN
| otherwise = x * y + z -- either x or y is Infinity or NaN
checkFMA :: (RealFloat a, Show a, Arbitrary a, Random a) => String -> (a -> a -> a -> a) -> [(a, a, a, a)] -> Spec
checkFMA name f cases = do
prop name $ forAllFloats3 $ \a b c -> do
-- sameFloatP は 0 の符号を区別し、NaN同士を区別しない比較関数
f a b c `sameFloatP` fusedMultiplyAdd_viaRational a b c
testSpecialValues name f cases -- 明示的に与えられた例についてテストする
spec :: Spec
spec = modifyMaxSuccess (* 1000) $ do
describe "Double" $ do
checkFMA "fusedMultiplyAdd (default)" fusedMultiplyAdd casesForDouble
checkFMA "fusedMultiplyAdd (via Rational)" fusedMultiplyAdd_viaRational casesForDouble
describe "Float" $ do
checkFMA "fusedMultiplyAdd (default)" fusedMultiplyAdd casesForFloat
checkFMA "fusedMultiplyAdd (via Rational)" fusedMultiplyAdd_viaRational casesForFloat
失敗の内容は、例えば次のようなものだ：
test/FMASpec.hs:87:3:
1) FMA.Float fusedMultiplyAdd (default)
Falsified (after 31878 tests):
0x1.189192p6
0x1.6e6c4p-5
0x1.6d6de2p5
0x1.86874ep5 =/= 0x1.8p5
To rerun use: --match "/FMA/Float/fusedMultiplyAdd (default)/"
Randomized with seed 324858585
test/FMASpec.hs:87:3:
1) FMA.Float fusedMultiplyAdd (via Rational)
Falsified (after 61282 tests):
0x1.17fb68p126
-0x1.797e3p-126
-0x1.00a96cp5
-0x1.08p5 =/= -0x1.0d9046p5
To rerun use: --match "/FMA/Float/fusedMultiplyAdd (via Rational)/"
Randomized with seed 783038111
単体で動くコードとして抜き出すと
#!/usr/bin/env cabal
{- cabal:
build-depends: base, QuickCheck
-}
{-# LANGUAGE TypeApplications #-}
import Numeric
import Test.QuickCheck
isFinite :: RealFloat a => a -> Bool
isFinite x = not (isNaN x) && not (isInfinite x)
-- 十六進で印字するためのnewtype wrapper
newtype HexFloat a = HexFloat { getHexFloat :: a }
deriving (Eq, Ord)
instance RealFloat a => Show (HexFloat a) where
showsPrec _ (HexFloat x) = showHFloat x
-- 0の符号や無限大、NaNの扱いを無視した単純なFMA
fusedMultiplyAdd_simple :: (RealFloat a) => a -> a -> a -> a
fusedMultiplyAdd_simple x y z
= fromRational (toRational x * toRational y + toRational z)
-- 与えられた関数が fusedMultiplyAdd_simple と一致するかテストする
{-# NOINLINE prop_fma_simple #-}
prop_fma_simple :: (RealFloat a, Show a) => (a -> a -> a -> a) -> HexFloat a -> HexFloat a -> HexFloat a -> Property
prop_fma_simple f (HexFloat x) (HexFloat y) (HexFloat z)
= HexFloat (f x y z) === HexFloat (fusedMultiplyAdd_simple x y z)
main :: IO ()
main = do
quickCheck $ withMaxSuccess 1000000 $
forAll (HexFloat <$> (arbitrary `suchThat` isFinite)) $ \x ->
forAll (HexFloat <$> (arbitrary `suchThat` isFinite)) $ \y ->
forAll (HexFloat <$> (arbitrary `suchThat` isFinite)) $ \z ->
prop_fma_simple @Float fusedMultiplyAdd_simple x y z
となるだろうか。同じ関数（fusedMultiplyAdd_simple）を同じ引数に適用してもこの問題は発生する。
こちらのコードは
$ cabal v2-run --with-compiler=/opt/ghc-9.0.1-alpha1/bin/ghc Test.hs
という感じで実行できて、出力は
*** Failed! Falsified (after 215548 tests):
0x1.3fa1c2p4
0x1.32a8cep3
0x1.36cf48p4
0x1.a5bc06p7 /= 0x1.a4p7
*** Failed! Falsified (after 54278 tests):
0x1.dd4246p5
0x1.49e178p5
-0x1.785d22p5
0x1.2d9dccp11 /= 0x1.2d8p11
という風になる。
観察
まずは条件を変えて色々と観察してみる。再現条件をなるべく特定したい。うまくいけばこの段階で原因を推測できる。
こういう観察と推測は探偵っぽくて楽しい。
結果は以下のようになった：
- FMAで発生する。FMAの実装は何種類か用意したけど関係なさそう。浮動小数点数の計算に頼らない「途中の計算を
Rationalで行う」実装でも誤った答えを返すことがある。 - 二項の足し算や掛け算では再現しない。FMA程度の複雑さが必要なのだろうか。
- 再現性がない。乱数のシードを揃えても違うエラーになる。
-
Doubleでは発生せず、Floatで発生する。 -
Floatの下位ビットが切り捨てられるような挙動をしている。 - 最適化を切ると発生しない。
- 特殊化・インライン化はしない方が再現しやすい気がする。
- LLVMバックエンドでは発生しない気がする。
- 古いGHC (8.10.2)では再現しない。
- 元々macOSでテストしていたが、同じ物理マシンで動かしている仮想マシン内のUbuntuでも再現した。別マシンのWindowsではstackによるGHC 9.0.1-alpha1のインストールに失敗したのでテストを諦めた。
- ラズパイ（AArch64）では発生しない気がする。→気のせいだった
この中には、結果的には間違っていた推測も含まれている。まあ確率的に発生する事象だから仕方がない。
推測
原因となる部分は、大まかにいくつか考えられる。
- 自作ライブラリー
- GHC付属ライブラリー（baseやghc-bignum）
- GHCのコード生成器、特にネイティブコード生成器 (NCG)
自作ライブラリーのバグは、GHC 9.0.1-alpha1でのテストですでに1個発見していた。しかし、今回のケースではコードを切り出して単体のプログラム（前掲）として動かしても再現するので、自前のコードのバグという可能性は低い。
となるとGHC側のバグの可能性が濃厚ということになる。GHC 9.0では線形型の導入が話題だが、フロントエンド側の変化でこういうメモリ破壊系のバグが起きる可能性は低い。どちらかというとバックエンド寄りの問題だろう。
GHC 9.0でのバックエンド寄りでは、Integer の新実装であるghc-bignumが最大の変化だろう。新しいライブラリーにはバグがつきものである。
コード生成器（NCG）の浮動小数点数周りでは最近「x86ターゲットで常にSSE2を使う（コンパイラー内部にあったx87 FPU向けのコードを一掃する）」という変化があったが、それはGHC 8.10ですでに実装済みなので、GHC 9.0で新たに問題を起こしている可能性は低い。
こう書くと「ghc-bignumがバグってそう」と考えるのが自然かもしれないが、筆者はこの時点でコード生成器のバグの可能性を捨てきれずに、「NCGがバグってたら相当嫌だな（デバッグしたくない）」と怯えていた。
デバッグ
まずはコード生成器のバグの可能性を検証したい。ただ、 -ddump-*** 系のオプションで生成される中間コードやアセンブリを見るのはしんどい。気休めにGHC本体のデバッグ用のオプション -dcore-lint -dstg-lint -dcmm-lint -fllvm-fill-undef-with-garbage を付けてみたが、特に効果はなかった。GHCに git bisect をかけるのは、GHCのビルドに時間がかかるのでやりたくない。
結局、コード生成器のバグの可能性を追求するのは諦め、コード生成器のバグの可能性と付属ライブラリー(baseやghc-bignum)のバグの可能性の両方を考慮しながらデバッグを進めることにした。
ここで、もう一度基本（観察）に立ち返る。「Float の下位ビットが切り捨てられている」と書いたが、切り捨てられるビット数を数えてみるとどうも中途半端（13ビットとか）である。メモリ破壊ならバイト単位（8ビット単位）で発生するのが自然だし、レジスターの誤使用であれば指数部も含めてガッツリ変な値になるはずだ。そうすると、浮動小数点数周りのコード生成バグという可能性は低いのではないか。
あとは Float 型に関する関数のバグの可能性を考える。テストコードで使っている Float の関数は toRational と fromRational くらいで、怪しいのは fromRational と、その内部で使われている encodeFloat くらいだ。
というわけで Float 型を自前でラップした独自の型を作って fromRational と encodeFloat のをフックし、sanity checkを仕込む。具体的には次のような感じだ：
newtype FloatX = FloatX Float
deriving (Eq, Ord, Show, Arbitrary, Random, Num, Floating, Real, RealFrac)
instance Fractional FloatX where
{-# NOINLINE fromRational #-}
fromRational x = let FloatX a = fromRat x
b = fromRational x
in if a == b then
FloatX a
else
error $ "fromRational/FloatX: mismatch " ++ showHFloat a (" " ++ showHFloat b "")
(/) = coerce ((/) @Float)
recip = coerce (recip @Float)
-- Cで実装したやつ
foreign import ccall unsafe my_encodeFloat :: Int32 -> Int32 -> Float
foreign import ccall unsafe my_encodeFloatU :: Int32 -> Int32 -> Word32
instance RealFloat FloatX where
floatRadix = coerce (floatRadix @Float)
floatDigits = coerce (floatDigits @Float)
floatRange = coerce (floatRange @Float)
decodeFloat = coerce (decodeFloat @Float)
{-# NOINLINE encodeFloat #-}
encodeFloat m e = if abs m > 2^24 then
error "encodeFloat/FloatX: Mantissa too large"
else if integerCheck m then
let r = my_encodeFloat (fromInteger m) (fromIntegral e)
w = my_encodeFloatU (fromInteger m) (fromIntegral e)
r' = castWord32ToFloat w
in if r == r' then
if toRational m * 2^^e == toRational r then
FloatX r
else
error $ "encodeFloat/FloatX: Wrong value: " ++ showHFloat r " " ++ show (m,e)
else
error $ "encodeFloat/FloatX: FFI mismatch: " ++ showHFloat r " vs " ++ showHFloat r' ""
else
error "encodeFloat/FloatX: Invalid integer"
isNaN = coerce (isNaN @Float)
isInfinite = coerce (isInfinite @Float)
isDenormalized = coerce (isDenormalized @Float)
isNegativeZero = coerce (isNegativeZero @Float)
isIEEE = coerce (isIEEE @Float)
atan2 = coerce (atan2 @Float)
しかし、埋め込んだチェックに引っかかることはなかった。
とはいえ、これで浮動小数点数周りの問題という可能性を小さくできた。問題があるとすれば、 encodeFloat や fromRational に値が渡ってくる前、有理数や多倍長整数の演算でバグってる可能性が濃厚である。（やっぱりghc-bignum関係じゃん！）
コードを改良して fromRational に渡される有理数の値を表示させてみると、ビンゴ。明らかに分母分子が小さくなっている：
-- テスト失敗時に fromRational への引数を表示させる
prop_fma_simple f (HexFloat x) (HexFloat y) (HexFloat z)
= let r = toRational x * toRational y + toRational z
in counterexample (show r) $
HexFloat (f x y z) === HexFloat (fromRational r)
{- 実行例：
*** Failed! Falsified (after 238018 tests):
0x1.9c55dp3
0x1.3fdddap3
0x1.d63dccp3
143 % 1 ←このFMAの計算結果がきっかり 143 になるはずがない！
0x1.1efdfep7 /= 0x1.1ep7
-}
さらに調べると、 toRational x * toRational y と toRational z はそれっぽい値になっているのに、 x * y + z は明らかにおかしい値が計算されている。
prop_fma_simple f (HexFloat x) (HexFloat y) (HexFloat z)
= let xy' = toRational x * toRational y
z' = toRational z
r = xy' + z'
in counterexample ("xy' = " ++ show xy') $
counterexample ("z' = " ++ show z') $
counterexample ("r = " ++ show r) $
HexFloat (f x y z) === HexFloat (fromRational r)
{- 実行例：
*** Failed! Falsified (after 68474 tests):
0x1.680692p5
-0x1.12d7aap2
-0x1.0f6604p5
xy' = (-106247109426877) % 549755813888 ←それっぽい値
z' = (-4446593) % 131072 ←それっぽい値
r = (-227) % 1 ←計算結果がこんな単純な値になるはずがない！
-0x1.c65fd6p7 /= -0x1.c6p7
-}
つまり、 Rational の (+) がバグっている。
GHCでは Rational 型はbaseパッケージの GHC.Real モジュールで定義されている。GHC.Real から有理数演算の実装をコピーして、さらに検証したところ、多倍長整数の足し算や掛け算は問題なさそうだが、分母分子をそのGCDで割る部分 reduce がおかしい、ということが判明した。
prop_fma_simple f (HexFloat x) (HexFloat y) (HexFloat z)
= let xy'@(n :% d) = toRational x * toRational y
z'@(n' :% d') = toRational z
n'' = n * d' + n' * d -- xy' + z' を既約分数にする前の分子
d'' = d * d' -- xy' + z' を既約分数にする前の分母
r = reduce n'' d'' -- xy' + z'
in counterexample ("xy' = " ++ show xy') $
counterexample ("z' = " ++ show z') $
counterexample ("n'' = " ++ show n'') $
counterexample ("d'' = " ++ show d'') $
counterexample ("r = " ++ show r) $
HexFloat (f x y z) === HexFloat (fromRational r)
{- 実行例：
*** Failed! Falsified (after 33464 tests):
-0x1.691a2cp2
-0x1.e6ff24p4
0x1.1ed1b8p3
xy' = 47205871654947 % 274877906944
z' = 2349623 % 262144
n'' = 13020595471461908480 ← それっぽい値
d'' = 72057594037927936 ← それっぽい値（これは 2^56）
r = 180 % 1 ← おかしい（これは n'' を d'' で割って小数部分を切り捨てた値）
0x1.6964e6p7 /= 0x1.68p7
ちなみに、 gcd 13020595471461908480 72057594037927936 の値は 262144 で、
r の正しい値は 49669629941795 % 274877906944
-}
GCDで不正に大きな値が計算されているということであれば、「浮動小数点数の下位ビットが切り捨てられているように見える」という現象を説明できる。（今回のケースでは分母は2の冪乗である）
補足： Integer や Rational の内部表現について
ここで Integer や Rational の内部表現を確認しておこう。ただし、1ワードは64ビットとする。
Integer の実体は、GHC付属ライブラリーのghc-bignumで定義されている。GHC 8.10まではinteger-gmpやinteger-simpleで定義されていたが、GHC 9.0ではghc-bignumに統合される。
Integer 等の多倍長整数型の定義は次のようになる。ghc-bignumの導入によりデータ構築子の名前が多少変わったが、基本的にはinteger-gmpのものと同等である。
-- GHC.Num.Integer
data Integer = IS !Int# -- 値が Int で表現できる場合、つまり -2^63 以上 2^63 未満の場合
| IP !BigNat# -- 値が 2^63 以上の場合
| IN !BigNat# -- 値が -2^63 よりも小さい場合
-- GHC.Num.Natural
data Natural = NS !Word# -- 値が Word で表現できる場合、つまり 2^64 未満の場合
| NS !BigNat# -- 値が 2^64 以上の場合
-- GHC.Num.BigNat
type BigNat# = WordArray#
-- GHC.Num.WordArray
type WordArray# = ByteArray# -- ByteArray# はGHCでバイト列を表すプリミティブな型。GHC.Exts参照
BigNat# というのはいかにも多倍長といった感じの型で、二進数で表した整数を64ビット単位の配列としてヒープ上に保持している。GMPの内部表現と互換性があり、GMPの関数をそのまま適用できるようになっている。
Integer や Natural は多倍長と固定長のハイブリッドといった感じで、値が1ワードに収まる場合はデータ構築子に値をそのまま埋め込み、収まらない場合は BigNat# を使って表現する。
有理数を表す Rational 型だが、これはLanguage Reportにもあるように Ratio Integer のエイリアスである。Ratio というのはLanguage Reportでは Data.Ratio モジュールで定義された型だが、GHCにおいて実際の定義は GHC.Real にある。
-- GHC.Real
data Ratio a = !a :% !a -- 不変条件：既約分数として保持する。分母は正
type Rational = Ratio Integer
Data.Ratio では Ratio のデータ構築子を公開しておらず、分子分母にアクセスするには (%) や numerator, denominator 等の関数を使う必要があるが、 GHC.Real からは普通にデータ構築子 (:%) が公開されている。
GHC.Real には reduce :: Integral a => a -> a -> Ratio a という内部的な関数があって、分子および分母（正）を与えると既約分数化（GCDで割る）して返す。(%) とほぼ同じだが、 (%) の引数は分母が負であっても良いという点が異なる。
reduce :: Integral a => a -> a -> Ratio a
reduce _ 0 = error "divide by zero"
reduce x y = (x `quot` d) :% (y `quot` d)
where d = gcd x y
有理数の足し算 (+) は reduce を使って
instance Integral a => Num (Ratio a) where
(x :% y) + (x' :% y') = reduce (x * y' + x' * y) (y * y')
という風に定義されている。
integerGcd のデバッグ
reduce では標準の gcd 関数を使ってGCDを計算しているが、 gcd :: Integral a => a -> a -> a は a = Integer の時に書き換え規則で integerGcd :: Integer -> Integer -> Integer という専用の関数に書き換えられるようになっている（GHC.Real 参照）。
integerGcd がバグっているとすれば、「最適化を切ると再現しなくなる」という現象に納得がいく。最適化が無効の状態では書き換え規則は発動しないからだ。
というわけで、 integerGcd を追求する。
#!/usr/bin/env cabal
{- cabal:
build-depends: base, ghc-bignum, QuickCheck
-}
import GHC.Num.Integer (integerGcd)
import Test.QuickCheck
naiveGcd :: Integer -> Integer -> Integer
naiveGcd x y = gcd' (abs x) (abs y)
where
gcd' x 0 = x
gcd' x y = gcd' y (x `rem` y)
main = do
quickCheck $ withMaxSuccess 1000000 $
forAll (choose (2^62, 2^65)) $ \x ->
let y = 16 -- y の値は重要ではなさそうなので適当に固定する
in integerGcd x y === naiveGcd x y
{- 実行例：
*** Failed! Falsified (after 647 tests):
15296509433926188867 ←これは63ビット整数（2^63 以上 2^64 未満）
16 /= 1
-}
どうやら、一方が63ビット整数でもう一方が63ビット未満の小さい数の時に、GCDを計算せずに小さい方の値をそのまま返してしまうことがあるようだ。
整数の範囲が条件、ということであれば「Float で発生するが Double では発生しない」「FMAでは発生するが足し算や掛け算では発生しない」という現象を説明できる。Double では分母が100ビット以上となるし、Float の足し算や掛け算では分子も分母も48ビット程度となるからだ。
ここで integerGcd 関数のソースコードを確認しよう。
integerGcd ではいくつかの自明なケースを除外したあと、
- 両方が
Intで表せる範囲（-2^{63} 2^{63}-1 gcdWord# - 両方が
Intで表せない場合：bigNatGcd - 片方が
Intで表せてもう片方がそうでない場合：bigNatGcdWord#
に処理を投げている。今回のケースでは bigNatGcdWord# が怪しい、ということになる。というわけで、 bigNatGcdWord# のソースコードを見てみた。
bigNatGcdWord# ではいくつかの自明なケースを除いたあと、バックエンドの bignat_gcd_word を呼び出している。GMPバックエンドの場合は bignat_gcd_word はCで実装された integer_gmp_mpn_gcd_1 を呼び出している。integer_gmp_mpn_gcd_1 の実体は gmp_wrappers.c にあり、GMPの mpn_gcd_1 関数を呼び出している。
GMPの mpn_gcd_1 のバグを疑ったが、問題はなさそうだった。呼び出す側の gmp_wrappers.c の integer_gmp_mpn_gcd_1 関数はここ数年変更がなさそうなので、シロだろう。
ソースコード上で関数呼び出しを追いかけてCのコードまで辿り着いたが、どうも深追いしすぎたようだ。問題はCコードを呼び出す、Haskell側の問題だと考えられる。
問題の切り分けのために、 integerGcd, naturalGcd, gcdWord, bigNatGcdWord# のそれぞれをテストしたところ、 integerGcd と bigNatGcdWord# が失敗した。
#!/usr/bin/env cabal
{- cabal:
build-depends: base, ghc-bignum, QuickCheck
-}
{-# LANGUAGE MagicHash #-}
import GHC.Exts
import GHC.Num.BigNat (bigNatFromWord, bigNatGcdWord#, gcdWord)
import GHC.Num.Integer (integerGcd)
import GHC.Num.Natural (Natural, naturalGcd)
import Test.QuickCheck
naiveGcd :: Integer -> Integer -> Integer
naiveGcd x y = gcd' (abs x) (abs y)
where
gcd' x 0 = x
gcd' x y = gcd' y (x `rem` y)
main = do
-- integerGcdのテスト
quickCheck $ withMaxSuccess 1000000 $
forAll (choose (2^63, 2^64-1)) $ \x ->
\y' ->
let y = fromIntegral (y' :: Word)
in integerGcd x y === naiveGcd x y
-- naturalGcdのテスト
quickCheck $ withMaxSuccess 1000000 $
forAll (choose (2^63, 2^64-1)) $ \x ->
\y' ->
let x' = fromInteger x :: Natural
y = fromIntegral (y' :: Word)
in toInteger (naturalGcd x' (fromInteger y)) === naiveGcd x y
-- gcdWordのテスト
quickCheck $ withMaxSuccess 1000000 $
forAll (choose (2^63, 2^64-1)) $ \x ->
\y' -> y' >= 2 ==>
let x' = fromInteger x :: Word
y = fromIntegral y'
in toInteger (gcdWord x' y') === naiveGcd x y
-- bigNatGcdWord# のテスト
quickCheck $ withMaxSuccess 1000000 $
forAll (choose (2^63, 2^64-1)) $ \x ->
\y'@(W# y#) -> y' >= 2 ==>
let x' = bigNatFromWord (fromInteger x)
y = fromIntegral y'
in toInteger (W# (bigNatGcdWord# x' y#)) === naiveGcd x y
{- 実行例：
*** Failed! Falsified (after 1049 tests):
14286959948854507901
40
40 /= 1
+++ OK, passed 1000000 tests.
+++ OK, passed 1000000 tests; 442425 discarded.
*** Failed! Falsified (after 318297 tests):
11289483210199994043
11
11 /= 1
-}
bigNatGcdWord# がバグっているのは確実なようだ。もう一度ソースコードを眺めると、 bigNatCompareWord# という関数を呼び出しているのに気づいた。この関数が EQ を返した場合、GCDの計算は行わずに引数の b :: Word# をそのまま返している。
-- GHC 9.0.1-alpha1 の src/GHC/Num/BigNat.hs より抜粋：
bigNatGcdWord# :: BigNat# -> Word# -> Word#
bigNatGcdWord# a b
| ...略...
| True = case bigNatCompareWord# a b of
EQ -> b
_ -> bignat_gcd_word a b
もしも、何らかの理由で bigNatCompareWord# が誤って EQ を返しているとすれば、 bigNatGcdWord# のバグを説明できる。というわけで、 bigNatCompareWord# を見てみよう。
-- GHC 9.0.1-alpha1 の src/GHC/Num/BigNat.hs より抜粋：
bigNatCompareWord# :: BigNat# -> Word# -> Ordering
bigNatCompareWord# a b
| bigNatIsZero a = cmpW# 0## b
| isTrue# (wordArraySize# a ># 1#) = GT
| True
= cmpW# (indexWordArray# a 1#) b -- インデックスに注目！！！
はい。
wordArraySize# a が 1 の場合に配列の 1 番目（0 始まり）の要素にアクセスしています。
これは明らかにバグです。メモリの範囲外アクセスです。
この原因を突き止めるのには数日かかったが、結果として修正すべきなのはたった「1文字」だった。プログラミングでこういうことはよくあるものだろう。ソースコードは書かれる時間よりも読まれる時間の方が長い。
報告
適当な再現用コードをでっち上げてGHCのIssuesに報告した。厳密なことを言えば bigNatCompareWord# だけで良いのだが、経緯を考えて integerGcd のテストも含めた。このバグは再現性がないことが特徴なので、QuickCheckへの依存は外せなかった。
{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE MagicHash #-}
import GHC.Exts
import GHC.Num.BigNat (bigNatCompareWord#, bigNatFromWord#)
import GHC.Num.Integer (integerGcd)
import Test.QuickCheck
main :: IO ()
main = do
quickCheck $ withMaxSuccess 1000000 $
forAll (choose (2^63, 2^64-1)) $ \x ->
\y ->
integerGcd x (toInteger y) === toInteger (gcd (fromInteger x) y :: Word)
quickCheck $ withMaxSuccess 1000000 $
\x@(W# x#) -> let !x' = bigNatFromWord# x# in
\y@(W# y#) ->
bigNatCompareWord# x' y# === compare x y
確認
バグの原因と思われるものを発見したら、それを修正して実際に直ることも確かめなければならない。
GHCのチェックアウトで
$ ./boot
$ ./configure LLC=llc-mp-9.0 OPT=opt-mp-9.0 --with-intree-gmp
$ hadrian/build --flavour=quick -j3
という感じでビルドを行い、 _build/stage1/bin/ghc を使ってテストする。cabal v2-runでプロジェクトを作らずにやろうとすると
Could not find module ‘Prelude’
Perhaps you haven't installed the profiling libraries for package ‘base-4.15.0.0’?
Use -v (or `:set -v` in ghci) to see a list of the files searched for.
とか言われてウザかったがそれは別の話。
修正後のGHCでは無事、再現しなくなっていることを確認できた。（元々の自作ライブラリーでは確認していないが、今後出るであろうGHC 9.0.1のalpha2かbetaかrcでテストできると見込んでいる）
雑感
安定したGHCを使って普通にHaskellでプログラミングする分には、こんなにデバッグに苦労することはない。配列に普通にアクセスする際には範囲チェックが行われるし、自前のコードならassert等を挟むのも容易だ。今回はGHC付属ライブラリーのバグということで、運が悪かった。
特定の範囲の Integer で発生するバグといえば、昔reflectionパッケージで「値がunsafeCoerce で型を意図的に無視した結果、型の不変条件が破壊されるやつだった。
今回は運よく筆者の独力でバグの原因の特定まで至り、原因（実質、修正方法）を添えてのバグ報告という形になった。OSSへバグ報告する際はなるべく原因を詳細に特定してから報告する、というのは筆者が心掛けていることである。
しかし、エスパー能力に自信のない者はこういう厄介なバグに遭遇したら自力で原因を特定することは適当なところで諦めて、再現コードを添えて開発元へ報告することを優先すべきだろう。ソースコードについては大抵は開発元の方が詳しい。
教訓（追記）
今回の件から頑張って教訓を捻り出すとすれば、
- テストは大事。今回の件は普段からテストを書いていたおかげで発見できた。
- 新しいものにはバグがつきものなので、alpha版やbeta版は積極的にテストしよう。
となるだろうか。
Discussion