Class | Radiation |
In: |
physics/radiation.f90
|
モデルの放射過程を計算するためのパッケージ型モジュール 具体的には以下の項を計算するための関数を格納する.
* 一様冷却
Subroutine : | |
cfgfile : | character(*), intent(in) |
NAMELIST から放射強制の設定を取得
This procedure input/output NAMELIST#radiation .
subroutine Radiation_init(cfgfile) ! !NAMELIST から放射強制の設定を取得 ! !暗黙の型宣言禁止 implicit none !入力変数 character(*), intent(in) :: cfgfile real(8) :: RadHeightUp = 0.0d0 !放射強制を与える鉛直領域の上限 real(8) :: RadHeightDown = 0.0d0 !放射強制を与える鉛直領域の下限 real(8) :: RadHeightUp2 = 0.0d0 !放射強制を与える鉛直領域の上限 real(8) :: RadHeightDown2= 0.0d0 !放射強制を与える鉛直領域の下限 integer :: k !ループ変数 ! NAMELIST から情報を取得 NAMELIST /radiation/ RadHeatRate, RadHeightUp, RadHeightDown, RadCoolRate, RadHeight1, RadHeight2, RadHeightUp2, RadHeightDown2 open (10, FILE=cfgfile) read(10, NML=radiation) close(10) ! IntHeightUp の特定 do k = DimZMin, DimZMax if ( (RadHeightUp - DelZ <= s_Z(k)) .and. (s_Z(k) < RadHeightUp) ) then IntHeightUp = k else if ( (RadHeightUp2 - DelZ <= s_Z(k)) .and. (s_Z(k) < RadHeightUp2) ) then IntHeightUp2 = k end if end do allocate( xz_RadHeight(DimXMin:DimXMax, DimZMin:DimZMax) ) allocate( xz_RadHeight2(DimXMin:DimXMax, DimZMin:DimZMax) ) ! 放射強制が存在する領域の設定 ! RadHeightDown < s_Z < RadHeightUp の範囲では値が 1 となる ! 係数を用意する. ! xz_RadHeight = 1.0d0 xz_RadHeight2 = 1.0d0 do k = DimZMin, DimZMax if ( s_Z(k) <= RadHeightDown ) then xz_RadHeight(:,k) = 0.0d0 elseif( s_Z(k) >= RadHeightUp ) then xz_RadHeight(:,k) = 0.0d0 end if end do do k = DimZMin, DimZMax if ( s_Z(k) <= RadHeightDown2 ) then xz_RadHeight2(:,k) = 0.0d0 elseif( s_Z(k) >= RadHeightUp2 ) then xz_RadHeight2(:,k) = 0.0d0 end if end do if (cpurank == 0) then call MessageNotify( "M", "Radiation", "RadHeatRate = %f", d=(/RadHeatRate/)) call MessageNotify( "M", "Radiation", "RadHeightUp = %f", d=(/RadHeightUP/)) call MessageNotify( "M", "Radiation", "RadHeightDown= %f", d=(/RadHeightDown/)) end if end subroutine Radiation_init
Function : | |||
xz_NewtonCool(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8)
| ||
xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
|
温位の放射強制項. 地表面から RadHeight で指定された高度までの間で一様放射冷却を与える.
function xz_NewtonCool(xz_PotTemp) ! ! 温位の放射強制項. ! 地表面から RadHeight で指定された高度までの間で一様放射冷却を与える. ! !暗黙の型宣言を禁止 implicit none !変数定義 real(8), intent(in) :: xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax) !温位の擾乱場 real(8) :: xz_NewtonCool(DimXMin:DimXMax, DimZMin:DimZMax) !放射強制 real(8) :: xz_PotTempMean(DimXMin:DimXMax, DimZMin:DimZMax) !温位擾乱の東西平均 real(8) :: EFTime !放射緩和時間 integer :: k !放射緩和時間は 5 日 EFTime = 5.0d0 * DayTime do k = DimZMin, DimZMax xz_PotTempMean(:,k) = sum( xz_PotTemp(:,k) ) / real(DimXMax - DimXMin + 1) end do xz_NewtonCool = - xz_PotTempMean / EFTime call StorePotTempDamp( xz_NewtonCool ) end function xz_NewtonCool
Function : | |||
xz_RadHeatBalance(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8)
| ||
xz_ExnerBasicZ(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in) | ||
xz_PotTempBasicZ(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in) | ||
xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in) | ||
xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in) |
温位の放射強制項. 地表面から RadHeight で指定された高度までの間で空間的に一様な加熱, RadHeight から RadHeight2 までの間で空間的な一様な冷却を与える. 放射強制全体として加熱と冷却がバランスするように時々刻々一様加熱率 および一様冷却率を変化させる. 冷却の振幅を与え, 加熱の振幅を調節する.
function xz_RadHeatBalance(xz_ExnerBasicZ, xz_PotTempBasicZ, xz_Exner, xz_PotTemp) ! ! 温位の放射強制項. ! 地表面から RadHeight で指定された高度までの間で空間的に一様な加熱, ! RadHeight から RadHeight2 までの間で空間的な一様な冷却を与える. ! 放射強制全体として加熱と冷却がバランスするように時々刻々一様加熱率 ! および一様冷却率を変化させる. ! 冷却の振幅を与え, 加熱の振幅を調節する. !暗黙の型宣言を禁止 implicit none !変数定義 real(8), intent(in) :: xz_ExnerBasicZ(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(in) :: xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(in) :: xz_PotTempBasicZ(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(in) :: xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax) real(8) :: xz_RadHeatBalance(DimXMin:DimXMax, DimZMin:DimZMax) !放射強制 real(8) :: xz_DensSum(DimXMin:DimXMax, DimZMin:DimZMax) !密度(基本場成分と擾乱成分の和) real(8) :: HeatSum !ある時刻での単位時間当たりの加熱量 real(8) :: CoolSum !ある時刻での単位時間当たりの冷却量 integer :: i integer :: k ! 密度の算出 xz_DensSum = PressSfc * (xz_ExnerBasicZ + xz_Exner)**( CvDry /GasRDry ) / (GasRDry * (xz_PotTempBasicZ + xz_PotTemp ) ) HeatSum = 0.0d0 CoolSum = 0.0d0 !加熱するのは, 高度 0 - 1000 m => 格子点は 1 (RegZMin + 1) から IntHeightUp !冷却するのは, 高度 1000 - 15000 m => 格子点は IntHeightUp + 1 から IntHeightUp2 !*** 抜粋 *** ! IntHeightUp の特定 !do k = DimZMin, DimZMax ! if ( (RadHeightUp - DelZ <= s_Z(k)) .and. (s_Z(k) < RadHeightUp) ) then ! IntHeightUp = k ! else if ( (RadHeightUp2 - DelZ <= s_Z(k)) .and. (s_Z(k) < RadHeightUp2) ) then ! IntHeightUp2 = k ! end if !end do ! ! 下端 1000 ! |--*--|--*--|--*--|--*--|--*--|--*--|--*--| ! 1 5 ! ! 1000 Max ! |--*--|--*--|--*--|--*--|--*--|--*--|--*--| ! 6 I ! ある時刻での (加熱量/加熱率) を計算 do k = RegZmin + 1, IntHeightUp do i = RegXmin + 1, RegXMax HeatSum = HeatSum + DelZ * CpDry * xz_DensSum(i,k) end do end do ! ある時刻での (冷却量/冷却率) を計算 do k = IntHeightUp + 1, IntHeightUp2 do i = RegXmin + 1, RegXMax CoolSum = CoolSum + DelZ * CpDry * xz_DensSum(i,k) end do end do ! 加熱率の算出 ! RadCoolRate が負値であることに注意. RadHeatRate = - RadCoolRate * CoolSum / HeatSum xz_RadHeatBalance = xz_RadHeight * RadHeatRate / ( xz_ExnerBasicZ * DayTime ) + xz_RadHeight2 * RadCoolRate / ( xz_ExnerBasicZ * DayTime ) call StorePotTempRad( xz_RadHeatBalance ) end function xz_RadHeatBalance
Function : | |||
xz_RadHeatConst(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8)
| ||
xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
|
温位の放射強制項. 地表面から RadHeight で指定された高度までの間で一様放射冷却を与える.
function xz_RadHeatConst(xz_Exner) ! ! 温位の放射強制項. ! 地表面から RadHeight で指定された高度までの間で一様放射冷却を与える. ! !暗黙の型宣言を禁止 implicit none !変数定義 real(8), intent(in) :: xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) !エクスナー関数の擾乱場 real(8) :: xz_RadHeatConst(DimXMin:DimXMax, DimZMin:DimZMax) !放射強制 xz_RadHeatConst = xz_RadHeight * RadHeatRate / ( ( xz_ExnerBasicZ + xz_Exner ) * DayTime ) call StorePotTempRad( xz_RadHeatConst ) end function xz_RadHeatConst
Function : | |||
xz_RadHeatVary(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8)
| ||
xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) : | real(8), intent(in)
|
温位の放射強制項. 地表面から RadHeight1 で指定された高度までの間は2[K/day](RadHeatRate)の 一様な放射冷却を与える. RadHeight1 から RadHeight2 までは (RadHeight2 - s_Z(k)) / (RadHeight2 - RadHeight1) の割合で放射冷却の効果を減少させる (RadHeight2でちょうどゼロにする) Nakajima and Matsuno(1988),中島(1944)を参考にした
function xz_RadHeatVary(xz_Exner) ! ! 温位の放射強制項. ! 地表面から RadHeight1 で指定された高度までの間は2[K/day](RadHeatRate)の ! 一様な放射冷却を与える. ! RadHeight1 から RadHeight2 までは (RadHeight2 - s_Z(k)) / (RadHeight2 - RadHeight1) ! の割合で放射冷却の効果を減少させる (RadHeight2でちょうどゼロにする) ! Nakajima and Matsuno(1988),中島(1944)を参考にした ! !暗黙の型宣言を禁止 implicit none !変数定義 real(8), intent(in) :: xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) !エクスナー関数の擾乱場 real(8) :: xz_RadHeatVary(DimXMin:DimXMax, DimZMin:DimZMax) !放射強制 integer :: k do k = RegZMin+1, DimZMax if (s_Z(k).le.RadHeight1) then xz_RadHeatVary(:,k) = RadHeatRate / ( ( xz_ExnerBasicZ(:,k) + xz_Exner(:,k) ) * DayTime ) else if (s_Z(k).gt.RadHeight1.and.s_Z(k).le.RadHeight2) then xz_RadHeatVary(:,k) = RadHeatRate * (RadHeight2 - s_Z(k)) / ( (RadHeight2 - RadHeight1) * ( xz_ExnerBasicZ(:,k) + xz_Exner(:,k) ) * DayTime ) else if (s_Z(k).gt.RadHeight2) then xz_RadHeatVary(:,k) = 0.0d0 end if end do call StorePotTempRad( xz_RadHeatVary ) end function xz_RadHeatVary