Index: branches/2015/nemo_v3_6_STABLE/NEMOGCM/CONFIG/SHARED/field_def.xml
===================================================================
--- branches/2015/nemo_v3_6_STABLE/NEMOGCM/CONFIG/SHARED/field_def.xml	(revision 8625)
+++ branches/2015/nemo_v3_6_STABLE/NEMOGCM/CONFIG/SHARED/field_def.xml	(revision 8627)
@@ -508,10 +508,10 @@
 
         <!-- variables available with key_diaar5 -->   
-        <field id="w_masstr"     long_name="Upward Ocean Mass Transport"             standard_name="upward_ocean_mass_transport"             unit="kg/s"   />
-        <field id="w_masstr2"    long_name="square of vertical mass transport"   standard_name="square_of_upward_ocean_mass_transport"   unit="kg2/s2" />
+        <field id="w_masstr"     long_name="Upward Ocean Mass Transport"           standard_name="upward_ocean_mass_transport"             unit="kg/s"   />
+        <field id="w_masstr2"    long_name="square of vertical mass transport"     standard_name="square_of_upward_ocean_mass_transport"   unit="kg2/s2" />
 
         <!-- aht2d and  aht2d_eiv: available with key_traldf_eiv and key_traldf_c2d -->
-        <field id="aht2d"        long_name="lateral eddy diffusivity"       standard_name="ocean_tracer_xy_laplacian_diffusivity"      unit="m2/s"   grid_ref="grid_W_2D" />
-        <field id="aht2d_eiv"    long_name="EIV lateral eddy diffusivity"   standard_name="ocean_tracer_bolus_laplacian_diffusivity"   unit="m2/s"   grid_ref="grid_W_2D" />
+        <field id="aht2d"        long_name="lateral eddy diffusivity"              standard_name="ocean_tracer_xy_laplacian_diffusivity"      unit="m2/s"   grid_ref="grid_W_2D" />
+        <field id="aht2d_eiv"    long_name="EIV lateral eddy diffusivity"          standard_name="ocean_tracer_bolus_laplacian_diffusivity"   unit="m2/s"   grid_ref="grid_W_2D" />
 
          <!-- Variable from eosbn2 -->
@@ -541,12 +541,12 @@
 
       <field_group id="scalar"  domain_ref="1point" >
-         <field id="voltot"     long_name="Sea Water Volume"                             standard_name="sea_water_volume"                               unit="m3"   />
-         <field id="sshtot"     long_name="global mean ssh"                              standard_name="global_average_sea_level_change"                unit="m"    />
-         <field id="sshsteric"  long_name="global mean ssh steric"                       standard_name="global_average_steric_sea_level_change"         unit="m"    />
-         <field id="sshthster"  long_name="Global Average Thermosteric Sea Level Change" standard_name="global_average_thermosteric_sea_level_change"   unit="m"    />
-         <field id="masstot"    long_name="global total mass"                            standard_name="sea_water_mass"                                 unit="kg"   />
-         <field id="temptot"    long_name="Global Average Sea Water Potential Temperature"      standard_name="sea_water_potential_temperature"         unit="degree_C" />
-         <field id="saltot"     long_name="Global Average Sea Water Salinity"                   standard_name="sea_water_salinity"                      unit="0.001" />
-         <field id="fram_trans" long_name="Sea Ice Mass Transport Through Fram Strait"   standard_name="sea_ice_transport_across_line"                  unit="kg/s" />
+         <field id="voltot"     long_name="Sea Water Volume"                                standard_name="sea_water_volume"                               unit="m3"   />
+         <field id="sshtot"     long_name="global mean ssh"                                 standard_name="global_average_sea_level_change"                unit="m"    />
+         <field id="sshsteric"  long_name="global mean ssh steric"                          standard_name="global_average_steric_sea_level_change"         unit="m"    />
+         <field id="sshthster"  long_name="Global Average Thermosteric Sea Level Change"    standard_name="global_average_thermosteric_sea_level_change"   unit="m"    />
+         <field id="masstot"    long_name="global total mass"                               standard_name="sea_water_mass"                                 unit="kg"   />
+         <field id="temptot"    long_name="Global Average Sea Water Potential Temperature"  standard_name="sea_water_potential_temperature"         unit="degree_C" />
+         <field id="saltot"     long_name="Global Average Sea Water Salinity"               standard_name="sea_water_salinity"                      unit="0.001" />
+         <field id="fram_trans" long_name="Sea Ice Mass Transport Through Fram Strait"      standard_name="sea_ice_transport_across_line"                  unit="kg/s" />
 
       	 <!-- available with ln_diahsb -->
