%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
%表題   3 次元積雲モデル -- 数値計算法
%
%履歴   2003-07-25 高橋 こう子   新規作成
%       2003-09-03 高橋 こう子   骨子完成
%       2003-09-05 高橋 こう子   間違い訂正
%       2003-09-06 高橋 こう子   コード解説 (1 ステップ目の計算追加)
%       2003-09-08 高橋 こう子   間違い訂正
%       2003-09-18 高橋 こう子   湿潤飽和調節法追加
%       2003-09-26 高橋 こう子   Asselin の時間フィルター追加
%       2003-11-17 高橋 こう子   圧力方程式の離散化修正
%       2003-11-22 高橋 こう子   運動エネルギーの式修正
%       2003-11-22 高橋 こう子   鉛直運動方程式修正
%       2004-01-05 高橋 こう子   鉛直運動方程式修正
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\documentclass[a4j,12pt]{jarticle}

% 使用パッケージ
\usepackage{Dennou6}
\usepackage{ascmac}

% 字下げなし
\Dnoparindent

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% 分数の定義
% (分数の中に分数がある場合, 文字の大きさを揃えるため)
\def\dfrac#1#2{{\displaystyle\frac{#1}{#2}}}
\def\cfrac#1#2{\dfrac{\mathstrut #1}{#2}}
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% 節, 式, 図, 表の表示設定
%\renewcommand{\thesection}{\arabic{section}}
\renewcommand{\thesubsection}{\thesection.\arabic{subsection}}
\renewcommand{\thesubsubsection}{\thesubsection.\arabic{subsubsection}}
%\renewcommand{\theequation}{\thesection.\arabic{equation}}
%\renewcommand{\thefigure}{\thesection.\arabic{figure}}
%\renewcommand{\thetable}{\thesection.\arabic{table}}
%
% 節, 式, 図, 表, ページの番号設定
\setcounter{secnumdepth}{3} % 番号表示する節の深さ
\setcounter{tocdepth}{3}    % 目次で表示する節の深さ

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\Dtitle{3 次元積雲モデル}

\Dauthor[高橋, 杉山, 小高, 中島, 林]
        {高橋 こう子, 杉山 耕一朗, 小高 正嗣, 中島 健介, 林 祥介}

\Ddate{2004年01月05日}

%\Dpath{/GFD\_Dennou\_Club/dc-arch/deepconv/}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{document}

\maketitle

\tableofcontents

\begin{abstract}
3 次元, 圧縮性の地形を考慮した積雲モデルの数値計算法.
\end{abstract}

%--------------------------------------------------------------------
\Dparskip

%======================================================================
\newpage

\section{数値計算法}

%--------------------------------------------------------------------
%\pagebreak

\subsection{数値計算の概要}

 本雲モデルで用いる数値解法の概要を以下にまとめる.
 数値解法は CReSS に従う.
 数理モデル\footnote{この節は後で書き直し.今はよくわらんないので.}は

 \begin{itemize}
  \item 従属変数: 速度 $u, v, w,$ 温位偏差 $\theta^{\prime},$
	気圧偏差 $p^{\prime},$ 水蒸気混合比$q_{v},$ 雲水混合
	比 $q_{c},$ 雨水混合比 $q_{r}$
  \item 基本方程式: 圧縮系
	\begin{itemize}
	 \item 運動方程式の式
	 \item 圧力方程式
	 \item 熱力学方程式
	 \item 水蒸気混合比の式
	 \item 雲水混合比の式
	 \item	雨水混合比の式
	\end{itemize}
	
  \item 地形に沿った曲線座標
  \item スタッガード格子
	\begin{itemize}
	 \item 水平方向: Arakawa C グリッド
	 \item 鉛直方向: Lorenz グリッド
	\end{itemize}
  \item 運動方程式と圧力方程式の計算:

	音波に関する項は短い時間間隔 $\Delta \tau$ で解き,
	その他の項は長い時間間隔 $\Delta t$ で解く
	\begin{itemize}
	 \item 短い時間間隔で解く時間積分
	       \begin{itemize}
		\item 水平方向の運動方程式: 陽解法, 前進差分
		      \footnote{移流項が不安定を起こさないか.
		      これでいいのか.}, 1 次精度
		\item 鉛直方向の運動方程式: 陰解法,
		      Crank-Nicolson スキーム (圧力項のみ), 
		      前進差分 (圧力項以外), 1 次精度
		\item 圧力方程式: 陰解法, Crank-Nicolson
		      スキーム (圧力項のみ), 前進差分 (圧力項
		      以外), 1 次精度
		\item 長い時間間隔で解く時間積分は一定
	       \end{itemize}
	       
	 \item 長い時間間隔で解く時間積分
	       \begin{itemize}
		\item leap frog スキーム
		\item Asselin の時間フィルター
	       \end{itemize}
	       
	 \item 空間微分
	       \begin{itemize}
		\item 中心差分, 2 次精度
	       \end{itemize}
	 \item 2 次の中心差分で表される移流項の計算
	       不安定を解消するために, 移流項に人工
	       的な粘性項を追加
	\end{itemize}
  \item 熱力学方程式, 水蒸気と雲水と雨水の混合比の式の計算:

	長い時間間隔 $\Delta t$ で解く
	\begin{itemize}
	 \item 長い時間間隔で解く時間積分: 
	       \begin{itemize}
		\item leap-frog スキーム%, 2 次精度
		\item Asselin の時間フィルター
	       \end{itemize}
	 \item 空間微分
	       \begin{itemize}
		\item 中心差分, 2 次精度
	       \end{itemize}
	 \item 2 次の中心差分で表される移流項の計算
	       不安定を解消するために, 移流項に人工
	       的な粘性項を追加
	\end{itemize}
 \end{itemize}


%--------------------------------------------------------------------
\pagebreak

\subsection{モデル格子と変数の配置}

 本雲モデルでは空間の表現を格子点法を用いているので,
 空間微分は差分法を用いて表現される.
 ここでは本モデルで用いるモデル格子の構造と変数の配置
 についてまとめる.

 \Dfigref{変数の配置}に示したように, 水平・鉛直方向の
 変数をずらして互い違いに配置する.
 配置方法は, 水平方向が Arakawa C グリッド, 鉛直方向には
 Lorenz グリッドを用いている.
 \begin{figure}[htbp]
  \begin{center}
   \includegraphics{./images/grid.eps}
   \caption{ある格子における変数の配置図.}
   \Dfiglab{変数の配置}
  \end{center}
 \end{figure}

 すべてのスカラー量 ($\phi$ : 気圧偏差, 温位偏差, 水蒸気
 混合比, 雲水混合比, 雨水混合比) を格子ボックスの中心に
 配置する.
 ベクトル量 (速度 $u, v, w$) を中心から次のように半格子
 ずらして配置する.
 水平成分 $u, v$ は鉛直方向をスカラー量と同じレベルに
 して, $u$ を $x$ 方向に半格子ずらし, $v$ を $y$ 方向に
 半格子ずらす.
 鉛直成分 $w$ は, 水平方向をスカラー量と同じレベルにして,
 $z$方向に半格子ずらす.

 さらに, 本モデルでは地形に沿う座標系を用いているので,
 座標変換に伴うヤコビアン $G^{1/2}$ と座標変換の計量
 テンソルの逆テンソルの成分のうち 0 または 1 でないもの,
 つまり, $G^{13}, G^{23}$ についても考えなければならない.
 $G^{1/2}$ は格子ボックスの中心に配置する.
 $G^{13}$ は中心から $x$ 方向に半格子ずれた $u$ の定義点から
 さらに $z$ 方向に半格子ずらす.
 $G^{23}$ は中心から $y$ 方向に半格子ずれた $v$ の定義点から
 さらに $z$ 方向に半格子ずらす.

 座標 $x, y, z^{*}$ についても半格子ずらして配置し, ベクトル量
 の位置に対応している (\Dfigref{変数の配置}).

 物理領域(方程式を計算する領域)は外縁の 1 格子分を除いた領域である.
 外縁の 1 格子は境界条件のための領域とする.
 \begin{figure}[htbp]
  \begin{center}
   \includegraphics[width=14cm]{./images/allgrid.eps}
   \caption{全計算領域における変数の配置図.
            黒丸は半格子ずらしたベクトル量, ×印はスカラー量の
            位置を表す.}
   \Dfiglab{全領域の変数の配置}
  \end{center}
 \end{figure}

%--------------------------------------------------------------------
\pagebreak

\subsection{時間積分}

 音波に関係する項 $(u, v, w, \pi)$ を短い時間間隔
 $\Delta \tau$ で解き, その他の項を長い時間間隔
 $\Delta t$ で解く.
 これをモード別時間分解法またはモード分離時間法という
 (\Dfigref{モード別時間分解法}).
 ただし $\Delta t$ を $n$ 分割したものを $\Delta \tau$ とする:
 \begin{eqnarray}
  \Delta t = n \Delta \tau
  .
 \end{eqnarray}
 \begin{figure}[htbp]
  \begin{center}
   \includegraphics[width=14cm]{./images/timestep.eps}
   \caption{モード別時間分解法の概念図.}
   \Dfiglab{モード別時間分解法}
  \end{center}
 \end{figure}


%--------------------------------------------------------------------
\pagebreak

\subsection{平均化・微分オペレータの定義}
\label{オペレータの定義}

 基礎方程式系を, スタッガード格子系について有限差分法を用いて
 離散化し, 数値積分する.
 スタッガード格子を用いているのである格子点の値を求めるために
 平均化の操作が必要である.
 この平均化のためのオペレータを以下のように定義する.

 任意の物理量 $\phi$ を空間の $x$ 方向 (または $y, z^{*}$方向)
 に平均する.
 2 次精度の中心差分に現れる 2 点の平均は
 \begin{eqnarray}
  \overline{\phi}_{i, j, k}^{x}
     \equiv \frac{1}{2}
            \left(\phi_{i - \frac{1}{2}, j, k}
	          + \phi_{i + \frac{1}{2}, j, k}
	    \right)
 \end{eqnarray}
 とする.
 ただし下付添字 $i, j, k$ は側面からそれぞれ $x, y, z$ 方向に
 $i, j, k$ 番目の格子点であることを表す.

 平均化の添字が 2 つある場合, 例えば $x, z^{*}$ 方向の平均化は
 \begin{eqnarray}
  \overline{\phi}_{i, j, k}^{xz^{*}}
     \equiv \frac{1}{4}
            \left(\phi_{i - \frac{1}{2}, j, k - \frac{1}{2}}
	          + \phi_{i + \frac{1}{2}, j, k - \frac{1}{2}}
	          + \phi_{i - \frac{1}{2}, j, k + \frac{1}{2}}
	          + \phi_{i + \frac{1}{2}, j, k + \frac{1}{2}}
	    \right)
  .
 \end{eqnarray}

 微分オペレータを定義する.
 ある方向 $x$の格子間隔を $\Delta x$ とする.
 \begin{eqnarray}
  \left(\partial_{x} \phi\right)_{i, j, k}
     = \frac{1}{\Delta x}
       \left(\phi_{i + \frac{\Delta x}{2}, j, k}
             - \phi_{i - \frac{\Delta x}{2}, j, k}\right)
 \end{eqnarray}
 平均化オペレータと同様に, 下付添字 $i$ は側面から $i$ 番目
 の格子点であることを表す.

 2 次の微分オペレータは以下のように定義する:
 \begin{eqnarray}
  \left(\partial_{x}^{2} \phi\right)_{i, j, k}
      &=& \frac{1}{\Delta x}
          \left(\frac{\phi_{i + \Delta x, j, k}
	              - \phi_{i, j, k}}
		     {\Delta x}
		- \frac{\phi_{i, j, k}
	              - \phi_{i - \Delta x, j, k}}
		     {\Delta x}
	 \right)
  \nonumber \\
      &=& \frac{1}{\Delta x^{2}}
          \left(\phi_{i + \Delta x, j, k}
	        - 2 \phi_{i, j, k}
	        + \phi_{i - \Delta x, j, k}
          \right)
  .
 \end{eqnarray}

 密度 $\bar{\rho}$ と座標変換のヤコビアン $G^{1/2}$ は同じ
 スカラー点で定義されている.
 直交直線座標での鉛直速度 $w$ と地形に沿う座標系における
 鉛直速度 $\omega$ を同じ点で定義する:
 \begin{eqnarray}
  \omega = \overline{\overline{u}^{z^{*}} G^{13}}^{x}
           + \overline{\overline{v}^{z^{*}} G^{23}}^{y}
	   + \frac{w}{\overline{G^{1/2}}^{z^{*}}}
  .
 \end{eqnarray}

%--------------------------------------------------------------------
\pagebreak

\subsection{離散化前の基礎方程式}

 数理モデルで導出した基礎方程式をそのまま離散化して数値計算を
 行うと, 計算不安定などを起こす.
 そこで人工的な項を追加することによってこうした不具合を修正する.

 タイム・スプリット法, HE-VI 法を用いると, 不安定な音波
 モードが非線形計算によって励起され成長することがある
 (Skamarock and Klemp, 1992).
 そのため, 人工的な音波減衰項 $\alpha Div$を追加して発散を抑える:
 \begin{eqnarray}
  Div = \DP{u}{x} + \DP{v}{y}
        + G^{13} \DP{u}{z^{*}}
        + G^{23} \DP{v}{z^{*}}
	+ \frac{1}{G^{1/2}}\DP{w}{z^{*}}
 \end{eqnarray}

 また, 移流項を中心差分で離散化すると計算不安定を起こすため,
 移流項に人工的な粘性項を付加する.
 2 次の数値粘性は
 \begin{eqnarray}
  {\rm Diff}.\phi &=& \nu_{2h}
                    \left[\DP[2]{\phi}{x}
		          + \DP[2]{\phi}{y}
		    \right]
		    + \nu_{2v}
		      \DP[2]{\phi}{z^{*}}
  , \\
  \nu_{2h} &=& \frac{\alpha_{2h} \Delta_{h}^{2}}{\Delta t}
  , \\
  \nu_{2v} &=& \frac{\alpha_{2v} \Delta_{v}^{2}}{\Delta t}
  , \\
  \Delta_{h} &=& \left(\Delta x \Delta y\right)^{1/2}
  , \\
  \Delta_{v} &=& \Delta z^{*}
 \end{eqnarray}
 とする.
 $\phi$ は任意の物理量であり,
 $\phi$ は $u, v, w, q_{v}, q_{c}, q_{r}$ については
 初期値からの偏差, $\theta$ は基準状態からの偏差を表す.
 $\alpha_{2 h}, \alpha_{2v}$ は無次元量で安定計算のためには
 \begin{eqnarray}
  \alpha_{2h} &\le& \frac{1}{8}
  , \\
  \alpha_{2v} &\le& \frac{1}{8}
 \end{eqnarray}
 でなければならない
 \footnote{具体的な値は不明.}.

 これらの人工的に追加した項については後で詳述する.
 以下にこれら人工的な項を追加した基礎方程式を記述する.

 {\bf 運動方程式}
  \begin{eqnarray}
  \DP{u}{t} &=& - c_{p} \bar{\theta}_{v}
                \left\{\DP{\left(\pi - \alpha Div\right)}{x}
		       + G^{13}
		        \DP{\left(\pi - \alpha Div\right)}{z^{*}}
              \right\}
  \nonumber \\
            & & - {\rm Adv}.u + {\rm Cori}.u + D_{u}
  , \\
  \DP{v}{t} &=& - c_{p} \bar{\theta}_{v}
                \left\{\DP{\left(\pi - \alpha Div\right)}{y}
		       + G^{23}
		        \DP{\left(\pi - \alpha Div\right)}{z^{*}}
              \right\}
  \nonumber \\
            & & - {\rm Adv}.v - {\rm Cori}.v + D_{v}
  , \\
  \DP{w}{t} &=& - \frac{c_{p} \bar{\theta}_{v}}{G^{1/2}}
		  \DP{\left(\pi - \alpha Div\right)}{z^{*}}
                - {\rm Adv}.w
	        + {\rm Buoy}.w
                + D_{w}
  .
 \end{eqnarray}

 {\bf 圧力方程式} :
 \begin{eqnarray}
  \DP{\pi}{t} + \frac{\bar{c}^{2}}
                     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	        \frac{1}{G^{1/2}}
	        \DP{}{z^{*}}
		\left(\bar{\rho} \bar{\theta}_{v} w\right)
     = - \frac{\bar{c}^{2}}{c_{p} \bar{\theta}_{v}}
	 \Biggl(\DP{u}{x}
	   &+& \DP{v}{y} + G^{13} \DP{u}{z^{*}}
	    + G^{23} \DP{v}{z^{*}}
	 \Biggr)
  .
 \end{eqnarray}

 {\bf 熱力学方程式} :
 \begin{eqnarray}
  \DP{\theta}{t} &=& - {\rm Adv}.\theta
                     + M_{\theta} + D_{\theta}
  .
 \end{eqnarray}

 {\bf 水蒸気の混合比の方程式} :
 \begin{eqnarray}
  \DP{q_{v}}{t} &=& - {\rm Adv}.q_{v}
                    + M_{q_{v}} + D_{q_{v}}
  .
 \end{eqnarray}

 {\bf 雲水の混合比の方程式} :
 \begin{eqnarray}
  \DP{q_{c}}{t} &=& - {\rm Adv}.q_{v}
                    + M_{q_{c}} + D_{q_{c}}
  .
 \end{eqnarray}

 {\bf 雨水の混合比の方程式} :
 \begin{eqnarray}
  \DP{q_{r}}{t} &=& - {\rm Adv}.q_{v}
                    + M_{q_{r}} + D_{q_{r}}
  .
 \end{eqnarray}
 ただし
 \begin{eqnarray}
  {\rm Adv}.\phi &=& u \DP{\phi}{x} + v \DP{\phi}{y}
		     + \omega \DP{\phi}{z^{*}}
		     + {\rm Diff}.\phi
  , \\
  {\rm Cori}.u &=& f v
  , \\
  {\rm Cori}.v &=& f u
  , \\
  {\rm Buoy}.w &=& g \left[\frac{\theta}{\bar{\theta}}
	 	           - 1 + 0.61
		           \left(q_{v} - \bar{q}_{v}\right)
                           - q_{c} - q_{r}
		      \right]
 \end{eqnarray}
 である.

