<4D F736F F F696E74202D208D E9197BF288CF68A4A B8CDD8AB B83685D>

離散化手法とスキームの基礎

離散化手法とスキームの基礎

離散化手法とスキームの基礎

離散化手法とスキームの基礎

と選択法

と選択法

選択法

選択法

2007/1/16 宇宙航空研究開発機構 宇宙航空研究開発機構 情報・計算工学センター 嶋英志 嶋英志

本講習の目的

本講習の目的

本講習の目的

本講習の目的

基礎的な計算法の性質を述べ 各手法 y 基礎的な計算法の性質を述べ、各手法 の持つ長所短所を理解することによっ 手法 背景を理解 た 選択 て、手法の背景を理解した正しい選択 に近づくこと y 「クーラン数」「風上差分」等の広い 範囲のCFD技術に共通の概念について 範囲のCFD技術に共通の概念について 、その意味とイメージを把握すること

本講習の方針

本講習の方針

本講習の方針

本講習の方針

様々な流体方程式の基礎となる移流方 y 様々な流体方程式の基礎となる移流方 程式を用いて色々な計算法の特徴を計 算例を示しながら解説する 算例を示しながら解説する y 上記については十分理解できるように上記については十分理解できるように 丁寧に説明する 圧縮性と非圧縮性 差分法と有限体積 y 圧縮性と非圧縮性、差分法と有限体積 法と有限要素法、各々で扱いの違う非 線形性や多次 扱 は省略する 線形性や多次元の扱いは省略する

目次

目次

目次

目次

y 移流方程式とその離散化の基本 y 移流方程式とその離散化の基本 ◦ 移流方程式の位置付 差分法 有限体積法 有限要素法の関係 ◦ 差分法、有限体積法、有限要素法の関係 y スキームの安定性、振幅誤差、位相誤差 y 実用コードでの高度なスキーム ◦ 陰解法 ◦ 高次精度化と無振動化 x 風上法 x TVD法、ENO法、WENO法 x CIP法 その他の手法 y その他の手法

目次

目次

目次

目次

y 移流方程式とその離散化の基本 y 移流方程式とその離散化の基本 ◦ 移流方程式の位置付 ◦ 差分法、有限体積法、有限要素法の関係 y スキームの安定性、振幅誤差、位相誤差スキ ムの安定性、振幅誤差、位相誤差 y 実用コードでの高度なスキーム 陰解法 ◦ 陰解法 ◦ 高次精度化と無振動化 x 風上法 x TVD法、ENO法、WENO法 x CIP法

流体の方程式と移流方程式

流体の方程式と移流方程式

流体の方程式と移流方程式

流体の方程式と移流方程式

一次元線形移流方程式（_OLCE)

0

=

+

_x t

aq

q

q

=

f

(

x

−

at

)

線 ₎ 0 = + _x t E Q 一次元圧縮性_{Euler方程式}

0

=

+ q

a

q

係数行列の直行分解

0

=

+

_{k k x} t k

a

q

q

