プラズマ乱流解析の方法論を通じて他分野とつながる
前山伸也
核融合科学研究所 メタ階層ダイナミクスユニット

第9回研究部セミナー, 管理・福利棟４階第一会議室, 2023年11月24日15:00 - 16:00

Outline➢ 自己紹介 - 3 min➢ 手法1. 三つ組み相互作用の対称化と有向表現 - 10 min➢ 手法2. 射影演算子法に基づく統計的時系列データ解析 - 10 min➢ 手法3. マルチフィデリティデータ融合モデリング - 10 min➢ まとめと議論 - 7 min
2 / 40自己紹介
名前：前山伸也
所属：核融合科学研究所 メタ階層ダイナミクスユニット・准教授
経歴：生まれも育ちも青森県 → 東工大（工）→ 東工大（原子核）→ 原子力機構（核融合）→ 原子力機構（システム計算）→ 名古屋大（物理）→ 核融合研
専門：プラズマ乱流、ジャイロ運動論、高性能コンピューティング

https://www.p.phys.nagoya-u.ac.jp/gkv/
https://www.github.com/GKV-developers/

3 / 40これまでの研究：プラズマ乱流マルチスケール相互作用問題

従来：スケール分離仮定に基づく単一スケール解析

手法

電子・イオン系マルチスケール乱流シミュレーションを世界に先駆けて実現 ~4×107 core-hour / simulation

成果✓ マルチスケール相互作用が存在し、乱流輸送が大きく変化 [Maeyama, PRL (2015), Maeyama, NF (2017)]✓ 電子スケールが支配的なケースも。多様性と普遍性

イオンスケール
揺らぎ ρi~1cm

電子スケー
ル揺らぎ
ρe~0.1mm[Maeyama, PRL (2017)]✓ 核燃焼プラズマを見据えた高Te・多種プラズマ解析[Maeyama, Nature Commun. (2022)]✓ 日欧米共著によるマルチスケール乱流に関するオーバービュー [Maeyama, IAEA-FEC 2023]

米国のAlcator C-Mod, DIII-D実験もマルチスケール相互作用の重要性を支持 [Howard (2016), Holland (2017)]
欧州でも解析が進められる [Mariani (2021), Mantica (2021)]

4 / 40本日のお話

これまでのマルチスケール乱流研究やそこでの問題意識から開発した、解析手法・モデリング手法についてお話しします。
方法論を軸として、考えたい問題のもつ構造（階層）を俯瞰する（メタ視点）ことで、他分野とも共通点が見つかると期待。
・・・メタ階層的アプローチの一つの形

5 / 40Outline➢ 自己紹介 - 3 min➢ 手法1. 三つ組み相互作用の対称化と有向表現 - 10 min➢ 手法2. 射影演算子法に基づく統計的時系列データ解析 - 10 min➢ 手法3. マルチフィデリティデータ融合モデリング - 10 min➢ まとめと議論 - 7 min
6 / 40研究概要
問題 複雑な乱流非線形相互作用を如何にして抽出するか？
Fourier解析 [Nakata’12PoP]

SVD解析[Sasaki’20PPCF]

EAEBECED

周波数解析[Xu’09PoP]

手法 • 2次の非線形性を介したエネルギー伝達という一般的な問題→ 任意の直交関数展開に対する三つ組み相互作用解析として定式化。
成果 ➢ マルチスケール乱流相互作用の非等方・非局所な集団相互作用を解析[Maeyama, NF (2017)]➢ 特異値分解(SVD)に対するエネルギー伝達解析 [Maeyama, New J. Phys. (2021)]
7 / 40二次の非線形性による三つ組み相互作用• 二次の非線形性のモデル方程式 (Cf. 渦度方程式の移流項)
𝜕𝑓 𝑡, 𝒙𝑑 𝑓2= {𝜙, 𝑓}න 𝑑𝒙 = 0
エネルギー保存則
𝑑𝑡2𝜕𝑡∵ න𝑓{𝜙, 𝑓}𝑑𝒙 = 0• 揺動を直交展開によりモードに分ける (Cf. Fourier, SVD, …)
𝑓 𝑡, 𝒙 = ෍ 𝑓𝑘 𝑡, 𝒙