\pagebreak

 拡散項:
 \begin{eqnarray}
  D_{u} &=& 2 \DP{}{x}
            \left(K_{m} \DP{u}{x}
	          + K_{m} G^{13} \DP{u}{z^{*}}
            \right)
	    + 2 G^{13} \DP{}{z^{*}}
            \left(K_{m} \DP{u}{x}
	          + K_{m} G^{13} \DP{u}{z^{*}}
            \right)
  \nonumber \\
	& & + \DP{}{y}
              \left(K_{m} \DP{u}{y}
		    + K_{m} \DP{v}{x}
		    + K_{m} G^{23} \DP{u}{z^{*}}
		    + K_{m} G^{13} \DP{v}{z^{*}}
		   \right)
  \nonumber \\
        & & + G^{23} \DP{}{z^{*}}
              \left(K_{m} \DP{u}{y}
		    + K_{m} \DP{v}{x}
		    + K_{m} G^{23} \DP{u}{z^{*}}
		    + K_{m} G^{13} \DP{v}{z^{*}}
		   \right)
  \nonumber \\
	& & + \frac{1}{G^{1/2}}\DP{}{z^{*}}
	      \left(K_{m} \DP{w}{x}
		    + \frac{K_{m}}{G^{1/2}} \DP{u}{z^{*}}
		    + K_{m} G^{13} \DP{w}{z^{*}}
		   \right)
  \nonumber \\
        & & - \frac{2}{3 C_{m}^{2} l^{2}}
	      \left(\DP{K_{m}^{2}}{x}
	            + G^{13} \DP{K_{m}^{2}}{z^{*}}
		    + \DP{K_{m}^{2}}{y}
	            + G^{23} \DP{K_{m}^{2}}{z^{*}}
		    + \frac{1}{G^{1/2}} \DP{K_{m}^{2}}{z^{*}}
	      \right)
  , \\
  D_{v} &=& \DP{}{x}
	    \left(K_{m} \DP{v}{x}
		  + K_{m} \DP{u}{y}
		  + K_{m} G^{23} \DP{u}{z^{*}}
		  + K_{m} G^{13} \DP{v}{z^{*}}
	    \right)
  \nonumber \\
        & & + G^{13} \DP{}{z^{*}}
	      \left(K_{m} \DP{v}{x}
		    + K_{m} \DP{u}{y}
		    + K_{m} G^{23} \DP{u}{z^{*}}
		    + K_{m} G^{13} \DP{v}{z^{*}}
	      \right)
  \nonumber \\
	& & + 2 \DP{}{y}
	      \left(K_{m} \DP{v}{y}
	            + K_{m} G^{23} \DP{v}{z^{*}}
	      \right)
  \nonumber \\
	& & + 2 G^{23} \DP{}{z^{*}}
	      \left(K_{m} \DP{v}{y}
	            + K_{m} G^{23} \DP{v}{z^{*}}
	      \right)
  \nonumber \\
	& & + \frac{1}{G^{1/2}}\DP{}{z^{*}}
	     \left(K_{m} \DP{w}{y}
	           + \frac{K_{m}}{G^{1/2}} \DP{v}{z^{*}}
		   + K_{m} G^{23} \DP{w}{z^{*}}
	     \right)
  \nonumber \\
        & & - \frac{2}{3 C_{m}^{2} l^{2}}
	      \left(\DP{K_{m}^{2}}{x}
	            + G^{13} \DP{K_{m}^{2}}{z^{*}}
		    + \DP{K_{m}^{2}}{y}
	            + G^{23} \DP{K_{m}^{2}}{z^{*}}
		    + \frac{1}{G^{1/2}} \DP{K_{m}^{2}}{z^{*}}
	      \right)
  , \\
  D_{w} &=& \DP{}{x}
            \left(K_{m} \DP{w}{x}
	          + \frac{K_{m}}{G^{1/2}} \DP{u}{z^{*}}
		  + K_{m} G^{13} \DP{w}{z^{*}}
	    \right)
  \nonumber \\
        &+& G^{13} \DP{}{z^{*}}
            \left(K_{m} \DP{w}{x}
	          + \frac{K_{m}}{G^{1/2}} \DP{u}{z^{*}}
		  + K_{m} G^{13} \DP{w}{z^{*}}
	    \right)
  \nonumber \\
	& & + \DP{}{y}
 	      \left(K_{m} \DP{w}{y}
		    + \frac{K_{m}}{G^{1/2}} \DP{v}{z^{*}}
		    + K_{m} G^{23} \DP{w}{z^{*}}
	      \right)
  \nonumber \\
	& & + G^{23} \DP{}{z^{*}}
 	      \left(K_{m} \DP{w}{y}
		    + \frac{K_{m}}{G^{1/2}} \DP{v}{z^{*}}
		    + K_{m} G^{23} \DP{w}{z^{*}}
	      \right)
  \nonumber \\
	& & + \frac{2}{G^{1/2}}\DP{}{z^{*}}
              \left(\frac{K_{m}}{G^{1/2}}
	            \DP{w}{z^{*}}
              \right)
  \nonumber \\
        & & - \frac{2}{3 C_{m}^{2} l^{2}}
	      \left(\DP{K_{m}^{2}}{x}
	            + G^{13} \DP{K_{m}^{2}}{z^{*}}
		    + \DP{K_{m}^{2}}{y}
	            + G^{23} \DP{K_{m}^{2}}{z^{*}}
		    + \frac{1}{G^{1/2}} \DP{K_{m}^{2}}{z^{*}}
	      \right)
 \end{eqnarray}
\pagebreak
 \begin{eqnarray}
  D_{\theta} &=& 3 \DP{}{x}
                 \left(K_{m}\DP{\theta}{x}
		       + K_{m} G^{13}
		         \DP{\theta}{z^{*}}\right)
                 + 3 \DP{}{y}
                     \left(K_{m}\DP{\theta}{y}
		           + K_{m} G^{23}
			   \DP{\theta}{z^{*}}
		     \right)
  \nonumber \\
             & & + 3 G^{13} \DP{}{z^{*}}
                 \left(K_{m}\DP{\theta}{x}
		       + K_{m} G^{13}
		         \DP{\theta}{z^{*}}\right)
  \nonumber \\
             & & + 3 G^{23} \DP{}{z^{*}}
                     \left(K_{m}\DP{\theta}{y}
		           + K_{m} G^{23}
			   \DP{\theta}{z^{*}}
		     \right)
  \nonumber \\
             & & + \frac{3}{G^{1/2}} \DP{}{z^{*}}
	           \left(\frac{K_{m}}{G^{1/2}}
		         \DP{\theta}{z^{*}}\right)
  , \\
  D_{q_{v}} &=& 3 \DP{}{x}
                 \left(K_{m}\DP{q_{v}}{x}
		       + K_{m} G^{13}
		       \DP{q_{v}}{z^{*}}
		 \right)
                 + 3 \DP{}{y}
	             \left(K_{m}\DP{q_{v}}{y}
		           + K_{m} G^{23}
			   \DP{q_{v}}{z^{*}}
		     \right)
  \nonumber \\
             & & + 3 G^{13} \DP{}{z^{*}}
                 \left(K_{m}\DP{q_{v}}{x}
		       + K_{m} G^{13}
		       \DP{q_{v}}{z^{*}}
		 \right)
  \nonumber \\
             & & + 3  G^{23} \DP{}{z^{*}}
	             \left(K_{m}\DP{q_{v}}{y}
		           + K_{m} G^{23}
			   \DP{q_{v}}{z^{*}}
		     \right)
  \nonumber \\
             & & + \frac{3}{G^{1/2}} \DP{}{z^{*}}
	           \left(\frac{K_{m}}{G^{1/2}}
		         \DP{q_{v}}{z^{*}}\right)
  , \\
  D_{q_{c}} &=& 3 \DP{}{x}
                  \left(K_{m}\DP{q_{c}}{x}
		        + K_{m} G^{13}
		        \DP{q_{c}}{z^{*}}
		  \right)
                + 3 \DP{}{y}
                      \left(K_{m}\DP{q_{c}}{y}
		           + K_{m} G^{23}
			   \DP{q_{c}}{z^{*}}
		      \right)
  \nonumber \\
            & & + 3 G^{13} \DP{}{z^{*}}
                  \left(K_{m}\DP{q_{c}}{x}
		        + K_{m} G^{13}
		        \DP{q_{c}}{z^{*}}
		  \right)
  \nonumber \\
             & &  + 3 G^{23} \DP{}{z^{*}}
                      \left(K_{m}\DP{q_{c}}{y}
		           + K_{m} G^{23}
			   \DP{q_{c}}{z^{*}}
		      \right)
  \nonumber \\
             & & + \frac{3}{G^{1/2}} \DP{}{z^{*}}
	           \left(\frac{K_{m}}{G^{1/2}}
		         \DP{q_{c}}{z^{*}}\right)
  , \\
  D_{q_{r}} &=& 3 \DP{}{x}
                  \left(K_{m}\DP{q_{r}}{x}
		        + K_{m} G^{13}
		        \DP{q_{r}}{z^{*}}
		  \right)
                + 3 \DP{}{y}
                    \left(K_{m}\DP{q_{r}}{y}
		         + K_{m} G^{23}
			   \DP{q_{r}}{z^{*}}\right)
  \nonumber \\
             & & + 3 G^{13} \DP{}{z^{*}}
                  \left(K_{m}\DP{q_{r}}{x}
		        + K_{m} G^{13}
		        \DP{q_{r}}{z^{*}}
		  \right)
  \nonumber \\
             & & + 3 G^{23} \DP{}{z^{*}}
                    \left(K_{m}\DP{q_{r}}{y}
		         + K_{m} G^{23}
			   \DP{q_{r}}{z^{*}}\right)
  \nonumber \\
             & & + \frac{3}{G^{1/2}} \DP{}{z^{*}}
	           \left(\frac{K_{m}}{G^{1/2}}
		         \DP{q_{r}}{z^{*}}\right)
  .
 \end{eqnarray}
 ただし
 \begin{eqnarray}
  K_{h} &=& 3 K_{m}
  , \\
  K_{m} &=& C_{m} E^{1/2} l
 \end{eqnarray}
 を用いて書き換えている.

 生成・消滅項:
 \begin{eqnarray}
  M_{\theta} &=& - \gamma
                 \left(\DP{q_{vs}}{t} + {\rm Adv}.q_{vs}
		       + E_{r}\right)
  , \\
  M_{q_{v}} &=& \DP{q_{vs}}{t} + {\rm Adv}.q_{vs} + E_{r}
  , \\
  M_{q_{c}} &=& - \left(\DP{q_{vs}}{t}
		         + {\rm Adv}.q_{vs}\right)
		- A_{r} - C_{r}
  , \\
  M_{q_{r}} &=& F_{r} - E_{r} + A_{r} +C_{r}
  .
 \end{eqnarray}

 ただし, 
 \begin{eqnarray}
  {\rm Adv}.q_{vs} &=& u \DP{q_{vs}}{x}
		       + v \DP{q_{vs}}{y}
		       + \omega \DP{q_{vs}}{z^{*}}
		       + {\rm Diff}.q_{vs}
 , \\
  E_{r} &=& \frac{1}{\bar{\rho}}
	  \frac{\left(1 - \dfrac{q_{v}}{q_{vs}}\right)
	        C \left(\bar{\rho} q_{r}\right)^{0.525}}
	       {5.4 \times 10^{5}
	  + \dfrac{2.55 \times 10^{6}}{\bar{p} q_{vs}}}
  , \\
  C &=& 1.6
      + 124.9 \left(\bar{\rho}q_{r}\right)^{0.2046}
  , \\
  A_{r} &=& k_{1} \left(q_{c} - a\right)
  , \\
  C_{r} &=& k_{2} q_{c} q_{r}^{0.875}
  , \\
  F_{r} &=& \frac{1}{G^{1/2} \bar{\rho}}
	      \DP{}{z^{*}}
	      \left(\bar{\rho} V q_{r}\right)
  , \\
  V &=& 36.34 \left(\bar{\rho} q_{r}\right)^{0.1346}
        \left(\frac{\bar{\rho}}{\rho_{0}}\right)^{-\frac{1}{2}}
  , \\
  k_{1} &=& 0.001 ~ {\rm s^{-1}}
  , \\
  a &=& 0.001 ~ {\rm kg \, kg^{-1}}
  , \\
  k_{2} &=& 2.2 ~ {\rm s^{-1}}
  \
 \end{eqnarray}
 
\pagebreak

 サブグリッドスケールの運動エネルギーの方程式:
 \begin{eqnarray}
  \DP{E}{t} &=& - {\rm Adv}.E
  \nonumber \\
            & & + \frac{3 g C_{m} l E^{1/2}}{G^{1/2}}
		  \left(- A \DP{\theta_{e}}{z^{*}}
		         + \DP{q_{l}}{z^{*}}
		  \right)
  \nonumber \\
            & & + 2 C_{m} l E^{1/2}
	          \left\{\left(\DP{u}{x}
			         + G^{13} \DP{u}{z^{*}}
			   \right)^{2}
			   + \left(\DP{v}{y}
			         + G^{23} \DP{v}{z^{*}}
			     \right)^{2}
			   + \left(\frac{1}{G^{1/2}}
			           \DP{w}{z}
			     \right)^{2}
		   \right\}
  \nonumber \\
	      & & + C_{m} l E^{1/2}
	            \Biggl\{\left(\DP{v}{x}
			        + \DP{u}{y}
				+ G^{13} \DP{v}{z^{*}}
				+ G^{23} \DP{u}{z^{*}}
			  \right)^{2}
  \nonumber \\
	      & & \hspace*{4em}
			+ \left(\DP{w}{x}
				+ G^{13} \DP{w}{z^{*}}
			        + \frac{1}{G^{1/2}}
				  \DP{u}{z^{*}}
			  \right)^{2}
			+ \left(\DP{w}{y}
				+ G^{23} \DP{w}{z^{*}}
			        + \frac{1}{G^{1/2}}
				  \DP{v}{z^{*}}
			  \right)^{2}
	        \Biggr\}
  \nonumber \\
              & & + C_{m} l E^{1/2}
	            \Biggl\{
		           \left(\DP{}{x}
			         + G^{13} \DP{}{z^{*}}
			   \right)
		           \left(\DP{E}{x}
			         + G^{13} \DP{E}{z^{*}}
			   \right)
  \nonumber \\
	      & & \hspace*{4em}
		           + \left(\DP{}{y}
			         + G^{23} \DP{}{z^{*}}
			   \right)
		           \left(\DP{E}{y}
			         + G^{23} \DP{E}{z^{*}}
			   \right)
			   + \left(\frac{1}{G^{1/2}}
			           \DP{}{z^{*}}
			     \right)
			     \left(\frac{1}{G^{1/2}}
			           \DP{E}{z^{*}}
			     \right)
	            \Biggr\}
  \nonumber \\
              & & + 2 C_{m} l E^{1/2}
	            \Biggl\{\left(\DP{E^{1/2}}{x}
				  + G^{13} \DP{E^{1/2}}{z^{*}}
			    \right)^{2}
  \nonumber \\
              & & \hspace*{8em}
			    + \left(\DP{E^{1/2}}{y}
				    + G^{23} \DP{E^{1/2}}{z^{*}}
			      \right)^{2}
			    + \left(\frac{1}{G^{1/2}}
				    \DP{E^{1/2}}{z^{*}}
			      \right)^{2}
			     \Biggr\}
  \nonumber \\
	     & & - \frac{C_{\varepsilon}}{l} E^{3/2}
  \nonumber \\
  {\rm Adv}.E &=& u \DP{E}{x}
                  + v \DP{E}{y} + \omega \DP{E}{z^{*}}
                  + {\rm Diff}.E
 \end{eqnarray}


