Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynadv_cen2.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynadv_cen2.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynadv_cen2.F90	(revision 10946)
@@ -92,5 +92,5 @@
          zfu_uw(:,:,:) = puu(:,:,:,Krhs) - zfu_uw(:,:,:)
          zfv_vw(:,:,:) = pvv(:,:,:,Krhs) - zfv_vw(:,:,:)
-         CALL trd_dyn( zfu_uw, zfv_vw, jpdyn_keg, kt )
+         CALL trd_dyn( zfu_uw, zfv_vw, jpdyn_keg, kt, Kmm )
          zfu_t(:,:,:) = puu(:,:,:,Krhs)
          zfv_t(:,:,:) = pvv(:,:,:,Krhs)
@@ -138,5 +138,5 @@
          zfu_t(:,:,:) = puu(:,:,:,Krhs) - zfu_t(:,:,:)
          zfv_t(:,:,:) = pvv(:,:,:,Krhs) - zfv_t(:,:,:)
-         CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt )
+         CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt, Kmm )
       ENDIF
       !                                   ! Control print
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynadv_ubs.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynadv_ubs.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynadv_ubs.F90	(revision 10946)
@@ -183,5 +183,5 @@
          zfu_uw(:,:,:) = puu(:,:,:,Krhs) - zfu_uw(:,:,:)
          zfv_vw(:,:,:) = pvv(:,:,:,Krhs) - zfv_vw(:,:,:)
-         CALL trd_dyn( zfu_uw, zfv_vw, jpdyn_keg, kt )
+         CALL trd_dyn( zfu_uw, zfv_vw, jpdyn_keg, kt, Kmm )
          zfu_t(:,:,:) = puu(:,:,:,Krhs)
          zfv_t(:,:,:) = pvv(:,:,:,Krhs)
@@ -231,5 +231,5 @@
          zfu_t(:,:,:) = puu(:,:,:,Krhs) - zfu_t(:,:,:)
          zfv_t(:,:,:) = pvv(:,:,:,Krhs) - zfv_t(:,:,:)
-         CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt )
+         CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt, Kmm )
       ENDIF
       !                                         ! Control print
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynhpg.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynhpg.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynhpg.F90	(revision 10946)
@@ -118,5 +118,5 @@
          ztrdu(:,:,:) = puu(:,:,:,Krhs) - ztrdu(:,:,:)
          ztrdv(:,:,:) = pvv(:,:,:,Krhs) - ztrdv(:,:,:)
-         CALL trd_dyn( ztrdu, ztrdv, jpdyn_hpg, kt )
+         CALL trd_dyn( ztrdu, ztrdv, jpdyn_hpg, kt, Kmm )
          DEALLOCATE( ztrdu , ztrdv )
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynkeg.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynkeg.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynkeg.F90	(revision 10946)
@@ -183,5 +183,5 @@
          ztrdu(:,:,:) = puu(:,:,:,Krhs) - ztrdu(:,:,:)
          ztrdv(:,:,:) = pvv(:,:,:,Krhs) - ztrdv(:,:,:)
-         CALL trd_dyn( ztrdu, ztrdv, jpdyn_keg, kt )
+         CALL trd_dyn( ztrdu, ztrdv, jpdyn_keg, kt, Kmm )
          DEALLOCATE( ztrdu , ztrdv )
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynldf.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynldf.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynldf.F90	(revision 10946)
@@ -78,5 +78,5 @@
          ztrdu(:,:,:) = puu(:,:,:,Krhs) - ztrdu(:,:,:)
          ztrdv(:,:,:) = pvv(:,:,:,Krhs) - ztrdv(:,:,:)
-         CALL trd_dyn( ztrdu, ztrdv, jpdyn_ldf, kt )
+         CALL trd_dyn( ztrdu, ztrdv, jpdyn_ldf, kt, Kmm )
          DEALLOCATE ( ztrdu , ztrdv )
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynnxt.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynnxt.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynnxt.F90	(revision 10946)
@@ -64,5 +64,5 @@
 CONTAINS
 
-   SUBROUTINE dyn_nxt ( kt )
+   SUBROUTINE dyn_nxt ( kt, Kmm )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE dyn_nxt  ***
@@ -93,4 +93,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
+      INTEGER, INTENT( in ) ::   Kmm     ! time level index
       !
       INTEGER  ::   ji, jj, jk   ! dummy loop indices
@@ -157,5 +158,5 @@
          !
          !                                  ! Kinetic energy and Conversion
-         IF( ln_KE_trd  )   CALL trd_dyn( ua, va, jpdyn_ken, kt )
+         IF( ln_KE_trd  )   CALL trd_dyn( ua, va, jpdyn_ken, kt, Kmm )
          !
          IF( ln_dyn_trd ) THEN              ! 3D output: total momentum trends
@@ -365,5 +366,5 @@
          zua(:,:,:) = ( ub(:,:,:) - zua(:,:,:) ) * z1_2dt
          zva(:,:,:) = ( vb(:,:,:) - zva(:,:,:) ) * z1_2dt
-         CALL trd_dyn( zua, zva, jpdyn_atf, kt )
+         CALL trd_dyn( zua, zva, jpdyn_atf, kt, Kmm )
       ENDIF
       !
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynspg.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynspg.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynspg.F90	(revision 10946)
@@ -171,5 +171,5 @@
          ztrdu(:,:,:) = puu(:,:,:,Krhs) - ztrdu(:,:,:)
          ztrdv(:,:,:) = pvv(:,:,:,Krhs) - ztrdv(:,:,:)
-         CALL trd_dyn( ztrdu, ztrdv, jpdyn_spg, kt )
+         CALL trd_dyn( ztrdu, ztrdv, jpdyn_spg, kt, Kmm )
          DEALLOCATE( ztrdu , ztrdv ) 
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynvor.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynvor.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynvor.F90	(revision 10946)
@@ -135,5 +135,5 @@
          ztrdu(:,:,:) = puu(:,:,:,Krhs) - ztrdu(:,:,:)
          ztrdv(:,:,:) = pvv(:,:,:,Krhs) - ztrdv(:,:,:)
-         CALL trd_dyn( ztrdu, ztrdv, jpdyn_pvo, kt )
+         CALL trd_dyn( ztrdu, ztrdv, jpdyn_pvo, kt, Kmm )
          !
          IF( n_dynadv /= np_LIN_dyn ) THEN   !* relative vorticity or metric trend (only in non-linear case)
@@ -149,5 +149,5 @@
             ztrdu(:,:,:) = puu(:,:,:,Krhs) - ztrdu(:,:,:)
             ztrdv(:,:,:) = pvv(:,:,:,Krhs) - ztrdv(:,:,:)
-            CALL trd_dyn( ztrdu, ztrdv, jpdyn_rvo, kt )
+            CALL trd_dyn( ztrdu, ztrdv, jpdyn_rvo, kt, Kmm )
          ENDIF
          !
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynzad.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynzad.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynzad.F90	(revision 10946)
@@ -112,5 +112,5 @@
          ztrdu(:,:,:) = puu(:,:,:,Krhs) - ztrdu(:,:,:)
          ztrdv(:,:,:) = pvv(:,:,:,Krhs) - ztrdv(:,:,:)
-         CALL trd_dyn( ztrdu, ztrdv, jpdyn_zad, kt )
+         CALL trd_dyn( ztrdu, ztrdv, jpdyn_zad, kt, Kmm )
          DEALLOCATE( ztrdu, ztrdv ) 
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynzdf.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynzdf.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynzdf.F90	(revision 10946)
@@ -98,5 +98,5 @@
       !
       !                             !* explicit top/bottom drag case
-      IF( .NOT.ln_drgimp )   CALL zdf_drg_exp( kt, puu(:,:,:,Kbb), pvv(:,:,:,Kbb), puu(:,:,:,Krhs), pvv(:,:,:,Krhs) )  ! add top/bottom friction trend to (puu(:,:,:,Kaa),pvv(:,:,:,Kaa))
+      IF( .NOT.ln_drgimp )   CALL zdf_drg_exp( kt, Kmm, puu(:,:,:,Kbb), pvv(:,:,:,Kbb), puu(:,:,:,Krhs), pvv(:,:,:,Krhs) )  ! add top/bottom friction trend to (puu(Kaa),pvv(Kaa))
       !
       !
@@ -490,5 +490,5 @@
          ztrdu(:,:,:) = ( puu(:,:,:,Kaa) - puu(:,:,:,Kbb) ) / r2dt - ztrdu(:,:,:)
          ztrdv(:,:,:) = ( pvv(:,:,:,Kaa) - pvv(:,:,:,Kbb) ) / r2dt - ztrdv(:,:,:)
-         CALL trd_dyn( ztrdu, ztrdv, jpdyn_zdf, kt )
+         CALL trd_dyn( ztrdu, ztrdv, jpdyn_zdf, kt, Kmm )
          DEALLOCATE( ztrdu, ztrdv ) 
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/LDF/ldftra.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/LDF/ldftra.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/LDF/ldftra.F90	(revision 10946)
@@ -727,5 +727,5 @@
 
 
-   SUBROUTINE ldf_eiv_trp( kt, kit000, pu, pv, pw, cdtype, Kmm )
+   SUBROUTINE ldf_eiv_trp( kt, kit000, pu, pv, pw, cdtype, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE ldf_eiv_trp  ***
@@ -745,8 +745,8 @@
       !! ** Action  : pu, pv increased by the eiv transport
       !!----------------------------------------------------------------------
-      INTEGER                         , INTENT(in   ) ::   kt       ! ocean time-step index
-      INTEGER                         , INTENT(in   ) ::   kit000   ! first time step index
-      INTEGER                         , INTENT(in   ) ::   Kmm      ! ocean time level indices
-      CHARACTER(len=3)                , INTENT(in   ) ::   cdtype   ! =TRA or TRC (tracer indicator)
+      INTEGER                         , INTENT(in   ) ::   kt        ! ocean time-step index
+      INTEGER                         , INTENT(in   ) ::   kit000    ! first time step index
+      INTEGER                         , INTENT(in   ) ::   Kmm, Krhs ! ocean time level indices
+      CHARACTER(len=3)                , INTENT(in   ) ::   cdtype    ! =TRA or TRC (tracer indicator)
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pu      ! in : 3 ocean transport components   [m3/s]
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pv      ! out: 3 ocean transport components   [m3/s]
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/SBC/sbcmod.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/SBC/sbcmod.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/SBC/sbcmod.F90	(revision 10946)
@@ -422,5 +422,5 @@
       SELECT CASE( nsbc )                                ! Compute ocean surface boundary condition
       !                                                  ! (i.e. utau,vtau, qns, qsr, emp, sfx)
-      CASE( jp_usr   )     ;   CALL usrdef_sbc_oce( kt )                             ! user defined formulation 
+      CASE( jp_usr   )     ;   CALL usrdef_sbc_oce( kt, Kbb )                        ! user defined formulation 
       CASE( jp_flx     )   ;   CALL sbc_flx       ( kt )                             ! flux formulation
       CASE( jp_blk     )
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv.F90	(revision 10946)
@@ -127,5 +127,5 @@
       !
       IF( ln_ldfeiv .AND. .NOT. ln_traldf_triad )   &
-         &              CALL ldf_eiv_trp( kt, nit000, zuu, zvv, zww, 'TRA', Kmm )   ! add the eiv transport (if necessary)
+         &              CALL ldf_eiv_trp( kt, nit000, zuu, zvv, zww, 'TRA', Kmm, Krhs )   ! add the eiv transport (if necessary)
       !
       IF( ln_mle    )   CALL tra_mle_trp( kt, nit000, zuu, zvv, zww, 'TRA' )   ! add the mle transport (if necessary)
@@ -165,6 +165,6 @@
             ztrds(:,:,jk) = pts(:,:,jk,jp_sal,Krhs) - ztrds(:,:,jk)
          END DO
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_totad, ztrdt )
-         CALL trd_tra( kt, 'TRA', jp_sal, jptra_totad, ztrds )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_totad, ztrdt )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_totad, ztrds )
          DEALLOCATE( ztrdt, ztrds )
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_cen.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_cen.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_cen.F90	(revision 10946)
@@ -191,7 +191,7 @@
          !                             ! trend diagnostics
          IF( l_trd ) THEN
-            CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kmm) )
-            CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kmm) )
-            CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pW, pt(:,:,:,jn,Kmm) )
+            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kmm) )
+            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kmm) )
+            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zwz, pW, pt(:,:,:,jn,Kmm) )
          END IF
          !                                 ! "Poleward" heat and salt transports 
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_fct.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_fct.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_fct.F90	(revision 10946)
@@ -303,7 +303,7 @@
             !
             IF( l_trd ) THEN              ! trend diagnostics
-               CALL trd_tra( kt, cdtype, jn, jptra_xad, ztrdx, pU, pt(:,:,:,jn,Kmm) )
-               CALL trd_tra( kt, cdtype, jn, jptra_yad, ztrdy, pV, pt(:,:,:,jn,Kmm) )
-               CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pW, pt(:,:,:,jn,Kmm) )
+               CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztrdx, pU, pt(:,:,:,jn,Kmm) )
+               CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztrdy, pV, pt(:,:,:,jn,Kmm) )
+               CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, ztrdz, pW, pt(:,:,:,jn,Kmm) )
             ENDIF
             !                             ! heat/salt transport
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_mus.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_mus.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_mus.F90	(revision 10946)
@@ -201,6 +201,6 @@
          !                                ! trend diagnostics
          IF( l_trd )  THEN
-            CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kbb) )
-            CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kbb) )
+            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kbb) )
+            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kbb) )
          END IF
          !                                 ! "Poleward" heat and salt transports 
@@ -268,5 +268,5 @@
          END DO
          !                                ! send trends for diagnostic
-         IF( l_trd )  CALL trd_tra( kt, cdtype, jn, jptra_zad, zwx, pW, pt(:,:,:,jn,Kbb) )
+         IF( l_trd )  CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zwx, pW, pt(:,:,:,jn,Kbb) )
          !
       END DO                     ! end of tracer loop
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_qck.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_qck.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_qck.F90	(revision 10946)
@@ -220,5 +220,5 @@
          END DO
          !                                 ! trend diagnostics
