Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_BFM/cpp_GYRE_BFM.fcm
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_BFM/cpp_GYRE_BFM.fcm	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_BFM/cpp_GYRE_BFM.fcm	(revision 4617)
@@ -1,2 +1,2 @@
-bld::tool::fppkeys key_dynspg_flt key_ldfslp key_zdftke key_vectopt_loop key_top key_my_trc key_mpp_mpi key_iomput
+bld::tool::fppkeys key_dynspg_flt key_ldfslp key_zdftke key_top key_my_trc key_mpp_mpi key_iomput
 inc $BFMDIR/src/nemo/bfm.fcm
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_PISCES/EXP00/iodef.xml
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_PISCES/EXP00/iodef.xml	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_PISCES/EXP00/iodef.xml	(revision 4617)
@@ -74,5 +74,5 @@
         </file>
 
-	<file id="file5" name="_ptrc_T" description="lobster sms variables" >
+	<file id="file5" name_suffix="_ptrc_T" description="lobster sms variables" >
           <field field_ref="DET"      />
           <field field_ref="ZOO"      />
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_PISCES/EXP00/namelist_cfg
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_PISCES/EXP00/namelist_cfg	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_PISCES/EXP00/namelist_cfg	(revision 4617)
@@ -31,4 +31,8 @@
    ln_zco      = .true.    !  z-coordinate - full    steps   (T/F)      ("key_zco" may also be defined)
    ln_zps      = .false.   !  z-coordinate - partial steps   (T/F)
+/
+!-----------------------------------------------------------------------
+&namzgr_sco    !   s-coordinate or hybrid z-s-coordinate
+!-----------------------------------------------------------------------
 /
 !-----------------------------------------------------------------------
@@ -78,5 +82,5 @@
 !-----------------------------------------------------------------------
    nn_fsbc     = 1         !  frequency of surface boundary condition computation
-                           !     (also = the frequency of sea-ice model call)
+   !                       !     (also = the frequency of sea-ice model call)
    ln_ana      = .true.    !  analytical formulation                    (T => fill namsbc_ana )
    ln_blk_core = .false.   !  CORE bulk formulation                     (T => fill namsbc_core)
@@ -87,4 +91,30 @@
 /
 !-----------------------------------------------------------------------
+&namsbc_ana    !   analytical surface boundary condition
+!-----------------------------------------------------------------------
+   nn_tau000   =   100     !  gently increase the stress over the first ntau_rst time-steps
+   rn_utau0    =   0.1e0   !  uniform value for the i-stress
+/
+!-----------------------------------------------------------------------
+&namsbc_flx    !   surface boundary condition : flux formulation
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namsbc_clio   !   namsbc_clio  CLIO bulk formulae
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namsbc_core   !   namsbc_core  CORE bulk formulae
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namsbc_mfs   !   namsbc_mfs  MFS bulk formulae
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namsbc_cpl    !   coupled ocean/atmosphere model                       ("key_coupled")
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
 &namtra_qsr    !   penetrative solar radiation
 !-----------------------------------------------------------------------
@@ -93,5 +123,24 @@
    nn_chldta   =      0    !  RGB : Chl data (=1) or cst value (=0)
 /
-
+!-----------------------------------------------------------------------
+&namsbc_rnf    !   runoffs namelist surface boundary condition
+!-----------------------------------------------------------------------
+   ln_rnf_mouth = .false.   !  specific treatment at rivers mouths
+/
+!-----------------------------------------------------------------------
+&namsbc_apr    !   Atmospheric pressure used as ocean forcing or in bulk
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namsbc_ssr    !   surface boundary condition : sea surface restoring
+!-----------------------------------------------------------------------
+   nn_sssr     =     0     !  add a damping     term in the surface freshwater flux (=2)
+   rn_deds     =   -27.7   !  magnitude of the damping on salinity   [mm/day]
+   ln_sssr_bnd =   .false. !  flag to bound erp term (associated with nn_sssr=2)
+/
+!-----------------------------------------------------------------------
+&namsbc_alb    !   albedo parameters
+!-----------------------------------------------------------------------
+/
 !-----------------------------------------------------------------------
 &namberg       !   iceberg parameters
