subroutine initialdata_takemi2007_basic( z_Temp, z_Press, za_MolFr )
!
!== 概要
! * deepconv の地球用のテスト計算としてTakemi(2007)の再現計算を
! するための湿度の基本場を作成する
! * 基本場の温度の式が温位で与えられているため, 温度に変換する必要がある
!
implicit none
real(DP), intent(out):: z_Press(DimZMin:DimZMax)!圧力
real(DP), intent(out):: z_Temp(DimZMin:DimZMax) !温度
real(DP), intent(out):: za_MolFr(DimZMin:DimZMax, 1:SpcNum) !モル比
!
! ECCM_Takemiでのみ必要な変数
!
real(DP) :: z_TakemiRH(DimZMin:DimZMax) !高度 1.5km 以上の相対湿度
real(DP) :: z_TakemiTheta(DimZMin:DimZMax) !基本場の温位
real(DP) :: Theta_0 ! 地表面の温位
real(DP) :: Theta_tr ! 対流圏界面の温位
real(DP) :: TakemiZtr ! 対流圏界面の高度
real(DP) :: HumSfc ! 地表面の相対湿度
real(DP) :: a_MolFrSfc ! 地表面のモル比
real(DP) :: z_LLhumidity(DimZMin:DimZMax) ! 1.5 km までの相対湿度
!!--------------------------------------------------
integer :: k, s, i, m, b, T, L
z_TakemiRH = 0.0d0
z_TakemiTheta = 0.0d0
! 湿度場の選択
select case (initialdata_takemi2007_hum)
case (ID_initialdata_takemi2007_hum)
call HUM_Takemi2007(z_TakemiRH,TakemiZtr,T,L)
case (ID_initialdata_takemi2007_TDRY)
call HUM_Takemi2007_TDRY(z_TakemiRH,TakemiZtr,T,L)
case (ID_initialdata_takemi2007_MDRY)
call HUM_Takemi2007_MDRY(z_TakemiRH,TakemiZtr,T,L)
end select
! 温度・湿度・気圧を求める
!地表面の温位は300K
Theta_0 = 300.0d0
! 中緯度場と熱帯場の選択
! 対流圏界面の温位で区別している
select case (initialdata_takemi2007_*)
case (ID_initialdata_takemi2007_midlat)
Theta_tr = 343.0d0 ! 中緯度
case (ID_initialdata_takemi2007_tropic)
Theta_tr = 358.0d0 ! 熱帯
end select
! 大気最下層の温位
z_TakemiTheta(RegZMin+1) = Theta_0 + (Theta_tr - Theta_0)*(z_Z(RegZMin+1) / TakemiZtr)**1.25
z_Temp = 1.0d-60
z_Press = 1.0d-60
z_LLhumidity = 0.0d0
! 大気最下層の気温・気圧・水蒸気混合比・分子量の計算
! 地表面のモル比
! 10g/kg〜18g/kgの混合比で与えられているためモル比に直す
! なお熱帯場はQ18のみ使用
a_MolFrSfc = Qb * MolWtDry / MolWtWet(1)
! a_MolFrSfc = 0.0180d0 * MolWtDry / MolWtWet(1) ! Q18
! a_MolFrSfc = 0.0160d0 * MolWtDry / MolWtWet(1) ! Q16
! a_MolFrSfc = 0.0140d0 * MolWtDry / MolWtWet(1) ! Q14
! a_MolFrSfc = 0.0120d0 * MolWtDry / MolWtWet(1) ! Q12
! a_MolFrSfc = 0.0100d0 * MolWtDry / MolWtWet(1) ! Q10
! 地表面の相対湿度
! モル比から逆算. 気圧を求めるのに必要
HumSfc = a_MolFrSfc * PressSfc / SvapPress(6,TempSfc)
! 大気最下層の気圧
z_Press(RegZMin+1) = PressSfc - (Grav * PressSfc * DelZ) /((GasRUniv/ (MolWtDry + a_MolFrSfc * (MolWtWet(1) - MolWtDry))) *TempSfc)
! 大気最下層の温度. 温位を元に計算
! 地表面気圧は 1000hPa と仮定する
z_Temp(RegZMin+1) = z_TakemiTheta(RegZMin+1) * (z_Press(RegZMin+1) / 100000.0d0)**(GasRDry / CpDry)
! 大気最下層のモル比
! 熱帯場はQ18のみ使用
za_MolFr(RegZMin+1,1) = Qb * MolWtDry / MolWtWet(1)
! za_MolFr(RegZMin+1,1) = 0.0180d0 * MolWtDry / MolWtWet(1) ! Q18
! za_MolFr(RegZMin+1,1) = 0.0160d0 * MolWtDry / MolWtWet(1) ! Q16
! za_MolFr(RegZMin+1,1) = 0.0140d0 * MolWtDry / MolWtWet(1) ! Q14
! za_MolFr(RegZMin+1,1) = 0.0120d0 * MolWtDry / MolWtWet(1) ! Q12
! za_MolFr(RegZMin+1,1) = 0.0100d0 * MolWtDry / MolWtWet(1) ! Q10
! 大気最下層の相対湿度
! モル比から逆算. 気圧を求めるのに必要
z_LLhumidity(RegZMin+1) = za_MolFr(RegZMin+1,1) * z_Press(RegZMin+1) / SvapPress(6,z_Temp(RegZMin+1))
! 大気最下層の種々の物理量をもとに温度・湿度気圧の基本場を作成する
! 高度 0 - 1.5km と 1.5km - で湿度の与え方が異なり,
! さらに, 対流圏界面以上では温度の式が異なるので, 以下の 3パターンに分けて
! 基本場を求める
! 高度 0 - 1.5km, 1.5km - 対流圏界面(Ztr), Ztr以上
do k = RegZMin+1, DimZMax
!
! 基本場(I)
! 大気最下層から高度 1.5 km 以下(k が L-1 以下)の基本場
!
if (k.le.L-1) then
! モル比
! 高度 1.5km まではモル比が任意の値で固定されている
! ただし, 熱帯場は 18 g/kg のみ
za_MolFr(k+1,1) = 0.0180d0 * MolWtDry / MolWtWet(1)
! za_MolFr(k+1,1) = 0.0160d0 * MolWtDry / MolWtWet(1)
! za_MolFr(k+1,1) = 0.0140d0 * MolWtDry / MolWtWet(1)
! za_MolFr(k+1,1) = 0.0120d0 * MolWtDry / MolWtWet(1)
! za_MolFr(k+1,1) = 0.0100d0 * MolWtDry / MolWtWet(1)
! 気圧
z_Press(k+1) = z_Press(k)-(Grav*z_Press(k)*DelZ) /((GasRUniv/ (MolWtDry + za_MolFr(k,1) * (MolWtWet(1) - MolWtDry))) *z_Temp(k))
! 温位から温度を計算
z_TakemiTheta(k+1) = Theta_0 + (Theta_tr - Theta_0)*(z_Z(k+1) / TakemiZtr)**1.25
z_Temp(k+1) = z_TakemiTheta(k+1) * (z_Press(k+1) / 100000.0d0)**(GasRDry / CpDry)
! モル比から相対湿度を計算
z_LLhumidity(k+1) = za_MolFr(k+1,1) * z_Press(k+1) / SvapPress(6,z_Temp(k+1))
!
!== Takemi(2007)の基本場の使用上の注意
! * 低層で混合比を固定し続けると高度 1.5km までに相対湿度が100%を超える高度がある
! * 100% を超えてしまった高度の湿度は, その直下の高度の湿度にすると
! Takemi(2007)と同じ湿度の鉛直プロファイルが描ける
! * なぞ
! * MQ16 は 湿度 99.・・・% という値が出て, 100% を超えなくてもほぼ100%
! な値が出てしまう
! * Takemi(2007)ではこのような値があるようには見えない
! しょうがないのでMQ18 とMQ14 の 100% を超える直前の値の間を取った値を入れておく
! * 0.94511507230d0という怪しい数字が上記の値
if (z_LLhumidity(k+1).gt.1) then
z_LLhumidity(k+1) = z_LLhumidity(k)
! if (z_LLhumidity(k+1).gt.0.99) then !MQ16 用
! z_LLhumidity(k+1) = 0.94511507230d0 !MQ16 用
end if
!
! 基本場(II)
! 高度 1.5km より上で, 対流圏界面以下(L < k < T)の基本場
!
else if (k.le.T.and.k.gt.L) then
! k = L での相対湿度
! この if 文の初めの k は L+1 から始まり z_TakemiRH(L) が区間外なので値を与える
z_TakemiRH(L) = z_LLhumidity(L)
! 気圧の基本場
z_Press(k) = z_Press(k-1)-(Grav*z_Press(k-1)*DelZ) /((GasRUniv/ (MolWtDry + za_MolFr(k+1,1) * (MolWtWet(1) - MolWtDry))) *z_Temp(k-1))
! 温位から温度を計算
z_TakemiTheta(k) = Theta_0 + (Theta_tr - Theta_0)*(z_Z(k) / TakemiZtr)**1.25
z_Temp(k) = z_TakemiTheta(k) * (z_Press(k) / 100000.0d0)**(GasRDry / CpDry)
! モル比と大気の平均分子量
za_MolFr(k,1) = SvapPress( 6,z_Temp(k))* z_TakemiRH(k)/z_Press(k)
!
! 基本場(III)
! 対流圏界面より上 (k > T) の場
!
else if (k.gt.T) then
! 気圧の基本場
z_Press(k) = z_Press(k-1)-(Grav*z_Press(k-1)*DelZ) /((GasRUniv/ (MolWtDry + za_MolFr(k,1) * (MolWtWet(1) - MolWtDry))) *z_Temp(k-1))
! 温位から温度を計算
z_TakemiTheta(k) = Theta_tr * exp(Grav * (z_Z(k) - TakemiZtr) / (CpDry * z_Temp(T)))
z_Temp(k) = z_TakemiTheta(k) * (z_Press(k) / 100000.0d0)**(GasRDry / CpDry)
! モル比と大気の平均分子量
za_MolFr(k,1) = SvapPress( 6,z_Temp(k))* z_TakemiRH(k)/z_Press(k)
end if
end do
end subroutine Initialdata_takemi2007_basic