TABLE OF CONTENTS


ABINIT/m_rwwf [ Modules ]

[ Top ] [ Modules ]

NAME

  m_rwwf

FUNCTION

   Read/Write wavefunctions.

COPYRIGHT

  Copyright (C) 1998-2018 ABINIT group (DCA,XG,GMR,MVer,MB,MT)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_rwwf
28 
29  use defs_basis
30  use defs_abitypes
31  use m_errors
32  use m_wffile
33  use m_abicore
34  use m_xmpi
35 #if defined HAVE_MPI2
36  use mpi
37 #endif
38  use m_nctk
39 #ifdef HAVE_NETCDF
40  use netcdf
41 #endif
42 
43  use m_time,   only : timab
44 
45  implicit none
46 
47  private
48 
49 #if defined HAVE_MPI1
50  include 'mpif.h'
51 #endif

m_rwwf/readwf [ Functions ]

[ Top ] [ m_rwwf ] [ Functions ]

NAME

 readwf

FUNCTION

  This subroutine reads the block of records related to one k point, and one spin-polarization, that
  contains the wavefunctions (as well as the eigenvalues).
  The disk file unitwf should have been prepared outside of this routine, in order to read the correct records.

INPUTS

  formeig=format of the eigenvalues
    0 => vector of eigenvalues
    1 => hermitian matrix of eigenvalues
  icg=shift to be given to the location of the cg array
  ikpt=index of current k point (only needed for error message)
  isppol=spin polarization currently treated (only needed for error message)
  kg_k(3,optkg*npw)=k+g data  (only if optkg==1)
  mband=maximum number of bands (dimension of cg, eigen and occ)
  mcg=dimension of cg
  nband=number of bands actually in cg, eigen and occ
  nband_disk=number of bands on the disk file
  npw=number of plane waves (on current proc)
  nspinor=number of spinorial components of the wavefunctions (on current proc)
  option= 1 for reading cg, eigen and occ,
         -1 for reading/skipping,
         -2 for reading cg only
          3 for reading the eigenvalues only
         -4 for reading a file written with 4
  optkg= if 1 , read or write kg_k ; if 0, do not care about kg_k
  wff=struct info for wavefunction

SIDE EFFECTS

  Current kpt and spin updated
  cg(2,npw*nspinor*mband)=planewave coefficients of wavefunctions,
  eigen((2*mband)**formeig *mband)=array for holding eigenvalues (hartree)
  occ(mband)=array for holding electronic occupations

NOTES

  WARNING : occ is not read in the present status of this routine
  WARNING : skipping k-blocks is also done in the randac subroutine
  WARNING : reading the two first records is also done in the rdnpw routine

PARENTS

      rwwf

CHILDREN

      mpi_bcast,wffreadwrite_mpio,wffwritenpwrec,xderivewrecend
      xderivewrecinit,xderivewrite,xmpi_sum

SOURCE