-         IF( l_trd )   CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kmm) )
+         IF( l_trd )   CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kmm) )
          !
       END DO
@@ -336,5 +336,5 @@
          END DO
          !                                 ! trend diagnostics
-         IF( l_trd )   CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kmm) )
+         IF( l_trd )   CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kmm) )
          !                                 ! "Poleward" heat and salt transports (contribution of upstream fluxes)
          IF( l_ptr )   CALL dia_ptr_hst( jn, 'adv', zwy(:,:,:) )
@@ -395,5 +395,5 @@
          END DO
          !                                 ! Send trends for diagnostic
-         IF( l_trd )  CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pW, pt(:,:,:,jn,Kmm) )
+         IF( l_trd )  CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zwz, pW, pt(:,:,:,jn,Kmm) )
          !
       END DO
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_ubs.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_ubs.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_ubs.F90	(revision 10946)
@@ -177,6 +177,6 @@
          !                
          IF( l_trd ) THEN                  ! trend diagnostics
-             CALL trd_tra( kt, cdtype, jn, jptra_xad, ztu, pU, pt(:,:,:,jn,Kmm) )
-             CALL trd_tra( kt, cdtype, jn, jptra_yad, ztv, pV, pt(:,:,:,jn,Kmm) )
+             CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztu, pU, pt(:,:,:,jn,Kmm) )
+             CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztv, pV, pt(:,:,:,jn,Kmm) )
          END IF
          !     
@@ -273,5 +273,5 @@
                END DO
             END DO
-            CALL trd_tra( kt, cdtype, jn, jptra_zad, zltv )
+            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zltv )
          ENDIF
          !
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trabbc.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trabbc.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trabbc.F90	(revision 10946)
@@ -51,5 +51,5 @@
 CONTAINS
 
-   SUBROUTINE tra_bbc( kt )
+   SUBROUTINE tra_bbc( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE tra_bbc  ***
@@ -74,4 +74,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
+      INTEGER, INTENT(in) ::   Kmm, Krhs  ! time level indices
       !
       INTEGER  ::   ji, jj    ! dummy loop indices
@@ -96,5 +97,5 @@
       IF( l_trdtra ) THEN        ! Send the trend for diagnostics
          ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_bbc, ztrdt )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_bbc, ztrdt )
          DEALLOCATE( ztrdt )
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trabbl.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trabbl.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trabbl.F90	(revision 10946)
@@ -89,5 +89,5 @@
 
 
-   SUBROUTINE tra_bbl( kt )
+   SUBROUTINE tra_bbl( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE bbl  ***
@@ -102,4 +102,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT( in ) ::   kt   ! ocean time-step
+      INTEGER, INTENT( in ) ::   Kmm, Krhs  ! time level indices
       !
       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrdt, ztrds
@@ -145,6 +146,6 @@
          ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
          ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_bbl, ztrdt )
-         CALL trd_tra( kt, 'TRA', jp_sal, jptra_bbl, ztrds )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_bbl, ztrdt )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_bbl, ztrds )
          DEALLOCATE( ztrdt, ztrds )
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tradmp.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tradmp.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tradmp.F90	(revision 10946)
@@ -72,5 +72,5 @@
 
 
-   SUBROUTINE tra_dmp( kt )
+   SUBROUTINE tra_dmp( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                   ***  ROUTINE tra_dmp  ***
@@ -91,4 +91,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
+      INTEGER, INTENT(in) ::   Kmm, Krhs  ! time level indices
       !
       INTEGER ::   ji, jj, jk, jn   ! dummy loop indices
@@ -151,6 +152,6 @@
       IF( l_trdtra )   THEN       ! trend diagnostic
          ztrdts(:,:,:,:) = tsa(:,:,:,:) - ztrdts(:,:,:,:)
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_dmp, ztrdts(:,:,:,jp_tem) )
-         CALL trd_tra( kt, 'TRA', jp_sal, jptra_dmp, ztrdts(:,:,:,jp_sal) )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_dmp, ztrdts(:,:,:,jp_tem) )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_dmp, ztrdts(:,:,:,jp_sal) )
          DEALLOCATE( ztrdts ) 
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traldf.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traldf.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traldf.F90	(revision 10946)
@@ -47,5 +47,5 @@
 CONTAINS
 
-   SUBROUTINE tra_ldf( kt, Kmm )
+   SUBROUTINE tra_ldf( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE tra_ldf  ***
@@ -54,5 +54,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
-      INTEGER, INTENT( in ) ::   Kmm  ! ocean time level indices
+      INTEGER, INTENT( in ) ::   Kmm, Krhs  ! ocean time level indices
       !!
       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrdt, ztrds
@@ -81,6 +81,6 @@
          ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
          ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_ldf, ztrdt )
-         CALL trd_tra( kt, 'TRA', jp_sal, jptra_ldf, ztrds )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_ldf, ztrdt )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_ldf, ztrds )
          DEALLOCATE( ztrdt, ztrds ) 
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tranpc.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tranpc.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tranpc.F90	(revision 10946)
@@ -42,5 +42,5 @@
 CONTAINS
 
-   SUBROUTINE tra_npc( kt )
+   SUBROUTINE tra_npc( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE tranpc  ***
@@ -59,4 +59,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
+      INTEGER, INTENT(in) ::   Kmm, Krhs  ! time level indices
       !
       INTEGER  ::   ji, jj, jk   ! dummy loop indices
@@ -304,6 +305,6 @@
             ztrdt(:,:,:) = ( tsa(:,:,:,jp_tem) - ztrdt(:,:,:) ) * z1_r2dt
             ztrds(:,:,:) = ( tsa(:,:,:,jp_sal) - ztrds(:,:,:) ) * z1_r2dt
-            CALL trd_tra( kt, 'TRA', jp_tem, jptra_npc, ztrdt )
-            CALL trd_tra( kt, 'TRA', jp_sal, jptra_npc, ztrds )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_npc, ztrdt )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_npc, ztrds )
             DEALLOCATE( ztrdt, ztrds )
          ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tranxt.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tranxt.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tranxt.F90	(revision 10946)
@@ -64,5 +64,5 @@
 CONTAINS
 
-   SUBROUTINE tra_nxt( kt )
+   SUBROUTINE tra_nxt( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                   ***  ROUTINE tranxt  ***
@@ -87,4 +87,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt    ! ocean time-step index
+      INTEGER, INTENT(in) ::   Kmm, Krhs  ! time level indices
       !!
       INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
@@ -122,6 +123,6 @@
          ztrds(:,:,jpk) = 0._wp
          IF( ln_traldf_iso ) THEN              ! diagnose the "pure" Kz diffusive trend 
-            CALL trd_tra( kt, 'TRA', jp_tem, jptra_zdfp, ztrdt )
-            CALL trd_tra( kt, 'TRA', jp_sal, jptra_zdfp, ztrds )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_zdfp, ztrdt )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_zdfp, ztrds )
          ENDIF
          ! total trend for the non-time-filtered variables. 
@@ -132,6 +133,6 @@
             ztrds(:,:,jk) = ( tsa(:,:,jk,jp_sal)*e3t_a(:,:,jk) / e3t_n(:,:,jk) - tsn(:,:,jk,jp_sal)) * zfact
          END DO
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_tot, ztrdt )
-         CALL trd_tra( kt, 'TRA', jp_sal, jptra_tot, ztrds )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_tot, ztrdt )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_tot, ztrds )
          IF( ln_linssh ) THEN       ! linear sea surface height only
             ! Store now fields before applying the Asselin filter 
@@ -152,6 +153,6 @@
             ztrdt(:,:,:) = 0._wp
             ztrds(:,:,:) = 0._wp
-            CALL trd_tra( kt, 'TRA', jp_tem, jptra_atf, ztrdt )
-            CALL trd_tra( kt, 'TRA', jp_sal, jptra_atf, ztrds )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_atf, ztrdt )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_atf, ztrds )
          END IF
          !
@@ -159,5 +160,5 @@
          !
          IF( ln_linssh ) THEN   ;   CALL tra_nxt_fix( kt, nit000,      'TRA', tsb, tsn, tsa, jpts )  ! linear free surface 
-         ELSE                   ;   CALL tra_nxt_vvl( kt, nit000, rdt, 'TRA', tsb, tsn, tsa,   &
+         ELSE                   ;   CALL tra_nxt_vvl( kt, Kmm, Krhs, nit000, rdt, 'TRA', tsb, tsn, tsa,   &
            &                                                                sbc_tsc, sbc_tsc_b, jpts )  ! non-linear free surface
          ENDIF
@@ -175,6 +176,6 @@
             ztrds(:,:,jk) = ( tsb(:,:,jk,jp_sal) - ztrds(:,:,jk) ) * zfact
          END DO
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_atf, ztrdt )
-         CALL trd_tra( kt, 'TRA', jp_sal, jptra_atf, ztrds )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_atf, ztrdt )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_atf, ztrds )
       END IF
       IF( l_trdtra )   DEALLOCATE( ztrdt , ztrds )
@@ -238,5 +239,5 @@
 
 
-   SUBROUTINE tra_nxt_vvl( kt, kit000, p2dt, cdtype, ptb, ptn, pta, psbc_tc, psbc_tc_b, kjpt )
+   SUBROUTINE tra_nxt_vvl( kt, Kmm, Krhs, kit000, p2dt, cdtype, ptb, ptn, pta, psbc_tc, psbc_tc_b, kjpt )
       !!----------------------------------------------------------------------
       !!                   ***  ROUTINE tra_nxt_vvl  ***
@@ -254,4 +255,5 @@
       !!----------------------------------------------------------------------
       INTEGER                              , INTENT(in   ) ::  kt        ! ocean time-step index
+      INTEGER                              , INTENT(in   ) ::  Kmm, Krhs ! time level indices
       INTEGER                              , INTENT(in   ) ::  kit000    ! first time step index
       REAL(wp)                             , INTENT(in   ) ::  p2dt      ! time-step
@@ -367,10 +369,10 @@
       IF( ( l_trdtra .AND. cdtype == 'TRA' ) .OR. ( l_trdtrc .AND. cdtype == 'TRC' ) )   THEN
          IF( l_trdtra .AND. cdtype == 'TRA' ) THEN 
-            CALL trd_tra( kt, cdtype, jp_tem, jptra_atf, ztrd_atf(:,:,:,jp_tem) )
-            CALL trd_tra( kt, cdtype, jp_sal, jptra_atf, ztrd_atf(:,:,:,jp_sal) )
+            CALL trd_tra( kt, Kmm, Krhs, cdtype, jp_tem, jptra_atf, ztrd_atf(:,:,:,jp_tem) )
+            CALL trd_tra( kt, Kmm, Krhs, cdtype, jp_sal, jptra_atf, ztrd_atf(:,:,:,jp_sal) )
          ENDIF
          IF( l_trdtrc .AND. cdtype == 'TRC' ) THEN
             DO jn = 1, kjpt
-               CALL trd_tra( kt, cdtype, jn, jptra_atf, ztrd_atf(:,:,:,jn) )
+               CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_atf, ztrd_atf(:,:,:,jn) )
             END DO
          ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traqsr.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traqsr.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traqsr.F90	(revision 10946)
@@ -75,5 +75,5 @@
 CONTAINS
 
-   SUBROUTINE tra_qsr( kt )
+   SUBROUTINE tra_qsr( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE tra_qsr  ***
@@ -102,4 +102,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt     ! ocean time-step
+      INTEGER, INTENT(in) ::   Kmm, Krhs     ! time level indices
       !
       INTEGER  ::   ji, jj, jk               ! dummy loop indices
@@ -296,5 +297,5 @@
       IF( l_trdtra ) THEN     ! qsr tracers trends saved for diagnostics
          ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_qsr, ztrdt )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_qsr, ztrdt )
          DEALLOCATE( ztrdt ) 
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trasbc.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trasbc.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trasbc.F90	(revision 10946)
@@ -51,5 +51,5 @@
 CONTAINS
 
-   SUBROUTINE tra_sbc ( kt )
+   SUBROUTINE tra_sbc ( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE tra_sbc  ***
@@ -73,4 +73,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
+      INTEGER, INTENT(in) ::   Kmm, Krhs  ! time level indices
       !
       INTEGER  ::   ji, jj, jk, jn              ! dummy loop indices  
@@ -261,6 +262,6 @@
          ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
          ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_nsr, ztrdt )
-         CALL trd_tra( kt, 'TRA', jp_sal, jptra_nsr, ztrds )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_nsr, ztrdt )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_nsr, ztrds )
          DEALLOCATE( ztrdt , ztrds ) 
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trazdf.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trazdf.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trazdf.F90	(revision 10946)
@@ -96,6 +96,6 @@
          CALL lbc_lnk_multi( 'trazdf', ztrdt, 'T', 1. , ztrds, 'T', 1. )
 !!gm
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_zdf, ztrdt )
-         CALL trd_tra( kt, 'TRA', jp_sal, jptra_zdf, ztrds )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_zdf, ztrdt )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_zdf, ztrds )
          DEALLOCATE( ztrdt , ztrds )
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trddyn.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trddyn.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trddyn.F90	(revision 10946)
@@ -44,5 +44,5 @@
 CONTAINS
 
-   SUBROUTINE trd_dyn( putrd, pvtrd, ktrd, kt )
+   SUBROUTINE trd_dyn( putrd, pvtrd, ktrd, kt, Kmm )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE trd_mod  ***
@@ -55,4 +55,5 @@
       INTEGER                   , INTENT(in   ) ::   ktrd           ! trend index
       INTEGER                   , INTENT(in   ) ::   kt             ! time step
+      INTEGER                   , INTENT(in   ) ::   Kmm            ! time level index
       !!----------------------------------------------------------------------
       !
@@ -66,20 +67,20 @@
       !   3D output of momentum and/or tracers trends using IOM interface
       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
-      IF( ln_dyn_trd )   CALL trd_dyn_iom( putrd, pvtrd, ktrd, kt )
+      IF( ln_dyn_trd )   CALL trd_dyn_iom( putrd, pvtrd, ktrd, kt, Kmm )
          
       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
       !  Integral Constraints Properties for momentum and/or tracers trends
       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
