Class timeset
In: setup/timeset.f90

時刻管理

Time control

Note that Japanese and English are described in parallel.

時刻情報の管理を行います.

Information of time is controlled.

Variables List

TimeB :ステップ $ t - Delta t $ の時刻
TimeN :ステップ $ t $ の時刻
TimeA :ステップ $ t + Delta t $ の時刻
DelTime :$ Delta t $ [s] (実数型の数値)
RestartTime :リスタート開始時刻
EndTime :計算終了時刻
InitialDate :計算開始日時
———— :————
TimeB :Time of step $ t - Delta t $
TimeN :Time of step $ t $
TimeA :Time of step $ t + Delta t $
DelTime :$ Delta t $ [s] (numerical value of float type)
RestartTime :Restart time of calculation
EndTime :End time of calculation
InitialDate :Start date of calculation

Procedures List

TimesetInit :timeset モジュールの初期化
TimesetDelTimeHalf :Δt を一時的に半分に
TimesetProgress :時刻の進行
TimesetSetTimeN :TimeN の設定
TimesetSetInitialDate :計算開始日時の設定
TimesetClockStart :計算時間計測開始
TimesetClockStop :計算時間計一時停止
TimesetClose :timeset モジュールの終了処理
———— :————
TimesetInit :Initialize "timeset" module
TimesetDelTimeHalf :Reduce delta t to half temporarily
TimesetProgress :Progress time
TimesetSetTimeN :Set TimeN
TimesetSetInitialDate :Set start date
TimesetClockStart :Start measurement of computation time
TimesetClockStop :Pause measurement of computation time
TimesetClose :Terminate "timeset" module

NAMELIST

NAMELIST#timeset_nml

Methods

Included Modules

dc_types dc_message dc_clock dc_calendar namelist_util dc_iounit dc_string

Public Instance methods

DelTime
Variable :
DelTime :real(DP), save, public
: $ Delta t $ [s]
EndDate
Variable :
EndDate :type(DC_CAL_DATE), save, public
: 計算終了の日時
EndTime
Variable :
EndTime :real(DP), save, public
: 計算終了時刻. End time of calculation
InitialDate
Variable :
InitialDate :type(DC_CAL_DATE), save, public
: 計算開始日時. Start date of calculation
RestartTime
Variable :
RestartTime :real(DP), save, public
: リスタート開始時刻. Restart time of calculation
TimeA
Variable :
TimeA :real(DP), save, public
: ステップ $ t + Delta t $ の時刻. Time of step $ t + Delta t $.
TimeB
Variable :
TimeB :real(DP), save, public
: ステップ $ t - Delta t $ の時刻. Time of step $ t - Delta t $.
TimeN
Variable :
TimeN :real(DP), save, public
: ステップ $ t $ の時刻. Time of step $ t $.
Subroutine :
name :character(*), intent(in)
: モジュールの名称. Name of module

プログラム単位 (主にモジュールを想定) ごとの時間計測を開始します.

Start measurement of computation time by program unit (expected modules).

[Source]

  subroutine TimesetClockStart( name )
    !
    ! プログラム単位 (主にモジュールを想定) ごとの時間計測を開始します. 
    !
    ! Start measurement of computation time by program unit
    ! (expected modules). 

    ! モジュール引用 ; USE statements
    !

    ! CPU 時間計測
    ! CPU time monitor
    !
    use dc_clock, only: DCClockCreate, DCClockStart

    ! 宣言文 ; Declaration statements
    !
    implicit none
    character(*), intent(in):: name
                              ! モジュールの名称. 
                              ! Name of module

    ! 作業変数
    ! Work variables
    !
    integer:: i               ! clocks, clocks_name 用 DO ループ用作業変数
                              ! Work variables for DO loop for "clocks", "clocks_name"

    ! 実行文 ; Executable statement
    !

    if ( .not. CpuTimeMoniter ) return

    do i = 1, clk_proc_num
      if ( trim(clocks_name(i)) == trim(name) ) then
        call DCClockStart( clocks(i) ) ! (in)
        return
      end if
    end do

    call DCClockCreate( clocks(clk_proc_num + 1), name ) ! (in)
    call DCClockStart( clocks(clk_proc_num + 1) ) ! (in)
    clocks_name(clk_proc_num + 1) = name
    clk_proc_num = clk_proc_num + 1

  end subroutine TimesetClockStart
Subroutine :
name :character(*), intent(in)
: モジュールの名称. Name of module

プログラム単位 (主にモジュールを想定) ごとの時間計測を一時停止します.

Pause measurement of computation time by program unit (expected modules).

[Source]

  subroutine TimesetClockStop( name )
    !
    ! プログラム単位 (主にモジュールを想定) ごとの時間計測を一時停止します.
    !
    ! Pause measurement of computation time by program unit
    ! (expected modules). 

    ! モジュール引用 ; USE statements
    !

    ! CPU 時間計測
    ! CPU time monitor
    !
    use dc_clock, only: DCClockStop

    ! 宣言文 ; Declaration statements
    !
    implicit none
    character(*), intent(in):: name
                              ! モジュールの名称. 
                              ! Name of module

    ! 作業変数
    ! Work variables
    !
    integer:: i               ! clocks, clocks_name 用 DO ループ用作業変数
                              ! Work variables for DO loop for "clocks", "clocks_name"

    ! 実行文 ; Executable statement
    !

    if ( .not. CpuTimeMoniter ) return

    do i = 1, clk_proc_num
      if ( trim(clocks_name(i)) == trim(name) ) then
        call DCClockStop( clocks(i) ) ! (in)
        return
      end if
    end do

    call MessageNotify( 'W', module_name, ' name "%c" is not found in "TimesetClockStop"', c1 = trim(name) )

  end subroutine TimesetClockStop
