Index: branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist_cfg
===================================================================
--- branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist_cfg	(revision 6978)
+++ branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist_cfg	(revision 6979)
@@ -13,6 +13,8 @@
 &namcfg        !   parameters of the configuration
 !-----------------------------------------------------------------------
-   cp_cfg      =  "orca"               !  name of the configuration
-   jp_cfg      =       2               !  resolution of the configuration
+   ln_read_cfg = .true.   !  (=T) read the domain configuration in 'domain_cfg.nc" file
+   !                       !  (=F) user defined configuration  ==>>>  see usrdef(_...) modules
+   cp_cfg      = "orca" !  name of the configuration
+   jp_cfg      =      2    !  resolution of the configuration
 /
 !-----------------------------------------------------------------------
@@ -20,23 +22,22 @@
 !-----------------------------------------------------------------------
    ln_zps      = .true.    !  z-coordinate - partial steps
-   ln_linssh   = .true.    !  linear free surface
 /
 !-----------------------------------------------------------------------
 &namdom        !   space and time domain (bathymetry, mesh, timestep)
 !-----------------------------------------------------------------------
-   nn_msh      =    1      !  create (=1) a mesh file or not (=0)
+   ln_linssh   = .true.   !  =T  linear free surface  ==>>  model level are fixed in time
+   !
+   rn_rdt      = 21600.     !  time step for the dynamics (and tracer if nn_acc=0)
+   !
+
+
    rn_rdt      = 21600.    !  time step for the dynamics 
-                           !
-   ppsur       =   -4762.96143546300   !  ORCA r4, r2 and r05 coefficients
-   ppa0        =     255.58049070440   ! (default coefficients)
-   ppa1        =     245.58132232490   !
-   ppkth       =      21.43336197938   !
-   ppacr       =       3.0             !
-   ppdzmin     =  999999.              !  Minimum vertical spacing
-   pphmax      =  999999.              !  Maximum depth
-   ldbletanh   =  .FALSE.              !  Use/do not use double tanf function for vertical coordinates
-   ppa2        =  999999.              !  Double tanh function parameters
-   ppkth2      =  999999.              !
-   ppacr2      =  999999.              !
+/
+!-----------------------------------------------------------------------
+&namlbc        !   lateral momentum boundary condition
+!-----------------------------------------------------------------------
+   !                       !  free slip  !   partial slip  !   no slip   ! strong slip
+   rn_shlat    =    2.     !  shlat = 0  !  0 < shlat < 2  !  shlat = 2  !  2 < shlat
+   ln_vorlat   = .false.   !  consistency of vorticity boundary condition with analytical Eqs.
 /
 !-----------------------------------------------------------------------
