Index: /branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/DYN/dynzdf_imp.F90
===================================================================
--- /branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/DYN/dynzdf_imp.F90	(revision 4405)
+++ /branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/DYN/dynzdf_imp.F90	(revision 4406)
@@ -63,5 +63,4 @@
       USE oce     , ONLY:  zwd  => ta       , zws   => sa   ! (ta,sa) used as 3D workspace
       USE wrk_nemo, ONLY:   zwi => wrk_3d_3                 ! 3D workspace
-      USE arpdebugging, ONLY: dump_array
       !! DCSE_NEMO: need additional directives for renamed module variables
 !FTRANS zwd :I :I :z
@@ -90,6 +89,4 @@
       z1_p2dt = 1._wp / p2dt      ! inverse of the timestep
 
-      !CALL dump_array(kt, 'utau_pre_zdf',utau(:,:),withHalos=.TRUE.)
-      !CALL dump_array(kt, 'utaub_pre_zdf',utau_b(:,:),withHalos=.TRUE.)
 
       ! 1. Vertical diffusion on u
@@ -105,5 +102,5 @@
       DO jj = 2, jpjm1 
          DO ji = 2, jpim1
-            DO jk = 1, jpkm1
+            DO jk = 1, jpkfm1
                zcoef = - p2dt / fse3u(ji,jj,jk)
                zzwi          = zcoef * avmu (ji,jj,jk  ) / fse3uw(ji,jj,jk  )
@@ -119,5 +116,5 @@
       END DO
 #else
-      DO jk = 1, jpkm1        ! Matrix
+      DO jk = 1, jpkfm1        ! Matrix
          DO jj = 2, jpjm1 
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -160,5 +157,5 @@
             ua(ji,jj,1) = ub(ji,jj,1) + p2dt * (  ua(ji,jj,1) + 0.5_wp * ( utau_b(ji,jj) + utau(ji,jj) )   &
                &                                                       / ( fse3u(ji,jj,1) * rau0       )  )
-            DO jk = 2, jpkm1
+            DO jk = 2, jpkfm1
                zzwibd = zwi(ji,jj,jk) / zwd(ji,jj,jk-1)
                !==  First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1   (increasing k)  ==
@@ -169,10 +166,10 @@
             END DO
             !==  third recurrence : SOLk = ( Lk - Uk * Ek+1 ) / Dk  ==
-            ua(ji,jj,jpkm1) = ua(ji,jj,jpkm1) / zwd(ji,jj,jpkm1)
-            DO jk = jpk-2, 1, -1
+            ua(ji,jj,jpkfm1) = ua(ji,jj,jpkfm1) / zwd(ji,jj,jpkfm1)
+            DO jk = jpkf-2, 1, -1
                ua(ji,jj,jk) = ( ua(ji,jj,jk) - zws(ji,jj,jk) * ua(ji,jj,jk+1) ) / zwd(ji,jj,jk)
             END DO
             ! Normalization to obtain the general momentum trend ua
-            DO jk = 1, jpkm1
+            DO jk = 1, jpkfm1
                ua(ji,jj,jk) = ( ua(ji,jj,jk) - ub(ji,jj,jk) ) * z1_p2dt
             END DO
@@ -180,5 +177,5 @@
       END DO
 #else
-      DO jk = 2, jpkm1        !==  First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1   (increasing k)  ==
+      DO jk = 2, jpkfm1        !==  First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1   (increasing k)  ==
          DO jj = 2, jpjm1   
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -194,5 +191,5 @@
          END DO
       END DO
-      DO jk = 2, jpkm1
+      DO jk = 2, jpkfm1
          DO jj = 2, jpjm1   
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -205,8 +202,8 @@
       DO jj = 2, jpjm1        !==  third recurrence : SOLk = ( Lk - Uk * Ek+1 ) / Dk  ==
          DO ji = fs_2, fs_jpim1   ! vector opt.