-      IF( ln_glo_trd )   CALL trd_glo( putrd, pvtrd, ktrd, 'DYN', kt )
+      IF( ln_glo_trd )   CALL trd_glo( putrd, pvtrd, ktrd, 'DYN', kt, Kmm )
 
       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
       !  Kinetic Energy trends
       !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-      IF( ln_KE_trd  )   CALL trd_ken( putrd, pvtrd, ktrd, kt )
+      IF( ln_KE_trd  )   CALL trd_ken( putrd, pvtrd, ktrd, kt, Kmm )
 
       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
       !  Vorticity trends
       !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-      IF( ln_vor_trd )   CALL trd_vor( putrd, pvtrd, ktrd, kt )
+      IF( ln_vor_trd )   CALL trd_vor( putrd, pvtrd, ktrd, kt, Kmm )
 
       !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
@@ -91,5 +92,5 @@
 
 
-   SUBROUTINE trd_dyn_iom( putrd, pvtrd, ktrd, kt )
+   SUBROUTINE trd_dyn_iom( putrd, pvtrd, ktrd, kt, Kmm )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE trd_dyn_iom  ***
@@ -100,4 +101,5 @@
       INTEGER                   , INTENT(in   ) ::   ktrd           ! trend index
       INTEGER                   , INTENT(in   ) ::   kt             ! time step
+      INTEGER                   , INTENT(in   ) ::   Kmm            ! time level index
       !
       INTEGER ::   ji, jj, jk   ! dummy loop indices
@@ -124,6 +126,6 @@
                                  DO jj = 2, jpjm1
                                     DO ji = 2, jpim1
-                                       z3dx(ji,jj,jk) = un(ji,jj,jk) * ( un(ji+1,jj,jk) - un(ji-1,jj,jk) ) / ( 2._wp * e1u(ji,jj) )
-                                       z3dy(ji,jj,jk) = vn(ji,jj,jk) * ( vn(ji,jj+1,jk) - vn(ji,jj-1,jk) ) / ( 2._wp * e2v(ji,jj) )
+                                       z3dx(ji,jj,jk) = uu(ji,jj,jk,Kmm) * ( uu(ji+1,jj,jk,Kmm) - uu(ji-1,jj,jk,Kmm) ) / ( 2._wp * e1u(ji,jj) )
+                                       z3dy(ji,jj,jk) = vv(ji,jj,jk,Kmm) * ( vv(ji,jj+1,jk,Kmm) - vv(ji,jj-1,jk,Kmm) ) / ( 2._wp * e2v(ji,jj) )
                                     END DO
                                  END DO
@@ -142,6 +144,6 @@
                               !                                    ! wind stress trends
                               ALLOCATE( z2dx(jpi,jpj) , z2dy(jpi,jpj) )
-                              z2dx(:,:) = ( utau_b(:,:) + utau(:,:) ) / ( e3u_n(:,:,1) * rau0 )
-                              z2dy(:,:) = ( vtau_b(:,:) + vtau(:,:) ) / ( e3v_n(:,:,1) * rau0 )
+                              z2dx(:,:) = ( utau_b(:,:) + utau(:,:) ) / ( e3u(:,:,1,Kmm) * rau0 )
+                              z2dy(:,:) = ( vtau_b(:,:) + vtau(:,:) ) / ( e3v(:,:,1,Kmm) * rau0 )
                               CALL iom_put( "utrd_tau", z2dx )
                               CALL iom_put( "vtrd_tau", z2dy )
@@ -159,7 +161,7 @@
 !                                          ikbv = mbkv(ji,jj)
 !                                          z3dx(ji,jj,jk) = 0.5 * ( rCdU_bot(ji+1,jj) + rCdU_bot(ji,jj) ) & 
-!                                               &         * un(ji,jj,ikbu) / e3u_n(ji,jj,ikbu)
+!                                               &         * uu(ji,jj,ikbu,Kmm) / e3u(ji,jj,ikbu,Kmm)
 !                                          z3dy(ji,jj,jk) = 0.5 * ( rCdU_bot(ji,jj+1) + rCdU_bot(ji,jj) ) &
-!                                               &         * vn(ji,jj,ikbv) / e3v_n(ji,jj,ikbv)
+!                                               &         * vv(ji,jj,ikbv,Kmm) / e3v(ji,jj,ikbv,Kmm)
 !                                    END DO
 !                                 END DO
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdglo.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdglo.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdglo.F90	(revision 10946)
@@ -59,5 +59,5 @@
 CONTAINS
 
-   SUBROUTINE trd_glo( ptrdx, ptrdy, ktrd, ctype, kt )
+   SUBROUTINE trd_glo( ptrdx, ptrdy, ktrd, ctype, kt, Kmm )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE trd_glo  ***
@@ -72,4 +72,5 @@
       CHARACTER(len=3)          , INTENT(in   ) ::   ctype   ! momentum or tracers trends type (='DYN'/'TRA')
       INTEGER                   , INTENT(in   ) ::   kt      ! time step
+      INTEGER                   , INTENT(in   ) ::   Kmm     ! time level index
       !!
       INTEGER ::   ji, jj, jk      ! dummy loop indices
@@ -87,11 +88,11 @@
       	      DO jj = 1, jpj
       	         DO ji = 1, jpi        
-      	            zvm = e1e2t(ji,jj) * e3t_n(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj)
+      	            zvm = e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) * tmask_i(ji,jj)
       	            zvt = ptrdx(ji,jj,jk) * zvm
       	            zvs = ptrdy(ji,jj,jk) * zvm
          	         tmo(ktrd) = tmo(ktrd) + zvt   
          	         smo(ktrd) = smo(ktrd) + zvs
-         	         t2 (ktrd) = t2(ktrd)  + zvt * tsn(ji,jj,jk,jp_tem)
-         	         s2 (ktrd) = s2(ktrd)  + zvs * tsn(ji,jj,jk,jp_sal)
+         	         t2 (ktrd) = t2(ktrd)  + zvt * ts(ji,jj,jk,jp_tem,Kmm)
+         	         s2 (ktrd) = s2(ktrd)  + zvs * ts(ji,jj,jk,jp_sal,Kmm)
                   END DO
                END DO
@@ -99,8 +100,8 @@
             !                       ! linear free surface: diagnose advective flux trough the fixed k=1 w-surface
             IF( ln_linssh .AND. ktrd == jptra_zad ) THEN  
-               tmo(jptra_sad) = SUM( wn(:,:,1) * tsn(:,:,1,jp_tem) * e1e2t(:,:) * tmask_i(:,:) )
-               smo(jptra_sad) = SUM( wn(:,:,1) * tsn(:,:,1,jp_sal) * e1e2t(:,:) * tmask_i(:,:)  )
-               t2 (jptra_sad) = SUM( wn(:,:,1) * tsn(:,:,1,jp_tem) * tsn(:,:,1,jp_tem) * e1e2t(:,:) * tmask_i(:,:)  )
-               s2 (jptra_sad) = SUM( wn(:,:,1) * tsn(:,:,1,jp_sal) * tsn(:,:,1,jp_sal) * e1e2t(:,:) * tmask_i(:,:)  )
+               tmo(jptra_sad) = SUM( ww(:,:,1) * ts(:,:,1,jp_tem,Kmm) * e1e2t(:,:) * tmask_i(:,:) )
+               smo(jptra_sad) = SUM( ww(:,:,1) * ts(:,:,1,jp_sal,Kmm) * e1e2t(:,:) * tmask_i(:,:)  )
+               t2 (jptra_sad) = SUM( ww(:,:,1) * ts(:,:,1,jp_tem,Kmm) * ts(:,:,1,jp_tem,Kmm) * e1e2t(:,:) * tmask_i(:,:)  )
+               s2 (jptra_sad) = SUM( ww(:,:,1) * ts(:,:,1,jp_sal,Kmm) * ts(:,:,1,jp_sal,Kmm) * e1e2t(:,:) * tmask_i(:,:)  )
             ENDIF
             !
@@ -121,10 +122,10 @@
                   DO ji = 1, jpim1
                      zvt = ptrdx(ji,jj,jk) * tmask_i(ji+1,jj) * tmask_i(ji,jj) * umask(ji,jj,jk)   &
-                        &                                     * e1e2u  (ji,jj) * e3u_n(ji,jj,jk)
+                        &                                     * e1e2u  (ji,jj) * e3u(ji,jj,jk,Kmm)
                      zvs = ptrdy(ji,jj,jk) * tmask_i(ji,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)   &
-                        &                                     * e1e2v  (ji,jj) * e3u_n(ji,jj,jk)
+                        &                                     * e1e2v  (ji,jj) * e3u(ji,jj,jk,Kmm)
                      umo(ktrd) = umo(ktrd) + zvt
                      vmo(ktrd) = vmo(ktrd) + zvs
-                     hke(ktrd) = hke(ktrd) + un(ji,jj,jk) * zvt + vn(ji,jj,jk) * zvs
+                     hke(ktrd) = hke(ktrd) + uu(ji,jj,jk,Kmm) * zvt + vv(ji,jj,jk,Kmm) * zvs
                   END DO
                END DO
@@ -141,5 +142,5 @@
                      umo(jpdyn_tau) = umo(jpdyn_tau) + zvt
                      vmo(jpdyn_tau) = vmo(jpdyn_tau) + zvs
-                     hke(jpdyn_tau) = hke(jpdyn_tau) + un(ji,jj,1) * zvt + vn(ji,jj,1) * zvs
+                     hke(jpdyn_tau) = hke(jpdyn_tau) + uu(ji,jj,1,Kmm) * zvt + vv(ji,jj,1,Kmm) * zvs
                   END DO
                END DO
@@ -155,9 +156,9 @@
 !                        ikbu = mbku(ji,jj)                  ! deepest ocean u- & v-levels
 !                        ikbv = mbkv(ji,jj)
-!                        zvt = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * un(ji,jj,ikbu) * e1e2u(ji,jj)
-!                        zvs = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * vn(ji,jj,ikbv) * e1e2v(ji,jj)
+!                        zvt = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * uu(ji,jj,ikbu,Kmm) * e1e2u(ji,jj)
+!                        zvs = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * vv(ji,jj,ikbv,Kmm) * e1e2v(ji,jj)
 !                        umo(jpdyn_bfri) = umo(jpdyn_bfri) + zvt
 !                        vmo(jpdyn_bfri) = vmo(jpdyn_bfri) + zvs
-!                        hke(jpdyn_bfri) = hke(jpdyn_bfri) + un(ji,jj,ikbu) * zvt + vn(ji,jj,ikbv) * zvs
+!                        hke(jpdyn_bfri) = hke(jpdyn_bfri) + uu(ji,jj,ikbu,Kmm) * zvt + vv(ji,jj,ikbv,Kmm) * zvs
 !                     END DO
 !                  END DO
@@ -183,5 +184,5 @@
 
 
-   SUBROUTINE glo_dyn_wri( kt )
+   SUBROUTINE glo_dyn_wri( kt, Kmm )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE glo_dyn_wri  ***
@@ -190,4 +191,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
+      INTEGER, INTENT(in) ::   Kmm  ! time level index
       !
       INTEGER  ::   ji, jj, jk   ! dummy loop indices
@@ -209,10 +211,10 @@
          zkepe(:,:,:) = 0._wp
    
-         CALL eos( tsn, rhd, rhop )       ! now potential density
+         CALL eos( ts(:,:,:,:,Kmm), rhd, rhop )       ! now potential density
 
          zcof = 0.5_wp / rau0             ! Density flux at w-point
          zkz(:,:,1) = 0._wp
          DO jk = 2, jpk
-            zkz(:,:,jk) = zcof * e1e2t(:,:) * wn(:,:,jk) * ( rhop(:,:,jk) + rhop(:,:,jk-1) ) * tmask_i(:,:)
+            zkz(:,:,jk) = zcof * e1e2t(:,:) * ww(:,:,jk) * ( rhop(:,:,jk) + rhop(:,:,jk-1) ) * tmask_i(:,:)
          END DO
          
@@ -221,6 +223,6 @@
             DO jj = 1, jpjm1
                DO ji = 1, jpim1
-                  zkx(ji,jj,jk) = zcof * e2u(ji,jj) * e3u_n(ji,jj,jk) * un(ji,jj,jk) * ( rhop(ji,jj,jk) + rhop(ji+1,jj,jk) )
-                  zky(ji,jj,jk) = zcof * e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk) * ( rhop(ji,jj,jk) + rhop(ji,jj+1,jk) )
+                  zkx(ji,jj,jk) = zcof * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * uu(ji,jj,jk,Kmm) * ( rhop(ji,jj,jk) + rhop(ji+1,jj,jk) )
+                  zky(ji,jj,jk) = zcof * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) * vv(ji,jj,jk,Kmm) * ( rhop(ji,jj,jk) + rhop(ji,jj+1,jk) )
                END DO
             END DO
@@ -233,5 +235,5 @@
                      &                 + zkx(ji,jj,jk) - zkx(ji-1,jj  ,jk  )               &
                      &                 + zky(ji,jj,jk) - zky(ji  ,jj-1,jk  )   )           &
-                     &              / ( e1e2t(ji,jj) * e3t_n(ji,jj,jk) ) * tmask(ji,jj,jk) * tmask_i(ji,jj)
+                     &              / ( e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) ) * tmask(ji,jj,jk) * tmask_i(ji,jj)
                END DO
             END DO
@@ -242,5 +244,5 @@
          peke = 0._wp
          DO jk = 1, jpkm1
-            peke = peke + SUM( zkepe(:,:,jk) * gdept_n(:,:,jk) * e1e2t(:,:) * e3t_n(:,:,jk) )
+            peke = peke + SUM( zkepe(:,:,jk) * gdept(:,:,jk,Kmm) * e1e2t(:,:) * e3t(:,:,jk,Kmm) )
          END DO
          peke = grav * peke
@@ -506,5 +508,5 @@
 
 
-   SUBROUTINE trd_glo_init
+   SUBROUTINE trd_glo_init( Kmm )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE trd_glo_init  ***
@@ -512,4 +514,5 @@
       !! ** Purpose :   Read the namtrd namelist
       !!----------------------------------------------------------------------