335 subroutine readwf(cg,eigen,formeig,headform,icg,ikpt,isppol,kg_k,mband,mcg,mpi_enreg,&
336 &                 nband,nband_disk,npw,nspinor,occ,option,optkg,wff)
337 
338 
339 !This section has been created automatically by the script Abilint (TD).
340 !Do not modify the following lines by hand.
341 #undef ABI_FUNC
342 #define ABI_FUNC 'readwf'
343 !End of the abilint section
344 
345  implicit none
346 
347 !Arguments ------------------------------------
348  integer,intent(in) :: formeig,headform,icg,ikpt,isppol,mband,mcg,nband,npw
349  integer,intent(in) :: nspinor,option,optkg
350  integer,intent(inout) :: nband_disk
351  integer,intent(inout),target :: kg_k(3,optkg*npw)
352  real(dp),intent(inout),target :: cg(2,mcg),eigen((2*mband)**formeig*mband),occ(mband)
353  type(MPI_type),intent(in) :: mpi_enreg
354  type(wffile_type),intent(inout) :: wff
355 
356 !Local variables-------------------------------
357  integer :: iband,indxx,ios,ipw,ispinor_index,nband1,ncid_hdr
358  integer :: npw1,npwso,npwso1,npwsotot,npwtot,nrec,nspinor1,nspinortot,unitwf,use_f90
359  integer :: band1,band2,ierr
360  character(len=500) :: msg
361  character(len=fnlen) :: fname
362  integer :: ikpt_this_proc,ispinor
363  integer,allocatable :: ind_cg_mpi_to_seq(:)
364 #ifdef HAVE_NETCDF
365  integer :: kg_varid,eig_varid,occ_varid,cg_varid,h1_varid,mband_varid,ncerr,ii,mband_file
366  integer,allocatable :: gall(:,:)
367  real(dp),allocatable :: cg_all(:,:,:),h1mat(:,:,:)
368 #endif
369 
370 ! *********************************************************************
371 
372  !write(std_out,*)"rwwf with ikpt, isppol, option, etsf",ikpt, isppol, option, wff%iomode == IO_MODE_ETSF
373 
374 !Check the options
375  if ( ALL(option /= [1,3,-1,-2,-4])) then
376    write(msg,'(a,i0,a)')'The argument option should be -4, -2, -1, 1 or 3.  However, option=',option,'.'
377    MSG_BUG(msg)
378  end if
379 
380  npwtot=npw; npwso=npw*nspinor
381  npwsotot=npwso
382  nspinortot=min(2,(1+mpi_enreg%paral_spinor)*nspinor)
383 
384  unitwf=wff%unwff
385  ncid_hdr=unitwf
386 
387  use_f90=0
388  if (wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) use_f90=1
389 
390  if (option==1.or.option==-2) then
391 
392    ! Compute mapping my_gtable --> sequential gtable.
393    if (any(wff%iomode==[IO_MODE_MPI, IO_MODE_ETSF]).and.nband>0) then
394      call xmpi_sum(npwsotot,wff%spaceComm_mpiio,ios)
395      npwtot=npwsotot/nspinortot
396      ABI_ALLOCATE(ind_cg_mpi_to_seq,(npwso))
397      if (allocated(mpi_enreg%my_kgtab)) then
398        ikpt_this_proc=mpi_enreg%my_kpttab(ikpt)
399        if ( ikpt_this_proc <= 0  ) then
400          MSG_BUG("rwwf: ikpt_this_proc <= 0")
401        end if
402        do ispinor=1,nspinor
403          ispinor_index=ispinor
404          if (mpi_enreg%nproc_spinor>1) ispinor_index=mpi_enreg%me_spinor + 1
405          ind_cg_mpi_to_seq(1+npw*(ispinor-1):npw*ispinor)=npwtot*(ispinor_index-1) &
406 &         + mpi_enreg%my_kgtab(1:npw,ikpt_this_proc)
407        end do
408      else
409        ind_cg_mpi_to_seq(1:npwso) = (/(ipw,ipw=1,npwso)/)
410      end if
411    end if
412  end if
413 
414 !---------------------------------------------------------------------------
415 !Read the first record: npw, nspinor, nband_disk
416 !---------------------------------------------------------------------------
417 
418  if (headform>=40.or.headform==0) then ! headform==0 refers to the current headform
419    call WffReadNpwRec(ios,ikpt,isppol,nband_disk,npw1,nspinor1,wff)
420    npwso1=npw1*nspinor1
421    if(ios/=0)then
422      inquire (unit=unitwf, NAME=fname)
423      write(msg,'(3a,i4,2a,i4,4a)') &
424 &     'Reading option of rwwf. Trying to read',ch10,&
425 &     'the (npw,nspinor,nband) record of a wf file, unit=',unitwf,ch10,&
426 &     'gave iostat=',ios,'. Your file is likely not correct.',ch10,&
427 &     'Action: check your input wf file:',trim(fname)
428      MSG_ERROR(msg)
429    end if
430 
431  else
432    ! Old format
433    if(use_f90==1)then
434      read (unitwf,iostat=ios) npwso1,nband_disk
435    else if(wff%iomode==IO_MODE_MPI)then
436      call xderiveRRecInit(wff,ios)
437      call xderiveRead(wff,npwso1,ios)
438      call xderiveRead(wff,nband_disk,ios)
439      call xderiveRRecEnd(wff,ios)
440    end if
441    if(ios/=0)then
442      inquire (unit=unitwf, NAME=fname)
443      write(msg,'(3a,i4,2a,i4,4a)') &
444 &     'Reading option of rwwf. Trying to read',ch10,&
445 &     'the (npw,nband) record of a wf file with headform <40 , unit=',unitwf,ch10,&
446 &     'gave iostat=',ios,'. Your file is likely not correct.',ch10,&
447 &     'Action: check your input wf file:',trim(fname)
448      MSG_ERROR(msg)
449    end if
450  end if ! headform
451 
452  if (option==1.or.option==-2) then !  Will read the wavefunction and/or kg data, so check npw and nspinor
453 
454    if (headform>=40.or.headform==0) then ! New format. headform==0 refers to the current headform
455      if (npwtot/=npw1) then
456        write(msg,'(3a,i0,a,i0,a)') &
457 &       'Reading option of rwwf. One should have npwtot=npw1',ch10,&
458 &       'However, npwtot= ',npwtot,', and npw1= ',npw1,'.'
459        MSG_BUG(msg)
460      end if
461      if(nspinortot/=nspinor1)then
462        write(msg,'(3a,i0,a,i0,a)') &
463 &       'Reading option of rwwf. One should have nspinor=nspinor1',ch10,&
464 &       'However, nspinortot= ',nspinortot,', and nspinor1= ',nspinor1,'.'
465        MSG_BUG(msg)
466      end if
467    else ! Treat the Old format.
468      if(npwsotot/=npwso1)then
469        write(msg,'(3a,i0,a,i0,a)') &
470 &       'Reading option of rwwf. One should have npwso=npwso1',ch10,&
471 &       'However, npwsotot= ',npwsotot,', and npwso1= ',npwso1,'.'
472        MSG_BUG(msg)
473      end if
474    end if ! headform
475 !
476  end if ! option==1.or.option==2
477 
478 !---------------------------------------------------------------------------
479 !Read the second record: (k+G) vectors
480 !---------------------------------------------------------------------------
481 
482  if (headform>=40.or.headform==0) then ! headform==0 refers to the current headform
483 
484    if ((option==1.or.option==-2.or.option==3).and.optkg/=0 )then
485 
486      if(use_f90==1)then
487        read(unitwf,iostat=ios) kg_k(1:3,1:npw)
488 
489      else if(wff%iomode==IO_MODE_MPI)then
490        call xderiveRRecInit(wff,ios)
491        if (allocated(mpi_enreg%my_kgtab)) then
492          ikpt_this_proc=mpi_enreg%my_kpttab(ikpt)
493          call xderiveRead(wff,kg_k(1:3,1:npw),3,npw,wff%spaceComm_mpiio,&
494 &         mpi_enreg%my_kgtab(1:npw,ikpt_this_proc),ios)
495        else
496 !        MG The call below uses MPI_SCAN but here we want to read the full set of G.
497          call xderiveRead(wff,kg_k(1:3,1:npw),3,npw,wff%spaceComm_mpiio,ios)
498          call xmpi_read_int2d(wff,kg_k(1:3,1:npw),wff%spaceComm_mpiio,xmpio_collective,ios)
499        end if
500        call xderiveRRecEnd(wff,ios)
501 
502 #ifdef HAVE_NETCDF
503      else if (wff%iomode == IO_MODE_ETSF) then
504        ! Read reduced_coordinates_of_plane_waves for this k point (npw1 is npw_disk).
505        ! TODO: spinor parallelism
506        NCF_CHECK(nf90_inq_varid(wff%unwff, "reduced_coordinates_of_plane_waves", kg_varid))
507        if (npw == npw1) then
508          ncerr = nf90_get_var(wff%unwff, kg_varid, kg_k, start=[1,1,ikpt], count=[3,npw1,1])
509          NCF_CHECK_MSG(ncerr, "getting kg_k")
510        else
511          write(std_out,*)"ETSF Reading distributed kg_k"
512          ABI_MALLOC(gall, (3, npw1))
513          ncerr = nf90_get_var(wff%unwff, kg_varid, gall, start=[1,1,ikpt], count=[3,npw1,1])
514          NCF_CHECK_MSG(ncerr, "getting kg_k")
515          do ipw=1,npw
516            kg_k(:,ipw) = gall(:,ind_cg_mpi_to_seq(ipw))
517          end do
518          ABI_FREE(gall)
519        end if
520 #endif
521      end if
522 
523    else ! option
524      call WffReadSkipRec(ios,1,wff) ! Skip the record
525    end if
526 
527    if(ios/=0)then
528      inquire (unit=unitwf, NAME=fname)
529      write(msg,'(3a,i4,2a,i4,4a)')  &
530 &     'Reading option of rwwf. Trying to read',ch10,&
531 &     'the k+g record of a wf file, unit=',unitwf,ch10,&
532 &     'gave iostat=',ios,'. Your file is likely not correct.',ch10,&
533 &     'Action: check your input wf file:',trim(fname)
534      MSG_ERROR(msg)
535    end if
536  end if ! headform
537 
538 !---------------------------------------------------------------------------
539 !Read the third record: eigenvalues
540 !---------------------------------------------------------------------------
541 !The reading of occ should be enabled, BUT taking into account
542 !of headform of the disk file : occ was NOT present in the disk files with headform=22
543 
544  nband1 = min(nband,nband_disk)
545 
546 !===== Case formeig=0: read eigenvalues =====
547  if (formeig==0) then
548 
549    if (option==1.or.option==3.or.option==-4) then
550      if(use_f90==1)then
551        read (unitwf,iostat=ios) eigen(1:nband1)
552      else if(wff%iomode==IO_MODE_MPI)then
553        call xderiveRRecInit(wff,ios)
554        call xderiveRead(wff,eigen,nband1,xmpi_comm_self,ios)
555        call xderiveRRecEnd(wff,ios)
556      end if
557    else
558      call WffReadSkipRec(ios,1,wff)
559    end if ! option
560 
561    if(ios/=0)then
562      inquire (unit=unitwf, NAME=fname)
563      write(msg,'(3a,i4,2a,i4,4a)') &
564 &     'Reading option of rwwf. Trying to read',ch10,&
565 &     'an eigenvalue record of a wf file, unit=',unitwf,ch10,&
566 &     'gave iostat=',ios,'. Your file is likely not correct.',ch10,&
567 &     'Action: check your input wf file.',trim(fname)
568      MSG_ERROR(msg)
569    end if
570 
571 #ifdef HAVE_NETCDF
572    if (wff%iomode == IO_MODE_ETSF) then
573      ! get eigenvalues and occupations
574      NCF_CHECK(nf90_inq_varid(wff%unwff, "eigenvalues", eig_varid))
575      ncerr = nf90_get_var(wff%unwff, eig_varid, eigen, start=[1,ikpt,isppol], count=[nband1,1,1])
576      NCF_CHECK_MSG(ncerr, "getting eig_k")
577 
578      NCF_CHECK(nf90_inq_varid(wff%unwff, "occupations", occ_varid))
579      ncerr = nf90_get_var(wff%unwff, occ_varid, occ, start=[1,ikpt,isppol], count=[nband1,1,1])
580      NCF_CHECK_MSG(ncerr, "getting occ_k")
581    end if
582 #endif
583 
584 !  ===== Case formeig=1: read matrix of eigenvalues =====
585 !  Will be written later (together with wave-functions)
586  else if(formeig==1)then
587  end if ! formeig
588 
589 !---------------------------------------------------------------------------
590 !Read the wave-function coefficients
591 !---------------------------------------------------------------------------
592 
593 !Select bands
594  nband1=min(nband,nband_disk)
595  if(nband1>0.and.option/=-1)then
596 
597 !  ===== Case formeig=0: read only wave-functions =====
598    if (formeig==0) then
599 
600      if (option==1.or.option==-2) then
601 
602        if (use_f90==1) then
603          do iband=1,nband1
604            ipw=(iband-1)*npwso+icg
605            read(unitwf,iostat=ios) cg(1:2,ipw+1:ipw+npwso)
606            if (ios/=0) exit
607          end do
608        else if (wff%iomode==IO_MODE_MPI) then
609          call WffReadWrite_mpio(wff,1,cg,mcg,icg,nband1,npwso,npwsotot,ind_cg_mpi_to_seq,ios)
610        end if
611 
612      else if (option/=-4) then
613        do iband=1,nband1
614          call WffReadSkipRec(ios,1,wff) ! Skip the record
615        end do
616      end if ! option
617 
618      if(ios/=0)then
619        inquire (unit=unitwf, NAME=fname)
620        write(msg,'(3a,i4,2a,i4,4a)') &
621 &       'Reading option of rwwf. Trying to read',ch10,&
622 &       'a RF wf record of a wf file, unit=',unitwf,ch10,&
623 &       'gave iostat=',ios,'. Your file is likely not correct.',ch10,&
624 &       'Action: check your input wf file.',trim(fname)
625        MSG_ERROR(msg)
626      end if
627 
628 #ifdef HAVE_NETCDF
629      if (wff%iomode == IO_MODE_ETSF) then
630        ! The coefficients_of_wavefunctions on file have shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
631        NCF_CHECK(nf90_inq_varid(wff%unwff, "coefficients_of_wavefunctions", cg_varid))
632        if (npw == npw1) then
633          ncerr = nf90_get_var(wff%unwff, cg_varid, cg(:, icg+1:), start=[1,1,1,1,ikpt,isppol], &
634          count=[2,npw,nspinor,nband1,1,1])
635          NCF_CHECK_MSG(ncerr, "getting cg_k")
636        else
637          write(std_out,*)"ETSF Reading distributed cg"
638          ABI_STAT_MALLOC(cg_all, (2, npw1*nspinor, nband1), ierr)
639          ABI_CHECK(ierr==0, "oom cg_all")
640          ncerr = nf90_get_var(wff%unwff, cg_varid, cg_all, start=[1,1,1,1,ikpt,isppol], &
641          count=[2,npw1,nspinor,nband1,1,1])
642          NCF_CHECK_MSG(ncerr, "getting cg_k")
643          ii = icg
644          do iband=1,nband1
645            do ipw=1,npw
646              ii = ii + 1
647              cg(:,ii) = cg_all(:,ind_cg_mpi_to_seq(ipw),iband)
648            end do
649          end do
650          ABI_FREE(cg_all)
651        end if
652      end if
653 #endif
654 
655 !    ===== Case formeig=1: read eigenvalues, occupations and wave-functions =====
656    else if(formeig==1)then
657      ! write(std_out,*)"nband1",nband1
658      ! ABI_CHECK(nband1==nband_disk,"nband != nband_disk")
659 
660      if (wff%iomode /= IO_MODE_ETSF) then
661        indxx=0
662        do iband=1,nband1
663 
664          if (option==1.or.option==3.or.option==-4) then
665            if(use_f90==1)then
666              read (unitwf,iostat=ios) eigen(1+indxx:2*nband1+indxx)
667            else if(wff%iomode==IO_MODE_MPI)then
668 !            Should use an access with a "view"
669              call xderiveRRecInit(wff,ios)
670              call xderiveRead(wff,eigen(1+indxx:2*nband1+indxx),2*nband1,xmpi_comm_self,ios)
671              call xderiveRRecEnd(wff,ios)
672            end if
673            indxx=indxx+2*nband1
674          else
675            call WffReadSkipRec(ios,1,wff) ! Skip the record
676          end if
677 
678          if(ios/=0)then
679            inquire (unit=unitwf, NAME=fname)
680            write(msg,'(3a,i4,2a,i4,4a)') &
681 &           'Reading option of rwwf. Trying to read',ch10,&
682 &           'a RF eigenvalue record of a wf file, unit=',unitwf,ch10,&
683 &           'gave iostat=',ios,'. Your file is likely not correct.',ch10,&
684 &           'Action: check your input wf file.',trim(fname)
685            MSG_ERROR(msg)
686          end if
687 
688          if(option==1.or.option==-2)then
689            ipw=(iband-1)*npwso+icg
690            if(use_f90==1)then
691              ipw=(iband-1)*npwso+icg
692              read(unitwf,iostat=ios) cg(1:2,ipw+1:ipw+npwso)
693            else if(wff%iomode==IO_MODE_MPI)then
694 !            Should use an access with a "view"
695              call xderiveRRecInit(wff,ios)
696              call xderiveRead(wff,cg(1:2,ipw+1:ipw+npwso),2,npwso,wff%spaceComm_mpiio,ind_cg_mpi_to_seq,ios)
697              call xderiveRRecEnd(wff,ios)
698            end if
699          else if (option/=-4) then
700            call WffReadSkipRec(ios,1,wff) ! Skip the record
701          end if ! option
702 
703          if(ios/=0)then
704            inquire (unit=unitwf, NAME=fname)
705            write(msg,'(3a,i4,2a,i4,4a)') &
706 &           'Reading option of rwwf. Trying to read',ch10,&
707 &           'a RF wf record of a wf file, unit=',unitwf,ch10,&
708 &           'gave iostat=',ios,'. Your file is likely not correct.',ch10,&
709 &           'Action: check your input wf file.',trim(fname)
710            MSG_ERROR(msg)
711          end if
712        end do ! iband
713 
714      else
715 #ifdef HAVE_NETCDF
716         ! ETSF-IO
717        if (any(option == [1, 3, -4])) then
718           ! Read eigen. Remember that the matrix on file has shape [2, mband, mband, nkpt, nspin]
719           ! whereas the eigen array used in Abinit is packed. Read the full matrix first, then pack data.
720          NCF_CHECK(nf90_inq_dimid(wff%unwff, "max_number_of_states", mband_varid))
721          NCF_CHECK(nf90_inquire_dimension(wff%unwff, mband_varid, len=mband_file))
722          h1_varid = nctk_idname(wff%unwff, "h1_matrix_elements")
723 
724          ABI_MALLOC(h1mat, (2, mband_file, mband_file))
725          ncerr = nf90_get_var(wff%unwff, h1_varid, h1mat, start=[1,1,1,ikpt,isppol])
726          NCF_CHECK_MSG(ncerr, "getting h1_matrix_elements")
727 
728          indxx=1
729          do band2=1,nband1
730            do band1=1,nband1
731              eigen(indxx:indxx+1) = h1mat(:,band1,band2)
732              indxx = indxx + 2
733            end do
734          end do
735          ABI_FREE(h1mat)
736        end if
737 
738        if (any(option == [1, -2])) then
739           ! Read wavefunctions.
740           ! The coefficients_of_wavefunctions on file have shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
741          ABI_CHECK(npw == npw1, "npw != npw1 not coded")
742 
743          cg_varid = nctk_idname(wff%unwff, "coefficients_of_wavefunctions")
744          ncerr = nf90_get_var(wff%unwff, cg_varid, cg(:, icg+1:), start=[1,1,1,1,ikpt,isppol], &
745          count=[2,npw,nspinor,nband1,1,1])
746          NCF_CHECK_MSG(ncerr, "getting cg_k")
747        end if
748 #endif
749      end if
750 
751    end if ! formeig == 1
752 
753  end if ! nband >0
754 
755 !If fewer than all bands were read wind disk file forward to end of bands for this k point.
756 !Will have to fill the non-filled bands outside of this routine ...
757  if (nband<nband_disk .or. option==-1) then
758    nrec=(formeig+1)*(nband_disk-nband)
759    if(option==-1)nrec=(formeig+1)*nband_disk
760    call WffReadSkipRec(ios,nrec,wff)
761  end if
762 
763 !---------------------------------------------------------------------------
764 ! Free memory
765 !---------------------------------------------------------------------------
766  if (allocated(ind_cg_mpi_to_seq)) then
767    ABI_DEALLOCATE(ind_cg_mpi_to_seq)
768  end if
769 
770 end subroutine readwf