%--------------------------------------------------------------------
\pagebreak

\subsection{水平方向の運動方程式の離散化}

 \ref{オペレータの定義} 節で定義した表記を用いて, $(i, j, k)$ 
 番目の格子について音波を解に含む基礎方程式を離散化する.

 まず, 短い時間間隔 $\tau$ で計算する項について離散化したものを
 書き下す.

 水平方向の運動方程式は陽解法で離散化する.

 {\bf x 方向の運動方程式}:
 \begin{eqnarray}
  \frac{u_{i + \frac{\Delta x}{2}, j, k}^{\tau + \Delta \tau}
      - u_{i + \frac{\Delta x}{2}, j, k}^{\tau}}{\Delta \tau}
      &=& - c_{p} \bar{\theta}_{v}
            \Biggl[\left\{
		         \left(\partial_{x}\pi\right)
			 _{i + \frac{\Delta x}{2}, j, k}
			 - \alpha
			 \left(\partial_{x} Div\right)
			 _{i + \frac{\Delta x}{2}, j, k}
	          \right\}^{\tau}
  \nonumber \\
      & &  \vspace*{4em}
                  + \overline{G^{13}}^{z^{*}}
		    \left\{
		           \left(\partial_{z^{*}}
			         \overline{\pi}^{xz^{*}}
			   \right)_{i + \frac{\Delta x}{2}, j, k}
			   - \alpha
			   \left(\partial_{z^{*}}
			         \overline{Div}^{xz^{*}}
			   \right)_{i + \frac{\Delta x}{2}, j, k}
		    \right\}^{\tau}
	    \Biggr]
  \nonumber \\
      & & - {\rm Adv}.u_{i + \frac{\Delta x}{2}, j, k}^{t}
          + {\rm Cori}.u_{i + \frac{\Delta x}{2}, j, k}^{t}
	  + \left(D_{u}\right)_{i + \frac{\Delta x}{2}, j, k}
	                      ^{t - \Delta t}
  \Deqlab{離散水平運動u}
 \end{eqnarray}
 {\bf y 方向の運動方程式}:
 \begin{eqnarray}
  \frac{v_{i, j + \frac{\Delta y}{2}, k}^{\tau + \Delta \tau}
      - v_{i, j + \frac{\Delta y}{2}, k}^{\tau}}{\Delta \tau}
      &=& - c_{p} \bar{\theta}_{v}
            \Biggl[\left\{
		         \left(\partial_{y}\pi\right)
			 _{i, j + \frac{\Delta y}{2}, k}
			 - \alpha
			 \left(\partial_{x} Div\right)
			 _{i, j + \frac{\Delta y}{2}, k}
	          \right\}^{\tau}
  \nonumber \\
      & &  \vspace*{4em}
                  + \overline{G^{23}}^{z^{*}}
		    \left\{
		           \left(\partial_{z^{*}}
			         \overline{\pi}^{yz^{*}}
			   \right)_{i, j + \frac{\Delta y}{2}, k}
			   - \alpha
			   \left(\partial_{z^{*}}
			         \overline{Div}^{yz^{*}}
			   \right)_{i, j + \frac{\Delta y}{2}, k}
		    \right\}^{\tau}
	    \Biggr]
  \nonumber \\
      & & - {\rm Adv}.v_{i, j + \frac{\Delta y}{2}, k}^{t}
          + {\rm Cori}.v_{i, j + \frac{\Delta y}{2}, k}^{t}
	  + \left(D_{v}\right)_{i, j + \frac{\Delta y}{2}, k}
	                      ^{t - \Delta t}
  \Deqlab{離散水平運動v}
  .
 \end{eqnarray}
 $Div$ を含む項を書き下すと以下のようになる.
 ただし上付添字 $\tau$ は省略した.
 \begin{eqnarray}
  Div_{i, j, k} &=& \left(\partial_{x} u\right)_{i, j, k}
                    + \left(\partial_{y} v\right)_{i, j, k}
          + \overline{G^{13}}_{i, j, k}^{xz^{*}}
	    \left(\partial_{z^{*}}
	          \overline{u}^{xz^{*}}
	    \right)_{i, j, k}
  \nonumber \\
      & & + \overline{G^{23}}_{i, j, k}^{yz^{*}}
	    \left(\partial_{z^{*}}
	          \overline{v}^{xz^{*}}
	    \right)_{i, j, k}
          + \frac{1}{G_{i, j, k}^{1/2}}
	    \left(\partial_{z^{*}} w\right)_{i, j, k}
  , \\
  \left(\partial_{x} Div\right)_{i + \frac{\Delta x}{2}, j, k}
      &=& \frac{Div_{i + \Delta x, j, k}- Div_{i, j, k}}
               {\Delta x}
  , \\
  \left(\partial_{z^{*}} \overline{Div}^{xz^{*}}
  \right)_{i + \frac{\Delta x}{2}, j, k}
      &=& \frac{\overline{Div}^{xz^{*}}
                _{i + \frac{\Delta x}{2}, j, 
		  k+ \frac{\Delta z^{*}}{2}}
                - \overline{Div}^{xz^{*}}
	       _{i + \frac{\Delta x}{2}, j,
		k - \frac{\Delta z^{*}}{2}}}
	       {\Delta z^{*}}
  , \\
  \left(\partial_{y} Div\right)_{i, j + \frac{\Delta y}{2}, k}
      &=& \frac{Div_{i + \Delta x, j, k}- Div_{i, j, k}}
               {\Delta y}
  , \\
  \left(\partial_{z^{*}} \overline{Div}^{yz^{*}}
  \right)_{i, j + \frac{\Delta y}{2}, k}
      &=& \frac{\overline{Div}^{yz^{*}}
                _{i, j + \frac{\Delta y}{2},
		  k+ \frac{\Delta z^{*}}{2}}
                - \overline{Div}^{yz^{*}}
	       _{i, j + \frac{\Delta y}{2},
		k - \frac{\Delta z^{*}}{2}}}
	       {\Delta z^{*}}
 \end{eqnarray}

 次に, 長い時間間隔 $t$ で計算する項 (各方程式の上付添字
 $t$ がついた項)について離散化したものを書き下す.

 移流項:
 \begin{eqnarray}
  {\rm Adv}.u^{t}_{i + \frac{\Delta x}{2}, j, k}
     &=& \left(u \partial_{x}
	       \overline{u}^{x}
	 \right)_{i + \frac{\Delta x}{2}, j, k}^{t}
         + \left(\overline{v}^{xy}
	         \partial_{y} \overline{u}^{y}
	   \right)_{i + \frac{\Delta x}{2}, j, k}^{t}
  \nonumber \\
     & & + \left(\overline{\omega}^{xz^{*}}
	         \partial_{z^{*}} \overline{u}^{z^{*}}
	   \right)_{i + \frac{\Delta x}{2}, j, k}^{t}
         + {\rm Diff}.u_{i + \frac{\Delta x}{2}, j, k}^{t}
  , \\
  {\rm Adv}.v^{t}_{i, j + \frac{\Delta y}{2}, k}
     &=& \left(\overline{u}^{xy}
	       \partial_{x} \overline{v}^{x}
	 \right)_{i, j + \frac{\Delta x}{2}, k}^{t}
         + \left(v \partial_{y}
	         \overline{v}^{y}
	   \right)_{i, j + \frac{\Delta x}{2}, k}^{t}
  \nonumber \\
     & & + \left(\overline{\omega}^{yz^{*}}
	         \partial_{z^{*}} \overline{v}^{z^{*}}
	   \right)_{i, j + \frac{\Delta y}{2}, k}^{t}
         + {\rm Diff}.v_{i, j + \frac{\Delta x}{2}, k}^{t}
  .
 \end{eqnarray}

 数値粘性項:
 \begin{eqnarray}
  {\rm Diff}.u_{i + \frac{\Delta x}{2}, j, k}^{t}
      &=& \nu_{2h}
          \left[\partial_{x}^{2} u
	        + \partial_{y}^{2} u
	  \right]_{i + \frac{\Delta x}{2}, j, k}^{t}
	  + \nu_{2v} 
	    \left[\partial_{z^{*}}^{2} u
	    \right]_{i + \frac{\Delta x}{2}, j, k}^{t}
  , \\
  {\rm Diff}.v_{i, j + \frac{\Delta y}{2}, k}^{t}
      &=& \nu_{2h}
          \left[\partial_{x}^{2} v
	       + \partial_{y}^{2} v
	  \right]_{i, j + \frac{\Delta y}{2}, k}^{t}
	  + \nu_{2v} 
	    \left[\partial_{z^{*}}^{2} v
	    \right]_{i, j + \frac{\Delta y}{2}, k}^{t}
 \end{eqnarray}

 コリオリ項:
 \begin{eqnarray}
  {\rm Cori}.u_{i + \frac{\Delta x}{2}, j, k}^{t}
      &=& \left[\overline{f \overline{v}^{y}}^{x}
	  \right]_{i + \frac{\Delta x}{2}, j, k}^{t}
  , \\
  {\rm Cori}.v_{i, j + \frac{\Delta y}{2}, k}^{t}
      &=& \left[\overline{f \overline{u}^{x}}^{y}
	  \right]_{i, j + \frac{\Delta y}{2}, k}^{t}
 \end{eqnarray}

\pagebreak
 拡散項:
 (右辺の下付添字と上付添字 $t - \Delta t$ は省略した)
 \begin{eqnarray}
  D_{u.i + \frac{1}{2}, j, k}
   &=& 2 \partial_{x}
            \left\{K_{m}
	         \left(\partial_{x} u\right)
	          + K_{m} \overline{G^{13}}^{xz^{*}}
		    \left(\partial_{z^{*}}
		    \overline{u}^{z^{*}}\right)
            \right\}
  \nonumber \\
	& & + \partial_{y}
              \Biggl\{\overline{K_{m}}^{xy}
	            \left(\partial_{y} u\right)
		    + \overline{K_{m}}^{xy}
		      \left(\partial_{x} v\right)
		    + \overline{K_{m}}^{xy}
		                \overline{G^{23}}^{xz^{*}}
		      \left(\partial_{z^{*}}
		            \overline{u}^{yz{*}}
		      \right)
  \nonumber \\
	& & \hspace*{4em}
	    + \overline{K_{m}}^{xy}
	      \overline{G^{13}}^{yz^{*}}
		      \left(\partial_{z^{*}}
		            \overline{v}^{xz^{*}}
		      \right)
	      \Biggr\}
  \nonumber \\
	& & + 2 \overline{G^{13}}^{z^{*}}
	      \partial_{z^{*}}
              \left\{\overline{K_{m}}^{xz^{*}}
	            \left(\partial_{x}
		    \overline{u}^{xz^{*}}\right)
	          + \overline{K_{m}}^{xz^{*}} G^{13}
		    \left(\partial_{z^{*}}
		          u\right)
            \right\}
  \nonumber \\
        & & + \overline{G^{23}}^{xyz^{*}}
	      \partial_{z^{*}}
              \Biggl\{\overline{K_{m}}^{xz^{*}}
	             \left(\partial_{y}
		           \overline{u}^{yz^{*}}
		     \right)
		    + \overline{K_{m}}^{xz^{*}}
		      \partial_{x} \overline{v}^{yz^{*}}
  \nonumber \\
	& & \hspace*{7em}
	 + \overline{K_{m}}^{xz^{*}}
		      \overline{G^{23}}^{xy}
		      \left(\partial_{z^{*}} u\right)
		    + \overline{K_{m}}^{xz^{*}}
		      \overline{G^{13}}^{z^{*}}
		      \left(\partial_{z^{*}}
		            \overline{v}^{xy}
		      \right)
		   \Biggr\}
  \nonumber \\
	& & + \frac{1}{\overline{G^{1/2}}^{x}}
	      \partial_{z^{*}}
	      \left\{\overline{K_{m}}^{xz^{*}}
	             \left(\partial_{x} w\right)
		    + \frac{\overline{K_{m}}^{xz^{*}}}
		           {\overline{G^{1/2}}^{xz^{*}}}
		      \left(\partial_{z^{*}} u\right)
		    + \overline{K_{m}}^{xz^{*}}
		      G^{13}
		      \partial_{z^{*}} \overline{w}^{xz^{*}}
	      \right\}
  \nonumber \\
        & & - \frac{2}{3 C_{m}^{2} l^{2}}
	      \Biggl\{\left(\partial_{x} K_{m}^{2}\right)
	             + \overline{G^{13}}^{z^{*}}
		       \left(\partial_{z^{*}}
			     \overline{K_{m}^{2}}^{xz^{*}}
		       \right)
  \nonumber \\
	& & \hspace*{6em}
	             + \left(\partial_{y}
			     \overline{K_{m}^{2}}^{xy}
		       \right)
		     + \overline{G^{23}}^{xz^{*}}
		       \left(\partial_{z^{*}}
			     \overline{K_{m}^{2}}^{yz^{*}}
		       \right)
  \nonumber \\
	& & \hspace*{6em}
		     + \frac{1}{\overline{G^{1/2}}^{x}}
		       \left(\partial_{z^{*}}
			     \overline{K_{m}^{2}}^{xz^{*}}
		       \right)
              \Biggr\}
  , \\
  D_{v.i, j + \frac{1}{2}, k}
       &=& \partial_{x}
	    \Biggl\{\overline{K_{m}}^{xy}
	          \left(\partial_{x} v\right)
		  + \overline{K_{m}}^{xy}
		    \left(\partial_{y} u\right)
  \nonumber \\
       & & \hspace*{4em}
		  + \overline{K_{m}}^{xy}
		    \overline{G^{23}}^{xz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{u}^{yz^{*}}
		    \right)
		  + \overline{K_{m}}^{xy}
		    \overline{G^{13}}^{yz^{*}}
		    \partial_{z^{*}}
		    \left(\overline{v}^{xz^{*}}\right)
	    \Biggr\}
  \nonumber \\
	& & + 2 \partial_{y}
	      \left\{K_{m} \left(\partial_{y} v\right)
	            + K_{m} \overline{G^{23}}^{yz^{*}}
		      \left(\partial_{z^{*}}
		            \overline{v}^{yz^{*}}\right)
	      \right\}
  \nonumber \\
        & & + \overline{G^{13}}^{xyz^{*}}
	      \partial_{z^{*}}
	      \Biggl\{\overline{K_{m}}^{yz^{*}}
	              \left(\overline{\partial_{x} v}^{xz^{*}}\right)
		    + \overline{K_{m}}^{yz^{*}}
		      \left(\partial_{y} \overline{u}^{xz^{*}}\right)
  \nonumber \\
        & & \hspace*{8em}
		    + \overline{K_{m}}^{yz^{*}}
		      G^{23} \left(\partial_{z^{*}}\right)
		      \overline{u}^{xy}
		    + \overline{K_{m}}^{yz^{*}}
		      \overline{G^{13}}^{xy}
		      \left(\partial_{z^{*}} v\right)
	      \Biggr\}
  \nonumber \\
	& & + 2 \overline{G^{23}}^{z^{*}}
	      \partial_{z^{*}}
	      \left\{\overline{K_{m}}^{yz^{*}}
	              \left(\partial_{y} \overline{v}^{yz^{*}}\right)
	              + \overline{K_{m}}^{yz^{*}}
		        G^{23} \left(\partial_{z^{*}} v\right)
	      \right\}
  \nonumber \\
        & & - \frac{2}{3 C_{m}^{2} l^{2}}
	      \Biggl\{\left(\partial_{x}
			     \overline{K_{m}^{2}}^{xy}
		       \right)
	             + \overline{G^{13}}^{yz^{*}}
		       \left(\partial_{z^{*}}
			     \overline{K_{m}^{2}}^{yz^{*}}
		       \right)
  \nonumber \\
	& & \hspace*{6em}
	             + \left(\partial_{y} K_{m}^{2}\right)
		     + \overline{G^{23}}^{z^{*}}
		       \left(\partial_{z^{*}}
			     \overline{K_{m}^{2}}^{yz^{*}}
		       \right)
  \nonumber \\
	& & \hspace*{6em}
		     + \frac{1}{\overline{G^{1/2}}^{y}}
		       \left(\partial_{z^{*}}
			     \overline{K_{m}^{2}}^{yz^{*}}
		       \right)
              \Biggr\}
 \end{eqnarray}