+      INTEGER, INTENT(in) ::   Kmm   ! time level index
       INTEGER  ::   ji, jj, jk   ! dummy loop indices
       !!----------------------------------------------------------------------
@@ -524,5 +527,5 @@
       tvolt = 0._wp
       DO jk = 1, jpkm1
-         tvolt = tvolt + SUM( e1e2t(:,:) * e3t_n(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:) )
+         tvolt = tvolt + SUM( e1e2t(:,:) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) * tmask_i(:,:) )
       END DO
       CALL mpp_sum( 'trdglo', tvolt )   ! sum over the global domain
@@ -541,6 +544,6 @@
          DO jj = 2, jpjm1
             DO ji = fs_2, fs_jpim1   ! vector opt.
-               tvolu = tvolu + e1u(ji,jj) * e2u(ji,jj) * e3u_n(ji,jj,jk) * tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk)
-               tvolv = tvolv + e1v(ji,jj) * e2v(ji,jj) * e3v_n(ji,jj,jk) * tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)
+               tvolu = tvolu + e1u(ji,jj) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk)
+               tvolv = tvolv + e1v(ji,jj) * e2v(ji,jj) * e3v(ji,jj,jk,Kmm) * tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)
             END DO
          END DO
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdini.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdini.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdini.F90	(revision 10946)
@@ -34,5 +34,5 @@
 CONTAINS
 
-   SUBROUTINE trd_init
+   SUBROUTINE trd_init( Kmm )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE trd_init  ***
@@ -40,4 +40,5 @@
       !! ** Purpose :   Initialization of trend diagnostics
       !!----------------------------------------------------------------------
+      INTEGER, INTENT(in) ::   Kmm  ! time level index
       INTEGER ::   ios   ! local integer
       !!
@@ -96,5 +97,5 @@
 
       !                             ! diagnostic initialization  
-      IF( ln_glo_trd )   CALL trd_glo_init      ! global domain averaged trends
+      IF( ln_glo_trd )   CALL trd_glo_init( Kmm )      ! global domain averaged trends
       IF( ln_tra_mxl )   CALL trd_mxl_init      ! mixed-layer          trends  
       IF( ln_vor_trd )   CALL trd_vor_init      ! barotropic vorticity trends
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdken.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdken.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdken.F90	(revision 10946)
@@ -59,5 +59,5 @@
 
 
-   SUBROUTINE trd_ken( putrd, pvtrd, ktrd, kt )
+   SUBROUTINE trd_ken( putrd, pvtrd, ktrd, kt, Kmm )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE trd_ken  ***
@@ -67,5 +67,5 @@
       !! ** Method  : - apply lbc to the input masked velocity trends 
       !!              - compute the associated KE trend:
-      !!          zke = 0.5 * (  mi-1[ un * putrd * bu ] + mj-1[ vn * pvtrd * bv]  ) / bt
+      !!          zke = 0.5 * (  mi-1[ uu(Kmm) * putrd * bu ] + mj-1[ vv(Kmm) * pvtrd * bv]  ) / bt
       !!      where bu, bv, bt are the volume of u-, v- and t-boxes. 
       !!              - vertical diffusion case (jpdyn_zdf): 
@@ -80,4 +80,5 @@
       INTEGER                   , INTENT(in   ) ::   ktrd           ! trend index
       INTEGER                   , INTENT(in   ) ::   kt             ! time step
+      INTEGER                   , INTENT(in   ) ::   Kmm            ! time level index
       !
       INTEGER ::   ji, jj, jk       ! dummy loop indices
@@ -92,7 +93,7 @@
       nkstp = kt
       DO jk = 1, jpkm1
-         bu   (:,:,jk) =    e1e2u(:,:) * e3u_n(:,:,jk)
-         bv   (:,:,jk) =    e1e2v(:,:) * e3v_n(:,:,jk)
-         r1_bt(:,:,jk) = r1_e1e2t(:,:) / e3t_n(:,:,jk) * tmask(:,:,jk)
+         bu   (:,:,jk) =    e1e2u(:,:) * e3u(:,:,jk,Kmm)
+         bv   (:,:,jk) =    e1e2v(:,:) * e3v(:,:,jk,Kmm)
+         r1_bt(:,:,jk) = r1_e1e2t(:,:) / e3t(:,:,jk,Kmm) * tmask(:,:,jk)
       END DO
       !
@@ -103,8 +104,8 @@
          DO jj = 2, jpj
             DO ji = 2, jpi
-               zke(ji,jj,jk) = 0.5_wp * rau0 *( un(ji  ,jj,jk) * putrd(ji  ,jj,jk) * bu(ji  ,jj,jk)  &
-                  &                           + un(ji-1,jj,jk) * putrd(ji-1,jj,jk) * bu(ji-1,jj,jk)  &
-                  &                           + vn(ji,jj  ,jk) * pvtrd(ji,jj  ,jk) * bv(ji,jj  ,jk)  &
-                  &                           + vn(ji,jj-1,jk) * pvtrd(ji,jj-1,jk) * bv(ji,jj-1,jk)  ) * r1_bt(ji,jj,jk)
+               zke(ji,jj,jk) = 0.5_wp * rau0 *( uu(ji  ,jj,jk,Kmm) * putrd(ji  ,jj,jk) * bu(ji  ,jj,jk)  &
+                  &                           + uu(ji-1,jj,jk,Kmm) * putrd(ji-1,jj,jk) * bu(ji-1,jj,jk)  &
+                  &                           + vv(ji,jj  ,jk,Kmm) * pvtrd(ji,jj  ,jk) * bv(ji,jj  ,jk)  &
+                  &                           + vv(ji,jj-1,jk,Kmm) * pvtrd(ji,jj-1,jk) * bv(ji,jj-1,jk)  ) * r1_bt(ji,jj,jk)
             END DO
          END DO
@@ -122,6 +123,6 @@
          !                   !                                          ! wind stress trends
                                  ALLOCATE( z2dx(jpi,jpj) , z2dy(jpi,jpj) , zke2d(jpi,jpj) )
-                           z2dx(:,:) = un(:,:,1) * ( utau_b(:,:) + utau(:,:) ) * e1e2u(:,:) * umask(:,:,1)
-                           z2dy(:,:) = vn(:,:,1) * ( vtau_b(:,:) + vtau(:,:) ) * e1e2v(:,:) * vmask(:,:,1)
+                           z2dx(:,:) = uu(:,:,1,Kmm) * ( utau_b(:,:) + utau(:,:) ) * e1e2u(:,:) * umask(:,:,1)
+                           z2dy(:,:) = vv(:,:,1,Kmm) * ( vtau_b(:,:) + vtau(:,:) ) * e1e2v(:,:) * vmask(:,:,1)
                            zke2d(1,:) = 0._wp   ;   zke2d(:,1) = 0._wp
                            DO jj = 2, jpj
@@ -141,6 +142,6 @@
 !                  ikbu = mbku(ji,jj)         ! deepest ocean u- & v-levels
 !                  ikbv = mbkv(ji,jj)   
-!                  z2dx(ji,jj) = un(ji,jj,ikbu) * bfrua(ji,jj) * un(ji,jj,ikbu)
-!                  z2dy(ji,jj) = vn(ji,jj,ikbu) * bfrva(ji,jj) * vn(ji,jj,ikbv)
+!                  z2dx(ji,jj) = uu(ji,jj,ikbu,Kmm) * bfrua(ji,jj) * uu(ji,jj,ikbu,Kmm)
+!                  z2dy(ji,jj) = vv(ji,jj,ikbu,Kmm) * bfrva(ji,jj) * vv(ji,jj,ikbv,Kmm)
 !               END DO
 !            END DO
@@ -157,5 +158,5 @@
          CASE( jpdyn_atf )   ;   CALL iom_put( "ketrd_atf"   , zke )    ! asselin filter trends 
 !! a faire !!!!  idee changer dynnxt pour avoir un appel a jpdyn_bfr avant le swap !!!
-!! reflechir a une possible sauvegarde du "vrai" un,vn pour le calcul de atf....
+!! reflechir a une possible sauvegarde du "vrai" uu(Kmm),vv(Kmm) pour le calcul de atf....
 !
 !         IF( ln_drgimp ) THEN                                          ! bottom friction (implicit case)
@@ -164,6 +165,6 @@
 !                  ikbu = mbku(ji,jj)          ! deepest ocean u- & v-levels
 !                  ikbv = mbkv(ji,jj)
-!                  z2dx(ji,jj) = un(ji,jj,ikbu) * bfrua(ji,jj) * un(ji,jj,ikbu) / e3u_n(ji,jj,ikbu)
-!                  z2dy(ji,jj) = un(ji,jj,ikbu) * bfrva(ji,jj) * vn(ji,jj,ikbv) / e3v_n(ji,jj,ikbv)
+!                  z2dx(ji,jj) = uu(ji,jj,ikbu,Kmm) * bfrua(ji,jj) * uu(ji,jj,ikbu,Kmm) / e3u(ji,jj,ikbu,Kmm)
+!                  z2dy(ji,jj) = uu(ji,jj,ikbu,Kmm) * bfrva(ji,jj) * vv(ji,jj,ikbv,Kmm) / e3v(ji,jj,ikbv,Kmm)
 !               END DO
 !            END DO
@@ -179,9 +180,9 @@
         CASE( jpdyn_ken )   ;   ! kinetic energy
                     ! called in dynnxt.F90 before asselin time filter
-                    ! with putrd=ua and pvtrd=va
+                    ! with putrd=uu(Krhs) and pvtrd=vv(Krhs)
                     zke(:,:,:) = 0.5_wp * zke(:,:,:)
                     CALL iom_put( "KE", zke )
                     !
-                    CALL ken_p2k( kt , zke )
+                    CALL ken_p2k( kt , zke, Kmm )
                       CALL iom_put( "ketrd_convP2K", zke )     ! conversion -rau*g*w
          !
@@ -191,5 +192,5 @@
 
 
-   SUBROUTINE ken_p2k( kt , pconv )
+   SUBROUTINE ken_p2k( kt , pconv, Kmm )
       !!---------------------------------------------------------------------
       !!                 ***  ROUTINE ken_p2k  ***
@@ -202,4 +203,5 @@
       !!---------------------------------------------------------------------- 
       INTEGER                   , INTENT(in   ) ::   kt      ! ocean time-step index
+      INTEGER                   , INTENT(in   ) ::   Kmm     ! time level index
       REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pconv   ! 
       !
@@ -214,9 +216,9 @@
       
       !  Surface value (also valid in partial step case)
-      zconv(:,:,1) = zcoef * ( 2._wp * rhd(:,:,1) ) * wn(:,:,1) * e3w_n(:,:,1)
+      zconv(:,:,1) = zcoef * ( 2._wp * rhd(:,:,1) ) * ww(:,:,1) * e3w(:,:,1,Kmm)
 
       ! interior value (2=<jk=<jpkm1)
       DO jk = 2, jpk
-         zconv(:,:,jk) = zcoef * ( rhd(:,:,jk) + rhd(:,:,jk-1) ) * wn(:,:,jk) * e3w_n(:,:,jk)
+         zconv(:,:,jk) = zcoef * ( rhd(:,:,jk) + rhd(:,:,jk-1) ) * ww(:,:,jk) * e3w(:,:,jk,Kmm)
       END DO
 
@@ -225,5 +227,5 @@
          DO jj = 1, jpj
             DO ji = 1, jpi
-               zcoef = 0.5_wp / e3t_n(ji,jj,jk)
+               zcoef = 0.5_wp / e3t(ji,jj,jk,Kmm)
                pconv(ji,jj,jk) = zcoef * ( zconv(ji,jj,jk) + zconv(ji,jj,jk+1) ) * tmask(ji,jj,jk)
             END DO
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdmxl.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdmxl.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdmxl.F90	(revision 10946)
@@ -86,5 +86,5 @@
 
 
-   SUBROUTINE trd_tra_mxl( ptrdx, ptrdy, ktrd, kt, p2dt, kmxln )
+   SUBROUTINE trd_tra_mxl( ptrdx, ptrdy, ktrd, kt, p2dt, kmxln, Kmm )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE trd_tra_mng  ***
@@ -98,4 +98,5 @@
       INTEGER                   , INTENT(in   ) ::   ktrd    ! tracer trend index
       INTEGER                   , INTENT(in   ) ::   kt      ! time step index
+      INTEGER                   , INTENT(in   ) ::   Kmm     ! time level index
       REAL(wp)                  , INTENT(in   ) ::   p2dt    ! time step  [s]
       REAL(wp), DIMENSION(:,:)  , INTENT(in   ) ::   kmxln   ! number of t-box for the vertical average 
@@ -119,5 +120,5 @@
             DO jj = 1,jpj
                DO ji = 1,jpi
-                  IF( jk - kmxln(ji,jj) < 0 )   wkx(ji,jj,jk) = e3t_n(ji,jj,jk) * tmask(ji,jj,jk)
+                  IF( jk - kmxln(ji,jj) < 0 )   wkx(ji,jj,jk) = e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk)
                END DO
             END DO
@@ -136,6 +137,6 @@
          tml(:,:) = 0._wp   ;   sml(:,:) = 0._wp
          DO jk = 1, jpktrd
-            tml(:,:) = tml(:,:) + wkx(:,:,jk) * tsn(:,:,jk,jp_tem)
-            sml(:,:) = sml(:,:) + wkx(:,:,jk) * tsn(:,:,jk,jp_sal)
+            tml(:,:) = tml(:,:) + wkx(:,:,jk) * ts(:,:,jk,jp_tem,Kmm)
+            sml(:,:) = sml(:,:) + wkx(:,:,jk) * ts(:,:,jk,jp_sal,Kmm)
          END DO
          !
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdpen.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdpen.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdpen.F90	(revision 10946)
@@ -55,5 +55,5 @@
 
 
-   SUBROUTINE trd_pen( ptrdx, ptrdy, ktrd, kt, pdt )
+   SUBROUTINE trd_pen( ptrdx, ptrdy, ktrd, kt, pdt, Kmm )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE trd_tra_mng  ***
@@ -66,4 +66,5 @@
       INTEGER                   , INTENT(in) ::   ktrd           ! tracer trend index
       INTEGER                   , INTENT(in) ::   kt             ! time step index
