TABLE OF CONTENTS


ABINIT/m_mpiotk [ Modules ]

[ Top ] [ Modules ]

NAME

  m_mpiotk

FUNCTION

  This module provides helper functions for MPI-IO operations.

COPYRIGHT

 Copyright (C) 2009-2018 ABINIT group (MG)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt.

PARENTS

CHILDREN

SOURCE

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 
27 #include "abi_common.h"
28 
29 MODULE m_mpiotk
30 
31  use defs_basis
32  use m_abicore
33  use m_errors
34  use m_xmpi
35 #if defined HAVE_MPI2 && defined HAVE_MPI_IO
36  use mpi
37 #endif
38 
39  implicit none
40 
41 #if defined HAVE_MPI1 && defined HAVE_MPI_IO
42  include 'mpif.h'
43 #endif
44 
45  private
46 
47 !public procedures.
48 #ifdef HAVE_MPI_IO
49  public :: mpiotk_read_fsuba_dp2D     ! (individual|collective) read of a 2D sub-matrix
50  public :: mpiotk_write_fsuba_dp2D    ! (individual|collective) write of a 2D sub-matrix
51 
52  public :: mpiotk_read_fsuba_dpc3D    ! (individual|collective) read of a 3D sub-matrix
53  !public :: mpiotk_write_fsuba_dpc3D  ! (individual|collective) write of a 3D sub-matrix
54 
55  public :: mpiotk_read_fsuba_dpc4D    ! (individual|collective) read of a 4D sub-matrix
56  !public :: mpiotk_write_fsuba_dpc4D  ! (individual|collective) write of a 4D sub-matrix
57 #else
58  public :: no_mpiotk
59 #endif

m_mpiotk/mpiotk_read_fsuba_dp2D [ Functions ]

[ Top ] [ m_mpiotk ] [ Functions ]

NAME

  mpiotk_read_fsuba_dp2D

FUNCTION

  Read a block of contiguous data stored in a 2D matrix.
  Data is placed within Fortran records. Target: complex data stored in a real array.

NOTES

  The value of ierr should always be checked by the caller

INPUTS

  fh = MPI-IO file handler.
  offset =
  sizes(2)
  subsizes(2)
  starts(2)
  bufsz = dimension of buffer (takes into accout both real and imaginary part)
  chunk_bsize =
  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading.
  comm = MPI communicator

 OUTPUTS
  buffer(bufsz)
  ierr=status error. A non zero value indicates that chunk_bsize is smaller that the fortran record
    and therefore bufsz has not been read.

PARENTS

      m_wfk

CHILDREN

      mpi_file_read,mpi_file_read_all,mpi_file_set_view,mpi_type_free
      xmpi_max,xmpi_min,xmpio_create_fsubarray_4d

SOURCE