%--------------------------------------------------------------------
\pagebreak

\subsection{鉛直方向の運動方程式と圧力方程式の離散化}

 鉛直方向の運動方程式と圧力方程式は陰解法で離散化する.
 $w, \pi$ は重み $\beta$ をつけて $\tau$ と
 $\tau + \Delta \tau$ の平均をしたものを用いる. 

 {\bf 鉛直方向の運動方程式}:
 \begin{eqnarray}
  \frac{w_{i, j, k + \frac{1}{2}}^{\tau + \Delta \tau}
     - w_{i, j, k + \frac{1}{2}}^{\tau}}{\Delta \tau}
     &=& \alpha
           \left(\frac{c_{p} \bar{\theta}_{v}}
	              {\overline{G^{1/2}}^{z^{*}}}
	   \right)_{i, j, k + \frac{1}{2}}
	   \left(\partial_{z^{*}} Div \right)
	   _{i, j, k + \frac{1}{2}}^{\tau}
  \nonumber \\
       & &  - \left(\frac{c_{p} \bar{\theta}_{v}}
	                 {\overline{G^{1/2}}^{z^{*}}}
	      \right)_{i, j,k + \frac{1}{2}}
	    \Biggl[\beta
	                  \left(\partial_{z^{*}} \pi \right)
			  _{i, j, k + \frac{1}{2}}
			  ^{\tau + \Delta \tau}
		 + \left(1 - \beta\right)
	                  \left(\partial_{z^{*}} \pi \right)
			  _{i, j, k + \frac{1}{2}}^{\tau}
	    \Biggr]
  \nonumber \\
     & & - {\rm Adv}.w_{i, j, k + \frac{1}{2}}^{t}
	 + {\rm Buoy}.w_{i, j, k + \frac{1}{2}}^{t}
	  + \left(D_{w}\right)_{i, j, k + \frac{1}{2}}
	                      ^{t - \Delta t}
  ,
  \Deqlab{離散鉛直運動}
 \end{eqnarray}
 {\bf 圧力方程式}:
 \begin{eqnarray}
  \frac{\pi_{i, j, k}^{\tau + \Delta \tau}
      - \pi_{i, j, k}^{\tau}}{\Delta \tau}
      &+& \frac{1}{G^{1/2}}
          \frac{\bar{c}^{2}}
               {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \Biggl\{\beta
	          \left(\partial_{z^{*}}
		       \bar{\rho} \bar{\theta}_{v} w
		  \right)_{i, j, k}^{\tau + \Delta \tau}
		  + \left(1 - \beta\right)
		    \left(\partial_{z^{*}}
		          \bar{\rho} \bar{\theta}_{v} w
		    \right)_{i, j, k}^{\tau}
	  \Biggr\}
  \nonumber \\
     &=& - \frac{\bar{c}^{2}}
	        {c_{p} \bar{\theta}_{v}}
	   \Biggl[\left(\partial_{x} u\right)_{i, j, k}
		  + \left(\partial_{y}v\right)_{i, j, k}
  \nonumber \\
     & & 
     \hspace*{4em}
	          + \overline{G^{13}}_{i, j, k}^{xz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{u}^{x z^{*}}
		    \right)_{i, j, k}
		  + \overline{G^{23}}_{i, j, k}^{yz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{v}^{y z^{*}}
		    \right)_{i, j, k}
	   \Biggr]^{\tau + \Delta \tau}
  \Deqlab{離散圧力}
  .
 \end{eqnarray}
 ここで $\beta \ne 0$ である.
 \Deqref{離散鉛直運動}と\Deqref{離散圧力}は
 $w^{\tau + \Delta \tau}$ と $\pi^{\tau + \Delta \tau}$ の
 連立方程式となり, Crank-Nicolson スキームにより鉛直について
 インプリシットに解くことができる.
 この連立方程式の解き方は後で述べる.

 $Div$ を含む項を書き下すと以下のようになる.
 ただし上付添字 $\tau$ は省略した.
 \begin{eqnarray}
  \left(\partial_{z^{*}} Div\right)_{i, j, k + \frac{\Delta z^{*}}{2}}
      &=& \frac{Div_{i, j, k + \Delta z^{*}}- Div_{i, j, k}}
               {\Delta z^{*}}
 \end{eqnarray}

 長い時間間隔 $t$ で計算する項 (各方程式の上付添字
 $t$ がついた項)について離散化したものを書き下す.

 移流項:
 \begin{eqnarray}
  {\rm Adv}.w^{t}_{i, j, k + \frac{\Delta z^{*}}{2}}
     &=& \left(\overline{u}^{xz^{*}}
               \partial_{x} \overline{w}^{x}
	 \right)_{i, j, k + \frac{\Delta z^{*}}{2}}^{t}
         + \left(\overline{v}^{yz^{*}}
	         \partial_{y} \overline{w}^{y}
	   \right)_{i, j, k + \frac{\Delta z^{*}}{2}}^{t}
  \nonumber \\
     & & + \left(\omega \partial_{z^{*}}
	        \overline{w}^{z^{*}}
	   \right)_{i, j, k + \frac{\Delta z^{*}}{2}}^{t}
         + {\rm Diff}.w_{i, j, k + \frac{\Delta z^{*}}{2}}^{t}
  .
 \end{eqnarray}

 数値粘性項:
 \begin{eqnarray}
  {\rm Diff}.w_{i, j, k + \frac{\Delta z^{*}}{2}}^{t}
      &=& \nu_{2h}
          \left[\partial_{x}^{2} w
	        + \partial_{y}^{2} w
	  \right]_{i, j, k + \frac{\Delta z^{*}}{2}}^{t}
	  + \nu_{2v}
	    \left[\partial_{z^{*}}^{2} w
	    \right]_{i, j, k + \frac{\Delta z^{*}}{2}}^{t}
  , \\
 \end{eqnarray}

 浮力項:
 \begin{eqnarray}
  {\rm Buoy}.w_{i, j, k + \frac{\Delta z^{*}}{2}}^{t}
     = g \left[\overline{\frac{\theta}
	                {\bar{\theta}}}^{z^{*}}
	       - 1 + 0.61
	       \left(\overline{q_{v}}^{z^{*}}
		- \overline{\bar{q}_{v}}^{z^{*}}\right)
	       - \overline{q_{c}}^{z^{*}}
	       - \overline{q_{r}}^{z^{*}}
	  \right]_{i, j, k + \frac{\Delta z^{*}}{2}}^{t}
 \end{eqnarray}

 拡散項:
 (右辺の下付添字と上付添字 $t - \Delta t$ は省略した)
 \begin{eqnarray}
  D_{w. i, j, k + \frac{1}{2}}
        &=& \partial_{x}
            \left\{\overline{K_{m}}^{xz^{*}}
	           \left(\partial_{x} w\right)
	           + \frac{\overline{K_{m}}^{xz^{*}}}
		          {\overline{G^{1/2}}^{xz^{*}}}
		     \left(\partial_{z^{*}}
		           \overline{u}^{xz^{*}}
		     \right)
		   + \overline{K_{m}}^{xz^{*}}
		     G^{13}
		     \left(\partial_{z^{*}}
		           \overline{w}^{xz^{*}}
		     \right)
	    \right\}
  \nonumber \\
	& & + \partial_{y}
 	      \left\{\overline{K_{m}}^{yz^{*}}
	             \left(\partial_{y} w\right)
		     + \frac{\overline{K_{m}}^{yz^{*}}}
		            {\overline{G^{1/2}}^{yz^{*}}}
		       \left(\partial_{z^{*}} v\right)
		     + \overline{K_{m}}^{yz^{*}}
		       G^{23}
		       \left(\partial_{z^{*}}
			     \overline{w}^{yz^{*}}
		       \right)
	      \right\}
  \nonumber \\
        & & + \overline{G^{13}}^{x}
	      \partial_{z^{*}}
            \left\{K_{m}
	           \left(\partial_{x}
		         \overline{w}^{xz^{*}}
		   \right)
	           + \frac{K_{m}}{G^{1/2}}
		     \left(\partial_{z^{*}}
		           \overline{u}^{xz^{*}}
		     \right)
		   + K_{m} \overline{G^{13}}^{xz^{*}}
		     \left(\partial_{z^{*}} w\right)
	    \right\}
  \nonumber \\
	& & + \overline{G^{23}}^{y}
	      \partial_{z^{*}}
 	      \left\{K_{m}
	             \left(\partial_{y}
		           \overline{w}^{yz^{*}}
		     \right)
		     + \frac{K_{m}}{G^{1/2}}
		       \left(\partial_{z^{*}}
			     \overline{v}^{yz^{*}}
		       \right)
		     + K_{m} \overline{G^{23}}^{yz^{*}}
		       \left(\partial_{z^{*}} w\right)
	      \right\}
  \nonumber \\
	& & + \frac{2}{\overline{G^{1/2}}^{z^{*}}}
	      \partial_{z^{*}}
              \left\{\frac{K_{m}}{G^{1/2}}
	            \left(\partial_{z^{*}} w\right)
              \right\}
  \nonumber \\
        & & - \frac{2}{3 C_{m}^{2} l^{2}}
	      \Biggl\{\left(\partial_{x}
			     \overline{K_{m}^{2}}^{xz^{*}}
		       \right)
	             + \overline{G^{13}}^{x}
		       \left(\partial_{z^{*}} K_{m}^{2}\right)
  \nonumber \\
	& & \hspace*{4em}
	             + \left(\partial_{y}
			     \overline{K_{m}^{2}}^{yz^{*}}
		       \right)
		     + \overline{G^{23}}^{y}
		       \left(\partial_{z^{*}} K_{m}^{2}\right)
		     + \frac{1}{\overline{G^{1/2}}^{z^{*}}}
		       \left(\partial_{z^{*}} K_{m}^{2}\right)
              \Biggr\}
  .
  \end{eqnarray}

 水平方向の運動方程式 \Deqref{離散水平運動u}と
 \Deqref{離散水平運動v}を計算して
 $u^{\tau + \Delta \tau}, v^{\tau + \Delta \tau}$ を求める.
 次に\Deqref{離散鉛直運動}と\Deqref{離散圧力}から
  $w^{\tau + \Delta \tau}$を消去して$\pi^{\tau + \Delta \tau}$ 
 を計算して求める.
 最後に $w^{\tau + \Delta \tau}$ を求める.

\pagebreak
 鉛直方向の運動方程式を離散化した式\Deqref{離散鉛直運動}より
 \begin{eqnarray}
  w_{i, j, k + \frac{1}{2}}^{\tau + \Delta \tau}
     &=& w_{i, j, k + \frac{1}{2}}^{\tau}
         - \left(\frac{c_{p} \bar{\theta}_{v} \Delta \tau}
		      {\overline{G^{1/2}}^{z^{*}}}
	   \right)_{i, j,k + \frac{1}{2}}
	   \left\{- \alpha
	          \left(\partial_{z^{*}} Div \right)
	          _{i, j, k + \frac{1}{2}}^{\tau}
		  + \left(1 - \beta\right)
		  \left(\partial_{z^{*}} \pi \right)
		  _{i, j, k + \frac{1}{2}}^{\tau}
	   \right\}
  \nonumber \\
       & & + F.w_{i, j, k + \frac{1}{2}}^{t} \Delta \tau
	   - \beta
	   \left(\frac{c_{p} \bar{\theta}_{v} \Delta \tau}
	              {\overline{G^{1/2}}^{z^{*}}}
	   \right)_{i, j,k + \frac{1}{2}}
	   \left(\partial_{z^{*}} \pi \right)
	   _{i, j, k + \frac{1}{2}}
	   ^{\tau + \Delta \tau}
 \Deqlab{離散鉛直運動2}
 \end{eqnarray}
 ただし
 \begin{eqnarray}
  F.w_{i, j, k + \frac{1}{2}}^{t}
          = - {\rm Adv}.w_{i, j, k + \frac{1}{2}}^{t}
            + {\rm Buoy}.w_{i, j, k + \frac{1}{2}}^{t}
	    + \left(D_{w}\right)_{i, j, k + \frac{1}{2}}
	      ^{t - \Delta t}
 \end{eqnarray}
 である.

 離散化した圧力方程式\Deqref{離散圧力}より,
 \begin{eqnarray}
  \pi_{i,j,k}^{\tau + \Delta \tau}
      &+& \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,k}
	  \left\{\partial_{z^{*}}
	         \left(\bar{\rho} \bar{\theta}_{v} w\right)
	  \right\}_{i,j,k}^{\tau + \Delta \tau}
  \nonumber \\
      &=& \pi_{i,j,k}^{\tau} 
         - \left(1 - \beta\right)
	   \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,k}
	  \left\{\partial_{z^{*}}
	         \left(\bar{\rho} \bar{\theta}_{v} w\right)
	  \right\}_{i,j,k}^{\tau}
  \nonumber \\
     & & - \frac{\bar{c}^{2} \Delta \tau}
	        {c_{p} \bar{\theta}_{v}}
	   \Biggl[\left(\partial_{x} u\right)_{i,j,k}
		  + \left(\partial_{y}v\right)_{i,j,k}
	          + \overline{G^{13}}_{i,j,k}^{xz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{u}^{x z^{*}}
		    \right)_{i,j,k}
		  + \overline{G^{23}}_{i,j,k}^{yz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{v}^{y z^{*}}
		    \right)_{i,j,k}
	   \Biggr]^{\tau + \Delta \tau}
  .
 \end{eqnarray}
 \Deqref{離散鉛直運動2}を代入して $w^{\tau + \Delta \tau}$ を
 消去すると, 
 \begin{eqnarray}
  \pi_{i,j,k}^{\tau + \Delta \tau}
      &-& \beta^{2}
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau^{2}}
		     {\bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,k}
	  \left[\partial_{z^{*}}
	        \left\{\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                          {G^{1/2}}
			\right)
		        \left(\partial_{z^{*}} \pi\right)
		\right\}
	  \right]_{i,j,k}^{\tau + \Delta \tau} 
 \nonumber \\
      &=& \pi_{i,j,k}^{\tau} 
         - \left(1 - \beta\right)
	   \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,k}
	  \left\{\partial_{z^{*}}
	         \left(\bar{\rho} \bar{\theta}_{v} w\right)
	  \right\}_{i,j,k}^{\tau}
  \nonumber \\
     & & - \frac{\bar{c}^{2} \Delta \tau}
	        {c_{p} \bar{\theta}_{v}}
	   \Biggl[\left(\partial_{x} u\right)_{i,j,k}
		  + \left(\partial_{y}v\right)_{i,j,k}
	          + \overline{G^{13}}_{i,j,k}^{xz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{u}^{x z^{*}}
		    \right)_{i,j,k}
		  + \overline{G^{23}}_{i,j,k}^{yz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{v}^{y z^{*}}
		    \right)_{i,j,k}
	   \Biggr]^{\tau + \Delta \tau}
  \nonumber \\
      & & - \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,k}
 \nonumber \\
      & & \hspace*{.5em}
          \times
	  \Biggl[\partial_{z^{*}}
	        \bar{\rho} \bar{\theta}_{v}
	        \Biggl\{
	        w^{\tau}
         - \left(\frac{c_{p} \bar{\theta}_{v} \Delta \tau}
		      {\overline{G^{1/2}}^{z^{*}}}
	   \right)
	   \left\{- \alpha
	          \left(\partial_{z^{*}} Div \right)^{\tau}
		  + \left(1 - \beta\right)
		  \left(\partial_{z^{*}} \pi \right)^{\tau}
	   \right\}
		- F.w^{t} \Delta \tau
		  \Biggr\}
	  \Biggr]_{i,j,k}
  \Deqlab{圧力の式}
  .
 \end{eqnarray}
 \pagebreak
 さらに, 左辺を $z^{*}$ 方向の微分をオペレータなしで書き下すと,
 \begin{eqnarray}
  {\rm (左辺)} &=& \pi_{i,j,k}^{\tau + \Delta \tau}
  \nonumber \\
         & & - \beta^{2}
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau^{2}}
		     {\bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,k}
  \nonumber \\
       & & \hspace*{2em}
	   \times
	   \frac{1}{\Delta z^{*}}
	  \Biggl[\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                   {G^{1/2}}
		\right)_{i,j,k}
	        \frac{1}{\Delta z^{*}}
		\left(\pi_{i,j,k+1} - \pi_{i,j,k}\right)
  \nonumber \\
       & & \hspace*{6em}
		- \left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                      {G^{1/2}}
		  \right)_{i,j,k-1}
	          \frac{1}{\Delta z^{*}}
		  \left(\pi_{i,j,k} - \pi_{i,j,k-1}\right)
	  \Biggr]^{\tau + \Delta \tau} 
  \nonumber \\
        &=& \pi_{i,j,k}^{\tau + \Delta \tau}
  \nonumber \\
         & & - \beta^{2}
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau^{2}}
		     {\bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,k}
	   \frac{1}{\Delta z^{* 2}}
  \nonumber \\
       & & \hspace*{1em}
	   \times
	  \Biggl[\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                   {G^{1/2}}
		\right)_{i,j,k}
	        \pi_{i,j,k+1}
		- \left\{\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                   {G^{1/2}}
		       \right)_{i,j,k}
		       + \left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                            {G^{1/2}}
			 \right)_{i,j,k-1}
		\right\}
		\pi_{i,j,k}
  \nonumber \\
       & & \hspace*{8em}
		+ \left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                      {G^{1/2}}
		  \right)_{i,j,k-1}
		  \pi_{i,j,k-1}
	  \Biggr]^{\tau + \Delta \tau} 
  .
 \end{eqnarray}
 
 鉛直陰解法における上部と下部の境界条件は, それぞれ
 \begin{eqnarray}
  w &=& 0 \hspace{7em} (k = 0)
  , \\
  \omega &=& u \DP{z_{s}}{x} + v \DP{z_{s}}{y}
    \hspace{2em} (k = km + 1)
 \end{eqnarray}
 である.
 上部境界における鉛直方向の運動方程式を離散化したものは
 \begin{eqnarray}
	   \beta
	   \left(\partial_{z^{*}} \pi \right)
	   _{i, j, w(km+1)}
	   ^{\tau + \Delta \tau}
     &=& \alpha
         \left(\partial_{z^{*}} Div \right)
	 _{i, j, w(km+1)}^{\tau}
	 - \left(1 - \beta\right)
	   \left(\partial_{z^{*}} \pi \right)
	   _{i, j, w(km+1)}^{\tau}
  \nonumber \\
      & &  + \left(\frac{G^{1/2}}
	                {c_{p} \bar{\theta}_{v}}
			F.w^{t}
	     \right)_{i, j, w(km+1)}
 \end{eqnarray}
 となり, 圧力方程式は