Subroutine :

計算時間の総計を表示します.

Total computation time is printed.

[Source]

  subroutine TimesetClose
    !
    ! 計算時間の総計を表示します. 
    !
    ! Total computation time is printed. 

    ! モジュール引用 ; USE statements
    !

    ! CPU 時間計測
    ! CPU time monitor
    !
    use dc_clock, only: DCClockStop, DCClockResult, DCClockSetName, operator(+), operator(-)

    ! 宣言文 ; Declaration statements
    !
    implicit none
    integer:: i               ! clocks, clocks_name 用 DO ループ用作業変数
                              ! Work variables for DO loop for "clocks", "clocks_name"

    ! 作業変数
    ! Work variables
    !

    ! 実行文 ; Executable statement
    !

    if ( .not. CpuTimeMoniter ) return

    ! CPU 時間計測終了 (モデル全体)
    ! Stop CPU time monitoring (for entire model)
    !
    call DCClockStop( clocks(1) ) ! (in)

    ! 「その他」の CPU 時間を算出
    ! Calculate CPU time of "Others"
    !
    clocks(clk_proc_num + 1) = clocks(1)
    do i = 2, clk_proc_num
      clocks(clk_proc_num + 1) = clocks(clk_proc_num + 1) - clocks(i)
    end do
    call DCClockSetName( clocks(clk_proc_num + 1), 'others' )

    ! CPU 時間の総計を表示
    ! Print total CPU time
    !
    call DCClockResult( clocks(2:clk_proc_num + 1), total_auto = .true. ) ! (in)

  end subroutine TimesetClose
Subroutine :

計算の初回だけはオイラー法を用いるため, 一時的に Δt を半分にします. TimesetProgress が呼ばれた段階で Δt は元に戻ります.

Delta t is reduced to half temporarily in order to use Euler method at initial step. Delta t is returned to default, when "TimesetProgress" is called.

[Source]

  subroutine TimesetDelTimeHalf
    ! 
    ! 計算の初回だけはオイラー法を用いるため, 
    ! 一時的に Δt を半分にします. 
    ! TimesetProgress が呼ばれた段階で Δt は元に戻ります. 
    ! 
    ! Delta t is reduced to half temporarily 
    ! in order to use Euler method at initial step. 
    ! Delta t is returned to default, when "TimesetProgress" is called. 
    !

    ! 宣言文 ; Declaration statements
    !
    implicit none

    ! 作業変数
    ! Work variables
    !

    ! 実行文 ; Executable statement
    !
    if ( flag_half ) return
    DelTime = DelTime / 2.0_DP
    flag_half = .true.

  end subroutine TimesetDelTimeHalf
Subroutine :

timeset モジュールの初期化を行います. NAMELIST#timeset_nml の読み込みはこの手続きで行われます.

"timeset" module is initialized. NAMELIST#timeset_nml is loaded in this procedure.

This procedure input/output NAMELIST#timeset_nml .

[Source]

  subroutine TimesetInit
    !
    ! timeset モジュールの初期化を行います. 
    ! NAMELIST#timeset_nml の読み込みはこの手続きで行われます. 
    !
    ! "timeset" module is initialized. 
    ! NAMELIST#timeset_nml is loaded in this procedure. 
    !

    ! モジュール引用 ; USE statements
    !

    ! NAMELIST ファイル入力に関するユーティリティ
    ! Utilities for NAMELIST file input
    !
    use namelist_util, only: namelist_filename, NmlutilMsg

    ! ファイル入出力補助
    ! File I/O support
    !
    use dc_iounit, only: FileOpen

    ! 種別型パラメタ
    ! Kind type parameter
    !
    use dc_types, only: STDOUT ! 標準出力の装置番号. Unit number of standard output

    ! 暦と日時の取り扱い
    ! Calendar and Date handler
    !
    use dc_calendar, only: DCCalCreate, DCCalConvertByUnit, DCCalToChar, DCCalInquire, DCCalDateCreate, DCCalDateDifference, DCCalDateToChar, DCCalDateInquire, DCCalDefault

    ! CPU 時間計測
    ! CPU time monitor
    !
    use dc_clock, only: DCClockCreate, DCClockStart

    ! 文字列操作
    ! Character handling
    !
    use dc_string, only: toChar

    ! 宣言文 ; Declaration statements
    !
    implicit none

    integer:: InitialYear, InitialMonth, InitialDay, InitialHour, InitialMin
                              ! 計算開始日時 (年月日時分)
                              ! Start date of calculation (year, month, day, hour, minute)
    real(DP):: InitialSec
                              ! 計算開始日時 (秒)
                              ! Start date of calculation (second)
    real(DP):: IntegPeriodValue
                              ! 積分時間. 
                              ! Integral time. 
    character(TOKEN):: IntegPeriodUnit
                              ! "IntegPeriodValue" の単位. 
                              ! Unit of "IntegPeriodValue". 
    integer:: EndYear, EndMonth, EndDay, EndHour, EndMin
                              ! 計算終了日時 (年月日時分). 
                              ! "IntegPeriodValue" が負の場合にこちらが使用される. 
                              ! 
                              ! End date of calculation (year, month, day, hour, minute)
                              ! These are used when "IntegPeriodValue" is negative
    real(DP):: EndSec
                              ! 計算終了日時 (秒). 
                              ! "IntegPeriodValue" が負の場合にこちらが使用される. 
                              ! 
                              ! End date of calculation (second)
                              ! These are used when "IntegPeriodValue" is negative

    ! 印字用作業変数
    ! Work variables for print
    !
    character(TOKEN):: cal_type_print
    type(DC_CAL):: cal_print
    real(DP):: EndTimeValue_print
    character(TOKEN):: date_print
    
    ! 作業変数
    ! Work variables
    !
    integer:: unit_nml        ! NAMELIST ファイルオープン用装置番号. 
                              ! Unit number for NAMELIST file open
    integer:: iostat_nml      ! NAMELIST 読み込み時の IOSTAT. 
                              ! IOSTAT of NAMELIST read

    ! NAMELIST 変数群
    ! NAMELIST group name
    !
    namelist /timeset_nml/ DelTimeValue,    DelTimeUnit, RestartTimeValue,  RestartTimeUnit, InitialYear, InitialMonth, InitialDay, InitialHour, InitialMin, InitialSec, EndYear, EndMonth, EndDay, EndHour, EndMin, EndSec, IntegPeriodValue, IntegPeriodUnit, PredictIntValue,  PredictIntUnit, CpuTimeMoniter
          !
          ! デフォルト値については初期化手続 "timeset#TimesetInit" 
          ! のソースコードを参照のこと. 
          !
          ! Refer to source codes in the initialization procedure
          ! "timeset#TimesetInit" for the default values. 
          !

    ! 実行文 ; Executable statement
    !

    if ( timeset_inited ) return
    call InitCheck

    ! 暦の情報を設定 (サンプルとしていくつかの case を記載)
    !

    ! case 1:
    ! dc_calendar で用意される典型的な暦をデフォルトの暦として設定
    !
    call DCCalCreate( cal_type = 'gregorian' ) ! (in)  ! グレゴリオ暦
