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_nmlTimeIntegScheme を変更することで, 重力波項をエクスプリシット法によって解くことも可能です.

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 for extension of $ Delta t $ . Explicit scheme can be applied to gravitational terms by changing "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml".

Procedures List

Dynamics :力学計算
———— :————
Dynamics :Calculate dynamics

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

gridset dc_types constants timeset wa_mpi_module wa_module gtool_historyauto dc_string axesset lumatrix w_module namelist_util gtool_history dc_iounit dc_date_types dc_date dc_message 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
xyz_QVapB(0:imax-1, 1:jmax, 1:kmax) :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
xyz_QVapN(0:imax-1, 1:jmax, 1:kmax) :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)
xyz_DQVapDtPhy(0:imax-1, 1:jmax, 1:kmax) :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
xyz_QVapA(0:imax-1, 1:jmax, 1:kmax) :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

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

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

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

Calculate dynamical processes. Eastward wind (xyz_UB, xyz_UN), northward wind (xyz_VB, xyz_VN), temperature (xyz_TempB, xyz_TempN), specific humidity (xyz_QVapB, xyz_QVapN), surface pressure (xyz_PsB, xyz_PsN) at $ t-\Delta t $ and $ t $, and surface height (xy_SurfHeight) are input, and eastward wind (xyz_UA), northward wind (xyz_VA), temperature (xyz_TempA), specific humidity (xyz_QVapA), surface pressure (xyz_PsA) are returned.

In order to add tendencies of vorticity, divergence, temperature and specific humidity by another physical process to tendencies of this dynamical process and to calculate the values at next time step, give these tendencies to "xyz_DUDtPhy", "xyz_DVDtPhy", "xyz_DTempDtPhy", "xyz_DQVapDtPhy"

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 Dynamics( xyz_UB,      xyz_VB,      xyz_TempB,      xyz_QVapB,   xy_PsB, xyz_UN,      xyz_VN,      xyz_TempN,      xyz_QVapN,   xy_PsN, xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyz_DQVapDtPhy, xy_SurfHeight, xyz_UA,      xyz_VA,      xyz_TempA,      xyz_QVapA,   xy_PsA )
    !
    ! 力学過程の演算を行い, 与えられた $ t-\Delta t $ および $ t $ の
    ! 東西風速 (xyz_UB, xyz_UN), 南北風速 (xyz_VB, xyz_VN), 
    ! 温度 (xyz_TempB, xyz_TempN), 比湿 (xyz_QVapB, xyz_QVapN), 
    ! 地表面気圧 (xyz_PsB, xyz_PsN), および地表面高度 (xy_SurfHeight) から, 
    ! $ t+\Delta t $ の
    ! 東西風速 (xyz_UA), 南北風速 (xyz_VA), 温度 (xyz_TempA), 
    ! 比湿 (xyz_QVapA), 地表面気圧 (xyz_PsA) を返します. 
    !
    ! 別の物理プロセスによる渦度, 発散, 温度, 比湿の変化を, 
    ! 力学過程の変化に足して次のステップを計算する場合には, 
    ! それらの変化を xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyz_DQVapDtPhy
    ! に与えてください. 
    !
    ! 時間積分法にはリープフロッグスキームを用いています.
    ! デフォルトでは, $ \Delta t $ を大きくとるために, 重力波項に 
    ! セミインプリシット法を適用しています.
    ! NAMELIST#dynamics_hspl_vas83_nml の TimeIntegScheme を変更することで, 
    ! 重力波項をエクスプリシット法によって解くことも可能です.
    !
    ! Calculate dynamical processes. 
    ! Eastward wind (xyz_UB, xyz_UN), 
    ! northward wind (xyz_VB, xyz_VN), 
    ! temperature (xyz_TempB, xyz_TempN), 
    ! specific humidity (xyz_QVapB, xyz_QVapN), 
    ! surface pressure (xyz_PsB, xyz_PsN) at $ t-\Delta t $ and $ t $, 
    ! and surface height (xy_SurfHeight) are input, 
    ! and 
    ! eastward wind (xyz_UA), northward wind (xyz_VA), 
    ! temperature (xyz_TempA), 
    ! specific humidity (xyz_QVapA), surface pressure (xyz_PsA) 
    ! are returned.
    !
    ! In order to add tendencies of vorticity, divergence, 
    ! temperature and specific humidity by another physical process to 
    ! tendencies of this dynamical process 
    ! and to calculate the values at next time step, 
    ! give these tendencies to 
    ! "xyz_DUDtPhy", "xyz_DVDtPhy", "xyz_DTempDtPhy", "xyz_DQVapDtPhy"
    !
    ! 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: imax, jmax, kmax    ! 鉛直層数. 
                               ! Number of vertical level

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

#ifdef LIB_MPI
    ! MPI 版 SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応) 
    ! SPMODEL library MPI version, problems on sphere are solved with spherical harmonics (multi layer is supported)
    !
    use wa_mpi_module, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, xy_GradLon_w   => xv_GradLon_w, xy_GradLat_w   => xv_GradLat_w, xya_GradLon_wa => xva_GradLon_wa, xya_GradLat_wa => xva_GradLat_wa, wa_Div_xya_xya => wa_Div_xva_xva, w_xy   => w_xv,   xy_w   => xv_w, wa_xya => wa_xva, xya_wa => xva_wa
