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.
sedrst.F90 in NEMO/branches/2020/dev_r13327_KERNEL-06_2_techene_e3/src/TOP/PISCES/SED – NEMO

source: NEMO/branches/2020/dev_r13327_KERNEL-06_2_techene_e3/src/TOP/PISCES/SED/sedrst.F90 @ 14018

Last change on this file since 14018 was 14018, checked in by techene, 3 years ago

#2385 branch updated with trunk 13970

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