masscondense_threshold-masstest.f90

Path: physics/masscondense_threshold-masstest.f90
Last Update: Thu Mar 03 13:31:20 +0900 2011

    Copyright (C) GFD Dennou Club, 2004. All rights reserved.

begin

Subroutine MassCondense

  * Developer: KITAMORI Taichi, YAMASHITA Tatsuya
  * Version: $
  * Tag Name: $Name:  $
  * Change History:

Overview

凝結/蒸発量を計算する.

Error Handling

Known Bugs

Note

拡散によって雲粒が成長する場合における成長方程式を解く. 蒸発量が存在する雲の量よりも多い場合, 蒸発量を存在する雲の量にする. 雲密度を時間発展(凝結部分)を計算.

Future Plans

end

Required files

Methods

Included Modules

dc_trace gridset cloudset basicset ChemData StoreDensCloud

Public Instance methods

Subroutine :
DelTime :real(8), intent(in)
: 時間ステップ [s]
xz_LatHeatPerMassNl(DimXMin:DimXMax, DimZMin:DimZMax) :real(8), intent(in)
: 単位質量あたりの凝結熱 [J/kg]
xz_SatRatioNl(DimXMin:DimXMax, DimZMin:DimZMax) :real(8), intent(in)
: 飽和比 [1]
xz_PotTempNl(DimXMin:DimXMax, DimZMin:DimZMax) :real(8), intent(in)
: 温位 [K]
xz_ExnerNl(DimXMin:DimXMax, DimZMin:DimZMax) :real(8), intent(in)
: 無次元圧力 [1]
xz_DensCloudNl(DimXMin:DimXMax, DimZMin:DimZMax) :real(8), intent(inout)
: 雲の密度 [kg/m^3]
xz_MassCondNl(DimXMin:DimXMax, DimZMin:DimZMax) :real(8), intent(out)
: 凝結量 [kg/m^3 s]

(out) 凝結量

begin

Dependency

[Source]

subroutine MassCondense( DelTime, xz_LatHeatPerMassNl, xz_SatRatioNl, xz_PotTempNl, xz_ExnerNl, xz_DensCloudNl, xz_MassCondNl)           !(out) 凝結量
                                                                 !=begin
  !==Dependency
  use dc_trace, only: BeginSub, EndSub, DbgMessage
  use gridset,  only: DimXMin, DimXMax, DimZMin, DimZMax
  use cloudset, only: DensIce, NumAerosol, RadiAerosol, Kd, SatRatioCr        ! 臨界飽和比
!  use basicset, only: ss_CpBasicZ, &       ! 定圧比熱
!    &                 ss_PotTempBasicZ,  & ! 温位基本場
!    &                 ss_ExnerBasicZ,  &   ! 無次元圧力
!    &                 ss_DensBasicZ    &   ! 密度
!    &                 GasRDry              ! 気体定数
!  use BasicSet, only: xz_PotTempBasicZ,  & ! 温位基本場
  use basicset, only: xz_PotTempBasicZ, xz_ExnerBasicZ, xz_DensBasicZ, GasRDry, PressSfc, CpDry                 ! 定圧比熱
  use ChemData, only: ChemData_Cp, ChemData_SVapPress_AntoineA, ChemData_SVapPress_AntoineB    ! Antoine の式の係数 B
  use StoreDensCloud, only: StoreDensCloudCond

                                                                 !=end
  !==暗黙の型宣言を禁止
  implicit none
  
  !==Input
  real(8), intent(in)  :: DelTime       ! 時間ステップ [s]
  real(8), intent(in)   :: xz_LatHeatPerMassNl(DimXMin:DimXMax, DimZMin:DimZMax) 
                                        ! 単位質量あたりの凝結熱 [J/kg]
  real(8), intent(in)   :: xz_SatRatioNl(DimXMin:DimXMax, DimZMin:DimZMax)
                                        ! 飽和比 [1]
  real(8), intent(in)   :: xz_PotTempNl(DimXMin:DimXMax, DimZMin:DimZMax) 
                                        ! 温位 [K]
  real(8), intent(in)   :: xz_ExnerNl(DimXMin:DimXMax, DimZMin:DimZMax) 
                                        ! 無次元圧力 [1]
  real(8), intent(inout) :: xz_DensCloudNl(DimXMin:DimXMax, DimZMin:DimZMax)   
                                        ! 雲の密度   [kg/m^3]


  !==Output
  real(8), intent(out)  :: xz_MassCondNl(DimXMin:DimXMax, DimZMin:DimZMax)
                                        ! 凝結量 [kg/m^3 s]
   
  !==Work
  real(8)               :: xz_Rh(DimXMin:DimXMax, DimZMin:DimZMax)
                                        ! 熱拡散による抵抗係数
  real(8)               :: xz_RadiCloud(DimXMin:DimXMax, DimZMin:DimZMax)
                                        ! 雲粒の半径 [m]
  real, parameter       :: Pi = 3.1415926535897932385d0
                                        ! 円周率
