GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

2020/2/1 GFD オンラインセミナー


Wigner変換を用いた流体波動の局所分離解析

話題提供者: 大貫陽平 (九州大学 応用力学研究所)


論文: Onuki 2020, Quasi-local method of wave decomposition in a slowly varying medium, Journal of Fluid Mechanics, 883, A56

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

目次

1. 導入: 地球流体における線形波動論の基礎

2. Wigner変換と擬微分作用素

3. Wigner変換を用いた流体方程式の理論解析

4. シミュレーションデータ分析への利用

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

1. 導入: 地球流体における線形波動論の基礎

  • Fourier変換に基づく分散関係, 偏光関係の解析

  • 空間に依存した系におけるWKB近似

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

はじめに:

  • 回転成層流体である, 大気や海洋には多種多様な波動が存在する.

    • 例: 音波, 慣性重力波, Rossby波, Kelvin波.
  • 大気/海洋に見られる流体現象は波動力学の観点で説明されることが多い.

    • 代表的なものとしては, 波動-平均流相互作用による惑星規模循環の形成.

    • 波動間相互作用(海洋における内部重力波, 大気Rossby波)も重要なテーマ.

  • 相互作用を考える際には, 非線形性が本質的に重要.

    • 大抵の場合は摂動法に基づいて議論される.

    • 線形波動論がすべての基礎になる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

GFDにおける線形波動論の基礎

回転系の浅水方程式:

流速を(u,v)(u, v), 水面変位をη\eta, 重力加速度をgg, コリオリパラメータをff, 静止状態での水深をHHとして,

tufv=gxη\partial_t u - f v = -g \partial_x \eta

tv+fu=gyη\partial_t v + f u = -g \partial_y \eta

tη+x(Hu)+y(Hv)=0\partial_t \eta + \partial_x (H u) + \partial_y (H v) = 0

  • 境界条件: 無限遠で, (u,v,η)0(u, v, \eta) \to 0

この方程式にはどのような性質の波動が含まれるか?

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

初等的な考察: (f,H)(f, H) が定数の場合

指数関数解: (u,v,η)=(u~,v~,η~)ei(kx+lyωt)(u, v, \eta) = (\tilde{u}, \tilde{v}, \tilde{\eta}) e^{{\rm i} (kx + ly - \omega t)} を方程式に代入, 微分作用素を波数と周波数に置き換える:

iωu~fv~=igkη~- {\rm i} \omega \tilde{u} - f \tilde{v} = - {\rm i} g k \tilde{\eta}

iωv~+fu~=iglη~- {\rm i} \omega \tilde{v} + f \tilde{u} = - {\rm i} g l \tilde{\eta}

iωη~+ikHu~+ilHv~=0- {\rm i} \omega \tilde{\eta} + {\rm i} k H \tilde{u} + {\rm i} l H \tilde{v} = 0

  • 固有値問題を解いて周波数が定まる: ω=0,±f2+gH(k2+l2)\omega = 0, \pm \sqrt{f^2 + gH(k^2 + l^2)} ... 分散関係

  • それぞれの分散関係に対応して, 従属変数の比である(u~,v~,η~)(\tilde{u}, \tilde{v}, \tilde{\eta})が定まる ... 偏光関係

独立した運動形態である, 慣性重力波モードと定常渦モードが得られる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

なぜ指数関数解を仮定して良い?

根底にあるのはFourier変換:

(u,v,η)=12π(u~(k,l),v~(k,l),η~(k,l))ei(kx+ly)dkdl(u, v, \eta) = \frac{1}{2\pi} \int (\tilde{u}(k,l), \tilde{v}(k,l), \tilde{\eta}(k,l)) e^{{\rm i} (kx + ly)} dk dl

  • それぞれの波数成分は独立に振る舞う.

  • あらゆる状態が波数成分の重ね合わせで表現できる.

  • 与えられた(u,v,η)(u, v, \eta)に対し, Fourier変換 \to 波数空間で成分分解 \to Fourier逆変換 という手順を踏むことで, 独立したモード(慣性重力波と定常渦)を分離することができる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(f,H)(f, H) が定数でない(x,yx, yの関数となる)場合

  • 波長に比べて f,Hf, Hの変化の空間スケールが十分大きいと仮定.

  • 長い空間スケールを導入し, μ(x,y)(x,y)\mu(x, y) \to (x, y), μ1\mu \ll 1と置き換える.

    • 空間微分が(x,y)μ(x,y)(\partial_x, \partial_y) \to \mu (\partial_x, \partial_y), と変更される.
  • 解の形式を (u,v,η)=(u~(x,y),v~(x,y),η~(x,y))ei(θ(x,y)/μσt)(u, v, \eta) = (\tilde{u}(x, y), \tilde{v}(x, y), \tilde{\eta}(x, y)) e^{{\rm i}(\theta(x, y) / \mu - \sigma t)} と仮定して方程式に代入 (WKB近似), O(μ)O(\mu)の項まで考慮して, (u~,v~,η~,θ,σ)(\tilde{u}, \tilde{v}, \tilde{\eta}, \nabla\theta, \sigma)の関係式を得る.

  • 局所的な波数kθ\boldsymbol{k} \equiv \nabla \thetaと定義すると, 分散関係式と偏光関係式が波数と空間位置に依存した形式で, σ=ω(x,k)\sigma = \omega(\boldsymbol{x}, \boldsymbol{k}), (u~,v~,η~)u~=u(x,k)(\tilde{u}, \tilde{v}, \tilde{\eta}) \equiv \tilde{\bm{u}} = \bm{u}(\boldsymbol{x}, \boldsymbol{k})と得られる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

指数関数解とWKB近似解の比較

解の形式 分散関係式 偏光関係式 波動分離
指数関数解 u~ei(kxσt)\tilde{\bm{u}} e^{{\rm i} (\boldsymbol{k} \cdot \boldsymbol{x} - \sigma t)} σ=ω(k)\sigma = \omega(\boldsymbol{k}) u~=u(k)\tilde{\bm{u}} = \bm{u}(\boldsymbol{k}) Fourier変換法
WKB近似解 u~ei(θ(x)/μσt)\tilde{\bm{u}} e^{{\rm i} (\theta(\boldsymbol{x}) / \mu - \sigma t)} σ=ω(x,k)\sigma = \omega(\boldsymbol{x}, \boldsymbol{k}) u~=u(x,k)\tilde{\bm{u}} = \bm{u}(\boldsymbol{x}, \boldsymbol{k}) ?

