J. Plasma Fusion Res. Vol.97, No.2 (2021)７
９‐８５

小特集 磁場閉じ込め核融合プラズマにおけるデータ駆動アプローチによる物理モデリングの進展

４．データ駆動アプローチを用いた動的乱流現象の解析
4. Data-Driven Analyses of Dynamical Turbulence Phenomena
佐 々 木 真１，２），河 原 吉 伸３，４），草 場
１，
２）

SASAKI Makoto

彰１）

３，
４）

, KAWAHARA Yoshinobu

and KUSABA Akira１）

１）

九州大学応用力学研究所，２）九州大学極限プラズマ連携研究センター，

３）

九州大学マス・フォア・インダストリ研究所，４）理化学研究所 革新知能統合研究センター（原稿受付：２
０
２
０年１
０月１
４日）

本章ではデータ駆動科学的手法として近年注目されている動的モード分解について，その理論的背景からプ
ラズマ乱流への応用例まで紹介する．動的モード分解は流体現象の時空間構造を観測／計測データから抽出するための数値解析手法の一つとして提案されたものであり，現在では様々なシステムへの応用がなされている．動的モード分解を用いると，時空間発展データから統計的に重要な空間構造とその構造に対応する複素周波数を得ることができ，より少数の自由度でシステムを理解することが可能となる．ここでは，（１）物理機構の理解をめざす例として，動的モード分解をリミットサイクル振動を伴う乱流に適用し，構造間の因果関係を抽出する試みについて紹介する．また，（２）乱流の時空間発展の予測に応用する試みについても紹介する．
Keywords:
plasma turbulence, abrupt phenomena, mode decomposition, DMD, causality, prediction

４．
１ 背景

観測データに基づき，システムの時空間発展を特徴付ける

磁場閉じ込めプラズマや宇宙のプラズマ現象はエネル

ことができる新たな方法である［５］．システムの時間発展

ギーの流入出がある非平衡開放系であり，自律的なリズムを考慮したモード分解であるため，系の動的挙動の理解やや自己組織的な空間パターンの形成が起こるため，システ

予測に適している．近年，流体力学分野をはじめとして，

ムが本質的に動的である．特に，極限的に空間不均一な系

様々な分野における時空間ダイナミクスに適用され，そのである磁場閉じ込めプラズマでは，プラズマ装置に対し有用性が認識されてきている．プラズマ研究への適用も始て，ミクロスケールの乱流が発達し，その非線形ダイナミまっており，例えば［６‐１０］などがある．ここでは，動的クスによってプラズマ閉じ込めが規定されている．ミクロ

モード分解の基礎を概観し，プラズマ乱流研究への応用例

スケール乱流が駆動する装置スケールのマクロ構造や中間を紹介し，物理的理解に向けた試みと定量的予測の試みにスケールのメゾスケール構造の重要性が認識されているついて紹介する．

［１］．そのため，プラズマの動的挙動を理解するためには，

４．
２ 動的モード分解

乱流とその非線形性によって駆動される構造間の相互作用を理解することが必要である．

動的モード分解は，当初，流体現象の時空間構造を観

磁場閉じ込めプラズマにおける重要な動的挙動には周辺

測／計測データから抽出するための数値解析手法の一つと局在モード（edge localized mode，ELM）［２］や閉じ込め遷して提案されたものである［１１，１２］．しかしその後，応用数移［３］等の自発的なものから，燃料供給時の外部摂動に対学分野などで議論される力学系のクープマン解析との関連する応答［４］などがある．これらの挙動の物理機構の理解も知られるようになり，流体分野に限らず広く科学・工学

やその予測は重要な問題である．これらの物理機構の解明

分野へ応用されるようになっている．

をめざした研究が理論・シミュレーション・実験の側面か流体分野における類似したモード分解としては，固有直ら精力的に進められており，理解が大きく進んできている交分解（POD）も古くからよく用いられている［１３］．固有

ものの，これらの現象の定量的な予測や制御には至ってい直交分解では，空間方向を一つの列としてまとめた行列ないのが現状である．近年，データ駆動科学的手法が急速（つまり行方向は時間に対応）に対して特異値分解を計算に発達してきており，この新たな手法を第４の科学的手法することで，時間発展に寄与する特徴的な空間構造の抽出として利用する方向性の研究も必要となっている．

