Index: trunk/NEMOGCM/NEMO/LIM_SRC_3/ice.F90
===================================================================
--- trunk/NEMOGCM/NEMO/LIM_SRC_3/ice.F90	(revision 6387)
+++ trunk/NEMOGCM/NEMO/LIM_SRC_3/ice.F90	(revision 6403)
@@ -301,4 +301,5 @@
 
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   ftr_ice   !: transmitted solar radiation under ice
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   rn_amax_2d  !: maximum ice concentration 2d array
 
    !!--------------------------------------------------------------------------
@@ -378,5 +379,6 @@
    LOGICAL          , PUBLIC ::   ln_limdyn       !: flag for ice dynamics (T) or not (F)
    LOGICAL          , PUBLIC ::   ln_icectl       !: flag for sea-ice points output (T) or not (F)
-   REAL(wp)         , PUBLIC ::   rn_amax         !: maximum ice concentration
+   REAL(wp)         , PUBLIC ::   rn_amax_n       !: maximum ice concentration Northern hemisphere
+   REAL(wp)         , PUBLIC ::   rn_amax_s       !: maximum ice concentration Southern hemisphere
    INTEGER          , PUBLIC ::   iiceprt         !: debug i-point
    INTEGER          , PUBLIC ::   jiceprt         !: debug j-point
@@ -438,4 +440,5 @@
          &      afx_tot(jpi,jpj) , afx_thd(jpi,jpj),  afx_dyn(jpi,jpj) ,                        &
          &      fhtur  (jpi,jpj) , ftr_ice(jpi,jpj,jpl), qlead  (jpi,jpj) ,                     &
+         &      rn_amax_2d(jpi,jpj),                                                            &
          &      sfx_res(jpi,jpj) , sfx_bri(jpi,jpj) , sfx_dyn(jpi,jpj) ,                        &
          &      sfx_bog(jpi,jpj) , sfx_bom(jpi,jpj) , sfx_sum(jpi,jpj) , sfx_sni(jpi,jpj) , sfx_opw(jpi,jpj) ,    &
Index: trunk/NEMOGCM/NEMO/LIM_SRC_3/limcons.F90
===================================================================
--- trunk/NEMOGCM/NEMO/LIM_SRC_3/limcons.F90	(revision 6387)
+++ trunk/NEMOGCM/NEMO/LIM_SRC_3/limcons.F90	(revision 6403)
@@ -256,5 +256,6 @@
             ENDIF
             IF (     zvmin   < -epsi10 ) WRITE(numout,*) 'violation v_i<0  [m]          (',cd_routine,') = ',zvmin
-            IF (     zamax   > rn_amax+epsi10 .AND. cd_routine /= 'limtrp' .AND. cd_routine /= 'limitd_me' ) THEN
+            IF (     zamax   > MAX( rn_amax_n, rn_amax_s ) + epsi10 .AND. &
+               &                         cd_routine /= 'limtrp' .AND. cd_routine /= 'limitd_me' ) THEN
                                          WRITE(numout,*) 'violation a_i>amax            (',cd_routine,') = ',zamax
             ENDIF
Index: trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90
===================================================================
--- trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90	(revision 6387)
+++ trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90	(revision 6403)
@@ -297,13 +297,14 @@
          END DO
 