Index: branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OFF_SRC/domrea.F90
===================================================================
--- branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OFF_SRC/domrea.F90	(revision 6978)
+++ 	(revision )
@@ -1,833 +1,0 @@
-MODULE domrea
-   !!==============================================================================
-   !!                       ***  MODULE domrea   ***
-   !! Ocean initialization : domain initialization
-   !!==============================================================================
-   !! History :  OPA  ! 1990-10  (C. Levy - G. Madec)  Original code
-   !!                 ! 1992-01  (M. Imbard) insert time step initialization
-   !!                 ! 1996-06  (G. Madec) generalized vertical coordinate 
-   !!                 ! 1997-02  (G. Madec) creation of domwri.F
-   !!                 ! 2001-05  (E.Durand - G. Madec) insert closed sea
-   !!  NEMO      1.0  ! 2002-08  (G. Madec)  F90: Free form and module
-   !!----------------------------------------------------------------------
-
-   !!----------------------------------------------------------------------
-   !!   dom_init       : initialize the space and time domain
-   !!   dom_nam        : read and contral domain namelists
-   !!   dom_ctl        : control print for the ocean domain
-   !!----------------------------------------------------------------------
-   USE oce             ! 
-   USE trc_oce         ! shared ocean/biogeochemical variables
-   USE dom_oce         ! ocean space and time domain
-   USE phycst          ! physical constants
-   USE domstp          ! domain: set the time-step
-   !
-   USE in_out_manager  ! I/O manager
-   USE lib_mpp         ! distributed memory computing library
-   USE lbclnk          ! lateral boundary condition - MPP exchanges
-   USE wrk_nemo  
-   
-   IMPLICIT NONE
-   PRIVATE
-
-   PUBLIC   dom_rea    ! called by nemogcm.F90
-
-   !! * Substitutions
-#  include "vectopt_loop_substitute.h90"
-   !!----------------------------------------------------------------------
-   !! NEMO/OFF 3.7 , NEMO Consortium (2015)
-   !! $Id$
-   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
-   !!----------------------------------------------------------------------
-CONTAINS
-
-   SUBROUTINE dom_rea
-      !!----------------------------------------------------------------------
-      !!                  ***  ROUTINE dom_rea  ***
-      !!                    
-      !! ** Purpose :   Domain initialization. Call the routines that are 
-      !!      required to create the arrays which define the space and time
-      !!      domain of the ocean model.
-      !!
-      !! ** Method  :
-      !!      - dom_stp: defined the model time step
-      !!      - dom_rea: read the meshmask file if nmsh=1
-      !!----------------------------------------------------------------------
-      INTEGER ::   jk          ! dummy loop index
-      INTEGER ::   iconf = 0   ! local integers
-      !!----------------------------------------------------------------------
-      !
-      IF(lwp) THEN
-         WRITE(numout,*)
-         WRITE(numout,*) 'dom_init : domain initialization'
-         WRITE(numout,*) '~~~~~~~~'
-      ENDIF
-      !
-      CALL dom_nam      ! read namelist ( namrun, namdom )
-      CALL dom_zgr      ! Vertical mesh and bathymetry option
-      CALL dom_grd      ! Create a domain file
-      !
-      !                                      ! associated horizontal metrics
-      !
-      r1_e1t(:,:) = 1._wp / e1t(:,:)   ;   r1_e2t (:,:) = 1._wp / e2t(:,:)
-      r1_e1u(:,:) = 1._wp / e1u(:,:)   ;   r1_e2u (:,:) = 1._wp / e2u(:,:)
-      r1_e1v(:,:) = 1._wp / e1v(:,:)   ;   r1_e2v (:,:) = 1._wp / e2v(:,:)
-      r1_e1f(:,:) = 1._wp / e1f(:,:)   ;   r1_e2f (:,:) = 1._wp / e2f(:,:)
-      !
-!!gm BUG if scale factor reduction !!!!
-      e1e2t (:,:) = e1t(:,:) * e2t(:,:)   ;   r1_e1e2t(:,:) = 1._wp / e1e2t(:,:)
-      e1e2u (:,:) = e1u(:,:) * e2u(:,:)   ;   r1_e1e2u(:,:) = 1._wp / e1e2u(:,:)
-      e1e2v (:,:) = e1v(:,:) * e2v(:,:)   ;   r1_e1e2v(:,:) = 1._wp / e1e2v(:,:)
-      e1e2f (:,:) = e1f(:,:) * e2f(:,:)   ;   r1_e1e2f(:,:) = 1._wp / e1e2f(:,:)
-      !   
-      e2_e1u(:,:) = e2u(:,:) / e1u(:,:)
-      e1_e2v(:,:) = e1v(:,:) / e2v(:,:)
-      !
-      hu_n(:,:) = e3u_n(:,:,1) * umask(:,:,1)     ! Ocean depth at U- and V-points
-      hv_n(:,:) = e3v_n(:,:,1) * vmask(:,:,1)
-      DO jk = 2, jpk
-         hu_n(:,:) = hu_n(:,:) + e3u_n(:,:,jk) * umask(:,:,jk)
-         hv_n(:,:) = hv_n(:,:) + e3v_n(:,:,jk) * vmask(:,:,jk)
-      END DO
-      !                                        ! Inverse of the local depth
-      r1_hu_n(:,:) = 1._wp / ( hu_n(:,:) + 1._wp - umask(:,:,1) ) * umask(:,:,1)
-      r1_hv_n(:,:) = 1._wp / ( hv_n(:,:) + 1._wp - vmask(:,:,1) ) * vmask(:,:,1)
-      !
-      CALL dom_stp      ! Time step
-      CALL dom_msk      ! Masks
-      CALL dom_ctl      ! Domain control
-      !
-   END SUBROUTINE dom_rea
-
-
-   SUBROUTINE dom_nam
-      !!----------------------------------------------------------------------
-      !!                     ***  ROUTINE dom_nam  ***
-      !!                    
-      !! ** Purpose :   read domaine namelists and print the variables.
-      !!
-      !! ** input   : - namrun namelist
-      !!              - namdom namelist
-      !!----------------------------------------------------------------------
-      USE ioipsl
-      INTEGER  ::   ios   ! Local integer output status for namelist read
-      !
-      NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list,                         &
-         &             nn_no   , cn_exp    , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl,             &
-         &             nn_it000, nn_itend  , nn_date0    , nn_time0, nn_leapy     , nn_istate , nn_stock ,   &
-         &             nn_write, ln_iscpl  , ln_mskland  , ln_cfmeta    , ln_clobber, nn_chunksz, nn_euler
-      NAMELIST/namdom/ nn_bathy , rn_bathy , rn_e3zps_min, rn_e3zps_rat , nn_msh    , rn_hmin   , rn_isfhmin,&
-         &             rn_atfp  , rn_rdt   , nn_baro     , nn_closea    , ln_crs    ,             &
-         &             ppsur, ppa0, ppa1, ppkth, ppacr, ppdzmin, pphmax, ldbletanh, &
-         &             ppa2, ppkth2, ppacr2
-#if defined key_netcdf4
-      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip
-#endif
-      !!----------------------------------------------------------------------
-
-      REWIND( numnam_ref )              ! Namelist namrun in reference namelist : Parameters of the run
-      READ  ( numnam_ref, namrun, IOSTAT = ios, ERR = 901)
-901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in reference namelist', lwp )
-
-      REWIND( numnam_cfg )              ! Namelist namrun in configuration namelist : Parameters of the run
-      READ  ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 )
-902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in configuration namelist', lwp )
-      IF(lwm) WRITE ( numond, namrun )
-      !
-      IF(lwp) THEN                  ! control print
-         WRITE(numout,*)
-         WRITE(numout,*) 'dom_nam  : domain initialization through namelist read'
-         WRITE(numout,*) '~~~~~~~ '
-         WRITE(numout,*) '   Namelist namrun'  
-         WRITE(numout,*) '      job number                      nn_no      = ', nn_no
-         WRITE(numout,*) '      experiment name for output      cn_exp     = ', cn_exp
-         WRITE(numout,*) '      restart logical                 ln_rstart  = ', ln_rstart
-         WRITE(numout,*) '      control of time step            nn_rstctl  = ', nn_rstctl
-         WRITE(numout,*) '      number of the first time step   nn_it000   = ', nn_it000
-         WRITE(numout,*) '      number of the last time step    nn_itend   = ', nn_itend
-         WRITE(numout,*) '      initial calendar date aammjj    nn_date0   = ', nn_date0
-         WRITE(numout,*) '      leap year calendar (0/1)        nn_leapy   = ', nn_leapy
-         WRITE(numout,*) '      initial state output            nn_istate  = ', nn_istate
-         WRITE(numout,*) '      frequency of restart file       nn_stock   = ', nn_stock
-         WRITE(numout,*) '      frequency of output file        nn_write   = ', nn_write
-         WRITE(numout,*) '      mask land points                ln_mskland = ', ln_mskland
-         WRITE(numout,*) '      additional CF standard metadata ln_cfmeta  = ', ln_cfmeta
-         WRITE(numout,*) '      overwrite an existing file      ln_clobber = ', ln_clobber
-         WRITE(numout,*) '      NetCDF chunksize (bytes)        nn_chunksz = ', nn_chunksz
-      ENDIF
-      no = nn_no                    ! conversion DOCTOR names into model names (this should disappear soon)
-      cexper = cn_exp
-      nrstdt = nn_rstctl
-      nit000 = nn_it000
-      nitend = nn_itend
-      ndate0 = nn_date0
-      nleapy = nn_leapy
-      ninist = nn_istate
-      nstock = nn_stock
-      nstocklist = nn_stocklist
-      nwrite = nn_write
-      !
-      !                             ! control of output frequency
-      IF ( nstock == 0 .OR. nstock > nitend ) THEN
-         WRITE(ctmp1,*) 'nstock = ', nstock, ' it is forced to ', nitend
-         CALL ctl_warn( ctmp1 )
-         nstock = nitend
-      ENDIF
-      IF ( nwrite == 0 ) THEN
-         WRITE(ctmp1,*) 'nwrite = ', nwrite, ' it is forced to ', nitend
-         CALL ctl_warn( ctmp1 )
-         nwrite = nitend
-      ENDIF
-
-      ! parameters correspondting to nit000 - 1 (as we start the step loop with a call to day)
-      ndastp = ndate0 - 1        ! ndate0 read in the namelist in dom_nam, we assume that we start run at 00:00
-      adatrj = ( REAL( nit000-1, wp ) * rdt ) / rday
-
-#if defined key_agrif
-      IF( Agrif_Root() ) THEN
-#endif
-      SELECT CASE ( nleapy )        ! Choose calendar for IOIPSL
-      CASE (  1 ) 
-         CALL ioconf_calendar('gregorian')
-         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "gregorian", i.e. leap year'
-      CASE (  0 )
-         CALL ioconf_calendar('noleap')
-         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "noleap", i.e. no leap year'
-      CASE ( 30 )
-         CALL ioconf_calendar('360d')
-         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "360d", i.e. 360 days in a year'
-      END SELECT
-#if defined key_agrif
-      ENDIF
-#endif
-
-      REWIND( numnam_ref )              ! Namelist namdom in reference namelist : space & time domain (bathymetry, mesh, timestep)
-      READ  ( numnam_ref, namdom, IOSTAT = ios, ERR = 903)
-903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in reference namelist', lwp )
-
-      REWIND( numnam_cfg )              ! Namelist namdom in configuration namelist : space & time domain (bathymetry, mesh, timestep)
-      READ  ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 )
-904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in configuration namelist', lwp )
-      IF(lwm) WRITE ( numond, namdom )
-
-      IF(lwp) THEN
-         WRITE(numout,*) 
-         WRITE(numout,*) '   Namelist namdom : space & time domain'
-         WRITE(numout,*) '      flag read/compute bathymetry      nn_bathy     = ', nn_bathy
-         WRITE(numout,*) '      Depth (if =0 bathy=jpkm1)         rn_bathy     = ', rn_bathy
-         WRITE(numout,*) '      min depth of the ocean    (>0) or    rn_hmin   = ', rn_hmin
-         WRITE(numout,*) '      minimum thickness of partial      rn_e3zps_min = ', rn_e3zps_min, ' (m)'
-         WRITE(numout,*) '         step level                     rn_e3zps_rat = ', rn_e3zps_rat
-         WRITE(numout,*) '      create mesh/mask file(s)          nn_msh       = ', nn_msh
-         WRITE(numout,*) '           = 0   no file created                 '
-         WRITE(numout,*) '           = 1   mesh_mask                       '
-         WRITE(numout,*) '           = 2   mesh and mask                   '
-         WRITE(numout,*) '           = 3   mesh_hgr, msh_zgr and mask      '
-         WRITE(numout,*) '      ocean time step                      rn_rdt    = ', rn_rdt
-         WRITE(numout,*) '      asselin time filter parameter        rn_atfp   = ', rn_atfp
-         WRITE(numout,*) '      time-splitting: nb of sub time-step  nn_baro   = ', nn_baro
-         WRITE(numout,*) '      suppression of closed seas (=0)      nn_closea = ', nn_closea
-
-         WRITE(numout,*) '      ORCA r4, r2 and r05 coefficients  ppsur           = ', ppsur
-         WRITE(numout,*) '                                        ppa0            = ', ppa0
-         WRITE(numout,*) '                                        ppa1            = ', ppa1
-         WRITE(numout,*) '                                        ppkth           = ', ppkth
-         WRITE(numout,*) '                                        ppacr           = ', ppacr
-         WRITE(numout,*) '      Minimum vertical spacing ppdzmin                  = ', ppdzmin
-         WRITE(numout,*) '      Maximum depth pphmax                              = ', pphmax
-         WRITE(numout,*) '      Use double tanf function for vertical coordinates ldbletanh = ', ldbletanh
-         WRITE(numout,*) '      Double tanh function parameters ppa2              = ', ppa2
-         WRITE(numout,*) '                                      ppkth2            = ', ppkth2
-         WRITE(numout,*) '                                      ppacr2            = ', ppacr2
-      ENDIF
-
-      ntopo     = nn_bathy          ! conversion DOCTOR names into model names (this should disappear soon)
-      e3zps_min = rn_e3zps_min
-      e3zps_rat = rn_e3zps_rat
-      nmsh      = nn_msh
-      atfp      = rn_atfp
-      rdt       = rn_rdt
-#if defined key_netcdf4
-      !                             ! NetCDF 4 case   ("key_netcdf4" defined)
-      REWIND( numnam_ref )              ! Namelist namnc4 in reference namelist : NETCDF
-      READ  ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907)
-907   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in reference namelist', lwp )
-
-      REWIND( numnam_cfg )              ! Namelist namnc4 in configuration namelist : NETCDF
-      READ  ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 )
-908   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in configuration namelist', lwp )
-      IF(lwm) WRITE( numond, namnc4 )
-      IF(lwp) THEN                        ! control print
-         WRITE(numout,*)
-         WRITE(numout,*) '   Namelist namnc4 - Netcdf4 chunking parameters'
-         WRITE(numout,*) '      number of chunks in i-dimension      nn_nchunks_i   = ', nn_nchunks_i
-         WRITE(numout,*) '      number of chunks in j-dimension      nn_nchunks_j   = ', nn_nchunks_j
-         WRITE(numout,*) '      number of chunks in k-dimension      nn_nchunks_k   = ', nn_nchunks_k
-         WRITE(numout,*) '      apply netcdf4/hdf5 chunking & compression ln_nc4zip = ', ln_nc4zip
-      ENDIF
-
-      ! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module)
-      ! Note the chunk size in the unlimited (time) dimension will be fixed at 1
-      snc4set%ni   = nn_nchunks_i
-      snc4set%nj   = nn_nchunks_j
-      snc4set%nk   = nn_nchunks_k
-      snc4set%luse = ln_nc4zip
-#else
-      snc4set%luse = .FALSE.        ! No NetCDF 4 case
-#endif
-      !
-   END SUBROUTINE dom_nam
-
-
-   SUBROUTINE dom_zgr
-      !!----------------------------------------------------------------------
-      !!                ***  ROUTINE dom_zgr  ***
-      !!                   
-      !! ** Purpose :  set the depth of model levels and the resulting 
-      !!      vertical scale factors.
-      !!
-      !! ** Method  : - reference 1D vertical coordinate (gdep._1d, e3._1d)
-      !!              - read/set ocean depth and ocean levels (bathy, mbathy)
-      !!              - vertical coordinate (gdep., e3.) depending on the 
-      !!                coordinate chosen :
-      !!                   ln_zco=T   z-coordinate  
-      !!                   ln_zps=T   z-coordinate with partial steps
-      !!                   ln_zco=T   s-coordinate 
-      !!
-      !! ** Action  :   define gdep., e3., mbathy and bathy
-      !!----------------------------------------------------------------------
-      INTEGER ::   ioptio = 0   ! temporary integer
-      INTEGER ::   ios
-      !!
-      NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco, ln_isfcav, ln_linssh
-      !!----------------------------------------------------------------------
-
-      REWIND( numnam_ref )              ! Namelist namzgr in reference namelist : Vertical coordinate
-      READ  ( numnam_ref, namzgr, IOSTAT = ios, ERR = 901 )
-901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzgr in reference namelist', lwp )
-
-      REWIND( numnam_cfg )              ! Namelist namzgr in configuration namelist : Vertical coordinate
-      READ  ( numnam_cfg, namzgr, IOSTAT = ios, ERR = 902 )
-902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzgr in configuration namelist', lwp )
-      IF(lwm) WRITE ( numond, namzgr )
-
-      IF(lwp) THEN                     ! Control print
-         WRITE(numout,*)
-         WRITE(numout,*) 'dom_zgr : vertical coordinate'
-         WRITE(numout,*) '~~~~~~~'
-         WRITE(numout,*) '          Namelist namzgr : set vertical coordinate'
-         WRITE(numout,*) '             z-coordinate - full steps      ln_zco    = ', ln_zco
-         WRITE(numout,*) '             z-coordinate - partial steps   ln_zps    = ', ln_zps
-         WRITE(numout,*) '             s- or hybrid z-s-coordinate    ln_sco    = ', ln_sco
-         WRITE(numout,*) '             ice shelf cavity               ln_isfcav = ', ln_isfcav
-         WRITE(numout,*) '             Linear free surface            ln_linssh = ', ln_linssh
-      ENDIF
-
-      ioptio = 0                       ! Check Vertical coordinate options
-      IF( ln_zco ) ioptio = ioptio + 1
-      IF( ln_zps ) ioptio = ioptio + 1
-      IF( ln_sco ) ioptio = ioptio + 1
-      IF( ln_isfcav ) ioptio = 33
-      IF ( ioptio /= 1  )   CALL ctl_stop( ' none or several vertical coordinate options used' )
-      IF ( ioptio == 33 )   CALL ctl_stop( ' isf cavity with off line module not yet done    ' )
-
-   END SUBROUTINE dom_zgr
-
-
-   SUBROUTINE dom_ctl
-      !!----------------------------------------------------------------------
-      !!                     ***  ROUTINE dom_ctl  ***
-      !!
-      !! ** Purpose :   Domain control.
-      !!
-      !! ** Method  :   compute and print extrema of masked scale factors
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER ::   iimi1, ijmi1, iimi2, ijmi2, iima1, ijma1, iima2, ijma2
-      INTEGER, DIMENSION(2) ::   iloc      ! 
-      REAL(wp) ::   ze1min, ze1max, ze2min, ze2max
-      !!----------------------------------------------------------------------
-
-      ! Extrema of the scale factors
-
-      IF(lwp)WRITE(numout,*)
-      IF(lwp)WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors'
-      IF(lwp)WRITE(numout,*) '~~~~~~~'
-
-      IF (lk_mpp) THEN
-         CALL mpp_minloc( e1t(:,:), tmask(:,:,1), ze1min, iimi1,ijmi1 )
-         CALL mpp_minloc( e2t(:,:), tmask(:,:,1), ze2min, iimi2,ijmi2 )
-         CALL mpp_maxloc( e1t(:,:), tmask(:,:,1), ze1max, iima1,ijma1 )
-         CALL mpp_maxloc( e2t(:,:), tmask(:,:,1), ze2max, iima2,ijma2 )
-      ELSE
-         ze1min = MINVAL( e1t(:,:), mask = tmask(:,:,1) == 1.e0 )    
-         ze2min = MINVAL( e2t(:,:), mask = tmask(:,:,1) == 1.e0 )    
-         ze1max = MAXVAL( e1t(:,:), mask = tmask(:,:,1) == 1.e0 )    
-         ze2max = MAXVAL( e2t(:,:), mask = tmask(:,:,1) == 1.e0 )    
-
-         iloc  = MINLOC( e1t(:,:), mask = tmask(:,:,1) == 1.e0 )
-         iimi1 = iloc(1) + nimpp - 1
-         ijmi1 = iloc(2) + njmpp - 1
-         iloc  = MINLOC( e2t(:,:), mask = tmask(:,:,1) == 1.e0 )
-         iimi2 = iloc(1) + nimpp - 1
-         ijmi2 = iloc(2) + njmpp - 1
-         iloc  = MAXLOC( e1t(:,:), mask = tmask(:,:,1) == 1.e0 )
-         iima1 = iloc(1) + nimpp - 1
-         ijma1 = iloc(2) + njmpp - 1
-         iloc  = MAXLOC( e2t(:,:), mask = tmask(:,:,1) == 1.e0 )
-         iima2 = iloc(1) + nimpp - 1
-         ijma2 = iloc(2) + njmpp - 1
-      ENDIF
-      !
-      IF(lwp) THEN
-         WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, iima1, ijma1
-         WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, iimi1, ijmi1
-         WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, iima2, ijma2
-         WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, iimi2, ijmi2
-      ENDIF
-      !
-   END SUBROUTINE dom_ctl
-
-
-   SUBROUTINE dom_grd
-      !!----------------------------------------------------------------------
-      !!                  ***  ROUTINE dom_grd  ***
-      !!                   
-      !! ** Purpose :  Read the NetCDF file(s) which contain(s) all the
-      !!      ocean domain informations (mesh and mask arrays). This (these)
-      !!      file(s) is (are) used for visualisation (SAXO software) and
-      !!      diagnostic computation.
-      !!
-      !! ** Method  :   Read in a file all the arrays generated in routines
-      !!      domhgr, domzgr, and dommsk. Note: the file contain depends on
-      !!      the vertical coord. used (z-coord, partial steps, s-coord)
-      !!                    nmsh = 1  :   'mesh_mask.nc' file
-      !!                         = 2  :   'mesh.nc' and mask.nc' files
-      !!                         = 3  :   'mesh_hgr.nc', 'mesh_zgr.nc' and
-      !!                                  'mask.nc' files
-      !!      For huge size domain, use option 2 or 3 depending on your 
-      !!      vertical coordinate.
-      !!
-      !! ** input file : 
-      !!      meshmask.nc  : domain size, horizontal grid-point position,
-      !!                     masks, depth and vertical scale factors
-      !!----------------------------------------------------------------------
-      USE iom
-      !!
-      INTEGER  ::   ji, jj, jk   ! dummy loop indices
-      INTEGER  ::   ik, inum0 , inum1 , inum2 , inum3 , inum4   ! local integers
-      REAL(wp) ::   zrefdep         ! local real
-      REAL(wp), POINTER, DIMENSION(:,:) :: zmbk, zprt, zprw
-      !!----------------------------------------------------------------------
-
-      IF(lwp) WRITE(numout,*)
-      IF(lwp) WRITE(numout,*) 'dom_rea : read NetCDF mesh and mask information file(s)'
-      IF(lwp) WRITE(numout,*) '~~~~~~~'
-
-      CALL wrk_alloc( jpi, jpj, zmbk, zprt, zprw )
-
-      zmbk(:,:) = 0._wp
-
-      SELECT CASE (nmsh)
-         !                                     ! ============================
-         CASE ( 1 )                            !  create 'mesh_mask.nc' file
-            !                                  ! ============================
-
-            IF(lwp) WRITE(numout,*) '          one file in "mesh_mask.nc" '
-            CALL iom_open( 'mesh_mask', inum0 )
-
-            inum2 = inum0                                            ! put all the informations
-            inum3 = inum0                                            ! in unit inum0
-            inum4 = inum0
-
-            !                                  ! ============================
-         CASE ( 2 )                            !  create 'mesh.nc' and 
-            !                                  !         'mask.nc' files
-            !                                  ! ============================
-
-            IF(lwp) WRITE(numout,*) '          two files in "mesh.nc" and "mask.nc" '
-            CALL iom_open( 'mesh', inum1 )
-            CALL iom_open( 'mask', inum2 )
-
-            inum3 = inum1                                            ! put mesh informations 
-            inum4 = inum1                                            ! in unit inum1 
-
-            !                                  ! ============================
-         CASE ( 3 )                            !  create 'mesh_hgr.nc'
-            !                                  !         'mesh_zgr.nc' and
-            !                                  !         'mask.nc'     files
-            !                                  ! ============================
-
-            IF(lwp) WRITE(numout,*) '          three files in "mesh_hgr.nc" , "mesh_zgr.nc" and "mask.nc" '
-            CALL iom_open( 'mesh_hgr', inum3 ) ! create 'mesh_hgr.nc'
-            CALL iom_open( 'mesh_zgr', inum4 ) ! create 'mesh_zgr.nc'
-            CALL iom_open( 'mask'    , inum2 ) ! create 'mask.nc'
-
-            !                                  ! ===========================
-         CASE DEFAULT                          ! return error 
-            !                                  ! mesh has to be provided
-            !                                  ! ===========================
-            CALL ctl_stop( ' OFFLINE mode requires the input mesh mask(s). ',   &
-            &                                 ' Invalid nn_msh value in the namelist (0 is not allowed)' )
-
-         END SELECT
-
-         !                                                         ! masks (inum2) 
-         CALL iom_get( inum2, jpdom_data, 'tmask', tmask )
-         CALL iom_get( inum2, jpdom_data, 'umask', umask )
-         CALL iom_get( inum2, jpdom_data, 'vmask', vmask )
-         CALL iom_get( inum2, jpdom_data, 'fmask', fmask )
-
-         CALL lbc_lnk( tmask, 'T', 1._wp )    ! Lateral boundary conditions
-         CALL lbc_lnk( umask, 'U', 1._wp )      
-         CALL lbc_lnk( vmask, 'V', 1._wp )
-         CALL lbc_lnk( fmask, 'F', 1._wp )
-
-#if defined key_c1d
-         ! set umask and vmask equal tmask in 1D configuration
-         IF(lwp) WRITE(numout,*)
-         IF(lwp) WRITE(numout,*) '**********  1D configuration : set umask and vmask equal tmask ********'
-         IF(lwp) WRITE(numout,*) '**********                                                     ********'
-
-         umask(:,:,:) = tmask(:,:,:)
-         vmask(:,:,:) = tmask(:,:,:)
-#endif
-
-#if defined key_degrad
-         CALL iom_get( inum2, jpdom_data, 'facvolt', facvol )
-#endif
-         !                                                         ! horizontal mesh (inum3)
-         CALL iom_get( inum3, jpdom_data, 'glamt', glamt )
-         CALL iom_get( inum3, jpdom_data, 'glamu', glamu )
-         CALL iom_get( inum3, jpdom_data, 'glamv', glamv )
-         CALL iom_get( inum3, jpdom_data, 'glamf', glamf )
-
-         CALL iom_get( inum3, jpdom_data, 'gphit', gphit )
-         CALL iom_get( inum3, jpdom_data, 'gphiu', gphiu )
-         CALL iom_get( inum3, jpdom_data, 'gphiv', gphiv )
-         CALL iom_get( inum3, jpdom_data, 'gphif', gphif )
-
-         CALL iom_get( inum3, jpdom_data, 'e1t', e1t )
-         CALL iom_get( inum3, jpdom_data, 'e1u', e1u )
-         CALL iom_get( inum3, jpdom_data, 'e1v', e1v )
-         
-         CALL iom_get( inum3, jpdom_data, 'e2t', e2t )
-         CALL iom_get( inum3, jpdom_data, 'e2u', e2u )
-         CALL iom_get( inum3, jpdom_data, 'e2v', e2v )
-
-         CALL iom_get( inum3, jpdom_data, 'ff', ff )
-
-         CALL iom_get( inum4, jpdom_data, 'mbathy', zmbk )              ! number of ocean t-points
-         mbathy (:,:) = INT( zmbk(:,:) )
-         misfdep(:,:) = 1                                               ! ice shelf case not yet done
-         
-         CALL zgr_bot_level                                             ! mbk. arrays (deepest ocean t-, u- & v-points
-         !
-         IF( ln_sco ) THEN                                         ! s-coordinate
-            CALL iom_get( inum4, jpdom_data, 'hbatt', hbatt )
-            CALL iom_get( inum4, jpdom_data, 'hbatu', hbatu )
-            CALL iom_get( inum4, jpdom_data, 'hbatv', hbatv )
-            CALL iom_get( inum4, jpdom_data, 'hbatf', hbatf )
-            
-            CALL iom_get( inum4, jpdom_unknown, 'gsigt', gsigt ) ! scaling coef.
-            CALL iom_get( inum4, jpdom_unknown, 'gsigw', gsigw )
-            CALL iom_get( inum4, jpdom_unknown, 'gsi3w', gsi3w ) 
-            CALL iom_get( inum4, jpdom_unknown, 'esigt', esigt )
-            CALL iom_get( inum4, jpdom_unknown, 'esigw', esigw )
-
-            CALL iom_get( inum4, jpdom_data, 'e3t_0', e3t_0(:,:,:) ) ! scale factors
-            CALL iom_get( inum4, jpdom_data, 'e3u_0', e3u_0(:,:,:) )
-            CALL iom_get( inum4, jpdom_data, 'e3v_0', e3v_0(:,:,:) )
-            CALL iom_get( inum4, jpdom_data, 'e3w_0', e3w_0(:,:,:) )
-
-            CALL iom_get( inum4, jpdom_unknown, 'gdept_1d', gdept_1d ) ! depth
-            CALL iom_get( inum4, jpdom_unknown, 'gdepw_1d', gdepw_1d )
-         ENDIF
-
- 
-         IF( ln_zps ) THEN                                           ! z-coordinate - partial steps
-            CALL iom_get( inum4, jpdom_unknown, 'gdept_1d', gdept_1d )  ! reference depth
-            CALL iom_get( inum4, jpdom_unknown, 'gdepw_1d', gdepw_1d )
-            CALL iom_get( inum4, jpdom_unknown, 'e3t_1d'  , e3t_1d   )    ! reference scale factors
-            CALL iom_get( inum4, jpdom_unknown, 'e3w_1d'  , e3w_1d   )
-            !
-            IF( nmsh <= 6 ) THEN                                        ! 3D vertical scale factors
-               CALL iom_get( inum4, jpdom_data, 'e3t_0', e3t_0(:,:,:) )
-               CALL iom_get( inum4, jpdom_data, 'e3u_0', e3u_0(:,:,:) )
-               CALL iom_get( inum4, jpdom_data, 'e3v_0', e3v_0(:,:,:) )
-               CALL iom_get( inum4, jpdom_data, 'e3w_0', e3w_0(:,:,:) )
-            ELSE                                                        ! 2D bottom scale factors
-               CALL iom_get( inum4, jpdom_data, 'e3t_ps', e3tp )
-               CALL iom_get( inum4, jpdom_data, 'e3w_ps', e3wp )
-               !                                                        ! deduces the 3D scale factors
-               DO jk = 1, jpk
-                  e3t_0(:,:,jk) = e3t_1d(jk)                                    ! set to the ref. factors
-                  e3u_0(:,:,jk) = e3t_1d(jk)
-                  e3v_0(:,:,jk) = e3t_1d(jk)
-                  e3w_0(:,:,jk) = e3w_1d(jk)
-               END DO
-               DO jj = 1,jpj                                                  ! adjust the deepest values
-                  DO ji = 1,jpi
-                     ik = mbkt(ji,jj)
-                     e3t_0(ji,jj,ik) = e3tp(ji,jj) * tmask(ji,jj,1) + e3t_1d(1) * ( 1._wp - tmask(ji,jj,1) )
-                     e3w_0(ji,jj,ik) = e3wp(ji,jj) * tmask(ji,jj,1) + e3w_1d(1) * ( 1._wp - tmask(ji,jj,1) )
-                  END DO
-               END DO
-               DO jk = 1,jpk                         ! Computed as the minimum of neighbooring scale factors
-                  DO jj = 1, jpjm1
-                     DO ji = 1, jpim1
-                        e3u_0(ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji+1,jj,jk) )
-                        e3v_0(ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji,jj+1,jk) )
-                     END DO
-                  END DO
-               END DO
-               CALL lbc_lnk( e3u_0(:,:,:) , 'U', 1._wp )   ;   CALL lbc_lnk( e3uw_0(:,:,:), 'U', 1._wp )   ! lateral boundary conditions
-               CALL lbc_lnk( e3v_0(:,:,:) , 'V', 1._wp )   ;   CALL lbc_lnk( e3vw_0(:,:,:), 'V', 1._wp )
-               !
-               DO jk = 1, jpk                        ! set to z-scale factor if zero (i.e. along closed boundaries)
-                  WHERE( e3u_0(:,:,jk) == 0._wp )   e3u_0(:,:,jk) = e3t_1d(jk)
-                  WHERE( e3v_0(:,:,jk) == 0._wp )   e3v_0(:,:,jk) = e3t_1d(jk)
-               END DO
-            END IF
-
-            IF( iom_varid( inum4, 'gdept_0', ldstop = .FALSE. ) > 0 ) THEN   ! 3D depth of t- and w-level
-               CALL iom_get( inum4, jpdom_data, 'gdept_0', gdept_0(:,:,:) )
-               CALL iom_get( inum4, jpdom_data, 'gdepw_0', gdepw_0(:,:,:) )
-            ELSE                                                           ! 2D bottom depth
-               CALL iom_get( inum4, jpdom_data, 'hdept', zprt )
-               CALL iom_get( inum4, jpdom_data, 'hdepw', zprw )
-               !
-               DO jk = 1, jpk                                              ! deduces the 3D depth
-                  gdept_0(:,:,jk) = gdept_1d(jk)
-                  gdepw_0(:,:,jk) = gdepw_1d(jk)
-               END DO
-               DO jj = 1, jpj
-                  DO ji = 1, jpi
-                     ik = mbkt(ji,jj)
-                     IF( ik > 0 ) THEN
-                        gdepw_0(ji,jj,ik+1) = zprw(ji,jj)
-                        gdept_0(ji,jj,ik  ) = zprt(ji,jj)
-                        gdept_0(ji,jj,ik+1) = gdept_0(ji,jj,ik) + e3t_0(ji,jj,ik)
-                     ENDIF
-                  END DO
-               END DO
-            ENDIF
-            !
-         ENDIF
-
-         IF( ln_zco ) THEN           ! Vertical coordinates and scales factors
-            CALL iom_get( inum4, jpdom_unknown, 'gdept_1d', gdept_1d ) ! depth
-            CALL iom_get( inum4, jpdom_unknown, 'gdepw_1d', gdepw_1d )
-            CALL iom_get( inum4, jpdom_unknown, 'e3t_1d'  , e3t_1d   )
-            CALL iom_get( inum4, jpdom_unknown, 'e3w_1d'  , e3w_1d   )
-            DO jk = 1, jpk
-               e3t_0(:,:,jk) = e3t_1d(jk)                              ! set to the ref. factors
-               e3u_0(:,:,jk) = e3t_1d(jk)
-               e3v_0(:,:,jk) = e3t_1d(jk)
-               e3w_0(:,:,jk) = e3w_1d(jk)
-               gdept_0(:,:,jk) = gdept_1d(jk)
-               gdepw_0(:,:,jk) = gdepw_1d(jk)
-            END DO
-         ENDIF
-
-      !
-      !              !==  time varying part of coordinate system  ==!
-      !
-      !       before        !          now          !       after         !
-      ;  gdept_b = gdept_0  ;   gdept_n = gdept_0   !        ---          ! depth of grid-points
-      ;  gdepw_b = gdepw_0  ;   gdepw_n = gdepw_0   !        ---          !
-      ;                     ;   gde3w_n = gde3w_0   !        ---          !
-      !
-      ;    e3t_b =   e3t_0  ;     e3t_n =   e3t_0   ;   e3t_a =  e3t_0    ! scale factors
-      ;    e3u_b =   e3u_0  ;     e3u_n =   e3u_0   ;   e3u_a =  e3u_0    !
-      ;    e3v_b =   e3v_0  ;     e3v_n =   e3v_0   ;   e3v_a =  e3v_0    !
-      ;                     ;     e3f_n =   e3f_0   !        ---          !
-      ;    e3w_b =   e3w_0  ;     e3w_n =   e3w_0   !        ---          !
-      ;   e3uw_b =  e3uw_0  ;    e3uw_n =  e3uw_0   !        ---          !
-      ;   e3vw_b =  e3vw_0  ;    e3vw_n =  e3vw_0   !        ---          !
-      !
-
-!!gm BUG in s-coordinate this does not work!
-      ! deepest/shallowest W level Above/Below ~10m
-      zrefdep = 10._wp - ( 0.1_wp * MINVAL(e3w_1d) )                 ! ref. depth with tolerance (10% of minimum layer thickness)
-      nlb10 = MINLOC( gdepw_1d, mask = gdepw_1d > zrefdep, dim = 1 ) ! shallowest W level Below ~10m
-      nla10 = nlb10 - 1                                              ! deepest    W level Above ~10m
-!!gm end bug
-
-      ! Control printing : Grid informations (if not restart)
-      ! ----------------
-
-      IF(lwp .AND. .NOT.ln_rstart ) THEN
-         WRITE(numout,*)
-         WRITE(numout,*) '          longitude and e1 scale factors'
-         WRITE(numout,*) '          ------------------------------'
-         WRITE(numout,9300) ( ji, glamt(ji,1), glamu(ji,1),   &
-            glamv(ji,1), glamf(ji,1),   &
-            e1t(ji,1), e1u(ji,1),   &
-            e1v(ji,1), ji = 1, jpi,10)
-9300     FORMAT( 1x, i4, f8.2,1x, f8.2,1x, f8.2,1x, f8.2, 1x,    &
-            f19.10, 1x, f19.10, 1x, f19.10 )
-
-         WRITE(numout,*)
-         WRITE(numout,*) '          latitude and e2 scale factors'
-         WRITE(numout,*) '          -----------------------------'
-         WRITE(numout,9300) ( jj, gphit(1,jj), gphiu(1,jj),   &
-            &                     gphiv(1,jj), gphif(1,jj),   &
-            &                     e2t  (1,jj), e2u  (1,jj),   &
-            &                     e2v  (1,jj), jj = 1, jpj, 10 )
-      ENDIF
-
-      IF(lwp) THEN
-         WRITE(numout,*)
-         WRITE(numout,*) '              Reference z-coordinate depth and scale factors:'
-         WRITE(numout, "(9x,' level   gdept    gdepw     e3t      e3w  ')" )
-         WRITE(numout, "(10x, i4, 4f9.2)" ) ( jk, gdept_1d(jk), gdepw_1d(jk), e3t_1d(jk), e3w_1d(jk), jk = 1, jpk )
-      ENDIF
-
-      DO jk = 1, jpk
-         IF( e3w_1d  (jk) <= 0._wp .OR. e3t_1d  (jk) <= 0._wp )   CALL ctl_stop( ' e3w_1d or e3t_1d =< 0 ' )
-         IF( gdepw_1d(jk) <  0._wp .OR. gdept_1d(jk) <  0._wp )   CALL ctl_stop( ' gdepw_1d or gdept_1d < 0 ' )
-      END DO
-      !                                     ! ============================
-      !                                     !        close the files 
-      !                                     ! ============================
-      SELECT CASE ( nmsh )
-         CASE ( 1 )                
-            CALL iom_close( inum0 )
-         CASE ( 2 )
-            CALL iom_close( inum1 )
-            CALL iom_close( inum2 )
-         CASE ( 3 )
-            CALL iom_close( inum2 )
-            CALL iom_close( inum3 )
-            CALL iom_close( inum4 )
-      END SELECT
-      !
-      CALL wrk_dealloc( jpi, jpj, zmbk, zprt, zprw )
-      !
-   END SUBROUTINE dom_grd
-
-
-   SUBROUTINE zgr_bot_level
-      !!----------------------------------------------------------------------
-      !!                    ***  ROUTINE zgr_bot_level  ***
-      !!
-      !! ** Purpose :   defines the vertical index of ocean bottom (mbk. arrays)
-      !!
-      !! ** Method  :   computes from mbathy with a minimum value of 1 over land
-      !!
-      !! ** Action  :   mbkt, mbku, mbkv :   vertical indices of the deeptest
-      !!                                     ocean level at t-, u- & v-points
-      !!                                     (min value = 1 over land)
-      !!----------------------------------------------------------------------
-      INTEGER ::   ji, jj   ! dummy loop indices
-      REAL(wp), POINTER, DIMENSION(:,:) :: zmbk
-      !!----------------------------------------------------------------------
-
-      !
-      IF(lwp) WRITE(numout,*)
-      IF(lwp) WRITE(numout,*) '    zgr_bot_level : ocean bottom k-index of T-, U-, V- and W-levels '
-      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~'
-      !
-      CALL wrk_alloc( jpi, jpj, zmbk )
-      !
-      mbkt(:,:) = MAX( mbathy(:,:) , 1 )    ! bottom k-index of T-level (=1 over land)
-      mikt(:,:) = 1 ; miku(:,:) = 1; mikv(:,:) = 1; ! top k-index of T-level (=1 over open ocean; >1 beneath ice shelf)
-      !                                     ! bottom k-index of W-level = mbkt+1
-      DO jj = 1, jpjm1                      ! bottom k-index of u- (v-) level
-         DO ji = 1, jpim1
-            mbku(ji,jj) = MIN(  mbkt(ji+1,jj  ) , mbkt(ji,jj)  )
-            mbkv(ji,jj) = MIN(  mbkt(ji  ,jj+1) , mbkt(ji,jj)  )
-         END DO
-      END DO
-      ! converte into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk 
-      zmbk(:,:) = REAL( mbku(:,:), wp )   ;   CALL lbc_lnk(zmbk,'U',1.)   ;   mbku  (:,:) = MAX( INT( zmbk(:,:) ), 1 )
-      zmbk(:,:) = REAL( mbkv(:,:), wp )   ;   CALL lbc_lnk(zmbk,'V',1.)   ;   mbkv  (:,:) = MAX( INT( zmbk(:,:) ), 1 )
-      !
-      CALL wrk_dealloc( jpi, jpj, zmbk )
-      !
-   END SUBROUTINE zgr_bot_level
-
-
-   SUBROUTINE dom_msk
-      !!---------------------------------------------------------------------
-      !!                 ***  ROUTINE dom_msk  ***
-      !!
-      !! ** Purpose :   Off-line case: defines the interior domain T-mask.
-      !!
-      !! ** Method  :   The interior ocean/land mask is computed from tmask
-      !!              setting to zero the duplicated row and lines due to
-      !!              MPP exchange halos, est-west cyclic and north fold
-      !!              boundary conditions.
-      !!
-      !! ** Action :   tmask_i  : interiorland/ocean mask at t-point
-      !!               tpol     : ???
-      !!----------------------------------------------------------------------
-      INTEGER  ::   ji, jj, jk           ! dummy loop indices
-      INTEGER  ::   iif, iil, ijf, ijl   ! local integers
-      INTEGER, POINTER, DIMENSION(:,:) ::  imsk 
-      !!---------------------------------------------------------------------
-      
-      CALL wrk_alloc( jpi, jpj, imsk )
-      !
-      ! Interior domain mask (used for global sum)
-      ! --------------------
-      ssmask(:,:)  = tmask(:,:,1)
-      tmask_i(:,:) = tmask(:,:,1)
-      iif = jpreci                        ! thickness of exchange halos in i-axis
-      iil = nlci - jpreci + 1
-      ijf = jprecj                        ! thickness of exchange halos in j-axis
-      ijl = nlcj - jprecj + 1
-      !
-      tmask_i( 1 :iif,   :   ) = 0._wp    ! first columns
-      tmask_i(iil:jpi,   :   ) = 0._wp    ! last  columns (including mpp extra columns)
-      tmask_i(   :   , 1 :ijf) = 0._wp    ! first rows
-      tmask_i(   :   ,ijl:jpj) = 0._wp    ! last  rows (including mpp extra rows)
-      !
-      !                                   ! north fold mask
-      tpol(1:jpiglo) = 1._wp
-      !                                
-      IF( jperio == 3 .OR. jperio == 4 )   tpol(jpiglo/2+1:jpiglo) = 0._wp    ! T-point pivot
-      IF( jperio == 5 .OR. jperio == 6 )   tpol(     1    :jpiglo) = 0._wp    ! F-point pivot
-      IF( jperio == 3 .OR. jperio == 4 ) THEN      ! T-point pivot: only half of the nlcj-1 row
-         IF( mjg(ijl-1) == jpjglo-1 ) THEN
-            DO ji = iif+1, iil-1
-               tmask_i(ji,ijl-1) = tmask_i(ji,ijl-1) * tpol(mig(ji))
-            END DO
-         ENDIF
-      ENDIF 
-      !
-      ! (ISF) MIN(1,SUM(umask)) is here to check if you have effectively at
-      ! least 1 wet u point
-      DO jj = 1, jpjm1
-         DO ji = 1, fs_jpim1   ! vector loop
-            ssumask(ji,jj)  = ssmask(ji,jj) * ssmask(ji+1,jj  )  * MIN(1._wp,SUM(umask(ji,jj,:)))
-            ssvmask(ji,jj)  = ssmask(ji,jj) * ssmask(ji  ,jj+1)  * MIN(1._wp,SUM(vmask(ji,jj,:)))
-         END DO
-         DO ji = 1, jpim1      ! NO vector opt.
-            ssfmask(ji,jj) =  ssmask(ji,jj  ) * ssmask(ji+1,jj  )   &
-               &            * ssmask(ji,jj+1) * ssmask(ji+1,jj+1) * MIN(1._wp,SUM(fmask(ji,jj,:)))
-         END DO
-      END DO
-      CALL lbc_lnk( ssumask, 'U', 1._wp )      ! Lateral boundary conditions
-      CALL lbc_lnk( ssvmask, 'V', 1._wp )
-      CALL lbc_lnk( ssfmask, 'F', 1._wp )
-
-      ! 3. Ocean/land mask at wu-, wv- and w points 
-      !----------------------------------------------
-      wmask (:,:,1) = tmask(:,:,1)        ! surface value
-      wumask(:,:,1) = umask(:,:,1) 
-      wvmask(:,:,1) = vmask(:,:,1)
-      DO jk = 2, jpk                      ! deeper value
-         wmask (:,:,jk) = tmask(:,:,jk) * tmask(:,:,jk-1)
-         wumask(:,:,jk) = umask(:,:,jk) * umask(:,:,jk-1)   
-         wvmask(:,:,jk) = vmask(:,:,jk) * vmask(:,:,jk-1)
-      END DO
-      !
-      CALL wrk_dealloc( jpi, jpj, imsk )
-      !
-   END SUBROUTINE dom_msk
-
-   !!======================================================================
-END MODULE domrea
-
Index: branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OFF_SRC/nemogcm.F90
===================================================================
--- branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OFF_SRC/nemogcm.F90	(revision 6978)
+++ branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OFF_SRC/nemogcm.F90	(revision 6979)
@@ -19,5 +19,5 @@
    USE domcfg          ! domain configuration               (dom_cfg routine)
    USE domain          ! domain initialization from coordinate & bathymetry (dom_init routine)
