Index: NEMO/trunk/src/OCE/ZDF/zdf_oce.F90
===================================================================
--- NEMO/trunk/src/OCE/ZDF/zdf_oce.F90	(revision 10364)
+++ NEMO/trunk/src/OCE/ZDF/zdf_oce.F90	(revision 10425)
@@ -69,5 +69,5 @@
          &      avmb(jpk)         , avtb(jpk)          , avtb_2d(jpi,jpj) , STAT = zdf_oce_alloc )
          !
-      IF( zdf_oce_alloc /= 0 )   CALL ctl_warn('zdf_oce_alloc: failed to allocate arrays')
+      IF( zdf_oce_alloc /= 0 )   CALL ctl_stop( 'STOP', 'zdf_oce_alloc: failed to allocate arrays' )
       !
    END FUNCTION zdf_oce_alloc
Index: NEMO/trunk/src/OCE/ZDF/zdfgls.F90
===================================================================
--- NEMO/trunk/src/OCE/ZDF/zdfgls.F90	(revision 10364)
+++ NEMO/trunk/src/OCE/ZDF/zdfgls.F90	(revision 10425)
@@ -119,6 +119,6 @@
          &      zwall (jpi,jpj,jpk) , ustar2_top (jpi,jpj) , ustar2_bot(jpi,jpj) , STAT= zdf_gls_alloc )
          !
-      IF( lk_mpp             )   CALL mpp_sum ( zdf_gls_alloc )
-      IF( zdf_gls_alloc /= 0 )   CALL ctl_warn('zdf_gls_alloc: failed to allocate arrays')
+      CALL mpp_sum ( 'zdfgls', zdf_gls_alloc )
+      IF( zdf_gls_alloc /= 0 )   CALL ctl_stop( 'STOP', 'zdf_gls_alloc: failed to allocate arrays' )
    END FUNCTION zdf_gls_alloc
 
Index: NEMO/trunk/src/OCE/ZDF/zdfiwm.F90
===================================================================
--- NEMO/trunk/src/OCE/ZDF/zdfiwm.F90	(revision 10364)
+++ NEMO/trunk/src/OCE/ZDF/zdfiwm.F90	(revision 10425)
@@ -64,6 +64,6 @@
       &         hbot_iwm(jpi,jpj),  hcri_iwm(jpi,jpj)                     , STAT=zdf_iwm_alloc )
       !
-      IF( lk_mpp             )   CALL mpp_sum ( zdf_iwm_alloc )
-      IF( zdf_iwm_alloc /= 0 )   CALL ctl_warn('zdf_iwm_alloc: failed to allocate arrays')
+      CALL mpp_sum ( 'zdfiwm', zdf_iwm_alloc )
+      IF( zdf_iwm_alloc /= 0 )   CALL ctl_stop( 'STOP', 'zdf_iwm_alloc: failed to allocate arrays' )
    END FUNCTION zdf_iwm_alloc
 
@@ -122,5 +122,5 @@
       !
       INTEGER  ::   ji, jj, jk   ! dummy loop indices
-      REAL(wp) ::   zztmp        ! scalar workspace
+      REAL(wp) ::   zztmp, ztmp1, ztmp2        ! scalar workspace
       REAL(wp), DIMENSION(jpi,jpj)     ::   zfact       ! Used for vertical structure
       REAL(wp), DIMENSION(jpi,jpj)     ::   zhdep       ! Ocean depth
@@ -157,8 +157,15 @@
 !!gm gde3w ==>>>  check for ssh taken into account.... seem OK gde3w_n=gdept_n - sshn
       DO jk = 2, jpkm1              ! complete with the level-dependent part
-         zemx_iwm(:,:,jk) = zfact(:,:) * (  EXP( ( gde3w_n(:,:,jk  ) - zhdep(:,:) ) / hcri_iwm(:,:) )                      &
-            &                             - EXP( ( gde3w_n(:,:,jk-1) - zhdep(:,:) ) / hcri_iwm(:,:) )  ) * wmask(:,:,jk)   &
-            &                          / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) )
-
+         DO jj = 1, jpj             
+            DO ji = 1, jpi
+               IF ( zfact(ji,jj) == 0._wp .OR. wmask(ji,jj,jk) == 0._wp ) THEN   ! optimization
+                  zemx_iwm(ji,jj,jk) = 0._wp
+               ELSE
+                  zemx_iwm(ji,jj,jk) = zfact(ji,jj) * (  EXP( ( gde3w_n(ji,jj,jk  ) - zhdep(ji,jj) ) / hcri_iwm(ji,jj) )     &
+                       &                               - EXP( ( gde3w_n(ji,jj,jk-1) - zhdep(ji,jj) ) / hcri_iwm(ji,jj) ) )   &
+                       &                            / ( gde3w_n(ji,jj,jk) - gde3w_n(ji,jj,jk-1) )
+               ENDIF
+            END DO
+         END DO
 !!gm delta(gde3w_n) = e3t_n  !!  Please verify the grid-point position w versus t-point
 !!gm it seems to me that only 1/hcri_iwm  is used ==>  compute it one for all