!  real(8)               :: DensCloudThreshold = 1.0d-5
!  real(8)               :: DensCloudThreshold = 5.0d-6
  real(8)               :: DensCloudThreshold = 5.0d-5
  integer               :: i,k          ! ループ変数

  call BeginSub("MassCondense", fmt="%c", c1="Calculate latent heat per unit mass.")

  
  call DbgMessage( fmt="*** Debug Message [MassCondense] *** xz_Rh ")

!write(*,*) xz_PotTempNl(1,:)

  xz_Rh = xz_LatHeatPerMassNl**2.0d0 / ( Kd * GasRDry * (xz_ExnerBasicZ + xz_ExnerNl)** 2.0d0 * (xz_PotTempBasicZ + xz_PotTempNl)**2.0d0 )
      ! 山下, 20080410, 
      ! L**2 / (k * R * Exner**2 * theta**2) という値を計算している. 
      ! 参考文献 : 北守修論の付録 A に導出過程が記述されている. 
      ! 但し, 北守修論本文の Rh の式 (3.4) には誤植があるので注意されたい. 
      ! 教科書として 1 次参照元はどれなのか調べる必要がありそうだ. 

  call DbgMessage( fmt="*** Debug Message [MassCondense] *** xz_RadiCloud ")


  xz_RadiCloud = ( 3.0d0 * xz_DensCloudNl / (4.0d0 * Pi * DensIce * NumAerosol * xz_DensBasicZ) + RadiAerosol**3.0d0 )**(1.0d0 / 3.0d0)

do k= DimZMin, DimZMax
   do i = DimXMin, DimXMax

      if (xz_SatRatioNl(i,k) - SatRatioCr > epsilon(0.0d0)) then
! 飽和比が定めた臨界飽和比よりも大きい場合
! CO2 の飽和比の式が地球大気で用いられる飽和比の式と等価なものなのか確認する必要がある. 

        xz_MassCondNl(i,k) = max( - xz_DensCloudNl(i,k) / DelTime, 4.0d0 * Pi * xz_RadiCloud(i,k) * NumAerosol * xz_DensBasicZ(i,k) / xz_Rh(i,k) * (xz_SatRatioNl(i,k) - 1.0d0) )

     else if (xz_DensCloudNl(i,k) > DensCloudThreshold ) then
! 臨界飽和比を超えていないが, 雲密度が閾値以上である場合
! 飽和比が 1 以上のときに凝結, 1 以下のときに蒸発. 
        xz_MassCondNl(i,k) = max( - xz_DensCloudNl(i,k) / DelTime, 4.0d0 * Pi * xz_RadiCloud(i,k) * NumAerosol * xz_DensBasicZ(i,k) / xz_Rh(i,k) * (xz_SatRatioNl(i,k) - 1.0d0) )
     else
! 臨界飽和比を超えておらず, 雲密度が閾値未満の場合         
! 雲密度がゼロではなく, 飽和比 1.0 未満のときに限って蒸発. 
! それ以外の場合には何も起きない. 
        xz_MassCondNl(i,k) = max( - xz_DensCloudNl(i,k) / DelTime, ( 0.5d0 * ( 1.0d0 + sign(1.0d0, xz_DensCloudNl(i,k) - epsilon(0.0d0)) ) ) * ( 0.5d0 * ( 1.0d0 - sign(1.0d0, xz_SatRatioNl(i,k) - 1.0d0) ) ) * 4.0d0 * Pi * xz_RadiCloud(i,k) * NumAerosol * xz_DensBasicZ(i,k) / xz_Rh(i,k) * (xz_SatRatioNl(i,k) - 1.0d0) )

     end if

     call StoreDensCloudCond( xz_MassCondNl(i,k) )

     xz_DensCloudNl(i,k) = xz_DensCloudNl(i,k) + DelTime * xz_MassCondNl(i,k)

  end do
end do

  call EndSub("MassCondense")

end subroutine MassCondense