LCOV - code coverage report
Current view: top level - shared/common/src/12_hide_mpi - m_xmpi.F90 (source / functions) Coverage Total Hit
Test: coverage.info Lines: 37.6 % 1245 468
Test Date: 2026-09-21 22:40:37 Functions: 53.8 % 91 49

            Line data    Source code
       1              : !!****m* ABINIT/m_xmpi
       2              : !! NAME
       3              : !!  m_xmpi
       4              : !!
       5              : !! FUNCTION
       6              : !!  This module provides MPI named constants, tools for inquiring the MPI environment
       7              : !!  and a set of generic interfaces wrapping the most commonly used MPI primitives.
       8              : !!
       9              : !! COPYRIGHT
      10              : !! Copyright (C) 2009-2026 ABINIT group (MG, MB, XG, YP, MT)
      11              : !! This file is distributed under the terms of the
      12              : !! GNU General Public License, see ~abinit/COPYING
      13              : !! or http://www.gnu.org/copyleft/gpl.txt .
      14              : !!
      15              : !! TODO
      16              : !!  Get rid of xmpi_paral. Sequential code is the **exception**. Developers should code parallel
      17              : !!  code or code that is compatible both with MPI and seq (thanks to the wrappers provided by this module)
      18              : !!
      19              : !! SOURCE
      20              : 
      21              : #if defined HAVE_CONFIG_H
      22              : #include "config.h"
      23              : #endif
      24              : 
      25              : #include "abi_common.h"
      26              : 
      27              : module m_xmpi
      28              : 
      29              :  use, intrinsic :: iso_c_binding
      30              : #ifdef HAVE_FC_ISO_FORTRAN_2008
      31              :  use ISO_FORTRAN_ENV, only : int16, int32, int64
      32              : #endif
      33              :  USE_MPI
      34              :  use defs_basis
      35              :  use m_profiling_abi
      36              : #ifdef FC_NAG
      37              :  use f90_unix_proc
      38              : #endif
      39              :  use m_clib
      40              : 
      41              :  implicit none
      42              : 
      43              :  private
      44              : !!***
      45              : 
      46              : #ifdef HAVE_MPI1
      47              :  include 'mpif.h'
      48              : #endif
      49              : 
      50              : #ifndef HAVE_FC_ISO_FORTRAN_2008
      51              :  integer,parameter :: int16=2,int32=4,int64=8
      52              : #endif
      53              : 
      54              : #ifdef HAVE_MPI
      55              :  ! MPI constants used in abinit. Make sure that a corresponding fake value is provided for the sequential version.
      56              :  integer,public,parameter :: xmpi_paral          = 1
      57              :  integer,public,parameter :: xmpi_world          = MPI_COMM_WORLD
      58              :  integer,public,parameter :: xmpi_comm_self      = MPI_COMM_SELF
      59              :  integer,public,parameter :: xmpi_undefined      = MPI_UNDEFINED
      60              :  integer,public,parameter :: xmpi_undefined_rank = MPI_UNDEFINED  ! MPI_UNDEFINED_RANK is not portable.
      61              :  integer,public,parameter :: xmpi_comm_null      = MPI_COMM_NULL
      62              :  integer,public,parameter :: xmpi_group_null     = MPI_GROUP_NULL
      63              :  integer,public,parameter :: xmpi_any_source     = MPI_ANY_SOURCE
      64              :  integer,public,parameter :: xmpi_request_null   = MPI_REQUEST_NULL
      65              :  integer,public,parameter :: xmpi_msg_len        = MPI_MAX_ERROR_STRING ! Length of fortran string used to store MPI error strings.
      66              :  integer,public,parameter :: xmpi_info_null      = MPI_INFO_NULL
      67              :  integer,public,parameter :: xmpi_success        = MPI_SUCCESS
      68              :  integer,public,parameter :: xmpi_max_processor_name = MPI_MAX_PROCESSOR_NAME
      69              :  integer,public,parameter :: XMPI_MODE_NOPRECEDE = MPI_MODE_NOPRECEDE
      70              :  integer,public,parameter :: XMPI_MODE_NOSTORE   = MPI_MODE_NOSTORE
      71              :  integer,public,parameter :: XMPI_MODE_NOPUT     = MPI_MODE_NOPUT
      72              :  integer,public,parameter :: XMPI_MODE_NOSUCCEED = MPI_MODE_NOSUCCEED
      73              : 
      74              : #else
      75              :  ! Fake replacements for the sequential version. Values are taken from
      76              :  ! http://www.mit.edu/course/13/13.715/sun-hpc-ct-8.2.1/Linux/sun/include/mpif-common.h
      77              :  ! Please use these conventions when adding new replacements in order to avoid collisions between values.
      78              :  integer,public,parameter :: xmpi_paral          = 0
      79              :  integer,public,parameter :: xmpi_world          = 0
      80              :  integer,public,parameter :: xmpi_comm_self      = 1
      81              :  integer,public,parameter :: xmpi_undefined      =-32766
      82              :  integer,public,parameter :: xmpi_undefined_rank =-32766
      83              :  integer,public,parameter :: xmpi_comm_null      = 2
      84              :  integer,public,parameter :: xmpi_group_null     = 0
      85              :  integer,public,parameter :: xmpi_any_source     = -1
      86              :  integer,public,parameter :: xmpi_request_null   = 0
      87              :  integer,public,parameter :: xmpi_msg_len        = 1000
      88              :  integer,public,parameter :: xmpi_info_null      = 0
      89              :  integer,public,parameter :: xmpi_success        = 0
      90              :  integer,public,parameter :: xmpi_max_processor_name = 128
      91              :  integer,public,parameter :: XMPI_MODE_NOPRECEDE = 1
      92              :  integer,public,parameter :: XMPI_MODE_NOSTORE   = 2
      93              :  integer,public,parameter :: XMPI_MODE_NOPUT     = 8
      94              :  integer,public,parameter :: XMPI_MODE_NOSUCCEED = 32
      95              : #endif
      96              : 
      97              : #ifdef HAVE_MPI
      98              :  integer,save,private  :: xmpi_tag_ub=32767
      99              :  ! The tag upper bound value must be at least 32767. An MPI implementation is free to make
     100              :  ! the value of MPI_TAG_UB larger than this hence xmpi_tag_ub is redefined when MPI is init in xmpi_init.
     101              : #endif
     102              : 
     103              :  ! Size in bytes of the entries used in MPI datatypes.
     104              :  integer,save, public ABI_PROTECTED:: xmpi_bsize_ch  = 0
     105              :  integer,save, public ABI_PROTECTED:: xmpi_bsize_int = 0
     106              :  integer,save, public ABI_PROTECTED:: xmpi_bsize_sp  = 0
     107              :  integer,save, public ABI_PROTECTED:: xmpi_bsize_dp  = 0
     108              :  integer,save, public ABI_PROTECTED:: xmpi_bsize_spc = 0
     109              :  integer,save, public ABI_PROTECTED:: xmpi_bsize_dpc = 0
     110              : 
     111              :  ! kind of the offset used for MPI-IO.
     112              : #ifdef HAVE_MPI_IO
     113              :  integer,public,parameter :: xmpi_offset_kind  = MPI_OFFSET_KIND
     114              :  integer,public,parameter :: xmpi_address_kind = MPI_ADDRESS_KIND
     115              :  integer,public,parameter :: xmpi_mpiio = 1
     116              : #else
     117              :  integer,public,parameter :: xmpi_offset_kind = i8b
     118              :  integer,public,parameter :: xmpi_address_kind = i8b
     119              :  integer,public,parameter :: xmpi_mpiio = 0
     120              : #endif
     121              : 
     122              :  ! The byte size and the MPI type of the Fortran record marker.
     123              :  ! These quantities are compiler-dependent and are initalized here
     124              :  ! for selected compilers or in xmpio_get_info_frm that is called by xmpi_init (only if MPI-IO is on).
     125              : #if defined HAVE_MPI && (defined FC_INTEL || defined FC_GNU || defined FC_IBM)
     126              :  integer,save,public ABI_PROTECTED :: xmpio_bsize_frm   = 4
     127              :  integer,save,public ABI_PROTECTED :: xmpio_mpi_type_frm= MPI_INTEGER4
     128              : #else
     129              :  integer,save,public ABI_PROTECTED :: xmpio_bsize_frm    = 0
     130              :  integer,save,public ABI_PROTECTED :: xmpio_mpi_type_frm = 0
     131              : #endif
     132              : 
     133              :  integer,save, public ABI_PROTECTED :: xmpio_info = xmpi_info_null
     134              :  ! Global variable used to pass hints to the MPI-IO routines.
     135              : 
     136              :  integer(XMPI_OFFSET_KIND),public,parameter :: xmpio_chunk_bsize = 2000 * (1024.0_dp**2)
     137              :  ! Defines the chunk size (in bytes) used to (read|write) data in a single MPI-IO call.
     138              :  ! MPI-IO, indeed, crashes if we try to do the IO of a large array with a single call.
     139              :  ! We use a value <= 2  Gb to avoid wraparound errors with standard integers.
     140              : 
     141              :  ! Options used for the MPI-IO wrappers used in abinit.
     142              :  integer,public,parameter :: xmpio_single     = 1  ! Individual IO.
     143              :  integer,public,parameter :: xmpio_collective = 2  ! Collective IO.
     144              : 
     145              :  integer,save, public ABI_PROTECTED :: xmpi_count_requests = 0
     146              :  ! Count number of requests (+1 for each call to non-blocking API, -1 for each call to xmpi_wait)
     147              :  ! This counter should be zero at the end of the run if all requests have been released.
     148              : 
     149              :  integer,save, public ABI_PROTECTED :: xmpi_count_wins = 0
     150              :  ! Count number of windows created
     151              :  ! This counter should be zero at the end of the run if all windows have been released.
     152              : 
     153              :  logical,save, private :: xmpi_use_inplace_operations = .True.
     154              :  ! Enable/disable usage of MPI_IN_PLACE in e.g. xmpi_sum
     155              : 
     156              :  ! For MPI < v4, collective communication routines accept only a 32bit integer as data count.
     157              :  ! To exchange more than 2^32 data we need to create specific user-defined datatypes
     158              :  ! For this, we need some parameters:
     159              :  integer(KIND=int32),public,parameter :: xmpi_maxint32 = huge(0_int32)
     160              :  integer(KIND=int64),public,parameter :: xmpi_maxint32_64 = int(xmpi_maxint32,kind=int64)
     161              :  ! Max. integer that can be represented with 32 bits
     162              :  integer(KIND=int64),public,save :: xmpi_largetype_size = 0
     163              :  ! Number of data to be used in user-defined operations related to user-defined "largetype" type
     164              : !!***
     165              : 
     166              : !----------------------------------------------------------------------
     167              : 
     168              : !!****t* m_xmpi/xcomm_t
     169              : !! NAME
     170              : !! xcomm_t
     171              : !!
     172              : !! FUNCTION
     173              : !!  A small object storing the MPI communicator, the rank of the process and the size of the communicator.
     174              : !!  Provides helper functions to perform typical operations and parallelize loops.
     175              : !!  The datatype is initialized with xmpi_comm_self.
     176              : !!
     177              : !! SOURCE
     178              : 
     179              :  type, public :: xcomm_t
     180              :    integer :: value = xmpi_comm_self
     181              :    integer :: nproc = 1
     182              :    integer :: me = 0
     183              :    integer,private :: can_use_shmem__ = -1
     184              :     ! -1 --> unitialized, 0 if ranks do not belong to a shared memory region else 1
     185              : 
     186              :  contains
     187              :    procedure :: skip => xcomm_skip                     ! Skip iteration according to rank
     188              :    procedure :: set_to_null => xcomm_set_to_null       ! Init object using xmpi_comm_null.
     189              :    procedure :: set_to_self => xcomm_set_to_self       ! Init object using xmpi_comm_self.
     190              :    procedure :: free => xcomm_free                     ! Free the communicator.
     191              :    procedure :: from_cart_sub => xcomm_from_cart_sub   ! Build sub-communicators in a Cartesian grid.
     192              :    procedure :: split_type => xcomm_split_type         ! Creates new communicators based on split types and keys
     193              :    procedure :: prep_gatherv => xcomm_prep_gatherv     ! Prepare a typical gatherv operation.
     194              :    procedure :: print_names => xcomm_print_names
     195              :    procedure :: can_use_shmem => xcomm_can_use_shmem   ! true if communicator can use shared memory.
     196              :    procedure :: allocate_shared_master => xcomm_allocate_shared_master  ! Allocate MPI shared memory
     197              :  end type xcomm_t
     198              : 
     199              :  public :: xcomm_from_mpi_int
     200              : !!***
     201              : 
     202              : !----------------------------------------------------------------------
     203              : 
     204              : !!****t* m_xmpi/xmpi_pool2d_t
     205              : !! NAME
     206              : !! xmpi_pool2d_t
     207              : !!
     208              : !! FUNCTION
     209              : !!  Pool of MPI processors operating on a 2D problem of shape (n1, n2).
     210              : !!  Each item in the (n1, n2) matrix is assigned to a single pool.
     211              : !!  Note that differerent pools do not necessarily have the same number of procs,
     212              : !!  thus a pool is more flexibile than a Cartesian grid although inter-pool communication becomes more complex.
     213              : !!
     214              : !! SOURCE
     215              : 
     216              :  type, public :: xmpi_pool2d_t
     217              : 
     218              :    integer :: n1 = -1, n2 = -1
     219              :    ! Dimensions of the 2d problem
     220              : 
     221              :    type(xcomm_t) :: comm
     222              :    ! MPI communicator.
     223              : 
     224              :    logical,allocatable :: treats(:,:)
     225              :    ! (n1, n2)
     226              :    ! True if this pool treats (i1, i2)
     227              : 
     228              :  contains
     229              :    procedure :: from_dims => pool2d_from_dims     ! Init pool from problem dims.
     230              :    procedure :: free => pool2d_free               ! Free memory.
     231              :  end type xmpi_pool2d_t
     232              : !!***
     233              : 
     234              : ! Public procedures.
     235              :  public :: xmpi_init                  ! Initialize the MPI environment.
     236              :  public :: xmpi_set_inplace_operations! Set internal flag to use MPI_IN_PLACE whenever possible.
     237              :  public :: xmpi_end                   ! Terminate the MPI environment.
     238              :  public :: xmpi_abort                 ! Hides MPI_ABORT from MPI library.
     239              :  public :: xmpi_show_info             ! Printout of the basic variables stored in this module (useful for debugging).
     240              :  public :: xmpi_group_free            ! Hides MPI_GROUP_FREE from MPI library.
     241              :  public :: xmpi_group_incl            ! Hides MPI_GROUP_INCL from MPI library.
     242              :  public :: xmpi_group_translate_ranks ! Hides MPI_GROUP_TRANSLATE_RANKS from MPI library.
     243              :  public :: xmpi_comm_create           ! Hides MPI_COMM_CREATE from MPI library.
     244              :  public :: xmpi_comm_rank             ! Hides MPI_COMM_RANK from MPI library.
     245              :  public :: xmpi_comm_size             ! Hides MPI_COMM_SIZE from MPI library.
     246              :  public :: xmpi_comm_free             ! Hides MPI_COMM_FREE from MPI library.
     247              :  public :: xmpi_comm_dup              ! Hides MPI_COMM_DUP from MPI library.
     248              :  public :: xmpi_comm_group            ! Hides MPI_COMM_GROUP from MPI library.
     249              :  public :: xmpi_comm_translate_ranks  ! Hides MPI_GROUP_TRANSLATE_RANKS from MPI library.
     250              :  public :: xmpi_comm_translate_rank   ! Translate one rank
     251              :  public :: xmpi_comm_split            ! Hides MPI_COMM_SPLIT from MPI library.
     252              :  public :: xmpi_subcomm               ! Creates a sub-communicator from an input communicator.
     253              :  public :: xmpi_comm_multiple_of      ! Creates sub-communicator with number of procs multiple of a certain number.
     254              :  public :: xmpi_barrier               ! Hides MPI_BARRIER from MPI library.
     255              :  public :: xmpi_name                  ! Hides MPI_NAME from MPI library.
     256              :  public :: xmpi_iprobe                ! Hides MPI_IPROBE from MPI library.
     257              :  public :: xmpi_wait                  ! Hides MPI_WAIT from MPI library.
     258              :  public :: xmpi_waitall               ! Hides MPI_WAITALL from MPI library.
     259              :  public :: xmpi_request_free          ! Hides MPI_REQUEST_FREE from MPI library.
     260              :  public :: xmpi_requests_add          ! Increase/decrement xmpi_count_requests internal counter
     261              :  public :: xmpi_comm_set_errhandler   ! Hides MPI_COMM_SET_ERRHANDLER from MPI library.
     262              :  public :: xmpi_error_string          ! Return a string describing the error from ierr.
     263              :  public :: xmpi_split_work            ! Splits tasks inside communicator using blocks
     264              :  public :: xmpi_split_block           ! Splits tasks inside communicator using block distribution.
     265              :  public :: xmpi_split_cyclic          ! Splits tasks inside communicator using cyclic distribution.
     266              :  public :: xmpi_split_list            ! Splits list of indices inside communicator using block distribution.
     267              :  public :: xmpi_distab                ! Fill table defining the distribution of the tasks according to the # of processors
     268              :  public :: xmpi_distrib_with_replicas ! Distribute tasks among MPI ranks (replicas are allowed)
     269              :  public :: xmpi_distrib_2d            ! Try to optimally distribute nprocs in a 2d grid of shape (n1, n2)
     270              :  public :: xmpi_split_nsppol          ! Distribute spins. Also create and return indirect mapping to spin index.
     271              : 
     272              : ! Private procedures.
     273              :  private :: xmpi_largetype_create      ! Build a large-count contiguous datatype (to handle a very large # of data)
     274              :  private :: xmpi_largetype_free        ! Release a large-count contiguous datatype
     275              : 
     276              :  interface xmpi_comm_free
     277              :    module procedure xmpi_comm_free_0D
     278              :    module procedure xmpi_comm_free_1D
     279              :    module procedure xmpi_comm_free_2D
     280              :    module procedure xmpi_comm_free_3D
     281              :  end interface xmpi_comm_free
     282              : 
     283              :  interface xmpi_waitall
     284              :    module procedure xmpi_waitall_1d
     285              :    module procedure xmpi_waitall_2d
     286              :  end interface xmpi_waitall
     287              : 
     288              :  interface xmpi_split_work
     289              :    module procedure xmpi_split_work_i4b
     290              :  end interface xmpi_split_work
     291              : 
     292              :  public :: xmpi_split_work2_i4b
     293              :  public :: xmpi_split_work2_i8b
     294              :  !public :: xmpi_split_work2
     295              :  !
     296              :  ! g95@green v0.93 is not able to resolve the interface.
     297              :  ! For the time being, this generic interface has been disabled.
     298              :  !interface xmpi_split_work2
     299              :  !  module procedure xmpi_split_work2_i4b
     300              :  !  module procedure xmpi_split_work2_i8b
     301              :  !end interface xmpi_split_work2
     302              : 
     303              :  interface xmpi_distab
     304              :    module procedure xmpi_distab_4D
     305              :  end interface xmpi_distab
     306              : 
     307              :  ! MPI generic interfaces.
     308              :  public :: xmpi_allgather
     309              :  public :: xmpi_iallgather
     310              :  public :: xmpi_allgatherv
     311              :  public :: xmpi_alltoall
     312              :  public :: xmpi_ialltoall
     313              :  public :: xmpi_alltoallv
     314              :  public :: xmpi_ialltoallv
     315              :  public :: xmpi_bcast
     316              :  public :: xmpi_ibcast
     317              :  public :: xmpi_exch
     318              :  public :: xmpi_gather
     319              :  public :: xmpi_gatherv
     320              :  public :: xmpi_max
     321              :  public :: xmpi_max_ip
     322              :  public :: xmpi_min         ! Out-of-place version
     323              :  public :: xmpi_min_ip      ! In-place version
     324              :  public :: xmpi_recv
     325              :  public :: xmpi_irecv
     326              :  public :: xmpi_scatterv
     327              :  public :: xmpi_send
     328              :  public :: xmpi_isend
     329              :  public :: xmpi_sum_master
     330              :  public :: xmpi_sum
     331              :  public :: xmpi_isum
     332              :  public :: xmpi_isum_ip
     333              :  public :: xmpi_land              ! allreduce with MPI_LAND
     334              :  public :: xmpi_lor               ! allreduce with MPI_LOR
     335              : 
     336              :  public :: xmpi_win_fence
     337              :  public :: xmpi_win_free
     338              : 
     339              : #ifdef HAVE_MPI_IO
     340              :  public :: xmpio_max_address      !  Returns .TRUE. if offset cannot be stored in integer(kind=XMPI_ADDRESS_KIND).
     341              :  public :: xmpio_type_struct
     342              :  public :: xmpio_get_info_frm
     343              :  public :: xmpio_check_frmarkers
     344              :  public :: xmpio_read_frm
     345              :  public :: xmpio_read_int
     346              :  public :: xmpio_read_dp
     347              :  public :: xmpio_write_frm
     348              :  public :: xmpio_write_frmarkers
     349              : 
     350              :  public :: xmpio_create_fstripes
     351              :  public :: xmpio_create_fsubarray_2D
     352              :  public :: xmpio_create_fsubarray_3D
     353              :  public :: xmpio_create_fsubarray_4D
     354              :  public :: xmpio_create_fherm_packed
     355              :  public :: xmpio_create_coldistr_from_fpacked
     356              :  public :: xmpio_create_coldistr_from_fp3blocks
     357              : #endif
     358              : 
     359              : !----------------------------------------------------------------------
     360              : 
     361              : interface xmpi_allgather
     362              :   module procedure xmpi_allgather_int
     363              :   module procedure xmpi_allgather_char
     364              :   module procedure xmpi_allgather_int1d_1b
     365              :   module procedure xmpi_allgather_int1d
     366              :   module procedure xmpi_allgather_int2d
     367              :   module procedure xmpi_allgather_dp1d
     368              :   module procedure xmpi_allgather_dp2d
     369              :   module procedure xmpi_allgather_dp3d
     370              :   module procedure xmpi_allgather_dp4d
     371              : end interface xmpi_allgather
     372              : 
     373              : !----------------------------------------------------------------------
     374              : 
     375              : ! non-blocking version (requires MPI3)
     376              : ! Prototype:
     377              : !
     378              : !   call xmpi_iallgather(SENDBUF, SENDCOUNT, SENDTYPE, RECVBUF, RECVCOUNT, RECVTYPE, COMM, REQUEST, IERROR)
     379              : !
     380              : ! If the MPI library does not provide ialltoall, we call the blocking version and
     381              : ! we return xmpi_request_null (see xmpi_iallgather.finc)
     382              : ! Client code should always test/wait the request so that code semantics is preserved.
     383              : 
     384              : interface xmpi_iallgather
     385              :   module procedure xmpi_iallgather_dp4d
     386              : end interface xmpi_iallgather
     387              : 
     388              : interface xmpi_allgatherv
     389              :   module procedure xmpi_allgatherv_int2d
     390              :   module procedure xmpi_allgatherv_int
     391              :   module procedure xmpi_allgatherv_int1_dp1
     392              :   module procedure xmpi_allgatherv_dp
     393              :   module procedure xmpi_allgatherv_dp2d
     394              :   module procedure xmpi_allgatherv_dp3d
     395              :   module procedure xmpi_allgatherv_dp4d
     396              :   module procedure xmpi_allgatherv_dp5d
     397              :   module procedure xmpi_allgatherv_dp6d
     398              :   module procedure xmpi_allgatherv_coeff2d
     399              :   module procedure xmpi_allgatherv_coeff2d_indx
     400              :   module procedure xmpi_allgatherv_dc
     401              : end interface xmpi_allgatherv
     402              : 
     403              : !----------------------------------------------------------------------
     404              : 
     405              : ! blocking
     406              : interface xmpi_alltoall
     407              :   module procedure xmpi_alltoall_int
     408              :   module procedure xmpi_alltoall_dp2d
     409              :   module procedure xmpi_alltoall_dp4d
     410              : end interface xmpi_alltoall
     411              : 
     412              : ! non-blocking version (requires MPI3)
     413              : ! Prototype:
     414              : !
     415              : !   call xmpi_ialltoall(xval, sendsize, recvbuf, recvsize, comm, request)
     416              : !
     417              : ! If the MPI library does not provide ialltoall, we call the blocking version and
     418              : ! we return xmpi_request_null (see xmpi_ialltoall.finc)
     419              : ! Client code should always test/wait the request so that code semantics is preserved.
     420              : 
     421              : interface xmpi_ialltoall
     422              :   module procedure xmpi_ialltoall_dp4d
     423              : end interface xmpi_ialltoall
     424              : 
     425              : !----------------------------------------------------------------------
     426              : 
     427              : interface xmpi_alltoallv
     428              :   module procedure xmpi_alltoallv_dp2d
     429              :   module procedure xmpi_alltoallv_int2d
     430              :   module procedure xmpi_alltoallv_dp1d
     431              :   module procedure xmpi_alltoallv_dp1d2
     432              : end interface xmpi_alltoallv
     433              : 
     434              : !----------------------------------------------------------------------
     435              : 
     436              : ! non-blocking version (requires MPI3)
     437              : ! Prototype:
     438              : !
     439              : !   call xmpi_ialltoallv(xval, sendcnts, sdispls, recvbuf, recvcnts, rdispls, comm, request)
     440              : !
     441              : ! If the MPI library does not provide ialltoallv, we call the blocking version and
     442              : ! we return xmpi_request_null (see xmpi_ialltoallv.finc)
     443              : ! Client code should always test/wait the request so that code semantics is preserved.
     444              : 
     445              : interface xmpi_ialltoallv
     446              :   module procedure xmpi_ialltoallv_dp2d
     447              :   module procedure xmpi_ialltoallv_int2d
     448              :   module procedure xmpi_ialltoallv_dp1d2
     449              : end interface xmpi_ialltoallv
     450              : 
     451              : !----------------------------------------------------------------------
     452              : 
     453              : interface xmpi_bcast
     454              :   module procedure xmpi_bcast_intv
     455              :   module procedure xmpi_bcast_int1d
     456              :   module procedure xmpi_bcast_int2d
     457              :   module procedure xmpi_bcast_int3d
     458              :   module procedure xmpi_bcast_int4d
     459              :   module procedure xmpi_bcast_dpv
     460              :   module procedure xmpi_bcast_dp1d
     461              :   module procedure xmpi_bcast_dp2d
     462              :   module procedure xmpi_bcast_dp3d
     463              :   module procedure xmpi_bcast_dp4d
     464              :   module procedure xmpi_bcast_dp5d
     465              :   module procedure xmpi_bcast_dp6d
     466              :   module procedure xmpi_bcast_spv
     467              :   module procedure xmpi_bcast_sp1d
     468              :   module procedure xmpi_bcast_sp2d
     469              :   module procedure xmpi_bcast_sp3d
     470              :   module procedure xmpi_bcast_sp4d
     471              :   module procedure xmpi_bcast_cplxv
     472              :   module procedure xmpi_bcast_cplx1d
     473              :   module procedure xmpi_bcast_cplx2d
     474              :   module procedure xmpi_bcast_cplx3d
     475              :   module procedure xmpi_bcast_cplx4d
     476              :   module procedure xmpi_bcast_dcv
     477              :   module procedure xmpi_bcast_dc1d
     478              :   module procedure xmpi_bcast_dc2d
     479              :   module procedure xmpi_bcast_dc3d
     480              :   module procedure xmpi_bcast_dc4d
     481              :   module procedure xmpi_bcast_ch0d
     482              :   module procedure xmpi_bcast_ch1d
     483              :   module procedure xmpi_bcast_log0d
     484              :   module procedure xmpi_bcast_coeffi2_1d
     485              :   module procedure xmpi_bcast_coeff2_1d
     486              : end interface xmpi_bcast
     487              : 
     488              : !----------------------------------------------------------------------
     489              : 
     490              : interface xmpi_ibcast
     491              :   module procedure xmpi_ibcast_int1d
     492              :   module procedure xmpi_ibcast_int4d
     493              :   module procedure xmpi_ibcast_dp1d
     494              :   module procedure xmpi_ibcast_dp2d
     495              :   module procedure xmpi_ibcast_dp3d
     496              :   module procedure xmpi_ibcast_dp4d
     497              :   module procedure xmpi_ibcast_dpc2d
     498              :   module procedure xmpi_ibcast_spc2d
     499              : end interface xmpi_ibcast
     500              : 
     501              : !----------------------------------------------------------------------
     502              : 
     503              : interface xmpi_exch
     504              :   module procedure xmpi_exch_int1d
     505              :   module procedure xmpi_exch_int2d
     506              :   module procedure xmpi_exch_dp1d
     507              :   module procedure xmpi_exch_dp2d
     508              :   module procedure xmpi_exch_dp3d
     509              :   module procedure xmpi_exch_dp4d
     510              :   module procedure xmpi_exch_dp5d
     511              :   module procedure xmpi_exch_spc1d
     512              :   module procedure xmpi_exch_dpc1d
     513              :   module procedure xmpi_exch_dpc2d
     514              : end interface xmpi_exch
     515              : 
     516              : !----------------------------------------------------------------------
     517              : 
     518              : interface xmpi_gather
     519              :   module procedure xmpi_gather_int
     520              :   module procedure xmpi_gather_int2d
     521              :   module procedure xmpi_gather_dp
     522              :   module procedure xmpi_gather_dp2d
     523              :   module procedure xmpi_gather_dp3d
     524              :   module procedure xmpi_gather_dp4d
     525              : end interface xmpi_gather
     526              : 
     527              : !----------------------------------------------------------------------
     528              : 
     529              : interface xmpi_gatherv
     530              :   module procedure xmpi_gatherv_int
     531              :   module procedure xmpi_gatherv_int1_dp1
     532              :   module procedure xmpi_gatherv_int2d
     533              :   module procedure xmpi_gatherv_dp
     534              :   module procedure xmpi_gatherv_dp2d
     535              :   module procedure xmpi_gatherv_dp3d
     536              :   module procedure xmpi_gatherv_dp4d
     537              :   module procedure xmpi_gatherv_dp5d
     538              :   module procedure xmpi_gatherv_dp6d
     539              :   module procedure xmpi_gatherv_dc
     540              : end interface xmpi_gatherv
     541              : 
     542              : !----------------------------------------------------------------------
     543              : 
     544              : interface xmpi_max
     545              :   module procedure xmpi_max_int0d_i4b
     546              :   module procedure xmpi_max_int0d_i8b
     547              :   module procedure xmpi_max_int
     548              :   module procedure xmpi_max_dpv
     549              :   module procedure xmpi_max_dp0d_ip
     550              :   module procedure xmpi_max_dp1d_ip
     551              : end interface xmpi_max
     552              : 
     553              : interface xmpi_max_ip
     554              :   module procedure xmpi_max_int0d_ip
     555              :   module procedure xmpi_max_int1d_ip
     556              :   module procedure xmpi_max_dp0d_ip
     557              :   module procedure xmpi_max_dp1d_ip
     558              : end interface xmpi_max_ip
     559              : 
     560              : !----------------------------------------------------------------------
     561              : 
     562              : interface xmpi_min
     563              :   module procedure xmpi_min_intv
     564              :   module procedure xmpi_min_dpv
     565              :   module procedure xmpi_min_dp
     566              : end interface xmpi_min
     567              : 
     568              : ! In-place version of xmpi_min
     569              : interface xmpi_min_ip
     570              :   module procedure xmpi_min_int1d
     571              :   module procedure xmpi_min_dp
     572              : end interface xmpi_min_ip
     573              : 
     574              : !----------------------------------------------------------------------
     575              : 
     576              : !interface xmpi_min_max
     577              : !  module procedure xmpi_min_max_int0d_i4b
     578              : !end interface xmpi_min_max
     579              : 
     580              : !----------------------------------------------------------------------
     581              : 
     582              : interface xmpi_recv
     583              :   module procedure xmpi_recv_char
     584              :   module procedure xmpi_recv_intv
     585              :   module procedure xmpi_recv_int1d
     586              :   module procedure xmpi_recv_int2d
     587              :   module procedure xmpi_recv_int3d
     588              :   module procedure xmpi_recv_dp
     589              :   module procedure xmpi_recv_dp1d
     590              :   module procedure xmpi_recv_dp2d
     591              :   module procedure xmpi_recv_dp3d
     592              :   module procedure xmpi_recv_dp4d
     593              : end interface xmpi_recv
     594              : 
     595              : !----------------------------------------------------------------------
     596              : 
     597              : interface xmpi_irecv
     598              :   module procedure xmpi_irecv_intv
     599              :   module procedure xmpi_irecv_int1d
     600              :   module procedure xmpi_irecv_dp1d
     601              :   module procedure xmpi_irecv_dp2d
     602              :   module procedure xmpi_irecv_dp3d
     603              : end interface xmpi_irecv
     604              : 
     605              : !----------------------------------------------------------------------
     606              : 
     607              : interface xmpi_scatterv
     608              :   module procedure xmpi_scatterv_int
     609              :   module procedure xmpi_scatterv_int2d
     610              :   module procedure xmpi_scatterv_dp
     611              :   module procedure xmpi_scatterv_dp2d
     612              :   module procedure xmpi_scatterv_dp3d
     613              :   module procedure xmpi_scatterv_dp4d
     614              : end interface xmpi_scatterv
     615              : 
     616              : !----------------------------------------------------------------------
     617              : 
     618              : interface xmpi_isend
     619              :   module procedure xmpi_isend_int1d
     620              :   module procedure xmpi_isend_dp1d
     621              :   module procedure xmpi_isend_dp2d
     622              :   module procedure xmpi_isend_dp3d
     623              : end interface xmpi_isend
     624              : 
     625              : !----------------------------------------------------------------------
     626              : 
     627              : interface xmpi_send
     628              :   module procedure xmpi_send_char
     629              :   module procedure xmpi_send_intv
     630              :   module procedure xmpi_send_int1d
     631              :   module procedure xmpi_send_int2d
     632              :   module procedure xmpi_send_int3d
     633              :   module procedure xmpi_send_dp
     634              :   module procedure xmpi_send_dp1d
     635              :   module procedure xmpi_send_dp2d
     636              :   module procedure xmpi_send_dp3d
     637              :   module procedure xmpi_send_dp4d
     638              : end interface xmpi_send
     639              : 
     640              : !----------------------------------------------------------------------
     641              : 
     642              : interface xmpi_sum_master
     643              :   module procedure xmpi_sum_master_int
     644              :   module procedure xmpi_sum_master_int2d
     645              :   module procedure xmpi_sum_master_int4d
     646              :   module procedure xmpi_sum_master_dp
     647              :   module procedure xmpi_sum_master_dp1d
     648              :   module procedure xmpi_sum_master_dp2d
     649              :   module procedure xmpi_sum_master_dp3d
     650              :   module procedure xmpi_sum_master_dp4d
     651              :   module procedure xmpi_sum_master_dp5d
     652              :   module procedure xmpi_sum_master_dp6d
     653              :   module procedure xmpi_sum_master_dp7d
     654              :   module procedure xmpi_sum_master_c1cplx
     655              :   module procedure xmpi_sum_master_c2cplx
     656              :   module procedure xmpi_sum_master_c3cplx
     657              :   module procedure xmpi_sum_master_c4cplx
     658              :   module procedure xmpi_sum_master_c5cplx
     659              :   module procedure xmpi_sum_master_c1dpc
     660              :   module procedure xmpi_sum_master_c2dpc
     661              :   module procedure xmpi_sum_master_c3dpc
     662              :   module procedure xmpi_sum_master_c4dpc
     663              :   module procedure xmpi_sum_master_c5dpc
     664              : end interface xmpi_sum_master
     665              : 
     666              : !----------------------------------------------------------------------
     667              : 
     668              : !MG:TODO procedure marked with !? are considered obsolete.
     669              : !   and will be removed in future versions.
     670              : !   Please use interfaces where array dimensions are not passed explicitly.
     671              : !   Rationale: The array descriptor is already passed to the routine
     672              : !   so it does not make sense to pass the dimension explicitly.
     673              : 
     674              : interface xmpi_sum
     675              :   module procedure xmpi_sum_int
     676              :   module procedure xmpi_sum_intv
     677              :   module procedure xmpi_sum_intv2
     678              :   module procedure xmpi_sum_intn   !?
     679              :   module procedure xmpi_sum_int2t  !?
     680              :   module procedure xmpi_sum_int2d
     681              :   module procedure xmpi_sum_int3d
     682              :   module procedure xmpi_sum_int4d
     683              :   module procedure xmpi_sum_dp
     684              :   module procedure xmpi_sum_dpvt
     685              :   module procedure xmpi_sum_dpv
     686              :   module procedure xmpi_sum_dpn    !?
     687              :   module procedure xmpi_sum_sp2d
     688              :   module procedure xmpi_sum_sp3d
     689              :   module procedure xmpi_sum_sp4d
     690              :   module procedure xmpi_sum_sp5d
     691              :   module procedure xmpi_sum_sp6d
     692              :   module procedure xmpi_sum_sp7d
     693              :   module procedure xmpi_sum_dp2d
     694              :   module procedure xmpi_sum_dp3d
     695              :   module procedure xmpi_sum_dp4d
     696              :   module procedure xmpi_sum_dp5d
     697              :   module procedure xmpi_sum_dp6d
     698              :   module procedure xmpi_sum_dp7d
     699              :   module procedure xmpi_sum_dp2t   !?
     700              :   module procedure xmpi_sum_dp2d2t
     701              :   module procedure xmpi_sum_dp3d2t !?
     702              :   module procedure xmpi_sum_dp4d2t !?
     703              :   module procedure xmpi_sum_c0dc
     704              :   module procedure xmpi_sum_c0sc
     705              :   module procedure xmpi_sum_c1dc
     706              :   module procedure xmpi_sum_c2dc
     707              :   module procedure xmpi_sum_c3dc
     708              :   module procedure xmpi_sum_c4dc
     709              :   module procedure xmpi_sum_c5dc
     710              :   module procedure xmpi_sum_c6dc
     711              :   module procedure xmpi_sum_c7dc
     712              :   module procedure xmpi_sum_c1cplx
     713              :   module procedure xmpi_sum_c2cplx
     714              :   module procedure xmpi_sum_c3cplx
     715              :   module procedure xmpi_sum_c4cplx
     716              :   module procedure xmpi_sum_c5cplx
     717              :   module procedure xmpi_sum_c6cplx
     718              :   module procedure xmpi_sum_coeff5d1
     719              :  end interface xmpi_sum
     720              : !!***
     721              : 
     722              : ! Non-blocking version
     723              : interface xmpi_isum
     724              :   module procedure xmpi_isum_int0d
     725              : end interface xmpi_isum
     726              : !!***
     727              : 
     728              : ! Non-blocking in-place version
     729              : interface xmpi_isum_ip
     730              :   module procedure xmpi_isum_ip_dp2d
     731              :   module procedure xmpi_isum_ip_dp3d
     732              :   module procedure xmpi_isum_ip_dp4d
     733              :   module procedure xmpi_isum_ip_spc1d
     734              :   module procedure xmpi_isum_ip_dpc1d
     735              :   module procedure xmpi_isum_ip_spc2d
     736              :   module procedure xmpi_isum_ip_dpc2d
     737              :   module procedure xmpi_isum_ip_spc3d
     738              :   module procedure xmpi_isum_ip_dpc3d
     739              : end interface xmpi_isum_ip
     740              : !!***
     741              : 
     742              : interface xmpi_land
     743              :   module procedure xmpi_land_log0d
     744              :   module procedure xmpi_land_log1d
     745              : end interface xmpi_land
     746              : !!***
     747              : 
     748              : interface xmpi_lor
     749              :   module procedure xmpi_lor_log1d
     750              :   module procedure xmpi_lor_log2d
     751              :   module procedure xmpi_lor_log3d
     752              : end interface xmpi_lor
     753              : !!!***
     754              : 
     755              : !----------------------------------------------------------------------
     756              : 
     757              : CONTAINS  !===========================================================
     758              : !!***
     759              : 
     760              : !!****f* m_xmpi/xmpi_init
     761              : !! NAME
     762              : !!  xmpi_init
     763              : !!
     764              : !! FUNCTION
     765              : !!  Hides MPI_INIT from MPI library. Perform the initialization of some basic variables
     766              : !!  used by the MPI routines employed in abinit.
     767              : !!
     768              : !! SOURCE
     769              : 
     770         3576 : subroutine xmpi_init()
     771              : 
     772              : !Local variables-------------------
     773              :  integer :: mpierr, ierr, unt
     774              :  integer(c_long) :: rlim_cur, rlim_max
     775              :  logical :: exists
     776              : #ifdef HAVE_MPI
     777              :  integer :: attribute_val
     778              :  logical :: lflag
     779              : #ifdef HAVE_OPENMP
     780              :  integer :: required,provided
     781              : #endif
     782              : #endif
     783              : ! *************************************************************************
     784              : 
     785              :  call set_num_threads_if_undef()
     786              : 
     787         1788 :  mpierr=0
     788              : #ifdef HAVE_MPI
     789              : 
     790              : #ifndef HAVE_OPENMP
     791         1788 :  call MPI_INIT(mpierr)
     792              : #else
     793              :  required = MPI_THREAD_SINGLE
     794              :  !required = MPI_THREAD_FUNNELED
     795              :  !required = MPI_THREAD_SERIALIZED
     796              :  !required = MPI_THREAD_MULTIPLE
     797              :  call MPI_INIT_THREAD(required,provided,mpierr)
     798              :  if (provided /= required) call xmpi_abort(msg="MPI_INIT_THREADS: provided /= required")
     799              : #endif
     800              : 
     801              :  !%comm_world = xmpi_world ! Needed to bypass a bug in some OMPI implementations (intent(inout))
     802              :  !%call xmpi_comm_set_errhandler(comm_world, MPI_ERRORS_RETURN, err_handler_sav, mpierr)
     803              : 
     804              :  ! Deprecated in MPI2 but not all MPI2 implementations provide MPI_Comm_get_attr !
     805         1788 :  call MPI_ATTR_GET(xmpi_world, MPI_TAG_UB, attribute_val, lflag, mpierr)
     806              :  !call MPI_Comm_get_attr(xmpi_world, MPI_TAG_UB, attribute_val, lflag, mpierr)
     807              : 
     808         1788 :  if (lflag) xmpi_tag_ub = attribute_val
     809              : 
     810              :  ! Define type values.
     811         1788 :  call MPI_TYPE_SIZE(MPI_CHARACTER, xmpi_bsize_ch, mpierr)
     812         1788 :  call MPI_TYPE_SIZE(MPI_INTEGER, xmpi_bsize_int, mpierr)
     813         1788 :  call MPI_TYPE_SIZE(MPI_REAL, xmpi_bsize_sp, mpierr)
     814         1788 :  call MPI_TYPE_SIZE(MPI_DOUBLE_PRECISION, xmpi_bsize_dp, mpierr)
     815         1788 :  call MPI_TYPE_SIZE(MPI_COMPLEX, xmpi_bsize_spc, mpierr)
     816         1788 :  call MPI_TYPE_SIZE(MPI_DOUBLE_COMPLEX, xmpi_bsize_dpc, mpierr)
     817              : 
     818              :  ! Find the byte size of Fortran record marker used in MPI-IO routines.
     819         1788 :  if (xmpio_bsize_frm == 0) call xmpio_get_info_frm(xmpio_bsize_frm, xmpio_mpi_type_frm, xmpi_world)
     820              : #endif
     821              : 
     822              :  ! Try to increase stack size.
     823         1788 :  call clib_ulimit_stack(rlim_cur, rlim_max, ierr)
     824              : 
     825         1788 :  if (xmpi_comm_rank(xmpi_world) == 0) then
     826              : 
     827         1486 :    if (ierr /= 0) then
     828            0 :      write(std_out, "(a)")" WARNING: cannot increase stack size limit. "
     829              :    end if
     830              :    ! rlim_cur/rlim_max are in bytes; -1 means RLIM_INFINITY (unlimited). Printed
     831              :    ! unconditionally, not just when ierr /= 0: clib_ulimit_stack's fallback to
     832              :    ! "soft = current hard limit" (see rlimit.c) can succeed (ierr == 0) at a
     833              :    ! hard limit far below what a large automatic/local array actually needs,
     834              :    ! which then only ever surfaces later as an unexplained OOM-kill rather
     835              :    ! than this warning.
     836              :    !write(std_out, "(2(a,i0),a)") &
     837              :    !  "- clib_ulimit_stack: rlim_cur= ", rlim_cur, " bytes, rlim_max= ", rlim_max, " bytes (-1 = RLIM_INFINITY)"
     838              : 
     839              :    ! Master Removes the ABI_MPIABORTFILE if present so that we start with a clean environment.
     840         1486 :    inquire(file=ABI_MPIABORTFILE, exist=exists)
     841         1486 :    if (exists) then
     842              :      ! Get free unit (emulate F2008 newunit for portability reasons)
     843            0 :      unt = xmpi_get_unit()
     844            0 :      if (unt == -1) call xmpi_abort(msg="Cannot find free unit!!")
     845            0 :      open(unit=unt, file=trim(ABI_MPIABORTFILE), status="old", iostat=ierr)
     846            0 :      if (ierr == 0) close(unit=unt, status="delete", iostat=ierr)
     847            0 :      if (ierr /= 0) call xmpi_abort(msg="Cannot remove ABI_MPIABORTFILE")
     848              :    end if
     849              : 
     850              :    ! If MPI interfaces are buggy, MPI_IN_PLACE is not allowed
     851              : #if defined HAVE_MPI2_INPLACE && defined HAVE_MPI_BUGGY_INTERFACES
     852              :    write(std_out, "(a)")"ERROR: Cannot use MPI_IN_PLACE with this buggy MPI version!"
     853              :    write(ab_out , "(a)")"ERROR: Cannot use MPI_IN_PLACE with this buggy MPI version!"
     854              :    call xmpi_abort(msg="Stopping here!")
     855              : #endif
     856              : 
     857              :  end if
     858              : 
     859         1788 : end subroutine xmpi_init
     860              : !!***
     861              : 
     862              : !----------------------------------------------------------------------
     863              : 
     864              : !!****f* m_xmpi/set_num_threads_if_undef
     865              : !! NAME
     866              : !!  set_num_threads_if_undef
     867              : !!
     868              : !! FUNCTION
     869              : !!  sets OMP_NUM_THREADS to 1 is the env variable is undefined.
     870              : !!
     871              : !! SOURCE
     872              : 
     873              : subroutine set_num_threads_if_undef()
     874              : 
     875              : #ifdef HAVE_OPENMP
     876              : !Local variables-------------------
     877              :  integer :: ierr
     878              :  character(len=100) :: omp_num_threads
     879              : ! *************************************************************************
     880              : 
     881              :  ! Get the value of OMP_NUM_THREADS environment variable
     882              :  call get_environment_variable('OMP_NUM_THREADS', omp_num_threads, status=ierr)
     883              : 
     884              :  ! If OMP_NUM_THREADS is not defined (ierr != 0), set it to 1
     885              :  if (ierr /= 0) then
     886              :    ierr = clib_setenv('OMP_NUM_THREADS', '1', 1)
     887              :    if (ierr == 0) then
     888              :      write(std_out,"(a)")'- OMP_NUM_THREADS was not defined. It has been set to 1.'
     889              :    else
     890              :      write(std_out,"(a)")'- WARNING: Failed to set OMP_NUM_THREADS.'
     891              :    end if
     892              :  else
     893              :    !write(std_out,*)'- OMP_NUM_THREADS is already set to: ', trim(omp_num_threads)
     894              :  end if
     895              : #endif
     896              : 
     897              : end subroutine set_num_threads_if_undef
     898              : !!***
     899              : 
     900              : !----------------------------------------------------------------------
     901              : 
     902              : !!****f* m_xmpi/xmpi_set_inplace_operations
     903              : !! NAME
     904              : !!  xmpi_set_inplace_operations
     905              : !!
     906              : !! FUNCTION
     907              : !!  Set internal flag to use MPI_IN_PLACE whenever possible.
     908              : !!
     909              : !! SOURCE
     910              : 
     911            0 : subroutine xmpi_set_inplace_operations(bool)
     912              : 
     913              : !Local variables-------------------
     914              :  logical,intent(in) :: bool
     915              : ! *************************************************************************
     916              : 
     917            0 :  xmpi_use_inplace_operations = bool
     918              : 
     919            0 : end subroutine xmpi_set_inplace_operations
     920              : !!***
     921              : 
     922              : !----------------------------------------------------------------------
     923              : 
     924              : !!****f* m_xmpi/xmpi_get_unit
     925              : !! NAME
     926              : !!  xmpi_get_unit
     927              : !!
     928              : !! FUNCTION
     929              : !! Get free unit (emulate F2008 newunit for portability reasons)
     930              : !! Return -1 if no unit is found.
     931              : !!
     932              : !! SOURCE
     933              : 
     934            0 : integer function xmpi_get_unit() result(unt)
     935              : 
     936              : !Local variables-------------------
     937              :  logical :: isopen
     938              : ! *************************************************************************
     939              : 
     940            0 :  do unt=1024,-1,-1
     941            0 :    inquire(unit=unt, opened=isopen)
     942            0 :    if (.not.isopen) exit
     943              :  end do
     944              : 
     945            0 : end function xmpi_get_unit
     946              : !!***
     947              : 
     948              : !----------------------------------------------------------------------
     949              : 
     950              : !!****f* m_xmpi/xmpi_end
     951              : !! NAME
     952              : !!  xmpi_end
     953              : !!
     954              : !! FUNCTION
     955              : !!  Hides MPI_FINALIZE from MPI library.
     956              : !!
     957              : !! SOURCE
     958              : 
     959         1772 : subroutine xmpi_end()
     960              : 
     961              : !Local variables-------------------
     962              :  integer :: mpierr
     963              : ! *************************************************************************
     964              : 
     965         1772 :  mpierr=0
     966              : #ifdef HAVE_MPI
     967         1772 :  call MPI_BARRIER(MPI_COMM_WORLD,mpierr)  !  Needed by some HPC architectures (MT, 20110315)
     968         1772 :  call MPI_FINALIZE(mpierr)
     969              : #endif
     970              : 
     971              : #ifndef FC_IBM
     972              :  ! IBM8 returns 260. 320 ...
     973              :  call sys_exit(0)
     974              : #endif
     975              : 
     976              : end subroutine xmpi_end
     977              : !!***
     978              : 
     979              : !----------------------------------------------------------------------
     980              : 
     981              : !!****f* m_xmpi/xmpi_abort
     982              : !! NAME
     983              : !!  xmpi_abort
     984              : !!
     985              : !! FUNCTION
     986              : !!  Hides MPI_ABORT from MPI library.
     987              : !!
     988              : !! INPUTS
     989              : !!  [comm]=communicator of tasks to abort.
     990              : !!  [mpierr]=Error code to return to invoking environment.
     991              : !!  [msg]=User message
     992              : !!  [exit_status]=optional, shell return code, default 1
     993              : !!
     994              : !! SOURCE
     995              : 
     996           16 : subroutine xmpi_abort(comm, mpierr, msg, exit_status)
     997              : 
     998              : !Arguments-------------------------
     999              :  integer,optional,intent(in) :: comm,mpierr,exit_status
    1000              :  character(len=*),optional,intent(in) :: msg
    1001              : 
    1002              : !Local variables-------------------
    1003              :  integer :: ierr,my_comm,my_errorcode,ilen,ierr2
    1004              :  logical :: testopen
    1005              :  character(len=xmpi_msg_len) :: mpi_msg_error
    1006              : ! *************************************************************************
    1007              : 
    1008           16 :  ierr=0
    1009           16 :  my_comm = xmpi_world; if (PRESENT(comm)) my_comm = comm
    1010              : 
    1011           16 :  if (PRESENT(msg)) then
    1012            0 :    write(std_out,'(2a)')"User message: ",TRIM(msg)
    1013              :  end if
    1014              : 
    1015              :  ! Close std_out and ab_out and flush units.
    1016              :  ! Note that flush does not guarantee that the data is committed to disk.
    1017              :  ! This is rather annoying because we may end up with incomplete log files
    1018              :  ! that cannot be parsed by Abipy
    1019              :  ! For a possible approach based on fsync, see
    1020              :  ! https://gcc.gnu.org/onlinedocs/gcc-4.7.4/gfortran/FLUSH.html
    1021              : 
    1022           16 :  inquire(std_out, opened=testopen)
    1023           16 :  if (testopen) then
    1024              : #if defined HAVE_FC_FLUSH
    1025           16 :    call flush(std_out)
    1026              : #endif
    1027           16 :    close(std_out)
    1028              :  end if
    1029              : 
    1030           16 :  inquire(ab_out,opened=testopen)
    1031           16 :  if (testopen) then
    1032              : #if defined HAVE_FC_FLUSH
    1033           16 :    call flush(ab_out)
    1034              : #endif
    1035           16 :    close(ab_out)
    1036              :  end if
    1037              : 
    1038              :  ! exit_status=0 signals a deliberate, non-error termination (e.g. testkgrid
    1039              :  ! stopping after dumping k-grids for prtkpt/=0), not a real abort. MPI_ABORT
    1040              :  ! never returns, so it would discard exit_status: under some MPI/job-scheduler
    1041              :  ! combinations (e.g. Slurm srun) the step then gets reported as killed by
    1042              :  ! SIGKILL (retcode 137) instead of exiting 0. This code path is reached
    1043              :  ! identically by every rank, so exiting the calling process directly is safe.
    1044           16 :  if (present(exit_status)) then
    1045           14 :    if (exit_status == 0) call sys_exit(0)
    1046              :  end if
    1047              : 
    1048              : #ifdef HAVE_MPI
    1049            2 :  my_errorcode=MPI_ERR_UNKNOWN; if (PRESENT(mpierr)) my_errorcode=mpierr
    1050              : 
    1051            2 :  call MPI_ERROR_STRING(my_errorcode, mpi_msg_error, ilen, ierr2)
    1052              : 
    1053              :  !if (ilen>xmpi_msg_len) write(std_out,*)" WARNING: MPI message has been truncated!"
    1054              :  !if (ierr2/=MPI_SUCCESS) then
    1055              :  !  write(std_out,'(a,i0)')" WARNING: MPI_ERROR_STRING returned ierr2= ",ierr2
    1056              :  !else
    1057              :  !  write(std_out,'(2a)')" MPI_ERROR_STRING: ",TRIM(mpi_msg_error)
    1058              :  !end if
    1059              : 
    1060            2 :  call MPI_ABORT(my_comm, my_errorcode, ierr)
    1061              : #endif
    1062              : 
    1063            0 :  if (present(exit_status)) then
    1064              :    call sys_exit(exit_status)
    1065              :  else
    1066              :    call sys_exit(1)
    1067              :  end if
    1068              : 
    1069              : end subroutine xmpi_abort
    1070              : !!***
    1071              : 
    1072              : !----------------------------------------------------------------------
    1073              : 
    1074              : !!****f* m_xmpi/sys_exit
    1075              : !! NAME
    1076              : !! sys_exit
    1077              : !!
    1078              : !! FUNCTION
    1079              : !! Routine for clean exit of f90 code by one processor
    1080              : !!
    1081              : !! INPUTS
    1082              : !!   exit_status: return code.
    1083              : !!
    1084              : !! NOTES
    1085              : !!  By default, it uses "call exit(1)", that is not completely portable.
    1086              : !!
    1087              : !! SOURCE
    1088              : 
    1089         1786 : subroutine sys_exit(exit_status)
    1090              : 
    1091              : !Arguments ------------------------------------
    1092              : !scalars
    1093              :  integer,intent(in) :: exit_status
    1094              : ! **********************************************************************
    1095              : 
    1096              : #if defined FC_NAG
    1097              :  call exit(exit_status)
    1098              : #elif defined HAVE_FC_EXIT
    1099           14 :  call exit(exit_status)
    1100              : #else
    1101              :  ! stop with exit_status
    1102              :  ! MT 06-2013:stop function only accept parameters!
    1103              :  if (exit_status== 0) stop  "0"
    1104              :  if (exit_status== 1) stop  "1"
    1105              :  if (exit_status==-1) stop "-1"
    1106              : #endif
    1107              :  stop 1
    1108              : 
    1109              : end subroutine sys_exit
    1110              : !!***
    1111              : 
    1112              : !----------------------------------------------------------------------
    1113              : 
    1114              : !!****f* m_xmpi/xmpi_show_info
    1115              : !! NAME
    1116              : !!  xmpi_show_info
    1117              : !!
    1118              : !! FUNCTION
    1119              : !!  Printout of the most important variables stored in this module (useful for debugging).
    1120              : !!
    1121              : !! INPUTS
    1122              : !!  unt=Unit number for formatted output.
    1123              : !!
    1124              : !! SOURCE
    1125              : 
    1126         1150 : subroutine xmpi_show_info(unit)
    1127              : 
    1128              : !Arguments-------------------------
    1129              :  integer,optional,intent(in) :: unit
    1130              : 
    1131              : !Local variables-------------------
    1132              :  integer :: my_unt
    1133              : ! *************************************************************************
    1134              : 
    1135         1150 :  my_unt = std_out; if (PRESENT(unit)) my_unt=unit
    1136              : 
    1137              : #ifdef HAVE_MPI1
    1138              :   write(my_unt,*)" ==== Using MPI-1 specifications ==== "
    1139              : #endif
    1140              : #ifdef HAVE_MPI2
    1141         1150 :   write(my_unt,*)" ==== Using MPI-2 specifications ==== "
    1142              : #endif
    1143              : 
    1144              : #ifdef HAVE_MPI_IO
    1145         1150 :   write(my_unt,*)" MPI-IO support is ON"
    1146              : #else
    1147              :   write(my_unt,*)" MPI-IO support is OFF"
    1148              : #endif
    1149              : 
    1150              : #ifdef HAVE_MPI
    1151         1150 :  write(my_unt,*)" xmpi_tag_ub ................ ",xmpi_tag_ub
    1152         1150 :  write(my_unt,*)" xmpi_bsize_ch .............. ",xmpi_bsize_ch
    1153         1150 :  write(my_unt,*)" xmpi_bsize_int ............. ",xmpi_bsize_int
    1154         1150 :  write(my_unt,*)" xmpi_bsize_sp .............. ",xmpi_bsize_sp
    1155         1150 :  write(my_unt,*)" xmpi_bsize_dp .............. ",xmpi_bsize_dp
    1156         1150 :  write(my_unt,*)" xmpi_bsize_spc ............. ",xmpi_bsize_spc
    1157         1150 :  write(my_unt,*)" xmpi_bsize_dpc ............. ",xmpi_bsize_dpc
    1158         1150 :  write(my_unt,*)" xmpio_bsize_frm ............ ",xmpio_bsize_frm
    1159         1150 :  write(my_unt,*)" xmpi_address_kind .......... ",xmpi_address_kind
    1160         1150 :  write(my_unt,*)" xmpi_offset_kind ........... ",xmpi_offset_kind
    1161         1150 :  write(my_unt,*)" MPI_WTICK .................. ",MPI_WTICK()
    1162              : #endif
    1163              : 
    1164         1150 : end subroutine xmpi_show_info
    1165              : !!***
    1166              : 
    1167              : !----------------------------------------------------------------------
    1168              : 
    1169              : !!****f* m_xmpi/xmpi_comm_rank
    1170              : !! NAME
    1171              : !!  xmpi_comm_rank
    1172              : !!
    1173              : !! FUNCTION
    1174              : !!  Hides MPI_COMM_RANK from MPI library.
    1175              : !!
    1176              : !! INPUTS
    1177              : !!  comm=MPI communicator.
    1178              : !!
    1179              : !! OUTPUT
    1180              : !!  xmpi_comm_rank=The rank of the node inside comm
    1181              : !!
    1182              : !! SOURCE
    1183              : 
    1184     10381177 : integer function xmpi_comm_rank(comm)
    1185              : 
    1186              : !Arguments-------------------------
    1187              :  integer,intent(in) :: comm
    1188              : 
    1189              : !Local variables-------------------
    1190              :  integer :: mpierr
    1191              : ! *************************************************************************
    1192              : 
    1193     10381177 :  mpierr=0
    1194              : #ifdef HAVE_MPI
    1195     10381177 :  xmpi_comm_rank=-1  ! Return non-sense value if the proc does not belong to the comm
    1196     10381177 :  if (comm/=xmpi_comm_null) then
    1197     10381177 :    call MPI_COMM_RANK(comm,xmpi_comm_rank,mpierr)
    1198              :  end if
    1199              : #else
    1200              :  xmpi_comm_rank=0
    1201              : #endif
    1202              : 
    1203     10381177 : end function xmpi_comm_rank
    1204              : !!***
    1205              : 
    1206              : !----------------------------------------------------------------------
    1207              : 
    1208              : !!****f* m_xmpi/xmpi_comm_size
    1209              : !! NAME
    1210              : !!  xmpi_comm_size
    1211              : !!
    1212              : !! FUNCTION
    1213              : !!  Hides MPI_COMM_SIZE from MPI library.
    1214              : !!
    1215              : !! INPUTS
    1216              : !!  comm=MPI communicator.
    1217              : !!
    1218              : !! OUTPUT
    1219              : !!  xmpi_comm_size=The number of processors inside comm. Return 0 if comm = xmpi_comm_null
    1220              : !!
    1221              : !! SOURCE
    1222              : 
    1223    213868291 : integer function xmpi_comm_size(comm)
    1224              : 
    1225              : !Arguments-------------------------
    1226              :  integer,intent(in) :: comm
    1227              : 
    1228              : !Local variables-------------------------------
    1229              :  integer :: mpierr
    1230              : ! *************************************************************************
    1231              : 
    1232    213868291 :  mpierr=0; xmpi_comm_size=1
    1233              : #ifdef HAVE_MPI
    1234    213868291 :  xmpi_comm_size = 0
    1235    213868291 :  if (comm /= xmpi_comm_null) call MPI_COMM_SIZE(comm,xmpi_comm_size,mpierr)
    1236              : #endif
    1237              : 
    1238    213868291 : end function xmpi_comm_size
    1239              : !!***
    1240              : 
    1241              : !----------------------------------------------------------------------
    1242              : 
    1243              : !!****f* m_xmpi/xmpi_comm_free_0D
    1244              : !! NAME
    1245              : !!  xmpi_comm_free_0D
    1246              : !!
    1247              : !! FUNCTION
    1248              : !!  Hides MPI_COMM_FREE from MPI library.
    1249              : !!  Does not abort MPI in case of an invalid communicator
    1250              : !!
    1251              : !! INPUTS
    1252              : !!  comm=MPI communicator.
    1253              : !!
    1254              : !! SOURCE
    1255              : 
    1256        79673 : subroutine xmpi_comm_free_0D(comm)
    1257              : 
    1258              : !Arguments-------------------------
    1259              :  integer,intent(inout) :: comm
    1260              : 
    1261              : !Local variables-------------------------------
    1262              : #ifdef HAVE_MPI
    1263              :  integer :: comm_world,err_handler_dum,err_handler_sav,ierr,mpierr,mpierr_class
    1264              : ! *************************************************************************
    1265              : 
    1266        79673 :  if (comm/=xmpi_comm_null.and.comm/=xmpi_world.and.comm/=xmpi_comm_self) then
    1267              : 
    1268         7820 :    comm_world=xmpi_world ! Needed to bypass a bug in some OMPI implementations (intent(inout))
    1269         7820 :    call xmpi_comm_set_errhandler(comm_world,MPI_ERRORS_RETURN,err_handler_sav,ierr)
    1270         7820 :    call MPI_COMM_FREE(comm,mpierr)
    1271         7820 :    call xmpi_comm_set_errhandler(comm_world,err_handler_sav,err_handler_dum,ierr)
    1272              : 
    1273         7820 :    if (mpierr/=MPI_SUCCESS) then
    1274          396 :      call MPI_ERROR_CLASS(mpierr,mpierr_class,ierr)
    1275          396 :      if (mpierr_class/=MPI_ERR_COMM) then
    1276            0 :        write(std_out,*)" WARNING: MPI_COMM_FREE returned ierr= ",mpierr
    1277              :      end if
    1278              :    end if
    1279              : 
    1280              :  end if
    1281              : 
    1282              : #else
    1283              :  if (.false.) write(std_out,*) comm
    1284              : #endif
    1285              : 
    1286        79673 : end subroutine xmpi_comm_free_0D
    1287              : !!***
    1288              : 
    1289              : !----------------------------------------------------------------------
    1290              : 
    1291              : !!****f* m_xmpi/xmpi_comm_free_1D
    1292              : !! NAME
    1293              : !!  xmpi_comm_free_1D
    1294              : !!
    1295              : !! FUNCTION
    1296              : !!  Hides MPI_COMM_FREE from MPI library. Target 1D arrays
    1297              : !!  Does not abort MPI in case of an invalid communicator
    1298              : !!
    1299              : !! INPUTS
    1300              : !!  comms(:)=MPI communicators
    1301              : !!
    1302              : !! SOURCE
    1303              : 
    1304           24 : subroutine xmpi_comm_free_1D(comms)
    1305              : 
    1306              : !Arguments-------------------------
    1307              :  integer,intent(inout) :: comms(:)
    1308              : 
    1309              : !Local variables-------------------------------
    1310              : !scalars
    1311              : #ifdef HAVE_MPI
    1312              :  integer :: comm_world,err_handler_dum,err_handler_sav,ii,mpierr
    1313              : ! *************************************************************************
    1314              : 
    1315           24 :  comm_world=xmpi_world ! Needed to bypass a bug in some OMPI implementations (intent(inout))
    1316           24 :  call xmpi_comm_set_errhandler(comm_world,MPI_ERRORS_RETURN,err_handler_sav,mpierr)
    1317              : 
    1318           72 :  do ii=LBOUND(comms,DIM=1),UBOUND(comms,DIM=1)
    1319           48 :    if (comms(ii)/=xmpi_comm_null.and.comms(ii)/=xmpi_world.and.comms(ii)/=xmpi_comm_self) then
    1320           24 :      call MPI_COMM_FREE(comms(ii),mpierr)
    1321              :    end if
    1322              :  end do
    1323              : 
    1324           24 :  call xmpi_comm_set_errhandler(comm_world,err_handler_sav,err_handler_dum,mpierr)
    1325              : 
    1326              : #else
    1327              :  if (.false.) write(std_out,*) comms(1)
    1328              : #endif
    1329              : 
    1330           24 : end subroutine xmpi_comm_free_1D
    1331              : !!***
    1332              : 
    1333              : !----------------------------------------------------------------------
    1334              : 
    1335              : !!****f* m_xmpi/xmpi_comm_free_2D
    1336              : !! NAME
    1337              : !!  xmpi_comm_free_2D
    1338              : !!
    1339              : !! FUNCTION
    1340              : !!  Hides MPI_COMM_FREE from MPI library. Target 2D arrays
    1341              : !!  Does not abort MPI in case of an invalid communicator
    1342              : !!
    1343              : !! INPUTS
    1344              : !!  comms=MPI communicator.
    1345              : !!
    1346              : !! SOURCE
    1347              : 
    1348            0 : subroutine xmpi_comm_free_2D(comms)
    1349              : 
    1350              : !Arguments-------------------------
    1351              :  integer,intent(inout) :: comms(:,:)
    1352              : 
    1353              : !Local variables-------------------------------
    1354              : !scalars
    1355              : #ifdef HAVE_MPI
    1356              :  integer :: comm_world,err_handler_dum,err_handler_sav,ii,jj,mpierr
    1357              : ! *************************************************************************
    1358              : 
    1359            0 :  comm_world=xmpi_world ! Needed to bypass a bug in some OMPI implementations (intent(inout))
    1360            0 :  call xmpi_comm_set_errhandler(comm_world,MPI_ERRORS_RETURN,err_handler_sav,mpierr)
    1361              : 
    1362            0 :  do jj=LBOUND(comms,DIM=2),UBOUND(comms,DIM=2)
    1363            0 :    do ii=LBOUND(comms,DIM=1),UBOUND(comms,DIM=1)
    1364            0 :      if (comms(ii,jj)/=xmpi_comm_null.and.comms(ii,jj)/=xmpi_world.and. comms(ii,jj)/=xmpi_comm_self) then
    1365            0 :        call MPI_COMM_FREE(comms(ii,jj),mpierr)
    1366              :      end if
    1367              :    end do
    1368              :  end do
    1369              : 
    1370            0 :  call xmpi_comm_set_errhandler(comm_world,err_handler_sav,err_handler_dum,mpierr)
    1371              : 
    1372              : #else
    1373              :  if (.false.) write(std_out,*) comms(1,1)
    1374              : #endif
    1375              : 
    1376            0 : end subroutine xmpi_comm_free_2D
    1377              : !!***
    1378              : 
    1379              : !----------------------------------------------------------------------
    1380              : 
    1381              : !!****f* m_xmpi/xmpi_comm_free_3D
    1382              : !! NAME
    1383              : !!  xmpi_comm_free_3D
    1384              : !!
    1385              : !! FUNCTION
    1386              : !!  Hides MPI_COMM_FREE from MPI library. Target 3D arrays
    1387              : !!  Does not abort MPI in case of an invalid communicator
    1388              : !!
    1389              : !! INPUTS
    1390              : !!  comms=MPI communicator.
    1391              : !!
    1392              : !! SOURCE
    1393              : 
    1394            0 : subroutine xmpi_comm_free_3D(comms)
    1395              : 
    1396              : !Arguments-------------------------
    1397              :  integer,intent(inout) :: comms(:,:,:)
    1398              : 
    1399              : !Local variables-------------------------------
    1400              : !scalars
    1401              : #ifdef HAVE_MPI
    1402              :  integer :: comm_world,err_handler_dum,err_handler_sav,ii,jj,kk,mpierr
    1403              : ! *************************************************************************
    1404              : 
    1405            0 :  comm_world=xmpi_world ! Needed to bypass a bug in some OMPI implementations (intent(inout))
    1406            0 :  call xmpi_comm_set_errhandler(comm_world,MPI_ERRORS_RETURN,err_handler_sav,mpierr)
    1407              : 
    1408            0 :  do kk=LBOUND(comms,DIM=3),UBOUND(comms,DIM=3)
    1409            0 :    do jj=LBOUND(comms,DIM=2),UBOUND(comms,DIM=2)
    1410            0 :      do ii=LBOUND(comms,DIM=1),UBOUND(comms,DIM=1)
    1411            0 :        if (comms(ii,jj,kk)/=xmpi_comm_null.and.comms(ii,jj,kk)/=xmpi_world.and. comms(ii,jj,kk)/=xmpi_comm_self) then
    1412            0 :          call MPI_COMM_FREE(comms(ii,jj,kk),mpierr)
    1413              :        end if
    1414              :      end do
    1415              :    end do
    1416              :  end do
    1417              : 
    1418            0 :  call xmpi_comm_set_errhandler(comm_world,err_handler_sav,err_handler_dum,mpierr)
    1419              : 
    1420              : #else
    1421              :  if (.false.) write(std_out,*) comms(1,1,1)
    1422              : #endif
    1423              : 
    1424            0 : end subroutine xmpi_comm_free_3D
    1425              : !!***
    1426              : 
    1427              : !----------------------------------------------------------------------
    1428              : 
    1429              : !!****f* m_xmpi/xmpi_group_free
    1430              : !! NAME
    1431              : !!  xmpi_group_free
    1432              : !!
    1433              : !! FUNCTION
    1434              : !!  Hides MPI_GROUP_FREE from MPI library.
    1435              : !!  Does not abort MPI in case of an invalid group
    1436              : !!
    1437              : !! INPUTS
    1438              : !!  spaceGroup=MPI group
    1439              : !!
    1440              : !! SOURCE
    1441              : 
    1442         1166 : subroutine xmpi_group_free(spaceGroup)
    1443              : 
    1444              : !Arguments-------------------------
    1445              :  integer,intent(inout) :: spaceGroup
    1446              : 
    1447              : !Local variables-------------------------------
    1448              : !scalars
    1449              : #ifdef HAVE_MPI
    1450              :  integer :: comm_world,err_handler_dum,err_handler_sav,ierr,mpierr,mpierr_class
    1451              : ! *************************************************************************
    1452              : 
    1453         1166 :  if (spaceGroup/=xmpi_group_null) then
    1454              : 
    1455         1166 :    comm_world=xmpi_world ! Needed to bypass a bug in some OMPI implementations (intent(inout))
    1456         1166 :    call xmpi_comm_set_errhandler(comm_world,MPI_ERRORS_RETURN,err_handler_sav,ierr)
    1457         1166 :    call MPI_GROUP_FREE(spaceGroup,mpierr)
    1458         1166 :    call xmpi_comm_set_errhandler(comm_world,err_handler_sav,err_handler_dum,ierr)
    1459              : 
    1460         1166 :    if (mpierr/=MPI_SUCCESS) then
    1461            0 :      call MPI_ERROR_CLASS(mpierr,mpierr_class,ierr)
    1462            0 :      if (mpierr_class/=MPI_ERR_GROUP) write(std_out,*)" WARNING: MPI_GROUP_FREE returned ierr= ",mpierr
    1463              :    end if
    1464              : 
    1465              :  end if
    1466              : 
    1467              : #else
    1468              :  if (.false.) write(std_out,*) spaceGroup
    1469              : #endif
    1470              : 
    1471         1166 : end subroutine xmpi_group_free
    1472              : !!***
    1473              : 
    1474              : !----------------------------------------------------------------------
    1475              : 
    1476              : !!****f* m_xmpi/xmpi_group_incl
    1477              : !! NAME
    1478              : !!  xmpi_group_incl
    1479              : !!
    1480              : !! FUNCTION
    1481              : !!  Hides MPI_GROUP_INCL from MPI library.
    1482              : !!
    1483              : !! INPUTS
    1484              : !!  group=input group
    1485              : !!  nrank=number of elements in array ranks (size of newgroup)
    1486              : !!  ranks=ranks of processes in group to appear in newgroup
    1487              : !!
    1488              : !! OUTPUT
    1489              : !!  newgroup= new group derived from above, in the order defined by ranks
    1490              : !!
    1491              : !! SOURCE
    1492              : 
    1493            0 : subroutine xmpi_group_incl(group,nranks,ranks,newgroup,mpierr)
    1494              : 
    1495              : !Arguments-------------------------
    1496              : !scalars
    1497              :  integer,intent(in) :: group,nranks
    1498              :  integer,intent(out) :: mpierr
    1499              :  integer,intent(inout) :: newgroup
    1500              : !arrays
    1501              :  integer,intent(in) :: ranks(nranks)
    1502              : ! *************************************************************************
    1503              : 
    1504            0 :  mpierr=0 ; newgroup=xmpi_group_null
    1505              : #ifdef HAVE_MPI
    1506            0 :  if (group/=xmpi_group_null) then
    1507            0 :    call MPI_GROUP_INCL(group,nranks,ranks,newgroup,mpierr)
    1508              :  end if
    1509              : #endif
    1510              : 
    1511            0 : end subroutine xmpi_group_incl
    1512              : !!***
    1513              : 
    1514              : !----------------------------------------------------------------------
    1515              : 
    1516              : !!****f* m_xmpi/xmpi_comm_create
    1517              : !! NAME
    1518              : !!  xmpi_comm_create
    1519              : !!
    1520              : !! FUNCTION
    1521              : !!  Hides MPI_COMM_CREATE from MPI library.
    1522              : !!
    1523              : !! INPUTS
    1524              : !!  comm=communicator
    1525              : !!  group=group, which is a subset of the group of comm
    1526              : !!
    1527              : !! OUTPUT
    1528              : !!  newcomm=new communicator
    1529              : !!
    1530              : !! SOURCE
    1531              : 
    1532            0 : subroutine xmpi_comm_create(comm,group,newcomm,mpierr)
    1533              : 
    1534              : !Arguments-------------------------
    1535              : !scalars
    1536              :  integer,intent(in) :: comm,group
    1537              :  integer,intent(out) :: mpierr
    1538              :  integer,intent(inout) :: newcomm
    1539              : ! *************************************************************************
    1540              : 
    1541            0 :  mpierr=0
    1542              : #ifdef HAVE_MPI
    1543            0 :  if (group/=xmpi_group_null) then
    1544            0 :    call MPI_comm_create(comm,group,newcomm,mpierr)
    1545              :  else
    1546            0 :    newcomm=xmpi_comm_null
    1547              :  end if
    1548              : #else
    1549              :   newcomm=xmpi_comm_self
    1550              : #endif
    1551              : 
    1552            0 : end subroutine xmpi_comm_create
    1553              : !!***
    1554              : 
    1555              : !----------------------------------------------------------------------
    1556              : 
    1557              : !!****f* m_xmpi/xmpi_subcomm
    1558              : !! NAME
    1559              : !!  xmpi_subcomm
    1560              : !!
    1561              : !! FUNCTION
    1562              : !!  Return a sub-communicator from an input communicator and a given proc. ranks set.
    1563              : !!  (hides subgroup creation/destruction)
    1564              : !!
    1565              : !! INPUTS
    1566              : !!  comm=input communicator
    1567              : !!  nrank=number of elements in array ranks (size of subcomm)
    1568              : !!  ranks=ranks of processes in group to appear in subcomm
    1569              : !!
    1570              : !! OUTPUT
    1571              : !!  [my_rank_in_group]=optional: my rank in the group of new sub-communicator
    1572              : !!  xmpi_subcomm=new (sub-)communicator
    1573              : !!
    1574              : !! SOURCE
    1575              : 
    1576          148 : integer function xmpi_subcomm(comm,nranks,ranks,my_rank_in_group)
    1577              : 
    1578              : !Arguments-------------------------
    1579              : !scalars
    1580              :  integer,intent(in) :: comm,nranks
    1581              :  integer,intent(out),optional :: my_rank_in_group
    1582              : !arrays
    1583              :  integer,intent(in) :: ranks(nranks)
    1584              : 
    1585              : !Local variables-------------------------------
    1586              : #ifdef HAVE_MPI
    1587              :  integer :: group,ierr,subgroup
    1588              : #endif
    1589              : ! *************************************************************************
    1590              : 
    1591          148 :  xmpi_subcomm=xmpi_comm_null
    1592          148 :  if (present(my_rank_in_group)) my_rank_in_group=xmpi_undefined
    1593              : 
    1594              : #ifdef HAVE_MPI
    1595          148 :  if (comm/=xmpi_comm_null.and.nranks>=0) then
    1596          148 :    call MPI_COMM_GROUP(comm,group,ierr)
    1597          148 :    call MPI_GROUP_INCL(group,nranks,ranks,subgroup,ierr)
    1598          148 :    call MPI_COMM_CREATE(comm,subgroup,xmpi_subcomm,ierr)
    1599          148 :    if ( nranks == 0 )xmpi_subcomm=xmpi_comm_self
    1600          148 :    if (present(my_rank_in_group)) then
    1601           24 :      call MPI_Group_rank(subgroup,my_rank_in_group,ierr)
    1602              :    end if
    1603          148 :    call MPI_GROUP_FREE(subgroup,ierr)
    1604          148 :    call MPI_GROUP_FREE(group,ierr)
    1605              :  end if
    1606              : #else
    1607              :  if (nranks>0) then
    1608              :    if (ranks(1)==0) then
    1609              :      xmpi_subcomm=xmpi_comm_self
    1610              :      if (present(my_rank_in_group)) my_rank_in_group=0
    1611              :    end if
    1612              :  end if
    1613              : #endif
    1614              : 
    1615          148 : end function xmpi_subcomm
    1616              : !!***
    1617              : 
    1618              : !----------------------------------------------------------------------
    1619              : 
    1620              : !!****f* m_xmpi/xmpi_comm_multiple
    1621              : !! NAME
    1622              : !!  xmpi_comm_multiple
    1623              : !!
    1624              : !! FUNCTION
    1625              : !!  Given an input communicator `input_comm`, create a new communicator
    1626              : !!  with number of procs multiple of a certain number `ntasks`.
    1627              : !!  Use all procs if ntasks >= input_nprocs.
    1628              : !!
    1629              : !! INPUTS
    1630              : !!  ntasks=Number of tasks.
    1631              : !!  comm=input communicator
    1632              : !!
    1633              : !! OUTPUT
    1634              : !!  idle_proc=True if this proc is idle. In this case, output_comm contains all the idle procs.
    1635              : !!  output_comm=Output communicator
    1636              : !!
    1637              : !! SOURCE
    1638              : 
    1639            0 : subroutine xmpi_comm_multiple_of(ntasks, input_comm, idle_proc, output_comm)
    1640              : 
    1641              : !Arguments-------------------------
    1642              : !scalars
    1643              :  integer,intent(in) :: ntasks, input_comm
    1644              :  integer,intent(out) :: output_comm
    1645              :  logical,intent(out) :: idle_proc
    1646              : 
    1647              : !Local variables-------------------------------
    1648              :  integer :: color, my_rank, ierr, input_nproc
    1649              : ! *************************************************************************
    1650              : 
    1651            0 :  my_rank = xmpi_comm_rank(input_comm)
    1652            0 :  input_nproc = xmpi_comm_size(input_comm)
    1653              : 
    1654            0 :  if (input_nproc <= ntasks) then
    1655              :    ! Use all procs in input comm.
    1656            0 :    idle_proc = .False.; output_comm = input_comm
    1657              : #ifdef HAVE_MPI
    1658            0 :    call MPI_Comm_dup(input_comm, output_comm, ierr)
    1659              : #endif
    1660              :  else
    1661            0 :    color = merge(0, 1, my_rank + 1 <= (ntasks / input_nproc) * input_nproc)
    1662            0 :    idle_proc = color == 1
    1663            0 :    call xmpi_comm_split(input_comm, color, my_rank, output_comm, ierr)
    1664              :  end if
    1665              : 
    1666            0 : end subroutine xmpi_comm_multiple_of
    1667              : !!***
    1668              : 
    1669              : !----------------------------------------------------------------------
    1670              : 
    1671              : !!****f* m_xmpi/xmpi_comm_dup
    1672              : !! NAME
    1673              : !!  xmpi_comm_dup
    1674              : !!
    1675              : !! FUNCTION
    1676              : !!  Hides MPI_COMM_DUP from MPI library.
    1677              : !!
    1678              : !! INPUTS
    1679              : !!  in_comm=input MPI communicator.
    1680              : !!
    1681              : !! OUTPUT
    1682              : !!  out_comm=Output MPI communicator.
    1683              : !!  mpierr=error code returned
    1684              : !!
    1685              : !! SOURCE
    1686              : 
    1687            0 : subroutine xmpi_comm_dup(in_comm, out_comm, mpierr)
    1688              : 
    1689              : !Arguments-------------------------
    1690              :  integer,intent(in) :: in_comm
    1691              :  integer,intent(out) :: out_comm, mpierr
    1692              : 
    1693              : !----------------------------------------------------------------------
    1694              : 
    1695              : #ifdef HAVE_MPI
    1696            0 :   call MPI_Comm_dup(in_comm, out_comm, mpierr)
    1697              : #else
    1698              :   out_comm = in_comm
    1699              : #endif
    1700              : 
    1701            0 : end subroutine xmpi_comm_dup
    1702              : !!***
    1703              : 
    1704              : !----------------------------------------------------------------------
    1705              : 
    1706              : !!****f* m_xmpi/xmpi_comm_group
    1707              : !! NAME
    1708              : !!  xmpi_comm_group
    1709              : !!
    1710              : !! FUNCTION
    1711              : !!  Hides MPI_COMM_GROUP from MPI library.
    1712              : !!
    1713              : !! INPUTS
    1714              : !!  comm=MPI communicator.
    1715              : !!
    1716              : !! OUTPUT
    1717              : !!  spaceGroup=The group associated to comm.
    1718              : !!  mpierr=error code returned
    1719              : !!
    1720              : !! SOURCE
    1721              : 
    1722         1166 : subroutine xmpi_comm_group(comm,spaceGroup,mpierr)
    1723              : 
    1724              : !Arguments-------------------------
    1725              :  integer,intent(in) :: comm
    1726              :  integer,intent(out) :: mpierr,spaceGroup
    1727              : ! *************************************************************************
    1728              : 
    1729         1166 :  mpierr=0; spaceGroup=xmpi_group_null
    1730              : #ifdef HAVE_MPI
    1731         1166 :  if (comm/=xmpi_comm_null) then
    1732         1166 :    call MPI_COMM_GROUP(comm,spaceGroup,mpierr)
    1733              :  end if
    1734              : #endif
    1735              : 
    1736         1166 : end subroutine xmpi_comm_group
    1737              : !!***
    1738              : 
    1739              : !----------------------------------------------------------------------
    1740              : 
    1741              : !!****f* m_xmpi/xmpi_comm_split
    1742              : !! NAME
    1743              : !!  xmpi_comm_split
    1744              : !!
    1745              : !! FUNCTION
    1746              : !!  Hides MPI_COMM_SPLIT from MPI library.
    1747              : !!
    1748              : !! INPUTS
    1749              : !!  input_comm=Input MPI communicator (to be splitted)
    1750              : !!  color=Control of subset assignment (nonnegative integer).
    1751              : !!        Processes with the same color are in the same new communicator
    1752              : !!  key=Control of rank assigment (integer)
    1753              : !!
    1754              : !! OUTPUT
    1755              : !!  mpierr=error code returned
    1756              : !!  output_comm=new splitted communicator
    1757              : !!
    1758              : !! SOURCE
    1759              : 
    1760          713 : subroutine xmpi_comm_split(input_comm, color, key, output_comm, mpierr)
    1761              : 
    1762              : !Arguments-------------------------
    1763              :  integer,intent(in) :: color,input_comm,key
    1764              :  integer,intent(out) :: mpierr,output_comm
    1765              : ! *************************************************************************
    1766              : 
    1767          713 :  mpierr=0; output_comm=input_comm
    1768              : #ifdef HAVE_MPI
    1769          713 :  if (input_comm/=xmpi_comm_null.and.input_comm/=xmpi_comm_self) then
    1770          702 :    call MPI_COMM_SPLIT(input_comm,color,key,output_comm,mpierr)
    1771              :  end if
    1772              : #endif
    1773              : 
    1774          713 : end subroutine xmpi_comm_split
    1775              : !!***
    1776              : 
    1777              : !----------------------------------------------------------------------
    1778              : 
    1779              : !!****f* m_xmpi/xmpi_group_translate_ranks
    1780              : !! NAME
    1781              : !!  xmpi_group_translate_ranks
    1782              : !!
    1783              : !! FUNCTION
    1784              : !!  Hides MPI_GROUP_TRANSLATE_RANKS from MPI library.
    1785              : !!
    1786              : !! INPUTS
    1787              : !!  nrank=number of ranks in ranks1 and ranks2 arrays
    1788              : !!  ranks1(nrank)=array of zero or more valid ranks in group1
    1789              : !!  spaceGroup1=group1
    1790              : !!  spaceGroup2=group2
    1791              : !!
    1792              : !! OUTPUT
    1793              : !!  mpierr=error code returned
    1794              : !!  ranks2(nrank)=array of corresponding ranks in group2,
    1795              : !!                xmpi_undefined when no correspondence exists
    1796              : !!
    1797              : !! SOURCE
    1798              : 
    1799          583 : subroutine xmpi_group_translate_ranks(spaceGroup1,nrank,ranks1,&
    1800          583 : &                                     spaceGroup2,ranks2,mpierr)
    1801              : 
    1802              : !Arguments-------------------------
    1803              : !scalars
    1804              :  integer,intent(in) :: nrank,spaceGroup1,spaceGroup2
    1805              :  integer,intent(out) :: mpierr
    1806              : !arrays
    1807              :  integer,intent(in) :: ranks1(nrank)
    1808              :  integer,intent(out) :: ranks2(nrank)
    1809              : ! *************************************************************************
    1810              : 
    1811         1173 :  mpierr=0; ranks2(:)=xmpi_undefined
    1812              : #ifdef HAVE_MPI
    1813          583 :  if (spaceGroup1/=xmpi_group_null.and.spaceGroup2/=xmpi_group_null) then
    1814          583 :    call MPI_GROUP_TRANSLATE_RANKS(spaceGroup1,nrank,ranks1, spaceGroup2,ranks2,mpierr)
    1815              :  end if
    1816              : #else
    1817              :  ranks2(1)=0
    1818              : #endif
    1819              : 
    1820          583 : end subroutine xmpi_group_translate_ranks
    1821              : !!***
    1822              : 
    1823              : !----------------------------------------------------------------------
    1824              : 
    1825              : !!****f* m_xmpi/xmpi_comm_translate_ranks
    1826              : !! NAME
    1827              : !!  xmpi_comm_translate_ranks
    1828              : !!
    1829              : !! FUNCTION
    1830              : !!  Helper function that translate the ranks from a communicator to another one.
    1831              : !!  Wraps xmpi_group_translate_ranks but provides a more user-friendly interface
    1832              : !!
    1833              : !! INPUTS
    1834              : !!  from_comm=MPI communicator where from_ranks are defined.
    1835              : !!  nrank=number of ranks in from_ranks and to_ranks arrays
    1836              : !!  from_ranks(nrank)=array of zero or more valid ranks in from_comm
    1837              : !!
    1838              : !! OUTPUT
    1839              : !!  to_ranks(nrank)=array of corresponding ranks in to_comm
    1840              : !!                xmpi_undefined when no correspondence exists
    1841              : !!
    1842              : !! SOURCE
    1843              : 
    1844           75 : subroutine xmpi_comm_translate_ranks(from_comm, nrank, from_ranks, to_comm, to_ranks)
    1845              : 
    1846              : !Arguments-------------------------
    1847              : !scalars
    1848              :  integer,intent(in) :: nrank,from_comm,to_comm
    1849              : !arrays
    1850              :  integer,intent(in) :: from_ranks(nrank)
    1851              :  integer,intent(out) :: to_ranks(nrank)
    1852              : 
    1853              : !Local variables-------------------------------
    1854              : !scalars
    1855              :  integer :: ierr,from_group,to_group
    1856              : ! *************************************************************************
    1857              : 
    1858              :  ! Get the groups
    1859           25 :  call xmpi_comm_group(from_comm,from_group,ierr)
    1860           25 :  call xmpi_comm_group(to_comm,to_group,ierr)
    1861              : 
    1862           25 :  call xmpi_group_translate_ranks(from_group,nrank,from_ranks,to_group,to_ranks,ierr)
    1863              : 
    1864              :  ! Release the groups
    1865           25 :  call xmpi_group_free(from_group)
    1866           25 :  call xmpi_group_free(to_group)
    1867              : 
    1868           25 : end subroutine xmpi_comm_translate_ranks
    1869              : !!***
    1870              : 
    1871              : !----------------------------------------------------------------------
    1872              : 
    1873              : !!****f* m_xmpi/xmpi_comm_translate_rank
    1874              : !! NAME
    1875              : !!  xmpi_comm_translate_rank
    1876              : !!
    1877              : !! FUNCTION
    1878              : !!  Helper function to translate a single rank `from_rank` in communicator `from_rank` to
    1879              : !!  the rank in communicator `to_comm`.
    1880              : 
    1881            0 : integer function xmpi_comm_translate_rank(from_comm, from_rank, to_comm) result(to_rank)
    1882              : 
    1883              : !Arguments ------------------------------------
    1884              :  integer,intent(in) :: from_comm, from_rank, to_comm
    1885              : 
    1886              : !Local variables-------------------------------
    1887              :  integer :: from_ranks(1), to_ranks(1)
    1888              : ! *************************************************************************
    1889              : 
    1890            0 :  from_ranks(1) = from_rank
    1891            0 :  call xmpi_comm_translate_ranks(from_comm, 1, from_ranks, to_comm, to_ranks)
    1892            0 :  to_rank = to_ranks(1)
    1893              : 
    1894            0 : end function xmpi_comm_translate_rank
    1895              : !!***
    1896              : 
    1897              : !----------------------------------------------------------------------
    1898              : 
    1899              : !!****f* m_xmpi/xmpi_barrier
    1900              : !! NAME
    1901              : !!  xmpi_barrier
    1902              : !!
    1903              : !! FUNCTION
    1904              : !!  Hides MPI_BARRIER from MPI library.
    1905              : !!
    1906              : !! INPUTS
    1907              : !!  comm=MPI communicator
    1908              : !!
    1909              : !! SOURCE
    1910              : 
    1911       104271 : subroutine xmpi_barrier(comm)
    1912              : 
    1913              : !Arguments-------------------------
    1914              :  integer,intent(in) :: comm
    1915              : 
    1916              : !Local variables-------------------
    1917              :  integer :: ier
    1918              : #ifdef HAVE_MPI
    1919              :  integer :: nprocs
    1920              : #endif
    1921              : ! *************************************************************************
    1922              : 
    1923       104271 :  ier = 0
    1924              : #ifdef HAVE_MPI
    1925       104271 :  if (comm/=xmpi_comm_null) then
    1926       104271 :    call MPI_COMM_SIZE(comm,nprocs,ier)
    1927       104271 :    if(nprocs>1) call MPI_BARRIER(comm,ier)
    1928              :  end if
    1929              : #endif
    1930              : 
    1931       104271 : end subroutine xmpi_barrier
    1932              : !!***
    1933              : 
    1934              : !----------------------------------------------------------------------
    1935              : 
    1936              : !!****f* m_xmpi/xmpi_name
    1937              : !! NAME
    1938              : !!  xmpi_name
    1939              : !!
    1940              : !! FUNCTION
    1941              : !!  Returns the name of the processor
    1942              : !!  Hides MPI_GET_PROCESSOR_NAME from MPI library.
    1943              : !!
    1944              : !! For the MPI standard:
    1945              : !!    The name returned should identify a particular piece of hardware; the exact format is implementation defined.
    1946              : !!    This name may or may not be the same as might be returned by gethostname, uname, or sysinfo.
    1947              : !!
    1948              : !! SOURCE
    1949              : 
    1950            0 : subroutine xmpi_name(name_ch, ierr)
    1951              : 
    1952              : !Arguments-------------------------
    1953              :  character(xmpi_max_processor_name),intent(out) :: name_ch
    1954              :  integer,intent(out) ::  ierr
    1955              : 
    1956              : !Local variables-------------------
    1957              :  integer :: len
    1958              : ! *************************************************************************
    1959              : 
    1960              :  ! Get the name of this processor (usually the hostname)
    1961            0 :  ierr = 0
    1962              : #ifdef HAVE_MPI
    1963            0 :  call MPI_GET_PROCESSOR_NAME(name_ch, len, ierr)
    1964            0 :  name_ch = trim(name_ch(1:len))
    1965              : 
    1966              : #else
    1967              :  name_ch = '0'
    1968              : #endif
    1969              : 
    1970            0 : end subroutine xmpi_name
    1971              : !!***
    1972              : 
    1973              : !----------------------------------------------------------------------
    1974              : 
    1975              : !!****f* m_xmpi/xmpi_iprobe
    1976              : !! NAME
    1977              : !!  xmpi_iprobe
    1978              : !!
    1979              : !! FUNCTION
    1980              : !!  Hides MPI_IPROBE from MPI library.
    1981              : !!  Nonblocking test for a message.
    1982              : !!
    1983              : !! INPUTS
    1984              : !!  source= source processes
    1985              : !!  tag= tag value
    1986              : !!  mpicomm= communicator
    1987              : !!
    1988              : !! OUTPUT
    1989              : !!  flag= True if a message with the specified source, tag, and communicator is available
    1990              : !!  mpierr= status error
    1991              : !!
    1992              : !! SOURCE
    1993              : 
    1994         4554 : subroutine xmpi_iprobe(source,tag,mpicomm,flag,mpierr)
    1995              : 
    1996              : !Arguments-------------------------
    1997              :  integer,intent(in) :: mpicomm,source,tag
    1998              :  integer,intent(out) :: mpierr
    1999              :  logical,intent(out) :: flag
    2000              : 
    2001              : !Local variables-------------------
    2002              : #ifdef HAVE_MPI
    2003              :  integer :: ier,status(MPI_STATUS_SIZE)
    2004              : #endif
    2005              : ! *************************************************************************
    2006              : 
    2007         4554 :  mpierr = 0
    2008              : #ifdef HAVE_MPI
    2009         4554 :   call MPI_IPROBE(source,tag,mpicomm,flag,status,ier)
    2010         4554 :   mpierr=ier
    2011              : #endif
    2012              : 
    2013         4554 : end subroutine xmpi_iprobe
    2014              : !!***
    2015              : 
    2016              : !----------------------------------------------------------------------
    2017              : 
    2018              : !!****f* m_xmpi/xmpi_wait
    2019              : !! NAME
    2020              : !!  xmpi_wait
    2021              : !!
    2022              : !! FUNCTION
    2023              : !!  Hides MPI_WAIT from MPI library.
    2024              : !!  Waits for an MPI request to complete.
    2025              : !!
    2026              : !! INPUTS
    2027              : !!  request= MPI request handle to wait for
    2028              : !!
    2029              : !! OUTPUT
    2030              : !!  mpierr= status error
    2031              : !!
    2032              : !! SOURCE
    2033              : 
    2034      1594023 : subroutine xmpi_wait(request, mpierr)
    2035              : 
    2036              : !Arguments-------------------------
    2037              :  integer,intent(inout) :: request
    2038              :  integer,intent(out) :: mpierr
    2039              : 
    2040              : !Local variables-------------------
    2041              : #ifdef HAVE_MPI
    2042              :  integer :: ier,status(MPI_STATUS_SIZE)
    2043              : #endif
    2044              : ! *************************************************************************
    2045              : 
    2046      1594023 :  mpierr = 0
    2047              : #ifdef HAVE_MPI
    2048      1594023 :  if (request /= xmpi_request_null) xmpi_count_requests = xmpi_count_requests - 1
    2049      1594023 :  call MPI_WAIT(request,status,ier)
    2050      1594023 :  mpierr=ier
    2051              : #endif
    2052              : 
    2053      1594023 : end subroutine xmpi_wait
    2054              : !!***
    2055              : 
    2056              : !----------------------------------------------------------------------
    2057              : 
    2058              : !!****f* m_xmpi/xmpi_waitall_1d
    2059              : !! NAME
    2060              : !!  xmpi_waitall_1d
    2061              : !!
    2062              : !! FUNCTION
    2063              : !!  Hides MPI_WAITALL from MPI library.
    2064              : !!  Waits for all given MPI Requests to complete.
    2065              : !!
    2066              : !! INPUTS
    2067              : !!  array_of_requests= array of request handles
    2068              : !!
    2069              : !! OUTPUT
    2070              : !!  mpierr= status error
    2071              : !!
    2072              : !! SOURCE
    2073              : 
    2074         9016 : subroutine xmpi_waitall_1d(array_of_requests, mpierr)
    2075              : 
    2076              : !Arguments-------------------------
    2077              :  integer,intent(inout) :: array_of_requests(:)
    2078              :  integer,intent(out) :: mpierr
    2079              : 
    2080              : !Local variables-------------------
    2081              : #ifdef HAVE_MPI
    2082        18032 :  integer :: ier,status(MPI_STATUS_SIZE,size(array_of_requests))
    2083              : #endif
    2084              : ! *************************************************************************
    2085              : 
    2086         9016 :  mpierr = 0
    2087              : #ifdef HAVE_MPI
    2088        42442 :  xmpi_count_requests = xmpi_count_requests - count(array_of_requests /= xmpi_request_null)
    2089         9016 :  call MPI_WAITALL(size(array_of_requests), array_of_requests, status, ier)
    2090         9016 :  mpierr=ier
    2091              : #endif
    2092              : 
    2093         9016 : end subroutine xmpi_waitall_1d
    2094              : !!***
    2095              : 
    2096              : !----------------------------------------------------------------------
    2097              : 
    2098              : !!****f* m_xmpi/xmpi_waitall_2d
    2099              : !! NAME
    2100              : !!  xmpi_waitall_2d
    2101              : !!
    2102              : !! FUNCTION
    2103              : !!  Hides MPI_WAITALL from MPI library.
    2104              : !!  Waits for all given MPI Requests to complete.
    2105              : !!
    2106              : !! INPUTS
    2107              : !!  array_of_requests= array of request handles
    2108              : !!
    2109              : !! OUTPUT
    2110              : !!  mpierr= status error
    2111              : !!
    2112              : !! SOURCE
    2113              : 
    2114         2050 : subroutine xmpi_waitall_2d(array_of_requests, mpierr)
    2115              : 
    2116              : !Arguments-------------------------
    2117              :  integer,intent(inout) :: array_of_requests(:,:)
    2118              :  integer,intent(out) :: mpierr
    2119              : 
    2120              : !Local variables-------------------
    2121         6150 :  integer :: flat_requests(product(shape(array_of_requests)))
    2122              : ! *************************************************************************
    2123              : 
    2124              :  ! MPI_WAITALL is a Fortran interface so cannot pass count and base address a la C
    2125              :  ! so flat 2d array and copy in-out. See https://github.com/open-mpi/ompi/issues/587
    2126         2050 :  flat_requests = pack(array_of_requests, mask=.True.)
    2127         2050 :  call xmpi_waitall_1d(flat_requests, mpierr)
    2128         6150 :  array_of_requests = reshape(flat_requests, shape(array_of_requests))
    2129              : 
    2130         2050 : end subroutine xmpi_waitall_2d
    2131              : !!***
    2132              : 
    2133              : !----------------------------------------------------------------------
    2134              : 
    2135              : !!****f* m_xmpi/xmpi_request_free
    2136              : !! NAME
    2137              : !!  xmpi_request_free
    2138              : !!
    2139              : !! FUNCTION
    2140              : !!  Hides MPI_REQUEST_FREE from MPI library.
    2141              : !!  Frees an array of communication request objects.
    2142              : !!
    2143              : !! INPUTS
    2144              : !!  requests(:)= communication request array (array of handles)
    2145              : !!
    2146              : !! OUTPUT
    2147              : !!  mpierr= status error
    2148              : !!
    2149              : !! SOURCE
    2150              : 
    2151            0 : subroutine xmpi_request_free(requests,mpierr)
    2152              : 
    2153              : !Arguments-------------------------
    2154              :  integer,intent(inout) :: requests(:)
    2155              :  integer,intent(out)  :: mpierr
    2156              : 
    2157              : !Local variables-------------------
    2158              : #ifdef HAVE_MPI
    2159              :  integer :: ier,ii
    2160              : #endif
    2161              : ! *************************************************************************
    2162              : 
    2163            0 :  mpierr = 0
    2164              : #ifdef HAVE_MPI
    2165            0 :  do ii=1,size(requests)
    2166            0 :    if (requests(ii) /= xmpi_request_null) xmpi_count_requests = xmpi_count_requests - 1
    2167            0 :    call MPI_REQUEST_FREE(requests(ii),ier)
    2168              :  end do
    2169            0 :  mpierr=ier
    2170              : #endif
    2171              : 
    2172            0 : end subroutine xmpi_request_free
    2173              : !!***
    2174              : 
    2175              : !!****f* m_xmpi/xmpi_requests_add
    2176              : !! NAME
    2177              : !!  xmpi_requests_add
    2178              : !!
    2179              : !! FUNCTION
    2180              : !!  Increase/decrement xmpi_count_requests internal counter
    2181              : !!
    2182              : !! SOURCE
    2183              : 
    2184            0 : subroutine xmpi_requests_add(count)
    2185              : 
    2186              : !Arguments-------------------------
    2187              :  integer,intent(in) :: count
    2188              : ! *************************************************************************
    2189              : 
    2190            0 :  xmpi_count_requests = xmpi_count_requests + count
    2191              : 
    2192            0 : end subroutine xmpi_requests_add
    2193              : !!***
    2194              : 
    2195              : !----------------------------------------------------------------------
    2196              : 
    2197              : !!****f* m_xmpi/xmpi_error_string
    2198              : !! NAME
    2199              : !!  xmpi_error_string
    2200              : !!
    2201              : !! FUNCTION
    2202              : !!  Hides MPI_ERROR_STRING from MPI library.
    2203              : !!
    2204              : !! INPUTS
    2205              : !!
    2206              : !! OUTPUT
    2207              : !!
    2208              : !! SOURCE
    2209              : 
    2210            0 : subroutine xmpi_error_string(mpierr,err_string,ilen,ierror)
    2211              : 
    2212              : !Arguments-------------------------
    2213              :  integer,intent(in) :: mpierr
    2214              :  integer,intent(out) :: ilen,ierror
    2215              :  character(len=*),intent(out) :: err_string
    2216              : ! *************************************************************************
    2217              : 
    2218            0 :  ilen=0
    2219              : #ifdef HAVE_MPI
    2220            0 :  call MPI_Error_string(mpierr,err_string,ilen,ierror)
    2221              : #else
    2222              :  ierror=1
    2223              :  err_string="Sorry, no MPI_Error_string routine is available to interpret the error message"
    2224              : #endif
    2225              : 
    2226            0 : end subroutine xmpi_error_string
    2227              : !!***
    2228              : 
    2229              : !----------------------------------------------------------------------
    2230              : 
    2231              : !!****f* m_xmpi/xmpi_comm_set_errhandler
    2232              : !! NAME
    2233              : !!  xmpi_set_errhandler
    2234              : !!
    2235              : !! FUNCTION
    2236              : !!  Hides MPI_COMM_SET_ERRHANDLER from MPI library.
    2237              : !!
    2238              : !! INPUTS
    2239              : !!  new_err_handler= new error handler
    2240              : !!
    2241              : !! OUTPUT
    2242              : !!  ierror=error code
    2243              : !!  old_err_handler= old error handler
    2244              : !!
    2245              : !! SIZE EFFECTS
    2246              : !!  comm= communicator (should be intent(in) but is intent(inout) in some
    2247              : !!             OMPI implementation ; known as a bug)
    2248              : !!
    2249              : !! SOURCE
    2250              : 
    2251        18020 : subroutine xmpi_comm_set_errhandler(comm,new_err_handler,old_err_handler,ierror)
    2252              : 
    2253              : !Arguments-------------------------
    2254              :  integer,intent(in) :: new_err_handler
    2255              :  integer,intent(in) :: comm
    2256              :  integer,intent(out) :: ierror,old_err_handler
    2257              : 
    2258              : !Local variables-------------------------
    2259              :  integer :: mpierr1,mpierr2,my_comm
    2260              : ! *************************************************************************
    2261              : 
    2262        18020 :  ierror=0
    2263        18020 :  my_comm = comm  !should be intent(in) but is intent(inout) in some OMPI implementation ; known as a bug)
    2264              : 
    2265              : #if defined HAVE_MPI
    2266              : 
    2267        18020 :  mpierr1=MPI_SUCCESS; mpierr2=MPI_SUCCESS
    2268              : 
    2269              : #if defined HAVE_MPI1
    2270              :    call MPI_Errhandler_get(my_comm,old_err_handler,mpierr1)
    2271              :    call MPI_Errhandler_set(my_comm,new_err_handler,mpierr2)
    2272              : #endif
    2273              : #if defined HAVE_MPI2
    2274        18020 :    call MPI_comm_get_Errhandler(my_comm,old_err_handler,mpierr1)
    2275        18020 :    call MPI_comm_set_Errhandler(my_comm,new_err_handler,mpierr2)
    2276              : #endif
    2277              : 
    2278        18020 :  ierror=MPI_SUCCESS
    2279        18020 :  if (mpierr1/=MPI_SUCCESS) then
    2280            0 :    ierror=mpierr1
    2281        18020 :  else if (mpierr2/=MPI_SUCCESS) then
    2282            0 :    ierror=mpierr2
    2283              :  end if
    2284              : #endif
    2285              : 
    2286        18020 : end subroutine xmpi_comm_set_errhandler
    2287              : !!***
    2288              : 
    2289              : !----------------------------------------------------------------------
    2290              : 
    2291              : !!****f* m_xmpi/xmpi_split_work_i4b
    2292              : !! NAME
    2293              : !!  xmpi_split_work_i4b
    2294              : !!
    2295              : !! FUNCTION
    2296              : !!  Splits the number of tasks, ntasks, among nprocs processors.
    2297              : !!  Used for the MPI parallelization of simple loops.
    2298              : !!
    2299              : !! INPUTS
    2300              : !!  ntasks=number of tasks
    2301              : !!  comm=MPI communicator.
    2302              : !!
    2303              : !! OUTPUT
    2304              : !!  my_start,my_stop= indices defining the initial and final task for this processor
    2305              : !!
    2306              : !! NOTES
    2307              : !!  If nprocs > ntasks then:
    2308              : !!
    2309              : !!    my_start = ntasks + 1
    2310              : !!    my_stop = ntask
    2311              : !!
    2312              : !!  In this particular case, loops of the form
    2313              : !!
    2314              : !!  do ii=my_start,my_stop
    2315              : !!   ...
    2316              : !!  end do
    2317              : !!
    2318              : !!  are not executed. Moreover allocation such as foo(my_start:my_stop) will generate a zero-sized array.
    2319              : !!
    2320              : !! SOURCE
    2321              : 
    2322        10513 : subroutine xmpi_split_work_i4b(ntasks, comm, my_start, my_stop)
    2323              : 
    2324              : !Arguments ------------------------------------
    2325              :  integer,intent(in)  :: ntasks,comm
    2326              :  integer,intent(out) :: my_start, my_stop
    2327              : 
    2328              : !Local variables-------------------------------
    2329              :  integer :: res,nprocs,my_rank,block_p1,block
    2330              : ! *************************************************************************
    2331              : 
    2332        10513 :  nprocs  = xmpi_comm_size(comm); my_rank = xmpi_comm_rank(comm)
    2333              : 
    2334        10513 :  block   = ntasks / nprocs
    2335        10513 :  res     = MOD(ntasks, nprocs)
    2336        10513 :  block_p1= block + 1
    2337              : 
    2338        10513 :  if (my_rank < res) then
    2339           60 :    my_start =  my_rank   *block_p1+1
    2340           60 :    my_stop  = (my_rank+1)*block_p1
    2341              :  else
    2342        10453 :    my_start = res*block_p1 + (my_rank-res  )*block + 1
    2343        10453 :    my_stop  = res*block_p1 + (my_rank-res+1)*block
    2344              :  end if
    2345              : 
    2346        10513 : end subroutine xmpi_split_work_i4b
    2347              : !!***
    2348              : 
    2349              : !----------------------------------------------------------------------
    2350              : 
    2351              : !!****f* m_xmpi/xmpi_split_block
    2352              : !! NAME
    2353              : !!  xmpi_split_block
    2354              : !!
    2355              : !! FUNCTION
    2356              : !!  Splits tasks inside communicator using block distribution. Used for the MPI parallelization of simple loops.
    2357              : !!
    2358              : !! INPUTS
    2359              : !!  ntasks: number of tasks
    2360              : !!  comm: MPI communicator.
    2361              : !!
    2362              : !! OUTPUT
    2363              : !!  my_ntasks: Number of tasks received by this rank. May be zero if ntasks > nprocs.
    2364              : !!  my_inds(my_ntasks): List of tasks treated by this rank. Allocated by the routine. May be zero-sized.
    2365              : !!
    2366              : !! SOURCE
    2367              : 
    2368           96 : subroutine xmpi_split_block(ntasks, comm, my_ntasks, my_inds)
    2369              : 
    2370              : !Arguments ------------------------------------
    2371              :  integer,intent(in)  :: ntasks, comm
    2372              :  integer,intent(out) :: my_ntasks
    2373              :  integer,allocatable,intent(out) :: my_inds(:)
    2374              : 
    2375              : !Local variables-------------------------------
    2376              :  integer :: ii, istart, istop
    2377              : ! *************************************************************************
    2378              : 
    2379           96 :  call xmpi_split_work(ntasks, comm, istart, istop)
    2380           96 :  my_ntasks = istop - istart + 1
    2381          288 :  ABI_MALLOC(my_inds, (my_ntasks))
    2382        12507 :  if (my_ntasks > 0) my_inds = [(istart + (ii - 1), ii=1, my_ntasks)]
    2383              : 
    2384           96 : end subroutine xmpi_split_block
    2385              : !!***
    2386              : 
    2387              : !----------------------------------------------------------------------
    2388              : 
    2389              : !!****f* m_xmpi/xmpi_split_cyclic
    2390              : !! NAME
    2391              : !!  xmpi_split_cyclic
    2392              : !!
    2393              : !! FUNCTION
    2394              : !!  Splits tasks inside communicator using cyclic distribution.
    2395              : !!  Used for the MPI parallelization of simple loops.
    2396              : !!
    2397              : !! INPUTS
    2398              : !!  ntasks: number of tasks
    2399              : !!  comm: MPI communicator.
    2400              : !!
    2401              : !! OUTPUT
    2402              : !!  my_ntasks: Number of tasks received by this rank. May be zero if ntasks > nprocs.
    2403              : !!  my_inds(my_ntasks): List of tasks treated by this rank. Allocated by the routine. May be zero-sized.
    2404              : !!
    2405              : !! SOURCE
    2406              : 
    2407           68 : subroutine xmpi_split_cyclic(ntasks, comm, my_ntasks, my_inds)
    2408              : 
    2409              : !Arguments ------------------------------------
    2410              :  integer,intent(in)  :: ntasks, comm
    2411              :  integer,intent(out) :: my_ntasks
    2412              :  integer,allocatable,intent(out) :: my_inds(:)
    2413              : 
    2414              : !Local variables-------------------------------
    2415              :  integer :: ii, cnt, itask, my_rank, nprocs
    2416              : ! *************************************************************************
    2417              : 
    2418           68 :  nprocs = xmpi_comm_size(comm); my_rank = xmpi_comm_rank(comm)
    2419              : 
    2420          204 :  do ii=1,2
    2421          136 :    if (ii == 2) then
    2422          204 :      ABI_MALLOC(my_inds, (my_ntasks))
    2423              :    end if
    2424          136 :    cnt = 0
    2425          632 :    do itask=1,ntasks
    2426          632 :      if (mod(itask, nprocs) == my_rank) then
    2427          496 :        cnt = cnt + 1
    2428          496 :        if (ii == 2) my_inds(cnt) = itask
    2429              :      end if
    2430              :    end do
    2431          204 :    if (ii == 1) my_ntasks = cnt
    2432              :  end do
    2433              : 
    2434           68 : end subroutine xmpi_split_cyclic
    2435              : !!***
    2436              : 
    2437              : !----------------------------------------------------------------------
    2438              : 
    2439              : !!****f* m_xmpi/xmpi_split_list
    2440              : !! NAME
    2441              : !!  xmpi_split_list
    2442              : !!
    2443              : !! FUNCTION
    2444              : !!  Splits list of items inside communicator using block distribution.
    2445              : !!  Used for the MPI parallelization of simple loops.
    2446              : !!
    2447              : !! INPUTS
    2448              : !!  ntasks:Number of items in list (global)
    2449              : !!  list(ntasks): List of indices
    2450              : !!  comm: MPI communicator.
    2451              : !!
    2452              : !! OUTPUT
    2453              : !!  my_ntasks: Number of tasks received by this rank. May be zero if ntasks > nprocs.
    2454              : !!  my_inds(my_ntasks): List of tasks treated by this rank. Allocated by the routine. May be zero-sized.
    2455              : !!
    2456              : !! SOURCE
    2457              : 
    2458           61 : subroutine xmpi_split_list(ntasks, list, comm, my_ntasks, my_inds)
    2459              : 
    2460              : !Arguments ------------------------------------
    2461              :  integer,intent(in)  :: ntasks, comm
    2462              :  integer,intent(out) :: my_ntasks
    2463              :  integer,intent(in) :: list(ntasks)
    2464              :  integer,allocatable,intent(out) :: my_inds(:)
    2465              : 
    2466              : !Local variables-------------------------------
    2467              :  integer :: my_start, my_stop
    2468              : ! *************************************************************************
    2469              : 
    2470           61 :  call xmpi_split_work(ntasks, comm, my_start, my_stop)
    2471              : 
    2472           61 :  my_ntasks = my_stop - my_start + 1
    2473              : 
    2474           61 :  if (my_stop >= my_start) then
    2475          183 :    ABI_MALLOC(my_inds, (my_ntasks))
    2476         6326 :    my_inds = list(my_start:my_stop)
    2477              :  else
    2478            0 :    my_ntasks = 0
    2479            0 :    ABI_MALLOC(my_inds, (0))
    2480              :  end if
    2481              : 
    2482           61 : end subroutine xmpi_split_list
    2483              : !!***
    2484              : 
    2485              : !----------------------------------------------------------------------
    2486              : 
    2487              : !!****f* m_xmpi/xmpi_split_work2_i4b
    2488              : !! NAME
    2489              : !!  xmpi_split_work2_i4b
    2490              : !!
    2491              : !! FUNCTION
    2492              : !!  Splits a number of tasks, ntasks, among nprocs processors.
    2493              : !!  The output arrays istart(1:nprocs) and istop(1:nprocs)
    2494              : !!  report the starting and final task index for each CPU.
    2495              : !!  Namely CPU with rank ii has to perform all the tasks between
    2496              : !!  istart(ii+1) and istop(ii+1). Note the Fortran convention of using 1 as first index of the array.
    2497              : !!  Note, moreover, that if a proc has rank > ntasks then:
    2498              : !!
    2499              : !!   istart(rank+1)=ntasks+1
    2500              : !!   istop(rank+1)=ntask
    2501              : !!
    2502              : !!  In this particular case, loops of the form
    2503              : !!
    2504              : !!  do ii=istart(rank),istop(rank)
    2505              : !!   ...
    2506              : !!  end do
    2507              : !!
    2508              : !!  are not executed. Moreover allocation such as foo(istart(rank):istop(rank))
    2509              : !!  will generate a zero-sized array
    2510              : !!
    2511              : !! INPUTS
    2512              : !!  ntasks= number of tasks
    2513              : !!  nprocs=Number of processors.
    2514              : !!
    2515              : !! OUTPUT
    2516              : !!  istart(nprocs),istop(nprocs)= indices defining the initial and final task for each processor
    2517              : !!
    2518              : !! SOURCE
    2519              : 
    2520          480 : subroutine xmpi_split_work2_i4b(ntasks, nprocs, istart, istop)
    2521              : 
    2522              : !Arguments ------------------------------------
    2523              :  integer,intent(in) :: ntasks,nprocs
    2524              :  integer,intent(inout) :: istart(nprocs), istop(nprocs)
    2525              : 
    2526              : !Local variables-------------------------------
    2527              :  integer :: res,irank,block,block_tmp
    2528              : ! *************************************************************************
    2529              : 
    2530          480 :  block_tmp = ntasks/nprocs
    2531          480 :  res       = MOD(ntasks,nprocs)
    2532          480 :  block     = block_tmp+1
    2533              : 
    2534         1618 :  do irank=0,nprocs-1
    2535         1618 :    if (irank<res) then
    2536          380 :      istart(irank+1) = irank    *block+1
    2537          380 :      istop (irank+1) = (irank+1)*block
    2538              :    else
    2539          758 :      istart(irank+1) = res*block + (irank-res  )*block_tmp+1
    2540          758 :      istop (irank+1) = res*block + (irank-res+1)*block_tmp
    2541              :    end if
    2542              :  end do
    2543              : 
    2544          480 : end subroutine xmpi_split_work2_i4b
    2545              : !!***
    2546              : 
    2547              : !----------------------------------------------------------------------
    2548              : 
    2549              : !!****f* m_xmpi/xmpi_split_work2_i8b
    2550              : !! NAME
    2551              : !!  xmpi_split_work2_i8b
    2552              : !!
    2553              : !! FUNCTION
    2554              : !!  Same as xmpi_split_work2_i8b but accepts 8 bytes integer.
    2555              : !!
    2556              : !! INPUTS
    2557              : !!  ntasks= number of tasks
    2558              : !!  nprocs=Number of processors.
    2559              : !!
    2560              : !! OUTPUT
    2561              : !!  istart(nprocs),istop(nprocs)= indices defining the initial and final task for each processor
    2562              : !!
    2563              : !! SOURCE
    2564              : 
    2565           24 : subroutine xmpi_split_work2_i8b(ntasks,nprocs,istart,istop)
    2566              : 
    2567              : !Arguments ------------------------------------
    2568              :  integer,intent(in)  :: nprocs
    2569              :  integer(i8b),intent(in)  :: ntasks
    2570              :  integer(i8b),intent(inout) :: istart(nprocs),istop(nprocs)
    2571              : 
    2572              : !Local variables-------------------------------
    2573              :  integer(i8b) :: res,irank,block,block_tmp
    2574              : ! *************************************************************************
    2575              : 
    2576           24 :  block_tmp = ntasks/nprocs
    2577           24 :  res       = MOD(ntasks,INT(nprocs,KIND=i8b))
    2578           24 :  block     = block_tmp+1
    2579              : 
    2580           62 :  do irank=0,nprocs-1
    2581           62 :    if (irank<res) then
    2582            0 :      istart(irank+1)= irank   *block+1
    2583            0 :      istop (irank+1)=(irank+1)*block
    2584              :    else
    2585           38 :      istart(irank+1)=res*block+(irank-res  )*block_tmp+1
    2586           38 :      istop (irank+1)=res*block+(irank-res+1)*block_tmp
    2587              :    end if
    2588              :  end do
    2589              : 
    2590           24 : end subroutine xmpi_split_work2_i8b
    2591              : !!***
    2592              : 
    2593              : !----------------------------------------------------------------------
    2594              : 
    2595              : !!****f* m_xmpi/xmpi_distab_4D
    2596              : !! NAME
    2597              : !!  xmpi_distab_4D
    2598              : !!
    2599              : !! FUNCTION
    2600              : !!  Fill table defining the distribution of the tasks according to the number of processors involved in the
    2601              : !!  calculation. For each set of indices, the table contains the rank of the node in the MPI communicator.
    2602              : !!
    2603              : !! INPUTS
    2604              : !!  nprocs=The number of processors performing the calculation in parallel.
    2605              : !!
    2606              : !! OUTPUT
    2607              : !!  task_distrib(:,:,:,:) = Contains the rank of the node that is taking care of this particular set of loop indices.
    2608              : !!  Tasks are distributed across the nodes in column-major order.
    2609              : !!
    2610              : !! SOURCE
    2611              : 
    2612           27 : subroutine xmpi_distab_4D(nprocs, task_distrib)
    2613              : 
    2614              : !Arguments ------------------------------------
    2615              :  integer,intent(in) :: nprocs
    2616              : !arrays
    2617              :  integer,intent(inout) :: task_distrib(:,:,:,:)
    2618              : 
    2619              : !Local variables ------------------------------
    2620              : !scalars
    2621              :  integer :: ii,jj,n1,n2,n3,n4,ntasks,irank,remainder,ntpblock
    2622              :  integer,allocatable :: list(:)
    2623              : !************************************************************************
    2624              : 
    2625           27 :  n1= SIZE(task_distrib,DIM=1)
    2626           27 :  n2= SIZE(task_distrib,DIM=2)
    2627           27 :  n3= SIZE(task_distrib,DIM=3)
    2628           27 :  n4= SIZE(task_distrib,DIM=4)
    2629           27 :  ntasks = n1*n2*n3*n4
    2630              : 
    2631           81 :  ABI_MALLOC(list, (ntasks))
    2632        29915 :  list=-999
    2633              : 
    2634           27 :  ntpblock  = ntasks/nprocs
    2635           27 :  remainder = MOD(ntasks,nprocs)
    2636              : 
    2637           27 :  if (ntpblock==0) then ! nprocs > ntasks
    2638            0 :    do ii=1,ntasks
    2639            0 :      list(ii) = ii-1
    2640              :    end do
    2641              :  else
    2642           27 :    ii=1
    2643           68 :    do irank=nprocs-1,0,-1 ! If remainder/=0, master will get less tasks.
    2644           41 :      jj = ii+ntpblock-1
    2645           41 :      if (remainder>0) then
    2646            0 :        jj=jj+1
    2647            0 :        remainder = remainder-1
    2648              :      end if
    2649        29929 :      list(ii:jj)=irank
    2650           68 :      ii=jj+1
    2651              :    end do
    2652              :  end if
    2653              : 
    2654          135 :  task_distrib = RESHAPE(list, [n1,n2,n3,n4])
    2655              : 
    2656        35893 :  if (ANY(task_distrib==-999)) call xmpi_abort(msg="task_distrib == -999")
    2657              : 
    2658           27 :  ABI_FREE(list)
    2659              : 
    2660           27 : end subroutine xmpi_distab_4D
    2661              : !!***
    2662              : 
    2663              : !----------------------------------------------------------------------
    2664              : 
    2665              : !!****f* m_xmpi/xmpi_distrib_with_replicas
    2666              : !! NAME
    2667              : !!  xmpi_distrib_with_replicas
    2668              : !!
    2669              : !! FUNCTION
    2670              : !!  This function distributes the i-th task `itask` among `nprocs` inside a MPI communicator.
    2671              : !!  If nprocs > ntasks, multiple MPI ranks will be assigned to a given task.
    2672              : !!
    2673              : !! INPUTS
    2674              : !!  itask=Index of the task (must be <= ntasks)
    2675              : !!  ntasks= number of tasks
    2676              : !!  rank=MPI Rank of this processor in the MPI communicator.
    2677              : !!  nprocs=Number of processors in the MPI communicator.
    2678              : !!
    2679              : !! OUTPUT
    2680              : !!  True if this node will treat itask (replicas are possible if nprocs > ntasks)
    2681              : !!
    2682              : !! SOURCE
    2683              : 
    2684         5081 : pure logical function xmpi_distrib_with_replicas(itask, ntasks, rank, nprocs) result(bool)
    2685              : 
    2686              : !Arguments ------------------------------------
    2687              :  integer,intent(in) :: itask,rank,nprocs,ntasks
    2688              : 
    2689              : !Local variables-------------------------------
    2690              :  integer :: ii,mnp_pool,rk_base
    2691              : ! *************************************************************************
    2692              : 
    2693              :  ! If the number of processors is less than ntasks, we have max one processor per task
    2694              :  ! else we replicate the tasks inside a pool of max size mnp_pool
    2695         5081 :  if (nprocs <= ntasks) then
    2696         5081 :    bool = modulo(itask - 1, nprocs) == rank
    2697              :  else
    2698            0 :    mnp_pool = (nprocs / ntasks)
    2699              :    !write(std_out,*)"Will duplicate itask, mnp_pool", mnp_pool, "nprocs, ntasks", nprocs, ntasks
    2700              : 
    2701            0 :    rk_base = modulo(itask - 1, nprocs)
    2702            0 :    bool = .False.
    2703            0 :    do ii=1,mnp_pool+1
    2704            0 :      if (rank == rk_base + (ii - 1) * ntasks) then
    2705              :        bool = .True.; exit
    2706              :      end if
    2707              :    end do
    2708              :  end if
    2709              : 
    2710         5081 : end function xmpi_distrib_with_replicas
    2711              : !!***
    2712              : 
    2713              : !----------------------------------------------------------------------
    2714              : 
    2715              : !!****f* m_xmpi/xmpi_largetype_create
    2716              : !! NAME
    2717              : !!  xmpi_largetype_create
    2718              : !!
    2719              : !! FUNCTION
    2720              : !!  This function builds a large-count contiguous datatype made of "small" adjacent
    2721              : !!  chunks (of same original type). The new type can then be used in MPI
    2722              : !!  routines when the number of elements to communicate exceeds a 32bit integer.
    2723              : !!
    2724              : !! INPUTS
    2725              : !!  largecount= total number of elements expressed as a 64bit integer
    2726              : !!  inputtype= (INTEGER) input type (typically INTEGER, REAL(dp), ...)
    2727              : !!  op_type= type of operation that will be applied during collective comms
    2728              : !!           At present, MPI_SUM, MPI_LOR, MPI_LAND are implemented
    2729              : !!
    2730              : !! OUTPUT
    2731              : !!  largetype= (INTEGER) new MPI type made of a serie of adjacent chunks
    2732              : !!  largetype_op= (INTEGER) MPI user-defined operation associated to largetype type
    2733              : !!
    2734              : !! NOTE
    2735              : !!  This routine is partially inspired by https://github.com/jeffhammond/BigMPI
    2736              : !!  See: J.R. Hammond. A. Schafer, R. Latham,
    2737              : !!       "ToINT_MAX. . . and beyond. Exploring large-count support in MPI",
    2738              : !!       2014 Workshop on Exascale MPI at Supercomputing Conference
    2739              : !!       MIT License (MIT)
    2740              : !!       Permission is hereby granted, free of charge, to any person obtaining a copy
    2741              : !!       of this software and associated documentation files (the "Software"), to deal
    2742              : !!       in the Software without restriction, including without limitation the rights
    2743              : !!       to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
    2744              : !!       copies of the Software, and to permit persons to whom the Software is
    2745              : !!       furnished to do so.
    2746              : !!
    2747              : !!  From MPI4 specification, this routine is useless as large-count MPI communications
    2748              : !!    can be called with the use of the MPI_count datatype (instead of INTEGER).
    2749              : !!
    2750              : !! SOURCE
    2751              : 
    2752            0 : subroutine xmpi_largetype_create(largecount,inputtype,largetype,largetype_op,op_type)
    2753              : 
    2754              : !Arguments ------------------------------------
    2755              : !scalars
    2756              :  integer(KIND=int64),intent(in) :: largecount
    2757              :  integer,intent(in) :: inputtype,op_type
    2758              :  integer,intent(out) :: largetype,largetype_op
    2759              : 
    2760              : !Local variables-------------------------------
    2761              : #ifdef HAVE_MPI
    2762              : !scalars
    2763              :  integer,parameter :: INT_MAX=max(1,xmpi_maxint32/2)
    2764              :  integer(KIND=int32) :: cc,rr,ierr
    2765              :  integer(KIND=XMPI_ADDRESS_KIND) :: extent,lb,remdisp
    2766              :  integer :: chunks,remainder
    2767              : !arrays
    2768              :  integer(KIND=int32) :: blklens(2)
    2769              :  integer(KIND=XMPI_ADDRESS_KIND) :: disps(2)
    2770              :  integer :: types(2)
    2771              : #endif
    2772              : ! *************************************************************************
    2773              : 
    2774              : #ifdef HAVE_MPI
    2775              :  if (XMPI_ADDRESS_KIND<int64) call xmpi_abort(msg="Too much data to communicate for this architecture!")
    2776              : 
    2777              : !Divide data in chunks
    2778            0 :  cc=int(largecount/INT_MAX,kind=int32)
    2779            0 :  rr=int(largecount-cc*INT_MAX,kind=int32)
    2780              : 
    2781              : !Create user-defined datatype
    2782            0 :  if (rr==0) then
    2783            0 :    call MPI_TYPE_VECTOR(cc,INT_MAX,INT_MAX,inputtype,largetype,ierr)
    2784            0 :    if (ierr==0) call MPI_TYPE_COMMIT(largetype,ierr)
    2785              :  else
    2786            0 :    call MPI_TYPE_VECTOR(cc,INT_MAX,INT_MAX,inputtype,chunks,ierr)
    2787            0 :    call MPI_TYPE_CONTIGUOUS(rr,inputtype,remainder,ierr)
    2788            0 :    if (ierr==0) then
    2789            0 :      call MPI_TYPE_GET_EXTENT(inputtype,lb,extent,ierr)
    2790            0 :      remdisp=cc*INT_MAX*extent
    2791            0 :      blklens(1:2)=1
    2792            0 :      disps(1)=0;disps(2)=remdisp
    2793            0 :      types(1)=chunks;types(2)=remainder
    2794              : #ifdef HAVE_MPI_TYPE_CREATE_STRUCT
    2795            0 :      call MPI_TYPE_CREATE_STRUCT(2,blklens,disps,types,largetype,ierr)
    2796              : #else
    2797              :      call MPI_TYPE_STRUCT(2,blklens,disps,types,largetype,ierr)
    2798              : #endif
    2799            0 :      if (ierr==0) then
    2800            0 :        call MPI_TYPE_COMMIT(largetype,ierr)
    2801            0 :        call MPI_TYPE_FREE(chunks,ierr)
    2802            0 :        call MPI_TYPE_FREE(remainder,ierr)
    2803              :      end if
    2804              :    end if
    2805              :  end if
    2806            0 :  if (ierr/=0) call xmpi_abort(msg="Cannot remove ABI_MPIABORTFILE")
    2807              : 
    2808              : !Associate user-defined MPI operation
    2809            0 :  xmpi_largetype_size=largecount ; largetype_op=-1111
    2810            0 :  if (op_type==MPI_SUM) then
    2811            0 :    select case(inputtype)
    2812              :      case(MPI_INTEGER)
    2813            0 :       call MPI_OP_CREATE(largetype_sum_int  ,.true.,largetype_op,ierr)
    2814              :      case(MPI_REAL)
    2815            0 :       call MPI_OP_CREATE(largetype_sum_real ,.true.,largetype_op,ierr)
    2816              :      case(MPI_DOUBLE_PRECISION)
    2817            0 :       call MPI_OP_CREATE(largetype_sum_dble ,.true.,largetype_op,ierr)
    2818              :      case(MPI_COMPLEX)
    2819            0 :       call MPI_OP_CREATE(largetype_sum_cplx ,.true.,largetype_op,ierr)
    2820              :      case(MPI_DOUBLE_COMPLEX)
    2821            0 :       call MPI_OP_CREATE(largetype_sum_dcplx,.true.,largetype_op,ierr)
    2822              :    end select
    2823            0 :  else if (op_type==MPI_LOR) then
    2824            0 :    select case(inputtype)
    2825              :      case(MPI_LOGICAL)
    2826            0 :       call MPI_OP_CREATE(largetype_lor_log,.true.,largetype_op,ierr)
    2827              :    end select
    2828            0 :  else if (op_type==MPI_LAND) then
    2829            0 :    select case(inputtype)
    2830              :      case(MPI_LOGICAL)
    2831            0 :       call MPI_OP_CREATE(largetype_land_log,.true.,largetype_op,ierr)
    2832              :    end select
    2833              :  else if (op_type==MPI_OP_NULL) then
    2834              :    largetype_op=-1111
    2835              :  end if
    2836              : #else
    2837              :  ABI_UNUSED(largecount)
    2838              :  ABI_UNUSED(inputtype)
    2839              :  ABI_UNUSED(largetype)
    2840              :  ABI_UNUSED(largetype_op)
    2841              :  ABI_UNUSED(op_type)
    2842              : #endif
    2843              : 
    2844            0 : end subroutine xmpi_largetype_create
    2845              : !!***
    2846              : 
    2847              : !--------------------------------------
    2848              : 
    2849              : !!****f* m_xmpi/largetype_sum_int
    2850              : !! NAME
    2851              : !!  largetype_sum_int
    2852              : !!
    2853              : !! FUNCTION
    2854              : !!  Routine used to overload MPI_SUM for integers
    2855              : 
    2856            0 :  subroutine largetype_sum_int(invec,inoutvec,len,datatype)
    2857              :   integer :: len,datatype
    2858              :   integer :: invec(len*xmpi_largetype_size),inoutvec(len*xmpi_largetype_size)
    2859              :   integer(KIND=int64) :: ii,jj,kk
    2860            0 :   kk=0
    2861            0 :   do ii=1,len
    2862            0 :     do jj=1,xmpi_largetype_size
    2863            0 :       kk=kk+1
    2864            0 :       inoutvec(kk)=inoutvec(kk)+invec(kk)
    2865              :     end do
    2866              :   end do
    2867              :   ! this macro is being used befor m_errors is compiled, so work around it
    2868              :   ! ABI_UNUSED(datatype)
    2869              :   if (.FALSE.) write(std_out,*) datatype
    2870            0 :  end subroutine largetype_sum_int
    2871              : !!***
    2872              : 
    2873              : !--------------------------------------
    2874              : 
    2875              : !!****f* m_xmpi/largetype_sum_real
    2876              : !! NAME
    2877              : !!  largetype_sum_real
    2878              : !!
    2879              : !! FUNCTION
    2880              : !!  Routine used to overload MPI_SUM for reals
    2881              : 
    2882            0 :  subroutine largetype_sum_real(invec,inoutvec,len,datatype)
    2883              :   integer :: len,datatype
    2884              :   real(sp) :: invec(len*xmpi_largetype_size),inoutvec(len*xmpi_largetype_size)
    2885              :   integer(KIND=int64) :: ii,jj,kk
    2886            0 :   kk=0
    2887            0 :   do ii=1,len
    2888            0 :     do jj=1,xmpi_largetype_size
    2889            0 :       kk=kk+1
    2890            0 :       inoutvec(kk)=inoutvec(kk)+invec(kk)
    2891              :     end do
    2892              :   end do
    2893              :   ! this macro is being used befor m_errors is compiled, so work around it
    2894              :   ! ABI_UNUSED(datatype)
    2895              :   if (.FALSE.) write(std_out,*) datatype
    2896            0 :  end subroutine largetype_sum_real
    2897              : !!***
    2898              : 
    2899              : !--------------------------------------
    2900              : 
    2901              : !!****f* m_xmpi/largetype_sum_dble
    2902              : !! NAME
    2903              : !!  largetype_sum_dble
    2904              : !!
    2905              : !! FUNCTION
    2906              : !!  Routine used to overload MPI_SUM for double precision reals
    2907              : 
    2908            0 :  subroutine largetype_sum_dble(invec,inoutvec,len,datatype)
    2909              :   integer :: len,datatype
    2910              :   real(dp) :: invec(len*xmpi_largetype_size),inoutvec(len*xmpi_largetype_size)
    2911              :   integer(KIND=int64) :: ii,jj,kk
    2912            0 :   kk=0
    2913            0 :   do ii=1,len
    2914            0 :     do jj=1,xmpi_largetype_size
    2915            0 :       kk=kk+1
    2916            0 :       inoutvec(kk)=inoutvec(kk)+invec(kk)
    2917              :     end do
    2918              :   end do
    2919              :   ! this macro is being used befor m_errors is compiled, so work around it
    2920              :   ! ABI_UNUSED(datatype)
    2921              :   if (.FALSE.) write(std_out,*) datatype
    2922            0 :  end subroutine largetype_sum_dble
    2923              : !!***
    2924              : 
    2925              : !--------------------------------------
    2926              : 
    2927              : !!****f* m_xmpi/largetype_sum_cplx
    2928              : !! NAME
    2929              : !!  largetype_sum_cplx
    2930              : !!
    2931              : !! FUNCTION
    2932              : !!  Routine used to overload MPI_SUM for complex
    2933              : 
    2934            0 :  subroutine largetype_sum_cplx(invec,inoutvec,len,datatype)
    2935              :   integer :: len,datatype
    2936              :   complex(sp) :: invec(len*xmpi_largetype_size),inoutvec(len*xmpi_largetype_size)
    2937              :   integer(KIND=int64) :: ii,jj,kk
    2938            0 :   kk=0
    2939            0 :   do ii=1,len
    2940            0 :     do jj=1,xmpi_largetype_size
    2941            0 :       kk=kk+1
    2942            0 :       inoutvec(kk)=inoutvec(kk)+invec(kk)
    2943              :     end do
    2944              :   end do
    2945              :   ! this macro is being used befor m_errors is compiled, so work around it
    2946              :   ! ABI_UNUSED(datatype)
    2947              :   if (.FALSE.) write(std_out,*) datatype
    2948            0 :  end subroutine largetype_sum_cplx
    2949              : !!***
    2950              : 
    2951              : !--------------------------------------
    2952              : 
    2953              : !!****f* m_xmpi/largetype_sum_dcplx
    2954              : !! NAME
    2955              : !!  largetype_sum_dcplx
    2956              : !!
    2957              : !! FUNCTION
    2958              : !!  Routine used to overload MPI_SUM for double complex
    2959              : 
    2960            0 :  subroutine largetype_sum_dcplx(invec,inoutvec,len,datatype)
    2961              :   integer :: len,datatype
    2962              :   complex(dp) :: invec(len*xmpi_largetype_size),inoutvec(len*xmpi_largetype_size)
    2963              :   integer(KIND=int64) :: ii,jj,kk
    2964            0 :   kk=0
    2965            0 :   do ii=1,len
    2966            0 :     do jj=1,xmpi_largetype_size
    2967            0 :       kk=kk+1
    2968            0 :       inoutvec(kk)=inoutvec(kk)+invec(kk)
    2969              :     end do
    2970              :   end do
    2971              :   ! this macro is being used befor m_errors is compiled, so work around it
    2972              :   ! ABI_UNUSED(datatype)
    2973              :   if (.FALSE.) write(std_out,*) datatype
    2974            0 :  end subroutine largetype_sum_dcplx
    2975              : !!***
    2976              : 
    2977              : !--------------------------------------
    2978              : 
    2979              : !!****f* m_xmpi/largetype_lor_log
    2980              : !! NAME
    2981              : !!  largetype_lor_log
    2982              : !!
    2983              : !! FUNCTION
    2984              : !!  Routine used to overload MPI_LOR for logicals
    2985              : 
    2986            0 :  subroutine largetype_lor_log(invec,inoutvec,len,datatype)
    2987              :   integer :: len,datatype
    2988              :   logical :: invec(len*xmpi_largetype_size),inoutvec(len*xmpi_largetype_size)
    2989              :   integer(KIND=int64) :: ii,jj,kk
    2990            0 :   kk=0
    2991            0 :   do ii=1,len
    2992            0 :     do jj=1,xmpi_largetype_size
    2993            0 :       kk=kk+1
    2994            0 :       inoutvec(kk)=inoutvec(kk).or.invec(kk)
    2995              :     end do
    2996              :   end do
    2997              :   ! this macro is being used befor m_errors is compiled, so work around it
    2998              :   ! ABI_UNUSED(datatype)
    2999              :   if (.FALSE.) write(std_out,*) datatype
    3000            0 :  end subroutine largetype_lor_log
    3001              : !!***
    3002              : 
    3003              : !--------------------------------------
    3004              : 
    3005              : !!****f* m_xmpi/largetype_land_log
    3006              : !! NAME
    3007              : !!  largetype_land_log
    3008              : !!
    3009              : !! FUNCTION
    3010              : !!  Routine used to overload MPI_LANG for logicals
    3011              : 
    3012            0 :  subroutine largetype_land_log(invec,inoutvec,len,datatype)
    3013              :   integer :: len,datatype
    3014              :   logical :: invec(len*xmpi_largetype_size),inoutvec(len*xmpi_largetype_size)
    3015              :   integer(KIND=int64) :: ii,jj,kk
    3016            0 :   kk=0
    3017            0 :   do ii=1,len
    3018            0 :     do jj=1,xmpi_largetype_size
    3019            0 :       kk=kk+1
    3020            0 :       inoutvec(kk)=inoutvec(kk).and.invec(kk)
    3021              :     end do
    3022              :   end do
    3023              :   ! this macro is being used befor m_errors is compiled, so work around it
    3024              :   ! ABI_UNUSED(datatype)
    3025              :   if (.FALSE.) write(std_out,*) datatype
    3026            0 :  end subroutine largetype_land_log
    3027              : !!***
    3028              : 
    3029              : !----------------------------------------------------------------------
    3030              : 
    3031              : !!****f* m_xmpi/xmpi_largetype_free
    3032              : !! NAME
    3033              : !!  xmpi_largetype_free
    3034              : !!
    3035              : !! FUNCTION
    3036              : !!  This function release a large-count contiguous datatype.
    3037              : !!
    3038              : !! SIDE EFFECTS
    3039              : !!  largetype= (INTEGER) MPI type to release
    3040              : !!  largetype_op= (INTEGER) MPI user-defined operation associated to largetype type
    3041              : !!
    3042              : !! SOURCE
    3043              : 
    3044            0 : subroutine xmpi_largetype_free(largetype,largetype_op)
    3045              : 
    3046              : !Arguments ------------------------------------
    3047              : !scalars
    3048              :  integer,intent(inout) :: largetype,largetype_op
    3049              : !Local variables-------------------------------
    3050              : #ifdef HAVE_MPI
    3051              :  integer :: ierr
    3052              : #endif
    3053              : ! *************************************************************************
    3054              : 
    3055              : #ifdef HAVE_MPI
    3056            0 :    xmpi_largetype_size=0
    3057            0 :    if (largetype_op/=-1111) call MPI_OP_FREE(largetype_op,ierr)
    3058            0 :    call MPI_TYPE_FREE(largetype,ierr)
    3059              : #else
    3060              :  ABI_UNUSED(largetype)
    3061              :  ABI_UNUSED(largetype_op)
    3062              : #endif
    3063              : 
    3064            0 : end subroutine xmpi_largetype_free
    3065              : !!***
    3066              : 
    3067              : !----------------------------------------------------------------------
    3068              : 
    3069              : ! Include files providing wrappers for some of the most commonly used MPI primitives.
    3070              : 
    3071              : #include "xmpi_iallgather.finc"
    3072              : #include "xmpi_allgather.finc"
    3073              : #include "xmpi_allgatherv.finc"
    3074              : #include "xmpi_alltoall.finc"
    3075              : #include "xmpi_ialltoall.finc"
    3076              : #include "xmpi_alltoallv.finc"
    3077              : #include "xmpi_ialltoallv.finc"
    3078              : #include "xmpi_bcast.finc"
    3079              : #include "xmpi_ibcast.finc"
    3080              : #include "xmpi_exch.finc"
    3081              : #include "xmpi_gather.finc"
    3082              : #include "xmpi_gatherv.finc"
    3083              : #include "xmpi_max.finc"
    3084              : #include "xmpi_min.finc"
    3085              : #include "xmpi_recv.finc"
    3086              : #include "xmpi_irecv.finc"
    3087              : #include "xmpi_scatterv.finc"
    3088              : #include "xmpi_send.finc"
    3089              : #include "xmpi_isend.finc"
    3090              : #include "xmpi_sum_master.finc"
    3091              : #include "xmpi_sum.finc"
    3092              : #include "xmpi_isum.finc"
    3093              : #include "xmpi_land_lor.finc"
    3094              : 
    3095              : !------------------------------------------------------------------------------------
    3096              : 
    3097              : !!****f* m_xmpi/xmpio_type_struct
    3098              : !! NAME
    3099              : !!  xmpio_type_struct
    3100              : !!
    3101              : !! FUNCTION
    3102              : !!  Some highly non-standard MPI implementations support MPI-IO without
    3103              : !!  implementing the full set of MPI-2 extensions.
    3104              : !!  This wrapper will call the obsolete MPI_TYPE_STRUCT if MPI_TYPE_CREATE_STRUCT
    3105              : !!  is not supported. Note that MPI_TYPE_STRUCT requires the displacement arrays
    3106              : !!  to be an array of default integers whereas the argument block_displ is an array of kind XMPI_ADDRESS_KIND.
    3107              : !!  The routine will abort if the displacement cannot be represented with a default integer.
    3108              : !!
    3109              : !! INPUTS
    3110              : !! ncount= number of blocks (integer) --- also number of entries in arrays
    3111              : !!         array_of_types, array_of_displacements and array_of_blocklengths
    3112              : !! array_of_blocklength(ncount)=number of elements in each block (array of integer)
    3113              : !! array_of_displacements(ncount)=byte displacement of each block (array of integer)
    3114              : !! array_of_types(ncount)=type of elements in each block (array of handles to datatype objects)
    3115              : !!
    3116              : !! OUTPUT
    3117              : !! new_type=new datatype (handle)
    3118              : !! mpierr=MPI status error
    3119              : !!
    3120              : !! SOURCE
    3121              : 
    3122              : #ifdef HAVE_MPI_IO
    3123              : 
    3124          179 : subroutine xmpio_type_struct(ncount, block_length, block_displ, block_type, new_type, mpierr)
    3125              : 
    3126              : !Arguments ------------------------------------
    3127              : !scalars
    3128              :  integer,intent(in) :: ncount
    3129              :  integer,intent(out) :: new_type,mpierr
    3130              : !arrays
    3131              :  integer,intent(in) :: block_length(ncount),block_type(ncount)
    3132              :  integer(XMPI_ADDRESS_KIND),intent(in) :: block_displ(ncount)
    3133              : 
    3134              : !Local variables-------------------
    3135              : #ifndef HAVE_MPI_TYPE_CREATE_STRUCT
    3136              :  integer,allocatable :: tmp_displ(:)
    3137              : #endif
    3138              : !************************************************************************
    3139              : 
    3140              : #ifdef HAVE_MPI_TYPE_CREATE_STRUCT
    3141          174 :  call MPI_TYPE_CREATE_STRUCT(ncount,block_length,block_displ,block_type,new_type,mpierr)
    3142              : #else
    3143              : 
    3144              :  ABI_MALLOC(tmp_displ,(ncount))
    3145              :  tmp_displ = block_displ
    3146              :  if (ANY(block_displ > HUGE(tmp_displ(1)) ))then
    3147              :    call xmpi_abort(msg=" byte displacement cannot be represented with a default integer")
    3148              :  end if
    3149              : 
    3150              :  call MPI_TYPE_STRUCT(ncount,block_length,block_displ,block_type,new_type,mpierr)
    3151              :  ABI_FREE(tmp_displ)
    3152              : #endif
    3153              : 
    3154          174 : end subroutine xmpio_type_struct
    3155              : !!***
    3156              : 
    3157              : #endif
    3158              : 
    3159              : !----------------------------------------------------------------------
    3160              : 
    3161              : !!****f* m_xmpi/xmpio_get_info_frm
    3162              : !! NAME
    3163              : !!  xmpio_marker_info
    3164              : !!
    3165              : !! FUNCTION
    3166              : !!  Return the byte size of the Fortran record and its corresponding MPI_type (compiler-dependent).
    3167              : !!  These two values are needed to access sequential binary Fortran files with MPI/IO routines where
    3168              : !!  C-streams are used.
    3169              : !!
    3170              : !! INPUTS
    3171              : !! comm=MPI communicator. Only master will find the values for the record marker. The results
    3172              : !! are then broadcast to all the other nodes in comm.
    3173              : !!
    3174              : !! OUTPUT
    3175              : !!  bsize_frm=Byte size of the Fortran record marker.
    3176              : !!  mpi_type_frm=MPI type of the marker.
    3177              : !!
    3178              : !! SOURCE
    3179              : 
    3180            0 : subroutine xmpio_get_info_frm(bsize_frm, mpi_type_frm, comm)
    3181              : 
    3182              : !Arguments ------------------------------------
    3183              : !scalars
    3184              :  integer,intent(in) :: comm
    3185              :  integer,intent(out) :: mpi_type_frm,bsize_frm
    3186              : 
    3187              : !Local variables-------------------------------
    3188              :  integer :: my_rank
    3189              : #ifdef HAVE_MPI_IO
    3190              : !scalars
    3191              :  integer,parameter :: master=0
    3192              :  integer :: spt,ept,ii
    3193              :  integer :: f90_unt,iimax,mpio_fh,bsize_int,mpierr
    3194              :  integer(XMPI_OFFSET_KIND) :: offset,rml
    3195              :  character(len=fnlen) :: fname
    3196              :  character(len=500) :: errmsg
    3197              :  logical :: file_exists
    3198              : !arrays
    3199              :  integer :: xvals(2),ivals(100),read_5ivals(5),ref_5ivals(5)
    3200              :  integer :: rm_lengths(4)=(/4,8,2,16/)
    3201              :  integer :: statux(MPI_STATUS_SIZE)
    3202              :  real(dp) :: xrand(fnlen)
    3203              : #endif
    3204              : !************************************************************************
    3205              : 
    3206            0 :  bsize_frm=0; mpi_type_frm=0
    3207              : 
    3208            0 :  my_rank = xmpi_comm_rank(comm) !; RETURN
    3209              : 
    3210              : #ifdef HAVE_MPI_IO
    3211            0 :  if ( my_rank == master ) then
    3212              :    ! Fortran scratch files cannot have a name so have to generate a random one.
    3213              :    ! cannot use pick_aname since it is higher level.
    3214            0 :    fname = "__MPI_IO_FRM__"
    3215            0 :    spt=LEN(trim(fname))+1; ept=spt
    3216              : 
    3217            0 :    inquire(file=trim(fname),exist=file_exists)
    3218              : 
    3219            0 :    do while (file_exists)
    3220            0 :      call RANDOM_NUMBER(xrand(spt:ept))
    3221            0 :      xrand(spt:ept) = 64+xrand(spt:ept)*26
    3222            0 :      do ii=spt,ept
    3223            0 :        fname(ii:ii) = ACHAR(NINT(xrand(ii)))
    3224              :      end do
    3225            0 :      ept = MIN(ept+1,fnlen)
    3226            0 :      inquire(file=trim(fname),exist=file_exists)
    3227              :    end do
    3228              :    !
    3229              :    ! Write five integers on the binary file open in Fortran mode, then try
    3230              :    ! to reread the values with MPI-IO using different offsets for the record marker.
    3231              :    !
    3232            0 :    f90_unt = xmpi_get_unit()
    3233            0 :    if (f90_unt == -1) call xmpi_abort(msg="Cannot find free unit!!")
    3234              :    ! MT dec 2013: suppress the new attribute: often cause unwanted errors
    3235              :    !              and theoretically useless because of the previous inquire
    3236            0 :    open(unit=f90_unt,file=trim(fname),form="unformatted",err=10, iomsg=errmsg)
    3237              : 
    3238            0 :    ref_5ivals = (/(ii, ii=5,9)/)
    3239            0 :    ivals = HUGE(1); ivals(5:9)=ref_5ivals
    3240            0 :    write(f90_unt, err=10, iomsg=errmsg) ivals
    3241            0 :    close(f90_unt, err=10, iomsg=errmsg)
    3242              : 
    3243            0 :    call MPI_FILE_OPEN(xmpi_comm_self, trim(fname), MPI_MODE_RDONLY, MPI_INFO_NULL, mpio_fh,mpierr)
    3244              : 
    3245              :    iimax=3 ! Define number of INTEGER types to be tested
    3246              : #ifdef HAVE_FC_INT_QUAD
    3247            0 :    iimax=4
    3248              : #endif
    3249              :    !
    3250              :    ! Try to read ivals(5:9) from file.
    3251            0 :    ii=0; bsize_frm=-1
    3252            0 :    call MPI_TYPE_SIZE(MPI_INTEGER,bsize_int,mpierr)
    3253              : 
    3254            0 :    do while (bsize_frm<=0 .and. ii<iimax)
    3255            0 :      ii=ii+1
    3256            0 :      rml = rm_lengths(ii)
    3257            0 :      offset = rml + 4 * bsize_int
    3258            0 :      call MPI_FILE_READ_AT(mpio_fh,offset,read_5ivals,5,MPI_INTEGER,statux,mpierr)
    3259              :      !write(std_out,*)read_5ivals
    3260            0 :      if (mpierr==MPI_SUCCESS .and. ALL(read_5ivals==ref_5ivals) ) bsize_frm=rml
    3261              :    end do
    3262              : 
    3263            0 :    if (ii==iimax.and.bsize_frm<=0) then
    3264              :      write(std_out,'(7a)') &
    3265            0 :        'Error during FORTRAN file record marker detection:',ch10,&
    3266            0 :        'It was not possible to read/write a small file!',ch10,&
    3267            0 :        'ACTION: check your access permissions to the file system.',ch10,&
    3268            0 :        'Common sources of this problem: quota limit exceeded, R/W incorrect permissions, ...'
    3269            0 :      call xmpi_abort()
    3270              :    else
    3271              :      !write(std_out,'(a,i0)')' Detected FORTRAN record mark length: ',bsize_frm
    3272              :    end if
    3273              : 
    3274            0 :    call MPI_FILE_CLOSE(mpio_fh, mpierr)
    3275              :    !
    3276              :    ! Select MPI datatype corresponding to the Fortran marker.
    3277            0 :    SELECT CASE (bsize_frm)
    3278              :    CASE (4)
    3279            0 :      mpi_type_frm=MPI_INTEGER4
    3280              :    CASE (8)
    3281            0 :      mpi_type_frm=MPI_INTEGER8
    3282              : #if defined HAVE_FC_INT_QUAD && defined HAVE_MPI_INTEGER16
    3283              :    CASE (16)
    3284            0 :      mpi_type_frm=MPI_INTEGER16
    3285              : #endif
    3286              :    CASE (2)
    3287            0 :      mpi_type_frm=MPI_INTEGER2
    3288              :    CASE DEFAULT
    3289            0 :      write(std_out,'(a,i0)')" Wrong bsize_frm: ",bsize_frm
    3290            0 :      call xmpi_abort()
    3291              :    END SELECT
    3292              : 
    3293            0 :    open(unit=f90_unt,file=trim(fname), err=10, iomsg=errmsg)
    3294            0 :    close(f90_unt,status="delete", err=10, iomsg=errmsg)
    3295              :  end if
    3296              :  !
    3297              :  ! Broadcast data.
    3298            0 :  xvals = (/bsize_frm,mpi_type_frm/)
    3299            0 :  call xmpi_bcast(xvals,master,comm,mpierr)
    3300              : 
    3301            0 :  bsize_frm    = xvals(1)
    3302            0 :  mpi_type_frm = xvals(2)
    3303              : 
    3304            0 :  return
    3305              : 
    3306              : !HANDLE IO ERROR
    3307              : 10 continue
    3308            0 :  call xmpi_abort(msg=errmsg)
    3309              : #endif
    3310              : 
    3311              : end subroutine xmpio_get_info_frm
    3312              : !!***
    3313              : 
    3314              : !----------------------------------------------------------------------
    3315              : 
    3316              : !!****f* m_wffile/xmpio_read_frm
    3317              : !! NAME
    3318              : !!  xmpio_read_frm
    3319              : !!
    3320              : !! FUNCTION
    3321              : !!  Read the content of a single record marker in a FORTRAN file at a given offset using MPI-IO.
    3322              : !!  the file pointer is modified according to the value of advance.
    3323              : !!
    3324              : !! INPUTS
    3325              : !!  fh=MPI-IO file handler.
    3326              : !!  sc_mode=
    3327              : !!    xmpio_single     ==> for reading by current proc.
    3328              : !!    xmpio_collective ==> for collective reading.
    3329              : !!  offset=MPI/IO file pointer
    3330              : !!  [advance]=By default the routine will move the file pointer to the next record.
    3331              : !!    advance=.FALSE. can be used so that the next read will continue picking information
    3332              : !!    off of the currect record.
    3333              : !!
    3334              : !! OUTPUT
    3335              : !!  fmarker=Content of the Fortran record marker.
    3336              : !!  mpierr= MPI error code
    3337              : !!
    3338              : !! SIDE EFFECTS
    3339              : !!  offset=
    3340              : !!     input: file pointer used to access the Fortran marker.
    3341              : !!     output: new offset updated after the reading, depending on advance.
    3342              : !!
    3343              : !! SOURCE
    3344              : 
    3345              : #ifdef HAVE_MPI_IO
    3346              : 
    3347         1042 : subroutine xmpio_read_frm(fh, offset, sc_mode, fmarker, mpierr, advance)
    3348              : 
    3349              : !Arguments ------------------------------------
    3350              : !scalars
    3351              :  integer,intent(in) :: fh,sc_mode
    3352              :  integer(XMPI_OFFSET_KIND),intent(inout) :: offset
    3353              :  integer(XMPI_OFFSET_KIND),intent(out) :: fmarker
    3354              :  integer,intent(out) :: mpierr
    3355              :  logical,optional,intent(in) :: advance
    3356              : 
    3357              : !Local variables-------------------------------
    3358              : !scalars
    3359              :  integer :: bsize_frm,mpi_type_frm,myfh
    3360              :  integer(kind=int16) :: delim_record2(1)
    3361              :  integer(kind=int32) :: delim_record4(1)
    3362              :  integer(kind=int64) :: delim_record8(1)
    3363              : #if defined HAVE_FC_INT_QUAD
    3364              :  integer*16 :: delim_record16(1)
    3365              : #endif
    3366              :  character(len=500) :: msg
    3367              : !arrays
    3368              :  integer :: statux(MPI_STATUS_SIZE)
    3369              : !************************************************************************
    3370              : 
    3371              :  !Workaround for XLF.
    3372         1042 :  myfh = fh
    3373              : 
    3374         1042 :  bsize_frm    = xmpio_bsize_frm    ! Byte size of the Fortran record marker.
    3375         1042 :  mpi_type_frm = xmpio_mpi_type_frm ! MPI type of the record marker.
    3376              : 
    3377         2084 :  SELECT CASE (sc_mode)
    3378              : 
    3379              :  CASE (xmpio_single)
    3380              : 
    3381         1042 :    if (bsize_frm==4) then
    3382         1042 :      call MPI_FILE_READ_AT(myfh,offset,delim_record4,1,mpi_type_frm,statux,mpierr)
    3383         1042 :      fmarker = delim_record4(1)
    3384            0 :    else if (bsize_frm==8) then
    3385            0 :      call MPI_FILE_READ_AT(myfh,offset,delim_record8,1,mpi_type_frm,statux,mpierr)
    3386            0 :      fmarker = delim_record8(1)
    3387              : #if defined HAVE_FC_INT_QUAD
    3388            0 :    else if (bsize_frm==16) then
    3389            0 :      call MPI_FILE_READ_AT(myfh,offset,delim_record16,1,mpi_type_frm,statux,mpierr)
    3390            0 :      fmarker = delim_record16(1)
    3391              : #endif
    3392            0 :    else if (bsize_frm==2) then
    3393            0 :      call MPI_FILE_READ_AT(myfh,offset,delim_record2,1,mpi_type_frm,statux,mpierr)
    3394            0 :      fmarker = delim_record2(1)
    3395              :    else
    3396            0 :      call xmpi_abort(msg='Wrong record marker length!')
    3397              :    end if
    3398              : 
    3399              :  CASE (xmpio_collective)
    3400              : 
    3401            0 :    if (bsize_frm==4) then
    3402            0 :      call MPI_FILE_READ_AT_ALL(myfh,offset,delim_record4,1,mpi_type_frm,statux,mpierr)
    3403            0 :      fmarker = delim_record4(1)
    3404            0 :    else if (bsize_frm==8) then
    3405            0 :      call MPI_FILE_READ_AT_ALL(myfh,offset,delim_record8,1,mpi_type_frm,statux,mpierr)
    3406            0 :      fmarker = delim_record8(1)
    3407              : #if defined HAVE_FC_INT_QUAD
    3408            0 :    else if (bsize_frm==16) then
    3409            0 :      call MPI_FILE_READ_AT_ALL(myfh,offset,delim_record16,1,mpi_type_frm,statux,mpierr)
    3410            0 :      fmarker = delim_record16(1)
    3411              : #endif
    3412            0 :    else if (bsize_frm==2) then
    3413            0 :      call MPI_FILE_READ_AT_ALL(myfh,offset,delim_record2,1,mpi_type_frm,statux,mpierr)
    3414            0 :      fmarker = delim_record2(1)
    3415              :    else
    3416            0 :      call xmpi_abort(msg='Wrong record marker length!')
    3417              :    end if
    3418              : 
    3419              :  CASE DEFAULT
    3420            0 :    write(msg,"(a,i0)")" Wrong value for sc_mode: ",sc_mode
    3421         1042 :    call xmpi_abort(msg=msg)
    3422              :  END SELECT
    3423              : 
    3424         1042 :  if (PRESENT(advance)) then
    3425            0 :    if (advance) then
    3426            0 :      offset = offset + fmarker + 2*bsize_frm ! Move the file pointer to the next record.
    3427              :    else
    3428            0 :      offset = offset + bsize_frm  ! Move the pointer after the marker.
    3429              :    end if
    3430              :  else
    3431         1042 :    offset = offset + fmarker + 2*bsize_frm
    3432              :  end if
    3433              : 
    3434         1042 : end subroutine xmpio_read_frm
    3435              : !!***
    3436              : 
    3437              : #endif
    3438              : 
    3439              : !------------------------------------------------------------------------------------
    3440              : 
    3441              : !!****f* m_wffile/xmpio_write_frm
    3442              : !! NAME
    3443              : !!  xmpio_write_frm
    3444              : !!
    3445              : !! FUNCTION
    3446              : !!  Write a single record marker in a FORTRAN file at a given offset using MPI-IO.
    3447              : !!  The file pointer is modified according to the value of advance.
    3448              : !!
    3449              : !! INPUTS
    3450              : !!  fh=MPI-IO file handler.
    3451              : !!  sc_mode=
    3452              : !!         xmpio_single     ==> for reading by current proc.
    3453              : !!         xmpio_collective ==> for collective reading.
    3454              : !!  fmarker=The content of the Fortran marker i.e. the size of the record in bytes.
    3455              : !!  [advance]=By default the routine will move the file pointer to the next record.
    3456              : !!    advance=.FALSE. can be used so that the next write will continue writing data
    3457              : !!    on the currect record.
    3458              : !!
    3459              : !! OUTPUT
    3460              : !!  mpierr= error code
    3461              : !!
    3462              : !! SIDE EFFECTS
    3463              : !!  offset=
    3464              : !!     input: offset of  the Fortran marker.
    3465              : !!     output: new offset updated after the writing, depending on advance.
    3466              : !!
    3467              : !! SOURCE
    3468              : 
    3469              : #ifdef HAVE_MPI_IO
    3470              : 
    3471          111 : subroutine xmpio_write_frm(fh, offset, sc_mode, fmarker, mpierr, advance)
    3472              : 
    3473              : !Arguments ------------------------------------
    3474              : !scalars
    3475              :  integer,intent(in) :: fh,sc_mode
    3476              :  integer(XMPI_OFFSET_KIND),intent(in) :: fmarker
    3477              :  integer(XMPI_OFFSET_KIND),intent(inout) :: offset
    3478              :  integer,intent(out) :: mpierr
    3479              :  logical,optional,intent(in) :: advance
    3480              : 
    3481              : !Local variables-------------------------------
    3482              : !scalars
    3483              :  integer :: myfh,bsize_frm,mpi_type_frm
    3484              :  integer(XMPI_OFFSET_KIND) :: last
    3485              :  integer(kind=int16)  :: delim_record2
    3486              :  integer(kind=int32)  :: delim_record4
    3487              :  integer(kind=int64)  :: delim_record8
    3488              : #if defined HAVE_FC_INT_QUAD
    3489              :  integer*16 :: delim_record16
    3490              : #endif
    3491              :  character(len=500) :: msg
    3492              : !arrays
    3493              :  integer :: statux(MPI_STATUS_SIZE)
    3494              : !************************************************************************
    3495              : 
    3496              :  ! Workaround for XLF
    3497          111 :  myfh = fh
    3498              : 
    3499          111 :  bsize_frm    = xmpio_bsize_frm      ! Byte size of the Fortran record marker.
    3500          111 :  mpi_type_frm = xmpio_mpi_type_frm   ! MPI type of the record marker.
    3501          111 :  last = offset + bsize_frm + fmarker ! position of the end marker
    3502              : 
    3503          222 :  SELECT CASE (sc_mode)
    3504              : 
    3505              :  CASE (xmpio_single)
    3506          111 :    if (bsize_frm==4) then
    3507          111 :      delim_record4 = fmarker
    3508          222 :      call MPI_FILE_WRITE_AT(myfh,offset,[delim_record4],1,mpi_type_frm,statux,mpierr)
    3509          222 :      call MPI_FILE_WRITE_AT(myfh,last,[delim_record4],1,mpi_type_frm,statux,mpierr)
    3510              : 
    3511            0 :    else if (bsize_frm==8) then
    3512            0 :      delim_record8 = fmarker
    3513            0 :      call MPI_FILE_WRITE_AT(myfh,offset,[delim_record8],1,mpi_type_frm,statux,mpierr)
    3514            0 :      call MPI_FILE_WRITE_AT(myfh,last,[delim_record8],1,mpi_type_frm,statux,mpierr)
    3515              : #if defined HAVE_FC_INT_QUAD
    3516            0 :    else if (bsize_frm==16) then
    3517            0 :      delim_record16 = fmarker
    3518            0 :      call MPI_FILE_WRITE_AT(myfh,offset,[delim_record16],1,mpi_type_frm,statux,mpierr)
    3519            0 :      call MPI_FILE_WRITE_AT(myfh,last,[delim_record16],1,mpi_type_frm,statux,mpierr)
    3520              : #endif
    3521            0 :    else if (bsize_frm==2) then
    3522            0 :      delim_record2 = fmarker
    3523            0 :      call MPI_FILE_WRITE_AT(myfh,offset,[delim_record2], 1,mpi_type_frm,statux,mpierr)
    3524            0 :      call MPI_FILE_WRITE_AT(myfh,last,[delim_record2],1,mpi_type_frm,statux,mpierr)
    3525              :    else
    3526            0 :      call xmpi_abort(msg='Wrong record marker length!')
    3527              :    end if
    3528              : 
    3529              :  CASE (xmpio_collective)
    3530            0 :    if (bsize_frm==4) then
    3531            0 :      delim_record4 = fmarker
    3532            0 :      call MPI_FILE_WRITE_AT_ALL(myfh,offset,[delim_record4],1,mpi_type_frm,statux,mpierr)
    3533            0 :      call MPI_FILE_WRITE_AT_ALL(myfh,last,[delim_record4],1,mpi_type_frm,statux,mpierr)
    3534            0 :    else if (bsize_frm==8) then
    3535            0 :      delim_record8 = fmarker
    3536            0 :      call MPI_FILE_WRITE_AT_ALL(myfh,offset,[delim_record8],1,mpi_type_frm,statux,mpierr)
    3537            0 :      call MPI_FILE_WRITE_AT_ALL(myfh,last,[delim_record8],1,mpi_type_frm,statux,mpierr)
    3538              : #if defined HAVE_FC_INT_QUAD
    3539            0 :    else if (bsize_frm==16) then
    3540            0 :      delim_record16 = fmarker
    3541            0 :      call MPI_FILE_WRITE_AT_ALL(myfh,offset,[delim_record16],1,mpi_type_frm,statux,mpierr)
    3542            0 :      call MPI_FILE_WRITE_AT_ALL(myfh,last,[delim_record16],1,mpi_type_frm,statux,mpierr)
    3543              : #endif
    3544            0 :    else if (bsize_frm==2) then
    3545            0 :      delim_record2 = fmarker
    3546            0 :      call MPI_FILE_WRITE_AT_ALL(myfh,offset,[delim_record2],1,mpi_type_frm,statux,mpierr)
    3547            0 :      call MPI_FILE_WRITE_AT_ALL(myfh,last,[delim_record2],1,mpi_type_frm,statux,mpierr)
    3548              :    else
    3549            0 :      call xmpi_abort(msg='Wrong record marker length!')
    3550              :    end if
    3551              : 
    3552              :  CASE DEFAULT
    3553            0 :    write(msg,"(a,i0)")" Wrong value for sc_mode: ",sc_mode
    3554          111 :    call xmpi_abort(msg=msg)
    3555              :  END SELECT
    3556              : 
    3557          111 :  if (PRESENT(advance)) then
    3558            0 :    if (advance) then
    3559            0 :      offset = offset + fmarker + 2*bsize_frm  ! Move the file pointer to the next record.
    3560              :    else
    3561            0 :      offset = offset + bsize_frm              ! Move the pointer after the marker.
    3562              :    end if
    3563              :  else
    3564          111 :    offset = offset + fmarker + 2*bsize_frm
    3565              :  end if
    3566              : 
    3567          111 : end subroutine xmpio_write_frm
    3568              : !!***
    3569              : #endif
    3570              : 
    3571              : !------------------------------------------------------------------------------------
    3572              : 
    3573              : !!****f* m_xmpi/xmpio_create_fstripes
    3574              : !! NAME
    3575              : !!  xmpio_create_fstripes
    3576              : !!
    3577              : !! FUNCTION
    3578              : !!  Return a MPI type that can be used to (read|write) a set of interleaved Fortran records.
    3579              : !!
    3580              : !!  <FRM> type(1), type(1), ... <FRM>  ! size(1) elements
    3581              : !!  <FRM> type(2), type(2), ... <FRM>  ! size(2) elements
    3582              : !!  <FRM> type(1), type(1), ... <FRM>  ! size(1) elements
    3583              : !!  ....
    3584              : !!
    3585              : !! INPUTS
    3586              : !!  ncount = Number of records with elements of type types(1) to (read|write)
    3587              : !!  sizes(1:2) = Number of elements of each type in the two sets of record
    3588              : !!  type(1:2) = MPI Type of the elements in the first and in the second record.
    3589              : !!
    3590              : !! OUTPUT
    3591              : !!  my_offpad=Offset to be added to the file pointer giving the position of the first Fortran record
    3592              : !!    marker individuating the beginning of the matrix. (lets call it "base").
    3593              : !!    Each node should (read|write) using my_offset = base + my_offpad.
    3594              : !!    my_offpad is used so that one can safely change the way the fileview is generated (for example
    3595              : !!    to make it more efficient) without having to change the client code.
    3596              : !!  new_type=New MPI type.
    3597              : !!  mpierr= MPI error code
    3598              : !!
    3599              : !! SOURCE
    3600              : 
    3601              : #ifdef HAVE_MPI_IO
    3602              : 
    3603            0 : subroutine xmpio_create_fstripes(ncount, sizes, types, new_type, my_offpad, mpierr)
    3604              : 
    3605              : !Arguments ------------------------------------
    3606              : !scalars
    3607              :  integer,intent(in) :: ncount
    3608              :  integer(XMPI_OFFSET_KIND),intent(out) :: my_offpad
    3609              :  integer,intent(out) :: new_type,mpierr
    3610              : !arrays
    3611              :  integer,intent(in) :: types(2),sizes(2)
    3612              : 
    3613              : !Local variables-------------------------------
    3614              : !scalars
    3615              :  integer :: type_x,type_y,bsize_frm,bsize_x,bsize_y,nx,ny,column_type
    3616              :  integer(MPI_ADDRESS_KIND) :: stride
    3617              : !************************************************************************
    3618              : 
    3619              :  ! Byte size of the Fortran record marker.
    3620            0 :  bsize_frm = xmpio_bsize_frm
    3621              : 
    3622              :  ! Number of elements in the two stripes.
    3623            0 :  nx = sizes(1)
    3624            0 :  ny = sizes(2)
    3625              : 
    3626            0 :  type_x = types(1)
    3627            0 :  type_y = types(2)
    3628              : 
    3629              :  ! Byte size of type_x and type_y
    3630            0 :  call MPI_TYPE_SIZE(type_x,bsize_x,mpierr)
    3631            0 :  ABI_HANDLE_MPIERR(mpierr)
    3632              : 
    3633            0 :  call MPI_TYPE_SIZE(type_y,bsize_y,mpierr)
    3634            0 :  ABI_HANDLE_MPIERR(mpierr)
    3635              : 
    3636              :  ! The view starts at the first element of the first stripe.
    3637            0 :  my_offpad = xmpio_bsize_frm
    3638              : 
    3639            0 :  call MPI_Type_contiguous(nx,type_x,column_type,mpierr)
    3640            0 :  ABI_HANDLE_MPIERR(mpierr)
    3641              : 
    3642              :  ! Byte size of the Fortran record + the two markers.
    3643            0 :  stride = nx*bsize_x + 2*bsize_frm  + ny*bsize_y + 2*bsize_frm
    3644              : 
    3645              :  ! ncount colum_type separated by stride bytes
    3646            0 :  if (ncount>0) then
    3647            0 :    call MPI_Type_create_hvector(ncount,1,stride,column_type,new_type,mpierr)
    3648              :  else
    3649            0 :    call MPI_Type_create_hvector(1,1,stride,column_type,new_type,mpierr)
    3650              :  end if
    3651            0 :  ABI_HANDLE_MPIERR(mpierr)
    3652              : 
    3653            0 :  call MPI_TYPE_COMMIT(new_type,mpierr)
    3654            0 :  ABI_HANDLE_MPIERR(mpierr)
    3655              : 
    3656            0 :  call MPI_TYPE_FREE(column_type,mpierr)
    3657            0 :  ABI_HANDLE_MPIERR(mpierr)
    3658              : 
    3659              : end subroutine xmpio_create_fstripes
    3660              : !!***
    3661              : #endif
    3662              : 
    3663              : !------------------------------------------------------------------------------------
    3664              : 
    3665              : !!****f* m_xmpi/xmpio_create_fsubarray_2D
    3666              : !! NAME
    3667              : !!  xmpio_create_fsubarray_2D
    3668              : !!
    3669              : !! FUNCTION
    3670              : !!  Return a MPI type that can be used to (read|write) a 2D matrix of elements of type old_type stored in a Fortran file.
    3671              : !!
    3672              : !! INPUTS
    3673              : !!  sizes(2)=number of elements of type old_type in each dimension of the full array (array of positive integers)
    3674              : !!  subsizes(2)=number of elements of type old_type in each dimension of the subarray (array of positive integers)
    3675              : !!  array_of_starts(2)=starting coordinates of the subarray in each dimension (array of nonnegative integers >=1, <=sizes)
    3676              : !!  old_type=Old MPI type.
    3677              : !!
    3678              : !! OUTPUT
    3679              : !!  my_offpad=Offset to be added to the file pointer giving the position of the first Fortran record
    3680              : !!    marker individuating the beginning of the matrix. (lets call it "base").
    3681              : !!    Each node should (read|write) using my_offset = base + my_offpad.
    3682              : !!    my_offpad is used so that one can safely change the way the fileview is generated (for example
    3683              : !!    to make it more efficient) without having to change the client code.
    3684              : !!  new_type=New MPI type.
    3685              : !!  mpierr= MPI error code
    3686              : !!
    3687              : !! SOURCE
    3688              : 
    3689              : #ifdef HAVE_MPI_IO
    3690              : 
    3691           30 : subroutine xmpio_create_fsubarray_2D(sizes, subsizes, array_of_starts, old_type, new_type, my_offpad, mpierr)
    3692              : 
    3693              : !Arguments ------------------------------------
    3694              : !scalars
    3695              :  integer,intent(in) :: old_type
    3696              :  integer(XMPI_OFFSET_KIND),intent(out) :: my_offpad
    3697              :  integer,intent(out) :: mpierr,new_type
    3698              : !arrays
    3699              :  integer,intent(in) :: sizes(2),subsizes(2),array_of_starts(2)
    3700              : !Local variables-------------------------------
    3701              : !scalars
    3702              :  integer :: bsize_frm,bsize_old,nx,ny,column_type,ldx
    3703              :  integer(XMPI_OFFSET_KIND) :: st_x,st_y
    3704              :  integer(MPI_ADDRESS_KIND) :: stride_x
    3705              :  !character(len=500) :: msg
    3706              : !************************************************************************
    3707              : 
    3708              :  ! Byte size of the Fortran record marker.
    3709           30 :  bsize_frm = xmpio_bsize_frm
    3710              : 
    3711              :  ! Byte size of old_type.
    3712           30 :  call MPI_TYPE_SIZE(old_type,bsize_old,mpierr)
    3713           30 :  ABI_HANDLE_MPIERR(mpierr)
    3714              :  !
    3715              :  ! Number of columns and rows of the submatrix.
    3716           30 :  nx = subsizes(1)
    3717           30 :  ny = subsizes(2)
    3718              : 
    3719           30 :  ldx = sizes(1)
    3720           30 :  st_x = array_of_starts(1)
    3721           30 :  st_y = array_of_starts(2)
    3722              : 
    3723              :  ! The view starts at the first element of the submatrix.
    3724           30 :  my_offpad = (st_x-1)*bsize_old + (st_y-1)*(ldx*bsize_old+2*xmpio_bsize_frm) + xmpio_bsize_frm
    3725              : 
    3726              :  ! Byte size of the Fortran record + the two markers.
    3727           30 :  stride_x = ldx*bsize_old + 2*bsize_frm
    3728              : 
    3729           30 :  call MPI_Type_contiguous(nx,old_type,column_type,mpierr)
    3730           30 :  ABI_HANDLE_MPIERR(mpierr)
    3731              : 
    3732           30 :  call MPI_Type_create_hvector(ny,1,stride_x,column_type,new_type,mpierr)
    3733           30 :  ABI_HANDLE_MPIERR(mpierr)
    3734              : 
    3735           30 :  call MPI_TYPE_COMMIT(new_type,mpierr)
    3736           30 :  ABI_HANDLE_MPIERR(mpierr)
    3737              : 
    3738           30 :  call MPI_TYPE_FREE(column_type, mpierr)
    3739           30 :  ABI_HANDLE_MPIERR(mpierr)
    3740              : 
    3741              : end subroutine xmpio_create_fsubarray_2D
    3742              : !!***
    3743              : #endif
    3744              : 
    3745              : !------------------------------------------------------------------------------------
    3746              : 
    3747              : !!****f* m_xmpi/xmpio_create_fsubarray_3D
    3748              : !! NAME
    3749              : !!  xmpio_create_fsubarray_3D
    3750              : !!
    3751              : !! FUNCTION
    3752              : !!  Return a MPI type that can be used to (read|write) a 3D matrix of elements of type old_type stored in a Fortran file.
    3753              : !!
    3754              : !! INPUTS
    3755              : !!  sizes(3)=number of elements of type old_type in each dimension of the full array (array of positive integers)
    3756              : !!  subsizes(3)=number of elements of type old_type in each dimension of the subarray (array of positive integers)
    3757              : !!  array_of_starts(3)=starting coordinates of the subarray in each dimension (array of nonnegative integers >=1, <=sizes)
    3758              : !!  old_type=Old MPI type.
    3759              : !!
    3760              : !! OUTPUT
    3761              : !!  my_offpad=Offset to be added to the file pointer giving the position of the first Fortran record
    3762              : !!    marker individuating the beginning of the matrix. (lets call it "base").
    3763              : !!    Each node should (read|write) using my_offset = base + my_offpad.
    3764              : !!    my_offpad is used so that one can safely change the way the fileview is generated (for example
    3765              : !!    to make it more efficient) without having to change the client code.
    3766              : !!  new_type=New MPI type.
    3767              : !!  mpierr= MPI error code
    3768              : !!
    3769              : !! SOURCE
    3770              : 
    3771              : #ifdef HAVE_MPI_IO
    3772              : 
    3773            0 : subroutine xmpio_create_fsubarray_3D(sizes, subsizes, array_of_starts, old_type, new_type, my_offpad, mpierr)
    3774              : 
    3775              : !Arguments ------------------------------------
    3776              : !scalars
    3777              :  integer,intent(in) :: old_type
    3778              :  integer,intent(out) :: mpierr,new_type
    3779              :  integer(XMPI_OFFSET_KIND),intent(out) :: my_offpad
    3780              : !arrays
    3781              :  integer,intent(in) :: sizes(3),subsizes(3),array_of_starts(3)
    3782              : !Local variables-------------------------------
    3783              : !scalars
    3784              :  integer :: bsize_frm,bsize_old,nx,ny,nz
    3785              :  integer :: column_type,plane_type,ldx,ldy,ldz
    3786              :  integer(XMPI_OFFSET_KIND) :: st_x,st_y,st_z
    3787              :  integer(MPI_ADDRESS_KIND) :: stride_x
    3788              :  !character(len=500) :: msg
    3789              : !************************************************************************
    3790              : 
    3791            0 :  bsize_frm = xmpio_bsize_frm    ! Byte size of the Fortran record marker.
    3792              : 
    3793              :  ! Byte size of old_type.
    3794            0 :  call MPI_TYPE_SIZE(old_type,bsize_old,mpierr)
    3795            0 :  ABI_HANDLE_MPIERR(mpierr)
    3796              :  !
    3797              :  ! Number of columns and rows of the submatrix.
    3798            0 :  nx = subsizes(1)
    3799            0 :  ny = subsizes(2)
    3800            0 :  nz = subsizes(3)
    3801              : 
    3802            0 :  ldx = sizes(1)
    3803            0 :  ldy = sizes(2)
    3804            0 :  ldz = sizes(3)
    3805              : 
    3806            0 :  st_x = array_of_starts(1)
    3807            0 :  st_y = array_of_starts(2)
    3808            0 :  st_z = array_of_starts(3)
    3809              : 
    3810              :  ! The view starts at the first element of the submatrix.
    3811              :  my_offpad = (st_x-1)*bsize_old + &
    3812              :              (st_y-1)*    (ldx*bsize_old+2*xmpio_bsize_frm) + &
    3813              :              (st_z-1)*ldy*(ldx*bsize_old+2*xmpio_bsize_frm) + &
    3814            0 :               xmpio_bsize_frm
    3815              : 
    3816              :  ! Byte size of the Fortran record + the two markers.
    3817            0 :  stride_x = ldx*bsize_old + 2*bsize_frm
    3818              : 
    3819            0 :  call MPI_Type_contiguous(nx,old_type,column_type,mpierr)
    3820            0 :  ABI_HANDLE_MPIERR(mpierr)
    3821              : 
    3822            0 :  call MPI_Type_create_hvector(ny,1,stride_x,column_type,plane_type,mpierr)
    3823            0 :  ABI_HANDLE_MPIERR(mpierr)
    3824              : 
    3825            0 :  call MPI_Type_create_hvector(nz,1,ldy*stride_x,plane_type,new_type,mpierr)
    3826            0 :  ABI_HANDLE_MPIERR(mpierr)
    3827              : 
    3828              :  ! Commit the datatype
    3829            0 :  call MPI_TYPE_COMMIT(new_type,mpierr)
    3830            0 :  ABI_HANDLE_MPIERR(mpierr)
    3831              : 
    3832              :  ! Free memory
    3833            0 :  call MPI_TYPE_FREE(plane_type, mpierr)
    3834            0 :  ABI_HANDLE_MPIERR(mpierr)
    3835              : 
    3836              : end subroutine xmpio_create_fsubarray_3D
    3837              : !!***
    3838              : #endif
    3839              : 
    3840              : !------------------------------------------------------------------------------------
    3841              : 
    3842              : !!****f* m_xmpi/xmpio_create_fsubarray_4D
    3843              : !! NAME
    3844              : !!  xmpio_create_fsubarray_4D
    3845              : !!
    3846              : !! FUNCTION
    3847              : !!  Return a MPI type that can be used to (read|write) a 2D matrix of elements of type old_type stored in a Fortran file.
    3848              : !!
    3849              : !! INPUTS
    3850              : !!  sizes(4)=number of elements of type old_type in each dimension of the full array (array of positive integers)
    3851              : !!  subsizes(4)=number of elements of type old_type in each dimension of the subarray (array of positive integers)
    3852              : !!  array_of_starts(4)=starting coordinates of the subarray in each dimension (array of nonnegative integers >=1, <=sizes)
    3853              : !!  old_type=Old MPI type.
    3854              : !!
    3855              : !! OUTPUT
    3856              : !!  my_offpad=Offset to be added to the file pointer giving the position of the first Fortran record
    3857              : !!    marker individuating the beginning of the matrix. (lets call it "base").
    3858              : !!    Each node should (read|write) using my_offset = base + my_offpad.
    3859              : !!    my_offpad is used so that one can safely change the way the fileview is generated (for example
    3860              : !!    to make it more efficient) without having to change the client code.
    3861              : !!  new_type=New MPI type.
    3862              : !!  mpierr= MPI error code
    3863              : !!
    3864              : !! SOURCE
    3865              : 
    3866              : #ifdef HAVE_MPI_IO
    3867              : 
    3868            0 : subroutine xmpio_create_fsubarray_4D(sizes, subsizes, array_of_starts, old_type, new_type, my_offpad, mpierr)
    3869              : 
    3870              : !Arguments ------------------------------------
    3871              : !scalars
    3872              :  integer,intent(in) :: old_type
    3873              :  integer,intent(out) :: mpierr,new_type
    3874              :  integer(XMPI_OFFSET_KIND),intent(out) :: my_offpad
    3875              : !arrays
    3876              :  integer,intent(in) :: sizes(4),subsizes(4),array_of_starts(4)
    3877              : 
    3878              : !Local variables-------------------------------
    3879              : !scalars
    3880              :  integer :: bsize_frm,bsize_old,nx,ny,nz,na
    3881              :  integer :: column_type,plane_type,ldx,ldy,ldz,lda,vol_type
    3882              :  integer(XMPI_OFFSET_KIND) :: st_x,st_y,st_z,st_a
    3883              :  integer(MPI_ADDRESS_KIND) :: stride_x
    3884              : !************************************************************************
    3885              : 
    3886            0 :  bsize_frm = xmpio_bsize_frm    ! Byte size of the Fortran record marker.
    3887              : 
    3888              :  ! Byte size of old_type.
    3889            0 :  call MPI_TYPE_SIZE(old_type,bsize_old,mpierr)
    3890            0 :  ABI_HANDLE_MPIERR(mpierr)
    3891              :  !
    3892              :  ! Number of columns and rows of the submatrix.
    3893            0 :  nx = subsizes(1)
    3894            0 :  ny = subsizes(2)
    3895            0 :  nz = subsizes(3)
    3896            0 :  na = subsizes(4)
    3897              : 
    3898            0 :  ldx = sizes(1)
    3899            0 :  ldy = sizes(2)
    3900            0 :  ldz = sizes(3)
    3901            0 :  lda = sizes(4)
    3902              : 
    3903            0 :  st_x = array_of_starts(1)
    3904            0 :  st_y = array_of_starts(2)
    3905            0 :  st_z = array_of_starts(3)
    3906            0 :  st_a = array_of_starts(4)
    3907              : 
    3908              :  ! The view starts at the first element of the submatrix.
    3909              :  my_offpad = (st_x-1)*bsize_old + &
    3910              :              (st_y-1)*        (ldx*bsize_old+2*xmpio_bsize_frm) + &
    3911              :              (st_z-1)*ldy*    (ldx*bsize_old+2*xmpio_bsize_frm) + &
    3912              :              (st_a-1)*lda*ldy*(ldx*bsize_old+2*xmpio_bsize_frm) + &
    3913            0 :              xmpio_bsize_frm
    3914              : 
    3915              :  ! Byte size of the Fortran record + the two markers.
    3916            0 :  stride_x = ldx*bsize_old + 2*bsize_frm
    3917              : 
    3918            0 :  call MPI_Type_contiguous(nx,old_type,column_type,mpierr)
    3919            0 :  ABI_HANDLE_MPIERR(mpierr)
    3920              : 
    3921            0 :  call MPI_Type_create_hvector(ny,1,stride_x,column_type,plane_type,mpierr)
    3922            0 :  ABI_HANDLE_MPIERR(mpierr)
    3923              : 
    3924            0 :  call MPI_Type_create_hvector(nz,1,ldy*stride_x,plane_type,vol_type,mpierr)
    3925            0 :  ABI_HANDLE_MPIERR(mpierr)
    3926              : 
    3927            0 :  call MPI_Type_create_hvector(na,1,ldz*ldy*stride_x,vol_type,new_type,mpierr)
    3928            0 :  ABI_HANDLE_MPIERR(mpierr)
    3929              : 
    3930              :  ! Commit the datatype
    3931            0 :  call MPI_TYPE_COMMIT(new_type,mpierr)
    3932            0 :  ABI_HANDLE_MPIERR(mpierr)
    3933              : 
    3934              :  ! Free memory
    3935            0 :  call MPI_TYPE_FREE(column_type, mpierr)
    3936            0 :  ABI_HANDLE_MPIERR(mpierr)
    3937              : 
    3938            0 :  call MPI_TYPE_FREE(plane_type, mpierr)
    3939            0 :  ABI_HANDLE_MPIERR(mpierr)
    3940              : 
    3941            0 :  call MPI_TYPE_FREE(vol_type, mpierr)
    3942            0 :  ABI_HANDLE_MPIERR(mpierr)
    3943              : 
    3944              : end subroutine xmpio_create_fsubarray_4D
    3945              : !!***
    3946              : #endif
    3947              : 
    3948              : !------------------------------------------------------------------------------------
    3949              : 
    3950              : !!****f* m_xmpi/xmpio_check_frmarkers
    3951              : !! NAME
    3952              : !!  xmpio_check_frmarkers
    3953              : !!
    3954              : !! FUNCTION
    3955              : !!  Check a set of Fortran record markers starting at a given offset using MPI-IO.
    3956              : !!
    3957              : !! INPUTS
    3958              : !!  fh=MPI-IO file handler.
    3959              : !!  offset=MPI-IO file pointer
    3960              : !!  sc_mode=Option for individual or collective reading.
    3961              : !!  nfrec=Number of Fortran records to be checked.
    3962              : !!  bsize_frecord(nfrec)=Byte size of the Fortran records (markers are NOT included)
    3963              : !!    These values will be compared with the markers reported in the file.
    3964              : !!
    3965              : !! OUTPUT
    3966              : !!  ierr=A non-zero error code signals failure.
    3967              : !!
    3968              : !! SOURCE
    3969              : 
    3970              : #ifdef HAVE_MPI_IO
    3971              : 
    3972            0 : subroutine xmpio_check_frmarkers(fh, offset, sc_mode, nfrec, bsize_frecord, ierr)
    3973              : 
    3974              : !Arguments ------------------------------------
    3975              : !scalars
    3976              :  integer,intent(in) :: fh,nfrec,sc_mode
    3977              :  integer(XMPI_OFFSET_KIND),intent(in) :: offset
    3978              :  integer,intent(out) :: ierr
    3979              : !arrays
    3980              :  integer(XMPI_OFFSET_KIND),intent(in) :: bsize_frecord(nfrec)
    3981              : 
    3982              : !Local variables-------------------------------
    3983              : !scalars
    3984              :  integer :: nb,irec,frmarkers_type,jj,bsize_frm,mpi_type_frm,mpierr,myfh
    3985              :  integer(XMPI_OFFSET_KIND) :: displ
    3986              : !arrays
    3987            0 :  integer(kind=int16),allocatable :: bufdelim2(:)
    3988            0 :  integer(kind=int32),allocatable :: bufdelim4(:)
    3989            0 :  integer(kind=int64),allocatable :: bufdelim8(:)
    3990              : #ifdef HAVE_FC_INT_QUAD
    3991            0 :  integer*16,allocatable :: bufdelim16(:)
    3992              : #endif
    3993              : !integer :: statux(MPI_STATUS_SIZE)
    3994            0 :  integer,allocatable :: block_length(:),block_type(:)
    3995            0 :  integer(XMPI_ADDRESS_KIND),allocatable :: block_displ(:)
    3996            0 :  integer(XMPI_OFFSET_KIND),allocatable :: delim_record(:)
    3997              : !************************************************************************
    3998              : 
    3999              :  ! Workaround for XLF
    4000            0 :  myfh = fh
    4001            0 :  ierr=0
    4002              : 
    4003            0 :  bsize_frm    = xmpio_bsize_frm     ! Byte size of the Fortran record marker.
    4004            0 :  mpi_type_frm = xmpio_mpi_type_frm  ! MPI type of the record marker.
    4005              : 
    4006              :  ! Define the view for the file.
    4007            0 :  nb=2*nfrec
    4008            0 :  ABI_MALLOC(block_length,(nb+2))
    4009            0 :  ABI_MALLOC(block_displ,(nb+2))
    4010            0 :  ABI_MALLOC(block_type,(nb+2))
    4011            0 :  block_length(1)=1
    4012            0 :  block_displ (1)=0
    4013            0 :  block_type  (1)=MPI_LB
    4014              : 
    4015            0 :  jj=2; displ=0
    4016            0 :  do irec=1,nfrec
    4017            0 :    block_type (jj:jj+1) =mpi_type_frm
    4018            0 :    block_length(jj:jj+1)=1
    4019            0 :    block_displ(jj  )     = displ
    4020            0 :    block_displ(jj+1)     = bsize_frm + displ + bsize_frecord(irec)
    4021            0 :    jj=jj+2
    4022            0 :    displ = displ + bsize_frecord(irec) + 2*bsize_frm ! Move to the beginning of the next column.
    4023            0 :    if (xmpio_max_address(displ)) ierr=-1  ! Check for wraparound.
    4024              :  end do
    4025              : 
    4026            0 :  block_length(nb+2)=1
    4027            0 :  block_displ (nb+2)=displ
    4028            0 :  block_type  (nb+2)=MPI_UB
    4029              : 
    4030            0 :  call xmpio_type_struct(nb+2,block_length,block_displ,block_type,frmarkers_type,mpierr)
    4031            0 :  ABI_FREE(block_length)
    4032            0 :  ABI_FREE(block_displ)
    4033            0 :  ABI_FREE(block_type)
    4034              : 
    4035            0 :  call MPI_TYPE_COMMIT(frmarkers_type,mpierr)
    4036            0 :  call MPI_FILE_SET_VIEW(myfh,offset,MPI_BYTE,frmarkers_type,"native",MPI_INFO_NULL,mpierr)
    4037              : 
    4038            0 :  jj=1
    4039            0 :  ABI_MALLOC(delim_record,(nb))
    4040            0 :  do irec=1,nfrec
    4041            0 :    delim_record(jj:jj+1)=bsize_frecord(irec)
    4042            0 :    jj=jj+2
    4043              :  end do
    4044              : 
    4045              :  ! Read markers according to the MPI type of the Fortran marker.
    4046            0 :  SELECT CASE (bsize_frm)
    4047              : 
    4048              :  CASE (4)
    4049            0 :    ABI_MALLOC(bufdelim4,(nb))
    4050            0 :    if (sc_mode==xmpio_single) then
    4051            0 :      call MPI_FILE_READ    (myfh,bufdelim4,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4052            0 :    else if (sc_mode==xmpio_collective) then
    4053            0 :      call MPI_FILE_READ_ALL(myfh,bufdelim4,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4054              :    else
    4055            0 :      ierr=2
    4056              :    end if
    4057            0 :    if (ANY(bufdelim4/=delim_record)) ierr=1
    4058            0 :    if (ierr==1) then
    4059            0 :      do irec=1,2*nfrec
    4060            0 :        write(std_out,*)"irec, bufdelim4, delim_record: ",irec,bufdelim4(irec),delim_record(irec)
    4061              :      end do
    4062              :    end if
    4063            0 :    ABI_FREE(bufdelim4)
    4064              : 
    4065              :  CASE (8)
    4066            0 :    ABI_MALLOC(bufdelim8,(nb))
    4067            0 :    if (sc_mode==xmpio_single) then
    4068            0 :      call MPI_FILE_READ    (myfh,bufdelim8,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4069            0 :    else if (sc_mode==xmpio_collective) then
    4070            0 :      call MPI_FILE_READ_ALL(myfh,bufdelim8,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4071              :    else
    4072            0 :      ierr=2
    4073              :    end if
    4074            0 :    if (ANY(bufdelim8/=delim_record)) ierr=1
    4075            0 :    ABI_FREE(bufdelim8)
    4076              : 
    4077              : #ifdef HAVE_FC_INT_QUAD
    4078              :  CASE (16)
    4079            0 :    ABI_MALLOC(bufdelim16,(nb))
    4080            0 :    if (sc_mode==xmpio_single) then
    4081            0 :      call MPI_FILE_READ    (myfh,bufdelim16,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4082            0 :    else if (sc_mode==xmpio_collective) then
    4083            0 :      call MPI_FILE_READ_ALL(myfh,bufdelim16,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4084              :    else
    4085            0 :      ierr=2
    4086              :    end if
    4087            0 :    if (ANY(bufdelim16/=delim_record)) ierr=1
    4088            0 :    ABI_FREE(bufdelim16)
    4089              : #endif
    4090              : 
    4091              :  CASE (2)
    4092            0 :    ABI_MALLOC(bufdelim2,(nb))
    4093            0 :    if (sc_mode==xmpio_single) then
    4094            0 :      call MPI_FILE_READ    (myfh,bufdelim2,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4095            0 :    else if (sc_mode==xmpio_collective) then
    4096            0 :      call MPI_FILE_READ_ALL(myfh,bufdelim2,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4097              :    else
    4098            0 :      ierr=2
    4099              :    end if
    4100            0 :    if (ANY(bufdelim2/=delim_record)) ierr=1
    4101            0 :    ABI_FREE(bufdelim2)
    4102              : 
    4103              :  CASE DEFAULT
    4104            0 :    ierr=-2
    4105              :  END SELECT
    4106              : 
    4107              :  ! Free memory
    4108            0 :  call MPI_TYPE_FREE(frmarkers_type,mpierr)
    4109            0 :  ABI_FREE(delim_record)
    4110              : 
    4111            0 : end subroutine xmpio_check_frmarkers
    4112              : !!***
    4113              : #endif
    4114              : 
    4115              : !----------------------------------------------------------------------
    4116              : 
    4117              : !!****f* m_xmpi/xmpio_read_int
    4118              : !! NAME
    4119              : !!  xmpio_read_int
    4120              : !!
    4121              : !! FUNCTION
    4122              : !!  Read the content of a single record marker in a FORTRAN file at a given offset using MPI-IO.
    4123              : !!  the file pointer is modified according to the value of advance.
    4124              : !!  target: integer array
    4125              : !!
    4126              : !! INPUTS
    4127              : !!  fh=MPI-IO file handler.
    4128              : !!  offset=MPI-IO file pointer
    4129              : !!  sc_mode=
    4130              : !!         xmpio_single     ==> for reading by current proc.
    4131              : !!         xmpio_collective ==> for collective reading.
    4132              : !!  ncount=Number of elements in the buffer
    4133              : !!  [advance]=By default the routine will move the file pointer to the next record.
    4134              : !!    advance=.FALSE. can be used so that the next read will continue picking information
    4135              : !!    off of the currect record.
    4136              : !!
    4137              : !! OUTPUT
    4138              : !!  buf(ncount)=array with the values read from file
    4139              : !!  fmarker=Content of the Fortran record marker.
    4140              : !!  mpierr= MPI error code
    4141              : !!
    4142              : !! SIDE EFFECTS
    4143              : !!  offset=
    4144              : !!     input: file pointer used to access the Fortran marker.
    4145              : !!     output: new offset updated after the reading, depending on advance.
    4146              : !!
    4147              : !! SOURCE
    4148              : 
    4149              : #ifdef HAVE_MPI_IO
    4150              : 
    4151            0 : subroutine xmpio_read_int(fh, offset, sc_mode, ncount, buf, fmarker, mpierr, advance)
    4152              : 
    4153              : !Arguments ------------------------------------
    4154              : !scalars
    4155              :  integer,intent(in) :: fh,sc_mode,ncount
    4156              :  integer(XMPI_OFFSET_KIND),intent(inout) :: offset
    4157              :  integer(XMPI_OFFSET_KIND),intent(out) :: fmarker
    4158              :  integer,intent(out) :: mpierr
    4159              :  logical,optional,intent(in) :: advance
    4160              : !arrays
    4161              :  integer,intent(out) :: buf(ncount)
    4162              : 
    4163              : !Local variables-------------------------------
    4164              : !scalars
    4165              :  integer :: myfh,bsize_frm
    4166              :  integer(XMPI_OFFSET_KIND) :: my_offset
    4167              :  character(len=500) :: msg
    4168              : !arrays
    4169              :  integer :: statux(MPI_STATUS_SIZE)
    4170              : !************************************************************************
    4171              : 
    4172              :  ! Workaround for XLF
    4173            0 :  myfh = fh
    4174              : 
    4175            0 :  my_offset = offset
    4176            0 :  bsize_frm = xmpio_bsize_frm  ! Byte size of the Fortran record marker.
    4177              : 
    4178            0 :  call xmpio_read_frm(myfh,my_offset,sc_mode,fmarker,mpierr,advance=.FALSE.)
    4179              : 
    4180            0 :  SELECT CASE (sc_mode)
    4181              :  CASE (xmpio_single)
    4182            0 :    call MPI_FILE_READ_AT(myfh, my_offset, buf, ncount, MPI_INTEGER, statux, mpierr)
    4183              : 
    4184              :  CASE (xmpio_collective)
    4185            0 :    call MPI_FILE_READ_AT_ALL(myfh, my_offset, buf, ncount, MPI_INTEGER, statux, mpierr)
    4186              : 
    4187              :  CASE DEFAULT
    4188            0 :    write(msg,"(a,i0)")" Wrong value for sc_mode: ",sc_mode
    4189            0 :    call xmpi_abort(msg=msg)
    4190              :  END SELECT
    4191              : 
    4192            0 :  if (PRESENT(advance)) then
    4193            0 :    if (advance) then
    4194            0 :      offset = offset + fmarker + 2*bsize_frm ! Move the file pointer to the next record.
    4195              :    else
    4196            0 :      offset = offset + bsize_frm  ! Move the pointer after the marker.
    4197              :    end if
    4198              :  else
    4199            0 :    offset = offset + fmarker + 2*bsize_frm
    4200              :  end if
    4201              : 
    4202            0 : end subroutine xmpio_read_int
    4203              : !!***
    4204              : #endif
    4205              : 
    4206              : !----------------------------------------------------------------------
    4207              : 
    4208              : !!****f* m_xmpi/xmpio_read_dp
    4209              : !! NAME
    4210              : !!  xmpio_read_dp
    4211              : !!
    4212              : !! FUNCTION
    4213              : !!  Read the content of a single record marker in a FORTRAN file at a given offset using MPI-IO.
    4214              : !!  the file pointer is modified according to the value of advance.
    4215              : !!  targer: double precision real array
    4216              : !!
    4217              : !! INPUTS
    4218              : !!  fh=MPI-IO file handler.
    4219              : !!  offset=MPI-IO file pointer
    4220              : !!  sc_mode=
    4221              : !!         xmpio_single     ==> for reading by current proc.
    4222              : !!         xmpio_collective ==> for collective reading.
    4223              : !!  ncount=Number of elements in the buffer
    4224              : !!  [advance]=By default the routine will move the file pointer to the next record.
    4225              : !!    advance=.FALSE. can be used so that the next read will continue picking information
    4226              : !!    off of the currect record.
    4227              : !!
    4228              : !! OUTPUT
    4229              : !!  buf(ncount)=array with the values read from file
    4230              : !!  fmarker=Content of the Fortran record marker.
    4231              : !!  mpierr= MPI error code
    4232              : !!
    4233              : !! SIDE EFFECTS
    4234              : !!  offset=
    4235              : !!     input: file pointer used to access the Fortran marker.
    4236              : !!     output: new offset updated after the reading, depending on advance.
    4237              : !!
    4238              : !! SOURCE
    4239              : 
    4240              : #ifdef HAVE_MPI_IO
    4241              : 
    4242            0 : subroutine xmpio_read_dp(fh, offset, sc_mode, ncount, buf, fmarker, mpierr, advance)
    4243              : 
    4244              : !Arguments ------------------------------------
    4245              : !scalars
    4246              :  integer,intent(in) :: fh,sc_mode,ncount
    4247              :  integer(XMPI_OFFSET_KIND),intent(inout) :: offset
    4248              :  integer(XMPI_OFFSET_KIND),intent(out) :: fmarker
    4249              :  integer,intent(out) :: mpierr
    4250              :  logical,optional,intent(in) :: advance
    4251              : !arrays
    4252              :  real(dp),intent(out) :: buf(ncount)
    4253              : 
    4254              : !Local variables-------------------------------
    4255              : !scalars
    4256              :  integer :: bsize_frm,myfh
    4257              :  integer(XMPI_OFFSET_KIND) :: my_offset
    4258              :  character(len=500) :: msg
    4259              : !arrays
    4260              :  integer :: statux(MPI_STATUS_SIZE)
    4261              : !************************************************************************
    4262              : 
    4263              :  ! Workaround for XLF
    4264            0 :  myfh = fh
    4265              : 
    4266            0 :  my_offset = offset
    4267            0 :  bsize_frm = xmpio_bsize_frm  ! Byte size of the Fortran record marker.
    4268              : 
    4269            0 :  call xmpio_read_frm(myfh,my_offset,sc_mode,fmarker,mpierr,advance=.FALSE.)
    4270              : 
    4271            0 :  SELECT CASE (sc_mode)
    4272              :  CASE (xmpio_single)
    4273            0 :    call MPI_FILE_READ_AT(myfh, my_offset, buf, ncount, MPI_DOUBLE_PRECISION, statux, mpierr)
    4274              : 
    4275              :  CASE (xmpio_collective)
    4276            0 :    call MPI_FILE_READ_AT_ALL(myfh, my_offset, buf, ncount, MPI_DOUBLE_PRECISION, statux, mpierr)
    4277              : 
    4278              :  CASE DEFAULT
    4279            0 :    write(msg,"(a,i0)")" Wrong value for sc_mode: ",sc_mode
    4280            0 :    call xmpi_abort(msg=msg)
    4281              :  END SELECT
    4282              : 
    4283            0 :  if (PRESENT(advance)) then
    4284            0 :    if (advance) then
    4285            0 :      offset = offset + fmarker + 2*bsize_frm ! Move the file pointer to the next record.
    4286              :    else
    4287            0 :      offset = offset + bsize_frm  ! Move the pointer after the marker.
    4288              :    end if
    4289              :  else
    4290            0 :    offset = offset + fmarker + 2*bsize_frm
    4291              :  end if
    4292              : 
    4293            0 : end subroutine xmpio_read_dp
    4294              : !!***
    4295              : #endif
    4296              : 
    4297              : !------------------------------------------------------------------------------------
    4298              : 
    4299              : !!****f* m_xmpi/xmpio_max_address
    4300              : !! NAME
    4301              : !!  xmpio_max_address
    4302              : !!
    4303              : !! FUNCTION
    4304              : !!  Returns .TRUE. if offset cannot be stored in a Fortran integer of kind XMPI_ADDRESS_KIND.
    4305              : !!
    4306              : !! SOURCE
    4307              : 
    4308              : #ifdef HAVE_MPI_IO
    4309              : 
    4310            0 : function xmpio_max_address(offset)
    4311              : 
    4312              : !Arguments ------------------------------------
    4313              : !scalars
    4314              :  logical :: xmpio_max_address
    4315              :  integer(XMPI_OFFSET_KIND),intent(in) :: offset
    4316              : !arrays
    4317              : 
    4318              : !Local variables-------------------------------
    4319              : !scalars
    4320              :  integer(XMPI_ADDRESS_KIND) :: address
    4321              :  integer(XMPI_OFFSET_KIND),parameter :: max_address=HUGE(address)-100
    4322              : !************************************************************************
    4323              : 
    4324            0 :  xmpio_max_address = (offset >= max_address)
    4325              : 
    4326            0 : end function xmpio_max_address
    4327              : !!***
    4328              : #endif
    4329              : 
    4330              : !------------------------------------------------------------------------------------
    4331              : 
    4332              : !!****f* m_xmpi/xmpio_write_frmarkers
    4333              : !! NAME
    4334              : !!  xmpio_write_frmarkers
    4335              : !!
    4336              : !! FUNCTION
    4337              : !!  Write a set of Fortran record markers starting at a given offset using MPI-IO.
    4338              : !!
    4339              : !! INPUTS
    4340              : !!  fh=MPI-IO file handler.
    4341              : !!  offset=MPI-IO file pointer
    4342              : !!  sc_mode=Option for individual or collective reading.
    4343              : !!  nfrec=Number of Fortran records to be written.
    4344              : !!  bsize_frecord(nfrec)=Byte size of the Fortran records to be written (markers are NOT included in the size)
    4345              : !!
    4346              : !! OUTPUT
    4347              : !!  ierr=A non-zero error code signals failure.
    4348              : !!
    4349              : !! SOURCE
    4350              : 
    4351              : #ifdef HAVE_MPI_IO
    4352              : 
    4353            5 : subroutine xmpio_write_frmarkers(fh, offset, sc_mode, nfrec, bsize_frecord, ierr)
    4354              : 
    4355              : !Arguments ------------------------------------
    4356              : !scalars
    4357              :  integer,intent(in) :: fh,nfrec,sc_mode
    4358              :  integer(XMPI_OFFSET_KIND),intent(in) :: offset
    4359              :  integer,intent(out) :: ierr
    4360              : !arrays
    4361              :  integer(XMPI_OFFSET_KIND),intent(in) :: bsize_frecord(nfrec)
    4362              : 
    4363              : !Local variables-------------------------------
    4364              : !scalars
    4365              :  integer :: nb,irec,frmarkers_type,jj,bsize_frm,mpi_type_frm,mpierr,myfh
    4366              :  integer(XMPI_OFFSET_KIND) :: displ
    4367              : !integer(XMPI_OFFSET_KIND) :: my_offset
    4368              : !character(len=500) :: msg
    4369              : !arrays
    4370            5 :  integer(kind=int16),allocatable :: bufdelim2(:)
    4371            5 :  integer(kind=int32),allocatable :: bufdelim4(:)
    4372            5 :  integer(kind=int64),allocatable :: bufdelim8(:)
    4373              : #ifdef HAVE_FC_INT_QUAD
    4374            5 :  integer*16,allocatable :: bufdelim16(:)
    4375              : #endif
    4376              : !integer :: statux(MPI_STATUS_SIZE)
    4377            5 :  integer,allocatable :: block_length(:),block_type(:)
    4378            5 :  integer(XMPI_ADDRESS_KIND),allocatable :: block_displ(:)
    4379            5 :  integer(XMPI_OFFSET_KIND),allocatable :: delim_record(:)
    4380              : !************************************************************************
    4381              : 
    4382              :  ! Workaround for XLF
    4383            5 :  myfh = fh; ierr=0
    4384              : 
    4385              :  !my_offset = offset
    4386              :  !do irec=1,nfrec
    4387              :  !  call xmpio_write_frm(myfh,my_offset,sc_mode,bsize_frecord(irec),mpierr)
    4388              :  !end do
    4389              :  !return
    4390              : 
    4391              :  ! FIXME: This is buggy
    4392            5 :  bsize_frm    = xmpio_bsize_frm     ! Byte size of the Fortran record marker.
    4393            5 :  mpi_type_frm = xmpio_mpi_type_frm  ! MPI type of the record marker.
    4394              : 
    4395              :  ! Define the view for the file
    4396            5 :  nb=2*nfrec
    4397           15 :  ABI_MALLOC(block_length,(nb+2))
    4398           15 :  ABI_MALLOC(block_displ,(nb+2))
    4399           10 :  ABI_MALLOC(block_type,(nb+2))
    4400            5 :  block_length(1)=1
    4401            5 :  block_displ (1)=0
    4402            5 :  block_type  (1)=MPI_LB
    4403              : 
    4404            5 :  jj=2; displ=0
    4405          335 :  do irec=1,nfrec
    4406          990 :    block_type (jj:jj+1)  = mpi_type_frm
    4407          990 :    block_length(jj:jj+1) = 1
    4408          330 :    block_displ(jj  )     = displ
    4409          330 :    block_displ(jj+1)     = displ + bsize_frm + bsize_frecord(irec)
    4410          330 :    jj=jj+2
    4411          330 :    displ = displ + bsize_frecord(irec) + 2*bsize_frm ! Move to the beginning of the next column.
    4412          335 :    if (xmpio_max_address(displ)) then ! Check for wraparound.
    4413            0 :       ierr = -1; return
    4414              :    end if
    4415              :  end do
    4416              : 
    4417            5 :  block_length(nb+2) = 1
    4418            5 :  block_displ (nb+2) = displ
    4419            5 :  block_type  (nb+2) = MPI_UB
    4420              : 
    4421            5 :  call xmpio_type_struct(nb+2,block_length,block_displ,block_type,frmarkers_type,mpierr)
    4422              : 
    4423            5 :  ABI_FREE(block_length)
    4424            5 :  ABI_FREE(block_displ)
    4425            5 :  ABI_FREE(block_type)
    4426              : 
    4427            5 :  call MPI_TYPE_COMMIT(frmarkers_type,mpierr)
    4428            5 :  call MPI_FILE_SET_VIEW(myfh,offset,MPI_BYTE,frmarkers_type,"native",MPI_INFO_NULL,mpierr)
    4429              : 
    4430            5 :  jj=1
    4431           15 :  ABI_MALLOC(delim_record,(nb))
    4432          335 :  do irec=1,nfrec
    4433          990 :    delim_record(jj:jj+1)=bsize_frecord(irec)
    4434          335 :    jj=jj+2
    4435              :  end do
    4436              : 
    4437              :  ! Write all markers according to the MPI type of the Fortran marker.
    4438            5 :  SELECT CASE (bsize_frm)
    4439              : 
    4440              :  CASE (4)
    4441           15 :    ABI_MALLOC(bufdelim4,(nb))
    4442          670 :    bufdelim4=delim_record
    4443            5 :    if (sc_mode==xmpio_single) then
    4444            0 :      call MPI_FILE_WRITE    (myfh,bufdelim4,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4445            5 :    else if (sc_mode==xmpio_collective) then
    4446            5 :      call MPI_FILE_WRITE_ALL(myfh,bufdelim4,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4447              :    else
    4448            0 :      ierr=2
    4449              :    end if
    4450            5 :    ABI_FREE(bufdelim4)
    4451              : 
    4452              :  CASE (8)
    4453            0 :    ABI_MALLOC(bufdelim8,(nb))
    4454            0 :    bufdelim8=delim_record
    4455            0 :    if (sc_mode==xmpio_single) then
    4456            0 :      call MPI_FILE_WRITE    (myfh,bufdelim8,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4457            0 :    else if (sc_mode==xmpio_collective) then
    4458            0 :      call MPI_FILE_WRITE_ALL(myfh,bufdelim8,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4459              :    else
    4460            0 :      ierr=2
    4461              :    end if
    4462            0 :    ABI_FREE(bufdelim8)
    4463              : 
    4464              : #ifdef HAVE_FC_INT_QUAD
    4465              :  CASE (16)
    4466            0 :    ABI_MALLOC(bufdelim16,(nb))
    4467            0 :    bufdelim16=delim_record
    4468            0 :    if (sc_mode==xmpio_single) then
    4469            0 :      call MPI_FILE_WRITE    (myfh,bufdelim16,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4470            0 :    else if (sc_mode==xmpio_collective) then
    4471            0 :      call MPI_FILE_WRITE_ALL(myfh,bufdelim16,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4472              :    else
    4473            0 :      ierr=2
    4474              :    end if
    4475            0 :    ABI_FREE(bufdelim16)
    4476              : #endif
    4477              : 
    4478              :  CASE (2)
    4479            0 :    ABI_MALLOC(bufdelim2,(nb))
    4480            0 :    bufdelim2=delim_record
    4481            0 :    if (sc_mode==xmpio_single) then
    4482            0 :      call MPI_FILE_WRITE    (myfh,bufdelim2,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4483            0 :    else if (sc_mode==xmpio_collective) then
    4484            0 :      call MPI_FILE_WRITE_ALL(myfh,bufdelim2,2*nfrec,mpi_type_frm,MPI_STATUS_IGNORE,mpierr)
    4485              :    else
    4486            0 :      ierr=2
    4487              :    end if
    4488            0 :    ABI_FREE(bufdelim2)
    4489              : 
    4490              :  CASE DEFAULT
    4491            5 :    ierr=-2
    4492              :  END SELECT
    4493              : 
    4494              :  ! Free memory
    4495            5 :  call MPI_TYPE_FREE(frmarkers_type,mpierr)
    4496            5 :  ABI_FREE(delim_record)
    4497              : 
    4498            5 : end subroutine xmpio_write_frmarkers
    4499              : #endif
    4500              : !!***
    4501              : 
    4502              : !------------------------------------------------------------------------------------
    4503              : 
    4504              : !!****f* m_xmpi/xmpio_create_fherm_packed
    4505              : !! NAME
    4506              : !!  xmpio_create_fherm_packed
    4507              : !!
    4508              : !! FUNCTION
    4509              : !!  Returns an MPI datatype that can be used to (read|write) with MPI-IO the columns of an
    4510              : !!  Hermitian matrix whose upper triangle is written on a Fortran binary file.
    4511              : !!  Note that the view assumes that the file pointer used to create the MPI-IO view
    4512              : !!  points to the first element of the first column. In other words,the first Fortran record marker
    4513              : !!  (if any) is not taken into account in the calculation of the displacements.
    4514              : !!
    4515              : !! INPUTS
    4516              : !!  array_of_starts(2)=starting coordinates in the global Hermitian matrix
    4517              : !!     (array of positive integers with jj>=ii, Fortran convention)
    4518              : !!  array_of_ends(2)=final coordinates in the global Hermitian matrix
    4519              : !!     (array of positive integers, jj>=ii, Fortran convention)
    4520              : !!  is_fortran_file=.FALSE. is C stream is used. .TRUE. for writing Fortran binary files.
    4521              : !!  old_type=MPI datatype of the elements of the matrix.
    4522              : !!
    4523              : !! OUTPUT
    4524              : !!  my_offset=Offset relative to the beginning of the matrix in the file.
    4525              : !!  hmat_type=New MPI type.
    4526              : !!  offset_err= error code
    4527              : !!
    4528              : !! NOTES
    4529              : !!  The matrix on file is written in the following FORTRAN format (let us assume a 3x3 matrix for simplicity)
    4530              : !!
    4531              : !!    m (1,1)             m
    4532              : !!    m (1,2) (2,2)       m
    4533              : !!    m (1,3) (2,3) (3,3) m
    4534              : !!
    4535              : !!  each Fortran record stores a column of the packed Hermitian matrix, "m" denotes the Fortran
    4536              : !!  record marker that introduces holes in the MPI-IO file view.
    4537              : !!  To read the columns from (1,2) up to (2,2) one should use array_of_starts=(1,2) and array_of_ends=(2,2).
    4538              : !!  The MPI-IO file view should be created by moving the file pointer so that it points to the elements (1,2).
    4539              : !!
    4540              : !!  File views for C-streams is not optimal since one can use a single slice of contigous data.
    4541              : !!
    4542              : !! SOURCE
    4543              : 
    4544              : #ifdef HAVE_MPI_IO
    4545              : 
    4546            0 : subroutine xmpio_create_fherm_packed(array_of_starts,array_of_ends,is_fortran_file,my_offset,old_type,hmat_type,offset_err)
    4547              : 
    4548              : !Arguments ------------------------------------
    4549              : !scalars
    4550              :  integer,intent(in) :: old_type
    4551              :  integer,intent(out) :: offset_err,hmat_type
    4552              :  integer(XMPI_OFFSET_KIND),intent(out) :: my_offset
    4553              :  logical,intent(in) :: is_fortran_file
    4554              : !arrays
    4555              :  integer,intent(in) :: array_of_starts(2),array_of_ends(2)
    4556              : 
    4557              : !Local variables-------------------------------
    4558              : !scalars
    4559              :  integer :: nrow,my_ncol,ii,bsize_old,col,jj_glob,bsize_frm,prev_col,mpierr
    4560              :  integer(XMPI_OFFSET_KIND) :: col_displ
    4561              : !arrays
    4562            0 :  integer,allocatable :: col_type(:),block_length(:),block_type(:)
    4563            0 :  integer(XMPI_ADDRESS_KIND),allocatable :: block_displ(:)
    4564              : !************************************************************************
    4565              : 
    4566            0 :  offset_err=0
    4567              : 
    4568              :  ! Byte size of old_type.
    4569            0 :  call MPI_TYPE_SIZE(old_type,bsize_old,mpierr)
    4570              : 
    4571            0 :  bsize_frm=0; if (is_fortran_file) bsize_frm = xmpio_bsize_frm
    4572              : 
    4573            0 :  my_ncol = array_of_ends(2) - array_of_starts(2) + 1
    4574              :  !
    4575              :  ! Calculate my offset relative to the beginning of the matrix in the file.
    4576            0 :  prev_col = array_of_starts(2)-1
    4577            0 :  my_offset = (prev_col*(prev_col+1)/2)*bsize_old + (array_of_starts(1)-1)*bsize_old + 2*prev_col*bsize_frm + bsize_frm
    4578              :  !
    4579              :  ! col_type(col) describes the col-th column of the packed matrix.
    4580              :  ! block_displ(col+1) stores its displacement taking into account the Fortran marker.
    4581            0 :  ABI_MALLOC(col_type,(my_ncol))
    4582            0 :  ABI_MALLOC(block_displ,(my_ncol+2))
    4583              : 
    4584            0 :  if (my_ncol>1) then
    4585              :    col_displ=0
    4586            0 :    do col=1,my_ncol
    4587            0 :     jj_glob = (col-1) + array_of_starts(2)
    4588            0 :     nrow = jj_glob
    4589            0 :     if (jj_glob==array_of_starts(2)) nrow = jj_glob - array_of_starts(1) + 1 ! First column treated by me.
    4590            0 :     if (jj_glob==array_of_ends(2))   nrow = array_of_ends(1)                 ! Last column treated by me.
    4591            0 :     call MPI_Type_contiguous(nrow,old_type,col_type(col),mpierr)
    4592              :     !
    4593            0 :     if (xmpio_max_address(col_displ)) offset_err=1  ! Test for wraparounds
    4594            0 :     block_displ(col+1) = col_displ
    4595            0 :     col_displ = col_displ + nrow * bsize_old + 2 * bsize_frm  ! Move to the next column.
    4596              :    end do
    4597              : 
    4598            0 :  else if (my_ncol==1) then  ! The case of a single column is treated separately.
    4599            0 :     block_displ(2) = 0
    4600            0 :     nrow = array_of_ends(1) - array_of_starts(1) + 1
    4601            0 :     call MPI_Type_contiguous(nrow,old_type,col_type(2),mpierr)
    4602            0 :     col_displ= nrow*bsize_old
    4603            0 :     if (xmpio_max_address(col_displ)) offset_err=1  ! Test for wraparounds
    4604              :  else
    4605            0 :    call xmpi_abort(msg="my_ncol cannot be negative!")
    4606              :  end if
    4607              : 
    4608            0 :  ABI_MALLOC(block_length,(my_ncol+2))
    4609            0 :  ABI_MALLOC(block_type,(my_ncol+2))
    4610              : 
    4611            0 :  block_length(1)=1
    4612            0 :  block_displ (1)=0
    4613            0 :  block_type  (1)=MPI_LB
    4614              : 
    4615            0 :  do ii=2,my_ncol+1
    4616            0 :    block_length(ii)=1
    4617            0 :    block_type(ii)  =col_type(ii-1)
    4618              :    !write(std_out,*)" ii-1, depl, length, type: ",ii-1,block_displ(ii),block_length(ii),block_type(ii)
    4619              :  end do
    4620              : 
    4621            0 :  block_length(my_ncol+2)= 1
    4622            0 :  block_displ (my_ncol+2)= col_displ
    4623            0 :  block_type  (my_ncol+2)= MPI_UB
    4624              : 
    4625            0 :  call xmpio_type_struct(my_ncol+2,block_length,block_displ,block_type,hmat_type,mpierr)
    4626              : 
    4627            0 :  call MPI_TYPE_COMMIT(hmat_type,mpierr)
    4628              : 
    4629            0 :  ABI_FREE(block_length)
    4630            0 :  ABI_FREE(block_displ)
    4631            0 :  ABI_FREE(block_type)
    4632              : 
    4633            0 :  do col=1,my_ncol
    4634            0 :    call MPI_TYPE_FREE(col_type(col),mpierr)
    4635              :  end do
    4636              : 
    4637            0 :  ABI_FREE(col_type)
    4638              : 
    4639            0 : end subroutine xmpio_create_fherm_packed
    4640              : !!***
    4641              : #endif
    4642              : 
    4643              : !------------------------------------------------------------------------------------
    4644              : 
    4645              : !!****f* m_xmpi/xmpio_create_coldistr_from_fpacked
    4646              : !! NAME
    4647              : !!  xmpio_create_coldistr_from_fpacked
    4648              : !!
    4649              : !! FUNCTION
    4650              : !!  Returns an MPI datatype that can be used to MPI-IO (read|write) the columns of an
    4651              : !!  (Hermitian|Symmetric) matrix whose upper triangle is written on a Fortran binary file.
    4652              : !!  Note that the view assumes that the file pointer used to instanciate the MPI-IO view
    4653              : !!  points to the first element of the first column. In other words,the first Fortran record marker
    4654              : !!  (if any) is not taken into account in the calculation of the displacements.
    4655              : !!
    4656              : !! INPUTS
    4657              : !!  sizes(2)=Number of elements of type old_type in each dimension of the full array (array of positive integers)
    4658              : !!  my_cols(2)=initial and final column to (read|write). Array of positive integers, Fortran convention.
    4659              : !!  old_type=MPI datatype of the elements of the matrix.
    4660              : !!
    4661              : !! OUTPUT
    4662              : !!  new_type=New MPI type that can be used to instanciate the MPI-IO view for the Fortran file.
    4663              : !!  my_offpad=Offset to be added to the file pointer giving the position of the first Fortran record
    4664              : !!    marker (lets call it "base"). Each node should (read|write) using my_offset = base + my_offpad.
    4665              : !!    my_offpad is used so that one can safely change the way the fileview is generated (for example
    4666              : !!    to make it more efficient) without having to change the client code.
    4667              : !!  offset_err=Error code. A non-zero returned value signals that the global matrix is tool large
    4668              : !!    for a single MPI-IO access (see notes below).
    4669              : !!
    4670              : !! NOTES
    4671              : !!  1) The matrix on file is written in the following FORTRAN format (let us assume a 3x3 matrix for simplicity)
    4672              : !!
    4673              : !!      m (1,1)             m
    4674              : !!      m (1,2) (2,2)       m
    4675              : !!      m (1,3) (2,3) (3,3) m
    4676              : !!
    4677              : !!     each Fortran record stores a column of the packed matrix, "m" denotes the Fortran
    4678              : !!     record marker that introduces holes in the file view.
    4679              : !!
    4680              : !!  2) With (signed) Fortran integers, the maximum size of the file that
    4681              : !!     that can be read in one-shot is around 2Gb when etype is set to byte.
    4682              : !!     Using a larger etype might create portability problems (real data on machines using
    4683              : !!     integer*16 for the marker) since etype must be a multiple of the Fortran record marker
    4684              : !!     Due to the above reason, block_displ is given in bytes but it has to be defined as Fortran
    4685              : !!     integer. If the displacement cannot be stored in a Fortran integer, the routine returns
    4686              : !!     offset_err=1 so that the caller will know that several MPI-IO reads are nedded to
    4687              : !!     read the file.
    4688              : !!
    4689              : !! SOURCE
    4690              : 
    4691              : #ifdef HAVE_MPI_IO
    4692              : 
    4693            0 : subroutine xmpio_create_coldistr_from_fpacked(sizes,my_cols,old_type,new_type,my_offpad,offset_err)
    4694              : 
    4695              : !Arguments ------------------------------------
    4696              : !scalars
    4697              :  integer,intent(in) :: old_type
    4698              :  integer,intent(out) :: new_type,offset_err
    4699              :  integer(XMPI_OFFSET_KIND),intent(out) :: my_offpad
    4700              : !arrays
    4701              :  integer,intent(in) :: sizes(2),my_cols(2)
    4702              : 
    4703              : !Local variables-------------------------------
    4704              : !scalars
    4705              :  integer :: my_ncol,bsize_old,my_col
    4706              :  integer :: my_nels,my_el,row_glob,ii_hpk,jj_hpk,col_glob,bsize_frm,mpierr
    4707              :  integer(XMPI_OFFSET_KIND) :: my_offset,ijp_glob
    4708              :  !character(len=500) :: msg
    4709              : !arrays
    4710            0 :  integer,allocatable :: block_length(:),block_type(:)
    4711            0 :  integer(XMPI_ADDRESS_KIND),allocatable :: block_displ(:)
    4712              : !************************************************************************
    4713              : 
    4714              :  ! Byte size of the Fortran record marker.
    4715            0 :  bsize_frm = xmpio_bsize_frm
    4716              : 
    4717              :  ! Byte size of old_type.
    4718            0 :  call MPI_TYPE_SIZE(old_type,bsize_old,mpierr)
    4719              : 
    4720              :  ! my number of columns and total numer of elements to be read.
    4721            0 :  my_ncol = my_cols(2) - my_cols(1) + 1
    4722            0 :  my_nels = my_ncol*sizes(1)
    4723              :  !
    4724              :  ! block_displ(el+1) stores the displacement of the local element el taking into account the Fortran marker.
    4725            0 :  ABI_MALLOC(block_displ,(my_nels+2))
    4726            0 :  ABI_MALLOC(block_length,(my_nels+2))
    4727            0 :  ABI_MALLOC(block_type,(my_nels+2))
    4728              : 
    4729            0 :  block_length(1)=1
    4730            0 :  block_displ (1)=0
    4731            0 :  block_type  (1)=MPI_LB
    4732              :  !
    4733              :  ! * the view assumes that the file pointer used to instanciate the MPI-IO view
    4734              :  !   points to the first element of the first column. In other words,the first Fortran record marker
    4735              :  !   is not taken into account in the calculation of the displacements.
    4736            0 :  my_offpad=xmpio_bsize_frm
    4737              : 
    4738              :  ! * Some matrix elements are read twice. This part has to be tested.
    4739            0 :  offset_err=0; my_el=0
    4740            0 :  do my_col=1,my_ncol
    4741            0 :    col_glob = (my_col-1) + my_cols(1)
    4742            0 :    do row_glob=1,sizes(1)
    4743            0 :      if (col_glob>=row_glob) then
    4744            0 :        ii_hpk = row_glob
    4745            0 :        jj_hpk = col_glob
    4746            0 :        ijp_glob = row_glob + col_glob*(col_glob-1)/2  ! Index for packed form
    4747              :      else ! Exchange the indices as (jj,ii) will be read.
    4748            0 :        ii_hpk = col_glob
    4749            0 :        jj_hpk = row_glob
    4750            0 :        ijp_glob = col_glob + row_glob*(row_glob-1)/2  ! Index for packed form
    4751              :      end if
    4752            0 :      my_el = my_el+1
    4753            0 :      my_offset = (ijp_glob-1)* bsize_old + (jj_hpk-1)*2*bsize_frm
    4754            0 :      if (xmpio_max_address(my_offset)) offset_err=1   ! Check for wraparounds.
    4755            0 :      block_displ (my_el+1)=my_offset
    4756            0 :      block_length(my_el+1)=1
    4757            0 :      block_type  (my_el+1)=old_type
    4758              :      !write(std_out,*)" my_el, displ: ",my_el,block_displ(my_el+1)
    4759              :    end do
    4760              :  end do
    4761              : 
    4762            0 :  block_length(my_nels+2)=1
    4763            0 :  block_displ (my_nels+2)=my_offset
    4764            0 :  block_type  (my_nels+2)=MPI_UB
    4765              : 
    4766            0 :  call xmpio_type_struct(my_nels+2,block_length,block_displ,block_type,new_type,mpierr)
    4767              : 
    4768            0 :  call MPI_TYPE_COMMIT(new_type,mpierr)
    4769              : 
    4770            0 :  ABI_FREE(block_length)
    4771            0 :  ABI_FREE(block_displ)
    4772            0 :  ABI_FREE(block_type)
    4773              : 
    4774            0 : end subroutine xmpio_create_coldistr_from_fpacked
    4775              : !!***
    4776              : #endif
    4777              : 
    4778              : !------------------------------------------------------------------------------------
    4779              : 
    4780              : !!****f* m_xmpi/xmpio_create_coldistr_from_fp3blocks
    4781              : !! NAME
    4782              : !!  xmpio_create_coldistr_from_fp3blocks
    4783              : !!
    4784              : !! FUNCTION
    4785              : !!  Returns an MPI datatype that can be used to MPI-IO (read|write) the columns of a
    4786              : !!  matrix of the form  M = (S1    F3)
    4787              : !!                          (F3^H  S2)
    4788              : !!  where S1 and S2 are square (symmetric|Hermitian) matrices whose upper triangle is stored on file
    4789              : !!  while F3 is a generic matrix (not necessarily square) stored in full mode.
    4790              : !!  The Fortran file contains the blocks in the following order.
    4791              : !!      upper(S1)
    4792              : !!      upper(S2)
    4793              : !!      F3
    4794              : !! INPUTS
    4795              : !!  sizes(2)=Number of elements of type old_type in each dimension of the full array M (array of positive integers)
    4796              : !!  my_cols(2)=initial and final column to (read|write). Array of positive integers, Fortran convention.
    4797              : !!  block_sizes(2,3)=The sizes of S1, S2, F.
    4798              : !!  old_type=MPI datatype of the elements of the matrix.
    4799              : !!
    4800              : !! OUTPUT
    4801              : !!  new_type=New MPI type that can be used to instanciate the MPI-IO view for the Fortran file.
    4802              : !!  my_offpad=Offset to be added to the file pointer giving the position of the first Fortran record
    4803              : !!    marker (lets call it "base"). Each node should (read|write) using my_offset = base + my_offpad.
    4804              : !!    my_offpad is used so that one can safely change the way the fileview is generated (for example
    4805              : !!    to make it more efficient) without having to change the client code.
    4806              : !!  offset_err=Error code. A non-zero returned value signals that the global matrix is tool large
    4807              : !!    for a single MPI-IO access (see notes below).
    4808              : !!
    4809              : !! NOTES
    4810              : !!  1) block_displ is given in bytes due to the presence of the marker.
    4811              : !!     If the displacement of an element is too large, the routine returns
    4812              : !!     offset_err=1 so that the caller knows that several MPI-IO reads are required to (read| write) the file.
    4813              : !!
    4814              : !! SOURCE
    4815              : 
    4816              : #ifdef HAVE_MPI_IO
    4817              : 
    4818            0 : subroutine xmpio_create_coldistr_from_fp3blocks(sizes,block_sizes,my_cols,old_type,new_type,my_offpad,offset_err)
    4819              : 
    4820              : !Arguments ------------------------------------
    4821              : !scalars
    4822              :  integer,intent(in) :: old_type
    4823              :  integer,intent(out) :: new_type,offset_err
    4824              :  integer(XMPI_OFFSET_KIND),intent(out) :: my_offpad
    4825              : !arrays
    4826              :  integer,intent(in) :: sizes(2),my_cols(2),block_sizes(2,3)
    4827              : 
    4828              : !Local variables-------------------------------
    4829              : !scalars
    4830              :  integer :: my_ncol,bsize_old,my_col,which_block,uplo,swap
    4831              :  integer :: my_nels,my_el,row_glob,ii_hpk,jj_hpk,ii,jj
    4832              :  integer :: col_glob,bsize_frm,mpierr,row_shift,col_shift,n1,n2
    4833              :  integer(XMPI_OFFSET_KIND) :: my_offset,ijp,bsize_tot,max_displ,min_displ
    4834              :  integer(XMPI_ADDRESS_KIND) :: address
    4835              : !arrays
    4836            0 :  integer,allocatable :: block_length(:),block_type(:)
    4837            0 :  integer(XMPI_ADDRESS_KIND),allocatable :: block_displ(:)
    4838              :  integer(XMPI_OFFSET_KIND) :: bsize_mat(2)
    4839              : !************************************************************************
    4840              : 
    4841            0 :  if (sizes(1) /= SUM(block_sizes(1,1:2)) .or. &
    4842              :      sizes(2) /= SUM(block_sizes(2,1:2)) ) then
    4843            0 :    write(std_out,*)" xmpio_create_coldistr_from_fp3blocks: Inconsistency between block_sizes ans sizes "
    4844            0 :    call xmpi_abort()
    4845              :  end if
    4846              : 
    4847            0 :  if (block_sizes(1,1) /= block_sizes(2,1) .or.&
    4848              :      block_sizes(1,2) /= block_sizes(2,2) ) then
    4849            0 :    write(std_out,*)" xmpio_create_coldistr_from_fp3blocks: first two blocks must be square"
    4850            0 :    call xmpi_abort()
    4851              :  end if
    4852              : 
    4853            0 :  if (block_sizes(2,3) /= block_sizes(2,2) .or.&
    4854              :      block_sizes(1,3) /= block_sizes(1,1) ) then
    4855            0 :    write(std_out,*)" xmpio_create_coldistr_from_fp3blocks: Full matrix must be square"
    4856            0 :    call xmpi_abort()
    4857              :  end if
    4858              : 
    4859            0 :  write(std_out,*)" xmpio_create_coldistr_from_fp3blocks is still under testing"
    4860              :  !call xmpi_abort()
    4861              : 
    4862              :  ! Byte size of the Fortran record marker.
    4863            0 :  bsize_frm = xmpio_bsize_frm
    4864              : 
    4865              :  ! Byte size of old_type.
    4866            0 :  call MPI_TYPE_SIZE(old_type,bsize_old,mpierr)
    4867              : 
    4868              :  ! my number of columns and total numer of elements to be read.
    4869            0 :  my_ncol = my_cols(2) - my_cols(1) + 1
    4870            0 :  my_nels = sizes(1)*my_ncol
    4871              :  !
    4872              :  ! block_displ(el+1) stores the displacement of the local element el taking into account the Fortran marker.
    4873            0 :  ABI_MALLOC(block_displ,(my_nels+2))
    4874            0 :  ABI_MALLOC(block_length,(my_nels+2))
    4875            0 :  ABI_MALLOC(block_type,(my_nels+2))
    4876              :  !
    4877              :  ! * the view assumes that the file pointer used to instanciate the MPI-IO view
    4878              :  !   points to the first element of the first column. In other words,the first Fortran record marker
    4879              :  !   is not taken into account in the calculation of the displacements.
    4880            0 :  my_offpad=xmpio_bsize_frm
    4881              :  !
    4882              :  ! Byte size of the first two blocks including the markers.
    4883            0 :  n1=block_sizes(1,1)
    4884            0 :  bsize_mat(1) = (n1*(n1+1)/2)*bsize_old + 2*n1*bsize_frm
    4885              : 
    4886            0 :  n2=block_sizes(1,2)
    4887            0 :  bsize_mat(2) = (n2*(n2+1)/2)*bsize_old + 2*n2*bsize_frm
    4888              : 
    4889            0 :  bsize_tot=SUM(bsize_mat) +  PRODUCT(block_sizes(:,3))*bsize_old + block_sizes(2,3)*2*bsize_frm - bsize_frm
    4890            0 :  write(std_out,*)"bsize_mat",bsize_mat,"bsize_tot",bsize_tot
    4891              :  !
    4892              :  ! * Some matrix elements are read twice. This part has to be tested.
    4893            0 :  offset_err=0; my_el=0; max_displ=0; min_displ=HUGE(address)
    4894            0 :  do my_col=1,my_ncol
    4895            0 :    col_glob = (my_col-1) + my_cols(1)
    4896            0 :    do row_glob=1,sizes(1)
    4897              :      !
    4898            0 :      which_block=3
    4899            0 :      if (row_glob<=block_sizes(1,1).and.col_glob<=block_sizes(2,1)) which_block=1
    4900            0 :      if (row_glob >block_sizes(1,1).and.col_glob >block_sizes(2,1)) which_block=2
    4901              : 
    4902            0 :      if ( ANY(which_block == (/1,2/)) ) then ! S1 or S2
    4903              :        !
    4904            0 :        row_shift=(which_block-1)*block_sizes(1,1)
    4905            0 :        col_shift=(which_block-1)*block_sizes(2,1)
    4906              : 
    4907            0 :        ii_hpk = row_glob - row_shift
    4908            0 :        jj_hpk = col_glob - col_shift
    4909            0 :        if (jj_hpk<ii_hpk) then ! Exchange the indices so that the symmetric is read.
    4910            0 :          swap   = jj_hpk
    4911            0 :          jj_hpk = ii_hpk
    4912            0 :          ii_hpk = swap
    4913              :        end if
    4914            0 :        ijp = ii_hpk + jj_hpk*(jj_hpk-1)/2  ! Index for packed form
    4915            0 :        my_offset = (ijp-1)*bsize_old + (jj_hpk-1)*2*bsize_frm
    4916            0 :        if (which_block==2) my_offset=my_offset+bsize_mat(1)    ! Shift the offset to account for S1.
    4917              :        !my_offset=4
    4918              :        !
    4919              :      else
    4920              :        ! The element belongs either to F3 of F3^H.
    4921              :        ! Now find whether it is the upper or the lower block since only F3 is stored on file.
    4922            0 :        uplo=1; if (row_glob>block_sizes(1,1)) uplo=2
    4923              : 
    4924              :        if (uplo==1) then
    4925            0 :          row_shift=0
    4926            0 :          col_shift=block_sizes(2,1)
    4927              :        else
    4928              :          row_shift=block_sizes(1,1)
    4929              :          col_shift=0
    4930              :        end if
    4931            0 :        ii = row_glob - row_shift
    4932            0 :        jj = col_glob - col_shift
    4933              : 
    4934            0 :        if (uplo==2) then ! Exchange the indices since the symmetric element will be read.
    4935            0 :          swap=jj
    4936            0 :          jj  =ii
    4937            0 :          ii  =swap
    4938              :        end if
    4939              : 
    4940            0 :        my_offset = (ii-1)*bsize_old + (jj-1)*block_sizes(1,3)*bsize_old + (jj-1)*2*bsize_frm
    4941            0 :        my_offset = my_offset + SUM(bsize_mat)
    4942              :        !if (uplo==1) my_offset=my_offset + bsize_mat(1)
    4943              :        !my_offset=0
    4944              :        !if (ii==1.and.jj==1) write(std_out,*)" (1,1) offset = ",my_offset
    4945              :        !if (ii==block_sizes(1,3).and.jj==block_sizes(2,3)) write(std_out,*)" (n,n) offset =", my_offset
    4946            0 :        if (my_offset>=bsize_tot-1*bsize_old) then
    4947            0 :          write(std_out,*)"WARNING (my_offset>bsize_tot-bsize_old),",ii,jj,my_offset,bsize_tot
    4948              :        end if
    4949              :      end if
    4950              : 
    4951            0 :      if (xmpio_max_address(my_offset)) offset_err=1   ! Check for wraparounds.
    4952            0 :      my_el = my_el+1
    4953            0 :      block_displ (my_el+1)=my_offset
    4954            0 :      block_length(my_el+1)=1
    4955            0 :      block_type  (my_el+1)=old_type
    4956            0 :      max_displ = MAX(max_displ,my_offset)
    4957            0 :      min_displ = MIN(min_displ,my_offset)
    4958              :      !if (which_block==3) write(std_out,*)" my_el, which, displ: ",my_el,which_block,block_displ(my_el+1)
    4959              :    end do
    4960              :  end do
    4961              : 
    4962              :  !write(std_out,*)" MAX displ = ",max_displ," my_nels = ",my_nels
    4963              :  !write(std_out,*)" MIN displ = ",MINVAL(block_displ(2:my_nels+1))
    4964              : 
    4965              :  !block_displ (1)=max_displ ! Do not change this value.
    4966              :  !if (min_displ>0) block_displ (1)=min_displ ! Do not change this value.
    4967              : 
    4968              :  block_displ (1)=min_displ
    4969            0 :  block_displ (1)=0
    4970            0 :  block_length(1)=0
    4971            0 :  block_type  (1)=MPI_LB
    4972              : 
    4973            0 :  block_length(my_nels+2)=0
    4974              :  !block_displ (my_nels+2)=bsize_tot
    4975            0 :  block_displ (my_nels+2)=max_displ
    4976            0 :  block_type  (my_nels+2)=MPI_UB
    4977              : 
    4978            0 :  call xmpio_type_struct(my_nels+2,block_length,block_displ,block_type,new_type,mpierr)
    4979              :  !call xmpio_type_struct(my_nels,block_length(2:),block_displ(2:),block_type(2:),new_type,mpierr)
    4980              : 
    4981              :  !call MPI_TYPE_CREATE_INDEXED_BLOCK(my_nels, block_length(2:), block_displ(2:), old_type, new_type, mpierr)
    4982              : 
    4983            0 :  call MPI_TYPE_COMMIT(new_type,mpierr)
    4984              : 
    4985            0 :  ABI_FREE(block_length)
    4986            0 :  ABI_FREE(block_displ)
    4987            0 :  ABI_FREE(block_type)
    4988              : 
    4989            0 : end subroutine xmpio_create_coldistr_from_fp3blocks
    4990              : !!***
    4991              : #endif
    4992              : 
    4993              : !!****f* m_xmpi/xmpi_distrib_2d
    4994              : !! NAME
    4995              : !! xmpi_distrib_2d
    4996              : !!
    4997              : !! FUNCTION
    4998              : !!  Try to optimally distribute nprocs in a 2d grid of shape (n1, n2) given a problem of dimension (size1, size2).
    4999              : !!  Use order string to define priorities:
    5000              : !!      "12" or "21" if both dimensions should be optimized (if not possibile the first one gets optimized).
    5001              : !!      "1" or "2" to optimize only one dimension.
    5002              : !!  Return: exit status in ierr.
    5003              : !!
    5004              : !! SOURCE
    5005              : 
    5006           26 : subroutine xmpi_distrib_2d(nprocs, order, size1, size2, n1, n2, ierr)
    5007              : 
    5008              : !Arguments ------------------------------------
    5009              :  integer,intent(in) :: nprocs, size1, size2
    5010              :  character(len=*),intent(in) :: order
    5011              :  integer,intent(out) :: n1, n2, ierr
    5012              : 
    5013              : !Local variables-------------------------------
    5014              :  integer :: ii
    5015              : !----------------------------------------------------------------------
    5016              : 
    5017           26 :  ierr = 1; n1 = -1; n2 = -1
    5018              : 
    5019           40 :  select case (order)
    5020              :  case ("12")
    5021           14 :    call balance_12()
    5022           14 :    if (ierr /= 0) call balance_1()
    5023              :  case ("21")
    5024            3 :    call balance_21()
    5025            3 :    if (ierr /= 0) call balance_2()
    5026              :  case ("1")
    5027            9 :    call balance_1()
    5028              :  case ("2")
    5029            0 :    call balance_2()
    5030              :  case default
    5031              :    ! Wrong order
    5032           26 :    ierr = -1
    5033              :  end select
    5034              : 
    5035              : contains
    5036              : 
    5037           14 : subroutine balance_12()
    5038              :  ! Try to find n1 x n2 = nprocs so that (size1, size2) are multiple of (n1, n2)
    5039           14 :  do ii=nprocs,1,-1
    5040           14 :    if (mod(size1, ii) == 0 .and. mod(nprocs, ii) == 0 .and. mod(size2, nprocs / ii) == 0) then
    5041           14 :      n1 = ii; n2 = nprocs / ii; ierr = 0; exit
    5042              :    end if
    5043              :  end do
    5044              : 
    5045           14 : end subroutine balance_12
    5046              : 
    5047            3 : subroutine balance_21()
    5048              :  ! Try to find n1 x n2 = nprocs so that (size1, size2) are multiple of (n1, n2)
    5049            3 :  do ii=nprocs,1,-1
    5050            3 :    if (mod(size2, ii) == 0 .and. mod(nprocs, ii) == 0 .and. mod(size1, nprocs / ii) == 0) then
    5051            3 :      n2 = ii; n1 = nprocs / ii; ierr = 0; exit
    5052              :    end if
    5053              :  end do
    5054            3 : end subroutine balance_21
    5055              : 
    5056            9 : subroutine balance_1()
    5057              :  integer :: imod1
    5058              :  ! Try to find n1 x n2 = nprocs so that only size1 is multiple of n1. Allow for some load imbalance.
    5059            9 :  do ii=nprocs,1,-1
    5060            9 :    imod1 = mod(size1, ii)
    5061            9 :    if ((imod1 == 0 .or. imod1 >= nprocs / 2) .and. mod(nprocs, ii) == 0 .and. size2 >= (nprocs/ii)) then
    5062            9 :      n1 = ii; n2 = nprocs / ii; ierr = 0; exit
    5063              :    end if
    5064              :  end do
    5065              : 
    5066            9 :  if (ierr /= 0 .and. nprocs <= size1) then
    5067            0 :    n1 = nprocs; n2 = 1; ierr = 0; return
    5068              :  end if
    5069              : end subroutine balance_1
    5070              : 
    5071            0 : subroutine balance_2()
    5072              :  integer :: imod2
    5073              :  ! Try to find n1 x n2 = nprocs so that only size2 is multiple of n2. Allow for some load imbalance.
    5074            0 :  do ii=nprocs,1,-1
    5075            0 :    imod2 = mod(size2, ii)
    5076            0 :    if ((imod2 == 0 .or. imod2 >= nprocs / 2) .and. mod(nprocs, ii) == 0 .and. size1 >= (nprocs/ii)) then
    5077            0 :      n2 = ii; n1 = nprocs / ii; ierr = 0; exit
    5078              :    end if
    5079              :  end do
    5080              : 
    5081            0 :  if (ierr /= 0 .and. nprocs <= size2) then
    5082            0 :    n2 = nprocs; n1 = 1; ierr = 0; return
    5083              :  end if
    5084              : end subroutine balance_2
    5085              : 
    5086              : end subroutine xmpi_distrib_2d
    5087              : !!***
    5088              : 
    5089              : !----------------------------------------------------------------------
    5090              : 
    5091              : !!****f* m_xmpi/xmpi_split_nsppol
    5092              : !! NAME
    5093              : !! xmpi_split_nsppol
    5094              : !!
    5095              : !! FUNCTION
    5096              : !!  Distribute collinear spins. Also create and return indirect mapping to spin index and init %brange_spin
    5097              : !!
    5098              : !! INPUTS
    5099              : !!  in_comm=Input communicator
    5100              : !!  nsppol=Number of spins
    5101              : !!
    5102              : !! OUTPUT
    5103              : !!  my_nspins=Number of spins treated by this MPI proc
    5104              : !!  my_spins(my_nspins)=Spin index
    5105              : !!  comm_my_is(my_nspins)=Spin communicator for each spin treated by this MPI proc.
    5106              : !!
    5107              : !! SOURCE
    5108              : 
    5109            2 : subroutine xmpi_split_nsppol(in_comm, nsppol, my_nspins, my_spins, comm_my_is)
    5110              : 
    5111              : !Arguments ------------------------------------
    5112              : !scalars
    5113              :  integer,intent(in) :: in_comm, nsppol
    5114              :  integer,intent(out) :: my_nspins
    5115              :  integer,allocatable,intent(out) :: my_spins(:)
    5116              :  type(xcomm_t),allocatable,intent(out) :: comm_my_is(:)
    5117              : 
    5118              : !Local variables-------------------------------
    5119              : !scalars
    5120              :  integer :: spin, my_rank, ierr, color, all_nprocs
    5121              : !arrays
    5122            4 :  integer :: buff_spin(nsppol), comm_spin(nsppol)
    5123              : !----------------------------------------------------------------------
    5124              : 
    5125            2 :  all_nprocs = xmpi_comm_size(in_comm); my_rank = xmpi_comm_rank(in_comm)
    5126              : 
    5127            2 :  my_nspins = 0
    5128            4 :  do spin=1,nsppol
    5129              :    ! NB: If MPI_UNDEFINED is passed as the colour value, the subgroup in which the calling MPI process will be placed is MPI_COMM_NULL
    5130            2 :    color = 1
    5131            2 :    if (nsppol == 2 .and. all_nprocs > 1) then
    5132            0 :      color = xmpi_undefined
    5133            0 :      if (spin == 1 .and. my_rank <= (all_nprocs - 1) / 2) color = 1
    5134            0 :      if (spin == 2 .and. my_rank > (all_nprocs - 1) / 2) color = 1
    5135              :    end if
    5136              : 
    5137            2 :    call xmpi_comm_split(in_comm, color, my_rank, comm_spin(spin), ierr)
    5138            4 :    if (comm_spin(spin) /= xmpi_comm_null) then
    5139            2 :      my_nspins = my_nspins + 1
    5140            2 :      buff_spin(my_nspins) = spin
    5141              :    end if
    5142              :  end do
    5143              : 
    5144            8 :  ABI_MALLOC(comm_my_is, (my_nspins))
    5145            2 :  my_nspins = 0
    5146            4 :  do spin=1,nsppol
    5147            4 :    if (comm_spin(spin) /= xmpi_comm_null) then
    5148            2 :      my_nspins = my_nspins + 1
    5149            2 :      comm_my_is(my_nspins) = xcomm_from_mpi_int(comm_spin(spin), free=.True.)
    5150              :    end if
    5151              :  end do
    5152              : 
    5153            6 :  ABI_MALLOC(my_spins, (my_nspins))
    5154            6 :  my_spins = buff_spin(1:my_nspins)
    5155              : 
    5156            2 : end subroutine xmpi_split_nsppol
    5157              : !!***
    5158              : 
    5159              : ! Init xcomm_t instance from MPI integer. Relase comm_int if optional argument free is set to .True.
    5160              : ! [root]: Rank of the proc treating iteration `iter`
    5161              : 
    5162           86 : type(xcomm_t) function xcomm_from_mpi_int(comm_int, free) result(new)
    5163              : 
    5164              : !Arguments ------------------------------------
    5165              :   integer,intent(inout) :: comm_int
    5166              :   integer :: new_comm, ierr
    5167              :   logical,optional,intent(in) :: free
    5168              : !----------------------------------------------------------------------
    5169              : 
    5170           86 :   new%value = comm_int; new%me = 0; new%nproc = 1
    5171              : #ifdef HAVE_MPI
    5172           86 :   call MPI_Comm_dup(comm_int, new_comm, ierr)
    5173           86 :   new%value = new_comm
    5174           86 :   new%nproc = xmpi_comm_size(new_comm)
    5175           86 :   new%me = xmpi_comm_rank(new_comm)
    5176           86 :   if (present(free)) then
    5177           32 :     if (free) call xmpi_comm_free(comm_int)
    5178              :   end if
    5179              : #endif
    5180           86 : end function xcomm_from_mpi_int
    5181              : 
    5182              : ! Skip iteration `iter` according to rank in xcomm.
    5183              : ! [root]: Rank of the proc treating iteration `iter`
    5184              : 
    5185       428171 : logical function xcomm_skip(xcomm, iter, root)
    5186              : 
    5187              : !Arguments ------------------------------------
    5188              :  class(xcomm_t),intent(in) :: xcomm
    5189              :  integer,intent(in) :: iter
    5190              :  integer,optional,intent(out) :: root
    5191              : 
    5192              : !Local variables-------------------------------
    5193              :  integer :: root__
    5194              : !----------------------------------------------------------------------
    5195              : 
    5196       428171 :  root__ = mod(iter, xcomm%nproc)
    5197       428171 :  xcomm_skip = root__ /= xcomm%me
    5198       428171 :  if (present(root)) root = root__
    5199       428171 : end function xcomm_skip
    5200              : 
    5201          177 : subroutine xcomm_set_to_self(xcomm)
    5202              :  class(xcomm_t),intent(inout) :: xcomm
    5203          177 :  call xcomm%free()
    5204          177 :  xcomm%value = xmpi_comm_self; xcomm%me = 0; xcomm%nproc = 1
    5205          177 : end subroutine xcomm_set_to_self
    5206              : 
    5207           58 : subroutine xcomm_set_to_null(xcomm)
    5208              : !Arguments ------------------------------------
    5209              :  class(xcomm_t),intent(inout) :: xcomm
    5210              : !----------------------------------------------------------------------
    5211              : 
    5212           58 :  call xcomm%free()
    5213           58 :  xcomm%value = xmpi_comm_null
    5214           58 : end subroutine xcomm_set_to_null
    5215              : 
    5216         1830 : subroutine xcomm_free(xcomm)
    5217              : !Arguments ------------------------------------
    5218              :  class(xcomm_t),intent(inout) :: xcomm
    5219              : !----------------------------------------------------------------------
    5220         1830 :  call xmpi_comm_free(xcomm%value)
    5221         1830 :  xcomm%me = -1; xcomm%nproc = 0
    5222         1830 : end subroutine xcomm_free
    5223              : 
    5224              : ! Build sub-communicators in a Cartesian grid.
    5225          652 : subroutine xcomm_from_cart_sub(xcomm, comm_cart, keepdim)
    5226              : 
    5227              : !Arguments ------------------------------------
    5228              :  class(xcomm_t),intent(out) :: xcomm
    5229              :  integer,intent(in) :: comm_cart
    5230              :  logical,intent(in) :: keepdim(:)
    5231              : 
    5232              : !Local variables-------------------------------
    5233              :  integer :: ierr
    5234              : !----------------------------------------------------------------------
    5235              : 
    5236              : #ifdef HAVE_MPI
    5237          652 :  call MPI_CART_SUB(comm_cart, keepdim, xcomm%value, ierr)
    5238              : #endif
    5239          652 :  xcomm%me = xmpi_comm_rank(xcomm%value)
    5240          652 :  xcomm%nproc = xmpi_comm_size(xcomm%value)
    5241              : 
    5242          652 : end subroutine xcomm_from_cart_sub
    5243              : 
    5244              : ! Creates new communicators based on split types and keys
    5245              : 
    5246           30 : type(xcomm_t) function xcomm_split_type(xcomm, split_type, key) result(out_xcomm)
    5247              : 
    5248              :  class(xcomm_t),intent(in) :: xcomm
    5249              :  integer,intent(in),optional :: split_type, key
    5250              : 
    5251              : !Local variables-------------------------------
    5252              :  integer :: split_type__, key__, shared_comm, ierr
    5253              : !----------------------------------------------------------------------
    5254              : 
    5255            0 :  key__ = 0; if (present(key)) key__ = key
    5256              : 
    5257              : #ifdef HAVE_MPI
    5258              :  ! Get node-level communicator
    5259           30 :  split_type__ = MPI_COMM_TYPE_SHARED; if (present(split_type)) split_type__ = split_type
    5260           30 :  call MPI_Comm_split_type(xcomm%value, split_type__, key__, MPI_INFO_NULL, shared_comm, ierr)
    5261           30 :  if (ierr /= MPI_SUCCESS) call xmpi_abort(msg="MPI_COMM_SPLIT_TYPE returned ierr /= 0")
    5262           30 :  out_xcomm = xcomm_from_mpi_int(shared_comm, free=.True.)
    5263              : #else
    5264              :  call out_xcomm%set_to_self()
    5265              : #endif
    5266              : 
    5267           30 : end function xcomm_split_type
    5268              : 
    5269              : ! Prepare a typical gatherv operation in which each MPI rank sends
    5270              : ! `nitems_per_rank(rank+1)` items and each item has length `nelem_per_item`.
    5271              : ! Final results are packed according to the rank of the processor.
    5272              : 
    5273            0 : subroutine xcomm_prep_gatherv(xcomm, nelem_per_item, nitems_per_rank, sendcount, recvcounts, displs)
    5274              : 
    5275              : !Arguments ------------------------------------
    5276              :  class(xcomm_t),intent(in) :: xcomm
    5277              :  integer,intent(in) :: nelem_per_item, nitems_per_rank(xcomm%nproc)
    5278              :  integer,intent(out) :: sendcount
    5279              :  integer, allocatable, intent(out) :: recvcounts(:), displs(:)
    5280              : 
    5281              : !Local variables-------------------
    5282              :  integer :: ii
    5283              : !----------------------------------------------------------------------
    5284              : 
    5285            0 :  ABI_MALLOC(recvcounts, (xcomm%nproc))
    5286            0 :  ABI_MALLOC(displs, (xcomm%nproc))
    5287            0 :  sendcount = nelem_per_item * nitems_per_rank(xcomm%me + 1)
    5288              : 
    5289            0 :  recvcounts(:) = nelem_per_item * nitems_per_rank
    5290            0 :  displs(1) = 0
    5291            0 :  do ii=2,xcomm%nproc
    5292            0 :    displs(ii) = nelem_per_item * sum(nitems_per_rank(1:ii-1))
    5293              :  end do
    5294            0 : end subroutine xcomm_prep_gatherv
    5295              : !!***
    5296              : 
    5297              : ! Debugging tool to print the hostname of the procs in the communicator
    5298            0 : subroutine xcomm_print_names(xcomm)
    5299              : 
    5300              : !Arguments ------------------------------------
    5301              :  class(xcomm_t),intent(in) :: xcomm
    5302              : 
    5303              : !Local variables-------------------
    5304              :  integer :: ip, ierr !, shared_comm
    5305            0 :  character(len=xmpi_max_processor_name) :: my_name, names(xcomm%nproc)
    5306              : !----------------------------------------------------------------------
    5307              : 
    5308            0 :  call xmpi_name(my_name, ierr)
    5309              :  ! FIXME
    5310              :  !call xmpi_allgather(my_name, names, xcomm%value, ierr)
    5311              : 
    5312              : ! ! Get node-level communicator.
    5313              : ! shared_comm = xmpi_comm_world
    5314              : !#ifdef HAVE_MPI
    5315              : ! call MPI_Comm_split_type(xcomm%value, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL, shared_comm, ierr)
    5316              : !#endif
    5317              : ! shared_rank = xmpi_comm_rank(shared_comm)
    5318              : ! shared_size = xmpi_comm_size(shared_comm)
    5319              : 
    5320            0 :  if (xcomm%me == 0) then
    5321            0 :    write(std_out, "(a5,2x,a20)")"rank", "hostname"
    5322            0 :    do ip=0,xcomm%nproc-1
    5323            0 :      write(std_out, "(i5,2x,a20)")ip, trim(names(ip+1))
    5324              :      !write(*,*) 'Global rank', xcomm%me, trim(name(:name_len)), 'Shared rank', shared_rank, 'Shared size', shared_size
    5325              :    end do
    5326              :  end if
    5327              : 
    5328            0 : end subroutine xcomm_print_names
    5329              : !!***
    5330              : 
    5331              : ! True if all procs in xcomm can create a shared memory region. Cache the result.
    5332          174 : logical function xcomm_can_use_shmem(xcomm) result(ok)
    5333              : 
    5334              : !Arguments ------------------------------------
    5335              :  class(xcomm_t),intent(inout) :: xcomm
    5336              : !Local variables-------------------
    5337              :  integer :: ierr, new_comm
    5338              : !----------------------------------------------------------------------
    5339              : 
    5340          174 :  ok = .False.
    5341              : #ifdef HAVE_MPI
    5342          174 :  if (xcomm%can_use_shmem__ == - 1) then
    5343              :    ! First call --> cache result
    5344          126 :    call MPI_COMM_SPLIT_TYPE(xcomm%value, MPI_COMM_TYPE_SHARED, xcomm%me, MPI_INFO_NULL, new_comm, ierr)
    5345          126 :    xcomm%can_use_shmem__ = merge(1, 0, xmpi_comm_size(new_comm) == xcomm%nproc)
    5346          126 :    call xmpi_comm_free(new_comm)
    5347              :  end if
    5348          174 :  ok = xcomm%can_use_shmem__ == 1
    5349              : #endif
    5350              : 
    5351          174 : end function xcomm_can_use_shmem
    5352              : !!***
    5353              : 
    5354           78 : subroutine xcomm_allocate_shared_master(xcomm, count, kind, info, baseptr, win)
    5355              : 
    5356              : !Arguments ------------------------------------
    5357              :  class(xcomm_t),intent(inout) :: xcomm
    5358              :  integer(kind=XMPI_ADDRESS_KIND), intent(in) :: count
    5359              :  integer,intent(in) :: kind, info
    5360              :  type(c_ptr),intent(out) :: baseptr
    5361              :  integer,intent(out) :: win
    5362              : 
    5363              : !Local variables-------------------
    5364              :  integer :: disp_unit, ierr
    5365              :  integer(kind=XMPI_ADDRESS_KIND) :: my_size
    5366              : !----------------------------------------------------------------------
    5367              : 
    5368           78 :  if (.not. xcomm%can_use_shmem()) call xmpi_abort(msg="MPI communicator does not support shared memory allocation!")
    5369              : 
    5370           78 :  select case (kind)
    5371              :  case (sp)
    5372            0 :    disp_unit = xmpi_bsize_sp
    5373              :  case (dp)
    5374           78 :    disp_unit = xmpi_bsize_dp
    5375              :  case default
    5376           78 :    call xmpi_abort(msg="Invalid kind")
    5377              :  end select
    5378              : 
    5379              : #ifdef HAVE_MPI_ALLOCATE_SHARED_CPTR
    5380              :  ! This call is problematic as the API with type(c_ptr) requires mpi_f08 else gcc complains with
    5381              :  ! Error: Type mismatch in argument 'baseptr' at (1); passed TYPE(c_ptr) to INTEGER(8)
    5382              :  ! See https://github.com/pmodels/mpich/issues/2659
    5383              :  ! Converting C_PTR to INTEGER(KIND=MPI_ADDRESS_KIND) with the trick below is not portable:
    5384              :  !address = transfer(baseptr, address)
    5385              : 
    5386           78 :  my_size = 0; if (xcomm%me == 0) my_size = count * disp_unit
    5387           78 :  call MPI_WIN_ALLOCATE_SHARED(my_size, disp_unit, info, xcomm%value, baseptr, win, ierr)
    5388           78 :  if (ierr /= MPI_SUCCESS) call xmpi_abort(msg="mpi_win_allocated_shared returned ierr /= 0")
    5389           78 :  xmpi_count_wins = xmpi_count_wins + 1
    5390              : 
    5391              :  ! Synchronize to ensure memory is allocated.
    5392           78 :  call MPI_Barrier(xcomm%value, ierr)
    5393              : 
    5394           78 :  if (xcomm%me /= 0) then
    5395           52 :    call MPI_WIN_SHARED_QUERY(win, 0, my_size, disp_unit, baseptr, ierr)
    5396           52 :    if (ierr /= MPI_SUCCESS) call xmpi_abort(msg="mpi_win_shared_query returned ierr /= 0")
    5397              :  end if
    5398              : 
    5399              :  ! No local operations prior to this epoch, so give an assertion
    5400           78 :  call MPI_Win_fence(MPI_MODE_NOPRECEDE, win, ierr)
    5401           78 :  if (ierr /= MPI_SUCCESS) call xmpi_abort(msg="mpi_win_shared_query returned ierr /= 0")
    5402              : 
    5403           78 :  call MPI_Barrier(xcomm%value, ierr)
    5404              : 
    5405              : #else
    5406              :  call xmpi_abort(msg="MPI_WIN_ALLOCATE_SHARED with C_PTR is not supported by your MPI library!")
    5407              : #endif
    5408              : 
    5409           78 : end subroutine xcomm_allocate_shared_master
    5410              : !!***
    5411              : 
    5412              : !!****f* m_xmpi/pool2d_from_dims
    5413              : !! NAME
    5414              : !!  pool2d_from_dims
    5415              : !!
    5416              : !! FUNCTION
    5417              : !!  Build pool of MPI procs to distribute (n1 x n2) tasks.
    5418              : !!
    5419              : !! INPUTS
    5420              : !!  n1, n2: dimensions of the problem
    5421              : !!  input_comm: Initial MPI communicator
    5422              : !!  with_pools: Set it to False to use just one pool.
    5423              : !!  [rectangular]: If True, change the number of procs in each pool so that it's possible to
    5424              : !!      create a rectangular grid. Useful for Scalapack algorithms in which 1d grid are not efficient.
    5425              : !!      Default: False.
    5426              : !!
    5427              : !! SOURCE
    5428              : 
    5429            0 : subroutine pool2d_from_dims(pool, n1, n2, input_comm, with_pools, rectangular)
    5430              : 
    5431              : !Arguments-------------------------
    5432              :  class(xmpi_pool2d_t),intent(out) :: pool
    5433              :  integer,intent(in) :: n1, n2, input_comm
    5434              :  logical,intent(in) :: with_pools
    5435              :  logical,optional,intent(in) :: rectangular
    5436              : 
    5437              : !Local variables-------------------
    5438              :  integer :: itask, ntasks, my_rank, nprocs, color, mpierr, jj, i1, i2, my_ntasks, new_comm
    5439              :  integer :: grid_dims(2) !, check(n1, n2)
    5440            0 :  integer,allocatable :: my_inds(:)
    5441              : !----------------------------------------------------------------------
    5442              : 
    5443            0 :  my_rank = xmpi_comm_rank(input_comm); nprocs = xmpi_comm_size(input_comm)
    5444              : 
    5445            0 :  pool%n1 = n1; pool%n2 = n2
    5446            0 :  ABI_MALLOC(pool%treats, (n1, n2))
    5447            0 :  pool%treats = .False.
    5448              : 
    5449            0 :  ntasks = n1 * n2; color = ntasks + 1
    5450              : 
    5451            0 :  if (.not. with_pools) then
    5452            0 :    pool%treats = .True.; color = 1
    5453              : 
    5454            0 :  else if (nprocs <= ntasks) then
    5455            0 :     color = my_rank
    5456            0 :     call xmpi_split_block(ntasks, input_comm, my_ntasks, my_inds)
    5457            0 :     do jj=1,size(my_inds)
    5458            0 :       itask = my_inds(jj) ! = i1 + (i2 - 1) * n1
    5459            0 :       i1 = mod(itask - 1, n1) + 1
    5460            0 :       i2 = 1 + (itask - i1) / n1
    5461            0 :       pool%treats(i1, i2) = .True.
    5462              :     end do
    5463            0 :     ABI_FREE(my_inds)
    5464              :  else
    5465            0 :    i2_loop: do i2=1,n2
    5466            0 :      do i1=1,n1
    5467            0 :        itask = i1 + (i2 - 1) * n1
    5468            0 :        if (xmpi_distrib_with_replicas(itask, ntasks, my_rank, nprocs)) then
    5469            0 :          pool%treats(i1, i2) = .True.; color = itask; exit i2_loop
    5470              :        end if
    5471              :      end do
    5472              :    end do i2_loop
    5473              :  end if
    5474              : 
    5475              :  !DEBUG
    5476              :  ! where (pool%treats)
    5477              :  !   check = 1
    5478              :  ! else where
    5479              :  !   check = 0
    5480              :  ! end where
    5481              :  ! call xmpi_sum(check, input_comm, mpierr)
    5482              :  ! if (any(check == 0)) then
    5483              :  !   write(std_out, *) check
    5484              :  !   call xmpi_abort(msg="Wrong distribution in pool2d_from_dims")
    5485              :  ! end if
    5486              :  !END_DEBUG
    5487              : 
    5488            0 :  call xmpi_comm_split(input_comm, color, my_rank, new_comm, mpierr)
    5489            0 :  pool%comm = xcomm_from_mpi_int(new_comm)
    5490            0 :  call xmpi_comm_free(new_comm)
    5491              : 
    5492            0 :  if (present(rectangular)) then
    5493            0 :    if (rectangular) then
    5494            0 :      if (pool%comm%nproc == 1 .or. is_rectangular_grid(pool%comm%nproc, grid_dims)) return
    5495              : 
    5496            0 :      do jj=pool%comm%nproc-1,1,-1
    5497            0 :        if (is_rectangular_grid(jj, grid_dims)) then
    5498            0 :          color = merge(1, 0, pool%comm%me < jj)
    5499            0 :          call xmpi_comm_split(pool%comm%value, color, pool%comm%me, new_comm, mpierr)
    5500            0 :          call pool%comm%free()
    5501            0 :          pool%comm = xcomm_from_mpi_int(new_comm)
    5502            0 :          call xmpi_comm_free(new_comm)
    5503            0 :          if (color == 0) pool%treats = .False.
    5504              :          exit
    5505              :        end if
    5506              :      end do
    5507              :    end if
    5508              :  end if
    5509              : 
    5510              : contains
    5511              : 
    5512            0 : logical function is_rectangular_grid(nproc, grid_dims) result (ans)
    5513              :  integer,intent(in) :: nproc
    5514              :  integer,intent(out) :: grid_dims(2)
    5515              : 
    5516              : !----------------------------------------------------------------------
    5517              :  integer :: i
    5518              :  ! Search for a rectangular grid of processors
    5519            0 :  i = INT(SQRT(float(nproc)))
    5520            0 :  do while (MOD(nproc,i) /= 0)
    5521            0 :    i = i - 1
    5522              :  end do
    5523            0 :  i = max(i, 1)
    5524              : 
    5525            0 :  grid_dims(1) = i
    5526            0 :  grid_dims(2) = int(nproc / i)
    5527            0 :  ans = grid_dims(1) > 1 .and. grid_dims(2) > 1
    5528              : 
    5529            0 : end function is_rectangular_grid
    5530              : 
    5531              : end subroutine pool2d_from_dims
    5532              : !!***
    5533              : 
    5534              : !!****f* m_xmpi/pool2d_free
    5535              : !! NAME
    5536              : !!  pool2d_free
    5537              : !!
    5538              : !! FUNCTION
    5539              : !!  Free memory
    5540              : 
    5541            0 : subroutine pool2d_free(pool)
    5542              : 
    5543              : !Arguments-------------------------
    5544              :  class(xmpi_pool2d_t),intent(inout) :: pool
    5545              : !----------------------------------------------------------------------
    5546              : 
    5547            0 :  ABI_SFREE(pool%treats)
    5548            0 :  call pool%comm%free()
    5549              : 
    5550            0 : end subroutine pool2d_free
    5551              : !!***
    5552              : 
    5553         6564 : subroutine xmpi_win_fence(assert, win, ierr)
    5554              : 
    5555              : !Arguments ------------------------------------
    5556              :  integer,intent(in) :: win, assert
    5557              :  integer,intent(out) :: ierr
    5558              : !----------------------------------------------------------------------
    5559              : 
    5560         6564 :  ierr = 0
    5561              : #ifdef HAVE_MPI
    5562         6564 :  call MPI_WIN_FENCE(assert, win, ierr)
    5563              : #endif
    5564              : 
    5565         6564 : end subroutine xmpi_win_fence
    5566              : 
    5567           78 : subroutine xmpi_win_free(win, ierr)
    5568              : 
    5569              : !Arguments-------------------------
    5570              :  integer,intent(inout) :: win
    5571              :  integer,intent(out) :: ierr
    5572              : !----------------------------------------------------------------------
    5573              : 
    5574           78 :  ierr = 0
    5575              : #ifdef HAVE_MPI
    5576           78 :  call MPI_WIN_FREE(win, ierr)
    5577           78 :  win = xmpi_undefined
    5578           78 :  xmpi_count_wins = xmpi_count_wins - 1
    5579              : #endif
    5580              : 
    5581           78 : end subroutine xmpi_win_free
    5582              : !!***
    5583              : 
    5584              : ! Return the number of nodes `num_nodes` in the `in_comm` communicator.
    5585              : 
    5586              : subroutine xmpi_get_nodes_in_comm(in_comm, num_nodes, nprocs_per_node)
    5587              : 
    5588              : !Arguments-------------------------
    5589              :  integer,intent(in) :: in_comm
    5590              :  integer,intent(out) :: num_nodes
    5591              :  integer,optional,allocatable,intent(out) :: nprocs_per_node(:)
    5592              : 
    5593              : !Local variables-------------------
    5594              :  integer :: ierr, in_rank, node_comm, node_rank !, masters_comm, color, np
    5595              : !----------------------------------------------------------------------
    5596              : 
    5597              : #ifndef HAVE_MPI
    5598              :  num_nodes = 1
    5599              :  if (present(nprocs_per_node)) then
    5600              :    ABI_MALLOC(nprocs_per_node, (num_nodes))
    5601              :    nprocs_per_node = 1
    5602              :  end if
    5603              : 
    5604              : #else
    5605              :  in_rank = xmpi_comm_rank(in_comm)
    5606              :  call MPI_COMM_SPLIT_TYPE(in_comm, MPI_COMM_TYPE_SHARED, in_rank, MPI_INFO_NULL, node_comm, ierr)
    5607              :  node_rank = xmpi_comm_rank(node_comm)
    5608              :  num_nodes = merge(1, 0, node_rank == 0)
    5609              :  call xmpi_sum(num_nodes, in_comm, ierr)
    5610              : 
    5611              :  !if (present(nprocs_per_node)) then
    5612              :  !  ABI_MALLOC(nprocs_per_node, (num_nodes))
    5613              :  !  color = merge(0, 1, node_rank == 0)
    5614              :  !  call xmpi_comm_split(in_comm, color, in_rank, masters_comm, ierr)
    5615              :  !  if (color == 0) then
    5616              :  !    np = xmpi_comm_size(node_comm)
    5617              :  !    call MPI_GATHER(np, 1, MPI_INT, nprocs_per_node, 1, MPI_INT, 0, masters_comm, ierr)
    5618              :  !  end if
    5619              :  !  call xmpi_comm_free(masters_comm)
    5620              :  !end if
    5621              : 
    5622              :  call xmpi_comm_free(node_comm)
    5623              : #endif
    5624              : 
    5625              : end subroutine xmpi_get_nodes_in_comm
    5626              : !!***
    5627              : 
    5628            0 : end module m_xmpi
    5629              : !!***
        

Generated by: LCOV version 2.3-1