where න𝑓𝑖 𝑓𝑗 𝑑𝒙 = 0 for 𝑖 ≠ 𝑗𝑘

• 各モードのエネルギーの時間発展は
𝑑𝐸𝑘= 𝑇𝑘𝑑𝑡
ここで𝑝,𝑞𝑇𝑘 = ෍ ෍ 𝐴𝑘𝑝𝑝,𝑞(Total transfer to the mode 𝑘)𝑞𝐴𝑘 = න𝑓𝑘 {𝜙𝑝 , 𝑓𝑞 } 𝑑𝒙

(Transfer via the coupling among 𝑘, 𝑝, 𝑞)
8 / 40

三つ組み伝達関数の対称化
非対称三つ組み伝達関数𝑝,𝑞𝐴𝑘 = න𝑓𝑘 {𝜙𝑝 , 𝑓𝑞 } 𝑑𝒙•𝑝,𝑞𝑝,𝑘𝑝,𝑞𝐴𝑘 = −𝐴𝑞 の性質から、 𝐴𝑘 を「モードpを介して、qからkへのエネルギー伝達と解釈」する提案 [Dar(2001)] 以降、プラズマ分野ではMHD・微視的乱流解析に広く用いられてきた [Mininni(2005),Holland(2007),Cerri(2014),Navarro(2014),Teaca(2017), …]

対称化三つ組み伝達関数
1 𝑝,𝑞1𝑝,𝑞𝑞,𝑝𝑆𝑘 = 𝐴𝑘 + 𝐴𝑘 = න𝑓𝑘 {𝜙𝑝 , 𝑓𝑞 } + {𝜙𝑞 , 𝑓𝑝 } 𝑑𝒙22•••••𝑝,𝑞𝑝,𝑞Energy consistency 𝑇𝑘 = σ𝑝 σ𝑞 𝐴𝑘 = σ𝑝 σ𝑞 𝑆𝑘
𝑝,𝑞𝑞,𝑝Symmetry 𝑆𝑘 = 𝑆𝑘
𝑝,𝑞𝑞,𝑘𝑘,𝑝Detailed balance 𝑆𝑘 + 𝑆𝑝 + 𝑆𝑞 = 0
𝑝,𝑞𝑝,𝑞
正味伝達に寄与しない𝐴𝑘 のp,q反対称成分を除いた。 𝑆𝑘 は「k,p,q三つ組み相互作用によりモードkが受け取る正味エネルギー伝達」
中性流体[Kraichnan(1952),Ohkitani(1992)]、ジャイロ運動論[Nakata(2012),Maeyama(2017)]
9 / 40Importance of symmetrization
An apparent example
ሚ

𝜕𝑓
ሚ
෨
෨ 𝑓ሚ = 𝑐 𝜙,෨ 𝜙෨ = 0.
Consider a no physical interaction case 𝑓 = 𝑐 𝜙, then= 𝜙,𝜕𝑡Nothing happens physically. Interactions among three modes 𝑓ሚ = 𝑓ሚ𝑘 + 𝑓ሚ𝑝 + 𝑓ሚ𝑞 ?

- Non-symmetrized transfer
𝑝,𝑞𝐴𝑘 = න𝑓𝑘 {𝜙𝑝 , 𝑓𝑞 } 𝑑𝒙 = 𝑐 2 න𝜙𝑘 {𝜙𝑝 , 𝜙𝑞 } 𝑑𝒙 ≠ 0Find newinteractions?
No, it’s fiction.
No fictitiousinteractions.

- Symmetrized transfer
1𝑐2𝑝,𝑞𝑆𝑘 = න𝑓𝑘 {𝜙𝑝 , 𝑓𝑞 } − {𝜙𝑝 , 𝑓𝑞 } 𝑑𝒙 = න𝜙𝑘 {𝜙𝑝 , 𝜙𝑞 } − {𝜙𝑝 , 𝜙𝑞 } 𝑑𝒙 = 02210 / 40

Difference among triad transfer functions
非対称伝達関数𝑝,𝑞• 𝐴𝑘 をq→kの1対1伝達として解釈。
• 三つ組み間の循環の任意性があるため、“見かけの伝達”が紛れ込みうる。