を行う．統計的にはこの計算は，時間方向の相関を無視す
ればいわゆる主成分分析と同様であり，変動の大きな空間

本章で焦点を当てる動的モード分解は，データ駆動的にResearch Institute for Applied Mechanics Kyushu University, Kasuga, FUKUOKA 816-8580,Japancorresponding author’s e-mail: sasaki@riam.kyushu-u.ac.jp

７９!2021 The Japan Society of Plasma
Science and Nuclear Fusion Research

Journal of Plasma and Fusion Research Vol.97, No.2 February 2021

構造，つまり物理的な解釈としてはエネルギーが大きい空間構造を取り出していることになる．そのように，固有直交分解は時間方向の周期性を陽に抽出している訳ではな
く，多くの場合に特徴的な周期性を持つ空間構造の抽出が困難であることが知られていた．しかし後に見るように，
動的モード分解は，時間方向の関係を陽に用いたモード分解を計算する方法であり対応する周波数（とその減衰・発展）も同時に得られる．

図１

クープマン作用素を用いた力学系の解析．

ここでは，まず動的モード分解の理論的背景について概説し，そしてその代表的な具体的手順を示しつつ，動的#'%%&'"*"'%&"*%+&%&-%&"*"" #&#"!$!"% &"!

モード分解の解釈について説明する．
４．
２．
１ クープマン解析（１）'個

ここでは動的モード分解の背景として，力学系のクープマン解析について概説する．動的モード分解は，クープマ

ただし，#&（$$）は固有関数 %&に対応する固有値を表す．

ン解析へのデータを用いた数値解析手法の一つとしても位

式（１）中の各コンポーネント中で時間発展に寄与する部分置付けられる．

は#&のみであるため，観測量の時間発展が#&で決まる固有

まず，状態空間"（一般には実ベクトル空間の部分集合）の周期性（#&の位相）と減衰率（#&の絶対値）を持つダイ$における（決定的な）力学系 ,（,$"）を考え
*
,&
*
&
'
"$
&
'
'

ナミクスに分解されていることがわかる．

る１．$: "→"は，系の軌跡を定める写像である．ここではクープマン作用素のスペクトル分解を用いた非線形力学

データを用いた解析を考えるため，次の離散時間での表現系の解析に関する応用として重要なものの一つに，縮約理も適宜用いる．

論との関係があげられる．特に，リミットサイクル振動子などの周期的な運動に近い挙動を示す系で用いられる位相

&
'!!'
"*

,
,
"
",
'!!"&
''
'!#
"*&

'
"*

$
,&
$
#$
&
'
'

縮約は，クープマン作用素の固有関数と直接的な関係を持つ［１７］．またもう一つの重要な応用としては，力学系から

ただし，,
'
"*
&
'であり，&"*: "→"は対応するフロー' ",得られる観測が流体現象のように空間的に広がりを持つ場である．これらの表現は状態の軌跡を記述することでダイ

合には，分解（１）は，時空間的なコヒーレンス構造の抽出とナミクスを表している．

して有用であることが知られている．つまり，式（１）の各一方で，ダイナミクスのもう一つの表現として観測量のコンポーネント中の固有値以外の部分，つまり +&%&にあた

時間発展を表す方法として，クープマン作用素を用いた記る部分は，各固有値で表される時間的変動に関する %への述が知られている［１５］．状態空間"上に定義されたなんら寄与として解釈できる．各観測量が各空間上の点におけるかの観測空間!内の関数 ': "→%（$!）（観測量）を考えなんらかの量（例えば，流体現象であれば流体速度や渦度たとき，クープマン作用素は次式により定義される［１６］など）に対応している場合には，この寄与の相対的な違い（図１も参照）．

を見ることで，各時間変動に対応する空間的な同調性を知ることができる．流体現象の解析においては，このような#%"%%&"*+#%&
,
,
"%&
''
'!!'

空間的なコヒーレンス構造（組織構造）を知ることは重要

ただし %は関数の合成である．定義から#の線形性は容易な課題として従来から盛んに研究されており，これが動的に示される．今，作用素#が点スペクトルのみを持ち固有

モード分解への注目のきっかけにもなっている．