-            ua(ji,jj,jpkm1) = ua(ji,jj,jpkm1) / zwd(ji,jj,jpkm1)
-         END DO
-      END DO
-      DO jk = jpk-2, 1, -1
+            ua(ji,jj,jpkfm1) = ua(ji,jj,jpkfm1) / zwd(ji,jj,jpkfm1)
+         END DO
+      END DO
+      DO jk = jpkf-2, 1, -1
          DO jj = 2, jpjm1   
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -216,5 +213,5 @@
       END DO
       ! Normalization to obtain the general momentum trend ua
-      DO jk = 1, jpkm1
+      DO jk = 1, jpkfm1
          DO jj = 2, jpjm1   
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -237,5 +234,5 @@
       DO jj = 2, jpjm1   
          DO ji = 2, jpim1
-            DO jk = 1, jpkm1        ! Matrix
+            DO jk = 1, jpkfm1        ! Matrix
                zcoef = -p2dt / fse3v(ji,jj,jk)
                zzwi          = zcoef * avmv (ji,jj,jk  ) / fse3vw(ji,jj,jk  )
@@ -251,5 +248,5 @@
       END DO
 #else
-      DO jk = 1, jpkm1        ! Matrix
+      DO jk = 1, jpkfm1        ! Matrix
          DO jj = 2, jpjm1   
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -292,5 +289,5 @@
             va(ji,jj,1) = vb(ji,jj,1) + p2dt * (  va(ji,jj,1) + 0.5_wp * ( vtau_b(ji,jj) + vtau(ji,jj) )   &
                &                                                       / ( fse3v(ji,jj,1) * rau0       )  )
-            DO jk = 2, jpkm1
+            DO jk = 2, jpkfm1
                zzwibd = zwi(ji,jj,jk) / zwd(ji,jj,jk-1)
                !==  First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1   (increasing k)  ==
@@ -301,10 +298,10 @@
             END DO
             !==  third recurrence : SOLk = ( Lk - Uk * SOLk+1 ) / Dk  ==
-            va(ji,jj,jpkm1) = va(ji,jj,jpkm1) / zwd(ji,jj,jpkm1)
-            DO jk = jpk-2, 1, -1
+            va(ji,jj,jpkfm1) = va(ji,jj,jpkfm1) / zwd(ji,jj,jpkfm1)
+            DO jk = jpkf-2, 1, -1
                va(ji,jj,jk) = ( va(ji,jj,jk) - zws(ji,jj,jk) * va(ji,jj,jk+1) ) / zwd(ji,jj,jk)
             END DO
             ! Normalization to obtain the general momentum trend va
-            DO jk = 1, jpkm1
+            DO jk = 1, jpkfm1
                va(ji,jj,jk) = ( va(ji,jj,jk) - vb(ji,jj,jk) ) * z1_p2dt
             END DO
@@ -312,5 +309,5 @@
       END DO
 #else
-      DO jk = 2, jpkm1        !==  First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1   (increasing k)  ==
+      DO jk = 2, jpkfm1        !==  First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1   (increasing k)  ==
          DO jj = 2, jpjm1   
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -326,5 +323,5 @@
          END DO
       END DO
-      DO jk = 2, jpkm1
+      DO jk = 2, jpkfm1
          DO jj = 2, jpjm1
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -337,8 +334,8 @@
       DO jj = 2, jpjm1        !==  third recurrence : SOLk = ( Lk - Uk * SOLk+1 ) / Dk  ==
          DO ji = fs_2, fs_jpim1   ! vector opt.
-            va(ji,jj,jpkm1) = va(ji,jj,jpkm1) / zwd(ji,jj,jpkm1)
-         END DO
-      END DO
-      DO jk = jpk-2, 1, -1
+            va(ji,jj,jpkfm1) = va(ji,jj,jpkfm1) / zwd(ji,jj,jpkfm1)
+         END DO
+      END DO
+      DO jk = jpkf-2, 1, -1
          DO jj = 2, jpjm1   
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -349,5 +346,5 @@
 
       ! Normalization to obtain the general momentum trend va
-      DO jk = 1, jpkm1
+      DO jk = 1, jpkfm1
          DO jj = 2, jpjm1   
             DO ji = fs_2, fs_jpim1   ! vector opt.