@@ -867,5 +867,9 @@
       <field id="ketrd_convP2K" long_name="ke-trend: conversion (potential to kinetic)"      unit="W/s^3"                        />
       <field id="KE"            long_name="kinetic energy: u(n)*u(n+1)/2"                    unit="W/s^2"                        />   
-      <field id="ketrd_eiv"     long_name="Kinetic energy trend from mesoscale eddy advection"   unit="W/s^2" grid_ref="grid_T_2D" />   
+      <!-- variables available when explicit mixing is used -->
+      <field id="dispkexyfo"    long_name="KE-trend: lateral  mixing induced dissipation"   standard_name="ocean_kinetic_energy_dissipation_per_unit_area_due_to_xy_friction"                   unit="W/m^2" grid_ref="grid_T_2D" />   
+      <field id="dispkevfo"     long_name="KE-trend: vertical mixing induced dissipation"   standard_name="ocean_kinetic_energy_dissipation_per_unit_area_due_to_vertical_friction"             unit="W/m^2" grid_ref="grid_T_2D" />   
+      <!-- variables available with key_traadv_eiv and key_diaeiv -->
+      <field id="eketrd_eiv"    long_name="EKE-trend due to parameterized eddy advection"   standard_name="tendency_of_ocean_eddy_kinetic_energy_content_due_to_parameterized_eddy_advection"   unit="W/m^2" grid_ref="grid_T_2D" />   
 
       <!-- variables available with ln_PE_trd -->
Index: branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilap.F90
===================================================================
--- branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilap.F90	(revision 8625)
+++ branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilap.F90	(revision 8627)
@@ -18,7 +18,9 @@
    USE oce             ! ocean dynamics and tracers
    USE dom_oce         ! ocean space and time domain
+   USE phycst          ! physical constants
    USE ldfdyn_oce      ! ocean dynamics: lateral physics
    !
    USE in_out_manager  ! I/O manager
+   USE iom             ! I/O library
    USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
    USE wrk_nemo        ! Memory Allocation
@@ -75,8 +77,9 @@
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
       !
-      INTEGER  ::   ji, jj, jk                  ! dummy loop indices
-      REAL(wp) ::   zua, zva, zbt, ze2u, ze2v   ! temporary scalar
-      REAL(wp), POINTER, DIMENSION(:,:  ) :: zcu, zcv
-      REAL(wp), POINTER, DIMENSION(:,:,:) :: zuf, zut, zlu, zlv
+      INTEGER  ::   ji, jj, jk                       ! dummy loop indices
+      REAL(wp) ::   zua, zva, zbt, ze2u, ze2v, zzz   ! local scalar
+      REAL(wp), POINTER, DIMENSION(:,:  ) ::   zcu, zcv
+      REAL(wp), POINTER, DIMENSION(:,:,:) ::   zuf, zut, zlu, zlv
+      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   z2d   ! 2D workspace 
       !!----------------------------------------------------------------------
       !
@@ -112,9 +115,9 @@
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
-                  zlu(ji,jj,jk) = - ( zuf(ji,jj,jk) - zuf(ji,jj-1,jk) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) )   &
-                     &         + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj,jk) ) / e1u(ji,jj)
+                  zlu(ji,jj,jk) = - (   zuf(ji  ,jj,jk) -   zuf(ji,jj-1,jk) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) )   &
+                     &            + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj  ,jk) ) /   e1u(ji,jj)
    
-                  zlv(ji,jj,jk) = + ( zuf(ji,jj,jk) - zuf(ji-1,jj,jk) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) )   &
-                     &         + ( hdivb(ji,jj+1,jk) - hdivb(ji,jj,jk) ) / e2v(ji,jj)
+                  zlv(ji,jj,jk) = + (   zuf(ji,jj  ,jk) -   zuf(ji-1,jj,jk) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) )   &
+                     &            + ( hdivb(ji,jj+1,jk) - hdivb(ji  ,jj,jk) ) /   e2v(ji,jj)
                END DO
             END DO
