Logo Search packages:      
Sourcecode: abinit version File versions  Download package


!!****f* ABINIT/newkpt
!! newkpt
!! This subroutine writes a starting guess for wave function (set 2)
!! It performs a "zero order" interpolation, ie simply
!! searches the nearest available k-point.
!! The data (set 1) associated with this point is either
!! read from a disk file (with a random access reading routine),
!! or input as argument.
!! Copyright (C) 1998-2007 ABINIT group (DCA, XG, GMR, ZL, AR, MB)
!! 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.
!!  ceksp2=if 1, center the sphere of pw on Gamma; if 0, on each k-point.
!!  doorth=1 to do orthogonalization
!!  debug=>0 for debugging output
!!  ecut1=kinetic energy cutoffs for basis sphere 1 (hartree)
!!  ecut2=kinetic energy cutoffs beyond which the coefficients of wf2 vanish (Ha)
!!  ecut2_eff=kinetic energy cut-off for basis sphere 2 (hartree)
!!  exchn2n3=if 1, n2 and n3 are exchanged
!!  fill=if 1, fill the supplementary bands ; if 0, reduce the number of bands
!!             Note : must have fill/=0 in the parallel execution
!!  formeig=if 0, GS format for wfs, eig and occ ; if 1, RF format.
!!  gmet1(3,3), gmet2(3,3)=reciprocal space metrics (bohr^-2)
!!  headform1=header format (might be needed to read the block of wfs)
!!  indkk(nkpt2*sppoldbl,6)=describe k point number of kptns1 that allows to
!!   generate wavefunctions closest to given kpt2 (and possibly isppol2=2)
!!   indkk(:,1)=k point number of kpt1
!!   indkk(:,2)=symmetry operation to be applied to kpt1, to give kpt1a
!!    (if 0, means no symmetry operation, equivalent to identity )
!!   indkk(:,3:5)=shift in reciprocal space to be given to kpt1a,
!!    to give kpt1b, that is the closest to ikpt2.
!!   indkk(:,6)=1 if time-reversal was used to generate kpt1a from kpt1, 0 otherwise
!!  iout=unit number for output file
!!  ireadwf=if 0, no reading of disk wavefunction file (random or 0.0 initialisation)
!!  istwfk1(nkpt1)=input parameter that describes the storage of wfs in set1
!!  istwfk2(nkpt2)=input parameter that describes the storage of wfs in set2
!!  kg2(3,mpw2*mkmem2)=dimensionless coords of G vecs in basis sphere at k point
!!  kptns1(3,nkpt1), kptns2(3,nkpt2)=k point sets (reduced coordinates)
!!  mband2= maximum number of bands of the output wavefunctions
!!  mcg=dimension of the cg array
!!   In case mkmem2/=0, all the output data must find their place in cg,
!!    so that mcg must be at least Sum(ikpt,isppol) [npw*nspinor*nband](ikpt,isppol)
!!    where these data are related to the output parameters
!!   In case mkmem1/=0, the same is true, for the input parameters,
!!    however, the maximum number of bands that will be read
!!    will be at most (mband2/nspinor2)*nspinor1
!!   In case mkmem1==0 and mkmem2==0, one must have at least mpw*nspinor*mband
!!    for BOTH the input and output parameters, taking into account the
!!    maximal number of band to be read, described above.
!!   In case mkmem1/=0 and mkmem2/=0, it is expected that the input cg array
!!    is organised using the output parameters nkpt2, nband2 ...
!!    This is needed, in order to use the same pointer.
!!  mkmem1= if 0, the input wf, eig, occ are available from disk
!!  mkmem2= if 0, the output wf, eig, occ must be written onto disk
!!  mpi_enreg=informations about MPI parallelization
!!  mpw1=maximum allowed number of planewaves at any k, for the input wf file
!!  mpw2=maximum allowed number of planewaves at any k, for the output wf file
!!  nband1(nkpt1*nsppol1)=number of bands, at each k point, on disk
!!  nband2(nkpt2*nsppol2)=desired number of bands at each k point
!!  ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/input_variables/vargs.htm#ngfft
!!  nkpt1, nkpt2=number of k points in each set
!!  npwarr1(nkpt1)=array holding npw for each k point (input wf file).
!!  npwarr2(nkpt2)=array holding npw for each k point (output wf file).
!!  nspinor1,nspinor2=number of spinorial components of the wavefunctions
!!   for each wf file (input or output)
!!  nsppol1=1 for unpolarized, 2 for spin-polarized, input wf file
!!  nsppol2=1 for unpolarized, 2 for spin-polarized, output wf file
!!  nsym=number of symmetry elements in space group
!!  optorth= 1 if the WFS have to be orthogonalized; 0 otherwise
!!  prtvol=control print volume and debugging
!!  restart= if 2, conversion between wavefunctions
!!           if 1, direct restart is allowed (see hdr_check.f)
!!  rprimd(3,3)=dimensional primitive translations for real space (bohr)
!!  sppoldbl= if 1, no doubling of the number if spins thanks to antiferromagn
!!    if 2, deduce nsppol=2 from nsppol=1, using Shubnikov symmetries
!!  symafm(nsym)=(anti)ferromagnetic part of symmetry operations
!!  symrel(3,3,nsym)=symmetry operations in real space in terms
!!   of primitive translations
!!  tnons(3,nsym)=nonsymmorphic translations for symmetry operations
!!  unkg2=unit number for storage of basis sphere data: stores indirect
!!   indexing array and integer coordinates for all planewaves in basis
!!   sphere for each k point being considered (kptns2 set)
!!  wffinp=structure info of input wf file unit number
!!  wffout=structure info of output wf file unit number
!!  (see side effects)
!!     The following arrays are input if mkmem1/=0, otherwise their input
!!     values are taken from disk, and are output if mkmem2/=0, otherwise
!!     their output values are written on disk.
!!     The location of the block for a given spin-k point at input MUST
!!     be the same as the location of the corresponding spin-k point at output.
!!  cg(2,mcg)=complex wf array
!!  eigen(mband2*(2*mband2)**formeig *nkpt2*nsppol2)=
!!    eigenvalues (input or init to large number for GS or init to 0.0 for RF), (Ha)
!!  occ(mband2*nkpt2*nsppol2)=occupation (input or init to 0.0)  NOT USED NOW
!! * When reading from disk, it is expected that the next record of
!! the wffinp%unwff disk unit is the first record of the first wavefunction block.
!! * When the data is input as argument, it is assumed that the
!! data for each spin- k wavefunction block is located at the proper
!! corresponding location of the output array (this is to be described).
!! * The information is pumped onto an fft box for the conversion.
!! This allows for changing the number of plane waves.
!! * In the present status of this routine, occ is not output.
!!      inwffil,newsp
!!      leave_new,leave_test,pareigocc,prmat,randac,rdnpw,rwwf,timab,wfconv
!!      wffreadskipk,wrtout,xcomm_init,xmaster_init,xme_init

