GRAPE プロジェクトが始まったのは、必然というよりは様々な偶然の結果であっ
た。 1985 年に大学院に進学した私は、天体物理学の世界的権威である
杉本大一郎の下で、球状星団の力学的進化の直接N体計算による研究と
いうテーマをやっていくことになった。球状星団というのは、古い星が典型的
には 100万個位集まって出来ているシステムで、我々の銀河系のなかには150
個ほどある。大きな銀河では数千個の球状星団を持つものもある。
球状星団のような恒星系がどのように進化していくかというのかというのは天体物理学的には大変重要な問題
であるが、それにもかかわらず当時は（実は現在でも）かなり大雑把なことし
かわかっていなかった。
一つの恒星については、フォン・ノイマンの意向もあって1940年代末に電子計 算機が利用可能になるとすぐ、その構造、進化の計算がマーティン・シュワル ツシルドらによって精力的に進められた。その結果、1960年代始めにはほぼ基本的 な理解が出来上がっていた。それに対し、重力でまとまっているという面では 一つの恒星と同じようにも見える恒星系については、1960年代始めに数値計算に よる研究が始まったものの、その理解は遅々として進んでいかなかった。
理解がなかなか進まなかった理由の大きなものは、「計算できなかった」とい うことである。星はガスが重力でかたまって出来ている。球対称を仮定して回 転等を考えなければ、解くべき方程式は1次元の流体力学の式であり、力学的 なつりあいや熱力学的な進化を数値計算することはそれほど難しく ない。このために、初期の電子計算機が理解の飛躍的な 進歩をもたらした。
ところが、恒星系では星は星団の中全体を他の星とぶつからないで動き回る。 このために、温度とか圧力といった熱力学的な概念がなりたたず、 星団の構造を記述するためには星の分布関数（位置と速度の6次元空間で）を 知る必要がある。またその進化を表現する方程式も6次元空間のなかでの分布 関数の偏微分方程式となる。分布関数をモンテカルロ積分するとか、ある いは角運動量を無視して次元を減らすとかの方法で近似的な計算が可能になっ たのが 1970年代から 1980年代にかけてであった。
近似的な計算でわかったことは、球状星団のような恒星系は「重力熱力学的
カタストロフ」という特異な現象を起こすはずであるということだった。恒星
系では局所的な温度という概念はないが、系全体としては熱力学第二法則にし
たがってエントロピーを生成し、熱平衡にむかって進化しようとする。これは、
平均的には星団の中心付近の速度の大きい星から、外側の速度の小さい星にむ
けてエネルギーが流れるということである。ところが、重力が働いているため
に、中心付近の速度の大きい星はエネルギーを失った結果軌道が縮まって、そ
の分重力が大きくなりますます速度が上がる。これに対し、外
側は逆に軌道が膨らんで速度
が下がる。このため、エントロピーを生成したのにいっそう熱平衡から遠ざか
るという不安定現象が起きるのである。この現象は「重力熱力学的カタストロ
フ」と名付けられた(例えば[7])。
この後の時間発展を近似的な計算で追跡すると、宇宙の 年齢よりも短い時間で恒星系の中心密度が無限大になってしまう。そのあとは どうなるかについては様々な理論が建てられた。ブラックホールが生成される、 あるいは核融合反応で水素からヘリウムが出来るように、星同士の３体重力散 乱で連星が出来てこれが星団を安定な定常状態にする、あるいはこの定常状態 がまた熱的に不安定で、「重力熱力学的振動」と名付けられた非線形振動を起 こすなどである。
球状星団の星の数は100万個程度である。従って、100万体問題を数値計算して しまえば、なにが起きるか原理的にはすべてわかる。私が大学院に入った時に 始めた研究は要するにそういうことをやろうとしたのである。しかし、1980 年代終りに至っても、当時のスーパーコンピュータを数百時間使ってさえ、扱 える星の数はわずか3000程度であり、これでは100万体問題でなにが起きるか についてはっきりしたことはいえなかった。
このように扱える星の数が少ないのは、基本的には時間積分のタイムステップ 数が多いからである。重力熱力学的不安定の成長するタイムスケールは、ほぼ 粒子数に比例して長くなるということがわかっている。このために、粒子数の 3乗に比例して計算量が増える。
並列計算や高速アルゴリズムに詳しい読者の方はなぜ Barnes-Hut ツリー [4]やFMM(高速多重極展開法) [5,3]を使わないのだと思われるであろう。 その理由は、積分時間が長いために高精度が必要で、そのためにこれらの近似 算法による利得があまり大きくないことが一つである。もう一つは、上に述べ たように中心の密度がどんどん上がっていくという特性のために、星ごとにば らばらにタイムステップを変える必要があり、その場合にはFMMのような方法 をベクトル計算機や並列計算機で効率良く使うのは難しいということである。