Class Thermo_Advanced_Routine
In: thermo_advanced_routine.f90

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

Methods

Included Modules

Thermo_Function Thermo_Routine Thermo_Advanced_Function analy

Public Instance methods

Subroutine :
nx :integer, intent(in)
: 第 1 要素数
ny :integer, intent(in)
: 第 2 要素数
nz :integer, intent(in)
: 第 3 要素数
dx :real, intent(in)
: x 方向の格子間隔 [m]
dy :real, intent(in)
: y 方向の格子間隔 [m]
dz :real, intent(in)
: z 方向の格子間隔 [m]
u(nx,ny,nz) :real, intent(in)
: 速度場の x 成分 [m/s]
v(nx,ny,nz) :real, intent(in)
: 速度場の y 成分 [m/s]
w(nx,ny,nz) :real, intent(in)
: 速度場の z 成分 [m/s]
rho(nx,ny,nz) :real, intent(in)
: 密度場 [kg/m^3]
pt(nx,ny,nz) :real, intent(in)
: 温位場 [K]
pv(nx,ny,nz) :real, intent(inout)
: PV [Km^2/kgs]
undeff :real, intent(in), optional

エルテルのポテンシャル渦度を計算する

[Source]

subroutine Ertel_PV( nx, ny, nz, dx, dy, dz, u, v, w, rho, pt, pv, undeff )
! エルテルのポテンシャル渦度を計算する
  use Thermo_Function
  use Thermo_Routine
  use analy 
  implicit none
  integer, intent(in) :: nx  ! 第 1 要素数
  integer, intent(in) :: ny  ! 第 2 要素数
  integer, intent(in) :: nz  ! 第 3 要素数
  real, intent(in) :: dx  ! x 方向の格子間隔 [m]
  real, intent(in) :: dy  ! y 方向の格子間隔 [m]
  real, intent(in) :: dz  ! z 方向の格子間隔 [m]
  real, intent(in) :: u(nx,ny,nz)  ! 速度場の x 成分 [m/s]
  real, intent(in) :: v(nx,ny,nz)  ! 速度場の y 成分 [m/s]
  real, intent(in) :: w(nx,ny,nz)  ! 速度場の z 成分 [m/s]
  real, intent(in) :: rho(nx,ny,nz)  ! 密度場 [kg/m^3]
  real, intent(in) :: pt(nx,ny,nz)  ! 温位場 [K]
  real, intent(inout) :: pv(nx,ny,nz)  ! PV [Km^2/kgs]
  real, intent(in), optional :: undeff
  real :: tmp1(nx,ny,nz)
  real :: tmp2(nx,ny,nz)
  real :: tmp3(nx,ny,nz)
  real :: tmp4(nx,ny,nz)
  real :: tmp5(nx,ny,nz)
  real :: tmp6(nx,ny,nz)
  real :: tmp7(nx,ny,nz)
  integer :: i, j, k

  if(present(undeff))then
!  温位の空間勾配を計算.
     call grad_vec_3d( nx, ny, nz, dx, dy, dz, pt, tmp1, tmp2, tmp3, undeff )
!  3 次元 rotation を計算.
     call rotate( nx, ny, nz, dx, dy, dz, u, v, w, tmp4, tmp5, tmp6, undeff )
!  omega と grad pt の内積を計算
     call dot_prod( tmp4, tmp5, tmp6, tmp1, tmp2, tmp3, tmp7, undeff )
  else
!  温位の空間勾配を計算.
     call grad_vec_3d( nx, ny, nz, dx, dy, dz, pt, tmp1, tmp2, tmp3 )
!  3 次元 rotation を計算.
     call rotate( nx, ny, nz, dx, dy, dz, u, v, w, tmp4, tmp5, tmp6 )
!  omega と grad pt の内積を計算
     call dot_prod( tmp4, tmp5, tmp6, tmp1, tmp2, tmp3, tmp7 )
  end if

  if(present(undeff))then
!$omp parallel default(shared)
!$omp do private(i,j,k)

     do k=1,nz
        do j=1,ny
           do i=1,nx
              if(tmp7(i,j,k)==undeff.or.rho(i,j,k)==undeff)then
                 pv(i,j,k)=undeff
              else
                 pv(i,j,k)=tmp7(i,j,k)/rho(i,j,k)
              end if
           end do
        end do
     end do

!$omp end do
!$omp end parallel

  else

!$omp parallel default(shared)
!$omp do private(i,j,k)

     do k=1,nz
        do j=1,ny
           do i=1,nx
              pv(i,j,k)=tmp7(i,j,k)/rho(i,j,k)
           end do
        end do
     end do

!$omp end do
!$omp end parallel

  end if