対称化伝達関数• 対称化により正味のエネルギー収支を抽出。循環の任意性を除外。
• 詳細つり合いは三つ組み間のエネルギー保存則。
𝑝,𝑞𝑞,𝑘𝑘,𝑝
𝑆𝑘 + 𝑆𝑝 + 𝑆𝑞 = 0

Not recommended.

Recommended!

有向表現（新提案）• 対称化伝達関数に基づき、「1つのモードはgainまたはlossのいずれかのみ」というルールの下で、1対1伝達関数を再構築。
• 循環の任意性を除外しつつ、伝達経路が解析可能。

11 / 40

解析の適用例： Sasaki, et al., PPCF (2020).
Flow patterns in Kelvin-Helmholtz
turbulence (extracted by SVD)

Time evolution of azimuthalenergy spectrum

Time evolution of each energy

EAECEBED

モードABCD間相互作用を三つ組み伝達関数により定量評価。
12 / 40

Consistency with its directional representation➢ 新提案の有向表現と対称化伝達との整合性を確認。
➢ 一方、従来の非対称伝達関数は、他とは異なる。また循環の任意性に起因して、評価式によりグラフが変わり得る。 → 物理過程の誤解釈の危険性。
非対称伝達関数（従来手法）

対称化伝達関数

有向表現（新提案手法）

B:帯状流

A: 背景
シア流

C:KelvinHelmholtz
モード

D: 間欠的ス
パイラル構造13 / 40

Contracted mode-to-mode transfer
有向表現は、対称化伝達では行えない伝達経路や縮約伝達解析などの利点。
𝑝𝐷𝑘←𝑞 = ෍ 𝐷𝑘←𝑞𝑝

非対称伝達とは異なり、見かけの伝達を除外。

14 / 40

オープンソースPythonパッケージ triadgraph
三つ組み相互作用の対称化・有向表現・ネットワーク可視化用Pythonパッケージ triadgraph として公開
https://github.com/smaeyama/triadgraph/t, A_kpq, A_qkp, A_pqk, …

[S. Maeyama, M. Sasaki, et al., New J. Phys. (2022)]

保存的な２次の非線形性による三つ組み相互作用
𝑝,𝑞を解析したい方、非対称伝達 𝐴𝑘 のデータがあれば関数呼び出しで試せます15 / 40

Outline➢ 自己紹介 - 3 min➢ 手法1. 三つ組み相互作用の対称化と有向表現 - 10 min➢ 手法2. 射影演算子法に基づく統計的時系列データ解析 - 10 min➢ 手法3. マルチフィデリティデータ融合モデリング - 10 min➢ まとめと議論 - 7 min
16 / 40

研究概要
問題 渦粘性モデルでは、自発励起する小スケール揺らぎの効果を表せない。
マルチスケール乱流相互作用

駆動
渦粘性モデルによる散逸

エネルギー

エネルギー

既存手法：渦粘性モデル

波数

手法 マルチスケール相互作用を、異なるスケール大スケーの揺らぎ間の統計相関として捉えなおす。
ル揺動
成果u(t)➢ 射影演算子法に基づく統計的時系列データ 小スケー解析を開発。物理解釈とモデリングに応用 ルからの寄与 f(t)[Maeyama, JPSJ (2020); 前山&三分一, 統計数理 (2023)]？

波数
時系列データアンサンブル…
時刻 t

時刻 t

時刻 t
17 / 40

Statistical analysis of time-series data
Set up of time-series data

A numberof samples• A response variable 𝑓 (𝑙) 𝑡 , which has dimensions of
the number of samples 𝑙 and time 𝑡.
• Explanatory variables 𝒈 𝑙 𝑡 can be a vector 𝒈 =𝑔0 , 𝑔1 , ⋯ , 𝑔𝑛−1 .
Response(𝑙)𝑙• We want to extract how 𝑓 𝑡 correlates with 𝒈 𝑡 . variable 𝒇 𝒕

Explanatory
variable 𝒈 𝒕
Time 𝒕

森-Zwanzig の射影演算子法
H. Mori, “Transport, Collective Motion, and Brownian Motion”, Prog. Theor. Phys. 33, 423 (1965).

• 非平衡統計物理で開発された、遅い集団運動とミクロな揺動力を分離する数理手法• 決定論的な運動方程式から確率的な一般化Langevin方程式を得る方法論
18 / 40

