Changeset 7180 for branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z/p5zrem.F90

Timestamp:

2016-11-03T16:41:10+01:00 (7 years ago)

Author:

aumont

Message:

various bug fixes on iron chemistry

File:

• branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z/p5zrem.F90 (modified) (5 diffs)

branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P5Z/p5zrem.F90

@@ -75 +75 @@
       REAL(wp) ::   zremik, zremikc, zremikn, zremikp, zsiremin 
       REAL(wp) ::   zsatur, zsatur2, znusil, znusil2, zdep, zdepmin, zfactdep
-      REAL(wp) ::   zbactfer, zolimit
+      REAL(wp) ::   zbactfer, zolimit, zonitr, zstep, zrfact2
       REAL(wp) ::   zosil, ztem, zdenitnh4, zolimic, zolimin, zolimip, zdenitrn, zdenitrp
-      REAL(wp) ::   zonitr, zstep, zrfact2
       CHARACTER (len=25) :: charout
       REAL(wp), POINTER, DIMENSION(:,:  ) :: ztempbac 
-      REAL(wp), POINTER, DIMENSION(:,:,:) :: zdepbac, zwork1, zdepprod
+      REAL(wp), POINTER, DIMENSION(:,:,:) :: zdepbac, zwork1, zdepprod, zfacsi, zfacsib
       !!---------------------------------------------------------------------
       !
@@ -87 +86 @@
       ! Allocate temporary workspace
       CALL wrk_alloc( jpi, jpj,      ztempbac                  )
-      CALL wrk_alloc( jpi, jpj, jpk, zdepbac, zdepprod, zwork1 )
+      CALL wrk_alloc( jpi, jpj, jpk, zdepbac, zdepprod, zwork1, zfacsi, zfacsib )
 
       ! Initialisation of temprary arrys
       zdepprod(:,:,:) = 1._wp
       ztempbac(:,:)   = 0._wp
+      zfacsib(:,:,:)  = xsilab / ( 1.0 - xsilab )
+      zfacsi(:,:,:)   = xsilab
 
       ! Computation of the mean phytoplankton concentration as
@@ -236 +237 @@
       ENDIF
 
+      ! Initialization of the array which contains the labile fraction
+      ! of bSi. Set to a constant in the upper ocean
+      ! ---------------------------------------------------------------
+
       DO jk = 1, jpkm1
          DO jj = 1, jpj
@@ -243 +248 @@
                zstep = zstep * facvol(ji,jj,jk)
 # endif
+               zdep     = MAX( hmld(ji,jj), heup_01(ji,jj) )
+               zsatur   = MAX( rtrn, ( sio3eq(ji,jj,jk) - trb(ji,jj,jk,jpsil) ) / ( sio3eq(ji,jj,jk) + rtrn ) )
+               zsatur2  = ( 1. + tsn(ji,jj,jk,jp_tem) / 400.)**37
+               znusil   = 0.225  * ( 1. + tsn(ji,jj,jk,jp_tem) / 15.) * zsatur + 0.775 * zsatur2 * zsatur**9.25
+
                ! Remineralization rate of BSi depedant on T and saturation
                ! ---------------------------------------------------------
-               zsatur   = ( sio3eq(ji,jj,jk) - trb(ji,jj,jk,jpsil) ) / ( sio3eq(ji,jj,jk) + rtrn )
-               zsatur   = MAX( rtrn, zsatur )
-               zsatur2  = ( 1. + tsn(ji,jj,jk,jp_tem) / 400.)**37
-               znusil   = 0.225  * ( 1. + tsn(ji,jj,jk,jp_tem) / 15.) * zsatur + 0.775 * zsatur2 * zsatur**9.25
-               znusil2  = 0.225  * ( 1. + tsn(ji,jj,1,jp_tem) / 15.) + 0.775 * zsatur2
-
-               ! Two classes of BSi are considered : a labile fraction and 
-               ! a more refractory one. The ratio between both fractions is
-               ! constant and specified in the namelist.
-               ! ----------------------------------------------------------
-               zdep     = MAX( hmld(ji,jj), heup_01(ji,jj) ) 
-               zdep     = MAX( 0., fsdept(ji,jj,jk) - zdep )
-               ztem     = MAX( tsn(ji,jj,1,jp_tem), 0. )
-               zfactdep = xsilab * EXP(-( xsiremlab - xsirem ) * znusil2 * zdep / wsbio2 ) * ztem / ( ztem + 10. )
-               zsiremin = ( xsiremlab * zfactdep + xsirem * ( 1. - zfactdep ) ) * zstep * znusil
+               IF ( fsdept(ji,jj,jk) > zdep ) THEN
+                  zfacsib(ji,jj,jk) = zfacsib(ji,jj,jk-1) * EXP( -0.5 * ( xsiremlab - xsirem )  &
+                  &                   * znusil * fse3t(ji,jj,jk) / wsbio4(ji,jj,jk) )
+                  zfacsi(ji,jj,jk)  = zfacsib(ji,jj,jk) / ( 1.0 + zfacsib(ji,jj,jk) )
+                  zfacsib(ji,jj,jk) = zfacsib(ji,jj,jk) * EXP( -0.5 * ( xsiremlab - xsirem )    &
+                  &                   * znusil * fse3t(ji,jj,jk) / wsbio4(ji,jj,jk) )
+               ENDIF
+               zsiremin = ( xsiremlab * zfacsi(ji,jj,jk) + xsirem * ( 1. - zfacsi(ji,jj,jk) ) ) * zstep * znusil
                zosil    = zsiremin * trb(ji,jj,jk,jpgsi)
                !
                tra(ji,jj,jk,jpgsi) = tra(ji,jj,jk,jpgsi) - zosil
                tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) + zosil
-               !
+               ENDIF
             END DO
          END DO
@@ -280 +283 @@
           CALL iom_put( "DENIT" , denitrc(:,:,:) * rdenit * rno3 * tmask(:,:,:) * zrfact2  )  ! Denitrification
       ENDIF
-
-      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
-         WRITE(charout, FMT="('rem4')")
-         CALL prt_ctl_trc_info(charout)
-         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
-      ENDIF
       !
       CALL wrk_dealloc( jpi, jpj,      ztempbac                  )
-      CALL wrk_dealloc( jpi, jpj, jpk, zdepbac, zdepprod, zwork1 )
+      CALL wrk_dealloc( jpi, jpj, jpk, zdepbac, zdepprod, zwork1, zfacsi, zfacsib )
       !
       IF( nn_timing == 1 )  CALL timing_stop('p5z_rem')