-         CALL tab_2d_1d( nbpac, qlead_1d  (1:nbpac)     , qlead  , jpi, jpj, npac(1:nbpac) )
-         CALL tab_2d_1d( nbpac, t_bo_1d   (1:nbpac)     , t_bo   , jpi, jpj, npac(1:nbpac) )
-         CALL tab_2d_1d( nbpac, sfx_opw_1d(1:nbpac)     , sfx_opw, jpi, jpj, npac(1:nbpac) )
-         CALL tab_2d_1d( nbpac, wfx_opw_1d(1:nbpac)     , wfx_opw, jpi, jpj, npac(1:nbpac) )
-         CALL tab_2d_1d( nbpac, hicol_1d  (1:nbpac)     , hicol  , jpi, jpj, npac(1:nbpac) )
-         CALL tab_2d_1d( nbpac, zvrel_1d  (1:nbpac)     , zvrel  , jpi, jpj, npac(1:nbpac) )
-
-         CALL tab_2d_1d( nbpac, hfx_thd_1d(1:nbpac)     , hfx_thd, jpi, jpj, npac(1:nbpac) )
-         CALL tab_2d_1d( nbpac, hfx_opw_1d(1:nbpac)     , hfx_opw, jpi, jpj, npac(1:nbpac) )
+         CALL tab_2d_1d( nbpac, qlead_1d  (1:nbpac)     , qlead     , jpi, jpj, npac(1:nbpac) )
+         CALL tab_2d_1d( nbpac, t_bo_1d   (1:nbpac)     , t_bo      , jpi, jpj, npac(1:nbpac) )
+         CALL tab_2d_1d( nbpac, sfx_opw_1d(1:nbpac)     , sfx_opw   , jpi, jpj, npac(1:nbpac) )
+         CALL tab_2d_1d( nbpac, wfx_opw_1d(1:nbpac)     , wfx_opw   , jpi, jpj, npac(1:nbpac) )
+         CALL tab_2d_1d( nbpac, hicol_1d  (1:nbpac)     , hicol     , jpi, jpj, npac(1:nbpac) )
+         CALL tab_2d_1d( nbpac, zvrel_1d  (1:nbpac)     , zvrel     , jpi, jpj, npac(1:nbpac) )
+
+         CALL tab_2d_1d( nbpac, hfx_thd_1d(1:nbpac)     , hfx_thd   , jpi, jpj, npac(1:nbpac) )
+         CALL tab_2d_1d( nbpac, hfx_opw_1d(1:nbpac)     , hfx_opw   , jpi, jpj, npac(1:nbpac) )
+         CALL tab_2d_1d( nbpac, rn_amax_1d(1:nbpac)     , rn_amax_2d, jpi, jpj, npac(1:nbpac) )
 
          !------------------------------------------------------------------------------!
@@ -409,6 +410,6 @@
          ! we keep the excessive volume in memory and attribute it later to bottom accretion
          DO ji = 1, nbpac
-            IF ( za_newice(ji) >  ( rn_amax - zat_i_1d(ji) ) ) THEN
-               zda_res(ji)   = za_newice(ji) - ( rn_amax - zat_i_1d(ji) )
+            IF ( za_newice(ji) >  ( rn_amax_1d(ji) - zat_i_1d(ji) ) ) THEN
+               zda_res(ji)   = za_newice(ji) - ( rn_amax_1d(ji) - zat_i_1d(ji) )
                zdv_res(ji)   = zda_res  (ji) * zh_newice(ji) 
                za_newice(ji) = za_newice(ji) - zda_res  (ji)
Index: trunk/NEMOGCM/NEMO/LIM_SRC_3/limtrp.F90
===================================================================
--- trunk/NEMOGCM/NEMO/LIM_SRC_3/limtrp.F90	(revision 6387)
+++ trunk/NEMOGCM/NEMO/LIM_SRC_3/limtrp.F90	(revision 6403)
@@ -422,5 +422,5 @@
             DO jj = 1, jpj
                DO ji = 1, jpi
-                  a_i(ji,jj,1)  = MIN( a_i(ji,jj,1), rn_amax )
+                  a_i(ji,jj,1)  = MIN( a_i(ji,jj,1), rn_amax_2d(ji,jj) )
                END DO
             END DO
Index: trunk/NEMOGCM/NEMO/LIM_SRC_3/limupdate1.F90
===================================================================
--- trunk/NEMOGCM/NEMO/LIM_SRC_3/limupdate1.F90	(revision 6387)
+++ trunk/NEMOGCM/NEMO/LIM_SRC_3/limupdate1.F90	(revision 6403)
@@ -80,7 +80,7 @@
          DO jj = 1, jpj
             DO ji = 1, jpi