Formulation of Mori-Zwanzig projection operator method (1)• Physical quantities 𝒖 𝒕 is governed by 𝑑𝒖/𝑑𝑡 = 𝑵 𝒖(𝑡) , with initial 𝑢𝑖 𝑡 = 0 = 𝑢𝑖 .
• Response and explanatory variables are represented by some functionals of 𝒖 𝒕 , i.e.,𝒇 𝑡 =𝑭 𝒖 𝑡 , 𝒈 𝑡 =𝑮 𝒖 𝑡 .
• Assume statistically steady state (𝑡 = 0 has no special meaning).

We introduce three operators 𝑃, 𝑄, Λ:• Projection operator 𝑃 and its complementary projection 𝑄 = 1 − 𝑃
𝑃𝑃 = 𝑃 (so that 𝑄𝑄 = 𝑄, 𝑃𝑄 = 𝑄𝑃 = 0). Not specified here.
•Liouville operator Λ =𝑑𝒖 𝑡𝑑𝑡⋅

𝜕= σ𝑖 𝑁𝑖
𝜕𝒖𝒖 𝑡

𝜕
𝜕𝑢𝑖

Time evolution of arbitrary functions 𝑓𝑖 𝑡 = 𝐹𝑖 𝒖 𝑡 can be rewritten as
𝑑𝑓𝑖 𝑡= Λ𝑓𝑖 𝑡 ,𝑓𝑖 𝑡 = 𝑒 𝑡Λ 𝑓𝑖𝑑𝑡
Time evolution operator 𝑒 𝑡Λ is also called as the Koopman operator.
19 / 40

Formulation of Mori-Zwanzig projection operator method (2)𝑓𝑖 𝑡 = 𝑒 𝑡Λ 𝑓𝑖= 𝑒 𝑡Λ 𝑃𝑓𝑖 + 𝑒 𝑡Λ 𝑄𝑓𝑖𝑡・・・Divide initial 𝑓𝑖 into projection and complement.

= 𝑒 𝑡Λ 𝑃𝑓𝑖 + න 𝑒 𝑡−𝑠 Λ 𝑃Λ𝑒 𝑠𝑄Λ 𝑄𝑓𝑖 𝑑𝑠 + 𝑒 𝑡𝑄Λ 𝑄𝑓𝑖0𝑡

where the Dyson formulae 𝑒 𝑡Λ = 𝑒 𝑡𝑄Λ + ‫׬‬0 𝑑𝑠𝑒 𝑡−𝑠 Λ 𝑃Λ𝑒 𝑠𝑄Λ is used.
Defining

𝑑𝑟𝑖 𝑡
𝑟𝑖 𝑡 = 𝑒 𝑄𝑓𝑖 ,⟺= 𝑄Λ𝑟𝑖 𝑡 ,𝑑𝑡
one yields the M-Z formalism𝑡𝑄Λ

𝑟𝑖 0 = 𝑄𝑓𝑖𝑡

𝑓𝑖 𝑡 = 𝑒 𝑡Λ 𝑃𝑓𝑖 + න 𝑒 𝑡−𝑠 Λ 𝑃Λ𝑟𝑖 𝑠 𝑑𝑠 + 𝑟𝑖 𝑡0(Term from projection of initial 𝑓𝑖 )+(Memory term)+(Uncorrelated term)
Projection of the uncorrelated term is always zero, 𝑃𝑟𝑖 𝑡 = 0 for any 𝑡 ≥ 0.
20 / 40射影演算子Pの選択
ここまでは射影と補射影の分解と記憶の獲得という一般的な形式論。
実用的な計算を進めるには、射影演算子Pを具体的に定義する必要がある。
ここでは、目的変数 f をある説明変数 g に射影するクラスの演算子を考えることにする。
■ 最適推定を与える射影演算子は、条件付き期待値関数Eである。
[Zwanzig(1980),Chorin(2002)](左) f,gのデータサンプルの散布図と条件付き期待値関数E(実線)、および、(右) f,gの同時確率密度分布

𝒫 nonlinear 𝑓𝑖 = E 𝑓𝑖 |𝒈
理想的には最適な射影だが、有限のデータセットから構築するのは困難。
■ 森の線形射影[Mori(1965)]
𝒫 linear 𝑓𝑖 = 𝑓𝑖 𝒈∗ ⋅ 𝒈𝒈∗ −1 ⋅ 𝒈は、条件付期待値に対する線形回帰と解釈できる。
21 / 40