m_rwwf/rwwf [ Functions ]

[ Top ] [ m_rwwf ] [ Functions ]

NAME

 rwwf

FUNCTION

  This subroutine reads (different options) or write (option=2) the block of records
  related to one k point, and one spin-polarization, that
  contains the wavefunctions (as well as the eigenvalues and occupations).
  If called with option -1, the records will be skipped.
  If called with option -2, only the wavefunctions are read.
  The disk file unitwf should have been prepared
  outside of this routine, in order to read or write the correct records.

INPUTS

  formeig=format of the eigenvalues
     0 => vector of eigenvalues
     1 => hermitian matrix of eigenvalues
  headform=format of the header of the wf file, also governing the k block format
    in case headform=0, use the default (current) format and headform
  icg=shift to be given to the location of the cg array
  ikpt=index of current k point (only needed for error message)
  isppol=spin polarization currently treated (only needed for error message)
  mband=maximum number of bands (dimension of cg, eigen and occ)
  mcg=dimention of cg
  nband=number of bands actually in cg, eigen and occ
   (if writing mode : must be larger or equal to nband_disk, only nband_disk bands are written ;
    if reading mode : can be equal, larger or smaller than nband_disk, but
     cg, eigen and occ will not be completely filled if nband>nband_disk)
  nband_disk=number of bands on the disk file
  npw=number of plane waves
  nspinor=number of spinorial components of the wavefunctions (on current proc)
  option= 2 for writing cg, eigen and occ,
          1 for reading cg and eigen,
         -1 for reading/skipping,
         -2 for reading cg only
          3 for reading the eigenvalues only
          4 for writing a file containing only eigenvalues and occupations (need to be read with option 4)
          5 for writing a file containing only eigenvalues and occupations, that can however be read with option 3)
         -4 for reading a file written with 4
                 (different from 3 which reads a normal option 2 file)
  optkg= if 1 , read or write kg_k ; if 0, do not care about kg_k
  tim_rwwf=timing code of the calling routine (set to 0 if not attributed)
  wff=struct info for wavefunction
   | unitwf=unit number for wavefunction