-   USE domrea          ! domain initialization from mesh_mask            (dom_init routine)
+!   USE domrea          ! domain initialization from mesh_mask            (dom_init routine)
    USE eosbn2          ! equation of state            (eos bn2 routine)
    !              ! ocean physics
@@ -54,4 +54,6 @@
    USE sbc_oce, ONLY: ln_rnf
    USE sbcrnf
+   USE usrdef_nam     ! user defined configuration
+
 
    IMPLICIT NONE
@@ -147,15 +149,17 @@
       INTEGER ::   ji            ! dummy loop indices
       INTEGER ::   ilocal_comm   ! local integer
-      INTEGER ::   ios
-      LOGICAL ::   llexist
-      CHARACTER(len=80), DIMENSION(16) ::   cltxt
+      INTEGER ::   ios, inum
+      REAL(wp) ::   ziglo, zjglo, zkglo, zperio   ! local scalars
+      CHARACTER(len=120), DIMENSION(30) ::   cltxt, cltxt2, clnam
       !!
       NAMELIST/namctl/ ln_ctl  , nn_print, nn_ictls, nn_ictle,   &
          &             nn_isplt, nn_jsplt, nn_jctls, nn_jctle,   &
-         &             nn_bench, nn_timing, nn_diacfl
-      NAMELIST/namcfg/ cp_cfg, cp_cfz, jp_cfg, jpidta, jpjdta, jpkdta, jpiglo, jpjglo, &
-         &             jperio, ln_use_jattr
-      !!----------------------------------------------------------------------
-      cltxt = ''
+         &             nn_timing, nn_diacfl
+
+      NAMELIST/namcfg/ ln_read_cfg, ln_write_cfg, cp_cfg, jp_cfg, ln_use_jattr
+      !!----------------------------------------------------------------------
+      cltxt  = ''
+      cltxt2 = ''
+      clnam  = ''  
       cxios_context = 'nemo'
       !
