!= MPI 関連ルーチン
!
!= MPI related routines
!
! Authors::   Yoshiyuki O. Takahashi, Shin-ichi Takehiro
! Version::   $Id: mpi_wrapper.F90,v 1.7 2025/12/14 18:00:00 takepiro Exp $
! Tag Name::  $Name:  $
! Copyright:: Copyright (C) GFD Dennou Club, 2008--2025. All rights reserved.
! License::   See COPYRIGHT[link:../../../COPYRIGHT]
!

module mpi_wrapper
  !
  != MPI 関連ルーチン
  !
  != MPI related routines
  !
  ! <b>Note that Japanese and English are described in parallel.</b>
  !
  ! MPI 関係の変数の管理と MPI 関係ラッパールーチンのモジュール. 
  !
  ! This is a module containing MPI-related variables and wrapper routines. 
  !
  !== Procedures List
  !
!!$  ! RadiationFluxDennouAGCM :: 放射フラックスの計算
!!$  ! RadiationFinalize       :: 終了処理 (モジュール内部の変数の割り付け解除)
!!$  ! ------------            :: ------------
!!$  ! RadiationFluxDennouAGCM :: Calculate radiation flux
!!$  ! RadiationFinalize       :: Termination (deallocate variables in this module)
  !
  !== NAMELIST
  !
!!$  ! NAMELIST#radiation_DennouAGCM_nml
  !

  ! モジュール引用 ; USE statements
  !

  ! 種別型パラメタ
  ! Kind type parameter
  !
  use dc_types, only: DP, &      ! 倍精度実数型. Double precision.
    &                 STRING, &  ! 文字列.       Strings.
    &                 TOKEN      ! キーワード.   Keywords.

#ifdef LIB_MPI
  ! MPI
  !
  use mpi
#endif

  ! 宣言文 ; Declaration statements
  !
  implicit none
  private

  ! 公開手続き
  ! Public procedure
  !
  public :: MPIWrapperInit
  public :: MPIWrapperFinalize
  public :: MPIWrapperISend
  public :: MPIWrapperIRecv
  public :: MPIWrapperWait
  public :: MPIWrapperFindMaxVal
  public :: MPIWrapperChkTrue

  ! 公開変数
  ! Public variables
  !
  integer, save, public :: nprocs
                           ! Number of MPI processes
  integer, save, public :: myrank
                           ! My rank
  logical, save, public :: flag_mpi_init
                           ! info for MPI_INIT 
  
  ! 非公開変数
  ! Private variables
  !


  interface MPIWrapperISend
    module procedure &
      MPIWrapperISend_logical_1d, &
      MPIWrapperISend_int_1d    , &
      MPIWrapperISend_dble_1d   , &
      MPIWrapperISend_dble_2d   , &
      MPIWrapperISend_dble_3d   , &
      MPIWrapperISend_dble_4d
  end interface

  interface MPIWrapperIRecv
    module procedure &
      MPIWrapperIRecv_logical_1d, &
      MPIWrapperIRecv_int_1d    , &
      MPIWrapperIRecv_dble_1d   , &
      MPIWrapperIRecv_dble_2d   , &
      MPIWrapperIRecv_dble_3d   , &
      MPIWrapperIRecv_dble_4d
  end interface

  interface MPIWrapperFindMaxVal
    module procedure &
      MPIWrapperFindMaxVal_dble_1d
  end interface

  interface MPIWrapperChkTrue
    module procedure &
      MPIWrapperChkTrue_1d
  end interface

  interface MPIWrapperAbort
    module procedure &
      MPIWrapperStop
  end interface


contains

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperInit
    !
    ! MPI の初期化
    !
    ! Initialization of MPI
    !

    ! モジュール引用 ; USE statements
    !


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr

#endif

    nprocs = 1
    myrank = 0
    flag_mpi_init = .false.
    
#ifdef LIB_MPI

    call mpi_init( ierr )
    call mpi_comm_size( mpi_comm_world, nprocs, ierr )
    call mpi_comm_rank( mpi_comm_world, myrank, ierr )

    flag_mpi_init = .true.
#endif

  end subroutine MPIWrapperInit

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperFinalize
    !
    ! MPI の終了処理
    !
    ! Finalization of MPI
    !

    ! モジュール引用 ; USE statements
    !


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr


    call mpi_finalize( ierr )
    flag_mpi_init = .false.

#endif

  end subroutine MPIWrapperFinalize

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperStop
    !
    ! MPI の異常終了処理
    !
    ! Abort of MPI
    !

    ! モジュール引用 ; USE statements
    !

#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: errorcode = 9
    integer :: ierr


    call mpi_abort( mpi_comm_world, errorcode, ierr )
    call MPIWrapperFinalize
    stop