(WKB近似における波数はk=θ\boldsymbol{k} = \nabla \theta として定義)

Q1. WKB近似解に対応する波動分離手法はないか?

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

空間に依存した系に特有の性質

浅水方程式系において, 係数の空間変化(f(y),H(x,y))(f(y), H(x,y))を考慮した場合:
惑星/地形性ベータ効果が生じて渦モードが有限の周波数を得る. (Rossby 波)

: df/dy=β,H=const.df / dy = \beta, H = const.とした時のRossby波:

ω=μβkk2+l2+f2/(gH)+O(μ2)\omega = \mu \frac{- \beta k}{k^2 + l^2 + f^2 / (gH)} + O(\mu^2)

浅水方程式で単に指数関数解を代入しただけでは導けない.

WKB近似を使って摂動項をゴリゴリ計算すれば導ける.
背景流の影響等を考慮するとかなり複雑な計算になる.

Q2. 分散関係式の摂動項を見通し良く求める方法はないか?

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

波動平均流相互作用研究の観点から

エネルギー/擬運動量フラックスの解析. (群速度則: F=ACg\bm{F} = A \bm{C}_g を満たすように定式化)

基本場に対する擾乱について線形化された方程式の理論解析が主な内容.

  • WKB近似を前提とするが, 摂動項を無視することで実質的には指数関数解を仮定.

  • Kinosita & Sato (2013) は, プリミティブ方程式をvvのみの式にしてから指数関数解を代入することで, Rossby波の分散関係式を導いている ... ややAdhocな印象.

問題設定をより複雑にして, 背景流の勾配とffの変化のスケールが同程度になるような場合(近慣性波やRossby波の分散関係が背景流シアを含む)には, 既存の方法論では太刀打ちできない.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

本研究: Wigner変換を用いて変数係数偏微分方程式の理論解析の見通しを良くする

複数種の波動を含む連立方程式を, 各波動成分が従う単独方程式に分解する.

Wigner_usage

Fourier変換に基づく波動分離では代数行列の固有値分解であったのに対し, Wigner変換に基づく波動分離では表象行列の対角化を行う.

ポイント: 物理量データを固有波動成分へ射影する作用素を導く.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

2. Wigner変換と擬微分作用素

  • 擬微分作用素と表象

  • スター積による計算

  • 表象行列に関するアルゴリズム

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

通常のフーリエ変換では, 微分作用素を波数に置き換える.

xik\partial_x \to {\rm i} k

Wigner変換は, 一般の線形作用素を関数に変換する.

代表的な例: 量子力学のHamiltonianから古典力学のHamiltonianへの対応,

22mx2+V(x)p22m+V(x)- \frac{\hbar^2}{2m} \partial_x^2 + V(x) \to \frac{p^2}{2m} + V(x)

上記の例では, 一座標xxを残したまま, 微分作用素ix- {\rm i} \hbar \partial_xを運動量ppに置き換えている. Wigner変換によって, 変数の自由度は2倍になる: x(x,p)x \to (x, p)

変換前の作用素と変換後の関数の関係は, 擬微分作用素とその表象の関係に対応.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

擬微分作用素(pseudo-diffential operator)は微分作用素を一般化した概念であり, フーリエ変換と逆変換を利用して, 積分変換の形式で

