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.
p5zprod.F90 in NEMO/branches/2019/dev_r11708_aumont_PISCES_QUOTA/src/TOP/PISCES/P4Z – NEMO

source: NEMO/branches/2019/dev_r11708_aumont_PISCES_QUOTA/src/TOP/PISCES/P4Z/p5zprod.F90 @ 14276

Last change on this file since 14276 was 14276, checked in by aumont, 3 years ago

numerous updates to PISCES, PISCES-QUOTA and the sediment module

  • Property svn:keywords set to Id
File size: 41.3 KB
Line 
1MODULE p5zprod
2   !!======================================================================
3   !!                         ***  MODULE p5zprod  ***
4   !! TOP :  Growth Rate of the three phytoplanktons groups
5   !!        PISCES-QUOTA version of the module
6   !!======================================================================
7   !! History :   1.0  !  2004     (O. Aumont) Original code
8   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
9   !!             3.4  !  2011-05  (O. Aumont, C. Ethe) New parameterization of light limitation
10   !!             3.6  !  2015-05  (O. Aumont) PISCES quota
11   !!----------------------------------------------------------------------
12   !!   p5z_prod       :   Compute the growth Rate of the two phytoplanktons groups
13   !!   p5z_prod_init  :   Initialization of the parameters for growth
14   !!   p5z_prod_alloc :   Allocate variables for growth
15   !!----------------------------------------------------------------------
16   USE oce_trc         !  shared variables between ocean and passive tracers
17   USE trc             !  passive tracers common variables
18   USE sms_pisces      !  PISCES Source Minus Sink variables
19   USE p4zlim
20   USE p5zlim          !  Co-limitations of differents nutrients
21   USE prtctl_trc      !  print control for debugging
22   USE iom             !  I/O manager
23
24   IMPLICIT NONE
25   PRIVATE
26
27   PUBLIC   p5z_prod         ! called in p5zbio.F90
28   PUBLIC   p5z_prod_init    ! called in trcsms_pisces.F90
29   PUBLIC   p5z_prod_alloc
30
31   !! * Shared module variables
32   REAL(wp), PUBLIC ::  pislopen        !: P-I slope of nanophytoplankton
33   REAL(wp), PUBLIC ::  pislopep        !: P-I slope of picophytoplankton
34   REAL(wp), PUBLIC ::  pisloped        !: P-I slope of diatoms
35   REAL(wp), PUBLIC ::  xadap           !: Adaptation factor to low light
36   REAL(wp), PUBLIC ::  excretn         !: Excretion ratio of nanophyto
37   REAL(wp), PUBLIC ::  excretp         !: Excretion ratio of picophyto
38   REAL(wp), PUBLIC ::  excretd         !: Excretion ratio of diatoms
39   REAL(wp), PUBLIC ::  bresp           !: Basal respiration rate
40   REAL(wp), PUBLIC ::  thetanpm        !: Maximum Chl/N ratio of picophyto
41   REAL(wp), PUBLIC ::  thetannm        !: Maximum Chl/N ratio of nanophyto
42   REAL(wp), PUBLIC ::  thetandm        !: Maximum Chl/N ratio of diatoms
43   REAL(wp), PUBLIC ::  chlcmin         !: Minimum Chl/C ratio of phytoplankton
44   REAL(wp), PUBLIC ::  grosip          !: Mean Si/C ratio of diatoms
45
46   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   zdaylen ! day length
47   
48   REAL(wp) :: r1_rday                !: 1 / rday
49   REAL(wp) :: texcretn               !: 1 - excretn
50   REAL(wp) :: texcretp               !: 1 - excretp
51   REAL(wp) :: texcretd               !: 1 - excretd       
52
53   !!----------------------------------------------------------------------
54   !! NEMO/TOP 4.0 , NEMO Consortium (2018)
55   !! $Id$
56   !! Software governed by the CeCILL license (see ./LICENSE)
57   !!----------------------------------------------------------------------
58CONTAINS
59
60   SUBROUTINE p5z_prod( kt , knt )
61      !!---------------------------------------------------------------------
62      !!                     ***  ROUTINE p5z_prod  ***
63      !!
64      !! ** Purpose :   Compute the phytoplankton production depending on
65      !!              light, temperature and nutrient availability
66      !!              Computes also the uptake of nutrients. PISCES-quota
67      !!              relies on a full quota formalism
68      !!---------------------------------------------------------------------
69      !
70      INTEGER, INTENT(in) :: kt, knt
71      !
72      INTEGER  ::   ji, jj, jk
73      REAL(wp) ::   zsilfac, znanotot, zpicotot, zdiattot, zconctemp, zconctemp2
74      REAL(wp) ::   zration, zratiop, zratiof, zmax, ztn, zadap
75      REAL(wp) ::   zpronmax, zpropmax, zprofmax, zratio
76      REAL(wp) ::   zlim, zsilfac2, zsiborn, zprod, zprontot, zproptot, zprodtot
77      REAL(wp) ::   zprnutmax, zdocprod, zprochln, zprochld, zprochlp
78      REAL(wp) ::   zpislopen, zpislopep, zpisloped
79      REAL(wp) ::   zrum, zcodel, zargu, zval, zfeup
80      REAL(wp) ::   zqfpmax, zqfnmax, zqfdmax
81      REAL(wp) ::   zfact, zrfact2, zmaxsi, zratiosi, zsizetmp, zlimfac, zsilim
82      CHARACTER (len=25) :: charout
83      REAL(wp), DIMENSION(jpi,jpj    ) :: zmixnano, zmixpico, zmixdiat
84      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpislopeadn, zpislopeadp, zpislopeadd
85      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprnut, zprmaxp, zprmaxn, zprmaxd
86      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprbio, zprpic, zprdia, zysopt
87      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprchln, zprchlp, zprchld
88      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprorcan, zprorcap, zprorcad 
89      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprofed, zprofep, zprofen
90      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpronewn, zpronewp, zpronewd
91      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zproregn, zproregp, zproregd
92      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpropo4n, zpropo4p, zpropo4d
93      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprodopn, zprodopp, zprodopd
94      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zrespn, zrespp, zrespd
95      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zcroissn, zcroissp, zcroissd
96      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmxl_fac, zmxl_chl
97      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpligprod1, zpligprod2
98      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d
99      REAL(wp), ALLOCATABLE, DIMENSION(:,:  ) :: zw2d
100      !!---------------------------------------------------------------------
101      !
102      IF( ln_timing )   CALL timing_start('p5z_prod')
103
104      ! Initialize the local arrays
105      zprorcan(:,:,:) = 0._wp ; zprorcap(:,:,:) = 0._wp ; zprorcad(:,:,:) = 0._wp
106      zprofed (:,:,:) = 0._wp ; zprofep (:,:,:) = 0._wp ; zprofen (:,:,:) = 0._wp
107      zpronewn(:,:,:) = 0._wp ; zpronewp(:,:,:) = 0._wp ; zpronewd(:,:,:) = 0._wp
108      zproregn(:,:,:) = 0._wp ; zproregp(:,:,:) = 0._wp ; zproregd(:,:,:) = 0._wp 
109      zpropo4n(:,:,:) = 0._wp ; zpropo4p(:,:,:) = 0._wp ; zpropo4d(:,:,:) = 0._wp
110      zprdia  (:,:,:) = 0._wp ; zprpic  (:,:,:) = 0._wp ; zprbio  (:,:,:) = 0._wp
111      zprodopn(:,:,:) = 0._wp ; zprodopp(:,:,:) = 0._wp ; zprodopd(:,:,:) = 0._wp
112      zysopt  (:,:,:) = 0._wp
113      zrespn  (:,:,:) = 0._wp ; zrespp  (:,:,:) = 0._wp ; zrespd  (:,:,:) = 0._wp 
114      zmxl_fac(:,:,:) = 0._wp ; zmxl_chl(:,:,:) = 0._wp
115      consfe3 (:,:,:) = 0._wp
116
117      ! Computation of the optimal production rates and nutrient uptake
118      ! rates. Based on a Q10 description of the thermal dependency.
119      zprnut (:,:,:) =  0.8_wp * r1_rday * tgfunc(:,:,:)
120      zprmaxn(:,:,:) =  0.8_wp * (1. + zpsino3 * qnnmax ) * r1_rday * tgfunc(:,:,:)
121      zprmaxd(:,:,:) =  0.8_wp * (1. + zpsino3 * qndmax ) * r1_rday * tgfunc(:,:,:)
122      zprmaxp(:,:,:) =  0.6_wp * (1. + zpsino3 * qnpmax ) * r1_rday * tgfunc(:,:,:)
123
124      ! Impact of the day duration and light intermittency on phytoplankton growth
125      ! Intermittency is supposed to have a similar effect on production as
126      ! day length (Shatwell et al., 2012). The correcting factor is zmxl_fac.
127      ! zmxl_chl is the fractional day length and is used to compute the mean
128      ! PAR during daytime. The effect of mixing is computed using the
129      ! absolute light level definition of the euphotic zone
130      ! -------------------------------------------------------------------------
131      DO jk = 1, jpkm1
132         DO jj = 1 ,jpj
133            DO ji = 1, jpi
134               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
135                  zval = MAX( 1., strn(ji,jj) )
136                  IF( gdepw_n(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
137                     zval = zval * MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn ))
138                  ENDIF
139                  zmxl_chl(ji,jj,jk) = zval / 24.
140                  zmxl_fac(ji,jj,jk) = 1.0 - exp( -0.26 * zval )
141               ENDIF
142            END DO
143         END DO
144      END DO
145
146      zprbio(:,:,:) = zprmaxn(:,:,:) * zmxl_fac(:,:,:)
147      zprdia(:,:,:) = zprmaxd(:,:,:) * zmxl_fac(:,:,:)
148      zprpic(:,:,:) = zprmaxp(:,:,:) * zmxl_fac(:,:,:)
149
150      ! Maximum light intensity
151      zdaylen(:,:) = MAX(1., strn(:,:)) / 24.
152
153      ! Computation of the P-I slope for nanos, picos and diatoms
154      ! The formulation proposed by Geider et al. (1997) has been used.
155      DO jk = 1, jpkm1
156         DO jj = 1, jpj
157            DO ji = 1, jpi
158               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
159                  ztn         = MAX( 0., tsn(ji,jj,jk,jp_tem) - 15. )
160                  zadap       = xadap * ztn / ( 2.+ ztn )
161                  ! Nanophytoplankton
162                  zpislopeadn(ji,jj,jk) = pislopen * trb(ji,jj,jk,jpnch)    &
163                  &                       /( trb(ji,jj,jk,jpphy) * 12. + rtrn)
164
165                  ! Picophytoplankton
166                  zpislopeadp(ji,jj,jk) = pislopep * ( 1. + zadap * EXP( -0.25 * epico(ji,jj,jk) ) )   &
167                  &                       * trb(ji,jj,jk,jppch) /( trb(ji,jj,jk,jppic) * 12. + rtrn)
168
169                  ! Diatoms
170                  zpislopeadd(ji,jj,jk) = pisloped * trb(ji,jj,jk,jpdch)    &
171                     &                    /( trb(ji,jj,jk,jpdia) * 12. + rtrn)
172                  !
173                  zpislopen = zpislopeadn(ji,jj,jk) / ( zprbio(ji,jj,jk) * rday * xlimphy(ji,jj,jk) + rtrn )
174                  zpislopep = zpislopeadp(ji,jj,jk) / ( zprpic(ji,jj,jk) * rday * xlimpic(ji,jj,jk) + rtrn )
175                  zpisloped = zpislopeadd(ji,jj,jk) / ( zprdia(ji,jj,jk) * rday * xlimdia(ji,jj,jk) + rtrn )
176
177                  ! Computation of production function for Carbon
178                  ! Actual light levels are used here
179                  !  ---------------------------------------------
180                  zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) / zmxl_fac(ji,jj,jk) )  )
181                  zprpic(ji,jj,jk) = zprpic(ji,jj,jk) * ( 1.- EXP( -zpislopep * epico(ji,jj,jk) / zmxl_fac(ji,jj,jk) )  )
182                  zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) / zmxl_fac(ji,jj,jk) )  )
183
184                  !  Computation of production function for Chlorophyll
185                  !  Mean light level in the mixed layer (when appropriate)
186                  !  is used here (acclimation is in general slower than
187                  !  the characteristic time scales of vertical mixing)
188                  !  ------------------------------------------------------
189                  zpislopen = zpislopen * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
190                  zpisloped = zpisloped * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
191                  zpislopep = zpislopep * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
192                  zprchln(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) )  )
193                  zprchlp(ji,jj,jk) = zprmaxp(ji,jj,jk) * ( 1.- EXP( -zpislopep * epicom(ji,jj,jk) )  )
194                  zprchld(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) )  )
195               ENDIF
196            END DO
197         END DO
198      END DO
199
200      DO jk = 1, jpkm1
201         DO jj = 1, jpj
202            DO ji = 1, jpi
203
204               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
205                  ! Si/C of diatoms
206                  ! ------------------------
207                  ! Si/C increases with iron stress and silicate availability (zsilfac)
208                  ! Si/C is arbitrariliy increased for very high Si concentrations
209                  ! to mimic the very high ratios observed in the Southern Ocean (zsilfac2)
210                  ! A parameterization derived from Flynn (2003) is used for the control
211                  ! when Si is not limiting which is similar to the parameterisation
212                  ! proposed by Gurney and Davidson (1999).
213                  ! -----------------------------------------------------------------------
214                  zlim  = trb(ji,jj,jk,jpsil) / ( trb(ji,jj,jk,jpsil) + xksi1 )
215                  zsilim = xlimdia(ji,jj,jk) * zprdia(ji,jj,jk) / ( zprmaxd(ji,jj,jk) + rtrn )
216                  zsiborn = trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil)
217                  IF (gphit(ji,jj) < -30 ) THEN
218                    zsilfac2 = 1. + 1. * zsiborn / ( zsiborn + xksi2**3 )
219                  ELSE
220                    zsilfac2 = 1.
221                  ENDIF
222                  zratiosi = 1.0 - trb(ji,jj,jk,jpdsi) / ( trb(ji,jj,jk,jpdia) + rtrn ) / ( zsilfac2 * grosip * 3.0 + rtrn )
223                  zratiosi = MAX(0., MIN(1.0, zratiosi) )
224                  zmaxsi  = (1.0 + 0.1**4) * zratiosi**4 / ( zratiosi**4 + 0.1**4 )
225                  IF ( xlimsi(ji,jj,jk) /= xlimdia(ji,jj,jk) ) THEN
226                     zysopt(ji,jj,jk) = zlim * zsilfac2 * grosip * 1.0 * zmaxsi
227                  ELSE
228                     zysopt(ji,jj,jk) = zlim * zsilfac2 * grosip * 1.0 * zsilim**0.7 * zmaxsi
229                  ENDIF
230               ENDIF
231            END DO
232         END DO
233      END DO
234
235      !  Sea-ice effect on production
236      ! No production is assumed below sea ice
237      ! --------------------------------------
238      DO jk = 1, jpkm1
239         DO jj = 1, jpj
240            DO ji = 1, jpi
241               zprbio(ji,jj,jk)  = zprbio(ji,jj,jk) * ( 1. - fr_i(ji,jj) )
242               zprpic(ji,jj,jk)  = zprpic(ji,jj,jk) * ( 1. - fr_i(ji,jj) ) 
243               zprdia(ji,jj,jk)  = zprdia(ji,jj,jk) * ( 1. - fr_i(ji,jj) ) 
244               zprnut(ji,jj,jk)  = zprnut(ji,jj,jk) * ( 1. - fr_i(ji,jj) )
245            END DO
246         END DO
247      END DO
248
249      ! Computation of the various production and uptake terms of nanophytoplankton
250      ! Interactions between N and P are modeled according to the Chain Model
251      ! of Pahlow et al. (2009). Iron uptake is modeled following traditional
252      ! Droop kinetics. When the quota is approaching the maximum achievable
253      ! quota, uptake is downregulated according to a sigmoidal function
254      ! (power 2), as proposed by Flynn (2003)
255      ! ---------------------------------------------------------------------------
256      DO jk = 1, jpkm1
257         DO jj = 1, jpj
258            DO ji = 1, jpi
259               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
260                  !  production terms for nanophyto.
261                  zprorcan(ji,jj,jk) = zprbio(ji,jj,jk)  * xlimphy(ji,jj,jk) * trb(ji,jj,jk,jpphy) * rfact2
262
263                  ! Size computation
264                  ! Size is made a function of the limitation of of phytoplankton growth
265                  ! Strongly limited cells are supposed to be smaller. sizena is the
266                  ! size at time step t+1 and is thus updated at the end of the
267                  ! current time step
268                  ! --------------------------------------------------------------------
269                  zlimfac = xlimphys(ji,jj,jk) * zprchln(ji,jj,jk) / ( zprmaxn(ji,jj,jk) + rtrn )
270                  zsizetmp = 1.0 + 1.3 * ( xsizern - 1.0 ) * zlimfac**3/(0.3 + zlimfac**3)
271                  sizena(ji,jj,jk) = MIN(xsizern, MAX( sizena(ji,jj,jk), zsizetmp ) )
272                  ! Maximum potential uptake rate
273                  zration = trb(ji,jj,jk,jpnph) / ( trb(ji,jj,jk,jpphy) + rtrn )
274                  zratiop = trb(ji,jj,jk,jppph) / ( trb(ji,jj,jk,jpphy) + rtrn )
275                  zratiof = trb(ji,jj,jk,jpnfe) / ( trb(ji,jj,jk,jpphy) + rtrn )
276                  zprnutmax = zprnut(ji,jj,jk) * fvnuptk(ji,jj,jk) / rno3 * trb(ji,jj,jk,jpphy) * rfact2
277                  ! Uptake of nitrogen
278                  zratio = 1.0 - MIN( 1., zration / (xqnnmax(ji,jj,jk) + rtrn) )
279                  zmax = MAX(0., MIN(1., zratio**2 / (0.05**2 + zratio**2) ) )
280                  zpronmax = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqpnmin(ji,jj,jk) )   &
281                  &          / ( xqpnmax(ji,jj,jk) - xqpnmin(ji,jj,jk) + rtrn ), xlimnfe(ji,jj,jk) ) )
282                  zpronmax = zpronmax * xqnnmin(ji,jj,jk) / qnnmin
283                  zpronewn(ji,jj,jk) = zpronmax * xnanono3(ji,jj,jk)
284                  zproregn(ji,jj,jk) = zpronmax * xnanonh4(ji,jj,jk)
285                  ! Uptake of phosphorus and DOP
286                  zratio = 1.0 - MIN( 1., zratiop / (xqpnmax(ji,jj,jk) + rtrn) )
287                  zmax = MAX(0., MIN(1., zratio**2 / (0.05**2 + zratio**2) ) )
288                  zpropmax = zprnutmax * zmax * xlimnfe(ji,jj,jk)
289                  zpropo4n(ji,jj,jk) = zpropmax * xnanopo4(ji,jj,jk)
290                  zprodopn(ji,jj,jk) = zpropmax * xnanodop(ji,jj,jk)
291                  ! Uptake of iron
292                  zqfnmax = xqfuncfecn(ji,jj,jk) + ( qfnmax - xqfuncfecn(ji,jj,jk) ) * xlimnpn(ji,jj,jk)
293                  zratio = 1.0 - MIN( 1., zratiof / zqfnmax )
294                  zmax = MAX(0., MIN(1., zratio**2/ (0.05**2 + zratio**2) ) )
295                  zprofmax = zprnutmax * zqfnmax * zmax 
296                  zprofen(ji,jj,jk) = zprofmax * xnanofer(ji,jj,jk)    &
297                  &          * (1. + 0.8 * xnanono3(ji,jj,jk) / ( rtrn  &
298                  &          + xnanono3(ji,jj,jk) + xnanonh4(ji,jj,jk) ) * (1. - xnanofer(ji,jj,jk) ) )
299               ENDIF
300            END DO
301         END DO
302      END DO
303
304      ! Computation of the various production and uptake terms of picophytoplankton
305      ! Interactions between N and P are modeled according to the Chain Model
306      ! of Pahlow et al. (2009). Iron uptake is modeled following traditional
307      ! Droop kinetics. When the quota is approaching the maximum achievable
308      ! quota, uptake is downregulated according to a sigmoidal function
309      ! (power 2), as proposed by Flynn (2003)
310      ! ---------------------------------------------------------------------------
311      DO jk = 1, jpkm1
312         DO jj = 1, jpj
313            DO ji = 1, jpi
314               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
315                  !  production terms for picophyto.
316                  zprorcap(ji,jj,jk) = zprpic(ji,jj,jk)  * xlimpic(ji,jj,jk) * trb(ji,jj,jk,jppic) * rfact2
317                  ! Size computation
318                  ! Size is made a function of the limitation of of phytoplankton growth
319                  ! Strongly limited cells are supposed to be smaller. sizepa is
320                  ! size at time step t+1 and is thus updated at the end of the
321                  ! current time step
322                  ! --------------------------------------------------------------------
323                  zlimfac = zprchlp(ji,jj,jk)  * xlimpics(ji,jj,jk) / ( zprmaxp(ji,jj,jk) + rtrn )
324                  zsizetmp = 1.0 + 1.3 * ( xsizerp - 1.0 ) * zlimfac**3/(0.3 + zlimfac**3)
325                  sizepa(ji,jj,jk) = min(xsizerp, max( sizepa(ji,jj,jk), zsizetmp ) )
326                  ! Maximum potential uptake rate of nutrients
327                  zration = trb(ji,jj,jk,jpnpi) / ( trb(ji,jj,jk,jppic) + rtrn )
328                  zratiop = trb(ji,jj,jk,jpppi) / ( trb(ji,jj,jk,jppic) + rtrn )
329                  zratiof = trb(ji,jj,jk,jppfe) / ( trb(ji,jj,jk,jppic) + rtrn )
330                  zprnutmax = zprnut(ji,jj,jk) * fvpuptk(ji,jj,jk) / rno3 * trb(ji,jj,jk,jppic) * rfact2
331                  ! Uptake of nitrogen
332                  zratio = 1.0 - MIN( 1., zration / (xqnpmax(ji,jj,jk) + rtrn) )
333                  zmax = MAX(0., MIN(1., zratio**2/ (0.05**2 + zratio**2) ) )
334                  zpronmax = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqppmin(ji,jj,jk) )   &
335                  &          / ( xqppmax(ji,jj,jk) - xqppmin(ji,jj,jk) + rtrn ), xlimpfe(ji,jj,jk) ) )
336                  zpronmax = zpronmax * xqnpmin(ji,jj,jk) / qnnmin
337                  zpronewp(ji,jj,jk) = zpronmax * xpicono3(ji,jj,jk) 
338                  zproregp(ji,jj,jk) = zpronmax * xpiconh4(ji,jj,jk)
339                  ! Uptake of phosphorus
340                  zratio = 1.0 - MIN( 1., zratiop / (xqppmax(ji,jj,jk) + rtrn) )
341                  zmax = MAX(0., MIN(1., zratio**2 / (0.05**2 + zratio**2) ) )
342                  zpropmax = zprnutmax * zmax * xlimpfe(ji,jj,jk) 
343                  zpropo4p(ji,jj,jk) = zpropmax * xpicopo4(ji,jj,jk)
344                  zprodopp(ji,jj,jk) = zpropmax * xpicodop(ji,jj,jk)
345                  ! Uptake of iron
346                  zqfpmax = xqfuncfecp(ji,jj,jk) + ( qfpmax - xqfuncfecp(ji,jj,jk) ) * xlimnpp(ji,jj,jk)
347                  zratio = 1.0 - MIN( 1., zratiof / zqfpmax )
348                  zmax = MAX(0., MIN(1., zratio**2 / (0.05**2 + zratio**2) ) )
349                  zprofmax = zprnutmax * zqfpmax * zmax
350                  zprofep(ji,jj,jk) = zprofmax * xpicofer(ji,jj,jk)  &
351                  &          * (1. + 0.8 * xpicono3(ji,jj,jk) / ( rtrn   &
352                  &          + xpicono3(ji,jj,jk) + xpiconh4(ji,jj,jk) ) * (1. - xpicofer(ji,jj,jk) ) )
353               ENDIF
354            END DO
355         END DO
356      END DO
357
358      ! Computation of the various production and uptake terms of diatoms
359      ! Interactions between N and P are modeled according to the Chain Model
360      ! of Pahlow et al. (2009). Iron uptake is modeled following traditional
361      ! Droop kinetics. When the quota is approaching the maximum achievable
362      ! quota, uptake is downregulated according to a sigmoidal function
363      ! (power 2), as proposed by Flynn (2003)
364      ! ---------------------------------------------------------------------------
365      DO jk = 1, jpkm1
366         DO jj = 1, jpj
367            DO ji = 1, jpi
368               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
369                  !  production terms for diatomees
370                  zprorcad(ji,jj,jk) = zprdia(ji,jj,jk) * xlimdia(ji,jj,jk) * trb(ji,jj,jk,jpdia) * rfact2
371                  ! Size computation
372                  ! Size is made a function of the limitation of of phytoplankton growth
373                  ! Strongly limited cells are supposed to be smaller. sizeda is
374                  ! size at time step t+1 and is thus updated at the end of the
375                  ! current time step.
376                  ! --------------------------------------------------------------------
377                  zlimfac = zprchld(ji,jj,jk) * xlimdias(ji,jj,jk) / ( zprmaxd(ji,jj,jk) + rtrn )
378                  zsizetmp = 1.0 + 1.3 * ( xsizerd - 1.0 ) * zlimfac**3/(0.3 + zlimfac**3)
379                  sizeda(ji,jj,jk) = min(xsizerd, max( sizeda(ji,jj,jk), zsizetmp ) )
380                  ! Maximum potential uptake rate of nutrients
381                  zration = trb(ji,jj,jk,jpndi) / ( trb(ji,jj,jk,jpdia) + rtrn )
382                  zratiop = trb(ji,jj,jk,jppdi) / ( trb(ji,jj,jk,jpdia) + rtrn )
383                  zratiof = trb(ji,jj,jk,jpdfe) / ( trb(ji,jj,jk,jpdia) + rtrn )
384                  zprnutmax = zprnut(ji,jj,jk) * fvduptk(ji,jj,jk) / rno3 * trb(ji,jj,jk,jpdia) * rfact2
385                  ! Uptake of nitrogen
386                  zratio = 1.0 - MIN( 1., zration / (xqndmax(ji,jj,jk) + rtrn) )
387                  zmax = MAX(0., MIN(1., zratio**2 / (0.05**2 + zratio**2) ) )
388                  zpronmax = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqpdmin(ji,jj,jk) )   &
389                  &          / ( xqpdmax(ji,jj,jk) - xqpdmin(ji,jj,jk) + rtrn ), xlimdfe(ji,jj,jk) ) )
390                  zpronmax = zpronmax * xqndmin(ji,jj,jk) / qnnmin
391                  zpronewd(ji,jj,jk) = zpronmax * xdiatno3(ji,jj,jk)
392                  zproregd(ji,jj,jk) = zpronmax * xdiatnh4(ji,jj,jk)
393                  ! Uptake of phosphorus
394                  zratio = 1.0 - MIN( 1., zratiop / (xqpdmax(ji,jj,jk) + rtrn) )
395                  zmax = MAX(0., MIN(1., zratio**2/ (0.05**2 + zratio**2) ) )
396                  zpropmax = zprnutmax * zmax * xlimdfe(ji,jj,jk)
397                  zpropo4d(ji,jj,jk) = zpropmax * xdiatpo4(ji,jj,jk)
398                  zprodopd(ji,jj,jk) = zpropmax * xdiatdop(ji,jj,jk)
399                  ! Uptake of iron
400                  zqfdmax = xqfuncfecd(ji,jj,jk) + ( qfdmax - xqfuncfecd(ji,jj,jk) ) * xlimnpd(ji,jj,jk)
401                  zratio = 1.0 - MIN( 1., zratiof / zqfdmax )
402                  zmax = MAX(0., MIN(1., zratio**2 / (0.05**2 + zratio**2) ) )
403                  zprofmax = zprnutmax * zqfdmax * zmax
404                  zprofed(ji,jj,jk) = zprofmax * xdiatfer(ji,jj,jk)    &
405                  &          * (1. + 0.8 * xdiatno3(ji,jj,jk) / ( rtrn   &
406                  &          + xdiatno3(ji,jj,jk) + xdiatnh4(ji,jj,jk) ) * (1. - xdiatfer(ji,jj,jk) ) )
407               ENDIF
408            END DO
409         END DO
410      END DO
411
412      ! Production of Chlorophyll. The formulation proposed by Geider et al.
413      ! is adopted here.
414      ! --------------------------------------------------------------------
415      DO jk = 1, jpkm1
416         DO jj = 1, jpj
417            DO ji = 1, jpi
418               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
419                     !  production terms for nanophyto. ( chlorophyll )
420                  znanotot = enanom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
421                  zprod = rday * (zpronewn(ji,jj,jk) + zproregn(ji,jj,jk)) * zprchln(ji,jj,jk) * xlimphy(ji,jj,jk)
422                  zprochln = thetannm * zprod / ( zpislopeadn(ji,jj,jk) * znanotot + rtrn )
423                  zprochln = MAX(zprochln, chlcmin * 12. * zprorcan (ji,jj,jk) )
424                     !  production terms for picophyto. ( chlorophyll )
425                  zpicotot = epicom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
426                  zprod = rday * (zpronewp(ji,jj,jk) + zproregp(ji,jj,jk)) * zprchlp(ji,jj,jk) * xlimpic(ji,jj,jk)
427                  zprochlp = thetanpm * zprod / ( zpislopeadp(ji,jj,jk) * zpicotot + rtrn )
428                  zprochlp = MAX(zprochlp, chlcmin * 12. * zprorcap(ji,jj,jk) )
429                  !  production terms for diatoms ( chlorophyll )
430                  zdiattot = ediatm(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
431                  zprod = rday * (zpronewd(ji,jj,jk) + zproregd(ji,jj,jk)) * zprchld(ji,jj,jk) * xlimdia(ji,jj,jk)
432                  zprochld = thetandm * zprod / ( zpislopeadd(ji,jj,jk) * zdiattot + rtrn )
433                  zprochld = MAX(zprochld, chlcmin * 12. * zprorcad(ji,jj,jk) )
434                  !   Update the arrays TRA which contain the Chla sources and sinks
435                  tra(ji,jj,jk,jpnch) = tra(ji,jj,jk,jpnch) + zprochln * texcretn
436                  tra(ji,jj,jk,jpdch) = tra(ji,jj,jk,jpdch) + zprochld * texcretd
437                  tra(ji,jj,jk,jppch) = tra(ji,jj,jk,jppch) + zprochlp * texcretp
438               ENDIF
439            END DO
440         END DO
441      END DO
442
443      !   Update the arrays TRA which contain the biological sources and sinks
444      DO jk = 1, jpkm1
445         DO jj = 1, jpj
446           DO ji =1 ,jpi
447              zprontot = zpronewn(ji,jj,jk) + zproregn(ji,jj,jk)
448              zproptot = zpronewp(ji,jj,jk) + zproregp(ji,jj,jk)
449              zprodtot = zpronewd(ji,jj,jk) + zproregd(ji,jj,jk)
450              zdocprod = excretd * zprorcad(ji,jj,jk) + excretn * zprorcan(ji,jj,jk)  &
451              &          + excretp * zprorcap(ji,jj,jk)
452              tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) - zpropo4n(ji,jj,jk) - zpropo4d(ji,jj,jk)  &
453              &                     - zpropo4p(ji,jj,jk)
454              tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) - zpronewn(ji,jj,jk) - zpronewd(ji,jj,jk)  &
455              &                     - zpronewp(ji,jj,jk)
456              tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) - zproregn(ji,jj,jk) - zproregd(ji,jj,jk)  &
457              &                     - zproregp(ji,jj,jk)
458              tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) + zprorcan(ji,jj,jk) * texcretn    &
459                 &                  - zpsino3 * zpronewn(ji,jj,jk) - zpsinh4 * zproregn(ji,jj,jk)   &
460                 &                  - zrespn(ji,jj,jk) 
461              zcroissn(ji,jj,jk) = tra(ji,jj,jk,jpphy) / rfact2/ (trb(ji,jj,jk,jpphy) + rtrn)
462              tra(ji,jj,jk,jpnph) = tra(ji,jj,jk,jpnph) + zprontot * texcretn
463              tra(ji,jj,jk,jppph) = tra(ji,jj,jk,jppph) + zpropo4n(ji,jj,jk) * texcretn   &
464              &                     + zprodopn(ji,jj,jk) * texcretn
465              tra(ji,jj,jk,jpnfe) = tra(ji,jj,jk,jpnfe) + zprofen(ji,jj,jk) * texcretn
466              tra(ji,jj,jk,jppic) = tra(ji,jj,jk,jppic) + zprorcap(ji,jj,jk) * texcretp     &
467                 &                  - zpsino3 * zpronewp(ji,jj,jk) - zpsinh4 * zproregp(ji,jj,jk)   &
468                 &                  - zrespp(ji,jj,jk) 
469              zcroissp(ji,jj,jk) = tra(ji,jj,jk,jppic) / rfact2/ (trb(ji,jj,jk,jppic) + rtrn)
470              tra(ji,jj,jk,jpnpi) = tra(ji,jj,jk,jpnpi) + zproptot * texcretp
471              tra(ji,jj,jk,jpppi) = tra(ji,jj,jk,jpppi) + zpropo4p(ji,jj,jk) * texcretp   &
472              &                     + zprodopp(ji,jj,jk) * texcretp
473              tra(ji,jj,jk,jppfe) = tra(ji,jj,jk,jppfe) + zprofep(ji,jj,jk) * texcretp
474              tra(ji,jj,jk,jpdia) = tra(ji,jj,jk,jpdia) + zprorcad(ji,jj,jk) * texcretd   &
475                 &                  - zpsino3 * zpronewd(ji,jj,jk) - zpsinh4 * zproregd(ji,jj,jk)   &
476                 &                  - zrespd(ji,jj,jk) 
477              zcroissd(ji,jj,jk) = tra(ji,jj,jk,jpdia) / rfact2 / (trb(ji,jj,jk,jpdia) + rtrn)
478              tra(ji,jj,jk,jpndi) = tra(ji,jj,jk,jpndi) + zprodtot * texcretd
479              tra(ji,jj,jk,jppdi) = tra(ji,jj,jk,jppdi) + zpropo4d(ji,jj,jk) * texcretd   &
480              &                     + zprodopd(ji,jj,jk) * texcretd
481              tra(ji,jj,jk,jpdfe) = tra(ji,jj,jk,jpdfe) + zprofed(ji,jj,jk) * texcretd
482              tra(ji,jj,jk,jpdsi) = tra(ji,jj,jk,jpdsi) + zprmaxd(ji,jj,jk) * zysopt(ji,jj,jk) * rfact2 * trb(ji,jj,jk,jpdia)
483              tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + excretd * zprorcad(ji,jj,jk) + excretn * zprorcan(ji,jj,jk)  &
484              &                     + excretp * zprorcap(ji,jj,jk)
485              tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) + excretd * zprodtot + excretn * zprontot   &
486              &                     + excretp * zproptot
487              tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) + excretd * zpropo4d(ji,jj,jk) + excretn * zpropo4n(ji,jj,jk)   &
488              &    - texcretn * zprodopn(ji,jj,jk) - texcretd * zprodopd(ji,jj,jk) + excretp * zpropo4p(ji,jj,jk)     &
489              &    - texcretp * zprodopp(ji,jj,jk)
490              tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + o2ut * ( zproregn(ji,jj,jk) + zproregd(ji,jj,jk)   &
491                 &                + zproregp(ji,jj,jk) ) + ( o2ut + o2nit ) * ( zpronewn(ji,jj,jk)           &
492                 &                + zpronewd(ji,jj,jk) + zpronewp(ji,jj,jk) )   &
493                 &                - o2ut * ( zrespn(ji,jj,jk) + zrespp(ji,jj,jk) + zrespd(ji,jj,jk) )
494              zfeup = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk) + texcretp * zprofep(ji,jj,jk)
495              tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zfeup
496              consfe3(ji,jj,jk)   = zfeup * 75.0 / ( rtrn + ( plig(ji,jj,jk) + 75.0 * (1.0 - plig(ji,jj,jk) ) )   &
497              &                     * trb(ji,jj,jk,jpfer) ) / rfact2
498              tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) - zprmaxd(ji,jj,jk) * zysopt(ji,jj,jk) * rfact2 * trb(ji,jj,jk,jpdia)
499              tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) - zprorcan(ji,jj,jk) - zprorcad(ji,jj,jk) - zprorcap(ji,jj,jk)  &
500              &                     + zpsino3 * zpronewn(ji,jj,jk) + zpsinh4 * zproregn(ji,jj,jk)   &
501              &                     + zpsino3 * zpronewp(ji,jj,jk) + zpsinh4 * zproregp(ji,jj,jk)   &
502              &                     + zpsino3 * zpronewd(ji,jj,jk) + zpsinh4 * zproregd(ji,jj,jk)  &
503              &                     + zrespn(ji,jj,jk) + zrespd(ji,jj,jk) + zrespp(ji,jj,jk) 
504              tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zpronewn(ji,jj,jk) + zpronewd(ji,jj,jk)  &
505              &                     + zpronewp(ji,jj,jk) ) - rno3 * ( zproregn(ji,jj,jk) + zproregd(ji,jj,jk)     &
506              &                     + zproregp(ji,jj,jk) ) 
507          END DO
508        END DO
509     END DO
510     
511     ! Production and uptake of ligands by phytoplankton. This part is activated
512     ! when ln_ligand is set to .true. in the namelist. Ligand uptake is small
513     ! and based on the FeL model by Morel et al. (2008) and on the study of
514     ! Shaked and Lis (2012)
515     ! -------------------------------------------------------------------------
516     IF( ln_ligand ) THEN
517         zpligprod1(:,:,:) = 0._wp    ;    zpligprod2(:,:,:) = 0._wp
518         DO jk = 1, jpkm1
519            DO jj = 1, jpj
520              DO ji =1 ,jpi
521                 zdocprod = excretd * zprorcad(ji,jj,jk) + excretn * zprorcan(ji,jj,jk) + excretp * zprorcap(ji,jj,jk)
522                 zfeup    = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk) + texcretp * zprofep(ji,jj,jk)
523                 tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zdocprod * ldocp    &
524                 &       - zfeup * plig(ji,jj,jk) / ( rtrn + plig(ji,jj,jk) + 75.0 * (1.0 - plig(ji,jj,jk) ) ) * lthet
525                 zpligprod1(ji,jj,jk) = zdocprod * ldocp
526                 zpligprod2(ji,jj,jk) = zfeup * plig(ji,jj,jk) / ( rtrn + plig(ji,jj,jk) &
527                 &                      + 75.0 * (1.0 - plig(ji,jj,jk) ) ) * lthet
528              END DO
529           END DO
530        END DO
531     ENDIF
532
533    ! Total primary production per year
534    IF( iom_use( "tintpp" ) .OR. ( ln_check_mass .AND. kt == nitend .AND. knt == nrdttrc )  )  &
535      & tpp = glob_sum( 'p5zprod', ( zprorcan(:,:,:) + zprorcad(:,:,:) + zprorcap(:,:,:) ) * cvol(:,:,:) )
536
537    IF( lk_iomput ) THEN
538       IF( knt == nrdttrc ) THEN
539          ALLOCATE( zw2d(jpi,jpj), zw3d(jpi,jpj,jpk) )
540          zfact = 1.e+3 * rfact2r  !  conversion from mol/l/kt to  mol/m3/s
541          !
542          IF( iom_use( "PPPHYN" ) .OR. iom_use( "PPPHYD" ) .OR. iom_use( "PPPHYP" ) )  THEN
543              zw3d(:,:,:) = zprorcan(:,:,:) * zfact * tmask(:,:,:)  ! primary production by nanophyto
544              CALL iom_put( "PPPHYN"  , zw3d )
545              !
546              zw3d(:,:,:) = zprorcap(:,:,:) * zfact * tmask(:,:,:)  ! primary production by picophyto
547              CALL iom_put( "PPPHYP"  , zw3d )
548              !
549              zw3d(:,:,:) = zprorcad(:,:,:) * zfact * tmask(:,:,:)  ! primary production by diatoms
550              CALL iom_put( "PPPHYD"  , zw3d )
551          ENDIF
552          IF( iom_use( "PPNEWN" ) .OR. iom_use( "PPNEWD" ) .OR. iom_use( "PPNEWP" ) )  THEN
553              zw3d(:,:,:) = zpronewn(:,:,:) * zfact * tmask(:,:,:)  ! new primary production by nanophyto
554              CALL iom_put( "PPNEWN"  , zw3d )
555              !
556              zw3d(:,:,:) = zpronewp(:,:,:) * zfact * tmask(:,:,:)  ! new primary production by picophyto
557              CALL iom_put( "PPNEWP"  , zw3d )
558              !
559              zw3d(:,:,:) = zpronewd(:,:,:) * zfact * tmask(:,:,:)  ! new primary production by diatoms
560              CALL iom_put( "PPNEWD"  , zw3d )
561          ENDIF
562          IF( iom_use( "PBSi" ) )  THEN
563              zw3d(:,:,:) = zprorcad(:,:,:) * zfact * tmask(:,:,:) * zysopt(:,:,:) ! biogenic silica production
564              CALL iom_put( "PBSi"  , zw3d )
565          ENDIF
566          IF( iom_use( "PFeN" ) .OR. iom_use( "PFeD" ) .OR. iom_use( "PFeP" ) )  THEN
567              zw3d(:,:,:) = zprofen(:,:,:) * zfact * tmask(:,:,:)  ! biogenic iron uptake by nanophyto
568              CALL iom_put( "PFeN"  , zw3d )
569              !
570              zw3d(:,:,:) = zprofep(:,:,:) * zfact * tmask(:,:,:)  ! biogenic iron uptake by picophyto
571              CALL iom_put( "PFeP"  , zw3d )
572              !
573              zw3d(:,:,:) = zprofed(:,:,:) * zfact * tmask(:,:,:)  ! biogenic iron uptake by  diatoms
574              CALL iom_put( "PFeD"  , zw3d )
575          ENDIF
576          IF( iom_use( "LPRODP" ) )  THEN
577              zw3d(:,:,:) = zpligprod1(:,:,:) * 1e9 * zfact * tmask(:,:,:)
578              CALL iom_put( "LPRODP"  , zw3d )  ! Ligand production by phytoplankton
579          ENDIF
580          IF( iom_use( "LDETP" ) )  THEN
581              zw3d(:,:,:) = zpligprod2(:,:,:) * 1e9 * zfact * tmask(:,:,:)
582              CALL iom_put( "LDETP"  , zw3d )  ! Uptake of ligands by phytoplankton
583          ENDIF
584          IF( iom_use( "Mumax" ) )  THEN
585              zw3d(:,:,:) = zprmaxn(:,:,:) * tmask(:,:,:)   ! Maximum growth rate
586              CALL iom_put( "Mumax"  , zw3d )
587          ENDIF
588          IF( iom_use( "MuN" ) .OR. iom_use( "MuD" ) .OR. iom_use( "MuP" ) )  THEN
589              zw3d(:,:,:) = zprbio(:,:,:) * xlimphy(:,:,:) * tmask(:,:,:)  ! Realized growth rate for nanophyto
590              CALL iom_put( "MuN"  , zw3d )
591              !
592              zw3d(:,:,:) = zprpic(:,:,:) * xlimpic(:,:,:) * tmask(:,:,:)  ! Realized growth rate for picophyto
593              CALL iom_put( "MuP"  , zw3d )
594              !
595              zw3d(:,:,:) =  zprdia(:,:,:) * xlimdia(:,:,:) * tmask(:,:,:)  ! Realized growth rate for diatoms
596              CALL iom_put( "MuD"  , zw3d )
597          ENDIF
598          IF( iom_use( "LNlight" ) .OR. iom_use( "LDlight" ) .OR. iom_use( "LPlight" ) )  THEN
599              zw3d(:,:,:) = zprbio (:,:,:) / (zprmaxn(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term of nanophytoplankton
600              CALL iom_put( "LNlight"  , zw3d )
601              !
602              zw3d(:,:,:) = zprpic (:,:,:) / (zprmaxp(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term of picophytoplankton
603              CALL iom_put( "LPlight"  , zw3d )
604              !
605              zw3d(:,:,:) =  zprdia (:,:,:) / (zprmaxd(:,:,:) + rtrn) * tmask(:,:,:)  ! light limitation term of diatoms
606              CALL iom_put( "LDlight"  , zw3d )
607          ENDIF
608          IF( iom_use( "MunetN" ) .OR. iom_use( "MunetD" ) .OR. iom_use( "MunetP" ) )  THEN
609              zw3d(:,:,:) = zcroissn(:,:,:) * tmask(:,:,:) ! ! Realized growth rate for nanophyto
610              CALL iom_put( "MunetN"  , zw3d )
611              !
612              zw3d(:,:,:) = zcroissp(:,:,:) * tmask(:,:,:) ! ! Realized growth rate for picophyto
613              CALL iom_put( "MunetP"  , zw3d )
614              !
615              zw3d(:,:,:) = zcroissd(:,:,:) * tmask(:,:,:) ! ! Realized growth rate for diatoms
616              CALL iom_put( "MunetD"  , zw3d )
617              !
618          ENDIF
619
620          IF( iom_use( "tintpp" ) )  CALL iom_put( "tintpp" , tpp * zfact )  !  global total integrated primary production molC/s
621          !
622          DEALLOCATE( zw2d, zw3d )
623       ENDIF
624     ENDIF
625
626      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
627         WRITE(charout, FMT="('prod')")
628         CALL prt_ctl_trc_info(charout)
629         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
630      ENDIF
631      !
632      IF( ln_timing )   CALL timing_stop('p5z_prod')
633      !
634   END SUBROUTINE p5z_prod
635
636
637   SUBROUTINE p5z_prod_init
638      !!----------------------------------------------------------------------
639      !!                  ***  ROUTINE p5z_prod_init  ***
640      !!
641      !! ** Purpose :   Initialization of phytoplankton production parameters
642      !!
643      !! ** Method  :   Read the namp5zprod namelist and check the parameters
644      !!      called at the first timestep (nittrc000)
645      !!
646      !! ** input   :   Namelist namp5zprod
647      !!----------------------------------------------------------------------
648      INTEGER :: ios    ! Local integer output status for namelist read
649      !!
650      NAMELIST/namp5zprod/ pislopen, pislopep, pisloped, excretn, excretp, excretd,     &
651         &                 thetannm, thetanpm, thetandm, chlcmin, grosip, bresp, xadap
652      !!----------------------------------------------------------------------
653
654      REWIND( numnatp_ref )              ! Namelist namp5zprod in reference namelist : Pisces phytoplankton production
655      READ  ( numnatp_ref, namp5zprod, IOSTAT = ios, ERR = 901)
656901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namp5zprod in reference namelist' )
657
658      REWIND( numnatp_cfg )              ! Namelist namp5zprod in configuration namelist : Pisces phytoplankton production
659      READ  ( numnatp_cfg, namp5zprod, IOSTAT = ios, ERR = 902 )
660902   IF( ios >  0 ) CALL ctl_nam ( ios , 'namp5zprod in configuration namelist' )
661      IF(lwm) WRITE ( numonp, namp5zprod )
662
663      IF(lwp) THEN                         ! control print
664         WRITE(numout,*) ' '
665         WRITE(numout,*) ' Namelist parameters for phytoplankton growth, namp5zprod'
666         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
667         WRITE(numout,*) '    mean Si/C ratio                           grosip       =', grosip
668         WRITE(numout,*) '    P-I slope                                 pislopen     =', pislopen
669         WRITE(numout,*) '    P-I slope  for diatoms                    pisloped     =', pisloped
670         WRITE(numout,*) '    P-I slope  for picophytoplankton          pislopep     =', pislopep
671         WRITE(numout,*) '    Acclimation factor to low light           xadap        =', xadap
672         WRITE(numout,*) '    excretion ratio of nanophytoplankton      excretn      =', excretn
673         WRITE(numout,*) '    excretion ratio of picophytoplankton      excretp      =', excretp
674         WRITE(numout,*) '    excretion ratio of diatoms                excretd      =', excretd
675         WRITE(numout,*) '    basal respiration in phytoplankton        bresp        =', bresp
676         WRITE(numout,*) '    Maximum Chl/C in phytoplankton            chlcmin      =', chlcmin
677         WRITE(numout,*) '    Minimum Chl/N in nanophytoplankton        thetannm     =', thetannm
678         WRITE(numout,*) '    Minimum Chl/N in picophytoplankton        thetanpm     =', thetanpm
679         WRITE(numout,*) '    Minimum Chl/N in diatoms                  thetandm     =', thetandm
680      ENDIF
681      !
682      r1_rday   = 1._wp / rday 
683      texcretn  = 1._wp - excretn
684      texcretp  = 1._wp - excretp
685      texcretd  = 1._wp - excretd
686      tpp       = 0._wp
687      !
688   END SUBROUTINE p5z_prod_init
689
690
691   INTEGER FUNCTION p5z_prod_alloc()
692      !!----------------------------------------------------------------------
693      !!                     ***  ROUTINE p5z_prod_alloc  ***
694      !!----------------------------------------------------------------------
695      ALLOCATE( zdaylen(jpi,jpj), STAT = p5z_prod_alloc )
696      !
697      IF( p5z_prod_alloc /= 0 ) CALL ctl_stop( 'STOP', 'p5z_prod_alloc : failed to allocate arrays.' )
698      !
699   END FUNCTION p5z_prod_alloc
700   !!======================================================================
701END MODULE p5zprod
Note: See TracBrowser for help on using the repository browser.