@@ -234,8 +241,15 @@
       zwkb(:,:,1) = zhdep(:,:) * wmask(:,:,1)
       !
-      zweight(:,:,:) = 0._wp
-      DO jk = 2, jpkm1
-         zweight(:,:,jk) = MAX( 0._wp, rn2(:,:,jk) ) * hbot_iwm(:,:) * wmask(:,:,jk)                    &
-            &   * (  EXP( -zwkb(:,:,jk) / hbot_iwm(:,:) ) - EXP( -zwkb(:,:,jk-1) / hbot_iwm(:,:) )  )
+      DO jk = 2, jpkm1
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF ( rn2(ji,jj,jk) <= 0._wp .OR. wmask(ji,jj,jk) == 0._wp ) THEN   ! optimization
+                  zweight(ji,jj,jk) = 0._wp
+               ELSE
+                  zweight(ji,jj,jk) = rn2(ji,jj,jk) * hbot_iwm(ji,jj)    &
+                     &   * (  EXP( -zwkb(ji,jj,jk) / hbot_iwm(ji,jj) ) - EXP( -zwkb(ji,jj,jk-1) / hbot_iwm(ji,jj) )  )
+               ENDIF
+            END DO
+         END DO
       END DO
       !
@@ -305,5 +319,5 @@
             END DO
          END DO
-         IF( lk_mpp )   CALL mpp_sum( zztmp )
+         CALL mpp_sum( 'zdfiwm', zztmp )
          zztmp = rau0 * zztmp ! Global integral of rauo * Kz * N^2 = power contributing to mixing 
          !
@@ -322,10 +336,14 @@
       !      
       IF( ln_tsdiff ) THEN          !* Option for differential mixing of salinity and temperature
+         ztmp1 = 0.505_wp + 0.495_wp * TANH( 0.92_wp * ( LOG10( 1.e-20_wp ) - 0.60_wp ) )
          DO jk = 2, jpkm1              ! Calculate S/T diffusivity ratio as a function of Reb
             DO jj = 1, jpj
                DO ji = 1, jpi
-                  zav_ratio(ji,jj,jk) = ( 0.505_wp + 0.495_wp *                                                                  &
-                      &   TANH(    0.92_wp * (   LOG10(  MAX( 1.e-20_wp, zReb(ji,jj,jk) * 5._wp * r1_6 )  ) - 0.60_wp   )    )   &
-                      &                 ) * wmask(ji,jj,jk)
+                  ztmp2 = zReb(ji,jj,jk) * 5._wp * r1_6
+                  IF ( ztmp2 > 1.e-20_wp .AND. wmask(ji,jj,jk) == 1._wp ) THEN
+                     zav_ratio(ji,jj,jk) = 0.505_wp + 0.495_wp * TANH( 0.92_wp * ( LOG10(ztmp2) - 0.60_wp ) )
+                  ELSE
+                     zav_ratio(ji,jj,jk) = ztmp1 * wmask(ji,jj,jk)
+                  ENDIF
                END DO
             END DO
@@ -463,7 +481,7 @@
       ecri_iwm(:,:) = ecri_iwm(:,:) * ssmask(:,:)
 
-      zbot = glob_sum( e1e2t(:,:) * ebot_iwm(:,:) )
-      zpyc = glob_sum( e1e2t(:,:) * epyc_iwm(:,:) )
-      zcri = glob_sum( e1e2t(:,:) * ecri_iwm(:,:) )
+      zbot = glob_sum( 'zdfiwm', e1e2t(:,:) * ebot_iwm(:,:) )
+      zpyc = glob_sum( 'zdfiwm', e1e2t(:,:) * epyc_iwm(:,:) )
+      zcri = glob_sum( 'zdfiwm', e1e2t(:,:) * ecri_iwm(:,:) )
       IF(lwp) THEN
          WRITE(numout,*) '      High-mode wave-breaking energy:             ', zbot * 1.e-12_wp, 'TW'
Index: NEMO/trunk/src/OCE/ZDF/zdfmxl.F90
===================================================================
--- NEMO/trunk/src/OCE/ZDF/zdfmxl.F90	(revision 10364)
+++ NEMO/trunk/src/OCE/ZDF/zdfmxl.F90	(revision 10425)
@@ -50,6 +50,6 @@
          ALLOCATE( nmln(jpi,jpj), hmld(jpi,jpj), hmlp(jpi,jpj), hmlpt(jpi,jpj), STAT= zdf_mxl_alloc )
          !