Index: /branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90
===================================================================
--- /branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90	(revision 4405)
+++ /branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90	(revision 4406)
@@ -285,5 +285,5 @@
             DO ji = 1, jpi
                zpelc(ji,jj,1) =  MAX( rn2b(ji,jj,1), 0._wp ) * fsdepw(ji,jj,1) * fse3w(ji,jj,1)
-               DO jk = 2, jpk
+               DO jk = 2, jpkf
                   zpelc(ji,jj,jk)  = zpelc(ji,jj,jk-1) + MAX( rn2b(ji,jj,jk), 0._wp ) * fsdepw(ji,jj,jk) * fse3w(ji,jj,jk)
                END DO
@@ -292,5 +292,5 @@
 #else
          zpelc(:,:,1) =  MAX( rn2b(:,:,1), 0._wp ) * fsdepw(:,:,1) * fse3w(:,:,1)
-         DO jk = 2, jpk
+         DO jk = 2, jpkf
             zpelc(:,:,jk)  = zpelc(:,:,jk-1) + MAX( rn2b(:,:,jk), 0._wp ) * fsdepw(:,:,jk) * fse3w(:,:,jk)
          END DO
@@ -303,7 +303,7 @@
             DO ji = 1, jpi               !      with us=0.016*wind(starting from jpk-1)
                zus  = zcof * taum(ji,jj)
-               DO jk = jpkm1, 2, -1
-#else
-         DO jk = jpkm1, 2, -1
+               DO jk = jpkfm1, 2, -1
+#else
+         DO jk = jpkfm1, 2, -1
             DO jj = 1, jpj               ! Last w-level at which zpelc>=0.5*us*us 
                DO ji = 1, jpi            !      with us=0.016*wind(starting from jpk-1)
@@ -330,8 +330,8 @@
             DO ji = 2, jpim1
                zus  = zcof * SQRT( taum(ji,jj) )           ! Stokes drift
-               DO jk = 2, jpkm1
-#else
-!CDIR NOVERRCHK
-         DO jk = 2, jpkm1         !* TKE Langmuir circulation source term added to en
+               DO jk = 2, jpkfm1
+#else
+!CDIR NOVERRCHK
+         DO jk = 2, jpkfm1         !* TKE Langmuir circulation source term added to en
 !CDIR NOVERRCHK
             DO jj = 2, jpjm1
@@ -363,7 +363,7 @@
       DO jj = 1, jpj
          DO ji = 1, jpi
-            DO jk = 2, jpkm1
-#else
-      DO jk = 2, jpkm1           !* Shear production at uw- and vw-points (energy conserving form)
+            DO jk = 2, jpkfm1
+#else
+      DO jk = 2, jpkfm1           !* Shear production at uw- and vw-points (energy conserving form)
          DO jj = 1, jpj                 ! here avmu, avmv used as workspace
             DO ji = 1, jpi
@@ -385,7 +385,7 @@
       DO jj = 2, jpjm1
          DO ji = 2, jpim1
-            DO jk = 2, jpkm1     !* Matrix and right hand side in en
-#else
-      DO jk = 2, jpkm1           !* Matrix and right hand side in en
+            DO jk = 2, jpkfm1     !* Matrix and right hand side in en
+#else
+      DO jk = 2, jpkfm1           !* Matrix and right hand side in en
          DO jj = 2, jpjm1
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -415,18 +415,18 @@
          DO ji = 2, jpim1
             ! First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1
-            DO jk = 3, jpkm1
+            DO jk = 3, jpkfm1
                zdiag(ji,jj,jk) = zdiag(ji,jj,jk) - zd_lw(ji,jj,jk) * zd_up(ji,jj,jk-1) / zdiag(ji,jj,jk-1)
             END DO
             ! Second recurrence : Lk = RHSk - Lk / Dk-1 * Lk-1
             zd_lw(ji,jj,2) = en(ji,jj,2) - zd_lw(ji,jj,2) * en(ji,jj,1)    ! Surface boudary conditions on tke
