再現性の規律¶

生成されたテーブルは科学的な成果物です。同じ入力からは同じビットが出なければなりません — ただし、それは下に書き出した範囲の中での話です。その範囲は仮定ではなく測定によって決めました。コードを 30 % 速くしても出力の末尾数ビットを動かしてしまう変更は改良ではありません — その変更より前に作られたテーブルがソースから再生成できなくなるので、新しいデータセット版を強いることになります。凍結されたアーカイブは「何を出荷したか」の記録としては有効なままです。そうした変更が壊すのは、そのアーカイブと現在のコードとのつながりです。

このページには方針を書きます。手順は CONTRIBUTING にあります。

なぜ末尾のビットが重要なのか¶

Temari が出荷するすべての値 — イオン化形状因子 F(s, E₀) や散乱因子 f_x(s) — は、浮動小数点の和・積・漸化式からなる長い連鎖の末端にあります。浮動小数点演算は 1 演算ごとに丸めるので、和の結果は項を足す順序に依存します。したがって、数学的には同一な 2 つのプログラムが、最後の桁だけ違う数を印字することがあり得ます。そしてコンパイラのフラグ 1 つや「無害な」ループの書き換え 1 つで、その順序は変わってしまいます。

浮動小数点の加算は結合的ではない

ごく普通の 3 つの数を、binary64 (Julia の Float64、Python の float) で 2 通りの順序で足します:

(0.1 + 0.2) + 0.3   →   0.6000000000000001
0.1 + (0.2 + 0.3)   →   0.6

どちらの答えも「15 桁で見れば 0.6」ですが、最後の桁で 1 単位だけ違います。どちらも間違いではありません — 同じ厳密な和に対する 2 通りの丸めです。これが問題のすべてを縮図で示しています: ビット同一な結果とは、出力のすべての値の 64 ビットのパターンが変わらないこと (Julia の Float64 に対する ===。+0.0 と -0.0 すら区別します) であり、それは途中のすべての丸めが同じ順序で起こるときにしか実現しません。

物理学者にとって、最後のビットは意味のあるどんな許容誤差よりもはるかに下にあります。ここでそれが重要になるのは別の理由からです: ビット同一性は、ある変更が検査したすべての出力値に手を触れなかったことを確実に示せる唯一のテストです — それがリファクタリングと物理の変更とを分けるものです。高速化がすべてのビットをそのまま残すなら、何も再検証する必要はありません。1 ビットでも動かしたなら、そこから生成されるテーブルは新しい世代であり、そのように再生成・検査・版付けをしなければなりません。

保証の範囲¶

3 つの命題を、強いものから弱いものへ並べます。強いものはコードが約束することで、残りの 2 つは観測されたことを、測定された数値とともに書いたものです。

命題	状態	何を測ったか
1 つのプロセスの中では、結果はスレッド数に依存しません。同じ入力、同じコード、同じ処理系なら、`-t 1` でも `-t 32` でも同じビットになります (F(s, E₀) エンジン)。	保証されます。	下の E8 の調査につながったフリート測定で、1 から 32 スレッドまで確認しました (`tools/e8_README.md`)。放出電子エネルギーのノード (ε) にわたるスレッド並列ループは互いに素なインデックスにだけ書き込み、ε についての和はその後で、完成したノードごとの値から 1 か所でとります。
2 つのプロセスの間では、自己無撞着場 (SCF) が別の反復で止まることがあります。	散発的に観測されています。リリースされたアーカイブのバイト列とその SHA-256 が正本の記録であり、「いま生成器が印字するもの」ではありません。	dataset-factors v1.0.0 (dt/16 格子) では、認証用の実行と出荷用の実行との間で 86 元素中 34 が別の SCF 反復で止まり、同じやり方で行った 2 回の本番実行の間では 85 中 6 でした。どの差も SCF の停止許容の内側にあり、リリースの許容をはるかに下回ります — データの factors の節を参照してください。
Windows 上の大規模な同時実行フリートで、高い割り当て率が続く状況での、実行間のバイト単位の一致。	保証しません。	まれに 1 行 (1 チャネルの 1 つの E₀ — 1 本の F(s) 曲線) が最後の桁で 1〜2 単位違うことがあります。下の E8 を参照してください。

最後の 2 行は、いい加減でよいという免罪符ではありません: フリート実行は QC 工程で監視されていますし、E8 の不一致は物理的な許容誤差より 6 桁下にあり、SCF の差は停止許容の内側にあります — 2 回の本番実行の間で f_x で見て最大でも 0.22 × B_scf です。ここで B_scf = 9.09 × 10⁻⁹ 電子は、f_x のリリースの許容のうち SCF の停止に割り当てた分です。これらはランタイムが現時点で約束できることについての記述であり、仮定ではなくフォレンジック (事後の追跡調査) によって得られたものです。