@@ -105,4 +154,28 @@
 !-----------------------------------------------------------------------
 &namcla        !   cross land advection
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namobc        !   open boundaries parameters                           ("key_obc")
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namagrif      !  AGRIF zoom                                            ("key_agrif")
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&nam_tide      !    tide parameters (#ifdef key_tide)
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&nambdy        !  unstructured open boundaries                          ("key_bdy")
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&nambdy_dta      !  open boundaries - external data           ("key_bdy")
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&nambdy_tide     ! tidal forcing at open boundaries
 !-----------------------------------------------------------------------
 /
@@ -119,4 +192,8 @@
 /
 !-----------------------------------------------------------------------
+&nambbl        !   bottom boundary layer scheme
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
 &nameos        !   ocean physical parameters
 !-----------------------------------------------------------------------
@@ -127,10 +204,45 @@
 !-----------------------------------------------------------------------
    ln_traadv_msc_ups=  .false.  !  use upstream scheme within muscl 
+/
+!-----------------------------------------------------------------------
+&namtra_adv_mle !  mixed layer eddy parametrisation (Fox-Kemper param)
+!-----------------------------------------------------------------------
 /
 !----------------------------------------------------------------------------------
 &namtra_ldf    !   lateral diffusion scheme for tracers
 !----------------------------------------------------------------------------------
-   rn_aeiv_0        =     0.    !  eddy induced velocity coefficient [m2/s]
-   rn_aht_0         =  1000.    !  horizontal eddy diffusivity for tracers [m2/s]
+   !                       !  Operator type:
+   ln_traldf_lap   =  .true.   !  laplacian operator
+   ln_traldf_blp   =  .false.  !  bilaplacian operator
+   !                       !  Direction of action:
+   ln_traldf_lev   =  .false.  !  iso-level
+   ln_traldf_hor   =  .false.  !  horizontal (geopotential)
+   ln_traldf_iso   =  .false.  !  iso-neutral
+   ln_traldf_triad =  .true.   !  iso-neutral with Griffies triads
+   !
+   !		       	         !  iso-neutral options         
+   ln_traldf_msc   =  .true.   !  Method of Stabilizing Correction (both operator)
+   rn_slpmax       =   0.01    !  slope limit                      (both operator)
+   ln_triad_iso    =  .false.  !  pure horizontal mixing in ML     (triad only)
+   ln_botmix_triad =  .false.  !  lateral mixing on bottom         (triad only)
+   !
+   !                       !  Coefficients
+   nn_aht_ijk_t    =  0        !  space/time variation of eddy coef
+   !                                   !   =0 constant ; =10 F(k) ; =20 F(i,j)=F(grid spacing) ; =30 F(i,j,k) 
+   !                                   !   =21 F(i,jt)=Treguier et al. JPO 1997 formulation
+   !                                   !   =-20 (=-30) read in eddy_induced_velocity_2D.nc (..._3D.nc) file
+   rn_aht_0        = 1000.     !  lateral eddy diffusivity   (lap. operator) [m2/s]
+   rn_bht_0        = 5.e+11    !  lateral eddy diffusivity (bilap. operator) [m4/s]
+/
+!----------------------------------------------------------------------------------
+&namtra_ldfeiv !   eddy induced velocity param.
+!----------------------------------------------------------------------------------
+   ln_ldfeiv     =.false.  ! use eddy induced velocity parameterization
+   ln_ldfeiv_dia =.false.  ! diagnose eiv stream function and velocities
+   rn_aeiv_0     = 1000.   ! eddy induced velocity coefficient   [m2/s]
+   nn_aei_ijk_t  =  0      ! space/time variation of the eiv coeficient
+   !                              !   =0 constant ; =10 F(k) ; =20 F(i,j) = F(grid spacing) ; =30 F(i,j,k) 
+   !                              !   =21 F(i,jt)=Treguier et al. JPO 1997 formulation
+   !                              !   =-20 (=-30) read in eddy_induced_velocity_2D.nc (..._3D.nc) file
 /
 !-----------------------------------------------------------------------
@@ -159,7 +271,24 @@
 /
 !-----------------------------------------------------------------------
+!namdyn_spg    !   surface pressure gradient   (CPP key only)
+!-----------------------------------------------------------------------
+
+!-----------------------------------------------------------------------
 &namdyn_ldf    !   lateral diffusion on momentum
 !-----------------------------------------------------------------------
-   rn_ahm_0_lap     = 100000.   !  horizontal laplacian eddy viscosity   [m2/s]
+   !                       !  Type of the operator :
+   ln_dynldf_lap =  .true.     !  laplacian operator
+   ln_dynldf_blp =  .false.    !  bilaplacian operator
+   !                       !  Direction of action  :
+   ln_dynldf_lev =  .false.    !  iso-level
+   ln_dynldf_hor =  .true.     !  horizontal (geopotential)
+   ln_dynldf_iso =  .false.    !  iso-neutral
+   !                       !  Coefficient
+   nn_ahm_ijk_t  = 0           !  space/time variation of eddy coef
+   !                                !   =0 constant ; =10 F(k) ; =20 F(i,j)=F(grid spacing) ; =30 F(i,j,k) 
+   !                                !   =-20 (=-30) read in eddy_induced_velocity_2D.nc (..._3D.nc) file
+   rn_ahm_0      = 100000.     !  horizontal laplacian eddy viscosity   [m2/s]
+   rn_ahm_b      =      0.     !  background eddy viscosity for ldf_iso [m2/s]
+   rn_bhm_0      = 1.e+12      !  horizontal bilaplacian eddy viscosity [m4/s]
 /
 !-----------------------------------------------------------------------
@@ -169,7 +298,28 @@
 /
 !-----------------------------------------------------------------------
+&namzdf_ric    !   richardson number dependent vertical diffusion       ("key_zdfric" )
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
 &namzdf_tke    !   turbulent eddy kinetic dependent vertical diffusion  ("key_zdftke")
 !-----------------------------------------------------------------------
    nn_etau     =   0       !  penetration of tke below the mixed layer (ML) due to internal & intertial waves
+/
+!------------------------------------------------------------------------
+&namzdf_kpp    !   K-Profile Parameterization dependent vertical mixing  ("key_zdfkpp", and optionally:
+!------------------------------------------------------------------------ "key_kppcustom" or "key_kpplktb")
+/
+!-----------------------------------------------------------------------
+&namzdf_gls                !   GLS vertical diffusion                   ("key_zdfgls")
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namzdf_ddm    !   double diffusive mixing parameterization             ("key_zdfddm")
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namzdf_tmx    !   tidal mixing parameterization                        ("key_zdftmx")
+!-----------------------------------------------------------------------
+   ln_tmx_itf  = .false.   !  ITF specific parameterisation
 /
 !-----------------------------------------------------------------------
@@ -189,4 +339,17 @@
 /
 !-----------------------------------------------------------------------
+&namnc4        !   netcdf4 chunking and compression settings            ("key_netcdf4")
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namtrd        !   diagnostics on dynamics and/or tracer trends         ("key_trddyn" and/or "key_trdtra")
+!              !       or mixed-layer trends or barotropic vorticity    ("key_trdmld" or     "key_trdvor")
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namflo       !   float parameters                                      ("key_float")
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
 &namptr       !   Poleward Transport Diagnostic
 !-----------------------------------------------------------------------
@@ -197,4 +360,23 @@
 /
 !-----------------------------------------------------------------------
+&namdct        ! transports through sections
+!-----------------------------------------------------------------------
+    nn_dct      = 60       !  time step frequency for transports computing
+    nn_dctwri   = 60       !  time step frequency for transports writing
+    nn_secdebug = 0        !      0 : no section to debug
+/
+!-----------------------------------------------------------------------
+&namobs       !  observation usage switch                               ('key_diaobs')
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&nam_asminc   !   assimilation increments                               ('key_asminc')
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
+&namsbc_wave   ! External fields from wave model
+!-----------------------------------------------------------------------
+/
+!-----------------------------------------------------------------------
 &namdyn_nept  !   Neptune effect (simplified: lateral and vertical diffusions removed)
 !-----------------------------------------------------------------------
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_PISCES/cpp_GYRE_PISCES.fcm
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_PISCES/cpp_GYRE_PISCES.fcm	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/CONFIG/GYRE_PISCES/cpp_GYRE_PISCES.fcm	(revision 4617)
@@ -1,1 +1,1 @@
-bld::tool::fppkeys  key_dynspg_flt key_ldfslp key_zdftke key_top key_pisces_reduced key_iomput key_mpp_mpi 
+bld::tool::fppkeys  key_dynspg_flt key_zdftke key_top key_pisces_reduced key_iomput key_mpp_mpi
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zbio.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zbio.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zbio.F90	(revision 4617)
@@ -66,5 +66,4 @@
    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
-
 CONTAINS
 
@@ -90,7 +89,6 @@
       !!              for passive tracers are saved for futher diagnostics.
       !!---------------------------------------------------------------------
-      !!
       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index      
-      !!
+      !
       INTEGER  ::   ji, jj, jk, jl
       REAL(wp) ::   zdet, zzoo, zphy, zno3, znh4, zdom      ! now concentrations
@@ -476,4 +474,5 @@
    END SUBROUTINE p2z_bio
 
+
    SUBROUTINE p2z_bio_init
       !!----------------------------------------------------------------------
@@ -485,4 +484,6 @@
       !!
       !!----------------------------------------------------------------------
+      INTEGER ::   ios   ! local integer 
+      !!
       NAMELIST/namlobphy/ tmumax, rgamma, fphylab, tmminp, aki
       NAMELIST/namlobnut/ akno3, aknh4, taunn, psinut
@@ -490,5 +491,4 @@
       NAMELIST/namlobdet/  taudn, fdetlab
       NAMELIST/namlobdom/ taudomn
-      INTEGER :: ios                 ! Local integer output status for namelist read
       !!----------------------------------------------------------------------
 
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zexp.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zexp.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zexp.F90	(revision 4617)
@@ -48,5 +48,4 @@
    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
-
 CONTAINS
 
@@ -64,7 +63,6 @@
       !!              COLUMN BELOW THE SURFACE LAYER.
       !!---------------------------------------------------------------------
-      !!
       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index      
-      !!
+      !
       INTEGER  ::   ji, jj, jk, jl, ikt
       REAL(wp) ::   zgeolpoc, zfact, zwork, ze3t, zsedpocd, zmaskt
@@ -77,7 +75,7 @@
       !
       IF( kt == nittrc000 )   CALL p2z_exp_init
-
+      !
       CALL wrk_alloc( jpi, jpj, zsedpoca ) 
-      zsedpoca(:,:) = 0.
+      zsedpoca(:,:) = 0._wp
 
       IF( l_trdtrc )  THEN
@@ -102,5 +100,4 @@
       ! Find the last level of the water column
       ! Compute fluxes due to sinking particles (slow)
-   
 
       zgeolpoc = 0.e0         !     Initialization
@@ -180,4 +177,5 @@
    END SUBROUTINE p2z_exp
 
+
    SUBROUTINE p2z_exp_init
       !!----------------------------------------------------------------------
@@ -272,4 +270,5 @@
    END SUBROUTINE p2z_exp_init
 
+
    INTEGER FUNCTION p2z_exp_alloc()
       !!----------------------------------------------------------------------
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zopt.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zopt.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zopt.F90	(revision 4617)
@@ -47,5 +47,4 @@
    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
-
 CONTAINS
 
@@ -64,7 +63,6 @@
       !!                neln   ???
       !!---------------------------------------------------------------------
-      !!
       INTEGER, INTENT( in ) ::   kt   ! index of the time stepping
-      !!
+      !
       INTEGER  ::   ji, jj, jk          ! dummy loop indices
       CHARACTER (len=25) ::   charout   ! temporary character
@@ -153,4 +151,5 @@
    END SUBROUTINE p2z_opt
 
+
    SUBROUTINE p2z_opt_init
       !!----------------------------------------------------------------------
@@ -164,16 +163,16 @@
       NAMELIST/namlobopt/ xkg0, xkr0, xkgp, xkrp, xlg, xlr, rpig
       NAMELIST/namlobrat/ rcchl, redf, reddom
-      INTEGER :: ios                 ! Local integer output status for namelist read
+      INTEGER ::   ios   ! local integer
       !!----------------------------------------------------------------------
-
+      !
       REWIND( numnatp_ref )              ! Namelist namlobopt in reference namelist : Lobster options
       READ  ( numnatp_ref, namlobopt, IOSTAT = ios, ERR = 901)
 901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namlobopt in reference namelist', lwp )
-
+      !
       REWIND( numnatp_cfg )              ! Namelist namlobopt in configuration namelist : Lobster options
       READ  ( numnatp_cfg, namlobopt, IOSTAT = ios, ERR = 902 )
 902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namlobopt in configuration namelist', lwp )
       WRITE ( numonp, namlobopt )
-
+      !
       IF(lwp) THEN
          WRITE(numout,*)
@@ -192,10 +191,10 @@
       READ  ( numnatp_ref, namlobrat, IOSTAT = ios, ERR = 903)
 903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namlobrat in reference namelist', lwp )
-
+      !
       REWIND( numnatp_cfg )              ! Namelist namlobrat in configuration namelist : Lobster ratios
       READ  ( numnatp_cfg, namlobrat, IOSTAT = ios, ERR = 904 )
 904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namlobrat in configuration namelist', lwp )
       WRITE ( numonp, namlobrat )
-
+      !
       IF(lwp) THEN
           WRITE(numout,*) ' Namelist namlobrat'
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zsed.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zsed.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zsed.F90	(revision 4617)
@@ -41,5 +41,4 @@
    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
-
 CONTAINS
 
@@ -62,5 +61,4 @@
       !!              trend of passive tracers is saved for futher diagnostics.
       !!---------------------------------------------------------------------
-      !!
       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index      
       !!
@@ -78,8 +76,7 @@
          IF(lwp) WRITE(numout,*) ' ~~~~~~~'
       ENDIF
-
-      ! Allocate temporary workspace
-                       CALL wrk_alloc( jpi, jpj, jpk, zwork, ztra )
-      IF( ln_diatrc )  CALL wrk_alloc( jpi, jpj, zw2d )
+      !
+      CALL wrk_alloc( jpi, jpj, jpk, zwork, ztra )
+      !
       IF( l_trdtrc ) THEN
          CALL wrk_alloc( jpi, jpj, jpk, ztrbio )
@@ -110,5 +107,7 @@
 
       IF( ln_diatrc ) THEN 
-         DO jk = 1, jpkm1
+         CALL wrk_alloc( jpi, jpj, zw2d )
+         zw2d(:,:) = ztra(:,:,1) * fse3t(:,:,1) * 86400.
+         DO jk = 2, jpkm1
             zw2d(:,:) = zw2d(:,:) + ztra(:,:,jk) * fse3t(:,:,jk) * 86400.
          END DO
@@ -121,5 +120,5 @@
       ENDIF
       !
-      IF( ln_diabio )  trbio(:,:,:,jp_pcs0_trd + 7) = ztra(:,:,:)
+      IF( ln_diabio .AND. .NOT.lk_iomput )  trbio(:,:,:,jp_pcs0_trd + 7) = ztra(:,:,:)
       CALL wrk_dealloc( jpi, jpj, jpk, zwork, ztra )
       !
@@ -130,5 +129,5 @@
          CALL wrk_dealloc( jpi, jpj, jpk, ztrbio )
       ENDIF
-
+      !
       IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
          WRITE(charout, FMT="('sed')")
@@ -141,4 +140,5 @@
    END SUBROUTINE p2z_sed
 
+
    SUBROUTINE p2z_sed_init
       !!----------------------------------------------------------------------
@@ -150,16 +150,18 @@
       !!
       !!----------------------------------------------------------------------
+      INTEGER ::   ios   ! Local integer
+      !
       NAMELIST/namlobsed/ sedlam, sedlostpoc, vsed, xhr
-      INTEGER :: ios                 ! Local integer output status for namelist read
-
+      !!----------------------------------------------------------------------
+      !
       REWIND( numnatp_ref )              ! Namelist namlobsed in reference namelist : Lobster sediments
       READ  ( numnatp_ref, namlobsed, IOSTAT = ios, ERR = 901)
 901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namlosed in reference namelist', lwp )
-
+      !
       REWIND( numnatp_cfg )              ! Namelist namlobsed in configuration namelist : Lobster sediments
       READ  ( numnatp_cfg, namlobsed, IOSTAT = ios, ERR = 902 )
 902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namlobsed in configuration namelist', lwp )
       WRITE ( numonp, namlobsed )
-
+      !
       IF(lwp) THEN
           WRITE(numout,*) ' Namelist namlobsed'
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zsms.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zsms.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P2Z/p2zsms.F90	(revision 4617)
@@ -58,5 +58,5 @@
       CALL p2z_sed( kt )      ! sedimentation model
       CALL p2z_exp( kt )      ! export
-
+      !
       IF( l_trdtrc ) THEN
          DO jn = jp_pcs0, jp_pcs1
@@ -64,9 +64,9 @@
          END DO
       END IF
-
+      !
       IF( lk_trdmld_trc )  CALL trd_mld_bio( kt )   ! trends: Mixed-layer
       !
-      IF ( kt == nittrc000 ) CALL FLUSH    ( numonp )     ! flush output namelist PISCES
-      IF( nn_timing == 1 )  CALL timing_stop('p2z_sms')
+      IF( kt == nittrc000 )   CALL FLUSH( numonp )       ! flush output namelist PISCES
+      IF( nn_timing == 1  )   CALL timing_stop('p2z_sms')
       !
    END SUBROUTINE p2z_sms
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zche.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zche.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zche.F90	(revision 4617)
@@ -195,7 +195,5 @@
       ! CHEMICAL CONSTANTS - SURFACE LAYER
       ! ----------------------------------
-!CDIR NOVERRCHK
       DO jj = 1, jpj
-!CDIR NOVERRCHK
          DO ji = 1, jpi
             !                             ! SET ABSOLUTE TEMPERATURE
@@ -227,9 +225,6 @@
       ! OXYGEN SOLUBILITY - DEEP OCEAN
       ! -------------------------------
-!CDIR NOVERRCHK
       DO jk = 1, jpk
-!CDIR NOVERRCHK
          DO jj = 1, jpj
-!CDIR NOVERRCHK
             DO ji = 1, jpi
               ztkel = tsn(ji,jj,jk,jp_tem) + 273.16
@@ -252,9 +247,6 @@
       ! CHEMICAL CONSTANTS - DEEP OCEAN
       ! -------------------------------
-!CDIR NOVERRCHK
       DO jk = 1, jpk
-!CDIR NOVERRCHK
          DO jj = 1, jpj
-!CDIR NOVERRCHK
             DO ji = 1, jpi
 
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zfechem.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zfechem.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zfechem.F90	(revision 4617)
@@ -113,9 +113,6 @@
          ! Chemistry is supposed to be fast enough to be at equilibrium
          ! ------------------------------------------------------------
-!CDIR NOVERRCHK
          DO jk = 1, jpkm1
-!CDIR NOVERRCHK
             DO jj = 1, jpj
-!CDIR NOVERRCHK
                DO ji = 1, jpi
                   ! Calculate ligand concentrations : assume 2/3rd of excess goes to
@@ -195,9 +192,6 @@
          ! Chemistry is supposed to be fast enough to be at equilibrium
          ! ------------------------------------------------------------
-!CDIR NOVERRCHK
          DO jk = 1, jpkm1
-!CDIR NOVERRCHK
             DO jj = 1, jpj
-!CDIR NOVERRCHK
                DO ji = 1, jpi
                   zTL1(ji,jj,jk) = ztotlig(ji,jj,jk)
@@ -217,10 +211,7 @@
       ENDIF
 
-      zdust = 0.         ! if no dust available
-!CDIR NOVERRCHK
+      zdust = 0._wp         ! if no dust available
       DO jk = 1, jpkm1
-!CDIR NOVERRCHK
          DO jj = 1, jpj
-!CDIR NOVERRCHK
             DO ji = 1, jpi
                zstep = xstep
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zflx.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zflx.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zflx.F90	(revision 4617)
@@ -122,7 +122,5 @@
 
       DO jm = 1, 10
-!CDIR NOVERRCHK
          DO jj = 1, jpj
-!CDIR NOVERRCHK
             DO ji = 1, jpi
 
@@ -155,7 +153,5 @@
       ! -------------------------------------------
 
-!CDIR NOVERRCHK
       DO jj = 1, jpj
-!CDIR NOVERRCHK
          DO ji = 1, jpi
             ztc  = MIN( 35., tsn(ji,jj,1,jp_tem) )
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zlys.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zlys.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zlys.F90	(revision 4617)
@@ -81,9 +81,6 @@
       DO jn = 1, 5                               !  BEGIN OF ITERATION
          !
-!CDIR NOVERRCHK
          DO jk = 1, jpkm1
-!CDIR NOVERRCHK
             DO jj = 1, jpj
-!CDIR NOVERRCHK
                DO ji = 1, jpi
                   zfact = rhop(ji,jj,jk) / 1000. + rtrn
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zopt.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zopt.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zopt.F90	(revision 4617)
@@ -90,10 +90,7 @@
       ze2(:,:,jpk) = 0._wp
       ze3(:,:,jpk) = 0._wp
-
       !                                        !* attenuation coef. function of Chlorophyll and wavelength (Red-Green-Blue)
       DO jk = 1, jpkm1                         !  --------------------------------------------------------
-!CDIR NOVERRCHK
-         DO jj = 1, jpj
-!CDIR NOVERRCHK
+         DO jj = 1, jpj
             DO ji = 1, jpi
                zchl = ( trn(ji,jj,jk,jpnch) + trn(ji,jj,jk,jpdch) + rtrn ) * 1.e6
@@ -107,9 +104,7 @@
          END DO
       END DO
-
-
+      !
       !                                        !* Photosynthetically Available Radiation (PAR)
       !                                        !  --------------------------------------
-
       IF( ln_varpar ) THEN
          ze1(:,:,1) = par_varsw(:,:) * qsr(:,:) * EXP( -0.5 * zekb(:,:,1) )
@@ -122,7 +117,5 @@
       ENDIF
 
-!CDIR NOVERRCHK
       DO jj = 1, jpj
-!CDIR NOVERRCHK
          DO ji = 1, jpi
             zc1 = ze1(ji,jj,1)
@@ -135,9 +128,6 @@
       END DO
 
-    
       DO jk = 2, nksrp      
-!CDIR NOVERRCHK
-         DO jj = 1, jpj
-!CDIR NOVERRCHK
+         DO jj = 1, jpj
             DO ji = 1, jpi
                zc1 = ze1(ji,jj,jk-1) * EXP( -0.5 * ( zekb(ji,jj,jk-1) + zekb(ji,jj,jk) ) )
@@ -175,7 +165,5 @@
          !
          DO jk = 2, nksrp + 1
-!CDIR NOVERRCHK
             DO jj = 1, jpj
-!CDIR NOVERRCHK
                DO ji = 1, jpi
                   zc0 = ze0(ji,jj,jk-1) * EXP( -fse3t(ji,jj,jk-1) * zxsi0r )
@@ -221,7 +209,5 @@
 
       DO jk = 1, nksrp
-!CDIR NOVERRCHK
-         DO jj = 1, jpj
-!CDIR NOVERRCHK
+         DO jj = 1, jpj
             DO ji = 1, jpi
                IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
@@ -238,7 +224,5 @@
       !
       DO jk = 1, nksrp
-!CDIR NOVERRCHK
-         DO jj = 1, jpj
-!CDIR NOVERRCHK
+         DO jj = 1, jpj
             DO ji = 1, jpi
                IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
@@ -273,4 +257,5 @@
    END SUBROUTINE p4z_opt
 
+
    SUBROUTINE p4z_optsbc( kt )
       !!----------------------------------------------------------------------
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zprod.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zprod.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zprod.F90	(revision 4617)
@@ -31,5 +31,4 @@
    PUBLIC   p4z_prod_alloc
 
-   !! * Shared module variables
    LOGICAL , PUBLIC ::  ln_newprod      !:
    REAL(wp), PUBLIC ::  pislope         !:
@@ -150,9 +149,6 @@
 
       IF( ln_newprod ) THEN
-!CDIR NOVERRCHK
          DO jk = 1, jpkm1
-!CDIR NOVERRCHK
             DO jj = 1, jpj
-!CDIR NOVERRCHK
                DO ji = 1, jpi
                   ! Computation of the P-I slope for nanos and diatoms
@@ -188,9 +184,6 @@
          END DO
       ELSE
-!CDIR NOVERRCHK
          DO jk = 1, jpkm1
-!CDIR NOVERRCHK
             DO jj = 1, jpj
-!CDIR NOVERRCHK
                DO ji = 1, jpi
 
@@ -231,9 +224,6 @@
       !  Computation of a proxy of the N/C ratio
       !  ---------------------------------------
-!CDIR NOVERRCHK
       DO jk = 1, jpkm1
-!CDIR NOVERRCHK
          DO jj = 1, jpj
-!CDIR NOVERRCHK
             DO ji = 1, jpi
                 zval = MIN( xnanopo4(ji,jj,jk), ( xnanonh4(ji,jj,jk) + xnanono3(ji,jj,jk) ) )   &
@@ -246,11 +236,9 @@
          END DO
       END DO
-
-
       DO jk = 1, jpkm1
          DO jj = 1, jpj
             DO ji = 1, jpi
-
-                IF( etot(ji,jj,jk) > 1.E-3 ) THEN
+               !
+               IF( etot(ji,jj,jk) > 1.E-3 ) THEN
                    !    Si/C of diatoms
                    !    ------------------------
@@ -268,5 +256,5 @@
                   ENDIF
                   zysopt(ji,jj,jk) = grosip * zlim * zsilfac * zsilfac2
-              ENDIF
+               ENDIF
             END DO
          END DO
@@ -296,9 +284,6 @@
 
       ! Computation of the various production terms 
-!CDIR NOVERRCHK
       DO jk = 1, jpkm1
-!CDIR NOVERRCHK
          DO jj = 1, jpj
-!CDIR NOVERRCHK
             DO ji = 1, jpi
                IF( etot(ji,jj,jk) > 1.E-3 ) THEN
@@ -331,9 +316,6 @@
 
       IF( ln_newprod ) THEN
-!CDIR NOVERRCHK
          DO jk = 1, jpkm1
-!CDIR NOVERRCHK
             DO jj = 1, jpj
-!CDIR NOVERRCHK
                DO ji = 1, jpi
                   IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
@@ -359,9 +341,6 @@
          END DO
       ELSE
-!CDIR NOVERRCHK
          DO jk = 1, jpkm1
-!CDIR NOVERRCHK
             DO jj = 1, jpj
-!CDIR NOVERRCHK
                DO ji = 1, jpi
                   IF( etot(ji,jj,jk) > 1.E-3 ) THEN
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90	(revision 4617)
@@ -15,5 +15,5 @@
    USE oce_trc         !  shared variables between ocean and passive tracers
    USE trc             !  passive tracers common variables 
-   USE trcdta
+   USE trcdta          !
    USE sms_pisces      !  PISCES Source Minus Sink variables
    USE p4zbio          !  Biological model
@@ -24,9 +24,10 @@
    USE p4zsed          !  Sedimentation
    USE p4zint          !  time interpolation
-   USE iom             !  I/O manager
    USE trdmod_oce      !  Ocean trends variables
    USE trdmod_trc      !  TOP trends variables
    USE sedmodel        !  Sediment model
+   !
    USE prtctl_trc      !  print control for debugging
+   USE iom             !  I/O manager
 
    IMPLICIT NONE
@@ -34,5 +35,5 @@
 
    PUBLIC   p4z_sms_init    ! called in p4zsms.F90
-   PUBLIC   p4z_sms    ! called in p4zsms.F90
+   PUBLIC   p4z_sms         ! called in p4zsms.F90
 
    REAL(wp) :: alkbudget, no3budget, silbudget, ferbudget
@@ -58,5 +59,4 @@
       !!              - ...
       !!---------------------------------------------------------------------
-      !
       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index      
       !!
@@ -73,5 +73,5 @@
             jl = jn + jp_pcs0 - 1
             ztrdpis(:,:,:,jn) = trn(:,:,:,jl)
-         ENDDO
+         END DO
       ENDIF
       !
@@ -91,16 +91,15 @@
          !
          ndayflxtr = nday_year
-
+         !
          IF(lwp) write(numout,*)
          IF(lwp) write(numout,*) ' New chemical constants and various rates for biogeochemistry at new day : ', nday_year
          IF(lwp) write(numout,*) '~~~~~~'
-
+         !
          CALL p4z_che              ! computation of chemical constants
          CALL p4z_int( kt )        ! computation of various rates for biogeochemistry
          !
       ENDIF
-
-      IF( ll_sbc ) CALL p4z_sbc( kt )   ! external sources of nutrients 
-
+      IF( ll_sbc )   CALL p4z_sbc( kt )   ! external sources of nutrients 
+      !
       DO jnt = 1, nrdttrc          ! Potential time splitting if requested
          !
@@ -110,5 +109,5 @@
          DO jn = jp_pcs0, jp_pcs1
             trb(:,:,:,jn) = trn(:,:,:,jn)
-         ENDDO
+         END DO
          !
       END DO
@@ -118,5 +117,5 @@
             jl = jn + jp_pcs0 - 1
             ztrdpis(:,:,:,jn) = ( ztrdpis(:,:,:,jn) - trn(:,:,:,jl) ) * rfact2r
-         ENDDO
+         END DO
       ENDIF
 
@@ -152,10 +151,10 @@
       !
       CALL p4z_chk_mass( kt ) ! Mass conservation checking
-
-      IF ( kt == nittrc000 ) CALL FLUSH    ( numonp )     ! flush output namelist PISCES
-      IF( nn_timing == 1 )  CALL timing_stop('p4z_sms')
-      !
+      !
+      IF( kt == nittrc000 )   CALL FLUSH    ( numonp )     ! flush output namelist PISCES
+      IF( nn_timing == 1  )   CALL timing_stop('p4z_sms')
       !
    END SUBROUTINE p4z_sms
+
 
    SUBROUTINE p4z_sms_init
@@ -260,6 +259,7 @@
          WRITE(numout,*) '    Flag to check mass conservation of NO3/Si/TALK ln_check_mass = ', ln_check_mass
       ENDIF
-
+      !
    END SUBROUTINE p4z_sms_init
+
 
    SUBROUTINE p4z_ph_ini
@@ -291,4 +291,5 @@
    END SUBROUTINE p4z_ph_ini
 
+
    SUBROUTINE p4z_rst( kt, cdrw )
       !!---------------------------------------------------------------------
@@ -341,4 +342,5 @@
    END SUBROUTINE p4z_rst
 
+
    SUBROUTINE p4z_dmp( kt )
       !!----------------------------------------------------------------------
@@ -347,6 +349,5 @@
       !! ** purpose  : Relaxation of some tracers
       !!----------------------------------------------------------------------
-      !
-      INTEGER, INTENT( in )  ::     kt ! time step
+      INTEGER, INTENT(in) ::   kt   ! time step
       !
       REAL(wp) ::  alkmean = 2426.     ! mean value of alkalinity ( Glodap ; for Goyet 2391. )
@@ -357,5 +358,4 @@
       REAL(wp) :: zarea, zalksum, zpo4sum, zno3sum, zsilsum
       !!---------------------------------------------------------------------
-
 
       IF(lwp)  WRITE(numout,*)
@@ -386,5 +386,5 @@
          !
       ENDIF
-
+      !
    END SUBROUTINE p4z_dmp
 
@@ -397,9 +397,7 @@
       !!
       !!---------------------------------------------------------------------
-      !
-      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index      
-      !!
-      !!---------------------------------------------------------------------
-
+      INTEGER, INTENT(in) ::   kt   ! ocean time-step index      
+      !!---------------------------------------------------------------------
+      !
       IF( kt == nittrc000 ) THEN 
          IF( ln_check_mass .AND. lwp) THEN      !   Open budget file of NO3, ALK, Si, Fer
@@ -434,6 +432,5 @@
             &                    + trn(:,:,:,jpzoo) * ferat3            &
             &                    + trn(:,:,:,jpmes) * ferat3            ) * cvol(:,:,:)  )