!!$    call DCCalCreate( cal_type = 'julian'    ) ! (in)  ! ユリウス暦
!!$    call DCCalCreate( cal_type = 'noleap'    ) ! (in)  ! 閏年無しの暦
!!$    call DCCalCreate( cal_type = '360day'    ) ! (in)  ! 1ヶ月が 30 日の暦
!!$    call DCCalCreate( cal_type = 'cyclic'    ) ! (in)  ! ある月の日数を 「30.6 × 月数 − 前月までの総日数」の小数点以下切捨とする暦 

!!$    ! case 2:
!!$    ! 決まり切った暦 (以下はユリウス暦) を設定
!!$    !
!!$    call DCCalCreate( cal_type = 'Julian' ) ! (in)
!!$
!!$    ! case 3:
!!$    ! 地球っぽい暦を利用者が明示的に指定
!!$    !
!!$    call DCCalCreate( month_in_year = 12 , &   ! (in)
!!$      &               day_in_month  = (/31, 28, 31, 30, 31, 30,   &
!!$      &                                 31, 31, 30, 31, 30, 31/), & ! (in) 
!!$      &               hour_in_day = 24,     &  ! (in) 
!!$      &               min_in_hour = 60 ,    &  ! (in) 
!!$      &               sec_in_min  = 60.0d0 )   ! (in) 
!!$
!!$    ! case 3:
!!$    ! 火星っぽい暦を利用者が明示的に指定. 
!!$    ! ※ 「暦に月がない」場合には, month_in_year = 1 とする. 
!!$    !
!!$    call DCCalCreate( month_in_year = 1 , &      ! (in)
!!$      &               day_in_month  = (/669/), & ! (in)
!!$      &               hour_in_day = 24, &        ! (in)
!!$      &               min_in_hour = 1 , &        ! (in)
!!$      &               sec_in_min  = 3694.0d0  )  ! (in)

    ! デフォルト値の設定
    ! Default values settings
    !
    DelTimeValue = 30.0_DP
    DelTimeUnit = 'min'
    DelTime = DCCalConvertByUnit( DelTimeValue, DelTimeUnit, 'sec' ) ! (in)
    DelTimeSave = DelTime
    flag_half = .false.

    InitialYear  = 2000
    InitialMonth = 1
    InitialDay   = 1
    InitialHour  = 0
    InitialMin   = 0
    InitialSec   = 0.0

!!$    InitialYear  = - 1 
!!$    InitialMonth = - 1
!!$    InitialDay   = - 1
!!$    InitialHour  = - 1
!!$    InitialMin   = - 1
!!$    InitialSec   = - 1.0

    IntegPeriodValue = -1.0  ! Negative value is invalid (Following EndXXX is used)
    IntegPeriodUnit  = 'sec'

    EndYear    = 2000
    EndMonth   = 2
    EndDay     = 1
    EndHour    = 0
    EndMin     = 0
    EndSec     = 0.0

