Index: NEMO/trunk/src/OCE/TRD/trdvor.F90
===================================================================
--- NEMO/trunk/src/OCE/TRD/trdvor.F90	(revision 11536)
+++ NEMO/trunk/src/OCE/TRD/trdvor.F90	(revision 12377)
@@ -56,5 +56,5 @@
 
    !! * Substitutions
-#  include "vectopt_loop_substitute.h90"
+#  include "do_loop_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/OCE 4.0 , NEMO Consortium (2018)
@@ -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,25 @@
 
       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 )
-            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 
+         DO_2D_00_00
+            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_2D
+         !
+         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 +121,5 @@
 
 
-   SUBROUTINE trd_vor_zint_2d( putrdvor, pvtrdvor, ktrd )
+   SUBROUTINE trd_vor_zint_2d( putrdvor, pvtrdvor, ktrd, Kmm )
       !!----------------------------------------------------------------------------
       !!                  ***  ROUTINE trd_vor_zint  ***
@@ -150,4 +149,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
@@ -171,22 +171,20 @@
       !
       CASE( jpvor_bfr )        ! bottom friction
-         DO jj = 2, jpjm1
-            DO ji = fs_2, fs_jpim1 
-               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)
-            END DO
-         END DO
+         DO_2D_00_00
+            ikbu = mbkv(ji,jj)
+            ikbv = mbkv(ji,jj)            
+            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_2D
          !
       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
 
       ! Average except for Beta.V
-      zudpvor(:,:) = zudpvor(:,:) * r1_hu_n(:,:)
-      zvdpvor(:,:) = zvdpvor(:,:) * r1_hv_n(:,:)
+      zudpvor(:,:) = zudpvor(:,:) * r1_hu(:,:,Kmm)
+      zvdpvor(:,:) = zvdpvor(:,:) * r1_hv(:,:,Kmm)
    
       ! Curl
@@ -207,5 +205,5 @@
 
 
-   SUBROUTINE trd_vor_zint_3d( putrdvor, pvtrdvor, ktrd )
+   SUBROUTINE trd_vor_zint_3d( putrdvor, pvtrdvor, ktrd , Kmm )
       !!----------------------------------------------------------------------------
       !!                  ***  ROUTINE trd_vor_zint  ***
@@ -236,4 +234,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 +256,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
 
@@ -273,10 +272,10 @@
          END DO
          ! Average of the Curl and Surface mask
-         vortrd(:,:,jpvor_bev) = vortrd(:,:,jpvor_bev) * r1_hu_n(:,:) * fmask(:,:,1)
+         vortrd(:,:,jpvor_bev) = vortrd(:,:,jpvor_bev) * r1_hu(:,:,Kmm) * fmask(:,:,1)
       ENDIF
       !
       ! Average 
-      zudpvor(:,:) = zudpvor(:,:) * r1_hu_n(:,:)
-      zvdpvor(:,:) = zvdpvor(:,:) * r1_hv_n(:,:)
+      zudpvor(:,:) = zudpvor(:,:) * r1_hu(:,:,Kmm)
+      zvdpvor(:,:) = zvdpvor(:,:) * r1_hv(:,:,Kmm)
       !
       ! Curl
@@ -298,5 +297,5 @@
 
 
-   SUBROUTINE trd_vor_iom( kt )
+   SUBROUTINE trd_vor_iom( kt , Kmm )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE trd_vor  ***
@@ -306,9 +305,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 +327,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 +334,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(:,:,Kmm)
+      zvv(:,:) = zvv(:,:) * r1_hv(:,:,Kmm)
 
       ! 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