-
-         !
+            !
          t_atm_co2_flx  = t_atm_co2_flx / glob_sum( e1e2t(:,:) )
          t_oce_co2_flx  = t_oce_co2_flx         * 12. / 1.e15 * (-1 )
@@ -452,8 +449,8 @@
          !
       ENDIF
-       !
- 9000  FORMAT(i8,f10.5,e18.10,f10.5,f10.5)
- 9500  FORMAT(i8,4e18.10)     
-       !
+      !
+ 9000 FORMAT(i8,f10.5,e18.10,f10.5,f10.5)
+ 9500 FORMAT(i8,4e18.10)     
+      !
    END SUBROUTINE p4z_chk_mass
 
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/TRP/trcadv.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/TRP/trcadv.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/TRP/trcadv.F90	(revision 4617)
@@ -23,6 +23,6 @@
    USE traadv_ubs      ! UBS      scheme           (tra_adv_ubs    routine)
    USE traadv_qck      ! QUICKEST scheme           (tra_adv_qck    routine)
-   USE traadv_eiv      ! eddy induced velocity     (tra_adv_eiv    routine)
-   USE ldftra_oce      ! lateral diffusion coefficient on tracers
+   USE ldftra          ! lateral diffusion: eddy diffusivity & EIV coeff.
+   USE ldfslp          ! Lateral diffusion: slopes of neutral surfaces
    USE prtctl_trc      ! Print control
 
@@ -41,5 +41,5 @@
 #  include "vectopt_loop_substitute.h90"
    !!----------------------------------------------------------------------
-   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
+   !! NEMO/TOP 3.7 , NEMO Consortium (2014)
    !! $Id$ 
    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
@@ -51,9 +51,6 @@
       !!                  ***  ROUTINE trc_adv_alloc  ***
       !!----------------------------------------------------------------------
-
       ALLOCATE( r2dt(jpk), STAT=trc_adv_alloc )
-
-      IF( trc_adv_alloc /= 0 ) CALL ctl_warn('trc_adv_alloc : failed to allocate array.')
-
+      IF( trc_adv_alloc /= 0 )   CALL ctl_warn('trc_adv_alloc : failed to allocate array.')
    END FUNCTION trc_adv_alloc
 
@@ -67,10 +64,9 @@
       !! ** Method  : - Update the tracer with the advection term following nadv
       !!----------------------------------------------------------------------
-      !!
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
       !
       INTEGER ::   jk 
       CHARACTER (len=22) ::   charout
-      REAL(wp), POINTER, DIMENSION(:,:,:) :: zun, zvn, zwn  ! effective velocity
+      REAL(wp), POINTER, DIMENSION(:,:,:) ::   zun, zvn, zwn   ! effective transport
       !!----------------------------------------------------------------------
       !
@@ -79,8 +75,7 @@
       CALL wrk_alloc( jpi, jpj, jpk, zun, zvn, zwn )
       !
-
       IF( kt == nittrc000 )   CALL trc_adv_ctl          ! initialisation & control of options
-
-      IF( ln_top_euler) THEN
+      !
+      IF( ln_top_euler ) THEN
          r2dt(:) =  rdttrc(:)              ! = rdttrc (use Euler time stepping)
       ELSE
@@ -92,18 +87,26 @@
       ENDIF
 
-      !                                                   ! effective transport
+      !                                               !==  effective transport  ==!
       DO jk = 1, jpkm1
-         !                                                ! eulerian transport only
-         zun(:,:,jk) = e2u  (:,:) * fse3u(:,:,jk) * un(:,:,jk)
+         zun(:,:,jk) = e2u  (:,:) * fse3u(:,:,jk) * un(:,:,jk)                   ! eulerian transport
          zvn(:,:,jk) = e1v  (:,:) * fse3v(:,:,jk) * vn(:,:,jk)
          zwn(:,:,jk) = e1e2t(:,:)                 * wn(:,:,jk)
-         !
       END DO
-      zwn(:,:,jpk) = 0.e0                                 ! no transport trough the bottom
-
-      IF( lk_traldf_eiv .AND. .NOT. ln_traldf_grif )   &  ! add the eiv transport (if necessary)
-         &              CALL tra_adv_eiv( kt, nittrc000, zun, zvn, zwn, 'TRC' )
+!!gm      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN        !!gm to be added like on T and S ????
+!!gm         zun(:,:,:) = zun(:,:,:) + un_td(:,:,:)
+!!gm         zvn(:,:,:) = zvn(:,:,:) + vn_td(:,:,:)
+!!gm      ENDIF
+      zun(:,:,jpk) = 0._wp                                                       ! no transport in the bottom
+      zvn(:,:,jpk) = 0._wp                                                       ! no transport in the bottom
+      zwn(:,:,jpk) = 0.e0                                                        ! no transport trough the bottom
+
+      IF( ln_ldfeiv .AND. .NOT. ln_traldf_triad )   &       
+         &              CALL ldf_eiv_trp( kt, nittrc000, zun, zvn, zwn, 'TRC' )  ! add the eiv transport (if necessary)
+      !
+!!gm  to be added like in T and S ?
+!!gm      IF( ln_mle    )   CALL tra_adv_mle( kt, nittrc000, zun, zvn, zwn, 'TRC' )  ! add the mle transport (if necessary)
       !
       SELECT CASE ( nadv )                            !==  compute advection trend and add it to general trend  ==!
+      !
       CASE ( 1 )   ;    CALL tra_adv_cen2  ( kt, nittrc000, 'TRC',       zun, zvn, zwn, trb, trn, tra, jptra )   !  2nd order centered
       CASE ( 2 )   ;    CALL tra_adv_tvd   ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra )   !  TVD 
@@ -132,5 +135,4 @@
          WRITE(charout, FMT="('adv6')")  ; CALL prt_ctl_trc_info(charout)
                                            CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd')
-         !
       END SELECT
 
@@ -157,5 +159,5 @@
       INTEGER ::   ioptio
       !!----------------------------------------------------------------------
-
+      !
       ioptio = 0                      ! Parameter control
       IF( ln_trcadv_cen2   )   ioptio = ioptio + 1
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90	(revision 4617)
@@ -4,27 +4,26 @@
    !! Ocean Passive tracers : lateral diffusive trends
    !!=====================================================================
-   !! History :  9.0  ! 2005-11 (G. Madec)  Original code
-   !!       NEMO 3.0  ! 2008-01  (C. Ethe, G. Madec)  merge TRC-TRA
+   !! History :  1.0  ! 2005-11  (G. Madec)  Original code
+   !!            3.0  ! 2008-01  (C. Ethe, G. Madec)  merge TRC-TRA
+   !!            3.7  ! 2014-03  (G. Madec)  LDF simplification
    !!----------------------------------------------------------------------
 #if defined key_top
    !!----------------------------------------------------------------------
    !!   'key_top'                                                TOP models
-   !!----------------------------------------------------------------------
    !!----------------------------------------------------------------------
    !!   trc_ldf     : update the tracer trend with the lateral diffusion
    !!       ldf_ctl : initialization, namelist read, and parameters control
    !!----------------------------------------------------------------------
-   USE oce_trc         ! ocean dynamics and active tracers
-   USE trc             ! ocean passive tracers variables
-   USE trcnam_trp      ! passive tracers transport namelist variables
-   USE ldftra_oce      ! lateral diffusion coefficient on tracers
-   USE ldfslp          ! ???
-   USE traldf_bilapg   ! lateral mixing            (tra_ldf_bilapg routine)
-   USE traldf_bilap    ! lateral mixing            (tra_ldf_bilap routine)
-   USE traldf_iso      ! lateral mixing            (tra_ldf_iso routine)
-   USE traldf_iso_grif ! lateral mixing          (tra_ldf_iso_grif routine)
-   USE traldf_lap      ! lateral mixing            (tra_ldf_lap routine)
+   USE trc              ! ocean passive tracers variables
+   USE oce_trc          ! ocean dynamics and active tracers
+   USE trcnam_trp       ! passive tracers transport namelist variables
+   USE ldftra           ! lateral diffusion: eddy diffusivity & EIV coeff.
+   USE ldfslp           ! lateral diffusion: iso-neutral slope 
+   USE traldf_lap       ! lateral diffusion: iso-level lap/blp (tra_ldf_lap/_blp  routine)
+   USE traldf_iso       ! lateral diffusion (Madec operator)   (tra_ldf_iso       routine)
+   USE traldf_iso_triad ! lateral diffusion (triad operator)   (tra_ldf_iso_triad routine)
    USE trdmod_oce
    USE trdtra
+   !
    USE prtctl_trc      ! Print control
 
@@ -40,9 +39,8 @@
 #  include "vectopt_loop_substitute.h90"
    !!----------------------------------------------------------------------
-   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
+   !! NEMO/TOP 3.7 , NEMO Consortium (2014)
    !! $Id$
    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
-
 CONTAINS
 
@@ -63,5 +61,8 @@
       IF( nn_timing == 1 )   CALL timing_start('trc_ldf')
       !
