| Class | WarmRainPrm |
| In: |
moist/warmrainprm.f90
|
暖かい雨のバルク法を用いた, 水蒸気と雨, 雲と雨の混合比の変換係数を求める.
* 中島健介 (1994) で利用した定式をそのまま利用.
| Subroutine : | |
| cfgfile : | character(*), intent(in) |
This procedure input/output NAMELIST#warmrainprm .
subroutine WarmRainPrm_Init( cfgfile )
!暗黙の型宣言禁止
implicit none
!変数定義
integer :: s
integer :: n1, n2
!入力変数
character(*), intent(in) :: cfgfile
!-----------------------------------------------------------
! NAMELIST から情報を取得
!-----------------------------------------------------------
!
NAMELIST /warmrainprm/ FactorJ, AutoConvTime, MixRt_AutoConvCr
open (10, FILE=cfgfile)
read(10, NML=warmrainprm)
close(10)
!-----------------------------------------------------------
! 雨粒と雲粒と気体の ID の組を作る
!-----------------------------------------------------------
!初期化
allocate( RainSW(SpcNum) )
LoopNum = 0
LoopNum2 = 0
RainSW = 0.0d0
!化学種の中から雨粒を作るものを選び, その配列添え字と分子量を保管.
SelectCloud: do s = 1, SpcNum
!'Cloud' という文字列が含まれるものの個数を数える
n1 = index(SpcWetSymbol(s), '-Cloud' )
if (n1 /= 0) then
LoopNum = LoopNum + 1
GasNum(LoopNum) = minloc(SpcWetID, 1, SpcWetID == ChemData_OneSpcID(SpcWetSymbol(s)(1:n1-3) // '-g'))
CloudNum(LoopNum) = s
end if
!'Rain' という文字列が含まれるものの個数を数える
n2 = index(SpcWetSymbol(s), '-Rain' )
if (n2 /= 0) then
LoopNum2 = LoopNum2 + 1
RainNum(LoopNum2) = s
RainSW(s) = 1.0d0
end if
! NH4SH が存在する場合は LoopNum を 1 つ減らす
if ( trim(SpcWetSymbol(s)) == 'NH4SH-s-Cloud' ) then
LoopNum = LoopNum - 1
end if
end do SelectCloud
!-----------------------------------------------------------
! 硫化アンモニウム, およびアンモニアと硫化水素の ID を取得
! 'Cloud' の方が 'Rain' よりも先に存在すると仮定している
!-----------------------------------------------------------
NH3Num = minloc(SpcWetID, 1, SpcWetID == ChemData_OneSpcID('NH3-g'))
H2SNum = minloc(SpcWetID, 1, SpcWetID == ChemData_OneSpcID('H2S-g'))
NH4SHCloudNum = minloc(SpcWetID, 1, SpcWetID == ChemData_OneSpcID('NH4SH-s'))
NH4SHRainNum = NH4SHCloudNum + 1
!-----------------------------------------------------------
! 確認
!-----------------------------------------------------------
if ( LoopNum == 0 ) then
write(*,*) "WarmRainPrm: CloudNum = 0, please comment out of WarmRainPrm"
! stop
end if
write(*,*) "WarmRainPrm_Init, LoopNum: ", LoopNum
write(*,*) "WarmRainPrm_Init, LoopNum2: ", LoopNum2
write(*,*) "WarmRainPrm_Init, GasNum: ", GasNum
write(*,*) "WarmRainPrm_Init, CloudNum: ", CloudNum
write(*,*) "WarmRainPrm_Init, RainNum: ", RainNum
write(*,*) "WarmRainPrm_Init, RainSW: ", RainSW
write(*,*) "WarmRainPrm_Init, NH3Num: ", NH3Num
write(*,*) "WarmRainPrm_Init, H2SNum: ", H2SNum
write(*,*) "WarmRainPrm_Init, NH4SHNum: ", NH4SHCloudNum
write(*,*) "WarmRainPrm_Init, NH4SHNum: ", NH4SHRainNum
end subroutine WarmRainPrm_Init
| Subroutine : | |
| M1 : | integer, intent(out) |
| MA : | integer, intent(out) |
| MB(10) : | integer, intent(out) |
| M2(10) : | integer, intent(out) |
| M3(10) : | integer, intent(out) |
| M4 : | integer, intent(out) |
| M5 : | integer, intent(out) |
| M6 : | integer, intent(out) |
subroutine WarmRainPrm_prm( M1, MA, MB, M2, M3, M4, M5, M6 )
implicit none
integer, intent(out) :: M1, MA, MB(10), M2(10), M3(10), M4, M5, M6
M1 = LoopNum
MA = LoopNum2
MB = RainNum
M2 = CloudNum
M3 = GasNum
M4 = NH3Num
M5 = H2SNum
M6 = NH4SHCloudNum
write(*,*) M1, MA, MB, M2, M3, M4, M5, M6
end subroutine WarmRainPrm_prm
| Function : | |||
| xz_FallRain(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8)
| ||
| xz_MixRt(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
| ||
| num : | integer, intent(in) |
雨粒の落下による移流を求める.
このルーチンの引数 xz_MixRt は 2 次元配列である. 引数に与えられた混合比に対し, 移流を計算する. もしも, xza_MixRt が雨粒の混合比を表すなら(RainSW = 1)なら その値を出力する. もしも, xza_MixRt が蒸気もしくは雲粒の混合比を表すなら(RainSW = 0) ならば, 落下による移流はゼロとする.
function xz_FallRain( xz_MixRt, num )
!
! 雨粒の落下による移流を求める.
!
! このルーチンの引数 xz_MixRt は 2 次元配列である.
! 引数に与えられた混合比に対し, 移流を計算する.
! もしも, xza_MixRt が雨粒の混合比を表すなら(RainSW = 1)なら
! その値を出力する.
! もしも, xza_MixRt が蒸気もしくは雲粒の混合比を表すなら(RainSW = 0)
! ならば, 落下による移流はゼロとする.
!
!暗黙の型宣言禁止
implicit none
!変数定義
integer, intent(in) :: num
real(8), intent(in) :: xz_MixRt(DimXMin:DimXMax, DimZMin:DimZMax)
!蒸気混合比(擾乱)
real(8) :: xz_MixRtAll(DimXMin:DimXMax, DimZMin:DimZMax)
!蒸気混合比(擾乱 + 平均場)
real(8) :: xz_FallRain(DimXMin:DimXMax, DimZMin:DimZMax)
!雨粒の落下効果
real(8) :: xz_VelZRain(DimXMin:DimXMax, DimZMin:DimZMax)
!雨粒落下速度
xz_MixRtAll = max( 0.0d0, xz_MixRt + xza_MixRtBasicZ(:,:,num) )
xz_FallRain = 0.0d0
xz_VelZRain = 0.0d0
!雨粒終端速度
xz_VelZRain = 12.2d0 * FactorJ * ( xz_MixRtAll ** 0.125d0 )
!落下による移流
! Dens の avr を取ってから割ると, ゼロ割が生じるので注意
xz_FallRain = xz_avr_xr( xr_dz_xz(xz_DensBasicZ * xz_VelZRain * xz_MixRtAll) ) / xz_DensBasicZ * RainSW(num)
! write(*,*) 'MixRt: ', minval( xz_MixRt ), maxval( xz_MixRt )
! write(*,*) 'Fall: ', minval( xz_FallRain ), maxval( xz_FallRain )
end function xz_FallRain
| Function : | |||
| xz_Rain2GasHeat(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8) | ||
| xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
| ||
| xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
| ||
| xza_DelMixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8), intent(in)
|
雨粒から蒸気への変換量を計算するためのルーチン
変換量および, 蒸気と雨粒の混合比は正の量なので, 計算の途中途中で 値が正になることを保証している. また, 元々存在する以上の雨粒が 蒸気に変換されないように, 元々の雨粒混合比を変換量の上限としている.
function xz_Rain2GasHeat(xz_PotTemp, xz_Exner, xza_DelMixRt)
!
! 雨粒から蒸気への変換量を計算するためのルーチン
!
! 変換量および, 蒸気と雨粒の混合比は正の量なので, 計算の途中途中で
! 値が正になることを保証している. また, 元々存在する以上の雨粒が
! 蒸気に変換されないように, 元々の雨粒混合比を変換量の上限としている.
!
!暗黙の型宣言禁止
implicit none
!変数定義
real(8), intent(in) :: xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax)
!温位の擾乱成分
real(8), intent(in) :: xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax)
!温度の擾乱成分
real(8), intent(in) :: xza_DelMixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!混合比の変化
real(8) :: xz_Rain2GasHeat(DimXMin:DimXMax, DimZMin:DimZMax)
real(8) :: xza_LatentHeat(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
real(8) :: xz_TempAll(DimXMin:DimXMax, DimZMin:DimZMax)
real(8) :: xz_ExnerAll(DimXMin:DimXMax, DimZMin:DimZMax)
integer :: s
!温度, 圧力, 混合比の全量を求める
!擾乱成分と平均成分の足し算
xz_ExnerAll = xz_Exner + xz_ExnerBasicZ
xz_TempAll = ( xz_PotTemp + xz_PotTempBasicZ ) * ( xz_Exner + xz_ExnerBasicZ )
xza_LatentHeat = 0.0d0
!雨から蒸気への相変化に伴う発熱
do s = 1, LoopNum
xza_LatentHeat(:,:,s) = xz_LatentHeat( SpcWetID(RainNum(s)), xz_TempAll ) * xza_DelMixRt(:,:,RainNum(s)) / (xz_ExnerAll * CpDry)
end do
xz_Rain2GasHeat = sum( xza_LatentHeat, 3 )
!値の保管
call StoreCond( xz_Rain2GasHeat )
end function xz_Rain2GasHeat
| Function : | |||
| xz_Rain2GasHeatNH4SH(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8) | ||
| xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
| ||
| xza_DelMixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8), intent(in)
|
雨粒から蒸気への変換量を計算するためのルーチン
変換量および, 蒸気と雨粒の混合比は正の量なので, 計算の途中途中で 値が正になることを保証している. また, 元々存在する以上の雨粒が 蒸気に変換されないように, 元々の雨粒混合比を変換量の上限としている.
function xz_Rain2GasHeatNH4SH(xz_Exner, xza_DelMixRt)
!
! 雨粒から蒸気への変換量を計算するためのルーチン
!
! 変換量および, 蒸気と雨粒の混合比は正の量なので, 計算の途中途中で
! 値が正になることを保証している. また, 元々存在する以上の雨粒が
! 蒸気に変換されないように, 元々の雨粒混合比を変換量の上限としている.
!
!暗黙の型宣言禁止
implicit none
!変数定義
real(8), intent(in) :: xza_DelMixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!混合比の変化量
real(8), intent(in) :: xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax)
!温度の擾乱成分
real(8) :: xz_Rain2GasHeatNH4SH(DimXMin:DimXMax, DimZMin:DimZMax)
!
real(8) :: xz_ExnerAll(DimXMin:DimXMax, DimZMin:DimZMax)
!圧力の全量を求める
!擾乱成分と平均成分の足し算
xz_ExnerAll = xz_Exner + xz_ExnerBasicZ
!雨から蒸気への相変化に伴う発熱
xz_Rain2GasHeatNH4SH = ReactHeatNH4SH * xza_DelMixRt(:,:,NH4SHRainNum) / (xz_ExnerAll * CpDry)
!値の保管
call StoreCond( xz_Rain2GasHeatNH4SH )
end function xz_Rain2GasHeatNH4SH
| Function : | |||
| xza_Cloud2Rain(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8)
| ||
| xza_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8), intent(in)
| ||
| DelTime : | real(8)
|
雲粒から雨粒への変換量を計算するためのルーチン 併合成長は Berry (1968) のパラメタリゼーションを利用し, 衝突合体成長は Kessler (1969) のパラメタリゼーションを利用する.
変換量および, 雲粒と雨粒の混合比は正の量なので, 計算の途中途中で 値が正になることを保証している. また, 元々存在する以上の雲粒が 雨粒に変換されないように, 元々の雲粒混合比を変換量の上限としている. 正の値を保証するために, 引数として時間刻みが必要となる. (AutoConv, Collect は時間刻み幅での積分値を計算)
このルーチンでは, 凝縮物質と反応生成物とを区別する必要が全くないので, ループを回す回数を LoopNum2 回としている.
function xza_Cloud2Rain( xza_MixRt, DelTime )
!
! 雲粒から雨粒への変換量を計算するためのルーチン
! 併合成長は Berry (1968) のパラメタリゼーションを利用し,
! 衝突合体成長は Kessler (1969) のパラメタリゼーションを利用する.
!
! 変換量および, 雲粒と雨粒の混合比は正の量なので, 計算の途中途中で
! 値が正になることを保証している. また, 元々存在する以上の雲粒が
! 雨粒に変換されないように, 元々の雲粒混合比を変換量の上限としている.
! 正の値を保証するために, 引数として時間刻みが必要となる.
! (AutoConv, Collect は時間刻み幅での積分値を計算)
!
! このルーチンでは, 凝縮物質と反応生成物とを区別する必要が全くないので,
! ループを回す回数を LoopNum2 回としている.
!
!暗黙の型宣言禁止
implicit none
!変数定義
real(8), intent(in) :: xza_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!混合比の擾乱成分
real(8) :: DelTime !時間刻み
real(8) :: xza_Cloud2Rain(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!雲から雨への変換量
real(8) :: xza_MixRtAll(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!混合比の擾乱成分 + 平均成分
real(8) :: xz_AutoConv(DimXMin:DimXMax, DimZMin:DimZMax)
!飽和混合比
real(8) :: xz_Collect(DimXMin:DimXMax, DimZMin:DimZMax)
!規格化された潜熱
real(8), parameter :: N0 = 5.0d7
real(8), parameter :: D0 = 3.66d-1
integer :: s
xza_Cloud2Rain = 0.0d0
!混合比は正の値を保証
!移流拡散で負になることもあり得るので.
xza_MixRtAll = max( 0.0d0, xza_MixRt + xza_MixRtBasicZ )
do s = 1, LoopNum2
xz_AutoConv = 0.0d0
xz_Collect = 0.0d0
!併合成長
! ! Kessler (1969) のパラメタリゼーション
! xz_AutoConv = &
! & DelTime / AutoConvTime &
! & * max( 0.0d0, ( xza_MixRtAll(:,:,CloudNum(s)) - MixRt_AutoConvCr) )
! Berry (1968) のパラメタリゼーション
xz_AutoConv = DelTime * xz_DensBasicZ * ( xza_MixRtAll(:,:,CloudNum(s)) ** 3.0d0 ) * 1.0d6 / ( 60.0d0 * ( 2.0d0 * xza_MixRtAll(:,:,CloudNum(s)) + 2.66d-8 * N0 / ( xz_DensBasicZ * D0 ) ) )
!衝突合体成長
! Kessler (1969) のパラメタリゼーション
xz_Collect = DelTime * 2.2d0 * FactorJ * xza_MixRtAll(:,:,CloudNum(s)) * ( xza_MixRtAll(:,:,RainNum(s)) * xz_DensBasicZ ) ** 0.875d0
!雲の変換量: 併合成長と合体衝突の和
! 元々の変化量を上限値として設定する. 負の値となる.
xza_Cloud2Rain(:,:,CloudNum(s)) = - min( xza_MixRtAll(:,:,CloudNum(s)), ( xz_AutoConv + xz_Collect ) )
!雨の変換量. 符号は雲の変換量とは反対.
xza_Cloud2Rain(:,:,RainNum(s)) = - xza_Cloud2Rain(:,:,CloudNum(s))
end do
! write(*,*) 'C2R: ', minval(xza_Cloud2Rain(:,:,1)), maxval(xza_Cloud2Rain(:,:,1))
! write(*,*) 'C2R: ', minval(xza_Cloud2Rain(:,:,2)), maxval(xza_Cloud2Rain(:,:,2))
! write(*,*) 'C2R: ', minval(xza_Cloud2Rain(:,:,3)), maxval(xza_Cloud2Rain(:,:,3))
end function xza_Cloud2Rain
| Function : | |||
| xza_Rain2Gas(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8) | ||
| xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
| ||
| xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
| ||
| xza_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8), intent(in)
| ||
| DelTime : | real(8)
|
雨粒から蒸気への変換量を計算するためのルーチン
変換量および, 蒸気と雨粒の混合比は正の量なので, 計算の途中途中で 値が正になることを保証している. また, 元々存在する以上の雨粒が 蒸気に変換されないように, 元々の雨粒混合比を変換量の上限としている.
木星の場合は, FactorJ で変換量を加速する.
function xza_Rain2Gas(xz_Exner, xz_PotTemp, xza_MixRt, DelTime)
!
! 雨粒から蒸気への変換量を計算するためのルーチン
!
! 変換量および, 蒸気と雨粒の混合比は正の量なので, 計算の途中途中で
! 値が正になることを保証している. また, 元々存在する以上の雨粒が
! 蒸気に変換されないように, 元々の雨粒混合比を変換量の上限としている.
!
! 木星の場合は, FactorJ で変換量を加速する.
!
!暗黙の型宣言禁止
implicit none
!変数定義
real(8), intent(in) :: xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax)
!温位の擾乱成分
real(8), intent(in) :: xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax)
!温度の擾乱成分
real(8), intent(in) :: xza_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!混合比の擾乱成分
real(8) :: DelTime !時間刻み
real(8) :: xza_Rain2Gas(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!
real(8) :: xz_TempAll(DimXMin:DimXMax, DimZMin:DimZMax)
!温度の擾乱成分 + 平均成分
real(8) :: xz_PressAll(DimXMin:DimXMax, DimZMin:DimZMax)
!全圧
real(8) :: xza_MixRtAll(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!混合比の擾乱成分 + 平均成分
real(8) :: xz_NonSaturate(DimXMin:DimXMax, DimZMin:DimZMax)
!未飽和度(飽和混合比と蒸気の混合比の差)
integer :: s
!温度, 圧力, 混合比の全量を求める
!擾乱成分と平均成分の足し算
xz_TempAll = ( xz_PotTemp + xz_PotTempBasicZ ) * ( xz_Exner + xz_ExnerBasicZ )
xz_PressAll = PressBasis * ((xz_Exner + xz_ExnerBasicZ) ** (CpDry / GasRDry))
xza_Rain2Gas = 0.0d0
!混合比は正の値であることを保証
!移流拡散計算で負になることがあり得るので.
xza_MixRtAll = max( 0.0d0, xza_MixRt + xza_MixRtBasicZ )
do s = 1, LoopNum
!飽和蒸気圧と混合比の差(飽和度)を計算.
! 雨から蒸気への変換量は飽和度に比例する.
xz_NonSaturate = max( 0.0d0, xz_SvapPress(SpcWetID(CloudNum(s)), xz_TempAll) * MolWtWet(CloudNum(s)) / ( MolWtDry * xz_PressAll) - xza_MixRtAll(:,:,GasNum(s)) )
!雨の変換量
! 元々の雨粒の混合比以上に蒸発が生じないように上限値を設定
xza_Rain2Gas(:,:,RainNum(s)) = - min( DelTime * 4.85d-2 * FactorJ * xz_NonSaturate * ( xza_MixRtAll(:,:,RainNum(s)) * xz_DensBasicZ )** 0.65d0, xza_MixRtAll(:,:,RainNum(s)) )
!蒸気の変換量
! 雨粒の変換量とは符号が逆となる
xza_Rain2Gas(:,:,GasNum(s)) = - xza_Rain2Gas(:,:,RainNum(s))
end do
end function xza_Rain2Gas
| Function : | |||
| xza_Rain2GasNH4SH(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8) | ||
| xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
| ||
| xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
| ||
| xza_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8), intent(in)
| ||
| DelTime : | real(8), intent(in)
|
雨粒から蒸気への変換量を計算するためのルーチン
変換量および, 蒸気と雨粒の混合比は正の量なので, 計算の途中途中で 値が正になることを保証している. また, 元々存在する以上の雨粒が 蒸気に変換されないように, 元々の雨粒混合比を変換量の上限としている.
function xza_Rain2GasNH4SH(xz_Exner, xz_PotTemp, xza_MixRt, DelTime)
!
! 雨粒から蒸気への変換量を計算するためのルーチン
!
! 変換量および, 蒸気と雨粒の混合比は正の量なので, 計算の途中途中で
! 値が正になることを保証している. また, 元々存在する以上の雨粒が
! 蒸気に変換されないように, 元々の雨粒混合比を変換量の上限としている.
!
!暗黙の型宣言禁止
implicit none
!変数定義
real(8), intent(in) :: xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax)
!温位の擾乱成分
real(8), intent(in) :: xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax)
!温度の擾乱成分
real(8), intent(in) :: xza_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!混合比の擾乱成分
real(8), intent(in) :: DelTime !時間刻み
real(8) :: xza_Rain2GasNH4SH(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!
real(8) :: xz_TempAll(DimXMin:DimXMax, DimZMin:DimZMax)
!温度の擾乱成分 + 平均成分
real(8) :: xz_PressAll(DimXMin:DimXMax, DimZMin:DimZMax)
!圧力の擾乱成分 + 平均成分
real(8) :: xza_MixRtAll(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
!混合比の擾乱成分 + 平均成分
real(8) :: xz_NonSaturate(DimXMin:DimXMax, DimZMin:DimZMax)
!未飽和度(飽和混合比と蒸気の混合比の差)
!温度, 圧力, 混合比の全量を求める
!擾乱成分と平均成分の足し算
xz_TempAll = ( xz_PotTemp + xz_PotTempBasicZ ) * ( xz_Exner + xz_ExnerBasicZ )
xz_PressAll = PressBasis * ((xz_Exner + xz_ExnerBasicZ) ** (CpDry / GasRDry))
xza_Rain2GasNH4SH = 0.0d0
!混合比は正の値であることを保証
!移流拡散計算で負になることがあり得るので.
xza_MixRtAll = max( 0.0d0, xza_MixRt + xza_MixRtBasicZ )
!飽和蒸気圧と混合比の差(飽和度)を計算.
! 雨から蒸気への変換量は飽和度に比例する.
! 未飽和度を求めたいので, マイナスをかけ算している
! (DelMixRtNH4SH は, NH4SH が増加する方向, すなわち飽和度を正としている)
xz_NonSaturate = max( 0.0d0, - xz_DelMixRtNH4SH( xz_TempAll, xz_PressAll, xza_MixRtAll(:,:,NH3Num), xza_MixRtAll(:,:,H2SNum), MolWtWet(NH3Num), MolWtWet(H2SNum) ) )
!雨の変換量
! 元々の雨粒の混合比以上に蒸発が生じないように上限値を設定
xza_Rain2GasNH4SH(:,:,NH4SHRainNum) = - min( DelTime * 4.85d-2 * FactorJ * xz_NonSaturate * ( xza_MixRtAll(:,:,NH4SHRainNum) * xz_DensBasicZ ) ** 0.65d0, xza_MixRtAll(:,:,NH4SHRainNum) )
!蒸気の変換量
! 雨粒の変換量とは符号が逆となる
xza_Rain2GasNH4SH(:,:,NH3Num) = - xza_Rain2GasNH4SH(:,:,NH4SHRainNum) * MolWtWet(NH3Num) / MolWtWet(NH4SHRainNum)
xza_Rain2GasNH4SH(:,:,H2SNum) = - xza_Rain2GasNH4SH(:,:,NH4SHRainNum) * MolWtWet(H2SNum) / MolWtWet(NH4SHRainNum)
end function xza_Rain2GasNH4SH