Index: NEMO/branches/UKMO/tools_r4.0-HEAD_dev_MEs/DOMAINcfg/namelist_cfg
===================================================================
--- NEMO/branches/UKMO/tools_r4.0-HEAD_dev_MEs/DOMAINcfg/namelist_cfg	(revision 15129)
+++ NEMO/branches/UKMO/tools_r4.0-HEAD_dev_MEs/DOMAINcfg/namelist_cfg	(revision 15131)
@@ -40,5 +40,5 @@
 &namzgr_mes    !   MEs-coordinate                     
 !-----------------------------------------------------------------------
-   ln_envl     =   .TRUE. , .TRUE. , .TRUE. , .TRUE., .TRUE.   ! (T/F) If the envelope is used
+   ln_envl     =   .TRUE. , .TRUE. , .TRUE. , .TRUE., .FALSE.  ! (T/F) If the envelope is used
    nn_strt     =     1    ,   1    ,    1   ,   1    ,   1     ! Stretch. funct.: Madec 1996 (0) or
                                                                ! Song & Haidvogel 1994 (1) or                    
Index: NEMO/branches/UKMO/tools_r4.0-HEAD_dev_MEs/DOMAINcfg/src/lib_mpp.f90
===================================================================
--- NEMO/branches/UKMO/tools_r4.0-HEAD_dev_MEs/DOMAINcfg/src/lib_mpp.f90	(revision 15129)
+++ NEMO/branches/UKMO/tools_r4.0-HEAD_dev_MEs/DOMAINcfg/src/lib_mpp.f90	(revision 15131)
@@ -76,5 +76,5 @@
    PUBLIC   mpp_lnk_2d_9 , mpp_lnk_2d_multiple 
    PUBLIC   mpp_lnk_sum_3d, mpp_lnk_sum_2d
-   PUBLIC   mppscatter, mppgather
+   PUBLIC   mppscatter, mppgather, mppbcast_a_real
    PUBLIC   mpp_ini_ice, mpp_ini_znl
    PUBLIC   mppsize
@@ -1847,4 +1847,27 @@
 
 
+   SUBROUTINE mppbcast_a_real( kvals, kno, kroot )
+      !!----------------------------------------------------------------------
+      !!                  ***  routine mppbcast_a_real  ***
+      !!
+      !! ** Purpose : Send array kvals to all processors
+      !!
+      !! ** Method  : MPI broadcast
+      !!
+      !!-----------------------------------------------------------------------
+      INTEGER                 , INTENT(in   ) :: kno     ! Number of elements in array
+      INTEGER                 , INTENT(in   ) :: kroot   ! Processor to send data
+      REAL(wp), DIMENSION(kno), INTENT(inout) :: kvals   ! Array to send on kroot, receive for non-kroot
+      !!
+      INTEGER                   ::   ierr    ! temporary integer
+      INTEGER                   ::   localcomm
+      !!-----------------------------------------------------------------------
+      !
+      localcomm = mpi_comm_opa
+      CALL mpi_bcast( kvals, kno, mpi_double_precision, kroot, localcomm, ierr )
+      !
+   END SUBROUTINE mppbcast_a_real
+
+
    SUBROUTINE mppmax_a_int( ktab, kdim, kcom )
       !!----------------------------------------------------------------------
Index: NEMO/branches/UKMO/tools_r4.0-HEAD_dev_MEs/DOMAINcfg/src/mes.F90
===================================================================
--- NEMO/branches/UKMO/tools_r4.0-HEAD_dev_MEs/DOMAINcfg/src/mes.F90	(revision 15129)
+++ NEMO/branches/UKMO/tools_r4.0-HEAD_dev_MEs/DOMAINcfg/src/mes.F90	(revision 15131)
@@ -35,5 +35,7 @@
    INTEGER            :: nn_strt(max_nn_env)     ! Array specifing the stretching function for each envelope:
                                                  ! Madec 1996 (0), Song and Haidvogel 1994 (1) 
