本研究はRT不安定性の線形及び非線形成長に対するHall項及びジャイロ粘性といった非理想 MHD効果の影響を解明するものである。RT不安定性は、第1章で述べたように、核融合プラズ マの圧力駆動型不安定性、バルーニング不安定性とも関連する現象である[20]。バルーニング不安 定性はトロイダル方向に対して最も不安定となっている。核融合プラズマのポロイダル断面の一部 に着目すると、圧力駆動型不安定性の駆動源となっている圧力勾配及び曲率効果の一部を一様重力 としてモデル化することでRT不安定性と類推することが可能である。(図6.1参照)
本論文の研究結果から、拡張MHDモデルを用いた高波数モードを高精度で表現する解析は核融 合プラズマでも必要であると考えられる。大型ヘリカル装置(LHD)では、圧力駆動型不安定性が 発生すると考えられる内寄せ配位の実験では、5%以上の高いβ値が得られている。従って、実験 結果は圧力駆動型不安定性には何らかの抑制・飽和機構が存在する事を示唆している。これまでに も実験結果を理論的に解明するために一流体モデルを用いた解析が行われているが、実験結果を
十分に説明する結果は得られていない。実験結果を十分に再現できていない事例としてMiura et al[56]及びSato[58]の事例を挙げる。
Sakakibara et alによる LHD内寄せ配位での実験[55]では、ある有理面での不安定性が存在す る間はβ値の上昇が制限されプラズマが不安定となるが、プラズマの崩壊までは至らない結果が得 られている。Miura and Nakajima[56]は一流体モデルを用いて実験と類似の初期平衡プロファイ ルによる非線形シミュレーションを行った。しかしながらシミュレーションでは、高波数モードが 早く成長しプラズマを不安定にした後に、遅れて成長してきた低波数モードによりプラズマコア部 が崩壊してしまう結果となった。この論文では中程度の(あるいはより短い)波長のモードに対し ては一流体モデルによる解析では十分に安定化されないことが指摘されている。
Sakakibara et alによる別の実験[57]では、プラズマ中心部では磁気井戸、周辺部では磁気丘と なる、周辺部でMHD不安定性が発生すると予測される磁場配位で実験が行われた。予測通り、周 辺部で圧力勾配駆動型不安定性と考えられるMHD不安定性が観測されたが、先の実験とは異な り中心部まで不安定性が影響を与えずに4.8%の高いβ 値が得られる結果が得られている。この 実験結果を数値シミュレーションで示すために、Sato[58]は一流体モデルを用いて、プラズマ中心 部では磁気井戸、周辺部では磁気丘となる磁場配位で数値シミュレーションを行った。このシミュ レーションにおいても、前述したシミュレーション研究と同様に、周辺部で発生したMHD不安定 性が中心部にまで影響しプラズマコア部が崩壊してしまう結果となった。このように、一流体モデ ルを用いた数値シミュレーションでは、LHDにおける実験結果を十分に再現出来ていないのが現 状である。一流体MHDモデルで圧力駆動型不安定性を上手く説明出来ていない理由として、圧力 駆動型不安定性の高波数モードに対しては、一流体MHDモデルでは考慮していないHall項及び ジャイロ粘性等の微視的効果がその成長に大きく影響することが指摘されている。従って、一流体 MHDモデルにこれらの微視的効果を取り入れた拡張MHDモデルを用いた3次元数値シミュレー ションを行う必要があるものと考えられる。
本研究を上記のようなトーラスプラズマに応用した場合には二次元の場合に比べて状況が複雑化 することが容易に想像出来る。バルーニング不安定性に応用した場合には高波数モードの固有関数 が動径方向に重なり合うことが予想される。一流体MHDモデルに対してジャイロ粘性を取り入 れた場合には、イオン温度勾配駆動型(ITG)モードが不安定化することが知られている。また、2 次元から3次元へと次元が増えることで非線形相互作用の組み合わせが増加し、形状効果(トロイ ダルカップリングなど)も起きると考えられる。しかしながら、以下で示されるような本研究で明 らかとなった基礎的な特性のいくつかは3次元の現象でも現れるものと考えられる。(いずれも、
Hall項及びジャイロ粘性を同時に加える事を考える。)
・Hall項及びジャイロ粘性が同時に存在することによる、高波数モードの成長率の大きな減少。
・高波数モードの安定化が低波数モードの成長に及ぼす影響。(低波数モードの線形段階の長期化)。
・密度変化が小さい安定性限界付近における、不安定性の発生する界面近傍での急峻な密度勾配の 形成及び二次不安定性(KH型不安定性)の成長。
・急峻な密度勾配、速度シアー層が発生する事に伴う、広い波数空間(高い数値解像度)を確保する
必要性。
この様な影響を明らかにするには、高波数モードを高精度で表現する必要があり、これにいくつか の困難が伴う事は第1章で述べた通りである。この論文では使わなかったが、核融合プラズマにお ける拡張MHDモデルを用いて高波数モードまで精度良く解析するための適応するための数値シ ミュレーションコードの開発を行った。このコード開発プロジェクトの企画書の概要については、
付録Bにまとめられている。この論文の発展として、この様な一般化座標系を用いた数値シミュ レーションにより、核融合プラズマの交換型モードなどに対するHall項、ジャイロ粘性の影響を 調べる事は、MHDに対する非理想的効果の解明のために有効であると考えられる。
図6.1 圧力駆動型不安定性とRT不安定性の類似性
謝辞
本研究を遂行し、学位論文をまとめるに当たり、多くの助言・御指導を頂きました、主任指導教 員である核融合科学研究所 三浦英昭准教授に心から感謝申し上げます。流体の数値シミュレー ションについて基礎的な内容から私に懇切丁寧に指導して頂きましてありがとうございました。不 安定性の線形解析手法についてご指導頂いた指導教員である核融合科学研究所 伊藤淳助教に感謝 申し上げます。研究内容や陰解法について数多くの助言をを頂きました、佐藤雅彦助教に感謝申し 上げます。審査委員として多くの助言・叱咤激励を頂きました、核融合科学研究所 榊原悟教授、
核融合科学研究所 藤堂泰教授、岡山理科大学 荒木圭典教授、日本原子力研究開発機構 相羽信 行博士には深く感謝致します。プラズマシミュレータでのプログラムの並列化について相談に乗っ て頂きました、株式会社日立製作所のエンジニアの皆様に感謝を申し上げます。また、研究を進 めるに当たって私生活面でも多くの助言・手助けを頂きました、同室の羽鳥智栄君、朝比祐一君、
Pianpanit Theerasarn君をはじめとする核融合科学研究所の学生の皆様、本当にありがとうござ
いました。
付録 A 線形安定性解析
第5章では、WKB近似による線形解析を行った。式(5.8)に現れる係数は、具体的には以下の ような表現をとる:
A1(x;ω, k) = β
2γp¯0+ ¯B02+δ
4(1−γ)τ βdi
L V∗
VAc
¯ p0V¯0′
B¯0
+λk¯di
L VAc
V∗ 1
¯ n0B¯0
B¯02+δ 4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0
Ω0−λ¯kdi
L VAc
V∗
¯ n′0
¯ n20B¯0
× [
B¯20 (¯n′0
¯
n0 − B¯0′ B¯0
) +β
2 (1−τ) (n¯′0
¯ n0
γp¯0−p¯′0 )]
(A.1)
B1(x;ω, k) = β 2γ
(
¯
p′0+kV¯0′ Ω0
¯ p0
) ( 1−δ
2 V∗
VAc
di
L τV¯0′
B¯0
)
+ 2 ¯B0B¯0′ +kV¯0′
Ω0
( B¯02+ δ
4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0
)
+ δ 4
V∗
VAc
di
L τ β B¯0
[(B¯0′ B¯0
V¯0′−V¯0′′
)
γp¯0+ ¯p′0V¯0′+ 2¯p0V¯0′′−B¯0′ B¯0
¯ p0V¯0′
]
+λ¯kdi
L VAc
V∗
1
¯ n0B¯0
1 Ω0−λk¯di
L VAc
V∗
¯ n′0
¯ n20B¯0
{(
B¯20+ δ 4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0 )
× {β
2γ(1−τ) 1
¯ n0
[(
−2¯n′0
¯
n0 −B¯0′ B¯0
+ kV¯0′ Ω0
)
¯
n′0p¯0+ ¯n′′0p¯0+ ¯n′0p¯′0 ]
−B¯0
(
B¯0′′+ kV¯0′ Ω0
B¯0′ )
+β
2 (1−τ) [(n¯′0
¯ n0
+B¯0′ B¯0 −kV¯0′
Ω0
)
¯ p′0−p¯′′0
]
+B¯02
¯ n0
[
¯ n′0
(
−2¯n′0
¯
n0 + 2 ¯B0′ B¯0
+ kV¯0′ Ω0
) + ¯n′′0
]}
+ {n¯′0
¯ n0
[β
2γ(1−τ) ¯p0+ ¯B02 ]
−B¯0B¯0′ −β
2 (1−τ) ¯p′0 }
× [
2 ¯B0B¯0′ + kV¯0′ Ω0
( B¯02+δ
4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0
)
+ δ 4τ β V∗
VAc di
L 1 B¯0
(
¯
p′0V¯0′+ ¯p0V¯0′′−B¯0′ B¯0
¯ p0V¯0′
)
+λ¯kdi
L VAc
V∗
1
¯ n20
B¯02+ δ 4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0 Ω0−λ¯kdi
L VAc
V∗
¯ n′0
¯ n20B¯0
× (kV¯0′
Ω0
¯
n′0B¯0+ ¯n′′0B¯0+ ¯n′0B¯0′ −2¯n′20
¯ n0
B¯0 )]}
(A.2)
C1(x;ω, k) = V∗2
VAc2 n¯0Ω20+ δ 4
V∗ VAc
di
Lτ βkΩ0 B¯0
(
¯ p′0−B¯′0
B¯0
¯ p0
)
+ β 2
(
¯
p′′0+kV¯0′ Ω0
¯ p′0
) ( 1−δ
2 V∗ VAc
di L
τV¯0′ B¯0
)
+ δ 4
V∗
VAc
di
Lτ βp¯′0 B¯0
(B¯0′ B¯0
V¯0′−V¯0′′
)
− Lg VAc2 n¯′0 + 1
B¯0
(
B¯0′′− B¯0′2 B¯0
) ( B¯02+δ
4τ β V∗
VAc di
L
¯ p0V¯0′
B¯0
)
+ B¯0′ B¯0
[kV¯0′ Ω0
( B¯20+ δ
4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0 )
+ 2 ¯B0B¯0′ +δ
4τ β V∗
VAc
di
L 1 B¯0
(
¯
p′0V¯0′+ ¯p0V¯0′′−B¯0′ B¯0
¯ p0V¯0′
)]
(A.3)
D1(x;ω, k) = δ 4τ β V∗
VAc
di
L
¯ p0
B¯0
Ω0 (A.4)
E1(x;ω, k) =−β 2kγp¯0
( 1− δ
2 V∗
VAc
di
L τV¯0′
B¯0
) +δ
4τ β V∗
VAc
di
L Ω0
B¯0
(
¯ p′0−B¯′0
B¯0
¯ p0
)
−k (
B¯02+δ 4τ β V∗
VAc di
L
¯ p0V¯0′
B¯0
)
−λk¯2VAc
V∗
di
L 1
¯ n0B¯0
B¯02+δ 4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0 Ω0−λk¯di
L VAc
V∗
¯ n′0
¯ n20B¯0
× [
B¯02 (n¯′0
¯
n0 −B¯0′ B¯0
) + β
2(1−τ) (n¯′0
¯ n0
γp¯0−p¯′0 )]
(A.5)
F1(x;ω, k) =−β 2γ
(
¯
p′0+kV¯0′ Ω0
¯ p0
) ( 1− δ
2 V∗ VAc
di L
τV¯0′ B¯0
)
+δ 4
V∗
VAc
di
Lτ βp¯0
B¯0
[
−kΩ0+γ (
V¯0′′− B¯0′ B¯0
V¯0′ )]
+ Lg VAc2 n¯0
−kV¯0′ Ω0
( B¯02+δ
4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0
)
− [
2 ¯B0B¯0′ +δ 4τ β V∗
VAc
di
L 1 B¯0
(
¯
p′0V¯0′+ ¯p0V¯0′′−B¯0′ B¯0
¯ p0V¯0′
)]
−λ¯kVAc
V∗
di
L 1
¯ n0B¯0
1 Ω0−λ¯kdi
L VAc
V∗
¯ n′0
¯ n20B¯0
{[(
B¯02+δ 4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0
)
×
−n¯′0
¯ n0
+2kV¯0′ Ω0
+λ¯kVAc
V∗
di
L B¯0
¯ n20
kV¯0′ Ω0
¯
n′0+ ¯n′0B¯′0 B¯0
+ ¯n′′0−2¯n′20
¯ n0
Ω0−λ¯kdi L
VAc V∗
¯ n′0
¯ n20B¯0
+ ¯B0B¯′0+δ 4τ β V∗
VAc di
L 1 B¯0
(
¯
p′0V¯0′+ ¯p0V¯0′′−2 ¯B′0 B¯0
¯ p0V¯0′
)]
× [
B¯02 (n¯′0
¯
n0 −B¯′0 B¯0
) + β
2 (1−τ) (n¯′0
¯ n0
γp¯0−p¯′0 )]
+ (
B¯20+ δ 4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0
) { B¯0B¯0′
(2¯n′0
¯
n0 − B¯0′ B¯0
)
−B¯0B¯′′0 +B¯02
¯ n0
(
¯
n′′0−n¯′20
¯ n0
) + β
2 (1−τ) [γp¯0
¯ n0
(
¯
n′′0− ¯n′20
¯ n0
) + ¯n′0
¯ n0
γp¯′0−p¯′′0 ]}}
(A.6)
A2(x;ω, k) = δ 4τ β V∗
VAc
di
L
¯ p0
B¯0
Ω0 (A.7)
B2(x;ω, k) = δ 4τ β V∗
VAc
di
L Ω0
B¯0
(
¯ p′0−B¯′0
B¯0
¯ p0
) + β
2kγp¯0
( 1 + δ
2 V∗
VAc
di
L τV¯0′
B¯0
)
+k (
B¯20− δ 4τ β V∗
VAc di
L
¯ p0V¯0′
B¯0
)
+λ¯k2di
L VAc
V∗ 1
¯ n0B¯0
B¯02− δ 4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0 Ω0−λ¯kdi
L VAc
V∗
¯ n′0
¯ n20B¯0
× [
B¯02 (n¯′0
¯
n0 −B¯0′ B¯0
) +β
2 (1−τ) (¯n′0
¯ n0
γp¯0−p¯′0 )]
(A.8)
C2(x;ω, k) = V∗2
VAc2 ¯n0V¯0′Ω0− δ 4τ β V∗
VAc
di
L
¯ p0
B¯0
¯k2Ω0 +β
2¯k¯p′0 (
1 + δ 2
V∗ VAc
di L
τV¯0′ B¯0
)
+ ¯kB¯′0 B¯0
( B¯02−δ
4τ β V∗ VAc
di L
¯ p0V¯0′
B¯0
)
(A.9)
F2(x;ω, k) = V∗2
VAc2 n¯0Ω20+ δ
4τ β¯k V∗ VAc
di L
Ω0 B¯0
(
¯ p′0− B¯0′
B¯0
¯ p0
)
−β 2¯k2γp¯0
( 1 +δ
2 V∗ VAc
di L
τV¯0′ B¯0
)
−¯k2 (
B¯02−δ 4τ β V∗
VAc
di L
¯ p0V¯0′
B¯0
)
−λ¯k3di
L VAc
V∗
1
¯ n0B¯0
B¯20− δ 4τ β V∗
VAc
di
L
¯ p0V¯0′
B¯0
Ω0−λ¯kdi
L VAc
V∗
¯ n′0
¯ n20B¯0
× [
B¯20 (n¯′0
¯
n0 − B¯0′ B¯0
) +β
2 (1−τ) (n¯′0
¯ n0
γp¯0−p¯′0 )]
(A.10) Ω0=ω−k¯V¯0
ここで以下の無次元化及び上付きの’は物理量のy方向への微分である。
B≡B∗B,¯ n≡n∗n,¯ p≡p∗p,¯ v≡VAcv,¯ x≡L¯x, ω = (V∗/L)ω, k=k/L これらの係数では以下の物理量を用いて表現されている。
β = 2p∗
B2∗/µ0
, VAc= B∗
√µ0n∗mi
, di=
√ mi
µ0n∗e2 又、線形解析では以下で表されるプロファイルを仮定している。
g=gey, (A.11)
g= const. (A.12)
B0=B0(y)ez, (A.13)
n0=n0(y), (A.14)
p0=p0(y), (A.15)
v0=V0(y)ex, (A.16)
pi0=τ p0, (A.17)
pe0= (1−τ)p0. (A.18)