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

source: NEMO/branches/2019/dev_r12072_MERGE_OPTION2_2019/src/TOP/PISCES/P4Z/p5zprod.F90 @ 12282

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

bugfix in dev_r12072_MERGE_OPTION2_2019 : initialise local array in PISCES, see ticket #2357

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