二次元_{(非圧縮）NS方程式} 高レイノル u p vu uu u_t + _x + _y − _x = + Δ Re 1 1

ρ

( 高レイノル ズ数では小

ρ

移流方程式の離散化（１）

移流方程式の離散化（１）

移流方程式の離散化（１）

移流方程式の離散化（１）

y 時空間の離散化 T(時間） y 時空間の離散化 Δt n n n q_in+1 q_i-1n+1 _q i+1n+1 時間発展 時間発展 x(空間) Δx q_in _q i+1n q_i-1n _時間発展 y 微分項の差分近似 x(空間) n i q x i t n q( Δ , Δ ) = 微分項の差分近似 x q q q_x i i Δ − ≈ + − 2 1 1 t q q q i n i n t _Δ − ≈ +1 x Δ 2 Δt

移流方程式の離散化（２）

移流方程式の離散化（２）

移流方程式の離散化（２）

移流方程式の離散化（２）

微分方程式の離散化 y 微分方程式の離散化 ◦ 空間中心差分+Euler陽解法 ) ( 6 2 4 2 1 1 _q x _{q O x} x q q xxx x i i = + Δ + Δ Δ − ₋ + _{：二次精度} 2 1 1 1 n i n i n i n i q q q q + = − λ + − − x t a Δ Δ = λ ：Courant数 2 Δx i i-1 i+1

移流方程式の離散化（３）

移流方程式の離散化（３）

移流方程式の離散化（３）

移流方程式の離散化（３）

微分方程式の離散化 y 微分方程式の離散化 ◦ 空間風上差分+Euler陽解法 ) ( 2 2 1 _q x _{q O x} x q q xx x i i = − Δ + Δ Δ − ₋ ：一次精度 i i-1 i+1

(

n

)

i n i n i n i q q q q +1 = − λ − ₋₁ 2 x Δ x q q q a x q q a x q q a i i i i i i i Δ + − − Δ − = Δ − _{− + − + − −} 2 2 2 1 1 1 1 1 1 2 2 2 1 1 1 1 1 n i n i n i n i n i n i n i q q q q q q q + = −λ + − − + λ + − + −

移流方程式の離散化（３）

移流方程式の離散化（３）

移流方程式の離散化（３）

移流方程式の離散化（３）

y 有限体積法と差分法の比較有限体積法と差分法の比較 fi 1/2 fi+1/2 ∫ Δ+ Δ ∫ + = + − t n t n x x qt aqx dxdt i i ) 1 ( 0 ) ( 2 / 1 2 / 1 時空間で積分 2 / 1 + i f i q q_i+1 1 − i q 2 / 1 − i f 0 ) ( 1 ) ( 1 2 / 1 2 / 1 1 _{− =} Δ + − Δ + − + i i n i n i q f f q 平均が正確なら下式は厳密 ) ( ) ( _{1/2 1/2} Δ Δ i i fi+ fi x q q t ∫ + − Δ = 1/2 2 / 1 1 i i x x qdx x q ∫ + Δ Δ Δ = n t t n aqdt t f 1 ( 1) 1 n n Δx i 1/2 Δt の近似次第で様々なスキームが定義可能 f 1 1 1 n i n i n n+ q ₊ − q ) ( 2 1 1 2 / 1 n i n i i a q q f ₊ = + ₊ 2 1 1 1 n _{i i} i n i q q q q + = − λ + − ：中心差分 ) ( 1 n i n i n i n i q q q q + = −λ − ：風上差分

{

( ) ( )

}

1 1 1 2 / 1 n n n n q q a q q a f =

{

( + )+ ( − )

}

qi qi λ(qi qi ) 2 1 1 2 / 1 i i i i i a q q a q q f ₊ + ₊ + ₊

移流方程式の離散化（４）

移流方程式の離散化（４）

移流方程式の離散化（４）

移流方程式の離散化（４）

y 有限要素法と差分法の比較 0 ) ( + = ∫ W_i q_t aq_x dx この式は厳密。局所的にのみ値を持つ W_iを用いて局所的性質を抽出する。 仮定１ 頂点の間を線形補完(例）で近似 1 1 + + + = q_iN_i q_i N_i q 仮定１：頂点の間を線形補完(例）で近似 仮定２：W_iにN_iと同じものを使用=ガラーキン法 ( ) ( ) ( ) ( 4 ) ( ) 0 6 1 1 1 1 1 + + + + Δ + − − = − t i t i i i i t q q x a q q q 仮定 _{i i} 同 を使用 ラ キ 法 t q q q i n i n t _Δ − ≈ +1 ( ) ( ) ( ) ( q_{t i−1} +4 q_{t i} + q_{t i+1})⇒( )q_{t i} 6 1 _{：質量集中行列} 2 1 1 1 _{+ −} + ₌ n _{− i} − _i i n i q q q q λ ：中心差分と一致

FDM FEM FVM

FDM FEM FVMの比較

の比較

FDM,FEM,FVM

FDM,FEM,FVMの比較

の比較

y FEM FVMのメリットは計算メッシュ y FEM、FVMのメリットは計算メッシュ に対する柔軟性 差分法用の格子 FVM、FEM用の格子

基本的離散化のまとめ

基本的離散化のまとめ

基本的離散化のまとめ

基本的離散化のまとめ

次元の簡単なスキ ムについては y 一次元の簡単なスキームについては、 差分法(FDM)、有限体積法(FVM)、有 限要素法 は結局同じようなもの 限要素法(FEM)は結局同じようなもの y したがって、今後は原則的にしたがって、今後は原則的にFDMを用を用 いて基礎的なスキームを説明する ただし より洗練されたスキ ムでは y ただし、より洗練されたスキームでは 相互の変換は一般には不可能 y FDM,FVM,FEM間で優れたスキームの相 互移植は研究テーマになりうる 互移植は研究テ マになりうる

目次

目次

目次

目次

y 移流方程式とその離散化の基本 y 移流方程式とその離散化の基本 ◦ 移流方程式の位置付 ◦ 差分法、有限体積法、有限要素法の関係 y スキームの安定性、振幅誤差、位相誤差スキ ムの安定性、振幅誤差、位相誤差 y 実用コードでの高度なスキーム 陰解法 ◦ 陰解法 ◦ 高次精度化と無振動化 x 風上法 x TVD法、ENO法、WENO法 x CIP法

スキームの安定性

スキームの安定性

(1)

(1)

スキームの安定性

スキームの安定性

(1)

(1)

y 周期的2重サイン波の伝搬問題 y 周期的2重サイン波の伝搬問題 1 1 1 n i n i n n+ λ q₊ −q +1 ( ) 中心差分+陽解法 _{風上差分+陽解法} 2 1 1 1 n i i i n i q q q q + = −λ + − q_in+1 =q_in −λ(q_in −q_i−1n) 無条件不安定 条件安定

スキームの安定性

スキームの安定性

(2)

(2)（矩形波）

（矩形波）

スキームの安定性

スキームの安定性

(2)

(2)（矩形波）

（矩形波）

中心差分+陽解法 ：無条件不安定 風上差分+陽解法 ：条件安定 矩形波やステップは最高波数を含むので高波数の誤差が顕著にわかる

スキームの安定性

スキームの安定性

(3)

(3)（安定性解析）

（安定性解析）

スキームの安定性

スキームの安定性

(3)

(3)（安定性解析）

（安定性解析）

y フォン・ノイマンの安定性解析 y フォン・ノイマンの安定性解析 ◦ サイン波の一つの波数に着目 ｋ：波数 0≦k≦ ) exp( jik x q_in = Δ j:虚数単位 ◦ ｋ：波数 0≦k≦π ◦ kΔx=πが扱える最高波数 Δx

スキームの安定性

スキームの安定性

(3)

(3)（安定性解析）

（安定性解析）

スキームの安定性

スキームの安定性

(3)

(3)（安定性解析）

（安定性解析）

y フォン・ノイマンの安定性解析 y フォン・ノイマンの安定性解析 ) exp( jik x q_in = Δ ◦ 厳密解 n n q A q_i +1 = A_exactq_i q = x k j x k x jk t jak A_exact Δ Δ = Δ = Δ − = λ λ λ ) cos sin exp( ) exp( ◦ 差分近似の_{A(増幅係数)} x k j x k x jk Δ = Δ − Δ −

= exp( λ ) cos λ sin λ

x k j A_c2 =1− λ sin Δ 中心差分+陽解法 常に|A_C2|>1 x k j x k A_u1 =1− λ(1− cos Δ ) − λ sin Δ 風上差分+陽解法 λ≦1で|A_U1|≦1

増幅係数からわかる振幅誤差

位相誤差

増幅係数からわかる振幅誤差

位相誤差

増幅係数からわかる振幅誤差、位相誤差

増幅係数からわかる振幅誤差、位相誤差

◦ 厳密解の厳密解のAA ) exp( )

exp( _{exact exact}

exact jk x r j A = − λΔ = θ 1 = r _tan(θ ₎ = _tan(−_kλΔ_x) ◦ 差分近似の_A k j A 1 λ i Δ 中心差分+陽解法 1 exact

r _tan(θ_exact) = _{tan( k}λΔ_x)

x k j A_c2 =1− λ sin Δ x k j x k A_u1 =1−λ(1− cos Δ ) − λ sin Δ 中心差分+陽解法 風上差分+陽解法 1 1.2 |A| Cn=0.4 1 5 2 t)) tan(θ)/tan(θexact)) 0.2 0.4 0.6 0.8 |A| FTCS UPWIND 0 0.5 1 1.5 θ)/tan( θexac t FTCS UPWIND 0 0 1 2 3 4 kΔx 0 0 1 2 3 4 tan( kΔｘ

時間積分法の変更による中心差分の安定化 時間積分法の変更による中心差分の安定化 時間積分法の変更による中心差分の安定化 時間積分法の変更による中心差分の安定化 y 予測子 修正子法（時間一次精度） y 予測子-修正子法（時間一次精度） 2 1 1 n i n i n i i q q q q = − λ + − − 2 2 1 1 1 + − + ₌ n _{− i} − _i i n i q q q q λ y 3段階ルンゲ・クッタ（時間三次精度） 2 2 3 1 _{i 1}n _{i 1}n n i i q q q q = − λ + − − 1 q q 2 2 1 ₊₁ − ₋₁ − = n i i i i q q q q λ 1 1 1 + n _i − _i n q q 2 1 1 1 _{+ −} + ₌ n _{− i i} i n i q q q q λ

予測子 予測子 修正子法でのクーラン数の効果修正子法でのクーラン数の効果 予測子 予測子--修正子法でのクーラン数の効果修正子法でのクーラン数の効果 C_N=0.1 C_N=0.5 C_N=0.9 C_N=1.1

数値例：スキームによる矩形波伝搬の違い 数値例：スキームによる矩形波伝搬の違い

中心差分+陽解法 中心差分+予測子 修正子法

中心差分+陽解法 中心差分+予測子-修正子法

数値例：振幅誤差（散逸性）と

数値例：振幅誤差（散逸性）と

位相誤差（分散性）

位相誤差（分散性）

中心差分+三段階R-K 風上差分+陽解法 高波数で位相誤差が顕著 高波数で位相誤差もあるが それより減衰が顕著

簡単なスキームからの改良方向

簡単なスキームからの改良方向

簡単なスキームからの改良方向

簡単なスキームからの改良方向

y ク ラン数の制限を外したい y クーラン数の制限を外したい ◦ 陰解法 y 高波数まで正しく解きたい ◦ ⇒高次精度法⇒高次精度法 y 振動を無くしたい 波数による移送速度の違いが振動の原因 ◦ 波数による移送速度の違いが振動の原因 ◦ 正しくない波を消す必要がある ◦ 現実の問題では振動により破綻すること が少なくない少 ◦ ⇒TVD、ENO、WENO、CIP

目次

目次

目次

目次

y 移流方程式とその離散化の基本 y 移流方程式とその離散化の基本 ◦ 移流方程式の位置付 ◦ 差分法、有限体積法、有限要素法の関係 y スキームの安定性、振幅誤差、位相誤差スキ ムの安定性、振幅誤差、位相誤差 y 実用コードでの高度なスキーム 陰解法 ◦ 陰解法 ◦ 高次精度化と無振動化 x 風上法 x TVD法、ENO法、WENO法 x CIP法

陰解法

陰解法

陰解法

陰解法

y 中心差分+Euler陰解法 1 1 2 1 1 1 1 1 + − + + + ₌ n _{− i} n − _i n i n i q q q q λ n n n n+1 +1 λ +1 λ n i n i n i n i q q q q ₊₁ +1 + +1 − ₋₁ +1 = 2 2 λ λ A 1 ：連立一次方程式を解く y クランク・ニコルソン（時間二次精度） x k j A Δ + = sin 1 λ ：無条件安定 y クランク・ニコルソン（時間二次精度） ⎟⎟ ⎞ ⎜⎜ ⎛ ₋ + − − = + − + + − + +1 _{1 1 1} 1 ₁ 1 n i n i n i n i n i n i q q q q q q λ _⎟ ⎠ ⎜ ⎝ 2 2 2 i i q q x k j A Δ − sin 2 1 1 λ 無条件安定（|A|≡1） x k j A Δ + = sin 2 1 1 2 λ ：無条件安定（|A|≡1）

数値例（陰解法）

数値例（陰解法）

数値例（陰解法）

数値例（陰解法）

中心差分+Euler陰解法 Cn=0.5 中心差分+Euler陰解法 Cn=2.0 中心差分+Euler陰解法 Cn=4.0

中心差分+C-N Cn=0.5 中心差分+C-N Cn=2.0 中心差分+C-N Cn=4.0

C

C--N

N法の振幅誤差、位相誤差

法の振幅誤差、位相誤差

((クーラン数

クーラン数

=0.4

=0.4の場合

の場合

))

|A| C =0 4 tan(θ)/tan(θexact)) 振幅誤差 位相誤差 0 8 1 1.2 |A| Cn=0.4 1.4 1.6 1.8 2 θexact)) tan(θ)/tan(θexact)) 0.2 0.4 0.6 0.8 |A| FTCS UPWIND C-N 0 2 0.4 0.6 0.8 1 1.2 tan( θ)/tan( θ FTCS UPWIND C-N 0 0 1 2 3 4 kΔx 0 0.2 0 1 2 3 4 kΔｘ 陰解法でクーラン数が大きくとれても、時間発展が正確に計算できる訳ではない

δδフォームの陰解法とその改良

フォームの陰解法とその改良

δδフォームの陰解法とその改良

フォームの陰解法とその改良

0 ) ( 1 = + Δ − + Q L t Q Qn n _{差分近似 Lは何かのオペレーター} Δt 0 ) (Q = L 定常解の条件 n n Q Q Q = +1 − δ 1 1 ∂L δQの導入 δフォームの導入 0 ) ( 1 ) ( 1 ₌ ∂ ∂ + + Δ ≈ + + Δ Q Q L Q L Q t Q Q L Q t δ δ δ δ L I ⎞ ⎛ ∂ ) (Q L Q Q L t I _{= −} ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ∂ ∂ + Δ δ Δtによらず定常解は不変∴時間精度は定常解に関係ない ) (Q L Q Q L t I _{= −} ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ∂ ∂ + Δ δ ) 右辺を近似しても定常解は不変 Q t _⎠ ⎝ Δ ∂ 実用的な陰解法はこのような手法を用いている