-               IF( at_i(ji,jj) > rn_amax .AND. a_i(ji,jj,jl) > 0._wp ) THEN
-                  a_i (ji,jj,jl) = a_i (ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax / at_i(ji,jj) ) )
-                  oa_i(ji,jj,jl) = oa_i(ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax / at_i(ji,jj) ) )
+               IF( at_i(ji,jj) > rn_amax_2d(ji,jj) .AND. a_i(ji,jj,jl) > 0._wp ) THEN
+                  a_i (ji,jj,jl) = a_i (ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax_2d(ji,jj) / at_i(ji,jj) ) )
+                  oa_i(ji,jj,jl) = oa_i(ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax_2d(ji,jj) / at_i(ji,jj) ) )
                ENDIF
             END DO
Index: trunk/NEMOGCM/NEMO/LIM_SRC_3/limupdate2.F90
===================================================================
--- trunk/NEMOGCM/NEMO/LIM_SRC_3/limupdate2.F90	(revision 6387)
+++ trunk/NEMOGCM/NEMO/LIM_SRC_3/limupdate2.F90	(revision 6403)
@@ -94,7 +94,7 @@
          DO jj = 1, jpj
             DO ji = 1, jpi
-               IF( at_i(ji,jj) > rn_amax .AND. a_i(ji,jj,jl) > 0._wp ) THEN
-                  a_i (ji,jj,jl) = a_i (ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax / at_i(ji,jj) ) )
-                  oa_i(ji,jj,jl) = oa_i(ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax / at_i(ji,jj) ) )
+               IF( at_i(ji,jj) > rn_amax_2d(ji,jj) .AND. a_i(ji,jj,jl) > 0._wp ) THEN
+                  a_i (ji,jj,jl) = a_i (ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax_2d(ji,jj) / at_i(ji,jj) ) )
+                  oa_i(ji,jj,jl) = oa_i(ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax_2d(ji,jj) / at_i(ji,jj) ) )
                ENDIF
             END DO
Index: trunk/NEMOGCM/NEMO/LIM_SRC_3/thd_ice.F90
===================================================================
--- trunk/NEMOGCM/NEMO/LIM_SRC_3/thd_ice.F90	(revision 6387)
+++ trunk/NEMOGCM/NEMO/LIM_SRC_3/thd_ice.F90	(revision 6403)
@@ -51,4 +51,5 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   hfx_dif_1d
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   hfx_opw_1d
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   rn_amax_1d
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   hfx_snw_1d
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   hfx_err_1d
@@ -144,4 +145,5 @@
          &      hfx_sum_1d(jpij) , hfx_bom_1d(jpij) ,hfx_bog_1d(jpij) ,    & 
          &      hfx_dif_1d(jpij) , hfx_opw_1d(jpij) ,                      &
+         &      rn_amax_1d(jpij) ,                                         &
          &      hfx_thd_1d(jpij) , hfx_spr_1d(jpij) ,                      &
          &      hfx_snw_1d(jpij) , hfx_sub_1d(jpij) , hfx_err_1d(jpij) ,   &
Index: trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
===================================================================
--- trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90	(revision 6387)
+++ trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90	(revision 6403)
@@ -104,5 +104,5 @@
       INTEGER, INTENT(in) ::   kblk    ! type of bulk (=3 CLIO, =4 CORE, =5 COUPLED)
       !!
-      INTEGER  ::   jl                 ! dummy loop index
+      INTEGER  ::    jl                 ! dummy loop index
       REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_os, zalb_cs  ! ice albedo under overcast/clear sky
       REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_ice          ! mean ice albedo (for coupled)
@@ -260,5 +260,5 @@
       !! ** purpose :   Allocate all the dynamic arrays of the LIM-3 modules
       !!----------------------------------------------------------------------
-      INTEGER :: ierr
+      INTEGER :: ji, jj, ierr
       !!----------------------------------------------------------------------
       IF(lwp) WRITE(numout,*)
