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.
p4zrem.F90 in branches/DEV_r1879_FCM/NEMOGCM/NEMO/TOP_SRC/PISCES – NEMO

source: branches/DEV_r1879_FCM/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zrem.F90 @ 2013

Last change on this file since 2013 was 2013, checked in by smasson, 14 years ago

remove propertie svn:executabe of fortran files in DEV_r1879_FCM

  • Property svn:keywords set to Id
File size: 18.0 KB
Line 
1MODULE p4zrem
2   !!======================================================================
3   !!                         ***  MODULE p4zrem  ***
4   !! TOP :   PISCES Compute remineralization/scavenging of organic compounds
5   !!======================================================================
6   !! History :   1.0  !  2004     (O. Aumont) Original code
7   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
8   !!----------------------------------------------------------------------
9#if defined key_pisces
10   !!----------------------------------------------------------------------
11   !!   'key_top'       and                                      TOP models
12   !!   'key_pisces'                                       PISCES bio-model
13   !!----------------------------------------------------------------------
14   !!   p4z_rem       :   Compute remineralization/scavenging of organic compounds
15   !!----------------------------------------------------------------------
16   USE trc
17   USE oce_trc         !
18   USE sms_pisces      !
19   USE prtctl_trc
20   USE p4zint
21   USE p4zopt
22   USE p4zmeso
23   USE p4zprod
24   USE p4zche
25
26   IMPLICIT NONE
27   PRIVATE
28
29   PUBLIC   p4z_rem    ! called in p4zbio.F90
30
31   !! * Shared module variables
32   REAL(wp), PUBLIC ::   &
33     xremik  = 0.3_wp      ,  & !:
34     xremip  = 0.025_wp    ,  & !:
35     nitrif  = 0.05_wp     ,  & !:
36     xsirem  = 0.015_wp    ,  & !:
37     xlam1   = 0.005_wp    ,  & !:
38     oxymin  = 1.e-6_wp         !:
39
40   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::    & !:
41     &                   denitr                     !: denitrification array
42
43   REAL(wp) ::   &
44     xstep            !: Time step duration for biology
45
46   !!* Substitution
47#  include "top_substitute.h90"
48   !!----------------------------------------------------------------------
49   !! NEMO/TOP 2.0 , LOCEAN-IPSL (2007)
50   !! $Id$
51   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
52   !!----------------------------------------------------------------------
53
54CONTAINS
55
56   SUBROUTINE p4z_rem(kt, jnt)
57      !!---------------------------------------------------------------------
58      !!                     ***  ROUTINE p4z_rem  ***
59      !!
60      !! ** Purpose :   Compute remineralization/scavenging of organic compounds
61      !!
62      !! ** Method  : - ???
63      !!---------------------------------------------------------------------
64      INTEGER, INTENT(in) ::   kt, jnt ! ocean time step
65      INTEGER  ::   ji, jj, jk
66      REAL(wp) ::   zremip, zremik , zlam1b
67      REAL(wp) ::   zkeq  , zfeequi, zsiremin
68      REAL(wp) ::   zsatur, zsatur2, znusil
69      REAL(wp) ::   zbactfer, zorem, zorem2, zofer
70      REAL(wp) ::   zosil, zdenom1, zscave, zaggdfe
71#if ! defined key_kriest
72      REAL(wp) ::   zofer2, zdenom, zdenom2
73#endif
74      REAL(wp) ::   zlamfac, zonitr
75      REAL(wp), DIMENSION(jpi,jpj)     ::   ztempbac
76      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zdepbac, zfesatur, zolimi
77      CHARACTER (len=25) :: charout
78
79      !!---------------------------------------------------------------------
80
81
82      IF( ( kt * jnt ) == nittrc000  )  THEN
83         CALL p4z_rem_init                ! Initialization (first time-step only)
84         xstep = rfact2 / rday            ! Time step duration for the biology
85         nitrfac(:,:,:) = 0.0
86         denitr (:,:,:) = 0.0 
87      ENDIF
88
89
90       ! Initialisation of temprary arrys
91       zdepbac (:,:,:) = 0.0
92       zfesatur(:,:,:) = 0.0
93       zolimi  (:,:,:) = 0.0
94       ztempbac(:,:)   = 0.0
95
96!      Computation of the mean phytoplankton concentration as
97!      a crude estimate of the bacterial biomass
98!      --------------------------------------------------
99
100      DO jk = 1, jpkm1
101         DO jj = 1, jpj
102            DO ji = 1, jpi
103               IF( fsdept(ji,jj,jk) < 120. ) THEN
104                  zdepbac(ji,jj,jk) = MIN( 0.7 * ( trn(ji,jj,jk,jpzoo) + 2.* trn(ji,jj,jk,jpmes) ), 4.e-6 )
105                  ztempbac(ji,jj)   = zdepbac(ji,jj,jk)
106               ELSE
107                  zdepbac(ji,jj,jk) = MIN( 1., 120./ fsdept(ji,jj,jk) ) * ztempbac(ji,jj)
108               ENDIF
109            END DO
110         END DO
111      END DO
112
113      DO jk = 1, jpkm1
114         DO jj = 1, jpj
115            DO ji = 1, jpi
116
117!    DENITRIFICATION FACTOR COMPUTED FROM O2 LEVELS
118!    ----------------------------------------------
119
120               nitrfac(ji,jj,jk) = MAX(  0.e0, 0.4 * ( 6.e-6  - trn(ji,jj,jk,jpoxy) )    &
121                  &                                / ( oxymin + trn(ji,jj,jk,jpoxy) )  )
122            END DO
123         END DO
124      END DO
125
126      nitrfac(:,:,:) = MIN( 1., nitrfac(:,:,:) )
127
128
129      DO jk = 1, jpkm1
130         DO jj = 1, jpj
131            DO ji = 1, jpi
132
133!     DOC ammonification. Depends on depth, phytoplankton biomass
134!     and a limitation term which is supposed to be a parameterization
135!     of the bacterial activity.
136!     ----------------------------------------------------------------
137               zremik = xremik * xstep / 1.e-6 * xlimbac(ji,jj,jk)         &
138# if defined key_off_degrad
139                  &            * facvol(ji,jj,jk)              &
140# endif
141                  &            * zdepbac(ji,jj,jk)
142               zremik = MAX( zremik, 5.5e-4 * xstep )
143
144!     Ammonification in oxic waters with oxygen consumption
145!     -----------------------------------------------------
146               zolimi(ji,jj,jk) = MIN(  ( trn(ji,jj,jk,jpoxy) - rtrn ) / o2ut,  &
147                  &                    zremik * ( 1.- nitrfac(ji,jj,jk) ) * trn(ji,jj,jk,jpdoc)  ) 
148
149!     Ammonification in suboxic waters with denitrification
150!     -------------------------------------------------------
151               denitr(ji,jj,jk) = MIN(  ( trn(ji,jj,jk,jpno3) - rtrn ) / rdenit,   &
152                  &                     zremik * nitrfac(ji,jj,jk) * trn(ji,jj,jk,jpdoc)  )
153            END DO
154         END DO
155      END DO
156
157      DO jk = 1, jpkm1
158         DO jj = 1, jpj
159            DO ji = 1, jpi
160               zolimi (ji,jj,jk) = MAX( 0.e0, zolimi (ji,jj,jk) )
161               denitr (ji,jj,jk) = MAX( 0.e0, denitr (ji,jj,jk) )
162            END DO
163         END DO
164      END DO
165
166      DO jk = 1, jpkm1
167         DO jj = 1, jpj
168            DO ji = 1, jpi
169
170!    NH4 nitrification to NO3. Ceased for oxygen concentrations
171!    below 2 umol/L. Inhibited at strong light
172!    ----------------------------------------------------------
173               zonitr  = nitrif * xstep * trn(ji,jj,jk,jpnh4) / ( 1.+ emoy(ji,jj,jk) )     &
174# if defined key_off_degrad
175                  &      * facvol(ji,jj,jk)              &
176# endif
177                  &      * ( 1.- nitrfac(ji,jj,jk) )
178
179!
180!   Update of the tracers trends
181!   ----------------------------
182
183              tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) - zonitr
184              tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) + zonitr
185              tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - o2nit * zonitr
186              tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) - rno3  * zonitr
187
188            END DO
189         END DO
190      END DO
191
192       IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
193         WRITE(charout, FMT="('rem1')")
194         CALL prt_ctl_trc_info(charout)
195         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
196       ENDIF
197
198      DO jk = 1, jpkm1
199         DO jj = 1, jpj
200            DO ji = 1, jpi
201
202!    Bacterial uptake of iron. No iron is available in DOC. So
203!    Bacteries are obliged to take up iron from the water. Some
204!    studies (especially at Papa) have shown this uptake to be
205!    significant
206!    ----------------------------------------------------------
207               zbactfer = 15.e-6 * rfact2 * 4.* 0.4 * prmax(ji,jj,jk)           &
208                  &               * ( xlimphy(ji,jj,jk) * zdepbac(ji,jj,jk))**2           &
209                  &                  / ( xkgraz2 + zdepbac(ji,jj,jk) )                    &
210                  &                  * ( 0.5 + SIGN( 0.5, trn(ji,jj,jk,jpfer) -2.e-11 )  )
211
212               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zbactfer
213#if defined key_kriest
214               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zbactfer
215#else
216               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zbactfer
217#endif
218
219            END DO
220         END DO
221      END DO
222
223       IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
224         WRITE(charout, FMT="('rem2')")
225         CALL prt_ctl_trc_info(charout)
226         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
227       ENDIF
228
229      DO jk = 1, jpkm1
230         DO jj = 1, jpj
231            DO ji = 1, jpi
232
233!    POC disaggregation by turbulence and bacterial activity.
234!    -------------------------------------------------------------
235               zremip = xremip * xstep * tgfunc(ji,jj,jk)   &
236# if defined key_off_degrad
237                  &            * facvol(ji,jj,jk)              &
238# endif
239                  &            * ( 1.- 0.5 * nitrfac(ji,jj,jk) )
240
241!    POC disaggregation rate is reduced in anoxic zone as shown by
242!    sediment traps data. In oxic area, the exponent of the martin s
243!    law is around -0.87. In anoxic zone, it is around -0.35. This
244!    means a disaggregation constant about 0.5 the value in oxic zones
245!    -----------------------------------------------------------------
246               zorem  = zremip * trn(ji,jj,jk,jppoc)
247               zofer  = zremip * trn(ji,jj,jk,jpsfe)
248#if ! defined key_kriest
249               zorem2 = zremip * trn(ji,jj,jk,jpgoc)
250               zofer2 = zremip * trn(ji,jj,jk,jpbfe)
251#else
252               zorem2 = zremip * trn(ji,jj,jk,jpnum)
253#endif
254
255!  Update the appropriate tracers trends
256!  -------------------------------------
257
258               tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zorem
259               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zofer
260#if defined key_kriest
261               tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - zorem
262               tra(ji,jj,jk,jpnum) = tra(ji,jj,jk,jpnum) - zorem2
263               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zofer
264#else
265               tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zorem2 - zorem
266               tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) - zorem2
267               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zofer2 - zofer
268               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) - zofer2
269#endif
270
271            END DO
272         END DO
273      END DO
274
275       IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
276         WRITE(charout, FMT="('rem3')")
277         CALL prt_ctl_trc_info(charout)
278         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
279       ENDIF
280
281      DO jk = 1, jpkm1
282         DO jj = 1, jpj
283            DO ji = 1, jpi
284
285!     Remineralization rate of BSi depedant on T and saturation
286!     ---------------------------------------------------------
287               zsatur  = ( sio3eq(ji,jj,jk) - trn(ji,jj,jk,jpsil) ) / ( sio3eq(ji,jj,jk) + rtrn )
288               zsatur  = MAX( rtrn, zsatur )
289               zsatur2 = zsatur * ( 1. + tn(ji,jj,jk) / 400.)**4
290               znusil  = 0.225  * ( 1. + tn(ji,jj,jk) / 15.) * zsatur + 0.775 * zsatur2**9
291#    if defined key_off_degrad
292               zsiremin = xsirem * xstep * znusil * facvol(ji,jj,jk)
293# else
294               zsiremin = xsirem * xstep * znusil
295#    endif
296               zosil = zsiremin * trn(ji,jj,jk,jpdsi)
297
298               tra(ji,jj,jk,jpdsi) = tra(ji,jj,jk,jpdsi) - zosil
299               tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) + zosil
300
301               !
302            END DO
303         END DO
304      END DO
305
306      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
307         WRITE(charout, FMT="('rem4')")
308         CALL prt_ctl_trc_info(charout)
309         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
310       ENDIF
311
312      zfesatur(:,:,:) = 0.6e-9
313!CDIR NOVERRCHK
314      DO jk = 1, jpkm1
315!CDIR NOVERRCHK
316         DO jj = 1, jpj
317!CDIR NOVERRCHK
318            DO ji = 1, jpi
319!
320!      Compute de different ratios for scavenging of iron
321!      --------------------------------------------------
322
323#if  defined key_kriest
324                zdenom1 = trn(ji,jj,jk,jppoc) / &
325           &           ( trn(ji,jj,jk,jppoc) + trn(ji,jj,jk,jpdsi) + trn(ji,jj,jk,jpcal) + rtrn )
326#else
327                zdenom = 1. / ( trn(ji,jj,jk,jppoc) + trn(ji,jj,jk,jpgoc)  &
328           &            + trn(ji,jj,jk,jpdsi) + trn(ji,jj,jk,jpcal) + rtrn )
329
330                zdenom1 = trn(ji,jj,jk,jppoc) * zdenom
331                zdenom2 = trn(ji,jj,jk,jpgoc) * zdenom
332#endif
333
334
335!     scavenging rate of iron. this scavenging rate depends on the
336!     load in particles on which they are adsorbed. The
337!     parameterization has been taken from studies on Th
338!     ------------------------------------------------------------
339               zkeq = fekeq(ji,jj,jk)
340               zfeequi = ( -( 1. + zfesatur(ji,jj,jk) * zkeq - zkeq * trn(ji,jj,jk,jpfer) )               &
341                  &        + SQRT( ( 1. + zfesatur(ji,jj,jk) * zkeq - zkeq * trn(ji,jj,jk,jpfer) )**2       &
342                  &               + 4. * trn(ji,jj,jk,jpfer) * zkeq) ) / ( 2. * zkeq )
343
344#if defined key_kriest
345               zlam1b = 3.e-5 + xlam1 * (  trn(ji,jj,jk,jppoc)                   &
346                  &                      + trn(ji,jj,jk,jpcal) + trn(ji,jj,jk,jpdsi)  ) * 1.e6
347#else
348               zlam1b = 3.e-5 + xlam1 * (  trn(ji,jj,jk,jppoc) + trn(ji,jj,jk,jpgoc)   &
349                  &                      + trn(ji,jj,jk,jpcal) + trn(ji,jj,jk,jpdsi)  ) * 1.e6
350#endif
351
352# if defined key_off_degrad
353               zscave = zfeequi * zlam1b * xstep  * facvol(ji,jj,jk)
354# else
355               zscave = zfeequi * zlam1b * xstep
356# endif
357
358!  Increased scavenging for very high iron concentrations
359!  found near the coasts due to increased lithogenic particles
360!  and let s say it unknown processes (precipitation, ...)
361!  -----------------------------------------------------------
362               zlamfac = MAX( 0.e0, ( gphit(ji,jj) + 55.) / 30. )
363               zlamfac = MIN( 1.  , zlamfac )
364#if ! defined key_kriest
365               zlam1b = (  80.* ( trn(ji,jj,jk,jpdoc) + 35.e-6 )                           &
366                  &     + 698.*   trn(ji,jj,jk,jppoc) + 1.05e4 * trn(ji,jj,jk,jpgoc)  )                    &
367                  &   * xdiss(ji,jj,jk) + 1E-4 * (1.-zlamfac)                &
368                  &   + xlam1 * MAX( 0.e0, ( trn(ji,jj,jk,jpfer) * 1.e9 - 1.)  )
369#else
370               zlam1b = (  80.* (trn(ji,jj,jk,jpdoc) + 35E-6)           &
371                  &     + 698.*  trn(ji,jj,jk,jppoc)  )                    &
372                  &   * xdiss(ji,jj,jk) + 1E-4 * (1.-zlamfac)           &
373                  &   + xlam1 * MAX( 0.e0, ( trn(ji,jj,jk,jpfer) * 1.e9 - 1.)  )
374#endif
375
376# if defined key_off_degrad
377               zaggdfe = zlam1b * xstep * 0.5 * ( trn(ji,jj,jk,jpfer) - zfeequi ) * facvol(ji,jj,jk)
378# else
379               zaggdfe = zlam1b * xstep * 0.5 * ( trn(ji,jj,jk,jpfer) - zfeequi )
380# endif
381
382               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zscave - zaggdfe
383
384#if defined key_kriest
385               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zscave * zdenom1
386#else
387               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zscave * zdenom1
388               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zscave * zdenom2
389#endif
390
391            END DO
392         END DO
393      END DO
394      !
395
396       IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
397         WRITE(charout, FMT="('rem5')")
398         CALL prt_ctl_trc_info(charout)
399         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
400       ENDIF
401
402!     Update the arrays TRA which contain the biological sources and sinks
403!     --------------------------------------------------------------------
404
405      DO jk = 1, jpkm1
406         tra(:,:,jk,jppo4) = tra(:,:,jk,jppo4) + zolimi(:,:,jk) + denitr(:,:,jk)
407         tra(:,:,jk,jpnh4) = tra(:,:,jk,jpnh4) + zolimi(:,:,jk) + denitr(:,:,jk)
408         tra(:,:,jk,jpno3) = tra(:,:,jk,jpno3) - denitr(:,:,jk) * rdenit
409         tra(:,:,jk,jpdoc) = tra(:,:,jk,jpdoc) - zolimi(:,:,jk) - denitr(:,:,jk)
410         tra(:,:,jk,jpoxy) = tra(:,:,jk,jpoxy) - zolimi(:,:,jk) * o2ut
411         tra(:,:,jk,jpdic) = tra(:,:,jk,jpdic) + zolimi(:,:,jk) + denitr(:,:,jk)
412         tra(:,:,jk,jptal) = tra(:,:,jk,jptal) + denitr(:,:,jk) * rno3 * rdenit
413     END DO
414
415       IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
416         WRITE(charout, FMT="('rem6')")
417         CALL prt_ctl_trc_info(charout)
418         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
419       ENDIF
420
421   END SUBROUTINE p4z_rem
422
423   SUBROUTINE p4z_rem_init
424
425      !!----------------------------------------------------------------------
426      !!                  ***  ROUTINE p4z_rem_init  ***
427      !!
428      !! ** Purpose :   Initialization of remineralization parameters
429      !!
430      !! ** Method  :   Read the nampisrem namelist and check the parameters
431      !!      called at the first timestep (nittrc000)
432      !!
433      !! ** input   :   Namelist nampisrem
434      !!
435      !!----------------------------------------------------------------------
436
437      NAMELIST/nampisrem/ xremik, xremip, nitrif, xsirem, xlam1, oxymin
438
439      REWIND( numnat )                     ! read numnat
440      READ  ( numnat, nampisrem )
441
442      IF(lwp) THEN                         ! control print
443         WRITE(numout,*) ' '
444         WRITE(numout,*) ' Namelist parameters for remineralization, nampisrem'
445         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
446         WRITE(numout,*) '    remineralisation rate of POC              xremip    =', xremip
447         WRITE(numout,*) '    remineralization rate of DOC              xremik    =', xremik
448         WRITE(numout,*) '    remineralization rate of Si               xsirem    =', xsirem
449         WRITE(numout,*) '    scavenging rate of Iron                   xlam1     =', xlam1
450         WRITE(numout,*) '    NH4 nitrification rate                    nitrif    =', nitrif
451         WRITE(numout,*) '    halk saturation constant for anoxia       oxymin    =', oxymin
452      ENDIF
453
454   END SUBROUTINE p4z_rem_init
455
456
457
458
459
460#else
461   !!======================================================================
462   !!  Dummy module :                                   No PISCES bio-model
463   !!======================================================================
464CONTAINS
465   SUBROUTINE p4z_rem                    ! Empty routine
466   END SUBROUTINE p4z_rem
467#endif 
468
469   !!======================================================================
470END MODULE p4zrem
Note: See TracBrowser for help on using the repository browser.