+      
+!!gm  this call should be put in trcini !
       IF( kt == nittrc000 )   CALL ldf_ctl          ! initialisation & control of options
+!!gm end
 
       rldf = rldf_rat
@@ -73,31 +74,28 @@
 
       SELECT CASE ( nldf )                       ! compute lateral mixing trend and add it to the general trend
-      CASE ( 0 )   ;   CALL tra_ldf_lap   ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            )  ! iso-level laplacian
-      CASE ( 1 )                                                                                            ! rotated laplacian
-                       IF( ln_traldf_grif ) THEN
-                          CALL tra_ldf_iso_grif( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 )
-                       ELSE
-                          CALL tra_ldf_iso     ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 )
-                       ENDIF
-      CASE ( 2 )   ;   CALL tra_ldf_bilap ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            )  ! iso-level bilaplacian
-      CASE ( 3 )   ;   CALL tra_ldf_bilapg( kt, nittrc000, 'TRC',             trb, tra, jptra            )  ! s-coord. horizontal bilaplacian
+      CASE ( 0 )   ;   CALL tra_ldf_lap     ( kt, nittrc000, 'TRC', gtru, gtrv, ahtu, ahtv, trb     , tra, jpts, 1)  ! iso-level laplacian
+      CASE ( 1 )                                                                                      ! rotated laplacian
+         IF( ln_traldf_triad ) THEN                                                          
+                       CALL tra_ldf_iso_triad( kt, nittrc000,'TRC', gtru, gtrv, ahtu, ahtv, trb, trb, tra, jpts, 1 )    ! triad operator
+         ELSE                                                                                
+                       CALL tra_ldf_iso      ( kt, nittrc000,'TRC', gtru, gtrv, ahtu, ahtv, trb, trb, tra, jpts, 1 )    ! Madec operator
+         ENDIF
+      CASE ( 2 )   ;   CALL tra_ldf_blp      ( kt, nittrc000,'TRC', gtru, gtrv, ahtu, ahtv, trb     , tra, jpts )    ! iso-level bilaplacian
+      CASE ( 3 )   ;   CALL tra_ldf_iso_blp  ( kt, nittrc000,'TRC', gtru, gtrv, ahtu, ahtv, trb     , tra, jpts )    ! rotated   bilaplacian
          !
-      CASE ( -1 )                                     ! esopa: test all possibility with control print
-         CALL tra_ldf_lap   ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            )
-         WRITE(charout, FMT="('ldf0 ')") ;  CALL prt_ctl_trc_info(charout)
-                                            CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
-         IF( ln_traldf_grif ) THEN
-            CALL tra_ldf_iso_grif( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 )
-         ELSE
-            CALL tra_ldf_iso     ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 )
-         ENDIF
-         WRITE(charout, FMT="('ldf1 ')") ;  CALL prt_ctl_trc_info(charout)
-                                            CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
-         CALL tra_ldf_bilap ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            )
-         WRITE(charout, FMT="('ldf2 ')") ;  CALL prt_ctl_trc_info(charout)
-                                            CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
-         CALL tra_ldf_bilapg( kt, nittrc000, 'TRC',             trb, tra, jptra            )
-         WRITE(charout, FMT="('ldf3 ')") ;  CALL prt_ctl_trc_info(charout)
-                                            CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
+      CASE ( -1 )                                ! esopa: test all possibility with control print
+         CALL tra_ldf_lap   ( kt, nittrc000, 'TRC', gtru, gtrv, ahtu, ahtv, trb, tra, jpts, 1 ) 
+         WRITE(charout, FMT="('ldf0')")   ;   CALL prt_ctl_trc_info(charout)
+                                              CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd')
+         IF( ln_traldf_triad ) THEN
+            CALL tra_ldf_iso_triad( kt, nittrc000, 'TRC', gtru, gtrv, ahtu, ahtv, trb, trb, tra, jpts, 1 )
+         ELSE
+            CALL tra_ldf_iso      ( kt, nittrc000, 'TRC', gtru, gtrv, ahtu, ahtv, trb, trb, tra, jpts, 1 )  
+         ENDIF
+         WRITE(charout, FMT="('ldf1')")   ;   CALL prt_ctl_trc_info(charout)
+                                              CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
+         CALL tra_ldf_blp   ( kt, nittrc000, 'TRC', gtru, gtrv, ahtu, ahtv, trb, tra, jpts        ) 
+         WRITE(charout, FMT="('ldf2')")   ;   CALL prt_ctl_trc_info(charout)
+                                              CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
       END SELECT
       !
@@ -111,6 +109,6 @@
       !                                          ! print mean trends (used for debugging)
       IF( ln_ctl )   THEN
-         WRITE(charout, FMT="('ldf ')") ;  CALL prt_ctl_trc_info(charout)
-                                           CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
+         WRITE(charout, FMT="('ldf ')")   ;   CALL prt_ctl_trc_info(charout)
+                                              CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
       ENDIF
       !
@@ -119,4 +117,5 @@
    END SUBROUTINE trc_ldf
 
+!!gm ldf_ctl should be called in trcini  so that l_ldfslp=T  cause the slope init and calculation
 
    SUBROUTINE ldf_ctl
@@ -155,4 +154,8 @@
       IF( ioptio >  1 )   CALL ctl_stop( '          use ONE or NONE of the 2 lap/bilap operator type on tracer' )
       IF( ioptio == 0 )   nldf = -2   ! No lateral diffusion
+      
+      IF( ln_traldf_lap .AND. ln_trcldf_bilap )   CALL ctl_stop( ' bilaplacian should be used on both TRC and TRA' )
+      IF( ln_traldf_blp .AND. ln_trcldf_lap   )   CALL ctl_stop( '   laplacian should be used on both TRC and TRA' )
+      
       ioptio = 0
       IF( ln_trcldf_level )   ioptio = ioptio + 1
@@ -202,12 +205,9 @@
       IF( ierr == 1 )   CALL ctl_stop( ' iso-level in z-coordinate - partial step, not allowed' )
       IF( ierr == 2 )   CALL ctl_stop( ' isoneutral bilaplacian operator does not exist' )
-      IF( lk_traldf_eiv .AND. .NOT.ln_trcldf_iso )   &
+      IF( ln_ldfeiv .AND. .NOT.ln_trcldf_iso )   &
            CALL ctl_stop( '          eddy induced velocity on tracers',   &
            &              ' the eddy induced velocity on tracers requires isopycnal laplacian diffusion' )
       IF( nldf == 1 .OR. nldf == 3 ) THEN      ! rotation
-         IF( .NOT.lk_ldfslp )   CALL ctl_stop( '          the rotation of the diffusive tensor require key_ldfslp' )
-#if defined key_offline
-         l_traldf_rot = .TRUE.                 ! needed for trazdf_imp
-#endif
+         IF( .NOT.l_ldfslp )   CALL ctl_stop( '          the rotation of the diffusive tensor require l_ldfslp' )
       ENDIF
 
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/TRP/trcnxt.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/TRP/trcnxt.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/TRP/trcnxt.F90	(revision 4617)
@@ -107,8 +107,4 @@
       END DO
 
-
-#if defined key_obc
-!!      CALL obc_trc( kt )               ! OBC open boundaries
-#endif
 #if defined key_bdy
 !!      CALL bdy_trc( kt )               ! BDY open boundaries
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/oce_trc.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/oce_trc.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/oce_trc.F90	(revision 4617)
@@ -11,6 +11,6 @@
    !!   'key_top'                                                TOP models
    !!----------------------------------------------------------------------
-
-   !* Domain size *
+   !
+   !                                            !* Domain size *
    USE par_oce , ONLY :   jpi      =>   jpi        !: first  dimension of grid --> i 
    USE par_oce , ONLY :   jpj      =>   jpj        !: second dimension of grid --> j  
@@ -24,59 +24,37 @@
    USE par_oce , ONLY :   jp_sal   =>   jp_sal     !: indice for salinity
 
-   !* IO manager *
-   USE in_out_manager    
- 
-   !* Memory Allocation *
-   USE wrk_nemo      
- 
-   !* Timing *
-   USE timing    
- 
-   !* MPP library                         
-   USE lib_mpp 
+   USE in_out_manager                           !* IO manager *
+   USE wrk_nemo                                 !* Memory Allocation *
+   USE timing                                   !* Timing * 
+   USE lib_mpp                                  !* MPP library                         
+   USE lib_fortran                              !* Fortran utilities                         
+   USE lbclnk                                   !* Lateral boundary conditions                         
+   USE phycst                                   !* physical constants *
+   USE c1d                                      !* 1D configuration
+   USE dom_oce                                  !* model domain *
 
-   !* Fortran utilities                         
-   USE lib_fortran
-
-   !* Lateral boundary conditions                         
-   USE lbclnk
-
-   !* physical constants *
-   USE phycst            
-
-   !* 1D configuration
-   USE c1d                                         
-
-   !* model domain *
-   USE dom_oce 
-
-
-   !* ocean fields: here now and after fields *
-   USE oce , ONLY :   ua      =>    ua      !: i-horizontal velocity (m s-1) 
-   USE oce , ONLY :   va      =>    va      !: j-horizontal velocity (m s-1)
-   USE oce , ONLY :   un      =>    un      !: i-horizontal velocity (m s-1) 
-   USE oce , ONLY :   vn      =>    vn      !: j-horizontal velocity (m s-1)
-   USE oce , ONLY :   wn      =>    wn      !: vertical velocity (m s-1)  
-   USE oce , ONLY :   tsn     =>    tsn     !: 4D array contaning ( tn, sn )
-   USE oce , ONLY :   tsb     =>    tsb     !: 4D array contaning ( tb, sb )
-   USE oce , ONLY :   tsa     =>    tsa     !: 4D array contaning ( ta, sa )
-   USE oce , ONLY :   rhop    =>    rhop    !: potential volumic mass (kg m-3) 
-   USE oce , ONLY :   rhd     =>    rhd     !: in situ density anomalie rhd=(rho-rau0)/rau0 (no units)
-   USE oce , ONLY :   hdivn   =>    hdivn   !: horizontal divergence (1/s)
-   USE oce , ONLY :   rotn    =>    rotn    !: relative vorticity    [s-1]
-   USE oce , ONLY :   hdivb   =>    hdivb   !: horizontal divergence (1/s)
-   USE oce , ONLY :   rotb    =>    rotb    !: relative vorticity    [s-1]
-   USE oce , ONLY :   sshn    =>    sshn    !: sea surface height at t-point [m]   
-   USE oce , ONLY :   sshb    =>    sshb    !: sea surface height at t-point [m]   
-   USE oce , ONLY :   ssha    =>    ssha    !: sea surface height at t-point [m]   
-   USE oce , ONLY :   l_traldf_rot => l_traldf_rot  !: rotated laplacian operator for lateral diffusion
+   !                                            !* ocean fields: here now and after fields *
+   USE oce     , ONLY :   ua      =>    ua         !: i-horizontal velocity (m s-1) 
+   USE oce     , ONLY :   va      =>    va         !: j-horizontal velocity (m s-1)
+   USE oce     , ONLY :   un      =>    un         !: i-horizontal velocity (m s-1) 
+   USE oce     , ONLY :   vn      =>    vn         !: j-horizontal velocity (m s-1)
+   USE oce     , ONLY :   wn      =>    wn         !: vertical velocity (m s-1)  
+   USE oce     , ONLY :   tsn     =>    tsn        !: 4D array contaning ( tn, sn )
+   USE oce     , ONLY :   tsb     =>    tsb        !: 4D array contaning ( tb, sb )
+   USE oce     , ONLY :   tsa     =>    tsa        !: 4D array contaning ( ta, sa )
+   USE oce     , ONLY :   rhop    =>    rhop       !: potential volumic mass (kg m-3) 
+   USE oce     , ONLY :   rhd     =>    rhd        !: in situ density anomalie rhd=(rho-rau0)/rau0 (no units)
+   USE oce     , ONLY :   hdivn   =>    hdivn      !: now horizontal divergence (1/s)
+   USE oce     , ONLY :   sshn    =>    sshn       !: sea surface height at t-point [m]   
+   USE oce     , ONLY :   sshb    =>    sshb       !: sea surface height at t-point [m]   
+   USE oce     , ONLY :   ssha    =>    ssha       !: sea surface height at t-point [m]   
 #if defined key_offline
-   USE oce , ONLY :   gtsu    =>    gtsu    !: t-, s- and rd horizontal gradient at u- and
-   USE oce , ONLY :   gtsv    =>    gtsv    !:
-   USE oce , ONLY :   gru     =>    gru     !:
-   USE oce , ONLY :   grv     =>    grv     !: 
+   USE oce     , ONLY :   gtsu    =>    gtsu       !: t-, s- and rd horizontal gradient at u- and
+   USE oce     , ONLY :   gtsv    =>    gtsv       !:
+   USE oce     , ONLY :   gru     =>    gru        !:
+   USE oce     , ONLY :   grv     =>    grv        !: 
 #endif
 
-   !* surface fluxes *
+   !                                            !* surface fluxes *
    USE sbc_oce , ONLY :   utau       =>    utau       !: i-surface stress component
    USE sbc_oce , ONLY :   vtau       =>    vtau       !: j-surface stress component
@@ -97,22 +75,22 @@
    USE sbcrnf  , ONLY :   h_rnf      =>    h_rnf      !: river runoff   [Kg/m2/s]
 
-   USE trc_oce
+   USE trc_oce                                  !* ??? *
 
-   !* lateral diffusivity (tracers) *
-   USE ldftra_oce , ONLY :  rldf     =>   rldf        !: multiplicative coef. for lateral diffusivity
-   USE ldftra_oce , ONLY :  rn_aht_0 =>   rn_aht_0    !: horizontal eddy diffusivity for tracers (m2/s)
-   USE ldftra_oce , ONLY :  aht0     =>   aht0        !: horizontal eddy diffusivity for tracers (m2/s)
-   USE ldftra_oce , ONLY :  ahtb0    =>   ahtb0       !: background eddy diffusivity for isopycnal diff. (m2/s)
-   USE ldftra_oce , ONLY :  ahtu     =>   ahtu        !: lateral diffusivity coef. at u-points 
-   USE ldftra_oce , ONLY :  ahtv     =>   ahtv        !: lateral diffusivity coef. at v-points 
-   USE ldftra_oce , ONLY :  ahtw     =>   ahtw        !: lateral diffusivity coef. at w-points 
-   USE ldftra_oce , ONLY :  ahtt     =>   ahtt        !: lateral diffusivity coef. at t-points
-   USE ldftra_oce , ONLY :  aeiv0    =>   aeiv0       !: eddy induced velocity coefficient (m2/s) 
-   USE ldftra_oce , ONLY :  aeiu     =>   aeiu        !: eddy induced velocity coef. at u-points (m2/s)   
-   USE ldftra_oce , ONLY :  aeiv     =>   aeiv        !: eddy induced velocity coef. at v-points (m2/s) 
-   USE ldftra_oce , ONLY :  aeiw     =>   aeiw        !: eddy induced velocity coef. at w-points (m2/s) 
-   USE ldftra_oce , ONLY :  lk_traldf_eiv  =>  lk_traldf_eiv     !: eddy induced velocity flag
+   !                                            !* lateral diffusivity (tracers) *
+!   USE ldftra , ONLY :  rn_aht_0 =>   rn_aht_0     !: horizontal eddy diffusivity for tracers (m2/s)
+!   USE ldftra , ONLY :  aht0     =>   aht0         !: horizontal eddy diffusivity for tracers (m2/s)
+!   USE ldftra , ONLY :  ahtb0    =>   ahtb0        !: background eddy diffusivity for isopycnal diff. (m2/s)
+!   USE ldftra , ONLY :  aeiv0    =>   aeiv0        !: eddy induced velocity coefficient (m2/s) 
+   !
+   USE ldftra , ONLY :  ln_ldfeiv      =>  ln_ldfeiv        !: eddy induced velocity flag
+   USE ldftra , ONLY :  ln_traldf_lap  =>  ln_traldf_lap    !: laplacian   flag
+   USE ldftra , ONLY :  ln_traldf_blp  =>  ln_traldf_blp    !: bilaplacian flag
+   USE ldftra , ONLY :  rldf     =>   rldf         !: multiplicative coef. for lateral diffusivity
+   USE ldftra , ONLY :  ahtu     =>   ahtu         !: lateral diffusivity coef. at u-points 
+   USE ldftra , ONLY :  ahtv     =>   ahtv         !: lateral diffusivity coef. at v-points 
+   USE ldftra , ONLY :  aeiu     =>   aeiu         !: eddy induced velocity coef. at u-points   
+   USE ldftra , ONLY :  aeiv     =>   aeiv         !: eddy induced velocity coef. at v-points
 
