出所：https://www.nagare.or.jp/download/noauth.html?d=30-2rensai2.pdf&dir=54

0.はじめに

機械学習において、主成分分析はデータの可視化や次元圧縮等に使用される。この主成分分析を流体の数値計算に適用すること（SnapshotPOD）で流れの空間モードの分解や数値計算の縮退化（高速化）が行える。
今回、SnapshotPODについて調べたので簡単にまとめてみる。

1.主成分分析（Principal Component Analysis）

機械学習の分野では、多次元データに対して変数間の分散が大きくなる方向にデータを射影することで低次元データに次元削減する目的で主成分分析が行われる。低次元のデータにできると、（例えば２次元、３次元に縮退できると）データの可視化が行えたり、各射影空間での軸が直行することが保証されるのでマルチコ（多重共線性）の回避ができる。
具体的には、下記のように１次元目に条件数を並べ、２次元目に説明変数を並べたデータ行列をX(N×M)とする。
このデータ行列の分散共分散行列に関する固有値問題を解くことで固有値と固有ベクトル行列が求まる。得られた固有ベクトル行列（固有値が大きい順に並び替え）とデータ行列とで行列積をとることで主成分得点が得られる。第１主成分はデータを最も分散が大きくなる軸に射影したデータ、第２主成分はそれに直行し、次に分散が大きくなる軸に射影したデータとなります。

pythonでは固有値問題は以下のように実装できる。この例では固有値問題を解き、第二主成分までを使用して次元削減を行っている。

import numpy as np
import numpy.linalg as LA
#X.T@Xの固有値、固有ベクトルを求める
sigma, U = LA.eig(X.T@X) 
#固有値が大きいものから並べ直し
indices = np.argsort(-sigma) 
sigma = sigma[indices]
 #Uもsigmaに合わせて並べ直す
U = U[:, indices] 
X_low = X_low @ U　
#指定した次元まで縮退化する
n_components = 2
X_new = X_low[:n_components, :]

上記の実装例では固有値問題を解いたが、特異値分解(SVD)によっても主成分分析ができる。例えば、scikit-learn.decompositionのPCAは特異値分解(SVD)が実装されており、下記のコードで主成文分析ができる。

scikit-learn.org

from sklearn.decomposition import PCA

pca = PCA(n_components=2)
X_new = pca.fit_transform(X)

この縮退化（次元削減）技術を流体の流れ場のモード分解に適用したものがSnapshotPODである。

2.SnapshotPOD

流体の計算では多数の計算格子を用意して支配方程式の時間発展を解かなければならないため、計算コストが一般的に高い。 SnapshotPODにより、流体の流れ場の時間発展データを時間軸方向に次元削減し、低次の空間モードを抽出することで計算コストを大幅に削減することができる。(後述) 以下では、２次元の流体計算を例にSnapshotPODを考えてみる。
まず、x座標、y座標の格子点での物理量を1つの縦ベクトルが1タイムステップを表すように変換する。この変換により、１次元目各時刻の計算格子全点での物理量（速度）、２次元目が時間になる。このデータ行列をXとして固有値問題を解くことで固有値、固有ベクトル行列が求まる。
ここで得られた固有値、固有ベクトル及びデータ行列Xを使用して空間モード（POD基底という）は以下のように求まる。 ※一般的にはデータ行列から時間平均流れ場を引き算したデータ行列の固有値問題を解く。 Pythonでは以下のように実装できる。

import numpy as np
import numpy.linalg as LA
#二次元解析データ(N1, N2, M)を(N1+N2, M)に変換
img_data = img_data.reshape(-1, N1*N2).T
#snapshot POD
sigma, U = LA.eig(img_data.T@img_data) #X.T@Xの固有値、固有ベクトルを求める
sigma = np.sqrt(sigma2)
sigma = np.nan_to_num(sigma) #Nanが発生していたので、0に置き換え
indices = np.argsort(-sigma) #固有値が大きいものから並べ直し
sigma = sigma[indices]
U = U[:, indices] #Uもsigmaに合わせて並べ直す
sigma_inv = np.diag(1/sigma) #∑-1を計算する
V = img_data_mat_l @ U @ sigma_inv #V = X @ U @ ∑-1で計算できる　
#空間モードを抽出する
mode = 1
V[mode, :]

POD基底は固有値の大きいモードほど、空間的に大きな構造をもつ。例えば、下の例では1次モード（上から順に固有値が大きい）が最も空間的に大きな渦構造をもつモードが抽出されてる。固有値が小さくなるほど、小さな渦構造を示している。
このようにSnapshotPODにより、時間発展する流れ場に対して空間モードが抽出でき、流れを支配する空間構造を把握できる。
さらに、SnapshotPODを使用することで流体解析の縮退化ができる。各空間モードに対して時間的に変化する係数をかけたものの足し合わせで流れ場を表現できると考えられる。(波の重ね合わせ、モードの重ね合わせ)
これによりR(<<M)までのモードと時間に依存する係数の和で流れを再構築できる。CFDで計算する場合、空間に離散化した(N点)の時間発展を解かなければならないが、SnapshotPODを使用すると、R(<<N)次元の各モードの係数の時間発展を解けばよくなり、大幅に計算コストを削減できる。

まさおのデータサイエンティストへの道

SnapshotPODについてまとめてみる

0.はじめに

1.主成分分析（Principal Component Analysis）

2.SnapshotPOD