-   REAL(wp)           :: rn_ebot_max(max_nn_env) ! maximum depth of envelopes (= envelopes depths) (>0) (m)
+   REAL(wp)           :: rn_ebot_max(max_nn_env) ! depth of envelopes in one point (for checking 1D transform.)
+   REAL(wp)           :: max_dep_env(max_nn_env) ! global maximum depth of envelopes
+   REAL(wp)           :: min_dep_env(max_nn_env) ! global minimum depth of envelopes
    INTEGER            :: nn_slev(max_nn_env)     ! Array specifing number of levels of each enveloped vertical zone
    REAL(wp)           :: rn_e_hc(max_nn_env)     ! Array specifing critical depth for transition to stretched 
@@ -161,6 +163,4 @@
       REAL(wp), POINTER, DIMENSION(:,:  ) :: env0, env1, env2, env3     ! for loops over envelopes
                                                                         ! for cubic splines
-      INTEGER                             :: gst_envl(max_nn_env)       ! Array to deal with a ghost last envelope 
-
       !!----------------------------------------------------------------------
       !
@@ -171,8 +171,4 @@
       CALL wrk_alloc( jpi, jpj, env_up, env_dw, env0, env1, env2, env3 )
       CALL wrk_alloc( jpk, z_gsigw1, z_gsigt1, z_esigw1, z_esigt1 )
-      !
-      ! Determining if there is a "ghost" envelope: 
-      gst_envl(:) = 0
-      IF ( nn_slev(tot_env) == 0 ) gst_envl(tot_env) = 1
       !
       ! Initializing to 0.0 some arrays:
@@ -227,35 +223,32 @@
 
          IF ( isodd(je) ) THEN
-            IF ( gst_envl(je) == 0 ) THEN
-               DO jk = 1, num_s
-                  ind = s_1st+jk-1
-                  ! W-GRID
-                  z_gsigw1(ind) = -s_coord(sigma(1,jk,num_s), coefa, &
-                                           coefb, alpha, num_s, nn_strt(je))
-                  ! T-GRID
-                  z_gsigt1(ind) = -s_coord(sigma(0,jk,num_s), coefa, &
-                                           coefb, alpha, num_s, nn_strt(je))
-                  zcoefw = z_gsigw1(ind)
-                  zcoeft = z_gsigt1(ind)
-                              
-                  IF (nn_strt(je) /= 2) THEN
-                     gdept_1d(ind) = z_mes(dep_up, dep_dw, zcoeft, &
-                                           -sigma(0,jk,num_s), env_hc)
-                     gdepw_1d(ind) = z_mes(dep_up, dep_dw, zcoefw, &
-                                           -sigma(1,jk,num_s), env_hc)
-                  ELSE
-                     gdept_1d(ind) = z_mes(dep_up, dep_dw, zcoeft, &
-                                           -sigma(0,jk,num_s), 0._wp)
-                     gdepw_1d(ind) = z_mes(dep_up, dep_dw, zcoefw, &
-                                           -sigma(1,jk,num_s), 0._wp)
-                  END IF
-               ENDDO
-            ENDIF
-
+            DO jk = 1, num_s
+               ind = s_1st+jk-1
+               ! W-GRID
+               z_gsigw1(ind) = -s_coord(sigma(1,jk,num_s), coefa, &
+                                        coefb, alpha, num_s, nn_strt(je))
+               ! T-GRID
+               z_gsigt1(ind) = -s_coord(sigma(0,jk,num_s), coefa, &
+                                        coefb, alpha, num_s, nn_strt(je))
+               zcoefw = z_gsigw1(ind)
+               zcoeft = z_gsigt1(ind)
+               !                  
+               IF (nn_strt(je) /= 2) THEN
+                   gdept_1d(ind) = z_mes(dep_up, dep_dw, zcoeft, &
+                                        -sigma(0,jk,num_s), env_hc)
+                   gdepw_1d(ind) = z_mes(dep_up, dep_dw, zcoefw, &
+                                        -sigma(1,jk,num_s), env_hc)
+               ELSE
+                   gdept_1d(ind) = z_mes(dep_up, dep_dw, zcoeft, &
+                                        -sigma(0,jk,num_s), 0._wp)
+                   gdepw_1d(ind) = z_mes(dep_up, dep_dw, zcoefw, &
+                                        -sigma(1,jk,num_s), 0._wp)
+               END IF
+            ENDDO
+            !
          ELSE
