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/UKMO/dev_r8183_ICEMODEL_svn_removed/NEMOGCM/NEMO/OPA_SRC/LBC – NEMO

source: branches/UKMO/dev_r8183_ICEMODEL_svn_removed/NEMOGCM/NEMO/OPA_SRC/LBC/lib_mpp.F90 @ 8738

Last change on this file since 8738 was 8738, checked in by dancopsey, 7 years ago

Merged in main ICEMODEL branch (branches/2017/dev_r8183_ICEMODEL) up to revision 8588

File size: 98.3 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   !!            3.5  !  2012  (S.Mocavero, I. Epicoco) Add mpp_lnk_bdy_3d/2d routines to optimize the BDY comm.
22   !!            3.5  !  2013  (C. Ethe, G. Madec)  message passing arrays as local variables
23   !!            3.5  !  2013  (S.Mocavero, I.Epicoco - CMCC) north fold optimizations
24   !!            3.6  !  2015  (O. Tintó and M. Castrillo - BSC) Added '_multiple' case for 2D lbc and max
25   !!            4.0  !  2017  (G. Madec) automatique allocation of array argument (use any 3rd dimension)
26   !!             -   !  2017  (G. Madec) create generic.h90 files to generate all lbc and north fold routines
27   !!----------------------------------------------------------------------
28
29   !!----------------------------------------------------------------------
30   !!   ctl_stop      : update momentum and tracer Kz from a tke scheme
31   !!   ctl_warn      : initialization, namelist read, and parameters control
32   !!   ctl_opn       : Open file and check if required file is available.
33   !!   ctl_nam       : Prints informations when an error occurs while reading a namelist
34   !!   get_unit      : give the index of an unused logical unit
35   !!----------------------------------------------------------------------
36#if   defined key_mpp_mpi
37   !!----------------------------------------------------------------------
38   !!   'key_mpp_mpi'             MPI massively parallel processing library
39   !!----------------------------------------------------------------------
40   !!   lib_mpp_alloc : allocate mpp arrays
41   !!   mynode        : indentify the processor unit
42   !!   mpp_lnk       : interface (defined in lbclnk) for message passing of 2d or 3d arrays (mpp_lnk_2d, mpp_lnk_3d)
43   !!   mpp_lnk_e     : interface (defined in lbclnk) for message passing of 2d array with extra halo (mpp_lnk_2d_e)
44   !!   mpp_lnk_icb   : interface for message passing of 2d arrays with extra halo for icebergs (mpp_lnk_2d_icb)
45   !!   mpprecv       :
46   !!   mppsend       :
47   !!   mppscatter    :
48   !!   mppgather     :
49   !!   mpp_min       : generic interface for mppmin_int , mppmin_a_int , mppmin_real, mppmin_a_real
50   !!   mpp_max       : generic interface for mppmax_int , mppmax_a_int , mppmax_real, mppmax_a_real
51   !!   mpp_sum       : generic interface for mppsum_int , mppsum_a_int , mppsum_real, mppsum_a_real
52   !!   mpp_minloc    :
53   !!   mpp_maxloc    :
54   !!   mppsync       :
55   !!   mppstop       :
56   !!   mpp_ini_north : initialisation of north fold
57!!gm   !!   mpp_lbc_north : north fold processors gathering
58   !!   mpp_lbc_north_e : variant of mpp_lbc_north for extra outer halo
59   !!   mpp_lbc_north_icb : variant of mpp_lbc_north for extra outer halo with icebergs
60   !!----------------------------------------------------------------------
61   USE dom_oce        ! ocean space and time domain
62   USE lbcnfd         ! north fold treatment
63   USE in_out_manager ! I/O manager
64   USE wrk_nemo       ! work arrays
65
66   IMPLICIT NONE
67   PRIVATE
68
69   INTERFACE mpp_nfd
70      MODULE PROCEDURE   mpp_nfd_2d      , mpp_nfd_3d      , mpp_nfd_4d
71      MODULE PROCEDURE   mpp_nfd_2d_ptr, mpp_nfd_3d_ptr, mpp_nfd_4d_ptr
72   END INTERFACE
73
74   ! Interface associated to the mpp_lnk_... routines is defined in lbclnk
75   PUBLIC   mpp_lnk_2d      , mpp_lnk_3d      , mpp_lnk_4d
76   PUBLIC   mpp_lnk_2d_ptr, mpp_lnk_3d_ptr, mpp_lnk_4d_ptr
77   PUBLIC   mpp_lnk_2d_e
78   !
79!!gm  this should be useless
80   PUBLIC   mpp_nfd_2d    , mpp_nfd_3d    , mpp_nfd_4d
81   PUBLIC   mpp_nfd_2d_ptr, mpp_nfd_3d_ptr, mpp_nfd_4d_ptr
82!!gm end
83   !
84   PUBLIC   ctl_stop, ctl_warn, get_unit, ctl_opn, ctl_nam
85   PUBLIC   mynode, mppstop, mppsync, mpp_comm_free
86   PUBLIC   mpp_ini_north, mpp_lbc_north_e
87!!gm   PUBLIC   mpp_ini_north, mpp_lbc_north, mpp_lbc_north_e
88   PUBLIC   mpp_lbc_north_icb, mpp_lnk_2d_icb
89   PUBLIC   mpp_min, mpp_max, mpp_sum, mpp_minloc, mpp_maxloc
90   PUBLIC   mpp_max_multiple
91!!gm   PUBLIC   mpp_lnk_2d_9
92!!gm   PUBLIC   mpp_lnk_sum_3d, mpp_lnk_sum_2d
93   PUBLIC   mppscatter, mppgather
94   PUBLIC   mpp_ini_ice, mpp_ini_znl
95   PUBLIC   mppsize
96   PUBLIC   mppsend, mpprecv                          ! needed by TAM and ICB routines
97   PUBLIC   mpp_lnk_bdy_2d, mpp_lnk_bdy_3d
98   PUBLIC   mpprank
99   
100   !! * Interfaces
101   !! define generic interface for these routine as they are called sometimes
102   !! with scalar arguments instead of array arguments, which causes problems
103   !! for the compilation on AIX system as well as NEC and SGI. Ok on COMPACQ
104   INTERFACE mpp_min
105      MODULE PROCEDURE mppmin_a_int, mppmin_int, mppmin_a_real, mppmin_real
106   END INTERFACE
107   INTERFACE mpp_max
108      MODULE PROCEDURE mppmax_a_int, mppmax_int, mppmax_a_real, mppmax_real
109   END INTERFACE
110   INTERFACE mpp_sum
111      MODULE PROCEDURE mppsum_a_int, mppsum_int, mppsum_a_real, mppsum_real,   &
112         &             mppsum_realdd, mppsum_a_realdd
113   END INTERFACE
114!!gm   INTERFACE mpp_lbc_north
115!!gm      MODULE PROCEDURE mpp_lbc_north_3d, mpp_lbc_north_2d
116!!gm   END INTERFACE
117   INTERFACE mpp_minloc
118      MODULE PROCEDURE mpp_minloc2d ,mpp_minloc3d
119   END INTERFACE
120   INTERFACE mpp_maxloc
121      MODULE PROCEDURE mpp_maxloc2d ,mpp_maxloc3d
122   END INTERFACE
123   INTERFACE mpp_max_multiple
124      MODULE PROCEDURE mppmax_real_multiple
125   END INTERFACE
126
127   !! ========================= !!
128   !!  MPI  variable definition !!
129   !! ========================= !!
130!$AGRIF_DO_NOT_TREAT
131   INCLUDE 'mpif.h'
132!$AGRIF_END_DO_NOT_TREAT
133
134   LOGICAL, PUBLIC, PARAMETER ::   lk_mpp = .TRUE.    !: mpp flag
135
136   INTEGER, PARAMETER         ::   nprocmax = 2**10   ! maximun dimension (required to be a power of 2)
137
138   INTEGER ::   mppsize        ! number of process
139   INTEGER ::   mpprank        ! process number  [ 0 - size-1 ]
140!$AGRIF_DO_NOT_TREAT
141   INTEGER, PUBLIC ::   mpi_comm_opa   ! opa local communicator
142!$AGRIF_END_DO_NOT_TREAT
143
144   INTEGER :: MPI_SUMDD
145
146   ! variables used in case of sea-ice
147   INTEGER, PUBLIC ::   ncomm_ice       !: communicator made by the processors with sea-ice (public so that it can be freed in icethd)
148   INTEGER         ::   ngrp_iworld     !  group ID for the world processors (for rheology)
149   INTEGER         ::   ngrp_ice        !  group ID for the ice processors (for rheology)
150   INTEGER         ::   ndim_rank_ice   !  number of 'ice' processors
151   INTEGER         ::   n_ice_root      !  number (in the comm_ice) of proc 0 in the ice comm
152   INTEGER, DIMENSION(:), ALLOCATABLE, SAVE ::   nrank_ice     ! dimension ndim_rank_ice
153
154   ! variables used for zonal integration
155   INTEGER, PUBLIC ::   ncomm_znl       !: communicator made by the processors on the same zonal average
156   LOGICAL, PUBLIC ::   l_znl_root      !: True on the 'left'most processor on the same row
157   INTEGER         ::   ngrp_znl        !  group ID for the znl processors
158   INTEGER         ::   ndim_rank_znl   !  number of processors on the same zonal average
159   INTEGER, DIMENSION(:), ALLOCATABLE, SAVE ::   nrank_znl  ! dimension ndim_rank_znl, number of the procs into the same znl domain
160
161   ! North fold condition in mpp_mpi with jpni > 1 (PUBLIC for TAM)
162   INTEGER, PUBLIC ::   ngrp_world        !: group ID for the world processors
163   INTEGER, PUBLIC ::   ngrp_opa          !: group ID for the opa processors
164   INTEGER, PUBLIC ::   ngrp_north        !: group ID for the northern processors (to be fold)
165   INTEGER, PUBLIC ::   ncomm_north       !: communicator made by the processors belonging to ngrp_north
166   INTEGER, PUBLIC ::   ndim_rank_north   !: number of 'sea' processor in the northern line (can be /= jpni !)
167   INTEGER, PUBLIC ::   njmppmax          !: value of njmpp for the processors of the northern line
168   INTEGER, PUBLIC ::   north_root        !: number (in the comm_opa) of proc 0 in the northern comm
169   INTEGER, PUBLIC, DIMENSION(:), ALLOCATABLE, SAVE ::   nrank_north   !: dimension ndim_rank_north
170
171   ! Type of send : standard, buffered, immediate
172   CHARACTER(len=1), PUBLIC ::   cn_mpi_send        !: type od mpi send/recieve (S=standard, B=bsend, I=isend)
173   LOGICAL         , PUBLIC ::   l_isend = .FALSE.  !: isend use indicator (T if cn_mpi_send='I')
174   INTEGER         , PUBLIC ::   nn_buffer          !: size of the buffer in case of mpi_bsend
175
176   REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE ::   tampon   ! buffer in case of bsend
177
178   LOGICAL, PUBLIC ::   ln_nnogather                !: namelist control of northfold comms
179   LOGICAL, PUBLIC ::   l_north_nogather = .FALSE.  !: internal control of northfold comms
180
181   !!----------------------------------------------------------------------
182   !! NEMO/OPA 4.0 , NEMO Consortium (2017)
183   !! $Id$
184   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
185   !!----------------------------------------------------------------------
186CONTAINS
187
188   FUNCTION mynode( ldtxt, ldname, kumnam_ref, kumnam_cfg, kumond, kstop, localComm )
189      !!----------------------------------------------------------------------
190      !!                  ***  routine mynode  ***
191      !!
192      !! ** Purpose :   Find processor unit
193      !!----------------------------------------------------------------------
194      CHARACTER(len=*),DIMENSION(:), INTENT(  out) ::   ldtxt        !
195      CHARACTER(len=*)             , INTENT(in   ) ::   ldname       !
196      INTEGER                      , INTENT(in   ) ::   kumnam_ref   ! logical unit for reference namelist
197      INTEGER                      , INTENT(in   ) ::   kumnam_cfg   ! logical unit for configuration namelist
198      INTEGER                      , INTENT(inout) ::   kumond       ! logical unit for namelist output
199      INTEGER                      , INTENT(inout) ::   kstop        ! stop indicator
200      INTEGER         , OPTIONAL   , INTENT(in   ) ::   localComm    !
201      !
202      INTEGER ::   mynode, ierr, code, ji, ii, ios
203      LOGICAL ::   mpi_was_called
204      !
205      NAMELIST/nammpp/ cn_mpi_send, nn_buffer, jpni, jpnj, jpnij, ln_nnogather
206      !!----------------------------------------------------------------------
207      !
208      ii = 1
209      WRITE(ldtxt(ii),*)                                                                  ;   ii = ii + 1
210      WRITE(ldtxt(ii),*) 'mynode : mpi initialisation'                                    ;   ii = ii + 1
211      WRITE(ldtxt(ii),*) '~~~~~~ '                                                        ;   ii = ii + 1
212      !
213      REWIND( kumnam_ref )              ! Namelist nammpp in reference namelist: mpi variables
214      READ  ( kumnam_ref, nammpp, IOSTAT = ios, ERR = 901)
215901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nammpp in reference namelist', lwp )
216      !
217      REWIND( kumnam_cfg )              ! Namelist nammpp in configuration namelist: mpi variables
218      READ  ( kumnam_cfg, nammpp, IOSTAT = ios, ERR = 902 )
219902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nammpp in configuration namelist', lwp )
220      !
221      !                              ! control print
222      WRITE(ldtxt(ii),*) '   Namelist nammpp'                                             ;   ii = ii + 1
223      WRITE(ldtxt(ii),*) '      mpi send type          cn_mpi_send = ', cn_mpi_send       ;   ii = ii + 1
224      WRITE(ldtxt(ii),*) '      size exported buffer   nn_buffer   = ', nn_buffer,' bytes';   ii = ii + 1
225      !
226#if defined key_agrif
227      IF( .NOT. Agrif_Root() ) THEN
228         jpni  = Agrif_Parent(jpni )
229         jpnj  = Agrif_Parent(jpnj )
230         jpnij = Agrif_Parent(jpnij)
231      ENDIF
232#endif
233      !
234      IF( jpnij < 1 ) THEN         ! If jpnij is not specified in namelist then we calculate it
235         jpnij = jpni * jpnj       ! this means there will be no land cutting out.
236      ENDIF
237
238      IF( jpni < 1 .OR. jpnj < 1  ) THEN
239         WRITE(ldtxt(ii),*) '      jpni, jpnj and jpnij will be calculated automatically' ;   ii = ii + 1
240      ELSE
241         WRITE(ldtxt(ii),*) '      processor grid extent in i         jpni = ',jpni       ;   ii = ii + 1
242         WRITE(ldtxt(ii),*) '      processor grid extent in j         jpnj = ',jpnj       ;   ii = ii + 1
243         WRITE(ldtxt(ii),*) '      number of local domains           jpnij = ',jpnij      ;   ii = ii + 1
244      ENDIF
245
246      WRITE(ldtxt(ii),*) '      avoid use of mpi_allgather at the north fold  ln_nnogather = ', ln_nnogather  ; ii = ii + 1
247
248      CALL mpi_initialized ( mpi_was_called, code )
249      IF( code /= MPI_SUCCESS ) THEN
250         DO ji = 1, SIZE(ldtxt)
251            IF( TRIM(ldtxt(ji)) /= '' )   WRITE(*,*) ldtxt(ji)      ! control print of mynode
252         END DO
253         WRITE(*, cform_err)
254         WRITE(*, *) 'lib_mpp: Error in routine mpi_initialized'
255         CALL mpi_abort( mpi_comm_world, code, ierr )
256      ENDIF
257
258      IF( mpi_was_called ) THEN
259         !
260         SELECT CASE ( cn_mpi_send )
261         CASE ( 'S' )                ! Standard mpi send (blocking)
262            WRITE(ldtxt(ii),*) '           Standard blocking mpi send (send)'             ;   ii = ii + 1
263         CASE ( 'B' )                ! Buffer mpi send (blocking)
264            WRITE(ldtxt(ii),*) '           Buffer blocking mpi send (bsend)'              ;   ii = ii + 1
265            IF( Agrif_Root() )   CALL mpi_init_opa( ldtxt, ii, ierr )
266         CASE ( 'I' )                ! Immediate mpi send (non-blocking send)
267            WRITE(ldtxt(ii),*) '           Immediate non-blocking send (isend)'           ;   ii = ii + 1
268            l_isend = .TRUE.
269         CASE DEFAULT
270            WRITE(ldtxt(ii),cform_err)                                                    ;   ii = ii + 1
271            WRITE(ldtxt(ii),*) '           bad value for cn_mpi_send = ', cn_mpi_send     ;   ii = ii + 1
272            kstop = kstop + 1
273         END SELECT
274         !
275      ELSEIF ( PRESENT(localComm) .AND. .NOT. mpi_was_called ) THEN
276         WRITE(ldtxt(ii),*) ' lib_mpp: You cannot provide a local communicator '          ;   ii = ii + 1
277         WRITE(ldtxt(ii),*) '          without calling MPI_Init before ! '                ;   ii = ii + 1
278         kstop = kstop + 1
279      ELSE
280         SELECT CASE ( cn_mpi_send )
281         CASE ( 'S' )                ! Standard mpi send (blocking)
282            WRITE(ldtxt(ii),*) '           Standard blocking mpi send (send)'             ;   ii = ii + 1
283            CALL mpi_init( ierr )
284         CASE ( 'B' )                ! Buffer mpi send (blocking)
285            WRITE(ldtxt(ii),*) '           Buffer blocking mpi send (bsend)'              ;   ii = ii + 1
286            IF( Agrif_Root() )   CALL mpi_init_opa( ldtxt, ii, ierr )
287         CASE ( 'I' )                ! Immediate mpi send (non-blocking send)
288            WRITE(ldtxt(ii),*) '           Immediate non-blocking send (isend)'           ;   ii = ii + 1
289            l_isend = .TRUE.
290            CALL mpi_init( ierr )
291         CASE DEFAULT
292            WRITE(ldtxt(ii),cform_err)                                                    ;   ii = ii + 1
293            WRITE(ldtxt(ii),*) '           bad value for cn_mpi_send = ', cn_mpi_send     ;   ii = ii + 1
294            kstop = kstop + 1
295         END SELECT
296         !
297      ENDIF
298
299      IF( PRESENT(localComm) ) THEN
300         IF( Agrif_Root() ) THEN
301            mpi_comm_opa = localComm
302         ENDIF
303      ELSE
304         CALL mpi_comm_dup( mpi_comm_world, mpi_comm_opa, code)
305         IF( code /= MPI_SUCCESS ) THEN
306            DO ji = 1, SIZE(ldtxt)
307               IF( TRIM(ldtxt(ji)) /= '' )   WRITE(*,*) ldtxt(ji)      ! control print of mynode
308            END DO
309            WRITE(*, cform_err)
310            WRITE(*, *) ' lib_mpp: Error in routine mpi_comm_dup'
311            CALL mpi_abort( mpi_comm_world, code, ierr )
312         ENDIF
313      ENDIF
314
315#if defined key_agrif
316      IF( Agrif_Root() ) THEN
317         CALL Agrif_MPI_Init(mpi_comm_opa)
318      ELSE
319         CALL Agrif_MPI_set_grid_comm(mpi_comm_opa)
320      ENDIF
321#endif
322
323      CALL mpi_comm_rank( mpi_comm_opa, mpprank, ierr )
324      CALL mpi_comm_size( mpi_comm_opa, mppsize, ierr )
325      mynode = mpprank
326
327      IF( mynode == 0 ) THEN
328         CALL ctl_opn( kumond, TRIM(ldname), 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. , 1 )
329         WRITE(kumond, nammpp)     
330      ENDIF
331      !
332      CALL MPI_OP_CREATE(DDPDD_MPI, .TRUE., MPI_SUMDD, ierr)
333      !
334   END FUNCTION mynode
335
336   !!----------------------------------------------------------------------
337   !!                   ***  routine mpp_lnk_(2,3,4)d  ***
338   !!
339   !!   * Argument : dummy argument use in mpp_lnk_... routines
340   !!                ptab   :   array or pointer of arrays on which the boundary condition is applied
341   !!                cd_nat :   nature of array grid-points
342   !!                psgn   :   sign used across the north fold boundary
343   !!                kfld   :   optional, number of pt3d arrays
344   !!                cd_mpp :   optional, fill the overlap area only
345   !!                pval   :   optional, background value (used at closed boundaries)
346   !!----------------------------------------------------------------------
347   !
348   !                       !==  2D array and array of 2D pointer  ==!
349   !
350#  define DIM_2d
351#     define ROUTINE_LNK           mpp_lnk_2d
352#     include "mpp_lnk_generic.h90"
353#     undef ROUTINE_LNK
354#     define MULTI
355#     define ROUTINE_LNK           mpp_lnk_2d_ptr
356#     include "mpp_lnk_generic.h90"
357#     undef ROUTINE_LNK
358#     undef MULTI
359#  undef DIM_2d
360   !
361   !                       !==  3D array and array of 3D pointer  ==!
362   !
363#  define DIM_3d
364#     define ROUTINE_LNK           mpp_lnk_3d
365#     include "mpp_lnk_generic.h90"
366#     undef ROUTINE_LNK
367#     define MULTI
368#     define ROUTINE_LNK           mpp_lnk_3d_ptr
369#     include "mpp_lnk_generic.h90"
370#     undef ROUTINE_LNK
371#     undef MULTI
372#  undef DIM_3d
373   !
374   !                       !==  4D array and array of 4D pointer  ==!
375   !
376#  define DIM_4d
377#     define ROUTINE_LNK           mpp_lnk_4d
378#     include "mpp_lnk_generic.h90"
379#     undef ROUTINE_LNK
380#     define MULTI
381#     define ROUTINE_LNK           mpp_lnk_4d_ptr
382#     include "mpp_lnk_generic.h90"
383#     undef ROUTINE_LNK
384#     undef MULTI
385#  undef DIM_4d
386
387   !!----------------------------------------------------------------------
388   !!                   ***  routine mpp_nfd_(2,3,4)d  ***
389   !!
390   !!   * Argument : dummy argument use in mpp_nfd_... routines
391   !!                ptab   :   array or pointer of arrays on which the boundary condition is applied
392   !!                cd_nat :   nature of array grid-points
393   !!                psgn   :   sign used across the north fold boundary
394   !!                kfld   :   optional, number of pt3d arrays
395   !!                cd_mpp :   optional, fill the overlap area only
396   !!                pval   :   optional, background value (used at closed boundaries)
397   !!----------------------------------------------------------------------
398   !
399   !                       !==  2D array and array of 2D pointer  ==!
400   !
401#  define DIM_2d
402#     define ROUTINE_NFD           mpp_nfd_2d
403#     include "mpp_nfd_generic.h90"
404#     undef ROUTINE_NFD
405#     define MULTI
406#     define ROUTINE_NFD           mpp_nfd_2d_ptr
407#     include "mpp_nfd_generic.h90"
408#     undef ROUTINE_NFD
409#     undef MULTI
410#  undef DIM_2d
411   !
412   !                       !==  3D array and array of 3D pointer  ==!
413   !
414#  define DIM_3d
415#     define ROUTINE_NFD           mpp_nfd_3d
416#     include "mpp_nfd_generic.h90"
417#     undef ROUTINE_NFD
418#     define MULTI
419#     define ROUTINE_NFD           mpp_nfd_3d_ptr
420#     include "mpp_nfd_generic.h90"
421#     undef ROUTINE_NFD
422#     undef MULTI
423#  undef DIM_3d
424   !
425   !                       !==  4D array and array of 4D pointer  ==!
426   !
427#  define DIM_4d
428#     define ROUTINE_NFD           mpp_nfd_4d
429#     include "mpp_nfd_generic.h90"
430#     undef ROUTINE_NFD
431#     define MULTI
432#     define ROUTINE_NFD           mpp_nfd_4d_ptr
433#     include "mpp_nfd_generic.h90"
434#     undef ROUTINE_NFD
435#     undef MULTI
436#  undef DIM_4d
437
438
439   !!----------------------------------------------------------------------
440   !!                   ***  routine mpp_lnk_bdy_(2,3,4)d  ***
441   !!
442   !!   * Argument : dummy argument use in mpp_lnk_... routines
443   !!                ptab   :   array or pointer of arrays on which the boundary condition is applied
444   !!                cd_nat :   nature of array grid-points
445   !!                psgn   :   sign used across the north fold boundary
446   !!                kb_bdy :   BDY boundary set
447   !!                kfld   :   optional, number of pt3d arrays
448   !!----------------------------------------------------------------------
449   !
450   !                       !==  2D array and array of 2D pointer  ==!
451   !
452#  define DIM_2d
453#     define ROUTINE_BDY           mpp_lnk_bdy_2d
454#     include "mpp_bdy_generic.h90"
455#     undef ROUTINE_BDY
456#     define MULTI
457#     define ROUTINE_BDY           mpp_lnk_bdy_2d_ptr
458#     include "mpp_bdy_generic.h90"
459#     undef ROUTINE_BDY
460#     undef MULTI
461#  undef DIM_2d
462   !
463   !                       !==  3D array and array of 3D pointer  ==!
464   !
465#  define DIM_3d
466#     define ROUTINE_BDY           mpp_lnk_bdy_3d
467#     include "mpp_bdy_generic.h90"
468#     undef ROUTINE_BDY
469#     define MULTI
470#     define ROUTINE_BDY           mpp_lnk_bdy_3d_ptr
471#     include "mpp_bdy_generic.h90"
472#     undef ROUTINE_BDY
473#     undef MULTI
474#  undef DIM_3d
475   !
476   !                       !==  4D array and array of 4D pointer  ==!
477   !
478!!#  define DIM_4d
479!!#     define ROUTINE_BDY           mpp_lnk_bdy_4d
480!!#     include "mpp_bdy_generic.h90"
481!!#     undef ROUTINE_BDY
482!!#     define MULTI
483!!#     define ROUTINE_BDY           mpp_lnk_bdy_4d_ptr
484!!#     include "mpp_bdy_generic.h90"
485!!#     undef ROUTINE_BDY
486!!#     undef MULTI
487!!#  undef DIM_4d
488
489   !!----------------------------------------------------------------------
490   !!
491   !!   load_array  &   mpp_lnk_2d_9    à generaliser a 3D et 4D
492   
493   
494   !!    mpp_lnk_2d_e     utilisé dans ICB
495
496
497   !!    mpp_lnk_sum_2d et 3D   ====>>>>>>   à virer du code !!!!
498   
499   
500   !!----------------------------------------------------------------------
501
502
503   SUBROUTINE mpp_lnk_2d_e( pt2d, cd_type, psgn, jpri, jprj )
504      !!----------------------------------------------------------------------
505      !!                  ***  routine mpp_lnk_2d_e  ***
506      !!
507      !! ** Purpose :   Message passing manadgement for 2d array (with halo)
508      !!
509      !! ** Method  :   Use mppsend and mpprecv function for passing mask
510      !!      between processors following neighboring subdomains.
511      !!            domain parameters
512      !!                    nlci   : first dimension of the local subdomain
513      !!                    nlcj   : second dimension of the local subdomain
514      !!                    jpri   : number of rows for extra outer halo
515      !!                    jprj   : number of columns for extra outer halo
516      !!                    nbondi : mark for "east-west local boundary"
517      !!                    nbondj : mark for "north-south local boundary"
518      !!                    noea   : number for local neighboring processors
519      !!                    nowe   : number for local neighboring processors
520      !!                    noso   : number for local neighboring processors
521      !!                    nono   : number for local neighboring processors
522      !!
523      !!----------------------------------------------------------------------
524      INTEGER                                             , INTENT(in   ) ::   jpri
525      INTEGER                                             , INTENT(in   ) ::   jprj
526      REAL(wp), DIMENSION(1-jpri:jpi+jpri,1-jprj:jpj+jprj), INTENT(inout) ::   pt2d     ! 2D array with extra halo
527      CHARACTER(len=1)                                    , INTENT(in   ) ::   cd_type  ! nature of ptab array grid-points
528      !                                                                                 ! = T , U , V , F , W and I points
529      REAL(wp)                                            , INTENT(in   ) ::   psgn     ! =-1 the sign change across the
530      !!                                                                                ! north boundary, =  1. otherwise
531      INTEGER  ::   jl   ! dummy loop indices
532      INTEGER  ::   imigr, iihom, ijhom        ! temporary integers
533      INTEGER  ::   ipreci, iprecj             ! temporary integers
534      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
535      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
536      !!
537      REAL(wp), DIMENSION(1-jpri:jpi+jpri,jprecj+jprj,2) :: r2dns
538      REAL(wp), DIMENSION(1-jpri:jpi+jpri,jprecj+jprj,2) :: r2dsn
539      REAL(wp), DIMENSION(1-jprj:jpj+jprj,jpreci+jpri,2) :: r2dwe
540      REAL(wp), DIMENSION(1-jprj:jpj+jprj,jpreci+jpri,2) :: r2dew
541      !!----------------------------------------------------------------------
542
543      ipreci = jpreci + jpri      ! take into account outer extra 2D overlap area
544      iprecj = jprecj + jprj
545
546
547      ! 1. standard boundary treatment   (CAUTION: the order matters Here !!!! )
548      ! ------------------------------
549      !                                !== North-South boundaries
550      !                                      !* cyclic
551      IF( nbondj == 2 .AND. jperio == 7 ) THEN
552         pt2d(:, 1-jprj:  1     ) = pt2d ( :, jpjm1-jprj:jpjm1 )
553         pt2d(:, jpj   :jpj+jprj) = pt2d ( :, 2         :2+jprj)
554      ELSE                                   !* closed
555         IF( .NOT. cd_type == 'F' )   pt2d(:,  1-jprj   :  jprecj  ) = 0._wp     ! south except at F-point
556                                      pt2d(:,nlcj-jprecj+1:jpj+jprj) = 0._wp     ! north
557      ENDIF
558      !                                !== East-West boundaries
559      !                                      !* Cyclic east-west
560      IF( nbondi == 2 .AND. (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
561         pt2d(1-jpri:     1    ,:) = pt2d(jpim1-jpri:  jpim1 ,:)              ! east
562         pt2d(   jpi  :jpi+jpri,:) = pt2d(     2      :2+jpri,:)              ! west
563      ELSE                                   !* closed
564         IF( .NOT. cd_type == 'F' )   pt2d(  1-jpri   :jpreci    ,:) = 0._wp  ! south except at F-point
565                                      pt2d(nlci-jpreci+1:jpi+jpri,:) = 0._wp  ! north
566      ENDIF
567      !
568      ! north fold treatment
569      ! --------------------
570      IF( npolj /= 0 ) THEN
571         !
572         SELECT CASE ( jpni )
573!!gm ERROR        CASE ( 1 )     ;   CALL lbc_nfd        ( pt2d(1:jpi,1:jpj+jprj), cd_type, psgn, pr2dj=jprj )
574!!gm ERROR         CASE DEFAULT   ;   CALL mpp_lbc_north_e( pt2d                  , cd_type, psgn             )
575         END SELECT
576         !
577      ENDIF
578
579      ! 2. East and west directions exchange
580      ! ------------------------------------
581      ! we play with the neigbours AND the row number because of the periodicity
582      !
583      SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
584      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
585         iihom = nlci-nreci-jpri
586         DO jl = 1, ipreci
587            r2dew(:,jl,1) = pt2d(jpreci+jl,:)
588            r2dwe(:,jl,1) = pt2d(iihom +jl,:)
589         END DO
590      END SELECT
591      !
592      !                           ! Migrations
593      imigr = ipreci * ( jpj + 2*jprj)
594      !
595      SELECT CASE ( nbondi )
596      CASE ( -1 )
597         CALL mppsend( 2, r2dwe(1-jprj,1,1), imigr, noea, ml_req1 )
598         CALL mpprecv( 1, r2dew(1-jprj,1,2), imigr, noea )
599         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
600      CASE ( 0 )
601         CALL mppsend( 1, r2dew(1-jprj,1,1), imigr, nowe, ml_req1 )
602         CALL mppsend( 2, r2dwe(1-jprj,1,1), imigr, noea, ml_req2 )
603         CALL mpprecv( 1, r2dew(1-jprj,1,2), imigr, noea )
604         CALL mpprecv( 2, r2dwe(1-jprj,1,2), imigr, nowe )
605         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
606         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
607      CASE ( 1 )
608         CALL mppsend( 1, r2dew(1-jprj,1,1), imigr, nowe, ml_req1 )
609         CALL mpprecv( 2, r2dwe(1-jprj,1,2), imigr, nowe )
610         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
611      END SELECT
612      !
613      !                           ! Write Dirichlet lateral conditions
614      iihom = nlci - jpreci
615      !
616      SELECT CASE ( nbondi )
617      CASE ( -1 )
618         DO jl = 1, ipreci
619            pt2d(iihom+jl,:) = r2dew(:,jl,2)
620         END DO
621      CASE ( 0 )
622         DO jl = 1, ipreci
623            pt2d(jl-jpri,:) = r2dwe(:,jl,2)
624            pt2d( iihom+jl,:) = r2dew(:,jl,2)
625         END DO
626      CASE ( 1 )
627         DO jl = 1, ipreci
628            pt2d(jl-jpri,:) = r2dwe(:,jl,2)
629         END DO
630      END SELECT
631
632      ! 3. North and south directions
633      ! -----------------------------
634      ! always closed : we play only with the neigbours
635      !
636      IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
637         ijhom = nlcj-nrecj-jprj
638         DO jl = 1, iprecj
639            r2dsn(:,jl,1) = pt2d(:,ijhom +jl)
640            r2dns(:,jl,1) = pt2d(:,jprecj+jl)
641         END DO
642      ENDIF
643      !
644      !                           ! Migrations
645      imigr = iprecj * ( jpi + 2*jpri )
646      !
647      SELECT CASE ( nbondj )
648      CASE ( -1 )
649         CALL mppsend( 4, r2dsn(1-jpri,1,1), imigr, nono, ml_req1 )
650         CALL mpprecv( 3, r2dns(1-jpri,1,2), imigr, nono )
651         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
652      CASE ( 0 )
653         CALL mppsend( 3, r2dns(1-jpri,1,1), imigr, noso, ml_req1 )
654         CALL mppsend( 4, r2dsn(1-jpri,1,1), imigr, nono, ml_req2 )
655         CALL mpprecv( 3, r2dns(1-jpri,1,2), imigr, nono )
656         CALL mpprecv( 4, r2dsn(1-jpri,1,2), imigr, noso )
657         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
658         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
659      CASE ( 1 )
660         CALL mppsend( 3, r2dns(1-jpri,1,1), imigr, noso, ml_req1 )
661         CALL mpprecv( 4, r2dsn(1-jpri,1,2), imigr, noso )
662         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
663      END SELECT
664      !
665      !                           ! Write Dirichlet lateral conditions
666      ijhom = nlcj - jprecj
667      !
668      SELECT CASE ( nbondj )
669      CASE ( -1 )
670         DO jl = 1, iprecj
671            pt2d(:,ijhom+jl) = r2dns(:,jl,2)
672         END DO
673      CASE ( 0 )
674         DO jl = 1, iprecj
675            pt2d(:,jl-jprj) = r2dsn(:,jl,2)
676            pt2d(:,ijhom+jl ) = r2dns(:,jl,2)
677         END DO
678      CASE ( 1 )
679         DO jl = 1, iprecj
680            pt2d(:,jl-jprj) = r2dsn(:,jl,2)
681         END DO
682      END SELECT
683      !
684   END SUBROUTINE mpp_lnk_2d_e
685
686
687   SUBROUTINE mppsend( ktyp, pmess, kbytes, kdest, md_req )
688      !!----------------------------------------------------------------------
689      !!                  ***  routine mppsend  ***
690      !!
691      !! ** Purpose :   Send messag passing array
692      !!
693      !!----------------------------------------------------------------------
694      REAL(wp), INTENT(inout) ::   pmess(*)   ! array of real
695      INTEGER , INTENT(in   ) ::   kbytes     ! size of the array pmess
696      INTEGER , INTENT(in   ) ::   kdest      ! receive process number
697      INTEGER , INTENT(in   ) ::   ktyp       ! tag of the message
698      INTEGER , INTENT(in   ) ::   md_req     ! argument for isend
699      !!
700      INTEGER ::   iflag
701      !!----------------------------------------------------------------------
702      !
703      SELECT CASE ( cn_mpi_send )
704      CASE ( 'S' )                ! Standard mpi send (blocking)
705         CALL mpi_send ( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_opa        , iflag )
706      CASE ( 'B' )                ! Buffer mpi send (blocking)
707         CALL mpi_bsend( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_opa        , iflag )
708      CASE ( 'I' )                ! Immediate mpi send (non-blocking send)
709         ! be carefull, one more argument here : the mpi request identifier..
710         CALL mpi_isend( pmess, kbytes, mpi_double_precision, kdest , ktyp, mpi_comm_opa, md_req, iflag )
711      END SELECT
712      !
713   END SUBROUTINE mppsend
714
715
716   SUBROUTINE mpprecv( ktyp, pmess, kbytes, ksource )
717      !!----------------------------------------------------------------------
718      !!                  ***  routine mpprecv  ***
719      !!
720      !! ** Purpose :   Receive messag passing array
721      !!
722      !!----------------------------------------------------------------------
723      REAL(wp), INTENT(inout) ::   pmess(*)   ! array of real
724      INTEGER , INTENT(in   ) ::   kbytes     ! suze of the array pmess
725      INTEGER , INTENT(in   ) ::   ktyp       ! Tag of the recevied message
726      INTEGER, OPTIONAL, INTENT(in) :: ksource    ! source process number
727      !!
728      INTEGER :: istatus(mpi_status_size)
729      INTEGER :: iflag
730      INTEGER :: use_source
731      !!----------------------------------------------------------------------
732      !
733      ! If a specific process number has been passed to the receive call,
734      ! use that one. Default is to use mpi_any_source
735      use_source = mpi_any_source
736      IF( PRESENT(ksource) )   use_source = ksource
737      !
738      CALL mpi_recv( pmess, kbytes, mpi_double_precision, use_source, ktyp, mpi_comm_opa, istatus, iflag )
739      !
740   END SUBROUTINE mpprecv
741
742
743   SUBROUTINE mppgather( ptab, kp, pio )
744      !!----------------------------------------------------------------------
745      !!                   ***  routine mppgather  ***
746      !!
747      !! ** Purpose :   Transfert between a local subdomain array and a work
748      !!     array which is distributed following the vertical level.
749      !!
750      !!----------------------------------------------------------------------
751      REAL(wp), DIMENSION(jpi,jpj)      , INTENT(in   ) ::   ptab   ! subdomain input array
752      INTEGER                           , INTENT(in   ) ::   kp     ! record length
753      REAL(wp), DIMENSION(jpi,jpj,jpnij), INTENT(  out) ::   pio    ! subdomain input array
754      !!
755      INTEGER :: itaille, ierror   ! temporary integer
756      !!---------------------------------------------------------------------
757      !
758      itaille = jpi * jpj
759      CALL mpi_gather( ptab, itaille, mpi_double_precision, pio, itaille     ,   &
760         &                            mpi_double_precision, kp , mpi_comm_opa, ierror )
761      !
762   END SUBROUTINE mppgather
763
764
765   SUBROUTINE mppscatter( pio, kp, ptab )
766      !!----------------------------------------------------------------------
767      !!                  ***  routine mppscatter  ***
768      !!
769      !! ** Purpose :   Transfert between awork array which is distributed
770      !!      following the vertical level and the local subdomain array.
771      !!
772      !!----------------------------------------------------------------------
773      REAL(wp), DIMENSION(jpi,jpj,jpnij)  ::   pio    ! output array
774      INTEGER                             ::   kp     ! Tag (not used with MPI
775      REAL(wp), DIMENSION(jpi,jpj)        ::   ptab   ! subdomain array input
776      !!
777      INTEGER :: itaille, ierror   ! temporary integer
778      !!---------------------------------------------------------------------
779      !
780      itaille = jpi * jpj
781      !
782      CALL mpi_scatter( pio, itaille, mpi_double_precision, ptab, itaille     ,   &
783         &                            mpi_double_precision, kp  , mpi_comm_opa, ierror )
784      !
785   END SUBROUTINE mppscatter
786
787   !!----------------------------------------------------------------------
788   !!    ***  mppmax_a_int, mppmax_int, mppmax_a_real, mppmax_real  ***
789   !!   
790   !!----------------------------------------------------------------------
791   !!
792   SUBROUTINE mppmax_a_int( ktab, kdim, kcom )
793      !!----------------------------------------------------------------------
794      INTEGER , INTENT(in   )                  ::   kdim   ! size of array
795      INTEGER , INTENT(inout), DIMENSION(kdim) ::   ktab   ! input array
796      INTEGER , INTENT(in   ), OPTIONAL        ::   kcom   !
797      INTEGER :: ierror, ilocalcomm   ! temporary integer
798      INTEGER, DIMENSION(kdim) ::   iwork
799      !!----------------------------------------------------------------------
800      ilocalcomm = mpi_comm_opa
801      IF( PRESENT(kcom) )   ilocalcomm = kcom
802      CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_max, ilocalcomm, ierror )
803      ktab(:) = iwork(:)
804   END SUBROUTINE mppmax_a_int
805   !!
806   SUBROUTINE mppmax_int( ktab, kcom )
807      !!----------------------------------------------------------------------
808      INTEGER, INTENT(inout)           ::   ktab   ! ???
809      INTEGER, INTENT(in   ), OPTIONAL ::   kcom   ! ???
810      INTEGER ::   ierror, iwork, ilocalcomm   ! temporary integer
811      !!----------------------------------------------------------------------
812      ilocalcomm = mpi_comm_opa
813      IF( PRESENT(kcom) )   ilocalcomm = kcom
814      CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_max, ilocalcomm, ierror )
815      ktab = iwork
816   END SUBROUTINE mppmax_int
817   !!
818   SUBROUTINE mppmax_a_real( ptab, kdim, kcom )
819      !!----------------------------------------------------------------------
820      REAL(wp), DIMENSION(kdim), INTENT(inout) ::   ptab
821      INTEGER                  , INTENT(in   ) ::   kdim
822      INTEGER , OPTIONAL       , INTENT(in   ) ::   kcom
823      INTEGER :: ierror, ilocalcomm
824      REAL(wp), DIMENSION(kdim) ::  zwork
825      !!----------------------------------------------------------------------
826      ilocalcomm = mpi_comm_opa
827      IF( PRESENT(kcom) )   ilocalcomm = kcom
828      CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_max, ilocalcomm, ierror )
829      ptab(:) = zwork(:)
830   END SUBROUTINE mppmax_a_real
831   !!
832   SUBROUTINE mppmax_real( ptab, kcom )
833      !!----------------------------------------------------------------------
834      REAL(wp), INTENT(inout)           ::   ptab   ! ???
835      INTEGER , INTENT(in   ), OPTIONAL ::   kcom   ! ???
836      INTEGER  ::   ierror, ilocalcomm
837      REAL(wp) ::   zwork
838      !!----------------------------------------------------------------------
839      ilocalcomm = mpi_comm_opa
840      IF( PRESENT(kcom) )   ilocalcomm = kcom!
841      CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_max, ilocalcomm, ierror )
842      ptab = zwork
843   END SUBROUTINE mppmax_real
844
845
846   !!----------------------------------------------------------------------
847   !!    ***  mppmin_a_int, mppmin_int, mppmin_a_real, mppmin_real  ***
848   !!   
849   !!----------------------------------------------------------------------
850   !!
851   SUBROUTINE mppmin_a_int( ktab, kdim, kcom )
852      !!----------------------------------------------------------------------
853      INTEGER , INTENT( in  )                  ::   kdim   ! size of array
854      INTEGER , INTENT(inout), DIMENSION(kdim) ::   ktab   ! input array
855      INTEGER , INTENT( in  ), OPTIONAL        ::   kcom   ! input array
856      !!
857      INTEGER ::   ierror, ilocalcomm   ! temporary integer
858      INTEGER, DIMENSION(kdim) ::   iwork
859      !!----------------------------------------------------------------------
860      ilocalcomm = mpi_comm_opa
861      IF( PRESENT(kcom) )   ilocalcomm = kcom
862      CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_min, ilocalcomm, ierror )
863      ktab(:) = iwork(:)
864   END SUBROUTINE mppmin_a_int
865   !!
866   SUBROUTINE mppmin_int( ktab, kcom )
867      !!----------------------------------------------------------------------
868      INTEGER, INTENT(inout) ::   ktab      ! ???
869      INTEGER , INTENT( in  ), OPTIONAL        ::   kcom        ! input array
870      !!
871      INTEGER ::  ierror, iwork, ilocalcomm
872      !!----------------------------------------------------------------------
873      ilocalcomm = mpi_comm_opa
874      IF( PRESENT(kcom) )   ilocalcomm = kcom
875      CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_min, ilocalcomm, ierror )
876      ktab = iwork
877   END SUBROUTINE mppmin_int
878   !!
879   SUBROUTINE mppmin_a_real( ptab, kdim, kcom )
880      !!----------------------------------------------------------------------
881      INTEGER , INTENT(in   )                  ::   kdim
882      REAL(wp), INTENT(inout), DIMENSION(kdim) ::   ptab
883      INTEGER , INTENT(in   ), OPTIONAL        ::   kcom
884      INTEGER :: ierror, ilocalcomm
885      REAL(wp), DIMENSION(kdim) ::   zwork
886      !!-----------------------------------------------------------------------
887      ilocalcomm = mpi_comm_opa
888      IF( PRESENT(kcom) )   ilocalcomm = kcom
889      CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_min, ilocalcomm, ierror )
890      ptab(:) = zwork(:)
891   END SUBROUTINE mppmin_a_real
892   !!
893   SUBROUTINE mppmin_real( ptab, kcom )
894      !!-----------------------------------------------------------------------
895      REAL(wp), INTENT(inout)           ::   ptab        !
896      INTEGER , INTENT(in   ), OPTIONAL :: kcom
897      INTEGER  ::   ierror, ilocalcomm
898      REAL(wp) ::   zwork
899      !!-----------------------------------------------------------------------
900      ilocalcomm = mpi_comm_opa
901      IF( PRESENT(kcom) )   ilocalcomm = kcom
902      CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_min, ilocalcomm, ierror )
903      ptab = zwork
904   END SUBROUTINE mppmin_real
905
906
907   !!----------------------------------------------------------------------
908   !!    ***  mppsum_a_int, mppsum_int, mppsum_a_real, mppsum_real  ***
909   !!   
910   !!   Global sum of 1D array or a variable (integer, real or complex)
911   !!----------------------------------------------------------------------
912   !!
913   SUBROUTINE mppsum_a_int( ktab, kdim )
914      !!----------------------------------------------------------------------
915      INTEGER, INTENT(in   )                   ::   kdim   ! ???
916      INTEGER, INTENT(inout), DIMENSION (kdim) ::   ktab   ! ???
917      INTEGER :: ierror
918      INTEGER, DIMENSION (kdim) ::  iwork
919      !!----------------------------------------------------------------------
920      CALL mpi_allreduce( ktab, iwork, kdim, mpi_integer, mpi_sum, mpi_comm_opa, ierror )
921      ktab(:) = iwork(:)
922   END SUBROUTINE mppsum_a_int
923   !!
924   SUBROUTINE mppsum_int( ktab )
925      !!----------------------------------------------------------------------
926      INTEGER, INTENT(inout) ::   ktab
927      INTEGER :: ierror, iwork
928      !!----------------------------------------------------------------------
929      CALL mpi_allreduce( ktab, iwork, 1, mpi_integer, mpi_sum, mpi_comm_opa, ierror )
930      ktab = iwork
931   END SUBROUTINE mppsum_int
932   !!
933   SUBROUTINE mppsum_a_real( ptab, kdim, kcom )
934      !!-----------------------------------------------------------------------
935      INTEGER                  , INTENT(in   ) ::   kdim   ! size of ptab
936      REAL(wp), DIMENSION(kdim), INTENT(inout) ::   ptab   ! input array
937      INTEGER , OPTIONAL       , INTENT(in   ) ::   kcom   ! specific communicator
938      INTEGER  ::   ierror, ilocalcomm    ! local integer
939      REAL(wp) ::   zwork(kdim)           ! local workspace
940      !!-----------------------------------------------------------------------
941      ilocalcomm = mpi_comm_opa
942      IF( PRESENT(kcom) )   ilocalcomm = kcom
943      CALL mpi_allreduce( ptab, zwork, kdim, mpi_double_precision, mpi_sum, ilocalcomm, ierror )
944      ptab(:) = zwork(:)
945   END SUBROUTINE mppsum_a_real
946   !!
947   SUBROUTINE mppsum_real( ptab, kcom )
948      !!-----------------------------------------------------------------------
949      REAL(wp)          , INTENT(inout)           ::   ptab   ! input scalar
950      INTEGER , OPTIONAL, INTENT(in   ) ::   kcom
951      INTEGER  ::   ierror, ilocalcomm
952      REAL(wp) ::   zwork
953      !!-----------------------------------------------------------------------
954      ilocalcomm = mpi_comm_opa
955      IF( PRESENT(kcom) )   ilocalcomm = kcom
956      CALL mpi_allreduce( ptab, zwork, 1, mpi_double_precision, mpi_sum, ilocalcomm, ierror )
957      ptab = zwork
958   END SUBROUTINE mppsum_real
959   !!
960   SUBROUTINE mppsum_realdd( ytab, kcom )
961      !!-----------------------------------------------------------------------
962      COMPLEX(wp)          , INTENT(inout) ::   ytab    ! input scalar
963      INTEGER    , OPTIONAL, INTENT(in   ) ::   kcom
964      INTEGER     ::   ierror, ilocalcomm
965      COMPLEX(wp) ::   zwork
966      !!-----------------------------------------------------------------------
967      ilocalcomm = mpi_comm_opa
968      IF( PRESENT(kcom) )   ilocalcomm = kcom
969      CALL MPI_ALLREDUCE( ytab, zwork, 1, MPI_DOUBLE_COMPLEX, MPI_SUMDD, ilocalcomm, ierror )
970      ytab = zwork
971   END SUBROUTINE mppsum_realdd
972   !!
973   SUBROUTINE mppsum_a_realdd( ytab, kdim, kcom )
974      !!----------------------------------------------------------------------
975      INTEGER                     , INTENT(in   ) ::   kdim   ! size of ytab
976      COMPLEX(wp), DIMENSION(kdim), INTENT(inout) ::   ytab   ! input array
977      INTEGER    , OPTIONAL       , INTENT(in   ) ::   kcom
978      INTEGER:: ierror, ilocalcomm    ! local integer
979      COMPLEX(wp), DIMENSION(kdim) :: zwork     ! temporary workspace
980      !!-----------------------------------------------------------------------
981      ilocalcomm = mpi_comm_opa
982      IF( PRESENT(kcom) )   ilocalcomm = kcom
983      CALL MPI_ALLREDUCE( ytab, zwork, kdim, MPI_DOUBLE_COMPLEX, MPI_SUMDD, ilocalcomm, ierror )
984      ytab(:) = zwork(:)
985   END SUBROUTINE mppsum_a_realdd
986   
987
988   SUBROUTINE mppmax_real_multiple( pt1d, kdim, kcom  )
989      !!----------------------------------------------------------------------
990      !!                  ***  routine mppmax_real  ***
991      !!
992      !! ** Purpose :   Maximum across processor of each element of a 1D arrays
993      !!
994      !!----------------------------------------------------------------------
995      REAL(wp), DIMENSION(kdim), INTENT(inout) ::   pt1d   ! 1D arrays
996      INTEGER                  , INTENT(in   ) ::   kdim
997      INTEGER , OPTIONAL       , INTENT(in   ) ::   kcom   ! local communicator
998      !!
999      INTEGER  ::   ierror, ilocalcomm
1000      REAL(wp), DIMENSION(kdim) ::  zwork
1001      !!----------------------------------------------------------------------
1002      ilocalcomm = mpi_comm_opa
1003      IF( PRESENT(kcom) )   ilocalcomm = kcom
1004      !
1005      CALL mpi_allreduce( pt1d, zwork, kdim, mpi_double_precision, mpi_max, ilocalcomm, ierror )
1006      pt1d(:) = zwork(:)
1007      !
1008   END SUBROUTINE mppmax_real_multiple
1009
1010
1011   SUBROUTINE mpp_minloc2d( ptab, pmask, pmin, ki,kj )
1012      !!------------------------------------------------------------------------
1013      !!             ***  routine mpp_minloc  ***
1014      !!
1015      !! ** Purpose :   Compute the global minimum of an array ptab
1016      !!              and also give its global position
1017      !!
1018      !! ** Method  :   Use MPI_ALLREDUCE with MPI_MINLOC
1019      !!
1020      !!--------------------------------------------------------------------------
1021      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   ptab    ! Local 2D array
1022      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   pmask   ! Local mask
1023      REAL(wp)                     , INTENT(  out) ::   pmin    ! Global minimum of ptab
1024      INTEGER                      , INTENT(  out) ::   ki, kj  ! index of minimum in global frame
1025      !
1026      INTEGER :: ierror
1027      INTEGER , DIMENSION(2)   ::   ilocs
1028      REAL(wp) ::   zmin   ! local minimum
1029      REAL(wp), DIMENSION(2,1) ::   zain, zaout
1030      !!-----------------------------------------------------------------------
1031      !
1032      zmin  = MINVAL( ptab(:,:) , mask= pmask == 1._wp )
1033      ilocs = MINLOC( ptab(:,:) , mask= pmask == 1._wp )
1034      !
1035      ki = ilocs(1) + nimpp - 1
1036      kj = ilocs(2) + njmpp - 1
1037      !
1038      zain(1,:)=zmin
1039      zain(2,:)=ki+10000.*kj
1040      !
1041      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MINLOC,MPI_COMM_OPA,ierror)
1042      !
1043      pmin = zaout(1,1)
1044      kj = INT(zaout(2,1)/10000.)
1045      ki = INT(zaout(2,1) - 10000.*kj )
1046      !
1047   END SUBROUTINE mpp_minloc2d
1048
1049
1050   SUBROUTINE mpp_minloc3d( ptab, pmask, pmin, ki, kj ,kk)
1051      !!------------------------------------------------------------------------
1052      !!             ***  routine mpp_minloc  ***
1053      !!
1054      !! ** Purpose :   Compute the global minimum of an array ptab
1055      !!              and also give its global position
1056      !!
1057      !! ** Method  :   Use MPI_ALLREDUCE with MPI_MINLOC
1058      !!
1059      !!--------------------------------------------------------------------------
1060      REAL(wp), DIMENSION(:,:,:), INTENT(in   ) ::   ptab         ! Local 2D array
1061      REAL(wp), DIMENSION(:,:,:), INTENT(in   ) ::   pmask        ! Local mask
1062      REAL(wp)                  , INTENT(  out) ::   pmin         ! Global minimum of ptab
1063      INTEGER                   , INTENT(  out) ::   ki, kj, kk   ! index of minimum in global frame
1064      !
1065      INTEGER  ::   ierror
1066      REAL(wp) ::   zmin     ! local minimum
1067      INTEGER , DIMENSION(3)   ::   ilocs
1068      REAL(wp), DIMENSION(2,1) ::   zain, zaout
1069      !!-----------------------------------------------------------------------
1070      !
1071      zmin  = MINVAL( ptab(:,:,:) , mask= pmask == 1._wp )
1072      ilocs = MINLOC( ptab(:,:,:) , mask= pmask == 1._wp )
1073      !
1074      ki = ilocs(1) + nimpp - 1
1075      kj = ilocs(2) + njmpp - 1
1076      kk = ilocs(3)
1077      !
1078      zain(1,:) = zmin
1079      zain(2,:) = ki + 10000.*kj + 100000000.*kk
1080      !
1081      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MINLOC,MPI_COMM_OPA,ierror)
1082      !
1083      pmin = zaout(1,1)
1084      kk   = INT( zaout(2,1) / 100000000. )
1085      kj   = INT( zaout(2,1) - kk * 100000000. ) / 10000
1086      ki   = INT( zaout(2,1) - kk * 100000000. -kj * 10000. )
1087      !
1088   END SUBROUTINE mpp_minloc3d
1089
1090
1091   SUBROUTINE mpp_maxloc2d( ptab, pmask, pmax, ki, kj )
1092      !!------------------------------------------------------------------------
1093      !!             ***  routine mpp_maxloc  ***
1094      !!
1095      !! ** Purpose :   Compute the global maximum of an array ptab
1096      !!              and also give its global position
1097      !!
1098      !! ** Method  :   Use MPI_ALLREDUCE with MPI_MINLOC
1099      !!
1100      !!--------------------------------------------------------------------------
1101      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   ptab     ! Local 2D array
1102      REAL(wp), DIMENSION (jpi,jpj), INTENT(in   ) ::   pmask    ! Local mask
1103      REAL(wp)                     , INTENT(  out) ::   pmax     ! Global maximum of ptab
1104      INTEGER                      , INTENT(  out) ::   ki, kj   ! index of maximum in global frame
1105      !!
1106      INTEGER  :: ierror
1107      INTEGER, DIMENSION (2)   ::   ilocs
1108      REAL(wp) :: zmax   ! local maximum
1109      REAL(wp), DIMENSION(2,1) ::   zain, zaout
1110      !!-----------------------------------------------------------------------
1111      !
1112      zmax  = MAXVAL( ptab(:,:) , mask= pmask == 1._wp )
1113      ilocs = MAXLOC( ptab(:,:) , mask= pmask == 1._wp )
1114      !
1115      ki = ilocs(1) + nimpp - 1
1116      kj = ilocs(2) + njmpp - 1
1117      !
1118      zain(1,:) = zmax
1119      zain(2,:) = ki + 10000. * kj
1120      !
1121      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MAXLOC,MPI_COMM_OPA,ierror)
1122      !
1123      pmax = zaout(1,1)
1124      kj   = INT( zaout(2,1) / 10000.     )
1125      ki   = INT( zaout(2,1) - 10000.* kj )
1126      !
1127   END SUBROUTINE mpp_maxloc2d
1128
1129
1130   SUBROUTINE mpp_maxloc3d( ptab, pmask, pmax, ki, kj, kk )
1131      !!------------------------------------------------------------------------
1132      !!             ***  routine mpp_maxloc  ***
1133      !!
1134      !! ** Purpose :  Compute the global maximum of an array ptab
1135      !!              and also give its global position
1136      !!
1137      !! ** Method : Use MPI_ALLREDUCE with MPI_MINLOC
1138      !!
1139      !!--------------------------------------------------------------------------
1140      REAL(wp), DIMENSION (:,:,:), INTENT(in   ) ::   ptab         ! Local 2D array
1141      REAL(wp), DIMENSION (:,:,:), INTENT(in   ) ::   pmask        ! Local mask
1142      REAL(wp)                   , INTENT(  out) ::   pmax         ! Global maximum of ptab
1143      INTEGER                    , INTENT(  out) ::   ki, kj, kk   ! index of maximum in global frame
1144      !
1145      INTEGER  ::   ierror   ! local integer
1146      REAL(wp) ::   zmax     ! local maximum
1147      REAL(wp), DIMENSION(2,1) ::   zain, zaout
1148      INTEGER , DIMENSION(3)   ::   ilocs
1149      !!-----------------------------------------------------------------------
1150      !
1151      zmax  = MAXVAL( ptab(:,:,:) , mask= pmask == 1._wp )
1152      ilocs = MAXLOC( ptab(:,:,:) , mask= pmask == 1._wp )
1153      !
1154      ki = ilocs(1) + nimpp - 1
1155      kj = ilocs(2) + njmpp - 1
1156      kk = ilocs(3)
1157      !
1158      zain(1,:) = zmax
1159      zain(2,:) = ki + 10000.*kj + 100000000.*kk
1160      !
1161      CALL MPI_ALLREDUCE( zain,zaout, 1, MPI_2DOUBLE_PRECISION,MPI_MAXLOC,MPI_COMM_OPA,ierror)
1162      !
1163      pmax = zaout(1,1)
1164      kk   = INT( zaout(2,1) / 100000000. )
1165      kj   = INT( zaout(2,1) - kk * 100000000. ) / 10000
1166      ki   = INT( zaout(2,1) - kk * 100000000. -kj * 10000. )
1167      !
1168   END SUBROUTINE mpp_maxloc3d
1169
1170
1171   SUBROUTINE mppsync()
1172      !!----------------------------------------------------------------------
1173      !!                  ***  routine mppsync  ***
1174      !!
1175      !! ** Purpose :   Massively parallel processors, synchroneous
1176      !!
1177      !!-----------------------------------------------------------------------
1178      INTEGER :: ierror
1179      !!-----------------------------------------------------------------------
1180      !
1181      CALL mpi_barrier( mpi_comm_opa, ierror )
1182      !
1183   END SUBROUTINE mppsync
1184
1185
1186   SUBROUTINE mppstop
1187      !!----------------------------------------------------------------------
1188      !!                  ***  routine mppstop  ***
1189      !!
1190      !! ** purpose :   Stop massively parallel processors method
1191      !!
1192      !!----------------------------------------------------------------------
1193      INTEGER ::   info
1194      !!----------------------------------------------------------------------
1195      !
1196      CALL mppsync
1197      CALL mpi_finalize( info )
1198      !
1199   END SUBROUTINE mppstop
1200
1201
1202   SUBROUTINE mpp_comm_free( kcom )
1203      !!----------------------------------------------------------------------
1204      INTEGER, INTENT(in) ::   kcom
1205      !!
1206      INTEGER :: ierr
1207      !!----------------------------------------------------------------------
1208      !
1209      CALL MPI_COMM_FREE(kcom, ierr)
1210      !
1211   END SUBROUTINE mpp_comm_free
1212
1213
1214   SUBROUTINE mpp_ini_ice( pindic, kumout )
1215      !!----------------------------------------------------------------------
1216      !!               ***  routine mpp_ini_ice  ***
1217      !!
1218      !! ** Purpose :   Initialize special communicator for ice areas
1219      !!      condition together with global variables needed in the ddmpp folding
1220      !!
1221      !! ** Method  : - Look for ice processors in ice routines
1222      !!              - Put their number in nrank_ice
1223      !!              - Create groups for the world processors and the ice processors
1224      !!              - Create a communicator for ice processors
1225      !!
1226      !! ** output
1227      !!      njmppmax = njmpp for northern procs
1228      !!      ndim_rank_ice = number of processors with ice
1229      !!      nrank_ice (ndim_rank_ice) = ice processors
1230      !!      ngrp_iworld = group ID for the world processors
1231      !!      ngrp_ice = group ID for the ice processors
1232      !!      ncomm_ice = communicator for the ice procs.
1233      !!      n_ice_root = number (in the world) of proc 0 in the ice comm.
1234      !!
1235      !!----------------------------------------------------------------------
1236      INTEGER, INTENT(in) ::   pindic
1237      INTEGER, INTENT(in) ::   kumout   ! ocean.output logical unit
1238      !!
1239      INTEGER :: jjproc
1240      INTEGER :: ii, ierr
1241      INTEGER, ALLOCATABLE, DIMENSION(:) ::   kice
1242      INTEGER, ALLOCATABLE, DIMENSION(:) ::   zwork
1243      !!----------------------------------------------------------------------
1244      !
1245      ! Since this is just an init routine and these arrays are of length jpnij
1246      ! then don't use wrk_nemo module - just allocate and deallocate.
1247      ALLOCATE( kice(jpnij), zwork(jpnij), STAT=ierr )
1248      IF( ierr /= 0 ) THEN
1249         WRITE(kumout, cform_err)
1250         WRITE(kumout,*) 'mpp_ini_ice : failed to allocate 2, 1D arrays (jpnij in length)'
1251         CALL mppstop
1252      ENDIF
1253
1254      ! Look for how many procs with sea-ice
1255      !
1256      kice = 0
1257      DO jjproc = 1, jpnij
1258         IF( jjproc == narea .AND. pindic .GT. 0 )   kice(jjproc) = 1
1259      END DO
1260      !
1261      zwork = 0
1262      CALL MPI_ALLREDUCE( kice, zwork, jpnij, mpi_integer, mpi_sum, mpi_comm_opa, ierr )
1263      ndim_rank_ice = SUM( zwork )
1264
1265      ! Allocate the right size to nrank_north
1266      IF( ALLOCATED ( nrank_ice ) )   DEALLOCATE( nrank_ice )
1267      ALLOCATE( nrank_ice(ndim_rank_ice) )
1268      !
1269      ii = 0
1270      nrank_ice = 0
1271      DO jjproc = 1, jpnij
1272         IF( zwork(jjproc) == 1) THEN
1273            ii = ii + 1
1274            nrank_ice(ii) = jjproc -1
1275         ENDIF
1276      END DO
1277
1278      ! Create the world group
1279      CALL MPI_COMM_GROUP( mpi_comm_opa, ngrp_iworld, ierr )
1280
1281      ! Create the ice group from the world group
1282      CALL MPI_GROUP_INCL( ngrp_iworld, ndim_rank_ice, nrank_ice, ngrp_ice, ierr )
1283
1284      ! Create the ice communicator , ie the pool of procs with sea-ice
1285      CALL MPI_COMM_CREATE( mpi_comm_opa, ngrp_ice, ncomm_ice, ierr )
1286
1287      ! Find proc number in the world of proc 0 in the north
1288      ! The following line seems to be useless, we just comment & keep it as reminder
1289      ! CALL MPI_GROUP_TRANSLATE_RANKS(ngrp_ice,1,0,ngrp_iworld,n_ice_root,ierr)
1290      !
1291      CALL MPI_GROUP_FREE(ngrp_ice, ierr)
1292      CALL MPI_GROUP_FREE(ngrp_iworld, ierr)
1293
1294      DEALLOCATE(kice, zwork)
1295      !
1296   END SUBROUTINE mpp_ini_ice
1297
1298
1299   SUBROUTINE mpp_ini_znl( kumout )
1300      !!----------------------------------------------------------------------
1301      !!               ***  routine mpp_ini_znl  ***
1302      !!
1303      !! ** Purpose :   Initialize special communicator for computing zonal sum
1304      !!
1305      !! ** Method  : - Look for processors in the same row
1306      !!              - Put their number in nrank_znl
1307      !!              - Create group for the znl processors
1308      !!              - Create a communicator for znl processors
1309      !!              - Determine if processor should write znl files
1310      !!
1311      !! ** output
1312      !!      ndim_rank_znl = number of processors on the same row
1313      !!      ngrp_znl = group ID for the znl processors
1314      !!      ncomm_znl = communicator for the ice procs.
1315      !!      n_znl_root = number (in the world) of proc 0 in the ice comm.
1316      !!
1317      !!----------------------------------------------------------------------
1318      INTEGER, INTENT(in) ::   kumout   ! ocean.output logical units
1319      !
1320      INTEGER :: jproc      ! dummy loop integer
1321      INTEGER :: ierr, ii   ! local integer
1322      INTEGER, ALLOCATABLE, DIMENSION(:) ::   kwork
1323      !!----------------------------------------------------------------------
1324      !-$$     WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ngrp_world     : ', ngrp_world
1325      !-$$     WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - mpi_comm_world : ', mpi_comm_world
1326      !-$$     WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - mpi_comm_opa   : ', mpi_comm_opa
1327      !
1328      ALLOCATE( kwork(jpnij), STAT=ierr )
1329      IF( ierr /= 0 ) THEN
1330         WRITE(kumout, cform_err)
1331         WRITE(kumout,*) 'mpp_ini_znl : failed to allocate 1D array of length jpnij'
1332         CALL mppstop
1333      ENDIF
1334
1335      IF( jpnj == 1 ) THEN
1336         ngrp_znl  = ngrp_world
1337         ncomm_znl = mpi_comm_opa
1338      ELSE
1339         !
1340         CALL MPI_ALLGATHER ( njmpp, 1, mpi_integer, kwork, 1, mpi_integer, mpi_comm_opa, ierr )
1341         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - kwork pour njmpp : ', kwork
1342         !-$$        CALL flush(numout)
1343         !
1344         ! Count number of processors on the same row
1345         ndim_rank_znl = 0
1346         DO jproc=1,jpnij
1347            IF ( kwork(jproc) == njmpp ) THEN
1348               ndim_rank_znl = ndim_rank_znl + 1
1349            ENDIF
1350         END DO
1351         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ndim_rank_znl : ', ndim_rank_znl
1352         !-$$        CALL flush(numout)
1353         ! Allocate the right size to nrank_znl
1354         IF (ALLOCATED (nrank_znl)) DEALLOCATE(nrank_znl)
1355         ALLOCATE(nrank_znl(ndim_rank_znl))
1356         ii = 0
1357         nrank_znl (:) = 0
1358         DO jproc=1,jpnij
1359            IF ( kwork(jproc) == njmpp) THEN
1360               ii = ii + 1
1361               nrank_znl(ii) = jproc -1
1362            ENDIF
1363         END DO
1364         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - nrank_znl : ', nrank_znl
1365         !-$$        CALL flush(numout)
1366
1367         ! Create the opa group
1368         CALL MPI_COMM_GROUP(mpi_comm_opa,ngrp_opa,ierr)
1369         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ngrp_opa : ', ngrp_opa
1370         !-$$        CALL flush(numout)
1371
1372         ! Create the znl group from the opa group
1373         CALL MPI_GROUP_INCL  ( ngrp_opa, ndim_rank_znl, nrank_znl, ngrp_znl, ierr )
1374         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ngrp_znl ', ngrp_znl
1375         !-$$        CALL flush(numout)
1376
1377         ! Create the znl communicator from the opa communicator, ie the pool of procs in the same row
1378         CALL MPI_COMM_CREATE ( mpi_comm_opa, ngrp_znl, ncomm_znl, ierr )
1379         !-$$        WRITE (numout,*) 'mpp_ini_znl ', nproc, ' - ncomm_znl ', ncomm_znl
1380         !-$$        CALL flush(numout)
1381         !
1382      END IF
1383
1384      ! Determines if processor if the first (starting from i=1) on the row
1385      IF ( jpni == 1 ) THEN
1386         l_znl_root = .TRUE.
1387      ELSE
1388         l_znl_root = .FALSE.
1389         kwork (1) = nimpp
1390         CALL mpp_min ( kwork(1), kcom = ncomm_znl)
1391         IF ( nimpp == kwork(1)) l_znl_root = .TRUE.
1392      END IF
1393
1394      DEALLOCATE(kwork)
1395
1396   END SUBROUTINE mpp_ini_znl
1397
1398
1399   SUBROUTINE mpp_ini_north
1400      !!----------------------------------------------------------------------
1401      !!               ***  routine mpp_ini_north  ***
1402      !!
1403      !! ** Purpose :   Initialize special communicator for north folding
1404      !!      condition together with global variables needed in the mpp folding
1405      !!
1406      !! ** Method  : - Look for northern processors
1407      !!              - Put their number in nrank_north
1408      !!              - Create groups for the world processors and the north processors
1409      !!              - Create a communicator for northern processors
1410      !!
1411      !! ** output
1412      !!      njmppmax = njmpp for northern procs
1413      !!      ndim_rank_north = number of processors in the northern line
1414      !!      nrank_north (ndim_rank_north) = number  of the northern procs.
1415      !!      ngrp_world = group ID for the world processors
1416      !!      ngrp_north = group ID for the northern processors
1417      !!      ncomm_north = communicator for the northern procs.
1418      !!      north_root = number (in the world) of proc 0 in the northern comm.
1419      !!
1420      !!----------------------------------------------------------------------
1421      INTEGER ::   ierr
1422      INTEGER ::   jjproc
1423      INTEGER ::   ii, ji
1424      !!----------------------------------------------------------------------
1425      !
1426      njmppmax = MAXVAL( njmppt )
1427      !
1428      ! Look for how many procs on the northern boundary
1429      ndim_rank_north = 0
1430      DO jjproc = 1, jpnij
1431         IF( njmppt(jjproc) == njmppmax )   ndim_rank_north = ndim_rank_north + 1
1432      END DO
1433      !
1434      ! Allocate the right size to nrank_north
1435      IF (ALLOCATED (nrank_north)) DEALLOCATE(nrank_north)
1436      ALLOCATE( nrank_north(ndim_rank_north) )
1437
1438      ! Fill the nrank_north array with proc. number of northern procs.
1439      ! Note : the rank start at 0 in MPI
1440      ii = 0
1441      DO ji = 1, jpnij
1442         IF ( njmppt(ji) == njmppmax   ) THEN
1443            ii=ii+1
1444            nrank_north(ii)=ji-1
1445         END IF
1446      END DO
1447      !
1448      ! create the world group
1449      CALL MPI_COMM_GROUP( mpi_comm_opa, ngrp_world, ierr )
1450      !
1451      ! Create the North group from the world group
1452      CALL MPI_GROUP_INCL( ngrp_world, ndim_rank_north, nrank_north, ngrp_north, ierr )
1453      !
1454      ! Create the North communicator , ie the pool of procs in the north group
1455      CALL MPI_COMM_CREATE( mpi_comm_opa, ngrp_north, ncomm_north, ierr )
1456      !
1457   END SUBROUTINE mpp_ini_north
1458
1459
1460   SUBROUTINE mpp_lbc_north_e( pt2d, cd_type, psgn)
1461      !!---------------------------------------------------------------------
1462      !!                   ***  routine mpp_lbc_north_2d  ***
1463      !!
1464      !! ** Purpose :   Ensure proper north fold horizontal bondary condition
1465      !!              in mpp configuration in case of jpn1 > 1 and for 2d
1466      !!              array with outer extra halo
1467      !!
1468      !! ** Method  :   North fold condition and mpp with more than one proc
1469      !!              in i-direction require a specific treatment. We gather
1470      !!              the 4+2*jpr2dj northern lines of the global domain on 1
1471      !!              processor and apply lbc north-fold on this sub array.
1472      !!              Then we scatter the north fold array back to the processors.
1473      !!
1474      !!----------------------------------------------------------------------
1475      REAL(wp), DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj), INTENT(inout) ::   pt2d     ! 2D array with extra halo
1476      CHARACTER(len=1)                                            , INTENT(in   ) ::   cd_type  ! nature of pt3d grid-points
1477      REAL(wp)                                                    , INTENT(in   ) ::   psgn     ! sign used across the north fold
1478      !
1479      INTEGER ::   ji, jj, jr
1480      INTEGER ::   ierr, itaille, ildi, ilei, iilb
1481      INTEGER ::   ijpj, ij, iproc
1482      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE  ::  ztab_e, znorthloc_e
1483      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE  ::  znorthgloio_e
1484      !!----------------------------------------------------------------------
1485      !
1486      ALLOCATE( ztab_e(jpiglo,4+2*jpr2dj), znorthloc_e(jpi,4+2*jpr2dj), znorthgloio_e(jpi,4+2*jpr2dj,jpni) )
1487      !
1488      ijpj=4
1489      ztab_e(:,:) = 0._wp
1490
1491      ij = 0
1492      ! put in znorthloc_e the last 4 jlines of pt2d
1493      DO jj = nlcj - ijpj + 1 - jpr2dj, nlcj +jpr2dj
1494         ij = ij + 1
1495         DO ji = 1, jpi
1496            znorthloc_e(ji,ij) = pt2d(ji,jj)
1497         END DO
1498      END DO
1499      !
1500      itaille = jpi * ( ijpj + 2 * jpr2dj )
1501      CALL MPI_ALLGATHER( znorthloc_e(1,1)    , itaille, MPI_DOUBLE_PRECISION,    &
1502         &                znorthgloio_e(1,1,1), itaille, MPI_DOUBLE_PRECISION, ncomm_north, ierr )
1503      !
1504      DO jr = 1, ndim_rank_north            ! recover the global north array
1505         iproc = nrank_north(jr) + 1
1506         ildi  = nldit (iproc)
1507         ilei  = nleit (iproc)
1508         iilb  = nimppt(iproc)
1509         DO jj = 1, ijpj+2*jpr2dj
1510            DO ji = ildi, ilei
1511               ztab_e(ji+iilb-1,jj) = znorthgloio_e(ji,jj,jr)
1512            END DO
1513         END DO
1514      END DO
1515
1516      ! 2. North-Fold boundary conditions
1517      ! ----------------------------------
1518!!gm ERROR      CALL lbc_nfd( ztab_e(:,:), cd_type, psgn, pr2dj = jpr2dj )
1519
1520      ij = jpr2dj
1521      !! Scatter back to pt2d
1522      DO jj = nlcj - ijpj + 1 , nlcj +jpr2dj
1523      ij  = ij +1
1524         DO ji= 1, nlci
1525            pt2d(ji,jj) = ztab_e(ji+nimpp-1,ij)
1526         END DO
1527      END DO
1528      !
1529      DEALLOCATE( ztab_e, znorthloc_e, znorthgloio_e )
1530      !
1531   END SUBROUTINE mpp_lbc_north_e
1532
1533
1534   SUBROUTINE mpi_init_opa( ldtxt, ksft, code )
1535      !!---------------------------------------------------------------------
1536      !!                   ***  routine mpp_init.opa  ***
1537      !!
1538      !! ** Purpose :: export and attach a MPI buffer for bsend
1539      !!
1540      !! ** Method  :: define buffer size in namelist, if 0 no buffer attachment
1541      !!            but classical mpi_init
1542      !!
1543      !! History :: 01/11 :: IDRIS initial version for IBM only
1544      !!            08/04 :: R. Benshila, generalisation
1545      !!---------------------------------------------------------------------
1546      CHARACTER(len=*),DIMENSION(:), INTENT(  out) ::   ldtxt
1547      INTEGER                      , INTENT(inout) ::   ksft
1548      INTEGER                      , INTENT(  out) ::   code
1549      INTEGER                                      ::   ierr, ji
1550      LOGICAL                                      ::   mpi_was_called
1551      !!---------------------------------------------------------------------
1552      !
1553      CALL mpi_initialized( mpi_was_called, code )      ! MPI initialization
1554      IF ( code /= MPI_SUCCESS ) THEN
1555         DO ji = 1, SIZE(ldtxt)
1556            IF( TRIM(ldtxt(ji)) /= '' )   WRITE(*,*) ldtxt(ji)      ! control print of mynode
1557         END DO
1558         WRITE(*, cform_err)
1559         WRITE(*, *) ' lib_mpp: Error in routine mpi_initialized'
1560         CALL mpi_abort( mpi_comm_world, code, ierr )
1561      ENDIF
1562      !
1563      IF( .NOT. mpi_was_called ) THEN
1564         CALL mpi_init( code )
1565         CALL mpi_comm_dup( mpi_comm_world, mpi_comm_opa, code )
1566         IF ( code /= MPI_SUCCESS ) THEN
1567            DO ji = 1, SIZE(ldtxt)
1568               IF( TRIM(ldtxt(ji)) /= '' )   WRITE(*,*) ldtxt(ji)      ! control print of mynode
1569            END DO
1570            WRITE(*, cform_err)
1571            WRITE(*, *) ' lib_mpp: Error in routine mpi_comm_dup'
1572            CALL mpi_abort( mpi_comm_world, code, ierr )
1573         ENDIF
1574      ENDIF
1575      !
1576      IF( nn_buffer > 0 ) THEN
1577         WRITE(ldtxt(ksft),*) 'mpi_bsend, buffer allocation of  : ', nn_buffer   ;   ksft = ksft + 1
1578         ! Buffer allocation and attachment
1579         ALLOCATE( tampon(nn_buffer), stat = ierr )
1580         IF( ierr /= 0 ) THEN
1581            DO ji = 1, SIZE(ldtxt)
1582               IF( TRIM(ldtxt(ji)) /= '' )   WRITE(*,*) ldtxt(ji)      ! control print of mynode
1583            END DO
1584            WRITE(*, cform_err)
1585            WRITE(*, *) ' lib_mpp: Error in ALLOCATE', ierr
1586            CALL mpi_abort( mpi_comm_world, code, ierr )
1587         END IF
1588         CALL mpi_buffer_attach( tampon, nn_buffer, code )
1589      ENDIF
1590      !
1591   END SUBROUTINE mpi_init_opa
1592
1593
1594   SUBROUTINE DDPDD_MPI( ydda, yddb, ilen, itype )
1595      !!---------------------------------------------------------------------
1596      !!   Routine DDPDD_MPI: used by reduction operator MPI_SUMDD
1597      !!
1598      !!   Modification of original codes written by David H. Bailey
1599      !!   This subroutine computes yddb(i) = ydda(i)+yddb(i)
1600      !!---------------------------------------------------------------------
1601      INTEGER                     , INTENT(in)    ::   ilen, itype
1602      COMPLEX(wp), DIMENSION(ilen), INTENT(in)    ::   ydda
1603      COMPLEX(wp), DIMENSION(ilen), INTENT(inout) ::   yddb
1604      !
1605      REAL(wp) :: zerr, zt1, zt2    ! local work variables
1606      INTEGER  :: ji, ztmp           ! local scalar
1607      !!---------------------------------------------------------------------
1608      !
1609      ztmp = itype   ! avoid compilation warning
1610      !
1611      DO ji=1,ilen
1612      ! Compute ydda + yddb using Knuth's trick.
1613         zt1  = real(ydda(ji)) + real(yddb(ji))
1614         zerr = zt1 - real(ydda(ji))
1615         zt2  = ((real(yddb(ji)) - zerr) + (real(ydda(ji)) - (zt1 - zerr))) &
1616                + aimag(ydda(ji)) + aimag(yddb(ji))
1617
1618         ! The result is zt1 + zt2, after normalization.
1619         yddb(ji) = cmplx ( zt1 + zt2, zt2 - ((zt1 + zt2) - zt1),wp )
1620      END DO
1621      !
1622   END SUBROUTINE DDPDD_MPI
1623
1624
1625   SUBROUTINE mpp_lbc_north_icb( pt2d, cd_type, psgn, pr2dj)
1626      !!---------------------------------------------------------------------
1627      !!                   ***  routine mpp_lbc_north_icb  ***
1628      !!
1629      !! ** Purpose :   Ensure proper north fold horizontal bondary condition
1630      !!              in mpp configuration in case of jpn1 > 1 and for 2d
1631      !!              array with outer extra halo
1632      !!
1633      !! ** Method  :   North fold condition and mpp with more than one proc
1634      !!              in i-direction require a specific treatment. We gather
1635      !!              the 4+2*jpr2dj northern lines of the global domain on 1
1636      !!              processor and apply lbc north-fold on this sub array.
1637      !!              Then we scatter the north fold array back to the processors.
1638      !!              This version accounts for an extra halo with icebergs.
1639      !!
1640      !!----------------------------------------------------------------------
1641      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   pt2d     ! 2D array with extra halo
1642      CHARACTER(len=1)        , INTENT(in   ) ::   cd_type  ! nature of pt3d grid-points
1643      !                                                     !   = T ,  U , V , F or W -points
1644      REAL(wp)                , INTENT(in   ) ::   psgn     ! = -1. the sign change across the
1645      !!                                                    ! north fold, =  1. otherwise
1646      INTEGER, OPTIONAL       , INTENT(in   ) ::   pr2dj
1647      !
1648      INTEGER ::   ji, jj, jr
1649      INTEGER ::   ierr, itaille, ildi, ilei, iilb
1650      INTEGER ::   ijpj, ij, iproc, ipr2dj
1651      !
1652      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE  ::  ztab_e, znorthloc_e
1653      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE  ::  znorthgloio_e
1654      !!----------------------------------------------------------------------
1655      !
1656      ijpj=4
1657      IF( PRESENT(pr2dj) ) THEN           ! use of additional halos
1658         ipr2dj = pr2dj
1659      ELSE
1660         ipr2dj = 0
1661      ENDIF
1662      ALLOCATE( ztab_e(jpiglo,4+2*ipr2dj), znorthloc_e(jpi,4+2*ipr2dj), znorthgloio_e(jpi,4+2*ipr2dj,jpni) )
1663      !
1664      ztab_e(:,:) = 0._wp
1665      !
1666      ij = 0
1667      ! put in znorthloc_e the last 4 jlines of pt2d
1668      DO jj = nlcj - ijpj + 1 - ipr2dj, nlcj +ipr2dj
1669         ij = ij + 1
1670         DO ji = 1, jpi
1671            znorthloc_e(ji,ij)=pt2d(ji,jj)
1672         END DO
1673      END DO
1674      !
1675      itaille = jpi * ( ijpj + 2 * ipr2dj )
1676      CALL MPI_ALLGATHER( znorthloc_e(1,1)  , itaille, MPI_DOUBLE_PRECISION,    &
1677         &                znorthgloio_e(1,1,1), itaille, MPI_DOUBLE_PRECISION, ncomm_north, ierr )
1678      !
1679      DO jr = 1, ndim_rank_north            ! recover the global north array
1680         iproc = nrank_north(jr) + 1
1681         ildi = nldit (iproc)
1682         ilei = nleit (iproc)
1683         iilb = nimppt(iproc)
1684         DO jj = 1, ijpj+2*ipr2dj
1685            DO ji = ildi, ilei
1686               ztab_e(ji+iilb-1,jj) = znorthgloio_e(ji,jj,jr)
1687            END DO
1688         END DO
1689      END DO
1690
1691      ! 2. North-Fold boundary conditions
1692      ! ----------------------------------
1693!!gm ERROR      CALL lbc_nfd( ztab_e(:,:), cd_type, psgn, pr2dj = ipr2dj )
1694
1695      ij = ipr2dj
1696      !! Scatter back to pt2d
1697      DO jj = nlcj - ijpj + 1 , nlcj +ipr2dj
1698      ij  = ij +1
1699         DO ji= 1, nlci
1700            pt2d(ji,jj) = ztab_e(ji+nimpp-1,ij)
1701         END DO
1702      END DO
1703      !
1704      DEALLOCATE( ztab_e, znorthloc_e, znorthgloio_e )
1705      !
1706   END SUBROUTINE mpp_lbc_north_icb
1707
1708
1709   SUBROUTINE mpp_lnk_2d_icb( pt2d, cd_type, psgn, jpri, jprj )
1710      !!----------------------------------------------------------------------
1711      !!                  ***  routine mpp_lnk_2d_icb  ***
1712      !!
1713      !! ** Purpose :   Message passing manadgement for 2d array (with extra halo and icebergs)
1714      !!
1715      !! ** Method  :   Use mppsend and mpprecv function for passing mask
1716      !!      between processors following neighboring subdomains.
1717      !!            domain parameters
1718      !!                    nlci   : first dimension of the local subdomain
1719      !!                    nlcj   : second dimension of the local subdomain
1720      !!                    jpri   : number of rows for extra outer halo
1721      !!                    jprj   : number of columns for extra outer halo
1722      !!                    nbondi : mark for "east-west local boundary"
1723      !!                    nbondj : mark for "north-south local boundary"
1724      !!                    noea   : number for local neighboring processors
1725      !!                    nowe   : number for local neighboring processors
1726      !!                    noso   : number for local neighboring processors
1727      !!                    nono   : number for local neighboring processors
1728      !!----------------------------------------------------------------------
1729      REAL(wp), DIMENSION(1-jpri:jpi+jpri,1-jprj:jpj+jprj), INTENT(inout) ::   pt2d     ! 2D array with extra halo
1730      CHARACTER(len=1)                                    , INTENT(in   ) ::   cd_type  ! nature of ptab array grid-points
1731      REAL(wp)                                            , INTENT(in   ) ::   psgn     ! sign used across the north fold
1732      INTEGER                                             , INTENT(in   ) ::   jpri
1733      INTEGER                                             , INTENT(in   ) ::   jprj
1734      !
1735      INTEGER  ::   jl   ! dummy loop indices
1736      INTEGER  ::   imigr, iihom, ijhom        ! local integers
1737      INTEGER  ::   ipreci, iprecj             !   -       -
1738      INTEGER  ::   ml_req1, ml_req2, ml_err   ! for key_mpi_isend
1739      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat   ! for key_mpi_isend
1740      !!
1741      REAL(wp), DIMENSION(1-jpri:jpi+jpri,jprecj+jprj,2) ::   r2dns, r2dsn
1742      REAL(wp), DIMENSION(1-jprj:jpj+jprj,jpreci+jpri,2) ::   r2dwe, r2dew
1743      !!----------------------------------------------------------------------
1744
1745      ipreci = jpreci + jpri      ! take into account outer extra 2D overlap area
1746      iprecj = jprecj + jprj
1747
1748
1749      ! 1. standard boundary treatment
1750      ! ------------------------------
1751      ! Order matters Here !!!!
1752      !
1753      !                                      ! East-West boundaries
1754      !                                           !* Cyclic east-west
1755      IF( nbondi == 2 .AND. (nperio == 1 .OR. nperio == 4 .OR. nperio == 6) ) THEN
1756         pt2d(1-jpri:     1    ,:) = pt2d(jpim1-jpri:  jpim1 ,:)       ! east
1757         pt2d(   jpi  :jpi+jpri,:) = pt2d(     2      :2+jpri,:)       ! west
1758         !
1759      ELSE                                        !* closed
1760         IF( .NOT. cd_type == 'F' )   pt2d(  1-jpri   :jpreci    ,:) = 0._wp    ! south except at F-point
1761                                      pt2d(nlci-jpreci+1:jpi+jpri,:) = 0._wp    ! north
1762      ENDIF
1763      !
1764
1765      ! north fold treatment
1766      ! -----------------------
1767      IF( npolj /= 0 ) THEN
1768         !
1769         SELECT CASE ( jpni )
1770!!gm ERROR         CASE ( 1 )     ;   CALL lbc_nfd        ( pt2d(1:jpi,1:jpj+jprj), cd_type, psgn, pr2dj=jprj )
1771!!gm ERROR         CASE DEFAULT   ;   CALL mpp_lbc_north_icb( pt2d(1:jpi,1:jpj+jprj)  , cd_type, psgn , pr2dj=jprj  )
1772         END SELECT
1773         !
1774      ENDIF
1775
1776      ! 2. East and west directions exchange
1777      ! ------------------------------------
1778      ! we play with the neigbours AND the row number because of the periodicity
1779      !
1780      SELECT CASE ( nbondi )      ! Read Dirichlet lateral conditions
1781      CASE ( -1, 0, 1 )                ! all exept 2 (i.e. close case)
1782         iihom = nlci-nreci-jpri
1783         DO jl = 1, ipreci
1784            r2dew(:,jl,1) = pt2d(jpreci+jl,:)
1785            r2dwe(:,jl,1) = pt2d(iihom +jl,:)
1786         END DO
1787      END SELECT
1788      !
1789      !                           ! Migrations
1790      imigr = ipreci * ( jpj + 2*jprj)
1791      !
1792      SELECT CASE ( nbondi )
1793      CASE ( -1 )
1794         CALL mppsend( 2, r2dwe(1-jprj,1,1), imigr, noea, ml_req1 )
1795         CALL mpprecv( 1, r2dew(1-jprj,1,2), imigr, noea )
1796         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1797      CASE ( 0 )
1798         CALL mppsend( 1, r2dew(1-jprj,1,1), imigr, nowe, ml_req1 )
1799         CALL mppsend( 2, r2dwe(1-jprj,1,1), imigr, noea, ml_req2 )
1800         CALL mpprecv( 1, r2dew(1-jprj,1,2), imigr, noea )
1801         CALL mpprecv( 2, r2dwe(1-jprj,1,2), imigr, nowe )
1802         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1803         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
1804      CASE ( 1 )
1805         CALL mppsend( 1, r2dew(1-jprj,1,1), imigr, nowe, ml_req1 )
1806         CALL mpprecv( 2, r2dwe(1-jprj,1,2), imigr, nowe )
1807         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1808      END SELECT
1809      !
1810      !                           ! Write Dirichlet lateral conditions
1811      iihom = nlci - jpreci
1812      !
1813      SELECT CASE ( nbondi )
1814      CASE ( -1 )
1815         DO jl = 1, ipreci
1816            pt2d(iihom+jl,:) = r2dew(:,jl,2)
1817         END DO
1818      CASE ( 0 )
1819         DO jl = 1, ipreci
1820            pt2d(jl-jpri,:) = r2dwe(:,jl,2)
1821            pt2d( iihom+jl,:) = r2dew(:,jl,2)
1822         END DO
1823      CASE ( 1 )
1824         DO jl = 1, ipreci
1825            pt2d(jl-jpri,:) = r2dwe(:,jl,2)
1826         END DO
1827      END SELECT
1828
1829
1830      ! 3. North and south directions
1831      ! -----------------------------
1832      ! always closed : we play only with the neigbours
1833      !
1834      IF( nbondj /= 2 ) THEN      ! Read Dirichlet lateral conditions
1835         ijhom = nlcj-nrecj-jprj
1836         DO jl = 1, iprecj
1837            r2dsn(:,jl,1) = pt2d(:,ijhom +jl)
1838            r2dns(:,jl,1) = pt2d(:,jprecj+jl)
1839         END DO
1840      ENDIF
1841      !
1842      !                           ! Migrations
1843      imigr = iprecj * ( jpi + 2*jpri )
1844      !
1845      SELECT CASE ( nbondj )
1846      CASE ( -1 )
1847         CALL mppsend( 4, r2dsn(1-jpri,1,1), imigr, nono, ml_req1 )
1848         CALL mpprecv( 3, r2dns(1-jpri,1,2), imigr, nono )
1849         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1850      CASE ( 0 )
1851         CALL mppsend( 3, r2dns(1-jpri,1,1), imigr, noso, ml_req1 )
1852         CALL mppsend( 4, r2dsn(1-jpri,1,1), imigr, nono, ml_req2 )
1853         CALL mpprecv( 3, r2dns(1-jpri,1,2), imigr, nono )
1854         CALL mpprecv( 4, r2dsn(1-jpri,1,2), imigr, noso )
1855         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1856         IF(l_isend) CALL mpi_wait(ml_req2,ml_stat,ml_err)
1857      CASE ( 1 )
1858         CALL mppsend( 3, r2dns(1-jpri,1,1), imigr, noso, ml_req1 )
1859         CALL mpprecv( 4, r2dsn(1-jpri,1,2), imigr, noso )
1860         IF(l_isend) CALL mpi_wait(ml_req1,ml_stat,ml_err)
1861      END SELECT
1862      !
1863      !                           ! Write Dirichlet lateral conditions
1864      ijhom = nlcj - jprecj
1865      !
1866      SELECT CASE ( nbondj )
1867      CASE ( -1 )
1868         DO jl = 1, iprecj
1869            pt2d(:,ijhom+jl) = r2dns(:,jl,2)
1870         END DO
1871      CASE ( 0 )
1872         DO jl = 1, iprecj
1873            pt2d(:,jl-jprj) = r2dsn(:,jl,2)
1874            pt2d(:,ijhom+jl ) = r2dns(:,jl,2)
1875         END DO
1876      CASE ( 1 )
1877         DO jl = 1, iprecj
1878            pt2d(:,jl-jprj) = r2dsn(:,jl,2)
1879         END DO
1880      END SELECT
1881      !
1882   END SUBROUTINE mpp_lnk_2d_icb
1883   
1884#else
1885   !!----------------------------------------------------------------------
1886   !!   Default case:            Dummy module        share memory computing
1887   !!----------------------------------------------------------------------
1888   USE in_out_manager
1889
1890   INTERFACE mpp_sum
1891      MODULE PROCEDURE mpp_sum_a2s, mpp_sum_as, mpp_sum_ai, mpp_sum_s, mpp_sum_i, mppsum_realdd, mppsum_a_realdd
1892   END INTERFACE
1893   INTERFACE mpp_max
1894      MODULE PROCEDURE mppmax_a_int, mppmax_int, mppmax_a_real, mppmax_real
1895   END INTERFACE
1896   INTERFACE mpp_min
1897      MODULE PROCEDURE mppmin_a_int, mppmin_int, mppmin_a_real, mppmin_real
1898   END INTERFACE
1899   INTERFACE mpp_minloc
1900      MODULE PROCEDURE mpp_minloc2d ,mpp_minloc3d
1901   END INTERFACE
1902   INTERFACE mpp_maxloc
1903      MODULE PROCEDURE mpp_maxloc2d ,mpp_maxloc3d
1904   END INTERFACE
1905   INTERFACE mpp_max_multiple
1906      MODULE PROCEDURE mppmax_real_multiple
1907   END INTERFACE
1908
1909   LOGICAL, PUBLIC, PARAMETER ::   lk_mpp = .FALSE.      !: mpp flag
1910   LOGICAL, PUBLIC            ::   ln_nnogather          !: namelist control of northfold comms (needed here in case "key_mpp_mpi" is not used)
1911   INTEGER :: ncomm_ice
1912   INTEGER, PUBLIC            ::   mpi_comm_opa          ! opa local communicator
1913   !!----------------------------------------------------------------------
1914CONTAINS
1915
1916   INTEGER FUNCTION lib_mpp_alloc(kumout)          ! Dummy function
1917      INTEGER, INTENT(in) ::   kumout
1918      lib_mpp_alloc = 0
1919   END FUNCTION lib_mpp_alloc
1920
1921   FUNCTION mynode( ldtxt, ldname, kumnam_ref, knumnam_cfg,  kumond , kstop, localComm ) RESULT (function_value)
1922      INTEGER, OPTIONAL            , INTENT(in   ) ::   localComm
1923      CHARACTER(len=*),DIMENSION(:) ::   ldtxt
1924      CHARACTER(len=*) ::   ldname
1925      INTEGER ::   kumnam_ref, knumnam_cfg , kumond , kstop
1926      IF( PRESENT( localComm ) ) mpi_comm_opa = localComm
1927      function_value = 0
1928      IF( .FALSE. )   ldtxt(:) = 'never done'
1929      CALL ctl_opn( kumond, TRIM(ldname), 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. , 1 )
1930   END FUNCTION mynode
1931
1932   SUBROUTINE mppsync                       ! Dummy routine
1933   END SUBROUTINE mppsync
1934
1935   SUBROUTINE mpp_sum_as( parr, kdim, kcom )      ! Dummy routine
1936      REAL   , DIMENSION(:) :: parr
1937      INTEGER               :: kdim
1938      INTEGER, OPTIONAL     :: kcom
1939      WRITE(*,*) 'mpp_sum_as: You should not have seen this print! error?', kdim, parr(1), kcom
1940   END SUBROUTINE mpp_sum_as
1941
1942   SUBROUTINE mpp_sum_a2s( parr, kdim, kcom )      ! Dummy routine
1943      REAL   , DIMENSION(:,:) :: parr
1944      INTEGER               :: kdim
1945      INTEGER, OPTIONAL     :: kcom
1946      WRITE(*,*) 'mpp_sum_a2s: You should not have seen this print! error?', kdim, parr(1,1), kcom
1947   END SUBROUTINE mpp_sum_a2s
1948
1949   SUBROUTINE mpp_sum_ai( karr, kdim, kcom )      ! Dummy routine
1950      INTEGER, DIMENSION(:) :: karr
1951      INTEGER               :: kdim
1952      INTEGER, OPTIONAL     :: kcom
1953      WRITE(*,*) 'mpp_sum_ai: You should not have seen this print! error?', kdim, karr(1), kcom
1954   END SUBROUTINE mpp_sum_ai
1955
1956   SUBROUTINE mpp_sum_s( psca, kcom )            ! Dummy routine
1957      REAL                  :: psca
1958      INTEGER, OPTIONAL     :: kcom
1959      WRITE(*,*) 'mpp_sum_s: You should not have seen this print! error?', psca, kcom
1960   END SUBROUTINE mpp_sum_s
1961
1962   SUBROUTINE mpp_sum_i( kint, kcom )            ! Dummy routine
1963      integer               :: kint
1964      INTEGER, OPTIONAL     :: kcom
1965      WRITE(*,*) 'mpp_sum_i: You should not have seen this print! error?', kint, kcom
1966   END SUBROUTINE mpp_sum_i
1967
1968   SUBROUTINE mppsum_realdd( ytab, kcom )
1969      COMPLEX(wp), INTENT(inout)         :: ytab    ! input scalar
1970      INTEGER , INTENT( in  ), OPTIONAL :: kcom
1971      WRITE(*,*) 'mppsum_realdd: You should not have seen this print! error?', ytab
1972   END SUBROUTINE mppsum_realdd
1973
1974   SUBROUTINE mppsum_a_realdd( ytab, kdim, kcom )
1975      INTEGER , INTENT( in )                        ::   kdim      ! size of ytab
1976      COMPLEX(wp), DIMENSION(kdim), INTENT( inout ) ::   ytab      ! input array
1977      INTEGER , INTENT( in  ), OPTIONAL :: kcom
1978      WRITE(*,*) 'mppsum_a_realdd: You should not have seen this print! error?', kdim, ytab(1), kcom
1979   END SUBROUTINE mppsum_a_realdd
1980
1981   SUBROUTINE mppmax_a_real( parr, kdim, kcom )
1982      REAL   , DIMENSION(:) :: parr
1983      INTEGER               :: kdim
1984      INTEGER, OPTIONAL     :: kcom
1985      WRITE(*,*) 'mppmax_a_real: You should not have seen this print! error?', kdim, parr(1), kcom
1986   END SUBROUTINE mppmax_a_real
1987
1988   SUBROUTINE mppmax_real( psca, kcom )
1989      REAL                  :: psca
1990      INTEGER, OPTIONAL     :: kcom
1991      WRITE(*,*) 'mppmax_real: You should not have seen this print! error?', psca, kcom
1992   END SUBROUTINE mppmax_real
1993
1994   SUBROUTINE mppmin_a_real( parr, kdim, kcom )
1995      REAL   , DIMENSION(:) :: parr
1996      INTEGER               :: kdim
1997      INTEGER, OPTIONAL     :: kcom
1998      WRITE(*,*) 'mppmin_a_real: You should not have seen this print! error?', kdim, parr(1), kcom
1999   END SUBROUTINE mppmin_a_real
2000
2001   SUBROUTINE mppmin_real( psca, kcom )
2002      REAL                  :: psca
2003      INTEGER, OPTIONAL     :: kcom
2004      WRITE(*,*) 'mppmin_real: You should not have seen this print! error?', psca, kcom
2005   END SUBROUTINE mppmin_real
2006
2007   SUBROUTINE mppmax_a_int( karr, kdim ,kcom)
2008      INTEGER, DIMENSION(:) :: karr
2009      INTEGER               :: kdim
2010      INTEGER, OPTIONAL     :: kcom
2011      WRITE(*,*) 'mppmax_a_int: You should not have seen this print! error?', kdim, karr(1), kcom
2012   END SUBROUTINE mppmax_a_int
2013
2014   SUBROUTINE mppmax_int( kint, kcom)
2015      INTEGER               :: kint
2016      INTEGER, OPTIONAL     :: kcom
2017      WRITE(*,*) 'mppmax_int: You should not have seen this print! error?', kint, kcom
2018   END SUBROUTINE mppmax_int
2019
2020   SUBROUTINE mppmin_a_int( karr, kdim, kcom )
2021      INTEGER, DIMENSION(:) :: karr
2022      INTEGER               :: kdim
2023      INTEGER, OPTIONAL     :: kcom
2024      WRITE(*,*) 'mppmin_a_int: You should not have seen this print! error?', kdim, karr(1), kcom
2025   END SUBROUTINE mppmin_a_int
2026
2027   SUBROUTINE mppmin_int( kint, kcom )
2028      INTEGER               :: kint
2029      INTEGER, OPTIONAL     :: kcom
2030      WRITE(*,*) 'mppmin_int: You should not have seen this print! error?', kint, kcom
2031   END SUBROUTINE mppmin_int
2032
2033   SUBROUTINE mpp_minloc2d( ptab, pmask, pmin, ki, kj )
2034      REAL                   :: pmin
2035      REAL , DIMENSION (:,:) :: ptab, pmask
2036      INTEGER :: ki, kj
2037      WRITE(*,*) 'mpp_minloc2d: You should not have seen this print! error?', pmin, ki, kj, ptab(1,1), pmask(1,1)
2038   END SUBROUTINE mpp_minloc2d
2039
2040   SUBROUTINE mpp_minloc3d( ptab, pmask, pmin, ki, kj, kk )
2041      REAL                     :: pmin
2042      REAL , DIMENSION (:,:,:) :: ptab, pmask
2043      INTEGER :: ki, kj, kk
2044      WRITE(*,*) 'mpp_minloc3d: You should not have seen this print! error?', pmin, ki, kj, kk, ptab(1,1,1), pmask(1,1,1)
2045   END SUBROUTINE mpp_minloc3d
2046
2047   SUBROUTINE mpp_maxloc2d( ptab, pmask, pmax, ki, kj )
2048      REAL                   :: pmax
2049      REAL , DIMENSION (:,:) :: ptab, pmask
2050      INTEGER :: ki, kj
2051      WRITE(*,*) 'mpp_maxloc2d: You should not have seen this print! error?', pmax, ki, kj, ptab(1,1), pmask(1,1)
2052   END SUBROUTINE mpp_maxloc2d
2053
2054   SUBROUTINE mpp_maxloc3d( ptab, pmask, pmax, ki, kj, kk )
2055      REAL                     :: pmax
2056      REAL , DIMENSION (:,:,:) :: ptab, pmask
2057      INTEGER :: ki, kj, kk
2058      WRITE(*,*) 'mpp_maxloc3d: You should not have seen this print! error?', pmax, ki, kj, kk, ptab(1,1,1), pmask(1,1,1)
2059   END SUBROUTINE mpp_maxloc3d
2060
2061   SUBROUTINE mppstop
2062      STOP      ! non MPP case, just stop the run
2063   END SUBROUTINE mppstop
2064
2065   SUBROUTINE mpp_ini_ice( kcom, knum )
2066      INTEGER :: kcom, knum
2067      WRITE(*,*) 'mpp_ini_ice: You should not have seen this print! error?', kcom, knum
2068   END SUBROUTINE mpp_ini_ice
2069
2070   SUBROUTINE mpp_ini_znl( knum )
2071      INTEGER :: knum
2072      WRITE(*,*) 'mpp_ini_znl: You should not have seen this print! error?', knum
2073   END SUBROUTINE mpp_ini_znl
2074
2075   SUBROUTINE mpp_comm_free( kcom )
2076      INTEGER :: kcom
2077      WRITE(*,*) 'mpp_comm_free: You should not have seen this print! error?', kcom
2078   END SUBROUTINE mpp_comm_free
2079   
2080   SUBROUTINE mppmax_real_multiple( ptab, kdim , kcom  )
2081      REAL, DIMENSION(:) ::   ptab   !
2082      INTEGER            ::   kdim   !
2083      INTEGER, OPTIONAL  ::   kcom   !
2084      WRITE(*,*) 'mppmax_real_multiple: You should not have seen this print! error?', ptab(1), kdim
2085   END SUBROUTINE mppmax_real_multiple
2086
2087#endif
2088
2089   !!----------------------------------------------------------------------
2090   !!   All cases:         ctl_stop, ctl_warn, get_unit, ctl_opn, ctl_nam   routines
2091   !!----------------------------------------------------------------------
2092
2093   SUBROUTINE ctl_stop( cd1, cd2, cd3, cd4, cd5 ,   &
2094      &                 cd6, cd7, cd8, cd9, cd10 )
2095      !!----------------------------------------------------------------------
2096      !!                  ***  ROUTINE  stop_opa  ***
2097      !!
2098      !! ** Purpose :   print in ocean.outpput file a error message and
2099      !!                increment the error number (nstop) by one.
2100      !!----------------------------------------------------------------------
2101      CHARACTER(len=*), INTENT(in), OPTIONAL ::  cd1, cd2, cd3, cd4, cd5
2102      CHARACTER(len=*), INTENT(in), OPTIONAL ::  cd6, cd7, cd8, cd9, cd10
2103      !!----------------------------------------------------------------------
2104      !
2105      nstop = nstop + 1
2106      IF(lwp) THEN
2107         WRITE(numout,cform_err)
2108         IF( PRESENT(cd1 ) )   WRITE(numout,*) cd1
2109         IF( PRESENT(cd2 ) )   WRITE(numout,*) cd2
2110         IF( PRESENT(cd3 ) )   WRITE(numout,*) cd3
2111         IF( PRESENT(cd4 ) )   WRITE(numout,*) cd4
2112         IF( PRESENT(cd5 ) )   WRITE(numout,*) cd5
2113         IF( PRESENT(cd6 ) )   WRITE(numout,*) cd6
2114         IF( PRESENT(cd7 ) )   WRITE(numout,*) cd7
2115         IF( PRESENT(cd8 ) )   WRITE(numout,*) cd8
2116         IF( PRESENT(cd9 ) )   WRITE(numout,*) cd9
2117         IF( PRESENT(cd10) )   WRITE(numout,*) cd10
2118      ENDIF
2119                               CALL FLUSH(numout    )
2120      IF( numstp     /= -1 )   CALL FLUSH(numstp    )
2121      IF( numrun     /= -1 )   CALL FLUSH(numrun    )
2122      IF( numevo_ice /= -1 )   CALL FLUSH(numevo_ice)
2123      !
2124      IF( cd1 == 'STOP' ) THEN
2125         IF(lwp) WRITE(numout,*)  'huge E-R-R-O-R : immediate stop'
2126         CALL mppstop()
2127      ENDIF
2128      !
2129   END SUBROUTINE ctl_stop
2130
2131
2132   SUBROUTINE ctl_warn( cd1, cd2, cd3, cd4, cd5,   &
2133      &                 cd6, cd7, cd8, cd9, cd10 )
2134      !!----------------------------------------------------------------------
2135      !!                  ***  ROUTINE  stop_warn  ***
2136      !!
2137      !! ** Purpose :   print in ocean.outpput file a error message and
2138      !!                increment the warning number (nwarn) by one.
2139      !!----------------------------------------------------------------------
2140      CHARACTER(len=*), INTENT(in), OPTIONAL ::  cd1, cd2, cd3, cd4, cd5
2141      CHARACTER(len=*), INTENT(in), OPTIONAL ::  cd6, cd7, cd8, cd9, cd10
2142      !!----------------------------------------------------------------------
2143      !
2144      nwarn = nwarn + 1
2145      IF(lwp) THEN
2146         WRITE(numout,cform_war)
2147         IF( PRESENT(cd1 ) ) WRITE(numout,*) cd1
2148         IF( PRESENT(cd2 ) ) WRITE(numout,*) cd2
2149         IF( PRESENT(cd3 ) ) WRITE(numout,*) cd3
2150         IF( PRESENT(cd4 ) ) WRITE(numout,*) cd4
2151         IF( PRESENT(cd5 ) ) WRITE(numout,*) cd5
2152         IF( PRESENT(cd6 ) ) WRITE(numout,*) cd6
2153         IF( PRESENT(cd7 ) ) WRITE(numout,*) cd7
2154         IF( PRESENT(cd8 ) ) WRITE(numout,*) cd8
2155         IF( PRESENT(cd9 ) ) WRITE(numout,*) cd9
2156         IF( PRESENT(cd10) ) WRITE(numout,*) cd10
2157      ENDIF
2158      CALL FLUSH(numout)
2159      !
2160   END SUBROUTINE ctl_warn
2161
2162
2163   SUBROUTINE ctl_opn( knum, cdfile, cdstat, cdform, cdacce, klengh, kout, ldwp, karea )
2164      !!----------------------------------------------------------------------
2165      !!                  ***  ROUTINE ctl_opn  ***
2166      !!
2167      !! ** Purpose :   Open file and check if required file is available.
2168      !!
2169      !! ** Method  :   Fortan open
2170      !!----------------------------------------------------------------------
2171      INTEGER          , INTENT(  out) ::   knum      ! logical unit to open
2172      CHARACTER(len=*) , INTENT(in   ) ::   cdfile    ! file name to open
2173      CHARACTER(len=*) , INTENT(in   ) ::   cdstat    ! disposition specifier
2174      CHARACTER(len=*) , INTENT(in   ) ::   cdform    ! formatting specifier
2175      CHARACTER(len=*) , INTENT(in   ) ::   cdacce    ! access specifier
2176      INTEGER          , INTENT(in   ) ::   klengh    ! record length
2177      INTEGER          , INTENT(in   ) ::   kout      ! number of logical units for write
2178      LOGICAL          , INTENT(in   ) ::   ldwp      ! boolean term for print
2179      INTEGER, OPTIONAL, INTENT(in   ) ::   karea     ! proc number
2180      !
2181      CHARACTER(len=80) ::   clfile
2182      INTEGER           ::   iost
2183      !!----------------------------------------------------------------------
2184      !
2185      ! adapt filename
2186      ! ----------------
2187      clfile = TRIM(cdfile)
2188      IF( PRESENT( karea ) ) THEN
2189         IF( karea > 1 )   WRITE(clfile, "(a,'_',i4.4)") TRIM(clfile), karea-1
2190      ENDIF
2191#if defined key_agrif
2192      IF( .NOT. Agrif_Root() )   clfile = TRIM(Agrif_CFixed())//'_'//TRIM(clfile)
2193      knum=Agrif_Get_Unit()
2194#else
2195      knum=get_unit()
2196#endif
2197      !
2198      iost=0
2199      IF( cdacce(1:6) == 'DIRECT' )  THEN
2200         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat, RECL=klengh, ERR=100, IOSTAT=iost )
2201      ELSE
2202         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat             , ERR=100, IOSTAT=iost )
2203      ENDIF
2204      IF( iost == 0 ) THEN
2205         IF(ldwp) THEN
2206            WRITE(kout,*) '     file   : ', clfile,' open ok'
2207            WRITE(kout,*) '     unit   = ', knum
2208            WRITE(kout,*) '     status = ', cdstat
2209            WRITE(kout,*) '     form   = ', cdform
2210            WRITE(kout,*) '     access = ', cdacce
2211            WRITE(kout,*)
2212         ENDIF
2213      ENDIF
2214100   CONTINUE
2215      IF( iost /= 0 ) THEN
2216         IF(ldwp) THEN
2217            WRITE(kout,*)
2218            WRITE(kout,*) ' ===>>>> : bad opening file: ', clfile
2219            WRITE(kout,*) ' =======   ===  '
2220            WRITE(kout,*) '           unit   = ', knum
2221            WRITE(kout,*) '           status = ', cdstat
2222            WRITE(kout,*) '           form   = ', cdform
2223            WRITE(kout,*) '           access = ', cdacce
2224            WRITE(kout,*) '           iostat = ', iost
2225            WRITE(kout,*) '           we stop. verify the file '
2226            WRITE(kout,*)
2227         ENDIF
2228         CALL FLUSH( kout ) 
2229         STOP 'ctl_opn bad opening'
2230      ENDIF
2231      !
2232   END SUBROUTINE ctl_opn
2233
2234
2235   SUBROUTINE ctl_nam ( kios, cdnam, ldwp )
2236      !!----------------------------------------------------------------------
2237      !!                  ***  ROUTINE ctl_nam  ***
2238      !!
2239      !! ** Purpose :   Informations when error while reading a namelist
2240      !!
2241      !! ** Method  :   Fortan open
2242      !!----------------------------------------------------------------------
2243      INTEGER         , INTENT(inout) ::   kios    ! IO status after reading the namelist
2244      CHARACTER(len=*), INTENT(in   ) ::   cdnam   ! group name of namelist for which error occurs
2245      CHARACTER(len=5)                ::   clios   ! string to convert iostat in character for print
2246      LOGICAL         , INTENT(in   ) ::   ldwp    ! boolean term for print
2247      !!----------------------------------------------------------------------
2248      !
2249      WRITE (clios, '(I5.0)')   kios
2250      IF( kios < 0 ) THEN         
2251         CALL ctl_warn( 'end of record or file while reading namelist '   &
2252            &           // TRIM(cdnam) // ' iostat = ' // TRIM(clios) )
2253      ENDIF
2254      !
2255      IF( kios > 0 ) THEN
2256         CALL ctl_stop( 'misspelled variable in namelist '   &
2257            &           // TRIM(cdnam) // ' iostat = ' // TRIM(clios) )
2258      ENDIF
2259      kios = 0
2260      RETURN
2261      !
2262   END SUBROUTINE ctl_nam
2263
2264
2265   INTEGER FUNCTION get_unit()
2266      !!----------------------------------------------------------------------
2267      !!                  ***  FUNCTION  get_unit  ***
2268      !!
2269      !! ** Purpose :   return the index of an unused logical unit
2270      !!----------------------------------------------------------------------
2271      LOGICAL :: llopn
2272      !!----------------------------------------------------------------------
2273      !
2274      get_unit = 15   ! choose a unit that is big enough then it is not already used in NEMO
2275      llopn = .TRUE.
2276      DO WHILE( (get_unit < 998) .AND. llopn )
2277         get_unit = get_unit + 1
2278         INQUIRE( unit = get_unit, opened = llopn )
2279      END DO
2280      IF( (get_unit == 999) .AND. llopn ) THEN
2281         CALL ctl_stop( 'get_unit: All logical units until 999 are used...' )
2282         get_unit = -1
2283      ENDIF
2284      !
2285   END FUNCTION get_unit
2286
2287   !!----------------------------------------------------------------------
2288END MODULE lib_mpp
Note: See TracBrowser for help on using the repository browser.