-   !* vertical diffusion *
+   !                                            !* vertical diffusion *
    USE zdf_oce , ONLY :   avt        =>   avt         !: vert. diffusivity coef. at w-point for temp  
 # if defined key_zdfddm
@@ -120,5 +98,5 @@
 # endif
 
-   !* mixing & mixed layer depth *
+   !                                            !* mixing & mixed layer depth *
    USE zdfmxl , ONLY :   nmln        =>   nmln        !: number of level in the mixed layer
    USE zdfmxl , ONLY :   hmld        =>   hmld        !: mixing layer depth (turbocline)
@@ -126,12 +104,10 @@
    USE zdfmxl , ONLY :   hmlpt       =>   hmlpt       !: mixed layer depth at t-points (m)
 
-   !* direction of lateral diffusion *
-   USE ldfslp , ONLY :   lk_ldfslp  =>  lk_ldfslp     !: slopes flag
-# if   defined key_ldfslp
+   !                                            !* direction of lateral diffusion *
+   USE ldfslp , ONLY :   l_ldfslp   =>   l_ldfslp     !: slopes flag
    USE ldfslp , ONLY :   uslp       =>   uslp         !: i-direction slope at u-, w-points
    USE ldfslp , ONLY :   vslp       =>   vslp         !: j-direction slope at v-, w-points
    USE ldfslp , ONLY :   wslpi      =>   wslpi        !: i-direction slope at u-, w-points
    USE ldfslp , ONLY :   wslpj      =>   wslpj        !: j-direction slope at v-, w-points
-# endif
 
 #else
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/top_substitute.h90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/top_substitute.h90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/top_substitute.h90	(revision 4617)
@@ -6,11 +6,10 @@
    !! History :   1.0  !  2004-03 (C. Ethe) Original code
    !!             2.0  !  2007-12 (C. Ethe, G. Madec) new architecture
+   !!             3.7  !  2014-03 (G. Madec)  LDF simplification 
    !!----------------------------------------------------------------------
 #  include "domzgr_substitute.h90"
-#  include "ldfeiv_substitute.h90"
-#  include "ldftra_substitute.h90"
 #  include "vectopt_loop_substitute.h90"
    !!----------------------------------------------------------------------
-   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
+   !! NEMO/TOP 3.7 , NEMO Consortium (2014)
    !! $Id$ 
    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trc.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trc.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trc.F90	(revision 4617)
@@ -119,10 +119,8 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  avs_tm      !: vertical double diffusivity coeff. at w-point   [m/s]
 # endif
-#if defined key_ldfslp
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  wslpi_tm    !: i-direction slope at u-, w-points
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  wslpj_tm    !: j-direction slope at u-, w-points
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  uslp_tm     !: j-direction slope at u-, w-points
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  vslp_tm     !: j-direction slope at u-, w-points
-#endif
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  wslpi_tm, wslpi_temp    !: i-direction slope at w-points
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  wslpj_tm, wslpj_temp    !: j-direction slope at w-points
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  uslp_tm , uslp_temp     !: i-direction slope at u-points
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  vslp_tm , vslp_temp     !: j-direction slope at v-points
 #if defined key_trabbl
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::  ahu_bbl_tm  !: u-, w-points
@@ -148,25 +146,6 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::  wndm_tm    !: 10m wind average [m]
    !
-#if defined key_traldf_c3d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  ahtt_tm, ahtu_tm, ahtv_tm, ahtw_tm   !: ** 3D coefficients ** at T-,U-,V-,W-points
-#elif defined key_traldf_c2d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::  ahtt_tm, ahtu_tm, ahtv_tm, ahtw_tm   !: ** 2D coefficients ** at T-,U-,V-,W-points
-#elif defined key_traldf_c1d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:)       ::  ahtt_tm, ahtu_tm, ahtv_tm, ahtw_tm   !: ** 1D coefficients ** at T-,U-,V-,W-points
-#else
-   REAL(wp), PUBLIC                                        ::  ahtt_tm, ahtu_tm, ahtv_tm, ahtw_tm   !: ** 0D coefficients ** at T-,U-,V-,W-points
-#endif
-   !
-#if defined key_traldf_eiv
-#  if defined key_traldf_c3d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  aeiu_tm , aeiv_tm , aeiw_tm   !: ** 3D coefficients **
-#  elif defined key_traldf_c2d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::  aeiu_tm , aeiv_tm , aeiw_tm   !: ** 2D coefficients **
-#  elif defined key_traldf_c1d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:)       ::  aeiu_tm , aeiv_tm, aeiw_tm   !: ** 1D coefficients **
-#  else
-   REAL(wp), PUBLIC                                        ::  aeiu_tm , aeiv_tm , aeiw_tm   !: ** 0D coefficients **
-#  endif
-#endif
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  ahtu_tm, ahtv_tm   !: ** 3D coefficients ** at T-,U-,V-,W-points
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  aeiu_tm, aeiv_tm   !: ** 3D coefficients **
 
    ! Temporary physical arrays for sub_stepping
@@ -188,34 +167,12 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::  ahu_bbl_temp, ahv_bbl_temp, utr_bbl_temp, vtr_bbl_temp !: hold current values 
 #endif
-   !
-#if defined key_ldfslp
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  wslpi_temp, wslpj_temp, uslp_temp, vslp_temp    !: hold current values 
-#endif
    ! 
 # if defined key_zdfddm
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  avs_temp      !: salinity vertical diffusivity coeff. at w-point   [m/s]
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   avs_temp      !: salinity vertical diffusivity coeff. at w-point   [m/s]
 # endif
    !
-#if defined key_traldf_c3d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  ahtt_temp, ahtu_temp, ahtv_temp, ahtw_temp   
-#elif defined key_traldf_c2d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::  ahtt_temp, ahtu_temp, ahtv_temp, ahtw_temp  
-#elif defined key_traldf_c1d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:)       ::  ahtt_temp, ahtu_temp, ahtv_temp, ahtw_temp 
-#else
-   REAL(wp), PUBLIC                                        ::  ahtt_temp, ahtu_temp, ahtv_temp, ahtw_temp
-#endif
-   !
-#if defined key_traldf_eiv
-# if defined key_traldf_c3d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  aeiu_temp , aeiv_temp , aeiw_temp   !: ** 3D coefficients **
-# elif defined key_traldf_c2d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::  aeiu_temp , aeiv_temp , aeiw_temp   !: ** 2D coefficients **
-# elif defined key_traldf_c1d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:)       ::  aeiu_temp , aeiv_temp, aeiw_temp   !: ** 1D coefficients **
-# else
-   REAL(wp), PUBLIC                                        ::  aeiu_temp , aeiv_temp , aeiw_temp   !: ** 0D coefficients **
-# endif
-# endif
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   ahtu_temp, ahtv_temp   
+   !
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   aeiu_temp , aeiv_temp , aeiw_temp   !: ** 3D coefficients **
 
    !!----------------------------------------------------------------------
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcsms.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcsms.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcsms.F90	(revision 4617)
@@ -52,5 +52,5 @@
       IF( lk_c14b    )   CALL trc_sms_c14b   ( kt )    ! surface fluxes of C14
       IF( lk_my_trc  )   CALL trc_sms_my_trc ( kt )    ! MY_TRC  tracers
-
+      !
       IF(ln_ctl) THEN      ! print mean trends (used for debugging)
          WRITE(charout, FMT="('sms ')")
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcstp.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcstp.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcstp.F90	(revision 4617)
@@ -10,9 +10,10 @@
    !!   trc_stp      : passive tracer system time-stepping
    !!----------------------------------------------------------------------
+   USE trc
    USE oce_trc          ! ocean dynamics and active tracers variables
    USE sbc_oce
-   USE trc
    USE trctrp           ! passive tracers transport
    USE trcsms           ! passive tracers sources and sinks
+   USE trcsub
    USE prtctl_trc       ! Print control for debbuging
    USE trcdia
@@ -21,7 +22,7 @@
    USE trdmod_trc_oce
    USE trdmld_trc
+   !
    USE iom
    USE in_out_manager
-   USE trcsub
 
    IMPLICIT NONE
@@ -50,8 +51,8 @@
       !!-------------------------------------------------------------------
       INTEGER, INTENT( in ) ::  kt      ! ocean time-step index
+      !
       INTEGER               ::  jk, jn  ! dummy loop indices
       REAL(wp)              ::  ztrai
       CHARACTER (len=25)    ::  charout 
-
       REAL(wp), DIMENSION(:,:), POINTER ::   zqsr_tmp ! save qsr during TOP time-step
       !!-------------------------------------------------------------------
@@ -79,4 +80,5 @@
       IF( nn_dttrc /= 1 )   CALL trc_sub_stp( kt )  ! averaging physical variables for sub-stepping
       !    
+
       IF( MOD( kt , nn_dttrc ) == 0 ) THEN      ! only every nn_dttrc time step
          !
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcsub.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcsub.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcsub.F90	(revision 4617)
@@ -25,10 +25,7 @@
    USE zdf_oce
    USE domvvl
-   USE divcur          ! hor. divergence and curl      (div & cur routines)
+   USE divhor          ! hor. divergence
    USE sbcrnf, ONLY: h_rnf, nk_rnf   ! River runoff 
    USE bdy_oce
-#if defined key_obc
-   USE obc_oce, ONLY: obctmsk
-#endif
 #if defined key_agrif
    USE agrif_opa_update
@@ -43,5 +40,4 @@
    PUBLIC   trc_sub_ssh      ! called by trc_stp to reset physics variables
 
-   !!* Module variables
    REAL(wp)  :: r1_ndttrc     !    1 /  nn_dttrc 
    REAL(wp)  :: r1_ndttrcp1   !    1 / (nn_dttrc+1) 
@@ -50,5 +46,5 @@
 #  include "top_substitute.h90"
    !!----------------------------------------------------------------------
-   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
+   !! NEMO/TOP 3.7 , NEMO Consortium (2014)
    !! $Id: trcstp.F90 2528 2010-12-27 17:33:53Z rblod $ 
    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
@@ -94,13 +90,13 @@
           tsn_tm  (:,:,:,jp_sal) = tsn_tm  (:,:,:,jp_sal) + tsn  (:,:,:,jp_sal) * fse3t(:,:,:)  
           rhop_tm (:,:,:)        = rhop_tm (:,:,:)        + rhop (:,:,:)        * fse3t(:,:,:)  
-          avt_tm  (:,:,:)        = avt_tm  (:,:,:)        + avt  (:,:,:)        * fse3w(:,:,:)  
+          avt_tm  (:,:,:)        = avt_tm  (:,:,:)        + avt  (:,:,:) 
 # if defined key_zdfddm
-          avs_tm  (:,:,:)        = avs_tm  (:,:,:)        + avs  (:,:,:)        * fse3w(:,:,:)  
+          avs_tm  (:,:,:)        = avs_tm  (:,:,:)        + avs  (:,:,:) 
 # endif
 #if defined key_ldfslp
-          wslpi_tm(:,:,:)        = wslpi_tm(:,:,:)        + wslpi(:,:,:)        * fse3w(:,:,:) 
-          wslpj_tm(:,:,:)        = wslpj_tm(:,:,:)        + wslpj(:,:,:)        * fse3w(:,:,:) 
-          uslp_tm (:,:,:)        = uslp_tm (:,:,:)        + uslp (:,:,:)        * fse3u(:,:,:) 
-          vslp_tm (:,:,:)        = vslp_tm (:,:,:)        + vslp (:,:,:)        * fse3v(:,:,:) 
+          wslpi_tm(:,:,:)        = wslpi_tm(:,:,:)        + wslpi(:,:,:) 
+          wslpj_tm(:,:,:)        = wslpj_tm(:,:,:)        + wslpj(:,:,:)
+          uslp_tm (:,:,:)        = uslp_tm (:,:,:)        + uslp (:,:,:) 
+          vslp_tm (:,:,:)        = vslp_tm (:,:,:)        + vslp (:,:,:)
 #endif
 # if defined key_trabbl
@@ -127,44 +123,10 @@
           wndm_tm  (:,:)         = wndm_tm  (:,:)         + wndm  (:,:)
           !