値分解可能であるとし，さらに観測量 'は作用素#の固有

４．
２．
２ 基本的な手順とその意味

関 数 %&: "→$を 用 い て 一 次 展 開 で き る，つ ま り

ここでは，一般によく用いられる動的モード分解（厳密

'"! #
&%
&の よ う に 表 す こ と が で き る と 仮 定 す る
&"!+

DMD
［１８］）を取り上げ，その手 順 を 説 明 し つ つ，動 的（+&$$は展開係数）．このとき，観測量 'の時間発展は次モード分解の性質について述べる．
動的モード分解は，未知の系に関するクープマン分解式のように展開することができる（クープマン・モード分（１）を，その系から得られたデータから近似的に推定する解とも呼ばれる）．

手順である．この際に，データは，状態変数を入力とする( 個の観測関数 '&: "→$（&"!
!
(）を通して得られて,!
１ 決定的でない場合は，例えば，%$"#（"）
を観測量としてとれば，それは確率的な力学系を扱うことに相当して，任意の$,に対して次のようにクープマン作用素を定義することができる
［１
４］
．
"
'
"#)"*&
#,
(
,
,
%&
,
%&
#'!!'
#' ",
,
#%&
''
'!!(
''
''
'!!"!)'*
ただし，)"*&"*における状態 ,,,'!!'"*における状' が与えられた際の時刻 &'!!(''は遷移状態確率である．つまり，ある時刻 '態の期待値を与える．

８０

Special Topic Article

4. Data-Driven Analyses of Dynamical Turbulence Phenomena

M. Sasaki et al.

いると仮定する．そのように，まずデータとしては，解析

コヒーレンス構造の情報を与えている．流体現象における，##!'したい未知の系 $"/に関して 0
0
0(として
##$
"(
"/'

適用では，こういった情報は，データから空間的な組織構生成された時系列データが与えられたとする造を抽出する目的などに利用できる［１３，１５，２４］．この動的T）．例えば，この系から生成したと想定さ（!: #*%"!
%, +,!

モードに関するコヒーレンスという視点は応用的にも重要

#
#
#
,!
れる有限長%の時系列データ #!!
"!
#!
% が与えられるであり，流体のみでなく，脳科学や感染症伝播，経済学な#場合には，（#（/#"!
%）のようにペアを構成し
,!
/）
!"!
/

ど，多くの分野でこの特徴を用いた解析が研究されているて用いる．このとき，厳密 DMD の手順は次のようになる．

［２５‐２８］．また，このコヒーレンス構造をデータの特徴として予測（教師あり／教師なし学習）に用いるための枠組と(#:#*
１．２つの行列(":#*
#
#
#
#
#
#
,!
,!
!!
"!
"!
#!
%!"+
%+

みも提案されており［２９，３０］，いくつかの応用場面においを構成する．

ても有用であることが確認されている［３１，３２］．

２．("の（打ち切り）特異値分解 ("%&$'H を計算する（!H はエルミート転置，採用ランクを .とする）．
/: #&H(#'$!"の固有値 #.
３．行列 !
+と固有ベクト ル

４．
３ プラズマ乱流の時空間ダイナミクス解析への応用/+を計算する（+#"
"
!
!
.）．
#
,!

４．
３．
１ 乱流構造間の因果関係同定への応用.
.
/+（動的モード）と対応する固有値 #
４．"
$!""
#'
+#(
+

乱流と背景場との動的な相互作用の理解に向けた試みとを得る．

して，動的モードを利用した構造間因果関係抽出方法の紹

この手順からわかるように，動的モード分解は比較的容易に実介をする．動的な相互作用の例として，大域的乱流シミュ装が可能である．なおこの手順は，)(*!&" )(#!!(")".
の最小二乗解 ! の固有値と固有ベクトルを計算している

レーションで得られたリミットサイクル振動に動的モード,",

分解を適用した成果について述べる［６］．
直線磁化プラズマにおける長谷川・若谷方程式の大域的ことに相当するが，主要な固有値のみが得られれば十分で.の固有値分解を計算する代わりに，これあるので，直接 !

乱流シミュレーションで得られた時系列データを解析に使/#&H !& の固有値分解を計
を & を用いて射影した行列 !

用した．本ミュレーションでは，密度・静電ポテンシャ

算している．なお，& は直交行列であるので，これらの固ル・平行速度場の空間３次元的分布の時間発展が出力さ有値分解で得られる最初の .個の固有値は同じものであれ，乱流データは１ケースおおよそ数十ギガバイト程度とる．一般に，DMD は高次元なデータに対して適用する必.も高次元になり固有値計算な
要があり２，そのため行列 !

なる．基礎プラズマ実験でしばしば観測されてる強い周方

どが困難になるため，上記のような計算がよく用いられる

流れの空間不均一性によって駆動されるケルビンヘルムホ

訳である．

ルツ不安定性の非線形発展を得た．渦度ソースの強度を増向流れ［３３］を渦度ソースを導入することで模擬し，周方向そして計算された動的モードとその固有値を用いて，例

加させると，定常的な非線形飽和状態からリミットサイク

えば時系列データの場合には，各サンプルは次のように分ル振動を伴うダイナミックな状態へ変化した［３４］
．ここで

解される．

は，リミットサイクル振動を伴うケースに動的モード分解.

を適用した結果に焦点を当てる．図２に乱流シミュレー

!"
/
.
.
#
+ 1
/%! #
+"
+

（２）

ションで得られた非線形飽和状態でのフーリエモード+#",!
-(の時間発展（,，-はそれぞれ周方向・軸方向モー
'

.
.: #*
T は # #'.
.
.
）
ただし，$:#*
（'
"
"
"
1
1
1
,!
,!
$
.+
.+
"!
#!
"!
#!
!

ド数）を示す．ここで背景場は'!
!
!(に対応している．軸方を満たす定数である．, 個の観測関数 )
（*
!
,）が，
,!
* #"

向に均一で周方向に有限の波数を持つ揺動が発達してい４．
２．
１節で述べた関数空間!内の要素であると仮定すると，る．これはケルビンヘルムホルツ揺動の特徴である．揺動

この分解（２）は，分解（１）の有限データによる近似を与えてエネルギーが大きく変動し，それと同期して背景場も変動いることが確認できる．実際，ある数学的仮定の下で，上している様子がわかる．また，図２下段で示しているス

記の動的モード分解により推定される固有値が，内在するナップショットを見ると，揺動の空間パターンも時間的に真のクープマン作用素のスペクトルの情報を捉えることも変動していることがわかる．リミットサイクルの周期は示されている［２３］．

%%$& $"!!程度であり，これは乱流の周期 %,-+'$"!に比

なお，動的モード分解で得られる分解（２）の各コンポー.
ネントは，推定した固有値 #+で決まる固有の周期性と減衰べ格段に長い（時間はイオンジャイロ周波数で規格化）．

率を持つ要素に分解されていることは，クープマン・モー

る．時刻 /における観測ベクトル #/を，そ の 時 刻 に お

ド分解（１）と同様である．そして，分解（２）における動的

け る 各 空 間 点 に お け る 密 度 #/'!
+(か ら 作 成 し，

このような動的乱流状態に動的モード分解を適用す

モード "+が，これらのダイナミクスがどのように各観測量

#
!
!
!
#/'
!
#/'
!
#/'
-!
+を用意する．そして，時
,(
/#*
"(
#(

%*へ相対的な重みで寄与するか，つまりある種の空間的な空間情報マトリックス ("#*
#
#
#
-!
!!
"!
%!"+と，観測時刻

２ 厳密 DMD において用いるデータの次元は，本質的には観測量 %*の数に相当する．分解（２）（または分解（１））では，%*の張る空間が十分大きいことを仮定しており（例えば［１８］），低次元なデータに対して適用する場合には，数理的に空間を広げる必要がある．具体的な手法もいくつか提案されており，例えば，基底関数を事前に準備する方法［１９］や，再生カーネルを用いた方法［２０］，ニューラルネットをデータから学習して用いる方法
［２
１，
２
２］
などが提案されている．

８１

Journal of Plasma and Fusion Research Vol.97, No.2 February 2021

図２

上段：乱流シミュレーションで得られた非線形飽和状態でのフーリエモード（m, n）の時間発展（m, n はそれぞれ周方向・軸方向
モード数）．背景流れ場は（0, 0）に対応している．下段：t ＝ 3050，3150，3200 における密度揺動のスナップショット［６］．

ルツ揺動，Mode4，5 はスパイラル構造に対応している．上述のように，解析している乱流状態は，プラズマの背景分
布，乱流揺動が乱流スケールに比べ長い時間でリミットサ
イクル振動しており，乱流の空間パターンも大きく変化する．この場合，各 DMD モードの時間発展を決める複素周
波数は指数関数的増加・減少しか記述することができない
ため，各モードの増減を表現するには工夫が必要である．
そこで，畳み込み積分を用いた固有関数の強度を評価する方法を提案する．この方法では，各時間ステップにおいて，得られた固有関数とオリジナルの乱流データとの間で畳み込み積分を行い，その空間最大値を固有関数（特徴的構造）の 強度 とする．すなわち，Mode-&の 強度 の瞬時
値 !&#($は，次のように計算される．
図３

支配的固有関数の空間パターン（Mode1~5）
［６］
．
!!

!
!'$
!!!'$
!'$
##(#
"&#
!
(
"&#
'!'" ##(#
'!'，
$
"!
（３）

を１ステップずらした %""%
"
"
"
(!
!!
"!
$ &とを準備し，前節で述べた特異値分解による低次元化した動的モード分解

!
!&#
(
"&#
!
(
)%
$
"(&
$
&
.

を行なった．ここで低次元化に用いたランクは '"%とし!

（４）

た．こうして，系の特徴量である固有値 #&とそれに対応す得られた !&#($を図４に示す．!&はうまく DMD モードのる固有関数 "&#!$を得ることができる．以下では，リミッ

増減を表現しており，プラズマ背景分布の変形とスパイラ

トサイクル振動とリミットサイクル周期に比べ高周波の乱ル構造の出現・消失まで含めた乱流パターンの変化を明確

流の挙動について，動的モード分解より得られた固有関数に捉えることができている．
DMD モードの 強度 !&を用いて，プラズマ分布や乱流（特徴的空間構造）間の因果関係を抽出する試みを紹介する．

構造間の因果関係を評価することができる．一例として，

図３に得られた固有関数（寄与率の大きな上位５つ）を

コ ヒ ー レ ン ト 揺 動（Mode3）と ス パ イ ラ ル 構 造

示す．ここでは，寄与率の最も大きなモードから順番に（Mode4）について，図４下段に示すように !#，!$の時間

Mode1 から Mode5 と呼ぶことにする．Mode１は背景密度発展からリサージュを書くことができる．コヒーレント揺分布，Mode2 は径方向に有限の波数を持ち周方向に一様な動の振幅が小さな時はコヒーレント揺動とスパイラル構造

帯状流成分，Mode３はコヒーレントなケルビンヘルムホの強度は共に増加していくが，ある程度コヒーレント揺動
８２

Special Topic Article

図４4. Data-Driven Analyses of Dynamical Turbulence Phenomena

M. Sasaki et al.

上段：各DMDモードの相関係数Cj の時間発展．下段：各リミットサイクル周期で描いた（C3, C4）間のリサージュ．ここでC3 はコヒーレントな KH 揺動，C4 はスパイラル構造に対応［６］．

が大きくなると，スパイラル揺動の強度が急増し，コヒーの次元%は時間ステップ数"よりも極端に小さくなること

レント揺動の強度が大きく減少することがわかる．この増が多い．これは設置できるセンサー数に限りがあるためで減を繰り返し，リミットサイクル振動が実現していることあり，数値シミュレーションとの大きな違いである．固有がわかる．また，スパイラル構造は，コヒーレント揺動の値分解に基づく標準的動的モード分解アルゴリズムでは，非線形飽和振幅を規定していることを示唆している．

モード数は最大でも % しか得られない．仮に，モード数が以上のように，DMD モードを用いて乱流構造間の因果""!%#以上あれば，観測ノイズまで含めてデータを正確関係を同定する方法を紹介した．リミットサイクル振動をに再構成できる．換言すると，回帰モデル #"!!#! は，縦伴うケルビンヘルムホルツ乱流に適用した結果を示した．

長のデータ行列には厳密にフィットし得るのに対し，極端リミットサイクルのように乱流の増減を繰り返す場合は，に横長のデータ行列には表現力が十分でない．別の見方を

畳み込み積分による構造 強度 の評価が有効である．得らすると，動的システムを構成する本質的なモードの数よりれた構造強度について，モード間のリサージュを描くことも少ないモード数でフィッティングすると，挙動の近しいで因果関係を抽出できることを示した．コヒーレント揺動

複数の本質的なモードが近似的に１つのモードでまとめての強度がスパイラル構造の出現に規定されていることを示表現される．さらに観測ノイズがあれば，ノイズに対応し唆する結果を得た．この方法は，他の様々な揺動観測に適て付加されるモードもまとめて表現される．こうして得ら用することが容易であり，幅広い波及効果を持つと期待しれたモードの重ね合わせでは，十分正確にデータを再構成ている．また，さらに拡張された発展的な動的モード分解

することができない．

手法もあり（例えば recursive DMD［３５］），今回の解析でこのような問題を解決するために，動的システムが有すは使用していない複素周波数の時間変化や位相ダイナミクるモード（本質的なモードとノイズ由来のモード）の数以スを抽出することで，乱流ダイナミクスの理解から輸送へ

上のモードを用いてフィッティングすることで必要なモーの理解への展開も期待される．また，ニューラルネット

ドを正確に推定したのちに，不要なモードを取り除くといワークを応用した構造間相互作用への発展も開発が進めらうアプローチを提案する．それぞれ，ハンケル行列作成とれている．

スパース動的モード分解によって実装することができる．

４．
３．
２ 乱流現象の定量的予測への応用

ハンケル行列は，元の時系列データから窓幅 'とスライド
幅 &で切り出したデータを $回積み上げて，次のように定動的モード分解をプラズマ乱流の観測データに適用する義される．

ことで，乱流挙動を予測する試みについて紹介する．ここでは，ハンケル行列とスパース動的モード分解を用いた方法について述べる．
一般に，実験から時系列データを取得する場合，データ
８３

Journal of Plasma and Fusion Research Vol.97, No.2 February 2021

"
&&
!
%!('
!
&
%
%&&
!
%!("'&&% & &+"(+'#&%!-(""'%&#: $%
%!
&
&&
!
%!("#
#
'
，
&
%
&
%
.
$
#
&&
(!
%'
!

:$+
.
!
!
!
-!
&&
)
!
)
"*
'
""!
"#!
"*,
)!
)
)
)

（５）

（６）

このハンケル行列に対して動的モード分解を適用すること
で，モード数は (+ だけ増加する．標準的動的モード分解
では，時刻 .における観測ベクトル !
.は行列 !の固有値分解から，式（２）のように展開された．スパース動的モード
分解では以下の手続きに よ り 寄 与 率 0*を 決 め 直す．まず，観測した時空間行列と再構成した時空間行列の差のFrobenius ノルムを最小化するための目的関数を $
"'と
&

図５

DMD モードの複素周波数．青がスパース DMD で除去された DMD モード，赤の DMD モードのみを解析に使用した．

図６

乱流時空間発展の学習と予測：上段から順にオリジナル計
測データ，ハンケル・スパース DMD を用いた再構成と予測，標準的 DMD を用いた再構成と予測［４０］．t < 3000［μs］
は学習領域，t > 3000［μs］
は予測領域．

し，L１正則化項を加えた次の最適化問題を解く．
"'
0
+)
,)
+)
0
'
&
"$' *
.
" $
**

（７）

*

最適化でゼロとなった寄与率 0*を拘束した条件下で，再び正則化項を外した $"'のみ最小化して，最終的な寄与率&0*を決定する．以上の手続きで，重要度の低いモードが取り除かれる．アルゴリズムの詳細は，文献［３６］を参照されたい．
我々は，動的モード分解による乱流データの解析を進めてきた［７］．解析に使用した乱流データは，九州大学にある直線磁化プラズマ PANTA 装置で観測されたイオン飽和
電流を用いた．周方向に６４チャンネル設置されているラン
グミュアプローブアレイによって１マイクロ秒ステップで計測された時空間発展を使用した．PANTA 装置では，放電条件によって帯状流やストリーマ，孤立波状態等の様々な乱流状態が存在するが［３７］，ストリーマ構造が支配的な条件で解析を行った［３８］．この放電では，5 kHz 程度のドリフト波乱流と 1 kHz 程度のストリーマ構造が共存している．このデータからハンケル行列を作成することは，上流側にプローブアレイを擬似的に増設していることに相当する．ここでは，-$"!!と ($%の設定で，時間区間%$3 ms測である［４０］．L1 正則化のパラメータ $は，DMD モードのデータで学習し，その先の予測を試みた．

の数が計測点数６４以下となるように調節した．標準的動的

この計測データにおいて，スパース動的モード分解の適モード分解では，時空間パターンは時間とともに減衰し用は，時空間データからのノイズ除去に相当する．ノイズて，細部も正確に再構成できていないのに対し，提案アプ除去の様子を，DMD モードの複素周波数の観点から図５

ローチではこの問題が大幅に改善されている．このようにに示す．ここで，複素平面における単位円付近の成分（赤精度良く再構成する動的モード分解を用いると，乱流の空丸）は動的システムの本質的なモードであり，単位円から間構造を含めた時間発展をより正確に予測できる．予測性離れた成分（青丸）はノイズ由来のモードと考えられる．

能を定量的に議論するため，以下の誤差指標"&.
'を導入す

実際，青丸で示した DMD モードは減衰が大きく，周波数る．

が高く，寄与率が小さい．ただし，ここで述べるノイズ除&$

去の方法は高周波なランダムノイズに対するもののみに適

"&
.
'
$
(

用可能なものであるため，計測器の基本的なノイズ対策は必要である［３９］．以降の結果は，赤丸で示した DMD モー#'' )$"(
/
.
.
'
!/
'
)
)&
)&
&$

#
' )$"/
.
'
)&

（８）

'
ここで，/
，計測
.
.
"%#.
'と /
'は，それぞれ，時刻 .
&
'
)&
)&

ドのみから得られるものである．

点 )における予測データと検証データである．誤差指標に

図６は計測データ（訓練データと検証データ），提案ア
プローチおよび標準的動的モード分解における再構成と予

対して，
８４

Special Topic Article

4. Data-Driven Analyses of Dynamical Turbulence Phenomena

M. Sasaki et al.

17H06089，16H02442，15H02155，15H02335，18K03578，
19J00871，九州大学応用力学研究所共同利用研究および
JST CREST JPMJCR1913 の助成をいただきました．R.O.
Dendy 教授，稲垣滋教授との有益な議論に感謝いたしま
す．

参考文献

図７

［１］P.H. Diamond et al., Plasma Phys. Control. Fusion 47, R35(2005).
［２］J.W. Connor et al., Phys. Plasmas 5, 2687 (1998).
［３］F. Wagner, Plasma Phys. Control. Fusion 49, B1 (2007).
［４］J.D. Callen and M.W. Kissick, Plasma Phys. Control. Fusion 39, B173 (1997).
［５］P.J. Schmid, J. Fluid Mech. 656, 5 (2010).
［６］M. Sasaki et al., Plasma Phys. Control. Fusion 61, 112001(2019).
［７］A. Kusaba et al., Plasma Fusion Res. 15, 1301001 (2020).
［８］A.A Kaptanoglu et al., Phys. Plasmas 27, 032108 (2020).
［９］R. Taylor et al., Rev. Sci. Instrum. 89, 053501 (2018).
［１０］H. Natsume et al., Phys. Plasmas 27, 042301 (2020).
［１１］C.W. Rowley et al., J. Fluid Mech. 641, 115 (2009).
［１２］P.J. Schmid, J. Fluid Mech. 656, 5 (2010).
［１３］K. Taira et al., Aiaa Journal 55, 4013 (2017).
［１４］S. Klus et al., J. Nonlinear Sci. 30, 283 (2020).
［１５］I. Mezi!, Annu. Rev. Fluid Mech. 45, 357 (2013).
［１６］B. Koopman, Proc. National Academy of Sciences of theUnited States of America 17, 315 (1931).
［１７］A. Mauroy et al., Physica D 261, 19 (2013).
［１８］J.H. Tu et al., J. Comput. Dyn. 1, 391 (2014).
［１９］M.O. Williams et al., J. Nonlinear Sci. 25, 1307 (2015).
［２０］Y. Kawahara, Neural Inf. Process. Syst. 29, 911 (2016).
［２１］N. Takeishi et al., Adv. Neural Inf. Process. Syst. 30, 1130
(2017).
［２２］B. Lusch et al., Nat. Commun 9, 4950 (2018).
［２３］M. Korda and I. Mezic, J. Nonlinear Sci. 28, 687 (2018).
［２４］B.J. Cantwell, Annu. Rev. Fluid Mech. 13, 457 (1981).
［２５］B.W. Brunton et al., J. Neurosci. Methods 258, 1 (2016).
［２６］J.L. Proctor and P.A. Eckhoff, Int. Health 7, 139 (2015).
［２７］J. Mann and J.N. Kutz, Finance 16, 643 (2016).
［２８］K. Fujii et al., Sci. Rep. 9, 16755 (2019).
［２９］K. Fujii et al., Proc. 2017 European Conf. on Machine
Learning and Principles and Practice of Knowledge Discovery in Databases (ECML-PKDD'17), pages 127 (2017).
［３０］T. Bito et al., Proc. 2019 Int'l Joint Conf. on Neural Networks (IJCNN'19), number 19278 (2019).
［３１］K. Fujii et al., Sci. Rep. 10, 3005 (2020).
［３２］H. Shiraishi et al., J. Neural Eng. 17, 036009 (2020).
［３３］B.M. Annaratone et al., Phys. Plasmas 18, 032108 (2011).
［３４］M. Sasaki et al., Phys. Plasmas 26, 042305 (2019).
［３５］B.R. Noack et al., J. Fluid Mech. 809, 843 (2016).
［３６］M. R. Jovanovi!et al., Phys. Fluids 26, 024103 (2014).
［３７］T. Kobayashi et al., Plasma Fusion Res. 12:1401019 (2017).
［３８］T. Yamada et al., Nat. Phys. 4, 721 (2008).
［３９］永 島 芳 彦，稲 垣 滋：プ ラ ズ マ・核 融 合 学 会 誌 96, 2(2020).
［４０］A. Kusaba et al., submitted to Plasma Fusion Res. (2020).
［４１］F. Kin et al., Phys. Plasmas 26, 042306 (2019).

予測誤差指標 E（t）の時間発展．t ＞ 3000［μs］
が予測領域．

!#
"
"
!
"
!$
#%$#
$
!!#
$

（９）

というフィッティングを行うことで，予測可能時間 "を評価する．図７は !#"$の対数プロットであり，"""!!!マイクロ秒程度であると評価できる．
現状では支配的乱流揺動（周期２００マイクロ秒程度）の数
周期程度の予測が空間構造を含めてできている．解析に用いたショットでは，ストリーマに起因するミリ秒スケールの空間局在した輸送が起こっている［４１］．今後，より長時間の予測が可能となれば，このような突発的・間欠的な輸送現象の空間構造を含めた予測にも応用できる可能性があ
る．

４．
４ まとめ
本章ではデータ駆動科学的手法として近年注目されている動的モード分解について，理論的基礎から応用例まで紹介した．動的モード分解は流体現象の時空間構造を観測／
計測データから抽出するための数値解析手法の一つとして提案されたものであり，ここではプラズマ乱流への応用について紹介した．動的モード分解を用いると，時空間発展
データから統計的に重要な空間構造とその構造に対応する
複素周波数を得ることができ，より少数の自由度でシステ
ムを理解することが可能となる．物理機構の理解をめざした応用例として，大域的乱流シミュレーションで得られた乱流場に適用した例を紹介した．ここでは動的モード分解
によって抽出された特徴的空間構造を元データとの畳み込
み積分を行うことで構造間の動的因果関係を抽出する方法を提案した．また，時空間発展予測へ向けた応用例として，
実験で観測された乱流の時空間発展について，ハンケル行
列とスパース動的モード分解を適用し，空間構造を含む乱流の発展を予測した試みを紹介した．今後，突発的輸送等の空間発展予測を可能にすると期待される．

謝辞

本研究は九州大学若手研究者グローバルリーダー育成型
Progress 100（NB８
０６４５０２８）および九州大学数理・データ
サイエンスに 関 する 教 育・研 究 支 援 プ ロ グ ラ ム の支援を受けました．また JSPS 科研費 JP18H03287，16K18335，

８５