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/dev_001_GM/NEMO/TOP_SRC – NEMO

source: branches/dev_001_GM/NEMO/TOP_SRC/trcrst.F90 @ 766

Last change on this file since 766 was 763, checked in by gm, 16 years ago

dev_001_GM - Style only addition in TOP F90 h90 routines

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 15.7 KB
Line 
1MODULE trcrst
2   !!======================================================================
3   !!                       ***  MODULE trcrst  ***
4   !! TOP :   create, write, read the restart files for passive tracers
5   !!======================================================================
6   !! History :   1.0  !  2007-02 (C. Ethe) adaptation from the ocean
7   !!----------------------------------------------------------------------
8#if defined key_passivetrc
9   !!----------------------------------------------------------------------
10   !!   'key_passivetrc'                                    Passive tracers
11   !!----------------------------------------------------------------------
12   !!   trc_rst_opn    : open  restart file
13   !!   trc_rst_read   : read  restart file
14   !!   trc_rst_wri    : write restart file
15   !!----------------------------------------------------------------------
16   USE oce_trc
17   USE trc
18   USE sms
19   USE trctrp_lec   
20   USE lib_mpp
21   USE iom
22   
23   IMPLICIT NONE
24   PRIVATE
25   
26   PUBLIC   trc_rst_opn       ! called by ???
27   PUBLIC   trc_rst_read      ! called by ???
28   PUBLIC   trc_rst_wri       ! called by ???
29   
30   LOGICAL, PUBLIC ::   lrst_trc         !: logical to control the trc restart write
31   INTEGER, PUBLIC ::   numrtr, numrtw   !: logical unit for trc restart (read and write)
32
33   !! * Substitutions
34#  include "passivetrc_substitute.h90"
35   !!----------------------------------------------------------------------
36   !! NEMO/TOP 1.0 , LOCEAN-IPSL (2005)
37   !! $Id:$
38   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
39   !!----------------------------------------------------------------------
40   
41CONTAINS
42   
43   SUBROUTINE trc_rst_opn( kt )
44      !!----------------------------------------------------------------------
45      !!                    ***  trc_rst_opn  ***
46      !!
47      !! ** purpose  :   output of sea-trc variable in a netcdf file
48      !!----------------------------------------------------------------------
49      INTEGER, INTENT(in) ::   kt       ! number of iteration
50      !
51      CHARACTER(LEN=20)   ::   clkt     ! ocean time-step define as a character
52      CHARACTER(LEN=50)   ::   clname   ! trc output restart file name
53      !!----------------------------------------------------------------------
54      !
55      IF( kt == nit000 )  THEN
56         lrst_trc = .FALSE.
57# if defined key_off_tra
58         nitrst = nitend  ! in online version, already done in rst_opn routine defined in restart.F90 module
59# endif
60      ENDIF
61     
62      IF( kt == nitrst - ndttrc .OR. nitend - nit000 + 1 < 2 * ndttrc ) THEN
63         ! beware if model runs less than 2*ndttrc time step
64         ! beware of the format used to write kt (default is i8.8, that should be large enough)
65         IF( nitrst > 1.0e9 ) THEN   ;   WRITE(clkt,*       ) nitrst
66         ELSE                        ;   WRITE(clkt,'(i8.8)') nitrst
67         ENDIF
68         ! create the file
69         IF(lwp) WRITE(numout,*)
70         clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_restart_trc"
71         IF(lwp) WRITE(numout,*) '             open trc restart.output NetCDF file: '//clname
72         CALL iom_open( clname, numrtw, ldwrt = .TRUE., kiolib = jprstlib )
73         lrst_trc = .TRUE.
74      ENDIF
75      !
76   END SUBROUTINE trc_rst_opn
77
78
79   SUBROUTINE trc_rst_read 
80      !!----------------------------------------------------------------------
81      !!                    ***  trc_rst_opn  ***
82      !!
83      !! ** purpose  :   read passive tracer fields in restart files
84      !!----------------------------------------------------------------------
85      INTEGER  ::   ji, jj, jk, jn 
86      INTEGER  ::   iarak0   
87      REAL(wp) ::   zkt, zarak0
88      REAL(wp) ::   caralk, bicarb, co3
89      REAL(wp) ::   ztrasum
90      !!----------------------------------------------------------------------
91
92      IF(lwp) WRITE(numout,*)
93      IF(lwp) WRITE(numout,*) 'trc_rst_read : read restart file of the passive tracers'
94      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
95
96      ztrasum = 0.e0
97      IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd ) THEN   ;   iarak0 = 1
98      ELSE                                           ;   iarak0 = 0
99      ENDIF
100
101      IF(lwp) WRITE(numout,*)
102      IF(lwp) WRITE(numout,*) ' the present run starts at the time step nit000 = ', nit000
103      IF(lwp .AND. iarak0 == 1 )   WRITE(numout,*) '   and needs previous fields for Arakawa sheme ??? '
104      ENDIF
105
106      ! Time domain : restart
107      ! -------------------------
108      IF(lwp) WRITE(numout,*)
109      IF(lwp) WRITE(numout,*) ' *** passive tracer restart option'
110      SELECT CASE ( nrsttr )
111      CASE ( 0 )
112         IF(lwp) WRITE(numout,*) ' nrsttr = 0 no control of nit000'
113      CASE ( 1 )
114         IF(lwp) WRITE(numout,*) ' nrsttr = 1 we control the date of nit000'
115      CASE ( 2 )
116         IF(lwp) WRITE(numout,*) ' nrsttr = 2 the date adatrj is read in restart file'
117      CASE DEFAULT
118         IF(lwp) WRITE(numout,*) '  ===>>>> nrsttr not equal 0, 1 or 2 : no control of the date'
119         IF(lwp) WRITE(numout,*) ' =======                   ========='
120      END SELECT
121
122      CALL iom_open ( 'restart.trc', numrtr, kiolib = jprstlib )
123
124      CALL iom_get( numrtr, 'kt'   , zkt    )
125      CALL iom_get( numrtr, 'arak0', zarak0 )
126
127      IF(lwp) WRITE(numout,*)
128      IF(lwp) WRITE(numout,*) ' Info on the restart file read : '
129      IF(lwp) WRITE(numout,*) '   time-step           : ', NINT( zkt    )
130      IF(lwp) WRITE(numout,*) '   arakawa option      : ', NINT( zarak0 )
131
132
133      IF( nittrc000 - NINT( zkt ) /= 1 .AND. nrsttr /= 0 )  &      ! control of date
134         &   CALL ctl_stop( ' ===>>>> : problem with nit000 for the restart', &
135         &                  ' verify the restart file or rerun with nrstdt = 0 (namelist)' )
136
137      IF( iarak0 /= NINT( zarak0 ) )   &                           ! Control of the scheme
138         & CALL ctl_stop( ' ===>>>> : problem with advection scheme', &
139         & ' it must be the same type for both restart and previous run', &
140         & ' centered or euler '  )
141
142
143      ! READ prognostic variables and computes diagnostic variable
144      DO jn = 1, jptra
145         CALL iom_get( numrtr, jpdom_local, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) ) 
146         CALL iom_get( numrtr, jpdom_local, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) ) 
147      END DO
148# if defined key_trc_lobster1
149      CALL iom_get( numrtr, jpdom_local, 'SEDB'//ctrcnm(jpdet), sedpocb(:,:) ) 
150      CALL iom_get( numrtr, jpdom_local, 'SEDN'//ctrcnm(jpdet), sedpocn(:,:) ) 
151# elif defined key_trc_pisces
152      CALL iom_get( numrtr, jpdom_local, 'Silicalim', xksi(:,:) ) 
153      xksimax = xksi
154# elif defined key_cfc
155      DO jn = 1, jptra
156         CALL iom_get( numrtr, jpdom_local, 'qint'//ctrcnm(jn), qint(:,:,jn) ) 
157         CALL iom_get( numrtr, jpdom_local, 'qtr'//ctrcnm(jn) , qtr( :,:,jn) ) 
158      END DO
159# endif
160
161# if defined key_trc_pisces 
162      !                                                         ! --------------------------- !
163      IF( cp_cfg == "orca" .AND. .NOT. lk_trccfg_1d ) THEN      ! ORCA condiguration (not 1D) !
164         !                                                      ! --------------------------- !
165         !                                                      ! set total alkalinity, phosphate, NO3 & silicate
166         !                          ! total alkalinity
167         ztrasum = 0.e0             ! ----------------
168         DO jk = 1, jpk
169            DO jj = 1, jpj
170               DO ji = 1, jpi
171                  ztrasum = ztrasum + trn(ji,jj,jk,jptal) * tmask(ji,jj,jk) * tmask_i(ji,jj)  &
172#  if defined key_off_degrad
173                     &              * facvol(ji,jj,jk)   &
174#  endif
175                     &              * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
176               END DO
177            END DO
178         END DO
179         IF( lk_mpp )   CALL mpp_sum( ztrasum )     ! sum over the global domain 
180
181         IF(lwp) WRITE(numout,*) 'TALK moyen ', ztrasum / areatot * 1.e6
182         ztrasum = ztrasum / areatot * 1.e6
183         trn(:,:,:,jptal) = trn(:,:,:,jptal) * 2391. / ztrasum
184
185         !                          ! phosphate
186         ztrasum = 0.e0             ! ---------
187         DO jk = 1, jpk
188            DO jj = 1, jpj
189               DO ji = 1, jpi
190                  ztrasum = ztrasum + trn(ji,jj,jk,jppo4) * tmask(ji,jj,jk) * tmask_i(ji,jj)  &
191#  if defined key_off_degrad
192                     &              * facvol(ji,jj,jk)   &
193#  endif
194                     &              * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
195               END DO
196            END DO
197         END DO
198         IF( lk_mpp )   CALL mpp_sum( ztrasum )     ! sum over the global domain 
199
200         IF(lwp) WRITE(numout,*) 'PO4 moyen ', ztrasum/areatot*1E6/122.
201         ztrasum = ztrasum / areatot * 1.e6 / 122.
202         trn(:,:,:,jppo4) = trn(:,:,:,jppo4) * 2.165 / ztrasum
203
204         !                          ! NO3
205         ztrasum = 0.e0             ! ---
206         DO jk = 1, jpk
207            DO jj = 1, jpj
208               DO ji = 1, jpi
209                  ztrasum = ztrasum + trn(ji,jj,jk,jpno3) * tmask(ji,jj,jk) * tmask_i(ji,jj)  &
210#  if defined key_off_degrad
211                     &              * facvol(ji,jj,jk)   &
212#  endif
213                     &              * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
214               END DO
215            END DO
216         END DO
217         IF( lk_mpp )   CALL mpp_sum( ztrasum )     ! sum over the global domain 
218
219         IF(lwp) WRITE(numout,*) 'NO3 moyen ', ztrasum / areatot*1.e6 / 7.6
220         ztrasum = ztrasum / areatot * 1.e6 / 7.6
221         trn(:,:,:,jpno3) = trn(:,:,:,jpno3) * 30.9 / ztrasum
222
223         !                          ! Silicate
224         ztrasum = 0.e0             ! --------
225         DO jk = 1, jpk
226            DO jj = 1, jpj
227               DO ji = 1, jpi
228                  ztrasum = ztrasum + trn(ji,jj,jk,jpsil) * tmask(ji,jj,jk) * tmask_i(ji,jj)   &
229#  if defined key_off_degrad
230                     &              * facvol(ji,jj,jk)   &
231#  endif
232                     &              * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
233               END DO
234            END DO
235         END DO
236         IF( lk_mpp )   CALL mpp_sum( ztrasum )     ! sum over the global domain 
237
238         IF(lwp) WRITE(numout,*) 'SiO3 moyen ', ztrasum / areatot * 1.e6
239         ztrasum = ztrasum / areatot * 1.e6
240         trn(:,:,:,jpsil) = MIN( 400.e-6,trn(:,:,:,jpsil) * 91.51 / ztrasum ) 
241         !
242      ENDIF
243
244!#if defined key_trc_kriest
245!      !! Initialize number of particles from a standart restart file
246!      !! The name of big organic particles jpgoc has been only change
247!      !! and replace by jpnum but the values here are concentration
248!      trn(:,:,:,jppoc) = trn(:,:,:,jppoc) + trn(:,:,:,jpnum)
249!      trn(:,:,:,jpnum) = trn(:,:,:,jppoc) / ( 6. * xkr_massp )
250!#endif
251      !!  Set hi (???) from  total alcalinity, borat (???), akb3 (???) and ak23 (???)
252      !!  ---------------------------------------------------------------------
253      DO jk = 1, jpk
254         DO jj = 1, jpj
255            DO ji = 1,jpi
256               caralk = trn(ji,jj,jk,jptal) - borat(ji,jj,jk) / ( 1. + 1.e-8 / ( rtrn + akb3(ji,jj,jk) ) )
257               co3    = ( caralk - trn(ji,jj,jk,jpdic) ) *       tmask(ji,jj,jk)   &
258                  &   +             0.5e-3               * ( 1.- tmask(ji,jj,jk) )
259               bicarb = 2.* trn(ji,jj,jk,jpdic) - caralk
260               hi(ji,jj,jk) = ( ak23(ji,jj,jk) * bicarb / co3 )   *       tmask(ji,jj,jk)
261                  &         +             1.0e-9                  * ( 1.- tmask(ji,jj,jk) )
262            END DO
263         END DO
264      END DO
265# endif
266      trb(:,:,:,:) = trn(:,:,:,:)
267
268      CALL iom_close( numrtr )
269      !
270   END SUBROUTINE trc_rst_read
271
272
273   SUBROUTINE trc_rst_wri( kt )
274      !!----------------------------------------------------------------------
275      !!                    ***  trc_rst_wri  ***
276      !!
277      !! ** purpose  :   write passive tracer fields in restart files
278      !!----------------------------------------------------------------------
279      INTEGER, INTENT( in ) ::    kt
280      !!
281      INTEGER  :: ji,jj,jk,jn
282      REAL(wp) :: zdiag_var, zdiag_varmin, zdiag_varmax, zdiag_tot
283      REAL(wp) :: zder
284      !!----------------------------------------------------------------------
285
286      IF( MOD( kt, nstock ) == 0 .OR. kt == nitend ) THEN
287
288         ! 0. initialisations
289         ! ------------------
290         IF(lwp) WRITE(numout,*)
291         IF(lwp) WRITE(numout,*) 'trc_wri : write the passive tracer restart file (NetCDF) ',   &
292            &                              'at it= ',kt,' date= ',ndastp
293         IF(lwp) WRITE(numout,*) '~~~~~~~~~'
294
295
296         CALL iom_rstput( kt, nitrst, numrtw, 'kt'   ,  REAL( kt, wp )  )
297
298         IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd ) THEN
299            CALL iom_rstput( kt, nitrst, numrtw, 'arak0', 1. )
300         ELSE
301            CALL iom_rstput( kt, nitrst, numrtw, 'arak0', 0. )
302         ENDIF
303
304
305         ! prognostic variables
306         ! --------------------
307         DO jn = 1, jptra
308            CALL iom_rstput( kt, nitrst, numrtw, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) )
309            CALL iom_rstput( kt, nitrst, numrtw, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) )
310         END DO
311
312#if defined key_trc_lobster1
313         CALL iom_rstput( kt, nitrst, numrtw, 'SEDB'//ctrcnm(jpdet), sedpocb(:,:) )
314         CALL iom_rstput( kt, nitrst, numrtw, 'SEDN'//ctrcnm(jpdet), sedpocn(:,:) )
315#elif defined key_trc_pisces
316         CALL iom_rstput( kt, nitrst, numrtw, 'Silicalim', xksi(:,:) )
317
318#elif defined key_cfc
319         DO jn = 1, jptra
320            CALL iom_rstput( kt, nitrst, numrtw, 'qint'//ctrcnm(jn), qint(:,:,jn) )
321            CALL iom_rstput( kt, nitrst, numrtw, 'qtr'//ctrcnm(jn) , qtr( :,:,jn) )
322         END DO
323#endif
324
325         IF(lwp) WRITE(numout,*) '----TRACER STAT----'
326
327         zdiag_tot = 0.e0
328         DO jn = 1, jptra
329            zdiag_var    = 0.e0
330            zdiag_varmin = 0.e0
331            zdiag_varmax = 0.e0
332            DO ji = 1, jpi
333               DO jj = 1, jpj
334                  DO jk = 1,jpk
335                     zdiag_var = zdiag_var + trn(ji,jj,jk,jn) * tmask(ji,jj,jk) * tmask_i(ji,jj)   &
336#if defined key_off_degrad
337                        &   * facvol(ji,jj,jk)   &
338#endif
339                        &   * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
340                  END DO
341               END DO
342            END DO
343
344            zdiag_varmin = MINVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) )
345            zdiag_varmax = MAXVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) )
346
347            IF( lk_mpp ) THEN
348               CALL mpp_min( zdiag_varmin )      ! min over the global domain
349               CALL mpp_max( zdiag_varmax )      ! max over the global domain
350               CALL mpp_sum( zdiag_var    )      ! sum over the global domain
351            END IF
352
353            zdiag_tot = zdiag_tot + zdiag_var
354            zdiag_var = zdiag_var / areatot
355
356            IF(lwp) WRITE(numout,*) 'MEAN NO ', jn, ctrcnm(jn), ' =', zdiag_var,   &
357               &                    'MIN= ', zdiag_varmin, 'MAX= ', zdiag_varmax
358         END DO
359
360         zdiag_tot=zdiag_tot
361         zder=((zdiag_tot-trai)/trai)*100._wp
362         IF (lwp) WRITE(numout,*) 'Integral of all tracers over the full domain  =',zdiag_tot
363         IF (lwp) WRITE(numout,*) 'Drift of the sum of all tracers =',zder, '%'
364
365         CALL iom_close(numrtw)
366         !
367      ENDIF
368      !
369   END SUBROUTINE trc_rst_wri
370
371#else
372   !!----------------------------------------------------------------------
373   !!  Dummy module :                                     No passive tracer
374   !!----------------------------------------------------------------------
375CONTAINS
376   SUBROUTINE trc_rst_read                      ! Empty routines
377   END SUBROUTINE trc_rst_read
378   SUBROUTINE trc_rst_wri( kt )
379      INTEGER, INTENT ( in ) :: kt
380      WRITE(*,*) 'trc_rst_wri: You should not have seen this print! error?', kt
381   END SUBROUTINE trc_rst_wri   
382#endif
383
384   !!======================================================================
385END MODULE trcrst
Note: See TracBrowser for help on using the repository browser.