SCF の観測は factors 生成器で行いました。2 つのデータセット族はどちらも同じ SCF 層 (l1_atomic.jl) で解いているので、リリースされたアーカイブ — バイト列、ファイルごとの SHA-256、MANIFEST.md — をどちらの族についても正本として扱ってください (F(s, E₀) 生成器では停止反復の差は観測されていません。注意をそちらへ広げるのは、共有している層からの推論であって測定ではありません)。そして再生成は、それに対して 突き合わせて検査するものとして扱い、等しいと仮定しないでください。検査器はこのページの末尾に列挙してあります。

計算コードで禁止していること¶

最初の 4 つは丸めの順序や回数を変え、したがってビットを変えます。最後の 1 つは同じリストに載っている性能の罠です:

総和 (reduction) に対する @simd。コンパイラに和の再結合 (足す順序の組み替え) を許してしまいます。
muladd と fma。乗算と加算を 1 回の丸めに融合します — 素直な式より丸めが 1 回少ないので、結果が変わります。
Base.sum() を手書きループに置き換えること、またはその逆。Base.sum() は内部でペアワイズ総和を使っており、左から右へ足すループとは別のアルゴリズムです。
総和順序のあらゆる変更。共有スカラーに違う順番で加算していく「無害な」ループの並べ替えも含みます。
ntuple に渡すクロージャに、ループ内で再代入される変数を捕まえさせること。変数が box 化され (Core.Box)、性能バグになります — このコードで実際に一度出荷されたことがあります。ビットは動かしませんが、善意の「掃除」が持ち込みがちな類のものなので、このリストに載せています。

許していること¶

キャッシュブロッキングは、タイルが同じスカラーにインデックスの昇順で加算していく限り許します。この制約は妥協ではありません: このプロジェクトで最大の高速化はこの制約のもとで得られました。現在のスタック — 8 レーン SIMD の球 Bessel カーネル、融合した角度パス、動径積分表のループ入れ替え、そして box 化の修正 — は、v3 データセットを生成したコードよりおよそ 11.7 倍速く、しかもその 1 段ごとがビット同一です (比較は 2026-08 に行いました。その後の v4 と v5 のデータセットは高速なコードで生成しています)。段階ごとの記録は性能にあります。

ビット同一性を必ず壊す最適化 (求積格子の変更、部分波和のより早い打ち切り) も却下はしません — 次のデータセット生成へ回すものとして保留し、そのマニフェストで宣言します。

正しさはビット互換性に優先する¶

帰結は 2 つあり、どちらも実際に起きた例があります:

本物のバグは、修正が値を動かすとしても直ちに直します。球 Bessel の 0/0 の修正は単独で出荷されましたが、それは閾値ガード (l0_numerics.jl の J0_MIN) として書かれていたからです: 壊れていた窓の外では命令列が変わらないので、旧コードがもともと間違っていなかった場所ではどこでもビット同一です。この主張自体もテストしました — ガードの閾値を 0.0 にした (つまりガードを無効化した) ビルドを作り、そのビルドが修正前のコードとビット同一であることを確認して、観測されたすべての差がガードが実際に発火したことに帰着できるようにしました (検証)。可能なら、修正はそのような形に整えてください。
非決定性はそれ自体が正しさのバグです。機械の負荷やスレッドの分割に依存する出力を、参照値として保護することは決してありません。バグを含む再現不能なビットを基準として凍結するのは無意味です。

処理系はデータセットの一部¶

Julia 自身のバージョンはデータセットの世代ごとにピン留めし、マニフェストに記録します — F(s, E₀) データセット (v3、v4、v5.0.0、そして現行の v7.0.0) には 1.11.9、 dataset-factors v1.0.0 には 1.12.6 です。コード生成や libm の実装はバージョン間で変わり、このリポジトリに一切変更が無くてもビット同一性を壊し得ます。したがって処理系の更新はテーブルの全再生成と同格に扱い、同じやり方で宣言します。

CI は selftest を Julia 1.11.9 と 1.12 で、Ubuntu と Windows の両方で走らせ、ずれを早期に捕まえます。リリースのゲートは、対象となるデータセット族のピン留めされたバージョンです。juliaup を使っているなら、julia +1.11 で F(s, E₀) のピンが選ばれます。

仮定せず、測る¶

上の方針は、禁止した最適化こそが価値のあるものだったなら高くつくはずです。しかしそうではありません。この作業負荷で測った結果:

変更	効果
内側の動径ループのキャッシュブロッキング	そのループで 2.4 倍 (ビット同一)、全体では 1.13 倍
転置 + `@simd`	ループでは 8.3 倍、しかしビット同一ではなく、全体では 1.27 倍にとどまる → 却下
除算を逆数の乗算に置き換え	1.03 倍 — レイテンシ律速の依存連鎖
`--heap-size-hint`	ゼロ。GC 回数 151 → 151。ライブセットは小さく、GC は割り当て率で駆動される
BigInt 階乗のルックアップテーブル	割り当て 0.777 → 0.775 GB — 仮説は外れ、害は無し

ルールを破ることになったはずの 1 つは、全体で 1.27 倍しか買えませんでした。ルールを守った側は 11.7 倍を買いました。