-            DO jk = 3, jpkm1
+            DO jk = 3, jpkfm1
                zd_lw(ji,jj,jk) = en(ji,jj,jk) - zd_lw(ji,jj,jk) / zdiag(ji,jj,jk-1) *zd_lw(ji,jj,jk-1)
             END DO
             ! Third recurrence : Ek = ( Lk - Uk * Ek+1 ) / Dk
-            en(ji,jj,jpkm1) = zd_lw(ji,jj,jpkm1) / zdiag(ji,jj,jpkm1)
-            DO jk = jpk-2, 2, -1
+            en(ji,jj,jpkfm1) = zd_lw(ji,jj,jpkfm1) / zdiag(ji,jj,jpkfm1)
+            DO jk = jpkf-2, 2, -1
                en(ji,jj,jk) = ( zd_lw(ji,jj,jk) - zd_up(ji,jj,jk) * en(ji,jj,jk+1) ) / zdiag(ji,jj,jk)
             END DO
-            DO jk = 2, jpkm1                       ! set the minimum value of tke
+            DO jk = 2, jpkfm1                       ! set the minimum value of tke
                en(ji,jj,jk) = MAX( en(ji,jj,jk), rn_emin ) * tmask(ji,jj,jk)
             END DO
@@ -434,5 +434,5 @@
       END DO
 #else
-      DO jk = 3, jpkm1                             ! First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1
+      DO jk = 3, jpkfm1                             ! First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1
          DO jj = 2, jpjm1
             DO ji = fs_2, fs_jpim1    ! vector opt.
@@ -446,5 +446,5 @@
          END DO
       END DO
-      DO jk = 3, jpkm1
+      DO jk = 3, jpkfm1
          DO jj = 2, jpjm1
             DO ji = fs_2, fs_jpim1    ! vector opt.
@@ -455,8 +455,8 @@
       DO jj = 2, jpjm1                             ! thrid recurrence : Ek = ( Lk - Uk * Ek+1 ) / Dk
          DO ji = fs_2, fs_jpim1    ! vector opt.
-            en(ji,jj,jpkm1) = zd_lw(ji,jj,jpkm1) / zdiag(ji,jj,jpkm1)
-         END DO
-      END DO
-      DO jk = jpk-2, 2, -1
+            en(ji,jj,jpkfm1) = zd_lw(ji,jj,jpkfm1) / zdiag(ji,jj,jpkfm1)
+         END DO
+      END DO
+      DO jk = jpkf-2, 2, -1
          DO jj = 2, jpjm1
             DO ji = fs_2, fs_jpim1    ! vector opt.
@@ -465,5 +465,5 @@
          END DO
       END DO
-      DO jk = 2, jpkm1                             ! set the minimum value of tke
+      DO jk = 2, jpkfm1                             ! set the minimum value of tke
          DO jj = 2, jpjm1
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -481,7 +481,7 @@
          DO jj = 2, jpjm1
             DO ji = 2, jpim1
-               DO jk = 2, jpkm1
-#else
-         DO jk = 2, jpkm1
+               DO jk = 2, jpkfm1
+#else
+         DO jk = 2, jpkfm1
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -505,5 +505,5 @@
          !!            unless we also make zdif a 2-d (jpi,jpj) array
 !CDIR NOVERRCHK
-         DO jk = 2, jpkm1
+         DO jk = 2, jpkfm1
 !CDIR NOVERRCHK
             DO jj = 2, jpjm1
@@ -574,4 +574,5 @@
       !!----------------------------------------------------------------------
       USE oce, ONLY:   zmpdl => ua , zmxlm => va , zmxld => ta   ! (ua,va,ta) used as workspace
+      USE arpdebugging, ONLY: dump_array
       !! DCSE_NEMO: need additional directives for renamed module variables
 !FTRANS zmpdl zmxlm zmxld :I :I :z
@@ -581,4 +582,5 @@
       REAL(wp) ::   zdku, zpdlr, zri, zsqen     !   -      -
       REAL(wp) ::   zdkv, zemxl, zemlm, zemlp   !   -      -
+      INTEGER, PARAMETER :: DUMP_LEVEL = 26 ! ARPDBG - level to dump to disk
       !!--------------------------------------------------------------------
 
