New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 10401 for NEMO/trunk/src/TOP – NEMO

Changeset 10401 for NEMO/trunk/src/TOP


Ignore:
Timestamp:
2018-12-17T16:57:34+01:00 (5 years ago)
Author:
cetlod
Message:

Remove obsolescence in PISCES

Location:
NEMO/trunk/src/TOP/PISCES/P4Z
Files:
2 edited

Legend:

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

    r10362 r10401  
    2525   PUBLIC   p4z_fechem_init   ! called in trcsms_pisces.F90 
    2626 
    27    LOGICAL          ::   ln_fechem    !: boolean for complex iron chemistry following Tagliabue and voelker 
    2827   LOGICAL          ::   ln_ligvar    !: boolean for variable ligand concentration following Tagliabue and voelker 
    2928   REAL(wp), PUBLIC ::   xlam1        !: scavenging rate of Iron  
     
    3231   REAL(wp), PUBLIC ::   kfep         !: rate constant for nanoparticle formation 
    3332 
    34    REAL(wp) :: kl1, kl2, kb1, kb2, ks, kpr, spd, con, kth      !!gm  <<<== DOCTOR names SVP !!! 
    35  
    3633   !!---------------------------------------------------------------------- 
    3734   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
     
    4744      !! ** Purpose :   Compute remineralization/scavenging of iron 
    4845      !! 
    49       !! ** Method  :   2 different chemistry models are available for iron 
    50       !!                (1) The simple chemistry model of Aumont and Bopp (2006) 
    51       !!                    based on one ligand and one inorganic form 
    52       !!                (2) The complex chemistry model of Tagliabue and  
    53       !!                    Voelker (2009) based on 2 ligands, 2 inorganic forms 
    54       !!                    and one particulate form (ln_fechem) 
     46      !! ** Method  :   A simple chemistry model of iron from Aumont and Bopp (2006) 
     47      !!                based on one ligand and one inorganic form 
    5548      !!--------------------------------------------------------------------- 
    5649      INTEGER, INTENT(in) ::   kt, knt   ! ocean time step 
     
    7467      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zTL1, zFe3, ztotlig, precip, zFeL1 
    7568      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zcoll3d, zscav3d, zlcoll3d 
    76       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zFeL2, zTL2, zFe2, zFeP 
    77       REAL(wp), ALLOCATABLE, DIMENSION(:,:  ) ::   zstrn, zstrn2 
    7869      !!--------------------------------------------------------------------- 
    7970      ! 
     
    8374      zFeL1(:,:,:) = 0. 
    8475      zTL1 (:,:,:) = 0. 
    85       IF( ln_fechem ) THEN 
    86          ALLOCATE( zstrn(jpi,jpj), zstrn2(jpi,jpj) ) 
    87          ALLOCATE( zFe2(jpi,jpj,jpk), zFeL2(jpi,jpj,jpk), zTL2(jpi,jpj,jpk), zFeP(jpi,jpj,jpk) ) 
    88          zFe2 (:,:,:) = 0. 
    89          zFeL2(:,:,:) = 0. 
    90          zTL2 (:,:,:) = 0. 
    91          zFeP (:,:,:) = 0. 
    92       ENDIF 
    9376 
    9477      ! Total ligand concentration : Ligands can be chosen to be constant or variable 
     
    10487      ENDIF 
    10588 
    106       IF( ln_fechem ) THEN 
    107          ! compute the day length depending on latitude and the day 
    108          zrum = REAL( nday_year - 80, wp ) / REAL( nyear_len(1), wp ) 
    109          zcodel = ASIN(  SIN( zrum * rpi * 2._wp ) * SIN( rad * 23.5_wp )  ) 
    110  
    111          ! day length in hours 
    112          zstrn(:,:) = 0. 
     89      ! ------------------------------------------------------------ 
     90      !  from Aumont and Bopp (2006) 
     91      ! This model is based on one ligand and Fe'  
     92      ! Chemistry is supposed to be fast enough to be at equilibrium 
     93      ! ------------------------------------------------------------ 
     94      DO jk = 1, jpkm1 
    11395         DO jj = 1, jpj 
    11496            DO ji = 1, jpi 
    115                zargu = TAN( zcodel ) * TAN( gphit(ji,jj) * rad ) 
    116                zargu = MAX( -1., MIN(  1., zargu ) ) 
    117                zstrn(ji,jj) = MAX( 0.0, 24. - 2. * ACOS( zargu ) / rad / 15. ) 
    118             END DO 
     97               zTL1(ji,jj,jk)  = ztotlig(ji,jj,jk) 
     98               zkeq            = fekeq(ji,jj,jk) 
     99               zfesatur        = zTL1(ji,jj,jk) * 1E-9 
     100               ztfe            = trb(ji,jj,jk,jpfer)  
     101               ! Fe' is the root of a 2nd order polynom 
     102               zFe3 (ji,jj,jk) = ( -( 1. + zfesatur * zkeq - zkeq * ztfe )               & 
     103                  &              + SQRT( ( 1. + zfesatur * zkeq - zkeq * ztfe )**2       & 
     104                  &              + 4. * ztfe * zkeq) ) / ( 2. * zkeq ) 
     105               zFe3 (ji,jj,jk) = zFe3(ji,jj,jk) * 1E9 
     106               zFeL1(ji,jj,jk) = MAX( 0., trb(ji,jj,jk,jpfer) * 1E9 - zFe3(ji,jj,jk) ) 
     107           END DO 
    119108         END DO 
    120  
    121          ! Maximum light intensity 
    122          zstrn2(:,:) = zstrn(:,:) / 24. 
    123          WHERE( zstrn(:,:) < 1.e0 ) zstrn(:,:) = 24. 
    124          zstrn(:,:) = 24. / zstrn(:,:) 
    125  
    126          ! ------------------------------------------------------------ 
    127          ! NEW FE CHEMISTRY ROUTINE from Tagliabue and Volker (2009) 
    128          ! This model is based on two ligands, Fe2+, Fe3+ and Fep 
    129          ! Chemistry is supposed to be fast enough to be at equilibrium 
    130          ! ------------------------------------------------------------ 
    131          DO jn = 1, 2 
    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 & 
    167                       & + zkappa1 * zkappa2 - ( zkappa1 + zkappa2 ) * ztfe / za 
    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. 
    181                      ELSE 
    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 
    195                      ENDIF 
    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 
    209                END DO 
    210             END DO 
    211          END DO 
    212       ELSE 
    213          ! ------------------------------------------------------------ 
    214          ! OLD FE CHEMISTRY ROUTINE from Aumont and Bopp (2006) 
    215          ! This model is based on one ligand and Fe'  
    216          ! Chemistry is supposed to be fast enough to be at equilibrium 
    217          ! ------------------------------------------------------------ 
    218          DO jk = 1, jpkm1 
    219             DO jj = 1, jpj 
    220                DO ji = 1, jpi 
    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)  
    225                   ! Fe' is the root of a 2nd order polynom 
    226                   zFe3 (ji,jj,jk) = ( -( 1. + zfesatur * zkeq - zkeq * ztfe )               & 
    227                      &              + SQRT( ( 1. + zfesatur * zkeq - zkeq * ztfe )**2       & 
    228                      &              + 4. * ztfe * zkeq) ) / ( 2. * zkeq ) 
    229                   zFe3 (ji,jj,jk) = zFe3(ji,jj,jk) * 1E9 
    230                   zFeL1(ji,jj,jk) = MAX( 0., trb(ji,jj,jk,jpfer) * 1E9 - zFe3(ji,jj,jk) ) 
    231               END DO 
    232             END DO 
    233          END DO 
    234          ! 
    235       ENDIF 
     109      END DO 
     110         ! 
    236111 
    237112      zdust = 0.         ! if no dust available 
     
    247122               &         + fesol(ji,jj,jk,3) * zhplus + fesol(ji,jj,jk,4)     & 
    248123               &         + fesol(ji,jj,jk,5) / zhplus ) 
    249                IF( ln_fechem ) THEN 
    250                   zfeequi = ( zFe3(ji,jj,jk) + zFe2(ji,jj,jk) + zFeP(ji,jj,jk) ) * 1E-9 
    251                   zfecoll = ( 0.3 * zFeL1(ji,jj,jk) + 0.5 * zFeL2(ji,jj,jk) ) * 1E-9 
    252                   precip(ji,jj,jk) = 0.0 
    253                ELSE 
    254                   zfeequi = zFe3(ji,jj,jk) * 1E-9 
    255                   zhplus  = max( rtrn, hi(ji,jj,jk) ) 
    256                   fe3sol  = fesol(ji,jj,jk,1) * ( zhplus**3 + fesol(ji,jj,jk,2) * zhplus**2  & 
     124               ! 
     125               zfeequi = zFe3(ji,jj,jk) * 1E-9 
     126               zhplus  = max( rtrn, hi(ji,jj,jk) ) 
     127               fe3sol  = fesol(ji,jj,jk,1) * ( zhplus**3 + fesol(ji,jj,jk,2) * zhplus**2  & 
    257128                  &         + fesol(ji,jj,jk,3) * zhplus + fesol(ji,jj,jk,4)     & 
    258129                  &         + 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 
    262                ENDIF 
     130               zfecoll = 0.5 * zFeL1(ji,jj,jk) * 1E-9 
     131               ! precipitation of Fe3+, creation of nanoparticles 
     132               precip(ji,jj,jk) = MAX( 0., ( zFe3(ji,jj,jk) * 1E-9 - fe3sol ) ) * kfep * xstep 
    263133               ! 
    264134               ztrc   = ( trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + trb(ji,jj,jk,jpcal) + trb(ji,jj,jk,jpgsi) ) * 1.e6  
     
    311181      !  Define the bioavailable fraction of iron 
    312182      !  ---------------------------------------- 
    313       IF( ln_fechem ) THEN  ;  biron(:,:,:) = MAX( 0., trb(:,:,:,jpfer) - zFeP(:,:,:) * 1E-9 ) 
    314       ELSE                  ;  biron(:,:,:) = trb(:,:,:,jpfer)  
    315       ENDIF 
     183      biron(:,:,:) = trb(:,:,:,jpfer)  
    316184      ! 
    317185      IF( ln_ligand ) THEN 
     
    334202         END DO 
    335203         ! 
    336          IF( .NOT.ln_fechem) THEN 
    337             plig(:,:,:) =  MAX( 0., ( ( zFeL1(:,:,:) * 1E-9 ) / ( trb(:,:,:,jpfer) +rtrn ) ) ) 
    338          ENDIF 
     204         plig(:,:,:) =  MAX( 0., ( ( zFeL1(:,:,:) * 1E-9 ) / ( trb(:,:,:,jpfer) +rtrn ) ) ) 
    339205         ! 
    340206      ENDIF 
     
    343209      IF( lk_iomput ) THEN 
    344210         IF( knt == nrdttrc ) THEN 
    345          zrfact2 = 1.e3 * rfact2r  ! conversion from mol/L/timestep into mol/m3/s 
    346          IF( iom_use("Fe3")    )  CALL iom_put("Fe3"    , zFe3   (:,:,:)       * tmask(:,:,:) )   ! Fe3+ 
    347          IF( iom_use("FeL1")   )  CALL iom_put("FeL1"   , zFeL1  (:,:,:)       * tmask(:,:,:) )   ! FeL1 
    348          IF( iom_use("TL1")    )  CALL iom_put("TL1"    , zTL1   (:,:,:)       * tmask(:,:,:) )   ! TL1 
    349          IF( iom_use("Totlig") )  CALL iom_put("Totlig" , ztotlig(:,:,:)       * tmask(:,:,:) )   ! TL 
    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 ) 
    354          IF( ln_fechem ) THEN 
    355             IF( iom_use("Fe2")  ) CALL iom_put("Fe2"    , zFe2   (:,:,:)       * tmask(:,:,:) )   ! Fe2+ 
    356             IF( iom_use("FeL2") ) CALL iom_put("FeL2"   , zFeL2  (:,:,:)       * tmask(:,:,:) )   ! FeL2 
    357             IF( iom_use("FeP")  ) CALL iom_put("FeP"    , zFeP   (:,:,:)       * tmask(:,:,:) )   ! FeP 
    358             IF( iom_use("TL2")  ) CALL iom_put("TL2"    , zTL2   (:,:,:)       * tmask(:,:,:) )   ! TL2 
    359          ENDIF 
     211            zrfact2 = 1.e3 * rfact2r  ! conversion from mol/L/timestep into mol/m3/s 
     212            IF( iom_use("Fe3")    )  CALL iom_put("Fe3"    , zFe3   (:,:,:)       * tmask(:,:,:) )   ! Fe3+ 
     213            IF( iom_use("FeL1")   )  CALL iom_put("FeL1"   , zFeL1  (:,:,:)       * tmask(:,:,:) )   ! FeL1 
     214            IF( iom_use("TL1")    )  CALL iom_put("TL1"    , zTL1   (:,:,:)       * tmask(:,:,:) )   ! TL1 
     215            IF( iom_use("Totlig") )  CALL iom_put("Totlig" , ztotlig(:,:,:)       * tmask(:,:,:) )   ! TL 
     216            IF( iom_use("Biron")  )  CALL iom_put("Biron"  , biron  (:,:,:)  * 1e9 * tmask(:,:,:) )   ! biron 
     217            IF( iom_use("FESCAV") )  CALL iom_put("FESCAV" , zscav3d(:,:,:)  * 1e9 * tmask(:,:,:) * zrfact2 ) 
     218            IF( iom_use("FECOLL") )  CALL iom_put("FECOLL" , zcoll3d(:,:,:)  * 1e9 * tmask(:,:,:) * zrfact2 ) 
     219            IF( iom_use("LGWCOLL"))  CALL iom_put("LGWCOLL", zlcoll3d(:,:,:) * 1e9 * tmask(:,:,:) * zrfact2 ) 
    360220         ENDIF 
    361221      ENDIF 
     
    367227      ENDIF 
    368228      ! 
    369       IF( ln_fechem )  THEN 
    370          DEALLOCATE( zstrn, zstrn2 ) 
    371          DEALLOCATE( zFe2, zFeL2, zTL2, zFeP ) 
    372       ENDIF 
    373       ! 
    374229      IF( ln_timing )   CALL timing_stop('p4z_fechem') 
    375230      ! 
     
    391246      INTEGER ::   ios   ! Local integer  
    392247      !! 
    393       NAMELIST/nampisfer/ ln_fechem, ln_ligvar, xlam1, xlamdust, ligand, kfep  
     248      NAMELIST/nampisfer/ ln_ligvar, xlam1, xlamdust, ligand, kfep  
    394249      !!---------------------------------------------------------------------- 
    395250      ! 
     
    410265      IF(lwp) THEN                     ! control print 
    411266         WRITE(numout,*) '   Namelist : nampisfer' 
    412          WRITE(numout,*) '      enable complex iron chemistry scheme      ln_fechem    =', ln_fechem 
    413267         WRITE(numout,*) '      variable concentration of ligand          ln_ligvar    =', ln_ligvar 
    414268         WRITE(numout,*) '      scavenging rate of Iron                   xlam1        =', xlam1 
     
    418272      ENDIF 
    419273      !  
    420       IF (ln_ligand .AND. ln_fechem) CALL ctl_stop( 'STOP', 'p4z_fechem_init: ln_ligand and ln_fechem are incompatible') 
    421       ! 
    422       IF( ln_fechem ) THEN             ! set some constants used by the complexe chemistry scheme 
    423          ! 
    424          spd = 3600. * 24. 
    425          con = 1.E9 
    426          ! LIGAND KINETICS (values from Witter et al. 2000) 
    427          ! Weak (L2) ligands 
    428          ! Phaeophytin 
    429          kl2 = 12.2E5  * spd / con 
    430          kb2 = 12.3E-6 * spd 
    431          ! Strong (L1) ligands 
    432          ! Saccharides 
    433          ! kl1 = 12.2E5  * spd / con 
    434          ! kb1 = 12.3E-6 * spd 
    435          ! DFOB- 
    436          kl1 = 19.6e5  * spd / con 
    437          kb1 = 1.5e-6  * spd 
    438          ! pcp and remin of Fe3p 
    439          ks  = 0.075 
    440          kpr = 0.05 
    441          ! thermal reduction of Fe3 
    442          kth = 0.0048 * 24. 
    443          ! 
    444       ENDIF 
    445       ! 
    446274   END SUBROUTINE p4z_fechem_init 
    447275    
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zprod.F90

    r10362 r10401  
    2626   PUBLIC   p4z_prod_alloc 
    2727 
    28    LOGICAL , PUBLIC ::   ln_newprod   !: 
    2928   REAL(wp), PUBLIC ::   pislopen     !: 
    3029   REAL(wp), PUBLIC ::   pisloped     !: 
     
    156155      END DO 
    157156 
    158       IF( ln_newprod ) THEN 
    159          DO jk = 1, jpkm1 
    160             DO jj = 1, jpj 
    161                DO ji = 1, jpi 
    162                   IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
    163                       ! Computation of production function for Carbon 
    164                       !  --------------------------------------------- 
    165                       zpislopen = zpislopeadn(ji,jj,jk) / ( ( r1_rday + bresp * r1_rday ) & 
    166                       &            * zmxl_fac(ji,jj,jk) * rday + rtrn) 
    167                       zpisloped = zpislopeadd(ji,jj,jk) / ( ( r1_rday + bresp * r1_rday ) & 
    168                       &            * zmxl_fac(ji,jj,jk) * rday + rtrn) 
    169                       zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) )  ) 
    170                       zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) )  ) 
    171                       !  Computation of production function for Chlorophyll 
    172                       !-------------------------------------------------- 
    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) ) ) 
    177                   ENDIF 
    178                END DO 
    179             END DO 
    180          END DO 
    181       ELSE 
    182          DO jk = 1, jpkm1 
    183             DO jj = 1, jpj 
    184                DO ji = 1, jpi 
    185                   IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
    186                       ! Computation of production function for Carbon 
    187                       !  --------------------------------------------- 
    188                       zpislopen = zpislopeadn(ji,jj,jk)  / ( zprbio(ji,jj,jk) * rday * xlimphy(ji,jj,jk) + rtrn ) 
    189                       zpisloped = zpislopeadd(ji,jj,jk) / ( zprdia(ji,jj,jk) * rday * xlimdia(ji,jj,jk) + rtrn ) 
    190                       zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) ) ) 
    191                       zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) ) ) 
    192                       !  Computation of production function for Chlorophyll 
    193                       !-------------------------------------------------- 
    194                       zpislopen = zpislopen * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    195                       zpisloped = zpisloped * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    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) ) ) 
    198                   ENDIF 
    199                END DO 
    200             END DO 
    201          END DO 
    202       ENDIF 
     157      DO jk = 1, jpkm1 
     158         DO jj = 1, jpj 
     159            DO ji = 1, jpi 
     160               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
     161                   ! Computation of production function for Carbon 
     162                   !  --------------------------------------------- 
     163                   zpislopen = zpislopeadn(ji,jj,jk) / ( ( r1_rday + bresp * r1_rday ) & 
     164                   &            * zmxl_fac(ji,jj,jk) * rday + rtrn) 
     165                   zpisloped = zpislopeadd(ji,jj,jk) / ( ( r1_rday + bresp * r1_rday ) & 
     166                   &            * zmxl_fac(ji,jj,jk) * rday + rtrn) 
     167                   zprbio(ji,jj,jk) = zprbio(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) )  ) 
     168                   zprdia(ji,jj,jk) = zprdia(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) )  ) 
     169                   !  Computation of production function for Chlorophyll 
     170                   !-------------------------------------------------- 
     171                   zpislopen = zpislopeadn(ji,jj,jk) / ( zprmaxn(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn ) 
     172                   zpisloped = zpislopeadd(ji,jj,jk) / ( zprmaxd(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn ) 
     173                   zprnch(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) ) ) 
     174                   zprdch(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) ) ) 
     175               ENDIF 
     176            END DO 
     177         END DO 
     178      END DO 
    203179 
    204180      !  Computation of a proxy of the N/C ratio 
     
    505481      INTEGER ::   ios   ! Local integer 
    506482      ! 
    507       NAMELIST/namp4zprod/ pislopen, pisloped, xadap, ln_newprod, bresp, excretn, excretd,  & 
     483      NAMELIST/namp4zprod/ pislopen, pisloped, xadap, bresp, excretn, excretd,  & 
    508484         &                 chlcnm, chlcdm, chlcmin, fecnm, fecdm, grosip 
    509485      !!---------------------------------------------------------------------- 
     
    525501      IF(lwp) THEN                         ! control print 
    526502         WRITE(numout,*) '   Namelist : namp4zprod' 
    527          WRITE(numout,*) '      Enable new parame. of production (T/F)   ln_newprod    =', ln_newprod 
    528503         WRITE(numout,*) '      mean Si/C ratio                           grosip       =', grosip 
    529504         WRITE(numout,*) '      P-I slope                                 pislopen     =', pislopen 
     
    531506         WRITE(numout,*) '      excretion ratio of nanophytoplankton      excretn      =', excretn 
    532507         WRITE(numout,*) '      excretion ratio of diatoms                excretd      =', excretd 
    533          IF( ln_newprod )  THEN 
    534             WRITE(numout,*) '      basal respiration in phytoplankton        bresp        =', bresp 
    535             WRITE(numout,*) '      Maximum Chl/C in phytoplankton            chlcmin      =', chlcmin 
    536          ENDIF 
     508         WRITE(numout,*) '      basal respiration in phytoplankton        bresp        =', bresp 
     509         WRITE(numout,*) '      Maximum Chl/C in phytoplankton            chlcmin      =', chlcmin 
    537510         WRITE(numout,*) '      P-I slope  for diatoms                    pisloped     =', pisloped 
    538511         WRITE(numout,*) '      Minimum Chl/C in nanophytoplankton        chlcnm       =', chlcnm 
Note: See TracChangeset for help on using the changeset viewer.