!!$    EndYear    = - 1  
!!$    EndMonth   = - 1  
!!$    EndDay     = - 1  
!!$    EndHour    = - 1  
!!$    EndMin     = - 1  
!!$    EndSec     = - 1.0

    PredictIntValue = 1.0_DP
    PredictIntUnit = 'day'

    PredictIntTime = DCCalConvertByUnit(PredictIntValue, PredictIntUnit, 'sec' )  ! (in)

    CpuTimeMoniter = .true.


    ! 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 = timeset_nml, iostat = iostat_nml ) ! (out)
      close( unit_nml )

      call NmlutilMsg( iostat_nml, module_name ) ! (in)
      if ( iostat_nml == 0 ) write( STDOUT, nml = timeset_nml )
    end if

    ! 時間変数 (秒) の設定
    ! Calculate time variable (sec)
    !
    DelTime   = DCCalConvertByUnit( DelTimeValue,    DelTimeUnit,  'sec' )  ! (in)
    RestartTime = DCCalConvertByUnit( RestartTimeValue, RestartTimeUnit, 'sec' )  ! (in)
    PredictIntTime = DCCalConvertByUnit( PredictIntValue, PredictIntUnit, 'sec' )  ! (in)

    ! 前回の終了予測日時表示時間の初期設定
    ! Initialize Time when predicted end time is output previously
    !
    PredictPrevTime = RestartTime - DelTime

    ! $ \Delta t $ [s] 算出
    ! Calculate $ \Delta t $ [s] 
    !
    DelTime = DCCalConvertByUnit( DelTimeValue, DelTimeUnit, 'sec' ) ! (in)
    DelTimeSave = DelTime

    ! 開始日時の情報を設定
    !
    ! 再計算時には restart_file_io モジュールによる TimesetSetInitialDate 
    ! 呼び出しによって上書きされることを想定. 
    !
    if ( .not. InitialYear  < 0 .and. .not. InitialMonth < 0 .and. .not. InitialDay   < 0 .and. .not. InitialHour  < 0 .and. .not. InitialMin   < 0 .and. .not. InitialSec   < 0.0 ) then

      call DCCalDateCreate( year  = InitialYear, month = InitialMonth, day   = InitialDay, hour  = InitialHour, min   = InitialMin, sec   = InitialSec, date  = InitialDate )  ! (out) optional
    end if

    ! 終了日時 (開始日時からの経過時間で表現) を設定
    ! 
    if ( IntegPeriodValue < 0.0 ) then

      IntegPeriod = IntegPeriodValue

      ! case 1
      ! 終了日時を使う場合
      !
      call DCCalDateCreate( year  = EndYear, month = EndMonth, day   = EndDay, hour  = EndHour, min   = EndMin, sec   = EndSec, date  = EndDate )    ! (out) optional

      EndTime = DCCalDateDifference( start_date = InitialDate, end_date   = EndDate )       ! (in)

      if ( .not. EndTime > 0.0 ) then
        call DCCalDateInquire( date_print, date = InitialDate )
        call MessageNotify('W', module_name, 'InitialDate=<%c>', c1 = trim(date_print) )
        call DCCalDateInquire( date_print, date = EndDate )
        call MessageNotify('W', module_name, 'EndDate=<%c>', c1 = trim(date_print) )
        call MessageNotify('E', module_name, '"EndTime" must be positive.')
      end if

    else
      ! case 2
      ! 積分時間を使う場合
      !
      IntegPeriod = DCCalConvertByUnit( IntegPeriodValue, IntegPeriodUnit, 'sec' )  ! (in)

      EndTime = RestartTime + IntegPeriod
    end if

    ! 時刻の正当性のチェック
    ! Check validation of time
    !
    call TimeValidCheck( RestartTime, EndTime, DelTime, PredictIntTime ) ! (in)

    ! 各タイムステップの時刻の設定
    ! Configure time on each time step
    !
    ! 再計算時には restart_file_io モジュールによる TimesetSetTimeN 
    ! 呼び出しによって上書きされることを想定. 
    ! 
    TimeN = RestartTime
    TimeB = TimeN - DelTime
    TimeA = TimeN + DelTime


    ! CPU 時間計測開始 (モデル全体)
    ! Start CPU time monitoring (for entire model)
    !
    if ( CpuTimeMoniter ) then
      call DCClockCreate( clocks(clk_proc_num + 1), 'total' ) ! (in)
      call DCClockStart( clocks(clk_proc_num + 1) ) ! (in)
      clocks_name(clk_proc_num + 1) = 'total'
      clk_proc_num = clk_proc_num + 1
    end if

    ! 印字 ; Print
    !
    call MessageNotify( 'M', module_name, '----- Initialization Messages -----' )
    call MessageNotify( 'M', module_name, '  RestartTime  = %f [%c]', d = (/ RestartTimeValue /), c1 = trim(RestartTimeUnit) )
    call DCCalDateInquire( date_print, date = InitialDate )
    call MessageNotify( 'M', module_name, '  InitialDate  = %c', c1 = trim(date_print) )
    if ( IntegPeriodValue < 0.0 ) then
      call DCCalDateInquire( date_print, date = EndDate )
      call MessageNotify( 'M', module_name, '  EndDate    = %c', c1 = trim(date_print) )
    end if
    EndTimeValue_print = DCCalConvertByUnit( EndTime, 'sec', RestartTimeUnit ) ! (in)
    call MessageNotify( 'M', module_name, '  EndTime    = %f [%c]', d = (/ EndTimeValue_print /), c1 = trim(RestartTimeUnit) )
    call MessageNotify( 'M', module_name, '  DelTime    = %f [%c]', d = (/ DelTimeValue /), c1 = trim(DelTimeUnit) )
    call MessageNotify( 'M', module_name, '             = %f [%c]', d = (/ DelTime /), c1 = 'sec' )
    call DCCalInquire( cal_type = cal_type_print ) ! (out) optional
    if ( cal_type_print /= 'user_defined' ) then
      call MessageNotify( 'M', module_name, '  Calendar   = %c', c1 = trim(cal_type_print) )
    else
      call DCCalDefault( cal_print ) ! (out) 
      call MessageNotify( 'M', module_name, '  Calendar   = %c', c1 = trim(DCCalToChar(cal_print)) )
    end if
    call MessageNotify( 'M', module_name, '  PredictInt = %f [%c]', d = (/ PredictIntValue /), c1 = trim(PredictIntUnit) )
    call MessageNotify( 'M', module_name, '  CpuTimeMoniter = %b', l = (/ CpuTimeMoniter /) )
    call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) )

    timeset_inited = .true.
  end subroutine TimesetInit