\pagebreak
 \begin{eqnarray}
  {\rm (左辺)} &=& \pi_{i,j,km}^{\tau + \Delta \tau}
  \nonumber \\
         & & + \beta^{2}
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau^{2}}
		     {\bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,km}
	   \frac{1}{\Delta z^{* 2}}
	  \Biggl[\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                      {G^{1/2}}
		  \right)_{i,j,w(km)}
		  \left(\pi_{i,j,km} - \pi_{i,j,km-1}\right)
	  \Biggr]^{\tau + \Delta \tau} 
  \nonumber \\
  {\rm (右辺)} &=& \pi_{i,j,km}^{\tau}
         - \left(1 - \beta\right)
	   \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,km}
	  \left\{\partial_{z^{*}}
	         \left(\bar{\rho} \bar{\theta}_{v} w\right)
	  \right\}_{i,j,km}^{\tau}
  \nonumber \\
     & & - \frac{\bar{c}^{2} \Delta \tau}
	        {c_{p} \bar{\theta}_{v}}
	   \Biggl[\left(\partial_{x} u\right)_{i,j,km}
		  + \left(\partial_{y}v\right)_{i,j,km}
 \nonumber \\
      & & \hspace*{5em}
	          + \overline{G^{13}}_{i,j,km}^{xz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{u}^{x z^{*}}
		    \right)_{i,j,km}
		  + \overline{G^{23}}_{i,j,km}^{yz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{v}^{y z^{*}}
		    \right)_{i,j,km}
	   \Biggr]^{\tau + \Delta \tau}
  \nonumber \\
      & & - \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,km}
 \nonumber \\
      & & \hspace*{.5em}
          \times
	  \Biggl[\partial_{z^{*}}
	        \bar{\rho} \bar{\theta}_{v}
	        \Biggl\{
	        w^{\tau}
         - \left(\frac{c_{p} \bar{\theta}_{v} \Delta \tau}
		      {\overline{G^{1/2}}^{z^{*}}}
	   \right)
	   \left\{- \alpha
	          \left(\partial_{z^{*}} Div \right)^{\tau}
		  + \left(1 - \beta\right)
		  \left(\partial_{z^{*}} \pi \right)^{\tau}
	   \right\}
		- F.w^{t} \Delta \tau
		  \Biggr\}
	  \Biggr]_{i,j,km}
  \nonumber \\
         & & + \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau^{2}}
		     {\bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,km}
	   \frac{1}{\Delta z^{*}}
	   \left(\frac{c_{p} \bar{\theta}_{v}^{2}}{G^{1/2}}
	   \right)_{i,j,w(km+1)}
  \nonumber \\
      & & \hspace*{.5em}
          \times
	   \Biggl\{\alpha
           \left(\partial_{z^{*}} Div \right)
	   _{i, j, w(km+1)}^{\tau}
	   - \left(1 - \beta\right)
	     \left(\partial_{z^{*}} \pi \right)
	     _{i, j, w(km+1)}^{\tau}
          + \left(\frac{G^{1/2}}
	                {c_{p} \bar{\theta}_{v}}
			F.w^{t}
	     \right)_{i, j, w(km+1)}
	   \Biggr\}
 \end{eqnarray}
 となる.

 また, 下部境界における鉛直方向の運動方程式を離散化したものは
 \begin{eqnarray}
	   \beta
	   \left(\partial_{z^{*}} \pi \right)
	   _{i, j, w(0)}
	   ^{\tau + \Delta \tau}
     &=& \alpha
         \left(\partial_{z^{*}} Div \right)
	 _{i, j, w(0)}^{\tau}
	 - \left(1 - \beta\right)
	   \left(\partial_{z^{*}} \pi \right)
	   _{i, j, w(0)}^{\tau}
  \nonumber \\
      & &  + \left(\frac{G^{1/2}}
	                {c_{p} \bar{\theta}_{v}}
			F.w^{t}
	     \right)_{i, j, w(0)}
 \end{eqnarray}
 となり, 圧力方程式は
\pagebreak
 \begin{eqnarray}
  {\rm (左辺)} &=& \pi_{i,j,0}^{\tau + \Delta \tau}
  \nonumber \\
         & & - \beta^{2}
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau^{2}}
		     {\bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,0}
	   \frac{1}{\Delta z^{* 2}}
	  \Biggl[\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                      {G^{1/2}}
		  \right)_{i,j,w(1)}
	  \left(\pi_{i,j,1} - \pi_{i,j,0}\right)
	  \Biggr]^{\tau + \Delta \tau} 
  \nonumber \\
  {\rm (右辺)} &=& \pi_{i,j,0}^{\tau} 
         - \left(1 - \beta\right)
	   \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,0}
	  \left\{\partial_{z^{*}}
	         \left(\bar{\rho} \bar{\theta}_{v} w\right)
	  \right\}_{i,j,0}^{\tau}
  \nonumber \\
     & & - \frac{\bar{c}^{2} \Delta \tau}
	        {c_{p} \bar{\theta}_{v}}
	   \Biggl[\left(\partial_{x} u\right)_{i,j,0}
		  + \left(\partial_{y}v\right)_{i,j,0}
 \nonumber \\
      & & \hspace*{5em}
	          + \overline{G^{13}}_{i,j,0}^{xz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{u}^{x z^{*}}
		    \right)_{i,j,1}
		  + \overline{G^{23}}_{i,j,0}^{yz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{v}^{y z^{*}}
		    \right)_{i,j,0}
	   \Biggr]^{\tau + \Delta \tau}
  \nonumber \\
      & & - \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,0}
 \nonumber \\
      & & \hspace*{.5em}
          \times
	  \Biggl[\partial_{z^{*}}
	        \bar{\rho} \bar{\theta}_{v}
	        \Biggl\{
	        w^{\tau}
         - \left(\frac{c_{p} \bar{\theta}_{v} \Delta \tau}
		      {\overline{G^{1/2}}^{z^{*}}}
	   \right)
	   \left\{- \alpha
	          \left(\partial_{z^{*}} Div \right)^{\tau}
		  + \left(1 - \beta\right)
		  \left(\partial_{z^{*}} \pi \right)^{\tau}
	   \right\}
		- F.w^{t} \Delta \tau
		  \Biggr\}
	  \Biggr]_{i,j,0}
  \nonumber \\
         & & - \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau^{2}}
		     {\bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,1}
	   \frac{1}{\Delta z^{*}}
	   \left(\frac{c_{p} \bar{\theta}_{v}^{2}}{G^{1/2}}
	   \right)_{i,j,w(0)}
  \nonumber \\
      & & \hspace*{.5em}
          \times
	   \Biggl\{\alpha
           \left(\partial_{z^{*}} Div \right)
	   _{i, j, w(0)}^{\tau}
	   - \left(1 - \beta\right)
	     \left(\partial_{z^{*}} \pi \right)
	     _{i, j, w(0)}^{\tau}
          + \left(\frac{G^{1/2}}
	                {c_{p} \bar{\theta}_{v}}
			F.w^{t}
	     \right)_{i, j, w(0)}
	   \Biggr\}
 \end{eqnarray}

 これらの式は以下のように行列を用いて書くことができる:
 \begin{eqnarray}
  \left(\begin{array}{cccc}
          A_{0} & B_{0}  &        &   0      \\ 
	  C_{0} & \ddots & \ddots &          \\
	        & \ddots & \ddots & B_{km-1}   \\
	    0   &        & C_{km-1} & A_{km} \\
	\end{array}
  \right)
  \left(\begin{array}{cccc}
          \pi_{00} & \pi_{02}  & \cdots & \pi_{0 \, im} \\
	  \pi_{20} & \ddots    & \ddots & \vdots        \\
	  \vdots   &           & \ddots & \vdots        \\
	  \pi_{km \, 0} & \cdots & \cdots & \pi_{km \, im} \\
	\end{array}
  \right)^{\tau + \Delta \tau}
  \nonumber \\
  = 
  \left(\begin{array}{cccc}
          D_{00} & D_{02}  & \cdots & D_{0 \, im} \\
	  D_{20} & \ddots    & \ddots & \vdots        \\
	  \vdots   &           & \ddots & \vdots        \\
	  D_{km \, 0} & \cdots & \cdots & D_{km \, im} \\
	\end{array}
  \right)^{\tau}
  .
 \end{eqnarray}
 \pagebreak
 ただし
 \begin{eqnarray}
  A_{k} &=& 1 + \beta^{2}
                \left(\frac{1}{G^{1/2}}
	              \frac{\bar{c}^{2} \Delta \tau^{2}}
		           {\bar{\rho} \bar{\theta}_{v}^{2}}
		\right)_{i,j,k}
		\frac{1}{\Delta z^{* 2}}
		\left\{\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                          {G^{1/2}}
		       \right)_{i,j,k}
		       + \left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                            {G^{1/2}}
			 \right)_{i,j,k-1}
		\right\}
  , \nonumber \\
        & & \hspace*{4em} (k = 1,3, \cdots, km-1)
  , \nonumber \\
  A_{0} &=& 1 + \beta^{2}
                \left(\frac{1}{G^{1/2}}
	              \frac{\bar{c}^{2} \Delta \tau^{2}}
		           {\bar{\rho} \bar{\theta}_{v}^{2}}
		\right)_{i,j,0}
		\frac{1}{\Delta z^{* 2}}
		\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                          {G^{1/2}}
		\right)_{i,j,w(1)}
  , \nonumber \\
  A_{km} &=& 1 + \beta^{2}
                \left(\frac{1}{G^{1/2}}
	              \frac{\bar{c}^{2} \Delta \tau^{2}}
		           {\bar{\rho} \bar{\theta}_{v}^{2}}
		\right)_{i,j,km}
		\frac{1}{\Delta z^{* 2}}
		\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                          {G^{1/2}}
		\right)_{i,j,w(km)}
  , \nonumber \\
  B_{k} &=& - \beta^{2}
            \left(\frac{1}{G^{1/2}}
	          \frac{\bar{c}^{2} \Delta \tau^{2}}
	               {\bar{\rho} \bar{\theta}_{v}^{2}}
	    \right)_{i,j,k}
	    \frac{1}{\Delta z^{* 2}}
		\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                          {G^{1/2}}
		\right)_{i,j,w(k+1)}
  , \nonumber \\
  C_{k} &=& - \beta^{2}
            \left(\frac{1}{G^{1/2}}
	          \frac{\bar{c}^{2} \Delta \tau^{2}}
	               {\bar{\rho} \bar{\theta}_{v}^{2}}
	    \right)_{i,j,k+1}
	    \frac{1}{\Delta z^{* 2}}
		\left(\frac{\bar{\rho} \overline{\theta}_{v}^{2}}
	                          {G^{1/2}}
		\right)_{i,j,w(k+1)}
  , \nonumber \\
  D_{ik}       &=& \pi_{i,j,k}^{\tau} 
         - \left(1 - \beta\right)
	   \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,k}
	  \left\{\partial_{z^{*}}
	         \left(\bar{\rho} \bar{\theta}_{v} w\right)
	  \right\}_{i,j,k}^{\tau}
  \nonumber \\
     & & - \frac{\bar{c}^{2} \Delta \tau}
	        {c_{p} \bar{\theta}_{v}}
	   \Biggl[\left(\partial_{x} u\right)_{i,j,k}
		  + \left(\partial_{y}v\right)_{i,j,k}
	          + \overline{G^{13}}_{i,j,k}^{xz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{u}^{x z^{*}}
		    \right)_{i,j,k}
		  + \overline{G^{23}}_{i,j,k}^{yz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{v}^{y z^{*}}
		    \right)_{i,j,k}
	   \Biggr]^{\tau + \Delta \tau}
  \nonumber \\
      & & - \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,k}
 \nonumber \\
      & & \hspace*{.5em}
          \times
	  \Biggl[\partial_{z^{*}}
	        \bar{\rho} \bar{\theta}_{v}
	        \Biggl\{
	        w^{\tau}
         - \left(\frac{c_{p} \bar{\theta}_{v} \Delta \tau}
		      {\overline{G^{1/2}}^{z^{*}}}
	   \right)
	   \left\{- \alpha
	          \left(\partial_{z^{*}} Div \right)^{\tau}
		  + \left(1 - \beta\right)
		  \left(\partial_{z^{*}} \pi \right)^{\tau}
	   \right\}
		- F.w^{t} \Delta \tau
		  \Biggr\}
	  \Biggr]_{i,j,k}
  , \nonumber \\
        & & \hspace*{4em} (k = 1,3, \cdots, km-3)
  , \nonumber \\
  D_{i0} &=& \pi_{i,j,0}^{\tau} 
         - \left(1 - \beta\right)
	   \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,0}
	  \left\{\partial_{z^{*}}
	         \left(\bar{\rho} \bar{\theta}_{v} w\right)
	  \right\}_{i,j,0}^{\tau}
  \nonumber \\
     & & - \frac{\bar{c}^{2} \Delta \tau}
	        {c_{p} \bar{\theta}_{v}}
	   \Biggl[\left(\partial_{x} u\right)_{i,j,0}
		  + \left(\partial_{y}v\right)_{i,j,0}
 \nonumber \\
      & & \hspace*{5em}
	          + \overline{G^{13}}_{i,j,0}^{xz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{u}^{x z^{*}}
		    \right)_{i,j,0}
		  + \overline{G^{23}}_{i,j,1}^{yz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{v}^{y z^{*}}
		    \right)_{i,j,0}
	   \Biggr]^{\tau + \Delta \tau}
  \nonumber \\
      & & - \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,0}
 \nonumber \\
      & & \hspace*{.5em}
          \times
	  \Biggl[\partial_{z^{*}}
	        \bar{\rho} \bar{\theta}_{v}
	        \Biggl\{
	        w^{\tau}
         - \left(\frac{c_{p} \bar{\theta}_{v} \Delta \tau}
		      {\overline{G^{1/2}}^{z^{*}}}
	   \right)
	   \left\{- \alpha
	          \left(\partial_{z^{*}} Div \right)^{\tau}
		  + \left(1 - \beta\right)
		  \left(\partial_{z^{*}} \pi \right)^{\tau}
	   \right\}
		- F.w^{t} \Delta \tau
		  \Biggr\}
	  \Biggr]_{i,j,0}
  \nonumber \\
         & & - \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau^{2}}
		     {\bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,0}
	   \frac{1}{\Delta z^{*}}
	   \left(\frac{c_{p} \bar{\theta}_{v}^{2}}{G^{1/2}}
	   \right)_{i,j,w(0)}
  \nonumber \\
      & & \hspace*{.5em}
          \times
	   \Biggl\{\alpha
           \left(\partial_{z^{*}} Div \right)
	   _{i, j, w(0)}^{\tau}
	   - \left(1 - \beta\right)
	     \left(\partial_{z^{*}} \pi \right)
	     _{i, j, w(0)}^{\tau}
          + \left(\frac{G^{1/2}}
	                {c_{p} \bar{\theta}_{v}}
			F.w^{t}
	     \right)_{i, j, w(0)}
	   \Biggr\}
  , \nonumber \\
  D_{ikm-2} &=& \pi_{i,j,km-2}^{\tau} 
         - \left(1 - \beta\right)
	   \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,km-2}
	  \left\{\partial_{z^{*}}
	         \left(\bar{\rho} \bar{\theta}_{v} w\right)
	  \right\}_{i,j,km-2}^{\tau}
  \nonumber \\
     & & - \frac{\bar{c}^{2} \Delta \tau}
	        {c_{p} \bar{\theta}_{v}}
	   \Biggl[\left(\partial_{x} u\right)_{i,j,km-2}
		  + \left(\partial_{y}v\right)_{i,j,km-2}
 \nonumber \\
      & & \hspace*{5em}
	          + \overline{G^{13}}_{i,j,km}^{xz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{u}^{x z^{*}}
		    \right)_{i,j,km-2}
		  + \overline{G^{23}}_{i,j,km}^{yz^{*}}
		    \left(\partial_{z^{*}}
		          \overline{v}^{y z^{*}}
		    \right)_{i,j,km-2}
	   \Biggr]^{\tau + \Delta \tau}
  \nonumber \\
      & & - \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau}
		     {c_{p} \bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,km-2}
 \nonumber \\
      & & \hspace*{.5em}
          \times
	  \Biggl[\partial_{z^{*}}
	        \bar{\rho} \bar{\theta}_{v}
	        \Biggl\{
	        w^{\tau}
         - \left(\frac{c_{p} \bar{\theta}_{v} \Delta \tau}
		      {\overline{G^{1/2}}^{z^{*}}}
	   \right)
	   \left\{- \alpha
	          \left(\partial_{z^{*}} Div \right)^{\tau}
		  + \left(1 - \beta\right)
		  \left(\partial_{z^{*}} \pi \right)^{\tau}
	   \right\}
		- F.w^{t} \Delta \tau
		  \Biggr\}
	  \Biggr]_{i,j,km-2}
  \nonumber \\
         & & + \beta
          \left(\frac{1}{G^{1/2}}
	        \frac{\bar{c}^{2} \Delta \tau^{2}}
		     {\bar{\rho} \bar{\theta}_{v}^{2}}
	  \right)_{i,j,km-2}
	   \frac{1}{\Delta z^{*}}
	   \left(\frac{c_{p} \bar{\theta}_{v}^{2}}{G^{1/2}}
	   \right)_{i,j,w(km-1)}
  \nonumber \\
      & & \hspace*{.5em}
          \times
	   \Biggl\{\alpha
           \left(\partial_{z^{*}} Div \right)
	   _{i, j, w(km-1)}^{\tau}
	   - \left(1 - \beta\right)
	     \left(\partial_{z^{*}} \pi \right)
	     _{i, j, w(km-1)}^{\tau}
          + \left(\frac{G^{1/2}}
	                {c_{p} \bar{\theta}_{v}}
			F.w^{t}
	     \right)_{i, j, w(km-1)}
	   \Biggr\}
 \end{eqnarray}

 この連立 1 次方程式の係数行列は 3 項行列なので, LU 分解を行って
 計算すれば大幅に計算量を減らすことができる.

 したがって, 求めた $\pi^{\tau + \Delta \tau}$ と鉛直方向の
 運動方程式から $w^{\tau + \Delta \tau}$ を求めることができる.

