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 5348 for branches/2015/dev_r5204_CNRS_PISCES_dcy/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsed.F90 – NEMO

Ignore:
Timestamp:
2015-06-04T13:10:38+02:00 (9 years ago)
Author:
cetlod
Message:

dev_r5204_CNRS_PISCES_dcy : Add a leapfrog temporal scheme in PISCES model

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5204_CNRS_PISCES_dcy/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsed.F90

    r5206 r5348  
    2121   USE p4zopt          !  optical model 
    2222   USE p4zlim          !  Co-limitations of differents nutrients 
    23    USE p4zrem          !  Remineralisation of organic matter 
    2423   USE p4zsbc          !  External source of nutrients  
    2524   USE p4zint          !  interpolation and computation of various fields 
     
    3029   PRIVATE 
    3130 
    32    PUBLIC   p4z_sed    
     31   PUBLIC   p4z_sed   
     32   PUBLIC   p4z_sed_alloc 
     33  
    3334 
    3435   !! * Module variables 
     36   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: nitrpot    !: Nitrogen fixation  
     37   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:  ) :: sdenit     !: Nitrate reduction in the sediments 
    3538   REAL(wp) :: r1_rday                  !: inverse of rday 
    36  
    37    INTEGER ::  numnit   
    38  
    3939 
    4040   !!* Substitution 
     
    4747CONTAINS 
    4848 
    49    SUBROUTINE p4z_sed( kt, jnt ) 
     49   SUBROUTINE p4z_sed( kt, knt ) 
    5050      !!--------------------------------------------------------------------- 
    5151      !!                     ***  ROUTINE p4z_sed  *** 
     
    5858      !!--------------------------------------------------------------------- 
    5959      ! 
    60       INTEGER, INTENT(in) ::   kt, jnt ! ocean time step 
     60      INTEGER, INTENT(in) ::   kt, knt ! ocean time step 
    6161      INTEGER  ::   ji, jj, jk, ikt 
    6262#if ! defined key_sed 
     
    6969      REAL(wp) ::  zsiloss, zcaloss, zws3, zws4, zwsc, zdep, zwstpoc 
    7070      REAL(wp) ::  ztrfer, ztrpo4, zwdust, zlight 
    71       REAL(wp) ::  zrdenittot, zsdenittot, znitrpottot 
    7271      ! 
    7372      CHARACTER (len=25) :: charout 
    74       REAL(wp), POINTER, DIMENSION(:,:  ) :: zpdep, zsidep, zwork1, zwork2, zwork3, zwork4 
     73      REAL(wp), POINTER, DIMENSION(:,:  ) :: zpdep, zsidep, zwork1, zwork2, zwork3 
    7574      REAL(wp), POINTER, DIMENSION(:,:  ) :: zdenit2d, zironice, zbureff 
    7675      REAL(wp), POINTER, DIMENSION(:,:  ) :: zwsbio3, zwsbio4, zwscal 
    77       REAL(wp), POINTER, DIMENSION(:,:,:) :: znitrpot, zirondep, zsoufer 
     76      REAL(wp), POINTER, DIMENSION(:,:,:) :: zirondep, zsoufer 
    7877      !!--------------------------------------------------------------------- 
    7978      ! 
    8079      IF( nn_timing == 1 )  CALL timing_start('p4z_sed') 
    8180      ! 
    82       IF( kt == nittrc000 .AND. jnt == 1 )  THEN 
    83          r1_rday  = 1. / rday 
    84          IF( ln_check_mass .AND. lwp)  & 
    85            &  CALL ctl_opn( numnit, 'nitrogen.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea ) 
    86       ENDIF 
     81      IF( kt == nittrc000 .AND. knt == 1 )   r1_rday  = 1. / rday 
    8782      ! 
    8883      ! Allocate temporary workspace 
    89       CALL wrk_alloc( jpi, jpj, zdenit2d, zwork1, zwork2, zwork3, zwork4, zbureff ) 
     84      CALL wrk_alloc( jpi, jpj, zdenit2d, zwork1, zwork2, zwork3, zbureff ) 
    9085      CALL wrk_alloc( jpi, jpj, zwsbio3, zwsbio4, zwscal ) 
    91       CALL wrk_alloc( jpi, jpj, jpk, znitrpot, zsoufer ) 
     86      CALL wrk_alloc( jpi, jpj, jpk, zsoufer ) 
    9287 
    9388      zdenit2d(:,:) = 0.e0 
     
    9691      zwork2  (:,:) = 0.e0 
    9792      zwork3  (:,:) = 0.e0 
    98       zwork4  (:,:) = 0.e0 
    9993 
    10094      ! Iron input/uptake due to sea ice : Crude parameterization based on Lancelot et al. 
     
    108102               zdep    = rfact2 / fse3t(ji,jj,1) 
    109103               zwflux  = fmmflx(ji,jj) / 1000._wp 
    110                zfminus = MIN( 0._wp, -zwflux ) * trn(ji,jj,1,jpfer) * zdep 
     104               zfminus = MIN( 0._wp, -zwflux ) * trb(ji,jj,1,jpfer) * zdep 
    111105               zfplus  = MAX( 0._wp, -zwflux ) * icefeinput * zdep 
    112106               zironice(ji,jj) =  zfplus + zfminus 
     
    114108         END DO 
    115109         ! 
    116          trn(:,:,1,jpfer) = trn(:,:,1,jpfer) + zironice(:,:)  
     110         tra(:,:,1,jpfer) = tra(:,:,1,jpfer) + zironice(:,:)  
    117111         !  
    118          IF( lk_iomput .AND. jnt == nrdttrc .AND. iom_use( "Ironice" ) )   & 
     112         IF( lk_iomput .AND. knt == nrdttrc .AND. iom_use( "Ironice" ) )   & 
    119113            &   CALL iom_put( "Ironice", zironice(:,:) * 1.e+3 * rfact2r * fse3t(:,:,1) * tmask(:,:,1) ) ! iron flux from ice 
    120114         ! 
     
    144138         END DO 
    145139         !                                              ! Iron solubilization of particles in the water column 
    146          trn(:,:,1,jppo4) = trn(:,:,1,jppo4) + zpdep   (:,:) 
    147          trn(:,:,1,jpsil) = trn(:,:,1,jpsil) + zsidep  (:,:) 
    148          trn(:,:,:,jpfer) = trn(:,:,:,jpfer) + zirondep(:,:,:)  
     140         tra(:,:,1,jppo4) = tra(:,:,1,jppo4) + zpdep   (:,:) 
     141         tra(:,:,1,jpsil) = tra(:,:,1,jpsil) + zsidep  (:,:) 
     142         tra(:,:,:,jpfer) = tra(:,:,:,jpfer) + zirondep(:,:,:)  
    149143         !  
    150144         IF( lk_iomput ) THEN 
    151             IF( jnt == nrdttrc ) THEN 
     145            IF( knt == nrdttrc ) THEN 
    152146                IF( iom_use( "Irondep" ) )   & 
    153147                &  CALL iom_put( "Irondep", zirondep(:,:,1) * 1.e+3 * rfact2r * fse3t(:,:,1) * tmask(:,:,1) ) ! surface downward dust depo of iron 
     
    170164            DO ji = 1, jpi 
    171165               DO jk = 1, nk_rnf(ji,jj) 
    172                   trn(ji,jj,jk,jppo4) = trn(ji,jj,jk,jppo4) +  rivdip(ji,jj) * rfact2 
    173                   trn(ji,jj,jk,jpno3) = trn(ji,jj,jk,jpno3) +  rivdin(ji,jj) * rfact2 
    174                   trn(ji,jj,jk,jpfer) = trn(ji,jj,jk,jpfer) +  rivdic(ji,jj) * 5.e-5 * rfact2 
    175                   trn(ji,jj,jk,jpsil) = trn(ji,jj,jk,jpsil) +  rivdsi(ji,jj) * rfact2 
    176                   trn(ji,jj,jk,jpdic) = trn(ji,jj,jk,jpdic) +  rivdic(ji,jj) * rfact2 
    177                   trn(ji,jj,jk,jptal) = trn(ji,jj,jk,jptal) +  ( rivalk(ji,jj) - rno3 * rivdin(ji,jj) ) * rfact2 
     166                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) +  rivdip(ji,jj) * rfact2 
     167                  tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) +  rivdin(ji,jj) * rfact2 
     168                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) +  rivdic(ji,jj) * 5.e-5 * rfact2 
     169                  tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) +  rivdsi(ji,jj) * rfact2 
     170                  tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) +  rivdic(ji,jj) * rfact2 
     171                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) +  ( rivalk(ji,jj) - rno3 * rivdin(ji,jj) ) * rfact2 
    178172               ENDDO 
    179173            ENDDO 
     
    184178      ! ---------------------------------------------------------- 
    185179      IF( ln_ndepo ) THEN 
    186          trn(:,:,1,jpno3) = trn(:,:,1,jpno3) + nitdep(:,:) * rfact2 
    187          trn(:,:,1,jptal) = trn(:,:,1,jptal) - rno3 * nitdep(:,:) * rfact2 
     180         tra(:,:,1,jpno3) = tra(:,:,1,jpno3) + nitdep(:,:) * rfact2 
     181         tra(:,:,1,jptal) = tra(:,:,1,jptal) - rno3 * nitdep(:,:) * rfact2 
    188182      ENDIF 
    189183 
     
    191185      ! ------------------------------------------------------ 
    192186      IF( ln_ironsed ) THEN 
    193          trn(:,:,:,jpfer) = trn(:,:,:,jpfer) + ironsed(:,:,:) * rfact2 
     187         tra(:,:,:,jpfer) = tra(:,:,:,jpfer) + ironsed(:,:,:) * rfact2 
    194188         ! 
    195          IF( lk_iomput .AND. jnt == nrdttrc .AND. iom_use( "Ironsed" ) )   & 
     189         IF( lk_iomput .AND. knt == nrdttrc .AND. iom_use( "Ironsed" ) )   & 
    196190            &   CALL iom_put( "Ironsed", ironsed(:,:,:) * 1.e+3 * tmask(:,:,:) ) ! iron inputs from sediments 
    197191      ENDIF 
     
    200194      ! ------------------------------------------------------ 
    201195      IF( ln_hydrofe ) THEN 
    202          trn(:,:,:,jpfer) = trn(:,:,:,jpfer) + hydrofe(:,:,:) * rfact2 
     196         tra(:,:,:,jpfer) = tra(:,:,:,jpfer) + hydrofe(:,:,:) * rfact2 
    203197         ! 
    204          IF( lk_iomput .AND. jnt == nrdttrc .AND. iom_use( "HYDR" ) )   & 
     198         IF( lk_iomput .AND. knt == nrdttrc .AND. iom_use( "HYDR" ) )   & 
    205199            &   CALL iom_put( "HYDR", hydrofe(:,:,:) * 1.e+3 * tmask(:,:,:) ) ! hydrothermal iron input 
    206200      ENDIF 
     
    228222              ikt = mbkt(ji,jj) 
    229223# if defined key_kriest 
    230               zflx =    trn(ji,jj,ikt,jppoc) * zwsbio3(ji,jj)    * 1E3 * 1E6 / 1E4 
     224              zflx =    trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj)    * 1E3 * 1E6 / 1E4 
    231225# else 
    232               zflx = (  trn(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj)   & 
    233                 &     + trn(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) )  * 1E3 * 1E6 / 1E4 
     226              zflx = (  trb(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj)   & 
     227                &     + trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) )  * 1E3 * 1E6 / 1E4 
    234228#endif 
    235229              zflx  = LOG10( MAX( 1E-3, zflx ) ) 
    236               zo2   = LOG10( MAX( 10. , trn(ji,jj,ikt,jpoxy) * 1E6 ) ) 
    237               zno3  = LOG10( MAX( 1.  , trn(ji,jj,ikt,jpno3) * 1E6 * rno3 ) ) 
     230              zo2   = LOG10( MAX( 10. , trb(ji,jj,ikt,jpoxy) * 1E6 ) ) 
     231              zno3  = LOG10( MAX( 1.  , trb(ji,jj,ikt,jpno3) * 1E6 * rno3 ) ) 
    238232              zdep  = LOG10( fsdepw(ji,jj,ikt+1) ) 
    239233              zdenit2d(ji,jj) = -2.2567 - 1.185 * zflx - 0.221 * zflx**2 - 0.3995 * zno3 * zo2 + 1.25 * zno3    & 
     
    241235              zdenit2d(ji,jj) = 10.0**( zdenit2d(ji,jj) ) 
    242236              ! 
    243               zflx = (  trn(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj)   & 
    244                 &     + trn(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) ) * 1E6 
     237              zflx = (  trb(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj)   & 
     238                &     + trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) ) * 1E6 
    245239              zbureff(ji,jj) = 0.013 + 0.53 * zflx**2 / ( 7.0 + zflx )**2 
    246240           ENDIF 
     
    257251               ikt = mbkt(ji,jj)  
    258252# if defined key_kriest 
    259                zwork1(ji,jj) = trn(ji,jj,ikt,jpgsi) * zwscal (ji,jj) 
    260                zwork2(ji,jj) = trn(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) 
     253               zwork1(ji,jj) = trb(ji,jj,ikt,jpgsi) * zwscal (ji,jj) 
     254               zwork2(ji,jj) = trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) 
    261255# else 
    262                zwork1(ji,jj) = trn(ji,jj,ikt,jpgsi) * zwsbio4(ji,jj) 
    263                zwork2(ji,jj) = trn(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj) + trn(ji,jj,ikt,jppoc) * zwsbio3(ji,jj)  
     256               zwork1(ji,jj) = trb(ji,jj,ikt,jpgsi) * zwsbio4(ji,jj) 
     257               zwork2(ji,jj) = trb(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj) + trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj)  
    264258# endif 
    265259               ! For calcite, burial efficiency is made a function of saturation 
    266260               zfactcal      = MIN( excess(ji,jj,ikt), 0.2 ) 
    267261               zfactcal      = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) ) 
    268                zwork3(ji,jj) = trn(ji,jj,ikt,jpcal) * zwscal(ji,jj) * 2.e0 * zfactcal 
     262               zwork3(ji,jj) = trb(ji,jj,ikt,jpcal) * zwscal(ji,jj) * 2.e0 * zfactcal 
    269263            ENDIF 
    270264         END DO 
     
    285279         DO ji = 1, jpi 
    286280            ikt  = mbkt(ji,jj) 
    287             zdep = xstep / fse3t(ji,jj,ikt) 
     281            zdep = xstep / fse3t(ji,jj,ikt)  
    288282            zws4 = zwsbio4(ji,jj) * zdep 
    289283            zwsc = zwscal (ji,jj) * zdep 
    290284# if defined key_kriest 
    291             zsiloss = trn(ji,jj,ikt,jpgsi) * zws4 
     285            zsiloss = trb(ji,jj,ikt,jpgsi) * zws4 
    292286# else 
    293             zsiloss = trn(ji,jj,ikt,jpgsi) * zwsc 
     287            zsiloss = trb(ji,jj,ikt,jpgsi) * zwsc 
    294288# endif 
    295             zcaloss = trn(ji,jj,ikt,jpcal) * zwsc 
     289            zcaloss = trb(ji,jj,ikt,jpcal) * zwsc 
    296290            ! 
    297             trn(ji,jj,ikt,jpgsi) = trn(ji,jj,ikt,jpgsi) - zsiloss 
    298             trn(ji,jj,ikt,jpcal) = trn(ji,jj,ikt,jpcal) - zcaloss 
     291            tra(ji,jj,ikt,jpgsi) = tra(ji,jj,ikt,jpgsi) - zsiloss 
     292            tra(ji,jj,ikt,jpcal) = tra(ji,jj,ikt,jpcal) - zcaloss 
    299293#if ! defined key_sed 
    300             trn(ji,jj,ikt,jpsil) = trn(ji,jj,ikt,jpsil) + zsiloss * zrivsil  
     294            tra(ji,jj,ikt,jpsil) = tra(ji,jj,ikt,jpsil) + zsiloss * zrivsil  
    301295            zfactcal = MIN( excess(ji,jj,ikt), 0.2 ) 
    302296            zfactcal = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) ) 
    303297            zrivalk  =  1._wp - ( rivalkinput * r1_ryyss ) * zfactcal / ( zsumsedcal + rtrn ) 
    304             trn(ji,jj,ikt,jptal) =  trn(ji,jj,ikt,jptal) + zcaloss * zrivalk * 2.0 
    305             trn(ji,jj,ikt,jpdic) =  trn(ji,jj,ikt,jpdic) + zcaloss * zrivalk 
     298            tra(ji,jj,ikt,jptal) =  tra(ji,jj,ikt,jptal) + zcaloss * zrivalk * 2.0 
     299            tra(ji,jj,ikt,jpdic) =  tra(ji,jj,ikt,jpdic) + zcaloss * zrivalk 
    306300#endif 
    307301         END DO 
     
    310304      DO jj = 1, jpj 
    311305         DO ji = 1, jpi 
    312             ikt     = mbkt(ji,jj) 
    313             zdep    = xstep / fse3t(ji,jj,ikt) 
     306            ikt  = mbkt(ji,jj) 
     307            zdep = xstep / fse3t(ji,jj,ikt)  
    314308            zws4 = zwsbio4(ji,jj) * zdep 
    315309            zws3 = zwsbio3(ji,jj) * zdep 
    316310            zrivno3 = 1. - zbureff(ji,jj) 
    317311# if ! defined key_kriest 
    318             trn(ji,jj,ikt,jpgoc) = trn(ji,jj,ikt,jpgoc) - trn(ji,jj,ikt,jpgoc) * zws4 
    319             trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) - trn(ji,jj,ikt,jppoc) * zws3 
    320             trn(ji,jj,ikt,jpbfe) = trn(ji,jj,ikt,jpbfe) - trn(ji,jj,ikt,jpbfe) * zws4 
    321             trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) - trn(ji,jj,ikt,jpsfe) * zws3 
    322             zwstpoc              =  trn(ji,jj,ikt,jpgoc) * zws4 + trn(ji,jj,ikt,jppoc) * zws3  
     312            tra(ji,jj,ikt,jpgoc) = tra(ji,jj,ikt,jpgoc) - trb(ji,jj,ikt,jpgoc) * zws4  
     313            tra(ji,jj,ikt,jppoc) = tra(ji,jj,ikt,jppoc) - trb(ji,jj,ikt,jppoc) * zws3 
     314            tra(ji,jj,ikt,jpbfe) = tra(ji,jj,ikt,jpbfe) - trb(ji,jj,ikt,jpbfe) * zws4 
     315            tra(ji,jj,ikt,jpsfe) = tra(ji,jj,ikt,jpsfe) - trb(ji,jj,ikt,jpsfe) * zws3 
     316            zwstpoc              = trb(ji,jj,ikt,jpgoc) * zws4 + trb(ji,jj,ikt,jppoc) * zws3 
    323317# else 
    324             trn(ji,jj,ikt,jpnum) = trn(ji,jj,ikt,jpnum) - trn(ji,jj,ikt,jpnum) * zws4 
    325             trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) - trn(ji,jj,ikt,jppoc) * zws3 
    326             trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) - trn(ji,jj,ikt,jpsfe) * zws3 
    327             zwstpoc = trn(ji,jj,ikt,jppoc) * zws3  
     318            tra(ji,jj,ikt,jpnum) = tra(ji,jj,ikt,jpnum) - trb(ji,jj,ikt,jpnum) * zws4  
     319            tra(ji,jj,ikt,jppoc) = tra(ji,jj,ikt,jppoc) - trb(ji,jj,ikt,jppoc) * zws3 
     320            tra(ji,jj,ikt,jpsfe) = tra(ji,jj,ikt,jpsfe) - trb(ji,jj,ikt,jpsfe) * zws3 
     321            zwstpoc = trb(ji,jj,ikt,jppoc) * zws3  
    328322# endif 
    329323 
     
    331325            ! The 0.5 factor in zpdenit and zdenitt is to avoid negative NO3 concentration after both denitrification 
    332326            ! in the sediments and just above the sediments. Not very clever, but simpliest option. 
    333             zpdenit  = MIN( 0.5 * ( trn(ji,jj,ikt,jpno3) - rtrn ) / rdenit, zdenit2d(ji,jj) * zwstpoc * zrivno3 ) 
     327            zpdenit  = MIN( 0.5 * ( trb(ji,jj,ikt,jpno3) - rtrn ) / rdenit, zdenit2d(ji,jj) * zwstpoc * zrivno3 ) 
    334328            z1pdenit = zwstpoc * zrivno3 - zpdenit 
    335             zolimit = MIN( ( trn(ji,jj,ikt,jpoxy) - rtrn ) / o2ut, z1pdenit * ( 1.- nitrfac(ji,jj,ikt) ) ) 
    336             zdenitt = MIN(  0.5 * ( trn(ji,jj,ikt,jpno3) - rtrn ) / rdenit, z1pdenit * nitrfac(ji,jj,ikt) ) 
    337             trn(ji,jj,ikt,jpdoc) = trn(ji,jj,ikt,jpdoc) + z1pdenit - zolimit - zdenitt 
    338             trn(ji,jj,ikt,jppo4) = trn(ji,jj,ikt,jppo4) + zpdenit + zolimit + zdenitt 
    339             trn(ji,jj,ikt,jpnh4) = trn(ji,jj,ikt,jpnh4) + zpdenit + zolimit + zdenitt 
    340             trn(ji,jj,ikt,jpno3) = trn(ji,jj,ikt,jpno3) - rdenit * (zpdenit + zdenitt) 
    341             trn(ji,jj,ikt,jpoxy) = trn(ji,jj,ikt,jpoxy) - zolimit * o2ut 
    342             trn(ji,jj,ikt,jptal) = trn(ji,jj,ikt,jptal) + rno3 * (zolimit + (1.+rdenit) * (zpdenit + zdenitt) ) 
    343             trn(ji,jj,ikt,jpdic) = trn(ji,jj,ikt,jpdic) + zpdenit + zolimit + zdenitt 
    344             zwork4(ji,jj) = rdenit * zpdenit * fse3t(ji,jj,ikt) 
     329            zolimit = MIN( ( trb(ji,jj,ikt,jpoxy) - rtrn ) / o2ut, z1pdenit * ( 1.- nitrfac(ji,jj,ikt) ) ) 
     330            zdenitt = MIN(  0.5 * ( trb(ji,jj,ikt,jpno3) - rtrn ) / rdenit, z1pdenit * nitrfac(ji,jj,ikt) ) 
     331            tra(ji,jj,ikt,jpdoc) = tra(ji,jj,ikt,jpdoc) + z1pdenit - zolimit - zdenitt 
     332            tra(ji,jj,ikt,jppo4) = tra(ji,jj,ikt,jppo4) + zpdenit + zolimit + zdenitt 
     333            tra(ji,jj,ikt,jpnh4) = tra(ji,jj,ikt,jpnh4) + zpdenit + zolimit + zdenitt 
     334            tra(ji,jj,ikt,jpno3) = tra(ji,jj,ikt,jpno3) - rdenit * (zpdenit + zdenitt) 
     335            tra(ji,jj,ikt,jpoxy) = tra(ji,jj,ikt,jpoxy) - zolimit * o2ut 
     336            tra(ji,jj,ikt,jptal) = tra(ji,jj,ikt,jptal) + rno3 * (zolimit + (1.+rdenit) * (zpdenit + zdenitt) ) 
     337            tra(ji,jj,ikt,jpdic) = tra(ji,jj,ikt,jpdic) + zpdenit + zolimit + zdenitt 
     338            sdenit(ji,jj) = rdenit * zpdenit * fse3t(ji,jj,ikt) 
    345339#endif 
    346340         END DO 
     
    362356#endif 
    363357               ztrfer = biron(ji,jj,jk)       / ( concfediaz + biron(ji,jj,jk)       ) 
    364                ztrpo4 = trn  (ji,jj,jk,jppo4) / ( concnnh4   + trn  (ji,jj,jk,jppo4) )  
     358               ztrpo4 = trb  (ji,jj,jk,jppo4) / ( concnnh4   + trb  (ji,jj,jk,jppo4) )  
    365359               zlight =  ( 1.- EXP( -etot_ndcy(ji,jj,jk) / diazolight ) )  
    366                znitrpot(ji,jj,jk) =  MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) * r1_rday )   & 
     360               nitrpot(ji,jj,jk) =  MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) * r1_rday )   & 
    367361                 &         *  zfact * MIN( ztrfer, ztrpo4 ) * zlight 
    368362               zsoufer(ji,jj,jk) = zlight * 2E-11 / (2E-11 + biron(ji,jj,jk)) 
     
    376370         DO jj = 1, jpj 
    377371            DO ji = 1, jpi 
    378                zfact = znitrpot(ji,jj,jk) * nitrfix 
    379                trn(ji,jj,jk,jpnh4) = trn(ji,jj,jk,jpnh4) +             zfact 
    380                trn(ji,jj,jk,jptal) = trn(ji,jj,jk,jptal) + rno3      * zfact 
    381                trn(ji,jj,jk,jpoxy) = trn(ji,jj,jk,jpoxy) + o2nit     * zfact  
    382                trn(ji,jj,jk,jppo4) = trn(ji,jj,jk,jppo4) + concdnh4 / ( concdnh4 + trn(ji,jj,jk,jppo4) ) & 
    383                &                     * 0.002 * trn(ji,jj,jk,jpdoc) * rfact2 / rday 
    384                trn(ji,jj,jk,jpfer) = trn(ji,jj,jk,jpfer) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * rfact2 / rday 
     372               zfact = nitrpot(ji,jj,jk) * nitrfix 
     373               tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) +             zfact 
     374               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3      * zfact 
     375               tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + o2nit     * zfact  
     376               tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + concdnh4 / ( concdnh4 + trb(ji,jj,jk,jppo4) ) & 
     377               &                     * 0.002 * trb(ji,jj,jk,jpdoc) * xstep 
     378               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * xstep 
    385379           END DO 
    386380         END DO  
    387381      END DO 
    388382 
    389       ! Global budget of N SMS : denitrification in the water column and in the sediment 
    390       !                          nitrogen fixation by the diazotrophs 
    391       ! -------------------------------------------------------------------------------- 
    392       zrdenittot   = glob_sum ( denitr(:,:,:) * rdenit * xnegtr(:,:,:) * cvol(:,:,:) ) 
    393       zsdenittot   = glob_sum ( zwork4(:,:)   * e1e2t(:,:) ) 
    394       znitrpottot  = glob_sum ( znitrpot(:,:,:) * nitrfix              * cvol(:,:,:) ) 
    395       zfact = 1.e+3 * rfact2r * rno3 * 14. / 1.e12 * ryyss   ! conversion molC/l/kt ----> TgN/m3/year 
    396       ! 
    397       IF( ln_check_mass .AND. ( kt == nitend .AND. jnt == nrdttrc ) .AND. ( lwp )  )  & 
    398          &  WRITE(numnit,9100) ndastp, znitrpottot * zfact  , & 
    399          &                             zrdenittot  * zfact  , & 
    400          &                             zsdenittot  * zfact 
    401       ! 
    402383      IF( lk_iomput ) THEN 
    403          IF( jnt == nrdttrc ) THEN 
     384         IF( knt == nrdttrc ) THEN 
    404385            zfact = 1.e+3 * rfact2r * rno3  !  conversion from molC/l/kt  to molN/m3/s 
    405             IF( iom_use("Nfix"   ) ) CALL iom_put( "Nfix"  , znitrpot(:,:,:) * nitrfix * zfact * tmask(:,:,:) )  ! nitrogen fixation  
    406             IF( iom_use("Sdenit" ) ) CALL iom_put( "Sdenit", zwork4(:,:) * zfact * tmask(:,:,1) )  ! Nitrate reduction in the sediments 
    407             IF( iom_use("tnfix"  ) ) CALL iom_put( "tnfix"  , znitrpottot * zfact  )               ! Global  nitrogen fixation 
    408             IF( iom_use("tdenit" ) ) CALL iom_put( "tdenit" , zrdenittot  * zfact  )               ! Total denitrification 
     386            IF( iom_use("Nfix"   ) ) CALL iom_put( "Nfix", nitrpot(:,:,:) * nitrfix * zfact * tmask(:,:,:) )  ! nitrogen fixation  
    409387            IF( iom_use("INTNFIX") ) THEN   ! nitrogen fixation rate in ocean ( vertically integrated ) 
    410388               zwork1(:,:) = 0. 
    411389               DO jk = 1, jpkm1 
    412                  zwork1(:,:) = zwork1(:,:) + znitrpot(:,:,jk) * nitrfix * zfact * fse3t(:,:,jk) * tmask(:,:,jk) 
     390                 zwork1(:,:) = zwork1(:,:) + nitrpot(:,:,jk) * nitrfix * zfact * fse3t(:,:,jk) * tmask(:,:,jk) 
    413391               ENDDO 
    414392               CALL iom_put( "INTNFIX" , zwork1 )  
     
    417395      ELSE 
    418396         IF( ln_diatrc )  & 
    419             &  trc2d(:,:,jp_pcs0_2d + 12) = znitrpot(:,:,1) * nitrfix * rno3 * 1.e+3 * rfact2r * fse3t(:,:,1) * tmask(:,:,1) 
     397            &  trc2d(:,:,jp_pcs0_2d + 12) = nitrpot(:,:,1) * nitrfix * rno3 * 1.e+3 * rfact2r * fse3t(:,:,1) * tmask(:,:,1) 
    420398      ENDIF 
    421399      ! 
     
    423401         WRITE(charout, fmt="('sed ')") 
    424402         CALL prt_ctl_trc_info(charout) 
    425          CALL prt_ctl_trc(tab4d=trn, mask=tmask, clinfo=ctrcnm) 
    426       ENDIF 
    427       ! 
    428       CALL wrk_dealloc( jpi, jpj, zdenit2d, zwork1, zwork2, zwork3, zwork4, zbureff ) 
     403         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm) 
     404      ENDIF 
     405      ! 
     406      CALL wrk_dealloc( jpi, jpj, zdenit2d, zwork1, zwork2, zwork3, zbureff ) 
    429407      CALL wrk_dealloc( jpi, jpj, zwsbio3, zwsbio4, zwscal ) 
    430       CALL wrk_dealloc( jpi, jpj, jpk, znitrpot, zsoufer ) 
     408      CALL wrk_dealloc( jpi, jpj, jpk, zsoufer ) 
    431409      ! 
    432410      IF( nn_timing == 1 )  CALL timing_stop('p4z_sed') 
     
    435413      ! 
    436414   END SUBROUTINE p4z_sed 
     415 
     416 
     417   INTEGER FUNCTION p4z_sed_alloc() 
     418      !!---------------------------------------------------------------------- 
     419      !!                     ***  ROUTINE p4z_sed_alloc  *** 
     420      !!---------------------------------------------------------------------- 
     421      ALLOCATE( nitrpot(jpi,jpj,jpk), sdenit(jpi,jpj), STAT=p4z_sed_alloc ) 
     422      ! 
     423      IF( p4z_sed_alloc /= 0 )   CALL ctl_warn('p4z_sed_alloc: failed to allocate arrays') 
     424      ! 
     425   END FUNCTION p4z_sed_alloc 
     426 
    437427 
    438428#else 
Note: See TracChangeset for help on using the changeset viewer.