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 @ 13233

Last change on this file since 13233 was 13233, checked in by aumont, 4 years ago

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