@@ -181,4 +185,24 @@
 904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcfg in configuration namelist', .TRUE. )   
 
+
+      !                             !--------------------------!
+      !                             !  Set global domain size  !   (control print return in cltxt2)
+      !                             !--------------------------!
+      IF( ln_read_cfg ) THEN              ! Read sizes in configuration "mesh_mask" file
+         CALL iom_open( 'domain_cfg', inum )
+         CALL iom_get( inum, 'jpiglo', ziglo  )   ;   jpiglo = INT( ziglo )
+         CALL iom_get( inum, 'jpjglo', zjglo  )   ;   jpjglo = INT( zjglo )
+         CALL iom_get( inum, 'jpkglo', zkglo  )   ;   jpkglo = INT( zkglo )
+         CALL iom_get( inum, 'jperio', zperio )   ;   jperio = INT( zperio )
+         CALL iom_close( inum )
+         WRITE(cltxt2(1),*) '~~~~~~~~~~   '      
+         WRITE(cltxt2(2),*) 'domain_cfg : domain size read in "domain_cfg" file : jp(i,j,k)glo = '
+         WRITE(cltxt2(3),*) '            ', jpiglo, jpjglo, jpkglo        
+         WRITE(cltxt2(1),*) '~~~~~~~~~~   lateral boudary type of the global domain jperio= ', jperio
+         !        
+      ELSE                                ! user-defined namelist
+         CALL usr_def_nam( cltxt2, clnam, jpiglo, jpjglo, jpkglo, jperio )
+      ENDIF
+      jpk    = jpkglo
       !
       !                             !--------------------------------------------!