Formulation of Mori-Zwanzig projection operator method (3)Now we employ the Mori’s linear projection to the initial values of the explanatory
variables 𝒈 [Mori’65Prog. Theor. Phys.],
𝑃𝑓𝑖 𝑡 = 𝑓𝑖 𝑡 𝒈∗ ⋅ 𝒈𝒈∗ −1 ⋅ 𝒈 = ෍ ෍ 𝑓𝑖 𝑡 𝑔𝑗∗ ⋅ 𝒈𝒈∗ −1 𝑗𝑘 ⋅ 𝑢𝑘𝑗𝑘

where 〈⋯ 〉 is the ensemble average. Then, 𝑃𝑟𝑖 𝑡 = 0 means 𝑟𝑖 𝑡 𝒈∗ = 0.
Calculate the terms in the M-Z formalism.

𝑒 𝑡Λ 𝑃𝑓𝑖 = 𝑒 𝑡Λ 𝑓𝑖 𝒈∗ ⋅ 𝒈𝒈∗ −1 ⋅ 𝒈 = 𝑓𝑖 𝒈∗ ⋅ 𝒈𝒈∗ −1 ⋅ 𝒈 𝑡 ≡ ෍ Ω𝑖𝑗 𝑔𝑗 𝑡𝑗𝑡𝑡

න 𝑒 𝑡−𝑠 Λ 𝑃Λ𝑟𝑖 𝑠 𝑑𝑠 = න 𝑒 𝑡−𝑠 Λ Λ𝑟𝑖 𝑠 𝒈∗ ⋅ 𝒈𝒈∗ −1 ⋅ 𝒖𝑑𝑠
00

Markovcorrelation Ω𝑑𝒈∗= − න 𝑟𝑖 𝑠⋅ 𝒈𝒈∗ −1 ⋅ 𝒈 𝑡 − 𝑠 𝑑𝑠
𝑑𝑡0𝑡𝑡

≡ − න ෍ Γ𝑖𝑗 (𝑠) 𝑔𝑗 𝑡 − 𝑠 𝑑𝑠
0𝑗

Memory function Γ22 / 40

Formulation of Mori-Zwanzig projection operator method (4)Finally, the (extended) generalized Langevin form of 𝒇(𝑡) is𝑡𝒇(𝑡) = Ω ⋅ 𝒈 𝑡 − න Γ 𝑠 ⋅ 𝒈(𝑡 − 𝑠)𝑑𝑠 + 𝒓 𝑡0

Markov coefficient: Ω = 𝒇𝒈∗ ⋅ 𝒈𝒈∗ −1 ,𝑑𝒈∗
Memory function: Γ 𝑡 = 𝒓 𝑡⋅ 𝒈𝒈∗ −1 ,𝑑𝑡
Uncorrelated term: 𝒓 𝑡 = 𝑒 𝑡𝑄Λ 𝑄𝒇which has the following properties:𝒓 𝑡𝑑𝒈∗
𝑑𝑡

= Γ 𝑡 ⋅ 𝒈𝒈∗

𝒓 𝑡 𝒈∗ = 0・・・ (Extended) generalized fluctuation-dissipation
theorem of the 2nd kind・・・ Uncorrelated term 𝒓 𝑡 has no corr. with 𝒈.

† This is NOT a theoretical approximation. MZ projection operator method splits a time-series data 𝒇 𝑡 into correlatedterms and an uncorrelated term 𝒓 𝑡 with regard to the explanatory variables 𝒈 𝑡 .

23 / 40

Computation of Ω and Γ 𝑡 from data set•Markov correlation ΩThe ensemble average is calculated by average over samples,𝑁sample −1 𝑖1
∗𝑖 ∗∗
∗ −1
σ𝒇𝒈 =
𝒇0
𝒈0,
andthenΩ=𝒇𝒈⋅𝒈𝒈
𝑖=0𝑁sample•Memory function Γ 𝑡𝑑𝒈∗
Multiply⋅
𝑑𝑡𝒈𝒈∗ −1 on generalized Langevin eq. and taking ensemble average,
one yields the memory equation𝑡Γ 𝑡 = 𝐸 𝑡 + න Γ 𝑠 ⋅ 𝐻 𝑡 − 𝑠 𝑑𝑠0where 𝐸 𝑡 = 𝒇 𝑡 − Ω ⋅ 𝒈 𝑡𝑑𝒈∗𝑑𝑡