@@ -123,8 +126,8 @@
                DO ji = fs_2, fs_jpim1   ! vector opt.
                   zlu(ji,jj,jk) = - ( rotb (ji  ,jj,jk) - rotb (ji,jj-1,jk) ) / e2u(ji,jj)   &
-                     &         + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj  ,jk) ) / e1u(ji,jj)
+                     &            + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj  ,jk) ) / e1u(ji,jj)
    
                   zlv(ji,jj,jk) = + ( rotb (ji,jj  ,jk) - rotb (ji-1,jj,jk) ) / e1v(ji,jj)   &
-                     &         + ( hdivb(ji,jj+1,jk) - hdivb(ji  ,jj,jk) ) / e2v(ji,jj)
+                     &            + ( hdivb(ji,jj+1,jk) - hdivb(ji  ,jj,jk) ) / e2v(ji,jj)
                END DO  
             END DO  
@@ -133,15 +136,40 @@
       CALL lbc_lnk( zlu, 'U', -1. )   ;   CALL lbc_lnk( zlv, 'V', -1. )   ! Boundary conditions
 
-         
+      IF( iom_use('dispkexyfo') ) THEN   ! ocean kinetic energy dissipation per unit area
+         !                               ! due to xy friction (xy=lateral) 
+         !   see NEMO_book appendix C, §C.7.2 (N.B. here averaged at t-points)
+         !   local dissipation of KE at t-point due to bilaplacian operator is given by :
+         !      + ahmu mi( zlu**2 ) + mj( ahmv zlv**2 )
+         !   Note that a sign + is used as in v3.6 ahm is negative for bilaplacian viscosity
+         !
+         ! NB: ahm is negative when bilaplacian is used
+         ALLOCATE( z2d(jpi,jpj) )
+         z2d(:,:) = 0._wp
+         DO jk = 1, jpkm1
+            DO jj = 2, jpjm1
+               DO ji = 2, jpim1
+                  z2d(ji,jj) = z2d(ji,jj)                                                                  &
+                     &  +  (  fsahmu(ji,jj,jk)*zlu(ji,jj,jk)**2 + fsahmu(ji-1,jj,jk)*zlu(ji-1,jj,jk)**2    &
+                     &      + fsahmv(ji,jj,jk)*zlv(ji,jj,jk)**2 + fsahmv(ji,jj-1,jk)*zlv(ji,jj-1,jk)**2  ) * tmask(ji,jj,jk)
+               END DO
+            END DO
+         END DO
+         zzz = 0.5_wp * rau0
+         z2d(:,:) = zzz * z2d(:,:) 
+         CALL lbc_lnk( z2d,'T', 1. )
+         CALL iom_put( 'dispkexyfo', z2d )
+         DEALLOCATE( z2d )
+      ENDIF
+      
+   
+      ! Third derivative
+      ! ----------------
+      !
       DO jk = 1, jpkm1
-   
-         ! Third derivative
-         ! ----------------
-         
+         !
          ! Multiply by the eddy viscosity coef. (at u- and v-points)
-         zlu(:,:,jk) = zlu(:,:,jk) * ( fsahmu(:,:,jk) * (1-nkahm_smag) + nkahm_smag)
-
-         zlv(:,:,jk) = zlv(:,:,jk) * ( fsahmv(:,:,jk) * (1-nkahm_smag) + nkahm_smag)
-         
+         zlu(:,:,jk) = zlu(:,:,jk) * fsahmu(:,:,jk) 
+         zlv(:,:,jk) = zlv(:,:,jk) * fsahmv(:,:,jk)
+         !         
          ! Contravariant "laplacian"
          zcu(:,:) = e1u(:,:) * zlu(:,:,jk)
@@ -152,6 +180,6 @@
             DO ji = 1, fs_jpim1   ! vector opt.
                zuf(ji,jj,jk) = fmask(ji,jj,jk) * (  zcv(ji+1,jj  ) - zcv(ji,jj)      &
-                  &                            - zcu(ji  ,jj+1) + zcu(ji,jj)  )   &
-                  &       * fse3f(ji,jj,jk) / ( e1f(ji,jj)*e2f(ji,jj) )
+                  &                               - zcu(ji  ,jj+1) + zcu(ji,jj)  )   &
+                  &       * fse3f(ji,jj,jk) * r1_e12f(ji,jj)
             END DO  
          END DO  
@@ -175,5 +203,4 @@
       END DO
 
-
       ! boundary conditions on the laplacian curl and div (zuf,zut)
 !!bug gm no need to do this 2 following lbc...