#if defined HAVE_CONFIG_H
#include "config.h"

subroutine newkpt(ceksp2,cg,debug,doorth,ecut1,ecut2,ecut2_eff,eigen,exchn2n3,fill,&
&                  formeig,gmet1,gmet2,headform1,indkk,iout,ireadwf,istwfk1,istwfk2,&
&                  kg2,kptns1,kptns2,mband2,mcg,mkmem1,mkmem2,mpi_enreg,mpw1,mpw2,&
&                  nband1,nband2,ngfft,nkpt1,nkpt2,npwarr1,npwarr2,nspinor1,nspinor2,&
&                  nsppol1,nsppol2,nsym,occ,optorth,prtvol,restart,rprimd,sppoldbl,symafm,&
&                  symrel,tnons,unkg2,wffinp,wffout)

 use defs_basis
 use defs_datatypes

!This section has been created automatically by the script Abilint (TD). Do not modify these by hand.
 use interfaces_01manage_mpi
 use interfaces_11util
 use interfaces_13io_mpi
 use interfaces_14iowfdenpot
 use interfaces_14wfs
 use interfaces_18seqpar, except_this_one => newkpt
 use interfaces_lib01hidempi
!End of the abilint section

 implicit none

!Arguments ------------------------------------
 integer, intent(in) :: ceksp2,debug,doorth,exchn2n3,fill,formeig,headform1,iout,ireadwf
 integer, intent(in) :: mband2,mcg,mkmem1,mkmem2,mpw1,mpw2,nkpt1,nkpt2,nspinor1
 integer, intent(inout) :: nspinor2
 integer, intent(in) :: nsppol1,nsppol2,nsym,optorth,prtvol,restart,sppoldbl,unkg2
 real(dp), intent(in) :: ecut1,ecut2,ecut2_eff
 type(MPI_type), intent(inout) :: mpi_enreg
 type(wffile_type), intent(inout) :: wffinp,wffout
 integer, intent(in) :: indkk(nkpt2*sppoldbl,6),istwfk1(nkpt1),istwfk2(nkpt2)
 integer, intent(in) :: kg2(3,mpw2*mkmem2),nband1(nkpt1*nsppol1),nband2(nkpt2*nsppol2)
 integer, intent(in) :: ngfft(18),npwarr1(nkpt1),npwarr2(nkpt2),symafm(nsym)
 integer, intent(in) :: symrel(3,3,nsym)
 real(dp), intent(in) :: gmet1(3,3),gmet2(3,3),kptns1(3,nkpt1)
 real(dp), intent(in) :: kptns2(3,nkpt2),rprimd(3,3),tnons(3,nsym)
 real(dp), intent(out) :: cg(2,mcg),eigen(mband2*(2*mband2)**formeig*nkpt2*nsppol2)
 real(dp), intent(out) :: occ(mband2*nkpt2*nsppol2)

