Class | MoistBuoyancy |
In: |
moist/moistbuoyancy.f90
|
Subroutine : |
subroutine MoistBuoy_Init( ) use ChemData, only: ChemData_OneSpcID !暗黙の型宣言禁止 implicit none !変数定義 integer :: s integer :: n1, n2 real(8), allocatable :: xza_MixRtBasicZPerMolWt(:,:,:) !基本場の混合比 / 分子量 !----------------------------------------------------------- ! 混合距離 !----------------------------------------------------------- MixLen = sqrt(DelX * DelZ) !----------------------------------------------------------- ! 雲粒と気体の ID の組を作る !----------------------------------------------------------- !化学種の中から雲粒を作るものを選び, その配列添え字と分子量を保管. LoopNum = 0 GasCount = 0 SelectCloud: do s = 1, SpcNum !'Cloud' という文字列が含まれるものの個数を数える n1 = index(SpcWetSymbol(s), '-Cloud' ) if (n1 /= 0) then LoopNum = LoopNum + 1 CloudNum(LoopNum)= s end if n2 = index(SpcWetSymbol(s), '-g' ) if (n2 /= 0) then GasCount = GasCount + 1 GasNum(GasCount)= s end if end do SelectCloud !----------------------------------------------------------- ! 確認 !----------------------------------------------------------- call MessageNotify( "M", "MoistBuoy_Init", "LoopNum = %d", i=(/LoopNum/) ) call MessageNotify( "M", "MoistBuoy_Init", "GasCount = %d", i=(/GasCount/)) call MessageNotify( "M", "MoistBuoy_Init", "GasNum = %d", i=(/GasNum/) ) call MessageNotify( "M", "MoistBuoy_Init", "CloudNum= % ad", i=(/CloudNum/)) if ( LoopNum /= GasCount ) then call MessageNotify( "E", "MoistBuoy_Init", "GasCount /= LoopNum" ) stop end if !----------------------------------------------------------- ! 作業配列の初期化. !----------------------------------------------------------- allocate( xza_MixRtBasicZPerMolWt(DimXMin:DimXMax, DimZMin:DimZMax, LoopNum), xr_MixRtBasicZPerMolWt(DimXMin:DimXMax, DimZMin:DimZMax), xz_MixRtBasicZPerMolWt(DimXMin:DimXMax, DimZMin:DimZMax), xr_MixRtBasicZ(DimXMin:DimXMax, DimZMin:DimZMax), xz_MixRtBasicZ(DimXMin:DimXMax, DimZMin:DimZMax) ) xza_MixRtBasicZPerMolWt = 0.0d0 xr_MixRtBasicZPerMolWt = 0.0d0 xz_MixRtBasicZPerMolWt = 0.0d0 xr_MixRtBasicZ = 0.0d0 xz_MixRtBasicZ = 0.0d0 call MessageNotify( "M", "MoistBuoy_Init", "MolWtWet = %f", d=(/pack(MolWtWet(:), .true.) /) ) do s = 1, LoopNum xza_MixRtBasicZPerMolWt(:,:,GasNum(s)) = xza_MixRtBasicZ(:,:,GasNum(s)) / MolWtWet(GasNum(s)) end do xz_MixRtBasicZPerMolWt = sum(xza_MixRtBasicZPerMolWt, 3) xr_MixRtBasicZPerMolWt = xr_avr_xz( sum(xza_MixRtBasicZPerMolWt, 3) ) xz_MixRtBasicZ = sum(xza_MixRtBasicZ, 3) xr_MixRtBasicZ = xr_avr_xz( sum(xza_MixRtBasicZ, 3) ) end subroutine MoistBuoy_Init
Function : | |||
xr_BuoyDrag(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8)
| ||
xza_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8), intent(in)
|
鉛直方向の運動方程式に現れる浮力項のうち, 引きづりのの効果だけを求める
function xr_BuoyDrag(xza_MixRt) ! ! 鉛直方向の運動方程式に現れる浮力項のうち, ! 引きづりのの効果だけを求める ! !暗黙の型宣言禁止 implicit none !変数定義 real(8), intent(in) :: xza_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) !凝縮成分の混合比 real(8) :: xr_BuoyDrag(DimXMin:DimXMax, DimZMin:DimZMax) !浮力項(分子量効果) !初期化 xr_BuoyDrag = 0.0d0 !浮力項の計算 xr_BuoyDrag = - Grav * xr_avr_xz( sum(xza_MixRt(:,:,LoopNum+1:SpcNum), 3) ) / ( 1.0d0 + xr_MixRtBasicZ ) end function xr_BuoyDrag
Function : | |||
xr_BuoyMolWt(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8)
| ||
xza_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8), intent(in)
|
鉛直方向の運動方程式に現れる浮力項のうち, 分子量の効果だけを求める
function xr_BuoyMolWt(xza_MixRt) ! ! 鉛直方向の運動方程式に現れる浮力項のうち, ! 分子量の効果だけを求める ! !暗黙の型宣言禁止 implicit none !変数定義 real(8), intent(in) :: xza_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) !凝縮成分の混合比 real(8) :: xr_BuoyMolWt(DimXMin:DimXMax, DimZMin:DimZMax) !浮力項(分子量効果) real(8) :: xza_MixRtPerMolWt(DimXMin:DimXMax, DimZMin:DimZMax, LoopNum) !混合比/分子量 integer :: s !初期化 xr_BuoyMolWt = 0.0d0 !作業配列の初期化. 気体のみ利用 do s = 1, LoopNum xza_MixRtPerMolWt(:,:,GasNum(s)) = xza_MixRt(:,:,GasNum(s)) / MolWtWet(GasNum(s)) end do !浮力項の計算 xr_BuoyMolWt = + Grav * xr_avr_xz( sum(xza_MixRtPerMolWt, 3) ) / ( 1.0d0 / MolWtDry + xr_MixRtBasicZPerMolWt ) - Grav * xr_avr_xz( sum(xza_MixRt(:,:,1:LoopNum), 3) ) / ( 1.0d0 + xr_MixRtBasicZ ) end function xr_BuoyMolWt
Function : | |||
xz_BuoyMoistKm(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_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) : | real(8), intent(in)
|
本来は NH4SH の凝結熱を考慮しないとならないのだが, ここでは とりあえず無視する.
function xz_BuoyMoistKm(xz_PotTemp, xz_Exner, xza_MixRt) ! !本来は NH4SH の凝結熱を考慮しないとならないのだが, ここでは !とりあえず無視する. ! !暗黙の型宣言禁止 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) :: xza_MixRtAll(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) !凝縮成分の混合比 real(8) :: xz_BuoyMoistKm(DimXMin:DimXMax, DimZMin:DimZMax) ! real(8) :: xza_LatentHeat(DimXMin:DimXMax, DimZMin:DimZMax, LoopNum) !潜熱 real(8) :: xz_TempAll(DimXMin:DimXMax, DimZMin:DimZMax) !温度 real(8) :: xz_EffHeat(DimXMin:DimXMax, DimZMin:DimZMax) ! real(8) :: xz_EffPotTemp(DimXMin:DimXMax, DimZMin:DimZMax) real(8) :: xz_EffMolWt(DimXMin:DimXMax, DimZMin:DimZMax) ! real(8) :: xza_MixRtPerMolWt(DimXMin:DimXMax, DimZMin:DimZMax, LoopNum) !混合比/分子量 integer :: s !温度, 圧力, 混合比の全量を求める !擾乱成分と平均成分の足し算 xz_TempAll = ( xz_PotTemp + xz_PotTempBasicZ ) * ( xz_Exner + xz_ExnerBasicZ ) xza_MixRtAll = xza_MixRtBasicZ + xza_MixRt xza_LatentHeat = 0.0d0 !作業配列の初期化. 気体のみ利用 do s = 1, LoopNum xza_MixRtPerMolWt(:,:,GasNum(s)) = xza_MixRt(:,:,GasNum(s)) / MolWtWet(GasNum(s)) end do !温度の効果 xz_EffPotTemp = xz_PotTemp / xz_PotTempBasicZ !分子量効果 + 引きづりの効果 xz_EffMolWt = + sum(xza_MixRtPerMolWt, 3) / ( 1.0d0 / MolWtDry + xz_MixRtBasicZPerMolWt ) - sum(xza_MixRt, 3) / ( 1.0d0 + xz_MixRtBasicZ ) !蒸気が蒸発する場合の潜熱を計算 ! 分子量の部分はいつでも効くが潜熱は飽和していないと効かないので, ! 雲の混合比がゼロの時には, 潜熱の寄与はゼロとなるように調節している do s = 1, LoopNum xza_LatentHeat(:,:,s) = xz_LatentHeat( SpcWetID(CloudNum(s)), xz_TempAll ) * xza_MixRtAll(:,:,GasNum(s)) * ( 5.0d-1 + sign( 5.0d-1, (xza_MixRtAll(:,:,CloudNum(s)) - 1.0d-4) ) ) end do xz_EffHeat = sum( xza_LatentHeat, 3 ) * xz_EffMolWtBasicZ / ( CpDry * xz_ExnerBasicZ ) !乱流拡散係数の時間発展式の浮力項を決める xz_BuoyMoistKm = - 3.0d0 * Grav * ( Cm ** 2.0d0 ) * ( MixLen ** 2.0d0 ) * xz_avr_xr( xr_dz_xz( xz_EffHeat + xz_PotTempBasicZ / xz_EffMolWtBasicZ * ( 1.0d0 + xz_EffPotTemp + xz_EffMolWt ) ) ) / ( 2.0d0 * xz_PotTempBasicZ / xz_EffMolWtBasicZ) end function xz_BuoyMoistKm