Index: utils/tools_AGRIF_CMEMS_2020/DOMAINcfg/src/lib_mpp.F90
===================================================================
--- utils/tools_AGRIF_CMEMS_2020/DOMAINcfg/src/lib_mpp.f90	(revision 10082)
+++ utils/tools_AGRIF_CMEMS_2020/DOMAINcfg/src/lib_mpp.F90	(revision 10727)
@@ -8,5 +8,5 @@
    !!            8.0  !  1998  (M. Imbard, J. Escobar, L. Colombet ) SHMEM and MPI
    !!                 !  1998  (J.M. Molines) Open boundary conditions
-   !!   NEMO     1.0  !  2003  (J.-M. Molines, G. Madec)  F90, free form
+   !!   NEMO     1.0  !  2003  (J.M. Molines, G. Madec)  F90, free form
    !!                 !  2003  (J.M. Molines) add mpp_ini_north(_3d,_2d)
    !!             -   !  2004  (R. Bourdalle Badie)  isend option in mpi
@@ -19,10 +19,10 @@
    !!            3.2  !  2009  (O. Marti)    add mpp_ini_znl
    !!            4.0  !  2011  (G. Madec)  move ctl_ routines from in_out_manager
-   !!            3.5  !  2012  (S.Mocavero, I. Epicoco) Add 'mpp_lnk_bdy_3d', 'mpp_lnk_obc_3d', 
-   !!                          'mpp_lnk_bdy_2d' and 'mpp_lnk_obc_2d' routines and update
-   !!                          the mppobc routine to optimize the BDY and OBC communications
-   !!            3.5  !  2013  ( C. Ethe, G. Madec ) message passing arrays as local variables 
+   !!            3.5  !  2012  (S.Mocavero, I. Epicoco) Add mpp_lnk_bdy_3d/2d routines to optimize the BDY comm.
+   !!            3.5  !  2013  (C. Ethe, G. Madec)  message passing arrays as local variables 
    !!            3.5  !  2013  (S.Mocavero, I.Epicoco - CMCC) north fold optimizations
-   !!            3.6  !  2015 (O. Tintó and M. Castrillo - BSC) Added 'mpp_lnk_2d_multiple', 'mpp_lbc_north_2d_multiple', 'mpp_max_multiple' 
+   !!            3.6  !  2015  (O. Tintó and M. Castrillo - BSC) Added '_multiple' case for 2D lbc and max
+   !!            4.0  !  2017  (G. Madec) automatique allocation of array argument (use any 3rd dimension)
+   !!             -   !  2017  (G. Madec) create generic.h90 files to generate all lbc and north fold routines
    !!----------------------------------------------------------------------
 
@@ -34,5 +34,5 @@
    !!   get_unit      : give the index of an unused logical unit
    !!----------------------------------------------------------------------
-
+#if   defined key_mpp_mpi
    !!----------------------------------------------------------------------
    !!   'key_mpp_mpi'             MPI massively parallel processing library
@@ -41,9 +41,7 @@
    !!   mynode        : indentify the processor unit
    !!   mpp_lnk       : interface (defined in lbclnk) for message passing of 2d or 3d arrays (mpp_lnk_2d, mpp_lnk_3d)
-   !!   mpp_lnk_3d_gather :  Message passing manadgement for two 3D arrays
-   !!   mpp_lnk_e     : interface (defined in lbclnk) for message passing of 2d array with extra halo (mpp_lnk_2d_e)
    !!   mpp_lnk_icb   : interface for message passing of 2d arrays with extra halo for icebergs (mpp_lnk_2d_icb)
    !!   mpprecv       :
-   !!   mppsend       :   SUBROUTINE mpp_ini_znl
+   !!   mppsend       :
    !!   mppscatter    :
    !!   mppgather     :
@@ -56,36 +54,38 @@
    !!   mppstop       :
    !!   mpp_ini_north : initialisation of north fold
-   !!   mpp_lbc_north : north fold processors gathering
-   !!   mpp_lbc_north_e : variant of mpp_lbc_north for extra outer halo
-   !!   mpp_lbc_north_icb : variant of mpp_lbc_north for extra outer halo with icebergs
+   !!   mpp_lbc_north_icb : alternative to mpp_nfd for extra outer halo with icebergs
    !!----------------------------------------------------------------------
    USE dom_oce        ! ocean space and time domain
    USE lbcnfd         ! north fold treatment
    USE in_out_manager ! I/O manager
-   USE wrk_nemo       ! work arrays
 
    IMPLICIT NONE
    PRIVATE
-   
+
+   INTERFACE mpp_nfd
+      MODULE PROCEDURE   mpp_nfd_2d    , mpp_nfd_3d    , mpp_nfd_4d
+      MODULE PROCEDURE   mpp_nfd_2d_ptr, mpp_nfd_3d_ptr, mpp_nfd_4d_ptr
+   END INTERFACE
+
+   ! Interface associated to the mpp_lnk_... routines is defined in lbclnk
+   PUBLIC   mpp_lnk_2d    , mpp_lnk_3d    , mpp_lnk_4d
+   PUBLIC   mpp_lnk_2d_ptr, mpp_lnk_3d_ptr, mpp_lnk_4d_ptr
+   !
+!!gm  this should be useless
+   PUBLIC   mpp_nfd_2d    , mpp_nfd_3d    , mpp_nfd_4d
+   PUBLIC   mpp_nfd_2d_ptr, mpp_nfd_3d_ptr, mpp_nfd_4d_ptr
+!!gm end
+   !
    PUBLIC   ctl_stop, ctl_warn, get_unit, ctl_opn, ctl_nam
    PUBLIC   mynode, mppstop, mppsync, mpp_comm_free
-   PUBLIC   mpp_ini_north, mpp_lbc_north, mpp_lbc_north_e
+   PUBLIC   mpp_ini_north
+   PUBLIC   mpp_lnk_2d_icb
+   PUBLIC   mpp_lbc_north_icb
    PUBLIC   mpp_min, mpp_max, mpp_sum, mpp_minloc, mpp_maxloc
-   PUBLIC   mpp_max_multiple
-   PUBLIC   mpp_lnk_3d, mpp_lnk_3d_gather, mpp_lnk_2d, mpp_lnk_2d_e
-   PUBLIC   mpp_lnk_2d_9 , mpp_lnk_2d_multiple 
-   PUBLIC   mpp_lnk_sum_3d, mpp_lnk_sum_2d
+   PUBLIC   mpp_delay_max, mpp_delay_sum, mpp_delay_rcv
    PUBLIC   mppscatter, mppgather
-   PUBLIC   mpp_ini_ice, mpp_ini_znl
-   PUBLIC   mppsize
+   PUBLIC   mpp_ini_znl
    PUBLIC   mppsend, mpprecv                          ! needed by TAM and ICB routines
-   PUBLIC   mpp_lnk_bdy_2d, mpp_lnk_bdy_3d
-   PUBLIC   mpp_lbc_north_icb, mpp_lnk_2d_icb
-   PUBLIC   mpprank
-
-   TYPE arrayptr
-      REAL , DIMENSION (:,:),  POINTER :: pt2d
-   END TYPE arrayptr
-   PUBLIC   arrayptr
+   PUBLIC   mpp_lnk_bdy_2d, mpp_lnk_bdy_3d, mpp_lnk_bdy_4d
    
    !! * Interfaces
@@ -101,8 +101,5 @@
    INTERFACE mpp_sum
       MODULE PROCEDURE mppsum_a_int, mppsum_int, mppsum_a_real, mppsum_real,   &
-                       mppsum_realdd, mppsum_a_realdd
-   END INTERFACE
-   INTERFACE mpp_lbc_north
-      MODULE PROCEDURE mpp_lbc_north_3d, mpp_lbc_north_2d
+         &             mppsum_realdd, mppsum_a_realdd
    END INTERFACE
    INTERFACE mpp_minloc
@@ -113,8 +110,4 @@
    END INTERFACE
 
-   INTERFACE mpp_max_multiple
-      MODULE PROCEDURE mppmax_real_multiple
-   END INTERFACE
-
    !! ========================= !!
    !!  MPI  variable definition !!
@@ -128,56 +121,79 @@
    INTEGER, PARAMETER         ::   nprocmax = 2**10   ! maximun dimension (required to be a power of 2)
 
-   INTEGER ::   mppsize        ! number of process
-   INTEGER ::   mpprank        ! process number  [ 0 - size-1 ]
+   INTEGER, PUBLIC ::   mppsize        ! number of process
+   INTEGER, PUBLIC ::   mpprank        ! process number  [ 0 - size-1 ]
 !$AGRIF_DO_NOT_TREAT
-   INTEGER, PUBLIC ::   mpi_comm_opa   ! opa local communicator
+   INTEGER, PUBLIC ::   mpi_comm_oce   ! opa local communicator
 !$AGRIF_END_DO_NOT_TREAT
 
    INTEGER :: MPI_SUMDD
-
-   ! variables used in case of sea-ice
-   INTEGER, PUBLIC ::   ncomm_ice       !: communicator made by the processors with sea-ice (public so that it can be freed in limthd)
-   INTEGER ::   ngrp_iworld     !  group ID for the world processors (for rheology)
-   INTEGER ::   ngrp_ice        !  group ID for the ice processors (for rheology)
-   INTEGER ::   ndim_rank_ice   !  number of 'ice' processors
-   INTEGER ::   n_ice_root      !  number (in the comm_ice) of proc 0 in the ice comm
-   INTEGER, DIMENSION(:), ALLOCATABLE, SAVE ::   nrank_ice     ! dimension ndim_rank_ice
 
    ! variables used for zonal integration
    INTEGER, PUBLIC ::   ncomm_znl       !: communicator made by the processors on the same zonal average
-   LOGICAL, PUBLIC ::   l_znl_root      ! True on the 'left'most processor on the same row
-   INTEGER ::   ngrp_znl        ! group ID for the znl processors
-   INTEGER ::   ndim_rank_znl   ! number of processors on the same zonal average
+   LOGICAL, PUBLIC ::   l_znl_root      !: True on the 'left'most processor on the same row
+   INTEGER         ::   ngrp_znl        !  group ID for the znl processors
+   INTEGER         ::   ndim_rank_znl   !  number of processors on the same zonal average
    INTEGER, DIMENSION(:), ALLOCATABLE, SAVE ::   nrank_znl  ! dimension ndim_rank_znl, number of the procs into the same znl domain
 
    ! North fold condition in mpp_mpi with jpni > 1 (PUBLIC for TAM)
-   INTEGER, PUBLIC ::   ngrp_world        ! group ID for the world processors
-   INTEGER, PUBLIC ::   ngrp_opa          ! group ID for the opa processors
-   INTEGER, PUBLIC ::   ngrp_north        ! group ID for the northern processors (to be fold)
-   INTEGER, PUBLIC ::   ncomm_north       ! communicator made by the processors belonging to ngrp_north
-   INTEGER, PUBLIC ::   ndim_rank_north   ! number of 'sea' processor in the northern line (can be /= jpni !)
-   INTEGER, PUBLIC ::   njmppmax          ! value of njmpp for the processors of the northern line
-   INTEGER, PUBLIC ::   north_root        ! number (in the comm_opa) of proc 0 in the northern comm
-   INTEGER, DIMENSION(:), ALLOCATABLE, SAVE, PUBLIC ::   nrank_north   ! dimension ndim_rank_north
+   INTEGER, PUBLIC ::   ngrp_world        !: group ID for the world processors
+   INTEGER, PUBLIC ::   ngrp_opa          !: group ID for the opa processors
+   INTEGER, PUBLIC ::   ngrp_north        !: group ID for the northern processors (to be fold)
+   INTEGER, PUBLIC ::   ncomm_north       !: communicator made by the processors belonging to ngrp_north
+   INTEGER, PUBLIC ::   ndim_rank_north   !: number of 'sea' processor in the northern line (can be /= jpni !)
+   INTEGER, PUBLIC ::   njmppmax          !: value of njmpp for the processors of the northern line
+   INTEGER, PUBLIC ::   north_root        !: number (in the comm_opa) of proc 0 in the northern comm
+   INTEGER, PUBLIC, DIMENSION(:), ALLOCATABLE, SAVE ::   nrank_north   !: dimension ndim_rank_north
 
    ! Type of send : standard, buffered, immediate
-   CHARACTER(len=1), PUBLIC ::   cn_mpi_send   ! type od mpi send/recieve (S=standard, B=bsend, I=isend)
-   LOGICAL, PUBLIC          ::   l_isend = .FALSE.   ! isend use indicator (T if cn_mpi_send='I')
-   INTEGER, PUBLIC          ::   nn_buffer     ! size of the buffer in case of mpi_bsend
-
-   REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE :: tampon  ! buffer in case of bsend
-
-   LOGICAL, PUBLIC                                  ::   ln_nnogather       ! namelist control of northfold comms
-   LOGICAL, PUBLIC                                  ::   l_north_nogather = .FALSE.  ! internal control of northfold comms
-   INTEGER, PUBLIC                                  ::   ityp
+   CHARACTER(len=1), PUBLIC ::   cn_mpi_send        !: type od mpi send/recieve (S=standard, B=bsend, I=isend)
+   LOGICAL         , PUBLIC ::   l_isend = .FALSE.  !: isend use indicator (T if cn_mpi_send='I')
+   INTEGER         , PUBLIC ::   nn_buffer          !: size of the buffer in case of mpi_bsend
+
+   ! Communications summary report
+   CHARACTER(len=400), DIMENSION(:), ALLOCATABLE ::   crname_lbc                   !: names of lbc_lnk calling routines
+   CHARACTER(len=400), DIMENSION(:), ALLOCATABLE ::   crname_glb                   !: names of global comm calling routines
+   CHARACTER(len=400), DIMENSION(:), ALLOCATABLE ::   crname_dlg                   !: names of delayed global comm calling routines
+   INTEGER, PUBLIC                               ::   ncom_stp = 0                 !: copy of time step # istp
+   INTEGER, PUBLIC                               ::   ncom_fsbc = 1                !: copy of sbc time step # nn_fsbc
+   INTEGER, PUBLIC                               ::   ncom_dttrc = 1               !: copy of top time step # nn_dttrc
+   INTEGER, PUBLIC                               ::   ncom_freq                    !: frequency of comm diagnostic
+   INTEGER, PUBLIC , DIMENSION(:,:), ALLOCATABLE ::   ncomm_sequence               !: size of communicated arrays (halos)
+   INTEGER, PARAMETER, PUBLIC                    ::   ncom_rec_max = 3000          !: max number of communication record
+   INTEGER, PUBLIC                               ::   n_sequence_lbc = 0           !: # of communicated arraysvia lbc
+   INTEGER, PUBLIC                               ::   n_sequence_glb = 0           !: # of global communications
+   INTEGER, PUBLIC                               ::   n_sequence_dlg = 0           !: # of delayed global communications
+   INTEGER, PUBLIC                               ::   numcom = -1                  !: logical unit for communicaton report
+   LOGICAL, PUBLIC                               ::   l_full_nf_update = .TRUE.    !: logical for a full (2lines) update of bc at North fold report
+   INTEGER,                    PARAMETER, PUBLIC ::   nbdelay = 2       !: number of delayed operations
+   !: name (used as id) of allreduce-delayed operations
+   ! Warning: we must use the same character length in an array constructor (at least for gcc compiler)
+   CHARACTER(len=32), DIMENSION(nbdelay), PUBLIC ::   c_delaylist = (/ 'cflice', 'fwb   ' /)
+   !: component name where the allreduce-delayed operation is performed
+   CHARACTER(len=3),  DIMENSION(nbdelay), PUBLIC ::   c_delaycpnt = (/ 'ICE'   , 'OCE' /)
+   TYPE, PUBLIC ::   DELAYARR
+      REAL(   wp), POINTER, DIMENSION(:) ::  z1d => NULL()
+      COMPLEX(wp), POINTER, DIMENSION(:) ::  y1d => NULL()
+   END TYPE DELAYARR
+   TYPE( DELAYARR ), DIMENSION(nbdelay), PUBLIC  ::   todelay              
+   INTEGER,          DIMENSION(nbdelay), PUBLIC  ::   ndelayid = -1     !: mpi request id of the delayed operations
+
+   ! timing summary report
+   REAL(wp), DIMENSION(2), PUBLIC ::  waiting_time = 0._wp
+   REAL(wp)              , PUBLIC ::  compute_time = 0._wp, elapsed_time = 0._wp
+   
+   REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE ::   tampon   ! buffer in case of bsend
+
+   LOGICAL, PUBLIC ::   ln_nnogather                !: namelist control of northfold comms
+   LOGICAL, PUBLIC ::   l_north_nogather = .FALSE.  !: internal control of northfold comms
+
    !!----------------------------------------------------------------------
    !! NEMO/OCE 4.0 , NEMO Consortium (2018)
-   !! $Id: lib_mpp.F90 6490 2016-04-20 14:55:58Z mcastril $
-   !! Software governed by the CeCILL licence     (./LICENSE)
+   !! $Id: lib_mpp.F90 10538 2019-01-17 10:41:10Z clem $
+   !! Software governed by the CeCILL license (see ./LICENSE)
    !!----------------------------------------------------------------------
 CONTAINS
 
-
-   FUNCTION mynode( ldtxt, ldname, kumnam_ref , kumnam_cfg , kumond , kstop, localComm )
+   FUNCTION mynode( ldtxt, ldname, kumnam_ref, kumnam_cfg, kumond, kstop, localComm )
       !!----------------------------------------------------------------------
       !!                  ***  routine mynode  ***
@@ -196,5 +212,5 @@
       LOGICAL ::   mpi_was_called
       !
-      NAMELIST/nammpp/ cn_mpi_send, nn_buffer, jpni, jpnj, jpnij, ln_nnogather
+      NAMELIST/nammpp/ cn_mpi_send, nn_buffer, jpni, jpnj, ln_nnogather
       !!----------------------------------------------------------------------
       !
@@ -204,39 +220,23 @@
       WRITE(ldtxt(ii),*) '~~~~~~ '                                                        ;   ii = ii + 1
       !
-
       REWIND( kumnam_ref )              ! Namelist nammpp in reference namelist: mpi variables
       READ  ( kumnam_ref, nammpp, IOSTAT = ios, ERR = 901)
-901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nammpp in reference namelist', lwp )
-
+901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nammpp in reference namelist', lwp )
+      !
       REWIND( kumnam_cfg )              ! Namelist nammpp in configuration namelist: mpi variables
       READ  ( kumnam_cfg, nammpp, IOSTAT = ios, ERR = 902 )
-902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nammpp in configuration namelist', lwp )
-
+902   IF( ios >  0 )   CALL ctl_nam ( ios , 'nammpp in configuration namelist', lwp )
+      !
       !                              ! control print
       WRITE(ldtxt(ii),*) '   Namelist nammpp'                                             ;   ii = ii + 1
       WRITE(ldtxt(ii),*) '      mpi send type          cn_mpi_send = ', cn_mpi_send       ;   ii = ii + 1
       WRITE(ldtxt(ii),*) '      size exported buffer   nn_buffer   = ', nn_buffer,' bytes';   ii = ii + 1
-
-
-
-
-
-
-
-
-
-      IF(jpnij < 1)THEN
-         ! If jpnij is not specified in namelist then we calculate it - this
-         ! means there will be no land cutting out.
-         jpnij = jpni * jpnj
-      END IF
-
-      IF( (jpni < 1) .OR. (jpnj < 1) )THEN
-         WRITE(ldtxt(ii),*) '      jpni, jpnj and jpnij will be calculated automatically' ;   ii = ii + 1
+      !
+      IF( jpni < 1 .OR. jpnj < 1  ) THEN
+         WRITE(ldtxt(ii),*) '      jpni and jpnj will be calculated automatically' ;   ii = ii + 1
       ELSE
          WRITE(ldtxt(ii),*) '      processor grid extent in i         jpni = ',jpni       ;   ii = ii + 1
          WRITE(ldtxt(ii),*) '      processor grid extent in j         jpnj = ',jpnj       ;   ii = ii + 1
-         WRITE(ldtxt(ii),*) '      number of local domains           jpnij = ',jpnij      ;   ii = ii + 1
-      END IF
+      ENDIF
 
       WRITE(ldtxt(ii),*) '      avoid use of mpi_allgather at the north fold  ln_nnogather = ', ln_nnogather  ; ii = ii + 1
@@ -259,5 +259,5 @@
          CASE ( 'B' )                ! Buffer mpi send (blocking)
             WRITE(ldtxt(ii),*) '           Buffer blocking mpi send (bsend)'              ;   ii = ii + 1
-            IF( Agrif_Root() )   CALL mpi_init_opa( ldtxt, ii, ierr )
+            IF( Agrif_Root() )   CALL mpi_init_oce( ldtxt, ii, ierr )
          CASE ( 'I' )                ! Immediate mpi send (non-blocking send)
             WRITE(ldtxt(ii),*) '           Immediate non-blocking send (isend)'           ;   ii = ii + 1
@@ -268,5 +268,7 @@
             kstop = kstop + 1
          END SELECT
-      ELSE IF ( PRESENT(localComm) .and. .not. mpi_was_called ) THEN
+         !
+      ELSEIF ( PRESENT(localComm) .AND. .NOT. mpi_was_called ) THEN
+         WRITE(ldtxt(ii),cform_err)                                                    ;   ii = ii + 1
          WRITE(ldtxt(ii),*) ' lib_mpp: You cannot provide a local communicator '          ;   ii = ii + 1
          WRITE(ldtxt(ii),*) '          without calling MPI_Init before ! '                ;   ii = ii + 1
@@ -279,5 +281,5 @@
          CASE ( 'B' )                ! Buffer mpi send (blocking)
             WRITE(ldtxt(ii),*) '           Buffer blocking mpi send (bsend)'              ;   ii = ii + 1
-            IF( Agrif_Root() )   CALL mpi_init_opa( ldtxt, ii, ierr )
+            IF( Agrif_Root() )   CALL mpi_init_oce( ldtxt, ii, ierr )
          CASE ( 'I' )                ! Immediate mpi send (non-blocking send)
             WRITE(ldtxt(ii),*) '           Immediate non-blocking send (isend)'           ;   ii = ii + 1
@@ -294,8 +296,8 @@
       IF( PRESENT(localComm) ) THEN
          IF( Agrif_Root() ) THEN
-            mpi_comm_opa = localComm
+            mpi_comm_oce = localComm
          ENDIF
       ELSE
-         CALL mpi_comm_dup( mpi_comm_world, mpi_comm_opa, code)
+         CALL mpi_comm_dup( mpi_comm_world, mpi_comm_oce, code)
          IF( code /= MPI_SUCCESS ) THEN
             DO ji = 1, SIZE(ldtxt)
@@ -308,14 +310,14 @@
       ENDIF
 
-
-
-
-
-
-
-
-
-      CALL mpi_comm_rank( mpi_comm_opa, mpprank, ierr )
-      CALL mpi_comm_size( mpi_comm_opa, mppsize, ierr )
+#if defined key_agrif
+      IF( Agrif_Root() ) THEN
+         CALL Agrif_MPI_Init(mpi_comm_oce)
+      ELSE
+         CALL Agrif_MPI_set_grid_comm(mpi_comm_oce)
+      ENDIF
+#endif
+
+      CALL mpi_comm_rank( mpi_comm_oce, mpprank, ierr )
+      CALL mpi_comm_size( mpi_comm_oce, mppsize, ierr )
       mynode = mpprank
 
@@ -329,1420 +331,153 @@
    END FUNCTION mynode
 