@@ -181,9 +208,5 @@
       CALL lbc_lnk( zut, 'T', 1. )
 
-      DO jk = 1, jpkm1      
-   
-         ! Bilaplacian
-         ! -----------
-
+      DO jk = 1, jpkm1               ! Bilaplacian
          DO jj = 2, jpjm1
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -197,9 +220,8 @@
                   &  + ( zut(ji,jj+1,jk) - zut(ji  ,jj,jk) ) / e2v(ji,jj)
                ! add it to the general momentum trends
-               ua(ji,jj,jk) = ua(ji,jj,jk) + zua * ( fsahmu(ji,jj,jk)*nkahm_smag +(1 -nkahm_smag ))
-               va(ji,jj,jk) = va(ji,jj,jk) + zva * ( fsahmv(ji,jj,jk)*nkahm_smag +(1 -nkahm_smag ))
-            END DO
-         END DO
-
+               ua(ji,jj,jk) = ua(ji,jj,jk) + zua
+               va(ji,jj,jk) = va(ji,jj,jk) + zva
+            END DO
+         END DO
          !                                             ! ===============
       END DO                                           !   End of slab
Index: branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_lap.F90
===================================================================
--- branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_lap.F90	(revision 8625)
+++ branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_lap.F90	(revision 8627)
@@ -17,8 +17,10 @@
    USE oce             ! ocean dynamics and tracers
    USE dom_oce         ! ocean space and time domain
+   USE phycst          ! physical constants
    USE ldfdyn_oce      ! ocean dynamics: lateral physics
    USE zdf_oce         ! ocean vertical physics
    !
    USE in_out_manager  ! I/O manager
+   USE iom             ! I/O library
    USE timing          ! Timing
 
@@ -62,4 +64,5 @@
       INTEGER  ::   ji, jj, jk             ! dummy loop indices
       REAL(wp) ::   zua, zva, ze2u, ze1v   ! local scalars
+      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   z2d   ! 2D workspace 
       !!----------------------------------------------------------------------
       !
@@ -93,4 +96,33 @@
       END DO                                           !   End of slab
       !                                                ! ===============
+      
+      IF( iom_use('dispkexyfo') ) THEN   ! ocean Kinetic Energy dissipation per unit area
+         !                               ! due to lateral friction (xy=lateral) 
+         !   see NEMO_book appendix C, §C.7.2 (N.B. here averaged at t-points)
+         !   local dissipation of KE at t-point due to laplacian operator is given by :
+         !      - ahmt hdivb**2 - mi( mj(ahmf rotb**2 e1f*e2f*e3t) ) / (e1e2t*e3t)
+         !
+         ALLOCATE( z2d(jpi,jpj) )
+         z2d(:,:) = 0._wp
+         DO jk = 1, jpkm1
+            DO jj = 2, jpjm1
+               DO ji = 2, jpim1
+                  z2d(ji,jj) = z2d(ji,jj)          -   (                                                         &
+                     &   hdivb(ji,jj,jk)**2 * fsahmt(ji,jj,jk) * fse3t_n(ji,jj,jk)                               &
+                     & + 0.25_wp * (                                                                             &
+                     &   rotb (ji  ,jj  ,jk)**2 * fsahmf(ji  ,jj  ,jk) * e12f(ji  ,jj  ) * fse3f(ji  ,jj  ,jk)   &
+                     & + rotb (ji-1,jj  ,jk)**2 * fsahmf(ji-1,jj  ,jk) * e12f(ji-1,jj  ) * fse3f(ji-1,jj  ,jk)   &
+                     & + rotb (ji  ,jj-1,jk)**2 * fsahmf(ji  ,jj-1,jk) * e12f(ji  ,jj-1) * fse3f(ji  ,jj-1,jk)   &
+                     & + rotb (ji-1,jj-1,jk)**2 * fsahmf(ji-1,jj-1,jk) * e12f(ji-1,jj-1) * fse3f(ji-1,jj-1,jk)   &
+                     &             ) * r1_e12t(ji,jj)  ) * tmask(ji,jj,jk)
+               END DO
+            END DO
+         END DO
+         z2d(:,:) = rau0 * z2d(:,:) 
+         CALL lbc_lnk( z2d,'T', 1. )
+         CALL iom_put( 'dispkexyfo', z2d )
+         DEALLOCATE( z2d )
+      ENDIF
+
       IF( nn_timing == 1 )  CALL timing_stop('dyn_ldf_lap')
       !
