3. 力学過程

3.1 力学過程の概要

本章では力学過程のコードの構造と計算手順を説明する. ここでいう力学過程と はスペクトルで計算される部分のことを指し, 移流項や水平拡散項, 圧力勾配項 などである. あるいは物理過程以外と言っても良いかもしれない. スペクトル とは球面調和関数展開の係数のことであり, スペクトルへの変換とは球関数展開 に他ならない.

実際の時間積分は力学過程の中で行われる. 物理過程で計算された時間変化 項は力学過程に引き渡され, 力学過程の時間変化項に加えられ, セミインプ リシット法による時間積分に組み込まれる. 但し, 降水過程などの調節過程は時 間積分後に適用される.

3.2 力学過程 Dynamics の流れ

Fig.3.1〜 Fig.3.3 に力学過程サブルーチン Dynamics の流れを示す. 四角枠内の太文字はサブルーチン名を表している. Dynamics とサブルーチン間にやりとりされる変数が示されている. 四角枠内が空白であるものは, 簡略に記述できるため Dynamics サブルーチン本体に直接コードが記述されているものである. これらのサブルーチンは, モジュール dynamics_hspl_vas83 内に ひとまとめにされている. DynamicsInitおよびSemiImplMatrix 内で設定されたパラメータや配列はモジュール内の共有変数として保存され, 各サブルーチンはその変数を参照する仕組みとなっているため, 原則的に Dynamics と各サブルーチンとでやり取りされるのは 時間変化する変数のみとなる.

Dynamics で呼ばれている NonLinearOnGrid は時間変化率の非重力波(NG)項を計算する. TimeIntegration はセミインプリシット時間積分 (3.3.1節参照)を行う. 以上の2つが Dynamics の主な下請ルーチンである.

NonLinearOnGrid 以前に呼ばれるサブルーチンは 必要な定数等の準備にあてられる. SemiImplMatrix までは原則的に時間に依存しない. TimeIntegration の前後で格子点値とスペクトル値の変換を行う. その後, 水平粘性に失われた運動エネルギーを補償する 摩擦熱を行い, Dynamics 外に出ない量の出力(DiagOutput)を行う. DiagOutput では診断量の計算を行うが, これは モニターのためであり, 時間積分には影響しない.

図 3.1: 力学過程 Dynamics の流れ1. 四角枠内の太文字はサブルーチン名であり, Dynamics が呼ぶサブルーチンと引数を示している. 引数は, 特に注記されていないものは格子点値である. 四角枠内が空白であるものは, 簡略に記述できるため Dynamics サブルーチン本体に直接コードが記述されているものである. ここでは Dynamics が引数として値を受け取り, 初期設定 (DynamicsInit, SemiImplMatrix) および 非線形項の計算 (NonLinearOnGrid) を行っている ところまでが記されている. これ以後の処理については Fig.3.2を参照されたい.

図 3.2: 力学過程 Dynamics の流れ2. Fig.3.1の続きである. 図の見方についてはFig.3.1を参照せよ. ここでは格子点データをスペクトルデータに変換し, それらを用いて, セミインプリシット時間積分 TimeIntegration を行っている. これ以後の処理については Fig.3.3を参照されたい.

図 3.3: Fig.3.2の続きである. 図の見方についてはFig.3.1を参照せよ. ここではスペクトルデータを格子点データに戻し, 摩擦熱補正をかけた後, ヒストリデータを出力し, 最後に $t+\Delta t$の値をサブルーチンの外に返している.

3.3 セミインプリシット時間積分

3.3.1 セミインプリシット時間積分の概要

セミインプリシット法については, 支配方程式系とその離散化ドキュメント の 3.5 節に解説があるので 参照のこと. 時刻 $t-\Delta t$ から $t+\Delta t$ へのセミインプリシットの 計算手順は以下のようにまとめられる.

1. 発散項 $\Dvect{D}$ のみに関する単一のセミインプリシット方程式,

   $$\underline{M}\, \overline{\Dvect{D}}^{t} = \Dvect{f}$$ (3.1)

   
   
   を解く. ここで $$\overline{\Dvect{D}}^{t} \equiv \frac{1}{2} (\Dvect{D}^{t+\Delta t} + \Dvect{D}^{t-\Delta t})$$ である.