Subroutine :

timeset モジュール内の時刻を進めます. また, TimesetProgress#PredictIntStep で設定された値に応じて, 現在までの計算時間と計算終了予測時刻を表示します.

Progress time configured in "timeset" module. And, computation time until now and predicted end of computation time are printed according to configured "TimesetProgress#PredictIntStep"

[Source]

  subroutine TimesetProgress
    !
    ! timeset モジュール内の時刻を進めます. 
    ! また, TimesetProgress#PredictIntStep で設定された値に応じて, 
    ! 現在までの計算時間と計算終了予測時刻を表示します. 
    !
    ! Progress time configured in "timeset" module. 
    ! And, computation time until now and 
    ! predicted end of computation time are printed 
    ! according to configured "TimesetProgress#PredictIntStep"
    !

    ! モジュール引用 ; USE statements
    !

    ! CPU 時間計測
    ! CPU time monitor
    !
    use dc_clock, only: DCClockPredict, DCClockStop, DCClockClose, operator(+), operator(-)

    ! 宣言文 ; Declaration statements
    !
    implicit none

    ! 作業変数
    ! Work variables
    !
    type(CLOCK):: clock_tmp

    ! 実行文 ; Executable statement
    !

    ! Δt を元に戻す. 
    ! Delta t is returned to default
    !
    if ( flag_half ) then
      DelTime = DelTimeSave
      flag_half = .false.
    end if

    ! 終了予測日時表示
    ! Print predicted end time
    !
    if ( .not. TimeA - PredictPrevTime < PredictIntTime ) then
      PredictPrevTime = PredictPrevTime + PredictIntTime
      if ( CpuTimeMoniter ) then
        clock_tmp = clocks(1)
        call DCClockStop( clock_tmp ) ! (in)
        call DCClockPredict( clock_tmp, real( ( TimeA - RestartTime ) / ( ( EndTime + DelTimeSave ) - RestartTime ) ) ) ! (in)
        call DCClockClose( clock_tmp ) ! (in)
      end if
    end if

    ! 時刻の進行
    ! Progress time
    !
    TimeB = TimeB + DelTime
    TimeN = TimeN + DelTime
    TimeA = TimeA + DelTime
  end subroutine TimesetProgress
Subroutine :
cal_type :character(*), intent(in)
: 暦のタイプ.
month_in_year :integer, intent(in)
day_in_month(:) :integer, intent(in)
hour_in_day :integer, intent(in)
min_in_hour :integer, intent(in)
sec_in_min :real(DP), intent(in)
: 暦情報. Information of Calendar.

(in)

暦の再設定を行います.

TimesetInit が既に呼ばれることが前提です. TimesetInit が呼ばれる前にこのサブルーチンが呼ばれた場合, 何もせずにこのサブルーチンは終了します.

Calendar is reconfigured.

"TimesetInit" must be called before this subroutine is called. If "TimesetInit" is not called previously, this subroutine is finished with no changes.

[Source]

  subroutine TimesetSetCalendar( cal_type, month_in_year, day_in_month, hour_in_day, min_in_hour, sec_in_min )       ! (in)
    !
    ! 暦の再設定を行います. 
    !
    ! TimesetInit が既に呼ばれることが前提です. 
    ! TimesetInit が呼ばれる前にこのサブルーチンが呼ばれた場合,
    ! 何もせずにこのサブルーチンは終了します. 
    !
    ! Calendar is reconfigured. 
    !
    ! "TimesetInit" must be called before this subroutine is called. 
    ! If "TimesetInit" is not called previously, this subroutine 
    ! is finished with no changes. 
    !

    ! モジュール引用 ; USE statements
    !

    ! 暦と日時の取り扱い
    ! Calendar and Date handler
    !
    use dc_calendar, only: DCCalCreate, DCCalInquire, DCCalDefault, DCCalToChar

    ! 文字列操作
    ! Character handling
    !
    use dc_string, only: toChar

    ! 宣言文 ; Declaration statements
    !
    implicit none

    character(*), intent(in):: cal_type
                              ! 暦のタイプ. 
    integer, intent(in):: month_in_year, day_in_month(:), hour_in_day, min_in_hour
    real(DP), intent(in):: sec_in_min
                              ! 暦情報. 
                              ! Information of Calendar. 

    ! 作業変数
    ! Work variables
    !
    character(TOKEN):: cal_type_print
    type(DC_CAL):: cal_print

    ! 実行文 ; Executable statement
    !

    if ( .not. timeset_inited ) return

    ! 開始日時の情報を設定
    !
    if ( cal_type /= 'user_defined' ) then
      call DCCalCreate( cal_type = cal_type )  ! (in)
    else
      call DCCalCreate( month_in_year, day_in_month , hour_in_day, min_in_hour , sec_in_min )     ! (in) 
    end if

    ! 印字 ; Print
    !
    call MessageNotify( 'M', module_name, '----- Reconfigure Calendar Messages -----' )
    call DCCalInquire( cal_type = cal_type_print ) ! (out) optional
    if ( cal_type_print /= 'user_defined' ) then
      call MessageNotify( 'M', module_name, '  Calendar   = %c', c1 = trim(cal_type_print) )
    else
      call DCCalDefault( cal_print ) ! (out) 
      call MessageNotify( 'M', module_name, '  Calendar   = %c', c1 = trim(DCCalToChar(cal_print)) )
    end if

  end subroutine TimesetSetCalendar
