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

source: NEMO/trunk/src/TOP/PISCES/SED/sedrst.F90

Last change on this file was 15450, checked in by cetlod, 2 years ago

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