-#if defined key_traldf_c3d
-          ahtt_tm  (:,:,:)       = ahtt_tm  (:,:,:)       + ahtt(:,:,:)         * fse3t(:,:,:)
-          ahtu_tm  (:,:,:)       = ahtu_tm  (:,:,:)       + ahtu(:,:,:)         * fse3u(:,:,:)
-          ahtv_tm  (:,:,:)       = ahtv_tm  (:,:,:)       + ahtv(:,:,:)         * fse3v(:,:,:)
-          ahtw_tm  (:,:,:)       = ahtw_tm  (:,:,:)       + ahtw(:,:,:)         * fse3w(:,:,:)
-#elif defined key_traldf_c2d
-          ahtt_tm  (:,:)         = ahtt_tm  (:,:)         + ahtt(:,:)
-          ahtu_tm  (:,:)         = ahtu_tm  (:,:)         + ahtu(:,:)
-          ahtv_tm  (:,:)         = ahtv_tm  (:,:)         + ahtv(:,:)
-          ahtw_tm  (:,:)         = ahtw_tm  (:,:)         + ahtw(:,:)
-#elif defined key_traldf_c1d
-          ahtt_tm  (:)           = ahtt_tm  (:,:)         + ahtt(:)
-          ahtu_tm  (:)           = ahtu_tm  (:,:)         + ahtu(:)
-          ahtv_tm  (:)           = ahtv_tm  (:,:)         + ahtv(:)
-          ahtw_tm  (:)           = ahtw_tm  (:,:)         + ahtw(:)
-#else
-          ahtt_tm                = ahtt_tm                + ahtt
-          ahtu_tm                = ahtu_tm                + ahtu
-          ahtv_tm                = ahtv_tm                + ahtv
-          ahtw_tm                = ahtw_tm                + ahtw
-#endif
-#if defined key_traldf_eiv
-#  if defined key_traldf_c3d
-          aeiu_tm  (:,:,:)       = aeiu_tm  (:,:,:)       + aeiu(:,:,:)         * fse3u(:,:,:)
-          aeiv_tm  (:,:,:)       = aeiv_tm  (:,:,:)       + aeiv(:,:,:)         * fse3v(:,:,:)
-          aeiw_tm  (:,:,:)       = aeiw_tm  (:,:,:)       + aeiw(:,:,:)         * fse3w(:,:,:)
-#  elif defined key_traldf_c2d
-          aeiu_tm  (:,:)         = aeiu_tm  (:,:)         + aeiu(:,:)
-          aeiv_tm  (:,:)         = aeiv_tm  (:,:)         + aeiv(:,:)
-          aeiw_tm  (:,:)         = aeiw_tm  (:,:)         + aeiw(:,:)
-#  elif defined key_traldf_c1d
-          aeiu_tm  (:)           = aeiu_tm  (:,:)         + aeiu(:)
-          aeiv_tm  (:)           = aeiv_tm  (:,:)         + aeiv(:)
-          aeiw_tm  (:)           = aeiw_tm  (:,:)         + aeiw(:)
-#  else
-          aeiu_tm                = aeiu_tm                + aeiu
-          aeiv_tm                = aeiv_tm                + aeiv
-          aeiw_tm                = aeiw_tm                + aeiw
-#  endif
-#endif
+          ahtu_tm  (:,:,:)       = ahtu_tm  (:,:,:)       + ahtu(:,:,:)
+          ahtv_tm  (:,:,:)       = ahtv_tm  (:,:,:)       + ahtv(:,:,:)
+          IF( ln_ldfeiv ) THEN
+            aeiu_tm  (:,:,:)     = aeiu_tm  (:,:,:)       + aeiu(:,:,:)
+            aeiv_tm  (:,:,:)     = aeiv_tm  (:,:,:)       + aeiv(:,:,:)
+          ENDIF
 
       ELSE                           !  It is time to substep 
@@ -179,10 +141,10 @@
          avs_temp   (:,:,:)      = avs   (:,:,:)
 # endif
-#if defined key_ldfslp
-         wslpi_temp (:,:,:)      = wslpi (:,:,:)
-         wslpj_temp (:,:,:)      = wslpj (:,:,:)
-         uslp_temp  (:,:,:)      = uslp  (:,:,:)
-         vslp_temp  (:,:,:)      = vslp  (:,:,:)
-#endif
+         IF( l_ldfslp ) THEN
+            uslp_temp  (:,:,:)      = uslp  (:,:,:)
+            vslp_temp  (:,:,:)      = vslp  (:,:,:)
+            wslpi_temp (:,:,:)      = wslpi (:,:,:)
+            wslpj_temp (:,:,:)      = wslpj (:,:,:)
+         ENDIF
 # if defined key_trabbl
           IF( nn_bbl_ldf == 1 ) THEN
@@ -214,50 +176,13 @@
          qsr_temp   (:,:)        = qsr   (:,:)
          wndm_temp  (:,:)        = wndm  (:,:)
-#if defined key_traldf_c3d
+         !
          ahtu_temp  (:,:,:)      = ahtu  (:,:,:)
          ahtv_temp  (:,:,:)      = ahtv  (:,:,:)
-         ahtw_temp  (:,:,:)      = ahtw  (:,:,:)
-         ahtt_temp  (:,:,:)      = ahtt  (:,:,:)
-#elif defined key_traldf_c2d
-         ahtu_temp  (:,:)        = ahtu  (:,:)
-         ahtv_temp  (:,:)        = ahtv  (:,:)
-         ahtw_temp  (:,:)        = ahtw  (:,:)
-         ahtt_temp  (:,:)        = ahtt  (:,:)
-#elif defined key_traldf_c1d
-         ahtu_temp  (:)          = ahtu  (:)
-         ahtv_temp  (:)          = ahtv  (:)
-         ahtw_temp  (:)          = ahtw  (:)
-         ahtt_temp  (:)          = ahtt  (:)
-#else
-         ahtu_temp               = ahtu
-         ahtv_temp               = ahtv
-         ahtw_temp               = ahtw
-         ahtt_temp               = ahtt
-#endif
-
-#if defined key_traldf_eiv
-# if defined key_traldf_c3d
-         aeiu_temp  (:,:,:)      = aeiu(:,:,:)
-         aeiv_temp  (:,:,:)      = aeiv(:,:,:)
-         aeiw_temp  (:,:,:)      = aeiw(:,:,:)
-# elif defined key_traldf_c2d
-         aeiu_temp  (:,:)        = aeiu(:,:)
-         aeiv_temp  (:,:)        = aeiv(:,:)
-         aeiw_temp  (:,:)        = aeiw(:,:)
-# elif defined key_traldf_c1d
-         aeiu_temp  (:)          = aeiu(:)
-         aeiv_temp  (:)          = aeiv(:)
-         aeiw_temp  (:)          = aeiw(:)
-# else
-         aeiu_temp               = aeiu
-         aeiv_temp               = aeiv
-         aeiw_temp               = aeiw
-# endif
-#endif
+         IF( ln_ldfeiv ) THEN
+            aeiu_temp  (:,:,:)   = aeiu(:,:,:)
+            aeiv_temp  (:,:,:)   = aeiv(:,:,:)
+         ENDIF
          !                                    !  Variables reset in trc_sub_ssh
-         rotn_temp  (:,:,:)      = rotn  (:,:,:)
          hdivn_temp (:,:,:)      = hdivn (:,:,:)
-         rotb_temp  (:,:,:)      = rotb  (:,:,:)
-         hdivb_temp (:,:,:)      = hdivb (:,:,:)
          hu_temp    (:,:)        = hu    (:,:)
          hv_temp    (:,:)        = hv    (:,:)
@@ -287,14 +212,14 @@
          tsn_tm   (:,:,:,jp_sal) = tsn_tm  (:,:,:,jp_sal) + tsn  (:,:,:,jp_sal) * e3t_temp(:,:,:)  
          rhop_tm (:,:,:)         = rhop_tm (:,:,:)        + rhop (:,:,:)        * e3t_temp(:,:,:)  
-         avt_tm   (:,:,:)        = avt_tm  (:,:,:)        + avt  (:,:,:)        * e3w_temp(:,:,:)  
+         avt_tm   (:,:,:)        = avt_tm  (:,:,:)        + avt  (:,:,:)  
 # if defined key_zdfddm
-         avs_tm   (:,:,:)        = avs_tm  (:,:,:)        + avs  (:,:,:)        * e3w_temp(:,:,:)  
-# endif
-#if defined key_ldfslp
-         wslpi_tm (:,:,:)        = wslpi_tm(:,:,:)        + wslpi(:,:,:)        * e3w_temp(:,:,:) 
-         wslpj_tm (:,:,:)        = wslpj_tm(:,:,:)        + wslpj(:,:,:)        * e3w_temp(:,:,:) 
-         uslp_tm  (:,:,:)        = uslp_tm (:,:,:)        + uslp (:,:,:)        * e3u_temp(:,:,:) 
-         vslp_tm  (:,:,:)        = vslp_tm (:,:,:)        + vslp (:,:,:)        * e3v_temp(:,:,:) 
-#endif
+         avs_tm   (:,:,:)        = avs_tm  (:,:,:)        + avs  (:,:,:)  
+# endif
+         IF( l_ldfslp ) THEN
+            wslpi_tm (:,:,:)        = wslpi_tm(:,:,:)        + wslpi(:,:,:)
+            wslpj_tm (:,:,:)        = wslpj_tm(:,:,:)        + wslpj(:,:,:)
+            uslp_tm  (:,:,:)        = uslp_tm (:,:,:)        + uslp (:,:,:)
+            vslp_tm  (:,:,:)        = vslp_tm (:,:,:)        + vslp (:,:,:) 
+         ENDIF
 # if defined key_trabbl
           IF( nn_bbl_ldf == 1 ) THEN
@@ -377,96 +302,32 @@
                   tsn  (ji,jj,jk,jp_sal) = tsn_tm  (ji,jj,jk,jp_sal) * z1_ne3t
                   rhop (ji,jj,jk)        = rhop_tm (ji,jj,jk)        * z1_ne3t
-                  avt  (ji,jj,jk)        = avt_tm  (ji,jj,jk)        * z1_ne3w
+                  avt  (ji,jj,jk)        = avt_tm  (ji,jj,jk)        * r1_ndttrcp1
 # if defined key_zdfddm
-                  avs  (ji,jj,jk)        = avs_tm  (ji,jj,jk)        * z1_ne3w
-# endif
-#if defined key_ldfslp
-                  wslpi(ji,jj,jk)        = wslpi_tm(ji,jj,jk)        * z1_ne3w 
-                  wslpj(ji,jj,jk)        = wslpj_tm(ji,jj,jk)        * z1_ne3w 
-                  uslp (ji,jj,jk)        = uslp_tm (ji,jj,jk)        * z1_ne3u 
-                  vslp (ji,jj,jk)        = vslp_tm (ji,jj,jk)        * z1_ne3v 
-#endif
-               ENDDO
-            ENDDO
-         ENDDO
-
-#if defined key_traldf_c3d
-          ahtt_tm  (:,:,:)       = ahtt_tm  (:,:,:)       + ahtt(:,:,:)         * e3t_temp(:,:,:)
-          ahtu_tm  (:,:,:)       = ahtu_tm  (:,:,:)       + ahtu(:,:,:)         * e3u_temp(:,:,:)
-          ahtv_tm  (:,:,:)       = ahtv_tm  (:,:,:)       + ahtv(:,:,:)         * e3v_temp(:,:,:)
-          ahtw_tm  (:,:,:)       = ahtw_tm  (:,:,:)       + ahtw(:,:,:)         * e3w_temp(:,:,:)
-          !
-          ahtt     (:,:,:)       = ahtt_tm  (:,:,:) * r1_ndttrcp1  / fse3t(:,:,:)
-          ahtu     (:,:,:)       = ahtu_tm  (:,:,:) * r1_ndttrcp1  / fse3u(:,:,:)
-          ahtv     (:,:,:)       = ahtv_tm  (:,:,:) * r1_ndttrcp1  / fse3v(:,:,:)
-          ahtw     (:,:,:)       = ahtw_tm  (:,:,:) * r1_ndttrcp1  / fse3w(:,:,:)
-#elif defined key_traldf_c2d
-          ahtt_tm  (:,:)         = ahtt_tm  (:,:)         + ahtt(:,:)
-          ahtu_tm  (:,:)         = ahtu_tm  (:,:)         + ahtu(:,:)
-          ahtv_tm  (:,:)         = ahtv_tm  (:,:)         + ahtv(:,:)
-          ahtw_tm  (:,:)         = ahtw_tm  (:,:)         + ahtw(:,:)
-          !
-          ahtt     (:,:)         = ahtt_tm  (:,:)   * r1_ndttrcp1
-          ahtu     (:,:)         = ahtu_tm  (:,:)   * r1_ndttrcp1
-          ahtv     (:,:)         = ahtv_tm  (:,:)   * r1_ndttrcp1
-          ahtw     (:,:)         = ahtw_tm  (:,:)   * r1_ndttrcp1
-#elif defined key_traldf_c1d
-          ahtt_tm  (:)           = ahtt_tm  (:,:)         + ahtt(:)
-          ahtu_tm  (:)           = ahtu_tm  (:,:)         + ahtu(:)
-          ahtv_tm  (:)           = ahtv_tm  (:,:)         + ahtv(:)
-          ahtw_tm  (:)           = ahtw_tm  (:,:)         + ahtw(:)
-          !
-          ahtt     (:)           = ahtt_tm  (:)     * r1_ndttrcp1
-          ahtu     (:)           = ahtu_tm  (:)     * r1_ndttrcp1
-          ahtv     (:)           = ahtv_tm  (:)     * r1_ndttrcp1
-          ahtw     (:)           = ahtw_tm  (:)     * r1_ndttrcp1
-#else
-          ahtt_tm                = ahtt_tm                + ahtt
-          ahtu_tm                = ahtu_tm                + ahtu
-          ahtv_tm                = ahtv_tm                + ahtv
-          ahtw_tm                = ahtw_tm                + ahtw
-          !
-          ahtt                   = ahtt_tm          * r1_ndttrcp1
-          ahtu                   = ahtu_tm          * r1_ndttrcp1
-          ahtv                   = ahtv_tm          * r1_ndttrcp1
-          ahtw                   = ahtw_tm          * r1_ndttrcp1
-#endif
-
-#if defined key_traldf_eiv
-# if defined key_traldf_c3d
-          aeiu_tm  (:,:,:)       = aeiu_tm  (:,:,:)       + aeiu(:,:,:)         * e3u_temp(:,:,:)
-          aeiv_tm  (:,:,:)       = aeiv_tm  (:,:,:)       + aeiv(:,:,:)         * e3v_temp(:,:,:)
-          aeiw_tm  (:,:,:)       = aeiw_tm  (:,:,:)       + aeiw(:,:,:)         * e3w_temp(:,:,:)
-          !
-          aeiu     (:,:,:)       = aeiu_tm  (:,:,:) * r1_ndttrcp1  / fse3u(:,:,:)
-          aeiv     (:,:,:)       = aeiv_tm  (:,:,:) * r1_ndttrcp1  / fse3v(:,:,:)
-          aeiw     (:,:,:)       = aeiw_tm  (:,:,:) * r1_ndttrcp1  / fse3w(:,:,:)
-# elif defined key_traldf_c2d
-          aeiu_tm  (:,:)         = aeiu_tm  (:,:)         + aeiu(:,:)
-          aeiv_tm  (:,:)         = aeiv_tm  (:,:)         + aeiv(:,:)
-          aeiw_tm  (:,:)         = aeiw_tm  (:,:)         + aeiw(:,:)
-          !
-          aeiu     (:,:)         = aeiu_tm  (:,:)   * r1_ndttrcp1
-          aeiv     (:,:)         = aeiv_tm  (:,:)   * r1_ndttrcp1
-          aeiw     (:,:)         = aeiw_tm  (:,:)   * r1_ndttrcp1
-# elif defined key_traldf_c1d
-          aeiu_tm  (:)           = aeiu_tm  (:,:)         + aeiu(:)
-          aeiv_tm  (:)           = aeiv_tm  (:,:)         + aeiv(:)
-          aeiw_tm  (:)           = aeiw_tm  (:,:)         + aeiw(:)
-          !
-          aeiu     (:)           = aeiu_tm  (:)     * r1_ndttrcp1
-          aeiv     (:)           = aeiv_tm  (:)     * r1_ndttrcp1
-          aeiw     (:)           = aeiw_tm  (:)     * r1_ndttrcp1
-# else
-          aeiu_tm                = aeiu_tm                + aeiu
-          aeiv_tm                = aeiv_tm                + aeiv
-          aeiw_tm                = aeiw_tm                + aeiw
-          !
-          aeiu                   = aeiu_tm          * r1_ndttrcp1
-          aeiv                   = aeiv_tm          * r1_ndttrcp1
-          aeiw                   = aeiw_tm          * r1_ndttrcp1
-# endif
-#endif
-
+                  avs  (ji,jj,jk)        = avs_tm  (ji,jj,jk)        * r1_ndttrcp1
+# endif
+               END DO
+            END DO
+         END DO
+         IF( l_ldfslp ) THEN
+            wslpi(:,:,:)        = wslpi_tm(:,:,:)        * r1_ndttrcp1 
+            wslpj(:,:,:)        = wslpj_tm(:,:,:)        * r1_ndttrcp1 
+            uslp (:,:,:)        = uslp_tm (:,:,:)        * r1_ndttrcp1 
+            vslp (:,:,:)        = vslp_tm (:,:,:)        * r1_ndttrcp1 
+         ENDIF
+         !
+         ahtu_tm(:,:,:) = ahtu_tm(:,:,:) + ahtu(:,:,:)
+         ahtv_tm(:,:,:) = ahtv_tm(:,:,:) + ahtv(:,:,:)
+         !
+         ahtu   (:,:,:) = ahtu_tm(:,:,:) * r1_ndttrcp1
+         ahtv   (:,:,:) = ahtv_tm(:,:,:) * r1_ndttrcp1
+         IF( ln_ldfeiv ) THEN
+            aeiu_tm(:,:,:) = aeiu_tm(:,:,:) + aeiu(:,:,:)
+            aeiv_tm(:,:,:) = aeiv_tm(:,:,:) + aeiv(:,:,:)
+            !
+            aeiu   (:,:,:) = aeiu_tm(:,:,:) * r1_ndttrcp1
+            aeiv   (:,:,:) = aeiv_tm(:,:,:) * r1_ndttrcp1
+         ENDIF
+
+!!gm  All those lbc_lnk are useless since the averaging is done over full arrays !
          CALL lbc_lnk( un    (:,:,:)       , 'U',-1. ) 
          CALL lbc_lnk( vn    (:,:,:)       , 'V',-1. ) 