-
             ! Ensure that for splines we do not use SF12
             !IF (nn_strt(je) > 1) nn_strt(je) = 1
-
+            !
             ! Interval's endpoints TAU are 0 and 1.
             tau(1) = 0._wp
@@ -269,6 +262,6 @@
             dep1   = dep_up
             dep2   = dep_dw
-            dep3   = rn_ebot_max(je+1)
-
+            IF ( je < tot_env ) dep3 = rn_ebot_max(je+1)
+            !
             ! Computing derivative analytically
             !
@@ -297,5 +290,5 @@
             !
             ! 2. Derivative for deeper envelope
-            IF ( gst_envl(je+1) == 0 ) THEN
+            IF ( je < tot_env ) THEN
                ibcend = 1 ! Bound. cond for the 1st derivative          
                IF (nn_strt(je+1) == 2) THEN
@@ -405,5 +398,4 @@
       jk = 1
       e3w_1d(jk) = 2.0_wp * (gdept_1d(1) - gdepw_1d(1))
-
  
       ! Control print
@@ -463,62 +455,60 @@
 
          IF ( isodd(je) ) THEN
-
-            IF ( gst_envl(je) == 0 ) THEN
-               DO jk = 1, num_s
-                  ind = s_1st+jk-1
-                  DO jj = 1,jpj
-                     DO ji = 1,jpi
-                        IF (env_dw(ji,jj) < env_hc) THEN   
-                           ! shallow water, uniform sigma
-                           zcoefw = -sigma(1, jk, num_s) 
-                           zcoeft = -sigma(0, jk, num_s)
-                           !IF (nn_strt(je) == 2) THEN
-                           !   ! shallow water, z coordinates
-                           !   ! SF12 in AMM* configurations uses 
-                           !   ! this (ln_sigcrit=.false.)
-                           !   zcoefw = zcoefw*(env_hc/(env_dw(ji,jj)-env_up(ji,jj)))
-                           !   zcoeft = zcoeft*(env_hc/(env_dw(ji,jj)-env_up(ji,jj)))
-                           !ENDIF
-                        ELSE                               ! deep water, stretched s
-                           IF (nn_strt(je) == 2) THEN
-                              coefa = rn_e_th(je) /  (env_dw(ji,jj)-env_up(ji,jj))
-                              coefb = rn_e_bb(je) + ((env_dw(ji,jj)-env_up(ji,jj))*rn_e_ba(je))
-                              coefb = 1.0_wp-(coefb/(env_dw(ji,jj)-env_up(ji,jj)))
-                              ! W-GRID
-                              zcoefw = -s_coord( sigma(1,jk,num_s), coefa, &
-                                                 coefb, alpha, num_s, nn_strt(je) )
-                              ! T-GRID
-                              zcoeft = -s_coord( sigma(0,jk,num_s), coefa, &
-                                                 coefb, alpha, num_s, nn_strt(je) )
-                           ELSE
-                              zcoefw = z_gsigw1(ind)
-                              zcoeft = z_gsigt1(ind)
-                           END IF
-                        ENDIF
-
-                        IF (nn_strt(je) /= 2) THEN
-                           gdept_0(ji,jj,ind) = z_mes(env_up(ji,jj), env_dw(ji,jj), zcoeft, &
-                                                      -sigma(0, jk, num_s), env_hc)
+            !
+            DO jk = 1, num_s
+               ind = s_1st+jk-1
+               DO jj = 1,jpj
+                  DO ji = 1,jpi
+                     IF (env_dw(ji,jj) < env_hc) THEN   
+                        ! shallow water, uniform sigma
+                        zcoefw = -sigma(1, jk, num_s) 
+                        zcoeft = -sigma(0, jk, num_s)
+                        !IF (nn_strt(je) == 2) THEN
+                        !   ! shallow water, z coordinates
+                        !   ! SF12 in AMM* configurations uses 
+                        !   ! this (ln_sigcrit=.false.)
+                        !   zcoefw = zcoefw*(env_hc/(env_dw(ji,jj)-env_up(ji,jj)))
+                        !   zcoeft = zcoeft*(env_hc/(env_dw(ji,jj)-env_up(ji,jj)))
+                        !ENDIF
+                     ELSE                               ! deep water, stretched s
+                        IF (nn_strt(je) == 2) THEN
+                           coefa = rn_e_th(je) /  (env_dw(ji,jj)-env_up(ji,jj))
+                           coefb = rn_e_bb(je) + ((env_dw(ji,jj)-env_up(ji,jj))*rn_e_ba(je))
+                           coefb = 1.0_wp-(coefb/(env_dw(ji,jj)-env_up(ji,jj)))
+                           ! W-GRID
+                           zcoefw = -s_coord( sigma(1,jk,num_s), coefa, &
+                                              coefb, alpha, num_s, nn_strt(je) )
+                           ! T-GRID
+                           zcoeft = -s_coord( sigma(0,jk,num_s), coefa, &
+                                              coefb, alpha, num_s, nn_strt(je) )
+                        ELSE
+                           zcoefw = z_gsigw1(ind)
+                           zcoeft = z_gsigt1(ind)
+                        END IF
+                     ENDIF
+                     !
+                     IF (nn_strt(je) /= 2) THEN
+                        gdept_0(ji,jj,ind) = z_mes(env_up(ji,jj), env_dw(ji,jj), zcoeft, &
+                                                   -sigma(0, jk, num_s), env_hc)
                                                 
