精度記録型量子化学計算

4P110（北教大・釧路）小原繁

Accuracy-Recording Quantum-Chemical Calculations(Hokkaido University of Education, Kushiro Campus)

1

Shigeru Obara【序】

まうことが分った。別の言い方をすると、より正確な精コンピュータの発達に伴い量子化学計算の研究対象 度算出には「精度が増す要因」を探す必要があることは大型化して来ている。その一方で計算そのものは依 が分った。
然として倍精度実数を使って行われており、精度不足「確実に精度が増す要因」を四則演算だけから探しの発生が懸念される。精度不足発生の有無そのものは、 出すことはできない。考え得るものは「確率的に精度別途精度の高い計算を行って比較することにより確認 が増す要因」だけである。この要因として、本研究ではすることはできる。しかし、このような別途計算を常 以下の三つの演算を取り上げ「精度追加」を行った。に行うことは不便である。それに、きめ細かな解析を
1. 減算精度落ち時追加 減算精度落ちにより n ビッ
行うことも容易でない。通常の量子化学計算を行いなト残ったとする。この場合、精度は n ビットへ減っがら計算精度も常に記録していくことができるとこれたと考えるのではなく、さらに下位の 1 ビット（0らの不便が解消する。この観点から、どのように精度か 1 かのいずれか）が確率的にある程度正しいのを記録していけば良いのか等を検討して本研究室で開で、(n + α1 ) ビットの精度へ減ったと考える。この発した高精度量子化学計算プログラムに精度記録機能α1 が「精度追加」である。
を追加した。これについて報告する。
2. 加算繰り上り時追加 例として有効ビット数が 32 【有効ビット数の増減要因】
ビットの下記 2 進数の加算を考える。
理系の数値計算では有効数字という概念が用いられ、101その桁数不足発生の有無が常に注視される。この場合の有効数字は 10 進数のことであるが、コンピュータで(1)+ 100
は 2 進数が使用されるので「2 進数の有効数字」という1001意味で「有効ビット」と表記し、また、「その桁数」という意味で「有効ビット数」と今後表記する。
和の値は繰り上りにより 4 ビットになっている。こ有効ビット数が減少する主要因は近接数同士の減算の数の有効ビット数は元々の 3 ビットよりも多くが行われることである。より正確には、絶対値が近接なったと考えてよいだろう。つまり、加算結果が繰する同符号二数の減算、あるいは、異符号二数の加算り上る時は有効ビット数が α2 ビットだけ増すと考が行われることである。
える。

−

3. 丸め繰り上り時追加 例として計算結果が下記の値だったとして、これを小数第 1 位で丸める場合を考える。
111.1 ⇒ 1000
(2)

これらの演算が行われると、有効ビットの上位が削ら丸めた値は 4 ビットになる。この数の有効ビット数れて下位だけが残り有効ビット数が減少する（上図）。
が元々の 3 ビットよりも多くなると考える。つまこの減少要因以外に有効ビット数が減少する要因を、り、丸め時に繰り上る時は有効ビット数が α3 ビッ通常、省られることはない。一方、有効ビット数が増トだけ増すと考える。なお、この繰り上りの起き加する要因については、通常、全く考慮されることがる頻度はかなり低い。しかし、一応、有効ビットない。
数が増す要因の一つとして掲げておく。
計算精度を決める以上の要因、つまり、「減算だけによる精度減少」により (18 × 18) 対称行列の対角化計算 3 【計算方法】（16 倍精度計算）を行い、計算結果の各要素の有効ビッ
上記の「精度追加」α1 , α2 , α3 によって有効ビット数ト数を算出した。また、この算出結果の精度を確認するを適切に算出することができるかどうかを確認するたために、高精度計算（32 倍精度計算）を別途行って得めに、Rn 原子の量子化学計算を行ない、重なり積分のた結果から精度を算出した。両者を比較すると、極めHouseholder 法による三重対角化における有効ビットて残念なことに、差が極端に大きかった。実際には、20数の算出を行なった。
ビット強しか精度落ちしていないのに、40 ビット強の量子化学計算には本研究室において Java 言語を使用
精度落ちがあると算出され（表 1 の 2 行目と最終行）、して開発した量子化学計算プログラムを使用し、16 倍「減算だけによる精度減少」では精度を低く算出してし精度（仮数部 496 ビット）の計算を解析対象データと

表 1: Rn 原子 16 倍精度計算の重なり積分を House-holder 法で三重対角化したときの最低精度要素の有効ビット数の算出値と実際の値
表 2:

