% 表題  移流方程式の差分解法: FCT 法  オリジナルFCTの問題点
%
% 履歴  1998/01/19 小高正嗣
%       1998/05/19 小高正嗣
%       2006/02/13 小高正嗣: ps ファイルへの相対パスを修正
%
\section{いくつかの問題}

    オリジナルの FCT には, 種々のスキームに対する応用面, 補正フラック
    スの与え方に関して問題点がある( Zalesak, 1979 ). 以下順を追ってこ
    れらを見ていくことにする.

    \subsection{応用性}

    Boris and Book(1973) では\Deqref{eq: 4}または\Deqref{eq: 5}のスキー
    ムに対して FCT を施したが, 拡散を含むような他のスキームに対しても
    応用することができる. 例えば上流差分に FCT を応用すると以下のよう
    な手順で計算を行なう(Book {\it et al.}, 1975).

    \begin{enumerate}
      \item 普通に上流差分で計算する.
        \[
          \tilde{\rho }_{i} = \rho _{i}^{n}-\frac{\Delta t}{\Delta x}
          [\rho _{i}^{n}u_{i}^{n}-\rho _{i-1}^{n}u_{i-1}^{n}].
        \]

      \item 各格子点間で差分を計算する.
        \[
          \Delta _{i+\frac{1}{2}} = \tilde{\rho }_{i+1}-\tilde{\rho }_{i}.
        \]

      \item 補正フラックスを計算する.
        \[
        f^{c}_{i+\frac{1}{2}} = S_{i+\frac{1}{2}} \mbox{max} \{0,
        \mbox{min}(S_{i+\frac{1}{2}}\Delta _{i-\frac{1}{2}}, \;
        \eta _{i+\frac{1}{2}}|\Delta _{i+\frac{1}{2}}|, \; 
        S_{i+\frac{1}{2}}\Delta _{i+\frac{3}{2}})\},
        \]
        \[
        S = \left\{
          \begin{array}{lcl} 
            +1 & \mbox{if} & \Delta _{i+\frac{1}{2}} \geq 0, \\
            -1 & \mbox{if} & \Delta f_{i+\frac{1}{2}} < 0.
          \end{array}
        \right. 
        \]
        ここで,
        \[
        \eta _{i+\frac{1}{2}}=\frac{1}{2}\epsilon _{i+\frac{1}{2}}
        (1-\epsilon _{i+\frac{1}{2}}),
        \]
        \[
        \epsilon _{i+\frac{1}{2}} = \frac{\Delta t}{\Delta x}
        \mbox{max}[u_{i}^{n}, u_{i+1}^{n}].
        \]

      \item 補正を行なう.
        \[
        \rho _{i}^{n+1}=\tilde{\rho }_{i} - (f_{i+\frac{1}{2}}^{c}-
        f_{i-\frac{1}{2}}^{c}).
        \]
    \end{enumerate}

    この場合 antidiffusive フラックスが $\eta _{i+\frac{1}{2}}|\Delta _
    {i+\frac{1}{2}}|$ となっているだけで, 手順は先に示した SHASTA の場
    合と同様である. つまり用いるスキームの antidiffusive フラックスが
    わかれば, FCT は簡単にそのスキームに適用することができる.

    しかし当然のことながらスキームごとに antidiffusive フラックスは異
    なる. よって適用するスキームを変更すると一から考察しなおさなければ
    ならない. これは結構面倒である. 実際に計算するときには一からコード
    を組み直さなければならなくなるからである. 何らかのより一般的な FCT
    の定義法があれば, スキーム全体の見通しの点や実際に計算を行なう上で
    便利である.

    \subsection{補正フラックスの問題(1): clipping}

    Boris and Book(1973) で与えられた補正フラックスは\Deqref{eq:
    10}$\sim $ \Deqref{eq: 12}式で与えられる. 先に解説したようにこれは
    新たに正負の極値を作ったり, すでに存在する極値をさらに強めないよう
    な工夫がなされている. 

    しかしこの工夫が災いして, 本当はピークを持つような分布が保存されな
    いという欠点が生じてしまう. これは ``clipping'' と呼ばれている.
    具体的に\Dfigref{fig: 5}を見ながら考えることにする. 簡単のため速度
    は 0 とし数値拡散だけを考慮する. 最初に(a)のような分布であったとす
    ると拡散のため(b)のようになる. これに対して\Deqref{eq: 10}$\sim $
    \Deqref{eq: 12}を基に補正フラックスを与えても, $i$ 点に存在したピー
    クは再現されない. 一方 $i-1$, $i+1$ 点での値は補正を受けるので $i$
    点の値に次第に近付く. やがて3点はほぼ同じ値をとるようになる
    ((c)). こうなると拡散も補正も効かなくなってしまいこの状態が維持さ
    れ, あたかも当初存在した山をはさみで切り取ってしまったかのようにな
    る.

    この問題は補正フラックスを決める際に1ステップ前の $\Delta _
    {i-\frac{1}{2}}, \Delta _{i+\frac{1}{2}}$ を参照することである程度
    解消されるが, 基本的には回避できない. しかし鋭いピークを持つような
    成分が重要かどうかは計算の対象となる現象に依存するので, clipping
    が問題とならない場合もある.
    
    \begin{figure}[b]
      \begin{center}
        \Depsf[][]{ps-fig/fig5.ps}
        \caption{``clipping'' の模式図(Book {\it et al.}(1973),図7を基に作成)}
        \Dfiglab{fig: 5}
      \end{center}
    \end{figure}

    \subsection{補正フラックスの問題(2): 多次元への拡張性}

    多次元に拡張した場合に\Deqref{eq: 10}$\sim $ \Deqref{eq: 12}式をそ
    のまま応用すると clipping とは別の問題が発生する. ここでは 2次元の
    場合についてその問題を考えることにする. 

    今1ステップの輸送過程の計算が終了して\Dfigref{fig: 6}のような分布
    が得られたとする.  このとき $(x,y)=(i\Delta x, j\Delta y)$ に働く
    antidiffusionフラックスの向きは通常は図の矢印の向きになる. $y$方向
    のantidiffusion フラックスは極値を強める方向に働くので, 補正フラッ
    クスは 0 になる. しかし $x$ 方向のフラックスには制限がつかないので
    補正の結果 $(i\Delta x, j\Delta y)$ での値は増加することになる.
    よって $y$ 方向の断面で見るとピークが成長していることになってしま
    う. これは FCT の基本思想と反する結果であるので, 多次元に拡張する
    際には何らかの改良を加えなければならないことになる.
    
    \begin{figure}[h]
      \begin{center}
        \Depsf[][]{ps-fig/fig6.ps}
        \caption{2次元の場合の補正フラックス
        (Zalesak (1979),図2を基に作成)}
        \Dfiglab{fig: 6}
      \end{center}
    \end{figure}