BGK
モデルを用いた相対論的な運動論的方程式の
分散関係と相対論における散逸過程の研究
高本亮
京都大学大学院理学研究科
物理学・宇宙物理学専攻物理学第二分野
天体核研究室
2009
年
2
月
16
日
概要
近年のX線やガンマ線等の観測技術の発展に刺激され、パルサーを初めとしてガンマ 線バーストやAGNジェットなど高エネルギーの天体現象に対して興味が高まっている。 これらの現象の多くは相対論的な流体現象としての取り扱いが必要であり、実際世界中で そのような研究が精力的に行なわれている。しかしこれまでの研究では多くの場合散逸の ない理想流体として扱っており、散逸は適当な現象論的な形で導入するか、もしくは全く 考えないかのどちらかで、満足な形で取り扱われているとは言えない。散逸はガスの熱化 をもたらし観測される電磁波などのスペクトルと決定的に結びついているため、現状の理 論的解析では不十分である。 相対論的な流体に散逸をきちんと取り入れた研究が少ない理由は次の通りである。相 対論的な流体力学において散逸過程を取り扱うことを意図した最も簡単な定式化として Eckart の理論[1]とLandau による理論[2]がある。しかしこれらの理論に基づく時間発 展では因果律が破れ、非物理的な不安定モードが含まれる [3] 。従って相対論的な流体力 学において散逸の過程を取り扱うための基礎方程式としては使えない。これらの因果律の 破れと不安定モードの発生という問題を解決するために様々な研究が行なわれており、例 えば最も有名な理論としてはIsraelとStewartらによるGradのmoment展開を用いた 理論がある [4]。しかし相対論的な流体における実験検証が困難であることやボルツマン 方程式を近似なしで解くことが非常に困難であることから、これらの理論が実際どの程度 正確で適切なのかを確かめることは非常に困難である。 本研究ではこのような問題に対し、BGKモデルを用いたボルツマン方程式の線形摂動 を行い、得られた分散関係と比べることで散逸の理論がどの程度のスケールまで正しく再 現しているのかを判定する。BGKモデルは線形化された衝突項の全ての固有値(緩和時 間)を、その最も小さい固有値に置き換えるという近似であり、流体力学的なモードであ る減衰音波と熱伝導に関しては十分短波長までボルツマン方程式から得られる正しい結果 を導くと期待される。この結果を用いて実際の相対論的天体現象を解析する際に使える、 散逸項の実用的な取り扱い方法を吟味する。目次
第1章 イントロダクション 1 第2章 運動論 7 2.1 非相対論的運動論 . . . 7 2.1.1 分布関数 . . . 8 2.1.2 詳細釣り合いの原理 . . . 9 2.1.3 Boltzmann方程式の導出 . . . 10 2.1.4 H定理 . . . 12 2.1.5 平衡解 . . . 14 2.1.6 巨視的方程式 . . . 15 2.1.7 Chapman-Enskog展開 . . . 18 2.1.8 気体の熱伝導 . . . 20 2.1.9 気体の粘性. . . 22 2.1.10 moment展開 . . . 24 2.2 相対論的運動論 . . . 32 2.2.1 相対論的な分布関数と位置空間での物理量の定義 . . . 32 2.2.2 相対論的Boltzmann方程式の導出 . . . 33 2.2.3 H定理 . . . 36 2.2.4 相対論的な平衡解と熱平衡状態 . . . 38 2.2.5 Eckart分解とLandau-Lifshitz分解 . . . 43 2.2.6 巨視的方程式 . . . 47 2.2.7 Chapman-Enskog展開 . . . 51 2.2.8 moment展開 . . . 54 2.2.9 散逸の理論と安定性 . . . 59 第3章 運動論的方程式の線形摂動 61 3.1 非相対論的Boltzmann方程式 . . . 623.1.1 BGKモデル . . . 62 3.1.2 線形解析と分散関係の導出 . . . 63 3.1.3 Result . . . 66 分散関係 . . . 66 固有関数 . . . 70 3.1.4 Discussion . . . 72 流体力学的モードの長波長近似による導出 . . . 72 分散関係式から得られた解の考察 . . . 73 3.1.5 各項の大きさの比較 . . . 74 音波モード. . . 74 熱伝導モード . . . 75 その他のモード . . . 76 3.1.6 continuous spectrum . . . 77 3.2 相対論的Boltzmann方程式 − Marle流の取り扱い . . . 79 3.2.1 MarleのBGK model . . . 79 3.2.2 Marle流のBGKモデルの修正 . . . 80 3.2.3 ボルツマン方程式の線形摂動と分散関係 . . . 82 3.2.4 Result . . . 86 分散関係 . . . 86 3.2.5 Discussion . . . 93 3.3 相対論的Boltzmann方程式 − Anderson-Witting流の取り扱い . . . . 95 3.3.1 Anderson-WittingのBGK model . . . 95 3.3.2 ボルツマン方程式の線形摂動と分散関係 . . . 96 3.3.3 Result . . . 100 3.3.4 Discussion . . . 107 3.4 Israel-Stewart理論から得られた分散関係との比較 . . . 110 第4章 まとめと今後の展望 113 付録A 線形化された衝突項 115 付録B moment展開の双曲性と分散関係 123 付録C 位相速度の上限 125 付録D Fredholm型積分方程式 127
付録E 分散関係式の計算の詳細 129 E.1 非相対論的なBGKモデル . . . 129 E.2 修正したMarleのBGKモデル . . . 132 E.3 Anderson-WittingのBGKモデル . . . 135 付録F 解析接続 139 F.1 非相対論的なBGKモデルの解析接続 . . . 139 F.2 Anderson-WittingのBGKモデルの解析接続 . . . 141 F.3 修正したMarleのBGKモデルの解析接続 . . . 143
第
1
章
イントロダクション
宇宙における高エネルギー天体現象
20世紀初頭まで人類は可視光でのみ宇宙を眺めてきた。可視光領域は電磁波の波長の ごく限られた領域に過ぎなかったが、それでも人類はこの可視光を用いることで惑星の運 動や太陽の黒点活動、そして恒星や銀河などを観測してきた。この可視光のエネルギーは 1eV程度であり、約1万度の温度の現象を観測していることになる。そして20世紀中頃 から電波、赤外線、紫外線、そして特にX線とγ線での観測が始まったことにより、可視 光領域の観測では予想されていなかった様々な高エネルギーの天体現象が発見された。例 えばコンパクトオブジェクトとしては白色矮星、中性子星、ブラックホールなどそれまで 予想されていたとしても理論上の存在であると考えられていたような天体が実際に同定さ れた。X線バーストやAGNジェット、パルサー、高エネルギー宇宙線などはこれらの天 体の特徴である強重力場や強磁場などにより起こると考えられている。また中心エンジン は未だ不明であるが、γ線バーストとよばれる宇宙で最も激しく明るい現象もコンパクト オブジェクト由来であると考えられている。一方大きなスケールの現象では銀河団からの X線放射が観測されており、これは銀河団に大規模な高温プラズマが付随していると考え られている。 このように様々な興味深い観測結果に伴い高エネルギー天体の理論的な理解もすすんで きた。特に60年代から80年代はパルサーなど重要な天体の発見が相次いだこともあり、 観測結果を説明する様々な解析的モデルが立てられた。当時は観測結果がまだ少ないこと もあり、最も基本的な仮定のみを行なったモデルが数多く作られた。その中で観測から否 定されなかったモデルが現在の高エネルギー天体の理論的な基礎になっている。現在は観 測技術の向上に伴い精密な観測データが増えてきており、さらにコンピュータの数値演算 能力が飛躍的に発展したことから、理論的な研究は基礎的なモデルを立てる段階からより 正確な描像を突き詰める段階に移ったといえるであろう。さて理論的な研究に必要な道具自体についてはどうであろうか? まず低温の非相対論 的な現象はガスを流体とみなした場合、磁場や輻射を考慮した研究は、解析的なアプロー チも数値的なアプローチも盛んに行なわれていてかなり理解がすすんでいるといえるだろ う。ところが高温の相対論的な現象については話が変わってくる。まず一般相対論的な効 果を入れる場合は、重力場を含めて解析的に解くことはほとんど不可能なため数値相対論 を用いなければいけないが、新しい分野として発展が期待されているがいまだ十分進んで いるとはいえないであろう。一般相対論的な重力の効果を考えない場合でも輻射流体や磁 気流体としての扱いは難しく、まだ本格的な数値シミュレーションなどは始まったばかり である。このように流体に外場を入れると扱いが難しくなるのは当然であるが、相対論的 な流体の場合散逸を考慮するだけでも後に説明するように独特の困難が発生してしまう。 このため現在の解析では流体を完全流体として扱うか、衝撃波によるエネルギーの散逸 率のみを取り入れて計算する程度しか行なわれていない。宇宙流体の場合一般に熱伝導 や粘性などよりも輻射や磁場のほうが大きな影響を及ぼすことが多いが、散逸を入れる ことでMRI(磁気回転不安定性)由来の乱流のsaturation rateが定まったり、衝撃波や
reconnectionなど狭い領域で大量のエネルギーを散逸する現象をきちんと考慮すること が出来るなど、解析に大きな影響を与える可能性は大きいと期待される。 上記の相対論的な流体に散逸を考慮することで起きる困難は微分方程式の性質に主に由 来している。非相対論の場合の散逸を考慮した流体方程式であるNavier-Stokes方程式は 方程式系が放物型になっており、一般に因果律を守らない。非相対論の場合は速度に上限 が無かったため因果律の問題は発生しなかったが、光速の有限性が本質的である相対論的 力学ではこの因果律を守らないことは致命的であり、これに関係する非物理的で病的な不 安定性を誘起してしまう(セクション2.2.9を参照)。この問題を防ぐため近似を更に進 めた散逸項を用いる研究がなされているが、残念ながら実験の困難などからどの理論がど の程度正しいのかなどはわかっていない。そこで本研究では流体近似のもとの式に相当す る運動論的な方程式を研究することでこの困難の解決への指針を与える。
非相対論的な希薄気体の研究の歩み
運動論的な方程式の研究内容は広く知られてはいないので、運動論的方程式の研究の歩 みを簡単に説明しよう。 流体現象は通常Navie-Stokes方程式を用いて解析されるが、この方程式は現象の特性 長が平均自由行程より十分に短い場合には適用できなくなる。そのため宇宙物理の分野や プラズマ物理でしばしば興味をもたれる衝撃波の内部構造や波長が平均自由行程程度の音 波 [21] [22] [23]についてはNavier-Stokes方程式は全く予言能力を持たない。 このような現象を解析するには運動論的な解析が必要で、Boltzmann 方程式を解く必要がある。しかしこの方程式を一般的に解くことは非常に困難で、様々な近似法が研究 されてきた。その中で最も一般的に用いられているものはChapmanとEnskog により 開発されたもので、これは0次近似をボルツマン分布として分布関数をKnudsen数(平 均自由行程を特性長で割ったもの)で展開し逐次的に解いていくという手法である [24]。 この方法は 0次でEuler方程式を、1 次でNavier-Stokes方程式を再現し、それに加え 現象論的には求めることが出来なかった散逸係数の表式を与え、方程式のミクロな基礎 付けとなった。さらに高次まで展開することにより、より短波長領域の現象を記述する ことを期待されたが、Burnett [25] [26] などの研究によると、残念ながらKnudsen 数 が1に比べて十分小さくないような現象に関しての精度はほとんど上がらない。これは
Chapman-Enskog近似がKnudsen数について展開するためで、Knudsen 数が1に比べ
て十分小さくない場合は展開の収束性が悪くなるためである。また Burnett近似により 求まった項はとても複雑で非物理的な不安定性を含むため実用的ではない。分布関数をボ ルツマン分布のまわりで展開する近似としては、Chen [27] [28]らによる衝突項をBGK 近似 [18]を用い逐次近似を用いない方法が近年提唱されており、Navier-Stokes方程式よ りもKnudsen数が大きい領域でさらに良い近似になっているが、Knudsen数が1に近い か大きい現象に関しては精度は悪い。
別の手法としてはGrad [8] [9]により導入されたmoment methodがある。この手法は 分布関数をKnudsen数について展開するのではなくHilbert展開を用いて各momentの 発展方程式を解くというもので、momentとして密度、速度、温度、粘性、熱伝導を考慮
した 13 moment expansionではleading order でNavier-Stokes方程式を再現する。平
衡分布が5個のmoment しか持たないボルツマン分布であることから、一般に無限に現 れるmomentのうち平衡分布の5個以外は小さいと通常は期待しており、Knudsen数が 小さくない領域でも使えるための補正としてはmoment の数を増やすという方向が考え られ多くの研究がなされているが [29] [30] [31] [32]、やはりKnudsen数が1 付近より大 きい領域では正しい結果を導かない。またGradとは別の方法でmoment展開したもの としてはLevermore [33] [34] の方法があるが、一般的な場合に使えるかどうかは不明で ある。
相対論的な希薄気体の運動学の研究の歩み
非相対論的なボルツマン方程式の研究内容をそのまま応用する形で相対論的な運動論的 方程式の研究はすすんできた。運動論的方程式自体も理論的には興味深いが、ここでは特 に流体力学との関わりを中心に研究の歩みを説明しよう。相対論的な流体方程式の研究は Eckartの理論 [1]やLandau [2]らにより独立に行なわれた。完全流体を考慮した場合導 かれた方程式系は全て一致していたが、散逸を考慮した場合は四元速度の定義がユニークにならないなど、相対論的流体力学特有の性質が明らかになった。しかし最も深刻な性質 は、非相対論の場合の散逸項をそのまま相対論的な場合に適用した場合に発生する非物理 的で病的な不安定性の存在である [3]。この問題を解決するために流体力学の散逸項を現 象論的に求めるのではなく、流体近似のない運動論的な方程式から求めようと考えるの は自然であろう。まず非相対論的な散逸項の拡張は、そのままChapman-Enskogの方法 を相対論的な運動論的方程式に適用することで得られた。しかしChapman-Enskogの方 法は一次で一般に方程式系を放物型にし、近似を上げると微分の階数が増えてしまい扱 いにくい微分方程式になってしまう。そこで IsraelとStewart はGradのmoment展開 を行なうことで問題の解決を図った [4]。moment展開は一般に方程式系を双曲型にする (付録B 参照)ため、因果律を守り安定な方程式系が得られることが期待されており、確 かにパラメータを適切に選ぶことで安定で因果律を守ることがわかっている。ところが moment展開はその近似法からどの程度まで良い近似なのかがはっきりしなく、実験結果 も無いため理論の正しさがわかっていない。
本論分の目的と構成
私はこの修士論文で運動論的な方程式の線形摂動を行い分散関係を得て、この運動論的 に正しい分散関係式と得られている相対論的な散逸理論を比較することで散逸理論の精度 を判定する。ただしボルツマン方程式は散逸項の複雑さから分散関係を得ることが難しい ため、BGK近似を用いたボルツマン方程式の摂動方程式の厳密解を求め、それを用いて 流体方程式系の解析を行う。この方法は厳密解を用いているためBGK近似が成り立つ精 度で厳密に正しい分散関係が得られる。 BGK近似がどれだけ正しいかについて非相対論の場合には、波長依存性についての 減衰波の実験 [21]との対比により一定の精度が得られることが分かっている [18]。また BGK近似の導出過程から高次のmoment については近似が悪くなることも分かってい る [19]。BGK近似の衝突項の意味は物理的で明白であり、また相対論へ拡張されたモデ ルも非相対論の場合と大きく変わっていないため非相対論の場合と同様の近似の精度であ ることが期待される。また得られた分散関係とボルツマン方程式の衝突項自身が持つ性質 から相対論的なBGK近似の精度についても議論する。 本文の構成は、まず第二章で非相対論的な運動論と相対論的な運動論の、特に散逸項の 導出を中心としたreviewを行なう。そして第三章でBGK近似を用いたボルツマン方程 式の線形摂動の方程式の一般解を導出し、それを用いて流体方程式系の分散関係を導出す る。この章ではまず最初にBGK近似の特徴と散逸の性質を知るために非相対論の場合を 解析する。そして続いて相対論的な場合のBGKモデルを用いて粒子流速の静止系と熱流 の静止系でそれぞれ分散関係式を求め、この結果とIsrael-Stewart理論の分散関係を比較し議論をする。最後に第四章でまとめを行い、実用的な相対論的な散逸の理論について議 論を行なう。計算の詳細や発展的な内容は付録に載せた。
第
2
章
運動論
この章では運動論のreview をする。まず非相対論的な場合でのBoltzmann 方程式を 導出し、その平衡解を求める。続いて巨視的方程式に移行し、散逸項を運動論的に求め る。次に相対論的な運動論について reviewを行い、同様にして巨視的方程式を導出し、 そのときの散逸がどのように求まるかを見る。最後に求まった散逸項がどのように問題が 表れるかをみる。 この修士論文全体では電気的中性の原子あるいは分子からなる気体の運動論について扱 う。このため粒子同士の相互作用は常に直接衝突の間の非常に短い時間でのみ働くとし、 それ以外の時間は自由運動していると仮定する。物理的に表すと、分子間の平均距離を ¯ r ∼ N−1/3 (N :単位体積中の分子数)とし、自身の固有の大きさ、つまり分子間の力の 作用半径をdとしたときに、気体度パラメータN d3 ∼ (d/¯r)3 が1より十分小さいと仮定 する。この仮定は気体が十分希薄であることを意味し、粒子間の衝突として2体衝突のみ を考えるという仮定に相当する。また考える気体粒子は原子もしくは単原子分子とする。 またこの修士論文では一貫して古典論のみを扱うが、エントロピーの定義などでプランク 定数が必要な部分はプランク定数を1とする自然単位系を用いている。またボルツマン定 数もkB = 1とする。2.1
非相対論的運動論
この章では非相対論的な運動論のreview をする。この review は主にランダウ-リフ シッツの物理的運動学 [5]に従う。また同じレベルの教科書としてはMihalas らの教科 書 [6]などがある。2.1.1
分布関数
力学において物理系の記述は、原理的には全ての粒子の運動方程式を解くことによって 求められる。しかし気体運動論においては考えている粒子数が多いため、通常このような ことはせずに統計的な記述が用いられる。気体の統計的な記述には、気体分子の位相空間 での分布関数f (t, q, p)が使われる。ここでq は分子の一般化座標(まとめて qで表す) で、pは分子の一般化運動量である。また分子の位相空間の体積要素をdµ = dqdpとする と、f dµは与えられた体積要素内にある分子の平均数を表す。 定義 1. f (t, q, p):分布関数 q:分子の一般化座標(まとめてqで表す) p:分子の一般化運動量 f dµ = f dqdp:与えられた位相体積要素内の分子の平均数 f の引数は一般化座標と一般化運動量なのでそれらの選び方は任意だが、この修士論文 ではqは分子の慣性中心の座標rを、pは分子の並進運動量pを意味するものとする。 また通常用いられる座標空間での粒子密度は次のように定義される。 定義 2. Z f (t, r, p)d3p = n(t, r):気体分子の空間分布の平均密度 (2.1) この定義からndV は体積要素dV 内の平均分子数を表すが、ここで Boltzmann方程 式を扱う場合無限小体積要素dV とは次のような意味を持つ。 定義 3. dV:無限小体積要素 d3¿ dV ¿ L3 d:分子間の力の作用半径 L:系の特性的な大きさ これは気体運動論において無限小の点が数学的な無限小ではなく物理的な大きさを持つ ことを表し、これにより気体運動論の予言能力に制限が課されることになる。つまり気体 運動論においてこの下限のスケール未満の現象は追うことが出来ず、分子の位置をせいぜ いその大きさの程度の精度で定義することを意味する。また時間スケールに関しても粒子 の平均速度でこの下限のスケールを横切る時間未満は分解できず、これにより時間スケー ルの下限を定めることになる。 時間と空間に物理的な下限を設けることは次のことを意味する。古典力学では気体粒子 の古典的軌道が定められると、2つの粒子の衝突は完全に決定され、決定論的に時間発展が求まることになる。しかしここで物理的なスケールの下限があるとこのスケール未満の 現象が追えなくなり、結果として衝突の結果も不定となり、その結果の確率についてのみ が予言できることになる。 ここまでで分子の大きさdと平均分子間距離 ¯r、系の特性長Lの相対関係は定まった が、平均分子間距離r¯と系の特性長Lの相対関係は任意であった。しかし分布関数によ り定義される密度nの特性はこの比により差が生じる。この修士論文では巨視的な量に 興味があるためr ¿ L¯ という場合のみを考えることにする。もしこのような比でなかっ た場合、考えている特性スケール内に十分な数の粒子数がいないことになり、nの揺らぎ がnの時間平均と同程度になってしまいnが巨視的な量にはならなくなってしまうので ある。
2.1.2
詳細釣り合いの原理
一方のpの値が区間dpの中にあり、もう一方が区間dp1の中にあるような分子の衝突 を考える。この衝突の結果これらの粒子がそれぞれ区間dp0 、dp01にあるpの値に移った とする。(以降遷移p, p1 → p0, p01 と表す)この場合単位時間、単位体積当たりに起こる遷 移p, p1 → p0, p01 の衝突数は次のようになる。 w(p0, p01; p, p1) f dp f1dp1dp0dp01 (2.2) ここでw(p0, p01; p, p1)は遷移p, p1 → p0, p01 を伴う遷移確率を表し、力学の問題を解く ことにより得られる。この式は単位体積、単位時間当たりの衝突の数が衝突する粒子の数 f dpと衝突確率w、遷移後の区間dp0に比例することから得られる。なおf ≡ f (t, r, p)、 f1 ≡ f (t, r, p1)としてある。 またこのwを用いて衝突の有効断面積は次のように定義される。 定義 4. dσ:衝突の有効断面積 dσ = w(p 0, p0 1; p, p1) |v − v1| dp 0dp0 1 (2.3) wは次の一つの一般的な関係を満たす。 Z w(p0, p0 1; p, p1)dp0dp01 = Z w(p, p1; p0, p01)dp0dp01 (2.4) この関係は詳細釣り合いの原理と呼ばれる。 この原理は量子力学の言葉を用いると関係の導出が分かりやすくなる。量子力学では衝 突の個々の過程の確率振幅はユニタリー行列Sˆにより表される。ここでユニタリー性の 条件は次のようになる。X n Sin† Snk = δik (2.5) ここで特にi = k の場合 X n |Sni|2 = 1 (2.6) となる。|Sni|2 は遷移i → nなる衝突の確率なので、上式は与えられた初期状態から全て の可能な状態への遷移確率の和をとると1になる、という確率の規格化を意味する。 ところでユニタリー性の条件は次の形に書くことも出来る。 X n SinSnk† = δik (2.7) この場合i = kでは X n |Sin|2 = 1 (2.8) となる。この場合は与えられた終状態への全ての可能な遷移の確率の和が1に規格化され ていることを示している。 以上の2つからn = iという状態を変えない遷移を除くと次の等式が成り立つ。 X n6=i |Sni|2 = X n6=i |Sin|2 (2.9) この式をwを用いて表すと詳細釣り合いの原理の式になる。
2.1.3 Boltzmann
方程式の導出
気体運動論において系の状態は分布関数f により定められる。つまり系の時間発展を 追うにはこのf の時間発展方程式を求めればよい。なおこの論文では『衝突』という言葉 を通常の力学的な散乱の意味ではなく、運動論的な衝突という意味で用いる。運動論的な 衝突の意味は後ほど定義する。 まず衝突が全く無視できるような場合は各粒子が閉じた部分系になり、各々力学方程式 に従い発展することになる。位相空間内で粒子数は保存されるので、次の連続の式が成り 立つ。 ∂f ∂t + ∂ ∂q (f ˙q) + ∂ ∂p(f ˙p) = 0 (2.10) 微分を展開すると次のようになる。 ∂f ∂t + µ ˙q∂f ∂q + ˙p ∂f ∂p ¶ + f µ ∂ ˙q ∂q + ∂ ˙p ∂p ¶ = 0 (2.11)ここで衝突が無視できることから、Hamilton方程式を考える。 ˙q = ∂H ∂p, ˙p = ∂H ∂q (2.12) この式より次のような関係が成り立つ。 ∂ ˙q ∂q = ∂2H ∂q∂p = − ∂ ˙p ∂p (2.13) これを用いると式(2.11)の第三項が0になり、次の式が得られる。 D Dtf = ∂f ∂t + µ ˙q∂f ∂q + ˙p ∂f ∂p ¶ = 0 (2.14) この式はリウヴィルの定理と呼ばれていて、粒子数が位相空間上の軌道に沿って一定で あることを示している。微分D/Dtは位相空間上の軌道に沿って取っていると考える。 衝突がある場合はリウヴィルの定理は成り立たず、衝突により位相空間上での軌道から 外れる粒子達が現れる。これらの寄与を考慮して輸送方程式を次のように書く。 ∂f ∂t + v · ∇f + ˙p · ∂f ∂p = µ ∂f ∂t ¶ coll (2.15) ここで(∂f /∂t)collは衝突により分布関数が変化する速さである。以降この項を衝突積 分と呼ぶ。f の発展が解けるようにするためにこの衝突積分の具体的な形を求める。 衝突によるf の変化とは、衝突により考えている位相空間上の点 (x, v)から出て行く (転出事象)か入ってくる(転入事象)の2つである。今考えている理論の精度では衝突 において位置は変化せず変化するのはpなので、運動量空間での転出・転入事象を考えれ ば十分である。 まず転出事象において遷移p, p1 → p0, p01 を伴う衝突を考える。pが与えられp1、p0、 p0 1 が可能な全ての値を取るような、単位時間に体積dV で起こる衝突の総数は dV dp Z w(p0, p0 1; p, p1)f f1d3p1d3p0d3p01 (2.16) に等しい。 同様に転入事象として遷移p0, p01 → p, p1 を伴う衝突を考えればよい。するとこのよう な衝突の総数は次のようになる。 dV dp Z w(p, p1; p0, p01)f0f10d3p1d3p0d3p01 (2.17) 転出事象数から転入事象数を差し引くと、全衝突の結果考えている位相空間上の点の粒 子数は毎秒 dV dp Z (w0f0f10 − wf f1)d3p1d3p0d3p01 (2.18)
だけ増加する。ここで簡単のため次のように置きなおしてある。 w ≡ w(p0, p01; p, p1), w0 ≡ w(p, p1; p0, p01) (2.19) 以上により衝突積分として次の式が得られる。 µ ∂f ∂t ¶ coll = Z (w0f0f0 1− wf f1)d3p1d3p0d3p01 (2.20) 第2項に対して詳細釣り合いの原理(2.4)を用いて変形すると次のような形になる。 µ ∂f ∂t ¶ coll = Z w0(f0f10 − f f1)d3p1d3p0d3p01 (2.21) 以上で衝突積分が求まったので輸送方程式が求まった。この修士論文では中性粒子のみ を扱うので外場による運動量の変化を0にする( ˙p = 0)と、次のボルツマン方程式が得ら れる。 ∂f ∂t + v · ∇f = Z w0(f0f10 − f f1)d3p1d3p0d3p01 (2.22) また衝突積分は式(2.3)を用いて変形することができ、次の形に表すことも出来る。 µ ∂f ∂t ¶ coll = Z vrel(f0f10 − f f1)dσd3p1 (2.23) 分布関数のところで述べたように、この方程式は分子間の力の作用半径d程度のスケー ル未満は分解できず、このスケール以上の現象についてのみ予言能力がある。
2.1.4 H
定理
考えている気体が孤立系の場合、エントロピーが最大の状態へと移行して平衡分布とな る。先に導出したボルツマン方程式が物理的な理論であるためにはこの時間発展を示す必 要がある。この章ではこれを示す。この章は一般的なボルツマン方程式の議論をしたいの で外場がある場合を考える。 非平衡状態の気体のエントロピーは分布関数f を用いて次のように表される。 S = Z f ln e f dV d 3p (2.24) この式を時間微分するとリウヴィルの定理より次のようになる。 dS dt = Z ∂ ∂t µ f ln e f ¶ dV d3p = − Z ln f∂f ∂tdV d 3p (2.25) ボルツマン方程式の導出過程から自由運動部分は力学的に求めているため時間反転対象 性を保つと考えられ、エントロピー増大による時間反転の非対称性は衝突項に結びついて いると考えられる。これは次のようにして示せる。ボルツマン方程式より式(2.25)に含まれる時間微分は次のように表せる。 ∂f ∂t = −v · ∇f − F · ∂f ∂p + µ ∂f ∂t ¶ coll (2.26) 自由運動に当たる最初の二項を式(2.25)に代入すると次のようになる。 − Z ln f · −v · ∇f − F · ∂f ∂p ¸ dV d3p = Z · v · ∇ + F · ∂ ∂p ¸ µ f lnf e ¶ dV d3p (2.27) この式はガウス積分により表面項になる。孤立系でエネルギーは有限な気体を仮定すると 位置空間と運動量空間の無限遠でf = 0となり、上の積分は0となる。以上により自由運 動部分はエントロピーを増大させないことが示された。 以上の議論によりエントロピーの変化は次のようになる。 dS dt = − Z ln f µ ∂f ∂t ¶ coll dV d3p (2.28) この積分を考えるにあたり、次のような更に一般の形を用いて変形していこう。セク ション 2.1.1のように位相空間上の座標をq、pで代表させる。次の形の積分を考える。 Z φ(p) µ ∂f ∂t ¶ coll dp (2.29) φはpの任意関数である。式(2.20)を代入すると次のようになる。 Z φ(p) µ ∂f ∂t ¶ coll dp = Z φ w(p, p1; p0, p01)f0f10d4p − Z φ w(p0, p01; p, p1)f f1d4p (2.30) 簡単のため d4p = dpdp1dp0dp01 としてある。4種類の pは全て積分変数なので積分の値 は変えずに変数の任意の置き換えができる。第2項でp、p1 とp0、p01 を入れ替えると次 の式になる。 Z φ(p) µ ∂f ∂t ¶ coll dp = Z (φ − φ0)w0f0f0 1d4p (2.31) 同様にして第2項でp、p0とp1、p01の入れ替えを行い、得られた積分の和を1/2倍し、 関数wの対称性を用いると次のようになる。 Z φ(p) µ ∂f ∂t ¶ coll dp = 1 2 Z (φ + φ1− φ0− φ01)w0f0f10d4p (2.32) この φにm、p、(v − u)2/2を代入すると保存則より0になることがわかる。φ = 1を 代入した式は次のように書き換えられる。 Z µ ∂f ∂t ¶ coll dp = Z w0(f0f10 − f f1)d4p = 0 (2.33)
式(2.32)を用いて式(2.28)を変形すると次のようになる。 dS dt = 1 2 Z w0f0f10lnf 0f0 1 f f1 d4pdV = 1 2 Z w0f f1x ln x d4pdV (2.34) ここでも先のd4pの表記を用いている。またx = f0f10/f f1 である。この式から0である 式(2.33)の1/2倍を差し引くと次のようになる。 dS dt = 1 2 Z w0f f1(x ln x − x + 1) d4pdV (2.35) この被積分関数の括弧内はx > 0で負にならず、x = 1のときにのみ0になる。またw0、 f、f1 は正なので結局次の式が示された。 dS dt ≥ 0 (2.36) これはエントロピー増大則そのものであり、H定理と呼ばれる。 注意すべきこととしては、式(2.35)は被積分関数が正であるためd4pについて積分し ても正になる。つまり衝突は気体のエントロピーを各体積要素ごとに増大させることを意 味する。なおこれはあくまで衝突の効果のみを考えた場合の話であり、自由運動まで考慮 すると自由運動により気体のエントロピーが値はそのままで別の場所に運ばれるというこ とが起き得るため、一般にはエントロピーは各気体要素ごとに増大するとは限らない。 なおエントロピー増大則が出てきたのは、本質的には衝突項を式(2.21)の形にしたこ とに起因している。これはボルツマン方程式の更に基礎的な理論であるリウヴィルの定理 から動力学的にボルツマン方程式を導出する際にmolecular chaos仮定を入れた結果であ る。この仮定により衝突現象がそれまでの物理過程によらず f f1vreldσd3p1d3pdV dtと 記述でき、この『それまでの物理過程を無視した』ということから時間反転の非対称性が 入ったのである。詳しい動力学的導出はランダウの教科書 [5]などを参照。
2.1.5
平衡解
式(2.22)の平衡解を求める。左辺には時間と空間微分しか含まれていないので、運動 量pのみに依存する任意の分布関数f (p)は左辺を恒等的に0にする。また右辺を恒等的 に0にするにはこの分布関数f (p)が次の関係を満たせばよい。 f (p0)f (p0 1) − f (p)f (p1) = 0 (2.37) この式のlogを取ると次のようになる。 ln f (p0) + ln f (p0 1) = ln f (p) + ln f (p1) (2.38)この式はln f (p)が2体衝突の運動の保存量の線形結合で書けることを示している。こ れらのうち運動量と角運動量は系の平行移動、回転により一般に消すことが出来るので、 ln f (p)は次の形に書くことが出来る。 ln f (p) = α + β² (2.39) α、β は定数で²は系のエネルギーで、第一項が粒子数保存を、第二項がエネルギー保 存を表している。 このf を式(2.1)を満たすように規格化し、平衡分布がエントロピーを最大にすること からβを定めると、ボルツマン方程式の平衡解は次のようになることが分かる。 f0 = ρ m(2πRT )3/2 exp h −² T i (2.40) R = 1 m (2.41) ここでmは粒子の質量、ρ = mnは質量密度、T は気体の温度を表している。この式は 統計力学で良く知られたMaxwell-Boltzmann分布であり、粒子同士の衝突によりf は Maxwell-Boltzmann分布へと緩和していくことをBoltzmann 方程式は示している。ま たこの解は式(2.35)の右辺を0にする解でもあり、物理的な平衡解にもなっていること が分かる。 この平衡分布について詳細に研究する理論としては統計力学がある。統計力学は膨大な 分野に応用されているためこのreviewでは説明は割愛する。統計力学の参考文献として はランダウ-リフシッツの物理学教程 [7]などがある。
2.1.6
巨視的方程式
先に得られたBoltzmann方程式は気体の状態の時間発展の微視的記述を与える。流体 力学方程式はその粗い巨視的な記述のはずであり、Boltzmann方程式から導出できるは ずである。この章ではそれがどのように行なわれるかを見ていく。 まず先に導入したスケールのほかに、次のようなスケールを導入する。 定義 5. 平均自由行程l:分子が相次ぐ二つの衝突の間に通過する平均距離 l ∼ 1 nσ (2.42) この量は、定義自身が気体中のどのような輸送現象を考えるかに依存するので、定性的 な意味しか持たないことに注意するべきである。衝突断面積は分子のdを用いてσ ∼ d2と書ける。n ∼ ¯r−3とすると、定義式よりlは 次のように書きなおせる。 l ∼ ¯r³ ¯r d ´2 (2.43) 気体中ではr À d¯ なので、l À ¯rとなる。 一般に流体力学近似は粒子同士が十分衝突できて、分布関数が平衡分布に近い場合に良 い近似になる。言い換えると気体の巨視的性質である温度、密度、速度などが変化する特 徴的な長さLが、気体の平均自由行程lに対してL À lを満たす程度に十分長い場合に 流体近似は良い近似になる。このことを見ていこう。 式(2.1)において気体の数密度を定義したが、同様にして次の二つの量を定義する。 定義 6. Z d3vvf (t, r, v) = u :気体分子の平均速度 (2.44) Z d3v1 2(v − u) 2f (t, r, v) = 3 2ρRT ≡ ¯² :気体分子の平均エネルギー (2.45) 巨視的な量n、u、¯²の発展方程式はBoltzmann方程式 ∂f ∂t + ∇ · (vf ) = µ ∂f ∂t ¶ coll (2.46) にm、p、(v − u)2/2をかけてd3vで積分することにより得られる。ここで衝突は、衝 突に関係する粒子数、全エネルギー、全運動量を変えないはずなので、巨視的な量を変 化させないはずである。実際式(2.32)においてφにm、p、(v − u)2/2 を代入すると保 存則より0になることがわかる。このことを用ると、Boltzmann方程式の両辺にm、p、 (v − u)2/2をかけてpで積分したものは次のようになる。 ∂ρ ∂t + ∇ · (ρu) = 0 (2.47) ∂ ∂tρuα + ∂Παβ ∂xβ = 0 (2.48) ∂ ∂tn¯² + ∇ · q = 0 (2.49) ここで Παβ は運動量流束密度テンソル、qはエネルギー流束を表し、それぞれ次のよう に定義される。 Παβ = Z mvαvβf d3p (2.50) q = Z 1 2(v − u) 2 vf d3p (2.51)
これらの式から流体力学の方程式を導出しよう。そのためには分布関数f の形を定め、 Παβ、qを巨視的な量で表さなければならない。そこで平均自由行程をl、系の特徴的な スケールをLとして、先に述べた通りl ¿ Lを仮定する。この時特徴的な時間スケール L/¯vでは衝突は十分な回数が起こり、系は平衡に十分近いと考えられる。そこで第0近似 として、気体は全体としては平衡にはなっていないが、局所的には熱平衡になっていると 考えられる。つまり分布関数として次の局所的平衡関数f0(t, x)を考える。 f0(t, x, v) = ρ(t, x) m(2πRT (t, x))3/2exp · − 1 2RT (t, x) (v − u(t, x)) 2 ¸ (2.52) ここで質量密度ρ、平均速度u、温度T は真のf と同じものをとることにする。局所的平 衡関数はBoltzmann方程式の右辺を恒等的に0にするが左辺は一般に0にならず、正確 な意味でBoltzmann方程式の解にはなっていないことに注意するべきである。このこと は次の散逸の取り扱いで詳しく述べる。 まずΠαβ、qから巨視的な速度uを引き抜こう。気体と共に運動する座標系K0 系の速 度をv0 とすると、v = v0+ uとなる。するとまずΠαβ は
Παβ = ρhvαvβi = ρh(uα+ vα0)(uβ + v0β)i = ρ(uαuβ + hvα0 vβ0i) (2.53)
となる。K0 系は気体静止系で局所平衡分布では等方的であることを用いた。またこれ より hv0αv0βi = 1 3hv 02iδ αβ = T mδαβ (2.54) となる。ここで熱速度の自乗の平均がhv02i = 3T /mであることを用いた。これをΠαβ の式に代入し、nT = P とするとΠαβ は巨視的な量から次のように表される。 Παβ = ρuαuβ+ δαβP (2.55) この式はP を圧力とみなせばオイラー方程式のと同等である。つまり気体の圧力は以上 のように定義される。Boltzmann方程式の仮定として十分希薄であると仮定したことか ら理想気体の状態方程式が得られたことに注意するべきである。 次にエネルギー流速qについても同様のことを行なう。気体と共に運動する座標系K0 系の速度をv0とすると次のような関係が得られる。 1 2mv 2 = 1 2mv 02 + mu · v0+ 1 2mu 2 (2.56) この式をqの定義式に代入すると次のようになる。 q = nu · mu2 2 + m 3 hv 02i + ¯²0 ¸ = u µ ρu2 2 + h ¶ (2.57)
この式は理想流体のエネルギー流速の式である。ここで気体の単位堆積あたりのエンタル ピーh = P + n¯²0 を用いた。 以上により、オイラー方程式を運動論的に導出することが出来た。このオイラー方程式 や次の章で導出する散逸項まで含めたNavier-Stokes方程式を扱う理論は流体力学と呼ば れる。この reviewでは流体力学は扱わないが参考文献としてはランダウの流体力学 [2] を挙げておく。
2.1.7 Chapman-Enskog
展開
Navier-Stokes方程式を導出するには、f のf0 からのずれを考慮に入れなければならな い。そこでまず分布関数f をl/Lで展開することを考える。f を局所平衡関数f0 の周り で展開すると次のようになる。 f = f0+ f1 l L + f2 µ l L ¶2 + · · · (2.58) ここで局所平衡関数に関する付加条件を考慮する。局所平衡関数は式の形は式(2.52) で定義されるが、この中のρ、u、T は任意の関数でも衝突項を0にする。この式を意味 のあるものにするためには、これらの任意関数を物理的に決めてやらなくてはならない。 今考えている系は少なくとも局所的には、つまり特徴的なスケール Lよりも十分小さい スケールでは平衡状態になっている。つまりこの系ではf ∼ f0となっていて、このf0か らのずれを全て無視していることになる。言い換えると次の matching conditionと呼ば れる条件からf0 を決めていることになる。 Z m(f − f0)d3p = 0, Z vα(f − f0)d3p = 0, Z 1 2mv 2(f − f 0)d3p = 0 (2.59) 密度、平均速度に関しては定義通りだが、温度T に関しては注意が必要である。本来温 度とは平衡分布であるMaxwell-Boltzmann分布のときにのみ存在する概念であり、非平 衡分布では定義されていなかった。このためmatching conditionは非平衡気体での温度 の定義になっており、逆に言うと非平衡気体で温度を定義するためにはこの付加条件が必 要とされる。 展開したf をBoltzmann方程式に代入しよう。この際式の近似の次数はl/Lによって 決まるが、Boltzmann 方程式の左辺に含まれる微分を無次元化するために両辺に平均自 由行程lを掛けておく。すると左辺は全て微分が入っているので、l/Lの次数は1次から 始まることになる。また右辺の衝突項は式 (2.23)より、f vrel dσ ∼ 1/lが打ち消され、 衝突項の l/Lの次数はそのまま展開したf の次数に等しくなる。よって次数別に式を立てると次のようになる。 0 = µ ∂f0 ∂t ¶ coll (2.60) l D Dtf0 = µ ∂f1 ∂t ¶ coll l L (2.61) l D Dtfn−1 µ l L ¶n−1 = µ ∂fn ∂t ¶ coll µ l L ¶n (2.62) D/Dtは位相空間上の軌道に沿った微分を表している。 以上によりこの連立微分方程式をmatching conditionを課して解いていけばいい。こ の式はf0 を求めればそれ以降は漸化式として解いていけばよく、l/Lが十分小さければ 収束する数列をなしている。またfnがf0 により決まるということは、fnに関係する物 理量は全てf0 の物理量ρ、u、T の勾配と関係付けられることになるということに注意す るべきである。以上のようなBoltzmann方程式の解の求め方はChapman-Enskogの方 法と呼ばれる。(有限のfnで打ち切った場合はChapman-Enskog近似と呼ばれる)具体 例と実用的な例として以降ではf1について考えていく。 まずf1 を次のように書く。 f = f0+ f1, f1 = −∂f0 ∂² χ(p) = 1 Tf0χ (2.63) ここで運動量をpで代表させる。 このχはmatching conditionより、次の条件を満たさなければならない。 Z mχd3p = 0, Z vαχd3p = 0, Z 1 2mv 2χd3p = 0 (2.64) 分布関数はcollisionによりまず局所的に局所平衡f0に収束するべきだが、このmatching conditionを課さなければf0がユニークに定まらず、f0 からρ、u、T のみが異なる別の f0 0 への緩和を記述する解が含まれてしまうことを注意するべきである。 ではChapman-Enskogの方法に従いf1を求める式を導出しよう。まず0次の方程式 0 = µ ∂f0 ∂t ¶ coll (2.65) は局所平衡関数(2.52)が解になる。 ではこのf0 を1次の方程式の左辺に代入しよう。この際簡単のためガリレイ変換を行 い、u = 0となる流体のある一点で考える。このあと求める散逸係数はガリレイ変換では 変化を受けないはずなので、一般にこのようにしても散逸係数は正しく求まる。すると少 し長い計算の後に次の式が得られる。 f0 T ½ ² − cpT T v · ∇T + · vαvβ R − δαβ ² cv ¸ Vαβ ¾ ≡ f0 T I(χ) (2.66)
ここで次の量を定義した。 ² ≡ v 2 2R, Vαβ ≡ 1 2 µ ∂uα ∂xβ + ∂uβ ∂xα ¶ (2.67) またI(χ)は衝突積分から出てくる項で、次のように定義される。 µ ∂f ∂t ¶ coll = f0 T I(χ) (2.68) I(χ) = Z w0f01(χ0+ χ01− χ − χ1)dp1 dp0 dp01 (2.69) また計算の途中で次の関係式を用いて時間微分を消し、T とVの勾配のみになるように した。 ∂ρ ∂t + ∇ · (ρu) = 0 (2.70) ∂u ∂t + u · ∇u = − ∇P ρ (2.71) ∂s ∂t + u · ∇s = 0 (2.72) cvdT − P ρ2dρ = T ds (2.73) P = ρRT (2.74) 上から連続方程式、オイラー方程式、熱輸送の一般式、熱力学第一法則、状態方程式であ る。式(2.66)をl/Lの1次の式の左辺に代入すれば原理的にf1が求まることになり、気 体の熱伝導係数と粘性係数を求めることが出来る。なお流体方程式において散逸項は完全 流体の項に比べて微小量にあたるため含まれていない。同様にして例えばf2 を求める場 合はf1から得られるNavier-Stokes方程式を代入すればよい。
2.1.8
気体の熱伝導
温度勾配と速度勾配が線形独立だと仮定して、式(2.66)の第一項を考えよう。 ² − cpT T v · ∇T = I(χ) (2.75) この方程式の解を次の形で求める。 χ = g · ∇T (2.76) ここでgは運動量のみの関数である。このようなχを考えると、式(2.75)は∇T の係数 が等しいことを示す方程式になり、次のようになる。 ² − cpT T v = I(g) (2.77)ここでmatching condition(2.64) について考えよう。まず理論には特別な方向はない と考えられるので、ベクトルパラメータを与える積分R f0gdpと R f0²gdpは0になるこ とが分かる。結局matching conditionより得られる付加条件は次のようなものになる。 Z f0vgdp = 0 (2.78) これ以上具体的にχを求めることはせずに、熱伝導率の表式について考えていこう。巨 視的運動がない場合、式(2.51)より熱伝導はf = f0+ δf を代入すると次のようになる。 q = 1 T Z 1 2(v − u) 2vf 0χdp = 1 T Z 1 2(v − u) 2vf 0(g · ∇T )dp (2.79) ここで平衡分布は等方的であることを用いた。この式を成分で書くと次のようになる。 qα = −καβ T xβ (2.80) καβ = −1 T Z 1 2(v − u) 2f 0vαgβdp (2.81) ここで平衡気体の等方性のために καβ は単位テンソル δαβ を用いて次のようなスカ ラーで書ける。 καβ = κδαβ, κ = καα 3 (2.82) 以上より熱伝導は次のようになる。 q = −κ∇T (2.83) ここでスカラー熱伝導率は次のようになる。 κ = − 1 3T Z 1 2(v − u) 2f 0v · gdp (2.84) この量は輸送方程式の性質からも正である事が示せる。これは現象論ではエントロピーの 増大側から要請されていた条件で、ボルツマン方程式ではH定理を満たすことから自然 にこの条件が要請されることに注意しよう。 最後にf1 に比例する項g から求まった熱伝導率κが、Chapman-Enskogの方法の要 請通り平均自由行程lに比例することを示そう。まず衝突項は次のように評価できる。 µ f t ¶ coll ∼ −v¯ l(f − f0) = − ¯ vf0 l T χ (2.85) ここで式(2.66)を用いた。この式と式(2.68)を比較すると結局次の関係が得られる。 I ∼ −v¯ lχ (2.86)
この式と式(2.75)を比べると、次の関係が得られる。 g ∼ l (2.87) この式と式(2.66)からf1 ∼ f0l/Lが得られ、これはChapman-Enskogの方法の要請を 満たしていることが分かる。このgを用いてκを評価すると次のようになる。 κ ∼ cnl¯v (2.88) ここでcは分子一個あたりの気体の熱容量である。この表式より熱伝導は平均自由行程l に比例し、平衡に到達する前のlが大きいような状況において大きくなることが分かる。
2.1.9
気体の粘性
先と同様に温度勾配と速度勾配が線形独立だと仮定して、今回は式(2.66)の第二項を 考えよう。 まず流体力学より、粘性は運動量流速密度テンソルに含まれ、次のように二種類の粘性 が定義される。 Παβ = P δαβ + ρuαuβ− σαβ0 (2.89) σ0 αβ = 2η µ Vαβ − 1 3δαβ∇ · u ¶ + ζδαβ∇ · u (2.90) 式(2.90)の第一項は流体の体積変化に関係しない変形に由来する粘性で、shear viscosity と呼ばれる。第二項目は逆に体積変化に由来する粘性で、bulk viscosityと呼ばれる。こ れらを別々に計算するために式(2.66)を次のように変形しよう。 mvαvβ µ Vαβ − 1 3δαβ∇ · u ¶ + µ mv2 3 − ² cv ¶ ∇ · u = I(χ) (2.91) まずはshear viscosityを計算しよう。そのために式(2.91)で∇ · u = 0と考える。す ると次のような式が得られる。 m µ vαvβ − 1 3δαβv 2 ¶ Vαβ = I(χ) (2.92) この式の左辺のテンソルは両方ともtraceが0であることに注意しよう。 この方程式の解を次の形で求めよう。 χ = gαβVαβ (2.93) ここで gαβ は運動量のみの関数である。またVαβ のtraceが0であることから、δαβ に 比例する項を加えることでχを変えずにgαβ のtraceを0にしてある。このようなχを考えると、式(2.92)はVαβ の係数が等しいことを示す方程式になり、次のようになる。 m µ vαvβ − 1 3δαβv 2 ¶ = I(gαβ) (2.94) このgαβ はスカラー、ベクトル部分は含まれていないのでmatching conditionは自動的 に満たされる。 熱伝導の場合と同様にしてgαβ の具体的な値を求めることはせずに、粘性係数の表式が どのように与えられるかを見ていこう。式(2.50) より、運動量流速密度テンソルのうち 散逸に関係する部分は次のように書ける。 σ0αβ = −m T Z vαvβf0χdp = ηαβγδVγδ (2.95) ηαβγδ = −m T Z vαvβf0gγδdp (2.96) 四階テンソルηαβγδ は添え字対{α, β}、{γ, δ}についてそれぞれ対称で、対{γ, δ}は 縮約すると0になる。気体の等方性よりこのテンソルは単位テンソルδαβ のみで表され るはずで、それは次のようになる。 ηαβγδ = η · δαγδβδ + δαδδβγ − 2 3δαβδγδ ¸ (2.97) η = − m 10T Z vαvβf0gαβdp (2.98) この表式を用いるとσαβ0 = 2ηVαβ となり、η は求めるスカラー粘性係数であることが分 かる。 ではこの粘性係数η を評価しよう。熱伝導のときと同様な評価を行なうと次のように なる。 η ∼ m¯vnl (2.99) この表式より粘性係数と熱伝導係数は同程度の大きさであることがわかる。具体的には次 のように書ける。 κ nc ∼ η mn ∼ ¯vl (2.100) bulk viscosityの粘性係数を求めるには式(2.91)の第二項が0でない場合を考えればよ く、次の式を考えることになる。 µ mv2 3 − ² cv ¶ ∇ · u = I(χ) (2.101) ここで次のようなχを考える。 χ = g∇ · u (2.102)
shear viscosityの場合と同様に式(2.101)は∇ · uの係数が等しいことを示す方程式にな り、次のようになる。 mv2 3 − ² cv = I(g) (2.103) 応力テンソルσ0αβ をshear viscosityの場合と同様に計算し、ζδαβ∇ · uと比べると、粘 性係数は次のように求まる。 ζ = −m 3T Z v2f0gdp (2.104) ここで単原子気体の場合を考えよう。この場合 ² = mv2/2、cv = 3/2 なので、式 (2.101) の左辺は 0になる。すると式 (2.103)はI(g) = 0 となり、g = 0 が得られる。 これを式(2.104)に代入するとζ = 0となる。つまり非相対論の場合は単原子気体では bulk viscosityの粘性係数は0になることが分かる。 このreviewでは散逸係数を具体的には求めなかったが、散逸係数の具体的な形と実際 の求め方についてはランダウの教科書 [5] やMihalas [6]を参照。
2.1.10 moment
展開
ボルツマン方程式の解を求める方法として、前の章ではChapman-Enskogの方法を説 明した。この方法は分布関数f を平均自由行程l と系の特徴的長さスケールLの比 l/L で展開するという手法で、流体近似のような長波長領域では良い近似になっていたが、比 が小さくなくなるような短波長領域ではl/Lによる展開の収束が悪くなり近似が悪くな る。この問題はl/Lで展開する限り避けられない問題である。 この章ではl/Lをあらわには用いない別のf の展開を説明する。これは分布関数f をHermite polynomialsを用いて展開するもので、Grad [8] [9]らにより導入された。
まず局所平衡分布(2.52)を考えよう。 f0(t, x, v) = ρ(t, x) m(2πRT (t, x))3/2 exp · − c 2 2RT (t, x) ¸ (2.105) c ≡ v − u (2.106) この関数は matching condition により、分布関数f の五個の局所的な量で定義されて いる。 では分布関数f (t, x, v)をこの局所平衡分布関数f0(t, x, v)を重み関数として、ρ(t, x)、 u(t, x)、T (t, x)、a(n)(t, x)の無限この組に展開することを考えよう。この展開によりボ ルツマン方程式は ρ、u、T、a(n) の無限個のmoment の連立準線形微分方程式になる。 全ての微分方程式で時間微分は一階の線形微分となっている。また保存則に対応するρ、
u、T の時間微分方程式以外は衝突項からa(n) の二次形式が現れる。一般にこれらの方程 式は単一のmoment では記述されず、厳密な解を得るには全ての方程式を解かなければ ならない。そこで実用的な解法としては適当なa(m)までを考え、残りはa(n) = 0として 解く。一般にこの解法が良い近似になっている保証はなく、またこのmoment展開の収 束性に関してはわからないが、実験値との比較によれば長波長では比較的少ないmoment 数で良い近似になることがわかっている。 展開を考える前に物理量を無次元化しておこう。 ξ = √c RT (2.107) g = (RT ) 3/2 n f (2.108) N、u、T の定義より、gには次のような条件が課される。 Z gdξ = 1 (2.109) Z ξgdξ = 0 (2.110) Z ξ2gdξ = 3 (2.111) また無次元化された局所平衡関数は次のようになる。 g0 = 1 (2π)3/2 exp · −ξ 2 2 ¸ (2.112) gをHermite polynomialsを用いて次のように展開しよう。 g(ξ, x, t) = g0(ξ) ∞ X n=0 1 n!a (n) i (x, t)H (n) i (ξ) (2.113) = g0(ξ) µ a(0)i H(0)+ a(1)i H(1)i + 1 2!a (2) ij H (2) ij + · · · ¶ f に直した場合は次のようになる。 f (ξ, x, t) = f0(ξ) ∞ X n=0 1 n!a (n) i (x, t)H (n) i (ξ) (2.114) Hermite polynomialsにおいてHi(n)はξ空間のn次のテンソルであり、a(n)i もn次の テンソルになっている。i = (i1· · · in)は成分添え字で、今後は簡単のためCartesian座 標で考える。またHermite polynomialsは直交多項式で、次のような直交関係を満たす。 a(n)i = Z gH(n)i d3ξ = 1 n Z f H(n)i d3v (2.115)
最初のいくつかのH(n)i は次のようになっている。 H(0) = 1 (2.116) H(1)i = ξi H(2)ij = ξiξj− δij Hijk(3) = ξiξjξk− (ξiδjk+ ξjδik+ ξkδij) H(4)ijkl = ξiξjξkξl− (ξiξjδkl + ξiξkδjl+ ξiξlδjk + ξjξkδil+ ξjξlδik+ ξkξlδij) +(δijδkl + δikδjl+ δilδjk) H(n)i の特徴としては、H(m)j とはn = mでかつ添え字j が添え字iを入れ替えものであ る場合を除けば互いに独立で直交する。またξi の係数としては±1以外は現れない。 a(n)i については最初のいくつかは次のように定義される。 a(0) = 1 (2.117) a(1)i = 0 a(2)ij = pij p a(3)ijk = Sijk p√RT a(4)ijkl = Qijkl pRT − 1 p(pijδkl + pikδjl+ pilδjk + pjkδil+ pjlδik+ pklδij) +(δijδkl + δikδjl+ δilδjk) これらは定義によりa(0)を除けばf = f0の場合に全て0になる。 a(3) の添え字を縮約することで次のように熱伝導を得る。 a(3)i = a(3)irr = 2qi p√RT (2.118) このようにして Hermite 展開の係数で最初の意味のある momentとして粘性テンソル pij と熱流ベクトルqi が得られた。注意するべきことは、これらを特に取り上げたのは流 体近似で我々が特に慣れ親しんでいるからであり、これら以外の a(n)i も物理的に同等の 価値を持っているのである。 式(2.113)を用いてボルツマン方程式を展開するとa(n)i についての無限階の連立微分
方程式になる。例えばa(2)ij については次のようになる。 ∂a(2)ij ∂t + ur ∂a(2)ij ∂xr + a (2) ir ∂uj ∂xr + a (2) rj ∂ui ∂xr + (a (2) ij + δij) 1 RT DRT Dt (2.119) +√RT∂a (3) ijr ∂xr + √ RT n a (3) ijr ∂N ∂xr + 3 2RTa (3) ijr ∂RT ∂xr + ∂ui ∂xj + ∂uj ∂xi = J (2) ij D/Dtはラグランジュ微分である。またJij(n) は衝突項からの寄与であり、次のようにし て定義される。 Ji(n) = N 2 ∞ X r,s=0 βijk(nrs)a(r)j a(s)k (2.120) i = (i1· · · in), j = (j1· · · jr), k = (k1· · · ks) βijk(nrs) = 1 r!s! Z 1 mB µ θ,|ξ√1− ξ| RT ¶ g0(ξ)g0(ξ1)H(r)j (ξ)H (s) k (ξ1)[H(n)i (ξ)]dθ d² dξ dξ1 (2.121) 上の式で次の定義を用いた。 [φ] = φ0+ φ01− φ − φ1 (2.122) Bは相互作用に関係する関数である。 式(2.120)により J(n) は一般にa(n)i の二次形式の無限和になっていることがわかる。 そして式(2.119)よりa(2)i の時間微分方程式は全てのa(n)i が含まれており、a(2)i の時間 発展を解くには全ての連立微分方程式を解く必要があることが分かる。このことは一般の a(n)i の時間微分方程式にも成り立つ。なお相互作用がポテンシャルU ∝ r−4 に従う場合 は特別で、この場合は式(2.120)は有限個の和で閉じることが知られている。(このよう な相互作用に従う粒子はMaxwell粒子と呼ばれる) このように一般の相互作用の粒子の場合、これらの方程式を解くためには分布関数f はある有限個のa(n)i までに依存しないと近似し、無限階連立微分方程式を有限階にする 必要がある。単原子気体の流体近似を考える場合に良く用いられるのは、moment とし て熱流、粘性までの13個のmoment を残す13 moment expansionである。なお単原子 でない場合や相対論的な粒子の場合はbulk viscosityが残るため、これを加えた14個の momentが考慮される。このように近似するとf は次のようになる。 f = f0 ½ 1 + pij 2pRTvivj − qi pRTvi µ 1 − v 2 5RT ¶¾ (2.123)
これにより粘性テンソルpij と熱流ベクトルqiの方程式として次の式が得られる。 ∂pij ∂t + ∂ ∂xr (urpij) + 2 5 µ ∂qi ∂xj + ∂qj ∂xi − 2 3δij ∂qr ∂xr ¶ (2.124) + pir ∂uj ∂xr + pjr ∂ui ∂xr − 2 3δijprs ∂ur ∂xs + p µ ∂ui ∂xj + ∂uj ∂xi − 2 3δij ∂ur ∂xr ¶ + βnpij = 0 ∂qi ∂t + ∂ ∂xr(urqi) + 7 5qr ∂ur ∂xj + 2 5qr ∂ur ∂xi + 2 5qi ∂ur ∂xr (2.125) + RT∂pir ∂xr + 7 2pir ∂RT ∂xr − pir N ∂Prs ∂xs + 5 2p ∂RT ∂xi + 2 3βnqi = 0 ここで Pij は応力テンソルである。この式は衝突のパラメータは次で与えられるβ しか 入っていないことが分かる。 β(RT ) = 2 5 √ 2π Z ∞ 0 x6e−x2/2 ½Z 1 mB(θ, x √ 2RT ) sin2θ cos2θ dθ ¾ dx (2.126) この β の解釈について考えてみよう。式(2.124)、(2.125)において全ての量の位置の 依存性を無視すると次のような式になる。 ∂pij ∂t + βnpij = 0 (2.127) ∂qi ∂t + 2 3βnqi = 0 (2.128) これらの式より、もし温度や平均速度の勾配がなかった場合は粘性 pij や熱流 qi は 1/βN 程度の時間で指数関数的に減衰することが分かる。ここでβ は定義式(2.120)から β ∼ vrelσ であることが分かる。このことから1/βnは自由走行時間程度であることが分 かる。 また式(2.124)、(2.125)をChapman-Enskog展開しよう。セクション2.1.7のf1 に対 応する項まで残すことにすると次のようになる。 p µ ∂ui ∂xj + ∂uj ∂xi − 2 3δij ∂ur ∂xr ¶ + βN p(1)ij = 0 (2.129) 5 2p ∂RT ∂xi + 2 3βN q (1) i = 0 (2.130) f1 までの精度という意味でp(1)ij 、q (1) i とした。これは一次のChapman-Enskog展開から 得られる粘性と熱伝導と同等であることがわかる。この式から散逸係数を計算すると次の