-                           gdepw_0(ji,jj,ind) = z_mes(env_up(ji,jj), env_dw(ji,jj), zcoefw, &
-                                                      -sigma(1, jk, num_s), env_hc)
-                        ELSE
-                           gdept_0(ji,jj,ind) = z_mes(env_up(ji,jj), env_dw(ji,jj), zcoeft, &
-                                                      -sigma(0, jk, num_s), 0._wp)
-
-                           gdepw_0(ji,jj,ind) = z_mes(env_up(ji,jj), env_dw(ji,jj), zcoefw, &
-                                                      -sigma(1, jk, num_s), 0._wp)
-                        END IF
-                     END DO
+                        gdepw_0(ji,jj,ind) = z_mes(env_up(ji,jj), env_dw(ji,jj), zcoefw, &
+                                                   -sigma(1, jk, num_s), env_hc)
+                     ELSE
+                        gdept_0(ji,jj,ind) = z_mes(env_up(ji,jj), env_dw(ji,jj), zcoeft, &
+                                                   -sigma(0, jk, num_s), 0._wp)
+ 
+                        gdepw_0(ji,jj,ind) = z_mes(env_up(ji,jj), env_dw(ji,jj), zcoefw, &
+                                                   -sigma(1, jk, num_s), 0._wp)
+                     END IF
                   END DO
                END DO
-            END IF
-
+            END DO
+            !
          ELSE
-
+            !
             ! Interval's endpoints TAU are 0 and 1.
             tau(1) = 0._wp
             tau(n) = 1._wp
-
+            !
             IF ( je == 2 ) THEN
                env0 = scosrf
@@ -528,6 +518,6 @@
             env1   = env_up
             env2   = env_dw
-            env3   = envlt(:,:,je+1)
-
+            IF ( je < tot_env ) env3 = envlt(:,:,je+1)
+            !            
             DO jj = 1,jpj
                DO ji = 1,jpi