@@ -478,10 +339,10 @@
           CALL lbc_lnk( avs  (:,:,:)       , 'W', 1. ) 
 # endif
-#if defined key_ldfslp
-         CALL lbc_lnk( uslp  (:,:,:)       , 'U',-1. ) 
-         CALL lbc_lnk( vslp  (:,:,:)       , 'V',-1. ) 
-         CALL lbc_lnk( wslpi (:,:,:)       , 'W',-1. ) 
-         CALL lbc_lnk( wslpj (:,:,:)       , 'W',-1. ) 
-#endif
+         IF( l_ldfslp ) THEN
+            CALL lbc_lnk( uslp  (:,:,:)       , 'U',-1. ) 
+            CALL lbc_lnk( vslp  (:,:,:)       , 'V',-1. ) 
+            CALL lbc_lnk( wslpi (:,:,:)       , 'W',-1. ) 
+            CALL lbc_lnk( wslpj (:,:,:)       , 'W',-1. ) 
+         ENDIF
          CALL lbc_lnk( sshn  (:,:)         , 'T', 1. ) 
 !!Z~     CALL lbc_lnk( sshu_n(:,:)         , 'U', 1. ) 
@@ -513,34 +374,15 @@
          ENDIF
 # endif
-#if defined key_traldf_c3d
-         CALL lbc_lnk( ahtt  (:,:,:)       , 'T', 1. ) 
          CALL lbc_lnk( ahtu  (:,:,:)       , 'U', 1. ) 
          CALL lbc_lnk( ahtv  (:,:,:)       , 'V', 1. ) 
-         CALL lbc_lnk( ahtw  (:,:,:)       , 'W', 1. ) 
-#elif defined key_traldf_c2d
-         CALL lbc_lnk( ahtt  (:,:)         , 'T', 1. ) 
-         CALL lbc_lnk( ahtu  (:,:)         , 'U', 1. ) 
-         CALL lbc_lnk( ahtv  (:,:)         , 'V', 1. ) 
-         CALL lbc_lnk( ahtw  (:,:)         , 'W', 1. ) 
-#endif
-#if defined key_traldf_eiv
-#if defined key_traldf_c3d
-         CALL lbc_lnk( aeiu  (:,:,:)       , 'U', 1. ) 
-         CALL lbc_lnk( aeiv  (:,:,:)       , 'V', 1. ) 
-         CALL lbc_lnk( aeiw  (:,:,:)       , 'W', 1. ) 
-#elif defined key_traldf_c2d
-         CALL lbc_lnk( aeiu  (:,:)         , 'U', 1. ) 
-         CALL lbc_lnk( aeiv  (:,:)         , 'V', 1. ) 
-         CALL lbc_lnk( aeiw  (:,:)         , 'W', 1. ) 
-#endif
-#endif
+         IF( ln_ldfeiv ) THEN
+            CALL lbc_lnk( aeiu  (:,:,:)    , 'U', 1. ) 
+            CALL lbc_lnk( aeiv  (:,:,:)    , 'V', 1. ) 
+         ENDIF
          !
          CALL trc_sub_ssh( kt )         ! after ssh & vertical velocity
          !
          CALL lbc_lnk( wn    (:,:,:)       , 'W',-1. ) 
-         CALL lbc_lnk( rotn  (:,:,:)       , 'F', 1. ) 
          CALL lbc_lnk( hdivn (:,:,:)       , 'T', 1. ) 
-         CALL lbc_lnk( rotb  (:,:,:)       , 'F', 1. ) 
-         CALL lbc_lnk( hdivb (:,:,:)       , 'T', 1. ) 
          CALL lbc_lnk( hu    (:,:)         , 'U', 1. ) 
          CALL lbc_lnk( hv    (:,:)         , 'V', 1. ) 
@@ -553,12 +395,12 @@
    END SUBROUTINE trc_sub_stp
 
+
    SUBROUTINE trc_sub_ini
       !!-------------------------------------------------------------------
       !!                     ***  ROUTINE trc_sub_ini  ***
       !!                      
-      !! ** Purpose : Initialize variables needed for sub-stepping passive tracers
+      !! ** Purpose :   Initialize variables needed for sub-stepping passive tracers
       !! 
-      !! ** Method  : 
-      !!              Compute the averages for sub-stepping
+      !! ** Method  :   Compute the averages for sub-stepping
       !!-------------------------------------------------------------------
       INTEGER ::   ierr
@@ -580,14 +422,14 @@
       tsn_tm  (:,:,:,jp_sal) = tsn  (:,:,:,jp_sal) * fse3t(:,:,:)  
       rhop_tm (:,:,:)        = rhop (:,:,:)        * fse3t(:,:,:)  
-      avt_tm  (:,:,:)        = avt  (:,:,:)        * fse3w(:,:,:)  
+      avt_tm  (:,:,:)        = avt  (:,:,:) 
 # if defined key_zdfddm
-      avs_tm  (:,:,:)        = avs  (:,:,:)        * fse3w(:,:,:)  
-# endif
-#if defined key_ldfslp
-      wslpi_tm(:,:,:)        = wslpi(:,:,:)        * fse3w(:,:,:) 
-      wslpj_tm(:,:,:)        = wslpj(:,:,:)        * fse3w(:,:,:) 
-      uslp_tm (:,:,:)        = uslp (:,:,:)        * fse3u(:,:,:) 
-      vslp_tm (:,:,:)        = vslp (:,:,:)        * fse3v(:,:,:) 
-#endif
+      avs_tm  (:,:,:)        = avs  (:,:,:)  
+# endif
+      IF( l_ldfslp ) THEN
+         wslpi_tm(:,:,:)        = wslpi(:,:,:)
+         wslpj_tm(:,:,:)        = wslpj(:,:,:)
+         uslp_tm (:,:,:)        = uslp (:,:,:)
+         vslp_tm (:,:,:)        = vslp (:,:,:)
+      ENDIF
       sshn_tm  (:,:) = sshn  (:,:) 
 !!Z~  sshu_n_tm(:,:) = sshu_n(:,:) 
@@ -614,49 +456,16 @@
 # endif
       !
-#if defined key_traldf_c3d
-      ahtt_tm(:,:,:) = ahtt(:,:,:) * fse3t(:,:,:)
-      ahtu_tm(:,:,:) = ahtu(:,:,:) * fse3u(:,:,:)
-      ahtv_tm(:,:,:) = ahtv(:,:,:) * fse3v(:,:,:)
-      ahtw_tm(:,:,:) = ahtw(:,:,:) * fse3w(:,:,:)
-#elif defined key_traldf_c2d
-      ahtt_tm(:,:)   = ahtt(:,:)
-      ahtu_tm(:,:)   = ahtu(:,:)
-      ahtv_tm(:,:)   = ahtv(:,:)
-      ahtw_tm(:,:)   = ahtw(:,:)
-#elif defined key_traldf_c1d
-      ahtt_tm(:)     = ahtt(:)
-      ahtu_tm(:)     = ahtu(:)
-      ahtv_tm(:)     = ahtv(:)
-      ahtw_tm(:)     = ahtw(:)
-#else
-      ahtt_tm        = ahtt
-      ahtu_tm        = ahtu
-      ahtv_tm        = ahtv
-      ahtw_tm        = ahtw
-#endif
-      !
-#if defined key_traldf_eiv
-#  if defined key_traldf_c3d
-      aeiu_tm(:,:,:) = aeiu(:,:,:) * fse3u(:,:,:)
-      aeiv_tm(:,:,:) = aeiv(:,:,:) * fse3v(:,:,:)
-      aeiw_tm(:,:,:) = aeiw(:,:,:) * fse3w(:,:,:)
-#  elif defined key_traldf_c2d
-      aeiu_tm(:,:)   = aeiu(:,:)
-      aeiv_tm(:,:)   = aeiv(:,:)
-      aeiw_tm(:,:)   = aeiw(:,:)
-#  elif defined key_traldf_c1d
-      aeiu_tm(:)     = aeiu(:)
-      aeiv_tm(:)     = aeiv(:)
-      aeiw_tm(:)     = aeiw(:)
-#  else
-      aeiu_tm        = aeiu
-      aeiv_tm        = aeiv
-      aeiw_tm        = aeiw
-#  endif
-#endif
+      ahtu_tm(:,:,:) = ahtu(:,:,:)
+      ahtv_tm(:,:,:) = ahtv(:,:,:)
+      !
+      IF( ln_ldfeiv ) THEN
+         aeiu_tm(:,:,:) = aeiu(:,:,:)
+         aeiv_tm(:,:,:) = aeiv(:,:,:)
+      ENDIF
       !
       IF( nn_timing == 1 )  CALL timing_stop('trc_sub_ini')
       !
    END SUBROUTINE trc_sub_ini
+
 
    SUBROUTINE trc_sub_reset( kt )
@@ -685,10 +494,10 @@
       avs   (:,:,:)   =  avs_temp   (:,:,:)
 # endif
-#if defined key_ldfslp
-      wslpi (:,:,:)   =  wslpi_temp (:,:,:)
-      wslpj (:,:,:)   =  wslpj_temp (:,:,:)
-      uslp  (:,:,:)   =  uslp_temp  (:,:,:)
-      vslp  (:,:,:)   =  vslp_temp  (:,:,:)
-#endif
+      IF( l_ldfslp ) THEN
+         wslpi (:,:,:)   =  wslpi_temp (:,:,:)
+         wslpj (:,:,:)   =  wslpj_temp (:,:,:)
+         uslp  (:,:,:)   =  uslp_temp  (:,:,:)
+         vslp  (:,:,:)   =  vslp_temp  (:,:,:)
+      ENDIF
       sshn  (:,:)     =  sshn_temp  (:,:)
       sshb  (:,:)     =  sshb_temp  (:,:)
@@ -722,49 +531,12 @@
 # endif
       !
-#if defined key_traldf_c3d
       ahtu  (:,:,:)   =  ahtu_temp  (:,:,:)
       ahtv  (:,:,:)   =  ahtv_temp  (:,:,:)
-      ahtw  (:,:,:)   =  ahtw_temp  (:,:,:)
-      ahtt  (:,:,:)   =  ahtt_temp  (:,:,:)
-#elif defined key_traldf_c2d
-      ahtu  (:,:)     =  ahtu_temp  (:,:)
-      ahtv  (:,:)     =  ahtv_temp  (:,:)
-      ahtw  (:,:)     =  ahtw_temp  (:,:)
-      ahtt  (:,:)     =  ahtt_temp  (:,:)
-#elif defined key_traldf_c1d
-      ahtu  (:)       =  ahtu_temp  (:)
-      ahtv  (:)       =  ahtv_temp  (:)
-      ahtw  (:)       =  ahtw_temp  (:)
-      ahtt  (:)       =  ahtt_temp  (:)
-#else
-      ahtu            =  ahtu_temp
-      ahtv            =  ahtv_temp
-      ahtw            =  ahtw_temp
-      ahtt            =  ahtt_temp
-#endif
-      !
-#if defined key_traldf_eiv
-#if defined key_traldf_c3d
-      aeiu  (:,:,:)  =  aeiu_temp(:,:,:)
-      aeiv  (:,:,:)  =  aeiv_temp(:,:,:)
-      aeiw  (:,:,:)  =  aeiw_temp(:,:,:)
-#elif defined key_traldf_c2d
-      aeiu  (:,:)    =  aeiu_temp(:,:)
-      aeiv  (:,:)    =  aeiv_temp(:,:)
-      aeiw  (:,:)    =  aeiw_temp(:,:)
-#elif defined key_traldf_c1d
-      aeiu  (:)      =  aeiu_temp(:)
-      aeiv  (:)      =  aeiv_temp(:)
-      aeiw  (:)      =  aeiw_temp(:)
-#else
-      aeiu           =  aeiu_temp
-      aeiv           =  aeiv_temp
-      aeiw           =  aeiw_temp
-#endif
-#endif 
+      !
+      IF( ln_ldfeiv ) THEN
+         aeiu  (:,:,:)  =  aeiu_temp(:,:,:)
+         aeiv  (:,:,:)  =  aeiv_temp(:,:,:)
+      ENDIF
       hdivn (:,:,:)   =  hdivn_temp (:,:,:)
-      rotn  (:,:,:)   =  rotn_temp  (:,:,:)
-      hdivb (:,:,:)   =  hdivb_temp (:,:,:)
-      rotb  (:,:,:)   =  rotb_temp  (:,:,:)
       hu    (:,:)     =  hu_temp    (:,:)
       hv    (:,:)     =  hv_temp    (:,:)
@@ -773,12 +545,12 @@
       !                                      
       DO jk = 1, jpk