#endif

  end subroutine MPIWrapperstop

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperWait( ireq )
    !
    ! MPI 通信終了まで待機
    !
    ! Wait finishing MPI transfer
    !

    ! モジュール引用 ; USE statements
    !


    integer, intent(inout) :: ireq
                               ! request number


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: istatus( MPI_STATUS_SIZE )


    call mpi_wait( ireq, istatus, ierr )

#endif

  end subroutine MPIWrapperWait

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperISend_logical_1d(  &
    & idest, im,                      & ! (in)
    & buf,                            & ! (in)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 1D 倍精度配列の非ブロッキング通信(送信)
    !
    ! Non-blocking transfer (send) of real(8) 1D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idest
                              ! Process number of destination
    integer , intent(in ) :: im
                              ! Size of 1st dimension of sent data
    logical , intent(in ) :: buf( im )
                              ! Array to be sent
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag

    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_isend( buf, isize, &
      mpi_logical, idest, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperISend_logical_1d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperIRecv_logical_1d(  &
    & idep, im,                       & ! (in)
    & buf,                            & ! (out)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 1D 倍精度配列の非ブロッキング通信(受信)
    !
    ! Non-blocking transfer (receive) of real(8) 1D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idep
                              ! Process number of departure
    integer , intent(in ) :: im
                              ! Size of 1st dimension of received data
    logical , intent(out) :: buf( im )
                              ! Array to be received
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_irecv( buf, isize, &
      mpi_logical, idep, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperIRecv_logical_1d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperISend_int_1d(  &
    & idest, im,                      & ! (in)
    & buf,                            & ! (in)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 1D 倍精度配列の非ブロッキング通信(送信)
    !
    ! Non-blocking transfer (send) of real(8) 1D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idest
                              ! Process number of destination
    integer , intent(in ) :: im
                              ! Size of 1st dimension of sent data
    integer , intent(in ) :: buf( im )
                              ! Array to be sent
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_isend( buf, isize, &
      mpi_integer, idest, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperISend_int_1d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperIRecv_int_1d(  &
    & idep, im,                       & ! (in)
    & buf,                            & ! (out)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 1D 倍精度配列の非ブロッキング通信(受信)
    !
    ! Non-blocking transfer (receive) of real(8) 1D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idep
                              ! Process number of departure
    integer , intent(in ) :: im
                              ! Size of 1st dimension of received data
    integer , intent(out) :: buf( im )
                              ! Array to be received
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_irecv( buf, isize, &
      mpi_integer, idep, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperIRecv_int_1d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperISend_dble_1d( &
    & idest, im,                      & ! (in)
    & buf,                            & ! (in)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 1D 倍精度配列の非ブロッキング通信(送信)
    !
    ! Non-blocking transfer (send) of real(8) 1D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idest
                              ! Process number of destination
    integer , intent(in ) :: im
                              ! Size of 1st dimension of sent data
    real(DP), intent(in ) :: buf( im )
                              ! Array to be sent
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_isend( buf, isize, &
      mpi_double_precision, idest, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperISend_dble_1d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperIRecv_dble_1d( &
    & idep, im,                       & ! (in)
    & buf,                            & ! (out)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 1D 倍精度配列の非ブロッキング通信(受信)
    !
    ! Non-blocking transfer (receive) of real(8) 1D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idep
                              ! Process number of departure
    integer , intent(in ) :: im
                              ! Size of 1st dimension of received data
    real(DP), intent(out) :: buf( im )
                              ! Array to be received
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_irecv( buf, isize, &
      mpi_double_precision, idep, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperIRecv_dble_1d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperISend_dble_2d( &
    & idest, im, jm,                  & ! (in)
    & buf,                            & ! (in)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 2D 倍精度配列の非ブロッキング通信(送信)
    !
    ! Non-blocking transfer (send) of real(8) 2D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idest
                              ! Process number of destination
    integer , intent(in ) :: im
                              ! Size of 1st dimension of sent data
    integer , intent(in ) :: jm
                              ! Size of 2nd dimension of sent data
    real(DP), intent(in ) :: buf( im, jm )
                              ! Array to be sent
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_isend( buf, isize, &
      mpi_double_precision, idest, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperISend_dble_2d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperIRecv_dble_2d( &
    & idep, im, jm,                   & ! (in)
    & buf,                            & ! (out)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 2D 倍精度配列の非ブロッキング通信(受信)
    !
    ! Non-blocking transfer (receive) of real(8) 2D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idep
                              ! Process number of destination
    integer , intent(in ) :: im
                              ! Size of 1st dimension of received data
    integer , intent(in ) :: jm
                              ! Size of 2nd dimension of received data
    real(DP), intent(out) :: buf( im, jm )
                              ! Array to be received
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_irecv( buf, isize, &
      mpi_double_precision, idep, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperIRecv_dble_2d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperISend_dble_3d( &
    & idest, im, jm, km,              & ! (in)
    & buf,                            & ! (in)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 3D 倍精度配列の非ブロッキング通信(送信)
    !
    ! Non-blocking transfer (send) of real(8) 3D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idest
                              ! Process number of destination
    integer , intent(in ) :: im
                              ! Size of 1st dimension of sent data
    integer , intent(in ) :: jm
                              ! Size of 2nd dimension of sent data
    integer , intent(in ) :: km
                              ! Size of 3rd dimension of sent data
    real(DP), intent(in ) :: buf( im, jm, km )
                              ! Array to be sent
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_isend( buf, isize, &
      mpi_double_precision, idest, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperISend_dble_3d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperIRecv_dble_3d( &
    & idep, im, jm, km,               & ! (in)
    & buf,                            & ! (out)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 3D 倍精度配列の非ブロッキング通信(受信)
    !
    ! Non-blocking transfer (receive) of real(8) 3D array
    !

    ! モジュール引用 ; USE statements
    !

    integer , intent(in ) :: idep
                              ! Process number of departure
    integer , intent(in ) :: im
                              ! Size of 1st dimension of received data
    integer , intent(in ) :: jm
                              ! Size of 2nd dimension of received data
    integer , intent(in ) :: km
                              ! Size of 3rd dimension of received data
    real(DP), intent(out) :: buf( im, jm, km )
                              ! Array to be received
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_irecv( buf, isize, &
      mpi_double_precision, idep, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperIRecv_dble_3d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperISend_dble_4d( &
    & idest, im, jm, km, lm,          & ! (in)
    & buf,                            & ! (in)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 4D 倍精度配列の非ブロッキング通信(送信)
    !
    ! Non-blocking transfer (send) of real(8) 4D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idest
                              ! Process number of destination
    integer , intent(in ) :: im
                              ! Size of 1st dimension of sent data
    integer , intent(in ) :: jm
                              ! Size of 2nd dimension of sent data
    integer , intent(in ) :: km
                              ! Size of 3rd dimension of sent data
    integer , intent(in ) :: lm
                              ! Size of 4th dimension of sent data
    real(DP), intent(in ) :: buf( im, jm, km, lm )
                              ! Array to be sent
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_isend( buf, isize, &
      mpi_double_precision, idest, tag, mpi_comm_world, &
      ireq, ierr )

#endif

  end subroutine MPIWrapperISend_dble_4d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperIRecv_dble_4d( &
    & idep, im, jm, km, lm,           & ! (in)
    & buf,                            & ! (out)
    & ireq,                           & ! (out)
    & itag                            & ! (in) optional
    & )
    !
    ! 4D 倍精度配列の非ブロッキング通信(受信)
    !
    ! Non-blocking transfer (receive) of real(8) 4D array
    !

    ! モジュール引用 ; USE statements
    !


    integer , intent(in ) :: idep
                              ! Process number of departure
    integer , intent(in ) :: im
                              ! Size of 1st dimension of received data
    integer , intent(in ) :: jm
                              ! Size of 2nd dimension of received data
    integer , intent(in ) :: km
                              ! Size of 3rd dimension of received data
    integer , intent(in ) :: lm
                              ! Size of 4th dimension of received data
    real(DP), intent(out) :: buf( im, jm, km, lm )
                              ! Array to be received
    integer , intent(out) :: ireq
                              ! Request number
    integer , intent(in ), optional :: itag
                              ! Size of 1st dimension of sent data


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer :: ierr
    integer :: isize
    integer :: tag


    isize = size( buf )

    if ( present( itag ) ) then
      tag = itag
    else
      tag = 1
    end if

    call mpi_irecv( buf, isize, &
      mpi_double_precision, idep, tag, mpi_comm_world, &
      ireq, ierr )

#endif


  end subroutine MPIWrapperIRecv_dble_4d

  !--------------------------------------------------------------------------------------

  subroutine MPIWrapperFindMaxVal_dble_1d( &
    & lmax, a_LocalMax,                    & ! (in)
    & a_GlobalMax                          & ! (out)
    & )
    !
    ! 全球の最大値を探す
    !
    ! Find the maximum of the globe
    !

    ! モジュール引用 ; USE statements
    !

    integer , intent(in ) :: lmax
    real(DP), intent(in ) :: a_LocalMax (1:lmax)
    real(DP), intent(out) :: a_GlobalMax(1:lmax)


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer  :: idep          ! Process number of departure
    integer  :: idest         ! Process number of destination
    real(DP) :: a_Buf(1:lmax) ! Array to be received
    integer  :: ireq          ! Request number
    real(DP), allocatable :: aa_LocalMax(:,:)
    integer               :: l
    integer               :: n


    if ( myrank == 0 ) then

      allocate( aa_LocalMax( 1:lmax, 0:nprocs-1 ) )

      aa_LocalMax(:,0) = a_LocalMax

      do n = 1, nprocs-1
        idep = n
        call MPIWrapperIRecv(  &
          & idep, lmax,        & ! (in)
          & aa_LocalMax(:,n),  & ! (out)
          & ireq               & ! (out)
          & )
        call MPIWrapperWait( ireq )
      end do

      do l = 1, lmax
        n = 0
        a_GlobalMax(l) = aa_LocalMax(l,n)
        do n = 1, nprocs-1
          if ( a_GlobalMax(l) < aa_LocalMax(l,n) ) then
            a_GlobalMax(l) = aa_LocalMax(l,n)
          end if
        end do
      end do

      a_Buf = a_GlobalMax
      do n = 1, nprocs-1
        idest = n
        call MPIWrapperISend(  &
          & idest, lmax,       & ! (in)
          & a_Buf,             & ! (in)
          & ireq               & ! (out)
          & )
        call MPIWrapperWait( ireq )
      end do

      deallocate( aa_LocalMax )

    else

      idest = 0
      a_Buf = a_LocalMax
      call MPIWrapperISend(  &
        & idest, lmax,       & ! (in)
        & a_Buf,             & ! (in)
        & ireq               & ! (out)
        & )
      call MPIWrapperWait( ireq )

      idep = 0
      call MPIWrapperIRecv(  &
        & idep, lmax,        & ! (in)
        & a_Buf,             & ! (out)
        & ireq               & ! (out)
        & )
      call MPIWrapperWait( ireq )

      a_GlobalMax = a_Buf

    end if

#else

    a_GlobalMax = a_LocalMax

#endif


  end subroutine MPIWrapperFindMaxVal_dble_1d

  !--------------------------------------------------------------------------------------
  !
  ! A value of a_GlobalLogical(k) is true, if a_LocalLogical(k) is true 
  ! at least in a process.
  !
  subroutine MPIWrapperChkTrue_1d( &
    & lmax, a_LocalLogical,        & ! (in)
    & a_GlobalLogical              & ! (out)
    & )
    !
    ! 全球の最大値を探す
    !
    ! Find the maximum of the globe
    !

    ! モジュール引用 ; USE statements
    !

    integer, intent(in ) :: lmax
    logical, intent(in ) :: a_LocalLogical (1:lmax)
    logical, intent(out) :: a_GlobalLogical(1:lmax)


#ifdef LIB_MPI

    ! 作業変数
    ! Work variables
    !
    integer  :: idep          ! Process number of departure
    integer  :: idest         ! Process number of destination
    logical  :: a_Buf(1:lmax) ! Array to be received
    integer  :: ireq          ! Request number
    logical, allocatable :: aa_LocalLogical(:,:)
    integer               :: l
    integer               :: n


    if ( myrank == 0 ) then

      allocate( aa_LocalLogical( 1:lmax, 0:nprocs-1 ) )

      aa_LocalLogical(:,0) = a_LocalLogical

      do n = 1, nprocs-1
        idep = n
        call MPIWrapperIRecv(     &
          & idep, lmax,           & ! (in)
          & aa_LocalLogical(:,n), & ! (out)
          & ireq                  & ! (out)
          & )
        call MPIWrapperWait( ireq )
      end do

      do l = 1, lmax
        n = 0
        a_GlobalLogical(l) = aa_LocalLogical(l,n)
        do n = 1, nprocs-1
          if (  aa_LocalLogical(l,n) ) then
            a_GlobalLogical(l) = aa_LocalLogical(l,n)
          end if
        end do
      end do

      a_Buf = a_GlobalLogical
      do n = 1, nprocs-1
        idest = n
        call MPIWrapperISend(  &
          & idest, lmax,       & ! (in)
          & a_Buf,             & ! (in)
          & ireq               & ! (out)
          & )
        call MPIWrapperWait( ireq )
      end do

      deallocate( aa_LocalLogical )

    else

      idest = 0
      a_Buf = a_LocalLogical
      call MPIWrapperISend(  &
        & idest, lmax,       & ! (in)
        & a_Buf,             & ! (in)
        & ireq               & ! (out)
        & )
      call MPIWrapperWait( ireq )

      idep = 0
      call MPIWrapperIRecv(  &
        & idep, lmax,        & ! (in)
        & a_Buf,             & ! (out)
        & ireq               & ! (out)
        & )
      call MPIWrapperWait( ireq )

      a_GlobalLogical = a_Buf

    end if

#else

    a_GlobalLogical = a_LocalLogical

#endif


  end subroutine MPIWrapperChkTrue_1d

  !--------------------------------------------------------------------------------------

end module mpi_wrapper
