New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
lib_mpp.F90 in branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/LBC – NEMO

source: branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/LBC/lib_mpp.F90 @ 3432

Last change on this file since 3432 was 3432, checked in by trackstand2, 12 years ago

Merge branch 'ksection_partition'

  • Property svn:keywords set to Id
File size: 128.6 KB
Line 
1MODULE lib_mpp
2   !!======================================================================
3   !!                       ***  MODULE  lib_mpp  ***
4   !! Ocean numerics:  massively parallel processing library
5   !!=====================================================================
6   !! History :  OPA  !  1994  (M. Guyon, J. Escobar, M. Imbard)  Original code
7   !!            7.0  !  1997  (A.M. Treguier)  SHMEM additions
8   !!            8.0  !  1998  (M. Imbard, J. Escobar, L. Colombet ) SHMEM and MPI
9   !!                 !  1998  (J.M. Molines) Open boundary conditions
10   !!   NEMO     1.0  !  2003  (J.-M. Molines, G. Madec)  F90, free form
11   !!                 !  2003  (J.M. Molines) add mpp_ini_north(_3d,_2d)
12   !!             -   !  2004  (R. Bourdalle Badie)  isend option in mpi
13   !!                 !  2004  (J.M. Molines) minloc, maxloc
14   !!             -   !  2005  (G. Madec, S. Masson)  npolj=5,6 F-point & ice cases
15   !!             -   !  2005  (R. Redler) Replacement of MPI_COMM_WORLD except for MPI_Abort
16   !!             -   !  2005  (R. Benshila, G. Madec)  add extra halo case
17   !!             -   !  2008  (R. Benshila) add mpp_ini_ice
18   !!            3.2  !  2009  (R. Benshila) SHMEM suppression, north fold in lbc_nfd
19   !!            3.2  !  2009  (O. Marti)    add mpp_ini_znl
20   !!            4.0  !  2011  (G. Madec)  move ctl_ routines from in_out_manager
21   !!----------------------------------------------------------------------
22
23   !!----------------------------------------------------------------------
24   !!   ctl_stop   : update momentum and tracer Kz from a tke scheme
25   !!   ctl_warn   : initialization, namelist read, and parameters control
26   !!   ctl_opn    : Open file and check if required file is available.
27   !!   get_unit    : give the index of an unused logical unit
28   !!----------------------------------------------------------------------
29#if   defined key_mpp_mpi 
30   !!----------------------------------------------------------------------
31   !!   'key_mpp_mpi'             MPI massively parallel processing library
32   !!----------------------------------------------------------------------
33   !!   lib_mpp_alloc : allocate mpp arrays
34   !!   mynode        : indentify the processor unit
35   !!   mpp_lnk       : interface (defined in lbclnk) for message passing of 2d or 3d arrays (mpp_lnk_2d, mpp_lnk_3d)
36   !!   mpp_lnk_3d_gather :  Message passing manadgement for two 3D arrays
37   !!   mpp_lnk_e     : interface (defined in lbclnk) for message passing of 2d array with extra halo (mpp_lnk_2d_e)
38   !!   mpprecv         :
39   !!   mppsend       :   SUBROUTINE mpp_ini_znl
40   !!   mppscatter    :
41   !!   mppgather     :
42   !!   mpp_min       : generic interface for mppmin_int , mppmin_a_int , mppmin_real, mppmin_a_real
43   !!   mpp_max       : generic interface for mppmax_int , mppmax_a_int , mppmax_real, mppmax_a_real
44   !!   mpp_sum       : generic interface for mppsum_int , mppsum_a_int , mppsum_real, mppsum_a_real
45   !!   mpp_minloc    :
46   !!   mpp_maxloc    :
47   !!   mppsync       :
48   !!   mppstop       :
49   !!   mppobc        : variant of mpp_lnk for open boundary condition
50   !!   mpp_ini_north : initialisation of north fold
51   !!   mpp_lbc_north : north fold processors gathering
52   !!   mpp_lbc_north_e : variant of mpp_lbc_north for extra outer halo
53   !!----------------------------------------------------------------------
54   USE dom_oce        ! ocean space and time domain
55   USE lbcnfd         ! north fold treatment
56   USE in_out_manager ! I/O manager
57
58   IMPLICIT NONE
59   PRIVATE
60   
61   PUBLIC   ctl_stop, ctl_warn, get_unit, ctl_opn
62   PUBLIC   mynode, mppstop, mppsync, mpp_comm_free
63   PUBLIC   mpp_ini_north, mpp_lbc_north, mpp_lbc_north_e
64   PUBLIC   mpp_min, mpp_max, mpp_sum, mpp_minloc, mpp_maxloc
65   PUBLIC   mpp_lnk_3d, mpp_lnk_3d_gather, mpp_lnk_2d, mpp_lnk_2d_e
66   PUBLIC   mppobc, mpp_ini_ice, mpp_ini_znl
67   PUBLIC   mppsize, MAX_FACTORS, nxfactors, xfactors, nyfactors, yfactors
68   PUBLIC   lib_mpp_alloc   ! Called in nemogcm.F90
69
70   !! * Interfaces
71   !! define generic interface for these routine as they are called sometimes
72   !! with scalar arguments instead of array arguments, which causes problems
73   !! for the compilation on AIX system as well as NEC and SGI. Ok on COMPACQ
74   INTERFACE mpp_min
75      MODULE PROCEDURE mppmin_a_int, mppmin_int, mppmin_a_real, mppmin_real
76   END INTERFACE
77   INTERFACE mpp_max
78      MODULE PROCEDURE mppmax_a_int, mppmax_int, mppmax_a_real, mppmax_real
79   END INTERFACE
80   INTERFACE mpp_sum
81# if defined key_mpp_rep
82      MODULE PROCEDURE mppsum_a_int, mppsum_int, mppsum_a_real, mppsum_real, &
83                       mppsum_realdd, mppsum_a_realdd
84# else
85      MODULE PROCEDURE mppsum_a_int, mppsum_int, mppsum_a_real, mppsum_real
86# endif
87   END INTERFACE
88   INTERFACE mpp_lbc_north
89      MODULE PROCEDURE mpp_lbc_north_3d, mpp_lbc_north_2d 
90   END INTERFACE
91   INTERFACE mpp_minloc
92      MODULE PROCEDURE mpp_minloc2d ,mpp_minloc3d
93   END INTERFACE
94   INTERFACE mpp_maxloc
95      MODULE PROCEDURE mpp_maxloc2d ,mpp_maxloc3d
96   END INTERFACE
97   
98   !! ========================= !!
99   !!  MPI  variable definition !!
100   !! ========================= !!
101!$AGRIF_DO_NOT_TREAT
102   INCLUDE 'mpif.h'
103!$AGRIF_END_DO_NOT_TREAT
104   
105   LOGICAL, PUBLIC, PARAMETER ::   lk_mpp = .TRUE.    !: mpp flag
106
107   INTEGER, PARAMETER         ::   nprocmax = 2**10   ! maximun dimension (required to be a power of 2)
108   
109   INTEGER ::   mppsize        ! number of process
110   INTEGER ::   mpprank        ! process number  [ 0 - size-1 ]
111!$AGRIF_DO_NOT_TREAT
112   INTEGER, PUBLIC ::   mpi_comm_opa   ! opa local communicator
113!$AGRIF_END_DO_NOT_TREAT
114
115# if defined key_mpp_rep
116   INTEGER :: MPI_SUMDD
117# endif
118
119   ! variables used in case of sea-ice
120   INTEGER, PUBLIC ::   ncomm_ice       !: communicator made by the processors with sea-ice
121   INTEGER ::   ngrp_ice        !  group ID for the ice processors (for rheology)
122   INTEGER ::   ndim_rank_ice   !  number of 'ice' processors
123   INTEGER ::   n_ice_root      !  number (in the comm_ice) of proc 0 in the ice comm
124   INTEGER, DIMENSION(:), ALLOCATABLE, SAVE ::   nrank_ice     ! dimension ndim_rank_ice
125
126   ! variables used for zonal integration
127   INTEGER, PUBLIC ::   ncomm_znl       !: communicator made by the processors on the same zonal average
128   LOGICAL, PUBLIC ::   l_znl_root      ! True on the 'left'most processor on the same row
129   INTEGER ::   ngrp_znl        ! group ID for the znl processors
130   INTEGER ::   ndim_rank_znl   ! number of processors on the same zonal average
131   INTEGER, DIMENSION(:), ALLOCATABLE, SAVE ::   nrank_znl  ! dimension ndim_rank_znl,
132   !                                                        ! number of the procs into the same znl domain
133   
134   ! North fold condition in mpp_mpi with jpni > 1
135   INTEGER ::   ngrp_world        ! group ID for the world processors
136   INTEGER ::   ngrp_opa          ! group ID for the opa processors
137   INTEGER ::   ngrp_north        ! group ID for the northern processors (to be fold)
138   INTEGER ::   ncomm_north       ! communicator made by the processors belonging to ngrp_north
139   INTEGER ::   ndim_rank_north   ! number of 'sea' processor in the northern line (can be /= jpni !)
140   INTEGER ::   njmppmax          ! value of njmpp for the processors of the northern line
141   INTEGER ::   north_root        ! number (in the comm_opa) of proc 0 in the northern comm
142   INTEGER, DIMENSION(:), ALLOCATABLE, SAVE ::   nrank_north   ! dimension ndim_rank_north
143
144   ! Type of send : standard, buffered, immediate
145   CHARACTER(len=1) ::   cn_mpi_send = 'S'    ! type od mpi send/recieve (S=standard, B=bsend, I=isend)
146   LOGICAL          ::   l_isend = .FALSE.   ! isend use indicator (T if cn_mpi_send='I')
147   INTEGER          ::   nn_buffer = 0       ! size of the buffer in case of mpi_bsend
148   CHARACTER(len=256) :: nn_xfactors = ''    ! String holding factors to use for PE grid in x direction     
149   CHARACTER(len=256) :: nn_yfactors = ''    ! String holding factors to use for PE grid in y direction     
150   INTEGER, PARAMETER :: MAX_FACTORS = 20    ! Maximum no. of factors factor() can return
151   ! Arrays to hold specific factorisation of the processor grid specified
152   ! in the namelist. Set to -1 if no specific factorisation requested.
153   INTEGER, SAVE, DIMENSION(MAX_FACTORS) :: xfactors, yfactors
154   INTEGER, SAVE                         :: nxfactors, nyfactors
155   LOGICAL,       PUBLIC                 :: nn_pttrim = .FALSE. ! Whether to trim dry
156                                                                ! land from PE domains
157   INTEGER, SAVE, PUBLIC                 :: nn_cpnode = 4 ! Number of cores per
158                                                          ! compute node on current computer
159 
160   REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE :: tampon  ! buffer in case of bsend
161
162   ! message passing arrays
163   REAL(wp), DIMENSION(:,:,:,:,:), ALLOCATABLE, SAVE ::   t4ns, t4sn   ! 2 x 3d for north-south & south-north
164   REAL(wp), DIMENSION(:,:,:,:,:), ALLOCATABLE, SAVE ::   t4ew, t4we   ! 2 x 3d for east-west & west-east
165   REAL(wp), DIMENSION(:,:,:,:,:), ALLOCATABLE, SAVE ::   t4p1, t4p2   ! 2 x 3d for north fold
166   REAL(wp), DIMENSION(:,:,:,:)  , ALLOCATABLE, SAVE ::   t3ns, t3sn   ! 3d for north-south & south-north
167   REAL(wp), DIMENSION(:,:,:,:)  , ALLOCATABLE, SAVE ::   t3ew, t3we   ! 3d for east-west & west-east
168   REAL(wp), DIMENSION(:,:,:,:)  , ALLOCATABLE, SAVE ::   t3p1, t3p2   ! 3d for north fold
169   REAL(wp), DIMENSION(:,:,:)    , ALLOCATABLE, SAVE ::   t2ns, t2sn   ! 2d for north-south & south-north
170   REAL(wp), DIMENSION(:,:,:)    , ALLOCATABLE, SAVE ::   t2ew, t2we   ! 2d for east-west & west-east
171   REAL(wp), DIMENSION(:,:,:)    , ALLOCATABLE, SAVE ::   t2p1, t2p2   ! 2d for north fold
172   REAL(wp), DIMENSION(:,:,:)    , ALLOCATABLE, SAVE ::   tr2ns, tr2sn ! 2d for north-south & south-north + extra outer halo
173   REAL(wp), DIMENSION(:,:,:)    , ALLOCATABLE, SAVE ::   tr2ew, tr2we ! 2d for east-west   & west-east   + extra outer halo
174
175   ! Arrays used in mpp_lbc_north_3d()
176   REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE, SAVE   ::   ztab, znorthloc
177   REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE, SAVE   ::   znorthgloio
178
179   ! Arrays used in mpp_lbc_north_2d()
180   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, SAVE    ::   ztab_2d, znorthloc_2d
181   REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, SAVE    ::   znorthgloio_2d
182
183   ! Arrays used in mpp_lbc_north_e()
184   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, SAVE    ::   ztab_e, znorthloc_e
185   REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, SAVE    ::   znorthgloio_e
186
187   !! * Control permutation of array indices
188#  include "dom_oce_ftrans.h90"
189!! These arrays are all private to the module
190!FTRANS t4ns :I :I :z :I :I
191!FTRANS t4sn :I :I :z :I :I
192!FTRANS t4ew :I :I :z :I :I
193!FTRANS t4we :I :I :z :I :I
194!FTRANS t3ns :I :I :z :I
195!FTRANS t3sn :I :I :z :I
196!FTRANS t3ew :I :I :z :I
197!FTRANS t3we :I :I :z :I
198!FTRANS ztab :I :I :z
199!FTRANS znorthloc :I :I :z
200!FTRANS znorthgloio :I :I :z :I
201
202   !!----------------------------------------------------------------------
203   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
204   !! $Id$
205   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
206   !!----------------------------------------------------------------------
207CONTAINS
208
209   INTEGER FUNCTION lib_mpp_alloc( kumout )
210      !!----------------------------------------------------------------------
211      !!              ***  routine lib_mpp_alloc  ***
212      !!----------------------------------------------------------------------
213      INTEGER, INTENT(in) ::   kumout   ! ocean.output logical unit
214      !!----------------------------------------------------------------------
215      !
216      ALLOCATE( &
217#if !defined key_mpp_rkpart
218                t4ns(jpi,jprecj,jpk,2,2) , t4sn(jpi,jprecj,jpk,2,2) ,                                            &
219         &      t4ew(jpj,jpreci,jpk,2,2) , t4we(jpj,jpreci,jpk,2,2) ,                                            &
220         &      t4p1(jpi,jprecj,jpk,2,2) , t4p2(jpi,jprecj,jpk,2,2) ,                                            &
221         &      t3ns(jpi,jprecj,jpk,2)   , t3sn(jpi,jprecj,jpk,2)   ,                                            &
222         &      t3ew(jpj,jpreci,jpk,2)   , t3we(jpj,jpreci,jpk,2)   ,                                            &
223         &      t3p1(jpi,jprecj,jpk,2)   , t3p2(jpi,jprecj,jpk,2)   ,                                            &
224         &      t2ns(jpi,jprecj    ,2)   , t2sn(jpi,jprecj    ,2)   ,                                            &
225         &      t2ew(jpj,jpreci    ,2)   , t2we(jpj,jpreci    ,2)   ,                                            &
226         &      t2p1(jpi,jprecj    ,2)   , t2p2(jpi,jprecj    ,2)   ,                                            &
227#endif
228         !
229         &      tr2ns(1-jpr2di:jpi+jpr2di,jprecj+jpr2dj,2) ,                                                     &
230         &      tr2sn(1-jpr2di:jpi+jpr2di,jprecj+jpr2dj,2) ,                                                     &
231         &      tr2ew(1-jpr2dj:jpj+jpr2dj,jpreci+jpr2di,2) ,                                                     &
232         &      tr2we(1-jpr2dj:jpj+jpr2dj,jpreci+jpr2di,2) ,                                                     &
233         !
234         &      ztab(jpiglo,4,jpk) , znorthloc(jpi,4,jpk) , znorthgloio(jpi,4,jpk,jpni) ,                        &
235         !
236         &      ztab_2d(jpiglo,4)  , znorthloc_2d(jpi,4)  , znorthgloio_2d(jpi,4,jpni)  ,                        &
237         !
238         &      ztab_e(jpiglo,4+2*jpr2dj) , znorthloc_e(jpi,4+2*jpr2dj) , znorthgloio_e(jpi,4+2*jpr2dj,jpni) ,   &
239         !
240         &      STAT=lib_mpp_alloc )
241         !
242      IF( lib_mpp_alloc /= 0 ) THEN
243         WRITE(kumout,cform_war)
244         WRITE(kumout,*) 'lib_mpp_alloc : failed to allocate arrays'
245      ENDIF
246      !
247   END FUNCTION lib_mpp_alloc
248
249
250   FUNCTION mynode( ldtxt, kumnam, kstop, localComm )
251      !!----------------------------------------------------------------------
252      !!                  ***  routine mynode  ***
253      !!                   
254      !! ** Purpose :   Find processor unit
255      !!----------------------------------------------------------------------
256      CHARACTER(len=*),DIMENSION(:), INTENT(  out) ::   ldtxt 
257      INTEGER                      , INTENT(in   ) ::   kumnam       ! namelist logical unit
258      INTEGER                      , INTENT(inout) ::   kstop        ! stop indicator
259      INTEGER, OPTIONAL            , INTENT(in   ) ::   localComm
260      !
261      INTEGER ::   mynode, ierr, code, ji, ii
262      LOGICAL ::   mpi_was_called
263      !
264      NAMELIST/nammpp/ cn_mpi_send, nn_buffer, jpni, jpnj, jpnij, &
265                       nn_xfactors, nn_yfactors, nn_pttrim, nn_cpnode
266      !!----------------------------------------------------------------------
267      !
268      ii = 1
269      WRITE(ldtxt(ii),*)                                                                          ;   ii = ii + 1
270      WRITE(ldtxt(ii),*) 'mynode : mpi initialisation'                                            ;   ii = ii + 1
271      WRITE(ldtxt(ii),*) '~~~~~~ '                                                                ;   ii = ii + 1
272      !
273      jpni = -1; jpnj = -1; jpnij = -1
274      REWIND( kumnam )               ! Namelist namrun : parameters of the run
275      READ  ( kumnam, nammpp )
276      !                              ! control print
277      WRITE(ldtxt(ii),*) '   Namelist nammpp'                                                     ;   ii = ii + 1
278      WRITE(ldtxt(ii),*) '      mpi send type                      cn_mpi_send = ', cn_mpi_send   ;   ii = ii + 1
279      WRITE(ldtxt(ii),*) '      size in bytes of exported buffer   nn_buffer   = ', nn_buffer     ;   ii = ii + 1
280      WRITE(ldtxt(ii),*) '      whether to trim dry points         nn_pttrim   = ', nn_pttrim     ;   ii = ii + 1
281      WRITE(ldtxt(ii),*) '      number of cores per compute node   nn_cpn      = ', nn_cpnode     ;   ii = ii + 1
282#if defined key_agrif
283      IF( .NOT. Agrif_Root() ) THEN
284         jpni  = Agrif_Parent(jpni ) 
285         jpnj  = Agrif_Parent(jpnj )
286         jpnij = Agrif_Parent(jpnij)
287      ENDIF
288#endif
289
290      IF(jpnij < 1)THEN
291         ! If jpnij is not specified in namelist then we calculate it - this
292         ! means there will be no land cutting out.
293         jpnij = jpni * jpnj
294      END IF
295
296      IF( (jpni < 1) .OR. (jpnj < 1) )THEN
297         WRITE(ldtxt(ii),*) '      jpni, jpnj and jpnij will be calculated automatically'; ii = ii + 1
298      ELSE
299         WRITE(ldtxt(ii),*) '      processor grid extent in i         jpni = ',jpni; ii = ii + 1
300         WRITE(ldtxt(ii),*) '      processor grid extent in j         jpnj = ',jpnj; ii = ii + 1
301         WRITE(ldtxt(ii),*) '      number of local domains           jpnij = ',jpnij; ii = ii +1
302      END IF
303
304      ! Check to see whether a specific factorisation of the number of
305      ! processors has been specified in the namelist file
306      nxfactors = 0; xfactors(:) = -1
307      nyfactors = 0; yfactors(:) = -1
308      IF((LEN_TRIM(nn_xfactors) > 0) .OR. (LEN_TRIM(nn_yfactors) > 0))THEN
309
310         IF( (VERIFY(TRIM(nn_xfactors),'0123456789,') > 0) .OR. &
311             (VERIFY(TRIM(nn_yfactors),'0123456789,') > 0) )THEN
312            WRITE(ldtxt(ii),*)'Invalid character in nn_xfactors/nn_yfactors namelist string. '; ii = ii + 1
313            WRITE(ldtxt(ii),*)'Will ignore requested factorisation.' ; ii = ii + 1
314         ELSE
315            READ(nn_xfactors, *,end=80,err=100) xfactors
316         80 CONTINUE
317            READ(nn_yfactors, *,end=90,err=100) yfactors
318         90 CONTINUE
319
320            trim_xarray: DO ji=MAX_FACTORS,1,-1
321               IF (xfactors(ji) .GE. 0) THEN
322                  nxfactors = ji
323                  EXIT trim_xarray
324               ENDIF
325            ENDDO trim_xarray
326            trim_yarray: DO ji=MAX_FACTORS,1,-1
327               IF (yfactors(ji) .GE. 0) THEN
328                  nyfactors = ji
329                  EXIT trim_yarray
330               ENDIF
331            ENDDO trim_yarray
332
333        100 CONTINUE
334            WRITE (*,*) 'ARPDBG: n{x,y}factors = ',nxfactors,nyfactors
335            IF(nxfactors < 1 .AND. nyfactors < 1)THEN
336               WRITE(ldtxt(ii),*)'Failed to parse factorisation string'    ; ii = ii + 1
337               WRITE(ldtxt(ii),*)' - will ignore requested factorisation.' ; ii = ii + 1
338            ELSE
339               WRITE(ldtxt(ii),*)'      automatic factorisation overridden'   ; ii = ii + 1
340               WRITE(ldtxt(ii),*)'      factors:', xfactors(1:nxfactors), &
341                                 '-',yfactors(1:nyfactors)
342               ii = ii + 1
343            END IF
344         ENDIF
345
346      END IF
347
348      CALL mpi_initialized ( mpi_was_called, code )
349      IF( code /= MPI_SUCCESS ) THEN
350         DO ji = 1, SIZE(ldtxt) 
351            IF( TRIM(ldtxt(ji)) /= '' )   WRITE(*,*) ldtxt(ji)      ! control print of mynode
352         END DO         
353         WRITE(*, cform_err)
354         WRITE(*, *) 'lib_mpp: Error in routine mpi_initialized'
355         CALL mpi_abort( mpi_comm_world, code, ierr )
356      ENDIF
357
358      IF( mpi_was_called ) THEN
359         !
360         SELECT CASE ( cn_mpi_send )
361         CASE ( 'S' )                ! Standard mpi send (blocking)
362            WRITE(ldtxt(ii),*) '           Standard blocking mpi send (send)'                     ;   ii = ii + 1
363         CASE ( 'B' )                ! Buffer mpi send (blocking)
364            WRITE(ldtxt(ii),*) '           Buffer blocking mpi send (bsend)'                      ;   ii = ii + 1
365            IF( Agrif_Root() )   CALL mpi_init_opa( ldtxt, ii, ierr ) 
366         CASE ( 'I' )                ! Immediate mpi send (non-blocking send)
367            WRITE(ldtxt(ii),*) '           Immediate non-blocking send (isend)'                   ;   ii = ii + 1
368            l_isend = .TRUE.
369         CASE DEFAULT
370            WRITE(ldtxt(ii),cform_err)                                                            ;   ii = ii + 1
371            WRITE(ldtxt(ii),*) '           bad value for cn_mpi_send = ', cn_mpi_send             ;   ii = ii + 1
372            kstop = kstop + 1
373         END SELECT
374      ELSE IF ( PRESENT(localComm) .and. .not. mpi_was_called ) THEN
375         WRITE(ldtxt(ii),*) ' lib_mpp: You cannot provide a local communicator '                  ;   ii = ii + 1
376         WRITE(ldtxt(ii),*) '          without calling MPI_Init before ! '                        ;   ii = ii + 1
377         kstop = kstop + 1
378      ELSE
379         SELECT CASE ( cn_mpi_send )
380         CASE ( 'S' )                ! Standard mpi send (blocking)
381            WRITE(ldtxt(ii),*) '           Standard blocking mpi send (send)'                     ;   ii = ii + 1
382            CALL mpi_init( ierr )
383         CASE ( 'B' )                ! Buffer mpi send (blocking)
384            WRITE(ldtxt(ii),*) '           Buffer blocking mpi send (bsend)'                      ;   ii = ii + 1
385            IF( Agrif_Root() )   CALL mpi_init_opa( ldtxt, ii, ierr )
386         CASE ( 'I' )                ! Immediate mpi send (non-blocking send)
387            WRITE(ldtxt(ii),*) '           Immediate non-blocking send (isend)'                   ;   ii = ii + 1
388            l_isend = .TRUE.
389            CALL mpi_init( ierr )
390         CASE DEFAULT
391            WRITE(ldtxt(ii),cform_err)                                                            ;   ii = ii + 1
392            WRITE(ldtxt(ii),*) '           bad value for cn_mpi_send = ', cn_mpi_send             ;   ii = ii + 1
393            kstop = kstop + 1
394         END SELECT
395         !
396      ENDIF
397
398      IF( PRESENT(localComm) ) THEN
399         IF( Agrif_Root() ) THEN
400            mpi_comm_opa = localComm
401         ENDIF
402      ELSE
403         CALL mpi_comm_dup( mpi_comm_world, mpi_comm_opa, code)
404         IF( code /= MPI_SUCCESS ) THEN
405            DO ji = 1, SIZE(ldtxt) 
406               IF( TRIM(ldtxt(ji)) /= '' )   WRITE(*,*) ldtxt(ji)      ! control print of mynode
407            END DO
408            WRITE(*, cform_err)
409            WRITE(*, *) ' lib_mpp: Error in routine mpi_comm_dup'
410            CALL mpi_abort( mpi_comm_world, code, ierr )
411         ENDIF
412      ENDIF
413
414      CALL mpi_comm_rank( mpi_comm_opa, mpprank, ierr )
415      CALL mpi_comm_size( mpi_comm_opa, mppsize, ierr )
416      mynode = mpprank
417      !
418      IF( (jpni < 1) .OR. (jpnj < 1) )THEN
419         IF(mynode == 0)WRITE(ldtxt(ii),*) '      Running on ',mppsize,' MPI processes'; ii = ii + 1
420      END IF
421
422#if defined key_mpp_rep
423      CALL MPI_OP_CREATE(DDPDD_MPI, .TRUE., MPI_SUMDD, ierr)
424#endif
425      !
426   END FUNCTION mynode
427
428
429   SUBROUTINE mpp_lnk_3d( ptab3d, cd_type, psgn, cd_mpp, pval, lzero )
430      !!----------------------------------------------------------------------
431      !!                  ***  routine mpp_lnk_3d  ***
432      !!
433      !! ** Purpose :   Message passing manadgement
434      !!
435      !! ** Method  :   Use mppsend and mpprecv function for passing mask
436      !!      between processors following neighboring subdomains.
437      !!            domain parameters
438      !!                    nlci   : first dimension of the local subdomain
439      !!                    nlcj   : second dimension of the local subdomain
440      !!                    nbondi : mark for "east-west local boundary"
441      !!                    nbondj : mark for "north-south local boundary"
442      !!                    noea   : number for local neighboring processors
443      !!                    nowe   : number for local neighboring processors
444      !!                    noso   : number for local neighboring processors
445      !!                    nono   : number for local neighboring processors
446      !!
447      !! ** Action  :   ptab with update value at its periphery
448      !!
449      !!----------------------------------------------------------------------
450!FTRANS ptab3d :I :I :z
451!! DCSE_NEMO: work around a deficiency in ftrans
452!     REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptab3d   ! 3D array on which the boundary condition is applied
453      REAL(wp),                         INTENT(inout) ::   ptab3d(jpi,jpj,jpk)
454      CHARACTER(len=1)                , INTENT(in   ) ::   cd_type  ! define the nature of ptab array grid-points
455      !                                                             ! = T , U , V , F , W points
456      REAL(wp)                        , INTENT(in   ) ::   psgn     ! =-1 the sign change across the north fold boundary
457      !                                                             ! =  1. , the sign is kept
458      CHARACTER(len=3), OPTIONAL      , INTENT(in   ) ::   cd_mpp   ! fill the overlap area only
459      REAL(wp)        , OPTIONAL      , INTENT(in   ) ::   pval     ! background value (used at closed boundaries)
460      LOGICAL         , OPTIONAL      , INTENT(in   ) ::   lzero    ! Whether to zero field at closed boundaries
461      !!
462      INTEGER  ::   ji, jj, jk, jl             ! dummy loop indices
463      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
464      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
465      REAL(wp) ::   zland
466      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
467      LOGICAL  ::   lzeroarg
468      !!----------------------------------------------------------------------
469
470#if defined key_mpp_rkpart
471      CALL ctl_stop('mpp_lnk_3d: should not have been called when key_mpp_rkpart defined!')
472      RETURN
473#endif
474      ! Deal with optional routine arguments
475      lzeroarg = .TRUE.
476      IF( PRESENT(lzero) ) lzeroarg = lzero
477
478      IF( PRESENT( pval ) ) THEN   ;   zland = pval      ! set land value
479      ELSE                         ;   zland = 0.e0      ! zero by default
480      ENDIF
481
482      ! 1. standard boundary treatment
483      ! ------------------------------
484      IF( PRESENT( cd_mpp ) ) THEN      ! only fill added line/raw with existing values
485         !
486         ! WARNING ptab3d is defined only between nld and nle
487#if defined key_z_first
488         DO jj = nlcj+1, jpj                    ! added line(s)   (inner only)
489            DO jk = 1, jpk
490#else
491         DO jk = 1, jpk
492            DO jj = nlcj+1, jpj                 ! added line(s)   (inner only)
493#endif
494               ptab3d(nldi  :nlei  , jj          ,jk) = ptab3d(nldi:nlei,     nlej,jk)   
495               ptab3d(1     :nldi-1, jj          ,jk) = ptab3d(nldi     ,     nlej,jk)
496               ptab3d(nlei+1:nlci  , jj          ,jk) = ptab3d(     nlei,     nlej,jk)
497            END DO
498#if defined key_z_first
499         END DO
500         DO ji = nlci+1, jpi                    ! added column(s) (full)
501            DO jk = 1, jpk
502#else
503            DO ji = nlci+1, jpi                 ! added column(s) (full)
504#endif
505               ptab3d(ji           ,nldj  :nlej  ,jk) = ptab3d(     nlei,nldj:nlej,jk)
506               ptab3d(ji           ,1     :nldj-1,jk) = ptab3d(     nlei,nldj     ,jk)
507               ptab3d(ji           ,nlej+1:jpj   ,jk) = ptab3d(     nlei,     nlej,jk)
508            END DO
509         END DO
510         !
511      ELSE                              ! standard close or cyclic treatment
512         !
513         !                                   ! East-West boundaries
514         !                                        !* Cyclic east-west
515         IF( nbondi == 2 .AND. (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
516            ptab3d( 1 ,:,:) = ptab3d(jpim1,:,:)
517            ptab3d(jpi,:,:) = ptab3d(  2  ,:,:)
518         ELSE                                     !* closed
519            IF( lzeroarg )THEN
520               IF( .NOT. cd_type == 'F' )   ptab3d(     1       :jpreci,:,:) = zland    ! south except F-point
521                                            ptab3d(nlci-jpreci+1:jpi   ,:,:) = zland    ! north
522            END IF
523         ENDIF
524         !                                   ! North-South boundaries (always closed)
525         IF( lzeroarg )THEN
526            IF( .NOT. cd_type == 'F' )   ptab3d(:,     1       :jprecj,:) = zland       ! south except F-point
527                                         ptab3d(:,nlcj-jprecj+1:jpj   ,:) = zland       ! north
528         END IF
529         !
530      ENDIF
531
532      ! 2. East and west directions exchange
533      ! ------------------------------------
534      ! we play with the neigbours AND the row number because of the periodicity
535      !
536      SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
537      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
538         iihom = nlci-nreci
539         DO jl = 1, jpreci
540            t3ew(:,jl,:,1) = ptab3d(jpreci+jl,:,:)
541            t3we(:,jl,:,1) = ptab3d(iihom +jl,:,:)
542         END DO
543      END SELECT 
544      !
545      !                           ! Migrations
546      imigr = jpreci * jpj * jpk
547      !
548      SELECT CASE ( nbondi ) 
549      CASE ( -1 )
550         CALL mppsend( 2, t3we(1,1,1,1), imigr, noea, ml_req1 )
551         CALL mpprecv( 1, t3ew(1,1,1,2), imigr )
552         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
553      CASE ( 0 )
554         CALL mppsend( 1, t3ew(1,1,1,1), imigr, nowe, ml_req1 )
555         CALL mppsend( 2, t3we(1,1,1,1), imigr, noea, ml_req2 )
556         CALL mpprecv( 1, t3ew(1,1,1,2), imigr )
557         CALL mpprecv( 2, t3we(1,1,1,2), imigr )
558         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
559         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
560      CASE ( 1 )
561         CALL mppsend( 1, t3ew(1,1,1,1), imigr, nowe, ml_req1 )
562         CALL mpprecv( 2, t3we(1,1,1,2), imigr )
563         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
564      END SELECT
565      !
566      !                           ! Write Dirichlet lateral conditions
567      iihom = nlci-jpreci
568      !
569      SELECT CASE ( nbondi )
570      CASE ( -1 )
571         DO jl = 1, jpreci
572            ptab3d(iihom+jl,:,:) = t3ew(:,jl,:,2)
573         END DO
574      CASE ( 0 ) 
575         DO jl = 1, jpreci
576            ptab3d(jl      ,:,:) = t3we(:,jl,:,2)
577            ptab3d(iihom+jl,:,:) = t3ew(:,jl,:,2)
578         END DO
579      CASE ( 1 )
580         DO jl = 1, jpreci
581            ptab3d(jl      ,:,:) = t3we(:,jl,:,2)
582         END DO
583      END SELECT
584
585
586      ! 3. North and south directions
587      ! -----------------------------
588      ! always closed : we play only with the neigbours
589      !
590      IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
591         ijhom = nlcj-nrecj
592         DO jl = 1, jprecj
593            t3sn(:,jl,:,1) = ptab3d(:,ijhom +jl,:)
594            t3ns(:,jl,:,1) = ptab3d(:,jprecj+jl,:)
595         END DO
596      ENDIF
597      !
598      !                           ! Migrations
599      imigr = jprecj * jpi * jpk
600      !
601      SELECT CASE ( nbondj )     
602      CASE ( -1 )
603         CALL mppsend( 4, t3sn(1,1,1,1), imigr, nono, ml_req1 )
604         CALL mpprecv( 3, t3ns(1,1,1,2), imigr )
605         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
606      CASE ( 0 )
607         CALL mppsend( 3, t3ns(1,1,1,1), imigr, noso, ml_req1 )
608         CALL mppsend( 4, t3sn(1,1,1,1), imigr, nono, ml_req2 )
609         CALL mpprecv( 3, t3ns(1,1,1,2), imigr )
610         CALL mpprecv( 4, t3sn(1,1,1,2), imigr )
611         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
612         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
613      CASE ( 1 ) 
614         CALL mppsend( 3, t3ns(1,1,1,1), imigr, noso, ml_req1 )
615         CALL mpprecv( 4, t3sn(1,1,1,2), imigr )
616         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
617      END SELECT
618      !
619      !                           ! Write Dirichlet lateral conditions
620      ijhom = nlcj-jprecj
621      !
622      SELECT CASE ( nbondj )
623      CASE ( -1 )
624         DO jl = 1, jprecj
625            ptab3d(:,ijhom+jl,:) = t3ns(:,jl,:,2)
626         END DO
627      CASE ( 0 ) 
628         DO jl = 1, jprecj
629            ptab3d(:,jl      ,:) = t3sn(:,jl,:,2)
630            ptab3d(:,ijhom+jl,:) = t3ns(:,jl,:,2)
631         END DO
632      CASE ( 1 )
633         DO jl = 1, jprecj
634            ptab3d(:,jl,:) = t3sn(:,jl,:,2)
635         END DO
636      END SELECT
637
638
639      ! 4. north fold treatment
640      ! -----------------------
641      !
642      IF( npolj /= 0 .AND. .NOT. PRESENT(cd_mpp) ) THEN
643         !
644         SELECT CASE ( jpni )
645         CASE ( 1 )     ;   CALL lbc_nfd      ( ptab3d, cd_type, psgn )   ! only 1 northern proc, no mpp
646         CASE DEFAULT   ;   CALL mpp_lbc_north( ptab3d, cd_type, psgn )   ! for all northern procs.
647         END SELECT
648         !
649      ENDIF
650      !
651   END SUBROUTINE mpp_lnk_3d
652
653
654   SUBROUTINE mpp_lnk_2d( pt2d, cd_type, psgn, cd_mpp, pval, lzero )
655      !!----------------------------------------------------------------------
656      !!                  ***  routine mpp_lnk_2d  ***
657      !!                 
658      !! ** Purpose :   Message passing manadgement for 2d array
659      !!
660      !! ** Method  :   Use mppsend and mpprecv function for passing mask
661      !!      between processors following neighboring subdomains.
662      !!            domain parameters
663      !!                    nlci   : first dimension of the local subdomain
664      !!                    nlcj   : second dimension of the local subdomain
665      !!                    nbondi : mark for "east-west local boundary"
666      !!                    nbondj : mark for "north-south local boundary"
667      !!                    noea   : number for local neighboring processors
668      !!                    nowe   : number for local neighboring processors
669      !!                    noso   : number for local neighboring processors
670      !!                    nono   : number for local neighboring processors
671      !!
672      !!----------------------------------------------------------------------
673      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   pt2d     ! 2D array on which the boundary condition is applied
674      CHARACTER(len=1)            , INTENT(in   ) ::   cd_type  ! define the nature of ptab array grid-points
675      !                                                         ! = T , U , V , F , W and I points
676      REAL(wp)                    , INTENT(in   ) ::   psgn     ! =-1 the sign change across the north fold boundary
677      !                                                         ! =  1. , the sign is kept
678      CHARACTER(len=3), OPTIONAL  , INTENT(in   ) ::   cd_mpp   ! fill the overlap area only
679      REAL(wp)        , OPTIONAL  , INTENT(in   ) ::   pval     ! background value (used at closed boundaries)
680      LOGICAL         , OPTIONAL  , INTENT(in   ) ::   lzero    ! Whether to zero field at closed boundaries
681      !!
682      INTEGER  ::   ji, jj, jl   ! dummy loop indices
683      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
684      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
685      REAL(wp) ::   zland
686      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
687      LOGICAL  ::   lzeroarg
688      !!----------------------------------------------------------------------
689
690#if defined key_mpp_rkpart
691      CALL ctl_stop('mpp_lnk_3d: should not have been called when key_mpp_rkpart defined!')
692      RETURN
693#endif
694
695      ! Deal with optional routine arguments
696      lzeroarg = .TRUE.
697      IF( PRESENT(lzero) ) lzeroarg = lzero
698
699      IF( PRESENT( pval ) ) THEN   ;   zland = pval      ! set land value
700      ELSE                         ;   zland = 0.e0      ! zero by default
701      ENDIF
702
703      ! 1. standard boundary treatment
704      ! ------------------------------
705      !
706      IF( PRESENT( cd_mpp ) ) THEN      ! only fill added line/raw with existing values
707         !
708         ! WARNING pt2d is defined only between nld and nle
709         DO jj = nlcj+1, jpj                 ! added line(s)   (inner only)
710            pt2d(nldi  :nlei  , jj          ) = pt2d(nldi:nlei,     nlej)   
711            pt2d(1     :nldi-1, jj          ) = pt2d(nldi     ,     nlej)
712            pt2d(nlei+1:nlci  , jj          ) = pt2d(     nlei,     nlej)
713         END DO
714         DO ji = nlci+1, jpi                 ! added column(s) (full)
715            pt2d(ji           ,nldj  :nlej  ) = pt2d(     nlei,nldj:nlej)
716            pt2d(ji           ,1     :nldj-1) = pt2d(     nlei,nldj     )
717            pt2d(ji           ,nlej+1:jpj   ) = pt2d(     nlei,     nlej)
718         END DO
719         !
720      ELSE                              ! standard close or cyclic treatment
721         !
722         !                                   ! East-West boundaries
723         IF( nbondi == 2 .AND.   &                ! Cyclic east-west
724            &    (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
725            pt2d( 1 ,:) = pt2d(jpim1,:)                                    ! west
726            pt2d(jpi,:) = pt2d(  2  ,:)                                    ! east
727         ELSE                                     ! closed
728            IF( lzeroarg )THEN
729               IF( .NOT. cd_type == 'F' )   pt2d(     1       :jpreci,:) = zland    ! south except F-point
730                                            pt2d(nlci-jpreci+1:jpi   ,:) = zland    ! north
731            END IF
732         ENDIF
733         !                                   ! North-South boundaries (always closed)
734            IF( lzeroarg )THEN
735               IF( .NOT. cd_type == 'F' )   pt2d(:,     1       :jprecj) = zland    !south except F-point
736                                            pt2d(:,nlcj-jprecj+1:jpj   ) = zland    ! north
737            END IF
738         !
739      ENDIF
740
741      ! 2. East and west directions exchange
742      ! ------------------------------------
743      ! we play with the neigbours AND the row number because of the periodicity
744      !
745      SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
746      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
747         iihom = nlci-nreci
748         DO jl = 1, jpreci
749            t2ew(:,jl,1) = pt2d(jpreci+jl,:)
750            t2we(:,jl,1) = pt2d(iihom +jl,:)
751         END DO
752      END SELECT
753      !
754      !                           ! Migrations
755      imigr = jpreci * jpj
756      !
757      SELECT CASE ( nbondi )
758      CASE ( -1 )
759         CALL mppsend( 2, t2we(1,1,1), imigr, noea, ml_req1 )
760         CALL mpprecv( 1, t2ew(1,1,2), imigr )
761         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
762      CASE ( 0 )
763         CALL mppsend( 1, t2ew(1,1,1), imigr, nowe, ml_req1 )
764         CALL mppsend( 2, t2we(1,1,1), imigr, noea, ml_req2 )
765         CALL mpprecv( 1, t2ew(1,1,2), imigr )
766         CALL mpprecv( 2, t2we(1,1,2), imigr )
767         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
768         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
769      CASE ( 1 )
770         CALL mppsend( 1, t2ew(1,1,1), imigr, nowe, ml_req1 )
771         CALL mpprecv( 2, t2we(1,1,2), imigr )
772         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
773      END SELECT
774      !
775      !                           ! Write Dirichlet lateral conditions
776      iihom = nlci - jpreci
777      !
778      SELECT CASE ( nbondi )
779      CASE ( -1 )
780         DO jl = 1, jpreci
781            pt2d(iihom+jl,:) = t2ew(:,jl,2)
782         END DO
783      CASE ( 0 )
784         DO jl = 1, jpreci
785            pt2d(jl      ,:) = t2we(:,jl,2)
786            pt2d(iihom+jl,:) = t2ew(:,jl,2)
787         END DO
788      CASE ( 1 )
789         DO jl = 1, jpreci
790            pt2d(jl      ,:) = t2we(:,jl,2)
791         END DO
792      END SELECT
793
794
795      ! 3. North and south directions
796      ! -----------------------------
797      ! always closed : we play only with the neigbours
798      !
799      IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
800         ijhom = nlcj-nrecj
801         DO jl = 1, jprecj
802            t2sn(:,jl,1) = pt2d(:,ijhom +jl)
803            t2ns(:,jl,1) = pt2d(:,jprecj+jl)
804         END DO
805      ENDIF
806      !
807      !                           ! Migrations
808      imigr = jprecj * jpi
809      !
810      SELECT CASE ( nbondj )
811      CASE ( -1 )
812         CALL mppsend( 4, t2sn(1,1,1), imigr, nono, ml_req1 )
813         CALL mpprecv( 3, t2ns(1,1,2), imigr )
814         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
815      CASE ( 0 )
816         CALL mppsend( 3, t2ns(1,1,1), imigr, noso, ml_req1 )
817         CALL mppsend( 4, t2sn(1,1,1), imigr, nono, ml_req2 )
818         CALL mpprecv( 3, t2ns(1,1,2), imigr )
819         CALL mpprecv( 4, t2sn(1,1,2), imigr )
820         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
821         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
822      CASE ( 1 )
823         CALL mppsend( 3, t2ns(1,1,1), imigr, noso, ml_req1 )
824         CALL mpprecv( 4, t2sn(1,1,2), imigr )
825         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
826      END SELECT
827      !
828      !                           ! Write Dirichlet lateral conditions
829      ijhom = nlcj - jprecj
830      !
831      SELECT CASE ( nbondj )
832      CASE ( -1 )
833         DO jl = 1, jprecj
834            pt2d(:,ijhom+jl) = t2ns(:,jl,2)
835         END DO
836      CASE ( 0 )
837         DO jl = 1, jprecj
838            pt2d(:,jl      ) = t2sn(:,jl,2)
839            pt2d(:,ijhom+jl) = t2ns(:,jl,2)
840         END DO
841      CASE ( 1 ) 
842         DO jl = 1, jprecj
843            pt2d(:,jl      ) = t2sn(:,jl,2)
844         END DO
845      END SELECT
846
847
848      ! 4. north fold treatment
849      ! -----------------------
850      !
851      IF( npolj /= 0 .AND. .NOT. PRESENT(cd_mpp) ) THEN
852         !
853         SELECT CASE ( jpni )
854         CASE ( 1 )     ;   CALL lbc_nfd      ( pt2d, cd_type, psgn )   ! only 1 northern proc, no mpp
855         CASE DEFAULT   ;   CALL mpp_lbc_north( pt2d, cd_type, psgn )   ! for all northern procs.
856         END SELECT
857         !
858      ENDIF
859      !
860   END SUBROUTINE mpp_lnk_2d
861
862
863   SUBROUTINE mpp_lnk_3d_gather( ptab1, cd_type1, ptab2, cd_type2, psgn )
864      !!----------------------------------------------------------------------
865      !!                  ***  routine mpp_lnk_3d_gather  ***
866      !!
867      !! ** Purpose :   Message passing management for two 3D arrays
868      !!
869      !! ** Method  :   Use mppsend and mpprecv function for passing mask
870      !!      between processors following neighboring subdomains.
871      !!            domain parameters
872      !!                    nlci   : first dimension of the local subdomain
873      !!                    nlcj   : second dimension of the local subdomain
874      !!                    nbondi : mark for "east-west local boundary"
875      !!                    nbondj : mark for "north-south local boundary"
876      !!                    noea   : number for local neighboring processors
877      !!                    nowe   : number for local neighboring processors
878      !!                    noso   : number for local neighboring processors
879      !!                    nono   : number for local neighboring processors
880      !!
881      !! ** Action  :   ptab1 and ptab2  with update value at its periphery
882      !!
883      !!----------------------------------------------------------------------
884!FTRANS ptab1 :I :I :z
885!FTRANS ptab2 :I :I :z
886!! DCSE_NEMO: work around a deficiency in ftrans
887!     REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptab1     ! first and second 3D array on which
888!     REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptab2     ! the boundary condition is applied
889      REAL(wp),                         INTENT(inout) ::   ptab1(jpi,jpj,jpk)
890      REAL(wp),                         INTENT(inout) ::   ptab2(jpi,jpj,jpk)
891      CHARACTER(len=1)                , INTENT(in   ) ::   cd_type1  ! nature of ptab1 and ptab2 arrays
892      CHARACTER(len=1)                , INTENT(in   ) ::   cd_type2  ! i.e. grid-points = T , U , V , F or W points
893      REAL(wp)                        , INTENT(in   ) ::   psgn      ! =-1 the sign change across the north fold boundary
894      !!                                                             ! =  1. , the sign is kept
895      INTEGER  ::   jl   ! dummy loop indices
896      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
897      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
898      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
899      !!----------------------------------------------------------------------
900
901      ! 1. standard boundary treatment
902      ! ------------------------------
903      !                                      ! East-West boundaries
904      !                                           !* Cyclic east-west
905      IF( nbondi == 2 .AND. (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
906         ptab1( 1 ,:,:) = ptab1(jpim1,:,:)
907         ptab1(jpi,:,:) = ptab1(  2  ,:,:)
908         ptab2( 1 ,:,:) = ptab2(jpim1,:,:)
909         ptab2(jpi,:,:) = ptab2(  2  ,:,:)
910      ELSE                                        !* closed
911         IF( .NOT. cd_type1 == 'F' )   ptab1(     1       :jpreci,:,:) = 0.e0    ! south except at F-point
912         IF( .NOT. cd_type2 == 'F' )   ptab2(     1       :jpreci,:,:) = 0.e0
913                                       ptab1(nlci-jpreci+1:jpi   ,:,:) = 0.e0    ! north
914                                       ptab2(nlci-jpreci+1:jpi   ,:,:) = 0.e0
915      ENDIF
916
917     
918      !                                      ! North-South boundaries
919      IF( .NOT. cd_type1 == 'F' )   ptab1(:,     1       :jprecj,:) = 0.e0    ! south except at F-point
920      IF( .NOT. cd_type2 == 'F' )   ptab2(:,     1       :jprecj,:) = 0.e0
921                                    ptab1(:,nlcj-jprecj+1:jpj   ,:) = 0.e0    ! north
922                                    ptab2(:,nlcj-jprecj+1:jpj   ,:) = 0.e0
923
924
925      ! 2. East and west directions exchange
926      ! ------------------------------------
927      ! we play with the neigbours AND the row number because of the periodicity
928      !
929      SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
930      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
931         iihom = nlci-nreci
932         DO jl = 1, jpreci
933            t4ew(:,jl,:,1,1) = ptab1(jpreci+jl,:,:)
934            t4we(:,jl,:,1,1) = ptab1(iihom +jl,:,:)
935            t4ew(:,jl,:,2,1) = ptab2(jpreci+jl,:,:)
936            t4we(:,jl,:,2,1) = ptab2(iihom +jl,:,:)
937         END DO
938      END SELECT
939      !
940      !                           ! Migrations
941      imigr = jpreci * jpj * jpk *2
942      !
943      SELECT CASE ( nbondi ) 
944      CASE ( -1 )
945         CALL mppsend( 2, t4we(1,1,1,1,1), imigr, noea, ml_req1 )
946         CALL mpprecv( 1, t4ew(1,1,1,1,2), imigr )
947         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
948      CASE ( 0 )
949         CALL mppsend( 1, t4ew(1,1,1,1,1), imigr, nowe, ml_req1 )
950         CALL mppsend( 2, t4we(1,1,1,1,1), imigr, noea, ml_req2 )
951         CALL mpprecv( 1, t4ew(1,1,1,1,2), imigr )
952         CALL mpprecv( 2, t4we(1,1,1,1,2), imigr )
953         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
954         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
955      CASE ( 1 )
956         CALL mppsend( 1, t4ew(1,1,1,1,1), imigr, nowe, ml_req1 )
957         CALL mpprecv( 2, t4we(1,1,1,1,2), imigr )
958         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
959      END SELECT
960      !
961      !                           ! Write Dirichlet lateral conditions
962      iihom = nlci - jpreci
963      !
964      SELECT CASE ( nbondi )
965      CASE ( -1 )
966         DO jl = 1, jpreci
967            ptab1(iihom+jl,:,:) = t4ew(:,jl,:,1,2)
968            ptab2(iihom+jl,:,:) = t4ew(:,jl,:,2,2)
969         END DO
970      CASE ( 0 ) 
971         DO jl = 1, jpreci
972            ptab1(jl      ,:,:) = t4we(:,jl,:,1,2)
973            ptab1(iihom+jl,:,:) = t4ew(:,jl,:,1,2)
974            ptab2(jl      ,:,:) = t4we(:,jl,:,2,2)
975            ptab2(iihom+jl,:,:) = t4ew(:,jl,:,2,2)
976         END DO
977      CASE ( 1 )
978         DO jl = 1, jpreci
979            ptab1(jl      ,:,:) = t4we(:,jl,:,1,2)
980            ptab2(jl      ,:,:) = t4we(:,jl,:,2,2)
981         END DO
982      END SELECT
983
984
985      ! 3. North and south directions
986      ! -----------------------------
987      ! always closed : we play only with the neigbours
988      !
989      IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
990         ijhom = nlcj - nrecj
991         DO jl = 1, jprecj
992            t4sn(:,jl,:,1,1) = ptab1(:,ijhom +jl,:)
993            t4ns(:,jl,:,1,1) = ptab1(:,jprecj+jl,:)
994            t4sn(:,jl,:,2,1) = ptab2(:,ijhom +jl,:)
995            t4ns(:,jl,:,2,1) = ptab2(:,jprecj+jl,:)
996         END DO
997      ENDIF
998      !
999      !                           ! Migrations
1000      imigr = jprecj * jpi * jpk * 2
1001      !
1002      SELECT CASE ( nbondj )     
1003      CASE ( -1 )
1004         CALL mppsend( 4, t4sn(1,1,1,1,1), imigr, nono, ml_req1 )
1005         CALL mpprecv( 3, t4ns(1,1,1,1,2), imigr )
1006         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
1007      CASE ( 0 )
1008         CALL mppsend( 3, t4ns(1,1,1,1,1), imigr, noso, ml_req1 )
1009         CALL mppsend( 4, t4sn(1,1,1,1,1), imigr, nono, ml_req2 )
1010         CALL mpprecv( 3, t4ns(1,1,1,1,2), imigr )
1011         CALL mpprecv( 4, t4sn(1,1,1,1,2), imigr )
1012         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
1013         IF(l_isend) CALL mpi_wait(ml_req2, ml_stat, ml_err)
1014      CASE ( 1 ) 
1015         CALL mppsend( 3, t4ns(1,1,1,1,1), imigr, noso, ml_req1 )
1016         CALL mpprecv( 4, t4sn(1,1,1,1,2), imigr )
1017         IF(l_isend) CALL mpi_wait(ml_req1, ml_stat, ml_err)
1018      END SELECT
1019      !
1020      !                           ! Write Dirichlet lateral conditions
1021      ijhom = nlcj - jprecj
1022      !
1023      SELECT CASE ( nbondj )
1024      CASE ( -1 )
1025         DO jl = 1, jprecj
1026            ptab1(:,ijhom+jl,:) = t4ns(:,jl,:,1,2)
1027            ptab2(:,ijhom+jl,:) = t4ns(:,jl,:,2,2)
1028         END DO
1029      CASE ( 0 ) 
1030         DO jl = 1, jprecj
1031            ptab1(:,jl      ,:) = t4sn(:,jl,:,1,2)
1032            ptab1(:,ijhom+jl,:) = t4ns(:,jl,:,1,2)
1033            ptab2(:,jl      ,:) = t4sn(:,jl,:,2,2)
1034            ptab2(:,ijhom+jl,:) = t4ns(:,jl,:,2,2)
1035         END DO
1036      CASE ( 1 )
1037         DO jl = 1, jprecj
1038            ptab1(:,jl,:) = t4sn(:,jl,:,1,2)
1039            ptab2(:,jl,:) = t4sn(:,jl,:,2,2)
1040         END DO
1041      END SELECT
1042
1043
1044      ! 4. north fold treatment
1045      ! -----------------------
1046      IF( npolj /= 0 ) THEN
1047         !
1048         SELECT CASE ( jpni )
1049         CASE ( 1 )                                           
1050            CALL lbc_nfd      ( ptab1, cd_type1, psgn )   ! only for northern procs.
1051            CALL lbc_nfd      ( ptab2, cd_type2, psgn )
1052         CASE DEFAULT
1053            CALL mpp_lbc_north( ptab1, cd_type1, psgn )   ! for all northern procs.
1054            CALL mpp_lbc_north (ptab2, cd_type2, psgn)
1055         END SELECT 
1056         !
1057      ENDIF
1058      !
1059   END SUBROUTINE mpp_lnk_3d_gather
1060
1061
1062   SUBROUTINE mpp_lnk_2d_e( pt2d, cd_type, psgn )
1063      !!----------------------------------------------------------------------
1064      !!                  ***  routine mpp_lnk_2d_e  ***
1065      !!                 
1066      !! ** Purpose :   Message passing manadgement for 2d array (with halo)
1067      !!
1068      !! ** Method  :   Use mppsend and mpprecv function for passing mask
1069      !!      between processors following neighboring subdomains.
1070      !!            domain parameters
1071      !!                    nlci   : first dimension of the local subdomain
1072      !!                    nlcj   : second dimension of the local subdomain
1073      !!                    jpr2di : number of rows for extra outer halo
1074      !!                    jpr2dj : number of columns for extra outer halo
1075      !!                    nbondi : mark for "east-west local boundary"
1076      !!                    nbondj : mark for "north-south local boundary"
1077      !!                    noea   : number for local neighboring processors
1078      !!                    nowe   : number for local neighboring processors
1079      !!                    noso   : number for local neighboring processors
1080      !!                    nono   : number for local neighboring processors
1081      !!
1082      !!----------------------------------------------------------------------
1083      REAL(wp), DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj), INTENT(inout) ::   pt2d     ! 2D array with extra halo
1084      CHARACTER(len=1)                                            , INTENT(in   ) ::   cd_type  ! nature of ptab array grid-points
1085      !                                                                                         ! = T , U , V , F , W and I points
1086      REAL(wp)                                                    , INTENT(in   ) ::   psgn     ! =-1 the sign change across the
1087      !!                                                                                        ! north boundary, =  1. otherwise
1088      INTEGER  ::   jl   ! dummy loop indices
1089      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
1090      INTEGER  ::   ipreci, iprecj             ! temporary integers
1091      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
1092      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
1093      !!----------------------------------------------------------------------
1094
1095      ipreci = jpreci + jpr2di      ! take into account outer extra 2D overlap area
1096      iprecj = jprecj + jpr2dj
1097
1098
1099      ! 1. standard boundary treatment
1100      ! ------------------------------
1101      ! Order matters Here !!!!
1102      !
1103      !                                      !* North-South boundaries (always colsed)
1104      IF( .NOT. cd_type == 'F' )   pt2d(:,  1-jpr2dj   :  jprecj  ) = 0.e0    ! south except at F-point
1105                                   pt2d(:,nlcj-jprecj+1:jpj+jpr2dj) = 0.e0    ! north
1106                               
1107      !                                      ! East-West boundaries
1108      !                                           !* Cyclic east-west
1109      IF( nbondi == 2 .AND. (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
1110         pt2d(1-jpr2di:     1    ,:) = pt2d(jpim1-jpr2di:  jpim1 ,:)       ! east
1111         pt2d(   jpi  :jpi+jpr2di,:) = pt2d(     2      :2+jpr2di,:)       ! west
1112         !
1113      ELSE                                        !* closed
1114         IF( .NOT. cd_type == 'F' )   pt2d(  1-jpr2di   :jpreci    ,:) = 0.e0    ! south except at F-point
1115                                      pt2d(nlci-jpreci+1:jpi+jpr2di,:) = 0.e0    ! north
1116      ENDIF
1117      !
1118
1119      ! north fold treatment
1120      ! -----------------------
1121      IF( npolj /= 0 ) THEN
1122         !
1123         SELECT CASE ( jpni )
1124         CASE ( 1 )     ;   CALL lbc_nfd        ( pt2d(1:jpi,1:jpj+jpr2dj), cd_type, psgn, pr2dj=jpr2dj )
1125         CASE DEFAULT   ;   CALL mpp_lbc_north_e( pt2d                    , cd_type, psgn               )
1126         END SELECT 
1127         !
1128      ENDIF
1129
1130      ! 2. East and west directions exchange
1131      ! ------------------------------------
1132      ! we play with the neigbours AND the row number because of the periodicity
1133      !
1134      SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
1135      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
1136         iihom = nlci-nreci-jpr2di
1137         DO jl = 1, ipreci
1138            tr2ew(:,jl,1) = pt2d(jpreci+jl,:)
1139            tr2we(:,jl,1) = pt2d(iihom +jl,:)
1140         END DO
1141      END SELECT
1142      !
1143      !                           ! Migrations
1144      imigr = ipreci * ( jpj + 2*jpr2dj)
1145      !
1146      SELECT CASE ( nbondi )
1147      CASE ( -1 )
1148         CALL mppsend( 2, tr2we(1-jpr2dj,1,1), imigr, noea, ml_req1 )
1149         CALL mpprecv( 1, tr2ew(1-jpr2dj,1,2), imigr )
1150         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1151      CASE ( 0 )
1152         CALL mppsend( 1, tr2ew(1-jpr2dj,1,1), imigr, nowe, ml_req1 )
1153         CALL mppsend( 2, tr2we(1-jpr2dj,1,1), imigr, noea, ml_req2 )
1154         CALL mpprecv( 1, tr2ew(1-jpr2dj,1,2), imigr )
1155         CALL mpprecv( 2, tr2we(1-jpr2dj,1,2), imigr )
1156         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1157         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
1158      CASE ( 1 )
1159         CALL mppsend( 1, tr2ew(1-jpr2dj,1,1), imigr, nowe, ml_req1 )
1160         CALL mpprecv( 2, tr2we(1-jpr2dj,1,2), imigr )
1161         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1162      END SELECT
1163      !
1164      !                           ! Write Dirichlet lateral conditions
1165      iihom = nlci - jpreci
1166      !
1167      SELECT CASE ( nbondi )
1168      CASE ( -1 )
1169         DO jl = 1, ipreci
1170            pt2d(iihom+jl,:) = tr2ew(:,jl,2)
1171         END DO
1172      CASE ( 0 )
1173         DO jl = 1, ipreci
1174            pt2d(jl-jpr2di,:) = tr2we(:,jl,2)
1175            pt2d( iihom+jl,:) = tr2ew(:,jl,2)
1176         END DO
1177      CASE ( 1 )
1178         DO jl = 1, ipreci
1179            pt2d(jl-jpr2di,:) = tr2we(:,jl,2)
1180         END DO
1181      END SELECT
1182
1183
1184      ! 3. North and south directions
1185      ! -----------------------------
1186      ! always closed : we play only with the neigbours
1187      !
1188      IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
1189         ijhom = nlcj-nrecj-jpr2dj
1190         DO jl = 1, iprecj
1191            tr2sn(:,jl,1) = pt2d(:,ijhom +jl)
1192            tr2ns(:,jl,1) = pt2d(:,jprecj+jl)
1193         END DO
1194      ENDIF
1195      !
1196      !                           ! Migrations
1197      imigr = iprecj * ( jpi + 2*jpr2di )
1198      !
1199      SELECT CASE ( nbondj )
1200      CASE ( -1 )
1201         CALL mppsend( 4, tr2sn(1-jpr2di,1,1), imigr, nono, ml_req1 )
1202         CALL mpprecv( 3, tr2ns(1-jpr2di,1,2), imigr )
1203         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1204      CASE ( 0 )
1205         CALL mppsend( 3, tr2ns(1-jpr2di,1,1), imigr, noso, ml_req1 )
1206         CALL mppsend( 4, tr2sn(1-jpr2di,1,1), imigr, nono, ml_req2 )
1207         CALL mpprecv( 3, tr2ns(1-jpr2di,1,2), imigr )
1208         CALL mpprecv( 4, tr2sn(1-jpr2di,1,2), imigr )
1209         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1210         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
1211      CASE ( 1 )
1212         CALL mppsend( 3, tr2ns(1-jpr2di,1,1), imigr, noso, ml_req1 )
1213         CALL mpprecv( 4, tr2sn(1-jpr2di,1,2), imigr )
1214         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1215      END SELECT
1216      !
1217      !                           ! Write Dirichlet lateral conditions
1218      ijhom = nlcj - jprecj 
1219      !
1220      SELECT CASE ( nbondj )
1221      CASE ( -1 )
1222         DO jl = 1, iprecj
1223            pt2d(:,ijhom+jl) = tr2ns(:,jl,2)
1224         END DO
1225      CASE ( 0 )
1226         DO jl = 1, iprecj
1227            pt2d(:,jl-jpr2dj) = tr2sn(:,jl,2)
1228            pt2d(:,ijhom+jl ) = tr2ns(:,jl,2)
1229         END DO
1230      CASE ( 1 ) 
1231         DO jl = 1, iprecj
1232            pt2d(:,jl-jpr2dj) = tr2sn(:,jl,2)
1233         END DO
1234      END SELECT
1235
1236   END SUBROUTINE mpp_lnk_2d_e
1237
1238
1239   SUBROUTINE mppsend( ktyp, pmess, kbytes, kdest, md_req )
1240      !!----------------------------------------------------------------------
1241      !!                  ***  routine mppsend  ***
1242      !!                   
1243      !! ** Purpose :   Send message passing array
1244      !!
1245      !!----------------------------------------------------------------------
1246      REAL(wp), INTENT(inout) ::   pmess(*)   ! array of real
1247      INTEGER , INTENT(in   ) ::   kbytes     ! size of the array pmess
1248      INTEGER , INTENT(in   ) ::   kdest      ! receive process number
1249      INTEGER , INTENT(in   ) ::   ktyp       ! tag of the message
1250      INTEGER , INTENT(in   ) ::   md_req     ! argument for isend
1251      !!
1252      INTEGER ::   iflag
1253      !!----------------------------------------------------------------------
1254      !
1255      SELECT CASE ( cn_mpi_send )
1256      CASE ( 'S' )                ! Standard mpi send (blocking)
1257         CALL mpi_send ( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_opa        , iflag )
1258      CASE ( 'B' )                ! Buffer mpi send (blocking)
1259         CALL mpi_bsend( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_opa        , iflag )
1260      CASE ( 'I' )                ! Immediate mpi send (non-blocking send)
1261         ! be carefull, one more argument here : the mpi request identifier..
1262         CALL mpi_isend( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_opa, md_req, iflag )
1263      END SELECT
1264      !
1265   END SUBROUTINE mppsend
1266
1267
1268   SUBROUTINE mpprecv( ktyp, pmess, kbytes )
1269      !!----------------------------------------------------------------------
1270      !!                  ***  routine mpprecv  ***
1271      !!
1272      !! ** Purpose :   Receive messag passing array
1273      !!
1274      !!----------------------------------------------------------------------
1275      REAL(wp), INTENT(inout) ::   pmess(*)   ! array of real
1276      INTEGER , INTENT(in   ) ::   kbytes     ! suze of the array pmess
1277      INTEGER , INTENT(in   ) ::   ktyp       ! Tag of the recevied message
1278      !!
1279      INTEGER :: istatus(mpi_status_size)
1280      INTEGER :: iflag
1281      !!----------------------------------------------------------------------
1282      !
1283      CALL mpi_recv( pmess, kbytes, mpi_double_precision, mpi_any_source, ktyp, mpi_comm_opa, istatus, iflag )
1284      !
1285   END SUBROUTINE mpprecv
1286
1287
1288   SUBROUTINE mppgather( ptab, kp, pio )
1289      !!----------------------------------------------------------------------
1290      !!                   ***  routine mppgather  ***
1291      !!                   
1292      !! ** Purpose :   Transfert between a local subdomain array and a work
1293      !!     array which is distributed following the vertical level.
1294      !!
1295      !!----------------------------------------------------------------------
1296      REAL(wp), DIMENSION(jpi,jpj),       INTENT(in   ) ::   ptab   ! subdomain input array
1297      INTEGER ,                           INTENT(in   ) ::   kp     ! record length
1298      REAL(wp), DIMENSION(jpi,jpj,jpnij), INTENT(  out) ::   pio    ! subdomain input array
1299      !!
1300      INTEGER :: itaille, ierror   ! temporary integer
1301      !!---------------------------------------------------------------------
1302      !
1303      itaille = jpi * jpj
1304      CALL mpi_gather( ptab, itaille, mpi_double_precision, pio, itaille     ,   &
1305         &                            mpi_double_precision, kp , mpi_comm_opa, ierror ) 
1306      !
1307   END SUBROUTINE mppgather
1308
1309
1310   SUBROUTINE mppscatter( pio, kp, ptab )
1311      !!----------------------------------------------------------------------
1312      !!                  ***  routine mppscatter  ***
1313      !!
1314      !! ** Purpose :   Transfert between awork array which is distributed
1315      !!      following the vertical level and the local subdomain array.
1316      !!
1317      !!----------------------------------------------------------------------
1318      REAL(wp), DIMENSION(jpi,jpj,jpnij)  ::  pio        ! output array
1319      INTEGER                             ::   kp        ! Tag (not used with MPI
1320      REAL(wp), DIMENSION(jpi,jpj)        ::  ptab       ! subdomain array input
1321      !!
1322      INTEGER :: itaille, ierror   ! temporary integer
1323      !!---------------------------------------------------------------------
1324      !
1325      itaille=jpi*jpj
1326      !
1327      CALL mpi_scatter( pio, itaille, mpi_double_precision, ptab, itaille     ,   &
1328         &                            mpi_double_precision, kp  , mpi_comm_opa, ierror )
1329      !
1330   END SUBROUTINE mppscatter
1331
1332
1333   SUBROUTINE mppmax_a_int( ktab, kdim, kcom )
1334      !!----------------------------------------------------------------------
1335      !!                  ***  routine mppmax_a_int  ***
1336      !!
1337      !! ** Purpose :   Find maximum value in an integer layout array
1338      !!
1339      !!----------------------------------------------------------------------
1340      INTEGER , INTENT(in   )                  ::   kdim   ! size of array
1341      INTEGER , INTENT(inout), DIMENSION(kdim) ::   ktab   ! input array
1342      INTEGER , INTENT(in   ), OPTIONAL        ::   kcom   !
1343      !!
1344      INTEGER :: ierror, localcomm   ! temporary integer
1345      INTEGER, DIMENSION(kdim) ::   iwork
1346      !!----------------------------------------------------------------------
1347      !
1348      localcomm = mpi_comm_opa
1349      IF( PRESENT(kcom) )   localcomm = kcom
1350      !
1351      CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_max, localcomm, ierror )
1352      !
1353      ktab(:) = iwork(:)
1354      !
1355   END SUBROUTINE mppmax_a_int
1356
1357
1358   SUBROUTINE mppmax_int( ktab, kcom )
1359      !!----------------------------------------------------------------------
1360      !!                  ***  routine mppmax_int  ***
1361      !!
1362      !! ** Purpose :   Find maximum value in an integer layout array
1363      !!
1364      !!----------------------------------------------------------------------
1365      INTEGER, INTENT(inout)           ::   ktab      ! ???
1366      INTEGER, INTENT(in   ), OPTIONAL ::   kcom      ! ???
1367      !!
1368      INTEGER ::   ierror, iwork, localcomm   ! temporary integer
1369      !!----------------------------------------------------------------------
1370      !
1371      localcomm = mpi_comm_opa 
1372      IF( PRESENT(kcom) )   localcomm = kcom
1373      !
1374      CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_max, localcomm, ierror)
1375      !
1376      ktab = iwork
1377      !
1378   END SUBROUTINE mppmax_int
1379
1380
1381   SUBROUTINE mppmin_a_int( ktab, kdim, kcom )
1382      !!----------------------------------------------------------------------
1383      !!                  ***  routine mppmin_a_int  ***
1384      !!
1385      !! ** Purpose :   Find minimum value in an integer layout array
1386      !!
1387      !!----------------------------------------------------------------------
1388      INTEGER , INTENT( in  )                  ::   kdim        ! size of array
1389      INTEGER , INTENT(inout), DIMENSION(kdim) ::   ktab        ! input array
1390      INTEGER , INTENT( in  ), OPTIONAL        ::   kcom        ! input array
1391      !!
1392      INTEGER ::   ierror, localcomm   ! temporary integer
1393      INTEGER, DIMENSION(kdim) ::   iwork
1394      !!----------------------------------------------------------------------
1395      !
1396      localcomm = mpi_comm_opa
1397      IF( PRESENT(kcom) )   localcomm = kcom
1398      !
1399      CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_min, localcomm, ierror )
1400      !
1401      ktab(:) = iwork(:)
1402      !
1403   END SUBROUTINE mppmin_a_int
1404
1405
1406   SUBROUTINE mppmin_int( ktab, kcom )
1407      !!----------------------------------------------------------------------
1408      !!                  ***  routine mppmin_int  ***
1409      !!
1410      !! ** Purpose :   Find minimum value in an integer layout array
1411      !!
1412      !!----------------------------------------------------------------------
1413      INTEGER, INTENT(inout) ::   ktab      ! ???
1414      INTEGER , INTENT( in  ), OPTIONAL        ::   kcom        ! input array
1415      !!
1416      INTEGER ::  ierror, iwork, localcomm
1417      !!----------------------------------------------------------------------
1418      !
1419      localcomm = mpi_comm_opa
1420      IF( PRESENT(kcom) )   localcomm = kcom
1421      !
1422     CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_min, localcomm, ierror )
1423      !
1424      ktab = iwork
1425      !
1426   END SUBROUTINE mppmin_int
1427
1428
1429   SUBROUTINE mppsum_a_int( ktab, kdim )
1430      !!----------------------------------------------------------------------
1431      !!                  ***  routine mppsum_a_int  ***
1432      !!                   
1433      !! ** Purpose :   Global integer sum, 1D array case
1434      !!
1435      !!----------------------------------------------------------------------
1436      INTEGER, INTENT(in   )                   ::   kdim      ! ???
1437      INTEGER, INTENT(inout), DIMENSION (kdim) ::   ktab      ! ???
1438      !!
1439      INTEGER :: ierror
1440      INTEGER, DIMENSION (kdim) ::  iwork
1441      !!----------------------------------------------------------------------
1442      !
1443      CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_sum, mpi_comm_opa, ierror )
1444      !
1445      ktab(:) = iwork(:)
1446      !
1447   END SUBROUTINE mppsum_a_int
1448
1449
1450   SUBROUTINE mppsum_int( ktab )
1451      !!----------------------------------------------------------------------
1452      !!                 ***  routine mppsum_int  ***
1453      !!                 
1454      !! ** Purpose :   Global integer sum
1455      !!
1456      !!----------------------------------------------------------------------
1457      INTEGER, INTENT(inout) ::   ktab
1458      !!
1459      INTEGER :: ierror, iwork
1460      !!----------------------------------------------------------------------
1461      !
1462      CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_sum, mpi_comm_opa, ierror )
1463      !
1464      ktab = iwork
1465      !
1466   END SUBROUTINE mppsum_int
1467
1468
1469   SUBROUTINE mppmax_a_real( ptab, kdim, kcom )
1470      !!----------------------------------------------------------------------
1471      !!                 ***  routine mppmax_a_real  ***
1472      !!                 
1473      !! ** Purpose :   Maximum
1474      !!
1475      !!----------------------------------------------------------------------
1476      INTEGER , INTENT(in   )                  ::   kdim
1477      REAL(wp), INTENT(inout), DIMENSION(kdim) ::   ptab
1478      INTEGER , INTENT(in   ), OPTIONAL        ::   kcom
1479      !!
1480      INTEGER :: ierror, localcomm
1481      REAL(wp), DIMENSION(kdim) ::  zwork
1482      !!----------------------------------------------------------------------
1483      !
1484      localcomm = mpi_comm_opa
1485      IF( PRESENT(kcom) ) localcomm = kcom
1486      !
1487      CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_max, localcomm, ierror )
1488      ptab(:) = zwork(:)
1489      !
1490   END SUBROUTINE mppmax_a_real
1491
1492
1493   SUBROUTINE mppmax_real( ptab, kcom )
1494      !!----------------------------------------------------------------------
1495      !!                  ***  routine mppmax_real  ***
1496      !!                   
1497      !! ** Purpose :   Maximum
1498      !!
1499      !!----------------------------------------------------------------------
1500      REAL(wp), INTENT(inout)           ::   ptab   ! ???
1501      INTEGER , INTENT(in   ), OPTIONAL ::   kcom   ! ???
1502      !!
1503      INTEGER  ::   ierror, localcomm
1504      REAL(wp) ::   zwork
1505      !!----------------------------------------------------------------------
1506      !
1507      localcomm = mpi_comm_opa 
1508      IF( PRESENT(kcom) )   localcomm = kcom
1509      !
1510      CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_max, localcomm, ierror )
1511      ptab = zwork
1512      !
1513   END SUBROUTINE mppmax_real
1514
1515
1516   SUBROUTINE mppmin_a_real( ptab, kdim, kcom )
1517      !!----------------------------------------------------------------------
1518      !!                 ***  routine mppmin_a_real  ***
1519      !!                 
1520      !! ** Purpose :   Minimum of REAL, array case
1521      !!
1522      !!-----------------------------------------------------------------------
1523      INTEGER , INTENT(in   )                  ::   kdim
1524      REAL(wp), INTENT(inout), DIMENSION(kdim) ::   ptab
1525      INTEGER , INTENT(in   ), OPTIONAL        ::   kcom
1526      !!
1527      INTEGER :: ierror, localcomm
1528      REAL(wp), DIMENSION(kdim) ::   zwork
1529      !!-----------------------------------------------------------------------
1530      !
1531      localcomm = mpi_comm_opa 
1532      IF( PRESENT(kcom) ) localcomm = kcom
1533      !
1534      CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_min, localcomm, ierror )
1535      ptab(:) = zwork(:)
1536      !
1537   END SUBROUTINE mppmin_a_real
1538
1539
1540   SUBROUTINE mppmin_real( ptab, kcom )
1541      !!----------------------------------------------------------------------
1542      !!                  ***  routine mppmin_real  ***
1543      !!
1544      !! ** Purpose :   minimum of REAL, scalar case
1545      !!
1546      !!-----------------------------------------------------------------------
1547      REAL(wp), INTENT(inout)           ::   ptab        !
1548      INTEGER , INTENT(in   ), OPTIONAL :: kcom
1549      !!
1550      INTEGER  ::   ierror
1551      REAL(wp) ::   zwork
1552      INTEGER :: localcomm
1553      !!-----------------------------------------------------------------------
1554      !
1555      localcomm = mpi_comm_opa 
1556      IF( PRESENT(kcom) )   localcomm = kcom
1557      !
1558      CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_min, localcomm, ierror )
1559      ptab = zwork
1560      !
1561   END SUBROUTINE mppmin_real
1562
1563
1564   SUBROUTINE mppsum_a_real( ptab, kdim, kcom )
1565      !!----------------------------------------------------------------------
1566      !!                  ***  routine mppsum_a_real  ***
1567      !!
1568      !! ** Purpose :   global sum, REAL ARRAY argument case
1569      !!
1570      !!-----------------------------------------------------------------------
1571      INTEGER , INTENT( in )                     ::   kdim      ! size of ptab
1572      REAL(wp), DIMENSION(kdim), INTENT( inout ) ::   ptab      ! input array
1573      INTEGER , INTENT( in ), OPTIONAL           :: kcom
1574      !!
1575      INTEGER                   ::   ierror    ! temporary integer
1576      INTEGER                   ::   localcomm 
1577      REAL(wp), DIMENSION(kdim) ::   zwork     ! temporary workspace
1578      !!-----------------------------------------------------------------------
1579      !
1580      localcomm = mpi_comm_opa 
1581      IF( PRESENT(kcom) )   localcomm = kcom
1582      !
1583      CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_sum, localcomm, ierror )
1584      ptab(:) = zwork(:)
1585      !
1586   END SUBROUTINE mppsum_a_real
1587
1588
1589   SUBROUTINE mppsum_real( ptab, kcom )
1590      !!----------------------------------------------------------------------
1591      !!                  ***  routine mppsum_real  ***
1592      !!             
1593      !! ** Purpose :   global sum, SCALAR argument case
1594      !!
1595      !!-----------------------------------------------------------------------
1596      REAL(wp), INTENT(inout)           ::   ptab   ! input scalar
1597      INTEGER , INTENT(in   ), OPTIONAL ::   kcom
1598      !!
1599      INTEGER  ::   ierror, localcomm 
1600      REAL(wp) ::   zwork
1601      !!-----------------------------------------------------------------------
1602      !
1603      localcomm = mpi_comm_opa 
1604      IF( PRESENT(kcom) ) localcomm = kcom
1605      !
1606      CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_sum, localcomm, ierror )
1607      ptab = zwork
1608      !
1609   END SUBROUTINE mppsum_real
1610
1611# if defined key_mpp_rep
1612   SUBROUTINE mppsum_realdd( ytab, kcom )
1613      !!----------------------------------------------------------------------
1614      !!                  ***  routine mppsum_realdd ***
1615      !!
1616      !! ** Purpose :   global sum in Massively Parallel Processing
1617      !!                SCALAR argument case for double-double precision
1618      !!
1619      !!-----------------------------------------------------------------------
1620      COMPLEX(wp), INTENT(inout)         :: ytab    ! input scalar
1621      INTEGER , INTENT( in  ), OPTIONAL :: kcom
1622
1623      !! * Local variables   (MPI version)
1624      INTEGER  ::    ierror
1625      INTEGER  ::   localcomm
1626      COMPLEX(wp) :: zwork
1627
1628      localcomm = mpi_comm_opa
1629      IF( PRESENT(kcom) ) localcomm = kcom
1630
1631      ! reduce local sums into global sum
1632      CALL MPI_ALLREDUCE (ytab, zwork, 1, MPI_DOUBLE_COMPLEX, &
1633                       MPI_SUMDD,localcomm,ierror)
1634      ytab = zwork
1635
1636   END SUBROUTINE mppsum_realdd
1637 
1638 
1639   SUBROUTINE mppsum_a_realdd( ytab, kdim, kcom )
1640      !!----------------------------------------------------------------------
1641      !!                  ***  routine mppsum_a_realdd  ***
1642      !!
1643      !! ** Purpose :   global sum in Massively Parallel Processing
1644      !!                COMPLEX ARRAY case for double-double precision
1645      !!
1646      !!-----------------------------------------------------------------------
1647      INTEGER , INTENT( in )                        ::   kdim      ! size of ytab
1648      COMPLEX(wp), DIMENSION(kdim), INTENT( inout ) ::   ytab      ! input array
1649      INTEGER , INTENT( in  ), OPTIONAL :: kcom
1650
1651      !! * Local variables   (MPI version)
1652      INTEGER                      :: ierror    ! temporary integer
1653      INTEGER                      ::   localcomm
1654      COMPLEX(wp), DIMENSION(kdim) :: zwork     ! temporary workspace
1655
1656      localcomm = mpi_comm_opa
1657      IF( PRESENT(kcom) ) localcomm = kcom
1658
1659      CALL MPI_ALLREDUCE (ytab, zwork, kdim, MPI_DOUBLE_COMPLEX, &
1660                       MPI_SUMDD,localcomm,ierror)
1661      ytab(:) = zwork(:)
1662
1663   END SUBROUTINE mppsum_a_realdd
1664# endif   
1665   
1666   SUBROUTINE mpp_minloc2d( ptab, pmask, pmin, ki,kj )
1667      !!------------------------------------------------------------------------
1668      !!             ***  routine mpp_minloc  ***
1669      !!
1670      !! ** Purpose :   Compute the global minimum of an array ptab
1671      !!              and also give its global position
1672      !!
1673      !! ** Method  :   Use MPI_ALLREDUCE with MPI_MINLOC
1674      !!
1675      !!--------------------------------------------------------------------------
1676      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   ptab    ! Local 2D array
1677      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   pmask   ! Local mask
1678      REAL(wp)                     , INTENT(  out) ::   pmin    ! Global minimum of ptab
1679      INTEGER                      , INTENT(  out) ::   ki, kj   ! index of minimum in global frame
1680      !!
1681      INTEGER , DIMENSION(2)   ::   ilocs
1682      INTEGER :: ierror
1683      REAL(wp) ::   zmin   ! local minimum
1684      REAL(wp), DIMENSION(2,1) ::   zain, zaout
1685      !!-----------------------------------------------------------------------
1686      !
1687      zmin  = MINVAL( ptab(:,:) , mask= pmask == 1.e0 )
1688      ilocs = MINLOC( ptab(:,:) , mask= pmask == 1.e0 )
1689      !
1690      ki = ilocs(1) + nimpp - 1
1691      kj = ilocs(2) + njmpp - 1
1692      !
1693      zain(1,:)=zmin
1694      zain(2,:)=ki+10000.*kj
1695      !
1696      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MINLOC,MPI_COMM_OPA,ierror)
1697      !
1698      pmin = zaout(1,1)
1699      kj = INT(zaout(2,1)/10000.)
1700      ki = INT(zaout(2,1) - 10000.*kj )
1701      !
1702   END SUBROUTINE mpp_minloc2d
1703
1704
1705   SUBROUTINE mpp_minloc3d( ptab3d, pmask3d, pmin, ki, kj ,kk)
1706      !!------------------------------------------------------------------------
1707      !!             ***  routine mpp_minloc  ***
1708      !!
1709      !! ** Purpose :   Compute the global minimum of an array ptab
1710      !!              and also give its global position
1711      !!
1712      !! ** Method  :   Use MPI_ALLREDUCE with MPI_MINLOC
1713      !!
1714      !!--------------------------------------------------------------------------
1715!FTRANS ptab3d  :I :I :z
1716!FTRANS pmask3d :I :I :z
1717!! DCSE_NEMO: work around a deficiency in ftrans
1718!     REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in   ) ::   ptab3d       ! Local 3D array
1719!     REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in   ) ::   pmask3d      ! Local mask
1720      REAL(wp), INTENT(in   )                          ::   ptab3d(jpi,jpj,jpk)
1721      REAL(wp), INTENT(in   )                          ::   pmask3d(jpi,jpj,jpk)
1722      REAL(wp)                         , INTENT(  out) ::   pmin         ! Global minimum of ptab
1723      INTEGER                          , INTENT(  out) ::   ki, kj, kk   ! index of minimum in global frame
1724      !!
1725      INTEGER  ::   ierror
1726      REAL(wp) ::   zmin     ! local minimum
1727      INTEGER , DIMENSION(3)   ::   ilocs
1728      REAL(wp), DIMENSION(2,1) ::   zain, zaout
1729      !!-----------------------------------------------------------------------
1730      !
1731      zmin  = MINVAL( ptab3d(:,:,:) , mask= pmask3d == 1.e0 )
1732      ilocs = MINLOC( ptab3d(:,:,:) , mask= pmask3d == 1.e0 )
1733      !
1734!! DCSE_NEMO: Attention!
1735#if defined key_z_first
1736      ki = ilocs(2) + nimpp - 1
1737      kj = ilocs(3) + njmpp - 1
1738      kk = ilocs(1)
1739#else
1740      ki = ilocs(1) + nimpp - 1
1741      kj = ilocs(2) + njmpp - 1
1742      kk = ilocs(3)
1743#endif
1744      !
1745      zain(1,:)=zmin
1746      zain(2,:)=ki+10000.*kj+100000000.*kk
1747      !
1748      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MINLOC,MPI_COMM_OPA,ierror)
1749      !
1750      pmin = zaout(1,1)
1751      kk   = INT( zaout(2,1) / 100000000. )
1752      kj   = INT( zaout(2,1) - kk * 100000000. ) / 10000
1753      ki   = INT( zaout(2,1) - kk * 100000000. -kj * 10000. )
1754      !
1755   END SUBROUTINE mpp_minloc3d
1756
1757
1758   SUBROUTINE mpp_maxloc2d( ptab, pmask, pmax, ki, kj )
1759      !!------------------------------------------------------------------------
1760      !!             ***  routine mpp_maxloc  ***
1761      !!
1762      !! ** Purpose :   Compute the global maximum of an array ptab
1763      !!              and also give its global position
1764      !!
1765      !! ** Method  :   Use MPI_ALLREDUCE with MPI_MINLOC
1766      !!
1767      !!--------------------------------------------------------------------------
1768      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   ptab     ! Local 2D array
1769      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   pmask    ! Local mask
1770      REAL(wp)                     , INTENT(  out) ::   pmax     ! Global maximum of ptab
1771      INTEGER                      , INTENT(  out) ::   ki, kj   ! index of maximum in global frame
1772      !! 
1773      INTEGER  :: ierror
1774      INTEGER, DIMENSION (2)   ::   ilocs
1775      REAL(wp) :: zmax   ! local maximum
1776      REAL(wp), DIMENSION(2,1) ::   zain, zaout
1777      !!-----------------------------------------------------------------------
1778      !
1779      zmax  = MAXVAL( ptab(:,:) , mask= pmask == 1.e0 )
1780      ilocs = MAXLOC( ptab(:,:) , mask= pmask == 1.e0 )
1781      !
1782      ki = ilocs(1) + nimpp - 1
1783      kj = ilocs(2) + njmpp - 1
1784      !
1785      zain(1,:) = zmax
1786      zain(2,:) = ki + 10000. * kj
1787      !
1788      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MAXLOC,MPI_COMM_OPA,ierror)
1789      !
1790      pmax = zaout(1,1)
1791      kj   = INT( zaout(2,1) / 10000.     )
1792      ki   = INT( zaout(2,1) - 10000.* kj )
1793      !
1794   END SUBROUTINE mpp_maxloc2d
1795
1796
1797   SUBROUTINE mpp_maxloc3d( ptab3d, pmask3d, pmax, ki, kj, kk )
1798      !!------------------------------------------------------------------------
1799      !!             ***  routine mpp_maxloc  ***
1800      !!
1801      !! ** Purpose :  Compute the global maximum of an array ptab
1802      !!              and also give its global position
1803      !!
1804      !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC
1805      !!
1806      !!--------------------------------------------------------------------------
1807!FTRANS ptab3d  :I :I :z
1808!FTRANS pmask3d :I :I :z
1809!! DCSE_NEMO: work around a deficiency in ftrans
1810!     REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in   ) ::   ptab3d       ! Local 2D array
1811!     REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in   ) ::   pmask3d      ! Local mask
1812      REAL(wp), INTENT(in   )                          ::   ptab3d(jpi,jpj,jpk)       ! Local 2D array
1813      REAL(wp), INTENT(in   )                          ::   pmask3d(jpi,jpj,jpk)      ! Local mask
1814      REAL(wp)                         , INTENT(  out) ::   pmax         ! Global maximum of ptab
1815      INTEGER                          , INTENT(  out) ::   ki, kj, kk   ! index of maximum in global frame
1816      !!   
1817      REAL(wp) :: zmax   ! local maximum
1818      REAL(wp), DIMENSION(2,1) ::   zain, zaout
1819      INTEGER , DIMENSION(3)   ::   ilocs
1820      INTEGER :: ierror
1821      !!-----------------------------------------------------------------------
1822      !
1823      zmax  = MAXVAL( ptab3d(:,:,:) , mask= pmask3d == 1.e0 )
1824      ilocs = MAXLOC( ptab3d(:,:,:) , mask= pmask3d == 1.e0 )
1825      !
1826!! DCSE_NEMO: Attention!
1827#if defined key_z_first
1828      ki = ilocs(2) + nimpp - 1
1829      kj = ilocs(3) + njmpp - 1
1830      kk = ilocs(1)
1831#else
1832      ki = ilocs(1) + nimpp - 1
1833      kj = ilocs(2) + njmpp - 1
1834      kk = ilocs(3)
1835#endif
1836      !
1837      zain(1,:)=zmax
1838      zain(2,:)=ki+10000.*kj+100000000.*kk
1839      !
1840      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MAXLOC,MPI_COMM_OPA,ierror)
1841      !
1842      pmax = zaout(1,1)
1843      kk   = INT( zaout(2,1) / 100000000. )
1844      kj   = INT( zaout(2,1) - kk * 100000000. ) / 10000
1845      ki   = INT( zaout(2,1) - kk * 100000000. -kj * 10000. )
1846      !
1847   END SUBROUTINE mpp_maxloc3d
1848
1849
1850   SUBROUTINE mppsync()
1851      !!----------------------------------------------------------------------
1852      !!                  ***  routine mppsync  ***
1853      !!                   
1854      !! ** Purpose :   Massively parallel processors, synchroneous
1855      !!
1856      !!-----------------------------------------------------------------------
1857      INTEGER :: ierror
1858      !!-----------------------------------------------------------------------
1859      !
1860      CALL mpi_barrier( mpi_comm_opa, ierror )
1861      !
1862   END SUBROUTINE mppsync
1863
1864
1865   SUBROUTINE mppstop
1866      !!----------------------------------------------------------------------
1867      !!                  ***  routine mppstop  ***
1868      !!                   
1869      !! ** purpose :   Stop massilively parallel processors method
1870      !!
1871      !!----------------------------------------------------------------------
1872      INTEGER ::   info
1873      !!----------------------------------------------------------------------
1874      !
1875      CALL mppsync
1876      CALL mpi_finalize( info )
1877      STOP
1878      !
1879   END SUBROUTINE mppstop
1880
1881
1882   SUBROUTINE mppobc( ptab, kd1, kd2, kl, kk, ktype, kij , kumout)
1883      !!----------------------------------------------------------------------
1884      !!                  ***  routine mppobc  ***
1885      !!
1886      !! ** Purpose :   Message passing management for open boundary
1887      !!     conditions array
1888      !!
1889      !! ** Method  :   Use mppsend and mpprecv function for passing mask
1890      !!       between processors following neighboring subdomains.
1891      !!       domain parameters
1892      !!                    nlci   : first dimension of the local subdomain
1893      !!                    nlcj   : second dimension of the local subdomain
1894      !!                    nbondi : mark for "east-west local boundary"
1895      !!                    nbondj : mark for "north-south local boundary"
1896      !!                    noea   : number for local neighboring processors
1897      !!                    nowe   : number for local neighboring processors
1898      !!                    noso   : number for local neighboring processors
1899      !!                    nono   : number for local neighboring processors
1900      !!
1901      !!----------------------------------------------------------------------
1902      USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released
1903!! DCSE_NEMO: Warning! ztab is also a lib_mpp module variable
1904!     USE wrk_nemo, ONLY:   ztab => wrk_2d_1
1905      USE wrk_nemo, ONLY:   ztab2d => wrk_2d_1
1906      !
1907      INTEGER , INTENT(in   )                     ::   kd1, kd2   ! starting and ending indices
1908      INTEGER , INTENT(in   )                     ::   kl         ! index of open boundary
1909      INTEGER , INTENT(in   )                     ::   kk         ! vertical dimension
1910      INTEGER , INTENT(in   )                     ::   ktype      ! define north/south or east/west cdt
1911      !                                                           !  = 1  north/south  ;  = 2  east/west
1912      INTEGER , INTENT(in   )                     ::   kij        ! horizontal dimension
1913      INTEGER , INTENT(in   )                     ::   kumout     ! ocean.output logical unit
1914      REAL(wp), INTENT(inout), DIMENSION(kij,kk)  ::   ptab       ! variable array
1915      !
1916      INTEGER ::   ji, jj, jk, jl        ! dummy loop indices
1917      INTEGER ::   iipt0, iipt1, ilpt1   ! local integers
1918      INTEGER ::   ijpt0, ijpt1          !   -       -
1919      INTEGER ::   imigr, iihom, ijhom   !   -       -
1920      INTEGER ::   ml_req1, ml_req2, ml_err    ! for key_mpi_isend
1921      INTEGER ::   ml_stat(MPI_STATUS_SIZE)    ! for key_mpi_isend
1922      !!----------------------------------------------------------------------
1923
1924      IF( wrk_in_use(2, 1) ) THEN
1925         WRITE(kumout, cform_err)
1926         WRITE(kumout,*) 'mppobc : requested workspace array unavailable'
1927         CALL mppstop
1928      ENDIF
1929
1930      ! boundary condition initialization
1931      ! ---------------------------------
1932      ztab2d(:,:) = 0.e0
1933      !
1934      IF( ktype==1 ) THEN                                  ! north/south boundaries
1935         iipt0 = MAX( 1, MIN(kd1 - nimpp+1, nlci     ) )
1936         iipt1 = MAX( 0, MIN(kd2 - nimpp+1, nlci - 1 ) )
1937         ilpt1 = MAX( 1, MIN(kd2 - nimpp+1, nlci     ) )
1938         ijpt0 = MAX( 1, MIN(kl  - njmpp+1, nlcj     ) )
1939         ijpt1 = MAX( 0, MIN(kl  - njmpp+1, nlcj - 1 ) )
1940      ELSEIF( ktype==2 ) THEN                              ! east/west boundaries
1941         iipt0 = MAX( 1, MIN(kl  - nimpp+1, nlci     ) )
1942         iipt1 = MAX( 0, MIN(kl  - nimpp+1, nlci - 1 ) )
1943         ijpt0 = MAX( 1, MIN(kd1 - njmpp+1, nlcj     ) )
1944         ijpt1 = MAX( 0, MIN(kd2 - njmpp+1, nlcj - 1 ) )
1945         ilpt1 = MAX( 1, MIN(kd2 - njmpp+1, nlcj     ) )
1946      ELSE
1947         WRITE(kumout, cform_err)
1948         WRITE(kumout,*) 'mppobc : bad ktype'
1949         CALL mppstop
1950      ENDIF
1951     
1952      ! Communication level by level
1953      ! ----------------------------
1954!!gm Remark : this is very time consumming!!!
1955      !                                         ! ------------------------ !
1956      DO jk = 1, kk                             !   Loop over the levels   !
1957         !                                      ! ------------------------ !
1958         !
1959         IF( ktype == 1 ) THEN                               ! north/south boundaries
1960            DO jj = ijpt0, ijpt1
1961               DO ji = iipt0, iipt1
1962                  ztab2d(ji,jj) = ptab(ji,jk)
1963               END DO
1964            END DO
1965         ELSEIF( ktype == 2 ) THEN                           ! east/west boundaries
1966            DO jj = ijpt0, ijpt1
1967               DO ji = iipt0, iipt1
1968                  ztab2d(ji,jj) = ptab(jj,jk)
1969               END DO
1970            END DO
1971         ENDIF
1972
1973
1974         ! 1. East and west directions
1975         ! ---------------------------
1976         !
1977         IF( nbondi /= 2 ) THEN         ! Read Dirichlet lateral conditions
1978            iihom = nlci-nreci
1979            DO jl = 1, jpreci
1980               t2ew(:,jl,1) = ztab2d(jpreci+jl,:)
1981               t2we(:,jl,1) = ztab2d(iihom +jl,:)
1982            END DO
1983         ENDIF
1984         !
1985         !                              ! Migrations
1986         imigr=jpreci*jpj
1987         !
1988         IF( nbondi == -1 ) THEN
1989            CALL mppsend( 2, t2we(1,1,1), imigr, noea, ml_req1 )
1990            CALL mpprecv( 1, t2ew(1,1,2), imigr )
1991            IF(l_isend)   CALL mpi_wait( ml_req1, ml_stat, ml_err )
1992         ELSEIF( nbondi == 0 ) THEN
1993            CALL mppsend( 1, t2ew(1,1,1), imigr, nowe, ml_req1 )
1994            CALL mppsend( 2, t2we(1,1,1), imigr, noea, ml_req2 )
1995            CALL mpprecv( 1, t2ew(1,1,2), imigr )
1996            CALL mpprecv( 2, t2we(1,1,2), imigr )
1997            IF(l_isend)   CALL mpi_wait( ml_req1, ml_stat, ml_err )
1998            IF(l_isend)   CALL mpi_wait( ml_req2, ml_stat, ml_err )
1999         ELSEIF( nbondi == 1 ) THEN
2000            CALL mppsend( 1, t2ew(1,1,1), imigr, nowe, ml_req1 )
2001            CALL mpprecv( 2, t2we(1,1,2), imigr )
2002            IF(l_isend) CALL mpi_wait( ml_req1, ml_stat, ml_err )
2003         ENDIF
2004         !
2005         !                              ! Write Dirichlet lateral conditions
2006         iihom = nlci-jpreci
2007         !
2008         IF( nbondi == 0 .OR. nbondi == 1 ) THEN
2009            DO jl = 1, jpreci
2010               ztab2d(jl,:) = t2we(:,jl,2)
2011            END DO
2012         ENDIF
2013         IF( nbondi == -1 .OR. nbondi == 0 ) THEN
2014            DO jl = 1, jpreci
2015               ztab2d(iihom+jl,:) = t2ew(:,jl,2)
2016            END DO
2017         ENDIF
2018
2019
2020         ! 2. North and south directions
2021         ! -----------------------------
2022         !
2023         IF( nbondj /= 2 ) THEN         ! Read Dirichlet lateral conditions
2024            ijhom = nlcj-nrecj
2025            DO jl = 1, jprecj
2026               t2sn(:,jl,1) = ztab2d(:,ijhom +jl)
2027               t2ns(:,jl,1) = ztab2d(:,jprecj+jl)
2028            END DO
2029         ENDIF
2030         !
2031         !                              ! Migrations
2032         imigr = jprecj * jpi
2033         !
2034         IF( nbondj == -1 ) THEN
2035            CALL mppsend( 4, t2sn(1,1,1), imigr, nono, ml_req1 )
2036            CALL mpprecv( 3, t2ns(1,1,2), imigr )
2037            IF(l_isend) CALL mpi_wait( ml_req1, ml_stat, ml_err )
2038         ELSEIF( nbondj == 0 ) THEN
2039            CALL mppsend( 3, t2ns(1,1,1), imigr, noso, ml_req1 )
2040            CALL mppsend( 4, t2sn(1,1,1), imigr, nono, ml_req2 )
2041            CALL mpprecv( 3, t2ns(1,1,2), imigr )
2042            CALL mpprecv( 4, t2sn(1,1,2), imigr )
2043            IF( l_isend )   CALL mpi_wait( ml_req1, ml_stat, ml_err )
2044            IF( l_isend )   CALL mpi_wait( ml_req2, ml_stat, ml_err )
2045         ELSEIF( nbondj == 1 ) THEN
2046            CALL mppsend( 3, t2ns(1,1,1), imigr, noso, ml_req1 )
2047            CALL mpprecv( 4, t2sn(1,1,2), imigr)
2048            IF( l_isend )   CALL mpi_wait( ml_req1, ml_stat, ml_err )
2049         ENDIF
2050         !
2051         !                              ! Write Dirichlet lateral conditions
2052         ijhom = nlcj - jprecj
2053         IF( nbondj == 0 .OR. nbondj == 1 ) THEN
2054            DO jl = 1, jprecj
2055               ztab2d(:,jl) = t2sn(:,jl,2)
2056            END DO
2057         ENDIF
2058         IF( nbondj == 0 .OR. nbondj == -1 ) THEN
2059            DO jl = 1, jprecj
2060               ztab2d(:,ijhom+jl) = t2ns(:,jl,2)
2061            END DO
2062         ENDIF
2063         IF( ktype==1 .AND. kd1 <= jpi+nimpp-1 .AND. nimpp <= kd2 ) THEN
2064            DO jj = ijpt0, ijpt1            ! north/south boundaries
2065               DO ji = iipt0,ilpt1
2066                  ptab(ji,jk) = ztab2d(ji,jj) 
2067               END DO
2068            END DO
2069         ELSEIF( ktype==2 .AND. kd1 <= jpj+njmpp-1 .AND. njmpp <= kd2 ) THEN
2070            DO jj = ijpt0, ilpt1            ! east/west boundaries
2071               DO ji = iipt0,iipt1
2072                  ptab(jj,jk) = ztab2d(ji,jj) 
2073               END DO
2074            END DO
2075         ENDIF
2076         !
2077      END DO
2078      !
2079      IF( wrk_not_released(2, 1) ) THEN
2080         WRITE(kumout, cform_err)
2081         WRITE(kumout,*) 'mppobc : failed to release workspace array'
2082         CALL mppstop
2083      ENDIF
2084      !
2085   END SUBROUTINE mppobc
2086   
2087
2088   SUBROUTINE mpp_comm_free( kcom )
2089      !!----------------------------------------------------------------------
2090      !!----------------------------------------------------------------------
2091      INTEGER, INTENT(in) ::   kcom
2092      !!
2093      INTEGER :: ierr
2094      !!----------------------------------------------------------------------
2095      !
2096      CALL MPI_COMM_FREE(kcom, ierr)
2097      !
2098   END SUBROUTINE mpp_comm_free
2099
2100
2101   SUBROUTINE mpp_ini_ice( pindic, kumout )
2102      !!----------------------------------------------------------------------
2103      !!               ***  routine mpp_ini_ice  ***
2104      !!
2105      !! ** Purpose :   Initialize special communicator for ice areas
2106      !!      condition together with global variables needed in the ddmpp folding
2107      !!
2108      !! ** Method  : - Look for ice processors in ice routines
2109      !!              - Put their number in nrank_ice
2110      !!              - Create groups for the world processors and the ice processors
2111      !!              - Create a communicator for ice processors
2112      !!
2113      !! ** output
2114      !!      njmppmax = njmpp for northern procs
2115      !!      ndim_rank_ice = number of processors with ice
2116      !!      nrank_ice (ndim_rank_ice) = ice processors
2117      !!      ngrp_world = group ID for the world processors
2118      !!      ngrp_ice = group ID for the ice processors
2119      !!      ncomm_ice = communicator for the ice procs.
2120      !!      n_ice_root = number (in the world) of proc 0 in the ice comm.
2121      !!
2122      !!----------------------------------------------------------------------
2123      INTEGER, INTENT(in) ::   pindic
2124      INTEGER, INTENT(in) ::   kumout   ! ocean.output logical unit
2125      !!
2126      INTEGER :: jjproc
2127      INTEGER :: ii, ierr
2128      INTEGER, ALLOCATABLE, DIMENSION(:) ::   kice
2129      INTEGER, ALLOCATABLE, DIMENSION(:) ::   zwork
2130      !!----------------------------------------------------------------------
2131      !
2132      ! Since this is just an init routine and these arrays are of length jpnij
2133      ! then don't use wrk_nemo module - just allocate and deallocate.
2134      ALLOCATE( kice(jpnij), zwork(jpnij), STAT=ierr )
2135      IF( ierr /= 0 ) THEN
2136         WRITE(kumout, cform_err)
2137         WRITE(kumout,*) 'mpp_ini_ice : failed to allocate 2, 1D arrays (jpnij in length)'
2138         CALL mppstop
2139      ENDIF
2140
2141      ! Look for how many procs with sea-ice
2142      !
2143      kice = 0
2144      DO jjproc = 1, jpnij
2145         IF( jjproc == narea .AND. pindic .GT. 0 )   kice(jjproc) = 1   
2146      END DO
2147      !
2148      zwork = 0
2149      CALL MPI_ALLREDUCE( kice, zwork, jpnij, mpi_integer, mpi_sum, mpi_comm_opa, ierr )
2150      ndim_rank_ice = SUM( zwork )         
2151
2152      ! Allocate the right size to nrank_north
2153      IF( ALLOCATED ( nrank_ice ) )   DEALLOCATE( nrank_ice )
2154      ALLOCATE( nrank_ice(ndim_rank_ice) )
2155      !
2156      ii = 0     
2157      nrank_ice = 0
2158      DO jjproc = 1, jpnij
2159         IF( zwork(jjproc) == 1) THEN
2160            ii = ii + 1
2161            nrank_ice(ii) = jjproc -1 
2162         ENDIF
2163      END DO
2164
2165      ! Create the world group
2166      CALL MPI_COMM_GROUP( mpi_comm_opa, ngrp_world, ierr )
2167
2168      ! Create the ice group from the world group
2169      CALL MPI_GROUP_INCL( ngrp_world, ndim_rank_ice, nrank_ice, ngrp_ice, ierr )
2170
2171      ! Create the ice communicator , ie the pool of procs with sea-ice
2172      CALL MPI_COMM_CREATE( mpi_comm_opa, ngrp_ice, ncomm_ice, ierr )
2173
2174      ! Find proc number in the world of proc 0 in the north
2175      ! The following line seems to be useless, we just comment & keep it as reminder
2176      ! CALL MPI_GROUP_TRANSLATE_RANKS(ngrp_ice,1,0,ngrp_world,n_ice_root,ierr)
2177      !
2178      DEALLOCATE(kice, zwork)
2179      !
2180   END SUBROUTINE mpp_ini_ice
2181
2182
2183   SUBROUTINE mpp_ini_znl( kumout )
2184      !!----------------------------------------------------------------------
2185      !!               ***  routine mpp_ini_znl  ***
2186      !!
2187      !! ** Purpose :   Initialize special communicator for computing zonal sum
2188      !!
2189      !! ** Method  : - Look for processors in the same row
2190      !!              - Put their number in nrank_znl
2191      !!              - Create group for the znl processors
2192      !!              - Create a communicator for znl processors
2193      !!              - Determine if processor should write znl files
2194      !!
2195      !! ** output
2196      !!      ndim_rank_znl = number of processors on the same row
2197      !!      ngrp_znl = group ID for the znl processors
2198      !!      ncomm_znl = communicator for the ice procs.
2199      !!      n_znl_root = number (in the world) of proc 0 in the ice comm.
2200      !!
2201      !!----------------------------------------------------------------------
2202      INTEGER, INTENT(in) ::   kumout   ! ocean.output logical units
2203      !
2204      INTEGER :: jproc      ! dummy loop integer
2205      INTEGER :: ierr, ii   ! local integer
2206      INTEGER, ALLOCATABLE, DIMENSION(:) ::   kwork
2207      !!----------------------------------------------------------------------
2208      !-$$     WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ngrp_world     : ', ngrp_world
2209      !-$$     WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - mpi_comm_world : ', mpi_comm_world
2210      !-$$     WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - mpi_comm_opa   : ', mpi_comm_opa
2211      !
2212      ALLOCATE( kwork(jpnij), STAT=ierr )
2213      IF( ierr /= 0 ) THEN
2214         WRITE(kumout, cform_err)
2215         WRITE(kumout,*) 'mpp_ini_znl : failed to allocate 1D array of length jpnij'
2216         CALL mppstop
2217      ENDIF
2218
2219      IF( jpnj == 1 ) THEN
2220         ngrp_znl  = ngrp_world
2221         ncomm_znl = mpi_comm_opa
2222      ELSE
2223         !
2224         CALL MPI_ALLGATHER ( njmpp, 1, mpi_integer, kwork, 1, mpi_integer, mpi_comm_opa, ierr )
2225         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - kwork pour njmpp : ', kwork
2226         !-$$        CALL flush(numout)
2227         !
2228         ! Count number of processors on the same row
2229         ndim_rank_znl = 0
2230         DO jproc=1,jpnij
2231            IF ( kwork(jproc) == njmpp ) THEN
2232               ndim_rank_znl = ndim_rank_znl + 1
2233            ENDIF
2234         END DO
2235         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ndim_rank_znl : ', ndim_rank_znl
2236         !-$$        CALL flush(numout)
2237         ! Allocate the right size to nrank_znl
2238         IF (ALLOCATED (nrank_znl)) DEALLOCATE(nrank_znl)
2239         ALLOCATE(nrank_znl(ndim_rank_znl))
2240         ii = 0     
2241         nrank_znl (:) = 0
2242         DO jproc=1,jpnij
2243            IF ( kwork(jproc) == njmpp) THEN
2244               ii = ii + 1
2245               nrank_znl(ii) = jproc -1 
2246            ENDIF
2247         END DO
2248         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - nrank_znl : ', nrank_znl
2249         !-$$        CALL flush(numout)
2250
2251         ! Create the opa group
2252         CALL MPI_COMM_GROUP(mpi_comm_opa,ngrp_opa,ierr)
2253         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ngrp_opa : ', ngrp_opa
2254         !-$$        CALL flush(numout)
2255
2256         ! Create the znl group from the opa group
2257         CALL MPI_GROUP_INCL  ( ngrp_opa, ndim_rank_znl, nrank_znl, ngrp_znl, ierr )
2258         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ngrp_znl ', ngrp_znl
2259         !-$$        CALL flush(numout)
2260
2261         ! Create the znl communicator from the opa communicator, ie the pool of procs in the same row
2262         CALL MPI_COMM_CREATE ( mpi_comm_opa, ngrp_znl, ncomm_znl, ierr )
2263         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ncomm_znl ', ncomm_znl
2264         !-$$        CALL flush(numout)
2265         !
2266      END IF
2267
2268      ! Determines if processor if the first (starting from i=1) on the row
2269      IF ( jpni == 1 ) THEN
2270         l_znl_root = .TRUE.
2271      ELSE
2272         l_znl_root = .FALSE.
2273         kwork (1) = nimpp
2274         CALL mpp_min ( kwork(1), kcom = ncomm_znl)
2275         IF ( nimpp == kwork(1)) l_znl_root = .TRUE.
2276      END IF
2277
2278      DEALLOCATE(kwork)
2279
2280   END SUBROUTINE mpp_ini_znl
2281
2282
2283   SUBROUTINE mpp_ini_north
2284      !!----------------------------------------------------------------------
2285      !!               ***  routine mpp_ini_north  ***
2286      !!
2287      !! ** Purpose :   Initialize special communicator for north folding
2288      !!      condition together with global variables needed in the mpp folding
2289      !!
2290      !! ** Method  : - Look for northern processors
2291      !!              - Put their number in nrank_north
2292      !!              - Create groups for the world processors and the north processors
2293      !!              - Create a communicator for northern processors
2294      !!
2295      !! ** output
2296      !!      njmppmax = njmpp for northern procs
2297      !!      ndim_rank_north = number of processors in the northern line
2298      !!      nrank_north (ndim_rank_north) = number  of the northern procs.
2299      !!      ngrp_world = group ID for the world processors
2300      !!      ngrp_north = group ID for the northern processors
2301      !!      ncomm_north = communicator for the northern procs.
2302      !!      north_root = number (in the world) of proc 0 in the northern comm.
2303      !!
2304      !!----------------------------------------------------------------------
2305      INTEGER ::   ierr
2306      INTEGER ::   jjproc
2307      INTEGER ::   ii, ji
2308      !!----------------------------------------------------------------------
2309      !
2310      njmppmax = MAXVAL( njmppt )
2311      !
2312      ! Look for how many procs on the northern boundary
2313      ndim_rank_north = 0
2314      DO jjproc = 1, jpnij
2315         IF( njmppt(jjproc) == njmppmax )   ndim_rank_north = ndim_rank_north + 1
2316      END DO
2317      !
2318      ! Allocate the right size to nrank_north
2319      IF (ALLOCATED (nrank_north)) DEALLOCATE(nrank_north)
2320      ALLOCATE( nrank_north(ndim_rank_north) )
2321
2322      ! Fill the nrank_north array with proc. number of northern procs.
2323      ! Note : the rank start at 0 in MPI
2324      ii = 0
2325      DO ji = 1, jpnij
2326         IF ( njmppt(ji) == njmppmax   ) THEN
2327            ii=ii+1
2328            nrank_north(ii)=ji-1
2329         END IF
2330      END DO
2331      !
2332      ! create the world group
2333      CALL MPI_COMM_GROUP( mpi_comm_opa, ngrp_world, ierr )
2334      !
2335      ! Create the North group from the world group
2336      CALL MPI_GROUP_INCL( ngrp_world, ndim_rank_north, nrank_north, ngrp_north, ierr )
2337      !
2338      ! Create the North communicator , ie the pool of procs in the north group
2339      CALL MPI_COMM_CREATE( mpi_comm_opa, ngrp_north, ncomm_north, ierr )
2340      !
2341   END SUBROUTINE mpp_ini_north
2342
2343
2344   SUBROUTINE mpp_lbc_north_3d( pt3d, cd_type, psgn )
2345      !!---------------------------------------------------------------------
2346      !!                   ***  routine mpp_lbc_north_3d  ***
2347      !!
2348      !! ** Purpose :   Ensure proper north fold horizontal bondary condition
2349      !!              in mpp configuration in case of jpn1 > 1
2350      !!
2351      !! ** Method  :   North fold condition and mpp with more than one proc
2352      !!              in i-direction require a specific treatment. We gather
2353      !!              the 4 northern lines of the global domain on 1 processor
2354      !!              and apply lbc north-fold on this sub array. Then we
2355      !!              scatter the north fold array back to the processors.
2356      !!
2357      !!----------------------------------------------------------------------
2358!FTRANS pt3d :I :I :z
2359!! DCSE_NEMO: work around a deficiency in ftrans
2360!     REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pt3d      ! 3D array on which the b.c. is applied
2361      REAL(wp),                         INTENT(inout) ::   pt3d(jpi,jpj,jpk)
2362      CHARACTER(len=1)                , INTENT(in   ) ::   cd_type   ! nature of pt3d grid-points
2363      !                                                              !   = T ,  U , V , F or W  gridpoints
2364      REAL(wp)                        , INTENT(in   ) ::   psgn      ! = -1. the sign change across the north fold
2365      !!                                                             ! =  1. , the sign is kept
2366      INTEGER ::   ji, jj, jr
2367      INTEGER ::   ierr, itaille, ildi, ilei, iilb
2368      INTEGER ::   ijpj, ijpjm1, ij, iproc
2369      !!----------------------------------------------------------------------
2370      !   
2371      ijpj   = 4
2372      ijpjm1 = 3
2373      ztab(:,:,:) = 0.e0
2374      !
2375      DO jj = nlcj - ijpj +1, nlcj          ! put in znorthloc the last 4 jlines of pt3d
2376         ij = jj - nlcj + ijpj
2377         znorthloc(:,ij,:) = pt3d(:,jj,:)
2378      END DO
2379      !
2380      !                                     ! Build in procs of ncomm_north the znorthgloio
2381      itaille = jpi * jpk * ijpj
2382      CALL MPI_ALLGATHER( znorthloc  , itaille, MPI_DOUBLE_PRECISION,                &
2383         &                znorthgloio, itaille, MPI_DOUBLE_PRECISION, ncomm_north, ierr )
2384      !
2385      !                                     ! recover the global north array
2386      DO jr = 1, ndim_rank_north
2387         iproc = nrank_north(jr) + 1
2388         ildi  = nldit (iproc)
2389         ilei  = nleit (iproc)
2390         iilb  = nimppt(iproc)
2391         DO jj = 1, 4
2392            DO ji = ildi, ilei
2393               ztab(ji+iilb-1,jj,:) = znorthgloio(ji,jj,:,jr)
2394            END DO
2395         END DO
2396      END DO
2397      !
2398      CALL lbc_nfd( ztab, cd_type, psgn )   ! North fold boundary condition
2399      !
2400      DO jj = nlcj-ijpj+1, nlcj             ! Scatter back to pt3d
2401         ij = jj - nlcj + ijpj
2402         DO ji= 1, nlci
2403            pt3d(ji,jj,:) = ztab(ji+nimpp-1,ij,:)
2404         END DO
2405      END DO
2406      !
2407   END SUBROUTINE mpp_lbc_north_3d
2408
2409
2410   SUBROUTINE mpp_lbc_north_2d( pt2d, cd_type, psgn)
2411      !!---------------------------------------------------------------------
2412      !!                   ***  routine mpp_lbc_north_2d  ***
2413      !!
2414      !! ** Purpose :   Ensure proper north fold horizontal bondary condition
2415      !!              in mpp configuration in case of jpn1 > 1 (for 2d array )
2416      !!
2417      !! ** Method  :   North fold condition and mpp with more than one proc
2418      !!              in i-direction require a specific treatment. We gather
2419      !!              the 4 northern lines of the global domain on 1 processor
2420      !!              and apply lbc north-fold on this sub array. Then we
2421      !!              scatter the north fold array back to the processors.
2422      !!
2423      !!----------------------------------------------------------------------
2424      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   pt2d      ! 3D array on which the b.c. is applied
2425      CHARACTER(len=1)            , INTENT(in   ) ::   cd_type   ! nature of pt3d grid-points
2426      !                                                          !   = T ,  U , V , F or W  gridpoints
2427      REAL(wp)                    , INTENT(in   ) ::   psgn      ! = -1. the sign change across the north fold
2428      !!                                                             ! =  1. , the sign is kept
2429      INTEGER ::   ji, jj, jr
2430      INTEGER ::   ierr, itaille, ildi, ilei, iilb
2431      INTEGER ::   ijpj, ijpjm1, ij, iproc
2432      !!----------------------------------------------------------------------
2433      !
2434      ijpj   = 4
2435      ijpjm1 = 3
2436      ztab_2d(:,:) = 0.e0
2437      !
2438      DO jj = nlcj-ijpj+1, nlcj             ! put in znorthloc_2d the last 4 jlines of pt2d
2439         ij = jj - nlcj + ijpj
2440         znorthloc_2d(:,ij) = pt2d(:,jj)
2441      END DO
2442
2443      !                                     ! Build in procs of ncomm_north the znorthgloio_2d
2444      itaille = jpi * ijpj
2445      CALL MPI_ALLGATHER( znorthloc_2d  , itaille, MPI_DOUBLE_PRECISION,        &
2446         &                znorthgloio_2d, itaille, MPI_DOUBLE_PRECISION, ncomm_north, ierr )
2447      !
2448      DO jr = 1, ndim_rank_north            ! recover the global north array
2449         iproc = nrank_north(jr) + 1
2450         ildi=nldit (iproc)
2451         ilei=nleit (iproc)
2452         iilb=nimppt(iproc)
2453         DO jj = 1, 4
2454            DO ji = ildi, ilei
2455               ztab_2d(ji+iilb-1,jj) = znorthgloio_2d(ji,jj,jr)
2456            END DO
2457         END DO
2458      END DO
2459      !
2460      CALL lbc_nfd( ztab_2d, cd_type, psgn )   ! North fold boundary condition
2461      !
2462      !
2463      DO jj = nlcj-ijpj+1, nlcj             ! Scatter back to pt2d
2464         ij = jj - nlcj + ijpj
2465         DO ji = 1, nlci
2466            pt2d(ji,jj) = ztab_2d(ji+nimpp-1,ij)
2467         END DO
2468      END DO
2469      !
2470   END SUBROUTINE mpp_lbc_north_2d
2471
2472
2473   SUBROUTINE mpp_lbc_north_e( pt2d, cd_type, psgn)
2474      !!---------------------------------------------------------------------
2475      !!                   ***  routine mpp_lbc_north_2d  ***
2476      !!
2477      !! ** Purpose :   Ensure proper north fold horizontal bondary condition
2478      !!              in mpp configuration in case of jpn1 > 1 and for 2d
2479      !!              array with outer extra halo
2480      !!
2481      !! ** Method  :   North fold condition and mpp with more than one proc
2482      !!              in i-direction require a specific treatment. We gather
2483      !!              the 4+2*jpr2dj northern lines of the global domain on 1
2484      !!              processor and apply lbc north-fold on this sub array.
2485      !!              Then we scatter the north fold array back to the processors.
2486      !!
2487      !!----------------------------------------------------------------------
2488      REAL(wp), DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj), INTENT(inout) ::   pt2d     ! 2D array with extra halo
2489      CHARACTER(len=1)                                            , INTENT(in   ) ::   cd_type  ! nature of pt3d grid-points
2490      !                                                                                         !   = T ,  U , V , F or W -points
2491      REAL(wp)                                                    , INTENT(in   ) ::   psgn     ! = -1. the sign change across the 
2492      !!                                                                                        ! north fold, =  1. otherwise
2493      INTEGER ::   ji, jj, jr
2494      INTEGER ::   ierr, itaille, ildi, ilei, iilb
2495      INTEGER ::   ijpj, ij, iproc
2496      !!----------------------------------------------------------------------
2497      !
2498      ijpj=4
2499      ztab_e(:,:) = 0.e0
2500
2501      ij=0
2502      ! put in znorthloc_e the last 4 jlines of pt2d
2503      DO jj = nlcj - ijpj + 1 - jpr2dj, nlcj +jpr2dj
2504         ij = ij + 1
2505         DO ji = 1, jpi
2506            znorthloc_e(ji,ij)=pt2d(ji,jj)
2507         END DO
2508      END DO
2509      !
2510      itaille = jpi * ( ijpj + 2 * jpr2dj )
2511      CALL MPI_ALLGATHER( znorthloc_e(1,1)  , itaille, MPI_DOUBLE_PRECISION,    &
2512         &                znorthgloio_e(1,1,1), itaille, MPI_DOUBLE_PRECISION, ncomm_north, ierr )
2513      !
2514      DO jr = 1, ndim_rank_north            ! recover the global north array
2515         iproc = nrank_north(jr) + 1
2516         ildi = nldit (iproc)
2517         ilei = nleit (iproc)
2518         iilb = nimppt(iproc)
2519         DO jj = 1, ijpj+2*jpr2dj
2520            DO ji = ildi, ilei
2521               ztab_e(ji+iilb-1,jj) = znorthgloio_e(ji,jj,jr)
2522            END DO
2523         END DO
2524      END DO
2525
2526
2527      ! 2. North-Fold boundary conditions
2528      ! ----------------------------------
2529      CALL lbc_nfd( ztab_e(:,:), cd_type, psgn, pr2dj = jpr2dj )
2530
2531      ij = jpr2dj
2532      !! Scatter back to pt2d
2533      DO jj = nlcj - ijpj + 1 , nlcj +jpr2dj
2534      ij  = ij +1 
2535         DO ji= 1, nlci
2536            pt2d(ji,jj) = ztab_e(ji+nimpp-1,ij)
2537         END DO
2538      END DO
2539      !
2540   END SUBROUTINE mpp_lbc_north_e
2541
2542
2543   SUBROUTINE mpi_init_opa( ldtxt, ksft, code )
2544      !!---------------------------------------------------------------------
2545      !!                   ***  routine mpp_init.opa  ***
2546      !!
2547      !! ** Purpose :: export and attach a MPI buffer for bsend
2548      !!
2549      !! ** Method  :: define buffer size in namelist, if 0 no buffer attachment
2550      !!            but classical mpi_init
2551      !!
2552      !! History :: 01/11 :: IDRIS initial version for IBM only 
2553      !!            08/04 :: R. Benshila, generalisation
2554      !!---------------------------------------------------------------------
2555      CHARACTER(len=*),DIMENSION(:), INTENT(  out) ::   ldtxt 
2556      INTEGER                      , INTENT(inout) ::   ksft
2557      INTEGER                      , INTENT(  out) ::   code
2558      INTEGER                                      ::   ierr, ji
2559      LOGICAL                                      ::   mpi_was_called
2560      !!---------------------------------------------------------------------
2561      !
2562      CALL mpi_initialized( mpi_was_called, code )      ! MPI initialization
2563      IF ( code /= MPI_SUCCESS ) THEN
2564         DO ji = 1, SIZE(ldtxt) 
2565            IF( TRIM(ldtxt(ji)) /= '' )   WRITE(*,*) ldtxt(ji)      ! control print of mynode
2566         END DO         
2567         WRITE(*, cform_err)
2568         WRITE(*, *) ' lib_mpp: Error in routine mpi_initialized'
2569         CALL mpi_abort( mpi_comm_world, code, ierr )
2570      ENDIF
2571      !
2572      IF( .NOT. mpi_was_called ) THEN
2573         CALL mpi_init( code )
2574         CALL mpi_comm_dup( mpi_comm_world, mpi_comm_opa, code )
2575         IF ( code /= MPI_SUCCESS ) THEN
2576            DO ji = 1, SIZE(ldtxt) 
2577               IF( TRIM(ldtxt(ji)) /= '' )   WRITE(*,*) ldtxt(ji)      ! control print of mynode
2578            END DO
2579            WRITE(*, cform_err)
2580            WRITE(*, *) ' lib_mpp: Error in routine mpi_comm_dup'
2581            CALL mpi_abort( mpi_comm_world, code, ierr )
2582         ENDIF
2583      ENDIF
2584      !
2585      IF( nn_buffer > 0 ) THEN
2586         WRITE(ldtxt(ksft),*) 'mpi_bsend, buffer allocation of  : ', nn_buffer   ;   ksft = ksft + 1
2587         ! Buffer allocation and attachment
2588         ALLOCATE( tampon(nn_buffer), stat = ierr )
2589         IF( ierr /= 0 ) THEN
2590            DO ji = 1, SIZE(ldtxt) 
2591               IF( TRIM(ldtxt(ji)) /= '' )   WRITE(*,*) ldtxt(ji)      ! control print of mynode
2592            END DO
2593            WRITE(*, cform_err)
2594            WRITE(*, *) ' lib_mpp: Error in ALLOCATE', ierr
2595            CALL mpi_abort( mpi_comm_world, code, ierr )
2596         END IF
2597         CALL mpi_buffer_attach( tampon, nn_buffer, code )
2598      ENDIF
2599      !
2600   END SUBROUTINE mpi_init_opa
2601
2602#if defined key_mpp_rep
2603   SUBROUTINE DDPDD_MPI (ydda, yddb, ilen, itype)
2604      !!---------------------------------------------------------------------
2605      !!   Routine DDPDD_MPI: used by reduction operator MPI_SUMDD
2606      !!
2607      !!   Modification of original codes written by David H. Bailey
2608      !!   This subroutine computes yddb(i) = ydda(i)+yddb(i)
2609      !!---------------------------------------------------------------------
2610      INTEGER, INTENT(in)                         :: ilen, itype
2611      COMPLEX(wp), DIMENSION(ilen), INTENT(in)     :: ydda
2612      COMPLEX(wp), DIMENSION(ilen), INTENT(inout)  :: yddb
2613      !
2614      REAL(wp) :: zerr, zt1, zt2    ! local work variables
2615      INTEGER :: ji, ztmp           ! local scalar
2616
2617      ztmp = itype   ! avoid compilation warning
2618
2619      DO ji=1,ilen
2620      ! Compute ydda + yddb using Knuth's trick.
2621         zt1  = real(ydda(ji)) + real(yddb(ji))
2622         zerr = zt1 - real(ydda(ji))
2623         zt2  = ((real(yddb(ji)) - zerr) + (real(ydda(ji)) - (zt1 - zerr))) &
2624                + aimag(ydda(ji)) + aimag(yddb(ji))
2625
2626         ! The result is zt1 + zt2, after normalization.
2627         yddb(ji) = cmplx ( zt1 + zt2, zt2 - ((zt1 + zt2) - zt1),wp )
2628      END DO
2629
2630   END SUBROUTINE DDPDD_MPI
2631#endif
2632
2633#else
2634   !!----------------------------------------------------------------------
2635   !!   Default case:            Dummy module        share memory computing
2636   !!----------------------------------------------------------------------
2637   USE in_out_manager
2638
2639   INTERFACE mpp_sum
2640      MODULE PROCEDURE mpp_sum_a2s, mpp_sum_as, mpp_sum_ai, mpp_sum_s, mpp_sum_i
2641   END INTERFACE
2642   INTERFACE mpp_max
2643      MODULE PROCEDURE mppmax_a_int, mppmax_int, mppmax_a_real, mppmax_real
2644   END INTERFACE
2645   INTERFACE mpp_min
2646      MODULE PROCEDURE mppmin_a_int, mppmin_int, mppmin_a_real, mppmin_real
2647   END INTERFACE
2648   INTERFACE mppobc
2649      MODULE PROCEDURE mppobc_1d, mppobc_2d, mppobc_3d, mppobc_4d
2650   END INTERFACE
2651   INTERFACE mpp_minloc
2652      MODULE PROCEDURE mpp_minloc2d ,mpp_minloc3d
2653   END INTERFACE
2654   INTERFACE mpp_maxloc
2655      MODULE PROCEDURE mpp_maxloc2d ,mpp_maxloc3d
2656   END INTERFACE
2657
2658   LOGICAL, PUBLIC, PARAMETER ::   lk_mpp = .FALSE.      !: mpp flag
2659   INTEGER :: ncomm_ice
2660   !!----------------------------------------------------------------------
2661CONTAINS
2662
2663   INTEGER FUNCTION lib_mpp_alloc(kumout)          ! Dummy function
2664      INTEGER, INTENT(in) ::   kumout
2665      lib_mpp_alloc = 0
2666   END FUNCTION lib_mpp_alloc
2667
2668   FUNCTION mynode( ldtxt, kumnam, kstop, localComm ) RESULT (function_value)
2669      INTEGER, OPTIONAL            , INTENT(in   ) ::   localComm
2670      CHARACTER(len=*),DIMENSION(:) ::   ldtxt 
2671      INTEGER ::   kumnam, kstop
2672      IF( PRESENT( localComm ) .OR. .NOT.PRESENT( localComm ) )   function_value = 0
2673      IF( .FALSE. )   ldtxt(:) = 'never done'
2674   END FUNCTION mynode
2675
2676   SUBROUTINE mppsync                       ! Dummy routine
2677   END SUBROUTINE mppsync
2678
2679   SUBROUTINE mpp_sum_as( parr, kdim, kcom )      ! Dummy routine
2680      REAL   , DIMENSION(:) :: parr
2681      INTEGER               :: kdim
2682      INTEGER, OPTIONAL     :: kcom 
2683      WRITE(*,*) 'mpp_sum_as: You should not have seen this print! error?', kdim, parr(1), kcom
2684   END SUBROUTINE mpp_sum_as
2685
2686   SUBROUTINE mpp_sum_a2s( parr, kdim, kcom )      ! Dummy routine
2687      REAL   , DIMENSION(:,:) :: parr
2688      INTEGER               :: kdim
2689      INTEGER, OPTIONAL     :: kcom 
2690      WRITE(*,*) 'mpp_sum_a2s: You should not have seen this print! error?', kdim, parr(1,1), kcom
2691   END SUBROUTINE mpp_sum_a2s
2692
2693   SUBROUTINE mpp_sum_ai( karr, kdim, kcom )      ! Dummy routine
2694      INTEGER, DIMENSION(:) :: karr
2695      INTEGER               :: kdim
2696      INTEGER, OPTIONAL     :: kcom 
2697      WRITE(*,*) 'mpp_sum_ai: You should not have seen this print! error?', kdim, karr(1), kcom
2698   END SUBROUTINE mpp_sum_ai
2699
2700   SUBROUTINE mpp_sum_s( psca, kcom )            ! Dummy routine
2701      REAL                  :: psca
2702      INTEGER, OPTIONAL     :: kcom 
2703      WRITE(*,*) 'mpp_sum_s: You should not have seen this print! error?', psca, kcom
2704   END SUBROUTINE mpp_sum_s
2705
2706   SUBROUTINE mpp_sum_i( kint, kcom )            ! Dummy routine
2707      integer               :: kint
2708      INTEGER, OPTIONAL     :: kcom 
2709      WRITE(*,*) 'mpp_sum_i: You should not have seen this print! error?', kint, kcom
2710   END SUBROUTINE mpp_sum_i
2711
2712   SUBROUTINE mppmax_a_real( parr, kdim, kcom )
2713      REAL   , DIMENSION(:) :: parr
2714      INTEGER               :: kdim
2715      INTEGER, OPTIONAL     :: kcom 
2716      WRITE(*,*) 'mppmax_a_real: You should not have seen this print! error?', kdim, parr(1), kcom
2717   END SUBROUTINE mppmax_a_real
2718
2719   SUBROUTINE mppmax_real( psca, kcom )
2720      REAL                  :: psca
2721      INTEGER, OPTIONAL     :: kcom 
2722      WRITE(*,*) 'mppmax_real: You should not have seen this print! error?', psca, kcom
2723   END SUBROUTINE mppmax_real
2724
2725   SUBROUTINE mppmin_a_real( parr, kdim, kcom )
2726      REAL   , DIMENSION(:) :: parr
2727      INTEGER               :: kdim
2728      INTEGER, OPTIONAL     :: kcom 
2729      WRITE(*,*) 'mppmin_a_real: You should not have seen this print! error?', kdim, parr(1), kcom
2730   END SUBROUTINE mppmin_a_real
2731
2732   SUBROUTINE mppmin_real( psca, kcom )
2733      REAL                  :: psca
2734      INTEGER, OPTIONAL     :: kcom 
2735      WRITE(*,*) 'mppmin_real: You should not have seen this print! error?', psca, kcom
2736   END SUBROUTINE mppmin_real
2737
2738   SUBROUTINE mppmax_a_int( karr, kdim ,kcom)
2739      INTEGER, DIMENSION(:) :: karr
2740      INTEGER               :: kdim
2741      INTEGER, OPTIONAL     :: kcom 
2742      WRITE(*,*) 'mppmax_a_int: You should not have seen this print! error?', kdim, karr(1), kcom
2743   END SUBROUTINE mppmax_a_int
2744
2745   SUBROUTINE mppmax_int( kint, kcom)
2746      INTEGER               :: kint
2747      INTEGER, OPTIONAL     :: kcom 
2748      WRITE(*,*) 'mppmax_int: You should not have seen this print! error?', kint, kcom
2749   END SUBROUTINE mppmax_int
2750
2751   SUBROUTINE mppmin_a_int( karr, kdim, kcom )
2752      INTEGER, DIMENSION(:) :: karr
2753      INTEGER               :: kdim
2754      INTEGER, OPTIONAL     :: kcom 
2755      WRITE(*,*) 'mppmin_a_int: You should not have seen this print! error?', kdim, karr(1), kcom
2756   END SUBROUTINE mppmin_a_int
2757
2758   SUBROUTINE mppmin_int( kint, kcom )
2759      INTEGER               :: kint
2760      INTEGER, OPTIONAL     :: kcom 
2761      WRITE(*,*) 'mppmin_int: You should not have seen this print! error?', kint, kcom
2762   END SUBROUTINE mppmin_int
2763
2764   SUBROUTINE mppobc_1d( parr, kd1, kd2, kl, kk, ktype, kij, knum )
2765      INTEGER  ::   kd1, kd2, kl , kk, ktype, kij, knum
2766      REAL, DIMENSION(:) ::   parr           ! variable array
2767      WRITE(*,*) 'mppobc: You should not have seen this print! error?', parr(1), kd1, kd2, kl, kk, ktype, kij, knum
2768   END SUBROUTINE mppobc_1d
2769
2770   SUBROUTINE mppobc_2d( parr, kd1, kd2, kl, kk, ktype, kij, knum )
2771      INTEGER  ::   kd1, kd2, kl , kk, ktype, kij, knum
2772      REAL, DIMENSION(:,:) ::   parr           ! variable array
2773      WRITE(*,*) 'mppobc: You should not have seen this print! error?', parr(1,1), kd1, kd2, kl, kk, ktype, kij, knum
2774   END SUBROUTINE mppobc_2d
2775
2776   SUBROUTINE mppobc_3d( parr, kd1, kd2, kl, kk, ktype, kij, knum )
2777      INTEGER  ::   kd1, kd2, kl , kk, ktype, kij, knum
2778      REAL, DIMENSION(:,:,:) ::   parr           ! variable array
2779      WRITE(*,*) 'mppobc: You should not have seen this print! error?', parr(1,1,1), kd1, kd2, kl, kk, ktype, kij, knum
2780   END SUBROUTINE mppobc_3d
2781
2782   SUBROUTINE mppobc_4d( parr, kd1, kd2, kl, kk, ktype, kij, knum )
2783      INTEGER  ::   kd1, kd2, kl , kk, ktype, kij, knum
2784      REAL, DIMENSION(:,:,:,:) ::   parr           ! variable array
2785      WRITE(*,*) 'mppobc: You should not have seen this print! error?', parr(1,1,1,1), kd1, kd2, kl, kk, ktype, kij, knum
2786   END SUBROUTINE mppobc_4d
2787
2788   SUBROUTINE mpp_minloc2d( ptab, pmask, pmin, ki, kj )
2789      REAL                   :: pmin
2790      REAL , DIMENSION (:,:) :: ptab, pmask
2791      INTEGER :: ki, kj
2792      WRITE(*,*) 'mpp_minloc2d: You should not have seen this print! error?', pmin, ki, kj, ptab(1,1), pmask(1,1)
2793   END SUBROUTINE mpp_minloc2d
2794
2795   SUBROUTINE mpp_minloc3d( ptab, pmask, pmin, ki, kj, kk )
2796      REAL                     :: pmin
2797      REAL , DIMENSION (:,:,:) :: ptab, pmask
2798      INTEGER :: ki, kj, kk
2799      WRITE(*,*) 'mpp_minloc3d: You should not have seen this print! error?', pmin, ki, kj, kk, ptab(1,1,1), pmask(1,1,1)
2800   END SUBROUTINE mpp_minloc3d
2801
2802   SUBROUTINE mpp_maxloc2d( ptab, pmask, pmax, ki, kj )
2803      REAL                   :: pmax
2804      REAL , DIMENSION (:,:) :: ptab, pmask
2805      INTEGER :: ki, kj
2806      WRITE(*,*) 'mpp_maxloc2d: You should not have seen this print! error?', pmax, ki, kj, ptab(1,1), pmask(1,1)
2807   END SUBROUTINE mpp_maxloc2d
2808
2809   SUBROUTINE mpp_maxloc3d( ptab, pmask, pmax, ki, kj, kk )
2810      REAL                     :: pmax
2811      REAL , DIMENSION (:,:,:) :: ptab, pmask
2812      INTEGER :: ki, kj, kk
2813      WRITE(*,*) 'mpp_maxloc3d: You should not have seen this print! error?', pmax, ki, kj, kk, ptab(1,1,1), pmask(1,1,1)
2814   END SUBROUTINE mpp_maxloc3d
2815
2816   SUBROUTINE mppstop
2817      WRITE(*,*) 'mppstop: You should not have seen this print! error?'
2818   END SUBROUTINE mppstop
2819
2820   SUBROUTINE mpp_ini_ice( kcom, knum )
2821      INTEGER :: kcom, knum
2822      WRITE(*,*) 'mpp_ini_ice: You should not have seen this print! error?', kcom, knum
2823   END SUBROUTINE mpp_ini_ice
2824
2825   SUBROUTINE mpp_ini_znl( knum )
2826      INTEGER :: knum
2827      WRITE(*,*) 'mpp_ini_znl: You should not have seen this print! error?', knum
2828   END SUBROUTINE mpp_ini_znl
2829
2830   SUBROUTINE mpp_comm_free( kcom )
2831      INTEGER :: kcom
2832      WRITE(*,*) 'mpp_comm_free: You should not have seen this print! error?', kcom
2833   END SUBROUTINE mpp_comm_free
2834#endif
2835
2836   !!----------------------------------------------------------------------
2837   !!   All cases:         ctl_stop, ctl_warn, get_unit, ctl_opn   routines
2838   !!----------------------------------------------------------------------
2839
2840   SUBROUTINE ctl_stop( cd1, cd2, cd3, cd4, cd5 ,   &
2841      &                 cd6, cd7, cd8, cd9, cd10 )
2842      !!----------------------------------------------------------------------
2843      !!                  ***  ROUTINE  stop_opa  ***
2844      !!
2845      !! ** Purpose :   print in ocean.outpput file a error message and
2846      !!                increment the error number (nstop) by one.
2847      !!----------------------------------------------------------------------
2848      CHARACTER(len=*), INTENT(in), OPTIONAL ::  cd1, cd2, cd3, cd4, cd5
2849      CHARACTER(len=*), INTENT(in), OPTIONAL ::  cd6, cd7, cd8, cd9, cd10
2850      !!----------------------------------------------------------------------
2851      !
2852      nstop = nstop + 1 
2853      IF(lwp) THEN
2854         WRITE(numout,cform_err)
2855         IF( PRESENT(cd1 ) )   WRITE(numout,*) cd1
2856         IF( PRESENT(cd2 ) )   WRITE(numout,*) cd2
2857         IF( PRESENT(cd3 ) )   WRITE(numout,*) cd3
2858         IF( PRESENT(cd4 ) )   WRITE(numout,*) cd4
2859         IF( PRESENT(cd5 ) )   WRITE(numout,*) cd5
2860         IF( PRESENT(cd6 ) )   WRITE(numout,*) cd6
2861         IF( PRESENT(cd7 ) )   WRITE(numout,*) cd7
2862         IF( PRESENT(cd8 ) )   WRITE(numout,*) cd8
2863         IF( PRESENT(cd9 ) )   WRITE(numout,*) cd9
2864         IF( PRESENT(cd10) )   WRITE(numout,*) cd10
2865      ENDIF
2866                               CALL FLUSH(numout    )
2867      IF( numstp     /= -1 )   CALL FLUSH(numstp    )
2868      IF( numsol     /= -1 )   CALL FLUSH(numsol    )
2869      IF( numevo_ice /= -1 )   CALL FLUSH(numevo_ice)
2870      !
2871      IF( cd1 == 'STOP' ) THEN
2872         IF(lwp) WRITE(numout,*)  'huge E-R-R-O-R : immediate stop'
2873         CALL mppstop()
2874      ENDIF
2875      !
2876   END SUBROUTINE ctl_stop
2877
2878
2879   SUBROUTINE ctl_warn( cd1, cd2, cd3, cd4, cd5,   &
2880      &                 cd6, cd7, cd8, cd9, cd10 )
2881      !!----------------------------------------------------------------------
2882      !!                  ***  ROUTINE  stop_warn  ***
2883      !!
2884      !! ** Purpose :   print in ocean.outpput file a error message and
2885      !!                increment the warning number (nwarn) by one.
2886      !!----------------------------------------------------------------------
2887      CHARACTER(len=*), INTENT(in), OPTIONAL ::  cd1, cd2, cd3, cd4, cd5
2888      CHARACTER(len=*), INTENT(in), OPTIONAL ::  cd6, cd7, cd8, cd9, cd10
2889      !!----------------------------------------------------------------------
2890      !
2891      nwarn = nwarn + 1 
2892      IF(lwp) THEN
2893         WRITE(numout,cform_war)
2894         IF( PRESENT(cd1 ) ) WRITE(numout,*) cd1
2895         IF( PRESENT(cd2 ) ) WRITE(numout,*) cd2
2896         IF( PRESENT(cd3 ) ) WRITE(numout,*) cd3
2897         IF( PRESENT(cd4 ) ) WRITE(numout,*) cd4
2898         IF( PRESENT(cd5 ) ) WRITE(numout,*) cd5
2899         IF( PRESENT(cd6 ) ) WRITE(numout,*) cd6
2900         IF( PRESENT(cd7 ) ) WRITE(numout,*) cd7
2901         IF( PRESENT(cd8 ) ) WRITE(numout,*) cd8
2902         IF( PRESENT(cd9 ) ) WRITE(numout,*) cd9
2903         IF( PRESENT(cd10) ) WRITE(numout,*) cd10
2904      ENDIF
2905      CALL FLUSH(numout)
2906      !
2907   END SUBROUTINE ctl_warn
2908
2909
2910   SUBROUTINE ctl_opn( knum, cdfile, cdstat, cdform, cdacce, klengh, kout, ldwp, karea )
2911      !!----------------------------------------------------------------------
2912      !!                  ***  ROUTINE ctl_opn  ***
2913      !!
2914      !! ** Purpose :   Open file and check if required file is available.
2915      !!
2916      !! ** Method  :   Fortan open
2917      !!----------------------------------------------------------------------
2918      INTEGER          , INTENT(  out) ::   knum      ! logical unit to open
2919      CHARACTER(len=*) , INTENT(in   ) ::   cdfile    ! file name to open
2920      CHARACTER(len=*) , INTENT(in   ) ::   cdstat    ! disposition specifier
2921      CHARACTER(len=*) , INTENT(in   ) ::   cdform    ! formatting specifier
2922      CHARACTER(len=*) , INTENT(in   ) ::   cdacce    ! access specifier
2923      INTEGER          , INTENT(in   ) ::   klengh    ! record length
2924      INTEGER          , INTENT(in   ) ::   kout      ! number of logical units for write
2925      LOGICAL          , INTENT(in   ) ::   ldwp      ! boolean term for print
2926      INTEGER, OPTIONAL, INTENT(in   ) ::   karea     ! proc number
2927      !!
2928      CHARACTER(len=80) ::   clfile
2929      INTEGER           ::   iost
2930      !!----------------------------------------------------------------------
2931
2932      ! adapt filename
2933      ! ----------------
2934      clfile = TRIM(cdfile)
2935      IF( PRESENT( karea ) ) THEN
2936         IF( karea > 1 )   WRITE(clfile, "(a,'_',i4.4)") TRIM(clfile), karea-1
2937      ENDIF
2938#if defined key_agrif
2939      IF( .NOT. Agrif_Root() )   clfile = TRIM(Agrif_CFixed())//'_'//TRIM(clfile)
2940      knum=Agrif_Get_Unit()
2941#else
2942      knum=get_unit()
2943#endif
2944
2945      iost=0
2946      IF( cdacce(1:6) == 'DIRECT' )  THEN
2947         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat, RECL=klengh, ERR=100, IOSTAT=iost )
2948      ELSE
2949         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat             , ERR=100, IOSTAT=iost )
2950      ENDIF
2951      IF( iost == 0 ) THEN
2952         IF(ldwp) THEN
2953            WRITE(kout,*) '     file   : ', clfile,' open ok'
2954            WRITE(kout,*) '     unit   = ', knum
2955            WRITE(kout,*) '     status = ', cdstat
2956            WRITE(kout,*) '     form   = ', cdform
2957            WRITE(kout,*) '     access = ', cdacce
2958            WRITE(kout,*)
2959         ENDIF
2960      ENDIF
2961100   CONTINUE
2962      IF( iost /= 0 ) THEN
2963         IF(ldwp) THEN
2964            WRITE(kout,*)
2965            WRITE(kout,*) ' ===>>>> : bad opening file: ', clfile
2966            WRITE(kout,*) ' =======   ===  '
2967            WRITE(kout,*) '           unit   = ', knum
2968            WRITE(kout,*) '           status = ', cdstat
2969            WRITE(kout,*) '           form   = ', cdform
2970            WRITE(kout,*) '           access = ', cdacce
2971            WRITE(kout,*) '           iostat = ', iost
2972            WRITE(kout,*) '           we stop. verify the file '
2973            WRITE(kout,*)
2974         ENDIF
2975         STOP 'ctl_opn bad opening'
2976      ENDIF
2977     
2978   END SUBROUTINE ctl_opn
2979
2980
2981   INTEGER FUNCTION get_unit()
2982      !!----------------------------------------------------------------------
2983      !!                  ***  FUNCTION  get_unit  ***
2984      !!
2985      !! ** Purpose :   return the index of an unused logical unit
2986      !!----------------------------------------------------------------------
2987      LOGICAL :: llopn 
2988      !!----------------------------------------------------------------------
2989      !
2990      get_unit = 15   ! choose a unit that is big enough then it is not already used in NEMO
2991      llopn = .TRUE.
2992      DO WHILE( (get_unit < 998) .AND. llopn )
2993         get_unit = get_unit + 1
2994         INQUIRE( unit = get_unit, opened = llopn )
2995      END DO
2996      IF( (get_unit == 999) .AND. llopn ) THEN
2997         CALL ctl_stop( 'get_unit: All logical units until 999 are used...' )
2998         get_unit = -1
2999      ENDIF
3000      !
3001   END FUNCTION get_unit
3002
3003   !!----------------------------------------------------------------------
3004END MODULE lib_mpp
Note: See TracBrowser for help on using the repository browser.