Subroutine :
InitialYear :integer, intent(in)
: 計算開始の年月日時分
InitialMonth :integer, intent(in)
: 計算開始の年月日時分
InitialDay :integer, intent(in)
: 計算開始の年月日時分
InitialHour :integer, intent(in)
: 計算開始の年月日時分
InitialMin :integer, intent(in)
: 計算開始の年月日時分
InitialSec :real(DP), intent(in)
: 計算開始の秒

開始日時の再設定を行います.

TimesetInit が既に呼ばれることが前提です. TimesetInit が呼ばれる前にこのサブルーチンが呼ばれた場合, 何もせずにこのサブルーチンは終了します.

Start date is reconfigured.

"TimesetInit" must be called before this subroutine is called. If "TimesetInit" is not called previously, this subroutine is finished with no changes.

[Source]

  subroutine TimesetSetInitialDate( InitialYear, InitialMonth, InitialDay, InitialHour, InitialMin, InitialSec )
    !
    ! 開始日時の再設定を行います. 
    !
    ! TimesetInit が既に呼ばれることが前提です. 
    ! TimesetInit が呼ばれる前にこのサブルーチンが呼ばれた場合,
    ! 何もせずにこのサブルーチンは終了します. 
    !
    ! Start date is reconfigured. 
    !
    ! "TimesetInit" must be called before this subroutine is called. 
    ! If "TimesetInit" is not called previously, this subroutine 
    ! is finished with no changes. 
    !

    ! モジュール引用 ; USE statements
    !

    ! 暦と日時の取り扱い
    ! Calendar and Date handler
    !
    use dc_calendar, only: DCCalCreate, DCCalConvertByUnit, DCCalToChar, DCCalInquire, DCCalDateCreate, DCCalDateDifference, DCCalDateToChar, DCCalDateInquire

    ! 文字列操作
    ! Character handling
    !
    use dc_string, only: toChar

    ! 宣言文 ; Declaration statements
    !
    implicit none

    integer, intent(in):: InitialYear, InitialMonth, InitialDay, InitialHour, InitialMin
                              ! 計算開始の年月日時分
    real(DP), intent(in):: InitialSec
                              ! 計算開始の秒

    ! 作業変数
    ! Work variables
    !
    real(DP):: EndTimeValue_print
    character(TOKEN):: date_print

    ! 実行文 ; Executable statement
    !

    if ( .not. timeset_inited ) return

    ! "InitialDate" を再設定
    ! Reconfigure "InitialDate"
    !
    call DCCalDateCreate( year  = InitialYear, month = InitialMonth, day   = InitialDay, hour  = InitialHour, min   = InitialMin, sec   = InitialSec, date  = InitialDate )  ! (out) optional

    ! "EndTime" の再設定
    ! Reconfigure "EndTime"
    !
    if ( IntegPeriod < 0.0 ) then
      EndTime = DCCalDateDifference( start_date = InitialDate, end_date = EndDate )         ! (in)

      if ( .not. EndTime > 0.0 ) then
        call DCCalDateInquire( date_print, date = InitialDate )
        call MessageNotify('W', module_name, 'InitialDate=<%c>', c1 = trim(date_print) )
        call DCCalDateInquire( date_print, date = EndDate )
        call MessageNotify('W', module_name, 'EndDate=<%c>', c1 = trim(date_print) )
        call MessageNotify('E', module_name, '"EndTime" must be positive.')
      end if

    end if

    ! 印字 ; Print
    !
    call MessageNotify( 'M', module_name, '----- Reconfigure InitialDate Messages -----' )
    call DCCalDateInquire( date_print, date = InitialDate )
    call MessageNotify( 'M', module_name, '  InitialDate  = %c', c1 = trim(date_print) )
    if ( IntegPeriod < 0.0 ) then
      call DCCalDateInquire( date_print, date = EndDate )
      call MessageNotify( 'M', module_name, '  EndDate    = %c', c1 = trim(date_print) )

      EndTimeValue_print = DCCalConvertByUnit( EndTime, 'sec', RestartTimeUnit ) ! (in)
      call MessageNotify( 'M', module_name, '  EndTime    = %f [%c]', d = (/ EndTimeValue_print /), c1 = trim(RestartTimeUnit) )
    end if

  end subroutine TimesetSetInitialDate
Subroutine :
TimeNSet :real(DP), intent(in)
: ステップ $ t $ の時刻. Time of step $ t $.

TimeN の再設定を行います. 自動的に TimeB, TimeA, EndTime についても 再設定を行います.

TimesetInit が既に呼ばれることが前提です. TimesetInit が呼ばれる前にこのサブルーチンが呼ばれた場合, 何もせずにこのサブルーチンは終了します.

"TimeN" is reconfigured. "TimeB", "TimeA", "EndTime" are reconfigured automatically.

"TimesetInit" must be called before this subroutine is called. If "TimesetInit" is not called previously, this subroutine is finished with no changes.

