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/2019/dev_r11708_aumont_PISCES_QUOTA/src/TOP/PISCES/P4Z – NEMO

source: NEMO/branches/2019/dev_r11708_aumont_PISCES_QUOTA/src/TOP/PISCES/P4Z/p4zprod.F90 @ 12537

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

Comments in routines have been revised and significantly augmented
  • Property svn:keywords set to Id
File size: 29.7 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_trc      ! 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     !:
29   REAL(wp), PUBLIC ::   pisloped     !:
30   REAL(wp), PUBLIC ::   xadap        !:
31   REAL(wp), PUBLIC ::   excretn      !:
32   REAL(wp), PUBLIC ::   excretd      !:
33   REAL(wp), PUBLIC ::   bresp        !:
34   REAL(wp), PUBLIC ::   chlcnm       !:
35   REAL(wp), PUBLIC ::   chlcdm       !:
36   REAL(wp), PUBLIC ::   chlcmin      !:
37   REAL(wp), PUBLIC ::   fecnm        !:
38   REAL(wp), PUBLIC ::   fecdm        !:
39   REAL(wp), PUBLIC ::   grosip       !:
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   !!----------------------------------------------------------------------
49   !! NEMO/TOP 4.0 , NEMO Consortium (2018)
50   !! $Id$
51   !! Software governed by the CeCILL license (see ./LICENSE)
52   !!----------------------------------------------------------------------
53CONTAINS
54
55   SUBROUTINE p4z_prod( kt , knt )
56      !!---------------------------------------------------------------------
57      !!                     ***  ROUTINE p4z_prod  ***
58      !!
59      !! ** Purpose :   Computes phytoplankton production depending on
60      !!                light, temperature and nutrient availability
61      !!                Computes also the uptake of Iron and Si as well
62      !!                as the chlorophyll content of the cells
63      !!---------------------------------------------------------------------
64      INTEGER, INTENT(in) ::   kt, knt   !
65      !
66      INTEGER  ::   ji, jj, jk
67      REAL(wp) ::   zsilfac, znanotot, zdiattot, zconctemp, zconctemp2
68      REAL(wp) ::   zratio, zmax, zsilim, ztn, zadap, zlim, zsilfac2, zsiborn
69      REAL(wp) ::   zprod, zproreg, zproreg2, zprochln, zprochld
70      REAL(wp) ::   zmaxday, zdocprod, zpislopen, zpisloped
71      REAL(wp) ::   zmxltst, zmxlday
72      REAL(wp) ::   zrum, zcodel, zargu, zval, zfeup, chlcnm_n, chlcdm_n
73      REAL(wp) ::   zfact
74      CHARACTER (len=25) :: charout
75      REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zw2d
76      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d
77      REAL(wp), DIMENSION(jpi,jpj    ) :: zstrn, zmixnano, zmixdiat
78      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprmaxn,zprmaxd
79      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpislopeadn, zpislopeadd, zysopt 
80      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprdia, zprbio, zprchld, zprchln   
81      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprorcan, zprorcad, zprofed, zprofen
82      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpronewn, zpronewd
83      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmxl_fac, zmxl_chl
84      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpligprod1, zpligprod2
85      !!---------------------------------------------------------------------
86      !
87      IF( ln_timing )   CALL timing_start('p4z_prod')
88      !
89      !  Allocate temporary workspace
90      !
91      zprorcan(:,:,:) = 0._wp ; zprorcad(:,:,:) = 0._wp ; zprofed (:,:,:) = 0._wp
92      zprofen (:,:,:) = 0._wp ; zysopt  (:,:,:) = 0._wp
93      zpronewn(:,:,:) = 0._wp ; zpronewd(:,:,:) = 0._wp ; zprdia  (:,:,:) = 0._wp
94      zprbio  (:,:,:) = 0._wp ; zprchld (:,:,:) = 0._wp ; zprchln (:,:,:) = 0._wp 
95      zmxl_fac(:,:,:) = 0._wp ; zmxl_chl(:,:,:) = 0._wp 
96
97      ! Computation of the maximimum production
98      ! Parameters are taken from Bissinger et al. (2008)
99      zprmaxn(:,:,:) = 0.8_wp * r1_rday * tgfunc(:,:,:)
100      zprmaxd(:,:,:) = zprmaxn(:,:,:)
101
102      ! compute the day length depending on latitude and the day
103      zrum = REAL( nday_year - 80, wp ) / REAL( nyear_len(1), wp )
104      zcodel = ASIN(  SIN( zrum * rpi * 2._wp ) * SIN( rad * 23.5_wp )  )
105
106      ! day length in hours
107      zstrn(:,:) = 0.
108      DO jj = 1, jpj
109         DO ji = 1, jpi
110            zargu = TAN( zcodel ) * TAN( gphit(ji,jj) * rad )
111            zargu = MAX( -1., MIN(  1., zargu ) )
112            zstrn(ji,jj) = MAX( 0.0, 24. - 2. * ACOS( zargu ) / rad / 15. )
113         END DO
114      END DO
115
116      ! Impact of the day duration and light intermittency on phytoplankton growth
117      ! Intermittency is supposed to have a similar effect on production as
118      ! day length. The correcting factor is zmxl_fac. zmxl_chl is the fractional
119      ! day length and is used to compute the mean PAR during daytime.
120      ! Formulation for the impact of day length on PP is from Thompson (1999)
121      ! -------------------------------------------------------------------------
122      DO jk = 1, jpkm1
123         DO jj = 1 ,jpj
124            DO ji = 1, jpi
125               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
126                  zval = MAX( 1., zstrn(ji,jj) )
127                  IF( gdepw_n(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
128                     zval = zval * MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn ))
129                  ENDIF
130                  zmxl_chl(ji,jj,jk) = zval / 24.
131                  zmxl_fac(ji,jj,jk) = 1.5 * zval / ( 12. + zval )
132               ENDIF
133            END DO
134         END DO
135      END DO
136
137      zprbio(:,:,:) = zprmaxn(:,:,:) * zmxl_fac(:,:,:)
138      zprdia(:,:,:) = zprmaxd(:,:,:) * zmxl_fac(:,:,:)
139
140      WHERE( zstrn(:,:) < 1.e0 ) zstrn(:,:) = 24.
141
142      ! Computation of the P-I slope for nanos and diatoms
143      ! The formulation proposed by Geider et al. (1997) has been modified
144      ! to exclude the effect of nutrient limitation and temperature in the PI
145      ! curve following Vichi et al. (2007)
146      ! -----------------------------------------------------------------------
147      DO jk = 1, jpkm1
148         DO jj = 1, jpj
149            DO ji = 1, jpi
150               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
151                  ztn         = MAX( 0., tsn(ji,jj,jk,jp_tem) - 15. )
152                  zadap       = xadap * ztn / ( 2.+ ztn )
153                  zconctemp   = MAX( 0.e0 , trb(ji,jj,jk,jpdia) - xsizedia )
154                  zconctemp2  = trb(ji,jj,jk,jpdia) - zconctemp
155                  !
156                  ! The initial slope of the PI curve can be increased for nano
157                  ! to account for photadaptation, for instance in the DCM
158                  zpislopeadn(ji,jj,jk) = pislopen * ( 1.+ zadap  * EXP( -0.25 * enano(ji,jj,jk) ) )  &
159                  &                   * trb(ji,jj,jk,jpnch) /( trb(ji,jj,jk,jpphy) * 12. + rtrn)
160                  !
161                  zpislopeadd(ji,jj,jk) = (pislopen * zconctemp2 + pisloped * zconctemp) / ( trb(ji,jj,jk,jpdia) + rtrn )   &
162                  &                   * trb(ji,jj,jk,jpdch) /( trb(ji,jj,jk,jpdia) * 12. + rtrn)
163               ENDIF
164            END DO
165         END DO
166      END DO
167
168      DO jk = 1, jpkm1
169         DO jj = 1, jpj
170            DO ji = 1, jpi
171               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
172                   ! Computation of production function for Carbon
173                   !  ---------------------------------------------
174                   zpislopen = zpislopeadn(ji,jj,jk) / ( ( r1_rday + bresp * r1_rday ) &
175                   &            * zmxl_fac(ji,jj,jk) * rday + rtrn)
176                   zpisloped = zpislopeadd(ji,jj,jk) / ( ( r1_rday + bresp * r1_rday ) &
177                   &            * zmxl_fac(ji,jj,jk) * rday + rtrn)
178                   zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) )  )
179                   zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) )  )
180                   !  Computation of production function for Chlorophyll
181                   !--------------------------------------------------
182                   zpislopen = zpislopeadn(ji,jj,jk) / ( zprmaxn(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn )
183                   zpisloped = zpislopeadd(ji,jj,jk) / ( zprmaxd(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn )
184                   zprchln(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) ) )
185                   zprchld(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) ) )
186               ENDIF
187            END DO
188         END DO
189      END DO
190
191      !  Computation of a proxy of the N/C ratio
192      !  Steady state is assumed
193      !  ---------------------------------------
194      DO jk = 1, jpkm1
195         DO jj = 1, jpj
196            DO ji = 1, jpi
197                zval = MIN( xnanopo4(ji,jj,jk), ( xnanonh4(ji,jj,jk) + xnanono3(ji,jj,jk) ) )   &
198                &      * zprmaxn(ji,jj,jk) / ( zprbio(ji,jj,jk) + rtrn )
199                quotan(ji,jj,jk) = MIN( 1., 0.2 + 0.8 * zval )
200                zval = MIN( xdiatpo4(ji,jj,jk), ( xdiatnh4(ji,jj,jk) + xdiatno3(ji,jj,jk) ) )   &
201                &      * zprmaxd(ji,jj,jk) / ( zprdia(ji,jj,jk) + rtrn )
202                quotad(ji,jj,jk) = MIN( 1., 0.2 + 0.8 * zval )
203            END DO
204         END DO
205      END DO
206
207
208      DO jk = 1, jpkm1
209         DO jj = 1, jpj
210            DO ji = 1, jpi
211
212                IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
213                   ! Si/C of diatoms
214                   ! ------------------------
215                   ! Si/C increases with iron stress and silicate availability (zsilfac)
216                   ! Si/C is arbitrariliy increased for very high Si concentrations
217                   ! to mimic the very high ratios observed in the Southern Ocean (zsilfac2)
218                   ! -----------------------------------------------------------------------
219                  zlim  = trb(ji,jj,jk,jpsil) / ( trb(ji,jj,jk,jpsil) + xksi1 )
220                  zsilim = MIN( zprdia(ji,jj,jk) / ( zprmaxd(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) )
221                  zsilfac = 4.4 * EXP( -4.23 * zsilim ) * MAX( 0.e0, MIN( 1., 2.2 * ( zlim - 0.5 ) )  ) + 1.e0
222                  zsiborn = trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil)
223                  IF (gphit(ji,jj) < -30 ) THEN
224                    zsilfac2 = 1. + 2. * zsiborn / ( zsiborn + xksi2**3 )
225                  ELSE
226                    zsilfac2 = 1. +      zsiborn / ( zsiborn + xksi2**3 )
227                  ENDIF
228                  zysopt(ji,jj,jk) = grosip * zlim * zsilfac * zsilfac2
229              ENDIF
230            END DO
231         END DO
232      END DO
233
234      ! Sea-ice effect on production
235      ! No production is assumed below sea ice
236      ! --------------------------------------
237      DO jk = 1, jpkm1
238         DO jj = 1, jpj
239            DO ji = 1, jpi
240               zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1. - fr_i(ji,jj) )
241               zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1. - fr_i(ji,jj) )
242            END DO
243         END DO
244      END DO
245
246      ! Computation of the various production  and nutrient uptake terms
247      ! ---------------------------------------------------------------
248      DO jk = 1, jpkm1
249         DO jj = 1, jpj
250            DO ji = 1, jpi
251               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
252                  !  production term of nanophyto. (C)
253                  zprorcan(ji,jj,jk) = zprbio(ji,jj,jk)  * xlimphy(ji,jj,jk) * trb(ji,jj,jk,jpphy) * rfact2
254
255                  !  New production (uptake of NO3)
256                  zpronewn(ji,jj,jk)  = zprorcan(ji,jj,jk)* xnanono3(ji,jj,jk) / ( xnanono3(ji,jj,jk) + xnanonh4(ji,jj,jk) + rtrn )
257                  !
258                  zratio = trb(ji,jj,jk,jpnfe) / ( trb(ji,jj,jk,jpphy) * fecnm + rtrn )
259
260                  !  Iron uptake rates of nanophytoplankton. Upregulation 
261                  !  is parameterized at low iron concentrations. Typical
262                  !  formulation used in quota formulations. Uptake is downregulated
263                  !  when the quota is close to the maximum quota
264                  zmax   = MAX( 0., ( 1. - zratio ) / ABS( 1.05 - zratio ) ) 
265                  zprofen(ji,jj,jk) = fecnm * zprmaxn(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  &
266                  &             * ( 4. - 4.5 * xlimnfe(ji,jj,jk) / ( xlimnfe(ji,jj,jk) + 0.5 ) )    &
267                  &             * biron(ji,jj,jk) / ( biron(ji,jj,jk) + concnfe(ji,jj,jk) )  &
268                  &             * zmax * trb(ji,jj,jk,jpphy) * rfact2
269                  !  production terms of diatoms (C)
270                  zprorcad(ji,jj,jk) = zprdia(ji,jj,jk) * xlimdia(ji,jj,jk) * trb(ji,jj,jk,jpdia) * rfact2
271
272                  ! New production (uptake of NO3)
273                  zpronewd(ji,jj,jk) = zprorcad(ji,jj,jk) * xdiatno3(ji,jj,jk) / ( xdiatno3(ji,jj,jk) + xdiatnh4(ji,jj,jk) + rtrn )
274
275                  !  Iron uptake rates of nanophytoplankton. Upregulation 
276                  !  is parameterized at low iron concentrations. Typical
277                  !  formulation used in quota formulations. Uptake is downregulated
278                  !  when the quota is close to the maximum quota
279                  zratio = trb(ji,jj,jk,jpdfe) / ( trb(ji,jj,jk,jpdia) * fecdm + rtrn )
280                  zmax   = MAX( 0., ( 1. - zratio ) / ABS( 1.05 - zratio ) ) 
281                  zprofed(ji,jj,jk) = fecdm * zprmaxd(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  &
282                  &             * ( 4. - 4.5 * xlimdfe(ji,jj,jk) / ( xlimdfe(ji,jj,jk) + 0.5 ) )    &
283                  &             * biron(ji,jj,jk) / ( biron(ji,jj,jk) + concdfe(ji,jj,jk) )  &
284                  &             * zmax * trb(ji,jj,jk,jpdia) * rfact2
285               ENDIF
286            END DO
287         END DO
288      END DO
289
290      ! Computation of the chlorophyll production terms
291      ! The parameterization is taken from Geider et al. (1997)
292      ! -------------------------------------------------------
293      DO jk = 1, jpkm1
294         DO jj = 1, jpj
295            DO ji = 1, jpi
296               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
297                  !  production term for nanophyto. ( chlorophyll )
298                  znanotot = enanom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
299                  zprod    = rday * zprorcan(ji,jj,jk) * zprchln(ji,jj,jk) * xlimphy(ji,jj,jk)
300                  zprochln = chlcmin * 12. * zprorcan (ji,jj,jk)
301                  ! The maximum reachable Chl quota is modulated by temperature
302                  ! following Geider (1987)
303                  chlcnm_n   = MIN ( chlcnm, ( chlcnm / (1. - 1.14 / 43.4 *tsn(ji,jj,jk,jp_tem))) * (1. - 1.14 / 43.4 * 20.))
304                  zprochln = zprochln + (chlcnm_n-chlcmin) * 12. * zprod / &
305                                        & (  zpislopeadn(ji,jj,jk) * znanotot +rtrn)
306                  !  production terms of diatoms ( chlorophyll )
307                  zdiattot = ediatm(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn )
308                  zprod    = rday * zprorcad(ji,jj,jk) * zprchld(ji,jj,jk) * xlimdia(ji,jj,jk)
309                  zprochld = chlcmin * 12. * zprorcad(ji,jj,jk)
310                  ! The maximum reachable Chl quota is modulated by temperature
311                  ! following Geider (1987)
312                  chlcdm_n   = MIN ( chlcdm, ( chlcdm / (1. - 1.14 / 43.4 * tsn(ji,jj,jk,jp_tem))) * (1. - 1.14 / 43.4 * 20.))
313                  zprochld = zprochld + (chlcdm_n-chlcmin) * 12. * zprod / &
314                                        & ( zpislopeadd(ji,jj,jk) * zdiattot +rtrn )
315                  !   Update the arrays TRA which contain the Chla sources and sinks
316                  tra(ji,jj,jk,jpnch) = tra(ji,jj,jk,jpnch) + zprochln * texcretn
317                  tra(ji,jj,jk,jpdch) = tra(ji,jj,jk,jpdch) + zprochld * texcretd
318               ENDIF
319            END DO
320         END DO
321      END DO
322
323      !   Update the arrays TRA which contain the biological sources and sinks
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                 zproreg  = zprorcan(ji,jj,jk) - zpronewn(ji,jj,jk)
329                 zproreg2 = zprorcad(ji,jj,jk) - zpronewd(ji,jj,jk)
330                 zdocprod = excretd * zprorcad(ji,jj,jk) + excretn * zprorcan(ji,jj,jk)
331                 tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) - zprorcan(ji,jj,jk) - zprorcad(ji,jj,jk)
332                 tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) - zpronewn(ji,jj,jk) - zpronewd(ji,jj,jk)
333                 tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) - zproreg - zproreg2
334                 tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) + zprorcan(ji,jj,jk) * texcretn
335                 tra(ji,jj,jk,jpnfe) = tra(ji,jj,jk,jpnfe) + zprofen(ji,jj,jk) * texcretn
336                 tra(ji,jj,jk,jpdia) = tra(ji,jj,jk,jpdia) + zprorcad(ji,jj,jk) * texcretd
337                 tra(ji,jj,jk,jpdfe) = tra(ji,jj,jk,jpdfe) + zprofed(ji,jj,jk) * texcretd
338                 tra(ji,jj,jk,jpdsi) = tra(ji,jj,jk,jpdsi) + zprorcad(ji,jj,jk) * zysopt(ji,jj,jk) * texcretd
339                 tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zdocprod
340                 tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + o2ut * ( zproreg + zproreg2) &
341                 &                   + ( o2ut + o2nit ) * ( zpronewn(ji,jj,jk) + zpronewd(ji,jj,jk) )
342                 !
343                 zfeup = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk)
344                 tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zfeup
345                 tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) - texcretd * zprorcad(ji,jj,jk) * zysopt(ji,jj,jk)
346                 tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) - zprorcan(ji,jj,jk) - zprorcad(ji,jj,jk)
347                 tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zpronewn(ji,jj,jk) + zpronewd(ji,jj,jk) ) &
348                 &                                         - rno3 * ( zproreg + zproreg2 )
349              ENDIF
350           END DO
351        END DO
352     END DO
353     !
354     ! Production and uptake of ligands by phytoplankton. This part is activated
355     ! when ln_ligand is set to .true. in the namelist. Ligand uptake is small
356     ! and based on the FeL model by Morel et al. (2008) and on the study of
357     ! Shaked and Lis (2012)
358     ! -------------------------------------------------------------------------
359     IF( ln_ligand ) THEN
360         zpligprod1(:,:,:) = 0._wp    ;    zpligprod2(:,:,:) = 0._wp
361         DO jk = 1, jpkm1
362            DO jj = 1, jpj
363              DO ji =1 ,jpi
364                 IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN
365                    zdocprod = excretd * zprorcad(ji,jj,jk) + excretn * zprorcan(ji,jj,jk)
366                    zfeup    = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk)
367                    tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zdocprod * ldocp    &
368                    &       - zfeup * plig(ji,jj,jk) / ( rtrn + plig(ji,jj,jk) + 2.E3 * (1.0 - plig(ji,jj,jk) ) ) * lthet
369                    zpligprod1(ji,jj,jk) = zdocprod * ldocp
370                    zpligprod2(ji,jj,jk) = zfeup * plig(ji,jj,jk) / ( rtrn + plig(ji,jj,jk) &
371                    &                      + 2.E3 * (1.0 - plig(ji,jj,jk) ) ) * lthet
372                 ENDIF
373              END DO
374           END DO
375        END DO
376     ENDIF
377
378
379    ! Total primary production per year
380    IF( iom_use( "tintpp" ) .OR. ( ln_check_mass .AND. kt == nitend .AND. knt == nrdttrc )  )  &
381         & tpp = glob_sum( 'p4zprod', ( zprorcan(:,:,:) + zprorcad(:,:,:) ) * cvol(:,:,:) )
382
383    IF( lk_iomput ) THEN
384       IF( knt == nrdttrc ) THEN
385          ALLOCATE( zw2d(jpi,jpj), zw3d(jpi,jpj,jpk) )
386          zfact = 1.e+3 * rfact2r  !  conversion from mol/l/kt to  mol/m3/s
387          !
388          IF( iom_use( "PPPHYN" ) .OR. iom_use( "PPPHYD" ) )  THEN
389              zw3d(:,:,:) = zprorcan(:,:,:) * zfact * tmask(:,:,:)  ! primary production by nanophyto
390              CALL iom_put( "PPPHYN"  , zw3d )
391              !
392              zw3d(:,:,:) = zprorcad(:,:,:) * zfact * tmask(:,:,:)  ! primary production by diatomes
393              CALL iom_put( "PPPHYD"  , zw3d )
394          ENDIF
395          IF( iom_use( "PPNEWN" ) .OR. iom_use( "PPNEWD" ) )  THEN
396              zw3d(:,:,:) = zpronewn(:,:,:) * zfact * tmask(:,:,:)  ! new primary production by nanophyto
397              CALL iom_put( "PPNEWN"  , zw3d )
398              !
399              zw3d(:,:,:) = zpronewd(:,:,:) * zfact * tmask(:,:,:)  ! new primary production by diatomes
400              CALL iom_put( "PPNEWD"  , zw3d )
401          ENDIF
402          IF( iom_use( "PBSi" ) )  THEN
403              zw3d(:,:,:) = zprorcad(:,:,:) * zfact * tmask(:,:,:) * zysopt(:,:,:) ! biogenic silica production
404              CALL iom_put( "PBSi"  , zw3d )
405          ENDIF
406          IF( iom_use( "PFeN" ) .OR. iom_use( "PFeD" ) )  THEN
407              zw3d(:,:,:) = zprofen(:,:,:) * zfact * tmask(:,:,:)  ! biogenic iron production by nanophyto
408              CALL iom_put( "PFeN"  , zw3d )
409              !
410              zw3d(:,:,:) = zprofed(:,:,:) * zfact * tmask(:,:,:)  ! biogenic iron production by  diatomes
411              CALL iom_put( "PFeD"  , zw3d )
412          ENDIF
413          IF( iom_use( "LPRODP" ) )  THEN
414              zw3d(:,:,:) = zpligprod1(:,:,:) * 1e9 * zfact * tmask(:,:,:)  ! Ligand production by phytoplankton
415              CALL iom_put( "LPRODP"  , zw3d )
416          ENDIF
417          IF( iom_use( "LDETP" ) )  THEN
418              zw3d(:,:,:) = zpligprod2(:,:,:) * 1e9 * zfact * tmask(:,:,:)  ! Uptake of ligands by phytoplankton
419              CALL iom_put( "LDETP"  , zw3d )
420          ENDIF
421          IF( iom_use( "Mumax" ) )  THEN
422              zw3d(:,:,:) = zprmaxn(:,:,:) * tmask(:,:,:)   ! Maximum growth rate
423              CALL iom_put( "Mumax"  , zw3d )
424          ENDIF
425          IF( iom_use( "MuN" ) .OR. iom_use( "MuD" ) )  THEN
426              zw3d(:,:,:) = zprbio(:,:,:) * xlimphy(:,:,:) * tmask(:,:,:)  ! Realized growth rate for nanophyto
427              CALL iom_put( "MuN"  , zw3d )
428              !
429              zw3d(:,:,:) =  zprdia(:,:,:) * xlimdia(:,:,:) * tmask(:,:,:)  ! Realized growth rate for diatoms
430              CALL iom_put( "MuD"  , zw3d )
431          ENDIF
432          IF( iom_use( "LNlight" ) .OR. iom_use( "LDlight" ) )  THEN
433              zw3d(:,:,:) = zprbio (:,:,:) / (zprmaxn(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term of nanophytoplankton
434              CALL iom_put( "LNlight"  , zw3d )
435              !
436              zw3d(:,:,:) = zprdia (:,:,:) / (zprmaxd(:,:,:) + rtrn) * tmask(:,:,:)  ! light limitation term of diatoms
437              CALL iom_put( "LDlight"  , zw3d )
438          ENDIF
439          IF( iom_use( "TPP" ) )  THEN
440              zw3d(:,:,:) = ( zprorcan(:,:,:) + zprorcad(:,:,:) ) * zfact * tmask(:,:,:)  ! total primary production
441              CALL iom_put( "TPP"  , zw3d )
442          ENDIF
443          IF( iom_use( "TPNEW" ) )  THEN
444              zw3d(:,:,:) = ( zpronewn(:,:,:) + zpronewd(:,:,:) ) * zfact * tmask(:,:,:)  ! total new production
445              CALL iom_put( "TPNEW"  , zw3d )
446          ENDIF
447          IF( iom_use( "TPBFE" ) )  THEN
448              zw3d(:,:,:) = ( zprofen(:,:,:) + zprofed(:,:,:) ) * zfact * tmask(:,:,:)  ! total biogenic iron production
449              CALL iom_put( "TPBFE"  , zw3d )
450          ENDIF
451          IF( iom_use( "INTPPPHYN" ) .OR. iom_use( "INTPPPHYD" ) ) THEN 
452             zw2d(:,:) = 0.
453             DO jk = 1, jpkm1
454               zw2d(:,:) = zw2d(:,:) + zprorcan(:,:,jk) * e3t_n(:,:,jk) * zfact * tmask(:,:,jk)  ! vert. integrated  primary produc. by nano
455             ENDDO
456             CALL iom_put( "INTPPPHYN" , zw2d )
457             !
458             zw2d(:,:) = 0.
459             DO jk = 1, jpkm1
460                zw2d(:,:) = zw2d(:,:) + zprorcad(:,:,jk) * e3t_n(:,:,jk) * zfact * tmask(:,:,jk) ! vert. integrated  primary produc. by diatom
461             ENDDO
462             CALL iom_put( "INTPPPHYD" , zw2d )
463          ENDIF
464          IF( iom_use( "INTPP" ) ) THEN   
465             zw2d(:,:) = 0.
466             DO jk = 1, jpkm1
467                zw2d(:,:) = zw2d(:,:) + ( zprorcan(:,:,jk) + zprorcad(:,:,jk) ) * e3t_n(:,:,jk) * zfact * tmask(:,:,jk) ! vert. integrated pp
468             ENDDO
469             CALL iom_put( "INTPP" , zw2d )
470          ENDIF
471          IF( iom_use( "INTPNEW" ) ) THEN   
472             zw2d(:,:) = 0.
473             DO jk = 1, jpkm1
474                zw2d(:,:) = zw2d(:,:) + ( zpronewn(:,:,jk) + zpronewd(:,:,jk) ) * e3t_n(:,:,jk) * zfact * tmask(:,:,jk)  ! vert. integrated new prod
475             ENDDO
476             CALL iom_put( "INTPNEW" , zw2d )
477          ENDIF
478          IF( iom_use( "INTPBFE" ) ) THEN           !   total biogenic iron production  ( vertically integrated )
479             zw2d(:,:) = 0.
480             DO jk = 1, jpkm1
481                zw2d(:,:) = zw2d(:,:) + ( zprofen(:,:,jk) + zprofed(:,:,jk) ) * e3t_n(:,:,jk) * zfact * tmask(:,:,jk) ! vert integr. bfe prod
482             ENDDO
483            CALL iom_put( "INTPBFE" , zw2d )
484          ENDIF
485          IF( iom_use( "INTPBSI" ) ) THEN           !   total biogenic silica production  ( vertically integrated )
486             zw2d(:,:) = 0.
487             DO jk = 1, jpkm1
488                zw2d(:,:) = zw2d(:,:) + zprorcad(:,:,jk) * zysopt(:,:,jk) * e3t_n(:,:,jk) * zfact * tmask(:,:,jk)  ! vert integr. bsi prod
489             ENDDO
490             CALL iom_put( "INTPBSI" , zw2d )
491          ENDIF
492          IF( iom_use( "tintpp" ) )  CALL iom_put( "tintpp" , tpp * zfact )  !  global total integrated primary production molC/s
493          !
494          DEALLOCATE( zw2d, zw3d )
495       ENDIF
496     ENDIF
497
498     IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
499         WRITE(charout, FMT="('prod')")
500         CALL prt_ctl_trc_info(charout)
501         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
502     ENDIF
503      !
504      IF( ln_timing )  CALL timing_stop('p4z_prod')
505      !
506   END SUBROUTINE p4z_prod
507
508
509   SUBROUTINE p4z_prod_init
510      !!----------------------------------------------------------------------
511      !!                  ***  ROUTINE p4z_prod_init  ***
512      !!
513      !! ** Purpose :   Initialization of phytoplankton production parameters
514      !!
515      !! ** Method  :   Read the namp4zprod namelist and check the parameters
516      !!      called at the first timestep (nittrc000)
517      !!
518      !! ** input   :   Namelist namp4zprod
519      !!----------------------------------------------------------------------
520      INTEGER ::   ios   ! Local integer
521      !
522      NAMELIST/namp4zprod/ pislopen, pisloped, xadap, bresp, excretn, excretd,  &
523         &                 chlcnm, chlcdm, chlcmin, fecnm, fecdm, grosip
524      !!----------------------------------------------------------------------
525      !
526      IF(lwp) THEN                         ! control print
527         WRITE(numout,*)
528         WRITE(numout,*) 'p4z_prod_init : phytoplankton growth'
529         WRITE(numout,*) '~~~~~~~~~~~~~'
530      ENDIF
531      !
532      REWIND( numnatp_ref )              ! Namelist namp4zprod in reference namelist : Pisces phytoplankton production
533      READ  ( numnatp_ref, namp4zprod, IOSTAT = ios, ERR = 901)
534901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namp4zprod in reference namelist' )
535      REWIND( numnatp_cfg )              ! Namelist namp4zprod in configuration namelist : Pisces phytoplankton production
536      READ  ( numnatp_cfg, namp4zprod, IOSTAT = ios, ERR = 902 )
537902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namp4zprod in configuration namelist' )
538      IF(lwm) WRITE( numonp, namp4zprod )
539
540      IF(lwp) THEN                         ! control print
541         WRITE(numout,*) '   Namelist : namp4zprod'
542         WRITE(numout,*) '      mean Si/C ratio                           grosip       =', grosip
543         WRITE(numout,*) '      P-I slope                                 pislopen     =', pislopen
544         WRITE(numout,*) '      Acclimation factor to low light           xadap        =', xadap
545         WRITE(numout,*) '      excretion ratio of nanophytoplankton      excretn      =', excretn
546         WRITE(numout,*) '      excretion ratio of diatoms                excretd      =', excretd
547         WRITE(numout,*) '      basal respiration in phytoplankton        bresp        =', bresp
548         WRITE(numout,*) '      Maximum Chl/C in phytoplankton            chlcmin      =', chlcmin
549         WRITE(numout,*) '      P-I slope  for diatoms                    pisloped     =', pisloped
550         WRITE(numout,*) '      Minimum Chl/C in nanophytoplankton        chlcnm       =', chlcnm
551         WRITE(numout,*) '      Minimum Chl/C in diatoms                  chlcdm       =', chlcdm
552         WRITE(numout,*) '      Maximum Fe/C in nanophytoplankton         fecnm        =', fecnm
553         WRITE(numout,*) '      Minimum Fe/C in diatoms                   fecdm        =', fecdm
554      ENDIF
555      !
556      r1_rday   = 1._wp / rday 
557      texcretn  = 1._wp - excretn
558      texcretd  = 1._wp - excretd
559      tpp       = 0._wp
560      !
561   END SUBROUTINE p4z_prod_init
562
563
564   INTEGER FUNCTION p4z_prod_alloc()
565      !!----------------------------------------------------------------------
566      !!                     ***  ROUTINE p4z_prod_alloc  ***
567      !!----------------------------------------------------------------------
568      ALLOCATE( quotan(jpi,jpj,jpk), quotad(jpi,jpj,jpk), STAT = p4z_prod_alloc )
569      !
570      IF( p4z_prod_alloc /= 0 ) CALL ctl_stop( 'STOP', 'p4z_prod_alloc : failed to allocate arrays.' )
571      !
572   END FUNCTION p4z_prod_alloc
573
574   !!======================================================================
575END MODULE p4zprod
Note: See TracBrowser for help on using the repository browser.