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.
trcrst.F90 in branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC – NEMO

source: branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/trcrst.F90 @ 6227

Last change on this file since 6227 was 6227, checked in by jpalmier, 8 years ago

JPALM -- 08-01-2016 -- MEDUSA debugg
File size: 32.1 KB
Line 
1MODULE trcrst
2   !!======================================================================
3   !!                         ***  MODULE trcrst  ***
4   !! TOP :   Manage the passive tracer restart
5   !!======================================================================
6   !! History :    -   !  1991-03  ()  original code
7   !!             1.0  !  2005-03 (O. Aumont, A. El Moussaoui) F90
8   !!              -   !  2005-10 (C. Ethe) print control
9   !!             2.0  !  2005-10 (C. Ethe, G. Madec) revised architecture
10   !!----------------------------------------------------------------------
11#if defined key_top
12   !!----------------------------------------------------------------------
13   !!   'key_top'                                                TOP models
14   !!----------------------------------------------------------------------
15   !!----------------------------------------------------------------------
16   !!   trc_rst :   Restart for passive tracer
17   !!----------------------------------------------------------------------
18   !!----------------------------------------------------------------------
19   !!   'key_top'                                                TOP models
20   !!----------------------------------------------------------------------
21   !!   trc_rst_opn    : open  restart file
22   !!   trc_rst_read   : read  restart file
23   !!   trc_rst_wri    : write restart file
24   !!----------------------------------------------------------------------
25   USE oce_trc
26   USE trc
27   USE trcnam_trp
28   USE iom
29   USE ioipsl, ONLY : ju2ymds    ! for calendar
30   USE daymod
31   !! AXY (05/11/13): need these for MEDUSA to input/output benthic reservoirs
32   USE sms_medusa
33   USE trcsms_medusa
34   !!
35#if defined key_idtra
36   USE trcsms_idtra
37#endif
38   !!
39#if defined key_cfc
40   USE trcsms_cfc
41#endif
42   IMPLICIT NONE
43   PRIVATE
44
45   PUBLIC   trc_rst_opn       ! called by ???
46   PUBLIC   trc_rst_read      ! called by ???
47   PUBLIC   trc_rst_wri       ! called by ???
48   PUBLIC   trc_rst_cal
49   PUBLIC   trc_rst_stat
50   PUBLIC   trc_rst_dia_stat
51
52   !! * Substitutions
53#  include "top_substitute.h90"
54   
55CONTAINS
56   
57   SUBROUTINE trc_rst_opn( kt )
58      !!----------------------------------------------------------------------
59      !!                    ***  trc_rst_opn  ***
60      !!
61      !! ** purpose  :   output of sea-trc variable in a netcdf file
62      !!----------------------------------------------------------------------
63      INTEGER, INTENT(in) ::   kt       ! number of iteration
64      INTEGER             ::   iyear, imonth, iday
65      REAL (wp)           ::   zsec
66      !
67      CHARACTER(LEN=20)   ::   clkt     ! ocean time-step define as a character
68      CHARACTER(LEN=50)   ::   clname   ! trc output restart file name
69      CHARACTER(LEN=256)  ::   clpath   ! full path to ocean output restart file
70      !!----------------------------------------------------------------------
71      !
72      IF( lk_offline ) THEN
73         IF( kt == nittrc000 ) THEN
74            lrst_trc = .FALSE.
75            IF( ln_rst_list ) THEN
76               nrst_lst = 1
77               nitrst = nstocklist( nrst_lst )
78            ELSE
79               nitrst = nitend
80            ENDIF
81         ENDIF
82
83         IF( .NOT. ln_rst_list .AND. MOD( kt - 1, nstock ) == 0 ) THEN
84            ! we use kt - 1 and not kt - nittrc000 to keep the same periodicity from the beginning of the experiment
85            nitrst = kt + nstock - 1                  ! define the next value of nitrst for restart writing
86            IF( nitrst > nitend )   nitrst = nitend   ! make sure we write a restart at the end of the run
87         ENDIF
88      ELSE
89         IF( kt == nittrc000 ) lrst_trc = .FALSE.
90      ENDIF
91
92      ! to get better performances with NetCDF format:
93      ! we open and define the tracer restart file one tracer time step before writing the data (-> at nitrst - 2*nn_dttrc + 1)
94      ! except if we write tracer restart files every tracer time step or if a tracer restart file was writen at nitend - 2*nn_dttrc + 1
95      IF( kt == nitrst - 2*nn_dttrc .OR. nstock == nn_dttrc .OR. ( kt == nitend - nn_dttrc .AND. .NOT. lrst_trc ) ) THEN
96         IF ( ln_rstdate ) THEN
97            !! JPALM -- 22-12-2015 -- modif to get the good date on restart trc file name
98            !!                     -- the condition to open the rst file is not the same than for the dynamic rst.
99            !!                     -- here it - for an obscure reason - is open 2 time-step before the restart writing process
100            !!                     instead of 1.
101            !!                     -- i am not sure if someone forgot +1 in the if loop condition as
102            !!                     it is writen in all comments nitrst -2*nn_dttrc + 1 and the condition is
103            !!                     nitrst - 2*nn_dttrc
104            !!                     -- nevertheless we didn't wanted to broke something already working
105            !!                     and just adapted the part we added.
106            !!                     -- So instead of calling ju2ymds( fjulday + (rdttra(1))
107            !!                     we call ju2ymds( fjulday + (2*rdttra(1))
108            !!--------------------------------------------------------------------     
109            CALL ju2ymds( fjulday + (2*rdttra(1)) / rday, iyear, imonth, iday, zsec )
110            WRITE(clkt, '(i4.4,2i2.2)') iyear, imonth, iday
111         ELSE
112            ! beware of the format used to write kt (default is i8.8, that should be large enough)
113            IF( nitrst > 1.0e9 ) THEN   ;   WRITE(clkt,*       ) nitrst
114            ELSE                        ;   WRITE(clkt,'(i8.8)') nitrst
115            ENDIF
116         ENDIF
117         ! create the file
118         IF(lwp) WRITE(numout,*)
119         clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_"//TRIM(cn_trcrst_out)
120         clpath = TRIM(cn_trcrst_outdir)
121         IF( clpath(LEN_TRIM(clpath):) /= '/' ) clpath = TRIM(clpath) // '/'
122         IF(lwp) WRITE(numout,*) &
123             '             open trc restart.output NetCDF file: ',TRIM(clpath)//clname
124         CALL iom_open( TRIM(clpath)//TRIM(clname), numrtw, ldwrt = .TRUE., kiolib = jprstlib )
125         lrst_trc = .TRUE.
126      ENDIF
127      !
128   END SUBROUTINE trc_rst_opn
129
130   SUBROUTINE trc_rst_read
131      !!----------------------------------------------------------------------
132      !!                    ***  trc_rst_opn  ***
133      !!
134      !! ** purpose  :   read passive tracer fields in restart files
135      !!----------------------------------------------------------------------
136      INTEGER  ::  jn, jl     
137      !! AXY (05/11/13): temporary variables
138      REAL(wp) ::    fq0,fq1,fq2
139
140      !!----------------------------------------------------------------------
141      !
142      IF(lwp) WRITE(numout,*)
143      IF(lwp) WRITE(numout,*) 'trc_rst_read : read data in the TOP restart file'
144      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
145
146      ! READ prognostic variables and computes diagnostic variable
147      DO jn = 1, jptra
148         CALL iom_get( numrtr, jpdom_autoglo, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) )
149      END DO
150
151      DO jn = 1, jptra
152         CALL iom_get( numrtr, jpdom_autoglo, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) )
153      END DO
154
155      !! AXY (09/06/14): the ARCHER version of MEDUSA-2 does not include an equivalent
156      !!                 call to MEDUSA-2 at this point; this suggests that the FCM
157      !!                 version of NEMO date significantly earlier than the current
158      !!                 version
159
160#if defined key_medusa
161      !! AXY (13/01/12): check if the restart contains sediment fields;
162      !!                 this is only relevant for simulations that include
163      !!                 biogeochemistry and are restarted from earlier runs
164      !!                 in which there was no sediment component
165      !!
166      IF( iom_varid( numrtr, 'B_SED_N', ldstop = .FALSE. ) > 0 ) THEN
167         !! YES; in which case read them
168         !!
169         IF(lwp) WRITE(numout,*) ' MEDUSA sediment fields present - reading in ...'
170         CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_N',  zb_sed_n(:,:)  )
171         CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_N',  zn_sed_n(:,:)  )
172         CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_FE', zb_sed_fe(:,:) )
173         CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_FE', zn_sed_fe(:,:) )
174         CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_SI', zb_sed_si(:,:) )
175         CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_SI', zn_sed_si(:,:) )
176         CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_C',  zb_sed_c(:,:)  )
177         CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_C',  zn_sed_c(:,:)  )
178         CALL iom_get( numrtr, jpdom_autoglo, 'B_SED_CA', zb_sed_ca(:,:) )
179         CALL iom_get( numrtr, jpdom_autoglo, 'N_SED_CA', zn_sed_ca(:,:) )
180      ELSE
181         !! NO; in which case set them to zero
182         !!
183         IF(lwp) WRITE(numout,*) ' MEDUSA sediment fields absent - setting to zero ...'
184         zb_sed_n(:,:)  = 0.0   !! organic N
185         zn_sed_n(:,:)  = 0.0
186         zb_sed_fe(:,:) = 0.0   !! organic Fe
187         zn_sed_fe(:,:) = 0.0
188         zb_sed_si(:,:) = 0.0   !! inorganic Si
189         zn_sed_si(:,:) = 0.0
190         zb_sed_c(:,:)  = 0.0   !! organic C
191         zn_sed_c(:,:)  = 0.0
192         zb_sed_ca(:,:) = 0.0   !! inorganic C
193         zn_sed_ca(:,:) = 0.0
194      ENDIF
195      !!
196      !! calculate stats on these fields
197      IF(lwp) WRITE(numout,*) ' MEDUSA sediment field stats (min, max, sum) ...'
198      fq0 = MINVAL(zn_sed_n(:,:))
199      fq1 = MAXVAL(zn_sed_n(:,:))
200      fq2 = SUM(zn_sed_n(:,:))
201      if (lwp) write (numout,'(a,3f15.5)') 'Sediment  N ', &
202         &        fq0, fq1, fq2
203      fq0 = MINVAL(zn_sed_fe(:,:))
204      fq1 = MAXVAL(zn_sed_fe(:,:))
205      fq2 = SUM(zn_sed_fe(:,:))
206      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Fe ', &
207         &        fq0, fq1, fq2
208      fq0 = MINVAL(zn_sed_si(:,:))
209      fq1 = MAXVAL(zn_sed_si(:,:))
210      fq2 = SUM(zn_sed_si(:,:))
211      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Si ', &
212         &        fq0, fq1, fq2
213      fq0 = MINVAL(zn_sed_c(:,:))
214      fq1 = MAXVAL(zn_sed_c(:,:))
215      fq2 = SUM(zn_sed_c(:,:))
216      if (lwp) write (numout,'(a,3f15.5)') 'Sediment  C ', &
217         &        fq0, fq1, fq2
218      fq0 = MINVAL(zn_sed_ca(:,:))
219      fq1 = MAXVAL(zn_sed_ca(:,:))
220      fq2 = SUM(zn_sed_ca(:,:))
221      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Ca ', &
222         &        fq0, fq1, fq2
223      !!
224      !! AXY (07/07/15): read in temporally averaged fields for DMS
225      !!                 calculations
226      !!
227      IF( iom_varid( numrtr, 'B_DMS_CHN', ldstop = .FALSE. ) > 0 ) THEN
228         !! YES; in which case read them
229         !!
230         IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS present - reading in ...'
231         CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_CHN',  zb_dms_chn(:,:)  )
232         CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_CHN',  zn_dms_chn(:,:)  )
233         CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_CHD',  zb_dms_chd(:,:)  )
234         CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_CHD',  zn_dms_chd(:,:)  )
235         CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_MLD',  zb_dms_mld(:,:)  )
236         CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_MLD',  zn_dms_mld(:,:)  )
237         CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_QSR',  zb_dms_qsr(:,:)  )
238         CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_QSR',  zn_dms_qsr(:,:)  )
239         CALL iom_get( numrtr, jpdom_autoglo, 'B_DMS_DIN',  zb_dms_din(:,:)  )
240         CALL iom_get( numrtr, jpdom_autoglo, 'N_DMS_DIN',  zn_dms_din(:,:)  )
241      ELSE
242         !! NO; in which case set them to zero
243         !!
244         IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS absent - setting to zero ...'
245         zb_dms_chn(:,:)  = 0.0   !! CHN
246         zn_dms_chn(:,:)  = 0.0
247         zb_dms_chd(:,:)  = 0.0   !! CHD
248         zn_dms_chd(:,:)  = 0.0
249         zb_dms_mld(:,:)  = 0.0   !! MLD
250         zn_dms_mld(:,:)  = 0.0
251         zb_dms_qsr(:,:)  = 0.0   !! QSR
252         zn_dms_qsr(:,:)  = 0.0
253         zb_dms_din(:,:)  = 0.0   !! DIN
254         zn_dms_din(:,:)  = 0.0
255      ENDIF
256      !!
257      !! calculate stats on these fields
258      IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS stats (min, max, sum) ...'
259      fq0 = MINVAL(zn_dms_chn(:,:))
260      fq1 = MAXVAL(zn_dms_chn(:,:))
261      fq2 = SUM(zn_dms_chn(:,:))
262      if (lwp) write (numout,'(a,3f15.5)') 'DMS, CHN ', fq0, fq1, fq2
263      fq0 = MINVAL(zn_dms_chd(:,:))
264      fq1 = MAXVAL(zn_dms_chd(:,:))
265      fq2 = SUM(zn_dms_chd(:,:))
266      if (lwp) write (numout,'(a,3f15.5)') 'DMS, CHD ', fq0, fq1, fq2
267      fq0 = MINVAL(zn_dms_mld(:,:))
268      fq1 = MAXVAL(zn_dms_mld(:,:))
269      fq2 = SUM(zn_dms_mld(:,:))
270      if (lwp) write (numout,'(a,3f15.5)') 'DMS, MLD ', fq0, fq1, fq2
271      fq0 = MINVAL(zn_dms_qsr(:,:))
272      fq1 = MAXVAL(zn_dms_qsr(:,:))
273      fq2 = SUM(zn_dms_qsr(:,:))
274      if (lwp) write (numout,'(a,3f15.5)') 'DMS, QSR ', fq0, fq1, fq2
275      fq0 = MINVAL(zn_dms_din(:,:))
276      fq1 = MAXVAL(zn_dms_din(:,:))
277      fq2 = SUM(zn_dms_din(:,:))
278      if (lwp) write (numout,'(a,3f15.5)') 'DMS, DIN ', fq0, fq1, fq2
279#endif
280      !
281#if defined key_idtra
282      !! JPALM -- 05-01-2016 -- mv idtra and CFC restart reading and
283      !!                        writting here undre their key.
284      !!                        problems in CFC restart, maybe because of this...
285      !!                        and pb in idtra diag or diad-restart writing.
286      !!----------------------------------------------------------------------
287      IF( iom_varid( numrtr, 'qint_idtra', ldstop = .FALSE. ) > 0 ) THEN
288         !! YES; in which case read them
289         !!
290         IF(lwp) WRITE(numout,*) ' IDTRA averaged properties present - reading in ...'
291         CALL iom_get( numrtr, jpdom_autoglo, 'qint_IDTRA',  qint_idtra(:,:,1)  )
292      ELSE
293         !! NO; in which case set them to zero
294         !!
295         IF(lwp) WRITE(numout,*) ' IDTRA averaged properties absent - setting to zero ...'
296         qint_idtra(:,:,1)  = 0.0   !! CHN
297      ENDIF
298      !!
299      !! calculate stats on these fields
300      call trc_rst_dia_stat(qint_idtra(:,:,1), 'qint_IDTRA')
301      IF(lwp) WRITE(numout,*) ' IDTRA averaged properties stats (min, max, sum) ...'
302      fq0 = MINVAL(qint_idtra(:,:,1))
303      fq1 = MAXVAL(qint_idtra(:,:,1))
304      fq2 = SUM(qint_idtra(:,:,1))
305      if (lwp) write (numout,'(a,3f15.5)') 'qint_IDTRA ', fq0, fq1, fq2
306#endif
307      !
308#if defined key_cfc
309      DO jl = 1, jp_cfc
310         jn = jp_cfc0 + jl - 1
311         IF( iom_varid( numrtr, 'qint_'//ctrcnm(jn), ldstop = .FALSE. ) > 0 ) THEN
312            !! YES; in which case read them
313            !!
314            IF(lwp) WRITE(numout,*) ' CFC averaged properties present - reading in ...'
315            CALL iom_get( numrtr, jpdom_autoglo, 'qint_'//ctrcnm(jn), qint_cfc(:,:,jl) )
316         ELSE
317            !! NO; in which case set them to zero
318            !!
319            IF(lwp) WRITE(numout,*) ' CFC averaged properties absent - setting to zero ...'
320            qint_cfc(:,:,jn)  = 0.0   !! CHN
321         ENDIF
322         !!
323         !! calculate stats on these fields
324         call trc_rst_dia_stat(qint_cfc(:,:,jl), 'qint_'//ctrcnm(jn))
325         IF(lwp) WRITE(numout,*) ' CFC averaged properties stats (min, max, sum) ...'
326         fq0 = MINVAL(qint_cfc(:,:,jl))
327         fq1 = MAXVAL(qint_cfc(:,:,jl))
328         fq2 = SUM(qint_cfc(:,:,jl))
329         if (lwp) write (numout,'(a,3f15.5)') 'qint_'//ctrcnm(jn), fq0, fq1, fq2
330      END DO
331#endif
332      !
333   END SUBROUTINE trc_rst_read
334
335   SUBROUTINE trc_rst_wri( kt )
336      !!----------------------------------------------------------------------
337      !!                    ***  trc_rst_wri  ***
338      !!
339      !! ** purpose  :   write passive tracer fields in restart files
340      !!----------------------------------------------------------------------
341      INTEGER, INTENT( in ) ::   kt    ! ocean time-step index
342      !!
343      INTEGER  :: jn, jl
344      REAL(wp) :: zarak0
345      !! AXY (05/11/13): temporary variables
346      REAL(wp) ::    fq0,fq1,fq2
347      !!----------------------------------------------------------------------
348      !
349      CALL iom_rstput( kt, nitrst, numrtw, 'rdttrc1', rdttrc(1) )   ! surface passive tracer time step
350      ! prognostic variables
351      ! --------------------
352      DO jn = 1, jptra
353         CALL iom_rstput( kt, nitrst, numrtw, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) )
354      END DO
355
356      DO jn = 1, jptra
357         CALL iom_rstput( kt, nitrst, numrtw, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) )
358      END DO
359
360      !! AXY (09/06/14): the ARCHER version of MEDUSA-2 does not include an equivalent
361      !!                 call to MEDUSA-2 at this point; this suggests that the FCM
362      !!                 version of NEMO date significantly earlier than the current
363      !!                 version
364
365#if defined key_medusa
366      !! AXY (13/01/12): write out "before" and "now" state of seafloor
367      !!                 sediment pools into restart; this happens
368      !!                 whether or not the pools are to be used by
369      !!                 MEDUSA (which is controlled by a switch in the
370      !!                 namelist_top file)
371      !!
372      IF(lwp) WRITE(numout,*) ' MEDUSA sediment fields - writing out ...'
373      CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_N',  zb_sed_n(:,:)  )
374      CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_N',  zn_sed_n(:,:)  )
375      CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_FE', zb_sed_fe(:,:) )
376      CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_FE', zn_sed_fe(:,:) )
377      CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_SI', zb_sed_si(:,:) )
378      CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_SI', zn_sed_si(:,:) )
379      CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_C',  zb_sed_c(:,:)  )
380      CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_C',  zn_sed_c(:,:)  )
381      CALL iom_rstput( kt, nitrst, numrtw, 'B_SED_CA', zb_sed_ca(:,:) )
382      CALL iom_rstput( kt, nitrst, numrtw, 'N_SED_CA', zn_sed_ca(:,:) )
383      !!
384      !! calculate stats on these fields
385      IF(lwp) WRITE(numout,*) ' MEDUSA sediment field stats (min, max, sum) ...'
386      fq0 = MINVAL(zn_sed_n(:,:))
387      fq1 = MAXVAL(zn_sed_n(:,:))
388      fq2 = SUM(zn_sed_n(:,:))
389      if (lwp) write (numout,'(a,3f15.5)') 'Sediment  N ', &
390         &        fq0, fq1, fq2
391      fq0 = MINVAL(zn_sed_fe(:,:))
392      fq1 = MAXVAL(zn_sed_fe(:,:))
393      fq2 = SUM(zn_sed_fe(:,:))
394      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Fe ', &
395         &        fq0, fq1, fq2
396      fq0 = MINVAL(zn_sed_si(:,:))
397      fq1 = MAXVAL(zn_sed_si(:,:))
398      fq2 = SUM(zn_sed_si(:,:))
399      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Si ', &
400         &        fq0, fq1, fq2
401      fq0 = MINVAL(zn_sed_c(:,:))
402      fq1 = MAXVAL(zn_sed_c(:,:))
403      fq2 = SUM(zn_sed_c(:,:))
404      if (lwp) write (numout,'(a,3f15.5)') 'Sediment  C ', &
405         &        fq0, fq1, fq2
406      fq0 = MINVAL(zn_sed_ca(:,:))
407      fq1 = MAXVAL(zn_sed_ca(:,:))
408      fq2 = SUM(zn_sed_ca(:,:))
409      if (lwp) write (numout,'(a,3f15.5)') 'Sediment Ca ', &
410         &        fq0, fq1, fq2
411      !!
412      !! AXY (07/07/15): write out temporally averaged fields for DMS
413      !!                 calculations
414      !!
415      IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS - writing out ...'
416      CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_CHN',  zb_dms_chn(:,:)  )
417      CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_CHN',  zn_dms_chn(:,:)  )
418      CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_CHD',  zb_dms_chd(:,:)  )
419      CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_CHD',  zn_dms_chd(:,:)  )
420      CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_MLD',  zb_dms_mld(:,:)  )
421      CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_MLD',  zn_dms_mld(:,:)  )
422      CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_QSR',  zb_dms_qsr(:,:)  )
423      CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_QSR',  zn_dms_qsr(:,:)  )
424      CALL iom_rstput( kt, nitrst, numrtw, 'B_DMS_DIN',  zb_dms_din(:,:)  )
425      CALL iom_rstput( kt, nitrst, numrtw, 'N_DMS_DIN',  zn_dms_din(:,:)  )
426      !!
427      !! calculate stats on these fields
428      IF(lwp) WRITE(numout,*) ' MEDUSA averaged properties for DMS stats (min, max, sum) ...'
429      fq0 = MINVAL(zn_dms_chn(:,:))
430      fq1 = MAXVAL(zn_dms_chn(:,:))
431      fq2 = SUM(zn_dms_chn(:,:))
432      if (lwp) write (numout,'(a,3f15.5)') 'DMS, CHN ', fq0, fq1, fq2
433      fq0 = MINVAL(zn_dms_chd(:,:))
434      fq1 = MAXVAL(zn_dms_chd(:,:))
435      fq2 = SUM(zn_dms_chd(:,:))
436      if (lwp) write (numout,'(a,3f15.5)') 'DMS, CHD ', fq0, fq1, fq2
437      fq0 = MINVAL(zn_dms_mld(:,:))
438      fq1 = MAXVAL(zn_dms_mld(:,:))
439      fq2 = SUM(zn_dms_mld(:,:))
440      if (lwp) write (numout,'(a,3f15.5)') 'DMS, MLD ', fq0, fq1, fq2
441      fq0 = MINVAL(zn_dms_qsr(:,:))
442      fq1 = MAXVAL(zn_dms_qsr(:,:))
443      fq2 = SUM(zn_dms_qsr(:,:))
444      if (lwp) write (numout,'(a,3f15.5)') 'DMS, QSR ', fq0, fq1, fq2
445      fq0 = MINVAL(zn_dms_din(:,:))
446      fq1 = MAXVAL(zn_dms_din(:,:))
447      fq2 = SUM(zn_dms_din(:,:))
448      if (lwp) write (numout,'(a,3f15.5)') 'DMS, DIN ', fq0, fq1, fq2
449      !!
450#endif
451      !
452#if defined key_idtra
453      !! JPALM -- 05-01-2016 -- mv idtra and CFC restart reading and
454      !!                        writting here undre their key.
455      !!                        problems in CFC restart, maybe because of this...
456      !!                        and pb in idtra diag or diad-restart writing.
457      !!----------------------------------------------------------------------
458      IF(lwp) WRITE(numout,*) ' IDTRA averaged properties - writing out ...'
459      CALL iom_rstput( kt, nitrst, numrtw, 'qint_IDTRA',  qint_idtra(:,:,1) )
460      !!
461      !! calculate stats on these fields
462      IF(lwp) WRITE(numout,*) ' IDTRA averaged properties stats (min, max, sum) ...'
463      call trc_rst_dia_stat(qint_idtra(:,:,1), 'qint_IDTRA')
464      fq0 = MINVAL(qint_idtra(:,:,1))
465      fq1 = MAXVAL(qint_idtra(:,:,1))
466      fq2 = SUM(qint_idtra(:,:,1))
467      if (lwp) write (numout,'(a,3f15.5)') 'qint_IDTRA ', fq0, fq1, fq2
468#endif
469      !
470#if defined key_cfc
471      DO jl = 1, jp_cfc
472         jn = jp_cfc0 + jl - 1
473         IF(lwp) WRITE(numout,*) ' CFC averaged properties - writing out ...'
474         CALL iom_rstput( kt, nitrst, numrtw, 'qint_'//ctrcnm(jn), qint_cfc(:,:,jl) )
475         !!
476         !! calculate stats on these fields
477         call trc_rst_dia_stat(qint_cfc(:,:,jl), 'qint_'//ctrcnm(jn))
478         IF(lwp) WRITE(numout,*) ' CFC averaged properties stats (min, max, sum) ...'
479         fq0 = MINVAL(qint_cfc(:,:,jl))
480         fq1 = MAXVAL(qint_cfc(:,:,jl))
481         fq2 = SUM(qint_cfc(:,:,jl))
482         if (lwp) write (numout,'(a,3f15.5)') 'qint_'//ctrcnm(jn), fq0, fq1, fq2
483      END DO
484#endif
485      !
486
487      IF( kt == nitrst ) THEN
488          CALL trc_rst_stat            ! statistics
489          CALL iom_close( numrtw )     ! close the restart file (only at last time step)
490#if ! defined key_trdmxl_trc
491          lrst_trc = .FALSE.
492#endif
493          IF( lk_offline .AND. ln_rst_list ) THEN
494             nrst_lst = nrst_lst + 1
495             nitrst = nstocklist( nrst_lst )
496          ENDIF
497      ENDIF
498      !
499   END SUBROUTINE trc_rst_wri 
500
501
502   SUBROUTINE trc_rst_cal( kt, cdrw )
503      !!---------------------------------------------------------------------
504      !!                   ***  ROUTINE trc_rst_cal  ***
505      !!
506      !!  ** Purpose : Read or write calendar in restart file:
507      !!
508      !!  WRITE(READ) mode:
509      !!       kt        : number of time step since the begining of the experiment at the
510      !!                   end of the current(previous) run
511      !!       adatrj(0) : number of elapsed days since the begining of the experiment at the
512      !!                   end of the current(previous) run (REAL -> keep fractions of day)
513      !!       ndastp    : date at the end of the current(previous) run (coded as yyyymmdd integer)
514      !!
515      !!   According to namelist parameter nrstdt,
516      !!       nn_rsttr = 0  no control on the date (nittrc000 is  arbitrary).
517      !!       nn_rsttr = 1  we verify that nittrc000 is equal to the last
518      !!                   time step of previous run + 1.
519      !!       In both those options, the  exact duration of the experiment
520      !!       since the beginning (cumulated duration of all previous restart runs)
521      !!       is not stored in the restart and is assumed to be (nittrc000-1)*rdt.
522      !!       This is valid is the time step has remained constant.
523      !!
524      !!       nn_rsttr = 2  the duration of the experiment in days (adatrj)
525      !!                    has been stored in the restart file.
526      !!----------------------------------------------------------------------
527      INTEGER         , INTENT(in) ::   kt         ! ocean time-step
528      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag
529      !
530      INTEGER  ::  jlibalt = jprstlib
531      LOGICAL  ::  llok
532      REAL(wp) ::  zkt, zrdttrc1
533      REAL(wp) ::  zndastp
534
535      ! Time domain : restart
536      ! ---------------------
537
538      IF( TRIM(cdrw) == 'READ' ) THEN
539
540         IF(lwp) WRITE(numout,*)
541         IF(lwp) WRITE(numout,*) 'trc_rst_cal : read the TOP restart file for calendar'
542         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
543
544         IF ( jprstlib == jprstdimg ) THEN
545           ! eventually read netcdf file (monobloc)  for restarting on different number of processors
546           ! if {cn_trcrst_in}.nc exists, then set jlibalt to jpnf90
547           INQUIRE( FILE = TRIM(cn_trcrst_indir)//'/'//TRIM(cn_trcrst_in)//'.nc', EXIST = llok )
548           IF ( llok ) THEN ; jlibalt = jpnf90  ; ELSE ; jlibalt = jprstlib ; ENDIF
549         ENDIF
550
551         IF( ln_rsttr ) THEN
552            CALL iom_open( TRIM(cn_trcrst_indir)//'/'//cn_trcrst_in, numrtr, kiolib = jlibalt )
553            CALL iom_get ( numrtr, 'kt', zkt )   ! last time-step of previous run
554
555            IF(lwp) THEN
556               WRITE(numout,*) ' *** Info read in restart : '
557               WRITE(numout,*) '   previous time-step                               : ', NINT( zkt )
558               WRITE(numout,*) ' *** restart option'
559               SELECT CASE ( nn_rsttr )
560               CASE ( 0 )   ;   WRITE(numout,*) ' nn_rsttr = 0 : no control of nittrc000'
561               CASE ( 1 )   ;   WRITE(numout,*) ' nn_rsttr = 1 : no control the date at nittrc000 (use ndate0 read in the namelist)'
562               CASE ( 2 )   ;   WRITE(numout,*) ' nn_rsttr = 2 : calendar parameters read in restart'
563               END SELECT
564               WRITE(numout,*)
565            ENDIF
566            ! Control of date
567            IF( nittrc000  - NINT( zkt ) /= nn_dttrc .AND.  nn_rsttr /= 0 )                                  &
568               &   CALL ctl_stop( ' ===>>>> : problem with nittrc000 for the restart',                 &
569               &                  ' verify the restart file or rerun with nn_rsttr = 0 (namelist)' )
570         ENDIF
571         !
572         IF( lk_offline ) THEN   
573            !                                          ! set the date in offline mode
574            IF( ln_rsttr .AND. nn_rsttr == 2 ) THEN
575               CALL iom_get( numrtr, 'ndastp', zndastp ) 
576               ndastp = NINT( zndastp )
577               CALL iom_get( numrtr, 'adatrj', adatrj  )
578             ELSE
579               ndastp = ndate0 - 1     ! ndate0 read in the namelist in dom_nam
580               adatrj = ( REAL( nittrc000-1, wp ) * rdttra(1) ) / rday
581               ! note this is wrong if time step has changed during run
582            ENDIF
583            !
584            IF(lwp) THEN
585              WRITE(numout,*) ' *** Info used values : '
586              WRITE(numout,*) '   date ndastp                                      : ', ndastp
587              WRITE(numout,*) '   number of elapsed days since the begining of run : ', adatrj
588              WRITE(numout,*)
589            ENDIF
590            !
591            IF( ln_rsttr )  THEN   ;    neuler = 1
592            ELSE                   ;    neuler = 0
593            ENDIF
594            !
595            CALL day_init          ! compute calendar
596            !
597         ENDIF
598         !
599      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN
600         !
601         IF(  kt == nitrst ) THEN
602            IF(lwp) WRITE(numout,*)
603            IF(lwp) WRITE(numout,*) 'trc_wri : write the TOP restart file (NetCDF) at it= ', kt, ' date= ', ndastp
604            IF(lwp) WRITE(numout,*) '~~~~~~~'
605         ENDIF
606         CALL iom_rstput( kt, nitrst, numrtw, 'kt'     , REAL( kt    , wp) )   ! time-step
607         CALL iom_rstput( kt, nitrst, numrtw, 'ndastp' , REAL( ndastp, wp) )   ! date
608         CALL iom_rstput( kt, nitrst, numrtw, 'adatrj' , adatrj            )   ! number of elapsed days since
609         !                                                                     ! the begining of the run [s]
610      ENDIF
611
612   END SUBROUTINE trc_rst_cal
613
614
615   SUBROUTINE trc_rst_stat
616      !!----------------------------------------------------------------------
617      !!                    ***  trc_rst_stat  ***
618      !!
619      !! ** purpose  :   Compute tracers statistics
620      !!----------------------------------------------------------------------
621      INTEGER  :: jk, jn
622      REAL(wp) :: ztraf, zmin, zmax, zmean, zdrift
623      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zvol
624      !!----------------------------------------------------------------------
625
626      IF( lwp ) THEN
627         WRITE(numout,*) 
628         WRITE(numout,*) '           ----TRACER STAT----             '
629         WRITE(numout,*) 
630      ENDIF
631      !
632      DO jk = 1, jpk
633         zvol(:,:,jk) = e1e2t(:,:) * fse3t_a(:,:,jk) * tmask(:,:,jk)
634      END DO
635      !
636      DO jn = 1, jptra
637         ztraf = glob_sum( trn(:,:,:,jn) * zvol(:,:,:) )
638         zmin  = MINVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) )
639         zmax  = MAXVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) )
640         IF( lk_mpp ) THEN
641            CALL mpp_min( zmin )      ! min over the global domain
642            CALL mpp_max( zmax )      ! max over the global domain
643         END IF
644         zmean  = ztraf / areatot
645         zdrift = ( ( ztraf - trai(jn) ) / ( trai(jn) + 1.e-12 )  ) * 100._wp
646         IF(lwp) WRITE(numout,9000) jn, TRIM( ctrcnm(jn) ), zmean, zmin, zmax, zdrift
647      END DO
648      IF(lwp) WRITE(numout,*) 
6499000  FORMAT(' tracer nb :',i2,'    name :',a10,'    mean :',e18.10,'    min :',e18.10, &
650      &      '    max :',e18.10,'    drift :',e18.10, ' %')
651      !
652   END SUBROUTINE trc_rst_stat
653
654
655   SUBROUTINE trc_rst_dia_stat( dgtr, names)
656      !!----------------------------------------------------------------------
657      !!                    ***  trc_rst_dia_stat  ***
658      !!
659      !! ** purpose  :   Compute tracers statistics
660      !!----------------------------------------------------------------------
661      REAL(wp), DIMENSION(jpi,jpj) , INTENT(in) ::   dgtr      ! 2D diag var
662      CHARACTER(len=*)             , INTENT(in) ::   names     ! 2D diag name
663      !!---------------------------------------------------------------------
664      INTEGER  :: jk, jn
665      REAL(wp) :: ztraf, zmin, zmax, zmean, areasf
666      REAL(wp), DIMENSION(jpi,jpj) :: zvol
667      !!----------------------------------------------------------------------
668
669      IF( lwp )  WRITE(numout,*) 'STAT- ', names
670      !
671      zvol(:,:) = e1e2t(:,:) * fse3t_a(:,:,1) * tmask(:,:,1)
672      ztraf = glob_sum( dgtr(:,:) * zvol(:,:) )
673      areasf = glob_sum(e1e2t(:,:) * tmask(:,:,1) )
674      zmin  = MINVAL( dgtr(:,:), mask= ((tmask(:,:,1).NE.0.)) )
675      zmax  = MAXVAL( dgtr(:,:), mask= ((tmask(:,:,1).NE.0.)) )
676      IF( lk_mpp ) THEN
677         CALL mpp_min( zmin )      ! min over the global domain
678         CALL mpp_max( zmax )      ! max over the global domain
679      END IF
680      zmean  = ztraf / areatot
681      IF(lwp) WRITE(numout,9001) TRIM( names ), zmean, zmin, zmax
682      !
683      IF(lwp) WRITE(numout,*)
6849001  FORMAT(' tracer name :',a10,'    mean :',e18.10,'    min :',e18.10, &
685      &      '    max :',e18.10 )
686      !
687   END SUBROUTINE trc_rst_dia_stat
688
689
690#else
691   !!----------------------------------------------------------------------
692   !!  Dummy module :                                     No passive tracer
693   !!----------------------------------------------------------------------
694CONTAINS
695   SUBROUTINE trc_rst_read                      ! Empty routines
696   END SUBROUTINE trc_rst_read
697   SUBROUTINE trc_rst_wri( kt )
698      INTEGER, INTENT ( in ) :: kt
699      WRITE(*,*) 'trc_rst_wri: You should not have seen this print! error?', kt
700   END SUBROUTINE trc_rst_wri   
701#endif
702
703   !!----------------------------------------------------------------------
704   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
705   !! $Id$
706   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
707   !!======================================================================
708END MODULE trcrst
Note: See TracBrowser for help on using the repository browser.