[Source]

  subroutine TimesetSetTimeN(TimeNSet)
    !
    ! TimeN の再設定を行います. 
    ! 自動的に TimeB, TimeA, EndTime についても
    ! 再設定を行います. 
    !
    ! TimesetInit が既に呼ばれることが前提です. 
    ! TimesetInit が呼ばれる前にこのサブルーチンが呼ばれた場合,
    ! 何もせずにこのサブルーチンは終了します. 
    !
    ! "TimeN" is reconfigured. 
    ! "TimeB", "TimeA", "EndTime" are reconfigured automatically. 
    !
    ! "TimesetInit" must be called before this subroutine is called. 
    ! If "TimesetInit" is not called previously, this subroutine 
    ! is finished with no changes. 
    !

    ! モジュール引用 ; USE statements
    !

    ! 暦と日時の取り扱い
    ! Calendar and Date handler
    !
    use dc_calendar, only: DCCalConvertByUnit

    ! 文字列操作
    ! Character handling
    !
    use dc_string, only: toChar

    ! 宣言文 ; Declaration statements
    !
    implicit none

    real(DP), intent(in):: TimeNSet
                              ! ステップ $ t $ の時刻. 
                              ! Time of step $ t $. 
    
    ! 作業変数
    ! Work variables
    !
    real(DP):: EndTimeValue_print

    ! 実行文 ; Executable statement
    !

    if ( .not. timeset_inited ) return

    ! TimeN が負の場合、エラーを発生
    ! Cause an error, if TimeN is negative
    !
    if ( TimeN < 0.0 ) then
      call MessageNotify( 'E', module_name, 'TimeN=<%f [sec]> must be positive', d = (/ TimeNSet /) )
    end if

    ! 各タイムステップの時刻の設定
    ! Reconfigure time on each time step
    !
    TimeN = TimeNSet
    TimeB = TimeN - DelTime
    TimeA = TimeN + DelTime

    ! 前回の終了予測日時表示時間の再設定
    ! Reconfigure Time when predicted end time is output previously
    !
    PredictPrevTime = TimeN - DelTime

    ! 印字 ; Print
    !
    call MessageNotify( 'M', module_name, '----- Reconfigure Time Messages -----' )
    call MessageNotify( 'M', module_name, '  TimeB  = %f [sec]', d = (/ TimeB /) )
    call MessageNotify( 'M', module_name, '  TimeN  = %f [sec]', d = (/ TimeN /) )
    call MessageNotify( 'M', module_name, '  TimeA  = %f [sec]', d = (/ TimeA /) )
    if ( IntegPeriod > 0.0 ) then
      EndTimeValue_print = DCCalConvertByUnit( EndTime, 'sec', RestartTimeUnit ) ! (in)
      call MessageNotify( 'M', module_name, '  EndTime    = %f [%c]', d = (/ EndTimeValue_print /), c1 = trim(RestartTimeUnit) )
    end if

  end subroutine TimesetSetTimeN
timeset_inited
Variable :
timeset_inited = .false. :logical, save, public
: 初期設定フラグ. Initialization flag

Private Instance methods

CpuTimeMoniter
Variable :
CpuTimeMoniter :logical, save
: CPU 時間計測のオンオフ On/off of CPU time monitoring
DCdiffDelTime
Variable :
DCdiffDelTime :real(DP), save
: $ Delta t $
Date
Variable :
Date(8) :integer, save
: 計算開始日時. (年月日時分秒, タイムゾーン時分) Start date of calculation. (year, month, day, hour, minute, second, and hour, minute of time zone)
DelTimeSave
Variable :
DelTimeSave :real(DP), save
: $ Delta t $ [s] のデフォルト値. ("TimesetDelTimeHalf" で使用される)

Default value of $ Delta t $ [s]. (for "TimesetDelTimeHalf")

DelTimeUnit
Variable :
DelTimeUnit :character(TOKEN), save
: $ Delta t $ の単位. Unit of $ Delta t $
DelTimeValue
Variable :
DelTimeValue :real(DP), save
: $ Delta t $ . 単位は DelTimeUnit にて指定. Unit is specified by "DelTimeUnit".
Subroutine :

依存モジュールの初期化チェック

Check initialization of dependency modules

[Source]

  subroutine InitCheck
    !
    ! 依存モジュールの初期化チェック
    !
    ! Check initialization of dependency modules

    ! モジュール引用 ; USE statements
    !

    ! NAMELIST ファイル入力に関するユーティリティ
    ! Utilities for NAMELIST file input
    !
    use namelist_util, only: namelist_util_inited

    ! 実行文 ; Executable statement
    !

    if ( .not. namelist_util_inited ) call MessageNotify( 'E', module_name, '"namelist_util" module is not initialized.' )

  end subroutine InitCheck
IntegPeriod
Variable :
IntegPeriod :real(DP), save
: 積分時間. 単位は秒. Unit is second.
PredictIntTime
Variable :
PredictIntTime :real(DP), save
: 終了予測日時表示時間間隔. Interval time of predicted end time output
PredictIntUnit
Variable :
PredictIntUnit :character(TOKEN), save
: 終了予測日時表示間隔 (単位). Unit for interval of predicted end time output
PredictIntValue
Variable :
PredictIntValue :real(DP), save
: 終了予測日時表示間隔. Interval of predicted end time output
PredictPrevTime
Variable :
PredictPrevTime :real(DP), save
: 前回の終了予測日時表示時間. Time when predicted end time is output previously
RestartTimeUnit
Variable :
RestartTimeUnit :character(TOKEN), save
: リスタート開始時刻の単位. Unit of restart time of calculation
RestartTimeValue
Variable :
RestartTimeValue :real(DP), save
: リスタート開始時刻. Restart time of calculation
Subroutine :
RestartTimeChk :real(DP), intent(in)
: リスタート開始時刻. Restart time of calculation
EndTimeChk :real(DP), intent(in)
: 計算終了時刻. End time of calculation
DelTimeChk :real(DP), intent(in)
: $ Delta t $ [s]
PredictIntTimeChk :real(DP), intent(in)
: 終了予測日時表示間隔. Interval of predicted end time output