SIDE EFFECTS

  cg(2,npw*nspinor*mband)=planewave coefficients of wavefunctions,
    input if option=2; output if option=1 or -2
  eigen((2*mband)**formeig *mband)=array for holding eigenvalues (hartree)
    input if option=2 or 4 or 5; output if option=1
  kg_k(3,optkg*npw)=k+g data  (only if optkg==1)
    input if option=2; output if option=1 or -2
  nband_disk=number of bands on disk
    input if option=2 or 4 or 5; output in the other cases
  occ(mband)=array for holding eigenvalues (hartree)
    input if option=2 or 4 or 5; output if option=1
    no meaning if frmeig/=0

NOTES

  WARNING : occ is not read in the present status of this routine
  WARNING : skipping k-blocks is also done in the randac subroutine
  WARNING : reading the two first records is also done in the rdnpw routine
  WARNING : writing the two first records is also done in the dfpt_vtowfk routine

 TODO!     if (mpi_enreg%flag_ind_kg_mpi_to_seq==1 ) then
  Some arguments are contained in the wff datastructure, and should be eliminated.
  option 3 should be called -3 (reading -> negative option) and others (-1,1) re-shuffled.

PARENTS

      WffReadEigK,WffReadSkipK,initwf,m_iowf,m_wfk,newkpt,uderiv