@@ -317,4 +317,12 @@
       tn_ice(:,:,:) = t_su(:,:,:)       ! initialisation of surface temp for coupled simu
       !
+      DO jj = 1, jpj
+         DO ji = 1, jpi
+            IF( gphit(ji,jj) > 0._wp ) THEN  ;  rn_amax_2d(ji,jj) = rn_amax_n  ! NH
+            ELSE                             ;  rn_amax_2d(ji,jj) = rn_amax_s  ! SH
+            ENDIF
+        ENDDO
+      ENDDO 
+      !
       nstart = numit  + nn_fsbc      
       nitrun = nitend - nit000 + 1 
@@ -339,5 +347,5 @@
       INTEGER  ::   ios                 ! Local integer output status for namelist read
       NAMELIST/namicerun/ jpl, nlay_i, nlay_s, cn_icerst_in, cn_icerst_indir, cn_icerst_out, cn_icerst_outdir,  &
-         &                ln_limdyn, rn_amax, ln_limdiahsb, ln_limdiaout, ln_icectl, iiceprt, jiceprt  
+         &                ln_limdyn, rn_amax_n, rn_amax_s, ln_limdiahsb, ln_limdiaout, ln_icectl, iiceprt, jiceprt  
       !!-------------------------------------------------------------------
       !                    
@@ -359,5 +367,6 @@
          WRITE(numout,*) '   number of snow layers                                   = ', nlay_s
          WRITE(numout,*) '   switch for ice dynamics (1) or not (0)      ln_limdyn   = ', ln_limdyn
-         WRITE(numout,*) '   maximum ice concentration                               = ', rn_amax 
+         WRITE(numout,*) '   maximum ice concentration for NH                        = ', rn_amax_n 
+         WRITE(numout,*) '   maximum ice concentration for SH                        = ', rn_amax_s
          WRITE(numout,*) '   Diagnose heat/salt budget or not          ln_limdiahsb  = ', ln_limdiahsb
          WRITE(numout,*) '   Output   heat/salt budget or not          ln_limdiaout  = ', ln_limdiaout
Index: trunk/NEMOGCM/NEMO/OPA_SRC/TRA/traqsr.F90
===================================================================
--- trunk/NEMOGCM/NEMO/OPA_SRC/TRA/traqsr.F90	(revision 6387)
+++ trunk/NEMOGCM/NEMO/OPA_SRC/TRA/traqsr.F90	(revision 6403)
@@ -11,4 +11,5 @@
    !!            3.2  !  2009-04  (G. Madec & NEMO team) 
    !!            3.6  !  2012-05  (C. Rousset) store attenuation coef for use in ice model 
+   !!            3.6  !  2015-12  (O. Aumont, J. Jouanno, C. Ethe) use vertical profile of chlorophyll
    !!            3.7  !  2015-11  (G. Madec, A. Coward)  remove optimisation for fix volume 
    !!----------------------------------------------------------------------
@@ -100,4 +101,5 @@
       !! Reference  : Jerlov, N. G., 1968 Optical Oceanography, Elsevier, 194pp.
       !!              Lengaigne et al. 2007, Clim. Dyn., V28, 5, 503-516.
+      !!              Morel, A. et Berthon, JF, 1989, Limnol Oceanogr 34(8), 1545-1562
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt     ! ocean time-step
@@ -109,7 +111,9 @@
       REAL(wp) ::   zzc0, zzc1, zzc2, zzc3   !    -         -
       REAL(wp) ::   zz0 , zz1                !    -         -
+      REAL(wp) ::   zCb, zCmax, zze, zpsi, zpsimax, zdelpsi, zCtot, zCze
+      REAL(wp) ::   zlogc, zlogc2, zlogc3 
       REAL(wp), POINTER, DIMENSION(:,:)   :: zekb, zekg, zekr
       REAL(wp), POINTER, DIMENSION(:,:,:) :: ze0, ze1, ze2, ze3, zea, ztrdt