F^ϕ(x)=12πμf(x+x2,p)ϕ(x)eip(xx)/μdxdp\hat{F} \phi(x) = \frac{1}{2\pi\mu} \int f \left(\frac{x + x'}{2}, p \right)\phi(x') e^{{\rm i} p(x-x')/\mu} dx' dp

といった形で定義される. (上式でμ\muは適当な定数であり, 量子論ではμ=\mu = \hbarとする)

  • 被積分関数に含まれるf(x,p)f(x, p)は, 擬微分作用素F^\hat{F}表象(symbol)と呼ばれる.

  • 微分作用素も擬微分作用素の一種: f=(ip/μ)nf = ({\rm i} p / \mu)^n \leftrightarrow F^=xn\hat{F} = \partial_x^n.

  • 表象をf=(ip/μ)α,αRf = ({\rm i} p / \mu)^\alpha, \alpha \in \mathbb{R}という形に設定することで、逆ラプラシアンや分数階微分なども擬微分作用素として定義できる.

  • Wigner変換とは, 擬微分作用素F^\hat{F}をその表象f(x,p)f(x, p)へ対応させる変換である.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

Wigner変換の具体的手順:

F^ϕ(x)=F(x,x)ϕ(x)dx\hat{F} \bm{\phi}(x) = \int F (x, x') \bm{\phi}(x') dx'

  • 核関数を相対位置座標についてFourier変換したものが作用素の表象になる:

f(x,p)=F(x+x2,xx2)eipx/μdxf (x, p) = \int F \left(x + \frac{x'}{2}, x - \frac{x'}{2} \right) {\rm e}^{- {\rm i} p x' / \mu} dx'

例: 移流作用素 u(x)xu(x) \partial_x を考えると, 核関数は F(x,x)=u(x)δ(xx)F(x, x') = u(x) \delta'(x - x')であり, これをフーリエ変換してf(x,p)=i(p/μ)u+xu/2f(x,p) = - {\rm i} (p/\mu) u + \partial_x u / 2を得る. 単に微分作用素を波数へ置き換えた場合に比べて追加の項が生じる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析
  • 作用素の線形演算(足し引きやスカラー倍)は, 表象の演算に自然に対応する.
    作用素の積は単純には対応付けられない.
作用素1 作用素2 作用素の線形演算 作用素の演算
作用素 F^\hat{F} G^\hat{G} αF^+βG^\alpha \hat{F} + \beta \hat{G} F^G^\hat{F} \hat{G}
表象 f(x,p)f(x, p) g(x,p)g(x, p) αf(x,p)+βg(x,p)\alpha f(x, p) + \beta g(x, p) f(x,p)g(x,p)f(x, p) \star g(x, p) (スター積)

積演算の定義が特殊なのは, 擬微分作用素が一般に非可換であるため

  • 基本となる2つの作用素をp^iμx,x^x\hat{p} \equiv - {\rm i} \mu \partial_x, \hat{x} \equiv x と定義したとき, x^p^p^x^\hat{x}\hat{p} \neq \hat{p}\hat{x}.

WKB近似ではμ1\mu \ll 1を仮定. 最も粗い近似では, 作用素を可換とみなす. 作用素の非可換性によって摂動が生じる. このことがRossby波の発生に関与する. (ffy\partial_yが非可換)

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

擬微分作用素の積に対応する表象の演算: スター積

擬微分作用素の定義式(p15)に立ち返って考察すると

f(x,p)g(x,p)f(x, p) \star g(x, p)

=1(πμ)2 ⁣ ⁣ ⁣ ⁣f(x+x,p+p)g(x+x,p+p)e2i(pxpx)/μdxdxdpdp= \frac{1}{(\pi \mu)^2} \iint\!\!\!\!\iint f(x+x', p+p') g(x+x'', p+p'') e^{2 {\rm i} (p'' x' - p' x'') / \mu} dx' dx'' dp' dp''

が得られる. 被積分関数のffgg(x,p)(x, p)の周りでテイラー展開し, 積分を実行すると,

f(x,p)g(x,p)=m=0,n=0(1)nm!n!(iμ2)(m+n)m+nfxmpnm+ngxnpmf(x, p) \star g(x, p) = \sum_{m=0,n=0}^\infty \frac{(-1)^n}{m!n!} \left( \frac{{\rm i} \mu}{2} \right)^{(m+n)} \frac{\partial^{m+n} f}{\partial x^m \partial p^n} \frac{\partial^{m+n} g}{\partial x^n \partial p^m}

となる. 結果として, パラメータμ\muのベキ級数展開で表現される. 以降は μ1\mu \ll 1を仮定し, 高次の項を摂動としてみなして話を進める.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

高次元系への拡張

一般の多変量関数 ϕ(x):RdCn\bm{\phi}(\bm{x}): \mathbb{R}^d \to \mathbb{C}^nに対して定義される擬微分作用素:

F^ϕ(x)=1(2πμ)df(x+x2,p)ϕ(x)eip(xx)/μdxdp\hat{\bm{F}}\bm{\phi}(\bm{x}) = \frac{1}{(2 \pi \mu)^d} \int \bm{f}\left( \frac{\bm{x} + \bm{x}'}{2}, \bm{p} \right) \bm{\phi}(\bm{x}') e^{i\bm{p}' \cdot (\bm{x} - \bm{x}') / \mu} d\bm{x}'d\bm{p}

ここで, f\bm{f}R2d\mathbb{R}^{2d}で定義されたn×nn \times nの行列値関数. スター積は

f(x,p)g(x,p)=m,n(1)nm!n!(iμ2)m+nxmpnf(x,p)xnpmg(x,p) \bm{f}(\bm{x},\bm{p}) \star \bm{g}(\bm{x},\bm{p}) = \sum_{\bm{m},\bm{n}}\frac{(-1)^{|\bm{n}|}}{\bm{m}!\bm{n}!} \left(\frac{{\rm i} \mu}{2}\right)^{|\bm{m}|+|\bm{n}|} \nabla^{\bm{m}}_x \nabla^{\bm{n}}_p \bm{f}(\bm{x},\bm{p}) \nabla^{\bm{n}}_x \nabla^{\bm{m}}_p \bm{g}(\bm{x},\bm{p})

(m\bm{m}n\bm{n}には多重指数記法を用いている)

以降は簡略化のため, 表象f(x,p)\bm{f}(\bm{x}, \bm{p})に対応する擬微分作用素をF^=f(x^,p^)\hat{\bm{F}} = \bm{f}(\hat{\bm{x}}, \hat{\bm{p}})等と表記する.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

Wigner変換の視覚的表現

width:1000px

  • 核関数の相対位置についてのFourier変換.

  • p\bm{p}準局所的に定義された波数.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

Wigner変換のご利益

作用素の性質や微積演算を, 表象行列の性質や演算へ置きかえる.

  • エルミート作用素 ↔ エルミート行列

  • 作用素の分解 ↔ 表象行列の分解

  • 逆作用素の計算 ↔ スター積についての逆要素の計算

  • 作用素の対角化 ↔ 表象行列の対角化

表象行列の対角化により, 流体方程式が含む波動成分を分離できる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

逆作用素の計算

ある作用素L^\hat{\bm{L}}に対し, 以下の式を満たす逆作用素L^1\hat{\bm{L}}^{-1}を求めたいとする:

L^L^1=I\hat{\bm{L}} \hat{\bm{L}}^{-1} = \bm{I}

(I^\hat{\bm{I}}は単位行列). ここで L^\hat{\bm{L}}L^1\hat{\bm{L}}^{-1}をWigner変換し, それぞれの表象をl(x,p)\bm{l}(\bm{x}, \bm{p}), lI(x,p)\bm{l}^{I}(\bm{x}, \bm{p})とすると, 解くべき方程式は l\bm{l} を既知として,

llI=I\bm{l} \star \bm{l}^{I} = \bm{I}

となる. この式を解くため, 表象を

l=l0+μl1+μ2l2+\bm{l} = \bm{l}_0 + \mu \bm{l}_1 + \mu^2 \bm{l}_2 + \ldots

lI=l0I+μl1I+μ2l2I+\bm{l}^I = \bm{l}^I_0 + \mu \bm{l}^I_1 + \mu^2 \bm{l}^I_2 + \ldots

と展開して代入する.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(逆作用素続き) 方程式の各項をμ\muのオーダーで分けると,

μ0: l0l0I=I\mu^0: \ \bm{l}_0 \bm{l}_0^I = \bm{I}

μ1: l1l0I+l0l1I+l0l0I1=0\mu^1: \ \bm{l}_1 \bm{l}_0^I + \bm{l}_0 \bm{l}_1^I + \overline{\bm{l}_0 \star \bm{l}_0^I}^1 = 0

\vdots

となる (表象a\bm{a}に対してa1\overline{\bm{a}}^1μ\muに比例した項を取り出すことを意味する). μ0\mu^0の式は通常の逆行列の計算により, l0I=l01\bm{l}^I_0 = \bm{l}^{-1}_0と解くことができ, 高次の式は

l1I=l01(l1l0I+l0l0I1)\bm{l}^I_1 = - \bm{l}_0^{-1} \left( \bm{l}_1 \bm{l}_0^I + \overline{\bm{l}_0 \star \bm{l}_0^I}^1 \right)

といった形で順次求めることができる. こうして得られたlI\bm{l}^Iを表象にもつ擬微分作用素が所期の L1\bm{L}^{-1}になる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

作用素の対角化

与えられた作用素行列H=h(x^,p^)\bm{H} = \bm{h}(\hat{\bm{x}}, \hat{\bm{p}})を, 対角形式に変換する問題を考える:

U^H^U^1=Ω^=diag(Ω^1,Ω^2,,Ω^n)\hat{\bm{U}} \hat{\bm{H}} \hat{\bm{U}}^{-1} = \hat{\bm{\Omega}} = {\rm diag}(\hat{\Omega}_1, \hat{\Omega}_2, \cdots, \hat{\Omega}_n)

表象行列に変換すると,

uhuI=ω=(ω1(x,p)000ω2(x,p)000ωn(x,p))\bm{u} \star \bm{h} \star \bm{u}^I = \bm{\omega} = \left( \begin{array}{cccc} \omega_1(\bm{x}, \bm{p}) & 0 & \ldots & 0 \\ 0 & \omega_2(\bm{x}, \bm{p}) & & \vdots \\ \vdots & & \ddots & 0 \\ 0 & \ldots & 0 & \omega_n(\bm{x}, \bm{p}) \end{array} \right)

となる. 通常の行列の固有値問題に似ているが, 積がスター積に置き換わっている.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(対角化続き1)

表象 h,ω,u,uI\bm{h}, \bm{\omega}, \bm{u}, \bm{u}^Iμ\muについて級数展開し, u\bm{u}uI\bm{u}^Iを以下のように書き換える:

u=(Iμ2u0u011+μα)u0+O(μ2)\bm{u} = \left( \bm{I} - \frac{\mu}{2} \overline{\bm{u}_0 \star \bm{u}^{-1}_0}^1 + \mu \bm{\alpha} \right) \bm{u}_0 + O(\mu^2)

uI=u01(Iμ2u0u011μα)+O(μ2)\bm{u}^I = \bm{u}^{-1}_0 \left( \bm{I} - \frac{\mu}{2} \overline{\bm{u}_0 \star \bm{u}^{-1}_0}^1 - \mu \bm{\alpha} \right) + O(\mu^2)

解くべき方程式からμ0\mu^0μ1\mu^1に比例した項を取り出すと, それぞれ

u0h0u01=ω0\bm{u}_0 \bm{h}_0 \bm{u}_0^{-1} = \bm{\omega}_0

u0hu01112u0u011ω012ω0u0u011+αω0ω0α=ω1 \overline{\bm{u}_0 \star \bm{h} \star \bm{u}_0^{-1}}^1 - \frac{1}{2} \overline{\bm{u}_0 \star \bm{u}^{-1}_0}^1 \bm{\omega}_0 - \frac{1}{2} \bm{\omega}_0 \overline{\bm{u}_0 \star \bm{u}^{-1}_0}^1 + \bm{\alpha} \bm{\omega}_0 - \bm{\omega}_0 \bm{\alpha} = \bm{\omega}_1

となる. 1つめの式は通常の行列の固有値問題であり、これは解けると仮定する.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(対角化続き2)

μ1\mu^1の式の性質は, μ0\mu^0の方程式の結果によって異なる:

  • (1) ω0\bm{\omega}_0のすべての対角成分 (ω0,1,ω0,2,)(\omega_{0,1}, \omega_{0,2}, \ldots) が異なる場合 (縮退がない)

  • (2) ω0\bm{\omega}_0の対角成分 (ω0,1,ω0,2,)(\omega_{0,1}, \omega_{0,2}, \ldots) に重複がある場合 (縮退がある)

(1) 縮退のない場合: ω\bm{\omega}u\bm{u}μ1\mu^1の補正項が次のように求まる:

ω1(k,k)=[u0hu011](k,k)ω0,k[u0u011](k,k) \bm{\omega}_{1(k,k)} = \left[ \overline{\bm{u}_0 \star \bm{h} \star \bm{u}_0^{-1}}^1 \right]_{(k,k)} - \omega_{0,k} \left[ \overline{ \bm{u}_0 \star \bm{u}_0^{-1}}^1 \right]_{(k,k)}

α(k,l)=1ω0,kω0,l[u0hu011](k,l)12ω0,k+ω0,lω0,kω0,l[u0u011](k,l) \bm{\alpha}_{(k,l)} = \frac{1}{\omega_{0,k}-\omega_{0,l}} \left[ \overline{ \bm{u}_0 \star \bm{h} \star \bm{u}^{-1}_0}^1 \right]_{(k,l)} - \frac{1}{2}\frac{\omega_{0,k} + \omega_{0,l}} {\omega_{0,k}-\omega_{0,l}} \left[ \overline{ \bm{u}_0 \star \bm{u}_0^{-1} }^1 \right]_{(k,l)}

(括弧書きの下付き添字は行列内における要素の位置を表す)

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(対角化続き3)

(2) 縮退がある場合:

  • μ1\mu^1の式は一般には計算できない.
    周波数が縮退した成分は分離することが出来ない.

特別な状況として縮退した固有値が0の場合 (ω0,i1=ω0,i2==0\omega_{0, i_1} = \omega_{0, i_2} = \ldots = 0):

  • μ1\mu^1の式を計算することができ, 縮退が解ける.
    → 周波数0に縮退していた成分が空間の不均質性により, 固有の周波数を獲得する.

[具体的な計算手順については論文(Onuki 2020)を参照]

回転成層流体では, β\beta効果により定常渦モードがRossby波として周波数を獲得することに相当する.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

ここまでのまとめ:

  • Wigner変換は, ベクトル値関数に対する一般の作用素を, 表象行列に変換する.

  • 作用素の積は表象行列のスター積に対応する.

  • 逆作用素の計算や作用素の対角化等の演算は, 対応する表象行列を級数展開することで, 代数的な計算アルゴリズムに帰着される.

変数係数偏微分方程式の厄介な計算を, Wigner変換を通して見通しよく行うことができる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

3. Wigner変換を用いた流体方程式の理論解析

  • 流体方程式の非正準ハミルトン形式

  • 線形化方程式の擬微分作用素表現

  • 単独方程式への分離

  • Wigner分布関数と輸送方程式

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

流体方程式の一般的表現

非粘性の流体運動を支配する方程式は, 一般に次の非正準ハミルトン形式で書かれる:

tυ=J^(υ)δHδυ   ()\partial_t \bm{\upsilon} = \hat{\bm{J}}(\bm{\upsilon}) \frac{\delta \mathcal{H}}{\delta \bm{\upsilon}} \ \ \ \cdots (\star)

  • υ(x,t)\bm{\upsilon}(\bm{x}, t)は状態ベクトルであり, 流速や密度等を表す.

  • H[v]\mathcal{H}[\bm{v}]はハミルトニアン汎関数で, エネルギーと適当なCasimirの和である.

  • J^\hat{\bm{J}}υ\bm{\upsilon}の関数として定義された歪エルミート作用素 (J^=J^\hat{\bm{J}}^\dag = -\hat{\bm{J}}を満たす).

地球流体力学における非正準ハミルトン形式の解説はShepherd (1990) あるいは Salmon (1998) を参照.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

方程式の線形化とエネルギーの定義

H[υ]\mathcal{H}[\bm{\upsilon}]υ=Υ\bm{\upsilon} = \bm{\Upsilon}で極小となると仮定, そこからの擾乱をυ=Υ+υ\bm{\upsilon} = \bm{\Upsilon} + \bm{\upsilon}'と定義して方程式()(\star)に代入, A^δ2H/δυ2,B^=iA^J^(Υ)A^\hat{\bm{A}} \equiv \delta^2 \mathcal{H} / \delta\bm{\upsilon}^2, \hat{\bm{B}} = {\rm i} \hat{\bm{A}} \hat{\bm{J}}(\bm{\Upsilon}) \hat{\bm{A}}として方程式を線形化すると

iA^tυ=B^υ  (){\rm i} \hat{\bm{A}} \partial_t \bm{\upsilon}' = \hat{\bm{B}} \bm{\upsilon}' \ \cdots \ (\star\star)

が得られる. ここでA^\hat{\bm{A}}は正値エルミート作用素, B^\hat{\bm{B}}はエルミート作用素.

この式から(擬)エネルギーの保存則が導かれる: E=(1/2)υA^υdx\mathcal{E} = (1/2)\int \bm{\upsilon}'^\dag \hat{\bm{A}} \bm{\upsilon}' d\bm{x}として,

E0;  E˙=0\mathcal{E} \geq 0; \ \ \dot{\mathcal{E}} = 0

が成立.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

作用素の分解

前のページで定義した正値エルミート作用素A^\hat{\bm{A}}を, Wigner変換と表彰行列の演算(論文参照)によって, 以下のように分解する:

A^=L^L^\hat{\bm{A}} = \hat{\bm{L}} \hat{\bm{L}}^\dag

ここで, L^\hat{\bm{L}}は可逆作用素であり, L^1\hat{\bm{L}}^{-1}が存在. 状態ベクトルをψ=L^υ/2\bm{\psi} = \hat{\bm{L}}^\dag \bm{\upsilon}' / \sqrt{2}と再定義すると, 方程式()(\star \star)

itψ=H^ψ  (){\rm i} \partial_t \bm{\psi} = \hat{\bm{H}} \bm{\psi} \ \cdots \ (\star\star\star)

と変形できる. ここでH^=L^1B^L^1\hat{\bm{H}} = \hat{{\bm{L}}}^{-1}\hat{\bm{B}}\hat{\bm{L}}^{\dag-1}はエルミート作用素. 系のエネルギーは, ψ\bm{\psi}のノルムとして, E=ψ2dx\mathcal{E} = \int |\bm{\psi}|^2 d\bm{x}と書かれる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

方程式の分離

Wigner変換によってエルミート作用素に対応する表象行列 h(x,p)\bm{h}(\bm{x}, \bm{p}) を定義, 固有値方程式uhuI=ω\bm{u} \star \bm{h} \star \bm{u}^I = \bm{\omega}を解き対角行列 ω=diag(ω1,ω2,)\bm{\omega} = diag(\omega_1, \omega_2, \ldots) を求める. 状態変数をψ=u(x^,p^)ψ\bm{\psi}' = \bm{u}(\hat{\bm{x}}, \hat{\bm{p}}) \bm{\psi}と再定義すると, 方程式(\star\star\star)は

itψ1=ω1(x^,p^)ψ1{\rm i} \partial_t \psi'_1 = \omega_1(\hat{\bm{x}}, \hat{\bm{p}}) \psi'_1

itψ2=ω2(x^,p^)ψ2{\rm i} \partial_t \psi'_2 = \omega_2(\hat{\bm{x}}, \hat{\bm{p}}) \psi'_2

\vdots

と分離できる. ここで, ωi(x^,p^)\omega_i(\hat{\bm{x}}, \hat{\bm{p}})はそれぞれエルミート作用素であり, u(x^,p^)\bm{u}(\hat{\bm{x}}, \hat{\bm{p}})は, ユニタリ作用素である.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

分散関係と偏光関係

  • 前のページで得られた表象行列ω(x,p)\bm{\omega}(\bm{x}, \bm{p})の対角成分(ω1,ω2,)(\omega_1, \omega_2, \ldots)はそれぞれの波動成分の分散関係を表す.

  • 対角化に用いた表象行列u\bm{u}の各列ベクトルが偏光関係を表す.

  • いずれの表象も, 背景パラメータの空間勾配(xf,xH)(\nabla_x f, \nabla_x H)波束の大きさの有限性を考慮に入れて, 準局所的に定義された物理量と言える.

  • ユニタリ作用素u(x^,p^)\bm{u}(\hat{\bm{x}}, \hat{\bm{p}})は, 状態ベクトルψ\bm{\psi}を固有空間へ写す射影作用素である.

以上の手順により, 空間に依存した系における波動分離を準局所的に行うことができる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

密度作用素とWigner分布関数

ある波動成分の状態ベクトルψ(x)\psi(\bm{x})に対し, 次のように定まる作用素W^\hat{W}を定義する:

W^ϕ(x)ψ(x)ψ(x)ϕ(x)dx\hat{W} \phi(\bm{x}) \equiv \psi(\bm{x}) \int \psi^\dag(\bm{x}') \phi(\bm{x}') d\bm{x}'

(ϕ(x)\phi(\bm{x})は任意の関数). このように定義されたW^\hat{W}は, 量子力学では密度作用素と呼ばれる. さらに, 密度作用素をWigner変換して得られる表象w(x,p)w(\bm{x}, \bm{p})は, Wigner分布関数と呼ばれる:

w(x,p)=ψ(x+x2)ψ(xx2)eipx/μdxw (\bm{x}, \bm{p}) = \int \psi \left( \bm{x} + \frac{\bm{x}'}{2} \right) {\psi}^\dag \left( \bm{x} - \frac{\bm{x'}}{2} \right) {\rm e}^{-{\rm i} \bm{p} \cdot \bm{x'} / \mu} d\bm{x}'

Wigner分布関数は物理-波数空間内で定義されたエネルギー密度と解釈でき, 積分値は全エネルギーに一致する: (1/2πμ)dw(x,p)dxdp=E(1/2\pi\mu)^d\iint w(\bm{x}, \bm{p}) d\bm{x}d\bm{p} = \mathcal{E}.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

輸送方程式

状態ベクトルの時間発展が itψ=ω(x^,p^)ψ{\rm i} \partial_t \psi = \omega (\hat{\bm{x}}, \hat{\bm{p}}) \psiであるとき, Wigner分布関数の時間発展は

itw=ωwwω{\rm i} \partial_t w = \omega \star w - w \star \omega

となる. スター積をμ\muで展開し, O(μ3)O(\mu^3)の項を無視すると,

tw+μpωxwμxωpw=0 \partial_t w + \mu \nabla_p \omega \cdot \nabla_x w - \mu \nabla_x \omega \cdot \nabla_p w = 0

となる. この式は物理-波数空間内におけるエネルギー密度の伝達を記述するものであり, 輸送方程式(あるいは放射伝達方程式など)と呼ばれる.

  • 輸送方程式の特性曲線は, x˙=pω,  p˙=xω\dot{\bm{x}} = \nabla_p \omega, \ \ \dot{\bm{p}} = - \nabla_x \omega で定められる.
    これは通常のRay-tracing 方程式である.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

エネルギーフラックス

輸送方程式を物理空間に射影(波数 p\bm{p} で積分)すると, エネルギー密度の時間発展が

tE+xF=0\partial_t E + \nabla_x \cdot \bm{F} = 0

と得られる. ここで,

E(x)=1(2μπ)dw(x,p)dpE(\bm{x}) = \frac{1}{(2\mu\pi)^d} \int w(\bm{x}, \bm{p}) d\bm{p}

F(x)=μ(2μπ)dw(x,p)pωdp\bm{F}(\bm{x}) = \frac{\mu}{(2\mu\pi)^d} \int w(\bm{x}, \bm{p}) \nabla_p \omega d\bm{p}

である. エネルギーフラックス == エネルギー密度 ×\times 群速度という形式になっている.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

ここまでのまとめ

  • 一般の流体方程式に対し, Wigner変換を利用して波動成分を準局所的に分離する手順を構築.

  • Wigner分布関数を導入することで, 物理-波数空間内でのエネルギー伝達を記述.

  • 群速度に比例した形式のフラックスを計算できることを示した.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

4. シミュレーションデータ分析への利用

  • 浅水方程式の擬微分作用素表現

  • 数値モデルデータから波動成分を分離してエネルギーフラックスを描画

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

浅水方程式の擬微分作用素表現

流速を(u,v)(u, v), 水深をhh, 重力加速度をgg, コリオリパラメータをf(y)f(y), 水底の地形をD(x,y)D(x, y)とする. ハミルトニアン(全エネルギー)を

H=(hu2+hv22+g(h+D)22)dx\mathcal{H} = \int \left( \frac{hu^2 + hv^2}{2} + \frac{g ( h + D )^2 }{2} \right)d\bm{x}

と定義して, 浅水方程式は非正準ハミルトン形式で次のように書ける:

(tutvth)=(0f/hxf/h0yxy0)(δH/δuδH/δvδH/δh)\left( \begin{array}{c} \partial_t u \\ \partial_t v \\ \partial_t h \end{array} \right) = \left( \begin{array}{ccc} 0 & f / h & - \partial_x \\ - f / h & 0 & - \partial_y \\ - \partial_x & - \partial_y & 0 \end{array} \right) \left( \begin{array}{c} \delta \mathcal{H} / \delta u \\ \delta \mathcal{H} / \delta v \\ \delta \mathcal{H} / \delta h \end{array} \right)

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(続き) 定常解として静止状態

u=0, v=0, h=H(x)H0D(x)u = 0, \ v = 0, \ h = H(\bm{x}) \equiv H_0 - D(\bm{x})

を考え, 前述の手順(p31)に沿って擾乱について線形化された方程式を導く:

i(H000H000g)Hermite operator(tutvtη)=(0ifHigHxifH0igHyigxHigyH0)Hermite operator(uvη){\rm i} \underbrace{ \left( \begin{array}{ccc} H & 0 & 0 \\ 0 & H & 0 \\ 0 & 0 & g \end{array} \right) }_{\rm Hermite \ operator} \left( \begin{array}{c} \partial_t u' \\ \partial_t v' \\ \partial_t \eta' \end{array} \right) = \underbrace{ \left( \begin{array}{ccc} 0 & {\rm i} f H & - {\rm i} gH \partial_x \\ - {\rm i} f H & 0 & - {\rm i} gH \partial_y \\ - {\rm i} g \partial_x H & - {\rm i} g \partial_y H & 0 \end{array} \right) }_{\rm Hermite \ operator} \left( \begin{array}{c} u' \\ v' \\ \eta' \end{array} \right)

状態ベクトルをψ=(H/2u,H/2v,g/2η)\bm{\psi} = ( \sqrt{H/2}u', \sqrt{H/2}v', \sqrt{g/2}\eta')と定義すると, 上の式は

itψ=h(x^,p^)ψ{\rm i} \partial_t \bm{\psi} = h(\hat{\bm{x}}, \hat{\bm{p}}) \bm{\psi}

と書ける. (表象行列h\bm{h}の具体的な形は次のページ)

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(続き) 浅水方程式を特徴づけるエルミート作用素の表象行列h(x,p)\bm{h}(\bm{x}, \bm{p}):

h=(0ifcpxif0cpycpxcpy0)h0+μ(00icx200icy2icx2icy20)h1\bm{h} = \underbrace{ \left( \begin{array}{ccc} 0 & {\rm i} f & - c p_x \\ - {\rm i} f & 0 & - c p_y \\ - c p_x & - c p_y & 0 \end{array} \right) }_{\displaystyle\bm{h}_0} + \mu \underbrace{ \left( \begin{array}{ccc} 0 & 0 & \displaystyle\frac{ {\rm i} c_x}{2} \\[2ex] 0 & 0 & \displaystyle\frac{ {\rm i} c_y}{2} \\[2ex] - \displaystyle\frac{{\rm i} c_x}{2} & - \displaystyle\frac{{\rm i} c_y}{2} & 0 \end{array} \right) }_{\displaystyle\bm{h}_1}

ここで, p=(px,py)\bm{p} = (p_x, p_y), c=gHc = \sqrt{gH}, xc=(cx,cy)\nabla_x c = (c_x, c_y)とした.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

浅水方程式の対角化

定められた手順(p24-p27)に沿い, 方程式uh(x,p)uI=ω\bm{u} \star \bm{h}(\bm{x}, \bm{p}) \star \bm{u}^I = \bm{\omega}を解いて作用素を対角化する. 偏光関係を表す行列u\bm{u}

u=(θpxifpy2θpθpy+ifpx2θpcp2θθpxifpy2θpθpy+ifpx2θpcp2θicpyθicpxθfθ)+O(μ)\bm{u} = \left( \begin{array}{ccc} \displaystyle\frac{\theta p_x - {\rm i} f p_y}{\sqrt{2} \theta p} & \displaystyle\frac{\theta p_y + {\rm i} f p_x}{\sqrt{2} \theta p} & \displaystyle\frac{cp}{\sqrt{2} \theta} \\[2ex] \displaystyle\frac{ - \theta p_x - {\rm i} f p_y}{\sqrt{2} \theta p} & \displaystyle\frac{ - \theta p_y + {\rm i} f p_x}{\sqrt{2} \theta p} & \displaystyle\frac{cp}{\sqrt{2} \theta} \\[2ex] \displaystyle\frac{{\rm i} c p_y}{\theta} & \displaystyle\frac{- {\rm i} c p_x}{\theta} & \displaystyle\frac{f}{\theta} \\ \end{array} \right) + O(\mu) \\

となる. (\bigl( p=px2+py2p = \sqrt{p_x^2 + p_y^2}, θ=f2+c2p2\theta = \sqrt{f^2 + c^2 p^2}とした)\bigr)

アルゴリズムの実行にはMaximaを利用.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(続き) 分散関係を表す ω=diag(ωG+,ωG,ωR)\bm{\omega} = diag(\omega_{G+}, \omega_{G-}, \omega_R)は,

ωG±=±f2+c2p2μ[βpx(2f2+c2p2)2p2(f2+c2p2)+fc(pxcypycx)f2+c2p2]+O(μ2)\omega_{G\pm} = \pm \sqrt{ f^2 + c^2 p^2 } - \mu \left[ \frac{\beta p_x (2f^2 + c^2 p^2)}{2 p^2 (f^2 + c^2 p^2)} + \frac{f c (p_x c_y - p_y c_x)}{f^2 + c^2 p^2} \right] + O(\mu^2)

ωR=μc(βcpx+2fpxcy2fpycx)f2+c2p2+O(μ2)\omega_R = \mu \frac{c ( - \beta c p_x + 2f p_x c_y - 2f p_y c_x)}{f^2 + c^2 p^2} + O(\mu^2)

となる(βyf\beta \equiv \partial_y fとした). それぞれ慣性重力波とRossby波の周波数に対応する.

  • パラメータf,Hf, Hの空間勾配を摂動として生じるRossby波の分散関係を正確に導くことができた.

  • β\beta効果により, 慣性重力波は波数ベクトルの方向によって分散関係が変化する.
    (このことは最近 Perez et al. でも議論されている.)

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

数値モデルデータの分析

  • 浅水系の数値シミュレーション.

  • 二重周期境界条件の正方形領域で, 水底の形状を右上の図に設定.

  • 初期の水面に変位を与えて, 波束の運動を観察する.


境界条件: 二重周期
格子点数: 128×128128 \times 128
格子幅\sim変形半径/22

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

基本的な出力データ

左から: u,v,ηu, v, \eta

UVE_anime

慣性重力波の速い伝搬と, 地形性Rossby波の遅い伝搬が見られる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

波動分離解析

射影作用素 u(x^,p^)\bm{u}(\hat{\bm{x}}, \hat{\bm{p}})を用いて波動成分を分離する.

  • 左: 慣性重力波の信号 ψG\psi_G
  • 右: Rossby波の信号 ψR\psi_R
    (上が実部, 下が虚部)

(Rossby波の信号を解析する際, ωR\omega_Rが正の部分のみを取り出すフィルタ操作を行っている. 詳細は論文参照)

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

エネルギー解析

  • Wigner分布関数w(x,p)w(\bm{x}, \bm{p})を計算, 物理空間に射影してエネルギー密度E(x)E(\bm{x})とフラックスF(x)\bm{F}(\bm{x})を得る. [左:慣性重力波, 右:Rossby波]

width:950px

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

エネルギー解析

  • ある瞬間のスナップショット.
    それぞれの波の群速度に比例したフラックスを描画することができた.

width:950px

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

全体まとめ

  • 流体波動を記述する偏微分方程式に対し, Wigner変換を用いることで, 見通しよく方程式の変形操作を行うことができる.

  • Wigner変換により得られた表象行列の対角化により, 連立偏微分方程式に含まれる独立した波動成分が準局所的に分離される.

  • それぞれの波動の信号からWigner分布関数を定義し, 物理空間へ射影することで, 群速度に比例したフラックスが得られる.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(補足)

一般的な振動モードの分離

  • 閉じた領域での固有値解析

    • 離散モードとして得られる.

    • 各モードは領域全体に広がった構造.
      (右: Rossby-Haurwitz 波)

  • 今回のテーマ

    • 無限に広い領域.

    • 空間に局在した波束.

GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(補足) Maximaを用いた表象行列の対角化計算

define(th(x,y,px,py),sqrt(f^2+(c+d*y)^2*(px^2+py^2)))$
define(p(px,py),sqrt(px^2+py^2))$
define(h0(x,y,px,py),matrix([0,%i*f,(c+d*y)*px],
 [-%i*f,0,(c+d*y)*py],
 [(c+d*y)*px,(c+d*y)*py,0]
))$
define(h1(x,y,px,py),matrix([0,0,0],
 [0,0,%i*d/2],
 [0,-%i*d/2,0]
))$
define(u0(x,y,px,py),
matrix(
[(th(x,y,px,py)*px - %i*f*py)/(sqrt(2)*th(x,y,px,py)*p(px,py)),
 (th(x,y,px,py)*py + %i*f*px)/(sqrt(2)*th(x,y,px,py)*p(px,py)),
 (c+d*y)*p(px,py)/(sqrt(2)*th(x,y,px,py))],
[(th(x,y,px,py)*px + %i*f*py)/(sqrt(2)*th(x,y,px,py)*p(px,py)),
 (th(x,y,px,py)*py - %i*f*px)/(sqrt(2)*th(x,y,px,py)*p(px,py)),
 -(c+d*y)*p(px,py)/(sqrt(2)*th(x,y,px,py))],
[%i*(c+d*y)*py/th(x,y,px,py),
 -%i*(c+d*y)*px/th(x,y,px,py),
 f/th(x,y,px,py)]
))$
define(uc0(x,y,px,py),
matrix(
[(th(x,y,px,py)*px + %i*f*py)/(sqrt(2)*th(x,y,px,py)*p(px,py)),
 (th(x,y,px,py)*px - %i*f*py)/(sqrt(2)*th(x,y,px,py)*p(px,py)),
 -%i*(c+d*y)*py/th(x,y,px,py)],
[(th(x,y,px,py)*py - %i*f*px)/(sqrt(2)*th(x,y,px,py)*p(px,py)),
 (th(x,y,px,py)*py + %i*f*px)/(sqrt(2)*th(x,y,px,py)*p(px,py)),
 %i*(c+d*y)*px/th(x,y,px,py)],
[ (c+d*y)*p(px,py)/(sqrt(2)*th(x,y,px,py)),
 -(c+d*y)*p(px,py)/(sqrt(2)*th(x,y,px,py)),
 f/th(x,y,px,py)]
))$
define(om(x,y,px,py),
matrix(
[th(x,y,px,py), 0 , 0],
[0, -th(x,y,px,py), 0],
[0, 0, 0]
))$
GFD オンラインセミナー: Wigner変換を用いた流体波動の局所分離解析

(続き)

define(du0dy(x,y,px,py)   , diff(u0(x,y,px,py)  , y ) )$
define(du0dpy(x,y,px,py)  , diff(u0(x,y,px,py)  , py) )$
define(duc0dy(x,y,px,py)  , diff(uc0(x,y,px,py) , y ) )$
define(duc0dpy(x,y,px,py) , diff(uc0(x,y,px,py) , py) )$
define(dh0dy(x,y,px,py)    , diff(h0(x,y,px,py)   , y ) )$
define(dh0dpy(x,y,px,py)   , diff(h0(x,y,px,py)   , py) )$
define(domdy(x,y,px,py)   , diff(om(x,y,px,py)   , y) )$
define(domdpy(x,y,px,py)  , diff(om(x,y,px,py)   , py) )$
define(uhu1(x,px,py),
u0(x,0,px,py).h1(x,0,px,py).uc0(x,0,px,py) +
%i/2*(du0dy(x,0,px,py).dh0dpy(x,0,px,py).uc0(x,0,px,py)
    - du0dpy(x,0,px,py).dh0dy(x,0,px,py).uc0(x,0,px,py)
    + u0(x,0,px,py).dh0dy(x,0,px,py).duc0dpy(x,0,px,py)
    - u0(x,0,px,py).dh0dpy(x,0,px,py).duc0dy(x,0,px,py)
    + du0dy(x,0,px,py).h0(x,0,px,py).duc0dpy(x,0,px,py)
    - du0dpy(x,0,px,py).h0(x,0,px,py).duc0dy(x,0,px,py) )
)$
define(uu1(x,px,py),
%i/2*(du0dy(x,0,px,py).duc0dpy(x,0,px,py)
                     - du0dpy(x,0,px,py).duc0dy(x,0,px,py))
)$
define(a(x,px,py),
uhu1(x,px,py)- om(x,0,px,py).uu1(x,px,py)
)$
factor(ratsimp(a(x,px,py)));
 
define(lm(x,px,py),matrix([0,0,0],[0,0,0],[0,0,0]))$
for k:1 thru 3 do
for l:1 thru 3 do
if k = l then lm(x,px,py)[k,l]:0
else lm(x,px,py)[k,l]:
 - 1/(om(x,0,px,py)[k,k] - om(x,0,px,py)[l,l]) * uhu1(x,px,py)[k,l]
 + 1/2 * (om(x,0,px,py)[k,k] + om(x,0,px,py)[l,l]) / (om(x,0,px,py)[k,k] - om(x,0,px,py)[l,l])
 * uu1(x,px,py)[k,l]$
define(u1(x,px,py),
(- uu1(x,px,py)/2 + lm(x,px,py)).u0(x,0,px,py)
)$
factor(ratsimp(u1(x,px,py)));

分散関係と偏光関係は, 振動モードの固有値と固有ベクトルに対応

- 基本場は定常(支配方程式は時間変化を含まない)としている

# 本題: 指数関数解とWKB近似のギャップを埋め, 空間的に変動する媒質を伝わる波動現象に対する理解を深める ## キーワード: Wigner変換, 擬微分作用素, 表象, スター積 ---

具体的には$\hat{x}\hat{p} - \hat{p}\hat{x} = {\rm i} \mu$ となる.

## 逆作用素の計算 (続き)

- $\hat{\bm{H}}$をハミルトニアンとすると, $(\star\star\star)$はSchrödinger方程式に一致.

なお, $\omega(\bm{x}, \bm{p})$を系のハミルトニアンとみなすとLiouville方程式に相当する.

- この時点では惑星$\beta$効果は表面的に見えない.

--- # Wigner変換を用いた作用素の分解と逆作用素の計算 作用素$\hat{\bm{A}}$をWigner変換して表象$\bm{a}(\bm{x}, \bm{p})$を求める - ここで, エルミート作用素に対応する表象はエルミート行列になる. 以下の方程式を満たす表象$\bm{l}(\bm{x}, \bm{p})$を求める: $$ \bm{a} = \bm{l} \star \bm{l}^\dag $$ 表象$\bm{a}$と$\bm{l}$が次のように漸近展開できるとする: $$ \bm{a} = \bm{a}_0 + \mu \bm{a}_1 + ...$$ --- 中立安定な線形波動は, 一般に以下の形式の方程式で記述される: $$\bm{n}(\hat{\bm{x}}, \hat{\bm{p}})\frac{\partial \bm{\psi}}{\partial t} = \bm{m}(\hat{\bm{x}}, \hat{\bm{p}}) \bm{\psi}$$ ここで$\bm{h}(\bm{x}, \bm{p})$はエルミート行列であり, ほとんどすべての$(\bm{x}, \bm{p}) \in \mathbb{R}^2$において正則であるとする. - エルミート行列を表象にもつ擬微分作用素は**エルミート作用素**である --- ## 擬微分作用素の ---