\section{氷相を含む Le Treut and Li (1991) の非対流性凝結過程 (大規模凝結)}
\label{sec:ll91ice}

Le Treut and Li (1991) の非対流性凝結 (大規模凝結) 過程に
氷相を考慮して定式化する. なお, 氷相がない場合においても, 
扱う保存量等において Le Treut and Li (1991) とは異なることに注意
\footnote{何が異なるのか忘れてしまった (yot, 2016/05/11)}. 


\subsection{数理表現}

\subsubsection{支配方程式}

凝結に関わる物理量は下の方程式を満たす. 
%
\begin{align}
  \DD{}{t} \int \rho \left( C_p T + L q_v - L_i q_i \right) dV &= 0
\\
  \DD{}{t} \int \rho q_t dV &= 0
\\
  q_t = q_v + q_c &= q_v + q_l + q_i
\end{align}

しかし, 上記の式のみでは未知数の数に対して方程式の数が足らない.
そこで, 下のように仮定する.
これは, 格子内で, 全水量が $\Delta q$ の大きさの揺らぎを一様分布で
持つことを意味する. 
%
\begin{align}
  q_l &= \alpha\left(T\right) q_c
\\
  q_i &= \left\{ 1 - \alpha\left(T\right) \right\} q_c
       = \beta\left(T\right) q_c
\\
  q_c &= 
    \begin{cases}
      0
        & (q_v^* > q_t + \Delta q)
    \\
      \frac{ \left( q_t + \Delta q - q_v^* \right)^2 }{4 \Delta q}
        & (q_t - \Delta q \le q_v^* \le q_t + \Delta q)
    \\
      q_t - q_v^*
        & (q_v^* < q_t - \Delta q)
    \end{cases}
\\
  \Delta q &= \gamma q_t
\\
  \gamma &= 
    \begin{cases}
      0
        & (q_v^* > q_t + \Delta q)
    \\
      \frac{ q_t + \Delta q - q_v^* }{2 \Delta q}
        & (q_t - \Delta q \le q_v^* \le q_t + \Delta q)
    \\
      1
        & (q_v^* < q_t - \Delta q)
    \end{cases}
\end{align}


\subsection{離散表現}

\subsubsection{支配方程式}

前節の支配方程式は下のように離散化される. 
%
\begin{align}
  C_p \hat{T} + L \hat{q}_v - L_i \hat{q}_i &= C_p T + L q_v - L_i q_i 
\\
  \hat{q}_t &= q_t
\\
  q_t &= q_v + q_c = q_v + q_l + q_i
%\\
%  q_v &= q_t - q_c
%\\
%  q_c &= q_l + q_i
\\
  q_l &= \alpha\left(T\right) q_c
\\
  q_i &= \left\{ 1 - \alpha\left(T\right) \right\} q_c
       = \beta\left(T\right) q_c
\\
  q_v^* &= q_v^*(T, p)
\\
  q_c &= 
    \begin{cases}
      0
        & (q_v^* > q_t + \Delta q)
    \\
      \frac{ \left( q_t + \Delta q - q_v^* \right)^2 }{4 \Delta q}
        & (q_t - \Delta q \le q_v^* \le q_t + \Delta q)
    \\
      q_t - q_v^*
        & (q_v^* < q_t - \Delta q)
    \end{cases}
\\
  \Delta q &= \gamma q_t
\\
  \gamma &= 
    \begin{cases}
      0
        & (q_v^* > q_t + \Delta q)
    \\
      \frac{ q_t + \Delta q - q_v^* }{2 \Delta q}
        & (q_t - \Delta q \le q_v^* \le q_t + \Delta q)
    \\
      1
        & (q_v^* < q_t - \Delta q)
    \end{cases}
\end{align}
%
ここで, $\hat{}$ は調節前の値を表す.


\subsubsection{数値解法}

上式を解けばよいが, 非線形であるため解析的には解けない. 
そこで, 繰り返し法によって解を求める. 

$q_v^*$ を, テイラー展開して一次の項までとると,
%
\begin{align}
  q_v^*(T, p)
%
    &=   q_v^*(\hat{T}, p)
      + \left. \frac{\partial q_v^*}{\partial T}\right|_{T=\hat{T}} \Delta T
  \Deqlab{ll91ice:math:lscond_qsat_taylor}
\\
  \Delta T &= T - \hat{T}
\end{align}
%
となることを用い, また, 非線形になることを防ぐため, 
$\alpha = \alpha\left(\hat{T}\right)$, 
$\beta = \beta\left(\hat{T}\right)$ として整理すると, 
%
\begin{align}
  \Delta T &=
    \frac{ 
           - L \hat{q}_l
           - \left( L + L_i \right) \hat{q}_i
           + \left\{ L + L_i \beta\left(\hat{T}\right) \right\} A 
         }
         {
         C_p - \left\{ L + L_i \beta\left(\hat{T}\right) \right\} B
         }
\\
  q_v &= q_t - q_c
\\
  q_l &= \alpha\left( \hat{T} \right) q_c
\\
  q_i &= \beta\left( \hat{T} \right) q_c
\\
  q_c &= A + B \Delta T
\\
  A &=
    \begin{cases}
      0
        & (q_v^* > q_t + \Delta q)
    \\
      \frac{ \left( \hat{q}_t + \Delta q - q_v^*(\hat{T}) \right)^2 }{4 \Delta q}
        & (q_t - \Delta q \le q_v^* \le q_t + \Delta q)
    \\
      \hat{q}_t - q_v^*(\hat{T})
        & (q_v^* < q_t - \Delta q)
    \end{cases}
\\
  B &=
    \begin{cases}
      0
        & (q_v^* > q_t + \Delta q)
    \\
      - \frac{ \hat{q}_t + \Delta q - q_v^*(\hat{T}) }{2 \Delta q}
        \left. \frac{\partial q_v^*}{\partial T}\right|_{T=\hat{T}} \Delta T
        & (q_t - \Delta q \le q_v^* \le q_t + \Delta q)
    \\
      - \left. \frac{\partial q_v^*}{\partial T}\right|_{T=\hat{T}} \Delta T
        & (q_v^* < q_t - \Delta q)
    \end{cases}
\end{align}
%
となる.

ただし, \Deqref{ll91ice:math:lscond_qsat_taylor} で $q_{v,k}$ をテイラー
展開で近似しているため, 上記の結果は近似値である.
したがって, 上記の計算を繰り返し行う. 繰り返しの回数を $l$ とすると, 値を
%
\begin{eqnarray}
  (T_{k})_{l+1} &=& (T_{k})_l + \Delta T_{k}
\\
  (q_{v,k})_{l+1} &=& (q_{v,k})_l + \Delta q_{v,k}
\\
  (q_{l,k})_{l+1} &=& (q_{l,k})_l + \Delta q_{l,k}
\\
  (q_{i,k})_{l+1} &=& (q_{i,k})_l + \Delta q_{i,k}
\end{eqnarray}
%
のように更新することで解を求める.