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/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/SED – NEMO

source: NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/SED/sedrst.F90 @ 15075

Last change on this file since 15075 was 15075, checked in by aumont, 3 years ago

major update of the sediment module

  • Property svn:keywords set to Id
File size: 15.5 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
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 = nn_stocklist( nrst_lst )
52            ELSE
53               nitrst = nitend
54            ENDIF
55         ENDIF
56         IF( .NOT. ln_rst_list .AND. MOD( kt - 1, nn_stock ) == 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 + nn_stock - 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      IF( .NOT. ln_rst_list .AND. nn_stock == -1 )   RETURN   ! we will never do any restart
66      ! to get better performances with NetCDF format:
67      ! we open and define the tracer restart file one tracer time step before writing the data (-> at nitrst - 2*nn_dttrc + 1)
68      ! 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
69      IF( kt == nitrst - 1 .OR. nn_stock == 1 .OR. ( kt == nitend - 1 .AND. .NOT. lrst_sed ) ) THEN
70         ! beware of the format used to write kt (default is i8.8, that should be large enough)
71         IF( nitrst > 1.0e9 ) THEN   ;   WRITE(clkt,*       ) nitrst
72         ELSE                        ;   WRITE(clkt,'(i8.8)') nitrst
73         ENDIF
74         ! create the file
75         IF(lwp) WRITE(numsed,*)
76         clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_"//TRIM(cn_sedrst_out)
77         clpath = TRIM(cn_sedrst_outdir)
78         IF( clpath(LEN_TRIM(clpath):) /= '/' ) clpath = TRIM(clpath) // '/'
79         IF(lwp) WRITE(numsed,*) &
80             '             open sed restart.output NetCDF file: ',TRIM(clpath)//clname
81         CALL iom_open( TRIM(clpath)//TRIM(clname), numrsw, ldwrt = .TRUE., kdlev = jpksed )
82         lrst_sed = .TRUE.
83      ENDIF
84      !
85   END SUBROUTINE sed_rst_opn
86
87   SUBROUTINE sed_rst_read 
88      !!----------------------------------------------------------------------
89      !!                   ***  ROUTINE sed_rst_read  ***
90      !!
91      !! ** Purpose :  Initialization of sediment module
92      !!               - sets initial sediment composition
93      !!                 ( only clay or reading restart file )
94      !!
95      !!   History :
96      !!        !  06-07  (C. Ethe)  original
97      !!----------------------------------------------------------------------
98
99      !! * local declarations
100      INTEGER :: ji, jj, jk, jn 
101      REAL(wp), DIMENSION(jpi,jpj,jpksed,jptrased) :: zdta
102      REAL(wp), DIMENSION(jpi,jpj,jpksed,2)        :: zdta1 
103      REAL(wp), DIMENSION(jpi,jpj,jpksed)          :: zdta2
104      REAL(wp), DIMENSION(jpoce,jpksed)            :: zhipor
105      REAL(wp) :: zkt
106      CHARACTER(len = 20) ::   cltra
107      CHARACTER(LEN=20)   ::   name1
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_auto, 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_auto, 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_auto, 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_auto, 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      IF( ln_timing )  CALL timing_stop('sed_rst_read')
189     
190   END SUBROUTINE sed_rst_read
191
192   SUBROUTINE sed_rst_wri( kt )
193      !!----------------------------------------------------------------------
194      !!                   ***  ROUTINE sed_rst_wri  ***
195      !!
196      !! ** Purpose :  save field which are necessary for sediment restart
197      !!
198      !!   History :
199      !!        !  06-07  (C. Ethe)  original
200      !!----------------------------------------------------------------------
201      !!* Modules used
202      INTEGER, INTENT(in) ::   kt       ! number of iteration
203      !! * local declarations
204      INTEGER  :: ji, jj, jk, jn
205      REAL(wp), DIMENSION(1) ::  zinfo
206      CHARACTER(len=50) :: clname
207      CHARACTER(len=20) :: cltra, name1 
208      REAL(wp), DIMENSION(jpoce,jpksed)   :: zdta   
209      REAL(wp), DIMENSION(jpi,jpj,jpksed) :: zdta2
210      !! -----------------------------------------------------------------------
211
212      IF( ln_timing )  CALL timing_start('sed_rst_wri')
213
214         !! 0. initialisations
215         !! ------------------
216         
217      IF(lwp) WRITE(numsed,*) ' '
218      IF(lwp) WRITE(numsed,*) 'sed_rst_write : write the sediment restart file in NetCDF format ',   &
219            'at it= ',kt
220      IF(lwp) WRITE(numsed,*) '~~~~~~~~~'
221
222      zdta(:,:)          = 1.0
223      zdta2(:,:,:)       = 0.0
224         
225      !! 1. WRITE in nutwrs
226      !! ------------------
227      zinfo(1) = REAL( kt)
228      CALL iom_rstput( kt, nitrst, numrsw, 'kt', zinfo  )
229
230      ! Back to 2D geometry
231      DO jn = 1, jpsol
232         CALL unpack_arr( jpoce, zdta2(1:jpi,1:jpj,1:jpksed) , iarroce(1:jpoce), &
233         &                       solcp(1:jpoce,1:jpksed,jn ) )
234         cltra = TRIM(sedtrcd(jn))
235         CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), zdta2(:,:,:) )
236      END DO
237
238      DO jn = 1, jpwat
239         CALL unpack_arr( jpoce, zdta2(1:jpi,1:jpj,1:jpksed) , iarroce(1:jpoce), &
240         &                       pwcp(1:jpoce,1:jpksed,jn  )  )
241         cltra = TRIM(sedtrcd(jpsol+jn))
242         CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), zdta2(:,:,:) )
243      END DO
244      ! pH
245      DO jk = 1, jpksed
246         DO ji = 1, jpoce
247            zdta(ji,jk) = -LOG10( hipor(ji,jk) / ( densSW(ji) + rtrn ) + rtrn )
248         ENDDO
249      ENDDO
250
251      CALL unpack_arr( jpoce, zdta2(1:jpi,1:jpj,1:jpksed)  , iarroce(1:jpoce), &
252      &                   zdta(1:jpoce,1:jpksed)  )
253      cltra = TRIM(seddia3d(1))
254      CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), zdta2(:,:,:) )
255         
256      CALL unpack_arr( jpoce, zdta2(1:jpi,1:jpj,1:jpksed)  , iarroce(1:jpoce), &
257      &                   co3por(1:jpoce,1:jpksed)  )
258      cltra = TRIM(seddia3d(2))
259      CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), zdta2(:,:,:) )
260
261      ! prognostic variables
262      ! --------------------
263
264      CALL unpack_arr( jpoce, zdta2(1:jpi,1:jpj,1:jpksed)  , iarroce(1:jpoce), &
265      &                   db(1:jpoce,1:jpksed)  )
266
267      cltra = "dbioturb"
268      CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), zdta2(:,:,:) )
269
270      CALL unpack_arr( jpoce, zdta2(1:jpi,1:jpj,1:jpksed)  , iarroce(1:jpoce), &
271      &                   irrig(1:jpoce,1:jpksed)  )
272
273      cltra = "irrig"
274      CALL iom_rstput( kt, nitrst, numrsw, TRIM(cltra), zdta2(:,:,:) )
275
276      IF( kt == nitrst ) THEN
277          CALL iom_close( numrsw )     ! close the restart file (only at last time step)
278          IF( l_offline .AND. ln_rst_list ) THEN
279             nrst_lst = nrst_lst + 1
280             nitrst = nn_stocklist( nrst_lst )
281          ENDIF
282      ENDIF
283
284      IF( ln_timing )  CALL timing_stop('sed_rst_wri')
285         
286   END SUBROUTINE sed_rst_wri
287
288
289   SUBROUTINE sed_rst_cal( kt, cdrw )
290      !!---------------------------------------------------------------------
291      !!                   ***  ROUTINE sed_rst_cal  ***
292      !!
293      !!  ** Purpose : Read or write calendar in restart file:
294      !!
295      !!  WRITE(READ) mode:
296      !!       kt        : number of time step since the begining of the experiment at the
297      !!                   end of the current(previous) run
298      !!       adatrj(0) : number of elapsed days since the begining of the experiment at the
299      !!                   end of the current(previous) run (REAL -> keep fractions of day)
300      !!       ndastp    : date at the end of the current(previous) run (coded as yyyymmdd integer)
301      !!
302      !!   According to namelist parameter nrstdt,
303      !!       nn_rsttr = 0  no control on the date (nittrc000 is  arbitrary).
304      !!       nn_rsttr = 1  we verify that nittrc000 is equal to the last
305      !!                   time step of previous run + 1.
306      !!       In both those options, the  exact duration of the experiment
307      !!       since the beginning (cumulated duration of all previous restart runs)
308      !!       is not stored in the restart and is assumed to be (nittrc000-1)*rdt.
309      !!       This is valid is the time step has remained constant.
310      !!
311      !!       nn_rsttr = 2  the duration of the experiment in days (adatrj)
312      !!                    has been stored in the restart file.
313      !!----------------------------------------------------------------------
314      INTEGER         , INTENT(in) ::   kt         ! ocean time-step
315      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag
316      !
317      LOGICAL  ::  llok
318      REAL(wp) ::  zkt, zrdttrc1
319      REAL(wp) ::  zndastp
320
321      ! Time domain : restart
322      ! ---------------------
323
324      IF( TRIM(cdrw) == 'READ' ) THEN
325
326         IF(lwp) WRITE(numsed,*)
327         IF(lwp) WRITE(numsed,*) 'sed_rst_cal : read the SED restart file for calendar'
328         IF(lwp) WRITE(numsed,*) '~~~~~~~~~~~~'
329
330         IF( ln_rst_sed ) THEN
331            CALL iom_open( TRIM(cn_sedrst_indir)//'/'//cn_sedrst_in, numrsr )
332            CALL iom_get ( numrsr, 'kt', zkt )   ! last time-step of previous run
333
334            IF(lwp) THEN
335               WRITE(numsed,*) ' *** Info read in restart : '
336               WRITE(numsed,*) '   previous time-step                               : ', NINT( zkt )
337               WRITE(numsed,*) ' *** restart option'
338               SELECT CASE ( nn_rstsed )
339               CASE ( 0 )   ;   WRITE(numsed,*) ' nn_rstsed = 0 : no control of nittrc000'
340               CASE ( 1 )   ;   WRITE(numsed,*) ' nn_rstsed = 1 : no control the date at nittrc000 (use ndate0 read in the namelist)'
341               CASE ( 2 )   ;   WRITE(numsed,*) ' nn_rstsed = 2 : calendar parameters read in restart'
342               END SELECT
343               WRITE(numsed,*)
344            ENDIF
345            ! Control of date
346            IF( nittrc000  - NINT( zkt ) /= 1 .AND.  nn_rstsed /= 0 )                                  &
347               &   CALL ctl_stop( ' ===>>>> : problem with nittrc000 for the restart',                 &
348               &                  ' verify the restart file or rerun with nn_rsttr = 0 (namelist)' )
349         ENDIF
350         !
351         IF( l_offline ) THEN
352            !                                          ! set the date in offline mode
353            IF( ln_rst_sed .AND. nn_rstsed == 2 ) THEN
354               CALL iom_get( numrsr, 'ndastp', zndastp )
355               ndastp = NINT( zndastp )
356               CALL iom_get( numrsr, 'adatrj', adatrj  )
357             ELSE
358               ndastp = ndate0 - 1     ! ndate0 read in the namelist in dom_nam
359               adatrj = ( REAL( nittrc000-1, wp ) * rdt ) / rday
360               ! note this is wrong if time step has changed during run
361            ENDIF
362            !
363            IF(lwp) THEN
364              WRITE(numsed,*) ' *** Info used values : '
365              WRITE(numsed,*) '   date ndastp                                      : ', ndastp
366              WRITE(numsed,*) '   number of elapsed days since the begining of run : ', adatrj
367              WRITE(numsed,*)
368            ENDIF
369            !
370            CALL day_init          ! compute calendar
371            !
372         ENDIF
373         !
374      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN
375         !
376         IF(  kt == nitrst ) THEN
377            IF(lwp) WRITE(numsed,*)
378            IF(lwp) WRITE(numsed,*) 'trc_wri : write the TOP restart file (NetCDF) at it= ', kt, ' date= ', ndastp
379            IF(lwp) WRITE(numsed,*) '~~~~~~~'
380         ENDIF
381         CALL iom_rstput( kt, nitrst, numrsw, 'kt'     , REAL( kt    , wp) )   ! time-step
382         CALL iom_rstput( kt, nitrst, numrsw, 'ndastp' , REAL( ndastp, wp) )   ! date
383         CALL iom_rstput( kt, nitrst, numrsw, 'adatrj' , adatrj            )   ! number of elapsed days since
384         !                                                                     ! the begining of the run [s]
385      ENDIF
386
387   END SUBROUTINE sed_rst_cal
388
389END MODULE sedrst
Note: See TracBrowser for help on using the repository browser.