2. $\overline{\Dvect{D}}^{t}$ と時刻 $t$ における時間変化の非重力波 (NG)項から各変数の $t+\Delta t$ での値を求める.

以上の手順の説明は 3.3.2 節で行う.

3.3.2 セミインプリシット時間積分の詳細

3.3.1 節で述べたように, セミインプリシット時間積分 は2つのステップに分けられる. これを詳しく書くと以下のようになる.

● 第1段階: $\overline{\Dvect{D}}^{t}$ の計算

$$\underline{M}\, \overline{\Dvect{D}}^{t} = \Dvect{f}$$ (3.2)

という形の式 を解くのが第1段階 である. ここでは支配方程式系とその離散化ドキュメント の 3.5 節と同様, 太字で鉛直に離散化したベクトルを表 し, 下線で行列を表す. 水平方向に関しては, 各スペクトル係数が独立に計算され るので離散化したことが表向きは現れない. $\Dvect{D}$ は発散で,

$$\overline{\Dvect{D}}^{t} \equiv \frac{1}{2} (\Dvect{D}^{t+\Delta t} + \Dvect{D}^{t-\Delta t})$$ (3.3)

である. 行列 $\underline{M}$ は

$$\underline{M} \equiv ( 1-2\Delta t {\cal D}_H )( 1-2\Delta t {\ca... ...ne{h} + ( 1-2\Delta t {\cal D}_H )\Dvect{G} \Dvect{C}^{T} ] \nabla^{2}_{\sigma}$$ (3.4)

であり, ベクトル $\Dvect{f}$ は,

$$\Dvect{f} \equiv ( 1-2\Delta t {\cal D}_H )( 1-\Delta t {\cal D}_M ) \Dvect{D}^{t-... ... t} + ( 1-2\Delta t {\cal D}_H ) \Delta t \left( \DP{\Dvect{D}}{t} \right)_{NG}$$

$$- \Delta t \nabla^{2}_{\sigma} \left\{ ( 1-2\Delta t {\cal D}_H )... ...}^{t-\Delta t} + \Delta t \left( \DP{\Dvect{T}}{t} \right)_{NG} \right] \right.$$

$$\left. \hspace*{20mm} + ( 1-2\Delta t {\cal D}_H ) \Dvect{G} \left[ \pi^{t-\Delta t} + \Delta t \left( \DP{\pi}{t} \right)_{NG} \right] \right\} .$$ (3.5)

である.

$\underline{M}$ はサブルーチン SemiImplMatrix で計算され, LU分解される. $\underline{M}$ は時間刻みの大きさが変らない限り設定し直す必要がない. $\Dvect{f}$ の計算と (3.2) を解く作業 はサブルーチン TimeIntegration で行われる.

● 第2段階: 時間積分
時間積分も TimeIntegration が行う. 第1段階で $\overline{\Dvect{D}}^{t}$ が 求まったので発散の時間積分は容易である:

$$\Dvect{D}^{t+\Delta t} = 2 \overline{\Dvect{D}}^{t} - \Dvect{D}^{t-\Delta t} .$$ (3.6)

他の物理量についても, $\overline{\Dvect{D}}^{t}$ と時刻 $t$ における時間 変化率の非重力波(NG)項から, $t+\Delta t$ での値が求まる:

$${\cal X}^{t+\Delta t} = \gamma_{\cal X}^{-1} \left\{ {\cal X}^{t-... ...ight)_{NG} + {\cal G}_{{\cal X}} \overline{\Dvect{D}}^{t} \right] \right\} \, .$$ (3.7)

ここで,

$$\begin{array}{cccccc} \displaystyle {\cal X} & = & \displayst... ...tyle -\underline{h}, & \displaystyle \underline{0}. \end{array}$$ (3.8)

以下の各節で, 具体的な計算手順とプログラムソースとの対応を述べる.

3.4 モデル変数表

ここでは, これ以後の節で述べられるモデル変数の表記法を述べる. モデル変数は    $$ のようなフォントで示され, その入出力に関連して表 3.1 のように分類される.

表 3.1: 数式と対応させて書く場合の モデル変数の表記法. 本文中で定数の下に付した文字列は, 値を設定したサブルーチン名を示す. 緑で示される定数は, 物理定数管理, 座標データ管理モジュールなどから参照するか, もしくはそれらの値を用いて DynamicsInit や SemiImplMatrix で時間ループの初回に設定される.

