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 branches/CNRS/dev_r4826_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z – NEMO

source: branches/CNRS/dev_r4826_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p5zprod.F90 @ 5266

Last change on this file since 5266 was 5266, checked in by cetlod, 9 years ago

PISCES_QUOTA : First commits, see ticket #1516
  • Property svn:executable set to *
File size: 35.3 KB
Line 
1MODULE p5zprod
2   !!======================================================================
3   !!                         ***  MODULE p5zprod  ***
4   !! TOP :  Growth Rate of the two phytoplanktons groups
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#if defined key_pisces_quota
11   !!----------------------------------------------------------------------
12   !!   'key_pisces_quota'     PISCES bio-model with variable stoichiometry
13   !!----------------------------------------------------------------------
14   !!   p5z_prod       :   Compute the growth Rate of the two phytoplanktons groups
15   !!   p5z_prod_init  :   Initialization of the parameters for growth
16   !!   p5z_prod_alloc :   Allocate variables for growth
17   !!----------------------------------------------------------------------
18   USE oce_trc         !  shared variables between ocean and passive tracers
19   USE trc             !  passive tracers common variables
20   USE sms_pisces      !  PISCES Source Minus Sink variables
21   USE p4zopt          !  optical model
22   USE p5zlim          !  Co-limitations of differents nutrients
23   USE prtctl_trc      !  print control for debugging
24   USE iom             !  I/O manager
25
26   IMPLICIT NONE
27   PRIVATE
28
29   PUBLIC   p5z_prod         ! called in p5zbio.F90
30   PUBLIC   p5z_prod_init    ! called in trcsms_pisces.F90
31   PUBLIC   p5z_prod_alloc
32
33   !! * Shared module variables
34   REAL(wp), PUBLIC ::  pislope         !:
35   REAL(wp), PUBLIC ::  pislopep        !:
36   REAL(wp), PUBLIC ::  pislope2        !:
37   REAL(wp), PUBLIC ::  xadap           !:
38   REAL(wp), PUBLIC ::  excret          !:
39   REAL(wp), PUBLIC ::  excretp         !:
40   REAL(wp), PUBLIC ::  excret2         !:
41   REAL(wp), PUBLIC ::  bresp           !:
42   REAL(wp), PUBLIC ::  chlcpm          !:
43   REAL(wp), PUBLIC ::  chlcnm          !:
44   REAL(wp), PUBLIC ::  chlcdm          !:
45   REAL(wp), PUBLIC ::  chlcmin         !:
46   REAL(wp), PUBLIC ::  grosip          !:
47
48   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmaxn   !: optimal production = f(temperature)
49   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmaxp   !: optimal production = f(temperature)
50   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmaxd   !: optimal production = f(temperature)
51   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   zdaylen
52   
53   REAL(wp) :: r1_rday                !: 1 / rday
54   REAL(wp) :: texcret                !: 1 - excret
55   REAL(wp) :: texcretp               !: 1 - excretp
56   REAL(wp) :: texcret2               !: 1 - excret2       
57
58
59   !!* Substitution
60#  include "top_substitute.h90"
61   !!----------------------------------------------------------------------
62   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
63   !! $Id: p4zprod.F90 3160 2011-11-20 14:27:18Z cetlod $
64   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
65   !!----------------------------------------------------------------------
66CONTAINS
67
68   SUBROUTINE p5z_prod( kt , jnt )
69      !!---------------------------------------------------------------------
70      !!                     ***  ROUTINE p5z_prod  ***
71      !!
72      !! ** Purpose :   Compute the phytoplankton production depending on
73      !!              light, temperature and nutrient availability
74      !!
75      !! ** Method  : - ???
76      !!---------------------------------------------------------------------
77      !
78      INTEGER, INTENT(in) :: kt, jnt
79      !
80      INTEGER  ::   ji, jj, jk
81      REAL(wp) ::   zsilfac, znanotot, zpicotot, zdiattot, zconctemp, zconctemp2
82      REAL(wp) ::   zration, zratiop, zratiof, zmax, zmax2, zsilim, ztn, zadap
83      REAL(wp) ::   zpronmax, zpropmax, zprofmax, zrat
84      REAL(wp) ::   zlim, zsilfac2, zsiborn, zprod, zprontot, zproptot, zprodtot
85      REAL(wp) ::   zmxltst, zmxlday, zprnutmax
86      REAL(wp) ::   zpislopen, zpislopep, zpislope2n, zpsichl
87      REAL(wp) ::   zrum, zcodel, zargu, zval, thetannm, thetanpm, thetandm
88      REAL(wp) ::   zrfact2
89      CHARACTER (len=25) :: charout
90      REAL(wp), POINTER, DIMENSION(:,:  ) :: zmixnano, zmixpico, zmixdiat, zstrn
91      REAL(wp), POINTER, DIMENSION(:,:,:) :: zpislopead, zpislopeadp, zpislopead2
92      REAL(wp), POINTER, DIMENSION(:,:,:) :: zprbio, zprpic, zprdia, zysopt
93      REAL(wp), POINTER, DIMENSION(:,:,:) :: zprchln, zprchlp, zprchld
94      REAL(wp), POINTER, DIMENSION(:,:,:) :: zprorcan, zprorcap, zprorcad 
95      REAL(wp), POINTER, DIMENSION(:,:,:) :: zprofed, zprofep, zprofen
96      REAL(wp), POINTER, DIMENSION(:,:,:) :: zprochln, zprochlp, zprochld
97      REAL(wp), POINTER, DIMENSION(:,:,:) :: zpronew, zpronewp, zpronewd
98      REAL(wp), POINTER, DIMENSION(:,:,:) :: zproreg, zproregp, zproregd
99      REAL(wp), POINTER, DIMENSION(:,:,:) :: zpropo4, zpropo4p, zpropo4d
100      REAL(wp), POINTER, DIMENSION(:,:,:) :: zprodopn, zprodopp, zprodopd
101      REAL(wp), POINTER, DIMENSION(:,:,:) :: zrespn, zrespp, zrespd, zprnut
102      REAL(wp), POINTER, DIMENSION(:,:,:) :: zcroissn, zcroissp, zcroissd
103      !!---------------------------------------------------------------------
104      !
105      IF( nn_timing == 1 )  CALL timing_start('p5z_prod')
106      !
107      !  Allocate temporary workspace
108      CALL wrk_alloc( jpi, jpj,      zmixnano, zmixpico, zmixdiat, zstrn )
109      CALL wrk_alloc( jpi, jpj, jpk, zpislopead, zpislopeadp, zpislopead2, zysopt ) 
110      CALL wrk_alloc( jpi, jpj, jpk, zprdia, zprpic, zprbio, zprorcan, zprorcap, zprorcad )
111      CALL wrk_alloc( jpi, jpj, jpk, zprofed, zprofep, zprofen, zprochln, zprochlp, zprochld )
112      CALL wrk_alloc( jpi, jpj, jpk, zpronew, zpronewp, zpronewd, zproreg, zproregp, zproregd )
113      CALL wrk_alloc( jpi, jpj, jpk, zpropo4, zpropo4p, zpropo4d, zrespn, zrespp, zrespd, zprnut )
114      CALL wrk_alloc( jpi, jpj, jpk, zprchln, zprchlp, zprchld, zprodopn, zprodopp, zprodopd )
115      CALL wrk_alloc( jpi, jpj, jpk, zcroissp, zcroissn, zcroissd )
116      !
117      zprorcan(:,:,:) = 0._wp
118      zprorcap(:,:,:) = 0._wp
119      zprorcad(:,:,:) = 0._wp
120      zprofed (:,:,:) = 0._wp
121      zprofep (:,:,:) = 0._wp
122      zprofen (:,:,:) = 0._wp
123      zprochln(:,:,:) = 0._wp
124      zprochlp(:,:,:) = 0._wp
125      zprochld(:,:,:) = 0._wp
126      zpronew (:,:,:) = 0._wp
127      zpronewp(:,:,:) = 0._wp
128      zpronewd(:,:,:) = 0._wp
129      zproreg (:,:,:) = 0._wp
130      zproregp(:,:,:) = 0._wp
131      zproregd(:,:,:) = 0._wp
132      zpropo4 (:,:,:) = 0._wp
133      zpropo4p(:,:,:) = 0._wp
134      zpropo4d(:,:,:) = 0._wp
135      zprdia  (:,:,:) = 0._wp
136      zprpic  (:,:,:) = 0._wp
137      zprbio  (:,:,:) = 0._wp
138      zysopt  (:,:,:) = 0._wp
139      zrespn   (:,:,:) = 0._wp
140      zrespp   (:,:,:) = 0._wp
141      zrespd   (:,:,:) = 0._wp
142
143      ! Computation of the optimal production
144      zpsichl = 0.
145      thetannm = 0.25
146      thetanpm = 0.25
147      thetandm = 0.35
148
149      prmaxn(:,:,:) = ( 0.8_wp * (1. + zpsino3 * qnpmax + zpsichl ) ) * r1_rday * tgfunc(:,:,:)
150      prmaxp(:,:,:) = 0.6 / 0.8 * prmaxn(:,:,:) 
151      prmaxd(:,:,:) = 0.8 / 0.8 * prmaxn(:,:,:) 
152      zprnut(:,:,:) = 0.8_wp * r1_rday * tgfunc(:,:,:)
153
154      IF( lk_degrad )  THEN
155         prmaxn(:,:,:) = prmaxn(:,:,:) * facvol(:,:,:) 
156         prmaxp(:,:,:) = prmaxp(:,:,:) * facvol(:,:,:)
157         prmaxd(:,:,:) = prmaxd(:,:,:) * facvol(:,:,:)
158      ENDIF
159
160      ! compute the day length depending on latitude and the day
161      zrum = REAL( nday_year - 80, wp ) / REAL( nyear_len(1), wp )
162      zcodel = ASIN(  SIN( zrum * rpi * 2._wp ) * SIN( rad * 23.5_wp )  )
163
164      ! day length in hours
165      zstrn(:,:) = 0.
166      DO jj = 1, jpj
167         DO ji = 1, jpi
168            zargu = TAN( zcodel ) * TAN( gphit(ji,jj) * rad )
169            zargu = MAX( -1., MIN(  1., zargu ) )
170            zstrn(ji,jj) = MAX( 0.0, 24. - 2. * ACOS( zargu ) / rad / 15. )
171         END DO
172      END DO
173
174         ! Impact of the day duration on phytoplankton growth
175      DO jk = 1, jpkm1
176         DO jj = 1 ,jpj
177            DO ji = 1, jpi
178               IF( etot(ji,jj,jk) > 1.E-3 ) THEN
179                  zval = MAX( 1., zstrn(ji,jj) )
180                  zval = 1.5 * zval / (12. + zval)
181                  zprbio(ji,jj,jk) = prmaxn(ji,jj,jk) * zval
182                  zprpic(ji,jj,jk) = prmaxp(ji,jj,jk) * zval
183                  zprdia(ji,jj,jk) = prmaxd(ji,jj,jk) * zval
184               ENDIF
185            END DO
186         END DO
187      END DO
188
189      ! Maximum light intensity
190      zdaylen(:,:) = MAX(1., zstrn(:,:)) / 24.
191      WHERE( zstrn(:,:) < 1.e0 ) zstrn(:,:) = 24.
192      zstrn(:,:) = 24. / zstrn(:,:)
193
194      DO jk = 1, jpkm1
195!CDIR NOVERRCHK
196         DO jj = 1, jpj
197!CDIR NOVERRCHK
198            DO ji = 1, jpi
199                  ! Computation of the P-I slope for nanos and diatoms
200               IF( etot(ji,jj,jk) > 1.E-3 ) THEN
201                  ztn         = MAX( 0., tsn(ji,jj,jk,jp_tem) - 15. )
202                  zadap       = xadap * ztn / ( 2.+ ztn )
203                  znanotot    = enano(ji,jj,jk) * zstrn(ji,jj)
204                  zpicotot    = epico(ji,jj,jk) * zstrn(ji,jj)
205                  zdiattot    = ediat(ji,jj,jk) * zstrn(ji,jj)
206                  !
207                  zpislopead (ji,jj,jk) = pislope * trn(ji,jj,jk,jpnch)    &
208                  &                       /( trn(ji,jj,jk,jpphy) * 12. + rtrn)
209                  zpislopeadp(ji,jj,jk) = pislopep * ( 1. + zadap * EXP( -0.5 * zpicotot ) )   &
210                  &                       * trn(ji,jj,jk,jppch) /( trn(ji,jj,jk,jppic) * 12. + rtrn)
211                  zpislopead2(ji,jj,jk) = pislope2 * trn(ji,jj,jk,jpdch)    &
212                     &                    /( trn(ji,jj,jk,jpdia) * 12. + rtrn)
213                  zpislopen   = zpislopead (ji,jj,jk) / ( prmaxn(ji,jj,jk) * rday * xlimphy(ji,jj,jk) + rtrn )
214                  zpislopep   = zpislopeadp(ji,jj,jk) / ( prmaxp(ji,jj,jk) * rday * xlimpic(ji,jj,jk) + rtrn )
215                  zpislope2n  = zpislopead2(ji,jj,jk) / ( prmaxd(ji,jj,jk) * rday * xlimdia(ji,jj,jk) + rtrn )
216
217                  ! Computation of production function for Carbon
218                  !  ---------------------------------------------
219                  zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen  * znanotot )  )
220                  zprpic(ji,jj,jk) = zprpic(ji,jj,jk) * ( 1.- EXP( -zpislopep  * zpicotot )  )
221                  zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpislope2n * zdiattot )  )
222
223                  ! Computation of production function for Chlorophyll
224                  !  -------------------------------------------------
225                  zprchln(ji,jj,jk) = prmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen  * enano(ji,jj,jk) )  )
226                  zprchlp(ji,jj,jk) = prmaxp(ji,jj,jk) * ( 1.- EXP( -zpislopep  * epico(ji,jj,jk) )  )
227                  zprchld(ji,jj,jk) = prmaxd(ji,jj,jk) * ( 1.- EXP( -zpislope2n * ediat(ji,jj,jk) )  )
228               ENDIF
229            END DO
230         END DO
231      END DO
232
233      DO jk = 1, jpkm1
234         DO jj = 1, jpj
235            DO ji = 1, jpi
236              IF( etot(ji,jj,jk) > 1.E-3 ) THEN
237                  !    Si/C of diatoms
238                  !    ------------------------
239                  !    Si/C increases with iron stress and silicate availability
240                  !    Si/C is arbitrariliy increased for very high Si concentrations
241                  !    to mimic the very high ratios observed in the Southern Ocean (silpot2)
242                  zlim  = trn(ji,jj,jk,jpsil) / ( trn(ji,jj,jk,jpsil) + xksi1 )
243                  zsilim = MIN( zprdia(ji,jj,jk) / ( prmaxd(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) )
244                  zsilfac = 4.4 * EXP( -4.23 * zsilim ) * MAX( 0.e0, MIN( 1., 2.2 * ( zlim - 0.5 ) )  ) + 1.e0
245                  zsiborn = trn(ji,jj,jk,jpsil) * trn(ji,jj,jk,jpsil) * trn(ji,jj,jk,jpsil)
246                  IF (gphit(ji,jj) < -30 ) THEN
247                    zsilfac2 = 1. + 2. * zsiborn / ( zsiborn + xksi2**3 )
248                  ELSE
249                    zsilfac2 = 1. +      zsiborn / ( zsiborn + xksi2**3 )
250                  ENDIF
251                  zysopt(ji,jj,jk) = grosip * zlim * zsilfac * zsilfac2
252              ENDIF
253            END DO
254         END DO
255      END DO
256
257      !  Computation of the limitation term due to a mixed layer deeper than the euphotic depth
258      DO jj = 1, jpj
259         DO ji = 1, jpi
260            zmxltst = MAX( 0.e0, hmld(ji,jj) - heup(ji,jj) )
261            zmxlday = zmxltst * zmxltst * r1_rday
262            zmixnano(ji,jj) = 1. - zmxlday / ( 2. + zmxlday )
263            zmixpico(ji,jj) = 1. - zmxlday / ( 2. + zmxlday )
264            zmixdiat(ji,jj) = 1. - zmxlday / ( 4. + zmxlday )
265         END DO
266      END DO
267 
268      !  Mixed-layer effect on production                                                                               
269      DO jk = 1, jpkm1
270         DO jj = 1, jpj
271            DO ji = 1, jpi
272               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
273                  zprbio(ji,jj,jk)  = zprbio(ji,jj,jk) * zmixnano(ji,jj)
274                  zprpic(ji,jj,jk)  = zprpic(ji,jj,jk) * zmixpico(ji,jj)
275                  zprdia(ji,jj,jk)  = zprdia(ji,jj,jk) * zmixdiat(ji,jj)
276                  zprchln(ji,jj,jk) = zprchln(ji,jj,jk) * zmixnano(ji,jj)
277                  zprchlp(ji,jj,jk) = zprchlp(ji,jj,jk) * zmixpico(ji,jj)
278                  zprchld(ji,jj,jk) = zprchld(ji,jj,jk) * zmixdiat(ji,jj)
279               ENDIF
280            END DO
281         END DO
282      END DO
283
284      ! Computation of the various production terms of nanophytoplankton
285      DO jk = 1, jpkm1
286         DO jj = 1, jpj
287            DO ji = 1, jpi
288               IF( etot(ji,jj,jk) > 1.E-3 ) THEN
289                  !  production terms for nanophyto.
290                  zprorcan(ji,jj,jk) = zprbio(ji,jj,jk)  * xlimphy(ji,jj,jk) * trn(ji,jj,jk,jpphy) * rfact2
291                  !
292                  zrespn(ji,jj,jk) = zprorcan(ji,jj,jk) * zpsichl
293                  zration = trn(ji,jj,jk,jpnph) / ( trn(ji,jj,jk,jpphy) + rtrn )
294                  zratiop = trn(ji,jj,jk,jppph) / ( trn(ji,jj,jk,jpphy) + rtrn )
295                  zratiof = trn(ji,jj,jk,jpnfe) / ( trn(ji,jj,jk,jpphy) + rtrn )
296                  zprnutmax = zprnut(ji,jj,jk) * fvnuptk(ji,jj,jk) / rno3 * trn(ji,jj,jk,jpphy) * rfact2
297                  ! Uptake of nitrogen
298                  zrat = MIN( 1., zration / (xqnnmax(ji,jj,jk) + rtrn) ) 
299                  zmax = MAX(0., MIN(1., (1. - zrat)/ (1.05 - zrat) * 1.05))
300                  zpronmax = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqpnmin(ji,jj,jk) )   &
301                  &          / ( xqpnmax(ji,jj,jk) - xqpnmin(ji,jj,jk) ), xlimnfe(ji,jj,jk) ) )
302                  zpronew(ji,jj,jk) = zpronmax * zdaylen(ji,jj) * xnanono3(ji,jj,jk)
303                  zproreg(ji,jj,jk) = zpronmax * xnanonh4(ji,jj,jk)
304                  ! Uptake of phosphorus
305                  zrat = MIN( 1., zratiop / (xqpnmax(ji,jj,jk) + rtrn) )
306                  zmax = MAX(0., MIN(1., (1. - zrat)/ (1.05 - zrat) * 1.05))
307                  zpropmax = zprnutmax * zmax * xlimnfe(ji,jj,jk)
308                  zpropo4(ji,jj,jk) = zpropmax * xnanopo4(ji,jj,jk)
309                  zprodopn(ji,jj,jk) = zpropmax * xnanodop(ji,jj,jk)
310                  ! Uptake of iron
311                  zrat = MIN( 1., zratiof / qfnmax )
312                  zmax = MAX(0., MIN(1., (1. - zrat)/ (1.05 - zrat) * 1.05))
313                  zprofmax = zprnutmax * qfnmax * zmax
314                  zprofen(ji,jj,jk) = zprofmax * xnanofer(ji,jj,jk) * ( 4. - 4.5 * xlimnfe(ji,jj,jk)    &
315                  &          / ( xlimnfe(ji,jj,jk) + 0.2 ) ) * (1. + 0.8 * xnanono3(ji,jj,jk) / ( rtrn  &
316                  &          + xnanono3(ji,jj,jk) + xnanonh4(ji,jj,jk) ) * (1. - xnanofer(ji,jj,jk) ) )
317               ENDIF
318            END DO
319         END DO
320      END DO
321
322      ! Computation of the various production terms of picophytoplankton
323      DO jk = 1, jpkm1
324         DO jj = 1, jpj
325            DO ji = 1, jpi
326               IF( etot(ji,jj,jk) > 1.E-3 ) THEN
327                  !  production terms for picophyto.
328                  zprorcap(ji,jj,jk) = zprpic(ji,jj,jk)  * xlimpic(ji,jj,jk) * trn(ji,jj,jk,jppic) * rfact2
329                  !
330                  zrespp(ji,jj,jk) = zprorcap(ji,jj,jk) * zpsichl
331                  zration = trn(ji,jj,jk,jpnpi) / ( trn(ji,jj,jk,jppic) + rtrn )
332                  zratiop = trn(ji,jj,jk,jpppi) / ( trn(ji,jj,jk,jppic) + rtrn )
333                  zratiof = trn(ji,jj,jk,jppfe) / ( trn(ji,jj,jk,jppic) + rtrn )
334                  zprnutmax = zprnut(ji,jj,jk) * fvpuptk(ji,jj,jk) / rno3 * trn(ji,jj,jk,jppic) * rfact2
335                  ! Uptake of nitrogen
336                  zrat = MIN( 1., zration / (xqnpmax(ji,jj,jk) + rtrn) )
337                  zmax = MAX(0., MIN(1., (1. - zrat)/ (1.05 - zrat) * 1.05))
338                  zpronmax = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqppmin(ji,jj,jk) )   &
339                  &          / ( xqppmax(ji,jj,jk) - xqppmin(ji,jj,jk) ), xlimpfe(ji,jj,jk) ) )
340                  zpronewp(ji,jj,jk) = zpronmax * zdaylen(ji,jj) * xpicono3(ji,jj,jk) 
341                  zproregp(ji,jj,jk) = zpronmax * xpiconh4(ji,jj,jk)
342                  ! Uptake of phosphorus
343                  zrat = MIN( 1., zratiop / (xqppmax(ji,jj,jk) + rtrn) )
344                  zmax = MAX(0., MIN(1., (1. - zrat)/ (1.05 - zrat) * 1.05))
345                  zpropmax = zprnutmax * zmax * xlimpfe(ji,jj,jk)
346                  zpropo4p(ji,jj,jk) = zpropmax * xpicopo4(ji,jj,jk)
347                  zprodopp(ji,jj,jk) = zpropmax * xpicodop(ji,jj,jk)
348                  ! Uptake of iron
349                  zrat = MIN( 1., zratiof / qfpmax )
350                  zmax = MAX(0., MIN(1., (1. - zrat)/ (1.05 - zrat) * 1.05))
351                  zprofmax = zprnutmax * qfpmax * zmax
352                  zprofep(ji,jj,jk) = zprofmax * xpicofer(ji,jj,jk) * ( 4. - 4.5 * xlimpfe(ji,jj,jk)   &
353                  &          / ( xlimpfe(ji,jj,jk) + 0.2 ) ) * (1. + 0.8 * xpicono3(ji,jj,jk) / ( rtrn   &
354                  &          + xpicono3(ji,jj,jk) + xpiconh4(ji,jj,jk) ) * (1. - xpicofer(ji,jj,jk) ) )
355               ENDIF
356            END DO
357         END DO
358      END DO
359
360      ! Computation of the various production terms of diatoms
361      DO jk = 1, jpkm1
362         DO jj = 1, jpj
363            DO ji = 1, jpi
364               IF( etot(ji,jj,jk) > 1.E-3 ) THEN
365                  !  production terms for diatomees
366                  zprorcad(ji,jj,jk) = zprdia(ji,jj,jk) * xlimdia(ji,jj,jk) * trn(ji,jj,jk,jpdia) * rfact2
367                  ! Computation of the respiration term according to pahlow
368                  ! & oschlies (2013)
369                  zrespd(ji,jj,jk) = zprorcad(ji,jj,jk) * zpsichl
370                  !
371                  zration = trn(ji,jj,jk,jpndi) / ( trn(ji,jj,jk,jpdia) + rtrn )
372                  zratiop = trn(ji,jj,jk,jppdi) / ( trn(ji,jj,jk,jpdia) + rtrn )
373                  zratiof = trn(ji,jj,jk,jpdfe) / ( trn(ji,jj,jk,jpdia) + rtrn )
374                  zprnutmax = zprnut(ji,jj,jk) * fvduptk(ji,jj,jk) / rno3 * trn(ji,jj,jk,jpdia) * rfact2
375                  ! Uptake of nitrogen
376                  zrat = MIN( 1., zration / (xqndmax(ji,jj,jk) + rtrn) )
377                  zmax = MAX(0., MIN(1., (1. - zrat)/ (1.05 - zrat) * 1.05)) 
378                  zpronmax = zprnutmax * zmax * MAX(0., MIN(1., ( zratiop - xqpdmin(ji,jj,jk) )   &
379                  &          / ( xqpdmax(ji,jj,jk) - xqpdmin(ji,jj,jk) ), xlimdfe(ji,jj,jk) ) )
380                  zpronewd(ji,jj,jk) = zpronmax * zdaylen(ji,jj) * xdiatno3(ji,jj,jk)
381                  zproregd(ji,jj,jk) = zpronmax * xdiatnh4(ji,jj,jk)
382                  ! Uptake of phosphorus
383                  zrat = MIN( 1., zratiop / (xqpdmax(ji,jj,jk) + rtrn) )
384                  zmax = MAX(0., MIN(1., (1. - zrat)/ (1.05 - zrat) * 1.05)) 
385                  zpropmax = zprnutmax * zmax * xlimdfe(ji,jj,jk)
386                  zpropo4d(ji,jj,jk) = zpropmax * xdiatpo4(ji,jj,jk)
387                  zprodopd(ji,jj,jk) = zpropmax * xdiatdop(ji,jj,jk)
388                  ! Uptake of iron
389                  zrat = MIN( 1., zratiof / qfdmax )
390                  zmax = MAX(0., MIN(1., (1. - zrat)/ (1.05 - zrat) * 1.05))
391                  zprofmax = zprnutmax * qfdmax * zmax
392                  zprofed(ji,jj,jk) = zprofmax * xdiatfer(ji,jj,jk) * ( 4. - 4.5 * xlimdfe(ji,jj,jk)     &
393                  &          / ( xlimdfe(ji,jj,jk) + 0.2 ) ) * (1. + 0.8 * xdiatno3(ji,jj,jk) / ( rtrn   &
394                  &          + xdiatno3(ji,jj,jk) + xdiatnh4(ji,jj,jk) ) * (1. - xdiatfer(ji,jj,jk) ) )
395               ENDIF
396            END DO
397         END DO
398      END DO
399
400      DO jk = 1, jpkm1
401         DO jj = 1, jpj
402            DO ji = 1, jpi
403               IF( etot(ji,jj,jk) > 1.E-3 ) THEN
404                     !  production terms for nanophyto. ( chlorophyll )
405                  zprod = rday * (zpronew(ji,jj,jk) + zproreg(ji,jj,jk)) * zprchln(ji,jj,jk) * xlimphy(ji,jj,jk)
406                  zprochln(ji,jj,jk) = thetannm * zprod / ( zpislopead(ji,jj,jk) * enano(ji,jj,jk) +rtrn )
407                  zprochln(ji,jj,jk) = MAX(zprochln(ji,jj,jk), chlcmin * 12. * zprorcan (ji,jj,jk) )
408                     !  production terms for picophyto. ( chlorophyll )
409                  zprod = rday * (zpronewp(ji,jj,jk) + zproregp(ji,jj,jk)) * zprchlp(ji,jj,jk) * xlimpic(ji,jj,jk)
410                  zprochlp(ji,jj,jk) = thetanpm * zprod / ( zpislopeadp(ji,jj,jk) * epico(ji,jj,jk) +rtrn )
411                  zprochlp(ji,jj,jk) = MAX(zprochlp(ji,jj,jk), chlcmin * 12. * zprorcap(ji,jj,jk) )
412                  !  production terms for diatomees ( chlorophyll )
413                  zprod = rday * (zpronewd(ji,jj,jk) + zproregd(ji,jj,jk)) * zprchld(ji,jj,jk) * xlimdia(ji,jj,jk)
414                  zprochld(ji,jj,jk) = thetandm * zprod / ( zpislopead2(ji,jj,jk) * ediat(ji,jj,jk) +rtrn )
415                  zprochld(ji,jj,jk) = MAX(zprochld(ji,jj,jk), chlcmin * 12. * zprorcad(ji,jj,jk) )
416               ENDIF
417            END DO
418         END DO
419      END DO
420
421      !   Update the arrays TRA which contain the biological sources and sinks
422      DO jk = 1, jpkm1
423         DO jj = 1, jpj
424           DO ji =1 ,jpi
425              zprontot = zpronew(ji,jj,jk) + zproreg(ji,jj,jk)
426              zproptot = zpronewp(ji,jj,jk) + zproregp(ji,jj,jk)
427              zprodtot = zpronewd(ji,jj,jk) + zproregd(ji,jj,jk)
428              tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) - zpropo4(ji,jj,jk) - zpropo4d(ji,jj,jk)  &
429              &                     - zpropo4p(ji,jj,jk)
430              tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) - zpronew(ji,jj,jk) - zpronewd(ji,jj,jk)  &
431              &                     - zpronewp(ji,jj,jk)
432              tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) - zproreg(ji,jj,jk) - zproregd(ji,jj,jk)  &
433              &                     - zproregp(ji,jj,jk)
434              tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) + zprorcan(ji,jj,jk) * texcret     &
435                 &                  - zpsino3 * zpronew(ji,jj,jk) - zpsinh4 * zproreg(ji,jj,jk)   &
436                 &                  - zrespn(ji,jj,jk) 
437              zcroissn(ji,jj,jk) = tra(ji,jj,jk,jpphy) / rfact2/ (trn(ji,jj,jk,jpphy) + rtrn)
438              tra(ji,jj,jk,jpnph) = tra(ji,jj,jk,jpnph) + zprontot * texcret
439              tra(ji,jj,jk,jppph) = tra(ji,jj,jk,jppph) + zpropo4(ji,jj,jk) * texcret   &
440              &                     + zprodopn(ji,jj,jk) * texcret
441              tra(ji,jj,jk,jpnch) = tra(ji,jj,jk,jpnch) + zprochln(ji,jj,jk) * texcret
442              tra(ji,jj,jk,jpnfe) = tra(ji,jj,jk,jpnfe) + zprofen(ji,jj,jk) * texcret
443              tra(ji,jj,jk,jppic) = tra(ji,jj,jk,jppic) + zprorcap(ji,jj,jk) * texcretp     &
444                 &                  - zpsino3 * zpronewp(ji,jj,jk) - zpsinh4 * zproregp(ji,jj,jk)   &
445                 &                  - zrespp(ji,jj,jk) 
446              zcroissp(ji,jj,jk) = tra(ji,jj,jk,jppic) / rfact2/ (trn(ji,jj,jk,jppic) + rtrn)
447              tra(ji,jj,jk,jpnpi) = tra(ji,jj,jk,jpnpi) + zproptot * texcretp
448              tra(ji,jj,jk,jpppi) = tra(ji,jj,jk,jpppi) + zpropo4p(ji,jj,jk) * texcretp   &
449              &                     + zprodopp(ji,jj,jk) * texcretp
450              tra(ji,jj,jk,jppch) = tra(ji,jj,jk,jppch) + zprochlp(ji,jj,jk) * texcretp
451              tra(ji,jj,jk,jppfe) = tra(ji,jj,jk,jppfe) + zprofep(ji,jj,jk) * texcretp
452              tra(ji,jj,jk,jpdia) = tra(ji,jj,jk,jpdia) + zprorcad(ji,jj,jk) * texcret2   &
453                 &                  - zpsino3 * zpronewd(ji,jj,jk) - zpsinh4 * zproregd(ji,jj,jk)   &
454                 &                  - zrespd(ji,jj,jk) 
455              zcroissd(ji,jj,jk) = tra(ji,jj,jk,jpdia) / rfact2 / (trn(ji,jj,jk,jpdia) + rtrn)
456              tra(ji,jj,jk,jpndi) = tra(ji,jj,jk,jpndi) + zprodtot * texcret2
457              tra(ji,jj,jk,jppdi) = tra(ji,jj,jk,jppdi) + zpropo4d(ji,jj,jk) * texcret2   &
458              &                     + zprodopd(ji,jj,jk) * texcret2
459              tra(ji,jj,jk,jpdch) = tra(ji,jj,jk,jpdch) + zprochld(ji,jj,jk) * texcret2
460              tra(ji,jj,jk,jpdfe) = tra(ji,jj,jk,jpdfe) + zprofed(ji,jj,jk) * texcret2
461              tra(ji,jj,jk,jpdsi) = tra(ji,jj,jk,jpdsi) + zprorcad(ji,jj,jk) * zysopt(ji,jj,jk) * texcret2
462              tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + excret2 * zprorcad(ji,jj,jk) + excret * zprorcan(ji,jj,jk)  &
463              &                     + excretp * zprorcap(ji,jj,jk)
464              tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) + excret2 * zprodtot + excret * zprontot   &
465              &                     + excretp * zproptot
466              tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) + excret2 * zpropo4d(ji,jj,jk) + excret * zpropo4(ji,jj,jk)   &
467              &    - texcret * zprodopn(ji,jj,jk) - texcret2 * zprodopd(ji,jj,jk) + excretp * zpropo4p(ji,jj,jk)     &
468              &    - texcretp * zprodopp(ji,jj,jk)
469              tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + o2ut * ( zproreg(ji,jj,jk) + zproregd(ji,jj,jk)   &
470                 &                + zproregp(ji,jj,jk) ) + ( o2ut + o2nit ) * ( zpronew(ji,jj,jk)           &
471                 &                + zpronewd(ji,jj,jk) + zpronewp(ji,jj,jk) )   &
472                 &                - o2ut * ( zrespn(ji,jj,jk) + zrespp(ji,jj,jk) + zrespd(ji,jj,jk) )
473              tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - texcret * zprofen(ji,jj,jk) - texcret2 * zprofed(ji,jj,jk)  &
474              &                     - texcretp * zprofep(ji,jj,jk)
475              tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) - texcret2 * zprorcad(ji,jj,jk) * zysopt(ji,jj,jk)
476              tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) - zprorcan(ji,jj,jk) - zprorcad(ji,jj,jk) - zprorcap(ji,jj,jk)  &
477              &                     + zpsino3 * zpronew(ji,jj,jk) + zpsinh4 * zproreg(ji,jj,jk)   &
478              &                     + zpsino3 * zpronewp(ji,jj,jk) + zpsinh4 * zproregp(ji,jj,jk)   &
479              &                     + zpsino3 * zpronewd(ji,jj,jk) + zpsinh4 * zproregd(ji,jj,jk)  &
480              &                     + zrespn(ji,jj,jk) + zrespd(ji,jj,jk) + zrespp(ji,jj,jk) 
481              tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zpronew(ji,jj,jk) + zpronewd(ji,jj,jk)  &
482              &                     + zpronewp(ji,jj,jk) ) - rno3 * ( zproreg(ji,jj,jk) + zproregd(ji,jj,jk)     &
483              &                     + zproregp(ji,jj,jk) ) 
484          END DO
485        END DO
486     END DO
487
488     ! Total primary production per year
489     tpp = tpp + glob_sum( ( zprorcan(:,:,:) + zprorcad(:,:,:) + zprorcap(:,:,:) ) * cvol(:,:,:) )
490
491     IF( ln_diatrc ) THEN
492         !
493         zrfact2 = 1.e3 * rfact2r  ! conversion from mol/L/timestep into mol/m3/s
494         IF( lk_iomput ) THEN
495           IF( jnt == nrdttrc ) THEN
496              CALL iom_put( "PPPHY"  , zprorcan(:,:,:) * zrfact2 * tmask(:,:,:) )  ! primary production by nanophyto
497              CALL iom_put( "PPPHYP" , zprorcap(:,:,:) * zrfact2 * tmask(:,:,:) )  ! primary production by nanophyto
498              CALL iom_put( "PPPHY2" , zprorcad(:,:,:) * zrfact2 * tmask(:,:,:) )  ! primary production by diatom
499              CALL iom_put( "PPNEWN" , zpronew (:,:,:) * zrfact2 * tmask(:,:,:) )  ! new primary production by nanophyto
500              CALL iom_put( "PPNEWP" , zpronewp(:,:,:) * zrfact2 * tmask(:,:,:) )  ! new primary production by nanophyto
501              CALL iom_put( "PPNEWD" , zpronewd(:,:,:) * zrfact2 * tmask(:,:,:) )  ! new primary production by diatom
502              CALL iom_put( "PBSi"   , zprorcad(:,:,:) * zrfact2 * tmask(:,:,:) * zysopt(:,:,:) ) ! biogenic silica production
503              CALL iom_put( "PFeD"   , zprofed (:,:,:) * zrfact2 * tmask(:,:,:) )  ! biogenic iron production by diatom
504              CALL iom_put( "PFeP"   , zprofep (:,:,:) * zrfact2 * tmask(:,:,:) )  ! biogenic iron production by diatom
505              CALL iom_put( "PFeN"   , zprofen (:,:,:) * zrfact2 * tmask(:,:,:) )  ! biogenic iron production by nanophyto
506              CALL iom_put( "Mumax"  , prmaxn(:,:,:) * tmask(:,:,:) )  ! Maximum growth rate
507              CALL iom_put( "MuN"    , zprbio(:,:,:) * xlimphy(:,:,:) * tmask(:,:,:) )  ! Realized growth rate for nanophyto
508              CALL iom_put( "MuP"    , zprpic(:,:,:) * xlimpic(:,:,:) * tmask(:,:,:) )  ! Realized growth rate for nanophyto
509              CALL iom_put( "MuD"    , zprdia(:,:,:) * xlimdia(:,:,:) * tmask(:,:,:) )  ! Realized growth rate for diatoms
510              CALL iom_put( "MunetN"    , zcroissn(:,:,:) * tmask(:,:,:) )  ! Realized growth rate for nanophyto
511              CALL iom_put( "MunetP"    , zcroissp(:,:,:) * tmask(:,:,:) )  ! Realized growth rate for nanophyto
512              CALL iom_put( "MunetD"    , zcroissd(:,:,:) * tmask(:,:,:) )  ! Realized growth rate for diatoms
513              CALL iom_put( "LNlight", zprbio (:,:,:) / (prmaxn(:,:,:) + rtrn) * tmask(:,:,:) )  ! light limitation term
514              CALL iom_put( "LPlight", zprpic (:,:,:) / (prmaxp(:,:,:) + rtrn) * tmask(:,:,:) )  ! light limitation term
515              CALL iom_put( "LDlight", zprdia (:,:,:) / (prmaxd(:,:,:) + rtrn) * tmask(:,:,:) )  ! light limitation term
516           ENDIF
517         ELSE
518              trc3d(:,:,:,jp_pcs0_3d + 4)  = zprorcan(:,:,:) * zrfact2 * tmask(:,:,:)
519              trc3d(:,:,:,jp_pcs0_3d + 5)  = zprorcad(:,:,:) * zrfact2 * tmask(:,:,:)
520              trc3d(:,:,:,jp_pcs0_3d + 6)  = zpronew (:,:,:) * zrfact2 * tmask(:,:,:)
521              trc3d(:,:,:,jp_pcs0_3d + 7)  = zpronewd(:,:,:) * zrfact2 * tmask(:,:,:)
522              trc3d(:,:,:,jp_pcs0_3d + 8)  = zprorcad(:,:,:) * zrfact2 * tmask(:,:,:) * zysopt(:,:,:)
523              trc3d(:,:,:,jp_pcs0_3d + 9)  = zprofed (:,:,:) * zrfact2 * tmask(:,:,:)
524#  if ! defined key_kriest
525              trc3d(:,:,:,jp_pcs0_3d + 10) = zprofen (:,:,:) * zrfact2 * tmask(:,:,:)
526#  endif
527         ENDIF
528         !
529      ENDIF
530
531      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
532         WRITE(charout, FMT="('prod')")
533         CALL prt_ctl_trc_info(charout)
534         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
535      ENDIF
536      !
537      CALL wrk_dealloc( jpi, jpj,      zmixnano, zmixpico, zmixdiat, zstrn )
538      CALL wrk_dealloc( jpi, jpj, jpk, zpislopead, zpislopeadp, zpislopead2, zysopt )                           
539      CALL wrk_dealloc( jpi, jpj, jpk, zprdia, zprpic, zprbio, zprorcan, zprorcap, zprorcad )
540      CALL wrk_dealloc( jpi, jpj, jpk, zprofed, zprofep, zprofen, zprochln, zprochlp, zprochld ) 
541      CALL wrk_dealloc( jpi, jpj, jpk, zpronew, zpronewp, zpronewd, zproreg, zproregp, zproregd )
542      CALL wrk_dealloc( jpi, jpj, jpk, zpropo4, zpropo4p, zpropo4d, zrespn, zrespp, zrespd, zprnut )
543      CALL wrk_dealloc( jpi, jpj, jpk, zprchln, zprchlp, zprchld, zprodopn, zprodopp, zprodopd )
544      CALL wrk_dealloc( jpi, jpj, jpk, zcroissp, zcroissn, zcroissd )
545      !
546      IF( nn_timing == 1 )  CALL timing_stop('p5z_prod')
547      !
548   END SUBROUTINE p5z_prod
549
550
551   SUBROUTINE p5z_prod_init
552      !!----------------------------------------------------------------------
553      !!                  ***  ROUTINE p5z_prod_init  ***
554      !!
555      !! ** Purpose :   Initialization of phytoplankton production parameters
556      !!
557      !! ** Method  :   Read the nampisprod namelist and check the parameters
558      !!      called at the first timestep (nittrc000)
559      !!
560      !! ** input   :   Namelist nampisprod
561      !!----------------------------------------------------------------------
562      !
563      NAMELIST/nampisprod/ pislope, pislopep, pislope2, xadap, bresp, excret, excretp, excret2,  &
564         &                 chlcnm, chlcpm, chlcdm, chlcmin, grosip
565
566      INTEGER :: ios                 ! Local integer output status for namelist read
567      !!----------------------------------------------------------------------
568
569      REWIND( numnatp_ref )              ! Namelist nampisprod in reference namelist : Pisces phytoplankton production
570      READ  ( numnatp_ref, nampisprod, IOSTAT = ios, ERR = 901)
571901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisprod in reference namelist', lwp )
572
573      REWIND( numnatp_cfg )              ! Namelist nampisprod in configuration namelist : Pisces phytoplankton production
574      READ  ( numnatp_cfg, nampisprod, IOSTAT = ios, ERR = 902 )
575902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisprod in configuration namelist', lwp )
576      IF(lwm) WRITE ( numonp, nampisprod )
577
578      IF(lwp) THEN                         ! control print
579         WRITE(numout,*) ' '
580         WRITE(numout,*) ' Namelist parameters for phytoplankton growth, nampisprod'
581         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
582         WRITE(numout,*) '    mean Si/C ratio                           grosip       =', grosip
583         WRITE(numout,*) '    P-I slope                                 pislope      =', pislope
584         WRITE(numout,*) '    Acclimation factor to low light           xadap       =', xadap
585         WRITE(numout,*) '    excretion ratio of nanophytoplankton      excret       =', excret
586         WRITE(numout,*) '    excretion ratio of picophytoplankton      excretp      =', excretp
587         WRITE(numout,*) '    excretion ratio of diatoms                excret2      =', excret2
588         WRITE(numout,*) '    basal respiration in phytoplankton        bresp        =', bresp
589         WRITE(numout,*) '    Maximum Chl/C in phytoplankton            chlcmin      =', chlcmin
590         WRITE(numout,*) '    P-I slope  for diatoms                    pislope2     =', pislope2
591         WRITE(numout,*) '    P-I slope  for picophytoplankton          pislopep     =', pislopep
592         WRITE(numout,*) '    Minimum Chl/C in nanophytoplankton        chlcnm       =', chlcnm
593         WRITE(numout,*) '    Minimum Chl/C in picophytoplankton        chlcpm       =', chlcpm
594         WRITE(numout,*) '    Minimum Chl/C in diatoms                  chlcdm       =', chlcdm
595      ENDIF
596      !
597      r1_rday   = 1._wp / rday 
598      texcret   = 1._wp - excret
599      texcretp  = 1._wp - excretp
600      texcret2  = 1._wp - excret2
601      tpp       = 0._wp
602      !
603   END SUBROUTINE p5z_prod_init
604
605
606   INTEGER FUNCTION p5z_prod_alloc()
607      !!----------------------------------------------------------------------
608      !!                     ***  ROUTINE p5z_prod_alloc  ***
609      !!----------------------------------------------------------------------
610      ALLOCATE( prmaxn(jpi,jpj,jpk), prmaxp(jpi,jpj,jpk), prmaxd(jpi,jpj,jpk),  &
611      &          zdaylen(jpi,jpj), STAT = p5z_prod_alloc )
612      !
613      IF( p5z_prod_alloc /= 0 ) CALL ctl_warn('p5z_prod_alloc : failed to allocate arrays.')
614      !
615   END FUNCTION p5z_prod_alloc
616
617#else
618   !!======================================================================
619   !!  Dummy module :                                   No PISCES bio-model
620   !!======================================================================
621CONTAINS
622   SUBROUTINE p5z_prod                    ! Empty routine
623   END SUBROUTINE p5z_prod
624#endif 
625
626   !!======================================================================
627END MODULE  p5zprod
Note: See TracBrowser for help on using the repository browser.