Index: branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynzdf_imp.F90
===================================================================
--- branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynzdf_imp.F90	(revision 8625)
+++ branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynzdf_imp.F90	(revision 8627)
@@ -20,11 +20,13 @@
    USE zdf_oce         ! ocean vertical physics
    USE phycst          ! physical constants
+   USE dynadv          ! dynamics: vector invariant versus flux form
+   USE dynspg_oce, ONLY: lk_dynspg_ts
+   USE zdfbfr          ! Bottom friction setup
+   !
    USE in_out_manager  ! I/O manager
+   USE iom             ! I/O library
    USE lib_mpp         ! MPP library
-   USE zdfbfr          ! Bottom friction setup
    USE wrk_nemo        ! Memory Allocation
    USE timing          ! Timing
-   USE dynadv          ! dynamics: vector invariant versus flux form
-   USE dynspg_oce, ONLY: lk_dynspg_ts
 
    IMPLICIT NONE
@@ -69,5 +71,6 @@
       INTEGER  ::   ikbu, ikbv   ! local integers
       REAL(wp) ::   z1_p2dt, zcoef, zzwi, zzws, zrhs   ! local scalars
-      REAL(wp) ::   ze3ua, ze3va
+      REAL(wp) ::   ze3ua, ze3va, zzz
+      REAL(wp), POINTER, DIMENSION(:,:)   ::  z2d
       REAL(wp), POINTER, DIMENSION(:,:,:) ::  zwi, zwd, zws
       !!----------------------------------------------------------------------
@@ -257,15 +260,4 @@
       END DO
 
-#if ! defined key_dynspg_ts
-      ! Normalization to obtain the general momentum trend ua
-      DO jk = 1, jpkm1
-         DO jj = 2, jpjm1   
-            DO ji = fs_2, fs_jpim1   ! vector opt.
-               ua(ji,jj,jk) = ( ua(ji,jj,jk) - ub(ji,jj,jk) ) * z1_p2dt
-            END DO
-         END DO
-      END DO
-#endif
-
       ! 3. Vertical diffusion on v
       ! ---------------------------
@@ -357,9 +349,38 @@
       END DO
 
-      ! Normalization to obtain the general momentum trend va
+      IF( iom_use( 'dispkevfo' ) ) THEN   ! ocean kinetic energy dissipation per unit area
+         !                                ! due to v friction (v=vertical) 
+         !                                ! see NEMO_book appendix C, §C.8 (N.B. here averaged at t-points)
+         !                                ! Note that formally, in a Leap-Frog environment, the shear**2 should be the product of 
+         !                                ! now by before shears, i.e. the source term of TKE (local positivity is not ensured).
+         CALL wrk_alloc(jpi,jpj,   z2d )
+         z2d(:,:) = 0._wp
+         DO jk = 1, jpkm1
+            DO jj = 2, jpjm1
+               DO ji = 2, jpim1
+                  z2d(ji,jj) = z2d(ji,jj)  +  (                                                                                  &
+                     &   avmu(ji  ,jj,jk) * ( ua(ji  ,jj,jk-1) - ua(ji  ,jj,jk) )**2 / fse3uw(ji  ,jj,jk) * wumask(ji  ,jj,jk)   &
+                     & + avmu(ji-1,jj,jk) * ( ua(ji-1,jj,jk-1) - ua(ji-1,jj,jk) )**2 / fse3uw(ji-1,jj,jk) * wumask(ji-1,jj,jk)   &
+                     & + avmv(ji,jj  ,jk) * ( va(ji,jj  ,jk-1) - va(ji,jj  ,jk) )**2 / fse3vw(ji,jj  ,jk) * wvmask(ji,jj  ,jk)   &
+                     & + avmv(ji,jj-1,jk) * ( va(ji,jj-1,jk-1) - va(ji,jj-1,jk) )**2 / fse3vw(ji,jj-1,jk) * wvmask(ji,jj-1,jk)   &
+                     &                        )
+               END DO
+            END DO
+         END DO
+         zzz= - 0.5_wp* rau0           ! caution sign minus here
+         z2d(:,:) = zzz * z2d(:,:) 
+         CALL lbc_lnk( z2d,'T', 1. )
+         CALL iom_put( 'dispkevfo', z2d )
+         CALL wrk_dealloc(jpi,jpj,   z2d )
+      ENDIF
+
 #if ! defined key_dynspg_ts
