Class Thermo_Routine
In: thermo_routine.f90

熱力学関係変数を用いた処理

Methods

CTTemp   hydro_pres  

Included Modules

Algebra Phys_Const

Public Instance methods

Subroutine :
nx :integer, intent(in)
: x 方向の配列数
ny :integer, intent(in)
: y 方向の配列数
nz :integer, intent(in)
: z 方向の配列数
Q :real, intent(in), dimension(nx,ny,nz)
: 凝結物の混合比 [kg/kg]
hnum :integer, intent(inout), dimension(nx,ny)
: 雲頂に相当する鉛直格子点番号
thres :real, intent(in), optional
: 雲内と判断するための閾値 [kg/kg]
undeff :real, intent(in), optional
: 未定義

計算領域内の雲頂温度を計算するルーチン. 厳密には温度を計算するのではなく, 雲頂に相当する鉛直格子番号を 水平面に渡って計算するルーチン (さもないと引数に 3 次元温度の情報が 必要になるので, スタックの無駄遣い). 計算時には, 鉛直上端から下に下ろし, 閾値を越えたところを雲頂と判断する. 閾値はオプションで変更可能. デフォルトでは 5.e-5 [kg/kg] で設定. 雲として勘定するのは任意であり, サブルーチンを呼ぶ前に足し合わせておく. [重要] 雲がまったくなかった場合, nz の配列を 1 つ足したものが返されるので, 親ルーチンで, それを undef と認識させる処理をしておく必要あり.

[Source]

subroutine CTTemp( nx, ny, nz, Q, hnum, thres, undeff )
! 計算領域内の雲頂温度を計算するルーチン.
! 厳密には温度を計算するのではなく, 雲頂に相当する鉛直格子番号を
! 水平面に渡って計算するルーチン (さもないと引数に 3 次元温度の情報が
! 必要になるので, スタックの無駄遣い).
! 計算時には, 鉛直上端から下に下ろし, 閾値を越えたところを雲頂と判断する.
! 閾値はオプションで変更可能. デフォルトでは 5.e-5 [kg/kg] で設定.
! 雲として勘定するのは任意であり, サブルーチンを呼ぶ前に足し合わせておく.
! [重要]
! 雲がまったくなかった場合, nz の配列を 1 つ足したものが返されるので,
! 親ルーチンで, それを undef と認識させる処理をしておく必要あり.
  implicit none
  integer, intent(in) :: nx  ! x 方向の配列数
  integer, intent(in) :: ny  ! y 方向の配列数
  integer, intent(in) :: nz  ! z 方向の配列数
  real, intent(in), dimension(nx,ny,nz) :: Q  ! 凝結物の混合比 [kg/kg]
  integer, intent(inout), dimension(nx,ny) :: hnum  ! 雲頂に相当する鉛直格子点番号
  real, intent(in), optional :: thres  ! 雲内と判断するための閾値 [kg/kg]
  real, intent(in), optional :: undeff  ! 未定義
  real :: defun
  integer :: i, j, k
  real :: threshold

  if(present(thres))then
     threshold=thres
  else
     threshold=5.e-6
  end if

  if(present(undeff))then
     defun=undeff
  else
     defun=999.0
  end if

  do i=1,nx
     do j=1,ny
        do k=nz,1,-1
           if(Q(i,j,k)>threshold.and.Q(i,j,k)/=undeff)then
              hnum(i,j)=k
              exit
           else
              if(k==1)then
                 hnum(i,j)=nz+1
              end if
           end if
        end do
     end do
  end do

end subroutine
Subroutine :
z(:) :real, intent(in)
: r 方向の位置座標 [m]
rho(size(z)) :real, intent(in)
: 密度 [kg/m^3]
z_ref :real, intent(in)
: 積分定数となる位置座標 [m]
p_ref :real, intent(in)
: z_ref での気圧 (積分定数) [Pa]
pres(size(z)) :real, intent(inout)
: 静力学平衡場での気圧 [Pa]
 鉛直 1 次元のサウンディングデータから静力学平衡場を満たす気圧場を計算する.

[Source]

subroutine hydro_pres( z, rho, z_ref, p_ref, pres )
!  鉛直 1 次元のサウンディングデータから静力学平衡場を満たす気圧場を計算する.
  use Algebra
  use Phys_Const
  implicit none
  real, intent(in) :: z(:)  ! r 方向の位置座標 [m]
  real, intent(in) :: rho(size(z))  ! 密度 [kg/m^3]
  real, intent(in) :: z_ref  ! 積分定数となる位置座標 [m]
  real, intent(in) :: p_ref  ! z_ref での気圧 (積分定数) [Pa]
  real, intent(inout) :: pres(size(z))  ! 静力学平衡場での気圧 [Pa]
  integer :: i, nz

  nz=size(z)

  do i=1,nz
     if(z(i)<z_ref)then
        call rectangle_int( z, rho, z(i), z_ref, pres(i) )
        pres(i)=p_ref+pres(i)*g
     else
        if(z(i)>z_ref)then
           call rectangle_int( z, rho, z_ref, z(i), pres(i) )
           pres(i)=p_ref-pres(i)*g
        else
           pres(i)=p_ref
        end if
     end if
  end do

end subroutine