program arare
!
! 非静力学モデル deepconv/arare 湿潤大気対流計算用主プログラム (三次元版)
!
! モジュール引用 ; use statement
!
! gt4f90io 関連
! gt4f90io modules
!
use dc_types, only: STRING, DP
use dc_message, only: MessageNotify
! 初期設定モジュール
! Initialize module
!
use argset, only: argset_init
use fileset_3d, only: fileset_init, InitFile
use timeset, only: timeset_init, DelTimeLong, DelTimeShort, NstepLong, NstepShort, NstepDisp
use gridset_3d, only: gridset_init, DimXMin, DimXMax, DimYMin, DimYMax, DimZMin, DimZMax, SpcNum
use basicset_3d, only: basicset_init, xyza_MixRtBasicZ, xyz_DensBasicZ, xyz_PotTempBasicZ, xyz_VelSoundBasicZ
! 化学量計算モジュール
! Chemical calculation modules
!
use ChemCalc_3d, only: ChemCalc_init
use chemdata, only: chemdata_init
! 力学過程計算用関数モジュール
! Dynamical processes module
!
use DynFunc_3d, only: xyz_AdvScalar, xyz_AdvKm, pyz_AdvVelX, xqz_AdvVelY, xyr_Buoy, xyr_AdvVelZ, pyz_GradPi, xqz_GradPi
use DynImpFunc_3d, only: xyz_Exner_init, xyz_Exner, xyr_GradPi
! 乱流拡散計算用モジュール
! Turbulent diffusion module
!
use Turbulence_3d, only: Turbulence_Init, xyz_TurbScalar, xyza_TurbScalar, pyz_TurbVelX, xyr_TurbVelZ, xyz_ShearKm, xyz_DispKm, xyz_DispHeat, xqz_TurbVelY, EddyViscosity, xyz_BuoyKm
! 境界からのフラックス計算用モジュール
! Surface flux module
!
use HeatFlux_3d, only: xyz_HeatFluxBulk, pyz_MomFluxBulk, xqz_MomFluxBulk
! 放射強制計算用モジュール
! Radiative forceing module
!
use Radiation_3d, only: Radiation_init, xyz_RadHeatConst, xyz_NewtonCool
! 湿潤過程計算用モジュール
! Moist processes modules
!
use moistset, only: moistset_init
! 数値拡散/摩擦計算用モジュール
! Numerical diffussion /dumping module
!
use NumDiffusion_3d, only: NumDiffusion_Init, xyz_NumDiffScalar, xyz_NumDiffKm, pyz_NumDiffVelX, xqz_NumDiffVelY, xyr_NumDiffVelZ
use damping_3d, only: damping_init, xyz_DampSponge, pyz_DampSponge, xyr_DampSponge, xqz_DampSponge
! 積算値管理モジュール
! Monitor variables setup modules
!
use StorePotTemp_3d, only: StorePotTemp_init, StorePotTempClean, StorePotTempCond
! ファイル入出力モジュール
! File I/O module
!
use RestartFileIO_3d_dry,only: ReStartFile_Open, ReStartFile_OutPut, ReStartFile_Close, ReStartFile_Get
use HistoryFileIO_3d_dry,only: HistoryFile_Open, HistoryFile_OutPut, HistoryFile_Close
! 下請けモジュール
! Utility modules
!
use cflcheck_3d, only: CFLCheckTimeShort, CFLCheckTimeLongVelX, CFLCheckTimeLongVelY, CFLCheckTimeLongVelZ
use timefilter_3d, only: AsselinFilter
use xyz_bc_module, only: BoundaryXCyc_xyz, BoundaryYCyc_xyz, BoundaryZSym_xyz, BoundaryZCyc_xyz, BoundaryXCyc_pyz, BoundaryYCyc_pyz, BoundaryZSym_pyz, BoundaryZCyc_pyz, BoundaryXCyc_xqz, BoundaryYCyc_xqz, BoundaryZSym_xqz, BoundaryZCyc_xqz, BoundaryXCyc_xyr, BoundaryYCyc_xyr, BoundaryZSym_xyr, BoundaryZCyc_xyr, BoundaryZAsym_xyr
implicit none
! 内部変数
! Internal variables
!
character(80) :: cfgfile
! NAMELIST ファイル名 ; NAMELIST fine name
real(DP), allocatable :: pyz_VelXBl(:,:,:)
! $ u (t-\Delta t) $ 東西風 ; zonal wind
real(DP), allocatable :: pyz_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 :: 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 :: 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_PotTempBl(:,:,:)
! $ \theta (t-\Delta t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_PotTempNl(:,:,:)
! $ \theta (t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_PotTempAl(:,:,:)
! $ \theta (t+\Delta t) $ 温位 ; Potential temp.
real(DP), allocatable :: xyz_PotTempWork(:,:,:)
! 温位 $ \theta $ の作業配列 ; Work array
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 :: xyza_MixRtBl(:,:,:,:)
! $ q (t-\Delta t) $ 湿潤量の混合比
! Mixing ratio of moist variables.
real(DP), allocatable :: xyza_MixRtNl(:,:,:,:)
! $ q (t) $ 湿潤量の混合比
! Mixing ratio of moist variables
real(DP), allocatable :: xyza_MixRtAl(:,:,:,:) !
! $ q (t+\Delta t) $ 湿潤量の混合比
!Mixing ratio of moist variables
real(DP), allocatable :: pyz_AccelVelXNl(:,:,:)
! 気圧傾度力を除いた $u$ の変化率
! Tendency of $u$ except for pressure gradient term
real(DP), allocatable :: xqz_AccelVelYNl(:,:,:)
! 気圧傾度力を除いた $v$ の変化率
! Tendency of $v$ except for pressure gradient term
real(DP), allocatable :: xyr_AccelVelZNl(:,:,:)
! 気圧傾度力を除いた $w$ の変化率
! Tendency of $w$ except for pressure gradient term
real(DP) :: Time ! 時刻 ; Time
real(DP) :: ReStartTime(2) ! リスタートファイル出力時刻用配列
! Output time array for restart file
real(DP), allocatable :: DelTimeLFrog(:)
! リープフロッグスキーム用時間格子間隔
! Time interval for Leap-frog scheme
real(DP) :: DelTimeEuler ! オイラースキーム用時間格子
! Time interval for Eular scheme
integer :: NStepLFrog ! リープフロッグスキーム用時間ステップ数
! The number of time step for Leap-frog scheme
integer, allocatable :: NStepEuler(:)
! オイラースキーム用時間ステップ数
! The number of time step for Eular scheme
integer :: t, tau, t1, t2, s ! do ループ変数 ; do loop variable
! 初期化手続き ; Initialize procedure
!
! NAMELIST ファイル名の読み込み
! Loading NAMELIST file.
!
call argset_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( )
! 内部変数の初期化
! 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.
!
call MessageNotify( "M", "main", "Initial value setup." )
if (trim(InitFile) /= '') then
call ReStartFile_Get( ReStartTime, xyz_PotTempBl, xyz_ExnerBl, pyz_VelXBl, xqz_VelYBl, xyr_VelZBl, xyz_KmBl, xyz_KhBl, xyz_PotTempNl, xyz_ExnerNl, pyz_VelXNl, xqz_VelYBl, xyr_VelZNl, xyz_KmNl, xyz_KhNl )
else
call BasicEnv_3d()
call DisturbEnv_3d( cfgfile, xyz_PotTempBl, xyz_ExnerBl, pyz_VelXBl, xqz_VelYBl, xyr_VelZBl, xyza_MixRtBl, xyz_KmBl, xyz_KhBl )
call BoundaryXCyc_pyz( pyz_VelXBl )
call BoundaryYCyc_pyz( pyz_VelXBl )
call BoundaryZSym_pyz( pyz_VelXBl )
call BoundaryXCyc_xqz( xqz_VelYBl )
call BoundaryYCyc_xqz( xqz_VelYBl )
call BoundaryZSym_xqz( xqz_VelYBl )
call BoundaryXCyc_xyr( xyr_VelZBl )
call BoundaryYCyc_xyr( xyr_VelZBl )
call BoundaryZAsym_xyr( xyr_VelZBl )
call BoundaryXCyc_xyz( xyz_PotTempBl )
call BoundaryYCyc_xyz( xyz_PotTempBl )
call BoundaryZSym_xyz( xyz_PotTempBl )
do s = 1, SpcNum
call BoundaryXCyc_xyz( xyza_MixRtBl(:,:,:,s) )
call BoundaryYCyc_xyz( xyza_MixRtBl(:,:,:,s) )
call BoundaryZSym_xyz( xyza_MixRtBl(:,:,:,s) )
end do
! 時刻 $ t $ の変数値の初期値の設定
! * 1 ループ目だけは $ t $ の値を $ t-\Delta t$ の値と同じにする.
! 1 ステップ目はオイラー法で解く必要があるが, 1 ステップ目と
! それ以外のステップを別々にコーディングしたくない為
!
! Set up initial value of time = "t" variables.
!
xyz_PotTempNl = xyz_PotTempBl
xyz_ExnerNl = xyz_ExnerBl
pyz_VelXNl = pyz_VelXBl
xqz_VelYNl = xqz_VelYBl
xyr_VelZNl = xyr_VelZBl
xyza_MixRtNl = xyza_MixRtBl
xyz_KmNl = xyz_KmBl
xyz_KhNl = xyz_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 radiative forcing.
!
call Radiation_Init( cfgfile )
! 圧力計算用係数行列の初期化
! Initialization of coefficient matrix for exner function calculation.
!
call xyz_Exner_Init()
! 時刻とループ回数の初期化
! Initialization of time integration parameters.
!
NstepLFrog = NstepLong
NstepEuler = NstepShort
DelTimeLFrog = DelTimeLong * 2.0d0
DelTimeEuler = 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
NstepEuler(1) = NstepShort /2
DelTimeLFrog(1) = DelTimeLong
end if
call MessageNotify( "M", "main", "NstepLong= %d", i=(/NstepLong/) )
call MessageNotify( "M", "main", "NstepEuler= %d", i=(/NstepEuler/) )
call MessageNotify( "M", "main", "DelTimeLFrog= %f", d=(/DelTimeLFrog/) )
call MessageNotify( "M", "main", "DelTimeEuler= %f", d=(/DelTimeEuler/) )
! ヒストリーファイルへの出力
! Out put to history file.
!
call HistoryFile_Open( )
if ( trim(InitFile) /= '') then
call HistoryFile_Output( ReStartTime(2), xyz_PotTempNl, xyz_ExnerNl, pyz_VelXNl, xqz_VelYNl, xyr_VelZNl, xyz_KmNl, xyz_KhNl )
end if
if ( Time == 0 ) then
call HistoryFile_Output( Time, xyz_PotTempNl, xyz_ExnerNl, pyz_VelXNl, xqz_VelYNl, xyr_VelZNl, xyz_KmNl, xyz_KhNl )
end if
! 音波に対する CFL 条件のチェック
! CFL condtion check for sound wave.
!
call CFLCheckTimeShort( xyz_VelSoundBasicZ )
! 時間積分 ; time integration
!
call MessageNotify( "M", "main", "Time Integration Start" )
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.
!
! 1 次クロージャの場合
!
!call EddyViscosity( &
! & pyz_VelXNl, & ! (in)
! & xqz_VelYNl, & ! (in)
! & xyr_VelZNl, & ! (in)
! & xyz_PotTempNl, & ! (in)
! & xyz_KmNl, & ! (out)
! & xyz_KmAl & ! (out)
! & )
xyz_KmAl = xyz_KmBl + DelTimeLFrog(t) * ( + xyz_AdvKm(xyz_KmNl, pyz_VelXNl, xqz_VelYNl, xyr_VelZNl) + xyz_BuoyKm(xyz_PotTempBl) + xyz_ShearKm(xyz_KmBl, pyz_VelXBl, xqz_VelYBl, xyr_VelZBl) + xyz_NumDiffKm(xyz_KmBl) + xyz_DispKm(xyz_KmBl) )
!値の上限下限の設定
! * 値は正になることを保証する
! * 値の上限は 800 とする. 中島健介(1994, 学位論文)参照
!
! Upper and lower bound value are specified.
!
xyz_KmAl = max( 0.0d0, min( xyz_KmAl, 800.0d0 ) )
! 境界条件 ; Boundary condition
!
call BoundaryXCyc_xyz( xyz_KmAl ) ! (inout)
call BoundaryYCyc_xyz( xyz_KmAl ) ! (inout)
call BoundaryZSym_xyz( xyz_KmAl ) ! (inout)
! スカラーに対する渦拡散係数の計算
! Specify turbulent diffusion coefficient for scalar variables.
!
xyz_KhAl = 3.0d0 * xyz_KmAl
! 温位の移流拡散の計算
! Advection and diffusion of potential temperature.
!
xyz_PotTempAl = xyz_PotTempBl + DelTimeLFrog(t) * ( + xyz_AdvScalar( xyz_PotTempNl, pyz_VelXNl, xqz_VelYNl, xyr_VelZNl) + xyz_AdvScalar( xyz_PotTempBasicZ, pyz_VelXNl, xqz_VelYNl, xyr_VelZNl) + xyz_NumDiffScalar(xyz_PotTempBl) + xyz_TurbScalar(xyz_PotTempBl, xyz_KhBl) + xyz_TurbScalar(xyz_PotTempBasicZ, xyz_KhBl) + xyz_DispHeat( xyz_KmBl ) )
! 境界条件 ; Boundary condition
!
call BoundaryXCyc_xyz( xyz_PotTempAl ) ! (inout)
call BoundaryYCyc_xyz( xyz_PotTempAl ) ! (inout)
call BoundaryZSym_xyz( xyz_PotTempAl ) ! (inout)
! 速度の移流拡散.
! Advection and diffusion of velocity components.
!
pyz_AccelVelXNl = + pyz_AdvVelX(pyz_VelXNl, xqz_VelYNl, xyr_VelZNl) + pyz_NumDiffVelX(pyz_VelXBl) + pyz_TurbVelX(xyz_KmBl, pyz_VelXBl, xqz_VelYBl, xyr_VelZBl)
! & + pyz_MomFluxBulk( pyz_VelXBl )
xqz_AccelVelYNl = + xqz_AdvVelY(pyz_VelXNl, xqz_VelYNl, xyr_VelZNl) + xqz_NumDiffVelY(xqz_VelYBl) + xqz_TurbVelY(xyz_KmBl, pyz_VelXBl, xqz_VelYBl, xyr_VelZBl)
! & + xqz_MomFluxBulk( xqz_VelYBl )
xyr_AccelVelZNl = + xyr_AdvVelZ(pyz_VelXNl, xqz_VelYNl, xyr_VelZNl) + xyr_Buoy(xyz_PotTempNl) + xyr_NumDiffVelZ(xyr_VelZBl) + xyr_TurbVelZ(xyz_KmBl, pyz_VelXBl, xqz_VelYBl, xyr_VelZBl)
! 短い時間ステップの初期値作成.
! 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
! 短い時間ステップの初期値作成.
! Initial values set up for time integration with short time step.
!
Euler: do tau = 1, NstepEuler(t)
! 速度 u の計算.
! Time integration horizontal velocity.
!
pyz_VelXAs = pyz_VelXNs + DelTimeEuler * ( - pyz_GradPi(xyz_ExnerNs, pyz_VelXNs, xqz_VelYNs, xyr_VelZNs) + pyz_AccelVelXNl )
! 境界条件 ; Boundary condition
!
call BoundaryXCyc_pyz( pyz_VelXAs ) ! (inout)
call BoundaryYCyc_pyz( pyz_VelXAs ) ! (inout)
call BoundaryZSym_pyz( pyz_VelXAs ) ! (inout)
! 速度 v の計算.
! Time integration horizontal velocity.
!
xqz_VelYAs = xqz_VelYNs + DelTimeEuler * ( - xqz_GradPi(xyz_ExnerNs, pyz_VelXNs, xqz_VelYNs, xyr_VelZNs) + xqz_AccelVelYNl )
! 境界条件 ; Boundary condition
!
call BoundaryXCyc_xqz( xqz_VelYAs ) ! (inout)
call BoundaryYCyc_xqz( xqz_VelYAs ) ! (inout)
call BoundaryZSym_xqz( xqz_VelYAs ) ! (inout)
! エクスナー関数の計算.
! Time integration exner function.
!
xyz_ExnerAs = xyz_Exner( xyr_AccelVelZNl, pyz_VelXNs, pyz_VelXAs, xqz_VelYNs, xqz_VelYAs, xyr_VelZNs, xyz_ExnerNs)
! 境界条件 ; Boundary condition
!
call BoundaryXCyc_xyz( xyz_ExnerAs ) ! (iuout)
call BoundaryYCyc_xyz( xyz_ExnerAs ) ! (iuout)
call BoundaryZSym_xyz( xyz_ExnerAs ) ! (iuout)
! 速度 w の計算
! Time integration vertical velocity.
!
xyr_VelZAs = xyr_VelZNs + DelTimeEuler * ( - xyr_GradPi(xyz_ExnerAs,xyz_ExnerNs, pyz_VelXNs,xqz_VelYNs, xyr_VelZNs) + xyr_AccelVelZNl )
! 境界条件 ; Boundary condition
!
call BoundaryXCyc_xyr( xyr_VelZAs ) ! (inout)
call BoundaryYCyc_xyr( xyr_VelZAs ) ! (inout)
call BoundaryZAsym_xyr( xyr_VelZAs )! (inout)
! 短い時間ステップのループを回すための処置
! 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
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
! 時間フィルタ.
! Time filter.
!
call AsselinFilter( xyz_ExnerAl, xyz_ExnerNl, xyz_ExnerBl )
call AsselinFilter( pyz_VelXAl, pyz_VelXNl, pyz_VelXBl )
call AsselinFilter( xqz_VelYAl, xqz_VelYNl, xqz_VelYBl )
call AsselinFilter( xyr_VelZAl, xyr_VelZNl, xyr_VelZBl )
call AsselinFilter( xyz_PotTempAl, xyz_PotTempNl, xyz_PotTempBl )
call AsselinFilter( xyz_KmAl, xyz_KmNl, xyz_KmBl )
! スポンジ層
! Numerical dumping.
!
! pyz_VelXAl = pyz_DampSponge( pyz_VelXAl, pyz_VelXBl, DelTimeLFrog(t) )
! xqz_VelYAl = xqz_DampSponge( xqz_VelYAl, xqz_VelYBl, DelTimeLFrog(t) )
! xyr_VelZAl = xyr_DampSponge( xyr_VelZAl, xyr_VelZBl, DelTimeLFrog(t) )
! 長い時間ステップのループを回すための処置.
! Renew prognostic variables for next long time step integration.
!
xyz_PotTempBl = xyz_PotTempNl
xyz_ExnerBl = xyz_ExnerNl
pyz_VelXBl = pyz_VelXNl
xqz_VelYBl = xqz_VelYNl
xyr_VelZBl = xyr_VelZNl
xyz_KmBl = xyz_KmNl
xyz_KhBl = xyz_KhNl
xyz_PotTempNl = xyz_PotTempAl
xyz_ExnerNl = xyz_ExnerAl
pyz_VelXNl = pyz_VelXAl
xqz_VelYNl = xqz_VelYAl
xyr_VelZNl = xyr_VelZAl
xyz_KmNl = xyz_KmAl
xyz_KhNl = xyz_KhAl
end do
! ヒストリーファイルへの出力.
! Out put to history file.
!
call MessageNotify( "M", "main", "Time = %f", d=(/Time/) )
! 移流に対する CFL 条件のチェック
! CFL condtion check for advection
!
call CFLCheckTimeLongVelX( pyz_VelXNl )
call CFLCheckTimeLongVelY( xqz_VelYNl )
call CFLCheckTimeLongVelZ( xyr_VelZNl )
!ヒストリファイル出力.
! Out put to history file.
!
call HistoryFile_Output( Time, xyz_PotTempNl, xyz_ExnerNl, pyz_VelXNl, xqz_VelYNl, xyr_VelZNl, xyz_KmNl, xyz_KhNl )
! 積算値のクリア.
! Clear monitor variables.
!
call StorePotTempClean
end do
call MessageNotify( "M", "main", "Time Integration is finished." )
! 出力ファイルのクローズ
! Close out put files.
!
call HistoryFile_Close
! リスタートファイルの作成
! Out put to restart file.
!
call ReStartFile_Open( )
call ReStartFile_OutPut( Time - DelTimeLong, xyz_PotTempBl, xyz_ExnerBl, pyz_VelXBl, xqz_VelYBl, xyr_VelZBl, xyz_KmBl, xyz_KhBl )
call ReStartFile_OutPut( Time, xyz_PotTempNl, xyz_ExnerNl, pyz_VelXNl, xqz_VelYNl, xyr_VelZNl, xyz_KmNl, xyz_KhNl )
call ReStartFile_Close
contains
!-----------------------------------------------------------------------
subroutine VariableAllocate
!
!初期化として, 配列を定義し, 値としてゼロを代入する.
!
!暗黙の型宣言禁止
implicit none
!配列割り当て
allocate( pyz_VelXBl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), pyz_VelXNl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), pyz_VelXAl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), pyz_VelXNs(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), pyz_VelXAs(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xqz_VelYBl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xqz_VelYNl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xqz_VelYAl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xqz_VelYNs(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xqz_VelYAs(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyr_VelZBl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyr_VelZNl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyr_VelZAl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyr_VelZNs(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyr_VelZAs(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_ExnerBl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_ExnerNl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_ExnerAl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_ExnerNs(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_ExnerAs(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_PotTempWork(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_PotTempBl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_PotTempNl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_PotTempAl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_KmBl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_KmNl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_KmAl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_KhBl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_KhNl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyz_KhAl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyza_MixRtBl (DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax,SpcNum), xyza_MixRtNl (DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax,SpcNum), xyza_MixRtAl (DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax,SpcNum), pyz_AccelVelXNl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xqz_AccelVelYNl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), xyr_AccelVelZNl(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax), DelTimeLFrog(NstepLong), NStepEuler(NStepLong) )
pyz_VelXNs = 0.0d0
pyz_VelXAs = 0.0d0
xqz_VelYNs = 0.0d0
xqz_VelYAs = 0.0d0
xyr_VelZNs = 0.0d0
xyr_VelZAs = 0.0d0
xyz_ExnerNs = 0.0d0
xyz_ExnerAs = 0.0d0
pyz_VelXBl = 0.0d0
pyz_VelXNl = 0.0d0
pyz_VelXAl = 0.0d0
xqz_VelYBl = 0.0d0
xqz_VelYNl = 0.0d0
xqz_VelYAl = 0.0d0
xyr_VelZBl = 0.0d0
xyr_VelZNl = 0.0d0
xyr_VelZAl = 0.0d0
xyz_ExnerBl = 0.0d0
xyz_ExnerNl = 0.0d0
xyz_ExnerAl = 0.0d0
xyz_KmBl = 0.0d0
xyz_KmNl = 0.0d0
xyz_KmAl = 0.0d0
xyz_KhBl = 0.0d0
xyz_KhNl = 0.0d0
xyz_KhAl = 0.0d0
xyz_PotTempBl = 0.0d0
xyz_PotTempNl = 0.0d0
xyz_PotTempAl = 0.0d0
xyza_MixRtBl = 0.0d0
xyza_MixRtNl = 0.0d0
xyza_MixRtAl = 0.0d0
pyz_AccelVelXNl = 0.0d0
xqz_AccelVelYNl = 0.0d0
xyr_AccelVelZNl = 0.0d0
DelTimeEuler = 0.0d0
DelTimeLFrog = 0.0d0
NStepEuler = 0.0d0
NstepLFrog = 0.0d0
end subroutine VariableAllocate
!-----------------------------------------------------------------------
end program arare