program arare
  !
  ! 非静力学モデル deepconv/arare 木星大気対流計算用主プログラム 
  !
  ! モジュール引用 use statement 
  !
  ! gtool5 関連 
  ! gtool5 modules
  !
  use dc_types,       only : STRING
  use dc_string,      only : StoA
  use dc_message,     only: MessageNotify
  ! 初期設定モジュール
  ! Initialize module
  !
  use argset,        only: argset_init  
  use fileset,       only: fileset_init, InitFile
  use debugset,      only: debugset_init
  use timeset,       only: timeset_init, DelTimeLong, DelTimeShort, NstepDisp, NstepLong, NstepShort
  use gridset,       only: gridset_init, DimXMin, DimXMax, DimZMin, DimZMax, SpcNum
  use basicset,      only: basicset_init, xza_MixRtBasicZ, xz_DensBasicZ, xz_PotTempBasicZ, xz_VelSoundBasicZ
  ! 化学量計算モジュール 
  ! Chemical calculation modules
  !
  use ChemCalc,      only: ChemCalc_init
  use chemdata,      only: chemdata_init
  ! 力学過程計算用関数モジュール
  ! Dynamical processes module
  !
  use DynFunc,       only: xz_AdvScalar, xz_AdvKm, xza_AdvScalar, pz_AdvVelX, xr_Buoy, xr_AdvVelZ, pz_GradPi
  use DynImpFunc,    only: xz_Exner_init, xz_Exner, xr_GradPi
  ! 乱流拡散計算用モジュール
  ! Turbulent diffusion module
  !
  use Turbulence,    only: Turbulence_Init, xz_TurbScalar, xza_TurbScalar, pz_TurbVelX, xr_TurbVelZ, xz_ShearKm, xz_DispKm, xz_DispHeat
  ! 境界からのフラックス計算用モジュール
  ! Surface flux module
  !
  use HeatFlux,      only: xz_HeatFluxDiff, xza_MixRtFluxDiff
  ! 放射強制計算用モジュール
  ! Radiative forceing module
  !
  use Radiation,     only: Radiation_init, xz_RadHeatConst  
  ! 湿潤過程計算用モジュール
  ! Moist processes modules
  !
  use moistset,      only: moistset_init
  use MoistAdjust,   only: MoistAdjustSvapPress, MoistAdjustNH4SH
  use WarmRainPrm,   only: WarmRainPrm_Init, xz_Rain2GasHeat, xza_Rain2Gas, xza_Rain2GasNH4SH, xz_Rain2GasHeatNH4SH, xza_Cloud2Rain, xza_FallRain
  use MoistBuoyancy, only: MoistBuoy_Init, xz_BuoyMoistKm, xr_BuoyMolWt, xr_BuoyDrag
  use fillnegative,  only: FillNegative_init, xza_FillNegative_xza
  ! 安定度の計算
  ! Static stability calculation module
  !
  use ECCM, only: ECCM_Stab
  ! 数値拡散/摩擦計算用モジュール
  ! Numerical diffussion /dumping module
  !
  use NumDiffusion,  only: NumDiffusion_Init, xz_NumDiffScalar, xz_NumDiffKm, xza_NumDiffScalar, pz_NumDiffVelX, xr_NumDiffVelZ
  use damping,       only: damping_init, DampSponge_xz, DampSponge_pz, DampSponge_xr
  ! 積算値管理モジュール
  ! Monitor variables setup modules
  !
  use StorePotTemp,  only: StorePotTemp_init, StorePotTempClean, StorePotTempCond
  use StoreMixRt,    only: StoreMixRt_init, StoreMixRtClean, StoreMixRtCond, StoreMixRtFill1, StoreMixRtFill2
  use StoreBuoy,     only: StoreBuoy_init, StoreBuoyClean
  use StoreStab,     only: StoreStab_init, StoreStabClean
  ! ファイル入出力モジュール
  ! File I/O module
  !
  use RestartFileIO, only: ReStartFile_Open, ReStartFile_OutPut, ReStartFile_Close, ReStartFile_Get
  use HistoryFileIO, only: HistoryFile_Open, HistoryFile_OutPut, HistoryFile_Close
  ! 下請けモジュール
  ! Utility modules
  !
  use cflcheck,      only: CFLCheckTimeShort, CFLCheckTimeLongVelX, CFLCheckTimeLongVelZ
  use timefilter,    only: AsselinFilter_xz, AsselinFilter_xr, AsselinFilter_pz, AsselinFilter_xza
  use boundary,      only: BoundaryXCyc_xz, BoundaryZSym_xz, BoundaryXCyc_xza, BoundaryZSym_xza, BoundaryXCyc_pz, BoundaryZSym_pz, BoundaryXCyc_xr, BoundaryZAntiSym_xr
  ! 変数宣言 Variables definition 
  !
  ! 内部変数
  ! Internal variables
  !
  character(80) :: cfgfile               ! NAMELIST ファイル名 NAMELIST file name
  real(8), allocatable :: pz_VelXBl(:,:) 
                            ! $ u(t-\Delta t) $ 水平風 Horizontal wind
  real(8), allocatable :: pz_VelXNl(:,:) 
                            ! $ u(t) $ 水平風 Horizontal wind
  real(8), allocatable :: pz_VelXAl(:,:) 
                            ! $ u (t+\Delta t)a $ 水平風 Horizontal wind
  real(8), allocatable :: pz_VelXNs(:,:) 
                            ! $ u (\tau) $ 水平風 Horizontal wind
  real(8), allocatable :: pz_VelXAs(:,:) 
                            ! $ u (\tau +\Delta \tau) 水平風 Horizontal wind
  real(8), allocatable :: xr_VelZBl(:,:) 
                            ! $ w (t-\Delta t) 鉛直風 Vertical wind
  real(8), allocatable :: xr_VelZNl(:,:) 
                            ! $ w (t) 鉛直風 Vertical wind
  real(8), allocatable :: xr_VelZAl(:,:) 
                            ! $ w (t+\Delta t) 鉛直風 Vertical wind
  real(8), allocatable :: xr_VelZNs(:,:) 
                            ! $ w (\tau) 鉛直風 Vertical wind
  real(8), allocatable :: xr_VelZAs(:,:) 
                            ! $ w (\tau +\Delta \tau) $ 鉛直風 Vertical wind
  real(8), allocatable :: xz_ExnerBl(:,:)
                            ! $ \pi (t-\Delta t) $ 圧力関数 Exner function
  real(8), allocatable :: xz_ExnerNl(:,:)
                            ! $ \pi (t) $ 圧力関数 Exner function
  real(8), allocatable :: xz_ExnerAl(:,:)
                            ! $ \pi (t+\Delta t) $ 圧力関数 Exner function
  real(8), allocatable :: xz_ExnerNs(:,:)
                            ! $ \pi (\tau) $ 圧力関数 Exner function
  real(8), allocatable :: xz_ExnerAs(:,:)
                            ! $ \pi (\tau +\Delta \tau) $ 圧力関数 Exner function
  real(8), allocatable :: xz_PotTempBl(:,:) 
                            ! $ \theta (t-\Delta t) $ 温位 Potential temp.
  real(8), allocatable :: xz_PotTempNl(:,:) 
                            ! $ \theta (t) $ 温位 Potential temp.
  real(8), allocatable :: xz_PotTempAl(:,:) 
                            ! $ \theta (t+\Delta t) $ 温位 Potential temp.
  real(8), allocatable :: xz_PotTempWork(:,:)
                            ! 温位の作業配列 Work array for potential temp.
  real(8), allocatable :: xz_KmBl(:,:)
                            ! $ Km (t-\Delta t) $ 乱流拡散係数 
                            ! Turbulent diffusion coeff. 
  real(8), allocatable :: xz_KmNl(:,:)         
                            ! $ K_m (t) 乱流拡散係数 
                            ! Turbulent diffusion coeff. 
  real(8), allocatable :: xz_KmAl(:,:)
                            ! $ K_m (t+\Delta t) $ 乱流拡散係数 
                            ! Turbulent diffusion coeff. 
  real(8), allocatable :: xz_KhBl(:,:)
                            ! $ K_h (t-\Delta t) $ 乱流拡散係数 
                            ! Turbulent diffusion coeff. 
  real(8), allocatable :: xz_KhNl(:,:)
                            ! $ K_h (t) $ 乱流拡散係数 
                            ! Turbulent diffusion coeff. 
  real(8), allocatable :: xz_KhAl(:,:)
                            ! $ K_h (t+\Delta t) $ 乱流拡散係数 
                            ! Turbulent diffusion coeff.
  real(8), allocatable :: xza_MixRtBl(:,:,:)
                            ! $ q (t-\Delta t) $ 湿潤量の混合比 
                            ! Mixing ratio of moist variables.
  real(8), allocatable :: xza_MixRtNl(:,:,:) 
                            ! $ q (t) $ 湿潤量の混合比 
                            ! Mixing ratio of moist variables.
  real(8), allocatable :: xza_MixRtAl(:,:,:)
                            ! $ q (t+\Delta t) $ 湿潤量の混合比 
                            ! Mixing ratio of moist variables.
  real(8), allocatable :: xza_MixRtWork(:,:,:) 
                            ! 湿潤量の作業配列
                            ! Work array for mixing ratio.
  real(8), allocatable :: pz_AccelVelXNl(:,:)  
                            ! 気圧傾度力を除いた $u$ の変化率
                            ! Tendency of $u$ except for pressure gradient term
  real(8), allocatable :: xr_AccelVelZNl(:,:)  
                            ! 気圧傾度力を除いた $w$ の変化率
                            ! Tendency of $w$ except for pressure gradient term
  real(8), allocatable :: xza_DelMixRt(:,:,:)  
                            ! 湿潤量の混合比 $ q $ の増分
                            ! Mixing ratio variation.
  real(8) :: Time           ! 時刻 Time 
  real(8) :: ReStartTime(2) ! リスタートファイル出力時刻用配列
                            ! Output time array for restart file
  real(8), allocatable :: DelTimeLFrog(:)
                            ! リープフロッグスキーム用時間格子間隔
                            ! Time interval for Leap-frog scheme
  real(8) :: DelTimeEular   ! オイラースキーム用時間格子
                            ! Time interval for Eular scheme
  integer :: NStepLFrog     ! リープフロッグスキーム用時間ステップ数
                            ! The number of time step for Leap-frog scheme
  integer, allocatable :: NStepEular(:)
                            ! オイラースキーム用時間ステップ数
                            ! The number of time step for Eular scheme
  integer :: t, tau, t1, t2       ! do ループ変数 do loop variable 
  ! 初期化手続き 
  ! Initialize procedure 
  !
  call MainInit
  ! 時間積分 time integration 
  !
  do t1 = 1, NstepLFrog / NstepDisp
    do t2 = 1, NstepDisp
    
      ! 時刻の設定 
      ! Time setting.
      !
      Time = Time + DelTimeLong
      t = (t1 - 1) * NstepDisp + t2
      ! 渦拡散係数の移流拡散
      ! Advection and diffusion of turbulent diffusion coefficient.
      !
      xz_KmAl = xz_KmBl + DelTimeLFrog(t) * ( + xz_AdvKm(xz_KmNl, pz_VelXNl, xr_VelZNl) + xz_BuoyMoistKm(xz_PotTempBl, xz_ExnerBl, xza_MixRtBl) + xz_ShearKm(xz_KmBl, pz_VelXBl, xr_VelZBl) + xz_NumDiffKm(xz_KmBl) + xz_DispKm(xz_KmBl) )
      
      ! 値の上限下限の設定
      !  * 値は正になることを保証する
      !  * 値の上限は 800 とする. 中島健介(1994, 学位論文)参照
      ! 
      ! Upper and lower bound value are specified.
      !
      xz_KmAl = max( 0.0d0, min( xz_KmAl, 800.0d0 ) )
      ! 境界条件 Boundary condition
      call BoundaryXCyc_xz( xz_KmAl ) ! (inout)
      call BoundaryZSym_xz( xz_KmAl ) ! (inout)
      
      ! スカラーに対する渦拡散係数の計算 
      ! Specify turbulent diffusion coefficient for scalar variables.
      !
      xz_KhAl = 3.0d0 * xz_KmAl
      
      ! 温位の移流拡散の計算 
      ! Advection and diffusion of potential temperature.
      !
      xz_PotTempAl = xz_PotTempBl + DelTimeLFrog(t) * ( + xz_AdvScalar( xz_PotTempNl,     pz_VelXNl, xr_VelZNl) + xz_AdvScalar( xz_PotTempBasicZ, pz_VelXNl, xr_VelZNl) + xz_TurbScalar(xz_PotTempBl,     xz_KhBl) + xz_TurbScalar(xz_PotTempBasicZ, xz_KhBl) + xz_NumDiffScalar(xz_PotTempBl) + xz_DispHeat( xz_KmBl ) + xz_RadHeatConst( xz_ExnerBl ) + xz_HeatFluxDiff( xz_PotTempNl ) )
      ! 境界条件 Boundary condition
      !
      call BoundaryXCyc_xz( xz_PotTempAl ) ! (inout)
      call BoundaryZSym_xz( xz_PotTempAl ) ! (inout)
      
      ! 凝縮成分混合比の移流拡散 
      ! Advection and diffusion of vapor, cloud and rain mixing ratios.
      !
      xza_MixRtAl = xza_MixRtBl + DelTimeLFrog(t) * ( + xza_AdvScalar(xza_MixRtNl,     pz_VelXNl, xr_VelZNl) + xza_AdvScalar(xza_MixRtBasicZ, pz_VelXNl, xr_VelZNl) + xza_TurbScalar(xza_MixRtBl,    xz_KhBl) + xza_TurbScalar(xza_MixRtBasicZ,xz_KhBl) + xza_NumDiffScalar(xza_MixRtBl) + xza_FallRain(xza_MixRtBl) + xza_MixRtFluxDiff(xza_MixRtNl) )
      ! 移流によって負になった部分を埋める
      ! Negative values due to advection are corrected.
      !
      xza_MixRtWork = xza_MixRtAl
      xza_MixRtAl = xza_FillNegative_xza( xza_MixRtWork ) 
      
      ! 埋めた/削った量を保管
      ! Correction value is stored.
      !
      call StoreMixRtFill1( (xza_MixRtAl - xza_MixRtWork) / DelTimeLFrog(t) )    
      ! 境界条件 Boundary condition      
      !
      call BoundaryXCyc_xza( xza_MixRtAl ) ! (inout)
      call BoundaryZSym_xza( xza_MixRtAl ) ! (inout)
      
      ! 暖かい雨のパラメタリゼーション.
      ! * 雲<-->雨 の変換を行う.
      !
      ! Warm rain parameterization.
      ! * Conversion from cloud to rain.
      ! これまでの値を作業配列に保管
      ! Previous values are stored to work area.
      !
      xza_MixRtWork = xza_MixRtAl
      
      ! 雨への変化量を計算
      ! Conversion values are calculated.
      !
      xza_MixRtAl   = xza_MixRtWork + xza_Cloud2Rain( xza_MixRtAl, DelTimeLFrog(t) )
      
      ! 雲から雨への変換量を保管
      ! Conversion values are sotred.
      !
      call StoreMixRtCond( (xza_MixRtAl - xza_MixRtWork) / DelTimeLFrog(t) ) 
      
      ! 境界条件 Boundary condition
      call BoundaryXCyc_xza( xza_MixRtAl ) ! (inout)
      call BoundaryZSym_xza( xza_MixRtAl ) ! (inout)
      
      ! 湿潤飽和調節
      ! * 水蒸気<-->雲の変換を行う.
      !
      ! Moist adjustment.
      ! * Conversion from vapor to cloud.
      ! これまでの値を作業配列に保管
      ! Previous values are stored to work area.
      !
      xz_PotTempWork = xz_PotTempAl
      xza_MixRtWork  = xza_MixRtAl
      
      ! 湿潤調節法を適用
      ! Moist adjustment is applied.
      !
      call MoistAdjustSvapPress(  xz_ExnerNl, xz_PotTempAl, xza_MixRtAl )
      call MoistAdjustNH4SH( xz_ExnerNl, xz_PotTempAl, xza_MixRtAl )
      
      ! 湿潤調節法による温位と混合比の変化量を保管
      ! Adjustment values are stored.
      !
      call StorePotTempCond( (xz_PotTempAl - xz_PotTempWork) / DelTimeLFrog(t) ) 
      call StoreMixRtCond( (xza_MixRtAl - xza_MixRtWork) / DelTimeLFrog(t) ) 
      
      ! 境界条件 Boundary condition
      !
      call BoundaryXCyc_xza( xza_MixRtAl ) ! (inout)
      call BoundaryZSym_xza( xza_MixRtAl ) ! (inout)
      ! 暖かい雨のパラメタリゼーション.
      ! * 蒸気<-->雨 の変換を行う
      !
      ! Warm rain parameterization.
      ! * Conversion from rain to vapor.
      ! これまでの値を作業配列に保管
      ! Previous values are stored to work area.
      !
      xz_PotTempWork = xz_PotTempAl
      xza_MixRtWork = xza_MixRtAl
      
      ! 雨から蒸気への混合比変化を求める
      ! * 温位の計算において, 混合比変化が必要となるため, 
      !   混合比変化を 1 つの配列として用意する.
      ! Conversion values are calculated.
      !
      xza_DelMixRt = 0.0d0
      xza_DelMixRt = ( + xza_Rain2Gas( xz_ExnerNl, xz_PotTempAl, xza_MixRtAl, DelTimeLFrog(t) ) + xza_Rain2GasNH4SH( xz_ExnerNl, xz_PotTempAl, xza_MixRtAl, DelTimeLFrog(t) ) )    
      
      ! 温位の計算. 雨から蒸気への変換に伴う潜熱・反応熱を追加.
      !
      xz_PotTempAl = xz_PotTempWork + ( + xz_Rain2GasHeat( xz_PotTempAl, xz_ExnerNl, xza_DelMixRt ) + xz_Rain2GasHeatNH4SH( xz_ExnerNl, xza_DelMixRt ) )
      
      ! 混合比の計算. 雨から蒸気への変換分を追加
      !
      xza_MixRtAl   = xza_MixRtWork + xza_DelMixRt
      
      ! 雨から蒸気の値の保管
      !
      call StorePotTempCond( (xz_PotTempAl - xz_PotTempWork) / DelTimeLFrog(t) ) 
      call StoreMixRtCond( xza_DelMixRt / DelTimeLFrog(t) ) 
      
      ! 境界条件 Boundary condition
      !
      call BoundaryXCyc_xz( xz_PotTempAl ) ! (inout)
      call BoundaryZSym_xz( xz_PotTempAl ) ! (inout)
      call BoundaryXCyc_xza( xza_MixRtAl ) ! (inout)
      call BoundaryZSym_xza( xza_MixRtAl ) ! (inout)
      
      ! 速度の移流拡散.
      ! Advection and diffusion of velocity components.
      !
      pz_AccelVelXNl = + pz_AdvVelX(pz_VelXNl, xr_VelZNl) + pz_TurbVelX(xz_KmBl,   pz_VelXBl, xr_VelZBl) + pz_NumDiffVelX(pz_VelXBl)
      
      xr_AccelVelZNl = + xr_AdvVelZ(xr_VelZNl, pz_VelXNl) + xr_Buoy(xz_PotTempNl) + xr_BuoyMolWt(xza_MixRtNl) + xr_BuoyDrag(xza_MixRtNl) + xr_TurbVelZ(xz_KmBl,   pz_VelXBl, xr_VelZBl) + xr_NumDiffVelZ(xr_VelZBl)                        
      ! 短い時間ステップの初期値作成.
      ! Initial values set up for time integration with short time step.
      !
      pz_VelXNs  = pz_VelXBl
      xr_VelZNs  = xr_VelZBl
      xz_ExnerNs = xz_ExnerBl
      
      ! 短い時間ステップの時間積分. オイラー法を利用.
      ! Time integration with short time step.
      !
      Eular: do tau = 1, NstepEular(t)
        
        ! 速度 u の計算.
        ! Time integration horizontal velocity.
        !
        pz_VelXAs = pz_VelXNs + DelTimeEular * ( - pz_GradPi(xz_ExnerNs, pz_VelXNs, xr_VelZNs) + pz_AccelVelXNl )
        
        ! 境界条件 Boundary conditions.
        ! 
        call BoundaryXCyc_pz( pz_VelXAs ) ! (inout)
        call BoundaryZSym_pz( pz_VelXAs ) ! (inout)
        
        ! エクスナー関数の計算.
        ! Time integration exner function.
        !
        xz_ExnerAs = xz_Exner( xr_AccelVelZNl, pz_VelXNs, pz_VelXAs, xr_VelZNs, xz_ExnerNs)
        
        ! 境界条件  Boundary condition
        !
        call BoundaryXCyc_xz( xz_ExnerAs ) ! (inout)
        call BoundaryZSym_xz( xz_ExnerAs ) ! (inout)
        
        ! 速度 w の計算
        ! Time integration vertical velocity.
        !
        xr_VelZAs = xr_VelZNs + DelTimeEular * ( - xr_GradPi(xz_ExnerAs,xz_ExnerNs,pz_VelXNs,xr_VelZNs) + xr_AccelVelZNl )
        
        ! 境界条件 ; Boundary condition
        !
        call BoundaryXCyc_xr( xr_VelZAs )     ! (inout)
        call BoundaryZAntiSym_xr( xr_VelZAs ) ! (inout)
        
        ! 短い時間ステップのループを回すための処置
        ! Renew prognostic variables for next short time step integration.
        !
        xz_ExnerNs  = xz_ExnerAs
        pz_VelXNs   = pz_VelXAs
        xr_VelZNs   = xr_VelZAs
        
      end do Eular
      
      ! 最終的な短い時間ステップでの値を長い時間ステップでの値とみなす
      ! Renew prognostic variables for next long time step integration.
      !
      xz_ExnerAl  = xz_ExnerAs
      pz_VelXAl   = pz_VelXAs
      xr_VelZAl   = xr_VelZAs
      
      ! 時間フィルタ. 
      ! Time filter. 
      !
      call AsselinFilter_xz( xz_ExnerAl, xz_ExnerNl, xz_ExnerBl )         
      call AsselinFilter_pz( pz_VelXAl, pz_VelXNl, pz_VelXBl )
      call AsselinFilter_xr( xr_VelZAl, xr_VelZNl, xr_VelZBl )
      call AsselinFilter_xz( xz_PotTempAl, xz_PotTempNl, xz_PotTempBl )
      call AsselinFilter_xz( xz_KmAl, xz_KmNl, xz_KmBl )
      call AsselinFilter_xza( xza_MixRtAl, xza_MixRtNl, xza_MixRtBl )
    
      ! スポンジ層.
      ! Numerical dumping.
      !
      call DampSponge_pz( pz_VelXAl, pz_VelXBl, DelTimeLFrog(t) )
      call DampSponge_xr( xr_VelZAl, xr_VelZBl, DelTimeLFrog(t) )
      ! 混合比がゼロ以下にならないための処置.
      ! Negative values mixing ratios are corrected. 
      !
      xza_MixRtWork = xza_MixRtAl 
      xza_MixRtAl = max( - xza_MixRtBasicZ, xza_MixRtWork )
      call StoreMixRtFill2( (xza_MixRtAl - xza_MixRtWork) / DelTimeLFrog(t) )
      ! 安定度の計算.
      ! Calculation static stability.
      !
      call ECCM_Stab( xz_PotTempAl, xz_ExnerAl, xza_MixRtAl )
      
      ! 長い時間ステップのループを回すための処置.
      ! Renew prognostic variables for next long time step integration.
      !
      xz_PotTempBl = xz_PotTempNl
      xza_MixRtBl  = xza_MixRtNl
      xz_ExnerBl   = xz_ExnerNl
      pz_VelXBl    = pz_VelXNl
      xr_VelZBl    = xr_VelZNl
      xz_KmBl      = xz_KmNl
      xz_KhBl      = xz_KhNl
      
      xz_PotTempNl = xz_PotTempAl
      xza_MixRtNl  = xza_MixRtAl
      xz_ExnerNl   = xz_ExnerAl
      pz_VelXNl    = pz_VelXAl
      xr_VelZNl    = xr_VelZAl
      xz_KmNl      = xz_KmAl
      xz_KhNl      = xz_KhAl
      
    end do
    
    ! ヒストリーファイルへの出力.
    ! Out put to history file.
    !
    call MessageNotify( "M", "main", "Time = %f", d=(/Time/) )
    ! 移流に対する CFL 条件のチェック
    ! CFL condtion check for advection
    call CFLCheckTimeLongVelX( pz_VelXNl )
    call CFLCheckTimeLongVelZ( xr_VelZNl )
    ! ヒストリファイル出力
    ! Out put to history file.
    !
    call HistoryFile_Output( Time, xz_PotTempNl, xz_ExnerNl, pz_VelXNl, xr_VelZNl, xza_MixRtNl, xz_KmNl, xz_KhNl )
    
    ! 積算値のクリア
    ! Clear monitor variables.
    !
    call StorePotTempClean
    call StoreMixRtClean
    call StoreBuoyClean
    call StoreStabClean
    
  end do
  
  ! 出力ファイルのクローズ
  ! Close out put files.
  !
  call HistoryFile_Close
  ! リスタートファイルの作成
  ! Make restartfile.
  !
  call ReStartFile_Open( )
  call ReStartFile_OutPut( Time - DelTimeLong, xz_PotTempBl, xz_ExnerBl, pz_VelXBl, xr_VelZBl, xza_MixRtBl, xz_KmBl, xz_KhBl )
  call ReStartFile_OutPut( Time, xz_PotTempNl, xz_ExnerNl, pz_VelXNl, xr_VelZNl, xza_MixRtNl, xz_KmNl, xz_KhNl )
  call ReStartFile_Close
  