δδフォームの具体例

フォームの具体例

δδフォームの具体例

フォームの具体例

2 1 1 1 1 1 + − + + + ₌ n _{− i} n − _i n i n i q q q q λ 中心差分+陰解法： 2 i i q q 2 2 1 1 1 1 n i n i i i i q q q q q +λ δ + −δ − = −λ + − − δ n n q q 中 分 陰解法 δフォームの導入 左辺を風上差分に ( ) 2 1 1 1 i i i i i q q q q q + λ δ −δ ₋ = −λ + − − δ ( ) 2 1 1 1 1 n i n i i i q q q q − ₋ = − + − − + λ δ λδ λ を風 分 整理 ( ) 2 1 i i q q 一方向スィープで解け、メリット大 中心差分+Euler陰解法 Cn=2.0 中心差分+Euler風上陰解法 Cn=2.0

陰解法

陰解法 まとめと注意

まとめと注意

陰解法

陰解法 まとめと注意

まとめと注意

y 陰解法はクーラン数を大きく取れ 特に定 y 陰解法はクーラン数を大きく取れ、特に定 常解を求めるのには有利 δフォ ムの陰解法では左辺側に近似形式 y δフォームの陰解法では左辺側に近似形式 の可能 数を きく れ が 時 時 y クーラン数を大きく取れるが、同時に時間 精度も急速に低下することに注意（特に時 間 次精度 場合） 間一次精度の場合） y 一次元線形の場合には適用容易だが、非線 形性や多次元の場合の適用法は自明ではな く、様々な工夫が必要