時刻情報についての有効性をチェックします.

Check validation about infomation time

[Source]

  subroutine TimeValidCheck( RestartTimeChk, EndTimeChk, DelTimeChk, PredictIntTimeChk )
    !
    ! 時刻情報についての有効性をチェックします. 
    !
    ! Check validation about infomation time
    !

    ! モジュール引用 ; USE statements
    !

    ! 暦と日時の取り扱い
    ! Calendar and Date handler
    !
    use dc_calendar, only: DCCalConvertByUnit

    implicit none
    real(DP), intent(in):: RestartTimeChk
                              ! リスタート開始時刻. 
                              ! Restart time of calculation
    real(DP), intent(in):: EndTimeChk
                              ! 計算終了時刻. 
                              ! End time of calculation
    real(DP), intent(in):: DelTimeChk
                              ! $ \Delta t $ [s]
    real(DP), intent(in):: PredictIntTimeChk
                              ! 終了予測日時表示間隔. 
                              ! Interval of predicted end time output

    ! 印字用作業変数
    ! Work variables for print
    !
    real(DP):: RestartTimeValue_print
    real(DP):: EndTimeValue_print
    real(DP):: DelTimeValue_print
    real(DP):: PredictIntValue_print

    ! 作業変数
    ! Work variables
    !

    ! 実行文 ; Executable statement
    !

    if ( .not. 0.0_DP < ( EndTimeChk - RestartTimeChk ) ) then
      RestartTimeValue_print = DCCalConvertByUnit( RestartTimeChk, 'sec', RestartTimeUnit ) ! (in)
      EndTimeValue_print = DCCalConvertByUnit( EndTimeChk, 'sec', RestartTimeUnit ) ! (in)
      
      call MessageNotify( 'E', module_name, 'RestartTime=<%f[%c]> is later than EndTime=<%f[%c]>', d = (/ RestartTimeValue_print, EndTimeValue_print /), c1 = trim(RestartTimeUnit), c2 = trim(RestartTimeUnit) )
    end if

    if ( DelTimeChk > ( EndTimeChk - RestartTimeChk ) ) then
      RestartTimeValue_print = DCCalConvertByUnit( RestartTimeChk, 'sec', RestartTimeUnit ) ! (in)
      EndTimeValue_print = DCCalConvertByUnit( EndTimeChk, 'sec', RestartTimeUnit ) ! (in)
      DelTimeValue_print = DCCalConvertByUnit( DelTimeChk, 'sec', DelTimeUnit ) ! (in)
      call MessageNotify( 'E', module_name, 'DelTime=<%f[%c]> is larger than ' // 'EndTime=<%f[%c]> - RestartTime=<%f[%c]>.', d = (/ DelTimeValue_print, EndTimeValue_print, RestartTimeValue_print /), c1 = trim(DelTimeUnit), c2 = trim(RestartTimeUnit), c3 = trim(RestartTimeUnit) )
    end if

    if ( .not. DelTimeChk > 0.0_DP ) then
      DelTimeValue_print = DCCalConvertByUnit( DelTimeChk, 'sec', DelTimeUnit ) ! (in)
      call MessageNotify( 'E', module_name, 'DelTime=<%f[%c]> must be more than 0.', d = (/ DelTimeValue_print /), c1 = trim(DelTimeUnit) )
    end if

    if ( .not. PredictIntTimeChk > 0.0_DP ) then
      PredictIntValue_print = DCCalConvertByUnit( PredictIntTimeChk, 'sec', PredictIntUnit ) ! (in)
      call MessageNotify( 'E', module_name, 'PredictInt=<%f[%c]> must be more than 0.', d = (/ PredictIntValue_print /), c1 = trim(PredictIntUnit) )
    end if

  end subroutine TimeValidCheck
clk_proc_num
Variable :
clk_proc_num = 0 :integer, save
: CPU 時間計測を行っているプロセスの数. Number of processes monitored CPU time
clkmax
Constant :
clkmax = 64 :integer, parameter
: CPU 時間計測を行うプロセスの最大数. Maximum number of processes monitored CPU time
clocks
Variable :
clocks(1:clkmax) :type(CLOCK), save
: CPU 時間計測用構造体 Derived type for monitoring CPU time
clocks_name
Variable :
clocks_name(1:clkmax) = ’’ :character(TOKEN), save
: CPU 時間計測を行っているプロセスの名称 Names of processes monitored CPU time
flag_half
Variable :
flag_half :logical, save
: TimesetDelTimeHalf によって $ Delta t $ が 半分になっていることを示すフラグ.

Flag that shows $ Delta t $ is reduced to half by "TimesetDelTimeHalf"

module_name
Constant :
module_name = ‘timeset :character(*), parameter
: モジュールの名称. Module name
version
Constant :
version = ’$Name: $’ // ’$Id: timeset.f90,v 1.1.1.1 2010-08-17 05:24:51 takepiro Exp $’ :character(*), parameter
: モジュールのバージョン Module version