Class Thermo_Advanced_Function
In: thermo_advanced_function.f90

基礎ルーチン, 関数集を use して複雑な熱力学関数を計算するモジュール

Methods

Louis   Rich   cm   ust  

Included Modules

Thermo_Function Thermo_Const Thermo_Routine Math_Const Phys_Const

Public Instance methods

Function :
Louis :real
z :real, intent(in)
: cm を求める高度 [m]
z0m :real, intent(in)
: モデルで計算される粗度高度 [m]
richard :real, intent(in)
: バルクリチャードソン数

Louis(1980) で提案されている大気の不安定度を考慮したバルク係数の補正係数を計算する関数

[Source]

real function Louis( z, z0m, richard )
! Louis(1980) で提案されている大気の不安定度を考慮したバルク係数の補正係数を計算する関数
  use Thermo_Const
  implicit none
  real, intent(in) :: z  ! cm を求める高度 [m]
  real, intent(in) :: z0m  ! モデルで計算される粗度高度 [m]
  real, intent(in) :: richard  ! バルクリチャードソン数
  real, parameter :: b=5.0, c=5.0
  real :: cm_tmp, zratio

  cm_tmp=(kalm/(log(z)-log(z0m)))**2
  zratio=z/z0m

  if(richard<0.0)then
     Louis=1.0-((2.0*b*richard)/(1.0+3.0*b*c*cm_tmp*sqrt(-richard*zratio)))
  else
     Louis=1.0/(1.0+2.0*b*richard*sqrt(1.0+c*richard))
  end if

  return
end function
Function :
Rich :real
za :real, intent(in)
: リチャードソン数を計算する高度 [m]
pta :real, intent(in)
: za での仮温位 [K]
ptg :real, intent(in)
: 地表面での温位 [K]
va :real, intent(in)
: 高度 za での水平風速の絶対値 [m/s]
qva :real, intent(in)
: za での混合比 [kg/kg]
qvs :real, intent(in)
: 地表面での飽和混合比 [kg/kg]

バルクリチャードソン数を計算する関数

[Source]

real function Rich( za, pta, ptg, va, qva, qvs )
! バルクリチャードソン数を計算する関数
  use Phys_Const
  use Thermo_Function
  implicit none
  real, intent(in) :: za  ! リチャードソン数を計算する高度 [m]
  real, intent(in) :: pta  ! za での仮温位 [K]
  real, intent(in) :: ptg  ! 地表面での温位 [K]
  real, intent(in) :: va  ! 高度 za での水平風速の絶対値 [m/s]
  real, intent(in) :: qva  ! za での混合比 [kg/kg]
  real, intent(in) :: qvs  ! 地表面での飽和混合比 [kg/kg]
  real :: ptvg, ptva, dpt

  ptvg=ptg*((1.0+eps_rdrv*qvs)/(1.0+qvs))
  ptva=pta*((1.0+eps_rdrv*qva)/(1.0+qva))
  dpt=ptva-ptvg
  Rich=(g*za*dpt)/(ptva*(va**2))

  return
end function
Function :
cm :real
z :real, intent(in)
: cm を求める高度 [m]
z0m :real, intent(in)
: モデルで計算される粗度高度 [m]
richard :real, intent(in), optional
: Louis (1980) のスキームで計算する場合のバルクリチャードソン数

運動量に関するバルク係数を計算する関数

[Source]

real function cm( z, z0m, richard )
! 運動量に関するバルク係数を計算する関数
  use Thermo_Const
  implicit none
  real, intent(in) :: z  ! cm を求める高度 [m]
  real, intent(in) :: z0m  ! モデルで計算される粗度高度 [m]
  real, intent(in), optional :: richard  ! Louis (1980) のスキームで計算する場合のバルクリチャードソン数

  if(present(richard))then
     cm=(kalm/(log(z)-log(z0m)))**2
     cm=cm*Louis( z, z0m, richard )
  else
     cm=(kalm/(log(z)-log(z0m)))**2
  end if

  return
end function
Function :
ust :real
cm :real, intent(in)
: 高度 za でのバルク係数
va :real, intent(in)
: 高度 za での水平風の絶対値 [m/s]

摩擦速度 u_* を計算する関数

[Source]

real function ust( cm, va )
! 摩擦速度 u_* を計算する関数
  implicit none
  real, intent(in) :: cm  ! 高度 za でのバルク係数
  real, intent(in) :: va  ! 高度 za での水平風の絶対値 [m/s]

  ust=va*sqrt(cm)

end function