+!!gm this can be removed if tranxt is changed like in the trunk so that implicit outcome with 
+!!gm the after velocity, not a trend
+      ! Normalization to obtain the general momentum trend ua
       DO jk = 1, jpkm1
          DO jj = 2, jpjm1   
             DO ji = fs_2, fs_jpim1   ! vector opt.
+               ua(ji,jj,jk) = ( ua(ji,jj,jk) - ub(ji,jj,jk) ) * z1_p2dt
                va(ji,jj,jk) = ( va(ji,jj,jk) - vb(ji,jj,jk) ) * z1_p2dt
             END DO
@@ -393,4 +414,5 @@
       CALL wrk_dealloc( jpi,jpj,jpk, zwi, zwd, zws) 
       !
+      !
       IF( nn_timing == 1 )  CALL timing_stop('dyn_zdf_imp')
       !
Index: branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_eiv.F90
===================================================================
--- branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_eiv.F90	(revision 8625)
+++ branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_eiv.F90	(revision 8627)
@@ -165,13 +165,15 @@
          CALL iom_put( "voce_eiv", v_eiv )    ! j-eiv current
          CALL iom_put( "woce_eiv", w_eiv )    ! vert. eiv current
+         !
          IF( iom_use('weiv_masstr') ) THEN   ! vertical mass transport & its square value
-           z2d(:,:) = rau0 * e12t(:,:)
-           DO jk = 1, jpk
-              z3d(:,:,jk) = w_eiv(:,:,jk) * z2d(:,:)
-           END DO
-           CALL iom_put( "weiv_masstr" , z3d )  
-         ENDIF
+            z2d(:,:) = rau0 * e12t(:,:)
+            DO jk = 1, jpk
+               z3d(:,:,jk) = w_eiv(:,:,jk) * z2d(:,:)
+            END DO
+            CALL iom_put( "weiv_masstr" , z3d )  
+         ENDIF
+         !
          IF( iom_use("ueiv_masstr") .OR. iom_use("ueiv_heattr") .OR. iom_use('ueiv_heattr3d')        &
-                                    .OR. iom_use("ueiv_salttr") .OR. iom_use('ueiv_salttr3d') ) THEN
+            &                       .OR. iom_use("ueiv_salttr") .OR. iom_use('ueiv_salttr3d') ) THEN
             z3d(:,:,jpk) = 0.e0
             z2d(:,:) = 0.e0
@@ -182,5 +184,5 @@
             CALL iom_put( "ueiv_masstr", z3d )                  ! mass transport in i-direction
          ENDIF
-
+         !
          IF( iom_use('ueiv_heattr') .OR. iom_use('ueiv_heattr3d') ) THEN
             zztmp = 0.5 * rcp 
@@ -204,5 +206,5 @@
             ENDIF
          ENDIF
-
+         !
          IF( iom_use('ueiv_salttr') .OR. iom_use('ueiv_salttr3d') ) THEN
             zztmp = 0.5 * 0.001
@@ -226,5 +228,5 @@
             ENDIF
          ENDIF
-
+         !
          IF( iom_use("veiv_masstr") .OR. iom_use("veiv_heattr") .OR. iom_use('veiv_heattr3d')       &
                                     .OR. iom_use("veiv_salttr") .OR. iom_use('veiv_salttr3d') ) THEN
@@ -235,5 +237,5 @@
             CALL iom_put( "veiv_masstr", z3d )                  ! mass transport in j-direction
          ENDIF
-            
+         !   
          IF( iom_use('veiv_heattr') .OR. iom_use('veiv_heattr3d') ) THEN
             zztmp = 0.5 * rcp 
@@ -257,5 +259,5 @@
             ENDIF
          ENDIF
-
+         !
          IF( iom_use('veiv_salttr') .OR. iom_use('veiv_salttr3d') ) THEN
             zztmp = 0.5 * 0.001
@@ -279,5 +281,5 @@
             ENDIF
          ENDIF
-
+         !
          IF( iom_use('weiv_masstr') .OR. iom_use('weiv_heattr3d') .OR. iom_use('weiv_salttr3d')) THEN   ! vertical mass transport & its square value
            z2d(:,:) = rau0 * e12t(:,:)
