富阪幸治（新潟大）

１．はじめに

「宇宙における流体力学」の特徴といったらなんだろうか？それはなんといっても自己重力の存在であろう。我々の周りの流体たとえば気象や海洋を考えるときにも重力は確かに重要である、しかし、これらは流体の流れとは無関係に存在する地球の重力であり、外場である。それに対して宇宙スケールでは、流れとともにそれ自身の万有引力によって重力場も変化する（この性質を強調して自己重力という言葉を使う）。これが宇宙の流体力学を考える上でエッセンスだと筆者は考えている。

この自己重力の性質によって、広がった気体から星や銀河の形成が起こるのだが、ここでは広がった気体から星への収縮を取り上げて、数値流体力学として自己重力の働く現象をどのように取り扱うかを述べてみたい。

２．自己重力の性質

自己重力の特徴は簡単な系でもわかるので、説明のために１次元球対称系（独立変数は半径rと時間tのみ）を考えよう。流体に働く力には圧力によるものとともに重力があるので、流体の運動方程式は、密度ρ、速度v、圧力p、重力定数Gなどを使って以下のように書ける。

ρdv/dt=-∂p/∂r-GρM(r)/r２

ここで、rの点での重力に効くのは半径rより内側に含まれる質量M(r)で、収縮が球対称のまま起これば不変な量である事に注意していただきたい。式の左辺と右辺第２項の比較から、重力が効いて収縮している時の収縮の時間スケール［自由落下時間］が見積もれる。（M〜ρr３を仮定すると）t〜r/v〜1/(Gρ)１／２ である。また自己重力が効いている構造のサイズは、右辺の第１、２項の比較から等温音速cｓ=(p/ρ)１／２を使って、r〜cｓ/(Gρ)１／２［特性的長さ］と書ける。ここでは、球対称１次元系についての運動方程式から考えたが、次元解析を行えば、代表的な量としてこの自由落下時間と特性的長さが出てくるのはあきらかであろう。

さて、収縮が進む（密度が上がる）につれて自由落下時間が短くなるということは、時間がたつにつれて収縮が急速に進むという事を示しており、特性的長さが密度の上昇に応じて小さくなるということは、収縮が進めば非常に小さい構造が出現することを示している。また、収縮の時間スケールは密度の1/2乗に反比例するのだから、場所によって密度が異なっているような場合にはもっとも密度の高い部分の収縮がまず起こり、密度の低い部分はその収縮に取り残されることがわかる。

星は、密度が約100cm−３で大きさが〜10１９cm程度の星間雲から生まれる。その星間雲の内部の〜3x10１７cm程度の大きさのさらに密度の高いコアと呼ばれるいくつかの領域で星形成が起こることが観測からわかっている。しかし自由落下時間t〜1/(Gρ)１／２で収縮が進むといっても、平均密度100cm−３の星間雲がその自由落下時間t〜10７年ですべて星になってしまったなら、星が過剰に形成されすぎることになる。現在標準的と考えられているモデルでは、ある程度の星間雲やそれが含むコアは収縮しているのではなく釣り合いの状態にあると考えている。その場合、力学的な釣り合いは、引力としての重力と反発する力である熱圧力、磁気力、遠心力などの間のバランスで保たれている。しかし星間雲の内部では電離度が低く、ゆっくりとではあるが、中心部から磁場が逃げ出して行く。そのため磁場のローレンツ力の支えを失って星間雲ははゆっくりと収縮し新しい釣合を保とうとする（準静的な進化）。しかし、これには限度があるので、もはや力の釣り合いを保てなくなると、星間雲は、自由落下時間程度で急速に収縮すると考えられている。この動的な収縮過程について多重格子法を用いた磁気流体力学シミュレーションで調べた結果を次の節でみてみよう。

３．多重格子法と星間雲の収縮

さて、前の節で述べた自己重力の性質から、星間雲の収縮を数値流体力学の問題として解くために必要なのは、急速に収縮する密度の高い中心部の細かな構造とゆっくりと進化する周辺部の密度の低い大きな構造を同時に計算できることであるといえる。もちろん中心部で細かくなるように格子を張った非構造格子を用いることも可能であるが、中心部の収縮に応じて、空間分解能を増してゆく、動的な格子形成が必要になる。ここでは、構造系の格子でこの要請を実現する方法として多重格子法（multigridあるいはnested grid scheme）を紹介する。