当該サブルーチンからの出力変数 (入力に上書きする場合を含む) の出力時	$^{\mbox{$\bullet$}}$
関数からの返り値	〃
当該サブルーチン内の作業変数または, 入(出)力変数が入(出)力時の値と違う 場合	$$
当該サブルーチンへの入力変数	$$
関数への入力変数	〃
当該サブルーチンへの入力に上書きして出力される変数の入力時	$^{\mbox{$\bullet$}}$
定数 (座標データや物理定数など)	$$
関数	$$

3.5 非重力波項の計算 (サブルーチン NonLinearOnGrid)

3.5.1 NonLinearOnGrid前の準備

NonLinearOnGridの直前で計算される, 地表面気圧の空間変化と渦度発散について記す. 灰色文字で示されるのは spml ライブラリの関数である.

$$\underbrace{\left(\nabla \pi (t)\right)_{\lambda}}_{ \mbox{{\cmssbx xy\_GradLonPiN}}^{\mbox{$\bullet$}} } = \underbrace{\Dinv{a}}_{ \mbox{{\cmssbx\textcolor{PineGreen}{RPl... ...\mbox{{\cmssbx\textcolor{Gray}{\qquad\qquad xy\_GradLon\_w( w\_xy( ... ) )}}} }$$ (3.9)

$$\underbrace{\left(\nabla \pi (t)\right)_{\varphi}}_{ \mbox{{\cmssbx xy\_GradLatPiN}}^{\mbox{$\bullet$}} } = \underbrace{\Dinv{a}}_{ \mbox{{\cmssbx\textcolor{PineGreen}{RPl... ...{{\cmssbx\textcolor{Gray}{\qquad\qquad\quad xy\_GradLat\_w( w\_xy( ... ) )}}} }$$ (3.10)

$$\underbrace{\zeta (t)}_{ \mbox{{\cmssbx xyz\_VorN}}^{\mbox{$\bullet$}} } = \underbrace{\Dinv{a}}_{ \mbox{{\cmssbx\textcolor{PineGreen}{RPl... ...\} }_{ \mbox{{\cmssbx\textcolor{Gray}{xya\_wa( wa\_Div\_xya\_xya( ... ) )}}} },$$ (3.11)

$$\underbrace{D (t)}_{ \mbox{{\cmssbx xyz\_DivN}}^{\mbox{$\bullet$}} } = \underbrace{\Dinv{a}}_{ \mbox{{\cmssbx\textcolor{PineGreen}{RPl... ...\} }_{ \mbox{{\cmssbx\textcolor{Gray}{xya\_wa( wa\_Div\_xya\_xya( ... ) )}}} }.$$ (3.12)

3.5.2 NonLinearOnGridでの計算

非重力波(NG)項はサブルーチン NonLinearOnGrid で計算される.

NG 項の計算に関する式の各項の下にモデル中の変数名を記す. 以下で, 数式と対応 させて書く場合のモデル変数の表記法は Table 3.1 に従って統一する.

$$\left( \DP{\pi}{t} \right)^{\rm NG} = \underbrace{ Z_{k}, }_{ \mbox{{\cmssbx xy\_DPiDtN}}^{\mbox{$\bullet$}} }$$ (3.13)

$$\begin{align*}\begin{split}\left( \DP{D_{k}}{t} \right)^{\rm NG} &= \Dinv{a \cos... ..._{ \mbox{{\cmssbx xyz\_KinEngyN}}^{\mbox{$\bullet$}} } , \end{split}\end{align*}$$ (3.14)

$$\begin{align*}\begin{split}\left( \DP{\zeta_{k}}{t} \right)^{\rm NG} &= \Dinv{a ... ...!\!\!\!\!\!\!\!\!\!\!} \cos \varphi)}{\varphi} \biggl\}, \end{split}\end{align*}$$ (3.15)

$$\begin{align*}\begin{split}\left( \DP{T_{k}}{t} \right)^{\rm NG} &= \Dinv{a \cos... ...{\mbox{$\bullet$}} \!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!} , \end{split}\end{align*}$$ (3.16)

$$\begin{align*}\begin{split}\left( \DP{q_{k}}{t} \right)^{\rm NG} &= \Dinv{a \cos... ...\!\! \mbox{{\cmssbx xyz\_DQVapDtN}}^{\mbox{$\bullet$}} } \end{split}\end{align*}$$ (3.17)

NonLinearOnGrid の出力変数は, 以上の10個の時間変化項と, 以下の変数である.

$$^{\mbox{$\bullet$}} \quad \equiv \quad \dot{\sigma}.$$ (3.18)

NonLinearOnGrid の入力変数は時刻 $t$ でのグリッド値    $$,    $$ ($\ast$ = 変数名) である.

次に時間変化項の具体的な中身とモデル変数の対応を見ることにする.

地表面気圧の変化に関連する項の計算

$$\underbrace{ Z_{k} }_{ \mbox{{\cmssbx xy\_DPiDtN}}^{\mbox{$\bullet$}} } = - \underbrace{ \sum_{k=1}^{K} \overbrace{\Dvect{v}_{k} \cdot \n... ...tcolor{PineGreen}{z\_DelSigma}}} } }_{ \mbox{{\cmssbx xyz\_PiAdvSum$_{k=1}$}} }$$ (3.19)

ここで,

$$\underbrace{\Dvect{v}_{k} \cdot \nabla \pi}_{ \mbox{{\cmssbx xyz\_PiAdv}} } = \underbrace{u_k}_{ \mbox{{\cmssbx\textcolor{blue}{xyz\_UN}}} } ... ...\!\!\!\!\!\!\!\!\!\!\!\!\!\! \mbox{{\cmssbx\textcolor{blue}{xy\_GradLatPiN}}} }$$ (3.20)

渦度, 発散の変化に関連する項の計算

$$\begin{align*}\begin{split}\overbrace{{u_A}_{,k=1}}^{\!\!\!\!\!\!\!\! \mbox{{\cm... ...\cos \varphi} \DP{\pi}{\lambda} + {\cal F}_{\lambda, k}, \end{split}\end{align*}$$ (3.21)

$$\begin{align*}\begin{split}\overbrace{{v_A}_{,k=1}}^{\!\!\!\!\!\!\!\! \mbox{{\cm... ... T_{v,k}'}{a} \DP{\pi}{\varphi} + {\cal F}_{\varphi, k}, \end{split}\end{align*}$$ (3.22)

$$\underbrace{ (\mbox{\sl KE})_{k} + \sum_{l=1}^{K}W_{kl} (T_{v,l}-... ...}_{\!\!\!\!\!\! \mbox{\footnotesize {\tt HydroGrid}を呼ぶと $W_{kl}$がかかる} }$$ (3.23)

ここで,

$$\begin{align*}\begin{split}\underbrace{\dot{\sigma}_{k-1/2}}_{ \mbox{{\cmssbx xy... ...\ \dot{\sigma}_{1/2} &= 0, \\ \dot{\sigma}_{K+1/2} &= 0, \end{split}\end{align*}$$ (3.24)

$$\underbrace{T_{v,k}^{\prime}}_{ \mbox{{\cmssbx xyz\_TempVirEdd}} } = \overbrace{T_{v,k}}^{ \mbox{{\cmss... ...{T}_k}_{\!\!\!\!\!\!\!\!\! \mbox{{\cmssbx\textcolor{PineGreen}{z\_RefTemp}}} },$$ (3.25)

$$\underbrace{T_{v,k}}_{ \mbox{{\cmssbx xyz\_TempVir}} } = \underbrace{T_k}_{ \mbox{{\cmssbx\textcolor{blue}{xyz\_TempN}}}... ...\mbox{{\cmssbx\textcolor{blue}{\!\!\!\!\!\!\!\!\!\!\!\!xyz\_QVapN}}} } \bigg\},$$ (3.26)

$$\begin{align*}\begin{split}\underbrace{\hat{\kappa}_k}_{\!\!\!\!\!\!\!\!\!\!\!\!... ...\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!} } { \Delta \sigma_k }, \end{split}\end{align*}$$ (3.27)

$$\begin{align*}\begin{split}\underbrace{\alpha_k}_{\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!... ...lor{PineGreen}{Kappa}}} \!\!\!\!\!\!\!\!\!\!\!\!\!}} -1, \end{split}\end{align*}$$ (3.28)

$$\begin{align*}\begin{split}\underbrace{\beta_k}_{\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\... ...( \frac{ \sigma_{k+1/2} } { \sigma_k } \biggl)^{\kappa}, \end{split}\end{align*}$$ (3.29)

$$\begin{align*}\begin{split}\underbrace{\kappa}_{\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!... ...{PineGreen}{CpDry}}} \!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!}}. \end{split}\end{align*}$$ (3.30)

温度の変化に関連する項の計算

$$^{\mbox{$\bullet$}} \equiv \underbrace{u_k}_{ \mbox{{\cmssbx\textcolor{blue}{xyz\_UN}... ...mssbx xyz\_TempEdd}} \ \mbox{\scriptsize\Deqref{code-nongrav:温度の擾乱}参照} }$$ (3.31)

$$^{\mbox{$\bullet$}} \equiv \underbrace{v_k}_{ \mbox{{\cmssbx\textcolor{blue}{xyz\_VN}}} } \!\!\!\! T'_k,$$ (3.32)

$$\begin{align*}\begin{split}\underbrace{ H_k }_{ \mbox{{\cmssbx xyz\_DTempDt}}^{\... ...lta \sigma_{K} + T'_{v,K} D_K \Delta \sigma_{K} \right]. \end{split}\end{align*}$$ (3.33)

ここで,

$$\underbrace{T_{k}^{\prime}}_{ \mbox{{\cmssbx xyz\_TempEdd}} } = \overbrace{T_{k}}^{ \mbox{{\cmssbx\textcolor{blue}{xyz\_TempN}}... ...{T}_k}_{\!\!\!\!\!\!\!\!\! \mbox{{\cmssbx\textcolor{PineGreen}{z\_RefTemp}}} },$$ (3.34)

$$\underbrace{\hat{T}_{k-1/2}^{\prime}}_{ \mbox{{\cmssbx xyr\_TempE... ...} & \text{($k = 2, \cdots, K$)} \\ 0 , & \text{($k = K+1$)} \end{array} \right.$$ (3.35)

$$\underbrace{\hat{\overline{T}}_{k-1/2}}_{ \begin{array}{r} \mbox{... ..., & \text{($k = 2, \cdots, K$)} \\ 0 , & \text{($k = K+1$)} \end{array} \right.$$ (3.36)

$$\underbrace{a_k}_{\begin{array}{r} \mbox{{\cmssbx\textcolor{PineGreen}{z\_TInpCoefA}}} \\ \mbox{{\footnotesize\tt DynamicsInit}} \end{array}} = \overbrace{\alpha_k}^{ \mbox{{\cmssbx\textcolor{PineGreen}{z\_H... ...eft[ 1- \left( \frac{ \sigma_k }{ \sigma_{k-1} } \right)^{\kappa} \right]^{-1},$$ (3.37)

$$\underbrace{b_k}_{\begin{array}{r} \mbox{{\cmssbx\textcolor{PineGreen}{z\_TInpCoefB}}} \\ \mbox{{\footnotesize\tt DynamicsInit}} \end{array}} = \overbrace{\beta_k}^{ \mbox{{\cmssbx\textcolor{PineGreen}{z\_Hy... ...ft[ \left( \frac{ \sigma_k }{ \sigma_{k+1} } \right)^{\kappa} - 1 \right]^{-1},$$ (3.38)

$$\underbrace{ \dot{\sigma}^{\rm NG}_{k-1/2}}_{ \mbox{{\cmssbx xyr\... ...{l} \cdot \nabla \pi \Delta \sigma_{l}}_{ \mbox{{\cmssbx xyz\_PiAdvSum}} } \, .$$ (3.39)

比湿の変化に関連する項の計算

$$^{\mbox{$\bullet$}} \equiv \underbrace{u_k}_{ \mbox{{\cmssbx\textcolor{blue}{xyz\_UN}... ...!\!\!\!\!\! \overbrace{q_{k},}^{ \mbox{{\cmssbx\textcolor{blue}{xyz\_QVapN}}} }$$ (3.40)

$$^{\mbox{$\bullet$}} \equiv \underbrace{v_k}_{ \mbox{{\cmssbx\textcolor{blue}{xyz\_VN}}} } \!\!\!\! q_k,$$ (3.41)

$$\begin{align*}\begin{split}\overbrace{R_1}^{ \mbox{{\cmssbx xyz\_DQVapDtN}}^{\mb... ...Delta \sigma_K} \dot{\sigma}_{K-1/2} ( q_{K-1} - q_K ) . \end{split}\end{align*}$$ (3.42)

3.6 格子点データからスペクトルデータへ

サブルーチン NonLinearOnGrid から得られた 非重力波項    $,$   $$ や, 物理過程の時間変化項    $,$   $$, $t-\Delta t$の予報変数    $,$   $$を スペクトルデータ    $^{\mbox{$\bullet$}} , \mbox{{\cmssbx wz\_D*DtN}}^{\mbox{$\bullet$}} , \mbox{{\cmssbx w\_*B}}^{\mbox{$\bullet$}} , \mbox{{\cmssbx wz\_*B}}^{\mbox{$\bullet$}}$ へ変換する.

3.7 セミインプリシット時間積分
(サブルーチン TimeIntegration)

3.3.2節 述べた時間積分の2つの手続きは TimeIntegration で行われる. ここでは TimeIntegration 内の計算の流れを, 物理量とモデ ル変数の対応を示しながら説明する.

3.7.1 (3.2) の $\Dvect{f}$ の計算

時間積分の第1段階は (3.2) を解くことである. 時間に依存しない行列 $\underline{M}$ は SemiImplMatrix で計算され, LU分解されて    $$として TimeIntegration で使用されるので, 最初にするべきことは $\Dvect{f}$ を求めることである. (3.5) を再掲し, 各 項の下にモデル内の変数名を記す.

$$\begin{align*}\begin{split}\underbrace{\Dvect{f}}_{ \mbox{{\cmssbx wz\_siDivAvrT... ...r.2ex\hbox{]}} }_{ \mbox{{\cmssbx w\_siPi}} } \Biggr\} . \end{split}\end{align*}$$ (3.43)

3.7.2 (3.2) を解く

サブルーチン LUSolve により解かれる. 結果は    $$ に上 書きされる. その他の LUSolve の引数は $\underline{M}$ とそのピボット (    $$ と    $$) である.

3.7.3 $t+\Delta t$の値の算出

(3.6) と (3.7) を解く. 以下に式とモデル変数との対応を記す.

$$\underbrace{\pi^{t+\Delta t}}_{ \mbox{{\cmssbx w\_PiA}}^{\mbox{$\bullet$}} } = \underbrace{\pi^{t-\Delta t}}_{ \mbox{{\cmssbx\textcolor{blue}{... ...{ \mbox{{\cmssbx wz\_siDivAvrTime}} } \Bigg] }_{ \mbox{{\cmssbx w\_siDPiDt}} },$$ (3.44)

$$\underbrace{\Dvect{\zeta}^{t+\Delta t}}_{ \mbox{{\cmssbx wz\_VorA}}^{\mbox{$\bullet$}} } = \Dinv{1-2\Delta t {\cal D}_M} \Biggl\{ \underbrace{\Dvect{\zeta... ...}}{t}^{\rm NG}}_{ \mbox{{\cmssbx\textcolor{blue}{wz\_DVorDtN}}} } \Biggr\} \, ,$$ (3.45)

$$\underbrace{\Dvect{D}^{t+\Delta t}}_{ \mbox{{\cmssbx wz\_DivA}}^{\mbox{$\bullet$}} } = 2 \overline{\Dvect{D}}^{t} - \underbrace{\Dvect{D}^{t-\Delta t}}_{ \mbox{{\cmssbx\textcolor{blue}{wz\_DivB}}} },$$ (3.46)

$$\underbrace{\Dvect{T}^{t+\Delta t}}_{ \mbox{{\cmssbx wz\_TempA}}^{\mbox{$\bullet$}} } = \Dinv{1-2\Delta t {\cal D}_H} \Biggl\{ \underbrace{\Dvect{T}^{t... ...rline{\Dvect{D}}^{t} \Biggr] }_{ \mbox{{\cmssbx wz\_siDTempDt}} } \Biggr\} \, ,$$ (3.47)

$$\underbrace{\Dvect{q}^{t+\Delta t}}_{ \mbox{{\cmssbx wz\_QVapA}}^{\mbox{$\bullet$}} } = \Dinv{1-2\Delta t {\cal D}_H} \Biggl\{ \underbrace{\Dvect{q}^{t... ...\! \mbox{{\cmssbx\textcolor{blue}{wz\_DQVapDtN}}} \!\!\!\!\!\!\!\!\!} \Biggr\}.$$ (3.48)