⋅ 𝒈𝒈

∗ −1

and 𝐻 𝑡 − 𝑠 = 𝒈 𝑡 − 𝑠𝑑𝒈∗
𝑑𝑡

⋅ 𝒈𝒈∗ −1 areevaluated from data. One can construct Γ(𝑡) successively from Γ 0 = 𝐸 0 .
•Uncorrelated term 𝑟(𝑡)𝑡𝒓 𝑡 = 𝒇 𝑡 − Ω ⋅ 𝒈 𝑡 + න Γ 𝑠 ⋅ 𝒈(𝑡 − 𝑠)𝑑𝑠
024 / 40

Example: Kuramoto-Sivashinsky turbulence [Maeyama’20JPSJ]
𝜕𝑢(𝑥, 𝑡)𝜕𝑢 𝜕 2 𝑢 𝜕 4 𝑢+𝑢++=0𝜕𝑡
𝜕𝑥 𝜕𝑥 2 𝜕𝑥 4

A simple 1D turbulence modelwith small-scale instability

Procedure for the analysis
1. From K-S turbulence simulation data, evaluate time evolution of low-k fluctuation 𝒈 𝑡 = 𝑢ത 𝑘 andhigh-k contribution term 𝒇 𝑡 = − 𝑢𝜕𝑥 𝑢 − 𝑢ത 𝜕𝑥 𝑢ത 𝑘 .
𝑡
2. Applying M-Z projection of high-k contribution term to low-k, 𝑓𝑘 (𝑡) = Ω𝑘 𝑔𝑘 𝑡 − ‫׬‬0 Γ𝑘 𝑠 𝑔𝑘 𝑡 − 𝑠 𝑑𝑠 +𝑟𝑘 𝑡 , extract physical information from evaluated Ω𝑘 , Γ𝑘 𝑠 and 𝑟𝑘 𝑡 .
3. If possible, model high-k contribution term and test reproducibility.
25 / 40

Example: Kuramoto-Sivashinsky turbulence [Maeyama’20JPSJ]
Results✓ The memory term of high-k contribution has
damping effects on low-k in K-S turbulence.

