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 14786 – NEMO

Changeset 14786


Ignore:
Timestamp:
2021-05-05T09:09:10+02:00 (3 years ago)
Author:
aumont
Message:

Various bug fixes in PISCES

Location:
NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP
Files:
20 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zbc.F90

    r14416 r14786  
    108108             ! surface downward dust depo of iron 
    109109             jl = n_trc_indsbc(jpfer) 
    110              CALL iom_put( "Irondep", ( rf_trsfac(jl) * sf_trcsbc(jl)%fnow(:,:,1) / rn_sbc_time ) * 1.e+3 * rfactr * tmask(:,:,1) ) 
     110             CALL iom_put( "Irondep", ( rf_trsfac(jl) * sf_trcsbc(jl)%fnow(:,:,1) / rn_sbc_time ) * 1.e+3 * tmask(:,:,1) ) 
    111111 
    112112             ! dust concentration at surface 
    113              CALL iom_put( "pdust"  , dust(:,:) / ( wdust * rday ) * tmask(:,:,1) ) ! dust concentration at surface 
     113             CALL iom_put( "pdust"  , dust(:,:) / ( wdust / rday ) * tmask(:,:,1) ) ! dust concentration at surface 
    114114         ENDIF 
    115115      ENDIF 
     
    143143            DO_2D( 1, 1, 1, 1 ) 
    144144               zndep = rf_trsfac(jl) * sf_trcsbc(jl)%fnow(ji,jj,1) / e3t(ji,jj,1,Kmm) / rn_sbc_time 
    145                tr(ji,jj,1,jptal,Krhs) = tr(ji,jj,1,jptal,Krhs) - rno3 * zndep * rfact 
     145               tr(ji,jj,1,jptal,Krhs) = tr(ji,jj,1,jptal,Krhs) + rno3 * zndep * rfact 
    146146            END_2D 
    147147         ENDIF 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zche.F90

    r14086 r14786  
    195195         !                             ! LN(K0) OF SOLUBILITY OF CO2 (EQ. 12, WEISS, 1980) 
    196196         !                             !     AND FOR THE ATMOSPHERE FOR NON IDEAL GAS 
    197          zcek1 = 9345.17/ztkel - 60.2409 + 23.3585 * LOG(zt) + zsal*(0.023517 - 0.00023656*ztkel    & 
    198          &       + 0.0047036e-4*ztkel**2) 
    199          chemc(ji,jj,1) = EXP( zcek1 ) * 1E-6 * rhop(ji,jj,1) / 1000. ! mol/(L atm) 
     197         zcek1 = 9050.69/ztkel - 58.0931 + 22.2940 * LOG(zt) + zsal*(0.027766 - 0.00025888*ztkel    & 
     198         &       + 0.0050578e-4*ztkel**2) 
     199         chemc(ji,jj,1) = EXP( zcek1 ) * 1E-6   ! mol/(L atm) 
    200200         chemc(ji,jj,2) = -1636.75 + 12.0408*ztkel - 0.0327957*ztkel**2 + 0.0000316528*ztkel**3 
    201201         chemc(ji,jj,3) = 57.7 - 0.118*ztkel 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zfechem.F90

    r14416 r14786  
    120120         &         + fesol(ji,jj,jk,5) / zhplus ) 
    121121         ! 
    122          zfecoll = 0.5 * zFeL1(ji,jj,jk) * 1E-9 
     122         zfecoll = 0.5 * zFeL1(ji,jj,jk) 
    123123         ! precipitation of Fe3+, creation of nanoparticles 
    124124         zprecip = MAX( 0., ( zFe3(ji,jj,jk) - fe3sol ) ) * kfep * xstep * ( 1.0 - nitrfac(ji,jj,jk) )  
     
    198198             zlcoll3d(ji,jj,jk)  = zaggliga + zaggligb 
    199199         END_3D 
    200          plig(:,:,:) =  MAX( 0., ( ( zFeL1(:,:,:) * 1E-9 ) / ( tr(:,:,:,jpfer,Kbb) +rtrn ) ) ) 
    201          ! 
    202       ENDIF 
     200      ENDIF 
     201 
    203202      !  Output of some diagnostics variables 
    204203      !     --------------------------------- 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zflx.F90

    r14385 r14786  
    8383      REAL(wp) ::   zyr_dec, zdco2dt 
    8484      CHARACTER (len=25) ::   charout 
    85       REAL(wp), DIMENSION(jpi,jpj) ::   zkgco2, zkgo2, zh2co3, zoflx,  zpco2atm   
     85      REAL(wp), DIMENSION(jpi,jpj) ::   zkgco2, zkgo2, zh2co3, zoflx,  zpco2atm, zpco2oce   
    8686      !!--------------------------------------------------------------------- 
    8787      ! 
     
    158158         zfld = zfco2 * chemc(ji,jj,1) * zkgco2(ji,jj)  ! (mol/L) * (m/s) 
    159159         zflu = zh2co3(ji,jj) * zkgco2(ji,jj)                                   ! (mol/L) (m/s) ? 
    160          oce_co2(ji,jj) = ( zfld - zflu ) * tmask(ji,jj,1)  
     160         zpco2oce(ji,jj) = zh2co3(ji,jj) / ( chemc(ji,jj,1) * zfugcoeff + rtrn ) 
     161         oce_co2(ji,jj)  = ( zfld - zflu ) * tmask(ji,jj,1)  
    161162         ! compute the trend 
    162163         tr(ji,jj,1,jpdic,Krhs) = tr(ji,jj,1,jpdic,Krhs) + oce_co2(ji,jj) * rfact2 / e3t(ji,jj,1,Kmm) 
     
    187188         CALL iom_put( "Oflx"    , zoflx(:,:) * 1000.  ) 
    188189         CALL iom_put( "Kg"      , zkgco2(:,:) * tmask(:,:,1)  ) 
    189          CALL iom_put( "Dpco2"   , ( zpco2atm(:,:) - zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) ) 
    190          CALL iom_put( "pCO2sea" , ( zh2co3(:,:) / ( chemc(:,:,1) + rtrn ) ) * tmask(:,:,1) ) 
     190         CALL iom_put( "Dpco2"   , ( zpco2atm(:,:) - zpco2oce(:,:) ) * tmask(:,:,1) ) 
     191         CALL iom_put( "pCO2sea" , zpco2oce(:,:) * tmask(:,:,1) ) 
    191192         CALL iom_put( "Dpo2"    , ( atcox * patm(:,:) - atcox * tr(:,:,1,jpoxy,Kbb) / ( chemo2(:,:,1) + rtrn ) ) * tmask(:,:,1) ) 
    192193         CALL iom_put( "tcflx"   , t_oce_co2_flx     )   ! molC/s 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zligand.F90

    r13295 r14786  
    2525   REAL(wp), PUBLIC ::  rlig     !: Remin ligand production 
    2626   REAL(wp), PUBLIC ::  prlgw    !: Photochemical of weak ligand 
     27   REAL(wp), PUBLIC ::  xklig    !: 1/2 saturation constant of photolysis 
    2728 
    2829   !! * Substitutions 
     
    6566         zlgwr = 1. / zlgwr * tgfunc(ji,jj,jk) * ( xstep / nyear_len(1) ) * blim(ji,jj,jk) * tr(ji,jj,jk,jplgw,Kbb) 
    6667         ! photochem loss of weak ligand 
    67          zlgwpr = prlgw * xstep * etot(ji,jj,jk) * tr(ji,jj,jk,jplgw,Kbb) * (1. - fr_i(ji,jj)) 
     68         zlgwpr = prlgw * xstep * etot(ji,jj,jk) * tr(ji,jj,jk,jplgw,Kbb)**3 * (1. - fr_i(ji,jj))   & 
     69         &        / ( tr(ji,jj,jk,jplgw,Kbb)**2 + (xklig)**2) 
    6870         tr(ji,jj,jk,jplgw,Krhs) = tr(ji,jj,jk,jplgw,Krhs) + zlgwp - zlgwr - zlgwpr 
    6971         zligrem(ji,jj,jk)   = zlgwr 
    7072         zligpr(ji,jj,jk)    = zlgwpr 
    71          zligprod(ji,jj,jk) = zlgwp 
     73         zligprod(ji,jj,jk)  = zlgwp 
    7274         ! 
    7375      END_3D 
     
    110112      INTEGER ::   ios   ! Local integer  
    111113      ! 
    112       NAMELIST/nampislig/ rlgw, prlgw, rlgs, rlig 
     114      NAMELIST/nampislig/ rlgw, prlgw, rlgs, rlig, xklig 
    113115      !!---------------------------------------------------------------------- 
    114116      ! 
     
    130132         WRITE(numout,*) '      Photolysis of weak ligand                    prlgw =', prlgw 
    131133         WRITE(numout,*) '      Lifetime (years) of strong ligands           rlgs  =', rlgs 
     134         WRITE(numout,*) '      1/2 saturation for photolysis                xklig =', xklig 
    132135      ENDIF 
    133136      ! 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zlim.F90

    r14416 r14786  
    4242   REAL(wp), PUBLIC ::  qnfelim     !:  optimal Fe quota for nanophyto 
    4343   REAL(wp), PUBLIC ::  qdfelim     !:  optimal Fe quota for diatoms 
     44   REAL(wp), PUBLIC ::  caco3r      !:  mean rainratio  
    4445 
    4546   !!* Phytoplankton limitation terms 
     
    6566   ! Coefficient for iron limitation following Flynn and Hipkin (1999) 
    6667   REAL(wp) ::  xcoef1   = 0.0016  / 55.85   
    67    REAL(wp) ::  xcoef2   = 1.21E-5 * 14. / 55.85 / 7.625 * 0.5 * 1.5 
    68    REAL(wp) ::  xcoef3   = 1.15E-4 * 14. / 55.85 / 7.625 * 0.5  
     68   REAL(wp) ::  xcoef2   = 1.21E-5 * 14. / 55.85 / 7.3125 * 0.5 * 1.5 
     69   REAL(wp) ::  xcoef3   = 1.15E-4 * 14. / 55.85 / 7.3125 * 0.5  
    6970 
    7071   !! * Substitutions 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zlys.F90

    r13295 r14786  
    6060      INTEGER  ::   ji, jj, jk, jn 
    6161      REAL(wp) ::   zdispot, zfact, zcalcon 
    62       REAL(wp) ::   zomegaca, zexcess, zexcess0 
     62      REAL(wp) ::   zomegaca, zexcess, zexcess0, zkd 
    6363      CHARACTER (len=25) ::   charout 
    6464      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zco3, zcaldiss, zhinit, zhi, zco3sat 
     
    6868      ! 
    6969      zhinit  (:,:,:) = hi(:,:,:) * 1000. / ( rhop(:,:,:) + rtrn ) 
     70      zco3    (:,:,:) = 0._wp     ;   zco3sat(:,:,:) = 0._wp 
     71      zcaldiss(:,:,:) = 0._wp 
    7072      ! 
    7173      !     ------------------------------------------- 
     
    99101         excess(ji,jj,jk) = 1._wp - zomegaca 
    100102         zexcess0 = MAX( 0., excess(ji,jj,jk) ) 
    101          zexcess  = zexcess0**nca 
    102103 
    103          ! AMOUNT CACO3 (12C) THAT RE-ENTERS SOLUTION 
    104          !       (ACCORDING TO THIS FORMULATION ALSO SOME PARTICULATE 
    105          !       CACO3 GETS DISSOLVED EVEN IN THE CASE OF OVERSATURATION) 
    106          zdispot = kdca * zexcess * tr(ji,jj,jk,jpcal,Kbb) 
     104         IF (zomegaca < 0.8) THEN 
     105            zexcess  = zexcess0**nca 
     106            ! AMOUNT CACO3 THAT RE-ENTERS SOLUTION 
     107            zdispot = kdca * zexcess * tr(ji,jj,jk,jpcal,Kbb) 
     108         ELSE 
     109 
     110            zkd = kdca * 0.1**(nca - 0.4) 
     111            zexcess  = zexcess0**0.4 
     112            zdispot = zkd * zexcess * tr(ji,jj,jk,jpcal,Kbb) 
     113        ENDIF 
     114 
    107115        !  CHANGE OF [CO3--] , [ALK], PARTICULATE [CACO3], 
    108116        !       AND [SUM(CO2)] DUE TO CACO3 DISSOLUTION/PRECIPITATION 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zmeso.F90

    r14385 r14786  
    141141         ! small cells 
    142142         ! ------------------------------------------------------------------------------- 
    143          zcompaph  = MAX( ( tr(ji,jj,jk,jpphy,Kbb) - xthresh2phy ), 0.e0 ) & 
    144             &      * MIN(1., MAX( 0., ( quotan(ji,jj,jk) - 0.2) / 0.3 ) ) 
    145  
     143         zcompaph  = MAX( ( tr(ji,jj,jk,jpphy,Kbb) - xthresh2phy ), 0.e0 ) 
    146144 
    147145         ! Mesozooplankton grazing 
     
    214212         ! Compute the proportion of filter feeders. It is assumed steady state. 
    215213         ! ---------------------------------------------------------------------   
    216          zproport  = (zgrazffep + zgrazffeg)/(rtrn + zgraztotc) 
     214         zproport  = 0._wp 
     215         IF( gdepw(ji,jj,jk+1,Kmm) > MAX(hmld(ji,jj), heup_01(ji,jj) ) ) THEN 
     216            zproport  = (zgrazffep + zgrazffeg)/(rtrn + zgraztotc) 
     217         ENDIF 
    217218 
    218219         ! Compute fractionation of aggregates. It is assumed that  
     
    228229         zfrac     = zproport * grazflux  * xstep * wsbio4(ji,jj,jk)      & 
    229230         &          * tr(ji,jj,jk,jpgoc,Kbb) * tr(ji,jj,jk,jpmes,Kbb)          & 
    230          &          * ( 0.2 + 3.8 * zratio2 / ( 1.**2 + zratio2 ) ) 
     231         &          * ( 0.4 + 3.6 * zratio2 / ( 1.**2 + zratio2 ) ) 
    231232         zfracfe   = zfrac * tr(ji,jj,jk,jpbfe,Kbb) / (tr(ji,jj,jk,jpgoc,Kbb) + rtrn) 
    232233 
     
    236237         zgrazfffp = zproport * zgrazfffp 
    237238         zgrazfffg = zproport * zgrazfffg 
     239         zgrazdc   = (1.0 - zproport) * zgrazdc 
     240         zgraznc   = (1.0 - zproport) * zgraznc 
     241         zgrazz    = (1.0 - zproport) * zgrazz 
     242         zgrazpoc  = (1.0 - zproport) * zgrazpoc 
     243         zgrazdf   = (1.0 - zproport) * zgrazdf 
     244         zgraznf   = (1.0 - zproport) * zgraznf 
     245         zgrazpof  = (1.0 - zproport) * zgrazpof 
     246 
    238247 
    239248         ! Total ingestion rates in C, N, Fe 
     
    280289         tr(ji,jj,jk,jpdch,Krhs) = tr(ji,jj,jk,jpdch,Krhs) - zgrazdc * tr(ji,jj,jk,jpdch,Kbb) / ( tr(ji,jj,jk,jpdia,Kbb) + rtrn ) 
    281290         tr(ji,jj,jk,jpdsi,Krhs) = tr(ji,jj,jk,jpdsi,Krhs) - zgrazdc * tr(ji,jj,jk,jpdsi,Kbb) / ( tr(ji,jj,jk,jpdia,Kbb) + rtrn ) 
    282          tr(ji,jj,jk,jpgsi,Krhs) = tr(ji,jj,jk,jpgsi,Krhs) + zgrazdc * tr(ji,jj,jk,jpdsi,Kbb) / ( tr(ji,jj,jk,jpdia,Kbb) + rtrn ) 
    283291         zgrabsi(ji,jj,jk)       = zgrazdc * tr(ji,jj,jk,jpdsi,Kbb) / ( tr(ji,jj,jk,jpdia,Kbb) + rtrn ) 
    284292         ! 
     
    308316         zprcaca = part2 * zprcaca 
    309317         tr(ji,jj,jk,jpdic,Krhs) = tr(ji,jj,jk,jpdic,Krhs) + zgrazcal - zprcaca 
    310          tr(ji,jj,jk,jptal,Krhs) = tr(ji,jj,jk,jptal,Krhs) - 2. * ( zgrazcal - zprcaca ) 
     318         tr(ji,jj,jk,jptal,Krhs) = tr(ji,jj,jk,jptal,Krhs) + 2. * ( zgrazcal - zprcaca ) 
    311319         tr(ji,jj,jk,jpcal,Krhs) = tr(ji,jj,jk,jpcal,Krhs) - zgrazcal + zprcaca 
    312320 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zmort.F90

    r14385 r14786  
    7676      prodcal(:,:,:) = 0._wp   ! calcite production variable set to zero 
    7777      DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
    78          zcompaph = MAX( ( tr(ji,jj,jk,jpphy,Kbb) - 1e-8 ), 0.e0 ) 
     78         zcompaph = MAX( ( tr(ji,jj,jk,jpphy,Kbb) - 1e-9 ), 0.e0 ) 
    7979 
    8080         ! Quadratic mortality of nano due to aggregation during 
     
    110110         tr(ji,jj,jk,jptal,Krhs) = tr(ji,jj,jk,jptal,Krhs) - 2. * zprcaca 
    111111         tr(ji,jj,jk,jpcal,Krhs) = tr(ji,jj,jk,jpcal,Krhs) + zprcaca 
    112          tr(ji,jj,jk,jpgoc,Krhs) = tr(ji,jj,jk,jpgoc,Krhs) + zfracal * zmortp 
    113          tr(ji,jj,jk,jppoc,Krhs) = tr(ji,jj,jk,jppoc,Krhs) + ( 1. - zfracal ) * zmortp 
    114          prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + ( 1. - zfracal ) * zmortp 
    115          prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zfracal * zmortp 
     112         tr(ji,jj,jk,jppoc,Krhs) = tr(ji,jj,jk,jppoc,Krhs) + zmortp 
     113         prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + zmortp 
    116114 
    117115         ! Update the arrays TRA which contains the biological sources and sinks 
    118          tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) + ( 1. - zfracal ) * zmortp * zfactfe 
    119          tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + zfracal * zmortp * zfactfe 
     116         tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) + zmortp * zfactfe 
    120117         ! 
    121118      END_3D 
     
    189186         ! Half of the linear mortality term is routed to big particles 
    190187         ! becaue of the ballasting effect 
    191          tr(ji,jj,jk,jpgoc,Krhs) = tr(ji,jj,jk,jpgoc,Krhs) + zrespp2 + 0.5 * ztortp2 
    192          tr(ji,jj,jk,jppoc,Krhs) = tr(ji,jj,jk,jppoc,Krhs) + 0.5 * ztortp2 
    193          prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + 0.5 * ztortp2 
    194          prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zrespp2 + 0.5 * ztortp2 
    195          tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) + 0.5 * ztortp2 * zfactfe 
    196          tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + ( zrespp2 + 0.5 * ztortp2 ) * zfactfe 
     188         tr(ji,jj,jk,jpgoc,Krhs) = tr(ji,jj,jk,jpgoc,Krhs) + zrespp2 
     189         tr(ji,jj,jk,jppoc,Krhs) = tr(ji,jj,jk,jppoc,Krhs) + ztortp2 
     190         prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + ztortp2 
     191         prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zrespp2 
     192         tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) + ztortp2 * zfactfe 
     193         tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + zrespp2 * zfactfe 
    197194      END_3D 
    198195      ! 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zopt.F90

    r14416 r14786  
    109109      ! -------------------------------------------- 
    110110      IF( l_trcdm2dc ) THEN                     !  diurnal cycle 
    111          ! 
    112           ! 
    113          ! SW over the ice free zone of the grid cell. This assumes that 
    114          ! SW is zero below sea ice which is a very crude assumption that is  
    115          ! not fully correct with LIM3 and SI3 but no information is  
    116          ! currently available to do a better job. SHould be improved in the  
    117          ! (near) future. 
    118          zqsr_corr(:,:) = qsr_mean(:,:) / ( 1.-fr_i(:,:) + rtrn ) 
    119          ! 
    120          CALL p4z_opt_par( kt, Kmm, zqsr_corr, ze1, ze2, ze3, pqsr100 = zqsr100 )  
    121          ! 
    122          ! Used PAR is computed for each phytoplankton species 
    123          ! etot_ndcy is PAR at level jk averaged over 24h. 
    124          ! Due to their size, they have different light absorption characteristics 
    125          DO jk = 1, nksr       
    126             etot_ndcy(:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk) 
    127             enano    (:,:,jk) =  1.85 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.46 * ze3(:,:,jk) 
    128             ediat    (:,:,jk) =  1.62 * ze1(:,:,jk) + 0.74 * ze2(:,:,jk) + 0.63 * ze3(:,:,jk) 
    129          END DO 
    130          IF( ln_p5z ) THEN 
     111         IF ( ln_p4z_dcyc ) THEN   ! Diurnal cycle in PISCES 
     112            ! 
     113            ! 
     114            ! SW over the ice free zone of the grid cell. This assumes that 
     115            ! SW is zero below sea ice which is a very crude assumption that is  
     116            ! not fully correct with LIM3 and SI3 but no information is  
     117            ! currently available to do a better job. SHould be improved in the  
     118            ! (near) future. 
     119            zqsr_corr(:,:) = qsr_mean(:,:) / ( 1.-fr_i(:,:) + rtrn ) 
     120            ! 
     121            CALL p4z_opt_par( kt, Kmm, zqsr_corr, ze1, ze2, ze3, pqsr100 = zqsr100 ) 
     122            ! 
     123            ! Used PAR is computed for each phytoplankton species 
     124            ! etot_ndcy is PAR at level jk averaged over 24h. 
     125            ! Due to their size, they have different light absorption characteristics 
     126            DO jk = 1, nksr 
     127               etot_ndcy(:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk) 
     128            END DO 
     129            ! 
     130            ! SW over the ice free zone of the grid cell. This assumes that 
     131            ! SW is zero below sea ice which is a very crude assumption that is  
     132            ! not fully correct with LIM3 and SI3 but no information is  
     133            ! currently available to do a better job. SHould be improved in the  
     134            ! (near) future. 
     135            zqsr_corr(:,:) = qsr(:,:) / ( 1.-fr_i(:,:) + rtrn ) 
     136            ! 
     137            CALL p4z_opt_par( kt, Kmm, zqsr_corr, ze1, ze2, ze3 ) 
     138            ! 
     139            ! Total PAR computation at level jk that includes the diurnal cycle 
     140            DO jk = 1, nksr 
     141               etot (:,:,jk) =  ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk) 
     142               enano(:,:,jk) =  1.85 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.46 * ze3(:,:,jk) 
     143               ediat(:,:,jk) =  1.62 * ze1(:,:,jk) + 0.74 * ze2(:,:,jk) + 0.63 * ze3(:,:,jk) 
     144            END DO 
     145            IF( ln_p5z ) THEN 
     146               DO jk = 1, nksr 
     147                  epico  (:,:,jk) =  1.94 * ze1(:,:,jk) + 0.66 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
     148               END DO 
     149            ENDIF 
     150 
     151         ELSE ! No diurnal cycle in PISCES 
     152 
     153            ! 
     154            ! 
     155            ! SW over the ice free zone of the grid cell. This assumes that 
     156            ! SW is zero below sea ice which is a very crude assumption that is  
     157            ! not fully correct with LIM3 and SI3 but no information is  
     158            ! currently available to do a better job. SHould be improved in the  
     159            ! (near) future. 
     160            zqsr_corr(:,:) = qsr_mean(:,:) / ( 1.-fr_i(:,:) + rtrn ) 
     161            ! 
     162            CALL p4z_opt_par( kt, Kmm, zqsr_corr, ze1, ze2, ze3, pqsr100 = zqsr100 )  
     163            ! 
     164            ! Used PAR is computed for each phytoplankton species 
     165            ! etot_ndcy is PAR at level jk averaged over 24h. 
     166            ! Due to their size, they have different light absorption characteristics 
    131167            DO jk = 1, nksr       
    132               epico  (:,:,jk) =  1.94 * ze1(:,:,jk) + 0.66 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
     168               etot_ndcy(:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk) 
     169               enano    (:,:,jk) =  1.85 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.46 * ze3(:,:,jk) 
     170               ediat    (:,:,jk) =  1.62 * ze1(:,:,jk) + 0.74 * ze2(:,:,jk) + 0.63 * ze3(:,:,jk) 
    133171            END DO 
    134          ENDIF 
    135          ! 
    136          ! SW over the ice free zone of the grid cell. This assumes that 
    137          ! SW is zero below sea ice which is a very crude assumption that is  
    138          ! not fully correct with LIM3 and SI3 but no information is  
    139          ! currently available to do a better job. SHould be improved in the  
    140          ! (near) future. 
    141  
    142          zqsr_corr(:,:) = qsr(:,:) / ( 1.-fr_i(:,:) + rtrn ) 
    143          ! 
    144          CALL p4z_opt_par( kt, Kmm, zqsr_corr, ze1, ze2, ze3 )  
    145          ! 
    146          ! Total PAR computation at level jk that includes the diurnal cycle 
    147          DO jk = 1, nksr       
    148             etot(:,:,jk) =  ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk) 
    149          END DO 
     172            IF( ln_p5z ) THEN 
     173               DO jk = 1, nksr       
     174                  epico  (:,:,jk) =  1.94 * ze1(:,:,jk) + 0.66 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
     175               END DO 
     176            ENDIF 
     177            ! 
     178            ! SW over the ice free zone of the grid cell. This assumes that 
     179            ! SW is zero below sea ice which is a very crude assumption that is  
     180            ! not fully correct with LIM3 and SI3 but no information is  
     181            ! currently available to do a better job. SHould be improved in the  
     182            ! (near) future. 
     183            zqsr_corr(:,:) = qsr(:,:) / ( 1.-fr_i(:,:) + rtrn ) 
     184            ! 
     185            CALL p4z_opt_par( kt, Kmm, zqsr_corr, ze1, ze2, ze3 )  
     186            ! 
     187            ! Total PAR computation at level jk that includes the diurnal cycle 
     188            DO jk = 1, nksr       
     189               etot(:,:,jk) =  ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk) 
     190            END DO 
     191         ENDIF 
    150192         ! 
    151193      ELSE   ! no diurnal cycle 
     
    401443      TYPE(FLD_N) ::   sn_par                ! informations about the fields to be read 
    402444      ! 
    403       NAMELIST/nampisopt/cn_dir, sn_par, ln_varpar, parlux 
     445      NAMELIST/nampisopt/cn_dir, sn_par, ln_varpar, parlux, ln_p4z_dcyc 
    404446      !!---------------------------------------------------------------------- 
    405447      IF(lwp) THEN 
     
    416458      IF(lwp) THEN 
    417459         WRITE(numout,*) '   Namelist : nampisopt ' 
    418          WRITE(numout,*) '      PAR as a variable fraction of SW     ln_varpar      = ', ln_varpar 
    419          WRITE(numout,*) '      Default value for the PAR fraction   parlux         = ', parlux 
     460         WRITE(numout,*) '      PAR as a variable fraction of SW       ln_varpar      = ', ln_varpar 
     461         WRITE(numout,*) '      Default value for the PAR fraction     parlux         = ', parlux 
     462         WRITE(numout,*) '      Activate the diurnal cycle in PISCES   ln_p4z_dcyc    = ', ln_p4z_dcyc 
    420463      ENDIF 
    421464      ! 
    422465      xparsw = parlux / 3.0 
    423466      xsi0r  = 1.e0 / rn_si0 
     467 
     468      ! Warning : activate the diurnal cycle with no diurnal cycle in the forcing fields makes no sense 
     469      ! That does not produce a bug because the model does not use the flag but a warning is necessary 
     470      ! ---------------------------------------------------------------------------------------------- 
     471      IF ( ln_p4z_dcyc .AND. l_trcdm2dc ) THEN 
     472         IF (lwp) WRITE(numout,*) 'No diurnal cycle in the PAR forcing field ' 
     473         IF (lwp) WRITE(numout,*) 'Activating the diurnal cycle in PISCES has no effect' 
     474      ENDIF 
    424475      ! 
    425476      ! Variable PAR at the surface of the ocean 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zprod.F90

    r14416 r14786  
    100100      ! of the thermal dependency 
    101101      ! Parameters are taken from Bissinger et al. (2008) 
    102       zprmaxn(:,:,:) = 0.8_wp * r1_rday * tgfunc(:,:,:) 
     102      zprmaxn(:,:,:) = 0.65_wp * r1_rday * tgfunc(:,:,:) 
    103103      zprmaxd(:,:,:) = zprmaxn(:,:,:) 
    104104 
     
    109109      ! absolute light level definition of the euphotic zone 
    110110      ! -------------------------------------------------------------------------  
    111       DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
    112          IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
    113             zval = MAX( 1., strn(ji,jj) ) 
    114             IF( gdepw(ji,jj,jk+1,Kmm) <= hmld(ji,jj) ) THEN 
    115                zval = zval * MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn )) 
     111      IF ( ln_p4z_dcyc ) THEN    ! Diurnal cycle in PISCES 
     112 
     113         DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
     114            IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
     115               IF( gdepw(ji,jj,jk+1,Kmm) <= hmld(ji,jj) ) THEN 
     116                  zval = MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn )) 
     117               ENDIF 
     118               zmxl_chl(ji,jj,jk) = zval / 24. 
     119               zmxl_fac(ji,jj,jk) = 1.0 - exp( -0.26 * zval ) 
    116120            ENDIF 
    117             zmxl_chl(ji,jj,jk) = zval / 24. 
    118             zmxl_fac(ji,jj,jk) = 1.0 - exp( -0.26 * zval ) 
    119          ENDIF 
    120       END_3D 
     121         END_3D 
     122  
     123      ELSE ! No diurnal cycle in PISCES 
     124 
     125         DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
     126            IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
     127               zval = MAX( 1., strn(ji,jj) ) 
     128               IF( gdepw(ji,jj,jk+1,Kmm) <= hmld(ji,jj) ) THEN 
     129                  zval = zval * MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn )) 
     130               ENDIF 
     131               zmxl_chl(ji,jj,jk) = zval / 24. 
     132               zmxl_fac(ji,jj,jk) = 1.0 - exp( -0.26 * zval ) 
     133            ENDIF 
     134         END_3D 
     135 
     136      ENDIF 
    121137 
    122138      zprbio(:,:,:) = zprmaxn(:,:,:) * zmxl_fac(:,:,:) 
     
    179195          zval = MIN( xnanopo4(ji,jj,jk), ( xnanonh4(ji,jj,jk) + xnanono3(ji,jj,jk) ) )   & 
    180196          &      * zprmaxn(ji,jj,jk) / ( zprbio(ji,jj,jk) + rtrn ) 
    181           quotan(ji,jj,jk) = MIN( 1., 0.2 + 0.8 * zval ) 
     197          quotan(ji,jj,jk) = MIN( 1., 0.3 + 0.7 * zval ) 
    182198          zval = MIN( xdiatpo4(ji,jj,jk), ( xdiatnh4(ji,jj,jk) + xdiatno3(ji,jj,jk) ) )   & 
    183199          &      * zprmaxd(ji,jj,jk) / ( zprdia(ji,jj,jk) + rtrn ) 
    184           quotad(ji,jj,jk) = MIN( 1., 0.2 + 0.8 * zval ) 
     200          quotad(ji,jj,jk) = MIN( 1., 0.3 + 0.7 * zval ) 
    185201      END_3D 
    186202 
     
    256272            &          + xnanonh4(ji,jj,jk) ) * (1. - xnanofer(ji,jj,jk) ) )   & 
    257273            &          * xnanofer(ji,jj,jk) * zmax * tr(ji,jj,jk,jpphy,Kbb) * rfact2 
    258  
    259274            ! production terms of diatoms (C) 
    260275            zprorcad(ji,jj,jk) = zprdia(ji,jj,jk) * xlimdia(ji,jj,jk) * tr(ji,jj,jk,jpdia,Kbb) * rfact2 
     
    337352           consfe3(ji,jj,jk)   = zprodfer * 75.0 / ( rtrn + ( plig(ji,jj,jk) + 75.0 * (1.0 - plig(ji,jj,jk) ) )   & 
    338353           &                   * tr(ji,jj,jk,jpfer,Kbb) ) / rfact2 
    339             
    340354           ! 
    341355           tr(ji,jj,jk,jpdic,Krhs) = tr(ji,jj,jk,jpdic,Krhs) - zpptot 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zrem.F90

    r14385 r14786  
    9797         zdep = MAX( hmld(ji,jj), heup_01(ji,jj) ) 
    9898         zdepbac(ji,jj,jk) = 0.6 * ( MAX(0.0, tr(ji,jj,jk,jpzoo,Kbb) + tr(ji,jj,jk,jpmes,Kbb) ) * 1.0E6 )**0.6 * 1.E-6 
    99          IF( gdept(ji,jj,jk,Kmm) < zdep ) THEN 
     99         IF( gdept(ji,jj,jk,Kmm) >= zdep ) THEN 
    100100            zdepmin = MIN( 1., zdep / gdept(ji,jj,jk,Kmm) ) 
    101101            zdepeff(ji,jj,jk) = zdepeff(ji,jj,jk) * zdepmin**0.3 
     
    114114         ! ----------------------------------------------------- 
    115115         zolimic = zremikc * ( 1.- nitrfac(ji,jj,jk) ) * tr(ji,jj,jk,jpdoc,Kbb)  
    116          zolimic = MIN( ( tr(ji,jj,jk,jpoxy,Kbb) - rtrn ) / o2ut, zolimic )  
     116         zolimic = MAX(0., MIN( ( tr(ji,jj,jk,jpoxy,Kbb) - rtrn ) / o2ut, zolimic ) )  
    117117         zolimi(ji,jj,jk) = zolimic 
    118118 
     
    155155         &         / ( 1.+ emoy(ji,jj,jk) ) * ( 1. + fr_i(ji,jj) * emoy(ji,jj,jk) )  
    156156         zdenitnh4 = nitrif * xstep * tr(ji,jj,jk,jpnh4,Kbb) * nitrfac(ji,jj,jk) 
    157          zdenitnh4 = MIN(  ( tr(ji,jj,jk,jpno3,Kbb) - rtrn ) / rdenita, zdenitnh4 )  
     157         zdenitnh4 = MAX(0., MIN(  ( tr(ji,jj,jk,jpno3,Kbb) - rtrn ) / rdenita, zdenitnh4 ) ) 
    158158         ! Update of the tracers trends 
    159159         ! ---------------------------- 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p4zsed.F90

    r14385 r14786  
    119119              zflx = (  tr(ji,jj,ikt,jpgoc,Kbb) * zwsbio4(ji,jj)   & 
    120120                &     + tr(ji,jj,ikt,jppoc,Kbb) * zwsbio3(ji,jj) ) * 1E6 
    121               zbureff(ji,jj) = 0.013 + 0.53 * zflx**2 / ( 7.0 + zflx )**2 
     121              zbureff(ji,jj) = 0.013 + 0.13 * zflx**2 / ( 7.0 + zflx )**2 
    122122           ENDIF 
    123123         END_2D 
     
    157157            tr(ji,jj,ikt,jpsil,Krhs) = tr(ji,jj,ikt,jpsil,Krhs) + zsiloss * zrivsil  
    158158            ! 
    159             zfactcal = MIN( excess(ji,jj,ikt), 0.2 ) 
    160             zfactcal = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) ) 
     159            zfactcal = MAX(-0.1, MIN( excess(ji,jj,ikt), 0.2 ) ) 
     160            zfactcal = 0.3 + 0.7 * MIN( 1., (0.1 + zfactcal) / ( 0.5 - zfactcal ) ) 
    161161            zrivalk  = sedcalfrac * zfactcal 
    162162            tr(ji,jj,ikt,jptal,Krhs) =  tr(ji,jj,ikt,jptal,Krhs) + zcaloss * zrivalk * 2.0 
     
    257257            ! Nitrogen fixation limitation by PO4 and Fe             
    258258            ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) ) 
    259             ztrpo4 = tr(ji,jj,jk,jppo4,Kbb) / ( 1E-6 + tr(ji,jj,jk,jppo4,Kbb) ) 
     259            ztrpo4 = tr(ji,jj,jk,jppo4,Kbb) / ( 1.E-6 + tr(ji,jj,jk,jppo4,Kbb) ) 
    260260            zlight =  ( 1.- EXP( -etot_ndcy(ji,jj,jk) / diazolight ) ) * ( 1. - fr_i(ji,jj) )  
    261261            nitrpot(ji,jj,jk) =  zmudia * r1_rday * zfact * MIN( ztrfer, ztrpo4 ) * zlight 
     
    271271            xdiano3 = tr(ji,jj,jk,jpno3,Kbb) / ( concnno3 + tr(ji,jj,jk,jpno3,Kbb) ) * (1. - xdianh4) 
    272272            zlim = ( 1.- xdiano3 - xdianh4 ) 
    273             IF( zlim <= 0.1 )   zlim = 0.01 
    274273            zfact = zlim * rfact2 
    275274            ! Nitrogen fixation limitation by PO4/DOP and Fe             
     
    300299            ! Fe/c of diazotrophs is assumed to be 30umol Fe/mol C at max             
    301300            zlight  =  ( 1.- EXP( -etot_ndcy(ji,jj,jk) / diazolight ) ) * ( 1. - fr_i(ji,jj) )  
    302             zsoufer = zlight * 2E-11 / ( 2E-11 + biron(ji,jj,jk) ) 
    303             tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) - 30E-6 * zfact * 1.0 / 3.0 
    304             tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) + 30E-6 * zfact * 1.0 / 3.0 * 2.0 / 3.0 
    305             tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + 30E-6 * zfact * 1.0 / 3.0 * 1.0 / 3.0 
    306             tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + 0.002 * 4E-10 * zsoufer * rfact2 / rday 
     301            zsoufer = zlight * 1E-10 / ( 1E-10 + biron(ji,jj,jk) ) 
     302            ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) ) 
     303            tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) - 30E-6 * ztrfer * zfact * 1.0 / 3.0 
     304            tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) + 30E-6 * ztrfer * zfact * 1.0 / 3.0 * 2.0 / 3.0 
     305            tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + 30E-6 * ztrfer * zfact * 1.0 / 3.0 * 1.0 / 3.0 
     306            tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + 0.005 * 4E-10 * zsoufer * rfact2 / rday 
    307307            tr(ji,jj,jk,jppo4,Krhs) = tr(ji,jj,jk,jppo4,Krhs) + concdnh4 / ( concdnh4 + tr(ji,jj,jk,jppo4,Kbb) ) & 
    308308            &                     * 0.001 * tr(ji,jj,jk,jpdoc,Kbb) * xstep 
     
    320320            ztrpo4 = tr(ji,jj,jk,jppo4,Kbb) / ( 1E-6 + tr(ji,jj,jk,jppo4,Kbb) ) 
    321321            ztrdop = tr(ji,jj,jk,jpdop,Kbb) / ( 1E-6 + tr(ji,jj,jk,jpdop,Kbb) ) * (1. - ztrpo4) 
     322            ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) ) 
    322323         
    323324            tr(ji,jj,jk,jppo4,Krhs) = tr(ji,jj,jk,jppo4,Krhs) - 16.0 / 46.0 * zfact * ( 1.0 - 1.0 / 3.0 ) & 
     
    336337            ! Fe/c of diazotrophs is assumed to be 30umol Fe/mol C at max            
    337338            zlight  =  ( 1.- EXP( -etot_ndcy(ji,jj,jk) / diazolight ) ) * ( 1. - fr_i(ji,jj) )  
    338             zsoufer = zlight * 2E-11 / ( 2E-11 + biron(ji,jj,jk) ) 
    339             tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) - 30E-6 * zfact * 1.0 / 3.0  
    340             tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) + 30E-6 * zfact * 1.0 / 3.0 * 2.0 / 3.0 
    341             tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + 30E-6 * zfact * 1.0 / 3.0 * 1.0 / 3.0 
     339            zsoufer = zlight * 1E-10 / ( 1E-10 + biron(ji,jj,jk) ) 
     340            tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) - 30E-6 * ztrfer * zfact * 1.0 / 3.0  
     341            tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) + 30E-6 * ztrfer * zfact * 1.0 / 3.0 * 2.0 / 3.0 
     342            tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + 30E-6 * ztrfer * zfact * 1.0 / 3.0 * 1.0 / 3.0 
    342343            tr(ji,jj,jk,jpfer,Krhs) = tr(ji,jj,jk,jpfer,Krhs) + 0.002 * 4E-10 * zsoufer * rfact2 / rday 
    343344         END_3D 
     
    401402      ! 
    402403      sedsilfrac = 0.03     ! percentage of silica loss in the sediments 
    403       sedcalfrac = 0.6      ! percentage of calcite loss in the sediments 
     404      sedcalfrac = 0.99      ! percentage of calcite loss in the sediments 
    404405      ! 
    405406      lk_sed = ln_sediment .AND. ln_sed_2way  
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p5zlim.F90

    r14416 r14786  
    137137      REAL(wp) ::   zqfemn, zqfemp, zqfemd, zbactno3, zbactnh4, zbiron 
    138138      REAL(wp) ::   znutlimtot, zlimno3, zlimnh4, zlim1f, zsizetmp 
    139       REAL(wp) ::   zrtp, zrtn 
    140139      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zrassn, zrassp, zrassd 
    141140      !!--------------------------------------------------------------------- 
     
    248247         ! Michaelis-Menten Limitation term for nutrients Small flagellates 
    249248         ! ----------------------------------------------- 
    250          zrtn    = tr(ji,jj,jk,jpnh4,Kbb) + tr(ji,jj,jk,jpno3,Kbb) 
    251          zrtp    = tr(ji,jj,jk,jppo4,Kbb) + tr(ji,jj,jk,jpdop,Kbb)  
     249         ztrn    = tr(ji,jj,jk,jpnh4,Kbb) + tr(ji,jj,jk,jpno3,Kbb) 
     250         ztrp    = tr(ji,jj,jk,jppo4,Kbb) + tr(ji,jj,jk,jpdop,Kbb)  
    252251         ! 
    253252         ! Limitation of N based nutrients uptake (NO3 and NH4) 
     
    255254         zlimnh4 = tr(ji,jj,jk,jpnh4,Kbb) / ( zconc0n + tr(ji,jj,jk,jpnh4,Kbb) ) 
    256255         zlimno3 = tr(ji,jj,jk,jpno3,Kbb) / ( zconc0n + tr(ji,jj,jk,jpno3,Kbb) ) 
    257          znutlimtot = (1. - fanano) * zrtn  / ( zfalim * zconc0n + ztrn ) 
     256         znutlimtot = (1. - fanano) * ztrn  / ( zfalim * zconc0n + ztrn ) 
    258257         xnanonh4(ji,jj,jk) = znutlimtot * 5.0 * zlimnh4 / ( zlimno3 + 5.0 * zlimnh4 + rtrn ) 
    259258         xnanono3(ji,jj,jk) = znutlimtot * zlimno3 / ( zlimno3 + 5.0 * zlimnh4 + rtrn ) 
     
    297296         xlimphys(ji,jj,jk) = MIN( 1., zlim1/( zlim1f + rtrn ), zlim3 ) 
    298297         xlimnpn (ji,jj,jk) = MIN( 1., zlim1) 
    299  
    300  
    301298         ! 
    302299         ! Michaelis-Menten Limitation term for nutrients picophytoplankton 
     
    306303         zlimnh4 = tr(ji,jj,jk,jpnh4,Kbb) / ( zconc0p + tr(ji,jj,jk,jpnh4,Kbb) ) 
    307304         zlimno3 = tr(ji,jj,jk,jpno3,Kbb) / ( zconc0p + tr(ji,jj,jk,jpno3,Kbb) ) 
    308          znutlimtot = (1. - fapico) * zrtn / ( zfalim * zconc0p + ztrn ) 
     305         znutlimtot = (1. - fapico) * ztrn / ( zfalim * zconc0p + ztrn ) 
    309306         xpiconh4(ji,jj,jk) = znutlimtot * 5.0 * zlimnh4 / ( zlimno3 + 5.0 * zlimnh4 + rtrn ) 
    310307         xpicono3(ji,jj,jk) = znutlimtot * zlimno3 / ( zlimno3 + 5.0 * zlimnh4 + rtrn ) 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p5zmeso.F90

    r14385 r14786  
    242242         ! Compute the proportion of filter feeders. It is assumed steady state. 
    243243         ! ---------------------------------------------------------------------   
    244          zproport  = (zgrazffep + zgrazffeg)/(rtrn + zgraztotc) 
     244         zproport  = 0._wp 
     245         IF( gdepw(ji,jj,jk+1,Kmm) > MAX(hmld(ji,jj), heup_01(ji,jj) ) ) THEN 
     246            zproport  = (zgrazffep + zgrazffeg)/(rtrn + zgraztotc) 
     247         ENDIF 
    245248 
    246249         !   Compute fractionation of aggregates. It is assumed that  
     
    252255         zfracc    = zproport * grazflux  * xstep * wsbio4(ji,jj,jk)      & 
    253256         &          * tr(ji,jj,jk,jpgoc,Kbb) * tr(ji,jj,jk,jpmes,Kbb)          & 
    254          &          * ( 0.2 + 3.8 * zratio2 / ( 1.**2 + zratio2 ) ) 
     257         &          * ( 0.4 + 3.6 * zratio2 / ( 1.**2 + zratio2 ) ) 
    255258         zfracfe   = zfracc * tr(ji,jj,jk,jpbfe,Kbb) / (tr(ji,jj,jk,jpgoc,Kbb) + rtrn) 
    256259         zfracn    = zfracc * tr(ji,jj,jk,jpgon,Kbb) / (tr(ji,jj,jk,jpgoc,Kbb) + rtrn) 
     
    258261 
    259262         ! Flux feeding is multiplied by the fractional biomass of flux feeders 
    260          zgrazffep = zproport * zgrazffep   ;   zgrazffeg = zproport * zgrazffeg 
    261          zgrazfffp = zproport * zgrazfffp   ;   zgrazfffg = zproport * zgrazfffg 
    262          zgrazffnp = zproport * zgrazffnp   ;   zgrazffng = zproport * zgrazffng 
    263          zgrazffpp = zproport * zgrazffpp   ;   zgrazffpg = zproport * zgrazffpg 
     263         zgrazffep = zproport * zgrazffep         ;   zgrazffeg = zproport * zgrazffeg 
     264         zgrazfffp = zproport * zgrazfffp         ;   zgrazfffg = zproport * zgrazfffg 
     265         zgrazffnp = zproport * zgrazffnp         ;   zgrazffng = zproport * zgrazffng 
     266         zgrazffpp = zproport * zgrazffpp         ;   zgrazffpg = zproport * zgrazffpg 
     267         zgrazdc   = (1.0 - zproport) * zgrazdc   ;   zgraznc   = (1.0 - zproport) * zgraznc 
     268         zgrazz    = (1.0 - zproport) * zgrazz    ;   zgrazpoc  = (1.0 - zproport) * zgrazpoc 
     269         zgrazm    = (1.0 - zproport) * zgrazm 
     270         zgrazdn   = (1.0 - zproport) * zgrazdn   ;   zgraznn   = (1.0 - zproport) * zgraznn 
     271         zgrazpon  = (1.0 - zproport) * zgrazpon 
     272         zgrazdp   = (1.0 - zproport) * zgrazdp   ;   zgraznp   = (1.0 - zproport) * zgraznp 
     273         zgrazpop  = (1.0 - zproport) * zgrazpop 
     274         zgrazdf   = (1.0 - zproport) * zgrazdf   ;   zgraznf   = (1.0 - zproport) * zgraznf 
     275         zgrazpof  = (1.0 - zproport) * zgrazpof 
    264276 
    265277         zgraztotc  = zgrazdc + zgrazz + zgraznc + zgrazm + zgrazpoc + zgrazffep + zgrazffeg 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/P4Z/p5zprod.F90

    r14416 r14786  
    130130      ! absolute light level definition of the euphotic zone 
    131131      ! -------------------------------------------------------------------------  
    132       DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
    133          IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
    134             zval = MAX( 1., strn(ji,jj) ) 
    135             IF( gdepw(ji,jj,jk+1,Kmm) <= hmld(ji,jj) ) THEN 
    136                zval = zval * MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn )) 
     132 
     133      IF ( ln_p4z_dcyc ) THEN    ! Diurnal cycle in PISCES 
     134 
     135         DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
     136            IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
     137               IF( gdepw(ji,jj,jk+1,Kmm) <= hmld(ji,jj) ) THEN 
     138                  zval = MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn )) 
     139               ENDIF 
     140               zmxl_chl(ji,jj,jk) = zval / 24. 
     141               zmxl_fac(ji,jj,jk) = 1.0 - exp( -0.26 * zval ) 
    137142            ENDIF 
    138             zmxl_chl(ji,jj,jk) = zval / 24. 
    139             zmxl_fac(ji,jj,jk) = 1.0 - exp( -0.26 * zval ) 
    140          ENDIF 
    141       END_3D 
     143         END_3D 
     144 
     145      ELSE ! No diurnal cycle in PISCES 
     146 
     147         DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 
     148            IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
     149               zval = MAX( 1., strn(ji,jj) ) 
     150               IF( gdepw(ji,jj,jk+1,Kmm) <= hmld(ji,jj) ) THEN 
     151                  zval = zval * MIN(1., heup_01(ji,jj) / ( hmld(ji,jj) + rtrn )) 
     152               ENDIF 
     153               zmxl_chl(ji,jj,jk) = zval / 24. 
     154               zmxl_fac(ji,jj,jk) = 1.0 - exp( -0.26 * zval ) 
     155            ENDIF 
     156         END_3D 
     157 
     158      ENDIF 
    142159 
    143160      zprbio(:,:,:) = zprmaxn(:,:,:) * zmxl_fac(:,:,:) 
     
    433450        zprodfer = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk) + texcretp * zprofep(ji,jj,jk) 
    434451        zresptot = zrespn(ji,jj,jk) + zrespp(ji,jj,jk) + zrespd(ji,jj,jk)  
    435   
    436452        ! 
    437453        tr(ji,jj,jk,jppo4,Krhs) = tr(ji,jj,jk,jppo4,Krhs) - zpptot 
     
    470486        tr(ji,jj,jk,jppdi,Krhs) = tr(ji,jj,jk,jppdi,Krhs) + ( zpropo4d(ji,jj,jk) + zprodopd(ji,jj,jk) ) * texcretd 
    471487        tr(ji,jj,jk,jpdfe,Krhs) = tr(ji,jj,jk,jpdfe,Krhs) + zprofed(ji,jj,jk) * texcretd 
    472         tr(ji,jj,jk,jpdsi,Krhs) = tr(ji,jj,jk,jpdsi,Krhs) + zprorcad(ji,jj,jk) * zysopt(ji,jj,jk) * texcretd 
    473  
     488        tr(ji,jj,jk,jpdsi,Krhs) = tr(ji,jj,jk,jpdsi,Krhs) + zprmaxd(ji,jj,jk) * zysopt(ji,jj,jk) * rfact2 * tr(ji,jj,jk,jpdia,Kbb) 
    474489        tr(ji,jj,jk,jpdoc,Krhs) = tr(ji,jj,jk,jpdoc,Krhs) + zproddoc 
    475490        tr(ji,jj,jk,jpdon,Krhs) = tr(ji,jj,jk,jpdon,Krhs) + zproddon                                         
     
    482497        consfe3(ji,jj,jk)       = zprodfer * 75.0 / ( rtrn + ( plig(ji,jj,jk) + 75.0 * (1.0 - plig(ji,jj,jk) ) )   & 
    483498           &                   * tr(ji,jj,jk,jpfer,Kbb) ) / rfact2 
    484          
    485         tr(ji,jj,jk,jpsil,Krhs) = tr(ji,jj,jk,jpsil,Krhs) - texcretd * zprorcad(ji,jj,jk) * zysopt(ji,jj,jk)  
     499        tr(ji,jj,jk,jpsil,Krhs) = tr(ji,jj,jk,jpsil,Krhs) - zprmaxd(ji,jj,jk) * zysopt(ji,jj,jk) * rfact2 * tr(ji,jj,jk,jpdia,Kbb)  
    486500 
    487501        tr(ji,jj,jk,jpdic,Krhs) = tr(ji,jj,jk,jpdic,Krhs) - zpptot  & 
     
    522536       CALL iom_put( "PPNEWN"  , zpronewn(:,:,:) * zfact * tmask(:,:,:)    ) ! new primary production by nanophyto 
    523537       CALL iom_put( "PPNEWD"  , zpronewd(:,:,:) * zfact * tmask(:,:,:)   ) ! new primary production by diatomes 
    524        CALL iom_put( "PBSi"    , zprorcad(:,:,:) * zfact * tmask(:,:,:) * zysopt(:,:,:)  ) ! biogenic silica production 
    525        CALL iom_put( "PFeP"    , zprofep(:,:,:) * zfact * tmask(:,:,:)  ) ! biogenic iron production by picophyto 
     538       CALL iom_put( "PBSi"    , zprmaxd (:,:,:) * zfact * tmask(:,:,:) * zysopt(:,:,:)  ) ! biogenic silica production 
     539       CALL iom_put( "PFeP"    , zprofep (:,:,:) * zfact * tmask(:,:,:)  ) ! biogenic iron production by picophyto 
    526540       CALL iom_put( "PFeN"    , zprofen(:,:,:) * zfact * tmask(:,:,:)  ) ! biogenic iron production by nanophyto 
    527541       CALL iom_put( "PFeD"    , zprofed(:,:,:) * zfact * tmask(:,:,:)  ) ! biogenic iron production by  diatomes 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/SED/sedarr.F90

    r14385 r14786  
    5151         jid        = MOD( tab_ind(jn) - 1, jpi ) + 1 
    5252         jjd        = ( tab_ind(jn) - 1 ) / jpi + 1 
    53 !         IF ( mig(jid) == 112 .and. mjg(jjd) == 25) write(0,*) 'plante ',jn,ndim1d 
    5453         tab1d(jn)  = tab2d(jid, jjd) 
    5554      END DO  
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/sms_pisces.F90

    r14416 r14786  
    4242   REAL(wp) ::   wsbio2max        !: Maximum sinking speed of the largest particles 
    4343   REAL(wp) ::   wsbio2scale      !: Length scale for the variations of wsbio2 
     44   REAL(wp) ::   oxymin           !:  half saturation constant for anoxia 
    4445   REAL(wp) ::   xkmort           !: Mortality half-saturation constant 
    4546   REAL(wp) ::   feratz           !: Fe/C in microzooplankton 
     
    5051   REAL(wp) ::   no3rat3          !: C/N ratio of zooplankton 
    5152   REAL(wp) ::   po4rat3          !: C/P ratio of zooplankton 
    52    REAL(wp) ::  oxymin            !:  half saturation constant for anoxia 
    53    REAL(wp) ::  caco3r            !:  mean rainratio 
    5453 
    5554   !!*  diagnostic parameters  
     
    6766   LOGICAL  ::  ln_check_mass     !: Flag to check mass conservation 
    6867   LOGICAL, PUBLIC ::   ln_ironice   !: boolean for Fe input from sea ice 
     68 
     69   !!* Diurnal cycle in PISCES 
     70   LOGICAL  ::  ln_p4z_dcyc       !: Flag to activate diurnal cycle in PISCES 
    6971 
    7072   !!*  Biological fluxes for light : variables shared by pisces & lobster 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/PISCES/trcini_pisces.F90

    r14416 r14786  
    107107      bioma0 =  1.000e-8_wp   
    108108      silic1 =  91.51e-6_wp   
    109       no3    =  30.9e-6_wp * 7.625_wp 
     109      no3    =  30.9e-6_wp * 7.3125_wp 
    110110      ! 
    111111      ! Allocate PISCES arrays 
     
    187187      ! Set biological ratios 
    188188      ! --------------------- 
    189       rno3    =  16._wp / 122._wp   ! C/N 
    190       po4r    =   1._wp / 122._wp   ! C/P 
    191       o2nit   =  32._wp / 122._wp   ! O2/C for nitrification 
    192       o2ut    = 133._wp / 122._wp   ! O2/C for ammonification 
     189      rno3    =  16._wp / 117._wp   ! C/N 
     190      po4r    =   1._wp / 117._wp   ! C/P 
     191      o2nit   =  32._wp / 117._wp   ! O2/C for nitrification 
     192      o2ut    = 138._wp / 117._wp   ! O2/C for ammonification 
    193193      rdenit  =  ( ( o2ut + o2nit ) * 0.80 - rno3 - rno3 * 0.60 ) / rno3  ! Denitrification 
    194194      rdenita =   3._wp /  5._wp    ! Denitrification 
  • NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/TOP/trcini.F90

    r14227 r14786  
    244244      ! 
    245245      IF( ln_trcdta )   CALL trc_dta_ini( jptra )           ! set initial tracers values 
    246       ! 
    247       IF( ln_trcbc .AND. lltrcbc )  THEN  
    248         CALL trc_bc_ini ( jptra, Kmm  )            ! set tracers Boundary Conditions 
    249         CALL trc_bc     ( nit000, Kmm, tr, Kaa )   ! tracers: surface and lateral Boundary Conditions 
    250       ENDIF 
    251       ! 
    252       IF( ln_trcais ) CALL trc_ais_ini   ! set tracers from Antarctic Ice Sheet 
    253246      ! 
    254247      IF( ln_rsttr ) THEN              ! restart from a file 
     
    273266      ! 
    274267      tr(:,:,:,:,Kaa) = 0._wp 
     268      ! 
     269      IF( ln_trcbc .AND. lltrcbc )  THEN 
     270        CALL trc_bc_ini ( jptra, Kmm  )            ! set tracers Boundary Conditions 
     271        CALL trc_bc     ( nit000, Kmm, tr, Kaa )   ! tracers: surface and lateral Boundary Conditions 
     272      ENDIF 
     273      ! 
     274      IF( ln_trcais ) CALL trc_ais_ini   ! set tracers from Antarctic Ice Sheet 
    275275      !                                                         ! Partial top/bottom cell: GRADh(tr(Kmm)) 
    276276   END SUBROUTINE trc_ini_state 
Note: See TracChangeset for help on using the changeset viewer.