Assumption of lattice vectors

[atbug]

I was trying to parse WAVECAR myself and was trying to read both your program and Feenstra's program (http://www.andrew.cmu.edu/user/feenstra/wavetrans/source/WaveTrans.f90). It seems the following code assumes the lattice vectors are orthogonal: https://github.com/QijingZheng/VaspBandUnfolding/blob/8cd0a45a531bbafb1ab095df2ab387e65c590d3a/vaspwfc.py#L178-L180

If the reciprocal lattice vectors are b1=[-0.99, 0.01, 0.0], b2=[0.99, 0.01, 0.0], b3=[0, 0, 1.0] and the wave vector corresponding to the energy cutoff is G=10, the norm of n*b1+n*b2 is smaller than G as long as n<500. But the above code would estimate n<=100.

[QijingZheng]

Hi

In my code, there is no such assumption that the lattice should be orthogonal. In fact, the snippet you showed is just a rewrite of the Fortran version in VASP. Please refer to the VASP source code, e.g. vasp.5.4.4 main.F ln966-968. Moreover, the number of g-vectors generated by my code and VASP are consistent, otherwise there will be an error and the code will not proceed.

In your example, your estimation is correct about n <= 500. However, please note that Anorm is the norm of the real-space basis vectors other than the reciprocal-space ones.

Since you mentioned WaveTrans.f90, I have to tell you that the code calculates the real-space wavefunction in a very stupid way, i.e. the author just added all the plane-waves one by one, which is inefficient and takes a long time. On the other hand, the more efficient way is to use the FFT to do the transformation. I have seen someone used the wavetrans.f90 algorithm to plot the wavefunction and it took about 1 hour to get just one single KS orbital, while my code get the same orbital in just 1 second.

http://www.andrew.cmu.edu/user/feenstra/wavetrans/source/WaveTransPlot.f90

do iz=0,2*nb3max
   z=dble(iz)/dble(1+2*nb3max)
   xyz(3)=z
   csum=cmplx(0.,0.)
   do iplane=1,nplane
      ig=igall(:,iplane)
      wkpg=wk+ig
      csum=csum+coeff(iplane)* &
           cdexp(2.*pi*cmplx(0.,1.)*dot_product(wkpg,xyz))
   enddo
   csum=csum/dsqrt(Vcell)
!!$ output z*a3mag for units of Angstroms
   write(6,*) sngl(z),sngl(real(csum)),sngl(dimag(csum))
enddo

[atbug]

Thank you.

I am trying to port the code to julia. I now realize your code is indeed correct. Sorry for the interuption.