!Local variables-------------------------------
 integer,parameter :: fform=2,init_random=-5,nkpt_max=50,tobox=1,tosph=-1,wr=2
 integer :: aux_stor,band_index,fftalg,i1,i2,iband,icg,icg_aux,idum,ierr
 integer :: igs1,igs2,ii,ikg2,ikl,ikpt1,ikpt10,ikpt2,ikptsp_prev,index,inplace
 integer :: iproc,ipw,ispinor,isppol1,isppol2,istwf10_k,istwf1_k,istwf2_k
 integer :: localrdwf,master,mband1,mband_rd,mband_rw,mcg_aux,me,mgfft,n1,n2,n3
 integer :: n4,n5,n6,nb_band,nbd1,nbd1_rd,nbd2,nbremn,ngb1,ngb2,nkpt_eff,npw1
 integer :: npw2,nrecwf,spaceComm,test_cycle,tim_rwwf,posoffloc,posoff
 character(len=500) :: message
 integer,allocatable :: kg1(:,:),kg2_k(:,:),kg_dum(:,:)
 real(dp) :: cgdum(2,1),dummy(1),kpoint(3),kptsph(3),tsec(2)
 real(dp),allocatable :: cfft(:,:,:,:),cg_aux(:,:),eig_k(:),occ_k(:)
 real(dp),allocatable :: wavef1(:,:),wavef2(:,:)

! *************************************************************************

!write(6,*)' newkpt : enter'


! write(6,*)' newkpt : paralbd=',mpi_enreg%paralbd
! if(restart==2)stop
!write(6,*)' newkpt : ireadwf,nkpt2= ',ireadwf,nkpt2
!write(6,*)' newkpt : npwarr2(:)=',npwarr2(:)
!write(6,*)' newkpt : mkmem1=',mkmem1
!write(6,*)' newkpt : set debug=1'
!write(6,*)' newkpt : cg array='
!do ikpt2=1,nkpt2
! nbd2=nband2(ikpt2)
! npw2=npwarr2(ikpt2)
! do iband=1,nbd2
!  write(6,*)' new band, icg=',icg
!  do ipw=1,1
!   write(6, '(3i4,2es20.10)' )ikpt2,iband,ipw,cg(:,icg+ipw)
!  end do
!  icg=icg+npw2
! end do
!end do

!Init me
 call xme_init(mpi_enreg,me)
!Init master
 call xmaster_init(mpi_enreg,master)
