source: trunk/NEMO/TOP_SRC/PISCES/p4zrem.F90 @ 935

Last change on this file since 935 was 935, checked in by cetlod, 13 years ago

adding modules for PISCES SMS model, see ticket 141

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