| Class | sltt |
| In: |
sltt/sltt.F90
|
Note that Japanese and English are described in parallel.
物質移流を非保存型のセミラグランジュ法で演算するモジュールです. 上流点探索には Williamson and Rasch (1989, MWR) を 補間には Enomoto (2008) を応用した方法を用いています。 すなわちスペクトルから求めた１階微分の値を利用した５次精度の変則エルミート補間です。 非負を保証するために arcsine 変換フィルタを用いています。 スペクトル変換・高精度補間に由来する人工的な短波を除去するために Sun et al. (1996) の 単調フィルタを応用したものを部分的に用いている。
This is a tracer transport module. Semi-Lagrangian method (Enomoto 2008 modified) Arcsine transformation filter is used to avoid negative values. Monotonicity filter (Sun et al 1996) is partly used.
| SLTTMain : | 移流計算 |
| SLTTInit : | 初期化 |
| SLTTTest : | 移流テスト用 |
| ——————— : | ———— |
| SLTTMain : | Main subroutine for SLTT |
| SLTTInit : | Initialization for SLTT |
| SLTTTest : | Generate velocity for SLTT Test |
NAMELIST#
| Subroutine : |
セミラグランジュ法の初期化処理 Initialization for Semi-Lagrangian method
This procedure input/output NAMELIST#sltt_nml .
subroutine SLTTInit
! セミラグランジュ法の初期化処理
! Initialization for Semi-Lagrangian method
use axesset, only : x_Lon, y_Lat
! 座標データ設定
! Axes data settings
!
use axesset, only: r_Sigma, z_Sigma ! $ \sigma $ レベル (整数).
! Full $ \sigma $ level
use sltt_const , only : SLTTConstInit
use sltt_extarr, only : SLTTExtArrInit
! NAMELIST ファイル入力に関するユーティリティ
! Utilities for NAMELIST file input
!
use namelist_util, only: namelist_filename, NmlutilMsg
! 種別型パラメタ
! Kind type parameter
!
use dc_types, only: STDOUT, STRING ! 文字列. Strings.
! ファイル入出力補助
! File I/O support
!
use dc_iounit, only: FileOpen
use sltt_const , only : iexmin, iexmax, jexmins, jexmaxs, jexminn, jexmaxn
!
! local variables
!
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:: unit_nml ! NAMELIST ファイルオープン用装置番号.
! Unit number for NAMELIST file open
integer:: iostat_nml ! NAMELIST 読み込み時の IOSTAT.
! IOSTAT of NAMELIST read
! NAMELIST 変数群
! NAMELIST group name
!
namelist /sltt_nml/ FlagSLTTArcsine, SLTTIntHor, SLTTIntVer
! 実行文 ; Executable statement
!
if ( sltt_inited ) return
if ( mod( jmax, 2 ) /= 0 ) then
stop 'jmax cannot be divided by 2.'
end if
call SLTTConstInit
! デフォルト値の設定
! Default values settings
!
FlagSLTTArcsine = .true.
SLTTIntHor = "HQ"
SLTTIntVer = "HQ"
! 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 = sltt_nml, iostat = iostat_nml ) ! (out)
close( unit_nml )
call NmlutilMsg( iostat_nml, module_name ) ! (in)
if ( iostat_nml == 0 ) write( STDOUT, nml = sltt_nml )
end if
allocate( x_LonS (0:imax-1) )
allocate( x_SinLonS(0:imax-1) )
allocate( x_CosLonS(0:imax-1) )
allocate( y_latS (1:jmax/2) )
allocate( y_SinLatS(1:jmax/2) )
allocate( y_CosLatS(1:jmax/2) )
do i = 0, imax-1
x_LonS (i) = x_Lon(i)
x_SinLonS(i) = sin( x_LonS(i) )
x_CosLonS(i) = cos( x_LonS(i) )
end do
do j = 1, jmax/2
y_LatS (j) = y_Lat(j)
y_SinLatS(j) = sin( y_LatS(j) )
y_CosLatS(j) = cos( y_LatS(j) )
end do
allocate( x_LonN (0:imax-1) )
allocate( x_SinLonN(0:imax-1) )
allocate( x_CosLonN(0:imax-1) )
allocate( y_latN (1:jmax/2) )
allocate( y_SinLatN(1:jmax/2) )
allocate( y_CosLatN(1:jmax/2) )
do i = 0, imax-1
x_LonN (i) = x_Lon(i)
x_SinLonN(i) = sin( x_LonN(i) )
x_CosLonN(i) = cos( x_LonN(i) )
end do
do j = 1, jmax/2
y_LatN (j) = y_Lat(j+jmax/2)
y_SinLatN(j) = sin( y_LatN(j) )
y_CosLatN(j) = cos( y_LatN(j) )
end do
allocate( x_ExtLonS( iexmin:iexmax ) )
allocate( x_ExtLonN( iexmin:iexmax ) )
allocate( y_ExtLatS( jexmins:jexmaxs ) )
allocate( y_ExtLatN( jexminn:jexmaxn ) )
call SLTTExtArrInit( x_LonS, y_LatS, x_LonN, y_LatN, x_ExtLonS, y_ExtLatS, x_ExtLonN, y_ExtLatN )
sltt_inited = .true.
end subroutine SLTTInit
| Subroutine : | |||
| xyz_UN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in )
| ||
| xyz_VN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in )
| ||
| xyr_SigDotN(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in )
| ||
| xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(inout)
|
セミラグランジュ法による物質移流計算を行う。 Calculates tracer transports by Semi-Lagrangian method
subroutine SLTTMain( xyz_UN, xyz_VN, xyr_SigDotN, xyzf_QMixA )
! セミラグランジュ法による物質移流計算を行う。
! Calculates tracer transports by Semi-Lagrangian method
real(DP), intent(in ) :: xyz_UN (0:imax-1, 1:jmax, 1:kmax)
! 東西風速
! Zonal Wind
real(DP), intent(in ) :: xyz_VN (0:imax-1, 1:jmax, 1:kmax)
! 南北風速
! Meridional Wind
real(DP), intent(in ) :: xyr_SigDotN(0:imax-1, 1:jmax, 0:kmax)
! 鉛直流速（SigmaDot）
real(DP), intent(inout) :: xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! 物質混合比
! Mix ratio of the tracers
! 作業変数
! Work variables
!
real(DP) :: f_QMixMax(1:ncmax)
! 各物質混合比の最大値
! Maximum of each mix ratio of the tracers
real(DP) :: f_QMixProcMax(1:ncmax)
! 各物質混合比のプロセス内最大値
! Maximum of each mix ratio of the tracers in each process
integer:: n ! 組成方向に回る DO ループ用作業変数
! Work variables for DO loop in dimension of constituents
real(DP) :: xyz_UTest (0:imax-1, 1:jmax, 1:kmax)
! 東西風速（テスト用）
! Zonal Wind (for test)
real(DP) :: xyz_VTest (0:imax-1, 1:jmax, 1:kmax)
! 南北風速（テスト用）
! Meridional Wind (for test)
real(DP) :: xyr_SigDotTest(0:imax-1, 1:jmax, 0:kmax)
! 鉛直流速（テスト用）;SigmaDot (for test)
if (FlagSLTTArcsine) then
! 非負を保証するための arcsine変換フィルタ
! Arcsine transformation for non-negative filter
do n = 1, ncmax
f_QMixProcMax(n) = maxval( xyzf_QMixA(:,:,:,n) )
end do
call MPIWrapperFindMaxVal( ncmax, f_QMixProcMax, f_QMixMax )
f_QMixMax = f_QMixMax * (1.05_DP) + 1.0D-14
do n = 1, ncmax
xyzf_QMixA(:,:,:,n) = 0.5_DP*(asin(2.0_DP*xyzf_QMixA(:,:,:,n)/f_QMixMax(n) - 1.0_DP))
end do
end if
xyzf_QMixA = SLTTHorAdv( xyzf_QMixA, xyz_UN, xyz_VN ) ! 水平セミラグ
! Horizontal
xyzf_QMixA = SLTTVerAdv( xyr_SigDotN, xyzf_QMixA ) ! 鉛直セミラグ
! Vertical
! 移流テスト
! call SLTTTest(xyz_UTest, xyz_VTest, xyr_SigDotTest)
! xyzf_QMixA = SLTTHorAdv( xyzf_QMixA, xyz_UTest, xyz_VTest ) ! 水平セミラグ
! xyzf_QMixA = SLTTVerAdv( xyr_SigDotTest, xyzf_QMixA ) ! 鉛直セミラグ
if (FlagSLTTArcsine) then
! 非負を保証するための arcsine変換フィルタ（逆変換）
! Arcsine transformation for non-negative filter
do n = 1, ncmax
xyzf_QMixA(:,:,:,n) = f_QMixMax(n)*(0.5_DP)*(sin(2.0_DP*xyzf_QMixA(:,:,:,n))+1.0_DP)
enddo
endif
end subroutine SLTTMain
| Variable : | |||
| FlagSLTTArcsine : | logical, save
|
| Function : | |||
| xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP)
| ||
| xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in )
| ||
| xyz_U(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in )
| ||
| xyz_V(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in )
|
セミラグランジュ法による水平移流の計算 Calculates tracer transports by Semi-Lagrangian method for horizontal direction
function SLTTHorAdv( xyzf_QMix, xyz_U, xyz_V ) result( xyzf_QMixA )
! セミラグランジュ法による水平移流の計算
! Calculates tracer transports by Semi-Lagrangian method for horizontal direction
use timeset , only : DelTime
! $\Delta t$
use axesset , only : x_Lon, y_Lat
! $\lambda, \varphai$ lon and lat
use sltt_const , only : dtjw, iexmin, iexmax, jexmins, jexmaxs, jexminn, jexmaxn
use sltt_extarr, only : SLTTExtArrExt, SLTTExtArrExt2
! 配列拡張ルーチン
! Expansion of arrays
use sltt_dp , only : SLTTDPHor
! 水平上流点探索
! Finding departure point in horizontal
use sltt_lagint, only : SLTTIrrHerIntK13
! 水平２次元の補間
! 2D Interpolation in horizontal
! SPMODEL ライブラリ, 球面上の問題を球面調和函数変換により解く(多層対応)
! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
!
#ifdef LIB_MPI
use wa_mpi_module, only: wa_xya => wa_xva, xya_wa => xva_wa, wa_DLon_wa, xya_GradLat_wa => xva_GradLat_wa
#elif AXISYMMETRY
use wa_zonal_module, only: wa_xya, xya_wa, wa_DLon_wa, xya_GradLat_wa
#elif SJPACK
use wa_module_sjpack, only: wa_xya, xya_wa, wa_DLon_wa, xya_GradLat_wa
#elif AXISYMMETRY_SJPACK
use wa_zonal_module_sjpack, only: wa_xya, xya_wa, wa_DLon_wa, xya_GradLat_wa
#else
use wa_module, only: wa_xya, xya_wa , wa_DLon_wa, xya_GradLat_wa
#endif
real(DP), intent(in ) :: xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! 現在時刻の物質混合比
! Present mix ratio of the tracers
real(DP), intent(in ) :: xyz_U (0:imax-1, 1:jmax, 1:kmax)
! 東西風速
! Zonal Wind
real(DP), intent(in ) :: xyz_V (0:imax-1, 1:jmax, 1:kmax)
! 南北風速
! Meridional Wind
real(DP) :: xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! 次ステップの物質混合比
! Next mix ratio of the tracers
!
! local variables
!
real(DP) :: xyzf_ExtQMixS(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax)
! 現在時刻の物質混合比の拡張配列（南半球）
! Extended array (SH) of present mix ratio of the tracers.
real(DP) :: xyzf_ExtQMixN(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax)
! 現在時刻の物質混合比の拡張配列（北半球）
! Extended array (NH) of present mix ratio of the tracers.
real(DP) :: xyz_ExtUS (iexmin:iexmax, jexmins:jexmaxs, 1:kmax)
! 東西風速の拡張配列（南半球）
! Extended array (SH) of Zonal Wind
real(DP) :: xyz_ExtUN (iexmin:iexmax, jexminn:jexmaxn, 1:kmax)
! 東西風速の拡張配列（北半球）
! Extended array (NH) of Zonal Wind
real(DP) :: xyz_ExtVS (iexmin:iexmax, jexmins:jexmaxs, 1:kmax)
! 南北風速の拡張配列（南半球）
! Extended array (SH) of Meridional Wind
real(DP) :: xyz_ExtVN (iexmin:iexmax, jexminn:jexmaxn, 1:kmax)
! 南北風速の拡張配列（北半球）
! Extended array (NH) of Meridional Wind
integer:: i, ii ! 東西方向に回る 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 ! 組成方向に回る DO ループ用作業変数
! Work variables for DO loop in dimension of constituents
real(DP) :: xyz_DPLonS(0:imax-1, 1:jmax/2, 1:kmax)
! 上流点経度（南半球）
! Lon of the departure point (SH)
real(DP) :: xyz_DPLonN(0:imax-1, 1:jmax/2, 1:kmax)
! 上流点経度（北半球）
! Lon of the departure point (NH)
real(DP) :: xyz_DPLatS(0:imax-1, 1:jmax/2, 1:kmax)
! 上流点緯度（南半球）
! Lat of the departure point (SH)
real(DP) :: xyz_DPLatN(0:imax-1, 1:jmax/2, 1:kmax)
! 上流点緯度（北半球）
! Lat of the departure point (NH)
real(DP) :: xyzf_QMixAS(0:imax-1, 1:jmax/2, 1:kmax, 1:ncmax)
! 次ステップの物質混合比（南半球）
! Next mix ratio of the tracers (SH)
real(DP) :: xyzf_QMixAN(0:imax-1, 1:jmax/2, 1:kmax, 1:ncmax)
! 次ステップの物質混合比（北半球）
! Next mix ratio of the tracers (NH)
!---fx, fy, fxy
real(DP) :: xyzf_QMix_dlon(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! 物質混合比の経度微分（グリッド）
! Zonal derivative of the mix ratio (on grid)
real(DP) :: xyzf_QMix_dlat(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! 物質混合比の緯度微分（グリッド）
! Meridional derivative of the mix ratio (on grid)
real(DP) :: xyzf_QMix_dlonlat(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! 物質混合比の緯度経度微分（グリッド）
! Zonal and meridional derivative of the mix ratio (on grid)
real(DP) :: xyzf_ExtQMixS_dlon(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax)
! 物質混合比の経度微分の拡張配列（南半球）
! Extended array (SH) of zonal derivative of the mix ratio
real(DP) :: xyzf_ExtQMixN_dlon(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax)
! 物質混合比の経度微分の拡張配列（北半球）
! Extended array (NH) of zonal derivative of the mix ratio
real(DP) :: xyzf_ExtQMixS_dlat(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax)
! 物質混合比の緯度微分の拡張配列（南半球）
! Extended array (SH) of meridional derivative of the mix ratio
real(DP) :: xyzf_ExtQMixN_dlat(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax)
! 物質混合比の緯度微分の拡張配列（北半球）
! Extended array (NH) of meridional derivative of the mix ratio
real(DP) :: xyzf_ExtQMixS_dlonlat(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax)
! 物質混合比の緯度経度微分の拡張配列（南半球）
! Extended array (SH) of zonal and meridional derivative of the mix ratio
real(DP) :: xyzf_ExtQMixN_dlonlat(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax)
! 物質混合比の緯度経度微分の拡張配列（北半球）
! Extended array (NH) of zonal and meridional derivative of the mix ratio
real(DP) :: wzf_QMix(1:lmax, 1:kmax, 1:ncmax)
! 物質混合比の経度微分（スペクトル）
! Zonal derivative of the mix ratio (on grid)
real(DP) :: wzf_QMix_dlon(1:lmax, 1:kmax, 1:ncmax)
! 物質混合比の経度微分（スペクトル）
! Zonal derivative of the mix ratio (on grid)
real(DP) :: pm ! 配列拡張する際、極ごえ後に符号が変わる場合は -1.0を与える。そうでない場合は1.0を与える。
! Sign change flag for array extension; -1.0 for sign change over the pole, 1.0 for no sign change
!---fxx, fyy, fxxyy
! real(DP) :: xyzf_QMix_dlon2(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_QMix_dlat2(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_QMix_dlon2lat2(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_ExtQMixS_dlon2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_ExtQMixN_dlon2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_ExtQMixS_dlat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_ExtQMixN_dlat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_ExtQMixS_dlon2lat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_ExtQMixN_dlon2lat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!----fxxy
! real(DP) :: xyzf_QMix_dlon2lat(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_ExtQMixS_dlon2lat(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_ExtQMixN_dlon2lat(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!----fxyy
! real(DP) :: xyzf_QMix_dlonlat2(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_ExtQMixS_dlonlat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
! real(DP) :: xyzf_ExtQMixN_dlonlat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!----
! real(DP) :: wzf_QMix_dlon2(1:lmax, 1:kmax, 1:ncmax)
! 実行文 ; Executable statement
!
! 初期化確認
! Initialization check
!
if ( .not. sltt_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
! QMixの微分計算（スペクトル変換利用）
! Derivatives of QMix
do n = 1, ncmax
wzf_QMix(:,:,n) = wa_xya(xyzf_QMix(:,:,:,n)) ! グリッド→スペクトル
! grid -> spectrum
xyzf_QMix_dlat(:,:,:,n) = xya_GradLat_wa(wzf_QMix(:,:,n)) ! スペクトル→グリッド緯度微分
! spectrum -> grid (dQ/dlat)
wzf_QMix_dlon(:,:,n) = wa_Dlon_wa(wzf_QMix(:,:,n)) ! スペクトル→スペクトル経度微分
! spectrum -> spectrum (dQ/dlon)
xyzf_QMix_dlon(:,:,:,n) = xya_wa(wzf_QMix_dlon(:,:,n)) ! スペクトル経度微分→グリッド経度微分
! spectrum (dQ/dlon) -> grid (dQ/dlon)
xyzf_QMix_dlonlat(:,:,:,n) = xya_GradLat_wa(wzf_QMix_dlon(:,:,n))! スペクトル経度微分→グリッド緯度経度微分
! spectrum (dQ/dlon) -> grid (d^2Q/dlon dlat)
!---fxx, fyy, fxxy, fxyy, fxxyy を計算
!xyzf_QMix_dlon2(:,:,:,n) = xya_wa(wa_Dlon_wa(wzf_QMix_dlon(:,:,n)))
!xyzf_QMix_dlat2(:,:,:,n) = xya_GradLat_wa(wa_xya(xyzf_QMix_dlat(:,:,:,n)))
!xyzf_QMix_dlon2lat(:,:,:,n) = xya_GradLat_wa(wa_xya(xyzf_QMix_dlon2(:,:,:,n)))
!xyzf_QMix_dlonlat2(:,:,:,n) = xya_GradLat_wa(wa_xya(xyzf_QMix_dlonlat(:,:,:,n)))
!xyzf_QMix_dlon2lat2(:,:,:,n) = xya_GradLat_wa(wa_xya(xyzf_QMix_dlon2lat(:,:,:,n)))
enddo
! 配列の分割と拡張
! Division and extension of arrays
!
! 配列の分割と拡張
! Division and extension of arrays
pm = -1.0_DP ! 配列拡張する際、極ごえ後に符号が変わる場合は -1.0を与える。そうでない場合は1.0を与える。
! -1.0 if the sign of value changes over the poles; if not 1.0.
!!$ call SLTTExtArrExt2( &
!!$ & xyzf_QMix_dlon, pm, & ! (in)
!!$ & xyzf_ExtQMixS_dlon, xyzf_ExtQMixN_dlon & ! (out)
!!$ & )
call SLTTExtArrExt2( x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, xyzf_QMix_dlon, pm, xyzf_ExtQMixS_dlon, xyzf_ExtQMixN_dlon, "Wave1" )
pm = -1.0_DP ! 配列拡張する際、極ごえ後に符号が変わる場合は -1.0を与える。そうでない場合は1.0を与える。
! -1.0 if the sign of value changes over the poles; if not 1.0.
!!$ call SLTTExtArrExt2( &
!!$ & xyzf_QMix_dlat, pm, & ! (in)
!!$ & xyzf_ExtQMixS_dlat, xyzf_ExtQMixN_dlat & ! (out)
!!$ & )
call SLTTExtArrExt2( x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, xyzf_QMix_dlat, pm, xyzf_ExtQMixS_dlat, xyzf_ExtQMixN_dlat, "Wave1" )
pm = +1.0_DP ! 配列拡張する際、極ごえ後に符号が変わる場合は -1.0を与える。そうでない場合は1.0を与える。
! -1.0 if the sign of value changes over the poles; if not 1.0.
!!$ call SLTTExtArrExt2( &
!!$ & xyzf_QMix_dlonlat, pm, & ! (in)
!!$ & xyzf_ExtQMixS_dlonlat, xyzf_ExtQMixN_dlonlat & ! (out)
!!$ & )
call SLTTExtArrExt2( x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, xyzf_QMix_dlonlat, pm, xyzf_ExtQMixS_dlonlat, xyzf_ExtQMixN_dlonlat, "Wave1" )
!-----fxx, fyy, fxxy, fxyy, fxxyy の配列拡張
! pm = +1.0_DP ! 配列拡張する際、極ごえ後に符号が変わる場合は -1.0を与える。そうでない場合は1.0を与える。
! -1.0 if the sign of value changes over the poles; if not 1.0.
! call SLTTExtArrExt2( &
! & xyzf_QMix_dlon2, pm, & ! (in)
! & xyzf_ExtQMixS_dlon2, xyzf_ExtQMixN_dlon2 & ! (out)
! & )
! pm = +1.0_DP ! 配列拡張する際、極ごえ後に符号が変わる場合は -1.0を与える。そうでない場合は1.0を与える。
! -1.0 if the sign of value changes over the poles; if not 1.0.
! call SLTTExtArrExt2( &
! & xyzf_QMix_dlat2, pm, & ! (in)
! & xyzf_ExtQMixS_dlat2, xyzf_ExtQMixN_dlat2 & ! (out)
! & )
! pm = -1.0_DP ! 配列拡張する際、極ごえ後に符号が変わる場合は -1.0を与える。そうでない場合は1.0を与える。
! -1.0 if the sign of value changes over the poles; if not 1.0.
! call SLTTExtArrExt2( &
! & xyzf_QMix_dlon2lat, pm, & ! (in)
! & xyzf_ExtQMixS_dlon2lat, xyzf_ExtQMixN_dlon2lat & ! (out)
! & )
! pm = -1.0_DP ! 配列拡張する際、極ごえ後に符号が変わる場合は -1.0を与える。そうでない場合は1.0を与える。
! -1.0 if the sign of value changes over the poles; if not 1.0.
! call SLTTExtArrExt2( &
! & xyzf_QMix_dlonlat2, pm, & ! (in)
! & xyzf_ExtQMixS_dlonlat2, xyzf_ExtQMixN_dlonlat2 & ! (out)
! & )
! pm = +1.0_DP ! 配列拡張する際、極ごえ後に符号が変わる場合は -1.0を与える。そうでない場合は1.0を与える。
! -1.0 if the sign of value changes over the poles; if not 1.0.
! call SLTTExtArrExt2( &
! & xyzf_QMix_dlon2lat2, pm, & ! (in)
! & xyzf_ExtQMixS_dlon2lat2, xyzf_ExtQMixN_dlon2lat2 & ! (out)
! & )
!!$ call SLTTExtArrExt( &
!!$ & x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, & ! (in)
!!$ & xyzf_QMix, xyz_U, xyz_V, & ! (in)
!!$ & xyzf_ExtQMixS, xyzf_ExtQMixN, & ! (out)
!!$ & xyz_ExtUS, xyz_ExtUN, & ! (out)
!!$ & xyz_ExtVS, xyz_ExtVN & ! (out)
!!$ & )
call SLTTExtArrExt( y_ExtLatS, y_ExtLatN, x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, xyzf_QMix, xyz_U, xyz_V, xyzf_ExtQMixS_dlat, xyzf_ExtQMixN_dlat, xyzf_ExtQMixS, xyzf_ExtQMixN, xyz_ExtUS, xyz_ExtUN, xyz_ExtVS, xyz_ExtVN )
! 上流点の計算
! estimation of departure point
! 南半球
! south array
call SLTTDPHor( DelTime, x_LonS, y_LatS, y_SinLatS, y_CosLatS, iexmin, iexmax, jexmins, jexmaxs, x_ExtLonS, y_ExtLatS, xyz_ExtUS, xyz_ExtVS, xyz_DPLonS, xyz_DPLatS )
! 北半球
! north array
call SLTTDPHor( DelTime, x_LonN, y_LatN, y_SinLatN, y_CosLatN, iexmin, iexmax, jexminn, jexmaxn, x_ExtLonN, y_ExtLatN, xyz_ExtUN, xyz_ExtVN, xyz_DPLonN, xyz_DPLatN )
! 補間
! Interpolation
! do n = 1, ncmax
call SLTTIrrHerIntK13( iexmin, iexmax, jexmins, jexmaxs, x_ExtLonS, y_ExtLatS, xyz_DPLonS, xyz_DPLatS, xyzf_ExtQMixS(:,:,:,:), xyzf_ExtQMixS_dlon(:,:,:,:), xyzf_ExtQMixS_dlat(:,:,:,:), xyzf_ExtQMixS_dlonlat(:,:,:,:), SLTTIntHor, xyzf_QMixAS(:,:,:,:) )
call SLTTIrrHerIntK13( iexmin, iexmax, jexminn, jexmaxn, x_ExtLonN, y_ExtLatN, xyz_DPLonN, xyz_DPLatN, xyzf_ExtQMixN(:,:,:,:), xyzf_ExtQMixN_dlon(:,:,:,:), xyzf_ExtQMixN_dlat(:,:,:,:), xyzf_ExtQMixN_dlonlat(:,:,:,:), SLTTIntHor, xyzf_QMixAN(:,:,:,:) )
! enddo
! 南北半球の配列の結合
! joint of each array
xyzf_QMixA(:,1:jmax/2,:,:) = xyzf_QMixAS(:,1:jmax/2,:,:)
xyzf_QMixA(:,jmax/2+1:jmax,:,:) = xyzf_QMixAN(:,1:jmax/2,:,:)
end function SLTTHorAdv
| Variable : | |||
| SLTTIntHor : | character(TOKEN), save
|
| Variable : | |||
| SLTTIntVer : | character(TOKEN), save
|
| Subroutine : | |||
| xyz_UTest(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyz_VTest(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyr_SigDotTest(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(out)
|
——-セミラグのテスト用の流速分布を与える——— Gives a velocity for Test. Only used for debug.
subroutine SLTTTest( xyz_UTest, xyz_VTest, xyr_SigDotTest )
!-------セミラグのテスト用の流速分布を与える--------
!Gives a velocity for Test. Only used for debug.
use constants0, only : PI
use axesset , only : x_Lon, y_Lat, r_Sigma
use constants, only: RPlanet
! $ a $ [m].
! 惑星半径.
! Radius of planet
use timeset, only: DelTime, TimeN ! ステップ $ t $ の時刻. Time of step $ t $.
real(DP), intent(out) :: xyz_UTest (0:imax-1, 1:jmax, 1:kmax)
! 東西風速
! Zonal Wind
real(DP), intent(out) :: xyz_VTest (0:imax-1, 1:jmax, 1:kmax)
! 南北風速
! Meridional Wind
real(DP), intent(out) :: xyr_SigDotTest(0:imax-1, 1:jmax, 0:kmax)
! 鉛直流速（SigmaDot）
! 作業変数
! Work variables
!
real(DP) :: u0, t, shape
real(DP), parameter :: lat0 = 0.5_8*PI, lon0 = 0.0_8*PI
real(DP), parameter :: tau = 345600.0_DP, p0 = 100000.0_DP
real(DP), parameter :: omega0 = PI*30000.0_DP/tau
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
! 水平流速分布を与える
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_UTest(i,j,k) = u0*(cos(y_Lat(j)) * cos(lat0) + cos(x_Lon(i)) * sin(y_Lat(j)) * sin(lat0))
xyz_VTest(i,j,k) = -u0*(sin(x_Lon(i))*sin(lat0))
end do
end do
end do
!鉛直流速分布を与える
t = TimeN
do k = 0, kmax
shape = min(1.0_DP, 0.5_DP*(sin((r_Sigma(k)-r_Sigma(kmax))/(1.0_DP - r_Sigma(kmax))*PI)) )
do j = 1, jmax
do i = 0, imax-1
xyr_SigDotTest(i,j,k) = -omega0/p0*cos(2.0_DP*PI/tau*t)*sin(shape*PI*0.5_DP)
enddo
enddo
enddo
end subroutine SLTTTest
| Function : | |||
| xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP)
| ||
| xyr_SigmaDot(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in)
| ||
| xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
|
セミラグランジュ法による鉛直移流の計算 Calculates tracer transports by Semi-Lagrangian method for vertical direction
function SLTTVerAdv( xyr_SigmaDot, xyzf_QMix ) result( xyzf_QMixA )
! セミラグランジュ法による鉛直移流の計算
! Calculates tracer transports by Semi-Lagrangian method for vertical direction
use axesset, only : z_Sigma ! 鉛直座標; Sigma coordinate
use timeset, only : DelTime ! $\Delta t$
use sltt_dp, only : SLTTDPVer ! 鉛直上流点探索; Finding departure point in vertical
use sltt_lagint, only : SLTTIrrHerIntQui1DNonUni, SLTTHerIntCub1D
real(DP), intent(in) :: xyr_SigmaDot(0:imax-1, 1:jmax, 0:kmax)
! 鉛直流速（SigmaDot）
real(DP), intent(in) :: xyzf_QMix (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! 現在時刻の物質混合比
! Present mix ratio of the tracers
real(DP) :: xyzf_QMixA (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! 次ステップの物質混合比
! Next mix ratio of the tracers
!
! local variables
!
real(DP) :: xyz_DPSigma(0:imax-1, 1:jmax, 1:kmax)
! 上流点高度
! Sigma of the departure point
integer:: i ! 東西方向に回る DO ループ用作業変数
! Work variables for DO loop in zonal direction
integer:: j ! 南北方向に回る DO ループ用作業変数
! Work variables for DO loop in meridional direction
integer:: k, kk ! 鉛直方向に回る DO ループ用作業変数
! Work variables for DO loop in vertical direction
integer:: n ! 組成方向に回る DO ループ用作業変数
! Work variables for DO loop in dimension of constituents
integer:: xy_kk(0:imax-1, 1:jmax)
! 上流点の上下のグリッドを探索するための作業変数
! Work variable for finding the grid just above the departure point
real(DP) :: xyzf_QMix_dz(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! 物質混合比の鉛直微分
! Vertical derivative of the mix ratio
real(DP) :: xyzf_ExtQMix(0:imax-1, 1:jmax, 1-2:kmax+2, 1:ncmax)
! 物質混合比の拡張配列
! Extended array of the mix ratio
real(DP) :: z_ExtSigma(1-2:kmax+2)
! σ座標の拡張配列
! Extended array of the sigma coordinate
real(DP) :: xyf_F11(0:imax-1, 1:jmax, 1:ncmax)
! 微分計算時に用いる作業変数
! work variable for the derivative calculation
real(DP) :: xyf_F22(0:imax-1, 1:jmax, 1:ncmax)
! 微分計算時に用いる作業変数
! work variable for the derivative calculation
real(DP) :: xyf_F12(0:imax-1, 1:jmax, 1:ncmax)
! 微分計算時に用いる作業変数
! work variable for the derivative calculation
real(DP) :: xyf_F21(0:imax-1, 1:jmax, 1:ncmax)
! 微分計算時に用いる作業変数
! work variable for the derivative calculation
real(DP) :: s1, t1, s2, t2, r1, r2
! 微分計算時に用いる作業変数
! work variable for the derivative calculation
! 実行文 ; Executable statement
!
! 初期化確認
! Initialization check
!
if ( .not. sltt_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
! 上流点探索
! estimation of departure point
!
call SLTTDPVer( DelTime, xyr_SigmaDot, xyz_DPSigma )
! 配列拡張（z_Sigma）
! Array extension for z_Sigma
z_ExtSigma(-1) = 2.0_DP - z_Sigma(2)
z_ExtSigma(0) = 2.0_DP - z_Sigma(1)
z_ExtSigma(1:kmax) = z_Sigma(1:kmax)
z_ExtSigma(kmax+1) = -z_Sigma(kmax)
z_ExtSigma(kmax+2) = -z_Sigma(kmax-1)
! 配列拡張（xyzf_QMix）
! Array extension for Q_Mix
xyzf_ExtQMix(:,:,-1,:) = xyzf_QMix(:,:,2,:)
xyzf_ExtQMix(:,:,0,:) = xyzf_QMix(:,:,1,:)
xyzf_ExtQMix(:,:,1:kmax,:) = xyzf_QMix(:,:,1:kmax,:)
xyzf_ExtQMix(:,:,kmax+1,:) = xyzf_QMix(:,:,kmax,:)
xyzf_ExtQMix(:,:,kmax+2,:) = xyzf_QMix(:,:,kmax-1,:)
! xyzf_QMix_dz（微分）を求める
! calculate xyzf_QMix_dz
do k = 1 , kmax
s1 = z_ExtSigma(k) - z_ExtSigma(k-1)
t1 = z_ExtSigma(k+1) - z_ExtSigma(k)
s2 = z_ExtSigma(k) - z_ExtSigma(k-2)
t2 = z_ExtSigma(k+2) - z_ExtSigma(k)
if (s1 == t1 .and. s2 == t2 .and. s1 + s1 == s2) then
! 格子が等間隔の場合
! Uniform depth
! 4次精度
! 4th order
xyzf_QMix_dz(:,:,k,:) = ( 8.0_DP*( xyzf_ExtQMix(:,:,k+1,:) - xyzf_ExtQMix(:,:,k-1,:)) - ( xyzf_ExtQMix(:,:,k+2,:) - xyzf_ExtQMix(:,:,k-2,:) ) )/12.0_DP
else
! 格子が不当間隔の場合
! Non-uniform depth
xyf_F11 = (s1*s1*xyzf_ExtQMix(:,:,k+1,:) +(t1*t1 - s1*s1)*xyzf_ExtQMix(:,:,k,:) - t1*t1*xyzf_ExtQMix(:,:,k-1,:)) /(s1*t1*(s1+t1))
xyf_F22 = (s2*s2*xyzf_ExtQMix(:,:,k+2,:) +(t2*t2 - s2*s2)*xyzf_ExtQMix(:,:,k,:) - t2*t2*xyzf_ExtQMix(:,:,k-2,:)) /(s2*t2*(s2+t2))
xyf_F21 = (s2*s2*xyzf_ExtQMix(:,:,k+1,:) +(t1*t1 - s2*s2)*xyzf_ExtQMix(:,:,k,:) - t1*t1*xyzf_ExtQMix(:,:,k-2,:)) /(s2*t1*(s2+t1))
xyf_F12 = (s1*s1*xyzf_ExtQMix(:,:,k+2,:) +(t2*t2 - s1*s1)*xyzf_ExtQMix(:,:,k,:) - t2*t2*xyzf_ExtQMix(:,:,k-1,:)) /(s1*t2*(s1+t2))
r1 = t1 - s1 - t2 + s2
r2 = t1 - s2 - t2 + s1
!４次精度
! 4th order
xyzf_QMix_dz(:,:,k,:) = ( (xyf_F11*s2*t2 - xyf_F22*s1*t1)*r2 - (xyf_F21*s1*t2 - xyf_F12*s2*t1)*r1 ) / ( (s2*t2-s1*t1)*r2 - (s1*t2-s2*t1)*r1 )
!3次精度
! 3rd order
! xyzf_QMix_dz(:,:,k,:) = (xyf_F11*s2*t2 - xyf_F22(:,:,:)*s1*t1)/(s2*t2 - s1*t1)
!2次精度
! 2nd order
! xyzf_QMix_dz(:,:,k,:) = xyf_F11
endif
enddo
xy_kk = 2
do k = 1, kmax
do j = 1, jmax
do i = 0, imax-1
if ( xyz_DPSigma(i,j,k) >= z_Sigma(1) ) then ! DPが z_Sigma(1) と 地表面(sigma = 1.0)の間の場合
! if DP is between z_Sigma(1) and the ground (sigma = 1.0)
xyzf_QMixA(i,j,k,:) = xyzf_QMix(i,j,1,:) ! Q_1で一定とする。
! use Q_1 for interpolated value
elseif (xyz_DPSigma(i,j,k) <= z_Sigma(kmax)) then! DPが z_Sigma(kmax) と 大気上端(sigma = 0.0)の間
! if DP is between z_Sigma(kmax) and the upper boundary (sigma = 0.0)
xyzf_QMixA(i,j,k,:) = xyzf_QMix(i,j,kmax,:) ! Q_kmaxで一定とする。
! use Q_kmax for interpolated value
else
do kk = xy_kk(i,j), kmax
if ( xyz_DPSigma(i,j,k) > z_Sigma(kk) ) then
select case (SLTTIntVer)
case("HQ") ! 変則エルミート５次補間; Irregular Hermite Quintic interpolation
do n = 1, ncmax
xyzf_QMixA(i,j,k,n) = SLTTIrrHerIntQui1DNonUni(xyzf_ExtQMix(i,j,kk-2,n), xyzf_ExtQMix(i,j,kk-1,n), xyzf_ExtQMix(i,j,kk,n), xyzf_ExtQMix(i,j,kk+1,n), xyzf_QMix_dz(i,j,kk-1,n), xyzf_QMix_dz(i,j,kk,n), z_ExtSigma(kk-2)-z_ExtSigma(kk-1), z_ExtSigma(kk)-z_ExtSigma(kk-1), z_ExtSigma(kk+1)-z_ExtSigma(kk-1), xyz_DPSigma(i,j,k)-z_ExtSigma(kk-1))
enddo
case("HC") ! エルミート３次補間; Hermitian Cubic interpolation
do n = 1, ncmax
xyzf_QMixA(i,j,k,n) = SLTTHerIntCub1D( xyzf_ExtQMix(i,j,kk-1,n), xyzf_ExtQMix(i,j,kk,n), xyzf_QMix_dz(i,j,kk-1,n), xyzf_QMix_dz(i,j,kk,n), z_ExtSigma(kk)-z_ExtSigma(kk-1), xyz_DPSigma(i,j,k)-z_ExtSigma(kk-1))
enddo
case default
write( 6, * ) "ERROR : GIVE CORRECT KEYWORD FOR <SLTTIntVer> IN NAMELIST"
stop
end select
xy_kk(i,j) = kk
exit
endif
enddo
endif
enddo
enddo
enddo
end function SLTTVerAdv
| Variable : | |||
| y_ExtLatN(:) : | real(DP) , save, allocatable
|
| Variable : | |||
| y_ExtLatS(:) : | real(DP) , save, allocatable
|