!Init mpi_comm
 call xcomm_init(mpi_enreg,spaceComm)

 if( (nsppol1==2 .and. nspinor2==2) .or. &
&    (nspinor1==2.and. nsppol2==2 )       )then
! This is not yet possible. See later for a message about where
! to make the needed modifs.
  write(message, '(a,a,a,a,a,a,a,a,i2,a,i2,a,a,i2,a,i2,a,a,a,a)' )ch10,&
&  ' newkpt : ERROR -',ch10,&
&  '  The wavefunction translator is (still) unable to interchange',ch10,&
&  '  spin-polarized wfs and spinor wfs. However,',ch10,&
&  '  the input  variables are nsppol1=',nsppol1,', and nspinor1=',nspinor1,ch10,&
&  '  the output variables are nsppol2=',nsppol2,', and nspinor2=',nspinor2,ch10,&
&  '  Action : use a non-spin-polarized wf to start a spinor wf,',ch10,&
&  '           and a non-spinor wf to start a spin-polarized wf.'
  call wrtout(6,message,'COLL')
  call leave_new('COLL')
 end if


 if(mkmem1==0 .and. mkmem2==0)then
   write(message, '(a,a,a,a,i9,a,a,a,i9)' )ch10,&
&   ' newkpt : BUG -',ch10,&
&   '  The dimension mcg=',mcg,', should be larger than',ch10,&
&   '  mband_rd=',mband_rd
   call wrtout(6,message,'COLL')
   call leave_new('COLL')
  end if
   write(message, '(a,a,a,a,i9,a,a,a,i4,a,i6,a,i2)' )ch10,&
&   ' newkpt : BUG -',ch10,&
&   '  The dimension mcg=',mcg,', should be larger than',ch10,&
&   '  the product of mband2=',mband2,', mpw2=',mpw2,', and nspinor2=',nspinor2
   call wrtout(6,message,'COLL')
   call leave_new('COLL')
  end if
 end if

 ikpt10 = 0

 if( (prtvol==0.or.prtvol==1) .and. nkpt_eff>nkpt_max ) nkpt_eff=nkpt_max


 if (debug>0 .and. me==0) then
  write(6, '(a)' ) ' newkpt:  kptns1'
  call prmat (kptns1, 3, nkpt1, 3)
  write(6, '(a)' ) ' newkpt:  kptns2'
  call prmat (kptns2, 3, nkpt2, 3)
 end if

!write(6,*)' newkpt :me,restart before loop ',me,restart


!Do outer loop over spins
 do isppol2=1,nsppol2

  if (nsppol2==2 .and. me==0) then
   write(06, '(a,i5)' ) ' newkpt: spin channel isppol2 = ',isppol2
  end if

  if (restart==1 .and. mkmem2==0) rewind (unkg2)

! Do loop over new k point set
  do ikpt2=1,nkpt2



!   Announce the treatment of k point ikpt
!    This message might be overwritten in parallel
     write(message, '(a,i6,a,i8,a,i4)' ) &
&     'P newkpt: treating ',nbd2,' bands with npw=',npw2,&
&     ' for ikpt=',ikpt2
!     This message might be overwritten in parallel
       do iproc=0,mpi_enreg%nproc-1
        do iband=1,nbd2
         if(mpi_enreg%proc_distrb(ikpt2,iband,isppol2) == iproc)nb_band=nb_band+1
        end do
         write(message, '(a,i6,a,i8,a,i4,a,i4)' ) &
&         'P newkpt: treating ',nb_band,' bands with npw=',npw2,&
&         ' for ikpt=',ikpt2,' by node ',iproc
        end if
       end do
      end if
      if((mpi_enreg%paralbd==0) .or. (mpi_enreg%paralbd>1))then
       if (mpi_enreg%parareel == 0) then
        write(message, '(a,i6,a,i8,a,i4,a,i4)' )&
&        'P newkpt: treating ',nbd2,' bands with npw=',npw2,&
&        ' for ikpt=',ikpt2,' by node ',mpi_enreg%proc_distrb(ikpt2,1,isppol2)
        write(message, '(a,i6,a,i8,a,i4,a,i4)' )&