end subroutine
Subroutine :
z :real, intent(in)
: cm を求める高度 [m]
z0m :real, intent(in), dimension(:,:)
: モデルで計算される粗度高度 [m]
richard :real, intent(in), dimension(size(z0m,1),size(z0m,2))
: バルクリチャードソン数
Lo :real, intent(inout), dimension(size(z0m,1),size(z0m,2))
: 補正係数

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

[Source]

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

  nx=size(z0m,1)
  ny=size(z0m,2)

!$omp parallel default(shared)
!$omp do private(i,j)
  do j=1,ny
     do i=1,nx
        Lo(i,j)=Louis( z, z0m(i,j), richard(i,j) )
     end do
  end do
!$omp end do
!$omp end parallel

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

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

[Source]

subroutine Rich_horizon( za, pta, ptg, va, qva, qvs, Ri )
! バルクリチャードソン数を計算するルーチン
  use Thermo_Advanced_Function
  implicit none
  real, intent(in) :: za  ! リチャードソン数を計算する高度 [m]
  real, intent(in), dimension(:,:) :: pta  ! za での仮温位 [K]
  real, intent(in), dimension(size(pta,1),size(pta,2)) :: ptg  ! 地表面での温位 [K]
  real, intent(in), dimension(size(pta,1),size(pta,2)) :: va  ! 高度 za での水平風速の絶対値 [m/s]
  real, intent(in), dimension(size(pta,1),size(pta,2)) :: qva  ! za での混合比 [kg/kg]
  real, intent(in), dimension(size(pta,1),size(pta,2)) :: qvs  ! 地表面での飽和混合比 [kg/kg]
  real, intent(inout), dimension(size(pta,1),size(pta,2)) :: Ri  ! 求めるリチャードソン数
  real, dimension(size(pta,1),size(pta,2)) :: ptvg, ptva, dpt
  integer :: i, j, nx, ny

  nx=size(pta,1)
  ny=size(pta,2)

!$omp parallel default(shared)
!$omp do private(i,j)
  do j=1,ny
     do i=1,nx
        Ri(i,j)=Rich( za, pta(i,j), ptg(i,j), va(i,j), qva(i,j), qvs(i,j) )
     end do
  end do
!$omp end do
!$omp end parallel

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

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

[Source]

subroutine cm_horizon( z, z0m, coem, richard )
! 運動量に関するバルク係数を計算するルーチン
  use Thermo_Advanced_Function
  implicit none
  real, intent(in) :: z  ! cm を求める高度 [m]
  real, intent(in), dimension(:,:) :: z0m  ! モデルで計算される粗度高度 [m]
  real, intent(inout), dimension(size(z0m,1),size(z0m,2)) :: coem  ! バルク係数
  real, intent(in), dimension(size(z0m,1),size(z0m,2)), optional :: richard  ! Louis (1980) のスキームで計算する場合のバルクリチャードソン数
  integer :: i, j, nx, ny

  nx=size(z0m,1)
  ny=size(z0m,2)

  if(present(richard))then
!$omp parallel default(shared)
!$omp do private(i,j)
     do j=1,ny
        do i=1,nx
           coem(i,j)=cm( z, z0m(i,j), richard(i,j) )
        end do
     end do
!$omp end do
!$omp end parallel
  else
!$omp parallel default(shared)
!$omp do private(i,j)
     do j=1,ny
        do i=1,nx
           coem(i,j)=cm( z, z0m(i,j), richard(i,j) )
        end do
     end do
!$omp end do
!$omp end parallel
  end if

end subroutine
Subroutine :
cm :real, intent(in), dimension(:,:)
: 高度 za でのバルク係数
va :real, intent(in), dimension(size(cm,1),size(cm,2))
: 高度 za での水平風の絶対値 [m/s]
velst :real, intent(inout), dimension(size(cm,1),size(cm,2))
: 摩擦速度 [m/s]

摩擦速度 u_* を計算するルーチン

[Source]

subroutine ust_horizon( cm, va, velst )
! 摩擦速度 u_* を計算するルーチン
  use Thermo_Advanced_Function
  implicit none
  real, intent(in), dimension(:,:) :: cm  ! 高度 za でのバルク係数
  real, intent(in), dimension(size(cm,1),size(cm,2)) :: va  ! 高度 za での水平風の絶対値 [m/s]
  real, intent(inout), dimension(size(cm,1),size(cm,2)) :: velst  ! 摩擦速度 [m/s]
  integer :: i, j, nx, ny

!$omp parallel default(shared)
!$omp do private(i,j)
     do j=1,ny
        do i=1,nx
           velst(i,j)=ust( cm(i,j), va(i,j) )
        end do
     end do
!$omp end do
!$omp end parallel

end subroutine