Memory function 𝚪 𝒕✓ Low-resolution simulation with the M-Zbased sub-grid scale (SGS） model
qualitatively reproduces DNS result.
The uncorrelated term is crucial forstochastic forcing of low-k fluctuation.

Evaluated coefficients ofcorrelated part 𝛀𝒌 , −𝜸𝒌𝑑𝑢Since 𝑑𝑡𝑘 ≃ Ω𝑘 − 𝛾𝑘 𝑢𝑘 + 𝑟𝑘 ,Re Ω𝑘 − 𝛾𝑘 < 0 means damping.

Test of reproducibility by M-Z-based modeling k-ω space energy spectrum E(k,ω)
Low-res.
+ SGS model

26 / 40

オープンソースPythonパッケージ mzprojection
1. Download mzprojection.py from GitHub:
https://github.com/smaeyama/mzprojection† pip install mzprojection or conda install –c smaeyama mzprojection may also work.

2. Prepare data set.
f[nsample,nperiod,nf] # Response variablesg[nsample,nperiod,ng] # Explanatory variables

3. Import mzprojection module and call function.
from mzprojection import mzprojection_multivariate_discrete_time
omega, memoryf = mzprojection_multivariate_discrete_time(g, f)

4. That’s all! You get,Markov coefficient Ω:Memory function Γ 𝑡 :omega[nf,ng]
memoryf[nperiod,nf,ng]

See an example: https://github.com/smaeyama/mzprojection/blob/main/tests/Demo_Oct2021_multivariate.ipynb

27 / 40

Outline➢ 自己紹介 - 3 min➢ 手法1. 三つ組み相互作用の対称化と有向表現 - 10 min➢ 手法2. 射影演算子法に基づく統計的時系列データ解析 - 10 min➢ 手法3. マルチフィデリティデータ融合モデリング - 10 min➢ まとめと議論 - 7 min
28 / 40

研究概要
問題 核融合炉開発では、乱流について粗視化した統合輸送コードが用いられ、その信頼性は輸送モデルに依存。
（局所プラズマパラメータ→輸送フラックスの関係が知りたい）

しかし、マルチスケール乱流シミュレーションは大規模すぎてパラメータ調査できず、ニューラルネットワークなどの代替モデル構築も難しい。
手法 着想：マルチフィデリティデータ融合手法の乱流輸送モデリングへの適用• Low-fidelity data: 低精度・多数のデータ• High-fidelity data: 高精度・少数の限られたデータ→ Low-high fidelityデータ間の相関を取り入れることで予測精度を改善する。
29 / 40非線形自己回帰Gauss過程回帰 (NARGP) [Perdikaris (2017)]
入力𝑥𝑙 と出力𝑦𝑙 のデータセッ
トが、異なるフィデリティ(𝑙 = 0,1, ⋯ )毎にあるとする。
𝑦𝑙 = 𝑓𝑙 𝑥𝑙𝑥0[[𝑦0,𝑥1], [𝑦1,𝑥high], …, …, [𝑦high,]]➢ 最低次フィデリティ(l=0)・・・平均0,共分散カーネル関数𝑘0 のガウス過程として表す。
𝑓0 ~GP 𝑓0 |0, 𝑘0 𝑥, 𝑥 ′ ; 𝜃0➢ 高次フィデリティ(l>0)・・・l-1次の事後分布𝑓∗𝑙−1 を取り入れたガウス過程として表す。
𝑓𝑙 𝑥 = 𝑔𝑙 𝑥, 𝑓∗𝑙−1 𝑥𝑔𝑙 ~GP 𝑓𝑙 |0, 𝑘𝑙𝑔 (𝑥, 𝑓∗𝑙−1 𝑥 ), (𝑥 ′ , 𝑓∗𝑙−1 𝑥′ ; 𝜃𝑘𝑙𝑔 = 𝑘𝑙𝜌 𝑥, 𝑥 ′ ; 𝜃𝑙𝜌 𝑘𝑙𝑓 𝑓∗𝑙−1 𝑥 , 𝑓∗𝑙−1 (𝑥 ′ ); 𝜃𝑙𝑓 + 𝑘𝑙𝛿 𝑥, 𝑥 ′ ; 𝜃𝑙𝛿

NARGPは、多段のガウス過程回帰により Low-to-high fidelity データ間の非線形相関と入力xに対する非一様性を表現しつつ、計算コストは通常のガウス過程回帰と同オーダー。
30 / 40

ケース1: 1次元関数➢ Low-fidelity データ
正弦波, データ点数50点➢ High-fidelity データ
Low-fidelity関数 𝑓𝑙 を非一様(x依存)かつ非線形(2乗)に変換。データ点数15点（low-fidelityの50か所の評価位置xの内から15か所）→目標は、これらのデータからなるべく正確にHigh-fidelity関数（右図青線）を推定したい。
31 / 40

ケース1: 1次元関数✓ NARGPでは低フィデリティデータと高フィデリティデータ間の相関を学習することで、通常のGauss過程回帰に比べて少ないデータ点数でも良い予測を構築。
✓ 加えて、分散により予測精度の指標を与える。
Single GP (通常のGauss過程回帰)

NARGPによるマルチフィデリティ回帰

32 / 40

ケース2: 2次元問題(Branin関数) 結果のみ
入力xが2次元かつ、フィデリティも3段階としたテスト問題。
✓ 多次元・多段についてもNARGPの実装が正しく動作していることを確認。
✓ この例についても、通常のガウス過程回帰よりもNARGPの優位性を確認。

33 / 40

ケース3: プラズマ乱流輸送データ
E. Narita, et al., Contrib. Plasma Phys., e202200152 (2023)のJT-60Uデータ• 入力とする局所プラズマパラメータ(12個):'R/Ln', 'R/LTe', 'R/LTi', 'ni/ne', 'Te/Ti', 'beta', 'nu_ee', 'q', 's', 'epsilon', 'kappa', 'delta’• 低フィデリティデータ（最大9749点）: 半経験的準線形拡散係数 DQL• 高フィデリティデータ（最大170点）: 実験値に基づく拡散係数 Dexp
→ 半分は学習用データとして利用し、残り半分は検証用データとして未参照で残す

34 / 40

ケース3: プラズマ乱流輸送データ✓ Single GPに比べて、NARGPにより予測精度改善。(1)プラズマパラメータと準線形拡散係数（線形成長率や揺動波数を反映）の関係を広いパラメータで学習し、(2)少数のTrainデータの範囲ながら準線形拡散係数と実験値との相関を学習しておくことで、未知のTestデータに対しても予測性が向上。

Single GP

NARGP

MSLE=0.089
R2=0.84

MSLE=0.036
R2=0.94(参考)準線形モデル単体による予測精度

35 / 40

Outline➢ 自己紹介 - 3 min➢ 手法1. 三つ組み相互作用の対称化と有向表現 - 10 min➢ 手法2. 射影演算子法に基づく統計的時系列データ解析 - 10 min➢ 手法3. マルチフィデリティデータ融合モデリング - 10 min➢ まとめと議論 - 7 min
36 / 40

まとめと議論手法1. 三つ組み相互作用の対称化と有向表現• 二次の非線形性と揺らぎの直交関数展開 → 三つ組み相互作用• 対称化の重要性。有向表現を提案。ネットワーク可視化。
他の研究とのつながり• プラズマに限らず三つ組み相互作用の有向表現は適用可能•中性流体でも1対1伝達解析の提案[Obuse, JPSJ (2020)]があるが、比較・議論を進めて手法を標準化したい[前山,流体力学会予稿(2021)]•

回転乱流における局所伝達に基づく解析[Yokoyama, JFM(2021)]に対し、本手法では非等方・非局所伝達も表現可能。

•マルチフィールドSVD[Yatomi, PPCF (2023)]に対しても、対称化伝達関数の部分和をとるよりも、ネットワーク可視化や有向表現に基づいて解析することをお勧め37 / 40

まとめと議論手法2. 射影演算子法に基づく統計的時系列データ解析• 時系列データ間の相関を記憶と無相関揺らぎに分離するデータ展開手法• 解析：一般化Langevin描像としての物理的解釈を与える• モデリング：データに基づく統計的モデリングにも応用できる
他の研究とのつながり• 統計数理分野で用いる自己回帰分析との類似と相違 [前山&三分一, 統計数理(2023)]•

射影演算子法の非線形拡張[Lin (2023)]。一般化第二種揺動散逸定理の定式化や連続時間系アルゴリズムなど、情報科学的な理論拡張の余地。

•

射影演算子法は分子動力学計算のモデリングにも用いられる[Kasahara, J. Chem.
Phys. (2016)]。さきがけ“複雑流動”の枠組みで、理論↔データ解析の連携研究。

•プラズマ乱流のみならず、中性流体乱流へも適用性。
38 / 40

まとめと議論手法3. マルチフィデリティデータ融合モデリング• 低精度だが多数の低フィデリティデータと高精度だが少数の高フィデリティ
データの相関を利用することで、全体としての予測精度の向上を図る方法論。
他の研究とのつながり• 輸送モデリングを行っている研究者らとデータ提供を通じて議論中。
•マルチフィデリティデータ融合は、低精度の予測を高精度の実測で較正する手法の一般化。実験計測の較正などのデータにも適用してみたい。

•

将来装置の性能を理論・シミュレーションで予測（低フィデリティ）する際に、既存実験（高フィデリティだが将来装置パラメータの情報は持たない）の情報を活用することで、将来性能予測を高精度化し原型炉開発に貢献。

39 / 40

全体のまとめ
これまでのマルチスケール乱流研究やそこでの問題意識から開発した、解析手法・モデリング手法について説明した。
手法1. 三つ組み相互作用の対称化と有向表現 [Maeyama, New J Phys. (2021)]
手法2. 射影演算子法に基づく統計的時系列データ解析 [前山&三分一, 統計数理(2023)]
手法3. マルチフィデリティデータ融合モデリング
方法論を軸として、考えたい問題のもつ（階層）構造を（メタ視点で）俯瞰することで、他分野とも共通した議論ができると期待。
オープンソースコードとして整備し公開 https://github.com/smaeyama/
共同研究の提案は大歓迎。
40 / 40