CHILDREN

      mpi_bcast,wffreadwrite_mpio,wffwritenpwrec,xderivewrecend
      xderivewrecinit,xderivewrite,xmpi_sum

SOURCE

139 subroutine rwwf(cg,eigen,formeig,headform,icg,ikpt,isppol,kg_k,mband,mcg,mpi_enreg,&
140 &               nband,nband_disk,npw,nspinor,occ,option,optkg,tim_rwwf,wff)
141 
142 
143 !This section has been created automatically by the script Abilint (TD).
144 !Do not modify the following lines by hand.
145 #undef ABI_FUNC
146 #define ABI_FUNC 'rwwf'
147 !End of the abilint section
148 
149  implicit none
150 
151 !Arguments ------------------------------------
152  integer,intent(in) :: formeig,headform,icg,ikpt,isppol,mband,mcg,nband,npw
153  integer,intent(inout) :: nband_disk
154  integer,intent(in) :: nspinor,option,optkg,tim_rwwf
155  integer,intent(inout),target :: kg_k(3,optkg*npw)
156  real(dp),intent(inout),target :: cg(2,mcg),eigen((2*mband)**formeig*mband),occ(mband)
157  type(wffile_type),intent(inout) :: wff
158  type(MPI_type), intent(in) :: mpi_enreg
159 
160 !Local variables-------------------------------
161 !scalars
162  character(len=500) :: msg
163 !arrays
164  real(dp) :: tsec(2)
165 
166 ! *************************************************************************
167 
168  call timab(270+tim_rwwf,1,tsec)
169 
170 !Might check that icg+npw*nband*nspinor is smaller than mcg
171 
172 !Check that nband is smaller than mband, if one will not skip the records.
173  if (nband>mband .and. option/=-1)then
174    write(msg,'(a,i0,a,i0,a)')' One should have nband<=mband. However, nband=',nband,', and mband=',mband,'.'
175    MSG_BUG(msg)
176  end if
177 
178 !Check that formeig is 0 or 1.
179  if ( ALL(formeig/=(/0,1/)) ) then
180    write(msg,'(a,i0,a)')' The argument formeig should be 0 or 1. However, formeig=',formeig,'.'
181    MSG_BUG(msg)
182  end if
183 
184 !Check the value of option
185  if ( ALL( option /= (/1,2,3,4,5,-1,-2,-4/) )) then
186    write(msg,'(a,i0,a)')' The argument option should be between 1 and 5, or -1, -2, -4. However, option=',option,'.'
187    MSG_BUG(msg)
188  end if
189 
190  if (option/=2.and.option/=4 .and. option/=5 ) then ! read
191    call readwf(cg,eigen,formeig,headform,icg,ikpt,isppol,kg_k,mband,mcg,mpi_enreg,&
192 &   nband,nband_disk,npw,nspinor,occ,option,optkg,wff)
193  else                                ! write
194    call writewf(cg,eigen,formeig,icg,ikpt,isppol,kg_k,mband,mcg,mpi_enreg,&
195 &   nband,nband_disk,npw,nspinor,occ,option,optkg,wff)
196  end if
197 
198  call timab(270+tim_rwwf,2,tsec)
199 
200 end subroutine rwwf

m_rwwf/WffReadSkipK [ Functions ]

[ Top ] [ m_rwwf ] [ Functions ]

NAME

 WffReadSkipK

FUNCTION

  (Wavefunction file, read action : skip one k-point blok)
  This subroutine skips the block of records
  related to one k point, and one spin-polarization, that
  contains the wavefunctions as well as the eigenvalues and occupations,
  in a wavefunction file that has been already initialized.