&        'P newkpt: treating ',nbd2,' bands with npw=',npw2,&
&        ' for ikpt=',ikpt2,' by node ',mpi_enreg%proc_distrb_para(mpi_enreg%ipara,ikpt2)
       end if
      end if
     end if
     call wrtout(iout,message,'COLL')
    end if

!   Cut the writing if the limit is reached
     write(message, '(a)' ) &
&     '  newkpt: prtvol=0 or 1, do not print more k-points.'
     call wrtout(iout,message,'COLL')
    end if

!  End of restart==1
   end if

    if (mpi_enreg%parareel == 0) then
     end if
     end if
    end if
      eigen(1+band_index : nbd2+band_index) = zero
!     occ(1+band_index : nbd2+band_index) = zero
      eigen(1+band_index : 2*nbd2**2+band_index) = 0.0_dp
     end if
!    In the case this k point does not belong to me, cycle
     if ((mkmem1==0) .and. (ireadwf==1) .and. (mpi_enreg%paralbd==1))then
      call WffReadSkipK(formeig,headform1,ikpt2,isppol2,wffinp)
     end if
    end if

   end if


    else if(mkmem2==0)then
!    Read the first line of a block and performs some checks on the unkg file.
     call rdnpw(ikpt2,isppol2,nbd2,npw2,nspinor2,0,unkg2)
!    Read k+g data
     read (unkg2) kg2_k(1:3,1:npw2)
    end if

   end if

!  Get ikpt1, the closest k from original set, from indkk
   if(sppoldbl==2 .and. isppol2==2)ikpt1=indkk(ikpt2+nkpt2,1)


!  Determine the spin polarization of the input data
   if(nsppol2==2 .and. nsppol1==1)isppol1=1

     write(message, '(a,i4,i8,a,i4,i8)' ) &
&     '- newkpt: read input wf with ikpt,npw=',ikpt1,npw1,&
&     ', make ikpt,npw=',ikpt2,npw2
     call wrtout(6,message,'PERS')
     if(iout/=6 .and. me==0)then
      call wrtout(iout,message,'PERS')
     end if
    else if(ikpt2==nkpt_eff+1)then
     write(message, '(a)' ) &
&     '- newkpt : prtvol=0 or 1, do not print more k-points.'
     call wrtout(6,message,'PERS')
     if(iout/=6 .and. me==0)then
      call wrtout(iout,message,'PERS')
     end if
    end if
   end if

!  Set up the number of bands to be read

!  Check that number of bands is not being increased if fill==0 --if so
!  print warning and reset new wf file nband2 to only allowed number
   if ( nbd2/nspinor2 > nbd1/nspinor1 .and. fill==0) then
     write(message, '(a,i8,a,i8,a,i8)' ) &
&     ' newkpt: nband2=',nbd2,' < nband1=',nbd1,&
&     ' => reset nband2 to ',nbd1
     call wrtout(6,message,'PERS')
    end if
   end if

!  write(6,*)' newkpt : before randac, nbd1=',nbd1
!  stop

!  Prepare the reading of the wavefunctions: the correct record is selected
!  WARNING : works only for GS - for RF the number of record differs
   if(restart==2 .and. mkmem1==0)then

     write(message, '(a,a,a,a,i5,a,i5,a,a,i5,a,i5)' ) ch10,&
&     ' newkpt : about to call randac',ch10,&
&     '  for ikpt1=',ikpt1,', ikpt2=',ikpt2,ch10,&
&     '  and isppol1=',isppol1,', isppol2=',isppol2
     call wrtout(6,message,'PERS')
    end if

    call randac(debug,headform1,ikptsp_prev,ikpt1,isppol1,&
&    mband2,nband1,nkpt1,nsppol1,wffinp)

   end if