#else
    ! SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応) 
    ! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
    !
    use wa_module, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, xy_GradLon_w, xy_GradLat_w, xya_GradLon_wa, xya_GradLat_wa, wa_Div_xya_xya, 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. 

    ! 宣言文 ; 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):: xyz_QVapB (0:imax-1, 1:jmax, 1:kmax)
                              ! $ 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):: xyz_QVapN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ 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):: xyz_DQVapDtPhy (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \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):: xyz_QVapA (0:imax-1, 1:jmax, 1:kmax)
                              ! $ q (t+\Delta t) $ .   比湿. Specific humidity
    real(DP), intent(out):: xy_PsA    (0:imax-1, 1:jmax)
                              ! $ p_s (t+\Delta t) $ . 地表面気圧. Surface pressure

    ! 作業変数
    ! Work variables
    !
    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):: w_PiN ((nmax+1)**2)
                              ! $ \pi $ スペクトル
    real(DP):: xy_GradLonPiN (0:imax-1, 1:jmax)
                              ! $ \frac{1}{\cos \varphi} \DP{\pi}{\lambda} $
    real(DP):: xy_GradLatPiN (0:imax-1, 1:jmax)
                              ! $ \DP{\pi}{\varphi} $

    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_DTempDtN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ H (t) $ . 温度時間変化項. 
                              ! Temperature tendency
    real(DP):: xyz_DQVapDtN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ 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_DPiDtN (0:imax-1, 1:jmax)
                              ! $ Z $ . 地表面気圧時間変化項. 
                              ! Surface pressure tendency
    real(DP):: xyz_QVapUAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ uq (t) $ . 比湿東西移流項. 
                              ! Eastward advection of specific humidity
    real(DP):: xyz_QVapVAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ vq (t) $ . 比湿南北移流項. 
                              ! Northward advection of specific humidity
    real(DP):: wz_DVorDtN ((nmax+1)**2, 1:kmax)
                              ! $ \DD{\zeta}{t} (t) $ . 渦度変化 (スペクトル). 
                              ! Vorticity tendency (spectral)
    real(DP):: wz_DDivDtN ((nmax+1)**2, 1:kmax)
                              ! $ \DD{D}{t} (t) $ . 発散変化 (スペクトル). 
                              ! Divergence tendency (spectral)
    real(DP):: wz_DTempDtN ((nmax+1)**2, 1:kmax)
                              ! $ \DD{T}{t} (t) $ . 温度変化 (スペクトル). 
                              ! Temperature tendency (spectral)
    real(DP):: wz_DQVapDtN ((nmax+1)**2, 1:kmax)
                              ! $ \DD{q}{t} (t) $ . 比湿変化 (スペクトル). 
                              ! Specific humidity tendency (spectral)
    real(DP):: w_DPiDtN ((nmax+1)**2)
                              ! $ \DD{p_s}{t} (t) $ . 地表面気圧変化 (スペクトル). 
                              ! Surface pressure tendency (spectral)
    real(DP):: wz_VorB ((nmax+1)**2, 1:kmax)
                              ! $ \zeta (t-\Delta t) $ . 渦度 (スペクトル). 
                              ! Vorticity (spectral)
    real(DP):: wz_DivB ((nmax+1)**2, 1:kmax)
                              ! $ D (t-\Delta t) $ . 発散 (スペクトル). 
                              ! Divergence (spectral)
    real(DP):: wz_TempB ((nmax+1)**2, 1:kmax)
                              ! $ T (t-\Delta t) $ . 温度 (スペクトル). 
                              ! Temperature (spectral)
    real(DP):: w_PiB ((nmax+1)**2)
                              ! $ \pi = \ln p_s (t-\Delta t) $ . 地表面気圧 (スペクトル). 
                              ! Surface pressure (spectral)
    real(DP):: wz_QVapB ((nmax+1)**2, 1:kmax)
                              ! $ q (t-\Delta t) $ . 比湿 (スペクトル). 
                              ! Specific humidity (spectral)
    real(DP):: wz_VorA ((nmax+1)**2, 1:kmax)
                              ! $ \zeta (t+\Delta t) $ . 渦度 (スペクトル). 
                              ! Vorticity (spectral)
    real(DP):: wz_DivA ((nmax+1)**2, 1:kmax)
                              ! $ D (t+\Delta t) $ . 発散 (スペクトル). 
                              ! Divergence (spectral)
    real(DP):: wz_TempA ((nmax+1)**2, 1:kmax)
                              ! $ T (t+\Delta t) $ . 温度 (スペクトル). 
                              ! Temperature (spectral)
    real(DP):: w_PiA ((nmax+1)**2)
                              ! $ \pi = \ln p_s (t+\Delta t) $ . 地表面気圧 (スペクトル). 
                              ! Surface pressure (spectral)
    real(DP):: wz_QVapA ((nmax+1)**2, 1:kmax)
                              ! $ q (t+\Delta t) $ . 比湿 (スペクトル). 
                              ! Specific humidity (spectral)

    real(DP):: wz_Psi ((nmax+1)**2, 1:kmax)
                              ! $ \psi $ . 流線関数. Streamline function 
    real(DP):: wz_Chi ((nmax+1)**2, 1:kmax)
                              ! $ \chi $ . ポテンシャル. Potential

    real(DP):: wz_VorDiffA ((nmax+1)**2, 1:kmax)
                              ! $ \mathscr{D}(\zeta) $ . 
                              ! 運動量水平拡散による渦度変化 (スペクトル). 
                              ! Vorticity tendency by 
                              ! horizontal momentum diffusion (spectral)
    real(DP):: wz_DivDiffA ((nmax+1)**2, 1:kmax)
                              ! $ \mathscr{D}(D) $ . 
                              ! 運動量水平拡散による発散変化 (スペクトル). 
                              ! Divergence tendency by 
                              ! horizontal momentum diffusion (spectral)
    real(DP):: wz_PsiDiff ((nmax+1)**2, 1:kmax)
                              ! 運動量水平拡散による流線関数 $ \psi $ 変化
                              ! Streamline function tendency by 
                              ! horizontal momentum diffusion
    real(DP):: wz_ChiDiff ((nmax+1)**2, 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 ((nmax+1)**2)
                              ! $ \Phi_s $ . 地表ジオポテンシャル. 
                              ! Surface geo-potential

    ! エクスプリシット法のための作業変数
    ! Work variables for explicit scheme
    !
    real(DP):: xyz_exWTGPi (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \underline{W} \Dvect{T} + \Dvect{G} \pi $ .
    real(DP):: xyz_exWT (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \underline{W} \Dvect{T} $ .
    real(DP):: xyz_exGPi (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \Dvect{G} \pi $ .

    real(DP):: xyz_exHDiv (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \underline{h} \Dvect{D} $ .


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

    ! 実行文 ; Executable statement
    !

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

    ! 初期化
    ! Initialization
    !
    if ( .not. dynamics_hspl_vas83_inited ) call DynamicsInit

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

    ! 地表面気圧の空間変化項の計算
    ! Calculate spatial surface pressure tendency
    !
    xy_PiN = log( xy_PsN )
    w_PiN = w_xy( xy_PiN )
    xy_GradLonPiN = xy_GradLon_w( w_PiN ) / RPlanet
    xy_GradLatPiN = xy_GradLat_w( w_PiN ) / RPlanet

    ! 風速から渦度発散の計算
    ! Calculate vorticity and divergence from wind velocity
    ! 
    xyz_VorN = xya_wa( wa_Div_xya_xya( xyz_VN , - xyz_UN ) ) / RPlanet
    xyz_DivN = xya_wa( wa_Div_xya_xya( xyz_UN ,   xyz_VN ) ) / RPlanet

    ! 格子点上での非線形力学項の計算
    ! Calculate non-linear dynamical terms on grid points
    !
    call NonLinearOnGrid( xyz_UN,        xyz_VN, xyz_VorN,      xyz_DivN, xyz_TempN,     xyz_QVapN, xy_GradLonPiN, xy_GradLatPiN, xyz_UAdvN,     xyz_VAdvN, xyz_DTempDtN,  xyz_DQVapDtN, xyz_KinEngyN, xyz_TempUAdvN, xyz_TempVAdvN, xyr_SigDotN,   xy_DPiDtN, xyz_QVapUAdvN, xyz_QVapVAdvN )  ! (out)


    ! スペクトル時間変化項の計算
    ! Calculate spectral tendency terms
    !

    ! 渦度の時間変化 (スペクトル) の計算
    ! Calculate vorticity tendency (spectral)
    !
    wz_DVorDtN = wa_Div_xya_xya( xyz_VAdvN, - xyz_UAdvN ) / RPlanet

    ! 発散の時間変化 (スペクトル) の計算
    ! Calculate divergence tendency (spectral)
    !
    wz_DDivDtN = wa_Div_xya_xya( xyz_UAdvN, xyz_VAdvN ) / RPlanet - wa_Lapla_wa( wa_xya(xyz_KinEngyN) ) / RPlanet**2

    ! 温度の時間変化 (スペクトル) の計算
    ! Calculate temperature tendency (spectral)
    !
    wz_DTempDtN = - wa_Div_xya_xya( xyz_TempUAdvN, xyz_TempVAdvN ) / RPlanet + wa_xya( xyz_DTempDtN )

    ! 地表面気圧の時間変化 (スペクトル) の計算
    ! Calculate surface pressure tendency (spectral)
    !
    w_DPiDtN = w_xy( xy_DPiDtN )

    ! 比湿の時間変化 (スペクトル) の計算
    ! Calculate specific humidity tendency (spectral)
    !
    wz_DQVapDtN = - wa_Div_xya_xya( xyz_QVapUAdvN, xyz_QVapVAdvN ) / RPlanet + wa_xya( xyz_DQVapDtN )

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

    ! エクスプリシット法を用いる際の計算
    ! Calculate for explicit scheme
    !
    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('explicit')

      ! $ \underline{W} \Dvect{T} + \Dvect{G} \pi $ の格子点値の計算
      ! Calculate $ \underline{W} \Dvect{T} + \Dvect{G} \pi $ on grid
      !

      ! $ \underline{W} \Dvect{T} $ の計算
      ! Calculate $ \underline{W} \Dvect{T} $
      !
      call HydroGrid( xyz_TempN, xyz_exWT )   ! (out)

      ! $ \Dvect{G} \pi $ の計算
      ! Calculate $ \Dvect{G} \pi $
      !
      do k = 1, kmax
        xyz_exGPi(:,:,k) = CpDry * z_TInpCoefK(k) * z_RefTemp(k) * xy_PiN
      enddo

      ! $ \underline{W} \Dvect{T} + \Dvect{G} \pi $ の計算
      ! Calculate $ \underline{W} \Dvect{T} + \Dvect{G} \pi $
      !
      xyz_exWTGPi = xyz_exWT + xyz_exGPi


      ! $ \underline{h} \Dvect{D} $ の格子点値の計算
      ! Calculate $ \underline{h} \Dvect{D} $ on grid
      !
      xyz_exHDiv = 0.

      do k = 1, kmax
        do kk = 1, kmax
          xyz_exHDiv (:,:,k) = xyz_exHDiv (:,:,k) + zz_siMtxH (k,kk) * xyz_DivN (:,:,kk)
        enddo
      enddo


      ! 発散, 温度 (スペクトル) の時間変化項の修正
      ! Modify divergence and temperature tencency (spectral)
      !

      ! 発散の時間変化 (スペクトル) の計算
      ! Calculate divergence tendency (spectral)
      !
      wz_DDivDtN = wz_DDivDtN - wa_Lapla_wa( wa_xya( xyz_exWTGPi ) ) / RPlanet**2

      do k = 1, kmax
        wz_DDivDtN(:,k) = wz_DDivDtN(:,k) - w_Lapla_w( w_SurfGeoPot ) / RPlanet**2
      end do

      ! 温度の時間変化 (スペクトル) の計算
      ! Calculate temperature tendency (spectral)
      !
      wz_DTempDtN = wz_DTempDtN - wa_xya( xyz_exHDiv )

    end select


    ! 格子点値をスペクトル値へ ( $ t-\Delta t$ )
    ! Exchange grid values to spectral values ( $ t-\Delta t$ )
    !
    wz_VorB  = wa_Div_xya_xya( xyz_VB, - xyz_UB ) / RPlanet
    wz_DivB  = wa_Div_xya_xya( xyz_UB,   xyz_VB ) / RPlanet
    wz_TempB = wa_xya( xyz_TempB )
    wz_QVapB = wa_xya( xyz_QVapB )
    w_PiB    = w_xy( log( xy_PsB ) )


    ! 外部プロセスによる時間変化格子データをスペクトルデータへ変換
    ! Convert tendency data on grid by external processes into spectral data
    !
    wz_DVorDtN  =   wz_DVorDtN + wa_Div_xya_xya( xyz_DVDtPhy, - xyz_DUDtPhy ) / RPlanet
    wz_DDivDtN  =   wz_DDivDtN + wa_Div_xya_xya( xyz_DUDtPhy,   xyz_DVDtPhy ) / RPlanet
    wz_DTempDtN = wz_DTempDtN + wa_xya( xyz_DTempDtPhy )
    wz_DQVapDtN = wz_DQVapDtN + wa_xya( xyz_DQVapDtPhy )

    ! 時間積分
    ! Time integration
    !
    call TimeIntegration( wz_DVorDtN, wz_DDivDtN, wz_DTempDtN, wz_DQVapDtN, w_DPiDtN, wz_VorB,    wz_DivB,    wz_TempB,    wz_QVapB,    w_PiB, w_SurfGeoPot, wz_VorA,    wz_DivA,    wz_TempA,    wz_QVapA,    w_PiA )     ! (out)

    ! スペクトル値を格子点値へ ( $ 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_UA = (   xya_GradLon_wa( wz_Chi ) - xya_GradLat_wa( wz_Psi )    ) / RPlanet

    xyz_VA = (   xya_GradLon_wa( wz_Psi ) + xya_GradLat_wa( wz_Chi )    ) / RPlanet

    xyz_TempA = xya_wa( wz_TempA )
    xyz_QVapA = xya_wa( wz_QVapA )
    xy_PsA = exp( xy_w( w_PiA ) )

    ! 拡散による補正
    ! Correction by diffusion
    !

    ! 運動量水平拡散による渦度発散の時間変化
    ! Vorticity and divergence tendency by horizontal diffusion of momentum
    !
    wz_VorDiffA = wz_VorA * wz_HDifCoefM
    wz_DivDiffA = wz_DivA * wz_HDifCoefM


    ! 運動量水平拡散による摩擦熱補正
    ! Frictional thermal correction by horizontal momentum diffusion
    !
    wz_PsiDiff = wa_LaplaInv_wa( wz_VorDiffA ) * RPlanet**2
    wz_ChiDiff = wa_LaplaInv_wa( wz_DivDiffA ) * RPlanet**2

    xyz_UDiff = (   xya_GradLon_wa( wz_ChiDiff ) - xya_GradLat_wa( wz_PsiDiff )    ) / RPlanet

    xyz_VDiff = (   xya_GradLon_wa( wz_PsiDiff ) + xya_GradLat_wa( wz_ChiDiff )    ) / RPlanet

    xyz_TempA = xyz_TempA - (   xyz_UA * xyz_UDiff + xyz_VA * xyz_VDiff   ) / CpDry * 2. * DelTime

    ! ヒストリデータ出力
    ! History data output
    !
    call DiagOutput( xyz_UN, xyz_VN, xyz_TempN, xyz_QVapN, xy_PsN, xyr_SigDotN, xy_DPiDtN )                        ! (in)


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

  end subroutine Dynamics
dynamics_hspl_vas83_inited
Variable :
dynamics_hspl_vas83_inited = .false. :logical, save, public
: 初期設定フラグ. Initialization flag

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.

Subroutine :
xyz_U(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ u $ . 東西風速. Eastward wind
xyz_V(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ v $ . 南北風速. Northward wind
xyz_Temp(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ T $ . 温度. Temperature
xyz_QVap(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in)
: $ q $ . 比湿. Specific humidity
xy_Ps(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ p_s $ . 地表面気圧. Surface pressure
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

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

Diagnostic variables are output.

[Source]

  subroutine DiagOutput( xyz_U, xyz_V, xyz_Temp, xyz_QVap, xy_Ps, xyr_SigDot, xy_DPiDt )
    !
    ! 診断量の出力を行います. 
    !
    ! 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: TimeN  ! ステップ $ t $ の時刻. Time of step $ t $. 

#ifdef LIB_MPI
    ! MPI 版 SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応) 
    ! SPMODEL library MPI version, problems on sphere are solved with spherical harmonics (multi layer is supported)
    !
    use wa_mpi_module, only: wa_Div_xya_xya => wa_Div_xva_xva, xya_wa => xva_wa
#else
    ! SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応) 
    ! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
    !
    use wa_module, only: wa_Div_xya_xya, xya_wa
#endif

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

    ! 宣言文 ; Declaration statements
    !
    implicit none
    real(DP), intent(in):: xyz_U    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ u $ .   東西風速. Eastward wind
    real(DP), intent(in):: xyz_V    (0:imax-1, 1:jmax, 1:kmax)
                              ! $ v $ .   南北風速. Northward wind
    real(DP), intent(in):: xyz_Temp (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T $ .   温度. Temperature
    real(DP), intent(in):: xyz_QVap (0:imax-1, 1:jmax, 1:kmax)
                              ! $ q $ .   比湿. Specific humidity
    real(DP), intent(in):: xy_Ps    (0:imax-1, 1:jmax)
                              ! $ p_s $ . 地表面気圧. Surface pressure
    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

    ! 作業変数
    ! Work variables
    !
    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):: 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

    ! 実行文 ; Executable statement
    !

    ! 鉛直流と地表面気圧時間変化項の出力
    ! 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
    ! 
    xyz_Vor = xya_wa( wa_Div_xya_xya( xyz_V , - xyz_U ) / RPlanet )
    xyz_Div = xya_wa( wa_Div_xya_xya( xyz_U ,   xyz_V ) / RPlanet )

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

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

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

    do k = 1, kmax
      xyz_KinEngy(:,:,k) = ( xyz_U(:,:,k) ** 2 + xyz_V(:,:,k) ** 2 ) / 2. * xy_Mass

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

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

      xyz_LatEngy(:,:,k) = LatentHeat * xyz_QVap(:,:,k) * 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 DiagOutput
Subroutine :

This procedure input/output NAMELIST#dynamics_hspl_vas83_nml .

[Source]

  subroutine DynamicsInit

    ! モジュール引用 ; 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

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


#ifdef LIB_MPI
    ! MPI 版 SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応) 
    ! SPMODEL library MPI version, problems on sphere are solved with spherical harmonics (multi layer is supported)
    !
    use wa_mpi_module, only: xy_Lat => xv_Lat, w_xy => w_xv
#else

    ! SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応) 
    ! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
    !
    use wa_module, only: xy_Lat, w_xy
#endif
    use w_module, only: rn

    ! 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_date_types, only: DC_DIFFTIME
                              ! 日時の差を表現するデータ型. 
                              ! Data type for difference about date and time
    use dc_date, only: DCDiffTimeCreate, EvalSec

    ! メッセージ出力
    ! Message output
    !
    use dc_message, only: MessageNotify

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

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

    ! 宣言文 ; Declaration statements
    !
    implicit none


    ! 基準温度の設定のための作業変数
    ! Work variable for reference temperature
    !
    real(DP):: RefTemp
                              ! $ \overline{T} $ . 基準温度. 
                              ! Reference temperature

    ! 水平拡散係数の設定のための作業変数
    ! Work variable for coefficient of horizontal diffusion
    !
    integer:: VisOrder        ! 超粘性の次数.  Order of hyper-viscosity
    real(DP):: VisCoef        ! 超粘性係数. Hyper-viscosity coefficient
    real(DP):: EFoldTimeValue ! 最大波数に対する e-folding time. 
                              ! 負の値を与えると, 水平拡散係数をゼロにします. 
                              ! 
                              ! E-folding time for maximum wavenumber. 
                              ! If negative value is given, 
                              ! coefficients of horizontal diffusion become zero. 
    character(TOKEN):: EFoldTimeUnit
                              ! 最大波数に対する e-folding time の単位. 
                              ! Unit of e-folding time for maximum wavenumber
    real(DP):: EFoldTime      ! 最大波数に対する e-folding time [単位: 秒]. 
                              ! E-folding time for maximum wavenumber [Unit: sec]
    type(DC_DIFFTIME):: dcdiff_efold

    ! 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:: 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, VisOrder, EFoldTimeValue, EFoldTimeUnit, RefTemp
          !
          ! デフォルト値については初期化手続 "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
    call InitCheck

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

    VisOrder = 8
    EFoldTimeValue = 8640.
    EFoldTimeUnit  = 'sec'

    RefTemp = 300.


    ! 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')
    case ('explicit')
    case default
      call MessageNotify( 'E', module_name, '"TimeIntegScheme" must be "Semi-Implicit" or "Explicit".' )
    end select


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

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

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


    ! 静水圧の式の係数 $ \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.

      z_HydroBeta(k) = 1. - ( 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. - (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. ) ** ( -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


    ! 水平拡散係数の設定
    ! Configure coefficient of horizontal diffusion
    !
    allocate( wz_LaplaCoef  ((nmax+1)**2, 1:kmax) )
    allocate( wz_HDifCoefM ((nmax+1)**2, 1:kmax) )
    allocate( wz_HDifCoefT  ((nmax+1)**2, 1:kmax) )

    do k = 1, kmax
      wz_LaplaCoef(:,k) = rn(:,1) / RPlanet**2
    enddo

    ! 粘性係数の計算 (最大波数の e-folding time が EFoldTime となるように)
    ! Calculate viscosity coefficient
    !
    call DCDiffTimeCreate( dcdiff_efold, EFoldTimeValue, EFoldTimeUnit )       ! (in)
    EFoldTime = EvalSec( dcdiff_efold )

    if ( EFoldTimeValue > 0. ) then
      VisCoef = ( (nmax*(nmax+1)) / RPlanet**2 )**(-VisOrder / 2) / EFoldTime
    else
      VisCoef = 0.
    end if

    wz_HDifCoefT = - VisCoef * ( ( - wz_LaplaCoef )**( VisOrder / 2 ) )

    wz_HDifCoefM = wz_HDifCoefT - VisCoef * ( - ( 2. / RPlanet**2 )**( VisOrder / 2 ) )

    ! ヒストリデータ出力のためのへの変数登録
    ! Register of variables for history data output
    !
    call HistoryAutoAddVariable( 'SigDot', (/ 'lon ', 'lat ', 'sigm', 'time' /), 'sigma-vertical velocity', '1 s-1' )
    call HistoryAutoAddVariable( 'DPiDt', (/ 'lon ', 'lat ', 'time' /), 'Pi (log Ps) tendency)', 'Pa s-1' )

    call HistoryAutoAddVariable( 'Vor', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'vorticity', 's-1' )
    call HistoryAutoAddVariable( 'Div', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'divergence', 's-1' )

    call HistoryAutoAddVariable( 'Mass', (/ 'lon ', 'lat ', 'time' /), 'mass', 'kg' )
    call HistoryAutoAddVariable( 'KinEngy', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'kinetic energy', 'J' )
    call HistoryAutoAddVariable( 'IntEngy', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'internal energy', 'J' )
    call HistoryAutoAddVariable( 'PotEngy', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'potential energy', 'J' )
    call HistoryAutoAddVariable( 'LatEngy', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'latent energy', 'J' )
    call HistoryAutoAddVariable( 'TotEngy', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'total energy', 'J' )
    call HistoryAutoAddVariable( 'Enstro', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'enstrophy', 'kg' )

    ! 印字 ; Print
    !
    call MessageNotify( 'M', module_name, '----- Initialization Messages -----' )
    call MessageNotify( 'M', module_name, '  TimeIntegScheme  = %c', c1 = trim( TimeIntegScheme ) )
    call MessageNotify( 'M', module_name, '  EFoldTime = %f [%c]', d = (/ EFoldTimeValue /), c1 = trim(EFoldTimeUnit) )
    call MessageNotify( 'M', module_name, '  VisOrder  = %d', i = (/ VisOrder /) )
    call MessageNotify( 'M', module_name, '  VisCoef   = %f', d = (/ VisCoef /) )
    call MessageNotify( 'M', module_name, '  RefTemp   = %f', d = (/ RefTemp /) )
    call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) )

    dynamics_hspl_vas83_inited = .true.
  end subroutine DynamicsInit
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

    ! 宣言文 ; 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
Subroutine :

依存モジュールの初期化チェック

Check initialization of dependency modules

[Source]

  subroutine InitCheck
    !
    ! 依存モジュールの初期化チェック
    !
    ! Check initialization of dependency modules

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

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

    ! 格子点設定
    ! Grid points settings
    !
    use gridset, only: gridset_inited

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: constants_inited

    ! 座標データ設定
    ! Axes data settings
    !
    use axesset, only: axesset_inited

    ! 時刻管理
    ! Time control
    !
    use timeset, only: timeset_inited

    ! メッセージ出力
    ! Message output
    !
    use dc_message, only: MessageNotify

    ! 実行文 ; Executable statement
    !

    if ( .not. namelist_util_inited ) call MessageNotify( 'E', module_name, '"namelist_util" module is not initialized.' )

    if ( .not. gridset_inited ) call MessageNotify( 'E', module_name, '"gridset" module is not initialized.' )

    if ( .not. constants_inited ) call MessageNotify( 'E', module_name, '"constants" module is not initialized.' )

    if ( .not. axesset_inited ) call MessageNotify( 'E', module_name, '"axesset" module is not initialized.' )

    if ( .not. timeset_inited ) call MessageNotify( 'E', module_name, '"timeset" module is not initialized.' )

  end subroutine InitCheck
Subroutine :
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_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_GradLonPiN(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ frac{1}{cos varphi} DP{pi}{lambda} $
xy_GradLatPiN(0:imax-1, 1:jmax) :real(DP), intent(in)
: $ DP{pi}{varphi} $
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_DTempDtN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ hat{H} (t) $ . 温度時間変化項. Temperature tendency
xyz_DQVapDtN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ R (t) $ . 比湿時間変化項. Specific humidity 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_DPiDtN(0:imax-1, 1:jmax) :real(DP), intent(out)
: $ Z (t) $ . 地表面気圧時間変化項. Surface pressure tendency
xyz_QVapUAdvN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ uq (t) $ . 比湿東西移流項. Eastward advection of specific humidity
xyz_QVapVAdvN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: $ vq (t) $ . 比湿南北移流項. Northward advection of specific humidity

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

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

[Source]

  subroutine NonLinearOnGrid( xyz_UN,        xyz_VN, xyz_VorN,      xyz_DivN, xyz_TempN,     xyz_QVapN, xy_GradLonPiN, xy_GradLatPiN, xyz_UAdvN,     xyz_VAdvN, xyz_DTempDtN,  xyz_DQVapDtN, xyz_KinEngyN, xyz_TempUAdvN, xyz_TempVAdvN, xyr_SigDotN, xy_DPiDtN, xyz_QVapUAdvN, xyz_QVapVAdvN )
    !
    ! 非線形項 (非重力波項) を格子点上で計算します.
    !
    ! 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: RPlanet, CpDry, EpsV
                              ! $ \epsilon_v $ . 
                              ! 水蒸気分子量比. 
                              ! Molecular weight of water vapor

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

    implicit none
    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_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_GradLonPiN (0:imax-1, 1:jmax)
                              ! $ \frac{1}{\cos \varphi} \DP{\pi}{\lambda} $
    real(DP), intent(in):: xy_GradLatPiN (0:imax-1, 1:jmax)
                              ! $ \DP{\pi}{\varphi} $
    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_DTempDtN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \hat{H} (t) $ . 温度時間変化項. 
                              ! Temperature tendency
    real(DP), intent(out):: xyz_DQVapDtN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ R (t) $ . 比湿時間変化項. 
                              ! Specific humidity 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_DPiDtN (0:imax-1, 1:jmax)
                              ! $ Z (t) $ . 地表面気圧時間変化項. 
                              ! Surface pressure tendency
    real(DP), intent(out):: xyz_QVapUAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ uq (t) $ . 比湿東西移流項. 
                              ! Eastward advection of specific humidity
    real(DP), intent(out):: xyz_QVapVAdvN (0:imax-1, 1:jmax, 1:kmax)
                              ! $ vq (t) $ . 比湿南北移流項. 
                              ! Northward advection of specific humidity

    !-----------------------------------
    !  作業変数
    !  Work variables
    real(DP):: xyz_PiAdv (0:imax-1, 1:jmax, 1:kmax)
                              ! $ \Dvect{v} \cdot \nabla \pi $ . 
                              ! $ \pi $ の移流. Advection of $ \pi $
    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_TempVir (0:imax-1, 1:jmax, 1:kmax)
                              ! $ T_v $ . 
                              ! 仮温度. Virtual temperature
    real(DP):: xyz_TempVirEdd (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

    ! 実行文 ; Executable statement
    !

    ! $ \pi $ の移流, $ \pi $ 移流の積下げ, 発散の積下げの計算
    ! Calculate advection of $ \pi $, integral downward of advection 
    !   of $ \pi $, integral downward of divergence
    !
    do k = 1, kmax
      xyz_PiAdv(:,:,k) = ( xyz_UN(:,:,k) * xy_GradLonPiN + xyz_VN(:,:,k) * xy_GradLatPiN )
    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 surface pressure tendency $ Z $
    !
    xy_DPiDtN = - xyz_PiAdvSum(:,:,1)

    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('explicit')

      xy_DPiDtN = xy_DPiDtN - xyz_DivSum(:,:,1)

    end select

    ! $ \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_TempVir(:,:,k) = xyz_TempN(:,:,k) * ( 1. + ((( 1. / EpsV ) - 1. ) * xyz_QVapN(:,:,k)) )

      xyz_TempEdd(:,:,k) = xyz_TempN(:,:,k) - z_RefTemp(k)
      xyz_TempVirEdd(:,:,k) = xyz_TempVir(:,:,k) - z_RefTemp(k)
    enddo

    ! 温度の擾乱 (半整数レベル) の計算
    ! Calculate temperature eddy (half level)
    !
    xyr_TempEdd(:,:,0) = 0.
    xyr_TempEdd(:,:,kmax) = 0.
    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_VN(:,:,1) - 1. / ( 2. * z_DelSigma(1) ) * xyr_SigDotN(:,:,1) * ( xyz_UN(:,:,1) - xyz_UN(:,:,2) ) - CpDry * z_TInpCoefK(1) * xyz_TempVirEdd(:,:,1) * xy_GradLonPiN

    do k = 2, kmax-1
      xyz_UAdvN(:,:,k) = ( xyz_VorN(:,:,k) + xy_Cori ) * xyz_VN(:,:,k) - 1. / ( 2. * z_DelSigma(k) ) * ( xyr_SigDotN(:,:,k-1) * ( xyz_UN(:,:,k-1) - xyz_UN(:,:,k) ) + xyr_SigDotN(:,:,k)   * ( xyz_UN(:,:,k)   - xyz_UN(:,:,k+1) ) ) - CpDry * z_TInpCoefK(k) * xyz_TempVirEdd(:,:,k) * xy_GradLonPiN
    end do

    xyz_UAdvN(:,:,kmax) = ( xyz_VorN(:,:,kmax) + xy_Cori ) * xyz_VN(:,:,kmax) - 1. / ( 2. * z_DelSigma(kmax) ) * xyr_SigDotN(:,:,kmax-1) * ( xyz_UN(:,:,kmax-1) - xyz_UN(:,:,kmax) ) - CpDry * z_TInpCoefK(kmax) * xyz_TempVirEdd(:,:,kmax) * xy_GradLonPiN


    ! 南北運動量移流項の計算
    ! Calculate advection of northward momentum
    !
    xyz_VAdvN(:,:,1) = - ( xyz_VorN(:,:,1) + xy_Cori ) * xyz_UN(:,:,1) - 1. / ( 2. * z_DelSigma(1) ) * xyr_SigDotN(:,:,1) * ( xyz_VN(:,:,1) - xyz_VN(:,:,2) ) - CpDry * z_TInpCoefK(1) * xyz_TempVirEdd(:,:,1) * xy_GradLatPiN

    do k = 2, kmax-1
      xyz_VAdvN(:,:,k) = - ( xyz_VorN(:,:,k) + xy_Cori ) * xyz_UN(:,:,k) - 1. / ( 2. * z_DelSigma(k) ) * ( xyr_SigDotN(:,:,k-1) * ( xyz_VN(:,:,k-1) - xyz_VN(:,:,k) ) + xyr_SigDotN(:,:,k)   * ( xyz_VN(:,:,k)   - xyz_VN(:,:,k+1) ) ) - CpDry * z_TInpCoefK(k) * xyz_TempVirEdd(:,:,k) * xy_GradLatPiN
    end do

    xyz_VAdvN(:,:,kmax) = - ( xyz_VorN(:,:,kmax) + xy_Cori ) * xyz_UN(:,:,kmax) - 1. / ( 2. * z_DelSigma(kmax) ) * xyr_SigDotN(:,:,kmax-1) * ( xyz_VN(:,:,kmax-1) - xyz_VN(:,:,kmax) ) - CpDry * z_TInpCoefK(kmax) * xyz_TempVirEdd(:,:,kmax) * xy_GradLatPiN


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

    xyz_KinEngyN = ( xyz_UN**2 + xyz_VN**2 ) / 2. + xyz_KinEngyN


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

    ! 温度の時間変化項 $ \hat{H} $ の計算
    ! Calculate temperature tendency term $ \hat{H} $
    !
    do k = 1, kmax-1
      xyz_DTempDtN(:,:,k) = xyz_TempEdd(:,:,k) * xyz_DivN(:,:,k) - 1. / z_DelSigma(k) * (   xyr_SigDotN(:,:,k-1) * ( xyr_TempEdd(:,:,k-1) - xyz_TempEdd(:,:,k) ) + xyr_SigDotN(:,:,k) * ( xyz_TempEdd(:,:,k)   - xyr_TempEdd(:,:,k) ) ) - 1. / 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_TempVir(:,:,k) * xyz_PiAdv(:,:,k) - z_HydroAlpha(k) / z_DelSigma(k) * (   xyz_TempVir(:,:,k)    * xyz_PiAdvSum(:,:,k) + xyz_TempVirEdd(:,:,k) * xyz_DivSum(:,:,k) ) - z_HydroBeta(k) / z_DelSigma(k) * (   xyz_TempVir(:,:,k)    * xyz_PiAdvSum(:,:,k+1) + xyz_TempVirEdd(:,:,k) * xyz_DivSum(:,:,k+1) )
    enddo

    xyz_DTempDtN(:,:,kmax) = xyz_TempEdd(:,:,kmax) * xyz_DivN(:,:,kmax) - 1. / z_DelSigma(kmax) * (   xyr_SigDotN(:,:,kmax-1) * ( xyr_TempEdd(:,:,kmax-1) - xyz_TempEdd(:,:,kmax) ) + xyr_SigDotN(:,:,kmax) * ( xyz_TempEdd(:,:,kmax)   - xyr_TempEdd(:,:,kmax) ) ) - 1. / 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_TempVir(:,:,kmax) * xyz_PiAdv(:,:,kmax) - z_HydroAlpha(kmax) / z_DelSigma(kmax) * (   xyz_TempVir(:,:,kmax)    * xyz_PiAdvSum(:,:,kmax) + xyz_TempVirEdd(:,:,kmax) * xyz_DivSum(:,:,kmax) )


    ! 比湿東西移流項, 比湿南北移流項の計算
    ! Calculate eastward and northward advection of specific humidity
    !
    do k = 1, kmax
      xyz_QVapUAdvN(:,:,k) =  xyz_UN(:,:,k) * xyz_QVapN(:,:,k)
      xyz_QVapVAdvN(:,:,k) =  xyz_VN(:,:,k) * xyz_QVapN(:,:,k)
    end do

    ! 比湿時間変化項 $ R $ の計算
    ! Calculate specific humidity tendency $ R $
    !
    xyz_DQVapDtN(:,:,1) = xyz_QVapN(:,:,1) * xyz_DivN(:,:,1) - 1. / ( 2. * z_DelSigma(1) ) * xyr_SigDotN(:,:,1) * ( xyz_QVapN(:,:,1) - xyz_QVapN(:,:,2) )

    do k = 2, kmax - 1
      xyz_DQVapDtN(:,:,k) = xyz_QVapN(:,:,k) * xyz_DivN(:,:,k) - 1. / ( 2. * z_DelSigma(k) ) * ( xyr_SigDotN(:,:,k-1) * ( xyz_QVapN(:,:,k-1) - xyz_QVapN(:,:,k)   ) + xyr_SigDotN(:,:,k) * ( xyz_QVapN(:,:,k)   - xyz_QVapN(:,:,k+1) )   )
    end do

    xyz_DQVapDtN(:,:,kmax) = xyz_QVapN(:,:,kmax) * xyz_DivN(:,:,kmax) - 1. / ( 2. * z_DelSigma(kmax) ) * xyr_SigDotN(:,:,kmax-1) * ( xyz_QVapN(:,:,kmax-1) - xyz_QVapN(:,:,kmax) )

  end subroutine NonLinearOnGrid
Subroutine :

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

Configure coefficients for "TimeIntegration". Setting values that are set last time are used, except when first time and Δ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 Δ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

    ! 座標データ設定
    ! 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)
    !
    use wa_module, only: l_nm

    ! LU 分解法により連立 1 次方程式を解くための関数
    ! Functions to solve linear simultaneous equation by LU decomposition
    !
    use lumatrix, only: LUDecomp

    ! メッセージ出力
    ! Message output
    !
    use dc_message, only: MessageNotify

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

    ! 宣言文 ; Declaration statements
    !
    implicit none

    ! TimeIntegration 等で使用する係数の設定のための作業変数
    ! Work variable for coefficients for "TimeIntegration", etc.
    !
    real(DP), allocatable:: zz_siMtxW (:,:)
                              ! $ W $ . 
                              ! 発散の式での線形重力波項の効果による温度補正の係数. 
                              ! Coefficient for correction of temperature 
                              ! by effort of linear gravitational terms

    real(DP), allocatable:: zz_siMtxQ (:,:)
                              ! $ Q = \DD{T}{\sigma} $
    real(DP), allocatable:: zz_siMtxS (:,:)
                              ! $ S = \DD{\sigma}{t} $
    real(DP), allocatable:: zz_siMtxQS (:,:)
                              ! $ QS $ . 
                              ! この変数は $ \sigma $ 移流の効果に相当. 
                              ! This variable  corresponds to effort of $ \sigma $ advection
    real(DP), allocatable:: zz_siMtxR (:,:)
                              ! $ R $ . 
                              ! 発散と掛け合わせることで気圧変化の効果となる.
                              ! Product R and divergence become effort of 
                              ! surface pressure tendency.
                              ! $ RD = \kappa T 
                              ! (\DD{\pi}{t} + \Dinv{\sigma}\DD{\sigma}{t}) $ .

    integer, allocatable:: nmo (:,:,:)
                              ! スペクトルの添字順番. 
                              ! Spectral subscript expression
    real(DP), allocatable:: zz_siMtxM (:,:)
                              ! 行列 $ \underline{M} $. 
                              ! Matrix $ \underline{M} $
    integer, allocatable:: z_siMtxPivWork(:)
                              ! 行列のピボット作成の作業変数. 
                              ! Work variable for pivot
    real(DP):: flapla         ! $ \nabla_{\sigma}^{2} $


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

    integer:: mmax            ! 最大東西波数. 
                              ! Maximum truncated eastward wavenumber
    integer:: lmax            ! 最大南北波数. 
                              ! Maximum truncated northward wavenumber
    integer:: n, m, nm, mxnm
                              ! 波数方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in wavenumber direction

    logical:: lmax_err        ! 最大南北波数に関するエラーフラグ. 
                              ! Error flag for maximum truncated northward wavenumber

    ! 実行文 ; 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_siMtxG   ) )  allocate( z_siMtxG    (1:kmax) )
    if ( .not. allocated( zz_siMtxH  ) )  allocate( zz_siMtxH   (1:kmax, 1:kmax) )
    if ( .not. allocated( zz_siMtxWH ) )  allocate( zz_siMtxWH  (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_siMtxQS ) )  allocate( zz_siMtxQS (1:kmax, 1:kmax) )
    if ( .not. allocated( zz_siMtxR  ) )  allocate( zz_siMtxR  (1:kmax, 1:kmax) )


    z_siMtxG = CpDry * z_TInpCoefK * z_RefTemp

    do k = 1, kmax
      do l = 1, kmax
        zz_siMtxGCt(k,l) = z_siMtxG(k) * z_DelSigma(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_siMtxQS = 0.
    zz_siMtxQS = matmul(zz_siMtxQ, zz_siMtxS)

    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 = 0.
    zz_siMtxH = zz_siMtxQS - zz_siMtxR

    zz_siMtxWH = 0.
    zz_siMtxWH = matmul(zz_siMtxW, zz_siMtxH)

    if ( .not. allocated(nmo           ) )  allocate( nmo           (1:2, 0:nmax, 0:nmax) )
    if ( .not. allocated(zz_siMtxM     ) )  allocate( zz_siMtxM     (1:kmax, 1:kmax) )
    if ( .not. allocated(z_siMtxPivWork) )  allocate( z_siMtxPivWork(1:kmax) )
    if ( .not. allocated(wzz_siMtxLU   ) )  allocate( wzz_siMtxLU   ((nmax+1)**2, 1:kmax, 1:kmax) )
    if ( .not. allocated(wz_siMtxPiv   ) )  allocate( wz_siMtxPiv   ((nmax+1)**2, 1:kmax) )

    mmax = nmax
    lmax = nmax
    mxnm = 0

    ! スペクトル添字順序の取り出し
    ! Fetch spectral subscript expression
    !
    nmo = 0
    do l = 0, lmax
      do m = 0, min(mmax, nmax-l)
        nmo(1,m,l) = l_nm(m+l,  m)
        nmo(2,m,l) = l_nm(m+l, -m)
      enddo
    enddo

    Loop_n: do n = 0, nmax
      flapla = - real(n) * real(n+1) / RPlanet**2

      ! スペクトル添字順序の取り出し
      ! Fetch spectral subscript expression
      !
      lmax_err = .true.
      do m = 0, min(n,mmax)
        if ( n-m <= lmax ) then
          nm = nmo(1,m,n-m)
          lmax_err = .false.
          exit
        endif
      end do
      if ( lmax_err ) then
        call MessageNotify( 'E', module_name, 'n-m=<%d> must be less than or equal to lmax=<%d>.', i = (/ n-m, lmax /) )
      end if

      ! 行列 $ \underline{M} $ の計算
      ! Calculate matrix $ \underline{M} $
      !
      do k = 1, kmax
        do kk = 1, kmax
          zz_siMtxM ( k,kk ) = - DelTime**2 * flapla * (   zz_siMtxWH( k,kk ) + zz_siMtxGCt( k,kk ) * ( 1. - 2. * DelTime * wz_HDifCoefT(nm,1) )  )
          if ( k == kk ) then
            zz_siMtxM ( k,kk ) = zz_siMtxM ( k,kk ) +   ( 1. - 2. * DelTime * wz_HDifCoefM(nm,1) ) * ( 1. - 2. * DelTime * wz_HDifCoefT(nm,1) )
          endif
        end do
      end do

      ! LU 行列計算
      ! LU matrix calculation
      !
      call LUDecomp( zz_siMtxM, z_siMtxPivWork ) ! (out)

      ! ダミー値の代入. (位置 kmax はまだ未定義なため).
      ! Dummy value is subtituted. (Because position kmax is undefined yet).
      !
      z_siMtxPivWork(kmax) = 0

      ! 配列の詰め替え
      ! Repack matrices
      !
      do m = 0, mmax
        l = n - m
        if ( ( l >= 0 ) .and. ( l <= lmax ) ) then
          do k = 1, kmax
            do kk = 1, kmax
              wzz_siMtxLU ( nmo(1,m,l),k,kk ) = zz_siMtxM ( k,kk )
              wzz_siMtxLU ( nmo(2,m,l),k,kk ) = zz_siMtxM ( k,kk )
            end do
            wz_siMtxPiv ( nmo(1,m,l),k ) = z_siMtxPivWork ( k )
            wz_siMtxPiv ( nmo(2,m,l),k ) = z_siMtxPivWork ( k )
            mxnm = max( mxnm, nmo(1,m,l) )
            mxnm = max( mxnm, nmo(2,m,l) )
          end do
        endif
      end do

    end do Loop_n

    do nm = mxnm+1, (nmax+1)**2
      do k = 1, kmax
        do kk = 1, kmax
          wzz_siMtxLU ( nm,k,kk ) = zz_siMtxM ( k,kk )
        end do
        wz_siMtxPiv ( nm,k ) = z_siMtxPivWork ( k )
      end do
    end do

  end subroutine SemiImplMatrix
TimeIntegScheme
Variable :
TimeIntegScheme :character(TOKEN), save
: 時間積分法. 以下の方法を選択可能.

Time integration scheme. Available schemes are as follows.

  • "Semi-implicit"
  • "Explicit"
Subroutine :
wz_DVorDtN((nmax+1)**2, 1:kmax) :real(DP), intent(in)
: $ DD{zeta}{t} (t) $ . 渦度変化 (スペクトル). Vorticity tendency (spectral)
wz_DDivDtN((nmax+1)**2, 1:kmax) :real(DP), intent(in)
: $ DD{D}{t} (t) $ . 発散変化 (スペクトル). Divergence tendency (spectral)
wz_DTempDtN((nmax+1)**2, 1:kmax) :real(DP), intent(in)
: $ DD{T}{t} (t) $ . 温度変化 (スペクトル). Temperature tendency (spectral)
wz_DQVapDtN((nmax+1)**2, 1:kmax) :real(DP), intent(in)
: $ DD{q}{t} (t) $ . 比湿変化 (スペクトル). Specific humidity tendency (spectral)
w_DPiDtN((nmax+1)**2) :real(DP), intent(in)
: $ DD{p_s}{t} (t) $ . 地表面気圧変化 (スペクトル). Surface pressure tendency (spectral)
wz_VorB((nmax+1)**2, 1:kmax) :real(DP), intent(in)
: $ zeta (t-\Delta t) $ . 渦度 (スペクトル). Vorticity (spectral)
wz_DivB((nmax+1)**2, 1:kmax) :real(DP), intent(in)
: $ D (t-\Delta t) $ . 発散 (スペクトル). Divergence (spectral)
wz_TempB((nmax+1)**2, 1:kmax) :real(DP), intent(in)
: $ T (t-\Delta t) $ . 温度 (スペクトル). Temperature (spectral)
wz_QVapB((nmax+1)**2, 1:kmax) :real(DP), intent(in)
: $ q (t-\Delta t) $ . 比湿 (スペクトル). Specific humidity (spectral)
w_PiB((nmax+1)**2) :real(DP), intent(in)
: $ pi = ln p_s (t-\Delta t) $ . 地表面気圧 (スペクトル). Surface pressure (spectral)
w_SurfGeoPot((nmax+1)**2) :real(DP), intent(in)
: $ Phi_s $ . 地表ジオポテンシャル. Surface geo-potential
wz_VorA((nmax+1)**2, 1:kmax) :real(DP), intent(out)
: $ zeta (t+Delta t) $ . 渦度 (スペクトル). Vorticity (spectral)
wz_DivA((nmax+1)**2, 1:kmax) :real(DP), intent(out)
: $ D (t+Delta t) $ . 発散 (スペクトル). Divergence (spectral)
wz_TempA((nmax+1)**2, 1:kmax) :real(DP), intent(out)
: $ T (t+Delta t) $ . 温度 (スペクトル). Temperature (spectral)
wz_QVapA((nmax+1)**2, 1:kmax) :real(DP), intent(out)
: $ q (t+Delta t) $ . 比湿 (スペクトル). Specific humidity (spectral)
w_PiA((nmax+1)**2) :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_nmlTimeIntegScheme を変更することで, 重力波項をエクスプリシット法によって解くことも可能です.

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

Leap-frog scheme is used as time integration scheme. And forward difference is used to diffusion terms. 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 TimeIntegration( wz_DVorDtN, wz_DDivDtN, wz_DTempDtN, wz_DQVapDtN, w_DPiDtN, wz_VorB,    wz_DivB,    wz_TempB,    wz_QVapB,    w_PiB, w_SurfGeoPot, wz_VorA,    wz_DivA,    wz_TempA,    wz_QVapA,    w_PiA )
    !
    ! 時間積分を行い, 
    ! 時刻 $ t $ の物理量の時間変化と $ t-\Delta t$ の物理量から
    ! 時刻 $ t+\Delta t $ の物理量を計算します.
    !
    ! 時間積分法にはリープフロッグスキームを用いています. 
    ! ただし拡散項には前方差分を用いています. 
    ! デフォルトでは, $ \Delta t $ を大きくとるために, 重力波項に 
    ! セミインプリシット法を適用しています.
    ! NAMELIST#dynamics_hspl_vas83_nml の TimeIntegScheme を変更することで, 
    ! 重力波項をエクスプリシット法によって解くことも可能です.
    ! 
    ! With time integration, calculate physical values at $ t+\Delta t $
    ! from tendency at $ t $ and physical values at $ t-\Delta t $ .
    !
    ! Leap-frog scheme is used as time integration scheme. 
    ! And forward difference is used to diffusion terms.
    ! 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
    !

    ! 物理定数設定
    ! 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 次方程式を解くための関数
    ! Functions to solve linear simultaneous equation by LU decomposition
    !
    use lumatrix, only: LUSolve

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

    implicit none
    real(DP), intent(in):: wz_DVorDtN ((nmax+1)**2, 1:kmax)
                              ! $ \DD{\zeta}{t} (t) $ . 渦度変化 (スペクトル). 
                              ! Vorticity tendency (spectral)
    real(DP), intent(in):: wz_DDivDtN ((nmax+1)**2, 1:kmax)
                              ! $ \DD{D}{t} (t) $ . 発散変化 (スペクトル). 
                              ! Divergence tendency (spectral)
    real(DP), intent(in):: wz_DTempDtN ((nmax+1)**2, 1:kmax)
                              ! $ \DD{T}{t} (t) $ . 温度変化 (スペクトル). 
                              ! Temperature tendency (spectral)
    real(DP), intent(in):: wz_DQVapDtN ((nmax+1)**2, 1:kmax)
                              ! $ \DD{q}{t} (t) $ . 比湿変化 (スペクトル). 
                              ! Specific humidity tendency (spectral)
    real(DP), intent(in):: w_DPiDtN ((nmax+1)**2)
                              ! $ \DD{p_s}{t} (t) $ . 地表面気圧変化 (スペクトル). 
                              ! Surface pressure tendency (spectral)
    real(DP), intent(in):: wz_VorB ((nmax+1)**2, 1:kmax)
                              ! $ \zeta (t-\Delta t) $ . 渦度 (スペクトル). 
                              ! Vorticity (spectral)
    real(DP), intent(in):: wz_DivB ((nmax+1)**2, 1:kmax)
                              ! $ D (t-\Delta t) $ . 発散 (スペクトル). 
                              ! Divergence (spectral)
    real(DP), intent(in):: wz_TempB ((nmax+1)**2, 1:kmax)
                              ! $ T (t-\Delta t) $ . 温度 (スペクトル). 
                              ! Temperature (spectral)
    real(DP), intent(in):: w_PiB ((nmax+1)**2)
                              ! $ \pi = \ln p_s (t-\Delta t) $ . 地表面気圧 (スペクトル). 
                              ! Surface pressure (spectral)
    real(DP), intent(in):: wz_QVapB ((nmax+1)**2, 1:kmax)
                              ! $ q (t-\Delta t) $ . 比湿 (スペクトル). 
                              ! Specific humidity (spectral)
    real(DP), intent(in):: w_SurfGeoPot ((nmax+1)**2)
                              ! $ \Phi_s $ . 地表ジオポテンシャル. 
                              ! Surface geo-potential

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

    ! 作業変数
    ! Work variables
    !
    real(DP):: wz_siTemp ((nmax+1)**2, 1:kmax)
                              ! 温度 (セミインプリシット法のための作業変数). 
                              ! Temperature (work variable for semi-implicit scheme)
    real(DP):: wz_siDTempDt ((nmax+1)**2, 1:kmax)
                              ! $ \DD{T}{t} (t) $ . 温度変化 (スペクトル) の作業変数. 
                              ! Temperature tendency (spectral) work variable

    real(DP):: w_siPi ((nmax+1)**2)
                              ! $ \pi $ (セミインプリシット法のための作業変数). 
                              ! $ \pi $ (work variable for semi-implicit scheme)
    real(DP):: w_siDPiDt ((nmax+1)**2)
                              ! $ \DD{p_s}{t} (t) $ . 地表面気圧変化 (スペクトル). 
                              ! Surface pressure tendency (spectral)
    real(DP):: wz_siPhi ((nmax+1)**2, 1:kmax)
                              ! $ \Phi = \underline{W} \overline{ \Dvect{T} }^{t}$ .
                              ! (セミインプリシット法のための作業変数). 
                              ! (Work variable for semi-implicit scheme)
    real(DP):: wz_siDivAvrTime ((nmax+1)**2, 1:kmax)
                              ! $ \Dvect{f} $ . 
                              ! 時間平均の $ \Dvect{D} $ (セミインプリシット法のための作業変数). 
                              ! Time average $ \Dvect{D} $ (work variable for semi-implicit scheme)

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

    ! 実行文 ; Executable statement
    !

    ! 非重力波項の仮積分
    ! Integration non gravitational terms temporarily
    !
    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('semi-implicit')

      wz_siTemp = ( 1. - DelTime * wz_HDifCoefT ) * wz_TempB + DelTime * wz_DTempDtN

      w_siPi = w_PiB + DelTime * w_DPiDtN

    end select

    ! ジオポテンシャルの計算
    ! Calculate geopotential
    !
    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('semi-implicit')

      wz_siPhi (:,1) = CpDry * z_HydroAlpha(1) * wz_siTemp (:,1)

      do k = 2, kmax
        wz_siPhi (:,k) = wz_siPhi(:,k-1) + CpDry * z_HydroAlpha(k) * wz_siTemp (:,k) + CpDry * z_HydroBeta (k-1) * wz_siTemp (:,k-1)
      end do

    case ('explicit')

    end select

    ! 発散方程式の右辺の計算
    ! Calculate right side of divergence equation
    !
    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('semi-implicit')

      do k = 1, kmax
        wz_siDivAvrTime(:,k) = ( 1. - 2. * DelTime * wz_HDifCoefT(:,k)  ) * ( 1. -      DelTime * wz_HDifCoefM(:,k) ) * wz_DivB(:,k) + ( 1. - 2. * DelTime * wz_HDifCoefT(:,k)  ) * DelTime * wz_DDivDtN(:,k) - DelTime * wz_LaplaCoef(:,k) * (   ( 1. - 2. * DelTime * wz_HDifCoefT(:,k) ) * w_SurfGeoPot + wz_siPhi(:,k) + ( 1. - 2. * DelTime * wz_HDifCoefT(:,k) ) * z_siMtxG(k) * w_siPi                  )
      end do

    end select

    ! 時間平均の $ \Dvect{D} $ を LU 行列で解く
    ! Solve time average $ \Dvect{D} $ with LU matrix
    !
    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('semi-implicit')

      wz_siDivAvrTime = LUSolve( wzz_siMtxLU, wz_siMtxPiv, wz_siDivAvrTime )

    end select

    ! 温度, 地表気圧の時間変化項の計算
    ! Calculate tendency of temperature and surface pressure
    !
    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('semi-implicit')

      wz_siDTempDt = wz_DTempDtN
      do k = 1, kmax
        do kk = 1, kmax
          wz_siDTempDt(:,k) = wz_siDTempDt(:,k) - zz_siMtxH(k,kk) * wz_siDivAvrTime(:,kk)
        end do
      end do

      w_siDPiDt = w_DPiDtN
      do kk = 1, kmax
        w_siDPiDt = w_siDPiDt - z_DelSigma(kk) * wz_siDivAvrTime(:,kk)
      end do

    end select

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

    ! 渦度 ; Vorticity
    !
    wz_VorA = ( 1. / ( 1. - 2. * DelTime * wz_HDifCoefM ) ) * ( wz_VorB + 2. * DelTime * wz_DVorDtN )

    ! 発散 ; Divergence
    !
    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('semi-implicit')

      wz_DivA  = 2. * wz_siDivAvrTime - wz_DivB

    case ('explicit')

      wz_DivA  = ( 1. / ( 1. - 2. * DelTime * wz_HDifCoefM ) ) * ( wz_DivB + 2. * DelTime * wz_DDivDtN )

    end select

    ! 温度 ; Temperature
    !
    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('semi-implicit')

      wz_TempA = ( 1.  / ( 1. - 2. * DelTime * wz_HDifCoefT ) ) * ( wz_TempB + wz_siDTempDt * 2. * DelTime )

    case ('explicit')

      wz_TempA = ( 1.  / ( 1. - 2. * DelTime * wz_HDifCoefT ) ) * ( wz_TempB + wz_DTempDtN * 2. * DelTime )

    end select

    ! 比湿 ; Specific humidity
    !
    wz_QVapA = ( 1. / ( 1. - 2. * DelTime * wz_HDifCoefT ) ) * ( wz_QVapB + wz_DQVapDtN * 2. * DelTime )

    ! 地表面気圧 ; Surface pressure
    !
    select case ( LChar( trim( TimeIntegScheme ) ) )
    case ('semi-implicit')

      w_PiA  = w_PiB + 2. * DelTime * w_siDPiDt

    case ('explicit')

      w_PiA  = w_PiB + 2. * DelTime * w_DPiDtN

    end select

  end subroutine TimeIntegration
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: dcpam5-20090225-2 $’ // ’$Id: dynamics_hspl_vas83.F90,v 1.30 2009-02-24 22:26:43 morikawa Exp $’ :character(*), parameter
: モジュールのバージョン Module version
wz_HDifCoefM
Variable :
wz_HDifCoefM(:,:) :real(DP), allocatable
: $ -K_{HD} [(-1)^{N_D/2}nabla^{N_D}- (2/a^2)^{N_D/2}] $ . 運動量水平拡散係数. Coefficient of horizontal momentum diffusion
wz_HDifCoefT
Variable :
wz_HDifCoefT(:,:) :real(DP), allocatable
: $ -(-1)^{N_D/2}K_{HD}nabla^{N_D} $ . 熱, 水水平拡散係数. Coefficient of horizontal thermal and water diffusion
wz_LaplaCoef
Variable :
wz_LaplaCoef(:,:) :real(DP), allocatable
: $ nabla^2_{sigma} = -n times (n+1) $ . ラプラシアンの係数. Laplacian coefficient
wz_siMtxPiv
Variable :
wz_siMtxPiv(:,:) :integer, allocatable
: セミインプリシット行列のピボット. Pivot for semi-implicit matrix
wzz_siMtxLU
Variable :
wzz_siMtxLU(:,:,:) :real(DP), allocatable
: セミインプリシット行列の LU 分解. LU decomposition for semi-implicit matrix
xy_Cori
Variable :
xy_Cori(:,:) :real(DP), allocatable
: $ f\equiv 2\Omega\sin\varphi $ . コリオリパラメータ. Coriolis parameter
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_siMtxG
Variable :
z_siMtxG(:) :real(DP), allocatable
: $ G = C_p kappa T $
zz_siMtxGCt
Variable :
zz_siMtxGCt(:,:) :real(DP), allocatable
: $ G C^{T} $ ( $ C = Delta sigma $ )
zz_siMtxH
Variable :
zz_siMtxH(:,:) :real(DP), allocatable
: $ h = QS - R $
zz_siMtxWH
Variable :
zz_siMtxWH(:,:) :real(DP), allocatable
: $ W h $

[Validate]