INPUTS

  formeig=format of the eigenvalues
   0 => vector of eigenvalues (for Ground-State files)
   1 => hermitian matrix of eigenvalues (for Response-Function files)
  headform=format of the header of the wf file, also governing the k block format
   in case headform=0, use the default (current) format and headform
  ikpt=index of current k point (only needed for error message)
  isppol=spin polarization currently treated (only needed for error message)
  mpi_enreg=information about MPI parallelization
  wff=structured info for wavefunction file

OUTPUT

NOTES

PARENTS

      initwf,newkpt,wfsinp

CHILDREN

      rwwf

SOURCE

237 subroutine WffReadSkipK(formeig,headform,ikpt,isppol,mpi_enreg,wff)
238 
239 
240 !This section has been created automatically by the script Abilint (TD).
241 !Do not modify the following lines by hand.
242 #undef ABI_FUNC
243 #define ABI_FUNC 'WffReadSkipK'
244 !End of the abilint section
245 
246  implicit none
247 
248 !Arguments ------------------------------------
249 !scalars
250  integer,intent(in) :: formeig,headform,ikpt,isppol
251  type(MPI_type),intent(in) :: mpi_enreg
252  type(wffile_type),intent(inout) :: wff
253 
254 !Local variables-------------------------------
255 !scalars
256  integer :: icg,mband,mcg,nband,nband_disk,option,optkg,tim_rwwf
257  integer,parameter :: nspinor1=1,npw1=1
258 !arrays
259  integer,allocatable :: kg_dum(:,:)
260  real(dp) :: cg_dum(2,1),occ_dum(1)
261  real(dp),allocatable :: eig_dum(:)
262 
263 ! *************************************************************************
264 
265  ! No need to skip if netcdf
266  if (wff%iomode == IO_MODE_ETSF) return
267 
268  option=-1
269  tim_rwwf=0 ; mcg=1 ; mband=1 ; icg=0 ; optkg=0 ; nband=0
270  ABI_ALLOCATE(eig_dum,(2**formeig))
271  ABI_ALLOCATE(kg_dum,(3,optkg*npw1))
272 
273  call rwwf(cg_dum,eig_dum,formeig,headform,icg,ikpt,isppol,kg_dum,mband,mcg,mpi_enreg,nband,&
274 & nband_disk,npw1,nspinor1,occ_dum,option,optkg,tim_rwwf,wff)
275 
276  ABI_DEALLOCATE(eig_dum)
277  ABI_DEALLOCATE(kg_dum)
278 
279 end subroutine WffReadSkipK

m_rwwf/writewf [ Functions ]

[ Top ] [ m_rwwf ] [ Functions ]

NAME

 writewf

FUNCTION

  This subroutine writes the block of records
  related to one k point, and one spin-polarization, that
  contains the wavefunctions (as well as the eigenvalues and occupations).
  The disk file unitwf should have been prepared
  outside of this routine, in order to write the correct records.

INPUTS

  cg(2,npw*nspinor*mband)=planewave coefficients of wavefunctions,
  eigen((2*mband)**formeig *mband)=array for holding eigenvalues (hartree)
  formeig=format of the eigenvalues
   0 => vector of eigenvalues
   1 => hermitian matrix of eigenvalues
  icg=shift to be given to the location of the cg array
  ikpt=index of current k point (only needed for error message)
  isppol=spin polarization currently treated (only needed for error message)
  kg_k(3,optkg*npw)=k+g data  (only if optkg==1)
  mband=maximum number of bands (dimension of cg, eigen and occ)
  mcg=dimension of cg
  nband=number of bands actually in cg, eigen and occ
   (must be larger or equal to nband_disk, only nband_disk bands are written)
  nband_disk=number of bands on the disk file
  npw=number of plane waves
  nspinor=number of spinorial components of the wavefunctions (on current proc)
  occ(mband)=array for holding electronic occupations
  option= 2 for writing cg, eigen and occ,
          4 for writing a file containing only eigenvalues and occupations
  optkg= if 1 , read or write kg_k ; if 0, do not care about kg_k
  wff=struct info for wavefunction

OUTPUT

 (none, only writing)

SIDE EFFECTS

NOTES

  WARNING : skipping k-blocks is also done in the randac subroutine
  WARNING : writing the two first records is also done in the dfpt_vtowfk routine

PARENTS

      rwwf

CHILDREN

      mpi_bcast,wffreadwrite_mpio,wffwritenpwrec,xderivewrecend
      xderivewrecinit,xderivewrite,xmpi_sum

