Index: NEMO/branches/2020/dev_r12558_HPC-08_epico_Extra_Halo/src/OCE/LBC/mppini.F90
===================================================================
--- NEMO/branches/2020/dev_r12558_HPC-08_epico_Extra_Halo/src/OCE/LBC/mppini.F90	(revision 13174)
+++ NEMO/branches/2020/dev_r12558_HPC-08_epico_Extra_Halo/src/OCE/LBC/mppini.F90	(revision 13176)
@@ -32,6 +32,9 @@
    PRIVATE
 
-   PUBLIC   mpp_init   ! called by opa.F90
-
+   PUBLIC   mpp_init       ! called by nemogcm.F90
+   PUBLIC   mpp_getnum     ! called by prtctl
+   PUBLIC   mpp_basesplit  ! called by prtctl
+   PUBLIC   mpp_is_ocean   ! called by prtctl
+   
    INTEGER ::   numbot = -1   ! 'bottom_level' local logical unit
    INTEGER ::   numbdy = -1   ! 'bdy_msk'      local logical unit
@@ -76,5 +79,6 @@
       jpnj   = 1
       jpnij  = jpni*jpnj
-      nimpp  = 1           ! 
+      nn_hls = 1
+      nimpp  = 1
       njmpp  = 1
       nbondi = 2
@@ -137,7 +141,6 @@
       INTEGER ::   ji, jj, jn, jproc, jarea   ! dummy loop indices
       INTEGER ::   inijmin
-      INTEGER ::   i2add
       INTEGER ::   inum                       ! local logical unit
-      INTEGER ::   idir, ifreq, icont         ! local integers
+      INTEGER ::   idir, ifreq                ! local integers
       INTEGER ::   ii, il1, ili, imil         !   -       -
       INTEGER ::   ij, il2, ilj, ijm1         !   -       -
@@ -186,4 +189,5 @@
          ENDIF
             WRITE(numout,*) '      avoid use of mpi_allgather at the north fold  ln_nnogather = ', ln_nnogather
+            WRITE(numout,*) '      halo width (applies to both rows and columns)       nn_hls = ', nn_hls
       ENDIF
       !
@@ -225,7 +229,7 @@
          CALL bestpartition( mppsize, inbi, inbj, icnt2 )    ! best mpi decomposition for mppsize mpi processes
          ! largest subdomain size for mpi decoposition jpni*jpnj given in the namelist
-         CALL basic_decomposition( jpiglo, jpjglo, nn_hls, jpni, jpnj, jpimax, jpjmax )
+         CALL mpp_basesplit( jpiglo, jpjglo, nn_hls, jpni, jpnj, jpimax, jpjmax )
          ! largest subdomain size for mpi decoposition inbi*inbj given by bestpartition
-         CALL basic_decomposition( jpiglo, jpjglo, nn_hls, inbi, inbj,  iimax,  ijmax )
+         CALL mpp_basesplit( jpiglo, jpjglo, nn_hls, inbi, inbj,  iimax,  ijmax )
          icnt1 = jpni*jpnj - mppsize   ! number of land subdomains that should be removed to use mppsize mpi processes
          IF(lwp) THEN
@@ -258,5 +262,5 @@
       ! look for land mpi subdomains...
       ALLOCATE( llisoce(jpni,jpnj) )
-      CALL is_ocean( jpni, jpnj, llisoce )
+      CALL mpp_is_ocean( llisoce )
       inijmin = COUNT( llisoce )   ! number of oce subdomains
 
@@ -316,7 +320,4 @@
 9003  FORMAT (a, i5)
 