@@ -287,5 +289,5 @@
            CALL iom_put( "weiv_masstr" , z3d )                  ! mass transport in k-direction
          ENDIF
-
+         !
          IF( iom_use('weiv_heattr3d') ) THEN
             zztmp = 0.5 * rcp 
@@ -300,5 +302,5 @@
             CALL iom_put( "weiv_heattr3d", zztmp * z3d_T )                 ! 3D heat transport in k-direction
          ENDIF
-
+         !
          IF( iom_use('weiv_salttr3d') ) THEN
             zztmp = 0.5 * 0.001 
@@ -313,31 +315,50 @@
             CALL iom_put( "weiv_salttr3d", zztmp * z3d_T )                 ! 3D salt transport in k-direction
          ENDIF
-
-    END IF
-!
-    IF( ln_diaptr .AND. cdtype == 'TRA' ) THEN
-       z3d(:,:,:) = 0._wp
-       DO jk = 1, jpkm1
-          DO jj = 2, jpjm1
-             DO ji = fs_2, fs_jpim1   ! vector opt.
-                z3d(ji,jj,jk) = v_eiv(ji,jj,jk) * 0.5 * (tsn(ji,jj,jk,jp_tem)+tsn(ji,jj+1,jk,jp_tem)) &
-                &             * e1v(ji,jj) * fse3v(ji,jj,jk)
-             END DO
-          END DO
-       END DO
-       CALL dia_ptr_ohst_components( jp_tem, 'eiv', z3d )
-       z3d(:,:,:) = 0._wp
-       DO jk = 1, jpkm1
-          DO jj = 2, jpjm1
-             DO ji = fs_2, fs_jpim1   ! vector opt.
-                z3d(ji,jj,jk) = v_eiv(ji,jj,jk) * 0.5 * (tsn(ji,jj,jk,jp_sal)+tsn(ji,jj+1,jk,jp_sal)) &
-                &             * e1v(ji,jj) * fse3v(ji,jj,jk)
-             END DO
-          END DO
-       END DO
-       CALL dia_ptr_ohst_components( jp_sal, 'eiv', z3d )
-    ENDIF
-
-    IF( ln_KE_trd ) CALL trd_dyn(u_eiv, v_eiv, jpdyn_eivke, kt )
+         !
+         IF( ln_diaptr ) THEN
+            z3d(:,:,:) = 0._wp
+            DO jk = 1, jpkm1
+               DO jj = 2, jpjm1
+                  DO ji = fs_2, fs_jpim1   ! vector opt.
+                     z3d(ji,jj,jk) = v_eiv(ji,jj,jk) * 0.5 * (tsn(ji,jj,jk,jp_tem)+tsn(ji,jj+1,jk,jp_tem)) &
+                        &             * e1v(ji,jj) * fse3v(ji,jj,jk)
+                  END DO
+               END DO
+            END DO
+            CALL dia_ptr_ohst_components( jp_tem, 'eiv', z3d )
+            z3d(:,:,:) = 0._wp
+            DO jk = 1, jpkm1
+               DO jj = 2, jpjm1
+                  DO ji = fs_2, fs_jpim1   ! vector opt.
+                     z3d(ji,jj,jk) = v_eiv(ji,jj,jk) * 0.5 * (tsn(ji,jj,jk,jp_sal)+tsn(ji,jj+1,jk,jp_sal)) &
+                     &             * e1v(ji,jj) * fse3v(ji,jj,jk)
+                  END DO
+               END DO
+            END DO
+            CALL dia_ptr_ohst_components( jp_sal, 'eiv', z3d )
+         ENDIF
+         !
+!!gm add CMIP6 diag here instead of been done in trdken.F90
+         !
+         IF( iom_use('eketrd_eiv') ) THEN     ! tendency of EKE from parameterized eddy advection
+            ! CMIP6 diagnostic tknebto = tendency of EKE from parameterized mesoscale eddy advection
+            ! = vertical_integral( k (N S)^2 ) rho dz   where rho = rau0 and S = isoneutral slope.
+            z2d(:,:) = 0._wp
+            DO jk = 1, jpkm1
+               DO ji = 1, jpi
+                  DO jj = 1,jpj
+                     z2d(ji,jj) = z2d(ji,jj) + rau0 * fsaeiw(ji,jj,jk)                 &
+                        &                    * rn2b(ji,jj,jk) * fse3w(ji,jj,jk)        &
+                        &                    * (  wslpi(ji,jj,jk) * wslpi(ji,jj,jk)    &
+                        &                       + wslpj(ji,jj,jk) * wslpj(ji,jj,jk)  ) * wmask(ji,jj,jk)
+                  END DO
+               END DO
+            END DO
+            CALL iom_put( "eketrd_eiv", z2d )
+         ENDIF
+         !
+!!gm  removed from trdken.F90    IF( ln_KE_trd )   CALL trd_dyn(u_eiv, v_eiv, jpdyn_eivke, kt )
+         !
+      ENDIF
 # endif  
 
