program deepconv_arare
!
! 非静力学モデル deepconv/arare 湿潤大気対流計算用主プログラム (三次元版)
!
! モジュール引用 use statement
!
! gtool5 関連
! gtool5 modules
!
use dc_types, only: STRING, DP
use dc_message, only: MessageNotify
use gtool_history, only: HistoryPut
use gtool_historyauto, only: HistoryAutoPut
! 初期設定モジュール
! Initialize module
!
use mpi_wrapper, only : MPIWrapperInit, MPIWrapperFinalize, myrank
use argset, only: argset_init
use clockset, only: ClocksetInit, ClocksetClose, ClocksetPredict, ClockSetPreStart, ClockSetLoopStart, ClockSetPreStop, ClockSetLoopStop
use gridset, only: gridset_init, imin, imax, jmin, jmax, kmin, kmax, nx, ny, nz, ncmax
use timeset, only: timeset_init, TimesetDelTimeHalf, TimesetProgress, TimeA, TimeN, TimeB, NstepShort, NstepOutput, EndTime, RestartTime, DelTimeLong, DelTimeShort, DelTimeOutput, FlagInitialRun
use axesset, only: axesset_init, xyz_avr_pyz, xyz_avr_xqz, xyz_avr_xyr
use constants,only: constants_init
use composition, only: composition_init, SpcWetSymbol
use fileset, only: fileset_init
use basicset, only: xyzf_QMixBZ, xyz_DensBZ, xyz_EffMolWtBZ, xyz_PTempBZ, xyz_TempBZ, xyz_PressBZ, xyz_VelSoundBZ, xyz_ExnerBZ
use ChemCalc, only: ChemCalc_init
use namelist_util, only: NmlutilInit, namelist_filename
! 下請けモジュール
! Utility modules
!
use cflcheck, only : CFLCheckTimeShort, CFLCheckTimeLongVelX, CFLCheckTimeLongVelY, CFLCheckTimeLongVelZ
use damping, only: Damping_init, SpongeLayer_forcing
! 力学過程計算用関数モジュール
! Dynamical processes module
!
use DynamicsHEVI, only: Dynamics_Init, Dynamics_Long_forcing, Dynamics_Short_forcing, Dynamics_Km_forcing
use fillnegative,only: FillNegative_init
! 乱流拡散計算用モジュール
! Turbulent diffusion module
!
use Turbulence_kw1978, only: Turbulence_kw1978_Init, Turbulence_KW1978_forcing
! 放射強制計算用モジュール
! Radiative forceing module
!
use Radiation_simple, only: Radiation_simple_init, xyz_PTempRadConst, xyz_PTempRadVary, xyz_ExnerRadConst, xyz_ExnerRadVary
use radiation_heatbalance, only: Radiation_heatbalance_init, Radiation_heatbalance_forcing
! 境界からのフラックス計算用モジュール
! Surface flux module
!
use Surfaceflux_Diff, only: Surfaceflux_Diff_init, Surfaceflux_Diff_forcing
! use Surfaceflux_Bulk_l1982, only: Surfaceflux_Bulk_init, Surfaceflux_Bulk_forcing
use Surfaceflux_Bulk, only: Surfaceflux_Bulk_init, Surfaceflux_Bulk_forcing
use Surfaceflux_Const, only: Surfaceflux_Const_init, Surfaceflux_Const_forcing
! 湿潤過程計算用モジュール
! Moist processes modules
!
use Cloudphys_K1969, only: Cloudphys_K1969_Init, Cloudphys_K1969_forcing, Cloudphys_K1969_FallRain
use Cloudphys_marscond, only: cloudphys_marscond_Init, cloudphys_marscond_forcing
! ファイル入出力モジュール
! File I/O module
!
use RestartFileIO, only : ReStartFileio_init, ReStartFileio_Finalize, ReStartFileio_BZ_Get, ReStartFileio_Var_Get, rstat
use HistoryFileIO, only: HistoryFileio_init, HistoryFileio_Finalize
! 暗黙の型宣言禁止
!
implicit none
! 内部変数
! Internal variables
!
character(*), parameter:: prog_name = 'arare'
! 主プログラム名.
! Main program name
real(DP), allocatable :: pyz_VelXBl(:,:,:)
! $ u (t-\Delta t) $ 東西風 ; zonal wind
real(DP), allocatable :: pyz_VelXNl(:,:,:)
! $ u (t) $ 東西風 ; zonal wind
real(DP), allocatable :: xyz_VelXNl(:,:,:)
! $ u (t) $ 東西風 ; zonal wind
real(DP), allocatable :: pyz_VelXAl(:,:,:)
! $ u (t+\Delta t) $ 東西風 ; zonal wind
real(DP), allocatable :: pyz_VelXNs(:,:,:)
! $ u (\tau) $ 東西風 ; zonal wind
real(DP), allocatable :: pyz_VelXAs(:,:,:)
! $ u (\tau +\Delta \tau) $ 東西風 ; zonal wind
real(DP), allocatable :: xqz_VelYBl(:,:,:)
! $ v (t-\Delta t) $ 南北風 ; meridonal wind
real(DP), allocatable :: xqz_VelYNl(:,:,:)
! $ v (t) $ 南北風 ; meridonal wind
real(DP), allocatable :: xyz_VelYNl(:,:,:)
! $ v (t) $ 南北風 ; meridonal wind
real(DP), allocatable :: xqz_VelYAl(:,:,:)
! $ v (t+\Delta t) $ 南北風 ; meridonal wind
real(DP), allocatable :: xqz_VelYNs(:,:,:)
! $ v (\tau -\tau) $ 南北風 ; meridonal wind
real(DP), allocatable :: xqz_VelYAs(:,:,:)
! $ v (\tau) $ 南北風 ; meridonal wind
real(DP), allocatable :: xyr_VelZBl(:,:,:)
! $ w (t-\Delta t) $ 鉛直風 ; vertical wind
real(DP), allocatable :: xyr_VelZNl(:,:,:)
! $ w (t) $ 鉛直風 ; vertical wind
real(DP), allocatable :: xyz_VelZNl(:,:,:)
! $ w (t) $ 鉛直風 ; vertical wind
real(DP), allocatable :: xyr_VelZAl(:,:,:)
! $ w (t+\Delta t) $ 鉛直風 ; vertical wind
real(DP), allocatable :: xyr_VelZNs(:,:,:)
! $ w (\tau) $ 鉛直風 ; vertical wind
real(DP), allocatable :: xyr_VelZAs(:,:,:)
! $ w (\tau +\Delta \tau) 鉛直風 ; vertical wind
real(DP), allocatable :: xyz_ExnerBl(:,:,:)
! $ \pi (t-\Delta t) $ 圧力関数 ; Exner function
real(DP), allocatable :: xyz_ExnerNl(:,:,:)
! $ \pi (t) $ 圧力関数 ; Exner function
real(DP), allocatable :: xyz_ExnerAl(:,:,:)
! $ \pi (t+\Delta t) $ 圧力関数 ; Exner function
real(DP), allocatable :: xyz_ExnerNs(:,:,:)
! $ \pi (\tau -\Delta \tau) $ 圧力関数 ; Exner function
real(DP), allocatable :: xyz_ExnerAs(:,:,:)
! $ \pi (\tau) $ 圧力関数 ; Exner function
real(DP), allocatable :: xyz_PTempBl(:,:,:)
! $ \theta (t-\Delta t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_PTempNl(:,:,:)
! $ \theta (t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_PTempAl(:,:,:)
! $ \theta (t+\Delta t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_PTempNs(:,:,:)
! $ \theta (t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_PTempAs(:,:,:)
! $ \theta (t+\Delta t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_CDensBl(:,:,:)
! $ \theta (t-\Delta t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_CDensNl(:,:,:)
! $ \theta (t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_CDensAl(:,:,:)
! $ \theta (t+\Delta t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_CDensNs(:,:,:)
! $ \theta (t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_CDensAs(:,:,:)
! $ \theta (t+\Delta t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_KmBl(:,:,:)
! $ Km (t-\Delta t) $ 乱流拡散係数
! Turbulent diffusion coeff.
real(DP), allocatable :: xyz_KmNl(:,:,:)
! $ K_m (t) $ 乱流拡散係数
! Turbulent diffusion coeff.
real(DP), allocatable :: xyz_KmAl(:,:,:)
! $ K_m (t+\Delta t) $ 乱流拡散係数
! Turbulent diffusion coeff.
real(DP), allocatable :: xyz_KhBl(:,:,:)
! $ K_h (t-\Delta t) $ 乱流拡散係数
! Turbulent diffusion coeff.
real(DP), allocatable :: xyz_KhNl(:,:,:)
! $ K_h (t) $ 乱流拡散係数
! Turbulent diffusion coeff.
real(DP), allocatable :: xyz_KhAl(:,:,:)
! $ K_h (t+\Delta t) $ 乱流拡散係数
! Turbulent diffusion coeff.
real(DP), allocatable :: xyzf_QMixBl(:,:,:,:)
! $ q (t-\Delta t) $ 湿潤量の混合比
! Mixing ratio of moist variables.
real(DP), allocatable :: xyzf_QMixNl(:,:,:,:)
! $ q (t) $ 湿潤量の混合比
! Mixing ratio of moist variables
real(DP), allocatable :: xyzf_QMixAl(:,:,:,:) !
! $ q (t+\Delta t) $ 湿潤量の混合比
!Mixing ratio of moist variables
real(DP), allocatable :: pyz_DVelXDtNl(:,:,:)
real(DP), allocatable :: xqz_DVelYDtNl(:,:,:)
real(DP), allocatable :: xyr_DVelZDtNl(:,:,:)
real(DP), allocatable :: xyz_DPTempDtNl(:,:,:)
real(DP), allocatable :: xyz_DExnerDtNl(:,:,:)
real(DP), allocatable :: xyz_DExnerDtNs(:,:,:)
real(DP), allocatable :: xyzf_DQMixDtNl(:,:,:,:)
real(DP), allocatable :: xyz_DKmDtNl(:,:,:)
real(DP), allocatable :: xyz_DCDensDtNl(:,:,:)
integer :: s, t, tau ! do ループ変数 ; do loop variable
integer :: IDTurbMethod = 0
integer, parameter:: IDTurbMethodKW1978 = 2
integer :: IDRadMethod = 0
integer, parameter:: IDRadMethodHeatConst = 1
integer, parameter:: IDRadMethodHeatVary = 2
integer, parameter:: IDRadMethodHeatBalance = 3
integer :: IDSurfaceMethod = 0
integer, parameter:: IDSurfaceMethodDiff = 1
integer, parameter:: IDSurfaceMethodBulk = 2
integer, parameter:: IDSurfaceMethodConst = 3
integer :: IDCloudMethod = 0
integer, parameter:: IDCloudMethodK1969 = 1
integer, parameter:: IDCloudMethodMarsCond = 2
!------------------------------------------
! 初期化手続き ; Initialize procedure
!
call MainInit
!------------------------------------------
! 時間積分 time integration
!
if (myrank == 0) call MessageNotify( "M", "main", "Time Integration Start" )
! CPU 時間計測開始
! Start CPU time counting
!
call ClocksetLoopStart
! 時間発展ループのスタート
!
do while (TimeN <= EndTime )
!-------------------------------
! 物理過程: 乱流
!
select case ( IDTurbMethod )
case ( IDTurbMethodKW1978 )
! Km の移流
!
call Dynamics_Km_forcing( pyz_VelXNl, xqz_VelYNl, xyr_VelZNl, xyz_KmBl, xyz_KmNl, xyz_DKmDtNl )
call turbulence_KW1978_forcing( pyz_VelXBl, xqz_VelYBl, xyr_VelZBl, xyz_PTempBl, xyz_ExnerBl, xyzf_QMixBl, xyz_KmBl, xyz_KhBl, xyz_CDensBl, pyz_DVelXDtNl, xqz_DVelYDtNl, xyr_DVelZDtNl, xyz_DPTempDtNl,xyz_DExnerDtNl, xyzf_DQMixDtNl, xyz_DKmDtNl, xyz_DCDensDtNl, xyz_KmAl, xyz_KhAl )
end select
!-------------------------------
! 物理過程: 放射
!
select case (IDRadMethod)
case (IDRadMethodHeatConst)
xyz_DPTempDtNl = xyz_DPTempDtNl + xyz_PTempRadConst
xyz_DExnerDtNl = xyz_DExnerDtNl + xyz_ExnerRadConst
call HistoryAutoPut(TimeN, 'PTempRad', xyz_PTempRadConst(1:nx,1:ny,1:nz))
call HistoryAutoPut(TimeN, 'ExnerRad', xyz_ExnerRadConst(1:nx,1:ny,1:nz))
case (IDRadMethodHeatVary)
xyz_DPTempDtNl = xyz_DPTempDtNl + xyz_PTempRadVary
xyz_DExnerDtNl = xyz_DExnerDtNl + xyz_ExnerRadVary
call HistoryAutoPut(TimeN, 'PTempRad', xyz_PTempRadVary(1:nx,1:ny,1:nz))
call HistoryAutoPut(TimeN, 'ExnerRad', xyz_ExnerRadVary(1:nx,1:ny,1:nz))
case (IDRadMethodHeatBalance)
call radiation_heatbalance_forcing( xyz_ExnerNl, xyz_PTempNl, xyz_DPTempDtNl, xyz_DExnerDtNl )
end select
!--------------------------------
! 境界からの熱・運動量輸送
!
select case (IDSurfaceMethod)
case (IDSurfaceMethodDiff)
call Surfaceflux_Diff_forcing( xyz_PTempNl, xyzf_QMixNl, xyz_DPTempDtNl, xyz_DExnerDtNl, xyzf_DQMixDtNl )
case (IDSurfaceMethodBulk)
call Surfaceflux_Bulk_forcing( pyz_VelXNl, xqz_VelYNl, xyz_PTempNl, xyz_ExnerNl, xyzf_QMixNl, pyz_DVelXDtNl, xqz_DVelYDtNl, xyz_DPTempDtNl, xyz_DExnerDtNl, xyzf_DQMixDtNl )
case (IDSurfaceMethodConst)
call Surfaceflux_Const_forcing( pyz_DVelXDtNl, xqz_DVelYDtNl, xyz_DPTempDtNl, xyz_DExnerDtNl, xyzf_DQMixDtNl )
end select
!-----------------------------------------
! 凝結過程
!
!
select case (IDCloudMethod)
case (IDCloudMethodK1969)
call CloudPhys_K1969_FallRain( xyzf_QMixNl, xyz_DExnerDtNl, xyzf_DQMixDtNl )
end select
!------------------------------------------
! スポンジ層; sponge layer
!
!!!!!!!!!!!!!!!!!!!!!!!
!
! test
!
!!!!!!!!!!!!!!!!!!!!!!!
call SpongeLayer_forcing( pyz_VelXBl, xqz_VelYBl, xyr_VelZBl, xyz_PTempBl, xyz_ExnerBl, pyz_DVelXDtNl, xqz_DVelYDtNl, xyr_DVelZDtNl, xyz_DPTempDtNl, xyz_DExnerDtNl )
!-----------------------------------------
! 移流拡散.
! Advection and diffusion
!
call Dynamics_Long_forcing( pyz_VelXBl, pyz_VelXNl, xqz_VelYBl, xqz_VelYNl, xyr_VelZBl, xyr_VelZNl, xyz_PTempBl, xyz_PTempNl, xyzf_QMixBl, xyzf_QMixNl, xyz_CDensBl, xyz_CDensNl, pyz_DVelXDtNl, xqz_DVelYDtNl, xyr_DVelZDtNl, xyz_DPTempDtNl, xyzf_DQMixDtNl, xyz_DCDensDtNl, xyz_PTempAl, xyzf_QMixAl )
!------------------------------------------
! 凝結過程
!
select case (IDCloudMethod)
case (IDCloudMethodK1969)
call Cloudphys_K1969_forcing( xyz_ExnerNl, xyz_DExnerDtNl, xyz_PTempAl, xyzf_QMixAl )
end select
! 短い時間ステップの初期値作成.
! Initial values set up for time integration with short time step.
!
pyz_VelXNs = pyz_VelXBl
xqz_VelYNs = xqz_VelYBl
xyr_VelZNs = xyr_VelZBl
xyz_ExnerNs = xyz_ExnerBl
xyz_PTempNs = xyz_PTempBl
xyz_CDensNs = xyz_CDensBl
! 短い時間ステップの時間積分. オイラー法を利用.
! Time integration with short time step.
!
! 圧力関数への強制は 0 にする
!
! xyz_DExnerDtNl = 0.0d0
Euler: do tau = 1, NstepShort
! 火星計算の場合. 凝結量の評価はここで行う.
!
select case (IDCloudMethod)
case (IDCloudMethodMarsCond)
call cloudphys_marscond_forcing( xyz_PTempNs, xyz_ExnerNs, xyz_CDensNs, xyz_DPTempDtNl, xyz_DExnerDtNl, xyz_DCDensDtNl, xyz_PTempAs, xyz_CDensAs, xyz_DExnerDtNs )
end select
! 陽解法: 速度 u, v の計算.
! Time integration horizontal velocity (u).
!
call Dynamics_Short_forcing( pyz_VelXNs, xqz_VelYNs, xyr_VelZNs, xyz_ExnerNs, pyz_DVelXDtNl, xqz_DVelYDtNl, xyr_DVelZDtNl, xyz_DExnerDtNl, xyz_DExnerDtNs, pyz_VelXAs, xqz_VelYAs, xyr_VelZAs, xyz_ExnerAs )
! 短い時間ステップのループを回すための処置
! Renew prognostic variables for next short time step integration.
!
xyz_ExnerNs = xyz_ExnerAs
pyz_VelXNs = pyz_VelXAs
xqz_VelYNs = xqz_VelYAs
xyr_VelZNs = xyr_VelZAs
xyz_PTempNs = xyz_PTempAs
xyz_CDensNs = xyz_CDensAs
!
! clear tendency
!
xyz_DExnerDtNs = 0.0d0
end do Euler
! 最終的な短い時間ステップでの値を長い時間ステップでの値とみなす
! Renew prognostic variables for next long time step integration.
!
xyz_ExnerAl = xyz_ExnerAs
pyz_VelXAl = pyz_VelXAs
xqz_VelYAl = xqz_VelYAs
xyr_VelZAl = xyr_VelZAs
select case (IDCloudMethod)
case (IDCloudMethodMarsCond)
xyz_PTempAl = xyz_PTempAs
xyz_CDensAl = xyz_CDensAs
end select
!
! clear tendency
!
pyz_DVelXDtNl = 0.0d0
xqz_DVelYDtNl = 0.0d0
xyr_DVelZDtNl = 0.0d0
xyz_DKmDtNl = 0.0d0
xyz_DPTempDtNl = 0.0d0
xyz_DExnerDtNl = 0.0d0
xyz_DCDensDtNl = 0.0d0
xyzf_DQMixDtNl = 0.0d0
! 時間フィルタ.
! Time filter.
!
call AsselinTimeFilter
!------------------------------------------
! 移流に対する CFL 条件のチェック
! CFL condtion check for advection
!
! call CFLCheckTimeLongVelX( pyz_VelXNl ) !(in)
! call CFLCheckTimeLongVelY( xqz_VelYNl ) !(in)
! call CFLCheckTimeLongVelZ( xyr_VelZNl ) !(in)
! ヒストリファイル出力. 時間フィルタをかけた後の数値を出力.
! Out put to history file.
!
xyz_VelXNl = xyz_avr_pyz(pyz_VelXNl)
xyz_VelYNl = xyz_avr_xqz(xqz_VelYNl)
xyz_VelZNl = xyz_avr_xyr(xyr_VelZNl)
call HistoryAutoPut(TimeN, 'PTemp', xyz_PTempNl(1:nx, 1:ny, 1:nz))
call HistoryAutoPut(TimeN, 'PTempAll', xyz_PTempNl(1:nx, 1:ny, 1:nz) + xyz_PTempBZ(1:nx, 1:ny, 1:nz))
call HistoryAutoPut(TimeN, 'Exner', xyz_ExnerNl(1:nx, 1:ny, 1:nz))
call HistoryAutoPut(TimeN, 'ExnerAll', xyz_ExnerNl(1:nx, 1:ny, 1:nz) + xyz_ExnerBZ(1:nx, 1:ny, 1:nz))
call HistoryAutoPut(TimeN, 'VelX', xyz_VelXNl(1:nx, 1:ny, 1:nz))
call HistoryAutoPut(TimeN, 'VelY', xyz_VelYNl(1:nx, 1:ny, 1:nz))
call HistoryAutoPut(TimeN, 'VelZ', xyz_VelZNl(1:nx, 1:ny, 1:nz))
call HistoryAutoPut(TimeN, 'Km', xyz_KmNl(1:nx, 1:ny, 1:nz))
call HistoryAutoPut(TimeN, 'Kh', xyz_KhNl(1:nx, 1:ny, 1:nz))
call HistoryAutoPut(TimeN, 'CDens', xyz_CDensNl(1:nx, 1:ny, 1:nz))
do s = 1, ncmax
call HistoryAutoPut(TimeN, trim(SpcWetSymbol(s)), xyzf_QMixNl(1:nx, 1:ny, 1:nz, s))
call HistoryAutoPut(TimeN, trim(SpcWetSymbol(s))//'All', xyzf_QMixNl(1:nx, 1:ny, 1:nz, s) + xyzf_QMixBZ(1:nx, 1:ny, 1:nz, s))
end do
!----------------------------------------------------
! リスタートファイルの作成
! Make restartfile.
!
if (mod(t, NstepOutput) == 0) then
! ファイル出力
!
call HistoryPut( 't', TimeN, rstat)
call HistoryPut( 'VelX', pyz_VelXNl, rstat)
call HistoryPut( 'VelY', xqz_VelYNl, rstat)
call HistoryPut( 'VelZ', xyr_VelZNl, rstat)
call HistoryPut( 'Exner', xyz_ExnerNl, rstat)
call HistoryPut( 'PTemp', xyz_PTempNl, rstat)
call HistoryPut( 'Km', xyz_KmNl, rstat)
call HistoryPut( 'Kh', xyz_KhNl, rstat)
call HistoryPut( 'CDens', xyz_CDensNl, rstat)
call HistoryPut( 'QMix', xyzf_QMixNl, rstat)
call HistoryPut( 't', TimeA, rstat)
call HistoryPut( 'VelX', pyz_VelXAl, rstat)
call HistoryPut( 'VelY', xqz_VelYAl, rstat)
call HistoryPut( 'VelZ', xyr_VelZAl, rstat)
call HistoryPut( 'Exner', xyz_ExnerAl, rstat)
call HistoryPut( 'PTemp', xyz_PTempAl, rstat)
call HistoryPut( 'Km', xyz_KmAl, rstat)
call HistoryPut( 'Kh', xyz_KhAl, rstat)
call HistoryPut( 'CDens', xyz_CDensAl, rstat)
call HistoryPut( 'QMix', xyzf_QMixAl, rstat)
! Show CPU time
!
call ClocksetPredict
end if
! 長い時間ステップのループを回すための処置.
! Renew prognostic variables for next long time step integration.
!
pyz_VelXBl = pyz_VelXNl
pyz_VelXNl = pyz_VelXAl
xqz_VelYBl = xqz_VelYNl
xqz_VelYNl = xqz_VelYAl
xyr_VelZBl = xyr_VelZNl
xyr_VelZNl = xyr_VelZAl
xyz_PTempBl = xyz_PTempNl
xyz_PTempNl = xyz_PTempAl
xyz_ExnerBl = xyz_ExnerNl
xyz_ExnerNl = xyz_ExnerAl
xyz_KmBl = xyz_KmNl
xyz_KmNl = xyz_KmAl
xyz_KhBl = xyz_KhNl
xyz_KhNl = xyz_KhAl
xyz_CDensBl = xyz_CDensNl
xyz_CDensNl = xyz_CDensAl
xyzf_QMixBl = xyzf_QMixNl
xyzf_QMixNl = xyzf_QMixAl
t = t + 1
! 時刻の進行
! Progress time
!
call TimesetProgress
end do
!----------------------------------------------------
! 出力ファイルのクローズ
! Close out put files.
!
call HistoryFileio_finalize
call ReStartFileio_finalize
!----------------------------------------------------
! MPI END
!
call MPIWrapperFinalize
!----------------------------------------------------
! CPU Time
!
call ClocksetLoopStop
call ClocksetClose
contains
!-----------------------------------------------------------------------
subroutine VariableAllocate
!
!初期化として, 配列を定義し, 値としてゼロを代入する.
!
!暗黙の型宣言禁止
implicit none
!配列割り当て
allocate( pyz_VelXBl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( pyz_VelXNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_VelXNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( pyz_VelXAl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( pyz_VelXNs(imin:imax,jmin:jmax,kmin:kmax) )
allocate( pyz_VelXAs(imin:imax,jmin:jmax,kmin:kmax) )
pyz_VelXBl = 0.0d0
pyz_VelXNl = 0.0d0
pyz_VelXAl = 0.0d0
pyz_VelXNs = 0.0d0
pyz_VelXAs = 0.0d0
xyz_VelXNl = 0.0d0
allocate( xqz_VelYBl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xqz_VelYNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_VelYNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xqz_VelYAl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xqz_VelYNs(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xqz_VelYAs(imin:imax,jmin:jmax,kmin:kmax) )
xqz_VelYBl = 0.0d0
xqz_VelYNl = 0.0d0
xqz_VelYAl = 0.0d0
xqz_VelYNs = 0.0d0
xqz_VelYAs = 0.0d0
xyz_VelYNl = 0.0d0
allocate( xyr_VelZBl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyr_VelZNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_VelZNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyr_VelZAl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyr_VelZNs(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyr_VelZAs(imin:imax,jmin:jmax,kmin:kmax) )
xyr_VelZBl = 0.0d0
xyr_VelZNl = 0.0d0
xyr_VelZAl = 0.0d0
xyr_VelZNs = 0.0d0
xyr_VelZAs = 0.0d0
xyz_VelZNl = 0.0d0
allocate( xyz_ExnerBl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_ExnerNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_ExnerAl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_ExnerNs(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_ExnerAs(imin:imax,jmin:jmax,kmin:kmax) )
xyz_ExnerBl = 0.0d0
xyz_ExnerNl = 0.0d0
xyz_ExnerAl = 0.0d0
xyz_ExnerNs = 0.0d0
xyz_ExnerAs = 0.0d0
allocate( xyz_PTempBl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_PTempNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_PTempAl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_PTempNs(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_PTempAs(imin:imax,jmin:jmax,kmin:kmax) )
xyz_PTempBl = 0.0d0
xyz_PTempNl = 0.0d0
xyz_PTempAl = 0.0d0
xyz_PTempNs = 0.0d0
xyz_PTempAs = 0.0d0
allocate( xyz_CDensBl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_CDensNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_CDensAl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_CDensNs(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_CDensAs(imin:imax,jmin:jmax,kmin:kmax) )
xyz_CDensBl = 0.0d0
xyz_CDensNl = 0.0d0
xyz_CDensAl = 0.0d0
xyz_CDensNs = 0.0d0
xyz_CDensAs = 0.0d0
allocate( xyz_KmBl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_KmNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_KmAl(imin:imax,jmin:jmax,kmin:kmax) )
xyz_KmBl = 0.0d0
xyz_KmNl = 0.0d0
xyz_KmAl = 0.0d0
allocate( xyz_KhBl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_KhNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_KhAl(imin:imax,jmin:jmax,kmin:kmax) )
xyz_KhBl = 0.0d0
xyz_KhNl = 0.0d0
xyz_KhAl = 0.0d0
allocate( xyzf_QMixBl(imin:imax,jmin:jmax,kmin:kmax,ncmax) )
allocate( xyzf_QMixNl(imin:imax,jmin:jmax,kmin:kmax,ncmax) )
allocate( xyzf_QMixAl(imin:imax,jmin:jmax,kmin:kmax,ncmax) )
xyzf_QMixBl = 0.0d0
xyzf_QMixNl = 0.0d0
xyzf_QMixAl = 0.0d0
allocate( pyz_DVelXDtNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xqz_DVelYDtNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyr_DVelZDtNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_DKmDtNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_DPTempDtNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_DExnerDtNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_DExnerDtNs(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyz_DCDensDtNl(imin:imax,jmin:jmax,kmin:kmax) )
allocate( xyzf_DQMixDtNl(imin:imax,jmin:jmax,kmin:kmax,ncmax) )
pyz_DVelXDtNl = 0.0d0
xqz_DVelYDtNl = 0.0d0
xyr_DVelZDtNl = 0.0d0
xyz_DKmDtNl = 0.0d0
xyz_DPTempDtNl = 0.0d0
xyz_DExnerDtNl = 0.0d0
xyz_DExnerDtNs = 0.0d0
xyz_DCDensDtNl = 0.0d0
xyzf_DQMixDtNl = 0.0d0
! 時間ループの初期化
!
t = 1
end subroutine VariableAllocate
!-----------------------------------------------------------------------
subroutine AsselinTimeFilter
!
! 時間フィルター; Asselin のタイムフィルターを利用
! t = 0.0 の場合には tfil = 0.0d0, それ以外は tfil = 1.0d-1
! (t = 0 の時はオイラー法で, それ以外はリープフロッグ法で積分するため)
!
use TimeSet, only: tfil
implicit none
real(DP) :: aaa_VarNl(imin:imax,jmin:jmax,kmin:kmax)
real(DP) :: aaaa_VarNl(imin:imax,jmin:jmax,kmin:kmax,ncmax)
aaa_VarNl = pyz_VelXNl + tfil * (pyz_VelXAl - 2.0d0 * pyz_VelXNl + pyz_VelXBl)
pyz_VelXNl = aaa_VarNl
aaa_VarNl = xqz_VelYNl + tfil * (xqz_VelYAl - 2.0d0 * xqz_VelYNl + xqz_VelYBl)
xqz_VelYNl = aaa_VarNl
aaa_VarNl = xyr_VelZNl + tfil * (xyr_VelZAl - 2.0d0 * xyr_VelZNl + xyr_VelZBl)
xyr_VelZNl = aaa_VarNl
aaa_VarNl = xyz_PTempNl + tfil * (xyz_PTempAl - 2.0d0 * xyz_PTempNl + xyz_PTempBl)
xyz_PTempNl= aaa_VarNl
aaa_VarNl = xyz_ExnerNl + tfil * (xyz_ExnerAl - 2.0d0 * xyz_ExnerNl + xyz_ExnerBl)
xyz_ExnerNl= aaa_VarNl
aaa_VarNl = xyz_KmNl + tfil * (xyz_KmAl - 2.0d0 * xyz_KmNl + xyz_KmBl)
xyz_KmNl = aaa_VarNl
aaa_VarNl = xyz_KhNl + tfil * (xyz_KhAl - 2.0d0 * xyz_KhNl + xyz_KhBl)
xyz_KhNl = aaa_VarNl
aaa_VarNl = xyz_CDensNl + tfil * (xyz_CDensAl - 2.0d0 * xyz_CDensNl + xyz_CDensBl)
xyz_CDensNl = aaa_VarNl
aaaa_VarNl = xyzf_QMixNl + tfil * (xyzf_QMixAl - 2.0d0 * xyzf_QMixNl + xyzf_QMixBl)
xyzf_QMixNl = aaaa_VarNl
end subroutine AsselinTimeFilter
!-----------------------------------------------------------------------
subroutine MainInit
implicit none
character(STRING) :: cfgfile
! NAMELIST ファイル名 ; NAMELIST fine name
! MPI
!
call MPIWrapperInit
! 時間計測
!
call ClocksetInit
call ClocksetPreStart
! NAMELIST ファイル名の読み込み
! Loading NAMELIST file.
!
call argset_init( cfgfile ) !(out)
! NAMELIST ファイル名のモジュールへの登録
! Loading NAMELIST file.
!
call NmlutilInit( cfgfile ) !(in)
! 時間積分の初期化
! Initialization of time integration.
!
call timeset_init
! 格子点情報の初期化
! Initialization of grid arrangement.
!
call gridset_init
! 化学計算ルーチンの初期化
! Initialization of chemical routines.
!
call chemcalc_init
! 定数の情報の初期化
! Initialization of constant variables.
!
call constants_init
! 軸の情報の初期化
! Initialization of axis variables.
!
call axesset_init
! I/O ファイル名の初期化
! Initialization of output file name.
!
call fileset_init
! 湿潤過程共有変数の初期化
! Initialization of common variables for moist process.
!
call composition_init
! ヒストリーファイル・リスタートファイルの初期化
! Initialize restart & history files.
!
call HistoryFileio_init
call ReStartFileio_init
! 数値摩擦係数の初期化
! Initialization of numerical friction coefficient.
!
call Damping_Init
! 内部変数の初期化
! Initialization of internal variables.
!
call VariableAllocate
! 初期値の代入
! * ReStartFile が設定されている場合にはファイルを読み込む.
! 設定されていない場合にはデフォルトの基本場と擾乱場を作る.
!
! Initial value set up.
! * Read restartfile if it is specified. If not, make default basic
! state and disturbance.
!
if (myrank == 0) call MessageNotify( "M", "main", "Initial value setup." )
! 基本場, 擾乱場の初期値を netCDF ファイルから取得する.
!
call ReStartFileio_BZ_Get
call ReStartFileio_Var_Get( pyz_VelXBl, pyz_VelXNl, xqz_VelYBl, xqz_VelYNl, xyr_VelZBl, xyr_VelZNl, xyz_PTempBl, xyz_PTempNl, xyz_ExnerBl, xyz_ExnerNl, xyzf_QMixBl, xyzf_QMixNl, xyz_KmBl, xyz_KmNl, xyz_KhBl, xyz_KhNl, xyz_CDensBl, xyz_CDensNl )
! 力学モジュールの初期化
! Initialization of dynamical modules
!
call Dynamics_init
! 負の湿潤量の補填計算の初期化
! Initialization of negative moist value correction.
!
call FillNegative_Init
! 物理過程の設定
!
call deepconv_main
!-------------------------------------------------------------
! 基本場のファイル出力
!
call HistoryPut( 'DensBZ', xyz_DensBZ , rstat)
call HistoryPut( 'ExnerBZ', xyz_ExnerBZ , rstat)
call HistoryPut( 'PTempBZ', xyz_PTempBZ , rstat)
call HistoryPut( 'VelSoundBZ', xyz_VelSoundBZ , rstat)
call HistoryPut( 'TempBZ', xyz_TempBZ , rstat)
call HistoryPut( 'PressBZ', xyz_PressBZ , rstat)
call HistoryPut( 'QMixBZ', xyzf_QMixBZ , rstat)
call HistoryPut( 'EffMolWtBZ', xyz_EffMolWtBZ, rstat)
! 音波に対する CFL 条件のチェック
! CFL condtion check for sound wave.
!
call CFLCheckTimeShort( xyz_VelSoundBZ )
! CPU time
!
call ClocksetPreStop
! 最初の一回目の場合はオイラー法で回す.
!
if ( FlagInitialRun ) call TimesetDelTimeHalf
end subroutine MainInit
subroutine deepconv_main
use dc_message, only: MessageNotify
use dc_iounit, only : FileOpen
implicit none
integer :: unit !装置番号
character(STRING) :: FlagTurbMethod = ""
character(STRING) :: FlagRadMethod = ""
character(STRING) :: FlagCloudMethod = ""
character(STRING) :: FlagSurfaceMethod = ""
NAMELIST /deepconv_main_nml / FlagTurbMethod, FlagRadMethod, FlagCloudMethod, FlagSurfaceMethod
! デフォルト値の設定
! Default values settings
!
FlagTurbMethod = "Nothing"
!!$ FlagTurbMethod = "KW1978"
FlagRadMethod = "Nothing"
!!$ FlagRadMethod = "HeatConst"
!!$ FlagRadMethod = "HeatVary"
!!$ FlagRadMethod = "HeatBalance"
FlagSurfaceMethod = "Nothing"
!!$ FlagSurfaceMethod = "Diff"
!!$ FlagSurfaceMethod = "Bulk"
!!$ FlagSurfaceMethod = "Const"
FlagCloudMethod = "Nothing"
!!$ FlagCloudMethod = "K1969"
!!$ FlagCloudMethod = "MarsCond"
call FileOpen(unit, file=namelist_filename, mode='r')
read(unit, NML=deepconv_main_nml)
close(unit)
select case ( FlagTurbMethod )
case ( "Nothing" )
case ( "KW1978" )
IDTurbMethod = IDTurbMethodKW1978
call turbulence_kw1978_init
case default
call MessageNotify( 'E', prog_name, 'FlagTurbMethod=<%c> is not supported.', c1 = trim(FlagTurbMethod) )
end select
select case ( FlagRadMethod )
case ( "Nothing" )
case ( "HeatConst" )
IDRadMethod = IDRadMethodHeatConst
call Radiation_Simple_init
case ( "HeatVary" )
IDRadMethod = IDRadMethodHeatVary
call Radiation_Simple_init
case ( "HeatBalance" )
IDRadMethod = IDRadMethodHeatBalance
call radiation_heatbalance_init
case default
call MessageNotify( 'E', prog_name, 'FlagRadMethod=<%c> is not supported.', c1 = trim(FlagRadMethod) )
end select
select case ( FlagSurfaceMethod )
case ( "Nothing" )
case ( "Diff" )
IDSurfaceMethod = IDSurfaceMethodDiff
call Surfaceflux_Diff_init
case ( "Bulk" )
IDSurfaceMethod = IDSurfaceMethodBulk
call Surfaceflux_Bulk_init
case ( "Const" )
IDSurfaceMethod = IDSurfaceMethodConst
call Surfaceflux_Const_init
case default
call MessageNotify( 'E', prog_name, 'FlagSurfaceMethod=<%c> is not supported.', c1 = trim(FlagSurfaceMethod) )
end select
select case ( FlagCloudMethod )
case ( "Nothing" )
case ( "K1969" )
IDCloudMethod = IDCloudMethodK1969
call cloudphys_K1969_init
case ( "MarsCond" )
IDCloudMethod = IDCloudMethodMarsCond
call cloudphys_marscond_init
case default
call MessageNotify( 'E', prog_name, 'FlagCloudMethod=<%c> is not supported.', c1 = trim(FlagCloudMethod) )
end select
if (myrank == 0) then
call MessageNotify( "M", "main", "FlagTurbMethod = %c", c1=trim(FlagTurbMethod))
call MessageNotify( "M", "main", "IDTurbMethod = %d", i=(/IDTurbMethod/))
call MessageNotify( "M", "main", "FlagRadMethod = %c", c1=trim(FlagRadMethod))
call MessageNotify( "M", "main", "IDRadMethod = %d", i=(/IDRadMethod/))
call MessageNotify( "M", "main", "FlagSurfaceMethod = %c", c1=trim(FlagSurfaceMethod))
call MessageNotify( "M", "main", "IDSurfaceMethod = %d", i=(/IDSurfaceMethod/))
call MessageNotify( "M", "main", "FlagCloudMethod = %c", c1=trim(FlagCloudMethod))
call MessageNotify( "M", "main", "IDCloudMethod = %d", i=(/IDCloudMethod/))
end if
end subroutine deepconv_main
end program deepconv_arare