E8 — 負荷依存の ULP 反転¶

記録しておく価値があります。調査が思いがけないところに行き着いたからです。(ULP = unit in the last place は、隣り合う 2 つの binary64 の数の間隔 — あり得る最小の差です。)

症状。フリート実行で、まれに 1 行が、同じ行を計算し直した結果と最後の桁で 1〜2 単位違いました。

調査。計装した張り込み (stakeout) がその事象を捕捉し、468 件のサイドカー比較と 240 回の単一ノード再生試行を使って、エンジン内部の候補となる機構をひとつずつ検査しました: 総和 (reduction)、LAPACK 固有値ソルバの呼び出し、アライメント依存のループ剥がし (loop peeling) です。すべて容疑が晴れました。

結論。残った整合的な説明は、ランタイムの並行スイープ (concurrent-sweep) ガベージコレクタによる一過性の擾乱です — トラブルシューティングに書いた Windows の GC クラッシュの静かな親戚です。頻度は、-t 4 で GC スレッドが既定のとき 1 行あたり約 10⁻²、-t 2 かつ --gcthreads=1 では 1704 行中 0 でした。振幅は、ただ 1 つの中間値の下位ビットで、10⁻¹⁰ の物理的な許容誤差より 6 桁下です。

とった措置: エンジンには何もしていません。明白な「修正」 — 総和を固定のインデックス順に書き直すこと — では観測された反転を 1 つも防げなかったはずで、しかもビット互換性を犠牲にしていたはずです。代わりに、サイドカーの計装を休止状態のまま出荷コードに残し (環境変数 E8_SIDECAR で目覚め、未設定ならコストはゼロ)、以後の事象が自然に捕捉されるようにし、ビット同一性の保証の範囲を上のとおり書き出しました。

一般原則はこの一件を経ても無傷です — 証拠はこのコードから離れる方向を指していますが、残るのはフォレンジックからの推論であって、証明ではありません。

実務チェックリスト¶

計算コードに手を触れる前に、「前」のスナップショットを取ってください。数チャネルを完全精度でダンプするので、素の diff がそのまま === 比較になります。変更後に diff が空なら、ビットは 1 つも動いていません。

# Snapshot BOTH prescriptions before touching anything. The plain form runs the
# v3 prescription (five channels); --v4 runs the shipping v4 prescription
# (seven channels, including M1 = 3s and M5 = 3d, so the l_init = 2 angular
# path and the kappa-resolved Dirac continuum are actually exercised).
julia +1.11 -t 4 tools/bitident_snapshot.jl      before.txt    # FIRST
julia +1.11 -t 4 tools/bitident_snapshot.jl --v4 before4.txt   # FIRST
# ... change ...
julia +1.11 -t 4 tools/bitident_snapshot.jl      after.txt
julia +1.11 -t 4 tools/bitident_snapshot.jl --v4 after4.txt
diff before.txt after.txt && diff before4.txt after4.txt

julia +1.11 -t auto src/ionization.jl selftest
julia +1.11 -t auto src/ionization.jl refcheck
julia +1.11 -t 1 tools/verify_simd_bessel.jl        # 8-lane Bessel kernel vs scalar
julia +1.11 -t 1 tools/verify_e5_qlane.jl           # radial-integral q lane vs reference
julia +1.11 -t 1 tools/verify_e5_qlane_dirac.jl     # the same on the Dirac (v4 shipping) path
julia +1.11 -t 1 tools/verify_angular_pack.jl       # v4 angular fast path vs the oracle

「前」のスナップショットを先に取る

変更の後からは再構成できません — 旧コードはもう無いからです。

変更が生成済みの F(s, E₀) テーブルに触れたなら、データセット検査器も走らせてください:

julia +1.11 -t auto tools/check_tables.jl <prod_dir> [--eb]

散乱因子データセットに対する同等の検査器は tools/check_factor_tables.jl (F1〜F10 を検査します: 元素集合、メタデータの一様性、s 格子の SHA-256、値の構造、 Mott–Bethe の恒等式、ゲート台帳、loader 規約。--certify-dir を付けると tight な参照に対する SCF 停止誤差 (F8) も、--golden を付けると言語横断の golden ベクトル、すなわち Python loader と Julia loader の 1e-12 での比較 (F9) も検査します。F10 は丸めの寄与、最悪のゲート値、SCF の秒数を記録するだけです) と、Python で書かれた実行可能な仕様です。後者の --negative モードは、各検査が捕まえるべき欠陥を実際に検知することも実演します:

julia tools/check_factor_tables.jl <factors_dir>
python tools/temari_factors_contract.py <factors_dir> --negative

上の範囲の表を思い出してください: 再生成した factors テーブルはアーカイブとバイト同一であるとは保証されません。検査器は値をその許容に対して比較し、アーカイブの SHA-256 がどのバイトが出荷されたかを教えてくれます。

そして、変更が値を変えることを意図したものであるときは、無効化した変種をビルドして、それが旧コードとビット同一であることを確認してください。検証を参照してください。