@@ -564,5 +554,5 @@
                   !
                   ! 2. Derivative for lower envelope
-                  IF ( gst_envl(je+1) == 0 ) THEN
+                  IF ( je < tot_env ) THEN
                      ibcend = 1 ! Bound. cond for the 1st derivative
                      IF (nn_strt(je+1) == 2) THEN
@@ -742,14 +732,10 @@
             s_1st  = s_1st + num_s - 1
             num_s  = nn_slev(je) + 1
-            max_env_up = MAXVAL(envlt(:,:,je-1))
-            min_env_up = MINVAL(envlt(:,:,je-1))
-            IF( lk_mpp ) CALL mpp_max( max_env_up )
-            IF( lk_mpp ) CALL mpp_min( min_env_up )
+            max_env_up = max_dep_env(je-1)
+            min_env_up = min_dep_env(je-1)
          ENDIF
 
-         max_env_dw  = MAXVAL(envlt(:,:,je))
-         min_env_dw  = MINVAL(envlt(:,:,je))
-         IF( lk_mpp ) CALL mpp_max( max_env_dw )
-         IF( lk_mpp ) CALL mpp_min( min_env_dw )
+         max_env_dw  = max_dep_env(je)
+         min_env_dw  = min_dep_env(je)
 
          IF ( max_env_up == min_env_up .AND. max_env_dw == min_env_dw ) THEN
@@ -891,9 +877,10 @@
    SUBROUTINE mes_init
       CHARACTER(lc)                :: env_name   ! name of the externally defined envelope
-      INTEGER                      :: ji, jj, je, iiemax, ijemax
+      INTEGER                      :: ji, jj, je, eii, eij
       INTEGER                      :: cor_env, ios, inum
       REAL(wp)                     :: rn_bot_min ! minimum depth of ocean bottom (>0) (m)
       REAL(wp)                     :: rn_bot_max ! maximum depth of ocean bottom (= ocean depth) (>0) (m)
-      REAL(wp), DIMENSION(jpi,jpj) :: pmsk       ! for loops over envelopes
+      REAL(wp)                     :: mindep
+      REAL(wp), DIMENSION(jpi,jpj) :: pmsk, ewrk ! working arrays
 
       NAMELIST/namzgr_mes/rn_bot_min , rn_bot_max , ln_envl, &
@@ -906,5 +893,16 @@
       !
       CALL wrk_alloc( jpi, jpj, max_nn_env , envlt )
-      pmsk(:,:) = 1._wp
+      pmsk(:,:) = 0.0
+      !
+      ! Initialising some arrays
+      rn_ebot_max(:) = 0.0
+      max_dep_env(:) = 0.0
+      min_dep_env(:) = 0.0
+      !
+      IF(lwp) THEN                           ! control print
+         WRITE(numout,*) ''
+         WRITE(numout,*) 'mes_init: Setting a Multi-Envelope s-ccordinate system (Bruciaferri et al. 2018)'
+         WRITE(numout,*) '~~~~~~~~'
+      END IF
       !
       ! Namelist namzgr_mes in reference namelist : envelopes and
@@ -931,4 +929,8 @@
       WRITE(ctlmes,*) 'num. of REQUESTED env. and num. of CORRECTLY defined env. are DIFFERENT'
       IF ( tot_env /= cor_env ) CALL ctl_stop( ctlmes )
+      IF(lwp) THEN
+        WRITE(numout,*) '          Num. or requested envelopes: ', tot_env
+        WRITE(numout,*) ''
+      END IF
       !
       ! 2) Checking consistency of user defined parameters
@@ -938,4 +940,9 @@
       ! 3) Reading Bathymetry and envelopes
       IF( ntopo == 1 ) THEN
+        IF(lwp) THEN
+           WRITE(numout,*) '          Reading bathymetry and envelopes ...'
+           WRITE(numout,*) ''
+        END IF
+
         CALL iom_open ( 'bathy_meter.nc', inum )
         CALL iom_get  ( inum, jpdom_data, 'Bathymetry' , bathy, lrowattr=ln_use_jattr  )