@@ -599,11 +601,11 @@
       DO jj = 2, jpjm1
          DO ji = 2, jpim1
-            zmxlm(ji,jj,jpk) = rmxl_min     ! last level set to the interior minium value
-            DO jk = 2, jpkm1        ! interior value : l=sqrt(2*e/n^2)
-#else
-      zmxlm(:,:,jpk)  = rmxl_min    ! last level set to the interior minium value
-      !
-!CDIR NOVERRCHK
-      DO jk = 2, jpkm1              ! interior value : l=sqrt(2*e/n^2)
+            zmxlm(ji,jj,jpkf) = rmxl_min     ! last level set to the interior minium value
+            DO jk = 2, jpkfm1        ! interior value : l=sqrt(2*e/n^2)
+#else
+      zmxlm(:,:,jpkf)  = rmxl_min    ! last level set to the interior minium value
+      !
+!CDIR NOVERRCHK
+      DO jk = 2, jpkfm1              ! interior value : l=sqrt(2*e/n^2)
 !CDIR NOVERRCHK
          DO jj = 2, jpjm1
@@ -619,6 +621,6 @@
       !                     !* Physical limits for the mixing length
       !
-      zmxld(:,:, 1 ) = zmxlm(:,:,1)   ! surface set to the zmxlm   value
-      zmxld(:,:,jpk) = rmxl_min       ! last level  set to the minimum value
+      zmxld(:,:, 1  ) = zmxlm(:,:,1)   ! surface set to the zmxlm   value
+      zmxld(:,:,jpkf) = rmxl_min       ! last level  set to the minimum value
       !
       SELECT CASE ( nn_mxl )
@@ -628,7 +630,7 @@
          DO jj = 2, jpjm1
             DO ji = 2, jpim1
-               DO jk = 2, jpkm1
-#else
-         DO jk = 2, jpkm1
+               DO jk = 2, jpkfm1
+#else
+         DO jk = 2, jpkfm1
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -646,7 +648,7 @@
          DO jj = 2, jpjm1
             DO ji = 2, jpim1
-               DO jk = 2, jpkm1
-#else
-         DO jk = 2, jpkm1
+               DO jk = 2, jpkfm1
+#else
+         DO jk = 2, jpkfm1
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -663,8 +665,8 @@
          DO jj = 2, jpjm1
             DO ji = 2, jpim1
-               DO jk = 2, jpkm1   ! from the surface to the bottom :
+               DO jk = 2, jpkfm1   ! from the surface to the bottom :
                   zmxlm(ji,jj,jk) = MIN( zmxlm(ji,jj,jk-1) + fse3t(ji,jj,jk-1), zmxlm(ji,jj,jk) )
                END DO
-               DO jk = jpkm1, 2, -1     ! from the bottom to the surface :
+               DO jk = jpkfm1, 2, -1     ! from the bottom to the surface :
                   zemxl = MIN( zmxlm(ji,jj,jk+1) + fse3t(ji,jj,jk+1), zmxlm(ji,jj,jk) )
                   zmxlm(ji,jj,jk) = zemxl
@@ -674,5 +676,5 @@
          END DO
 #else
-         DO jk = 2, jpkm1         ! from the surface to the bottom :
+         DO jk = 2, jpkfm1         ! from the surface to the bottom :
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -681,5 +683,5 @@
             END DO
          END DO
-         DO jk = jpkm1, 2, -1     ! from the bottom to the surface :
+         DO jk = jpkfm1, 2, -1     ! from the bottom to the surface :
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -696,11 +698,11 @@
          DO jj = 2, jpjm1
             DO ji = 2, jpim1
-               DO jk = 2, jpkm1         ! from the surface to the bottom : lup
+               DO jk = 2, jpkfm1         ! from the surface to the bottom : lup
                   zmxld(ji,jj,jk) = MIN( zmxld(ji,jj,jk-1) + fse3t(ji,jj,jk-1), zmxlm(ji,jj,jk) )
                END DO
