Advection

簡単のため1次元の波で考える 時刻\(t\)における座標 \(x\) の値を \(f(x-ut)\) とする 例えば，正弦波が速度\(u\)でx軸正方向に進むなら，\(f(x-ut) = sin(x-ut)\) この時，\(\frac{\partial f}{\partial t} + u \frac{\partial f}{\partial x} =0 \)が成り立つ．これが移流方程式 1次風上差分法 直線で近似する 数値解析なので，\(x\)は\(\dots,x_{i-1},x_i,x_{i+1},\dots\)のように離散的 \(x = x_i\)における\(f\)の近似を\(F_i^n(x)\)とする この手法は直線で近似するので，\(F_i^n(x)=a(x-x_i) + f_i^n, a = \frac{f_i^n - f_{i-1}^n}{\Delta x}\) ただし，\(f_i^n\)は\(x = x_i\)において時間ステップを\(n\)回進めた時の値 1時間ステップあたり時間は\(\Delta t\)だけ進む 時間ステップを進めると，$f_i^{n+1} = F_i^n(x_i-u\Delta t) = a(-u \Delta t) + f_i^n= - \frac{u \Delta t}{\Delta x}(f_i^n - f_{i-1}^n) + f_i^n$ プログラムでは，$K=u \frac{\Delta t}{\Delta x}$として以下を繰り返すだけ 図は，矩形波を与えたもの 風下もある for(i=1; i<99; i++) f_new[i] = -K * (f[i]-f[i-1]) + f[i]; for(i=1; i<99; i++) f[i] = f_new[i]; Lax Wendroff法 二次関数で近似する $F_i^n(x) = a(x-x_i)^2 + b(x-x_i) + c$で近似する $a,b,c$3つの定数が必要なので，$x$に適当に値を入れて， $F_i^n(x_{i-1}) = a \Delta x ^2 - b \Delta x + c = f_{i-1}^n$ $F_i^n(x_{i}) = c = f_{i-1}^n$ $F_i^n(x_{i+1}) = a \Delta x ^2 + b \Delta x + c = f_{i+1}^n$ ただし，$x_{i+1} - x_{i} = x_{i} - x_{i-1}=\Delta x$ これを解いて， $$ F_{i}^n(x) = a(x-x_i)^2 + b(x-x_i) + f_i^n, \\ a = \frac{f_{i+1}^n - 2f_{i}^n + f_{i-1}^n}{2\Delta x^2}, b = \frac{f_{i+1}^n - f_{i-1}^n}{2\Delta x} $$ 時間ステップを進めると， $$ f_i^{n+1} = F_i^n(x_i-u\Delta t) = f_i^n - \frac{u\Delta t}{2 \Delta x}(f_{i+1}^n - f_{i-1}^n) + \frac{(u\Delta t)^2}{2 \Delta x ^2}(f_{i+1}^n -2 f_i^n- f_{i-1}^n) $$ プログラムでは，風上差分と同様に$K=u \frac{\Delta t}{\Delta x}$として以下を繰り返すだけ for(i=1; i<99; i++) f_new[i] = f[i] - K*K*(f[i+1]-f[i-1])/2....