!= deepconv/arare 湿潤大気対流計算用主プログラム (三次元版) ! != deepconv/arare main program for moist atmospheric convection (3D) ! ! Authors:: SUGIYAMA Ko-ichiro, ODAKA Masatsugu ! Version:: $Id: ararempi_3d.f90,v 1.1 2011-03-30 13:53:08 sugiyama Exp $ ! Tag Name:: $Name: arare4-20120911 $ ! Copyright:: Copyright (C) GFD Dennou Club, 2007. All rights reserved. ! License:: See COPYRIGHT[link:../../COPYRIGHT] ! program ararempi ! ! 非静力学モデル deepconv/arare 湿潤大気対流計算用主プログラム (三次元版) ! ! モジュール引用 use statement ! ! gtool5 関連 ! gtool5 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, xyz_DensBasicZ, xyza_MixRtBasicZ, & & xyz_PotTempBasicZ, xyz_VelSoundBasicZ ! 化学量計算モジュール ! Chemical calculation modules ! use ChemCalc_3d, only: ChemCalc_init use chemdata, only: chemdata_init ! MPI の初期化 ! use mpiset, only : mpii_init, mpii_end, myrank use mpi_wrapper, only : MPIWrapperInit, MPIWrapperFinalize, myrank ! 力学過程計算用関数モジュール ! Dynamical processes module ! use DynFunc_3d, only: xyz_AdvScalar, xyza_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_BuoyKm, & & xyz_TurbScalar, xyza_TurbScalar, pyz_TurbVelX, & & xyr_TurbVelZ , xyz_ShearKm , xyz_DispKm, & & xyz_DispHeat , xqz_TurbVelY , EddyViscosity ! 境界からのフラックス計算用モジュール ! 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_3d, only: moistset_init use MoistAdjust_3d, only: MoistAdjustSvapPress, MoistAdjustNH4SH use WarmRainPrm_3d, only: WarmRainPrm_Init, xyz_Rain2GasHeat, xyza_Rain2Gas,& & xyza_Rain2GasNH4SH, xyz_Rain2GasHeatNH4SH, & & xyza_Cloud2Rain, xyza_FallRain use MoistBuoyancy_3d,only: MoistBuoy_Init, xyz_BuoyMoistKm, xyr_BuoyMolWt, & & xyr_BuoyDrag use fillnegative_3d,only: FillNegative_init, xyza_FillNegative_xyza ! 安定度の計算 ! Static stability calculation module ! use ECCM_3D, only: ECCM_Stab ! 数値拡散/摩擦計算用モジュール ! Numerical diffussion /dumping module ! use NumDiffusion_3d,only: NumDiffusion_Init, xyz_NumDiffScalar, & & xyz_NumDiffKm, xyza_NumDiffScalar, & & pyz_NumDiffVelX, xqz_NumDiffVelY, & & xyr_NumDiffVelZ use damping_3d, only: damping_init, & ! & xyz_DampSponge, pyz_DampSponge, xyr_DampSponge, xqz_DampSponge & DampSponge_xyz, DampSponge_xyr, DampSponge_pyz,DampSponge_xqz ! 積算値管理モジュール ! Monitor variables setup modules ! use StorePotTemp_3d,only : StorePotTemp_init, StorePotTempClean, & & StorePotTempCond use StoreMixRt_3d, only : StoreMixRt_init, StoreMixRtClean, & & StoreMixRtCond, StoreMixRtFill1, StoreMixRtFill2 use StoreBuoy_3d, only : StoreBuoy_init, StoreBuoyClean use StoreStab_3d, only : StoreStab_init, StoreStabClean ! ファイル入出力モジュール ! File I/O module ! use RestartFileIO_3d, only : ReStartFile_Open, ReStartFile_OutPut, & & ReStartFile_Close, ReStartFile_Get use HistoryFileIO_3d, 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_mpi, 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 :: xyza_MixRtWork(:,:,:,:) ! 湿潤量の作業配列 ! Work array for mixing ratio. 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), allocatable :: xyza_DelMixRt(:,:,:,:) ! 湿潤量の混合比 $ q $ の増分 ! Mixing ratio variation. 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, & ! do ループ変数 ; do loop variable & t1, & ! do ループ変数 ; do loop variable & t2, & ! do ループ変数 ; do loop variable & s ! do ループ変数 ; do loop variable ! 初期化手続き ; Initialize procedure ! !MPI 並列 ! call mpii_init call MPIWrapperInit ! NAMELIST ファイル名の読み込み ! Loading NAMELIST file. ! call argset_init( & & cfgfile & ! (out) & ) ! 化学定数の初期化 ! Initialization of chemical constatns. ! call chemdata_init( ) ! 時間積分の初期化 ! Initialization of time integration. ! call timeset_init( & & cfgfile & ! (in) & ) ! 格子点情報の初期化 ! Initialization of grid arrangement. ! call gridset_init( & & cfgfile & ! (in) & ) ! 化学計算ルーチンの初期化 ! Initialization of chemical routines. ! call chemcalc_init( ) ! 基本場の情報の初期化 ! Initialization of basic state variables. ! call basicset_init( & & cfgfile & ! (in) & ) ! I/O ファイル名の初期化 ! Initialization of output file name. ! call fileset_init( & & cfgfile, myrank & ! (in) & ) ! 湿潤過程共有変数の初期化 ! 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. ! 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 !基本場, 擾乱場の初期値を netCDF ファイルから取得する. call ReStartFile_Get( & & ReStartTime, & ! (out) & xyz_PotTempBl, & ! (out) & xyz_ExnerBl, & ! (out) & pyz_VelXBl, & ! (out) & xqz_VelYBl, & ! (out) & xyr_VelZBl, & ! (out) & xyza_MixRtBl, & ! (out) & xyz_KmBl, & ! (out) & xyz_KhBl, & ! (out) & xyz_PotTempNl, & ! (out) & xyz_ExnerNl, & ! (out) & pyz_VelXNl, & ! (out) & xqz_VelYNl, & ! (out) & xyr_VelZNl, & ! (out) & xyza_MixRtNl, & ! (out) & xyz_KmNl, & ! (out) & xyz_KhNl & ! (out) & ) else call BasicEnv_3d() call MessageNotify( "M", "main", "Time Integration Start 1-1!" ) call DisturbEnv_3d( & & cfgfile, & ! (in) & xyz_PotTempBl, & ! (out) & xyz_ExnerBl, & ! (out) & pyz_VelXBl, & ! (out) & xqz_VelYBl, & ! (out) & xyr_VelZBl, & ! (out) & xyza_MixRtBl, & ! (out) & xyz_KmBl, & ! (out) & xyz_KhBl & ! (out) & ) call MessageNotify( "M", "main", "Time Integration Start 1-2!" ) call BoundaryXCyc_pyz( pyz_VelXBl ) ! (inout) call BoundaryYCyc_pyz( pyz_VelXBl ) ! (inout) call BoundaryZSym_pyz( pyz_VelXBl ) ! (inout) call BoundaryXCyc_xqz( xqz_VelYBl ) ! (inout) call BoundaryYCyc_xqz( xqz_VelYBl ) ! (inout) call BoundaryZSym_xqz( xqz_VelYBl ) ! (inout) call BoundaryXCyc_xyr( xyr_VelZBl ) ! (inout) call BoundaryYCyc_xyr( xyr_VelZBl ) ! (inout) call BoundaryZAsym_xyr( xyr_VelZBl ) ! (inout) call BoundaryXCyc_xyz( xyz_PotTempBl ) ! (inout) call BoundaryYCyc_xyz( xyz_PotTempBl ) ! (inout) call BoundaryZSym_xyz( xyz_PotTempBl ) ! (inout) do s = 1, SpcNum call BoundaryXCyc_xyz( xyza_MixRtBl(:,:,:,s) ) ! (inout) call BoundaryYCyc_xyz( xyza_MixRtBl(:,:,:,s) ) ! (inout) call BoundaryZSym_xyz( xyza_MixRtBl(:,:,:,s) ) ! (inout) 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 call MessageNotify( "M", "main", "Time Integration Start 2!" ) ! 数値摩擦係数の初期化 ! Initialization of numerical friction coefficient. ! call Damping_Init( & & cfgfile & ! (in) & ) ! 数値拡散項の初期化 ! Initialization of numerical diffusion term. ! call NumDiffusion_Init( ) ! 乱流拡散項の初期化 ! Initialization of turbulent diffusion term. ! call Turbulence_Init() ! 暖かい雨のパラメタリゼーションの初期化 ! Initialization of warmrain parameterization. ! call WarmRainPrm_Init( & & cfgfile & ! (in) & ) ! 負の湿潤量の補填計算の初期化 ! Initialization of negative moist value correction. ! call FillNegative_Init( & & xyza_MixRtBasicZ, & ! (in) & xyz_DensBasicZ & ! (in) & ) ! 放射強制の初期化 ! Initialization of radiative forcing. ! call Radiation_Init( & & cfgfile & ! (in) & ) ! 湿潤気塊の浮力計算の初期化 ! Initialization of moist buoyancy calculation. ! call MoistBuoy_Init( ) ! 圧力計算用係数行列の初期化 ! Initialization of coefficient matrix for exner function calculation. ! call xyz_Exner_Init() !陰解法の初期化 ! 時刻とループ回数の初期化 ! Initialization of time integration. ! 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), & ! (in) & xyz_PotTempNl, & ! (in) & xyz_ExnerNl, & ! (in) & pyz_VelXNl, & ! (in) & xqz_VelYNl, & ! (in) & xyr_VelZNl, & ! (in) & xyza_MixRtNl, & ! (in) & xyz_KmNl, & ! (in) & xyz_KhNl & ! (in) & ) end if !時刻 0 の場合には出力 if ( Time == 0 ) then call HistoryFile_Output( & & Time, & ! (in) & xyz_PotTempNl, & ! (in) & xyz_ExnerNl, & ! (in) & pyz_VelXNl, & ! (in) & xqz_VelYNl, & ! (in) & xyr_VelZNl, & ! (in) & xyza_MixRtNl, & ! (in) & xyz_KmNl, & ! (in) & xyz_KhNl & ! (in) & ) end if ! 音波に対する CFL 条件のチェック ! CFL condtion check for sound wave. ! call CFLCheckTimeShort( & & xyz_VelSoundBasicZ & ! (in) & ) ! 時間積分 time integration ! call MessageNotify( "M", "main", "Time Integration Start" ) do t1 = 1, NstepLFrog / NstepDisp do t2 = 1, NstepDisp ! 時刻の設定 ! Time setting. ! t = (t1 - 1) * NstepDisp + t2 Time = Time + DelTimeLong ! 渦拡散係数の移流拡散 ! 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 BoundaryZSym_xyz( xyz_KmAl ) ! (inout) call BoundaryYCyc_xyz( xyz_KmAl ) ! (inout) call BoundaryXCyc_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_TurbScalar(xyz_PotTempBl, xyz_KhBl) & & + xyz_TurbScalar(xyz_PotTempBasicZ, xyz_KhBl) & & + xyz_NumDiffScalar(xyz_PotTempBl) & & + xyz_DispHeat( xyz_KmBl ) & !! & + xyz_RadHeatConst( xyz_ExnerBl ) & !! & + xyz_HeatFluxDiff( xyz_PotTempNl ) & !! & + xyz_HeatFluxBulk( xyz_PotTempBl, pyz_VelXBl, xqz_VelYBl ) & !! & + xyz_NewtonCool( xyz_PotTempBl ) & & ) ! 境界条件 ; Boundary condition ! call BoundaryXCyc_xyz( xyz_PotTempAl ) ! (inout) call BoundaryYCyc_xyz( xyz_PotTempAl ) ! (inout) call BoundaryZSym_xyz( xyz_PotTempAl ) ! (inout) ! 凝縮成分混合比の移流拡散 ! Advection and diffusion of vapor, cloud and rain mixing ratios. ! xyza_MixRtAl = & & xyza_MixRtBl & & + DelTimeLFrog(t) & & * ( & & + xyza_AdvScalar (xyza_MixRtNl, pyz_VelXNl, xqz_VelYNl, xyr_VelZNl)& & + xyza_AdvScalar (xyza_MixRtBasicZ, pyz_VelXNl, xqz_VelYNl, xyr_VelZNl)& & + xyza_TurbScalar(xyza_MixRtBl, xyz_KhBl) & & + xyza_TurbScalar(xyza_MixRtBasicZ,xyz_KhBl) & & + xyza_NumDiffScalar(xyza_MixRtBl) & & + xyza_FallRain(xyza_MixRtBl) & !! & + xyza_MixRtFluxDiff(xyza_MixRtNl) & !! & + xyza_MixRtFluxBulk(xyza_MixRtNl) & & ) ! 境界条件 Boundary condition ! do s = 1, SpcNum call BoundaryXCyc_xyz( xyza_MixRtAl(:,:,:,s) ) ! (inout) call BoundaryYCyc_xyz( xyza_MixRtAl(:,:,:,s) ) ! (inout) call BoundaryZSym_xyz( xyza_MixRtAl(:,:,:,s) ) ! (inout) end do ! 移流によって負になった部分を埋める ! Negative values due to advection are corrected. ! xyza_MixRtWork = xyza_MixRtAl xyza_MixRtAl = xyza_FillNegative_xyza( xyza_MixRtWork ) ! 埋めた/削った量を保管 ! Correction value is stored. ! call StoreMixRtFill1( & & (xyza_MixRtAl - xyza_MixRtWork) / DelTimeLFrog(t) & ! (in) & ) ! 埋め切れなかった部分をゼロにする. Fill2 に保管 ! Negative values mixing ratios are corrected. xyza_MixRtWork = xyza_MixRtAl xyza_MixRtAl = max( - xyza_MixRtBasicZ, xyza_MixRtWork ) ! call StoreMixRtFill2( & ! & (xza_MixRtAl - xza_MixRtWork) / DelTimeLFrog(t) & ! (in) ! & ) ! 境界条件 Boundary condition ! do s = 1, SpcNum call BoundaryXCyc_xyz( xyza_MixRtAl(:,:,:,s) ) ! (inout) call BoundaryYCyc_xyz( xyza_MixRtAl(:,:,:,s) ) ! (inout) call BoundaryZSym_xyz( xyza_MixRtAl(:,:,:,s) ) ! (inout) end do ! 暖かい雨のパラメタリゼーション. ! * 雲<-->雨 の変換を行う. ! ! Warm rain parameterization. ! * Conversion from cloud to rain. !これまでの値を作業配列に保管 ! Previous values are stored to work area. ! xyza_MixRtWork = xyza_MixRtAl !雨への変化量を計算 ! Conversion values are calculated. ! xyza_MixRtAl = xyza_MixRtWork & & + xyza_Cloud2Rain( xyza_MixRtAl, DelTimeLFrog(t) ) ! 雲から雨への変換量を保管 ! Conversion values are sotred. ! call StoreMixRtCond( & & (xyza_MixRtAl - xyza_MixRtWork) / DelTimeLFrog(t) & ! (in) & ) ! 境界条件 Boundary condition ! do s = 1, SpcNum call BoundaryXCyc_xyz( xyza_MixRtAl(:,:,:,s) ) ! (inout) call BoundaryYCyc_xyz( xyza_MixRtAl(:,:,:,s) ) ! (inout) call BoundaryZSym_xyz( xyza_MixRtAl(:,:,:,s) ) ! (inout) end do ! 湿潤飽和調節 ! * 蒸気<-->雲の変換を行う. ! ! Moist adjustment. ! * Conversion from vapor to cloud. ! これまでの値を作業配列に保管 ! Previous values are stored to work area. ! xyz_PotTempWork = xyz_PotTempAl xyza_MixRtWork = xyza_MixRtAl ! 湿潤調節法を適用 ! Moist adjustment is applied. ! call MoistAdjustSvapPress( & & xyz_ExnerNl, & ! (in) & xyz_PotTempAl, & ! (inout) & xyza_MixRtAl & & ) call MoistAdjustNH4SH( & & xyz_ExnerNl, & !(in) & xyz_PotTempAl, & !(inout) & xyza_MixRtAl & !(inout) & ) ! 湿潤調節法による温位と混合比の変化量を保管 ! Adjustment values are stored. ! call StorePotTempCond( & & (xyz_PotTempAl - xyz_PotTempWork) / DelTimeLFrog(t) & ! (in) & ) call StoreMixRtCond( & & (xyza_MixRtAl - xyza_MixRtWork) / DelTimeLFrog(t) & ! (in) & ) ! 境界条件 Boundary condition ! do s = 1, SpcNum call BoundaryXCyc_xyz( xyza_MixRtAl(:,:,:,s) ) call BoundaryYCyc_xyz( xyza_MixRtAl(:,:,:,s) ) call BoundaryZSym_xyz( xyza_MixRtAl(:,:,:,s) ) end do ! 暖かい雨のパラメタリゼーション. ! * 蒸気<-->雨 の変換を行う ! ! Warm rain parameterization. ! * Conversion from rain to vapor. !これまでの値を作業配列に保管 ! Previous values are stored to work area. ! xyz_PotTempWork = xyz_PotTempAl xyza_MixRtWork = xyza_MixRtAl ! 雨から蒸気への混合比変化を求める ! * 温位の計算において, 混合比変化が必要となるため, ! 混合比変化を 1 つの配列として用意する. ! ! Conversion values are calculated. ! xyza_DelMixRt = 0.0d0 xyza_DelMixRt = & & ( & & + xyza_Rain2Gas( & & xyz_ExnerNl, xyz_PotTempAl, xyza_MixRtAl, DelTimeLFrog(t)& & ) & & + xyza_Rain2GasNH4SH( & & xyz_ExnerNl, xyz_PotTempAl, xyza_MixRtAl, DelTimeLFrog(t)& & ) & & ) ! 温位の計算. 雨から蒸気への変換に伴う潜熱・反応熱を追加. ! ! xyz_PotTempAl = & & xyz_PotTempWork & & + ( & & + xyz_Rain2GasHeat( xyz_PotTempAl, xyz_ExnerNl, xyza_DelMixRt ) & & + xyz_Rain2GasHeatNH4SH(xyz_ExnerNl, xyza_DelMixRt) & & ) ! 混合比の計算. 雨から蒸気への変換分を追加 ! ! xyza_MixRtAl = xyza_MixRtWork + xyza_DelMixRt !雨から蒸気の値の保管 ! ! call StorePotTempCond( & & (xyz_PotTempAl - xyz_PotTempWork) / DelTimeLFrog(t) & ! (in) & ) call StoreMixRtCond( & & xyza_DelMixRt / DelTimeLFrog(t) & ! (in) & ) ! 境界条件 Boundary condition ! call BoundaryXCyc_xyz( xyz_PotTempAl ) ! (inout) call BoundaryYCyc_xyz( xyz_PotTempAl ) ! (inout) call BoundaryZSym_xyz( xyz_PotTempAl ) ! (inout) do s = 1, SpcNum call BoundaryXCyc_xyz( xyza_MixRtAl(:,:,:,s) ) call BoundaryYCyc_xyz( xyza_MixRtAl(:,:,:,s) ) call BoundaryZSym_xyz( xyza_MixRtAl(:,:,:,s) ) end do ! 速度の移流拡散. ! 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_BuoyMolWt(xyza_MixRtNl) & & + xyr_BuoyDrag(xyza_MixRtNl) & & + 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 ! 短い時間ステップの時間積分. オイラー法を利用. ! Time integration with short time step. ! Euler: do tau = 1, NstepEuler(t) ! 速度 u の計算. ! Time integration horizontal velocity (u). ! 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 (v). ! 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 ) ! (inout) call BoundaryYCyc_xyz( xyz_ExnerAs ) ! (inout) call BoundaryZSym_xyz( xyz_ExnerAs ) ! (inout) ! 速度 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, & ! (in) & xyz_ExnerNl, & ! (inout) & xyz_ExnerBl & ! (in) & ) call AsselinFilter( & & pyz_VelXAl, & ! (in) & pyz_VelXNl, & ! (inout) & pyz_VelXBl & ! (in) & ) call AsselinFilter( & & xqz_VelYAl, & ! (in) & xqz_VelYNl, & ! (inout) & xqz_VelYBl & ! (in) & ) call AsselinFilter( & & xyr_VelZAl, & ! (in) & xyr_VelZNl, & ! (inout) & xyr_VelZBl & ! (in) & ) call AsselinFilter( & & xyz_PotTempAl, & ! (in) & xyz_PotTempNl, & ! (inout) & xyz_PotTempBl & ! (in) & ) call AsselinFilter( & & xyz_KmAl, & ! (in) & xyz_KmNl, & ! (inout) & xyz_KmBl & ! (in) & ) do s = 1, SpcNum call AsselinFilter(xyza_MixRtAl(:,:,:,s), & & xyza_MixRtNl(:,:,:,s), & & xyza_MixRtBl(:,:,:,s)) end do ! スポンジ層 ! Numerical dumping. ! call DampSponge_pyz( pyz_VelXAl, pyz_VelXBl, DelTimeLFrog(t) ) call DampSponge_xqz( xqz_VelYAl, xqz_VelYBl, DelTimeLFrog(t) ) call DampSponge_xyr( xyr_VelZAl, xyr_VelZBl, DelTimeLFrog(t) ) call DampSponge_xyz( xyz_ExnerAl, xyz_ExnerBl, DelTimeLFrog(t) ) call DampSponge_xyz( xyz_PotTempAl, xyz_PotTempBl, DelTimeLFrog(t) ) call DampSponge_xyz( xyz_KmAl, xyz_KmBl, DelTimeLFrog(t) ) ! 混合比がゼロ以下にならないための処置. ! Negative values mixing ratios are corrected. ! ! xyza_MixRtWork = xyza_MixRtAl ! xyza_MixRtAl = max( - xyza_MixRtBasicZ, xyza_MixRtWork ) ! call StoreMixRtFill2( & ! & (xyza_MixRtAl - xyza_MixRtWork) / DelTimeLFrog(t) & ! (in) ! & ) ! 安定度の計算. ! Calculation static stability. ! call ECCM_Stab( & & xyz_PotTempAl, & ! (in) & xyz_ExnerAl, & ! (in) & xyza_MixRtAl & ! (in) & ) ! 長い時間ステップのループを回すための処置. ! Renew prognostic variables for next long time step integration. ! xyz_PotTempBl = xyz_PotTempNl xyza_MixRtBl = xyza_MixRtNl 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 xyza_MixRtNl = xyza_MixRtAl 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 & ! (in) & ) call CFLCheckTimeLongVelY( & & xqz_VelYNl & ! (in) & ) call CFLCheckTimeLongVelZ( & & xyr_VelZNl & ! (in) & ) !ヒストリファイル出力 ! Out put to history file. ! call HistoryFile_Output( & & Time, & ! (in) & xyz_PotTempNl, & ! (in) & xyz_ExnerNl, & ! (in) & pyz_VelXNl, & ! (in) & xqz_VelYNl, & ! (in) & xyr_VelZNl, & ! (in) & xyza_MixRtNl, & ! (in) & xyz_KmNl, & ! (in) & xyz_KhNl & ! (in) & ) ! 積算値のクリア. ! 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, & ! (in) & xyz_PotTempBl, & ! (in) & xyz_ExnerBl, & ! (in) & pyz_VelXBl, & ! (in) & xqz_VelYBl, & ! (in) & xyr_VelZBl, & ! (in) & xyza_MixRtBl, & ! (in) & xyz_KmBl, & ! (in) & xyz_KhBl & ! (in) & ) call ReStartFile_OutPut( & ! (in) & Time, & ! (in) & xyz_PotTempNl, & ! (in) & xyz_ExnerNl, & ! (in) & pyz_VelXNl, & ! (in) & xqz_VelYNl, & ! (in) & xyr_VelZNl, & ! (in) & xyza_MixRtNl, & ! (in) & xyz_KmNl, & ! (in) & xyz_KhNl & ! (in) & ) call ReStartFile_Close !---------------------------------------------------- ! MPI END !---------------------------------------------------- ! call mpii_end call MPIWrapperFinalize 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_MixRtWork(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax,SpcNum), & & 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),& ! & xyza_DelMixRt(DimXMin:DimXMax,DimYMin:DimYMax,DimZMin:DimZMax,SpcNum), & ! & 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 xyz_PotTempWork = 0.0d0 xyza_MixRtWork = 0.0d0 pyz_AccelVelXNl = 0.0d0 xqz_AccelVelYNl = 0.0d0 xyr_AccelVelZNl = 0.0d0 xyza_DelMixRt = 0.0d0 DelTimeEuler = 0.0d0 DelTimeLFrog = 0.0d0; NStepEuler = 0.0d0 NstepLFrog = 0.0d0 end subroutine VariableAllocate !----------------------------------------------------------------------- end program ararempi