
% Mesinger and Arakawa
% 20110113 (担当者：川畑拓也)
% 1次元線形移流方程式を解く


\documentclass[a4j,12pt]{jreport}
\usepackage{Dennou6}
\usepackage{amsmath}
\usepackage{times}
\usepackage{url}
\usepackage{bm}
%\usepackage{booktabs}
%\usepackage{float}
\usepackage{ascmac}
\pagestyle{Dmyheadings}

\setlength{\parskip}{2ex}

\setcounter{chapter}{1}
\setcounter{section}{1}
\setcounter{subsection}{0}
\setcounter{figure}{0}
\setcounter{table}{0}
\setcounter{page}{6}
\setcounter{equation}{0}

\makeatletter  % 再定義する際に必要. documentclass の後ろに書く必要があ
       % る.
\renewcommand{\chapter}{%
  \thispagestyle{jpl@in}%
  \global\@topnum\z@
  \@afterindenttrue
  \secdef\@chapter\@schapter}
 % chapter の直前で改行しないように設定.
\makeatother

\def\thesection{付録 \arabic{chapter}}
\def\theequation{a. \arabic{equation}}
\def\thepage{\arabic{page}} % 数字
%\def\thepage{\roman{page}} % ローマ数字
\def\thetable{\arabic{chapter}.\arabic{section}.\arabic{table}}
\def\thefootnote{\arabic{footnote} )}
\def\thefigure{\arabic{section}.\arabic{figure}}

\Dtitle{GFDワークノート}{}
\Dauthor{川畑拓也}{}
\Ddate{2011/01/13}
\markright{第 10 章 1 次元線形移流方程式の数値解法}


\begin{document}
\section{1次元移流方程式の数値解法 (オイラースキーム)}
関数 $u=u(x,t)$ を考えると, 1次元線形移流方程式は \Deqref{adeq} で表される.
\begin{equation}
  \DP{u}{t} + c\DP{u}{x} = 0, \ \qquad \ c = \mathrm{const.} \Deqlab{adeq}
\end{equation}
ここでは, 時間方向に前進差分を, 空間方向に中央差分を用いて \Deqref{adeq} を数値
的に解く \footnote{この様な差分のとり方をするスキームを FTCS(Front Time Center 
Space) 法と言うこともある.
}. 
 ちなみに, \Deqref{adeq} の解析解は $u(x,t) = f(x-ct)$ という形をしている\footnote{
$u(x,t)= f(x-ct)$ が \Deqref{adeq} の解であることを確かめる. 実際に $f(x-ct)$ 
を\Deqref{adeq} に代入し合成関数の微分に注意して計算すれば,
\begin{equation*}
  \DP{u}{t} + c\DP{u}{x} = -c f^{\prime} (x-ct) + c f^{\prime}(x-ct) =0.
\end{equation*}
ゆえに, $f(x-ct)$ は\Deqref{adeq} の解である.
}. 

上記の差分法を用いて \Deqref{adeq} を差分化すると, 
\begin{equation*}
  u_j^{n+1} = u_j^n - c \frac{\Delta t}{2 \Delta x} (u_{j+1}^n - u_{j-1}^n)
\end{equation*}
となり, これが解くべき差分式である. 以下の図は計算結果を可視化したものである. 
格子間隔 $\Delta x$ は $\Delta x = 0.33$, 時間間隔 $\Delta t$ は $\Delta t = 0.01$ 
である. また, 境界条件は周期境界条件 $u_{-1} = u_{n}$ を用い, $c = 10$ とした. 


\begin{figure}[H]
  \begin{tabular}{cc}
    \begin{minipage}{0.5\hsize}
      \begin{center}
        \includegraphics[scale=0.3,angle=-90]{0001.ps}  
        \caption{初期状態 }
      \end{center}
      \label{t=1}
    \end{minipage}
    \begin{minipage}{0.5\hsize}
      \begin{center}
        \includegraphics[scale=0.3,angle=-90]{0010.ps}  
        \caption{1 次元移流方程式の解 (t = 0.1 のとき)}
      \end{center}
      \label{t = 100}
    \end{minipage}
  \end{tabular}
\end{figure}


\begin{figure}[H]
  \begin{tabular}{cc}
    \begin{minipage}{0.5\hsize}
      \begin{center}
        \includegraphics[scale=0.3,angle=-90]{0050.ps}  
        \caption{1 次元移流方程式の解 (t = 0.5 のとき) }
      \end{center}
      \label{t=500}
    \end{minipage}
    \begin{minipage}{0.5\hsize}
      \begin{center}
        \includegraphics[scale=0.3,angle=-90]{0100.ps}  
        \caption{1 次元移流方程式の解 (t= 1.0 のとき)}
      \end{center}
      \label{t = 1000}
    \end{minipage}
  \end{tabular}
\end{figure}

解析解であれば, $t=0$ での矩形が維持されたまま x 軸正方向へ伝播するはずである. 
しかし, 図を見てもらえばわかるように, 時間ステップの増加とともに波形は崩れている. 
これはフォンノイマンの安定性解析から考察すると, 次のような理由によるものだろうと
考えられる. すなわち, 解を$u_j^n = \mathrm{Re}[U(t)e^{ik\Delta x}]$ と仮定し, 
差分式に代入して整理すると,
\begin{equation*}
  U^{n+1} = \left( 1 - i c \frac{\Delta t}{\Delta x} \sin k \Delta x \right)U^n. 
\end{equation*}
増幅係数 $\lambda = U^{n+1} / U^n$ の大きさは, 
\begin{align*}
  \lvert \lambda \rvert &= \Dsqrt{\lambda \lambda^{\ast}} \\
  &= \Dsqrt{\left( 1- i c \frac{\Delta t}{\Delta x}\sin k \Delta x \right) \left( 1+ i c \frac{\Delta t}{\Delta x}\sin k \Delta x \right) }\\
  &= \Dsqrt{1 + \left( c\frac{\Delta t}{\Delta x} \right)^2 \sin^2 k \Delta x}
\end{align*}
ゆえに, 増幅係数が常に $\lvert \lambda \rvert > 1$ となるため不安定になると考えられる.
  

\end{document}
