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.
p4zpoc.F90 in branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z – NEMO

source: branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zpoc.F90 @ 6453

Last change on this file since 6453 was 6453, checked in by aumont, 8 years ago

New developments of PISCES (QUOTA, ligands, lability, ...)

File size: 23.7 KB
Line 
1MODULE p4zpoc
2   !!======================================================================
3   !!                         ***  MODULE p4zpoc  ***
4   !! TOP :   PISCES Compute remineralization of organic particles
5   !!=========================================================================
6   !! History :   1.0  !  2004     (O. Aumont) Original code
7   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
8   !!             3.4  !  2011-06  (O. Aumont, C. Ethe) Quota model for iron
9   !!             3.6  !  2016-03  (O. Aumont) Quota model and diverse
10   !!----------------------------------------------------------------------
11#if defined key_pisces
12   !!----------------------------------------------------------------------
13   !!   'key_top'       and                                      TOP models
14   !!   'key_pisces'                                       PISCES bio-model
15   !!----------------------------------------------------------------------
16   !!   p4z_poc       :  Compute remineralization/dissolution of organic compounds
17   !!   p4z_poc_init  :  Initialisation of parameters for remineralisation
18   !!----------------------------------------------------------------------
19   USE oce_trc         !  shared variables between ocean and passive tracers
20   USE trc             !  passive tracers common variables
21   USE sms_pisces      !  PISCES Source Minus Sink variables
22   USE p4zsink         !  Sedimentation of organic particles
23   USE prtctl_trc      !  print control for debugging
24   USE iom             !  I/O manager
25
26
27   IMPLICIT NONE
28   PRIVATE
29
30   PUBLIC   p4z_poc         ! called in p4zbio.F90
31   PUBLIC   p4z_poc_init    ! called in trcsms_pisces.F90
32
33   !! * Shared module variables
34   REAL(wp), PUBLIC ::  xremip     !: remineralisation rate of POC
35   INTEGER , PUBLIC ::  jcpoc      !: number of lability classes
36
37
38   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:)     ::   alphan, reminp
39   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: alphap
40
41
42   !!* Substitution
43#  include "top_substitute.h90"
44   !!----------------------------------------------------------------------
45   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
46   !! $Id: p4zrem.F90 3160 2011-11-20 14:27:18Z cetlod $
47   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
48   !!----------------------------------------------------------------------
49CONTAINS
50
51   SUBROUTINE p4z_poc( kt, jnt )
52      !!---------------------------------------------------------------------
53      !!                     ***  ROUTINE p4z_poc  ***
54      !!
55      !! ** Purpose :   Compute remineralization of organic particles
56      !!
57      !! ** Method  : - ???
58      !!---------------------------------------------------------------------
59      !
60      INTEGER, INTENT(in) ::   kt, jnt ! ocean time step
61      !
62      INTEGER  ::   ji, jj, jk, jn
63      REAL(wp) ::   zremip, zremig, zdep, zorem2, zofer
64      REAL(wp) ::   ztem, zsizek, zsizek1, alphat, remint
65      REAL(wp) ::   solgoc, zpoc, zpoctot, zremif
66#if ! defined key_kriest
67      REAL(wp) ::   zofer2, zofer3
68#endif
69      REAL(wp) ::   zstep, zrfact2
70      CHARACTER (len=25) :: charout
71      REAL(wp), POINTER, DIMENSION(:,:  ) :: totprod, totthick, totcons 
72      REAL(wp), POINTER, DIMENSION(:,:,:) :: zremipoc, zremigoc, zorem3
73      REAL(wp), POINTER, DIMENSION(:,:,:,:) :: alphag
74      !!---------------------------------------------------------------------
75      !
76      IF( nn_timing == 1 )  CALL timing_start('p4z_rem')
77      !
78      ! Allocate temporary workspace
79      CALL wrk_alloc( jpi, jpj,      totprod, totthick, totcons )
80      CALL wrk_alloc( jpi, jpj, jpk, zremipoc, zremigoc, zorem3 )
81      ALLOCATE( alphag(jpi,jpj,jpk,jcpoc) )
82
83      ! Initialization of local variables
84      ! ---------------------------------
85
86      ! Here we compute the GOC -> POC rate due to the shrinking
87      ! of the fecal pellets/aggregates as a result of bacterial
88      ! solubilization
89      ! This is based on a fractal dimension of 2.56 and a spectral
90      ! slope of -3.6 (identical to what is used in p4zsink to compute
91      ! aggregation
92      solgoc = 0.04/ 2.56 * 1./ ( 1.-50**(-0.04) )
93
94      ! Initialisation of temprary arrys
95      zremipoc(:,:,:) = xremip
96      zremigoc(:,:,:) = xremip
97      zorem3(:,:,:)   = 0.
98      orem  (:,:,:)   = 0.
99
100      DO jn = 1, jcpoc
101        alphag(:,:,:,jn) = alphan(jn)
102      END DO
103
104#if ! defined key_kriest
105! -----------------------------------------------------------------------
106! Lability parameterization. This is the big particles part (GOC)
107! This lability parameterization can be activated only with the standard
108! particle scheme. Does not work with Kriest parameterization.
109! -----------------------------------------------------------------------
110     DO jk = 2, jpkm1
111        DO jj = 1, jpj
112           DO ji = 1, jpi
113              IF (tmask(ji,jj,jk) == 1.) THEN
114                zdep = hmld(ji,jj)
115                !
116                ! In the case of GOC, lability is constant in the mixed layer
117                ! It is computed only below the mixed layer depth
118                ! ------------------------------------------------------------
119                !
120                IF( fsdept(ji,jj,jk) >= zdep ) THEN
121                  alphat = 0.
122                  remint = 0.
123                  !
124                  ! The first level just below the mixed layer needs a
125                  ! specific treatment because lability is supposed constant
126                  ! everywhere within the mixed layer. This means that
127                  ! change in lability in the bottom part of the previous cell
128                  ! should not be computed
129                  ! ----------------------------------------------------------
130                  !
131                  IF ( fsdept(ji,jj,jk-1) < zdep ) THEN
132                    DO jn = 1, jcpoc
133                       !
134                       ! Lagrangian based algorithm. The fraction of each
135                       ! lability class is computed starting from the previous
136                       ! level
137                       ! -----------------------------------------------------
138                       !
139                       zsizek  = zdep / (wsbio2 + rtrn) * tgfunc(ji,jj,jk-1)
140                       zsizek1 = fse3t(ji,jj,jk) / 2. / (wsbio4(ji,jj,jk) + rtrn) * tgfunc(ji,jj,jk)
141                       ! POC concentration is computed using the lagrangian
142                       ! framework. It is only used for the lability param
143                       zpoc = trb(ji,jj,jk-1,jpgoc) + consgoc(ji,jj,jk) * rday / rfact2               &
144                       &   * fse3t(ji,jj,jk) / 2. / (wsbio4(ji,jj,jk) + rtrn)
145                       zpoc = max(0., zpoc)
146                       ! the concentration of each lability class is calculated
147                       ! as the sum of the different sources and sinks
148                       ! Please note that production of new GOC experiences
149                       ! degradation
150                       alphag(ji,jj,jk,jn) = alphan(jn) / (reminp(jn) * tgfunc(ji,jj,jk-1) )       &
151                       &   * (1. - exp( -reminp(jn) * zsizek ) ) * exp( -reminp(jn) * zsizek1 )     &
152                       &   * zpoc + prodgoc(ji,jj,jk) * ( 1. - exp( -reminp(jn) * zsizek1 ) )        &
153                       &   * rday / rfact2 * alphan(jn) / reminp(jn) / tgfunc(ji,jj,jk)
154                       alphat = alphat + alphag(ji,jj,jk,jn)
155                       remint = remint + alphag(ji,jj,jk,jn) * reminp(jn)
156                    END DO
157                  ELSE
158                    !
159                    ! standard algorithm in the rest of the water column
160                    ! See the comments in the previous block.
161                    ! ---------------------------------------------------
162                    !
163                    DO jn = 1, jcpoc
164                       zsizek  = fse3t(ji,jj,jk-1) / 2. / (wsbio4(ji,jj,jk-1) + rtrn) * tgfunc(ji,jj,jk-1)
165                       zsizek1 = fse3t(ji,jj,jk) / 2. / (wsbio4(ji,jj,jk) + rtrn) * tgfunc(ji,jj,jk)
166                       zpoc = trb(ji,jj,jk-1,jpgoc) + consgoc(ji,jj,jk-1) * rday / rfact2               &
167                       &   * fse3t(ji,jj,jk-1) / 2. / (wsbio4(ji,jj,jk-1) + rtrn) + consgoc(ji,jj,jk)   &
168                       &   * rday / rfact2 * fse3t(ji,jj,jk) / 2. / (wsbio4(ji,jj,jk) + rtrn)
169                       zpoc = max(0., zpoc)
170                       alphag(ji,jj,jk,jn) = alphag(ji,jj,jk-1,jn) * exp( -reminp(jn)                &
171                       &   * ( zsizek + zsizek1 ) ) * zpoc + ( prodgoc(ji,jj,jk-1)  &
172                       &   / tgfunc(ji,jj,jk-1) * ( 1. - exp( -reminp(jn) * zsizek ) )                 &
173                       &   * exp( -reminp(jn) * zsizek1 ) + prodgoc(ji,jj,jk) / tgfunc(ji,jj,jk)       &
174                       &   * ( 1. - exp( -reminp(jn) * zsizek1 ) ) ) * rday / rfact2 * alphan(jn)     &
175                       &   / reminp(jn)
176                       alphat = alphat + alphag(ji,jj,jk,jn)
177                       remint = remint + alphag(ji,jj,jk,jn) * reminp(jn)
178                    END DO
179                  ENDIF
180                  DO jn = 1, jcpoc
181                     ! The contribution of each lability class at the current
182                     ! level is computed
183                     alphag(ji,jj,jk,jn) = alphag(ji,jj,jk,jn) / ( alphat + rtrn)
184                  END DO
185                  ! Computation of the mean remineralisation rate
186                  zremigoc(ji,jj,jk) = MIN(xremip, MAX(0., remint / ( alphat + rtrn) ))
187                ENDIF
188              ENDIF
189            END DO
190         END DO
191      END DO
192
193      DO jk = 1, jpkm1
194         DO jj = 1, jpj
195            DO ji = 1, jpi
196               zstep   = xstep
197# if defined key_degrad
198               zstep = zstep * facvol(ji,jj,jk)
199# endif
200               ! POC disaggregation by turbulence and bacterial activity.
201               ! --------------------------------------------------------
202               zremig = zremigoc(ji,jj,jk) * zstep * tgfunc(ji,jj,jk)
203               zorem2  = zremig * trb(ji,jj,jk,jpgoc)
204               zorem3(ji,jj,jk) = zremig * solgoc * trb(ji,jj,jk,jpgoc)
205               zofer2 = zremig * trb(ji,jj,jk,jpbfe)
206               zofer3 = zremig * solgoc * trb(ji,jj,jk,jpbfe)
207
208               ! Update the appropriate tracers trends
209               ! -------------------------------------
210
211               tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zorem3(ji,jj,jk)
212               tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) - zorem2 - zorem3(ji,jj,jk)
213               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zofer3
214               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) - zofer2 - zofer3
215               tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zorem2
216               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zofer2
217            END DO
218         END DO
219      END DO
220
221     IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
222        WRITE(charout, FMT="('poc1')")
223        CALL prt_ctl_trc_info(charout)
224        CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
225     ENDIF
226
227     ! ------------------------------------------------------------------
228     ! Lability parameterization for the small OM particles. This param
229     ! is based on the same theoretical background as the big particles.
230     ! However, because of its low sinking speed, lability is not supposed
231     ! to be equal to its initial value (the value of the freshly produced
232     ! organic matter). It is however uniform in the mixed layer.
233     ! -------------------------------------------------------------------
234     !
235     totprod(:,:) = 0.
236     totthick(:,:) = 0.
237     totcons(:,:) = 0.
238     ! intregrated production and consumption of POC in the mixed layer
239     ! ----------------------------------------------------------------
240     !
241     DO jk = 1, jpkm1
242        DO jj = 1, jpj
243           DO ji = 1, jpi
244              IF (tmask(ji,jj,jk) == 1.) THEN
245                zdep = hmld(ji,jj)
246                IF( fsdept(ji,jj,jk) <= zdep ) THEN
247                  totprod(ji,jj) = totprod(ji,jj) + prodpoc(ji,jj,jk) * fse3t(ji,jj,jk) * rday/ rfact2
248                  ! The temperature effect is included here
249                  totthick(ji,jj) = totthick(ji,jj) + fse3t(ji,jj,jk)* tgfunc(ji,jj,jk)
250                  totcons(ji,jj) = totcons(ji,jj) - conspoc(ji,jj,jk) * fse3t(ji,jj,jk) * rday/ rfact2    &
251                  &                / ( trb(ji,jj,jk,jppoc) + rtrn )
252                ENDIF
253              ENDIF
254           END DO
255        END DO
256     END DO
257
258     ! Computation of the lability spectrum in the mixed layer. In the mixed
259     ! layer, this spectrum is supposed to be uniform.
260     ! ---------------------------------------------------------------------
261     DO jk = 1, jpkm1
262        DO jj = 1, jpj
263           DO ji = 1, jpi
264              IF (tmask(ji,jj,jk) == 1.) THEN
265                zdep = hmld(ji,jj)
266                alphat = 0.0
267                remint = 0.0
268                IF( fsdept(ji,jj,jk) <= zdep ) THEN
269                   DO jn = 1, jcpoc
270                      ! For each lability class, the system is supposed to be
271                      ! at equilibrium: Prod - Sink - w alphap = 0.
272                      alphap(ji,jj,jk,jn) = totprod(ji,jj) * alphan(jn) / ( reminp(jn)    &
273                      &                     * totthick(ji,jj) + totcons(ji,jj) + wsbio + rtrn )
274                      alphat = alphat + alphap(ji,jj,jk,jn)
275                      remint = remint + alphap(ji,jj,jk,jn) * reminp(jn)
276                   END DO
277                   DO jn = 1, jcpoc
278                      alphap(ji,jj,jk,jn) = alphap(ji,jj,jk,jn) / ( alphat + rtrn)
279                   END DO
280                   ! Mean remineralization rate in the mixed layer
281                   zremipoc(ji,jj,jk) = MIN(xremip, MAX(0., remint / ( alphat + rtrn) ))
282                ENDIF
283              ENDIF
284           END DO
285        END DO
286     END DO
287     !
288     ! -----------------------------------------------------------------------
289     ! The lability parameterization is used here. The code is here
290     ! almost identical to what is done for big particles. The only difference
291     ! is that an additional source from GOC to POC is included. This means
292     ! that since we need the lability spectrum of GOC, GOC spectrum
293     ! should be determined before.
294     ! -----------------------------------------------------------------------
295     !
296     DO jk = 2, jpkm1
297        DO jj = 1, jpj
298           DO ji = 1, jpi
299              IF (tmask(ji,jj,jk) == 1.) THEN
300                zdep = hmld(ji,jj)
301                IF( fsdept(ji,jj,jk) > zdep ) THEN
302                  alphat = 0.
303                  remint = 0.
304                  !
305                  ! Special treatment of the level just below the MXL
306                  ! See the comments in the GOC section
307                  ! ---------------------------------------------------
308                  !
309                  IF ( fsdept(ji,jj,jk-1) <= zdep ) THEN
310                    DO jn = 1, jcpoc
311                       ! the scale factor is corrected with temperature
312                       zsizek1 = fse3t(ji,jj,jk) / 2. / (wsbio3(ji,jj,jk) + rtrn) * tgfunc(ji,jj,jk)
313                       ! Computation of the POC concentration using the
314                       ! lagrangian algorithm
315                       zpoc = trb(ji,jj,jk-1,jppoc) + conspoc(ji,jj,jk) * rday / rfact2               &
316                       &   * fse3t(ji,jj,jk) / 2. / (wsbio3(ji,jj,jk) + rtrn)
317                       zpoc = max(0., zpoc)
318                       ! computation of the lability spectrum applying the
319                       ! different sources and sinks
320                       alphap(ji,jj,jk,jn) = alphap(ji,jj,jk-1,jn) + zorem3(ji,jj,jk)                 &
321                       &   * ( 1. - exp( -reminp(jn) * zsizek1 ) ) * rday / rfact2                   &
322                       &   * alphag(ji,jj,jk,jn) / reminp(jn) / tgfunc(ji,jj,jk)
323                       alphat = alphat + alphap(ji,jj,jk,jn)
324                       remint = remint + alphap(ji,jj,jk,jn) * reminp(jn)
325                    END DO
326                  ELSE
327                    !
328                    ! Lability parameterization for the interior of the ocean
329                    ! This is very similar to what is done in the previous
330                    ! block
331                    ! --------------------------------------------------------
332                    !
333                    DO jn = 1, jcpoc
334                       zsizek  = fse3t(ji,jj,jk-1) / 2. / (wsbio3(ji,jj,jk-1) + rtrn) * tgfunc(ji,jj,jk-1)
335                       zsizek1 = fse3t(ji,jj,jk) / 2. / (wsbio3(ji,jj,jk) + rtrn) * tgfunc(ji,jj,jk)
336                       zpoc = trb(ji,jj,jk-1,jppoc) + conspoc(ji,jj,jk-1) * rday / rfact2               &
337                       &   * fse3t(ji,jj,jk-1) / 2. / (wsbio3(ji,jj,jk-1) + rtrn) + conspoc(ji,jj,jk)   &
338                       &   * rday / rfact2 * fse3t(ji,jj,jk) / 2. / (wsbio3(ji,jj,jk) + rtrn)
339                       zpoc = max(0., zpoc)
340                       alphap(ji,jj,jk,jn) = alphap(ji,jj,jk-1,jn) * exp( -reminp(jn)             &
341                       &   * ( zsizek + zsizek1 ) ) * zpoc + ( prodpoc(ji,jj,jk-1)  &
342                       &   / tgfunc(ji,jj,jk-1) * ( 1. - exp( -reminp(jn) * zsizek ) )              &
343                       &   * exp( -reminp(jn) * zsizek1 ) + prodpoc(ji,jj,jk) / tgfunc(ji,jj,jk)    &
344                       &   * ( 1. - exp( -reminp(jn) * zsizek1 ) ) ) * rday / rfact2 * alphan(jn)  &
345                       &   / reminp(jn)
346                       alphap(ji,jj,jk,jn) = alphap(ji,jj,jk,jn) + zorem3(ji,jj,jk-1)              &
347                       &  * alphag(ji,jj,jk-1,jn) / tgfunc(ji,jj,jk-1) * rday / rfact2 * ( 1.       &
348                       &  - exp( -reminp(jn) * zsizek ) ) * exp( -reminp(jn) * zsizek1 )                  &
349                       &  + zorem3(ji,jj,jk) * ( 1. - exp( -reminp(jn) * zsizek1 ) ) * rday         &
350                       &  / rfact2 * alphag(ji,jj,jk,jn) / reminp(jn) / tgfunc(ji,jj,jk)
351                       alphat = alphat + alphap(ji,jj,jk,jn)
352                       remint = remint + alphap(ji,jj,jk,jn) * reminp(jn)
353                    END DO
354                  ENDIF
355                  ! Normalization of the lability spectrum so that the
356                  ! integral is equal to 1
357                  DO jn = 1, jcpoc
358                     alphap(ji,jj,jk,jn) = alphap(ji,jj,jk,jn) / ( alphat + rtrn)
359                  END DO
360                   ! Mean remineralization rate in the water column
361                   zremipoc(ji,jj,jk) = MIN(xremip, MAX(0., remint / ( alphat + rtrn) ))
362                ENDIF
363              ENDIF
364            END DO
365         END DO
366      END DO
367#endif
368
369
370      DO jk = 1, jpkm1
371         DO jj = 1, jpj
372            DO ji = 1, jpi
373               IF (tmask(ji,jj,jk) == 1.) THEN
374                 zstep   = xstep
375# if defined key_degrad
376                 zstep = zstep * facvol(ji,jj,jk)
377# endif
378                 ! POC disaggregation by turbulence and bacterial activity.
379                 ! --------------------------------------------------------
380                 zremip          = zremipoc(ji,jj,jk) * zstep * tgfunc(ji,jj,jk)
381                 orem(ji,jj,jk)  = zremip * trb(ji,jj,jk,jppoc)
382                 zofer           = zremip * trb(ji,jj,jk,jpsfe)
383#if defined key_kriest
384                 zorem2          = zremip * trb(ji,jj,jk,jpnum)
385#endif
386
387                 ! Update the appropriate tracers trends
388                 ! -------------------------------------
389
390                 tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + orem(ji,jj,jk)
391                 tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zofer
392                 tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - orem(ji,jj,jk)
393#if defined key_kriest
394                 tra(ji,jj,jk,jpnum) = tra(ji,jj,jk,jpnum) - zorem2
395#endif
396                 tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zofer
397
398               ENDIF
399            END DO
400         END DO
401      END DO
402
403      IF( ln_diatrc .AND. lk_iomput .AND. jnt == nrdttrc ) THEN
404          zrfact2 = 1.e3 * rfact2r
405          CALL iom_put( "REMINP" , zremipoc(:,:,:) * tmask(:,:,:) )  ! Remineralisation rate
406          CALL iom_put( "REMING" , zremigoc(:,:,:) * tmask(:,:,:) )  ! Remineralisation rate
407      ENDIF
408
409      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
410         WRITE(charout, FMT="('poc2')")
411         CALL prt_ctl_trc_info(charout)
412         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
413      ENDIF
414      !
415      CALL wrk_dealloc( jpi, jpj,      totprod, totthick, totcons )
416      CALL wrk_dealloc( jpi, jpj, jpk, zremipoc, zremigoc, zorem3 )
417      DEALLOCATE( alphag )
418      !
419      IF( nn_timing == 1 )  CALL timing_stop('p4z_poc')
420      !
421   END SUBROUTINE p4z_poc
422
423
424   SUBROUTINE p4z_poc_init
425      !!----------------------------------------------------------------------
426      !!                  ***  ROUTINE p4z_poc_init  ***
427      !!
428      !! ** Purpose :   Initialization of remineralization parameters
429      !!
430      !! ** Method  :   Read the nampispoc namelist and check the parameters
431      !!      called at the first timestep
432      !!
433      !! ** input   :   Namelist nampispoc
434      !!
435      !!----------------------------------------------------------------------
436      INTEGER :: jn
437      REAL(wp) :: remindelta, reminup, remindown
438      NAMELIST/nampispoc/ xremip, jcpoc
439      INTEGER :: ios                 ! Local integer output status for namelist read
440
441      REWIND( numnatp_ref )              ! Namelist nampisrem in reference namelist : Pisces remineralization
442      READ  ( numnatp_ref, nampispoc, IOSTAT = ios, ERR = 901)
443901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampispoc in reference namelist', lwp )
444
445      REWIND( numnatp_cfg )              ! Namelist nampisrem in configuration namelist : Pisces remineralization
446      READ  ( numnatp_cfg, nampispoc, IOSTAT = ios, ERR = 902 )
447902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampispoc in configuration namelist', lwp )
448      IF(lwm) WRITE ( numonp, nampispoc )
449
450      IF(lwp) THEN                         ! control print
451         WRITE(numout,*) ' '
452         WRITE(numout,*) ' Namelist parameters for remineralization, nampispoc'
453         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
454         WRITE(numout,*) '    remineralisation rate of POC              xremip    =', xremip
455         WRITE(numout,*) '    Number of lability classes for POC        jcpoc     =', jcpoc
456      ENDIF
457      !
458      ! Discretization along the lability space
459      ! ---------------------------------------
460      !
461      ALLOCATE( alphan(jcpoc), reminp(jcpoc) )
462      ALLOCATE( alphap(jpi,jpj,jpk,jcpoc) )
463      !
464      IF (jcpoc > 1) THEN
465         remindelta = log(4. * 1000. ) / float(jcpoc-1)
466         reminup = xremip/400. * exp(remindelta)
467         alphan(1) = 1. - exp(-reminup/xremip)
468         reminp(1) = xremip - reminup * exp(-reminup/xremip) / alphan(1)
469         DO jn = 2, jcpoc-1
470            reminup = xremip/400. * exp(float(jn) * remindelta)
471            remindown = xremip / 400. * exp(float(jn-1) * remindelta)
472            alphan(jn) = exp(-remindown /xremip) - exp(-reminup/xremip)
473            reminp(jn) = xremip + (remindown * exp(-remindown /xremip)    &
474            &            - reminup * exp(-reminup/xremip) ) / alphan(jn)
475         END DO
476         remindown = xremip/400. * exp(float(jcpoc-1) * remindelta)
477         alphan(jcpoc) = exp(-remindown /xremip)
478         reminp(jcpoc) = xremip + remindown
479      ELSE
480         alphan(jcpoc) = 1.
481         reminp(jcpoc) = xremip
482      ENDIF
483
484      DO jn = 1, jcpoc
485         alphap(:,:,:,jn) = alphan(jn)
486      END DO
487
488   END SUBROUTINE p4z_poc_init
489
490#else
491   !!======================================================================
492   !!  Dummy module :                                   No PISCES bio-model
493   !!======================================================================
494CONTAINS
495   SUBROUTINE p4z_poc                    ! Empty routine
496   END SUBROUTINE p4z_poc
497#endif 
498
499   !!======================================================================
500END MODULE p4zpoc
Note: See TracBrowser for help on using the repository browser.