Previous | ToC | Next |
最近ぼちぼち色々な計算機で FP16 が使えるようになってきたので、 それで多体問題、少なくとも重力多体問題ができないかという雑談を 何人かでしていた時にでてきた方法について。
仮数長は 10ビットもあるので、GRAPE-1, 1A, 3, 5, 8 あたりに比べて長く、 相互作用計算にはいってからの精度には特に問題はありませんが、問題は
とはいえ、自分は原点近くにいて、遠くのほうは相対精度が落ちてよいとすれば、遠くは「桁を合わせない」という方法でなんかできる気がします。
積算は複数語のアキュムレータを使う方法が古くから知られています。
問題は途中でのオーバーフローです。距離のダイナミックレンジ自体は30ビッ ト程度はほしいので、途中で必要な距離の3乗の逆数を計算すると90ビットに なり、 FP16 で持つ5ビットの指数では表現できなくなります。
ここで、ツリー法を使い、質量の空間分布が2ないし3次元でほぼ一様な場合だ けをまず考えると、相互作用を計算する超粒子の質量は距離の 2 ないし3乗に 比例します。まず、色々上手くいきそうな 2次元の場合を考えると、 距離と質量の平方根が比例するので、距離を質量の平方根でわったものは ほぼ一定になります。従って、普通なら(dx 等は3次元ベクトルで内積があるとして)
dx = xj-xi r2= dx*dx r=sqrt(r2) mr3inv=mj/(r*r*r) fi += mr3inv*dxと書くところを
q=1/sqrt(mj) dxhq = q*(xjh-xih) dxlq = q*(xjh-xil) dxq= dxhq+dxlq r2q2= dxq*dxq rq=sqrt(r2q2) mr3invqinv=1/(rq*rq*rq) fi += mr3invqinv*dxqとしてやると、 dxhq や dxlq はダイナミックレンジが小さく、基本的に1の オーダーになっているので、指数のオーバーフローやアンダーフローが ほぼ起きなくなることが期待できます。最初の、質量の逆数の平方根をとると ころは、ツリー構築を作った時に1度だけすればいいので、計算量の増加は無視できます。
3次元の場合には、一様分布の場合遠くの超粒子からの力が近くの粒子からの 力より大きくなりますが、これは距離に比例する程度なのでやはり上のような 修正を行うことでダイナミックレンジを大幅に小さくできます。
3次元で TreePM とか使うなら距離にはカットオフがあるのでこれで問題ない でしょう。開放境界の3次元問題の時には、結局、トータルの力が決まるあた りの特徴的な距離スケールがあるので、そこでオーバーフローが起きないよう にスケールを決めることができるはずです。
Previous | ToC | Next |