TABLE OF CONTENTS
- ABINIT/m_wffile
- ABINIT/WffOffset
- ABINIT/xderiveRead_dp1d
- m_wffile/clsopn
- m_wffile/getRecordMarkerLength_wffile
- m_wffile/rwRecordMarker
- m_wffile/wff_ireadf90
- m_wffile/wff_usef90
- m_wffile/WffClose
- m_wffile/WffDelete
- m_wffile/wffile_type
- m_wffile/WffKg
- m_wffile/WffOpen
- m_wffile/WffReadDataRec_dp1d
- m_wffile/WffReadDataRec_dp2d
- m_wffile/WffReadNpwRec
- m_wffile/WffReadSkipRec
- m_wffile/WffReadWrite_mpio
- m_wffile/WffWriteDataRec_dp1d
- m_wffile/WffWriteDataRec_dp2d
- m_wffile/WffWriteDataRec_int2d
- m_wffile/WffWriteNpwRec
- m_wffile/xdefineOff
- m_wffile/xderiveRead_dp
- m_wffile/xderiveRead_dp2d
- m_wffile/xderiveRead_dp2d_displ
- m_wffile/xderiveRead_int
- m_wffile/xderiveRead_int1d
- m_wffile/xderiveRead_int2d
- m_wffile/xderiveRead_int2d_displ
- m_wffile/xderiveReadVal_char
- m_wffile/xderiveRRecEnd
- m_wffile/xderiveRRecInit
- m_wffile/xderiveWRecEnd
- m_wffile/xderiveWRecInit
- m_wffile/xderiveWrite_char
- m_wffile/xderiveWrite_dp
- m_wffile/xderiveWrite_dp1d
- m_wffile/xderiveWrite_dp2d
- m_wffile/xderiveWrite_dp2d_displ
- m_wffile/xderiveWrite_dp2d_seq
- m_wffile/xderiveWrite_int
- m_wffile/xderiveWrite_int1d
- m_wffile/xderiveWrite_int2d
- m_wffile/xderiveWrite_int2d_displ
- m_wffile/xmoveOff
- m_wffile/xmpi_read_int2d
- m_wffile/xnullifyOff
ABINIT/m_wffile [ Modules ]
NAME
m_wffile
FUNCTION
This module provides the definition of the wffile_type used to WF file data. As the type contains MPI-dependent fields, it has to be declared in a MPI-managed directory.
COPYRIGHT
Copyright (C) 2009-2024 ABINIT group (MT,MB,MVer,ZL,MD) 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.
NOTES
wffile_type: a handler for dealing with the IO of a wavefunction file
SOURCE
21 #if defined HAVE_CONFIG_H 22 #include "config.h" 23 #endif 24 25 #include "abi_common.h" 26 27 MODULE m_wffile 28 29 use defs_basis 30 use m_errors 31 use m_abicore 32 use m_xmpi 33 #ifdef HAVE_MPI2 34 use mpi 35 #endif 36 use m_nctk 37 #ifdef HAVE_NETCDF 38 use netcdf 39 #endif 40 41 use defs_abitypes, only : MPI_Type 42 use m_io_tools, only : mvrecord, open_file 43 use m_fstrings, only : toupper, endswith, sjoin 44 use iso_c_binding 45 46 implicit none 47 48 private 49 50 #ifdef HAVE_MPI1 51 include 'mpif.h' 52 #endif 53 54 #define DEV_DEBUG_THIS 0 55 56 !public procedures. 57 public :: WffOpen 58 public :: wffclose 59 public :: wffdelete 60 public :: wffkg 61 public :: wffoffset 62 public :: wffreaddatarec ! Generic subroutines to read data in one record of a wavefunction file 63 public :: wffreadnpwrec 64 public :: wffreadskiprec 65 public :: wffreadwrite_mpio 66 public :: wffwritedatarec ! Generic subroutines to write data in one record of a wavefunction file 67 public :: wffwritenpwrec 68 public :: xderiveread ! Generic subroutines to read wf files. 69 public :: xmpi_read_int2d 70 public :: xderivewrite 71 72 public :: getRecordMarkerLength_wffile 73 public :: xnullifyOff 74 public :: xmoveOff 75 public :: xderiveWRecEnd 76 public :: xderiveWRecInit 77 public :: xderiveRRecEnd 78 public :: xderiveRRecInit 79 #if defined HAVE_MPI_IO 80 public :: rwRecordMarker 81 #endif 82 public :: clsopn 83 public :: wff_usef90 84 public :: xdefineOff
ABINIT/WffOffset [ Functions ]
NAME
WffOffset
FUNCTION
Tool to manage WF file in the MPI/IO case : broadcast the offset of the first k-point data block
INPUTS
wff <type(wffile_type)> = structured info about the wavefunction file sender = id of the sender spaceComm = id of the space communicator handler
OUTPUT
ier = error code returned by the MPI call
SOURCE
1275 subroutine WffOffset(wff,sender,spaceComm,ier) 1276 1277 !Arguments ------------------------------------ 1278 type(wffile_type),intent(inout) :: wff 1279 integer ,intent(inout) :: sender 1280 integer ,intent(in) :: spaceComm 1281 integer ,intent(out) :: ier 1282 1283 !Local variables ------------------------------ 1284 #if defined HAVE_MPI_IO 1285 integer :: icom 1286 integer(kind=MPI_OFFSET_KIND) :: off(1) 1287 #endif 1288 1289 ! ********************************************************************* 1290 1291 #if defined HAVE_MPI_IO 1292 if (wff%iomode == IO_MODE_MPI) then 1293 call xmpi_max(sender,icom,spaceComm,ier) 1294 if (icom>=0)then 1295 off(1)=wff%offwff 1296 call MPI_BCAST(off,1,wff%offset_mpi_type,icom,spaceComm,ier) 1297 wff%offwff=off(1) 1298 end if 1299 end if ! iomode 1300 #else 1301 ier = 0 1302 ABI_UNUSED((/wff%iomode,sender,spaceComm/)) 1303 #endif 1304 1305 end subroutine WffOffset
ABINIT/xderiveRead_dp1d [ Functions ]
NAME
xderiveRead_dp1d
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: one-dimensional double precision arrays.
INPUTS
n1= first dimension of the array spaceComm= MPI communicator
OUTPUT
ierr= exit status, a non-zero value meaning there is an error xval= data buffer array
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2436 subroutine xderiveRead_dp1d(wff,xval,n1,spaceComm,ierr) 2437 2438 !Arguments ------------------------------------ 2439 type(wffile_type),intent(inout) :: wff 2440 integer,intent(in) :: n1,spaceComm 2441 integer,intent(out) :: ierr 2442 real(dp),intent(out) :: xval(:) 2443 2444 !Local variables------------------------------- 2445 #if defined HAVE_MPI_IO 2446 integer(kind=MPI_OFFSET_KIND) :: delim_record,nboct,dispoct(1),posit,totoct 2447 integer :: statux(MPI_STATUS_SIZE) 2448 #endif 2449 character(len=500) :: msg 2450 2451 !********************************************************************* 2452 2453 xval(:)=zero ; ierr=0 ! Initialization, for the compiler 2454 if(.false.)write(std_out,*)wff%me,n1,spaceComm 2455 2456 #if defined HAVE_MPI_IO 2457 nboct = wff%nbOct_dp * n1 2458 posit = wff%offwff 2459 delim_record = posit - wff%off_recs + wff%lght_recs - wff%nbOct_recMarker 2460 2461 if (delim_record >= nboct) then 2462 ! Compute offset for local part 2463 ! dispoct = sum (nboct, rank=0..me) 2464 if (spaceComm/=MPI_COMM_SELF) then 2465 call MPI_SCAN([nboct],dispoct,1,wff%offset_mpi_type,MPI_SUM,spaceComm,ierr) 2466 posit = posit + dispoct(1) - nboct 2467 end if 2468 2469 call MPI_FILE_READ_AT(wff%fhwff,posit,xval,n1,MPI_DOUBLE_PRECISION,statux,ierr) 2470 2471 ! get the total number of bits wrote by processors 2472 if (spaceComm/=MPI_COMM_SELF) then 2473 call xmpi_max(dispoct(1),totoct,spaceComm,ierr) 2474 !call MPI_ALLREDUCE(dispoct,totoct,1,wff%offset_mpi_type,MPI_MAX,spaceComm,ierr) 2475 else 2476 totoct=nboct 2477 end if 2478 else 2479 ierr = 1 2480 nboct = 0 2481 totoct = 0 2482 end if 2483 2484 !new offset 2485 wff%offwff=wff%offwff + totoct 2486 return 2487 #endif 2488 2489 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 2490 ABI_WARNING(msg) 2491 2492 end subroutine xderiveRead_dp1d
m_wffile/clsopn [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
clsopn
FUNCTION
Close wavefunction file (provided its access is standard F90 IO), then reopen the same. Uses fortran inquire statement to reopen with same characteristics.
INPUTS
wff=number of unit to which on which file is already opened.
SIDE EFFECTS
SOURCE
835 subroutine clsopn(wff) 836 837 !Arguments ------------------------------------ 838 !scalars 839 type(wffile_type),intent(inout) :: wff 840 841 !Local variables------------------------------- 842 !scalars 843 integer :: ios,unit 844 logical :: nmd,od 845 character(len=11) :: fm 846 character(len=500) :: message 847 character(len=fnlen) :: filnam 848 849 ! ************************************************************************* 850 851 if ( ANY(wff%iomode==(/IO_MODE_FORTRAN_MASTER,IO_MODE_FORTRAN/) ))then 852 853 unit=wff%unwff 854 inquire (unit=unit,iostat=ios,opened=od,name=filnam,form=fm,named=nmd) 855 856 ! ios is a status specifier. If an error condition exists, 857 ! ios is assigned a processor-dependent value > 0. 858 if (ios/=0) then 859 write(message, '(/,a,/,a,i8,a,i8,/,a,/,a,/,a)' ) & 860 & ' clsopn : ERROR -',& 861 & ' Attempt to inquire about unit=',unit,& 862 & ' indicates error condition iostat=',ios,& 863 & ' May be due to temporary problem with file, disks or network.',& 864 & ' Action: check whether there might be some external problem,',& 865 & ' then resubmit.' 866 ABI_ERROR(message) 867 868 ! od is a logical variable which is set to true if the specified 869 ! unit is connected to a file; otherwise it is set to false. 870 else if (.not.od) then 871 write(message, '(/,a,/,a,i8,/,a,/,a,/,a,/,a)' ) & 872 & ' clsopn : ERROR -',& 873 & ' Tried to inquire about unit',unit,& 874 & ' and found it not connected to a file.',& 875 & ' May be due to temporary problem with file, disks or network.',& 876 & ' Action: check whether there might be some external problem,',& 877 & ' then resubmit.' 878 ABI_ERROR(message) 879 880 ! nmd is a logical variable assigned the value true if the file 881 ! has a name; otherwise false. A scratch file is not named. 882 else if (.not.nmd) then 883 884 ! No action for the time being. Possibility to debug. 885 886 else 887 888 ! May now close the file and then reopen it 889 ! (file is already opened according to above checks) 890 close (unit=unit) 891 open (unit=unit,file=filnam,form=fm,status='old') !VALGRIND complains filnam is just a few thousand bytes inside a block of 8300 892 893 end if 894 895 else if (wff%iomode == IO_MODE_MPI) then 896 call xnullifyOff(wff) 897 else if (wff%iomode == IO_MODE_ETSF) then 898 ! We do nothing, ETSF access already not being sequential. 899 end if 900 901 end subroutine clsopn
m_wffile/getRecordMarkerLength_wffile [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
getRecordMarkerLength_wffile
FUNCTION
Get the record marker length of the FORTRAN header of a file to access it in MPI/IO. This routine assumes that the header has been written (and flushed) in the file.
SIDE EFFECTS
wff=<type(wffile_type)>=structured info for reading/writing the wavefunctions only%nbOct_recMarker is changed
SOURCE
242 subroutine getRecordMarkerLength_wffile(wff) 243 244 !Arguments ------------------------------------ 245 !scalars 246 type(wffile_type),intent(inout) :: wff 247 248 !Local variables------------------------------- 249 #if defined HAVE_MPI_IO 250 !scalars 251 integer :: ierr,ii,iimax 252 integer(kind=MPI_OFFSET_KIND) :: posit,rml 253 character(len=500) :: msg 254 !arrays 255 integer :: headform(1),statux(MPI_STATUS_SIZE) 256 integer(kind=MPI_OFFSET_KIND) :: off(1) 257 #endif 258 259 !************************************************************************ 260 261 #ifdef DEV_DEBUG_THIS 262 return 263 ! Already done in WffOpen 264 #endif 265 266 #if defined HAVE_MPI_IO 267 268 if (wff%nbOct_recMarker>0) return 269 270 !wff%nbOct_recMarker=4;return 271 !call flush(wff%unwff) 272 !call MPI_FILE_SYNC(wff%fhwff,ierr) 273 274 !Only master do that 275 ierr=0 276 if (wff%master==wff%me) then 277 278 ! Define number of INTEGER types to be tested 279 #if defined HAVE_FC_INT_QUAD 280 iimax=4 281 #else 282 iimax=3 283 #endif 284 285 ! Try to read headform 286 rml=-1;ii=0 287 do while (wff%nbOct_recMarker<=0.and.ii<iimax) 288 ii=ii+1 289 if (ii==1) rml=4 290 if (ii==2) rml=8 291 if (ii==3) rml=2 292 if (ii==4) rml=16 293 posit=rml+6*wff%nbOct_ch 294 call MPI_FILE_READ_AT(wff%fhwff,posit,headform,1,MPI_INTEGER,statux,ierr) 295 if (ierr==MPI_SUCCESS) then 296 if (headform(1)==wff%headform) wff%nbOct_recMarker=rml 297 end if 298 end do 299 300 if (ierr/=MPI_SUCCESS) then 301 ABI_BUG("Header problem") 302 end if 303 304 if (ii==iimax.and.wff%nbOct_recMarker<=0) then 305 ! if (iimax>=4) then 306 ! write(msg,'(3a)') & 307 ! & ' Your architecture is not able to handle 16, 8, 4 or 2-bytes FORTRAN file record markers !',ch10,& 308 ! & ' You cannot use ABINIT and MPI/IO.' 309 ! else 310 ! write(msg,'(3a)') & 311 ! & ' Your architecture is not able to handle 8, 4 or 2-bytes FORTRAN file record markers !',ch10,& 312 ! & ' You cannot use ABINIT and MPI/IO.' 313 ! end if 314 write(msg,'(13a)') & 315 & ' Error during FORTRAN file record marker detection:',ch10,& 316 & ' It was not possible to read/write a small file!',ch10,& 317 & ' ACTION: check your access permissions to the file system.',ch10,& 318 & ' Common sources of this problem:',ch10,& 319 & ' - Quota limit exceeded,',ch10,& 320 & ' - R/W incorrect permissions,',ch10,& 321 & ' - WFK file requested as input (irdwfk=1/getwfk=1) but not existing ...' 322 ABI_ERROR(msg) 323 else 324 write(msg,'(a,i0)') & 325 & ' MPI/IO accessing FORTRAN file header: detected record mark length=',wff%nbOct_recMarker 326 ABI_COMMENT(msg) 327 end if 328 329 end if ! me=master 330 331 !Broadcast record marker length 332 if (wff%spaceComm/=MPI_COMM_SELF) then 333 off(1)=wff%nbOct_recMarker 334 call MPI_BCAST(off,1,wff%offset_mpi_type,wff%master,wff%spaceComm,ierr) 335 wff%nbOct_recMarker=off(1) 336 end if 337 338 !Select MPI datatype for markers 339 if (wff%nbOct_recMarker==4) then 340 wff%marker_mpi_type=MPI_INTEGER4 341 else if (wff%nbOct_recMarker==8) then 342 wff%marker_mpi_type=MPI_INTEGER8 343 #if defined HAVE_FC_INT_QUAD && defined HAVE_MPI_INTEGER16 344 else if (wff%nbOct_recMarker==16) then 345 wff%marker_mpi_type=MPI_INTEGER16 346 #endif 347 else if (wff%nbOct_recMarker==2) then 348 wff%marker_mpi_type=MPI_INTEGER2 349 end if 350 351 #endif 352 353 RETURN 354 ABI_UNUSED(wff%me) 355 356 end subroutine getRecordMarkerLength_wffile
m_wffile/rwRecordMarker [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
rwRecordMarker
FUNCTION
Read/Write a record marker in a FORTRAN file at a given file pointer position. This is needed to access data in a FORTRAN file with MPI/IO.
INPUTS
option=1 for reading by current proc 2 for writing by current proc 3 for reading by all procs 4 for writing by all procs posit= position of the MPI/IO file pointer wff=<type(wffile_type)>=structured info for reading/writing Use here only: wff%fhwff= handle of the MPI/IO file wff%nbOct_recMarker= length of Fortran record markers
OUTPUT
ierr= error code
SIDE EFFECTS
posit= position of the MPI/IO file pointer updated after the reading (with the length of the record) recordmarker= content of the record marker
SOURCE
390 #if defined HAVE_MPI_IO 391 392 subroutine rwRecordMarker(option,posit,recordmarker,wff,ierr) 393 394 !Arguments ------------------------------------ 395 !scalars 396 integer,intent(in) :: option 397 integer(kind=MPI_OFFSET_KIND),intent(inout) :: posit,recordmarker 398 integer,intent(out) :: ierr 399 type(wffile_type),intent(inout) :: wff 400 401 !Local variables------------------------------- 402 !scalars 403 integer(kind=2) :: delim_record2(1) 404 integer(kind=4) :: delim_record4(1) 405 integer(kind=8) :: delim_record8(1) 406 #if defined HAVE_FC_INT_QUAD 407 integer(kind=16) :: delim_record16(1) 408 #endif 409 !character(len=500) :: msg 410 !arrays 411 integer :: statux(MPI_STATUS_SIZE) 412 413 !************************************************************************ 414 415 ierr=0 416 417 if (option==1) then 418 if (wff%nbOct_recMarker==4) then 419 call MPI_FILE_READ_AT(wff%fhwff,posit,delim_record4,1,wff%marker_mpi_type,statux,ierr) 420 recordmarker = delim_record4(1) 421 else if (wff%nbOct_recMarker==8) then 422 call MPI_FILE_READ_AT(wff%fhwff,posit,delim_record8,1,wff%marker_mpi_type,statux,ierr) 423 recordmarker = delim_record8(1) 424 #if defined HAVE_FC_INT_QUAD 425 else if (wff%nbOct_recMarker==16) then 426 call MPI_FILE_READ_AT(wff%fhwff,posit,delim_record16,1,wff%marker_mpi_type,statux,ierr) 427 recordmarker = delim_record16(1) 428 #endif 429 else if (wff%nbOct_recMarker==2) then 430 call MPI_FILE_READ_AT(wff%fhwff,posit,delim_record2,1,wff%marker_mpi_type,statux,ierr) 431 recordmarker = delim_record2(1) 432 else 433 ABI_BUG('Wrong record marker length!') 434 end if 435 436 else if (option==2) then 437 if (wff%nbOct_recMarker==4) then 438 delim_record4(1) = recordmarker 439 call MPI_FILE_WRITE_AT(wff%fhwff,posit,delim_record4,1,wff%marker_mpi_type,statux,ierr) 440 else if (wff%nbOct_recMarker==8) then 441 delim_record8(1) = recordmarker 442 call MPI_FILE_WRITE_AT(wff%fhwff,posit,delim_record8,1,wff%marker_mpi_type,statux,ierr) 443 #if defined HAVE_FC_INT_QUAD 444 else if (wff%nbOct_recMarker==16) then 445 delim_record16(1) = recordmarker 446 call MPI_FILE_WRITE_AT(wff%fhwff,posit,delim_record16,1,wff%marker_mpi_type,statux,ierr) 447 #endif 448 else if (wff%nbOct_recMarker==2) then 449 delim_record2(1) = recordmarker 450 call MPI_FILE_WRITE_AT(wff%fhwff,posit,delim_record2,1,wff%marker_mpi_type,statux,ierr) 451 else 452 ABI_BUG('Wrong record marker length!') 453 end if 454 455 else if (option==3) then 456 if (wff%nbOct_recMarker==4) then 457 call MPI_FILE_READ_AT_ALL(wff%fhwff,posit,delim_record4,1,wff%marker_mpi_type,statux,ierr) 458 recordmarker = delim_record4(1) 459 else if (wff%nbOct_recMarker==8) then 460 call MPI_FILE_READ_AT_ALL(wff%fhwff,posit,delim_record8,1,wff%marker_mpi_type,statux,ierr) 461 recordmarker = delim_record8(1) 462 #if defined HAVE_FC_INT_QUAD 463 else if (wff%nbOct_recMarker==16) then 464 call MPI_FILE_READ_AT_ALL(wff%fhwff,posit,delim_record16,1,wff%marker_mpi_type,statux,ierr) 465 recordmarker = delim_record16(1) 466 #endif 467 else if (wff%nbOct_recMarker==2) then 468 call MPI_FILE_READ_AT_ALL(wff%fhwff,posit,delim_record2,1,wff%marker_mpi_type,statux,ierr) 469 recordmarker = delim_record2(1) 470 else 471 ABI_BUG('Wrong record marker length !') 472 end if 473 474 else if (option==4) then 475 if (wff%nbOct_recMarker==4) then 476 delim_record4(1) = recordmarker 477 call MPI_FILE_WRITE_AT_ALL(wff%fhwff,posit,delim_record4,1,wff%marker_mpi_type,statux,ierr) 478 else if (wff%nbOct_recMarker==8) then 479 delim_record8(1) = recordmarker 480 call MPI_FILE_WRITE_AT_ALL(wff%fhwff,posit,delim_record8,1,wff%marker_mpi_type,statux,ierr) 481 #if defined HAVE_FC_INT_QUAD 482 else if (wff%nbOct_recMarker==16) then 483 delim_record16(1) = recordmarker 484 call MPI_FILE_WRITE_AT_ALL(wff%fhwff,posit,delim_record16,1,wff%marker_mpi_type,statux,ierr) 485 #endif 486 else if (wff%nbOct_recMarker==2) then 487 delim_record2(1) = recordmarker 488 call MPI_FILE_WRITE_AT_ALL(wff%fhwff,posit,delim_record2,1,wff%marker_mpi_type,statux,ierr) 489 else 490 ABI_BUG('Wrong record marker length!') 491 end if 492 493 else 494 ABI_BUG('Wrong value for option!') 495 end if 496 497 posit = posit + recordmarker + 2*wff%nbOct_recMarker 498 499 end subroutine rwRecordMarker 500 #endif
m_wffile/wff_ireadf90 [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
wff_ireadf90
FUNCTION
1 if a Fortran file is going to be read by this node, 0 otherwise.
INPUTS
SOURCE
945 function wff_ireadf90(wff) 946 947 !Arguments ------------------------------------ 948 !scalars 949 integer :: wff_ireadf90 950 type(wffile_type),intent(in) :: wff 951 952 ! ************************************************************************* 953 954 wff_ireadf90=0 955 if (wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) wff_ireadf90=1 956 957 end function wff_ireadf90
m_wffile/wff_usef90 [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
wff_usef90
FUNCTION
1 if a Fortran file is going to be read by this node, 0 otherwise.
INPUTS
SOURCE
917 function wff_usef90(wff) 918 919 !Arguments ------------------------------------ 920 !scalars 921 integer :: wff_usef90 922 type(wffile_type),intent(in) :: wff 923 924 ! ************************************************************************* 925 926 wff_usef90=0 927 if (wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode ==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) wff_usef90=1 928 929 end function wff_usef90
m_wffile/WffClose [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffClose
FUNCTION
This subroutine closes a Wf file.
INPUTS
wff= structured info about the wavefunction file
OUTPUT
ierr=error code
SOURCE
1136 subroutine WffClose(wff,ier) 1137 1138 !Arguments ------------------------------------ 1139 type(wffile_type), intent(inout) :: wff 1140 integer, intent(out) :: ier 1141 1142 ! ************************************************************************* 1143 1144 ier=0 1145 if(wff%iomode==IO_MODE_FORTRAN) then ! All processors see a local file 1146 close(unit=wff%unwff) 1147 1148 #ifdef HAVE_NETCDF 1149 else if(wff%iomode == IO_MODE_ETSF)then 1150 NCF_CHECK(nf90_close(wff%unwff)) 1151 #endif 1152 1153 else if(wff%iomode==IO_MODE_FORTRAN_MASTER)then ! Only the master processor see a local file 1154 if(wff%master==wff%me) close (unit=wff%unwff) ! VALGRIND complains buf points to uninitialized bytes 1155 1156 #if defined HAVE_MPI_IO 1157 else if(wff%iomode==IO_MODE_MPI)then 1158 call MPI_FILE_CLOSE(wff%fhwff,ier) 1159 if (wff%master==wff%me ) close(unit=wff%unwff) 1160 wff%offwff=0;wff%off_recs=0;wff%lght_recs=0 1161 wff%nbOct_recMarker=-1 1162 wff%kgwff=-1 1163 #endif 1164 1165 end if 1166 1167 end subroutine WffClose
m_wffile/WffDelete [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffDelete
FUNCTION
This subroutine closes a Wf file, and delete it.
INPUTS
wff= structured info about the wavefunction file
OUTPUT
ierr=error code
SOURCE
1187 subroutine WffDelete(wff,ier) 1188 1189 !Arguments ------------------------------------ 1190 type(wffile_type),intent(inout) :: wff 1191 integer, intent(out) :: ier 1192 1193 ! ************************************************************************* 1194 1195 ier=0 1196 if (wff%iomode==IO_MODE_FORTRAN) then ! All processors see a local file 1197 close(unit=wff%unwff,status='delete') 1198 1199 else if (wff%iomode==IO_MODE_FORTRAN_MASTER)then ! Only the master processor see a local file 1200 if (wff%master==wff%me) close (unit=wff%unwff,status='delete') 1201 1202 1203 else if (wff%iomode==IO_MODE_MPI)then 1204 #if defined HAVE_MPI_IO 1205 if ( wff%fhwff /= -1 )then 1206 call MPI_FILE_CLOSE(wff%fhwff,ier) 1207 end if 1208 if (wff%master==wff%me ) then 1209 close(unit=wff%unwff,status='delete') 1210 wff%fhwff = -1 1211 end if 1212 wff%offwff=0;wff%off_recs=0;wff%lght_recs=0 1213 wff%nbOct_recMarker=-1 1214 wff%kgwff=-1 1215 #endif 1216 end if 1217 1218 end subroutine WffDelete
m_wffile/wffile_type [ Types ]
[ Top ] [ m_wffile ] [ Types ]
NAME
wffile_type
FUNCTION
This structure datatype is a handler for dealing with the IO of a wavefunction file. It contains, among other things, the method of access to the file (standard F90 read/write, or NetCDF call, or MPI IO), the unit number if applicable, the filename, the information on the parallelism, etc ...
SOURCE
142 type, public :: wffile_type 143 144 ! WARNING : if you modify this datatype, please check there there is no creation/destruction/copy routine, 145 ! declared in another part of ABINIT, that might need to take into account your modification. 146 147 ! Integer scalar 148 integer :: unwff 149 ! unwff unit number of unformatted wavefunction disk file 150 151 integer :: iomode 152 ! Method to access the wavefunction file 153 ! IO_MODE_FORTRAN for usual Fortran IO routines 154 ! IO_MODE_FORTRAN_MASTER if usual Fortran IO routines, but only the master node in the parallel case 155 ! IO_MODE_MPI if MPI/IO routines (this access method is only available in parallel) 156 ! IO_MODE_NETCDF if NetCDF routines (obsolete, do not use) 157 ! IO_MODE_ETSF, NetCDF format read via etsf-io. 158 159 integer :: formwff 160 ! formwff=format of the eigenvalues 161 ! -1 => not used 162 ! 0 => vector of eigenvalues 163 ! 1 => hermitian matrix of eigenvalues 164 165 integer :: headform 166 ! headform=format of the header 167 168 integer :: kgwff 169 ! kgwff if 1 , read or write kg_k ; if 0, do not care about kg_k 170 171 ! Character 172 character(len=fnlen) :: fname 173 ! filename (if available) 174 175 ! In case of MPI parallel use 176 integer :: master 177 ! index of the processor master of the IO procedure when the WffOpen call is issued 178 179 integer :: me 180 ! index of my processor in the spaceComm communicator 181 182 integer :: me_mpiio 183 ! index of my processor in the spaceComm_mpiio communicator 184 185 integer :: nproc 186 ! number of processors that will have access to the file 187 188 integer :: spaceComm 189 ! space communicator for the standard FORTRAN access to the file 190 191 integer :: spaceComm_mpiio 192 ! space communicator for the MPI/IO access to the file 193 194 ! In case of MPI/IO : additional information 195 integer :: fhwff 196 ! file handle used to access the file with MPI/IO. 197 198 integer(kind=XMPI_OFFSET_KIND) :: nbOct_int,nbOct_dp,nbOct_ch 199 ! nbOct_int byte number of int value 200 ! nbOct_dp byte number of dp value 201 ! nbOct_ch byte number of character value 202 203 integer(kind=XMPI_OFFSET_KIND) :: nbOct_recMarker 204 ! byte number of Fortran file record markers 205 206 integer(kind=XMPI_OFFSET_KIND) :: lght_recs 207 ! length of record 208 209 integer :: marker_mpi_type 210 ! MPI Datatype for Fortran record markers 211 212 integer(kind=XMPI_OFFSET_KIND) :: offwff,off_recs 213 ! offwff offset position of unformatted wavefunction disk file 214 ! off_recs offset position of start record 215 ! (used in parallel MPI-IO) 216 217 integer :: offset_mpi_type 218 ! MPI Datatype for INTEGER(kind=MPI_OFFSET_KIND) 219 220 end type wffile_type 221 222 223 CONTAINS
m_wffile/WffKg [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffKg
FUNCTION
Check kgwff to manage WF file in the MPI/IO case
INPUTS
wff <type(wffile_type)> = structured info about the wavefunction file optkg= if 1 , read or write kg_k ; if 0,do not care about kg_k in rwwf
OUTPUT
SOURCE
1238 subroutine WffKg(wff,optkg) 1239 1240 !Arguments ------------------------------------ 1241 type(wffile_type),intent(inout) :: wff 1242 integer,intent(in) :: optkg 1243 1244 ! ********************************************************************* 1245 1246 #if defined HAVE_MPI_IO 1247 if (wff%iomode == IO_MODE_MPI) wff%kgwff=optkg 1248 #else 1249 ABI_UNUSED((/wff%iomode,optkg/)) 1250 #endif 1251 1252 end subroutine WffKg
m_wffile/WffOpen [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffOpen
FUNCTION
This subroutine opens a Wf file. It might be accessed by different mechanisms (usual F90 IO routines, MPI I/O, or, in the future, NetCDF). The routine provides a file handler, wff (a data structure containing all needed information).
INPUTS
iomode=access mode (0 means all procs access using usual F90 routines ; -1 means only the master proc access, using usual F90 routines ; 1 means MPI I/O; 2 means netcdf I/O) filename=name of the file master=the number of the master proc (only needed in parallel) me=my number (only needed in parallel) spaceComm= the space communicator handler (only needed in MPI parallel I/O) spaceWorld= the space communicator for the whole set of procs unwff=the file unit number
OUTPUT
ier=error code wff= structured info about the wavefunction file
SOURCE
990 subroutine WffOpen(iomode,spaceComm,filename,ier,wff,master,me,unwff,& 991 & spaceComm_mpiio) ! optional argument 992 993 !Arguments ------------------------------------ 994 integer, intent(in) :: iomode,spaceComm,master,me,unwff 995 integer, intent(in),optional :: spaceComm_mpiio 996 integer, intent(out) :: ier 997 character(len=fnlen), intent(in) :: filename 998 type(wffile_type), intent(inout) :: wff !vz_i 999 1000 !Local variables------------------------------- 1001 character(len=500) :: message 1002 character(len=fnlen) :: fildata 1003 #ifdef HAVE_MPI_IO 1004 integer :: isize 1005 #endif 1006 1007 ! ************************************************************************* 1008 1009 !Initialize the mandatory data of the wff datastructure 1010 wff%unwff =unwff 1011 wff%iomode =iomode; if (endswith(filename, ".nc")) wff%iomode = IO_MODE_ETSF 1012 if (filename/=wff%fname) wff%fname=filename 1013 1014 !Initialize info useful for parallel use 1015 wff%nproc =1 1016 wff%master =master 1017 wff%me =me 1018 wff%me_mpiio =0 1019 wff%spaceComm=spaceComm 1020 wff%spaceComm_mpiio=xmpi_comm_self 1021 1022 #if defined HAVE_MPI 1023 ! This case occurs when wff is connected to a DENSITY file 1024 ! abinit_comm_output is generally equal to MPI_COMM_WORLD (except if paral. over images) 1025 if (spaceComm==MPI_COMM_SELF) wff%spaceComm=abinit_comm_output 1026 ! if (spaceComm==MPI_COMM_SELF) wff%spaceComm=MPI_COMM_WORLD 1027 call MPI_COMM_SIZE(wff%spaceComm,wff%nproc,ier) 1028 ! Redefine the default MPIIO communicator if MPI, although MPIIO features should not be used unless 1029 ! present(spaceComm_mpiio).and.wff%iomode==1 1030 wff%spaceComm_mpiio=wff%spaceComm 1031 wff%me_mpiio=wff%me 1032 #endif 1033 1034 if (present(spaceComm_mpiio).and.any(wff%iomode==[IO_MODE_MPI, IO_MODE_ETSF])) wff%spaceComm_mpiio=spaceComm_mpiio 1035 #if defined HAVE_MPI 1036 call MPI_COMM_RANK(wff%spaceComm_mpiio,wff%me_mpiio,ier) 1037 #endif 1038 1039 ier=0 1040 if (wff%iomode==IO_MODE_FORTRAN) then ! All processors see a local file 1041 if (open_file(filename, message, unit=unwff, form="unformatted") /= 0) then 1042 ABI_ERROR(message) 1043 end if 1044 rewind(unwff) 1045 1046 else if (wff%iomode==IO_MODE_FORTRAN_MASTER) then ! Only the master processor see a local file 1047 if(master==me)then 1048 if (open_file(filename, message, unit=unwff, form="unformatted") /= 0) then 1049 ABI_ERROR(message) 1050 end if 1051 rewind(unwff) 1052 end if 1053 1054 #if defined HAVE_MPI_IO 1055 else if (wff%iomode==IO_MODE_MPI)then ! In the parallel case, only the master open filename file 1056 if(master==me)then 1057 if (open_file(filename, message, unit=unwff, form="unformatted") /= 0) then 1058 ABI_ERROR(message) 1059 end if 1060 rewind(unwff) 1061 end if 1062 ! MG: Great! These barriers lead to a deadlock if prtded hence MPI_FILE_OPEN is not called by all the processors! 1063 !call xmpi_barrier(wff%spaceComm) 1064 !call xmpi_barrier(wff%spaceComm_mpiio) 1065 1066 call MPI_FILE_OPEN(wff%spaceComm,filename,MPI_MODE_CREATE + MPI_MODE_RDWR,MPI_INFO_NULL,wff%fhwff,ier) 1067 ABI_CHECK_MPI(ier,sjoin("WffOpen:", filename)) 1068 1069 ! Define all type values 1070 call MPI_Type_size(MPI_INTEGER,isize,ier) 1071 wff%nbOct_int=isize 1072 call MPI_Type_size(MPI_DOUBLE_PRECISION,isize,ier) 1073 wff%nbOct_dp=isize 1074 call MPI_Type_size(MPI_CHARACTER,isize,ier) 1075 wff%nbOct_ch=isize 1076 wff%nbOct_recMarker=-1;wff%kgwff=-1;wff%formwff=-1 1077 wff%offwff=0;wff%off_recs=0;wff%lght_recs=0 1078 wff%marker_mpi_type=MPI_INTEGER ! Default value 1079 1080 #ifdef DEV_DEBUG_THIS 1081 wff%nbOct_recMarker=xmpio_bsize_frm 1082 wff%marker_mpi_type=xmpio_mpi_type_frm 1083 #endif 1084 1085 if (MPI_OFFSET_KIND==4) then 1086 wff%offset_mpi_type=MPI_INTEGER4 1087 else if (MPI_OFFSET_KIND==8) then 1088 wff%offset_mpi_type=MPI_INTEGER8 1089 #if defined HAVE_FC_INT_QUAD && defined HAVE_MPI_INTEGER16 1090 else if (MPI_OFFSET_KIND==16) then 1091 wff%offset_mpi_type=MPI_INTEGER16 1092 #endif 1093 else if (MPI_OFFSET_KIND==2) then 1094 wff%offset_mpi_type=MPI_INTEGER2 1095 end if 1096 #endif 1097 1098 #ifdef HAVE_NETCDF 1099 else if (wff%iomode==IO_MODE_ETSF)then 1100 fildata = nctk_ncify(filename) 1101 NCF_CHECK(nctk_open_modify(wff%unwff, fildata, xmpi_comm_self)) 1102 wff%fname = fildata 1103 !write(message,'(3A,I0)')'WffOpen: opening ', trim(wff%fname)," on unit ", wff%unwff 1104 !call wrtout(std_out, message, 'COLL') 1105 #endif 1106 else 1107 write(message, '(7a,i0,3a)' )& 1108 & 'For the time being the input variable iomode is restricted ',ch10,& 1109 & 'to 0 (all cases), 1 (in case MPI is enabled),',ch10,& 1110 & 'or 3 (only sequential, and if the NetCDF and ETSF_IO libraries have been enabled).',ch10,& 1111 & 'Its value is iomode= ',wff%iomode,'.',ch10,& 1112 & 'Action: change iomode or use ABINIT in parallel or enable NetCDF and/or ETSF_IO.' 1113 ABI_ERROR(message) 1114 end if 1115 1116 end subroutine WffOpen
m_wffile/WffReadDataRec_dp1d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffReadDataRec_dp1d
FUNCTION
Subroutine to read data in one record of a wavefunction file Handles double precision 1D arrays
INPUTS
ndp=size of the double precision array to be read wff= structured info about the wavefunction file
OUTPUT
dparray=array of double precision numbers ierr=error code
SIDE EFFECTS
SOURCE
1330 subroutine WffReadDataRec_dp1d(dparray,ierr,ndp,wff) 1331 1332 !Arguments ------------------------------------ 1333 type(wffile_type),intent(inout) :: wff 1334 integer,intent(in) :: ndp 1335 integer,intent(out) :: ierr 1336 real(dp),intent(out) :: dparray(ndp) 1337 1338 !Local variables------------------------------- 1339 character(len=500) :: msg 1340 1341 ! ************************************************************************* 1342 1343 ierr=0 1344 if (wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) then 1345 read (wff%unwff,iostat=ierr) dparray(1:ndp) 1346 1347 else if(wff%iomode==IO_MODE_MPI)then 1348 #if defined HAVE_MPI_IO 1349 call xderiveRRecInit(wff,ierr) 1350 call xderiveRead(wff,dparray,ndp,MPI_COMM_SELF,ierr) 1351 call xderiveRRecEnd(wff,ierr) 1352 #endif 1353 else 1354 write(msg,'(a,i0)')"Wrong iomode: ",wff%iomode 1355 ABI_ERROR(msg) 1356 end if 1357 1358 end subroutine WffReadDataRec_dp1d
m_wffile/WffReadDataRec_dp2d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffReadDataRec_dp2d
FUNCTION
Subroutine to read data in one record of a wavefunction file Handles double precision 2D arrays
INPUTS
n1,n2=sizes of the double precision array to be read wff= structured info about the wavefunction file
OUTPUT
dparray=array of double precision numbers ierr=error code
SIDE EFFECTS
SOURCE
1384 subroutine WffReadDataRec_dp2d(dparray,ierr,n1,n2,wff) 1385 1386 !Arguments ------------------------------------ 1387 type(wffile_type),intent(inout) :: wff 1388 integer,intent(in) :: n1,n2 1389 integer,intent(out) :: ierr 1390 real(dp),intent(out) :: dparray(n1,n2) 1391 1392 !Local variables------------------------------- 1393 character(len=500) :: msg 1394 1395 ! ************************************************************************* 1396 1397 ierr=0 1398 if (wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) then 1399 read (wff%unwff,iostat=ierr) dparray(1:n1,1:n2) 1400 1401 else if(wff%iomode==IO_MODE_MPI)then 1402 #if defined HAVE_MPI_IO 1403 call xderiveRRecInit(wff,ierr) 1404 call xderiveRead(wff,dparray,n1,n2,MPI_COMM_SELF,ierr) 1405 call xderiveRRecEnd(wff,ierr) 1406 #endif 1407 else 1408 write(msg,'(a,i0)')"Wrong iomode: ",wff%iomode 1409 ABI_ERROR(msg) 1410 end if 1411 1412 end subroutine WffReadDataRec_dp2d
m_wffile/WffReadNpwRec [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffReadNpwRec
FUNCTION
This subroutine read the npw record of a wavefunction file
INPUTS
wff= structured info about the wavefunction file wff%access == -1 and wf%master == Wff%me: read binary data wff%iomode == 0: read binary data wff%iomode == 1: use MPI/IO routines (MPIO defined) wff%iomode == 2: read netcdf format (NETCDF defined) ikpt= the i-th kpoint. isppol= the given spin polarisation element.
OUTPUT
ierr=error code (iostat integer from read statement) nband_disk=number of bands npw=number of plane waves nspinor=number of spinorial components of the wavefunctions
SIDE EFFECTS
SOURCE
1447 subroutine WffReadNpwRec(ierr,ikpt,isppol,nband_disk,npw,nspinor,wff) 1448 1449 !Arguments ------------------------------------ 1450 type(wffile_type),intent(inout) :: wff 1451 integer,intent(in) :: ikpt, isppol 1452 integer,intent(out) :: ierr,nband_disk,npw,nspinor 1453 1454 !Local variables------------------------------- 1455 !character(len=500) :: msg 1456 #if defined HAVE_NETCDF 1457 integer :: vid 1458 #endif 1459 1460 ! ************************************************************************* 1461 1462 ierr=0 1463 1464 if (wff%iomode == IO_MODE_FORTRAN.or.(wff%iomode==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me) ) then 1465 read (wff%unwff,iostat=ierr) npw,nspinor,nband_disk 1466 1467 else if(wff%iomode==IO_MODE_MPI)then 1468 #if defined HAVE_MPI_IO 1469 call xderiveRRecInit(wff,ierr) 1470 call xderiveRead(wff,npw,ierr) 1471 call xderiveRead(wff,nspinor,ierr) 1472 call xderiveRead(wff,nband_disk,ierr) 1473 call xderiveRRecEnd(wff,ierr) 1474 #endif 1475 1476 else if (wff%iomode == IO_MODE_ETSF) then 1477 1478 #if defined HAVE_NETCDF 1479 !write(std_out,*)"readnpwrec: ikpt, spin", ikpt, spin 1480 NCF_CHECK(nctk_get_dim(wff%unwff, "number_of_spinor_components", nspinor)) 1481 vid = nctk_idname(wff%unwff, "number_of_coefficients") 1482 NCF_CHECK(nf90_get_var(wff%unwff, vid, npw, start=[ikpt])) 1483 vid = nctk_idname(wff%unwff, "number_of_states") 1484 NCF_CHECK(nf90_get_var(wff%unwff, vid, nband_disk, start=[ikpt, isppol])) 1485 #endif 1486 1487 else 1488 ! MG: I don't understand why we have to use this ugly code!!!!!!!! 1489 ! Only master knows npw,nspinor,nband_disk in IO_MODE_FORTRAN_MASTE mode 1490 ! To the person who wrote this stuff: 1491 ! Have you ever heard about the "IF" statement of Fortran and the typical construct 1492 ! 1493 ! if (rank==master) call mpifoo_seq() 1494 1495 ABI_WARNING("Skipping read in WffReadNpwRec. Keep fingers crossed") 1496 ! MG: Must initialze these values somehow to avoid overflows. 1497 npw = 0; nspinor = 0; nband_disk = 0 1498 end if 1499 1500 !write(std_out,*)"nband_disk,npw,nspinor",nband_disk,npw,nspinor 1501 ABI_CHECK(ierr==0,"ierr!=0") 1502 1503 end subroutine WffReadNpwRec
m_wffile/WffReadSkipRec [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffReadSkipRec
FUNCTION
This subroutine move forward or backward in a Wf file by nrec records.
INPUTS
nrec=number of records wff= structured info about the wavefunction file
OUTPUT
ierr=error code
TODO
For the future : one should treat the possible errors of backspace
SOURCE
1528 subroutine WffReadSkipRec(ierr,nrec,wff) 1529 1530 !Arguments ------------------------------------ 1531 integer,intent(in) :: nrec 1532 integer,intent(out) :: ierr 1533 type(wffile_type),intent(inout) :: wff 1534 1535 !Local variables------------------------------- 1536 #if defined HAVE_MPI_IO 1537 integer :: irec 1538 integer(kind=MPI_OFFSET_KIND) :: delim_record,offset 1539 #endif 1540 1541 ! ************************************************************************* 1542 1543 ierr=0 1544 if( wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) then 1545 1546 call mvrecord(wff%unwff,nrec,ierr) 1547 ABI_CHECK(ierr==0,"error in mvrecord") 1548 1549 1550 else if(wff%iomode==IO_MODE_MPI)then 1551 #if defined HAVE_MPI_IO 1552 if (nrec>0) then ! Move forward nrec records 1553 do irec=1,nrec 1554 wff%off_recs = wff%offwff 1555 call rwRecordMarker(1,wff%offwff,delim_record,wff,ierr) 1556 wff%lght_recs = delim_record 1557 end do 1558 else ! Move backward -nrec records 1559 do irec=1,-nrec 1560 offset = wff%offwff-wff%nbOct_recMarker 1561 call rwRecordMarker(1,offset,delim_record,wff,ierr) 1562 wff%lght_recs = delim_record 1563 wff%offwff = wff%offwff - delim_record - 2*wff%nbOct_recMarker 1564 wff%off_recs = wff%offwff 1565 end do 1566 end if 1567 #endif 1568 end if ! wff%iomode==0,1 or -1 1569 1570 end subroutine WffReadSkipRec
m_wffile/WffReadWrite_mpio [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffReadWrite_mpio
FUNCTION
This procedure read or write cg in the file _WFK using MPI_IO when cg are dispatched amoung commcart communicator
INPUTS
wff=struct info for wavefunction nband_disk=number of bands on disk files to be write icg=shift to be given to the location of the cg array mcg=second dimention of cg mpi_enreg=information about parallelisation depl_mpi_to_seq=for each proc, index of cg in sequential mode npwso=npw*nspinor number of plane waves treated by this node. npwsotot=npwtot*nspinor Total number of planewaves Used to calculate the size of data to be written. rdwr=1 if reading, 2 if writing
OUTPUT
ierr=error status
SIDE EFFECTS
cg(2,mcg)=planewave coefficients of wavefunctions,
NOTES
cg is written like the following: BeginMarker cg ( iband = 1 ) EndMarker BeginMarker cg ( iband = 2 ) EndMarker ... BeginMarker cg( iband = nband_disk ) EndMarker BeginMarker and EndMarker give the value of the total length of cg for one band For MPI-IO library the performance is improved by the use a "view" of the file for each proc.
SOURCE
1613 subroutine WffReadWrite_mpio(wff,rdwr,cg,mcg,icg,nband_disk,npwso,npwsotot,depl_mpi_to_seq,ierr) 1614 1615 !Arguments ------------------------------------ 1616 !scalars 1617 integer,intent(in) :: icg,mcg,nband_disk,npwso,npwsotot,rdwr 1618 integer,intent(out) :: ierr 1619 type(wffile_type),intent(inout) :: wff 1620 !arrays 1621 integer,intent(in) :: depl_mpi_to_seq(npwso) 1622 real(dp),intent(inout) :: cg(2,mcg) 1623 1624 !Local variables------------------------------- 1625 !scalars 1626 #if defined HAVE_MPI_IO 1627 integer,parameter :: MAXBAND=500, check_markers=1 1628 integer :: filetype,iband,ibandmax,ibandmin,iblock,ii,iloc,ipw,jerr,jj 1629 integer :: nb,nband_block,step,totsize1bandcg,wfftempo 1630 integer(kind=MPI_OFFSET_KIND) :: delim_record,loc_depl_band,offset,totsize1bandByte 1631 character(len=500) :: msg 1632 !arrays 1633 integer,allocatable :: BlockLength(:),BlockType(:),map(:),tempo_map(:) 1634 integer(kind=MPI_OFFSET_KIND),allocatable :: BlockDepl(:) 1635 integer(kind=2),allocatable :: bufdelim2(:) 1636 integer(kind=4),allocatable :: bufdelim4(:) 1637 integer(kind=8),allocatable :: bufdelim8(:) 1638 real(dp),allocatable :: buf(:),tempo_buf(:) 1639 #if defined HAVE_FC_INT_QUAD 1640 integer(kind=16),allocatable :: bufdelim16(:) 1641 #endif 1642 #endif 1643 1644 ! ********************************************************************* 1645 1646 ierr=0 1647 1648 #if defined HAVE_MPI_IO 1649 !---------------------------------------------- 1650 !! Prepare WF data 1651 !---------------------------------------------- 1652 !Init offset of record 1653 wff%off_recs = wff%offwff 1654 1655 !Total size to be written (in number of bands and in bytes) 1656 totsize1bandcg=2*npwsotot 1657 !call xmpi_sum(totsize1bandcg,wff%spaceComm_mpiio,ierr) 1658 1659 totsize1bandByte=totsize1bandcg*wff%nbOct_dp+2*wff%nbOct_recMarker 1660 1661 !Check file size 1662 offset=wff%offwff+nband_disk*totsize1bandByte 1663 if (offset>Huge(offset)) then 1664 msg='File is too large for MPI-IO specifications !' 1665 ABI_ERROR(msg) 1666 end if 1667 1668 !Open file 1669 call MPI_FILE_OPEN(wff%spaceComm_mpiio,wff%fname,MPI_MODE_RDWR,MPI_INFO_NULL,wfftempo,ierr) 1670 ABI_CHECK_MPI(ierr, sjoin("MPI_FILE_OPEN:", wff%fname)) 1671 1672 !---------------------------------------------------------- 1673 !Loop blocks of bands (to decrease offsets inside the file) 1674 !---------------------------------------------------------- 1675 ibandmax=0;ibandmin=1 1676 ii=huge(check_markers)/totsize1bandByte;step=min(ii,MAXBAND,nband_disk) 1677 do iblock=1,nband_disk/step+1 1678 ibandmax=min(ibandmin+step-1,nband_disk) 1679 nband_block=ibandmax-ibandmin+1 1680 offset=wff%offwff+(ibandmin-1)*totsize1bandByte 1681 1682 ! ---------------------------------------------- 1683 ! Read/Write bands 1684 ! ---------------------------------------------- 1685 1686 ! Build map; for better performance, map must be in increasing order 1687 ABI_MALLOC_OR_DIE(map,(2*npwso*nband_block), ierr) 1688 1689 ABI_STAT_MALLOC(buf,(2*npwso*nband_block), ierr) 1690 ABI_CHECK(ierr==0, "out of memory in wavefunction buffer. Try to decrease MAXBAND in WffReadWrite_mpio") 1691 1692 if (rdwr==1) then 1693 ! If reading, only build map 1694 nb=0;loc_depl_band=0 1695 ABI_MALLOC(tempo_map,(2*npwso)) 1696 do iband=ibandmin,ibandmax 1697 tempo_map(1:2*npwso)=-1 1698 jj=1;ipw=(iband-1)*npwso+icg 1699 do ii=1,npwso 1700 iloc=loc_depl_band+wff%nbOct_recMarker+2*(depl_mpi_to_seq(ii)-1)*wff%nbOct_dp 1701 tempo_map(jj )=iloc ! Real part 1702 tempo_map(jj+1)=iloc+wff%nbOct_dp ! Imag part 1703 jj=jj+2 1704 end do 1705 do ii=1,2*npwso ! Now, elimate holes 1706 if (tempo_map(ii)/=-1) then 1707 nb=nb+1 1708 map(nb)=tempo_map(ii) 1709 end if 1710 end do 1711 loc_depl_band=loc_depl_band+totsize1bandByte ! Location in bytes 1712 end do 1713 else if (rdwr==2) then 1714 ! If writing, build map and store cg in a buffer 1715 nb=0;loc_depl_band=0 1716 ABI_MALLOC(tempo_map,(2*npwso)) 1717 ABI_MALLOC(tempo_buf,(2*npwso)) 1718 do iband=ibandmin,ibandmax 1719 tempo_map(1:2*npwso)=-1 1720 jj=1;ipw=(iband-1)*npwso+icg 1721 do ii=1,npwso 1722 iloc=loc_depl_band+wff%nbOct_recMarker+2*(depl_mpi_to_seq(ii)-1)*wff%nbOct_dp 1723 tempo_map(jj )=iloc ! Real part 1724 tempo_map(jj+1)=iloc+wff%nbOct_dp ! Imag part 1725 tempo_buf(jj:jj+1)=cg(1:2,ipw+ii) 1726 jj=jj+2 1727 end do 1728 do ii=1,2*npwso ! Now, elimate holes 1729 if (tempo_map(ii)/=-1) then 1730 nb=nb+1 1731 map(nb)=tempo_map(ii) 1732 buf(nb)=tempo_buf(ii) 1733 end if 1734 end do 1735 loc_depl_band=loc_depl_band+totsize1bandByte ! Location in bytes 1736 end do 1737 ABI_FREE(tempo_map) 1738 ABI_FREE(tempo_buf) 1739 end if ! rdwr 1740 1741 ! Build and commit MPI datatype 1742 ABI_MALLOC(BlockLength,(nb+2)) 1743 ABI_MALLOC(BlockDepl,(nb+2)) 1744 ABI_MALLOC(BlockType,(nb+2)) 1745 BlockLength(1)=1;BlockDepl(1)=0;BlockType(1)=MPI_LB 1746 do ii=2,nb+1 1747 BlockLength(ii)=1 1748 BlockDepl(ii)=map(ii-1) 1749 BlockType(ii)=MPI_DOUBLE_PRECISION 1750 end do 1751 BlockLength(nb+2)=1;BlockDepl(nb+2)=totsize1bandByte*nband_block;BlockType(nb+2)=MPI_UB 1752 call xmpio_type_struct(nb+2,BlockLength,BlockDepl,BlockType,filetype,ierr) 1753 call MPI_TYPE_COMMIT(filetype,ierr) 1754 ABI_FREE(BlockLength) 1755 ABI_FREE(BlockDepl) 1756 ABI_FREE(BlockType) 1757 1758 ! Read/Write data on disk 1759 call MPI_FILE_SET_VIEW(wfftempo,offset,MPI_BYTE,filetype,"native",MPI_INFO_NULL,ierr) 1760 if (rdwr==1) then 1761 call MPI_FILE_READ_ALL (wfftempo,buf,nb,MPI_DOUBLE_PRECISION,MPI_STATUS_IGNORE,ierr) 1762 else 1763 call MPI_FILE_WRITE_ALL(wfftempo,buf,nb,MPI_DOUBLE_PRECISION,MPI_STATUS_IGNORE,ierr) 1764 end if 1765 1766 ! In case of reading, retrieve cg 1767 if (rdwr==1) then 1768 nb=0;loc_depl_band=0 1769 ABI_MALLOC(tempo_buf,(2*npwso)) 1770 do iband=ibandmin,ibandmax 1771 do ii=1,2*npwso ! Now, elimate holes 1772 if (tempo_map(ii)/=-1) then 1773 nb=nb+1;tempo_buf(ii)=buf(nb) 1774 end if 1775 end do 1776 jj=1;ipw=(iband-1)*npwso+icg 1777 do ii=1,npwso 1778 iloc=loc_depl_band+wff%nbOct_recMarker+2*(depl_mpi_to_seq(ii)-1)*wff%nbOct_dp 1779 cg(1:2,ipw+ii)=tempo_buf(jj:jj+1) 1780 jj=jj+2 1781 end do 1782 loc_depl_band=loc_depl_band+totsize1bandByte ! Location in bytes 1783 end do 1784 ABI_FREE(tempo_map) 1785 ABI_FREE(tempo_buf) 1786 end if ! rdwr 1787 1788 ! Free memory 1789 ABI_FREE(map) 1790 ABI_FREE(buf) 1791 call MPI_TYPE_FREE(filetype,ierr) 1792 1793 ! ---------------------------------------------- 1794 ! Check/Write record markers (only master proc) 1795 ! ---------------------------------------------- 1796 if ((rdwr==1.and.check_markers==1).or.(rdwr==2)) then 1797 1798 ! Define view for the file 1799 nb=2*nband_block 1800 ABI_MALLOC(BlockLength,(nb+2)) 1801 ABI_MALLOC(BlockDepl,(nb+2)) 1802 ABI_MALLOC(BlockType,(nb+2)) 1803 BlockLength(1)=1;BlockDepl(1)=0;BlockType(1)=MPI_LB 1804 jj=2 1805 do ii=1,nband_block 1806 BlockType(jj:jj+1) =wff%marker_mpi_type 1807 BlockLength(jj:jj+1)=1 1808 BlockDepl(jj )=(ii-1)*totsize1bandByte 1809 BlockDepl(jj+1)= ii *totsize1bandByte-wff%nbOct_recMarker 1810 jj=jj+2 1811 end do 1812 BlockLength(nb+2)=1;BlockDepl(nb+2)=nband_block*totsize1bandByte;BlockType(nb+2)=MPI_UB 1813 call xmpio_type_struct(nb+2,BlockLength,BlockDepl,BlockType,filetype,ierr) 1814 call MPI_TYPE_COMMIT(filetype,ierr) 1815 call MPI_FILE_SET_VIEW(wfftempo,offset,MPI_BYTE,filetype,"native",MPI_INFO_NULL,ierr) 1816 ABI_FREE(BlockLength) 1817 ABI_FREE(BlockDepl) 1818 ABI_FREE(BlockType) 1819 1820 ! Read/Write all markers (depend on Fortran marker MPI type) 1821 if (wff%me_mpiio==0) then 1822 jerr=0;delim_record=totsize1bandByte-2*wff%nbOct_recMarker 1823 if (wff%nbOct_recMarker==4) then 1824 ABI_MALLOC(bufdelim4,(nb)) 1825 if (rdwr==2) bufdelim4(:)=delim_record 1826 if (rdwr==1) then 1827 call MPI_FILE_READ (wfftempo,bufdelim4,2*nband_block,wff%marker_mpi_type,MPI_STATUS_IGNORE,ierr) 1828 if (any(bufdelim4(:)/=delim_record)) jerr=1 1829 else 1830 call MPI_FILE_WRITE(wfftempo,bufdelim4,2*nband_block,wff%marker_mpi_type,MPI_STATUS_IGNORE,ierr) 1831 end if 1832 ABI_FREE(bufdelim4) 1833 else if (wff%nbOct_recMarker==8) then 1834 ABI_MALLOC(bufdelim8,(nb)) 1835 if (rdwr==2) bufdelim8(:)=delim_record 1836 if (rdwr==1) then 1837 call MPI_FILE_READ (wfftempo,bufdelim8,2*nband_block,wff%marker_mpi_type,MPI_STATUS_IGNORE,ierr) 1838 if (any(bufdelim8(:)/=delim_record)) jerr=1 1839 else 1840 call MPI_FILE_WRITE(wfftempo,bufdelim8,2*nband_block,wff%marker_mpi_type,MPI_STATUS_IGNORE,ierr) 1841 end if 1842 ABI_FREE(bufdelim8) 1843 #if defined HAVE_FC_INT_QUAD 1844 else if (wff%nbOct_recMarker==16) then 1845 ABI_MALLOC(bufdelim16,(nb)) 1846 if (rdwr==2) bufdelim16(:)=delim_record 1847 if (rdwr==1) then 1848 call MPI_FILE_READ (wfftempo,bufdelim16,2*nband_block,wff%marker_mpi_type,MPI_STATUS_IGNORE,ierr) 1849 if (any(bufdelim16(:)/=delim_record)) jerr=1 1850 else 1851 call MPI_FILE_WRITE(wfftempo,bufdelim16,2*nband_block,wff%marker_mpi_type,MPI_STATUS_IGNORE,ierr) 1852 end if 1853 ABI_FREE(bufdelim16) 1854 #endif 1855 else if (wff%nbOct_recMarker==2) then 1856 ABI_MALLOC(bufdelim2,(nb)) 1857 if (rdwr==2) bufdelim2(:)=delim_record 1858 if (rdwr==1) then 1859 call MPI_FILE_READ (wfftempo,bufdelim2,2*nband_block,wff%marker_mpi_type,MPI_STATUS_IGNORE,ierr) 1860 if (any(bufdelim2(:)/=delim_record)) jerr=1 1861 else 1862 call MPI_FILE_WRITE(wfftempo,bufdelim2,2*nband_block,wff%marker_mpi_type,MPI_STATUS_IGNORE,ierr) 1863 end if 1864 ABI_FREE(bufdelim2) 1865 end if 1866 if (rdwr==1.and.jerr==1) then 1867 write(unit=msg,fmt='(2a)') 'Error when reading record markers of file ',trim(wff%fname) 1868 ABI_ERROR(msg) 1869 end if 1870 end if ! me_mpiio=0 1871 1872 ! Free memory 1873 call MPI_TYPE_FREE(filetype,ierr) 1874 1875 end if ! rdwr 1876 1877 ! ----------------------------------------- 1878 ! End loop on blocks of bands 1879 ! ----------------------------------------- 1880 if (ibandmax>=nband_disk) exit 1881 ibandmin=ibandmax+1 1882 end do 1883 1884 !----------------------------------------- 1885 !End statements 1886 !----------------------------------------- 1887 !Close file 1888 call MPI_FILE_CLOSE(wfftempo,ierr) 1889 1890 !Update offset 1891 wff%offwff=wff%offwff+totsize1bandByte*nband_disk 1892 #endif 1893 1894 #if !defined HAVE_MPI_IO 1895 !Dummy check to avoid warning from compilers. 1896 ABI_UNUSED((/wff%iomode,rdwr,size(cg),mcg,icg,nband_disk,npwso,depl_mpi_to_seq(1),npwsotot/)) 1897 #endif 1898 1899 end subroutine WffReadWrite_mpio
m_wffile/WffWriteDataRec_dp1d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffWriteDataRec_dp1d
FUNCTION
Subroutine to write data in one record of a wavefunction file Handles double precision 1D arrays
INPUTS
dparray=array of double precision numbers ndp=size of the double precision array to be written wff= structured info about the wavefunction file
OUTPUT
ierr=error code
SIDE EFFECTS
SOURCE
1978 subroutine WffWriteDataRec_dp1d(dparray,ierr,ndp,wff) 1979 1980 !Arguments ------------------------------------ 1981 type(wffile_type),intent(inout) :: wff 1982 integer,intent(in) :: ndp 1983 integer,intent(out) :: ierr 1984 real(dp),intent(in) :: dparray(ndp) 1985 1986 !Local variables------------------------------- 1987 character(len=500) :: msg 1988 1989 ! ************************************************************************* 1990 1991 ierr=0 1992 if (wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) then 1993 write(wff%unwff,iostat=ierr) dparray(1:ndp) 1994 1995 else if(wff%iomode==IO_MODE_MPI)then 1996 #if defined HAVE_MPI_IO 1997 call xderiveWRecInit(wff,ierr) 1998 call xderiveWrite(wff,dparray,ndp,MPI_COMM_SELF,ierr) 1999 call xderiveWRecEnd(wff,ierr) 2000 #endif 2001 else 2002 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 2003 ABI_WARNING(msg) 2004 end if 2005 2006 end subroutine WffWriteDataRec_dp1d
m_wffile/WffWriteDataRec_dp2d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffWriteDataRec_dp2d
FUNCTION
Subroutine to write data in one record of a wavefunction file Handles double precision 2D arrays
INPUTS
dparray=array of double precision numbers n1,n2=sizes of the double precision array to be written wff= structured info about the wavefunction file
OUTPUT
ierr=error code
SIDE EFFECTS
SOURCE
2033 subroutine WffWriteDataRec_dp2d(dparray,ierr,n1,n2,wff) 2034 2035 !Arguments ------------------------------------ 2036 type(wffile_type),intent(inout) :: wff 2037 integer,intent(in) :: n1,n2 2038 integer,intent(out) :: ierr 2039 real(dp),intent(in) :: dparray(n1,n2) 2040 2041 !Local variables------------------------------- 2042 character(len=500) :: msg 2043 2044 ! ************************************************************************* 2045 2046 ierr=0 2047 if (wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) then 2048 write(wff%unwff,iostat=ierr) dparray(1:n1,1:n2) 2049 2050 else if(wff%iomode==IO_MODE_MPI)then 2051 #if defined HAVE_MPI_IO 2052 call xderiveWRecInit(wff,ierr) 2053 call xderiveWrite(wff,dparray,n1,n2,MPI_COMM_SELF,ierr) 2054 call xderiveWRecEnd(wff,ierr) 2055 #endif 2056 else 2057 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 2058 ABI_WARNING(msg) 2059 end if 2060 2061 end subroutine WffWriteDataRec_dp2d
m_wffile/WffWriteDataRec_int2d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffWriteDataRec_int2d
FUNCTION
Subroutine to write data in one record of a wavefunction file Handles integer 2D arrays
INPUTS
intarray=array of integer numbers n1,n2=sizes of the integer array to be written wff= structured info about the wavefunction file
OUTPUT
ierr=error code
SIDE EFFECTS
SOURCE
1924 subroutine WffWriteDataRec_int2d(intarray,ierr,n1,n2,wff) 1925 1926 !Arguments ------------------------------------ 1927 type(wffile_type),intent(inout) :: wff 1928 integer,intent(in) :: n1,n2 1929 integer,intent(out) :: ierr 1930 integer,intent(in) :: intarray(n1,n2) 1931 1932 !Local variables------------------------------- 1933 character(len=500) :: msg 1934 1935 ! ************************************************************************* 1936 1937 ierr=0 1938 if (wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) then 1939 write(wff%unwff,iostat=ierr) intarray(1:n1,1:n2) 1940 1941 else if(wff%iomode==IO_MODE_MPI)then 1942 #if defined HAVE_MPI_IO 1943 call xderiveWRecInit(wff,ierr) 1944 call xderiveWrite(wff,intarray,n1,n2,MPI_COMM_SELF,ierr) 1945 call xderiveWRecEnd(wff,ierr) 1946 #endif 1947 else 1948 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 1949 ABI_WARNING(msg) 1950 end if 1951 1952 end subroutine WffWriteDataRec_int2d
m_wffile/WffWriteNpwRec [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
WffWriteNpwRec
FUNCTION
This subroutine writes the npw record of a wavefunction file
INPUTS
wff= structured info about the wavefunction file nband_disk=number of bands npw=number of plane waves nspinor=number of spinorial components of the wavefunctions opt_paral=(optional argument, default=1, only used for MPI-IO) 1: all procs in the communicator write the data 2: only master in the communicator writes the data
OUTPUT
ierr=error code
SIDE EFFECTS
SOURCE
2091 subroutine WffWriteNpwRec(ierr,nband_disk,npw,nspinor,wff,& 2092 & opt_paral) ! optional argument 2093 2094 !Arguments ------------------------------------ 2095 type(wffile_type),intent(inout) :: wff 2096 integer,intent(in) :: nband_disk,npw,nspinor 2097 integer,intent(in),optional :: opt_paral 2098 integer,intent(out) :: ierr 2099 2100 !Local variables------------------------------- 2101 integer :: opt_paral_ 2102 character(len=500) :: msg 2103 #if defined HAVE_MPI_IO 2104 integer :: me 2105 integer(kind=MPI_OFFSET_KIND) :: off(1) 2106 #endif 2107 2108 ! ************************************************************************* 2109 2110 ierr=0 2111 opt_paral_=1;if (present(opt_paral)) opt_paral_=opt_paral 2112 2113 if (wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode ==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) then 2114 write(wff%unwff,iostat=ierr) npw,nspinor,nband_disk 2115 2116 else if(wff%iomode==IO_MODE_MPI)then 2117 #if defined HAVE_MPI_IO 2118 me=-1;if (opt_paral_==2) me=wff%me_mpiio 2119 if ((me==-1.and.opt_paral_==1).or.(me==0.and.opt_paral_==2)) then 2120 call xderiveWRecInit(wff,ierr) 2121 call xderiveWrite(wff,npw,ierr) 2122 call xderiveWrite(wff,nspinor,ierr) 2123 call xderiveWrite(wff,nband_disk,ierr) 2124 call xderiveWRecEnd(wff,ierr) 2125 end if 2126 if (opt_paral_==2.and.wff%spaceComm_mpiio/=MPI_COMM_SELF) then 2127 call xmpi_barrier(wff%spaceComm_mpiio) 2128 off(1)=wff%offwff 2129 call MPI_BCAST(off,1,wff%offset_mpi_type,0,wff%spaceComm_mpiio,ierr) 2130 wff%offwff=off(1) 2131 end if 2132 #endif 2133 else 2134 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 2135 ABI_WARNING(msg) 2136 end if 2137 2138 end subroutine WffWriteNpwRec
m_wffile/xdefineOff [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xdefineOff
FUNCTION
In case of MPI I/O, defines the offset for each processor
INPUTS
formeig option (format of the eigenvalues and occupations) : 0 => ground-state format (initialisation of eigenvectors with random numbers, vector of eigenvalues, occupations are present) 1 => respfn format (initialisation of eigenvectors with 0 s, hermitian matrix of eigenvalues) nkpt = number of k points nspinor = total number of spinorial components of the wavefunctions nsppol = number of channels for spin-polarization (1 or 2) nband(nkpt*nsppol) = number of bands at each k point, for each polarization npwarr(nkpt) = number of planewaves at each k point mpi_enreg <type(MPI_type)> = information about MPI parallelization
OUTPUT
(no output)
SIDE EFFECTS
wff <type(wffile_type)> =
SOURCE
3686 subroutine xdefineOff(formeig,wff,mpi_enreg,nband,npwarr,nspinor,nsppol,nkpt) 3687 3688 !Arguments ------------------------------------ 3689 integer, intent(in) :: nsppol,nkpt,nspinor,formeig 3690 integer, intent(in) :: nband(nkpt*nsppol),npwarr(nkpt) 3691 type(wffile_type),intent(inout) :: wff 3692 type(MPI_type),intent(in) :: mpi_enreg 3693 3694 !Local variables------------------------------- 3695 #if defined HAVE_MPI_IO 3696 !scalars 3697 integer :: comm,iproc 3698 integer :: nband_k,npw_k,nproc,me,ipp 3699 integer :: nbrec,isppol,ikpt,nbint,nbreal,nbd,ippband 3700 integer :: nrecnpw,nreckg 3701 integer(kind=XMPI_OFFSET_KIND) :: pos_start 3702 !arrays 3703 integer(kind=XMPI_OFFSET_KIND),allocatable :: offproc(:) 3704 #endif 3705 3706 ! ************************************************************************* 3707 !nbOct_int octet number of int value 3708 !nbOct_dp octet number of dp value 3709 !nbOct_ch octet number of character value 3710 !lght_recs length of record 3711 3712 if(.false.)write(std_out,*)wff%me,mpi_enreg%nproc,formeig,nband,npwarr,nspinor,nkpt 3713 #if defined HAVE_MPI_IO 3714 if(wff%iomode==IO_MODE_MPI)then 3715 3716 comm=mpi_enreg%comm_cell 3717 me=xmpi_comm_rank(comm) 3718 nproc=xmpi_comm_size(comm) 3719 pos_start=wff%offwff 3720 3721 ABI_MALLOC(offproc,(0:nproc)) 3722 offproc = 0 3723 nbrec =2 3724 nrecnpw=3+nbrec 3725 3726 do isppol=1,nsppol 3727 do ikpt=1,nkpt 3728 nband_k=nband(ikpt+(isppol-1)*nkpt) 3729 npw_k=npwarr(ikpt) 3730 iproc=mpi_enreg%proc_distrb(ikpt,1,isppol) 3731 if (mpi_enreg%paralbd==1) iproc=mpi_enreg%proc_distrb(ikpt,1,isppol) 3732 ! record kg 3733 nreckg=nbrec+ wff%kgwff*3*npw_k 3734 3735 ! Record npw,nspinor,nband, Record kg 3736 offproc(iproc) = offproc(iproc) + wff%nbOct_int*(nrecnpw+nreckg) 3737 3738 if (formeig == 0) then 3739 ! Records eigen,occ 3740 nbint=nbrec 3741 nbreal = 2 *nband_k 3742 offproc(iproc) = offproc(iproc) + (wff%nbOct_int*nbint+wff%nbOct_dp*nbreal) 3743 3744 ! Records cg 3745 offproc(iproc) = offproc(iproc) & 3746 & + (wff%nbOct_int*nbrec+wff%nbOct_dp*2*npw_k*nspinor)*nband_k 3747 3748 ippband=iproc 3749 do nbd=1,nband_k 3750 ipp=mpi_enreg%proc_distrb(ikpt,nbd,isppol) 3751 if (ipp /= ippband ) then 3752 ippband=ipp 3753 offproc(ippband)=offproc(ippband)+ wff%nbOct_int*(nrecnpw+nreckg) 3754 offproc(ippband) = offproc(ippband) + (wff%nbOct_int*nbint & 3755 & +wff%nbOct_dp*nbreal) 3756 offproc(ippband) = offproc(ippband) + (wff%nbOct_int*nbrec & 3757 & + wff%nbOct_dp*2*npw_k*nspinor)*nband_k 3758 end if 3759 end do 3760 else if (formeig == 1) then 3761 ! record eigen 3762 offproc(iproc) = offproc(iproc) + (wff%nbOct_int*2*nbrec & 3763 & + wff%nbOct_dp*2*npw_k*nspinor & 3764 & + wff%nbOct_dp*2*nband_k)*nband_k 3765 ippband=iproc 3766 do nbd=1,nband_k 3767 ipp=mpi_enreg%proc_distrb(ikpt,nbd,isppol) 3768 if (ipp /= ippband) then 3769 ippband=ipp 3770 offproc(ippband)=offproc(ippband)+ wff%nbOct_int*(nrecnpw+nreckg) 3771 offproc(ippband) = offproc(ippband) + (wff%nbOct_int*2*nbrec & 3772 & + wff%nbOct_dp*2*npw_k*nspinor & 3773 & + wff%nbOct_dp*2*nband_k)*nband_k 3774 end if 3775 end do 3776 end if ! formeig 3777 end do ! ikpt 3778 3779 end do ! isppol 3780 3781 ! pos_start=wff%offwff 3782 ! wff%offwff = pos_start 3783 3784 if (me/=0)then 3785 do iproc=0,me-1 3786 wff%offwff=wff%offwff+offproc(iproc) 3787 end do 3788 end if 3789 ABI_FREE(offproc) 3790 3791 end if ! iomode 3792 #endif 3793 3794 end subroutine xdefineOff
m_wffile/xderiveRead_dp [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveRead_dp
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: double precision scalar.
INPUTS
(none)
OUTPUT
xval= data buffer ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2381 subroutine xderiveRead_dp(wff,xval,ierr) 2382 2383 !Arguments ------------------------------------ 2384 type(wffile_type),intent(inout) :: wff 2385 integer,intent(out) :: ierr 2386 real(dp),intent(out) :: xval 2387 2388 !Local variables------------------------------- 2389 #if defined HAVE_MPI_IO 2390 integer :: statux(MPI_STATUS_SIZE) 2391 real(dp) :: arr_xval(1) 2392 #endif 2393 character(len=500) :: msg 2394 2395 ! ********************************************************************* 2396 2397 xval=zero ; ierr=0 2398 if(.false.)write(std_out,*)wff%me 2399 #if defined HAVE_MPI_IO 2400 call MPI_FILE_READ_AT(wff%fhwff,wff%offwff,arr_xval,1,MPI_DOUBLE_PRECISION,statux,ierr) 2401 xval=arr_xval(1) 2402 wff%offwff = wff%offwff + wff%nbOct_dp 2403 return 2404 #endif 2405 2406 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 2407 ABI_WARNING(msg) 2408 2409 end subroutine xderiveRead_dp
m_wffile/xderiveRead_dp2d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveRead_dp2d
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: double precision two-dimensional arrays.
INPUTS
n1= first dimension of the array n2= second dimension of the array spaceComm= MPI communicator
OUTPUT
ierr= exit status, a non-zero value meaning there is an error xval= data buffer array
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2519 subroutine xderiveRead_dp2d(wff,xval,n1,n2,spaceComm,ierr) 2520 2521 !Arguments ------------------------------------ 2522 type(wffile_type),intent(inout) :: wff 2523 integer,intent(in) :: n1,n2,spaceComm 2524 integer,intent(out) :: ierr 2525 real(dp),intent(out) :: xval(:,:) 2526 2527 !Local variables------------------------------- 2528 #if defined HAVE_MPI_IO 2529 integer(kind=MPI_OFFSET_KIND) :: delim_record,nboct,dispoct(1),posit,totoct 2530 integer :: statux(MPI_STATUS_SIZE) 2531 #endif 2532 character(len=500) :: msg 2533 2534 ! ********************************************************************* 2535 2536 xval(:,:)=zero ; ierr=0 ! Initialization, for the compiler 2537 if(.false.)write(std_out,*)wff%me,n1,n2,spaceComm 2538 2539 #if defined HAVE_MPI_IO 2540 nboct = wff%nbOct_dp * n1 *n2 2541 posit = wff%offwff 2542 delim_record = posit - wff%off_recs + wff%lght_recs - wff%nbOct_recMarker 2543 2544 if (delim_record >= nboct) then 2545 ! Compute offset for local part 2546 ! dispoct = sum (nboct, rank=0..me) 2547 if (spaceComm/=MPI_COMM_SELF) then 2548 call MPI_SCAN([nboct],dispoct,1,wff%offset_mpi_type,MPI_SUM,spaceComm,ierr) 2549 posit = posit + dispoct(1) - nboct 2550 end if 2551 call MPI_FILE_READ_AT(wff%fhwff,posit,xval,n1*n2,MPI_DOUBLE_PRECISION,statux,ierr) 2552 2553 ! get the total number of bits wrote by processors 2554 if (spaceComm/=MPI_COMM_SELF) then 2555 call xmpi_max(dispoct(1),totoct,spaceComm,ierr) 2556 !call MPI_ALLREDUCE(dispoct,totoct,1,wff%offset_mpi_type,MPI_MAX,spaceComm,ierr) 2557 else 2558 totoct=nboct 2559 end if 2560 else 2561 ierr = 1 2562 nboct = 0 2563 totoct = 0 2564 end if 2565 2566 !new offset 2567 wff%offwff=wff%offwff + totoct 2568 return 2569 #endif 2570 2571 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 2572 ABI_WARNING(msg) 2573 2574 end subroutine xderiveRead_dp2d
m_wffile/xderiveRead_dp2d_displ [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveRead_dp2d_displ
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: double precision two-dimensional arrays.
INPUTS
n1= first dimension of the array n2= second dimension of the array spaceComm= MPI communicator displace= number of elements for the offset
OUTPUT
ierr= exit status, a non-zero value meaning there is an error xval= data buffer array
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2736 subroutine xderiveRead_dp2d_displ(wff,xval,n1,n2,spaceComm,displace,ierr) 2737 2738 !Arguments ------------------------------------ 2739 type(wffile_type),intent(inout) :: wff 2740 integer,intent(in) :: n1,n2,spaceComm 2741 integer,intent(out) :: ierr 2742 real(dp),intent(out):: xval(:,:) 2743 integer,intent(in):: displace(:) 2744 2745 !Local variables------------------------------- 2746 #if defined HAVE_MPI_IO 2747 !scalars 2748 integer :: filetype,i1,i2,ipos,nb,nbval,totsize,wfftempo 2749 !arrays 2750 integer :: statux(MPI_STATUS_SIZE) 2751 integer,allocatable :: length1(:),type1(:) 2752 integer(kind=MPI_OFFSET_KIND),allocatable :: depl(:),depl1(:),depl_val(:) 2753 real(dp), allocatable :: buf_val(:),val(:) 2754 #endif 2755 2756 ! ********************************************************************* 2757 2758 xval(:,:)=zero ; ierr=0 2759 if(.false.)write(std_out,*)wff%me,n1,n2,displace,spaceComm 2760 2761 #if defined HAVE_MPI_IO 2762 nb=n1*n2 2763 call xmpi_sum(nb,totsize,spaceComm,ierr) 2764 ABI_MALLOC(depl_val,(0:totsize-1)) 2765 ABI_MALLOC(depl,(nb)) 2766 ABI_MALLOC(buf_val,(0:totsize-1)) 2767 ABI_MALLOC(val,(nb)) 2768 2769 !Map displacements 2770 depl_val(0:totsize-1)=-1 2771 do i2=1,n2 2772 do i1=1,n1 2773 ipos=(displace(i2)-1)*n1 + i1-1 2774 depl_val(ipos)=ipos 2775 end do 2776 end do 2777 !To save time, the location describe by array map must be in increasing order 2778 nbval=0 2779 do i1=0,totsize-1 2780 if (depl_val(i1)/=-1) then 2781 nbval=nbval+1 2782 depl(nbval)=depl_val(i1) 2783 end if 2784 end do 2785 2786 !Build MPI datatype for view 2787 ABI_MALLOC(length1,(nbval+2)) 2788 ABI_MALLOC(depl1,(nbval+2)) 2789 ABI_MALLOC(type1,(nbval+2)) 2790 length1(1)=1;depl1(1)=0;type1(1)=MPI_LB 2791 do i1=2,nbval+1 2792 length1(i1) = 1 2793 depl1(i1)= depl(i1-1)*wff%nbOct_dp 2794 type1(i1)= MPI_DOUBLE_PRECISION 2795 end do 2796 length1(nbval+2)=1;depl1(nbval+2)=totsize*wff%nbOct_dp;type1(nbval+2)=MPI_UB 2797 call xmpio_type_struct(nbval+2,length1,depl1,type1,filetype,ierr) 2798 call MPI_TYPE_COMMIT(filetype,ierr) 2799 ABI_FREE(length1) 2800 ABI_FREE(depl1) 2801 ABI_FREE(type1) 2802 2803 !Write data 2804 call MPI_FILE_OPEN(spaceComm,wff%fname,MPI_MODE_RDWR,MPI_INFO_NULL,wfftempo,ierr) 2805 ABI_CHECK_MPI(ierr, sjoin("MPI_FILE_OPEN:", wff%fname)) 2806 call MPI_FILE_SET_VIEW(wfftempo,wff%offwff,MPI_BYTE,filetype,"native",MPI_INFO_NULL,ierr) 2807 call MPI_FILE_READ_ALL(wfftempo,val,nbval,MPI_DOUBLE_PRECISION,statux,ierr) 2808 call MPI_FILE_CLOSE(wfftempo,ierr) 2809 2810 !Retrieve xval 2811 nbval=0 2812 do i1=0,totsize-1 2813 if (depl_val(i1)/=-1) then 2814 nbval=nbval+1 2815 buf_val(i1)=val(nbval) 2816 end if 2817 end do 2818 do i2=1,n2 2819 do i1=1,n1 2820 ipos=(displace(i2)-1)*n1 + i1-1 2821 xval(i1,i2)=buf_val(ipos) 2822 end do 2823 end do 2824 2825 !Update offset 2826 wff%offwff = wff%offwff + totsize*wff%nbOct_dp 2827 2828 !Free memory 2829 call MPI_TYPE_FREE(filetype,ierr) 2830 ABI_FREE(depl) 2831 ABI_FREE(depl_val) 2832 ABI_FREE(buf_val) 2833 ABI_FREE(val) 2834 #endif 2835 2836 end subroutine xderiveRead_dp2d_displ
m_wffile/xderiveRead_int [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveRead_int
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: integer scalar.
INPUTS
(none)
OUTPUT
xval= data buffer ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2161 subroutine xderiveRead_int(wff,xval,ierr) 2162 2163 !Arguments ------------------------------------ 2164 type(wffile_type),intent(inout) :: wff 2165 integer,intent(out) :: xval 2166 integer,intent(out) :: ierr 2167 2168 !Local variables------------------------------- 2169 #if defined HAVE_MPI_IO 2170 integer :: arr_xval(1),statux(MPI_STATUS_SIZE) 2171 #endif 2172 character(len=500) :: msg 2173 2174 ! ********************************************************************* 2175 2176 xval=0; ierr=0 2177 2178 #if defined HAVE_MPI_IO 2179 call MPI_FILE_READ_AT(wff%fhwff,wff%offwff,arr_xval,1,MPI_INTEGER,statux,ierr) 2180 xval=arr_xval(1) 2181 wff%offwff = wff%offwff + wff%nbOct_int 2182 RETURN 2183 #endif 2184 2185 ABI_UNUSED(wff%me) 2186 2187 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 2188 ABI_WARNING(msg) 2189 2190 end subroutine xderiveRead_int
m_wffile/xderiveRead_int1d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveRead_int1d
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: one-dimensional integer arrays.
INPUTS
n1= first dimension of the array spaceComm= MPI communicator
OUTPUT
xval= data buffer array ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2217 subroutine xderiveRead_int1d(wff,xval,n1,spaceComm,ierr) 2218 2219 !Arguments ------------------------------------ 2220 type(wffile_type),intent(inout) :: wff 2221 integer,intent(out) :: xval(:) 2222 integer,intent(in) :: n1,spaceComm 2223 integer,intent(out) :: ierr 2224 2225 !Local variables------------------------------- 2226 #if defined HAVE_MPI_IO 2227 integer(kind=MPI_OFFSET_KIND) :: delim_record,nboct,dispoct(1),posit,totoct 2228 integer :: statux(MPI_STATUS_SIZE) 2229 #endif 2230 character(len=500) :: msg 2231 2232 ! ********************************************************************* 2233 2234 xval(:)=0 ; ierr=0 ! Initialization, for the compiler 2235 if(.false.)write(std_out,*)wff%me,n1,spaceComm 2236 2237 #if defined HAVE_MPI_IO 2238 nboct = wff%nbOct_int * n1 2239 posit = wff%offwff 2240 delim_record = posit - wff%off_recs + wff%lght_recs - wff%nbOct_recMarker 2241 2242 if (delim_record >= nboct) then 2243 ! Compute offset for local part 2244 ! dispoct = sum (nboct, rank=0..me) 2245 if (spaceComm/=MPI_COMM_SELF) then 2246 call MPI_SCAN([nboct],dispoct,1,wff%offset_mpi_type,MPI_SUM,spaceComm,ierr) 2247 posit = posit + dispoct(1) - nboct 2248 end if 2249 call MPI_FILE_READ_AT(wff%fhwff,posit,xval,n1,MPI_INTEGER,statux,ierr) 2250 2251 ! get the total number of bits wrote by processors 2252 if (spaceComm/=MPI_COMM_SELF) then 2253 call xmpi_max(dispoct(1),totoct,spaceComm,ierr) 2254 !call MPI_ALLREDUCE(dispoct,totoct,1,wff%offset_mpi_type,MPI_MAX,spaceComm,ierr) 2255 else 2256 totoct=nboct 2257 end if 2258 else 2259 ierr = 1 2260 nboct = 0 2261 totoct = 0 2262 end if 2263 2264 !new offset 2265 wff%offwff = wff%offwff + totoct 2266 return 2267 #endif 2268 2269 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 2270 ABI_WARNING(msg) 2271 2272 end subroutine xderiveRead_int1d
m_wffile/xderiveRead_int2d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveRead_int2d
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: two-dimensional integer arrays.
INPUTS
n1= first dimension of the array n2= second dimension of the array spaceComm= MPI communicator
OUTPUT
xval= data buffer array ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2300 subroutine xderiveRead_int2d(wff,xval,n1,n2,spaceComm,ierr) 2301 2302 !Arguments ------------------------------------ 2303 type(wffile_type),intent(inout) :: wff 2304 integer,intent(out) :: xval(:,:) 2305 integer,intent(in) :: n1,n2,spaceComm 2306 integer,intent(out) :: ierr 2307 2308 !Local variables------------------------------- 2309 #if defined HAVE_MPI_IO 2310 integer(kind=MPI_OFFSET_KIND) :: delim_record,nboct,dispoct(1),posit,totoct 2311 integer :: statux(MPI_STATUS_SIZE) 2312 #endif 2313 character(len=500) :: msg 2314 2315 ! ********************************************************************* 2316 2317 xval(:,:)=0 ; ierr=0 ! Initialization, for the compiler 2318 if(.false.)write(std_out,*)wff%me,n1,n2,spaceComm 2319 2320 #if defined HAVE_MPI_IO 2321 nboct = wff%nbOct_int * n1 * n2 2322 posit = wff%offwff 2323 delim_record = posit - wff%off_recs + wff%lght_recs - wff%nbOct_recMarker 2324 2325 if (delim_record >= nboct) then 2326 ! Compute offset for local part 2327 ! dispoct = sum (nboct, rank=0..me) 2328 if (spaceComm/=MPI_COMM_SELF) then 2329 call MPI_SCAN([nboct],dispoct,1,wff%offset_mpi_type,MPI_SUM,spaceComm,ierr) 2330 posit = posit + dispoct(1) - nboct 2331 end if 2332 call MPI_FILE_READ_AT(wff%fhwff,posit,xval,n1*n2,MPI_INTEGER,statux,ierr) 2333 2334 ! get the total number of bits wrote by processors 2335 if (spaceComm/=MPI_COMM_SELF) then 2336 call xmpi_max(dispoct(1),totoct,spaceComm,ierr) 2337 !call MPI_ALLREDUCE(dispoct,totoct,1,wff%offset_mpi_type,MPI_MAX,spaceComm,ierr) 2338 else 2339 totoct=nboct 2340 end if 2341 else 2342 ierr = 1 2343 nboct = 0 2344 totoct = 0 2345 end if 2346 2347 !new offset 2348 wff%offwff=wff%offwff + totoct 2349 return 2350 #endif 2351 2352 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 2353 ABI_WARNING(msg) 2354 2355 end subroutine xderiveRead_int2d
m_wffile/xderiveRead_int2d_displ [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveRead_int2d_displ
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: two-dimensional integer arrays.
INPUTS
n1= first dimension of the array n2= second dimension of the array spaceComm= MPI communicator displace= number of elements for the offset
OUTPUT
ierr= exit status, a non-zero value meaning there is an error xval= data buffer array
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2603 subroutine xderiveRead_int2d_displ(wff,xval,n1,n2,spaceComm,displace,ierr) 2604 2605 !Arguments ------------------------------------ 2606 type(wffile_type),intent(inout) :: wff 2607 integer,intent(in) :: n1,n2,spaceComm 2608 integer,intent(out) :: ierr 2609 integer,intent(out):: xval(:,:) 2610 integer,intent(in):: displace(:) 2611 2612 !Local variables------------------------------- 2613 #if defined HAVE_MPI_IO 2614 !scalars 2615 integer :: filetype,i1,i2,ipos,nb,nbval,totsize,wfftempo 2616 !arrays 2617 integer :: statux(MPI_STATUS_SIZE) 2618 integer,allocatable :: buf_val(:),length1(:),type1(:),val(:) 2619 integer(kind=MPI_OFFSET_KIND),allocatable :: depl(:),depl1(:),depl_val(:) 2620 #endif 2621 character(len=500) :: msg 2622 2623 ! ********************************************************************* 2624 2625 xval(:,:)=0 ; ierr=0 2626 if(.false.)write(std_out,*)wff%me,n1,n2,spaceComm,displace 2627 2628 #if defined HAVE_MPI_IO 2629 nb=n1*n2 2630 call xmpi_sum(nb,totsize,spaceComm,ierr) 2631 ABI_MALLOC(depl_val,(0:totsize-1)) 2632 ABI_MALLOC(depl,(nb)) 2633 ABI_MALLOC(buf_val,(0:totsize-1)) 2634 ABI_MALLOC(val,(nb)) 2635 2636 !Map displacements 2637 depl_val(0:totsize-1)=-1 2638 do i2=1,n2 2639 do i1=1,n1 2640 ipos=(displace(i2)-1)*n1 + i1-1 2641 depl_val(ipos)=ipos 2642 end do 2643 end do 2644 !To save time, the location describe by array map must be in increasing order 2645 nbval=0 2646 do i1=0,totsize-1 2647 if (depl_val(i1)/=-1) then 2648 nbval=nbval+1 2649 depl(nbval)=depl_val(i1) 2650 end if 2651 end do 2652 2653 !Build MPI datatype for view 2654 ABI_MALLOC(length1,(nbval+2)) 2655 ABI_MALLOC(depl1,(nbval+2)) 2656 ABI_MALLOC(type1,(nbval+2)) 2657 length1(1)=1;depl1(1)=0;type1(1)=MPI_LB 2658 do i1=2,nbval+1 2659 length1(i1) = 1 2660 depl1(i1)= depl(i1-1)*wff%nbOct_int 2661 type1(i1)= MPI_INTEGER 2662 end do 2663 length1(nbval+2)=1;depl1(nbval+2)=totsize*wff%nbOct_int;type1(nbval+2)=MPI_UB 2664 call xmpio_type_struct(nbval+2,length1,depl1,type1,filetype,ierr) 2665 call MPI_TYPE_COMMIT(filetype,ierr) 2666 ABI_FREE(length1) 2667 ABI_FREE(depl1) 2668 ABI_FREE(type1) 2669 2670 !Write data 2671 call MPI_FILE_OPEN(spaceComm,wff%fname,MPI_MODE_RDWR,MPI_INFO_NULL,wfftempo,ierr) 2672 ABI_CHECK_MPI(ierr, sjoin("MPI_FILE_OPEN:", wff%fname)) 2673 call MPI_FILE_SET_VIEW(wfftempo,wff%offwff,MPI_BYTE,filetype,"native",MPI_INFO_NULL,ierr) 2674 call MPI_FILE_READ_ALL(wfftempo,val,nbval,MPI_INTEGER,statux,ierr) 2675 call MPI_FILE_CLOSE(wfftempo,ierr) 2676 2677 !Retrieve xval 2678 nbval=0 2679 do i1=0,totsize-1 2680 if (depl_val(i1)/=-1) then 2681 nbval=nbval+1 2682 buf_val(i1)=val(nbval) 2683 end if 2684 end do 2685 do i2=1,n2 2686 do i1=1,n1 2687 ipos=(displace(i2)-1)*n1 + i1-1 2688 xval(i1,i2)=buf_val(ipos) 2689 end do 2690 end do 2691 2692 !Update offset 2693 wff%offwff = wff%offwff + totsize*wff%nbOct_int 2694 2695 !Free memory 2696 call MPI_TYPE_FREE(filetype,ierr) 2697 ABI_FREE(depl) 2698 ABI_FREE(depl_val) 2699 ABI_FREE(buf_val) 2700 ABI_FREE(val) 2701 return 2702 #endif 2703 2704 write(msg,'(a,i0,a)')' The value of wff%iomode=',wff%iomode,' is not allowed.' 2705 ABI_WARNING(msg) 2706 2707 2708 end subroutine xderiveRead_int2d_displ
m_wffile/xderiveReadVal_char [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveReadVal_char
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: character string.
INPUTS
n= number of elements in the array
OUTPUT
xval= data buffer array ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2861 subroutine xderiveReadVal_char(wff,xval,n,ierr) 2862 2863 !Arguments ------------------------------------ 2864 type(wffile_type),intent(inout) :: wff 2865 integer,intent(in) :: n 2866 integer,intent(out) :: ierr 2867 character(len=*),intent(out),target :: xval 2868 2869 !Local variables------------------------------- 2870 #if defined HAVE_MPI_IO 2871 integer :: statux(MPI_STATUS_SIZE) 2872 character,pointer :: arr_xval(:) 2873 type(c_ptr) :: cptr 2874 #endif 2875 2876 ! ********************************************************************* 2877 2878 xval=' ' ; ierr=0 2879 if(.false.)write(std_out,*)wff%me,n 2880 2881 #if defined HAVE_MPI_IO 2882 cptr=c_loc(xval) ; call c_f_pointer(cptr,arr_xval,[n]) 2883 call MPI_FILE_READ_AT(wff%fhwff,wff%offwff,arr_xval,n,MPI_CHARACTER,statux,ierr) 2884 wff%offwff = wff%offwff + wff%nbOct_ch * n 2885 #endif 2886 2887 end subroutine xderiveReadVal_char
m_wffile/xderiveRRecEnd [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveRRecEnd
FUNCTION
Initializes the end-of-record offset for MPI/IO.
INPUTS
me_proc= (optional argument) index of current proc
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
NOTES
We assume that: wff%off_recs contains the position of the beginning of the record
SOURCE
742 subroutine xderiveRRecEnd(wff,ierr) 743 744 !Arguments ------------------------------------ 745 integer,intent(out) :: ierr 746 type(wffile_type),intent(inout) :: wff 747 748 ! ************************************************************************* 749 750 ierr=0 751 #if defined HAVE_MPI_IO 752 !Define offset end of record 753 wff%offwff = wff%off_recs + wff%lght_recs + 2*wff%nbOct_recMarker 754 #endif 755 756 RETURN 757 ABI_UNUSED(wff%me) 758 759 end subroutine xderiveRRecEnd
m_wffile/xderiveRRecInit [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveRRecInit
FUNCTION
Initializes the record length for MPI/IO.
INPUTS
me_proc= (optional argument) index of current proc
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
NOTES
We assume that: wff%offwff contains the position of the beginning of the record
SOURCE
786 subroutine xderiveRRecInit(wff,ierr) 787 788 !Arguments ------------------------------------ 789 type(wffile_type),intent(inout) :: wff 790 integer,intent(out) :: ierr 791 792 !Local variables------------------------------- 793 #if defined HAVE_MPI_IO 794 integer(kind=MPI_OFFSET_KIND) :: delim_record,posit 795 #endif 796 797 ! ************************************************************************* 798 799 ierr=0 800 801 #if defined HAVE_MPI_IO 802 wff%off_recs = wff%offwff 803 804 !Read the length of the record 805 posit=wff%off_recs 806 call rwRecordMarker(1,posit,delim_record,wff,ierr) 807 808 wff%lght_recs = delim_record 809 wff%offwff = wff%offwff + wff%nbOct_recMarker 810 #endif 811 812 RETURN 813 ABI_UNUSED(wff%me) 814 815 end subroutine xderiveRRecInit
m_wffile/xderiveWRecEnd [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWRecEnd
FUNCTION
Writes the first and last wavefunction block marker using MPI/IO
INPUTS
me_proc= (optional argument) index of current proc
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
NOTES
We assume that: wff%offwff contains the position of the end of the record wff%off_recs contains the position of the beginning of the record
SOURCE
611 subroutine xderiveWRecEnd(wff,ierr,me_proc) 612 613 !Arguments ------------------------------------ 614 type(wffile_type),intent(inout) :: wff 615 integer,intent(in),optional :: me_proc 616 integer,intent(out) :: ierr 617 618 !Local variables------------------------------- 619 #if defined HAVE_MPI_IO 620 !scalars 621 integer :: me 622 integer(kind=MPI_OFFSET_KIND) :: delim_record,posit 623 !arrays 624 #endif 625 626 ! ************************************************************************* 627 628 ierr=0 629 630 #if defined HAVE_MPI_IO 631 me=-1;if (present(me_proc)) me=me_proc 632 if (me==-1.or.me==0) then 633 634 delim_record=wff%offwff-wff%off_recs-wff%nbOct_recMarker 635 636 ! Write the first word of the record 637 posit=wff%off_recs 638 call rwRecordMarker(2,posit,delim_record,wff,ierr) 639 640 ! Write the last word of the record 641 posit=wff%offwff 642 call rwRecordMarker(2,posit,delim_record,wff,ierr) 643 644 end if 645 646 wff%offwff = wff%offwff + wff%nbOct_recMarker 647 #endif 648 649 RETURN 650 ABI_UNUSED((/wff%me,me_proc/)) 651 652 end subroutine xderiveWRecEnd
m_wffile/xderiveWRecInit [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWRecInit
FUNCTION
Writes the first wavefunction block marker using MPI/IO.
INPUTS
me_proc= (optional argument) index of current proc
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
NOTES
We assume that: wff%offwff contains the position of the beginning of the record
SOURCE
679 subroutine xderiveWRecInit(wff,ierr,me_proc) 680 681 !Arguments ------------------------------------ 682 type(wffile_type),intent(inout) :: wff 683 integer,intent(in),optional :: me_proc 684 integer,intent(out) :: ierr 685 686 !Local variables------------------------------- 687 #if defined HAVE_MPI_IO 688 !scalars 689 integer :: me 690 integer(kind=MPI_OFFSET_KIND) :: delim_record,posit 691 !arrays 692 #endif 693 694 ! ************************************************************************* 695 696 ierr=0 697 698 #if defined HAVE_MPI_IO 699 me=-1;if (present(me_proc)) me=me_proc 700 if (me==-1.or.me==0) then 701 702 ! Write the first word of the record 703 posit=wff%offwff;delim_record=0 704 call rwRecordMarker(2,posit,delim_record,wff,ierr) 705 706 end if 707 708 wff%off_recs = wff%offwff 709 wff%offwff = wff%offwff + wff%nbOct_recMarker 710 #endif 711 712 RETURN 713 ABI_UNUSED((/wff%me,me_proc/)) 714 715 end subroutine xderiveWRecInit
m_wffile/xderiveWrite_char [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWrite_char
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: character string.
INPUTS
xval= data buffer array n= number of elements in the string
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
3630 subroutine xderiveWrite_char(wff,xval,n,ierr) 3631 3632 !Arguments ------------------------------------ 3633 type(wffile_type),intent(inout) :: wff 3634 integer,intent(in) :: n 3635 integer,intent(out) :: ierr 3636 character(len=*),intent(in),target :: xval 3637 3638 !Local variables------------------------------- 3639 #if defined HAVE_MPI_IO 3640 integer :: statux(MPI_STATUS_SIZE) 3641 character,pointer :: buf(:) 3642 type(c_ptr) :: cptr 3643 #endif 3644 ! ********************************************************************* 3645 3646 ierr=0 3647 if(.false.)write(std_out,*)wff%me,xval,n 3648 3649 #if defined HAVE_MPI_IO 3650 cptr=c_loc(xval) ; call c_f_pointer(cptr,buf,[n]) 3651 call MPI_FILE_WRITE_AT(wff%fhwff,wff%offwff,buf,n,MPI_CHARACTER,statux,ierr) 3652 wff%offwff = wff%offwff + wff%nbOct_ch * n 3653 #endif 3654 3655 end subroutine xderiveWrite_char
m_wffile/xderiveWrite_dp [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWrite_dp
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: double precision scalar.
INPUTS
xval= data buffer
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
3176 subroutine xderiveWrite_dp(wff,xval,ierr) 3177 3178 !Arguments ------------------------------------ 3179 integer,intent(out) :: ierr 3180 real(dp),intent(in):: xval 3181 type(wffile_type),intent(inout) :: wff 3182 3183 !Local variables------------------------------- 3184 #if defined HAVE_MPI_IO 3185 integer :: statux(MPI_STATUS_SIZE) 3186 #endif 3187 ! ********************************************************************* 3188 3189 ierr=0 3190 if(.false.)write(std_out,*)wff%me,xval 3191 #if defined HAVE_MPI_IO 3192 call MPI_FILE_WRITE_AT(wff%fhwff,wff%offwff,[xval],1,MPI_DOUBLE_PRECISION,statux,ierr) 3193 wff%offwff = wff%offwff+wff%nbOct_dp 3194 #endif 3195 3196 end subroutine xderiveWrite_dp
m_wffile/xderiveWrite_dp1d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWrite_dp1d
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: one-dimensional double precision arrays.
INPUTS
n1= first dimension of the array spaceComm= MPI communicator xval= data buffer array
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
3223 subroutine xderiveWrite_dp1d(wff,xval,n1,spaceComm,ierr) 3224 3225 !Arguments ------------------------------------ 3226 integer,intent(in) :: n1,spaceComm 3227 integer,intent(out) :: ierr 3228 real(dp),intent(in):: xval(:) 3229 type(wffile_type),intent(inout) :: wff 3230 3231 !Local variables------------------------------- 3232 #if defined HAVE_MPI_IO 3233 integer(kind=MPI_OFFSET_KIND) :: nboct,dispoct(1),posit,totoct 3234 integer :: statux(MPI_STATUS_SIZE) 3235 #endif 3236 3237 ! ********************************************************************* 3238 3239 ierr=0 3240 if(.false.)write(std_out,*)wff%me,n1,spaceComm,xval 3241 #if defined HAVE_MPI_IO 3242 nboct = n1*wff%nbOct_dp 3243 posit = wff%offwff 3244 !dispoct = sum (nboct, rank = 0..me) 3245 if (spaceComm/=MPI_COMM_SELF) then 3246 call MPI_SCAN([nboct],dispoct,1,wff%offset_mpi_type,MPI_SUM,spaceComm,ierr) 3247 posit = posit + dispoct(1) - nboct 3248 end if 3249 call MPI_FILE_WRITE_AT(wff%fhwff,posit,xval,n1,MPI_DOUBLE_PRECISION,statux,ierr) 3250 !Gather the biggest offset 3251 if (spaceComm/=MPI_COMM_SELF) then 3252 call xmpi_max(dispoct(1),totoct,spaceComm,ierr) 3253 !call MPI_ALLREDUCE(dispoct,totoct,1,wff%offset_mpi_type,MPI_MAX,spaceComm,ierr) 3254 else 3255 totoct=nboct 3256 end if 3257 wff%offwff = wff%offwff + totoct 3258 #endif 3259 3260 end subroutine xderiveWrite_dp1d
m_wffile/xderiveWrite_dp2d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWrite_dp2d
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: double precision two-dimensional arrays.
INPUTS
n1= first dimension of the array n2= second dimension of the array spaceComm= MPI communicator xval= data buffer array
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
3287 subroutine xderiveWrite_dp2d(wff,xval,n1,n2,spaceComm,ierr) 3288 3289 !Arguments ------------------------------------ 3290 integer,intent(in) :: n1,n2,spaceComm 3291 integer,intent(out) :: ierr 3292 real(dp),intent(in):: xval(:,:) 3293 type(wffile_type),intent(inout) :: wff 3294 3295 !Local variables------------------------------- 3296 #if defined HAVE_MPI_IO 3297 integer(kind=MPI_OFFSET_KIND) :: nboct,dispoct(1),posit,totoct 3298 integer :: statux(MPI_STATUS_SIZE) 3299 #endif 3300 3301 ! ********************************************************************* 3302 3303 ierr=0 3304 if(.false.)write(std_out,*)wff%me,xval,n1,n2,spaceComm 3305 3306 #if defined HAVE_MPI_IO 3307 nboct = n1*n2*wff%nbOct_dp 3308 posit = wff%offwff 3309 !dispoct = sum(nboct, rank=0..me) 3310 if (spaceComm/=MPI_COMM_SELF) then 3311 call MPI_SCAN([nboct],dispoct,1,wff%offset_mpi_type,MPI_SUM,spaceComm,ierr) 3312 posit = posit + dispoct(1) - nboct 3313 end if 3314 call MPI_FILE_WRITE_AT(wff%fhwff,posit,xval,n1*n2,MPI_DOUBLE_PRECISION,statux,ierr) 3315 posit = posit + nboct 3316 !gather the biggest offset 3317 if (spaceComm/=MPI_COMM_SELF) then 3318 call xmpi_max(dispoct(1),totoct,spaceComm,ierr) 3319 !call MPI_ALLREDUCE(dispoct,totoct,1,wff%offset_mpi_type,MPI_MAX,spaceComm,ierr) 3320 else 3321 totoct=nboct 3322 end if 3323 wff%offwff = wff%offwff + totoct 3324 #endif 3325 3326 end subroutine xderiveWrite_dp2d
m_wffile/xderiveWrite_dp2d_displ [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWrite_dp2d_displ
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: two-dimensional double precision arrays.
INPUTS
n1= first dimension of the array n2= second dimension of the array spaceComm= MPI communicator xval= data buffer array displace= number of elements for the offset
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
3516 subroutine xderiveWrite_dp2d_displ(wff,xval,n1,n2,spaceComm,displace,ierr) 3517 3518 !Arguments ------------------------------------ 3519 integer,intent(in) :: n1,n2,spaceComm 3520 integer,intent(out) :: ierr 3521 integer,intent(in):: displace(:) 3522 real(dp),intent(in) :: xval(:,:) 3523 type(wffile_type),intent(inout) :: wff 3524 3525 !Local variables------------------------------- 3526 #if defined HAVE_MPI_IO 3527 !scalars 3528 integer :: filetype,i1,i2,ipos,nb,nbval,totsize,wfftempo 3529 !arrays 3530 integer :: statux(MPI_STATUS_SIZE) 3531 integer, allocatable :: length1(:),type1(:) 3532 integer(kind=XMPI_ADDRESS_KIND),allocatable :: depl(:),depl1(:),depl_val(:) 3533 real(dp),allocatable :: buf_val(:),val(:) 3534 #endif 3535 3536 ! ********************************************************************* 3537 3538 ierr=0 3539 if(.false.)write(std_out,*)wff%me,xval,n1,n2,spaceComm,displace 3540 3541 #if defined HAVE_MPI_IO 3542 nb = n1*n2 3543 call xmpi_sum(nb,totsize,spaceComm,ierr) 3544 ABI_MALLOC(depl_val,(0:totsize-1)) 3545 ABI_MALLOC(depl,(nb)) 3546 ABI_MALLOC(buf_val,(0:totsize-1)) 3547 ABI_MALLOC(val,(nb)) 3548 3549 !Map displacements 3550 !Put xval in a buffer at its position 3551 depl_val(0:totsize-1)=-1 3552 do i2=1,n2 3553 do i1=1,n1 3554 ! ipos location of xval(i1,i2) in the array associated with record to be written 3555 ipos=(displace(i2)-1)*n1 + i1-1 3556 buf_val(ipos) = xval(i1,i2) 3557 depl_val(ipos) = ipos 3558 end do 3559 end do 3560 !To save time, the location describe by array map must be in increasing order 3561 nbval=0 3562 do i1=0,totsize-1 3563 if (depl_val(i1)/=-1) then 3564 nbval=nbval+1 3565 val(nbval)=buf_val(i1) 3566 depl(nbval)=depl_val(i1) 3567 end if 3568 end do 3569 3570 !Build MPI datatype for view 3571 ABI_MALLOC(length1,(nbval+2)) 3572 ABI_MALLOC(depl1,(nbval+2)) 3573 ABI_MALLOC(type1,(nbval+2)) 3574 length1(1)=1;depl1(1)=0;type1(1)=MPI_LB 3575 do i1=2,nbval+1 3576 length1(i1) = 1 3577 depl1(i1)= depl(i1-1)*wff%nbOct_dp 3578 type1(i1)= MPI_DOUBLE_PRECISION 3579 end do 3580 length1(nbval+2)=1;depl1(nbval+2)=totsize*wff%nbOct_dp;type1(nbval+2)=MPI_UB 3581 call xmpio_type_struct(nbval+2,length1,depl1,type1,filetype,ierr) 3582 call MPI_TYPE_COMMIT(filetype,ierr) 3583 ABI_FREE(length1) 3584 ABI_FREE(depl1) 3585 ABI_FREE(type1) 3586 3587 !Write data 3588 call MPI_FILE_OPEN(spaceComm,wff%fname,MPI_MODE_RDWR,MPI_INFO_NULL,wfftempo,ierr) 3589 ABI_CHECK_MPI(ierr, sjoin("MPI_FILE_OPEN:", wff%fname)) 3590 call MPI_FILE_SET_VIEW(wfftempo,wff%offwff,MPI_BYTE,filetype,"native",MPI_INFO_NULL,ierr) 3591 call MPI_FILE_WRITE_ALL(wfftempo,val,nbval,MPI_DOUBLE_PRECISION,statux,ierr) 3592 call MPI_FILE_CLOSE(wfftempo,ierr) 3593 3594 wff%offwff = wff%offwff + totsize*wff%nbOct_dp 3595 3596 !Free memory 3597 call MPI_TYPE_FREE(filetype,ierr) 3598 ABI_FREE(depl) 3599 ABI_FREE(depl_val) 3600 ABI_FREE(buf_val) 3601 ABI_FREE(val) 3602 #endif 3603 3604 end subroutine xderiveWrite_dp2d_displ
m_wffile/xderiveWrite_dp2d_seq [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWrite_dp2d_seq
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: double precision two-dimensional arrays.
INPUTS
xval= data buffer array
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
3350 subroutine xderiveWrite_dp2d_seq(wff,xval,ierr) 3351 3352 !Arguments ------------------------------------ 3353 integer,intent(out) :: ierr 3354 real(dp),intent(in):: xval(:,:) 3355 type(wffile_type),intent(inout) :: wff 3356 3357 !Local variables------------------------------- 3358 #if defined HAVE_MPI_IO 3359 integer :: n1,n2 3360 integer :: statux(MPI_STATUS_SIZE) 3361 #endif 3362 ! ********************************************************************* 3363 3364 ierr=0 3365 if(.false.)write(std_out,*)wff%me,xval 3366 #if defined HAVE_MPI_IO 3367 n1=size(xval,1);n2=size(xval,2) 3368 call MPI_FILE_WRITE_AT(wff%fhwff,wff%offwff,xval,n1*n2,MPI_DOUBLE_PRECISION,statux,ierr) 3369 wff%offwff = wff%offwff+wff%nbOct_dp*n1*n2 3370 #endif 3371 3372 end subroutine xderiveWrite_dp2d_seq
m_wffile/xderiveWrite_int [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWrite_int
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: integer scalar.
INPUTS
xval= data buffer
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2999 subroutine xderiveWrite_int(wff,xval,ierr) 3000 3001 !Arguments ------------------------------------ 3002 integer,intent(out) :: ierr 3003 integer,intent(in):: xval 3004 type(wffile_type),intent(inout) :: wff 3005 3006 !Local variables------------------------------- 3007 #if defined HAVE_MPI_IO 3008 integer :: statux(MPI_STATUS_SIZE) 3009 #endif 3010 ! ********************************************************************* 3011 3012 ierr=0 3013 if(.false.)write(std_out,*)wff%me,xval 3014 #if defined HAVE_MPI_IO 3015 call MPI_FILE_WRITE_AT(wff%fhwff,wff%offwff,[xval],1,MPI_INTEGER,statux,ierr) 3016 wff%offwff = wff%offwff+wff%nbOct_int 3017 #endif 3018 3019 end subroutine xderiveWrite_int
m_wffile/xderiveWrite_int1d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWrite_int1d
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: one-dimensional integer arrays.
INPUTS
n1= first dimension of the array spaceComm= MPI communicator xval= data buffer array
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
3046 subroutine xderiveWrite_int1d(wff,xval,n1,spaceComm,ierr) 3047 3048 !Arguments ------------------------------------ 3049 integer,intent(in) :: n1,spaceComm 3050 integer,intent(out) :: ierr 3051 integer,intent(in):: xval(:) 3052 type(wffile_type),intent(inout) :: wff 3053 3054 !Local variables------------------------------- 3055 #if defined HAVE_MPI_IO 3056 integer(kind=MPI_OFFSET_KIND) :: nboct,dispoct(1),posit,totoct 3057 integer :: statux(MPI_STATUS_SIZE) 3058 #endif 3059 ! ********************************************************************* 3060 3061 ierr=0 3062 if(.false.)write(std_out,*)wff%me,n1,spaceComm,xval 3063 #if defined HAVE_MPI_IO 3064 nboct = n1*wff%nbOct_int 3065 posit = wff%offwff 3066 3067 !dispoct = sum (nboct, rank=0..me) 3068 if (spaceComm/=MPI_COMM_SELF) then 3069 call MPI_SCAN([nboct],dispoct,1,wff%offset_mpi_type,MPI_SUM,spaceComm,ierr) 3070 posit = posit + dispoct(1) - nboct 3071 end if 3072 call MPI_FILE_WRITE_AT(wff%fhwff,posit,xval,n1,MPI_INTEGER,statux,ierr) 3073 !gather the bigest offset 3074 3075 if (spaceComm/=MPI_COMM_SELF) then 3076 call xmpi_max(dispoct(1),totoct,spaceComm,ierr) 3077 !call MPI_ALLREDUCE(dispoct,totoct,1,wff%offset_mpi_type,MPI_MAX,spaceComm,ierr) 3078 else 3079 totoct=nboct 3080 end if 3081 wff%offwff = wff%offwff + totoct 3082 3083 !Disable old code 3084 #endif 3085 3086 end subroutine xderiveWrite_int1d
m_wffile/xderiveWrite_int2d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWrite_int2d
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: two-dimensional integer arrays.
INPUTS
n1= first dimension of the array n2= second dimension of the array spaceComm= MPI communicator xval= data buffer array
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
3113 subroutine xderiveWrite_int2d(wff,xval,n1,n2,spaceComm,ierr) 3114 3115 !Arguments ------------------------------------ 3116 integer,intent(in) :: n1,n2,spaceComm 3117 integer,intent(out) :: ierr 3118 integer,intent(in):: xval(:,:) 3119 type(wffile_type),intent(inout) :: wff 3120 3121 !Local variables------------------------------- 3122 #if defined HAVE_MPI_IO 3123 integer(kind=MPI_OFFSET_KIND) :: nboct,dispoct(1),posit,totoct 3124 integer :: statux(MPI_STATUS_SIZE) 3125 #endif 3126 3127 ! ********************************************************************* 3128 3129 ierr=0 3130 if(.false.)write(std_out,*)wff%me,n1,n2,spaceComm,xval 3131 #if defined HAVE_MPI_IO 3132 nboct = n1*n2*wff%nbOct_int 3133 posit = wff%offwff 3134 3135 !dispoct = sum(nboct, rank=0..me) 3136 if (spaceComm/=MPI_COMM_SELF) then 3137 call MPI_SCAN([nboct],dispoct,1,wff%offset_mpi_type,MPI_SUM,spaceComm,ierr) 3138 posit = posit + dispoct(1) - nboct 3139 end if 3140 call MPI_FILE_WRITE_AT(wff%fhwff,posit,xval,n1*n2,MPI_INTEGER,statux,ierr) 3141 !gather the biggest offset 3142 if (spaceComm/=MPI_COMM_SELF) then 3143 call xmpi_max(dispoct(1),totoct,spaceComm,ierr) 3144 !call MPI_ALLREDUCE(dispoct,totoct,1,wff%offset_mpi_type,MPI_MAX,spaceComm,ierr) 3145 else 3146 totoct=nboct 3147 end if 3148 wff%offwff = wff%offwff + totoct 3149 #endif 3150 3151 end subroutine xderiveWrite_int2d
m_wffile/xderiveWrite_int2d_displ [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xderiveWrite_int2d_displ
FUNCTION
Generic routine to read/write wf files with MPI I/O. Target: two-dimensional integer arrays.
INPUTS
n1= first dimension of the array n2= second dimension of the array spaceComm= MPI communicator xval= data buffer array displace= number of elements for the offset
OUTPUT
ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
3401 subroutine xderiveWrite_int2d_displ(wff,xval,n1,n2,spaceComm,displace,ierr) 3402 3403 !Arguments ------------------------------------ 3404 integer,intent(in) :: n1,n2,spaceComm 3405 integer,intent(out) :: ierr 3406 integer,intent(in):: displace(:),xval(:,:) 3407 type(wffile_type),intent(inout) :: wff 3408 3409 !Local variables------------------------------- 3410 #if defined HAVE_MPI_IO 3411 !scalars 3412 integer :: filetype,i1,i2,ipos,nb,nbval,totsize,wfftempo 3413 !arrays 3414 integer :: statux(MPI_STATUS_SIZE) 3415 integer, allocatable :: buf_val(:),length1(:),type1(:),val(:) 3416 integer(kind=MPI_OFFSET_KIND),allocatable :: depl(:),depl1(:),depl_val(:) 3417 #endif 3418 3419 ! ********************************************************************* 3420 3421 ierr=0 3422 if(.false.)write(std_out,*)wff%me,xval,n1,n2,spaceComm,displace 3423 3424 #if defined HAVE_MPI_IO 3425 nb = n1*n2 3426 call xmpi_sum(nb,totsize,spaceComm,ierr) 3427 ABI_MALLOC(depl_val,(0:totsize-1)) 3428 ABI_MALLOC(depl,(nb)) 3429 ABI_MALLOC(buf_val,(0:totsize-1)) 3430 ABI_MALLOC(val,(nb)) 3431 3432 !Map displacements 3433 !Put xval in a buffer at its position 3434 depl_val(0:totsize-1)=-1 3435 do i2=1,n2 3436 do i1=1,n1 3437 ! ipos location of xval(i1,i2) in the array associated with record to be written 3438 ipos=(displace(i2)-1)*n1 + i1-1 3439 buf_val(ipos) = xval(i1,i2) 3440 depl_val(ipos) = ipos 3441 end do 3442 end do 3443 !To save time, the location describe by array map must be in increasing order 3444 nbval=0 3445 do i1=0,totsize-1 3446 if (depl_val(i1)/=-1) then 3447 nbval=nbval+1 3448 val(nbval)=buf_val(i1) 3449 depl(nbval)=depl_val(i1) 3450 end if 3451 end do 3452 3453 !Build MPI datatype for view 3454 ABI_MALLOC(length1,(nbval+2)) 3455 ABI_MALLOC(depl1,(nbval+2)) 3456 ABI_MALLOC(type1,(nbval+2)) 3457 length1(1)=1;depl1(1)=0;type1(1)=MPI_LB 3458 do i1=2,nbval+1 3459 length1(i1) = 1 3460 depl1(i1)= depl(i1-1)*wff%nbOct_int 3461 type1(i1)= MPI_INTEGER 3462 end do 3463 length1(nbval+2)=1;depl1(nbval+2)=totsize*wff%nbOct_int;type1(nbval+2)=MPI_UB 3464 call xmpio_type_struct(nbval+2,length1,depl1,type1,filetype,ierr) 3465 call MPI_TYPE_COMMIT(filetype,ierr) 3466 ABI_FREE(length1) 3467 ABI_FREE(depl1) 3468 ABI_FREE(type1) 3469 3470 !Write data 3471 call MPI_FILE_OPEN(spaceComm,wff%fname,MPI_MODE_RDWR,MPI_INFO_NULL,wfftempo,ierr) 3472 ABI_CHECK_MPI(ierr, sjoin("MPI_FILE_OPEN:", wff%fname)) 3473 call MPI_FILE_SET_VIEW(wfftempo,wff%offwff,MPI_BYTE,filetype,"native",MPI_INFO_NULL,ierr) 3474 call MPI_FILE_WRITE_ALL(wfftempo,val,nbval,MPI_INTEGER,statux,ierr) 3475 call MPI_FILE_CLOSE(wfftempo,ierr) 3476 3477 !Update offset 3478 wff%offwff = wff%offwff + totsize*wff%nbOct_int 3479 3480 !Free memory 3481 call MPI_TYPE_FREE(filetype,ierr) 3482 ABI_FREE(depl) 3483 ABI_FREE(depl_val) 3484 ABI_FREE(buf_val) 3485 ABI_FREE(val) 3486 #endif 3487 3488 end subroutine xderiveWrite_int2d_displ
m_wffile/xmoveOff [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xmoveOff
FUNCTION
In case of MPI I/O, move the offset of a WF file
INPUTS
[n_int] = number if integers to skip [n_dp] = number if double precision reals to skip [n_ch] = number if characters to skip [n_mark]= number if record markers to skip
OUTPUT
SIDE EFFECTS
wff=<type(wffile_type)>=structured info for reading/writing
SOURCE
562 subroutine xmoveOff(wff,n_int,n_dp,n_ch,n_mark) 563 564 !Arguments ------------------------------------ 565 integer,intent(in),optional :: n_int,n_dp,n_ch,n_mark 566 type(wffile_type),intent(inout) :: wff 567 568 ! ************************************************************************* 569 570 #if defined HAVE_MPI_IO 571 if (present(n_int) ) wff%offwff=wff%offwff+n_int *wff%nbOct_int 572 if (present(n_dp) ) wff%offwff=wff%offwff+n_dp *wff%nbOct_dp 573 if (present(n_ch) ) wff%offwff=wff%offwff+n_ch *wff%nbOct_ch 574 if (present(n_mark)) wff%offwff=wff%offwff+n_mark*wff%nbOct_recMarker 575 #else 576 !This section should not be used... 577 if (present(n_int) .and.(.false.)) write(std_out,*) n_int 578 if (present(n_dp) .and.(.false.)) write(std_out,*) n_dp 579 if (present(n_ch) .and.(.false.)) write(std_out,*) n_ch 580 if (present(n_mark).and.(.false.)) write(std_out,*) n_mark 581 #endif 582 583 end subroutine xmoveOff
m_wffile/xmpi_read_int2d [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xmpi_read_int2d
FUNCTION
Generic routine to read arrays with MPI I/O. Target: integer two-dimensional arrays.
INPUTS
sc_mode= xmpio_single ==> Local reading. xmpio_collective ==> Collective reading. spaceComm= MPI communicator
OUTPUT
xval= data buffer array ierr= exit status, a non-zero value meaning there is an error
SIDE EFFECTS
wff= structured info for reading/writing the wavefunctions
SOURCE
2915 subroutine xmpi_read_int2d(wff,xval,spaceComm,sc_mode,ierr) 2916 2917 !Arguments ------------------------------------ 2918 !scalars 2919 integer,intent(in) :: spaceComm,sc_mode 2920 integer,intent(out) :: ierr 2921 type(wffile_type),intent(inout) :: wff 2922 !array 2923 integer,intent(out) :: xval(:,:) 2924 2925 !Local variables------------------------------- 2926 integer :: n1,n2 2927 #ifdef HAVE_MPI_IO 2928 integer(kind=MPI_OFFSET_KIND) :: delim_record,nboct,posit,totoct 2929 character(len=500) :: msg 2930 !arrays 2931 integer :: statux(MPI_STATUS_SIZE) 2932 #endif 2933 2934 ! ********************************************************************* 2935 2936 ierr=0 2937 n1 = SIZE(xval,DIM=1) 2938 n2 = SIZE(xval,DIM=2) 2939 2940 #ifdef HAVE_MPI_IO 2941 nboct = wff%nbOct_int * n1 *n2 2942 posit = wff%offwff 2943 delim_record = posit - wff%off_recs + wff%lght_recs - wff%nbOct_recMarker 2944 2945 if (delim_record >= nboct) then 2946 2947 select case (sc_mode) 2948 case (xmpio_single) 2949 call MPI_FILE_READ_AT(wff%fhwff,posit,xval,n1*n2,MPI_INTEGER,statux,ierr) 2950 2951 case (xmpio_collective) 2952 call MPI_FILE_READ_AT_ALL(wff%fhwff,posit,xval,n1*n2,MPI_INTEGER,statux,ierr) 2953 2954 case default 2955 write(msg,('(a,i0)'))" Wrong value for sc_mode: ",sc_mode 2956 ABI_ERROR(msg) 2957 end select 2958 2959 totoct=nboct 2960 else 2961 write(msg,('(a,2(i0,1x))'))" delim_record < nboct: ",delim_record,nboct 2962 ABI_WARNING(msg) 2963 ierr=MPI_ERR_UNKNOWN 2964 totoct=0 2965 end if 2966 ! 2967 !Increment the offset. 2968 wff%offwff=wff%offwff + totoct 2969 #endif 2970 2971 RETURN 2972 ABI_UNUSED(xval(1,1)) 2973 ABI_UNUSED((/wff%me,spaceComm,sc_mode/)) 2974 2975 end subroutine xmpi_read_int2d
m_wffile/xnullifyOff [ Functions ]
[ Top ] [ m_wffile ] [ Functions ]
NAME
xnullifyOff
FUNCTION
In case of MPI I/O, nullify the offset of a WF file
INPUTS
OUTPUT
SIDE EFFECTS
wff=<type(wffile_type)>=structured info for reading/writing
SOURCE
521 subroutine xnullifyOff(wff) 522 523 !Arguments ------------------------------------ 524 type(wffile_type),intent(inout) :: wff 525 526 ! ************************************************************************* 527 528 #if defined HAVE_MPI_IO 529 wff%offwff = 0 530 wff%off_recs = 0 531 wff%lght_recs = 0 532 #endif 533 534 RETURN 535 ABI_UNUSED(wff%me) 536 537 end subroutine xnullifyOff