SOURCE

 827 subroutine writewf(cg,eigen,formeig,icg,ikpt,isppol,kg_k,mband,mcg,mpi_enreg,&
 828 &                  nband,nband_disk,npw,nspinor,occ,option,optkg,wff)
 829 
 830 
 831 !This section has been created automatically by the script Abilint (TD).
 832 !Do not modify the following lines by hand.
 833 #undef ABI_FUNC
 834 #define ABI_FUNC 'writewf'
 835 !End of the abilint section
 836 
 837  implicit none
 838 
 839 !Arguments ------------------------------------
 840  integer,intent(in) :: formeig,icg,ikpt,isppol,mband,mcg,nband,nband_disk,npw,nspinor,option,optkg
 841  integer,intent(in),target :: kg_k(3,optkg*npw)
 842  real(dp),intent(in),target :: cg(2,mcg),eigen((2*mband)**formeig*mband),occ(mband)
 843  type(MPI_type),intent(in) :: mpi_enreg
 844  type(wffile_type),intent(inout) :: wff
 845 
 846 !Local variables-------------------------------
 847  integer :: iband,ii,ios,ipw,ispinor_index,nband2,ncid_hdr,npwso,npwsotot,npwtot,nspinortot
 848  integer :: unitwf,use_f90
 849  character(len=500) :: msg
 850  integer :: ikpt_this_proc,ispinor,me_cart_3d
 851  integer,allocatable :: ind_cg_mpi_to_seq(:)
 852  real(dp),ABI_CONTIGUOUS pointer :: cg_ptr(:,:)
 853 #ifdef HAVE_NETCDF
 854  integer :: kg_varid,eig_varid,occ_varid,cg_varid,ncerr
 855  character(len=nctk_slen) :: kdep
 856 #endif
 857 
 858 ! *********************************************************************
 859 
 860 !Check the options
 861  if ( ALL(option /= (/2,4,5/)) ) then
 862    write(msg,'(a,i0)')' The argument option should be 2, 4 or 5. However, option=',option
 863    MSG_BUG(msg)
 864  end if
 865 
 866  if(wff%iomode==IO_MODE_MPI)then
 867    if(option==5)then
 868      write(msg,'(a,i0)')' With MPI-IO activated, the argument option should be 2 or 4. However, option=',option
 869      MSG_BUG(msg)
 870    end if
 871  end if
 872 
 873  if (wff%iomode==IO_MODE_ETSF) then
 874    ! outkss still calls this routine!
 875    MSG_WARNING("You should use outwff or ncwrite_cg to write data in netcdf format!")
 876  end if
 877 
 878 !Check that nband_disk is not larger than nband (only for writing)
 879  if (nband<nband_disk) then
 880    write(msg,'(3a,i5,a,i5,a)') &
 881 &   'Writing option of rwwf. One should have nband<=nband_disk',ch10,&
 882 &   'However, nband= ',nband,', and nband_disk= ',nband_disk,'.'
 883    MSG_BUG(msg)
 884  end if
 885 
 886  npwtot=npw; npwso=npw*nspinor
 887  unitwf=wff%unwff; ncid_hdr=unitwf
 888  npwsotot=npwso
 889  nspinortot=min(2,(1+mpi_enreg%paral_spinor)*nspinor)
 890 
 891  use_f90=0
 892  if (wff%iomode==IO_MODE_FORTRAN.or.(wff%iomode ==IO_MODE_FORTRAN_MASTER.and.wff%master==wff%me)) use_f90=1
 893 
 894  if (wff%iomode==IO_MODE_MPI) then
 895 
 896    call xmpi_sum(npwsotot,wff%spaceComm_mpiio,ios)
 897    npwtot=npwsotot/nspinortot
 898 
 899    if (option/=4) then
 900      ABI_ALLOCATE(ind_cg_mpi_to_seq,(npwso))
 901      if (allocated(mpi_enreg%my_kgtab)) then
 902        ikpt_this_proc=mpi_enreg%my_kpttab(ikpt)
 903        do ispinor=1,nspinor
 904          ispinor_index=ispinor
 905          if (mpi_enreg%nproc_spinor > 1) ispinor_index = mpi_enreg%me_spinor + 1
 906          ind_cg_mpi_to_seq(1+npw*(ispinor-1):npw*ispinor)=npwtot*(ispinor_index-1) &
 907 &         + mpi_enreg%my_kgtab(1:npw,ikpt_this_proc)
 908        end do
 909      else
 910        ind_cg_mpi_to_seq(1:npwso) = (/(ipw,ipw=1,npwso)/)
 911      end if
 912      !write(std_out,*)"MPI-IO ind_cg_mpi_to_seq", ind_cg_mpi_to_seq(1:5)
 913    end if
 914  end if
 915 
 916 !---------------------------------------------------------------------------
 917 !Write the first record: npw, nspinor, nband_disk
 918 !---------------------------------------------------------------------------
 919 !Not modified for netCDF: no need to add writing of nband_disk,npw,nspinor
 920 
 921  call WffWriteNpwRec(ios,nband_disk,npwtot,nspinortot,wff,opt_paral=2)
 922 
 923 !---------------------------------------------------------------------------
 924 !Write the second record: (k+G) vectors
 925 !---------------------------------------------------------------------------
 926 
 927  if (optkg/=0.and.option/=4) then
 928    if(use_f90==1)then
 929      write(unitwf) kg_k(1:3,1:optkg*npw)
 930    else if (wff%iomode==IO_MODE_MPI) then
 931 
 932      if (allocated(mpi_enreg%my_kgtab)) then
 933        me_cart_3d=xmpi_comm_rank(mpi_enreg%comm_bandspinorfft)
 934        ikpt_this_proc=mpi_enreg%my_kpttab(ikpt)
 935        call xderiveWRecInit(wff,ios,me_cart_3d)
 936        if (mpi_enreg%me_spinor==0) then
 937          call xderiveWrite(wff,kg_k,3,npw,mpi_enreg%comm_bandfft, &
 938 &         mpi_enreg%my_kgtab(1:npw,ikpt_this_proc),ios)
 939        end if
 940        call xderiveWRecEnd(wff,ios,me_cart_3d)
 941      else
 942 !      MG does it work if we are not using FFT distribution ?
 943        call xderiveWRecInit(wff,ios )
 944        if (mpi_enreg%me_spinor==0) then
 945          call xderiveWrite(wff,kg_k,3,optkg*npw,Wff%spaceComm_mpiio,ios)
 946        end if
 947        call xderiveWRecEnd(wff,ios)
 948      end if
 949 
 950 #ifdef HAVE_NETCDF
 951    else if (wff%iomode == IO_MODE_ETSF) then
 952      ! Write the reduced_coordinates_of_plane_waves for this k point.
 953      NCF_CHECK(nf90_inq_varid(wff%unwff, "reduced_coordinates_of_plane_waves", kg_varid))
 954      NCF_CHECK(nf90_get_att(wff%unwff, kg_varid, "k_dependent", kdep))
 955      if (kdep == "no") then
 956        ncerr = nf90_put_var(wff%unwff, kg_varid, kg_k, start=[1,1], count=[3,npw])
 957      else
 958        ncerr = nf90_put_var(wff%unwff, kg_varid, kg_k, start=[1,1,ikpt], count=[3,npw,1])
 959      end if
 960      NCF_CHECK_MSG(ncerr, "putting kg_k")
 961 #endif
 962    end if ! end if wff%iomode
 963  else ! Still skip the record
 964    if (use_f90==1) then
 965      write(unitwf)
 966    else if (wff%iomode==IO_MODE_MPI) then
 967      call xderiveWRecInit(wff,wff%spaceComm_mpiio,ios)
 968      call xderiveWRecEnd(wff,wff%spaceComm_mpiio,ios)
 969    end if
 970  end if
 971 
 972 !---------------------------------------------------------------------------
 973 !Write the third record: eigenvalues and occupations
 974 !---------------------------------------------------------------------------
 975 
 976 !===== Case formeig=0: write eigenvalues and occupations =====
 977  if (formeig==0) then
 978    if (use_f90==1) then
 979      write(unitwf) (eigen(iband),iband=1,nband_disk),(occ(iband),iband=1,nband_disk)
 980    else if(wff%iomode==IO_MODE_MPI) then
 981      if (wff%me_mpiio==0) then
 982        call xderiveWRecInit(wff,ios)
 983        call xderiveWrite(wff,eigen,nband_disk,xmpi_comm_self,ios)
 984        call xderiveWrite(wff,occ,nband_disk,xmpi_comm_self,ios)
 985        call xderiveWRecEnd(wff,ios)
 986      end if
 987 #ifdef HAVE_MPI
 988      call MPI_BCAST(wff%offwff,1,wff%offset_mpi_type,0,wff%spaceComm_mpiio,ios)
 989 #endif
 990 #ifdef HAVE_NETCDF
 991    else if (wff%iomode == IO_MODE_ETSF) then
 992      ! Write eigenvalues and occupation factors.
 993      NCF_CHECK(nf90_inq_varid(wff%unwff, "eigenvalues", eig_varid))
 994      ncerr = nf90_put_var(wff%unwff, eig_varid, eigen, start=[1,ikpt,isppol], count=[mband,1,1])
 995      NCF_CHECK_MSG(ncerr, "putting eig_k")
 996 
 997      NCF_CHECK(nf90_inq_varid(wff%unwff, "occupations", occ_varid))
 998      ncerr = nf90_put_var(wff%unwff, occ_varid, occ, start=[1,ikpt,isppol], count=[mband,1,1])
 999      NCF_CHECK_MSG(ncerr, "putting occ_k")
