%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % %表題 2 次元非静力学モデル -- 数理モデル % %履歴 2003-07-25 高橋 こう子 新規作成 % 2003-09-03 高橋 こう子 骨子完成 % 2003-09-08 高橋 こう子 参考文献追加 % 2003-09-09 高橋 こう子 Appendix 追加 % 2003-11-17 高橋 こう子 修正 % 2003-11-22 高橋 こう子 運動エネルギーの式修正 % 2004-07-25 小高正嗣 2 次元乾燥大気版へ再構成 % 2004-11-30 小高正嗣 % 2005-01-31 杉山耕一朗 % 2005-08-25 小高正嗣 熱力学の式に散逸加熱項を追加 % 2005-08-26 小高正嗣 熱力学の式に基本場の拡散項を追加 % 2006-03-10 杉山耕一朗 % 2006-03-12 杉山耕一朗 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %\chapter{座標系} %空間座標系は 2 次元直線直交座標系を用いる. 水平方向の座標変数を $x$, %鉛直方向の座標変数を $z$ と表す. 時間方向の変数は $t$ とする. \chapter{基礎方程式系} 本数値モデルは水平・鉛直の 2 次元モデルである. 水平方向の座標変数を $x$, 鉛直方向の座標変数を $z$ と表し, 時間方向の変数は $t$ と表す. %\section{予報変数の定義} % % モデルの独立変数は空間方向の座標変数 $x, z$ と時間方向の変数 $t$ である. % モデルの予報変数はこれらの関数として定義される % \begin{center} % \begin{tabular}[tb]{lll} % $u$ & : & 速度の $x$ 成分 \\ % $w$ & : & 速度の $z$ 成分 \\ % $\theta$ & : & 温位 \\ % $\pi$ & : & 無次元圧力関数 % \end{tabular} % \end{center} % である. % 無次元圧力関数 (エクスナー関数) $\pi$ は圧力 $p$ を用いて以下のように定義され % る. % \begin{equation} % \Pi \equiv \left(\frac{p}{p_{0}}\right)^{R_{d}/c_{p}} % \Deqlab{エクスナー関数} % \end{equation} % $p_{0}$ は地表面での気圧である. %温位 $\Theta$ は % \begin{eqnarray} % \Theta &\equiv& T \left(\frac{p_{0}}{p}\right)^{R_{d}/c_{p}} % \nonumber \\ % &=& \frac{T}{\Pi} % \Deqlab{温位} % \end{eqnarray} % で与えられる. % ここで $T$ は温度, $c_{p}$ は % 定圧比熱, $R_{d}$ は乾燥空気の気体定数である. %-------------------------------------------------------------------- \section{運動方程式・圧力方程式・熱の式・混合比の保存式} 力学的な枠組みは, 準圧縮方程式系(Klemp and Wilhelmson,1978)を用いる. この方程式系では, 予報変数を水平一様な基本場とそこからのずれに分離し, 方程式の線形化を行っている. 方程式中の変数は Appendix \ref{appendix:hensu} に示す. 以下に準圧縮方程式系の時間発展方程式を一覧する. 密度の式では乾燥成分と湿潤成分の分子量の差を考慮するが, 熱の式では考慮しない. また圧力方程式では非断熱加熱による大気の膨張と, 凝縮に伴う圧力変化を無視している. % \begin{description} \item[運動方程式] % \begin{eqnarray} \DP{u}{t} &=& - \left( u \DP{u}{x} + w \DP{u}{z} \right) - {c_{p}}_{d} \bar{\theta_{v}} \DP{\pi}{x} + Turb.u \\ %% \DP{w}{t} &=& - \left( u \DP{u}{x} + w \DP{u}{z} \right) - {c_{p}}_{d} \bar{\theta_{v}} \DP{\pi}{z} + Turb.w \nonumber \\ && + \left( \frac{\theta}{\bar{\theta}} + \frac{\sum q_{v}/M_{v}}{1/M_{d} + \sum \bar{q_{v}}/M_{v}} - \frac{\sum q_{v} + \sum q_{c} + \sum q_{r}} {1 + \sum \bar{q_{v}}} \right) g \end{eqnarray} \item[圧力方程式] % \begin{eqnarray} % \DP{\pi}{t} = - \frac{\overline{{C_{s}}^{2}}}{{c_{p}}_{d} \bar{\rho} \bar{\theta_{v}}^{2}} \DP{}{x_{j}}(\bar{\rho} \bar{\theta_{v}} u_{j}) \end{eqnarray} \item[熱の式] % \begin{eqnarray} \DP{\theta}{t} &=& - \left( u \DP{\theta}{x} + w \DP{\theta}{z} \right) - w\DP{\bar{\theta}}{x} + \Dinv{\bar{\pi}} \left(Q_{cnd} + Q_{rad} + Q_{dis}\right) % + \frac{L}{{c_{p}}_{d} \bar{\pi}} \left( CN_{vc} - EV_{cv} - EV_{rv} \right) \nonumber \\ && + Turb.\bar{\theta} + Turb.\theta \Deqlab{def_theta} \end{eqnarray} \item[混合比の保存式] \begin{eqnarray} \DP{q_{v}}{t} &=& - \left( u \DP{q_{v}}{x} + w \DP{q_{v}}{z} \right) - w\DP{\bar{q_{v}}}{x} % - \left( CN_{vc} - EV_{cv} - EV_{rv} \right) % + Turb.q_{v} + Turb.\bar{q_{v}}, + Src.q_{v} + Turb.q_{v} + Turb.\bar{q_{v}}, \Deqlab{def_qv} \\ % \DP{q_{c}}{t} &=& - \left( u \DP{q_{c}}{x} + w \DP{q_{c}}{z} \right) % + ( CN_{vc} - EV_{cv} - CN_{cr} - CL_{cr}) + Turb.q_{c}, + Src.q_{c} + Turb.q_{c} \Deqlab{def_qc} \\ % \DP{q_{r}}{t} &=& - \left( u \DP{q_{c}}{x} + w \DP{q_{c}}{z} \right) % + (CN_{cr} + CL_{cr} - EV_{rv} ) + PR_{r} + Turb.q_{r} + Src.q_{r} + Fall.q_{r} + Turb.q_{r} \Deqlab{def_qr} % \end{eqnarray} % \end{description} ただし, $\bar{~}$ の付いた変数は水平一様な基本場であることを示し, 上付き添え字 $c$ は個々の凝縮成分を示す. % \begin{description} \item[エクスナー関数 $\pi$] % \begin{equation} \pi \equiv \left(\frac{p}{p_{0}}\right)^{R_{d}/{c_{p}}_{d}} \Deqlab{エクスナー関数の定義} \end{equation} % \item[温位 $\theta$] % \begin{eqnarray} \theta \equiv T \left(\frac{p_{0}}{p}\right)^{R_{d}/{c_{p}}_{d}} = \frac{T}{\pi} \Deqlab{温位の定義} \end{eqnarray} % \item[密度 $\rho$] % \begin{eqnarray} \rho &=& \frac{p}{R_{d}T} \left( \frac{1/{M_{d}}} {1/M_{d} + \sum{{q_{v}}/{M_{v}}}} \right) \left( 1 + \sum q_{v} + \sum q_{c} + \sum q_{r} \right) \nonumber \\ &=& \frac{p}{R_{d}T_{v} } = \frac{p_{0} \pi^{{c_{v}}_{d}/R_{d}}}{R_{d} \theta_{v}} \Deqlab{密度の定義} \end{eqnarray} \item[仮温位 $\theta_{v}$] % \begin{eqnarray} \theta_{v} = \frac{\theta}{ \left( \frac{1/M_{d}} {1/M_{d} + \sum{{q_{v}}/{M_{v}}}} \right) \left( 1 + \sum q_{v} + \sum q_{c} + \sum q_{r} \right) } \end{eqnarray} \item[音波速度 ${C_{s}}^{2}$] % \begin{eqnarray} && {C_{s}}^{2} = \frac{{c_{p}}_{d}}{{c_{v}}_{d}} R_{d} \pi \theta_{v} \end{eqnarray} % \end{description} %-------------------------------------------------------------------- %\pagebreak \section{雲微物理過程のパラメタリゼーション} 方程式系に含まれる凝縮による加熱項 $Q_{cnd}$, 生成項 $Src$, 落下項 $Fall$ の評価は, 中島(1998)で用いられた Kessler (1969) のパラメタリゼーションに従う. 暖かい雨のバルク法のパラメタリゼーションでは, 気相と凝縮相を以下の 3 つのカテゴリーに分ける. \\ % \begin{tabular}{lll}\\ \hline 記号 & 意味 & 内容 \\ \hline $q_{v}$ & 気相の混合比 & 気体の状態で大気中に存在する水 \\ $q_{c}$ & 雲水混合比 & 落下速度がゼロな液体の粒子で, 実際の大気中の 雲粒に対応する.\\ & & 通常 100 $\mu$m 以下の微小な流体粒子である. \\ $q_{r}$ & 雨水混合比 & 有意な落下速度を持つ液体の粒子で, 実際の 大気中の雨粒に対応する. \\ \hline \end{tabular} そして, 微物理素過程として以下を考慮する. ただし, これらの量は全て正の値として定義され, 水蒸気が直接雨水に凝結する過程は無視されている. \\ % \begin{tabular}{ll}\\ \hline 記号 & 内容 \\ \hline $CN_{vc}$ & 凝結による水蒸気から雲水への変換 (condensation) \\ $EV_{cv}$ & 蒸発による雲水から水蒸気への変換 (evaporation) \\ $EV_{rv}$ & 蒸発による雨水から水蒸気への変換 (evaporation) \\ $CN_{cr}$ & 併合成長による雲水から雨水への変換. \\ & 併合や水蒸気拡散により, 雲粒子が雨粒の大きさにまで成長する (autocondensation) \\ $CL_{cr}$ & 衝突併合による雲水から雨水への変換. 大水滴が小水滴を衝突併 合する (collection) \\ $PR_{r}$ & 雨水の重力落下に伴う雨水混合比の変化率 (Precipitation)\\ \hline \end{tabular} この微物理素過程を用いて \Deqref{def_qv} -- \Deqref{def_qr} 式を書き直すと, 以下のようになる. % \begin{eqnarray} \DP{\theta}{t} &=& - \left( u \DP{\theta}{x} + w \DP{\theta}{z} \right) - w\DP{\bar{\theta}}{x} + \frac{L}{{c_{p}}_{d} \bar{\pi}} \left( CN_{vc} - EV_{cv} - EV_{rv} \right) \nonumber \\ && + \Dinv{\bar{\pi}} \left(Q_{rad} + Q_{dis}\right) + Turb.\bar{\theta} + Turb.\theta \\ \DP{q_{v}}{t} &=& - \left( u \DP{q_{v}}{x} + w \DP{q_{v}}{z} \right) - w\DP{\bar{q_{v}}}{x} - \left( CN_{vc} - EV_{cv} - EV_{rv} \right) \nonumber \\ && + Turb.q_{v} + Turb.\bar{q_{v}}, \Deqlab{def_qv2} \\ % \DP{q_{c}}{t} &=& - \left( u \DP{q_{c}}{x} + w \DP{q_{c}}{z} \right) + ( CN_{vc} - EV_{cv} - CN_{cr} - CL_{cr}) + Turb.q_{c}, \Deqlab{def_qc2} \\ % \DP{q_{r}}{t} &=& - \left( u \DP{q_{c}}{x} + w \DP{q_{c}}{z} \right) + (CN_{cr} + CL_{cr} - EV_{rv} ) + PR_{r} + Turb.q_{r} \Deqlab{def_qr2} % \end{eqnarray} % ここで, $\gamma = L_{v}/ ({c_{p}}_{d} \pi)$ であり, $L_{v}$ は水の蒸発の潜熱[J K$^{-1}$ kg$^{-1}$], ${c_{p}}_{d}$ は乾燥大気の定圧比熱[J K kg$^{-1}$], $\pi$ はエクスナー関数である. 微物理素過程は以下のように定式化する. % %CReSS(坪木と榊原, 2001)と Klemp and Wilhelmson (1978) を基に %雲の微物理素過程を定式化する. CReSSユーザーマニュアル(坪木と榊原, 2001) %と Klemp and Wilhelmson (1978) は同じ定式化を用いているが, %飽和蒸気圧 $q_{vsw}$ と雨水の終端速度 $U_{r}$ の式の係数がそれぞれ異なる. %その理由は両者の用いる単位系が互いに異なるためで, %Klemp and Wilhelmson (1978) では圧力として mbar, %速度として cm s$^{-1}$ を用いている. \begin{description} \item[水蒸気と雲水の間の変換: $-CN_{vc} + EV_{cv}$]~\\ % 雲水は粒が小さく, 水蒸気との間で瞬間的に飽和調節が起こるもの とする. すなわち, 移流などの項を計算した後の温度と水蒸気量が 過飽和状態となっている場合には, ちょうど飽和になる量の水蒸気 を凝縮させる. 一方, 移流などの項を計算した後に, 雲水が存在す るにも拘わらず未飽和になっている場所では, ちょうど飽和になる 量の雲水を蒸発させる. % \begin{eqnarray} -CN_{vc} + EV_{cv} = \DD{q_{vsw}}{t}. \Deqlab{teishiki:def_CNvc+EVcv} \end{eqnarray} % \item[飽和混合比: $q_{vsw}$] ~\\ % % 飽和蒸気圧 $p_{vsw}$ は Antoine の式を用いて以下のように表現 % する. % % \begin{eqnarray} % p_{vsw} % = A - \frac{B}{C + T} % = A - \frac{B}{C + (\theta \pi)} % \end{eqnarray} %% % 水の場合, Antoine 式中の係数の値は, %% % \begin{eqnarray} % A = 23.1964 \hspace{1em} %% B = 3816.44 \hspace{1em} % C = -46.13 % \end{eqnarray} %% % である. 値は化学工学便覧を参照した. 飽和蒸気圧より飽和混合比を % 求めると以下となる. %% % \begin{eqnarray} % q_{vsw} % = \varepsilon \frac{p_{vsw}}{p - p_{vsw}} % = \varepsilon \frac{p_{vsw}}{p_{0} \pi^{{c_{p}}_{v}/R_{d}} - p_{vsw}} % \end{eqnarray} %% % ここで $\varepsilon$ は水蒸気と乾燥空気の分子量比であり, % 地球の場合には $\varepsilon = 0.622$ である. % この定式化では凝縮成分の体積は気体成分の体積に比べ無視できるくらい % 小さいと仮定し, 凝縮成分混合比の圧力に与える影響を無視している. \item[雲水の併合成長: $CN_{cr}$] ~\\ % Berry(1968)に従って以下で与える. % \begin{eqnarray} && CN_{cr} = \frac{10^{6} \bar{\rho} q_{c}^{3}} {60 (2 q_{c} + 2.66 \times 10^{-8} N_{0}/ \bar{\rho} D_{0})} \Deqlab{def_CNcr} \\ && \hspace{2em} N_{0} = 5.0 \times 10^{7}, \hspace{1em} D_{0} = 0.366 \nonumber \end{eqnarray} % これを選択した理由は, 第一には Kessler (1969) のもともとの式 が閾値を含んでいるために生じる不自然さ(長時間積分で, 領域の 上部に雲水が溜まり, ついには全域が雲に覆われてしまう)を嫌っ たためである. 第二には, Berry(1969) の式がまがりなりにも, 海洋上の雲物理を扱えると主張しているからである. \item[雲水の衝突併合: $CL_{cr}$] ~\\ % Kessler (1969) に従って, 以下で定式化する. % \begin{eqnarray} && CL_{cr} = 2.2 q_{c} (\bar{\rho} q_{r})^{0.875} . \Deqlab{def_CLcr} \end{eqnarray} % \item[雨水の蒸発: $EV_{rv}$]~\\ % \begin{eqnarray} EV_{rv} = 4.85 \times 10^{-2} (q_{vsw} - q_{v}) (\bar{\rho} q_{r})^{0.65} \end{eqnarray} \item[雨水のフラックス: $PR_{r}$ ]~\\ % 雨水の重力落下による混合比の変化率は, % \begin{eqnarray} PR_{r} = \Dinv{\bar{\rho}} \DP{}{z}(\bar{\rho} U_{r} q_{r}). \Deqlab{teishiki:def_PRr} \end{eqnarray} % であり, 雨水の終端落下速度 $U_{r}$ [m s$^{-1}$] は % \begin{eqnarray} U_{r} = 12.2 (q_{r})^{0.125} \end{eqnarray} % で与える. \end{description} %-------------------------------------------------------------------- %\pagebreak \section{簡略化した放射冷却} 本モデルにおいて, 放射強制は高度のみに依存するパラメタとして与える. 温位の式の右辺に, 放射強制 $Q_{rad}$ を考慮する. %-------------------------------------------------------------------- %\pagebreak \section{乱流混合のパラメタリゼーション} \subsection{乱流運動エネルギーの式} Klemp and Wilhelmson (1978) および CReSS (坪木と榊原, 2001) と同様に, 1.5 次のクロージャーを用いることで, 乱流エネルギーの時間発展方程式は以 下ように書ける. % \begin{eqnarray} \DP{K_{m}}{t} &=& - \left( u \DP{K_{m}}{x} + w \DP{K_{m}}{z} \right) - \frac{3 g C_{m}^{2} l^{2}}{ 2 \overline{\theta}} \left(\DP{\theta_{el}}{z} \right) \nonumber \\ && + \left( C_{m}^{2} l^{2} \right) \left\{ \left( \DP{u}{x} \right)^{2} + \left( \DP{w}{z} \right)^{2} \right\} \nonumber \\ && + \frac{ C_{m}^{2} l^{2} }{2} \left( \DP{u}{z} + \DP{w}{x}\right)^{2} - \frac{K_{m}}{3} \left( \DP{u}{x} + \DP{w}{z} \right) \nonumber \\ && + \Dinv{2} \left(\DP[2]{K_{m}^{2}}{x} + \DP[2]{K_{m}^{2}}{z} \right) + \left(\DP{K_{m}}{x}\right)^{2} + \left(\DP{K_{m}}{z}\right)^{2} \nonumber \\ && - \Dinv{2 l^{2}} K_{m}^{2} \Deqlab{kiso:TurbE} \end{eqnarray} % ここで $C_{\varepsilon} = C_{m} = 0.2$, 混合距離 $l = \left(\Delta x \Delta z \right)^{1/2}$ とする. $\theta_{el}$ は以下のように定義する \footnote{これで正しいのか要チェック. 中島さんの書き方とは整合的だが... 何を見るべきだか....}. % \begin{eqnarray} && \theta_{el} = \theta + \frac{L q_{v}}{{c_p}_{d} \bar{\pi}} + \bar{\theta} \left\{ \frac{\sum q_{v}/M_{v}}{1/M_{d} + \sum \bar{q_{v}}/M_{v}} - \frac{\sum q_{v} + \sum q_{c} + \sum q_{r}} {1 + \sum \bar{q_{v}}} \right\} \end{eqnarray} \subsection{運動方程式中の拡散項} Klemp and Wilhelmson (1978) および CReSS (坪木と榊原, 2001) と同様に, 1.5 次のクロージャーを用いることで粘性拡散項は以下のように書ける. % \begin{eqnarray} Turb.{u_{i}} &=& - \DP{}{x_{j}} \overline{(u_{i}^{\prime} u_{j}^{\prime})} \nonumber \\ &=& - \DP{}{x_{j}} \left[ - K_{m} \left(\DP{u_{i}}{x_{j}} + \DP{u_{j}}{x_{i}}\right) + \frac{2}{3} \delta_{ij} E \right]. \end{eqnarray} ここで $K_{m}$ は運動量に対する乱流拡散係数であり, $E$ は サブグリッドスケールの乱流運動エネルギー \begin{eqnarray} E = \frac{1}{2} \overline{(u^{\prime})^{2} + (w^{\prime})^{2}} = \frac{{K_{m}}^{2}}{C_{m} l^{2}} \end{eqnarray} である. \subsection{熱力学の式の拡散項} Klemp and Wilhelmson (1978) および CReSS (坪木と榊原, 2001) と同様に, 1.5 次のクロージャーを用いることで温位の粘性拡散項は以下のように書ける. % \begin{eqnarray} Turb.{\theta} &=& - \DP{}{x_{j}} \overline{u_{i}^{\prime} \theta} \nonumber \\ &=& - \DP{}{x_{j}} \left(K_{h}\DP{\theta}{x_{j}}\right) . \end{eqnarray} ここで $K_{h}$ は温位に対する乱流拡散係数である. \subsection{散逸加熱項の表現} 散逸加熱項 $Q_{dis}$ は, 乱流運動エ ネルギーの散逸項をもとに, 以下のように与える. \begin{equation} Q_{dis} = \frac{1}{\overline{c_{p}}}\frac{C_{\varepsilon}}{l} \frac{K_{m}^{3}}{(C_{m}l)^{3}}. \end{equation} ここで $l=(\Delta x\Delta z)^{1/2}$ である. %--------------------------------------------------------------------