%--------------------------------------------------------------------
\pagebreak

\subsection{熱力学方程式と各混合比の方程式の離散化}

 温位偏差 $(\theta)$と水蒸気混合比 $(q_{v})$, 雲水の
 混合比 $(q_{c})$, 雨水の混合比 $(q_{r})$について,
 長い時間間隔 $\Delta t$ で離散化する.
 過飽和がない場合を考えているので, 予報方程式から得られた
 水蒸気の混合比と飽和混合比の関係が
 $q_{v}^{t + \Delta t} > q_{vs}^{t + \Delta t}$  の場合は
 湿潤飽和調節法\footnote{\ref{湿潤飽和調節法}参照.}を用いる.

 これらの物理量は格子点 $(i, j, k)$ についてなので,
 ここでは物理量の下付添字 $i, j, k$ は省略する.

 {\bf 熱力学方程式}:
 \begin{eqnarray}
  \frac{\theta^{t + \Delta t}
      - \theta^{t - \Delta t}}{2 \Delta t}
      = - {\rm Adv}.\theta^{t}
        + M_{\theta}^{t}
        + D_{\theta}^{t - \Delta t}
 \end{eqnarray}

 {\bf 水蒸気の混合比の方程式}:
 \begin{eqnarray}
  \frac{q_{v}^{t + \Delta t}
      - q_{v}^{t - \Delta t}}{2 \Delta t}
      = - {\rm Adv}.q_{v}^{t}
        + M_{q_{v}}^{t}
	+ D_{q_{v}}^{t - \Delta t}
 \end{eqnarray}

 {\bf 雲水の混合比の方程式}:
 \begin{eqnarray}
  \frac{q_{c}^{t + \Delta t}
      - q_{c}^{t - \Delta t}}{2 \Delta t}
      = - {\rm Adv}.q_{c}^{t}
        + M_{q_{c}}^{t}
	+ D_{q_{c}}^{t - \Delta t}
 \end{eqnarray}

 {\bf 雨水の混合比の方程式}:
 \begin{eqnarray}
  \frac{q_{r}^{t + \Delta t}
      - q_{r}^{t - \Delta t}}{2 \Delta t}
      = - {\rm Adv}.q_{r}^{t}
        + M_{q_{r}}^{t}
	+ D_{q_{r}}^{t - \Delta t}
 \end{eqnarray}

 長い時間間隔における Leap-frog 法による積分の場合,
 それに伴う計算モードを抑えるため Asselin の時間フィルタ
 (Asselin, 1972) を用いる\footnote{\ref{asselin}節参照.}.
 時間 $t + \Delta t$ の物理量が得られた後, 時間 $t$ の物理
 量を以下のようにフィルタリングする:
 \begin{eqnarray}
  \phi^{t} &=& \phi^{* t}
               + \mu_{a}
	         \left(\phi^{* t + \Delta t} - 2 \phi^{* t}
		       + \phi^{t - \Delta t}
		 \right)
  .
 \end{eqnarray}
 ただし $\mu_{a}$ はフィルターの係数で, $\mu_{a} = 0.1$ が
 標準的な値である.

 移流項:
 \begin{eqnarray}
  {\rm Adv}.\theta^{t}
      &=& \left[\overline{u \partial_{x} \theta}^{x} 
	        + \overline{v \partial_{y} \theta}^{y}
		+ \overline{\omega \partial_{z^{*}}
		  \theta}^{z^{*}}
	  \right]^{t}
	  + {\rm Diff}.\theta^{t}
  , \\
  {\rm Adv}.q_{v}^{t}
      &=& \left[\overline{u \partial_{x} q_{v}}^{x} 
	        + \overline{v \partial_{y} q_{v}}^{y}
		+ \overline{\omega \partial_{z^{*}}
		  q_{v}}^{z^{*}}
	  \right]^{t}
	  + {\rm Diff}.q_{v}^{t}
  , \\
  {\rm Adv}.q_{c}^{t}
      &=& \left[\overline{u \partial_{x} q_{c}}^{x} 
	        + \overline{v \partial_{y} q_{c}}^{y}
		+ \overline{\omega \partial_{z^{*}}
		  q_{c}}^{z^{*}}
	  \right]^{t}
          + {\rm Diff}.q_{c}^{t}
  , \\
  {\rm Adv}.q_{r}^{t}
      &=& \left[\overline{u \partial_{x} q_{r}}^{x} 
	        + \overline{v \partial_{y} q_{r}}^{y}
		+ \overline{\omega \partial_{z^{*}}
		  q_{r}}^{z^{*}}
	  \right]^{t}
          + {\rm Diff}.q_{r}^{t}
  .
 \end{eqnarray}

 数値粘性項:
 \begin{eqnarray}
  {\rm Diff}.\theta^{t}
      &=& \nu_{2h}
          \left[\partial_{x}^{2} \theta
	        + \partial_{y}^{2} \theta
	  \right]^{t}
	  + \nu_{2v}
	     \left[\partial_{z^{*}}^{2} \theta\right]^{t}
  , \\
  {\rm Diff}.q_{v}^{t}
      &=& \nu_{2h}
          \left[\partial_{x}^{2} q_{v}
	        + \partial_{y}^{2} q_{v}
	  \right]^{t}
	  + \nu_{2v}
	     \left[\partial_{z^{*}}^{2} q_{v}\right]^{t}
  , \\
  {\rm Diff}.q_{c}^{t}
      &=& \nu_{2h}
          \left[\partial_{x}^{2} q_{c}
	        + \partial_{y}^{2} q_{c}
	  \right]^{t}
	  + \nu_{2v}
	     \left[\partial_{z^{*}}^{2} q_{c}\right]^{t}
  , \\
  {\rm Diff}.q_{r}^{t}
      &=& \nu_{2h}
          \left[\partial_{x}^{2} q_{r}
	        + \partial_{y}^{2} q_{r}
	  \right]^{t}
	  + \nu_{2v}
	     \left[\partial_{z^{*}}^{2} q_{r}\right]^{t}
 \end{eqnarray}

 生成・消滅項:
 \begin{eqnarray}
  M_{\theta}^{t} &=& - \gamma
                    \left(\frac{q_{vs}^{t + \Delta t}
	                  - q_{vs}^{t - \Delta t}}
			  {2 \Delta t}
			  + {\rm Adv}.q_{vs}^{t}
			  + E_{r}^{t}
		    \right)
  , \\
  M_{q_{v}}^{t} &=& \frac{q_{vs}^{t + \Delta t}
		          - q_{vs}^{t - \Delta t}}
			 {2 \Delta t}
		    + {\rm Adv}.q_{vs}^{t - \Delta t}
		    + E_{r}^{t}
  , \\
  M_{q_{c}}^{t} &=& - \frac{q_{vs}^{t + \Delta t}
		            - q_{vs}^{t - \Delta t}}
			   {2 \Delta t}
		    - {\rm Adv}.q_{vs}^{t}
		    - A_{r}^{t}
		    - C_{r}^{t}
  , \\
  M_{q_{r}}^{t} &=& F_{r}^{t}
                    - E_{r}^{t}
                    + A_{r}^{t}
		    + C_{r}^{t}
 \end{eqnarray}
 ただし
 \begin{eqnarray}
  {\rm Adv}.q_{vs}^{t}
      &=& \left[\overline{u \partial_{x} q_{vs}}^{x} 
	        + \overline{v \partial_{y} q_{vs}}^{y}
		+ \overline{\omega \partial_{z^{*}}
		q_{vs}}^{z^{*}}
	  \right]^{t}
	+ {\rm Diff}.q_{vs}^{t}
  , \\
  q_{vs} &=& \frac{3.8}{\bar{p}}
           \exp
	   \left(17.27
	    \frac{\bar{\Pi} \theta -273}
	    {\bar{\Pi} \theta -36}\right)
  \nonumber \\ 
	 &\simeq& q_{vs}^{*}
	          \left[1
		        + \frac{4093 a \bar{\Pi}}
			       {\left(\bar{\Pi} \theta^{*}
				- 36\right)^{2}}
			  \left(\theta^{t + \Delta t}
			        - \theta^{*}
			  \right)
		  \right]
  ,
 \end{eqnarray}
 $q_{vs}$ は Tetens の式より $\theta = \theta^{*}$ について
 線形化を行った.
 $q_{vs}^{*}$ は $\theta = \theta^{*}$ のときの $q_{vs}$ の
 値である.
 \begin{eqnarray}
  F_{r}^{t}
       &=& \left[\frac{1}{G^{1/2} \bar{\rho}}
	         \partial_{z^{*}}
		 \left(\bar{\rho} V q_{r}\right)
	   \right]^{t}
  , \\
  V &=& 36.34 \left(\bar{\rho} q_{r}\right)^{0.1346}
        \left(\frac{\bar{\rho}}{\rho_{0}}\right)^{-\frac{1}{2}}
 .
 \end{eqnarray}