264 subroutine mpiotk_read_fsuba_dp2D(fh,offset,sizes,subsizes,starts,bufsz,buffer,chunk_bsize,sc_mode,comm,ierr)
265 
266 
267 !This section has been created automatically by the script Abilint (TD).
268 !Do not modify the following lines by hand.
269 #undef ABI_FUNC
270 #define ABI_FUNC 'mpiotk_read_fsuba_dp2D'
271 !End of the abilint section
272 
273  implicit none
274 
275 !Arguments ------------------------------------
276 !scalars
277  integer,intent(in) :: fh,comm,bufsz,sc_mode
278  integer,intent(out) :: ierr
279  integer(XMPI_OFFSET_KIND),intent(in) :: offset,chunk_bsize
280 !arrays
281  integer,intent(in) :: sizes(2),subsizes(2),starts(2)
282  real(dp),intent(out) :: buffer(bufsz)
283 
284 !Local variables ------------------------------
285 !scalars
286  integer :: mpierr,ptr,ncount,myfh
287  integer :: fsub_type,my_ncalls,icall,ncalls,subs_x
288  integer(XMPI_OFFSET_KIND) :: my_offset,my_offpad
289  !character(len=500) :: msg
290 !arrays
291  integer :: call_subsizes(2),call_starts(2)
292  integer,allocatable :: my_basead(:),my_subsizes(:,:),my_starts(:,:)
293  real(dp),allocatable :: dummy_buf(:,:)
294 
295 !************************************************************************
296 
297  ! Workaround for XLF
298  myfh = fh
299 
300  if (bufsz < 2 * PRODUCT(subsizes) ) then
301    MSG_ERROR("bufsz is too small")
302  end if
303 
304  if (sc_mode==xmpio_single) then
305    ! This makes the automatic tests fail (cut3d)
306    !MSG_WARNING("comm != xmpi_comm_self")
307  else if (sc_mode==xmpio_collective) then
308    continue
309  else
310    MSG_ERROR("Wrong sc_mode")
311  end if
312 
313  subs_x  = subsizes(1)
314 
315  call setup_fsuba_dp2D(sizes,subsizes,starts,chunk_bsize,my_basead,my_subsizes,my_starts,my_ncalls,ncalls,comm,ierr)
316  if (ierr/=0) RETURN
317 
318  do icall=1,ncalls
319 
320    if (icall <= my_ncalls) then
321      call_subsizes = my_subsizes(:,icall)
322      call_starts   = my_starts(:,icall)
323      ptr           = my_basead(icall)
324    else
325      ! Fake values needed to call read_all collectively.
326      call_subsizes = (/subs_x, 1/)
327      call_starts   = starts
328    end if
329    ncount = PRODUCT(call_subsizes)
330    !write(std_out,*)"  icall,ptr, ncount, ",icall,ptr,ncount
331    !write(std_out,*)"  call_starts",call_starts
332    !write(std_out,*)"  call_subsizes",call_subsizes
333 
334    ! Create subarry file view.
335    call xmpio_create_fsubarray_2D(sizes,call_subsizes,call_starts,MPI_DOUBLE_COMPLEX,fsub_type,my_offpad,mpierr)
336    ABI_CHECK_MPI(mpierr,"fsubarray_2D")
337 
338    ! Update the offset.
339    my_offset = offset + my_offpad
340 
341    call MPI_FILE_SET_VIEW(myfh, my_offset, MPI_BYTE, fsub_type, 'native', MPI_INFO_NULL, mpierr)
342    ABI_CHECK_MPI(mpierr,"SET_VIEW")
343 
344    call MPI_TYPE_FREE(fsub_type, mpierr)
345    ABI_CHECK_MPI(mpierr,"MPI_TYPE_FREE")
346 
347    if (sc_mode==xmpio_collective) then
348      ! Collective read
349      if (icall <= my_ncalls) then
350        call MPI_FILE_READ_ALL(myfh, buffer(ptr), ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
351      else
352        ABI_MALLOC(dummy_buf,(2,subs_x))
353        call MPI_FILE_READ_ALL(myfh, dummy_buf, ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
354        ABI_FREE(dummy_buf)
355      end if
356      ABI_CHECK_MPI(mpierr,"FILE_READ_ALL")
357 
358    else
359      ! Individual read.
360      call MPI_FILE_READ(myfh, buffer(ptr), ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
361      ABI_CHECK_MPI(mpierr,"FILE_READ")
362    end if
363 
364  end do
365 
366  ABI_FREE(my_subsizes)
367  ABI_FREE(my_starts)
368  ABI_FREE(my_basead)
369 
370 end subroutine mpiotk_read_fsuba_dp2D

m_mpiotk/mpiotk_read_fsuba_dpc3D [ Functions ]

[ Top ] [ m_mpiotk ] [ Functions ]

NAME

  mpiotk_read_fsuba_dpc3D

FUNCTION

  Read of a block of contiguous data stored in a 3D matrix.
  Data is placed within Fortran records. Target: complex data stored in a complex array.

NOTES

  The value of ierr should always be checked by the caller

INPUTS

  fh = MPI-IO file handler.
  offset =
  sizes(3)
  subsizes(3)
  starts(3)
  bufsz = dimension of cbuffer
  chunk_bsize =
  comm = MPI communicator
  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading.

 OUTPUTS
  cbuffer(bufsz)
  ierr=status error. A non zero value indicates that chunk_bsize is smaller that the fortran record
    and therefore bufsz has not been read.

PARENTS

      m_io_screening

CHILDREN

      mpi_file_read,mpi_file_read_all,mpi_file_set_view,mpi_type_free
      xmpi_max,xmpi_min,xmpio_create_fsubarray_4d

SOURCE

560 subroutine mpiotk_read_fsuba_dpc3D(fh,offset,sizes,subsizes,starts,bufsz,cbuffer,chunk_bsize,sc_mode,comm,ierr)
561 
562 
563 !This section has been created automatically by the script Abilint (TD).
564 !Do not modify the following lines by hand.
565 #undef ABI_FUNC
566 #define ABI_FUNC 'mpiotk_read_fsuba_dpc3D'
567 !End of the abilint section
568 
569  implicit none
570 
571 !Arguments ------------------------------------
572 !scalars
573  integer,intent(in) :: fh,comm,bufsz,sc_mode
574  integer,intent(out) :: ierr
575  integer(XMPI_OFFSET_KIND),intent(in) :: offset,chunk_bsize
576 !arrays
577  integer,intent(in) :: sizes(3),subsizes(3),starts(3)
578  complex(dpc),intent(out) :: cbuffer(bufsz)
579 
580 !Local variables-------------------------------
581 !scalars
582  integer :: mpierr,nz2read,nz_chunk,ptr,ncount,myfh
583  integer :: fsub_type,my_ncalls,icall,zrest,ncalls
584  integer :: size_x,size_y,subs_x,subs_y,subs_xy,start_x,start_y,start_z,stop_z
585  integer(XMPI_OFFSET_KIND) :: my_offset,my_offpad
586  !character(len=500) :: msg
587 !arrays
588  integer :: call_subsizes(3),call_starts(3)
589  integer,allocatable :: my_basead(:),my_subsizes(:,:),my_starts(:,:)
590  complex(dpc),allocatable :: dummy_cbuf(:)
591 
592 !************************************************************************
593 
594  ! Workaround for XLF
595  myfh = fh
596 
597  if (bufsz < PRODUCT(subsizes) ) then
598    MSG_ERROR("bufsz is too small")
599  end if
600 
601  if (sc_mode==xmpio_single) then
602   ABI_CHECK(comm==xmpi_comm_self,"comm != xmpi_comm_self")
603  else if (sc_mode==xmpio_collective) then
604    continue
605  else
606    MSG_ERROR("Wrong sc_mode")
607  end if
608 
609  size_x  = sizes(1)
610  size_y  = sizes(2)
611  subs_x  = subsizes(1)
612  subs_y  = subsizes(2)
613  subs_xy = subs_x * subs_y
614  start_x = starts(1)
615  start_y = starts(2)
616  start_z = starts(3)
617  stop_z  = start_z + subsizes(3)-1 ! last column to read
618 
619  ! Read rows in blocks of size nz_chunk:
620  ! MPI-IO crashes if we try to read data > 2Gb in a single call.
621  nz2read = subsizes(3)
622  nz_chunk = nz2read
623  if ( (subs_xy*nz2read*xmpi_bsize_dpc) > chunk_bsize) then
624    nz_chunk = chunk_bsize / (subs_xy*xmpi_bsize_dpc)
625    !if (nz_chunk == 0) nz_chunk = 50
626  end if
627 
628  call xmpi_min(nz_chunk,ierr,comm,mpierr)
629  if (ierr == 0) then
630    ierr = 1
631    RETURN
632  end if
633  ierr = 0
634 
635  ! my_ncalls : number of read needed to fill my cbuffer.
636  ! ncalls    : max number of read in comm (needed for collective operations).
637  my_ncalls = nz2read / nz_chunk
638  zrest = MOD(nz2read, nz_chunk)
639  if (zrest /= 0) my_ncalls = my_ncalls + 1
640 
641  call xmpi_max(my_ncalls,ncalls,comm,mpierr)
642 
643  ! Compute arrays used to define the file view.
644  ABI_MALLOC(my_subsizes,(3,ncalls))
645  ABI_MALLOC(my_starts,(3,ncalls))
646  ABI_MALLOC(my_basead,(ncalls))
647 
648  do icall=1,my_ncalls
649    if (icall*nz_chunk <= nz2read) then
650      my_subsizes(:,icall) = (/subs_x, subs_y, nz_chunk/)
651      my_starts(:,icall) = (/start_x, start_y, (icall-1) * nz_chunk + start_z/)
652      my_basead(icall) = 1 + (icall-1)*subs_xy*nz_chunk
653    else
654      ! Two cases:
655      ! 1) nz2read > nz_chunk and not divisible by nz2read
656      ! 2) nz2read < nz_chunk
657      my_subsizes(:,icall) = (/subs_x, subs_y, zrest/)
658      if (nz2read >= nz_chunk) then
659        my_starts(:,icall) = (/start_x, start_y, stop_z-zrest+1/)
660        my_basead(icall) = 1 + (nz2read-zrest) * subs_xy
661      else
662        my_starts(:,icall) = starts
663        my_basead(icall) = 1
664      end if
665    end if
666  end do
667  !write(std_out,*)" >>>> my_ncalls, ncalls, nz2read, nz_chunk ",my_ncalls,ncalls,nz2read,nz_chunk
668 
669  do icall=1,ncalls
670 
671    if (icall <= my_ncalls) then
672      call_subsizes = my_subsizes(:,icall)
673      call_starts   = my_starts(:,icall)
674      ptr           = my_basead(icall)
675    else
676      ! Fake values needed to call read_all collectively.
677      call_subsizes = (/subs_x, 1, 1/)
678      call_starts   = starts
679    end if
680    ncount = PRODUCT(call_subsizes)
681    !write(std_out,*)"  icall,ptr, ncount, ",icall,ptr,ncount
682    !write(std_out,*)"  call_starts",call_starts
683    !write(std_out,*)"  call_subsizes",call_subsizes
684 
685    ! Create subarry file view.
686    call xmpio_create_fsubarray_3D(sizes,call_subsizes,call_starts,MPI_DOUBLE_COMPLEX,&
687 &    fsub_type,my_offpad,mpierr)
688    ABI_CHECK_MPI(mpierr,"fsubarray_3D")
689 
690    ! Update the offset.
691    my_offset = offset + my_offpad
692 
693    call MPI_FILE_SET_VIEW(myfh, my_offset, MPI_BYTE, fsub_type, 'native', MPI_INFO_NULL, mpierr)
694    ABI_CHECK_MPI(mpierr,"SET_VIEW")
695 
696    call MPI_TYPE_FREE(fsub_type, mpierr)
697    ABI_CHECK_MPI(mpierr,"MPI_TYPE_FREE")
698 
699    if (sc_mode==xmpio_collective) then
700      ! Collective read
701      if (icall <= my_ncalls) then
702        call MPI_FILE_READ_ALL(myfh, cbuffer(ptr), ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
703      else
704        ABI_MALLOC(dummy_cbuf,(subs_x))
705        call MPI_FILE_READ_ALL(myfh, dummy_cbuf, ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
706        ABI_FREE(dummy_cbuf)
707      end if
708      ABI_CHECK_MPI(mpierr,"FILE_READ_ALL")
709 
710    else
711      ! Individual read.
712      call MPI_FILE_READ(myfh, cbuffer(ptr), ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
713      ABI_CHECK_MPI(mpierr,"FILE_READ")
714    end if
715 
716  end do
717 
718  ABI_FREE(my_subsizes)
719  ABI_FREE(my_starts)
720  ABI_FREE(my_basead)
721 
722 end subroutine mpiotk_read_fsuba_dpc3D

m_mpiotk/mpiotk_read_fsuba_dpc4D [ Functions ]

[ Top ] [ m_mpiotk ] [ Functions ]

NAME

  mpiotk_read_fsuba_dpc4D

FUNCTION

  Reading a block of contiguous data stored in a 4D matrix.
  Data is placed within Fortran records. Target: complex data stored in a complex array.

NOTES

  The value of ierr should always be checked by the caller

INPUTS

  fh = MPI-IO file handler.
  offset =
  sizes(4)
  subsizes(4)
  starts(4)
  bufsz = dimension of cbuffer
  chunk_bsize =
  comm = MPI communicator
  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading.

 OUTPUTS
  cbuffer(bufsz)
  ierr=status error. A non zero value indicates that chunk_bsize is smaller that the fortran record
    and therefore bufsz has not been read.

PARENTS

      m_io_screening

CHILDREN

      mpi_file_read,mpi_file_read_all,mpi_file_set_view,mpi_type_free
      xmpi_max,xmpi_min,xmpio_create_fsubarray_4d

SOURCE

765 subroutine mpiotk_read_fsuba_dpc4D(fh,offset,sizes,subsizes,starts,bufsz,cbuffer,chunk_bsize,sc_mode,comm,ierr)
766 
767 
768 !This section has been created automatically by the script Abilint (TD).
769 !Do not modify the following lines by hand.
770 #undef ABI_FUNC
771 #define ABI_FUNC 'mpiotk_read_fsuba_dpc4D'
772 !End of the abilint section
773 
774  implicit none
775 
776 !Arguments ------------------------------------
777 !scalars
778  integer,intent(in) :: fh,comm,bufsz,sc_mode
779  integer,intent(out) :: ierr
780  integer(XMPI_OFFSET_KIND),intent(in) :: offset,chunk_bsize
781 !arrays
782  integer,intent(in) :: sizes(4),subsizes(4),starts(4)
783  complex(dpc),intent(out) :: cbuffer(bufsz)
784 
785 !Local variables-------------------------------
786 !scalars
787  integer :: mpierr,na2read,na_chunk,ptr,ncount,myfh
788  integer :: fsub_type,my_ncalls,icall,arest,ncalls
789  integer :: size_x,size_y,size_z,subs_x,subs_y,subs_z,subs_xyz,start_x,start_y,start_z,start_a,stop_a
790  integer(XMPI_OFFSET_KIND) :: my_offset,my_offpad
791  !character(len=500) :: msg
792 !arrays
793  integer :: call_subsizes(4),call_starts(4)
794  integer,allocatable :: my_basead(:),my_subsizes(:,:),my_starts(:,:)
795  complex(dpc),allocatable :: dummy_cbuf(:)
796 
797 !************************************************************************
798 
799  ! Workaround for XLF
800  myfh = fh
801 
802  if (bufsz < PRODUCT(subsizes) ) then
803    MSG_ERROR("bufsz is too small")
804  end if
805 
806  if (sc_mode==xmpio_single) then
807   ABI_CHECK(comm==xmpi_comm_self,"comm != xmpi_comm_self")
808  else if (sc_mode==xmpio_collective) then
809    continue
810  else
811    MSG_ERROR("Wrong sc_mode")
812  end if
813 
814  size_x  = sizes(1)
815  size_y  = sizes(2)
816  size_z  = sizes(3)
817 
818  subs_x  = subsizes(1)
819  subs_y  = subsizes(2)
820  subs_z  = subsizes(3)
821  subs_xyz = subs_x * subs_y * subs_z
822 
823  start_x = starts(1)
824  start_y = starts(2)
825  start_z = starts(3)
826  start_a = starts(4)
827  stop_a  = start_a + subsizes(4)-1 ! last column to read
828 
829  ! Read rows in blocks of size na_chunk:
830  ! MPI-IO crashes if we try to read data > 2Gb in a single call.
831  na2read = subsizes(4)
832  na_chunk = na2read
833  if ( (subs_xyz*na2read*xmpi_bsize_dpc) > chunk_bsize) then
834    na_chunk = chunk_bsize / (subs_xyz*xmpi_bsize_dpc)
835  end if
836 
837  call xmpi_min(na_chunk,ierr,comm,mpierr)
838  if (ierr == 0) then
839    ierr = 1
840    RETURN
841  end if
842  ierr = 0
843 
844  ! my_ncalls : number of read needed to fill my cbuffer.
845  ! ncalls    : max number of read in comm (needed for collective operations).
846  my_ncalls = na2read / na_chunk
847  arest = MOD(na2read, na_chunk)
848  if (arest /= 0) my_ncalls = my_ncalls + 1
849 
850  call xmpi_max(my_ncalls,ncalls,comm,mpierr)
851 
852  ! Compute arrays used to define the file view.
853  ABI_MALLOC(my_subsizes,(4,ncalls))
854  ABI_MALLOC(my_starts,(4,ncalls))
855  ABI_MALLOC(my_basead,(ncalls))
856 
857  do icall=1,my_ncalls
858 
859    if (icall*na_chunk <= na2read) then
860      my_subsizes(:,icall) = (/subs_x, subs_y, subs_z, na_chunk/)
861      my_starts(:,icall) = (/start_x, start_y, start_z, (icall-1) * na_chunk + start_a/)
862      my_basead(icall) = 1 + (icall-1)*subs_xyz*na_chunk
863    else
864      ! Two cases:
865      ! 1) na2read > na_chunk and not divisible by na2read
866      ! 2) na2read < na_chunk
867      my_subsizes(:,icall) = (/subs_x, subs_y, subs_z, arest/)
868      if (na2read >= na_chunk) then
869        my_starts(:,icall) = (/start_x, start_y, start_z, stop_a-arest+1/)
870        my_basead(icall) = 1 + (na2read-arest) * subs_xyz
871      else
872        my_starts(:,icall) = starts
873        my_basead(icall) = 1
874      end if
875    end if
876  end do
877  !write(std_out,*)" >>>> my_ncalls, ncalls, na2read, na_chunk ",my_ncalls,ncalls,na2read,na_chunk
878 
879  do icall=1,ncalls
880 
881    if (icall <= my_ncalls) then
882      call_subsizes = my_subsizes(:,icall)
883      call_starts   = my_starts(:,icall)
884      ptr           = my_basead(icall)
885    else
886      ! Fake values needed to call read_all collectively.
887      call_subsizes = (/subs_x, 1, 1, 1/)
888      call_starts   = starts
889    end if
890    ncount = PRODUCT(call_subsizes)
891    !
892    !write(std_out,*)"  icall,ptr, ncount, ",icall,ptr,ncount
893    !write(std_out,*)"  call_starts",call_starts
894    !write(std_out,*)"  call_subsizes",call_subsizes
895    !
896    ! Create subarry file view.
897    call xmpio_create_fsubarray_4D(sizes,call_subsizes,call_starts,MPI_DOUBLE_COMPLEX,&
898 &    fsub_type,my_offpad,mpierr)
899    ABI_CHECK_MPI(mpierr,"fsubarray_4D")
900 
901    ! Update the offset.
902    my_offset = offset + my_offpad
903 
904    call MPI_FILE_SET_VIEW(myfh, my_offset, MPI_BYTE, fsub_type, 'native', MPI_INFO_NULL, mpierr)
905    ABI_CHECK_MPI(mpierr,"SET_VIEW")
906 
907    call MPI_TYPE_FREE(fsub_type, mpierr)
908    ABI_CHECK_MPI(mpierr,"MPI_TYPE_FREE")
909 
910    if (sc_mode==xmpio_collective) then
911      ! Collective read
912      if (icall <= my_ncalls) then
913        call MPI_FILE_READ_ALL(myfh, cbuffer(ptr), ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
914      else
915        ABI_MALLOC(dummy_cbuf,(subs_x))
916        call MPI_FILE_READ_ALL(myfh, dummy_cbuf, ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
917        ABI_FREE(dummy_cbuf)
918      end if
919      ABI_CHECK_MPI(mpierr,"FILE_READ_ALL")
920    else
921      ! Individual read.
922      call MPI_FILE_READ(myfh, cbuffer(ptr), ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
923      ABI_CHECK_MPI(mpierr,"FILE_READ")
924    end if
925 
926  end do
927 
928  ABI_FREE(my_subsizes)
929  ABI_FREE(my_starts)
930  ABI_FREE(my_basead)
931 
932 end subroutine mpiotk_read_fsuba_dpc4D

m_mpiotk/mpiotk_write_fsuba_dp2D [ Functions ]

[ Top ] [ m_mpiotk ] [ Functions ]

NAME

  mpiotk_write_fsuba_dp2D

FUNCTION

  Write a block of contiguous data stored in a 2D matrix.
  Data is placed within Fortran records. Target: complex data stored in a real array.

NOTES

  The value of ierr should always be checked by the caller

INPUTS

  fh = MPI-IO file handler.
  offset =
  sizes(2)
  subsizes(2)
  starts(2)
  bufsz = dimension of buffer (takes into accout both real and imaginary part)
  chunk_bsize =
  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading.
  comm = MPI communicator

 OUTPUTS
  buffer(bufsz)
  ierr=status error. A non zero value indicates that chunk_bsize is smaller that the fortran record
    and therefore bufsz has not been written.

PARENTS

      m_wfk

CHILDREN

      mpi_file_read,mpi_file_read_all,mpi_file_set_view,mpi_type_free
      xmpi_max,xmpi_min,xmpio_create_fsubarray_4d

SOURCE

413 subroutine mpiotk_write_fsuba_dp2D(fh,offset,sizes,subsizes,starts,bufsz,buffer,chunk_bsize,sc_mode,comm,ierr)
414 
415 
416 !This section has been created automatically by the script Abilint (TD).
417 !Do not modify the following lines by hand.
418 #undef ABI_FUNC
419 #define ABI_FUNC 'mpiotk_write_fsuba_dp2D'
420 !End of the abilint section
421 
422  implicit none
423 
424 !Arguments ------------------------------------
425 !scalars
426  integer,intent(in) :: fh,comm,bufsz,sc_mode
427  integer,intent(out) :: ierr
428  integer(XMPI_OFFSET_KIND),intent(in) :: offset,chunk_bsize
429 !arrays
430  integer,intent(in) :: sizes(2),subsizes(2),starts(2)
431  real(dp),intent(in) :: buffer(bufsz)
432 
433 !Local variables ------------------------------
434 !scalars
435  integer :: mpierr,ptr,ncount,myfh
436  integer :: fsub_type,my_ncalls,icall,ncalls,subs_x
437  integer(XMPI_OFFSET_KIND) :: my_offset,my_offpad
438  !character(len=500) :: msg
439 !arrays
440  integer :: call_subsizes(2),call_starts(2)
441  integer,allocatable :: my_basead(:),my_subsizes(:,:),my_starts(:,:)
442 
443 !************************************************************************
444 
445  DBG_ENTER("COLL")
446 
447  ! Workaround for XLF
448  myfh = fh
449 
450  if (bufsz < 2 * PRODUCT(subsizes) ) then
451    MSG_ERROR("bufsz is too small")
452  end if
453 
454  if (sc_mode==xmpio_single) then
455    call setup_fsuba_dp2D(sizes,subsizes,starts,chunk_bsize,my_basead,my_subsizes,my_starts,my_ncalls,ncalls,xmpi_comm_self,ierr)
456  else if (sc_mode==xmpio_collective) then
457    call setup_fsuba_dp2D(sizes,subsizes,starts,chunk_bsize,my_basead,my_subsizes,my_starts,my_ncalls,ncalls,comm,ierr)
458  else
459    MSG_ERROR("Wrong sc_mode")
460  end if
461  if (ierr/=0) RETURN
462 
463  subs_x = subsizes(1)
464  do icall=1,ncalls
465 
466    if (icall <= my_ncalls) then
467      call_subsizes = my_subsizes(:,icall)
468      call_starts   = my_starts(:,icall)
469      ptr           = my_basead(icall)
470    else
471      ! Fake values needed to call write_all collectively.
472      call_subsizes = (/subs_x, 1/)
473      call_starts   = starts
474    end if
475    ncount = PRODUCT(call_subsizes)
476    !write(std_out,*)"  icall,ptr, ncount, ",icall,ptr,ncount
477    !write(std_out,*)"  call_starts",call_starts
478    !write(std_out,*)"  call_subsizes",call_subsizes
479 
480    ! Create subarry file view.
481    call xmpio_create_fsubarray_2D(sizes,call_subsizes,call_starts,MPI_DOUBLE_COMPLEX,fsub_type,my_offpad,mpierr)
482    ABI_CHECK_MPI(mpierr,"fsubarray_2D")
483 
484    ! Update the offset.
485    my_offset = offset + my_offpad
486 
487    call MPI_FILE_SET_VIEW(myfh, my_offset, MPI_BYTE, fsub_type, 'native', MPI_INFO_NULL, mpierr)
488    ABI_CHECK_MPI(mpierr,"SET_VIEW")
489 
490    call MPI_TYPE_FREE(fsub_type, mpierr)
491    ABI_CHECK_MPI(mpierr,"MPI_TYPE_FREE")
492 
493    if (sc_mode==xmpio_collective) then
494      ! Collective write
495      if (icall <= my_ncalls) then
496        call MPI_FILE_WRITE_ALL(myfh, buffer(ptr), ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
497      else
498        ! Re-write my first chunk of data.
499        call MPI_FILE_WRITE_ALL(myfh, buffer(1), ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
500      end if
501      ABI_CHECK_MPI(mpierr,"FILE_WRITE_ALL")
502 
503    else
504      ! Individual write.
505      call MPI_FILE_WRITE(myfh, buffer(ptr), ncount, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
506      ABI_CHECK_MPI(mpierr,"FILE_WRITE")
507    end if
508 
509  end do
510 
511  ABI_FREE(my_subsizes)
512  ABI_FREE(my_starts)
513  ABI_FREE(my_basead)
514 
515  DBG_EXIT("COLL")
516 
517 end subroutine mpiotk_write_fsuba_dp2D

m_mpiotk/no_mpiotk [ Functions ]

[ Top ] [ m_mpiotk ] [ Functions ]

NAME

  no_mpiotk

FUNCTION

   Empty placeholder.

PARENTS

CHILDREN

      mpi_file_read,mpi_file_read_all,mpi_file_set_view,mpi_type_free
      xmpi_max,xmpi_min,xmpio_create_fsubarray_4d

SOURCE

84 subroutine no_mpiotk()
85 
86 
87 !This section has been created automatically by the script Abilint (TD).
88 !Do not modify the following lines by hand.
89 #undef ABI_FUNC
90 #define ABI_FUNC 'no_mpiotk'
91 !End of the abilint section
92 
93  implicit none
94 
95 ! *************************************************************************
96 
97 end subroutine no_mpiotk

m_mpiotk/setup_fsuba_dp2D [ Functions ]

[ Top ] [ m_mpiotk ] [ Functions ]

NAME

  setup_fsuba_dp2D

FUNCTION

  Setup tables used in (read|write) a 2D array stored in a Fortran file

NOTES

  The value of ierr should always be checked by the caller

INPUTS

  sizes(2)
  subsizes(2)
  starts(2)
  chunk_bsize =
  comm = MPI communicator

 OUTPUTS
  my_basead(:)
  my_subsizes(:,:)
  my_starts(:,:)
  ierr=status error. A non zero value indicates that chunk_bsize is smaller that the fortran record
    and therefore bufsz has not been read.

PARENTS

      m_mpiotk

CHILDREN

      mpi_file_read,mpi_file_read_all,mpi_file_set_view,mpi_type_free
      xmpi_max,xmpi_min,xmpio_create_fsubarray_4d

SOURCE

135 subroutine setup_fsuba_dp2D(sizes,subsizes,starts,chunk_bsize,&
136 &  my_basead,my_subsizes,my_starts,my_ncalls,ncalls,comm,ierr)
137 
138 
139 !This section has been created automatically by the script Abilint (TD).
140 !Do not modify the following lines by hand.
141 #undef ABI_FUNC
142 #define ABI_FUNC 'setup_fsuba_dp2D'
143 !End of the abilint section
144 
145  implicit none
146 
147 !Arguments ------------------------------------
148 !scalars
149  integer,intent(in) :: comm
150  integer,intent(out) :: ncalls,my_ncalls,ierr
151  integer(XMPI_OFFSET_KIND),intent(in) :: chunk_bsize
152 !arrays
153  integer,intent(in) :: sizes(2),subsizes(2),starts(2)
154  integer,allocatable,intent(out) :: my_basead(:),my_subsizes(:,:),my_starts(:,:)
155 
156 !Local variables ------------------------------
157 !scalars
158  integer :: mpierr,ny2read,ny_chunk,icall,yrest
159  integer :: size_x,subs_x,start_x,start_y,stop_y
160  !character(len=500) :: msg
161 
162 !************************************************************************
163 
164  size_x  = sizes(1)
165  subs_x  = subsizes(1)
166  start_x = starts(1)
167  start_y = starts(2)
168  stop_y  = start_y + subsizes(2)-1 ! last column to read
169  !
170  ! Read rows in blocks of size ny_chunk:
171  ! MPI-IO crashes if we try to read data > 2Gb in a single call.
172  ny2read = subsizes(2)
173  ny_chunk = ny2read
174  if ((2*subs_x*ny2read*xmpi_bsize_dp) > chunk_bsize) then
175    ny_chunk = chunk_bsize / (2*subs_x*xmpi_bsize_dp)
176    !if (ny_chunk == 0) ny_chunk = 50
177  end if
178 
179  call xmpi_min(ny_chunk,ierr,comm,mpierr)
180  if (ierr == 0) then
181    ierr = 1
182    RETURN
183  end if
184  ierr = 0
185  !
186  ! my_ncalls : number of read needed to fill my buffer.
187  ! ncalls    : max number of read in comm (needed for collective operations).
188  my_ncalls = ny2read / ny_chunk
189  yrest = MOD(ny2read, ny_chunk)
190  if (yrest /= 0) my_ncalls = my_ncalls + 1
191 
192  call xmpi_max(my_ncalls,ncalls,comm,mpierr)
193  !
194  ! Compute arrays used to define the file view.
195  ABI_MALLOC(my_subsizes,(2,ncalls))
196  ABI_MALLOC(my_starts,(2,ncalls))
197  ABI_MALLOC(my_basead,(ncalls))
198 
199  do icall=1,my_ncalls
200 
201    if (icall*ny_chunk <= ny2read) then
202      my_subsizes(:,icall) = (/subs_x, ny_chunk/)
203      my_starts(:,icall) = (/start_x, (icall-1) * ny_chunk + start_y/)
204      my_basead(icall) = 1 + 2*(icall-1)*ny_chunk*subs_x ! 2 accounts for real and imag part.
205    else
206      ! Two cases:
207      ! 1) ny2read > ny_chunk and not divisible by ny2read
208      ! 2) ny2read < ny_chunk
209      my_subsizes(:,icall) = (/subs_x, yrest/)
210      if (ny2read >= ny_chunk) then
211        my_starts(:,icall) = (/start_x, stop_y-yrest+1/)
212        my_basead(icall) = 1 + 2 * (ny2read-yrest) * subs_x ! 2 accounts for real and imag part.
213      else
214        my_starts(:,icall) = starts
215        my_basead(icall) = 1
216      end if
217    end if
218  end do
219  !write(std_out,*)" >>>> my_ncalls, ncalls, ny2read, ny_chunk ",my_ncalls,ncalls,ny2read,ny_chunk
220 
221 end subroutine setup_fsuba_dp2D