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 trunk/NEMOGCM/NEMO/TOP_SRC/PISCES – NEMO

source: trunk/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zrem.F90 @ 2730

Last change on this file since 2730 was 2715, checked in by rblod, 13 years ago

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