-               DO jk = jpkm1, 2, -1     ! from the bottom to the surface : ldown
+               DO jk = jpkfm1, 2, -1     ! from the bottom to the surface : ldown
                   zmxlm(ji,jj,jk) = MIN( zmxlm(ji,jj,jk+1) + fse3t(ji,jj,jk+1), zmxlm(ji,jj,jk) )
                END DO
-               DO jk = 2, jpkm1
+               DO jk = 2, jpkfm1
                   zemlm = MIN ( zmxld(ji,jj,jk),  zmxlm(ji,jj,jk) )
                   zemlp = SQRT( zmxld(ji,jj,jk) * zmxlm(ji,jj,jk) )
@@ -711,5 +713,5 @@
          END DO
 #else
-         DO jk = 2, jpkm1         ! from the surface to the bottom : lup
+         DO jk = 2, jpkfm1         ! from the surface to the bottom : lup
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -718,5 +720,5 @@
             END DO
          END DO
-         DO jk = jpkm1, 2, -1     ! from the bottom to the surface : ldown
+         DO jk = jpkfm1, 2, -1     ! from the bottom to the surface : ldown
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -726,5 +728,5 @@
          END DO
 !CDIR NOVERRCHK
-         DO jk = 2, jpkm1
+         DO jk = 2, jpkfm1
 !CDIR NOVERRCHK
             DO jj = 2, jpjm1
@@ -753,8 +755,8 @@
       DO jj = 2, jpjm1
          DO ji = 2, jpim1
-            DO jk = 1, jpkm1      !* vertical eddy viscosity & diffivity at w-points
-#else
-!CDIR NOVERRCHK
-      DO jk = 1, jpkm1            !* vertical eddy viscosity & diffivity at w-points
+            DO jk = 1, jpkfm1      !* vertical eddy viscosity & diffivity at w-points
+#else
+!CDIR NOVERRCHK
+      DO jk = 1, jpkfm1            !* vertical eddy viscosity & diffivity at w-points
 !CDIR NOVERRCHK
          DO jj = 2, jpjm1
@@ -775,7 +777,7 @@
       DO jj = 2, jpjm1
          DO ji = 2, jpim1
-            DO jk = 2, jpkm1      !* vertical eddy viscosity at u- and v-points
-#else
-      DO jk = 2, jpkm1            !* vertical eddy viscosity at u- and v-points
+            DO jk = 2, jpkfm1      !* vertical eddy viscosity at u- and v-points
+#else
+      DO jk = 2, jpkfm1            !* vertical eddy viscosity at u- and v-points
          DO jj = 2, jpjm1
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -786,13 +788,17 @@
          END DO
       END DO
+
       CALL lbc_lnk( avmu, 'U', 1. )   ;   CALL lbc_lnk( avmv, 'V', 1. )      ! Lateral boundary conditions
       !
+      CALL dump_array(1, 'avmu_tke',avmu(:,:,DUMP_LEVEL), &
+                      withHalos=.TRUE.,atStep=1)
+
       IF( nn_pdl == 1 ) THEN      !* Prandtl number case: update avt
 #if defined key_z_first
          DO jj = 2, jpjm1
             DO ji = 2, jpim1
-               DO jk = 2, jpkm1
-#else
-         DO jk = 2, jpkm1
+               DO jk = 2, jpkfm1
+#else
+         DO jk = 2, jpkfm1
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -920,8 +926,8 @@
       DO jj = 1, jpj
          DO ji = 1, jpi
-             avt (ji,jj,:) = avtb(:) * tmask(ji,jj,:)
-             avm (ji,jj,:) = avmb(:) * tmask(ji,jj,:)
-             avmu(ji,jj,:) = avmb(:) * umask(ji,jj,:)
-             avmv(ji,jj,:) = avmb(:) * vmask(ji,jj,:)
+             avt (ji,jj,:jpk) = avtb(:jpk) * tmask(ji,jj,:jpk)
+             avm (ji,jj,:jpk) = avmb(:jpk) * tmask(ji,jj,:jpk)
+             avmu(ji,jj,:jpk) = avmb(:jpk) * umask(ji,jj,:jpk)
+             avmv(ji,jj,:jpk) = avmb(:jpk) * vmask(ji,jj,:jpk)
          END DO
       END DO