図１にここで用いた格子系の概念図を示した。一様間隔の格子で空間分解能の異なる複数の格子系を準備し、中心部分は細かな格子が、周辺部分は荒い格子が覆うようにする。ここで、おのおのの格子系（Ln：n番目のレベル）の格子間隔は、一つあらい格子系(Ln-1：n-1番目のレベル)の1/2となるようにとる。このようにとると、Ln-1の流体素片はLnの（２次元の場合）４個と重なることになる。一つの格子系はそれぞれ同じ数の格子点で分割する（ここでは64x64）。こうすることによって、すべてのレベルの計算は同じデータ構造とプログラムを使うことができる。すなわち、同一のプログラムを異なるレベルの物理量に対して適用することができる。15段（L0からL14まで）の格子系を使えば１次元あたり2６＋１４≒10６点分の空間分解能を中心部で持つことができる。CFL条件からLnの時間ステップΔtｎはLn-1のそれの1/2程度と予想されるので、Ln-1を１回時間推進するごとにLnは２回計算しなければならないことになる。従って、計算量は、L0のそれの1+2+4+...+2１４≒2１５倍が必要になることになるが、もし多重格子法を用いないで2２０x2２０の格子点を使ったならば、計算量は［データが（２次元の場合）2２８倍でCFL条件からこれを2１４回時間推進する必要がある。計算量がデータ量に線形にしか比例しないとしても］2３ｘ１４倍にもなるのだから、多重格子法がこのような自己重力収縮系を数値流体力学として取り扱う場合の非常に優れた方法だということがわかる。

粗い格子系で細かい格子にも覆われている部分は計算してもしなくても良いが、対応する細かい格子での値が求まった後、これを平均した量で置き換える。（後で使われない部分を計算しても、凹型の領域を掃くよりループ構造が単純になるのでベクトル型スーパーコンピュータで実行するならたいした損にはならない。）全体を覆っていない細かい格子系の境界条件は、その外側の一つ荒い格子系の量から決める。この二つのプロセスで分解能の異なる格子系は情報をやりとりし、全体としてつじつまのあった解を与える。

さて、実際の軸対称２次元の多重格子法を用いた磁気流体力学シミュレーションで星間雲の重力収縮を調べた結果についてみてみよう。力学的釣合状態の磁場に貫かれた等温星間雲の構造（図2aはそのときの構造を示す）を前もって解き（これは核融合プラズマの平衡配位を求める問題と自己重力の効果をのぞいて同じなのだが）、これに中心部の密度を10％あげるような小さな擾乱を与えて釣り合いを破って重力収縮を開始させた。これはぎりぎり支えられた状態から磁場が抜けて動的な重力収縮が始まるという前節のシナリオに沿って、その後の進化を調べたことにあたっている。

次の図は赤道（z=0）面での密度（実線）、磁場の強さ（破線）、動径方向の収縮速度（点線）の分布を示している。中心からの距離で10万倍、密度比では10１０倍のダイナミックレンジをを持っていることがわかる。それぞれの線は異なる時刻に対応する。これをみると、周辺部分（r>0.1）の密度分布はほとんど変化せず、収縮速度も遅い。特徴的なことは、高密度のコア（密度の一定の部分をここではコアと呼ぶ）の大きさは時間とともに、あるいはコアの密度上昇とともに減少する。また、中心密度がほぼ10倍になるたびに線を引いているが、その間の時間間隔が次第に短くなっていることもわかるだろう。これらはいずれも前節で述べた自己重力の性質である。もう一つの特徴は、密度も磁場の強さも、ある指数分布ρ=Arーα、BＺ=r−βに漸近しているということである（α≒β/2）。これらは、等温の自己重力収縮系で普遍的にみられる性質で、このまま、状態方程式が等温で収縮が続くならば有限の時間でr=0の点の密度が（数学的には）発散する。コアの質量は減少しながら［Mｃｏｒｅ〜ρｃｏｒｅ｛cｓ/(Gρｃｏｒｅ)１／２}３∝ρｃｏｒｅ−１／２だから］密度が際限なく上昇する過程で、これは暴走的な収縮（runaway collapse）と呼ばれている。ちなみにこの収縮は回転している等温星間雲においても起こることが知られている。角運動量の保存からr−３に比例する遠心力バリアが働いても起こるところが一見パラドックスめいているが、遠心力バリアで外部の物質は収縮とともに支えられてゆき、（数学的には無限小の）中心部の質量が暴走的に収縮するだけなので、矛盾ではない。この中心密度が（数学的には）発散するまでを「前期の解」と呼んでいる。