-      REAL(wp), POINTER, DIMENSION(:,:,:) :: zetot
+      REAL(wp), POINTER, DIMENSION(:,:,:) :: zetot, zchl3d
       !!----------------------------------------------------------------------
       !
@@ -158,28 +162,35 @@
          !
          CALL wrk_alloc( jpi,jpj,       zekb, zekg, zekr        ) 
-         CALL wrk_alloc( jpi,jpj,jpk,   ze0, ze1, ze2, ze3, zea ) 
+         CALL wrk_alloc( jpi,jpj,jpk,   ze0, ze1, ze2, ze3, zea, zchl3d ) 
          !
          IF( nqsr == np_RGBc ) THEN          !*  Variable Chlorophyll
             CALL fld_read( kt, 1, sf_chl )         ! Read Chl data and provides it at the current time step
-            DO jj = 2, jpjm1                       ! Separation in R-G-B depending of the surface Chl
-               DO ji = fs_2, fs_jpim1
-                  zchl = MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) )
-                  irgb = NINT( 41 + 20.*LOG10(zchl) + 1.e-15 )
-                  zekb(ji,jj) = rkrgb(1,irgb)
-                  zekg(ji,jj) = rkrgb(2,irgb)
-                  zekr(ji,jj) = rkrgb(3,irgb)
+            DO jk = 1, nksr + 1
+               DO jj = 2, jpjm1                       ! Separation in R-G-B depending of the surface Chl
+                  DO ji = fs_2, fs_jpim1
+                     zchl    = sf_chl(1)%fnow(ji,jj,1)
+                     zCtot   = 40.6  * zchl**0.459
+                     zze     = 568.2 * zCtot**(-0.746)
+                     IF( zze > 102. ) zze = 200.0 * zCtot**(-0.293)
+                     zpsi    = gdepw_n(ji,jj,jk) / zze
+                     !
+                     zlogc   = LOG( zchl )
+                     zlogc2  = zlogc * zlogc
+                     zlogc3  = zlogc * zlogc * zlogc
+                     zCb     = 0.768 + 0.087 * zlogc - 0.179 * zlogc2 - 0.025 * zlogc3
+                     zCmax   = 0.299 - 0.289 * zlogc + 0.579 * zlogc2
+                     zpsimax = 0.6   - 0.640 * zlogc + 0.021 * zlogc2 + 0.115 * zlogc3
+                     zdelpsi = 0.710 + 0.159 * zlogc + 0.021 * zlogc2
+                     zCze    = 1.12  * (zchl)**0.803 
+                     !
+                     zchl3d(ji,jj,jk) = zCze * ( zCb + zCmax * EXP( -( (zpsi - zpsimax) / zdelpsi )**2 ) )
+                  END DO
+                  !
                END DO
             END DO
          ELSE                                !* constant chrlorophyll
-            zchl = 0.05                            ! constant chlorophyll
-            !                                      ! Separation in R-G-B depending of the chlorophyll
-            irgb = NINT( 41 + 20.*LOG10( zchl ) + 1.e-15 )
-            DO jj = 2, jpjm1
-               DO ji = fs_2, fs_jpim1
-                  zekb(ji,jj) = rkrgb(1,irgb)                      
-                  zekg(ji,jj) = rkrgb(2,irgb)
-                  zekr(ji,jj) = rkrgb(3,irgb)
-               END DO
-            END DO
+           DO jk = 1, nksr + 1
+              zchl3d(:,:,jk) = 0.05 
+            ENDDO
          ENDIF
          !
@@ -195,5 +206,15 @@
          END DO
          !
