source: NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/SED/sedrst.F90 @ 10345

Last change on this file since 10345 was 10345, checked in by smasson, 2 years ago

dev_r10164_HPC09_ESIWACE_PREP_MERGE: merge with trunk@10344, see #2133

  • Property svn:keywords set to Id
File size: 16.1 KB
Line 
1MODULE sedrst
2   !!======================================================================
3   !!                       *** MODULE sedrst ***
4   !!   Read and write the restart files for sediment
5   !!======================================================================
6
7   !!----------------------------------------------------------------------
8   !! * Modules used
9   !! ==============
10   USE sed
11   USE sedarr
12   USE trc_oce, ONLY : l_offline, nn_dttrc
13   USE phycst , ONLY : rday
14   USE iom
15   USE daymod
16   USE lib_mpp         ! distribued memory computing library
17
18
19   !! * Accessibility
20   IMPLICIT NONE
21   PRIVATE
22
23   !! * Accessibility
24   PUBLIC sed_rst_opn       ! called by ???
25   PUBLIC sed_rst_read
26   PUBLIC sed_rst_wri
27   PUBLIC sed_rst_cal
28
29   !! $Id$
30CONTAINS
31
32
33   SUBROUTINE sed_rst_opn( kt )
34      !!----------------------------------------------------------------------
35      !!                    ***  sed_rst_opn  ***
36      !!
37      !! ** purpose  :   output of sed-trc variable in a netcdf file
38      !!----------------------------------------------------------------------
39      INTEGER, INTENT(in) ::   kt       ! number of iteration
40      !
41      CHARACTER(LEN=20)   ::   clkt     ! ocean time-step define as a character
42      CHARACTER(LEN=50)   ::   clname   ! trc output restart file name
43      CHARACTER(LEN=256)  ::   clpath   ! full path to ocean output restart file
44      !!----------------------------------------------------------------------
45      !
46      IF( l_offline ) THEN
47         IF( kt == nittrc000 ) THEN
48            lrst_sed = .FALSE.
49            IF( ln_rst_list ) THEN
50               nrst_lst = 1
51               nitrst = nstocklist( nrst_lst )
52            ELSE
53               nitrst = nitend
54            ENDIF
55         ENDIF
56         IF( .NOT. ln_rst_list .AND. MOD( kt - 1, nstock ) == 0 ) THEN
57            ! we use kt - 1 and not kt - nittrc000 to keep the same periodicity from the beginning of the experiment
58            nitrst = kt + nstock - 1                  ! define the next value of nitrst for restart writing
59            IF( nitrst > nitend )   nitrst = nitend   ! make sure we write a restart at the end of the run
60         ENDIF
61      ELSE
62         IF( kt == nittrc000 ) lrst_sed = .FALSE.
63      ENDIF
64
65      ! to get better performances with NetCDF format:
66      ! we open and define the tracer restart file one tracer time step before writing the data (-> at nitrst - 2*nn_dttrc + 1)
67      ! 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
68      IF( kt == nitrst - 2*nn_dtsed .OR. nstock == nn_dtsed .OR. ( kt == nitend - nn_dtsed .AND. .NOT. lrst_sed ) ) THEN
69         ! beware of the format used to write kt (default is i8.8, that should be large enough)
70         IF( nitrst > 1.0e9 ) THEN   ;   WRITE(clkt,*       ) nitrst
71         ELSE                        ;   WRITE(clkt,'(i8.8)') nitrst
72         ENDIF
73         ! create the file
74         IF(lwp) WRITE(numsed,*)
75         clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_"//TRIM(cn_sedrst_out)
76         clpath = TRIM(cn_sedrst_outdir)
77         IF( clpath(LEN_TRIM(clpath):) /= '/' ) clpath = TRIM(clpath) // '/'
78         IF(lwp) WRITE(numsed,*) &
79             '             open sed restart.output NetCDF file: ',TRIM(clpath)//clname
80         CALL iom_open( TRIM(clpath)//TRIM(clname), numrsw, ldwrt = .TRUE., kiolib = jprstlib, kdlev = jpksed )
81         lrst_sed = .TRUE.
82      ENDIF
83      !
84   END SUBROUTINE sed_rst_opn
85
86   SUBROUTINE sed_rst_read 
87      !!----------------------------------------------------------------------
88      !!                   ***  ROUTINE sed_rst_read  ***
89      !!
90      !! ** Purpose :  Initialization of sediment module
91      !!               - sets initial sediment composition
92      !!                 ( only clay or reading restart file )
93      !!
94      !!   History :
95      !!        !  06-07  (C. Ethe)  original
96      !!----------------------------------------------------------------------
97
98      !! * local declarations
99      INTEGER :: ji, jj, jk, jn 
100      REAL(wp), DIMENSION(jpi,jpj,jpksed,jptrased) :: zdta
101      REAL(wp), DIMENSION(jpi,jpj,jpksed,2)        :: zdta1 
102      REAL(wp), DIMENSION(jpi,jpj,jpksed)          :: zdta2
103      REAL(wp), DIMENSION(jpoce,jpksed)            :: zhipor
104      REAL(wp) :: zkt
105      CHARACTER(len = 20) ::   cltra
106      CHARACTER(LEN=20)   ::   name1
107      INTEGER             ::   jlibalt = jprstlib
108      LOGICAL             ::   llok
109      !--------------------------------------------------------------------
110
111      IF( ln_timing )  CALL timing_start('sed_rst_read')
112
113      IF (lwp) WRITE(numsed,*) ' '     
114      IF (lwp) WRITE(numsed,*) ' Initilization of Sediment components from restart'
115      IF (lwp) WRITE(numsed,*) ' '
116
117      zdta  = 1.
118      zdta1 = 1.
119      zdta2 = 0.
120
121      DO jn = 1, jptrased
122         cltra = TRIM(sedtrcd(jn))
123         IF( iom_varid( numrsr, TRIM(cltra) , ldstop = .FALSE. ) > 0 ) THEN
124            CALL iom_get( numrsr, jpdom_autoglo, TRIM(cltra), zdta(:,:,:,jn) )
125         ELSE
126            zdta(:,:,:,jn) = 0.0
127         ENDIF
128      ENDDO
129
130      DO jn = 1, jpsol
131         CALL pack_arr( jpoce, solcp(1:jpoce,1:jpksed,jn), &
132         &              zdta(1:jpi,1:jpj,1:jpksed,jn), iarroce(1:jpoce) )
133      END DO
134
135      DO jn = 1, jpwat
136         CALL pack_arr( jpoce, pwcp(1:jpoce,1:jpksed,jn), &
137         &              zdta(1:jpi,1:jpj,1:jpksed,jpsol+jn), iarroce(1:jpoce) )
138      END DO
139
140      DO jn = 1, 2
141         cltra = TRIM(seddia3d(jn))
142         IF( iom_varid( numrsr, TRIM(cltra) , ldstop = .FALSE. ) > 0 ) THEN
143            CALL iom_get( numrsr, jpdom_autoglo, TRIM(cltra), zdta1(:,:,:,jn) )
144         ELSE
145            zdta1(:,:,:,jn) = 0.0
146         ENDIF
147      ENDDO
148
149      zhipor(:,:) = 0.
150      CALL pack_arr( jpoce, zhipor(1:jpoce,1:jpksed), &
151         &             zdta1(1:jpi,1:jpj,1:jpksed,1), iarroce(1:jpoce) )
152
153      ! Initialization of [h+] in mol/kg
154      DO jk = 1, jpksed
155         DO ji = 1, jpoce
156            hipor (ji,jk) = 10.**( -1. * zhipor(ji,jk) )
157         ENDDO
158      ENDDO
159
160      CALL pack_arr( jpoce, co3por(1:jpoce,1:jpksed), &
161         &             zdta1(1:jpi,1:jpj,1:jpksed,2), iarroce(1:jpoce) )
162
163      ! Initialization of sediment composant only ie jk=2 to jk=jpksed
164      ! ( nothing in jk=1)
165      solcp(1:jpoce,1,:) = 0.
166      pwcp (1:jpoce,1,:) = 0.
167
168      cltra = "dbioturb"
169      IF( iom_varid( numrsr, TRIM(cltra) , ldstop = .FALSE. ) > 0 ) THEN
170         CALL iom_get( numrsr, jpdom_autoglo, TRIM(cltra), zdta2(:,:,:) )
171      ELSE
172         zdta2(:,:,:) = 0.0
173      ENDIF
174
175      CALL pack_arr( jpoce, db(1:jpoce,1:jpksed), &
176         &             zdta2(1:jpi,1:jpj,1:jpksed), iarroce(1:jpoce) )
177
178      cltra = "irrig"
179      IF( iom_varid( numrsr, TRIM(cltra) , ldstop = .FALSE. ) > 0 ) THEN
180         CALL iom_get( numrsr, jpdom_autoglo, TRIM(cltra), zdta2(:,:,:) )
181      ELSE
182         zdta2(:,:,:) = 0.0
183      ENDIF
184
185      CALL pack_arr( jpoce, irrig(1:jpoce,1:jpksed), &
186         &             zdta2(1:jpi,1:jpj,1:jpksed), iarroce(1:jpoce) )
187
188      cltra = "sedligand"
189      IF( iom_varid( numrsr, TRIM(cltra) , ldstop = .FALSE. ) > 0 ) THEN
190         CALL iom_get( numrsr, jpdom_autoglo, TRIM(cltra), zdta2(:,:,:) )
191      ELSE
192         zdta2(:,:,:) = 0.0
193      ENDIF
194
195      CALL pack_arr( jpoce, sedligand(1:jpoce,1:jpksed), &
196         &             zdta2(1:jpi,1:jpj,1:jpksed), iarroce(1:jpoce) )
197
198      IF( ln_timing )  CALL timing_stop('sed_rst_read')
199     
200   END SUBROUTINE sed_rst_read
201
202   SUBROUTINE sed_rst_wri( kt )
203      !!----------------------------------------------------------------------
204      !!                   ***  ROUTINE sed_rst_wri  ***
205      !!
206      !! ** Purpose :  save field which are necessary for sediment restart
207      !!
208      !!   History :
209      !!        !  06-07  (C. Ethe)  original
210      !!----------------------------------------------------------------------
211      !!* Modules used
212      INTEGER, INTENT(in) ::   kt       ! number of iteration
213      !! * local declarations
214      INTEGER  :: ji, jj, jk, jn
215      REAL(wp), DIMENSION(1) ::  zinfo
216      CHARACTER(len=50) :: clname
217      CHARACTER(len=20) :: cltra, name1 
218      REAL(wp), DIMENSION(jpoce,jpksed)   :: zdta   
219      REAL(wp), DIMENSION(jpi,jpj,jpksed) :: zdta2
220      !! -----------------------------------------------------------------------
221
222      IF( ln_timing )  CALL timing_start('sed_rst_wri')
223
224         !! 0. initialisations
225         !! ------------------
226         
227      IF(lwp) WRITE(numsed,*) ' '
228      IF(lwp) WRITE(numsed,*) 'sed_rst_write : write the sediment restart file in NetCDF format ',   &
229            'at it= ',kt
230      IF(lwp) WRITE(numsed,*) '~~~~~~~~~'
231
232
233      trcsedi(:,:,:,:)   = 0.0
234      flxsedi3d(:,:,:,:) = 0.0
235      zdta(:,:)          = 1.0
236      zdta2(:,:,:)       = 0.0
237
238         
239      !! 1. WRITE in nutwrs
240      !! ------------------
241
242      zinfo(1) = REAL( kt)
243      CALL iom_rstput( kt, nitrst, numrsw, 'kt', zinfo  )
244
245      ! Back to 2D geometry
246      DO jn = 1, jpsol
247         CALL unpack_arr( jpoce, trcsedi(1:jpi,1:jpj,1:jpksed,jn) , iarroce(1:jpoce), &
248         &                       solcp(1:jpoce,1:jpksed,jn ) )
249      END DO
250
251      DO jn = 1, jpwat
252         CALL unpack_arr( jpoce, trcsedi(1:jpi,1:jpj,1:jpksed,jpsol+jn) , iarroce(1:jpoce), &
253         &                       pwcp(1:jpoce,1:jpksed,jn  )  )
254      END DO
255      ! pH
256      DO jk = 1, jpksed
257         DO ji = 1, jpoce
258            zdta(ji,jk) = -LOG10( hipor(ji,jk) / ( densSW(ji) + rtrn ) + rtrn )
259         ENDDO
260      ENDDO
261
262      CALL unpack_arr( jpoce, flxsedi3d(1:jpi,1:jpj,1:jpksed,1)  , iarroce(1:jpoce), &
263      &                   zdta(1:jpoce,1:jpksed)  )
264         
265      CALL unpack_arr( jpoce, flxsedi3d(1:jpi,1:jpj,1:jpksed,2)  , iarroce(1:jpoce), &
266      &                   co3por(1:jpoce,1:jpksed)  )
267
268      ! prognostic variables
269      ! --------------------
270
271      DO jn = 1, jptrased
272         cltra = TRIM(sedtrcd(jn))
273         CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), trcsedi(:,:,:,jn) )
274      ENDDO
275
276      DO jn = 1, 2
277         cltra = TRIM(seddia3d(jn))
278         CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), flxsedi3d(:,:,:,jn) )
279      ENDDO
280
281      CALL unpack_arr( jpoce, zdta2(1:jpi,1:jpj,1:jpksed)  , iarroce(1:jpoce), &
282      &                   db(1:jpoce,1:jpksed)  )
283
284      cltra = "dbioturb"
285      CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), zdta2(:,:,:) )
286
287      CALL unpack_arr( jpoce, zdta2(1:jpi,1:jpj,1:jpksed)  , iarroce(1:jpoce), &
288      &                   irrig(1:jpoce,1:jpksed)  )
289
290      cltra = "irrig"
291      CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), zdta2(:,:,:) )
292
293      CALL unpack_arr( jpoce, zdta2(1:jpi,1:jpj,1:jpksed)  , iarroce(1:jpoce), &
294      &                   sedligand(1:jpoce,1:jpksed)  )
295
296      cltra = "sedligand"
297      CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), zdta2(:,:,:) )
298
299      IF( kt == nitrst ) THEN
300          CALL iom_close( numrsw )     ! close the restart file (only at last time step)
301          IF( l_offline .AND. ln_rst_list ) THEN
302             nrst_lst = nrst_lst + 1
303             nitrst = nstocklist( nrst_lst )
304          ENDIF
305      ENDIF
306
307      IF( ln_timing )  CALL timing_stop('sed_rst_wri')
308         
309   END SUBROUTINE sed_rst_wri
310
311
312   SUBROUTINE sed_rst_cal( kt, cdrw )
313      !!---------------------------------------------------------------------
314      !!                   ***  ROUTINE sed_rst_cal  ***
315      !!
316      !!  ** Purpose : Read or write calendar in restart file:
317      !!
318      !!  WRITE(READ) mode:
319      !!       kt        : number of time step since the begining of the experiment at the
320      !!                   end of the current(previous) run
321      !!       adatrj(0) : number of elapsed days since the begining of the experiment at the
322      !!                   end of the current(previous) run (REAL -> keep fractions of day)
323      !!       ndastp    : date at the end of the current(previous) run (coded as yyyymmdd integer)
324      !!
325      !!   According to namelist parameter nrstdt,
326      !!       nn_rsttr = 0  no control on the date (nittrc000 is  arbitrary).
327      !!       nn_rsttr = 1  we verify that nittrc000 is equal to the last
328      !!                   time step of previous run + 1.
329      !!       In both those options, the  exact duration of the experiment
330      !!       since the beginning (cumulated duration of all previous restart runs)
331      !!       is not stored in the restart and is assumed to be (nittrc000-1)*rdt.
332      !!       This is valid is the time step has remained constant.
333      !!
334      !!       nn_rsttr = 2  the duration of the experiment in days (adatrj)
335      !!                    has been stored in the restart file.
336      !!----------------------------------------------------------------------
337      INTEGER         , INTENT(in) ::   kt         ! ocean time-step
338      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag
339      !
340      INTEGER  ::  jlibalt = jprstlib
341      LOGICAL  ::  llok
342      REAL(wp) ::  zkt, zrdttrc1
343      REAL(wp) ::  zndastp
344
345      ! Time domain : restart
346      ! ---------------------
347
348      IF( TRIM(cdrw) == 'READ' ) THEN
349
350         IF(lwp) WRITE(numsed,*)
351         IF(lwp) WRITE(numsed,*) 'sed_rst_cal : read the SED restart file for calendar'
352         IF(lwp) WRITE(numsed,*) '~~~~~~~~~~~~'
353
354         IF( ln_rst_sed ) THEN
355            CALL iom_open( TRIM(cn_sedrst_indir)//'/'//cn_sedrst_in, numrsr, kiolib = jlibalt )
356            CALL iom_get ( numrsr, 'kt', zkt )   ! last time-step of previous run
357
358            IF(lwp) THEN
359               WRITE(numsed,*) ' *** Info read in restart : '
360               WRITE(numsed,*) '   previous time-step                               : ', NINT( zkt )
361               WRITE(numsed,*) ' *** restart option'
362               SELECT CASE ( nn_rstsed )
363               CASE ( 0 )   ;   WRITE(numsed,*) ' nn_rstsed = 0 : no control of nittrc000'
364               CASE ( 1 )   ;   WRITE(numsed,*) ' nn_rstsed = 1 : no control the date at nittrc000 (use ndate0 read in the namelist)'
365               CASE ( 2 )   ;   WRITE(numsed,*) ' nn_rstsed = 2 : calendar parameters read in restart'
366               END SELECT
367               WRITE(numsed,*)
368            ENDIF
369            ! Control of date
370            IF( nittrc000  - NINT( zkt ) /= nn_dtsed .AND.  nn_rstsed /= 0 )                                  &
371               &   CALL ctl_stop( ' ===>>>> : problem with nittrc000 for the restart',                 &
372               &                  ' verify the restart file or rerun with nn_rsttr = 0 (namelist)' )
373         ENDIF
374         !
375         IF( l_offline ) THEN
376            !                                          ! set the date in offline mode
377            IF( ln_rst_sed .AND. nn_rstsed == 2 ) THEN
378               CALL iom_get( numrsr, 'ndastp', zndastp )
379               ndastp = NINT( zndastp )
380               CALL iom_get( numrsr, 'adatrj', adatrj  )
381             ELSE
382               ndastp = ndate0 - 1     ! ndate0 read in the namelist in dom_nam
383               adatrj = ( REAL( nittrc000-1, wp ) * rdt ) / rday
384               ! note this is wrong if time step has changed during run
385            ENDIF
386            !
387            IF(lwp) THEN
388              WRITE(numsed,*) ' *** Info used values : '
389              WRITE(numsed,*) '   date ndastp                                      : ', ndastp
390              WRITE(numsed,*) '   number of elapsed days since the begining of run : ', adatrj
391              WRITE(numsed,*)
392            ENDIF
393            !
394            CALL day_init          ! compute calendar
395            !
396         ENDIF
397         !
398      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN
399         !
400         IF(  kt == nitrst ) THEN
401            IF(lwp) WRITE(numsed,*)
402            IF(lwp) WRITE(numsed,*) 'trc_wri : write the TOP restart file (NetCDF) at it= ', kt, ' date= ', ndastp
403            IF(lwp) WRITE(numsed,*) '~~~~~~~'
404         ENDIF
405         CALL iom_rstput( kt, nitrst, numrsw, 'kt'     , REAL( kt    , wp) )   ! time-step
406         CALL iom_rstput( kt, nitrst, numrsw, 'ndastp' , REAL( ndastp, wp) )   ! date
407         CALL iom_rstput( kt, nitrst, numrsw, 'adatrj' , adatrj            )   ! number of elapsed days since
408         !                                                                     ! the begining of the run [s]
409      ENDIF
410
411   END SUBROUTINE sed_rst_cal
412
413END MODULE sedrst
Note: See TracBrowser for help on using the repository browser.