@@ -956,12 +963,29 @@
          IF (MAXVAL(envlt(:,:,je+1)) < MAXVAL(envlt(:,:,je))) CALL ctl_stop( ctlmes )
       ENDDO
-      ! 5) Computing max depth of envelopes in the deepest
-      !    point of the domain for first check of monotonicity
+      ! 5) Computing max depth of envelopes in the shallowest
+      !    point of the deeper envelope for first check of monotonicity
       !    of transformation and computing 1D MEs-levels depths
       !    (used by diawri.F90)
-      CALL mpp_maxloc( envlt(:,:,tot_env), pmsk(:,:), rn_ebot_max(tot_env), iiemax, ijemax )
+      WHERE ( bathy > 0.0 ) pmsk = 1.0
+      ewrk(:,:) = envlt(:,:,tot_env)
+      CALL mpp_minloc( ewrk(:,:), pmsk(:,:), mindep, eii, eij )
+      IF ((mi0(eii)>1 .AND. mi0(eii)<jpi) .AND. (mj0(eij)>1 .AND. mj0(eij)<jpj)) THEN        
+         inum = mpprank
+         DO je = 1, tot_env
+            rn_ebot_max(je) = envlt(mi0(eii),mj0(eij),je)
+         END DO
+      ELSE
+         inum = -1
+      END IF
+      IF( lk_mpp ) CALL mppsync
+      IF( lk_mpp ) CALL mpp_max(inum)
+      IF( lk_mpp ) CALL mppbcast_a_real(rn_ebot_max, max_nn_env, inum )
       DO je = 1, tot_env
-         rn_ebot_max(je) = envlt(iiemax,ijemax,je)
+         max_dep_env(je) = MAXVAL(envlt(:,:,je))
+         min_dep_env(je) = MINVAL(envlt(:,:,je))
+         IF( lk_mpp ) CALL mpp_max( max_dep_env(je) )
+         IF( lk_mpp ) CALL mpp_min( min_dep_env(je) )
       END DO
+      IF( lk_mpp ) CALL mppsync
       !
       ! 6) Set maximum and minimum ocean depth
@@ -975,41 +999,48 @@
       IF(lwp) THEN                           ! control print
          WRITE(numout,*)
-         WRITE(numout,*) 'domzgr_mes : Multi Enveloped S-coordinate (Bruciaferri, Shapiro and Wobus 2017)'
-         WRITE(numout,*) '~~~~~~~~~~~'
-         WRITE(numout,*) '   Namelist namzgr_mes'
+         WRITE(numout,*) '          GEOMETRICAL FEATURES OF THE REQUESTED MEs GRID:'
+         WRITE(numout,*) '          -----------------------------------------------'
+         WRITE(numout,*) '          Minimum depth of the ocean   rn_bot_min = ', rn_bot_min
+         WRITE(numout,*) '          Maximum depth of the ocean   rn_bot_max = ', rn_bot_max
          WRITE(numout,*) ''
-         WRITE(numout,*) '   Minimum depth of the ocean   rn_bot_min, ', rn_bot_min
-         WRITE(numout,*) '   Maximum depth of the ocean   rn_bot_max, ', rn_bot_max
+         WRITE(numout,*) '          ENVELOPES:'
+         WRITE(numout,*) '          ========= '
          WRITE(numout,*) ''
-         WRITE(numout,*) '-------------------------------------------------------------------------------'
-         DO je = 1, max_nn_env
-            WRITE(numout,*) 'SUBDOMAIN ', je,':'
+         DO je = 1, tot_env
+            WRITE(numout,*) '             * envelope ', je, ':   min depth = ', min_dep_env(je)
+            WRITE(numout,*) '                               max depth = ', max_dep_env(je)
+            WRITE(numout,*) ''
+         END DO
+         WRITE(numout,*) ''
+         WRITE(numout,*) '          SUBDOMAINS:'
+         WRITE(numout,*) '          ========== '
+         WRITE(numout,*) ''
+         DO je = 1, tot_env
+            WRITE(numout,*) '             * subdomain ', je,':'
             IF ( je == 1) THEN
