| Class | major_comp_phase_change |
| In: |
saturate/major_comp_phase_change.f90
|
Note that Japanese and English are described in parallel.
| !$ ! DryConvAdjust : | 乾燥対流調節 |
| !$ ! ———— : | ———— |
| !$ ! DryConvAdjust : | Dry convective adjustment |
| Subroutine : | |||
| xyr_Press(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in )
| ||
| xyz_Press(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in )
| ||
| xy_Ps(0:imax-1, 1:jmax) : | real(DP), intent(inout)
| ||
| xyz_Temp(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout)
| ||
| xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(inout) | ||
| xyz_U(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout) | ||
| xyz_V(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout) | ||
| xy_SurfMajCompIce(0:imax-1, 1:jmax) : | real(DP), intent(inout)
|
CO2 相変化
CO2 phase change
subroutine MajorCompPhaseChangeInAtm( xyr_Press, xyz_Press, xy_Ps, xyz_Temp, xyzf_QMix, xyz_U, xyz_V, xy_SurfMajCompIce )
!
! CO2 相変化
!
! CO2 phase change
!
! モジュール引用 ; USE statements
!
! 時刻管理
! Time control
!
use timeset, only: DelTime, TimeN, TimesetClockStart, TimesetClockStop
! ヒストリデータ出力
! History data output
!
use gtool_historyauto, only: HistoryAutoPut
! 組成に関わる配列の設定
! Settings of array for atmospheric composition
!
use composition, only : ncmax, IndexTKE
! 物理定数設定
! Physical constants settings
!
use constants, only: Grav, CpDry
! $ C_p $ [J kg-1 K-1].
! 乾燥大気の定圧比熱.
! Specific heat of air at constant pressure
! 温度の半整数σレベルの補間, 気圧と高度の算出
! Interpolate temperature on half sigma level,
! and calculate pressure and height
!
use auxiliary, only: AuxVars
! 主成分相変化
! Phase change of atmospheric major component
!
use saturate_major_comp, only : SaturateMajorCompCalcCondTemp, SaturateMajorCompInqLatentHeat
! 質量の補正
! Mass fixer
!
use mass_fixer, only: MassFixerColumn
! 宣言文 ; Declaration statements
!
implicit none
real(DP), intent(in ):: xyr_Press (0:imax-1, 1:jmax, 0:kmax)
! $ \hat{p} $ . 気圧 (半整数レベル).
! Air pressure (half level)
real(DP), intent(in ):: xyz_Press (0:imax-1, 1:jmax, 1:kmax)
! $ p $ . 気圧 (整数レベル).
! Air pressure (full level)
real(DP), intent(inout):: xy_Ps (0:imax-1, 1:jmax)
! $ T $ . 温度. Temperature
real(DP), intent(inout):: xyz_Temp (0:imax-1, 1:jmax, 1:kmax)
! $ T $ . 温度. Temperature
real(DP), intent(inout):: xyzf_QMix (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
real(DP), intent(inout):: xyz_U (0:imax-1, 1:jmax, 1:kmax)
real(DP), intent(inout):: xyz_V (0:imax-1, 1:jmax, 1:kmax)
real(DP), intent(inout):: xy_SurfMajCompIce(0:imax-1, 1:jmax)
!
! Surface major component ice amount
! 作業変数
! Work variables
!
real(DP):: LatentHeatMajCompSubl
real(DP):: xy_PsB (0:imax-1, 1:jmax)
real(DP):: xy_PsA (0:imax-1, 1:jmax)
real(DP):: xyr_PressB (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyr_PressA (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyz_TempB (0:imax-1, 1:jmax, 1:kmax)
! 調節前の温度.
! Temperature before adjustment
real(DP):: xyzf_QMixB (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
real(DP):: xyz_DelAtmMass (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xy_TempTmp (0:imax-1, 1:jmax)
real(DP):: xy_DTempDtSubl (0:imax-1, 1:jmax)
real(DP):: xy_DTempDtCond (0:imax-1, 1:jmax)
real(DP):: xyr_AtmMassFallFlux (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyz_TempTmp (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyz_QH2OVapTmp (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyz_DTempDt (0:imax-1, 1:jmax, 1:kmax)
! 温度変化率.
! Temperature tendency
real(DP):: xy_DSurfMajCompIceDt(0:imax-1, 1:jmax)
!
! Surface major component ice tendency
real(DP):: xy_DPsDt (0:imax-1, 1:jmax)
real(DP):: xy_DSurfMajCompIceDtOneLayer(0:imax-1, 1:jmax)
real(DP):: xyr_DPressDt (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyz_TempCond (0:imax-1, 1:jmax, 1:kmax)
integer :: mmax
real(DP):: xyza_Array (0:imax-1, 1:jmax, 1:kmax, 1:ncmax+1+1)
logical :: a_FlagSurfaceSink (1:ncmax+1+1)
real(DP):: xyra_MassFlow (0:imax-1, 1:jmax, 0:kmax, 1:ncmax+1+1)
real(DP):: xyrf_MassFlow (0:imax-1, 1:jmax, 0:kmax, 1:ncmax)
real(DP):: xyr_MomXFlow (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyr_MomYFlow (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyz_DelAtmMassB (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyz_DelAtmMassA (0:imax-1, 1:jmax, 1:kmax)
integer:: i ! 経度方向に回る DO ループ用作業変数
! Work variables for DO loop in longitude
integer:: j ! 緯度方向に回る DO ループ用作業変数
! Work variables for DO loop in latitude
integer:: k ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
integer:: m
integer:: n
logical :: FlagCheckPs
! 実行文 ; Executable statement
!
! 初期化
! Initialization
!
if ( .not. major_comp_phase_change_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
if ( .not. FlagMajCompPhaseChange ) return
! 計算時間計測開始
! Start measurement of computation time
!
call TimesetClockStart( module_name )
! Set latent heat
LatentHeatMajCompSubl = SaturateMajorCompInqLatentHeat()
FlagCheckPs = .false.
do j = 1, jmax
do i = 0, imax-1
if ( xyr_Press(i,j,0) > 1.0e4_DP ) then
FlagCheckPs = .true.
end if
end do
end do
if ( FlagCheckPs ) then
call MessageNotify( 'W', module_name, 'Surface pressure is greater than 10000 Pa.' )
end if
! Store variables
!
xyz_TempB = xyz_Temp
xyzf_QMixB = xyzf_QMix
call SaturateMajorCompCalcCondTemp( xyz_Press, xyz_TempCond )
do k = 1, kmax
xyz_DelAtmMass(:,:,k) = ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav
end do
xyr_AtmMassFallFlux(:,:,kmax) = 0.0_DP
do k = kmax, 1, -1
! sublimation of falling condensate
!!$ xy_DTempDtSubl = &
!!$ & - LatentHeatMajCompSubl * xyr_AtmMassFallFlux(:,:,k) &
!!$ & / ( CpDry * xyz_DelAtmMass(:,:,k) )
!!$ xyz_Temp(:,:,k) = xyz_Temp(:,:,k) + xy_DTempDtSubl * ( 2.0_DP * DelTime )
!!$ xyr_AtmMassFallFlux(:,:,k) = 0.0_DP
! condensation
xy_TempTmp = xyz_Temp(:,:,k)
xyz_Temp(:,:,k) = max( xyz_TempCond(:,:,k), xyz_Temp(:,:,k) )
xy_DTempDtCond = ( xyz_Temp(:,:,k) - xy_TempTmp ) / ( 2.0_DP * DelTime )
xy_DSurfMajCompIceDtOneLayer = CpDry * xyz_DelAtmMass(:,:,k) * xy_DTempDtCond / LatentHeatMajCompSubl
xyr_AtmMassFallFlux(:,:,k-1) = xyr_AtmMassFallFlux(:,:,k) + xy_DSurfMajCompIceDtOneLayer
end do
xyr_DPressDt = - xyr_AtmMassFallFlux * Grav
xy_DPsDt = xyr_DPressDt(:,:,0)
xy_DSurfMajCompIceDt = - xy_DPsDt / Grav
! packing
if ( FlagModMom ) then
mmax = ncmax + 1 + 1
else
mmax = ncmax
end if
do m = 1, ncmax
n = m
xyza_Array(:,:,:,m) = xyzf_QMix(:,:,:,n)
if ( n == IndexTKE ) then
a_FlagSurfaceSink(m) = .true.
else
a_FlagSurfaceSink(m) = .false.
end if
end do
if ( FlagModMom ) then
m = ncmax
m = m + 1
xyza_Array(:,:,:,m) = xyz_U
a_FlagSurfaceSink(m) = .true.
m = m + 1
xyza_Array(:,:,:,m) = xyz_V
a_FlagSurfaceSink(m) = .true.
end if
call MajorCompPhaseChangeCalcFlow( xyr_Press, xyr_DPressDt, mmax, a_FlagSurfaceSink(1:mmax), xyza_Array(:,:,:,1:mmax), xyra_MassFlow(:,:,:,1:mmax) )
! unpacking
do m = 1, ncmax
xyrf_MassFlow(:,:,:,m) = xyra_MassFlow(:,:,:,m)
end do
if ( FlagModMom ) then
m = ncmax
m = m + 1
xyr_MomXFlow = xyra_MassFlow(:,:,:,m)
m = m + 1
xyr_MomYFlow = xyra_MassFlow(:,:,:,m)
end if
! Adjustment
! preparation
xy_PsB = xy_Ps
xy_PsA = xy_PsB + xy_DPsDt * ( 2.0_DP * DelTime )
! 温度の半整数σレベルの補間, 気圧と高度の算出
! Interpolate temperature on half sigma level,
! and calculate pressure and height
!
xyz_TempTmp = 300.0_DP
xyz_QH2OVapTmp = 0.0_DP
call AuxVars( xy_PsB, xyz_TempTmp, xyz_QH2OVapTmp, xyr_Press = xyr_PressB )
call AuxVars( xy_PsA, xyz_TempTmp, xyz_QH2OVapTmp, xyr_Press = xyr_PressA )
do k = 1, kmax
xyz_DelAtmMassB(:,:,k) = ( xyr_PressB(:,:,k-1) - xyr_PressB(:,:,k) ) / Grav
xyz_DelAtmMassA(:,:,k) = ( xyr_PressA(:,:,k-1) - xyr_PressA(:,:,k) ) / Grav
end do
! Atmospheric composition
do n = 1, ncmax
do k = 1, kmax
xyzf_QMix(:,:,k,n) = ( xyz_DelAtmMassB(:,:,k) * xyzf_QMix(:,:,k,n) - ( xyrf_MassFlow(:,:,k,n) - xyrf_MassFlow(:,:,k-1,n) ) ) / xyz_DelAtmMassA(:,:,k)
end do
end do
if ( FlagModMom ) then
do k = 1, kmax
! Zonal wind
xyz_U(:,:,k) = ( xyz_DelAtmMassB(:,:,k) * xyz_U(:,:,k) - ( xyr_MomXFlow(:,:,k) - xyr_MomXFlow(:,:,k-1) ) ) / xyz_DelAtmMassA(:,:,k)
! Meridional wind
xyz_V(:,:,k) = ( xyz_DelAtmMassB(:,:,k) * xyz_V(:,:,k) - ( xyr_MomYFlow(:,:,k) - xyr_MomYFlow(:,:,k-1) ) ) / xyz_DelAtmMassA(:,:,k)
end do
end if
! Surface major component ice adjustment
xy_SurfMajCompIce = xy_SurfMajCompIce + xy_DSurfMajCompIceDt * ( 2.0_DP * DelTime )
! Surface pressure adjustment
xy_Ps = xy_PsA
call AuxVars( xy_PsA, xyz_TempTmp, xyz_QH2OVapTmp, xyr_Press = xyr_PressA )
! 成分の質量の補正
! Fix masses of constituents
!
call MassFixerColumn( xyr_PressA, xyzf_QMix )
! Check
call MajorCompPhaseChangeConsChk( a_FlagSurfaceSink, xyz_DelAtmMassB, xyzf_QMixB, xyz_DelAtmMassA, xyzf_QMix )
! ヒストリデータ出力
! History data output
!
call HistoryAutoPut( TimeN, 'DTempDtMajCompPhaseChange', xyz_DTempDt )
! 計算時間計測一時停止
! Pause measurement of computation time
!
call TimesetClockStop( module_name )
end subroutine MajorCompPhaseChangeInAtm
| Subroutine : | |
| ArgFlagMajCompPhaseChange : | logical , intent(in) |
| CondMajCompName : | character(*), intent(in) |
major_comp_phase_change モジュールの初期化を行います. NAMELIST#major_comp_phase_change_nml の読み込みはこの手続きで行われます.
"major_comp_phase_change" module is initialized. "NAMELIST#major_comp_phase_change_nml" is loaded in this procedure.
This procedure input/output NAMELIST#major_comp_phase_change_nml .
subroutine MajorCompPhaseChangeInit( ArgFlagMajCompPhaseChange, CondMajCompName )
!
! major_comp_phase_change モジュールの初期化を行います.
! NAMELIST#major_comp_phase_change_nml の読み込みはこの手続きで行われます.
!
! "major_comp_phase_change" module is initialized.
! "NAMELIST#major_comp_phase_change_nml" is loaded in this procedure.
!
! モジュール引用 ; USE statements
!
! NAMELIST ファイル入力に関するユーティリティ
! Utilities for NAMELIST file input
!
use namelist_util, only: namelist_filename, NmlutilMsg, NmlutilAryValid
! ファイル入出力補助
! File I/O support
!
use dc_iounit, only: FileOpen
! 種別型パラメタ
! Kind type parameter
!
use dc_types, only: STDOUT ! 標準出力の装置番号. Unit number of standard output
! 文字列操作
! Character handling
!
use dc_string, only: StoA
! ヒストリデータ出力
! History data output
!
use gtool_historyauto, only: HistoryAutoAddVariable
! 補助的な変数を計算するサブルーチン・関数群
! Subroutines and functions for calculating auxiliary variables
!
use auxiliary, only : AuxVarsInit
! 主成分相変化
! Phase change of atmospheric major component
!
use saturate_major_comp, only : SaturateMajorCompInit
! 質量の補正
! Mass fixer
!
use mass_fixer, only : MassFixerInit
! 宣言文 ; Declaration statements
!
implicit none
logical , intent(in) :: ArgFlagMajCompPhaseChange
character(*), intent(in) :: CondMajCompName
integer:: unit_nml ! NAMELIST ファイルオープン用装置番号.
! Unit number for NAMELIST file open
integer:: iostat_nml ! NAMELIST 読み込み時の IOSTAT.
! IOSTAT of NAMELIST read
! NAMELIST 変数群
! NAMELIST group name
!
namelist /major_comp_phase_change_nml/ FlagModMom
! デフォルト値については初期化手続 "major_comp_phase_change#MajorCompPhaseChangeInit"
! のソースコードを参照のこと.
!
! Refer to source codes in the initialization procedure
! "major_comp_phase_change#MajorCompPhaseChangeInit" for the default values.
!
! 実行文 ; Executable statement
!
if ( major_comp_phase_change_inited ) return
FlagMajCompPhaseChange = ArgFlagMajCompPhaseChange
! デフォルト値の設定
! Default values settings
!
FlagModMom = .false.
! NAMELIST の読み込み
! NAMELIST is input
!
if ( trim(namelist_filename) /= '' ) then
call FileOpen( unit_nml, namelist_filename, mode = 'r' ) ! (in)
rewind( unit_nml )
read( unit_nml, nml = major_comp_phase_change_nml, iostat = iostat_nml ) ! (out)
close( unit_nml )
call NmlutilMsg( iostat_nml, module_name ) ! (in)
! if ( iostat_nml == 0 ) write( STDOUT, nml = cumulus_adjust_nml )
end if
if ( FlagMajCompPhaseChange ) then
! 主成分相変化
! Phase change of atmospheric major component
!
call SaturateMajorCompInit( CondMajCompName )
end if
! 補助的な変数を計算するサブルーチン・関数群
! Subroutines and functions for calculating auxiliary variables
!
call AuxVarsInit
! 質量の補正
! Mass fixer
!
call MassFixerInit
! ヒストリデータ出力のためのへの変数登録
! Register of variables for history data output
!
call HistoryAutoAddVariable( 'DSurfTempDtMajCompPhaseChange', (/ 'lon ', 'lat ', 'time' /), 'heating by major component phase change', 'K s-1' )
call HistoryAutoAddVariable( 'DTempDtMajCompPhaseChange', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'heating by major component phase change', 'K s-1' )
! 印字 ; Print
!
call MessageNotify( 'M', module_name, '----- Initialization Messages -----' )
call MessageNotify( 'M', module_name, ' FlagModMom = %b', l = (/ FlagModMom /) )
call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) )
major_comp_phase_change_inited = .true.
end subroutine MajorCompPhaseChangeInit
| Subroutine : | |||
| xy_DPsDt(0:imax-1, 1:jmax) : | real(DP), intent(in ) | ||
| xy_DSurfMajCompIceDt(0:imax-1, 1:jmax) : | real(DP), intent(in ) | ||
| xy_Ps(0:imax-1, 1:jmax) : | real(DP), intent(inout) | ||
| xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(inout) | ||
| xyz_U(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout) | ||
| xyz_V(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout) | ||
| xy_SurfMajCompIce(0:imax-1, 1:jmax) : | real(DP), intent(inout)
|
CO2 相変化
CO2 phase change
subroutine MajorCompPhaseChangeOnGround( xy_DPsDt, xy_DSurfMajCompIceDt, xy_Ps, xyzf_QMix, xyz_U, xyz_V, xy_SurfMajCompIce )
!
! CO2 相変化
!
! CO2 phase change
!
! モジュール引用 ; USE statements
!
! 時刻管理
! Time control
!
use timeset, only: DelTime, TimeN, TimesetClockStart, TimesetClockStop
! ヒストリデータ出力
! History data output
!
use gtool_historyauto, only: HistoryAutoPut
! 組成に関わる配列の設定
! Settings of array for atmospheric composition
!
use composition, only : ncmax, IndexTKE
! 物理定数設定
! Physical constants settings
!
use constants, only: Grav, CpDry
! $ C_p $ [J kg-1 K-1].
! 乾燥大気の定圧比熱.
! Specific heat of air at constant pressure
! 温度の半整数σレベルの補間, 気圧と高度の算出
! Interpolate temperature on half sigma level,
! and calculate pressure and height
!
use auxiliary, only: AuxVars
! 質量の補正
! Mass fixer
!
use mass_fixer, only: MassFixerColumn
! 宣言文 ; Declaration statements
!
implicit none
real(DP), intent(in ):: xy_DPsDt (0:imax-1, 1:jmax)
real(DP), intent(in ):: xy_DSurfMajCompIceDt(0:imax-1, 1:jmax)
real(DP), intent(inout):: xy_Ps (0:imax-1, 1:jmax)
real(DP), intent(inout):: xyzf_QMix (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
real(DP), intent(inout):: xyz_U (0:imax-1, 1:jmax, 1:kmax)
real(DP), intent(inout):: xyz_V (0:imax-1, 1:jmax, 1:kmax)
real(DP), intent(inout):: xy_SurfMajCompIce (0:imax-1, 1:jmax)
!
! Surface major component ice amount
! 作業変数
! Work variables
!
real(DP):: xy_PsB (0:imax-1, 1:jmax)
real(DP):: xy_PsA (0:imax-1, 1:jmax)
real(DP):: xyr_PressB (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyr_PressA (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyzf_QMixB (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
real(DP):: xyz_TempTmp (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyz_QH2OVapTmp (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyr_DPressDt (0:imax-1, 1:jmax, 0:kmax)
integer :: mmax
real(DP):: xyza_Array (0:imax-1, 1:jmax, 1:kmax, 1:ncmax+1+1)
logical :: a_FlagSurfaceSink (1:ncmax+1+1)
real(DP):: xyra_MassFlow (0:imax-1, 1:jmax, 0:kmax, 1:ncmax+1+1)
real(DP):: xyrf_MassFlow (0:imax-1, 1:jmax, 0:kmax, 1:ncmax)
real(DP):: xyr_MomXFlow (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyr_MomYFlow (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyz_DelAtmMassB (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyz_DelAtmMassA (0:imax-1, 1:jmax, 1:kmax)
integer:: i ! 経度方向に回る DO ループ用作業変数
! Work variables for DO loop in longitude
integer:: j ! 緯度方向に回る DO ループ用作業変数
! Work variables for DO loop in latitude
integer:: k ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
integer:: m
integer:: n
! 実行文 ; Executable statement
!
! 初期化
! Initialization
!
if ( .not. major_comp_phase_change_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
if ( .not. FlagMajCompPhaseChange ) return
! 計算時間計測開始
! Start measurement of computation time
!
call TimesetClockStart( module_name )
xy_PsB = xy_Ps
xy_PsA = xy_PsB + xy_DPsDt * ( 2.0_DP * DelTime )
xyzf_QMixB = xyzf_QMix
xyz_TempTmp = 300.0_DP
xyz_QH2OVapTmp = 0.0_DP
! 温度の半整数σレベルの補間, 気圧と高度の算出
! Interpolate temperature on half sigma level,
! and calculate pressure and height
!
call AuxVars( xy_PsB, xyz_TempTmp, xyz_QH2OVapTmp, xyr_Press = xyr_PressB )
call AuxVars( xy_PsA, xyz_TempTmp, xyz_QH2OVapTmp, xyr_Press = xyr_PressA )
xyr_DPressDt = ( xyr_PressA - xyr_PressB ) / ( 2.0_DP * Deltime )
! packing
if ( FlagModMom ) then
mmax = ncmax + 1 + 1
else
mmax = ncmax
end if
do m = 1, ncmax
n = m
xyza_Array(:,:,:,m) = xyzf_QMix(:,:,:,n)
if ( n == IndexTKE ) then
a_FlagSurfaceSink(m) = .true.
else
a_FlagSurfaceSink(m) = .false.
end if
end do
if ( FlagModMom ) then
m = ncmax
m = m + 1
xyza_Array(:,:,:,m) = xyz_U
a_FlagSurfaceSink(m) = .true.
m = m + 1
xyza_Array(:,:,:,m) = xyz_V
a_FlagSurfaceSink(m) = .true.
end if
call MajorCompPhaseChangeCalcFlow( xyr_PressB, xyr_DPressDt, mmax, a_FlagSurfaceSink(1:mmax), xyza_Array(:,:,:,1:mmax), xyra_MassFlow(:,:,:,1:mmax) )
! unpacking
do m = 1, ncmax
xyrf_MassFlow(:,:,:,m) = xyra_MassFlow(:,:,:,m)
end do
if ( FlagModMom ) then
m = ncmax
m = m + 1
xyr_MomXFlow = xyra_MassFlow(:,:,:,m)
m = m + 1
xyr_MomYFlow = xyra_MassFlow(:,:,:,m)
end if
! Adjustment
! preparation
do k = 1, kmax
xyz_DelAtmMassB(:,:,k) = ( xyr_PressB(:,:,k-1) - xyr_PressB(:,:,k) ) / Grav
xyz_DelAtmMassA(:,:,k) = ( xyr_PressA(:,:,k-1) - xyr_PressA(:,:,k) ) / Grav
end do
! Atmospheric composition
do n = 1, ncmax
do k = 1, kmax
xyzf_QMix(:,:,k,n) = ( xyz_DelAtmMassB(:,:,k) * xyzf_QMix(:,:,k,n) - ( xyrf_MassFlow(:,:,k,n) - xyrf_MassFlow(:,:,k-1,n) ) ) / xyz_DelAtmMassA(:,:,k)
end do
end do
if ( FlagModMom ) then
do k = 1, kmax
! Zonal wind
xyz_U(:,:,k) = ( xyz_DelAtmMassB(:,:,k) * xyz_U(:,:,k) - ( xyr_MomXFlow(:,:,k) - xyr_MomXFlow(:,:,k-1) ) ) / xyz_DelAtmMassA(:,:,k)
! Meridional wind
xyz_V(:,:,k) = ( xyz_DelAtmMassB(:,:,k) * xyz_V(:,:,k) - ( xyr_MomYFlow(:,:,k) - xyr_MomYFlow(:,:,k-1) ) ) / xyz_DelAtmMassA(:,:,k)
end do
end if
! Surface major component ice
xy_SurfMajCompIce = xy_SurfMajCompIce + xy_DSurfMajCompIceDt * ( 2.0_DP * DelTime )
! Surface pressure
xy_Ps = xy_PsA
call AuxVars( xy_PsA, xyz_TempTmp, xyz_QH2OVapTmp, xyr_Press = xyr_PressA )
! 成分の質量の補正
! Fix masses of constituents
!
call MassFixerColumn( xyr_PressA, xyzf_QMix )
! Check
call MajorCompPhaseChangeConsChk( a_FlagSurfaceSink, xyz_DelAtmMassB, xyzf_QMixB, xyz_DelAtmMassA, xyzf_QMix )
! 計算時間計測一時停止
! Pause measurement of computation time
!
call TimesetClockStop( module_name )
end subroutine MajorCompPhaseChangeOnGround
| Variable : | |||
| major_comp_phase_change_inited = .false. : | logical, save, public
|
| Subroutine : | |||
| xyr_Press(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in )
| ||
| xyr_DPressDt(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in ) | ||
| mmax : | integer , intent(in ) | ||
| a_FlagSurfaceSink(1:mmax) : | logical , intent(in ) | ||
| xyza_Array(0:imax-1, 1:jmax, 1:kmax, 1:mmax) : | real(DP), intent(in ) | ||
| xyra_MassFlow(0:imax-1, 1:jmax, 0:kmax, 1:mmax) : | real(DP), intent(out) |
CO2 相変化
CO2 phase change
subroutine MajorCompPhaseChangeCalcFlow( xyr_Press, xyr_DPressDt, mmax, a_FlagSurfaceSink, xyza_Array, xyra_MassFlow )
!
! CO2 相変化
!
! CO2 phase change
!
! モジュール引用 ; USE statements
!
! 時刻管理
! Time control
!
use timeset, only: DelTime ! $ \Delta t $
! 物理定数設定
! Physical constants settings
!
use constants, only: Grav ! $ g $ [m s-2].
! 重力加速度.
! Gravitational acceleration
! 宣言文 ; Declaration statements
!
implicit none
real(DP), intent(in ):: xyr_Press (0:imax-1, 1:jmax, 0:kmax)
! pressure
real(DP), intent(in ):: xyr_DPressDt (0:imax-1, 1:jmax, 0:kmax)
integer , intent(in ):: mmax
logical , intent(in ):: a_FlagSurfaceSink(1:mmax)
real(DP), intent(in ):: xyza_Array (0:imax-1, 1:jmax, 1:kmax, 1:mmax)
real(DP), intent(out):: xyra_MassFlow(0:imax-1, 1:jmax, 0:kmax, 1:mmax)
! 作業変数
! Work variables
!
real(DP):: xyr_DPPress (0:imax-1, 1:jmax, 0:kmax)
! pressure at departure point
real(DP):: DelAtmMass
integer:: i ! 経度方向に回る DO ループ用作業変数
! Work variables for DO loop in longitude
integer:: j ! 緯度方向に回る DO ループ用作業変数
! Work variables for DO loop in latitude
integer:: k ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
integer:: k2 ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
integer:: m
! 実行文 ; Executable statement
!
! 初期化
! Initialization
!
if ( .not. major_comp_phase_change_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
if ( .not. FlagMajCompPhaseChange ) return
xyr_DPPress = xyr_Press + xyr_DPressDt * ( 2.0_DP * DelTime )
! check
do k = 1, kmax
do j = 1, jmax
do i = 0, imax-1
if ( xyr_DPPress(i,j,k-1) < xyr_DPPress(i,j,k) ) then
call MessageNotify( 'E', module_name, 'Order of departure points are inappropriate, P(k=%d)=%f < P(k=%d)=%f.', i = (/ k-1, k /), d = (/ xyr_DPPress(i,j,k-1), xyr_DPPress(i,j,k) /) )
end if
end do
end do
end do
xyra_MassFlow = 0.0_DP
do k = 0, kmax-1
do j = 1, jmax
do i = 0, imax-1
!!$ if ( xyr_DPressDt(i,j,k) >= 0.0_DP ) then
if ( xyr_DPPress(i,j,k) > xyr_Press(i,j,k) ) then
sum_upward_mass_transport : do k2 = k, 1, -1
if ( xyr_DPPress(i,j,k) > xyr_Press(i,j,k2-1) ) then
DelAtmMass = ( xyr_Press(i,j,k2-1) - xyr_Press(i,j,k2) ) / Grav
do m = 1, mmax
xyra_MassFlow(i,j,k,m) = xyra_MassFlow(i,j,k,m) + xyza_Array(i,j,k2,m) * DelAtmMass
end do
else
DelAtmMass = ( xyr_DPPress(i,j,k) - xyr_Press(i,j,k2) ) / Grav
do m = 1, mmax
xyra_MassFlow(i,j,k,m) = xyra_MassFlow(i,j,k,m) + xyza_Array(i,j,k2,m) * DelAtmMass
end do
exit sum_upward_mass_transport
end if
end do sum_upward_mass_transport
else
sum_downward_mass_transport : do k2 = k+1, kmax
if ( xyr_DPPress(i,j,k) < xyr_Press(i,j,k2 ) ) then
DelAtmMass = ( xyr_Press(i,j,k2-1) - xyr_Press(i,j,k2) ) / Grav
do m = 1, mmax
xyra_MassFlow(i,j,k,m) = xyra_MassFlow(i,j,k,m) - xyza_Array(i,j,k2,m) * DelAtmMass
end do
else
DelAtmMass = ( xyr_Press(i,j,k2-1) - xyr_DPPress(i,j,k) ) / Grav
do m = 1, mmax
xyra_MassFlow(i,j,k,m) = xyra_MassFlow(i,j,k,m) - xyza_Array(i,j,k2,m) * DelAtmMass
end do
exit sum_downward_mass_transport
end if
end do sum_downward_mass_transport
end if
end do
end do
end do
do k = 0, 0
do j = 1, jmax
do i = 0, imax-1
! not surface sink
if ( xyr_DPressDt(i,j,k) <= 0.0_DP ) then
do m = 1, mmax
if ( .not. a_FlagSurfaceSink(m) ) then
xyra_MassFlow(i,j,k,m) = 0.0_DP
end if
end do
end if
end do
end do
end do
end subroutine MajorCompPhaseChangeCalcFlow
| Subroutine : | |
| a_FlagSurfaceSink(1:ncmax) : | logical , intent(in) |
| xyz_DelAtmMassB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in) |
| xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, ncmax) : | real(DP), intent(in) |
| xyz_DelAtmMassA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in) |
| xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, ncmax) : | real(DP), intent(in) |
subroutine MajorCompPhaseChangeConsChk( a_FlagSurfaceSink, xyz_DelAtmMassB, xyzf_QMixB, xyz_DelAtmMassA, xyzf_QMixA )
! 組成に関わる配列の設定
! Settings of array for atmospheric composition
!
use composition, only : ncmax
! 物理定数設定
! Physical constants settings
!
use constants, only: Grav, CpDry
! $ C_p $ [J kg-1 K-1].
! 乾燥大気の定圧比熱.
! Specific heat of air at constant pressure
logical , intent(in) :: a_FlagSurfaceSink(1:ncmax)
real(DP), intent(in) :: xyz_DelAtmMassB(0:imax-1, 1:jmax, 1:kmax)
real(DP), intent(in) :: xyzf_QMixB (0:imax-1, 1:jmax, 1:kmax, ncmax)
real(DP), intent(in) :: xyz_DelAtmMassA(0:imax-1, 1:jmax, 1:kmax)
real(DP), intent(in) :: xyzf_QMixA (0:imax-1, 1:jmax, 1:kmax, ncmax)
! Local variables
!
real(DP) :: ValB
real(DP) :: ValA
real(DP) :: xyf_SumB(0:imax-1, 1:jmax, 1:ncmax)
real(DP) :: xyf_SumA(0:imax-1, 1:jmax, 1:ncmax)
real(DP) :: Ratio
integer :: i
integer :: j
integer :: k
integer :: n
xyf_SumB = 0.0_DP
xyf_SumA = 0.0_DP
do n = 1, ncmax
do k = kmax, 1, -1
xyf_SumB(:,:,n) = xyf_SumB(:,:,n) + xyz_DelAtmMassB(:,:,k) * xyzf_QMixB(:,:,k,n)
xyf_SumA(:,:,n) = xyf_SumA(:,:,n) + xyz_DelAtmMassA(:,:,k) * xyzf_QMixA(:,:,k,n)
end do
end do
do n = 1, ncmax
if ( .not. a_FlagSurfaceSink(n) ) then
do j = 1, jmax
do i = 0, imax-1
ValB = xyf_SumB(i,j,n)
ValA = xyf_SumA(i,j,n)
Ratio = ( ValA - ValB ) / ( ValA + 1.0d-100 )
if ( abs( Ratio ) > 1.0d-10 ) then
call MessageNotify( 'M', module_name, 'Mass No. %d is not conserved, %f.', i = (/ n /), d = (/ Ratio /) )
end if
end do
end do
end if
end do
end subroutine MajorCompPhaseChangeConsChk
| Subroutine : | |||
| xyr_Press(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in )
| ||
| xyz_Press(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in )
| ||
| xy_Ps(0:imax-1, 1:jmax) : | real(DP), intent(inout)
| ||
| xyz_Temp(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout)
| ||
| xy_SurfMajCompIce(0:imax-1, 1:jmax) : | real(DP), intent(inout)
|
CO2 相変化
CO2 phase change
subroutine MajorCompPhaseChangeInAtmBK( xyr_Press, xyz_Press, xy_Ps, xyz_Temp, xy_SurfMajCompIce )
!
! CO2 相変化
!
! CO2 phase change
!
! モジュール引用 ; USE statements
!
! 時刻管理
! Time control
!
use timeset, only: DelTime, TimeN, TimesetClockStart, TimesetClockStop
! ヒストリデータ出力
! History data output
!
use gtool_historyauto, only: HistoryAutoPut
! 物理定数設定
! Physical constants settings
!
use constants, only: Grav, CpDry
! $ C_p $ [J kg-1 K-1].
! 乾燥大気の定圧比熱.
! Specific heat of air at constant pressure
! 主成分相変化
! Phase change of atmospheric major component
!
use saturate_major_comp, only : SaturateMajorCompCalcCondTemp, SaturateMajorCompInqLatentHeat
! 宣言文 ; Declaration statements
!
implicit none
real(DP), intent(in ):: xyr_Press (0:imax-1, 1:jmax, 0:kmax)
! $ \hat{p} $ . 気圧 (半整数レベル).
! Air pressure (half level)
real(DP), intent(in ):: xyz_Press (0:imax-1, 1:jmax, 1:kmax)
! $ p $ . 気圧 (整数レベル).
! Air pressure (full level)
real(DP), intent(inout):: xy_Ps (0:imax-1, 1:jmax)
! $ T $ . 温度. Temperature
real(DP), intent(inout):: xyz_Temp (0:imax-1, 1:jmax, 1:kmax)
! $ T $ . 温度. Temperature
real(DP), intent(inout):: xy_SurfMajCompIce(0:imax-1, 1:jmax)
!
! Surface major component ice amount
! 作業変数
! Work variables
!
real(DP):: xyz_TempB (0:imax-1, 1:jmax, 1:kmax)
! 調節前の温度.
! Temperature before adjustment
real(DP):: xyz_DTempDt (0:imax-1, 1:jmax, 1:kmax)
! 温度変化率.
! Temperature tendency
real(DP):: xy_DSurfMajCompIceDt(0:imax-1, 1:jmax)
!
! Surface major component ice tendency
real(DP):: xyz_TempCond (0:imax-1, 1:jmax, 1:kmax)
real(DP):: LatentHeatMajCompSubl
integer:: i ! 経度方向に回る DO ループ用作業変数
! Work variables for DO loop in longitude
integer:: j ! 緯度方向に回る DO ループ用作業変数
! Work variables for DO loop in latitude
integer:: k ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
logical :: FlagCheckPs
! 実行文 ; Executable statement
!
! 初期化
! Initialization
!
if ( .not. major_comp_phase_change_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
if ( .not. FlagMajCompPhaseChange ) return
! 計算時間計測開始
! Start measurement of computation time
!
call TimesetClockStart( module_name )
! Set latent heat
LatentHeatMajCompSubl = SaturateMajorCompInqLatentHeat()
FlagCheckPs = .false.
do j = 1, jmax
do i = 0, imax-1
if ( xyr_Press(i,j,0) > 1.0e4_DP ) then
FlagCheckPs = .true.
end if
end do
end do
if ( FlagCheckPs ) then
call MessageNotify( 'W', module_name, 'Surface pressure is greater than 10000 Pa.' )
end if
! 調節前 "Temp" の保存
! Store "Temp" before adjustment
!
xyz_TempB = xyz_Temp
call SaturateMajorCompCalcCondTemp( xyz_Press, xyz_TempCond )
do k = 1, kmax
do j = 1, jmax
do i = 0, imax-1
if ( xyz_Temp(i,j,k) < xyz_TempCond(i,j,k) ) then
xyz_Temp(i,j,k) = xyz_TempCond(i,j,k)
end if
end do
end do
end do
! 温度変化率
! Calculate temperature tendency
!
xyz_DTempDt = ( xyz_Temp - xyz_TempB ) / ( 2.0_DP * DelTime )
!
! Surface major component ice adjustment
!
xy_DSurfMajCompIceDt = 0.0_DP
do k = kmax, 1, -1
xy_DSurfMajCompIceDt = xy_DSurfMajCompIceDt + CpDry * xyz_DTempDt(:,:,k) * ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav / LatentHeatMajCompSubl
end do
xy_SurfMajCompIce = xy_SurfMajCompIce + xy_DSurfMajCompIceDt * ( 2.0_DP * DelTime )
!
! Surface pressure adjustment
!
xy_Ps = xy_Ps - xy_DSurfMajCompIceDt * Grav * ( 2.0_DP * DelTime )
! ヒストリデータ出力
! History data output
!
call HistoryAutoPut( TimeN, 'DTempDtMajCompPhaseChange', xyz_DTempDt )
! 計算時間計測一時停止
! Pause measurement of computation time
!
call TimesetClockStop( module_name )
end subroutine MajorCompPhaseChangeInAtmBK
| Subroutine : | |||
| xyr_Press(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in )
| ||
| xyz_Press(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in )
| ||
| xy_Ps(0:imax-1, 1:jmax) : | real(DP), intent(inout)
| ||
| xyz_Temp(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout)
| ||
| xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(inout) | ||
| xy_SurfMajCompIce(0:imax-1, 1:jmax) : | real(DP), intent(inout)
|
CO2 相変化
CO2 phase change
subroutine MajorCompPhaseChangeInAtmBK2( xyr_Press, xyz_Press, xy_Ps, xyz_Temp, xyzf_QMix, xy_SurfMajCompIce )
!
! CO2 相変化
!
! CO2 phase change
!
! モジュール引用 ; USE statements
!
! 時刻管理
! Time control
!
use timeset, only: DelTime, TimeN, TimesetClockStart, TimesetClockStop
! ヒストリデータ出力
! History data output
!
use gtool_historyauto, only: HistoryAutoPut
! 組成に関わる配列の設定
! Settings of array for atmospheric composition
!
use composition, only : ncmax, IndexTKE
! 物理定数設定
! Physical constants settings
!
use constants, only: Grav, CpDry
! $ C_p $ [J kg-1 K-1].
! 乾燥大気の定圧比熱.
! Specific heat of air at constant pressure
! 温度の半整数σレベルの補間, 気圧と高度の算出
! Interpolate temperature on half sigma level,
! and calculate pressure and height
!
use auxiliary, only: AuxVars
! 主成分相変化
! Phase change of atmospheric major component
!
use saturate_major_comp, only : SaturateMajorCompCalcCondTemp, SaturateMajorCompInqLatentHeat
! 質量の補正
! Mass fixer
!
use mass_fixer, only: MassFixerColumn
! 宣言文 ; Declaration statements
!
implicit none
real(DP), intent(in ):: xyr_Press (0:imax-1, 1:jmax, 0:kmax)
! $ \hat{p} $ . 気圧 (半整数レベル).
! Air pressure (half level)
real(DP), intent(in ):: xyz_Press (0:imax-1, 1:jmax, 1:kmax)
! $ p $ . 気圧 (整数レベル).
! Air pressure (full level)
real(DP), intent(inout):: xy_Ps (0:imax-1, 1:jmax)
! $ T $ . 温度. Temperature
real(DP), intent(inout):: xyz_Temp (0:imax-1, 1:jmax, 1:kmax)
! $ T $ . 温度. Temperature
real(DP), intent(inout):: xyzf_QMix (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
real(DP), intent(inout):: xy_SurfMajCompIce(0:imax-1, 1:jmax)
!
! Surface major component ice amount
! 作業変数
! Work variables
!
real(DP):: LatentHeatMajCompSubl
real(DP):: xy_PsB (0:imax-1, 1:jmax)
real(DP):: xy_PsA (0:imax-1, 1:jmax)
real(DP):: xyr_PressB (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyr_PressA (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyz_TempB (0:imax-1, 1:jmax, 1:kmax)
! 調節前の温度.
! Temperature before adjustment
real(DP):: xyzf_QMixB (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
real(DP):: xyz_TempTmp (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyz_QH2OVapTmp (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyz_DTempDt (0:imax-1, 1:jmax, 1:kmax)
! 温度変化率.
! Temperature tendency
real(DP):: xy_DSurfMajCompIceDt(0:imax-1, 1:jmax)
!
! Surface major component ice tendency
real(DP):: xy_DPsDt (0:imax-1, 1:jmax)
real(DP):: xy_DSurfMajCompIceDtOneLayer(0:imax-1, 1:jmax)
real(DP):: xyr_DPressDt (0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyz_TempCond (0:imax-1, 1:jmax, 1:kmax)
integer :: mmax
real(DP):: xyza_Array (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
logical :: a_FlagSurfaceSink(1:ncmax)
real(DP):: xyra_MassFlow (0:imax-1, 1:jmax, 0:kmax, 1:ncmax)
real(DP):: xyrf_MassFlow (0:imax-1, 1:jmax, 0:kmax, 1:ncmax)
real(DP):: xyz_DelAtmMassB (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyz_DelAtmMassA (0:imax-1, 1:jmax, 1:kmax)
integer:: i ! 経度方向に回る DO ループ用作業変数
! Work variables for DO loop in longitude
integer:: j ! 緯度方向に回る DO ループ用作業変数
! Work variables for DO loop in latitude
integer:: k ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
integer:: m
integer:: n
logical :: FlagCheckPs
! 実行文 ; Executable statement
!
! 初期化
! Initialization
!
if ( .not. major_comp_phase_change_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
if ( .not. FlagMajCompPhaseChange ) return
! 計算時間計測開始
! Start measurement of computation time
!
call TimesetClockStart( module_name )
! Set latent heat
LatentHeatMajCompSubl = SaturateMajorCompInqLatentHeat()
FlagCheckPs = .false.
do j = 1, jmax
do i = 0, imax-1
if ( xyr_Press(i,j,0) > 1.0e4_DP ) then
FlagCheckPs = .true.
end if
end do
end do
if ( FlagCheckPs ) then
call MessageNotify( 'W', module_name, 'Surface pressure is greater than 10000 Pa.' )
end if
! Store variables
!
xyz_TempB = xyz_Temp
xyzf_QMixB = xyzf_QMix
call SaturateMajorCompCalcCondTemp( xyz_Press, xyz_TempCond )
do k = 1, kmax
do j = 1, jmax
do i = 0, imax-1
if ( xyz_Temp(i,j,k) < xyz_TempCond(i,j,k) ) then
xyz_Temp(i,j,k) = xyz_TempCond(i,j,k)
end if
end do
end do
end do
! 温度変化率
! Calculate temperature tendency
!
xyz_DTempDt = ( xyz_Temp - xyz_TempB ) / ( 2.0_DP * DelTime )
xyr_DPressDt(:,:,kmax) = 0.0_DP
do k = kmax, 1, -1
xy_DSurfMajCompIceDtOneLayer = CpDry * xyz_DTempDt(:,:,k) * ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav / LatentHeatMajCompSubl
xyr_DPressDt(:,:,k-1) = xyr_DPressDt(:,:,k) - xy_DSurfMajCompIceDtOneLayer * Grav
end do
xy_DPsDt = xyr_DPressDt(:,:,0)
xy_DSurfMajCompIceDt = - xy_DPsDt / Grav
! packing
mmax = ncmax
do m = 1, mmax
n = m
xyza_Array(:,:,:,m) = xyzf_QMix(:,:,:,n)
if ( n == IndexTKE ) then
a_FlagSurfaceSink(m) = .true.
else
a_FlagSurfaceSink(m) = .false.
end if
end do
call MajorCompPhaseChangeCalcFlow( xyr_Press, xyr_DPressDt, mmax, a_FlagSurfaceSink, xyza_Array, xyra_MassFlow )
! unpacking
do m = 1, mmax
xyrf_MassFlow(:,:,:,m) = xyra_MassFlow(:,:,:,m)
end do
! Adjustment
! preparation
xy_PsB = xy_Ps
xy_PsA = xy_PsB + xy_DPsDt * ( 2.0_DP * DelTime )
! 温度の半整数σレベルの補間, 気圧と高度の算出
! Interpolate temperature on half sigma level,
! and calculate pressure and height
!
xyz_TempTmp = 300.0_DP
xyz_QH2OVapTmp = 0.0_DP
call AuxVars( xy_PsB, xyz_TempTmp, xyz_QH2OVapTmp, xyr_Press = xyr_PressB )
call AuxVars( xy_PsA, xyz_TempTmp, xyz_QH2OVapTmp, xyr_Press = xyr_PressA )
do k = 1, kmax
xyz_DelAtmMassB(:,:,k) = ( xyr_PressB(:,:,k-1) - xyr_PressB(:,:,k) ) / Grav
xyz_DelAtmMassA(:,:,k) = ( xyr_PressA(:,:,k-1) - xyr_PressA(:,:,k) ) / Grav
end do
! Atmospheric composition
do n = 1, ncmax
do k = 1, kmax
xyzf_QMix(:,:,k,n) = ( xyz_DelAtmMassB(:,:,k) * xyzf_QMix(:,:,k,n) - ( xyrf_MassFlow(:,:,k,n) - xyrf_MassFlow(:,:,k-1,n) ) ) / xyz_DelAtmMassA(:,:,k)
end do
end do
! Surface major component ice adjustment
xy_SurfMajCompIce = xy_SurfMajCompIce + xy_DSurfMajCompIceDt * ( 2.0_DP * DelTime )
! Surface pressure adjustment
xy_Ps = xy_PsA
call AuxVars( xy_PsA, xyz_TempTmp, xyz_QH2OVapTmp, xyr_Press = xyr_PressA )
! 成分の質量の補正
! Fix masses of constituents
!
call MassFixerColumn( xyr_PressA, xyzf_QMix )
! Check
call MajorCompPhaseChangeConsChk( a_FlagSurfaceSink, xyz_DelAtmMassB, xyzf_QMixB, xyz_DelAtmMassA, xyzf_QMix )
! ヒストリデータ出力
! History data output
!
call HistoryAutoPut( TimeN, 'DTempDtMajCompPhaseChange', xyz_DTempDt )
! 計算時間計測一時停止
! Pause measurement of computation time
!
call TimesetClockStop( module_name )
end subroutine MajorCompPhaseChangeInAtmBK2
| Subroutine : | |||
| xyr_Press(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in )
| ||
| xyz_Press(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in )
| ||
| xyz_Height(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in ) | ||
| xy_Ps(0:imax-1, 1:jmax) : | real(DP), intent(inout)
| ||
| xyz_Temp(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout)
| ||
| xy_SurfMajCompIce(0:imax-1, 1:jmax) : | real(DP), intent(inout)
|
CO2 相変化
CO2 phase change
subroutine MajorCompPhaseChangeInAtmTest( xyr_Press, xyz_Press, xyz_Height, xy_Ps, xyz_Temp, xy_SurfMajCompIce )
!
! CO2 相変化
!
! CO2 phase change
!
! モジュール引用 ; USE statements
!
! 時刻管理
! Time control
!
use timeset, only: DelTime, TimeN, TimesetClockStart, TimesetClockStop
! ヒストリデータ出力
! History data output
!
use gtool_historyauto, only: HistoryAutoPut
! 物理定数設定
! Physical constants settings
!
use constants, only: Grav, CpDry
! $ C_p $ [J kg-1 K-1].
! 乾燥大気の定圧比熱.
! Specific heat of air at constant pressure
! 主成分相変化
! Phase change of atmospheric major component
!
use saturate_major_comp, only : SaturateMajorCompCalcCondTemp, SaturateMajorCompInqLatentHeat
! 宣言文 ; Declaration statements
!
implicit none
real(DP), intent(in ):: xyr_Press (0:imax-1, 1:jmax, 0:kmax)
! $ \hat{p} $ . 気圧 (半整数レベル).
! Air pressure (half level)
real(DP), intent(in ):: xyz_Press (0:imax-1, 1:jmax, 1:kmax)
! $ p $ . 気圧 (整数レベル).
! Air pressure (full level)
real(DP), intent(in ):: xyz_Height(0:imax-1, 1:jmax, 1:kmax)
!
!
real(DP), intent(inout):: xy_Ps (0:imax-1, 1:jmax)
! $ T $ . 温度. Temperature
real(DP), intent(inout):: xyz_Temp (0:imax-1, 1:jmax, 1:kmax)
! $ T $ . 温度. Temperature
real(DP), intent(inout):: xy_SurfMajCompIce(0:imax-1, 1:jmax)
!
! Surface major component ice amount
! 作業変数
! Work variables
!
real(DP):: xyz_TempB (0:imax-1, 1:jmax, 1:kmax)
! 調節前の温度.
! Temperature before adjustment
real(DP):: xyz_DelAtmMass (0:imax-1, 1:jmax, 1:kmax)
!
! Atmospheric mass in a layer
real(DP):: xy_FallingIce (0:imax-1, 1:jmax)
!
!
real(DP):: xyz_DelMajCompIce (0:imax-1, 1:jmax, 1:kmax)
!
!
real(DP):: xy_DelSurfMajCompIce(0:imax-1, 1:jmax)
!
!
real(DP):: xyz_DTempDt (0:imax-1, 1:jmax, 1:kmax)
! 温度変化率.
! Temperature tendency
real(DP):: xy_DSurfMajCompIceDt(0:imax-1, 1:jmax)
!
! Surface major component ice tendency
real(DP):: xyz_TempCond (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xy_SurfTempCond(0:imax-1, 1:jmax)
real(DP):: SpecHeatCO2Ice
real(DP):: LatentHeatMajCompSubl
integer:: i ! 経度方向に回る DO ループ用作業変数
! Work variables for DO loop in longitude
integer:: j ! 緯度方向に回る DO ループ用作業変数
! Work variables for DO loop in latitude
integer:: k ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
logical :: FlagCheckPs
! 実行文 ; Executable statement
!
! 初期化
! Initialization
!
if ( .not. major_comp_phase_change_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
if ( .not. FlagMajCompPhaseChange ) return
! 計算時間計測開始
! Start measurement of computation time
!
call TimesetClockStart( module_name )
! Set latent heat
LatentHeatMajCompSubl = SaturateMajorCompInqLatentHeat()
FlagCheckPs = .false.
do j = 1, jmax
do i = 0, imax-1
if ( xyr_Press(i,j,0) > 1.0e4_DP ) then
FlagCheckPs = .true.
end if
end do
end do
if ( FlagCheckPs ) then
call MessageNotify( 'W', module_name, 'Surface pressure is greater than 10000 Pa.' )
end if
! 調節前 "Temp" の保存
! Store "Temp" before adjustment
!
xyz_TempB = xyz_Temp
call SaturateMajorCompCalcCondTemp( xyz_Press, xyz_TempCond )
do k = 1, kmax
xyz_DelAtmMass(:,:,k) = ( xyr_Press(i,j,k-1) - xyr_Press(i,j,k) ) / Grav
end do
k = kmax
do j = 1, jmax
do i = 0, imax-1
if ( xyz_Temp(i,j,k) < xyz_TempCond(i,j,k) ) then
xyz_DelMajCompIce(i,j,k) = CpDry * ( xyz_TempCond(i,j,k) - xyz_Temp(i,j,k) ) * xyz_DelAtmMass(i,j,k) / LatentHeatMajCompSubl
xyz_Temp(i,j,k) = xyz_TempCond(i,j,k)
else
xyz_DelMajCompIce(i,j,k) = 0.0_DP
end if
end do
end do
!
xy_FallingIce = 0.0_DP
do k = kmax-1, 1, -1
xy_FallingIce = xy_FallingIce + xyz_DelMajCompIce(:,:,k+1)
do j = 1, jmax
do i = 0, imax-1
SpecHeatCO2Ice = 349.0_DP + 4.8_DP * xyz_TempCond(i,j,k)
! Forget et al. (1998)
xyz_DelMajCompIce(i,j,k) = CpDry * ( xyz_TempCond(i,j,k) - xyz_Temp(i,j,k) ) * xyz_DelAtmMass(i,j,k) / LatentHeatMajCompSubl - ( Grav * ( xyz_Height(i,j,k+1) - xyz_Height(i,j,k) ) + SpecHeatCO2Ice * ( xyz_TempCond(i,j,k+1) - xyz_TempCond(i,j,k) ) ) / LatentHeatMajCompSubl * xy_FallingIce(i,j)
if ( ( xy_FallingIce(i,j) + xyz_DelMajCompIce(i,j,k) ) >= 0.0_DP ) then
xyz_Temp(i,j,k) = xyz_TempCond(i,j,k)
else
xyz_DelMajCompIce(i,j,k) = - xy_FallingIce(i,j)
xyz_Temp(i,j,k) = xyz_Temp(i,j,k) + ( - LatentHeatMajCompSubl + Grav * ( xyz_Height(i,j,k+1) - xyz_Height(i,j,k) ) + SpecHeatCO2Ice * ( xyz_TempCond(i,j,k+1) - xyz_TempCond(i,j,k) ) ) / ( CpDry * xyz_DelAtmMass(i,j,k) ) * xy_FallingIce(i,j)
end if
end do
end do
end do
! Ice falling on the surface
! This may result in supersaturation in the lowest level.
!
xy_FallingIce = xy_FallingIce + xyz_DelMajCompIce(:,:,1)
k = 1
do j = 1, jmax
do i = 0, imax-1
xy_DelSurfMajCompIce(i,j) = - ( Grav * ( xyz_Height(i,j,1) - 0.0_DP ) ) / LatentHeatMajCompSubl * xy_FallingIce(i,j)
SpecHeatCO2Ice = 349.0_DP + 4.8_DP * xy_SurfTempCond(i,j)
! Forget et al. (1998)
xy_DelSurfMajCompIce(i,j) = - ( Grav * ( xyz_Height(i,j,1) - 0.0_DP ) + SpecHeatCO2Ice * ( xyz_TempCond(i,j,1) - xy_SurfTempCond(i,j) ) ) / LatentHeatMajCompSubl * xy_FallingIce(i,j)
if ( ( xy_FallingIce(i,j) + xy_DelSurfMajCompIce(i,j) ) >= 0.0_DP ) then
! Part of ice sublimes.
! NOTE: In this case, temperature in the lowest layer should be
! condensation temperature. So, actually, the set of temperature is
! meaningless.
xyz_Temp(i,j,k) = xyz_TempCond(i,j,k)
else
! All falling ice sublimes.
! NOTE: The formulation below is different from that by Forget et al.
! (1998). The latent heat and heat by potential energy release and
! heating ice is distributed in the lowest layer in this model, not
! to the soil.
xy_DelSurfMajCompIce(i,j) = - xy_FallingIce(i,j)
xyz_Temp(i,j,k) = xyz_Temp(i,j,k) + ( - LatentHeatMajCompSubl + Grav * ( xyz_Height(i,j,1) - 0.0_DP ) + SpecHeatCO2Ice * ( xyz_TempCond(i,j,1) - xy_SurfTempCond(i,j) ) ) / ( CpDry * xyz_DelAtmMass(i,j,k) ) * xy_FallingIce(i,j)
end if
end do
end do
! 温度変化率
! Calculate temperature tendency
!
xyz_DTempDt = ( xyz_Temp - xyz_TempB ) / ( 2.0_DP * DelTime )
!
! Surface major component ice adjustment
!
xy_DSurfMajCompIceDt = 0.0_DP
do k = kmax, 1, -1
xy_DSurfMajCompIceDt = xy_DSurfMajCompIceDt + xyz_DelMajCompIce(:,:,k)
end do
xy_DSurfMajCompIceDt = xy_DSurfMajCompIceDt + xy_DelSurfMajCompIce(i,j)
!
xy_SurfMajCompIce = xy_SurfMajCompIce + xy_DSurfMajCompIceDt * ( 2.0_DP * DelTime )
!
! Surface pressure adjustment
!
xy_Ps = xy_Ps - xy_DSurfMajCompIceDt * Grav * ( 2.0_DP * DelTime )
! ヒストリデータ出力
! History data output
!
call HistoryAutoPut( TimeN, 'DTempDtMajCompPhaseChange', xyz_DTempDt )
! 計算時間計測一時停止
! Pause measurement of computation time
!
call TimesetClockStop( module_name )
end subroutine MajorCompPhaseChangeInAtmTest
| Constant : | |||
| module_name = ‘major_comp_phase_change‘ : | character(*), parameter
|
| Constant : | |||
| version = ’$Name: $’ // ’$Id: major_comp_phase_change.f90,v 1.3 2014/05/07 09:39:21 murashin Exp $’ : | character(*), parameter
|