Class dynamics_hspl_vas83
In: dynamics/dynamics_hspl_vas83.F90

力学過程 (スペクトル法, Arakawa and Suarez (1983))

Dynamical Core (Spectral method, Arakawa and Suarez (1983))

Note that Japanese and English are described in parallel.

力学過程を演算するモジュールです. 水平離散化にスペクトル法 (Bourke, 1988) を, 鉛直離散化には Arakawa and Suarez (1983) を用いています.

時間積分法にはリープフロッグスキームを用いています. デフォルトでは, $ Delta t $ を大きくとるために, 重力波項に セミインプリシット法を適用しています. NAMELIST#dynamics_hspl_vas83_nml の TimeIntegScheme を変更することで, 重力波項をエクスプリシット法によって解くことも可能です.

This is a dynamical core module. Spectral method (Bourke, 1988) (for horizontal) and Arakawa and Suarez (1983) method (for vertical) are used.

Leap-frog scheme is used as time integration. By default, semi-implicit scheme is applied to gravitational terms in order to enlarge the value of $ Delta t $ . Explicit scheme can be applied to gravitational terms by changing "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml".

Procedures List

DynamicsHSplVAS83 :力学計算
DynamicsHSplVAS83Init :初期化
DynamicsHSplVAS83Finalize :終了処理 (モジュール内部の変数の割り付け解除)
——————— :————
DynamicsHSplVAS83 :Calculate dynamics
DynamicsHSplVAS83Init :Initialization
DynamicsHSplVAS83Finalize :Termination (deallocate variables in this module)

NAMELIST

NAMELIST#dynamics_hspl_vas83_nml

References

  • Bourke, W. P., 1988: Spectral methods in global climate and weather prediction models. Physically Based Modelling and Simulation of Climates and Climatic Change. Part I., M.E. Schlesinger (ed.), Kluwer Academic Publishers, Dordrecht, 169--220.
  • Arakawa, A., Suarez, M. J., 1983: Vertical differencing of the primitive equations in sigma coordinates. Mon. Wea. Rev., 111, 34--35.

Methods

Included Modules

dc_types dc_message constants gridset composition timeset axesset wa_mpi_module_sjpack wa_mpi_module wa_zonal_module wa_module_sjpack wa_zonal_module_sjpack wa_module gtool_historyauto dc_string intavr_operate mass_fixer sltt adv_test lumatrix namelist_util gtool_history dc_iounit dc_calendar dc_present dc_trace

Public Instance methods

