導入: なぜ「15桁書いたのに合わない」のか
数値計算の現場で、レガシーなFortranコード（F77など）をメンテナンスしている際、「ソースコード上では15桁の精度で数値を記述しているのに、計算結果の精度が7桁程度までしか出ていない」という現象に遭遇したことはありませんか。これは、数値計算エンジニアが必ず一度は通る「定数代入の罠」です。本記事では、コンパイラが定数をどのように評価し、なぜそこで精度が失われるのか、そのメカニズムと正しい記述方法を解説します。
基礎知識: 単精度と倍精度の内部表現
コンピュータ上で実数を扱う際、主に「単精度（Single Precision）」と「倍精度（Double Precision）」の2種類が使われます。
・単精度: 32ビットで表現。有効桁数は約7桁。
・倍精度: 64ビットで表現。有効桁数は約15～17桁。
問題は、ソースコードに記述された数値定数が、コンパイラによって「どの型として解釈されるか」という点です。多くのコンパイラは、指数表記（E）や末尾の指定がない小数を「単精度」として評価します。そのため、倍精度の変数に単精度の定数を代入すると、代入前に一度「精度が7桁に丸められ」、その後に倍精度へ拡張されるため、失われた精度は二度と戻りません。
実装/解決策: 指数指定子の正しい使用
この問題を解決する最も確実な方法は、定数の末尾に倍精度であることを明示する「指数指定子」を付与することです。Fortranにおいて、倍精度を明示するには `D` を使用します。これにより、コンパイラはその数値をコンパイル時から倍精度として正しく処理します。
サンプルプログラム
以下のコードは、精度の欠落が発生するケースと、正しく記述するケースを比較した例です。
! プログラム名: precision_test.f90
PROGRAM precision_test
IMPLICIT NONE
REAL(8) :: x_bad, x_good
! 誤った代入: 単精度として解釈され、丸め誤差が発生する
! 3.141592653589793 という値が内部で単精度に変換されてから代入される
x_bad = 3.141592653589793
! 正しい代入: D0 を付与し、倍精度定数として定義する
x_good = 3.141592653589793D0
PRINT , "単精度として代入 (誤差あり):", x_bad
PRINT , "倍精度として代入 (精度維持):", x_good
END PROGRAM precision_test
応用・注意点: 現場で役立つチェックリスト
1. 定数の末尾を確認する: 数値計算を行う変数に定数を代入する際は、必ず末尾に `D0` を付ける習慣をつけましょう。特に円周率や物理定数など、コード内に直書きされている定数は要注意です。
2. コンパイルオプションの活用: 一部のコンパイラでは、オプション（例: gfortranの -fdefault-real-8 など）によってデフォルトの実数型を倍精度に変更できます。しかし、これはコード全体の挙動を変えるため、既存のライブラリや単精度が意図されている箇所との整合性に注意が必要です。
3. 警告メッセージを無視しない: 最新のコンパイラであれば、精度が低下するような代入に対して警告（Warning）を出力するものもあります。ビルドログに警告が出ていないか、定期的に確認することをお勧めします。
これらの基本ルールを守るだけで、数値計算コードの信頼性は劇的に向上します。特にシミュレーションの初期値や収束判定値など、小さな誤差が積み重なって大きな問題になる箇所では、この「D0」の有無が運命を分けます。
コメント