Class | sltt_dp |
In: |
sltt/sltt_dp.f90
|
Note that Japanese and English are described in parallel.
SLTTDPHor : | 水平方向の上流点探索 |
SLTTDPVer : | 鉛直方向の上流点探索 |
——————— : | ———— |
SLTTDPHor : | Finding DP in Horizontal |
SLTTDPVer : | Finding DP in Vertical |
NAMELIST#
Subroutine : | |||
DelTime : | real(DP), intent(in )
| ||
x_Lon(0:imax-1) : | real(DP), intent(in )
| ||
y_Lat(1:jmax/2) : | real(DP), intent(in )
| ||
y_SinLat(1:jmax/2) : | real(DP), intent(in )
| ||
y_CosLat(1:jmax/2) : | real(DP), intent(in )
| ||
iexmin : | integer , intent(in ) | ||
iexmax : | integer , intent(in ) | ||
jexmin : | integer , intent(in ) | ||
jexmax : | integer , intent(in ) | ||
x_ExtLon(iexmin:iexmax) : | real(DP), intent(in )
| ||
y_ExtLat(jexmin:jexmax) : | real(DP), intent(in )
| ||
xyz_ExtU(iexmin:iexmax, jexmin:jexmax, 1:kmax) : | real(DP), intent(in )
| ||
xyz_ExtV(iexmin:iexmax, jexmin:jexmax, 1:kmax) : | real(DP), intent(in )
| ||
xyz_DPLon(0:imax-1, 1:jmax/2, 1:kmax) : | real(DP), intent(out)
| ||
xyz_DPLat(0:imax-1, 1:jmax/2, 1:kmax) : | real(DP), intent(out)
|
水平方向の上流点探索 Finding DP in Horizontal
subroutine SLTTDPHor( DelTime, x_Lon, y_Lat, y_SinLat, y_CosLat, iexmin, iexmax, jexmin, jexmax, x_ExtLon, y_ExtLat, xyz_ExtU, xyz_ExtV, xyz_DPLon, xyz_DPLat ) ! 水平方向の上流点探索 ! Finding DP in Horizontal use sltt_const, only : PIx2, dtjw, nloop_dp_h real(DP), intent(in ) :: DelTime ! $\Delta t$ real(DP), intent(in ) :: x_Lon(0:imax-1) ! $\lambda$ longitude real(DP), intent(in ) :: y_Lat(1:jmax/2) ! $\varphi$ latitude real(DP), intent(in ) :: y_SinLat(1:jmax/2) ! $\sin\varphi$ real(DP), intent(in ) :: y_CosLat(1:jmax/2) ! $\cos\varphi$ integer , intent(in ) :: iexmin integer , intent(in ) :: iexmax integer , intent(in ) :: jexmin integer , intent(in ) :: jexmax real(DP), intent(in ) :: x_ExtLon(iexmin:iexmax) ! 経度の拡張配列 ! Extended array of Lon real(DP), intent(in ) :: y_ExtLat(jexmin:jexmax) ! 緯度の拡張配列 ! Extended array of Lat real(DP), intent(in ) :: xyz_ExtU(iexmin:iexmax, jexmin:jexmax, 1:kmax) ! 東西風速の拡張配列 ! Extended array of zonal velocity real(DP), intent(in ) :: xyz_ExtV(iexmin:iexmax, jexmin:jexmax, 1:kmax) ! 南北風速の拡張配列 ! Extended array of meridional velocity real(DP), intent(out) :: xyz_DPLon(0:imax-1, 1:jmax/2, 1:kmax) ! 上流点の経度 ! Lon of departure point real(DP), intent(out) :: xyz_DPLat(0:imax-1, 1:jmax/2, 1:kmax) ! 上流点の緯度 ! Lat of departure point ! ! local variables ! real(DP) :: xyz_MPLon (0:imax-1, 1:jmax/2, 1:kmax) ! 中間点の経度 ! Lon of mid-point real(DP) :: xyz_MPLat (0:imax-1, 1:jmax/2, 1:kmax) ! 中間点の緯度 ! Lat of mid-point real(DP) :: xyz_MPLonN(0:imax-1, 1:jmax/2, 1:kmax) ! 中間点の経度(現在推定値) ! Lon of mid-point (temporal estimation) real(DP) :: xyz_MPLatN(0:imax-1, 1:jmax/2, 1:kmax) ! 中間点の緯度(現在推定値) ! Lat of mid-point (now) logical :: FlagWindInt ! 初期化フラッグ ! Flag for initialization logical :: xyz_FlagEstDP(0:imax-1, 1:jmax/2, 1:kmax) 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 :: t ! 推定回数の DO ループ用作業変数 ! Work variables for DO loop for DP estimation ! 初期化処理 ! initialization ! xyz_FlagEstDP = .true. do k = 1, kmax do j = 1, jmax/2 do i = 0, imax-1 xyz_MPLonN(i,j,k) = x_Lon(i) xyz_MPLatN(i,j,k) = y_Lat(j) end do end do end do FlagWindInt = .false. mp_search_loop : do t = 1, nloop_dp_h if( mod( t, 5 ) .eq. 0 ) write( 6, * ) "### Loop for searching departure points : ", t xyz_MPLon = xyz_MPLonN xyz_MPLat = xyz_MPLatN call SLTTDPHorCore( DelTime, x_Lon, y_Lat, y_SinLat, y_CosLat, iexmin, iexmax, jexmin, jexmax, x_ExtLon, y_ExtLat, xyz_ExtU, xyz_ExtV, xyz_MPLon, xyz_MPLat, xyz_MPLonN, xyz_MPLatN, xyz_FlagEstDP, FlagWindInt ) FlagWindInt = .true. end do mp_search_loop xyz_FlagEstDP = .true. FlagWindInt = .true. call SLTTDPHorCore( 2.0_DP * DelTime, x_Lon, y_Lat, y_SinLat, y_CosLat, iexmin, iexmax, jexmin, jexmax, x_ExtLon, y_ExtLat, xyz_ExtU, xyz_ExtV, xyz_MPLon, xyz_MPLat, xyz_DPLon, xyz_DPLat, xyz_FlagEstDP, FlagWindInt ) end subroutine SLTTDPHor
Subroutine : | |||
DelTime : | real(DP), intent(in )
| ||
xyr_SigmaDot(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in )
| ||
xyz_DPSigma(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
|
鉛直方向の上流点探索 Finding DP in Vertical
subroutine SLTTDPVer( DelTime, xyr_SigmaDot, xyz_DPSigma ) ! 鉛直方向の上流点探索 ! Finding DP in Vertical ! 座標データ設定 ! Axes data settings ! use axesset, only : r_Sigma, z_Sigma use sltt_const, only : nloop_dp_v real(DP), intent(in ) :: DelTime ! $\Delta t$ real(DP), intent(in ) :: xyr_SigmaDot(0:imax-1, 1:jmax, 0:kmax) ! 鉛直流速(SigmaDot) real(DP), intent(out) :: xyz_DPSigma (0:imax-1, 1:jmax, 1:kmax) ! 上流点高度 ! Sigma of the departure point ! ! local variables ! real(DP) :: xyz_MPSigma (0:imax-1, 1:jmax, 1:kmax) ! 中間点の高度 ! Sigma of mid-point real(DP) :: xyz_MPSigmaN (0:imax-1, 1:jmax, 1:kmax) ! 中間点の高度(現在推定値) ! Lat of mid-point (temporal estimation) real(DP) :: xyz_MPSigmaDot(0:imax-1, 1:jmax, 1:kmax) ! 中間点のSigmaDot ! SigmaDot of mid-point real(DP) :: xyza_lcifz(0:imax-1, 1:jmax, 1:kmax, -1:2) ! ラグランジュ補間のための作業変数 ! Work variable for Lagrange interpolation integer :: xyz_kk(0:imax-1, 1:jmax, 1:kmax) ! ラグランジュ補間のための作業変数 ! Work variable for Lagrange interpolation 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, k2 ! 鉛直方向に回る DO ループ用作業変数 ! Work variables for DO loop in vertical direction integer :: t ! 推定回数の DO ループ用作業変数 ! Work variables for DO loop for DP estimation ! 初期化処理 ! initialization ! do k = 1, kmax xyz_MPSigmaN(:,:,k) = z_Sigma(k) end do ! 上流点推定のループ ! Loop for finding departure point mp_search_loop : do t = 1, nloop_dp_v xyz_MPSigma = xyz_MPSigmaN ! 中間点の鉛直流速をラグランジュ3次補間で求める ! vertical wind velocity at mid-point is estimated by using Lagrange cubic ! interpolation ! do k = 1, kmax do j = 1, jmax do i = 0, imax-1 ! ! Routine for dcpam ! ! Departure points, xyz_MPSigma(:,:,k), must be located between ! r_Sigma(kk) > xyz_MPSigma(k) > r_Sigma(kk+1). ! Further, 1 <= kk <= kmax-2. ! ! ! economical method ! if( xyz_MPSigma(i,j,k) > r_Sigma(k) ) then k_search_1 : do k2 = k, 1, -1 if( r_Sigma(k2) > xyz_MPSigma(i,j,k) ) exit k_search_1 end do k_search_1 xyz_kk(i,j,k) = k2 else k_search_2 : do k2 = min( k+1, kmax ), kmax if( r_Sigma(k2) < xyz_MPSigma(i,j,k) ) exit k_search_2 end do k_search_2 xyz_kk(i,j,k) = k2 - 1 end if xyz_kk(i,j,k) = min( max( xyz_kk(i,j,k), 1 ), kmax-2 ) end do end do end do ! do k = 1, kmax do j = 1, jmax do i = 0, imax-1 kk = xyz_kk(i,j,k) xyza_lcifz(i,j,k,-1) = ( xyz_MPSigma(i,j,k) - r_Sigma(kk ) ) * ( xyz_MPSigma(i,j,k) - r_Sigma(kk+1) ) * ( xyz_MPSigma(i,j,k) - r_Sigma(kk+2) ) / ( ( r_Sigma(kk-1) - r_Sigma(kk ) ) * ( r_Sigma(kk-1) - r_Sigma(kk+1) ) * ( r_Sigma(kk-1) - r_Sigma(kk+2) ) ) xyza_lcifz(i,j,k, 0) = ( xyz_MPSigma(i,j,k) - r_Sigma(kk-1) ) * ( xyz_MPSigma(i,j,k) - r_Sigma(kk+1) ) * ( xyz_MPSigma(i,j,k) - r_Sigma(kk+2) ) / ( ( r_Sigma(kk ) - r_Sigma(kk-1) ) * ( r_Sigma(kk ) - r_Sigma(kk+1) ) * ( r_Sigma(kk ) - r_Sigma(kk+2) ) ) xyza_lcifz(i,j,k, 1) = ( xyz_MPSigma(i,j,k) - r_Sigma(kk-1) ) * ( xyz_MPSigma(i,j,k) - r_Sigma(kk ) ) * ( xyz_MPSigma(i,j,k) - r_Sigma(kk+2) ) / ( ( r_Sigma(kk+1) - r_Sigma(kk-1) ) * ( r_Sigma(kk+1) - r_Sigma(kk ) ) * ( r_Sigma(kk+1) - r_Sigma(kk+2) ) ) xyza_lcifz(i,j,k, 2) = ( xyz_MPSigma(i,j,k) - r_Sigma(kk-1) ) * ( xyz_MPSigma(i,j,k) - r_Sigma(kk ) ) * ( xyz_MPSigma(i,j,k) - r_Sigma(kk+1) ) / ( ( r_Sigma(kk+2) - r_Sigma(kk-1) ) * ( r_Sigma(kk+2) - r_Sigma(kk ) ) * ( r_Sigma(kk+2) - r_Sigma(kk+1) ) ) end do end do end do ! do k = 1, kmax do j = 1, jmax do i = 0, imax-1 kk = xyz_kk(i,j,k) xyz_MPSigmaDot(i,j,k) = xyza_lcifz(i,j,k,-1) * xyr_SigmaDot(i,j,kk-1) + xyza_lcifz(i,j,k, 0) * xyr_SigmaDot(i,j,kk ) + xyza_lcifz(i,j,k, 1) * xyr_SigmaDot(i,j,kk+1) + xyza_lcifz(i,j,k, 2) * xyr_SigmaDot(i,j,kk+2) end do end do end do do k = 1, kmax xyz_MPSigmaN(:,:,k) = z_Sigma(k) - xyz_MPSigmaDot(:,:,k) * DelTime end do xyz_MPSigmaN = min( max( xyz_MPSigmaN, r_Sigma(kmax) ), r_Sigma(0) ) end do mp_search_loop ! 上流点高度の推定 ! estimating departure point ! do k = 1, kmax xyz_DPSigma(:,:,k) = z_Sigma(k) - xyz_MPSigmaDot(:,:,k) * 2.0_DP * DelTime end do end subroutine SLTTDPVer
Subroutine : | |||
DelTime : | real(DP), intent(in )
| ||
x_APLon(0:imax-1) : | real(DP), intent(in )
| ||
y_APLat(1:jmax/2 ) : | real(DP), intent(in )
| ||
y_APSinLat(1:jmax/2) : | real(DP), intent(in )
| ||
y_APCosLat(1:jmax/2) : | real(DP), intent(in )
| ||
iexmin : | integer , intent(in ) | ||
iexmax : | integer , intent(in ) | ||
jexmin : | integer , intent(in ) | ||
jexmax : | integer , intent(in ) | ||
x_ExtLon(iexmin:iexmax) : | real(DP), intent(in )
| ||
y_ExtLat(jexmin:jexmax) : | real(DP), intent(in )
| ||
xyz_ExtU(iexmin:iexmax, jexmin:jexmax, 1:kmax) : | real(DP), intent(in )
| ||
xyz_ExtV(iexmin:iexmax, jexmin:jexmax, 1:kmax) : | real(DP), intent(in )
| ||
xyz_MPLon(0:imax-1, 1:jmax/2, 1:kmax) : | real(DP), intent(in )
| ||
xyz_MPLat(0:imax-1, 1:jmax/2, 1:kmax) : | real(DP), intent(in )
| ||
xyz_DPLon(0:imax-1, 1:jmax/2, 1:kmax) : | real(DP), intent(inout)
| ||
xyz_DPLat(0:imax-1, 1:jmax/2, 1:kmax) : | real(DP), intent(inout)
| ||
xyz_FlagEstDP(0:imax-1, 1:jmax/2, 1:kmax) : | logical , intent(in )
| ||
FlagWindInt : | logical , intent(in )
|
水平方向の上流点探索のコア部分 Finding DP in Horizontal (Core)
subroutine SLTTDPHorCore( DelTime, x_APLon, y_APLat, y_APSinLat, y_APCosLat, iexmin, iexmax, jexmin, jexmax, x_ExtLon, y_ExtLat, xyz_ExtU, xyz_ExtV, xyz_MPLon, xyz_MPLat, xyz_DPLon, xyz_DPLat, xyz_FlagEstDP, FlagWindInt ) ! 水平方向の上流点探索のコア部分 ! Finding DP in Horizontal (Core) use constants0 , only : PI use constants , only : RPlanet use mpi_wrapper, only : myrank use sltt_const , only : PIx2, PIH use sltt_lagint, only : SLTTLagIntCubCalcFactHor , SLTTLagIntCubIntHor real(DP), intent(in ) :: DelTime ! $\Delta t$ real(DP), intent(in ) :: x_APLon(0:imax-1) ! 到着点(グリッド上)の経度 ! Lon of arrival point (which is on grid) real(DP), intent(in ) :: y_APLat(1:jmax/2 ) ! 到着点(グリッド上)の緯度 ! Lat of arrival point (which is on grid) real(DP), intent(in ) :: y_APSinLat(1:jmax/2) ! 到着点の sin(緯度) ! sin(lat) of arrival point real(DP), intent(in ) :: y_APCosLat(1:jmax/2) ! 到着点の cos(緯度) ! cos(lat) of arrival point integer , intent(in ) :: iexmin integer , intent(in ) :: iexmax integer , intent(in ) :: jexmin integer , intent(in ) :: jexmax real(DP), intent(in ) :: x_ExtLon(iexmin:iexmax) ! 経度の拡張配列 ! Extended array of Lon real(DP), intent(in ) :: y_ExtLat(jexmin:jexmax) ! 緯度の拡張配列 ! Extended array of Lat real(DP), intent(in ) :: xyz_ExtU(iexmin:iexmax, jexmin:jexmax, 1:kmax) ! 東西風速の拡張配列 ! Extended array of zonal wind real(DP), intent(in ) :: xyz_ExtV(iexmin:iexmax, jexmin:jexmax, 1:kmax) ! 南北風速の拡張配列 ! Extended array of meridional wind real(DP), intent(in ) :: xyz_MPLon(0:imax-1, 1:jmax/2, 1:kmax) ! 中間点の経度 ! Lon of mid-point real(DP), intent(in ) :: xyz_MPLat(0:imax-1, 1:jmax/2, 1:kmax) ! 中間点の緯度 ! Lat of mid-point real(DP), intent(inout) :: xyz_DPLon(0:imax-1, 1:jmax/2, 1:kmax) ! 上流点の経度 ! Lon of departure point real(DP), intent(inout) :: xyz_DPLat(0:imax-1, 1:jmax/2, 1:kmax) ! 上流点の緯度 ! Lat of departure point logical , intent(in ) :: xyz_FlagEstDP(0:imax-1, 1:jmax/2, 1:kmax) ! 上流点探索のフラグ ! Flag for finding departure point logical , intent(in ) :: FlagWindInt ! 初回フラグ ! Initial flag ! ! local variables ! real(DP) :: xyz_MPU(0:imax-1, 1:jmax/2, 1:kmax) ! 中間点の東西風速 ! Zonal wind at mid-point real(DP) :: xyz_MPV(0:imax-1, 1:jmax/2, 1:kmax) ! 中間点の南北風速 ! Meridional wind at mid-point real(DP) :: MPCosLat ! 中間点のcos(lat) ! cos(lat) at mid-point real(DP) :: xyza_lcifx(0:imax-1, 1:jmax/2, 1:kmax, -1:2) real(DP) :: xyza_lcify(0:imax-1, 1:jmax/2, 1:kmax, -1:2) integer :: xyz_ii(0:imax-1, 1:jmax/2, 1:kmax) integer :: xyz_jj(0:imax-1, 1:jmax/2, 1:kmax) ! ラグランジュ補間のための作業変数 ! Working variables for Lagrange interpolation 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 ループ用作業変数 real(DP) :: InvRPlanet InvRPlanet = 1.0_DP/RPlanet if( FlagWindInt ) then ! ラグランジュ3次補間 ! Cubic interpolation call SLTTLagIntCubCalcFactHor( iexmin, iexmax, jexmin, jexmax, x_ExtLon, y_ExtLat, xyz_MPLon, xyz_MPLat, xyz_ii, xyz_jj, xyza_lcifx, xyza_lcify ) ! 中間点の風速推定 ! estimating wind velocities at mid-point ! call SLTTLagIntCubIntHor( iexmin, iexmax, jexmin, jexmax, xyz_ii, xyz_jj, xyza_lcifx, xyza_lcify, xyz_ExtU, xyz_MPU ) call SLTTLagIntCubIntHor( iexmin, iexmax, jexmin, jexmax, xyz_ii, xyz_jj, xyza_lcifx, xyza_lcify, xyz_ExtV, xyz_MPV ) else ! ! 推定風速の更新 ! xyz_MPU(0:imax-1,1:jmax/2,:) = xyz_ExtU(0:imax-1,1:jmax/2,:) xyz_MPV(0:imax-1,1:jmax/2,:) = xyz_ExtV(0:imax-1,1:jmax/2,:) end if do k = 1, kmax do j = 1, jmax/2 do i = 0, imax-1 if( xyz_FlagEstDP(i,j,k) ) then ! MPCosLat = cos( xyz_MPLat(i,j,k) ) xyz_DPLon(i,j,k) = x_APLon(i) - DelTime * xyz_MPU(i,j,k) / cos( xyz_MPLat(i,j,k) ) * InvRPlanet !/ RPlanet xyz_DPLat(i,j,k) = y_APLat(j) - DelTime * xyz_MPV(i,j,k)* InvRPlanet ! / RPlanet end if end do end do end do ! 普通の緯度経度の範囲に直す do k = 1, kmax do j = 1, jmax/2 do i = 0, imax-1 if( xyz_DPLat(i,j,k) < -PIH ) then xyz_DPLon(i,j,k) = xyz_DPLon(i,j,k) + PI xyz_DPLat(i,j,k) = -PIH + ( -PIH - xyz_DPLat(i,j,k) ) else if( xyz_DPLat(i,j,k) > PIH ) then xyz_DPLon(i,j,k) = xyz_DPLon(i,j,k) + PI xyz_DPLat(i,j,k) = PIH - ( xyz_DPLat(i,j,k) - PIH ) end if xyz_DPLon(i,j,k) = mod( xyz_DPLon(i,j,k) + PIx2, PIx2 ) end do end do end do end subroutine SLTTDPHorCore