@@ -206,4 +230,9 @@
          WRITE( numond, namctl )
          WRITE( numond, namcfg )
+         IF( .NOT.ln_read_cfg ) THEN
+            DO ji = 1, SIZE(clnam)
+               IF( TRIM(clnam (ji)) /= '' )   WRITE(numond, * ) clnam(ji)    ! namusr_def print
+            END DO
+         ENDIF
       ENDIF
 
@@ -225,5 +254,4 @@
       jpi = ( jpiglo-2*jpreci + (jpni-1) ) / jpni + 2*jpreci   ! first  dim.
       jpj = ( jpjglo-2*jprecj + (jpnj-1) ) / jpnj + 2*jprecj   ! second dim.
-      jpk = jpkdta                                             ! third dim
       jpim1 = jpi-1                                            ! inner domain indices
       jpjm1 = jpj-1                                            !   "           "
@@ -274,11 +302,8 @@
                             CALL     eos_init   ! Equation of state
       IF( lk_c1d        )   CALL     c1d_init   ! 1D column configuration
+
                             CALL     dom_cfg    ! Domain configuration
-      !
-      INQUIRE( FILE='coordinates.nc', EXIST = llexist )   ! Check if coordinate file exist
-      !
-      IF( llexist )  THEN  ;  CALL  dom_init   !  compute the grid from coordinates and bathymetry
-      ELSE                 ;  CALL  dom_rea    !  read grid from the meskmask
-      ENDIF
+                            CALL     dom_init   ! Domain
+
                             CALL  istate_init   ! ocean initial state (Dynamics and tracers)
 