1000 #endif
1001    end if
1002 
1003 !  ===== Case formeig=1: write matrix of eigenvalues =====
1004 !  Will be written later (together with wave-functions)
1005  else if(formeig==1)then
1006  end if ! formeig
1007 
1008 !---------------------------------------------------------------------------
1009 !Write the wave-function coefficients
1010 !---------------------------------------------------------------------------
1011 
1012 !===== Case formeig=0: write only wave-functions =====
1013  if (formeig==0) then
1014 !  If option=4, do not write wave functions
1015    if (option/=4) then
1016      if (use_f90==1) then
1017        do iband=1,nband_disk
1018          ipw=(iband-1)*npwso+icg
1019          if(option/=5)then
1020            write(unitwf) cg(1:2,ipw+1:ipw+npwso) ! VALGRIND complains some elements of cg are not initialized, but written
1021          else
1022            write(unitwf)
1023          end if
1024        end do
1025      else if(wff%iomode==IO_MODE_MPI)then
1026        cg_ptr => cg ! Need pointer to bypass "inout" intent attribute
1027        call WffReadWrite_mpio(wff,2,cg_ptr,mcg,icg,nband_disk,npwso,npwsotot,ind_cg_mpi_to_seq,ios)
1028        nullify(cg_ptr)
1029 #ifdef HAVE_NETCDF
1030      else if (wff%iomode == IO_MODE_ETSF .and. option/=5) then
1031 
1032        NCF_CHECK(nf90_inq_varid(wff%unwff, "reduced_coordinates_of_plane_waves", kg_varid))
1033        NCF_CHECK(nf90_get_att(wff%unwff, kg_varid, "k_dependent", kdep))
1034 
1035        ! The coefficients_of_wavefunctions on file have shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
1036        NCF_CHECK(nf90_inq_varid(wff%unwff, "coefficients_of_wavefunctions", cg_varid))
1037        !write(std_out,*)"put cg, count: ",[2,npw,nspinor,nband,1,1]
1038        ncerr = nf90_put_var(wff%unwff, cg_varid, cg(:, icg+1:), start=[1,1,1,1,ikpt,isppol], &
1039        count=[2,npw,nspinor,nband,1,1])
1040        NCF_CHECK_MSG(ncerr, "putting cg_k")
1041        !write(std_out,*)"after cg"
1042 #endif
1043      end if
1044    end if ! option/=4
1045 
1046 !  ===== Case formeig=1: write eigenvalues and wave-functions =====
1047  else if(formeig==1)then
1048 
1049 !  Not available for NETCDF and ETSF_IO
1050    ABI_CHECK(wff%iomode /= IO_MODE_ETSF, "ETSF-write-eigen1 not coded!")
1051 
1052 !  ABI_CHECK(nband_disk==nband,"nband_disk!=nband")
1053 
1054    nband2=2*nband_disk
1055    do iband=1,nband_disk
1056      ipw=(iband-1)*npwso+icg
1057      ii=(iband-1)*nband2
1058      if(use_f90==1)then
1059        write(unitwf) eigen(1+ii:nband2+ii)
1060        if (option/=5) then
1061          write(unitwf) cg(1:2,1+ipw:npwso+ipw)
1062        else
1063          write(unitwf)
1064        end if
1065      else if(wff%iomode==IO_MODE_MPI)then
1066 !      Should use an access with a "view"
1067        call xderiveWRecInit(wff,ios)
1068        call xderiveWrite(wff,eigen(1+ii:ii+nband2),nband2,wff%spaceComm_mpiio,ios)
1069        call xderiveWRecEnd(wff,ios)
1070        if (option/=4) then
1071          call xderiveWRecInit(wff,ios)
1072          call xderiveWrite(wff,cg(1:2,ipw+1:ipw+npwso),2,npwso,wff%spaceComm_mpiio,ios)
1073          call xderiveWRecEnd(wff,ios)
1074        end if
1075      end if
1076    end do
1077 
1078  end if ! formeig
1079 
1080 !---------------------------------------------------------------------------
1081 !Final statements
1082 !---------------------------------------------------------------------------
1083 
1084  if (allocated(ind_cg_mpi_to_seq)) then
1085    ABI_DEALLOCATE(ind_cg_mpi_to_seq)
1086  end if
1087 
1088  RETURN
1089 
1090 !Silence compiler warning.
1091  ABI_UNUSED((/ii,mpi_enreg%me/))
1092 
1093 end subroutine writewf