-      IF( numbot /= -1 )   CALL iom_close( numbot )
-      IF( numbdy /= -1 )   CALL iom_close( numbdy )
-      
       ALLOCATE(  nfimpp(jpni ) , nfproc(jpni ) ,   nfjpi(jpni ) ,                     &
          &       nimppt(jpnij) , ibonit(jpnij) ,  jpiall(jpnij) ,  jpjall(jpnij) ,    &
@@ -346,5 +347,8 @@
       ! -----------------------------------
       !
-      CALL basic_decomposition( jpiglo, jpjglo, nn_hls, jpni, jpnj, jpimax, jpjmax, iimppt, ijmppt, ijpi, ijpj )
+      CALL mpp_basesplit( jpiglo, jpjglo, nn_hls, jpni, jpnj, jpimax, jpjmax, iimppt, ijmppt, ijpi, ijpj )
+      CALL mpp_getnum( llisoce, ipproc, iin, ijn )
+      !
+      nfproc(:) = ipproc(:,jpnj)
       nfimpp(:) = iimppt(:,jpnj)
       nfjpi (:) =   ijpi(:,jpnj)
@@ -357,4 +361,5 @@
          WRITE(numout,*) '      jpni = ', jpni  
          WRITE(numout,*) '      jpnj = ', jpnj
+         WRITE(numout,*) '     jpnij = ', jpnij
          WRITE(numout,*)
          WRITE(numout,*) '      sum ijpi(i,1) = ', sum(ijpi(:,1)), ' jpiglo = ', jpiglo
@@ -431,34 +436,4 @@
       ! ----------------------------
       !
-      ! specify which subdomains are oce subdomains; other are land subdomains
-      ipproc(:,:) = -1
-      icont = -1
-      DO jarea = 1, jpni*jpnj
-         iarea0 = jarea - 1
-         ii = 1 + MOD(iarea0,jpni)
-         ij = 1 +     iarea0/jpni
-         IF( llisoce(ii,ij) ) THEN
-            icont = icont + 1
-            ipproc(ii,ij) = icont
-            iin(icont+1) = ii
-            ijn(icont+1) = ij
-         ENDIF
-      END DO
-      ! if needed add some land subdomains to reach jpnij active subdomains
-      i2add = jpnij - inijmin
-      DO jarea = 1, jpni*jpnj
-         iarea0 = jarea - 1
-         ii = 1 + MOD(iarea0,jpni)
-         ij = 1 +     iarea0/jpni
-         IF( .NOT. llisoce(ii,ij) .AND. i2add > 0 ) THEN
-            icont = icont + 1
-            ipproc(ii,ij) = icont
-            iin(icont+1) = ii
-            ijn(icont+1) = ij
-            i2add = i2add - 1
-         ENDIF
-      END DO
-      nfproc(:) = ipproc(:,jpnj)
-
       ! neighbour treatment: change ibondi, ibondj if next to a land zone
       DO jarea = 1, jpni*jpnj
@@ -655,5 +630,4 @@
          WRITE(numout,*) '      nimpp  = ', nimpp
          WRITE(numout,*) '      njmpp  = ', njmpp
-         WRITE(numout,*) '      nn_hls = ', nn_hls 
       ENDIF
 
@@ -700,7 +674,7 @@
 
 
-    SUBROUTINE basic_decomposition( kiglo, kjglo, khls, knbi, knbj, kimax, kjmax, kimppt, kjmppt, klci, klcj)
-      !!----------------------------------------------------------------------
-      !!                  ***  ROUTINE basic_decomposition  ***
+    SUBROUTINE mpp_basesplit( kiglo, kjglo, khls, knbi, knbj, kimax, kjmax, kimppt, kjmppt, klci, klcj)
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE mpp_basesplit  ***
       !!                    
       !! ** Purpose :   Lay out the global domain over processors.
@@ -757,5 +731,5 @@
       klci(iresti+1:knbi ,:) = kimax-1
       IF( MINVAL(klci) < 2*i2hls ) THEN
-         WRITE(ctmp1,*) '   basic_decomposition: minimum value of jpi must be >= ', 2*i2hls
+         WRITE(ctmp1,*) '   mpp_basesplit: minimum value of jpi must be >= ', 2*i2hls
          WRITE(ctmp2,*) '   We have ', MINVAL(klci)
         CALL ctl_stop( 'STOP', ctmp1, ctmp2 )
@@ -775,5 +749,5 @@
       klcj(:,1:irestj) = kjmax
       IF( MINVAL(klcj) < 2*i2hls ) THEN
-         WRITE(ctmp1,*) '   basic_decomposition: minimum value of jpj must be >= ', 2*i2hls
+         WRITE(ctmp1,*) '   mpp_basesplit: minimum value of jpj must be >= ', 2*i2hls
          WRITE(ctmp2,*) '   We have ', MINVAL(klcj)
          CALL ctl_stop( 'STOP', ctmp1, ctmp2 )
@@ -802,5 +776,5 @@
       ENDIF
       
-   END SUBROUTINE basic_decomposition
+   END SUBROUTINE mpp_basesplit
 
 
@@ -909,5 +883,5 @@
       iszij1(:) = iszi1(:) * iszj1(:)
 
-      ! if therr is no land and no print
+      ! if there is no land and no print
       IF( .NOT. llist .AND. numbot == -1 .AND. numbdy == -1 ) THEN
          ! get the smaller partition which gives the smallest subdomain size
@@ -957,10 +931,10 @@
          ji = isz0   ! initialization with the largest value
          ALLOCATE( llisoce(inbi0(ji), inbj0(ji)) )
-         CALL is_ocean( inbi0(ji), inbj0(ji), llisoce ) ! Warning: must be call by all cores (call mpp_sum)
+         CALL mpp_is_ocean( llisoce )   ! Warning: must be call by all cores (call mpp_sum)
          inbijold = COUNT(llisoce)
          DEALLOCATE( llisoce )
          DO ji =isz0-1,1,-1
             ALLOCATE( llisoce(inbi0(ji), inbj0(ji)) )
-            CALL is_ocean( inbi0(ji), inbj0(ji), llisoce ) ! Warning: must be call by all cores (call mpp_sum)
+            CALL mpp_is_ocean( llisoce )   ! Warning: must be call by all cores (call mpp_sum)
             inbij = COUNT(llisoce)
             DEALLOCATE( llisoce )
@@ -988,5 +962,5 @@
          ii = ii -1 
          ALLOCATE( llisoce(inbi0(ii), inbj0(ii)) )
-         CALL is_ocean( inbi0(ii), inbj0(ii), llisoce )            ! must be done by all core
+         CALL mpp_is_ocean( llisoce )            ! must be done by all core
          inbij = COUNT(llisoce)
          DEALLOCATE( llisoce )
@@ -1052,9 +1026,9 @@
    
    
-   SUBROUTINE is_ocean( knbi, knbj, ldisoce )
-      !!----------------------------------------------------------------------
-      !!                  ***  ROUTINE mpp_init_nboce  ***
-      !!
-      !! ** Purpose : Check for a mpi domain decomposition knbi x knbj which
+   SUBROUTINE mpp_is_ocean( ldisoce )
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE mpp_is_ocean  ***
+      !!
+      !! ** Purpose : Check for a mpi domain decomposition inbi x inbj which
       !!              subdomains, including 1 halo (even if nn_hls>1), contain
       !!              at least 1 ocean point.
@@ -1065,17 +1039,16 @@
       !!              a subdomain with a closed boundary.
       !!
-      !! ** Method  : read knbj strips (of length Ni0glo) of the land-sea mask
-      !!----------------------------------------------------------------------
-      INTEGER,                       INTENT(in   ) ::   knbi, knbj     ! domain decomposition
-      LOGICAL, DIMENSION(knbi,knbj), INTENT(  out) ::   ldisoce        ! .true. if a sub domain constains 1 ocean point 
-      !
-      INTEGER, DIMENSION(knbi,knbj) ::   inboce                        ! number oce oce pint in each mpi subdomain
-      INTEGER, DIMENSION(knbi*knbj) ::   inboce_1d
+      !! ** Method  : read inbj strips (of length Ni0glo) of the land-sea mask
+      !!----------------------------------------------------------------------
+      LOGICAL, DIMENSION(:,:), INTENT(  out) ::   ldisoce        ! .true. if a sub domain constains 1 ocean point 
+      !
       INTEGER :: idiv, iimax, ijmax, iarea
-      INTEGER :: inx, iny, inry, isty
+      INTEGER :: inbi, inbj, inx, iny, inry, isty
       INTEGER :: ji, jn
-      LOGICAL, ALLOCATABLE, DIMENSION(:,:) ::   lloce                  ! lloce(i,j) = .true. if the point (i,j) is ocean 
+      INTEGER, ALLOCATABLE, DIMENSION(:,:) ::   inboce           ! number oce oce pint in each mpi subdomain
+      INTEGER, ALLOCATABLE, DIMENSION(:  ) ::   inboce_1d
       INTEGER, ALLOCATABLE, DIMENSION(:,:) ::   iimppt, ijpi
       INTEGER, ALLOCATABLE, DIMENSION(:,:) ::   ijmppt, ijpj
+      LOGICAL, ALLOCATABLE, DIMENSION(:,:) ::   lloce            ! lloce(i,j) = .true. if the point (i,j) is ocean 
       !!----------------------------------------------------------------------
       ! do nothing if there is no land-sea mask
@@ -1084,23 +1057,28 @@
          RETURN
       ENDIF
-
-      ! we want to read knbj strips of the land-sea mask. -> pick up knbj processes every idiv processes starting at 1
-      IF           ( knbj == 1 ) THEN   ;   idiv = mppsize
-      ELSE IF ( mppsize < knbj ) THEN   ;   idiv = 1
-      ELSE                              ;   idiv = ( mppsize - 1 ) / ( knbj - 1 )
-      ENDIF
+      !
+      inbi = SIZE( ldisoce, dim = 1 )
+      inbj = SIZE( ldisoce, dim = 2 )
+      !
+      ! we want to read inbj strips of the land-sea mask. -> pick up inbj processes every idiv processes starting at 1
+      IF           ( inbj == 1 ) THEN   ;   idiv = mppsize
+      ELSE IF ( mppsize < inbj ) THEN   ;   idiv = 1
+      ELSE                              ;   idiv = ( mppsize - 1 ) / ( inbj - 1 )
+      ENDIF
+      !
+      ALLOCATE( inboce(inbi,inbj), inboce_1d(inbi*inbj) )
       inboce(:,:) = 0          ! default no ocean point found
-
-      DO jn = 0, (knbj-1)/mppsize   ! if mppsize < knbj : more strips than mpi processes (because of potential land domains)
+      !
+      DO jn = 0, (inbj-1)/mppsize   ! if mppsize < inbj : more strips than mpi processes (because of potential land domains)
          !
          iarea = (narea-1)/idiv + jn * mppsize + 1                     ! involed process number (starting counting at 1)
-         IF( MOD( narea-1, idiv ) == 0 .AND. iarea <= knbj ) THEN      ! beware idiv can be = to 1
+         IF( MOD( narea-1, idiv ) == 0 .AND. iarea <= inbj ) THEN      ! beware idiv can be = to 1
             !
-            ALLOCATE( iimppt(knbi,knbj), ijmppt(knbi,knbj), ijpi(knbi,knbj), ijpj(knbi,knbj) )
-            CALL basic_decomposition( Ni0glo, Nj0glo, 0, knbi, knbj, iimax, ijmax, iimppt, ijmppt, ijpi, ijpj )
+            ALLOCATE( iimppt(inbi,inbj), ijmppt(inbi,inbj), ijpi(inbi,inbj), ijpj(inbi,inbj) )
+            CALL mpp_basesplit( Ni0glo, Nj0glo, 0, inbi, inbj, iimax, ijmax, iimppt, ijmppt, ijpi, ijpj )
             !
             inx = Ni0glo + 2   ;   iny = ijpj(1,iarea) + 2             ! strip size + 1 halo on each direction (even if nn_hls>1)
             ALLOCATE( lloce(inx, iny) )                                ! allocate the strip
-            inry = iny - COUNT( (/ iarea == 1, iarea == knbj /) )      ! number of point to read in y-direction
+            inry = iny - COUNT( (/ iarea == 1, iarea == inbj /) )      ! number of point to read in y-direction
             isty = 1 + COUNT( (/ iarea == 1 /) )                       ! read from the first or the second line?
             CALL readbot_strip( ijmppt(1,iarea) - 2 + isty, inry, lloce(2:inx-1, isty:inry+isty-1) )   ! read the strip
@@ -1113,5 +1091,5 @@
                ENDIF
             ENDIF
-            IF( iarea == knbj ) THEN                                   ! the last line was not read
+            IF( iarea == inbj ) THEN                                   ! the last line was not read
                IF( jperio == 2 .OR. jperio == 7 ) THEN                 !   north-south periodocity
                   CALL readbot_strip( 1, 1, lloce(2:inx-1,iny) )       !   read the first line -> last line of lloce
@@ -1127,5 +1105,5 @@
             ENDIF
             !
-            DO  ji = 1, knbi
+            DO  ji = 1, inbi
                inboce(ji,iarea) = COUNT( lloce(iimppt(ji,1):iimppt(ji,1)+ijpi(ji,1)+1,:) )   ! lloce as 2 points more than Ni0glo
             END DO
@@ -1137,10 +1115,11 @@
       END DO
    
-      inboce_1d = RESHAPE(inboce, (/ knbi*knbj /))
+      inboce_1d = RESHAPE(inboce, (/ inbi*inbj /))
       CALL mpp_sum( 'mppini', inboce_1d )
-      inboce = RESHAPE(inboce_1d, (/knbi, knbj/))
+      inboce = RESHAPE(inboce_1d, (/inbi, inbj/))
       ldisoce(:,:) = inboce(:,:) /= 0
-      !
-   END SUBROUTINE is_ocean
+      DEALLOCATE(inboce, inboce_1d)
+      !
+   END SUBROUTINE mpp_is_ocean
    
    
@@ -1155,7 +1134,7 @@
       !! ** Method  : read stipe of size (Ni0glo,...)
       !!----------------------------------------------------------------------
-      INTEGER                         , INTENT(in   ) :: kjstr       ! starting j position of the reading
-      INTEGER                         , INTENT(in   ) :: kjcnt       ! number of lines to read
-      LOGICAL, DIMENSION(Ni0glo,kjcnt), INTENT(  out) :: ldoce       ! ldoce(i,j) = .true. if the point (i,j) is ocean 
+      INTEGER                         , INTENT(in   ) ::   kjstr       ! starting j position of the reading
+      INTEGER                         , INTENT(in   ) ::   kjcnt       ! number of lines to read
+      LOGICAL, DIMENSION(Ni0glo,kjcnt), INTENT(  out) ::   ldoce       ! ldoce(i,j) = .true. if the point (i,j) is ocean 
       !
       INTEGER                           ::   inumsave                ! local logical unit
@@ -1180,4 +1159,57 @@
       !
    END SUBROUTINE readbot_strip
+
+
+   SUBROUTINE mpp_getnum( ldisoce, kproc, kipos, kjpos )
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE mpp_getnum  ***
+      !!
+      !! ** Purpose : give a number to each MPI subdomains (starting at 0)
+      !!
+      !! ** Method  : start from bottom left. First skip land subdomain, and finally use them if needed
+      !!----------------------------------------------------------------------
+      LOGICAL, DIMENSION(:,:), INTENT(in   ) ::   ldisoce     ! F if land process
+      INTEGER, DIMENSION(:,:), INTENT(  out) ::   kproc       ! subdomain number (-1 if supressed, starting at 0)
+      INTEGER, DIMENSION(  :), INTENT(  out) ::   kipos       ! i-position of the subdomain (from 1 to jpni)
+      INTEGER, DIMENSION(  :), INTENT(  out) ::   kjpos       ! j-position of the subdomain (from 1 to jpnj)
+      !
+      INTEGER :: ii, ij, jarea, iarea0
+      INTEGER :: icont, i2add , ini, inj, inij
+      !!----------------------------------------------------------------------
+      !
+      ini = SIZE(ldisoce, dim = 1)
+      inj = SIZE(ldisoce, dim = 2)
+      inij = SIZE(kipos)
+      !
+      ! specify which subdomains are oce subdomains; other are land subdomains
+      kproc(:,:) = -1
+      icont = -1
+      DO jarea = 1, ini*inj
+         iarea0 = jarea - 1
+         ii = 1 + MOD(iarea0,ini)
+         ij = 1 +     iarea0/ini
+         IF( ldisoce(ii,ij) ) THEN
+            icont = icont + 1
+            kproc(ii,ij) = icont
+            kipos(icont+1) = ii
+            kjpos(icont+1) = ij
+         ENDIF
+      END DO
+      ! if needed add some land subdomains to reach inij active subdomains
+      i2add = inij - COUNT( ldisoce )
+      DO jarea = 1, ini*inj
+         iarea0 = jarea - 1
+         ii = 1 + MOD(iarea0,ini)
+         ij = 1 +     iarea0/ini
+         IF( .NOT. ldisoce(ii,ij) .AND. i2add > 0 ) THEN
+            icont = icont + 1
+            kproc(ii,ij) = icont
+            kipos(icont+1) = ii
+            kjpos(icont+1) = ij
+            i2add = i2add - 1
+         ENDIF
+      END DO
+      !
+   END SUBROUTINE mpp_getnum
 
 