!  Read the data for nbd2 bands at this k point
!  Must decide whether an auxiliary storage is needed
!  When mkmem1==0 and mkmem2==0 , the cg array should be large enough ...
!  When mkmem1==0 and mkmem2/=0 , each k-point block in cg might not be large enough
!   however, will read at most (nbd2/nspinor2)*nspinor1 bands from disk
!  When mkmem1/=0 , it is supposed that each input k-point block is smaller
!   than the corresponding output k-point block, so that the input data
!   have been placed already in cg, at the k-point location where they are needed
   if(mkmem2/=0 .and. mkmem1==0)then
    if( mcg_aux > npw2*nspinor2*nbd2 )then
     aux_stor=1 ; icg_aux=0
    end if
   end if


   if(mkmem1/=0 .and. ireadwf==1)then
!   Checks that nbd1 and nbd1_rd are equal if eig and occ are input
     write(message, '(a,a,a,a,a,a,i6,a,i6)' )ch10,&
&     ' newkpt: BUG -',ch10,&
&     '  When mkmem1/=0, one must have nbd1=nbd1_rd, while',ch10,&
&     '  nbd1=',nbd1,', and nbd1_rd=',nbd1_rd
     call wrtout(6,message,'PERS')
     call leave_new('PERS')
    end if
!   Need to put eigenvalues in eig_k, same for occ
!   Note use of band_index, since it is assumed that eigen and occ
!   already have spin-k point location structure than output.
     eig_k(1:nbd1_rd)=eigen(1+band_index : nbd1_rd+band_index)
!    occ_k(1:nbd1_rd)=occ(1+band_index : nbd1_rd+band_index)
    else if(formeig==1)then
!    The matrix of eigenvalues has size nbd1 ,  that must be equal
!    to nbd1_rd in the case mkmem1/=0)
&     eigen(1+band_index : 2*nbd1_rd**2+band_index)
    end if
   end if

!  Must read the wavefunctions if they are not yet in place
   if(mkmem1==0 .and. ireadwf==1)then

    if (debug>0 .and. restart==2) then
     write(message,'(a,i5,a,a,i5,a,i5,a)' ) &
&     ' newkpt : about to call rwwf with ikpt1=',ikpt1,ch10,&
&     ' and nband(ikpt1)=',nband1(ikpt1),' nbd2=',nbd2,'.'
     call wrtout(6,message,'PERS')
    end if

     call rwwf(cg,eig_k,formeig,headform1,icg,ikpt1,isppol1,kg_dum,mband_rw,mcg,nbd1_rd,nbd1,&
&     npw1,nspinor1,occ_k,1,0,tim_rwwf,wffinp)
     call rwwf(cg_aux,eig_k,formeig,headform1,icg_aux,ikpt1,isppol1,kg_dum,mband_rw,mcg_aux,nbd1_rd,nbd1,&
&     npw1,nspinor1,occ_k,1,0,tim_rwwf,wffinp)
    end if
   end if

   if(formeig==1 .and. nbd2/=nbd1_rd .and. ireadwf==1)then
!   Change the storage of eig_k
     do iband=nbd1_rd,1,-1
!     The factor of two is for complex eigenvalues
      do ii=2*nbd2,2*nbd1_rd+1,-1
      end do
      do ii=2*nbd1_rd,1,-1
      end do
     end do
    else if(nbd1_rd>nbd2)then
     do iband=1,nbd2
!     The factor of two is for complex eigenvalues
      do ii=1,2*nbd2
      end do
     end do
    end if
   end if

!  If change nsppol, must adapt the occupation numbers
!  if(nsppol1/=nsppol2)then
!   occ_k(1:nbd2)=occ_k(1:nbd2)*nsppol1/dbl(nsppol2)
!  then

!  In case nsppol1=2 and nspinor2=2, one should read
!  the other spin-component, and form a spinor wf here, before calling
!  wfconv. One should treat eig_k and occ_k as well.
!  A similar operation is to be performed when nspino1=2 and nsppol2=2

!  write(6,*)' newkpt: before wfconv'
!  write(6,*)' newkpt: mkmem2=',mkmem2
!  stop