contains
!-----------------------------------------------------------------------
  subroutine MainInit
    ! NAMELIST ファイル名の読み込み
    ! Loading NAMELIST file.
    !
    call argset_init( cfgfile )
    
    ! デバッグ設定の初期化
    ! Initialization of debug output control.
    !
    call debugset_init( cfgfile )
    
    ! 化学定数の初期化
    ! Initialization of chemical constatns.
    !
    call chemdata_init()
    
    ! 時間積分の初期化
    ! Initialization of time integration.
    !
    call timeset_init( cfgfile )
    
    ! 格子点情報の初期化
    ! Initialization of grid arrangement.
    !
    call gridset_init( cfgfile )
    
    ! 化学計算ルーチンの初期化
    ! Initialization of chemical routines.
    !
    call chemcalc_init()
    
    ! 基本場変数の初期化
    ! Initialization of basic state variables.
    !
    call basicset_init( cfgfile )
    
    ! I/O ファイル名の初期化
    ! Initialization of output file name. 
    !
    call fileset_init( cfgfile )
    
    ! 湿潤過程共有変数の初期化
    ! Initialization of common variables for moist process.
    !
    call moistset_init( )
    
    ! 積算値保管変数の初期化
    ! Initialization of monitor variables.
    !
    call StorePotTemp_init( )
    call StoreMixRt_init( )
    call StoreBuoy_init( )  
    call StoreStab_init( )  
    
    ! 内部変数の初期化
    ! Initialization of internal variables.
    !
    allocate( pz_VelXBl(DimXMin:DimXMax, DimZMin:DimZMax), pz_VelXNl(DimXMin:DimXMax, DimZMin:DimZMax), pz_VelXAl(DimXMin:DimXMax, DimZMin:DimZMax), pz_VelXNs(DimXMin:DimXMax, DimZMin:DimZMax), pz_VelXAs(DimXMin:DimXMax, DimZMin:DimZMax), xr_VelZBl(DimXMin:DimXMax, DimZMin:DimZMax), xr_VelZNl(DimXMin:DimXMax, DimZMin:DimZMax), xr_VelZAl(DimXMin:DimXMax, DimZMin:DimZMax), xr_VelZNs(DimXMin:DimXMax, DimZMin:DimZMax), xr_VelZAs(DimXMin:DimXMax, DimZMin:DimZMax), xz_ExnerBl(DimXMin:DimXMax, DimZMin:DimZMax), xz_ExnerNl(DimXMin:DimXMax, DimZMin:DimZMax), xz_ExnerAl(DimXMin:DimXMax, DimZMin:DimZMax), xz_ExnerNs(DimXMin:DimXMax, DimZMin:DimZMax), xz_ExnerAs(DimXMin:DimXMax, DimZMin:DimZMax), xz_PotTempWork(DimXMin:DimXMax, DimZMin:DimZMax), xz_PotTempBl(DimXMin:DimXMax, DimZMin:DimZMax), xz_PotTempNl(DimXMin:DimXMax, DimZMin:DimZMax), xz_PotTempAl(DimXMin:DimXMax, DimZMin:DimZMax), xz_KmBl(DimXMin:DimXMax, DimZMin:DimZMax), xz_KmNl(DimXMin:DimXMax, DimZMin:DimZMax), xz_KmAl(DimXMin:DimXMax, DimZMin:DimZMax), xz_KhBl(DimXMin:DimXMax, DimZMin:DimZMax), xz_KhNl(DimXMin:DimXMax, DimZMin:DimZMax), xz_KhAl(DimXMin:DimXMax, DimZMin:DimZMax), xza_MixRtWork(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum), xza_MixRtBl(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum), xza_MixRtNl(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum), xza_MixRtAl(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum), pz_AccelVelXNl(DimXMin:DimXMax, DimZMin:DimZMax), xr_AccelVelZNl(DimXMin:DimXMax, DimZMin:DimZMax), xza_DelMixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum), DelTimeLFrog(NstepLong), NStepEular(NStepLong) )
    
    pz_VelXNs = 0.0d0
         pz_VelXAs = 0.0d0    
    xr_VelZNs = 0.0d0
         xr_VelZAs = 0.0d0
    xz_ExnerNs = 0.0d0
        xz_ExnerAs = 0.0d0
    pz_VelXBl = 0.0d0
         pz_VelXNl = 0.0d0
         pz_VelXAl = 0.0d0
    xr_VelZBl = 0.0d0
         xr_VelZNl = 0.0d0
         xr_VelZAl = 0.0d0
    xz_ExnerBl = 0.0d0
        xz_ExnerNl = 0.0d0
        xz_ExnerAl = 0.0d0
    xz_KmBl = 0.0d0
           xz_KmNl = 0.0d0
           xz_KmAl = 0.0d0
    xz_KhBl = 0.0d0
           xz_KhNl = 0.0d0
           xz_KhAl = 0.0d0
    xz_PotTempBl = 0.0d0
      xz_PotTempNl = 0.0d0
      xz_PotTempAl = 0.0d0
    xza_MixRtBl = 0.0d0
        xza_MixRtNl = 0.0d0
        xza_MixRtAl = 0.0d0
    pz_AccelVelXNl = 0.0d0
    xr_AccelVelZNl = 0.0d0
    xza_DelMixRt = 0.0d0
    DelTimeEular = 0.0d0
    DelTimeLFrog = 0.0d0
     
    NStepEular = 0.0d0
    NstepLFrog = 0.0d0
    
    ! 初期値の代入 
    ! * ReStartFile が設定されている場合にはファイルを読み込む. 
    !   設定されていない場合にはデフォルトの基本場と擾乱場を作る. 
    !
    ! Initial value set up.
    ! * Read restartfile if it is specified. If not, make default basic
    !   state and disturbance.
    !
    call MessageNotify( "M", "main", "Initial value setup." )
    
    if (trim(InitFile) /= '') then    
       
       call MessageNotify( "M", "main", "Restart file is %c", c1=trim(Initfile) )
       
       call ReStartFile_Get( ReStartTime, xz_PotTempBl, xz_ExnerBl, pz_VelXBl, xr_VelZBl, xza_MixRtBl, xz_KmBl, xz_KhBl, xz_PotTempNl, xz_ExnerNl, pz_VelXNl, xr_VelZNl, xza_MixRtNl, xz_KmNl, xz_KhNl )
    else
       call MessageNotify( "W", "main", "Restart file is not specified." )
       call BasicEnv( )
       call DisturbEnv( cfgfile, xz_PotTempBl, xz_ExnerBl, pz_VelXBl, xr_VelZBl, xza_MixRtBl, xz_KmBl, xz_KhBl )
       
       call BoundaryXCyc_pz( pz_VelXBl )     ! (inout)
       call BoundaryZSym_pz( pz_VelXBl )     ! (inout)
       
       call BoundaryXCyc_xr( xr_VelZBl )     ! (inout)
       call BoundaryZAntiSym_xr( xr_VelZBl ) ! (inout)
       
       call BoundaryXCyc_xz( xz_PotTempBl )  ! (inout)
       call BoundaryZSym_xz( xz_PotTempBl )  ! (inout)
       
       call BoundaryXCyc_xza( xza_MixRtBl )  ! (inout)
       call BoundaryZSym_xza( xza_MixRtBl )  ! (inout)
       
       
       ! 時刻 $ t $ の変数値の初期値の設定
       ! * 1 ループ目だけは $ t $ の値を $ t-\Delta t$ の値と同じにする. 
       !   1 ステップ目はオイラー法で解く必要があるが, 1 ステップ目と
       !   それ以外のステップを別々にコーディングしたくない為
       !
       ! Set up initial value of time = "t" variables.
       !  
       xz_PotTempNl = xz_PotTempBl
       xz_ExnerNl   = xz_ExnerBl
       pz_VelXNl    = pz_VelXBl
       xr_VelZNl    = xr_VelZBl
       xza_MixRtNl  = xza_MixRtBl
       xz_KmNl      = xz_KmBl
       xz_KhNl      = xz_KhBl
       
    end if
    
    ! 数値摩擦係数の初期化
    ! Initialization of numerical friction coefficient.
    !
    call Damping_Init( cfgfile )      
    
    ! 数値拡散項の初期化
    ! Initialization of numerical diffusion term.
    !
    call NumDiffusion_Init( )          
    
    ! 乱流拡散項の初期化
    ! Initialization of turbulent diffusion term.
    !
    call Turbulence_Init( )            
    
    ! 暖かい雨のパラメタリゼーションの初期化 
    ! Initialization of warmrain parameterization.
    !
    call WarmRainPrm_Init( cfgfile )  
    
    ! 負の湿潤量の補填計算の初期化
    ! Initialization of negative moist value correction.
    !
    call FillNegative_Init( xza_MixRtBasicZ, xz_DensBasicZ )
    
    ! 放射強制の初期化
    !  Initialization of radiative forcing.
    !
    call Radiation_Init( cfgfile )
    
    ! 湿潤気塊の浮力計算の初期化
    ! Initialization of moist buoyancy calculation.
    !
    call MoistBuoy_Init( ) 
    
    ! 圧力計算用係数行列の初期化 
    ! Initialization of coefficient matrix for exner function calculation.
    !
    call xz_Exner_Init( )              
    
    ! 時刻とループ回数の初期化
    ! Initialization of time integration.
    !
    NstepLFrog   = NstepLong 
    NstepEular   = NstepShort 
    DelTimeLFrog = DelTimeLong * 2.0d0 
    DelTimeEular = DelTimeShort
    
    ! 計算開始時刻と時間格子間隔の初期化
    ! * ReStartFile が設定されている場合, ファイルから読み込んだ値を指定.
    ! * ReStartFile が設定されてない場合
    !   * 開始時刻は 0.0
    !   * 1 ステップ目の時間格子間隔を陽に指定
    !
    ! Setup restart time and time interval. 
    ! * Read restartfile if it is specified.
    ! * If not, 
    !   * "t" is set to be 0.
    !   * Time intervals for 1st step are specified explicitly.
    !
    if ( trim(InitFile) /= '') then    
       Time = ReStartTime(2)               
    else
       Time = 0.0d0                    
       NstepEular(1)   = NstepShort /2 
       DelTimeLFrog(1) = DelTimeLong   
    end if
    
    ! ヒストリーファイルへの出力
    ! Out put to history file.
    !
    call HistoryFile_Open( )
    
    if ( trim(InitFile) /= '') then    
       call HistoryFile_Output( ReStartTime(2), xz_PotTempNl, xz_ExnerNl, pz_VelXNl, xr_VelZNl, xza_MixRtNl, xz_KmNl, xz_KhNl )
    end if
    
    if ( Time == 0 ) then
       call HistoryFile_Output( Time, xz_PotTempNl, xz_ExnerNl, pz_VelXNl, xr_VelZNl, xza_MixRtNl, xz_KmNl, xz_KhNl )
    end if
    
    ! 音波に対する CFL 条件のチェック
    ! CFL condtion check for sound wave.
    !
    call CFLCheckTimeShort( xz_VelSoundBasicZ )
    
end subroutine MainInit
!-----------------------------------------------------------------------
end program arare