Class | vertical_filter |
In: |
util/vertical_filter.f90
|
Note that Japanese and English are described in parallel.
鉛直モードを抑えるための鉛直フィルター (Ishiwatari et al., 2002) を適用するためのモジュールです.
すなわち begin{align}
T_{B,k} &\equiv \frac{T_{k+1/2} + T_{k-1/2}}{2} \\
end{align} と定義した上で, begin{align}
\bar{T}_{k-1} &= T_{k-1} + F_T ( T_{B, k-1} - T_{k-1} ) - \frac{ \Sigma_{k'=k-1}^{k+1} F_T (T_{B,k'} - T_{k'}) \delta p_{k'} }{\Sigma_{k'=k-1}^{k+1} \delta p_{k'}}, \\ \bar{T}_{k} &= T_{k} + F_T ( T_{B, k} - T_{k} ) - \frac{ \Sigma_{k'=k-1}^{k+1} F_T (T_{B,k'} - T_{k'}) \delta p_{k'} }{\Sigma_{k'=k-1}^{k+1} \delta p_{k'}}, \\ \bar{T}_{k+1} &= T_{k+1} + F_T ( T_{B, k+1} - T_{k+1} ) - \frac{ \Sigma_{k'=k-1}^{k+1} F_T (T_{B,k'} - T_{k'}) \delta p_{k'} }{\Sigma_{k'=k-1}^{k+1} \delta p_{k'}}
end{align} という操作を下の層から順に適用することで $ bar{T} $ を求めます. ここで添え字 $ k $ は層の番号です. 鉛直フィルターが適用された項には上線 ( $ bar{quad} $ ) が付いています.
鉛直フィルターの係数 $F_T$ は Create の際に設定します.
++
Vertical filter (Ishiwatari et al., 2002) for suppression of vertical mode is applied.
Concretely, $ bar{A}^{t} $ is derived as follows. Reference Temperature is defined as begin{align}
T_{B,k} &\equiv \frac{T_{k+1/2} + T_{k-1/2}}{2}.
end{align} Operations as follows are applied from lower layer to higher layer. begin{align}
\bar{T}_{k-1} &= T_{k-1} + F_T ( T_{B, k-1} - T_{k-1} ) - \frac{ \Sigma_{k'=k-1}^{k+1} F_T (T_{B,k'} - T_{k'}) \delta p_{k'} }{\Sigma_{k'=k-1}^{k+1} \delta p_{k'}}, \\ \bar{T}_{k} &= T_{k} + F_T ( T_{B, k} - T_{k} ) - \frac{ \Sigma_{k'=k-1}^{k+1} F_T (T_{B,k'} - T_{k'}) \delta p_{k'} }{\Sigma_{k'=k-1}^{k+1} \delta p_{k'}}, \\ \bar{T}_{k+1} &= T_{k+1} + F_T ( T_{B, k+1} - T_{k+1} ) - \frac{ \Sigma_{k'=k-1}^{k+1} F_T (T_{B,k'} - T_{k'}) \delta p_{k'} }{\Sigma_{k'=k-1}^{k+1} \delta p_{k'}}
end{align} where a suffix $ k $ is the number of layer. Over-bar ( $ bar{quad} $ ) represents the terms applied vertical filter.
Vertical filter coefficient $ F_T $ is specified by "Create".
VerticalFilter : | 鉛直フィルターの適用 |
—————- : | ——————— |
VerticalFilter : | Apply vertical filter |
Subroutine : | |||
xyz_UA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout)
| ||
xyz_VA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout)
| ||
xyz_TempA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout)
| ||
xyr_Press(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in)
| ||
xyr_Temp(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in)
|
予報変数に鉛直フィルターを適用します.
Vertical filter is applied to predictional variables.
subroutine VerticalFilter( xyz_UA, xyz_VA, xyz_TempA, xyr_Press, xyr_Temp ) ! ! 予報変数に鉛直フィルターを適用します. ! ! Vertical filter is applied to predictional variables. ! ! モジュール引用 ; USE statements ! ! 時刻管理 ! Time control ! use timeset, only: TimesetClockStart, TimesetClockStop ! 宣言文 ; Declaration statements ! implicit none real(DP), intent(inout):: xyz_UA (0:imax-1, 1:jmax, 1:kmax) ! $ u (t+\Delta t) $ . 東西風速. Eastward wind real(DP), intent(inout):: xyz_VA (0:imax-1, 1:jmax, 1:kmax) ! $ v (t+\Delta t) $ . 南北風速. Northward wind real(DP), intent(inout):: xyz_TempA (0:imax-1, 1:jmax, 1:kmax) ! $ T (t+\Delta t) $ . 温度. Temperature real(DP), intent(in):: xyr_Press (0:imax-1, 1:jmax, 0:kmax) ! $ p_{1/2} (t+\Delta t) $ . 気圧. Pressure real(DP), intent(in):: xyr_Temp (0:imax-1, 1:jmax, 0:kmax) ! $ T_{1/2} (t+\Delta t) $ . 温度. Temperature ! 作業変数 ! Work variables ! integer:: i ! 経度方向に回る DO ループ用作業変数 ! Work variables for DO loop in longitude integer:: j ! 緯度方向に回る DO ループ用作業変数 ! Work variables for DO loop in latitude integer:: k ! 鉛直方向に回る DO ループ用作業変数 ! Work variables for DO loop in vertical direction integer:: k1 ! 鉛直方向に回る DO ループ用作業変数 ! Work variables for DO loop in vertical direction real(DP):: Temp_ref(0:imax-1, 1:jmax, 1:kmax), Temp_corr(0:imax-1, 1:jmax) ! 実行文 ; Executable statement ! ! 計算時間計測開始 ! Start measurement of computation time ! call TimesetClockStart( module_name ) ! 初期化 ! Initialization ! if ( .not. vertical_filter_inited ) call VerticalFilterInit !!$ ! ステップのチェック !!$ ! Check Step !!$ ! !!$ CurStep = CurStep + 1 !!$ if ( CurStep - PrevStep < IntStep ) return !!$ PrevStep = CurStep ! FilterParameter の有効性のチェック ! Check validity of "FilterParameter" ! if ( FilterParameter == 0.0_DP ) return ! 鉛直フィルターを適用 ! Apply vertical filter ! do k=1,kmax Temp_ref(:,:,k) = 0.5_DP * ( xyr_Temp(:,:,k-1) + xyr_Temp(:,:,k) ) end do do k=2,kmax-1 Temp_corr = 0.0_DP do k1=-1,1 Temp_corr = Temp_corr + FilterParameter * ( Temp_ref(:,:,k+k1) - xyz_TempA(:,:,k+k1) ) * ( xyr_Press(:,:,k+k1-1) - xyr_Press(:,:,k+k1) ) end do Temp_corr = Temp_corr / ( xyr_Press(:,:,k-2) - xyr_Press(:,:,k+1) ) do k1=-1,1 xyz_TempA(:,:,k+k1) = xyz_TempA(:,:,k+k1) + FilterParameter * ( Temp_ref(:,:,k+k1) - xyz_TempA(:,:,k+k1) ) - Temp_corr end do end do ! 計算時間計測一時停止 ! Pause measurement of computation time ! call TimesetClockStop( module_name ) end subroutine VerticalFilter
Variable : | |||
vertical_filter_inited = .false. : | logical, save, public
|
Variable : | |||||||||
FilterParameter : | real(DP), save
|
Subroutine : |
vertical_filter モジュールの初期化を行います. NAMELIST#vertical_filter_nml の読み込みはこの手続きで行われます.
"vertical_filter" module is initialized. "NAMELIST#vertical_filter_nml" is loaded in this procedure.
This procedure input/output NAMELIST#vertical_filter_nml .
subroutine VerticalFilterInit ! ! vertical_filter モジュールの初期化を行います. ! NAMELIST#vertical_filter_nml の読み込みはこの手続きで行われます. ! ! "vertical_filter" module is initialized. ! "NAMELIST#vertical_filter_nml" is loaded in this procedure. ! ! モジュール引用 ; USE statements ! ! NAMELIST ファイル入力に関するユーティリティ ! Utilities for NAMELIST file input ! use namelist_util, only: namelist_filename, NmlutilMsg, NmlutilAryValid ! ファイル入出力補助 ! File I/O support ! use dc_iounit, only: FileOpen ! 種別型パラメタ ! Kind type parameter ! use dc_types, only: STDOUT ! 標準出力の装置番号. Unit number of standard output ! 宣言文 ; Declaration statements ! implicit none integer:: unit_nml ! NAMELIST ファイルオープン用装置番号. ! Unit number for NAMELIST file open integer:: iostat_nml ! NAMELIST 読み込み時の IOSTAT. ! IOSTAT of NAMELIST read ! NAMELIST 変数群 ! NAMELIST group name ! namelist /vertical_filter_nml/ FilterParameter !, IntStep ! ! デフォルト値については初期化手続 "vertical_filter#VerticalFilterInit" ! のソースコードを参照のこと. ! ! Refer to source codes in the initialization procedure ! "vertical_filter#VerticalFilterInit" for the default values. ! ! 実行文 ; Executable statement ! if ( vertical_filter_inited ) return ! デフォルト値の設定 ! Default values settings ! FilterParameter = 0.1_DP !!$ IntStep = 1 !!$ CurStep = 0 ! 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 = vertical_filter_nml, iostat = iostat_nml ) ! (out) close( unit_nml ) call NmlutilMsg( iostat_nml, module_name ) ! (in) if ( iostat_nml == 0 ) write( STDOUT, nml = vertical_filter_nml ) end if !!$ ! IntStep のチェック !!$ ! Check "IntStep" !!$ ! !!$ if ( IntStep < 1 ) then !!$ call MessageNotify( 'E', module_name, & !!$ & 'IntStep=<%d> must be greater than 0', & !!$ & i = (/ IntStep /) ) !!$ end if !!$ !!$ ! PrevStep の設定 !!$ ! Configure "PrevStep" !!$ ! !!$ PrevStep = - IntStep ! 印字 ; Print ! call MessageNotify( 'M', module_name, '----- Initialization Messages -----' ) call MessageNotify( 'M', module_name, ' FilterParameter = %f', d = (/ FilterParameter /) ) !!$ call MessageNotify( 'M', module_name, ' IntStep = %d', i = (/ IntStep /) ) call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) ) vertical_filter_inited = .true. end subroutine VerticalFilterInit
Constant : | |||
module_name = ‘vertical_filter‘ : | character(*), parameter
|
Constant : | |||
version = ’$Name: dcpam5-20120921 $’ // ’$Id: vertical_filter.f90,v 1.2 2011-06-19 10:56:28 yot Exp $’ : | character(*), parameter
|