+      INTEGER                   , INTENT(in) ::   Kmm            ! time level index
       REAL(wp)                  , INTENT(in) ::   pdt            ! time step [s]
       !
@@ -77,5 +78,5 @@
       IF( kt /= nkstp ) THEN     ! full eos: set partial derivatives at the 1st call of kt time step
          nkstp = kt
-         CALL eos_pen( tsn, rab_PE, zpe )
+         CALL eos_pen( ts(:,:,:,:,Kmm), rab_PE, zpe )
          CALL iom_put( "alphaPE", rab_pe(:,:,:,jp_tem) )
          CALL iom_put( "betaPE" , rab_pe(:,:,:,jp_sal) )
@@ -95,8 +96,8 @@
                                 IF( ln_linssh ) THEN                   ! cst volume : adv flux through z=0 surface
                                    ALLOCATE( z2d(jpi,jpj) )
-                                   z2d(:,:) = wn(:,:,1) * ( &
-                                     &   - ( rab_n(:,:,1,jp_tem) + rab_pe(:,:,1,jp_tem) ) * tsn(:,:,1,jp_tem)    &
-                                     &   + ( rab_n(:,:,1,jp_sal) + rab_pe(:,:,1,jp_sal) ) * tsn(:,:,1,jp_sal)    &
-                                     & ) / e3t_n(:,:,1)
+                                   z2d(:,:) = ww(:,:,1) * ( &
+                                     &   - ( rab_n(:,:,1,jp_tem) + rab_pe(:,:,1,jp_tem) ) * ts(:,:,1,jp_tem,Kmm)    &
+                                     &   + ( rab_n(:,:,1,jp_sal) + rab_pe(:,:,1,jp_sal) ) * ts(:,:,1,jp_sal,Kmm)    &
+                                     & ) / e3t(:,:,1,Kmm)
                                    CALL iom_put( "petrd_sad" , z2d )
                                    DEALLOCATE( z2d )
@@ -112,12 +113,4 @@
       CASE ( jptra_bbc  )   ;   CALL iom_put( "petrd_bbc" , zpe )   ! bottom bound cond (geoth flux)
       CASE ( jptra_atf  )   ;   CALL iom_put( "petrd_atf" , zpe )   ! asselin time filter (last trend)
-                                !IF( ln_linssh ) THEN                   ! cst volume : ssh term (otherwise include in e3t variation)
-                                !   ALLOCATE( z2d(jpi,jpj) )
-                                !   z2d(:,:) = ( ssha(:,:) - sshb(:,:) )                 &
-                                !      &     * (   dPE_dt(:,:,1) * tsn(:,:,1,jp_tem)    &
-                                !      &         + dPE_ds(:,:,1) * tsn(:,:,1,jp_sal)  ) / ( e3t_n(:,:,1) * pdt )
-                                !   CALL iom_put( "petrd_sad" , z2d )
-                                !   DEALLOCATE( z2d )
-                                !ENDIF
          !
       END SELECT
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdtra.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdtra.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdtra.F90	(revision 10946)
@@ -60,5 +60,5 @@
 
 
-   SUBROUTINE trd_tra( kt, ctype, ktra, ktrd, ptrd, pun, ptra )
+   SUBROUTINE trd_tra( kt, Kmm, Krhs, ctype, ktra, ktrd, ptrd, pu, ptra )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE trd_tra  ***
@@ -77,6 +77,7 @@
       INTEGER                         , INTENT(in)           ::   ktra    ! tracer index
       INTEGER                         , INTENT(in)           ::   ktrd    ! tracer trend index
+      INTEGER                         , INTENT(in)           ::   Kmm, Krhs ! time level indices
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)           ::   ptrd    ! tracer trend  or flux
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::   pun     ! now velocity 
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::   pu      ! now velocity 
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::   ptra    ! now tracer variable
       !
@@ -94,11 +95,11 @@
          SELECT CASE( ktrd )
          !                            ! advection: transform the advective flux into a trend
-         CASE( jptra_xad )   ;   CALL trd_tra_adv( ptrd, pun, ptra, 'X', trdtx ) 
-         CASE( jptra_yad )   ;   CALL trd_tra_adv( ptrd, pun, ptra, 'Y', trdty ) 
-         CASE( jptra_zad )   ;   CALL trd_tra_adv( ptrd, pun, ptra, 'Z', trdt  ) 
+         CASE( jptra_xad )   ;   CALL trd_tra_adv( ptrd, pu, ptra, 'X', trdtx, Kmm ) 
+         CASE( jptra_yad )   ;   CALL trd_tra_adv( ptrd, pu, ptra, 'Y', trdty, Kmm ) 
+         CASE( jptra_zad )   ;   CALL trd_tra_adv( ptrd, pu, ptra, 'Z', trdt, Kmm )
          CASE( jptra_bbc,    &        ! qsr, bbc: on temperature only, send to trd_tra_mng
             &  jptra_qsr )   ;   trdt(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)
                                  ztrds(:,:,:) = 0._wp
-                                 CALL trd_tra_mng( trdt, ztrds, ktrd, kt )
+                                 CALL trd_tra_mng( trdt, ztrds, ktrd, kt, Kmm )
  !!gm Gurvan, verify the jptra_evd trend please !
          CASE( jptra_evd )   ;   avt_evd(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)
@@ -114,10 +115,10 @@
          !                            ! advection: transform the advective flux into a trend
          !                            !            and send T & S trends to trd_tra_mng
-         CASE( jptra_xad  )   ;   CALL trd_tra_adv( ptrd , pun  , ptra, 'X'  , ztrds ) 
-                                  CALL trd_tra_mng( trdtx, ztrds, ktrd, kt   )
-         CASE( jptra_yad  )   ;   CALL trd_tra_adv( ptrd , pun  , ptra, 'Y'  , ztrds ) 
-                                  CALL trd_tra_mng( trdty, ztrds, ktrd, kt   )
-         CASE( jptra_zad  )   ;   CALL trd_tra_adv( ptrd , pun  , ptra, 'Z'  , ztrds ) 
-                                  CALL trd_tra_mng( trdt , ztrds, ktrd, kt   )
+         CASE( jptra_xad  )   ;   CALL trd_tra_adv( ptrd , pu  , ptra, 'X'  , ztrds, Kmm ) 
+                                  CALL trd_tra_mng( trdtx, ztrds, ktrd, kt, Kmm   )
+         CASE( jptra_yad  )   ;   CALL trd_tra_adv( ptrd , pu  , ptra, 'Y'  , ztrds, Kmm ) 
+                                  CALL trd_tra_mng( trdty, ztrds, ktrd, kt, Kmm   )
+         CASE( jptra_zad  )   ;   CALL trd_tra_adv( ptrd , pu  , ptra, 'Z'  , ztrds, Kmm ) 
+                                  CALL trd_tra_mng( trdt , ztrds, ktrd, kt, Kmm   )
          CASE( jptra_zdfp )           ! diagnose the "PURE" Kz trend (here: just before the swap)
             !                         ! iso-neutral diffusion case otherwise jptra_zdf is "PURE"
@@ -127,28 +128,28 @@
             zwt(:,:,jpk) = 0._wp   ;   zws(:,:,jpk) = 0._wp
             DO jk = 2, jpk
-               zwt(:,:,jk) = avt(:,:,jk) * ( tsa(:,:,jk-1,jp_tem) - tsa(:,:,jk,jp_tem) ) / e3w_n(:,:,jk) * tmask(:,:,jk)
-               zws(:,:,jk) = avs(:,:,jk) * ( tsa(:,:,jk-1,jp_sal) - tsa(:,:,jk,jp_sal) ) / e3w_n(:,:,jk) * tmask(:,:,jk)
+               zwt(:,:,jk) = avt(:,:,jk) * ( ts(:,:,jk-1,jp_tem,Krhs) - ts(:,:,jk,jp_tem,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk)
+               zws(:,:,jk) = avs(:,:,jk) * ( ts(:,:,jk-1,jp_sal,Krhs) - ts(:,:,jk,jp_sal,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk)
             END DO
             !
             ztrdt(:,:,jpk) = 0._wp   ;   ztrds(:,:,jpk) = 0._wp
             DO jk = 1, jpkm1
-               ztrdt(:,:,jk) = ( zwt(:,:,jk) - zwt(:,:,jk+1) ) / e3t_n(:,:,jk)
-               ztrds(:,:,jk) = ( zws(:,:,jk) - zws(:,:,jk+1) ) / e3t_n(:,:,jk) 
-            END DO
-            CALL trd_tra_mng( ztrdt, ztrds, jptra_zdfp, kt )  
+               ztrdt(:,:,jk) = ( zwt(:,:,jk) - zwt(:,:,jk+1) ) / e3t(:,:,jk,Kmm)
+               ztrds(:,:,jk) = ( zws(:,:,jk) - zws(:,:,jk+1) ) / e3t(:,:,jk,Kmm) 
+            END DO
+            CALL trd_tra_mng( ztrdt, ztrds, jptra_zdfp, kt, Kmm )  
             !
             !                         ! Also calculate EVD trend at this point. 
             zwt(:,:,:) = 0._wp   ;   zws(:,:,:) = 0._wp            ! vertical diffusive fluxes
             DO jk = 2, jpk
-               zwt(:,:,jk) = avt_evd(:,:,jk) * ( tsa(:,:,jk-1,jp_tem) - tsa(:,:,jk,jp_tem) ) / e3w_n(:,:,jk) * tmask(:,:,jk)
-               zws(:,:,jk) = avt_evd(:,:,jk) * ( tsa(:,:,jk-1,jp_sal) - tsa(:,:,jk,jp_sal) ) / e3w_n(:,:,jk) * tmask(:,:,jk)
+               zwt(:,:,jk) = avt_evd(:,:,jk) * ( ts(:,:,jk-1,jp_tem,Krhs) - ts(:,:,jk,jp_tem,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk)
+               zws(:,:,jk) = avt_evd(:,:,jk) * ( ts(:,:,jk-1,jp_sal,Krhs) - ts(:,:,jk,jp_sal,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk)
             END DO
             !
             ztrdt(:,:,jpk) = 0._wp   ;   ztrds(:,:,jpk) = 0._wp
             DO jk = 1, jpkm1
-               ztrdt(:,:,jk) = ( zwt(:,:,jk) - zwt(:,:,jk+1) ) / e3t_n(:,:,jk)
-               ztrds(:,:,jk) = ( zws(:,:,jk) - zws(:,:,jk+1) ) / e3t_n(:,:,jk) 
-            END DO
-            CALL trd_tra_mng( ztrdt, ztrds, jptra_evd, kt )  
+               ztrdt(:,:,jk) = ( zwt(:,:,jk) - zwt(:,:,jk+1) ) / e3t(:,:,jk,Kmm)
+               ztrds(:,:,jk) = ( zws(:,:,jk) - zws(:,:,jk+1) ) / e3t(:,:,jk,Kmm) 
+            END DO
+            CALL trd_tra_mng( ztrdt, ztrds, jptra_evd, kt, Kmm )  
             !
             DEALLOCATE( zwt, zws, ztrdt )
@@ -156,5 +157,5 @@
          CASE DEFAULT                 ! other trends: mask and send T & S trends to trd_tra_mng
             ztrds(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)
-            CALL trd_tra_mng( trdt, ztrds, ktrd, kt )  
+            CALL trd_tra_mng( trdt, ztrds, ktrd, kt, Kmm )  
          END SELECT
       ENDIF
@@ -164,7 +165,7 @@
          SELECT CASE( ktrd )
          !                            ! advection: transform the advective flux into a masked trend
-         CASE( jptra_xad )   ;   CALL trd_tra_adv( ptrd , pun , ptra, 'X', ztrds ) 
-         CASE( jptra_yad )   ;   CALL trd_tra_adv( ptrd , pun , ptra, 'Y', ztrds ) 
-         CASE( jptra_zad )   ;   CALL trd_tra_adv( ptrd , pun , ptra, 'Z', ztrds ) 
+         CASE( jptra_xad )   ;   CALL trd_tra_adv( ptrd , pu , ptra, 'X', ztrds, Kmm ) 
+         CASE( jptra_yad )   ;   CALL trd_tra_adv( ptrd , pu , ptra, 'Y', ztrds, Kmm ) 
+         CASE( jptra_zad )   ;   CALL trd_tra_adv( ptrd , pu , ptra, 'Z', ztrds, Kmm ) 
          CASE DEFAULT                 ! other trends: just masked 
                                  ztrds(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)
@@ -178,5 +179,5 @@
 
 
-   SUBROUTINE trd_tra_adv( pf, pun, ptn, cdir, ptrd )
+   SUBROUTINE trd_tra_adv( pf, pu, pt, cdir, ptrd, Kmm )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE trd_tra_adv  ***
@@ -191,8 +192,9 @@
       !!----------------------------------------------------------------------
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   pf      ! advective flux in one direction
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   pun     ! now velocity   in one direction
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   ptn     ! now or before tracer 
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   pu      ! now velocity   in one direction
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   pt      ! now or before tracer 
       CHARACTER(len=1)                , INTENT(in   ) ::   cdir    ! X/Y/Z direction
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(  out) ::   ptrd    ! advective trend in one direction
+      INTEGER,  INTENT(in)                            ::   Kmm     ! time level index
       !
       INTEGER  ::   ji, jj, jk   ! dummy loop indices
@@ -215,6 +217,6 @@
             DO ji = fs_2, fs_jpim1   ! vector opt.
                ptrd(ji,jj,jk) = - (     pf (ji,jj,jk) - pf (ji-ii,jj-ij,jk-ik)                        &
-                 &                  - ( pun(ji,jj,jk) - pun(ji-ii,jj-ij,jk-ik) ) * ptn(ji,jj,jk)  )   &
-                 &              * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) * tmask(ji,jj,jk)
+                 &                  - ( pu(ji,jj,jk) - pu(ji-ii,jj-ij,jk-ik) ) * pt(ji,jj,jk)  )   &
+                 &              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk)
             END DO
          END DO
@@ -224,5 +226,5 @@
 
 
-   SUBROUTINE trd_tra_mng( ptrdx, ptrdy, ktrd, kt )
+   SUBROUTINE trd_tra_mng( ptrdx, ptrdy, ktrd, kt, Kmm )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE trd_tra_mng  ***
@@ -236,4 +238,5 @@
       INTEGER                   , INTENT(in   ) ::   ktrd    ! tracer trend index
       INTEGER                   , INTENT(in   ) ::   kt      ! time step
+      INTEGER                   , INTENT(in   ) ::   Kmm     ! time level index
       !!----------------------------------------------------------------------
 
@@ -243,11 +246,11 @@
 
       !                   ! 3D output of tracers trends using IOM interface
-      IF( ln_tra_trd )   CALL trd_tra_iom ( ptrdx, ptrdy, ktrd, kt )
+      IF( ln_tra_trd )   CALL trd_tra_iom ( ptrdx, ptrdy, ktrd, kt, Kmm )
 
       !                   ! Integral Constraints Properties for tracers trends                                       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
-      IF( ln_glo_trd )   CALL trd_glo( ptrdx, ptrdy, ktrd, 'TRA', kt )
+      IF( ln_glo_trd )   CALL trd_glo( ptrdx, ptrdy, ktrd, 'TRA', kt, Kmm )
 
       !                   ! Potential ENergy trends
-      IF( ln_PE_trd  )   CALL trd_pen( ptrdx, ptrdy, ktrd, kt, r2dt )
+      IF( ln_PE_trd  )   CALL trd_pen( ptrdx, ptrdy, ktrd, kt, r2dt, Kmm )
 
       !                   ! Mixed layer trends for active tracers
@@ -290,5 +293,5 @@
 
 
-   SUBROUTINE trd_tra_iom( ptrdx, ptrdy, ktrd, kt )
+   SUBROUTINE trd_tra_iom( ptrdx, ptrdy, ktrd, kt, Kmm )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE trd_tra_iom  ***
@@ -300,4 +303,5 @@
       INTEGER                   , INTENT(in   ) ::   ktrd    ! tracer trend index
       INTEGER                   , INTENT(in   ) ::   kt      ! time step
+      INTEGER                   , INTENT(in   ) ::   Kmm     ! time level index
       !!
       INTEGER ::   ji, jj, jk   ! dummy loop indices
@@ -326,6 +330,6 @@
                                   IF( ln_linssh ) THEN                   ! cst volume : adv flux through z=0 surface
                                      ALLOCATE( z2dx(jpi,jpj), z2dy(jpi,jpj) )
-                                     z2dx(:,:) = wn(:,:,1) * tsn(:,:,1,jp_tem) / e3t_n(:,:,1)
-                                     z2dy(:,:) = wn(:,:,1) * tsn(:,:,1,jp_sal) / e3t_n(:,:,1)
+                                     z2dx(:,:) = ww(:,:,1) * ts(:,:,1,jp_tem,Kmm) / e3t(:,:,1,Kmm)
+                                     z2dy(:,:) = ww(:,:,1) * ts(:,:,1,jp_sal,Kmm) / e3t(:,:,1,Kmm)
                                      CALL iom_put( "ttrd_sad", z2dx )
                                      CALL iom_put( "strd_sad", z2dy )
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdvor.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdvor.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdvor.F90	(revision 10946)
@@ -78,5 +78,5 @@
 
 
-   SUBROUTINE trd_vor( putrd, pvtrd, ktrd, kt )
+   SUBROUTINE trd_vor( putrd, pvtrd, ktrd, kt, Kmm )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE trd_vor  ***
@@ -88,4 +88,5 @@
       INTEGER                   , INTENT(in   ) ::   ktrd           ! trend index
       INTEGER                   , INTENT(in   ) ::   kt             ! time step
+      INTEGER                   , INTENT(in   ) ::   Kmm            ! time level index
       !
       INTEGER ::   ji, jj   ! dummy loop indices
@@ -94,27 +95,27 @@
 
       SELECT CASE( ktrd ) 
-      CASE( jpdyn_hpg )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_prg )   ! Hydrostatique Pressure Gradient 
-      CASE( jpdyn_keg )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_keg )   ! KE Gradient 
-      CASE( jpdyn_rvo )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_rvo )   ! Relative Vorticity 
-      CASE( jpdyn_pvo )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_pvo )   ! Planetary Vorticity Term 
-      CASE( jpdyn_ldf )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_ldf )   ! Horizontal Diffusion 
-      CASE( jpdyn_zad )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_zad )   ! Vertical Advection 
-      CASE( jpdyn_spg )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_spg )   ! Surface Pressure Grad. 
+      CASE( jpdyn_hpg )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_prg, Kmm )   ! Hydrostatique Pressure Gradient 
+      CASE( jpdyn_keg )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_keg, Kmm )   ! KE Gradient 
+      CASE( jpdyn_rvo )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_rvo, Kmm )   ! Relative Vorticity 
+      CASE( jpdyn_pvo )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_pvo, Kmm )   ! Planetary Vorticity Term 
+      CASE( jpdyn_ldf )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_ldf, Kmm )   ! Horizontal Diffusion 
+      CASE( jpdyn_zad )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_zad, Kmm )   ! Vertical Advection 
+      CASE( jpdyn_spg )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_spg, Kmm )   ! Surface Pressure Grad. 
       CASE( jpdyn_zdf )                                                      ! Vertical Diffusion 
          ztswu(:,:) = 0.e0   ;   ztswv(:,:) = 0.e0
          DO jj = 2, jpjm1                                                             ! wind stress trends
             DO ji = fs_2, fs_jpim1   ! vector opt.
