Changeset 10362


Ignore:
Timestamp:
2018-11-30T16:38:17+01:00 (22 months ago)
Author:
aumont
Message:

Various bug fixes and improvements

Location:
NEMO/trunk/src/TOP/PISCES
Files:
22 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zfechem.F90

    r10069 r10362  
    2727   LOGICAL          ::   ln_fechem    !: boolean for complex iron chemistry following Tagliabue and voelker 
    2828   LOGICAL          ::   ln_ligvar    !: boolean for variable ligand concentration following Tagliabue and voelker 
    29    LOGICAL          ::   ln_fecolloid !: boolean for variable colloidal fraction 
    3029   REAL(wp), PUBLIC ::   xlam1        !: scavenging rate of Iron  
    3130   REAL(wp), PUBLIC ::   xlamdust     !: scavenging rate of Iron by dust  
     
    6867      REAL(wp) ::   za, zb, zc, zkappa1, zkappa2, za0, za1, za2 
    6968      REAL(wp) ::   zxs, zfunc, zp, zq, zd, zr, zphi, zfff, zp3, zq2 
    70       REAL(wp) ::   ztfe, zoxy, zhplus 
     69      REAL(wp) ::   ztfe, zoxy, zhplus, zxlam 
    7170      REAL(wp) ::   zaggliga, zaggligb 
    7271      REAL(wp) ::   dissol, zligco 
     72      REAL(wp) :: zrfact2 
    7373      CHARACTER (len=25) :: charout 
    7474      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zTL1, zFe3, ztotlig, precip, zFeL1 
     75      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zcoll3d, zscav3d, zlcoll3d 
    7576      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zFeL2, zTL2, zFe2, zFeP 
    7677      REAL(wp), ALLOCATABLE, DIMENSION(:,:  ) ::   zstrn, zstrn2 
     
    129130         ! ------------------------------------------------------------ 
    130131         DO jn = 1, 2 
    131           DO jk = 1, jpkm1 
    132             DO jj = 1, jpj 
    133                DO ji = 1, jpi 
    134                   zlight = etot(ji,jj,jk) * zstrn(ji,jj) * REAL( 2-jn, wp ) 
    135                   zzstrn2 = zstrn2(ji,jj) * REAL( 2-jn, wp ) + (1. - zstrn2(ji,jj) ) * REAL( jn-1, wp ) 
    136                   ! Calculate ligand concentrations : assume 2/3rd of excess goes to 
    137                   ! strong ligands (L1) and 1/3rd to weak ligands (L2) 
    138                   ztligand       = ztotlig(ji,jj,jk) - ligand * 1E9 
    139                   zTL1(ji,jj,jk) =                0.000001 + 0.67 * ztligand 
    140                   zTL2(ji,jj,jk) = ligand * 1E9 - 0.000001 + 0.33 * ztligand 
    141                   ! ionic strength from Millero et al. 1987 
    142                   zph    = -LOG10( MAX( hi(ji,jj,jk), rtrn) ) 
    143                   zoxy   = trb(ji,jj,jk,jpoxy) 
    144                   ! Fe2+ oxydation rate from Santana-Casiano et al. (2005) 
    145                   zkox   = 35.407 - 6.7109 * zph + 0.5342 * zph * zph - 5362.6 / ( tempis(ji,jj,jk) + 273.15 )  & 
    146                     &    - 0.04406 * SQRT( salinprac(ji,jj,jk) ) - 0.002847 * salinprac(ji,jj,jk) 
    147                   zkox   = ( 10.** zkox ) * spd 
    148                   zkox   = zkox * MAX( 1.e-6, zoxy) / ( chemo2(ji,jj,jk) + rtrn ) 
    149                   ! PHOTOREDUCTION of complexed iron : Tagliabue and Arrigo (2006) 
    150                   zkph2 = MAX( 0., 15. * zlight / ( zlight + 2. ) ) * (1. - fr_i(ji,jj)) 
    151                   zkph1 = zkph2 / 5. 
    152                   ! pass the dfe concentration from PISCES 
    153                   ztfe = trb(ji,jj,jk,jpfer) * 1e9 
    154                   ! ---------------------------------------------------------- 
    155                   ! ANALYTICAL SOLUTION OF ROOTS OF THE FE3+ EQUATION 
    156                   ! As shown in Tagliabue and Voelker (2009), Fe3+ is the root of a 3rd order polynom.  
    157                   ! ---------------------------------------------------------- 
    158                   ! calculate some parameters 
    159                   za = 1 + ks / kpr 
    160                   zb = 1 + ( zkph1 + kth ) / ( zkox + rtrn ) 
    161                   zc = 1 + zkph2 / ( zkox + rtrn ) 
    162                   zkappa1 = ( kb1 + zkph1 + kth ) / kl1 
    163                   zkappa2 = ( kb2 + zkph2 ) / kl2 
    164                   za2 = zTL1(ji,jj,jk) * zb / za + zTL2(ji,jj,jk) * zc / za + zkappa1 + zkappa2 - ztfe / za 
    165                   za1 = zkappa2 * zTL1(ji,jj,jk) * zb / za + zkappa1 * zTL2(ji,jj,jk) * zc / za & 
     132            DO jk = 1, jpkm1 
     133               DO jj = 1, jpj 
     134                  DO ji = 1, jpi 
     135                     zlight = etot(ji,jj,jk) * zstrn(ji,jj) * REAL( 2-jn, wp ) 
     136                     zzstrn2 = zstrn2(ji,jj) * REAL( 2-jn, wp ) + (1. - zstrn2(ji,jj) ) * REAL( jn-1, wp ) 
     137                     ! Calculate ligand concentrations : assume 2/3rd of excess goes to 
     138                     ! strong ligands (L1) and 1/3rd to weak ligands (L2) 
     139                     ztligand       = ztotlig(ji,jj,jk) - ligand * 1E9 
     140                     zTL1(ji,jj,jk) =                0.000001 + 0.67 * ztligand 
     141                     zTL2(ji,jj,jk) = ligand * 1E9 - 0.000001 + 0.33 * ztligand 
     142                     ! ionic strength from Millero et al. 1987 
     143                     zph    = -LOG10( MAX( hi(ji,jj,jk), rtrn) ) 
     144                     zoxy   = trb(ji,jj,jk,jpoxy) 
     145                     ! Fe2+ oxydation rate from Santana-Casiano et al. (2005) 
     146                     zkox   = 35.407 - 6.7109 * zph + 0.5342 * zph * zph - 5362.6 / ( tempis(ji,jj,jk) + 273.15 )  & 
     147                     &        - 0.04406 * SQRT( salinprac(ji,jj,jk) ) - 0.002847 * salinprac(ji,jj,jk) 
     148                     zkox   = ( 10.** zkox ) * spd 
     149                     zkox   = zkox * MAX( 1.e-6, zoxy) / ( chemo2(ji,jj,jk) + rtrn ) 
     150                     ! PHOTOREDUCTION of complexed iron : Tagliabue and Arrigo (2006) 
     151                     zkph2 = MAX( 0., 15. * zlight / ( zlight + 2. ) ) * (1. - fr_i(ji,jj)) 
     152                     zkph1 = zkph2 / 5. 
     153                     ! pass the dfe concentration from PISCES 
     154                     ztfe = trb(ji,jj,jk,jpfer) * 1e9 
     155                     ! ---------------------------------------------------------- 
     156                     ! ANALYTICAL SOLUTION OF ROOTS OF THE FE3+ EQUATION 
     157                     ! As shown in Tagliabue and Voelker (2009), Fe3+ is the root of a 3rd order polynom.  
     158                     ! ---------------------------------------------------------- 
     159                     ! calculate some parameters 
     160                     za = 1.0 + ks / kpr 
     161                     zb = 1.0 + zkph2 / ( zkox ) 
     162                     zc = 1.0 + ( zkph1 + kth ) / ( zkox ) 
     163                     zkappa1 = ( kb1 + zkph1 + kth ) / kl1 
     164                     zkappa2 = ( kb2 + zkph2 ) / kl2 
     165                     za2 = zTL2(ji,jj,jk) * zb / za + zTL2(ji,jj,jk) * zc / za + zkappa1 + zkappa2 - ztfe / za 
     166                     za1 = zkappa1 * zTL2(ji,jj,jk) * zb / za + zkappa2 * zTL1(ji,jj,jk) * zc / za & 
    166167                      & + zkappa1 * zkappa2 - ( zkappa1 + zkappa2 ) * ztfe / za 
    167                   za0 = -zkappa1 * zkappa2 * ztfe / za 
    168                   zp  = za1 - za2 * za2 / 3. 
    169                   zq  = za2 * za2 * za2 * 2. / 27. - za2 * za1 / 3. + za0 
    170                   zp3 = zp / 3. 
    171                   zq2 = zq / 2. 
    172                   zd  = zp3 * zp3 * zp3 + zq2 * zq2 
    173                   zr  = zq / ABS( zq ) * SQRT( ABS( zp ) / 3. ) 
    174                   ! compute the roots 
    175                   IF( zp > 0.) THEN 
    176                      ! zphi = ASINH( zq / ( 2. * zr * zr * zr ) ) 
    177                      zphi =  zq / ( 2. * zr * zr * zr )  
    178                      zphi = LOG( zphi + SQRT( zphi * zphi + 1 ) )  ! asinh(x) = log(x + sqrt(x^2+1)) 
    179                      zxs  = -2. * zr * SINH( zphi / 3. ) - za1 / 3. 
    180                   ELSE 
    181                      IF( zd > 0. ) THEN 
    182                         zfff = MAX( 1., zq / ( 2. * zr * zr * zr ) ) 
    183                         ! zphi = ACOSH( zfff ) 
    184                         zphi = LOG( zfff + SQRT( zfff * zfff - 1 ) )  ! acosh(x) = log(x + sqrt(x^2-1)) 
    185                         zxs = -2. * zr * COSH( zphi / 3. ) - za1 / 3. 
     168                     za0 = -zkappa1 * zkappa2 * ztfe / za 
     169                     zp  = za1 - za2 * za2 / 3. 
     170                     zq  = za2 * za2 * za2 * 2. / 27. - za2 * za1 / 3. + za0 
     171                     zp3 = zp / 3. 
     172                     zq2 = zq / 2. 
     173                     zd  = zp3 * zp3 * zp3 + zq2 * zq2 
     174                     zr  = zq / ABS( zq ) * SQRT( ABS( zp ) / 3. ) 
     175                     ! compute the roots 
     176                     IF( zp > 0.) THEN 
     177                        ! zphi = ASINH( zq / ( 2. * zr * zr * zr ) ) 
     178                        zphi =  zq / ( 2. * zr * zr * zr )  
     179                        zphi = LOG( zphi + SQRT( zphi * zphi + 1 ) )  ! asinh(x) = log(x + sqrt(x^2+1)) 
     180                        zxs  = -2. * zr * SINH( zphi / 3. ) - za1 / 3. 
    186181                     ELSE 
    187                         zfff = MIN( 1., zq / ( 2. * zr * zr * zr ) ) 
    188                         zphi = ACOS( zfff ) 
    189                         DO jic = 1, 3 
    190                            zfunc = -2 * zr * COS( zphi / 3. + 2. * REAL( jic - 1, wp ) * rpi / 3. ) - za2 / 3. 
    191                            IF( zfunc > 0. .AND. zfunc <= ztfe)  zxs = zfunc 
    192                         END DO 
     182                        IF( zd > 0. ) THEN 
     183                           zfff = MAX( 1., zq / ( 2. * zr * zr * zr ) ) 
     184                           ! zphi = ACOSH( zfff ) 
     185                           zphi = LOG( zfff + SQRT( zfff * zfff - 1 ) )  ! acosh(x) = log(x + sqrt(x^2-1)) 
     186                           zxs = -2. * zr * COSH( zphi / 3. ) - za1 / 3. 
     187                        ELSE 
     188                           zfff = MIN( 1., zq / ( 2. * zr * zr * zr ) ) 
     189                           zphi = ACOS( zfff ) 
     190                           DO jic = 1, 3 
     191                              zfunc = -2 * zr * COS( zphi / 3. + 2. * REAL( jic - 1, wp ) * rpi / 3. ) - za2 / 3. 
     192                              IF( zfunc > 0. .AND. zfunc <= ztfe)  zxs = zfunc 
     193                           END DO 
     194                        ENDIF 
    193195                     ENDIF 
    194                   ENDIF 
    195                   ! solve for the other Fe species 
    196                   zzFe3 = MAX( 0., zxs ) 
    197                   zzFep = MAX( 0., ( ks * zzFe3 / kpr ) ) 
    198                   zkappa2 = ( kb2 + zkph2 ) / kl2 
    199                   zzFeL2 = MAX( 0., ( zzFe3 * zTL2(ji,jj,jk) ) / ( zkappa2 + zzFe3 ) ) 
    200                   zzFeL1 = MAX( 0., ( ztfe / zb - za / zb * zzFe3 - zc / zb * zzFeL2 ) ) 
    201                   zzFe2  = MAX( 0., ( ( zkph1 * zzFeL1 + zkph2 * zzFeL2 ) / zkox ) ) 
    202                   zFe3(ji,jj,jk)  = zFe3(ji,jj,jk)  + zzFe3 * zzstrn2 
    203                   zFe2(ji,jj,jk)  = zFe2(ji,jj,jk)  + zzFe2 * zzstrn2 
    204                   zFeL2(ji,jj,jk) = zFeL2(ji,jj,jk) + zzFeL2 * zzstrn2 
    205                   zFeL1(ji,jj,jk) = zFeL1(ji,jj,jk) + zzFeL1 * zzstrn2 
    206                   zFeP(ji,jj,jk)  = zFeP(ji,jj,jk)  + zzFeP * zzstrn2 
     196                     ! solve for the other Fe species 
     197                     zzFe3  = MAX( 0., zxs ) 
     198                     zzFep  = MAX( 0., ( ks * zzFe3 / kpr ) ) 
     199                     zzFeL1 = MAX( 0., ( zzFe3 * zTL1(ji,jj,jk) ) / ( zkappa1 + zzFe3 ) ) 
     200                     zzFeL2 = (ztfe - za * zzFe3 - zc * zzFeL1 ) / zb 
     201                     zzFe2  = MAX( 0., ( ( ( zkph1 + kth ) * zzFeL1 + zkph2 * zzFeL2 ) / zkox ) ) 
     202                     zzFep  = ztfe - zzFe3 - zzFe2 - zzFeL1 - zzFeL2 
     203                     zFe3(ji,jj,jk)  = zFe3(ji,jj,jk)  + zzFe3 * zzstrn2 
     204                     zFe2(ji,jj,jk)  = zFe2(ji,jj,jk)  + zzFe2 * zzstrn2 
     205                     zFeL2(ji,jj,jk) = zFeL2(ji,jj,jk) + zzFeL2 * zzstrn2 
     206                     zFeL1(ji,jj,jk) = zFeL1(ji,jj,jk) + zzFeL1 * zzstrn2 
     207                     zFeP(ji,jj,jk)  = zFeP(ji,jj,jk)  + zzFeP * zzstrn2 
     208                  END DO 
    207209               END DO 
    208210            END DO 
    209          END DO 
    210211         END DO 
    211212      ELSE 
     
    218219            DO jj = 1, jpj 
    219220               DO ji = 1, jpi 
    220                   zTL1(ji,jj,jk) = ztotlig(ji,jj,jk) 
    221                   zkeq           = fekeq(ji,jj,jk) 
    222                   zfesatur       = zTL1(ji,jj,jk) * 1E-9 
    223                   ztfe           = trb(ji,jj,jk,jpfer)  
     221                  zTL1(ji,jj,jk)  = ztotlig(ji,jj,jk) 
     222                  zkeq            = fekeq(ji,jj,jk) 
     223                  zfesatur        = zTL1(ji,jj,jk) * 1E-9 
     224                  ztfe            = trb(ji,jj,jk,jpfer)  
    224225                  ! Fe' is the root of a 2nd order polynom 
    225226                  zFe3 (ji,jj,jk) = ( -( 1. + zfesatur * zkeq - zkeq * ztfe )               & 
    226                      &             + SQRT( ( 1. + zfesatur * zkeq - zkeq * ztfe )**2       & 
    227                      &               + 4. * ztfe * zkeq) ) / ( 2. * zkeq ) 
     227                     &              + SQRT( ( 1. + zfesatur * zkeq - zkeq * ztfe )**2       & 
     228                     &              + 4. * ztfe * zkeq) ) / ( 2. * zkeq ) 
    228229                  zFe3 (ji,jj,jk) = zFe3(ji,jj,jk) * 1E9 
    229230                  zFeL1(ji,jj,jk) = MAX( 0., trb(ji,jj,jk,jpfer) * 1E9 - zFe3(ji,jj,jk) ) 
     
    242243               ! Scavenging onto dust is also included as evidenced from the DUNE experiments. 
    243244               ! -------------------------------------------------------------------------------------- 
     245               zhplus  = max( rtrn, hi(ji,jj,jk) ) 
     246               fe3sol  = fesol(ji,jj,jk,1) * ( zhplus**3 + fesol(ji,jj,jk,2) * zhplus**2  & 
     247               &         + fesol(ji,jj,jk,3) * zhplus + fesol(ji,jj,jk,4)     & 
     248               &         + fesol(ji,jj,jk,5) / zhplus ) 
    244249               IF( ln_fechem ) THEN 
    245250                  zfeequi = ( zFe3(ji,jj,jk) + zFe2(ji,jj,jk) + zFeP(ji,jj,jk) ) * 1E-9 
    246251                  zfecoll = ( 0.3 * zFeL1(ji,jj,jk) + 0.5 * zFeL2(ji,jj,jk) ) * 1E-9 
     252                  precip(ji,jj,jk) = 0.0 
    247253               ELSE 
    248254                  zfeequi = zFe3(ji,jj,jk) * 1E-9 
    249                   IF (ln_fecolloid) THEN 
    250                      zhplus   = max( rtrn, hi(ji,jj,jk) ) 
    251                      fe3sol  = fesol(ji,jj,jk,1) * ( zhplus**3 + fesol(ji,jj,jk,2) * zhplus**2  & 
    252                      &         + fesol(ji,jj,jk,3) * zhplus + fesol(ji,jj,jk,4)     & 
    253                      &         + fesol(ji,jj,jk,5) / zhplus ) 
    254                      zfecoll = max( ( 0.1 * zFeL1(ji,jj,jk) * 1E-9 ), ( zFeL1(ji,jj,jk) * 1E-9 -fe3sol ) ) 
    255                   ELSE 
    256                      zfecoll = 0.5 * zFeL1(ji,jj,jk) * 1E-9 
    257                      fe3sol  = 0. 
    258                   ENDIF 
     255                  zhplus  = max( rtrn, hi(ji,jj,jk) ) 
     256                  fe3sol  = fesol(ji,jj,jk,1) * ( zhplus**3 + fesol(ji,jj,jk,2) * zhplus**2  & 
     257                  &         + fesol(ji,jj,jk,3) * zhplus + fesol(ji,jj,jk,4)     & 
     258                  &         + fesol(ji,jj,jk,5) / zhplus ) 
     259                  zfecoll = 0.5 * zFeL1(ji,jj,jk) * 1E-9 
     260                  ! precipitation of Fe3+, creation of nanoparticles 
     261                  precip(ji,jj,jk) = MAX( 0., ( zFe3(ji,jj,jk) * 1E-9 - fe3sol ) ) * kfep * xstep 
    259262               ENDIF 
    260263               ! 
    261264               ztrc   = ( trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + trb(ji,jj,jk,jpcal) + trb(ji,jj,jk,jpgsi) ) * 1.e6  
    262                IF( ln_dust )  zdust  = dust(ji,jj) / ( wdust / rday ) * tmask(ji,jj,jk) ! dust in kg/m2/s 
    263                zlam1b = 3.e-5 + xlamdust * zdust + xlam1 * ztrc 
     265               IF( ln_dust )  zdust  = dust(ji,jj) / ( wdust / rday ) * tmask(ji,jj,jk) & 
     266               &  * EXP( -gdept_n(ji,jj,jk) / 540. ) 
     267               IF (ln_ligand) THEN 
     268                  zxlam  = xlam1 * MAX( 1.E-3, EXP(-2 * etot(ji,jj,jk) / 10. ) * (1. - EXP(-2 * trb(ji,jj,jk,jpoxy) / 100.E-6 ) )) 
     269               ELSE 
     270                  zxlam  = xlam1 * 1.0 
     271               ENDIF 
     272               zlam1b = 3.e-5 + xlamdust * zdust + zxlam * ztrc 
    264273               zscave = zfeequi * zlam1b * xstep 
    265274 
     
    267276               ! to later allocate scavenged iron to the different organic pools 
    268277               ! --------------------------------------------------------- 
    269                zdenom1 = xlam1 * trb(ji,jj,jk,jppoc) / zlam1b 
    270                zdenom2 = xlam1 * trb(ji,jj,jk,jpgoc) / zlam1b 
     278               zdenom1 = zxlam * trb(ji,jj,jk,jppoc) / zlam1b 
     279               zdenom2 = zxlam * trb(ji,jj,jk,jpgoc) / zlam1b 
    271280 
    272281               !  Increased scavenging for very high iron concentrations found near the coasts  
     
    276285               zlamfac = MIN( 1.  , zlamfac ) 
    277286               zdep    = MIN( 1., 1000. / gdept_n(ji,jj,jk) ) 
    278                zlam1b  = xlam1 * MAX( 0.e0, ( trb(ji,jj,jk,jpfer) * 1.e9 - ztotlig(ji,jj,jk) ) ) 
    279                zcoag   = zfeequi * zlam1b * xstep + 1E-4 * ( 1. - zlamfac ) * zdep * xstep * trb(ji,jj,jk,jpfer) 
     287               zcoag   = 1E-4 * ( 1. - zlamfac ) * zdep * xstep * trb(ji,jj,jk,jpfer) 
    280288 
    281289               !  Compute the coagulation of colloidal iron. This parameterization  
     
    283291               !  It requires certainly some more work as it is very poorly constrained. 
    284292               !  ---------------------------------------------------------------- 
    285                zlam1a  = ( 0.369  * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4  * trb(ji,jj,jk,jppoc) ) * xdiss(ji,jj,jk)    & 
    286                    &   + ( 114.   * 0.3 * trb(ji,jj,jk,jpdoc) ) 
     293               zlam1a   = ( 0.369  * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4  * trb(ji,jj,jk,jppoc) ) * xdiss(ji,jj,jk)    & 
     294                   &      + ( 114.   * 0.3 * trb(ji,jj,jk,jpdoc) ) 
    287295               zaggdfea = zlam1a * xstep * zfecoll 
    288296               ! 
    289                zlam1b = 3.53E3 *  trb(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk) 
     297               zlam1b   = 3.53E3 * trb(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk) 
    290298               zaggdfeb = zlam1b * xstep * zfecoll 
    291                ! 
    292                ! precipitation of Fe3+, creation of nanoparticles 
    293                precip(ji,jj,jk) = MAX( 0., ( zfeequi - fe3sol ) ) * kfep * xstep 
    294299               ! 
    295300               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zscave - zaggdfea - zaggdfeb & 
     
    297302               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zscave * zdenom1 + zaggdfea 
    298303               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zscave * zdenom2 + zaggdfeb 
     304               zscav3d(ji,jj,jk)   = zscave 
     305               zcoll3d(ji,jj,jk)   = zaggdfea + zaggdfeb 
    299306               ! 
    300307            END DO 
     
    317324                  ! 
    318325                  zlam1b   = 3.53E3 *   trb(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk) 
    319                   zligco   = MAX( ( 0.1 * trb(ji,jj,jk,jplgw) ), ( trb(ji,jj,jk,jplgw) - fe3sol ) ) 
     326                  zligco   = 0.5 * trn(ji,jj,jk,jplgw) 
    320327                  zaggliga = zlam1a * xstep * zligco 
    321328                  zaggligb = zlam1b * xstep * zligco 
    322329                  tra(ji,jj,jk,jpfep) = tra(ji,jj,jk,jpfep) + precip(ji,jj,jk) 
    323330                  tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) - zaggliga - zaggligb 
     331                  zlcoll3d(ji,jj,jk)  = zaggliga + zaggligb 
    324332               END DO 
    325333            END DO 
     
    328336         IF( .NOT.ln_fechem) THEN 
    329337            plig(:,:,:) =  MAX( 0., ( ( zFeL1(:,:,:) * 1E-9 ) / ( trb(:,:,:,jpfer) +rtrn ) ) ) 
    330             plig(:,:,:) =  MAX( 0. , plig(:,:,:) ) 
    331338         ENDIF 
    332339         ! 
     
    336343      IF( lk_iomput ) THEN 
    337344         IF( knt == nrdttrc ) THEN 
     345         zrfact2 = 1.e3 * rfact2r  ! conversion from mol/L/timestep into mol/m3/s 
    338346         IF( iom_use("Fe3")    )  CALL iom_put("Fe3"    , zFe3   (:,:,:)       * tmask(:,:,:) )   ! Fe3+ 
    339347         IF( iom_use("FeL1")   )  CALL iom_put("FeL1"   , zFeL1  (:,:,:)       * tmask(:,:,:) )   ! FeL1 
    340348         IF( iom_use("TL1")    )  CALL iom_put("TL1"    , zTL1   (:,:,:)       * tmask(:,:,:) )   ! TL1 
    341349         IF( iom_use("Totlig") )  CALL iom_put("Totlig" , ztotlig(:,:,:)       * tmask(:,:,:) )   ! TL 
    342          IF( iom_use("Biron")  )  CALL iom_put("Biron"  , biron  (:,:,:) * 1e9 * tmask(:,:,:) )   ! biron 
     350         IF( iom_use("Biron")  )  CALL iom_put("Biron"  , biron  (:,:,:)  * 1e9 * tmask(:,:,:) )   ! biron 
     351         IF( iom_use("FESCAV") )  CALL iom_put("FESCAV" , zscav3d(:,:,:)  * 1e9 * tmask(:,:,:) * zrfact2 ) 
     352         IF( iom_use("FECOLL") )  CALL iom_put("FECOLL" , zcoll3d(:,:,:)  * 1e9 * tmask(:,:,:) * zrfact2 ) 
     353         IF( iom_use("LGWCOLL"))  CALL iom_put("LGWCOLL", zlcoll3d(:,:,:) * 1e9 * tmask(:,:,:) * zrfact2 ) 
    343354         IF( ln_fechem ) THEN 
    344355            IF( iom_use("Fe2")  ) CALL iom_put("Fe2"    , zFe2   (:,:,:)       * tmask(:,:,:) )   ! Fe2+ 
     
    380391      INTEGER ::   ios   ! Local integer  
    381392      !! 
    382       NAMELIST/nampisfer/ ln_fechem, ln_ligvar, ln_fecolloid, xlam1, xlamdust, ligand, kfep  
     393      NAMELIST/nampisfer/ ln_fechem, ln_ligvar, xlam1, xlamdust, ligand, kfep  
    383394      !!---------------------------------------------------------------------- 
    384395      ! 
     
    401412         WRITE(numout,*) '      enable complex iron chemistry scheme      ln_fechem    =', ln_fechem 
    402413         WRITE(numout,*) '      variable concentration of ligand          ln_ligvar    =', ln_ligvar 
    403          WRITE(numout,*) '      Variable colloidal fraction of Fe3+       ln_fecolloid =', ln_fecolloid 
    404414         WRITE(numout,*) '      scavenging rate of Iron                   xlam1        =', xlam1 
    405415         WRITE(numout,*) '      scavenging rate of Iron by dust           xlamdust     =', xlamdust 
     
    407417         WRITE(numout,*) '      rate constant for nanoparticle formation  kfep         =', kfep 
    408418      ENDIF 
     419      !  
     420      IF (ln_ligand .AND. ln_fechem) CALL ctl_stop( 'STOP', 'p4z_fechem_init: ln_ligand and ln_fechem are incompatible') 
    409421      ! 
    410422      IF( ln_fechem ) THEN             ! set some constants used by the complexe chemistry scheme 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zligand.F90

    r10069 r10362  
    1313   USE sms_pisces      ! PISCES Source Minus Sink variables 
    1414   USE prtctl_trc      ! print control for debugging 
     15   USE iom             !  I/O manager 
    1516 
    1617   IMPLICIT NONE 
     
    4344      INTEGER  ::   ji, jj, jk 
    4445      REAL(wp) ::   zlgwp, zlgwpr, zlgwr, zlablgw, zrfepa, zfepr 
     46      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zligrem, zligpr, zrligprod 
     47      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zw3d 
    4548      CHARACTER (len=25) ::   charout 
    4649      !!--------------------------------------------------------------------- 
     
    5659               ! ------------------------------------------------------------------ 
    5760               ! production from remineralisation of organic matter 
    58                zlgwp  = orem(ji,jj,jk) * rlig 
     61               zlgwp = orem(ji,jj,jk) * rlig 
    5962               ! decay of weak ligand 
    6063               ! This is based on the idea that as LGW is lower 
    6164               ! there is a larger fraction of refractory OM 
    6265               zlgwr = max( rlgs , rlgw * exp( -2 * (trb(ji,jj,jk,jplgw)*1e9) ) ) ! years 
    63                zlgwr = 1. / zlgwr * tgfunc(ji,jj,jk) * ( xstep / nyear_len(1) ) * trb(ji,jj,jk,jplgw) 
     66               zlgwr = 1. / zlgwr * tgfunc(ji,jj,jk) * ( xstep / nyear_len(1) ) * blim(ji,jj,jk) * trb(ji,jj,jk,jplgw) 
    6467               ! photochem loss of weak ligand 
    6568               zlgwpr = prlgw * xstep * etot(ji,jj,jk) * trb(ji,jj,jk,jplgw) * (1. - fr_i(ji,jj)) 
    6669               tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zlgwp - zlgwr - zlgwpr 
     70               zligrem(ji,jj,jk)   = zlgwr 
     71               zligpr(ji,jj,jk)    = zlgwpr 
     72               zrligprod(ji,jj,jk) = zlgwp 
    6773               ! 
    6874               ! ---------------------------------------------------------- 
     
    7480               zrfepa = rfep * ( 1. - EXP( -1. * etot(ji,jj,jk) / 25. ) ) * (1.- fr_i(ji,jj)) 
    7581               zrfepa = MAX( (zrfepa / 10.0), zrfepa ) ! min of 10 days lifetime 
    76                zfepr = rfep * xstep * trb(ji,jj,jk,jpfep) 
     82               zfepr  = rfep * xstep * trb(ji,jj,jk,jpfep) 
    7783               tra(ji,jj,jk,jpfep) = tra(ji,jj,jk,jpfep) - zfepr 
    7884               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zfepr 
     
    8187         END DO 
    8288      END DO 
     89      ! 
     90      !  Output of some diagnostics variables 
     91      !     --------------------------------- 
     92      IF( lk_iomput .AND. knt == nrdttrc ) THEN 
     93         ALLOCATE( zw3d(jpi,jpj,jpk) ) 
     94         IF( iom_use( "LIGREM" ) ) THEN 
     95            zw3d(:,:,:) = zligrem(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) 
     96            CALL iom_put( "LIGREM", zw3d ) 
     97         ENDIF 
     98         IF( iom_use( "LIGPR" ) ) THEN 
     99            zw3d(:,:,:) = zligpr(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)  
     100            CALL iom_put( "LIGPR", zw3d ) 
     101         ENDIF 
     102         IF( iom_use( "LPRODR" ) ) THEN 
     103            zw3d(:,:,:) = zrligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)  
     104            CALL iom_put( "LPRODR", zw3d ) 
     105         ENDIF 
     106         DEALLOCATE( zw3d ) 
     107      ENDIF 
    83108      ! 
    84109      IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zmeso.F90

    r10222 r10362  
    2525 
    2626   REAL(wp), PUBLIC ::  part2        !: part of calcite not dissolved in mesozoo guts 
    27    REAL(wp), PUBLIC ::  xprefc       !: mesozoo preference for POC  
    28    REAL(wp), PUBLIC ::  xprefp       !: mesozoo preference for nanophyto 
    29    REAL(wp), PUBLIC ::  xprefz       !: mesozoo preference for diatoms 
    30    REAL(wp), PUBLIC ::  xprefpoc     !: mesozoo preference for POC  
     27   REAL(wp), PUBLIC ::  xpref2d      !: mesozoo preference for diatoms 
     28   REAL(wp), PUBLIC ::  xpref2n      !: mesozoo preference for nanophyto 
     29   REAL(wp), PUBLIC ::  xpref2z      !: mesozoo preference for microzooplankton 
     30   REAL(wp), PUBLIC ::  xpref2c      !: mesozoo preference for POC  
    3131   REAL(wp), PUBLIC ::  xthresh2zoo  !: zoo feeding threshold for mesozooplankton  
    3232   REAL(wp), PUBLIC ::  xthresh2dia  !: diatoms feeding threshold for mesozooplankton  
     
    4040   REAL(wp), PUBLIC ::  unass2       !: Efficicency of mesozoo growth  
    4141   REAL(wp), PUBLIC ::  sigma2       !: Fraction of mesozoo excretion as DOM  
    42    REAL(wp), PUBLIC ::  epsher2      !: half sturation constant for grazing 2 
     42   REAL(wp), PUBLIC ::  epsher2      !: growth efficiency 
     43   REAL(wp), PUBLIC ::  epsher2min   !: minimum growth efficiency at high food for grazing 2 
    4344   REAL(wp), PUBLIC ::  grazflux     !: mesozoo flux feeding rate 
    4445 
     
    6364      REAL(wp) :: zcompadi, zcompaph, zcompapoc, zcompaz, zcompam 
    6465      REAL(wp) :: zgraze2 , zdenom, zdenom2 
    65       REAL(wp) :: zfact   , zfood, zfoodlim, zproport 
     66      REAL(wp) :: zfact   , zfood, zfoodlim, zproport, zbeta 
    6667      REAL(wp) :: zmortzgoc, zfrac, zfracfe, zratio, zratio2, zfracal, zgrazcal 
    67       REAL(wp) :: zepshert, zepsherv, zgrarsig, zgraztot, zgraztotn, zgraztotf 
    68       REAL(wp) :: zgrarem2, zgrafer2, zgrapoc2, zprcaca, zmortz2, zgrasrat, zgrasratn 
    69       REAL(wp) :: zrespz2, ztortz2, zgrazd, zgrazz, zgrazpof 
     68      REAL(wp) :: zepsherf, zepshert, zepsherv, zgrarsig, zgraztotc, zgraztotn, zgraztotf 
     69      REAL(wp) :: zgrarem2, zgrafer2, zgrapoc2, zprcaca, zmortz, zgrasrat, zgrasratn 
     70      REAL(wp) :: zrespz, ztortz, zgrazd, zgrazz, zgrazpof 
    7071      REAL(wp) :: zgrazn, zgrazpoc, zgraznf, zgrazf 
    7172      REAL(wp) :: zgrazfffp, zgrazfffg, zgrazffep, zgrazffeg 
    7273      CHARACTER (len=25) :: charout 
    73       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing 
    74       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zw3d 
     74      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing, zfezoo2 
     75      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zw3d, zz2ligprod 
    7576      !!--------------------------------------------------------------------- 
    7677      ! 
     
    7879      ! 
    7980      zgrazing(:,:,:) = 0._wp 
    80  
     81      zfezoo2 (:,:,:) = 0._wp 
     82      ! 
     83      IF (ln_ligand) THEN 
     84         ALLOCATE( zz2ligprod(jpi,jpj,jpk) ) 
     85         zz2ligprod(:,:,:) = 0._wp 
     86      ENDIF 
     87      ! 
    8188      DO jk = 1, jpkm1 
    8289         DO jj = 1, jpj 
     
    8794               !  Respiration rates of both zooplankton 
    8895               !  ------------------------------------- 
    89                zrespz2   = resrat2 * zfact * trb(ji,jj,jk,jpmes) / ( xkmort + trb(ji,jj,jk,jpmes) )  & 
    90                   &      + resrat2 * zfact * 3. * nitrfac(ji,jj,jk) 
     96               zrespz    = resrat2 * zfact * ( trb(ji,jj,jk,jpmes) / ( xkmort + trb(ji,jj,jk,jpmes) )  & 
     97               &           + 3. * nitrfac(ji,jj,jk) ) 
    9198 
    9299               !  Zooplankton mortality. A square function has been selected with 
    93100               !  no real reason except that it seems to be more stable and may mimic predation 
    94101               !  --------------------------------------------------------------- 
    95                ztortz2   = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes)  * (1. - nitrfac(ji,jj,jk) ) 
     102               ztortz    = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes)  * (1. - nitrfac(ji,jj,jk) ) 
    96103               ! 
    97104               zcompadi  = MAX( ( trb(ji,jj,jk,jpdia) - xthresh2dia ), 0.e0 ) 
    98105               zcompaz   = MAX( ( trb(ji,jj,jk,jpzoo) - xthresh2zoo ), 0.e0 ) 
     106               zcompapoc = MAX( ( trb(ji,jj,jk,jppoc) - xthresh2poc ), 0.e0 ) 
    99107               ! Size effect of nanophytoplankton on grazing : the smaller it is, the less prone 
    100108               ! it is to predation by mesozooplankton 
     
    102110               zcompaph  = MAX( ( trb(ji,jj,jk,jpphy) - xthresh2phy ), 0.e0 ) & 
    103111                  &      * MIN(1., MAX( 0., ( quotan(ji,jj,jk) - 0.2) / 0.3 ) ) 
    104                zcompapoc = MAX( ( trb(ji,jj,jk,jppoc) - xthresh2poc ), 0.e0 ) 
    105  
    106                zfood     = xprefc * zcompadi + xprefz * zcompaz + xprefp * zcompaph + xprefpoc * zcompapoc  
     112 
     113               !   Mesozooplankton grazing 
     114               !   ------------------------ 
     115               zfood     = xpref2d * zcompadi + xpref2z * zcompaz + xpref2n * zcompaph + xpref2c * zcompapoc  
    107116               zfoodlim  = MAX( 0., zfood - MIN( 0.5 * zfood, xthresh2 ) ) 
    108117               zdenom    = zfoodlim / ( xkgraz2 + zfoodlim ) 
     
    110119               zgraze2   = grazrat2 * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpmes) * (1. - nitrfac(ji,jj,jk))  
    111120 
    112                zgrazd    = zgraze2  * xprefc   * zcompadi  * zdenom2  
    113                zgrazz    = zgraze2  * xprefz   * zcompaz   * zdenom2  
    114                zgrazn    = zgraze2  * xprefp   * zcompaph  * zdenom2  
    115                zgrazpoc  = zgraze2  * xprefpoc * zcompapoc * zdenom2  
     121               zgrazd    = zgraze2  * xpref2d  * zcompadi  * zdenom2  
     122               zgrazz    = zgraze2  * xpref2z  * zcompaz   * zdenom2  
     123               zgrazn    = zgraze2  * xpref2n  * zcompaph  * zdenom2  
     124               zgrazpoc  = zgraze2  * xpref2c * zcompapoc * zdenom2  
    116125 
    117126               zgraznf   = zgrazn   * trb(ji,jj,jk,jpnfe) / ( trb(ji,jj,jk,jpphy) + rtrn) 
     
    129138               &           * (1. - nitrfac(ji,jj,jk)) 
    130139               zgrazfffp = zgrazffep * trb(ji,jj,jk,jpsfe) / (trb(ji,jj,jk,jppoc) + rtrn) 
    131               ! 
    132               zgraztot = zgrazd + zgrazz + zgrazn + zgrazpoc + zgrazffep + zgrazffeg 
    133               ! Compute the proportion of filter feeders 
    134               zproport  = (zgrazffep + zgrazffeg)/(rtrn + zgraztot) 
    135               ! Compute fractionation of aggregates. It is assumed that  
    136               ! diatoms based aggregates are more prone to fractionation 
    137               ! since they are more porous (marine snow instead of fecal pellets) 
    138               zratio    = trb(ji,jj,jk,jpgsi) / ( trb(ji,jj,jk,jpgoc) + rtrn ) 
    139               zratio2   = zratio * zratio 
    140               zfrac     = zproport * grazflux  * xstep * wsbio4(ji,jj,jk)      & 
     140               ! 
     141               zgraztotc = zgrazd + zgrazz + zgrazn + zgrazpoc + zgrazffep + zgrazffeg 
     142               ! Compute the proportion of filter feeders 
     143               zproport  = (zgrazffep + zgrazffeg)/(rtrn + zgraztotc) 
     144               ! Compute fractionation of aggregates. It is assumed that  
     145               ! diatoms based aggregates are more prone to fractionation 
     146               ! since they are more porous (marine snow instead of fecal pellets) 
     147               zratio    = trb(ji,jj,jk,jpgsi) / ( trb(ji,jj,jk,jpgoc) + rtrn ) 
     148               zratio2   = zratio * zratio 
     149               zfrac     = zproport * grazflux  * xstep * wsbio4(ji,jj,jk)      & 
    141150               &          * trb(ji,jj,jk,jpgoc) * trb(ji,jj,jk,jpmes)          & 
    142151               &          * ( 0.2 + 3.8 * zratio2 / ( 1.**2 + zratio2 ) ) 
    143               zfracfe   = zfrac * trb(ji,jj,jk,jpbfe) / (trb(ji,jj,jk,jpgoc) + rtrn) 
    144  
    145               zgrazffep = zproport * zgrazffep 
    146               zgrazffeg = zproport * zgrazffeg 
    147               zgrazfffp = zproport * zgrazfffp 
    148               zgrazfffg = zproport * zgrazfffg 
    149               zgraztot = zgrazd + zgrazz + zgrazn + zgrazpoc + zgrazffep + zgrazffeg 
    150               zgraztotn = zgrazd * quotad(ji,jj,jk) + zgrazz + zgrazn * quotan(ji,jj,jk)   & 
    151               &   + zgrazpoc + zgrazffep + zgrazffeg 
    152               zgraztotf = zgrazf + zgraznf + zgrazz * ferat3 + zgrazpof + zgrazfffp + zgrazfffg 
    153  
    154               ! Total grazing ( grazing by microzoo is already computed in p4zmicro ) 
    155               zgrazing(ji,jj,jk) = zgraztot 
    156  
    157               !    Mesozooplankton efficiency 
    158               !    -------------------------- 
    159                zgrasrat  =  ( zgraztotf +rtrn )/ ( zgraztot + rtrn ) 
    160                zgrasratn =  ( zgraztotn +rtrn )/ ( zgraztot + rtrn ) 
     152               zfracfe   = zfrac * trb(ji,jj,jk,jpbfe) / (trb(ji,jj,jk,jpgoc) + rtrn) 
     153 
     154               zgrazffep = zproport * zgrazffep 
     155               zgrazffeg = zproport * zgrazffeg 
     156               zgrazfffp = zproport * zgrazfffp 
     157               zgrazfffg = zproport * zgrazfffg 
     158               zgraztotc = zgrazd + zgrazz + zgrazn + zgrazpoc + zgrazffep + zgrazffeg 
     159               zgraztotn = zgrazd * quotad(ji,jj,jk) + zgrazz + zgrazn * quotan(ji,jj,jk)   & 
     160               &   + zgrazpoc + zgrazffep + zgrazffeg 
     161               zgraztotf = zgrazf + zgraznf + zgrazz * ferat3 + zgrazpof + zgrazfffp + zgrazfffg 
     162 
     163               ! Total grazing ( grazing by microzoo is already computed in p4zmicro ) 
     164               zgrazing(ji,jj,jk) = zgraztotc 
     165 
     166               !    Mesozooplankton efficiency 
     167               !    -------------------------- 
     168               zgrasrat  =  ( zgraztotf + rtrn )/ ( zgraztotc + rtrn ) 
     169               zgrasratn =  ( zgraztotn + rtrn )/ ( zgraztotc + rtrn ) 
    161170               zepshert  = MIN( 1., zgrasratn, zgrasrat / ferat3) 
    162                zepsherv  = zepshert * MIN( epsher2, (1. - unass2) * zgrasrat / ferat3, (1. - unass2) * zgrasratn ) 
    163                zgrarem2  = zgraztot * ( 1. - zepsherv - unass2 ) & 
    164                 &       + ( 1. - epsher2 - unass2 ) / ( 1. - epsher2 ) * ztortz2 
    165                zgrafer2  = zgraztot * MAX( 0. , ( 1. - unass2 ) * zgrasrat - ferat3 * zepsherv )    & 
    166                 &       + ferat3 * ( ( 1. - epsher2 - unass2 ) /( 1. - epsher2 ) * ztortz2 ) 
    167                zgrapoc2  = zgraztot * unass2 
     171               zbeta     = MAX(0., (epsher2 - epsher2min) ) 
     172               zepsherf  = epsher2min + zbeta / ( 1.0 + 0.04E6 * 12. * zfood * zbeta )  
     173               zepsherv  = zepsherf * zepshert  
     174 
     175               zgrarem2  = zgraztotc * ( 1. - zepsherv - unass2 ) & 
     176               &         + ( 1. - epsher2 - unass2 ) / ( 1. - epsher2 ) * ztortz 
     177               zgrafer2  = zgraztotc * MAX( 0. , ( 1. - unass2 ) * zgrasrat - ferat3 * zepsherv )    & 
     178               &         + ferat3 * ( ( 1. - epsher2 - unass2 ) /( 1. - epsher2 ) * ztortz ) 
     179               zgrapoc2  = zgraztotc * unass2 
    168180 
    169181               !   Update the arrays TRA which contain the biological sources and sinks 
     
    173185               tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zgrarem2 - zgrarsig 
    174186               ! 
    175                IF( ln_ligand ) tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + (zgrarem2 - zgrarsig) * ldocz 
     187               IF( ln_ligand ) THEN  
     188                  tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + (zgrarem2 - zgrarsig) * ldocz 
     189                  zz2ligprod(ji,jj,jk) = (zgrarem2 - zgrarsig) * ldocz 
     190               ENDIF 
    176191               ! 
    177192               tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - o2ut * zgrarsig 
    178193               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zgrafer2 
     194               zfezoo2(ji,jj,jk)   = zgrafer2 
    179195               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrarsig 
    180196               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * zgrarsig               
    181197 
    182                zmortz2 = ztortz2 + zrespz2 
    183                zmortzgoc = unass2 / ( 1. - epsher2 ) * ztortz2 + zrespz2 
    184                tra(ji,jj,jk,jpmes) = tra(ji,jj,jk,jpmes) - zmortz2 + zepsherv * zgraztot  
     198               zmortz = ztortz + zrespz 
     199               zmortzgoc = unass2 / ( 1. - epsher2 ) * ztortz + zrespz 
     200               tra(ji,jj,jk,jpmes) = tra(ji,jj,jk,jpmes) - zmortz + zepsherv * zgraztotc  
    185201               tra(ji,jj,jk,jpdia) = tra(ji,jj,jk,jpdia) - zgrazd 
    186202               tra(ji,jj,jk,jpzoo) = tra(ji,jj,jk,jpzoo) - zgrazz 
     
    193209               tra(ji,jj,jk,jpdfe) = tra(ji,jj,jk,jpdfe) - zgrazf 
    194210 
    195               tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - zgrazpoc - zgrazffep + zfrac 
    196               prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + zfrac 
    197               conspoc(ji,jj,jk) = conspoc(ji,jj,jk) - zgrazpoc - zgrazffep 
    198               tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zmortzgoc - zgrazffeg + zgrapoc2 - zfrac 
    199               prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zmortzgoc + zgrapoc2 
    200               consgoc(ji,jj,jk) = consgoc(ji,jj,jk) - zgrazffeg - zfrac 
    201               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zgrazpof - zgrazfffp + zfracfe 
    202               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + ferat3 * zmortzgoc - zgrazfffg     & 
     211               tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - zgrazpoc - zgrazffep + zfrac 
     212               prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + zfrac 
     213               conspoc(ji,jj,jk) = conspoc(ji,jj,jk) - zgrazpoc - zgrazffep 
     214               tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zmortzgoc - zgrazffeg + zgrapoc2 - zfrac 
     215               prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zmortzgoc + zgrapoc2 
     216               consgoc(ji,jj,jk) = consgoc(ji,jj,jk) - zgrazffeg - zfrac 
     217               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zgrazpof - zgrazfffp + zfracfe 
     218               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + ferat3 * zmortzgoc - zgrazfffg     & 
    203219                 &                + zgraztotf * unass2 - zfracfe 
    204               zfracal = trb(ji,jj,jk,jpcal) / (trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + rtrn ) 
    205               zgrazcal = (zgrazffeg + zgrazpoc) * (1. - part2) * zfracal 
    206               ! calcite production 
    207               zprcaca = xfracal(ji,jj,jk) * zgrazn 
    208               prodcal(ji,jj,jk) = prodcal(ji,jj,jk) + zprcaca  ! prodcal=prodcal(nanophy)+prodcal(microzoo)+prodcal(mesozoo) 
    209               ! 
    210               zprcaca = part2 * zprcaca 
    211               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrazcal - zprcaca 
    212               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) - 2. * ( zgrazcal + zprcaca ) 
    213               tra(ji,jj,jk,jpcal) = tra(ji,jj,jk,jpcal) - zgrazcal + zprcaca 
     220               zfracal = trb(ji,jj,jk,jpcal) / (trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + rtrn ) 
     221               zgrazcal = (zgrazffeg + zgrazpoc) * (1. - part2) * zfracal 
     222               ! calcite production 
     223               zprcaca = xfracal(ji,jj,jk) * zgrazn 
     224               prodcal(ji,jj,jk) = prodcal(ji,jj,jk) + zprcaca  ! prodcal=prodcal(nanophy)+prodcal(microzoo)+prodcal(mesozoo) 
     225               ! 
     226               zprcaca = part2 * zprcaca 
     227               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrazcal - zprcaca 
     228               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) - 2. * ( zgrazcal + zprcaca ) 
     229               tra(ji,jj,jk,jpcal) = tra(ji,jj,jk,jpcal) - zgrazcal + zprcaca 
    214230            END DO 
    215231         END DO 
     
    226242            CALL iom_put( "PCAL", zw3d )   
    227243         ENDIF 
     244         IF( iom_use( "FEZOO2" ) ) THEN 
     245            zw3d(:,:,:) = zfezoo2(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)   ! 
     246            CALL iom_put( "FEZOO2", zw3d ) 
     247         ENDIF 
     248         IF( iom_use( "LPRODZ2" ) .AND. ln_ligand )  THEN 
     249            zw3d(:,:,:) = zz2ligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) 
     250            CALL iom_put( "LPRODZ2"  , zw3d ) 
     251         ENDIF 
    228252         DEALLOCATE( zw3d ) 
    229253      ENDIF 
     
    253277      INTEGER ::   ios   ! Local integer 
    254278      ! 
    255       NAMELIST/namp4zmes/ part2, grazrat2, resrat2, mzrat2, xprefc, xprefp, xprefz,   & 
    256          &                xprefpoc, xthresh2dia, xthresh2phy, xthresh2zoo, xthresh2poc, & 
    257          &                xthresh2, xkgraz2, epsher2, sigma2, unass2, grazflux 
     279      NAMELIST/namp4zmes/ part2, grazrat2, resrat2, mzrat2, xpref2n, xpref2d, xpref2z,   & 
     280         &                xpref2c, xthresh2dia, xthresh2phy, xthresh2zoo, xthresh2poc, & 
     281         &                xthresh2, xkgraz2, epsher2, epsher2min, sigma2, unass2, grazflux 
    258282      !!---------------------------------------------------------------------- 
    259283      ! 
     
    275299         WRITE(numout,*) '   Namelist : namp4zmes' 
    276300         WRITE(numout,*) '      part of calcite not dissolved in mesozoo guts  part2        =', part2 
    277          WRITE(numout,*) '      mesozoo preference for phyto                   xprefc       =', xprefc 
    278          WRITE(numout,*) '      mesozoo preference for POC                     xprefp       =', xprefp 
    279          WRITE(numout,*) '      mesozoo preference for zoo                     xprefz       =', xprefz 
    280          WRITE(numout,*) '      mesozoo preference for poc                     xprefpoc     =', xprefpoc 
     301         WRITE(numout,*) '      mesozoo preference for phyto                   xpref2n      =', xpref2n 
     302         WRITE(numout,*) '      mesozoo preference for diatoms                 xpref2d      =', xpref2d 
     303         WRITE(numout,*) '      mesozoo preference for zoo                     xpref2z      =', xpref2z 
     304         WRITE(numout,*) '      mesozoo preference for poc                     xpref2c      =', xpref2c 
    281305         WRITE(numout,*) '      microzoo feeding threshold  for mesozoo        xthresh2zoo  =', xthresh2zoo 
    282306         WRITE(numout,*) '      diatoms feeding threshold  for mesozoo         xthresh2dia  =', xthresh2dia 
     
    289313         WRITE(numout,*) '      mesozoo flux feeding rate                      grazflux     =', grazflux 
    290314         WRITE(numout,*) '      non assimilated fraction of P by mesozoo       unass2       =', unass2 
    291          WRITE(numout,*) '      Efficicency of Mesozoo growth                  epsher2      =', epsher2 
     315         WRITE(numout,*) '      Efficiency of Mesozoo growth                   epsher2      =', epsher2 
     316         WRITE(numout,*) '      Minimum Efficiency of Mesozoo growth           epsher2min  =', epsher2min 
    292317         WRITE(numout,*) '      Fraction of mesozoo excretion as DOM           sigma2       =', sigma2 
    293318         WRITE(numout,*) '      half sturation constant for grazing 2          xkgraz2      =', xkgraz2 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zmicro.F90

    r10227 r10362  
    2626 
    2727   REAL(wp), PUBLIC ::   part        !: part of calcite not dissolved in microzoo guts 
    28    REAL(wp), PUBLIC ::   xpref2c     !: microzoo preference for POC  
    29    REAL(wp), PUBLIC ::   xpref2p     !: microzoo preference for nanophyto 
    30    REAL(wp), PUBLIC ::   xpref2d     !: microzoo preference for diatoms 
     28   REAL(wp), PUBLIC ::   xprefc      !: microzoo preference for POC  
     29   REAL(wp), PUBLIC ::   xprefn      !: microzoo preference for nanophyto 
     30   REAL(wp), PUBLIC ::   xprefd      !: microzoo preference for diatoms 
    3131   REAL(wp), PUBLIC ::   xthreshdia  !: diatoms feeding threshold for microzooplankton  
    3232   REAL(wp), PUBLIC ::   xthreshphy  !: nanophyto threshold for microzooplankton  
     
    3636   REAL(wp), PUBLIC ::   mzrat       !: microzooplankton mortality rate  
    3737   REAL(wp), PUBLIC ::   grazrat     !: maximal microzoo grazing rate 
    38    REAL(wp), PUBLIC ::   xkgraz      !: non assimilated fraction of P by microzoo  
    39    REAL(wp), PUBLIC ::   unass       !: Efficicency of microzoo growth  
     38   REAL(wp), PUBLIC ::   xkgraz      !: Half-saturation constant of assimilation 
     39   REAL(wp), PUBLIC ::   unass       !: Non-assimilated part of food 
    4040   REAL(wp), PUBLIC ::   sigma1      !: Fraction of microzoo excretion as DOM  
    41    REAL(wp), PUBLIC ::   epsher      !: half sturation constant for grazing 1  
     41   REAL(wp), PUBLIC ::   epsher      !: growth efficiency for grazing 1  
     42   REAL(wp), PUBLIC ::   epshermin   !: minimum growth efficiency for grazing 1 
    4243 
    4344   !!---------------------------------------------------------------------- 
     
    6263      REAL(wp) :: zcompadi, zcompaz , zcompaph, zcompapoc 
    6364      REAL(wp) :: zgraze  , zdenom, zdenom2 
    64       REAL(wp) :: zfact   , zfood, zfoodlim 
    65       REAL(wp) :: zepshert, zepsherv, zgrarsig, zgraztot, zgraztotn, zgraztotf 
     65      REAL(wp) :: zfact   , zfood, zfoodlim, zbeta 
     66      REAL(wp) :: zepsherf, zepshert, zepsherv, zgrarsig, zgraztotc, zgraztotn, zgraztotf 
    6667      REAL(wp) :: zgrarem, zgrafer, zgrapoc, zprcaca, zmortz 
    6768      REAL(wp) :: zrespz, ztortz, zgrasrat, zgrasratn 
    6869      REAL(wp) :: zgrazp, zgrazm, zgrazsd 
    6970      REAL(wp) :: zgrazmf, zgrazsf, zgrazpf 
    70       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing  
    71       REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zw3d 
     71      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing, zfezoo 
     72      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zw3d, zzligprod 
    7273      CHARACTER (len=25) :: charout 
    7374      !!--------------------------------------------------------------------- 
    7475      ! 
    7576      IF( ln_timing )   CALL timing_start('p4z_micro') 
     77      ! 
     78      IF (ln_ligand) THEN 
     79         ALLOCATE( zzligprod(jpi,jpj,jpk) ) 
     80         zzligprod(:,:,:) = 0._wp 
     81      ENDIF 
    7682      ! 
    7783      DO jk = 1, jpkm1 
     
    97103               !     Microzooplankton grazing 
    98104               !     ------------------------ 
    99                zfood     = xpref2p * zcompaph + xpref2c * zcompapoc + xpref2d * zcompadi 
     105               zfood     = xprefn * zcompaph + xprefc * zcompapoc + xprefd * zcompadi 
    100106               zfoodlim  = MAX( 0. , zfood - min(xthresh,0.5*zfood) ) 
    101107               zdenom    = zfoodlim / ( xkgraz + zfoodlim ) 
     
    103109               zgraze    = grazrat * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpzoo) * (1. - nitrfac(ji,jj,jk)) 
    104110 
    105                zgrazp    = zgraze  * xpref2p * zcompaph  * zdenom2  
    106                zgrazm    = zgraze  * xpref2c * zcompapoc * zdenom2  
    107                zgrazsd   = zgraze  * xpref2d * zcompadi  * zdenom2  
     111               zgrazp    = zgraze  * xprefn * zcompaph  * zdenom2  
     112               zgrazm    = zgraze  * xprefc * zcompapoc * zdenom2  
     113               zgrazsd   = zgraze  * xprefd * zcompadi  * zdenom2  
    108114 
    109115               zgrazpf   = zgrazp  * trb(ji,jj,jk,jpnfe) / (trb(ji,jj,jk,jpphy) + rtrn) 
     
    111117               zgrazsf   = zgrazsd * trb(ji,jj,jk,jpdfe) / (trb(ji,jj,jk,jpdia) + rtrn) 
    112118               ! 
    113                zgraztot  = zgrazp  + zgrazm  + zgrazsd  
     119               zgraztotc = zgrazp  + zgrazm  + zgrazsd  
    114120               zgraztotf = zgrazpf + zgrazsf + zgrazmf  
    115121               zgraztotn = zgrazp * quotan(ji,jj,jk) + zgrazm + zgrazsd * quotad(ji,jj,jk) 
    116122 
    117123               ! Grazing by microzooplankton 
    118                zgrazing(ji,jj,jk) = zgraztot 
     124               zgrazing(ji,jj,jk) = zgraztotc 
    119125 
    120126               !    Various remineralization and excretion terms 
    121127               !    -------------------------------------------- 
    122                zgrasrat  = ( zgraztotf + rtrn ) / ( zgraztot + rtrn ) 
    123                zgrasratn = ( zgraztotn + rtrn ) / ( zgraztot + rtrn ) 
     128               zgrasrat  = ( zgraztotf + rtrn ) / ( zgraztotc + rtrn ) 
     129               zgrasratn = ( zgraztotn + rtrn ) / ( zgraztotc + rtrn ) 
    124130               zepshert  =  MIN( 1., zgrasratn, zgrasrat / ferat3) 
    125                zepsherv  = zepshert * MIN( epsher, (1. - unass) * zgrasrat / ferat3, (1. - unass) * zgrasratn ) 
    126                zgrafer   = zgraztot * MAX( 0. , ( 1. - unass ) * zgrasrat - ferat3 * zepsherv )  
    127                zgrarem   = zgraztot * ( 1. - zepsherv - unass ) 
    128                zgrapoc   = zgraztot * unass 
     131               zbeta     = MAX(0., (epsher - epshermin) ) 
     132               zepsherf  = epshermin + zbeta / ( 1.0 + 0.04E6 * 12. * zfood * zbeta ) 
     133               zepsherv  = zepsherf * zepshert  
     134 
     135               zgrafer   = zgraztotc * MAX( 0. , ( 1. - unass ) * zgrasrat - ferat3 * zepsherv )  
     136               zgrarem   = zgraztotc * ( 1. - zepsherv - unass ) 
     137               zgrapoc   = zgraztotc * unass 
    129138 
    130139               !  Update of the TRA arrays 
     
    135144               tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zgrarem - zgrarsig 
    136145               ! 
    137                IF( ln_ligand ) tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + (zgrarem - zgrarsig) * ldocz 
     146               IF( ln_ligand ) THEN 
     147                  tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + (zgrarem - zgrarsig) * ldocz 
     148                  zzligprod(ji,jj,jk) = (zgrarem - zgrarsig) * ldocz 
     149               ENDIF 
    138150               ! 
    139151               tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - o2ut * zgrarsig 
    140152               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zgrafer 
     153               zfezoo(ji,jj,jk)    = zgrafer 
    141154               tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zgrapoc 
    142                prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + zgrapoc 
     155               prodpoc(ji,jj,jk)   = prodpoc(ji,jj,jk) + zgrapoc 
    143156               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zgraztotf * unass 
    144157               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrarsig 
     
    147160               !   -------------------------------------------------------------------- 
    148161               zmortz = ztortz + zrespz 
    149                tra(ji,jj,jk,jpzoo) = tra(ji,jj,jk,jpzoo) - zmortz + zepsherv * zgraztot  
     162               tra(ji,jj,jk,jpzoo) = tra(ji,jj,jk,jpzoo) - zmortz + zepsherv * zgraztotc  
    150163               tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) - zgrazp 
    151164               tra(ji,jj,jk,jpdia) = tra(ji,jj,jk,jpdia) - zgrazsd 
     
    180193              CALL iom_put( "GRAZ1", zw3d ) 
    181194           ENDIF 
     195           IF( iom_use( "FEZOO" ) ) THEN 
     196              zw3d(:,:,:) = zfezoo(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)   ! 
     197              CALL iom_put( "FEZOO", zw3d ) 
     198           ENDIF 
     199           IF( iom_use( "LPRODZ" ) .AND. ln_ligand )  THEN 
     200              zw3d(:,:,:) = zzligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) 
     201              CALL iom_put( "LPRODZ"  , zw3d ) 
     202           ENDIF 
    182203           DEALLOCATE( zw3d ) 
    183204         ENDIF 
     
    209230      INTEGER ::   ios   ! Local integer 
    210231      ! 
    211       NAMELIST/namp4zzoo/ part, grazrat, resrat, mzrat, xpref2c, xpref2p, & 
    212          &                xpref2d,  xthreshdia,  xthreshphy,  xthreshpoc, & 
    213          &                xthresh, xkgraz, epsher, sigma1, unass 
     232      NAMELIST/namp4zzoo/ part, grazrat, resrat, mzrat, xprefn, xprefc, & 
     233         &                xprefd,  xthreshdia,  xthreshphy,  xthreshpoc, & 
     234         &                xthresh, xkgraz, epsher, epshermin, sigma1, unass 
    214235      !!---------------------------------------------------------------------- 
    215236      ! 
     
    231252         WRITE(numout,*) '   Namelist : namp4zzoo' 
    232253         WRITE(numout,*) '      part of calcite not dissolved in microzoo guts  part        =', part 
    233          WRITE(numout,*) '      microzoo preference for POC                     xpref2c     =', xpref2c 
    234          WRITE(numout,*) '      microzoo preference for nano                    xpref2p     =', xpref2p 
    235          WRITE(numout,*) '      microzoo preference for diatoms                 xpref2d     =', xpref2d 
     254         WRITE(numout,*) '      microzoo preference for POC                     xprefc      =', xprefc 
     255         WRITE(numout,*) '      microzoo preference for nano                    xprefn      =', xprefn 
     256         WRITE(numout,*) '      microzoo preference for diatoms                 xprefd      =', xprefd 
    236257         WRITE(numout,*) '      diatoms feeding threshold  for microzoo         xthreshdia  =', xthreshdia 
    237258         WRITE(numout,*) '      nanophyto feeding threshold for microzoo        xthreshphy  =', xthreshphy 
     
    243264         WRITE(numout,*) '      non assimilated fraction of P by microzoo       unass       =', unass 
    244265         WRITE(numout,*) '      Efficicency of microzoo growth                  epsher      =', epsher 
     266         WRITE(numout,*) '      Minimum efficicency of microzoo growth          epshermin   =', epshermin 
    245267         WRITE(numout,*) '      Fraction of microzoo excretion as DOM           sigma1      =', sigma1 
    246268         WRITE(numout,*) '      half sturation constant for grazing 1           xkgraz      =', xkgraz 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zopt.F90

    r10069 r10362  
    109109         DO jk = 1, nksrp       
    110110            etot_ndcy(:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk) 
    111             enano    (:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
    112             ediat    (:,:,jk) =  1.6 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.7 * ze3(:,:,jk) 
     111            enano    (:,:,jk) =  1.85 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.46 * ze3(:,:,jk) 
     112            ediat    (:,:,jk) =  1.62 * ze1(:,:,jk) + 0.74 * ze2(:,:,jk) + 0.63 * ze3(:,:,jk) 
    113113         END DO 
    114114         IF( ln_p5z ) THEN 
    115115            DO jk = 1, nksrp       
    116               epico  (:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
     116              epico  (:,:,jk) =  1.94 * ze1(:,:,jk) + 0.66 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
    117117            END DO 
    118118         ENDIF 
     
    134134         DO jk = 1, nksrp       
    135135            etot (:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk) 
    136             enano(:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
    137             ediat(:,:,jk) =  1.6 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.7 * ze3(:,:,jk) 
     136            enano(:,:,jk) =  1.85 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.46 * ze3(:,:,jk) 
     137            ediat(:,:,jk) =  1.62 * ze1(:,:,jk) + 0.74 * ze2(:,:,jk) + 0.63 * ze3(:,:,jk) 
    138138         END DO 
    139139         IF( ln_p5z ) THEN 
    140140            DO jk = 1, nksrp       
    141               epico(:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
     141              epico(:,:,jk) =  1.94 * ze1(:,:,jk) + 0.66 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
    142142            END DO 
    143143         ENDIF 
     
    182182      zetmp1 (:,:)   = 0.e0 
    183183      zetmp2 (:,:)   = 0.e0 
    184       zetmp3 (:,:)   = 0.e0 
    185       zetmp4 (:,:)   = 0.e0 
    186184 
    187185      DO jk = 1, nksrp 
     
    191189                  zetmp1 (ji,jj) = zetmp1 (ji,jj) + etot     (ji,jj,jk) * e3t_n(ji,jj,jk) ! remineralisation 
    192190                  zetmp2 (ji,jj) = zetmp2 (ji,jj) + etot_ndcy(ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
    193                   zetmp3 (ji,jj) = zetmp3 (ji,jj) + enano    (ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
    194                   zetmp4 (ji,jj) = zetmp4 (ji,jj) + ediat    (ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
    195191                  zdepmoy(ji,jj) = zdepmoy(ji,jj) +                       e3t_n(ji,jj,jk) 
    196192               ENDIF 
     
    209205                  emoy (ji,jj,jk) = zetmp1(ji,jj) * z1_dep 
    210206                  zpar (ji,jj,jk) = zetmp2(ji,jj) * z1_dep 
    211                   enano(ji,jj,jk) = zetmp3(ji,jj) * z1_dep 
    212                   ediat(ji,jj,jk) = zetmp4(ji,jj) * z1_dep 
     207               ENDIF 
     208            END DO 
     209         END DO 
     210      END DO 
     211      ! 
     212      zdepmoy(:,:)   = 0.e0 
     213      zetmp3 (:,:)   = 0.e0 
     214      zetmp4 (:,:)   = 0.e0 
     215      ! 
     216      DO jk = 1, nksrp 
     217         DO jj = 1, jpj 
     218            DO ji = 1, jpi 
     219               IF( gdepw_n(ji,jj,jk+1) <= MIN(hmld(ji,jj), heup_01(ji,jj)) ) THEN 
     220                  zetmp3 (ji,jj) = zetmp3 (ji,jj) + enano    (ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
     221                  zetmp4 (ji,jj) = zetmp4 (ji,jj) + ediat    (ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
     222                  zdepmoy(ji,jj) = zdepmoy(ji,jj) +                       e3t_n(ji,jj,jk) 
     223               ENDIF 
     224            END DO 
     225         END DO 
     226      END DO 
     227      enanom(:,:,:) = enano(:,:,:) 
     228      ediatm(:,:,:) = ediat(:,:,:) 
     229      ! 
     230      DO jk = 1, nksrp 
     231         DO jj = 1, jpj 
     232            DO ji = 1, jpi 
     233               IF( gdepw_n(ji,jj,jk+1) <= hmld(ji,jj) ) THEN 
     234                  z1_dep = 1. / ( zdepmoy(ji,jj) + rtrn ) 
     235                  enanom(ji,jj,jk) = zetmp3(ji,jj) * z1_dep 
     236                  ediatm(ji,jj,jk) = zetmp4(ji,jj) * z1_dep 
    213237               ENDIF 
    214238            END DO 
     
    221245            DO jj = 1, jpj 
    222246               DO ji = 1, jpi 
    223                   IF( gdepw_n(ji,jj,jk+1) <= hmld(ji,jj) ) THEN  
     247                  IF( gdepw_n(ji,jj,jk+1) <= MIN(hmld(ji,jj), heup_01(ji,jj)) ) THEN 
     248                     zetmp5(ji,jj)  = zetmp5 (ji,jj) + epico(ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
     249                  ENDIF 
     250               END DO 
     251            END DO 
     252         END DO 
     253         ! 
     254         epicom(:,:,:) = epico(:,:,:) 
     255         ! 
     256         DO jk = 1, nksrp 
     257            DO jj = 1, jpj 
     258               DO ji = 1, jpi 
     259                  IF( gdepw_n(ji,jj,jk+1) <= hmld(ji,jj) ) THEN 
    224260                     z1_dep = 1. / ( zdepmoy(ji,jj) + rtrn ) 
    225                      zetmp5(ji,jj)  = zetmp5 (ji,jj) + epico(ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
    226                      epico(ji,jj,jk) = zetmp5(ji,jj) * z1_dep 
     261                     epicom(ji,jj,jk) = zetmp5(ji,jj) * z1_dep 
    227262                  ENDIF 
    228263               END DO 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zpoc.F90

    r10069 r10362  
    6262      CHARACTER (len=25) :: charout 
    6363      REAL(wp), DIMENSION(jpi,jpj  )   :: totprod, totthick, totcons  
    64       REAL(wp), DIMENSION(jpi,jpj,jpk)   :: zremipoc, zremigoc, zorem3, ztremint 
     64      REAL(wp), DIMENSION(jpi,jpj,jpk)   :: zremipoc, zremigoc, zorem3, ztremint, zfolimi 
    6565      REAL(wp), DIMENSION(jpi,jpj,jpk,jcpoc) :: alphag 
    6666      !!--------------------------------------------------------------------- 
     
    9090      orem  (:,:,:)   = 0. 
    9191      ztremint(:,:,:) = 0. 
     92      zfolimi (:,:,:) = 0. 
    9293 
    9394      DO jn = 1, jcpoc 
     
    209210                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zorem2 
    210211                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zofer2 
     212                  zfolimi(ji,jj,jk)   = zofer2 
    211213               END DO 
    212214            END DO 
     
    239241                  tra(ji,jj,jk,jpgop) = tra(ji,jj,jk,jpgop) - zopop2 * (1. + solgoc) 
    240242                  tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) - zofer2 * (1. + solgoc) 
     243                  zfolimi(ji,jj,jk)   = zofer2 
    241244               END DO 
    242245            END DO 
     
    414417                    tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - zorem 
    415418                    tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zofer 
     419                    zfolimi(ji,jj,jk)   = zfolimi(ji,jj,jk) + zofer 
    416420                  ENDIF 
    417421               END DO 
     
    439443                tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) + zopop  
    440444                tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zofer  
     445                zfolimi(ji,jj,jk)   = zfolimi(ji,jj,jk) + zofer 
    441446             END DO 
    442447           END DO 
     
    449454          CALL iom_put( "REMINP" , zremipoc(:,:,:) * tmask(:,:,:) )  ! Remineralisation rate 
    450455          CALL iom_put( "REMING" , zremigoc(:,:,:) * tmask(:,:,:) )  ! Remineralisation rate 
     456          CALL iom_put( "REMINF" , zfolimi(:,:,:)  * tmask(:,:,:)  * 1.e+9 * zrfact2 )  ! Remineralisation rate 
    451457        ENDIF 
    452458     ENDIF 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zprod.F90

    r10227 r10362  
    4040   REAL(wp), PUBLIC ::   grosip       !: 
    4141 
    42    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmax    !: optimal production = f(temperature) 
    4342   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   quotan   !: proxy of N quota in Nanophyto 
    4443   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   quotad   !: proxy of N quota in diatomee 
     
    7877      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
    7978      REAL(wp), DIMENSION(jpi,jpj    ) :: zstrn, zmixnano, zmixdiat 
     79      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprmaxn,zprmaxd 
    8080      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpislopeadn, zpislopeadd, zysopt   
    8181      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprdia, zprbio, zprdch, zprnch    
     
    8383      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpronewn, zpronewd 
    8484      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmxl_fac, zmxl_chl 
     85      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpligprod1, zpligprod2 
    8586      !!--------------------------------------------------------------------- 
    8687      ! 
     
    9697 
    9798      ! Computation of the optimal production 
    98       prmax(:,:,:) = 0.8_wp * r1_rday * tgfunc(:,:,:)  
     99      zprmaxn(:,:,:) = 0.8_wp * r1_rday * tgfunc(:,:,:) 
     100      zprmaxd(:,:,:) = zprmaxn(:,:,:) 
    99101 
    100102      ! compute the day length depending on latitude and the day 
     
    128130      END DO 
    129131 
    130       zprbio(:,:,:) = prmax(:,:,:) * zmxl_fac(:,:,:) 
    131       zprdia(:,:,:) = zprbio(:,:,:) 
     132      zprbio(:,:,:) = zprmaxn(:,:,:) * zmxl_fac(:,:,:) 
     133      zprdia(:,:,:) = zprmaxd(:,:,:) * zmxl_fac(:,:,:) 
    132134 
    133135      ! Maximum light intensity 
     
    169171                      !  Computation of production function for Chlorophyll 
    170172                      !-------------------------------------------------- 
    171                       zpislopen = zpislopeadn(ji,jj,jk) / ( prmax(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn ) 
    172                       zpisloped = zpislopeadd(ji,jj,jk) / ( prmax(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn ) 
    173                       zprnch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) ) ) 
    174                       zprdch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) ) ) 
     173                      zpislopen = zpislopeadn(ji,jj,jk) / ( zprmaxn(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn ) 
     174                      zpisloped = zpislopeadd(ji,jj,jk) / ( zprmaxd(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn ) 
     175                      zprnch(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) ) ) 
     176                      zprdch(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) ) ) 
    175177                  ENDIF 
    176178               END DO 
     
    192194                      zpislopen = zpislopen * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    193195                      zpisloped = zpisloped * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    194                       zprnch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) ) ) 
    195                       zprdch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) ) ) 
     196                      zprnch(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) ) ) 
     197                      zprdch(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) ) ) 
    196198                  ENDIF 
    197199               END DO 
     
    206208            DO ji = 1, jpi 
    207209                zval = MIN( xnanopo4(ji,jj,jk), ( xnanonh4(ji,jj,jk) + xnanono3(ji,jj,jk) ) )   & 
    208                 &      * prmax(ji,jj,jk) / ( zprbio(ji,jj,jk) + rtrn ) 
     210                &      * zprmaxn(ji,jj,jk) / ( zprbio(ji,jj,jk) + rtrn ) 
    209211                quotan(ji,jj,jk) = MIN( 1., 0.2 + 0.8 * zval ) 
    210212                zval = MIN( xdiatpo4(ji,jj,jk), ( xdiatnh4(ji,jj,jk) + xdiatno3(ji,jj,jk) ) )   & 
    211                 &      * prmax(ji,jj,jk) / ( zprdia(ji,jj,jk) + rtrn ) 
     213                &      * zprmaxd(ji,jj,jk) / ( zprdia(ji,jj,jk) + rtrn ) 
    212214                quotad(ji,jj,jk) = MIN( 1., 0.2 + 0.8 * zval ) 
    213215            END DO 
     
    227229                   !    to mimic the very high ratios observed in the Southern Ocean (silpot2) 
    228230                  zlim  = trb(ji,jj,jk,jpsil) / ( trb(ji,jj,jk,jpsil) + xksi1 ) 
    229                   zsilim = MIN( zprdia(ji,jj,jk) / ( prmax(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) ) 
     231                  zsilim = MIN( zprdia(ji,jj,jk) / ( zprmaxd(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) ) 
    230232                  zsilfac = 4.4 * EXP( -4.23 * zsilim ) * MAX( 0.e0, MIN( 1., 2.2 * ( zlim - 0.5 ) )  ) + 1.e0 
    231233                  zsiborn = trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) 
     
    266268                  zratio = trb(ji,jj,jk,jpnfe) / ( trb(ji,jj,jk,jpphy) * fecnm + rtrn ) 
    267269                  zmax   = MAX( 0., ( 1. - zratio ) / ABS( 1.05 - zratio ) )  
    268                   zprofen(ji,jj,jk) = fecnm * prmax(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  & 
     270                  zprofen(ji,jj,jk) = fecnm * zprmaxn(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  & 
    269271                  &             * ( 4. - 4.5 * xlimnfe(ji,jj,jk) / ( xlimnfe(ji,jj,jk) + 0.5 ) )    & 
    270272                  &             * biron(ji,jj,jk) / ( biron(ji,jj,jk) + concnfe(ji,jj,jk) )  & 
     
    276278                  zratio = trb(ji,jj,jk,jpdfe) / ( trb(ji,jj,jk,jpdia) * fecdm + rtrn ) 
    277279                  zmax   = MAX( 0., ( 1. - zratio ) / ABS( 1.05 - zratio ) )  
    278                   zprofed(ji,jj,jk) = fecdm * prmax(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  & 
     280                  zprofed(ji,jj,jk) = fecdm * zprmaxd(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  & 
    279281                  &             * ( 4. - 4.5 * xlimdfe(ji,jj,jk) / ( xlimdfe(ji,jj,jk) + 0.5 ) )    & 
    280282                  &             * biron(ji,jj,jk) / ( biron(ji,jj,jk) + concdfe(ji,jj,jk) )  & 
     
    291293               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
    292294                  !  production terms for nanophyto. ( chlorophyll ) 
    293                   znanotot = enano(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
     295                  znanotot = enanom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    294296                  zprod    = rday * zprorcan(ji,jj,jk) * zprnch(ji,jj,jk) * xlimphy(ji,jj,jk) 
    295297                  zprochln = chlcmin * 12. * zprorcan (ji,jj,jk) 
     
    298300                                        & (  zpislopeadn(ji,jj,jk) * znanotot +rtrn) 
    299301                  !  production terms for diatoms ( chlorophyll ) 
    300                   zdiattot = ediat(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
     302                  zdiattot = ediatm(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    301303                  zprod    = rday * zprorcad(ji,jj,jk) * zprdch(ji,jj,jk) * xlimdia(ji,jj,jk) 
    302304                  zprochld = chlcmin * 12. * zprorcad(ji,jj,jk) 
     
    351353                    zfeup    = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk) 
    352354                    tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zdocprod * ldocp - zfeup * plig(ji,jj,jk) * lthet 
     355                    zpligprod1(ji,jj,jk) = zdocprod * ldocp 
     356                    zpligprod2(ji,jj,jk) = zfeup * plig(ji,jj,jk) * lthet 
    353357                 ENDIF 
    354358              END DO 
     
    392396              CALL iom_put( "PFeD"  , zw3d ) 
    393397          ENDIF 
     398          IF( iom_use( "LPRODP" ) )  THEN 
     399              zw3d(:,:,:) = zpligprod1(:,:,:) * 1e9 * zfact * tmask(:,:,:) 
     400              CALL iom_put( "LPRODP"  , zw3d ) 
     401          ENDIF 
     402          IF( iom_use( "LDETP" ) )  THEN 
     403              zw3d(:,:,:) = zpligprod2(:,:,:) * 1e9 * zfact * tmask(:,:,:) 
     404              CALL iom_put( "LDETP"  , zw3d ) 
     405          ENDIF 
    394406          IF( iom_use( "Mumax" ) )  THEN 
    395               zw3d(:,:,:) = prmax(:,:,:) * tmask(:,:,:)   ! Maximum growth rate 
     407              zw3d(:,:,:) = zprmaxn(:,:,:) * tmask(:,:,:)   ! Maximum growth rate 
    396408              CALL iom_put( "Mumax"  , zw3d ) 
    397409          ENDIF 
     
    404416          ENDIF 
    405417          IF( iom_use( "LNlight" ) .OR. iom_use( "LDlight" ) )  THEN 
    406               zw3d(:,:,:) = zprbio (:,:,:) / (prmax(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
     418              zw3d(:,:,:) = zprbio (:,:,:) / (zprmaxn(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
    407419              CALL iom_put( "LNlight"  , zw3d ) 
    408420              ! 
    409               zw3d(:,:,:) =  zprdia (:,:,:) / (prmax(:,:,:) + rtrn) * tmask(:,:,:)  ! light limitation term 
     421              zw3d(:,:,:) = zprdia (:,:,:) / (zprmaxd(:,:,:) + rtrn) * tmask(:,:,:)  ! light limitation term 
    410422              CALL iom_put( "LDlight"  , zw3d ) 
    411423          ENDIF 
     
    542554      !!                     ***  ROUTINE p4z_prod_alloc  *** 
    543555      !!---------------------------------------------------------------------- 
    544       ALLOCATE( prmax(jpi,jpj,jpk), quotan(jpi,jpj,jpk), quotad(jpi,jpj,jpk), STAT = p4z_prod_alloc ) 
     556      ALLOCATE( quotan(jpi,jpj,jpk), quotad(jpi,jpj,jpk), STAT = p4z_prod_alloc ) 
    545557      ! 
    546558      IF( p4z_prod_alloc /= 0 ) CALL ctl_warn('p4z_prod_alloc : failed to allocate arrays.') 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zrem.F90

    r10227 r10362  
    6363      REAL(wp) ::   zsatur, zsatur2, znusil, znusil2, zdep, zdepmin, zfactdep 
    6464      REAL(wp) ::   zbactfer, zolimit, zonitr, zrfact2 
    65       REAL(wp) ::   zammonic, zoxyrem 
     65      REAL(wp) ::   zammonic, zoxyremc, zoxyremn, zoxyremp 
    6666      REAL(wp) ::   zosil, ztem, zdenitnh4, zolimic, zolimin, zolimip, zdenitrn, zdenitrp 
    6767      CHARACTER (len=25) :: charout 
    6868      REAL(wp), DIMENSION(jpi,jpj    ) :: ztempbac 
    69       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdepbac, zolimi, zdepprod, zfacsi, zfacsib 
     69      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdepbac, zolimi, zdepprod, zfacsi, zfacsib, zdepeff, zfebact 
    7070      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
    7171      !!--------------------------------------------------------------------- 
     
    7575      ! Initialisation of arrys 
    7676      zdepprod(:,:,:) = 1._wp 
     77      zdepeff (:,:,:) = 0.3_wp 
    7778      ztempbac(:,:)   = 0._wp 
    7879      zfacsib(:,:,:)  = xsilab / ( 1.0 - xsilab ) 
     80      zfebact(:,:,:)  = 0._wp 
    7981      zfacsi(:,:,:)   = xsilab 
    8082 
     
    9597                  zdepbac (ji,jj,jk) = zdepmin**0.683 * ztempbac(ji,jj) 
    9698                  zdepprod(ji,jj,jk) = zdepmin**0.273 
     99                  zdepeff (ji,jj,jk) = zdepeff(ji,jj,jk) * zdepmin**0.3 
    97100               ENDIF 
    98101            END DO 
     
    117120                  denitr(ji,jj,jk)  = zammonic * ( 1. - nitrfac2(ji,jj,jk) ) 
    118121                  denitr(ji,jj,jk)  = MIN( ( trb(ji,jj,jk,jpno3) - rtrn ) / rdenit, denitr(ji,jj,jk) ) 
    119                   zoxyrem           = zammonic - denitr(ji,jj,jk) 
     122                  zoxyremc          = zammonic - denitr(ji,jj,jk) 
    120123                  ! 
    121124                  zolimi (ji,jj,jk) = MAX( 0.e0, zolimi (ji,jj,jk) ) 
    122125                  denitr (ji,jj,jk) = MAX( 0.e0, denitr (ji,jj,jk) ) 
    123                   zoxyrem           = MAX( 0.e0, zoxyrem ) 
     126                  zoxyremc          = MAX( 0.e0, zoxyremc ) 
    124127 
    125128                  ! 
    126                   tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyrem 
    127                   tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyrem 
     129                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyremc 
     130                  tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyremc 
    128131                  tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) - denitr (ji,jj,jk) * rdenit 
    129                   tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zolimi (ji,jj,jk) - denitr(ji,jj,jk) - zoxyrem 
     132                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zolimi (ji,jj,jk) - denitr(ji,jj,jk) - zoxyremc 
    130133                  tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - zolimi (ji,jj,jk) * o2ut 
    131                   tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyrem 
    132                   tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zolimi(ji,jj,jk) + zoxyrem    & 
     134                  tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyremc 
     135                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zolimi(ji,jj,jk) + zoxyremc    & 
    133136                  &                     + ( rdenit + 1.) * denitr(ji,jj,jk) ) 
    134137               END DO 
     
    159162                  ! Ammonification in suboxic waters with denitrification 
    160163                  ! ------------------------------------------------------- 
    161                   zolimit = zremikc * nitrfac(ji,jj,jk) * trb(ji,jj,jk,jpdoc) 
    162                   denitr(ji,jj,jk)  = MIN(  ( trb(ji,jj,jk,jpno3) - rtrn ) / rdenit, zolimit ) 
    163                   denitr(ji,jj,jk) = MAX( 0.e0, denitr(ji,jj,jk) ) 
     164                  zammonic = zremikc * nitrfac(ji,jj,jk) * trb(ji,jj,jk,jpdoc) 
     165                  denitr(ji,jj,jk)  = zammonic * ( 1. - nitrfac2(ji,jj,jk) ) 
     166                  denitr(ji,jj,jk)  = MAX(0., MIN(  ( trb(ji,jj,jk,jpno3) - rtrn ) / rdenit, denitr(ji,jj,jk) ) ) 
     167                  zoxyremc          = MAX(0., zammonic - denitr(ji,jj,jk)) 
    164168                  zdenitrn  = zremikn * denitr(ji,jj,jk) * trb(ji,jj,jk,jpdon) / ( trb(ji,jj,jk,jpdoc) + rtrn ) 
    165169                  zdenitrp  = zremikp * denitr(ji,jj,jk) * trb(ji,jj,jk,jpdop) / ( trb(ji,jj,jk,jpdoc) + rtrn ) 
    166  
    167                   tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + zolimip + zdenitrp 
    168                   tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zolimin + zdenitrn 
     170                  zoxyremn  = zremikn * zoxyremc * trb(ji,jj,jk,jpdon) / ( trb(ji,jj,jk,jpdoc) + rtrn ) 
     171                  zoxyremp  = zremikp * zoxyremc * trb(ji,jj,jk,jpdop) / ( trb(ji,jj,jk,jpdoc) + rtrn ) 
     172 
     173                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + zolimip + zdenitrp + zoxyremp 
     174                  tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zolimin + zdenitrn + zoxyremn 
    169175                  tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) - denitr(ji,jj,jk) * rdenit 
    170                   tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zolimic - denitr(ji,jj,jk) 
    171                   tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) - zolimin - zdenitrn 
    172                   tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) - zolimip - zdenitrp 
     176                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zolimic - denitr(ji,jj,jk) - zoxyremc 
     177                  tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) - zolimin - zdenitrn - zoxyremn 
     178                  tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) - zolimip - zdenitrp - zoxyremp 
    173179                  tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - zolimic * o2ut 
    174                   tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zolimic + denitr(ji,jj,jk) 
    175                   tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zolimin + ( rdenit + 1.) * zdenitrn ) 
     180                  tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zolimic + denitr(ji,jj,jk) + zoxyremc 
     181                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zolimin + zoxyremn + ( rdenit + 1.) * zdenitrn ) 
    176182               END DO 
    177183            END DO 
     
    215221               ! studies (especially at Papa) have shown this uptake to be significant 
    216222               ! ---------------------------------------------------------- 
    217                zbactfer = feratb *  rfact2 * prmax(ji,jj,jk) * xlimbacl(ji,jj,jk)             & 
     223               zbactfer = feratb *  rfact2 * 0.6_wp / rday * tgfunc(ji,jj,jk) * xlimbacl(ji,jj,jk)     & 
    218224                  &              * trb(ji,jj,jk,jpfer) / ( xkferb + trb(ji,jj,jk,jpfer) )    & 
    219                   &              * zdepprod(ji,jj,jk) * zdepbac(ji,jj,jk) 
    220                tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zbactfer*0.16 
    221                tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zbactfer*0.12 
    222                tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zbactfer*0.04 
     225                  &              * zdepprod(ji,jj,jk) * zdepeff(ji,jj,jk) * zdepbac(ji,jj,jk) 
     226               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zbactfer*0.33 
     227               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zbactfer*0.25 
     228               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zbactfer*0.08 
     229               zfebact(ji,jj,jk)   = zbactfer * 0.33 
     230               blim(ji,jj,jk)      = xlimbacl(ji,jj,jk)  * zdepbac(ji,jj,jk) / 1.e-6 * zdepprod(ji,jj,jk) 
    223231            END DO 
    224232         END DO 
     
    256264               tra(ji,jj,jk,jpgsi) = tra(ji,jj,jk,jpgsi) - zosil 
    257265               tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) + zosil 
    258                ! 
    259266            END DO 
    260267         END DO 
     
    268275 
    269276      IF( knt == nrdttrc ) THEN 
     277          zrfact2 = 1.e3 * rfact2r 
    270278          ALLOCATE( zw3d(jpi,jpj,jpk) ) 
    271279          zfact = 1.e+3 * rfact2r  !  conversion from mol/l/kt to  mol/m3/s 
     
    278286              zw3d(:,:,:) = denitr(:,:,:) * rdenit * rno3 * tmask(:,:,:) * zfact ! Denitrification 
    279287              CALL iom_put( "DENIT"  , zw3d ) 
     288          ENDIF 
     289          IF( iom_use( "BACT" ) )  THEN 
     290               zw3d(:,:,:) = zdepbac(:,:,:) * 1.E6 * tmask(:,:,:)  ! Bacterial biomass 
     291               CALL iom_put( "BACT", zw3d ) 
     292          ENDIF 
     293          IF( iom_use( "FEBACT" ) )  THEN 
     294               zw3d(:,:,:) = zfebact(:,:,:) * 1E9 * tmask(:,:,:) * zrfact2   ! Bacterial iron consumption 
     295               CALL iom_put( "FEBACT" , zw3d ) 
    280296          ENDIF 
    281297          ! 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zsbc.F90

    r10127 r10362  
    3131   REAL(wp), PUBLIC ::   dustsolub    !: Solubility of the dust 
    3232   REAL(wp), PUBLIC ::   mfrac        !: Mineral Content of the dust 
     33   REAL(wp), PUBLIC ::   rdustfep     !: Fraction of dust that is dissolvable 
    3334   REAL(wp), PUBLIC ::   icefeinput   !: Iron concentration in sea ice 
    3435   REAL(wp), PUBLIC ::   wdust        !: Sinking speed of the dust  
     
    134135                  DO ji = 1, jpi 
    135136                     zcoef = ryyss * e1e2t(ji,jj) * h_rnf(ji,jj)  
    136                      rivalk(ji,jj) =   sf_river(jr_dic)%fnow(ji,jj,1)                                    & 
     137                     rivalk(ji,jj) =   sf_river(jr_dic)%fnow(ji,jj,1)  & 
    137138                        &              * 1.E3        / ( 12. * zcoef + rtrn ) 
    138                      rivdic(ji,jj) = ( sf_river(jr_dic)%fnow(ji,jj,1) + sf_river(jr_doc)%fnow(ji,jj,1) ) & 
     139                     rivdic(ji,jj) =   sf_river(jr_dic)%fnow(ji,jj,1) & 
    139140                        &              * 1.E3         / ( 12. * zcoef + rtrn ) 
    140                      rivdin(ji,jj) = ( sf_river(jr_din)%fnow(ji,jj,1) + sf_river(jr_don)%fnow(ji,jj,1) ) & 
     141                     rivdin(ji,jj) =   sf_river(jr_din)%fnow(ji,jj,1) & 
    141142                        &              * 1.E3 / rno3 / ( 14. * zcoef + rtrn ) 
    142                      rivdip(ji,jj) = ( sf_river(jr_dip)%fnow(ji,jj,1) + sf_river(jr_dop)%fnow(ji,jj,1) ) & 
     143                     rivdip(ji,jj) =   sf_river(jr_dip)%fnow(ji,jj,1) & 
    143144                        &              * 1.E3 / po4r / ( 31. * zcoef + rtrn ) 
    144                      rivdsi(ji,jj) =   sf_river(jr_dsi)%fnow(ji,jj,1)                                    & 
     145                     rivdsi(ji,jj) =   sf_river(jr_dsi)%fnow(ji,jj,1)  & 
    145146                        &              * 1.E3        / ( 28.1 * zcoef + rtrn ) 
     147                     rivdoc(ji,jj) =   sf_river(jr_doc)%fnow(ji,jj,1)  & 
     148                        &              * 1.E3        / ( 12. * zcoef + rtrn )  
    146149                  END DO 
    147150               END DO 
     
    158161                     rivdip(ji,jj) = ( sf_river(jr_dip)%fnow(ji,jj,1) ) & 
    159162                        &              * 1.E3 / po4r / ( 31. * zcoef + rtrn ) * tmask(ji,jj,1) 
    160                      rivdoc(ji,jj) = ( sf_river(jr_doc)%fnow(ji,jj,1) ) & 
    161                         &              * 1.E3 / ( 12. * zcoef + rtrn ) * tmask(ji,jj,1) 
    162163                     rivdon(ji,jj) = ( sf_river(jr_don)%fnow(ji,jj,1) ) & 
    163164                        &              * 1.E3 / rno3 / ( 14. * zcoef + rtrn ) * tmask(ji,jj,1) 
    164165                     rivdop(ji,jj) = ( sf_river(jr_dop)%fnow(ji,jj,1) ) & 
    165166                        &              * 1.E3 / po4r / ( 31. * zcoef + rtrn ) * tmask(ji,jj,1) 
     167                     rivdsi(ji,jj) =   sf_river(jr_dsi)%fnow(ji,jj,1)  & 
     168                        &              * 1.E3        / ( 28.1 * zcoef + rtrn ) 
     169                     rivdoc(ji,jj) =   sf_river(jr_doc)%fnow(ji,jj,1)  & 
     170                        &              * 1.E3        / ( 12. * zcoef + rtrn ) 
    166171                  END DO 
    167172               END DO 
     
    223228        &                ln_dust, ln_solub, ln_river, ln_ndepo, ln_ironsed, ln_ironice, ln_hydrofe,    & 
    224229        &                sedfeinput, distcoast, dustsolub, icefeinput, wdust, mfrac, nitrfix, diazolight, concfediaz, & 
    225         &                hratio, fep_rats, fep_rath, lgw_rath 
     230        &                hratio, fep_rats, fep_rath, rdustfep, lgw_rath 
    226231      !!---------------------------------------------------------------------- 
    227232      ! 
     
    262267            WRITE(numout,*) '      Fep/Fer ratio from sed sources            fep_rats   = ', fep_rats 
    263268            WRITE(numout,*) '      Fep/Fer ratio from sed hydro sources      fep_rath   = ', fep_rath 
     269            WRITE(numout,*) '      Fraction of dust that is dissolvable      rdustfep   = ', rdustfep 
    264270            WRITE(numout,*) '      Weak ligand ratio from sed hydro sources  lgw_rath   = ', lgw_rath 
    265271         ENDIF 
     
    343349         slf_river(jr_dsi) = sn_riverdsi   
    344350         ! 
    345          ALLOCATE( rivdic(jpi,jpj), rivalk(jpi,jpj), rivdin(jpi,jpj), rivdip(jpi,jpj), rivdsi(jpi,jpj) )  
    346          IF( ln_p5z )  ALLOCATE( rivdon(jpi,jpj), rivdop(jpi,jpj), rivdoc(jpi,jpj) ) 
     351         ALLOCATE( rivdic(jpi,jpj), rivalk(jpi,jpj), rivdin(jpi,jpj), rivdip(jpi,jpj), rivdsi(jpi,jpj), rivdoc(jpi,jpj) )  
     352         IF( ln_p5z )  ALLOCATE( rivdon(jpi,jpj), rivdop(jpi,jpj) ) 
    347353         ! 
    348354         ALLOCATE( sf_river(jpriv), rivinput(jpriv), STAT=ierr1 )    !* allocate and fill sf_river (forcing structure) with sn_river_ 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zsed.F90

    r10227 r10362  
    6464      CHARACTER (len=25) :: charout 
    6565      REAL(wp), DIMENSION(jpi,jpj    ) :: zdenit2d, zbureff, zwork 
    66       REAL(wp), DIMENSION(jpi,jpj    ) :: zwsbio3, zwsbio4, zwscal 
     66      REAL(wp), DIMENSION(jpi,jpj    ) :: zwsbio3, zwsbio4 
    6767      REAL(wp), DIMENSION(jpi,jpj    ) :: zsedcal, zsedsi, zsedc 
    6868      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zsoufer, zlight 
     
    8585      ! 
    8686      ! Allocate temporary workspace 
    87       IF( ln_p5z )    ALLOCATE( ztrpo4(jpi,jpj,jpk), ztrdop(jpi,jpj,jpk) ) 
     87      ALLOCATE( ztrpo4(jpi,jpj,jpk) ) 
     88      IF( ln_p5z )    ALLOCATE( ztrdop(jpi,jpj,jpk) ) 
    8889      IF( ln_ligand ) ALLOCATE( zwsfep(jpi,jpj) ) 
    89  
    9090 
    9191      zdenit2d(:,:) = 0.e0 
     
    9595      zsedcal (:,:) = 0.e0 
    9696      zsedc   (:,:) = 0.e0 
    97  
    9897 
    9998      ! Iron input/uptake due to sea ice : Crude parameterization based on Lancelot et al. 
     
    132131         ELSE 
    133132            zirondep(:,:,1) = dustsolub  * dust(:,:) * mfrac * rfact2 / e3t_n(:,:,1) / 55.85 + 3.e-10 * r1_ryyss  
     133         ENDIF 
     134         IF ( ln_ligand ) THEN 
     135            IF( ln_solub ) THEN 
     136               tra(:,:,1,jpfep) = tra(:,:,1,jpfep) + rdustfep * (1.0 - solub(:,:)) * dust(:,:) * mfrac * rfact2 / e3t_n(:,:,1) / 55.85 
     137            ELSE 
     138               tra(:,:,1,jpfep) = tra(:,:,1,jpfep) + rdustfep * (1.0 - dustsolub) * dust(:,:) * mfrac * rfact2 / e3t_n(:,:,1) / 55.85 
     139            ENDIF 
    134140         ENDIF 
    135141         zsidep(:,:)   = 8.8 * 0.075 * dust(:,:) * mfrac * rfact2 / e3t_n(:,:,1) / 28.1  
     
    173179                  tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) +  rivdic(ji,jj) * rfact2 
    174180                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) +  ( rivalk(ji,jj) - rno3 * rivdin(ji,jj) ) * rfact2 
     181                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) +  rivdoc(ji,jj) * rfact2 
    175182               ENDDO 
    176183            ENDDO 
    177184         ENDDO 
     185         IF (ln_ligand) THEN 
     186            DO jj = 1, jpj 
     187               DO ji = 1, jpi 
     188                  DO jk = 1, nk_rnf(ji,jj) 
     189                     tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) +  rivdic(ji,jj) * 5.e-5 * rfact2 
     190                  ENDDO 
     191               ENDDO 
     192            ENDDO 
     193         ENDIF 
    178194         IF( ln_p5z ) THEN 
    179195            DO jj = 1, jpj 
     
    182198                     tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) + rivdop(ji,jj) * rfact2 
    183199                     tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) + rivdon(ji,jj) * rfact2 
    184                      tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + rivdoc(ji,jj) * rfact2 
    185200                  ENDDO 
    186201               ENDDO 
     
    216231            zdep = e3t_n(ji,jj,ikt) / xstep 
    217232            zwsbio4(ji,jj) = MIN( 0.99 * zdep, wsbio4(ji,jj,ikt) ) 
    218             zwscal (ji,jj) = MIN( 0.99 * zdep, wscal (ji,jj,ikt) ) 
    219233            zwsbio3(ji,jj) = MIN( 0.99 * zdep, wsbio3(ji,jj,ikt) ) 
    220234         END DO 
     
    278292            ikt  = mbkt(ji,jj) 
    279293            zdep = xstep / e3t_n(ji,jj,ikt)  
    280             zwsc = zwscal (ji,jj) * zdep 
     294            zwsc = zwsbio4(ji,jj) * zdep 
    281295            zsiloss = trb(ji,jj,ikt,jpgsi) * zwsc 
    282296            zcaloss = trb(ji,jj,ikt,jpcal) * zwsc 
     
    292306               ikt  = mbkt(ji,jj) 
    293307               zdep = xstep / e3t_n(ji,jj,ikt)  
    294                zwsc = zwscal (ji,jj) * zdep 
     308               zwsc = zwsbio4(ji,jj) * zdep 
    295309               zsiloss = trb(ji,jj,ikt,jpgsi) * zwsc 
    296310               zcaloss = trb(ji,jj,ikt,jpcal) * zwsc 
     
    393407               DO ji = 1, jpi 
    394408                  !                      ! Potential nitrogen fixation dependant on temperature and iron 
    395                   zlim = ( 1.- xnanono3(ji,jj,jk) - xnanonh4(ji,jj,jk) ) 
    396                   IF( zlim <= 0.2 )   zlim = 0.01 
     409                  ztemp = tsn(ji,jj,jk,jp_tem) 
     410                  zmudia = MAX( 0.,-0.001096*ztemp**2 + 0.057*ztemp -0.637 ) * 7.625 
     411                  !       Potential nitrogen fixation dependant on temperature and iron 
     412                  xdianh4 = trb(ji,jj,jk,jpnh4) / ( concnnh4 + trb(ji,jj,jk,jpnh4) ) 
     413                  xdiano3 = trb(ji,jj,jk,jpno3) / ( concnno3 + trb(ji,jj,jk,jpno3) ) * (1. - xdianh4) 
     414                  zlim = ( 1.- xdiano3 - xdianh4 ) 
     415                  IF( zlim <= 0.1 )   zlim = 0.01 
    397416                  zfact = zlim * rfact2 
    398  
    399                   ztrfer  = biron(ji,jj,jk)       / ( concfediaz + biron(ji,jj,jk)       ) 
    400                   ztrpo4s = trb  (ji,jj,jk,jppo4) / ( concnnh4   + trb  (ji,jj,jk,jppo4) )  
    401                   nitrpot(ji,jj,jk) =  MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) * r1_rday ) & 
    402                     &                *  zfact * MIN( ztrfer, ztrpo4s ) * zlight(ji,jj,jk) 
     417                  ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) ) 
     418                  ztrpo4(ji,jj,jk) = trb(ji,jj,jk,jppo4) / ( 1E-6 + trb(ji,jj,jk,jppo4) ) 
     419                  ztrdp = ztrpo4(ji,jj,jk) 
     420                  nitrpot(ji,jj,jk) =  zmudia * r1_rday * zfact * MIN( ztrfer, ztrdp ) * zlight(ji,jj,jk) 
    403421               END DO 
    404422            END DO 
     
    434452               DO ji = 1, jpi 
    435453                  zfact = nitrpot(ji,jj,jk) * nitrfix 
    436                   tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) +             zfact 
    437                   tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3      * zfact 
    438                   tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + o2nit     * zfact  
     454                  tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zfact / 3.0 
     455                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * zfact / 3.0 
     456                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) - zfact * 2.0 / 3.0 
     457                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zfact * 1.0 / 3.0 
     458                  tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zfact * 1.0 / 3.0 * 2.0 / 3.0 
     459                  tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zfact * 1.0 / 3.0 * 1.0 / 3.0 
     460                  tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + ( o2ut + o2nit ) * zfact * 2.0 / 3.0 + o2nit * zfact / 3.0 
     461                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - 30E-6 * zfact * 1.0 / 3.0 
     462                  tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + 30E-6 * zfact * 1.0 / 3.0 * 2.0 / 3.0 
     463                  tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + 30E-6 * zfact * 1.0 / 3.0 * 1.0 / 3.0 
     464                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * rfact2 / rday 
    439465                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + concdnh4 / ( concdnh4 + trb(ji,jj,jk,jppo4) ) & 
    440                   &                     * 0.002 * trb(ji,jj,jk,jpdoc) * xstep 
    441                   tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * xstep 
     466                  &                     * 0.001 * trb(ji,jj,jk,jpdoc) * xstep 
    442467              END DO 
    443468            END DO  
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zsink.F90

    r10069 r10362  
    143143      END DO 
    144144 
    145       wscal (:,:,:) = wsbio4(:,:,:) 
    146  
    147145      !  Initializa to zero all the sinking arrays  
    148146      !   ----------------------------------------- 
     
    165163        CALL p4z_sink2( wsbio4, sinkfer2, jpbfe, iiter2 ) 
    166164        CALL p4z_sink2( wsbio4, sinksil , jpgsi, iiter2 ) 
    167         CALL p4z_sink2( wscal , sinkcal , jpcal, iiter2 ) 
     165        CALL p4z_sink2( wsbio4, sinkcal , jpcal, iiter2 ) 
    168166      END DO 
    169167 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p5zlim.F90

    r10227 r10362  
    1414   USE oce_trc         ! Shared ocean-passive tracers variables 
    1515   USE trc             ! Tracers defined 
     16   USE p4zlim 
    1617   USE sms_pisces      ! PISCES variables 
    1718   USE iom             !  I/O manager 
     
    2526 
    2627   !! * Shared module variables 
    27    REAL(wp), PUBLIC ::  concnno3    !:  NO3, PO4 half saturation    
    2828   REAL(wp), PUBLIC ::  concpno3    !:  NO3, PO4 half saturation    
    29    REAL(wp), PUBLIC ::  concdno3    !:  Phosphate half saturation for diatoms   
    30    REAL(wp), PUBLIC ::  concnnh4    !:  NH4 half saturation for phyto   
    3129   REAL(wp), PUBLIC ::  concpnh4    !:  NH4 half saturation for phyto   
    32    REAL(wp), PUBLIC ::  concdnh4    !:  NH4 half saturation for diatoms 
    3330   REAL(wp), PUBLIC ::  concnpo4    !:  NH4 half saturation for diatoms 
    3431   REAL(wp), PUBLIC ::  concppo4    !:  NH4 half saturation for diatoms 
    3532   REAL(wp), PUBLIC ::  concdpo4    !:  NH4 half saturation for diatoms 
    36    REAL(wp), PUBLIC ::  concnfer    !:  Iron half saturation for nanophyto  
    3733   REAL(wp), PUBLIC ::  concpfer    !:  Iron half saturation for nanophyto  
    38    REAL(wp), PUBLIC ::  concdfer    !:  Iron half saturation for diatoms   
    39    REAL(wp), PUBLIC ::  concbno3    !:  NO3 half saturation  for bacteria  
    40    REAL(wp), PUBLIC ::  concbnh4    !:  NH4 half saturation for bacteria 
    4134   REAL(wp), PUBLIC ::  concbpo4    !:  PO4 half saturation for bacteria 
    42    REAL(wp), PUBLIC ::  xsizedia    !:  Minimum size criteria for diatoms 
    4335   REAL(wp), PUBLIC ::  xsizepic    !:  Minimum size criteria for diatoms 
    44    REAL(wp), PUBLIC ::  xsizephy    !:  Minimum size criteria for nanophyto 
    45    REAL(wp), PUBLIC ::  xsizern     !:  Size ratio for nanophytoplankton 
    4636   REAL(wp), PUBLIC ::  xsizerp     !:  Size ratio for nanophytoplankton 
    47    REAL(wp), PUBLIC ::  xsizerd     !:  Size ratio for diatoms 
    48    REAL(wp), PUBLIC ::  xksi1       !:  half saturation constant for Si uptake  
    49    REAL(wp), PUBLIC ::  xksi2       !:  half saturation constant for Si/C  
    50    REAL(wp), PUBLIC ::  xkdoc       !:  2nd half-sat. of DOC remineralization   
    51    REAL(wp), PUBLIC ::  concbfe     !:  Fe half saturation for bacteria  
    52    REAL(wp), PUBLIC ::  oxymin      !:  half saturation constant for anoxia 
    5337   REAL(wp), PUBLIC ::  qfnopt      !:  optimal Fe quota for nanophyto 
    5438   REAL(wp), PUBLIC ::  qfpopt      !:  optimal Fe quota for nanophyto 
    5539   REAL(wp), PUBLIC ::  qfdopt      !:  optimal Fe quota for diatoms 
    56    REAL(wp), PUBLIC ::  caco3r      !:  mean rainratio  
    5740   REAL(wp), PUBLIC ::  qnnmin      !:  optimal Fe quota for diatoms 
    5841   REAL(wp), PUBLIC ::  qnnmax      !:  optimal Fe quota for diatoms 
     
    8972 
    9073   !!* Phytoplankton limitation terms 
    91    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xnanono3   !: ??? 
    92    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xdiatno3   !: ??? 
    93    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xnanonh4   !: ??? 
    94    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xdiatnh4   !: ??? 
    95    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xnanopo4   !: ??? 
    96    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xdiatpo4   !: ??? 
    97    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimphy    !: ??? 
    98    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimdia    !: ??? 
    99    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimnfe    !: ??? 
    100    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimdfe    !: ??? 
    101    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimsi     !: ??? 
    102    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimbac    !: ?? 
    103    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimbacl   !: ?? 
    10474   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xpicono3   !: ??? 
    10575   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xpiconh4   !: ??? 
     
    574544         &      xpicopo4(jpi,jpj,jpk), xpicodop(jpi,jpj,jpk),       & 
    575545         &      xnanodop(jpi,jpj,jpk), xdiatdop(jpi,jpj,jpk),       & 
    576          &      xnanono3(jpi,jpj,jpk), xdiatno3(jpi,jpj,jpk),       & 
    577          &      xnanonh4(jpi,jpj,jpk), xdiatnh4(jpi,jpj,jpk),       & 
    578          &      xnanopo4(jpi,jpj,jpk), xdiatpo4(jpi,jpj,jpk),       & 
    579          &      xlimphy (jpi,jpj,jpk), xlimdia (jpi,jpj,jpk),       & 
    580          &      xlimnfe (jpi,jpj,jpk), xlimdfe (jpi,jpj,jpk),       & 
    581          &      xlimbac (jpi,jpj,jpk), xlimbacl(jpi,jpj,jpk),       & 
    582546         &      xnanofer(jpi,jpj,jpk), xdiatfer(jpi,jpj,jpk),       & 
    583547         &      xpicofer(jpi,jpj,jpk), xlimpfe (jpi,jpj,jpk),       & 
    584548         &      fvnuptk (jpi,jpj,jpk), fvduptk (jpi,jpj,jpk),       & 
    585          &      fvpuptk (jpi,jpj,jpk), xlimpic (jpi,jpj,jpk),       & 
    586          &      xlimsi  (jpi,jpj,jpk), STAT=ierr(1) ) 
     549         &      fvpuptk (jpi,jpj,jpk), xlimpic (jpi,jpj,jpk),    STAT=ierr(1) ) 
    587550         ! 
    588551      !*  Minimum/maximum quotas of phytoplankton 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p5zmeso.F90

    r10070 r10362  
    2727   REAL(wp), PUBLIC ::  part2        !: part of calcite not dissolved in mesozoo guts 
    2828   REAL(wp), PUBLIC ::  xpref2c      !: mesozoo preference for POC  
    29    REAL(wp), PUBLIC ::  xpref2p      !: mesozoo preference for nanophyto 
     29   REAL(wp), PUBLIC ::  xpref2n      !: mesozoo preference for nanophyto 
    3030   REAL(wp), PUBLIC ::  xpref2z      !: mesozoo preference for zooplankton 
    3131   REAL(wp), PUBLIC ::  xpref2d      !: mesozoo preference for Diatoms  
     
    4040   REAL(wp), PUBLIC ::  mzrat2       !: microzooplankton mortality rate  
    4141   REAL(wp), PUBLIC ::  grazrat2     !: maximal mesozoo grazing rate 
    42    REAL(wp), PUBLIC ::  xkgraz2      !: non assimilated fraction of P by mesozoo 
    43    REAL(wp), PUBLIC ::  unass2c      !: Efficicency of mesozoo growth  
    44    REAL(wp), PUBLIC ::  unass2n      !: Efficicency of mesozoo growth  
    45    REAL(wp), PUBLIC ::  unass2p      !: Efficicency of mesozoo growth  
    46    REAL(wp), PUBLIC ::  epsher2      !: half sturation constant for grazing 2 
     42   REAL(wp), PUBLIC ::  xkgraz2      !: Half-saturation constant of assimilation 
     43   REAL(wp), PUBLIC ::  unass2c      !: Non-assimilated fraction of food 
     44   REAL(wp), PUBLIC ::  unass2n      !: Non-assimilated fraction of food 
     45   REAL(wp), PUBLIC ::  unass2p      !: Non-assimilated fraction of food 
     46   REAL(wp), PUBLIC ::  epsher2      !: Growth efficiency of mesozoo 
     47   REAL(wp), PUBLIC ::  epsher2min   !: Minimum growth efficiency of mesozoo 
    4748   REAL(wp), PUBLIC ::  ssigma2      !: Fraction excreted as semi-labile DOM 
    4849   REAL(wp), PUBLIC ::  srespir2     !: Active respiration 
     
    8586      CHARACTER (len=25) :: charout 
    8687      REAL(wp) :: zrfact2, zmetexcess 
    87       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing 
    88       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
     88      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing, zfezoo2 
     89      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d, zz2ligprod 
    8990 
    9091      !!--------------------------------------------------------------------- 
     
    9293      IF( ln_timing )   CALL timing_start('p5z_meso') 
    9394      ! 
     95 
    9496      zgrazing(:,:,:) = 0._wp 
     97      zfezoo2 (:,:,:) = 0._wp 
     98      ! 
     99      IF (ln_ligand) THEN 
     100         ALLOCATE( zz2ligprod(jpi,jpj,jpk) ) 
     101         zz2ligprod(:,:,:) = 0._wp 
     102      ENDIF 
    95103 
    96104      zmetexcess = 0.0 
     
    111119               !   no real reason except that it seems to be more stable and may mimic predation 
    112120               !   --------------------------------------------------------------- 
    113                ztortz   = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes) 
     121               ztortz   = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes) * (1. - nitrfac(ji,jj,jk)) 
    114122 
    115123               !   Computation of the abundance of the preys 
     
    124132               !   Mesozooplankton grazing 
    125133               !   ------------------------ 
    126                zfood     = xpref2d * zcompadi + xpref2z * zcompaz + xpref2p * zcompaph + xpref2c * zcompapoc   & 
     134               zfood     = xpref2d * zcompadi + xpref2z * zcompaz + xpref2n * zcompaph + xpref2c * zcompapoc   & 
    127135               &           + xpref2m * zcompames  
    128136               zfoodlim  = MAX( 0., zfood - MIN( 0.5 * zfood, xthresh2 ) ) 
    129137               zdenom    = zfoodlim / ( xkgraz2 + zfoodlim ) 
    130                zgraze2   = grazrat2 * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpmes)  
     138               zgraze2   = grazrat2 * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpmes) * (1. - nitrfac(ji,jj,jk))  
    131139 
    132140               !   An active switching parameterization is used here. 
     
    138146               !   most abundant species 
    139147               !   ------------------------------------------------------------   
    140                ztmp1 = xpref2p * zcompaph**1.5 
     148               ztmp1 = xpref2n * zcompaph**1.5 
    141149               ztmp2 = xpref2m * zcompames**1.5 
    142150               ztmp3 = xpref2c * zcompapoc**1.5 
     
    170178               !   ---------------------------------- 
    171179               zgrazffeg = grazflux  * xstep * wsbio4(ji,jj,jk)      & 
    172                &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpgoc) * trb(ji,jj,jk,jpmes) 
     180               &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpgoc) * trb(ji,jj,jk,jpmes)  & 
     181               &           * (1. - nitrfac(ji,jj,jk)) 
    173182               zgrazfffg = zgrazffeg * trb(ji,jj,jk,jpbfe) / (trb(ji,jj,jk,jpgoc) + rtrn) 
    174183               zgrazffng = zgrazffeg * trb(ji,jj,jk,jpgon) / (trb(ji,jj,jk,jpgoc) + rtrn) 
    175184               zgrazffpg = zgrazffeg * trb(ji,jj,jk,jpgop) / (trb(ji,jj,jk,jpgoc) + rtrn) 
    176185               zgrazffep = grazflux  * xstep *  wsbio3(ji,jj,jk)     & 
    177                &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpmes) 
     186               &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpmes)   & 
     187               &           * (1. - nitrfac(ji,jj,jk)) 
    178188               zgrazfffp = zgrazffep * trb(ji,jj,jk,jpsfe) / (trb(ji,jj,jk,jppoc) + rtrn) 
    179189               zgrazffnp = zgrazffep * trb(ji,jj,jk,jppon) / (trb(ji,jj,jk,jppoc) + rtrn) 
     
    226236               !   --------------------------------------------------- 
    227237               zepshert  = MIN( 1., zgrasratn/ no3rat3, zgrasratp/ po4rat3, zgrasratf / ferat3) 
    228                zbeta = 1./ (epsher2 - 0.2) 
    229                zepsherf = 0.2 + 1./ (zbeta + 0.04 * 12. * zfood *1E6 ) 
     238               zbeta     = MAX(0., (epsher2 - epsher2min) ) 
     239               zepsherf  = epsher2min + zbeta / ( 1.0 + 0.04E6 * 12. * zfood * zbeta ) 
    230240               zepsherv  = zepsherf * zepshert 
    231241 
     
    290300               tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zgradoc 
    291301               ! 
    292                IF( ln_ligand ) tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zgradoc * ldocz 
     302               IF( ln_ligand ) THEN 
     303                  tra(ji,jj,jk,jplgw)  = tra(ji,jj,jk,jplgw) + zgradoc * ldocz 
     304                  zz2ligprod(ji,jj,jk) = zgradoc * ldocz 
     305               ENDIF 
    293306               ! 
    294307               tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) + zgradon 
     
    296309               tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - o2ut * zgrarem 
    297310               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zgraref 
     311               zfezoo2(ji,jj,jk)   = zgraref 
    298312               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrarem 
    299313               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * zgraren 
     
    351365            CALL iom_put( "PCAL", zw3d ) 
    352366         ENDIF 
     367         IF( iom_use( "FEZOO2" ) ) THEN 
     368            zw3d(:,:,:) = zfezoo2(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)   ! 
     369            CALL iom_put( "FEZOO2", zw3d ) 
     370         ENDIF 
     371         IF( iom_use( "LPRODZ2" ) .AND. ln_ligand )  THEN 
     372            zw3d(:,:,:) = zz2ligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) 
     373            CALL iom_put( "LPRODZ2"  , zw3d ) 
     374         ENDIF 
    353375         DEALLOCATE( zw3d ) 
    354376      ENDIF 
     
    379401      INTEGER :: ios                 ! Local integer output status for namelist read 
    380402      !! 
    381       NAMELIST/namp5zmes/part2, bmetexc2, grazrat2, resrat2, mzrat2, xpref2c, xpref2p, xpref2z, & 
     403      NAMELIST/namp5zmes/part2, bmetexc2, grazrat2, resrat2, mzrat2, xpref2c, xpref2n, xpref2z, & 
    382404         &                xpref2m, xpref2d, xthresh2dia, xthresh2phy, xthresh2zoo, xthresh2poc, & 
    383          &                xthresh2mes, xthresh2, xkgraz2, epsher2, ssigma2, unass2c, & 
     405         &                xthresh2mes, xthresh2, xkgraz2, epsher2, epsher2min, ssigma2, unass2c, & 
    384406         &                unass2n, unass2p, srespir2, grazflux 
    385407      !!---------------------------------------------------------------------- 
     
    399421         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~' 
    400422         WRITE(numout,*) '    part of calcite not dissolved in mesozoo guts  part2       = ', part2 
    401          WRITE(numout,*) '    mesozoo preference for nano.                   xpref2p     = ', xpref2p 
     423         WRITE(numout,*) '    mesozoo preference for nano.                   xpref2n     = ', xpref2n 
    402424         WRITE(numout,*) '    mesozoo preference for diatoms                 xpref2d     = ', xpref2d 
    403425         WRITE(numout,*) '    mesozoo preference for zoo                     xpref2z     = ', xpref2z 
     
    415437         WRITE(numout,*) '    mesozoo flux feeding rate                      grazflux    = ', grazflux 
    416438         WRITE(numout,*) '    C egested fraction of food by mesozoo          unass2c     = ', unass2c 
    417          WRITE(numout,*) '    C egested fraction of food by mesozoo          unass2n     = ', unass2n 
    418          WRITE(numout,*) '    C egested fraction of food by mesozoo          unass2p     = ', unass2p 
     439         WRITE(numout,*) '    N egested fraction of food by mesozoo          unass2n     = ', unass2n 
     440         WRITE(numout,*) '    P egested fraction of food by mesozoo          unass2p     = ', unass2p 
    419441         WRITE(numout,*) '    Efficicency of Mesozoo growth                  epsher2     = ', epsher2 
     442         WRITE(numout,*) '    Minimum Efficiency of Mesozoo growth           epsher2min  =', epsher2min 
    420443         WRITE(numout,*) '    Fraction excreted as semi-labile DOM           ssigma2     = ', ssigma2 
    421444         WRITE(numout,*) '    Active respiration                             srespir2    = ', srespir2 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p5zmicro.F90

    r10227 r10362  
    1515   USE trc             !  passive tracers common variables  
    1616   USE sms_pisces      !  PISCES Source Minus Sink variables 
     17   USE p4zlim 
    1718   USE p5zlim          !  Phytoplankton limitation terms 
    1819   USE iom             !  I/O manager 
     
    4142   REAL(wp), PUBLIC ::  mzrat       !: microzooplankton mortality rate  
    4243   REAL(wp), PUBLIC ::  grazrat     !: maximal microzoo grazing rate 
    43    REAL(wp), PUBLIC ::  xkgraz      !: non assimilated fraction of P by microzoo 
    44    REAL(wp), PUBLIC ::  unassc      !: Efficicency of microzoo growth  
    45    REAL(wp), PUBLIC ::  unassn      !: Efficicency of microzoo growth  
    46    REAL(wp), PUBLIC ::  unassp      !: Efficicency of microzoo growth  
    47    REAL(wp), PUBLIC ::  epsher      !: half sturation constant for grazing 1  
     44   REAL(wp), PUBLIC ::  xkgraz      !: Half-saturation constant of assimilation 
     45   REAL(wp), PUBLIC ::  unassc      !: Non-assimilated part of food 
     46   REAL(wp), PUBLIC ::  unassn      !: Non-assimilated part of food 
     47   REAL(wp), PUBLIC ::  unassp      !: Non-assimilated part of food 
     48   REAL(wp), PUBLIC ::  epsher      !: Growth efficiency for microzoo 
     49   REAL(wp), PUBLIC ::  epshermin   !: Minimum growth efficiency for microzoo 
    4850   REAL(wp), PUBLIC ::  srespir     !: half sturation constant for grazing 1  
    4951   REAL(wp), PUBLIC ::  ssigma      !: Fraction excreted as semi-labile DOM 
     
    8284      REAL(wp) :: zgrazdc, zgrazdn, zgrazdp, zgrazdf, zgraznf, zgrazz 
    8385      REAL(wp) :: zgrazpc, zgrazpn, zgrazpp, zgrazpf, zbeta, zrfact2, zmetexcess 
    84       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing 
    85       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
     86      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing, zfezoo 
     87      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d, zzligprod 
    8688      CHARACTER (len=25) :: charout 
    8789      !!--------------------------------------------------------------------- 
    8890      ! 
    8991      IF( ln_timing )   CALL timing_start('p5z_micro') 
     92      ! 
     93      IF (ln_ligand) THEN 
     94         ALLOCATE( zzligprod(jpi,jpj,jpk) ) 
     95         zzligprod(:,:,:) = 0._wp 
     96      ENDIF 
    9097      ! 
    9198      zmetexcess = 0.0 
     
    106113               !   no real reason except that it seems to be more stable and may mimic predation. 
    107114               !   ------------------------------------------------------------------------------ 
    108                ztortz = mzrat * 1.e6 * zfact * trb(ji,jj,jk,jpzoo) 
     115               ztortz = mzrat * 1.e6 * zfact * trb(ji,jj,jk,jpzoo) * (1. - nitrfac(ji,jj,jk)) 
    109116 
    110117               !   Computation of the abundance of the preys 
     
    123130               zfoodlim  = MAX( 0. , zfood - min(xthresh,0.5*zfood) ) 
    124131               zdenom    = zfoodlim / ( xkgraz + zfoodlim ) 
    125                zgraze    = grazrat * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpzoo)  
     132               zgraze    = grazrat * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpzoo) * (1. - nitrfac(ji,jj,jk))  
    126133 
    127134               !   An active switching parameterization is used here. 
     
    183190               !   --------------------------------------------------- 
    184191               zepshert  = MIN( 1., zgrasratn/ no3rat3, zgrasratp/ po4rat3, zgrasratf / ferat3) 
    185                zbeta = 1./ (epsher - 0.2) 
    186                zepsherf = 0.2 + 1./ (zbeta + 0.04 * 12. * zfood * 1E6 ) 
     192               zbeta     = MAX( 0., (epsher - epshermin) ) 
     193               zepsherf  = epshermin + zbeta / ( 1.0 + 0.04E6 * 12. * zfood * zbeta ) 
    187194               zepsherv  = zepsherf * zepshert 
    188195 
     
    244251               tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zgradoc 
    245252               ! 
    246                IF( ln_ligand ) tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zgradoc * ldocz 
     253               IF( ln_ligand ) THEN  
     254                  tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zgradoc * ldocz 
     255                  zzligprod(ji,jj,jk) = zgradoc * ldocz 
     256               ENDIF 
    247257               ! 
    248258               tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) + zgradon 
     
    250260               tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - o2ut * zgrarem  
    251261               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zgraref 
     262               zfezoo(ji,jj,jk)    = zgraref 
    252263               tra(ji,jj,jk,jpzoo) = tra(ji,jj,jk,jpzoo) + zepsherv * zgraztotc - zrespirc - ztortz - zgrazz 
    253264               tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) - zgraznc 
     
    288299      END DO 
    289300      ! 
    290       IF( lk_iomput .AND. knt == nrdttrc ) THEN 
    291          ALLOCATE( zw3d(jpi,jpj,jpk) ) 
    292          IF( iom_use( "GRAZ1" ) ) THEN 
    293             zw3d(:,:,:) = zgrazing(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:)  !  Total grazing of phyto by zooplankton 
    294             CALL iom_put( "GRAZ1", zw3d ) 
     301      IF( lk_iomput ) THEN 
     302         IF( knt == nrdttrc ) THEN 
     303            ALLOCATE( zw3d(jpi,jpj,jpk) ) 
     304            IF( iom_use( "GRAZ1" ) ) THEN 
     305               zw3d(:,:,:) = zgrazing(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:)  !  Total grazing of phyto by zooplankton 
     306               CALL iom_put( "GRAZ1", zw3d ) 
     307            ENDIF 
     308            IF( iom_use( "FEZOO" ) ) THEN 
     309               zw3d(:,:,:) = zfezoo(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)   ! 
     310               CALL iom_put( "FEZOO", zw3d ) 
     311            ENDIF 
     312            IF( iom_use( "LPRODZ" ) .AND. ln_ligand )  THEN 
     313               zw3d(:,:,:) = zzligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) 
     314               CALL iom_put( "LPRODZ"  , zw3d ) 
     315            ENDIF 
     316            DEALLOCATE( zw3d ) 
    295317         ENDIF 
    296          DEALLOCATE( zw3d ) 
    297318      ENDIF 
    298319      ! 
     
    325346         &                xprefp, xprefd, xprefz, xthreshdia, xthreshphy, & 
    326347         &                xthreshpic, xthreshpoc, xthreshzoo, xthresh, xkgraz, & 
    327          &                epsher, ssigma, srespir, unassc, unassn, unassp 
     348         &                epsher, epshermin, ssigma, srespir, unassc, unassn, unassp 
    328349      !!---------------------------------------------------------------------- 
    329350      ! 
     
    360381         WRITE(numout,*) '    P egested fraction of fodd by microzoo          unassp      =', unassp 
    361382         WRITE(numout,*) '    Efficicency of microzoo growth                  epsher      =', epsher 
     383         WRITE(numout,*) '    Minimum Efficiency of Microzoo growth           epshermin   =', epshermin 
    362384         WRITE(numout,*) '    Fraction excreted as semi-labile DOM            ssigma      =', ssigma 
    363385         WRITE(numout,*) '    Active respiration                              srespir     =', srespir 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p5zmort.F90

    r10227 r10362  
    1414   USE trc             !  passive tracers common variables  
    1515   USE sms_pisces      !  PISCES Source Minus Sink variables 
     16   USE p4zlim 
    1617   USE p5zlim          !  Phytoplankton limitation terms 
    1718   USE prtctl_trc      !  print control for debugging 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p5zprod.F90

    r10227 r10362  
    1616   USE trc             !  passive tracers common variables  
    1717   USE sms_pisces      !  PISCES Source Minus Sink variables 
     18   USE p4zlim 
    1819   USE p5zlim          !  Co-limitations of differents nutrients 
    1920   USE prtctl_trc      !  print control for debugging 
     
    4243   REAL(wp), PUBLIC ::  grosip          !: 
    4344 
    44    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmaxn   !: optimal production = f(temperature) 
    45    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmaxp   !: optimal production = f(temperature) 
    46    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmaxd   !: optimal production = f(temperature) 
    4745   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   zdaylen 
    4846    
     
    8381      REAL(wp), DIMENSION(jpi,jpj    ) :: zmixnano, zmixpico, zmixdiat, zstrn 
    8482      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpislopeadn, zpislopeadp, zpislopeadd 
     83      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprnut, zprmaxp, zprmaxn, zprmaxd 
    8584      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprbio, zprpic, zprdia, zysopt 
    8685      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprchln, zprchlp, zprchld 
     
    9190      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpropo4n, zpropo4p, zpropo4d 
    9291      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprodopn, zprodopp, zprodopd 
    93       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zrespn, zrespp, zrespd, zprnut 
     92      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zrespn, zrespp, zrespd 
    9493      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zcroissn, zcroissp, zcroissd 
    9594      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmxl_fac, zmxl_chl 
     95      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpligprod1, zpligprod2 
    9696      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
    9797      REAL(wp), ALLOCATABLE, DIMENSION(:,:  ) :: zw2d 
     
    106106      zpropo4n(:,:,:) = 0._wp ; zpropo4p(:,:,:) = 0._wp ; zpropo4d(:,:,:) = 0._wp 
    107107      zprdia  (:,:,:) = 0._wp ; zprpic  (:,:,:) = 0._wp ; zprbio  (:,:,:) = 0._wp 
     108      zprodopn(:,:,:) = 0._wp ; zprodopp(:,:,:) = 0._wp ; zprodopd(:,:,:) = 0._wp 
    108109      zysopt  (:,:,:) = 0._wp 
    109110      zrespn  (:,:,:) = 0._wp ; zrespp  (:,:,:) = 0._wp ; zrespd  (:,:,:) = 0._wp  
    110111 
    111112      ! Computation of the optimal production 
    112       prmaxn(:,:,:) = ( 0.65_wp * (1. + zpsino3 * qnpmax ) ) * r1_rday * tgfunc(:,:,:) 
    113       prmaxp(:,:,:) = 0.5 / 0.65 * prmaxn(:,:,:)  
    114       prmaxd(:,:,:) = prmaxn(:,:,:)  
    115       zprnut(:,:,:) = 0.65_wp * r1_rday * tgfunc(:,:,:) 
     113      zprnut (:,:,:) = 0.65_wp * r1_rday * tgfunc(:,:,:) 
     114      zprmaxn(:,:,:) = ( 0.65_wp * (1. + zpsino3 * qnpmax ) ) * r1_rday * tgfunc(:,:,:) 
     115      zprmaxp(:,:,:) = 0.5 / 0.65 * zprmaxn(:,:,:)  
     116      zprmaxd(:,:,:) = zprmaxn(:,:,:)  
    116117 
    117118      ! compute the day length depending on latitude and the day 
     
    145146      END DO 
    146147 
    147       zprbio(:,:,:) = prmaxn(:,:,:) * zmxl_fac(:,:,:) 
    148       zprdia(:,:,:) = prmaxd(:,:,:) * zmxl_fac(:,:,:) 
    149       zprpic(:,:,:) = prmaxp(:,:,:) * zmxl_fac(:,:,:) 
     148      zprbio(:,:,:) = zprmaxn(:,:,:) * zmxl_fac(:,:,:) 
     149      zprdia(:,:,:) = zprmaxd(:,:,:) * zmxl_fac(:,:,:) 
     150      zprpic(:,:,:) = zprmaxp(:,:,:) * zmxl_fac(:,:,:) 
    150151 
    151152 
     
    184185                  zpisloped = zpisloped * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    185186                  zpislopep = zpislopep * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    186                   zprchln(ji,jj,jk) = prmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) )  ) 
    187                   zprchlp(ji,jj,jk) = prmaxp(ji,jj,jk) * ( 1.- EXP( -zpislopep * epico(ji,jj,jk) )  ) 
    188                   zprchld(ji,jj,jk) = prmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) )  ) 
     187                  zprchln(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) )  ) 
     188                  zprchlp(ji,jj,jk) = zprmaxp(ji,jj,jk) * ( 1.- EXP( -zpislopep * epicom(ji,jj,jk) )  ) 
     189                  zprchld(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) )  ) 
    189190               ENDIF 
    190191            END DO 
     
    203204                  !    to mimic the very high ratios observed in the Southern Ocean (silpot2) 
    204205                  zlim  = trb(ji,jj,jk,jpsil) / ( trb(ji,jj,jk,jpsil) + xksi1 ) 
    205                   zsilim = MIN( zprdia(ji,jj,jk) / ( prmaxd(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) ) 
     206                  zsilim = MIN( zprdia(ji,jj,jk) / ( zprmaxd(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) ) 
    206207                  zsilfac = 3.4 * EXP( -4.23 * zsilim ) * MAX( 0.e0, MIN( 1., 2.2 * ( zlim - 0.5 ) )  ) + 1.e0 
    207208                  zsiborn = trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) 
     
    347348               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
    348349                     !  production terms for nanophyto. ( chlorophyll ) 
    349                   znanotot = enano(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
     350                  znanotot = enanom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    350351                  zprod = rday * (zpronewn(ji,jj,jk) + zproregn(ji,jj,jk)) * zprchln(ji,jj,jk) * xlimphy(ji,jj,jk) 
    351352                  thetannm_n   = MIN ( thetannm, ( thetannm / (1. - 1.14 / 43.4 *tsn(ji,jj,jk,jp_tem)))   & 
     
    354355                  zprochln = MAX(zprochln, chlcmin * 12. * zprorcan (ji,jj,jk) ) 
    355356                     !  production terms for picophyto. ( chlorophyll ) 
    356                   zpicotot = epico(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
     357                  zpicotot = epicom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    357358                  zprod = rday * (zpronewp(ji,jj,jk) + zproregp(ji,jj,jk)) * zprchlp(ji,jj,jk) * xlimpic(ji,jj,jk) 
    358359                  thetanpm_n   = MIN ( thetanpm, ( thetanpm / (1. - 1.14 / 43.4 *tsn(ji,jj,jk,jp_tem)))   & 
     
    361362                  zprochlp = MAX(zprochlp, chlcmin * 12. * zprorcap(ji,jj,jk) ) 
    362363                  !  production terms for diatomees ( chlorophyll ) 
    363                   zdiattot = ediat(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
     364                  zdiattot = ediatm(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    364365                  zprod = rday * (zpronewd(ji,jj,jk) + zproregd(ji,jj,jk)) * zprchld(ji,jj,jk) * xlimdia(ji,jj,jk) 
    365366                  thetandm_n   = MIN ( thetandm, ( thetandm / (1. - 1.14 / 43.4 *tsn(ji,jj,jk,jp_tem)))   & 
     
    449450                 zfeup    = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk) + texcretp * zprofep(ji,jj,jk) 
    450451                 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 
    451454              END DO 
    452455           END DO 
     
    500503              CALL iom_put( "PFeD"  , zw3d ) 
    501504          ENDIF 
     505          IF( iom_use( "LPRODP" ) )  THEN 
     506              zw3d(:,:,:) = zpligprod1(:,:,:) * 1e9 * zfact * tmask(:,:,:) 
     507              CALL iom_put( "LPRODP"  , zw3d ) 
     508          ENDIF 
     509          IF( iom_use( "LDETP" ) )  THEN 
     510              zw3d(:,:,:) = zpligprod2(:,:,:) * 1e9 * zfact * tmask(:,:,:) 
     511              CALL iom_put( "LDETP"  , zw3d ) 
     512          ENDIF 
    502513          IF( iom_use( "Mumax" ) )  THEN 
    503               zw3d(:,:,:) = prmaxn(:,:,:) * tmask(:,:,:)   ! Maximum growth rate 
     514              zw3d(:,:,:) = zprmaxn(:,:,:) * tmask(:,:,:)   ! Maximum growth rate 
    504515              CALL iom_put( "Mumax"  , zw3d ) 
    505516          ENDIF 
     
    515526          ENDIF 
    516527          IF( iom_use( "LNlight" ) .OR. iom_use( "LDlight" ) .OR. iom_use( "LPlight" ) )  THEN 
    517               zw3d(:,:,:) = zprbio (:,:,:) / (prmaxn(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
     528              zw3d(:,:,:) = zprbio (:,:,:) / (zprmaxn(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
    518529              CALL iom_put( "LNlight"  , zw3d ) 
    519530              ! 
    520               zw3d(:,:,:) = zprpic (:,:,:) / (prmaxp(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
     531              zw3d(:,:,:) = zprpic (:,:,:) / (zprmaxp(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
    521532              CALL iom_put( "LPlight"  , zw3d ) 
    522533              ! 
    523               zw3d(:,:,:) =  zprdia (:,:,:) / (prmaxd(:,:,:) + rtrn) * tmask(:,:,:)  ! light limitation term 
     534              zw3d(:,:,:) =  zprdia (:,:,:) / (zprmaxd(:,:,:) + rtrn) * tmask(:,:,:)  ! light limitation term 
    524535              CALL iom_put( "LDlight"  , zw3d ) 
    525536          ENDIF 
     
    611622      !!                     ***  ROUTINE p5z_prod_alloc  *** 
    612623      !!---------------------------------------------------------------------- 
    613       ALLOCATE( prmaxn(jpi,jpj,jpk), prmaxp(jpi,jpj,jpk), prmaxd(jpi,jpj,jpk),  & 
    614       &          zdaylen(jpi,jpj), STAT = p5z_prod_alloc ) 
     624      ALLOCATE( zdaylen(jpi,jpj), STAT = p5z_prod_alloc ) 
    615625      ! 
    616626      IF( p5z_prod_alloc /= 0 ) CALL ctl_warn('p5z_prod_alloc : failed to allocate arrays.') 
  • NEMO/trunk/src/TOP/PISCES/SED/oce_sed.F90

    r10225 r10362  
    1111   USE par_trc 
    1212 
    13    USE dom_oce , ONLY :   nidom     =>   nidom          !: 
    1413   USE dom_oce , ONLY :   glamt     =>   glamt          !: longitude of t-point (degre) 
    1514   USE dom_oce , ONLY :   gphit     =>   gphit          !: latitude  of t-point (degre) 
     
    2120   USE dom_oce , ONLY :   rdt       =>   rdt            !: time step for the dynamics 
    2221   USE dom_oce , ONLY :   nyear     =>   nyear          !: Current year 
    23    USE dom_oce , ONLY :   nmonth    =>   nmonth         !: Current month 
    24    USE dom_oce , ONLY :   nday      =>   nday           !: Current day 
    2522   USE dom_oce , ONLY :   ndastp    =>   ndastp         !: time step date in year/month/day aammjj 
    26    USE dom_oce , ONLY :   nday_year =>   nday_year      !: curent day counted from jan 1st of the current year 
    2723   USE dom_oce , ONLY :   adatrj    =>   adatrj         !: number of elapsed days since the begining of the run 
    2824   USE trc     , ONLY :  nittrc000  =>   nittrc000 
     
    3531   USE sms_pisces, ONLY : wsbio4    =>   wsbio4          !: sinking flux for POC 
    3632   USE sms_pisces, ONLY : wsbio3    =>   wsbio3          !: sinking flux for GOC 
    37    USE sms_pisces, ONLY : wscal     =>   wscal           !: sinking flux for calcite 
    3833   USE sms_pisces, ONLY : wsbio2    =>   wsbio2           !: sinking flux for calcite 
    3934   USE sms_pisces, ONLY : wsbio     =>   wsbio           !: sinking flux for calcite 
     35   USE sms_pisces, ONLY : ln_p5z    =>   ln_p5z          !: PISCES-QUOTA flag 
    4036   USE p4zche, ONLY     : akb3      =>   akb3            !: Chemical constants   
    4137   USE sms_pisces, ONLY : ak13      =>   ak13            !: Chemical constants   
  • NEMO/trunk/src/TOP/PISCES/SED/seddsr.F90

    r10222 r10362  
    591591            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(4) 
    592592            zdelta = pwcp(ji,jk,jwpo4) - redfep * zsedtra(4) 
    593             zsedtra(4) = ( zsedtra(4) * zalpha ) / ( 0.25 * zsedtra(4) * ( exp( reac_fe2 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    594             &            + zalpha * exp( reac_fe2 * zalpha * dtsed2 / 2. ) + rtrn ) 
     593            IF ( zalpha == 0. ) THEN 
     594               zsedtra(4) = zsedtra(4) / ( 1.0 + zsedtra(4) * reac_fe2 * dtsed2 / 2.0 ) 
     595            ELSE 
     596               zsedtra(4) = ( zsedtra(4) * zalpha ) / ( 0.25 * zsedtra(4) * ( exp( reac_fe2 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     597               &            + zalpha * exp( reac_fe2 * zalpha * dtsed2 / 2. ) ) 
     598            ENDIF 
    595599            zsedtra(1) = zalpha + 0.25 * zsedtra(4)  
    596600            zsedtra(5) = zbeta  - zsedtra(4) 
     
    601605            zbeta  = pwcp(ji,jk,jwso4) + zsedtra(2) 
    602606            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(2) 
    603             zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_h2s * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    604             &            + zalpha * exp( reac_h2s * zalpha * dtsed2 / 2. ) + rtrn ) 
     607            IF ( zalpha == 0. ) THEN 
     608               zsedtra(2) = zsedtra(2) / ( 1.0 + zsedtra(2) * reac_h2s * dtsed2 / 2.0 ) 
     609            ELSE 
     610               zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_h2s * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     611               &            + zalpha * exp( reac_h2s * zalpha * dtsed2 / 2. ) ) 
     612            ENDIF 
    605613            zsedtra(1) = zalpha + 2.0 * zsedtra(2) 
    606614            pwcp(ji,jk,jwalk) = zgamma + 2.0 * zsedtra(2) 
     
    609617            zalpha = zsedtra(1) - 2.0 * zsedtra(3) 
    610618            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(3) 
    611             zsedtra(3) = ( zsedtra(3) * zalpha ) / ( 2.0 * zsedtra(3) * ( exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    612             &            + zalpha * exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 /2. ) + rtrn ) 
     619            IF ( zalpha == 0. ) THEN 
     620               zsedtra(3) = zsedtra(3) / ( 1.0 + zsedtra(3) * reac_nh4 * zadsnh4 * dtsed2 / 2.0 ) 
     621            ELSE 
     622               zsedtra(3) = ( zsedtra(3) * zalpha ) / ( 2.0 * zsedtra(3) * ( exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     623               &            + zalpha * exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 /2. ) ) 
     624            ENDIF 
    613625            zsedtra(1) = zalpha + 2.0 * zsedtra(3) 
    614626            pwcp(ji,jk,jwalk) = zgamma + 2.0 * zsedtra(3) 
     
    617629            zbeta  = zsedtra(4) + zsedtra(6) 
    618630            zgamma = pwcp(ji,jk,jwso4) + zsedtra(6) 
    619             zsedtra(6) = ( zsedtra(6) * zalpha ) / ( 2.0 * zsedtra(6) * ( exp( reac_feso * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    620             &            + zalpha * exp( reac_feso * zalpha * dtsed2 /2. ) + rtrn ) 
     631            IF ( zalpha == 0. ) THEN 
     632               zsedtra(6) = zsedtra(6) / ( 1.0 + zsedtra(6) * reac_feso * dtsed2 / 2.0 ) 
     633            ELSE 
     634               zsedtra(6) = ( zsedtra(6) * zalpha ) / ( 2.0 * zsedtra(6) * ( exp( reac_feso * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     635               &            + zalpha * exp( reac_feso * zalpha * dtsed2 /2. ) ) 
     636            ENDIF 
    621637            zsedtra(1) = zalpha + 2.0 * zsedtra(6) 
    622638            zsedtra(4) = zbeta  - zsedtra(6) 
     
    626642            zbeta  = zsedtra(4) + zsedtra(6) 
    627643            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(4) 
    628             zsedtra(4) = ( zsedtra(4) * zalpha ) / ( zsedtra(4) * ( exp( reac_fes * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    629             &            + zalpha * exp( reac_fes * zalpha * dtsed2 /2. ) + rtrn ) 
     644            IF ( zalpha == 0. ) THEN 
     645               zsedtra(4) = zsedtra(4) / ( 1.0 + zsedtra(4) * reac_fes * dtsed2 / 2.0 ) 
     646            ELSE 
     647               zsedtra(4) = ( zsedtra(4) * zalpha ) / ( zsedtra(4) * ( exp( reac_fes * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     648               &            + zalpha * exp( reac_fes * zalpha * dtsed2 /2. ) ) 
     649            ENDIF 
    630650            zsedtra(2) = zalpha + zsedtra(4) 
    631651            zsedtra(6) = zbeta  - zsedtra(4) 
     
    637657            zdelta = pwcp(ji,jk,jwso4) + zsedtra(2) 
    638658            zepsi  = pwcp(ji,jk,jwpo4) + redfep * zsedtra(5) 
    639             zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_feh2s * zalpha * dtsed2 ) - 1.0 )  & 
    640             &            + zalpha * exp( reac_feh2s * zalpha * dtsed2 ) + rtrn ) 
     659            IF ( zalpha == 0. ) THEN 
     660               zsedtra(2) = zsedtra(2) / ( 1.0 + zsedtra(2) * reac_feh2s * dtsed2 ) 
     661            ELSE 
     662               zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_feh2s * zalpha * dtsed2 ) - 1.0 )  & 
     663               &            + zalpha * exp( reac_feh2s * zalpha * dtsed2 ) ) 
     664            ENDIF 
    641665            zsedtra(5) = zalpha + 2.0 * zsedtra(2) 
    642666            zsedtra(4) = zbeta  - zsedtra(5) 
     
    648672            zbeta  = zsedtra(4) + zsedtra(6) 
    649673            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(4) 
    650             zsedtra(4) = ( zsedtra(4) * zalpha ) / ( zsedtra(4) * ( exp( reac_fes * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    651             &            + zalpha * exp( reac_fes * zalpha * dtsed2 /2. ) + rtrn ) 
     674            IF ( zalpha == 0. ) THEN 
     675               zsedtra(4) = zsedtra(4) / ( 1.0 + zsedtra(4) * reac_fes * dtsed2 / 2.0 ) 
     676            ELSE 
     677               zsedtra(4) = ( zsedtra(4) * zalpha ) / ( zsedtra(4) * ( exp( reac_fes * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     678               &            + zalpha * exp( reac_fes * zalpha * dtsed2 /2. ) ) 
     679            ENDIF 
    652680            zsedtra(2) = zalpha + zsedtra(4) 
    653681            zsedtra(6) = zbeta  - zsedtra(4) 
     
    657685            zbeta  = zsedtra(4) + zsedtra(6) 
    658686            zgamma = pwcp(ji,jk,jwso4) + zsedtra(6) 
    659             zsedtra(6) = ( zsedtra(6) * zalpha ) / ( 2.0 * zsedtra(6) * ( exp( reac_feso * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    660             &            + zalpha * exp( reac_feso * zalpha * dtsed2 /2. ) + rtrn ) 
     687            IF (zalpha == 0.) THEN 
     688               zsedtra(6) = zsedtra(6) / ( 1.0 + zsedtra(6) * reac_feso * dtsed2 / 2. ) 
     689            ELSE 
     690               zsedtra(6) = ( zsedtra(6) * zalpha ) / ( 2.0 * zsedtra(6) * ( exp( reac_feso * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     691               &            + zalpha * exp( reac_feso * zalpha * dtsed2 /2. ) ) 
     692            ENDIF 
    661693            zsedtra(1) = zalpha + 2.0 * zsedtra(6) 
    662694            zsedtra(4) = zbeta  - zsedtra(6) 
     
    665697            zalpha = zsedtra(1) - 2.0 * zsedtra(3) 
    666698            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(3) 
    667             zsedtra(3) = ( zsedtra(3) * zalpha ) / ( 2.0 * zsedtra(3) * ( exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    668             &            + zalpha * exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 /2. ) + rtrn ) 
     699            IF (zalpha == 0.) THEN 
     700               zsedtra(3) = zsedtra(3) / ( 1.0 + zsedtra(3) * reac_nh4 * zadsnh4 * dtsed2 / 2.0)  
     701            ELSE 
     702               zsedtra(3) = ( zsedtra(3) * zalpha ) / ( 2.0 * zsedtra(3) * ( exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     703               &            + zalpha * exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 /2. ) ) 
     704            ENDIF 
    669705            zsedtra(1) = zalpha + 2.0 * zsedtra(3) 
    670706            pwcp(ji,jk,jwalk) = zgamma + 2.0 * zsedtra(3) 
     
    673709            zbeta  = pwcp(ji,jk,jwso4) + zsedtra(2) 
    674710            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(2) 
    675             zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_h2s * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    676             &            + zalpha * exp( reac_h2s * zalpha * dtsed2 / 2. ) + rtrn ) 
     711            IF ( zalpha == 0. ) THEN 
     712               zsedtra(2) = zsedtra(2) / ( 1.0 + zsedtra(2) * reac_h2s * dtsed2 / 2.0 ) 
     713            ELSE 
     714               zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_h2s * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     715               &            + zalpha * exp( reac_h2s * zalpha * dtsed2 / 2. ) ) 
     716            ENDIF 
    677717            zsedtra(1) = zalpha + 2.0 * zsedtra(2) 
    678718            pwcp(ji,jk,jwso4) = zbeta - zsedtra(2) 
     
    683723            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(4) 
    684724            zdelta = pwcp(ji,jk,jwpo4) - redfep * zsedtra(4) 
    685             zsedtra(4) = ( zsedtra(4) * zalpha ) / ( 0.25 * zsedtra(4) * ( exp( reac_fe2 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    686             &            + zalpha * exp( reac_fe2 * zalpha * dtsed2 / 2. ) + rtrn ) 
     725            IF ( zalpha == 0. ) THEN 
     726               zsedtra(4) = zsedtra(4) / ( 1.0 + zsedtra(4) * reac_fe2 * dtsed2 / 2.0 ) 
     727            ELSE 
     728               zsedtra(4) = ( zsedtra(4) * zalpha ) / ( 0.25 * zsedtra(4) * ( exp( reac_fe2 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     729               &            + zalpha * exp( reac_fe2 * zalpha * dtsed2 / 2. ) ) 
     730            ENDIF 
    687731            zsedtra(1) = zalpha + 0.25 * zsedtra(4) 
    688732            zsedtra(5) = zbeta  - zsedtra(4) 
  • NEMO/trunk/src/TOP/PISCES/SED/seddta.F90

    r10333 r10362  
    4949 
    5050      REAL(wp), DIMENSION(jpoce) :: zdtap, zdtag 
    51       REAL(wp), DIMENSION(jpi,jpj) :: zwsbio4, zwsbio3, zwscal 
     51      REAL(wp), DIMENSION(jpi,jpj) :: zwsbio4, zwsbio3 
    5252      REAL(wp) :: zf0, zf1, zf2, zkapp, zratio, zdep 
    5353 
     
    9797               zwsbio4(ji,jj) = wsbio2 / rday 
    9898               zwsbio3(ji,jj) = wsbio  / rday 
    99                zwscal (ji,jj) = wsbio2 / rday 
    10099            END DO 
    101100         END DO 
     
    106105               zdep = e3t_n(ji,jj,ikt) / r2dttrc 
    107106               zwsbio4(ji,jj) = MIN( 0.99 * zdep, wsbio4(ji,jj,ikt) / rday ) 
    108                zwscal (ji,jj) = MIN( 0.99 * zdep, wscal (ji,jj,ikt) / rday ) 
    109107               zwsbio3(ji,jj) = MIN( 0.99 * zdep, wsbio3(ji,jj,ikt) / rday ) 
    110108            END DO 
     
    130128               trc_data(ji,jj,12 ) = MIN(trb(ji,jj,ikt,jppoc), 1E-4) * zwsbio3(ji,jj) * 1E3 
    131129               trc_data(ji,jj,13 ) = MIN(trb(ji,jj,ikt,jpgoc), 1E-4) * zwsbio4(ji,jj) * 1E3 
    132                trc_data(ji,jj,14)  = MIN(trb(ji,jj,ikt,jpcal), 1E-4) * zwscal (ji,jj) * 1E3 
     130               trc_data(ji,jj,14)  = MIN(trb(ji,jj,ikt,jpcal), 1E-4) * zwsbio4(ji,jj) * 1E3 
    133131               trc_data(ji,jj,15)  = tsn(ji,jj,ikt,jp_tem) 
    134132               trc_data(ji,jj,16)  = tsn(ji,jj,ikt,jp_sal) 
  • NEMO/trunk/src/TOP/PISCES/SED/sedini.F90

    r10333 r10362  
    472472      ENDIF 
    473473 
     474      IF ( ln_p5z .AND. ln_sed_2way ) CALL ctl_stop( '2 ways coupling with sediment cannot be activated with PISCES-QUOTA' ) 
     475 
    474476      REWIND( numnamsed_ref )              ! Namelist nam_geom in reference namelist : Pisces variables 
    475477      READ  ( numnamsed_ref, nam_geom, IOSTAT = ios, ERR = 903) 
  • NEMO/trunk/src/TOP/PISCES/sms_pisces.F90

    r10222 r10362  
    7676   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  etot_ndcy      !: PAR over 24h in case of diurnal cycle 
    7777   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  enano, ediat   !: PAR for phyto, nano and diat  
     78   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  enanom, ediatm !: PAR for phyto, nano and diat  
    7879   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  epico          !: PAR for pico 
     80   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  epicom         !: PAR for pico 
    7981   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  emoy           !: averaged PAR in the mixed layer 
    8082   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::  heup_01 !: Absolute euphotic layer depth 
     
    8991   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wsbio3   !: POC sinking speed  
    9092   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wsbio4   !: GOC sinking speed 
    91    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wscal    !: Calcite and BSi sinking speeds 
    9293   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wsfep 
    9394 
     
    105106   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prodgoc    !: Calcite production 
    106107   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   consgoc    !: Calcite production 
    107  
     108   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   blim       !: bacterial production factor 
    108109   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sizen      !: size of diatoms  
    109110   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sizep      !: size of diatoms  
     
    147148         !*  Biological fluxes for light  
    148149         ALLOCATE(  enano(jpi,jpj,jpk)    , ediat(jpi,jpj,jpk) ,   & 
     150           &        enanom(jpi,jpj,jpk)   , ediatm(jpi,jpj,jpk),   & 
    149151           &        etot_ndcy(jpi,jpj,jpk), emoy(jpi,jpj,jpk)  ,  STAT=ierr(2) )  
    150152 
     
    157159            &      prodcal(jpi,jpj,jpk) , xdiss   (jpi,jpj,jpk),    & 
    158160            &      prodpoc(jpi,jpj,jpk) , conspoc(jpi,jpj,jpk) ,    & 
    159             &      prodgoc(jpi,jpj,jpk) , consgoc(jpi,jpj,jpk) ,  STAT=ierr(4) ) 
     161            &      prodgoc(jpi,jpj,jpk) , consgoc(jpi,jpj,jpk) ,    & 
     162            &      blim   (jpi,jpj,jpk) ,                         STAT=ierr(4) ) 
    160163 
    161164         !* Variable for chemistry of the CO2 cycle 
     
    170173         !* Sinkong speed 
    171174         ALLOCATE( wsbio3 (jpi,jpj,jpk) , wsbio4 (jpi,jpj,jpk),     & 
    172             &      wscal(jpi,jpj,jpk)                         ,   STAT=ierr(7) )    
     175            &                             STAT=ierr(7) )    
    173176         !  
    174177         IF( ln_ligand ) THEN 
     
    180183      IF( ln_p5z ) THEN 
    181184         !        
    182          ALLOCATE( epico(jpi,jpj,jpk)                         ,   STAT=ierr(9) )  
     185         ALLOCATE( epico(jpi,jpj,jpk)   , epicom(jpi,jpj,jpk) ,   STAT=ierr(9) )  
    183186 
    184187         !*  Size of phytoplankton cells 
  • NEMO/trunk/src/TOP/PISCES/trcini_pisces.F90

    r10227 r10362  
    109109      ierr = ierr +  p4z_flx_alloc() 
    110110      ierr = ierr +  p4z_sed_alloc() 
    111       ierr = ierr +  p4z_rem_alloc() 
     111      ierr = ierr +  p4z_lim_alloc() 
    112112      IF( ln_p4z ) THEN 
    113          ierr = ierr +  p4z_lim_alloc() 
    114113         ierr = ierr +  p4z_prod_alloc() 
    115114      ELSE 
     
    117116         ierr = ierr +  p5z_prod_alloc() 
    118117      ENDIF 
     118      ierr = ierr +  p4z_rem_alloc() 
    119119      ! 
    120120      IF( lk_mpp    )   CALL mpp_sum( ierr ) 
     
    165165        IF( cltra == 'PICN'     )   jpnpi = jn      !: Picophytoplankton N biomass 
    166166        IF( cltra == 'PICP'     )   jpppi = jn      !: Picophytoplankton P biomass 
     167        IF( cltra == 'PCHL'     )   jppch = jn      !: Diatoms Chlorophyll Concentration 
    167168        IF( cltra == 'PFe'      )   jppfe = jn      !: Picophytoplankton Fe biomass 
    168169        IF( cltra == 'LGW'      )   jplgw = jn      !: Weak ligands 
     
    238239         xksi(:,:)    = 2.e-6 
    239240         xksimax(:,:) = xksi(:,:) 
     241         IF( ln_p5z ) THEN 
     242            sized(:,:,:) = 1.0 
     243            sizen(:,:,:) = 1.0 
     244            sized(:,:,:) = 1.0 
     245         ENDIF 
    240246      END IF 
    241247 
Note: See TracChangeset for help on using the changeset viewer.