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.
p4zsms.F90 in branches/UKMO/dev_r10171_test_crs_AMM7/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z – NEMO

source: branches/UKMO/dev_r10171_test_crs_AMM7/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90 @ 10207

Last change on this file since 10207 was 10207, checked in by cmao, 6 years ago

remove svn keyword
File size: 26.1 KB
Line 
1MODULE p4zsms
2   !!======================================================================
3   !!                         ***  MODULE p4zsms  ***
4   !! TOP :   PISCES Source Minus Sink manager
5   !!======================================================================
6   !! History :   1.0  !  2004-03 (O. Aumont) Original code
7   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
8   !!----------------------------------------------------------------------
9#if defined key_pisces
10   !!----------------------------------------------------------------------
11   !!   'key_pisces'                                       PISCES bio-model
12   !!----------------------------------------------------------------------
13   !!   p4zsms         :  Time loop of passive tracers sms
14   !!----------------------------------------------------------------------
15   USE oce_trc         !  shared variables between ocean and passive tracers
16   USE trcdta
17   USE sms_pisces      !  PISCES Source Minus Sink variables
18   USE p4zbio          !  Biological model
19   USE p4zche          !  Chemical model
20   USE p4zlys          !  Calcite saturation
21   USE p4zflx          !  Gas exchange
22   USE p4zsbc          !  External source of nutrients
23   USE p4zsed          !  Sedimentation
24   USE p4zint          !  time interpolation
25   USE p4zrem          !  remineralisation
26   USE trd_oce         !  Ocean trends variables
27   USE trdtrc          !  TOP trends variables
28   USE sedmodel        !  Sediment model
29   USE prtctl_trc, ONLY:prt_ctl_trc_info,prt_ctl_trc     !  print control for debugging
30
31   IMPLICIT NONE
32   PRIVATE
33
34   PUBLIC   p4z_sms_init   ! called in p4zsms.F90
35   PUBLIC   p4z_sms        ! called in p4zsms.F90
36
37   REAL(wp) :: alkbudget, no3budget, silbudget, ferbudget, po4budget
38   REAL(wp) :: xfact1, xfact2, xfact3
39   INTEGER ::  numco2, numnut, numnit  !: logical unit for co2 budget
40
41   !!* Array used to indicate negative tracer values
42   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   xnegtr     !: ???
43
44
45   !! * Substitutions
46#  include "top_substitute.h90"
47   !!----------------------------------------------------------------------
48   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
49   !! $Id$
50   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
51   !!----------------------------------------------------------------------
52
53CONTAINS
54
55   SUBROUTINE p4z_sms( kt )
56      !!---------------------------------------------------------------------
57      !!                     ***  ROUTINE p4z_sms  ***
58      !!
59      !! ** Purpose :   Managment of the call to Biological sources and sinks
60      !!              routines of PISCES bio-model
61      !!
62      !! ** Method  : - at each new day ...
63      !!              - several calls of bio and sed ???
64      !!              - ...
65      !!---------------------------------------------------------------------
66      !
67      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index     
68      !!
69      INTEGER ::   ji, jj, jk, jnt, jn, jl
70      REAL(wp) ::  ztra
71#if defined key_kriest
72      REAL(wp) ::  zcoef1, zcoef2
73#endif
74      CHARACTER (len=25) :: charout
75      !!---------------------------------------------------------------------
76      !
77      IF( nn_timing == 1 )  CALL timing_start('p4z_sms')
78      !
79      IF( kt == nittrc000 ) THEN
80        !
81        ALLOCATE( xnegtr(jpi,jpj,jpk) )
82        !
83        CALL p4z_che                              ! initialize the chemical constants
84        !
85        IF( .NOT. ln_rsttr ) THEN  ;   CALL p4z_che_ahini( hi )   !  set PH at kt=nit000
86        ELSE                       ;   CALL p4z_rst( nittrc000, 'READ' )  !* read or initialize all required fields
87        ENDIF
88        !
89      ENDIF
90      !
91      IF( ln_pisdmp .AND. MOD( kt - nn_dttrc, nn_pisdmp ) == 0 )   CALL p4z_dmp( kt )      ! Relaxation of some tracers
92      !
93      !                                                                    !   set time step size (Euler/Leapfrog)
94      IF( ( neuler == 0 .AND. kt == nittrc000 ) .OR. ln_top_euler ) THEN   ;    rfact = rdttrc(1)     !  at nittrc000
95      ELSEIF( kt <= nittrc000 + nn_dttrc )                          THEN   ;    rfact = 2. * rdttrc(1)   ! at nittrc000 or nittrc000+nn_dttrc (Leapfrog)
96      ENDIF
97      !
98      IF( ( ln_top_euler .AND. kt == nittrc000 )  .OR. ( .NOT.ln_top_euler .AND. kt <= nittrc000 + nn_dttrc ) ) THEN
99         rfactr  = 1. / rfact
100         rfact2  = rfact / FLOAT( nrdttrc )
101         rfact2r = 1. / rfact2
102         xstep = rfact2 / rday         ! Time step duration for biology
103         IF(lwp) WRITE(numout,*) 
104         IF(lwp) WRITE(numout,*) '    Passive Tracer  time step    rfact  = ', rfact, ' rdt = ', rdttra(1)
105         IF(lwp) write(numout,*) '    PISCES  Biology time step    rfact2 = ', rfact2
106         IF(lwp) WRITE(numout,*)
107      ENDIF
108
109      IF( ( neuler == 0 .AND. kt == nittrc000 ) .OR. ln_top_euler ) THEN
110         DO jn = jp_pcs0, jp_pcs1              !   SMS on tracer without Asselin time-filter
111            trb(:,:,:,jn) = trn(:,:,:,jn)
112         END DO
113      ENDIF
114      !
115      IF( ndayflxtr /= nday_year ) THEN      ! New days
116         !
117         ndayflxtr = nday_year
118
119         IF(lwp) write(numout,*)
120         IF(lwp) write(numout,*) ' New chemical constants and various rates for biogeochemistry at new day : ', nday_year
121         IF(lwp) write(numout,*) '~~~~~~'
122
123         CALL p4z_che              ! computation of chemical constants
124         CALL p4z_int( kt )        ! computation of various rates for biogeochemistry
125         !
126      ENDIF
127
128      IF( ll_sbc ) CALL p4z_sbc( kt )   ! external sources of nutrients
129
130      DO jnt = 1, nrdttrc          ! Potential time splitting if requested
131         !
132         CALL p4z_bio( kt, jnt )   ! Biology
133         CALL p4z_lys( kt, jnt )   ! Compute CaCO3 saturation
134         CALL p4z_sed( kt, jnt )   ! Surface and Bottom boundary conditions
135         CALL p4z_flx( kt, jnt )   ! Compute surface fluxes
136         !
137         xnegtr(:,:,:) = 1.e0
138         DO jn = jp_pcs0, jp_pcs1
139            DO jk = 1, jpk
140               DO jj = 1, jpj
141                  DO ji = 1, jpi
142                     IF( ( trb(ji,jj,jk,jn) + tra(ji,jj,jk,jn) ) < 0.e0 ) THEN
143                        ztra             = ABS( trb(ji,jj,jk,jn) ) / ( ABS( tra(ji,jj,jk,jn) ) + rtrn )
144                        xnegtr(ji,jj,jk) = MIN( xnegtr(ji,jj,jk),  ztra )
145                     ENDIF
146                 END DO
147               END DO
148            END DO
149         END DO
150         !                                ! where at least 1 tracer concentration becomes negative
151         !                                !
152         DO jn = jp_pcs0, jp_pcs1
153           trb(:,:,:,jn) = trb(:,:,:,jn) + xnegtr(:,:,:) * tra(:,:,:,jn)
154         END DO
155        !
156         DO jn = jp_pcs0, jp_pcs1
157            tra(:,:,:,jn) = 0._wp
158         END DO
159         !
160         IF( ln_top_euler ) THEN
161            DO jn = jp_pcs0, jp_pcs1
162               trn(:,:,:,jn) = trb(:,:,:,jn)
163            END DO
164         ENDIF
165      END DO
166
167#if defined key_kriest
168      !
169      zcoef1 = 1.e0 / xkr_massp 
170      zcoef2 = 1.e0 / xkr_massp / 1.1
171      DO jk = 1,jpkm1
172         trb(:,:,jk,jpnum) = MAX(  trb(:,:,jk,jpnum), trb(:,:,jk,jppoc) * zcoef1 / xnumm(jk)  )
173         trb(:,:,jk,jpnum) = MIN(  trb(:,:,jk,jpnum), trb(:,:,jk,jppoc) * zcoef2              )
174      END DO
175      !
176#endif
177      !
178      !
179      IF( l_trdtrc ) THEN
180         DO jn = jp_pcs0, jp_pcs1
181           CALL trd_trc( tra(:,:,:,jn), jn, jptra_sms, kt )   ! save trends
182         END DO
183      END IF
184      !
185      IF( lk_sed ) THEN 
186         !
187         CALL sed_model( kt )     !  Main program of Sediment model
188         !
189         DO jn = jp_pcs0, jp_pcs1
190           CALL lbc_lnk( trb(:,:,:,jn), 'T', 1. )
191         END DO
192         !
193      ENDIF
194      !
195      IF( lrst_trc )  CALL p4z_rst( kt, 'WRITE' )  !* Write PISCES informations in restart file
196      !
197
198      IF( lk_iomput .OR. ln_check_mass )  CALL p4z_chk_mass( kt ) ! Mass conservation checking
199
200      IF ( lwm .AND. kt == nittrc000 ) CALL FLUSH    ( numonp )     ! flush output namelist PISCES
201      IF( nn_timing == 1 )  CALL timing_stop('p4z_sms')
202      !
203      !
204   END SUBROUTINE p4z_sms
205
206   SUBROUTINE p4z_sms_init
207      !!----------------------------------------------------------------------
208      !!                     ***  p4z_sms_init  *** 
209      !!
210      !! ** Purpose :   read PISCES namelist
211      !!
212      !! ** input   :   file 'namelist.trc.s' containing the following
213      !!             namelist: natext, natbio, natsms
214      !!                       natkriest ("key_kriest")
215      !!----------------------------------------------------------------------
216      NAMELIST/nampisbio/ nrdttrc, wsbio, xkmort, ferat3, wsbio2, niter1max, niter2max
217#if defined key_kriest
218      NAMELIST/nampiskrp/ xkr_eta, xkr_zeta, xkr_ncontent, xkr_mass_min, xkr_mass_max
219#endif
220      NAMELIST/nampisdmp/ ln_pisdmp, nn_pisdmp
221      NAMELIST/nampismass/ ln_check_mass
222      INTEGER :: ios                 ! Local integer output status for namelist read
223      !!----------------------------------------------------------------------
224
225      REWIND( numnatp_ref )              ! Namelist nampisbio in reference namelist : Pisces variables
226      READ  ( numnatp_ref, nampisbio, IOSTAT = ios, ERR = 901)
227901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisbio in reference namelist', lwp )
228
229      REWIND( numnatp_cfg )              ! Namelist nampisbio in configuration namelist : Pisces variables
230      READ  ( numnatp_cfg, nampisbio, IOSTAT = ios, ERR = 902 )
231902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisbio in configuration namelist', lwp )
232      IF(lwm) WRITE ( numonp, nampisbio )
233
234      IF(lwp) THEN                         ! control print
235         WRITE(numout,*) ' Namelist : nampisbio'
236         WRITE(numout,*) '    frequence pour la biologie                nrdttrc   =', nrdttrc
237         WRITE(numout,*) '    POC sinking speed                         wsbio     =', wsbio
238         WRITE(numout,*) '    half saturation constant for mortality    xkmort    =', xkmort
239         WRITE(numout,*) '    Fe/C in zooplankton                       ferat3    =', ferat3
240         WRITE(numout,*) '    Big particles sinking speed               wsbio2    =', wsbio2
241         WRITE(numout,*) '    Maximum number of iterations for POC      niter1max =', niter1max
242         WRITE(numout,*) '    Maximum number of iterations for GOC      niter2max =', niter2max
243      ENDIF
244
245#if defined key_kriest
246
247      !                               ! nampiskrp : kriest parameters
248      !                               ! -----------------------------
249      REWIND( numnatp_ref )              ! Namelist nampiskrp in reference namelist : Pisces Kriest
250      READ  ( numnatp_ref, nampiskrp, IOSTAT = ios, ERR = 903)
251903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampiskrp in reference namelist', lwp )
252
253      REWIND( numnatp_cfg )              ! Namelist nampiskrp in configuration namelist : Pisces Kriest
254      READ  ( numnatp_cfg, nampiskrp, IOSTAT = ios, ERR = 904 )
255904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampiskrp in configuration namelist', lwp )
256      IF(lwm) WRITE ( numonp, nampiskrp )
257
258      IF(lwp) THEN
259         WRITE(numout,*)
260         WRITE(numout,*) ' Namelist : nampiskrp'
261         WRITE(numout,*) '    Sinking  exponent                        xkr_eta      = ', xkr_eta
262         WRITE(numout,*) '    N content exponent                       xkr_zeta     = ', xkr_zeta
263         WRITE(numout,*) '    N content factor                         xkr_ncontent = ', xkr_ncontent
264         WRITE(numout,*) '    Minimum mass for Aggregates              xkr_mass_min = ', xkr_mass_min
265         WRITE(numout,*) '    Maximum mass for Aggregates              xkr_mass_max = ', xkr_mass_max
266         WRITE(numout,*)
267     ENDIF
268
269
270     ! Computation of some variables
271     xkr_massp = xkr_ncontent * 7.625 * xkr_mass_min**xkr_zeta
272
273#endif
274
275      REWIND( numnatp_ref )              ! Namelist nampisdmp in reference namelist : Pisces damping
276      READ  ( numnatp_ref, nampisdmp, IOSTAT = ios, ERR = 905)
277905   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisdmp in reference namelist', lwp )
278
279      REWIND( numnatp_cfg )              ! Namelist nampisdmp in configuration namelist : Pisces damping
280      READ  ( numnatp_cfg, nampisdmp, IOSTAT = ios, ERR = 906 )
281906   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisdmp in configuration namelist', lwp )
282      IF(lwm) WRITE ( numonp, nampisdmp )
283
284      IF(lwp) THEN                         ! control print
285         WRITE(numout,*)
286         WRITE(numout,*) ' Namelist : nampisdmp'
287         WRITE(numout,*) '    Relaxation of tracer to glodap mean value             ln_pisdmp      =', ln_pisdmp
288         WRITE(numout,*) '    Frequency of Relaxation                               nn_pisdmp      =', nn_pisdmp
289         WRITE(numout,*) ' '
290      ENDIF
291
292      REWIND( numnatp_ref )              ! Namelist nampismass in reference namelist : Pisces mass conservation check
293      READ  ( numnatp_ref, nampismass, IOSTAT = ios, ERR = 907)
294907   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampismass in reference namelist', lwp )
295
296      REWIND( numnatp_cfg )              ! Namelist nampismass in configuration namelist : Pisces mass conservation check
297      READ  ( numnatp_cfg, nampismass, IOSTAT = ios, ERR = 908 )
298908   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampismass in configuration namelist', lwp )
299      IF(lwm) WRITE ( numonp, nampismass )
300
301      IF(lwp) THEN                         ! control print
302         WRITE(numout,*) ' '
303         WRITE(numout,*) ' Namelist parameter for mass conservation checking'
304         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
305         WRITE(numout,*) '    Flag to check mass conservation of NO3/Si/TALK ln_check_mass = ', ln_check_mass
306      ENDIF
307
308   END SUBROUTINE p4z_sms_init
309
310   SUBROUTINE p4z_rst( kt, cdrw )
311      !!---------------------------------------------------------------------
312      !!                   ***  ROUTINE p4z_rst  ***
313      !!
314      !!  ** Purpose : Read or write variables in restart file:
315      !!
316      !!  WRITE(READ) mode:
317      !!       kt        : number of time step since the begining of the experiment at the
318      !!                   end of the current(previous) run
319      !!---------------------------------------------------------------------
320      INTEGER         , INTENT(in) ::   kt         ! ocean time-step
321      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag
322      !!---------------------------------------------------------------------
323
324      IF( TRIM(cdrw) == 'READ' ) THEN
325         !
326         IF(lwp) WRITE(numout,*)
327         IF(lwp) WRITE(numout,*) ' p4z_rst : Read specific variables from pisces model '
328         IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~'
329         !
330         IF( iom_varid( numrtr, 'PH', ldstop = .FALSE. ) > 0 ) THEN
331            CALL iom_get( numrtr, jpdom_autoglo, 'PH' , hi(:,:,:)  )
332         ELSE
333            CALL p4z_che_ahini( hi )
334         ENDIF
335         CALL iom_get( numrtr, jpdom_autoglo, 'Silicalim', xksi(:,:) )
336         IF( iom_varid( numrtr, 'Silicamax', ldstop = .FALSE. ) > 0 ) THEN
337            CALL iom_get( numrtr, jpdom_autoglo, 'Silicamax' , xksimax(:,:)  )
338         ELSE
339            xksimax(:,:) = xksi(:,:)
340         ENDIF
341         !
342         IF( iom_varid( numrtr, 'tcflxcum', ldstop = .FALSE. ) > 0 ) THEN  ! cumulative total flux of carbon
343            CALL iom_get( numrtr, 'tcflxcum' , t_oce_co2_flx_cum  )
344         ELSE
345            t_oce_co2_flx_cum = 0._wp
346         ENDIF
347         !
348      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN
349         IF( kt == nitrst ) THEN
350            IF(lwp) WRITE(numout,*)
351            IF(lwp) WRITE(numout,*) 'p4z_rst : write pisces restart file  kt =', kt
352            IF(lwp) WRITE(numout,*) '~~~~~~~'
353         ENDIF
354         CALL iom_rstput( kt, nitrst, numrtw, 'PH', hi(:,:,:) )
355         CALL iom_rstput( kt, nitrst, numrtw, 'Silicalim', xksi(:,:) )
356         CALL iom_rstput( kt, nitrst, numrtw, 'Silicamax', xksimax(:,:) )
357         CALL iom_rstput( kt, nitrst, numrtw, 'tcflxcum', t_oce_co2_flx_cum )
358      ENDIF
359      !
360   END SUBROUTINE p4z_rst
361
362   SUBROUTINE p4z_dmp( kt )
363      !!----------------------------------------------------------------------
364      !!                    ***  p4z_dmp  ***
365      !!
366      !! ** purpose  : Relaxation of some tracers
367      !!----------------------------------------------------------------------
368      !
369      INTEGER, INTENT( in )  ::     kt ! time step
370      !
371      REAL(wp) ::  alkmean = 2426.     ! mean value of alkalinity ( Glodap ; for Goyet 2391. )
372      REAL(wp) ::  po4mean = 2.165     ! mean value of phosphates
373      REAL(wp) ::  no3mean = 30.90     ! mean value of nitrate
374      REAL(wp) ::  silmean = 91.51     ! mean value of silicate
375      !
376      REAL(wp) :: zarea, zalksumn, zpo4sumn, zno3sumn, zsilsumn
377      REAL(wp) :: zalksumb, zpo4sumb, zno3sumb, zsilsumb
378      !!---------------------------------------------------------------------
379
380
381      IF(lwp)  WRITE(numout,*)
382      IF(lwp)  WRITE(numout,*) ' p4z_dmp : Restoring of nutrients at time-step kt = ', kt
383      IF(lwp)  WRITE(numout,*)
384
385      IF( cp_cfg == "orca" .AND. .NOT. lk_c1d ) THEN      ! ORCA configuration (not 1D) !
386         !                                                    ! --------------------------- !
387         ! set total alkalinity, phosphate, nitrate & silicate
388         zarea          = 1._wp / glob_sum( cvol(:,:,:) ) * 1e6             
389
390         zalksumn = glob_sum( trn(:,:,:,jptal) * cvol(:,:,:)  ) * zarea
391         zpo4sumn = glob_sum( trn(:,:,:,jppo4) * cvol(:,:,:)  ) * zarea * po4r
392         zno3sumn = glob_sum( trn(:,:,:,jpno3) * cvol(:,:,:)  ) * zarea * rno3
393         zsilsumn = glob_sum( trn(:,:,:,jpsil) * cvol(:,:,:)  ) * zarea
394 
395         IF(lwp) WRITE(numout,*) '       TALKN mean : ', zalksumn
396         trn(:,:,:,jptal) = trn(:,:,:,jptal) * alkmean / zalksumn
397
398         IF(lwp) WRITE(numout,*) '       PO4N  mean : ', zpo4sumn
399         trn(:,:,:,jppo4) = trn(:,:,:,jppo4) * po4mean / zpo4sumn
400
401         IF(lwp) WRITE(numout,*) '       NO3N  mean : ', zno3sumn
402         trn(:,:,:,jpno3) = trn(:,:,:,jpno3) * no3mean / zno3sumn
403
404         IF(lwp) WRITE(numout,*) '       SiO3N mean : ', zsilsumn
405         trn(:,:,:,jpsil) = MIN( 400.e-6,trn(:,:,:,jpsil) * silmean / zsilsumn )
406         !
407         !
408         IF( .NOT. ln_top_euler ) THEN
409            zalksumb = glob_sum( trb(:,:,:,jptal) * cvol(:,:,:)  ) * zarea
410            zpo4sumb = glob_sum( trb(:,:,:,jppo4) * cvol(:,:,:)  ) * zarea * po4r
411            zno3sumb = glob_sum( trb(:,:,:,jpno3) * cvol(:,:,:)  ) * zarea * rno3
412            zsilsumb = glob_sum( trb(:,:,:,jpsil) * cvol(:,:,:)  ) * zarea
413 
414            IF(lwp) WRITE(numout,*) ' '
415            IF(lwp) WRITE(numout,*) '       TALKB mean : ', zalksumb
416            trb(:,:,:,jptal) = trb(:,:,:,jptal) * alkmean / zalksumb
417
418            IF(lwp) WRITE(numout,*) '       PO4B  mean : ', zpo4sumb
419            trb(:,:,:,jppo4) = trb(:,:,:,jppo4) * po4mean / zpo4sumb
420
421            IF(lwp) WRITE(numout,*) '       NO3B  mean : ', zno3sumb
422            trb(:,:,:,jpno3) = trb(:,:,:,jpno3) * no3mean / zno3sumb
423
424            IF(lwp) WRITE(numout,*) '       SiO3B mean : ', zsilsumb
425            trb(:,:,:,jpsil) = MIN( 400.e-6,trb(:,:,:,jpsil) * silmean / zsilsumb )
426        ENDIF
427        !
428      ENDIF
429        !
430   END SUBROUTINE p4z_dmp
431
432
433   SUBROUTINE p4z_chk_mass( kt )
434      !!----------------------------------------------------------------------
435      !!                  ***  ROUTINE p4z_chk_mass  ***
436      !!
437      !! ** Purpose :  Mass conservation check
438      !!
439      !!---------------------------------------------------------------------
440      !
441      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index     
442      REAL(wp)             ::  zrdenittot, zsdenittot, znitrpottot
443      CHARACTER(LEN=100)   ::   cltxt
444      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zvol
445      INTEGER :: jk
446      !!----------------------------------------------------------------------
447
448      !
449      !!---------------------------------------------------------------------
450
451      IF( kt == nittrc000 ) THEN
452         xfact1 = rfact2r * 12. / 1.e15 * ryyss    ! conversion molC/kt --> PgC/yr
453         xfact2 = 1.e+3 * rno3 * 14. / 1.e12 * ryyss   ! conversion molC/l/s ----> TgN/m3/yr
454         xfact3 = 1.e+3 * rfact2r * rno3   ! conversion molC/l/kt ----> molN/m3/s
455         IF( ln_check_mass .AND. lwp) THEN      !   Open budget file of NO3, ALK, Si, Fer
456            CALL ctl_opn( numco2, 'carbon.budget'  , 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
457            CALL ctl_opn( numnut, 'nutrient.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
458            CALL ctl_opn( numnit, 'nitrogen.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
459            cltxt='time-step   Alkalinity        Nitrate        Phosphorus         Silicate           Iron'
460            IF( lwp ) WRITE(numnut,*)  TRIM(cltxt)
461            IF( lwp ) WRITE(numnut,*) 
462         ENDIF
463      ENDIF
464
465      !
466      IF( iom_use( "pno3tot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN
467         !   Compute the budget of NO3, ALK, Si, Fer
468         no3budget = glob_sum( (   trn(:,:,:,jpno3) + trn(:,:,:,jpnh4)  &
469            &                    + trn(:,:,:,jpphy) + trn(:,:,:,jpdia)  &
470            &                    + trn(:,:,:,jpzoo) + trn(:,:,:,jpmes)  &
471            &                    + trn(:,:,:,jppoc)                     &
472#if ! defined key_kriest
473            &                    + trn(:,:,:,jpgoc)                     &
474#endif
475            &                    + trn(:,:,:,jpdoc)                     ) * cvol(:,:,:)  )
476         !
477         no3budget = no3budget / areatot
478         CALL iom_put( "pno3tot", no3budget )
479      ENDIF
480      !
481      IF( iom_use( "ppo4tot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN
482         po4budget = glob_sum( (   trn(:,:,:,jppo4)                     &
483            &                    + trn(:,:,:,jpphy) + trn(:,:,:,jpdia)  &
484            &                    + trn(:,:,:,jpzoo) + trn(:,:,:,jpmes)  &
485            &                    + trn(:,:,:,jppoc)                     &
486#if ! defined key_kriest
487            &                    + trn(:,:,:,jpgoc)                     &
488#endif
489            &                    + trn(:,:,:,jpdoc)                     ) * cvol(:,:,:)  )
490         po4budget = po4budget / areatot
491         CALL iom_put( "ppo4tot", po4budget )
492      ENDIF
493      !
494      IF( iom_use( "psiltot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN
495         silbudget = glob_sum( (   trn(:,:,:,jpsil) + trn(:,:,:,jpgsi)  &
496            &                    + trn(:,:,:,jpdsi)                     ) * cvol(:,:,:)  )
497         !
498         silbudget = silbudget / areatot
499         CALL iom_put( "psiltot", silbudget )
500      ENDIF
501      !
502      IF( iom_use( "palktot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN
503         alkbudget = glob_sum( (   trn(:,:,:,jpno3) * rno3              &
504            &                    + trn(:,:,:,jptal)                     &
505            &                    + trn(:,:,:,jpcal) * 2.                ) * cvol(:,:,:)  )
506         !
507         alkbudget = alkbudget / areatot
508         CALL iom_put( "palktot", alkbudget )
509      ENDIF
510      !
511      IF( iom_use( "pfertot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN
512         ferbudget = glob_sum( (   trn(:,:,:,jpfer) + trn(:,:,:,jpnfe)  &
513            &                    + trn(:,:,:,jpdfe)                     &
514#if ! defined key_kriest
515            &                    + trn(:,:,:,jpbfe)                     &
516#endif
517            &                    + trn(:,:,:,jpsfe)                     &
518            &                    + trn(:,:,:,jpzoo) * ferat3            &
519            &                    + trn(:,:,:,jpmes) * ferat3            ) * cvol(:,:,:)  )
520         !
521         ferbudget = ferbudget / areatot
522         CALL iom_put( "pfertot", ferbudget )
523      ENDIF
524      !
525
526      ! Global budget of N SMS : denitrification in the water column and in the sediment
527      !                          nitrogen fixation by the diazotrophs
528      ! --------------------------------------------------------------------------------
529      IF( iom_use( "tnfix" ) .OR.  ( ln_check_mass .AND. kt == nitend )  ) THEN
530         znitrpottot  = glob_sum ( nitrpot(:,:,:) * nitrfix * cvol(:,:,:) )
531         CALL iom_put( "tnfix"  , znitrpottot * xfact3 )  ! Global  nitrogen fixation molC/l  to molN/m3
532      ENDIF
533      !
534      IF( iom_use( "tdenit" ) .OR.  ( ln_check_mass .AND. kt == nitend )  ) THEN
535         zrdenittot = glob_sum ( denitr(:,:,:) * rdenit * xnegtr(:,:,:) * cvol(:,:,:) ) ! denitrification in the water column
536         zsdenittot = glob_sum ( sdenit(:,:) * e1e2t(:,:) * tmask(:,:,1) )              ! denitrification in the sediments
537         CALL iom_put( "tdenit"  , ( zrdenittot + zsdenittot ) * xfact3 )               ! Total denitrification in molN/m3
538      ENDIF
539
540      IF( ln_check_mass .AND. kt == nitend ) THEN   ! Compute the budget of NO3, ALK, Si, Fer
541         t_atm_co2_flx  = t_atm_co2_flx / glob_sum( e1e2t(:,:) )
542         t_oce_co2_flx  = t_oce_co2_flx         * xfact1 * (-1 )
543         tpp            = tpp           * 1000. * xfact1
544         t_oce_co2_exp  = t_oce_co2_exp * 1000. * xfact1
545         IF( lwp ) WRITE(numco2,9000) ndastp, t_atm_co2_flx, t_oce_co2_flx, tpp, t_oce_co2_exp
546         IF( lwp ) WRITE(numnut,9100) ndastp, alkbudget        * 1.e+06, &
547             &                                no3budget * rno3 * 1.e+06, &
548             &                                po4budget * po4r * 1.e+06, &
549             &                                silbudget        * 1.e+06, &
550             &                                ferbudget        * 1.e+09
551         !
552         IF( lwp ) WRITE(numnit,9200) ndastp, znitrpottot * xfact2  , &
553         &                             zrdenittot  * xfact2  , &
554         &                             zsdenittot  * xfact2
555
556      ENDIF
557      !
558 9000  FORMAT(i8,f10.5,e18.10,f10.5,f10.5)
559 9100  FORMAT(i8,5e18.10)
560 9200  FORMAT(i8,3f10.5)
561
562       !
563   END SUBROUTINE p4z_chk_mass
564
565#else
566   !!======================================================================
567   !!  Dummy module :                                   No PISCES bio-model
568   !!======================================================================
569CONTAINS
570   SUBROUTINE p4z_sms( kt )                   ! Empty routine
571      INTEGER, INTENT( in ) ::   kt
572      WRITE(*,*) 'p4z_sms: You should not have seen this print! error?', kt
573   END SUBROUTINE p4z_sms
574#endif 
575
576   !!======================================================================
577END MODULE p4zsms 
Note: See TracBrowser for help on using the repository browser.