source: NEMO/branches/2020/dev_r12563_ASINTER-06_ABL_improvement/src/TOP/PISCES/P4Z/p4zsms.F90 @ 13159

Last change on this file since 13159 was 13159, checked in by gsamson, 3 months ago

merge trunk@r13136 into ASINTER-06 branch; pass all SETTE tests; results identical to trunk@r13136; ticket #2419

  • Property svn:keywords set to Id
File size: 27.2 KB
RevLine 
[3443]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   !!----------------------------------------------------------------------
[9169]9   !!   p4z_sms        : Time loop of passive tracers sms
[3443]10   !!----------------------------------------------------------------------
[9169]11   USE oce_trc         ! shared variables between ocean and passive tracers
12   USE trc             ! passive tracers common variables
13   USE trcdta          !
14   USE sms_pisces      ! PISCES Source Minus Sink variables
15   USE p4zbio          ! Biological model
16   USE p4zche          ! Chemical model
17   USE p4zlys          ! Calcite saturation
18   USE p4zflx          ! Gas exchange
[12377]19   USE p4zbc           ! External source of nutrients
[9169]20   USE p4zsed          ! Sedimentation
21   USE p4zint          ! time interpolation
22   USE p4zrem          ! remineralisation
23   USE iom             ! I/O manager
24   USE trd_oce         ! Ocean trends variables
25   USE trdtrc          ! TOP trends variables
26   USE sedmodel        ! Sediment model
27   USE prtctl_trc      ! print control for debugging
[3443]28
29   IMPLICIT NONE
30   PRIVATE
31
[4990]32   PUBLIC   p4z_sms_init   ! called in p4zsms.F90
33   PUBLIC   p4z_sms        ! called in p4zsms.F90
[3443]34
[9169]35   INTEGER ::    numco2, numnut, numnit      ! logical unit for co2 budget
36   REAL(wp) ::   alkbudget, no3budget, silbudget, ferbudget, po4budget
[12276]37   REAL(wp) ::   xfact, xfact1, xfact2, xfact3
[3443]38
[9169]39   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   xnegtr     ! Array used to indicate negative tracer values
[5385]40
[12377]41   !! * Substitutions
42#  include "do_loop_substitute.h90"
[3443]43   !!----------------------------------------------------------------------
[10067]44   !! NEMO/TOP 4.0 , NEMO Consortium (2018)
[10069]45   !! $Id$
[10068]46   !! Software governed by the CeCILL license (see ./LICENSE)
[3443]47   !!----------------------------------------------------------------------
48CONTAINS
49
[12377]50   SUBROUTINE p4z_sms( kt, Kbb, Kmm, Krhs )
[3443]51      !!---------------------------------------------------------------------
52      !!                     ***  ROUTINE p4z_sms  ***
53      !!
54      !! ** Purpose :   Managment of the call to Biological sources and sinks
55      !!              routines of PISCES bio-model
56      !!
57      !! ** Method  : - at each new day ...
58      !!              - several calls of bio and sed ???
59      !!              - ...
60      !!---------------------------------------------------------------------
61      !
[12377]62      INTEGER, INTENT( in ) ::   kt              ! ocean time-step index     
63      INTEGER, INTENT( in ) ::   Kbb, Kmm, Krhs  ! time level index
[3443]64      !!
[5385]65      INTEGER ::   ji, jj, jk, jnt, jn, jl
66      REAL(wp) ::  ztra
[3443]67      CHARACTER (len=25) :: charout
[12276]68      REAL(wp), ALLOCATABLE, DIMENSION(:,:    ) :: zw2d
69      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:  ) :: zw3d
70      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) ::   ztrdt   ! 4D workspace
71
[3443]72      !!---------------------------------------------------------------------
73      !
[9124]74      IF( ln_timing )   CALL timing_start('p4z_sms')
[3443]75      !
[4152]76      IF( kt == nittrc000 ) THEN
77        !
[5385]78        ALLOCATE( xnegtr(jpi,jpj,jpk) )
79        !
[9559]80        IF( .NOT. ln_rsttr ) THEN
[12377]81            CALL p4z_che( Kbb, Kmm )                  ! initialize the chemical constants
82            CALL ahini_for_at( hi, Kbb )              !  set PH at kt=nit000
[10382]83            t_oce_co2_flx_cum = 0._wp
[9559]84        ELSE
[12377]85            CALL p4z_rst( nittrc000, Kbb, Kmm,  'READ' )  !* read or initialize all required fields
[4152]86        ENDIF
87        !
88      ENDIF
89      !
[12377]90      IF( ln_pisdmp .AND. MOD( kt - 1, nn_pisdmp ) == 0 )   CALL p4z_dmp( kt, Kbb, Kmm )      ! Relaxation of some tracers
[3496]91      !
[12489]92      rfact = rDt_trc
[5385]93      !
[12276]94      ! trends computation initialisation
95      IF( l_trdtrc )  THEN
96         ALLOCATE( ztrdt(jpi,jpj,jpk,jp_pisces) )  !* store now fields before applying the Asselin filter
[12377]97         ztrdt(:,:,:,:)  = tr(:,:,:,:,Kmm)
[12276]98      ENDIF
99      !
100
[12377]101      IF( ( ln_top_euler .AND. kt == nittrc000 )  .OR. ( .NOT.ln_top_euler .AND. kt <= nittrc000 + 1 ) ) THEN
[5385]102         rfactr  = 1. / rfact
[7646]103         rfact2  = rfact / REAL( nrdttrc, wp )
[5385]104         rfact2r = 1. / rfact2
105         xstep = rfact2 / rday         ! Time step duration for biology
[12276]106         xfact = 1.e+3 * rfact2r
[5385]107         IF(lwp) WRITE(numout,*) 
[12489]108         IF(lwp) WRITE(numout,*) '    Passive Tracer  time step    rfact  = ', rfact, ' rn_Dt = ', rn_Dt
[5385]109         IF(lwp) write(numout,*) '    PISCES  Biology time step    rfact2 = ', rfact2
110         IF(lwp) WRITE(numout,*)
111      ENDIF
112
[12489]113      IF( l_1st_euler .OR. ln_top_euler ) THEN
[5385]114         DO jn = jp_pcs0, jp_pcs1              !   SMS on tracer without Asselin time-filter
[12377]115            tr(:,:,:,jn,Kbb) = tr(:,:,:,jn,Kmm)
[5385]116         END DO
117      ENDIF
118      !
[12377]119      IF( ll_bc )    CALL p4z_bc( kt, Kbb, Kmm, Krhs )   ! external sources of nutrients
[10222]120      !
121#if ! defined key_sed_off
[12377]122      CALL p4z_che(     Kbb, Kmm       ) ! computation of chemical constants
123      CALL p4z_int( kt, Kbb, Kmm       ) ! computation of various rates for biogeochemistry
[9559]124      !
[3443]125      DO jnt = 1, nrdttrc          ! Potential time splitting if requested
126         !
[12377]127         CALL p4z_bio( kt, jnt, Kbb, Kmm, Krhs )   ! Biology
128         CALL p4z_lys( kt, jnt, Kbb,      Krhs )   ! Compute CaCO3 saturation
129         CALL p4z_sed( kt, jnt, Kbb, Kmm, Krhs )   ! Surface and Bottom boundary conditions
130         CALL p4z_flx( kt, jnt, Kbb, Kmm, Krhs )   ! Compute surface fluxes
[3443]131         !
[7753]132         xnegtr(:,:,:) = 1.e0
[3443]133         DO jn = jp_pcs0, jp_pcs1
[12377]134            DO_3D_11_11( 1, jpk )
135               IF( ( tr(ji,jj,jk,jn,Kbb) + tr(ji,jj,jk,jn,Krhs) ) < 0.e0 ) THEN
136                  ztra             = ABS( tr(ji,jj,jk,jn,Kbb) ) / ( ABS( tr(ji,jj,jk,jn,Krhs) ) + rtrn )
137                  xnegtr(ji,jj,jk) = MIN( xnegtr(ji,jj,jk),  ztra )
138               ENDIF
139            END_3D
[5385]140         END DO
141         !                                ! where at least 1 tracer concentration becomes negative
142         !                                !
143         DO jn = jp_pcs0, jp_pcs1
[12377]144           tr(:,:,:,jn,Kbb) = tr(:,:,:,jn,Kbb) + xnegtr(:,:,:) * tr(:,:,:,jn,Krhs)
[5385]145         END DO
146        !
[12276]147        IF(  iom_use( 'INTdtAlk' ) .OR. iom_use( 'INTdtDIC' ) .OR. iom_use( 'INTdtFer' ) .OR.  &
148          &  iom_use( 'INTdtDIN' ) .OR. iom_use( 'INTdtDIP' ) .OR. iom_use( 'INTdtSil' ) )  THEN
149          !
150          ALLOCATE( zw3d(jpi,jpj,jpk), zw2d(jpi,jpj) )
151          zw3d(:,:,jpk) = 0.
152          DO jk = 1, jpkm1
[12377]153              zw3d(:,:,jk) = xnegtr(:,:,jk) * xfact * e3t(:,:,jk,Kmm) * tmask(:,:,jk)
[12276]154          ENDDO
155          !
156          zw2d(:,:) = 0.
157          DO jk = 1, jpkm1
[12377]158             zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tr(:,:,jk,jptal,Krhs)
[12276]159          ENDDO
160          CALL iom_put( 'INTdtAlk', zw2d )
161          !
162          zw2d(:,:) = 0.
163          DO jk = 1, jpkm1
[12377]164             zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tr(:,:,jk,jpdic,Krhs)
[12276]165          ENDDO
166          CALL iom_put( 'INTdtDIC', zw2d )
167          !
168          zw2d(:,:) = 0.
169          DO jk = 1, jpkm1
[12377]170             zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * rno3 * ( tr(:,:,jk,jpno3,Krhs) + tr(:,:,jk,jpnh4,Krhs) )
[12276]171          ENDDO
172          CALL iom_put( 'INTdtDIN', zw2d )
173          !
174          zw2d(:,:) = 0.
175          DO jk = 1, jpkm1
[12377]176             zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * po4r * tr(:,:,jk,jppo4,Krhs)
[12276]177          ENDDO
178          CALL iom_put( 'INTdtDIP', zw2d )
179          !
180          zw2d(:,:) = 0.
181          DO jk = 1, jpkm1
[12377]182             zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tr(:,:,jk,jpfer,Krhs)
[12276]183          ENDDO
184          CALL iom_put( 'INTdtFer', zw2d )
185          !
186          zw2d(:,:) = 0.
187          DO jk = 1, jpkm1
[12377]188             zw2d(:,:) = zw2d(:,:) + zw3d(:,:,jk) * tr(:,:,jk,jpsil,Krhs)
[12276]189          ENDDO
190          CALL iom_put( 'INTdtSil', zw2d )
191          !
192          DEALLOCATE( zw3d, zw2d )
193        ENDIF
194        !
[5385]195         DO jn = jp_pcs0, jp_pcs1
[12377]196            tr(:,:,:,jn,Krhs) = 0._wp
[5385]197         END DO
[3443]198         !
[5385]199         IF( ln_top_euler ) THEN
200            DO jn = jp_pcs0, jp_pcs1
[12377]201               tr(:,:,:,jn,Kmm) = tr(:,:,:,jn,Kbb)
[5385]202            END DO
203         ENDIF
[3443]204      END DO
205      !
[5385]206      IF( l_trdtrc ) THEN
207         DO jn = jp_pcs0, jp_pcs1
[13159]208           ztrdt(:,:,:,jn) = ( tr(:,:,:,jn,Kbb) - ztrdt(:,:,:,jn) ) * rfactr 
[12377]209           CALL trd_trc( tr(:,:,:,jn,Krhs), jn, jptra_sms, kt, Kmm )   ! save trends
[5385]210         END DO
[12276]211         DEALLOCATE( ztrdt ) 
[5385]212      END IF
[10222]213#endif
[5385]214      !
[10222]215      IF( ln_sediment ) THEN 
[3443]216         !
[12377]217         CALL sed_model( kt, Kbb, Kmm, Krhs )     !  Main program of Sediment model
[3443]218         !
[10222]219         IF( ln_top_euler ) THEN
220            DO jn = jp_pcs0, jp_pcs1
[12377]221               tr(:,:,:,jn,Kmm) = tr(:,:,:,jn,Kbb)
[10222]222            END DO
223         ENDIF
224         !
[3443]225      ENDIF
226      !
[12377]227      IF( lrst_trc )  CALL p4z_rst( kt, Kbb, Kmm,  'WRITE' )           !* Write PISCES informations in restart file
[3443]228      !
[5385]229
[12377]230      IF( lk_iomput .OR. ln_check_mass )  CALL p4z_chk_mass( kt, Kmm ) ! Mass conservation checking
[3481]231
[12377]232      IF( lwm .AND. kt == nittrc000    )  CALL FLUSH( numonp )         ! flush output namelist PISCES
[3443]233      !
[9124]234      IF( ln_timing )  CALL timing_stop('p4z_sms')
[3481]235      !
[3443]236   END SUBROUTINE p4z_sms
237
[9124]238
[3443]239   SUBROUTINE p4z_sms_init
240      !!----------------------------------------------------------------------
241      !!                     ***  p4z_sms_init  *** 
242      !!
243      !! ** Purpose :   read PISCES namelist
244      !!
245      !! ** input   :   file 'namelist.trc.s' containing the following
246      !!             namelist: natext, natbio, natsms
247      !!----------------------------------------------------------------------
[9124]248      INTEGER :: ios                 ! Local integer output status for namelist read
249      !!
[7646]250      NAMELIST/nampisbio/ nrdttrc, wsbio, xkmort, ferat3, wsbio2, wsbio2max, wsbio2scale,    &
[10416]251         &                   ldocp, ldocz, lthet, no3rat3, po4rat3
[9124]252         !
[4148]253      NAMELIST/nampisdmp/ ln_pisdmp, nn_pisdmp
[3531]254      NAMELIST/nampismass/ ln_check_mass
255      !!----------------------------------------------------------------------
[9169]256      !
257      IF(lwp) THEN
258         WRITE(numout,*)
259         WRITE(numout,*) 'p4z_sms_init : PISCES initialization'
260         WRITE(numout,*) '~~~~~~~~~~~~'
261      ENDIF
[3443]262
[4147]263      READ  ( numnatp_ref, nampisbio, IOSTAT = ios, ERR = 901)
[11536]264901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nampisbio in reference namelist' )
[4147]265      READ  ( numnatp_cfg, nampisbio, IOSTAT = ios, ERR = 902 )
[11536]266902   IF( ios >  0 )   CALL ctl_nam ( ios , 'nampisbio in configuration namelist' )
[9169]267      IF(lwm) WRITE( numonp, nampisbio )
268      !
[3443]269      IF(lwp) THEN                         ! control print
[9169]270         WRITE(numout,*) '   Namelist : nampisbio'
271         WRITE(numout,*) '      frequency for the biology                 nrdttrc     =', nrdttrc
272         WRITE(numout,*) '      POC sinking speed                         wsbio       =', wsbio
273         WRITE(numout,*) '      half saturation constant for mortality    xkmort      =', xkmort 
[7646]274         IF( ln_p5z ) THEN
[9169]275            WRITE(numout,*) '      N/C in zooplankton                        no3rat3     =', no3rat3
276            WRITE(numout,*) '      P/C in zooplankton                        po4rat3     =', po4rat3
[7646]277         ENDIF
[9169]278         WRITE(numout,*) '      Fe/C in zooplankton                       ferat3      =', ferat3
279         WRITE(numout,*) '      Big particles sinking speed               wsbio2      =', wsbio2
280         WRITE(numout,*) '      Big particles maximum sinking speed       wsbio2max   =', wsbio2max
281         WRITE(numout,*) '      Big particles sinking speed length scale  wsbio2scale =', wsbio2scale
[7646]282         IF( ln_ligand ) THEN
283            IF( ln_p4z ) THEN
[9169]284               WRITE(numout,*) '      Phyto ligand production per unit doc           ldocp  =', ldocp
285               WRITE(numout,*) '      Zoo ligand production per unit doc             ldocz  =', ldocz
286               WRITE(numout,*) '      Proportional loss of ligands due to Fe uptake  lthet  =', lthet
[7646]287            ENDIF
288         ENDIF
[3443]289      ENDIF
290
291
[4147]292      READ  ( numnatp_ref, nampisdmp, IOSTAT = ios, ERR = 905)
[11536]293905   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nampisdmp in reference namelist' )
[4147]294      READ  ( numnatp_cfg, nampisdmp, IOSTAT = ios, ERR = 906 )
[11536]295906   IF( ios >  0 )   CALL ctl_nam ( ios , 'nampisdmp in configuration namelist' )
[9169]296      IF(lwm) WRITE( numonp, nampisdmp )
297      !
[3443]298      IF(lwp) THEN                         ! control print
299         WRITE(numout,*)
[9169]300         WRITE(numout,*) '   Namelist : nampisdmp --- relaxation to GLODAP'
301         WRITE(numout,*) '      Relaxation of tracer to glodap mean value   ln_pisdmp =', ln_pisdmp
302         WRITE(numout,*) '      Frequency of Relaxation                     nn_pisdmp =', nn_pisdmp
[3443]303      ENDIF
304
[4147]305      READ  ( numnatp_ref, nampismass, IOSTAT = ios, ERR = 907)
[11536]306907   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nampismass in reference namelist' )
[4147]307      READ  ( numnatp_cfg, nampismass, IOSTAT = ios, ERR = 908 )
[11536]308908   IF( ios >  0 )   CALL ctl_nam ( ios , 'nampismass in configuration namelist' )
[9169]309      IF(lwm) WRITE( numonp, nampismass )
[4147]310
[3531]311      IF(lwp) THEN                         ! control print
[9169]312         WRITE(numout,*)
313         WRITE(numout,*) '   Namelist : nampismass  --- mass conservation checking'
314         WRITE(numout,*) '      Flag to check mass conservation of NO3/Si/TALK   ln_check_mass = ', ln_check_mass
[3531]315      ENDIF
[9124]316      !
[3443]317   END SUBROUTINE p4z_sms_init
318
[9124]319
[12377]320   SUBROUTINE p4z_rst( kt, Kbb, Kmm, cdrw )
[3443]321      !!---------------------------------------------------------------------
322      !!                   ***  ROUTINE p4z_rst  ***
323      !!
324      !!  ** Purpose : Read or write variables in restart file:
325      !!
326      !!  WRITE(READ) mode:
327      !!       kt        : number of time step since the begining of the experiment at the
328      !!                   end of the current(previous) run
329      !!---------------------------------------------------------------------
330      INTEGER         , INTENT(in) ::   kt         ! ocean time-step
[12377]331      INTEGER         , INTENT(in) ::   Kbb, Kmm   ! time level indices
[3443]332      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag
[9124]333      !!---------------------------------------------------------------------
[3443]334      !
335      IF( TRIM(cdrw) == 'READ' ) THEN
336         !
337         IF(lwp) WRITE(numout,*)
338         IF(lwp) WRITE(numout,*) ' p4z_rst : Read specific variables from pisces model '
339         IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~'
340         !
341         IF( iom_varid( numrtr, 'PH', ldstop = .FALSE. ) > 0 ) THEN
342            CALL iom_get( numrtr, jpdom_autoglo, 'PH' , hi(:,:,:)  )
343         ELSE
[12377]344            CALL p4z_che( Kbb, Kmm )                  ! initialize the chemical constants
345            CALL ahini_for_at( hi, Kbb )
[3443]346         ENDIF
347         CALL iom_get( numrtr, jpdom_autoglo, 'Silicalim', xksi(:,:) )
348         IF( iom_varid( numrtr, 'Silicamax', ldstop = .FALSE. ) > 0 ) THEN
349            CALL iom_get( numrtr, jpdom_autoglo, 'Silicamax' , xksimax(:,:)  )
350         ELSE
351            xksimax(:,:) = xksi(:,:)
352         ENDIF
353         !
[4996]354         IF( iom_varid( numrtr, 'tcflxcum', ldstop = .FALSE. ) > 0 ) THEN  ! cumulative total flux of carbon
355            CALL iom_get( numrtr, 'tcflxcum' , t_oce_co2_flx_cum  )
356         ELSE
357            t_oce_co2_flx_cum = 0._wp
358         ENDIF
359         !
[7646]360         IF( ln_p5z ) THEN
361            IF( iom_varid( numrtr, 'sized', ldstop = .FALSE. ) > 0 ) THEN
[9909]362               CALL iom_get( numrtr, jpdom_autoglo, 'sizep' , sizep(:,:,:)  )
363               CALL iom_get( numrtr, jpdom_autoglo, 'sizen' , sizen(:,:,:)  )
[7646]364               CALL iom_get( numrtr, jpdom_autoglo, 'sized' , sized(:,:,:)  )
365            ELSE
366               sizep(:,:,:) = 1.
367               sizen(:,:,:) = 1.
368               sized(:,:,:) = 1.
369            ENDIF
370        ENDIF
371        !
[3443]372      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN
373         IF( kt == nitrst ) THEN
374            IF(lwp) WRITE(numout,*)
375            IF(lwp) WRITE(numout,*) 'p4z_rst : write pisces restart file  kt =', kt
376            IF(lwp) WRITE(numout,*) '~~~~~~~'
377         ENDIF
378         CALL iom_rstput( kt, nitrst, numrtw, 'PH', hi(:,:,:) )
379         CALL iom_rstput( kt, nitrst, numrtw, 'Silicalim', xksi(:,:) )
380         CALL iom_rstput( kt, nitrst, numrtw, 'Silicamax', xksimax(:,:) )
[4996]381         CALL iom_rstput( kt, nitrst, numrtw, 'tcflxcum', t_oce_co2_flx_cum )
[7646]382         IF( ln_p5z ) THEN
[9909]383            CALL iom_rstput( kt, nitrst, numrtw, 'sizep', sizep(:,:,:) )
384            CALL iom_rstput( kt, nitrst, numrtw, 'sizen', sizen(:,:,:) )
[7646]385            CALL iom_rstput( kt, nitrst, numrtw, 'sized', sized(:,:,:) )
386         ENDIF
[3443]387      ENDIF
388      !
389   END SUBROUTINE p4z_rst
390
[9124]391
[12377]392   SUBROUTINE p4z_dmp( kt, Kbb, Kmm )
[3443]393      !!----------------------------------------------------------------------
394      !!                    ***  p4z_dmp  ***
395      !!
396      !! ** purpose  : Relaxation of some tracers
397      !!----------------------------------------------------------------------
398      !
[12377]399      INTEGER, INTENT( in )  ::     kt            ! time step
400      INTEGER, INTENT( in )  ::     Kbb, Kmm      ! time level indices
[3443]401      !
402      REAL(wp) ::  alkmean = 2426.     ! mean value of alkalinity ( Glodap ; for Goyet 2391. )
403      REAL(wp) ::  po4mean = 2.165     ! mean value of phosphates
404      REAL(wp) ::  no3mean = 30.90     ! mean value of nitrate
405      REAL(wp) ::  silmean = 91.51     ! mean value of silicate
406      !
[5385]407      REAL(wp) :: zarea, zalksumn, zpo4sumn, zno3sumn, zsilsumn
408      REAL(wp) :: zalksumb, zpo4sumb, zno3sumb, zsilsumb
[3443]409      !!---------------------------------------------------------------------
410
411      IF(lwp)  WRITE(numout,*)
[3557]412      IF(lwp)  WRITE(numout,*) ' p4z_dmp : Restoring of nutrients at time-step kt = ', kt
[3443]413      IF(lwp)  WRITE(numout,*)
414
[10222]415      IF( cn_cfg == "ORCA" .OR. cn_cfg == "orca") THEN
[10213]416         IF( .NOT. lk_c1d ) THEN      ! ORCA configuration (not 1D) !
417            !                                                ! --------------------------- !
418            ! set total alkalinity, phosphate, nitrate & silicate
[10425]419            zarea          = 1._wp / glob_sum( 'p4zsms', cvol(:,:,:) ) * 1e6             
[3443]420
[12377]421            zalksumn = glob_sum( 'p4zsms', tr(:,:,:,jptal,Kmm) * cvol(:,:,:)  ) * zarea
422            zpo4sumn = glob_sum( 'p4zsms', tr(:,:,:,jppo4,Kmm) * cvol(:,:,:)  ) * zarea * po4r
423            zno3sumn = glob_sum( 'p4zsms', tr(:,:,:,jpno3,Kmm) * cvol(:,:,:)  ) * zarea * rno3
424            zsilsumn = glob_sum( 'p4zsms', tr(:,:,:,jpsil,Kmm) * cvol(:,:,:)  ) * zarea
[7753]425 
[10213]426            IF(lwp) WRITE(numout,*) '       TALKN mean : ', zalksumn
[12377]427            tr(:,:,:,jptal,Kmm) = tr(:,:,:,jptal,Kmm) * alkmean / zalksumn
[3443]428
[10213]429            IF(lwp) WRITE(numout,*) '       PO4N  mean : ', zpo4sumn
[12377]430            tr(:,:,:,jppo4,Kmm) = tr(:,:,:,jppo4,Kmm) * po4mean / zpo4sumn
[3443]431
[10213]432            IF(lwp) WRITE(numout,*) '       NO3N  mean : ', zno3sumn
[12377]433            tr(:,:,:,jpno3,Kmm) = tr(:,:,:,jpno3,Kmm) * no3mean / zno3sumn
[3443]434
[10213]435            IF(lwp) WRITE(numout,*) '       SiO3N mean : ', zsilsumn
[12377]436            tr(:,:,:,jpsil,Kmm) = MIN( 400.e-6,tr(:,:,:,jpsil,Kmm) * silmean / zsilsumn )
[10213]437            !
438            !
439            IF( .NOT. ln_top_euler ) THEN
[12377]440               zalksumb = glob_sum( 'p4zsms', tr(:,:,:,jptal,Kbb) * cvol(:,:,:)  ) * zarea
441               zpo4sumb = glob_sum( 'p4zsms', tr(:,:,:,jppo4,Kbb) * cvol(:,:,:)  ) * zarea * po4r
442               zno3sumb = glob_sum( 'p4zsms', tr(:,:,:,jpno3,Kbb) * cvol(:,:,:)  ) * zarea * rno3
443               zsilsumb = glob_sum( 'p4zsms', tr(:,:,:,jpsil,Kbb) * cvol(:,:,:)  ) * zarea
[7753]444 
[10213]445               IF(lwp) WRITE(numout,*) ' '
446               IF(lwp) WRITE(numout,*) '       TALKB mean : ', zalksumb
[12377]447               tr(:,:,:,jptal,Kbb) = tr(:,:,:,jptal,Kbb) * alkmean / zalksumb
[5385]448
[10213]449               IF(lwp) WRITE(numout,*) '       PO4B  mean : ', zpo4sumb
[12377]450               tr(:,:,:,jppo4,Kbb) = tr(:,:,:,jppo4,Kbb) * po4mean / zpo4sumb
[5385]451
[10213]452               IF(lwp) WRITE(numout,*) '       NO3B  mean : ', zno3sumb
[12377]453               tr(:,:,:,jpno3,Kbb) = tr(:,:,:,jpno3,Kbb) * no3mean / zno3sumb
[5385]454
[10213]455               IF(lwp) WRITE(numout,*) '       SiO3B mean : ', zsilsumb
[12377]456               tr(:,:,:,jpsil,Kbb) = MIN( 400.e-6,tr(:,:,:,jpsil,Kbb) * silmean / zsilsumb )
[10213]457           ENDIF
[5385]458        ENDIF
459        !
[3443]460      ENDIF
[5385]461        !
[3443]462   END SUBROUTINE p4z_dmp
463
464
[12377]465   SUBROUTINE p4z_chk_mass( kt, Kmm )
[3443]466      !!----------------------------------------------------------------------
467      !!                  ***  ROUTINE p4z_chk_mass  ***
468      !!
469      !! ** Purpose :  Mass conservation check
470      !!
471      !!---------------------------------------------------------------------
[5547]472      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index     
[12377]473      INTEGER, INTENT( in ) ::   Kmm     ! time level indices
[5547]474      REAL(wp)             ::  zrdenittot, zsdenittot, znitrpottot
[5385]475      CHARACTER(LEN=100)   ::   cltxt
476      INTEGER :: jk
[9125]477      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwork
[5385]478      !!----------------------------------------------------------------------
479      !
[3443]480      IF( kt == nittrc000 ) THEN
[8533]481         xfact1 = rfact2r * 12. / 1.e15 * ryyss    ! conversion molC/kt --> PgC/yr
482         xfact2 = 1.e+3 * rno3 * 14. / 1.e12 * ryyss   ! conversion molC/l/s ----> TgN/m3/yr
483         xfact3 = 1.e+3 * rfact2r * rno3   ! conversion molC/l/kt ----> molN/m3/s
[3451]484         IF( ln_check_mass .AND. lwp) THEN      !   Open budget file of NO3, ALK, Si, Fer
[3496]485            CALL ctl_opn( numco2, 'carbon.budget'  , 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
486            CALL ctl_opn( numnut, 'nutrient.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
[5385]487            CALL ctl_opn( numnit, 'nitrogen.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
488            cltxt='time-step   Alkalinity        Nitrate        Phosphorus         Silicate           Iron'
489            IF( lwp ) WRITE(numnut,*)  TRIM(cltxt)
490            IF( lwp ) WRITE(numnut,*) 
[3443]491         ENDIF
492      ENDIF
493
[4996]494      IF( iom_use( "pno3tot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN
495         !   Compute the budget of NO3, ALK, Si, Fer
[7646]496         IF( ln_p4z ) THEN
[12377]497            zwork(:,:,:) =    tr(:,:,:,jpno3,Kmm) + tr(:,:,:,jpnh4,Kmm)                      &
498               &          +   tr(:,:,:,jpphy,Kmm) + tr(:,:,:,jpdia,Kmm)                      &
499               &          +   tr(:,:,:,jppoc,Kmm) + tr(:,:,:,jpgoc,Kmm)  + tr(:,:,:,jpdoc,Kmm)  &       
500               &          +   tr(:,:,:,jpzoo,Kmm) + tr(:,:,:,jpmes,Kmm) 
[7646]501        ELSE
[12377]502            zwork(:,:,:) =    tr(:,:,:,jpno3,Kmm) + tr(:,:,:,jpnh4,Kmm) + tr(:,:,:,jpnph,Kmm)   &
503               &          +   tr(:,:,:,jpndi,Kmm) + tr(:,:,:,jpnpi,Kmm)                      & 
504               &          +   tr(:,:,:,jppon,Kmm) + tr(:,:,:,jpgon,Kmm) + tr(:,:,:,jpdon,Kmm)   &
505               &          + ( tr(:,:,:,jpzoo,Kmm) + tr(:,:,:,jpmes,Kmm) ) * no3rat3 
[7646]506        ENDIF
507        !
[10425]508        no3budget = glob_sum( 'p4zsms', zwork(:,:,:) * cvol(:,:,:)  ) 
[7646]509        no3budget = no3budget / areatot
510        CALL iom_put( "pno3tot", no3budget )
[4996]511      ENDIF
512      !
[5385]513      IF( iom_use( "ppo4tot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN
[7646]514         IF( ln_p4z ) THEN
[12377]515            zwork(:,:,:) =    tr(:,:,:,jppo4,Kmm)                                         &
516               &          +   tr(:,:,:,jpphy,Kmm) + tr(:,:,:,jpdia,Kmm)                      &
517               &          +   tr(:,:,:,jppoc,Kmm) + tr(:,:,:,jpgoc,Kmm)  + tr(:,:,:,jpdoc,Kmm)  &       
518               &          +   tr(:,:,:,jpzoo,Kmm) + tr(:,:,:,jpmes,Kmm) 
[7646]519        ELSE
[12377]520            zwork(:,:,:) =    tr(:,:,:,jppo4,Kmm) + tr(:,:,:,jppph,Kmm)                      &
521               &          +   tr(:,:,:,jppdi,Kmm) + tr(:,:,:,jpppi,Kmm)                      & 
522               &          +   tr(:,:,:,jppop,Kmm) + tr(:,:,:,jpgop,Kmm) + tr(:,:,:,jpdop,Kmm)   &
523               &          + ( tr(:,:,:,jpzoo,Kmm) + tr(:,:,:,jpmes,Kmm) ) * po4rat3 
[7646]524        ENDIF
525        !
[10425]526        po4budget = glob_sum( 'p4zsms', zwork(:,:,:) * cvol(:,:,:)  ) 
[7646]527        po4budget = po4budget / areatot
528        CALL iom_put( "ppo4tot", po4budget )
[5385]529      ENDIF
530      !
531      IF( iom_use( "psiltot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN
[12377]532         zwork(:,:,:) =  tr(:,:,:,jpsil,Kmm) + tr(:,:,:,jpgsi,Kmm) + tr(:,:,:,jpdsi,Kmm) 
[4996]533         !
[10425]534         silbudget = glob_sum( 'p4zsms', zwork(:,:,:) * cvol(:,:,:)  ) 
[4996]535         silbudget = silbudget / areatot
536         CALL iom_put( "psiltot", silbudget )
537      ENDIF
538      !
[5385]539      IF( iom_use( "palktot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN
[12377]540         zwork(:,:,:) =  tr(:,:,:,jpno3,Kmm) * rno3 + tr(:,:,:,jptal,Kmm) + tr(:,:,:,jpcal,Kmm) * 2.             
[4996]541         !
[10425]542         alkbudget = glob_sum( 'p4zsms', zwork(:,:,:) * cvol(:,:,:)  )         !
[4996]543         alkbudget = alkbudget / areatot
544         CALL iom_put( "palktot", alkbudget )
545      ENDIF
546      !
[5385]547      IF( iom_use( "pfertot" ) .OR. ( ln_check_mass .AND. kt == nitend )  ) THEN
[12377]548         zwork(:,:,:) =   tr(:,:,:,jpfer,Kmm) + tr(:,:,:,jpnfe,Kmm) + tr(:,:,:,jpdfe,Kmm)   &
549            &         +   tr(:,:,:,jpbfe,Kmm) + tr(:,:,:,jpsfe,Kmm)                      &
550            &         + ( tr(:,:,:,jpzoo,Kmm) + tr(:,:,:,jpmes,Kmm) )  * ferat3   
[3496]551         !
[10425]552         ferbudget = glob_sum( 'p4zsms', zwork(:,:,:) * cvol(:,:,:)  ) 
[3496]553         ferbudget = ferbudget / areatot
[4996]554         CALL iom_put( "pfertot", ferbudget )
[3496]555      ENDIF
[4996]556      !
[5385]557      ! Global budget of N SMS : denitrification in the water column and in the sediment
558      !                          nitrogen fixation by the diazotrophs
559      ! --------------------------------------------------------------------------------
560      IF( iom_use( "tnfix" ) .OR.  ( ln_check_mass .AND. kt == nitend )  ) THEN
[10425]561         znitrpottot  = glob_sum ( 'p4zsms', nitrpot(:,:,:) * nitrfix * cvol(:,:,:) )
[8533]562         CALL iom_put( "tnfix"  , znitrpottot * xfact3 )  ! Global  nitrogen fixation molC/l  to molN/m3
[5385]563      ENDIF
564      !
565      IF( iom_use( "tdenit" ) .OR.  ( ln_check_mass .AND. kt == nitend )  ) THEN
[10425]566         zrdenittot = glob_sum ( 'p4zsms', denitr(:,:,:) * rdenit * xnegtr(:,:,:) * cvol(:,:,:) )
567         zsdenittot = glob_sum ( 'p4zsms', sdenit(:,:) * e1e2t(:,:) * tmask(:,:,1) )
[8533]568         CALL iom_put( "tdenit" , ( zrdenittot + zsdenittot ) * xfact3 )  ! Total denitrification molC/l to molN/m3
[5385]569      ENDIF
570      !
[4996]571      IF( ln_check_mass .AND. kt == nitend ) THEN   ! Compute the budget of NO3, ALK, Si, Fer
[10425]572         t_atm_co2_flx  = t_atm_co2_flx / glob_sum( 'p4zsms', e1e2t(:,:) )
[5385]573         t_oce_co2_flx  = t_oce_co2_flx         * xfact1 * (-1 )
574         tpp            = tpp           * 1000. * xfact1
575         t_oce_co2_exp  = t_oce_co2_exp * 1000. * xfact1
[4996]576         IF( lwp ) WRITE(numco2,9000) ndastp, t_atm_co2_flx, t_oce_co2_flx, tpp, t_oce_co2_exp
[5385]577         IF( lwp ) WRITE(numnut,9100) ndastp, alkbudget        * 1.e+06, &
[4996]578             &                                no3budget * rno3 * 1.e+06, &
[5385]579             &                                po4budget * po4r * 1.e+06, &
[4996]580             &                                silbudget        * 1.e+06, &
581             &                                ferbudget        * 1.e+09
[5385]582         !
583         IF( lwp ) WRITE(numnit,9200) ndastp, znitrpottot * xfact2  , &
[9124]584            &                             zrdenittot  * xfact2  , &
585            &                             zsdenittot  * xfact2
[4996]586      ENDIF
587      !
[3496]588 9000  FORMAT(i8,f10.5,e18.10,f10.5,f10.5)
[5385]589 9100  FORMAT(i8,5e18.10)
590 9200  FORMAT(i8,3f10.5)
[3443]591       !
592   END SUBROUTINE p4z_chk_mass
593
594   !!======================================================================
595END MODULE p4zsms 
Note: See TracBrowser for help on using the repository browser.