-
-   SUBROUTINE mpp_lnk_3d( ptab, cd_type, psgn, cd_mpp, pval )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mpp_lnk_3d  ***
-      !!
-      !! ** Purpose :   Message passing manadgement
-      !!
-      !! ** Method  :   Use mppsend and mpprecv function for passing mask
-      !!      between processors following neighboring subdomains.
-      !!            domain parameters
-      !!                    nlci   : first dimension of the local subdomain
-      !!                    nlcj   : second dimension of the local subdomain
-      !!                    nbondi : mark for "east-west local boundary"
-      !!                    nbondj : mark for "north-south local boundary"
-      !!                    noea   : number for local neighboring processors
-      !!                    nowe   : number for local neighboring processors
-      !!                    noso   : number for local neighboring processors
-      !!                    nono   : number for local neighboring processors
-      !!
-      !! ** Action  :   ptab with update value at its periphery
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptab     ! 3D array on which the boundary condition is applied
-      CHARACTER(len=1)                , INTENT(in   ) ::   cd_type  ! define the nature of ptab array grid-points
-      !                                                             ! = T , U , V , F , W points
-      REAL(wp)                        , INTENT(in   ) ::   psgn     ! =-1 the sign change across the north fold boundary
-      !                                                             ! =  1. , the sign is kept
-      CHARACTER(len=3), OPTIONAL      , INTENT(in   ) ::   cd_mpp   ! fill the overlap area only
-      REAL(wp)        , OPTIONAL      , INTENT(in   ) ::   pval     ! background value (used at closed boundaries)
-      !
-      INTEGER  ::   ji, jj, jk, jl             ! dummy loop indices
-      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
-      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
-      REAL(wp) ::   zland
-      INTEGER , DIMENSION(MPI_STATUS_SIZE)      ::   ml_stat        ! for key_mpi_isend
-      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::   zt3ns, zt3sn   ! 3d for north-south & south-north
-      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::   zt3ew, zt3we   ! 3d for east-west & west-east
-      !!----------------------------------------------------------------------
-      
-      ALLOCATE( zt3ns(jpi,jprecj,jpk,2), zt3sn(jpi,jprecj,jpk,2),   &
-         &      zt3ew(jpj,jpreci,jpk,2), zt3we(jpj,jpreci,jpk,2)  )
-
-      !
-      IF( PRESENT( pval ) ) THEN   ;   zland = pval      ! set land value
-      ELSE                         ;   zland = 0._wp     ! zero by default
-      ENDIF
-
-      ! 1. standard boundary treatment
-      ! ------------------------------
-      IF( PRESENT( cd_mpp ) ) THEN      ! only fill added line/raw with existing values
-         !
-         ! WARNING ptab is defined only between nld and nle
-         DO jk = 1, jpk
-            DO jj = nlcj+1, jpj                 ! added line(s)   (inner only)
-               ptab(nldi  :nlei  , jj          ,jk) = ptab(nldi:nlei,     nlej,jk)
-               ptab(1     :nldi-1, jj          ,jk) = ptab(nldi     ,     nlej,jk)
-               ptab(nlei+1:nlci  , jj          ,jk) = ptab(     nlei,     nlej,jk)
-            END DO
-            DO ji = nlci+1, jpi                 ! added column(s) (full)
-               ptab(ji           ,nldj  :nlej  ,jk) = ptab(     nlei,nldj:nlej,jk)
-               ptab(ji           ,1     :nldj-1,jk) = ptab(     nlei,nldj     ,jk)
-               ptab(ji           ,nlej+1:jpj   ,jk) = ptab(     nlei,     nlej,jk)
-            END DO
-         END DO
-         !
-      ELSE                              ! standard close or cyclic treatment
-         !
-         !                                   ! East-West boundaries
-         !                                        !* Cyclic east-west
-         IF( nbondi == 2 .AND. (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
-            ptab( 1 ,:,:) = ptab(jpim1,:,:)
-            ptab(jpi,:,:) = ptab(  2  ,:,:)
-         ELSE                                     !* closed
-            IF( .NOT. cd_type == 'F' )   ptab(     1       :jpreci,:,:) = zland    ! south except F-point
-                                         ptab(nlci-jpreci+1:jpi   ,:,:) = zland    ! north
-         ENDIF
-         !                                   ! North-South boundaries (always closed)
-         IF( .NOT. cd_type == 'F' )   ptab(:,     1       :jprecj,:) = zland       ! south except F-point
-                                      ptab(:,nlcj-jprecj+1:jpj   ,:) = zland       ! north
-         !
-      ENDIF
-
-      ! 2. East and west directions exchange
-      ! ------------------------------------
-      ! we play with the neigbours AND the row number because of the periodicity
-      !
-      SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
-      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
-         iihom = nlci-nreci
-         DO jl = 1, jpreci
-            zt3ew(:,jl,:,1) = ptab(jpreci+jl,:,:)
-            zt3we(:,jl,:,1) = ptab(iihom +jl,:,:)
-         END DO
-      END SELECT
-      !
-      !                           ! Migrations
-      imigr = jpreci * jpj * jpk
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         CALL mppsend( 2, zt3we(1,1,1,1), imigr, noea, ml_req1 )
-         CALL mpprecv( 1, zt3ew(1,1,1,2), imigr, noea )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      CASE ( 0 )
-         CALL mppsend( 1, zt3ew(1,1,1,1), imigr, nowe, ml_req1 )
-         CALL mppsend( 2, zt3we(1,1,1,1), imigr, noea, ml_req2 )
-         CALL mpprecv( 1, zt3ew(1,1,1,2), imigr, noea )
-         CALL mpprecv( 2, zt3we(1,1,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
-      CASE ( 1 )
-         CALL mppsend( 1, zt3ew(1,1,1,1), imigr, nowe, ml_req1 )
-         CALL mpprecv( 2, zt3we(1,1,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      iihom = nlci-jpreci
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         DO jl = 1, jpreci
-            ptab(iihom+jl,:,:) = zt3ew(:,jl,:,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jpreci
-            ptab(jl      ,:,:) = zt3we(:,jl,:,2)
-            ptab(iihom+jl,:,:) = zt3ew(:,jl,:,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jpreci
-            ptab(jl      ,:,:) = zt3we(:,jl,:,2)
-         END DO
-      END SELECT
-
-      ! 3. North and south directions
-      ! -----------------------------
-      ! always closed : we play only with the neigbours
-      !
-      IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
-         ijhom = nlcj-nrecj
-         DO jl = 1, jprecj
-            zt3sn(:,jl,:,1) = ptab(:,ijhom +jl,:)
-            zt3ns(:,jl,:,1) = ptab(:,jprecj+jl,:)
-         END DO
-      ENDIF
-      !
-      !                           ! Migrations
-      imigr = jprecj * jpi * jpk
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         CALL mppsend( 4, zt3sn(1,1,1,1), imigr, nono, ml_req1 )
-         CALL mpprecv( 3, zt3ns(1,1,1,2), imigr, nono )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      CASE ( 0 )
-         CALL mppsend( 3, zt3ns(1,1,1,1), imigr, noso, ml_req1 )
-         CALL mppsend( 4, zt3sn(1,1,1,1), imigr, nono, ml_req2 )
-         CALL mpprecv( 3, zt3ns(1,1,1,2), imigr, nono )
-         CALL mpprecv( 4, zt3sn(1,1,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
-      CASE ( 1 )
-         CALL mppsend( 3, zt3ns(1,1,1,1), imigr, noso, ml_req1 )
-         CALL mpprecv( 4, zt3sn(1,1,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      ijhom = nlcj-jprecj
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         DO jl = 1, jprecj
-            ptab(:,ijhom+jl,:) = zt3ns(:,jl,:,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jprecj
-            ptab(:,jl      ,:) = zt3sn(:,jl,:,2)
-            ptab(:,ijhom+jl,:) = zt3ns(:,jl,:,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jprecj
-            ptab(:,jl,:) = zt3sn(:,jl,:,2)
-         END DO
-      END SELECT
-
-      ! 4. north fold treatment
-      ! -----------------------
-      !
-      IF( npolj /= 0 .AND. .NOT. PRESENT(cd_mpp) ) THEN
-         !
-         SELECT CASE ( jpni )
-         CASE ( 1 )     ;   CALL lbc_nfd      ( ptab, cd_type, psgn )   ! only 1 northern proc, no mpp
-         CASE DEFAULT   ;   CALL mpp_lbc_north( ptab, cd_type, psgn )   ! for all northern procs.
-         END SELECT
-         !
-      ENDIF
-      !
-      DEALLOCATE( zt3ns, zt3sn, zt3ew, zt3we )
-      !
-   END SUBROUTINE mpp_lnk_3d
-
-
-   SUBROUTINE mpp_lnk_2d_multiple( pt2d_array , type_array , psgn_array , num_fields , cd_mpp, pval )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mpp_lnk_2d_multiple  ***
-      !!
-      !! ** Purpose :   Message passing management for multiple 2d arrays
-      !!
-      !! ** Method  :   Use mppsend and mpprecv function for passing mask
-      !!      between processors following neighboring subdomains.
-      !!            domain parameters
-      !!                    nlci   : first dimension of the local subdomain
-      !!                    nlcj   : second dimension of the local subdomain
-      !!                    nbondi : mark for "east-west local boundary"
-      !!                    nbondj : mark for "north-south local boundary"
-      !!                    noea   : number for local neighboring processors
-      !!                    nowe   : number for local neighboring processors
-      !!                    noso   : number for local neighboring processors
-      !!                    nono   : number for local neighboring processors
-      !!----------------------------------------------------------------------
-      CHARACTER(len=1), DIMENSION(:), INTENT(in   ) ::   type_array   ! define the nature of ptab array grid-points
-      !                                                               ! = T , U , V , F , W and I points
-      REAL(wp)        , DIMENSION(:), INTENT(in   ) ::   psgn_array   ! =-1 the sign change across the north fold boundary
-      !                                                               ! =  1. , the sign is kept
-      CHARACTER(len=3), OPTIONAL    , INTENT(in   ) ::   cd_mpp       ! fill the overlap area only
-      REAL(wp)        , OPTIONAL    , INTENT(in   ) ::   pval         ! background value (used at closed boundaries)
-      !!
-      INTEGER  ::   ji, jj, jl   ! dummy loop indices
-      INTEGER  ::   ii    !!MULTI SEND DUMMY LOOP INDICES
-      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
-      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
-      INTEGER :: num_fields
-      TYPE( arrayptr ), DIMENSION(:) :: pt2d_array
-      REAL(wp) ::   zland
-      INTEGER , DIMENSION(MPI_STATUS_SIZE)    ::   ml_stat       ! for key_mpi_isend
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  zt2ns, zt2sn   ! 2d for north-south & south-north
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  zt2ew, zt2we   ! 2d for east-west & west-east
-
-      !!----------------------------------------------------------------------
-      !
-      ALLOCATE( zt2ns(jpi,jprecj,2*num_fields), zt2sn(jpi,jprecj,2*num_fields),  &
-         &      zt2ew(jpj,jpreci,2*num_fields), zt2we(jpj,jpreci,2*num_fields)   )
-      !
-      IF( PRESENT( pval ) ) THEN   ;   zland = pval      ! set land value
-      ELSE                         ;   zland = 0._wp     ! zero by default
-      ENDIF
-
-      ! 1. standard boundary treatment
-      ! ------------------------------
-      !
-      !First Array
-      DO ii = 1 , num_fields
-         IF( PRESENT( cd_mpp ) ) THEN      ! only fill added line/raw with existing values
-            !
-            ! WARNING pt2d is defined only between nld and nle
-            DO jj = nlcj+1, jpj                 ! added line(s)   (inner only)
-               pt2d_array(ii)%pt2d(nldi  :nlei  , jj) = pt2d_array(ii)%pt2d(nldi:nlei, nlej)
-               pt2d_array(ii)%pt2d(1     :nldi-1, jj) = pt2d_array(ii)%pt2d(nldi     , nlej)
-               pt2d_array(ii)%pt2d(nlei+1:nlci  , jj) = pt2d_array(ii)%pt2d(     nlei, nlej) 
-            END DO
-            DO ji = nlci+1, jpi                 ! added column(s) (full)
-               pt2d_array(ii)%pt2d(ji, nldj  :nlej  ) = pt2d_array(ii)%pt2d(nlei, nldj:nlej)
-               pt2d_array(ii)%pt2d(ji, 1     :nldj-1) = pt2d_array(ii)%pt2d(nlei, nldj     )
-               pt2d_array(ii)%pt2d(ji, nlej+1:jpj   ) = pt2d_array(ii)%pt2d(nlei,      nlej)
-            END DO
-            !
-         ELSE                              ! standard close or cyclic treatment
-            !
-            !                                   ! East-West boundaries
-            IF( nbondi == 2 .AND.   &                ! Cyclic east-west
-               &    (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
-               pt2d_array(ii)%pt2d(  1  , : ) = pt2d_array(ii)%pt2d( jpim1, : )                                    ! west
-               pt2d_array(ii)%pt2d( jpi , : ) = pt2d_array(ii)%pt2d(   2  , : )                                    ! east
-            ELSE                                     ! closed
-               IF( .NOT. type_array(ii) == 'F' )   pt2d_array(ii)%pt2d(            1 : jpreci,:) = zland    ! south except F-point
-                                                   pt2d_array(ii)%pt2d(nlci-jpreci+1 : jpi   ,:) = zland    ! north
-            ENDIF
-            !                                   ! North-South boundaries (always closed)
-               IF( .NOT. type_array(ii) == 'F' )   pt2d_array(ii)%pt2d(:,             1:jprecj ) = zland    ! south except F-point
-                                                   pt2d_array(ii)%pt2d(:, nlcj-jprecj+1:jpj    ) = zland    ! north
-            !
-         ENDIF
-      END DO
-
-      ! 2. East and west directions exchange
-      ! ------------------------------------
-      ! we play with the neigbours AND the row number because of the periodicity
-      !
-      DO ii = 1 , num_fields
-         SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
-         CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
-            iihom = nlci-nreci
-            DO jl = 1, jpreci
-               zt2ew( : , jl , ii ) = pt2d_array(ii)%pt2d( jpreci+jl , : )
-               zt2we( : , jl , ii ) = pt2d_array(ii)%pt2d( iihom +jl , : )
-            END DO
-         END SELECT
-      END DO
-      !
-      !                           ! Migrations
-      imigr = jpreci * jpj
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         CALL mppsend( 2, zt2we(1,1,1), num_fields*imigr, noea, ml_req1 )
-         CALL mpprecv( 1, zt2ew(1,1,num_fields+1), num_fields*imigr, noea )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      CASE ( 0 )
-         CALL mppsend( 1, zt2ew(1,1,1), num_fields*imigr, nowe, ml_req1 )
-         CALL mppsend( 2, zt2we(1,1,1), num_fields*imigr, noea, ml_req2 )
-         CALL mpprecv( 1, zt2ew(1,1,num_fields+1), num_fields*imigr, noea )
-         CALL mpprecv( 2, zt2we(1,1,num_fields+1), num_fields*imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
-      CASE ( 1 )
-         CALL mppsend( 1, zt2ew(1,1,1), num_fields*imigr, nowe, ml_req1 )
-         CALL mpprecv( 2, zt2we(1,1,num_fields+1), num_fields*imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      iihom = nlci - jpreci
-      !
-
-      DO ii = 1 , num_fields
-         SELECT CASE ( nbondi )
-         CASE ( -1 )
-            DO jl = 1, jpreci
-               pt2d_array(ii)%pt2d( iihom+jl , : ) = zt2ew(:,jl,num_fields+ii)
-            END DO
-         CASE ( 0 )
-            DO jl = 1, jpreci
-               pt2d_array(ii)%pt2d( jl , : ) = zt2we(:,jl,num_fields+ii)
-               pt2d_array(ii)%pt2d( iihom+jl , : ) = zt2ew(:,jl,num_fields+ii)
-            END DO
-         CASE ( 1 )
-            DO jl = 1, jpreci
-               pt2d_array(ii)%pt2d( jl , : )= zt2we(:,jl,num_fields+ii)
-            END DO
-         END SELECT
-      END DO
-      
-      ! 3. North and south directions
-      ! -----------------------------
-      ! always closed : we play only with the neigbours
-      !
-      !First Array
-      DO ii = 1 , num_fields
-         IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
-            ijhom = nlcj-nrecj
-            DO jl = 1, jprecj
-               zt2sn(:,jl , ii) = pt2d_array(ii)%pt2d( : , ijhom +jl )
-               zt2ns(:,jl , ii) = pt2d_array(ii)%pt2d( : , jprecj+jl )
-            END DO
-         ENDIF
-      END DO
-      !
-      !                           ! Migrations
-      imigr = jprecj * jpi
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         CALL mppsend( 4, zt2sn(1,1,1), num_fields*imigr, nono, ml_req1 )
-         CALL mpprecv( 3, zt2ns(1,1,num_fields+1), num_fields*imigr, nono )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      CASE ( 0 )
-         CALL mppsend( 3, zt2ns(1,1,1), num_fields*imigr, noso, ml_req1 )
-         CALL mppsend( 4, zt2sn(1,1,1), num_fields*imigr, nono, ml_req2 )
-         CALL mpprecv( 3, zt2ns(1,1,num_fields+1), num_fields*imigr, nono )
-         CALL mpprecv( 4, zt2sn(1,1,num_fields+1), num_fields*imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
-      CASE ( 1 )
-         CALL mppsend( 3, zt2ns(1,1,1), num_fields*imigr, noso, ml_req1 )
-         CALL mpprecv( 4, zt2sn(1,1,num_fields+1), num_fields*imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      ijhom = nlcj - jprecj
-      !
-
-      DO ii = 1 , num_fields
-         !First Array
-         SELECT CASE ( nbondj )
-         CASE ( -1 )
-            DO jl = 1, jprecj
-               pt2d_array(ii)%pt2d( : , ijhom+jl ) = zt2ns( : , jl , num_fields+ii )
-            END DO
-         CASE ( 0 )
-            DO jl = 1, jprecj
-               pt2d_array(ii)%pt2d( : , jl ) = zt2sn( : , jl , num_fields + ii)
-               pt2d_array(ii)%pt2d( : , ijhom + jl ) = zt2ns( : , jl , num_fields + ii )
-            END DO
-         CASE ( 1 )
-            DO jl = 1, jprecj
-               pt2d_array(ii)%pt2d( : , jl ) = zt2sn( : , jl , num_fields + ii )
-            END DO
-         END SELECT
-      END DO
-      
-      ! 4. north fold treatment
-      ! -----------------------
-      !
-         !First Array
-      IF( npolj /= 0 .AND. .NOT. PRESENT(cd_mpp) ) THEN
-         !
-         SELECT CASE ( jpni )
-         CASE ( 1 )     ;   
-             DO ii = 1 , num_fields  
-                       CALL lbc_nfd      ( pt2d_array(ii)%pt2d( : , : ), type_array(ii) , psgn_array(ii) )   ! only 1 northern proc, no mpp
-             END DO
-         CASE DEFAULT   ;   CALL mpp_lbc_north_2d_multiple( pt2d_array, type_array, psgn_array, num_fields )   ! for all northern procs.
-         END SELECT
-         !
-      ENDIF
-        !
-      !
-      DEALLOCATE( zt2ns, zt2sn, zt2ew, zt2we )
-      !
-   END SUBROUTINE mpp_lnk_2d_multiple
-
-   
-   SUBROUTINE load_array( pt2d, cd_type, psgn, pt2d_array, type_array, psgn_array, num_fields )
-      !!---------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj), TARGET, INTENT(inout) ::   pt2d    ! Second 2D array on which the boundary condition is applied
-      CHARACTER(len=1)                    , INTENT(in   ) ::   cd_type ! define the nature of ptab array grid-points
-      REAL(wp)                            , INTENT(in   ) ::   psgn    ! =-1 the sign change across the north fold boundary
-      TYPE(arrayptr)   , DIMENSION(9) ::   pt2d_array
-      CHARACTER(len=1) , DIMENSION(9) ::   type_array    ! define the nature of ptab array grid-points
-      REAL(wp)         , DIMENSION(9) ::   psgn_array    ! =-1 the sign change across the north fold boundary
-      INTEGER                            , INTENT (inout) :: num_fields 
-      !!---------------------------------------------------------------------
-      num_fields = num_fields + 1
-      pt2d_array(num_fields)%pt2d => pt2d
-      type_array(num_fields)      =  cd_type
-      psgn_array(num_fields)      =  psgn
-   END SUBROUTINE load_array
+   !!----------------------------------------------------------------------
+   !!                   ***  routine mpp_lnk_(2,3,4)d  ***
+   !!
+   !!   * Argument : dummy argument use in mpp_lnk_... routines
+   !!                ptab   :   array or pointer of arrays on which the boundary condition is applied
+   !!                cd_nat :   nature of array grid-points
+   !!                psgn   :   sign used across the north fold boundary
+   !!                kfld   :   optional, number of pt3d arrays
+   !!                cd_mpp :   optional, fill the overlap area only
+   !!                pval   :   optional, background value (used at closed boundaries)
+   !!----------------------------------------------------------------------
+   !
+   !                       !==  2D array and array of 2D pointer  ==!
+   !
+#  define DIM_2d
+#     define ROUTINE_LNK           mpp_lnk_2d
+#     include "mpp_lnk_generic.h90"
+#     undef ROUTINE_LNK
+#     define MULTI
+#     define ROUTINE_LNK           mpp_lnk_2d_ptr
+#     include "mpp_lnk_generic.h90"
+#     undef ROUTINE_LNK
+#     undef MULTI
+#  undef DIM_2d
+   !
+   !                       !==  3D array and array of 3D pointer  ==!
+   !
+#  define DIM_3d
+#     define ROUTINE_LNK           mpp_lnk_3d
+#     include "mpp_lnk_generic.h90"
+#     undef ROUTINE_LNK
+#     define MULTI
+#     define ROUTINE_LNK           mpp_lnk_3d_ptr
+#     include "mpp_lnk_generic.h90"
+#     undef ROUTINE_LNK
+#     undef MULTI
+#  undef DIM_3d
+   !
+   !                       !==  4D array and array of 4D pointer  ==!
+   !
+#  define DIM_4d
+#     define ROUTINE_LNK           mpp_lnk_4d
+#     include "mpp_lnk_generic.h90"
+#     undef ROUTINE_LNK
+#     define MULTI
+#     define ROUTINE_LNK           mpp_lnk_4d_ptr
+#     include "mpp_lnk_generic.h90"
+#     undef ROUTINE_LNK
+#     undef MULTI
+#  undef DIM_4d
+
+   !!----------------------------------------------------------------------
+   !!                   ***  routine mpp_nfd_(2,3,4)d  ***
+   !!
+   !!   * Argument : dummy argument use in mpp_nfd_... routines
+   !!                ptab   :   array or pointer of arrays on which the boundary condition is applied
+   !!                cd_nat :   nature of array grid-points
+   !!                psgn   :   sign used across the north fold boundary
+   !!                kfld   :   optional, number of pt3d arrays
+   !!                cd_mpp :   optional, fill the overlap area only
+   !!                pval   :   optional, background value (used at closed boundaries)
+   !!----------------------------------------------------------------------
+   !
+   !                       !==  2D array and array of 2D pointer  ==!
+   !
+#  define DIM_2d
+#     define ROUTINE_NFD           mpp_nfd_2d
+#     include "mpp_nfd_generic.h90"
+#     undef ROUTINE_NFD
+#     define MULTI
+#     define ROUTINE_NFD           mpp_nfd_2d_ptr
+#     include "mpp_nfd_generic.h90"
+#     undef ROUTINE_NFD
+#     undef MULTI
+#  undef DIM_2d
+   !
+   !                       !==  3D array and array of 3D pointer  ==!
+   !
+#  define DIM_3d
+#     define ROUTINE_NFD           mpp_nfd_3d
+#     include "mpp_nfd_generic.h90"
+#     undef ROUTINE_NFD
+#     define MULTI
+#     define ROUTINE_NFD           mpp_nfd_3d_ptr
+#     include "mpp_nfd_generic.h90"
+#     undef ROUTINE_NFD
+#     undef MULTI
+#  undef DIM_3d
+   !
+   !                       !==  4D array and array of 4D pointer  ==!
+   !
+#  define DIM_4d
+#     define ROUTINE_NFD           mpp_nfd_4d
+#     include "mpp_nfd_generic.h90"
+#     undef ROUTINE_NFD
+#     define MULTI
+#     define ROUTINE_NFD           mpp_nfd_4d_ptr
+#     include "mpp_nfd_generic.h90"
+#     undef ROUTINE_NFD
+#     undef MULTI
+#  undef DIM_4d
+
+
+   !!----------------------------------------------------------------------
+   !!                   ***  routine mpp_lnk_bdy_(2,3,4)d  ***
+   !!
+   !!   * Argument : dummy argument use in mpp_lnk_... routines
+   !!                ptab   :   array or pointer of arrays on which the boundary condition is applied
+   !!                cd_nat :   nature of array grid-points
+   !!                psgn   :   sign used across the north fold boundary
+   !!                kb_bdy :   BDY boundary set
+   !!                kfld   :   optional, number of pt3d arrays
+   !!----------------------------------------------------------------------
+   !
+   !                       !==  2D array and array of 2D pointer  ==!
+   !
+#  define DIM_2d
+#     define ROUTINE_BDY           mpp_lnk_bdy_2d
+#     include "mpp_bdy_generic.h90"
+#     undef ROUTINE_BDY
+#  undef DIM_2d
+   !
+   !                       !==  3D array and array of 3D pointer  ==!
+   !
+#  define DIM_3d
+#     define ROUTINE_BDY           mpp_lnk_bdy_3d
+#     include "mpp_bdy_generic.h90"
+#     undef ROUTINE_BDY
+#  undef DIM_3d
+   !
+   !                       !==  4D array and array of 4D pointer  ==!
+   !
+#  define DIM_4d
+#     define ROUTINE_BDY           mpp_lnk_bdy_4d
+#     include "mpp_bdy_generic.h90"
+#     undef ROUTINE_BDY
+#  undef DIM_4d
+
+   !!----------------------------------------------------------------------
+   !!
+   !!   load_array  &   mpp_lnk_2d_9    à generaliser a 3D et 4D
    
    
-   SUBROUTINE mpp_lnk_2d_9( pt2dA, cd_typeA, psgnA, pt2dB, cd_typeB, psgnB, pt2dC, cd_typeC, psgnC   &
-      &                   , pt2dD, cd_typeD, psgnD, pt2dE, cd_typeE, psgnE, pt2dF, cd_typeF, psgnF   &
-      &                   , pt2dG, cd_typeG, psgnG, pt2dH, cd_typeH, psgnH, pt2dI, cd_typeI, psgnI, cd_mpp, pval)
-      !!---------------------------------------------------------------------
-      ! Second 2D array on which the boundary condition is applied
-      REAL(wp), DIMENSION(jpi,jpj), TARGET          , INTENT(inout) ::   pt2dA    
-      REAL(wp), DIMENSION(jpi,jpj), TARGET, OPTIONAL, INTENT(inout) ::   pt2dB , pt2dC , pt2dD , pt2dE
-      REAL(wp), DIMENSION(jpi,jpj), TARGET, OPTIONAL, INTENT(inout) ::   pt2dF , pt2dG , pt2dH , pt2dI 
-      ! define the nature of ptab array grid-points
-      CHARACTER(len=1)                              , INTENT(in   ) ::   cd_typeA
-      CHARACTER(len=1)                    , OPTIONAL, INTENT(in   ) ::   cd_typeB , cd_typeC , cd_typeD , cd_typeE
-      CHARACTER(len=1)                    , OPTIONAL, INTENT(in   ) ::   cd_typeF , cd_typeG , cd_typeH , cd_typeI
-      ! =-1 the sign change across the north fold boundary
-      REAL(wp)                                      , INTENT(in   ) ::   psgnA    
-      REAL(wp)                            , OPTIONAL, INTENT(in   ) ::   psgnB , psgnC , psgnD , psgnE
-      REAL(wp)                            , OPTIONAL, INTENT(in   ) ::   psgnF , psgnG , psgnH , psgnI   
-      CHARACTER(len=3)                    , OPTIONAL, INTENT(in   ) ::   cd_mpp   ! fill the overlap area only
-      REAL(wp)                            , OPTIONAL, INTENT(in   ) ::   pval     ! background value (used at closed boundaries)
-      !!
-      TYPE(arrayptr)   , DIMENSION(9) ::   pt2d_array 
-      CHARACTER(len=1) , DIMENSION(9) ::   type_array    ! define the nature of ptab array grid-points
-      !                                                         ! = T , U , V , F , W and I points
-      REAL(wp)         , DIMENSION(9) ::   psgn_array    ! =-1 the sign change across the north fold boundary
-      INTEGER :: num_fields
-      !!---------------------------------------------------------------------
-      !
-      num_fields = 0
-      !
-      ! Load the first array
-      CALL load_array( pt2dA, cd_typeA, psgnA, pt2d_array, type_array, psgn_array, num_fields )
-      !
-      ! Look if more arrays are added
-      IF( PRESENT(psgnB) )   CALL load_array(pt2dB,cd_typeB,psgnB,pt2d_array, type_array, psgn_array,num_fields)
-      IF( PRESENT(psgnC) )   CALL load_array(pt2dC,cd_typeC,psgnC,pt2d_array, type_array, psgn_array,num_fields)
-      IF( PRESENT(psgnD) )   CALL load_array(pt2dD,cd_typeD,psgnD,pt2d_array, type_array, psgn_array,num_fields)
-      IF( PRESENT(psgnE) )   CALL load_array(pt2dE,cd_typeE,psgnE,pt2d_array, type_array, psgn_array,num_fields)
-      IF( PRESENT(psgnF) )   CALL load_array(pt2dF,cd_typeF,psgnF,pt2d_array, type_array, psgn_array,num_fields)
-      IF( PRESENT(psgnG) )   CALL load_array(pt2dG,cd_typeG,psgnG,pt2d_array, type_array, psgn_array,num_fields)
-      IF( PRESENT(psgnH) )   CALL load_array(pt2dH,cd_typeH,psgnH,pt2d_array, type_array, psgn_array,num_fields)
-      IF( PRESENT(psgnI) )   CALL load_array(pt2dI,cd_typeI,psgnI,pt2d_array, type_array, psgn_array,num_fields)
-      !
-      CALL mpp_lnk_2d_multiple( pt2d_array, type_array, psgn_array, num_fields, cd_mpp,pval )
-      !
-   END SUBROUTINE mpp_lnk_2d_9
-
-
-   SUBROUTINE mpp_lnk_2d( pt2d, cd_type, psgn, cd_mpp, pval )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mpp_lnk_2d  ***
-      !!
-      !! ** Purpose :   Message passing manadgement for 2d array
-      !!
-      !! ** Method  :   Use mppsend and mpprecv function for passing mask
-      !!      between processors following neighboring subdomains.
-      !!            domain parameters
-      !!                    nlci   : first dimension of the local subdomain
-      !!                    nlcj   : second dimension of the local subdomain
-      !!                    nbondi : mark for "east-west local boundary"
-      !!                    nbondj : mark for "north-south local boundary"
-      !!                    noea   : number for local neighboring processors
-      !!                    nowe   : number for local neighboring processors
-      !!                    noso   : number for local neighboring processors
-      !!                    nono   : number for local neighboring processors
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   pt2d     ! 2D array on which the boundary condition is applied
-      CHARACTER(len=1)            , INTENT(in   ) ::   cd_type  ! define the nature of ptab array grid-points
-      !                                                         ! = T , U , V , F , W and I points
-      REAL(wp)                    , INTENT(in   ) ::   psgn     ! =-1 the sign change across the north fold boundary
-      !                                                         ! =  1. , the sign is kept
-      CHARACTER(len=3), OPTIONAL  , INTENT(in   ) ::   cd_mpp   ! fill the overlap area only
-      REAL(wp)        , OPTIONAL  , INTENT(in   ) ::   pval     ! background value (used at closed boundaries)
-      !!
-      INTEGER  ::   ji, jj, jl   ! dummy loop indices
-      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
-      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
-      REAL(wp) ::   zland
-      INTEGER, DIMENSION(MPI_STATUS_SIZE)     ::   ml_stat       ! for key_mpi_isend
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  zt2ns, zt2sn   ! 2d for north-south & south-north
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  zt2ew, zt2we   ! 2d for east-west & west-east
-      !!----------------------------------------------------------------------
-      !
-      ALLOCATE( zt2ns(jpi,jprecj,2), zt2sn(jpi,jprecj,2),  &
-         &      zt2ew(jpj,jpreci,2), zt2we(jpj,jpreci,2)   )
-      !
-      IF( PRESENT( pval ) ) THEN   ;   zland = pval      ! set land value
-      ELSE                         ;   zland = 0._wp     ! zero by default
-      ENDIF
-
-      ! 1. standard boundary treatment
-      ! ------------------------------
-      !
-      IF( PRESENT( cd_mpp ) ) THEN      ! only fill added line/raw with existing values
-         !
-         ! WARNING pt2d is defined only between nld and nle
-         DO jj = nlcj+1, jpj                 ! added line(s)   (inner only)
-            pt2d(nldi  :nlei  , jj          ) = pt2d(nldi:nlei,     nlej)
-            pt2d(1     :nldi-1, jj          ) = pt2d(nldi     ,     nlej)
-            pt2d(nlei+1:nlci  , jj          ) = pt2d(     nlei,     nlej)
-         END DO
-         DO ji = nlci+1, jpi                 ! added column(s) (full)
-            pt2d(ji           ,nldj  :nlej  ) = pt2d(     nlei,nldj:nlej)
-            pt2d(ji           ,1     :nldj-1) = pt2d(     nlei,nldj     )
-            pt2d(ji           ,nlej+1:jpj   ) = pt2d(     nlei,     nlej)
-         END DO
-         !
-      ELSE                              ! standard close or cyclic treatment
-         !
-         !                                   ! East-West boundaries
-         IF( nbondi == 2 .AND.   &                ! Cyclic east-west
-            &    (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
-            pt2d( 1 ,:) = pt2d(jpim1,:)                                    ! west
-            pt2d(jpi,:) = pt2d(  2  ,:)                                    ! east
-         ELSE                                     ! closed
-            IF( .NOT. cd_type == 'F' )   pt2d(     1       :jpreci,:) = zland    ! south except F-point
-                                         pt2d(nlci-jpreci+1:jpi   ,:) = zland    ! north
-         ENDIF
-         !                                   ! North-South boundaries (always closed)
-            IF( .NOT. cd_type == 'F' )   pt2d(:,     1       :jprecj) = zland    !south except F-point
-                                         pt2d(:,nlcj-jprecj+1:jpj   ) = zland    ! north
-         !
-      ENDIF
-
-      ! 2. East and west directions exchange
-      ! ------------------------------------
-      ! we play with the neigbours AND the row number because of the periodicity
-      !
-      SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
-      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
-         iihom = nlci-nreci
-         DO jl = 1, jpreci
-            zt2ew(:,jl,1) = pt2d(jpreci+jl,:)
-            zt2we(:,jl,1) = pt2d(iihom +jl,:)
-         END DO
-      END SELECT
-      !
-      !                           ! Migrations
-      imigr = jpreci * jpj
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         CALL mppsend( 2, zt2we(1,1,1), imigr, noea, ml_req1 )
-         CALL mpprecv( 1, zt2ew(1,1,2), imigr, noea )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      CASE ( 0 )
-         CALL mppsend( 1, zt2ew(1,1,1), imigr, nowe, ml_req1 )
-         CALL mppsend( 2, zt2we(1,1,1), imigr, noea, ml_req2 )
-         CALL mpprecv( 1, zt2ew(1,1,2), imigr, noea )
-         CALL mpprecv( 2, zt2we(1,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
-      CASE ( 1 )
-         CALL mppsend( 1, zt2ew(1,1,1), imigr, nowe, ml_req1 )
-         CALL mpprecv( 2, zt2we(1,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      iihom = nlci - jpreci
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         DO jl = 1, jpreci
-            pt2d(iihom+jl,:) = zt2ew(:,jl,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jpreci
-            pt2d(jl      ,:) = zt2we(:,jl,2)
-            pt2d(iihom+jl,:) = zt2ew(:,jl,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jpreci
-            pt2d(jl      ,:) = zt2we(:,jl,2)
-         END DO
-      END SELECT
-
-
-      ! 3. North and south directions
-      ! -----------------------------
-      ! always closed : we play only with the neigbours
-      !
-      IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
-         ijhom = nlcj-nrecj
-         DO jl = 1, jprecj
-            zt2sn(:,jl,1) = pt2d(:,ijhom +jl)
-            zt2ns(:,jl,1) = pt2d(:,jprecj+jl)
-         END DO
-      ENDIF
-      !
-      !                           ! Migrations
-      imigr = jprecj * jpi
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         CALL mppsend( 4, zt2sn(1,1,1), imigr, nono, ml_req1 )
-         CALL mpprecv( 3, zt2ns(1,1,2), imigr, nono )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      CASE ( 0 )
-         CALL mppsend( 3, zt2ns(1,1,1), imigr, noso, ml_req1 )
-         CALL mppsend( 4, zt2sn(1,1,1), imigr, nono, ml_req2 )
-         CALL mpprecv( 3, zt2ns(1,1,2), imigr, nono )
-         CALL mpprecv( 4, zt2sn(1,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
-      CASE ( 1 )
-         CALL mppsend( 3, zt2ns(1,1,1), imigr, noso, ml_req1 )
-         CALL mpprecv( 4, zt2sn(1,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      ijhom = nlcj - jprecj
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         DO jl = 1, jprecj
-            pt2d(:,ijhom+jl) = zt2ns(:,jl,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jprecj
-            pt2d(:,jl      ) = zt2sn(:,jl,2)
-            pt2d(:,ijhom+jl) = zt2ns(:,jl,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jprecj
-            pt2d(:,jl      ) = zt2sn(:,jl,2)
-         END DO
-      END SELECT
-
-
-      ! 4. north fold treatment
-      ! -----------------------
-      !
-      IF( npolj /= 0 .AND. .NOT. PRESENT(cd_mpp) ) THEN
-         !
-         SELECT CASE ( jpni )
-         CASE ( 1 )     ;   CALL lbc_nfd      ( pt2d, cd_type, psgn )   ! only 1 northern proc, no mpp
-         CASE DEFAULT   ;   CALL mpp_lbc_north( pt2d, cd_type, psgn )   ! for all northern procs.
-         END SELECT
-         !
-      ENDIF
-      !
-      DEALLOCATE( zt2ns, zt2sn, zt2ew, zt2we )
-      !
-   END SUBROUTINE mpp_lnk_2d
-
-
-   SUBROUTINE mpp_lnk_3d_gather( ptab1, cd_type1, ptab2, cd_type2, psgn )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mpp_lnk_3d_gather  ***
-      !!
-      !! ** Purpose :   Message passing manadgement for two 3D arrays
-      !!
-      !! ** Method  :   Use mppsend and mpprecv function for passing mask
-      !!      between processors following neighboring subdomains.
-      !!            domain parameters
-      !!                    nlci   : first dimension of the local subdomain
-      !!                    nlcj   : second dimension of the local subdomain
-      !!                    nbondi : mark for "east-west local boundary"
-      !!                    nbondj : mark for "north-south local boundary"
-      !!                    noea   : number for local neighboring processors
-      !!                    nowe   : number for local neighboring processors
-      !!                    noso   : number for local neighboring processors
-      !!                    nono   : number for local neighboring processors
-      !!
-      !! ** Action  :   ptab1 and ptab2  with update value at its periphery
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptab1     ! first and second 3D array on which
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptab2     ! the boundary condition is applied
-      CHARACTER(len=1)                , INTENT(in   ) ::   cd_type1  ! nature of ptab1 and ptab2 arrays
-      CHARACTER(len=1)                , INTENT(in   ) ::   cd_type2  ! i.e. grid-points = T , U , V , F or W points
-      REAL(wp)                        , INTENT(in   ) ::   psgn      ! =-1 the sign change across the north fold boundary
-      !!                                                             ! =  1. , the sign is kept
-      INTEGER  ::   jl   ! dummy loop indices
-      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
-      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
-      INTEGER , DIMENSION(MPI_STATUS_SIZE)        ::   ml_stat   ! for key_mpi_isend
-      REAL(wp), DIMENSION(:,:,:,:,:), ALLOCATABLE ::   zt4ns, zt4sn   ! 2 x 3d for north-south & south-north
-      REAL(wp), DIMENSION(:,:,:,:,:), ALLOCATABLE ::   zt4ew, zt4we   ! 2 x 3d for east-west & west-east
-      !!----------------------------------------------------------------------
-      !
-      ALLOCATE( zt4ns(jpi,jprecj,jpk,2,2), zt4sn(jpi,jprecj,jpk,2,2) ,    &
-         &      zt4ew(jpj,jpreci,jpk,2,2), zt4we(jpj,jpreci,jpk,2,2) )
-      !
-      ! 1. standard boundary treatment
-      ! ------------------------------
-      !                                      ! East-West boundaries
-      !                                           !* Cyclic east-west
-      IF( nbondi == 2 .AND. (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
-         ptab1( 1 ,:,:) = ptab1(jpim1,:,:)
-         ptab1(jpi,:,:) = ptab1(  2  ,:,:)
-         ptab2( 1 ,:,:) = ptab2(jpim1,:,:)
-         ptab2(jpi,:,:) = ptab2(  2  ,:,:)
-      ELSE                                        !* closed
-         IF( .NOT. cd_type1 == 'F' )   ptab1(     1       :jpreci,:,:) = 0.e0    ! south except at F-point
-         IF( .NOT. cd_type2 == 'F' )   ptab2(     1       :jpreci,:,:) = 0.e0
-                                       ptab1(nlci-jpreci+1:jpi   ,:,:) = 0.e0    ! north
-                                       ptab2(nlci-jpreci+1:jpi   ,:,:) = 0.e0
-      ENDIF
-
-
-      !                                      ! North-South boundaries
-      IF( .NOT. cd_type1 == 'F' )   ptab1(:,     1       :jprecj,:) = 0.e0    ! south except at F-point
-      IF( .NOT. cd_type2 == 'F' )   ptab2(:,     1       :jprecj,:) = 0.e0
-                                    ptab1(:,nlcj-jprecj+1:jpj   ,:) = 0.e0    ! north
-                                    ptab2(:,nlcj-jprecj+1:jpj   ,:) = 0.e0
-
-
-      ! 2. East and west directions exchange
-      ! ------------------------------------
-      ! we play with the neigbours AND the row number because of the periodicity
-      !
-      SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
-      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
-         iihom = nlci-nreci
-         DO jl = 1, jpreci
-            zt4ew(:,jl,:,1,1) = ptab1(jpreci+jl,:,:)
-            zt4we(:,jl,:,1,1) = ptab1(iihom +jl,:,:)
-            zt4ew(:,jl,:,2,1) = ptab2(jpreci+jl,:,:)
-            zt4we(:,jl,:,2,1) = ptab2(iihom +jl,:,:)
-         END DO
-      END SELECT
-      !
-      !                           ! Migrations
-      imigr = jpreci * jpj * jpk *2
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         CALL mppsend( 2, zt4we(1,1,1,1,1), imigr, noea, ml_req1 )
-         CALL mpprecv( 1, zt4ew(1,1,1,1,2), imigr, noea )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      CASE ( 0 )
-         CALL mppsend( 1, zt4ew(1,1,1,1,1), imigr, nowe, ml_req1 )
-         CALL mppsend( 2, zt4we(1,1,1,1,1), imigr, noea, ml_req2 )
-         CALL mpprecv( 1, zt4ew(1,1,1,1,2), imigr, noea )
-         CALL mpprecv( 2, zt4we(1,1,1,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
-      CASE ( 1 )
-         CALL mppsend( 1, zt4ew(1,1,1,1,1), imigr, nowe, ml_req1 )
-         CALL mpprecv( 2, zt4we(1,1,1,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      iihom = nlci - jpreci
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         DO jl = 1, jpreci
-            ptab1(iihom+jl,:,:) = zt4ew(:,jl,:,1,2)
-            ptab2(iihom+jl,:,:) = zt4ew(:,jl,:,2,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jpreci
-            ptab1(jl      ,:,:) = zt4we(:,jl,:,1,2)
-            ptab1(iihom+jl,:,:) = zt4ew(:,jl,:,1,2)
-            ptab2(jl      ,:,:) = zt4we(:,jl,:,2,2)
-            ptab2(iihom+jl,:,:) = zt4ew(:,jl,:,2,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jpreci
-            ptab1(jl      ,:,:) = zt4we(:,jl,:,1,2)
-            ptab2(jl      ,:,:) = zt4we(:,jl,:,2,2)
-         END DO
-      END SELECT
-
-
-      ! 3. North and south directions
-      ! -----------------------------
-      ! always closed : we play only with the neigbours
-      !
-      IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
-         ijhom = nlcj - nrecj
-         DO jl = 1, jprecj
-            zt4sn(:,jl,:,1,1) = ptab1(:,ijhom +jl,:)
-            zt4ns(:,jl,:,1,1) = ptab1(:,jprecj+jl,:)
-            zt4sn(:,jl,:,2,1) = ptab2(:,ijhom +jl,:)
-            zt4ns(:,jl,:,2,1) = ptab2(:,jprecj+jl,:)
-         END DO
-      ENDIF
-      !
-      !                           ! Migrations
-      imigr = jprecj * jpi * jpk * 2
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         CALL mppsend( 4, zt4sn(1,1,1,1,1), imigr, nono, ml_req1 )
-         CALL mpprecv( 3, zt4ns(1,1,1,1,2), imigr, nono )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      CASE ( 0 )
-         CALL mppsend( 3, zt4ns(1,1,1,1,1), imigr, noso, ml_req1 )
-         CALL mppsend( 4, zt4sn(1,1,1,1,1), imigr, nono, ml_req2 )
-         CALL mpprecv( 3, zt4ns(1,1,1,1,2), imigr, nono )
-         CALL mpprecv( 4, zt4sn(1,1,1,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
-      CASE ( 1 )
-         CALL mppsend( 3, zt4ns(1,1,1,1,1), imigr, noso, ml_req1 )
-         CALL mpprecv( 4, zt4sn(1,1,1,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      ijhom = nlcj - jprecj
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         DO jl = 1, jprecj
-            ptab1(:,ijhom+jl,:) = zt4ns(:,jl,:,1,2)
-            ptab2(:,ijhom+jl,:) = zt4ns(:,jl,:,2,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jprecj
-            ptab1(:,jl      ,:) = zt4sn(:,jl,:,1,2)
-            ptab1(:,ijhom+jl,:) = zt4ns(:,jl,:,1,2)
-            ptab2(:,jl      ,:) = zt4sn(:,jl,:,2,2)
-            ptab2(:,ijhom+jl,:) = zt4ns(:,jl,:,2,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jprecj
-            ptab1(:,jl,:) = zt4sn(:,jl,:,1,2)
-            ptab2(:,jl,:) = zt4sn(:,jl,:,2,2)
-         END DO
-      END SELECT
-
-
-      ! 4. north fold treatment
-      ! -----------------------
-      IF( npolj /= 0 ) THEN
-         !
-         SELECT CASE ( jpni )
-         CASE ( 1 )
-            CALL lbc_nfd      ( ptab1, cd_type1, psgn )   ! only for northern procs.
-            CALL lbc_nfd      ( ptab2, cd_type2, psgn )
-         CASE DEFAULT
-            CALL mpp_lbc_north( ptab1, cd_type1, psgn )   ! for all northern procs.
-            CALL mpp_lbc_north (ptab2, cd_type2, psgn)
-         END SELECT
-         !
-      ENDIF
-      !
-      DEALLOCATE( zt4ns, zt4sn, zt4ew, zt4we )
-      !
-   END SUBROUTINE mpp_lnk_3d_gather
-
-
-   SUBROUTINE mpp_lnk_2d_e( pt2d, cd_type, psgn, jpri, jprj )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mpp_lnk_2d_e  ***
-      !!
-      !! ** Purpose :   Message passing manadgement for 2d array (with halo)
-      !!
-      !! ** Method  :   Use mppsend and mpprecv function for passing mask
-      !!      between processors following neighboring subdomains.
-      !!            domain parameters
-      !!                    nlci   : first dimension of the local subdomain
-      !!                    nlcj   : second dimension of the local subdomain
-      !!                    jpri   : number of rows for extra outer halo
-      !!                    jprj   : number of columns for extra outer halo
-      !!                    nbondi : mark for "east-west local boundary"
-      !!                    nbondj : mark for "north-south local boundary"
-      !!                    noea   : number for local neighboring processors
-      !!                    nowe   : number for local neighboring processors
-      !!                    noso   : number for local neighboring processors
-      !!                    nono   : number for local neighboring processors
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER                                             , INTENT(in   ) ::   jpri
-      INTEGER                                             , INTENT(in   ) ::   jprj
-      REAL(wp), DIMENSION(1-jpri:jpi+jpri,1-jprj:jpj+jprj), INTENT(inout) ::   pt2d     ! 2D array with extra halo
-      CHARACTER(len=1)                                    , INTENT(in   ) ::   cd_type  ! nature of ptab array grid-points
-      !                                                                                 ! = T , U , V , F , W and I points
-      REAL(wp)                                            , INTENT(in   ) ::   psgn     ! =-1 the sign change across the
-      !!                                                                                ! north boundary, =  1. otherwise
-      INTEGER  ::   jl   ! dummy loop indices
-      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
-      INTEGER  ::   ipreci, iprecj             ! temporary integers
-      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
-      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
-      !!
-      REAL(wp), DIMENSION(1-jpri:jpi+jpri,jprecj+jprj,2) :: r2dns
-      REAL(wp), DIMENSION(1-jpri:jpi+jpri,jprecj+jprj,2) :: r2dsn
-      REAL(wp), DIMENSION(1-jprj:jpj+jprj,jpreci+jpri,2) :: r2dwe
-      REAL(wp), DIMENSION(1-jprj:jpj+jprj,jpreci+jpri,2) :: r2dew
-      !!----------------------------------------------------------------------
-
-      ipreci = jpreci + jpri      ! take into account outer extra 2D overlap area
-      iprecj = jprecj + jprj
-
-
-      ! 1. standard boundary treatment
-      ! ------------------------------
-      ! Order matters Here !!!!
-      !
-      !                                      !* North-South boundaries (always colsed)
-      IF( .NOT. cd_type == 'F' )   pt2d(:,  1-jprj   :  jprecj  ) = 0.e0    ! south except at F-point
-                                   pt2d(:,nlcj-jprecj+1:jpj+jprj) = 0.e0    ! north
-
-      !                                      ! East-West boundaries
-      !                                           !* Cyclic east-west
-      IF( nbondi == 2 .AND. (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
-         pt2d(1-jpri:     1    ,:) = pt2d(jpim1-jpri:  jpim1 ,:)       ! east
-         pt2d(   jpi  :jpi+jpri,:) = pt2d(     2      :2+jpri,:)       ! west
-         !
-      ELSE                                        !* closed
-         IF( .NOT. cd_type == 'F' )   pt2d(  1-jpri   :jpreci    ,:) = 0.e0    ! south except at F-point
-                                      pt2d(nlci-jpreci+1:jpi+jpri,:) = 0.e0    ! north
-      ENDIF
-      !
-
-      ! north fold treatment
-      ! -----------------------
-      IF( npolj /= 0 ) THEN
-         !
-         SELECT CASE ( jpni )
-         CASE ( 1 )     ;   CALL lbc_nfd        ( pt2d(1:jpi,1:jpj+jprj), cd_type, psgn, pr2dj=jprj )
-         CASE DEFAULT   ;   CALL mpp_lbc_north_e( pt2d                    , cd_type, psgn               )
-         END SELECT
-         !
-      ENDIF
-
-      ! 2. East and west directions exchange
-      ! ------------------------------------
-      ! we play with the neigbours AND the row number because of the periodicity
-      !
-      SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
-      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
-         iihom = nlci-nreci-jpri
-         DO jl = 1, ipreci
-            r2dew(:,jl,1) = pt2d(jpreci+jl,:)
-            r2dwe(:,jl,1) = pt2d(iihom +jl,:)
-         END DO
-      END SELECT
-      !
-      !                           ! Migrations
-      imigr = ipreci * ( jpj + 2*jprj)
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         CALL mppsend( 2, r2dwe(1-jprj,1,1), imigr, noea, ml_req1 )
-         CALL mpprecv( 1, r2dew(1-jprj,1,2), imigr, noea )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      CASE ( 0 )
-         CALL mppsend( 1, r2dew(1-jprj,1,1), imigr, nowe, ml_req1 )
-         CALL mppsend( 2, r2dwe(1-jprj,1,1), imigr, noea, ml_req2 )
-         CALL mpprecv( 1, r2dew(1-jprj,1,2), imigr, noea )
-         CALL mpprecv( 2, r2dwe(1-jprj,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
-      CASE ( 1 )
-         CALL mppsend( 1, r2dew(1-jprj,1,1), imigr, nowe, ml_req1 )
-         CALL mpprecv( 2, r2dwe(1-jprj,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      iihom = nlci - jpreci
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         DO jl = 1, ipreci
-            pt2d(iihom+jl,:) = r2dew(:,jl,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, ipreci
-            pt2d(jl-jpri,:) = r2dwe(:,jl,2)
-            pt2d( iihom+jl,:) = r2dew(:,jl,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, ipreci
-            pt2d(jl-jpri,:) = r2dwe(:,jl,2)
-         END DO
-      END SELECT
-
-
-      ! 3. North and south directions
-      ! -----------------------------
-      ! always closed : we play only with the neigbours
-      !
-      IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
-         ijhom = nlcj-nrecj-jprj
-         DO jl = 1, iprecj
-            r2dsn(:,jl,1) = pt2d(:,ijhom +jl)
-            r2dns(:,jl,1) = pt2d(:,jprecj+jl)
-         END DO
-      ENDIF
-      !
-      !                           ! Migrations
-      imigr = iprecj * ( jpi + 2*jpri )
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         CALL mppsend( 4, r2dsn(1-jpri,1,1), imigr, nono, ml_req1 )
-         CALL mpprecv( 3, r2dns(1-jpri,1,2), imigr, nono )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      CASE ( 0 )
-         CALL mppsend( 3, r2dns(1-jpri,1,1), imigr, noso, ml_req1 )
-         CALL mppsend( 4, r2dsn(1-jpri,1,1), imigr, nono, ml_req2 )
-         CALL mpprecv( 3, r2dns(1-jpri,1,2), imigr, nono )
-         CALL mpprecv( 4, r2dsn(1-jpri,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
-      CASE ( 1 )
-         CALL mppsend( 3, r2dns(1-jpri,1,1), imigr, noso, ml_req1 )
-         CALL mpprecv( 4, r2dsn(1-jpri,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      ijhom = nlcj - jprecj
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         DO jl = 1, iprecj
-            pt2d(:,ijhom+jl) = r2dns(:,jl,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, iprecj
-            pt2d(:,jl-jprj) = r2dsn(:,jl,2)
-            pt2d(:,ijhom+jl ) = r2dns(:,jl,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, iprecj
-            pt2d(:,jl-jprj) = r2dsn(:,jl,2)
-         END DO
-      END SELECT
-      !
-   END SUBROUTINE mpp_lnk_2d_e
-
-   SUBROUTINE mpp_lnk_sum_3d( ptab, cd_type, psgn, cd_mpp, pval )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mpp_lnk_sum_3d  ***
-      !!
-      !! ** Purpose :   Message passing manadgement (sum the overlap region)
-      !!
-      !! ** Method  :   Use mppsend and mpprecv function for passing mask
-      !!      between processors following neighboring subdomains.
-      !!            domain parameters
-      !!                    nlci   : first dimension of the local subdomain
-      !!                    nlcj   : second dimension of the local subdomain
-      !!                    nbondi : mark for "east-west local boundary"
-      !!                    nbondj : mark for "north-south local boundary"
-      !!                    noea   : number for local neighboring processors
-      !!                    nowe   : number for local neighboring processors
-      !!                    noso   : number for local neighboring processors
-      !!                    nono   : number for local neighboring processors
-      !!
-      !! ** Action  :   ptab with update value at its periphery
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptab     ! 3D array on which the boundary condition is applied
-      CHARACTER(len=1)                , INTENT(in   ) ::   cd_type  ! define the nature of ptab array grid-points
-      !                                                             ! = T , U , V , F , W points
-      REAL(wp)                        , INTENT(in   ) ::   psgn     ! =-1 the sign change across the north fold boundary
-      !                                                             ! =  1. , the sign is kept
-      CHARACTER(len=3), OPTIONAL      , INTENT(in   ) ::   cd_mpp   ! fill the overlap area only
-      REAL(wp)        , OPTIONAL      , INTENT(in   ) ::   pval     ! background value (used at closed boundaries)
-      !!
-      INTEGER  ::   ji, jj, jk, jl             ! dummy loop indices
-      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
-      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
-      REAL(wp) ::   zland
-      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
-      !
-      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::   zt3ns, zt3sn   ! 3d for north-south & south-north
-      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::   zt3ew, zt3we   ! 3d for east-west & west-east
-
-      !!----------------------------------------------------------------------
-      
-      ALLOCATE( zt3ns(jpi,jprecj,jpk,2), zt3sn(jpi,jprecj,jpk,2),   &
-         &      zt3ew(jpj,jpreci,jpk,2), zt3we(jpj,jpreci,jpk,2)  )
-
-      !
-      IF( PRESENT( pval ) ) THEN   ;   zland = pval      ! set land value
-      ELSE                         ;   zland = 0.e0      ! zero by default
-      ENDIF
-
-      ! 1. standard boundary treatment
-      ! ------------------------------
-      ! 2. East and west directions exchange
-      ! ------------------------------------
-      ! we play with the neigbours AND the row number because of the periodicity
-      !
-      SELECT CASE ( nbondi )      ! Read lateral conditions
-      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
-      iihom = nlci-jpreci
-         DO jl = 1, jpreci
-            zt3ew(:,jl,:,1) = ptab(jl      ,:,:) ; ptab(jl      ,:,:) = 0.0_wp
-            zt3we(:,jl,:,1) = ptab(iihom+jl,:,:) ; ptab(iihom+jl,:,:) = 0.0_wp 
-         END DO
-      END SELECT
-      !
-      !                           ! Migrations
-      imigr = jpreci * jpj * jpk
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         CALL mppsend( 2, zt3we(1,1,1,1), imigr, noea, ml_req1 )
-         CALL mpprecv( 1, zt3ew(1,1,1,2), imigr, noea )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      CASE ( 0 )
-         CALL mppsend( 1, zt3ew(1,1,1,1), imigr, nowe, ml_req1 )
-         CALL mppsend( 2, zt3we(1,1,1,1), imigr, noea, ml_req2 )
-         CALL mpprecv( 1, zt3ew(1,1,1,2), imigr, noea )
-         CALL mpprecv( 2, zt3we(1,1,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
-      CASE ( 1 )
-         CALL mppsend( 1, zt3ew(1,1,1,1), imigr, nowe, ml_req1 )
-         CALL mpprecv( 2, zt3we(1,1,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      END SELECT
-      !
-      !                           ! Write lateral conditions
-      iihom = nlci-nreci
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         DO jl = 1, jpreci
-            ptab(iihom+jl,:,:) = ptab(iihom+jl,:,:) + zt3ew(:,jl,:,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jpreci
-            ptab(jpreci+jl,:,:) = ptab(jpreci+jl,:,:) + zt3we(:,jl,:,2)
-            ptab(iihom +jl,:,:) = ptab(iihom +jl,:,:) + zt3ew(:,jl,:,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jpreci
-            ptab(jpreci+jl,:,:) = ptab(jpreci+jl,:,:) + zt3we(:,jl,:,2)
-         END DO
-      END SELECT
-
-
-      ! 3. North and south directions
-      ! -----------------------------
-      ! always closed : we play only with the neigbours
-      !
-      IF( nbondj /= 2 ) THEN      ! Read lateral conditions
-         ijhom = nlcj-jprecj
-         DO jl = 1, jprecj
-            zt3sn(:,jl,:,1) = ptab(:,ijhom+jl,:) ; ptab(:,ijhom+jl,:) = 0.0_wp
-            zt3ns(:,jl,:,1) = ptab(:,jl      ,:) ; ptab(:,jl      ,:) = 0.0_wp
-         END DO
-      ENDIF
-      !
-      !                           ! Migrations
-      imigr = jprecj * jpi * jpk
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         CALL mppsend( 4, zt3sn(1,1,1,1), imigr, nono, ml_req1 )
-         CALL mpprecv( 3, zt3ns(1,1,1,2), imigr, nono )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      CASE ( 0 )
-         CALL mppsend( 3, zt3ns(1,1,1,1), imigr, noso, ml_req1 )
-         CALL mppsend( 4, zt3sn(1,1,1,1), imigr, nono, ml_req2 )
-         CALL mpprecv( 3, zt3ns(1,1,1,2), imigr, nono )
-         CALL mpprecv( 4, zt3sn(1,1,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
-      CASE ( 1 )
-         CALL mppsend( 3, zt3ns(1,1,1,1), imigr, noso, ml_req1 )
-         CALL mpprecv( 4, zt3sn(1,1,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      END SELECT
-      !
-      !                           ! Write lateral conditions
-      ijhom = nlcj-nrecj
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         DO jl = 1, jprecj
-            ptab(:,ijhom+jl,:) = ptab(:,ijhom+jl,:) + zt3ns(:,jl,:,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jprecj
-            ptab(:,jprecj+jl,:) = ptab(:,jprecj+jl,:) + zt3sn(:,jl,:,2)
-            ptab(:,ijhom +jl,:) = ptab(:,ijhom +jl,:) + zt3ns(:,jl,:,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jprecj
-            ptab(:,jprecj+jl,:) = ptab(:,jprecj+jl,:) + zt3sn(:,jl   ,:,2)
-         END DO
-      END SELECT
-
-
-      ! 4. north fold treatment
-      ! -----------------------
-      !
-      IF( npolj /= 0 .AND. .NOT. PRESENT(cd_mpp) ) THEN
-         !
-         SELECT CASE ( jpni )
-         CASE ( 1 )     ;   CALL lbc_nfd      ( ptab, cd_type, psgn )   ! only 1 northern proc, no mpp
-         CASE DEFAULT   ;   CALL mpp_lbc_north( ptab, cd_type, psgn )   ! for all northern procs.
-         END SELECT
-         !
-      ENDIF
-      !
-      DEALLOCATE( zt3ns, zt3sn, zt3ew, zt3we )
-      !
-   END SUBROUTINE mpp_lnk_sum_3d
-
-   SUBROUTINE mpp_lnk_sum_2d( pt2d, cd_type, psgn, cd_mpp, pval )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mpp_lnk_sum_2d  ***
-      !!
-      !! ** Purpose :   Message passing manadgement for 2d array (sum the overlap region)
-      !!
-      !! ** Method  :   Use mppsend and mpprecv function for passing mask
-      !!      between processors following neighboring subdomains.
-      !!            domain parameters
-      !!                    nlci   : first dimension of the local subdomain
-      !!                    nlcj   : second dimension of the local subdomain
-      !!                    nbondi : mark for "east-west local boundary"
-      !!                    nbondj : mark for "north-south local boundary"
-      !!                    noea   : number for local neighboring processors
-      !!                    nowe   : number for local neighboring processors
-      !!                    noso   : number for local neighboring processors
-      !!                    nono   : number for local neighboring processors
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   pt2d     ! 2D array on which the boundary condition is applied
-      CHARACTER(len=1)            , INTENT(in   ) ::   cd_type  ! define the nature of ptab array grid-points
-      !                                                         ! = T , U , V , F , W and I points
-      REAL(wp)                    , INTENT(in   ) ::   psgn     ! =-1 the sign change across the north fold boundary
-      !                                                         ! =  1. , the sign is kept
-      CHARACTER(len=3), OPTIONAL  , INTENT(in   ) ::   cd_mpp   ! fill the overlap area only
-      REAL(wp)        , OPTIONAL  , INTENT(in   ) ::   pval     ! background value (used at closed boundaries)
-      !!
-      INTEGER  ::   ji, jj, jl   ! dummy loop indices
-      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
-      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
-      REAL(wp) ::   zland
-      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
-      !
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  zt2ns, zt2sn   ! 2d for north-south & south-north
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  zt2ew, zt2we   ! 2d for east-west & west-east
-
-      !!----------------------------------------------------------------------
-
-      ALLOCATE( zt2ns(jpi,jprecj,2), zt2sn(jpi,jprecj,2),  &
-         &      zt2ew(jpj,jpreci,2), zt2we(jpj,jpreci,2)   )
-
-      !
-      IF( PRESENT( pval ) ) THEN   ;   zland = pval      ! set land value
-      ELSE                         ;   zland = 0.e0      ! zero by default
-      ENDIF
-
-      ! 1. standard boundary treatment
-      ! ------------------------------
-      ! 2. East and west directions exchange
-      ! ------------------------------------
-      ! we play with the neigbours AND the row number because of the periodicity
-      !
-      SELECT CASE ( nbondi )      ! Read lateral conditions
-      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
-         iihom = nlci - jpreci
-         DO jl = 1, jpreci
-            zt2ew(:,jl,1) = pt2d(jl       ,:) ; pt2d(jl       ,:) = 0.0_wp
-            zt2we(:,jl,1) = pt2d(iihom +jl,:) ; pt2d(iihom +jl,:) = 0.0_wp
-         END DO
-      END SELECT
-      !
-      !                           ! Migrations
-      imigr = jpreci * jpj
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         CALL mppsend( 2, zt2we(1,1,1), imigr, noea, ml_req1 )
-         CALL mpprecv( 1, zt2ew(1,1,2), imigr, noea )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      CASE ( 0 )
-         CALL mppsend( 1, zt2ew(1,1,1), imigr, nowe, ml_req1 )
-         CALL mppsend( 2, zt2we(1,1,1), imigr, noea, ml_req2 )
-         CALL mpprecv( 1, zt2ew(1,1,2), imigr, noea )
-         CALL mpprecv( 2, zt2we(1,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
-      CASE ( 1 )
-         CALL mppsend( 1, zt2ew(1,1,1), imigr, nowe, ml_req1 )
-         CALL mpprecv( 2, zt2we(1,1,2), imigr, nowe )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      END SELECT
-      !
-      !                           ! Write lateral conditions
-      iihom = nlci-nreci
-      !
-      SELECT CASE ( nbondi )
-      CASE ( -1 )
-         DO jl = 1, jpreci
-            pt2d(iihom+jl,:) = pt2d(iihom+jl,:) + zt2ew(:,jl,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jpreci
-            pt2d(jpreci+jl,:) = pt2d(jpreci+jl,:) + zt2we(:,jl,2)
-            pt2d(iihom +jl,:) = pt2d(iihom +jl,:) + zt2ew(:,jl,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jpreci
-            pt2d(jpreci+jl,:) = pt2d(jpreci+jl,:) + zt2we(:,jl,2)
-         END DO
-      END SELECT
-
-
-      ! 3. North and south directions
-      ! -----------------------------
-      ! always closed : we play only with the neigbours
-      !
-      IF( nbondj /= 2 ) THEN      ! Read lateral conditions
-         ijhom = nlcj - jprecj
-         DO jl = 1, jprecj
-            zt2sn(:,jl,1) = pt2d(:,ijhom +jl) ; pt2d(:,ijhom +jl) = 0.0_wp
-            zt2ns(:,jl,1) = pt2d(:,jl       ) ; pt2d(:,jl       ) = 0.0_wp
-         END DO
-      ENDIF
-      !
-      !                           ! Migrations
-      imigr = jprecj * jpi
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         CALL mppsend( 4, zt2sn(1,1,1), imigr, nono, ml_req1 )
-         CALL mpprecv( 3, zt2ns(1,1,2), imigr, nono )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      CASE ( 0 )
-         CALL mppsend( 3, zt2ns(1,1,1), imigr, noso, ml_req1 )
-         CALL mppsend( 4, zt2sn(1,1,1), imigr, nono, ml_req2 )
-         CALL mpprecv( 3, zt2ns(1,1,2), imigr, nono )
-         CALL mpprecv( 4, zt2sn(1,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
-      CASE ( 1 )
-         CALL mppsend( 3, zt2ns(1,1,1), imigr, noso, ml_req1 )
-         CALL mpprecv( 4, zt2sn(1,1,2), imigr, noso )
-         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
-      END SELECT
-      !
-      !                           ! Write lateral conditions
-      ijhom = nlcj-nrecj
-      !
-      SELECT CASE ( nbondj )
-      CASE ( -1 )
-         DO jl = 1, jprecj
-            pt2d(:,ijhom+jl) = pt2d(:,ijhom+jl) + zt2ns(:,jl,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jprecj
-            pt2d(:,jprecj+jl) = pt2d(:,jprecj+jl) + zt2sn(:,jl,2)
-            pt2d(:,ijhom +jl) = pt2d(:,ijhom +jl) + zt2ns(:,jl,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jprecj
-            pt2d(:,jprecj+jl) = pt2d(:,jprecj+jl) + zt2sn(:,jl,2)
-         END DO
-      END SELECT
-
-
-      ! 4. north fold treatment
-      ! -----------------------
-      !
-      IF( npolj /= 0 .AND. .NOT. PRESENT(cd_mpp) ) THEN
-         !
-         SELECT CASE ( jpni )
-         CASE ( 1 )     ;   CALL lbc_nfd      ( pt2d, cd_type, psgn )   ! only 1 northern proc, no mpp
-         CASE DEFAULT   ;   CALL mpp_lbc_north( pt2d, cd_type, psgn )   ! for all northern procs.
-         END SELECT
-         !
-      ENDIF
-      !
-      DEALLOCATE( zt2ns, zt2sn, zt2ew, zt2we )
-      !
-   END SUBROUTINE mpp_lnk_sum_2d
+   !!    mpp_lnk_sum_2d et 3D   ====>>>>>>   à virer du code !!!!
+   
+   
+   !!----------------------------------------------------------------------
+
+
 
    SUBROUTINE mppsend( ktyp, pmess, kbytes, kdest, md_req )
@@ -1764,10 +499,10 @@
       SELECT CASE ( cn_mpi_send )
       CASE ( 'S' )                ! Standard mpi send (blocking)
-         CALL mpi_send ( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_opa        , iflag )
+         CALL mpi_send ( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_oce        , iflag )
       CASE ( 'B' )                ! Buffer mpi send (blocking)
-         CALL mpi_bsend( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_opa        , iflag )
+         CALL mpi_bsend( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_oce        , iflag )
       CASE ( 'I' )                ! Immediate mpi send (non-blocking send)
          ! be carefull, one more argument here : the mpi request identifier..
-         CALL mpi_isend( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_opa, md_req, iflag )
+         CALL mpi_isend( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_oce, md_req, iflag )
       END SELECT
       !
@@ -1797,5 +532,5 @@
       IF( PRESENT(ksource) )   use_source = ksource
       !
-      CALL mpi_recv( pmess, kbytes, mpi_double_precision, use_source, ktyp, mpi_comm_opa, istatus, iflag )
+      CALL mpi_recv( pmess, kbytes, mpi_double_precision, use_source, ktyp, mpi_comm_oce, istatus, iflag )
       !
    END SUBROUTINE mpprecv
@@ -1819,5 +554,5 @@
       itaille = jpi * jpj
       CALL mpi_gather( ptab, itaille, mpi_double_precision, pio, itaille     ,   &
-         &                            mpi_double_precision, kp , mpi_comm_opa, ierror )
+         &                            mpi_double_precision, kp , mpi_comm_oce, ierror )
       !
    END SUBROUTINE mppgather
@@ -1842,523 +577,306 @@
       !
       CALL mpi_scatter( pio, itaille, mpi_double_precision, ptab, itaille     ,   &
-         &                            mpi_double_precision, kp  , mpi_comm_opa, ierror )
+         &                            mpi_double_precision, kp  , mpi_comm_oce, ierror )
       !
    END SUBROUTINE mppscatter
 
-
-   SUBROUTINE mppmax_a_int( ktab, kdim, kcom )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppmax_a_int  ***
-      !!
-      !! ** Purpose :   Find maximum value in an integer layout array
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER , INTENT(in   )                  ::   kdim   ! size of array
-      INTEGER , INTENT(inout), DIMENSION(kdim) ::   ktab   ! input array
-      INTEGER , INTENT(in   ), OPTIONAL        ::   kcom   !
-      !
-      INTEGER :: ierror, localcomm   ! temporary integer
-      INTEGER, DIMENSION(kdim) ::   iwork
-      !!----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) )   localcomm = kcom
-      !
-      CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_max, localcomm, ierror )
-      !
-      ktab(:) = iwork(:)
-      !
-   END SUBROUTINE mppmax_a_int
-
-
-   SUBROUTINE mppmax_int( ktab, kcom )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppmax_int  ***
-      !!
-      !! ** Purpose :   Find maximum value in an integer layout array
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER, INTENT(inout)           ::   ktab   ! ???
-      INTEGER, INTENT(in   ), OPTIONAL ::   kcom   ! ???
-      !
-      INTEGER ::   ierror, iwork, localcomm   ! temporary integer
-      !!----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) )   localcomm = kcom
-      !
-      CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_max, localcomm, ierror )
-      !
-      ktab = iwork
-      !
-   END SUBROUTINE mppmax_int
-
-
-   SUBROUTINE mppmin_a_int( ktab, kdim, kcom )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppmin_a_int  ***
-      !!
-      !! ** Purpose :   Find minimum value in an integer layout array
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER , INTENT( in  )                  ::   kdim   ! size of array
-      INTEGER , INTENT(inout), DIMENSION(kdim) ::   ktab   ! input array
-      INTEGER , INTENT( in  ), OPTIONAL        ::   kcom   ! input array
-      !!
-      INTEGER ::   ierror, localcomm   ! temporary integer
-      INTEGER, DIMENSION(kdim) ::   iwork
-      !!----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) )   localcomm = kcom
-      !
-      CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_min, localcomm, ierror )
-      !
-      ktab(:) = iwork(:)
-      !
-   END SUBROUTINE mppmin_a_int
-
-
-   SUBROUTINE mppmin_int( ktab, kcom )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppmin_int  ***
-      !!
-      !! ** Purpose :   Find minimum value in an integer layout array
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER, INTENT(inout) ::   ktab      ! ???
-      INTEGER , INTENT( in  ), OPTIONAL        ::   kcom        ! input array
-      !!
-      INTEGER ::  ierror, iwork, localcomm
-      !!----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) )   localcomm = kcom
-      !
-      CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_min, localcomm, ierror )
-      !
-      ktab = iwork
-      !
-   END SUBROUTINE mppmin_int
-
-
-   SUBROUTINE mppsum_a_int( ktab, kdim )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppsum_a_int  ***
-      !!
-      !! ** Purpose :   Global integer sum, 1D array case
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER, INTENT(in   )                   ::   kdim   ! ???
-      INTEGER, INTENT(inout), DIMENSION (kdim) ::   ktab   ! ???
-      !
-      INTEGER :: ierror
-      INTEGER, DIMENSION (kdim) ::  iwork
-      !!----------------------------------------------------------------------
-      !
-      CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_sum, mpi_comm_opa, ierror )
-      !
-      ktab(:) = iwork(:)
-      !
-   END SUBROUTINE mppsum_a_int
-
-
-   SUBROUTINE mppsum_int( ktab )
-      !!----------------------------------------------------------------------
-      !!                 ***  routine mppsum_int  ***
-      !!
-      !! ** Purpose :   Global integer sum
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER, INTENT(inout) ::   ktab
-      !!
-      INTEGER :: ierror, iwork
-      !!----------------------------------------------------------------------
-      !
-      CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_sum, mpi_comm_opa, ierror )
-      !
-      ktab = iwork
-      !
-   END SUBROUTINE mppsum_int
-
-
-   SUBROUTINE mppmax_a_real( ptab, kdim, kcom )
-      !!----------------------------------------------------------------------
-      !!                 ***  routine mppmax_a_real  ***
-      !!
-      !! ** Purpose :   Maximum
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER , INTENT(in   )                  ::   kdim
-      REAL(wp), INTENT(inout), DIMENSION(kdim) ::   ptab
-      INTEGER , INTENT(in   ), OPTIONAL        ::   kcom
-      !
-      INTEGER :: ierror, localcomm
-      REAL(wp), DIMENSION(kdim) ::  zwork
-      !!----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) ) localcomm = kcom
-      !
-      CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_max, localcomm, ierror )
-      ptab(:) = zwork(:)
-      !
-   END SUBROUTINE mppmax_a_real
-
-
-   SUBROUTINE mppmax_real( ptab, kcom )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppmax_real  ***
-      !!
-      !! ** Purpose :   Maximum
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), INTENT(inout)           ::   ptab   ! ???
-      INTEGER , INTENT(in   ), OPTIONAL ::   kcom   ! ???
-      !!
-      INTEGER  ::   ierror, localcomm
-      REAL(wp) ::   zwork
-      !!----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) )   localcomm = kcom
-      !
-      CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_max, localcomm, ierror )
-      ptab = zwork
-      !
-   END SUBROUTINE mppmax_real
-
-   SUBROUTINE mppmax_real_multiple( ptab, NUM , kcom  )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppmax_real  ***
-      !!
-      !! ** Purpose :   Maximum
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(:) ,  INTENT(inout)           ::   ptab   ! ???
-      INTEGER , INTENT(in   )           ::   NUM
-      INTEGER , INTENT(in   ), OPTIONAL ::   kcom   ! ???
-      !!
-      INTEGER  ::   ierror, localcomm
-      REAL(wp) , POINTER , DIMENSION(:) ::   zwork
-      !!----------------------------------------------------------------------
-      !
-      CALL wrk_alloc(NUM , zwork)
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) )   localcomm = kcom
-      !
-      CALL mpi_allreduce( ptab, zwork, NUM, mpi_double_precision, mpi_max, localcomm, ierror )
-      ptab = zwork
-      CALL wrk_dealloc(NUM , zwork)
-      !
-   END SUBROUTINE mppmax_real_multiple
-
-
-   SUBROUTINE mppmin_a_real( ptab, kdim, kcom )
-      !!----------------------------------------------------------------------
-      !!                 ***  routine mppmin_a_real  ***
-      !!
-      !! ** Purpose :   Minimum of REAL, array case
-      !!
-      !!-----------------------------------------------------------------------
-      INTEGER , INTENT(in   )                  ::   kdim
-      REAL(wp), INTENT(inout), DIMENSION(kdim) ::   ptab
-      INTEGER , INTENT(in   ), OPTIONAL        ::   kcom
-      !!
-      INTEGER :: ierror, localcomm
-      REAL(wp), DIMENSION(kdim) ::   zwork
-      !!-----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) ) localcomm = kcom
-      !
-      CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_min, localcomm, ierror )
-      ptab(:) = zwork(:)
-      !
-   END SUBROUTINE mppmin_a_real
-
-
-   SUBROUTINE mppmin_real( ptab, kcom )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppmin_real  ***
-      !!
-      !! ** Purpose :   minimum of REAL, scalar case
-      !!
-      !!-----------------------------------------------------------------------
-      REAL(wp), INTENT(inout)           ::   ptab        !
-      INTEGER , INTENT(in   ), OPTIONAL :: kcom
-      !!
-      INTEGER  ::   ierror
-      REAL(wp) ::   zwork
-      INTEGER :: localcomm
-      !!-----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) )   localcomm = kcom
-      !
-      CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_min, localcomm, ierror )
-      ptab = zwork
-      !
-   END SUBROUTINE mppmin_real
-
-
-   SUBROUTINE mppsum_a_real( ptab, kdim, kcom )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppsum_a_real  ***
-      !!
-      !! ** Purpose :   global sum, REAL ARRAY argument case
-      !!
-      !!-----------------------------------------------------------------------
-      INTEGER , INTENT( in )                     ::   kdim      ! size of ptab
-      REAL(wp), DIMENSION(kdim), INTENT( inout ) ::   ptab      ! input array
-      INTEGER , INTENT( in ), OPTIONAL           :: kcom
-      !!
-      INTEGER                   ::   ierror    ! temporary integer
-      INTEGER                   ::   localcomm
-      REAL(wp), DIMENSION(kdim) ::   zwork     ! temporary workspace
-      !!-----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) )   localcomm = kcom
-      !
-      CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_sum, localcomm, ierror )
-      ptab(:) = zwork(:)
-      !
-   END SUBROUTINE mppsum_a_real
-
-
-   SUBROUTINE mppsum_real( ptab, kcom )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppsum_real  ***
-      !!
-      !! ** Purpose :   global sum, SCALAR argument case
-      !!
-      !!-----------------------------------------------------------------------
-      REAL(wp), INTENT(inout)           ::   ptab   ! input scalar
-      INTEGER , INTENT(in   ), OPTIONAL ::   kcom
-      !!
-      INTEGER  ::   ierror, localcomm
-      REAL(wp) ::   zwork
-      !!-----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) ) localcomm = kcom
-      !
-      CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_sum, localcomm, ierror )
-      ptab = zwork
-      !
-   END SUBROUTINE mppsum_real
-
-
-   SUBROUTINE mppsum_realdd( ytab, kcom )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppsum_realdd ***
-      !!
-      !! ** Purpose :   global sum in Massively Parallel Processing
-      !!                SCALAR argument case for double-double precision
-      !!
-      !!-----------------------------------------------------------------------
-      COMPLEX(wp), INTENT(inout)           ::   ytab    ! input scalar
-      INTEGER    , INTENT(in   ), OPTIONAL ::   kcom
-      !
-      INTEGER     ::   ierror
-      INTEGER     ::   localcomm
-      COMPLEX(wp) ::   zwork
-      !!-----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) )   localcomm = kcom
-      !
-      ! reduce local sums into global sum
-      CALL MPI_ALLREDUCE (ytab, zwork, 1, MPI_DOUBLE_COMPLEX, MPI_SUMDD, localcomm, ierror )
-      ytab = zwork
-      !
-   END SUBROUTINE mppsum_realdd
-
-
-   SUBROUTINE mppsum_a_realdd( ytab, kdim, kcom )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mppsum_a_realdd  ***
-      !!
-      !! ** Purpose :   global sum in Massively Parallel Processing
-      !!                COMPLEX ARRAY case for double-double precision
-      !!
-      !!-----------------------------------------------------------------------
-      INTEGER                     , INTENT(in   ) ::   kdim   ! size of ytab
-      COMPLEX(wp), DIMENSION(kdim), INTENT(inout) ::   ytab   ! input array
-      INTEGER    , OPTIONAL       , INTENT(in   ) ::   kcom
-      !
-      INTEGER:: ierror, localcomm    ! local integer
-      COMPLEX(wp), DIMENSION(kdim) :: zwork     ! temporary workspace
-      !!-----------------------------------------------------------------------
-      !
-      localcomm = mpi_comm_opa
-      IF( PRESENT(kcom) )   localcomm = kcom
-      !
-      CALL MPI_ALLREDUCE( ytab, zwork, kdim, MPI_DOUBLE_COMPLEX, MPI_SUMDD, localcomm, ierror )
-      ytab(:) = zwork(:)
-      !
-   END SUBROUTINE mppsum_a_realdd
-
-
-   SUBROUTINE mpp_minloc2d( ptab, pmask, pmin, ki,kj )
-      !!------------------------------------------------------------------------
-      !!             ***  routine mpp_minloc  ***
-      !!
-      !! ** Purpose :   Compute the global minimum of an array ptab
-      !!              and also give its global position
-      !!
-      !! ** Method  :   Use MPI_ALLREDUCE with MPI_MINLOC
-      !!
-      !!--------------------------------------------------------------------------
-      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   ptab    ! Local 2D array
-      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   pmask   ! Local mask
-      REAL(wp)                     , INTENT(  out) ::   pmin    ! Global minimum of ptab
-      INTEGER                      , INTENT(  out) ::   ki, kj   ! index of minimum in global frame
-      !
-      INTEGER :: ierror
-      INTEGER , DIMENSION(2)   ::   ilocs
-      REAL(wp) ::   zmin   ! local minimum
-      REAL(wp), DIMENSION(2,1) ::   zain, zaout
-      !!-----------------------------------------------------------------------
-      !
-      zmin  = MINVAL( ptab(:,:) , mask= pmask == 1.e0 )
-      ilocs = MINLOC( ptab(:,:) , mask= pmask == 1.e0 )
-      !
-      ki = ilocs(1) + nimpp - 1
-      kj = ilocs(2) + njmpp - 1
-      !
-      zain(1,:)=zmin
-      zain(2,:)=ki+10000.*kj
-      !
-      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MINLOC,MPI_COMM_OPA,ierror)
-      !
-      pmin = zaout(1,1)
-      kj = INT(zaout(2,1)/10000.)
-      ki = INT(zaout(2,1) - 10000.*kj )
-      !
-   END SUBROUTINE mpp_minloc2d
-
-
-   SUBROUTINE mpp_minloc3d( ptab, pmask, pmin, ki, kj ,kk)
-      !!------------------------------------------------------------------------
-      !!             ***  routine mpp_minloc  ***
-      !!
-      !! ** Purpose :   Compute the global minimum of an array ptab
-      !!              and also give its global position
-      !!
-      !! ** Method  :   Use MPI_ALLREDUCE with MPI_MINLOC
-      !!
-      !!--------------------------------------------------------------------------
-      REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in   ) ::   ptab         ! Local 2D array
-      REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in   ) ::   pmask        ! Local mask
-      REAL(wp)                         , INTENT(  out) ::   pmin         ! Global minimum of ptab
-      INTEGER                          , INTENT(  out) ::   ki, kj, kk   ! index of minimum in global frame
-      !!
-      INTEGER  ::   ierror
-      REAL(wp) ::   zmin     ! local minimum
-      INTEGER , DIMENSION(3)   ::   ilocs
-      REAL(wp), DIMENSION(2,1) ::   zain, zaout
-      !!-----------------------------------------------------------------------
-      !
-      zmin  = MINVAL( ptab(:,:,:) , mask= pmask == 1.e0 )
-      ilocs = MINLOC( ptab(:,:,:) , mask= pmask == 1.e0 )
-      !
-      ki = ilocs(1) + nimpp - 1
-      kj = ilocs(2) + njmpp - 1
-      kk = ilocs(3)
-      !
-      zain(1,:)=zmin
-      zain(2,:)=ki+10000.*kj+100000000.*kk
-      !
-      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MINLOC,MPI_COMM_OPA,ierror)
-      !
-      pmin = zaout(1,1)
-      kk   = INT( zaout(2,1) / 100000000. )
-      kj   = INT( zaout(2,1) - kk * 100000000. ) / 10000
-      ki   = INT( zaout(2,1) - kk * 100000000. -kj * 10000. )
-      !
-   END SUBROUTINE mpp_minloc3d
-
-
-   SUBROUTINE mpp_maxloc2d( ptab, pmask, pmax, ki, kj )
-      !!------------------------------------------------------------------------
-      !!             ***  routine mpp_maxloc  ***
-      !!
-      !! ** Purpose :   Compute the global maximum of an array ptab
-      !!              and also give its global position
-      !!
-      !! ** Method  :   Use MPI_ALLREDUCE with MPI_MINLOC
-      !!
-      !!--------------------------------------------------------------------------
-      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   ptab     ! Local 2D array
-      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   pmask    ! Local mask
-      REAL(wp)                     , INTENT(  out) ::   pmax     ! Global maximum of ptab
-      INTEGER                      , INTENT(  out) ::   ki, kj   ! index of maximum in global frame
-      !!
-      INTEGER  :: ierror
-      INTEGER, DIMENSION (2)   ::   ilocs
-      REAL(wp) :: zmax   ! local maximum
-      REAL(wp), DIMENSION(2,1) ::   zain, zaout
-      !!-----------------------------------------------------------------------
-      !
-      zmax  = MAXVAL( ptab(:,:) , mask= pmask == 1.e0 )
-      ilocs = MAXLOC( ptab(:,:) , mask= pmask == 1.e0 )
-      !
-      ki = ilocs(1) + nimpp - 1
-      kj = ilocs(2) + njmpp - 1
-      !
-      zain(1,:) = zmax
-      zain(2,:) = ki + 10000. * kj
-      !
-      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MAXLOC,MPI_COMM_OPA,ierror)
-      !
-      pmax = zaout(1,1)
-      kj   = INT( zaout(2,1) / 10000.     )
-      ki   = INT( zaout(2,1) - 10000.* kj )
-      !
-   END SUBROUTINE mpp_maxloc2d
-
-
-   SUBROUTINE mpp_maxloc3d( ptab, pmask, pmax, ki, kj, kk )
-      !!------------------------------------------------------------------------
-      !!             ***  routine mpp_maxloc  ***
-      !!
-      !! ** Purpose :  Compute the global maximum of an array ptab
-      !!              and also give its global position
-      !!
-      !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC
-      !!
-      !!--------------------------------------------------------------------------
-      REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in   ) ::   ptab         ! Local 2D array
-      REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in   ) ::   pmask        ! Local mask
-      REAL(wp)                         , INTENT(  out) ::   pmax         ! Global maximum of ptab
-      INTEGER                          , INTENT(  out) ::   ki, kj, kk   ! index of maximum in global frame
-      !!
-      REAL(wp) :: zmax   ! local maximum
-      REAL(wp), DIMENSION(2,1) ::   zain, zaout
-      INTEGER , DIMENSION(3)   ::   ilocs
-      INTEGER :: ierror
-      !!-----------------------------------------------------------------------
-      !
-      zmax  = MAXVAL( ptab(:,:,:) , mask= pmask == 1.e0 )
-      ilocs = MAXLOC( ptab(:,:,:) , mask= pmask == 1.e0 )
-      !
-      ki = ilocs(1) + nimpp - 1
-      kj = ilocs(2) + njmpp - 1
-      kk = ilocs(3)
-      !
-      zain(1,:)=zmax
-      zain(2,:)=ki+10000.*kj+100000000.*kk
-      !
-      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MAXLOC,MPI_COMM_OPA,ierror)
-      !
-      pmax = zaout(1,1)
-      kk   = INT( zaout(2,1) / 100000000. )
-      kj   = INT( zaout(2,1) - kk * 100000000. ) / 10000
-      ki   = INT( zaout(2,1) - kk * 100000000. -kj * 10000. )
-      !
-   END SUBROUTINE mpp_maxloc3d
-
+   
+   SUBROUTINE mpp_delay_sum( cdname, cdelay, y_in, pout, ldlast, kcom )
+     !!----------------------------------------------------------------------
+      !!                   ***  routine mpp_delay_sum  ***
+      !!
+      !! ** Purpose :   performed delayed mpp_sum, the result is received on next call
+      !!
+      !!----------------------------------------------------------------------
+      CHARACTER(len=*), INTENT(in   )               ::   cdname  ! name of the calling subroutine
+      CHARACTER(len=*), INTENT(in   )               ::   cdelay  ! name (used as id) of the delayed operation
+      COMPLEX(wp),      INTENT(in   ), DIMENSION(:) ::   y_in
+      REAL(wp),         INTENT(  out), DIMENSION(:) ::   pout
+      LOGICAL,          INTENT(in   )               ::   ldlast  ! true if this is the last time we call this routine
+      INTEGER,          INTENT(in   ), OPTIONAL     ::   kcom
+      !!
+      INTEGER ::   ji, isz
+      INTEGER ::   idvar
+      INTEGER ::   ierr, ilocalcomm
+      COMPLEX(wp), ALLOCATABLE, DIMENSION(:) ::   ytmp
+      !!----------------------------------------------------------------------
+      ilocalcomm = mpi_comm_oce
+      IF( PRESENT(kcom) )   ilocalcomm = kcom
+
+      isz = SIZE(y_in)
+      
+      IF( narea == 1 .AND. numcom == -1 ) CALL mpp_report( cdname, ld_dlg = .TRUE. )
+
+      idvar = -1
+      DO ji = 1, nbdelay
+         IF( TRIM(cdelay) == TRIM(c_delaylist(ji)) ) idvar = ji
+      END DO
+      IF ( idvar == -1 )   CALL ctl_stop( 'STOP',' mpp_delay_sum : please add a new delayed exchange for '//TRIM(cdname) )
+
+      IF ( ndelayid(idvar) == 0 ) THEN         ! first call    with restart: %z1d defined in iom_delay_rst
+         !                                       --------------------------
+         IF ( SIZE(todelay(idvar)%z1d) /= isz ) THEN                  ! Check dimension coherence
+            IF(lwp) WRITE(numout,*) ' WARNING: the nb of delayed variables in restart file is not the model one'
+            DEALLOCATE(todelay(idvar)%z1d)
+            ndelayid(idvar) = -1                                      ! do as if we had no restart
+         ELSE
+            ALLOCATE(todelay(idvar)%y1d(isz))
+            todelay(idvar)%y1d(:) = CMPLX(todelay(idvar)%z1d(:), 0., wp)   ! create %y1d, complex variable needed by mpi_sumdd
+         END IF
+      ENDIF
+      
+      IF( ndelayid(idvar) == -1 ) THEN         ! first call without restart: define %y1d and %z1d from y_in with blocking allreduce
+         !                                       --------------------------
+         ALLOCATE(todelay(idvar)%z1d(isz), todelay(idvar)%y1d(isz))   ! allocate also %z1d as used for the restart
+         CALL mpi_allreduce( y_in(:), todelay(idvar)%y1d(:), isz, MPI_DOUBLE_COMPLEX, mpi_sumdd, ilocalcomm, ierr )   ! get %y1d
+         todelay(idvar)%z1d(:) = REAL(todelay(idvar)%y1d(:), wp)      ! define %z1d from %y1d
+      ENDIF
+
+      IF( ndelayid(idvar) > 0 )   CALL mpp_delay_rcv( idvar )         ! make sure %z1d is received
+
+      ! send back pout from todelay(idvar)%z1d defined at previous call
+      pout(:) = todelay(idvar)%z1d(:)
+
+      ! send y_in into todelay(idvar)%y1d with a non-blocking communication
+#if defined key_mpi2
+      IF( ln_timing ) CALL tic_tac( .TRUE., ld_global = .TRUE.)
+      CALL  mpi_allreduce( y_in(:), todelay(idvar)%y1d(:), isz, MPI_DOUBLE_COMPLEX, mpi_sumdd, ilocalcomm, ndelayid(idvar), ierr )
+      IF( ln_timing ) CALL tic_tac(.FALSE., ld_global = .TRUE.)
+#else
+      CALL mpi_iallreduce( y_in(:), todelay(idvar)%y1d(:), isz, MPI_DOUBLE_COMPLEX, mpi_sumdd, ilocalcomm, ndelayid(idvar), ierr )
+#endif
+
+   END SUBROUTINE mpp_delay_sum
+
+   
+   SUBROUTINE mpp_delay_max( cdname, cdelay, p_in, pout, ldlast, kcom )
+      !!----------------------------------------------------------------------
+      !!                   ***  routine mpp_delay_max  ***
+      !!
+      !! ** Purpose :   performed delayed mpp_max, the result is received on next call
+      !!
+      !!----------------------------------------------------------------------
+      CHARACTER(len=*), INTENT(in   )                 ::   cdname  ! name of the calling subroutine
+      CHARACTER(len=*), INTENT(in   )                 ::   cdelay  ! name (used as id) of the delayed operation
+      REAL(wp),         INTENT(in   ), DIMENSION(:)   ::   p_in    ! 
+      REAL(wp),         INTENT(  out), DIMENSION(:)   ::   pout    ! 
+      LOGICAL,          INTENT(in   )                 ::   ldlast  ! true if this is the last time we call this routine
+      INTEGER,          INTENT(in   ), OPTIONAL       ::   kcom
+      !!
+      INTEGER ::   ji, isz
+      INTEGER ::   idvar
+      INTEGER ::   ierr, ilocalcomm
+      !!----------------------------------------------------------------------
+      ilocalcomm = mpi_comm_oce
+      IF( PRESENT(kcom) )   ilocalcomm = kcom
+
+      isz = SIZE(p_in)
+
+      IF( narea == 1 .AND. numcom == -1 ) CALL mpp_report( cdname, ld_dlg = .TRUE. )
+
+      idvar = -1
+      DO ji = 1, nbdelay
+         IF( TRIM(cdelay) == TRIM(c_delaylist(ji)) ) idvar = ji
+      END DO
+      IF ( idvar == -1 )   CALL ctl_stop( 'STOP',' mpp_delay_max : please add a new delayed exchange for '//TRIM(cdname) )
+
+      IF ( ndelayid(idvar) == 0 ) THEN         ! first call    with restart: %z1d defined in iom_delay_rst
+         !                                       --------------------------
+         IF ( SIZE(todelay(idvar)%z1d) /= isz ) THEN                  ! Check dimension coherence
+            IF(lwp) WRITE(numout,*) ' WARNING: the nb of delayed variables in restart file is not the model one'
+            DEALLOCATE(todelay(idvar)%z1d)
+            ndelayid(idvar) = -1                                      ! do as if we had no restart
+         END IF
+      ENDIF
+
+      IF( ndelayid(idvar) == -1 ) THEN         ! first call without restart: define %z1d from p_in with a blocking allreduce
+         !                                       --------------------------
+         ALLOCATE(todelay(idvar)%z1d(isz))
+         CALL mpi_allreduce( p_in(:), todelay(idvar)%z1d(:), isz, MPI_DOUBLE_PRECISION, mpi_max, ilocalcomm, ierr )   ! get %z1d
+      ENDIF
+
+      IF( ndelayid(idvar) > 0 )   CALL mpp_delay_rcv( idvar )         ! make sure %z1d is received
+
+      ! send back pout from todelay(idvar)%z1d defined at previous call
+      pout(:) = todelay(idvar)%z1d(:)
+
+      ! send p_in into todelay(idvar)%z1d with a non-blocking communication
+#if defined key_mpi2
+      IF( ln_timing ) CALL tic_tac( .TRUE., ld_global = .TRUE.)
+      CALL  mpi_allreduce( p_in(:), todelay(idvar)%z1d(:), isz, MPI_DOUBLE_PRECISION, mpi_max, ilocalcomm, ndelayid(idvar), ierr )
+      IF( ln_timing ) CALL tic_tac(.FALSE., ld_global = .TRUE.)
+#else
+      CALL mpi_iallreduce( p_in(:), todelay(idvar)%z1d(:), isz, MPI_DOUBLE_PRECISION, mpi_max, ilocalcomm, ndelayid(idvar), ierr )
+#endif
+
+   END SUBROUTINE mpp_delay_max
+
+   
+   SUBROUTINE mpp_delay_rcv( kid )
+      !!----------------------------------------------------------------------
+      !!                   ***  routine mpp_delay_rcv  ***
+      !!
+      !! ** Purpose :  force barrier for delayed mpp (needed for restart) 
+      !!
+      !!----------------------------------------------------------------------
+      INTEGER,INTENT(in   )      ::  kid 
+      INTEGER ::   ierr
+      !!----------------------------------------------------------------------
+      IF( ndelayid(kid) /= -2 ) THEN  
+#if ! defined key_mpi2
+         IF( ln_timing ) CALL tic_tac( .TRUE., ld_global = .TRUE.)
+         CALL mpi_wait( ndelayid(kid), MPI_STATUS_IGNORE, ierr )                        ! make sure todelay(kid) is received
+         IF( ln_timing ) CALL tic_tac(.FALSE., ld_global = .TRUE.)
+#endif
+         IF( ASSOCIATED(todelay(kid)%y1d) )   todelay(kid)%z1d(:) = REAL(todelay(kid)%y1d(:), wp)  ! define %z1d from %y1d
+         ndelayid(kid) = -2   ! add flag to know that mpi_wait was already called on kid
+      ENDIF
+   END SUBROUTINE mpp_delay_rcv
+
+   
+   !!----------------------------------------------------------------------
+   !!    ***  mppmax_a_int, mppmax_int, mppmax_a_real, mppmax_real  ***
+   !!   
+   !!----------------------------------------------------------------------
+   !!
+#  define OPERATION_MAX
+#  define INTEGER_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppmax_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppmax_a_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef INTEGER_TYPE
+!
+#  define REAL_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppmax_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppmax_a_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef REAL_TYPE
+#  undef OPERATION_MAX
+   !!----------------------------------------------------------------------
+   !!    ***  mppmin_a_int, mppmin_int, mppmin_a_real, mppmin_real  ***
+   !!   
+   !!----------------------------------------------------------------------
+   !!
+#  define OPERATION_MIN
+#  define INTEGER_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppmin_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppmin_a_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef INTEGER_TYPE
+!
+#  define REAL_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppmin_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppmin_a_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef REAL_TYPE
+#  undef OPERATION_MIN
+
+   !!----------------------------------------------------------------------
+   !!    ***  mppsum_a_int, mppsum_int, mppsum_a_real, mppsum_real  ***
+   !!   
+   !!   Global sum of 1D array or a variable (integer, real or complex)
+   !!----------------------------------------------------------------------
+   !!
+#  define OPERATION_SUM
+#  define INTEGER_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppsum_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppsum_a_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef INTEGER_TYPE
+!
+#  define REAL_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppsum_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppsum_a_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef REAL_TYPE
+#  undef OPERATION_SUM
+
+#  define OPERATION_SUM_DD
+#  define COMPLEX_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppsum_realdd
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppsum_a_realdd
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef COMPLEX_TYPE
+#  undef OPERATION_SUM_DD
+
+   !!----------------------------------------------------------------------
+   !!    ***  mpp_minloc2d, mpp_minloc3d, mpp_maxloc2d, mpp_maxloc3d
+   !!   
+   !!----------------------------------------------------------------------
+   !!
+#  define OPERATION_MINLOC
+#  define DIM_2d
+#     define ROUTINE_LOC           mpp_minloc2d
+#     include "mpp_loc_generic.h90"
+#     undef ROUTINE_LOC
+#  undef DIM_2d
+#  define DIM_3d
+#     define ROUTINE_LOC           mpp_minloc3d
+#     include "mpp_loc_generic.h90"
+#     undef ROUTINE_LOC
+#  undef DIM_3d
+#  undef OPERATION_MINLOC
+
+#  define OPERATION_MAXLOC
+#  define DIM_2d
+#     define ROUTINE_LOC           mpp_maxloc2d
+#     include "mpp_loc_generic.h90"
+#     undef ROUTINE_LOC
+#  undef DIM_2d
+#  define DIM_3d
+#     define ROUTINE_LOC           mpp_maxloc3d
+#     include "mpp_loc_generic.h90"
+#     undef ROUTINE_LOC
+#  undef DIM_3d
+#  undef OPERATION_MAXLOC
 
    SUBROUTINE mppsync()
@@ -2372,10 +890,10 @@
       !!-----------------------------------------------------------------------
       !
-      CALL mpi_barrier( mpi_comm_opa, ierror )
+      CALL mpi_barrier( mpi_comm_oce, ierror )
       !
    END SUBROUTINE mppsync
 
 
-   SUBROUTINE mppstop
+   SUBROUTINE mppstop( ldfinal, ld_force_abort ) 
       !!----------------------------------------------------------------------
       !!                  ***  routine mppstop  ***
@@ -2384,9 +902,21 @@
       !!
       !!----------------------------------------------------------------------
+      LOGICAL, OPTIONAL, INTENT(in) :: ldfinal    ! source process number
+      LOGICAL, OPTIONAL, INTENT(in) :: ld_force_abort    ! source process number
+      LOGICAL ::   llfinal, ll_force_abort
       INTEGER ::   info
       !!----------------------------------------------------------------------
-      !
-      CALL mppsync
-      CALL mpi_finalize( info )
+      llfinal = .FALSE.
+      IF( PRESENT(ldfinal) ) llfinal = ldfinal
+      ll_force_abort = .FALSE.
+      IF( PRESENT(ld_force_abort) ) ll_force_abort = ld_force_abort
+      !
+      IF(ll_force_abort) THEN
+         CALL mpi_abort( MPI_COMM_WORLD )
+      ELSE
+         CALL mppsync
+         CALL mpi_finalize( info )
+      ENDIF
+      IF( .NOT. llfinal ) STOP 123456
       !
    END SUBROUTINE mppstop
@@ -2395,5 +925,4 @@
    SUBROUTINE mpp_comm_free( kcom )
       !!----------------------------------------------------------------------
-      !!----------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kcom
       !!
@@ -2404,89 +933,4 @@
       !
    END SUBROUTINE mpp_comm_free
-
-
-   SUBROUTINE mpp_ini_ice( pindic, kumout )
-      !!----------------------------------------------------------------------
-      !!               ***  routine mpp_ini_ice  ***
-      !!
-      !! ** Purpose :   Initialize special communicator for ice areas
-      !!      condition together with global variables needed in the ddmpp folding
-      !!
-      !! ** Method  : - Look for ice processors in ice routines
-      !!              - Put their number in nrank_ice
-      !!              - Create groups for the world processors and the ice processors
-      !!              - Create a communicator for ice processors
-      !!
-      !! ** output
-      !!      njmppmax = njmpp for northern procs
-      !!      ndim_rank_ice = number of processors with ice
-      !!      nrank_ice (ndim_rank_ice) = ice processors
-      !!      ngrp_iworld = group ID for the world processors
-      !!      ngrp_ice = group ID for the ice processors
-      !!      ncomm_ice = communicator for the ice procs.
-      !!      n_ice_root = number (in the world) of proc 0 in the ice comm.
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER, INTENT(in) ::   pindic
-      INTEGER, INTENT(in) ::   kumout   ! ocean.output logical unit
-      !!
-      INTEGER :: jjproc
-      INTEGER :: ii, ierr
-      INTEGER, ALLOCATABLE, DIMENSION(:) ::   kice
-      INTEGER, ALLOCATABLE, DIMENSION(:) ::   zwork
-      !!----------------------------------------------------------------------
-      !
-      ! Since this is just an init routine and these arrays are of length jpnij
-      ! then don't use wrk_nemo module - just allocate and deallocate.
-      ALLOCATE( kice(jpnij), zwork(jpnij), STAT=ierr )
-      IF( ierr /= 0 ) THEN
-         WRITE(kumout, cform_err)
-         WRITE(kumout,*) 'mpp_ini_ice : failed to allocate 2, 1D arrays (jpnij in length)'
-         CALL mppstop
-      ENDIF
-
-      ! Look for how many procs with sea-ice
-      !
-      kice = 0
-      DO jjproc = 1, jpnij
-         IF( jjproc == narea .AND. pindic .GT. 0 )   kice(jjproc) = 1
-      END DO
-      !
-      zwork = 0
-      CALL MPI_ALLREDUCE( kice, zwork, jpnij, mpi_integer, mpi_sum, mpi_comm_opa, ierr )
-      ndim_rank_ice = SUM( zwork )
-
-      ! Allocate the right size to nrank_north
-      IF( ALLOCATED ( nrank_ice ) )   DEALLOCATE( nrank_ice )
-      ALLOCATE( nrank_ice(ndim_rank_ice) )
-      !
-      ii = 0
-      nrank_ice = 0
-      DO jjproc = 1, jpnij
-         IF( zwork(jjproc) == 1) THEN
-            ii = ii + 1
-            nrank_ice(ii) = jjproc -1
-         ENDIF
-      END DO
-
-      ! Create the world group
-      CALL MPI_COMM_GROUP( mpi_comm_opa, ngrp_iworld, ierr )
-
-      ! Create the ice group from the world group
-      CALL MPI_GROUP_INCL( ngrp_iworld, ndim_rank_ice, nrank_ice, ngrp_ice, ierr )
-
-      ! Create the ice communicator , ie the pool of procs with sea-ice
-      CALL MPI_COMM_CREATE( mpi_comm_opa, ngrp_ice, ncomm_ice, ierr )
-
-      ! Find proc number in the world of proc 0 in the north
-      ! The following line seems to be useless, we just comment & keep it as reminder
-      ! CALL MPI_GROUP_TRANSLATE_RANKS(ngrp_ice,1,0,ngrp_iworld,n_ice_root,ierr)
-      !
-      CALL MPI_GROUP_FREE(ngrp_ice, ierr)
-      CALL MPI_GROUP_FREE(ngrp_iworld, ierr)
-
-      DEALLOCATE(kice, zwork)
-      !
-   END SUBROUTINE mpp_ini_ice
 
 
@@ -2518,5 +962,5 @@
       !-$$     WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ngrp_world     : ', ngrp_world
       !-$$     WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - mpi_comm_world : ', mpi_comm_world
-      !-$$     WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - mpi_comm_opa   : ', mpi_comm_opa
+      !-$$     WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - mpi_comm_oce   : ', mpi_comm_oce
       !
       ALLOCATE( kwork(jpnij), STAT=ierr )
@@ -2529,8 +973,8 @@
       IF( jpnj == 1 ) THEN
          ngrp_znl  = ngrp_world
-         ncomm_znl = mpi_comm_opa
+         ncomm_znl = mpi_comm_oce
       ELSE
          !
-         CALL MPI_ALLGATHER ( njmpp, 1, mpi_integer, kwork, 1, mpi_integer, mpi_comm_opa, ierr )
+         CALL MPI_ALLGATHER ( njmpp, 1, mpi_integer, kwork, 1, mpi_integer, mpi_comm_oce, ierr )
          !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - kwork pour njmpp : ', kwork
          !-$$        CALL flush(numout)
@@ -2560,5 +1004,5 @@
 
          ! Create the opa group
-         CALL MPI_COMM_GROUP(mpi_comm_opa,ngrp_opa,ierr)
+         CALL MPI_COMM_GROUP(mpi_comm_oce,ngrp_opa,ierr)
          !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ngrp_opa : ', ngrp_opa
          !-$$        CALL flush(numout)
@@ -2570,5 +1014,5 @@
 
          ! Create the znl communicator from the opa communicator, ie the pool of procs in the same row
-         CALL MPI_COMM_CREATE ( mpi_comm_opa, ngrp_znl, ncomm_znl, ierr )
+         CALL MPI_COMM_CREATE ( mpi_comm_oce, ngrp_znl, ncomm_znl, ierr )
          !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ncomm_znl ', ncomm_znl
          !-$$        CALL flush(numout)
@@ -2582,5 +1026,5 @@
          l_znl_root = .FALSE.
          kwork (1) = nimpp
-         CALL mpp_min ( kwork(1), kcom = ncomm_znl)
+         CALL mpp_min ( 'lib_mpp', kwork(1), kcom = ncomm_znl)
          IF ( nimpp == kwork(1)) l_znl_root = .TRUE.
       END IF
@@ -2641,5 +1085,5 @@
       !
       ! create the world group
-      CALL MPI_COMM_GROUP( mpi_comm_opa, ngrp_world, ierr )
+      CALL MPI_COMM_GROUP( mpi_comm_oce, ngrp_world, ierr )
       !
       ! Create the North group from the world group
@@ -2647,976 +1091,10 @@
       !
       ! Create the North communicator , ie the pool of procs in the north group
-      CALL MPI_COMM_CREATE( mpi_comm_opa, ngrp_north, ncomm_north, ierr )
+      CALL MPI_COMM_CREATE( mpi_comm_oce, ngrp_north, ncomm_north, ierr )
       !
    END SUBROUTINE mpp_ini_north
 
 
-   SUBROUTINE mpp_lbc_north_3d( pt3d, cd_type, psgn )
-      !!---------------------------------------------------------------------
-      !!                   ***  routine mpp_lbc_north_3d  ***
-      !!
-      !! ** Purpose :   Ensure proper north fold horizontal bondary condition
-      !!              in mpp configuration in case of jpn1 > 1
-      !!
-      !! ** Method  :   North fold condition and mpp with more than one proc
-      !!              in i-direction require a specific treatment. We gather
-      !!              the 4 northern lines of the global domain on 1 processor
-      !!              and apply lbc north-fold on this sub array. Then we
-      !!              scatter the north fold array back to the processors.
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pt3d      ! 3D array on which the b.c. is applied
-      CHARACTER(len=1)                , INTENT(in   ) ::   cd_type   ! nature of pt3d grid-points
-      !                                                              !   = T ,  U , V , F or W  gridpoints
-      REAL(wp)                        , INTENT(in   ) ::   psgn      ! = -1. the sign change across the north fold 
-      !!                                                             ! =  1. , the sign is kept
-      INTEGER ::   ji, jj, jr, jk
-      INTEGER ::   ierr, itaille, ildi, ilei, iilb
-      INTEGER ::   ijpj, ijpjm1, ij, iproc
-      INTEGER, DIMENSION (jpmaxngh)          ::   ml_req_nf          !for mpi_isend when avoiding mpi_allgather
-      INTEGER                                ::   ml_err             ! for mpi_isend when avoiding mpi_allgather
-      INTEGER, DIMENSION(MPI_STATUS_SIZE)    ::   ml_stat            ! for mpi_isend when avoiding mpi_allgather
-      !                                                              ! Workspace for message transfers avoiding mpi_allgather
-      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE   :: ztab
-      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE   :: znorthloc, zfoldwk      
-      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE   :: znorthgloio
-      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE   :: ztabl, ztabr
-
-      INTEGER :: istatus(mpi_status_size)
-      INTEGER :: iflag
-      !!----------------------------------------------------------------------
-      !
-      ALLOCATE( ztab(jpiglo,4,jpk) , znorthloc(jpi,4,jpk), zfoldwk(jpi,4,jpk), znorthgloio(jpi,4,jpk,jpni) )
-      ALLOCATE( ztabl(jpi,4,jpk), ztabr(jpi*jpmaxngh, 4, jpk) ) 
-
-      ijpj   = 4
-      ijpjm1 = 3
-      !
-      znorthloc(:,:,:) = 0
-      DO jk = 1, jpk
-         DO jj = nlcj - ijpj +1, nlcj          ! put in xnorthloc the last 4 jlines of pt3d
-            ij = jj - nlcj + ijpj
-            znorthloc(:,ij,jk) = pt3d(:,jj,jk)
-         END DO
-      END DO
-      !
-      !                                     ! Build in procs of ncomm_north the znorthgloio
-      itaille = jpi * jpk * ijpj
-
-      IF ( l_north_nogather ) THEN
-         !
-        ztabr(:,:,:) = 0
-        ztabl(:,:,:) = 0
-
-        DO jk = 1, jpk
-           DO jj = nlcj-ijpj+1, nlcj          ! First put local values into the global array
-              ij = jj - nlcj + ijpj
-              DO ji = nfsloop, nfeloop
-                 ztabl(ji,ij,jk) = pt3d(ji,jj,jk)
-              END DO
-           END DO
-        END DO
-
-         DO jr = 1,nsndto
-            IF ((nfipproc(isendto(jr),jpnj) .ne. (narea-1)) .and. (nfipproc(isendto(jr),jpnj) .ne. -1)) THEN
-              CALL mppsend( 5, znorthloc, itaille, nfipproc(isendto(jr),jpnj), ml_req_nf(jr) )
-            ENDIF
-         END DO
-         DO jr = 1,nsndto
-            iproc = nfipproc(isendto(jr),jpnj)
-            IF(iproc .ne. -1) THEN
-               ilei = nleit (iproc+1)
-               ildi = nldit (iproc+1)
-               iilb = nfiimpp(isendto(jr),jpnj) - nfiimpp(isendto(1),jpnj)
-            ENDIF
-            IF((iproc .ne. (narea-1)) .and. (iproc .ne. -1)) THEN
-              CALL mpprecv(5, zfoldwk, itaille, iproc)
-              DO jk = 1, jpk
-                 DO jj = 1, ijpj
-                    DO ji = ildi, ilei
-                       ztabr(iilb+ji,jj,jk) = zfoldwk(ji,jj,jk)
-                    END DO
-                 END DO
-              END DO
-           ELSE IF (iproc .eq. (narea-1)) THEN
-              DO jk = 1, jpk
-                 DO jj = 1, ijpj
-                    DO ji = ildi, ilei
-                       ztabr(iilb+ji,jj,jk) = pt3d(ji,nlcj-ijpj+jj,jk)
-                    END DO
-                 END DO
-              END DO
-           ENDIF
-         END DO
-         IF (l_isend) THEN
-            DO jr = 1,nsndto
-               IF ((nfipproc(isendto(jr),jpnj) .ne. (narea-1)) .and. (nfipproc(isendto(jr),jpnj) .ne. -1)) THEN
-                  CALL mpi_wait(ml_req_nf(jr), ml_stat, ml_err)
-               ENDIF    
-            END DO
-         ENDIF
-         CALL mpp_lbc_nfd( ztabl, ztabr, cd_type, psgn )   ! North fold boundary condition
-         DO jk = 1, jpk
-            DO jj = nlcj-ijpj+1, nlcj             ! Scatter back to pt3d
-               ij = jj - nlcj + ijpj
-               DO ji= 1, nlci
-                  pt3d(ji,jj,jk) = ztabl(ji,ij,jk)
-               END DO
-            END DO
-         END DO
-         !
-
-      ELSE
-         CALL MPI_ALLGATHER( znorthloc  , itaille, MPI_DOUBLE_PRECISION,                &
-            &                znorthgloio, itaille, MPI_DOUBLE_PRECISION, ncomm_north, ierr )
-         !
-         ztab(:,:,:) = 0.e0
-         DO jr = 1, ndim_rank_north         ! recover the global north array
-            iproc = nrank_north(jr) + 1
-            ildi  = nldit (iproc)
-            ilei  = nleit (iproc)
-            iilb  = nimppt(iproc)
-            DO jk = 1, jpk
-               DO jj = 1, ijpj
-                  DO ji = ildi, ilei
-                    ztab(ji+iilb-1,jj,jk) = znorthgloio(ji,jj,jk,jr)
-                  END DO
-               END DO
-            END DO
-         END DO
-         CALL lbc_nfd( ztab, cd_type, psgn )   ! North fold boundary condition
-         !
-         DO jk = 1, jpk
-            DO jj = nlcj-ijpj+1, nlcj             ! Scatter back to pt3d
-               ij = jj - nlcj + ijpj
-               DO ji= 1, nlci
-                  pt3d(ji,jj,jk) = ztab(ji+nimpp-1,ij,jk)
-               END DO
-            END DO
-         END DO
-         !
-      ENDIF
-      !
-      ! The ztab array has been either:
-      !  a. Fully populated by the mpi_allgather operation or
-      !  b. Had the active points for this domain and northern neighbours populated
-      !     by peer to peer exchanges
-      ! Either way the array may be folded by lbc_nfd and the result for the span of
-      ! this domain will be identical.
-      !
-      DEALLOCATE( ztab, znorthloc, zfoldwk, znorthgloio )
-      DEALLOCATE( ztabl, ztabr ) 
-      !
-   END SUBROUTINE mpp_lbc_north_3d
-
-
-   SUBROUTINE mpp_lbc_north_2d( pt2d, cd_type, psgn)
-      !!---------------------------------------------------------------------
-      !!                   ***  routine mpp_lbc_north_2d  ***
-      !!
-      !! ** Purpose :   Ensure proper north fold horizontal bondary condition
-      !!              in mpp configuration in case of jpn1 > 1 (for 2d array )
-      !!
-      !! ** Method  :   North fold condition and mpp with more than one proc
-      !!              in i-direction require a specific treatment. We gather
-      !!              the 4 northern lines of the global domain on 1 processor
-      !!              and apply lbc north-fold on this sub array. Then we
-      !!              scatter the north fold array back to the processors.
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   pt2d      ! 2D array on which the b.c. is applied
-      CHARACTER(len=1)            , INTENT(in   ) ::   cd_type   ! nature of pt2d grid-points
-      !                                                          !   = T ,  U , V , F or W  gridpoints
-      REAL(wp)                    , INTENT(in   ) ::   psgn      ! = -1. the sign change across the north fold 
-      !!                                                             ! =  1. , the sign is kept
-      INTEGER ::   ji, jj, jr
-      INTEGER ::   ierr, itaille, ildi, ilei, iilb
-      INTEGER ::   ijpj, ijpjm1, ij, iproc
-      INTEGER, DIMENSION (jpmaxngh)      ::   ml_req_nf          !for mpi_isend when avoiding mpi_allgather
-      INTEGER                            ::   ml_err             ! for mpi_isend when avoiding mpi_allgather
-      INTEGER, DIMENSION(MPI_STATUS_SIZE)::   ml_stat            ! for mpi_isend when avoiding mpi_allgather
-      !                                                              ! Workspace for message transfers avoiding mpi_allgather
-      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE   :: ztab
-      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE   :: znorthloc, zfoldwk      
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE   :: znorthgloio
-      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE   :: ztabl, ztabr
-      INTEGER :: istatus(mpi_status_size)
-      INTEGER :: iflag
-      !!----------------------------------------------------------------------
-      !
-      ALLOCATE( ztab(jpiglo,4), znorthloc(jpi,4), zfoldwk(jpi,4), znorthgloio(jpi,4,jpni) )
-      ALLOCATE( ztabl(jpi,4), ztabr(jpi*jpmaxngh, 4) ) 
-      !
-      ijpj   = 4
-      ijpjm1 = 3
-      !
-      DO jj = nlcj-ijpj+1, nlcj             ! put in znorthloc the last 4 jlines of pt2d
-         ij = jj - nlcj + ijpj
-         znorthloc(:,ij) = pt2d(:,jj)
-      END DO
-
-      !                                     ! Build in procs of ncomm_north the znorthgloio
-      itaille = jpi * ijpj
-      IF ( l_north_nogather ) THEN
-         !
-         ! Avoid the use of mpi_allgather by exchanging only with the processes already identified 
-         ! (in nemo_northcomms) as being  involved in this process' northern boundary exchange
-         !
-         ztabr(:,:) = 0
-         ztabl(:,:) = 0
-
-         DO jj = nlcj-ijpj+1, nlcj          ! First put local values into the global array
-            ij = jj - nlcj + ijpj
-              DO ji = nfsloop, nfeloop
-               ztabl(ji,ij) = pt2d(ji,jj)
-            END DO
-         END DO
-
-         DO jr = 1,nsndto
-            IF ((nfipproc(isendto(jr),jpnj) .ne. (narea-1)) .and. (nfipproc(isendto(jr),jpnj) .ne. -1)) THEN
-               CALL mppsend(5, znorthloc, itaille, nfipproc(isendto(jr),jpnj), ml_req_nf(jr))
-            ENDIF
-         END DO
-         DO jr = 1,nsndto
-            iproc = nfipproc(isendto(jr),jpnj)
-            IF(iproc .ne. -1) THEN
-               ilei = nleit (iproc+1)
-               ildi = nldit (iproc+1)
-               iilb = nfiimpp(isendto(jr),jpnj) - nfiimpp(isendto(1),jpnj)
-            ENDIF
-            IF((iproc .ne. (narea-1)) .and. (iproc .ne. -1)) THEN
-              CALL mpprecv(5, zfoldwk, itaille, iproc)
-              DO jj = 1, ijpj
-                 DO ji = ildi, ilei
-                    ztabr(iilb+ji,jj) = zfoldwk(ji,jj)
-                 END DO
-              END DO
-            ELSE IF (iproc .eq. (narea-1)) THEN
-              DO jj = 1, ijpj
-                 DO ji = ildi, ilei
-                    ztabr(iilb+ji,jj) = pt2d(ji,nlcj-ijpj+jj)
-                 END DO
-              END DO
-            ENDIF
-         END DO
-         IF (l_isend) THEN
-            DO jr = 1,nsndto
-               IF ((nfipproc(isendto(jr),jpnj) .ne. (narea-1)) .and. (nfipproc(isendto(jr),jpnj) .ne. -1)) THEN
-                  CALL mpi_wait(ml_req_nf(jr), ml_stat, ml_err)
-               ENDIF
-            END DO
-         ENDIF
-         CALL mpp_lbc_nfd( ztabl, ztabr, cd_type, psgn )   ! North fold boundary condition
-         !
-         DO jj = nlcj-ijpj+1, nlcj             ! Scatter back to pt2d
-            ij = jj - nlcj + ijpj
-            DO ji = 1, nlci
-               pt2d(ji,jj) = ztabl(ji,ij)
-            END DO
-         END DO
-         !
-      ELSE
-         CALL MPI_ALLGATHER( znorthloc  , itaille, MPI_DOUBLE_PRECISION,        &
-            &                znorthgloio, itaille, MPI_DOUBLE_PRECISION, ncomm_north, ierr )
-         !
-         ztab(:,:) = 0.e0
-         DO jr = 1, ndim_rank_north            ! recover the global north array
-            iproc = nrank_north(jr) + 1
-            ildi = nldit (iproc)
-            ilei = nleit (iproc)
-            iilb = nimppt(iproc)
-            DO jj = 1, ijpj
-               DO ji = ildi, ilei
-                  ztab(ji+iilb-1,jj) = znorthgloio(ji,jj,jr)
-               END DO
-            END DO
-         END DO
-         CALL lbc_nfd( ztab, cd_type, psgn )   ! North fold boundary condition
-         !
-         DO jj = nlcj-ijpj+1, nlcj             ! Scatter back to pt2d
-            ij = jj - nlcj + ijpj
-            DO ji = 1, nlci
-               pt2d(ji,jj) = ztab(ji+nimpp-1,ij)
-            END DO
-         END DO
-         !
-      ENDIF
-      DEALLOCATE( ztab, znorthloc, zfoldwk, znorthgloio )
-      DEALLOCATE( ztabl, ztabr ) 
-      !
-   END SUBROUTINE mpp_lbc_north_2d
-
-   SUBROUTINE mpp_lbc_north_2d_multiple( pt2d_array, cd_type, psgn, num_fields)
-      !!---------------------------------------------------------------------
-      !!                   ***  routine mpp_lbc_north_2d  ***
-      !!
-      !! ** Purpose :   Ensure proper north fold horizontal bondary condition
-      !!              in mpp configuration in case of jpn1 > 1
-      !!              (for multiple 2d arrays )
-      !!
-      !! ** Method  :   North fold condition and mpp with more than one proc
-      !!              in i-direction require a specific treatment. We gather
-      !!              the 4 northern lines of the global domain on 1 processor
-      !!              and apply lbc north-fold on this sub array. Then we
-      !!              scatter the north fold array back to the processors.
-      !!
-      !!----------------------------------------------------------------------
-      INTEGER ,  INTENT (in   ) ::   num_fields  ! number of variables contained in pt2d
-      TYPE( arrayptr ), DIMENSION(:) :: pt2d_array
-      CHARACTER(len=1), DIMENSION(:), INTENT(in   ) ::   cd_type   ! nature of pt2d grid-points
-      !                                                          !   = T ,  U , V , F or W  gridpoints
-      REAL(wp), DIMENSION(:), INTENT(in   ) ::   psgn      ! = -1. the sign change across the north fold 
-      !!                                                             ! =  1. , the sign is kept
-      INTEGER ::   ji, jj, jr, jk
-      INTEGER ::   ierr, itaille, ildi, ilei, iilb
-      INTEGER ::   ijpj, ijpjm1, ij, iproc
-      INTEGER, DIMENSION (jpmaxngh)      ::   ml_req_nf          !for mpi_isend when avoiding mpi_allgather
-      INTEGER                            ::   ml_err             ! for mpi_isend when avoiding mpi_allgather
-      INTEGER, DIMENSION(MPI_STATUS_SIZE)::   ml_stat            ! for mpi_isend when avoiding mpi_allgather
-      !                                                              ! Workspace for message transfers avoiding mpi_allgather
-      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE   :: ztab
-      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE   :: znorthloc, zfoldwk
-      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE   :: znorthgloio
-      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE   :: ztabl, ztabr
-      INTEGER :: istatus(mpi_status_size)
-      INTEGER :: iflag
-      !!----------------------------------------------------------------------
-      !
-      ALLOCATE( ztab(jpiglo,4,num_fields), znorthloc(jpi,4,num_fields), zfoldwk(jpi,4,num_fields), znorthgloio(jpi,4,num_fields,jpni) )   ! expanded to 3 dimensions
-      ALLOCATE( ztabl(jpi,4,num_fields), ztabr(jpi*jpmaxngh, 4,num_fields) )
-      !
-      ijpj   = 4
-      ijpjm1 = 3
-      !
-      
-      DO jk = 1, num_fields
-         DO jj = nlcj-ijpj+1, nlcj             ! put in znorthloc the last 4 jlines of pt2d (for every variable)
-            ij = jj - nlcj + ijpj
-            znorthloc(:,ij,jk) = pt2d_array(jk)%pt2d(:,jj)
-         END DO
-      END DO
-      !                                     ! Build in procs of ncomm_north the znorthgloio
-      itaille = jpi * ijpj
-                                                                  
-      IF ( l_north_nogather ) THEN
-         !
-         ! Avoid the use of mpi_allgather by exchanging only with the processes already identified 
-         ! (in nemo_northcomms) as being  involved in this process' northern boundary exchange
-         !
-         ztabr(:,:,:) = 0
-         ztabl(:,:,:) = 0
-
-         DO jk = 1, num_fields
-            DO jj = nlcj-ijpj+1, nlcj          ! First put local values into the global array
-               ij = jj - nlcj + ijpj
-               DO ji = nfsloop, nfeloop
-                  ztabl(ji,ij,jk) = pt2d_array(jk)%pt2d(ji,jj)
-               END DO
-            END DO
-         END DO
-
-         DO jr = 1,nsndto
-            IF ((nfipproc(isendto(jr),jpnj) .ne. (narea-1)) .and. (nfipproc(isendto(jr),jpnj) .ne. -1)) THEN
-               CALL mppsend(5, znorthloc, itaille*num_fields, nfipproc(isendto(jr),jpnj), ml_req_nf(jr)) ! Buffer expanded "num_fields" times
-            ENDIF
-         END DO
-         DO jr = 1,nsndto
-            iproc = nfipproc(isendto(jr),jpnj)
-            IF(iproc .ne. -1) THEN
-               ilei = nleit (iproc+1)
-               ildi = nldit (iproc+1)
-               iilb = nfiimpp(isendto(jr),jpnj) - nfiimpp(isendto(1),jpnj)
-            ENDIF
-            IF((iproc .ne. (narea-1)) .and. (iproc .ne. -1)) THEN
-              CALL mpprecv(5, zfoldwk, itaille*num_fields, iproc) ! Buffer expanded "num_fields" times
-              DO jk = 1 , num_fields
-                 DO jj = 1, ijpj
-                    DO ji = ildi, ilei
-                       ztabr(iilb+ji,jj,jk) = zfoldwk(ji,jj,jk)       ! Modified to 3D
-                    END DO
-                 END DO
-              END DO
-            ELSE IF (iproc .eq. (narea-1)) THEN
-              DO jk = 1, num_fields
-                 DO jj = 1, ijpj
-                    DO ji = ildi, ilei
-                          ztabr(iilb+ji,jj,jk) = pt2d_array(jk)%pt2d(ji,nlcj-ijpj+jj)       ! Modified to 3D
-                    END DO
-                 END DO
-              END DO
-            ENDIF
-         END DO
-         IF (l_isend) THEN
-            DO jr = 1,nsndto
-               IF ((nfipproc(isendto(jr),jpnj) .ne. (narea-1)) .and. (nfipproc(isendto(jr),jpnj) .ne. -1)) THEN
-                  CALL mpi_wait(ml_req_nf(jr), ml_stat, ml_err)
-               ENDIF
-            END DO
-         ENDIF
-         !
-         DO ji = 1, num_fields     ! Loop to manage 3D variables
-            CALL mpp_lbc_nfd( ztabl(:,:,ji), ztabr(:,:,ji), cd_type(ji), psgn(ji) )  ! North fold boundary condition
-         END DO
-         !
-         DO jk = 1, num_fields
-            DO jj = nlcj-ijpj+1, nlcj             ! Scatter back to pt2d
-               ij = jj - nlcj + ijpj
-               DO ji = 1, nlci
-                  pt2d_array(jk)%pt2d(ji,jj) = ztabl(ji,ij,jk)       ! Modified to 3D
-               END DO
-            END DO
-         END DO
-         
-         !
-      ELSE
-         !
-         CALL MPI_ALLGATHER( znorthloc  , itaille*num_fields, MPI_DOUBLE_PRECISION,        &
-            &                znorthgloio, itaille*num_fields, MPI_DOUBLE_PRECISION, ncomm_north, ierr )
-         !
-         ztab(:,:,:) = 0.e0
-         DO jk = 1, num_fields
-            DO jr = 1, ndim_rank_north            ! recover the global north array
-               iproc = nrank_north(jr) + 1
-               ildi = nldit (iproc)
-               ilei = nleit (iproc)
-               iilb = nimppt(iproc)
-               DO jj = 1, ijpj
-                  DO ji = ildi, ilei
-                     ztab(ji+iilb-1,jj,jk) = znorthgloio(ji,jj,jk,jr)
-                  END DO
-               END DO
-            END DO
-         END DO
-         
-         DO ji = 1, num_fields
-            CALL lbc_nfd( ztab(:,:,ji), cd_type(ji), psgn(ji) )   ! North fold boundary condition
-         END DO
-         !
-         DO jk = 1, num_fields
-            DO jj = nlcj-ijpj+1, nlcj             ! Scatter back to pt2d
-               ij = jj - nlcj + ijpj
-               DO ji = 1, nlci
-                  pt2d_array(jk)%pt2d(ji,jj) = ztab(ji+nimpp-1,ij,jk)
-               END DO
-            END DO
-         END DO
-         !
-         !
-      ENDIF
-      DEALLOCATE( ztab, znorthloc, zfoldwk, znorthgloio )
-      DEALLOCATE( ztabl, ztabr )
-      !
-   END SUBROUTINE mpp_lbc_north_2d_multiple
-
-   SUBROUTINE mpp_lbc_north_e( pt2d, cd_type, psgn)
-      !!---------------------------------------------------------------------
-      !!                   ***  routine mpp_lbc_north_2d  ***
-      !!
-      !! ** Purpose :   Ensure proper north fold horizontal bondary condition
-      !!              in mpp configuration in case of jpn1 > 1 and for 2d
-      !!              array with outer extra halo
-      !!
-      !! ** Method  :   North fold condition and mpp with more than one proc
-      !!              in i-direction require a specific treatment. We gather
-      !!              the 4+2*jpr2dj northern lines of the global domain on 1
-      !!              processor and apply lbc north-fold on this sub array.
-      !!              Then we scatter the north fold array back to the processors.
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj), INTENT(inout) ::   pt2d     ! 2D array with extra halo
-      CHARACTER(len=1)                                            , INTENT(in   ) ::   cd_type  ! nature of pt3d grid-points
-      !                                                                                         !   = T ,  U , V , F or W -points
-      REAL(wp)                                                    , INTENT(in   ) ::   psgn     ! = -1. the sign change across the
-      !!                                                                                        ! north fold, =  1. otherwise
-      INTEGER ::   ji, jj, jr
-      INTEGER ::   ierr, itaille, ildi, ilei, iilb
-      INTEGER ::   ijpj, ij, iproc
-      !
-      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE  ::  ztab_e, znorthloc_e
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE  ::  znorthgloio_e
-
-      !!----------------------------------------------------------------------
-      !
-      ALLOCATE( ztab_e(jpiglo,4+2*jpr2dj), znorthloc_e(jpi,4+2*jpr2dj), znorthgloio_e(jpi,4+2*jpr2dj,jpni) )
-
-      !
-      ijpj=4
-      ztab_e(:,:) = 0.e0
-
-      ij=0
-      ! put in znorthloc_e the last 4 jlines of pt2d
-      DO jj = nlcj - ijpj + 1 - jpr2dj, nlcj +jpr2dj
-         ij = ij + 1
-         DO ji = 1, jpi
-            znorthloc_e(ji,ij)=pt2d(ji,jj)
-         END DO
-      END DO
-      !
-      itaille = jpi * ( ijpj + 2 * jpr2dj )
-      CALL MPI_ALLGATHER( znorthloc_e(1,1)  , itaille, MPI_DOUBLE_PRECISION,    &
-         &                znorthgloio_e(1,1,1), itaille, MPI_DOUBLE_PRECISION, ncomm_north, ierr )
-      !
-      DO jr = 1, ndim_rank_north            ! recover the global north array
-         iproc = nrank_north(jr) + 1
-         ildi = nldit (iproc)
-         ilei = nleit (iproc)
-         iilb = nimppt(iproc)
-         DO jj = 1, ijpj+2*jpr2dj
-            DO ji = ildi, ilei
-               ztab_e(ji+iilb-1,jj) = znorthgloio_e(ji,jj,jr)
-            END DO
-         END DO
-      END DO
-
-
-      ! 2. North-Fold boundary conditions
-      ! ----------------------------------
-      CALL lbc_nfd( ztab_e(:,:), cd_type, psgn, pr2dj = jpr2dj )
-
-      ij = jpr2dj
-      !! Scatter back to pt2d
-      DO jj = nlcj - ijpj + 1 , nlcj +jpr2dj
-      ij  = ij +1
-         DO ji= 1, nlci
-            pt2d(ji,jj) = ztab_e(ji+nimpp-1,ij)
-         END DO
-      END DO
-      !
-      DEALLOCATE( ztab_e, znorthloc_e, znorthgloio_e )
-      !
-   END SUBROUTINE mpp_lbc_north_e
-
-
-   SUBROUTINE mpp_lnk_bdy_3d( ptab, cd_type, psgn, ib_bdy )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mpp_lnk_bdy_3d  ***
-      !!
-      !! ** Purpose :   Message passing management
-      !!
-      !! ** Method  :   Use mppsend and mpprecv function for passing BDY boundaries 
-      !!      between processors following neighboring subdomains.
-      !!            domain parameters
-      !!                    nlci   : first dimension of the local subdomain
-      !!                    nlcj   : second dimension of the local subdomain
-      !!                    nbondi_bdy : mark for "east-west local boundary"
-      !!                    nbondj_bdy : mark for "north-south local boundary"
-      !!                    noea   : number for local neighboring processors 
-      !!                    nowe   : number for local neighboring processors
-      !!                    noso   : number for local neighboring processors
-      !!                    nono   : number for local neighboring processors
-      !!
-      !! ** Action  :   ptab with update value at its periphery
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptab     ! 3D array on which the boundary condition is applied
-      CHARACTER(len=1)                , INTENT(in   ) ::   cd_type  ! define the nature of ptab array grid-points
-      !                                                             ! = T , U , V , F , W points
-      REAL(wp)                        , INTENT(in   ) ::   psgn     ! =-1 the sign change across the north fold boundary
-      !                                                             ! =  1. , the sign is kept
-      INTEGER                         , INTENT(in   ) ::   ib_bdy   ! BDY boundary set
-      !
-      INTEGER  ::   ji, jj, jk, jl             ! dummy loop indices
-      INTEGER  ::   imigr, iihom, ijhom        ! local integers
-      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
-      REAL(wp) ::   zland                      ! local scalar
-      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
-      !
-      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::   zt3ns, zt3sn   ! 3d for north-south & south-north
-      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::   zt3ew, zt3we   ! 3d for east-west & west-east
-      !!----------------------------------------------------------------------
-      !
-      ALLOCATE( zt3ns(jpi,jprecj,jpk,2), zt3sn(jpi,jprecj,jpk,2),   &
-         &      zt3ew(jpj,jpreci,jpk,2), zt3we(jpj,jpreci,jpk,2)  )
-
-      zland = 0._wp
-
-      ! 1. standard boundary treatment
-      ! ------------------------------
-      !                                   ! East-West boundaries
-      !                                        !* Cyclic east-west
-      IF( nbondi == 2) THEN
-         IF( nperio == 1 .OR. nperio == 4 .OR. nperio == 6 ) THEN
-            ptab( 1 ,:,:) = ptab(jpim1,:,:)
-            ptab(jpi,:,:) = ptab(  2  ,:,:)
-         ELSE
-            IF( .NOT. cd_type == 'F' )   ptab(1:jpreci,:,:) = zland    ! south except F-point
-            ptab(nlci-jpreci+1:jpi,:,:) = zland    ! north
-         ENDIF
-      ELSEIF(nbondi == -1) THEN
-         IF( .NOT. cd_type == 'F' )   ptab(1:jpreci,:,:) = zland    ! south except F-point
-      ELSEIF(nbondi == 1) THEN
-         ptab(nlci-jpreci+1:jpi,:,:) = zland    ! north
-      ENDIF                                     !* closed
-
-      IF (nbondj == 2 .OR. nbondj == -1) THEN
-        IF( .NOT. cd_type == 'F' )   ptab(:,1:jprecj,:) = zland       ! south except F-point
-      ELSEIF (nbondj == 2 .OR. nbondj == 1) THEN
-        ptab(:,nlcj-jprecj+1:jpj,:) = zland       ! north
-      ENDIF
-      !
-      ! 2. East and west directions exchange
-      ! ------------------------------------
-      ! we play with the neigbours AND the row number because of the periodicity 
-      !
-      SELECT CASE ( nbondi_bdy(ib_bdy) )      ! Read Dirichlet lateral conditions
-      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
-         iihom = nlci-nreci
-         DO jl = 1, jpreci
-            zt3ew(:,jl,:,1) = ptab(jpreci+jl,:,:)
-            zt3we(:,jl,:,1) = ptab(iihom +jl,:,:)
-         END DO
-      END SELECT
-      !
-      !                           ! Migrations
-      imigr = jpreci * jpj * jpk
-      !
-      SELECT CASE ( nbondi_bdy(ib_bdy) )
-      CASE ( -1 )
-         CALL mppsend( 2, zt3we(1,1,1,1), imigr, noea, ml_req1 )
-      CASE ( 0 )
-         CALL mppsend( 1, zt3ew(1,1,1,1), imigr, nowe, ml_req1 )
-         CALL mppsend( 2, zt3we(1,1,1,1), imigr, noea, ml_req2 )
-      CASE ( 1 )
-         CALL mppsend( 1, zt3ew(1,1,1,1), imigr, nowe, ml_req1 )
-      END SELECT
-      !
-      SELECT CASE ( nbondi_bdy_b(ib_bdy) )
-      CASE ( -1 )
-         CALL mpprecv( 1, zt3ew(1,1,1,2), imigr, noea )
-      CASE ( 0 )
-         CALL mpprecv( 1, zt3ew(1,1,1,2), imigr, noea )
-         CALL mpprecv( 2, zt3we(1,1,1,2), imigr, nowe )
-      CASE ( 1 )
-         CALL mpprecv( 2, zt3we(1,1,1,2), imigr, nowe )
-      END SELECT
-      !
-      SELECT CASE ( nbondi_bdy(ib_bdy) )
-      CASE ( -1 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      CASE ( 0 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
-      CASE ( 1 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      iihom = nlci-jpreci
-      !
-      SELECT CASE ( nbondi_bdy_b(ib_bdy) )
-      CASE ( -1 )
-         DO jl = 1, jpreci
-            ptab(iihom+jl,:,:) = zt3ew(:,jl,:,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jpreci
-            ptab(      jl,:,:) = zt3we(:,jl,:,2)
-            ptab(iihom+jl,:,:) = zt3ew(:,jl,:,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jpreci
-            ptab(      jl,:,:) = zt3we(:,jl,:,2)
-         END DO
-      END SELECT
-
-
-      ! 3. North and south directions
-      ! -----------------------------
-      ! always closed : we play only with the neigbours
-      !
-      IF( nbondj_bdy(ib_bdy) /= 2 ) THEN      ! Read Dirichlet lateral conditions
-         ijhom = nlcj-nrecj
-         DO jl = 1, jprecj
-            zt3sn(:,jl,:,1) = ptab(:,ijhom +jl,:)
-            zt3ns(:,jl,:,1) = ptab(:,jprecj+jl,:)
-         END DO
-      ENDIF
-      !
-      !                           ! Migrations
-      imigr = jprecj * jpi * jpk
-      !
-      SELECT CASE ( nbondj_bdy(ib_bdy) )
-      CASE ( -1 )
-         CALL mppsend( 4, zt3sn(1,1,1,1), imigr, nono, ml_req1 )
-      CASE ( 0 )
-         CALL mppsend( 3, zt3ns(1,1,1,1), imigr, noso, ml_req1 )
-         CALL mppsend( 4, zt3sn(1,1,1,1), imigr, nono, ml_req2 )
-      CASE ( 1 )
-         CALL mppsend( 3, zt3ns(1,1,1,1), imigr, noso, ml_req1 )
-      END SELECT
-      !
-      SELECT CASE ( nbondj_bdy_b(ib_bdy) )
-      CASE ( -1 )
-         CALL mpprecv( 3, zt3ns(1,1,1,2), imigr, nono )
-      CASE ( 0 )
-         CALL mpprecv( 3, zt3ns(1,1,1,2), imigr, nono )
-         CALL mpprecv( 4, zt3sn(1,1,1,2), imigr, noso )
-      CASE ( 1 )
-         CALL mpprecv( 4, zt3sn(1,1,1,2), imigr, noso )
-      END SELECT
-      !
-      SELECT CASE ( nbondj_bdy(ib_bdy) )
-      CASE ( -1 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      CASE ( 0 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
-      CASE ( 1 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      ijhom = nlcj-jprecj
-      !
-      SELECT CASE ( nbondj_bdy_b(ib_bdy) )
-      CASE ( -1 )
-         DO jl = 1, jprecj
-            ptab(:,ijhom+jl,:) = zt3ns(:,jl,:,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jprecj
-            ptab(:,jl      ,:) = zt3sn(:,jl,:,2)
-            ptab(:,ijhom+jl,:) = zt3ns(:,jl,:,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jprecj
-            ptab(:,jl,:) = zt3sn(:,jl,:,2)
-         END DO
-      END SELECT
-
-
-      ! 4. north fold treatment
-      ! -----------------------
-      !
-      IF( npolj /= 0) THEN
-         !
-         SELECT CASE ( jpni )
-         CASE ( 1 )     ;   CALL lbc_nfd      ( ptab, cd_type, psgn )   ! only 1 northern proc, no mpp
-         CASE DEFAULT   ;   CALL mpp_lbc_north( ptab, cd_type, psgn )   ! for all northern procs.
-         END SELECT
-         !
-      ENDIF
-      !
-      DEALLOCATE( zt3ns, zt3sn, zt3ew, zt3we  )
-      !
-   END SUBROUTINE mpp_lnk_bdy_3d
-
-
-   SUBROUTINE mpp_lnk_bdy_2d( ptab, cd_type, psgn, ib_bdy )
-      !!----------------------------------------------------------------------
-      !!                  ***  routine mpp_lnk_bdy_2d  ***
-      !!
-      !! ** Purpose :   Message passing management
-      !!
-      !! ** Method  :   Use mppsend and mpprecv function for passing BDY boundaries 
-      !!      between processors following neighboring subdomains.
-      !!            domain parameters
-      !!                    nlci   : first dimension of the local subdomain
-      !!                    nlcj   : second dimension of the local subdomain
-      !!                    nbondi_bdy : mark for "east-west local boundary"
-      !!                    nbondj_bdy : mark for "north-south local boundary"
-      !!                    noea   : number for local neighboring processors 
-      !!                    nowe   : number for local neighboring processors
-      !!                    noso   : number for local neighboring processors
-      !!                    nono   : number for local neighboring processors
-      !!
-      !! ** Action  :   ptab with update value at its periphery
-      !!
-      !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   ptab     ! 3D array on which the boundary condition is applied
-      CHARACTER(len=1)            , INTENT(in   ) ::   cd_type  ! define the nature of ptab array grid-points
-      !                                                         ! = T , U , V , F , W points
-      REAL(wp)                    , INTENT(in   ) ::   psgn     ! =-1 the sign change across the north fold boundary
-      !                                                         ! =  1. , the sign is kept
-      INTEGER                     , INTENT(in   ) ::   ib_bdy   ! BDY boundary set
-      !
-      INTEGER  ::   ji, jj, jl             ! dummy loop indices
-      INTEGER  ::   imigr, iihom, ijhom        ! local integers
-      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
-      REAL(wp) ::   zland
-      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
-      !
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  zt2ns, zt2sn   ! 2d for north-south & south-north
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  zt2ew, zt2we   ! 2d for east-west & west-east
-      !!----------------------------------------------------------------------
-
-      ALLOCATE( zt2ns(jpi,jprecj,2), zt2sn(jpi,jprecj,2),  &
-         &      zt2ew(jpj,jpreci,2), zt2we(jpj,jpreci,2)   )
-
-      zland = 0._wp
-
-      ! 1. standard boundary treatment
-      ! ------------------------------
-      !                                   ! East-West boundaries
-      !                                      !* Cyclic east-west
-      IF( nbondi == 2 ) THEN
-         IF (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) THEN
-            ptab( 1 ,:) = ptab(jpim1,:)
-            ptab(jpi,:) = ptab(  2  ,:)
-         ELSE
-            IF(.NOT.cd_type == 'F' )  ptab(     1       :jpreci,:) = zland    ! south except F-point
-                                      ptab(nlci-jpreci+1:jpi   ,:) = zland    ! north
-         ENDIF
-      ELSEIF(nbondi == -1) THEN
-         IF( .NOT.cd_type == 'F' )    ptab(     1       :jpreci,:) = zland    ! south except F-point
-      ELSEIF(nbondi == 1) THEN
-                                      ptab(nlci-jpreci+1:jpi   ,:) = zland    ! north
-      ENDIF
-      !                                      !* closed
-      IF( nbondj == 2 .OR. nbondj == -1 ) THEN
-         IF( .NOT.cd_type == 'F' )    ptab(:,     1       :jprecj) = zland    ! south except F-point
-      ELSEIF (nbondj == 2 .OR. nbondj == 1) THEN
-                                      ptab(:,nlcj-jprecj+1:jpj   ) = zland    ! north
-      ENDIF
-      !
-      ! 2. East and west directions exchange
-      ! ------------------------------------
-      ! we play with the neigbours AND the row number because of the periodicity 
-      !
-      SELECT CASE ( nbondi_bdy(ib_bdy) )      ! Read Dirichlet lateral conditions
-      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
-         iihom = nlci-nreci
-         DO jl = 1, jpreci
-            zt2ew(:,jl,1) = ptab(jpreci+jl,:)
-            zt2we(:,jl,1) = ptab(iihom +jl,:)
-         END DO
-      END SELECT
-      !
-      !                           ! Migrations
-      imigr = jpreci * jpj
-      !
-      SELECT CASE ( nbondi_bdy(ib_bdy) )
-      CASE ( -1 )
-         CALL mppsend( 2, zt2we(1,1,1), imigr, noea, ml_req1 )
-      CASE ( 0 )
-         CALL mppsend( 1, zt2ew(1,1,1), imigr, nowe, ml_req1 )
-         CALL mppsend( 2, zt2we(1,1,1), imigr, noea, ml_req2 )
-      CASE ( 1 )
-         CALL mppsend( 1, zt2ew(1,1,1), imigr, nowe, ml_req1 )
-      END SELECT
-      !
-      SELECT CASE ( nbondi_bdy_b(ib_bdy) )
-      CASE ( -1 )
-         CALL mpprecv( 1, zt2ew(1,1,2), imigr, noea )
-      CASE ( 0 )
-         CALL mpprecv( 1, zt2ew(1,1,2), imigr, noea )
-         CALL mpprecv( 2, zt2we(1,1,2), imigr, nowe )
-      CASE ( 1 )
-         CALL mpprecv( 2, zt2we(1,1,2), imigr, nowe )
-      END SELECT
-      !
-      SELECT CASE ( nbondi_bdy(ib_bdy) )
-      CASE ( -1 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      CASE ( 0 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
-      CASE ( 1 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      iihom = nlci-jpreci
-      !
-      SELECT CASE ( nbondi_bdy_b(ib_bdy) )
-      CASE ( -1 )
-         DO jl = 1, jpreci
-            ptab(iihom+jl,:) = zt2ew(:,jl,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jpreci
-            ptab(jl      ,:) = zt2we(:,jl,2)
-            ptab(iihom+jl,:) = zt2ew(:,jl,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jpreci
-            ptab(jl      ,:) = zt2we(:,jl,2)
-         END DO
-      END SELECT
-
-
-      ! 3. North and south directions
-      ! -----------------------------
-      ! always closed : we play only with the neigbours
-      !
-      IF( nbondj_bdy(ib_bdy) /= 2 ) THEN      ! Read Dirichlet lateral conditions
-         ijhom = nlcj-nrecj
-         DO jl = 1, jprecj
-            zt2sn(:,jl,1) = ptab(:,ijhom +jl)
-            zt2ns(:,jl,1) = ptab(:,jprecj+jl)
-         END DO
-      ENDIF
-      !
-      !                           ! Migrations
-      imigr = jprecj * jpi
-      !
-      SELECT CASE ( nbondj_bdy(ib_bdy) )
-      CASE ( -1 )
-         CALL mppsend( 4, zt2sn(1,1,1), imigr, nono, ml_req1 )
-      CASE ( 0 )
-         CALL mppsend( 3, zt2ns(1,1,1), imigr, noso, ml_req1 )
-         CALL mppsend( 4, zt2sn(1,1,1), imigr, nono, ml_req2 )
-      CASE ( 1 )
-         CALL mppsend( 3, zt2ns(1,1,1), imigr, noso, ml_req1 )
-      END SELECT
-      !
-      SELECT CASE ( nbondj_bdy_b(ib_bdy) )
-      CASE ( -1 )
-         CALL mpprecv( 3, zt2ns(1,1,2), imigr, nono )
-      CASE ( 0 )
-         CALL mpprecv( 3, zt2ns(1,1,2), imigr, nono )
-         CALL mpprecv( 4, zt2sn(1,1,2), imigr, noso )
-      CASE ( 1 )
-         CALL mpprecv( 4, zt2sn(1,1,2), imigr, noso )
-      END SELECT
-      !
-      SELECT CASE ( nbondj_bdy(ib_bdy) )
-      CASE ( -1 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      CASE ( 0 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
-      CASE ( 1 )
-         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
-      END SELECT
-      !
-      !                           ! Write Dirichlet lateral conditions
-      ijhom = nlcj-jprecj
-      !
-      SELECT CASE ( nbondj_bdy_b(ib_bdy) )
-      CASE ( -1 )
-         DO jl = 1, jprecj
-            ptab(:,ijhom+jl) = zt2ns(:,jl,2)
-         END DO
-      CASE ( 0 )
-         DO jl = 1, jprecj
-            ptab(:,jl      ) = zt2sn(:,jl,2)
-            ptab(:,ijhom+jl) = zt2ns(:,jl,2)
-         END DO
-      CASE ( 1 )
-         DO jl = 1, jprecj
-            ptab(:,jl) = zt2sn(:,jl,2)
-         END DO
-      END SELECT
-
-
-      ! 4. north fold treatment
-      ! -----------------------
-      !
-      IF( npolj /= 0) THEN
-         !
-         SELECT CASE ( jpni )
-         CASE ( 1 )     ;   CALL lbc_nfd      ( ptab, cd_type, psgn )   ! only 1 northern proc, no mpp
-         CASE DEFAULT   ;   CALL mpp_lbc_north( ptab, cd_type, psgn )   ! for all northern procs.
-         END SELECT
-         !
-      ENDIF
-      !
-      DEALLOCATE( zt2ns, zt2sn, zt2ew, zt2we  )
-      !
-   END SUBROUTINE mpp_lnk_bdy_2d
-
-
-   SUBROUTINE mpi_init_opa( ldtxt, ksft, code )
+   SUBROUTINE mpi_init_oce( ldtxt, ksft, code )
       !!---------------------------------------------------------------------
       !!                   ***  routine mpp_init.opa  ***
@@ -3649,5 +1127,5 @@
       IF( .NOT. mpi_was_called ) THEN
          CALL mpi_init( code )
-         CALL mpi_comm_dup( mpi_comm_world, mpi_comm_opa, code )
+         CALL mpi_comm_dup( mpi_comm_world, mpi_comm_oce, code )
          IF ( code /= MPI_SUCCESS ) THEN
             DO ji = 1, SIZE(ldtxt)
@@ -3675,7 +1153,8 @@
       ENDIF
       !
-   END SUBROUTINE mpi_init_opa
-
-   SUBROUTINE DDPDD_MPI (ydda, yddb, ilen, itype)
+   END SUBROUTINE mpi_init_oce
+
+
+   SUBROUTINE DDPDD_MPI( ydda, yddb, ilen, itype )
       !!---------------------------------------------------------------------
       !!   Routine DDPDD_MPI: used by reduction operator MPI_SUMDD
@@ -3684,13 +1163,14 @@
       !!   This subroutine computes yddb(i) = ydda(i)+yddb(i)
       !!---------------------------------------------------------------------
-      INTEGER, INTENT(in)                         :: ilen, itype
-      COMPLEX(wp), DIMENSION(ilen), INTENT(in)     :: ydda
-      COMPLEX(wp), DIMENSION(ilen), INTENT(inout)  :: yddb
+      INTEGER                     , INTENT(in)    ::   ilen, itype
+      COMPLEX(wp), DIMENSION(ilen), INTENT(in)    ::   ydda
+      COMPLEX(wp), DIMENSION(ilen), INTENT(inout) ::   yddb
       !
       REAL(wp) :: zerr, zt1, zt2    ! local work variables
-      INTEGER :: ji, ztmp           ! local scalar
-
+      INTEGER  :: ji, ztmp           ! local scalar
+      !!---------------------------------------------------------------------
+      !
       ztmp = itype   ! avoid compilation warning
-
+      !
       DO ji=1,ilen
       ! Compute ydda + yddb using Knuth's trick.
@@ -3703,9 +1183,9 @@
          yddb(ji) = cmplx ( zt1 + zt2, zt2 - ((zt1 + zt2) - zt1),wp )
       END DO
-
+      !
    END SUBROUTINE DDPDD_MPI
 
 
-   SUBROUTINE mpp_lbc_north_icb( pt2d, cd_type, psgn, pr2dj)
+   SUBROUTINE mpp_lbc_north_icb( pt2d, cd_type, psgn, kextj)
       !!---------------------------------------------------------------------
       !!                   ***  routine mpp_lbc_north_icb  ***
@@ -3717,8 +1197,9 @@
       !! ** Method  :   North fold condition and mpp with more than one proc
       !!              in i-direction require a specific treatment. We gather
-      !!              the 4+2*jpr2dj northern lines of the global domain on 1
+      !!              the 4+kextj northern lines of the global domain on 1
       !!              processor and apply lbc north-fold on this sub array.
       !!              Then we scatter the north fold array back to the processors.
-      !!              This version accounts for an extra halo with icebergs.
+      !!              This routine accounts for an extra halo with icebergs
+      !!              and assumes ghost rows and columns have been suppressed.
       !!
       !!----------------------------------------------------------------------
@@ -3728,9 +1209,9 @@
       REAL(wp)                , INTENT(in   ) ::   psgn     ! = -1. the sign change across the
       !!                                                    ! north fold, =  1. otherwise
-      INTEGER, OPTIONAL       , INTENT(in   ) ::   pr2dj
+      INTEGER                 , INTENT(in   ) ::   kextj    ! Extra halo width at north fold
       !
       INTEGER ::   ji, jj, jr
       INTEGER ::   ierr, itaille, ildi, ilei, iilb
-      INTEGER ::   ijpj, ij, iproc, ipr2dj
+      INTEGER ::   ipj, ij, iproc
       !
       REAL(wp), DIMENSION(:,:)  , ALLOCATABLE  ::  ztab_e, znorthloc_e
@@ -3738,26 +1219,27 @@
       !!----------------------------------------------------------------------
       !
-      ijpj=4
-      IF( PRESENT(pr2dj) ) THEN           ! use of additional halos
-         ipr2dj = pr2dj
-      ELSE
-         ipr2dj = 0
-      ENDIF
-      ALLOCATE( ztab_e(jpiglo,4+2*ipr2dj), znorthloc_e(jpi,4+2*ipr2dj), znorthgloio_e(jpi,4+2*ipr2dj,jpni) )
-      !
-      ztab_e(:,:) = 0._wp
-      !
-      ij = 0
-      ! put in znorthloc_e the last 4 jlines of pt2d
-      DO jj = nlcj - ijpj + 1 - ipr2dj, nlcj +ipr2dj
+      ipj=4
+      ALLOCATE(        ztab_e(jpiglo, 1-kextj:ipj+kextj)       ,       &
+     &            znorthloc_e(jpimax, 1-kextj:ipj+kextj)       ,       &
+     &          znorthgloio_e(jpimax, 1-kextj:ipj+kextj,jpni)    )
+      !
+      ztab_e(:,:)      = 0._wp
+      znorthloc_e(:,:) = 0._wp
+      !
+      ij = 1 - kextj
+      ! put the last ipj+2*kextj lines of pt2d into znorthloc_e 
+      DO jj = jpj - ipj + 1 - kextj , jpj + kextj
+         znorthloc_e(1:jpi,ij)=pt2d(1:jpi,jj)
          ij = ij + 1
-         DO ji = 1, jpi
-            znorthloc_e(ji,ij)=pt2d(ji,jj)
-         END DO
       END DO
       !
-      itaille = jpi * ( ijpj + 2 * ipr2dj )
-      CALL MPI_ALLGATHER( znorthloc_e(1,1)  , itaille, MPI_DOUBLE_PRECISION,    &
-         &                znorthgloio_e(1,1,1), itaille, MPI_DOUBLE_PRECISION, ncomm_north, ierr )
+      itaille = jpimax * ( ipj + 2*kextj )
+      !
+      IF( ln_timing ) CALL tic_tac(.TRUE.)
+      CALL MPI_ALLGATHER( znorthloc_e(1,1-kextj)    , itaille, MPI_DOUBLE_PRECISION,    &
+         &                znorthgloio_e(1,1-kextj,1), itaille, MPI_DOUBLE_PRECISION,    &
+         &                ncomm_north, ierr )
+      !
+      IF( ln_timing ) CALL tic_tac(.FALSE.)
       !
       DO jr = 1, ndim_rank_north            ! recover the global north array
@@ -3766,5 +1248,5 @@
          ilei = nleit (iproc)
          iilb = nimppt(iproc)
-         DO jj = 1, ijpj+2*ipr2dj
+         DO jj = 1-kextj, ipj+kextj
             DO ji = ildi, ilei
                ztab_e(ji+iilb-1,jj) = znorthgloio_e(ji,jj,jr)
@@ -3773,16 +1255,15 @@
       END DO
 
-
       ! 2. North-Fold boundary conditions
       ! ----------------------------------
-      CALL lbc_nfd( ztab_e(:,:), cd_type, psgn, pr2dj = ipr2dj )
-
-      ij = ipr2dj
+      CALL lbc_nfd( ztab_e(:,1-kextj:ipj+kextj), cd_type, psgn, kextj )
+
+      ij = 1 - kextj
       !! Scatter back to pt2d
-      DO jj = nlcj - ijpj + 1 , nlcj +ipr2dj
-      ij  = ij +1
-         DO ji= 1, nlci
+      DO jj = jpj - ipj + 1 - kextj , jpj + kextj
+         DO ji= 1, jpi
             pt2d(ji,jj) = ztab_e(ji+nimpp-1,ij)
          END DO
+         ij  = ij +1
       END DO
       !
@@ -3792,17 +1273,19 @@
 
 
-   SUBROUTINE mpp_lnk_2d_icb( pt2d, cd_type, psgn, jpri, jprj )
+   SUBROUTINE mpp_lnk_2d_icb( cdname, pt2d, cd_type, psgn, kexti, kextj )
       !!----------------------------------------------------------------------
       !!                  ***  routine mpp_lnk_2d_icb  ***
       !!
-      !! ** Purpose :   Message passing manadgement for 2d array (with extra halo and icebergs)
+      !! ** Purpose :   Message passing management for 2d array (with extra halo for icebergs)
+      !!                This routine receives a (1-kexti:jpi+kexti,1-kexti:jpj+kextj)
+      !!                array (usually (0:jpi+1, 0:jpj+1)) from lbc_lnk_icb calls.
       !!
       !! ** Method  :   Use mppsend and mpprecv function for passing mask
       !!      between processors following neighboring subdomains.
       !!            domain parameters
-      !!                    nlci   : first dimension of the local subdomain
-      !!                    nlcj   : second dimension of the local subdomain
-      !!                    jpri   : number of rows for extra outer halo
-      !!                    jprj   : number of columns for extra outer halo
+      !!                    jpi    : first dimension of the local subdomain
+      !!                    jpj    : second dimension of the local subdomain
+      !!                    kexti  : number of columns for extra outer halo
+      !!                    kextj  : number of rows for extra outer halo
       !!                    nbondi : mark for "east-west local boundary"
       !!                    nbondj : mark for "north-south local boundary"
@@ -3812,26 +1295,25 @@
       !!                    nono   : number for local neighboring processors
       !!----------------------------------------------------------------------
-      INTEGER                                             , INTENT(in   ) ::   jpri
-      INTEGER                                             , INTENT(in   ) ::   jprj
-      REAL(wp), DIMENSION(1-jpri:jpi+jpri,1-jprj:jpj+jprj), INTENT(inout) ::   pt2d     ! 2D array with extra halo
-      CHARACTER(len=1)                                    , INTENT(in   ) ::   cd_type  ! nature of ptab array grid-points
-      !                                                                                 ! = T , U , V , F , W and I points
-      REAL(wp)                                            , INTENT(in   ) ::   psgn     ! =-1 the sign change across the
-      !!                                                                                ! north boundary, =  1. otherwise
+      CHARACTER(len=*)                                        , INTENT(in   ) ::   cdname      ! name of the calling subroutine
+      REAL(wp), DIMENSION(1-kexti:jpi+kexti,1-kextj:jpj+kextj), INTENT(inout) ::   pt2d     ! 2D array with extra halo
+      CHARACTER(len=1)                                        , INTENT(in   ) ::   cd_type  ! nature of ptab array grid-points
+      REAL(wp)                                                , INTENT(in   ) ::   psgn     ! sign used across the north fold
+      INTEGER                                                 , INTENT(in   ) ::   kexti    ! extra i-halo width
+      INTEGER                                                 , INTENT(in   ) ::   kextj    ! extra j-halo width
+      !
       INTEGER  ::   jl   ! dummy loop indices
-      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
-      INTEGER  ::   ipreci, iprecj             ! temporary integers
+      INTEGER  ::   imigr, iihom, ijhom        ! local integers
+      INTEGER  ::   ipreci, iprecj             !   -       -
       INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
       INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
       !!
-      REAL(wp), DIMENSION(1-jpri:jpi+jpri,jprecj+jprj,2) :: r2dns
-      REAL(wp), DIMENSION(1-jpri:jpi+jpri,jprecj+jprj,2) :: r2dsn
-      REAL(wp), DIMENSION(1-jprj:jpj+jprj,jpreci+jpri,2) :: r2dwe
-      REAL(wp), DIMENSION(1-jprj:jpj+jprj,jpreci+jpri,2) :: r2dew
-      !!----------------------------------------------------------------------
-
-      ipreci = jpreci + jpri      ! take into account outer extra 2D overlap area
-      iprecj = jprecj + jprj
-
+      REAL(wp), DIMENSION(1-kexti:jpi+kexti,nn_hls+kextj,2) ::   r2dns, r2dsn
+      REAL(wp), DIMENSION(1-kextj:jpj+kextj,nn_hls+kexti,2) ::   r2dwe, r2dew
+      !!----------------------------------------------------------------------
+
+      ipreci = nn_hls + kexti      ! take into account outer extra 2D overlap area
+      iprecj = nn_hls + kextj
+
+      IF( narea == 1 .AND. numcom == -1 ) CALL mpp_report( cdname, 1, 1, 1, ld_lbc = .TRUE. )
 
       ! 1. standard boundary treatment
@@ -3841,11 +1323,19 @@
       !                                      ! East-West boundaries
       !                                           !* Cyclic east-west
-      IF( nbondi == 2 .AND. (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
-         pt2d(1-jpri:     1    ,:) = pt2d(jpim1-jpri:  jpim1 ,:)       ! east
-         pt2d(   jpi  :jpi+jpri,:) = pt2d(     2      :2+jpri,:)       ! west
+      IF( l_Iperio ) THEN
+         pt2d(1-kexti:     1   ,:) = pt2d(jpim1-kexti: jpim1 ,:)       ! east
+         pt2d(  jpi  :jpi+kexti,:) = pt2d(     2     :2+kexti,:)       ! west
          !
       ELSE                                        !* closed
-         IF( .NOT. cd_type == 'F' )   pt2d(  1-jpri   :jpreci    ,:) = 0.e0    ! south except at F-point
-                                      pt2d(nlci-jpreci+1:jpi+jpri,:) = 0.e0    ! north
+         IF( .NOT. cd_type == 'F' )   pt2d(  1-kexti   :nn_hls   ,:) = 0._wp    ! east except at F-point
+                                      pt2d(jpi-nn_hls+1:jpi+kexti,:) = 0._wp    ! west
+      ENDIF
+      !                                      ! North-South boundaries
+      IF( l_Jperio ) THEN                         !* cyclic (only with no mpp j-split)
+         pt2d(:,1-kextj:     1   ) = pt2d(:,jpjm1-kextj:  jpjm1)       ! north
+         pt2d(:,  jpj  :jpj+kextj) = pt2d(:,     2     :2+kextj)       ! south
+      ELSE                                        !* closed
+         IF( .NOT. cd_type == 'F' )   pt2d(:,  1-kextj   :nn_hls   ) = 0._wp    ! north except at F-point
+                                      pt2d(:,jpj-nn_hls+1:jpj+kextj) = 0._wp    ! south
       ENDIF
       !
@@ -3856,6 +1346,6 @@
          !
          SELECT CASE ( jpni )
-         CASE ( 1 )     ;   CALL lbc_nfd        ( pt2d(1:jpi,1:jpj+jprj), cd_type, psgn, pr2dj=jprj )
-         CASE DEFAULT   ;   CALL mpp_lbc_north_icb( pt2d(1:jpi,1:jpj+jprj)  , cd_type, psgn , pr2dj=jprj  )
+                   CASE ( 1 )     ;   CALL lbc_nfd          ( pt2d(1:jpi,1:jpj+kextj), cd_type, psgn, kextj )
+                   CASE DEFAULT   ;   CALL mpp_lbc_north_icb( pt2d(1:jpi,1:jpj+kextj), cd_type, psgn, kextj )
          END SELECT
          !
@@ -3868,7 +1358,7 @@
       SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
       CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
-         iihom = nlci-nreci-jpri
+         iihom = jpi-nreci-kexti
          DO jl = 1, ipreci
-            r2dew(:,jl,1) = pt2d(jpreci+jl,:)
+            r2dew(:,jl,1) = pt2d(nn_hls+jl,:)
             r2dwe(:,jl,1) = pt2d(iihom +jl,:)
          END DO
@@ -3876,26 +1366,30 @@
       !
       !                           ! Migrations
-      imigr = ipreci * ( jpj + 2*jprj)
+      imigr = ipreci * ( jpj + 2*kextj )
+      !
+      IF( ln_timing ) CALL tic_tac(.TRUE.)
       !
       SELECT CASE ( nbondi )
       CASE ( -1 )
-         CALL mppsend( 2, r2dwe(1-jprj,1,1), imigr, noea, ml_req1 )
-         CALL mpprecv( 1, r2dew(1-jprj,1,2), imigr, noea )
+         CALL mppsend( 2, r2dwe(1-kextj,1,1), imigr, noea, ml_req1 )
+         CALL mpprecv( 1, r2dew(1-kextj,1,2), imigr, noea )
          IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
       CASE ( 0 )
-         CALL mppsend( 1, r2dew(1-jprj,1,1), imigr, nowe, ml_req1 )
-         CALL mppsend( 2, r2dwe(1-jprj,1,1), imigr, noea, ml_req2 )
-         CALL mpprecv( 1, r2dew(1-jprj,1,2), imigr, noea )
-         CALL mpprecv( 2, r2dwe(1-jprj,1,2), imigr, nowe )
+         CALL mppsend( 1, r2dew(1-kextj,1,1), imigr, nowe, ml_req1 )
+         CALL mppsend( 2, r2dwe(1-kextj,1,1), imigr, noea, ml_req2 )
+         CALL mpprecv( 1, r2dew(1-kextj,1,2), imigr, noea )
+         CALL mpprecv( 2, r2dwe(1-kextj,1,2), imigr, nowe )
          IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
          IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
       CASE ( 1 )
-         CALL mppsend( 1, r2dew(1-jprj,1,1), imigr, nowe, ml_req1 )
-         CALL mpprecv( 2, r2dwe(1-jprj,1,2), imigr, nowe )
+         CALL mppsend( 1, r2dew(1-kextj,1,1), imigr, nowe, ml_req1 )
+         CALL mpprecv( 2, r2dwe(1-kextj,1,2), imigr, nowe )
          IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
       END SELECT
       !
+      IF( ln_timing ) CALL tic_tac(.FALSE.)
+      !
       !                           ! Write Dirichlet lateral conditions
-      iihom = nlci - jpreci
+      iihom = jpi - nn_hls
       !
       SELECT CASE ( nbondi )
@@ -3906,10 +1400,10 @@
       CASE ( 0 )
          DO jl = 1, ipreci
-            pt2d(jl-jpri,:) = r2dwe(:,jl,2)
-            pt2d( iihom+jl,:) = r2dew(:,jl,2)
+            pt2d(jl-kexti,:) = r2dwe(:,jl,2)
+            pt2d(iihom+jl,:) = r2dew(:,jl,2)
          END DO
       CASE ( 1 )
          DO jl = 1, ipreci
-            pt2d(jl-jpri,:) = r2dwe(:,jl,2)
+            pt2d(jl-kexti,:) = r2dwe(:,jl,2)
          END DO
       END SELECT
@@ -3921,34 +1415,38 @@
       !
       IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
-         ijhom = nlcj-nrecj-jprj
+         ijhom = jpj-nrecj-kextj
          DO jl = 1, iprecj
             r2dsn(:,jl,1) = pt2d(:,ijhom +jl)
-            r2dns(:,jl,1) = pt2d(:,jprecj+jl)
+            r2dns(:,jl,1) = pt2d(:,nn_hls+jl)
          END DO
       ENDIF
       !
       !                           ! Migrations
-      imigr = iprecj * ( jpi + 2*jpri )
+      imigr = iprecj * ( jpi + 2*kexti )
+      !
+      IF( ln_timing ) CALL tic_tac(.TRUE.)
       !
       SELECT CASE ( nbondj )
       CASE ( -1 )
-         CALL mppsend( 4, r2dsn(1-jpri,1,1), imigr, nono, ml_req1 )
-         CALL mpprecv( 3, r2dns(1-jpri,1,2), imigr, nono )
+         CALL mppsend( 4, r2dsn(1-kexti,1,1), imigr, nono, ml_req1 )
+         CALL mpprecv( 3, r2dns(1-kexti,1,2), imigr, nono )
          IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
       CASE ( 0 )
-         CALL mppsend( 3, r2dns(1-jpri,1,1), imigr, noso, ml_req1 )
-         CALL mppsend( 4, r2dsn(1-jpri,1,1), imigr, nono, ml_req2 )
-         CALL mpprecv( 3, r2dns(1-jpri,1,2), imigr, nono )
-         CALL mpprecv( 4, r2dsn(1-jpri,1,2), imigr, noso )
+         CALL mppsend( 3, r2dns(1-kexti,1,1), imigr, noso, ml_req1 )
+         CALL mppsend( 4, r2dsn(1-kexti,1,1), imigr, nono, ml_req2 )
+         CALL mpprecv( 3, r2dns(1-kexti,1,2), imigr, nono )
+         CALL mpprecv( 4, r2dsn(1-kexti,1,2), imigr, noso )
          IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
          IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
       CASE ( 1 )
-         CALL mppsend( 3, r2dns(1-jpri,1,1), imigr, noso, ml_req1 )
-         CALL mpprecv( 4, r2dsn(1-jpri,1,2), imigr, noso )
+         CALL mppsend( 3, r2dns(1-kexti,1,1), imigr, noso, ml_req1 )
+         CALL mpprecv( 4, r2dsn(1-kexti,1,2), imigr, noso )
          IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
       END SELECT
       !
+      IF( ln_timing ) CALL tic_tac(.FALSE.)
+      !
       !                           ! Write Dirichlet lateral conditions
-      ijhom = nlcj - jprecj
+      ijhom = jpj - nn_hls
       !
       SELECT CASE ( nbondj )
@@ -3959,15 +1457,408 @@
       CASE ( 0 )
          DO jl = 1, iprecj
-            pt2d(:,jl-jprj) = r2dsn(:,jl,2)
-            pt2d(:,ijhom+jl ) = r2dns(:,jl,2)
+            pt2d(:,jl-kextj) = r2dsn(:,jl,2)
+            pt2d(:,ijhom+jl) = r2dns(:,jl,2)
          END DO
       CASE ( 1 )
          DO jl = 1, iprecj
-            pt2d(:,jl-jprj) = r2dsn(:,jl,2)
+            pt2d(:,jl-kextj) = r2dsn(:,jl,2)
          END DO
       END SELECT
-
+      !
    END SUBROUTINE mpp_lnk_2d_icb
+
+
+   SUBROUTINE mpp_report( cdname, kpk, kpl, kpf, ld_lbc, ld_glb, ld_dlg )
+      !!----------------------------------------------------------------------
+      !!                  ***  routine mpp_report  ***
+      !!
+      !! ** Purpose :   report use of mpp routines per time-setp
+      !!
+      !!----------------------------------------------------------------------
+      CHARACTER(len=*),           INTENT(in   ) ::   cdname      ! name of the calling subroutine
+      INTEGER         , OPTIONAL, INTENT(in   ) ::   kpk, kpl, kpf
+      LOGICAL         , OPTIONAL, INTENT(in   ) ::   ld_lbc, ld_glb, ld_dlg
+      !!
+      LOGICAL ::   ll_lbc, ll_glb, ll_dlg
+      INTEGER ::    ji,  jj,  jk,  jh, jf   ! dummy loop indices
+      !!----------------------------------------------------------------------
+      !
+      ll_lbc = .FALSE.
+      IF( PRESENT(ld_lbc) ) ll_lbc = ld_lbc
+      ll_glb = .FALSE.
+      IF( PRESENT(ld_glb) ) ll_glb = ld_glb
+      ll_dlg = .FALSE.
+      IF( PRESENT(ld_dlg) ) ll_dlg = ld_dlg
+      !
+      ! find the smallest common frequency: default = frequency product, if multiple, choose the larger of the 2 frequency
+      IF( ncom_dttrc /= 1 )   CALL ctl_stop( 'STOP', 'mpp_report, ncom_dttrc /= 1 not coded...' ) 
+      ncom_freq = ncom_fsbc
+      !
+      IF ( ncom_stp == nit000+ncom_freq ) THEN   ! avoid to count extra communications in potential initializations at nit000
+         IF( ll_lbc ) THEN
+            IF( .NOT. ALLOCATED(ncomm_sequence) ) ALLOCATE( ncomm_sequence(ncom_rec_max,2) )
+            IF( .NOT. ALLOCATED(    crname_lbc) ) ALLOCATE(     crname_lbc(ncom_rec_max  ) )
+            n_sequence_lbc = n_sequence_lbc + 1
+            IF( n_sequence_lbc > ncom_rec_max ) CALL ctl_stop( 'STOP', 'lib_mpp, increase ncom_rec_max' )   ! deadlock
+            crname_lbc(n_sequence_lbc) = cdname     ! keep the name of the calling routine
+            ncomm_sequence(n_sequence_lbc,1) = kpk*kpl   ! size of 3rd and 4th dimensions
+            ncomm_sequence(n_sequence_lbc,2) = kpf       ! number of arrays to be treated (multi)
+         ENDIF
+         IF( ll_glb ) THEN
+            IF( .NOT. ALLOCATED(crname_glb) ) ALLOCATE( crname_glb(ncom_rec_max) )
+            n_sequence_glb = n_sequence_glb + 1
+            IF( n_sequence_glb > ncom_rec_max ) CALL ctl_stop( 'STOP', 'lib_mpp, increase ncom_rec_max' )   ! deadlock
+            crname_glb(n_sequence_glb) = cdname     ! keep the name of the calling routine
+         ENDIF
+         IF( ll_dlg ) THEN
+            IF( .NOT. ALLOCATED(crname_dlg) ) ALLOCATE( crname_dlg(ncom_rec_max) )
+            n_sequence_dlg = n_sequence_dlg + 1
+            IF( n_sequence_dlg > ncom_rec_max ) CALL ctl_stop( 'STOP', 'lib_mpp, increase ncom_rec_max' )   ! deadlock
+            crname_dlg(n_sequence_dlg) = cdname     ! keep the name of the calling routine
+         ENDIF
+      ELSE IF ( ncom_stp == nit000+2*ncom_freq ) THEN
+         CALL ctl_opn( numcom, 'communication_report.txt', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE., narea )
+         WRITE(numcom,*) ' '
+         WRITE(numcom,*) ' ------------------------------------------------------------'
+         WRITE(numcom,*) ' Communication pattern report (second oce+sbc+top time step):'
+         WRITE(numcom,*) ' ------------------------------------------------------------'
+         WRITE(numcom,*) ' '
+         WRITE(numcom,'(A,I4)') ' Exchanged halos : ', n_sequence_lbc
+         jj = 0; jk = 0; jf = 0; jh = 0
+         DO ji = 1, n_sequence_lbc
+            IF ( ncomm_sequence(ji,1) .GT. 1 ) jk = jk + 1
+            IF ( ncomm_sequence(ji,2) .GT. 1 ) jf = jf + 1
+            IF ( ncomm_sequence(ji,1) .GT. 1 .AND. ncomm_sequence(ji,2) .GT. 1 ) jj = jj + 1
+            jh = MAX (jh, ncomm_sequence(ji,1)*ncomm_sequence(ji,2))
+         END DO
+         WRITE(numcom,'(A,I3)') ' 3D Exchanged halos : ', jk
+         WRITE(numcom,'(A,I3)') ' Multi arrays exchanged halos : ', jf
+         WRITE(numcom,'(A,I3)') '   from which 3D : ', jj
+         WRITE(numcom,'(A,I10)') ' Array max size : ', jh*jpi*jpj
+         WRITE(numcom,*) ' '
+         WRITE(numcom,*) ' lbc_lnk called'
+         jj = 1
+         DO ji = 2, n_sequence_lbc
+            IF( crname_lbc(ji-1) /= crname_lbc(ji) ) THEN
+               WRITE(numcom,'(A, I4, A, A)') ' - ', jj,' times by subroutine ', TRIM(crname_lbc(ji-1))
+               jj = 0
+            END IF
+            jj = jj + 1 
+         END DO
+         WRITE(numcom,'(A, I4, A, A)') ' - ', jj,' times by subroutine ', TRIM(crname_lbc(n_sequence_lbc))
+         WRITE(numcom,*) ' '
+         IF ( n_sequence_glb > 0 ) THEN
+            WRITE(numcom,'(A,I4)') ' Global communications : ', n_sequence_glb
+            jj = 1
+            DO ji = 2, n_sequence_glb
+               IF( crname_glb(ji-1) /= crname_glb(ji) ) THEN
+                  WRITE(numcom,'(A, I4, A, A)') ' - ', jj,' times by subroutine ', TRIM(crname_glb(ji-1))
+                  jj = 0
+               END IF
+               jj = jj + 1 
+            END DO
+            WRITE(numcom,'(A, I4, A, A)') ' - ', jj,' times by subroutine ', TRIM(crname_glb(n_sequence_glb))
+            DEALLOCATE(crname_glb)
+         ELSE
+            WRITE(numcom,*) ' No MPI global communication '
+         ENDIF
+         WRITE(numcom,*) ' '
+         IF ( n_sequence_dlg > 0 ) THEN
+            WRITE(numcom,'(A,I4)') ' Delayed global communications : ', n_sequence_dlg
+            jj = 1
+            DO ji = 2, n_sequence_dlg
+               IF( crname_dlg(ji-1) /= crname_dlg(ji) ) THEN
+                  WRITE(numcom,'(A, I4, A, A)') ' - ', jj,' times by subroutine ', TRIM(crname_dlg(ji-1))
+                  jj = 0
+               END IF
+               jj = jj + 1 
+            END DO
+            WRITE(numcom,'(A, I4, A, A)') ' - ', jj,' times by subroutine ', TRIM(crname_dlg(n_sequence_dlg))
+            DEALLOCATE(crname_dlg)
+         ELSE
+            WRITE(numcom,*) ' No MPI delayed global communication '
+         ENDIF
+         WRITE(numcom,*) ' '
+         WRITE(numcom,*) ' -----------------------------------------------'
+         WRITE(numcom,*) ' '
+         DEALLOCATE(ncomm_sequence)
+         DEALLOCATE(crname_lbc)
+      ENDIF
+   END SUBROUTINE mpp_report
+
    
+   SUBROUTINE tic_tac (ld_tic, ld_global)
+
+    LOGICAL,           INTENT(IN) :: ld_tic
+    LOGICAL, OPTIONAL, INTENT(IN) :: ld_global
+    REAL(wp), DIMENSION(2), SAVE :: tic_wt
+    REAL(wp),               SAVE :: tic_ct = 0._wp
+    INTEGER :: ii
+
+    IF( ncom_stp <= nit000 ) RETURN
+    IF( ncom_stp == nitend ) RETURN
+    ii = 1
+    IF( PRESENT( ld_global ) ) THEN
+       IF( ld_global ) ii = 2
+    END IF
+    
+    IF ( ld_tic ) THEN
+       tic_wt(ii) = MPI_Wtime()                                                    ! start count tic->tac (waiting time)
+       IF ( tic_ct > 0.0_wp ) compute_time = compute_time + MPI_Wtime() - tic_ct   ! cumulate count tac->tic
+    ELSE
+       waiting_time(ii) = waiting_time(ii) + MPI_Wtime() - tic_wt(ii)              ! cumulate count tic->tac
+       tic_ct = MPI_Wtime()                                                        ! start count tac->tic (waiting time)
+    ENDIF
+    
+   END SUBROUTINE tic_tac
+
+   
+#else
+   !!----------------------------------------------------------------------
+   !!   Default case:            Dummy module        share memory computing
+   !!----------------------------------------------------------------------
+   USE in_out_manager
+
+   INTERFACE mpp_sum
+      MODULE PROCEDURE mppsum_int, mppsum_a_int, mppsum_real, mppsum_a_real, mppsum_realdd, mppsum_a_realdd
+   END INTERFACE
+   INTERFACE mpp_max
+      MODULE PROCEDURE mppmax_a_int, mppmax_int, mppmax_a_real, mppmax_real
+   END INTERFACE
+   INTERFACE mpp_min
+      MODULE PROCEDURE mppmin_a_int, mppmin_int, mppmin_a_real, mppmin_real
+   END INTERFACE
+   INTERFACE mpp_minloc
+      MODULE PROCEDURE mpp_minloc2d ,mpp_minloc3d
+   END INTERFACE
+   INTERFACE mpp_maxloc
+      MODULE PROCEDURE mpp_maxloc2d ,mpp_maxloc3d
+   END INTERFACE
+
+   LOGICAL, PUBLIC, PARAMETER ::   lk_mpp = .FALSE.      !: mpp flag
+   LOGICAL, PUBLIC            ::   ln_nnogather          !: namelist control of northfold comms (needed here in case "key_mpp_mpi" is not used)
+   INTEGER, PUBLIC            ::   mpi_comm_oce          ! opa local communicator
+
+   INTEGER, PARAMETER, PUBLIC               ::   nbdelay = 0   ! make sure we don't enter loops: DO ji = 1, nbdelay
+   CHARACTER(len=32), DIMENSION(1), PUBLIC  ::   c_delaylist = 'empty'
+   CHARACTER(len=32), DIMENSION(1), PUBLIC  ::   c_delaycpnt = 'empty'
+   LOGICAL, PUBLIC                          ::   l_full_nf_update = .TRUE.
+   TYPE ::   DELAYARR
+      REAL(   wp), POINTER, DIMENSION(:) ::  z1d => NULL()
+      COMPLEX(wp), POINTER, DIMENSION(:) ::  y1d => NULL()
+   END TYPE DELAYARR
+   TYPE( DELAYARR ), DIMENSION(1), PUBLIC  ::   todelay              
+   INTEGER,  PUBLIC, DIMENSION(1)           ::   ndelayid = -1
+   !!----------------------------------------------------------------------
+CONTAINS
+
+   INTEGER FUNCTION lib_mpp_alloc(kumout)          ! Dummy function
+      INTEGER, INTENT(in) ::   kumout
+      lib_mpp_alloc = 0
+   END FUNCTION lib_mpp_alloc
+
+   FUNCTION mynode( ldtxt, ldname, kumnam_ref, knumnam_cfg,  kumond , kstop, localComm ) RESULT (function_value)
+      INTEGER, OPTIONAL            , INTENT(in   ) ::   localComm
+      CHARACTER(len=*),DIMENSION(:) ::   ldtxt
+      CHARACTER(len=*) ::   ldname
+      INTEGER ::   kumnam_ref, knumnam_cfg , kumond , kstop
+      IF( PRESENT( localComm ) ) mpi_comm_oce = localComm
+      function_value = 0
+      IF( .FALSE. )   ldtxt(:) = 'never done'
+      CALL ctl_opn( kumond, TRIM(ldname), 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. , 1 )
+   END FUNCTION mynode
+
+   SUBROUTINE mppsync                       ! Dummy routine
+   END SUBROUTINE mppsync
+
+   !!----------------------------------------------------------------------
+   !!    ***  mppmax_a_int, mppmax_int, mppmax_a_real, mppmax_real  ***
+   !!   
+   !!----------------------------------------------------------------------
+   !!
+#  define OPERATION_MAX
+#  define INTEGER_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppmax_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppmax_a_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef INTEGER_TYPE
+!
+#  define REAL_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppmax_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppmax_a_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef REAL_TYPE
+#  undef OPERATION_MAX
+   !!----------------------------------------------------------------------
+   !!    ***  mppmin_a_int, mppmin_int, mppmin_a_real, mppmin_real  ***
+   !!   
+   !!----------------------------------------------------------------------
+   !!
+#  define OPERATION_MIN
+#  define INTEGER_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppmin_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppmin_a_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef INTEGER_TYPE
+!
+#  define REAL_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppmin_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppmin_a_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef REAL_TYPE
+#  undef OPERATION_MIN
+
+   !!----------------------------------------------------------------------
+   !!    ***  mppsum_a_int, mppsum_int, mppsum_a_real, mppsum_real  ***
+   !!   
+   !!   Global sum of 1D array or a variable (integer, real or complex)
+   !!----------------------------------------------------------------------
+   !!
+#  define OPERATION_SUM
+#  define INTEGER_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppsum_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppsum_a_int
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef INTEGER_TYPE
+!
+#  define REAL_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppsum_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppsum_a_real
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef REAL_TYPE
+#  undef OPERATION_SUM
+
+#  define OPERATION_SUM_DD
+#  define COMPLEX_TYPE
+#  define DIM_0d
+#     define ROUTINE_ALLREDUCE           mppsum_realdd
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_0d
+#  define DIM_1d
+#     define ROUTINE_ALLREDUCE           mppsum_a_realdd
+#     include "mpp_allreduce_generic.h90"
+#     undef ROUTINE_ALLREDUCE
+#  undef DIM_1d
+#  undef COMPLEX_TYPE
+#  undef OPERATION_SUM_DD
+
+   !!----------------------------------------------------------------------
+   !!    ***  mpp_minloc2d, mpp_minloc3d, mpp_maxloc2d, mpp_maxloc3d
+   !!   
+   !!----------------------------------------------------------------------
+   !!
+#  define OPERATION_MINLOC
+#  define DIM_2d
+#     define ROUTINE_LOC           mpp_minloc2d
+#     include "mpp_loc_generic.h90"
+#     undef ROUTINE_LOC
+#  undef DIM_2d
+#  define DIM_3d
+#     define ROUTINE_LOC           mpp_minloc3d
+#     include "mpp_loc_generic.h90"
+#     undef ROUTINE_LOC
+#  undef DIM_3d
+#  undef OPERATION_MINLOC
+
+#  define OPERATION_MAXLOC
+#  define DIM_2d
+#     define ROUTINE_LOC           mpp_maxloc2d
+#     include "mpp_loc_generic.h90"
+#     undef ROUTINE_LOC
+#  undef DIM_2d
+#  define DIM_3d
+#     define ROUTINE_LOC           mpp_maxloc3d
+#     include "mpp_loc_generic.h90"
+#     undef ROUTINE_LOC
+#  undef DIM_3d
+#  undef OPERATION_MAXLOC
+
+   SUBROUTINE mpp_delay_sum( cdname, cdelay, y_in, pout, ldlast, kcom )
+      CHARACTER(len=*), INTENT(in   )               ::   cdname  ! name of the calling subroutine
+      CHARACTER(len=*), INTENT(in   )               ::   cdelay  ! name (used as id) of the delayed operation
+      COMPLEX(wp),      INTENT(in   ), DIMENSION(:) ::   y_in
+      REAL(wp),         INTENT(  out), DIMENSION(:) ::   pout
+      LOGICAL,          INTENT(in   )               ::   ldlast  ! true if this is the last time we call this routine
+      INTEGER,          INTENT(in   ), OPTIONAL     ::   kcom
+      !
+      pout(:) = REAL(y_in(:), wp)
+   END SUBROUTINE mpp_delay_sum
+
+   SUBROUTINE mpp_delay_max( cdname, cdelay, p_in, pout, ldlast, kcom )
+      CHARACTER(len=*), INTENT(in   )               ::   cdname  ! name of the calling subroutine
+      CHARACTER(len=*), INTENT(in   )               ::   cdelay  ! name (used as id) of the delayed operation
+      REAL(wp),         INTENT(in   ), DIMENSION(:) ::   p_in
+      REAL(wp),         INTENT(  out), DIMENSION(:) ::   pout
+      LOGICAL,          INTENT(in   )               ::   ldlast  ! true if this is the last time we call this routine
+      INTEGER,          INTENT(in   ), OPTIONAL     ::   kcom
+      !
+      pout(:) = p_in(:)
+   END SUBROUTINE mpp_delay_max
+
+   SUBROUTINE mpp_delay_rcv( kid )
+      INTEGER,INTENT(in   )      ::  kid 
+      WRITE(*,*) 'mpp_delay_rcv: You should not have seen this print! error?', kid
+   END SUBROUTINE mpp_delay_rcv
+   
+   SUBROUTINE mppstop( ldfinal, ld_force_abort )
+      LOGICAL, OPTIONAL, INTENT(in) :: ldfinal    ! source process number
+      LOGICAL, OPTIONAL, INTENT(in) :: ld_force_abort    ! source process number
+      STOP      ! non MPP case, just stop the run
+   END SUBROUTINE mppstop
+
+   SUBROUTINE mpp_ini_znl( knum )
+      INTEGER :: knum
+      WRITE(*,*) 'mpp_ini_znl: You should not have seen this print! error?', knum
+   END SUBROUTINE mpp_ini_znl
+
+   SUBROUTINE mpp_comm_free( kcom )
+      INTEGER :: kcom
+      WRITE(*,*) 'mpp_comm_free: You should not have seen this print! error?', kcom
+   END SUBROUTINE mpp_comm_free
+   
+#endif
 
    !!----------------------------------------------------------------------
@@ -3988,25 +1879,28 @@
       !
       nstop = nstop + 1
-      IF(lwp) THEN
-         WRITE(numout,cform_err)
-         IF( PRESENT(cd1 ) )   WRITE(numout,*) cd1
-         IF( PRESENT(cd2 ) )   WRITE(numout,*) cd2
-         IF( PRESENT(cd3 ) )   WRITE(numout,*) cd3
-         IF( PRESENT(cd4 ) )   WRITE(numout,*) cd4
-         IF( PRESENT(cd5 ) )   WRITE(numout,*) cd5
-         IF( PRESENT(cd6 ) )   WRITE(numout,*) cd6
-         IF( PRESENT(cd7 ) )   WRITE(numout,*) cd7
-         IF( PRESENT(cd8 ) )   WRITE(numout,*) cd8
-         IF( PRESENT(cd9 ) )   WRITE(numout,*) cd9
-         IF( PRESENT(cd10) )   WRITE(numout,*) cd10
-      ENDIF
+
+      ! force to open ocean.output file
+      IF( numout == 6 ) CALL ctl_opn( numout, 'ocean.output', 'APPEND', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. )
+       
+      WRITE(numout,cform_err)
+      IF( PRESENT(cd1 ) )   WRITE(numout,*) TRIM(cd1)
+      IF( PRESENT(cd2 ) )   WRITE(numout,*) TRIM(cd2)
+      IF( PRESENT(cd3 ) )   WRITE(numout,*) TRIM(cd3)
+      IF( PRESENT(cd4 ) )   WRITE(numout,*) TRIM(cd4)
+      IF( PRESENT(cd5 ) )   WRITE(numout,*) TRIM(cd5)
+      IF( PRESENT(cd6 ) )   WRITE(numout,*) TRIM(cd6)
+      IF( PRESENT(cd7 ) )   WRITE(numout,*) TRIM(cd7)
+      IF( PRESENT(cd8 ) )   WRITE(numout,*) TRIM(cd8)
+      IF( PRESENT(cd9 ) )   WRITE(numout,*) TRIM(cd9)
+      IF( PRESENT(cd10) )   WRITE(numout,*) TRIM(cd10)
+
                                CALL FLUSH(numout    )
       IF( numstp     /= -1 )   CALL FLUSH(numstp    )
-      IF( numsol     /= -1 )   CALL FLUSH(numsol    )
+      IF( numrun     /= -1 )   CALL FLUSH(numrun    )
       IF( numevo_ice /= -1 )   CALL FLUSH(numevo_ice)
       !
       IF( cd1 == 'STOP' ) THEN
-         IF(lwp) WRITE(numout,*)  'huge E-R-R-O-R : immediate stop'
-         CALL mppstop()
+         WRITE(numout,*)  'huge E-R-R-O-R : immediate stop'
+         CALL mppstop(ld_force_abort = .true.)
       ENDIF
       !
@@ -4029,14 +1923,14 @@
       IF(lwp) THEN
          WRITE(numout,cform_war)
-         IF( PRESENT(cd1 ) ) WRITE(numout,*) cd1
-         IF( PRESENT(cd2 ) ) WRITE(numout,*) cd2
-         IF( PRESENT(cd3 ) ) WRITE(numout,*) cd3
-         IF( PRESENT(cd4 ) ) WRITE(numout,*) cd4
-         IF( PRESENT(cd5 ) ) WRITE(numout,*) cd5
-         IF( PRESENT(cd6 ) ) WRITE(numout,*) cd6
-         IF( PRESENT(cd7 ) ) WRITE(numout,*) cd7
-         IF( PRESENT(cd8 ) ) WRITE(numout,*) cd8
-         IF( PRESENT(cd9 ) ) WRITE(numout,*) cd9
-         IF( PRESENT(cd10) ) WRITE(numout,*) cd10
+         IF( PRESENT(cd1 ) ) WRITE(numout,*) TRIM(cd1)
+         IF( PRESENT(cd2 ) ) WRITE(numout,*) TRIM(cd2)
+         IF( PRESENT(cd3 ) ) WRITE(numout,*) TRIM(cd3)
+         IF( PRESENT(cd4 ) ) WRITE(numout,*) TRIM(cd4)
+         IF( PRESENT(cd5 ) ) WRITE(numout,*) TRIM(cd5)
+         IF( PRESENT(cd6 ) ) WRITE(numout,*) TRIM(cd6)
+         IF( PRESENT(cd7 ) ) WRITE(numout,*) TRIM(cd7)
+         IF( PRESENT(cd8 ) ) WRITE(numout,*) TRIM(cd8)
+         IF( PRESENT(cd9 ) ) WRITE(numout,*) TRIM(cd9)
+         IF( PRESENT(cd10) ) WRITE(numout,*) TRIM(cd10)
       ENDIF
       CALL FLUSH(numout)
@@ -4073,15 +1967,25 @@
          IF( karea > 1 )   WRITE(clfile, "(a,'_',i4.4)") TRIM(clfile), karea-1
       ENDIF
+#if defined key_agrif
+      IF( .NOT. Agrif_Root() )   clfile = TRIM(Agrif_CFixed())//'_'//TRIM(clfile)
+      knum=Agrif_Get_Unit()
+#else
       knum=get_unit()
+#endif
+      IF( TRIM(cdfile) == '/dev/null' )   clfile = TRIM(cdfile)   ! force the use of /dev/null
       !
       iost=0
-      IF( cdacce(1:6) == 'DIRECT' )  THEN
-         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat, RECL=klengh, ERR=100, IOSTAT=iost )
+      IF( cdacce(1:6) == 'DIRECT' )  THEN         ! cdacce has always more than 6 characters
+         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat, RECL=klengh         , ERR=100, IOSTAT=iost )
+      ELSE IF( TRIM(cdstat) == 'APPEND' )  THEN   ! cdstat can have less than 6 characters
+         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS='UNKNOWN', POSITION='APPEND', ERR=100, IOSTAT=iost )
       ELSE
-         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat             , ERR=100, IOSTAT=iost )
-      ENDIF
+         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat                      , ERR=100, IOSTAT=iost )
+      ENDIF
+      IF( iost /= 0 .AND. TRIM(clfile) == '/dev/null' ) &   ! for windows
+         &  OPEN(UNIT=knum,FILE='NUL', FORM=cdform, ACCESS=cdacce, STATUS=cdstat                      , ERR=100, IOSTAT=iost )   
       IF( iost == 0 ) THEN
          IF(ldwp) THEN
-            WRITE(kout,*) '     file   : ', clfile,' open ok'
+            WRITE(kout,*) '     file   : ', TRIM(clfile),' open ok'
             WRITE(kout,*) '     unit   = ', knum
             WRITE(kout,*) '     status = ', cdstat
@@ -4095,5 +1999,5 @@
          IF(ldwp) THEN
             WRITE(kout,*)
-            WRITE(kout,*) ' ===>>>> : bad opening file: ', clfile
+            WRITE(kout,*) ' ===>>>> : bad opening file: ', TRIM(clfile)
             WRITE(kout,*) ' =======   ===  '
             WRITE(kout,*) '           unit   = ', knum
@@ -4104,5 +2008,17 @@
             WRITE(kout,*) '           we stop. verify the file '
             WRITE(kout,*)
+         ELSE  !!! Force writing to make sure we get the information - at least once - in this violent STOP!!
+            WRITE(*,*)
+            WRITE(*,*) ' ===>>>> : bad opening file: ', TRIM(clfile)
+            WRITE(*,*) ' =======   ===  '
+            WRITE(*,*) '           unit   = ', knum
+            WRITE(*,*) '           status = ', cdstat
+            WRITE(*,*) '           form   = ', cdform
+            WRITE(*,*) '           access = ', cdacce
+            WRITE(*,*) '           iostat = ', iost
+            WRITE(*,*) '           we stop. verify the file '
+            WRITE(*,*)
          ENDIF
+         CALL FLUSH( kout ) 
          STOP 'ctl_opn bad opening'
       ENDIF
@@ -4121,9 +2037,9 @@
       INTEGER         , INTENT(inout) ::   kios    ! IO status after reading the namelist
       CHARACTER(len=*), INTENT(in   ) ::   cdnam   ! group name of namelist for which error occurs
-      CHARACTER(len=4)                ::   clios   ! string to convert iostat in character for print
+      CHARACTER(len=5)                ::   clios   ! string to convert iostat in character for print
       LOGICAL         , INTENT(in   ) ::   ldwp    ! boolean term for print
       !!----------------------------------------------------------------------
       !
-      WRITE (clios, '(I4.0)')   kios
+      WRITE (clios, '(I5.0)')   kios
       IF( kios < 0 ) THEN         
          CALL ctl_warn( 'end of record or file while reading namelist '   &