図2cは中心密度が表面の10１０倍に達した状態のもっとも細かいL14でみた構造を示している。縮尺が全く異なるので直接比べられないが、系全体はほぼ球状で図2aとほとんど変化していないのだが、磁場に垂直な方向の流れは妨げられるという磁場の効果によって、密度の高い部分はディスクを形成していることがわかる。磁場に平行なz方向（左右）からの流れは磁場に沿っており、ディスク内部のみが半径方向に収縮していることがわかる。これは、ディスク内部は熱圧力が磁気圧に対して支配的であるのに対して、ディスクの外部では磁気圧が熱圧に勝っており磁場が流れをコントロールしているためである。実際の星間気体では、ここまで中心が収縮するすこし前に、それまで星間雲を等温に保ってきた（星間雲は、星間塵が気体から熱をもらって、それを熱輻射として外部に捨てることによって冷やされている）星間塵からの熱放射が星間雲内を自由に抜けられなくなり、次第に温度が上がってくると考えられている。

この後は、将来星に成長する中心の高密度部分に外からの流れが重力によって降り積もる降着期に入る（この解を「後期の解」と呼ぶ）。この過程でディスクはさらに薄くなり、ディスクの収縮速度は音速の数倍程度まで加速される。この収縮しているガスが、そろそろ電波干渉計観測で見つかり始めており、数値流体力学と観測の直接の比較も可能になりつつあるのが現在の状況である。

４．まとめ

現在天文学では、星間気体から星や、銀河間気体から銀河といった自己重力が重要な役割を果たしている構造形成過程が注目を集めている。これらの過程を数値流体力学の手法で研究するために、多重格子法が非常に適した方法であることを、その応用として星間雲の重力収縮の例を挙げて述べてきた。この際重要なことは、スキームを多重格子にすること（多重格子化）は、どのような計算法をとっているかに関わりなくできるといういう点である。ここで紹介した方法は、細かい格子系を中心部に次々に生成するということを外から与えていたいわばアドホックな計算法であったが、格子発生の条件を自動判別し、任意の形状の格子系を自動的に発生・消滅させる適応型（adaptive）スキームへ発展させることがいま強く望まれている。これによって、重力によるディスクの分裂など、より複雑な構造の解析が可能になるだろうと考えている。

(1)多重格子法の流体力学への応用については
Berger M.J., Oliger J. 1984, JCP 53, 484
Berger M.J., Colella P. 1989, JCP 82, 64
(2)ここで紹介した多重格子法による星間雲の収縮については
Tomisaka K. 1996, Publ. Astron. Soc. Japan 48, No.5

図の説明

図１：多重格子法の概念図。L0、L1、L2の三つのレベルからなる多重格子を示す。図はz軸が横軸、動径座標が縦軸であることに注意。実際には15レベルを用いて計算した。

図２：
(a)釣合の状態にある初期状態。L2での密度の等高線（ほぼ円形）と磁力線（z方向に走っている）を示す。横軸のz軸が対称軸。

(b)赤道面上（z=0）の密度（実線）、磁束密度（破線）、収縮速度（点線）の分布。数字は収縮開始からの時間tを表しており、磁束密度、収縮速度ともより上側の線が時間が経過した状態に対応している。密度と磁束密度は対数目盛、速度は内向きの値をそのまま示している。横軸が距離の対数目盛である点に注意。

(c)中心密度が表面の10１０倍に達した最終状態をL14格子でみたもの。図(a)と2１２=4000倍だけ縮尺が異なることに注意。密度等高線（主に縦方向に走っている）の込んでいる部分に等温斜め衝撃波が生じており、z方向からきた流れはここで減速、圧縮される。衝撃波面に挟まれて高密度のディスクが生じている。磁力線（主に横方向に走っている）がディスクの動径方向への収縮に引きずられていることがわかる。