-         fse3t(:,:,jk)= e3t_temp(:,:,jk) 
-         fse3u(:,:,jk)= e3u_temp(:,:,jk) 
-         fse3v(:,:,jk)= e3v_temp(:,:,jk) 
-         fse3w(:,:,jk)= e3w_temp(:,:,jk) 
+         fse3t(:,:,jk) = e3t_temp(:,:,jk) 
+         fse3u(:,:,jk) = e3u_temp(:,:,jk) 
+         fse3v(:,:,jk) = e3v_temp(:,:,jk) 
+         fse3w(:,:,jk) = e3w_temp(:,:,jk) 
       END DO
       !                                           !------------------------------!
       IF( lk_vvl ) THEN                           !  Update Now Vertical coord.  !   (only in vvl case)
-        !                                           !------------------------------!
+        !                                         !------------------------------!
          DO jk = 1, jpkm1
             fsdept(:,:,jk) = fsdept_n(:,:,jk)          ! now local depths stored in fsdep. arrays
@@ -803,14 +575,14 @@
          tsn_tm  (:,:,:,jp_sal) = tsn  (:,:,:,jp_sal) * fse3t(:,:,:)  
          rhop_tm (:,:,:)        = rhop (:,:,:)        * fse3t(:,:,:)  
-         avt_tm  (:,:,:)        = avt  (:,:,:)        * fse3w(:,:,:)  
+         avt_tm  (:,:,:)        = avt  (:,:,:) 
 # if defined key_zdfddm
-         avs_tm  (:,:,:)        = avs  (:,:,:)        * fse3w(:,:,:)  
-# endif
-#if defined key_ldfslp
-         wslpi_tm(:,:,:)        = wslpi(:,:,:)        * fse3w(:,:,:) 
-         wslpj_tm(:,:,:)        = wslpj(:,:,:)        * fse3w(:,:,:) 
-         uslp_tm (:,:,:)        = uslp (:,:,:)        * fse3u(:,:,:) 
-         vslp_tm (:,:,:)        = vslp (:,:,:)        * fse3v(:,:,:) 
-#endif
+         avs_tm  (:,:,:)        = avs  (:,:,:)
+# endif
+      IF( l_ldfslp ) THEN
+         wslpi_tm(:,:,:)        = wslpi(:,:,:)
+         wslpj_tm(:,:,:)        = wslpj(:,:,:)
+         uslp_tm (:,:,:)        = uslp (:,:,:)
+         vslp_tm (:,:,:)        = vslp (:,:,:)
+      ENDIF
       !
       sshb_hold  (:,:) = sshn  (:,:)
@@ -840,49 +612,12 @@
 # endif
       !
-#if defined key_traldf_c3d
       DO jk = 1, jpkm1
-         ahtt_tm(:,:,jk) = ahtt(:,:,jk) * fse3t(:,:,jk)
-         ahtu_tm(:,:,jk) = ahtu(:,:,jk) * fse3u(:,:,jk)
-         ahtv_tm(:,:,jk) = ahtv(:,:,jk) * fse3v(:,:,jk)
-         ahtw_tm(:,:,jk) = ahtw(:,:,jk) * fse3w(:,:,jk)
+         ahtu_tm(:,:,jk) = ahtu(:,:,jk)
+         ahtv_tm(:,:,jk) = ahtv(:,:,jk)
+         IF( ln_ldfeiv ) THEN
+            aeiu_tm(:,:,jk) = aeiu(:,:,jk) * fse3u(:,:,jk)
+            aeiv_tm(:,:,jk) = aeiv(:,:,jk) * fse3v(:,:,jk)
+         ENDIF
       END DO
-#elif defined key_traldf_c2d
-      ahtt_tm(:,:)   = ahtt(:,:)
-      ahtu_tm(:,:)   = ahtu(:,:)
-      ahtv_tm(:,:)   = ahtv(:,:)
-      ahtw_tm(:,:)   = ahtw(:,:)
-#elif defined key_traldf_c1d
-      ahtt_tm(:)     = ahtt(:)
-      ahtu_tm(:)     = ahtu(:)
-      ahtv_tm(:)     = ahtv(:)
-      ahtw_tm(:)     = ahtw(:)
-#else
-      ahtt_tm        = ahtt
-      ahtu_tm        = ahtu
-      ahtv_tm        = ahtv
-      ahtw_tm        = ahtw
-#endif
-      !
-#if defined key_traldf_eiv
-#  if defined key_traldf_c3d
-      DO jk = 1, jpk
-         aeiu_tm(:,:,jk) = aeiu(:,:,jk) * fse3u(:,:,jk)
-         aeiv_tm(:,:,jk) = aeiv(:,:,jk) * fse3v(:,:,jk)
-         aeiw_tm(:,:,jk) = aeiw(:,:,jk) * fse3w(:,:,jk)
-      END DO
-#  elif defined key_traldf_c2d
-      aeiu_tm(:,:)   = aeiu(:,:)
-      aeiv_tm(:,:)   = aeiv(:,:)
-      aeiw_tm(:,:)   = aeiw(:,:)
-#  elif defined key_traldf_c1d
-      aeiu_tm(:)     = aeiu(:)
-      aeiv_tm(:)     = aeiv(:)
-      aeiw_tm(:)     = aeiw(:)
-#  else
-      aeiu_tm        = aeiu
-      aeiv_tm        = aeiv
-      aeiw_tm        = aeiw
-#  endif
-#endif
       !
       IF( nn_timing == 1 )  CALL timing_stop('trc_sub_reset')
@@ -958,5 +693,6 @@
       ENDIF
       !
-      CALL div_cur( kt )                              ! Horizontal divergence & Relative vorticity
+!!gm BUG here  the recomputed divergence will includes the wrong river runoff.....
+      CALL div_hor( kt )                              ! Horizontal divergence 
       !
       z2dt = 2._wp * rdt                              ! set time step size (Euler/Leapfrog)
@@ -979,31 +715,22 @@
       CALL agrif_ssh( kt )
 #endif
-#if defined key_obc
-      IF( Agrif_Root() ) THEN 
-         ssha(:,:) = ssha(:,:) * obctmsk(:,:)
-         CALL lbc_lnk( ssha, 'T', 1. )                 ! absolutly compulsory !! (jmm)
-      ENDIF
-#endif
 #if defined key_bdy
       ssha(:,:) = ssha(:,:) * bdytmask(:,:)
       CALL lbc_lnk( ssha, 'T', 1. ) 
 #endif
-
+      !
       !                                                ! Sea Surface Height at u-,v- and f-points (vvl case only)
 !!Z~  IF( lk_vvl ) THEN                                ! (required only in key_vvl case)
 !!Z~     DO jj = 1, jpjm1
 !!Z~        DO ji = 1, jpim1      ! NO Vector Opt.
-!!Z~           sshu_a(ji,jj) = 0.5  * umask(ji,jj,1) / ( e1u(ji  ,jj) * e2u(ji  ,jj) )                   &
-!!Z~              &                                  * ( e1t(ji  ,jj) * e2t(ji  ,jj) * ssha(ji  ,jj)     &
-!!Z~              &                                    + e1t(ji+1,jj) * e2t(ji+1,jj) * ssha(ji+1,jj) )
-!!Z~           sshv_a(ji,jj) = 0.5  * vmask(ji,jj,1) / ( e1v(ji,jj  ) * e2v(ji,jj  ) )                   &
-!!Z~              &                                  * ( e1t(ji,jj  ) * e2t(ji,jj  ) * ssha(ji,jj  )     &
-!!Z~              &                                    + e1t(ji,jj+1) * e2t(ji,jj+1) * ssha(ji,jj+1) )
+!!Z~           sshu_a(ji,jj) = 0.5  * umask(ji,jj,1) * r1_e1e2u(ji  ,jj) * ( e1e2t(ji  ,jj) * ssha(ji  ,jj)   &
+!!Z~              &                                                        + e1e2t(ji+1,jj) * ssha(ji+1,jj) )
+!!Z~           sshv_a(ji,jj) = 0.5  * vmask(ji,jj,1) * r1_e1e2v(ji,jj  ) * ( e1e2t(ji,jj  ) * ssha(ji,jj  )     &
+!!Z~              &                                                        + e1e2t(ji,jj+1) * ssha(ji,jj+1) )
 !!Z~        END DO
 !!Z~     END DO
 !!Z~     CALL lbc_lnk( sshu_a, 'U', 1. )   ;   CALL lbc_lnk( sshv_a, 'V', 1. )      ! Boundaries conditions
 !!Z~  ENDIF
-      
-
+      !
       !                                           !------------------------------!
       !                                           !     Now Vertical Velocity    !
@@ -1019,5 +746,4 @@
 #endif
       END DO
-
       !
       CALL wrk_dealloc( jpi, jpj, zhdiv )
@@ -1045,11 +771,15 @@
          &      e3t_temp(jpi,jpj,jpk)       ,  e3u_temp(jpi,jpj,jpk),    &
          &      e3v_temp(jpi,jpj,jpk)       ,  e3w_temp(jpi,jpj,jpk),    &
-#if defined key_ldfslp
+!!gm here define these arrays only if l_ldfslp = True 
          &      wslpi_temp(jpi,jpj,jpk)     ,  wslpj_temp(jpi,jpj,jpk),  &
          &      uslp_temp(jpi,jpj,jpk)      ,  vslp_temp(jpi,jpj,jpk),   &
-#endif
+         &      wslpi_tm(jpi,jpj,jpk)       ,  wslpj_tm(jpi,jpj,jpk),    &
+         &      uslp_tm(jpi,jpj,jpk)        ,  vslp_tm(jpi,jpj,jpk),     &
+!!gm
 #if defined key_trabbl
          &      ahu_bbl_temp(jpi,jpj)       ,  ahv_bbl_temp(jpi,jpj),    &
          &      utr_bbl_temp(jpi,jpj)       ,  vtr_bbl_temp(jpi,jpj),    &
+         &      ahu_bbl_tm(jpi,jpj)         ,  ahv_bbl_tm(jpi,jpj),      &
+         &      utr_bbl_tm(jpi,jpj)         ,  vtr_bbl_tm(jpi,jpj),      &
 #endif
          &      rnf_temp(jpi,jpj)           ,  h_rnf_temp(jpi,jpj) ,     &
@@ -1062,37 +792,12 @@
          &      avs_tm(jpi,jpj,jpk)         ,  avs_temp(jpi,jpj,jpk) ,   &
 # endif
-#if defined key_traldf_c3d
-         &      ahtt_tm(jpi,jpj,jpk)        ,  ahtt_temp(jpi,jpj,jpk),   &
+!!gm here define these arrays only if l_ldftra_time = True 
          &      ahtu_tm(jpi,jpj,jpk)        ,  ahtu_temp(jpi,jpj,jpk),   &
          &      ahtv_tm(jpi,jpj,jpk)        ,  ahtv_temp(jpi,jpj,jpk),   &
-         &      ahtw_tm(jpi,jpj,jpk)        ,  ahtw_temp(jpi,jpj,jpk),   &
-#elif defined key_traldf_c2d
-         &      ahtt_tm(jpi,jpj)            ,  ahtt_temp(jpi,jpj),       &
-         &      ahtu_tm(jpi,jpj)            ,  ahtu_temp(jpi,jpj),       &
-         &      ahtv_tm(jpi,jpj)            ,  ahtv_temp(jpi,jpj),       &
-         &      ahtw_tm(jpi,jpj)            ,  ahtw_temp(jpi,jpj),       &
-#elif defined key_traldf_c1d
-         &      ahtt_tm(jpk)                ,  ahtt_temp(jpk),           &
-         &      ahtu_tm(jpk)                ,  ahtu_temp(jpk),           &
-         &      ahtv_tm(jpk)                ,  ahtv_temp(jpk),           &
-         &      ahtw_tm(jpk)                ,  ahtw_temp(jpk),           &
-#endif
-#if defined key_traldf_eiv
-# if defined key_traldf_c3d
+!!gm here define these arrays only if l_ldfeiv_time = True 
          &      aeiu_tm(jpi,jpj,jpk)        ,  aeiu_temp(jpi,jpj,jpk),   &
          &      aeiv_tm(jpi,jpj,jpk)        ,  aeiv_temp(jpi,jpj,jpk),   &
-         &      aeiw_tm(jpi,jpj,jpk)        ,  aeiw_temp(jpi,jpj,jpk),   &
-# elif defined key_traldf_c2d
-         &      aeiu_tm(jpi,jpj)            ,  aeiu_temp(jpi,jpj),       &
-         &      aeiv_tm(jpi,jpj)            ,  aeiv_temp(jpi,jpj),       &
-         &      aeiw_tm(jpi,jpj)            ,  aeiw_temp(jpi,jpj),       &
-# elif defined key_traldf_c1d
-         &      aeiu_tm(jpk)                ,  aeiu_temp(jpk),           &
-         &      aeiv_tm(jpk)                ,  aeiv_temp(jpk),           &
-         &      aeiw_tm(jpk)                ,  aeiw_temp(jpk),           &
-# endif
-# endif
-         &      hdivn_temp(jpi,jpj,jpk)     ,  hdivb_temp(jpi,jpj,jpk),  &
-         &      rotn_temp(jpi,jpj,jpk)      ,  rotb_temp(jpi,jpj,jpk),   &
+!!gm end
+         &      hdivn_temp(jpi,jpj,jpk)     ,                            &
          &      hu_temp(jpi,jpj)            ,  hv_temp(jpi,jpj),         &
          &      hur_temp(jpi,jpj)           ,  hvr_temp(jpi,jpj),        &
@@ -1106,12 +811,4 @@
          &      emp_b_hold(jpi,jpj)         ,                            &
          &      hmld_tm(jpi,jpj)            ,  qsr_tm(jpi,jpj) ,         &
-#if defined key_ldfslp
-         &      wslpi_tm(jpi,jpj,jpk)       ,  wslpj_tm(jpi,jpj,jpk),    &
-         &      uslp_tm(jpi,jpj,jpk)        ,  vslp_tm(jpi,jpj,jpk),     &
-#endif
-#if defined key_trabbl
-         &      ahu_bbl_tm(jpi,jpj)         ,  ahv_bbl_tm(jpi,jpj),      &
-         &      utr_bbl_tm(jpi,jpj)         ,  vtr_bbl_tm(jpi,jpj),      &
-#endif
          &      rnf_tm(jpi,jpj)             ,  h_rnf_tm(jpi,jpj) ,       &
          &                                    STAT=trc_sub_alloc )  
Index: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcwri.F90
===================================================================
--- branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcwri.F90	(revision 4616)
+++ branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/TOP_SRC/trcwri.F90	(revision 4617)
@@ -10,5 +10,5 @@
    !!   'key_top'                                           TOP models
    !!----------------------------------------------------------------------
-   !! trc_wri_trc   :  outputs of concentration fields
+   !! trc_wri      :  outputs of concentration fields
    !!----------------------------------------------------------------------
    USE dom_oce     ! ocean space and time domain variables
@@ -29,5 +29,9 @@
    !! * Substitutions
 #  include "top_substitute.h90"
-
+   !!----------------------------------------------------------------------
+   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
+   !! $Id$ 
+   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
+   !!----------------------------------------------------------------------
 CONTAINS
 
@@ -76,8 +80,4 @@
 #endif
 
-   !!----------------------------------------------------------------------
-   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
-   !! $Id$ 
-   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
    !!======================================================================
 END MODULE trcwri