-               WRITE(numout,*) '   Envelope up  : envelope 0, free surface'
-               WRITE(numout,*) '   Envelope down: envelope ',je,',ln_envl(',je,') = ',ln_envl(je)
+               WRITE(numout,*) '                  Shallower envelope: free surface'
             ELSE
-               WRITE(numout,*) '   Envelope up  : envelope ',je-1,',ln_envl(',je,') = ',ln_envl(je-1)
-               WRITE(numout,*) '   Envelope down: envelope ',je  ,',ln_envl(',je,') = ',ln_envl(je)
+               WRITE(numout,*) '                  Shallower envelope: envelope ',je-1
             END IF
-            WRITE(numout,*) '   max dep of envlp down rn_ebot_max(',je,') =',rn_ebot_max(je)
-            WRITE(numout,*) '   num. of MEs-lev.          nn_slev(',je,') = ',nn_slev(je)
+            WRITE(numout,*) '                  Deeper    envelope: envelope ',je
+            WRITE(numout,*) '                  Number of MEs-levs: nn_slev(',je,') = ',nn_slev(je)
             IF ( isodd(je) ) THEN
-               WRITE(numout,*) '   Stretched s-coordinates: '
+               WRITE(numout,*) '                  Stretched s-coordinates: '
             ELSE
-               WRITE(numout,*) '   Stretched CUBIC SPLINES: '
+               WRITE(numout,*) '                  Stretched CUBIC SPLINES: '
             END IF
-            IF (nn_strt(je) == 0) WRITE(numout,*) '     M96  stretching function'
-            IF (nn_strt(je) == 1) WRITE(numout,*) '     SH94 stretching function'
-            IF (nn_strt(je) == 2) WRITE(numout,*) '     SF12 stretching function'
-            WRITE(numout,*) '     critical depth        rn_e_hc(',je,') = ',rn_e_hc(je)
-            WRITE(numout,*) '     surface stretc. coef. rn_e_th(',je,') = ',rn_e_th(je)
+            IF (nn_strt(je) == 0) WRITE(numout,*) '                    M96  stretching function'
+            IF (nn_strt(je) == 1) WRITE(numout,*) '                    SH94 stretching function'
+            IF (nn_strt(je) == 2) WRITE(numout,*) '                    SF12 stretching function'
+            WRITE(numout,*) '                    critical depth        rn_e_hc(',je,') = ',rn_e_hc(je)
+            WRITE(numout,*) '                    surface stretc. coef. rn_e_th(',je,') = ',rn_e_th(je)
             IF (nn_strt(je) == 2) THEN
-               WRITE(numout,*) '     bottom  stretc. coef. rn_e_ba(',je,') = ',rn_e_ba(je)
+               WRITE(numout,*) '                    bottom  stretc. coef. rn_e_ba(',je,') = ',rn_e_ba(je)
             END IF
-            WRITE(numout,*) '     bottom  stretc. coef. rn_e_bb(',je,') = ',rn_e_bb(je)
+            WRITE(numout,*) '                    bottom  stretc. coef. rn_e_bb(',je,') = ',rn_e_bb(je)
             IF (nn_strt(je) == 2) THEN
-               WRITE(numout,*) '     bottom  stretc. coef. rn_e_al(',je,') = ',rn_e_al(je)
+               WRITE(numout,*) '                 bottom  stretc. coef. rn_e_al(',je,') = ',rn_e_al(je)
             END IF
-            WRITE(numout,*) '-------------------------------------------------------------------------------'
+            WRITE(numout,*) '          ------------------------------------------------------------------'
          ENDDO
       ENDIF