目次

目次

目次

目次

y 移流方程式とその離散化の基本 y 移流方程式とその離散化の基本 ◦ 移流方程式の位置付 ◦ 差分法、有限体積法、有限要素法の関係 y スキームの安定性、振幅誤差、位相誤差スキ ムの安定性、振幅誤差、位相誤差 y 実用コードでの高度なスキーム 陰解法 ◦ 陰解法 ◦ 高次精度化と無振動化 x 風上法 x TVD法、ENO法、WENO法 x CIP法

高波数の減衰の必要性

高波数の減衰の必要性

高波数の減衰の必要性

高波数の減衰の必要性

y バーガース方程式バ ガ ス方程式(非線形）の場合(非線形）の場合 0 = + _x t uu u x un = sin 1 x x x unt sin 2 2 1 cos sin = − − = x t x un+1 ≈ sin x−Δt 1 sin2x 非線形性から高波数が生じる u sin2 2 sin Δ ≈ 中心差分+予測子修正子法 風上差分+陽解法 非線形性から高波数が生じる 中心差分+予測子修正子法 風上差分+陽解法

分子粘性および数値粘性

分子粘性および数値粘性

分子粘性および数値粘性

分子粘性および数値粘性

◦ 一次元移流拡散方程式（粘性流体）次元移流拡散方 （粘 流体） xx x t aq q q + =ν x セルレイノルズ数 が十分小ならば安定化 x DNSはこの流儀に近い ν x aΔ DNSはこの流儀に近い x 高レイノルズ数の実用問題では望み薄 ◦ 数値粘性の付加 ◦ 数値粘性の付加 x 必要最小限度だけ粘性項を加える 数値粘性係数をどのように決めるかが問題 xxxx a xx a x t aq q q q + =ν (2) −ν (4) x 数値粘性係数をどのように決めるかが問題

風上差分の効果

風上差分の効果

風上差分の効果

風上差分の効果

風上差分に含まれる偶数次の誤差項が数 y 風上差分に含まれる偶数次の誤差項が数 値散逸として働く ◦ 一次風上法 ) ( 3 2 1 x x _O q q_{i i} Δ Δ Δ − ◦ 三次風上法 ) ( 6 2 3 1 _{q q q O x} x q q xxx xx x i i = − + + Δ Δ − ◦ 三次風上法 ) ( 12 6 6 3 2 _{1 1 2} 3 ₄ x O q x q q q q q xxxx x i i i i = + Δ + Δ Δ + − + _{− −} + y だが それだけでは十分ではない・・・ ) ( 12 6Δx qx qxxxx y だが、それだけでは十分ではない

高次精度無振動スキーム

高次精度無振動スキーム

高次精度無振動スキ

ム

高次精度無振動スキ

ム

y ゴドノフの定理 ◦ 「単調な線形スキ ムは高々一次精度である」 ◦ 「単調な線形スキームは高々一次精度である」 ◦ 一見、無振動高次精度スキームは不可能にみえ る る。 ◦ しかし「線形」の制限を外せば可能 3次風上+2段階R-K(線形スキーム) _{3次精度TVD(C-O)+2段階R-K(非線形スキーム)}

目次

目次

目次

目次

y 移流方程式とその離散化の基本 y 移流方程式とその離散化の基本 ◦ 移流方程式の位置付 ◦ 差分法、有限体積法、有限要素法の関係 y スキームの安定性、振幅誤差、位相誤差スキ ムの安定性、振幅誤差、位相誤差 y 実用コードでの高度なスキーム 陰解法 ◦ 陰解法 ◦ 高次精度化と無振動化 x 風上法 x TVD法、ENO法、WENO法 x CIP法

TVD

TVD法

法

TVD

TVD法

法

y TV安定(Total Variation Diminishing) y TV安定(Total Variation Diminishing)

n i n i q q ≤ Δ Δ + ∑ ∑ 1 n Δ 次ステ _Δ n y 十分条件 n i n i n i q q q = − Δ+ ₊ 1 Δq ni−1 i n q Δ 1 + Δq ni 1 − Δq ni i n q Δ 1 + Δq ni 次ステ ップ y 十分条件 1 2 / 1 2 / 1 1 − + − + + + Δ − Δ + = n _{i i i i} i n i q c q d q q 1 0 , 0 _1/2 2 / 1 ≥ − ≥ + i i d c 1 2 / 1 2 / 1 + ≤ + i i d c n i i n i i n i i i n i c d q c q d q q +1 = (1− _+1/2 − _−1/2) + _{+1/2 +1} + _{−1/2 −1} 1 − i n q i+1 n q i n q +1 y 一次精度風上法はTVDである 2 / 1 2 / 1 − + i i i n q q q

三次風上法の

三次風上法の

TVD

TVD化

化

三次風上法の

三次風上法の

TVD

TVD化

化

y 三次風上法のFVM表記 ) ( _{1/2 1/2} 1 − + + ₋ Δ Δ − = n _{i i} i n i f f x t q q 三次風上法のFVM表記 制限関数の導入 ( ) ( ) ⎭ ⎬ ⎫ ⎩ ⎨ ⎧ _{+ + Δ + − Δ} = + + ₋ +1/2 1 1 4 1 1 4 1 i i i i a q q q f φ φ ：Φ=1/3 y 制限関数の導入 ◦ Chakravarthy-Osher法 ⎫ ⎧ _{( ) ( )} ⎭ ⎬ ⎫ ⎩ ⎨ ⎧ _{+ + Δ + − Δ} = + + ₋ +1/2 1 ~ ~ 1 4 1 ~ 1 4 1 i i i i a q q q f φ φ ) ( minmod ~ q q q + + + _{= Δ Δ} Δ q_i−1 = minmod(Δ q_i−1,βΔ q_i) Δ~~+q = minmod(Δ+q βΔ+q ) Δ β Δ q_i = minmod(Δ q_i,βΔ q_i−1) ) , ( minmod < ≤ ⎪ ⎨ ⎧ = b a b a b a b a 1≤ β ≤ 3−φ 非線形性の導入 0 0 ) , ( minmod ≤ < ⎪ ⎩ ⎨ ab a b b b a φ β − ≤ ≤ 1 1

制限関数は何をしているのか？

制限関数は何をしているのか？

制限関数は何をしているのか？

制限関数は何をしているのか？

) ( _{1/2 1/2} 1 + − + ₋ Δ Δ + = n _{i i} i n i f f x t q q _{( )} 2 1 2 1 1 1 − − + _{Δ +Δ} Δ = Δ − ≈ i i _{i i} x q q x x q q q xx x n i xq x q q q = + +φ 2 2 2Δx Δx ) ( 1 2 1 2 2 1 1 − − + _{Δ −Δ} Δ = Δ + − ≈ i i _{i i} xx q q x x q q q q 2 / 1 + i f 1/2 ( ) (1 ) 1 4 1 1 4 1 − + + + = i + + Δ i + − Δ i i q q q q φ φ ≤ q q n i q ( ) qi ( ) qi qi + − + + _{+ − Δ ≤ Δ} Δ + 1 ₁ 4 1 1 4 1 _{φ φ} 1 2 / 1 + + ≤ i i q q 2 / 1 − i f 4 4 i i q q₋ + + _Δ − − ≤ Δ φ φ 1 3 1 ) 1 3 , min( ~ 1 1 1 i i i i q q q q₋ + ₋ + ₋ + + _Δ − − Δ = Δ ⇒ Δ φ φ ) 3 mod( min ~ q q q + + + _Δ − _Δ Δ φ ) 1 , mod( min ₁ 1 i i i q q q_{− −} Δ − Δ = Δ φ

TVD

TVD類似の手法：

類似の手法：

ENO WENO

ENO WENO

TVD

TVD類似の手法：

類似の手法：

ENO,WENO

ENO,WENO

y Essensially Non Oscilating Sheme y Essensially Non-Oscilating Sheme

2 / 1 + i f xx x n i xq x q q q = + +φ 2 n i q 2 / 1 − i f _{一次、二次の微係数から判断し最も変化の少ない} =滑らかな内挿関数を使う（三次精度の場合） ⇒ENO ⇒ENO 滑らかな場合にはすべての関数の重み付き 平均を使う。この場合最大5次精度 ⇒Weighted ENO

目次

目次

目次

目次

y 移流方程式とその離散化の基本 y 移流方程式とその離散化の基本 ◦ 移流方程式の位置付 ◦ 差分法、有限体積法、有限要素法の関係 y スキームの安定性、振幅誤差、位相誤差スキ ムの安定性、振幅誤差、位相誤差 y 実用コードでの高度なスキーム 陰解法 ◦ 陰解法 ◦ 高次精度化と無振動化 x 風上法 x TVD法、ENO法、WENO法 x CIP法

CIP

CIP法（１）

法（１）

CIP

CIP法（１）

法（１）

y Cubic Interpolated Polynomial y Cubic Interpolated Polynomial

0 = + _x t aq q q = g(x − at) ：移流方程式の解 q ↓この関係を使う 区分的3次関数 移流速度 N+1ステップの値と 微係数の取得 aΔt x ) , ( ) , ) 1 (( 1 t a x i t n g x i t n g qn+ i = + Δ Δ = Δ Δ − Δ

CIP

CIP（２）

（２）

CIP

CIP（２）

（２）

q 1 1 + + i n x i n q q 区分的3次関数 n i n q x q 区分的3次関数 i x q x 0 = + _x t aq q qn+1i = g((n +1)Δt,iΔx) = g(nΔt,iΔx − aΔt) 0 ) ( ) (q_{x t} + a q_{x x} = q 1i g'((n 1) t,i x) g'(n t,i x a t) n x = + Δ Δ = Δ Δ − Δ +

CIP

CIP（３）：非線形の場合

（３）：非線形の場合

CIP

CIP（３）：非線形の場合

（３）：非線形の場合

y 非線形輸送方程式 y 非線形輸送方程式

( )

= 0 + _x t ρu ρ

( )

( )

{

}

0 1 2 / 1 2 / 1 1 = − − = + + i i i n i n u u ρ ρ ρ ρ FVM的方法

( )

( )

{

_{1/2 1/2}

}

Δ i+ i− x ρ ρ ρ ρ CIP法法 ρ ρ ρ_t + u _x = −u_x

( )

ρ_{x t} + u

( )

ρ_{x x} = −u_xxρ − 2u_xρ_x 方程式 ρ ρ_t = −u_x

( )

ρ_{x t} = −u_xxρ − 2u_xρ_x ①非移流ステップ 0

( )

+

( )

0 0 = + _x t uρ ρ

( )

ρ_{x t} + u

( )

ρ_{x x} = 0 ②移流ステップ

数値例：ステップの伝搬

数値例：ステップの伝搬

数値例：ステップの伝搬

数値例：ステップの伝搬

一次風上 C-O3次+RK2 ENO3次+RK2 CIP ENO3次+RK2 CIP

数値例：二重サイン波の伝搬

数値例：二重サイン波の伝搬

数値例：二重サイン波の伝搬

数値例：二重サイン波の伝搬

一次風上 C-O3次+RK2 ENO3次+RK2 CIP ENO3次+RK2 CIP

高次精度無振動スキームのまとめ

高次精度無振動スキームのまとめ

と注意

と注意

y TVD ENO CIP等の非線形スキ ム y TVD、ENO、CIP等の非線形スキーム を用いることで高次精度化と無振動化 の両立が可能 ある の両立が可能である y 高波数の波の伝搬を扱うには高波数の波の伝搬を扱うにはENOやOや CIPの方が良い y 定常解に関しては２次精度以上のTVD は十分な精度を持っている

ご清聴ありがとうございました

ご清聴ありがとうございました

ご清聴ありがとうございました

ご清聴ありがとうございました

その他

その他

その他

その他

コンパクト差分（陰的差分）

コンパクト差分（陰的差分）

コンパクト差分（陰的差分）

コンパクト差分（陰的差分）

fとf’に関するTaylor展開 ) ( ) ( 6 1 ) ( 2 1 ) (l x f l x 2 f l x 3 f O x4 f f_i+l = + Δ ′ + Δ ′′ + Δ ′ ′′ + Δ ) ( ) ( 6 1 ) ( 2 1 ) ( 6 2 4 ) 4 ( 3 2 x O f x l f x l f x l f f_i′_+l = ′ + Δ ′′ + Δ ′ ′′ + Δ + Δ 6 2 より、例えば、次の関係が成立 ) ( ) ( 3 2 ) ( 12 1 ₄ 1 1 2 2 f f O x x f f x f_{i i i i} − _i + Δ Δ + − Δ − = ′ _{+ − + −} 4次精度差分公式

(

)

( ) 2 1 6 1 3 2 6 1 ₄ 1 1 1 1 f f O x x f f f_{i i i i} − _i + Δ Δ = ′ + ′ + ′_{− + + −} 陰的差分公式： 3 12Δx Δx 2 6 3 6 Δx 連立一次方程式を解くことでf’が求まる

コンパクト差分の評価（１）

コンパクト差分の評価（１）

コンパクト差分の評価（１）

コンパクト差分の評価（１）

fが三角関数であらわされる場合 ) exp( ) exp( jkx jik x f_i = = Δ jkf f ′ 安定性解析とは異なり時間積分は考えていないことに注意 2次精度中心差分の場合 jkf f = 2次精度中心差分の場合

(

)

f jk f x x jk x jk f f x i 1 i 1 2 e2 ) exp( ) exp( 2 1 = Δ Δ − − Δ = − Δ + − x k k_e Δ Δ = sin 2 x x 2 2Δ Δ 有効波数 x e Δ 2 x k x kΔ 1 sin 2 Δ sin 4 4次精度中心差分の場合 x x k x x k k_e Δ Δ − Δ Δ = sin 2 3 1 sin 3 4 4

コンパクト差分の評価（２）

コンパクト差分の評価（２）

コンパクト差分の評価（２）

コンパクト差分の評価（２）

⎞ ⎛ 4次精度陰的差分公式： f x x k j f x jk x jk jk_i Δ Δ = ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ _{− Δ + + Δ} sin ) exp( 6 1 3 2 ) exp( 6 1 4 sin 3 4 x k k_i = Δ ) cos 2 ( 4 x k k i Δ + 有効波数(L/Lexact)kΔxの比較 2 2.5 3 3.5 厳密 0.5 1 1.5 2 厳密 2次 4次 陰的4次 0 0 0.5 1 1.5 2 2.5 3 3.5