-               ztswu(ji,jj) = 0.5 * ( utau_b(ji,jj) + utau(ji,jj) ) / ( e3u_n(ji,jj,1) * rau0 )
-               ztswv(ji,jj) = 0.5 * ( vtau_b(ji,jj) + vtau(ji,jj) ) / ( e3v_n(ji,jj,1) * rau0 )
+               ztswu(ji,jj) = 0.5 * ( utau_b(ji,jj) + utau(ji,jj) ) / ( e3u(ji,jj,1,Kmm) * rau0 )
+               ztswv(ji,jj) = 0.5 * ( vtau_b(ji,jj) + vtau(ji,jj) ) / ( e3v(ji,jj,1,Kmm) * rau0 )
             END DO
          END DO
          !
-         CALL trd_vor_zint( putrd, pvtrd, jpvor_zdf )                             ! zdf trend including surf./bot. stresses 
-         CALL trd_vor_zint( ztswu, ztswv, jpvor_swf )                             ! surface wind stress 
+         CALL trd_vor_zint( putrd, pvtrd, jpvor_zdf, Kmm )                             ! zdf trend including surf./bot. stresses 
+         CALL trd_vor_zint( ztswu, ztswv, jpvor_swf, Kmm )                             ! surface wind stress 
       CASE( jpdyn_bfr )
-         CALL trd_vor_zint( putrd, pvtrd, jpvor_bfr )                             ! Bottom stress
+         CALL trd_vor_zint( putrd, pvtrd, jpvor_bfr, Kmm )                             ! Bottom stress
          !
       CASE( jpdyn_atf )       ! last trends: perform the output of 2D vorticity trends
-         CALL trd_vor_iom( kt )
+         CALL trd_vor_iom( kt, Kmm )
       END SELECT
       !
@@ -122,5 +123,5 @@
 
 
-   SUBROUTINE trd_vor_zint_2d( putrdvor, pvtrdvor, ktrd )
+   SUBROUTINE trd_vor_zint_2d( putrdvor, pvtrdvor, ktrd, Kmm )
       !!----------------------------------------------------------------------------
       !!                  ***  ROUTINE trd_vor_zint  ***
@@ -150,4 +151,5 @@
       !!----------------------------------------------------------------------
       INTEGER                     , INTENT(in   ) ::   ktrd       ! ocean trend index
+      INTEGER                     , INTENT(in   ) ::   Kmm        ! time level index
       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   putrdvor   ! u vorticity trend 
       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   pvtrdvor   ! v vorticity trend
@@ -175,12 +177,12 @@
                ikbu = mbkv(ji,jj)
                ikbv = mbkv(ji,jj)            
-               zudpvor(ji,jj) = putrdvor(ji,jj) * e3u_n(ji,jj,ikbu) * e1u(ji,jj) * umask(ji,jj,ikbu)
-               zvdpvor(ji,jj) = pvtrdvor(ji,jj) * e3v_n(ji,jj,ikbv) * e2v(ji,jj) * vmask(ji,jj,ikbv)
+               zudpvor(ji,jj) = putrdvor(ji,jj) * e3u(ji,jj,ikbu,Kmm) * e1u(ji,jj) * umask(ji,jj,ikbu)
+               zvdpvor(ji,jj) = pvtrdvor(ji,jj) * e3v(ji,jj,ikbv,Kmm) * e2v(ji,jj) * vmask(ji,jj,ikbv)
             END DO
          END DO
          !
       CASE( jpvor_swf )        ! wind stress
-         zudpvor(:,:) = putrdvor(:,:) * e3u_n(:,:,1) * e1u(:,:) * umask(:,:,1)
-         zvdpvor(:,:) = pvtrdvor(:,:) * e3v_n(:,:,1) * e2v(:,:) * vmask(:,:,1)
+         zudpvor(:,:) = putrdvor(:,:) * e3u(:,:,1,Kmm) * e1u(:,:) * umask(:,:,1)
+         zvdpvor(:,:) = pvtrdvor(:,:) * e3v(:,:,1,Kmm) * e2v(:,:) * vmask(:,:,1)
          !
       END SELECT
@@ -207,5 +209,5 @@
 
 
-   SUBROUTINE trd_vor_zint_3d( putrdvor, pvtrdvor, ktrd )
+   SUBROUTINE trd_vor_zint_3d( putrdvor, pvtrdvor, ktrd , Kmm )
       !!----------------------------------------------------------------------------
       !!                  ***  ROUTINE trd_vor_zint  ***
@@ -236,4 +238,5 @@
       !
       INTEGER                         , INTENT(in   ) ::   ktrd       ! ocean trend index
+      INTEGER                         , INTENT(in   ) ::   Kmm        ! time level index
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   putrdvor   ! u vorticity trend 
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pvtrdvor   ! v vorticity trend
@@ -257,6 +260,6 @@
       ! putrdvor and pvtrdvor terms
       DO jk = 1,jpk
-        zudpvor(:,:) = zudpvor(:,:) + putrdvor(:,:,jk) * e3u_n(:,:,jk) * e1u(:,:) * umask(:,:,jk)
-        zvdpvor(:,:) = zvdpvor(:,:) + pvtrdvor(:,:,jk) * e3v_n(:,:,jk) * e2v(:,:) * vmask(:,:,jk)
+        zudpvor(:,:) = zudpvor(:,:) + putrdvor(:,:,jk) * e3u(:,:,jk,Kmm) * e1u(:,:) * umask(:,:,jk)
+        zvdpvor(:,:) = zvdpvor(:,:) + pvtrdvor(:,:,jk) * e3v(:,:,jk,Kmm) * e2v(:,:) * vmask(:,:,jk)
       END DO
 
@@ -298,5 +301,5 @@
 
 
-   SUBROUTINE trd_vor_iom( kt )
+   SUBROUTINE trd_vor_iom( kt , Kmm )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE trd_vor  ***
@@ -306,9 +309,10 @@
       !!----------------------------------------------------------------------
       INTEGER                   , INTENT(in   ) ::   kt             ! time step
+      INTEGER                   , INTENT(in   ) ::   Kmm            ! time level index
       !
       INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices
       INTEGER  ::   it, itmod        ! local integers
       REAL(wp) ::   zmean            ! local scalars
-      REAL(wp), DIMENSION(jpi,jpj) :: zun, zvn
+      REAL(wp), DIMENSION(jpi,jpj) :: zuu, zvv
       !!----------------------------------------------------------------------
 
@@ -327,6 +331,6 @@
 
       vor_avr   (:,:) = 0._wp
-      zun       (:,:) = 0._wp
-      zvn       (:,:) = 0._wp
+      zuu       (:,:) = 0._wp
+      zvv       (:,:) = 0._wp
       vor_avrtot(:,:) = 0._wp
       vor_avrres(:,:) = 0._wp
@@ -334,16 +338,16 @@
       ! Vertically averaged velocity
       DO jk = 1, jpk - 1
-         zun(:,:) = zun(:,:) + e1u(:,:) * un(:,:,jk) * e3u_n(:,:,jk)
-         zvn(:,:) = zvn(:,:) + e2v(:,:) * vn(:,:,jk) * e3v_n(:,:,jk)
+         zuu(:,:) = zuu(:,:) + e1u(:,:) * uu(:,:,jk,Kmm) * e3u(:,:,jk,Kmm)
+         zvv(:,:) = zvv(:,:) + e2v(:,:) * vv(:,:,jk,Kmm) * e3v(:,:,jk,Kmm)
       END DO
  
-      zun(:,:) = zun(:,:) * r1_hu_n(:,:)
-      zvn(:,:) = zvn(:,:) * r1_hv_n(:,:)
+      zuu(:,:) = zuu(:,:) * r1_hu_n(:,:)
+      zvv(:,:) = zvv(:,:) * r1_hv_n(:,:)
 
       ! Curl
       DO ji = 1, jpim1
          DO jj = 1, jpjm1
