New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
daymod.2.F90 on Ticket #1037 – Attachment – NEMO

Ticket #1037: daymod.2.F90

File daymod.2.F90, 17.6 KB (added by poddo, 11 years ago)

more portable

Line 
1MODULE daymod
2   !!======================================================================
3   !!                       ***  MODULE  daymod  ***
4   !! Ocean        :  calendar
5   !!=====================================================================
6   !! History :  OPA  ! 1994-09  (M. Pontaud M. Imbard)  Original code
7   !!                 ! 1997-03  (O. Marti)
8   !!                 ! 1997-05  (G. Madec)
9   !!                 ! 1997-08  (M. Imbard)
10   !!   NEMO     1.0  ! 2003-09  (G. Madec)  F90 + nyear, nmonth, nday
11   !!                 ! 2004-01  (A.M. Treguier) new calculation based on adatrj
12   !!                 ! 2006-08  (G. Madec)  surface module major update
13   !!                 ! 2012-12  (P. Oddo) generalizzation of the algorithm
14   !!                       it is now possible to start at any time of the day
15   !!----------------------------------------------------------------------     
16
17   !!----------------------------------------------------------------------
18   !!   day        : calendar
19   !! 
20   !!           -------------------------------
21   !!           ----------- WARNING -----------
22   !!
23   !!   we suppose that the time step is deviding the number of second of in a day
24   !!             ---> MOD( rday, rdttra(1) ) == 0
25   !!
26   !!           ----------- WARNING -----------
27   !!           -------------------------------
28   !! 
29   !!----------------------------------------------------------------------
30   USE dom_oce         ! ocean space and time domain
31   USE phycst          ! physical constants
32   USE in_out_manager  ! I/O manager
33   USE iom             !
34   USE ioipsl, ONLY :   ymds2ju,ju2ymds   ! for calendar
35   USE prtctl          ! Print control
36   USE restart         !
37   USE trc_oce, ONLY : lk_offline ! offline flag
38   USE timing          ! Timing
39
40   IMPLICIT NONE
41   PRIVATE
42
43   PUBLIC   day        ! called by step.F90
44   PUBLIC   day_init   ! called by istate.F90
45
46   INTEGER ::   nsecd, nsecd05, ndt, ndt05
47
48   !!----------------------------------------------------------------------
49   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
50   !! $Id$
51   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
52   !!----------------------------------------------------------------------
53CONTAINS
54
55   SUBROUTINE day_init
56      !!----------------------------------------------------------------------
57      !!                   ***  ROUTINE day_init  ***
58      !!
59      !! ** Purpose :   Initialization of the calendar values to their values 1 time step before nit000
60      !!                because day will be called at the beginning of step
61      !!
62      !! ** Action  : - nyear        : current year
63      !!              - nmonth       : current month of the year nyear
64      !!              - nday         : current day of the month nmonth
65      !!              - nday_year    : current day of the year nyear
66      !!              - nsec_year    : current time step counted in second since 00h jan 1st of the current year
67      !!              - nsec_month   : current time step counted in second since 00h 1st day of the current month
68      !!              - nsec_day     : current time step counted in second since 00h of the current day
69      !!              - nsec1jan000  : second since Jan. 1st 00h of nit000 year and Jan. 1st 00h of the current year
70      !!              - nmonth_len, nyear_len, nmonth_half, nmonth_end through day_mth
71      !!----------------------------------------------------------------------
72      INTEGER  ::   inbday, idweek
73      REAL(wp) ::   zjul,fracday
74      !!----------------------------------------------------------------------
75      !
76      ! all calendar staff is based on the fact that MOD( rday, rdttra(1) ) == 0
77      IF( MOD( rday     , rdttra(1) ) /= 0. )   CALL ctl_stop( 'the time step must devide the number of second of in a day' )
78      IF( MOD( rday     , 2.        ) /= 0. )   CALL ctl_stop( 'the number of second of in a day must be an even number'    )
79      IF( MOD( rdttra(1), 2.        ) /= 0. )   CALL ctl_stop( 'the time step (in second) must be an even number'           )
80      nsecd   = NINT(rday           )
81      nsecd05 = NINT(0.5 * rday     )
82      ndt     = NINT(      rdttra(1))
83      ndt05   = NINT(0.5 * rdttra(1))
84
85      IF( .NOT. lk_offline ) CALL day_rst( nit000, 'READ' ) 
86
87      ! set the calandar from ndastp (read in restart file and namelist)
88
89      nyear   =   ndastp / 10000
90      nmonth  = ( ndastp - (nyear * 10000) ) / 100
91      nday    =   ndastp - (nyear * 10000) - ( nmonth * 100 ) 
92 
93      ! round adatrj to the to 0.1 sec
94      adatrj= NINT(adatrj*rday*10._wp,i8)/(rday*10._wp)
95
96      nsec_day= INT((adatrj-INT(adatrj,i8))*rday)
97      ! compute the julian of the present day
98      CALL ymds2ju( nyear, nmonth, nday, REAL(nsec_day), fjulday ) 
99
100      ! Compute date (year, month, day, second) 1/2 time-step before the integration start
101      CALL ju2ymds(fjulday-(ndt05/rday),nyear,nmonth,nday,fracday)
102      ! Compute the number of seconds elapsed from the start of the day(1/2step before)
103      nsec_day=INT(fracday)
104
105      ! Re-Compute the julian day at 1/2 time-step before the integration start
106      CALL ymds2ju( nyear, nmonth, nday, REAL(nsec_day), fjulday ) 
107
108      ! For long integration avoid truncation error at the start/end of a day
109      IF( ABS(fjulday - REAL(NINT(fjulday),wp)) < 0.1 / rday )   fjulday = REAL(NINT(fjulday),wp)   
110
111      nsec1jan000 = 0
112      CALL day_mth
113
114      IF (nmonth == 12) THEN
115         nsec1jan000 = nsec1jan000 - nsecd * nyear_len(0)
116      ENDIF 
117
118      ! day since january 1st
119      nday_year = nday + SUM( nmonth_len(1:nmonth - 1) )
120
121      !compute number of days between last monday and today     
122      CALL ymds2ju( 1900, 01, 01, 0.0, zjul )  ! compute julian day value of 01.01.1900 (our reference that was a Monday)
123      inbday = NINT(INT(fjulday) - zjul)       ! compute nb day between  01.01.1900 and current day 
124      idweek = MOD(inbday, 7)                  ! compute nb day between last monday and current day 
125
126      ! number of seconds since the beginning of current year/month/week/day at the middle of the time-step
127      nsec_week  = idweek    * nsecd + nsec_day 
128      nsec_month = (nday-1)  * nsecd + nsec_day 
129      nsec_year  = (nday_year-1) * nsecd + nsec_day 
130
131      ! control print
132      IF(lwp) WRITE(numout,'(a,i6,a,i2,a,i2,a,i6)')' ==============>> 1/2 time step before the start of the run DATE Y/M/D = ',   &
133           &                   nyear, '/', nmonth, '/', nday, '  nsec_day:', nsec_day, '  nsec_week:', nsec_week
134
135      ! Up to now, calendar parameters are related to the end of previous run (nit000-1)
136      ! call day to set the calendar parameters at the begining of the current simulaton. needed by iom_init
137      CALL day( nit000 )
138      !
139   END SUBROUTINE day_init
140
141
142   SUBROUTINE day_mth
143      !!----------------------------------------------------------------------
144      !!                   ***  ROUTINE day_init  ***
145      !!
146      !! ** Purpose :   calendar values related to the months
147      !!
148      !! ** Action  : - nmonth_len    : length in days of the months of the current year
149      !!              - nyear_len     : length in days of the previous/current year
150      !!              - nmonth_half   : second since the beginning of the year and the halft of the months
151      !!              - nmonth_end    : second since the beginning of the year and the end of the months
152      !!----------------------------------------------------------------------
153      INTEGER  ::   jm               ! dummy loop indice
154      !!----------------------------------------------------------------------
155
156      ! length of the month of the current year (from nleapy, read in namelist)
157      IF ( nleapy < 2 ) THEN
158         nmonth_len(:) = (/ 31, 31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31, 31 /)
159         nyear_len(:) = 365
160         IF ( nleapy == 1 ) THEN   ! we are using calandar with leap years
161            IF ( MOD(nyear-1, 4) == 0 .AND. ( MOD(nyear-1, 400) == 0 .OR. MOD(nyear-1, 100) /= 0 ) ) THEN
162               nyear_len(0) = 366
163            ENDIF
164            IF ( MOD(nyear, 4) == 0 .AND. ( MOD(nyear, 400) == 0 .OR. MOD(nyear, 100) /= 0 ) ) THEN
165               nmonth_len(2) = 29
166               nyear_len(1) = 366
167            ENDIF
168         ENDIF
169      ELSE
170         nmonth_len(:) = nleapy   ! all months with nleapy days per year
171         nyear_len(:) = 12 * nleapy
172      ENDIF
173
174      ! half month in second since the begining of the year:
175      ! time since Jan 1st   0     1     2    ...    11    12    13
176      !          ---------*--|--*--|--*--| ... |--*--|--*--|--*--|--------------------------------------
177      !                 <---> <---> <--->  ...  <---> <---> <--->       
178      ! month number      0     1     2    ...    11    12    13
179      !
180      ! nmonth_half(jm) = rday * REAL( 0.5 * nmonth_len(jm) + SUM(nmonth_len(1:jm-1)) )
181      nmonth_half(0) = - nsecd05 * nmonth_len(0)
182      DO jm = 1, 13
183         nmonth_half(jm) = nmonth_half(jm-1) + nsecd05 * ( nmonth_len(jm-1) + nmonth_len(jm) )
184      END DO
185
186      nmonth_end(0) = 0
187      DO jm = 1, 13
188         nmonth_end(jm) = nmonth_end(jm-1) + nsecd * nmonth_len(jm)
189      END DO
190      !           
191   END SUBROUTINE
192
193
194   SUBROUTINE day( kt )
195      !!----------------------------------------------------------------------
196      !!                      ***  ROUTINE day  ***
197      !!
198      !! ** Purpose :   Compute the date with a day iteration IF necessary.
199      !!
200      !! ** Method  : - ???
201      !!
202      !! ** Action  : - nyear     : current year
203      !!              - nmonth    : current month of the year nyear
204      !!              - nday      : current day of the month nmonth
205      !!              - nday_year : current day of the year nyear
206      !!              - ndastp    : = nyear*10000 + nmonth*100 + nday
207      !!              - adatrj    : date in days since the beginning of the run
208      !!              - nsec_year : current time of the year (in second since 00h, jan 1st)
209      !!----------------------------------------------------------------------     
210      INTEGER, INTENT(in) ::   kt        ! ocean time-step indices
211      !
212      CHARACTER (len=25) ::   charout
213      REAL(wp)           ::   zprec      ! fraction of day corresponding to 0.1 second
214      !!----------------------------------------------------------------------
215      !
216      IF( nn_timing == 1 )  CALL timing_start('day')
217      !
218      zprec = 0.1 / rday
219      !                                                 ! New time-step
220      nsec_year  = nsec_year  + ndt 
221      nsec_month = nsec_month + ndt                 
222      nsec_week  = nsec_week  + ndt
223      nsec_day   = nsec_day   + ndt               
224      adatrj  = adatrj  + rdttra(1) / rday
225      fjulday = fjulday + rdttra(1) / rday
226      IF( ABS(fjulday - REAL(NINT(fjulday),wp)) < zprec )   fjulday = REAL(NINT(fjulday),wp)   ! avoid truncation error
227      IF( ABS(adatrj  - REAL(NINT(adatrj ),wp)) < zprec )   adatrj  = REAL(NINT(adatrj ),wp)   ! avoid truncation error
228     
229      IF( nsec_day > nsecd ) THEN                       ! New day
230         !
231         nday      = nday + 1
232         nday_year = nday_year + 1
233         nsec_day  = ndt05
234         !
235         IF( nday == nmonth_len(nmonth) + 1 ) THEN      ! New month
236            nday   = 1
237            nmonth = nmonth + 1
238            nsec_month = ndt05
239            IF( nmonth == 13 ) THEN                     ! New year
240               nyear     = nyear + 1
241               nmonth    = 1
242               nday_year = 1
243               nsec_year = ndt05
244               nsec1jan000 = nsec1jan000 + nsecd * nyear_len(1)
245               IF( nleapy == 1 )   CALL day_mth
246            ENDIF
247         ENDIF
248         !
249         ndastp = nyear * 10000 + nmonth * 100 + nday   ! New date
250         !
251         !compute first day of the year in julian days
252         CALL ymds2ju( nyear, 01, 01, 0.0, fjulstartyear )
253         !
254         IF(lwp) WRITE(numout,'(a,i8,a,i4.4,a,i2.2,a,i2.2,a,i3.3)') '======>> time-step =', kt,   &
255              &   '      New day, DATE Y/M/D = ', nyear, '/', nmonth, '/', nday, '      nday_year = ', nday_year
256         IF(lwp) WRITE(numout,'(a,i8,a,i7,a,i5)') '         nsec_year = ', nsec_year,   &
257              &   '   nsec_month = ', nsec_month, '   nsec_day = ', nsec_day, '   nsec_week = ', nsec_week
258      ENDIF
259
260      IF( nsec_week > 7*nsecd )   nsec_week = ndt05     ! New week
261     
262      IF(ln_ctl) THEN
263         WRITE(charout,FMT="('kt =', I4,'  d/m/y =',I2,I2,I4)") kt, nday, nmonth, nyear
264         CALL prt_ctl_info(charout)
265      ENDIF
266
267      IF( .NOT. lk_offline ) CALL rst_opn( kt )               ! Open the restart file if needed and control lrst_oce
268      IF( lrst_oce         ) CALL day_rst( kt, 'WRITE' )      ! write day restart information
269      !
270      IF( nn_timing == 1 )  CALL timing_stop('day')
271      !
272   END SUBROUTINE day
273
274
275   SUBROUTINE day_rst( kt, cdrw )
276      !!---------------------------------------------------------------------
277      !!                   ***  ROUTINE ts_rst  ***
278      !!
279      !!  ** Purpose : Read or write calendar in restart file:
280      !!
281      !!  WRITE(READ) mode:
282      !!       kt        : number of time step since the begining of the experiment at the
283      !!                   end of the current(previous) run
284      !!       adatrj(0) : number of elapsed days since the begining of the experiment at the
285      !!                   end of the current(previous) run (REAL -> keep fractions of day)
286      !!       ndastp    : date at the end of the current(previous) run (coded as yyyymmdd integer)
287      !!
288      !!   According to namelist parameter nrstdt,
289      !!       nrstdt = 0  no control on the date (nit000 is  arbitrary).
290      !!       nrstdt = 1  we verify that nit000 is equal to the last
291      !!                   time step of previous run + 1.
292      !!       In both those options, the  exact duration of the experiment
293      !!       since the beginning (cumulated duration of all previous restart runs)
294      !!       is not stored in the restart and is assumed to be (nit000-1)*rdt.
295      !!       This is valid is the time step has remained constant.
296      !!
297      !!       nrstdt = 2  the duration of the experiment in days (adatrj)
298      !!                    has been stored in the restart file.
299      !!----------------------------------------------------------------------
300      INTEGER         , INTENT(in) ::   kt         ! ocean time-step
301      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag
302      !
303      REAL(wp) ::   zkt, zndastp
304      !!----------------------------------------------------------------------
305     
306      IF( TRIM(cdrw) == 'READ' ) THEN
307
308         IF( iom_varid( numror, 'kt', ldstop = .FALSE. ) > 0 ) THEN
309            ! Get Calendar informations
310            CALL iom_get( numror, 'kt', zkt )   ! last time-step of previous run
311            IF(lwp) THEN
312               WRITE(numout,*) ' *** Info read in restart : '
313               WRITE(numout,*) '   previous time-step                               : ', NINT( zkt )
314               WRITE(numout,*) ' *** restart option'
315               SELECT CASE ( nrstdt )
316               CASE ( 0 )   ;   WRITE(numout,*) ' nrstdt = 0 : no control of nit000'
317               CASE ( 1 )   ;   WRITE(numout,*) ' nrstdt = 1 : no control the date at nit000 (use ndate0 read in the namelist)'
318               CASE ( 2 )   ;   WRITE(numout,*) ' nrstdt = 2 : calendar parameters read in restart'
319               END SELECT
320               WRITE(numout,*)
321            ENDIF
322            ! Control of date
323            IF( nit000 - NINT( zkt ) /= 1 .AND. nrstdt /= 0 )                                         & 
324                 &   CALL ctl_stop( ' ===>>>> : problem with nit000 for the restart',                 & 
325                 &                  ' verify the restart file or rerun with nrstdt = 0 (namelist)' )
326            ! define ndastp and adatrj
327            IF ( nrstdt == 2 ) THEN 
328               ! read the parameters correspondting to nit000 - 1 (last time step of previous run)
329               CALL iom_get( numror, 'ndastp', zndastp )
330               ndastp = NINT( zndastp )
331               CALL iom_get( numror, 'adatrj', adatrj  )
332            ELSE 
333               ! parameters correspondting to nit000
334               ndastp = ndate0 
335               adatrj = ( REAL( nit000 -1, wp ) * rdttra(1) ) / rday 
336               ! note this is wrong if time step has changed during run
337            ENDIF
338         ELSE
339            ! parameters correspondting to nit000
340            ndastp = ndate0 
341            adatrj = ( REAL( nit000 -1, wp ) * rdttra(1) ) / rday 
342         ENDIF
343         IF( ABS(adatrj  - REAL(NINT(adatrj),wp)) < 0.1 / rday )   adatrj = REAL(NINT(adatrj),wp)   ! avoid truncation error
344         !
345         IF(lwp) THEN
346            WRITE(numout,*) ' *** Info used values : '
347            WRITE(numout,*) '   date ndastp                                      : ', ndastp
348            WRITE(numout,*) '   number of elapsed days since the begining of run : ', adatrj
349            WRITE(numout,*)
350         ENDIF
351         !
352      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN
353         !
354         IF( kt == nitrst ) THEN
355            IF(lwp) WRITE(numout,*)
356            IF(lwp) WRITE(numout,*) 'rst_write : write oce restart file  kt =', kt
357            IF(lwp) WRITE(numout,*) '~~~~~~~'         
358         ENDIF
359         ! calendar control
360         CALL iom_rstput( kt, nitrst, numrow, 'kt'     , REAL( kt    , wp) )   ! time-step
361         CALL iom_rstput( kt, nitrst, numrow, 'ndastp' , REAL( ndastp, wp) )   ! date
362         CALL iom_rstput( kt, nitrst, numrow, 'adatrj' , adatrj            )   ! number of elapsed days since
363         !                                                                     ! the begining of the run [s]
364      ENDIF
365      !
366   END SUBROUTINE day_rst
367
368   !!======================================================================
369END MODULE daymod