-         IF( lk_mpp             )   CALL mpp_sum ( zdf_mxl_alloc )
-         IF( zdf_mxl_alloc /= 0 )   CALL ctl_warn('zdf_mxl_alloc: failed to allocate arrays.')
+         CALL mpp_sum ( 'zdfmxl', zdf_mxl_alloc )
+         IF( zdf_mxl_alloc /= 0 )   CALL ctl_stop( 'STOP', 'zdf_mxl_alloc: failed to allocate arrays.' )
          !
       ENDIF
Index: NEMO/trunk/src/OCE/ZDF/zdfosm.F90
===================================================================
--- NEMO/trunk/src/OCE/ZDF/zdfosm.F90	(revision 10364)
+++ NEMO/trunk/src/OCE/ZDF/zdfosm.F90	(revision 10425)
@@ -118,5 +118,5 @@
           &   etmean(jpi,jpj,jpk),  STAT= zdf_osm_alloc )
      IF( zdf_osm_alloc /= 0 )   CALL ctl_warn('zdf_osm_alloc: failed to allocate zdf_osm arrays')
-     IF( lk_mpp             )   CALL mpp_sum ( zdf_osm_alloc )
+     CALL mpp_sum ( 'zdfosm', zdf_osm_alloc )
    END FUNCTION zdf_osm_alloc
 
@@ -1287,5 +1287,5 @@
 
        ! Lateral boundary conditions on zvicos (sign unchanged), needed to caclulate viscosities on u and v grids
-       CALL lbc_lnk( zviscos(:,:,:), 'W', 1. )
+       CALL lbc_lnk( 'zdfosm', zviscos(:,:,:), 'W', 1. )
 
        ! GN 25/8: need to change tmask --> wmask
@@ -1300,5 +1300,5 @@
      END DO
       ! Lateral boundary conditions on ghamu and ghamv, currently on W-grid  (sign unchanged), needed to caclulate gham[uv] on u and v grids