@@ -315,12 +340,12 @@
       !!                     ***  ROUTINE nemo_ctl  ***
       !!
-      !! ** Purpose :   control print setting 
+      !! ** Purpose :   control print setting
       !!
       !! ** Method  : - print namctl information and check some consistencies
       !!----------------------------------------------------------------------
       !
-      IF(lwp) THEN                  ! Parameter print
+      IF(lwp) THEN                  ! control print
          WRITE(numout,*)
-         WRITE(numout,*) 'nemo_flg: Control prints & Benchmark'
+         WRITE(numout,*) 'nemo_ctl: Control prints'
          WRITE(numout,*) '~~~~~~~ '
          WRITE(numout,*) '   Namelist namctl'
@@ -333,5 +358,5 @@
          WRITE(numout,*) '      number of proc. following i     nn_isplt   = ', nn_isplt
          WRITE(numout,*) '      number of proc. following j     nn_jsplt   = ', nn_jsplt
-         WRITE(numout,*) '      benchmark parameter (0/1)       nn_bench   = ', nn_bench
+         WRITE(numout,*) '      timing activated    (0/1)       nn_timing  = ', nn_timing
       ENDIF
       !
@@ -343,13 +368,13 @@
       isplt     = nn_isplt
       jsplt     = nn_jsplt
