| Class | adv_test |
| In: |
dynamics/adv_test.f90
|
Note that Japanese and English are described in parallel.
| Subroutine : |
This procedure input/output NAMELIST#adv_test_nml .
subroutine AdvTestInit
! ヒストリデータ出力
! History data output
!
use gtool_historyauto, only: HistoryAutoAddVariable
! 文字列操作
! Character handling
!
use dc_string, only: toChar
! NAMELIST ファイル入力に関するユーティリティ
! Utilities for NAMELIST file input
!
use namelist_util, only: namelist_filename, NmlutilMsg, NmlutilAryValid
! ファイル入出力補助
! File I/O support
!
use dc_iounit, only: FileOpen
! 温度の半整数σレベルの補間, 気圧と高度の算出
! Interpolate temperature on half sigma level,
! and calculate pressure and height
!
use auxiliary, only: AuxVarsInit
use constants , only: RPlanet, Grav, GasRDry
! $ R $ [J kg-1 K-1].
! 乾燥大気の気体定数.
! Gas constant of air
! 宣言文 ; Declaration statements
!
character(STRING) :: VelDist
real(DP) :: AlphaInDeg
!!$ real(DP) :: HCW0
!!$ real(DP) :: Temp0
integer:: unit_nml ! NAMELIST ファイルオープン用装置番号.
! Unit number for NAMELIST file open
integer:: iostat_nml ! NAMELIST 読み込み時の IOSTAT.
! IOSTAT of NAMELIST read
! NAMELIST 変数群
! NAMELIST group name
!
namelist /adv_test_nml/ VelDist, AlphaInDeg, Tau, Press0, Omega0, HCSigmaTop, HCNumCell, HCV0
!
! デフォルト値については初期化手続 "set_GATE_profile#SetGATEProfileInit"
! のソースコードを参照のこと.
!
! Refer to source codes in the initialization procedure
! "set_GATE_profile#SetGATEProfileInit" for the default values.
!
! デフォルト値の設定
! Default values settings
!
VelDist = "NCARASPSummerCol"
AlphaInDeg = 0.0_DP
Tau = 345600.0_DP
Press0 = 1000.0e2_DP
Omega0 = PI * 4.0e4_DP / Tau
HCSigmaTop = 0.1_DP
HCNumCell = 3
HCV0 = 10.0_DP
! NAMELIST の読み込み
! NAMELIST is input
!
if ( trim(namelist_filename) /= '' ) then
call FileOpen( unit_nml, namelist_filename, mode = 'r' ) ! (in)
rewind( unit_nml )
read( unit_nml, nml = adv_test_nml, iostat = iostat_nml ) ! (out)
close( unit_nml )
call NmlutilMsg( iostat_nml, module_name ) ! (in)
end if
select case ( VelDist )
case ( "NCARASPSummerCol" )
IDVelDist = IDVelDistNCARASPSummerCol
case ( "SymHadley" )
IDVelDist = IDVelDistSymHadley
case ( "AsymHadley" )
IDVelDist = IDVelDistAsymHadley
case default
call MessageNotify( 'E', module_name, 'VelDist of %c is not supported.', c1 = trim( VelDist ) )
end select
Alpha = AlphaInDeg * PI / 180.0_DP
! 補助的な変数を計算するサブルーチン・関数群
! Subroutines and functions for calculating auxiliary variables
!
call AuxVarsInit
! 印字 ; Print
!
call MessageNotify( 'M', module_name, '----- Initialization Messages -----' )
call MessageNotify( 'M', module_name, 'VelDist = %c', c1 = trim( VelDist ) )
call MessageNotify( 'M', module_name, 'Alpha = %f', d = (/ Alpha /) )
call MessageNotify( 'M', module_name, 'Tau = %f', d = (/ Tau /) )
call MessageNotify( 'M', module_name, 'Press0 = %f', d = (/ Press0 /) )
call MessageNotify( 'M', module_name, 'Omega0 = %f', d = (/ Omega0 /) )
call MessageNotify( 'M', module_name, 'HCSigmaTop = %f', d = (/ HCSigmaTop /) )
call MessageNotify( 'M', module_name, 'HCNumCell = %d', i = (/ HCNumCell /) )
!!$ call MessageNotify( 'M', module_name, 'InFileNameForcing = %c', c1 = trim(InFileNameForcing ) )
call MessageNotify( 'M', module_name, 'HCV0 = %f', d = (/ HCV0 /) )
call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) )
adv_test_inited = .true.
end subroutine AdvTestInit
| Subroutine : | |||
| xyz_UN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyz_VN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
|
——-セミラグのテスト用の流速分布を与える——— Gives a velocity for Test. Only used for debug.
subroutine AdvTestSetHorVels( xyz_UN, xyz_VN )
!-------セミラグのテスト用の流速分布を与える--------
!Gives a velocity for Test. Only used for debug.
use timeset, only: DelTime, TimeN ! ステップ $ t $ の時刻. Time of step $ t $.
use axesset , only : x_Lon, y_Lat, z_Sigma, r_Sigma
use constants , only: RPlanet, Grav
real(DP), intent(out) :: xyz_UN (0:imax-1, 1:jmax, 1:kmax)
! 東西風速
! Zonal Wind
real(DP), intent(out) :: xyz_VN (0:imax-1, 1:jmax, 1:kmax)
! 南北風速
! Meridional Wind
! 作業変数
! Work variables
!
real(DP) :: U0
real(DP) :: Time
real(DP) :: StreamFunc
integer:: i ! 東西方向に回る DO ループ用作業変数
! Work variables for DO loop in zonal direction
integer:: j ! 南北方向に回る DO ループ用作業変数
! Work variables for DO loop in meridional direction
integer:: k ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
! 初期化確認
! Initialization check
!
if ( .not. adv_test_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
select case ( IDVelDist )
case ( IDVelDistNCARASPSummerCol )
! 水平流速分布を与える
U0 = 2.0_DP * PI * RPlanet / ( 86400.0_DP * 12.0_DP )
do k = 1, kmax
do j = 1, jmax
do i = 0, imax-1
xyz_UN(i,j,k) = U0 * ( cos( y_Lat(j) ) * cos( Alpha ) + cos( x_Lon(i) ) * sin( y_Lat(j) ) * sin( Alpha ) )
xyz_VN(i,j,k) = - U0 * ( sin( x_Lon(i) ) * sin( Alpha ) )
end do
end do
end do
case ( IDVelDistSymHadley )
do k = 1, kmax
do j = 1, jmax
do i = 0, imax-1
if ( z_Sigma(k) < HCSigmaTop ) then
StreamFunc = 0.0_DP
xyz_UN(i,j,k) = 0.0_DP
xyz_VN(i,j,k) = 0.0_DP
else
!!$ StreamFunc = &
!!$ & StreamFunc0 &
!!$ & * sin( ( y_Lat(j) - 0.0_DP ) / ( HCLatMax - 0.0_DP ) * PI ) &
!!$ & * sin( (z_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI )
!!$
!!$ xyz_UN(i,j,k) = 0.0_DP
!!$ xyz_VN(i,j,k) = &
!!$ & StreamFunc0 * Grav / Press0 &
!!$ & * sin( ( y_Lat(j) - 0.0_DP ) / ( HCLatMax - 0.0_DP ) * PI ) &
!!$ & * PI/(r_Sigma(0)-HCSigmaTop) &
!!$ & * cos( (z_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI )
Time = TimeN
!!$ StreamFunc = &
!!$ & StreamFunc0 &
!!$ & * sin( HCNumCell * y_Lat(j) ) &
!!$ & * sin( (z_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI )
!!$
!!$ xyz_UN(i,j,k) = 0.0_DP
!!$ xyz_VN(i,j,k) = &
!!$ & StreamFunc0 * Grav / Press0 &
!!$ & * sin( HCNumCell * y_Lat(j) ) &
!!$ & * PI/(r_Sigma(0)-HCSigmaTop) &
!!$ & * cos( (z_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI ) &
!!$ & * cos( PI * Time / Tau )
!!$ StreamFunc = &
!!$ & PI * RPlanet**2 * Press0 / Grav * SigDot0 &
!!$ & * ( sin( dble( 2 * HCNumCell + 1 ) * y_Lat(j) ) &
!!$ & / dble( 2 * HCNumCell + 1 ) &
!!$ & + sin( dble( 2 * HCNumCell - 1 ) * y_Lat(j) ) &
!!$ & / dble( 2 * HCNumCell - 1 ) ) &
!!$ & * sin( (z_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI ) &
!!$ & * cos( PI * Time / Tau )
!!$
!!$ xyz_UN(i,j,k) = 0.0_DP
!!$ xyz_VN(i,j,k) = &
!!$ & - RPlanet / ( 2.0_DP * cos( y_Lat(j) ) ) * SigDot0 &
!!$ & * PI/(r_Sigma(0)-HCSigmaTop) &
!!$ & * ( sin( dble( 2 * HCNumCell + 1 ) * y_Lat(j) ) &
!!$ & / dble( 2 * HCNumCell + 1 ) &
!!$ & + sin( dble( 2 * HCNumCell - 1 ) * y_Lat(j) ) &
!!$ & / dble( 2 * HCNumCell - 1 ) ) &
!!$ & * cos( (z_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI ) &
!!$ & * cos( PI * Time / Tau )
StreamFunc = - 2.0_DP * PI * RPlanet * (r_Sigma(0)-HCSigmaTop) / PI * Press0 / Grav * HCV0 * sin( dble( 2 * HCNumCell ) * y_Lat(j) ) * cos( y_Lat(j) )**2 * sin( (z_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI ) * cos( PI * Time / Tau )
xyz_UN(i,j,k) = 0.0_DP
xyz_VN(i,j,k) = HCV0 * sin( dble( 2 * HCNumCell ) * y_Lat(j) ) * cos( y_Lat(j) ) * cos( (z_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI ) * cos( PI * Time / Tau )
end if
end do
end do
end do
case ( IDVelDistAsymHadley )
end select
end subroutine AdvTestSetHorVels
| Subroutine : | |||
| xy_Ps(0:imax-1, 1:jmax) : | real(DP), intent(in )
| ||
| xyz_Temp(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in )
| ||
| xyz_U(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyz_V(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(out)
|
——-セミラグのテスト用の流速分布を与える——— Gives a velocity for Test. Only used for debug.
subroutine AdvTestSetICs( xy_Ps, xyz_Temp, xyz_U, xyz_V, xyzf_QMix )
!-------セミラグのテスト用の流速分布を与える--------
!Gives a velocity for Test. Only used for debug.
use axesset , only : x_Lon, y_Lat, z_Sigma, r_Sigma
use constants , only: RPlanet
! $ a $ [m].
! 惑星半径.
! Radius of planet
use timeset, only: DelTime, TimeN ! ステップ $ t $ の時刻. Time of step $ t $.
! 温度の半整数σレベルの補間, 気圧と高度の算出
! Interpolate temperature on half sigma level,
! and calculate pressure and height
!
use auxiliary, only: AuxVars
real(DP), intent(in ) :: xy_Ps (0:imax-1, 1:jmax)
! 東西風速
! Zonal Wind
real(DP), intent(in ) :: xyz_Temp (0:imax-1, 1:jmax, 1:kmax)
! 温度
! Temperature
real(DP), intent(out) :: xyz_U (0:imax-1, 1:jmax, 1:kmax)
! 東西風速
! Zonal Wind
real(DP), intent(out) :: xyz_V (0:imax-1, 1:jmax, 1:kmax)
! 南北風速
! Meridional Wind
real(DP), intent(out) :: xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
!
! Mass mixing ratio
! 作業変数
! Work variables
!
real(DP) :: Lon0
real(DP) :: Lat0
real(DP) :: Height0
real(DP) :: Height1
real(DP) :: Height2
real(DP) :: RScale
real(DP) :: ZScale
real(DP) :: xyz_QVap (0:imax-1, 1:jmax, 1:kmax)
real(DP) :: xy_SurfHeight(0:imax-1, 1:jmax)
real(DP) :: xyz_Height (0:imax-1, 1:jmax, 1:kmax)
real(DP) :: xyz_HorDist(0:imax-1, 1:jmax, 1:kmax)
real(DP) :: xyz_VerDist(0:imax-1, 1:jmax, 1:kmax)
real(DP) :: D1
real(DP) :: D2
integer:: i ! 東西方向に回る DO ループ用作業変数
! Work variables for DO loop in zonal direction
integer:: j ! 南北方向に回る DO ループ用作業変数
! Work variables for DO loop in meridional direction
integer:: k ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
integer:: n
! 初期化確認
! Initialization check
!
if ( .not. adv_test_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
xyz_QVap = 0.0_DP
xy_SurfHeight = 0.0_DP
! 温度の半整数σレベルの補間, 気圧と高度の算出
! Interpolate temperature on half sigma level,
! and calculate pressure and height
!
call AuxVars( xy_Ps, xyz_Temp, xyz_QVap, xy_SurfHeight = xy_SurfHeight, xyz_Height = xyz_Height )
!
! Utility module for advection test
!
call AdvTestSetHorVels( xyz_U, xyz_V )
select case ( IDVelDist )
case ( IDVelDistNCARASPSummerCol )
Lon0 = 3.0_DP * PI / 2.0_DP
Lat0 = 0.0_DP
Height0 = 4.5e3_DP
RScale = RPlanet / 3.0_DP
ZScale = 1.0e3_DP
do k = 1, kmax
do j = 1, jmax
do i = 0, imax-1
xyz_HorDist(i,j,k) = RPlanet * acos( sin( Lat0 ) * sin( y_Lat(j) ) + cos( Lat0 ) * cos( y_Lat(j) ) * cos( x_Lon(i) - Lon0 ) )
end do
end do
end do
xyz_VerDist = xyz_Height - Height0
do n = 1, ncmax
do k = 1, kmax
do j = 1, jmax
do i = 0, imax-1
D2 = ( xyz_HorDist(i,j,k) / RScale )**2 + ( xyz_VerDist(i,j,k) / ZScale )**2
D1 = min( 1.0_DP, D2 )
if ( mod(n,2) == 1 ) then
! q5
xyzf_QMix(i,j,k,n) = 1.0_DP / 2.0_DP * ( 1.0_DP + cos( PI * D1 ) )
else if ( mod(n,2) == 0 ) then
! q6
if ( D2 <= 1.0_DP ) then
xyzf_QMix(i,j,k,n) = 1.0_DP
else
xyzf_QMix(i,j,k,n) = 0.0_DP
end if
if ( ( xyz_Height(i,j,k) > Height0 ) .and. ( ( Lat0 - 1.0_DP / 8.0_DP < y_Lat(j) ) .and. ( y_Lat(j) < Lat0 + 1.0_DP / 8.0_DP ) ) ) then
xyzf_QMix(i,j,k,n) = 0.0_DP
end if
else
xyzf_QMix(i,j,k,n) = 0.0_DP
end if
end do
end do
end do
end do
case ( IDVelDistSymHadley )
Height1 = 2.0e3_DP
Height2 = 5.0e3_DP
Height0 = ( Height1 + Height2 ) / 2.0_DP
xyz_VerDist = xyz_Height - Height0
do n = 1, ncmax
do k = 1, kmax
do j = 1, jmax
do i = 0, imax-1
if ( ( Height1 < xyz_Height(i,j,k) ) .and. ( xyz_Height(i,j,k) < Height2 ) ) then
xyzf_QMix(i,j,k,n) = ( 1.0_DP + cos( 2.0_DP * PI * xyz_VerDist(i,j,k) / ( Height2 - Height1 ) ) ) / 2.0_DP
else
xyzf_QMix(i,j,k,n) = 0.0_DP
end if
end do
end do
end do
end do
case ( IDVelDistAsymHadley )
end select
end subroutine AdvTestSetICs
| Subroutine : | |||
| xyr_SigDotN(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(out)
| ||
| xyr_StreamFunc(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(out), optional |
——-セミラグのテスト用の流速分布を与える——— Gives a velocity for Test. Only used for debug.
subroutine AdvTestSetVerVel( xyr_SigDotN, xyr_StreamFunc )
!-------セミラグのテスト用の流速分布を与える--------
!Gives a velocity for Test. Only used for debug.
use axesset, only : y_Lat, r_Sigma
use timeset, only: DelTime, TimeN ! ステップ $ t $ の時刻. Time of step $ t $.
use constants , only: RPlanet, Grav
real(DP), intent(out) :: xyr_SigDotN(0:imax-1, 1:jmax, 0:kmax)
! 鉛直流速（SigmaDot）
real(DP), intent(out), optional :: xyr_StreamFunc(0:imax-1, 1:jmax, 0:kmax)
! 作業変数
! Work variables
!
real(DP) :: Time
real(DP) :: Shape
real(DP) :: StreamFunc
integer:: i ! 東西方向に回る DO ループ用作業変数
! Work variables for DO loop in zonal direction
integer:: j ! 南北方向に回る DO ループ用作業変数
! Work variables for DO loop in meridional direction
integer:: k ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
! 初期化確認
! Initialization check
!
if ( .not. adv_test_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
select case ( IDVelDist )
case ( IDVelDistNCARASPSummerCol )
!鉛直流速分布を与える
Time = TimeN
do k = 0, kmax
Shape = min( 1.0_DP, 2.0_DP * sqrt( sin( (r_Sigma(k)-r_Sigma(kmax))/(r_Sigma(0)-r_Sigma(kmax)) * PI ) ) )
do j = 1, jmax
do i = 0, imax-1
xyr_SigDotN(i,j,k) = Omega0 / Press0 * cos( 2.0_DP*PI/Tau*Time )*sin( Shape*PI/2.0_DP )
end do
end do
end do
if ( present( xyr_StreamFunc ) ) then
xyr_StreamFunc = -999.0_DP
end if
case ( IDVelDistSymHadley )
do k = 0, kmax
do j = 1, jmax
do i = 0, imax-1
if ( r_Sigma(k) < HCSigmaTop ) then
StreamFunc = 0.0_DP
xyr_SigDotN(i,j,k) = 0.0_DP
else
!!$ StreamFunc = &
!!$ & StreamFunc0 &
!!$ & * sin( ( y_Lat(j) - 0.0_DP ) / ( HCLatMax - 0.0_DP ) * PI ) &
!!$ & * sin( (r_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI )
!!$
!!$ xyr_SigDotN(i,j,k) = &
!!$ & - StreamFunc0 * Grav / Press0 / ( RPlanet * cos( y_Lat(j) ) ) &
!!$ & * ( PI / ( HCLatMax - 0.0_DP ) &
!!$ & * cos( ( y_Lat(j) - 0.0_DP ) / ( HCLatMax - 0.0_DP ) * PI ) &
!!$ & * cos( y_Lat(j) ) &
!!$ & - sin( ( y_Lat(j) - 0.0_DP ) / ( HCLatMax - 0.0_DP ) * PI ) &
!!$ & * sin( y_Lat(j) ) ) &
!!$ & * sin( (r_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI )
Time = TimeN
!!$ StreamFunc = &
!!$ & StreamFunc0 &
!!$ & * sin( HCNumCell * y_Lat(j) ) &
!!$ & * sin( (r_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI )
!!$
!!$ xyr_SigDotN(i,j,k) = &
!!$ & - StreamFunc0 * Grav / Press0 / ( RPlanet * cos( y_Lat(j) ) ) &
!!$ & * ( HCNumCell &
!!$ & * cos( HCNumCell * y_Lat(j) ) * cos( y_Lat(j) ) &
!!$ & - sin( HCNumCell * y_Lat(j) ) * sin( y_Lat(j) ) ) &
!!$ & * sin( (r_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI ) &
!!$ & * cos( PI * Time / Tau )
!!$ StreamFunc = &
!!$ & PI * RPlanet**2 * Press0 / Grav * SigDot0 &
!!$ & * ( sin( dble( 2 * HCNumCell + 1 ) * y_Lat(j) ) &
!!$ & / dble( 2 * HCNumCell + 1 ) &
!!$ & + sin( dble( 2 * HCNumCell - 1 ) * y_Lat(j) ) &
!!$ & / dble( 2 * HCNumCell - 1 ) ) &
!!$ & * sin( (r_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI ) &
!!$ & * cos( PI * Time / Tau )
!!$
!!$ xyr_SigDotN(i,j,k) = &
!!$ & SigDot0 &
!!$ & * cos( dble( 2 * HCNumCell ) * y_Lat(j) ) &
!!$ & * sin( (r_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI ) &
!!$ & * cos( PI * Time / Tau )
StreamFunc = - 2.0_DP * PI * RPlanet * (r_Sigma(0)-HCSigmaTop) / PI * Press0 / Grav * HCV0 * sin( dble( 2 * HCNumCell ) * y_Lat(j) ) * cos( y_Lat(j) )**2 * sin( (r_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI ) * cos( PI * Time / Tau )
xyr_SigDotN(i,j,k) = - (r_Sigma(0)-HCSigmaTop) / ( PI * RPlanet ) * HCV0 * ( dble( 2 * HCNumCell ) * cos( dble( 2 * HCNumCell ) * y_Lat(j) ) * cos( y_Lat(j) ) - 2.0_DP * sin( dble( 2 * HCNumCell ) * y_Lat(j) ) * sin( y_Lat(j) ) ) * sin( (r_Sigma(k)-HCSigmaTop)/(r_Sigma(0)-HCSigmaTop) * PI ) * cos( PI * Time / Tau )
end if
if ( present( xyr_StreamFunc ) ) then
xyr_StreamFunc(i,j,k) = StreamFunc
end if
end do
end do
end do
xyr_SigDotN(:,:,0 ) = 0.0_DP
xyr_SigDotN(:,:,kmax) = 0.0_DP
case ( IDVelDistAsymHadley )
if ( present( xyr_StreamFunc ) ) then
xyr_StreamFunc = -999.0_DP
end if
end select
end subroutine AdvTestSetVerVel
| Variable : | |||
| adv_test_inited = .false. : | logical, save, public
|
| Constant : | |||
| version = ’$Name: $’ // ’$Id: dynamics_1d_utils.f90,v 1.1 2015/01/31 06:16:26 yot Exp $’ : | character(*), parameter
|