Subroutine : |
|
xz_DensCloudNs(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
|
DelTimeShort : | real(8), intent(in)
|
xz_TendDensCloudNl(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
: | 雲の密度の長い時間で評価する項 [kg/m^3/s]
|
|
xz_DensCloudWorkAs(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(out)
|
subroutine DensityCloud( xz_DensCloudNs, DelTimeShort, xz_TendDensCloudNl, xz_DensCloudWorkAs) !(out) 雲の密度(負の密度の調整前の値)
!=begin
!== Dependancy
use dc_trace, only: BeginSub, EndSub, DbgMessage
use gridset, only: DimXMin, DimXMax, DimZMin, DimZMax
use average, only: xr_avr_xz
use differentiate_center4, only: xz_dz_xr
use basicset, only: Grav, PressBasis, GasRDry, CpDry, CvDry, xz_ExnerBasicZ, xz_PotTempBasicZ, xz_DensBasicZ
use cloudset, only: DensIce, NumAerosol, RadiAerosol ! 凝結核半径
use StoreDensCloud, only: StoreDensCloudFall
!=end
!== 暗黙の型宣言禁止
implicit none
!=begin
!== Input
real(8), intent(in) :: DelTimeShort ! タイムステップ [s]
real(8), intent(in) :: xz_DensCloudNs(DimXMin:DimXMax, DimZMin:DimZMax)
! 雲の密度 [kg/m^3]
real(8), intent(in) :: xz_TendDensCloudNl(DimXMin:DimXMax, DimZMin:DimZMax)
! 雲の密度の長い時間で評価する項 [kg/m^3/s]
!== Output
real(8), intent(out) :: xz_DensCloudWorkAs(DimXMin:DimXMax, DimZMin:DimZMax)
! 雲の密度(負の密度調整前) [kg/m^3]
!=end
!== Work
real(8) :: xz_RadiCloudNs(DimXMin:DimXMax, DimZMin:DimZMax)
! 雲粒の半径 [m]
real(8) :: xz_FallDensCloud(DimXMin:DimXMax, DimZMin:DimZMax)
! 雲粒落下による雲密度の時間変化率 [kg/m^3 s]
real(8) :: xz_VelZTerm(DimXMin:DimXMax, DimZMin:DimZMax)
! 雲粒の終端速度
real(8) :: xz_VisCoffCO2(DimXMin:DimXMax, DimZMin:DimZMax)
! CO2 の粘性係数
real(8) :: xz_MeanFreePath(DimXMin:DimXMax, DimZMin:DimZMax)
! CO2 の平均自由行程
real(8) :: xz_KnudsenNum(DimXMin:DimXMax, DimZMin:DimZMax)
! 雲粒に対する Knudsen 数
real, parameter :: Pi = 3.1415926535897932385d0
! 円周率
real, parameter :: BoltzmannConst = 1.38065d-23
! ボルツマン定数
real, parameter :: CO2Diam = 3.30d-10
! CO2 分子の kinetic diameter
real, parameter :: SutherlandConst = 2.40d2
! CO2 分子のサザランド定数
real, parameter :: VisCoffRef = 1.47d-5
! 粘性係数の参照値
real, parameter :: TempRef = 2.93d2
! 粘性係数の参照温度
call BeginSub("DensityCloud", fmt="%c", c1="Time integration of density of cloud.")
! write(*,*) xz_DensCloudNl(1,:)
!=== 粘性係数の計算
xz_VisCoffCO2 = VisCoffRef * (TempRef + SutherlandConst) / ( xz_PotTempBasicZ * xz_ExnerBasicZ + SutherlandConst ) * ( xz_PotTempBasicZ * xz_ExnerBasicZ / TempRef )**1.5d0
!=== 平均自由行程の計算
xz_MeanFreePath = BoltzMannConst * xz_PotTempBasicZ / ( 2.0**0.5d0 * Pi * CO2Diam**2.0d0 * PressBasis * xz_ExnerBasicZ**(CvDry/GasRDry) )
!=== 雲粒半径の計算
xz_RadiCloudNs = ( 3.0d0 * xz_DensCloudNs / (4.0d0 * Pi * DensIce * NumAerosol * xz_DensBasicZ) + RadiAerosol**3.0d0 )**(1.0d0 / 3.0d0)
!=== Knudsen 数の計算
xz_KnudsenNum = xz_MeanFreePath / xz_RadiCloudNs
!=== 終端速度の計算
! Rossow(1978) の式
xz_VelZTerm = ( 1.0d0 + 4.0d0 * xz_KnudsenNum / 3.0d0 ) * 2.0d0 * xz_RadiCloudNs**2.0d0 * Grav * DensIce / (9.0d0 * xz_VisCoffCO2)
! Davies(1945) の式
! xz_VelZTerm = &
! & ( 1.0d0 + &
! & xz_KnudsenNum &
! & * (1.257d0 + 0.400d0 * exp(- 1.10d0 / xz_KnudsenNum) ) ) &
! & * 2.0d0 * xz_RadiCloudNs * xz_RadiCloudNs * Grav * DensIce &
! & / (9.0d0 * xz_VisCoffCO2)
!=== 落下項の計算. 4 次精度中心差分を用いて計算
xz_FallDensCloud = xz_dz_xr(xr_avr_xz(xz_VelZTerm * xz_DensCloudNs))
call StoreDensCloudFall( xz_FallDensCloud )
xz_DensCloudWorkAs = xz_DensCloudNs + DelTimeShort * ( xz_FallDensCloud + xz_TendDensCloudNl )
! !=== 境界条件
! call boundary(xz_BC, xz_DensCloud_out)
!
call EndSub("DensityCloud")
end subroutine DensityCloud