-         DO jk = 2, nksr+1                   !* interior equi-partition in R-G-B
+         DO jk = 2, nksr+1                   !* interior equi-partition in R-G-B depending of vertical profile of Chl
+            DO jj = 2, jpjm1
+               DO ji = fs_2, fs_jpim1
+                  zchl = MIN( 10. , MAX( 0.03, zchl3d(ji,jj,jk) ) )
+                  irgb = NINT( 41 + 20.*LOG10(zchl) + 1.e-15 )
+                  zekb(ji,jj) = rkrgb(1,irgb)
+                  zekg(ji,jj) = rkrgb(2,irgb)
+                  zekr(ji,jj) = rkrgb(3,irgb)
+               END DO
+            END DO
+
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1
@@ -220,5 +241,5 @@
          !
          CALL wrk_dealloc( jpi,jpj,        zekb, zekg, zekr        ) 
-         CALL wrk_dealloc( jpi,jpj,jpk,   ze0, ze1, ze2, ze3, zea ) 
+         CALL wrk_dealloc( jpi,jpj,jpk,   ze0, ze1, ze2, ze3, zea, zchl3d ) 
          !
       CASE( np_2BD  )            !==  2-bands fluxes  ==!
Index: trunk/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90
===================================================================
--- trunk/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90	(revision 6387)
+++ trunk/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90	(revision 6403)
@@ -40,4 +40,5 @@
    REAL(wp), PUBLIC ::   rn_ahtrc_0          !:   laplacian diffusivity coefficient for passive tracer [m2/s]
    REAL(wp), PUBLIC ::   rn_bhtrc_0          !: bilaplacian      -          --     -       -   [m4/s]
+   REAL(wp), PUBLIC ::   rn_fact_lap         !: Enhanced zonal diffusivity coefficent in the equatorial domain
    !
    !                      !!: ** lateral mixing namelist (nam_trcldf) **
@@ -64,5 +65,6 @@
       INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
       !
-      INTEGER            :: jn
+      INTEGER            :: ji, jj, jk, jn
+      REAL(wp)           :: zdep
       CHARACTER (len=22) :: charout
       REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zahu, zahv
@@ -76,10 +78,20 @@
          ztrtrd(:,:,:,:)  = tra(:,:,:,:)
       ENDIF
-      !
-      !                                        !* set the lateral diffusivity coef. for passive tracer      
+      !                                  !* set the lateral diffusivity coef. for passive tracer      
       CALL wrk_alloc( jpi,jpj,jpk,   zahu, zahv )
-      zahu(:,:,:) = rldf * ahtu(:,:,:)
+      zahu(:,:,:) = rldf * ahtu(:,:,:) 
       zahv(:,:,:) = rldf * ahtv(:,:,:)
-
+      !                                  !* Enhanced zonal diffusivity coefficent in the equatorial domain
+      DO jk= 1, jpk
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF( gdept_n(ji,jj,jk) > 200. .AND. gphit(ji,jj) < 5. .AND. gphit(ji,jj) > -5. ) THEN
+                  zdep = MAX( gdept_n(ji,jj,jk) - 1000., 0. ) / 1000.
+                  zahu(ji,jj,jk) = zahu(ji,jj,jk) * MAX( 1., rn_fact_lap * EXP( -zdep ) )
+               ENDIF
+            END DO
+         END DO
+      END DO
+      !
       SELECT CASE ( nldf )                     !* compute lateral mixing trend and add it to the general trend
       !
@@ -136,5 +148,5 @@
       NAMELIST/namtrc_ldf/ ln_trcldf_lap, ln_trcldf_blp,                                  &
          &                 ln_trcldf_lev, ln_trcldf_hor, ln_trcldf_iso, ln_trcldf_triad,  &
-         &                 rn_ahtrc_0   , rn_bhtrc_0
+         &                 rn_ahtrc_0   , rn_bhtrc_0, rn_fact_lap  
       !!----------------------------------------------------------------------
       !
@@ -164,4 +176,6 @@
          WRITE(numout,*) '           laplacian                 rn_ahtrc_0      = ', rn_ahtrc_0
          WRITE(numout,*) '         bilaplacian                 rn_bhtrc_0      = ', rn_bhtrc_0
+         WRITE(numout,*) '      enhanced zonal diffusivity     rn_fact_lap     = ', rn_fact_lap
+
       ENDIF
       !      