-     CALL lbc_lnk_multi( p_avt, 'W', 1. , p_avm, 'W', 1.,   &
+     CALL lbc_lnk_multi( 'zdfosm', p_avt, 'W', 1. , p_avm, 'W', 1.,   &
       &                  ghamu, 'W', 1. , ghamv, 'W', 1. )
        DO jk = 2, jpkm1
@@ -1318,5 +1318,5 @@
         ! Lateral boundary conditions on final outputs for gham[ts],  on W-grid  (sign unchanged)
         ! Lateral boundary conditions on final outputs for gham[uv],  on [UV]-grid  (sign unchanged)
-        CALL lbc_lnk_multi( ghamt, 'W', 1. , ghams, 'W', 1.,   &
+        CALL lbc_lnk_multi( 'zdfosm', ghamt, 'W', 1. , ghams, 'W', 1.,   &
          &                  ghamu, 'U', 1. , ghamv, 'V', 1. )
 
@@ -1359,5 +1359,5 @@
       END IF
       ! Lateral boundary conditions on p_avt  (sign unchanged)
-      CALL lbc_lnk( p_avt(:,:,:), 'W', 1. )
+      CALL lbc_lnk( 'zdfosm', p_avt(:,:,:), 'W', 1. )
       !
    END SUBROUTINE zdf_osm
Index: NEMO/trunk/src/OCE/ZDF/zdfphy.F90
===================================================================
--- NEMO/trunk/src/OCE/ZDF/zdfphy.F90	(revision 10364)
+++ NEMO/trunk/src/OCE/ZDF/zdfphy.F90	(revision 10425)
@@ -300,13 +300,13 @@
       !                                         !* Lateral boundary conditions (sign unchanged)
       IF( l_zdfsh2 ) THEN
-         CALL lbc_lnk_multi( avm_k, 'W', 1. , avt_k, 'W', 1.,   &
+         CALL lbc_lnk_multi( 'zdfphy', avm_k, 'W', 1. , avt_k, 'W', 1.,   &
             &                avm  , 'W', 1. , avt  , 'W', 1. , avs , 'W', 1. )
       ELSE
-         CALL lbc_lnk_multi( avm  , 'W', 1. , avt  , 'W', 1. , avs , 'W', 1. )
+         CALL lbc_lnk_multi( 'zdfphy', avm  , 'W', 1. , avt  , 'W', 1. , avs , 'W', 1. )
       ENDIF
       !
       IF( l_zdfdrg ) THEN     ! drag  have been updated (non-linear cases)
-         IF( ln_isfcav ) THEN   ;  CALL lbc_lnk_multi( rCdU_top, 'T', 1. , rCdU_bot, 'T', 1. )   ! top & bot drag
-         ELSE                   ;  CALL lbc_lnk      ( rCdU_bot, 'T', 1. )                       ! bottom drag only
+         IF( ln_isfcav ) THEN   ;  CALL lbc_lnk_multi( 'zdfphy', rCdU_top, 'T', 1. , rCdU_bot, 'T', 1. )   ! top & bot drag
+         ELSE                   ;  CALL lbc_lnk      ( 'zdfphy', rCdU_bot, 'T', 1. )                       ! bottom drag only
          ENDIF
       ENDIF
@@ -331,5 +331,5 @@
      ALLOCATE(     wi(jpi,jpj,jpk), Cu_adv(jpi,jpj,jpk),  STAT= zdf_phy_alloc )
      IF( zdf_phy_alloc /= 0 )   CALL ctl_warn('zdf_phy_alloc: failed to allocate ln_zad_Aimp=T required arrays')
-     IF( lk_mpp             )   CALL mpp_sum ( zdf_phy_alloc )
+     CALL mpp_sum ( 'zdfphy', zdf_phy_alloc )
    END FUNCTION zdf_phy_alloc
 
Index: NEMO/trunk/src/OCE/ZDF/zdftke.F90
===================================================================
--- NEMO/trunk/src/OCE/ZDF/zdftke.F90	(revision 10364)
+++ NEMO/trunk/src/OCE/ZDF/zdftke.F90	(revision 10425)
@@ -103,6 +103,6 @@
       ALLOCATE( htau(jpi,jpj) , dissl(jpi,jpj,jpk) , apdlr(jpi,jpj,jpk) ,   STAT= zdf_tke_alloc )
       !
-      IF( lk_mpp             )   CALL mpp_sum ( zdf_tke_alloc )
-      IF( zdf_tke_alloc /= 0 )   CALL ctl_warn('zdf_tke_alloc: failed to allocate arrays')
+      CALL mpp_sum ( 'zdftke', zdf_tke_alloc )
+      IF( zdf_tke_alloc /= 0 )   CALL ctl_stop( 'STOP', 'zdf_tke_alloc: failed to allocate arrays' )
       !
    END FUNCTION zdf_tke_alloc
@@ -202,5 +202,5 @@
       REAL(wp) ::   zzd_up, zzd_lw             !   -         -
       INTEGER , DIMENSION(jpi,jpj)     ::   imlc
-      REAL(wp), DIMENSION(jpi,jpj)     ::   zhlc
+      REAL(wp), DIMENSION(jpi,jpj)     ::   zhlc, zfr_i
       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zpelc, zdiag, zd_up, zd_lw
       !!--------------------------------------------------------------------
@@ -290,14 +290,23 @@
          END DO
          zcof = 0.016 / SQRT( zrhoa * zcdrag )
+         DO jj = 2, jpjm1
+            DO ji = fs_2, fs_jpim1   ! vector opt.
+               zus  = zcof * SQRT( taum(ji,jj) )           ! Stokes drift
+               zfr_i(ji,jj) = ( 1._wp - 4._wp * fr_i(ji,jj) ) * zus * zus * zus * tmask(ji,jj,1) ! zus > 0. ok
+               IF (zfr_i(ji,jj) < 0. ) zfr_i(ji,jj) = 0.
+            END DO
+         END DO         
          DO jk = 2, jpkm1         !* TKE Langmuir circulation source term added to en
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
-                  zus  = zcof * SQRT( taum(ji,jj) )           ! Stokes drift
-                  !                                           ! vertical velocity due to LC
-                  zind = 0.5 - SIGN( 0.5, pdepw(ji,jj,jk) - zhlc(ji,jj) )
-                  zwlc = zind * rn_lc * zus * SIN( rpi * pdepw(ji,jj,jk) / zhlc(ji,jj) )
-                  !                                           ! TKE Langmuir circulation source term
-                  en(ji,jj,jk) = en(ji,jj,jk) + rdt * MAX(0.,1._wp - 4.*fr_i(ji,jj) ) * ( zwlc * zwlc * zwlc )   &
-                     &                              / zhlc(ji,jj) * wmask(ji,jj,jk) * tmask(ji,jj,1)
+                  IF ( zfr_i(ji,jj) /= 0. ) THEN               
+                     ! vertical velocity due to LC   
+                     IF ( pdepw(ji,jj,jk) - zhlc(ji,jj) < 0 .AND. wmask(ji,jj,jk) /= 0. ) THEN
+                        !                                           ! vertical velocity due to LC
+                        zwlc = rn_lc * SIN( rpi * pdepw(ji,jj,jk) / zhlc(ji,jj) )   ! warning: optimization: zus^3 is in zfr_i
+                        !                                           ! TKE Langmuir circulation source term
+                        en(ji,jj,jk) = en(ji,jj,jk) + rdt * zfr_i(ji,jj) * ( zwlc * zwlc * zwlc ) / zhlc(ji,jj)
+                     ENDIF
+                  ENDIF
                END DO
             END DO
