Up 数値予報モデルの計算 作成: 2022-08-08
更新: 2022-08-12

気象庁の数値予報モデル

    田中 (2017), pp.201-203 より：
    気圧座標系 ( p系) で書いたプリミティブ方程式系
      \[ \begin{eqnarray} &状態方程式 \quad \quad & &： \quad& p \alpha = R_d T &\quad \quad (1)& \\ &静力学の式 \quad \quad & &： \quad& \frac{\partial \phi}{\partial p} = -\alpha &\quad \quad (2)& \\ &連続の式　 \quad \quad & &： \quad& \nabla \cdot {\bf v} + \frac{\partial \omega}{\partial p} = 0 &\quad \quad (3)& \\ &運動方程式 \quad \quad & &： \quad& \frac{d {\bf v}}{d t} = - \nabla \phi - f k \times {\bf v} + F &\quad \quad (4)& \\ &熱力学の式 \quad \quad & &： \quad& \frac{d c_p T}{d t} = \omega \alpha + Q &\quad \quad (5)& \end{eqnarray} \]
    を，実際に時間積分する方法：


    変数は， \( {\bf v} = (u, v),\ \omega,\ \alpha,\ T,\ \phi \) の６つで， 独立な方程式が６本あるので方程式系は閉じている。
    ここで， ラグランジュ微分は
      \[ \frac{d}{d t} = \frac{\partial}{\partial t} + {\bf v} \cdot \nabla + \omega \frac{\partial}{\partial p} \]
    である。

    運動方程式と熱力学の式が予報方程式で， ほかは診断方程式である。
    時間積分によって得られるコアとなる気象要素は， 水平風速 \( v \) と気温 \( T \) であり， これらが各格子点上で初期値として与えられているとする.
    ただし， p系の場合は，海面更正気圧 \( p_s \) (sea level pressure) の値も予報変数として与えられる必要がある.

    はじめに， 微分を含まない状態方程式 (1) を用いて，気温Tから比容 \( \alpha \) を計算する。
    気圧容 \( p \) は座標軸なので，これは簡単な計算で完了する。

    比容 \( \alpha \) が得られたなら， 静力学の式 (2) を用いて， 比容 \( \alpha \) からジオポテンシャル \( \phi \) を計算する。
    この式は鉛直1階微分なので， \( \phi \) の境界条件が1つ必要であるが，\( p = p_s \) の海面で \( \phi = 0 \) を用いればよい。
    微分を差分に近似して，\( \phi \) を下から上へ積み上げることで， すべての格子点上の \( \phi \) を計算する。

    次に連続の式 (3) を用いて， 水平風 \( {\bf v} \) から \( \omega \) を計算する.
    発散 \( \nabla \cdot v \) は周辺の格子点から中央差分で求め， 発散が求まったら， それを用いて \( \omega \) を計算する.
    このときの \( \omega \) の境界条件は， 気圧 \( p = 0 \) で \( \omega = 0 \) を用いればよい。
    微分を差分に近似して， \( \omega \) を上層から下層に向かつて計算する。
    鉛直積分の結果 \( p = p_s \) で \( \omega = \omega_s \) となる。

    以上の診断方程式を用いた計算により， 6つの状態変数 \( {\bf v} = (u, v),\ \omega,\ \alpha,\ T,\ \phi \) が，すべての格子点で得られた。

    次に， これらの気象要素を用いて粘性摩擦 \( F \) と非断熱加熱 \( Q \) を求める。

    このプロセスは物理過程のパラメタリゼーションとよばれ， 多くの手法が開発されている。
    力学系モデルの誤差は1%程度なので， 物理過程も同程度の誤差で定式化できればよいのであるが， 放射過程や降水過程などの物理過程の定式化はたいへん困難であり， 時には真の値との比較で， 誤差100%を超えることも多い。
    精度良くこれらを求めることが望まれるが， ここではその詳細は省略する。

    最も簡単なパラメタリゼーションとしては， レイリー摩擦 \( F = - \nu v \) とニュートン冷却 \( Q = - \epsilon (T -T_0) \) などがある。
    風があれは反対向きに摩擦力がはたらき，気温が平衡状態からずれたときには， 平衡状態に線形近似で緩和させるというものである.

    粘性摩擦 \( F \) と非断熱加熱 \( Q \) が計算できれば， 運動方程式 (4) と熱力学の式 (5) の右辺が，移流項も含めてすべて計算可能になる。
    緯度経度格子モデルの各点で空間微分があれば， 差分で近似計算を行う。
    その結果， (4), (5) の左辺の \( \frac{\partial {\bf v}}{\partial t},\ \frac{\partial T}{\partial t} \) の時間変化が既知となる。

    ただし，p系の場合， 海面更正気圧 \( p_s \) を以下の傾向方程式で， 上と同様に時間積分する必要がある。
      \[ \frac{\partial p_s}{\partial t} = \omega_s - {\bf v}_s \cdot \nabla p_s \]

    これらの時間変化項がわかれば， あとは時間積分である。
    力学系の数値積分で説明したように， 現在の値に変化量を可算することで，将来の値を計算する。
      \[ A( t + \Delta t ) = A( t ) + \frac{\partial A}{\partial t} \Delta t \]
    ここで， 状態変数ベクトルは \( A = (u,\ v,\ T,\ p_s)^T \) であり， すべての格子点の値である.

    ここでは， 最も簡単な前方差分による時間積分を例示したが， 中央差分を用いたリープブロッグによる方法や， 後方差分を前方差分の後で繰り返して行う松野スキームによる時間積分など， より安定した時間積分法が開発されている。

    これで， 時間軸上で \( \Delta t \) だけ将来の \( A = (u,\ v,\ T,\ p_s)^T \) が得られたので，アルゴリズムのループの最初に戻って同様の計算を繰り返すことで\将来の気象要素の値を計算することができる。

  • 引用/参考文献
    • 田中 博 (2017) :『地球大気の科学』, 共立出版, 2017.
        第12章 数値予報の原理, 12.3 緯度経度格子モデル