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.
p4zprod.F90 in NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z – NEMO

source: NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zprod.F90 @ 14416

Last change on this file since 14416 was 14416, checked in by cetlod, 4 years ago

dev_r14383_PISCES_NEWDEV_PISCO : minor improvments
  • Property svn:keywords set to Id
File size: 26.5 KB
Line 
1MODULE p4zprod
2   !!======================================================================
3   !!                         ***  MODULE p4zprod  ***
4   !! TOP :  Growth Rate of the two phytoplankton groups of PISCES
5   !!======================================================================
6   !! History :   1.0  !  2004     (O. Aumont) Original code
7   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
8   !!             3.4  !  2011-05  (O. Aumont, C. Ethe) New parameterization of light limitation
9   !!----------------------------------------------------------------------
10   !!   p4z_prod       : Compute the growth Rate of the two phytoplanktons groups
11   !!   p4z_prod_init  : Initialization of the parameters for growth
12   !!   p4z_prod_alloc : Allocate variables for growth
13   !!----------------------------------------------------------------------
14   USE oce_trc         ! shared variables between ocean and passive tracers
15   USE trc             ! passive tracers common variables
16   USE sms_pisces      ! PISCES Source Minus Sink variables
17   USE p4zlim          ! Co-limitations of differents nutrients
18   USE prtctl          ! print control for debugging
19   USE iom             ! I/O manager
20
21   IMPLICIT NONE
22   PRIVATE
23
24   PUBLIC   p4z_prod         ! called in p4zbio.F90
25   PUBLIC   p4z_prod_init    ! called in trcsms_pisces.F90
26   PUBLIC   p4z_prod_alloc   ! called in trcini_pisces.F90
27
28   REAL(wp), PUBLIC ::   pislopen     !:  P-I slope of nanophytoplankton
29   REAL(wp), PUBLIC ::   pisloped     !:  P-I slope of diatoms
30   REAL(wp), PUBLIC ::   xadap        !:  Adaptation factor to low light
31   REAL(wp), PUBLIC ::   excretn      !:  Excretion ratio of nanophyto
32   REAL(wp), PUBLIC ::   excretd      !:  Excretion ratio of diatoms
33   REAL(wp), PUBLIC ::   bresp        !:  Basal respiration rate
34   REAL(wp), PUBLIC ::   chlcnm       !:  Maximum Chl/C ratio of nano
35   REAL(wp), PUBLIC ::   chlcdm       !:  Maximum Chl/C ratio of diatoms
36   REAL(wp), PUBLIC ::   chlcmin      !:  Minimum Chl/C ratio of phytoplankton
37   REAL(wp), PUBLIC ::   fecnm        !:  Maximum Fe/C ratio of nano
38   REAL(wp), PUBLIC ::   fecdm        !:  Maximum Fe/C ratio of diatoms
39   REAL(wp), PUBLIC ::   grosip       !:  Mean Si/C ratio of diatoms
40
41   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   quotan   !: proxy of N quota in Nanophyto
42   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   quotad   !: proxy of N quota in diatoms
43   
44   REAL(wp) ::   r1_rday    ! 1 / rday
45   REAL(wp) ::   texcretn   ! 1 - excretn
46   REAL(wp) ::   texcretd   ! 1 - excretd       
47
48   !! * Substitutions
49#  include "do_loop_substitute.h90"
50#  include "domzgr_substitute.h90"
51   !!----------------------------------------------------------------------
52   !! NEMO/TOP 4.0 , NEMO Consortium (2018)
53   !! $Id$
54   !! Software governed by the CeCILL license (see ./LICENSE)
55   !!----------------------------------------------------------------------
56CONTAINS
57
58   SUBROUTINE p4z_prod( kt , knt, Kbb, Kmm, Krhs )
59      !!---------------------------------------------------------------------
60      !!                     ***  ROUTINE p4z_prod  ***
61      !!
62      !! ** Purpose :   Computes phytoplankton production depending on
63      !!                light, temperature and nutrient availability
64      !!                Computes also the uptake of Iron and Si as well
65      !!                as the chlorophyll content of the cells
66      !!                PISCES relies on a mixed Monod-Quota formalism
67      !!---------------------------------------------------------------------
68      INTEGER, INTENT(in) ::   kt, knt   !
69      INTEGER, INTENT(in) ::   Kbb, Kmm, Krhs  ! time level indices
70      !
71      INTEGER  ::   ji, jj, jk
72      REAL(wp) ::   zsilfac, znanotot, zdiattot, zconctemp, zconctemp2
73      REAL(wp) ::   zratio, zmax, zsilim, ztn, zadap, zlim, zsiborn
74      REAL(wp) ::   zpptot, zpnewtot, zpregtot, zprochln, zprochld
75      REAL(wp) ::   zproddoc, zprodsil, zprodfer, zprodlig
76      REAL(wp) ::   zpislopen, zpisloped, zfact
77      REAL(wp) ::   zratiosi, zmaxsi, zlimfac, zsizetmp, zfecnm, zfecdm
78      REAL(wp) ::   zprod, zval
79      CHARACTER (len=25) :: charout
80      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprmaxn,zprmaxd
81      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpislopeadn, zpislopeadd, zysopt 
82      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprdia, zprbio, zprchld, zprchln   
83      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprorcan, zprorcad, zprofed, zprofen
84      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpronewn, zpronewd
85      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmxl_fac, zmxl_chl
86      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zpligprod
87      !!---------------------------------------------------------------------
88      !
89      IF( ln_timing )   CALL timing_start('p4z_prod')
90      !
91      !  Allocate temporary workspace
92      !
93      zprorcan(:,:,:) = 0._wp ; zprorcad(:,:,:) = 0._wp ; zprofed(:,:,:) = 0._wp
94      zprofen (:,:,:) = 0._wp ; zysopt  (:,:,:) = 0._wp
95      zpronewn(:,:,:) = 0._wp ; zpronewd(:,:,:) = 0._wp ; zprdia(:,:,:) = 0._wp
96      zprbio  (:,:,:) = 0._wp ; zprchld (:,:,:) = 0._wp ; zprchln(:,:,:) = 0._wp 
97      zmxl_fac(:,:,:) = 0._wp ; zmxl_chl(:,:,:) = 0._wp 
98
99      ! Computation of the maximimum production. Based on a Q10 description
100      ! of the thermal dependency
101      ! Parameters are taken from Bissinger et al. (2008)
102      zprmaxn(:,:,:) = 0.8_wp * r1_rday * tgfunc(:,:,:)
103      zprmaxd(:,:,:) = zprmaxn(:,:,:)
104
105      ! Intermittency is supposed to have a similar effect on production as
106      ! day length (Shatwell et al., 2012). The correcting factor is zmxl_fac.
107      ! zmxl_chl is the fractional day length and is used to compute the mean
108      ! PAR during daytime. The effect of mixing is computed using the
109      ! absolute light level definition of the euphotic zone
110      ! -------------------------------------------------------------------------
111      DO_3D( 1, 1, 1, 1, 1, jpkm1 )
112         IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
113            zval = MAX( 1., strn(ji,jj) )
114            IF( gdepw(ji,jj,jk+1,Kmm) <= hmld(ji,jj) ) THEN
115               zval = zval * MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn ))
116            ENDIF
117            zmxl_chl(ji,jj,jk) = zval / 24.
118            zmxl_fac(ji,jj,jk) = 1.0 - exp( -0.26 * zval )
119         ENDIF
120      END_3D
121
122      zprbio(:,:,:) = zprmaxn(:,:,:) * zmxl_fac(:,:,:)
123      zprdia(:,:,:) = zprmaxd(:,:,:) * zmxl_fac(:,:,:)
124
125      ! The formulation proposed by Geider et al. (1997) has been modified
126      ! to exclude the effect of nutrient limitation and temperature in the PI
127      ! curve following Vichi et al. (2007)
128      ! -----------------------------------------------------------------------
129      DO_3D( 1, 1, 1, 1, 1, jpkm1 )
130         IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
131            ztn         = MAX( 0., ts(ji,jj,jk,jp_tem,Kmm) - 15. )
132            zadap       = xadap * ztn / ( 2.+ ztn )
133            zconctemp   = MAX( 0.e0 , tr(ji,jj,jk,jpdia,Kbb) - xsizedia )
134            zconctemp2  = tr(ji,jj,jk,jpdia,Kbb) - zconctemp
135            !
136            ! The initial slope of the PI curve can be increased for nano
137            ! to account for photadaptation, for instance in the DCM
138            ! This parameterization is adhoc and should be either
139            ! improved or removed in future versions of the model
140
141            ! Nanophytoplankton
142            zpislopeadn(ji,jj,jk) = pislopen * ( 1.+ zadap  * EXP( -0.25 * enano(ji,jj,jk) ) )  &
143            &                   * tr(ji,jj,jk,jpnch,Kbb) /( tr(ji,jj,jk,jpphy,Kbb) * 12. + rtrn)
144
145            ! Diatoms
146            zpislopeadd(ji,jj,jk) = (pislopen * zconctemp2 + pisloped * zconctemp) / ( tr(ji,jj,jk,jpdia,Kbb) + rtrn )   &
147            &                   * tr(ji,jj,jk,jpdch,Kbb) /( tr(ji,jj,jk,jpdia,Kbb) * 12. + rtrn)
148         ENDIF
149      END_3D
150
151      DO_3D( 1, 1, 1, 1, 1, jpkm1 )
152         IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
153             ! Computation of production function for Carbon
154             ! Actual light levels are used here
155             ! ----------------------------------------------
156             zpislopen = zpislopeadn(ji,jj,jk) / ( ( r1_rday + bresp * r1_rday ) &
157             &            * zmxl_fac(ji,jj,jk) * rday + rtrn)
158             zpisloped = zpislopeadd(ji,jj,jk) / ( ( r1_rday + bresp * r1_rday ) &
159             &            * zmxl_fac(ji,jj,jk) * rday + rtrn)
160             zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) )  )
161             zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) )  )
162
163             !  Computation of production function for Chlorophyll
164             !  Mean light level in the mixed layer (when appropriate)
165             !  is used here (acclimation is in general slower than
166             !  the characteristic time scales of vertical mixing)
167             !  ------------------------------------------------------
168             zpislopen = zpislopeadn(ji,jj,jk) / ( zprmaxn(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn )
169             zpisloped = zpislopeadd(ji,jj,jk) / ( zprmaxd(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn )
170             zprchln(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) ) )
171             zprchld(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) ) )
172         ENDIF
173      END_3D
174
175      !  Computation of a proxy of the N/C quota from nutrient limitation
176      !  and light limitation. Steady state is assumed to allow the computation
177      !  ----------------------------------------------------------------------
178      DO_3D( 1, 1, 1, 1, 1, jpkm1 )
179          zval = MIN( xnanopo4(ji,jj,jk), ( xnanonh4(ji,jj,jk) + xnanono3(ji,jj,jk) ) )   &
180          &      * zprmaxn(ji,jj,jk) / ( zprbio(ji,jj,jk) + rtrn )
181          quotan(ji,jj,jk) = MIN( 1., 0.2 + 0.8 * zval )
182          zval = MIN( xdiatpo4(ji,jj,jk), ( xdiatnh4(ji,jj,jk) + xdiatno3(ji,jj,jk) ) )   &
183          &      * zprmaxd(ji,jj,jk) / ( zprdia(ji,jj,jk) + rtrn )
184          quotad(ji,jj,jk) = MIN( 1., 0.2 + 0.8 * zval )
185      END_3D
186
187
188      DO_3D( 1, 1, 1, 1, 1, jpkm1 )
189
190          IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
191             ! Si/C of diatoms
192             ! ------------------------
193             ! Si/C increases with iron stress and silicate availability
194             ! Si/C is arbitrariliy increased for very high Si concentrations
195             ! to mimic the very high ratios observed in the Southern Ocean (zsilfac)
196             ! A parameterization derived from Flynn (2003) is used for the control
197             ! when Si is not limiting which is similar to the parameterisation
198             ! proposed by Gurney and Davidson (1999).
199             ! -----------------------------------------------------------------------
200            zlim  = tr(ji,jj,jk,jpsil,Kbb) / ( tr(ji,jj,jk,jpsil,Kbb) + xksi1 )
201            zsilim = xlimdia(ji,jj,jk) * zprdia(ji,jj,jk) / ( zprmaxd(ji,jj,jk) + rtrn )
202            zsiborn = tr(ji,jj,jk,jpsil,Kbb) * tr(ji,jj,jk,jpsil,Kbb) * tr(ji,jj,jk,jpsil,Kbb)
203            IF (gphit(ji,jj) < -30 ) THEN
204              zsilfac = 1. + 2. * zsiborn / ( zsiborn + xksi2**3 )
205            ELSE
206              zsilfac = 1. +      zsiborn / ( zsiborn + xksi2**3 )
207            ENDIF
208            zratiosi = 1.0 - tr(ji,jj,jk,jpdsi,Kbb) / ( tr(ji,jj,jk,jpdia,Kbb) + rtrn ) / ( zsilfac * grosip * 3.0 + rtrn )
209            zratiosi = MAX(0., MIN(1.0, zratiosi) )
210            zmaxsi  = (1.0 + 0.1**4) * zratiosi**4 / ( zratiosi**4 + 0.1**4 )
211            IF( xlimsi(ji,jj,jk) /= xlimdia(ji,jj,jk) ) THEN
212               zysopt(ji,jj,jk) = zlim * zsilfac * grosip * 1.0 * zmaxsi
213            ELSE
214               zysopt(ji,jj,jk) = zlim * zsilfac * grosip * 1.0 * zsilim**0.7 * zmaxsi
215            ENDIF
216        ENDIF
217      END_3D
218
219      ! Sea-ice effect on production
220      ! No production is assumed below sea ice
221      ! --------------------------------------
222      DO_3D( 1, 1, 1, 1, 1, jpkm1 )
223         zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1. - fr_i(ji,jj) )
224         zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1. - fr_i(ji,jj) )
225      END_3D
226
227      ! Computation of the various production  and nutrient uptake terms
228      ! ---------------------------------------------------------------
229      DO_3D( 1, 1, 1, 1, 1, jpkm1 )
230         IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
231            !  production terms for nanophyto. (C)
232            zprorcan(ji,jj,jk) = zprbio(ji,jj,jk)  * xlimphy(ji,jj,jk) * tr(ji,jj,jk,jpphy,Kbb) * rfact2
233
234            !  New production (uptake of NO3)
235            zpronewn(ji,jj,jk)  = zprorcan(ji,jj,jk)* xnanono3(ji,jj,jk) / ( xnanono3(ji,jj,jk) + xnanonh4(ji,jj,jk) + rtrn )
236            !
237            ! Size computation
238            ! Size is made a function of the limitation of of phytoplankton growth
239            ! Strongly limited cells are supposed to be smaller. sizena is the
240            ! size at time step t+1 and is thus updated at the end of the
241            ! current time step
242            ! --------------------------------------------------------------------
243            zlimfac = xlimphy(ji,jj,jk) * zprchln(ji,jj,jk) / ( zprmaxn(ji,jj,jk) + rtrn )
244            zsizetmp = 1.0 + 1.3 * ( xsizern - 1.0 ) * zlimfac**3/(0.3 + zlimfac**3)
245            sizena(ji,jj,jk) = min(xsizern, max( sizena(ji,jj,jk), zsizetmp ) )
246
247            ! Iron uptake rates of nanophytoplankton. Upregulation is 
248            ! not parameterized at low iron concentrations as observations
249            ! do not suggest it for accimated cells. Uptake is
250            ! downregulated when the quota is close to the maximum quota
251            zfecnm = xqfuncfecn(ji,jj,jk) + ( fecnm - xqfuncfecn(ji,jj,jk) ) * ( xnanono3(ji,jj,jk) + xnanonh4(ji,jj,jk) )
252            zratio = 1.0 - MIN(1.0,tr(ji,jj,jk,jpnfe,Kbb) / ( tr(ji,jj,jk,jpphy,Kbb) * zfecnm + rtrn ) )
253            zmax   = MAX( 0., MIN( 1.0, zratio**2/ (0.05**2+zratio**2) ) ) 
254            zprofen(ji,jj,jk) = zfecnm * zprmaxn(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  &
255            &          * (1. + 0.8 * xnanono3(ji,jj,jk) / ( rtrn + xnanono3(ji,jj,jk)  &
256            &          + xnanonh4(ji,jj,jk) ) * (1. - xnanofer(ji,jj,jk) ) )   &
257            &          * xnanofer(ji,jj,jk) * zmax * tr(ji,jj,jk,jpphy,Kbb) * rfact2
258
259            ! production terms of diatoms (C)
260            zprorcad(ji,jj,jk) = zprdia(ji,jj,jk) * xlimdia(ji,jj,jk) * tr(ji,jj,jk,jpdia,Kbb) * rfact2
261
262            ! New production (uptake of NO3)
263            zpronewd(ji,jj,jk) = zprorcad(ji,jj,jk) * xdiatno3(ji,jj,jk) / ( xdiatno3(ji,jj,jk) + xdiatnh4(ji,jj,jk) + rtrn )
264
265            ! Size computation
266            ! Size is made a function of the limitation of of phytoplankton growth
267            ! Strongly limited cells are supposed to be smaller. sizeda is
268            ! size at time step t+1 and is thus updated at the end of the
269            ! current time step.
270            ! --------------------------------------------------------------------
271            zlimfac = zprchld(ji,jj,jk) * xlimdia(ji,jj,jk) / ( zprmaxd(ji,jj,jk) + rtrn )
272            zsizetmp = 1.0 + 1.3 * ( xsizerd - 1.0 ) * zlimfac**3/(0.3 + zlimfac**3)
273            sizeda(ji,jj,jk) = min(xsizerd, max( sizeda(ji,jj,jk), zsizetmp ) )
274
275            ! Iron uptake rates of diatoms. Upregulation is 
276            ! not parameterized at low iron concentrations as observations
277            ! do not suggest it for accimated cells. Uptake is
278            ! downregulated when the quota is close to the maximum quota
279            zfecdm = xqfuncfecd(ji,jj,jk) + ( fecdm - xqfuncfecd(ji,jj,jk) ) * ( xdiatno3(ji,jj,jk) + xdiatnh4(ji,jj,jk) )
280            zratio = 1.0 - MIN(1.0, tr(ji,jj,jk,jpdfe,Kbb) / ( tr(ji,jj,jk,jpdia,Kbb) * zfecdm + rtrn ) )
281            zmax   = MAX( 0., MIN( 1.0, zratio**2/ (0.05**2+zratio**2) ) ) 
282            zprofed(ji,jj,jk) = zfecdm * zprmaxd(ji,jj,jk) * (1.0 - fr_i(ji,jj) )  &
283            &          * (1. + 0.8 * xdiatno3(ji,jj,jk) / ( rtrn + xdiatno3(ji,jj,jk)  &
284            &          + xdiatnh4(ji,jj,jk) ) * (1. - xdiatfer(ji,jj,jk) ) )   &
285            &          * xdiatfer(ji,jj,jk) * zmax * tr(ji,jj,jk,jpdia,Kbb) * rfact2
286         ENDIF
287      END_3D
288
289      ! Computation of the chlorophyll production terms
290      ! The parameterization is taken from Geider et al. (1997)
291      ! -------------------------------------------------------
292      DO_3D( 1, 1, 1, 1, 1, jpkm1 )
293         IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
294            !  production terms for nanophyto. ( chlorophyll )
295            znanotot = enanom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
296            zprod    = rday * zprorcan(ji,jj,jk) * zprchln(ji,jj,jk) * xlimphy(ji,jj,jk)
297            zprochln = chlcmin * 12. * zprorcan (ji,jj,jk)
298            zprochln = zprochln + (chlcnm - chlcmin) * 12. * zprod / &
299                                  & (  zpislopeadn(ji,jj,jk) * znanotot +rtrn)
300
301            !  production terms for diatoms ( chlorophyll )
302            zdiattot = ediatm(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
303            zprod    = rday * zprorcad(ji,jj,jk) * zprchld(ji,jj,jk) * xlimdia(ji,jj,jk)
304            zprochld = chlcmin * 12. * zprorcad(ji,jj,jk)
305            zprochld = zprochld + (chlcdm - chlcmin) * 12. * zprod / &
306                                  & ( zpislopeadd(ji,jj,jk) * zdiattot +rtrn )
307
308            !   Update the arrays TRA which contain the Chla sources and sinks
309            tr(ji,jj,jk,jpnch,Krhs) = tr(ji,jj,jk,jpnch,Krhs) + zprochln * texcretn
310            tr(ji,jj,jk,jpdch,Krhs) = tr(ji,jj,jk,jpdch,Krhs) + zprochld * texcretd
311         ENDIF
312      END_3D
313
314      !   Update the arrays TRA which contain the biological sources and sinks
315      DO_3D( 1, 1, 1, 1, 1, jpkm1 )
316        IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
317           zpptot   = zprorcan(ji,jj,jk) + zprorcad(ji,jj,jk)
318           zpnewtot = zpronewn(ji,jj,jk) + zpronewd(ji,jj,jk)
319           zpregtot = zpptot - zpnewtot
320           zprodsil  = zprmaxd(ji,jj,jk) * zysopt(ji,jj,jk) * rfact2 * tr(ji,jj,jk,jpdia,Kbb)
321           zproddoc  = excretd * zprorcad(ji,jj,jk) + excretn * zprorcan(ji,jj,jk)
322           zprodfer  = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk)
323           !
324           tr(ji,jj,jk,jppo4,Krhs) = tr(ji,jj,jk,jppo4,Krhs) - zpptot
325           tr(ji,jj,jk,jpno3,Krhs) = tr(ji,jj,jk,jpno3,Krhs) - zpnewtot
326           tr(ji,jj,jk,jpnh4,Krhs) = tr(ji,jj,jk,jpnh4,Krhs) - zpregtot
327           tr(ji,jj,jk,jpphy,Krhs) = tr(ji,jj,jk,jpphy,Krhs) + zprorcan(ji,jj,jk) * texcretn
328           tr(ji,jj,jk,jpnfe,Krhs) = tr(ji,jj,jk,jpnfe,Krhs) + zprofen(ji,jj,jk)  * texcretn
329           tr(ji,jj,jk,jpdia,Krhs) = tr(ji,jj,jk,jpdia,Krhs) + zprorcad(ji,jj,jk) * texcretd
330           tr(ji,jj,jk,jpdfe,Krhs) = tr(ji,jj,jk,jpdfe,Krhs) + zprofed(ji,jj,jk)  * texcretd
331           tr(ji,jj,jk,jpdsi,Krhs) = tr(ji,jj,jk,jpdsi,Krhs) + zprodsil
332           tr(ji,jj,jk,jpsil,Krhs) = tr(ji,jj,jk,jpsil,Krhs) - zprodsil
333           tr(ji,jj,jk,jpdoc,Krhs) = tr(ji,jj,jk,jpdoc,Krhs) + zproddoc
334           tr(ji,jj,jk,jpoxy,Krhs) = tr(ji,jj,jk,jpoxy,Krhs) + o2ut * zpregtot + ( o2ut + o2nit ) * zpnewtot
335           !
336           tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) - zprodfer
337           consfe3(ji,jj,jk)   = zprodfer * 75.0 / ( rtrn + ( plig(ji,jj,jk) + 75.0 * (1.0 - plig(ji,jj,jk) ) )   &
338           &                   * tr(ji,jj,jk,jpfer,Kbb) ) / rfact2
339           
340           !
341           tr(ji,jj,jk,jpdic,Krhs) = tr(ji,jj,jk,jpdic,Krhs) - zpptot
342           tr(ji,jj,jk,jptal,Krhs) = tr(ji,jj,jk,jptal,Krhs) + rno3 * ( zpnewtot - zpregtot )
343        ENDIF
344      END_3D
345
346     ! Production and uptake of ligands by phytoplankton. This part is activated
347     ! when ln_ligand is set to .true. in the namelist. Ligand uptake is small
348     ! and based on the FeL model by Morel et al. (2008) and on the study of
349     ! Shaked et al. (2020)
350     ! -------------------------------------------------------------------------
351     IF( ln_ligand ) THEN
352         DO_3D( 1, 1, 1, 1, 1, jpkm1 )
353           IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
354              zproddoc = excretd * zprorcad(ji,jj,jk) + excretn * zprorcan(ji,jj,jk)
355              zprodfer = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk)
356              zprodlig = plig(ji,jj,jk) / ( rtrn + plig(ji,jj,jk) + 75.0 * (1.0 - plig(ji,jj,jk) ) ) * lthet
357              !
358              tr(ji,jj,jk,jplgw,Krhs) = tr(ji,jj,jk,jplgw,Krhs) + zproddoc * ldocp - zprodfer * zprodlig
359           ENDIF
360         END_3D
361     ENDIF
362
363
364    ! Output of the diagnostics
365    ! Total primary production per year
366    IF( iom_use( "tintpp" ) .OR. ( ln_check_mass .AND. kt == nitend .AND. knt == nrdttrc )  )  &
367         & tpp = glob_sum( 'p4zprod', ( zprorcan(:,:,:) + zprorcad(:,:,:) ) * cvol(:,:,:) )
368
369    IF( lk_iomput .AND.  knt == nrdttrc ) THEN
370       zfact = 1.e+3 * rfact2r  !  conversion from mol/l/kt to  mol/m3/s
371       !
372       CALL iom_put( "PPPHYN"  , zprorcan(:,:,:) * zfact * tmask(:,:,:) )  ! primary production by nanophyto
373       CALL iom_put( "PPPHYD"  , zprorcad(:,:,:) * zfact * tmask(:,:,:)   ) ! primary production by diatomes
374       CALL iom_put( "PPNEWN"  , zpronewn(:,:,:) * zfact * tmask(:,:,:)    ) ! new primary production by nanophyto
375       CALL iom_put( "PPNEWD"  , zpronewd(:,:,:) * zfact * tmask(:,:,:)   ) ! new primary production by diatomes
376       CALL iom_put( "PBSi"    , zprorcad(:,:,:) * zfact * tmask(:,:,:) * zysopt(:,:,:)  ) ! biogenic silica production
377       CALL iom_put( "PFeN"    , zprofen(:,:,:)  * zfact * tmask(:,:,:)  ) ! biogenic iron production by nanophyto
378       CALL iom_put( "PFeD"    , zprofed(:,:,:)  * zfact * tmask(:,:,:)  ) ! biogenic iron production by  diatomes
379       IF( ln_ligand .AND. ( iom_use( "LPRODP" ) .OR. iom_use( "LDETP" ) ) ) THEN
380           ALLOCATE(  zpligprod(jpi,jpj,jpk) )
381           zpligprod(:,:,:) = excretd * zprorcad(:,:,:) + excretn * zprorcan(:,:,:)
382           CALL iom_put( "LPRODP"  , zpligprod(:,:,:) * ldocp * 1e9 * zfact * tmask(:,:,:) )
383           !
384           zpligprod(:,:,:) = ( texcretn * zprofen(:,:,:) + texcretd * zprofed(:,:,:) ) & 
385             &                  * plig(:,:,:) / ( rtrn + plig(:,:,:) + 75.0 * (1.0 - plig(:,:,:) ) )
386           CALL iom_put( "LDETP"   , zpligprod(:,:,:)  * lthet * 1e9 * zfact * tmask(:,:,:) )
387           DEALLOCATE(  zpligprod )
388       ENDIF
389       CALL iom_put( "Mumax"   , zprmaxn(:,:,:) * tmask(:,:,:)  ) ! Maximum growth rate
390       CALL iom_put( "MuN"     , zprbio(:,:,:) * xlimphy(:,:,:) * tmask(:,:,:) ) ! Realized growth rate for nanophyto
391       CALL iom_put( "MuD"     , zprdia(:,:,:) * xlimdia(:,:,:) * tmask(:,:,:) ) ! Realized growth rate for diatoms
392       CALL iom_put( "LNlight" , zprbio (:,:,:) / (zprmaxn(:,:,:) + rtrn) * tmask(:,:,:)  )  ! light limitation term
393       CALL iom_put( "LDlight" , zprdia (:,:,:) / (zprmaxd(:,:,:) + rtrn) * tmask(:,:,:)   )
394       CALL iom_put( "TPP"     , ( zprorcan(:,:,:) + zprorcad(:,:,:) ) * zfact * tmask(:,:,:)  )  ! total primary production
395       CALL iom_put( "TPNEW"   , ( zpronewn(:,:,:) + zpronewd(:,:,:) ) * zfact * tmask(:,:,:)  ) ! total new production
396       CALL iom_put( "TPBFE"   , ( zprofen(:,:,:) + zprofed(:,:,:) ) * zfact * tmask(:,:,:)  )  ! total biogenic iron production
397       CALL iom_put( "tintpp"  , tpp * zfact )  !  global total integrated primary production molC/s
398     ENDIF
399
400     IF(sn_cfctl%l_prttrc)   THEN  ! print mean trends (used for debugging)
401         WRITE(charout, FMT="('prod')")
402         CALL prt_ctl_info( charout, cdcomp = 'top' )
403         CALL prt_ctl(tab4d_1=tr(:,:,:,:,Krhs), mask1=tmask, clinfo=ctrcnm)
404     ENDIF
405      !
406      IF( ln_timing )  CALL timing_stop('p4z_prod')
407      !
408   END SUBROUTINE p4z_prod
409
410
411   SUBROUTINE p4z_prod_init
412      !!----------------------------------------------------------------------
413      !!                  ***  ROUTINE p4z_prod_init  ***
414      !!
415      !! ** Purpose :   Initialization of phytoplankton production parameters
416      !!
417      !! ** Method  :   Read the namp4zprod namelist and check the parameters
418      !!      called at the first timestep (nittrc000)
419      !!
420      !! ** input   :   Namelist namp4zprod
421      !!----------------------------------------------------------------------
422      INTEGER ::   ios   ! Local integer
423      !
424      ! Namelist block
425      NAMELIST/namp4zprod/ pislopen, pisloped, xadap, bresp, excretn, excretd,  &
426         &                 chlcnm, chlcdm, chlcmin, fecnm, fecdm, grosip
427      !!----------------------------------------------------------------------
428      !
429      IF(lwp) THEN                         ! control print
430         WRITE(numout,*)
431         WRITE(numout,*) 'p4z_prod_init : phytoplankton growth'
432         WRITE(numout,*) '~~~~~~~~~~~~~'
433      ENDIF
434      !
435      READ  ( numnatp_ref, namp4zprod, IOSTAT = ios, ERR = 901)
436901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namp4zprod in reference namelist' )
437
438      READ  ( numnatp_cfg, namp4zprod, IOSTAT = ios, ERR = 902 )
439902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namp4zprod in configuration namelist' )
440      IF(lwm) WRITE( numonp, namp4zprod )
441
442      IF(lwp) THEN                         ! control print
443         WRITE(numout,*) '   Namelist : namp4zprod'
444         WRITE(numout,*) '      mean Si/C ratio                           grosip       =', grosip
445         WRITE(numout,*) '      P-I slope                                 pislopen     =', pislopen
446         WRITE(numout,*) '      Acclimation factor to low light           xadap        =', xadap
447         WRITE(numout,*) '      excretion ratio of nanophytoplankton      excretn      =', excretn
448         WRITE(numout,*) '      excretion ratio of diatoms                excretd      =', excretd
449         WRITE(numout,*) '      basal respiration in phytoplankton        bresp        =', bresp
450         WRITE(numout,*) '      Maximum Chl/C in phytoplankton            chlcmin      =', chlcmin
451         WRITE(numout,*) '      P-I slope  for diatoms                    pisloped     =', pisloped
452         WRITE(numout,*) '      Minimum Chl/C in nanophytoplankton        chlcnm       =', chlcnm
453         WRITE(numout,*) '      Minimum Chl/C in diatoms                  chlcdm       =', chlcdm
454         WRITE(numout,*) '      Maximum Fe/C in nanophytoplankton         fecnm        =', fecnm
455         WRITE(numout,*) '      Minimum Fe/C in diatoms                   fecdm        =', fecdm
456      ENDIF
457      !
458      r1_rday   = 1._wp / rday 
459      texcretn  = 1._wp - excretn
460      texcretd  = 1._wp - excretd
461      tpp       = 0._wp
462      !
463   END SUBROUTINE p4z_prod_init
464
465
466   INTEGER FUNCTION p4z_prod_alloc()
467      !!----------------------------------------------------------------------
468      !!                     ***  ROUTINE p4z_prod_alloc  ***
469      !!----------------------------------------------------------------------
470      ALLOCATE( quotan(jpi,jpj,jpk), quotad(jpi,jpj,jpk), STAT = p4z_prod_alloc )
471      !
472      IF( p4z_prod_alloc /= 0 ) CALL ctl_stop( 'STOP', 'p4z_prod_alloc : failed to allocate arrays.' )
473      !
474   END FUNCTION p4z_prod_alloc
475
476   !!======================================================================
477END MODULE p4zprod
Note: See TracBrowser for help on using the repository browser.