!  Note the use of mband2, while mband is used inside
!  write(6,*) 'in newkpt,before wfconv,npw1,npw2',npw1,npw2
     call wfconv(ceksp2,cg,cg,debug,ecut1,ecut2,ecut2_eff,&
&     eig_k,eig_k,exchn2n3,formeig,gmet1,gmet2,icg,icg,idum,&
&     ikpt1,ikpt10,ikpt2,indkk,inplace,isppol2,istwfk1,istwfk2,&
&     kg1,kg2_k,kptns1,kptns2,mband2,mband2,mcg,mcg,mgfft,mpi_enreg,mpw1,mpw2,&
&     nbd1_rd,nbd2,&
&     ngfft,nkpt1,nkpt2,npw1,npw2,nspinor1,nspinor2,nsym,&
&     occ_k,occ_k,optorth,restart,rprimd,sppoldbl,symafm,symrel,tnons)
     call wfconv(ceksp2,cg_aux,cg_aux,debug,ecut1,ecut2,ecut2_eff,&
&     eig_k,eig_k,exchn2n3,formeig,gmet1,gmet2,icg_aux,icg_aux,idum,&
&     ikpt1,ikpt10,ikpt2,indkk,inplace,isppol2,istwfk1,istwfk2,&
&     kg1,kg2_k,kptns1,kptns2,mband2,mband2,mcg,mcg,mgfft,mpi_enreg,mpw1,mpw2,&
&     nbd1_rd,nbd2,&
&     ngfft,nkpt1,nkpt2,npw1,npw2,nspinor1,nspinor2,nsym,&
&     occ_k,occ_k,optorth,restart,rprimd,sppoldbl,symafm,symrel,tnons)
   end if

!  Finally write new wf to disk file or save in permanent file

!   Note that in this case, we are sure aux_stor==0
    call rwwf(cg,eig_k,formeig,0,0,ikpt2,isppol2,kg2_k,nbd2,mcg,nbd2,nbd2,npw2,nspinor2,&
&    occ_k,wr,1,tim_rwwf,wffout)

    end if


   end if

!   write(6,*)' newkpt : final content of cg,ikpt,icg=',ikpt2,icg
!   do iband=1,nbd2
!    do ipw=1,npw2
!     if(aux_stor==0)then
!      write(6, '(3i4,2es20.10)' )&
!&      ikpt2,iband,ipw,cg(1:2,ipw+(iband-1)*npw2+icg)
!     else
!      write(6, '(3i4,2es20.10)' )&
!&      ikpt2,iband,ipw,cg_aux(1:2,ipw+(iband-1)*npw2+icg)
!     end if
!    end do
!   end do

&   eig_k(1:nbd2*(2*nbd2)**formeig)
!  occ(1+band_index:nbd2+band_index)=occ_k(1:nbd2)

   else if(formeig==1)then
   end if


  end do ! ikpt2
 end do ! isppol2
  call timab(67,1,tsec)
  call leave_test(mpi_enreg)
  write(message, '(a)') ' newkpt: loop on k-points done in parallel'
  call wrtout(06,message,'COLL')

! Transmit eigenvalues (not yet occupation numbers)
! newkpt.f is not yet suited for RF format
! This routine works in both localrdwf=0 or 1 cases.
! However, in the present routine, localrdwf is to be considered
! as 1 always, since the transfer has been made in wfsinp .
  call pareigocc(eigen,formeig,localrdwf,mpi_enreg,mband2,nband2, &
&                nkpt2,nsppol2,occ)

  call timab(67,2,tsec)
 end if

!write(6,*)' newkpt : cg array='
!do ikpt2=1,nkpt2
! nbd2=nband2(ikpt2)
! npw2=npwarr2(ikpt2)
! kg2_k(:,1:npw2)=kg2(:,1+ikg2:npw2+ikg2)
! do iband=1,nbd2
!  write(6, '(a,4i4)' )' new band, ikpt2,iband,npw2,icg=',ikpt2,iband,npw2,icg
!  do ipw=1,npw2
!   write(6, '(4i4,2es20.10)' )ipw,kg2_k(:,ipw),cg(:,icg+ipw)
!  end do
!  icg=icg+npw2
! end do
! ikg2=ikg2+npw2
!end do
! write(6,*)' newkpt : exit'


end subroutine newkpt

Generated by  Doxygen 1.6.0   Back to index