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 7180 for branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z – NEMO

Ignore:
Timestamp:
2016-11-03T16:41:10+01:00 (8 years ago)
Author:
aumont
Message:

various bug fixes on iron chemistry

Location:
branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z
Files:
6 edited

Legend:

Unmodified
Added
Removed
  • branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z/p5zmeso.F90

    r6453 r7180  
    371371               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zgrazpof - zgrazfffp + zfracfe 
    372372               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) - zgrazfffg + zgrapof - zfracfe 
    373                zfracal = trb(ji,jj,jk,jpcal) / (trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + rtrn ) 
    374                zgrazcal = ( zgrazffeg + zgrazpoc ) * (1. - part2) * zfracal 
     373               zfracal = trb(ji,jj,jk,jpcal) / ( trb(ji,jj,jk,jpgoc) + rtrn ) 
     374               zgrazcal = zgrazffeg * (1. - part2) * zfracal 
    375375#endif 
    376376               !  calcite production 
  • branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z/p5zmicro.F90

    r6453 r7180  
    8484      REAL(wp) :: zgraztotc, zgraztotn, zgraztotp, zgraztotf, zbasresn, zbasresp, zbasresf 
    8585      REAL(wp) :: zgradoc, zgradon, zgradop, zgraref, zgradoct, zgradont, zgradopt, zgrareft 
    86       REAL(wp) :: zexcess, zgraren, zgrarep, zgrarem, zfracal, zgrazcal 
     86      REAL(wp) :: zexcess, zgraren, zgrarep, zgrarem 
    8787      REAL(wp) :: zgrapoc, zgrapon, zgrapop, zgrapof, zprcaca, zmortz 
    8888      REAL(wp) :: zrespz, ztortz, zgrasratf, zgrasratn, zgrasratp 
     
    172172               zgrazpon  = zgrazpoc * trb(ji,jj,jk,jppon) / ( trb(ji,jj,jk,jppoc) + rtrn ) 
    173173               zgrazpop  = zgrazpoc * trb(ji,jj,jk,jppop) / ( trb(ji,jj,jk,jppoc) + rtrn ) 
     174               zgrazpof  = zgrazpoc* trb(ji,jj,jk,jpsfe) / (trb(ji,jj,jk,jppoc) + rtrn) 
    174175               zgrazdc   = zgraze  * ztmp4  * zdenom 
    175176               zgrazdn   = zgrazdc * trb(ji,jj,jk,jpndi) / (trb(ji,jj,jk,jpdia) + rtrn) 
    176177               zgrazdp   = zgrazdc * trb(ji,jj,jk,jppdi) / (trb(ji,jj,jk,jpdia) + rtrn) 
    177178               zgrazdf   = zgrazdc * trb(ji,jj,jk,jpdfe) / (trb(ji,jj,jk,jpdia) + rtrn) 
    178                zgrazpof  = zgrazpoc* trb(ji,jj,jk,jpsfe) / (trb(ji,jj,jk,jppoc) + rtrn) 
    179179               ! 
    180180               zgraztotc = zgraznc + zgrazpoc + zgrazdc + zgrazz + zgrazpc 
     
    296296               prodcal(ji,jj,jk) = prodcal(ji,jj,jk) + zprcaca  ! prodcal=prodcal(nanophy)+prodcal(microzoo)+prodcal(mesozoo) 
    297297               ! 
    298                zfracal = trb(ji,jj,jk,jpcal) / (trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + rtrn ) 
    299                zgrazcal = zgrazpoc * (1. - part ) * zfracal  
    300298               zprcaca = part * zprcaca 
    301                tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrarem - zprcaca + zgrazcal 
    302                tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + 2. * ( zgrazcal - zprcaca )    & 
     299               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrarem - zprcaca 
     300               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) - 2. * zprcaca     & 
    303301               &                     + rno3 * zgraren 
    304                tra(ji,jj,jk,jpcal) = tra(ji,jj,jk,jpcal) + zprcaca - zgrazcal 
     302               tra(ji,jj,jk,jpcal) = tra(ji,jj,jk,jpcal) + zprcaca 
    305303#if defined key_kriest 
    306304               tra(ji,jj,jk,jpnum) = tra(ji,jj,jk,jpnum) + ztortz * xkr_dmicro & 
  • branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z/p5zmort.F90

    r6453 r7180  
    204204               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zmortp * zfactfe 
    205205               prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + zmortp 
    206                prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zmortp 
    207206#endif 
    208207            END DO 
     
    285284               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zmortp2 * zfactfe 
    286285#else 
    287                tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zrespp2 + 0.5 * ztortp2 
    288                tra(ji,jj,jk,jpgon) = tra(ji,jj,jk,jpgon) + (zrespp2 + 0.5 * ztortp2) * zfactn 
    289                tra(ji,jj,jk,jpgop) = tra(ji,jj,jk,jpgop) + (zrespp2 + 0.5 * ztortp2) * zfactp 
    290                tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + (zrespp2 + 0.5 * ztortp2) * zfactfe 
    291                tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + 0.5 * ztortp2 
    292                tra(ji,jj,jk,jppon) = tra(ji,jj,jk,jppon) + 0.5 * ztortp2 * zfactn 
    293                tra(ji,jj,jk,jppop) = tra(ji,jj,jk,jppop) + 0.5 * ztortp2 * zfactp 
    294                tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + 0.5 * ztortp2 * zfactfe 
    295                prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + 0.5 * ztortp2 
    296                prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zrespp2 + 0.5 * ztortp2 
     286               tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zrespp2  
     287               tra(ji,jj,jk,jpgon) = tra(ji,jj,jk,jpgon) + zrespp2 * zfactn 
     288               tra(ji,jj,jk,jpgop) = tra(ji,jj,jk,jpgop) + zrespp2 * zfactp 
     289               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zrespp2 * zfactfe 
     290               tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + ztortp2 
     291               tra(ji,jj,jk,jppon) = tra(ji,jj,jk,jppon) + ztortp2 * zfactn 
     292               tra(ji,jj,jk,jppop) = tra(ji,jj,jk,jppop) + ztortp2 * zfactp 
     293               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + ztortp2 * zfactfe 
     294               prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + ztortp2 
     295               prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zrespp2 
    297296#endif 
    298297            END DO 
  • branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z/p5zpoc.F90

    r6966 r7180  
    128128                  remint = 0. 
    129129                  ! 
     130                  zsizek1  = fse3t(ji,jj,jk-1) / 2. / (wsbio4(ji,jj,jk-1) + rtrn) * tgfunc(ji,jj,jk-1) 
     131                  zsizek = fse3t(ji,jj,jk) / 2. / (wsbio4(ji,jj,jk) + rtrn) * tgfunc(ji,jj,jk) 
     132                  ! 
    130133                  IF ( fsdept(ji,jj,jk-1) <= zdep ) THEN 
    131134                    !  
     
    150153                       ! ----------------------------------------------------- 
    151154                       ! 
    152                        zsizek  = zdep / (wsbio2 + rtrn) * tgfunc(ji,jj,jk-1) 
    153                        zsizek1 = fse3t(ji,jj,jk) / 2. / (wsbio4(ji,jj,jk) + rtrn) * tgfunc(ji,jj,jk) 
    154155                       ! the concentration of each lability class is calculated 
    155156                       ! as the sum of the different sources and sinks 
    156157                       ! Please note that production of new GOC experiences 
    157158                       ! degradation  
    158                        alphag(ji,jj,jk,jn) = alphan(jn) / (reminp(jn) * tgfunc(ji,jj,jk-1) )       & 
    159                        &   * (1. - exp( -reminp(jn) * zsizek ) ) * exp( -reminp(jn) * zsizek1 )     & 
    160                        &   * zpoc + prodgoc(ji,jj,jk) * ( 1. - exp( -reminp(jn) * zsizek1 ) )        & 
    161                        &   * rday / rfact2 * alphan(jn) / reminp(jn) / tgfunc(ji,jj,jk) 
     159                       alphag(ji,jj,jk,jn) = alphag(ji,jj,jk-1,jn) * exp( -reminp(jn) * zsizek ) * zpoc & 
     160                       &   + prodgoc(ji,jj,jk) * alphan(jn) / tgfunc(ji,jj,jk) / reminp(jn)             & 
     161                       &   * ( 1. - exp( -reminp(jn) * zsizek ) ) * rday / rfact2 
    162162                       alphat = alphat + alphag(ji,jj,jk,jn) 
    163163                       remint = remint + alphag(ji,jj,jk,jn) * reminp(jn) 
     
    175175                    ! 
    176176                    DO jn = 1, jcpoc 
    177                        zsizek  = fse3t(ji,jj,jk-1) / 2. / (wsbio4(ji,jj,jk-1) + rtrn) * tgfunc(ji,jj,jk-1) 
    178                        zsizek1 = fse3t(ji,jj,jk) / 2. / (wsbio4(ji,jj,jk) + rtrn) * tgfunc(ji,jj,jk) 
    179                        alphag(ji,jj,jk,jn) = alphag(ji,jj,jk-1,jn) * exp( -reminp(jn)                & 
    180                        &   * ( zsizek + zsizek1 ) ) * zpoc + ( prodgoc(ji,jj,jk-1)  & 
    181                        &   / tgfunc(ji,jj,jk-1) * ( 1. - exp( -reminp(jn) * zsizek ) )                 & 
    182                        &   * exp( -reminp(jn) * zsizek1 ) + prodgoc(ji,jj,jk) / tgfunc(ji,jj,jk)       & 
    183                        &   * ( 1. - exp( -reminp(jn) * zsizek1 ) ) ) * rday / rfact2 * alphan(jn)     & 
    184                        &   / reminp(jn) 
     177                       alphag(ji,jj,jk,jn) = alphag(ji,jj,jk-1,jn) * exp( -reminp(jn) * ( zsizek              & 
     178                       &   + zsizek1 ) ) * zpoc + ( prodgoc(ji,jj,jk-1) / tgfunc(ji,jj,jk-1) * ( 1.           & 
     179                       &   - exp( -reminp(jn) * zsizek1 ) ) * exp( -reminp(jn) * zsizek ) + prodgoc(ji,jj,jk) & 
     180                       &   / tgfunc(ji,jj,jk) * ( 1. - exp( -reminp(jn) * zsizek ) ) ) * rday / rfact2 / reminp(jn) 
    185181                       alphat = alphat + alphag(ji,jj,jk,jn) 
    186182                       remint = remint + alphag(ji,jj,jk,jn) * reminp(jn) 
     
    212208               zopoc2 = zremig  * trb(ji,jj,jk,jpgoc) 
    213209               orem(ji,jj,jk) = zopoc2 
     210               zorem3(ji,jj,jk) = zremig * solgoc * trb(ji,jj,jk,jpgoc) 
    214211               zopon2 = xremipn / xremipc * zremig * trb(ji,jj,jk,jpgon) 
    215212               zopop2 = xremipp / xremipc * zremig * trb(ji,jj,jk,jpgop) 
    216213               zofer2 = xremipn / xremipc * zremig * trb(ji,jj,jk,jpbfe) 
    217                zorem3(ji,jj,jk) = zremig * solgoc * trb(ji,jj,jk,jpgoc) 
    218214 
    219215               ! Update the appropriate tracers trends 
     
    258254        DO jj = 1, jpj 
    259255           DO ji = 1, jpi 
    260               IF (tmask(ji,jj,jk) == 1.) THEN 
    261                 zdep = hmld(ji,jj) 
    262                 IF( fsdept(ji,jj,jk) <= zdep ) THEN 
    263                   totprod(ji,jj) = totprod(ji,jj) + prodpoc(ji,jj,jk) * fse3t(ji,jj,jk) * rday/ rfact2 
    264                   ! The temperature effect is included here 
    265                   totthick(ji,jj) = totthick(ji,jj) + fse3t(ji,jj,jk)* tgfunc(ji,jj,jk) 
    266                   totcons(ji,jj) = totcons(ji,jj) - conspoc(ji,jj,jk) * fse3t(ji,jj,jk) * rday/ rfact2    & 
    267                      &                / ( trb(ji,jj,jk,jppoc) + rtrn ) 
    268                 ENDIF 
     256              zdep = hmld(ji,jj) 
     257              IF (tmask(ji,jj,jk) == 1. .AND. fsdept(ji,jj,jk) <= zdep ) THEN 
     258                 totprod(ji,jj) = totprod(ji,jj) + prodpoc(ji,jj,jk) * fse3t(ji,jj,jk) * rday/ rfact2 
     259                 ! The temperature effect is included here 
     260                 totthick(ji,jj) = totthick(ji,jj) + fse3t(ji,jj,jk)* tgfunc(ji,jj,jk) 
     261                 totcons(ji,jj) = totcons(ji,jj) - conspoc(ji,jj,jk) * fse3t(ji,jj,jk) * rday/ rfact2    & 
     262                 &                / ( trb(ji,jj,jk,jppoc) + rtrn ) 
    269263              ENDIF 
    270264           END DO 
     
    319313                  remint = 0. 
    320314                  ! 
     315                  ! the scale factors are corrected with temperature 
     316                  zsizek1  = fse3t(ji,jj,jk-1) / 2. / (wsbio3(ji,jj,jk-1) + rtrn) * tgfunc(ji,jj,jk-1) 
     317                  zsizek = fse3t(ji,jj,jk) / 2. / (wsbio3(ji,jj,jk) + rtrn) * tgfunc(ji,jj,jk) 
     318                  ! 
    321319                  ! Special treatment of the level just below the MXL 
    322320                  ! See the comments in the GOC section 
     
    332330                    ! 
    333331                    DO jn = 1, jcpoc 
    334                        ! the scale factor is corrected with temperature 
    335                        zsizek1 = fse3t(ji,jj,jk) / 2. / (wsbio3(ji,jj,jk) + rtrn) * tgfunc(ji,jj,jk) 
    336332                       ! computation of the lability spectrum applying the  
    337333                       ! different sources and sinks 
    338                        alphap(ji,jj,jk,jn) = alphap(ji,jj,jk-1,jn) + zorem3(ji,jj,jk)                 & 
    339                        &   * ( 1. - exp( -reminp(jn) * zsizek1 ) ) * rday / rfact2                   & 
    340                        &   * alphag(ji,jj,jk,jn) / reminp(jn) / tgfunc(ji,jj,jk) 
     334                      alphap(ji,jj,jk,jn) = alphap(ji,jj,jk-1,jn) * exp( -reminp(jn) * zsizek ) * zpoc  & 
     335                       &   + ( prodpoc(ji,jj,jk) * alphan(jn) + zorem3(ji,jj,jk) * alphag(ji,jj,jk,jn) ) & 
     336                       &   / tgfunc(ji,jj,jk) / reminp(jn) * rday / rfact2 * ( 1. - exp( -reminp(jn)     & 
     337                       &   * zsizek ) ) 
    341338                       alphap(ji,jj,jk,jn) = MAX( 0., alphap(ji,jj,jk,jn) ) 
    342339                       alphat = alphat + alphap(ji,jj,jk,jn) 
     
    355352                    ! 
    356353                    DO jn = 1, jcpoc 
    357                        zsizek  = fse3t(ji,jj,jk-1) / 2. / (wsbio3(ji,jj,jk-1) + rtrn) * tgfunc(ji,jj,jk-1) 
    358                        zsizek1 = fse3t(ji,jj,jk) / 2. / (wsbio3(ji,jj,jk) + rtrn) * tgfunc(ji,jj,jk) 
    359                        alphap(ji,jj,jk,jn) = alphap(ji,jj,jk-1,jn) * exp( -reminp(jn)             & 
    360                        &   * ( zsizek + zsizek1 ) ) * zpoc + ( prodpoc(ji,jj,jk-1)  & 
    361                        &   / tgfunc(ji,jj,jk-1) * ( 1. - exp( -reminp(jn) * zsizek ) )              & 
    362                        &   * exp( -reminp(jn) * zsizek1 ) + prodpoc(ji,jj,jk) / tgfunc(ji,jj,jk)    & 
    363                        &   * ( 1. - exp( -reminp(jn) * zsizek1 ) ) ) * rday / rfact2 * alphan(jn)  & 
    364                        &   / reminp(jn) 
    365                        alphap(ji,jj,jk,jn) = alphap(ji,jj,jk,jn) + zorem3(ji,jj,jk-1)              & 
    366                        &  * alphag(ji,jj,jk-1,jn) / tgfunc(ji,jj,jk-1) * rday / rfact2 * ( 1.       & 
    367                        &  - exp( -reminp(jn) * zsizek ) ) * exp( -reminp(jn) * zsizek1 )                  & 
    368                        &  + zorem3(ji,jj,jk) * ( 1. - exp( -reminp(jn) * zsizek1 ) ) * rday         & 
    369                        &  / rfact2 * alphag(ji,jj,jk,jn) / reminp(jn) / tgfunc(ji,jj,jk) 
     354                       alphap(ji,jj,jk,jn) = alphap(ji,jj,jk-1,jn) * exp( -reminp(jn)                       & 
     355                       &   * ( zsizek + zsizek1 ) ) * zpoc + ( prodpoc(ji,jj,jk-1) * alphan(jn)             & 
     356                       &   + zorem3(ji,jj,jk-1) * alphag(ji,jj,jk-1,jn) ) * rday / rfact2 / reminp(jn)      & 
     357                       &   / tgfunc(ji,jj,jk-1) * ( 1. - exp( -reminp(jn) * zsizek1 ) ) * exp( -reminp(jn)  & 
     358                       &   * zsizek ) + ( prodpoc(ji,jj,jk) * alphan(jn) + zorem3(ji,jj,jk)                 & 
     359                       &   * alphag(ji,jj,jk,jn) ) * rday / rfact2 / reminp(jn) / tgfunc(ji,jj,jk) * ( 1.   & 
     360                       &   - exp( -reminp(jn) * zsizek ) ) 
    370361                       alphap(ji,jj,jk,jn) = max(0., alphap(ji,jj,jk,jn) ) 
    371362                       alphat = alphat + alphap(ji,jj,jk,jn) 
  • branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z/p5zprod.F90

    r6966 r7180  
    8888      REAL(wp) ::   zlim, zsilfac2, zsiborn, zprod, zprontot, zproptot, zprodtot 
    8989      REAL(wp) ::   zprnutmax, zdocprod, zprochln, zprochld, zprochlp 
    90       REAL(wp) ::   zpislopen, zpislopep, zpisloped 
     90      REAL(wp) ::   zpislopen, zpislopep, zpisloped, thetannm_n, thetandm_n, thetanpm_n 
    9191      REAL(wp) ::   zrum, zcodel, zargu, zval, zfeup 
    9292      REAL(wp) ::   zrfact2 
     
    382382                  znanotot = enano(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    383383                  zprod = rday * (zpronewn(ji,jj,jk) + zproregn(ji,jj,jk)) * zprchln(ji,jj,jk) * xlimphy(ji,jj,jk) 
    384                   zprochln = thetannm * zprod / ( zpislopeadn(ji,jj,jk) * znanotot + rtrn ) 
     384                  thetannm_n   = MIN ( thetannm, ( thetannm / (1. - 1.14 / 43.4 *tsn(ji,jj,jk,jp_tem)))   & 
     385                  &               * (1. - 1.14 / 43.4 * 20.)) 
     386                  zprochln = thetannm_n * zprod / ( zpislopeadn(ji,jj,jk) * znanotot + rtrn ) 
    385387                  zprochln = MAX(zprochln, chlcmin * 12. * zprorcan (ji,jj,jk) ) 
    386388                     !  production terms for picophyto. ( chlorophyll ) 
    387389                  zpicotot = epico(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    388390                  zprod = rday * (zpronewp(ji,jj,jk) + zproregp(ji,jj,jk)) * zprchlp(ji,jj,jk) * xlimpic(ji,jj,jk) 
    389                   zprochlp = thetanpm * zprod / ( zpislopeadp(ji,jj,jk) * zpicotot + rtrn ) 
     391                  thetanpm_n   = MIN ( thetanpm, ( thetanpm / (1. - 1.14 / 43.4 *tsn(ji,jj,jk,jp_tem)))   & 
     392                  &               * (1. - 1.14 / 43.4 * 20.)) 
     393                  zprochlp = thetanpm_n * zprod / ( zpislopeadp(ji,jj,jk) * zpicotot + rtrn ) 
    390394                  zprochlp = MAX(zprochlp, chlcmin * 12. * zprorcap(ji,jj,jk) ) 
    391395                  !  production terms for diatomees ( chlorophyll ) 
    392396                  zdiattot = ediat(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    393397                  zprod = rday * (zpronewd(ji,jj,jk) + zproregd(ji,jj,jk)) * zprchld(ji,jj,jk) * xlimdia(ji,jj,jk) 
    394                   zprochld = thetandm * zprod / ( zpislopeadd(ji,jj,jk) * zdiattot + rtrn ) 
     398                  thetandm_n   = MIN ( thetandm, ( thetandm / (1. - 1.14 / 43.4 *tsn(ji,jj,jk,jp_tem)))   & 
     399                  &               * (1. - 1.14 / 43.4 * 20.)) 
     400                  zprochld = thetandm_n * zprod / ( zpislopeadd(ji,jj,jk) * zdiattot + rtrn ) 
    395401                  zprochld = MAX(zprochld, chlcmin * 12. * zprorcad(ji,jj,jk) ) 
    396402                  !   Update the arrays TRA which contain the Chla sources and sinks 
  • branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z/p5zrem.F90

    r6841 r7180  
    7575      REAL(wp) ::   zremik, zremikc, zremikn, zremikp, zsiremin  
    7676      REAL(wp) ::   zsatur, zsatur2, znusil, znusil2, zdep, zdepmin, zfactdep 
    77       REAL(wp) ::   zbactfer, zolimit 
     77      REAL(wp) ::   zbactfer, zolimit, zonitr, zstep, zrfact2 
    7878      REAL(wp) ::   zosil, ztem, zdenitnh4, zolimic, zolimin, zolimip, zdenitrn, zdenitrp 
    79       REAL(wp) ::   zonitr, zstep, zrfact2 
    8079      CHARACTER (len=25) :: charout 
    8180      REAL(wp), POINTER, DIMENSION(:,:  ) :: ztempbac  
    82       REAL(wp), POINTER, DIMENSION(:,:,:) :: zdepbac, zwork1, zdepprod 
     81      REAL(wp), POINTER, DIMENSION(:,:,:) :: zdepbac, zwork1, zdepprod, zfacsi, zfacsib 
    8382      !!--------------------------------------------------------------------- 
    8483      ! 
     
    8786      ! Allocate temporary workspace 
    8887      CALL wrk_alloc( jpi, jpj,      ztempbac                  ) 
    89       CALL wrk_alloc( jpi, jpj, jpk, zdepbac, zdepprod, zwork1 ) 
     88      CALL wrk_alloc( jpi, jpj, jpk, zdepbac, zdepprod, zwork1, zfacsi, zfacsib ) 
    9089 
    9190      ! Initialisation of temprary arrys 
    9291      zdepprod(:,:,:) = 1._wp 
    9392      ztempbac(:,:)   = 0._wp 
     93      zfacsib(:,:,:)  = xsilab / ( 1.0 - xsilab ) 
     94      zfacsi(:,:,:)   = xsilab 
    9495 
    9596      ! Computation of the mean phytoplankton concentration as 
     
    236237      ENDIF 
    237238 
     239      ! Initialization of the array which contains the labile fraction 
     240      ! of bSi. Set to a constant in the upper ocean 
     241      ! --------------------------------------------------------------- 
     242 
    238243      DO jk = 1, jpkm1 
    239244         DO jj = 1, jpj 
     
    243248               zstep = zstep * facvol(ji,jj,jk) 
    244249# endif 
     250               zdep     = MAX( hmld(ji,jj), heup_01(ji,jj) ) 
     251               zsatur   = MAX( rtrn, ( sio3eq(ji,jj,jk) - trb(ji,jj,jk,jpsil) ) / ( sio3eq(ji,jj,jk) + rtrn ) ) 
     252               zsatur2  = ( 1. + tsn(ji,jj,jk,jp_tem) / 400.)**37 
     253               znusil   = 0.225  * ( 1. + tsn(ji,jj,jk,jp_tem) / 15.) * zsatur + 0.775 * zsatur2 * zsatur**9.25 
     254 
    245255               ! Remineralization rate of BSi depedant on T and saturation 
    246256               ! --------------------------------------------------------- 
    247                zsatur   = ( sio3eq(ji,jj,jk) - trb(ji,jj,jk,jpsil) ) / ( sio3eq(ji,jj,jk) + rtrn ) 
    248                zsatur   = MAX( rtrn, zsatur ) 
    249                zsatur2  = ( 1. + tsn(ji,jj,jk,jp_tem) / 400.)**37 
    250                znusil   = 0.225  * ( 1. + tsn(ji,jj,jk,jp_tem) / 15.) * zsatur + 0.775 * zsatur2 * zsatur**9.25 
    251                znusil2  = 0.225  * ( 1. + tsn(ji,jj,1,jp_tem) / 15.) + 0.775 * zsatur2 
    252  
    253                ! Two classes of BSi are considered : a labile fraction and  
    254                ! a more refractory one. The ratio between both fractions is 
    255                ! constant and specified in the namelist. 
    256                ! ---------------------------------------------------------- 
    257                zdep     = MAX( hmld(ji,jj), heup_01(ji,jj) )  
    258                zdep     = MAX( 0., fsdept(ji,jj,jk) - zdep ) 
    259                ztem     = MAX( tsn(ji,jj,1,jp_tem), 0. ) 
    260                zfactdep = xsilab * EXP(-( xsiremlab - xsirem ) * znusil2 * zdep / wsbio2 ) * ztem / ( ztem + 10. ) 
    261                zsiremin = ( xsiremlab * zfactdep + xsirem * ( 1. - zfactdep ) ) * zstep * znusil 
     257               IF ( fsdept(ji,jj,jk) > zdep ) THEN 
     258                  zfacsib(ji,jj,jk) = zfacsib(ji,jj,jk-1) * EXP( -0.5 * ( xsiremlab - xsirem )  & 
     259                  &                   * znusil * fse3t(ji,jj,jk) / wsbio4(ji,jj,jk) ) 
     260                  zfacsi(ji,jj,jk)  = zfacsib(ji,jj,jk) / ( 1.0 + zfacsib(ji,jj,jk) ) 
     261                  zfacsib(ji,jj,jk) = zfacsib(ji,jj,jk) * EXP( -0.5 * ( xsiremlab - xsirem )    & 
     262                  &                   * znusil * fse3t(ji,jj,jk) / wsbio4(ji,jj,jk) ) 
     263               ENDIF 
     264               zsiremin = ( xsiremlab * zfacsi(ji,jj,jk) + xsirem * ( 1. - zfacsi(ji,jj,jk) ) ) * zstep * znusil 
    262265               zosil    = zsiremin * trb(ji,jj,jk,jpgsi) 
    263266               ! 
    264267               tra(ji,jj,jk,jpgsi) = tra(ji,jj,jk,jpgsi) - zosil 
    265268               tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) + zosil 
    266                ! 
     269               ENDIF 
    267270            END DO 
    268271         END DO 
     
    280283          CALL iom_put( "DENIT" , denitrc(:,:,:) * rdenit * rno3 * tmask(:,:,:) * zrfact2  )  ! Denitrification 
    281284      ENDIF 
    282  
    283       IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
    284          WRITE(charout, FMT="('rem4')") 
    285          CALL prt_ctl_trc_info(charout) 
    286          CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm) 
    287       ENDIF 
    288285      ! 
    289286      CALL wrk_dealloc( jpi, jpj,      ztempbac                  ) 
    290       CALL wrk_dealloc( jpi, jpj, jpk, zdepbac, zdepprod, zwork1 ) 
     287      CALL wrk_dealloc( jpi, jpj, jpk, zdepbac, zdepprod, zwork1, zfacsi, zfacsib ) 
    291288      ! 
    292289      IF( nn_timing == 1 )  CALL timing_stop('p5z_rem') 
Note: See TracChangeset for help on using the changeset viewer.