-            vor_avr(ji,jj) = (  ( zvn(ji+1,jj) - zvn(ji,jj) )    &
-               &              - ( zun(ji,jj+1) - zun(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) ) * fmask(ji,jj,1)
+            vor_avr(ji,jj) = (  ( zvv(ji+1,jj) - zvv(ji,jj) )    &
+               &              - ( zuu(ji,jj+1) - zuu(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) ) * fmask(ji,jj,1)
          END DO
       END DO
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/USR/usrdef_sbc.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/USR/usrdef_sbc.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/USR/usrdef_sbc.F90	(revision 10946)
@@ -39,5 +39,5 @@
 CONTAINS
 
-   SUBROUTINE usrdef_sbc_oce( kt )
+   SUBROUTINE usrdef_sbc_oce( kt, Kbb )
       !!---------------------------------------------------------------------
       !!                    ***  ROUTINE usrdef_sbc  ***
@@ -55,4 +55,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time step
+      INTEGER, INTENT(in) ::   Kbb  ! ocean time index
       !!
       INTEGER  ::   ji, jj                 ! dummy loop indices
@@ -119,5 +120,5 @@
             ! 23.5 deg : tropics
             qsr (ji,jj) =  230 * COS( 3.1415 * ( gphit(ji,jj) - 23.5 * zcos_sais1 ) / ( 0.9 * 180 ) )
-            qns (ji,jj) = ztrp * ( tsb(ji,jj,1,jp_tem) - t_star ) - qsr(ji,jj)
+            qns (ji,jj) = ztrp * ( ts(ji,jj,1,jp_tem,Kbb) - t_star ) - qsr(ji,jj)
             IF( gphit(ji,jj) >= 14.845 .AND. 37.2 >= gphit(ji,jj) ) THEN    ! zero at 37.8 deg, max at 24.6 deg
                emp  (ji,jj) =   zemp_S * zconv   &
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfdrg.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfdrg.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfdrg.F90	(revision 10946)
@@ -144,5 +144,5 @@
 
 
-   SUBROUTINE zdf_drg_exp( kt, pub, pvb, pua, pva )
+   SUBROUTINE zdf_drg_exp( kt, Kmm, pub, pvb, pua, pva )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE zdf_drg_exp  ***
@@ -157,4 +157,5 @@
       !!---------------------------------------------------------------------
       INTEGER                         , INTENT(in   ) ::   kt         ! ocean time-step index
+      INTEGER                         , INTENT(in   ) ::   Kmm        ! time level indices
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pub, pvb   ! the two components of the before velocity
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pua, pva   ! the two components of the velocity tendency
@@ -209,5 +210,5 @@
          ztrdu(:,:,:) = pua(:,:,:) - ztrdu(:,:,:)
          ztrdv(:,:,:) = pva(:,:,:) - ztrdv(:,:,:)
-         CALL trd_dyn( ztrdu(:,:,:), ztrdv(:,:,:), jpdyn_bfr, kt )
+         CALL trd_dyn( ztrdu(:,:,:), ztrdv(:,:,:), jpdyn_bfr, kt, Kmm )
          DEALLOCATE( ztrdu, ztrdv )
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfevd.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfevd.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfevd.F90	(revision 10946)
@@ -38,5 +38,5 @@
 CONTAINS
 
-   SUBROUTINE zdf_evd( kt, p_avm, p_avt )
+   SUBROUTINE zdf_evd( kt, Kmm, Krhs, p_avm, p_avt )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE zdf_evd  ***
@@ -56,4 +56,5 @@
       !!----------------------------------------------------------------------
       INTEGER                    , INTENT(in   ) ::   kt             ! ocean time-step indexocean time step
+      INTEGER                    , INTENT(in   ) ::   Kmm, Krhs      ! time level indices
       REAL(wp), DIMENSION(:,:,:) , INTENT(inout) ::   p_avm, p_avt   !  momentum and tracer Kz (w-points)
       !
@@ -117,5 +118,5 @@
       zavt_evd(:,:,:) = p_avt(:,:,:) - zavt_evd(:,:,:)   ! change in avt due to evd
       CALL iom_put( "avt_evd", zavt_evd )              ! output this change
-      IF( l_trdtra ) CALL trd_tra( kt, 'TRA', jp_tem, jptra_evd, zavt_evd )
+      IF( l_trdtra ) CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_evd, zavt_evd )
       !
    END SUBROUTINE zdf_evd
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfosm.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfosm.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfosm.F90	(revision 10946)
@@ -122,5 +122,5 @@
 
 
-   SUBROUTINE zdf_osm( kt, Kbb, Kmm, p_avm, p_avt )
+   SUBROUTINE zdf_osm( kt, Kbb, Kmm, Krhs, p_avm, p_avt )
       !!----------------------------------------------------------------------
       !!                   ***  ROUTINE zdf_osm  ***
@@ -158,5 +158,5 @@
       !!----------------------------------------------------------------------
       INTEGER                   , INTENT(in   ) ::  kt             ! ocean time step
-      INTEGER                   , INTENT(in   ) ::  Kbb, Kmm       ! ocean time level indices
+      INTEGER                   , INTENT(in   ) ::  Kbb, Kmm, Krhs ! ocean time level indices
       REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::  p_avm, p_avt   ! momentum and tracer Kz (w-points)
       !!
@@ -1661,6 +1661,6 @@
          ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs) - ztrds(:,:,:)
 !!bug gm jpttdzdf ==> jpttosm
-         CALL trd_tra( kt, 'TRA', jp_tem, jptra_zdf, ztrdt )
-         CALL trd_tra( kt, 'TRA', jp_sal, jptra_zdf, ztrds )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_zdf, ztrdt )
+         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_zdf, ztrds )
          DEALLOCATE( ztrdt )      ;     DEALLOCATE( ztrds )
       ENDIF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfphy.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfphy.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfphy.F90	(revision 10946)
@@ -220,5 +220,5 @@
 
 
-   SUBROUTINE zdf_phy( kt, Kbb, Kmm )
+   SUBROUTINE zdf_phy( kt, Kbb, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                     ***  ROUTINE zdf_phy  ***
@@ -233,5 +233,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt         ! ocean time-step index
-      INTEGER, INTENT(in) ::   Kbb, Kmm   ! ocean time level indices
+      INTEGER, INTENT(in) ::   Kbb, Kmm, Krhs   ! ocean time level indices
       !
       INTEGER ::   ji, jj, jk   ! dummy loop indice
@@ -264,5 +264,5 @@
       CASE( np_TKE )   ;   CALL zdf_tke( kt, Kbb, Kmm, zsh2, avm_k, avt_k )    ! TKE closure scheme for Kz
       CASE( np_GLS )   ;   CALL zdf_gls( kt, Kbb, Kmm, zsh2, avm_k, avt_k )    ! GLS closure scheme for Kz
-      CASE( np_OSM )   ;   CALL zdf_osm( kt, Kbb, Kmm      , avm_k, avt_k )    ! OSMOSIS closure scheme for Kz
+      CASE( np_OSM )   ;   CALL zdf_osm( kt, Kbb, Kmm, Krhs, avm_k, avt_k )    ! OSMOSIS closure scheme for Kz
 !     CASE( np_CST )                                  ! Constant Kz (reset avt, avm to the background value)
 !         ! avt_k and avm_k set one for all at initialisation phase
@@ -283,5 +283,5 @@
       ENDIF
       !
-      IF( ln_zdfevd )   CALL zdf_evd( kt, avm, avt )  !* convection: enhanced vertical eddy diffusivity
+      IF( ln_zdfevd )   CALL zdf_evd( kt, Kmm, Krhs, avm, avt )  !* convection: enhanced vertical eddy diffusivity
       !
       !                                         !* double diffusive mixing
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/nemogcm.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/nemogcm.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/nemogcm.F90	(revision 10946)
@@ -490,6 +490,6 @@
       IF( lk_diadct    )   CALL dia_dct_init    ! Sections tranports
                            CALL dia_hsb_init    ! heat content, salt content and volume budgets
-                           CALL     trd_init    ! Mixed-layer/Vorticity/Integral constraints trends
-                           CALL dia_obs_init( Nnn )        ! Initialize observational data
+                           CALL     trd_init( Nnn )    ! Mixed-layer/Vorticity/Integral constraints trends
+                           CALL dia_obs_init( Nnn )    ! Initialize observational data
                            CALL dia_tmb_init    ! TMB outputs
                            CALL dia_25h_init    ! 25h mean  outputs
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/step.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/step.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/step.F90	(revision 10946)
@@ -137,5 +137,5 @@
 
       !  VERTICAL PHYSICS
-                         CALL zdf_phy( kstp, Nbb, Nnn )   ! vertical physics update (top/bot drag, avt, avs, avm + MLD)
+                         CALL zdf_phy( kstp, Nbb, Nnn, Nrhs )   ! vertical physics update (top/bot drag, avt, avs, avm + MLD)
 
       !  LATERAL  PHYSICS
@@ -197,5 +197,5 @@
                          CALL dyn_ldf( kstp, Nbb, Nnn      , uu, vv, Nrhs )  ! lateral mixing
       IF( ln_zdfosm  )   CALL dyn_osm( kstp,      Nnn      , uu, vv, Nrhs )  ! OSMOSIS non-local velocity fluxes ==> RHS
-                         CALL dyn_hpg( kstp,      Nnn,       uu, vv, Nrhs )  ! horizontal gradient of Hydrostatic pressure
+                         CALL dyn_hpg( kstp,      Nnn      , uu, vv, Nrhs )  ! horizontal gradient of Hydrostatic pressure
                          CALL dyn_spg( kstp, Nbb, Nnn, Nrhs, uu, vv, ssh, uu_b, vv_b, Naa )  ! surface pressure gradient
 
@@ -242,9 +242,9 @@
       IF(  lk_asminc .AND. ln_asmiau .AND. &
          & ln_trainc )   CALL tra_asm_inc   ( kstp )  ! apply tracer assimilation increment
-                         CALL tra_sbc       ( kstp )  ! surface boundary condition
-      IF( ln_traqsr  )   CALL tra_qsr       ( kstp )  ! penetrative solar radiation qsr
-      IF( ln_trabbc  )   CALL tra_bbc       ( kstp )  ! bottom heat flux
-      IF( ln_trabbl  )   CALL tra_bbl       ( kstp )  ! advective (and/or diffusive) bottom boundary layer scheme
-      IF( ln_tradmp  )   CALL tra_dmp       ( kstp )  ! internal damping trends
+                         CALL tra_sbc       ( kstp, Nnn, Nrhs )  ! surface boundary condition
+      IF( ln_traqsr  )   CALL tra_qsr       ( kstp, Nnn, Nrhs )  ! penetrative solar radiation qsr
+      IF( ln_trabbc  )   CALL tra_bbc       ( kstp, Nnn, Nrhs )  ! bottom heat flux
+      IF( ln_trabbl  )   CALL tra_bbl       ( kstp, Nnn, Nrhs )  ! advective (and/or diffusive) bottom boundary layer scheme
+      IF( ln_tradmp  )   CALL tra_dmp       ( kstp, Nnn, Nrhs )  ! internal damping trends
       IF( ln_bdy     )   CALL bdy_tra_dmp   ( kstp )  ! bdy damping trends
 #if defined key_agrif
@@ -256,5 +256,5 @@
       IF( lrst_oce .AND. ln_zdfosm ) &
            &             CALL osm_rst( kstp, Nnn, 'WRITE' )! write OSMOSIS outputs + wn (so must do here) to restarts
-                         CALL tra_ldf       ( kstp, Nnn )  ! lateral mixing
+                         CALL tra_ldf( kstp, Nnn, Nrhs )   ! lateral mixing
 
 !!gm : why CALL to dia_ptr has been moved here??? (use trends info?)
@@ -262,5 +262,5 @@
 !!gm
                          CALL tra_zdf( kstp, Nbb, Nnn, Nrhs, ts, Naa  )  ! vert. mixing & after tracer	==> after
-      IF( ln_zdfnpc  )   CALL tra_npc       ( kstp )  ! update after fields by non-penetrative convection
+      IF( ln_zdfnpc  )   CALL tra_npc( kstp,      Nnn, Nrhs           )  ! update after fields by non-penetrative convection
 
       !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
@@ -281,6 +281,6 @@
 !! 
 !!jc2: dynnxt must be the latest call. e3t_b are indeed updated in that routine
-                         CALL tra_nxt       ( kstp )  ! finalize (bcs) tracer fields at next time step and swap
-                         CALL dyn_nxt       ( kstp )  ! finalize (bcs) velocities at next time step and swap (always called after tra_nxt)
+                         CALL tra_nxt       ( kstp, Nnn, Nrhs )  ! finalize (bcs) tracer fields at next time step and swap
+                         CALL dyn_nxt       ( kstp, Nnn )        ! finalize (bcs) velocities at next time step and swap (always called after tra_nxt)
                          CALL ssh_swp       ( kstp )  ! swap of sea surface height
       IF(.NOT.ln_linssh) CALL dom_vvl_sf_swp( kstp )  ! swap of vertical scale factors
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcadv.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcadv.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcadv.F90	(revision 10946)
@@ -116,5 +116,5 @@
          !
          IF( ln_ldfeiv .AND. .NOT. ln_traldf_triad )   & 
-            &              CALL ldf_eiv_trp( kt, nittrc000, zun, zvn, zwn, 'TRC', Kmm )  ! add the eiv transport
+            &              CALL ldf_eiv_trp( kt, nittrc000, zun, zvn, zwn, 'TRC', Kmm, Krhs )  ! add the eiv transport
          !
          IF( ln_mle    )   CALL tra_mle_trp( kt, nittrc000, zun, zvn, zwn, 'TRC' )  ! add the mle transport
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcbbl.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcbbl.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcbbl.F90	(revision 10946)
@@ -36,5 +36,5 @@
 CONTAINS
 
-   SUBROUTINE trc_bbl( kt )
+   SUBROUTINE trc_bbl( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE bbl  ***
@@ -46,4 +46,5 @@
       !!----------------------------------------------------------------------  
       INTEGER, INTENT( in ) ::   kt   ! ocean time-step 
+      INTEGER, INTENT( in ) ::   Kmm, Krhs  ! time level indices
       INTEGER :: jn                   ! loop index
       CHARACTER (len=22) :: charout
@@ -88,5 +89,5 @@
         DO jn = 1, jptra
            ztrtrd(:,:,:,jn) = tra(:,:,:,jn) - ztrtrd(:,:,:,jn)
-           CALL trd_tra( kt, 'TRC', jn, jptra_bbl, ztrtrd(:,:,:,jn) )
+           CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_bbl, ztrtrd(:,:,:,jn) )
         END DO
         DEALLOCATE( ztrtrd ) ! temporary save of trends
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcdmp.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcdmp.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcdmp.F90	(revision 10946)
@@ -63,5 +63,5 @@
 
 
-   SUBROUTINE trc_dmp( kt )
+   SUBROUTINE trc_dmp( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                   ***  ROUTINE trc_dmp  ***
@@ -83,4 +83,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
+      INTEGER, INTENT(in) ::   Kmm, Krhs  ! time level indices
       !
       INTEGER ::   ji, jj, jk, jn, jl   ! dummy loop indices
@@ -146,5 +147,5 @@
             IF( l_trdtrc ) THEN
                ztrtrd(:,:,:) = tra(:,:,:,jn) -  ztrtrd(:,:,:)
-               CALL trd_tra( kt, 'TRC', jn, jptra_dmp, ztrtrd )
+               CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_dmp, ztrtrd )
             END IF
             !                                                       ! ===========
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcldf.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcldf.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcldf.F90	(revision 10946)
@@ -51,5 +51,5 @@
 CONTAINS
 
-   SUBROUTINE trc_ldf( kt, Kmm )
+   SUBROUTINE trc_ldf( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE tra_ldf  ***
@@ -59,5 +59,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
-      INTEGER, INTENT( in ) ::   Kmm  ! ocean time-level index
+      INTEGER, INTENT( in ) ::   Kmm, Krhs  ! ocean time-level index
       !
       INTEGER            :: ji, jj, jk, jn
@@ -106,5 +106,5 @@
         DO jn = 1, jptra
            ztrtrd(:,:,:,jn) = tra(:,:,:,jn) - ztrtrd(:,:,:,jn)
-           CALL trd_tra( kt, 'TRC', jn, jptra_ldf, ztrtrd(:,:,:,jn) )
+           CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_ldf, ztrtrd(:,:,:,jn) )
         END DO
         DEALLOCATE( ztrtrd )
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcnxt.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcnxt.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcnxt.F90	(revision 10946)
@@ -54,5 +54,5 @@
 CONTAINS
 
-   SUBROUTINE trc_nxt( kt )
+   SUBROUTINE trc_nxt( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                   ***  ROUTINE trcnxt  ***
@@ -79,4 +79,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT( in ) ::   kt     ! ocean time-step index
+      INTEGER, INTENT( in ) ::   Kmm, Krhs  ! time level indices
       !
       INTEGER  ::   jk, jn   ! dummy loop indices
@@ -106,5 +107,5 @@
          IF( ln_traldf_iso ) THEN                       ! diagnose the "pure" Kz diffusive trend 
             DO jn = 1, jptra
-               CALL trd_tra( kt, 'TRC', jn, jptra_zdfp, ztrdt(:,:,:,jn) )
+               CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_zdfp, ztrdt(:,:,:,jn) )
             ENDDO
          ENDIF
@@ -128,5 +129,5 @@
          !
          DO jn = 1, jptra
-            CALL trd_tra( kt, 'TRC', jn, jptra_tot, ztrdt(:,:,:,jn) )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_tot, ztrdt(:,:,:,jn) )
          ENDDO
          !