α1 と α2 の種々の値の組み合わせによる Rnα値†

主対角要素

副対角要素

原 子 16 倍 精 度 計 算 の 重 な り 積 分 三 重 対 角 化 後 の減 0.00, 加 0.00, 丸 0.00

450459
467
472
476

451
461468473477

最低精度要素の有効ビット数の算出値と実際の値。

450456461464466451457462464466

いは、減算桁落ち時追加精度を大きくして加算繰り上

450450450450450451
451
451
451
451

α値†

主対角要素

副対角要素

減 0.00, 加 0.00

472.54‡473.58‡

450
456
461
464
466

451
457
462
464
466

459
466
469
472
473

461
467
470
472
473

467
471
475
476
478

468
472
475
477
478

472
476
478
479
480

473
476
478
480
481

476
478
480
481
482

477
479
481
482
483

472.54‡473.58‡

減 0.25, 加 0.00, 丸 0.00
減 0.50, 加 0.00, 丸 0.00
減 0.75, 加 0.00, 丸 0.00
減 1.00, 加 0.00, 丸 0.00
減 0.00, 加 0.00, 丸 0.00
減 0.00, 加 0.25, 丸 0.00
減 0.00, 加 0.50, 丸 0.00
減 0.00, 加 0.75, 丸 0.00
減 0.00, 加 1.00, 丸 0.00
減 0.00, 加 0.00, 丸 0.00
減 0.00, 加 0.00, 丸 0.25
減 0.00, 加 0.00, 丸 0.50
減 0.00, 加 0.00, 丸 0.75
減 0.00, 加 0.00, 丸 1.00†

この列の 3 つの数は α1 , α2 , α3 の値を表わす‡

32 倍精度計算から算出した実際の有効ビット数

実際の値に近い [α1 , α2 ] の値は [0.25,0.75](472,472)、[0.25,1.00](473,473)、および、[0.75,0.00](472,473) だった。これは、減算桁落ち時追加精度を小さくして加算繰り上がり時追加精度を大きくする（前者２つ）、あるがり時追加精度をなしにする、と実際に近い値を算出できると解釈できるが、より確実なものにするには一層の比較・検討が必要だろう。

減 0.00, 加 0.25
減 0.00, 加 0.50
減 0.00, 加 0.75
減 0.00, 加 1.00
減 0.25, 加 0.00
減 0.25, 加 0.25
減 0.25, 加 0.50
減 0.25, 加 0.75

した。また、精確な値を得るために 32 倍精度（仮数部992 ビット）の計算を行ない、これを真値として用いた。
なお、このとき数値計算上の誤差だけを解析に使用するため、解析対象データである 16 倍精度の重なり積分の各要素の仮数部下位に 496 ビットの 0 を追加することにより得た行列を 32 倍精度計算で用いる 32 倍精度重なり積分とした。

4【精度追加と有効ビット数】

それぞれの「精度追加」α1 , α2 , α3 の種々の値と、それを使用して算出された重なり積分三重対角化行列の最低精度の主対角要素と副対角要素の有効ビット数を表 1 に掲げる。
この表の最左列は α1 , α2 , α3 の値を表わす。値が 0.00
は追加する有効ビット数が零ビットであることを表わし、1.00 は追加する有効ビット数が 1 ビットであることを表わす。表の中央列と最右列はそれぞれ主対角要素と副対角要素の中で最低精度の有効ビット数の算出値である。また、最終行には 32 倍精度計算から得られた当該要素の実際の有効ビット数を掲げる。
表から α1 や α2 の値によって算出された有効ビット数が大きく変わることがわかる。一方、α3 では当初の予想通り算出有効ビット数に影響が現れない。
α1 と α2 の種々の値の組み合わせによって三重対角化行列の最低精度要素の有効ビット数の算出値がどう変わるをまとめたのが表 2 である。詳細は当日発表する。

減 0.25, 加 1.00
減 0.50, 加 0.00
減 0.50, 加 0.25
減 0.50, 加 0.50
減 0.50, 加 0.75
減 0.50, 加 1.00
減 0.75, 加 0.00
減 0.75, 加 0.25
減 0.75, 加 0.50
減 0.75, 加 0.75
減 0.75, 加 1.00
減 1.00, 加 0.00
減 1.00, 加 0.25
減 1.00, 加 0.50
減 1.00, 加 0.75
減 1.00, 加 1.00†

この列の 2 つの数は α1 と α2 の値を表わす‡

32 倍精度計算から算出した実際の有効ビット数