Subroutine :
xyz_UB(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ u (t-\Delta t) $ . 東西風速. Eastward wind
xyz_VB(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ v (t-\Delta t) $ . 南北風速. Northward wind
xyz_TempB(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ T (t-\Delta t) $ . 温度. Temperature
xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: $ q (t-\Delta t) $ . 比湿. Specific humidity
xy_PsB(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ p_s (t-\Delta t) $ . 地表面気圧. Surface pressure
xyz_UN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ u (t) $ . 東西風速. Eastward wind
xyz_VN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ v (t) $ . 南北風速. Northward wind
xyz_TempN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ T (t) $ . 温度. Temperature
xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: $ q (t) $ . 比湿. Specific humidity
xy_PsN(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ p_s (t) $ . 地表面気圧. Surface pressure
xyz_DUDtPhy(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ left(DP{u}{t}right)^{phy} $ . 外力項 (物理過程) による東西風速変化. Eastward wind tendency by external force terms (physical processes)
xyz_DVDtPhy(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ left(DP{v}{t}right)^{phy} $ . 外力項 (物理過程) による南北風速変化. Northward wind tendency by external force terms (physical processes)
xyz_DTempDtPhy(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ left(DP{T}{t}right)^{phy} $ . 外力項 (物理過程) による温度変化. Temperature tendency by external force terms (physical processes)
xyzf_DQMixDtPhy(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: $ left(DP{q}{t}right)^{phy} $ . 外力項 (物理過程) による比湿変化. Temperature tendency by external force terms (physical processes)
xy_SurfHeight(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ z_s $ . 地表面高度. Surface height.
xyz_UA(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ u (t+Delta t) $ . 東西風速. Eastward wind
xyz_VA(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ v (t+Delta t) $ . 南北風速. Northward wind
xyz_TempA(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ T (t+Delta t) $ . 温度. Temperature
xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(out)
: $ q (t+Delta t) $ . 比湿. Specific humidity
xy_PsA(0:imax-1, 1:jmax) :real(DP), intent(out)
: $ p_s (t+Delta t) $ . 地表面気圧. Surface pressure
xyz_ArgOMG(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out), optional
: OMEGA, DP/Dt

力学過程の演算を行い, 与えられた $ t-\Delta t $ および $ t $ の 東西風速 (xyz_UB, xyz_UN), 南北風速 (xyz_VB, xyz_VN), 温度 (xyz_TempB, xyz_TempN), 質量混合比 (xyzf_QMixB, xyzf_QMixN), 地表面気圧 (xyz_PsB, xyz_PsN), および地表面高度 (xy_SurfHeight) から, $ t+Delta t $ の 東西風速 (xyz_UA), 南北風速 (xyz_VA), 温度 (xyz_TempA), 質量混合比 (xyzf_QMixA), 地表面気圧 (xyz_PsA) を返します.

別の物理プロセスによる渦度, 発散, 温度, 比湿の変化を, 力学過程の変化に足して次のステップを計算する場合には, それらの変化を xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyzf_DQMixDtPhy に与えてください.

時間積分法にはリープフロッグスキームを用いています. デフォルトでは, $ Delta t $ を大きくとるために, 重力波項に セミインプリシット法を適用しています. NAMELIST#dynamics_hspl_vas83_nml の TimeIntegScheme を変更することで, 重力波項をエクスプリシット法によって解くことも可能です.

Calculate dynamical processes. Input data are Eastward wind (xyz_UB, xyz_UN), northward wind (xyz_VB, xyz_VN), temperature (xyz_TempB, xyz_TempN), mass mixing ratio (xyzf_QMixB, xyzf_QMixN), surface pressure (xyz_PsB, xyz_PsN) at $ t-\Delta t $ and $ t $, and surface height (xy_SurfHeight). Eastward wind (xyz_UA), northward wind (xyz_VA), temperature (xyz_TempA), mass mixing ratio (xyzf_QMixA), surface pressure (xyz_PsA) are returned.

In order to add tendencies of vorticity, divergence, temperature and specific humidity calculated by physical processes other than those by this dynamical process, these tendencies should be given to "xyz_DUDtPhy", "xyz_DVDtPhy", "xyz_DTempDtPhy", "xyzf_DQMixDtPhy"

Leap-frog scheme is used as time integration. By default, semi-implicit scheme is applied to gravitational terms for extension of $ Delta t $ . Explicit scheme can be applied to gravitational terms by changing "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml".

[Source]

  subroutine DynamicsHSplVAS83( xyz_UB,      xyz_VB,      xyz_TempB,      xyzf_QMixB,   xy_PsB, xyz_UN,      xyz_VN,      xyz_TempN,      xyzf_QMixN,   xy_PsN, xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyzf_DQMixDtPhy, xy_SurfHeight, xyz_UA,      xyz_VA,      xyz_TempA,      xyzf_QMixA,   xy_PsA, xyz_ArgOMG )
    !
    ! 力学過程の演算を行い, 与えられた $ t-\Delta t $ および $ t $ の
    ! 東西風速 (xyz_UB, xyz_UN), 南北風速 (xyz_VB, xyz_VN), 
    ! 温度 (xyz_TempB, xyz_TempN), 質量混合比 (xyzf_QMixB, xyzf_QMixN), 
    ! 地表面気圧 (xyz_PsB, xyz_PsN), および地表面高度 (xy_SurfHeight) から, 
    ! $ t+\Delta t $ の
    ! 東西風速 (xyz_UA), 南北風速 (xyz_VA), 温度 (xyz_TempA), 
    ! 質量混合比 (xyzf_QMixA), 地表面気圧 (xyz_PsA) を返します. 
    !
    ! 別の物理プロセスによる渦度, 発散, 温度, 比湿の変化を, 
    ! 力学過程の変化に足して次のステップを計算する場合には, 
    ! それらの変化を xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyzf_DQMixDtPhy
    ! に与えてください. 
    !
    ! 時間積分法にはリープフロッグスキームを用いています.
    ! デフォルトでは, $ \Delta t $ を大きくとるために, 重力波項に 
    ! セミインプリシット法を適用しています.
    ! NAMELIST#dynamics_hspl_vas83_nml の TimeIntegScheme を変更することで, 
    ! 重力波項をエクスプリシット法によって解くことも可能です.
    !
    ! Calculate dynamical processes. 
    ! Input data are
    ! Eastward wind (xyz_UB, xyz_UN), 
    ! northward wind (xyz_VB, xyz_VN), 
    ! temperature (xyz_TempB, xyz_TempN), 
    ! mass mixing ratio (xyzf_QMixB, xyzf_QMixN), 
    ! surface pressure (xyz_PsB, xyz_PsN) at $ t-\Delta t $ and $ t $, 
    ! and surface height (xy_SurfHeight).
    ! Eastward wind (xyz_UA), northward wind (xyz_VA), 
    ! temperature (xyz_TempA), 
    ! mass mixing ratio (xyzf_QMixA), surface pressure (xyz_PsA) 
    ! are returned.
    !
    ! In order to add tendencies of vorticity, divergence, 
    ! temperature and specific humidity calculated by
    ! physical processes other than those by
    ! this dynamical process, 
    ! these tendencies should be given to 
    ! "xyz_DUDtPhy", "xyz_DVDtPhy", "xyz_DTempDtPhy", "xyzf_DQMixDtPhy"
    !
    ! Leap-frog scheme is used as time integration. 
    ! By default, semi-implicit scheme is applied to gravitational terms 
    ! for extension of $ \Delta t $ . 
    ! Explicit scheme can be applied to gravitational terms by changing
    ! "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml". 
    !

    ! モジュール引用 ; USE statements
    !

    use constants, only: RPlanet, CpDry, Grav
                              ! $ g $ [m s-2]. 
                              ! 重力加速度. 
                              ! Gravitational acceleration

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: lmax, imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

    ! 組成に関わる配列の設定
    ! Settings of array for atmospheric composition
    !
    use composition, only: ncmax, IndexH2OVap, CompositionInqFlagAdv

    ! 時刻管理
    ! Time control
    !
    use timeset, only: DelTime, TimeN, TimesetClockStart, TimesetClockStop

    ! 座標データ設定
    ! Axes data settings
    !
    use axesset, only: r_Sigma, z_Sigma               ! $ \sigma $ レベル (整数). 
                              ! Full $ \sigma $ level

    ! SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応) 
    ! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
    !
#ifdef LIB_MPI
#ifdef SJPACK
    use wa_mpi_module_sjpack, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya  => wa_DivLambda_xva, wa_DivMu_xya      => wa_DivMu_xva, xy_GradLambda_w   => xv_GradLambda_w, xya_GradLambda_wa => xva_GradLambda_wa, xy_GradMu_w       => xv_GradMu_w, xya_GradMu_wa     => xva_GradMu_wa, w_xy              => w_xv, xy_w              => xv_w, wa_xya            => wa_xva, xya_wa            => xva_wa
#else
    use wa_mpi_module, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya  => wa_DivLambda_xva, wa_DivMu_xya      => wa_DivMu_xva, xy_GradLambda_w   => xv_GradLambda_w, xya_GradLambda_wa => xva_GradLambda_wa, xy_GradMu_w       => xv_GradMu_w, xya_GradMu_wa     => xva_GradMu_wa, w_xy              => w_xv, xy_w              => xv_w, wa_xya            => wa_xva, xya_wa            => xva_wa
#endif
#elif AXISYMMETRY
    use wa_zonal_module, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya, wa_DivMu_xya, xy_GradLambda_w, xya_GradLambda_wa, xy_GradMu_w, xya_GradMu_wa, w_xy, xy_w, wa_xya, xya_wa
#elif SJPACK
    use wa_module_sjpack, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya, wa_DivMu_xya, xy_GradLambda_w, xya_GradLambda_wa, xy_GradMu_w, xya_GradMu_wa, w_xy, xy_w, wa_xya, xya_wa
#elif AXISYMMETRY_SJPACK
    use wa_zonal_module_sjpack, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya, wa_DivMu_xya, xy_GradLambda_w, xya_GradLambda_wa, xy_GradMu_w, xya_GradMu_wa, w_xy, xy_w, wa_xya, xya_wa
#else
    use wa_module, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya, wa_DivMu_xya, xy_GradLambda_w, xya_GradLambda_wa, xy_GradMu_w, xya_GradMu_wa, w_xy, xy_w, wa_xya, xya_wa
#endif

    ! ヒストリデータ出力
    ! History data output
    !
    use gtool_historyauto, only: HistoryAutoPut

    ! 文字列操作
    ! Character handling
    !
    use dc_string, only: LChar

    ! 種別型パラメタ
    ! Kind type parameter
    !
    use dc_types, only: DP      ! 倍精度実数型. Double precision. 

    ! 積分と平均の操作
    ! Operation for integral and average
    !
    use intavr_operate, only: IntLonLat_xy

    ! 質量の補正
    ! Mass fixer
    !
    use mass_fixer, only: MassFixer, MassFixerColumn

    ! セミラグランジュ法による物質移流計算
    ! Semi-Lagrangian method for tracer transport
    !
    use sltt, only: SLTTMain

!!$    ! 物質移流 (有限体積法)
!!$    ! Tracer Transport (finite volume method)
!!$    !
!!$    use tt_fv, only: TTFVMain

    !
    ! Utility module for advection test
    !
    use adv_test, only: AdvTestSetHorVels, AdvTestSetVerVel


    ! 宣言文 ; Declaration statements
    !
    implicit none

    real(DP), intent(in):: xyz_UB    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ u (t-\Delta t) $ .   東西風速. Eastward wind
    real(DP), intent(in):: xyz_VB    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ v (t-\Delta t) $ .   南北風速. Northward wind
    real(DP), intent(in):: xyz_TempB (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T (t-\Delta t) $ .   温度. Temperature
    real(DP), intent(in):: xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ q (t-\Delta t) $ .   比湿. Specific humidity
    real(DP), intent(in):: xy_PsB    (0:imax-1, 1:jmax)
                              ! $ p_s (t-\Delta t) $ . 地表面気圧. Surface pressure

    real(DP), intent(in):: xyz_UN   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ u (t) $ .     東西風速. Eastward wind
    real(DP), intent(in):: xyz_VN   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ v (t) $ .     南北風速. Northward wind
    real(DP), intent(in):: xyz_TempN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T (t) $ .     温度. Temperature
    real(DP), intent(in):: xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ q (t) $ .     比湿. Specific humidity
    real(DP), intent(in):: xy_PsN    (0:imax-1, 1:jmax)
                              ! $ p_s (t) $ .   地表面気圧. Surface pressure

    real(DP), intent(in):: xyz_DUDtPhy    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \left(\DP{u}{t}\right)^{phy} $ . 
                              ! 外力項 (物理過程) による東西風速変化. 
                              ! Eastward wind tendency by external force terms (physical processes)
    real(DP), intent(in):: xyz_DVDtPhy    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \left(\DP{v}{t}\right)^{phy} $ . 
                              ! 外力項 (物理過程) による南北風速変化. 
                              ! Northward wind tendency by external force terms (physical processes)
    real(DP), intent(in):: xyz_DTempDtPhy (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \left(\DP{T}{t}\right)^{phy} $ . 
                              ! 外力項 (物理過程) による温度変化. 
                              ! Temperature tendency by external force terms (physical processes)
    real(DP), intent(in):: xyzf_DQMixDtPhy(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ \left(\DP{q}{t}\right)^{phy} $ . 
                              ! 外力項 (物理過程) による比湿変化. 
                              ! Temperature tendency by external force terms (physical processes)

    real(DP), intent(in):: xy_SurfHeight (0:imax-1, 1:jmax)
                              ! $ z_s $ . 地表面高度. 
                              ! Surface height. 

    real(DP), intent(out):: xyz_UA    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ u (t+\Delta t) $ .   東西風速. Eastward wind
    real(DP), intent(out):: xyz_VA    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ v (t+\Delta t) $ .   南北風速. Northward wind
    real(DP), intent(out):: xyz_TempA (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T (t+\Delta t) $ .   温度. Temperature
    real(DP), intent(out):: xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ q (t+\Delta t) $ .   比湿. Specific humidity
    real(DP), intent(out):: xy_PsA    (0:imax-1, 1:jmax)
                              ! $ p_s (t+\Delta t) $ . 地表面気圧. Surface pressure
    real(DP), intent(out), optional :: xyz_ArgOMG(0:imax-1, 1:jmax, 1:kmax)
                              !
                              ! OMEGA, DP/Dt

    ! 作業変数
    ! Work variables
    !
    real(DP):: xyz_UCosLatN   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ U (t) = u (t) \cos \varphi $ .
    real(DP):: xyz_VCosLatN   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ V (t) = v (t) \cos \varphi $ .
    real(DP):: xyz_VorN       (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \zeta (t) $ . 渦度. Vorticity
    real(DP):: xyz_DivN       (0:imax-1, 1:jmax, 1:kmax)
                              ! $ D (t) $ .     発散. Divergence

    real(DP):: xy_PiN         (0:imax-1, 1:jmax)
                              ! $ \pi = \ln p_s $
    real(DP):: wz_DivN        (lmax, 1:kmax)
                              ! $ D (t) $ . 発散 (スペクトル). 
                              ! Divergence (spectral)
    real(DP):: wz_TempN       (lmax, 1:kmax)
                              ! $ T (t) $ . 温度 (スペクトル). 
                              ! Temperature (spectral)
    real(DP):: w_PiN          (lmax)
                              ! $ \pi $ スペクトル
    real(DP):: xy_GradLambdaPiN (0:imax-1, 1:jmax)
                              ! $ \DP{\pi}{\lambda} $
    real(DP):: xy_GradMuPiN     (0:imax-1, 1:jmax)
                              ! $ (1-\mu^2) \DP{\pi}{\mu} $

    real(DP):: xyz_PiAdv (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \Dvect{v} \cdot \nabla \pi $ . 
                              ! $ \pi $ の移流. Advection of $ \pi $

    real(DP):: xyz_UAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ U_A (t) $ . 東西運動量移流項. 
                              ! Eastward advection of momentum
    real(DP):: xyz_VAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ V_A (t) $ . 南北運動量移流項. 
                              ! Northward advection of momentum
    real(DP):: xyz_TempNonLinearN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ H (t) $ . 温度時間変化項. 
                              ! Temperature tendency
    real(DP):: xyzf_QMixNonLinearN (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ R (t) $ . 比湿時間変化項. 
                              ! Specific humidity tendency
    real(DP):: xyz_KinEngyN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ KE (t) $ . 運動エネルギー項. 
                              ! Kinetic energy
    real(DP):: xyz_TempUAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ UT' (t) $ . 温度東西移流項. 
                              ! Eastward advection of temperature
    real(DP):: xyz_TempVAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ VT' (t) $ . 温度南北移流項. 
                              ! Northward advection of temperature
    real(DP):: xyr_SigDotN (0:imax-1, 1:jmax, 0:kmax)
                              ! $ \dot{\sigma} $ .
                              ! 鉛直流. Vertical flow

    real(DP):: xy_DPiDtNG (0:imax-1, 1:jmax)
                              ! $ Z $ . 地表面気圧時間変化非重力波項. 
                              ! Non-gravity wave component of surface pressure tendency

    real(DP):: xyzf_QMixUAdvN (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ Uq (t) $ . 比湿東西移流項. 
                              ! Eastward advection of specific humidity
    real(DP):: xyzf_QMixVAdvN (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ Vq (t) $ . 比湿南北移流項. 
                              ! Northward advection of specific humidity
    real(DP):: wz_DVorDtNG (lmax, 1:kmax)
                              ! $ \left( \DD{\zeta}{t} (t) \right)^{NG}$ . 渦度変化の非重力波成分 (スペクトル). 
                              ! Non-gravity wave component of vorticity tendency (spectral)
    real(DP):: wz_DDivDtNG (lmax, 1:kmax)
                              ! $ \left( \DD{D}{t} (t) \right)^{NG}$ . 発散変化の非重力波成分 (スペクトル). 
                              ! Non-gravity wave component of divergence tendency (spectral)
    real(DP):: wz_DTempDtNG (lmax, 1:kmax)
                              ! $ \left( \DD{T}{t} (t) \right)^{NG}$ . 温度変化の非重力波成分 (スペクトル). 
                              ! Non-gravity wave component of temperature tendency (spectral)
    real(DP):: wzf_DQMixDtN(lmax, 1:kmax, 1:ncmax)
                              ! $ \DD{q}{t} (t) $ . 比湿変化 (スペクトル). 
                              ! Specific humidity tendency (spectral)
    real(DP):: w_DPiDtNG(lmax)
                              ! $ \left( \DD{p_s}{t} (t) \right)^{NG}$ . 地表面気圧変化費重力波項 (スペクトル). 
                              ! Non-gravity wave component of surface pressure tendency (spectral)
    real(DP):: xyz_UCosLatB   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ U (t-\Delta t) = u (t-\Delta t) \cos \varphi $ .
    real(DP):: xyz_VCosLatB   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ V (t-\Delta t) = v (t-\Delta t) \cos \varphi $ .
    real(DP):: wz_VorB (lmax, 1:kmax)
                              ! $ \zeta (t-\Delta t) $ . 渦度 (スペクトル). 
                              ! Vorticity (spectral)
    real(DP):: wz_DivB (lmax, 1:kmax)
                              ! $ D (t-\Delta t) $ . 発散 (スペクトル). 
                              ! Divergence (spectral)
    real(DP):: wz_TempB (lmax, 1:kmax)
                              ! $ T (t-\Delta t) $ . 温度 (スペクトル). 
                              ! Temperature (spectral)
    real(DP):: w_PiB (lmax)
                              ! $ \pi = \ln p_s (t-\Delta t) $ . 地表面気圧 (スペクトル). 
                              ! Surface pressure (spectral)
    real(DP):: wzf_QMixB(lmax, 1:kmax, 1:ncmax)
                              ! $ q (t-\Delta t) $ . 比湿 (スペクトル). 
                              ! Specific humidity (spectral)
    real(DP):: xyz_DUDtPhyCosLat   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \left(\DP{u}{t}\right)^{phy} \cos \varphi $ .
    real(DP):: xyz_DVDtPhyCosLat   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \left(\DP{v}{t}\right)^{phy} \cos \varphi $ . 
    real(DP):: wz_VorA (lmax, 1:kmax)
                              ! $ \zeta (t+\Delta t) $ . 渦度 (スペクトル). 
                              ! Vorticity (spectral)
    real(DP):: wz_DivA (lmax, 1:kmax)
                              ! $ D (t+\Delta t) $ . 発散 (スペクトル). 
                              ! Divergence (spectral)
    real(DP):: wz_TempA (lmax, 1:kmax)
                              ! $ T (t+\Delta t) $ . 温度 (スペクトル). 
                              ! Temperature (spectral)
    real(DP):: w_PiA (lmax)
                              ! $ \pi = \ln p_s (t+\Delta t) $ . 地表面気圧 (スペクトル).
                              ! Surface pressure (spectral)
    real(DP):: wzf_QMixA(lmax, 1:kmax, 1:ncmax)
                              ! $ q (t+\Delta t) $ . 比湿 (スペクトル). 
                              ! Specific humidity (spectral)

    real(DP):: wz_Psi (lmax, 1:kmax)
                              ! $ \psi $ . 流線関数. Streamline function 
    real(DP):: wz_Chi (lmax, 1:kmax)
                              ! $ \chi $ . ポテンシャル. Potential
    real(DP):: xyz_UCosLatA   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ U (t+\Delta t) = u (t+\Delta t) \cos \varphi $ .
    real(DP):: xyz_VCosLatA   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ V (t+\Delta t) = v (t+\Delta t) \cos \varphi $ .

    real(DP):: wz_VorDiffA (lmax, 1:kmax)
                              ! $ \mathscr{D}(\zeta) $ . 
                              ! 運動量水平拡散による渦度変化 (スペクトル). 
                              ! Vorticity tendency by 
                              ! horizontal momentum diffusion (spectral)
    real(DP):: wz_DivDiffA (lmax, 1:kmax)
                              ! $ \mathscr{D}(D) $ . 
                              ! 運動量水平拡散による発散変化 (スペクトル). 
                              ! Divergence tendency by 
                              ! horizontal momentum diffusion (spectral)
    real(DP):: wz_PsiDiff (lmax, 1:kmax)
                              ! 運動量水平拡散による流線関数 $ \psi $ 変化
                              ! Streamline function tendency by 
                              ! horizontal momentum diffusion
    real(DP):: wz_ChiDiff (lmax, 1:kmax)
                              ! 運動量水平拡散によるポテンシャル $ \chi $ 変化
                              ! Potential tendency by 
                              ! horizontal momentum diffusion
    real(DP):: xyz_UDiff (0:imax-1, 1:jmax, 1:kmax)
                              ! 運動量水平拡散による東西風変化. 
                              ! Eastward wind tendency by 
                              ! horizontal momentum diffusion
    real(DP):: xyz_VDiff (0:imax-1, 1:jmax, 1:kmax)
                              ! 運動量水平拡散による南北風変化. 
                              ! Northward wind tendency by 
                              ! horizontal momentum diffusion

    ! 地表面高度関連
    ! Surface height etc. 
    !
    real(DP):: w_SurfGeoPot (lmax)
                              ! $ \Phi_s $ . 地表ジオポテンシャル. 
                              ! Surface geo-potential

    ! Variables for mass fixing
    !
    real(DP):: MassB
    real(DP):: MassA
    real(DP):: xyr_PressA(0:imax-1, 1:jmax, 0:kmax)
    real(DP):: xyr_PressB(0:imax-1, 1:jmax, 0:kmax)

    real(DP) :: xyz_PressA       (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xy_HDifPsA       (0:imax-1, 1:jmax)
    real(DP) :: xyz_HDifPressA   (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyzf_DQMixDtHDCor(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
    integer  :: kp, kn

    real(DP) :: xyz_OMG          (0:imax-1, 1:jmax, 1:kmax)
                              !
                              ! OMEGA, DP/Dt
    real(DP) :: xyzf_QMixTentative(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
    logical  :: FlagExistNonAdvComp

    real(DP) :: xyz_UTracerAdv(0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_VTracerAdv(0:imax-1, 1:jmax, 1:kmax)

    real(DP) :: xyr_AdvTestStreamFunc(0:imax-1, 1:jmax, 0:kmax)

    integer:: k               ! 鉛直方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in vertical direction
    integer:: n               ! 組成方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in dimension of constituents

    ! 実行文 ; Executable statement
    !

    ! 初期化確認
    ! Initialization check
    !
    if ( .not. dynamics_hspl_vas83_inited ) then
      call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
    end if


    ! 計算時間計測開始
    ! Start measurement of computation time
    !
    call TimesetClockStart( module_name )


    ! TimeIntegration で使用する係数の設定
    ! Configure coefficients for "TimeIntegration"
    !
    call SemiImplMatrix

    ! 格子点値をスペクトル値へ ( $ t-\Delta t$ )
    ! Exchange grid values to spectral values ( $ t-\Delta t$ )
    !
    !   渦度発散の計算
    !   Calculate vorticity and divergence
    do k = 1, kmax
      xyz_UCosLatB(:,:,k) = xyz_UB(:,:,k) * xy_CosLat
      xyz_VCosLatB(:,:,k) = xyz_VB(:,:,k) * xy_CosLat
    end do
    wz_VorB  =   (   wa_DivLambda_xya( xyz_VCosLatB ) - wa_DivMu_xya    ( xyz_UCosLatB ) ) / RPlanet
    wz_DivB  =   (   wa_DivLambda_xya( xyz_UCosLatB ) + wa_DivMu_xya    ( xyz_VCosLatB ) ) / RPlanet
    !
    wz_TempB = wa_xya( xyz_TempB )
    !
    w_PiB    = w_xy( log( xy_PsB ) )


    ! 格子点値をスペクトル値へ ( $ t $ )
    ! Exchange grid values to spectral values ( $ t $ )
    !
    !   渦度発散の計算
    !   Calculate vorticity and divergence
    ! 
    do k = 1, kmax
      xyz_UCosLatN(:,:,k) = xyz_UN(:,:,k) * xy_CosLat
      xyz_VCosLatN(:,:,k) = xyz_VN(:,:,k) * xy_CosLat
    end do
    xyz_VorN = xya_wa(   wa_DivLambda_xya( xyz_VCosLatN ) - wa_DivMu_xya    ( xyz_UCosLatN ) ) / RPlanet
    wz_DivN  =       (   wa_DivLambda_xya( xyz_UCosLatN ) + wa_DivMu_xya    ( xyz_VCosLatN ) ) / RPlanet
    xyz_DivN = xya_wa( wz_DivN )
    !
    select case ( IDTimeIntegScheme )
    case ( IDTimeIntegSchemeExplicit )
      wz_TempN = wa_xya( xyz_TempN )
    end select
    !
    xy_PiN = log( xy_PsN )
    w_PiN = w_xy( xy_PiN )
    xy_GradLambdaPiN = xy_GradLambda_w( w_PiN ) / RPlanet
    xy_GradMuPiN     = xy_GradMu_w    ( w_PiN ) / RPlanet


    ! 地表ジオポテンシャルの計算
    ! Calculate surface geo-potential
    !
    w_SurfGeoPot = w_xy( Grav * xy_SurfHeight )


    ! 格子点上での非線形力学項の計算
    ! Calculate non-linear dynamical terms on grid points
    !
    call NonLinearOnGrid( xyz_UCosLatN,     xyz_VCosLatN, xyz_VorN,         xyz_DivN, xyz_TempN,        xyzf_QMixN(:,:,:,IndexH2OVap), xy_GradLambdaPiN, xy_GradMuPiN, xyz_PiAdv, xyz_UAdvN,        xyz_VAdvN, xyz_TempNonLinearN, xyz_KinEngyN, xyz_TempUAdvN,    xyz_TempVAdvN, xyr_SigDotN,      xy_DPiDtNG )


    ! 非線形項と物理過程による時間変化項の和をスペクトル変換
    ! Spectral transformation of sum of non-linear term and tendencies by physical 
    ! processes 
    !
    do k = 1, kmax
      xyz_DUDtPhyCosLat(:,:,k) = xyz_DUDtPhy(:,:,k) * xy_CosLat
      xyz_DVDtPhyCosLat(:,:,k) = xyz_DVDtPhy(:,:,k) * xy_CosLat
    end do
    wz_DVorDtNG =   (   wa_DivLambda_xya( xyz_VAdvN + xyz_DVDtPhyCosLat ) - wa_DivMu_xya    ( xyz_UAdvN + xyz_DUDtPhyCosLat ) ) / RPlanet
    wz_DDivDtNG =   (   wa_DivLambda_xya( xyz_UAdvN + xyz_DUDtPhyCosLat ) + wa_DivMu_xya    ( xyz_VAdvN + xyz_DVDtPhyCosLat ) ) / RPlanet - wa_Lapla_wa( wa_xya(xyz_KinEngyN) ) / RPlanet**2
    wz_DTempDtNG = - (   wa_DivLambda_xya( xyz_TempUAdvN ) + wa_DivMu_xya    ( xyz_TempVAdvN ) ) / RPlanet + wa_xya( xyz_TempNonLinearN + xyz_DTempDtPhy )
    w_DPiDtNG = w_xy( xy_DPiDtNG )


    ! 時間積分
    ! Time integration
    !
    call TimeIntegration( w_SurfGeoPot, wz_DVorDtNG, wz_DDivDtNG, wz_DTempDtNG, w_DPiDtNG, wz_DivN,                  wz_TempN,     w_PiN, wz_VorB,     wz_DivB,     wz_TempB,     w_PiB, wz_VorA,     wz_DivA,     wz_TempA,     w_PiA )


    ! Divergence damping
    !
    call DivergenceDamping( wz_DivA )


    ! スペクトル値を格子点値へ ( $ t+\Delta t$ )
    ! Exchange spectral values to grid values ( $ t+\Delta t$ )
    !
    wz_Psi = wa_LaplaInv_wa( wz_VorA ) * RPlanet**2
    wz_Chi = wa_LaplaInv_wa( wz_DivA ) * RPlanet**2

    xyz_UCosLatA = (   xya_GradLambda_wa( wz_Chi ) - xya_GradMu_wa    ( wz_Psi ) ) / RPlanet

    xyz_VCosLatA = (   xya_GradLambda_wa( wz_Psi ) + xya_GradMu_wa    ( wz_Chi ) ) / RPlanet

    do k = 1, kmax
      xyz_UA(:,:,k) = xyz_UCosLatA(:,:,k) / xy_CosLat
      xyz_VA(:,:,k) = xyz_VCosLatA(:,:,k) / xy_CosLat
    end do
    xyz_TempA = xya_wa( wz_TempA )

    xy_PsA = exp( xy_w( w_PiA ) )

    ! 水平拡散とスポンジ層による散逸の補正
    ! Correction for dissipation of kinetic enery by horizontal diffusion and sponge 
    ! layer
    !

    ! 水平拡散とスポンジ層による渦度発散の時間変化
    ! Vorticity and divergence tendency by horizontal diffusion and damping in a sponge 
    ! layer
    !
    wz_VorDiffA = wz_VorA * wz_DisCoefM
    wz_DivDiffA = wz_DivA * wz_DisCoefM


    ! 水平拡散とスポンジ層による散逸の摩擦熱補正
    ! Correction for internal energy by horizontal diffusion and damping in a sponge 
    ! layer
    !
    wz_PsiDiff = wa_LaplaInv_wa( wz_VorDiffA ) * RPlanet**2
    wz_ChiDiff = wa_LaplaInv_wa( wz_DivDiffA ) * RPlanet**2

    xyz_UDiff = (   xya_GradLambda_wa( wz_ChiDiff ) - xya_GradMu_wa    ( wz_PsiDiff ) ) / RPlanet

    xyz_VDiff = (   xya_GradLambda_wa( wz_PsiDiff ) + xya_GradMu_wa    ( wz_ChiDiff ) ) / RPlanet

    do k = 1, kmax
      xyz_TempA(:,:,k) = xyz_TempA(:,:,k) - (   xyz_UA(:,:,k) * xyz_UDiff(:,:,k) + xyz_VA(:,:,k) * xyz_VDiff(:,:,k)   ) / xy_CosLat / CpDry * 2.0_DP * DelTime
    end do




    ! Mass fixer
    !   Total Mass
    !   NOTICE: It is assumed that the total atmospheric mass does not change.
    !
    if ( FlagTotMassFix ) then
      MassB = IntLonLat_xy( xy_PsB )
      MassA = IntLonLat_xy( xy_PsA )
      if ( MassA /= 0.0_DP ) then
        xy_PsA = MassB / MassA * xy_PsA
      end if
    end if


    ! Set velocities for tracer advection
    !
    xyz_UTracerAdv = xyz_UN
    xyz_VTracerAdv = xyz_VN

    ! 主に移流テストのために使う if 文
    !
    if ( FlagAdvTest ) then
      xyz_UA      = xyz_UB
      xyz_VA      = xyz_VB
      xyz_TempA   = xyz_TempB
      xy_PsA      = xy_PsB
      xyr_SigDotN = 0.0_DP

      !
      ! Utility module for advection test
      !
      call AdvTestSetHorVels( xyz_UTracerAdv, xyz_VTracerAdv )
      ! calculation of variables required in spectral advection calculation
      do k = 1, kmax
        xyz_UCosLatN(:,:,k) = xyz_UTracerAdv(:,:,k) * xy_CosLat
        xyz_VCosLatN(:,:,k) = xyz_VTracerAdv(:,:,k) * xy_CosLat
      end do
      wz_DivN  =       (   wa_DivLambda_xya( xyz_UCosLatN ) + wa_DivMu_xya    ( xyz_VCosLatN ) ) / RPlanet
      xyz_DivN = xya_wa( wz_DivN )

      call AdvTestSetVerVel( xyr_SigDotN, xyr_AdvTestStreamFunc )

      call HistoryAutoPut( TimeN, 'AdvTestStreamFunc', xyr_AdvTestStreamFunc )

    end if

    if ( FlagSLTT ) then

!!$    select case ( IDTTMethod )
!!$    case ( IDTTMethodSL )
      ! Semi-Lagrangian advection

      do k = 0, kmax
        xyr_PressA(:,:,k) = xy_PsA * r_Sigma(k)
        xyr_PressB(:,:,k) = xy_PsB * r_Sigma(k)
      end do
      ! Mass fixer is included in SLTTMain
      call SLTTMain( xyr_PressB, xyr_PressA, xyz_UTracerAdv, xyz_VTracerAdv, xyr_SigDotN, xyzf_DQMixDtPhy, xyzf_QMixB, xyzf_QMixA )

!!$    case ( IDTTMethodFV )
!!$      ! Finite volume method
!!$
!!$      do k = 0, kmax
!!$        xyr_PressA(:,:,k) = xy_PsA * r_Sigma(k)
!!$        xyr_PressB(:,:,k) = xy_PsB * r_Sigma(k)
!!$      end do
!!$      call TTFVMain(             &
!!$        & xyr_PressB, xyr_PressA, &
!!$        & xyz_UTracerAdv, xyz_VTracerAdv, xyr_SigDotN, & ! (in)
!!$        & xyzf_DQMixDtPhy,             & ! (in)
!!$        & xyzf_QMixB,                  & ! (in)
!!$        & xyzf_QMixA                   & ! (out)
!!$        & )

    else

!!$    case ( IDTTMethodSp )
      ! Spectral advection

      do n = 1, ncmax
        wzf_QMixB(:,:,n) = wa_xya( xyzf_QMixB(:,:,:,n) )
      end do
      !
      call NonLinearOnGridQMix( xyz_UCosLatN, xyz_VCosLatN, xyz_DivN, xyr_SigDotN, xyzf_QMixN, xyzf_QMixNonLinearN, xyzf_QMixUAdvN, xyzf_QMixVAdvN )

      do n = 1, ncmax
        wzf_DQMixDtN(:,:,n) = - (   wa_DivLambda_xya( xyzf_QMixUAdvN(:,:,:,n) ) + wa_DivMu_xya    ( xyzf_QMixVAdvN(:,:,:,n) ) ) / RPlanet + wa_xya(   xyzf_QMixNonLinearN(:,:,:,n) + xyzf_DQMixDtPhy    (:,:,:,n) )
      end do

      call TimeIntegrationQMix( wzf_DQMixDtN, wzf_QMixB, wzf_QMixA )

      do n = 1, ncmax
        xyzf_QMixA(:,:,:,n) = xya_wa( wzf_QMixA(:,:,n) )
      end do

      if ( FlagMassHorDifCor ) then

        ! Correction for horizontal diffusion of constituents.
        ! This is performed for the aim of correcting the diffusion in the viscinity of 
        ! steep terrain slope.
        !

        do k = 1, kmax
          xyz_PressA(:,:,k) = xy_PsA * z_Sigma(k)
        end do
        xy_HDifPsA = xy_w( w_DisCoefQ * w_xy( xy_PsA ) )
        do k = 1, kmax
          xyz_HDifPressA(:,:,k)  = z_Sigma(k) * xy_HDifPsA(:,:)
        end do

        do n = 1, ncmax
          do k = 1, kmax
            kp = k-1
            kn = k+1
            if( k == 1    ) kp = 1
            if( k == kmax ) kn = kmax
            xyzf_DQMixDtHDCor(:,:,k,n) = - ( xyzf_QMixA(:,:,kn,n) - xyzf_QMixA(:,:,kp,n) ) / ( xyz_PressA(:,:,kn)   - xyz_PressA(:,:,kp)   ) * xyz_HDifPressA(:,:,k)
          end do
        end do
        xyzf_QMixA = xyzf_QMixA + xyzf_DQMixDtHDCor * 2.0_DP * DelTime
      endif

      ! Mass fixer
      !   Constituents
      !
      do k = 0, kmax
        xyr_PressA(:,:,k) = xy_PsA * r_Sigma(k)
        xyr_PressB(:,:,k) = xy_PsB * r_Sigma(k)
      end do
      call MassFixer( xyr_PressA, xyzf_QMixA, xyr_PressRef = xyr_PressB, xyzf_QMixRef = ( xyzf_QMixB+xyzf_DQMixDtPhy*2.0_DP*DelTime ) )
!!$    end select

    end if

    ! 主に移流テストのために使う if 文
    !
!!$    if ( FlagAdvTest ) then
!!$      xyr_SigDotN = 0.0_DP
!!$    end if


    ! Calculation for non-advection case
    ! オイラー移流計算のとき、ピーキーな分布の物質を移流させない。計算安定の
    ! ため。
    !
    FlagExistNonAdvComp = .false.
    do n = 1, ncmax
      if ( .not. CompositionInqFlagAdv( n ) ) then
        xyzf_QMixA(:,:,:,n) = xyzf_QMixB(:,:,:,n) + xyzf_DQMixDtPhy(:,:,:,n) * 2.0_DP * DelTime
        xyzf_QMixTentative(:,:,:,n) = xyzf_QMixA(:,:,:,n)
        FlagExistNonAdvComp = .true.
      else
        xyzf_QMixTentative(:,:,:,n) = 0.0_DP
      end if
    end do
    if ( FlagExistNonAdvComp ) then
      ! Mass fixer
      !   Constituents
      call MassFixerColumn( xyr_PressB, xyzf_QMixTentative, xyr_PressRef = xyr_PressB, xyzf_QMixRef = ( xyzf_QMixB+xyzf_DQMixDtPhy*2.0_DP*DelTime ) )
      do n = 1, ncmax
        if ( .not. CompositionInqFlagAdv( n ) ) then
          xyzf_QMixA(:,:,:,n) = xyzf_QMixTentative(:,:,:,n)
        end if
      end do
    end if


    ! ヒストリデータ出力
    ! History data output
    !
    call OutputDiagnosedVariables( xyz_UB, xyz_VB, xyz_TempB, xyzf_QMixB, xy_PsB, xyz_UN, xyz_VN, xyz_TempN, xyzf_QMixN, xy_PsN, xyz_UA, xyz_VA, xyz_TempA, xyzf_QMixA, xy_PsA, xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyzf_DQMixDtPhy, xyr_SigDotN, xy_DPiDtNG, xyz_PiAdv, xyz_OMG )


    if ( present( xyz_ArgOMG ) ) then
      xyz_ArgOMG = xyz_OMG
    end if


    ! 計算時間計測一時停止
    ! Pause measurement of computation time
    !
    call TimesetClockStop( module_name )

  end subroutine DynamicsHSplVAS83
Subroutine :

モジュール内部の変数の割り付け解除を行います.

Deallocate variables in this module.

[Source]

  subroutine DynamicsHSplVAS83Finalize
    !
    ! モジュール内部の変数の割り付け解除を行います. 
    !
    ! Deallocate variables in this module. 
    !

    ! 宣言文 ; Declaration statements
    !
    implicit none

    ! 実行文 ; Executable statement
    !

    if ( .not. dynamics_hspl_vas83_inited ) return

    ! デフォルト値へ戻す
    ! Return to default values
    !
    DelTimeSave = - 1.0_DP

    ! 割り付け解除
    ! Deallocation
    !
    if ( allocated( xy_SinLat ) )    deallocate( xy_SinLat )
    if ( allocated( xy_CosLat ) )    deallocate( xy_CosLat )

    if ( allocated( xy_Cori ) )      deallocate( xy_Cori )

    if ( allocated( z_HydroAlpha ) ) deallocate( z_HydroAlpha )
    if ( allocated( z_HydroBeta  ) ) deallocate( z_HydroBeta  )

    if ( allocated( z_TInpCoefA ) )  deallocate( z_TInpCoefA )
    if ( allocated( z_TInpCoefB ) )  deallocate( z_TInpCoefB )
    if ( allocated( z_TInpCoefK ) )  deallocate( z_TInpCoefK )

    if ( allocated( z_RefTemp ) )    deallocate( z_RefTemp )
    if ( allocated( r_RefTemp ) )    deallocate( r_RefTemp )

    if ( allocated( w_LaplaEigVal ) ) deallocate( w_LaplaEigVal )

    if ( allocated( wz_DisCoefM   ) ) deallocate( wz_DisCoefM     )
    if ( allocated( wz_DisCoefH   ) ) deallocate( wz_DisCoefH     )
    if ( allocated( w_DisCoefQ    ) ) deallocate( w_DisCoefQ      )

    if ( allocated( z_siMtxC      ) )  deallocate( z_siMtxC       )
    if ( allocated( z_siMtxG      ) )  deallocate( z_siMtxG       )
    if ( allocated( zz_siMtxH     ) )  deallocate( zz_siMtxH      )
    if ( allocated( zz_siMtxDiH   ) )  deallocate( zz_siMtxDiH    )
    if ( allocated( wzz_siMtxWDiH ) )  deallocate( wzz_siMtxWDiH  )
    if ( allocated( zz_siMtxGCt   ) )  deallocate( zz_siMtxGCt    )

    if ( allocated( zz_siMtxW  ) )  deallocate( zz_siMtxW  )
    if ( allocated( zz_siMtxQ  ) )  deallocate( zz_siMtxQ  )
    if ( allocated( zz_siMtxS  ) )  deallocate( zz_siMtxS  )
    if ( allocated( zz_siMtxR  ) )  deallocate( zz_siMtxR  )

!!$    if ( allocated(nmo           ) )  deallocate( nmo            )
!!$    if ( allocated(wzz_siMtxM    ) )  deallocate( wzz_siMtxM     )
!!$    if ( allocated(z_siMtxPivWork) )  deallocate( z_siMtxPivWork )
    if ( allocated(wzz_siMtxLU   ) )  deallocate( wzz_siMtxLU    )
    if ( allocated(wz_siMtxPiv   ) )  deallocate( wz_siMtxPiv    )

    dynamics_hspl_vas83_inited = .false.

  end subroutine DynamicsHSplVAS83Finalize
Subroutine :

計算に必要なパラメタの設定や NAMELIST#dynamics_hspl_vas83_nml の読み込みを行います.

Configure parameters for calculation, and load "NAMELIST#dynamics_hspl_vas83_nml"

This procedure input/output NAMELIST#dynamics_hspl_vas83_nml .

[Source]

  subroutine DynamicsHSplVAS83Init
    !
    ! 計算に必要なパラメタの設定や NAMELIST#dynamics_hspl_vas83_nml
    ! の読み込みを行います. 
    ! 
    ! Configure parameters for calculation, 
    ! and load "NAMELIST#dynamics_hspl_vas83_nml"
    !

    ! モジュール引用 ; USE statements
    !

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: RPlanet, Omega, GasRDry, CpDry
                              ! $ C_p $ [J kg-1 K-1]. 
                              ! 乾燥大気の定圧比熱. 
                              ! Specific heat of air at constant pressure

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: nmax, lmax, imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

    ! 組成に関わる配列の設定
    ! Settings of array for atmospheric composition
    !
    use composition, only: ncmax, a_QMixName

    ! 座標データ設定
    ! Axes data settings
    !
    use axesset, only: z_Sigma, r_Sigma, z_DelSigma, AxNameX, AxNameY, AxNameZ, AxNameR, AxNameT

    ! SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応) 
    ! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
    !
#ifdef LIB_MPI
#ifdef SJPACK
    use wa_mpi_module_sjpack, only: xy_Lat => xv_Lat, w_xy   => w_xv, rn, nm_l
#else
    use wa_mpi_module, only: xy_Lat => xv_Lat, w_xy   => w_xv, rn, nm_l
#endif
#elif AXISYMMETRY
    use wa_zonal_module, only: xy_Lat, w_xy, rn, nm_l
#elif SJPACK
    use wa_module_sjpack, only: xy_Lat, w_xy, rn, nm_l
#elif AXISYMMETRY_SJPACK
    use wa_zonal_module_sjpack, only: xy_Lat, w_xy, rn, nm_l
#else
    use wa_module, only: xy_Lat, w_xy, rn, nm_l
#endif

    ! NAMELIST ファイル入力に関するユーティリティ
    ! Utilities for NAMELIST file input
    !
    use namelist_util, only: namelist_filename, NmlutilMsg

    ! ヒストリデータ出力
    ! History data output
    !
    use gtool_historyauto, only: HistoryAutoAddVariable

    ! gtool4 データ入力
    ! Gtool4 data input
    !
    use gtool_history, only: HistoryGet

    ! ファイル入出力補助
    ! File I/O support
    !
    use dc_iounit, only: FileOpen

    ! 種別型パラメタ
    ! Kind type parameter
    !
    use dc_types, only: STDOUT, STRING                ! 文字列.       Strings. 

    ! 日付および時刻の取り扱い
    ! Date and time handler
    !
    use dc_calendar, only: DCCalConvertByUnit

    ! 組み込み関数 PRESENT の拡張版関数
    ! Extended functions of intrinsic function "PRESENT"
    !
    use dc_present, only: present_and_not_empty

    ! 文字列操作
    ! Character handling
    !
    use dc_string, only: LChar

    ! 時刻管理
    ! Time control
    !
    use timeset, only: TimeN                 ! ステップ $ t $ の時刻. Time of step $ t $. 

    ! 質量の補正
    ! Mass fixer
    !
    use mass_fixer, only : MassFixerInit

    ! セミラグランジュ法による物質移流計算
    ! Semi-Lagrangian method for tracer transport
    use sltt, only: SLTTInit

!!$    ! 物質移流 (有限体積法)
!!$    ! Tracer Transport (finite volume method)
!!$    !
!!$    use tt_fv, only: TTFVInit

    !
    ! Utility module for advection test
    !
    use adv_test, only: AdvTestInit


    ! 宣言文 ; Declaration statements
    !
    implicit none


    ! 基準温度の設定のための作業変数
    ! Work variable for reference temperature
    !
    character(TOKEN) :: TimeIntegScheme
                              ! 時間積分法. 
                              ! 以下の方法を選択可能. 
                              !
                              ! Time integration scheme. 
                              ! Available schemes are as follows.
                              !
                              ! * "Semi-implicit"
                              ! * "Explicit"

    real(DP):: RefTemp
                              ! $ \overline{T} $ . 基準温度. 
                              ! Reference temperature

    ! 水平拡散, スポンジ層のための作業変数
    ! Work variable for horizontal diffusion and sponge layer
    !
    real(DP) :: w_HDifCoefM        (lmax)
                              ! $ {\cal D}_M = -K_{HD} [(-1)^{N_D/2}\nabla^{N_D}- (2/a^2)^{N_D/2}] $ . 
                              ! 運動量水平拡散係数. 
                              ! Coefficient of horizontal momentum diffusion
    real(DP) :: w_HDifCoefH        (lmax)
                              ! $ {\cal D}_H = -(-1)^{N_D/2}K_{HD}\nabla^{N_D} $ . 
                              ! 熱, 水水平拡散係数. 
                              ! Coefficient of horizontal thermal and water diffusion
    real(DP) :: wz_SpongeLayerCoefM(lmax, 1:kmax)
                              ! スポンジ層における運動量の減衰係数
                              ! Damping coefficient for momentum in a sponge layer
    real(DP) :: wz_SpongeLayerCoefH(lmax, 1:kmax)
                              ! スポンジ層における熱の減衰係数
                              ! Damping coefficient for temperature zonal perturbation 
                              ! in a sponge layer

    integer:: HDOrder           ! 超粘性の次数.  Order of hyper-viscosity
    real(DP):: VisCoef          ! 超粘性係数. Hyper-viscosity coefficient
    real(DP):: HDEFoldTimeValue ! 最大波数に対する e-folding time. 
                                ! 負の値を与えると, 水平拡散係数をゼロにします. 
                                ! 
                                ! E-folding time for maximum wavenumber. 
                                ! If negative value is given, 
                                ! coefficients of horizontal diffusion become zero. 
    character(TOKEN):: HDEFoldTimeUnit
                                ! 最大波数に対する e-folding time の単位. 
                                ! Unit of e-folding time for maximum wavenumber
    real(DP):: HDEFoldTime      ! 最大波数に対する e-folding time [単位: 秒]. 
                                ! E-folding time for maximum wavenumber [Unit: sec]

    logical          :: FlagSpongeLayer
    logical          :: FlagSpongeLayerforZonalMean
    logical          :: FlagSpongeLayerforHeat
    real(DP)         :: SLEFoldTimeValue
                                ! スポンジ層の最上層における減衰時定数
                                ! Damping time constant at the top model layer 
                                ! in a sponge layer
    character(TOKEN) :: SLEFoldTimeUnit
                                ! SLEFoldTimeValue の単位
                                ! Unit of SLEFoldTimeValue
    real(DP)         :: SLEFoldTime
                                ! スポンジ層の最上層における減衰時定数 (秒)
                                ! Damping time constant at the top model layer 
                                ! in a sponge layer (sec)
    integer          :: SLOrder
                                ! スポンジ層の減衰係数の sigma 依存性のオーダ
                                ! Sigma dependence of damping coefficient 
                                ! in a sponge layer
    integer          :: SLNumLayer
                                ! スポンジ層が適応されるモデル上端からの層の数
                                ! Number of layers which the sponge layer is applied.
    integer          :: a_DegOrd(2)


    real(DP)         :: DivDampPeriodValue
                                ! Period for divergence damping application
    character(TOKEN) :: DivDampPeriodUnit
                                ! Unit of DivDampPeriodValue

    character(TOKEN):: TracerTransportMethod

    ! NonLinearOnGrid 等で使用する係数の設定のための作業変数
    ! Work variable for coefficients for "NonLinearOnGrid", etc.
    !
    real(DP):: Kappa          ! $ \kappa = R/C_p $ .
                              ! 気体定数の定圧比熱に対する比. 
                              ! Ratio of gas constant to specific heat

    integer:: k               ! 鉛直方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in vertical direction
    integer:: n
    integer:: l               ! 波数方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in wavenumber

    integer:: unit_nml        ! NAMELIST ファイルオープン用装置番号. 
                              ! Unit number for NAMELIST file open
    integer:: iostat_nml      ! NAMELIST 読み込み時の IOSTAT. 
                              ! IOSTAT of NAMELIST read

    ! NAMELIST 変数群
    ! NAMELIST group name
    !
    namelist /dynamics_hspl_vas83_nml/ TimeIntegScheme, HDOrder, HDEFoldTimeValue, HDEFoldTimeUnit, FlagSpongeLayer, FlagSpongeLayerforZonalMean, FlagSpongeLayerforHeat, SLEFoldTimeValue, SLEFoldTimeUnit, SLOrder, SLNumLayer, RefTemp, FlagDivDamp, DivDampPeriodValue, DivDampPeriodUnit, DivDampTime0, FlagTotMassFix, FlagMassHorDifCor, FlagSLTT, FlagAdvTest
          !
          ! デフォルト値については初期化手続 "dynamics_hspl_vas83#DynamicsInit" 
          ! のソースコードを参照のこと. 
          !
          ! Refer to source codes in the initialization procedure
          ! "dynamics_hspl_vas83#DynamicsInit" for the default values. 
          !

    ! 実行文 ; Executable statement
    !

    if ( dynamics_hspl_vas83_inited ) return


    ! デフォルト値の設定
    ! Default values settings
    !
    TimeIntegScheme             = 'Semi-implicit'

    HDOrder                     = 8
    HDEFoldTimeValue            = 8640.0_DP
    HDEFoldTimeUnit             = 'sec'

    FlagSpongeLayer             = .false.
    FlagSpongeLayerforZonalMean = .false.
    FlagSpongeLayerforHeat      = .false.
    SLEFoldTimeValue            = 86400.0_DP
    SLEFoldTimeUnit             = 'sec'
    SLOrder                     = 1
    SLNumLayer                  = kmax

    RefTemp                     = 300.

    FlagDivDamp                 = .false.
    DivDampPeriodValue          = 2.0_DP
    DivDampPeriodUnit           = 'day'
    DivDampTime0                = 0.0_DP

    FlagTotMassFix              = .true.

    FlagMassHorDifCor           = .false.

!!$    TracerTransportMethod       = "spectrum"
    FlagSLTT                    = .false.
    FlagAdvTest                 = .false.

    ! NAMELIST の読み込み
    ! NAMELIST is input
    !
    if ( trim(namelist_filename) /= '' ) then
      call FileOpen( unit_nml, namelist_filename, mode = 'r' ) ! (in)

      rewind( unit_nml )
      read( unit_nml, nml = dynamics_hspl_vas83_nml, iostat = iostat_nml )        ! (out)
      close( unit_nml )

      call NmlutilMsg( iostat_nml, module_name ) ! (in)
      if ( iostat_nml == 0 ) write( STDOUT, nml = dynamics_hspl_vas83_nml )
    end if

    ! 時間積分法のチェック
    ! Check time integration scheme
    !
    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('semi-implicit')
      IDTimeIntegScheme = IDTimeIntegSchemeSemiImplicit
    case ('explicit')
      IDTimeIntegScheme = IDTimeIntegSchemeExplicit
    case default
      call MessageNotify( 'E', module_name, '"TimeIntegScheme" must be "Semi-Implicit" or "Explicit".' )
    end select


    !
    ! Check tracer transport method
    !
!!$    select case ( LChar( trim( TracerTransportMethod ) ) )
!!$    case ('spectrum')
!!$      IDTTMethod = IDTTMethodSp
!!$    case ('semi-lagrange')
!!$      IDTTMethod = IDTTMethodSL
!!$    case ('finite volume')
!!$      IDTTMethod = IDTTMethodFV
!!$    case default
!!$      call MessageNotify( 'E', module_name, &
!!$        & '"TracerTransportMethod" must be "spectrum", "semi-lagrange", or "finite volume".' )
!!$    end select


    ! SemiImplMatrix サブルーチン用に $ \Delta t $ の保存値に初期値を設定. 
    ! Configure initial value to saved value of $ \Delta t $ for a subroutine "SemiImplMatrix"
    !
    DelTimeSave = - 1.0_DP

    ! NonLinearOnGrid 等で使用する係数の設定
    ! Configure coefficients for "NonLinearOnGrid", etc.
    !

    ! $ \sin \varphi $ と $ \cos \varphi $ の計算
    ! Calculate $ \sin \varphi $ and $ \cos \varphi $
    !
    allocate( xy_SinLat (0:imax-1, 1:jmax) )
    allocate( xy_CosLat (0:imax-1, 1:jmax) )
    xy_SinLat = sin( xy_Lat )
    xy_CosLat = cos( xy_Lat )


    ! コリオリパラメータの計算
    ! Calculate Coriolis parameter
    !
    allocate( xy_Cori (0:imax-1, 1:jmax) )
    xy_Cori = 2.0_DP * Omega * xy_SinLat


    ! 静水圧の式の係数 $ \alpha $ , $ \beta $ の計算
    ! Calculate coefficients $ \alpha $ and $ \beta $ in hydrostatic equation.
    !
    allocate( z_HydroAlpha(1:kmax) )
    allocate( z_HydroBeta (1:kmax) )

    Kappa = GasRDry / CpDry

    do k = 1, kmax
      z_HydroAlpha(k) = ( r_Sigma(k-1) / z_Sigma(k) ) ** Kappa - 1.0_DP

      z_HydroBeta(k) = 1.0_DP - ( r_Sigma(k) / z_Sigma(k) ) ** Kappa
    enddo

    ! 温度鉛直補間の係数 $ \kappa $, $ a $ , $ b $ の計算
    ! Calculate coefficients $ \kappa $, $ a $ , $ b $ 
    !   for interpolation of temperature
    !
    allocate( z_TInpCoefA (1:kmax) )
    allocate( z_TInpCoefB (1:kmax) )
    allocate( z_TInpCoefK (1:kmax) )

    do k = 1, kmax
      z_TInpCoefK(k) = (   r_Sigma(k-1) * z_HydroAlpha(k) + r_Sigma(k  ) * z_HydroBeta (k) ) / z_DelSigma(k)
    enddo

    z_TInpCoefA = 0.
    do k = 2, kmax
      z_TInpCoefA(k) = z_HydroAlpha(k) * ( 1.0_DP - (z_Sigma(k) / z_Sigma(k-1)) ** Kappa ) ** ( -1 )
    end do

    z_TInpCoefB = 0.
    do k = 1, kmax - 1
      z_TInpCoefB(k) = z_HydroBeta(k) * ( (z_Sigma(k) / z_Sigma(k+1) ) ** Kappa - 1.0_DP ) ** ( -1 )
    end do

    ! 基準温度 (半整数レベル) の計算
    ! Calculate reference temperature on half levels
    !
    allocate( z_RefTemp(1:kmax) )
    allocate( r_RefTemp(0:kmax) )

    z_RefTemp       = RefTemp
    r_RefTemp(0)    = 0.
    r_RefTemp(kmax) = 0.

    do k = 1, kmax - 1
      r_RefTemp(k) = z_TInpCoefA(k+1) * z_RefTemp(k+1) + z_TInpCoefB( k ) * z_RefTemp( k )
    enddo


    ! ルジャンドル陪関数の固有値の設定
    ! Set eigen value of Associated Legendre Function
    !
    allocate( w_LaplaEigVal(lmax) )

    w_LaplaEigVal(:) = rn(:,1) / RPlanet**2


    ! 水平拡散係数の設定
    ! Configure coefficient of horizontal diffusion
    !
    ! 粘性係数の計算 (最大波数の e-folding time が HDEFoldTime となるように)
    ! Calculate viscosity coefficient
    !
    HDEFoldTime = DCCalConvertByUnit( HDEFoldTimeValue, HDEFoldTimeUnit, 'sec' ) ! (in)

    if ( HDEFoldTimeValue > 0.0_DP ) then
      VisCoef = ( (nmax*(nmax+1)) / RPlanet**2 )**(-HDOrder / 2) / HDEFoldTime
    else
      VisCoef = 0.0_DP
    end if

    w_HDifCoefH = - VisCoef * ( ( - w_LaplaEigVal )**( HDOrder / 2 ) )

    w_HDifCoefM = w_HDifCoefH - VisCoef * ( - ( 2.0_DP / RPlanet**2 )**( HDOrder / 2 ) )

    ! スポンジ層の減衰係数の設定
    ! Set damping coefficient in a sponge layer
    !
    if ( SLEFoldTimeValue <= 0.0_DP ) then
      call MessageNotify( 'E', module_name, 'SLEFoldTimeValue must be greater than zero, SLEFoldTimeValue=<%f>.', d = (/ SLEFoldTimeValue /) )
    end if
    if ( ( SLNumLayer <= 0 ) .or. ( SLNumLayer > kmax ) ) then
      call MessageNotify( 'E', module_name, 'SLNumLayer must be greater than zero and less than/equal to kmax, SLNumLayer=<%d>.', i = (/ SLNumLayer /) )
    end if

    ! スポンジ層の減衰係数の計算
    ! Calculate damping coefficient for momentum in a sponge layer
    !
    SLEFoldTime = DCCalConvertByUnit( SLEFoldTimeValue, SLEFoldTimeUnit, 'sec' ) ! (in)

    if ( FlagSpongeLayer ) then

      ! Sponge layer is not applied for model layers, k = 1, kmax-SLNumLayers.
      !
      do k = 1, kmax-SLNumLayer
        wz_SpongeLayerCoefM(:,k) = 0.0_DP
      end do
      do k = kmax-SLNumLayer+1, kmax
        wz_SpongeLayerCoefM(:,k) = 1.0_DP / ( SLEFoldTime * ( z_Sigma(k) / z_Sigma(kmax) )**SLOrder )
      end do

      if ( .not. FlagSpongeLayerforZonalMean ) then
        ! Sponge layer is not applied for zonal mean component. 
        ! i.e., sponge layer is applied for zonal wave component only. 
        !
        do k = kmax-SLNumLayer+1, kmax
          do l = 1, lmax
            a_DegOrd = nm_l( l )
            if ( a_DegOrd(2) == 0 ) then
              wz_SpongeLayerCoefM(l,k) = 0.0_DP
            end if
          end do
        end do
      end if

      if ( FlagSpongeLayerforHeat ) then

        wz_SpongeLayerCoefH = wz_SpongeLayerCoefM

        ! Sponge layer is not applied for zonal mean component. 
        ! i.e., sponge layer is applied for zonal wave component only. 
        !
        do k = kmax-SLNumLayer+1, kmax
          do l = 1, lmax
            a_DegOrd = nm_l( l )
            if ( a_DegOrd(2) == 0 ) then
              wz_SpongeLayerCoefH(l,k) = 0.0_DP
            end if
          end do
        end do
      else
        wz_SpongeLayerCoefH = 0.0_DP
      end if

    else
      ! Sponge layer is not applied. 
      !
      wz_SpongeLayerCoefM = 0.0_DP
      wz_SpongeLayerCoefH = 0.0_DP
    end if

    ! 運動量の水平拡散係数とスポンジ層減衰係数の和
    ! Damping coefficients by horizonatl diffusion and a sponge layer
    !
    allocate( wz_DisCoefM(lmax, 1:kmax) )
    allocate( wz_DisCoefH(lmax, 1:kmax) )
    allocate( w_DisCoefQ (lmax) )

    do k = 1, kmax
      wz_DisCoefM(:,k) = w_HDifCoefM - wz_SpongeLayerCoefM(:,k)
      wz_DisCoefH(:,k) = w_HDifCoefH - wz_SpongeLayerCoefH(:,k)
    end do
    w_DisCoefQ = w_HDifCoefH


    ! Calculation of divergence damping period
    !
    DivDampPeriod = DCCalConvertByUnit( DivDampPeriodValue, DivDampPeriodUnit, 'sec' )
    if ( DivDampTime0 > TimeN ) then
      call MessageNotify( 'E', module_name, 'DivDampTime0, %f, must be less than or equal to TimeN, %f.', d = (/ DivDampTime0, TimeN /) )
    end if

    ! ヒストリデータ出力のためのへの変数登録
    ! Register of variables for history data output
    !
    call HistoryAutoAddVariable( 'DUDtDyn', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'dynamical tendency of zonal wind', 'm s-2' )
    call HistoryAutoAddVariable( 'DVDtDyn', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'dynamical tendency of meridional wind', 'm s-2' )
    call HistoryAutoAddVariable( 'DTempDtDyn', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'dynamical tendency of temperature', 'K s-1' )
    do n = 1, ncmax
      call HistoryAutoAddVariable( 'D'//trim(a_QMixName(n))//'DtDyn', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'dynamical tendency of water vapor', 's-1' )
    end do
    call HistoryAutoAddVariable( 'DPsDtDyn', (/ AxNameX, AxNameY, AxNameT /), 'dynamical tendency of surface pressure', 'Ps s-1' )
    call HistoryAutoAddVariable( 'OMG', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'vertical velocity in pressure coordinate (omega, DP/Dt)', 'Pa s-1' )

    call HistoryAutoAddVariable( 'SigDot', (/ AxNameX, AxNameY, AxNameR, AxNameT /), 'sigma-vertical velocity', '1 s-1' )
    call HistoryAutoAddVariable( 'DPiDt', (/ AxNameX, AxNameY, AxNameT /), 'Pi (log Ps) tendency)', 'Pa s-1' )

    call HistoryAutoAddVariable( 'Vor', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'vorticity', 's-1' )
    call HistoryAutoAddVariable( 'Div', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'divergence', 's-1' )

    call HistoryAutoAddVariable( 'Mass', (/ AxNameX, AxNameY, AxNameT /), 'mass', 'kg' )
    call HistoryAutoAddVariable( 'KinEngy', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'kinetic energy', 'J m-2' )
    call HistoryAutoAddVariable( 'IntEngy', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'internal energy', 'J m-2' )
    call HistoryAutoAddVariable( 'PotEngy', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'potential energy', 'J m-2' )
    call HistoryAutoAddVariable( 'LatEngy', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'latent energy', 'J m-2' )
    call HistoryAutoAddVariable( 'TotEngy', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'total energy', 'J m-2' )
    call HistoryAutoAddVariable( 'Enstro', (/ AxNameX, AxNameY, AxNameZ, AxNameT /), 'enstrophy', 'kg m-2 s-2' )

    call HistoryAutoAddVariable( 'AdvTestStreamFunc', (/ AxNameX, AxNameY, AxNameR, AxNameT /), 'mass stream function', 'kg s-1' )


    ! Initialization of modules used in this module
    !

    ! 質量の補正
    ! Mass fixer
    !
    call MassFixerInit


    if ( FlagAdvTest ) then
      !
      ! Utility module for advection test
      !
      call AdvTestInit
    end if

    if ( FlagSLTT ) then

!!$    select case ( IDTTMethod )
!!$    case ( IDTTMethodSp )
!!$    case ( IDTTMethodSL )
      ! セミラグランジュ法による物質移流
      ! Semi-Lagrangian method for tracer transport
      !
      call SLTTInit
!!$    case ( IDTTMethodFV )
!!$      ! 物質移流 (有限体積法)
!!$      ! Tracer Transport (finite volume method)
!!$      !
!!$      call TTFVInit
!!$    end select

    end if

    ! 印字 ; Print
    !
    call MessageNotify( 'M', module_name, '----- Initialization Messages -----' )
    call MessageNotify( 'M', module_name, '  TimeIntegScheme             = %c', c1 = trim( TimeIntegScheme ) )
    call MessageNotify( 'M', module_name, '  HDEFoldTime                 = %f [%c]', d = (/ HDEFoldTimeValue /), c1 = trim(HDEFoldTimeUnit) )
    call MessageNotify( 'M', module_name, '  HDOrder                     = %d', i = (/ HDOrder /) )
    call MessageNotify( 'M', module_name, '  VisCoef                     = %f', d = (/ VisCoef /) )
    call MessageNotify( 'M', module_name, '  FlagSpongeLayer             = %b', l = (/ FlagSpongeLayer /) )
    call MessageNotify( 'M', module_name, '  FlagSpongeLayerforZonalMean = %b', l = (/ FlagSpongeLayerforZonalMean /) )
    call MessageNotify( 'M', module_name, '  FlagSpongeLayerforHeat      = %b', l = (/ FlagSpongeLayerforHeat /) )
    call MessageNotify( 'M', module_name, '  SLEFoldTime                 = %f [%c]', d = (/ SLEFoldTimeValue /), c1 = trim(SLEFoldTimeUnit) )
    call MessageNotify( 'M', module_name, '  SLOrder                     = %d', i = (/ SLOrder /) )
    call MessageNotify( 'M', module_name, '  SLNumLayer                  = %d', i = (/ SLNumLayer /) )
    call MessageNotify( 'M', module_name, '  RefTemp                     = %f', d = (/ RefTemp /) )
    call MessageNotify( 'M', module_name, '  FlagDivDamp                 = %b', l = (/ FlagDivDamp /) )
    call MessageNotify( 'M', module_name, '  DivDampPeriod               = %f', d = (/ DivDampPeriod /) )
    call MessageNotify( 'M', module_name, '  DivDampTime0                = %f', d = (/ DivDampTime0 /) )
    call MessageNotify( 'M', module_name, '  FlagTotMassFix              = %b', l = (/ FlagTotMassFix /) )
    call MessageNotify( 'M', module_name, '  FlagMassHorDifCor           = %b', l = (/ FlagMassHorDifCor /) )
    call MessageNotify( 'M', module_name, '  FlagSLTT                    = %b', l = (/ FlagSLTT /) )
!!$    call MessageNotify( 'M', module_name, '  TracerTransportMethod       = %c', c1 = trim( TracerTransportMethod ) )
    call MessageNotify( 'M', module_name, '  FlagAdvTest                 = %b', l = (/ FlagAdvTest /) )
    call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) )

    dynamics_hspl_vas83_inited = .true.

  end subroutine DynamicsHSplVAS83Init
wa_module_inited
Variable :
wa_module_inited = .false. :logical, save, public

Private Instance methods

DelTimeSave
Variable :
DelTimeSave :real(DP), save
: 前回の $ Delta t $ [s]. 陰解法のための係数の再生成の必要性の チェックに使用する.

$ Delta t $ [s] at previous step This is used to check necessity of regeneration of coefficients for implicit method.

DivDampPeriod
Variable :
DivDampPeriod :real(DP), save
: Period for divergence damping application in unit of second
DivDampTime0
Variable :
DivDampTime0 :real(DP), save
: Time for starting divergence damping in unit of second
Subroutine :
wz_DivA(lmax, 1:kmax) :real(DP), intent(inout)
: $ D (t+Delta t) $ . 発散 (スペクトル). Divergence (spectral)

発散の減衰項の付加 (場がつり合っていない計算初期の擾乱の減衰を想定)

Addition of divergence damping term

[Source]

  subroutine DivergenceDamping( wz_DivA )
    !
    ! 発散の減衰項の付加 (場がつり合っていない計算初期の擾乱の減衰を想定)
    !
    ! Addition of divergence damping term
    !

    ! モジュール引用 ; USE statements
    !

    ! 時刻管理
    ! Time control
    !
    use timeset, only: DelTime, TimeN                 ! ステップ $ t $ の時刻. Time of step $ t $. 

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: lmax, imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

    implicit none

    real(DP), intent(inout):: wz_DivA (lmax, 1:kmax)
                              ! $ D (t+\Delta t) $ . 発散 (スペクトル). 
                              ! Divergence (spectral)

    ! 作業変数
    ! Work variables
    !
    real(DP) :: DivDampCoef

    ! 実行文 ; Executable statement
    !

    if ( FlagDivDamp ) then

      DivDampCoef = 1.0_DP / DelTime * ( DivDampPeriod - ( TimeN - DivDampTime0 ) ) / DivDampPeriod
      if ( DivDampCoef < 0.0_DP ) DivDampCoef = 0.0_DP

      wz_DivA = wz_DivA / ( 1.0_DP + 2.0_DP * DelTime * DivDampCoef )

    end if


  end subroutine DivergenceDamping
FlagAdvTest
Variable :
FlagAdvTest :logical , save
: Flag for advection test. When this is true, U, V, T, and Ps are not calculated prognostically.
FlagDivDamp
Variable :
FlagDivDamp :logical , save
: Flag for divergence damping
FlagMassHorDifCor
Variable :
FlagMassHorDifCor :logical , save
: Flag for correction of horizontal diffusion of mass
FlagSLTT
Variable :
FlagSLTT :logical , save
FlagTotMassFix
Variable :
FlagTotMassFix :logical , save
: Flag for fixing total mass
Subroutine :
xyz_Temp(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ T $ . 温度. Temperature
xyz_Phi(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ Phi $ . ジオポテンシャル高度. Getpotential height

格子点データである温度 $ T $ から, 静水圧の式を用いて 格子点データのジオポテンシャル高度 $ Phi $ を求めます.

[Source]

  subroutine HydroGrid( xyz_Temp, xyz_Phi )
    !
    ! 格子点データである温度 $ T $ から, 静水圧の式を用いて
    ! 格子点データのジオポテンシャル高度 $ \Phi $ を求めます. 
    !

    ! モジュール引用 ; USE statements
    !

    use constants, only: CpDry
                              ! $ C_p $ [J kg-1 K-1]. 
                              ! 乾燥大気の定圧比熱. 
                              ! Specific heat of air at constant pressure

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

    ! 宣言文 ; Declaration statements
    !
    implicit none
    real(DP), intent(in):: xyz_Temp (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T $ .     温度. Temperature
    real(DP), intent(out):: xyz_Phi (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \Phi $ .  ジオポテンシャル高度. 
                              ! Getpotential height

    ! 作業変数
    ! Work variables
    !
    integer:: k               ! 鉛直方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in vertical direction

    ! 実行文 ; Executable statement
    !
    xyz_Phi(:,:,1) = CpDry * z_HydroAlpha(1) * xyz_Temp(:,:,1)

    do k = 2, kmax
      xyz_Phi(:,:,k) =   xyz_Phi(:,:,k-1) + CpDry * z_HydroAlpha(k  ) * xyz_Temp(:,:,k  ) + CpDry * z_HydroBeta (k-1) * xyz_Temp(:,:,k-1)
    enddo

  end subroutine HydroGrid
IDTimeIntegScheme
Variable :
IDTimeIntegScheme :integer , save
: ID for used time integration scheme
IDTimeIntegSchemeExplicit
Constant :
IDTimeIntegSchemeExplicit = 0 :integer , parameter
IDTimeIntegSchemeSemiIMplicit
Constant :
IDTimeIntegSchemeSemiIMplicit = 1 :integer , parameter
Subroutine :
xyz_UCosLatN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ U (t) = u (t) cos varphi $ .
xyz_VCosLatN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ V (t) = v (t) cos varphi $ .
xyz_VorN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ zeta (t) $ . 渦度. Vorticity
xyz_DivN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ D (t) $ . 発散. Divergence
xyz_TempN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ T (t) $ . 温度. Temperature
xyz_QVapN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ q (t) $ . 比湿. Specific humidity
xy_GradLambdaPiN(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ DP{pi}{lambda} $
xy_GradMuPiN(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ (1-\mu^2) DP{pi}{mu} $
xyz_PiAdv(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ Dvect{v} cdot nabla pi $ . $ pi $ の移流. Advection of $ pi $
xyz_UAdvN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ U_A (t) $ . 東西運動量移流項. Eastward advection of momentum
xyz_VAdvN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ V_A (t) $ . 南北運動量移流項. Northward advection of momentum
xyz_TempNonLinearN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ hat{H} (t) $ . 温度時間変化項. Temperature tendency
xyz_KinEngyN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ KE (t) $ . 運動エネルギー項. Kinetic energy
xyz_TempUAdvN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ UT’ (t) $ . 温度東西移流項. Eastward advection of temperature
xyz_TempVAdvN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ VT’ (t) $ . 温度南北移流項. Northward advection of temperature
xyr_SigDotN(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(out)
: $ dot{sigma} (t) $ . 鉛直流. Vertical flow
xy_DPiDtNG(0:imax-1, 1:jmax) :real(DP), intent(out)
: $ Z (t) $ . 地表面気圧時間変化非重力波項. Non-gravity wave component of surface pressure tendency

非線形項 (非重力波項) を格子点上で計算します.

Non-linear terms (non gravitational terms) are calculated on grid points

[Source]

  subroutine NonLinearOnGrid( xyz_UCosLatN,     xyz_VCosLatN, xyz_VorN,         xyz_DivN, xyz_TempN,        xyz_QVapN, xy_GradLambdaPiN, xy_GradMuPiN, xyz_PiAdv, xyz_UAdvN,        xyz_VAdvN, xyz_TempNonLinearN, xyz_KinEngyN, xyz_TempUAdvN,    xyz_TempVAdvN, xyr_SigDotN,      xy_DPiDtNG )
    !
    ! 非線形項 (非重力波項) を格子点上で計算します.
    !
    ! Non-linear terms (non gravitational terms) are calculated on
    ! grid points
    !

    ! モジュール引用 ; USE statements
    !

    ! 座標データ設定
    ! Axes data settings
    !
    use axesset, only: r_Sigma, z_DelSigma            ! $ \Delta \sigma $ (整数). 
                              ! $ \Delta \sigma $ (Full)

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: CpDry, EpsV
                              ! $ \epsilon_v $ . 
                              ! 水蒸気分子量比. 
                              ! Molecular weight of water vapor

    ! 文字列操作
    ! Character handling
    !
    use dc_string, only: LChar

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

    ! 組成に関わる配列の設定
    ! Settings of array for atmospheric composition
    !
    use composition, only: ncmax, IndexH2OVap

    implicit none
    real(DP), intent(in):: xyz_UCosLatN   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ U (t) = u (t) \cos \varphi $ .
    real(DP), intent(in):: xyz_VCosLatN   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ V (t) = v (t) \cos \varphi $ .
    real(DP), intent(in):: xyz_VorN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \zeta (t) $ . 渦度. Vorticity
    real(DP), intent(in):: xyz_DivN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ D (t) $ . 発散. Divergence
    real(DP), intent(in):: xyz_TempN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T (t) $ . 温度. Temperature
    real(DP), intent(in):: xyz_QVapN(0:imax-1, 1:jmax, 1:kmax)
                              ! $ q (t) $ . 比湿. Specific humidity
    real(DP), intent(in):: xy_GradLambdaPiN (0:imax-1, 1:jmax)
                              ! $ \DP{\pi}{\lambda} $
    real(DP), intent(in):: xy_GradMuPiN (0:imax-1, 1:jmax)
                              ! $ (1-\mu^2) \DP{\pi}{\mu} $
    real(DP), intent(out):: xyz_PiAdv (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \Dvect{v} \cdot \nabla \pi $ . 
                              ! $ \pi $ の移流. Advection of $ \pi $
    real(DP), intent(out):: xyz_UAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ U_A (t) $ . 東西運動量移流項.
                              ! Eastward advection of momentum
    real(DP), intent(out):: xyz_VAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ V_A (t) $ . 南北運動量移流項. 
                              ! Northward advection of momentum
    real(DP), intent(out):: xyz_TempNonLinearN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \hat{H} (t) $ . 温度時間変化項. 
                              ! Temperature tendency
    real(DP), intent(out):: xyz_KinEngyN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ KE (t) $ . 運動エネルギー項. 
                              ! Kinetic energy
    real(DP), intent(out):: xyz_TempUAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ UT' (t) $ . 温度東西移流項. 
                              ! Eastward advection of temperature
    real(DP), intent(out):: xyz_TempVAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ VT' (t) $ . 温度南北移流項. 
                              ! Northward advection of temperature
    real(DP), intent(out):: xyr_SigDotN (0:imax-1, 1:jmax, 0:kmax)
                              ! $ \dot{\sigma} (t) $ .
                              ! 鉛直流. Vertical flow
    real(DP), intent(out):: xy_DPiDtNG (0:imax-1, 1:jmax)
                              ! $ Z (t) $ . 地表面気圧時間変化非重力波項. 
                              ! Non-gravity wave component of surface pressure tendency


    !-----------------------------------
    !  作業変数
    !  Work variables
    !
    real(DP):: xyz_PiAdvSum (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \sum_k^K(\Dvect{v}\cdot\nabla\pi)\Delta\sigma $ . 
                              ! $ \pi $ 移流の積分値. Integral downward of advection of $ \pi $
    real(DP):: xyz_DivSum (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \sum_k^K D\Delta\sigma $ .
                              ! 発散の積分値. Integral downward of divergence
    real(DP):: xyr_SigDotNonG (0:imax-1, 1:jmax, 0:kmax)
                              ! $ \dot{\sigma} $ .
                              ! 鉛直流 (非重力波). Vertical flow (non gravitational)
    real(DP):: xyz_TempEdd (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T' = T - \bar{T} $ . 
                              ! 温度の擾乱 (整数レベル). Temperature eddy (full level)
    real(DP):: xyr_TempEdd (0:imax-1, 1:jmax, 0:kmax)
                              ! $ \hat{T}' $ . 
                              ! 温度の擾乱 (半整数レベル). Temperature eddy (half level)
    real(DP):: xyz_VirTemp (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T_v $ . 
                              ! 仮温度. Virtual temperature
    real(DP):: xyz_VirTempEdd (0:imax-1, 1:jmax, 1:kmax)
                              ! $ {T_v}' = T_v - \bar{T} $ . 
                              ! 仮温度の擾乱. Virtual temperature eddy

    integer:: k               ! 鉛直方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in vertical direction
    integer:: n               ! 組成方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in dimension of constituents

    ! 実行文 ; Executable statement
    !

    ! $ \pi $ の移流, $ \pi $ 移流と発散の大気上端から下向きの積分
    ! Calculate advection of $ \pi $, integral advection of $ \pi $ and divergence
    !   from the top of the model downward
    !
    do k = 1, kmax
      xyz_PiAdv(:,:,k) = (   xyz_UCosLatN(:,:,k) * xy_GradLambdaPiN + xyz_VCosLatN(:,:,k) * xy_GradMuPiN     ) / ( 1.0_DP - xy_SinLat**2 )
    enddo

    xyz_PiAdvSum(:,:,kmax) = xyz_PiAdv(:,:,kmax) * z_DelSigma(kmax)
    do k = kmax-1, 1, -1
      xyz_PiAdvSum(:,:,k) = xyz_PiAdvSum(:,:,k+1) + xyz_PiAdv(:,:,k) * z_DelSigma(k)
    enddo

    xyz_DivSum(:,:,kmax) = xyz_DivN(:,:,kmax) * z_DelSigma(kmax)
    do k = kmax-1, 1, -1
      xyz_DivSum(:,:,k) = xyz_DivSum(:,:,k+1) + xyz_DivN(:,:,k) * z_DelSigma(k)
    enddo

    ! 地表面気圧時間変化の非重力波項 $ Z $ の計算
    ! Calculate non-gravity wave component of surface pressure tendency $ Z $
    !
    xy_DPiDtNG = - xyz_PiAdvSum(:,:,1)

    ! $ \dot{\sigma} $ の計算
    ! Calculate $ \dot{\sigma} $
    !
    do k = 1, kmax-1
      xyr_SigDotN(:,:,k) = r_Sigma(k) * ( xyz_PiAdvSum(:,:,1) + xyz_DivSum(:,:,1) ) - ( xyz_PiAdvSum(:,:,k+1) + xyz_DivSum(:,:,k+1) )

      xyr_SigDotNonG(:,:,k) = r_Sigma(k) * xyz_PiAdvSum(:,:,1) - xyz_PiAdvSum(:,:,k+1)
    enddo

    ! $ \dot{\sigma} $ の上下境界値
    ! $ \dot{\sigma} $ on upper and lower boundary
    !
    xyr_SigDotN   (:,:,0   ) = 0.
    xyr_SigDotN   (:,:,kmax) = 0.
    xyr_SigDotNonG(:,:,0   ) = 0.
    xyr_SigDotNonG(:,:,kmax) = 0.

    ! 温度の擾乱 (整数レベル), 仮温度, 仮温度の擾乱の計算
    ! Calculate temperature eddy (full level), virtual temperature, 
    !   virtual temperature eddy
    !
    do k = 1, kmax
      xyz_VirTemp(:,:,k) = xyz_TempN(:,:,k) * ( 1.0_DP + ((( 1.0_DP / EpsV ) - 1.0_DP ) * xyz_QVapN(:,:,k)) )

      xyz_TempEdd(:,:,k) = xyz_TempN(:,:,k) - z_RefTemp(k)
      xyz_VirTempEdd(:,:,k) = xyz_VirTemp(:,:,k) - z_RefTemp(k)
    enddo

    ! 温度の擾乱 (半整数レベル) の計算
    ! Calculate temperature eddy (half level)
    !
    xyr_TempEdd(:,:,0   ) = 0.0_DP
    xyr_TempEdd(:,:,kmax) = 0.0_DP
    do k = 1, kmax-1
      xyr_TempEdd(:,:,k) = z_TInpCoefA(k+1) * xyz_TempN(:,:,k+1) + z_TInpCoefB(k  ) * xyz_TempN(:,:,k  ) - r_RefTemp(k)
    enddo

    ! 東西運動量移流項の計算
    ! Calculate advection of eastward momentum
    !
    xyz_UAdvN(:,:,1) = ( xyz_VorN(:,:,1) + xy_Cori ) * xyz_VCosLatN(:,:,1) - 1.0_DP / ( 2.0_DP * z_DelSigma(1) ) * xyr_SigDotN(:,:,1) * ( xyz_UCosLatN(:,:,1) - xyz_UCosLatN(:,:,2) ) - CpDry * z_TInpCoefK(1) * xyz_VirTempEdd(:,:,1) * xy_GradLambdaPiN

    do k = 2, kmax-1
      xyz_UAdvN(:,:,k) = ( xyz_VorN(:,:,k) + xy_Cori ) * xyz_VCosLatN(:,:,k) - 1.0_DP / ( 2.0_DP * z_DelSigma(k) ) * ( xyr_SigDotN(:,:,k-1) * ( xyz_UCosLatN(:,:,k-1) - xyz_UCosLatN(:,:,k) ) + xyr_SigDotN(:,:,k)   * ( xyz_UCosLatN(:,:,k)   - xyz_UCosLatN(:,:,k+1) ) ) - CpDry * z_TInpCoefK(k) * xyz_VirTempEdd(:,:,k) * xy_GradLambdaPiN
    end do

    xyz_UAdvN(:,:,kmax) = ( xyz_VorN(:,:,kmax) + xy_Cori ) * xyz_VCosLatN(:,:,kmax) - 1.0_DP / ( 2.0_DP * z_DelSigma(kmax) ) * xyr_SigDotN(:,:,kmax-1) * ( xyz_UCosLatN(:,:,kmax-1) - xyz_UCosLatN(:,:,kmax) ) - CpDry * z_TInpCoefK(kmax) * xyz_VirTempEdd(:,:,kmax) * xy_GradLambdaPiN


    ! 南北運動量移流項の計算
    ! Calculate advection of northward momentum
    !
    xyz_VAdvN(:,:,1) = - ( xyz_VorN(:,:,1) + xy_Cori ) * xyz_UCosLatN(:,:,1) - 1.0_DP / ( 2.0_DP * z_DelSigma(1) ) * xyr_SigDotN(:,:,1) * ( xyz_VCosLatN(:,:,1) - xyz_VCosLatN(:,:,2) ) - CpDry * z_TInpCoefK(1) * xyz_VirTempEdd(:,:,1) * xy_GradMuPiN

    do k = 2, kmax-1
      xyz_VAdvN(:,:,k) = - ( xyz_VorN(:,:,k) + xy_Cori ) * xyz_UCosLatN(:,:,k) - 1.0_DP / ( 2.0_DP * z_DelSigma(k) ) * ( xyr_SigDotN(:,:,k-1) * ( xyz_VCosLatN(:,:,k-1) - xyz_VCosLatN(:,:,k) ) + xyr_SigDotN(:,:,k)   * ( xyz_VCosLatN(:,:,k)   - xyz_VCosLatN(:,:,k+1) ) ) - CpDry * z_TInpCoefK(k) * xyz_VirTempEdd(:,:,k) * xy_GradMuPiN
    end do

    xyz_VAdvN(:,:,kmax) = - ( xyz_VorN(:,:,kmax) + xy_Cori ) * xyz_UCosLatN(:,:,kmax) - 1.0_DP / ( 2.0_DP * z_DelSigma(kmax) ) * xyr_SigDotN(:,:,kmax-1) * ( xyz_VCosLatN(:,:,kmax-1) - xyz_VCosLatN(:,:,kmax) ) - CpDry * z_TInpCoefK(kmax) * xyz_VirTempEdd(:,:,kmax) * xy_GradMuPiN


    ! 運動エネルギー項 (仮温度補正含む) の計算
    ! Calculate kinematic energy term 
    !   (including virtual temperature correction)
    !
    call HydroGrid( xyz_VirTemp - xyz_TempN, xyz_KinEngyN )             ! (out)

    do k = 1, kmax
      xyz_KinEngyN(:,:,k) = xyz_KinEngyN(:,:,k) + ( xyz_UCosLatN(:,:,k)**2 + xyz_VCosLatN(:,:,k)**2 ) / ( 2.0_DP * ( 1.0_DP - xy_SinLat**2 ) )
    end do


    ! 温度東西移流項, 温度南北移流項の計算
    ! Calculate eastward and northward advection of temperature
    !
    do k = 1, kmax
      xyz_TempUAdvN(:,:,k) =  xyz_UCosLatN(:,:,k) * xyz_TempEdd(:,:,k)
      xyz_TempVAdvN(:,:,k) =  xyz_VCosLatN(:,:,k) * xyz_TempEdd(:,:,k)
    end do

    ! 温度の時間変化項 $ \hat{H} $ の計算
    ! Calculate temperature tendency term $ \hat{H} $
    !
    do k = 1, kmax-1
      xyz_TempNonLinearN(:,:,k) = xyz_TempEdd(:,:,k) * xyz_DivN(:,:,k) - 1.0_DP / z_DelSigma(k) * (   xyr_SigDotN(:,:,k-1) * ( xyr_TempEdd(:,:,k-1) - xyz_TempEdd(:,:,k) ) + xyr_SigDotN(:,:,k) * ( xyz_TempEdd(:,:,k)   - xyr_TempEdd(:,:,k) ) ) - 1.0_DP / z_DelSigma(k) * (   xyr_SigDotNonG(:,:,k-1) * ( r_RefTemp(k-1) - z_RefTemp(k) ) + xyr_SigDotNonG(:,:,k) * ( z_RefTemp(k)   - r_RefTemp(k) ) ) + z_TInpCoefK(k) * xyz_VirTemp(:,:,k) * xyz_PiAdv(:,:,k) - z_HydroAlpha(k) / z_DelSigma(k) * (   xyz_VirTemp(:,:,k)    * xyz_PiAdvSum(:,:,k) + xyz_VirTempEdd(:,:,k) * xyz_DivSum(:,:,k) ) - z_HydroBeta(k) / z_DelSigma(k) * (   xyz_VirTemp(:,:,k)    * xyz_PiAdvSum(:,:,k+1) + xyz_VirTempEdd(:,:,k) * xyz_DivSum(:,:,k+1) )
    enddo

    xyz_TempNonLinearN(:,:,kmax) = xyz_TempEdd(:,:,kmax) * xyz_DivN(:,:,kmax) - 1.0_DP / z_DelSigma(kmax) * (   xyr_SigDotN(:,:,kmax-1) * ( xyr_TempEdd(:,:,kmax-1) - xyz_TempEdd(:,:,kmax) ) + xyr_SigDotN(:,:,kmax) * ( xyz_TempEdd(:,:,kmax)   - xyr_TempEdd(:,:,kmax) ) ) - 1.0_DP / z_DelSigma(kmax) * (   xyr_SigDotNonG(:,:,kmax-1) * ( r_RefTemp(kmax-1) - z_RefTemp(kmax) ) + xyr_SigDotNonG(:,:,kmax) * ( z_RefTemp(kmax)   - r_RefTemp(kmax) ) ) + z_TInpCoefK(kmax) * xyz_VirTemp(:,:,kmax) * xyz_PiAdv(:,:,kmax) - z_HydroAlpha(kmax) / z_DelSigma(kmax) * (   xyz_VirTemp(:,:,kmax)    * xyz_PiAdvSum(:,:,kmax) + xyz_VirTempEdd(:,:,kmax) * xyz_DivSum(:,:,kmax) )


  end subroutine NonLinearOnGrid
Subroutine :
xyz_UCosLatN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ U (t) = u (t) cos varphi $ .
xyz_VCosLatN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ V (t) = v (t) cos varphi $ .
xyz_DivN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ D (t) $ . 発散. Divergence
xyr_SigDotN(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(in )
: $ dot{sigma} (t) $ . 鉛直流. Vertical flow
xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: $ q (t) $ . 比湿. Specific humidity
xyzf_QMixNonLinearN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(out)
: $ R (t) $ . 比湿時間変化項. Specific humidity tendency
xyzf_QMixUAdvN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(out)
: $ Uq (t) $ . 比湿東西移流項. Eastward advection of specific humidity
xyzf_QMixVAdvN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(out)
: $ Vq (t) $ . 比湿南北移流項. Northward advection of specific humidity

非線形項 (非重力波項) を格子点上で計算します.

Non-linear terms (non gravitational terms) are calculated on grid points

[Source]

  subroutine NonLinearOnGridQMix( xyz_UCosLatN, xyz_VCosLatN, xyz_DivN, xyr_SigDotN, xyzf_QMixN, xyzf_QMixNonLinearN, xyzf_QMixUAdvN, xyzf_QMixVAdvN )
    !
    ! 非線形項 (非重力波項) を格子点上で計算します.
    !
    ! Non-linear terms (non gravitational terms) are calculated on
    ! grid points
    !

    ! モジュール引用 ; USE statements
    !

    ! 座標データ設定
    ! Axes data settings
    !
    use axesset, only: r_Sigma, z_DelSigma            ! $ \Delta \sigma $ (整数). 
                              ! $ \Delta \sigma $ (Full)

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

    ! 組成に関わる配列の設定
    ! Settings of array for atmospheric composition
    !
    use composition, only: ncmax
                               ! 成分の数
                               ! Number of composition

    implicit none
    real(DP), intent(in):: xyz_UCosLatN   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ U (t) = u (t) \cos \varphi $ .
    real(DP), intent(in):: xyz_VCosLatN   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ V (t) = v (t) \cos \varphi $ .
    real(DP), intent(in):: xyz_DivN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ D (t) $ . 発散. Divergence
    real(DP), intent(in ):: xyr_SigDotN (0:imax-1, 1:jmax, 0:kmax)
                              ! $ \dot{\sigma} (t) $ .
                              ! 鉛直流. Vertical flow
    real(DP), intent(in):: xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ q (t) $ . 比湿. Specific humidity
    real(DP), intent(out):: xyzf_QMixNonLinearN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ R (t) $ . 比湿時間変化項. 
                              ! Specific humidity tendency
    real(DP), intent(out):: xyzf_QMixUAdvN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ Uq (t) $ . 比湿東西移流項. 
                              ! Eastward advection of specific humidity
    real(DP), intent(out):: xyzf_QMixVAdvN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ Vq (t) $ . 比湿南北移流項. 
                              ! Northward advection of specific humidity

    !-----------------------------------
    !  作業変数
    !  Work variables
    !
    integer:: k               ! 鉛直方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in vertical direction
    integer:: n               ! 組成方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in dimension of constituents

    ! 実行文 ; Executable statement
    !


    ! 比湿東西移流項, 比湿南北移流項の計算
    ! Calculate eastward and northward advection of specific humidity
    !
    do n = 1, ncmax
      do k = 1, kmax
        xyzf_QMixUAdvN(:,:,k,n) =  xyz_UCosLatN(:,:,k) * xyzf_QMixN(:,:,k,n)
        xyzf_QMixVAdvN(:,:,k,n) =  xyz_VCosLatN(:,:,k) * xyzf_QMixN(:,:,k,n)
      end do
    end do

    ! 比湿時間変化項 $ R $ の計算
    ! Calculate specific humidity tendency $ R $
    !
    do n = 1, ncmax
      xyzf_QMixNonLinearN(:,:,1,n) = xyzf_QMixN(:,:,1,n) * xyz_DivN(:,:,1) - 1.0_DP / ( 2.0_DP * z_DelSigma(1) ) * xyr_SigDotN(:,:,1) * ( xyzf_QMixN(:,:,1,n) - xyzf_QMixN(:,:,2,n) )

      do k = 2, kmax - 1
        xyzf_QMixNonLinearN(:,:,k,n) = xyzf_QMixN(:,:,k,n) * xyz_DivN(:,:,k) - 1.0_DP / ( 2.0_DP * z_DelSigma(k) ) * ( xyr_SigDotN(:,:,k-1) * ( xyzf_QMixN(:,:,k-1,n) - xyzf_QMixN(:,:,k  ,n)   ) + xyr_SigDotN(:,:,k) * ( xyzf_QMixN(:,:,k  ,n) - xyzf_QMixN(:,:,k+1,n) ) )
      end do

      xyzf_QMixNonLinearN(:,:,kmax,n) = xyzf_QMixN(:,:,kmax,n) * xyz_DivN(:,:,kmax) - 1.0_DP / ( 2.0_DP * z_DelSigma(kmax) ) * xyr_SigDotN(:,:,kmax-1) * ( xyzf_QMixN(:,:,kmax-1,n) - xyzf_QMixN(:,:,kmax,n) )
    end do


  end subroutine NonLinearOnGridQMix
Subroutine :
xyz_UB(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ u (t-\Delta t) $ . 東西風速. Eastward wind
xyz_VB(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ v (t-\Delta t) $ . 南北風速. Northward wind
xyz_TempB(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ T (t-\Delta t) $ . 温度. Temperature
xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: $ q (t-\Delta t) $ . 比湿. Specific humidity
xy_PsB(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ p_s (t-\Delta t) $ . 地表面気圧. Surface pressure
xyz_UN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ u (t) $ . 東西風速. Eastward wind
xyz_VN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ v (t) $ . 南北風速. Northward wind
xyz_TempN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ T (t) $ . 温度. Temperature
xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: $ q (t) $ . 比湿. Specific humidity
xy_PsN(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ p_s (t) $ . 地表面気圧. Surface pressure
xyz_UA(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ u (t+Delta t) $ . 東西風速. Eastward wind
xyz_VA(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ v (t+Delta t) $ . 南北風速. Northward wind
xyz_TempA(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ T (t+Delta t) $ . 温度. Temperature
xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: $ q (t+Delta t) $ . 比湿. Specific humidity
xy_PsA(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ p_s (t+Delta t) $ . 地表面気圧. Surface pressure
xyz_DUDtPhy(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ left(DP{u}{t}right)^{phy} $ . 外力項 (物理過程) による東西風速変化. Eastward wind tendency by external force terms (physical processes)
xyz_DVDtPhy(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ left(DP{v}{t}right)^{phy} $ . 外力項 (物理過程) による南北風速変化. Northward wind tendency by external force terms (physical processes)
xyz_DTempDtPhy(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ left(DP{T}{t}right)^{phy} $ . 外力項 (物理過程) による温度変化. Temperature tendency by external force terms (physical processes)
xyzf_DQMixDtPhy(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: $ left(DP{q}{t}right)^{phy} $ . 外力項 (物理過程) による比湿変化. Temperature tendency by external force terms (physical processes)
xyr_SigDot(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(in)
: $ dot{sigma} $ . 鉛直流. Vertical flow
xy_DPiDt(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ Z $ . 地表面気圧時間変化項. Surface pressure tendency
xyz_PiAdv(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ Dvect{v} cdot nabla pi $ . $ pi $ の移流. Advection of $ pi $
xyz_OMG(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: OMEGA, DP/Dt

診断量の出力を行います.

Diagnostic variables are output.

[Source]

  subroutine OutputDiagnosedVariables( xyz_UB, xyz_VB, xyz_TempB, xyzf_QMixB, xy_PsB, xyz_UN, xyz_VN, xyz_TempN, xyzf_QMixN, xy_PsN, xyz_UA, xyz_VA, xyz_TempA, xyzf_QMixA, xy_PsA, xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyzf_DQMixDtPhy, xyr_SigDot, xy_DPiDt, xyz_PiAdv, xyz_OMG )
    !
    ! 診断量の出力を行います. 
    !
    ! Diagnostic variables are output. 
    !

    ! モジュール引用 ; USE statements
    !

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: RPlanet, Grav, CpDry, LatentHeat
                              ! $ L $ [J kg-1] . 
                              ! 凝結の潜熱. 
                              ! Latent heat of condensation

    ! 時刻管理
    ! Time control
    !
    use timeset, only: DelTime, TimeN                 ! ステップ $ t $ の時刻. Time of step $ t $. 

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: lmax, imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

    ! 組成に関わる配列の設定
    ! Settings of array for atmospheric composition
    !
    use composition, only: ncmax, a_QMixName, IndexH2OVap

    ! 座標データ設定
    ! Axes data settings
    !
    use axesset, only: z_Sigma, r_Sigma, z_DelSigma            ! $ \Delta \sigma $ (整数). 
                              ! $ \Delta \sigma $ (Full)

    ! SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応) 
    ! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
    !
#ifdef LIB_MPI
#ifdef SJPACK
    use wa_mpi_module_sjpack, only: w_xy => w_xv, xy_GradLon_w => xv_GradLon_w, xy_GradLat_w => xv_GradLat_w, wa_DivLambda_xya => wa_DivLambda_xva, wa_DivMu_xya => wa_DivMu_xva, xya_wa => xva_wa
#else
    use wa_mpi_module, only: w_xy => w_xv, xy_GradLon_w => xv_GradLon_w, xy_GradLat_w => xv_GradLat_w, wa_DivLambda_xya => wa_DivLambda_xva, wa_DivMu_xya => wa_DivMu_xva, xya_wa => xva_wa
#endif
#elif AXISYMMETRY
    use wa_zonal_module, only: w_xy, xy_GradLon_w, xy_GradLat_w, wa_DivLambda_xya, wa_DivMu_xya, xya_wa
#elif SJPACK
    use wa_module_sjpack, only: w_xy, xy_GradLon_w, xy_GradLat_w, wa_DivLambda_xya, wa_DivMu_xya, xya_wa
#elif AXISYMMETRY_SJPACK
    use wa_zonal_module_sjpack, only: w_xy, xy_GradLon_w, xy_GradLat_w, wa_DivLambda_xya, wa_DivMu_xya, xya_wa
#else
    use wa_module, only: w_xy, xy_GradLon_w, xy_GradLat_w, wa_DivLambda_xya, wa_DivMu_xya, xya_wa
#endif

    ! ヒストリデータ出力
    ! History data output
    !
    use gtool_historyauto, only: HistoryAutoPut

    ! 宣言文 ; Declaration statements
    !
    implicit none

    real(DP), intent(in):: xyz_UB    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ u (t-\Delta t) $ .   東西風速. Eastward wind
    real(DP), intent(in):: xyz_VB    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ v (t-\Delta t) $ .   南北風速. Northward wind
    real(DP), intent(in):: xyz_TempB (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T (t-\Delta t) $ .   温度. Temperature
    real(DP), intent(in):: xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ q (t-\Delta t) $ .   比湿. Specific humidity
    real(DP), intent(in):: xy_PsB    (0:imax-1, 1:jmax)
                              ! $ p_s (t-\Delta t) $ . 地表面気圧. Surface pressure

    real(DP), intent(in):: xyz_UN   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ u (t) $ .     東西風速. Eastward wind
    real(DP), intent(in):: xyz_VN   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ v (t) $ .     南北風速. Northward wind
    real(DP), intent(in):: xyz_TempN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T (t) $ .     温度. Temperature
    real(DP), intent(in):: xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ q (t) $ .     比湿. Specific humidity
    real(DP), intent(in):: xy_PsN    (0:imax-1, 1:jmax)
                              ! $ p_s (t) $ .   地表面気圧. Surface pressure

    real(DP), intent(in):: xyz_UA    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ u (t+\Delta t) $ .   東西風速. Eastward wind
    real(DP), intent(in):: xyz_VA    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ v (t+\Delta t) $ .   南北風速. Northward wind
    real(DP), intent(in):: xyz_TempA (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T (t+\Delta t) $ .   温度. Temperature
    real(DP), intent(in):: xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ q (t+\Delta t) $ .   比湿. Specific humidity
    real(DP), intent(in):: xy_PsA    (0:imax-1, 1:jmax)
                              ! $ p_s (t+\Delta t) $ . 地表面気圧. Surface pressure

    real(DP), intent(in):: xyz_DUDtPhy    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \left(\DP{u}{t}\right)^{phy} $ . 
                              ! 外力項 (物理過程) による東西風速変化. 
                              ! Eastward wind tendency by external force terms (physical processes)
    real(DP), intent(in):: xyz_DVDtPhy    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \left(\DP{v}{t}\right)^{phy} $ . 
                              ! 外力項 (物理過程) による南北風速変化. 
                              ! Northward wind tendency by external force terms (physical processes)
    real(DP), intent(in):: xyz_DTempDtPhy (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \left(\DP{T}{t}\right)^{phy} $ . 
                              ! 外力項 (物理過程) による温度変化. 
                              ! Temperature tendency by external force terms (physical processes)
    real(DP), intent(in):: xyzf_DQMixDtPhy(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! $ \left(\DP{q}{t}\right)^{phy} $ . 
                              ! 外力項 (物理過程) による比湿変化. 
                              ! Temperature tendency by external force terms (physical processes)
    real(DP), intent(in):: xyr_SigDot (0:imax-1, 1:jmax, 0:kmax)
                              ! $ \dot{\sigma} $ .
                              ! 鉛直流. Vertical flow
    real(DP), intent(in):: xy_DPiDt (0:imax-1, 1:jmax)
                              ! $ Z $ . 地表面気圧時間変化項. 
                              ! Surface pressure tendency

    real(DP), intent(in):: xyz_PiAdv (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \Dvect{v} \cdot \nabla \pi $ . 
                              ! $ \pi $ の移流. Advection of $ \pi $
    real(DP), intent(out):: xyz_OMG  (0:imax-1, 1:jmax, 1:kmax)
                              !
                              ! OMEGA, DP/Dt

    ! 作業変数
    ! Work variables
    !
    real(DP):: xyz_DUDtDyn    (0:imax-1, 1:jmax, 1:kmax)
    real(DP):: xyz_DVDtDyn    (0:imax-1, 1:jmax, 1:kmax)
    real(DP):: xyz_DTempDtDyn (0:imax-1, 1:jmax, 1:kmax)
    real(DP):: xyzf_DQMixDtDyn(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
    real(DP):: xy_DPsDtDyn    (0:imax-1, 1:jmax)

    real(DP):: xyz_UCosLat   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ U = u \cos \varphi $ .
    real(DP):: xyz_VCosLat   (0:imax-1, 1:jmax, 1:kmax)
                              ! $ V = v \cos \varphi $ .
    real(DP):: xyz_Vor       (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \zeta $ . 渦度. Vorticity
    real(DP):: xyz_Div       (0:imax-1, 1:jmax, 1:kmax)
                              ! $ D $ . 発散. Divergence

    real(DP):: xyr_PiAdvDivSum(0:imax-1, 1:jmax, 0:kmax)

    real(DP):: xy_Mass (0:imax-1, 1:jmax)
                              ! 質量. 
                              ! Mass
    real(DP):: xyz_KinEngy (0:imax-1, 1:jmax, 1:kmax)
                              ! $ KE $ . 運動エネルギー.
                              ! Kinetic energy
    real(DP):: xyz_IntEngy (0:imax-1, 1:jmax, 1:kmax)
                              ! $ IE $ . 内部エネルギー. 
                              ! Internal energy
    real(DP):: xyz_PotEngy (0:imax-1, 1:jmax, 1:kmax)
                              ! $ PE $ . ポテンシャルエネルギー. 
                              ! Potential energy
    real(DP):: xyz_LatEngy (0:imax-1, 1:jmax, 1:kmax)
                              ! $ LE $ . 潜熱エネルギー. 
                              ! Latent heat energy
    real(DP):: xyz_TotEngy (0:imax-1, 1:jmax, 1:kmax)
                              ! $ TE $ . 全エネルギー. 
                              ! Total energy
    real(DP):: xyz_Enstro (0:imax-1, 1:jmax, 1:kmax)
                              ! エンストロフィー. 
                              ! Enstrophy

    integer:: k               ! 鉛直方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in vertical direction
    integer:: n

    ! 実行文 ; Executable statement
    !

    ! Calculate tendencies
    !
    xyz_DUDtDyn     = ( xyz_UA     - xyz_UB     ) / ( 2.0_DP * DelTime ) - xyz_DUDtPhy
    xyz_DVDtDyn     = ( xyz_VA     - xyz_VB     ) / ( 2.0_DP * DelTime ) - xyz_DVDtPhy
    xyz_DTempDtDyn  = ( xyz_TempA  - xyz_TempB  ) / ( 2.0_DP * DelTime ) - xyz_DTempDtPhy
    xyzf_DQMixDtDyn = ( xyzf_QMixA - xyzf_QMixB ) / ( 2.0_DP * DelTime ) - xyzf_DQMixDtPhy
    xy_DPsDtDyn     = ( xy_PsA     - xy_PsB     ) / ( 2.0_DP * DelTime )

    call HistoryAutoPut( TimeN, 'DUDtDyn'   , xyz_DUDtDyn                        )
    call HistoryAutoPut( TimeN, 'DVDtDyn'   , xyz_DVDtDyn                        )
    call HistoryAutoPut( TimeN, 'DTempDtDyn', xyz_DTempDtDyn                     )
    do n = 1, ncmax
      call HistoryAutoPut( TimeN, 'D'//trim(a_QMixName(n))//'DtDyn', xyzf_DQMixDtDyn(:,:,:,n) )
    end do
    call HistoryAutoPut( TimeN, 'DPsDtDyn'  , xy_DPsDtDyn                        )

    ! 鉛直流と地表面気圧時間変化項の出力
    ! Output vertical flow and surface pressure tendency
    !
    call HistoryAutoPut( TimeN, 'SigDot', xyr_SigDot )
    call HistoryAutoPut( TimeN, 'DPiDt',  xy_DPiDt )

    ! 風速から渦度発散の計算
    ! Calculate vorticity and divergence from wind velocity
    !
    do k = 1, kmax
      xyz_UCosLat(:,:,k) = xyz_UN(:,:,k) * xy_CosLat
      xyz_VCosLat(:,:,k) = xyz_VN(:,:,k) * xy_CosLat
    end do
    xyz_Vor = xya_wa(   wa_DivLambda_xya( xyz_VCosLat ) - wa_DivMu_xya(     xyz_UCosLat ) ) / RPlanet
    xyz_Div = xya_wa(   wa_DivLambda_xya( xyz_UCosLat ) + wa_DivMu_xya(     xyz_VCosLat ) ) / RPlanet

    call HistoryAutoPut( TimeN, 'Vor', xyz_Vor )
    call HistoryAutoPut( TimeN, 'Div', xyz_Div )


    ! Calculation of Omega (DPressDt)
    ! It should be recognized that the value of xy_Ps here may have been modified 
    ! by mass fixer.
    !
    !   Integration from top of the model to k's layer upper interface
    xyr_PiAdvDivSum(:,:,kmax) = 0.0_DP
    do k = kmax-1, 0, -1
      xyr_PiAdvDivSum(:,:,k) = xyr_PiAdvDivSum(:,:,k+1) + ( xyz_PiAdv(:,:,k+1) + xyz_Div(:,:,k+1) ) * z_DelSigma(k+1)
    end do
    do k = 1, kmax
      xyz_OMG(:,:,k) = xy_PsN * ( z_Sigma(k) * xyz_PiAdv(:,:,k) - xyr_PiAdvDivSum(:,:,k) - ( xyz_PiAdv(:,:,k) + xyz_Div(:,:,k) ) * ( z_Sigma(k) - r_Sigma(k) ) )
    end do
    call HistoryAutoPut( TimeN, 'OMG', xyz_OMG )

    ! 質量の計算
    ! Calculate mass
    !
    xy_Mass = xy_PsN / Grav

    ! エネルギー, エンストロフィーの計算
    ! Calculate energy and enstrophy
    !
    call HydroGrid( xyz_TempN, xyz_PotEngy )

    do k = 1, kmax
      xyz_KinEngy(:,:,k) = ( xyz_UN(:,:,k) ** 2 + xyz_VN(:,:,k) ** 2 ) / 2.0_DP * xy_Mass

      xyz_IntEngy(:,:,k) = CpDry * xyz_TempN(:,:,k) * xy_Mass

      xyz_PotEngy(:,:,k) = xyz_PotEngy(:,:,k) * xy_Mass

      xyz_LatEngy(:,:,k) = LatentHeat * xyzf_QMixN(:,:,k,IndexH2OVap) * xy_Mass
    end do

    xyz_TotEngy = xyz_KinEngy + xyz_IntEngy + xyz_PotEngy + xyz_LatEngy

    do k = 1, kmax
      xyz_Enstro(:,:,k) = xyz_Vor(:,:,k) ** 2 * xy_Mass
    end do

    call HistoryAutoPut( TimeN, 'Mass',    xy_Mass     )
    call HistoryAutoPut( TimeN, 'KinEngy', xyz_KinEngy )
    call HistoryAutoPut( TimeN, 'IntEngy', xyz_IntEngy )
    call HistoryAutoPut( TimeN, 'PotEngy', xyz_PotEngy )
    call HistoryAutoPut( TimeN, 'LatEngy', xyz_LatEngy )
    call HistoryAutoPut( TimeN, 'TotEngy', xyz_TotEngy )
    call HistoryAutoPut( TimeN, 'Enstro',  xyz_Enstro  )

  end subroutine OutputDiagnosedVariables
Subroutine :

TimeIntegration で使用する係数の設定を行います. 初回および $ Delta t $ が変更された場合以外は, 前回に設定した値をそのまま用います.

Configure coefficients for "TimeIntegration". Setting values that are set last time are used, except when first time and $ Delta t $ are changed.

[Source]

  subroutine SemiImplMatrix
    !
    ! TimeIntegration で使用する係数の設定を行います. 
    ! 初回および $ \Delta t $ が変更された場合以外は, 
    ! 前回に設定した値をそのまま用います. 
    !
    ! Configure coefficients for "TimeIntegration". 
    ! Setting values that are set last time are used, 
    ! except when first time and $ \Delta t $ are changed. 
    !

    ! モジュール引用 ; USE statements
    !

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: RPlanet, CpDry
                              ! $ C_p $ [J kg-1 K-1]. 
                              ! 乾燥大気の定圧比熱. 
                              ! Specific heat of air at constant pressure

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: lmax, imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

    ! 座標データ設定
    ! Axes data settings
    !
    use axesset, only: r_Sigma, z_DelSigma            ! $ \Delta \sigma $ (整数). 
                              ! $ \Delta \sigma $ (Full)


    ! 時刻管理
    ! Time control
    !
    use timeset, only: DelTime ! $ \Delta t $ [s]

    ! SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応) 
    ! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
    !
#ifdef LIB_MPI
#ifdef SJPACK
    use wa_mpi_module_sjpack, only: l_nm
#else
    use wa_mpi_module, only: l_nm
#endif
#elif AXISYMMETRY
    use wa_zonal_module, only: l_nm
#elif SJPACK
    use wa_module_sjpack, only: l_nm
#elif AXISYMMETRY_SJPACK
    use wa_zonal_module_sjpack, only: l_nm
#else
    use wa_module, only: l_nm
#endif

    ! LU 分解法により連立 1 次方程式を解くための関数 (spml 同梱モジュール)
    ! Functions to solve linear simultaneous equation by LU decomposition
    ! (a module included in spml)
    !
    use lumatrix, only: LUDecomp

    ! デバッグ用ユーティリティ
    ! Utilities for debug
    !
    use dc_trace, only: DbgMessage, BeginSub, EndSub, Debug

    ! 宣言文 ; Declaration statements
    !
    implicit none

    ! TimeIntegration 等で使用する係数の設定のための作業変数
    ! Work variable for coefficients for "TimeIntegration", etc.
    !
    integer:: k, l, kk        ! 鉛直方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in vertical direction
    integer:: n

    ! 実行文 ; Executable statement
    !

    ! $ \Delta t $ [s] のチェック・保存
    ! Check and save $ \Delta t $ [s]
    !
    if ( DelTimeSave == DelTime ) return
    DelTimeSave = DelTime

    call DbgMessage( '%c: %c: (DelTime=%f [sec]) coefficients for "TimeIntegration" is generated. ', c1 = module_name, c2 = 'SemiImplMatrix', d = (/ DelTime /) )

    ! TimeIntegration で使用する係数の設定
    ! Configure coefficients for "TimeIntegration"
    !
    if ( .not. allocated( z_siMtxC      ) )  allocate( z_siMtxC     (1:kmax) )
    if ( .not. allocated( z_siMtxG      ) )  allocate( z_siMtxG     (1:kmax) )
    if ( .not. allocated( zz_siMtxH     ) )  allocate( zz_siMtxH    (1:kmax, 1:kmax) )
    if ( .not. allocated( wz_siMtxDi    ) )  allocate( wz_siMtxDi   (lmax,1:kmax) )
    if ( .not. allocated( zz_siMtxDiH   ) )  allocate( zz_siMtxDiH  (1:kmax, 1:kmax) )
    if ( .not. allocated( wzz_siMtxWDiH ) )  allocate( wzz_siMtxWDiH(lmax,1:kmax, 1:kmax) )
    if ( .not. allocated( zz_siMtxGCt   ) )  allocate( zz_siMtxGCt  (1:kmax, 1:kmax) )

    if ( .not. allocated( zz_siMtxW  ) )  allocate( zz_siMtxW  (1:kmax, 1:kmax) )
    if ( .not. allocated( zz_siMtxQ  ) )  allocate( zz_siMtxQ  (1:kmax, 1:kmax) )
    if ( .not. allocated( zz_siMtxS  ) )  allocate( zz_siMtxS  (1:kmax, 1:kmax) )
    if ( .not. allocated( zz_siMtxR  ) )  allocate( zz_siMtxR  (1:kmax, 1:kmax) )


    z_siMtxC = z_DelSigma
    z_siMtxG = CpDry * z_TInpCoefK * z_RefTemp

    do k = 1, kmax
      do l = 1, kmax
        zz_siMtxGCt(k,l) = z_siMtxG(k) * z_siMtxC(l)
      end do
    end do

    zz_siMtxW = 0.
    do k = 1, kmax
      do l = 1, k
        zz_siMtxW(k,l) = CpDry * z_HydroAlpha(l)
      enddo

      do l = 1, k-1
        zz_siMtxW(k,l) = zz_siMtxW(k,l) + CpDry * z_HydroBeta(l)
      enddo
    enddo

    zz_siMtxS = 0.
    do k = 1, kmax
      do l = 1, kmax
        zz_siMtxS(k,l) = r_Sigma(k-1) * z_DelSigma(l)
      enddo
      do l = k, kmax
        zz_siMtxS(k,l) = zz_siMtxS(k,l) - z_DelSigma(l)
      enddo
    enddo

    zz_siMtxQ = 0.
    do k = 1, kmax
      zz_siMtxQ(k,k) = ( r_RefTemp(k-1) - z_RefTemp(k) ) / z_DelSigma(k)
    enddo
    do k = 1, kmax-1
      zz_siMtxQ(k,k+1) = ( z_RefTemp(k) - r_RefTemp(k) ) / z_DelSigma(k)
    enddo

    zz_siMtxR = 0.
    do k = 1, kmax
      do l = k, kmax
        zz_siMtxR(k,l) = - z_HydroAlpha(k) / z_DelSigma(k) * z_DelSigma(l) * z_RefTemp(k)
      enddo
      do l = k + 1, kmax
        zz_siMtxR(k,l) = zz_siMtxR(k,l) - z_HydroBeta(k) / z_DelSigma(k) * z_DelSigma(l) * z_RefTemp(k)
      enddo
    enddo

    zz_siMtxH = matmul(zz_siMtxQ, zz_siMtxS) - zz_siMtxR

    ! Check of coefficients for horizontal diffusion and sponge layer
    ! A threshold value used below is arbitrary.
    !
    do n = 1, lmax
      do k = 1, kmax
        if ( abs( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefM(n,k) ) < 1.0e-10_DP ) then
          call MessageNotify( 'E', module_name, 'Dissipation coefficient for momentum is inappropriate.' )
        end if
        if ( abs( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefH(n,k) ) < 1.0e-10_DP ) then
          call MessageNotify( 'E', module_name, 'Dissipation coefficient for heat is inappropriate.' )
        end if
      end do
      if ( abs( 1.0_DP - 2.0_DP * DelTime * w_DisCoefQ(n) ) < 1.0e-10_DP ) then
        call MessageNotify( 'E', module_name, 'Dissipation coefficient for composition is inappropriate.' )
      end if
    end do

    wz_siMtxDi = 1.0_DP / ( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefH )

    do n = 1, lmax
      ! NOTE: wz_siMtiDi is a diagonal matrix.
      !
      do k = 1, kmax
        do kk = 1, kmax
          zz_siMtxDiH(k,kk) = wz_siMtxDi(n,k) * zz_siMtxH(k,kk)
        end do
      end do
      zz_siMtxDiH = matmul( zz_siMtxW, zz_siMtxDiH )
      do k = 1, kmax
        do kk = 1, kmax
          wzz_siMtxWDiH(n,k,kk) = zz_siMtxDiH(k,kk)
        end do
      end do
    end do

    if ( .not. allocated(wzz_siMtxLU) )  allocate( wzz_siMtxLU(lmax, 1:kmax, 1:kmax) )
    if ( .not. allocated(wz_siMtxPiv) )  allocate( wz_siMtxPiv(lmax, 1:kmax) )

    ! 行列 $ \underline{M} $ の計算
    ! Calculate matrix $ \underline{M} $
    !
    do k = 1, kmax
      do kk = 1, kmax
        wzz_siMtxLU ( :,k,kk ) = - DelTime**2 * w_LaplaEigVal(:) * ( wzz_siMtxWDiH(:,k,kk) + zz_siMtxGCt(k,kk) )
        if ( k == kk ) then
          wzz_siMtxLU ( :,k,kk ) = wzz_siMtxLU ( :,k,kk ) + ( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefM (:,k) )
        endif
      end do
    end do


    ! LU 行列計算
    ! LU matrix calculation
    !
    call LUDecomp( wzz_siMtxLU, wz_siMtxPiv )   ! (out)


  end subroutine SemiImplMatrix
Subroutine :
w_SurfGeoPot(lmax) :real(DP), intent(in)
: $ Phi_s $ . 地表ジオポテンシャル. Surface geo-potential
wz_DVorDtNG(lmax, 1:kmax) :real(DP), intent(in)
: $ left( DD{zeta}{t} (t) right)^{NG}$ . 渦度変化の非重力波成分 (スペクトル). Non-gravity wave component of vorticity tendency (spectral)
wz_DDivDtNG(lmax, 1:kmax) :real(DP), intent(in)
: $ DD{D}{t} (t) $ . 発散変化の非重力波成分 (スペクトル). Divergence tendency (spectral)
wz_DTempDtNG(lmax, 1:kmax) :real(DP), intent(in)
: $ left( DD{T}{t} (t) right)^{NG}$ . 温度変化の非重力波成分 (スペクトル). Temperature tendency (spectral)
w_DPiDtNG(lmax) :real(DP), intent(in)
: $ left( DD{p_s}{t} (t) right) $ . 地表面気圧変化の非重力波項 (スペクトル). Non-gravity wave component of surface pressure tendency (spectral)
wz_DivN(lmax, 1:kmax) :real(DP), intent(in)
: $ D (t) $ . 発散 (スペクトル). Divergence (spectral)
wz_TempN(lmax, 1:kmax) :real(DP), intent(in)
: $ T (t) $ . 温度 (スペクトル). Temperature (spectral)
w_PiN(lmax) :real(DP), intent(in)
: $ pi = ln p_s (t) $ . 地表面気圧 (スペクトル).
wz_VorB(lmax, 1:kmax) :real(DP), intent(in)
: $ zeta (t-\Delta t) $ . 渦度 (スペクトル). Vorticity (spectral)
wz_DivB(lmax, 1:kmax) :real(DP), intent(in)
: $ D (t-\Delta t) $ . 発散 (スペクトル). Divergence (spectral)
wz_TempB(lmax, 1:kmax) :real(DP), intent(in)
: $ T (t-\Delta t) $ . 温度 (スペクトル). Temperature (spectral)
w_PiB(lmax) :real(DP), intent(in)
: $ pi = ln p_s (t-\Delta t) $ . 地表面気圧 (スペクトル). Surface pressure (spectral)
wz_VorA(lmax, 1:kmax) :real(DP), intent(out)
: $ zeta (t+Delta t) $ . 渦度 (スペクトル). Vorticity (spectral)
wz_DivA(lmax, 1:kmax) :real(DP), intent(out)
: $ D (t+Delta t) $ . 発散 (スペクトル). Divergence (spectral)
wz_TempA(lmax, 1:kmax) :real(DP), intent(out)
: $ T (t+Delta t) $ . 温度 (スペクトル). Temperature (spectral)
w_PiA(lmax) :real(DP), intent(out)
: $ pi = ln p_s (t+Delta t) $ . 地表面気圧 (スペクトル). Surface pressure (spectral)

時間積分を行い, 時刻 $ t $ の物理量の時間変化と $ t-\Delta t$ の物理量から 時刻 $ t+Delta t $ の物理量を計算します.

時間積分法にはリープフロッグスキームを用いています. ただし拡散項には後方差分を用いています. デフォルトでは, $ Delta t $ を大きくとるために, 重力波項に セミインプリシット法を適用しています. NAMELIST#dynamics_hspl_vas83_nml の TimeIntegScheme を変更することで, 重力波項をエクスプリシット法によって解くことも可能です.

With time integration, physical values at $ t+Delta t $ is calculated from tendency at $ t $ and physical values at $ t-\Delta t $ .

Leap-frog scheme is used as time integration scheme. And backward scheme is used for diffusion terms. As default setting, semi-implicit scheme is applied to gravitational terms in order to enlarge the value of $ Delta t $ . Explicit scheme can be applied to gravitational terms by changing "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml".

[Source]

  subroutine TimeIntegration( w_SurfGeoPot, wz_DVorDtNG, wz_DDivDtNG, wz_DTempDtNG, w_DPiDtNG, wz_DivN,                  wz_TempN,     w_PiN, wz_VorB,     wz_DivB,     wz_TempB,     w_PiB, wz_VorA,     wz_DivA,     wz_TempA,     w_PiA )
    !
    ! 時間積分を行い, 
    ! 時刻 $ t $ の物理量の時間変化と $ t-\Delta t$ の物理量から
    ! 時刻 $ t+\Delta t $ の物理量を計算します.
    !
    ! 時間積分法にはリープフロッグスキームを用いています. 
    ! ただし拡散項には後方差分を用いています. 
    ! デフォルトでは, $ \Delta t $ を大きくとるために, 重力波項に 
    ! セミインプリシット法を適用しています.
    ! NAMELIST#dynamics_hspl_vas83_nml の TimeIntegScheme を変更することで, 
    ! 重力波項をエクスプリシット法によって解くことも可能です.
    ! 
    ! With time integration, physical values at $ t+\Delta t $ is calculated
    ! from tendency at $ t $ and physical values at $ t-\Delta t $ .
    !
    ! Leap-frog scheme is used as time integration scheme. 
    ! And backward scheme is used for diffusion terms.
    ! As default setting, semi-implicit scheme is applied to gravitational terms 
    ! in order to enlarge the value of $ \Delta t $ . 
    ! Explicit scheme can be applied to gravitational terms by changing
    ! "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml". 
    !

    ! モジュール引用 ; USE statements
    !

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: RPlanet, CpDry
                              ! $ C_p $ [J kg-1 K-1]. 
                              ! 乾燥大気の定圧比熱. 
                              ! Specific heat of air at constant pressure

    ! 時刻管理
    ! Time control
    !
    use timeset, only: DelTime ! $ \Delta t $ [s]

    ! 座標データ設定
    ! Axes data settings
    !
    use axesset, only: z_DelSigma            ! $ \Delta \sigma $ (整数). 
                              ! $ \Delta \sigma $ (Full)

    ! LU 分解法により連立 1 次方程式を解くための関数 (spml 同梱モジュール)
    ! Functions to solve linear simultaneous equation by LU decomposition
    ! (a module included in spml)
    !
    use lumatrix, only: LUSolve

    ! 文字列操作
    ! Character handling
    !
    use dc_string, only: LChar

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: lmax, imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

    ! 組成に関わる配列の設定
    ! Settings of array for atmospheric composition
    !
    use composition, only: ncmax, IndexH2OVap


    implicit none
    real(DP), intent(in):: w_SurfGeoPot (lmax)
                              ! $ \Phi_s $ . 地表ジオポテンシャル. 
                              ! Surface geo-potential

    real(DP), intent(in):: wz_DVorDtNG (lmax, 1:kmax)
                              ! $ \left( \DD{\zeta}{t} (t) \right)^{NG}$ . 渦度変化の非重力波成分 (スペクトル). 
                              ! Non-gravity wave component of vorticity tendency (spectral)
    real(DP), intent(in):: wz_DDivDtNG (lmax, 1:kmax)
                              ! $ \DD{D}{t} (t) $ . 発散変化の非重力波成分 (スペクトル). 
                              ! Divergence tendency (spectral)
    real(DP), intent(in):: wz_DTempDtNG(lmax, 1:kmax)
                              ! $ \left( \DD{T}{t} (t) \right)^{NG}$ . 温度変化の非重力波成分 (スペクトル). 
                              ! Temperature tendency (spectral)
    real(DP), intent(in):: w_DPiDtNG(lmax)
                              ! $ \left( \DD{p_s}{t} (t) \right) $ . 地表面気圧変化の非重力波項 (スペクトル). 
                              ! Non-gravity wave component of surface pressure tendency (spectral)
    real(DP), intent(in):: wz_DivN (lmax, 1:kmax)
                              ! $ D (t) $ . 発散 (スペクトル). 
                              ! Divergence (spectral)
    real(DP), intent(in):: wz_TempN (lmax, 1:kmax)
                              ! $ T (t) $ . 温度 (スペクトル). 
                              ! Temperature (spectral)
    real(DP), intent(in):: w_PiN   (lmax)
                              ! $ \pi = \ln p_s (t) $ . 地表面気圧 (スペクトル). 

    real(DP), intent(in):: wz_VorB (lmax, 1:kmax)
                              ! $ \zeta (t-\Delta t) $ . 渦度 (スペクトル). 
                              ! Vorticity (spectral)
    real(DP), intent(in):: wz_DivB (lmax, 1:kmax)
                              ! $ D (t-\Delta t) $ . 発散 (スペクトル). 
                              ! Divergence (spectral)
    real(DP), intent(in):: wz_TempB (lmax, 1:kmax)
                              ! $ T (t-\Delta t) $ . 温度 (スペクトル). 
                              ! Temperature (spectral)
    real(DP), intent(in):: w_PiB (lmax)
                              ! $ \pi = \ln p_s (t-\Delta t) $ . 地表面気圧 (スペクトル). 
                              ! Surface pressure (spectral)

    real(DP), intent(out):: wz_VorA (lmax, 1:kmax)
                              ! $ \zeta (t+\Delta t) $ . 渦度 (スペクトル). 
                              ! Vorticity (spectral)
    real(DP), intent(out):: wz_DivA (lmax, 1:kmax)
                              ! $ D (t+\Delta t) $ . 発散 (スペクトル). 
                              ! Divergence (spectral)
    real(DP), intent(out):: wz_TempA (lmax, 1:kmax)
                              ! $ T (t+\Delta t) $ . 温度 (スペクトル). 
                              ! Temperature (spectral)
    real(DP), intent(out):: w_PiA (lmax)
                              ! $ \pi = \ln p_s (t+\Delta t) $ . 地表面気圧 (スペクトル). 
                              ! Surface pressure (spectral)

    ! 作業変数
    ! Work variables
    !
    real(DP):: wz_HDiv (lmax, 1:kmax)
    real(DP):: w_CtDiv (lmax)
    real(DP):: wz_WT   (lmax, 1:kmax)

    real(DP):: wz_siTemp (lmax, 1:kmax)
                              ! 温度 (セミインプリシット法のための作業変数). 
                              ! Temperature (work variable for semi-implicit scheme)
    real(DP):: w_siPi (lmax)
                              ! $ \pi $ (セミインプリシット法のための作業変数). 
                              ! $ \pi $ (work variable for semi-implicit scheme)
    real(DP):: wz_siPhi (lmax, 1:kmax)
                              ! $ \Phi = \underline{W} ( 1 - 2 \Delta t \underline{D_H} ) \overline{ \Dvect{T} }^{t}$ .
                              ! (セミインプリシット法のための作業変数). 
                              ! (Work variable for semi-implicit scheme)
    real(DP):: wz_siVectF (lmax, 1:kmax)
                              ! $ \Dvect{f} $ . 
                              ! 発散項に関するセミインプリシット方程式の右辺. 
                              ! Right-hand side of a semi-implicit equation of a divergence term. 

    real(DP):: wz_siDivAvrTime (lmax, 1:kmax)
                              ! $ \overline{\Dvect{D}}^{t} $ . 
                              ! 時間平均の $ \Dvect{D} $ (セミインプリシット法のための作業変数). 
                              ! Time average $ \Dvect{D} $ (a work variable for semi-implicit scheme)

    integer:: k, kk           ! 鉛直方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in vertical direction
    integer:: n               ! 組成方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in dimension of constituents

    ! 実行文 ; Executable statement
    !

    ! 時間積分. 拡散は後方差分
    ! Time integration. Backward difference is applied to diffusion
    !

    ! 渦度 ; Vorticity
    !
    wz_VorA = ( 1.0_DP / ( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefM ) ) * ( wz_VorB + 2.0_DP * DelTime * wz_DVorDtNG )

    ! 発散 ; Divergence
    !
    select case ( IDTimeIntegScheme )
    case ( IDTimeIntegSchemeSemiImplicit )

      ! セミインプリシット法で用いる行列の計算に使われる項の計算
      ! Calculate terms used in making a matrix for semi-implicit method
      !
      wz_siTemp = ( 1.0_DP - DelTime * wz_DisCoefH ) * wz_TempB + DelTime * wz_DTempDtNG

      !     NOTE:
      !     The matrix wz_siMtxDi = ( 1 - 2 \Delta t D_H )^{-1} is diagonal matrix.
      !
      wz_siPhi (:,1) = CpDry * z_HydroAlpha(1) * wz_siMtxDi(:,1) * wz_siTemp(:,1)
      do k = 2, kmax
        wz_siPhi (:,k) = wz_siPhi(:,k-1) + CpDry * z_HydroAlpha(k  ) * wz_siMtxDi(:,k  ) * wz_siTemp(:,k  ) + CpDry * z_HydroBeta (k-1) * wz_siMtxDi(:,k-1) * wz_siTemp(:,k-1)
      end do

      w_siPi = w_PiB + DelTime * w_DPiDtNG


      ! 発散方程式の右辺の計算
      ! Calculate right side of divergence equation
      !
      do k = 1, kmax
        wz_siVectF(:,k) = ( 1.0_DP - DelTime * wz_DisCoefM (:,k) ) * wz_DivB(:,k) + DelTime * wz_DDivDtNG(:,k) - DelTime * w_LaplaEigVal(:) * ( w_SurfGeoPot + wz_siPhi(:,k) + z_siMtxG(k) * w_siPi )
      end do


      ! 時間平均の $ \Dvect{D} $ を LU 行列で解く
      ! Solve time average $ \Dvect{D} $ with LU matrix
      !
      wz_siDivAvrTime = LUSolve( wzz_siMtxLU, wz_siMtxPiv, wz_siVectF )


      wz_DivA = 2. * wz_siDivAvrTime - wz_DivB

    case ( IDTimeIntegSchemeExplicit )

      wz_WT = 0.0_DP
      do k = 1, kmax
        do kk = 1, kmax
          wz_WT(:,k) = wz_WT(:,k) + zz_siMtxW(k,kk) * wz_TempN(:,kk)
        end do
      end do

      do k = 1, kmax
        wz_DivA(:,k) = ( 1.0_DP / ( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefM(:,k) ) ) * ( wz_DivB(:,k) + 2.0_DP * DelTime * (   wz_DDivDtNG(:,k) - w_LaplaEigVal(:) * (   w_SurfGeoPot + wz_WT(:,k) + z_siMtxG(k) * w_PiN ) ) )
      end do

    end select

    ! 温度 ; Temperature
    !
    select case ( IDTimeIntegScheme )
    case ( IDTimeIntegSchemeSemiImplicit )
      wz_HDiv = 0.
      do k = 1, kmax
        do kk = 1, kmax
          wz_HDiv(:,k) = wz_HDiv(:,k) + zz_siMtxH(k,kk) * wz_siDivAvrTime(:,kk)
        end do
      end do
    case ( IDTimeIntegSchemeExplicit )
      wz_HDiv = 0.
      do k = 1, kmax
        do kk = 1, kmax
          wz_HDiv(:,k) = wz_HDiv(:,k) + zz_siMtxH(k,kk) * wz_DivN(:,kk)
        end do
      end do
    end select

    wz_TempA = ( 1.0_DP  / ( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefH ) ) * ( wz_TempB + 2.0_DP * DelTime * ( wz_DTempDtNG - wz_HDiv ) )


    ! 地表面気圧 ; Surface pressure
    !
    select case ( IDTimeIntegScheme )
    case ( IDTimeIntegSchemeSemiImplicit )
      w_CtDiv = 0.0_DP
      do k = 1, kmax
        w_CtDiv = w_CtDiv + z_siMtxC(k) * wz_siDivAvrTime(:,k)
      end do
    case ( IDTimeIntegSchemeExplicit )
      w_CtDiv = 0.0_DP
      do k = 1, kmax
        w_CtDiv = w_CtDiv + z_siMtxC(k) * wz_DivN(:,k)
      end do
    end select

    w_PiA  = w_PiB + 2.0_DP * DelTime * ( w_DPiDtNG - w_CtDiv )


  end subroutine TimeIntegration
Subroutine :
wzf_DQMixDtN(lmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: $ DD{q}{t} (t) $ . 比湿変化 (スペクトル). Specific humidity tendency (spectral)
wzf_QMixB(lmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: $ q (t-\Delta t) $ . 比湿 (スペクトル). Specific humidity (spectral)
wzf_QMixA(lmax, 1:kmax, 1:ncmax) :real(DP), intent(out)
: $ q (t+Delta t) $ . 比湿 (スペクトル). Specific humidity (spectral)

時間積分を行い, 時刻 $ t $ の物理量の時間変化と $ t-\Delta t$ の物理量から 時刻 $ t+Delta t $ の物理量を計算します.

時間積分法にはリープフロッグスキームを用いています. ただし拡散項には後方差分を用いています. デフォルトでは, $ Delta t $ を大きくとるために, 重力波項に セミインプリシット法を適用しています. NAMELIST#dynamics_hspl_vas83_nml の TimeIntegScheme を変更することで, 重力波項をエクスプリシット法によって解くことも可能です.

With time integration, physical values at $ t+Delta t $ is calculated from tendency at $ t $ and physical values at $ t-\Delta t $ .

Leap-frog scheme is used as time integration scheme. And backward scheme is used for diffusion terms. As default setting, semi-implicit scheme is applied to gravitational terms in order to enlarge the value of $ Delta t $ . Explicit scheme can be applied to gravitational terms by changing "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml".

[Source]

  subroutine TimeIntegrationQMix( wzf_DQMixDtN, wzf_QMixB, wzf_QMixA )
    !
    ! 時間積分を行い, 
    ! 時刻 $ t $ の物理量の時間変化と $ t-\Delta t$ の物理量から
    ! 時刻 $ t+\Delta t $ の物理量を計算します.
    !
    ! 時間積分法にはリープフロッグスキームを用いています. 
    ! ただし拡散項には後方差分を用いています. 
    ! デフォルトでは, $ \Delta t $ を大きくとるために, 重力波項に 
    ! セミインプリシット法を適用しています.
    ! NAMELIST#dynamics_hspl_vas83_nml の TimeIntegScheme を変更することで, 
    ! 重力波項をエクスプリシット法によって解くことも可能です.
    ! 
    ! With time integration, physical values at $ t+\Delta t $ is calculated
    ! from tendency at $ t $ and physical values at $ t-\Delta t $ .
    !
    ! Leap-frog scheme is used as time integration scheme. 
    ! And backward scheme is used for diffusion terms.
    ! As default setting, semi-implicit scheme is applied to gravitational terms 
    ! in order to enlarge the value of $ \Delta t $ . 
    ! Explicit scheme can be applied to gravitational terms by changing
    ! "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml". 
    !

    ! モジュール引用 ; USE statements
    !

    ! 時刻管理
    ! Time control
    !
    use timeset, only: DelTime ! $ \Delta t $ [s]

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: lmax, imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

    ! 組成に関わる配列の設定
    ! Settings of array for atmospheric composition
    !
    use composition, only: ncmax
                               ! 成分の数
                               ! Number of composition


    implicit none

    real(DP), intent(in):: wzf_DQMixDtN(lmax, 1:kmax, 1:ncmax)
                              ! $ \DD{q}{t} (t) $ . 比湿変化 (スペクトル). 
                              ! Specific humidity tendency (spectral)
    real(DP), intent(in):: wzf_QMixB(lmax, 1:kmax, 1:ncmax)
                              ! $ q (t-\Delta t) $ . 比湿 (スペクトル). 
                              ! Specific humidity (spectral)

    real(DP), intent(out):: wzf_QMixA(lmax, 1:kmax, 1:ncmax)
                              ! $ q (t+\Delta t) $ . 比湿 (スペクトル). 
                              ! Specific humidity (spectral)

    ! 作業変数
    ! Work variables
    !
    integer:: k               ! 鉛直方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in vertical direction
    integer:: n               ! 組成方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in dimension of constituents

    ! 実行文 ; Executable statement
    !

    ! 時間積分. 拡散は後方差分
    ! Time integration. Backward difference is applied to diffusion
    !

    ! 比湿 ; Specific humidity
    !
    do n = 1, ncmax
      do k = 1, kmax
        wzf_QMixA(:,k,n) = ( 1.0_DP / ( 1.0_DP - 2.0_DP * DelTime * w_DisCoefQ ) ) * ( wzf_QMixB(:,k,n) + 2.0_DP * DelTime * wzf_DQMixDtN(:,k,n) )
      end do
    end do


  end subroutine TimeIntegrationQMix
dynamics_hspl_vas83_inited
Variable :
dynamics_hspl_vas83_inited = .false. :logical, save
: 初期設定フラグ. Initialization flag
module_name
Constant :
module_name = ‘dynamics_hspl_vas83 :character(*), parameter
: モジュールの名称. Module name
r_RefTemp
Variable :
r_RefTemp(:) :real(DP), allocatable
: $ hat{overline{T}} $ . 基準温度 (半整数レベル). Reference temperature on half sigma level
version
Constant :
version = ’$Name: $’ // ’$Id: dynamics_hspl_vas83.F90,v 1.83 2014/05/07 09:39:17 murashin Exp $’ :character(*), parameter
: モジュールのバージョン Module version
w_DisCoefQ
Variable :
w_DisCoefQ(:) :real(DP), save, allocatable
: 成分の水平拡散係数 Damping coefficient for composition by horizonatl diffusion
w_LaplaEigVal
Variable :
w_LaplaEigVal(:) :real(DP), save, allocatable
: $ nabla^2_{sigma} = -n (n+1) $ . ラプラシアンの固有値. Eigen value of laplacian operator.
wz_DisCoefH
Variable :
wz_DisCoefH(:,:) :real(DP), save, allocatable
: 熱の水平拡散係数とスポンジ層減衰係数の和 Damping coefficient for heat by horizonatl diffusion and a sponge layer
wz_DisCoefM
Variable :
wz_DisCoefM(:,:) :real(DP), save, allocatable
: 運動量の水平拡散係数とスポンジ層減衰係数の和 Damping coefficient for momentum by horizonatl diffusion and a sponge layer
wz_siMtxDi
Variable :
wz_siMtxDi(:,:) :real(DP), save, allocatable
: $ D = ( 1 - 2 Delta t D_h )^{-1}$ NOTE: This is a diagonal matrix.
wz_siMtxPiv
Variable :
wz_siMtxPiv(:,:) :integer , save, allocatable
: セミインプリシット行列のピボット. Pivot for semi-implicit matrix
wzz_siMtxLU
Variable :
wzz_siMtxLU(:,:,:) :real(DP), save, allocatable
: セミインプリシット行列の LU 分解. LU decomposition for semi-implicit matrix
wzz_siMtxWDiH
Variable :
wzz_siMtxWDiH(:,:,:) :real(DP), save, allocatable
: $ W Di h $
xy_Cori
Variable :
xy_Cori(:,:) :real(DP), allocatable
: $ f\equiv 2\Omega\sin\varphi $ . コリオリパラメータ. Coriolis parameter
xy_CosLat
Variable :
xy_CosLat(:,:) :real(DP), allocatable
: $ cos varphi $ .
xy_SinLat
Variable :
xy_SinLat(:,:) :real(DP), allocatable
: $ sin varphi $ .
z_HydroAlpha
Variable :
z_HydroAlpha(:) :real(DP), allocatable
: $ alpha $ . 静水圧の式の係数. Coefficient in hydrostatic equation.
z_HydroBeta
Variable :
z_HydroBeta(:) :real(DP), allocatable
: $ beta $ . 静水圧の式の係数. Coefficient in hydrostatic equation.
z_RefTemp
Variable :
z_RefTemp(:) :real(DP), allocatable
: $ overline{T} $ . 基準温度 (整数レベル). Reference temperature on full sigma level
z_TInpCoefA
Variable :
z_TInpCoefA(:) :real(DP), allocatable
: $ a $ . 温度鉛直補間の係数. Coefficient for vertical interpolation of temperature
z_TInpCoefB
Variable :
z_TInpCoefB(:) :real(DP), allocatable
: $ b $ . 温度鉛直補間の係数. Coefficient for vertical interpolation of temperature
z_TInpCoefK
Variable :
z_TInpCoefK(:) :real(DP), allocatable
: $ kappa $ . 温度鉛直補間の係数. Coefficient for vertical interpolation of temperature
z_siMtxC
Variable :
z_siMtxC(:) :real(DP), save, allocatable
: $ C = sigma $
z_siMtxG
Variable :
z_siMtxG(:) :real(DP), save, allocatable
: $ G = C_p kappa T $
zz_siMtxDiH
Variable :
zz_siMtxDiH(:,:) :real(DP), save, allocatable
: $ Di h $
zz_siMtxGCt
Variable :
zz_siMtxGCt(:,:) :real(DP), save, allocatable
: $ G C^{T} $ ( $ C = Delta sigma $ )
zz_siMtxH
Variable :
zz_siMtxH(:,:) :real(DP), save, allocatable
: $ h = QS - R $
zz_siMtxQ
Variable :
zz_siMtxQ(:,:) :real(DP), save, allocatable
: $ Q = DD{T}{sigma} $
zz_siMtxR
Variable :
zz_siMtxR(:,:) :real(DP), save, allocatable
: $ R $ . 発散と掛け合わせることで気圧変化の効果となる. Product R and divergence become effort of surface pressure tendency. $ RD = kappa T (DD{pi}{t} + Dinv{sigma}DD{sigma}{t}) $ .
!$ integer, save, allocatable:nmo (:,:,:)

!$ ! スペクトルの添字順番. !$ ! Spectral subscript expression

!$ real(DP), save, allocatable:wzz_siMtxM (:,:,:)

!$ ! 行列 $ underline{M} $. !$ ! Matrix $ underline{M} $

!$ integer, save, allocatable:z_siMtxPivWork(:)

!$ ! 行列のピボット作成の作業変数. !$ ! Work variable for pivot

zz_siMtxS
Variable :
zz_siMtxS(:,:) :real(DP), save, allocatable
: $ S = DD{sigma}{t} $
zz_siMtxW
Variable :
zz_siMtxW(:,:) :real(DP), save, allocatable
: $ W $ . 発散の式での線形重力波項の効果による温度補正の係数. Coefficient for correction of temperature by effort of linear gravitational terms