Class | timefilter_asselin1972 |
In: |
util/timefilter_asselin1972.f90
|
Note that Japanese and English are described in parallel.
計算モードを抑えるためのタイムフィルター (Asselin, 1972) を適用するためのモジュールです.
すなわち \[
\bar{A}^{t} = {A}^{t} + \epsilon_f \left( \bar{A}^{t-\Delta t} - 2{A}^{t} + {A}^{t+\Delta t} \right) \] \[ \qquad = ( 1-2 \epsilon_f ) {A}^{t} + \epsilon_f \left( \bar{A}^{t-\Delta t} + {A}^{t+\Delta t} \right)
\] として $ bar{A}^{t} $ を求めます. 添え字 $ t-\Delta t, t, t+Delta t $ はそれぞれ時間ステップを表し, タイムフィルターが適用された項には上線 ( $ bar{quad} $ ) が付いています.
タイムフィルターの係数 $epsilon_f$ は Create の際に設定します.
Time filter (Asselin, 1972) for suppression of computational mode is applied.
Concretely, $ bar{A}^{t} $ is derived as follows. \[
\bar{A}^{t} = {A}^{t} + \epsilon_f \left( \bar{A}^{t-\Delta t} - 2{A}^{t} + {A}^{t+\Delta t} \right) \] \[ \qquad = ( 1-2 \epsilon_f ) {A}^{t} + \epsilon_f \left( \bar{A}^{t-\Delta t} + {A}^{t+\Delta t} \right),
\] where suffices $ t-\Delta t, t, t+Delta t $ represent time step, and over-bar $ bar{quad} $ represents the terms applied time filter.
Time filter coefficient $ epsilon_f $ is specified by "Create".
TimeFilter : | タイムフィルターの適用 |
———— : | ———— |
TimeFilter : | Apply time filter |
NAMELIST#timefilter_asselin1972_nml
Subroutine : | |||
xyz_UB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
xyz_VB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
xyz_TempB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
xy_PsB(0:imax-1, 1:jmax) : | real(DP), intent(in)
| ||
xyz_UN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout)
| ||
xyz_VN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout)
| ||
xyz_TempN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(inout)
| ||
xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(inout)
| ||
xy_PsN(0:imax-1, 1:jmax) : | real(DP), intent(inout)
| ||
xyz_UA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
xyz_VA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
xyz_TempA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
xy_PsA(0:imax-1, 1:jmax) : | real(DP), intent(in)
|
予報変数にタイムフィルターを適用します.
Time filter is applied to predictional variables.
subroutine TimeFilter( xyz_UB, xyz_VB, xyz_TempB, xyzf_QMixB, xy_PsB, xyz_UN, xyz_VN, xyz_TempN, xyzf_QMixN, xy_PsN, xyz_UA, xyz_VA, xyz_TempA, xyzf_QMixA, xy_PsA ) ! ! 予報変数にタイムフィルターを適用します. ! ! Time filter is applied to predictional variables. ! ! モジュール引用 ; USE statements ! ! 時刻管理 ! Time control ! use timeset, only: TimesetClockStart, TimesetClockStop ! 宣言文 ; Declaration statements ! implicit none real(DP), intent(in):: xyz_UB (0:imax-1, 1:jmax, 1:kmax) ! $ u (t-\Delta t) $ . 東西風速. Eastward wind real(DP), intent(in):: xyz_VB (0:imax-1, 1:jmax, 1:kmax) ! $ v (t-\Delta t) $ . 南北風速. Northward wind real(DP), intent(in):: xyz_TempB (0:imax-1, 1:jmax, 1:kmax) ! $ T (t-\Delta t) $ . 温度. Temperature real(DP), intent(in):: xyzf_QMixB (0:imax-1, 1:jmax, 1:kmax, 1:ncmax) ! $ q (t-\Delta t) $ . 混合比. Mass mixing ratio of constituents real(DP), intent(in):: xy_PsB (0:imax-1, 1:jmax) ! $ p_s (t-\Delta t) $ . 地表面気圧. Surface pressure real(DP), intent(inout):: xyz_UN (0:imax-1, 1:jmax, 1:kmax) ! $ u (t) $ . 東西風速. Eastward wind real(DP), intent(inout):: xyz_VN (0:imax-1, 1:jmax, 1:kmax) ! $ v (t) $ . 南北風速. Northward wind real(DP), intent(inout):: xyz_TempN (0:imax-1, 1:jmax, 1:kmax) ! $ T (t) $ . 温度. Temperature real(DP), intent(inout):: xyzf_QMixN (0:imax-1, 1:jmax, 1:kmax, 1:ncmax) ! $ q (t) $ . 混合比. Mass mixing ratio of constituents real(DP), intent(inout):: xy_PsN (0:imax-1, 1:jmax) ! $ p_s (t) $ . 地表面気圧. Surface pressure real(DP), intent(in):: xyz_UA (0:imax-1, 1:jmax, 1:kmax) ! $ u (t+\Delta t) $ . 東西風速. Eastward wind real(DP), intent(in):: xyz_VA (0:imax-1, 1:jmax, 1:kmax) ! $ v (t+\Delta t) $ . 南北風速. Northward wind real(DP), intent(in):: xyz_TempA (0:imax-1, 1:jmax, 1:kmax) ! $ T (t+\Delta t) $ . 温度. Temperature real(DP), intent(in):: xyzf_QMixA (0:imax-1, 1:jmax, 1:kmax, 1:ncmax) ! $ q (t+\Delta t) $ . 混合比. Mass mixing ratio of constituents real(DP), intent(in):: xy_PsA (0:imax-1, 1:jmax) ! $ p_s (t+\Delta t) $ . 地表面気圧. Surface pressure ! 作業変数 ! Work variables ! ! 実行文 ; Executable statement ! ! 計算時間計測開始 ! Start measurement of computation time ! call TimesetClockStart( module_name ) ! 初期化 ! Initialization ! if ( .not. timefilter_asselin1972_inited ) call TimeFiltInit !!$ ! ステップのチェック !!$ ! Check Step !!$ ! !!$ CurStep = CurStep + 1 !!$ if ( CurStep - PrevStep < IntStep ) return !!$ PrevStep = CurStep ! FilterParameter の有効性のチェック ! Check validity of "FilterParameter" ! if ( FilterParameter == 0.0_DP ) return ! タイムフィルターを適用 ! Apply time filter ! xyz_UN = xyz_UN + FilterParameter * ( xyz_UB - 2.0d0 * xyz_UN + xyz_UA ) xyz_VN = xyz_VN + FilterParameter * ( xyz_VB - 2.0d0 * xyz_VN + xyz_VA ) xyz_TempN = xyz_TempN + FilterParameter * ( xyz_TempB - 2.0d0 * xyz_TempN + xyz_TempA ) xy_PsN = xy_PsN + FilterParameter * ( xy_PsB - 2.0d0 * xy_PsN + xy_PsA ) xyzf_QMixN = xyzf_QMixN + FilterParameter * ( xyzf_QMixB - 2.0d0 * xyzf_QMixN + xyzf_QMixA ) ! 計算時間計測一時停止 ! Pause measurement of computation time ! call TimesetClockStop( module_name ) end subroutine TimeFilter
Subroutine : | |||
xy_SoilMoistB(0:imax-1, 1:jmax) : | real(DP), intent(in )
| ||
xy_SurfSnowB(0:imax-1, 1:jmax) : | real(DP), intent(in )
| ||
xy_SoilMoistN(0:imax-1, 1:jmax) : | real(DP), intent(inout)
| ||
xy_SurfSnowN(0:imax-1, 1:jmax) : | real(DP), intent(inout)
| ||
xy_SoilMoistA(0:imax-1, 1:jmax) : | real(DP), intent(in )
| ||
xy_SurfSnowA(0:imax-1, 1:jmax) : | real(DP), intent(in )
|
予報変数にタイムフィルターを適用します.
Time filter is applied to predictional variables.
subroutine TimeFilterSurfVars( xy_SoilMoistB, xy_SurfSnowB, xy_SoilMoistN, xy_SurfSnowN, xy_SoilMoistA, xy_SurfSnowA ) ! ! 予報変数にタイムフィルターを適用します. ! ! Time filter is applied to predictional variables. ! ! モジュール引用 ; USE statements ! ! 時刻管理 ! Time control ! use timeset, only: TimesetClockStart, TimesetClockStop ! 宣言文 ; Declaration statements ! implicit none real(DP), intent(in ):: xy_SoilMoistB (0:imax-1, 1:jmax) ! $ M_ws (t-\Delta t) $ . 地表水質量. ! Soil moisture real(DP), intent(in ):: xy_SurfSnowB (0:imax-1, 1:jmax) ! $ M_ss (t-\Delta t) $ . 積雪量. ! Surface snow amount real(DP), intent(inout):: xy_SoilMoistN (0:imax-1, 1:jmax) ! $ M_ws (t-\Delta t) $ . 地表水質量. ! Soil moisture real(DP), intent(inout):: xy_SurfSnowN (0:imax-1, 1:jmax) ! $ M_ss (t) $ . 積雪量. ! Surface snow amount real(DP), intent(in ):: xy_SoilMoistA (0:imax-1, 1:jmax) ! $ M_ws (t-\Delta t) $ . 地表水質量. ! Soil moisture real(DP), intent(in ):: xy_SurfSnowA (0:imax-1, 1:jmax) ! $ M_ss (t+\Delta t) $ . 積雪量. ! Surface snow amount ! 作業変数 ! Work variables ! !!$ integer:: i ! 経度方向に回る DO ループ用作業変数 !!$ ! Work variables for DO loop in longitude !!$ integer:: j ! 緯度方向に回る DO ループ用作業変数 !!$ ! Work variables for DO loop in latitude ! 実行文 ; Executable statement ! ! 計算時間計測開始 ! Start measurement of computation time ! call TimesetClockStart( module_name ) ! 初期化 ! Initialization ! if ( .not. timefilter_asselin1972_inited ) call TimeFiltInit !!$ ! ステップのチェック !!$ ! Check Step !!$ ! !!$ CurStep = CurStep + 1 !!$ if ( CurStep - PrevStep < IntStep ) return !!$ PrevStep = CurStep ! FilterParameter の有効性のチェック ! Check validity of "FilterParameter" ! if ( FilterParameter == 0.0_DP ) return ! タイムフィルターを適用 ! Apply time filter ! xy_SoilMoistN = xy_SoilMoistN + FilterParameter * ( xy_SoilMoistB - 2.0d0 * xy_SoilMoistN + xy_SoilMoistA ) xy_SurfSnowN = xy_SurfSnowN + FilterParameter * ( xy_SurfSnowB - 2.0d0 * xy_SurfSnowN + xy_SurfSnowA ) ! 計算時間計測一時停止 ! Pause measurement of computation time ! call TimesetClockStop( module_name ) end subroutine TimeFilterSurfVars
Variable : | |||
timefilter_asselin1972_inited = .false. : | logical, save, public
|
Variable : | |||||||||
FilterParameter : | real(DP), save
|
Subroutine : |
依存モジュールの初期化チェック
Check initialization of dependency modules
subroutine InitCheck ! ! 依存モジュールの初期化チェック ! ! Check initialization of dependency modules ! モジュール引用 ; USE statements ! ! NAMELIST ファイル入力に関するユーティリティ ! Utilities for NAMELIST file input ! use namelist_util, only: namelist_util_inited ! 格子点設定 ! Grid points settings ! use gridset, only: gridset_inited ! 物理定数設定 ! Physical constants settings ! use constants, only: constants_inited ! 座標データ設定 ! Axes data settings ! use axesset, only: axesset_inited ! 時刻管理 ! Time control ! use timeset, only: timeset_inited ! 実行文 ; Executable statement ! if ( .not. namelist_util_inited ) call MessageNotify( 'E', module_name, '"namelist_util" module is not initialized.' ) if ( .not. gridset_inited ) call MessageNotify( 'E', module_name, '"gridset" module is not initialized.' ) if ( .not. constants_inited ) call MessageNotify( 'E', module_name, '"constants" module is not initialized.' ) if ( .not. axesset_inited ) call MessageNotify( 'E', module_name, '"axesset" module is not initialized.' ) if ( .not. timeset_inited ) call MessageNotify( 'E', module_name, '"timeset" module is not initialized.' ) end subroutine InitCheck
Subroutine : |
timefilter_asselin1972 モジュールの初期化を行います. NAMELIST#timefilter_asselin1972_nml の読み込みはこの手続きで行われます.
"timefilter_asselin1972" module is initialized. "NAMELIST#timefilter_asselin1972_nml" is loaded in this procedure.
This procedure input/output NAMELIST#timefilter_asselin1972_nml .
subroutine TimeFiltInit ! ! timefilter_asselin1972 モジュールの初期化を行います. ! NAMELIST#timefilter_asselin1972_nml の読み込みはこの手続きで行われます. ! ! "timefilter_asselin1972" module is initialized. ! "NAMELIST#timefilter_asselin1972_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 /timefilter_asselin1972_nml/ FilterParameter !, IntStep ! ! デフォルト値については初期化手続 "timefilter_asselin1972#Hs94Init" ! のソースコードを参照のこと. ! ! Refer to source codes in the initialization procedure ! "timefilter_asselin1972#Hs94Init" for the default values. ! ! 実行文 ; Executable statement ! if ( timefilter_asselin1972_inited ) return call InitCheck ! デフォルト値の設定 ! Default values settings ! FilterParameter = 0.05_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 = timefilter_asselin1972_nml, iostat = iostat_nml ) ! (out) close( unit_nml ) call NmlutilMsg( iostat_nml, module_name ) ! (in) if ( iostat_nml == 0 ) write( STDOUT, nml = timefilter_asselin1972_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) ) timefilter_asselin1972_inited = .true. end subroutine TimeFiltInit
Constant : | |||
module_name = ‘timefilter_asselin1972‘ : | character(*), parameter
|
Constant : | |||
version = ’$Name: dcpam5-20100424 $’ // ’$Id: timefilter_asselin1972.f90,v 1.8 2009-10-19 09:17:12 yot Exp $’ : | character(*), parameter
|