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

source: branches/UKMO/dev_merge_2017_GC_couple_pkg/NEMOGCM/NEMO/OPA_SRC/LBC/lib_mpp.F90 @ 9679

Last change on this file since 9679 was 9679, checked in by dancopsey, 6 years ago

Merge in r8183 version of this branch (dev_r8183_GC_couple_pkg [8730:8734])

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