TABLE OF CONTENTS


ABINIT/m_wfk [ Modules ]

[ Top ] [ Modules ]

NAME

  m_wfk

FUNCTION

  This module defines the wfk_t object providing a high-level API
  to perform common IO operations on the WFK file produced by ABINIT.
  The API wraps thee different formats/io-libraries:

    1) binary Fortran files with sequential Fortran IO (read, write)
    2) binary Fortran files with MPI-IO primitives (C-steam + Fortran records)
    3) Netcdf files with parallel IO a.k.a HDF5

  and emulate random access when binary Fortran files are used in *read-only* mode.
  See notes below for more info.

COPYRIGHT

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

PARENTS

NOTES

  1) The wfk_t object supports random access also when plain Fortran-IO is used.
     One can easily *read* the block of wavefunctions with a given (kpt,spin)
     by simply passing the appropriate indices (ik_ibz,spin) to the wfk_read_ routines.
     Note however that the same feature is not available in write mode when Fortran IO
     is used. In this case indeed one should access the block of wavefunctions
     according to their (kpt,spin) indices in order to write the correct record markers
     MPI-IO and NETCDF do not have such limitation.

  2) MPI-IO read operations are done with file views even for contiguous data of the same type.
     I found, indeed, that mixing views with explicit offset calls causes
     wrong results unless the file is closed and re-open! Very strange since, according to
     the documentation, the two APIs do not interfere and can be mixed. Calls to MPI_FILE_SEEK
     to reset the pointer to the start of the file do not solve the problem. Don't know if it's
     a feature or a bug (the problem showed up with MPICH2, I haven't tested other MPI libraries)

CHILDREN

SOURCE

 47 #if defined HAVE_CONFIG_H
 48 #include "config.h"
 49 #endif
 50 
 51 #include "abi_common.h"
 52 
 53 MODULE m_wfk
 54 
 55  use defs_basis
 56  use m_abicore
 57  use m_build_info
 58  use m_errors
 59 #ifdef HAVE_MPI2
 60  use mpi
 61 #endif
 62  use m_xmpi
 63  use m_mpiotk
 64  use m_hdr
 65  use m_sort
 66  use m_crystal
 67  use m_crystal_io
 68  use m_pawtab
 69  use m_ebands
 70  use m_pawrhoij
 71  use m_wffile
 72  use m_nctk
 73 #ifdef HAVE_NETCDF
 74  use netcdf
 75 #endif
 76 
 77  use defs_abitypes,  only : hdr_type, dataset_type, MPI_type
 78  use defs_datatypes, only : pseudopotential_type, ebands_t
 79  use defs_wvltypes,  only : wvl_internal_type
 80  use m_geometry,     only : metric
 81  use m_time,         only : cwtime, asctime
 82  use m_fstrings,     only : sjoin, strcat, endswith, itoa
 83  use m_io_tools,     only : get_unit, mvrecord, iomode_from_fname, open_file, close_unit, delete_file, file_exists
 84  use m_numeric_tools,only : mask2blocks
 85  use m_cgtk,         only : cgtk_rotate
 86  use m_fftcore,      only : get_kg, ngfft_seq
 87  use m_distribfft,   only : init_distribfft_seq
 88  use m_mpinfo,       only : destroy_mpi_enreg, initmpi_seq
 89  use m_rwwf,         only : rwwf
 90 
 91  implicit none
 92 
 93  private
 94 
 95 #ifdef HAVE_MPI1
 96  include 'mpif.h'
 97 #endif
 98 
 99  integer,private,parameter :: WFK_NOMODE    = 0
100  integer,private,parameter :: WFK_READMODE  = 1
101  integer,private,parameter :: WFK_WRITEMODE = 2

m_wfk/fill_or_check [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  fill_or_check

FUNCTION

INPUTS

PARENTS

      m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

4925 subroutine fill_or_check(task,Hdr,Kvars,ik_ibz,spin,formeig,kg_k,cg_k,eig_k,occ_k,ierr)
4926 
4927 
4928 !This section has been created automatically by the script Abilint (TD).
4929 !Do not modify the following lines by hand.
4930 #undef ABI_FUNC
4931 #define ABI_FUNC 'fill_or_check'
4932 !End of the abilint section
4933 
4934  implicit none
4935 
4936 !Arguments ------------------------------------
4937 !scalars
4938  integer,intent(in) :: ik_ibz,spin,formeig
4939  integer,intent(out) :: ierr
4940  character(len=*),intent(in) :: task
4941 !arrays
4942  integer,intent(in) :: kg_k(:,:)
4943  real(dp),intent(inout) :: cg_k(:,:),eig_k(:),occ_k(:)
4944  type(hdr_type),intent(in) :: Hdr
4945  type(kvars_t),intent(in) :: Kvars
4946 
4947 !Local variables-------------------------------
4948 !scalars
4949  integer :: nkpt,nsppol,nspinor,nband_k,npw_k,band,ipw,kspad,ii,base,idx,mpw,eigsz
4950  character(len=500) :: msg
4951 !arrays
4952  integer,allocatable :: ref_kg_k(:,:)
4953  real(dp),allocatable :: ref_eig_k(:),ref_occ_k(:),ref_cg_k(:,:)
4954 !************************************************************************
4955 
4956  ierr = 0
4957 
4958  nkpt    = Hdr%nkpt
4959  nsppol  = Hdr%nsppol
4960  nspinor = Hdr%nspinor
4961 
4962  nband_k = Hdr%nband(ik_ibz + (spin-1)*nkpt)
4963  npw_k   = Hdr%npwarr(ik_ibz)
4964 
4965  ABI_MALLOC(ref_kg_k,(3,npw_k))
4966  ABI_MALLOC(ref_eig_k,((2*nband_k)**formeig*nband_k))
4967  ABI_MALLOC(ref_occ_k,(nband_k))
4968  ABI_MALLOC(ref_cg_k,(2,npw_k*nspinor*nband_k))
4969 
4970  ref_kg_k = Kvars%kg_k
4971 
4972  ! Pad values according to (k,s).
4973  kspad = (spin-1)*nkpt + (ik_ibz-1) * nband_k
4974 
4975  if (formeig==0) then
4976    eigsz = nband_k
4977    do band=1,nband_k
4978      ref_eig_k(band) = half * (kspad + band)
4979      ref_occ_k(band) = two  * (kspad + band)
4980    end do
4981  else if (formeig==1) then
4982    eigsz = 2*nband_k**2
4983    base=0
4984    do band=1,nband_k
4985      do ii=1,2*nband_k
4986        idx = base + ii
4987        ref_eig_k(idx) = ii*(kspad + band)
4988      end do
4989      base = base + 2*nband_k
4990    end do
4991  end if
4992 
4993  mpw = npw_k*nspinor*nband_k
4994  do ipw=1,mpw
4995    ref_cg_k(1,ipw) =  ipw + kspad
4996    ref_cg_k(2,ipw) = -ipw + kspad
4997  end do
4998 
4999  SELECT CASE (task)
5000  CASE ("fill")
5001    cg_k(:,1:mpw) = ref_cg_k(:,1:mpw)
5002    if (formeig==0) then
5003      eig_k(1:nband_k) = ref_eig_k
5004      occ_k(1:nband_k) = ref_occ_k
5005    else
5006      eig_k(1:2*nband_k**2) = ref_eig_k
5007    end if
5008 
5009  CASE ("check")
5010 
5011    if (ANY( ABS(cg_k(:,1:mpw) - ref_cg_k) > zero)) then
5012      ierr = ierr + 1
5013      MSG_WARNING("Difference in cg_k")
5014    end if
5015 
5016    if (ANY( ABS(kg_k - ref_kg_k) > zero)) then
5017      ierr = ierr + 2
5018      MSG_WARNING("Difference in kg_k")
5019      !write(std_out,*)"ref_kg_k",ref_kg_k
5020      !write(std_out,*)"kg_k",kg_k
5021    end if
5022 
5023    if (ANY( ABS(eig_k(1:eigsz) - ref_eig_k) > zero)) then
5024      ierr = ierr + 4
5025      MSG_WARNING("Difference in eig_k")
5026      !write(std_out,*)"ref_eig_k",ref_eig_k
5027      !write(std_out,*)"eig_k",eig_k
5028    end if
5029 
5030    if (formeig==0) then
5031      if (ANY( ABS(occ_k(1:nband_k) - ref_occ_k) > zero)) then
5032        ierr = ierr + 8
5033        MSG_WARNING("occ_k")
5034        !write(std_out,*)"ref_occ_k",ref_occ_k
5035        !write(std_out,*)"occ_k",occ_k
5036      end if
5037    end if
5038 
5039    write(msg,"(a,3(i0,2x))")" (ik_ibz, spin, ierr) ",ik_ibz,spin,ierr
5040    if (ierr/=0) then
5041      MSG_WARNING(TRIM(msg)//": FAILED")
5042    else
5043      call wrtout(std_out,TRIM(msg)//": OK","COLL")
5044    end if
5045 
5046  CASE DEFAULT
5047    MSG_ERROR("Wrong task")
5048  END SELECT
5049 
5050  ABI_FREE(ref_kg_k)
5051  ABI_FREE(ref_eig_k)
5052  ABI_FREE(ref_occ_k)
5053  ABI_FREE(ref_cg_k)
5054 
5055 end subroutine fill_or_check

m_wfk/kvars_t [ Types ]

[ Top ] [ m_wfk ] [ Types ]

NAME

FUNCTION

SOURCE

271  type,public :: kvars_t
272    integer,allocatable  :: kg_k(:,:)
273    real(dp),pointer :: occ_k(:)   => null()
274    real(dp),pointer :: eig_k(:)   => null()
275  end type kvars_t
276 
277 CONTAINS

m_wfk/mpio_read_eigocc_k [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  mpio_read_eigocc_k

FUNCTION

  Helper functions to read the eigenvalues and the occupations with MPI-IO

INPUTS

  fh
  offset
  nband_disk
  formeig
  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading.

 OUTPUTS
  buffer(:)
  mpierr=MPI error status.

PARENTS

      m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3784 #ifdef HAVE_MPI_IO
3785 
3786 subroutine mpio_read_eigocc_k(fh,offset,nband_disk,formeig,sc_mode,buffer,mpierr)
3787 
3788 
3789 !This section has been created automatically by the script Abilint (TD).
3790 !Do not modify the following lines by hand.
3791 #undef ABI_FUNC
3792 #define ABI_FUNC 'mpio_read_eigocc_k'
3793 !End of the abilint section
3794 
3795  implicit none
3796 
3797 !Arguments ------------------------------------
3798 !scalars
3799  integer,intent(in) :: fh,nband_disk,formeig,sc_mode
3800  integer(XMPI_OFFSET_KIND),intent(in) :: offset
3801  integer,intent(out) :: mpierr
3802 !arrays
3803  real(dp),pointer :: buffer(:)
3804 
3805 !Local variables-------------------------------
3806 !scalars
3807  integer :: myfh,bufsz,gkk_type,eneocc_type
3808  integer(XMPI_OFFSET_KIND) :: my_offset,my_offpad !,fmarker
3809 !arrays
3810  integer :: sizes(2),subsizes(2),starts(2)
3811 
3812 !************************************************************************
3813 
3814  ! Workaround for XLF
3815  myfh = fh
3816 
3817  SELECT CASE (formeig)
3818  CASE (0)
3819    !
3820    ! Read both eig and occ in buffer
3821    bufsz = 2*nband_disk
3822    my_offset = offset
3823    ABI_MALLOC(buffer, (bufsz))
3824 
3825    call MPI_TYPE_CONTIGUOUS(bufsz, MPI_DOUBLE_PRECISION, eneocc_type, mpierr)
3826    ABI_HANDLE_MPIERR(mpierr)
3827 
3828    call MPI_TYPE_COMMIT(eneocc_type,mpierr)
3829    ABI_HANDLE_MPIERR(mpierr)
3830 
3831    call MPI_FILE_SET_VIEW(myfh, my_offset, MPI_BYTE, eneocc_type, 'native', xmpio_info, mpierr)
3832    ABI_HANDLE_MPIERR(mpierr)
3833 
3834    call MPI_TYPE_FREE(eneocc_type,mpierr)
3835    ABI_HANDLE_MPIERR(mpierr)
3836 
3837    if (sc_mode==xmpio_collective) then
3838      call MPI_FILE_READ_ALL(myfh,buffer,bufsz,MPI_DOUBLE_PRECISION,MPI_STATUS_IGNORE,mpierr)
3839    else if (sc_mode==xmpio_single) then
3840      call MPI_FILE_READ(myfh,buffer,bufsz,MPI_DOUBLE_PRECISION,MPI_STATUS_IGNORE,mpierr)
3841    else
3842      MSG_ERROR("Wrong sc_mode")
3843    end if
3844    ABI_HANDLE_MPIERR(mpierr)
3845 
3846  CASE (1)
3847    ! Read the (nband_k,nband_k) matrix with the (complex) GKK matrix elements.
3848    bufsz    = (nband_disk**2)
3849    sizes    = [nband_disk, nband_disk]
3850    subsizes = [nband_disk, nband_disk]
3851    starts   = [1, 1]
3852 
3853    ABI_MALLOC(buffer, (2*bufsz))
3854 
3855    !my_offset = offset - xmpio_bsize_frm
3856    !call xmpio_read_dp(myfh,my_offset,sc_mode,2*nband_disk,buffer,fmarker,mpierr)
3857    !write(std_out,*)buffer(1:2*nband_disk)
3858    !MSG_ERROR("Done")
3859 
3860    call xmpio_create_fsubarray_2D(sizes,subsizes,starts,MPI_DOUBLE_COMPLEX,gkk_type,my_offpad,mpierr)
3861    ABI_HANDLE_MPIERR(mpierr)
3862 
3863    ! TODO: Rationalize the offsets
3864    my_offset = offset + my_offpad - xmpio_bsize_frm
3865 
3866    call MPI_FILE_SET_VIEW(myfh,my_offset,MPI_BYTE,gkk_type,'native',xmpio_info,mpierr)
3867    ABI_HANDLE_MPIERR(mpierr)
3868 
3869    call MPI_TYPE_FREE(gkk_type, mpierr)
3870    ABI_HANDLE_MPIERR(mpierr)
3871 
3872    if (sc_mode==xmpio_collective) then
3873      call MPI_FILE_READ_ALL(myfh,buffer,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
3874    else if (sc_mode==xmpio_single) then
3875      call MPI_FILE_READ(myfh,buffer,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
3876    else
3877      MSG_ERROR("Wrong sc_mode")
3878    end if
3879    ABI_HANDLE_MPIERR(mpierr)
3880 
3881  CASE DEFAULT
3882    MSG_ERROR("formeig not in [0,1]")
3883  END SELECT
3884 
3885 end subroutine mpio_read_eigocc_k
3886 #endif

m_wfk/mpio_read_kg_k [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  mpio_read_kg_k

FUNCTION

  Helper functions to read the G-vectors with MPI-IO

INPUTS

  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading.

 OUTPUTS
  kg_k=(3,npw_disk) = G-vectors
  mpierr=MPI error status (error check is delegated to the caller)

PARENTS

      m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3603 #ifdef HAVE_MPI_IO
3604 
3605 subroutine mpio_read_kg_k(fh,offset,npw_disk,sc_mode,kg_k,mpierr)
3606 
3607 
3608 !This section has been created automatically by the script Abilint (TD).
3609 !Do not modify the following lines by hand.
3610 #undef ABI_FUNC
3611 #define ABI_FUNC 'mpio_read_kg_k'
3612 !End of the abilint section
3613 
3614  implicit none
3615 
3616 !Arguments ------------------------------------
3617 !scalars
3618  integer,intent(in) :: fh,npw_disk,sc_mode
3619  integer(XMPI_OFFSET_KIND),intent(in) :: offset
3620  integer,intent(out) :: mpierr
3621 !arrays
3622  integer,intent(out) :: kg_k(3,npw_disk)
3623 
3624 !Local variables-------------------------------
3625 !scalars
3626  integer :: kg_k_type,ncount,myfh
3627  integer(XMPI_OFFSET_KIND) :: my_offset
3628 
3629 !************************************************************************
3630 
3631  ! Workarounds for XLF
3632  myfh      = fh
3633  ncount    = 3*npw_disk
3634  my_offset = offset
3635 
3636  call MPI_TYPE_CONTIGUOUS(ncount, MPI_INTEGER, kg_k_type, mpierr)
3637  ABI_HANDLE_MPIERR(mpierr)
3638 
3639  call MPI_TYPE_COMMIT(kg_k_type,mpierr)
3640  ABI_HANDLE_MPIERR(mpierr)
3641 
3642  call MPI_FILE_SET_VIEW(myfh,my_offset,MPI_BYTE,kg_k_type,'native',xmpio_info,mpierr)
3643  ABI_HANDLE_MPIERR(mpierr)
3644 
3645  if (sc_mode==xmpio_collective) then
3646    call MPI_FILE_READ_ALL(myfh,kg_k,ncount,MPI_INTEGER,MPI_STATUS_IGNORE,mpierr)
3647  else if (sc_mode==xmpio_single) then
3648    !call MPI_File_seek(myfh, 0, MPI_SEEK_SET,mpierr)
3649    call MPI_FILE_READ(myfh,kg_k,ncount,MPI_INTEGER, MPI_STATUS_IGNORE,mpierr)
3650  else
3651    MSG_ERROR("Wrong sc_mode")
3652  end if
3653 
3654  ABI_HANDLE_MPIERR(mpierr)
3655 
3656  call MPI_TYPE_FREE(kg_k_type,mpierr)
3657  ABI_HANDLE_MPIERR(mpierr)
3658 
3659 end subroutine mpio_read_kg_k
3660 #endif

m_wfk/mpio_write_eigocc_k [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  mpio_write_eigocc_k

FUNCTION

  Helper functions to write the eigenvalues and the occupations with MPI-IO

INPUTS

  fh
  offset
  nband_disk
  formeig
  sc_mode= MPI-IO option
    xmpio_single     ==> for writing  by current proc.
    xmpio_collective ==> for collective write.

 OUTPUTS
  buffer(:)
  mpierr=MPI error status.

PARENTS

      m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3920 #ifdef HAVE_MPI_IO
3921 
3922 subroutine mpio_write_eigocc_k(fh,offset,nband_disk,formeig,sc_mode,buffer,mpierr)
3923 
3924 
3925 !This section has been created automatically by the script Abilint (TD).
3926 !Do not modify the following lines by hand.
3927 #undef ABI_FUNC
3928 #define ABI_FUNC 'mpio_write_eigocc_k'
3929 !End of the abilint section
3930 
3931  implicit none
3932 
3933 !Arguments ------------------------------------
3934 !scalars
3935  integer,intent(in) :: fh,nband_disk,formeig,sc_mode
3936  integer(XMPI_OFFSET_KIND),intent(in) :: offset
3937  integer,intent(out) :: mpierr
3938 !arrays
3939  real(dp),intent(in) :: buffer(:)
3940 
3941 !Local variables-------------------------------
3942 !scalars
3943  integer :: bufsz,gkk_type,eneocc_type,myfh
3944  integer(XMPI_OFFSET_KIND) :: my_offset,my_offpad
3945 !arrays
3946  integer :: sizes(2),subsizes(2),starts(2)
3947 
3948 !************************************************************************
3949 
3950  ! Workaround for XLF
3951  myfh = fh
3952 
3953  SELECT CASE (formeig)
3954  CASE (0)
3955    !
3956    ! write both eig and occ in buffer
3957    my_offset = offset
3958 
3959    bufsz = 2*nband_disk
3960    ABI_CHECK(SIZE(buffer) >= bufsz, "buffer too small")
3961 
3962    call MPI_TYPE_CONTIGUOUS(bufsz, MPI_DOUBLE_PRECISION, eneocc_type, mpierr)
3963    ABI_HANDLE_MPIERR(mpierr)
3964 
3965    call MPI_TYPE_COMMIT(eneocc_type,mpierr)
3966    ABI_HANDLE_MPIERR(mpierr)
3967 
3968    call MPI_FILE_SET_VIEW(myfh,my_offset,MPI_BYTE,eneocc_type,'native',xmpio_info,mpierr)
3969    ABI_HANDLE_MPIERR(mpierr)
3970 
3971    call MPI_TYPE_FREE(eneocc_type,mpierr)
3972    ABI_HANDLE_MPIERR(mpierr)
3973 
3974    if (sc_mode==xmpio_collective) then
3975      call MPI_FILE_WRITE_ALL(myfh,buffer,bufsz,MPI_DOUBLE_PRECISION,MPI_STATUS_IGNORE,mpierr)
3976    else if (sc_mode==xmpio_single) then
3977      call MPI_FILE_WRITE(myfh,buffer,bufsz,MPI_DOUBLE_PRECISION,MPI_STATUS_IGNORE,mpierr)
3978    else
3979      MSG_ERROR("Wrong sc_mode")
3980    end if
3981    ABI_HANDLE_MPIERR(mpierr)
3982 
3983  CASE (1)
3984    !MSG_ERROR("formeig ==1 with MPI-IO not tested")
3985    ! write the (nband_k,nband_k) matrix with the (complex) GKK matrix elements.
3986    bufsz    = (nband_disk**2)
3987    sizes    = [nband_disk, nband_disk]
3988    subsizes = [nband_disk, nband_disk]
3989    starts   = [1, 1]
3990 
3991    ABI_CHECK(SIZE(buffer) >= bufsz, "buffer too small")
3992 
3993    call xmpio_create_fsubarray_2D(sizes,subsizes,starts,MPI_DOUBLE_COMPLEX,gkk_type,my_offpad,mpierr)
3994    ABI_HANDLE_MPIERR(mpierr)
3995 
3996    ! TODO: Rationalize the offsets
3997    my_offset = offset + my_offpad - xmpio_bsize_frm
3998 
3999    call MPI_FILE_SET_VIEW(myfh,my_offset,MPI_BYTE,gkk_type,'native',xmpio_info,mpierr)
4000    ABI_HANDLE_MPIERR(mpierr)
4001 
4002    call MPI_TYPE_FREE(gkk_type, mpierr)
4003    ABI_HANDLE_MPIERR(mpierr)
4004 
4005    if (sc_mode==xmpio_collective) then
4006      call MPI_FILE_WRITE_ALL(myfh,buffer,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
4007    else if (sc_mode==xmpio_single) then
4008      call MPI_FILE_WRITE(myfh,buffer,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
4009    else
4010      MSG_ERROR("Wrong sc_mode")
4011    end if
4012    ABI_HANDLE_MPIERR(mpierr)
4013 
4014  CASE DEFAULT
4015    MSG_ERROR("formeig not in [0,1]")
4016  END SELECT
4017 
4018 end subroutine mpio_write_eigocc_k

m_wfk/mpio_write_kg_k [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  mpio_write_kg_k

FUNCTION

  Helper function to write the G-vectors with MPI-IO

INPUTS

  sc_mode= MPI-IO option
    xmpio_single     ==> for writing by current proc.
    xmpio_collective ==> for collective write.
  kg_k=(3,npw_disk) = G-vectors

 OUTPUTS
  mpierr=MPI error status (error check is delegated to the caller)

PARENTS

      m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3690 #ifdef HAVE_MPI_IO
3691 
3692 subroutine mpio_write_kg_k(fh,offset,npw_disk,sc_mode,kg_k,mpierr)
3693 
3694 
3695 !This section has been created automatically by the script Abilint (TD).
3696 !Do not modify the following lines by hand.
3697 #undef ABI_FUNC
3698 #define ABI_FUNC 'mpio_write_kg_k'
3699 !End of the abilint section
3700 
3701  implicit none
3702 
3703 !Arguments ------------------------------------
3704 !scalars
3705  integer,intent(in) :: fh,npw_disk,sc_mode
3706  integer(XMPI_OFFSET_KIND),intent(in) :: offset
3707  integer,intent(out) :: mpierr
3708 !arrays
3709  integer,intent(in) :: kg_k(3,npw_disk)
3710 
3711 !Local variables-------------------------------
3712 !scalars
3713  integer :: myfh,kg_k_type,ncount
3714  integer(XMPI_OFFSET_KIND) :: my_offset
3715 
3716 !************************************************************************
3717 
3718  DBG_ENTER("COLL")
3719 
3720  ! Workarounds for XLF
3721  myfh      = fh
3722  ncount    = 3*npw_disk
3723  my_offset = offset
3724 
3725  call MPI_TYPE_CONTIGUOUS(ncount, MPI_INTEGER, kg_k_type, mpierr)
3726  ABI_HANDLE_MPIERR(mpierr)
3727 
3728  call MPI_TYPE_COMMIT(kg_k_type,mpierr)
3729  ABI_HANDLE_MPIERR(mpierr)
3730 
3731  call MPI_FILE_SET_VIEW(myfh,my_offset,MPI_BYTE,kg_k_type,'native',xmpio_info,mpierr)
3732  ABI_HANDLE_MPIERR(mpierr)
3733 
3734  call MPI_TYPE_FREE(kg_k_type,mpierr)
3735  ABI_HANDLE_MPIERR(mpierr)
3736 
3737  if (sc_mode==xmpio_collective) then
3738    call MPI_FILE_WRITE_ALL(myfh,kg_k,ncount,MPI_INTEGER,MPI_STATUS_IGNORE,mpierr)
3739  else if (sc_mode==xmpio_single) then
3740    call MPI_FILE_WRITE(myfh,kg_k,ncount,MPI_INTEGER,MPI_STATUS_IGNORE,mpierr)
3741  else
3742    MSG_ERROR("Wrong sc_mode")
3743  end if
3744 
3745  ABI_HANDLE_MPIERR(mpierr)
3746 
3747  DBG_EXIT("COLL")
3748 
3749 end subroutine mpio_write_kg_k
3750 #endif

m_wfk/wfk_check_wfkfile [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_check_wfkfile

FUNCTION

INPUTS

PARENTS

      ioprof

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

4773 subroutine wfk_check_wfkfile(wfk_fname,Hdr,iomode,method,formeig,Kvars,cwtimes,comm,ierr)
4774 
4775 
4776 !This section has been created automatically by the script Abilint (TD).
4777 !Do not modify the following lines by hand.
4778 #undef ABI_FUNC
4779 #define ABI_FUNC 'wfk_check_wfkfile'
4780 !End of the abilint section
4781 
4782  implicit none
4783 
4784 !Arguments ------------------------------------
4785 !scalars
4786  integer,intent(in) :: iomode,formeig,comm,method
4787  integer,intent(out) :: ierr
4788  character(len=*),intent(in) :: wfk_fname
4789 !arrays
4790  real(dp),intent(out) :: cwtimes(2)
4791  type(hdr_type),intent(in) :: Hdr
4792  type(kvars_t),intent(in) :: Kvars(Hdr%nkpt)
4793 
4794 !Local variables-------------------------------
4795 !scalars
4796  integer :: nkpt,nsppol,nspinor,ik_ibz,spin,funt,nband_k,npw_k,sc_mode
4797  integer :: my_ierr,restart,restartpaw,is,ik,ntests,test,mband
4798  real(dp) :: cpu,wall,gflops
4799  character(len=500) :: msg
4800  type(wfk_t) :: Wfk
4801 !arrays
4802  integer :: nband(Hdr%nkpt,Hdr%nsppol),spins(Hdr%nsppol),kindices(Hdr%nkpt)
4803  integer,allocatable :: kg_k(:,:)
4804  real(dp),allocatable :: cg_k(:,:),eig_k(:),occ_k(:)
4805  logical,allocatable :: bmask(:)
4806 
4807 !************************************************************************
4808 
4809  !write(msg,"(3a,i2)")"Checking file: ",TRIM(wfk_fname),", with iomode = ",iomode
4810  !call wrtout(std_out,msg,"COLL")
4811 
4812  ierr = 0
4813  cwtimes = zero
4814 
4815  nband   = RESHAPE(Hdr%nband, (/Hdr%nkpt,Hdr%nsppol/) )
4816  nkpt    = Hdr%nkpt
4817  nsppol  = Hdr%nsppol
4818  nspinor = Hdr%nspinor
4819 
4820  ! Open the file for writing.
4821  call cwtime(cpu,wall,gflops,"start")
4822  funt = get_unit()
4823 
4824  call wfk_open_read(Wfk,wfk_fname,formeig,iomode,funt,comm)
4825  mband = Wfk%mband
4826 
4827  call cwtime(cpu,wall,gflops,"stop")
4828  cwtimes = cwtimes + (/cpu,wall/)
4829 
4830  ntests = 2
4831 
4832  do test=1,ntests
4833    spins    = (/(spin, spin=1,Hdr%nsppol)/)
4834    kindices = (/(ik_ibz, ik_ibz=1,Hdr%nkpt)/)
4835 
4836    if (test==2) then ! Reverse the indices
4837      spins    = (/(spin, spin=Hdr%nsppol,1,-1)/)
4838      kindices = (/(ik_ibz, ik_ibz=Hdr%nkpt,1,-1)/)
4839    end if
4840    !
4841    do is=1,SIZE(spins)
4842      spin = spins(is)
4843      do ik=1,SIZE(kindices)
4844        ik_ibz = kindices(ik)
4845 
4846        if (Wfk%debug) then
4847          call hdr_check(Wfk%fform,Wfk%fform,Hdr,Wfk%Hdr,"COLL",restart,restartpaw)
4848        end if
4849 
4850        nband_k = nband(ik_ibz,spin)
4851        npw_k   = Hdr%npwarr(ik_ibz)
4852 
4853        ABI_MALLOC(kg_k, (3,npw_k))
4854        ABI_MALLOC(cg_k, (2,npw_k*nspinor*nband_k))
4855        ABI_MALLOC(eig_k, ((2*mband)**Wfk%formeig*mband) )
4856        ABI_MALLOC(occ_k, (mband))
4857        !
4858        !sc_mode = xmpio_collective
4859        sc_mode = xmpio_single
4860        ABI_MALLOC(bmask, (mband))
4861        bmask = .FALSE.
4862        bmask(1:nband_k) = .TRUE.
4863 
4864        call cwtime(cpu,wall,gflops,"start")
4865 
4866        if (method==0) then
4867          call wfk_read_band_block(Wfk,(/1,nband_k/),ik_ibz,spin,sc_mode,kg_k=kg_k,cg_k=cg_k,eig_k=eig_k,occ_k=occ_k)
4868        else if (method==1) then
4869          call wfk_read_bmask(Wfk,bmask,ik_ibz,spin,sc_mode,kg_k=kg_k,cg_k=cg_k,eig_k=eig_k,occ_k=occ_k)
4870        else
4871          MSG_ERROR("Wrong method")
4872        end if
4873 
4874        !call wfk_read_eigk(Wfk,ik_ibz,spin,sc_mode,eig_k)
4875        !write(std_out,*)"eig_k",eig_k
4876 
4877        call cwtime(cpu,wall,gflops,"stop")
4878        cwtimes = cwtimes + (/cpu,wall/)
4879 
4880        ! Check the correctness of the reading.
4881        call fill_or_check("check",Hdr,Kvars(ik_ibz),ik_ibz,spin,formeig,kg_k,cg_k,eig_k,occ_k,my_ierr)
4882 
4883        if (my_ierr/=0) then
4884          write(msg,"(a,i0)")"fill_or_check returned my_ierr: ",my_ierr
4885          ierr = my_ierr
4886          MSG_WARNING(msg)
4887        end if
4888 
4889        ABI_FREE(kg_k)
4890        ABI_FREE(cg_k)
4891        ABI_FREE(eig_k)
4892        ABI_FREE(occ_k)
4893 
4894        ABI_FREE(bmask)
4895      end do
4896    end do
4897    !
4898  end do ! test
4899 
4900  ! Close the file
4901  call wfk_close(Wfk)
4902 
4903 end subroutine wfk_check_wfkfile

m_wfk/wfk_close [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_close

FUNCTION

  Close the wavefunction file handler and release the memory allocated
  Delete the file if `delete` is True. Default: False

PARENTS

      conducti_nc,d2frnl,dfpt_nstdy,dfpt_nstpaw,dfpt_scfcv,fold2Bloch,initwf
      ioprof,m_cut3d,m_iowf,m_wfd,m_wfk,mrggkk,optic

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

671 subroutine wfk_close(Wfk, delete)
672 
673 
674 !This section has been created automatically by the script Abilint (TD).
675 !Do not modify the following lines by hand.
676 #undef ABI_FUNC
677 #define ABI_FUNC 'wfk_close'
678 !End of the abilint section
679 
680  implicit none
681 
682 !Arguments ------------------------------------
683 !scalars
684  type(wfk_t),intent(inout) :: Wfk
685  logical,optional,intent(in) :: delete
686 
687 !Local variables-------------------------------
688 !scalars
689  integer :: ierr
690  character(len=500) :: msg
691 #ifdef HAVE_MPI_IO
692  integer :: mpierr,nfrec
693  integer(XMPI_OFFSET_KIND),allocatable :: bsize_frecords(:)
694 #endif
695 
696 ! *************************************************************************
697 
698  DBG_ENTER("COLL")
699 
700  !@wfk_t
701 
702  ! Close the file only if it was open.
703  if (wfk%rw_mode /= WFK_NOMODE) then
704    Wfk%rw_mode = WFK_NOMODE
705 
706    select case (Wfk%iomode)
707    case (IO_MODE_FORTRAN)
708       ABI_FCLOSE(Wfk%fh, msg)
709 
710 #ifdef HAVE_MPI_IO
711    case (IO_MODE_MPI)
712      call MPI_FILE_CLOSE(Wfk%fh,mpierr)
713      ABI_CHECK_MPI(mpierr,"FILE_CLOSE!")
714 
715      if (Wfk%debug .and. Wfk%my_rank == Wfk%master) then
716        ! Check the fortran records.
717        call MPI_FILE_OPEN(xmpi_comm_self, Wfk%fname, MPI_MODE_RDONLY, xmpio_info, Wfk%fh, mpierr)
718        ABI_CHECK_MPI(mpierr,"MPI_FILE_OPEN")
719        call hdr_bsize_frecords(Wfk%Hdr,Wfk%formeig,nfrec,bsize_frecords)
720        call xmpio_check_frmarkers(Wfk%fh,Wfk%hdr_offset,xmpio_single,nfrec,bsize_frecords,ierr)
721        ABI_CHECK(ierr==0,"xmpio_check_frmarkers returned ierr!=0")
722        ABI_FREE(bsize_frecords)
723        call MPI_FILE_CLOSE(Wfk%fh,mpierr)
724        ABI_CHECK_MPI(mpierr,"FILE_CLOSE!")
725      end if
726 #endif
727 
728 #ifdef HAVE_NETCDF
729    case (IO_MODE_ETSF)
730      NCF_CHECK(nf90_close(wfk%fh))
731 #endif
732 
733    case default
734      MSG_ERROR(sjoin('Wrong/unsupported value of iomode: ', itoa(Wfk%iomode)))
735    end select
736  end if
737 
738  ! Free memory.
739  call hdr_free(Wfk%Hdr)
740 
741  if (allocated(Wfk%nband)) then
742    ABI_FREE(Wfk%nband)
743  end if
744  if (allocated(Wfk%recn_ks)) then
745    ABI_FREE(Wfk%recn_ks)
746  end if
747  if (allocated(Wfk%offset_ks)) then
748    ABI_FREE(Wfk%offset_ks)
749  end if
750 
751  if (present(delete)) then
752    if (delete) call delete_file(wfk%fname, ierr)
753  end if
754 
755  DBG_EXIT("COLL")
756 
757 end subroutine wfk_close

m_wfk/wfk_compare [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_compare

FUNCTION

  Test two wfk_t objects for consistency. Return non-zero value if test fails.

INPUTS

  wfk1, wfk1<type(wfk_t)> = WFK handlers to be compared

OUTPUT

  ierr

PARENTS

CHILDREN

SOURCE

1085 integer function wfk_compare(wfk1, wfk2) result(ierr)
1086 
1087 
1088 !This section has been created automatically by the script Abilint (TD).
1089 !Do not modify the following lines by hand.
1090 #undef ABI_FUNC
1091 #define ABI_FUNC 'wfk_compare'
1092 !End of the abilint section
1093 
1094  implicit none
1095 
1096 !Arguments ------------------------------------
1097 !scalars
1098  type(wfk_t),intent(in) :: wfk1, wfk2
1099 
1100 !Local variables-------------------------------
1101 !scalars
1102  integer :: restart,restartpaw
1103  !character(len=500) :: msg
1104 
1105 !************************************************************************
1106 
1107  ierr = 0
1108 
1109  ! Test basic dimensions
1110  if (wfk1%hdr%nsppol /= wfk2%hdr%nsppol) then
1111    ierr = ierr + 1; MSG_WARNING("Different nsppol")
1112  end if
1113  if (wfk1%hdr%nspinor /= wfk2%hdr%nspinor) then
1114    ierr = ierr + 1; MSG_WARNING("Different nspinor")
1115  end if
1116  if (wfk1%hdr%nspden /= wfk2%hdr%nspden) then
1117    ierr = ierr + 1; MSG_WARNING("Different nspden")
1118  end if
1119  if (wfk1%hdr%nkpt /= wfk2%hdr%nkpt) then
1120    ierr = ierr + 1; MSG_WARNING("Different nkpt")
1121  end if
1122  if (wfk1%formeig /= wfk2%formeig) then
1123    ierr = ierr + 1; MSG_WARNING("Different formeig")
1124  end if
1125  if (wfk1%hdr%usepaw /= wfk2%hdr%usepaw) then
1126    ierr = ierr + 1; MSG_WARNING("Different usepaw")
1127  end if
1128  if (wfk1%hdr%ntypat /= wfk2%hdr%ntypat) then
1129    ierr = ierr + 1; MSG_WARNING("Different ntypat")
1130  end if
1131  if (wfk1%hdr%natom /= wfk2%hdr%natom) then
1132    ierr = ierr + 1; MSG_WARNING("Different natom")
1133  end if
1134  !if (wfk1%hdr%fform /= wfk2%hdr%fform) then
1135  !  ierr = ierr + 1; MSG_WARNING("Different fform")
1136  !end if
1137 
1138  ! Return immediately if important dimensions are not equal.
1139  if (ierr /= 0) return
1140 
1141  ! Test important arrays (rprimd is not tested)
1142  if (any(wfk1%hdr%typat /= wfk2%hdr%typat)) then
1143    ierr = ierr + 1; MSG_WARNING("Different typat")
1144  end if
1145  if (any(wfk1%hdr%npwarr /= wfk2%hdr%npwarr)) then
1146    ierr = ierr + 1; MSG_WARNING("Different npwarr array")
1147  end if
1148  if (any(wfk1%nband /= wfk2%nband)) then
1149    ierr = ierr + 1; MSG_WARNING("Different nband array")
1150  end if
1151  if (any(abs(wfk1%hdr%kptns - wfk2%hdr%kptns) > tol6)) then
1152    ierr = ierr + 1; MSG_WARNING("Different kptns array")
1153  end if
1154 
1155  ! Call hdr_check to get a nice diff of the header but don't check restart and restartpaw.
1156  call hdr_check(wfk1%fform,wfk2%fform,wfk1%hdr,wfk2%hdr,"PERS",restart,restartpaw)
1157 
1158 end function wfk_compare

m_wfk/wfk_create_wfkfile [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_create_wfkfile

FUNCTION

INPUTS

PARENTS

      ioprof

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

4657 subroutine wfk_create_wfkfile(wfk_fname,Hdr,iomode,formeig,Kvars,cwtimes,comm)
4658 
4659 
4660 !This section has been created automatically by the script Abilint (TD).
4661 !Do not modify the following lines by hand.
4662 #undef ABI_FUNC
4663 #define ABI_FUNC 'wfk_create_wfkfile'
4664 !End of the abilint section
4665 
4666  implicit none
4667 
4668 !Arguments ------------------------------------
4669 !scalars
4670  integer,intent(in) :: iomode,formeig,comm
4671  character(len=*),intent(in) :: wfk_fname
4672 !arrays
4673  real(dp),intent(out) :: cwtimes(2)
4674  type(hdr_type),intent(in) :: Hdr
4675  type(kvars_t),target,intent(out) :: Kvars(Hdr%nkpt)
4676 
4677 !Local variables-------------------------------
4678 !scalars
4679  integer :: nkpt,nsppol,nspinor,ierr,sc_mode
4680  integer :: ik_ibz,spin,funt,nband_k,npw_k,istwfk_k
4681  real(dp) :: cpu,wall,gflops,ucvol
4682  type(wfk_t) :: Wfk
4683 !arrays
4684  integer :: nband(Hdr%nkpt,Hdr%nsppol)
4685  integer,pointer :: kg_k(:,:)
4686  real(dp) :: kpoint(3),gmet(3,3),gprimd(3,3),rmet(3,3)
4687  real(dp),allocatable :: cg_k(:,:),eig_k(:),occ_k(:)
4688 
4689 !************************************************************************
4690 
4691  cwtimes = zero
4692 
4693  nband   = RESHAPE(Hdr%nband, (/Hdr%nkpt,Hdr%nsppol/) )
4694  nkpt    = Hdr%nkpt
4695  nsppol  = Hdr%nsppol
4696  nspinor = Hdr%nspinor
4697 
4698  call metric(gmet,gprimd,dev_null,rmet,Hdr%rprimd,ucvol)
4699 
4700  ! Generate the G-vectors from input Hdr%ecut.
4701  do ik_ibz=1,nkpt
4702    kpoint   = Hdr%kptns(:,ik_ibz)
4703    istwfk_k = Hdr%istwfk(ik_ibz)
4704    call get_kg(kpoint,istwfk_k,Hdr%ecut,gmet,npw_k,Kvars(ik_ibz)%kg_k)
4705    ABI_CHECK(npw_k == Hdr%npwarr(ik_ibz),"npw_k != Hdr%npwarr(ik)")
4706  end do
4707  !
4708  ! Open the file for writing.
4709  sc_mode = xmpio_collective
4710 
4711  call cwtime(cpu,wall,gflops,"start")
4712  funt = get_unit()
4713 
4714  call wfk_open_write(Wfk,Hdr,wfk_fname,formeig,iomode,funt,comm,write_frm=.TRUE.)
4715 
4716  call cwtime(cpu,wall,gflops,"stop")
4717  cwtimes = cwtimes + (/cpu,wall/)
4718 
4719  do spin=1,nsppol
4720    do ik_ibz=1,nkpt
4721 
4722      nband_k = nband(ik_ibz,spin)
4723      npw_k   = Hdr%npwarr(ik_ibz)
4724      ABI_MALLOC(cg_k, (2,npw_k*nspinor*nband_k))
4725      ABI_MALLOC(eig_k, ((2*Wfk%mband)**formeig*Wfk%mband) )
4726      ABI_MALLOC(occ_k, (Wfk%mband))
4727 
4728      kg_k => Kvars(ik_ibz)%kg_k
4729      !
4730      ! Fill cg_k, eig_k, occ_k using a deterministic algorithm so that
4731      ! we can check the correctness of the reading.
4732      call fill_or_check("fill",Hdr,Kvars(ik_ibz),ik_ibz,spin,formeig,kg_k,cg_k,eig_k,occ_k,ierr)
4733      ABI_CHECK(ierr==0,"filling")
4734 
4735      call cwtime(cpu,wall,gflops,"start")
4736 
4737      call wfk_write_band_block(Wfk,(/1,nband_k/),ik_ibz,spin,sc_mode,kg_k=kg_k,cg_k=cg_k,eig_k=eig_k,occ_k=occ_k)
4738 
4739      call cwtime(cpu,wall,gflops,"stop")
4740      cwtimes = cwtimes + (/cpu,wall/)
4741 
4742      ABI_FREE(cg_k)
4743      ABI_FREE(eig_k)
4744      ABI_FREE(occ_k)
4745    end do
4746  end do
4747 
4748  ! Close the file
4749  call wfk_close(Wfk)
4750 
4751 end subroutine wfk_create_wfkfile

m_wfk/wfk_diff [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_diff

FUNCTION

  Compare two WFK file for binary equality

INPUTS

PARENTS

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

5077 subroutine wfk_diff(fname1,fname2,formeig,comm,ierr)
5078 
5079 
5080 !This section has been created automatically by the script Abilint (TD).
5081 !Do not modify the following lines by hand.
5082 #undef ABI_FUNC
5083 #define ABI_FUNC 'wfk_diff'
5084 !End of the abilint section
5085 
5086  implicit none
5087 
5088 !Arguments ------------------------------------
5089  integer,intent(in) :: formeig,comm
5090  integer,intent(out) :: ierr
5091  character(len=*),intent(in) :: fname1,fname2
5092 
5093 !Local variables-------------------------------
5094 !scalars
5095  integer,parameter :: master=0
5096  integer :: iomode1,iomode2,ik_ibz,spin,mband,nband_k
5097  integer :: npw_k,mcg,fform1,fform2,sc_mode,my_rank,nproc
5098  character(len=500) :: msg
5099  type(hdr_type) :: Hdr1,Hdr2
5100  type(Wfk_t) :: Wfk1,Wfk2
5101 !arrays
5102  integer,allocatable :: kg1_k(:,:),kg2_k(:,:)
5103  real(dp),allocatable :: eig1_k(:),cg1_k(:,:),occ1_k(:)
5104  real(dp),allocatable :: eig2_k(:),cg2_k(:,:),occ2_k(:)
5105 
5106 ! *************************************************************************
5107 
5108  call wrtout(std_out,ABI_FUNC//": comparing "//TRIM(fname1)//" "//TRIM(fname2),"COLL")
5109 
5110  my_rank = xmpi_comm_rank(comm); nproc   = xmpi_comm_size(comm)
5111  sc_mode = xmpio_collective
5112 
5113  call hdr_read_from_fname(Hdr1,fname1,fform1,comm)
5114  call hdr_read_from_fname(Hdr2,fname2,fform2,comm)
5115 
5116  ABI_CHECK(fform1==fform2,"fform1 != fform2")
5117  ABI_CHECK(Hdr1%nsppol==Hdr2%nsppol,"nsppol1 != nsppol2")
5118  ABI_CHECK(Hdr1%nspinor==Hdr2%nspinor,"nspinor1 != nspinor2")
5119  ABI_CHECK(Hdr1%nkpt==Hdr2%nkpt,"nkpt1 != nkpt2")
5120  !call hdr_check(fform,fform0,hdr1,hdr2,"COLL",restart,restartpaw)
5121 
5122  iomode1 = iomode_from_fname(fname1)
5123  iomode2 = iomode_from_fname(fname1)
5124  ABI_CHECK(iomode1==iomode2,"iomode1 != iomode2")
5125  !iomode1 = IO_MODE_FORTRAN
5126  !iomode2 = IO_MODE_MPI
5127 
5128  call wfk_open_read(Wfk1,fname1,formeig,iomode1,get_unit(),comm)
5129  call wfk_open_read(Wfk2,fname2,formeig,iomode2,get_unit(),comm)
5130 
5131  if (wfk_compare(wfk1, wfk2) /= 0) then
5132    MSG_ERROR("wfk files are not consistent. see above messages")
5133  end if
5134 
5135  mband = Wfk1%mband
5136  ABI_CHECK(mband==Wfk2%mband,"different mband")
5137  ABI_CHECK(all(Wfk1%nband==Wfk2%nband),"different nband")
5138  ABI_CHECK(all(Wfk1%hdr%npwarr==Wfk2%hdr%npwarr),"different npwarr")
5139 
5140  ierr = 0
5141  do spin=1,Wfk1%nsppol
5142    do ik_ibz=1,Wfk1%nkpt
5143      npw_k    = Wfk1%Hdr%npwarr(ik_ibz)
5144      nband_k  = Wfk1%nband(ik_ibz,spin)
5145      ABI_CHECK(npw_k  ==Wfk2%Hdr%npwarr(ik_ibz),"different npw_k")
5146      ABI_CHECK(nband_k==Wfk2%nband(ik_ibz,spin),"different nband_k")
5147 
5148      mcg = npw_k*Hdr1%nspinor*nband_k
5149 
5150      ABI_MALLOC(eig1_k,((2*mband)**formeig*mband))
5151      ABI_MALLOC(occ1_k,(mband))
5152      ABI_MALLOC(kg1_k,(3,npw_k))
5153      ABI_STAT_MALLOC(cg1_k,(2,mcg), ierr)
5154      ABI_CHECK(ierr==0, "out of memory in cg1_k")
5155 
5156      ABI_MALLOC(eig2_k,((2*mband)**formeig*mband))
5157      ABI_MALLOC(occ2_k,(mband))
5158      ABI_MALLOC(kg2_k,(3,npw_k))
5159      ABI_STAT_MALLOC(cg2_k,(2,mcg), ierr)
5160      ABI_CHECK(ierr==0, "out of memory in cg2_k")
5161 
5162      ! Read the block of bands for this (k,s).
5163      call wfk_read_band_block(Wfk1,(/1,nband_k/),ik_ibz,spin,sc_mode,kg_k=kg1_k,eig_k=eig1_k,occ_k=occ1_k) !, cg_k=cg1_k,
5164      call wfk_read_band_block(Wfk2,(/1,nband_k/),ik_ibz,spin,sc_mode,kg_k=kg2_k,eig_k=eig2_k,occ_k=occ2_k) !, cg_k=cg2_k,
5165 
5166      if (ANY( ABS(kg1_k - kg2_k) > zero)) then
5167        ierr = ierr + 2
5168        MSG_WARNING("Difference in kg_k")
5169        !write(std_out,*)"kg1_k",kg1_k
5170        !write(std_out,*)"kg2_k",kg2_k
5171      end if
5172 
5173      if (ANY( ABS(eig1_k - eig2_k) > zero)) then
5174        ierr = ierr + 4
5175        MSG_WARNING("Difference in eig_k")
5176        !write(std_out,*)"eig1_k",eig1_k
5177        !write(std_out,*)"eig2_k",eig2_k
5178      end if
5179 
5180      if (formeig==0) then
5181        if (ANY( ABS(occ1_k - occ2_k) > zero)) then
5182          ierr = ierr + 8
5183          MSG_WARNING("occ_k")
5184          write(std_out,*)"occ1_k",occ1_k
5185          write(std_out,*)"occ2_k",occ2_k
5186        end if
5187      end if
5188 
5189      if (ANY( ABS(cg1_k - cg2_k) > zero)) then
5190        ierr = ierr + 1
5191        MSG_WARNING("Difference in cg_k")
5192      end if
5193 
5194      write(msg,"(a,3(i0,2x))")" (ik_ibz, spin, ierr) ",ik_ibz,spin,ierr
5195      if (ierr/=0) then
5196        MSG_WARNING(TRIM(msg)//": FAILED")
5197      else
5198        call wrtout(std_out,TRIM(msg)//": OK","COLL")
5199      end if
5200 
5201      ABI_FREE(eig1_k)
5202      ABI_FREE(occ1_k)
5203      ABI_FREE(kg1_k)
5204      ABI_FREE(cg1_k)
5205 
5206      ABI_FREE(eig2_k)
5207      ABI_FREE(occ2_k)
5208      ABI_FREE(kg2_k)
5209      ABI_FREE(cg2_k)
5210    end do !ik_ibz
5211  end do !spin
5212 
5213  call wfk_close(Wfk1)
5214  call wfk_close(Wfk2)
5215 
5216  call hdr_free(Hdr1)
5217  call hdr_free(Hdr2)
5218 
5219 end subroutine wfk_diff

m_wfk/wfk_findk [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_findk

FUNCTION

  Find the index of the k-point in the WKF file. umklapp vectors are not allowed.
  Return -1 if not found.

INPUTS

  wfk<type(wfk_t)> = WFK handler initialized and set in read mode
  kpt(3)=k-point in reduced coordinates.
  [ktol]=Optional tolerance for k-point comparison.
         For each reduced direction the absolute difference between the coordinates must be less that ktol

PARENTS

CHILDREN

SOURCE

924 integer pure function wfk_findk(wfk, kpt, ktol) result(ikpt)
925 
926 
927 !This section has been created automatically by the script Abilint (TD).
928 !Do not modify the following lines by hand.
929 #undef ABI_FUNC
930 #define ABI_FUNC 'wfk_findk'
931 !End of the abilint section
932 
933  implicit none
934 
935 !Arguments ------------------------------------
936 !scalars
937  real(dp),optional,intent(in) :: ktol
938  type(wfk_t),intent(in) :: wfk
939 !arrays
940  real(dp),intent(in) :: kpt(3)
941 
942 !Local variables-------------------------------
943 !scalars
944  integer :: ik
945  real(dp) :: my_ktol
946 
947 ! *************************************************************************
948 
949  my_ktol = 0.0001_dp; if (present(ktol)) my_ktol = ktol
950 
951  ikpt = -1
952  do ik=1,wfk%hdr%nkpt
953    if (all(abs(wfk%hdr%kptns(:, ik) - kpt) < my_ktol)) then
954      ikpt = ik; exit
955    end if
956  end do
957 
958 end function wfk_findk

m_wfk/wfk_nc2fort [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_nc2fort

FUNCTION

  Convert a netcdf WFK file (nc_path) to a Fortran WFK file (fort_path).

PARENTS

      ioprof

NOTES

  - This routine should be called by a single processor.
  - Only GS WFK files are supported (formeig==0)

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

4366 subroutine wfk_nc2fort(nc_path, fort_path)
4367 
4368 
4369 !This section has been created automatically by the script Abilint (TD).
4370 !Do not modify the following lines by hand.
4371 #undef ABI_FUNC
4372 #define ABI_FUNC 'wfk_nc2fort'
4373 !End of the abilint section
4374 
4375  implicit none
4376 
4377 !Arguments ------------------------------------
4378 !scalars
4379  character(len=*),intent(in) :: nc_path,fort_path
4380 
4381 !Local variables-------------------------------
4382 !scalars
4383  integer :: ik,spin,mband,mpw,nband_k
4384  type(wfk_t) :: in_wfk,out_wfk
4385 !arrays
4386  integer,parameter :: formeig0=0
4387  integer,allocatable :: kg_k(:,:)
4388  real(dp),allocatable :: cg_k(:,:),eig_k(:),occ_k(:)
4389 
4390 ! *************************************************************************
4391 
4392  call wrtout(std_out, sjoin("Converting:", nc_path, "to", fort_path))
4393 
4394  ! Open input file, extract dimensions and allocate workspace arrays.
4395  call wfk_open_read(in_wfk,nc_path,formeig0,IO_MODE_ETSF,get_unit(),xmpi_comm_self)
4396 
4397  mpw = maxval(in_wfk%hdr%npwarr); mband = in_wfk%mband
4398  ABI_MALLOC(kg_k, (3, mpw))
4399  ABI_MALLOC(cg_k, (2, mpw*in_wfk%nspinor*mband))
4400  ABI_MALLOC(eig_k, ((2*mband)**in_wfk%formeig*mband) )
4401  ABI_MALLOC(occ_k, (mband))
4402 
4403  ! Open output file.
4404  call wfk_open_write(out_wfk,in_wfk%hdr,fort_path,formeig0,IO_MODE_FORTRAN,get_unit(),xmpi_comm_self)
4405 
4406  do spin=1,in_wfk%nsppol
4407    do ik=1,in_wfk%nkpt
4408      nband_k = in_wfk%nband(ik,spin)
4409 
4410      call wfk_read_band_block(in_wfk,[1,nband_k],ik,spin,xmpio_single,&
4411        kg_k=kg_k,cg_k=cg_k,eig_k=eig_k,occ_k=occ_k)
4412 
4413      call wfk_write_band_block(out_wfk,[1,nband_k],ik,spin,xmpio_single,&
4414        kg_k=kg_k,cg_k=cg_k,eig_k=eig_k,occ_k=occ_k)
4415    end do
4416  end do
4417 
4418  ABI_FREE(kg_k)
4419  ABI_FREE(cg_k)
4420  ABI_FREE(eig_k)
4421  ABI_FREE(occ_k)
4422 
4423  call wfk_close(in_wfk)
4424  call wfk_close(out_wfk)
4425 
4426 end subroutine wfk_nc2fort

m_wfk/wfk_ncdef_dims_vars [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_ncdef_dims_vars

FUNCTION

  Write the heeder, fform as well as the etsf-io dimensions and variables.

INPUTS

  ncid=netcdf file handler.
  hdr<hdr_tyep>=Abinit header
  fform=File type format of the header
  [write_hdr]=True if the header should be written (default)
  [iskss]=True if this is a KSS file (activate kdependent=No)

PARENTS

      m_io_kss,m_iowf,m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

 986 subroutine wfk_ncdef_dims_vars(ncid, hdr, fform, write_hdr, iskss)
 987 
 988 
 989 !This section has been created automatically by the script Abilint (TD).
 990 !Do not modify the following lines by hand.
 991 #undef ABI_FUNC
 992 #define ABI_FUNC 'wfk_ncdef_dims_vars'
 993 !End of the abilint section
 994 
 995  implicit none
 996 
 997 !Arguments ------------------------------------
 998 !scalars
 999  integer,intent(in) :: ncid,fform
1000  type(hdr_type),intent(in) :: hdr
1001  logical,optional,intent(in) :: write_hdr,iskss
1002 
1003 !Local variables-------------------------------
1004 !scalars
1005 #ifdef HAVE_NETCDF
1006  character(len=500) :: title,history
1007  logical :: do_write_hdr,my_iskss
1008  integer :: ivar,mpw,ncerr
1009 ! *************************************************************************
1010 
1011  do_write_hdr = .True.; if (present(write_hdr)) do_write_hdr = write_hdr
1012  my_iskss = .False.; if (present(iskss)) my_iskss = iskss
1013  if (do_write_hdr) then
1014    NCF_CHECK(hdr_ncwrite(hdr, ncid, fform, nc_define=.True.))
1015  end if
1016 
1017  ! Add the etsf header.
1018  title = sjoin("WFK file generated by Abinit, version: ",  abinit_version)
1019  if (my_iskss) title = sjoin("KSS file generated by Abinit, version: ",  abinit_version)
1020  history = sjoin("Generated on: ", asctime())
1021  NCF_CHECK(nctk_add_etsf_header(ncid, title=title, history=history))
1022 
1023  mpw = MAXVAL(hdr%npwarr)
1024  ncerr = nctk_def_dims(ncid, [&
1025    nctkdim_t("real_or_complex_coefficients", 2), nctkdim_t("max_number_of_coefficients", mpw)])
1026  NCF_CHECK(ncerr)
1027 
1028  ! Define kg_k
1029  ncerr = nctk_def_arrays(ncid, [&
1030    nctkarr_t("reduced_coordinates_of_plane_waves", "int", &
1031      "number_of_reduced_dimensions, max_number_of_coefficients, number_of_kpoints")&
1032  ])
1033  NCF_CHECK(ncerr)
1034 
1035  NCF_CHECK(nf90_inq_varid(ncid, "reduced_coordinates_of_plane_waves", ivar))
1036  if (my_iskss) then
1037    NCF_CHECK(nf90_put_att(ncid, ivar, "k_dependent", "no"))
1038  else
1039    NCF_CHECK(nf90_put_att(ncid, ivar, "k_dependent", "yes"))
1040  end if
1041 
1042  ncerr = nctk_def_arrays(ncid, [&
1043    nctkarr_t("eigenvalues", "dp", "max_number_of_states, number_of_kpoints, number_of_spins") &
1044  ])
1045  NCF_CHECK(ncerr)
1046  NCF_CHECK(nctk_set_atomic_units(ncid, "eigenvalues"))
1047 
1048  ncerr = nctk_def_arrays(ncid, nctkarr_t("coefficients_of_wavefunctions", "dp", &
1049    "real_or_complex_coefficients, max_number_of_coefficients, number_of_spinor_components, &
1050 &max_number_of_states, number_of_kpoints, number_of_spins"))
1051  NCF_CHECK(ncerr)
1052 
1053  !NF90_DEF_VAR_FILL(INTEGER NCID, INTEGER VARID, INTEGER NO_FILL, FILL_VALUE)
1054  !NCF_CHECK(nf90_inq_varid(ncid, "coefficients_of_wavefunctions", ivar))
1055  !NCF_CHECK(nf90_def_var_fill(ncid, ivar, 0, -one))
1056 
1057 #else
1058  MSG_ERROR("netcdf not available")
1059 #endif
1060 
1061 end subroutine wfk_ncdef_dims_vars

m_wfk/wfk_open_read [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_open_read

FUNCTION

  Open the WFK file in read mode.

INPUTS

  fname = Name of the file
  formeig = 0 for GS wavefunctions, 1 for RF wavefunctions.
  iomode = access mode
  funt = Fortran unit numer. Only used if iomode == IO_MODE_FORTRAN
  comm = MPI communicator (used for collective parallel IO)

OUTPUT

  Wfk<type(wfk_t)> = WFK handler initialized and set in read mode
  [Hdr_out]=Copy of the abinit header

PARENTS

      conducti_nc,d2frnl,dfpt_looppert,dfpt_nstdy,dfpt_nstpaw,fold2Bloch
      initwf,ioprof,m_cut3d,m_wfd,m_wfk,mrggkk,optic

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

310 subroutine wfk_open_read(Wfk,fname,formeig,iomode,funt,comm,Hdr_out)
311 
312 
313 !This section has been created automatically by the script Abilint (TD).
314 !Do not modify the following lines by hand.
315 #undef ABI_FUNC
316 #define ABI_FUNC 'wfk_open_read'
317 !End of the abilint section
318 
319  implicit none
320 
321 !Arguments ------------------------------------
322 !scalars
323  integer,intent(in) :: iomode,comm,formeig,funt
324  character(len=*),intent(in) :: fname
325  type(wfk_t),intent(inout) :: Wfk
326  type(hdr_type),optional,intent(inout) :: Hdr_out  ! should be intent(out), but psc miscompiles the call!
327 
328 !Local variables-------------------------------
329 !scalars
330  integer :: ierr,mpierr
331  character(len=500) :: msg
332 #ifdef HAVE_MPI_IO
333  integer :: fform !,ncerr
334 #endif
335 
336 !************************************************************************
337 
338  DBG_ENTER("COLL")
339 
340  !Initialize the mandatory data of the Wfk datastructure
341  !@wfk_t
342  Wfk%rw_mode     = WFK_READMODE
343  Wfk%chunk_bsize = WFK_CHUNK_BSIZE
344 
345  Wfk%fname     = fname
346  Wfk%formeig   = formeig
347  Wfk%iomode    = iomode; if (endswith(fname, ".nc")) wfk%iomode = IO_MODE_ETSF
348  ! This is to test the different versions.
349  !wfk%iomode    = IO_MODE_MPI
350  !if (.not. endswith(fname, ".nc") .and. xmpi_comm_size == 1) wfk%iomode == IO_MODE_FORTRAN
351 
352  Wfk%comm      = comm
353  Wfk%master    = 0
354  Wfk%my_rank   = xmpi_comm_rank(comm)
355  Wfk%nproc     = xmpi_comm_size(comm)
356 
357  ! Reads fform and the Header.
358  call hdr_read_from_fname(Wfk%Hdr,fname,Wfk%fform,comm)
359  ABI_CHECK(Wfk%fform/=0,"fform ==0")
360 
361  if (Wfk%debug) call hdr_echo(Wfk%Hdr,Wfk%fform,4,unit=std_out)
362 
363  ! Copy the header if required.
364  if (present(Hdr_out)) call hdr_copy(Wfk%Hdr,Hdr_out)
365 
366  ! Useful dimensions
367  Wfk%mband   = MAXVAL(Wfk%Hdr%nband)
368  Wfk%nkpt    = Wfk%Hdr%nkpt
369  Wfk%nsppol  = Wfk%Hdr%nsppol
370  Wfk%nspinor = Wfk%Hdr%nspinor
371 
372  ABI_MALLOC(Wfk%nband, (Wfk%nkpt,Wfk%nsppol))
373  Wfk%nband = RESHAPE(Wfk%Hdr%nband, (/Wfk%nkpt,Wfk%nsppol/))
374 
375  ierr=0
376  select case (wfk%iomode)
377  case (IO_MODE_FORTRAN)
378    ! All processors see a local Fortran binary file.
379    ! Each node opens the file, skip the header and set f90_fptr.
380    Wfk%fh = funt
381    if (open_file(Wfk%fname,msg,unit=Wfk%fh,form="unformatted", status="old", action="read") /= 0) then
382      MSG_ERROR(msg)
383    end if
384 
385    ! Precompute number of records for Fortran IO.
386    call wfk_compute_offsets(Wfk)
387 
388    call hdr_skip(Wfk%fh,ierr)
389    ABI_CHECK(ierr==0, "hdr_skip returned ierr! /= 0")
390    Wfk%f90_fptr = [1,1,REC_NPW]
391 
392 #ifdef HAVE_MPI_IO
393  case (IO_MODE_MPI)
394    call MPI_FILE_OPEN(Wfk%comm, Wfk%fname, MPI_MODE_RDONLY, xmpio_info, Wfk%fh, mpierr)
395    ABI_CHECK_MPI(mpierr,"MPI_FILE_OPEN")
396    !call MPI_FILE_SET_VIEW(Wfk%fh,origin,MPI_BYTE,MPI_BYTE,'native',xmpio_info,mpierr)
397 
398    call hdr_mpio_skip(Wfk%fh,fform,Wfk%hdr_offset)
399 
400    ! Precompute offsets for MPI-IO access
401    if (Wfk%hdr_offset > 0) then
402      call wfk_compute_offsets(Wfk)
403    else
404      MSG_ERROR("hdr_offset <=0")
405    end if
406 #endif
407 
408 #ifdef HAVE_NETCDF
409  case (IO_MODE_ETSF)
410    NCF_CHECK(nctk_open_read(wfk%fh, wfk%fname, wfk%comm))
411 #endif
412 
413  case default
414    MSG_ERROR(sjoin('Wrong or unsupported iomode:', itoa(wfk%iomode)))
415  end select
416 
417  DBG_EXIT("COLL")
418 
419 end subroutine wfk_open_read

m_wfk/wfk_open_write [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_open_write

FUNCTION

  Open the WFK file in write mode.

INPUTS

  fname = Name of the file
  formeig = 0 for GS wavefunctions, 1 for RF wavefunctions.
  iomode = access mode
  funt = Fortran unit numer for  Only used if iomode == IO_MODE_FORTRAN
  comm = MPI communicator (used for MPI-IO)
  [write_hdr]=True if the header should be written (default)
  [write_frm]=True if the fortran record markers should be written (default). Only if Fortran binary file.

OUTPUT

  Wfk<type(wfk_t)> = WFK handler initialized and set in read mode

PARENTS

      m_iowf,m_wfd,m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

452 subroutine wfk_open_write(Wfk,Hdr,fname,formeig,iomode,funt,comm,write_hdr,write_frm)
453 
454 
455 !This section has been created automatically by the script Abilint (TD).
456 !Do not modify the following lines by hand.
457 #undef ABI_FUNC
458 #define ABI_FUNC 'wfk_open_write'
459 !End of the abilint section
460 
461  implicit none
462 
463 !Arguments ------------------------------------
464 !scalars
465  integer,intent(in) :: iomode,comm,formeig,funt
466  character(len=*),intent(in) :: fname
467  logical,optional,intent(in) :: write_hdr,write_frm
468  type(hdr_type),intent(in) :: Hdr
469  type(wfk_t),intent(out) :: Wfk
470 
471 !Local variables-------------------------------
472 !scalars
473  integer :: mpierr,ierr
474  real(dp) :: cpu,wall,gflops
475  logical :: do_write_frm,do_write_hdr
476  character(len=500) :: msg
477 #ifdef HAVE_MPI_IO
478  integer :: fform,nfrec,sc_mode
479  integer(XMPI_OFFSET_KIND) :: offset
480  integer(XMPI_OFFSET_KIND),allocatable :: bsize_frecords(:)
481 #endif
482 #ifdef HAVE_NETCDF
483  integer :: ncerr
484 #endif
485 
486 !************************************************************************
487 
488  DBG_ENTER("COLL")
489 
490  do_write_hdr = .TRUE.; if (present(write_hdr)) do_write_hdr = write_hdr
491  do_write_frm = .TRUE.; if (present(write_frm)) do_write_frm = write_frm
492 
493  !Initialize mandatory data of the Wfk datastructure
494  !@wfk_t
495  Wfk%rw_mode     = WFK_WRITEMODE
496  Wfk%chunk_bsize = WFK_CHUNK_BSIZE
497 
498  Wfk%fname     = fname
499  Wfk%formeig   = formeig
500  Wfk%iomode    = iomode; if (endswith(fname, ".nc")) wfk%iomode = IO_MODE_ETSF
501  ! This is to test the different versions.
502  !wfk%iomode   = IO_MODE_MPI
503  !if (.not. endswith(fname, ".nc") .and. xmpi_comm_size == 1) wfk%iomode == IO_MODE_FORTRAN
504 
505  Wfk%comm      = comm
506  Wfk%master    = 0
507  Wfk%my_rank   = xmpi_comm_rank(comm)
508  Wfk%nproc     = xmpi_comm_size(comm)
509  Wfk%fform     = 2
510 
511  ! Copy the header
512  call hdr_copy(Hdr,Wfk%Hdr)
513 
514  ! Master writes fform and the Header (write it afterwards if IO_MODE_ETSF)
515  if (Wfk%my_rank==Wfk%master .and. do_write_hdr .and. iomode /= IO_MODE_ETSF) then
516    call hdr_write_to_fname(Wfk%Hdr,Wfk%fname,Wfk%fform)
517    if (Wfk%debug) call hdr_echo(Wfk%Hdr,Wfk%fform,4,unit=std_out)
518  end if
519  call xmpi_barrier(Wfk%comm)
520 
521  ! Useful dimensions
522  Wfk%mband   = MAXVAL(Wfk%Hdr%nband)
523  Wfk%nkpt    = Wfk%Hdr%nkpt
524  Wfk%nsppol  = Wfk%Hdr%nsppol
525  Wfk%nspinor = Wfk%Hdr%nspinor
526 
527  ABI_MALLOC(Wfk%nband, (Wfk%nkpt,Wfk%nsppol))
528  Wfk%nband = RESHAPE(Wfk%Hdr%nband, (/Wfk%nkpt,Wfk%nsppol/))
529 
530  ierr=0
531 
532  select case (wfk%iomode)
533  case (IO_MODE_FORTRAN)
534    ABI_CHECK(wfk%nproc == 1, "Cannot use Fortran-IO to write WFK file with nprocs > 1")
535    Wfk%fh = funt
536    if (open_file(Wfk%fname,msg,unit=Wfk%fh,form="unformatted", status="unknown", action="readwrite") /= 0) then
537      MSG_ERROR(msg)
538    end if
539 
540    ! Precompute number of records for Fortran IO.
541    call wfk_compute_offsets(Wfk)
542 
543    call hdr_skip(Wfk%fh,ierr)
544    Wfk%f90_fptr = (/1,1,REC_NPW/)
545 
546 #ifdef HAVE_MPI_IO
547  case (IO_MODE_MPI)
548 
549    call cwtime(cpu,wall,gflops,"start")
550 
551    ! FIXME: mode flags should be rationalized
552    !call MPI_FILE_OPEN(Wfk%comm, Wfk%fname, MPI_MODE_CREATE + MPI_MODE_WRONLY, xmpio_info, Wfk%fh, mpierr)
553    !call MPI_FILE_OPEN(Wfk%comm, Wfk%fname, MPI_MODE_CREATE + MPI_MODE_RDWR, xmpio_info, Wfk%fh, mpierr)
554    call MPI_FILE_OPEN(Wfk%comm, Wfk%fname,  MPI_MODE_RDWR, xmpio_info, Wfk%fh, mpierr)
555    ABI_CHECK_MPI(mpierr,"MPI_FILE_OPEN")
556 
557    !call MPI_FILE_SET_VIEW(Wfk%fh,origin,MPI_BYTE,MPI_BYTE,'native',xmpio_info,mpierr)
558 
559    ! TODO
560    !%% call MPI_File_set_size(Wfk%fh, MPI_Offset size, mpierr)
561    !ABI_CHECK_MPI(mpierr,"MPI_FILE_SET_SIZE")
562 
563    call hdr_mpio_skip(Wfk%fh,fform,Wfk%hdr_offset)
564    ABI_CHECK(fform == Wfk%fform,"fform != Wfk%fform")
565    !call hdr_echo(wfk%Hdr, wfk%fform, 4, unit=std_out)
566 
567    ! Precompute offsets for MPI-IO access
568    if (Wfk%hdr_offset > 0) then
569      call wfk_compute_offsets(Wfk)
570    else
571      MSG_ERROR("hdr_offset <=0")
572    end if
573 
574    call cwtime(cpu,wall,gflops,"stop")
575    write(msg,"(2(a,f8.2))")" FILE_OPEN cpu: ",cpu,", wall:",wall
576    call wrtout(std_out,msg,"PERS")
577 
578    ! Write Fortran record markers.
579    !if (.FALSE.) then
580    if (do_write_frm) then
581      call cwtime(cpu,wall,gflops,"start")
582      call hdr_bsize_frecords(Wfk%Hdr,Wfk%formeig,nfrec,bsize_frecords)
583 
584      sc_mode = xmpio_collective
585      offset = Wfk%hdr_offset
586 
587      if (sc_mode == xmpio_collective) then
588        call xmpio_write_frmarkers(Wfk%fh,offset,sc_mode,nfrec,bsize_frecords,ierr)
589      else
590        ierr = 0
591        if (Wfk%my_rank == Wfk%master) then
592          call xmpio_write_frmarkers(Wfk%fh,offset,xmpio_single,nfrec,bsize_frecords,ierr)
593        end if
594      end if
595      ABI_CHECK(ierr==0,"xmpio_write_frmarkers returned ierr!=0")
596 
597      !call MPI_FILE_SYNC(Wfk%fh,mpierr)
598      !ABI_CHECK_MPI(mpierr,"FILE_SYNC")
599 
600      if (Wfk%debug) then
601        call xmpio_check_frmarkers(Wfk%fh,offset,sc_mode,nfrec,bsize_frecords,ierr)
602        ABI_CHECK(ierr==0,"xmpio_check_frmarkers returned ierr!=0")
603      end if
604 
605      ABI_FREE(bsize_frecords)
606 
607      call cwtime(cpu,wall,gflops,"stop")
608      write(msg,"(2(a,f8.2))")" write_frmarkers cpu: ",cpu,", wall:",wall
609      call wrtout(std_out,msg,"PERS")
610    end if
611 #endif
612 
613 #ifdef HAVE_NETCDF
614  CASE (IO_MODE_ETSF)
615    !NCF_CHECK(nctk_open_modify(wfk%fh, wfk%fname, wfk%comm))
616 
617    if (nctk_has_mpiio) then
618 #ifdef HAVE_NETCDF_MPI
619      ncerr = nf90_create(wfk%fname, cmode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_write), &
620          comm=wfk%comm, info=xmpio_info, ncid=wfk%fh)
621 
622      NCF_CHECK_MSG(ncerr, sjoin("nf90_create: ", wfk%fname))
623 #else
624      MSG_ERROR("You should not be here")
625 #endif
626    else
627      if (wfk%nproc > 1) then
628        MSG_ERROR("Your netcdf library does not support MPI-IO. Cannot write WFK file with nprocs > 1")
629      end if
630 
631      ncerr = nf90_create(wfk%fname, nf90_write, wfk%fh)
632      NCF_CHECK_MSG(ncerr, sjoin("nf90_create: ", wfk%fname))
633    end if
634 
635    call wfk_ncdef_dims_vars(wfk%fh, hdr, wfk%fform, write_hdr=.True.)
636    NCF_CHECK(nctk_def_basedims(wfk%fh))
637 
638    ! Switch to data mode.
639    NCF_CHECK(nctk_set_datamode(wfk%fh))
640 #endif
641 
642  case default
643    MSG_ERROR(sjoin('Wrong/unsupported iomode: ', itoa(wfk%iomode)))
644  end select
645 
646  DBG_EXIT("COLL")
647 
648 end subroutine wfk_open_write

m_wfk/wfk_print [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_print

FUNCTION

  Print information on the object.

INPUTS

  wfk<type(wfk_t)> = WFK handler
  [header]=String to be printed as header for additional info.
  [unit]=Unit number for output. Defaults to std_out
  [prtvol]=Verbosity level

PARENTS

      ioprof

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

784 subroutine wfk_print(wfk,unit,header,prtvol)
785 
786 
787 !This section has been created automatically by the script Abilint (TD).
788 !Do not modify the following lines by hand.
789 #undef ABI_FUNC
790 #define ABI_FUNC 'wfk_print'
791 !End of the abilint section
792 
793  implicit none
794 
795 !Arguments ------------------------------------
796 !scalars
797  type(wfk_t),intent(inout) :: wfk
798  integer,optional,intent(in) :: unit,prtvol
799  character(len=*),optional,intent(in) :: header
800 
801 !Local variables-------------------------------
802  integer,parameter :: rdwr4=4
803  integer :: my_unt,my_prtvol
804  character(len=500) :: msg
805 
806 ! *************************************************************************
807 
808  my_unt = std_out; if (present(unit)) my_unt = unit
809  my_prtvol = 0; if (present(prtvol)) my_prtvol = prtvol
810 
811  msg=' ==== Info on the wfk_t object ==== '
812  if (present(header)) msg=' ==== '//trim(adjustl(header))//' ==== '
813  call wrtout(my_unt,msg)
814  call wrtout(my_unt, sjoin(" iomode = ",itoa(wfk%iomode)))
815 
816  call hdr_echo(wfk%hdr, wfk%fform, rdwr4 ,unit=my_unt)
817 
818 end subroutine wfk_print

m_wfk/wfk_prof [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_prof

FUNCTION

  Profiling tool for IO routines

INPUTS

  wfk_fname=Filename
  formeig=0 for GS file, 1 for DFPT file
  nband=Number of bands to read.
  comm=MPI communicator

PARENTS

      ioprof

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

4453 subroutine wfk_prof(wfk_fname, formeig, nband, comm)
4454 
4455 
4456 !This section has been created automatically by the script Abilint (TD).
4457 !Do not modify the following lines by hand.
4458 #undef ABI_FUNC
4459 #define ABI_FUNC 'wfk_prof'
4460 !End of the abilint section
4461 
4462  implicit none
4463 
4464 !Arguments ------------------------------------
4465  integer,intent(in) :: nband,formeig,comm
4466  character(len=*),intent(in) :: wfk_fname
4467 
4468 !Local variables-------------------------------
4469 !scalars
4470  integer,parameter :: rdwr1=1,master=0,optkg1=1,option1=1,tim_rwwf0=0,icg0=0,headform0=0
4471  integer :: iomode,wfk_unt,ik_ibz,spin,ierr,ii,option,mband
4472  integer :: npw_disk,nband_disk,mcg,fform,nband_read,sc_mode,my_rank,nproc
4473  real(dp) :: cpu,wall,gflops
4474  character(len=500) :: msg
4475  type(hdr_type) :: Hdr
4476  type(Wfk_t) :: Wfk
4477  type(wffile_type) :: wff
4478  type(MPI_type) :: MPI_enreg_seq
4479 !arrays
4480  !integer,parameter :: io_modes(2) = (/IO_MODE_FORTRAN, IO_MODE_MPI/)
4481  integer,parameter :: io_modes(1) = (/IO_MODE_MPI/)
4482  integer :: ngfft(18)
4483  logical,allocatable :: my_bmask(:)
4484  integer,allocatable :: kg_k(:,:)
4485  real(dp),allocatable :: eig_k(:),cg_k(:,:),occ_k(:)
4486 
4487 ! *************************************************************************
4488 
4489  my_rank = xmpi_comm_rank(comm); nproc = xmpi_comm_size(comm)
4490  sc_mode = xmpio_collective
4491 
4492  call hdr_read_from_fname(hdr,wfk_fname,fform,comm)
4493 
4494  ! nband_read is the max number of bands we can read from this file.
4495  nband_read = nband
4496  if (nband_read <= 0) then
4497    nband_read = minval(hdr%nband)
4498    call wrtout(std_out, sjoin("nband == 0 --> Setting nband_read to:",itoa(nband_read)))
4499  end if
4500  if (nband_read > minval(hdr%nband)) then
4501    nband_read = minval(hdr%nband)
4502    call wrtout(std_out, sjoin("nband > hdr%nband --> Setting nband_read to:",itoa(nband_read)))
4503  end if
4504 
4505  wfk_unt = get_unit()
4506 
4507  do ii=1,SIZE(io_modes)
4508    iomode = io_modes(ii)
4509    !do option=1,3
4510    do option=1,3,2
4511      write(std_out,*)"iomode, option",iomode,option
4512      call cwtime(cpu,wall,gflops,"start")
4513 
4514      select case (option)
4515      case (1)
4516        call wfk_open_read(Wfk,wfk_fname,formeig,iomode,wfk_unt,comm)
4517 
4518        do spin=1,Hdr%nsppol
4519          do ik_ibz=1,Hdr%nkpt
4520            npw_disk   = Hdr%npwarr(ik_ibz)
4521            nband_disk = Wfk%nband(ik_ibz,spin)
4522 
4523            mcg = npw_disk*Hdr%nspinor*nband_read
4524 
4525            ABI_MALLOC(eig_k,((2*Wfk%mband)**formeig*Wfk%mband))
4526            ABI_MALLOC(occ_k,(Wfk%mband))
4527 
4528            ABI_MALLOC(kg_k,(3,npw_disk))
4529            ABI_STAT_MALLOC(cg_k,(2,mcg), ierr)
4530            ABI_CHECK(ierr==0, "out of memory in cg_k")
4531 
4532            ! Read the block of bands for this (k,s).
4533            call wfk_read_band_block(Wfk,[1,nband_read],ik_ibz,spin,xmpio_collective,&
4534              kg_k=kg_k,cg_k=cg_k,eig_k=eig_k,occ_k=occ_k)
4535 
4536            ABI_FREE(eig_k)
4537            ABI_FREE(occ_k)
4538            ABI_FREE(kg_k)
4539            ABI_FREE(cg_k)
4540          end do !ik_ibz
4541        end do !spin
4542 
4543        call wfk_close(Wfk)
4544 
4545      case (2)
4546        call wfk_open_read(Wfk,wfk_fname,formeig,iomode,wfk_unt,comm)
4547 
4548        do spin=1,Hdr%nsppol
4549          do ik_ibz=1,Hdr%nkpt
4550            npw_disk   = Hdr%npwarr(ik_ibz)
4551            nband_disk = Hdr%nband(ik_ibz+(spin-1)*Hdr%nkpt)
4552 
4553            ABI_MALLOC(my_bmask,(MAXVAL(Hdr%nband)))
4554            my_bmask=.False.; my_bmask(1:nband_read) = .True.
4555 
4556            ABI_MALLOC(eig_k,((2*nband_disk)**formeig*nband_disk))
4557            ABI_MALLOC(kg_k,(3,npw_disk))
4558            ABI_MALLOC(occ_k,(nband_disk))
4559 
4560            mcg = npw_disk*Hdr%nspinor*COUNT(my_bmask)
4561            ABI_STAT_MALLOC(cg_k,(2,mcg), ierr)
4562            ABI_CHECK(ierr==0, "out of memory in cg_k")
4563 
4564            call wfk_read_bmask(Wfk,my_bmask,ik_ibz,spin,sc_mode,kg_k=kg_k,cg_k=cg_k,eig_k=eig_k,occ_k=occ_k)
4565            !call wfk_read_band_block(Wfk,(/1,nband_read/),ik_ibz,spin,sc_mode,kg_k=kg_k,cg_k=cg_k,eig_k=eig_k,occ_k=occ_k)
4566 
4567            ABI_FREE(my_bmask)
4568            ABI_FREE(eig_k)
4569            ABI_FREE(occ_k)
4570            ABI_FREE(kg_k)
4571            ABI_FREE(cg_k)
4572          end do !ik_ibz
4573        end do !spin
4574 
4575        call wfk_close(Wfk)
4576 
4577      case (3)
4578        !Fake MPI_type for the sequential part.
4579        ngfft(1:6) = (/12,12,12,13,13,13/)
4580        call initmpi_seq(MPI_enreg_seq)
4581        call init_distribfft_seq(MPI_enreg_seq%distribfft,'c',ngfft(2),ngfft(3),'all')
4582        call init_distribfft_seq(MPI_enreg_seq%distribfft,'f',ngfft(2),ngfft(3),'all')
4583 
4584        call WffOpen(iomode,comm,wfk_fname,ierr,wff,master,my_rank,wfk_unt) !,spaceComm_mpiio) ! optional argument
4585        ABI_CHECK(ierr==0,"ierr!=0")
4586 
4587        call hdr_free(Hdr)
4588        call hdr_io(fform,Hdr,1,wff)
4589        call WffKg(wff,optkg1)
4590 
4591        do spin=1,Hdr%nsppol
4592          do ik_ibz=1,Hdr%nkpt
4593 
4594            npw_disk   = Hdr%npwarr(ik_ibz)
4595            nband_disk = Hdr%nband(ik_ibz+(spin-1)*Hdr%nkpt)
4596 
4597            mband = MAXVAL(Hdr%nband)
4598            mcg = npw_disk*Hdr%nspinor*nband_read
4599 
4600            ABI_MALLOC(eig_k,((2*mband)**formeig*mband))
4601            ABI_MALLOC(occ_k,(mband))
4602 
4603            ABI_MALLOC(kg_k,(3,optkg1*npw_disk))
4604            ABI_STAT_MALLOC(cg_k,(2,mcg), ierr)
4605            ABI_CHECK(ierr==0, "out of memory in cg_k")
4606            !
4607            ! Read the block of bands for this (k,s).
4608            call rwwf(cg_k,eig_k,formeig,headform0,icg0,ik_ibz,spin,kg_k,mband,mcg,MPI_enreg_seq,nband_read,&
4609              nband_disk,npw_disk,Hdr%nspinor,occ_k,option1,optkg1,tim_rwwf0,Wff)
4610 
4611            ABI_FREE(eig_k)
4612            ABI_FREE(occ_k)
4613            ABI_FREE(kg_k)
4614            ABI_FREE(cg_k)
4615 
4616          end do !ik_ibz
4617        end do !spin
4618 
4619        call WffClose(wff,ierr)
4620        call destroy_mpi_enreg(MPI_enreg_seq)
4621 
4622      case default
4623        MSG_ERROR("Wrong method")
4624      end select
4625 
4626      call cwtime(cpu,wall,gflops,"stop")
4627      write(msg,'(3(a,i2),2(a,f8.2))')&
4628        " iomode: ",iomode,", nproc: ",nproc,", option: ",option,", cpu: ",cpu,", wall:",wall
4629      call wrtout(std_out,msg,"COLL")
4630    end do
4631  end do
4632 
4633  call hdr_free(Hdr)
4634 
4635 end subroutine wfk_prof

m_wfk/wfk_read_band_block [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_read_band_block

FUNCTION

  Read a block of contiguous bands at a given (k-point, spin)

INPUTS

  Wfk<type(wfk_t)>= WFK file handler object.
  band_block(2)=Initial and final band index.
  ik_ibz=Index of the k-point in the IBZ.
  spin=Spin index
  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading in wfk%comm (use it wisely!)

 OUTPUTS
  [kg_k=(:,:)] = G-vectors
  [eig_k(:)] = Eigenvectors
  [cg_k(:,:)] = Fourier coefficients

NOTES

  The output arrays eig_k and occ_k contain the *full* set of eigenvalues and occupation
  factors stored in the file and are dimensioned with wfk%mband.

PARENTS

      fold2Bloch,initwf,m_cut3d,m_wfd,m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

1196 subroutine wfk_read_band_block(Wfk,band_block,ik_ibz,spin,sc_mode,kg_k,cg_k,eig_k,occ_k)
1197 
1198 
1199 !This section has been created automatically by the script Abilint (TD).
1200 !Do not modify the following lines by hand.
1201 #undef ABI_FUNC
1202 #define ABI_FUNC 'wfk_read_band_block'
1203 !End of the abilint section
1204 
1205  implicit none
1206 
1207 !Arguments ------------------------------------
1208 !scalars
1209  integer,intent(in) :: ik_ibz,spin,sc_mode
1210  type(wfk_t),intent(inout) :: Wfk
1211 !arrays
1212  integer,intent(in) :: band_block(2)
1213  integer,intent(out), DEV_CONTARRD  optional :: kg_k(:,:) !(3,npw_k)
1214  real(dp),intent(out), DEV_CONTARRD optional :: cg_k(:,:) !(2,npw_k*nspinor*nband)
1215  real(dp),intent(inout),optional :: eig_k((2*Wfk%mband)**Wfk%formeig*Wfk%mband)
1216  real(dp),intent(out),optional :: occ_k(Wfk%mband)
1217 
1218 !Local variables-------------------------------
1219 !scalars
1220  integer :: ierr,npw_disk,nspinor_disk,nband_disk,band
1221  integer :: ipw,my_bcount,npwso,npw_tot,nb_block,base
1222  integer :: npw_read,nspinor_read,nband_read
1223  character(len=500) :: msg,errmsg
1224 !arrays
1225  real(dp),ABI_CONTIGUOUS pointer :: tmp_eigk(:),tmp_occk(:)
1226 #ifdef HAVE_MPI_IO
1227  integer :: mpierr,bufsz,gkk_type,cgblock_type
1228  integer(XMPI_OFFSET_KIND) :: my_offset,my_offpad
1229  integer :: sizes(2),subsizes(2),starts(2),types(2)
1230 #endif
1231 #ifdef HAVE_NETCDF
1232  integer :: kg_varid,eig_varid,occ_varid,cg_varid,ncerr,h1_varid,idx,ib1,ib2
1233  real(dp),allocatable :: h1mat(:,:,:)
1234 #endif
1235 
1236 !************************************************************************
1237 
1238  DBG_ENTER("COLL")
1239 
1240  ABI_CHECK(Wfk%rw_mode==WFK_READMODE, "Wfk must be in READMODE")
1241 
1242  if (wfk_validate_ks(wfk, ik_ibz, spin) /= 0) then
1243    MSG_ERROR("Wrong (ik_ibz, spin) args, Aborting now")
1244  end if
1245 
1246  ! Look before you leap.
1247  npw_disk     = Wfk%Hdr%npwarr(ik_ibz)
1248  nspinor_disk = Wfk%nspinor
1249  nband_disk   = Wfk%nband(ik_ibz,spin)
1250  nb_block     = (band_block(2) - band_block(1) + 1)
1251  ABI_CHECK(nb_block>0,"nband <=0")
1252  npw_tot      = npw_disk * nspinor_disk * nb_block
1253 
1254  if (present(kg_k)) then
1255    ABI_CHECK(SIZE(kg_k,DIM=2) >= npw_disk,"kg_k too small")
1256  end if
1257 
1258  if (present(cg_k)) then
1259    ABI_CHECK(SIZE(cg_k, DIM=2) >= npw_tot,"cg_k too small")
1260  end if
1261 
1262  if (present(eig_k)) then
1263    if (Wfk%formeig==0) then
1264       ABI_CHECK(SIZE(eig_k) >= nband_disk, "GS eig_k too small")
1265    else if (Wfk%formeig==1) then
1266       ABI_CHECK(SIZE(eig_k) >= 2*nband_disk**2, "DFPT eig_k too small")
1267    else
1268      MSG_ERROR("formeig != [0,1]")
1269    end if
1270  end if
1271 
1272  if (present(occ_k)) then
1273    ABI_CHECK(SIZE(occ_k) >= nband_disk, "GS eig_k too small")
1274    if (Wfk%formeig==1) then
1275      MSG_ERROR("occ_k cannot be used when formeig ==1")
1276    end if
1277  end if
1278 
1279  select case (Wfk%iomode)
1280  case (IO_MODE_FORTRAN)
1281 
1282    ! Rewind the file to have the correct (k,s) block (if needed)
1283    call wfk_seek(Wfk,ik_ibz,spin)
1284    !
1285    ! Read the first record: npw, nspinor, nband_disk
1286    read(Wfk%fh, err=10, iomsg=errmsg) npw_read, nspinor_read, nband_read
1287 
1288    if (any( [npw_read, nspinor_read, nband_read] /= [npw_disk, nspinor_disk, nband_disk])) then
1289      write(msg,"(a,6(i0,2x))")"Mismatch between (npw, nspinor, nband) read from WFK and those found in HDR ",&
1290 &      npw_read, nspinor_read, nband_read, npw_disk, nspinor_disk, nband_disk
1291      MSG_ERROR(msg)
1292    end if
1293 
1294    ! The second record: (k+G) vectors
1295    if (present(kg_k)) then
1296      read(Wfk%fh, err=10, iomsg=errmsg) kg_k(1:3,1:npw_disk)
1297    else
1298      read(Wfk%fh, err=10, iomsg=errmsg) ! kg_k(1:3,1:npw_disk)
1299    end if
1300 
1301    select case (Wfk%formeig)
1302    case (0)
1303      ! The third record: eigenvalues and occupation factors.
1304      ! write(unitwf) (eigen(iband),iband=1,nband_disk),(occ(iband),iband=1,nband_disk)
1305      if (present(eig_k) .or. present(occ_k)) then
1306 
1307        ABI_MALLOC(tmp_eigk, (nband_disk))
1308        ABI_MALLOC(tmp_occk, (nband_disk))
1309 
1310        read(Wfk%fh, err=10, iomsg=errmsg) tmp_eigk, tmp_occk
1311 
1312        if (present(eig_k)) eig_k = tmp_eigk
1313        if (present(occ_k)) occ_k = tmp_occk
1314 
1315        ABI_FREE(tmp_eigk)
1316        ABI_FREE(tmp_occk)
1317 
1318      else
1319        read(Wfk%fh, err=10, iomsg=errmsg) ! eig_k(1:nband_disk), occ_k(1:nband_k)
1320      end if
1321 
1322      ! Fourth record with the wave-functions.
1323      if (present(cg_k)) then
1324        npwso = npw_disk*nspinor_disk
1325        my_bcount = 0
1326        do band=1,nband_disk
1327          if (band >= band_block(1) .and. band <= band_block(2)) then
1328            ipw = my_bcount * npwso
1329            my_bcount = my_bcount + 1
1330            read(Wfk%fh, err=10, iomsg=errmsg) cg_k(1:2,ipw+1:ipw+npwso)
1331          else
1332            read(Wfk%fh, err=10, iomsg=errmsg) ! cg_k(1:2,ipw+1:ipw+npwso)
1333          end if
1334        end do
1335 
1336      else
1337        do band=1,nband_disk
1338          read(Wfk%fh, err=10, iomsg=errmsg) ! cg_k(1:2,ipw+1:ipw+npwso)
1339        end do
1340      end if
1341 
1342    case (1)
1343      npwso = npw_disk*nspinor_disk
1344      my_bcount = 0
1345      do band=1,nband_disk
1346 
1347        if (present(eig_k)) then
1348          ! Read column matrix of size (2*nband_k)
1349          base = 2*(band-1)*nband_disk
1350          read(Wfk%fh, err=10, iomsg=errmsg) eig_k(base+1:base+2*nband_disk)
1351        else
1352          read(Wfk%fh, err=10, iomsg=errmsg ) ! eig_k(2*nband_disk)
1353        end if
1354 
1355        if (present(cg_k) .and. (band >= band_block(1) .and. band <= band_block(2)) ) then
1356          ipw = my_bcount * npwso
1357          my_bcount = my_bcount + 1
1358          read(Wfk%fh, err=10, iomsg=errmsg) cg_k(1:2,ipw+1:ipw+npwso)
1359        else
1360          read(Wfk%fh, err=10, iomsg=errmsg) ! cg_k(1:2,ipw+1:ipw+npwso)
1361        end if
1362      end do
1363 
1364    case default
1365      MSG_ERROR("formeig != [0,1]")
1366    end select
1367 
1368    ! Reached the end of the (k,s) block. Update f90_fptr
1369    if (ik_ibz < Wfk%nkpt) then
1370      Wfk%f90_fptr = [ik_ibz+1,spin,REC_NPW]
1371    else
1372      ABI_CHECK(ik_ibz == wfk%nkpt, "ik_ibz != nkpt")
1373      if (spin==Wfk%nsppol) then
1374        Wfk%f90_fptr = FPTR_EOF ! EOF condition
1375      else
1376        Wfk%f90_fptr = [1,spin+1,REC_NPW]
1377      end if
1378    end if
1379 
1380 #ifdef HAVE_MPI_IO
1381  case (IO_MODE_MPI)
1382    if (present(kg_k)) then
1383      my_offset = Wfk%offset_ks(ik_ibz,spin,REC_KG) + xmpio_bsize_frm
1384 
1385      call mpio_read_kg_k(Wfk%fh,my_offset,npw_disk,sc_mode,kg_k,mpierr)
1386      ABI_CHECK_MPI(mpierr,"reading kg")
1387    end if
1388 
1389    ! formeig=0 =>  Read both eig and occ in tmp_eigk
1390    ! formeig=1 =>  Read (nband_k,nband_k) matrix of complex numbers.
1391    select case (Wfk%formeig)
1392    case (0)
1393      if (present(eig_k) .or. present(occ_k)) then
1394        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_EIG) + xmpio_bsize_frm
1395 
1396        call mpio_read_eigocc_k(Wfk%fh,my_offset,nband_disk,Wfk%formeig,sc_mode,tmp_eigk,mpierr)
1397        ABI_CHECK_MPI(mpierr,"reading eigocc")
1398 
1399        if (present(eig_k)) eig_k(1:nband_disk) = tmp_eigk(1:nband_disk)
1400        if (present(occ_k)) occ_k(1:nband_disk) = tmp_eigk(nband_disk+1:)
1401 
1402        ABI_FREE(tmp_eigk)
1403      end if
1404 
1405      if (present(cg_k)) then
1406        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_CG)
1407        sizes     = [npw_disk*nspinor_disk, nband_disk]
1408        subsizes  = [npw_disk*nspinor_disk, band_block(2)-band_block(1)+1]
1409        bufsz     = 2 * npw_disk * nspinor_disk * nb_block
1410        starts    = [1, band_block(1)]
1411 
1412        call mpiotk_read_fsuba_dp2D(Wfk%fh,my_offset,sizes,subsizes,starts,bufsz,cg_k,&
1413          wfk%chunk_bsize,sc_mode,Wfk%comm,ierr)
1414        ABI_CHECK(ierr==0,"Fortran record too big")
1415      end if
1416 
1417    case (1)
1418      if (present(eig_k)) then
1419        sizes = [nband_disk, npw_disk*nspinor_disk]
1420        types = [MPI_DOUBLE_COMPLEX, MPI_DOUBLE_COMPLEX]
1421 
1422        call xmpio_create_fstripes(nband_disk,sizes,types,gkk_type,my_offpad,mpierr)
1423        ABI_CHECK_MPI(mpierr,"xmpio_create_fstripes")
1424 
1425        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_EIG) + my_offpad
1426 
1427        call MPI_FILE_SET_VIEW(Wfk%fh,my_offset,MPI_BYTE,gkk_type,'native',xmpio_info,mpierr)
1428        ABI_CHECK_MPI(mpierr,"")
1429        call MPI_TYPE_FREE(gkk_type,mpierr)
1430        ABI_CHECK_MPI(mpierr,"")
1431 
1432        bufsz = nband_disk**2
1433 
1434        if (sc_mode==xmpio_collective) then
1435          call MPI_FILE_READ_ALL(Wfk%fh,eig_k,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1436        else if (sc_mode==xmpio_single) then
1437          call MPI_FILE_READ(Wfk%fh,eig_k,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1438        else
1439          MSG_ERROR("Wrong sc_mode")
1440        end if
1441        ABI_CHECK_MPI(mpierr,"FILE_READ")
1442      end if
1443 
1444      if (present(cg_k)) then
1445        types = [MPI_DOUBLE_COMPLEX, MPI_DOUBLE_COMPLEX]
1446        sizes = [npw_disk*nspinor_disk, nband_disk]
1447 
1448        call xmpio_create_fstripes(nb_block,sizes,types,cgblock_type,my_offpad,mpierr)
1449        ABI_CHECK_MPI(mpierr,"xmpio_create_fstripes")
1450 
1451        ! Increment my_offset to account for previous eigen and cg records if band_block(1) != 1
1452        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_CG) + my_offpad
1453        my_offset = my_offset + (band_block(1) - 1) * ( &
1454           (2*npw_disk*wfk%nspinor*xmpi_bsize_dp + 2*xmpio_bsize_frm) + &
1455           (2*nband_disk*xmpi_bsize_dp + 2*xmpio_bsize_frm) )
1456 
1457        call MPI_FILE_SET_VIEW(Wfk%fh,my_offset,MPI_BYTE,cgblock_type,'native',xmpio_info,mpierr)
1458        ABI_CHECK_MPI(mpierr,"SET_VIEW")
1459 
1460        call MPI_TYPE_FREE(cgblock_type,mpierr)
1461        ABI_CHECK_MPI(mpierr,"")
1462 
1463        bufsz = npw_disk * nspinor_disk * nb_block
1464 
1465        if (sc_mode==xmpio_collective) then
1466          call MPI_FILE_READ_ALL(Wfk%fh,cg_k,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1467        else if (sc_mode==xmpio_single) then
1468          call MPI_FILE_READ(Wfk%fh,cg_k,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1469        else
1470          MSG_ERROR("Wrong sc_mode")
1471        end if
1472        ABI_CHECK_MPI(mpierr,"FILE_READ")
1473      end if
1474 
1475    case default
1476      MSG_ERROR("formeig != [0,1]")
1477    end select
1478 #endif
1479 
1480 #ifdef HAVE_NETCDF
1481  case (IO_MODE_ETSF)
1482    if (present(kg_k)) then
1483      ! Read the reduced_coordinates_of_plane_waves for this k point.
1484      NCF_CHECK(nf90_inq_varid(wfk%fh, "reduced_coordinates_of_plane_waves", kg_varid))
1485      if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
1486        NCF_CHECK(nctk_set_collective(wfk%fh, kg_varid))
1487      end if
1488      ncerr = nf90_get_var(wfk%fh, kg_varid, kg_k, start=[1,1,ik_ibz], count=[3,npw_disk,1])
1489      NCF_CHECK(ncerr)
1490    end if
1491 
1492    if (Wfk%formeig==0) then
1493      ! Read eigenvalues and occupations.
1494      if (present(eig_k)) then
1495        NCF_CHECK(nf90_inq_varid(wfk%fh, "eigenvalues", eig_varid))
1496        if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
1497          NCF_CHECK(nctk_set_collective(wfk%fh, eig_varid))
1498        end if
1499        ncerr = nf90_get_var(wfk%fh, eig_varid, eig_k, start=[1,ik_ibz,spin], count=[nband_disk,1,1])
1500        NCF_CHECK(ncerr)
1501      end if
1502 
1503      if (present(occ_k)) then
1504        NCF_CHECK(nf90_inq_varid(wfk%fh, "occupations", occ_varid))
1505        if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
1506          NCF_CHECK(nctk_set_collective(wfk%fh, occ_varid))
1507        end if
1508        ncerr = nf90_get_var(wfk%fh, occ_varid, occ_k, start=[1,ik_ibz,spin], count=[nband_disk,1,1])
1509        NCF_CHECK_MSG(ncerr, "getting occ_k")
1510      end if
1511 
1512    else ! formeig == 1
1513 
1514      if (present(eig_k)) then
1515        ! Read h1 matrix elements. The netcdf array has shape:
1516        ! [complex, max_number_of_states, max_number_of_states, number_of_kpoints, number_of_spins]
1517        ABI_MALLOC(h1mat, (2, wfk%mband, wfk%mband))
1518        NCF_CHECK(nf90_inq_varid(wfk%fh, "h1_matrix_elements", h1_varid))
1519        if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
1520          NCF_CHECK(nctk_set_collective(wfk%fh, h1_varid))
1521        end if
1522        ncerr = nf90_get_var(wfk%fh, h1_varid, h1mat, start=[1,1,1,ik_ibz,spin], &
1523          count=[2, wfk%mband, wfk%mband, 1, 1])
1524        NCF_CHECK_MSG(ncerr, "getting h1mat_k")
1525 
1526        ! For legacy reasons, I have to pack nband_k**2 elements in the first positions in eig_k
1527        ! This is important only if nband(:) depends on k i.e. mband != nband_k
1528        idx=1
1529        do ib2=1,nband_disk
1530          do ib1=1,nband_disk
1531             eig_k(idx:idx+1) = h1mat(:2,ib1,ib2)
1532             idx=idx+2
1533           end do
1534        end do
1535        ABI_FREE(h1mat)
1536      end if
1537    end if
1538 
1539    if (present(cg_k)) then
1540      ! Read the nb_block bands starting from band_block(1)
1541      ! The coefficients_of_wavefunctions on file have shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
1542      NCF_CHECK(nf90_inq_varid(wfk%fh, "coefficients_of_wavefunctions", cg_varid))
1543      if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
1544        NCF_CHECK(nctk_set_collective(wfk%fh, cg_varid))
1545      end if
1546 
1547      ncerr = nf90_get_var(wfk%fh, cg_varid, cg_k, start=[1,1,1,band_block(1),ik_ibz,spin], &
1548        count=[2,npw_disk,wfk%nspinor,nb_block,1,1])
1549      NCF_CHECK_MSG(ncerr, "getting cg_k")
1550   end if
1551 #endif
1552 
1553  case default
1554    MSG_ERROR(sjoin('Wrong/unsupported iomode: ', itoa(Wfk%iomode)))
1555  end select
1556 
1557  DBG_EXIT("COLL")
1558 
1559  return
1560 
1561  ! Handle Fortran IO error
1562 10 continue
1563  MSG_ERROR(errmsg)
1564 
1565 end subroutine wfk_read_band_block

m_wfk/wfk_read_bks [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_read_bks

FUNCTION

  Read the wavefunction and, optionally, the *DFPT* matrix elements
  for a given (band, k-point, spin).

INPUTS

  Wfk<type(wfk_t)>=WFK file handler.
  band=Band index
  ik_ibz=Index of the k-point in the IBZ.
  spin=Spin index
  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading.

 OUTPUTS
  cg_bks(2,npw_k*nspinor) = Fourier coefficients of the wavefunction
  [eig1_bks(2*wfk%mband)] = Matrix elements of the DFPT H1 Hamiltonian at the specified (k, spin).

PARENTS

      d2frnl,dfpt_nstpaw,dfpt_nstwf,dfpt_vtowfk,rf2_init

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

1600 subroutine wfk_read_bks(wfk, band, ik_ibz, spin, sc_mode, cg_bks, eig1_bks)
1601 
1602 
1603 !This section has been created automatically by the script Abilint (TD).
1604 !Do not modify the following lines by hand.
1605 #undef ABI_FUNC
1606 #define ABI_FUNC 'wfk_read_bks'
1607 !End of the abilint section
1608 
1609  implicit none
1610 
1611 !Arguments ------------------------------------
1612 !scalars
1613  integer,intent(in) :: band,ik_ibz,spin,sc_mode
1614  type(wfk_t),intent(inout) :: wfk
1615 !arrays
1616  real(dp),DEV_CONTARRD intent(out) :: cg_bks(:,:)
1617  real(dp),optional,intent(inout) :: eig1_bks(2*wfk%mband)
1618 
1619 !Local variables-------------------------------
1620 !scalars
1621  integer :: start,ib1,nspinor_disk,npw_disk,nband_disk !ierr,
1622 #ifdef HAVE_MPI_IO
1623  integer :: mpierr,bufsz,gkk_type !,cg_type
1624  integer(XMPI_OFFSET_KIND) :: my_offset,my_offpad
1625  integer :: sizes(2),types(2)
1626 #endif
1627 #ifdef HAVE_NETCDF
1628  integer :: h1_varid,cg_varid,ncerr
1629 #endif
1630  character(len=500) :: errmsg
1631 !arrays
1632  real(dp),allocatable :: all_eigk(:)
1633 
1634 !************************************************************************
1635 
1636  if (wfk_validate_ks(wfk, ik_ibz, spin, band=band) /= 0) then
1637    MSG_ERROR("Wrong (ik_ibz, spin, band) args, Aborting now")
1638  end if
1639  npw_disk = wfk%Hdr%npwarr(ik_ibz)
1640  nband_disk = wfk%nband(ik_ibz, spin)
1641  nspinor_disk = wfk%nspinor
1642  !cg_bks = one; if (present(eig1_bks)) eig1_bks = zero ; return
1643 
1644  if (.not. present(eig1_bks)) then
1645    call wfk_read_band_block(wfk, [band, band], ik_ibz, spin, sc_mode, cg_k=cg_bks)
1646    return
1647 
1648  else
1649    ABI_CHECK(wfk%formeig==1, "formeig must be 1 if eig1_bks is present")
1650    ABI_CHECK(size(cg_bks, dim=2) >= npw_disk*wfk%nspinor,"cg_bks too small")
1651    ABI_CHECK(size(eig1_bks) >= 2*nband_disk, "eig1_bks too small")
1652 
1653  if (.False.) then
1654  !if (.True.) then
1655    ! Due to the API of wfk_read_band_block, we have to read the full set of eigenvalues
1656    ! and then extract the relevant band
1657    ! TODO: Should write another routine to avoid doing that.
1658    ABI_MALLOC(all_eigk, (2*wfk%mband**2))
1659    call wfk_read_band_block(wfk, [band, band], ik_ibz, spin, sc_mode, cg_k=cg_bks, eig_k=all_eigk)
1660 
1661    ! Find the index of the slice.
1662    ! Remember that data in all_eigk does not have a constant stride if nband_disk != mband
1663    start = (band-1)*2*nband_disk
1664 
1665    !write(std_out,*)size(eig1_bks), nband_disk
1666    eig1_bks(1:2*nband_disk) = all_eigk(start+1:start+2*nband_disk)
1667    ABI_FREE(all_eigk)
1668 
1669  else
1670    ! Improved version
1671    select case (wfk%iomode)
1672    case (IO_MODE_FORTRAN)
1673      ! This code is not optimal because I'm rewinding the (k, s) block
1674      ! at each call but I'm not gonna spend time on this because I should
1675      ! refactor a lot of stuff. Use MPI-IO or HDF5!
1676      call wfk_seek(wfk,ik_ibz,spin)
1677 
1678      ! Read the first record: npw, nspinor, nband_disk
1679      read(Wfk%fh, err=10, iomsg=errmsg) !npw_read, nspinor_read, nband_read
1680      ! The second record: (k+G) vectors
1681      read(wfk%fh, err=10, iomsg=errmsg) ! kg_k(1:3,1:npw_disk)
1682 
1683      do ib1=1,nband_disk
1684        if (ib1 /= band) then
1685          read(Wfk%fh, err=10, iomsg=errmsg) ! eig_k(base+1:base+2*nband_disk)
1686          read(Wfk%fh, err=10, iomsg=errmsg) ! cg_k(1:2,ipw+1:ipw+npwso)
1687        else
1688          if (present(eig1_bks)) then
1689            read(wfk%fh, err=10, iomsg=errmsg) eig1_bks(1:2*nband_disk)
1690          else
1691            read(wfk%fh, err=10, iomsg=errmsg)
1692          end if
1693          read(wfk%fh, err=10, iomsg=errmsg) cg_bks(:,:npw_disk*wfk%nspinor)
1694        end if
1695      end do
1696 
1697      ! Reached the end of the (k,s) block. Update f90_fptr
1698      call wfk_update_f90ptr(wfk, ik_ibz, spin)
1699 
1700 #ifdef HAVE_MPI_IO
1701    case (IO_MODE_MPI)
1702      ! MPI-IO operations are done with file views even for contiguous data of the same type.
1703      ! See NOTES at the beginning of this module.
1704 
1705      sizes = [nband_disk, npw_disk*nspinor_disk]
1706      types = [MPI_DOUBLE_COMPLEX, MPI_DOUBLE_COMPLEX]
1707 
1708      call xmpio_create_fstripes(1,sizes,types,gkk_type,my_offpad,mpierr)
1709      ABI_CHECK_MPI(mpierr,"xmpio_create_fstripes")
1710 
1711      !call MPI_TYPE_CONTIGUOUS(nband_disk,MPI_DOUBLE_COMPLEX,gkk_type,mpierr)
1712      !ABI_CHECK_MPI(mpierr,"type_contigous")
1713      !call MPI_TYPE_COMMIT(gkk_type,mpierr)
1714      !ABI_CHECK_MPI(mpierr,"mpi_commit")
1715      !my_offpad = 0
1716 
1717      ! Increment my_offset to account for previous (iband -1) bands.
1718      my_offset = wfk%offset_ks(ik_ibz,spin,REC_EIG) + my_offpad
1719      my_offset = my_offset + (band - 1) * ( &
1720         (2*npw_disk*wfk%nspinor*xmpi_bsize_dp + 2*xmpio_bsize_frm) + &
1721         (2*nband_disk*xmpi_bsize_dp + 2*xmpio_bsize_frm) )
1722 
1723 #if 1
1724      bufsz = nband_disk
1725      if (sc_mode==xmpio_collective) then
1726        call MPI_FILE_READ_AT_ALL(wfk%fh,my_offset,eig1_bks,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1727      else
1728        call MPI_FILE_READ_AT(wfk%fh,my_offset,eig1_bks,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1729      end if
1730 
1731      ! Read the cg_ks(G)
1732      ! Increment my_offset to account for the previous eigen and cg records of (iband-1) bands.
1733      my_offset = wfk%offset_ks(ik_ibz,spin,REC_CG) + xmpio_bsize_frm
1734      my_offset = my_offset + (band - 1) * ( &
1735         (2*npw_disk*wfk%nspinor*xmpi_bsize_dp + 2*xmpio_bsize_frm) + &
1736         (2*nband_disk*xmpi_bsize_dp + 2*xmpio_bsize_frm) )
1737 
1738      bufsz = npw_disk * nspinor_disk
1739      if (sc_mode==xmpio_collective) then
1740        call MPI_FILE_READ_AT_ALL(wfk%fh,my_offset,cg_bks,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1741      else
1742        call MPI_FILE_READ_AT(wfk%fh,my_offset,cg_bks,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1743      end if
1744 
1745 #else
1746      call MPI_FILE_SET_VIEW(wfk%fh,my_offset,MPI_BYTE,gkk_type,'native',xmpio_info,mpierr)
1747      ABI_CHECK_MPI(mpierr,"")
1748      call MPI_TYPE_FREE(gkk_type,mpierr)
1749      ABI_CHECK_MPI(mpierr,"")
1750 
1751      bufsz = nband_disk
1752      if (sc_mode==xmpio_collective) then
1753        call MPI_FILE_READ_ALL(wfk%fh,eig1_bks,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1754      else if (sc_mode==xmpio_single) then
1755        call MPI_FILE_READ(wfk%fh,eig1_bks,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1756      else
1757        MSG_ERROR("Wrong sc_mode")
1758      end if
1759      ABI_CHECK_MPI(mpierr,"FILE_READ")
1760 
1761      ! Read the cg_ks(G)
1762      ! Increment my_offset to account for the previous eigen and cg records of (iband-1) bands.
1763      my_offset = wfk%offset_ks(ik_ibz,spin,REC_CG) + xmpio_bsize_frm
1764      my_offset = my_offset + (band - 1) * ( &
1765         (2*npw_disk*wfk%nspinor*xmpi_bsize_dp + 2*xmpio_bsize_frm) + &
1766         (2*nband_disk*xmpi_bsize_dp + 2*xmpio_bsize_frm) )
1767 
1768      call MPI_TYPE_CONTIGUOUS(npw_disk*nspinor_disk,MPI_DOUBLE_COMPLEX,cg_type,mpierr)
1769      ABI_CHECK_MPI(mpierr,"type_contigous")
1770      call MPI_TYPE_COMMIT(cg_type,mpierr)
1771      ABI_CHECK_MPI(mpierr,"mpi_commit")
1772 
1773      call MPI_FILE_SET_VIEW(wfk%fh,my_offset,MPI_BYTE,cg_type,'native',xmpio_info,mpierr)
1774      ABI_CHECK_MPI(mpierr,"")
1775      call MPI_TYPE_FREE(cg_type, mpierr)
1776      ABI_CHECK_MPI(mpierr,"MPI_TYPE_FREE")
1777 
1778      bufsz = npw_disk * nspinor_disk
1779      if (sc_mode==xmpio_collective) then
1780        call MPI_FILE_READ_ALL(wfk%fh,cg_bks,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1781      else if (sc_mode==xmpio_single) then
1782        call MPI_FILE_READ(wfk%fh,cg_bks,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
1783      else
1784        MSG_ERROR("Wrong sc_mode")
1785      end if
1786      ABI_CHECK_MPI(mpierr,"FILE_READ")
1787 #endif
1788 #endif
1789 
1790 #ifdef HAVE_NETCDF
1791    case (IO_MODE_ETSF)
1792      ! Read h1 matrix elements. The netcdf array has shape:
1793      ! [complex, max_number_of_states, max_number_of_states, number_of_kpoints, number_of_spins]
1794      NCF_CHECK(nf90_inq_varid(wfk%fh, "h1_matrix_elements", h1_varid))
1795      if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
1796        NCF_CHECK(nctk_set_collective(wfk%fh, h1_varid))
1797      end if
1798      ncerr = nf90_get_var(wfk%fh, h1_varid, eig1_bks, start=[1,1,band,ik_ibz,spin], count=[2, nband_disk, 1, 1, 1])
1799      NCF_CHECK_MSG(ncerr, "getting h1mat_k")
1800 
1801      ! Read the wavefunction.  The coefficients_of_wavefunctions on file have shape:
1802      ! [cplex, mpw, nspinor, mband, nkpt, nsppol]
1803      NCF_CHECK(nf90_inq_varid(wfk%fh, "coefficients_of_wavefunctions", cg_varid))
1804      if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
1805        NCF_CHECK(nctk_set_collective(wfk%fh, cg_varid))
1806      end if
1807 
1808      ncerr = nf90_get_var(wfk%fh, cg_varid, cg_bks, start=[1,1,1,band,ik_ibz,spin], &
1809        count=[2,npw_disk,wfk%nspinor,1,1,1])
1810      NCF_CHECK_MSG(ncerr, "getting cg_k")
1811 #endif
1812 
1813   case default
1814     MSG_ERROR(sjoin('Wrong value for iomode:', itoa(Wfk%iomode)))
1815   end select
1816  end if
1817 
1818  end if
1819 
1820  return
1821 
1822  ! Handle Fortran IO error
1823 10 continue
1824  MSG_ERROR(errmsg)
1825 
1826 end subroutine wfk_read_bks

m_wfk/wfk_read_bmask [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_read_bmask

FUNCTION

  Read a set of bands at a given k-point, spin. The bands to be read
  are specified by the logical mask `bmask`.

INPUTS

  Wfk<type(wfk_t)>=
  ik_ibz=Index of the k-point in the IBZ.
  spin=Spin index
  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading.

 OUTPUTS
  [kg_k=(:,:)] = G-vectors
  [cg_k(:,:)]  = Fourier coefficients
  [eig_k(:)] = Eigenvectors
  [occ_k(:)] = Occupation

NOTES

  The output arrays eig_k and occ_k contain the *full* set of eigenvalues and occupation
  factors stored in the file and are dimensioned with wfk%mband.

PARENTS

      m_wfd,m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

2272 subroutine wfk_read_bmask(Wfk,bmask,ik_ibz,spin,sc_mode,kg_k,cg_k,eig_k,occ_k)
2273 
2274 
2275 !This section has been created automatically by the script Abilint (TD).
2276 !Do not modify the following lines by hand.
2277 #undef ABI_FUNC
2278 #define ABI_FUNC 'wfk_read_bmask'
2279 !End of the abilint section
2280 
2281  implicit none
2282 
2283 !Arguments ------------------------------------
2284 !scalars
2285  integer,intent(in) :: ik_ibz,spin,sc_mode
2286  type(wfk_t),intent(inout) :: Wfk
2287 !arrays
2288  logical,intent(in) :: bmask(Wfk%mband)
2289  integer,intent(out), DEV_CONTARRD optional :: kg_k(:,:)  !(3,npw_k)
2290  real(dp),intent(out), DEV_CONTARRD optional :: cg_k(:,:) !(2,npw_k*nspinor*nband)
2291  real(dp),intent(out),optional :: eig_k((2*Wfk%mband)**Wfk%formeig*Wfk%mband)
2292  real(dp),intent(out),optional :: occ_k(Wfk%mband)
2293 
2294 !Local variables-------------------------------
2295 !scalars
2296  integer :: ierr,npw_disk,nspinor_disk,nband_disk,ipw,my_bcount,cnt,npwso,npw_tot,pt1,pt2,band
2297  integer :: npw_read,nspinor_read,nband_read,nb_tot,ncount,my_bcnt,my_maxb,base,nb
2298  character(len=500) :: msg,errmsg
2299 !arrays
2300  real(dp),ABI_CONTIGUOUS pointer :: tmp_eigk(:),tmp_occk(:)
2301  integer :: mpierr,cgscatter_type,cg_type,method,block,nblocks,nbxblock
2302  integer :: bstart,bstop,bufsz,ugsz,brest,max_nband
2303  integer(XMPI_OFFSET_KIND) :: my_offset,base_ofs,my_offpad
2304  integer :: band_block(2),sizes(2),subsizes(2),starts(2),types(2)
2305  integer,allocatable :: block_length(:),block_type(:)
2306  integer(XMPI_ADDRESS_KIND),allocatable :: block_displ(:)
2307  real(dp),allocatable :: buffer(:,:)
2308 #ifdef HAVE_NETCDF
2309  integer :: kg_varid,eig_varid,occ_varid,cg_varid,ncerr
2310  integer,allocatable :: blocks(:,:)
2311 #endif
2312 
2313 !************************************************************************
2314 
2315  DBG_ENTER("COLL")
2316 
2317  ABI_CHECK(Wfk%rw_mode==WFK_READMODE, "Wfk must be in READMODE")
2318 
2319  !do band=1,wfk%mband
2320  !  if (.not. bmask(band)) continue
2321  !  if (wfk_validate_ks(wfk, ik_ibz, spin, band=band) /= 0) then
2322  !    MSG_ERROR("Wrong (ik_ibz, spin, band) args, Aborting now")
2323  !  end if
2324  !end if
2325 
2326  !
2327  ! Look before you leap.
2328  npw_disk = Wfk%Hdr%npwarr(ik_ibz)
2329  nspinor_disk = Wfk%nspinor
2330  nband_disk = Wfk%nband(ik_ibz,spin)
2331  nb_tot = COUNT(bmask)
2332  npw_tot = npw_disk * nspinor_disk * nb_tot
2333 
2334  if (present(kg_k)) then
2335    ABI_CHECK((SIZE(kg_k,DIM=2) >= npw_disk),"kg_k too small")
2336  end if
2337 
2338  if (present(cg_k)) then
2339    ABI_CHECK(SIZE(cg_k, DIM=2) >= npw_tot, "Too small cg_k")
2340  end if
2341 
2342  if (present(eig_k)) then
2343    if (Wfk%formeig==0) then
2344       ABI_CHECK(SIZE(eig_k) >= nband_disk, "GS eig_k too small")
2345    else if (Wfk%formeig==1) then
2346       ABI_CHECK(SIZE(eig_k) >= 2*nband_disk**2, "DFPT eig_k too small")
2347    else
2348      MSG_ERROR("formeig != [0,1]")
2349    end if
2350  end if
2351 
2352  if (present(occ_k)) then
2353    ABI_CHECK(Wfk%formeig==0,"occ_k with formeig != 0")
2354    ABI_CHECK(SIZE(occ_k) >= nband_disk, "GS eig_k too small")
2355  end if
2356 
2357  select case (Wfk%iomode)
2358  case (IO_MODE_FORTRAN)
2359 
2360    ! Rewind the file to have the correct (k,s) block (if needed)
2361    call wfk_seek(Wfk,ik_ibz,spin)
2362 
2363    ! Read the first record: npw, nspinor, nband_disk
2364    read(Wfk%fh, err=10, iomsg=errmsg) npw_read, nspinor_read, nband_read
2365 
2366    if (any([npw_read, nspinor_read, nband_read] /= [npw_disk, nspinor_disk, nband_disk])) then
2367      write(msg,"(a,6(i0,2x))")"Mismatch between (npw, nspinor, nband) read from WFK and those found in HDR ",&
2368 &      npw_read, nspinor_read, nband_read, npw_disk, nspinor_disk, nband_disk
2369      MSG_ERROR(msg)
2370    end if
2371 
2372    ! The second record: (k+G) vectors
2373    if (present(kg_k)) then
2374      read(Wfk%fh, err=10, iomsg=errmsg) kg_k(1:3,1:npw_disk)
2375    else
2376      read(Wfk%fh, err=10, iomsg=errmsg) ! kg_k(1:3,1:npw_disk)
2377    end if
2378 
2379    ! The third record: eigenvalues and occupation factors.
2380    if (Wfk%formeig==0) then
2381 
2382      if (present(eig_k) .or. present(occ_k)) then
2383        ABI_MALLOC(tmp_eigk, (nband_disk))
2384        ABI_MALLOC(tmp_occk, (nband_disk))
2385 
2386        read(Wfk%fh, err=10, iomsg=errmsg) tmp_eigk, tmp_occk
2387 
2388        if (present(eig_k)) eig_k = tmp_eigk
2389        if (present(occ_k)) occ_k = tmp_occk
2390 
2391        ABI_FREE(tmp_eigk)
2392        ABI_FREE(tmp_occk)
2393 
2394      else
2395        read(Wfk%fh, err=10, iomsg=errmsg) ! eig_k(1:nband_disk)
2396      end if
2397 
2398      ! The wave-functions.
2399      if (present(cg_k)) then
2400        npwso = npw_disk*nspinor_disk
2401        my_bcount = 0
2402        do band=1,nband_disk
2403          if (bmask(band)) then
2404            ipw = my_bcount * npwso
2405            my_bcount = my_bcount + 1
2406            read(Wfk%fh, err=10, iomsg=errmsg) cg_k(1:2,ipw+1:ipw+npwso)
2407          else
2408            read(Wfk%fh, err=10, iomsg=errmsg) ! cg_k(1:2,ipw+1:ipw+npwso)
2409          end if
2410        end do
2411      else
2412        do band=1,nband_disk
2413          read(Wfk%fh, err=10, iomsg=errmsg) ! cg_k(1:2,ipw+1:ipw+npwso)
2414        end do
2415      end if
2416 
2417    else if (Wfk%formeig==1) then
2418      ! Read matrix of size (2*nband_k**2)
2419      npwso = npw_disk*nspinor_disk
2420      my_bcount = 0
2421 
2422      do band=1,nband_disk
2423        base = 2*(band-1)*nband_disk
2424        if (present(eig_k)) then
2425          read(Wfk%fh, err=10, iomsg=errmsg) eig_k(base+1:base+2*nband_disk)
2426        else
2427          read(Wfk%fh, err=10, iomsg=errmsg) ! eig_k(base+1:base+2*nband_disk)
2428        end if
2429 
2430        if (bmask(band).and.present(cg_k)) then
2431          ipw = my_bcount * npwso
2432          my_bcount = my_bcount + 1
2433          read(Wfk%fh, err=10, iomsg=errmsg) cg_k(1:2,ipw+1:ipw+npwso)
2434        else
2435          read(Wfk%fh, err=10, iomsg=errmsg) ! cg_k(1:2,ipw+1:ipw+npwso)
2436        end if
2437      end do
2438 
2439    else
2440      MSG_ERROR("formeig != [0,1]")
2441    end if
2442 
2443    ! Reached the end of the (k,s) block. Update f90_fptr
2444    call wfk_update_f90ptr(wfk, ik_ibz, spin)
2445 
2446 #ifdef HAVE_MPI_IO
2447  case (IO_MODE_MPI)
2448    if (present(kg_k)) then
2449      my_offset = Wfk%offset_ks(ik_ibz,spin,REC_KG) + xmpio_bsize_frm
2450 
2451      call mpio_read_kg_k(Wfk%fh,my_offset,npw_disk,sc_mode,kg_k,mpierr)
2452      ABI_CHECK_MPI(mpierr,"mpio_read_kg_k")
2453    end if
2454 
2455    ! The third record: eigenvalues and occupation factors.
2456    if (present(eig_k) .or. present(occ_k)) then
2457      my_offset = Wfk%offset_ks(ik_ibz,spin,REC_EIG) + xmpio_bsize_frm
2458      !
2459      ! formeig=0 =>  Read both eig and occ in tmp_eigk.
2460      ! formeig=1 =>  Read (nband_k,nband_k) matrix of complex numbers.
2461      !
2462      call mpio_read_eigocc_k(Wfk%fh,my_offset,nband_disk,Wfk%formeig,sc_mode,tmp_eigk,mpierr)
2463      ABI_CHECK_MPI(mpierr,"mpio_read_eigocc")
2464 
2465      if (Wfk%formeig==0) then
2466        if (present(eig_k)) eig_k(1:nband_disk) = tmp_eigk(1:nband_disk)
2467        if (present(occ_k)) occ_k(1:nband_disk) = tmp_eigk(nband_disk+1:)
2468      else if (Wfk%formeig==1) then
2469        if (present(eig_k)) eig_k(1:2*nband_disk**2) = tmp_eigk(1:2*nband_disk**2)
2470      else
2471        MSG_ERROR("formeig not in [0,1]")
2472      end if
2473 
2474      ABI_FREE(tmp_eigk)
2475    end if
2476 
2477    if (present(cg_k)) then
2478      method = 0
2479 
2480      select case (method)
2481      case (0)
2482        ! DATA SIEVING:
2483        !   read max_nband states in chuncks of nbxblock, then extract my states according to bmask.
2484        !
2485        ! MAX number of bands read by the procs in the communicator
2486        my_maxb = nband_disk
2487        do band=nband_disk,1,-1
2488          if (bmask(band)) then
2489            my_maxb = band
2490            EXIT
2491          end if
2492        end do
2493        call xmpi_max(my_maxb,max_nband,Wfk%comm,mpierr)
2494        !max_nband = nband_disk
2495        !
2496        ! MPI-IO crashes if we try to read a large number of bands in a single call.
2497        nbxblock = max_nband
2498        if ((2*npw_disk*nspinor_disk*nbxblock*xmpi_bsize_dp) > Wfk%chunk_bsize) then
2499          nbxblock = Wfk%chunk_bsize / (2*npw_disk*nspinor_disk*xmpi_bsize_dp)
2500          if (nbxblock == 0) nbxblock = 50
2501        end if
2502        !nbxblock = 2
2503 
2504        nblocks = max_nband / nbxblock
2505        brest   = MOD(max_nband, nbxblock)
2506        if (brest /= 0) nblocks = nblocks + 1
2507        !write(std_out,*)"in buffered bmask: "nb_tot",nb_tot,"nblocks",nblocks,"nbxblock",nbxblock
2508 
2509        base_ofs = Wfk%offset_ks(ik_ibz,spin,REC_CG)
2510        sizes = [npw_disk*nspinor_disk, nband_disk]
2511 
2512        my_bcnt = 0  ! index of my band in cg_k
2513        do block=1,nblocks
2514          bstart = 1 + (block-1) * nbxblock
2515          bstop  = bstart + nbxblock - 1
2516          if (bstop > max_nband) bstop = max_nband
2517          nb = bstop - bstart + 1
2518          !
2519          ! Allocate and read the buffer
2520          band_block = [bstart, bstop]
2521          ugsz = npw_disk*nspinor_disk
2522          bufsz = 2*ugsz*(bstop-bstart+1)
2523          ABI_MALLOC(buffer, (2,bufsz))
2524          !write(std_out,*)"  bstart,bstop, ",band_block
2525 
2526          ! Read the cg_ks(G). different versions depending on formeig
2527          if (wfk%formeig == 0) then
2528            subsizes = (/npw_disk*nspinor_disk, band_block(2)-band_block(1)+1/)
2529            starts = [1, bstart]
2530 
2531            call mpiotk_read_fsuba_dp2D(Wfk%fh,base_ofs,sizes,subsizes,starts,&
2532               bufsz,buffer,Wfk%chunk_bsize,sc_mode,Wfk%comm,ierr)
2533            ABI_CHECK(ierr==0,"Fortran record too big")
2534 
2535          else if (wfk%formeig == 1) then
2536 
2537            ! Increment my_offset to account for the previous eigen and cg records of (iband-1) bands.
2538            my_offset = wfk%offset_ks(ik_ibz,spin,REC_CG) + xmpio_bsize_frm
2539            my_offset = my_offset + (bstart - 1) * ( &
2540               (2*npw_disk*wfk%nspinor*xmpi_bsize_dp + 2*xmpio_bsize_frm) + &
2541               (2*nband_disk*xmpi_bsize_dp + 2*xmpio_bsize_frm) )
2542 
2543            sizes = [npw_disk*nspinor_disk, nband_disk]
2544            types = [MPI_DOUBLE_COMPLEX, MPI_DOUBLE_COMPLEX]
2545 
2546            call xmpio_create_fstripes(nb,sizes,types,cg_type,my_offpad,mpierr)
2547            ABI_CHECK_MPI(mpierr,"xmpio_create_fstripes")
2548 
2549            call MPI_FILE_SET_VIEW(wfk%fh,my_offset,MPI_BYTE,cg_type,'native',xmpio_info,mpierr)
2550            ABI_CHECK_MPI(mpierr,"")
2551            call MPI_TYPE_FREE(cg_type,mpierr)
2552            ABI_CHECK_MPI(mpierr,"")
2553 
2554            if (sc_mode==xmpio_collective) then
2555              call MPI_FILE_READ_ALL(wfk%fh,buffer,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
2556            else if (sc_mode==xmpio_single) then
2557              call MPI_FILE_READ(wfk%fh,buffer,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
2558            else
2559              MSG_ERROR("Wrong sc_mode")
2560            end if
2561            ABI_CHECK_MPI(mpierr,"FILE_READ")
2562          end if
2563 
2564          ! Extract my bands from buffer.
2565          do band=bstart,bstop
2566            if (bmask(band)) then
2567              my_bcnt = my_bcnt + 1
2568              pt1 = 1 + (my_bcnt - 1) * ugsz
2569              pt2 = 1 + (band - bstart) * ugsz
2570              cg_k(:,pt1:pt1+ugsz-1) = buffer(:,pt2:pt2+ugsz-1)
2571            end if
2572          end do
2573 
2574          ABI_FREE(buffer)
2575        end do
2576 
2577      case (1,2)
2578        ABI_CHECK(wfk%formeig == 0, "formeig == 1 not coded")
2579        call MPI_TYPE_CONTIGUOUS(npw_disk*nspinor_disk,MPI_DOUBLE_COMPLEX,cg_type,mpierr)
2580        ABI_CHECK_MPI(mpierr,"type_contigous")
2581 
2582        if (method==1) then
2583          ncount = nb_tot
2584          ABI_MALLOC(block_length,(ncount+2))
2585          ABI_MALLOC(block_type, (ncount+2))
2586          ABI_MALLOC(block_displ,(ncount+2))
2587 
2588          block_length(1)=1
2589          block_displ (1)=0
2590          block_type  (1)=MPI_LB
2591 
2592          my_bcount = 1
2593          do band=1,Wfk%mband
2594            if (bmask(band)) then
2595              my_bcount = my_bcount + 1
2596              block_length(my_bcount) = 1
2597              block_type(my_bcount) = cg_type
2598              block_displ(my_bcount) = xmpio_bsize_frm + &
2599 &              (band-1) * (2*npw_disk*nspinor_disk*xmpi_bsize_dp + 2*xmpio_bsize_frm)
2600            end if
2601          end do
2602 
2603          block_length(ncount+2) = 1
2604          block_displ (ncount+2) = block_displ(my_bcount)
2605          block_type  (ncount+2) = MPI_UB
2606 
2607        else if (method==2) then
2608          ! this file view is not efficient but it's similar to the
2609          ! one used in wff_readwrite. Let's see if MPI-IO likes it!
2610          ncount = nb_tot* nspinor_disk * npw_disk
2611 
2612          ABI_MALLOC(block_length,(ncount+2))
2613          ABI_MALLOC(block_type, (ncount+2))
2614          ABI_MALLOC(block_displ,(ncount+2))
2615 
2616          block_length(1)=1
2617          block_displ (1)=0
2618          block_type  (1)=MPI_LB
2619          !
2620          ! The view starts at REC_CG
2621          cnt = 1
2622          do band=1,Wfk%mband
2623            if (bmask(band)) then
2624              base_ofs =  xmpio_bsize_frm + &
2625 &              (band-1) * (2*npw_disk*nspinor_disk*xmpi_bsize_dp + 2*xmpio_bsize_frm)
2626              do ipw=1,npw_disk*nspinor_disk
2627                cnt = cnt + 1
2628                block_length(cnt) = 1
2629                block_type(cnt)   = MPI_DOUBLE_COMPLEX
2630                block_displ(cnt)  = base_ofs + 2*(ipw-1)*xmpi_bsize_dp
2631              end do
2632            end if
2633          end do
2634 
2635          block_length(ncount+2) = 1
2636          block_displ (ncount+2) = block_displ(cnt)
2637          block_type  (ncount+2) = MPI_UB
2638        end if
2639 
2640        call xmpio_type_struct(ncount+2,block_length,block_displ,block_type,cgscatter_type,mpierr)
2641        ABI_CHECK_MPI(mpierr,"type_struct")
2642 
2643        ABI_FREE(block_length)
2644        ABI_FREE(block_type)
2645        ABI_FREE(block_displ)
2646 
2647        call MPI_TYPE_FREE(cg_type, mpierr)
2648        ABI_CHECK_MPI(mpierr,"MPI_TYPE_FREE")
2649 
2650        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_CG)
2651 
2652        call MPI_FILE_SET_VIEW(Wfk%fh, my_offset, MPI_BYTE, cgscatter_type, 'native', xmpio_info, mpierr)
2653        ABI_CHECK_MPI(mpierr,"SET_VIEW")
2654 
2655        call MPI_TYPE_FREE(cgscatter_type, mpierr)
2656        ABI_CHECK_MPI(mpierr,"MPI_TYPE_FREE")
2657 
2658        call MPI_FILE_READ_ALL(Wfk%fh, cg_k, npw_tot, MPI_DOUBLE_COMPLEX, MPI_STATUS_IGNORE, mpierr)
2659        ABI_CHECK_MPI(mpierr,"FILE_READ_ALL")
2660 
2661      case default
2662        MSG_ERROR("Wrong method")
2663      end select
2664    end if
2665 #endif
2666 
2667 #ifdef HAVE_NETCDF
2668  case (IO_MODE_ETSF)
2669    ABI_CHECK(wfk%formeig == 0, "formeig != 0 not coded")
2670    !write(std_out,*)"bmask: ",bmask
2671 
2672    ! TODO: extract routines (see other similar calls)
2673    if (present(kg_k)) then
2674      ! Read the reduced_coordinates_of_plane_waves for this k point.
2675      NCF_CHECK(nf90_inq_varid(wfk%fh, "reduced_coordinates_of_plane_waves", kg_varid))
2676      if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
2677        NCF_CHECK(nctk_set_collective(wfk%fh, kg_varid))
2678      end if
2679 
2680      ncerr = nf90_get_var(wfk%fh, kg_varid, kg_k, start=[1,1,ik_ibz], count=[3,npw_disk,1])
2681      NCF_CHECK(ncerr)
2682    end if
2683 
2684    ! Read eigenvalues and occupations.
2685    if (Wfk%formeig==0) then
2686      if (present(eig_k)) then
2687        NCF_CHECK(nf90_inq_varid(wfk%fh, "eigenvalues", eig_varid))
2688        if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
2689          NCF_CHECK(nctk_set_collective(wfk%fh, eig_varid))
2690        end if
2691 
2692        ncerr = nf90_get_var(wfk%fh, eig_varid, eig_k, start=[1,ik_ibz,spin], count=[nband_disk,1,1])
2693        NCF_CHECK(ncerr)
2694      end if
2695 
2696      if (present(occ_k)) then
2697        NCF_CHECK(nf90_inq_varid(wfk%fh, "occupations", occ_varid))
2698        if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
2699          NCF_CHECK(nctk_set_collective(wfk%fh, occ_varid))
2700        end if
2701 
2702        ncerr = nf90_get_var(wfk%fh, occ_varid, occ_k, start=[1,ik_ibz,spin], count=[nband_disk,1,1])
2703        NCF_CHECK_MSG(ncerr, "getting occ_k")
2704      end if
2705    else
2706      MSG_ERROR("formeig !=0 not compatible with ETSF-IO")
2707    end if
2708 
2709    if (present(cg_k)) then
2710      ! The coefficients_of_wavefunctions on file have shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
2711      NCF_CHECK(nf90_inq_varid(wfk%fh, "coefficients_of_wavefunctions", cg_varid))
2712      ! TODO: Collective
2713      !if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
2714      !  NCF_CHECK(nctk_set_collective(wfk%fh, cg_varid))
2715      !end if
2716 #if 0
2717      ! Simple and very inefficient version for debugging.
2718      ipw = 1
2719      do band=1,wfk%mband
2720        if (.not. bmask(band)) cycle
2721        ncerr = nf90_get_var(wfk%fh, cg_varid, cg_k(:,ipw:), start=[1,1,1,band,ik_ibz,spin], &
2722          count=[2,npw_disk,wfk%nspinor,1,1,1])
2723        NCF_CHECK_MSG(ncerr, "getting cg_k block")
2724        ipw = ipw + wfk%nspinor * npw_disk
2725      end do
2726 #else
2727      ! Read bands in blocks defined by bmask.
2728      ! be careful when in collective mode because processors may call the routine with nblocks==0
2729      ! The coefficients_of_wavefunctions on file have shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
2730      call mask2blocks(bmask, nblocks, blocks)
2731      !ABI_CHECK(nblocks /= 0, "nblocks==0")
2732 
2733      ipw = 1
2734      do block=1,nblocks
2735        band_block = blocks(:,block)
2736        nb = band_block(2) - band_block(1) + 1
2737        ncerr = nf90_get_var(wfk%fh, cg_varid, cg_k(:,ipw:), start=[1,1,1,band_block(1),ik_ibz,spin], &
2738         count=[2,npw_disk,wfk%nspinor,nb,1,1])
2739        NCF_CHECK_MSG(ncerr, "getting cg_k block")
2740        ipw = ipw + wfk%nspinor * npw_disk * nb
2741      end do
2742      ABI_FREE(blocks)
2743 #endif
2744 
2745      ! Prototype for collective version.
2746      !min_band = lfind(bmask); if min_
2747      !max_band = lfind(bmask, back=.True.) ! n+1
2748      !! TODO: xmpi_min_max
2749      !call xmpi_max(my_min_band, min_band, comm_cell, ierr)
2750      !call xmpi_min(my_min_band, max_band, comm_cell, ierr)
2751      !nb = max_band - min_band + 1
2752      !ncalls = nb /
2753 
2754      !NCF_CHECK(nf90_var_par_access(ncid, cg_varid, nf90_collective))
2755      !do block=1,nblocks
2756      !  ! The coefficients_of_wavefunctions on file have shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
2757      !  ncerr = nf90_get_var(wfk%fh, cg_varid, cg_k, start=[1,1,1,band_block(1),ik_ibz,spin], &
2758      !   count=[2,npw_disk,wfk%nspinor,band_block(2)-band_block(1)+1,1,1])
2759      !  NCF_CHECK_MSG(ncerr, "getting cg_k block")
2760      !   do band=band_start,band_end
2761      !     if (.not. bmask(band)) cycle
2762      !     cg_k(:,:) =
2763      !   end do
2764      !end do
2765    end if
2766 #endif
2767 
2768  case default
2769    MSG_ERROR(sjoin('Wrong/unsupported value of iomode: ', itoa(wfk%iomode)))
2770  end select
2771 
2772  DBG_EXIT("COLL")
2773 
2774  return
2775 
2776  ! Handle Fortran IO error
2777 10 continue
2778  MSG_ERROR(errmsg)
2779 
2780 end subroutine wfk_read_bmask

m_wfk/wfk_read_ebands [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_read_ebands

FUNCTION

  Read the GS eigenvalues and return ebands_t object.

INPUTS

  path=WFK file name
  comm=MPI communicator

 OUTPUTS
  ebands<ebands_t>=GS band-structure.

PARENTS

CHILDREN

SOURCE

2805 type(ebands_t) function wfk_read_ebands(path, comm) result(ebands)
2806 
2807 
2808 !This section has been created automatically by the script Abilint (TD).
2809 !Do not modify the following lines by hand.
2810 #undef ABI_FUNC
2811 #define ABI_FUNC 'wfk_read_ebands'
2812 !End of the abilint section
2813 
2814  implicit none
2815 
2816 !Arguments ------------------------------------
2817 !scalars
2818  character(len=*),intent(in) :: path
2819  integer,intent(in) :: comm
2820 
2821 !Local variables-------------------------------
2822 !scalars
2823  type(hdr_type) :: hdr
2824 !arrays
2825  real(dp),pointer :: eigen(:,:,:)
2826 
2827 !************************************************************************
2828 
2829  call wfk_read_eigenvalues(path,eigen,hdr,comm)
2830  ebands = ebands_from_hdr(hdr,maxval(hdr%nband),eigen)
2831 
2832  ABI_FREE(eigen)
2833  call hdr_free(hdr)
2834 
2835 end function wfk_read_ebands

m_wfk/wfk_read_eigenvalues [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_read_eigenvalues

FUNCTION

  Read all the GS eigenvalues stored in the WFK file fname.

INPUTS

  fname=Name of the file
  comm=MPI communicator.

 OUTPUTS
  eigen = In input: nullified pointer
          In output: eigen(mband,nkpt,nsppol) contains the GS eigevalues.
  Hdr_out<hdr_type>=The header of the file

PARENTS

      dfpt_looppert,eph,m_wfk,setup_bse,setup_bse_interp,setup_screening
      setup_sigma,wfk_analyze

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

2936 subroutine wfk_read_eigenvalues(fname,eigen,Hdr_out,comm,occ)
2937 
2938 
2939 !This section has been created automatically by the script Abilint (TD).
2940 !Do not modify the following lines by hand.
2941 #undef ABI_FUNC
2942 #define ABI_FUNC 'wfk_read_eigenvalues'
2943 !End of the abilint section
2944 
2945  implicit none
2946 
2947 !Arguments ------------------------------------
2948 !scalars
2949  integer,intent(in) :: comm
2950  character(len=*),intent(in) :: fname
2951  type(hdr_type),intent(out) :: Hdr_out
2952 !arrays
2953 !MGTODO: Replace pointers with allocatable.
2954  real(dp),pointer :: eigen(:,:,:)
2955  real(dp),pointer,optional :: occ(:,:,:)
2956 
2957 !Local variables-------------------------------
2958 !scalars
2959  integer,parameter :: master=0,formeig0=0
2960  integer :: ik_ibz,spin,my_rank,ierr,iomode,funt,sc_mode,mband
2961  type(wfk_t) :: Wfk
2962 
2963 !************************************************************************
2964 
2965  iomode = iomode_from_fname(fname)
2966  my_rank = xmpi_comm_rank(comm)
2967 
2968  if (my_rank==master) then
2969    ! Open the file.
2970    sc_mode = xmpio_single
2971    funt = get_unit()
2972    call wfk_open_read(Wfk,fname,formeig0,iomode,funt,xmpi_comm_self,Hdr_out=Hdr_out)
2973 
2974    ! Read the eigenvalues and optionally the occupation factors.
2975    ABI_MALLOC(eigen, (Wfk%mband,Wfk%nkpt,Wfk%nsppol))
2976    eigen = HUGE(zero)
2977    if (present(occ)) then
2978      ABI_MALLOC(occ, (Wfk%mband,Wfk%nkpt,Wfk%nsppol))
2979      occ = HUGE(zero)
2980    end if
2981 
2982    do spin=1,Wfk%nsppol
2983      do ik_ibz=1,Wfk%nkpt
2984        if (present(occ)) then
2985          call wfk_read_eigk(Wfk,ik_ibz,spin,sc_mode,eigen(:,ik_ibz,spin),occ_k=occ(:,ik_ibz,spin))
2986        else
2987          call wfk_read_eigk(Wfk,ik_ibz,spin,sc_mode,eigen(:,ik_ibz,spin))
2988        end if
2989      end do
2990    end do
2991 
2992    ! Close the file.
2993    call wfk_close(Wfk)
2994  end if
2995 
2996  ! Broadcast data
2997  if (xmpi_comm_size(comm) > 1) then
2998    call hdr_bcast(Hdr_out,master,my_rank,comm)
2999    mband = MAXVAL(Hdr_out%nband)
3000    if (my_rank/=master) then
3001      ABI_MALLOC(eigen, (mband,Hdr_out%nkpt,Hdr_out%nsppol))
3002      if (present(occ)) then
3003        ABI_MALLOC(occ, (mband,Hdr_out%nkpt,Hdr_out%nsppol))
3004      end if
3005    end if
3006    call xmpi_bcast(eigen,master,comm,ierr)
3007    if (present(occ)) then
3008      call xmpi_bcast(occ,master,comm,ierr)
3009    end if
3010  end if
3011 
3012 end subroutine wfk_read_eigenvalues

m_wfk/wfk_read_eigk [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_read_eigk

FUNCTION

  Helper function to read all the eigenvalues for a given (k-point,spin)

INPUTS

  Wfk<type(wfk_t)>= WFK file handler
  ik_ibz=Index of the k-point in the IBZ.
  spin=spin index
  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading.

 OUTPUTS
  eig_k(1:nband_k) = GS Eigenvalues for the given (k,s)
  occ_k(1:nband_k) = Occupation factors for the given (k,s)

NOTES

  The buffers eig_k and occ_k are dimensions with wfk%mband. The routine
  will fill the first nband_k positions with data read from file where
  nband_k is the number of bands on file i.e. wfk%nband(ik_ibz,spin)

PARENTS

      conducti_nc,m_wfk,mrggkk,optic

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

2873 subroutine wfk_read_eigk(Wfk,ik_ibz,spin,sc_mode,eig_k,occ_k)
2874 
2875 
2876 !This section has been created automatically by the script Abilint (TD).
2877 !Do not modify the following lines by hand.
2878 #undef ABI_FUNC
2879 #define ABI_FUNC 'wfk_read_eigk'
2880 !End of the abilint section
2881 
2882  implicit none
2883 
2884 !Arguments ------------------------------------
2885 !scalars
2886  integer,intent(in) :: ik_ibz,spin,sc_mode
2887  type(wfk_t),intent(inout) :: Wfk
2888 !arrays
2889  real(dp),intent(out) :: eig_k((2*Wfk%mband)**Wfk%formeig*Wfk%mband)
2890  real(dp),optional,intent(out) :: occ_k(Wfk%mband)
2891 
2892 !Local variables-------------------------------
2893 !scalars
2894  integer,parameter :: band_block00(2)=[0,0]
2895 
2896 !************************************************************************
2897 
2898  if (present(occ_k)) then
2899    ABI_CHECK(Wfk%formeig==0,"formeig !=0")
2900    call wfk_read_band_block(Wfk,band_block00,ik_ibz,spin,sc_mode,eig_k=eig_k,occ_k=occ_k)
2901  else
2902    call wfk_read_band_block(Wfk,band_block00,ik_ibz,spin,sc_mode,eig_k=eig_k)
2903  end if
2904 
2905 end subroutine wfk_read_eigk

m_wfk/wfk_read_h1mat [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_read_h1mat

FUNCTION

  Read all H1 matrix elements in the WFK file fname inside the MPI communicator comm.

INPUTS

  path=File name
  comm=MPI communicator.

 OUTPUTS
  eigen(2*hdr_out%mband**2*hdr_out%nkpt*hdr_out%nsppol)=Array with the matrix elements of H1
   packed in the first positions. The array is allocated by the procedure.

  Hdr_out<hdr_type>=The header of the file

PARENTS

      m_ddk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3102 subroutine wfk_read_h1mat(fname, eigen, hdr_out, comm)
3103 
3104 
3105 !This section has been created automatically by the script Abilint (TD).
3106 !Do not modify the following lines by hand.
3107 #undef ABI_FUNC
3108 #define ABI_FUNC 'wfk_read_h1mat'
3109 !End of the abilint section
3110 
3111  implicit none
3112 
3113 !Arguments ------------------------------------
3114 !scalars
3115  character(len=*),intent(in) :: fname
3116  integer,intent(in) :: comm
3117  type(hdr_type),intent(out) :: Hdr_out
3118 !arrays
3119  real(dp),allocatable,intent(out) :: eigen(:)
3120 
3121 !Local variables-------------------------------
3122 !scalars
3123  integer,parameter :: master=0,formeig1=1
3124  integer :: spin,ik_ibz,nband_k,ptr,ierr,iomode,mband,my_rank
3125  type(wfk_t) :: wfk
3126 !arrays
3127  integer,parameter :: band_block00(2)=[0,0]
3128 
3129 !************************************************************************
3130 
3131  my_rank = xmpi_comm_rank(comm)
3132 
3133  if (my_rank==master) then
3134    ! Open the file.
3135    iomode = iomode_from_fname(fname)
3136    call wfk_open_read(wfk, fname, formeig1, iomode, get_unit(), xmpi_comm_self, hdr_out=hdr_out)
3137 
3138    ! Read h1 mat and pack them in the first positions.
3139    ABI_MALLOC(eigen, (2*wfk%mband**2*wfk%nkpt*wfk%nsppol))
3140 
3141    ptr=1
3142    do spin=1,wfk%nsppol
3143      do ik_ibz=1,wfk%nkpt
3144        nband_k = wfk%nband(ik_ibz,spin)
3145        call wfk_read_band_block(wfk, band_block00, ik_ibz, spin, xmpio_single, eig_k=eigen(ptr:))
3146        ptr = ptr + 2*nband_k**2
3147      end do
3148    end do
3149 
3150    call wfk_close(wfk)
3151  end if
3152 
3153  ! Broadcast data
3154  if (xmpi_comm_size(comm) > 1) then
3155    call hdr_bcast(hdr_out,master,my_rank,comm)
3156 
3157    mband = maxval(Hdr_out%nband)
3158    if (my_rank/=master) then
3159      ABI_MALLOC(eigen, (2*mband**2*hdr_out%nkpt*hdr_out%nsppol))
3160    end if
3161    call xmpi_bcast(eigen,master,comm,ierr)
3162  end if
3163 
3164 end subroutine wfk_read_h1mat

m_wfk/wfk_rewind [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_rewind

FUNCTION

  Rewind the file, skip the header and modifies Wfk%f90_fptr $
  Mainly used for debugging purposes when IO_MODE_FORTRAN is used.

PARENTS

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3185 subroutine wfk_rewind(wfk)
3186 
3187 
3188 !This section has been created automatically by the script Abilint (TD).
3189 !Do not modify the following lines by hand.
3190 #undef ABI_FUNC
3191 #define ABI_FUNC 'wfk_rewind'
3192 !End of the abilint section
3193 
3194  implicit none
3195 
3196 !Arguments ------------------------------------
3197  type(Wfk_t),intent(inout) :: wfk
3198 
3199 !Local variables-------------------------------
3200  integer :: ierr
3201 
3202 ! *************************************************************************
3203 
3204  select case (wfk%iomode)
3205  case (IO_MODE_FORTRAN)
3206    rewind(wfk%fh)
3207    call hdr_skip(wfk%fh,ierr)
3208    ABI_CHECK(ierr==0, "hdr_skip returned ierr! /= 0")
3209    wfk%f90_fptr = [1,1,REC_NPW]
3210 
3211  case default
3212    MSG_ERROR("should not be called when wfk%iomode /= IO_MODE_FORTRAN")
3213  end select
3214 
3215 end subroutine wfk_rewind

m_wfk/wfk_seek [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_seek

FUNCTION

   Move the internal file pointer so that it points to the
   block (ik_ibz, spin). Needed only if iomode==IO_MODE_FORTRAN

INPUTS

   ik_ibz,spin = (k-point,spin) indices

SIDE EFFECTS

   Wfk<type(Wfk_t)> : modifies Wfk%f90_fptr and the internal F90 file pointer.

PARENTS

      m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3243 subroutine wfk_seek(Wfk,ik_ibz,spin)
3244 
3245 
3246 !This section has been created automatically by the script Abilint (TD).
3247 !Do not modify the following lines by hand.
3248 #undef ABI_FUNC
3249 #define ABI_FUNC 'wfk_seek'
3250 !End of the abilint section
3251 
3252  implicit none
3253 
3254 !Arguments ------------------------------------
3255  integer,intent(in)  :: ik_ibz,spin
3256  type(Wfk_t),intent(inout) :: Wfk
3257 
3258 !Local variables-------------------------------
3259  integer :: ierr,ik_fpt,spin_fpt,recn_wanted,recn_fpt,rec_type
3260  character(len=500) :: msg
3261 
3262 ! *************************************************************************
3263 
3264  select case (Wfk%iomode)
3265  case (IO_MODE_FORTRAN)
3266    !
3267    ! Find the position inside the file.
3268    if (ALL(Wfk%f90_fptr==FPTR_EOF)) then ! handle the EOF condition
3269      if (Wfk%debug) then
3270        call wrtout(std_out,"EOF condition","PERS")
3271      end if
3272      recn_fpt = Wfk%recn_eof
3273    else
3274      ik_fpt   = Wfk%f90_fptr(1)
3275      spin_fpt = Wfk%f90_fptr(2)
3276      rec_type = Wfk%f90_fptr(3)
3277      recn_fpt = Wfk%recn_ks(ik_fpt,spin_fpt, rec_type)
3278    end if
3279 
3280    recn_wanted = Wfk%recn_ks(ik_ibz,spin, REC_NPW)
3281 
3282    if (Wfk%debug) then
3283      write(msg,'(a,3(i0,2x))')"seeking ik_ibz, spin, recn_wanted-recn_fpt: ",ik_ibz,spin,recn_wanted - recn_fpt
3284      call wrtout(std_out,msg,"PERS")
3285    end if
3286 
3287    call mvrecord(Wfk%fh, (recn_wanted - recn_fpt) ,ierr)
3288    ABI_CHECK(ierr==0,"error in mvrecord")
3289 
3290    Wfk%f90_fptr = [ik_ibz, spin, REC_NPW]
3291 
3292  case default
3293    MSG_ERROR("should not be called when Wfk%iomode /= IO_MODE_FORTRAN")
3294  end select
3295 
3296 end subroutine wfk_seek

m_wfk/wfk_t [ Types ]

[ Top ] [ m_wfk ] [ Types ]

NAME

  wfk_t

FUNCTION

  File handler for the WFK file.

SOURCE

115  type,public :: wfk_t
116 
117   integer :: fh
118    !  unit number if IO_MODE_FORTRAN
119    !  MPI file handler if IO_MODE_MPI
120    !  Netcdf file handler if IO_MODE_ETSF
121 
122   integer :: iomode
123    ! Method used to access the WFK file:
124    !   IO_MODE_FORTRAN for usual Fortran IO routines
125    !   IO_MODE_MPI if MPI/IO routines.
126    !   IO_MODE_ETSF, NetCDF format read/written via etsf-io.
127 
128   integer :: mband
129   ! Max number of bands stored on file (MAX(Hdr%nband))
130 
131   integer :: nkpt
132   ! Number of k-points.
133 
134   integer :: nsppol
135   ! Number of spins
136 
137   integer :: nspinor
138   ! Number of spinor components.
139 
140   integer :: formeig
141    ! formeig=format of the eigenvalues
142    !    0 => vector of eigenvalues
143    !    1 => hermitian matrix of eigenvalues
144    ! TODO: this should be reported somewhere in the WFK file, at present is passed to wfk_open
145 
146   integer :: fform
147    ! File type format of the header
148 
149   integer :: rw_mode = WFK_NOMODE
150    ! (Read|Write) mode
151 
152   character(len=fnlen) :: fname = ABI_NOFILE
153    ! File name
154 
155   integer :: master
156    ! master node of the IO procedure
157 
158   integer :: my_rank
159    ! index of my processor in the MPI communicator comm
160 
161   integer :: nproc
162    ! number of processors in comm
163 
164   integer :: comm
165    ! MPI communicator
166 
167   integer :: recn_eof
168    ! EOF record number (used for Fortran IO)
169 
170   integer(XMPI_OFFSET_KIND) :: offset_eof
171   ! EOF offset (used for MPI-IO access)
172 
173   logical :: debug=.FALSE.
174   !logical :: debug=.TRUE.
175 
176   type(hdr_type) :: Hdr
177    ! Abinit header.
178 
179   integer,allocatable :: nband(:,:)
180   ! nband(nkpt,nsppol) = Number of bands at each (k,s)
181 
182   integer :: f90_fptr(3) = [0,0,0]
183   ! The position of the file pointer used for sequential access with Fortran-IO.
184   !  f90_fprt(1) = Index of the k-point associated to the block.
185   !  f90_fprt(2) = the spin associated to the block.
186   !  f90_fprt(3) = Record Type (see REC_* variables).
187   !  [0,0,0] corresponds to the beginning of the file.
188   !  FPTR_EOF signals the end of file
189 
190   integer,allocatable :: recn_ks(:,:,:)
191    ! recn_ks(k,s,1) : record number of  (npw, nspinor, nband_disk)
192    ! recn_ks(k,s,2) : record number of the (k+G) vectors.
193    ! recn_ks(k,s,3) : record number of the eigenvalues.
194    ! recn_ks(k,s,4) : record number of the first wavefunction in the wf coefficients block.
195 
196   integer(XMPI_OFFSET_KIND),allocatable :: offset_ks(:,:,:)
197    ! offset_ks(k,s,1) : offset of the record: npw, nspinor, nband_disk.
198    ! offset_ks(k,s,2) : offset of the Second record: (k+G) vectors.
199    ! offset_ks(k,s,3) : offset of the third record eigenvalues.
200    ! offset_ks(k,s,4) : offset of the fourth record (wavefunction coefficients).
201    !
202    ! **********************************************************************
203    ! NB: The offset point to the Fortran record marker and not to the data
204    ! **********************************************************************
205 
206   integer(XMPI_OFFSET_KIND) :: hdr_offset
207    ! offset of the header
208    ! TODO this should be the output of a hdr method!
209 
210   integer(XMPI_OFFSET_KIND) :: chunk_bsize
211    ! IO is performed in chunks of max size chunk_bsize [bytes]
212 
213  end type wfk_t
214 
215 !public procedures.
216  public :: wfk_open_read           ! Open the WFK file in read mode.
217  public :: wfk_open_write          ! Open the WFK file in write mode.
218  public :: wfk_close               ! Close the WFK file and release the memory allocated in wfk_t.
219  public :: wfk_print               ! Print info on the wfk_t object
220  public :: wfk_findk               ! Returns the index of the k-point in the WFK file.
221  public :: wfk_ncdef_dims_vars     ! Define basic dimensions for netcdf file format.
222  public :: wfk_compare             ! Test two wfk_t objects for consistency.
223  public :: wfk_read_band_block     ! Read a contiguous block of bands for a given (kpoint, spin)
224  public :: wfk_read_bks            ! Read the wavefunction and the eigenvalues for a given (band, k-point, spin)
225  public :: wfk_write_band_block    ! Write a contiguous block of bands for a given (kpoint, spin)
226  public :: wfk_read_bmask          ! Read a scattered set of bands for a given (kpoint, spin).
227  public :: wfk_read_ebands         ! Read the GS eigenvalues and return ebands_t object.
228  public :: wfk_read_eigk           ! Read the eigenvalues at a given (kpoint,spin).
229  public :: wfk_read_eigenvalues    ! Read all the GS eigenvalues stored in the WFK file.
230  public :: wfk_write_h1mat         ! Write all the H1 matrix elements.
231  public :: wfk_read_h1mat          ! Read all the H1 matrix elements.
232  public :: wfk_tofullbz            ! Generate a new WFK file with wavefunctions in the full BZ and istwfk==1
233                                    ! Mainly used to interface ABINIT with other codes that
234                                    ! cannot handle symmetries e.g. lobster
235  public :: wfk_nc2fort             ! Convert a netcdf WFK file to a Fortran WFK file.
236 
237  ! Profiling tools
238  public :: wfk_prof                ! Profiling tool.
239 
240  ! Unit tests
241  public :: wfk_diff                ! Compare two WFK file for binary equality.
242  public :: wfk_create_wfkfile      ! Create a FAKE WFK file.
243  public :: wfk_check_wfkfile       ! Read a FAKE WFK file and perform basic tests.

m_wfk/wfk_tofullbz [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_tofullbz

FUNCTION

 Generate a new WFK file with wavefunctions in the full BZ and istwfk==1
 Mainly used to interface ABINIT with other codes that cannot handle symmetries e.g. lobster

INPUTS

  in_path = Input WFK file
  dtset <dataset_type>=all input variables for this dataset
  psps <pseudopotential_type>=all the information about psps
  pawtab(ntypat*usepaw) <type(pawtab_type)>=paw tabulated starting data
  out_path = Output WFK file.

OUTPUT

  Output is written to file out_path

NOTES

  - This routine should be called by a single processor.
  - Only GS WFK files are supported (formeig==0)

PARENTS

      gstate,wfk_analyze

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

4055 subroutine wfk_tofullbz(in_path, dtset, psps, pawtab, out_path)
4056 
4057 
4058 !This section has been created automatically by the script Abilint (TD).
4059 !Do not modify the following lines by hand.
4060 #undef ABI_FUNC
4061 #define ABI_FUNC 'wfk_tofullbz'
4062 !End of the abilint section
4063 
4064  implicit none
4065 
4066 !Arguments ------------------------------------
4067 !scalars
4068  character(len=*),intent(in) :: in_path,out_path
4069  type(pseudopotential_type),intent(in) :: psps
4070  type(dataset_type),intent(in) :: dtset
4071 !arrays
4072  type(pawtab_type),intent(in) :: pawtab(dtset%ntypat*psps%usepaw)
4073 
4074 !Local variables-------------------------------
4075 !scalars
4076  integer,parameter :: formeig0=0,kptopt3=3
4077  integer :: spin,ikf,ik_ibz,nband_k,mpw_ki,mpw_kf,mband,nspinor,nkfull
4078  integer :: in_iomode,nsppol,nkibz,out_iomode,isym,itimrev
4079  integer :: npw_ki,npw_kf,istwf_ki,istwf_kf,ii,jj,iqst,nqst
4080  real(dp) :: ecut_eff,dksqmax,cpu,wall,gflops
4081  character(len=500) :: msg
4082  character(len=fnlen) :: my_inpath
4083  logical :: isirred_kf
4084  logical,parameter :: force_istwfk1=.True.
4085  type(wfk_t),target :: in_wfk
4086  type(wfk_t) :: out_wfk
4087  type(crystal_t) :: cryst
4088  type(hdr_type) :: hdr_kfull
4089  type(hdr_type),pointer :: ihdr
4090  type(ebands_t) :: ebands_ibz
4091  type(ebands_t),target :: ebands_full
4092  type(wvl_internal_type) :: dummy_wvl
4093 !arrays
4094  integer :: g0(3),work_ngfft(18),gmax_ki(3),gmax_kf(3),gmax(3)
4095  integer,allocatable :: bz2ibz(:,:),kg_ki(:,:),kg_kf(:,:),iperm(:),bz2ibz_sort(:)
4096  real(dp) :: kf(3),kibz(3)
4097  real(dp),allocatable :: cg_ki(:,:),cg_kf(:,:),eig_ki(:),occ_ki(:),work(:,:,:,:)
4098  real(dp), ABI_CONTIGUOUS pointer :: kfull(:,:)
4099 
4100 ! *************************************************************************
4101 
4102  if (all(dtset%kptrlatt == 0)) then
4103    write(msg,"(5a)")&
4104      "Cannot produce full WFK file because kptrlatt == 0",ch10,&
4105      "Please use nkgpt and shiftk to define a homogeneous k-mesh.",ch10,&
4106      "Returning to caller"
4107    MSG_WARNING(msg)
4108    return
4109  end if
4110 
4111  call cwtime(cpu,wall,gflops,"start")
4112 
4113  my_inpath = in_path
4114 
4115  if (nctk_try_fort_or_ncfile(my_inpath, msg) /= 0) then
4116    MSG_ERROR(msg)
4117  end if
4118 
4119  call wrtout(std_out, sjoin("Converting:", my_inpath, "to", out_path))
4120  in_iomode = iomode_from_fname(my_inpath)
4121 
4122  ebands_ibz = wfk_read_ebands(my_inpath, xmpi_comm_self)
4123 
4124  ! Open input file, extract dimensions and allocate workspace arrays.
4125  call wfk_open_read(in_wfk,my_inpath,formeig0,in_iomode,get_unit(),xmpi_comm_self)
4126  ihdr => in_wfk%hdr
4127 
4128  mband = in_wfk%mband; mpw_ki = maxval(in_wfk%Hdr%npwarr); nkibz = in_wfk%nkpt
4129  nsppol = in_wfk%nsppol; nspinor = in_wfk%nspinor
4130  ecut_eff = in_wfk%hdr%ecut_eff ! ecut * dilatmx**2
4131 
4132  ABI_MALLOC(kg_ki, (3, mpw_ki))
4133  ABI_MALLOC(cg_ki, (2, mpw_ki*nspinor*mband))
4134  ABI_MALLOC(eig_ki, ((2*mband)**in_wfk%formeig*mband) )
4135  ABI_MALLOC(occ_ki, (mband))
4136 
4137  call crystal_from_hdr(cryst, in_wfk%hdr, 2)
4138 
4139  ! Build new header for out_wfk. This is the most delicate part since all the arrays in hdr_full
4140  ! that depend on k-points must be consistent with kfull and nkfull.
4141  call ebands_expandk(ebands_ibz, cryst, ecut_eff, force_istwfk1, dksqmax, bz2ibz, ebands_full)
4142 
4143  if (dksqmax > tol12) then
4144    write(msg, '(3a,es16.6,4a)' )&
4145    'At least one of the k points could not be generated from a symmetrical one.',ch10,&
4146    'dksqmax=',dksqmax,ch10,&
4147    'Action: check your WFK file and k-point input variables',ch10,&
4148    '        (e.g. kptopt or shiftk might be wrong in the present dataset or the preparatory one.'
4149    MSG_ERROR(msg)
4150  end if
4151 
4152  nkfull = ebands_full%nkpt
4153  kfull => ebands_full%kptns
4154 
4155  ! Build new header and update pawrhoij.
4156  call hdr_init_lowlvl(hdr_kfull,ebands_full,psps,pawtab,dummy_wvl,abinit_version,&
4157    ihdr%pertcase,ihdr%natom,ihdr%nsym,ihdr%nspden,ihdr%ecut,dtset%pawecutdg,ihdr%ecutsm,dtset%dilatmx,&
4158    ihdr%intxc,ihdr%ixc,ihdr%stmbias,ihdr%usewvl,dtset%pawcpxocc,dtset%pawspnorb,dtset%ngfft,dtset%ngfftdg,ihdr%so_psp,&
4159    ihdr%qptn,cryst%rprimd,cryst%xred,ihdr%symrel,ihdr%tnons,ihdr%symafm,ihdr%typat,ihdr%amu,ihdr%icoulomb,&
4160    kptopt3,dtset%nelect,dtset%charge,dtset%kptrlatt_orig,dtset%kptrlatt,&
4161    dtset%nshiftk_orig,dtset%nshiftk,dtset%shiftk_orig,dtset%shiftk)
4162 
4163  if (psps%usepaw == 1) call pawrhoij_copy(in_wfk%hdr%pawrhoij, hdr_kfull%pawrhoij)
4164 
4165  out_iomode = iomode_from_fname(out_path)
4166  call wfk_open_write(out_wfk,hdr_kfull,out_path,in_wfk%formeig,out_iomode,get_unit(),xmpi_comm_self)
4167  call hdr_free(hdr_kfull)
4168 
4169  ! workspace array for BZ wavefunction block.
4170  mpw_kf = maxval(ebands_full%npwarr)
4171  ABI_MALLOC(cg_kf, (2,mpw_kf*nspinor*mband))
4172 
4173  if (out_iomode == IO_MODE_FORTRAN) then
4174    call wrtout(std_out,"Using (slow) Fortran IO version to generate full WFK file", do_flush=.True.)
4175 
4176    ! Fortran IO does not support random access hence the external loop is on the k-points in the full BZ.
4177    !
4178    !   For each point in the BZ:
4179    !     - Find symmetric k-point in the IBZ and read IBZ wavefunctions from in_wfk
4180    !     - Rotate wavefunctions in G-space and write kbz data
4181    !
4182    ! Inefficient since we are reading the same IBZ block several times.
4183    do spin=1,nsppol
4184      do ikf=1,nkfull
4185        ik_ibz = bz2ibz(ikf,1); isym = bz2ibz(ikf,2); itimrev = bz2ibz(ikf,6); g0 = bz2ibz(ikf,3:5) ! IS(k_ibz) + g0 = k_bz
4186        isirred_kf = (isym == 1 .and. itimrev == 0 .and. all(g0 == 0))
4187 
4188        nband_k = in_wfk%nband(ik_ibz,spin)
4189        kf = kfull(:,ikf)
4190        kibz = ebands_ibz%kptns(:,ik_ibz)
4191 
4192        istwf_ki = in_wfk%hdr%istwfk(ik_ibz)
4193        istwf_kf = out_wfk%hdr%istwfk(ikf)
4194        npw_ki = in_wfk%hdr%npwarr(ik_ibz)
4195 
4196        ! Read IBZ data.
4197        call wfk_read_band_block(in_wfk,[1,nband_k],ik_ibz,spin,xmpio_single,&
4198          kg_k=kg_ki,cg_k=cg_ki,eig_k=eig_ki,occ_k=occ_ki)
4199 
4200        ! The test on npwarr is needed because we may change istwfk e.g. gamma.
4201        if (isirred_kf .and. in_wfk%hdr%npwarr(ik_ibz) == out_wfk%hdr%npwarr(ikf)) then
4202 
4203          call wfk_write_band_block(out_wfk,[1,nband_k],ikf,spin,xmpio_single,&
4204            kg_k=kg_ki,cg_k=cg_ki,eig_k=eig_ki,occ_k=occ_ki)
4205 
4206        else
4207          ! Compute G-sphere centered on kf
4208          call get_kg(kf,istwf_kf,ecut_eff,cryst%gmet,npw_kf,kg_kf)
4209          ABI_CHECK(npw_kf == out_wfk%hdr%npwarr(ikf), "Wrong npw_kf")
4210 
4211          ! FFT box must enclose the two spheres centered on ki and kf
4212          gmax_ki = maxval(abs(kg_ki(:,1:npw_ki)), dim=2)
4213          gmax_kf = maxval(abs(kg_kf), dim=2)
4214          do ii=1,3
4215            gmax(ii) = max(gmax_ki(ii), gmax_kf(ii))
4216          end do
4217          gmax = 2*gmax + 1
4218          call ngfft_seq(work_ngfft, gmax)
4219          ABI_CALLOC(work, (2, work_ngfft(4),work_ngfft(5),work_ngfft(6)))
4220 
4221          ! Rotate nband_k wavefunctions (output in cg_kf)
4222          call cgtk_rotate(cryst,kibz,isym,itimrev,g0,nspinor,nband_k,&
4223            npw_ki,kg_ki,npw_kf,kg_kf,istwf_ki,istwf_kf,cg_ki,cg_kf,work_ngfft,work)
4224 
4225          ABI_FREE(work)
4226 
4227          ! Write data
4228          call wfk_write_band_block(out_wfk,[1,nband_k],ikf,spin,xmpio_single,&
4229            kg_k=kg_kf,cg_k=cg_kf,eig_k=eig_ki,occ_k=occ_ki)
4230 
4231          ABI_FREE(kg_kf)
4232        end if
4233      end do
4234    end do
4235 
4236  else
4237    !
4238    ! More efficienct algorithm based on random access IO:
4239    !   For each point in the IBZ:
4240    !     - Read wavefunctions from in_wfk
4241    !     - For each k-point in the star of kpt_ibz:
4242    !        - Rotate wavefunctions in G-space to get the k-point in the full BZ.
4243    !        - Write kbz data to file.
4244    if (out_iomode == IO_MODE_MPI) call wrtout(std_out,"Using MPI-IO to generate full WFK file", do_flush=.True.)
4245    if (out_iomode == IO_MODE_ETSF) call wrtout(std_out,"Using Netcdf-IO to generate full WFK file", do_flush=.True.)
4246 
4247    ! Construct sorted mapping BZ --> IBZ to speedup qbz search below.
4248    ABI_MALLOC(iperm, (nkfull))
4249    ABI_MALLOC(bz2ibz_sort, (nkfull))
4250    iperm = [(ii, ii=1,nkfull)]
4251    bz2ibz_sort = bz2ibz(:,1)
4252    call sort_int(nkfull, bz2ibz_sort, iperm)
4253 
4254    do spin=1,nsppol
4255      iqst = 0
4256      do ik_ibz=1,in_wfk%nkpt
4257        nband_k = in_wfk%nband(ik_ibz,spin)
4258        kibz = ebands_ibz%kptns(:,ik_ibz)
4259        istwf_ki = in_wfk%hdr%istwfk(ik_ibz)
4260        npw_ki = in_wfk%hdr%npwarr(ik_ibz)
4261 
4262        call wfk_read_band_block(in_wfk,[1,nband_k],ik_ibz,spin,xmpio_single,&
4263          kg_k=kg_ki,cg_k=cg_ki,eig_k=eig_ki,occ_k=occ_ki)
4264 
4265        ! Find number of symmetric q-ponts associated to ik_ibz
4266        nqst = 0
4267        do ii=iqst+1,nkfull
4268          if (bz2ibz_sort(ii) /= ik_ibz) exit
4269          nqst = nqst + 1
4270        end do
4271        ABI_CHECK(nqst > 0 .and. bz2ibz_sort(iqst+1) == ik_ibz, "Wrong iqst")
4272 
4273        do jj=1,nqst
4274          iqst = iqst + 1
4275          ikf = iperm(iqst)
4276          ABI_CHECK(ik_ibz == bz2ibz(ikf,1), "ik_ibz !/ ind qq(1)")
4277 
4278          isym = bz2ibz(ikf,2); itimrev = bz2ibz(ikf,6); g0 = bz2ibz(ikf,3:5) ! IS(k_ibz) + g0 = k_bz
4279          isirred_kf = (isym == 1 .and. itimrev == 0 .and. all(g0 == 0))
4280 
4281          kf = kfull(:,ikf)
4282          istwf_kf = out_wfk%hdr%istwfk(ikf)
4283 
4284          ! The test on npwarr is needed because we may change istwfk e.g. gamma.
4285          if (isirred_kf .and. in_wfk%hdr%npwarr(ik_ibz) == out_wfk%hdr%npwarr(ikf)) then
4286 
4287            call wfk_write_band_block(out_wfk,[1,nband_k],ikf,spin,xmpio_single,&
4288              kg_k=kg_ki,cg_k=cg_ki,eig_k=eig_ki,occ_k=occ_ki)
4289 
4290          else
4291            ! Compute G-sphere centered on kf
4292            call get_kg(kf,istwf_kf,ecut_eff,cryst%gmet,npw_kf,kg_kf)
4293            ABI_CHECK(npw_kf == out_wfk%hdr%npwarr(ikf), "Wrong npw_kf")
4294 
4295            ! FFT box must enclose the two spheres centered on ki and kf
4296            gmax_ki = maxval(abs(kg_ki(:,1:npw_ki)), dim=2)
4297            gmax_kf = maxval(abs(kg_kf), dim=2)
4298            do ii=1,3
4299              gmax(ii) = max(gmax_ki(ii), gmax_kf(ii))
4300            end do
4301            gmax = 2*gmax + 1
4302            call ngfft_seq(work_ngfft, gmax)
4303            ABI_CALLOC(work, (2, work_ngfft(4),work_ngfft(5),work_ngfft(6)))
4304 
4305            ! Rotate nband_k wavefunctions (output in cg_kf)
4306            call cgtk_rotate(cryst,kibz,isym,itimrev,g0,nspinor,nband_k,&
4307              npw_ki,kg_ki,npw_kf,kg_kf,istwf_ki,istwf_kf,cg_ki,cg_kf,work_ngfft,work)
4308 
4309            ABI_FREE(work)
4310 
4311            ! Write data
4312            call wfk_write_band_block(out_wfk,[1,nband_k],ikf,spin,xmpio_single,&
4313              kg_k=kg_kf,cg_k=cg_kf,eig_k=eig_ki,occ_k=occ_ki)
4314 
4315            ABI_FREE(kg_kf)
4316          end if
4317        end do
4318      end do
4319    end do
4320 
4321    ABI_FREE(iperm)
4322    ABI_FREE(bz2ibz_sort)
4323  end if
4324 
4325  call cwtime(cpu,wall,gflops,"stop")
4326  write(std_out,"(2(a,f8.2))")" FULL_WFK written to file.  cpu: ",cpu,", wall:",wall
4327 
4328  ABI_FREE(kg_ki)
4329  ABI_FREE(cg_ki)
4330  ABI_FREE(eig_ki)
4331  ABI_FREE(occ_ki)
4332  ABI_FREE(bz2ibz)
4333  ABI_FREE(cg_kf)
4334 
4335  call crystal_free(cryst)
4336  call ebands_free(ebands_ibz)
4337  call ebands_free(ebands_full)
4338  call wfk_close(in_wfk)
4339  call wfk_close(out_wfk)
4340 
4341 end subroutine wfk_tofullbz

m_wfk/wfk_update_f90ptr [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_update_f90ptr

FUNCTION

  Update wfk%f90_ptr. Used if wfk%iomode == IO_MODE_FORTRAN.

INPUTS

  ik_ibz=K-point index,
  spin=Spin index.

PARENTS

      m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3321 subroutine wfk_update_f90ptr(wfk, ik_ibz, spin)
3322 
3323 
3324 !This section has been created automatically by the script Abilint (TD).
3325 !Do not modify the following lines by hand.
3326 #undef ABI_FUNC
3327 #define ABI_FUNC 'wfk_update_f90ptr'
3328 !End of the abilint section
3329 
3330  implicit none
3331 
3332 !Arguments ------------------------------------
3333  type(Wfk_t),intent(inout) :: wfk
3334  integer,intent(in) :: ik_ibz,spin
3335 
3336 ! *************************************************************************
3337 
3338  if (ik_ibz < wfk%nkpt) then
3339    wfk%f90_fptr = [ik_ibz+1,spin,REC_NPW]
3340  else
3341    ABI_CHECK(ik_ibz == wfk%nkpt, "ik_ibz != nkpt")
3342    if (spin == wfk%nsppol) then
3343      wfk%f90_fptr = FPTR_EOF ! EOF condition
3344    else
3345      wfk%f90_fptr = [1,spin+1,REC_NPW]
3346    end if
3347  end if
3348 
3349 end subroutine wfk_update_f90ptr

m_wfk/wfk_validate_ks [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_validate_ks

FUNCTION

  Validate the k-point, the spin index and, optionally, the band index.
  Return non-zero value if error.

INPUTS

  wfk<type(wfk_t)> = WFK handler
  ik_ibz=k-point index.
  spin=Spin index.
  [band]=Band index.

PARENTS

CHILDREN

SOURCE

843 integer function wfk_validate_ks(wfk, ik_ibz, spin, band) result(ierr)
844 
845 
846 !This section has been created automatically by the script Abilint (TD).
847 !Do not modify the following lines by hand.
848 #undef ABI_FUNC
849 #define ABI_FUNC 'wfk_validate_ks'
850 !End of the abilint section
851 
852  implicit none
853 
854 !Arguments ------------------------------------
855 !scalars
856  integer,intent(in) :: ik_ibz, spin
857  integer,optional,intent(in) :: band
858  type(wfk_t),intent(in) :: wfk
859 
860 !Local variables-------------------------------
861 !scalars
862  character(len=500) :: msg
863 
864 ! *************************************************************************
865  ierr = 0
866 
867  if (ik_ibz <= 0 .or. ik_ibz > wfk%nkpt) then
868    ierr = ierr + 1
869    write(msg, '(2(a,i0))')'ik_ibz = ',ik_ibz,' whereas it should be between 1 and ',wfk%nkpt
870    MSG_WARNING(msg)
871  end if
872 
873  if (spin <= 0 .or. spin > wfk%nsppol) then
874    ierr = ierr + 1
875    write(msg, '(2(a,i0))')'spin = ',spin,' whereas it should be between 1 and ',wfk%nsppol
876    MSG_WARNING(msg)
877  end if
878 
879  if (present(band)) then
880    if (band <=0) then
881      ierr = ierr + 1
882      MSG_WARNING(sjoin('Negative band index: band = ',itoa(band)))
883    end if
884 
885    ! Don't touch nband array if wrong indices.
886    if (spin > 0 .and. spin <= wfk%nsppol .and. ik_ibz > 0 .and. ik_ibz <= wfk%nkpt) then
887       if (band > wfk%nband(ik_ibz, spin)) then
888         ierr = ierr + 1
889         write(msg, '(2(a,i0))')'band = ',band,' whereas it should be between 1 and ',wfk%nband(ik_ibz,spin)
890         MSG_WARNING(msg)
891       end if
892    end if
893  end if
894 
895  !if (ierr /= 0) then
896  !  MSG_ERROR("Wrong (ik_ibz, spin) args, Aborting now")
897  !end if
898 
899 end function wfk_validate_ks

m_wfk/wfk_write_band_block [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_write_band_block

FUNCTION

  Write a block of contigous bands.

INPUTS

  Wfk<type(wfk_t)>=
  band_block(2)=Initial and final band index.
  ik_ibz=Index of the k-point in the IBZ.
  spin=Spin index
  sc_mode= MPI-IO option
    xmpio_single     ==> for reading by current proc.
    xmpio_collective ==> for collective reading.

 OUTPUTS
  [kg_k=(:,:)] = G-vectors
  [cg_k(:,:)]  = Fourier coefficients
  [eig_k(:)] = Eigenvectors
  [occ_k(:)] = Eigenvectors

PARENTS

      m_iowf,m_wfd,m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

1860 subroutine wfk_write_band_block(Wfk,band_block,ik_ibz,spin,sc_mode,kg_k,cg_k,eig_k,occ_k)
1861 
1862 
1863 !This section has been created automatically by the script Abilint (TD).
1864 !Do not modify the following lines by hand.
1865 #undef ABI_FUNC
1866 #define ABI_FUNC 'wfk_write_band_block'
1867 !End of the abilint section
1868 
1869  implicit none
1870 
1871 !Arguments ------------------------------------
1872 !scalars
1873  integer,intent(in) :: ik_ibz,spin,sc_mode !,mband,rdcg,rdeig,npw_k,nband_k
1874  type(wfk_t),intent(inout) :: Wfk
1875 !arrays
1876  integer,intent(in) :: band_block(2)
1877  integer,intent(in),optional :: kg_k(:,:)  !(3,npw_k)
1878  real(dp),intent(in),optional :: cg_k(:,:) ! cg_k(2,rdcg*cgsize2) !(2,npw_k*nspinor*nband)
1879  real(dp),intent(in),optional :: eig_k((2*Wfk%mband)**Wfk%formeig*Wfk%mband)
1880  real(dp),intent(in),optional :: occ_k(Wfk%mband)
1881 
1882 !Local variables-------------------------------
1883 !scalars
1884  integer :: ierr,npw_disk,nspinor_disk,nband_disk,band
1885  integer :: ipw,my_bcount,npwso,npw_tot,nb_block,base
1886  character(len=500) :: errmsg !msg,
1887 !arrays
1888  real(dp),ABI_CONTIGUOUS pointer :: tmp_eigk(:)
1889 #ifdef HAVE_MPI_IO
1890  integer :: mpierr,bufsz,recnpw_type,gkk_type,cgblock_type
1891  integer(XMPI_OFFSET_KIND) :: my_offset,my_offpad
1892  integer :: sizes(2),subsizes(2),starts(2),dims(3),types(2)
1893  integer(XMPI_OFFSET_KIND) :: bsize_rec(1)
1894  integer(XMPI_OFFSET_KIND),allocatable :: bsize_frecords(:)
1895 #endif
1896 #ifdef HAVE_NETCDF
1897  integer :: kg_varid,eig_varid,occ_varid,cg_varid,ncerr
1898 #endif
1899 
1900 !************************************************************************
1901 
1902  DBG_ENTER("COLL")
1903 
1904  ABI_CHECK(Wfk%rw_mode==WFK_WRITEMODE, "Wfk must be in WRITEMODE")
1905 
1906  ! Look before you leap.
1907  npw_disk     = Wfk%Hdr%npwarr(ik_ibz)
1908  nspinor_disk = Wfk%nspinor
1909  nband_disk   = Wfk%nband(ik_ibz,spin)
1910  nb_block     = (band_block(2) - band_block(1) + 1)
1911  npw_tot      = npw_disk * nspinor_disk * nb_block
1912 
1913  if (PRESENT(kg_k)) then
1914    ABI_CHECK(SIZE(kg_k,DIM=2) >= npw_disk,"kg_k too small")
1915  end if
1916 
1917  if (PRESENT(cg_k)) then
1918    ABI_CHECK(SIZE(cg_k, DIM=2) >= npw_tot,"cg_k too small")
1919  end if
1920 
1921  if (PRESENT(eig_k)) then
1922    if (Wfk%formeig==0) then
1923       ABI_CHECK(SIZE(eig_k) >= nband_disk, "GS eig_k too small")
1924       ABI_CHECK(PRESENT(occ_k),"both eig_k and occ_k must be present")
1925    else if (Wfk%formeig==1) then
1926       ABI_CHECK(SIZE(eig_k) >= 2*nband_disk**2, "DFPT eig_k too small")
1927    else
1928      MSG_ERROR("formeig != [0,1]")
1929    end if
1930  end if
1931 
1932  if (PRESENT(occ_k)) then
1933    ABI_CHECK(SIZE(occ_k) >= nband_disk, "GS eig_k too small")
1934    ABI_CHECK(PRESENT(eig_k),"both eig_k and occ_k must be present")
1935    ABI_CHECK(Wfk%formeig == 0, "formeig /=0 with occ_k in input!")
1936  end if
1937 
1938  select case (Wfk%iomode)
1939  case (IO_MODE_FORTRAN)
1940 
1941    ! Rewind the file to have the correct (k,s) block (if needed)
1942    call wfk_seek(Wfk,ik_ibz,spin)
1943 
1944    ! Write the first record: npw, nspinor, nband_disk
1945    write(Wfk%fh, err=10, iomsg=errmsg) npw_disk, nspinor_disk, nband_disk
1946 
1947    ! The second record: (k+G) vectors
1948    if (PRESENT(kg_k)) then
1949      write(Wfk%fh, err=10, iomsg=errmsg) kg_k(1:3,1:npw_disk)
1950    else
1951      read(Wfk%fh, err=10, iomsg=errmsg) ! kg_k(1:3,1:npw_disk)
1952    end if
1953 
1954    ! The third record: eigenvalues and occupation factors.
1955    select case (Wfk%formeig)
1956    case (0)
1957      !write(unitwf) (eigen(iband),iband=1,nband_disk),(occ(iband),iband=1,nband_disk)
1958 
1959      if (present(eig_k) .and. present(occ_k)) then
1960        write(Wfk%fh, err=10, iomsg=errmsg) eig_k, occ_k
1961      else
1962        MSG_ERROR("Not coded")
1963        write(Wfk%fh, err=10, iomsg=errmsg) ! eig_k(1:nband_disk), occ_k(1:nband_k)
1964      end if
1965 
1966      ! The wave-functions.
1967      if (present(cg_k)) then
1968        npwso = npw_disk*nspinor_disk
1969        my_bcount = 0
1970        do band=1,nband_disk
1971          if (band >= band_block(1) .and. band <= band_block(2)) then
1972            ipw = my_bcount * npwso
1973            my_bcount = my_bcount + 1
1974            write(Wfk%fh, err=10, iomsg=errmsg) cg_k(1:2,ipw+1:ipw+npwso)
1975          else
1976            MSG_ERROR("Not coded")
1977            write(Wfk%fh, err=10, iomsg=errmsg) ! cg_k(1:2,ipw+1:ipw+npwso)
1978          end if
1979        end do
1980 
1981      else
1982        MSG_ERROR("Not coded")
1983        do band=1,nband_disk
1984          write(Wfk%fh, err=10, iomsg=errmsg) ! cg_k(1:2,ipw+1:ipw+npwso)
1985        end do
1986      end if
1987 
1988    case (1)
1989      ! Write matrix of size (2*nband_k**2)
1990      ! The wave-functions.
1991      npwso = npw_disk*nspinor_disk
1992      my_bcount = 0
1993 
1994      do band=1,nband_disk
1995        base = 2*(band-1)*nband_disk
1996        write(Wfk%fh, err=10, iomsg=errmsg) eig_k(base+1:base+2*nband_disk)
1997        if (band >= band_block(1) .and. band <= band_block(2)) then
1998          ipw = my_bcount * npwso
1999          my_bcount = my_bcount + 1
2000          write(Wfk%fh, err=10, iomsg=errmsg) cg_k(1:2,ipw+1:ipw+npwso)
2001        else
2002          MSG_ERROR("Not coded")
2003          write(Wfk%fh, err=10, iomsg=errmsg) ! cg_k(1:2,ipw+1:ipw+npwso)
2004        end if
2005      end do
2006 
2007    case default
2008      MSG_ERROR("formeig != [0,1]")
2009    end select
2010 
2011    ! Reached the end of the (k,s) block. Update f90_fptr
2012    call wfk_update_f90ptr(wfk, ik_ibz, spin)
2013 
2014 #ifdef HAVE_MPI_IO
2015  case (IO_MODE_MPI)
2016    my_offset = Wfk%offset_ks(ik_ibz,spin,REC_NPW)
2017 
2018    bsize_rec(1) = 3 * xmpi_bsize_int
2019    call xmpio_write_frmarkers(Wfk%fh,my_offset,sc_mode,1,bsize_rec,ierr)
2020    ABI_CHECK(ierr==0,"ierr!=0")
2021 
2022    my_offset = Wfk%offset_ks(ik_ibz,spin,REC_NPW) + xmpio_bsize_frm
2023 
2024    call MPI_TYPE_CONTIGUOUS(3, MPI_INTEGER, recnpw_type, mpierr)
2025    ABI_CHECK_MPI(mpierr,"writing REC_NPW")
2026 
2027    call MPI_TYPE_COMMIT(recnpw_type,mpierr)
2028    ABI_CHECK_MPI(mpierr,"writing REC_NPW")
2029 
2030    call MPI_FILE_SET_VIEW(Wfk%fh,my_offset,MPI_BYTE,recnpw_type,'native',xmpio_info,mpierr)
2031    ABI_CHECK_MPI(mpierr,"writing REC_NPW")
2032 
2033    call MPI_TYPE_FREE(recnpw_type,mpierr)
2034    ABI_CHECK_MPI(mpierr,"writing REC_NPW")
2035 
2036    dims = [npw_disk, nspinor_disk, nband_disk]
2037 
2038    if (sc_mode==xmpio_collective) then
2039      call MPI_FILE_WRITE_ALL(Wfk%fh,dims,SIZE(dims),MPI_INTEGER,MPI_STATUS_IGNORE,mpierr)
2040    else if (sc_mode==xmpio_single) then
2041      call MPI_FILE_WRITE(Wfk%fh,dims,SIZE(dims),MPI_INTEGER,MPI_STATUS_IGNORE,mpierr)
2042    else
2043      MSG_ERROR("Wrong sc_mode")
2044    end if
2045    ABI_CHECK_MPI(mpierr,"writing REC_NPW")
2046 
2047    if (present(kg_k)) then
2048      my_offset = Wfk%offset_ks(ik_ibz,spin,REC_KG)
2049 
2050      bsize_rec(1) = 3 * npw_disk * xmpi_bsize_int
2051      call xmpio_write_frmarkers(Wfk%fh,my_offset,sc_mode,1,bsize_rec,ierr)
2052 
2053      my_offset = Wfk%offset_ks(ik_ibz,spin,REC_KG) + xmpio_bsize_frm
2054 
2055      call mpio_write_kg_k(Wfk%fh,my_offset,npw_disk,sc_mode,kg_k,mpierr)
2056      ABI_CHECK_MPI(mpierr,"mpio_write_kg_k")
2057    end if
2058 
2059    if (Wfk%formeig==0) then
2060 
2061      if (present(eig_k) .and. present(occ_k)) then
2062 
2063        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_EIG)
2064 
2065        bsize_rec(1) = 2 * nband_disk * xmpi_bsize_dp
2066        call xmpio_write_frmarkers(Wfk%fh,my_offset,sc_mode,1,bsize_rec,ierr)
2067 
2068        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_EIG) + xmpio_bsize_frm
2069        !
2070        ! Write both eig and occ in tmp_eigk
2071        bufsz = 2*nband_disk
2072        ABI_MALLOC(tmp_eigk, (bufsz))
2073 
2074        tmp_eigk(1:nband_disk)  = eig_k(1:nband_disk)
2075        tmp_eigk(nband_disk+1:) = occ_k(1:nband_disk)
2076 
2077        call mpio_write_eigocc_k(Wfk%fh,my_offset,nband_disk,Wfk%formeig,sc_mode,tmp_eigk,mpierr)
2078        ABI_CHECK_MPI(mpierr,"mpio_write_eigocc_k")
2079 
2080        ABI_FREE(tmp_eigk)
2081      end if
2082 
2083      if (present(cg_k)) then
2084        ABI_MALLOC(bsize_frecords, (nb_block))
2085        bsize_frecords = 2 * npw_disk * nspinor_disk * xmpi_bsize_dp
2086        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_CG) + (band_block(1)-1) * (bsize_frecords(1) + 2*xmpio_bsize_frm)
2087        call xmpio_write_frmarkers(Wfk%fh,my_offset,sc_mode,nb_block,bsize_frecords,ierr)
2088        ABI_CHECK(ierr==0,"ierr!=0")
2089        ABI_FREE(bsize_frecords)
2090 
2091        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_CG)
2092        sizes    = (/npw_disk*nspinor_disk, nband_disk/)
2093        subsizes = (/npw_disk*nspinor_disk, band_block(2)-band_block(1)+1/)
2094        bufsz = 2 * npw_disk * nspinor_disk * nb_block
2095        starts = [1, band_block(1)]
2096 
2097        call mpiotk_write_fsuba_dp2D(Wfk%fh,my_offset,sizes,subsizes,starts,bufsz,cg_k,Wfk%chunk_bsize,sc_mode,Wfk%comm,ierr)
2098        ABI_CHECK(ierr==0,"ierr!=0")
2099      end if
2100 
2101    else if (Wfk%formeig==1) then
2102 
2103      if (present(eig_k)) then
2104        types = [MPI_DOUBLE_COMPLEX,MPI_DOUBLE_COMPLEX]
2105        sizes = [nband_disk,npw_disk*nspinor_disk]
2106 
2107        call xmpio_create_fstripes(nband_disk,sizes,types,gkk_type,my_offpad,mpierr)
2108        ABI_CHECK_MPI(mpierr,"xmpio_create_fstripes")
2109 
2110        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_EIG) + my_offpad
2111 
2112        call MPI_FILE_SET_VIEW(Wfk%fh,my_offset,MPI_BYTE,gkk_type,'native',xmpio_info,mpierr)
2113        ABI_CHECK_MPI(mpierr,"SET_VIEW")
2114 
2115        call MPI_TYPE_FREE(gkk_type,mpierr)
2116        ABI_CHECK_MPI(mpierr,"")
2117 
2118        bufsz = (nband_disk**2)
2119 
2120        if (sc_mode==xmpio_collective) then
2121          call MPI_FILE_WRITE_ALL(Wfk%fh,eig_k,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
2122        else if (sc_mode==xmpio_single) then
2123          call MPI_FILE_WRITE(Wfk%fh,eig_k,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
2124        else
2125          MSG_ERROR("Wrong sc_mode")
2126        end if
2127 
2128        ABI_CHECK_MPI(mpierr,"FILE_WRITE")
2129      end if
2130 
2131      if (present(cg_k)) then
2132        ABI_CHECK(band_block(1)==1,"band_block(1) !=1 not coded")
2133 
2134        types = [MPI_DOUBLE_COMPLEX, MPI_DOUBLE_COMPLEX]
2135        sizes = [npw_disk*nspinor_disk, nband_disk]
2136 
2137        call xmpio_create_fstripes(nb_block,sizes,types,cgblock_type,my_offpad,mpierr)
2138        ABI_CHECK_MPI(mpierr,"xmpio_create_fstripes")
2139 
2140        my_offset = Wfk%offset_ks(ik_ibz,spin,REC_CG) + my_offpad
2141 
2142        call MPI_FILE_SET_VIEW(Wfk%fh,my_offset,MPI_BYTE,cgblock_type,'native',xmpio_info,mpierr)
2143        ABI_CHECK_MPI(mpierr,"SET_VIEW")
2144 
2145        call MPI_TYPE_FREE(cgblock_type,mpierr)
2146        ABI_CHECK_MPI(mpierr,"")
2147 
2148        bufsz = npw_disk * nspinor_disk * nb_block
2149        if (sc_mode==xmpio_collective) then
2150          call MPI_FILE_WRITE_ALL(Wfk%fh,cg_k,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
2151        else if (sc_mode==xmpio_single) then
2152          call MPI_FILE_WRITE(Wfk%fh,cg_k,bufsz,MPI_DOUBLE_COMPLEX,MPI_STATUS_IGNORE,mpierr)
2153        else
2154          MSG_ERROR("Wrong sc_mode")
2155        end if
2156        ABI_CHECK_MPI(mpierr,"FILE_WRITE")
2157      end if
2158 
2159    else
2160      MSG_ERROR("formeig not in [0,1]")
2161    end if
2162 #endif
2163 
2164 #ifdef HAVE_NETCDF
2165  case (IO_MODE_ETSF)
2166    if (present(kg_k)) then
2167      ! Write the reduced_coordinates_of_plane_waves for this k point.
2168      NCF_CHECK(nf90_inq_varid(wfk%fh, "reduced_coordinates_of_plane_waves", kg_varid))
2169      if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
2170        NCF_CHECK(nctk_set_collective(wfk%fh, kg_varid))
2171      end if
2172      ncerr = nf90_put_var(wfk%fh, kg_varid, kg_k, start=[1,1,ik_ibz], count=[3,npw_disk,1])
2173      NCF_CHECK_MSG(ncerr, "putting kg_k")
2174    end if
2175 
2176    ! Write eigenvalues and occupation factors.
2177    if (Wfk%formeig==0) then
2178 
2179      if (present(eig_k)) then
2180        NCF_CHECK(nf90_inq_varid(wfk%fh, "eigenvalues", eig_varid))
2181        if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
2182          NCF_CHECK(nctk_set_collective(wfk%fh, eig_varid))
2183        end if
2184        ncerr = nf90_put_var(wfk%fh, eig_varid, eig_k, start=[1,ik_ibz,spin], count=[nband_disk,1,1])
2185        NCF_CHECK_MSG(ncerr, "putting eig_k")
2186      end if
2187 
2188      if (present(occ_k)) then
2189        NCF_CHECK(nf90_inq_varid(wfk%fh, "occupations", occ_varid))
2190        if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
2191          NCF_CHECK(nctk_set_collective(wfk%fh, occ_varid))
2192        end if
2193        ncerr = nf90_put_var(wfk%fh, occ_varid, occ_k, start=[1,ik_ibz,spin], count=[nband_disk,1,1])
2194        NCF_CHECK_MSG(ncerr, "putting occ_k")
2195      end if
2196 
2197    else if (Wfk%formeig==1) then
2198      if (present(eig_k) .or. present(occ_k)) then
2199        MSG_ERROR("Don't pass eig_k or occ_k when formeig==1 and ETSF-IO")
2200      end if
2201 
2202    else
2203      MSG_ERROR("formeig != [0,1]")
2204    end if
2205 
2206    if (present(cg_k)) then
2207      ! Write the nb_block bands starting from band_block(1)
2208      ! The coefficients_of_wavefunctions on file have shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
2209      NCF_CHECK(nf90_inq_varid(wfk%fh, "coefficients_of_wavefunctions", cg_varid))
2210     if (sc_mode == xmpio_collective .and. wfk%nproc > 1) then
2211        NCF_CHECK(nctk_set_collective(wfk%fh, cg_varid))
2212      end if
2213 
2214      ncerr = nf90_put_var(wfk%fh, cg_varid, cg_k, start=[1,1,1,band_block(1),ik_ibz,spin], &
2215                           count=[2,npw_disk,wfk%nspinor,nb_block,1,1])
2216      NCF_CHECK_MSG(ncerr, "putting cg_k")
2217   end if
2218 #endif
2219 
2220  case default
2221    MSG_ERROR(sjoin('Wrong value of iomode:', itoa(Wfk%iomode)))
2222  end select
2223 
2224  DBG_EXIT("COLL")
2225 
2226  return
2227 
2228  ! Handle Fortran IO error
2229 10 continue
2230  MSG_ERROR(errmsg)
2231 
2232 end subroutine wfk_write_band_block

m_wfk/wfk_write_h1mat [ Functions ]

[ Top ] [ m_wfk ] [ Functions ]

NAME

  wfk_write_h1mat

FUNCTION

  Write all H1 matrix elements in the WFK file fname.

INPUTS

 OUTPUTS

PARENTS

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3036 subroutine wfk_write_h1mat(Wfk,sc_mode,eigen)
3037 
3038 
3039 !This section has been created automatically by the script Abilint (TD).
3040 !Do not modify the following lines by hand.
3041 #undef ABI_FUNC
3042 #define ABI_FUNC 'wfk_write_h1mat'
3043 !End of the abilint section
3044 
3045  implicit none
3046 
3047 !Arguments ------------------------------------
3048 !scalars
3049  integer,intent(in) :: sc_mode
3050  type(wfk_t),intent(inout) :: Wfk
3051 !arrays
3052  real(dp),intent(in) :: eigen(2*Wfk%mband**2*Wfk%nkpt*Wfk%nsppol)
3053 
3054 !Local variables-------------------------------
3055 !scalars
3056  integer :: spin,ik_ibz,nband_k,ptr
3057 !arrays
3058  integer,parameter :: band_block00(2)=[0,0]
3059 
3060 !************************************************************************
3061 
3062  ptr=1
3063  do spin=1,Wfk%nsppol
3064    do ik_ibz=1,Wfk%nkpt
3065      nband_k = Wfk%nband(ik_ibz,spin)
3066      call wfk_write_band_block(Wfk,band_block00,ik_ibz,spin,sc_mode,eig_k=eigen(ptr:))
3067      ptr = ptr + 2*nband_k**2
3068    end do
3069  end do
3070 
3071 end subroutine wfk_write_h1mat

m_wfkfile/wfk_compute_offsets [ Functions ]

[ Top ] [ Functions ]

NAME

  wfk_compute_offsets

FUNCTION

  Compute the offsets corresponding to the different sections of the file (G-vectors, eigenvalues, u(G).
  Needed only for Fortran-IO or MPI-IO.

PARENTS

      m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3371 subroutine wfk_compute_offsets(Wfk)
3372 
3373 
3374 !This section has been created automatically by the script Abilint (TD).
3375 !Do not modify the following lines by hand.
3376 #undef ABI_FUNC
3377 #define ABI_FUNC 'wfk_compute_offsets'
3378 !End of the abilint section
3379 
3380  implicit none
3381 
3382 !Arguments ------------------------------------
3383  type(wfk_t),intent(inout) :: Wfk
3384 
3385 !Local variables-------------------------------
3386 !scalars
3387  integer :: spin,ik_ibz,npw_k,nband_k,bsize_frm,mpi_type_frm,base !,band
3388  integer(XMPI_OFFSET_KIND) :: offset
3389 
3390 ! *************************************************************************
3391 
3392  select case (Wfk%iomode)
3393  case (IO_MODE_FORTRAN)
3394    !
3395    ! Compute record number for Fortran IO
3396    ABI_MALLOC(Wfk%recn_ks,(Wfk%nkpt,Wfk%nsppol,REC_NUM))
3397 
3398    ! We start to count the number of Fortran records from the end of the Header
3399    ! Hence recn gives the relative position from the header, it's not an absolute position.
3400    base = 0
3401    do spin=1,Wfk%nsppol
3402      do ik_ibz=1,Wfk%nkpt
3403        nband_k = Wfk%nband(ik_ibz,spin)
3404        Wfk%recn_ks(ik_ibz,spin, REC_NPW) = base + 1
3405        Wfk%recn_ks(ik_ibz,spin, REC_KG)  = base + 2
3406        Wfk%recn_ks(ik_ibz,spin, REC_EIG) = base + 3
3407        Wfk%recn_ks(ik_ibz,spin, REC_CG)  = base + 4
3408        base = Wfk%recn_ks(ik_ibz,spin,REC_CG)
3409        if (Wfk%formeig==0) then
3410          base = base + (nband_k-1)
3411        else if (Wfk%formeig==1) then
3412          base = base + 2*(nband_k-1)
3413        else
3414          MSG_ERROR("formeig != [0,1]")
3415        end if
3416      end do
3417    end do
3418    !
3419    ! Save EOF position
3420    Wfk%recn_eof = base + 1
3421 
3422  case (IO_MODE_MPI)
3423    !
3424    ! Compute offsets for MPI-IO.
3425    ABI_MALLOC(Wfk%offset_ks,(Wfk%nkpt,Wfk%nsppol,REC_NUM))
3426 
3427    bsize_frm    = xmpio_bsize_frm    ! Byte length of the Fortran record marker.
3428    mpi_type_frm = xmpio_mpi_type_frm ! MPI type of the record marker.
3429 
3430    ! The offset of the Header. TODO
3431    ! hdr_offset(Hdr)
3432    offset = Wfk%hdr_offset
3433 
3434    do spin=1,Wfk%nsppol
3435      do ik_ibz=1,Wfk%nkpt
3436        npw_k   = Wfk%Hdr%npwarr(ik_ibz)
3437        nband_k = Wfk%nband(ik_ibz,spin)
3438        !
3439        !---------------------------------------------------------------------------
3440        ! First record: npw, nspinor, nband_disk
3441        !---------------------------------------------------------------------------
3442        Wfk%offset_ks(ik_ibz,spin,REC_NPW) = offset
3443 
3444        if (Wfk%Hdr%headform>=40) then
3445          ! npw, nspinor, nband_disk
3446          offset = offset +  3*xmpi_bsize_int + 2*bsize_frm
3447        else
3448          MSG_ERROR("Old headforms < 40 are not supported")
3449        end if
3450        Wfk%offset_ks(ik_ibz,spin,REC_KG) = offset
3451        !
3452        !---------------------------------------------------------------------------
3453        ! Second record: (k+G) vectors
3454        ! kg_k(1:3,1:npw_k)
3455        !---------------------------------------------------------------------------
3456        offset = offset + 3*npw_k*xmpi_bsize_int + 2*bsize_frm
3457        Wfk%offset_ks(ik_ibz,spin,REC_EIG) = offset
3458        !
3459        !---------------------------------------------------------------------------
3460        ! Third record: eigenvalues
3461        !---------------------------------------------------------------------------
3462        if (Wfk%formeig==0) then
3463          ! eigen(1:nband_k), occ(1:nband_k)
3464          offset = offset + 2*nband_k*xmpi_bsize_dp + 2*bsize_frm
3465          Wfk%offset_ks(ik_ibz,spin,REC_CG) = offset
3466          !
3467          ! Wavefunction coefficients
3468          ! do band=1,nband_k; write(unitwf) cg_k(1:2,npw_k*nspinor); end do
3469          offset = offset + nband_k * (2*npw_k*Wfk%nspinor*xmpi_bsize_dp + 2*bsize_frm)
3470 
3471        else if (Wfk%formeig==1) then
3472          ! read(unitwf) eigen(2*nband_k)
3473          Wfk%offset_ks(ik_ibz,spin,REC_CG) = offset + 2*nband_k*xmpi_bsize_dp + 2*bsize_frm
3474 
3475          offset = offset + &
3476 &          nband_k * (2*npw_k*Wfk%nspinor*xmpi_bsize_dp + 2*bsize_frm) + &
3477 &          nband_k * (2*nband_k*xmpi_bsize_dp + 2*bsize_frm)
3478 
3479        else
3480          MSG_ERROR("Wrong formeig")
3481        end if
3482 
3483      end do ! ik_ibz
3484    end do ! spin
3485    !
3486    ! Save EOF offset
3487    Wfk%offset_eof = offset
3488 
3489    ! Check for possible wraparound errors.
3490    if (ANY(Wfk%offset_ks <= 0) .or. Wfk%offset_eof < 0) then
3491      MSG_ERROR("Found negative offset. File too large for MPI-IO!!!")
3492    end if
3493  end select
3494 
3495  if (Wfk%debug) call wfk_show_offsets(Wfk)
3496 
3497 end subroutine wfk_compute_offsets

m_wfkfile/wfk_show_offsets [ Functions ]

[ Top ] [ Functions ]

NAME

  wfk_show_offsets

FUNCTION

  Print the offsets.

PARENTS

      m_wfk

CHILDREN

      hdr_free,hdr_read_from_fname,wfk_close,wfk_open_read
      wfk_read_band_block,wrtout

SOURCE

3518 subroutine wfk_show_offsets(Wfk)
3519 
3520 
3521 !This section has been created automatically by the script Abilint (TD).
3522 !Do not modify the following lines by hand.
3523 #undef ABI_FUNC
3524 #define ABI_FUNC 'wfk_show_offsets'
3525 !End of the abilint section
3526 
3527  implicit none
3528 
3529 !Arguments ------------------------------------
3530  type(wfk_t),intent(inout) :: Wfk
3531 
3532 !Local variables-------------------------------
3533 !scalars
3534  integer :: spin,ik_ibz
3535 
3536 ! *************************************************************************
3537 
3538  select case (Wfk%iomode)
3539 
3540  case (IO_MODE_FORTRAN)
3541    write(std_out,*)"Record number relative to the header."
3542    do spin=1,Wfk%nsppol
3543      do ik_ibz=1,Wfk%nkpt
3544        write(std_out,"(a,2(i0,2x),a,4(a,i0,a))")                   &
3545 &        "(ik_ibz, spin) ",ik_ibz,spin,ch10,                       &
3546 &        "  recn(REC_NPW): ",Wfk%recn_ks(ik_ibz,spin,REC_NPW),ch10,&
3547 &        "  recn(REC_KG) : ",Wfk%recn_ks(ik_ibz,spin,REC_KG), ch10,&
3548 &        "  recn(REC_EIG): ",Wfk%recn_ks(ik_ibz,spin,REC_EIG),ch10,&
3549 &        "  recn(REC_CG) : ",Wfk%recn_ks(ik_ibz,spin,REC_CG),ch10
3550      end do
3551    end do
3552 
3553    write(std_out,"(a,i0)")"EOS position: ",Wfk%recn_eof
3554 
3555  case (IO_MODE_MPI)
3556    write(std_out,"(a,i0)")"hdr_offset ",Wfk%hdr_offset
3557 
3558    do spin=1,Wfk%nsppol
3559      do ik_ibz=1,Wfk%nkpt
3560        write(std_out,"(a,2(i0,2x),a,4(a,i0,a))")                       &
3561 &        "(ik_ibz, spin) ",ik_ibz,spin,ch10,                           &
3562 &        "  offset(REC_NPW): ",Wfk%offset_ks(ik_ibz,spin,REC_NPW),ch10,&
3563 &        "  offset(REC_KG) : ",Wfk%offset_ks(ik_ibz,spin,REC_KG), ch10,&
3564 &        "  offset(REC_EIG): ",Wfk%offset_ks(ik_ibz,spin,REC_EIG),ch10,&
3565 &        "  offset(REC_CG) : ",Wfk%offset_ks(ik_ibz,spin,REC_CG),ch10
3566      end do ! ik_ibz
3567    end do ! spin
3568    !
3569    ! Write EOF position
3570    write(std_out,"(a,i0)")"offset_eof ",Wfk%offset_eof
3571  end select
3572 
3573 end subroutine wfk_show_offsets