-      nbench    = nn_bench
-     IF(lwp) THEN                  ! control print
+
+      IF(lwp) THEN                  ! control print
          WRITE(numout,*)
          WRITE(numout,*) 'namcfg  : configuration initialization through namelist read'
          WRITE(numout,*) '~~~~~~~ '
          WRITE(numout,*) '   Namelist namcfg'
-         WRITE(numout,*) '      configuration name              cp_cfg      = ', TRIM(cp_cfg)
-         WRITE(numout,*) '      configuration resolution        jp_cfg      = ', jp_cfg
-         WRITE(numout,*) '      lateral cond. type (between 0 and 6) jperio = ', jperio   
+         WRITE(numout,*) '      read configuration definition files          ln_read_cfg = ', ln_read_cfg
+         WRITE(numout,*) '      configuration name                               cp_cfg  = ', TRIM(cp_cfg)
+         WRITE(numout,*) '      configuration resolution                         jp_cfg  = ', jp_cfg
          WRITE(numout,*) '      use file attribute if exists as i/p j-start ln_use_jattr = ', ln_use_jattr
       ENDIF
@@ -357,6 +382,6 @@
       !
       IF( ln_ctl ) THEN                 ! sub-domain area indices for the control prints
-         IF( lk_mpp ) THEN
-            isplt = jpni   ;   jsplt = jpnj   ;   ijsplt = jpni*jpnj   ! the domain is forced to the real splitted domain
+         IF( lk_mpp .AND. jpnij > 1 ) THEN
+            isplt = jpni   ;   jsplt = jpnj   ;   ijsplt = jpni*jpnj   ! the domain is forced to the real split domain
          ELSE
             IF( isplt == 1 .AND. jsplt == 1  ) THEN
