next up previous
Next: 2 少ない粒子数での計算 Up: 重力多体系 Previous: 重力多体系

1 はじめに

一応、ここでの課題は、12ページで重力多体系の数値計算の全てを議論しろ、 というものだが、こんなのは本3冊くらいの分量を書かないと無理と思われ る。多少詳しい資料は例えば牧野の講義ノート1をみていただくなり、何故か「物性研究」にのった牧野の文章 [Mak01] や Binney and Tremaine[BT87] (まあ、数値計算のことはろくに 書いてないが)を見ていただければいいわけなので、ここでは適当なこと を思いつくままに書く。

一般に$N$体計算に限らず、数値的な研究というのは難しい。紙と鉛筆で仕事 をしている理論家な人の中には、「数値計算なんて、基礎方程式を計算機にプ ログラムすればいいんでしょ」くらいに思っている人が珍しくない。そういう 人に指導された学生(あるいはその本人)がなにもわからないで数値計算のよう なことをして無意味な結果を出してしまうこともまたありがちなことである。

数値計算による研究の水準を向上させるためには、そういった「底辺」な研究の レベルアップが、トップレベルの研究をより高度なものにしていくのと少なく とも同程度 には重要であるというのが理屈ではある。まあ、しかし、問題は、 大抵のこういった「底辺」な人は、数値計算というものが難しいということを そもそもまるで認識していないために底辺にとどまっているわけであり、自分 が認識していないという認識ももちろんないので例えばこういったところにな にか書いたとしてもそれが「底辺」な人の目に止まる可能性は事実上ないであ ろう。つまり、底辺を上げようという努力自体には意味がないとも限らないが、 ここで何か書くことは底辺の向上にはつながらない。

つまり、ここになにか書くことには、集録のページふさぎとか予算消化とか書 く人の自己満足とかいった程度の意味しかないと考えられる。

そういうわけで、こんな原稿を書いてる暇があればもっと生産的なことをした ほうがいいし、読んでいるあなたも著者がそんなつもりで書いているものは読 むのはやめにして次の原稿に進むほうがいいような気がしてきたであろう。ね、 してきたでしょ?

cm

というわけで、読む人がいなくなったところで本題に入る。

本題は、自己重力多体系の数値計算を、どうすれば「正しく」できるかである。

自己重力多体系の基礎方程式はもちろん、各粒子の運動方程式


\begin{displaymath}
{d^2 \mbox{\boldmath$x$}_i \over dt^2} = \sum_{j\ne i} Gm_j...
...er \vert\mbox{\boldmath$x$}_j - \mbox{\boldmath$x$}_i\vert^3},
\end{displaymath} (1)

である。ここで $\mbox{\boldmath$x$}_i$ は粒子 $i$ の位置、 $m_j$ は粒子 $j$ の質量、 $G$ は重力定数である。ガスとか磁場とか恒星が進化するとかなんとかかんと かを考えなければ、多体系のシミュレーションとは単にこの方程式をプログラ ムして数値積分するだけである。

流体シミュレーションや一般相対論的な自己重力系の計算に比べると、古典 ニュートン重力に従った自己重力多体系のシミュレーションには著しい特徴が ある。それは、シミュレーションプログラムが「正しい」かどうかは極めて精 密にチェックできるということである。

さて、ここで、シミュレーションが「正しい」というのはどういう意味かとい うことを明確にしておく必要がある。理論天文学においては、これには少なく とも 3 つのレベルがある。

第一は、物理過程として何を取り込み、何を落としているかということである。

つまり、例えば銀河系を自己重力多体系として扱うなら、星間ガスやそこでの 星形成の影響はもちろん入ってこない。解明したい問題に対して、そうい う扱いが適当かどうかというのは、当たり前なことではあるが検討されている 必要がある。というか、もちろん検討して駄目だとわかるならばそれは駄目な わけで、何かしら別のことを考えないといけないことになるのはいうまでもな い。

次に問題になるのは、数値的に扱うことにした系と、実際に計算できる系の違 いである。 例えば、銀河系を自己重力多体系として扱うなら、粒子数は $10^{10}$ 以上である。これだけの粒子を扱うのはまだしばらくは無理なので、 シミュレーションではもっと粒子数の少ない系を扱う。また、計算を容易にす るために重力ポテンシャルをいじったりもする。これらの影響はどういうもの かというのが第二の問題である。

粒子数や重力ポテンシャルを決めてしまえば、解くべき問題はハミルトン力学 系の初期値問題であり、相互作用、つまりハミルトニアンは完全に決まってい る。数値解法が問題になるのはここからである。つまり、ここまで来ると話は 単純な常微分方程式の初期値問題であるわけだが、これを正しく解くというこ とが第三の問題である。

このように、近似の段階をかなり明確に分離して議論できるために、数値計算 の信頼性については比較的きちんとした評価が可能である。特に重要なのは、 最後の、ハミルトン系になった後、それを常微分方程式系として数値積分する 時には、数値計算があっているか間違っているかは非常に確実に評価できるこ とである。

具体的には、特殊な計算法を使わないかぎり、数値誤差は普通はハミルトニア ン、つまり全エネルギーや、全角運動量のような保存量の誤差となって現れて くる。従って、保存量の精度に問題があれば、計算が破綻していることは確実 である。

もっとも、逆は真ではないことはいうまでもない。つまり、保存量が保存した からといって、数値解があっているとは限らない。しかし、特別に保存量を保 存するようにデザインされた数値解法でなければ、保存量が保存していてしか も計算が破綻しているということは極めてありそうにない。

なお、保存量の精度が関係するのは、あくまでもここでいう第三の問題だけで あって、第一と第二の問題とはなんの関係もないということは、当たり前のこ とではあるが必ずしも理解されてはいない場合もあるので一応注意しておこう。

もちろん、この3個の問題は完全に独立に議論できるわけではない。以下の議 論では最初の問題、つまり我々が扱いたい問題は本当に重力多体系なのかとい う問題には目をつぶって、後の2つだけを考える。この時、第二の問題がどの 程度深刻であるかは、端的には粒子数で決まる。粒子数無限大の(と思ってい い)系を有限粒子数で近似するなら、粒子数が大きいほど問題が少ないことは 間違いないからである。

ところが、計算機の能力は有限なので、粒子数を増やせばそれだけ計算時間が かかる。しかし、計算時間は計算法と要求する計算精度にも依存するので、初 期問題として解く時の要求精度を落とすことでより大きな粒子数が扱えるかも しれない。つまり、実際に「もっとも良い」計算結果を得るためには、上の第 二の問題と第三の問題の両方を理解し、制御できる必要がある。

このことを言い換えると、第三の、ハミルトニアンが与えられた時にそれをい かにして正確に解くかという問題は、単純に正確に解けばいいというのではな く、必要な精度をいかに少ない計算時間(場合によってはメモリ使用量)で実現 するかという問題、つまりは計算をどうやって速くするかという問題なのであ る。

というわけで、まず第二の問題はどういうものかということをもう少し詳しく みた後、第三の問題を、計算速度という観点から見て行こう。



Jun Makino
平成15年3月19日