@@ -150,5 +151,5 @@
             ztrdt(:,:,:,:) = 0._wp            
             DO jn = 1, jptra
-               CALL trd_tra( kt, 'TRC', jn, jptra_atf, ztrdt(:,:,:,jn) )
+               CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_atf, ztrdt(:,:,:,jn) )
             ENDDO
          END IF
@@ -157,5 +158,5 @@
          IF( .NOT. l_offline ) THEN ! Leap-Frog + Asselin filter time stepping
             IF( ln_linssh ) THEN   ;   CALL tra_nxt_fix( kt, nittrc000,         'TRC', trb, trn, tra, jptra )  !     linear ssh
-            ELSE                   ;   CALL tra_nxt_vvl( kt, nittrc000, rdttrc, 'TRC', trb, trn, tra,      &
+            ELSE                   ;   CALL tra_nxt_vvl( kt, Kmm, Krhs, nittrc000, rdttrc, 'TRC', trb, trn, tra,      &
               &                                                                   sbc_trc, sbc_trc_b, jptra )  ! non-linear ssh
             ENDIF
@@ -173,5 +174,5 @@
                ztrdt(:,:,jk,jn) = ( trb(:,:,jk,jn) - ztrdt(:,:,jk,jn) ) * zfact 
             END DO
-            CALL trd_tra( kt, 'TRC', jn, jptra_atf, ztrdt(:,:,:,jn) )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_atf, ztrdt(:,:,:,jn) )
          END DO
       END IF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcrad.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcrad.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcrad.F90	(revision 10946)
@@ -37,5 +37,5 @@
 CONTAINS
 
-   SUBROUTINE trc_rad( kt )
+   SUBROUTINE trc_rad( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE trc_rad  ***
@@ -53,4 +53,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
+      INTEGER, INTENT(in) ::   Kmm, Krhs  ! time level indices
       !
       CHARACTER (len=22) :: charout
@@ -59,9 +60,9 @@
       IF( ln_timing )   CALL timing_start('trc_rad')
       !
-      IF( ln_age     )   CALL trc_rad_sms( kt, trb, trn, jp_age , jp_age                )  !  AGE
-      IF( ll_cfc     )   CALL trc_rad_sms( kt, trb, trn, jp_cfc0, jp_cfc1               )  !  CFC model
-      IF( ln_c14     )   CALL trc_rad_sms( kt, trb, trn, jp_c14 , jp_c14                )  !  C14
-      IF( ln_pisces  )   CALL trc_rad_sms( kt, trb, trn, jp_pcs0, jp_pcs1, cpreserv='Y' )  !  PISCES model
-      IF( ln_my_trc  )   CALL trc_rad_sms( kt, trb, trn, jp_myt0, jp_myt1               )  !  MY_TRC model
+      IF( ln_age     )   CALL trc_rad_sms( kt, Kmm, Krhs, trb, trn, jp_age , jp_age                )  !  AGE
+      IF( ll_cfc     )   CALL trc_rad_sms( kt, Kmm, Krhs, trb, trn, jp_cfc0, jp_cfc1               )  !  CFC model
+      IF( ln_c14     )   CALL trc_rad_sms( kt, Kmm, Krhs, trb, trn, jp_c14 , jp_c14                )  !  C14
+      IF( ln_pisces  )   CALL trc_rad_sms( kt, Kmm, Krhs, trb, trn, jp_pcs0, jp_pcs1, cpreserv='Y' )  !  PISCES model
+      IF( ln_my_trc  )   CALL trc_rad_sms( kt, Kmm, Krhs, trb, trn, jp_myt0, jp_myt1               )  !  MY_TRC model
       !
       IF(ln_ctl) THEN      ! print mean trends (used for debugging)
@@ -113,5 +114,5 @@
 
 
-   SUBROUTINE trc_rad_sms( kt, ptrb, ptrn, jp_sms0, jp_sms1, cpreserv )
+   SUBROUTINE trc_rad_sms( kt, Kmm, Krhs, ptrb, ptrn, jp_sms0, jp_sms1, cpreserv )
       !!-----------------------------------------------------------------------------
       !!                  ***  ROUTINE trc_rad_sms  ***
@@ -130,4 +131,5 @@
       !!--------------------------------------------------------------------------------
       INTEGER                                , INTENT(in   ) ::   kt                 ! ocean time-step index
+      INTEGER                                , INTENT(in   ) ::   Kmm, Krhs          ! time level indices
       INTEGER                                , INTENT(in   ) ::   jp_sms0, jp_sms1   ! First & last index of the passive tracer model
       REAL(wp), DIMENSION (jpi,jpj,jpk,jptra), INTENT(inout) ::   ptrb    , ptrn     ! before and now traceur concentration
@@ -183,5 +185,5 @@
             IF( l_trdtrc ) THEN
                ztrtrd(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrd(:,:,:) ) * zs2rdt
-               CALL trd_tra( kt, 'TRC', jn, jptra_radb, ztrtrd )       ! Asselin-like trend handling
+               CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_radb, ztrtrd )       ! Asselin-like trend handling
             ENDIF
             !
@@ -233,5 +235,5 @@
             IF( l_trdtrc ) THEN
                ztrtrd(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrd(:,:,:) ) * zs2rdt
-               CALL trd_tra( kt, 'TRC', jn, jptra_radn, ztrtrd )       ! standard     trend handling
+               CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_radn, ztrtrd )       ! standard     trend handling
             ENDIF
             !
@@ -261,5 +263,5 @@
             IF( l_trdtrc ) THEN
                ztrtrd(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrd(:,:,:) ) * zs2rdt
-               CALL trd_tra( kt, 'TRC', jn, jptra_radb, ztrtrd )       ! Asselin-like trend handling
+               CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_radb, ztrtrd )       ! Asselin-like trend handling
             ENDIF
             !
@@ -270,5 +272,5 @@
             IF( l_trdtrc ) THEN
                ztrtrd(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrd(:,:,:) ) * zs2rdt
-               CALL trd_tra( kt, 'TRC', jn, jptra_radn, ztrtrd )       ! standard     trend handling
+               CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_radn, ztrtrd )       ! standard     trend handling
             ENDIF
             !
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcsbc.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcsbc.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcsbc.F90	(revision 10946)
@@ -37,5 +37,5 @@
 CONTAINS
 
-   SUBROUTINE trc_sbc ( kt )
+   SUBROUTINE trc_sbc ( kt, Kmm, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE trc_sbc  ***
@@ -59,4 +59,5 @@
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt   ! ocean time-step index
+      INTEGER, INTENT(in) ::   Kmm, Krhs  ! time level indices
       !
       INTEGER  ::   ji, jj, jn                      ! dummy loop indices
@@ -163,5 +164,5 @@
          IF( l_trdtrc ) THEN
             ztrtrd(:,:,:) = tra(:,:,:,jn) - ztrtrd(:,:,:)
-            CALL trd_tra( kt, 'TRC', jn, jptra_nsr, ztrtrd )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_nsr, ztrtrd )
          END IF
          !                                                       ! ===========
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trctrp.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trctrp.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trctrp.F90	(revision 10946)
@@ -61,10 +61,10 @@
       IF( .NOT. lk_c1d ) THEN
          !
-                                CALL trc_sbc    ( kt )      ! surface boundary condition
-         IF( ln_trabbl )        CALL trc_bbl    ( kt )      ! advective (and/or diffusive) bottom boundary layer scheme
-         IF( ln_trcdmp )        CALL trc_dmp    ( kt )      ! internal damping trends
-         IF( ln_bdy )           CALL trc_bdy_dmp( kt )      ! BDY damping trends
+                                CALL trc_sbc    ( kt, Kmm, Krhs )      ! surface boundary condition
+         IF( ln_trabbl )        CALL trc_bbl    ( kt, Kmm, Krhs )      ! advective (and/or diffusive) bottom boundary layer scheme
+         IF( ln_trcdmp )        CALL trc_dmp    ( kt, Kmm, Krhs )      ! internal damping trends
+         IF( ln_bdy )           CALL trc_bdy_dmp( kt )                 ! BDY damping trends
                                 CALL trc_adv    ( kt, Kbb, Kmm, tr, Krhs )      ! horizontal & vertical advection 
-         !                                                         ! Partial top/bottom cell: GRADh( trb )  
+         !                                                             ! Partial top/bottom cell: GRADh( trb )  
          IF( ln_zps ) THEN
            IF( ln_isfcav ) THEN ; CALL zps_hde_isf( kt, jptra, trb, pgtu=gtru, pgtv=gtrv, pgtui=gtrui, pgtvi=gtrvi )  ! both top & bottom
@@ -73,20 +73,20 @@
          ENDIF
          !                                                      
-                                CALL trc_ldf    ( kt, Kmm ) ! lateral mixing
+                                CALL trc_ldf    ( kt, Kmm, Krhs ) ! lateral mixing
 #if defined key_agrif
          IF(.NOT. Agrif_Root()) CALL Agrif_Sponge_trc       ! tracers sponge
 #endif
                                 CALL trc_zdf( kt, Kbb, Kmm, Krhs, tr, Kaa  )  ! vert. mixing & after tracer	==> after
-                                CALL trc_nxt    ( kt )      ! tracer fields at next time step     
-         IF( ln_trcrad )        CALL trc_rad    ( kt )      ! Correct artificial negative concentrations
-         IF( ln_trcdmp_clo )    CALL trc_dmp_clo( kt )      ! internal damping trends on closed seas only
+                                CALL trc_nxt    ( kt, Kmm, Krhs )        ! tracer fields at next time step     
+         IF( ln_trcrad )        CALL trc_rad    ( kt, Kmm, Krhs )        ! Correct artificial negative concentrations
+         IF( ln_trcdmp_clo )    CALL trc_dmp_clo( kt )                   ! internal damping trends on closed seas only
 
          !
       ELSE                                               ! 1D vertical configuration
-                                CALL trc_sbc( kt )            ! surface boundary condition
-         IF( ln_trcdmp )        CALL trc_dmp( kt )            ! internal damping trends
+                                CALL trc_sbc( kt, Kmm, Krhs )            ! surface boundary condition
+         IF( ln_trcdmp )        CALL trc_dmp( kt, Kmm, Krhs )            ! internal damping trends
                                 CALL trc_zdf( kt, Kbb, Kmm, Krhs, tr, Kaa  )  ! vert. mixing & after tracer	==> after
-                                CALL trc_nxt( kt )            ! tracer fields at next time step     
-          IF( ln_trcrad )       CALL trc_rad( kt )            ! Correct artificial negative concentrations
+                                CALL trc_nxt( kt, Kmm, Krhs )            ! tracer fields at next time step     
+          IF( ln_trcrad )       CALL trc_rad( kt, Kmm, Krhs )            ! Correct artificial negative concentrations
          !
       END IF
Index: /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trczdf.F90
===================================================================
--- /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trczdf.F90	(revision 10945)
+++ /NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trczdf.F90	(revision 10946)
@@ -63,5 +63,5 @@
                ztrtrd(:,:,jk,jn) = ( ( ptr(:,:,jk,jn,Kaa) - ptr(:,:,jk,jn,Kbb) ) / r2dttrc ) - ztrtrd(:,:,jk,jn)
             END DO
-            CALL trd_tra( kt, 'TRC', jn, jptra_zdf, ztrtrd(:,:,:,jn) )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_zdf, ztrtrd(:,:,:,jn) )
          END DO
       ENDIF