@@ -393,8 +418,4 @@
       ENDIF
       !
-      IF( lk_c1d .AND. .NOT.lk_iomput )   CALL ctl_stop( 'nemo_ctl: The 1D configuration must be used ',   &
-         &                                               'with the IOM Input/Output manager. '        ,   &
-         &                                               'Compile with key_iomput enabled' )
-      !
       IF( 1_wp /= SIGN(1._wp,-0._wp)  )   CALL ctl_stop( 'nemo_ctl: The intrinsec SIGN function follows ',  &
          &                                               'f2003 standard. '                              ,  &
@@ -419,4 +440,6 @@
       IF( numnam_cfg /= -1 )   CLOSE( numnam_cfg )   ! oce configuration namelist
       IF( numout     /=  6 )   CLOSE( numout     )   ! standard model output file
+      IF( lwm.AND.numond  /= -1 )   CLOSE( numond          )   ! oce output namelist
+
       numout = 6                                     ! redefine numout in case it is used after this point...
       !
Index: branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90
===================================================================
--- branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90	(revision 6978)
+++ branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90	(revision 6979)
@@ -69,5 +69,5 @@
       !!              - 1D configuration, move Coriolis, u and v at T-point
       !!----------------------------------------------------------------------
-      INTEGER ::   ji, jj, jk   ! dummy loop indices
+      INTEGER ::   ji, jj, jk, ik   ! dummy loop indices
       INTEGER ::   iconf = 0    ! local integers
       CHARACTER (len=64) ::   cform = "(A12, 3(A13, I7))" 
@@ -107,6 +107,7 @@
       !
       DO jj = 1, jpj                   ! depth of the iceshelves
-         DO ji = 1, jpj
-            risfdep(ji,jj) = gdepw_0(ji,jj,mikt(ji,jj))
+         DO ji = 1, jpi
+            ik = mikt(ji,jj)
+            risfdep(ji,jj) = gdepw_0(ji,jj,ik)
          END DO
       END DO
Index: branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OPA_SRC/USR/usrdef_fmask.F90
===================================================================
--- branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OPA_SRC/USR/usrdef_fmask.F90	(revision 6978)
+++ branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OPA_SRC/USR/usrdef_fmask.F90	(revision 6979)
@@ -51,5 +51,5 @@
       !!                       in some user defined straits
       !!----------------------------------------------------------------------
-      CHARACTER(len=1)          , INTENT(in   ) ::   cd_cfg   ! configuration name
+      CHARACTER(len=*)          , INTENT(in   ) ::   cd_cfg   ! configuration name
       INTEGER                   , INTENT(in   ) ::   kcfg     ! configuration identifier 
       REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pfmsk    ! Ocean/Land f-point mask including lateral boundary cond.
@@ -60,5 +60,5 @@
       !!----------------------------------------------------------------------
       !
-      IF( cd_cfg == "orca" ) THEN      !==  ORCA Configurations  ==!
+      IF( TRIM( cd_cfg ) == "orca" ) THEN      !==  ORCA Configurations  ==!
          !
          SELECT CASE ( kcfg )