\pagebreak
 拡散項:
 \begin{eqnarray}
  D_{\theta}^{t - \Delta t} &=& 3 \partial_{x}
                 \left\{\overline{K_{m}}^{x}
		        \left(\partial_{x} \theta\right)
		        + \overline{K_{m}}^{x}
			  \overline{G^{13}}^{z^{*}}
			  \left(\partial_{z^{*}}
			        \overline{\theta}^{xz^{*}}
			  \right)
	         \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \partial_{y}
                     \left\{\overline{K_{m}}^{y}
		            \left(\partial_{y} \theta\right)
		           + \overline{K_{m}}^{y}
			     \overline{G^{23}}^{z^{*}}
			     \left(\partial_{z^{*}}
			           \overline{\theta}^{yz^{*}}
			     \right)
		     \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \overline{G^{13}}^{xz^{*}}
	           \partial_{z^{*}}
                 \left\{\overline{K_{m}}^{z^{*}}
		        \left(\partial_{x}
			      \overline{\theta}^{xz^{*}}
		        \right)
		       + \overline{K_{m}}^{z^{*}}
		         \overline{G^{13}}^{x}
			 \left(\partial_{z^{*}} \theta\right)
		 \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \overline{G^{23}}^{yz^{*}}
	           \partial_{z^{*}}
		   \left\{\overline{K_{m}}^{z^{*}}
		          \left(\partial_{y}
			        \overline{\theta}^{yz^{*}}
			  \right)
		          + \overline{K_{m}}^{z^{*}}
			    \overline{G^{23}}^{y}
			    \left(\partial_{z^{*}} \theta\right)
		     \right\}^{t - \Delta t}
  \nonumber \\
             & & + \frac{3}{G^{1/2}} \partial_{z^{*}}
	           \left\{\frac{K_{m}}{\overline{G^{1/2}}^{z^{*}}}
		          \left(\partial_{z^{*}} \theta\right)
		   \right\}^{t - \Delta t}
  , \\
  D_{q_{v}}^{t - \Delta t} &=& 3 \partial_{x}
                 \left\{\overline{K_{m}}^{x}
		        \left(\partial_{x} q_{v}\right)
		        + \overline{K_{m}}^{x}
			  \overline{G^{13}}^{z^{*}}
			  \left(\partial_{z^{*}}
			        \overline{q_{v}}^{xz^{*}}
			  \right)
	         \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \partial_{y}
                     \left\{\overline{K_{m}}^{y}
		            \left(\partial_{y} q_{v}\right)
		           + \overline{K_{m}}^{y}
			     \overline{G^{23}}^{z^{*}}
			     \left(\partial_{z^{*}}
			           \overline{q_{v}}^{yz^{*}}
			     \right)
		     \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \overline{G^{13}}^{xz^{*}}
	           \partial_{z^{*}}
                 \left\{\overline{K_{m}}^{z^{*}}
		        \left(\partial_{x}
			      \overline{q_{v}}^{xz^{*}}
		        \right)
		       + \overline{K_{m}}^{z^{*}}
		         \overline{G^{13}}^{x}
			 \left(\partial_{z^{*}} q_{v}\right)
		 \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \overline{G^{23}}^{yz^{*}}
	           \partial_{z^{*}}
		   \left\{\overline{K_{m}}^{z^{*}}
		          \left(\partial_{y}
			        \overline{q_{v}}^{yz^{*}}
			  \right)
		          + \overline{K_{m}}^{z^{*}}
			    \overline{G^{23}}^{y}
			    \left(\partial_{z^{*}} q_{v}\right)
		     \right\}^{t - \Delta t}
  \nonumber \\
             & & + \frac{3}{G^{1/2}} \partial_{z^{*}}
	           \left\{\frac{K_{m}}{\overline{G^{1/2}}^{z^{*}}}
		          \left(\partial_{z^{*}} q_{v}\right)
		   \right\}^{t - \Delta t}
  , \\
  D_{q_{c}}^{t - \Delta t} &=& 3 \partial_{x}
                 \left\{\overline{K_{m}}^{x}
		        \left(\partial_{x} q_{c}\right)
		        + \overline{K_{m}}^{x}
			  \overline{G^{13}}^{z^{*}}
			  \left(\partial_{z^{*}}
			        \overline{q_{c}}^{xz^{*}}
			  \right)
	         \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \partial_{y}
                     \left\{\overline{K_{m}}^{y}
		            \left(\partial_{y} q_{c}\right)
		           + \overline{K_{m}}^{y}
			     \overline{G^{23}}^{z^{*}}
			     \left(\partial_{z^{*}}
			           \overline{q_{v}}^{yz^{*}}
			     \right)
		     \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \overline{G^{13}}^{xz^{*}}
	           \partial_{z^{*}}
                 \left\{\overline{K_{m}}^{z^{*}}
		        \left(\partial_{x}
			      \overline{q_{c}}^{xz^{*}}
		        \right)
		       + \overline{K_{m}}^{z^{*}}
		         \overline{G^{13}}^{x}
			 \left(\partial_{z^{*}} q_{c}\right)
		 \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \overline{G^{23}}^{yz^{*}}
	           \partial_{z^{*}}
		   \left\{\overline{K_{m}}^{z^{*}}
		          \left(\partial_{y}
			        \overline{q_{c}}^{yz^{*}}
			  \right)
		          + \overline{K_{m}}^{z^{*}}
			    \overline{G^{23}}^{y}
			    \left(\partial_{z^{*}} q_{c}\right)
		     \right\}^{t - \Delta t}
  \nonumber \\
             & & + \frac{3}{G^{1/2}} \partial_{z^{*}}
	           \left\{\frac{K_{m}}{\overline{G^{1/2}}^{z^{*}}}
		          \left(\partial_{z^{*}} q_{c}\right)
		   \right\}^{t - \Delta t}
  , \\
  D_{q_{r}}^{t - \Delta t} &=& 3 \partial_{x}
                 \left\{\overline{K_{m}}^{x}
		        \left(\partial_{x} q_{r}\right)
		        + \overline{K_{m}}^{x}
			  \overline{G^{13}}^{z^{*}}
			  \left(\partial_{z^{*}}
			        \overline{q_{r}}^{xz^{*}}
			  \right)
	         \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \partial_{y}
                     \left\{\overline{K_{m}}^{y}
		            \left(\partial_{y} q_{r}\right)
		           + \overline{K_{m}}^{y}
			     \overline{G^{23}}^{z^{*}}
			     \left(\partial_{z^{*}}
			           \overline{q_{r}}^{yz^{*}}
			     \right)
		     \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \overline{G^{13}}^{xz^{*}}
	           \partial_{z^{*}}
                 \left\{\overline{K_{m}}^{z^{*}}
		        \left(\partial_{x}
			      \overline{q_{r}}^{xz^{*}}
		        \right)
		       + \overline{K_{m}}^{z^{*}}
		         \overline{G^{13}}^{x}
			 \left(\partial_{z^{*}} q_{r}\right)
		 \right\}^{t - \Delta t}
  \nonumber \\
             & & + 3 \overline{G^{23}}^{yz^{*}}
	           \partial_{z^{*}}
		   \left\{\overline{K_{m}}^{z^{*}}
		          \left(\partial_{y}
			        \overline{q_{r}}^{yz^{*}}
			  \right)
		          + \overline{K_{m}}^{z^{*}}
			    \overline{G^{23}}^{y}
			    \left(\partial_{z^{*}} q_{r}\right)
		     \right\}^{t - \Delta t}
  \nonumber \\
             & & + \frac{3}{G^{1/2}} \partial_{z^{*}}
	           \left\{\frac{K_{m}}{\overline{G^{1/2}}^{z^{*}}}
		          \left(\partial_{z^{*}} q_{r}\right)
		   \right\}^{t - \Delta t}
  .
 \end{eqnarray}

%--------------------------------------------------------------------
\pagebreak

\subsection{サブグリッドスケールの運動エネルギーの方程式の離散化}

 サブグリッドスケールの運動エネルギーの方程式は長い時間間隔
 $t$ で計算する.
 これらの物理量は格子点 $(i, j, k)$ についてなので,
 ここでは物理量の下付添字 $i, j, k$ は省略する.

 \begin{eqnarray}
  \frac{E^{t + \Delta t} - E^{t - \Delta t}}{2 \Delta t}
      = &-& {\rm Adv}.E^{t}
  \nonumber \\
            & & + \frac{3 g C_{m} l E^{1/2}}{G^{1/2}}
		  \left(- A \partial_{z^{*}}
		          \overline{\theta_{e}}^{z^{*}}
		         + \partial_{z^{*}}
			   \overline{q_{l}}^{z^{*}}
		    \right)^{t}
  \nonumber \\
            & & + C_{m}^{2} l^{2}
	          \Bigg\{\left(\partial_{x} u
			        + \overline{G^{13}}^{xz^{*}}
				  \partial_{z^{*}}
				  \overline{u}^{xz^{*}}
			   \right)^{2}
  \nonumber \\
	      & & \hspace*{4em}
			   + \left(\partial_{y} v
			           + \overline{G^{23}}^{yz^{*}}
				     \partial_{z^{*}}
				     \overline{v}^{yz^{*}}
			     \right)^{2}
			   + \left(\frac{1}{G^{1/2}}
			           \partial_{z^{*}} w
			     \right)^{2}
		   \Biggr\}^{t}
  \nonumber \\
	      & & + \frac{C_{m}^{2} l^{2}}{2}
	            \Biggl\{\left(\partial_{x}
			          \overline{v}^{xy}
			        + \partial_{y}
			          \overline{u}^{xy}
				+ \overline{G^{13}}^{xz^{*}}
				  \partial_{z^{*}}
				  \overline{v}^{yz^{*}}
				+ \overline{G^{23}}^{yz^{*}}
				  \partial_{z^{*}}
				  \overline{u}^{xz^{*}}
			  \right)^{2}
  \nonumber \\
	      & & \hspace*{4em}
			+ \left(\partial_{x} w
				+ \overline{G^{13}}^{xz^{*}}
				  \partial_{z^{*}} w
			        + \frac{1}{G^{1/2}}
				  \partial_{z^{*}}
				  \overline{u}^{xz^{*}}
			  \right)^{2}
  \nonumber \\
	      & & \hspace*{4em}
			+ \left(\partial_{y}
			        \overline{w}^{yz^{*}}
				+ \overline{G^{23}}^{yz^{*}}
				  \partial_{z^{*}} w
			        + \frac{1}{G^{1/2}}
				  \partial_{z^{*}}
				  \overline{v}^{yz^{*}}
			  \right)^{2}
	        \Biggr\}^{t}
  \nonumber \\
              & & + \frac{C_{m}^{2} l^{2}}{2}
                   \Biggl\{\partial_{x}^{2} E
                          + \partial_{x}
                            \left(\overline{G^{13}}^{z^{*}}
                                  \partial_{z^{*}}
                                  \overline{E}^{xz^{*}}
                            \right)
  \nonumber \\
         & & \hspace*{6em}
                          + \partial_{z^{*}}
                            \left(\overline{G^{13}}^{x}
                                  \partial_{x}
                                  \overline{E}^{xz^{*}}
                            \right)
                          + \overline{G^{13}}^{xz^{*}}
                            \partial_{z^{*}}
                            \left(\overline{G^{13}}^{x}
                                  \partial_{z^{*}} E
                            \right)
  \nonumber \\
         & & \hspace*{4em}
                          + \partial_{y}^{2} E
                          + \partial_{y}
                            \left(\overline{G^{23}}^{z^{*}}
                                  \partial_{z^{*}}
                                  \overline{E}^{yz^{*}}
                            \right)
  \nonumber \\
         & & \hspace*{6em}
                          + \partial_{z^{*}}
                            \left(\overline{G^{23}}^{y}
                                  \partial_{y}
                                  \overline{E}^{yz^{*}}
                            \right)
                          + \overline{G^{23}}^{yz^{*}}
                            \partial_{z^{*}}
                            \left(\overline{G^{23}}^{y}
                                  \partial_{z^{*}} E
                            \right)
  \nonumber \\
        & & \hspace*{4em}
                          + \frac{1}{G^{1/2}}
                            \partial_{z^{*}}
                            \left(\overline{\frac{1}{G^{1/2}}}^{*}
                                  \partial_{z^{*}} E
                            \right)
             \Biggr]^{t}
  \nonumber \\
              & & + C_{m}^{2} l^{2}
	            \Biggl\{\left(\partial_{x}
			          \overline{E^{1/2}}^{x}
				  + \overline{G^{13}}^{xz^{*}}
				    \partial_{z^{*}}
				    \overline{E^{1/2}}^{z^{*}}
			    \right)^{2}
  \nonumber \\
              & & \hspace*{4em}
			    + \left(\partial_{y}
			            \overline{E^{1/2}}^{y}
				    + \overline{G^{23}}^{yz^{*}}
				      \partial_{z^{*}}
				      \overline{E^{1/2}}^{z^{*}}
			      \right)^{2}
			    + \left(\frac{1}{G^{1/2}}
				    \partial_{z^{*}}
				    \overline{E^{1/2}}^{z^{*}}
			      \right)^{2}
			     \Biggr\}^{t}
  \nonumber \\
	     & & - \frac{C_{\varepsilon}}{l}
	           \left(E^{3/2}\right)^{t}
  .
 \end{eqnarray}

 移流項:
 \begin{eqnarray}
  {\rm Adv}.E^{t}
      &=& \left[\overline{u \partial_{x} E}^{x} 
	        + \overline{v \partial_{y} E}^{y}
		+ \overline{\omega \partial_{z^{*}}
		  E}^{z^{*}}
	  \right]^{t}
	  + {\rm Diff}.E^{t}
 \end{eqnarray}

 数値粘性項:
 \begin{eqnarray}
    {\rm Diff}.E^{t}
      &=& \nu_{2h}
          \left[\partial_{x}^{2} E
	        + \partial_{y}^{2} E
	  \right]^{t}
	  + \nu_{2v}
	     \left[\partial_{z^{*}}^{2} E\right]^{t}
 \end{eqnarray}


%======================================================================
\newpage

\section{コード解説}

%--------------------------------------------------------------------
%\pagebreak

\subsection{計算手順}

 計算手順は以下のように行う:
 \begin{enumerate}
   \item 基本場を求める
	 \begin{enumerate}
	  \item 基本場の
		\begin{enumerate}
		 \item 地表面における圧力$\bar{p} (0) = p_{0}$ 
		 \item 温度の鉛直分布 $\bar{T} (z)$
		 \item 水蒸気の混合比の鉛直分布 $\bar{q}_{v} (z)$ 
		\end{enumerate}
		を与える

	  \item 地表面でのエクスナー関数を設定
		\begin{eqnarray*}
		 \bar{\Pi} (0) = 1
		\end{eqnarray*}
		
	  \item エクスナー関数より, 圧力 $\bar{p} (z)$ を求める
		\begin{eqnarray*}
		 \bar{p} (z) = \bar{p} (0) \times
		               \bar{\Pi}^{c_{p}/R_{d}} (z)
		\end{eqnarray*}
		\label{基本場導出始め}
		
	  \item 温位 $\bar{\theta} (z)$ を求める
		\begin{eqnarray*}
		 \bar{\theta} (z)
		     = \frac{\bar{T} (z)}
		            {\bar{\Pi} (z)}
		\end{eqnarray*}
		
	  \item 仮温位 $\bar{\theta}_{v} (z)$ を求める
		\begin{eqnarray*}
		 \bar{\theta}_{v} (z)
		     = \bar{\theta} (z) \times
		       \left(1 + 0.62 \bar{q}_{v} (z)\right)
		\end{eqnarray*}
		
	  \item 湿潤空気の状態方程式より,
		密度 $\bar{\rho} (z)$ を求める
		\begin{eqnarray*}
		 \bar{\rho} (z)
		     = \frac{\bar{p} (0)}{R_{d}}
		       \frac{\bar{\Pi} (z)}
		            {\bar{\theta}_{v}^{c_{v}/R_{d}} (z)}
		\end{eqnarray*}

	  \item 静水圧平衡の式より, エクスナー関数 (無次元量の圧力)
		の次の値 $\bar{\Pi} (z + dz)$ を求める
		\begin{eqnarray*}
		 \bar{\Pi} (z + dz)
		    = \bar{\pi} (z)
		      - dz \times
		        \frac{g}{c_{p} \bar{\theta}_{v} (z)}
		\end{eqnarray*}
		\label{基本場導出終わり}

	  \item \ref{基本場導出始め}から\ref{基本場導出終わり}
		までを繰り返して基本場を求める
	 \end{enumerate}

  \item 1 ステップ目の $\theta^{\Delta t}, q_{v}^{\Delta t},
	 q_{c}^{\Delta t}, q_{r}^{\Delta t}, E^{\Delta t}$
	 のみを前進差分で求める
	 \begin{itemize}
	  \item 飽和混合比の時間微分は 0 とする
	  \item $q_{v}^{t + \Delta t} > q_{vs}^{t + \Delta t}$
		ならば湿潤調節法を用いる
	 \end{itemize}

  \item (初期)値 $u^{\tau - n \Delta \tau}, v^{\tau - n \Delta \tau},
	w^{\tau - n \Delta \tau}, \pi^{\tau - n \Delta \tau}$ と
	$\theta^{t - \Delta t}, q_{v}^{t - \Delta t},
	 q_{c}^{t - \Delta t}, q_{r}^{t - \Delta t},
	 E^{t - \Delta t}$ を与える
	\label{導出始め}

  \item $\tau - n \Delta \tau (= t - \Delta t)$ から $\tau (= t)$
	 までの $u, v, w, \pi$ の値を $\Delta \tau$ 毎に求める
	 \label{速度圧力求める}
	 \begin{itemize}
	  \item $\theta, q_{v}, q_{c}, q_{r}, E$ は $t - \Delta t$
		の値を用いる
	 \end{itemize}
	 \begin{enumerate}
	  \item 水平方向の速度
		$u^{\tau - (n -1) \Delta \tau},
		 w^{\tau - (n -1) \Delta \tau}$
		を求める
		\label{水平速度求める}
	  \item 鉛直方向の運動方程式と圧力方程式より
		$\pi^{\tau - (n - 1) \Delta \tau}$ を消去して
		$w^{\tau - (n - 1) \Delta \tau}$ を求める
	  \item $\pi^{\tau - (n - 1) \Delta \tau}$ を求める
		\label{圧力求める}
	  \item \ref{水平速度求める}から\ref{圧力求める}の手順と
		同様にして,
		$\tau - (n - 2) \Delta \tau, \tau - (n - 3)
		 \Delta \tau, \cdots, \tau $ の $u, v, w.\pi$ を求める
	 \end{enumerate}

  \item \ref{速度圧力求める} の手順と同様にして,
	$\tau (= t)$ から $\tau + n \Delta \tau (= t + \Delta t)$
	 までの $u, v, w, \pi$ の値を $\Delta \tau$ 毎に求める
	 \begin{itemize}
	  \item $\theta, q_{v}, q_{c}, q_{r}, E$ は $t$
		の値を用いる
	 \end{itemize}

  \item $\theta^{t - \Delta t}, q_{v}^{t - \Delta t},
	 q_{c}^{t - \Delta t}, q_{r}^{t - \Delta t},
	 E^{t - \Delta t}$ と
	$\theta^{t}, q_{v}^{t}, q_{c}^{t}, q_{r}^{t},
	 E^{t}$ から
	 $\theta^{t + \Delta t}, q_{v}^{t + \Delta t},
	 q_{c}^{t + \Delta t}, q_{r}^{t + \Delta t},
	 E^{t + \Delta t}$ を求める
	 \begin{enumerate}
	  \item 熱力学方程式より $\theta$

	  \item 混合比の式より $q_{v}, q_{c}, q_{r}$
		(乾燥空気の場合 $q_{v} = q_{c} = q_{r} = 0$ )		
		\begin{enumerate}
	 	  \item $q_{v}^{t + \Delta t} > q_{vs}^{t + \Delta t}$
		ならば湿潤調節法を用いる
		\end{enumerate}

	  \item サブグリッドスケールの運動エネルギーの時間発展
		方程式より $E$
	 \end{enumerate}
	 \label{導出終わり}

  \item 得られた
	$u, v, w, \pi, \theta, q_{v}, q_{c}, q_{r}, E$
	を用いて, \ref{導出始め}から\ref{導出終わり}行ったと
	同様の手順でさらに次の値を求める
 \end{enumerate}

%======================================================================
\newpage
\appendix

\section{付録}

\subsection{湿潤飽和調節法}
\label{湿潤飽和調節法}

 飽和した状態の下での凝結-蒸発過程を数値的に評価するために,
 Soong and Ogura (1973) が提唱した湿潤飽和調節法を用いて
 熱力学方程式と混合比の方程式を解く.
 湿潤飽和調節法とは過飽和しないように温位と混合比を計算する
 方法である.

 どういう条件で湿潤飽和調節法を用いるのかが不明.
 調節法を用いるのは「雲ができたかどうか」で分けている.
 「雲ができた」と判別するには以下の 2 つの条件が考えられる:
 \begin{itemize}
  \item 水蒸気の混合比が飽和混合比より大きいかどうか
	(CReSS 方式):
	\begin{itemize}
	 \item 数値上過飽和になっている:
	       \begin{eqnarray}
	        q_{v}^{*} - q_{vs}^{*} > 0
	       \end{eqnarray}
	       $q_{c}^{t + \Delta t} > 0$ の場合
	       繰り返し計算.

	       $q_{c}^{t + \Delta t} < 0$ の場合不明.
	 \item 数値上過飽和になっていない:
	       \begin{eqnarray}
	        q_{v}^{*} - q_{vs}^{*} < 0
	       \end{eqnarray}
	       場合, $q_{c}^{t + \Delta t} = 0$ として計算.
	\end{itemize}
  \item 全凝結量の正負 (斉藤 (1999) 方式):
	\begin{eqnarray}
	 q_{c}^{*} + q_{v}^{*} - q_{vs}^{*} > 0
	\end{eqnarray}
	かつ $q_{c}^{t + \Delta t} > 0$ の場合繰り返し計算.
	
	かつ $q_{c}^{t + \Delta t} < 0$ の場合
	$q_{c}^{t + \Delta t} = 0$ として計算.

	$q_{c}^{*} + q_{v}^{*} - q_{vs}^{*} < 0$ の場合不明.
 \end{itemize}
 条件が微妙に違う気がする.
 分けた後の計算方法はどちらも同じようである.


 温位, 水蒸気と雲水の混合比の方程式は
 \begin{eqnarray}
  \DD{\phi}{t} &=& M_{\phi} + D_{\phi}
  \Deqlab{予報方程式}
 \end{eqnarray}
 と書ける.
 ただし $\phi = \theta, q_{v}, q_{c}$ であり, $D_{\phi}$ は
 拡散項, 微物理過程の生成・消滅項は
 \begin{eqnarray}
  M_{\theta} &=& - \gamma
                   \left\{\DD{q_{vs}}{t} + E_{r}\right\}
  , \\
  M_{q_{v}} &=& \DD{q_{vs}}{t} + E_{r}
  , \\
  M_{q_{c}} &=& - \DD{q_{vs}}{t} - A_{r} - C_{r}
  , \\
  q_{vs} &=& \frac{3.8}{\bar{p}}
           \exp
	   \left(17.27
	    \frac{\bar{\Pi} \theta -273}
	    {\bar{\Pi} \theta -36}\right)
 \end{eqnarray}
 である.
 $q_{vs}$ は飽和混合比, $E_{r}$ は蒸発による
 雨水から水蒸気への変換, $A_{r}, C_{r}$ は雲水から
 雨水への変換を表す.

 $(温位の式) + \gamma \times (水蒸気の混合比の式)$ より
 \begin{eqnarray}
  \DD{\theta}{t} + \gamma \DD{q_{v}}{t}
      = D_{\theta} + D_{q_{v}}
  \Deqlab{温位水蒸気}
 ,
 \end{eqnarray}
 水蒸気の混合比の式 + 雲水の混合比の式より
 \begin{eqnarray}
  \DD{q_{v}}{t} + \DD{q_{c}}{t}
      = D_{q_{v}} + D_{q_{c}}
        + E_{r} - A_{r} - C_{r}
  \Deqlab{水蒸気雲水}
 \end{eqnarray}
 である.

 予報方程式\Deqref{予報方程式}から計算された $t + \Delta t$ で
 の各値を仮に $\theta^{*}, q_{v}^{*}, q_{c}^{*}$ と表す.
 ここでは過飽和しないという条件で考えているので,
 $q_{v}^{*} > q_{vs}^{*}$ ならば
 湿潤飽和調節を行って $t + \Delta t$ での各値を求める.
 $q_{v}^{*} \le q_{vs}^{*}$ ならば予報方程式から得られた値を
 用いる.
 以下で湿潤飽和調節法を解説する.

 式 \Deqref{温位水蒸気} から
 \begin{eqnarray}
  \frac{\theta^{t + \Delta t} - \theta^{t - \Delta t}}{2 \Delta t}
     + \gamma
       \frac{q_{v}^{t + \Delta t} - q_{v}^{t - \Delta t}}{2 \Delta t}
     &=& D_{\theta}^{t - \Delta t} + D_{q_{v}}^{t - \Delta t}
  , \\
  \frac{\theta^{*} - \theta^{t - \Delta t}}{2 \Delta t}
     + \gamma
       \frac{q_{v}^{*} - q_{v}^{t - \Delta t}}{2 \Delta t}
     &=& D_{\theta}^{t - \Delta t} + D_{q_{v}}^{t - \Delta t}
  .
 \end{eqnarray}
 辺々引くと
 \begin{eqnarray}
  \theta^{t + \Delta t} + \gamma q_{v}^{t + \Delta t}
     = \theta^{*} + \gamma q_{v}^{*}
  \Deqlab{温位水蒸気2}
 \end{eqnarray}
 が成り立つ.

 過飽和しないので, 飽和に達した後の $t + \Delta t$ における
 水蒸気の混合比と飽和混合比は等しくなければならない:
 \begin{eqnarray}
  q_{v}^{t + \Delta t}
      &=& q_{vs}^{t + \Delta t}
  \nonumber \\
      &=& b \exp \left[\frac{a \left(\bar{\Pi}
				     \theta^{t + \Delta t} - 273
			       \right)}
                            {\bar{\Pi} \theta^{t + \Delta t} - 36}
		 \right]
  \nonumber \\
     &\approx& q_{vs}^{*}
               \left\{1 + \frac{273 a \bar{\Pi}}
		               {\left(\bar{\Pi}
				      \theta^{*} - 36
				 \right)^{2}}
		      \left(\theta^{t + \Delta t} - \theta^{*}\right)
	       \right\}
  .
 \end{eqnarray}
 ただし $a = 7.5 \ln 10, b = 3.8 / \bar{p}$ である.
 したがって, 式\Deqref{温位水蒸気2}は
 \begin{eqnarray}
  \theta^{t + \Delta t}
     &=& \theta^{*}
         + \gamma \left(q_{v}^{*} - q_{vs}^{t + \Delta t}\right)
  \nonumber \\
     &=& \theta^{*}
         + \gamma r_{1} \left(q_{v}^{*} - q_{vs}^{*}\right)
  \Deqlab{調節温位}
 \end{eqnarray}
 となる.
 ただし
 \begin{eqnarray}
    r_{1} &=& \left[1 + \gamma q_{vs}^{*}
	              \frac{273 a \bar{\Pi}}
		           {\left(\bar{\Pi} \theta^{*} - 36\right)^{2}}
	    \right]^{-1}
  .
 \end{eqnarray}

 式\Deqref{温位水蒸気2} と \Deqref{調節温位}から
 \begin{eqnarray}
  q_{v}^{t + \Delta t}
      = q_{v}^{*}
        + \frac{1}{\gamma}
	  \left(\theta^{*} - \theta^{t + \Delta t}\right)
  \Deqlab{調節水蒸気}
 \end{eqnarray}

 また, \Deqref{温位水蒸気2}を求めたのと同様に,
 \Deqref{水蒸気雲水}より
 \begin{eqnarray}
  q_{v}^{t + \Delta t} + q_{c}^{t + \Delta t}
   = q_{v}^{*} + q_{c}^{*}
  \Deqlab{調節雲水}
 \end{eqnarray}
 が得られる.

 式\Deqref{調節温位}, \Deqref{調節水蒸気}, \Deqref{調節雲水}より
 暫定的に $\theta, q_{v}, q_{c}$ の値を求める.
 求めた値について $q_{c}^{t + \Delta t} > 0$ ならば
 \footnote{雲水ができる場合.},
 式\Deqref{調節温位}, \Deqref{調節水蒸気}, \Deqref{調節雲水}を
 繰り返し用いて値が収束するまで計算する
 \footnote{凝結して周りの温度が上がるため, 蒸発が起こるのを繰り返す.
           この状態が落ち着くまで計算する.}.

 $q_{c}^{t + \Delta t} \le 0$ ならば
 \footnote{雲水ができない場合.},
 \begin{eqnarray}
  q_{c}^{t + \Delta t} &=& 0
  , \\
  q_{v}^{t + \Delta t} &=& q_{v}^{*} + q_{c}^{*}
  \Deqlab{調節水蒸気2}
  , \\
  \theta^{t + \Delta t} &=& \theta^{*} - \gamma q_{c}^{*}
  \Deqlab{調節雲水2}
 \end{eqnarray}
 を用いる.
 ここで, 式\Deqref{調節水蒸気2}は式\Deqref{調節雲水}より,
 式\Deqref{調節雲水2}は\Deqref{調節水蒸気}と
 \Deqref{調節水蒸気2}より得られる.


%--------------------------------------------------------------------
\pagebreak

\subsection{Asselin の時間フィルタ}
\label{asselin}

 長い時間間隔における Leap-frog 法による積分を計算する場合,
 微小な物理モードの擾乱で計算モードが強く減衰する
 (Asselin, 1972).
 この計算モードを抑えるため Asselin の時間フィルターを
 用いる.
 $\partial \phi / \partial t = F$ という方程式を解く場合,
 \begin{eqnarray}
  \phi^{* t + \Delta t} &=& \phi^{t - \Delta t}
                          + 2 \Delta t F^{* t}
  , \\
  \phi^{t} &=& \phi^{* t}
               + \mu_{a}
	         \left(\phi^{* t + \Delta t} - 2 \phi^{* t}
		       + \phi^{t - \Delta t}
		 \right)
 \end{eqnarray}
 ただし, アスタリスク (*) 付変数はフィルタをかけられていない
 量であり, $\phi = \theta, q_{v}, q_{c}, q_{r}$ である.
 また $\mu_{a}$ は 0.1 とした.
 このフィルターは各時間ステップ毎にかけられる

%======================================================================
%\newpage

%--------------------------------------------------------------------
%\pagebreak

%\subsection{計算例}


%======================================================================
%\newpage

%\section{}

%--------------------------------------------------------------------
%\pagebreak

%\subsection{}


%======================================================================
\newpage

\section{参考文献}

 \begin{description}
  \item  Asselin, R., 1972:
         Frequency filter for time integrations.
         {\it Mon. Weather Rev.}, {\bf 100}, 487-490. 

  \item Arakawa, A. and Lamb, V. R.,1977:
	Computational Design of the Basic Dynamical Processes
	of the UCLA General Circulation Model.
        Methods in Computational Physics,
        Academic Press, {\bf 17}, 174-267.

  \item Klemp J. B. and Robert B. Wilhelmson, 1978:
	The simulation of three-dimensional convective storm
	dynamics.
	{\it J. Atmos. Sci.}, {\bf 35}, 1070-1096.

  \item Klemp J. B., 1978:
	A  splitting procedure for numerical solution of
	the compressible equations of motion.
	???, ???.

  \item 中野満寿男, 2003:
        鉛直シアーが台風初期渦形成に及ぼす影響.
	九州大学大学院理学府 地球惑星科学専攻 修士論文.

  \item Ogura, Y., and N. A. Phillips, 1962:
        A scale analysis of deep and shallow convection in the
        atmosphere.
        {\it J. Atmos. Sci.}, {\bf 19}, 173-179.

  \item	Ogura, Y., and T. Takahashi, 1971:
	Numerical simulation of the life cycle of a thunderstorm
	cell.
	{\it Mon. Wea. Rev.}, {\bf 99}, 895-911.

  \item	Ogura, Y. and M. Yoshizaki, 1988:
	Numerical study of orographic-convective precipitation
	over the eastern Arabic Sea and the Ghats Mountains
	during the summer monsoon.
	{\it J. Atoms. Sci.}, {\bf 45}, 2097-2122.

  \item Wilhelmson, R., and Y. Ogura, 1972:
        The pressure perturbation and the numerical modeling
        of a cloud.
        {\it J. Atmos. Sci.}, {\bf 29}, 1295-1307.

  \item 大宮司久明, 三宅裕, 吉澤徴, 1998:
	乱流の数値流体力学---モデルと計算法.
	東京大学出版会, 652.

  \item	斉藤和雄 編, 1999: 非静力学モデル.
	気象研究ノート 第 196 号, 日本気象学会.

  \item Skamarock, Wiliam C. and Joseph B. Klemp, 1992:
        The stability of time-split methods for the hydrostatic
        and the nonhydrostatic elastic equation.
        {\it Mon. Wea. Rev.}, {\bf 120}, 2109-2127.

  \item Soong, S-T., and Y. Ogura, 1973:
	A comparison between axi-symmetric and slab-symmetric
	cumulus cloud models.
	{\it J. Atmos. Sci.}, {\bf 30}, 879-893

  \item Tokioka, T., 1978:
        Some considerations on vertical differencing.
	{\it J. Meteor. soc. Japan}, {\bf 56}, 98-111.

  \item 坪木和久, 榊原篤志, 2001: CReSS ユーザーズガイド 第 2 版. \\
         http://www.tokyo.rist.or.jp/CReSS\_Fujin/CReSS.top.html
 \end{description}

\vspace{2em}

%======================================================================
\newpage

\section*{謝辞}
\markright{謝辞}

 本資源は, 
 地球流体電脳倶楽部の
 インターネット上での学術知識の集積と活用の実験の一環として
 \begin{center}
  http://www.gfd-dennou.org/arch/deepconv/kaminari/
 \end{center}
 において公開されているものである.
 \copyright 高橋 こう子, 杉山 耕一朗, 小高 正嗣, 中島 健介, 林 祥介
  (K. Takahashi, K. Sugiyama, M. Odaka, K. Nakajima and Y.-Y. Hayashi)
 2003.
 本資源は, 
 著作者の諸権利に抵触しない(迷惑をかけない)限りにおいて
 自由に利用していただいて構わない. 
 なお, 利用する際には今一度自ら内容を確かめることをお願いする
 (無保証無責任原則). 

 本資源に含まれる元資源提供者(図等の版元等を含む)からは, 
 直接的な形での WEB 上での著作権または使用許諾を得ていない
 場合があるが, 
 勝手ながら, 「未来の教育」のための実験という
 学術目的であることをご理解いただけるものと信じ, 
 学術標準の引用手順を守ることで諸手続きを略させていただいている. 
 本資源の利用者には, この点を理解の上, 
 注意して扱っていただけるようお願いする. 
 万一, 不都合のある場合には
 \begin{center}
  deepconv@gfd-dennou.org
 \end{center}
 まで連絡していただければ幸いである. 

%======================================================================
\end{document}
