!= 鉛直フィルター (Ishiwatari et al., 2002)
!
!= Vertical filter (Ishiwatari et al., 2002)
!
! Authors:: Satoshi Noda
! Version:: $Id: vertical_filter_wt_VFonlyT.f90,v 1.1 2010-05-07 09:58:22 noda Exp $
! Tag Name:: $Name: dcpam5-20110221-2 $
! Copyright:: Copyright (C) GFD Dennou Club, 2010. All rights reserved.
! License:: See COPYRIGHT[link:../../../COPYRIGHT]
!
module vertical_filter
!
!= 鉛直フィルター (Ishiwatari et al., 2002)
!
!= Vertical filter (Ishiwatari et al., 2002)
!
! 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".
!
!== Procedures List
!
! VerticalFilter :: 鉛直フィルターの適用
! ---------------- :: ---------------------
! VerticalFilter :: Apply vertical filter
!
!== NAMELIST
!
! NAMELIST#vertical_filter_nml
!
!== References
!
! * Ishiwatari, M., Takehiro S., Nakajima, K. and Hayashi, Y. -Y. 2002:
! A Numerical Study on Appearance of the Runaway Greenhouse State of a
! Three-Dimensional Gray Atmosphere.
! J. Atmos. Sci., 59, 3223--3238.
!
! モジュール引用 ; USE statements
!
! 格子点設定
! Grid points settings
!
use gridset, only: imax, & ! 経度格子点数.
! Number of grid points in longitude
& jmax, & ! 緯度格子点数.
! Number of grid points in latitude
& kmax ! 鉛直層数.
! Number of vertical level
use composition, only: ncmax ! 微量成分の数
! Number of composition
use timeset, only: DelTime
! 種別型パラメタ
! Kind type parameter
!
use dc_types, only: DP, & ! 倍精度実数型. Double precision.
& STRING ! 文字列. Strings.
! メッセージ出力
! Message output
!
use dc_message, only: MessageNotify
! 宣言文 ; Declaration statements
!
implicit none
private
! 公開手続き
! Public procedure
!
public:: VerticalFilter
! 公開変数
! Public variables
!
logical, save, public:: vertical_filter_inited = .false.
! 初期設定フラグ.
! Initialization flag
! 非公開変数
! Private variables
!
real(DP), save:: FilterParameter
! 鉛直フィルターの係数.
! Vertical filter coefficient
!!$ integer, save:: IntStep ! 鉛直フィルターを適用するステップ間隔.
!!$ ! Step interval of vertical filter application
!!$ integer, save:: CurStep ! VerticalFilter が呼ばれた回数.
!!$ ! Number of times that "VerticalFilter" is called
!!$ integer, save:: PrevStep ! 前回鉛直フィルターを適用したステップ数.
!!$ ! Step number that is filtered previously
character(*), parameter:: module_name = 'vertical_filter'
! モジュールの名称.
! Module name
character(*), parameter:: version = &
& '$Name: dcpam5-20110221-2 $' // &
& '$Id: vertical_filter_wt_VFonlyT.f90,v 1.1 2010-05-07 09:58:22 noda Exp $'
! モジュールのバージョン
! Module version
! INTERFACE 文 ; INTERFACE statements
!
interface VerticalFilter
module procedure VerticalFilter
end interface
contains
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):: xyz_TempRef(0:imax-1, 1:jmax, 1:kmax), &
& Temp_corr(0:imax-1, 1:jmax)
integer :: NLDump, NLNewtonZonal, NLRayleigh
real(DP) :: CoefDump(1:kmax), CoefNewtonZonal(1:kmax), CoefRayleigh(1:kmax)
real(DP) :: TimeScaleDump
!
! testing (鉛直フィルタ (Ishiwatari et al., 2002))
!
! Ishiwatari et al. で使われていた定数 (F_T = 0.1)
NLDump = 7
TimeScaleDump = 10800.0_DP
do k= 1, kmax - NLDump
CoefDump(k) = 0.0_DP
print *, '*** DEBUG ***', k, CoefDump(k)
end do
do k= kmax - NLDump + 1, kmax
CoefDump(k) = DelTime / (TimeScaleDump * ( kmax - k + 1 ))
print *, '*** DEBUG ***', k, CoefDump(k)
end do
NLNewtonZonal = NLDump
CoefNewtonZonal = CoefDump
NLRayleigh = NLDump
CoefRayleigh = CoefDump
! 実行文 ; 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
xyz_TempRef(:,:,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 &
& * ( xyz_TempRef(:,:,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 * ( xyz_TempRef(:,:,k+k1) - xyz_TempA(:,:,k+k1) ) - Temp_corr
end do
end do
!
! testing (ニュートン冷却(帯状方向への平滑化, Ishiwatari et al. 2002))
!
xyz_TempRef = 0.0_DP
! 参照温度の計算
do k = kmax-NLNewtonZonal+1, kmax
! do j=1,jmax
do i=0,imax-1
xyz_TempRef(0,:,k) = xyz_TempRef(0,:,k) + xyz_TempA(i,:,k)
end do
xyz_TempRef(0,:,k) = xyz_TempRef(0,:,k) / imax
do i=1,imax-1
xyz_TempRef(i,:,k) = xyz_TempRef(0,:,k)
end do
! end do
end do
! 温度の補正
do k = kmax-NLNewtonZonal+1, kmax
xyz_TempA(:,:,k) = xyz_TempA(:,:,k) &
& - CoefNewtonZonal(k) * ( xyz_TempA(:,:,k) - xyz_TempRef(:,:,k) )
end do
!
! testing (レイリー摩擦)
!
do k= kmax-NLRayleigh+1,kmax
xyz_UA(:,:,k) = (1.0_DP - CoefRayleigh(k)) * xyz_UA(:,:,k)
xyz_VA(:,:,k) = (1.0_DP - CoefRayleigh(k)) * xyz_VA(:,:,k)
end do
! 計算時間計測一時停止
! Pause measurement of computation time
!
call TimesetClockStop( module_name )
end subroutine VerticalFilter
!-------------------------------------------------------------------
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
call InitCheck
! デフォルト値の設定
! Default values settings
!
FilterParameter = 0.1_DP
!!$ IntStep = 1
!!$ CurStep = 0
! NAMELIST の読み込み
! NAMELIST is input
!
if ( trim(namelist_filename) /= '' ) then
call FileOpen( unit_nml, & ! (out)
& namelist_filename, mode = 'r' ) ! (in)
rewind( unit_nml )
read( unit_nml, & ! (in)
& nml = vertical_filter_nml, & ! (out)
& 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
!-------------------------------------------------------------------
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
end module vertical_filter