Index: branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRD/trd_oce.F90
===================================================================
--- branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRD/trd_oce.F90	(revision 8625)
+++ branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRD/trd_oce.F90	(revision 8627)
@@ -76,5 +76,4 @@
    INTEGER, PUBLIC, PARAMETER ::   jpdyn_spgflt  = 14  !: filter contribution to surface pressure gradient (spg_flt)
    INTEGER, PUBLIC, PARAMETER ::   jpdyn_spgexp  = 15  !: explicit contribution to surface pressure gradient (spg_flt)
-   INTEGER, PUBLIC, PARAMETER ::   jpdyn_eivke   = 16  !: K.E trend from Gent McWilliams scheme
    !
    !!----------------------------------------------------------------------
Index: branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRD/trdken.F90
===================================================================
--- branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRD/trdken.F90	(revision 8625)
+++ branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/TRD/trdken.F90	(revision 8627)
@@ -195,27 +195,29 @@
                     CALL ken_p2k( kt , zke )
                       CALL iom_put( "ketrd_convP2K", zke )     ! conversion -rau*g*w
-        CASE( jpdyn_eivke )
-            ! CMIP6 diagnostic tknebto = tendency of KE from
-            ! parameterized mesoscale eddy advection
-            ! = vertical_integral( k (N S)^2 ) rho dz
-            ! rho = reference density
-            ! S = isoneutral slope.
-            ! Most terms are on W grid so work on this grid
-#ifdef key_traldf_eiv
-            CALL wrk_alloc( jpi, jpj, zke2d )
-            zke2d(:,:) = 0._wp
-            DO jk = 1,jpk
-               DO ji = 1,jpi
-                  DO jj = 1,jpj
-                     zke2d(ji,jj) = zke2d(ji,jj) +  rau0 * fsaeiw(ji, jj, jk)               &
-                          &                      * ( wslpi(ji, jj, jk) * wslpi(ji,jj,jk)    &
-                          &                      +   wslpj(ji, jj, jk) * wslpj(ji,jj,jk) )  &
-                          &                      *   rn2(ji,jj,jk) * fse3w(ji, jj, jk)
-                  ENDDO
-               ENDDO
-            ENDDO
-            CALL iom_put("ketrd_eiv", zke2d)
-            CALL wrk_dealloc( jpi, jpj, zke2d )
-#endif
+!!gm moved in traadv_eiv ===>>> diag becomes accessible without ln_trdtra=T
+!        CASE( jpdyn_eivke )
+!            ! CMIP6 diagnostic tknebto = tendency of EKE from
+!            ! parameterized mesoscale eddy advection
+!            ! = vertical_integral( k (N S)^2 ) rho dz
+!            ! rho = reference density
+!            ! S = isoneutral slope.
+!            ! Most terms are on W grid so work on this grid
+!#ifdef key_traldf_eiv
+!            CALL wrk_alloc( jpi, jpj, zke2d )
+!            zke2d(:,:) = 0._wp
+!            DO jk = 1,jpk
+!               DO ji = 1,jpi
+!                  DO jj = 1,jpj
+!                     zke2d(ji,jj) = zke2d(ji,jj) +  rau0 * fsaeiw(ji, jj, jk)               &
+!                          &                      * ( wslpi(ji, jj, jk) * wslpi(ji,jj,jk)    &
+!                          &                      +   wslpj(ji, jj, jk) * wslpj(ji,jj,jk) )  &
+!                          &                      *   rn2(ji,jj,jk) * fse3w(ji, jj, jk)
+!                  ENDDO
+!               ENDDO
+!            ENDDO
+!            CALL iom_put("ketrd_eiv", zke2d)
+!            CALL wrk_dealloc( jpi, jpj, zke2d )
+!#endif
+!!gm end
          !
       END SELECT
