mirror of https://gitlab.com/QEF/q-e.git
290 lines
7.6 KiB
Fortran
290 lines
7.6 KiB
Fortran
!
|
|
! Copyright (C) 2001 PWSCF group
|
|
! This file is distributed under the terms of the
|
|
! GNU General Public License. See the file 'License'
|
|
! in the root directory of the present distribution,
|
|
! or http://www.gnu.org/copyleft/gpl.txt .
|
|
!
|
|
!
|
|
#include "f_defs.h"
|
|
!---------------------------------------------------------------------
|
|
program fhi2upf
|
|
!---------------------------------------------------------------------
|
|
!
|
|
! Convert a pseudopotential written in the Fritz-Haber format
|
|
! (numerical format) to unified pseudopotential format
|
|
! Restrictions:
|
|
! - no semicore states
|
|
! Adapted from the converter written by Andrea Ferretti
|
|
!
|
|
implicit none
|
|
character(len=256) filein, fileout
|
|
!
|
|
!
|
|
call get_file ( filein )
|
|
open (unit = 1, file = filein, status = 'old', form = 'formatted')
|
|
call read_fhi(1)
|
|
close (1)
|
|
|
|
! convert variables read from FHI format into those needed
|
|
! by the upf format - add missing quantities
|
|
|
|
call convert_fhi
|
|
|
|
fileout=trim(filein)//'.UPF'
|
|
print '(''Output PP file in UPF format : '',a)', fileout
|
|
|
|
open(unit=2,file=fileout,status='unknown',form='formatted')
|
|
call write_upf(2)
|
|
close (unit=2)
|
|
|
|
stop
|
|
20 write (6,'("fhi2upf: error reading pseudopotential file name")')
|
|
stop
|
|
end program fhi2upf
|
|
|
|
module fhi
|
|
!
|
|
! All variables read from FHI file format
|
|
!
|
|
|
|
type angular_comp
|
|
real(8), pointer :: pot(:)
|
|
real(8), pointer :: wfc(:)
|
|
real(8), pointer :: grid(:)
|
|
real(8) :: amesh
|
|
integer :: nmesh
|
|
integer :: lcomp
|
|
end type angular_comp
|
|
|
|
!------------------------------
|
|
|
|
real(8) :: Zval ! valence charge
|
|
integer :: lmax_ ! max l-component used
|
|
|
|
logical :: nlcc_
|
|
real(8), allocatable :: rho_atc_(:) ! core charge
|
|
|
|
type (angular_comp), pointer :: comp(:) ! PP numerical info
|
|
! (wfc, grid, potentials...)
|
|
!------------------------------
|
|
|
|
end module fhi
|
|
!
|
|
! ----------------------------------------------------------
|
|
subroutine read_fhi(iunps)
|
|
! ----------------------------------------------------------
|
|
!
|
|
use fhi
|
|
implicit none
|
|
integer, parameter :: Nl=7 ! max number of l-components
|
|
integer :: iunps
|
|
real(8) :: r, rhoc, drhoc, d2rhoc
|
|
!
|
|
|
|
integer :: l, i, idum, mesh
|
|
|
|
! Starting file reading
|
|
|
|
read(iunps,*) Zval, lmax_
|
|
lmax_ = lmax_ - 1
|
|
|
|
if (lmax_+1 > Nl) then
|
|
write (6,'("read_fhi: too many l-components...stopping")')
|
|
stop
|
|
end if
|
|
|
|
do i=1,10
|
|
read(iunps,*) ! skipping 11 lines
|
|
end do
|
|
|
|
allocate( comp(0:lmax_) )
|
|
|
|
do l=0,lmax_
|
|
comp(l)%lcomp = l
|
|
read(iunps,*) comp(l)%nmesh, comp(l)%amesh
|
|
if ( l > 0) then
|
|
if (comp(l)%nmesh /= comp(0)%nmesh .or. &
|
|
comp(l)%amesh /= comp(0)%amesh ) then
|
|
write(6,'("read_fhi: different radial grids not allowed...stopping")')
|
|
stop
|
|
end if
|
|
end if
|
|
mesh = comp(l)%nmesh
|
|
allocate( comp(l)%wfc(mesh), & ! wave-functions
|
|
comp(l)%pot(mesh), & ! potentials
|
|
comp(l)%grid(mesh) ) ! real space radial grid
|
|
! read the above quantities
|
|
do i=1,mesh
|
|
read(iunps,*) idum, comp(l)%grid(i), &
|
|
comp(l)%wfc(i), &
|
|
comp(l)%pot(i)
|
|
end do
|
|
end do
|
|
|
|
nlcc_ =.false.
|
|
allocate(rho_atc_(comp(0)%nmesh))
|
|
mesh = comp(0)%nmesh
|
|
do i=1,mesh
|
|
read(iunps,*,end=10, err=20) r, rho_atc_(i), drhoc, d2rhoc
|
|
if ( abs( r - comp(0)%grid(i) ) > 1.d-6 ) then
|
|
write(6,'("read_fhi: radial grid for core charge? stopping")')
|
|
stop
|
|
end if
|
|
end do
|
|
nlcc_ = .true.
|
|
! ----------------------------------------------------------
|
|
write (6,'(a)') 'Pseudopotential with NLCC successfully read'
|
|
! ----------------------------------------------------------
|
|
return
|
|
10 continue
|
|
! ----------------------------------------------------------
|
|
write (6,'(a)') 'Pseudopotential without NLCC successfully read'
|
|
! ----------------------------------------------------------
|
|
return
|
|
!
|
|
20 write(6,'("read_fhi: error reading core charge")')
|
|
stop
|
|
!
|
|
100 write(6,'("read_fhi: error reading pseudopotential file")')
|
|
stop
|
|
|
|
end subroutine read_fhi
|
|
|
|
! ----------------------------------------------------------
|
|
subroutine convert_fhi
|
|
! ----------------------------------------------------------
|
|
!
|
|
use fhi
|
|
use upf
|
|
use funct, ONLY : set_dft_from_name, get_iexch, get_icorr, get_igcx, get_igcc
|
|
implicit none
|
|
real(8), parameter :: rmax = 10.0
|
|
real(8), allocatable :: aux(:)
|
|
real(8) :: vll
|
|
character (len=20):: dft
|
|
integer :: lloc, kkbeta
|
|
integer :: l, i, ir, iv
|
|
!
|
|
print '("Atom name > ",$)'
|
|
read (5,'(a)') psd
|
|
print '("l local (max: ",i1,") > ",$)', lmax_
|
|
read (5,*) lloc
|
|
print '("DFT > ",$)'
|
|
read (5,'(a)') dft
|
|
|
|
write(generated, '("Generated using Fritz-Haber code")')
|
|
write(date_author,'("Author: unknown Generation date: as well")')
|
|
comment = 'Info: automatically converted from FHI format'
|
|
! reasonable assumption
|
|
rel = 1
|
|
rcloc = 0.0
|
|
nwfs = lmax_+1
|
|
allocate( els(nwfs), oc(nwfs), epseu(nwfs))
|
|
allocate(lchi(nwfs), nns(nwfs) )
|
|
allocate(rcut (nwfs), rcutus (nwfs))
|
|
do i=1, nwfs
|
|
print '("Wavefunction # ",i1,": label, occupancy > ",$)', i
|
|
read (5,*) els(i), oc(i)
|
|
nns (i) = 0
|
|
lchi(i) = i-1
|
|
rcut(i) = 0.0
|
|
rcutus(i)= 0.0
|
|
epseu(i) = 0.0
|
|
end do
|
|
|
|
pseudotype = 'NC'
|
|
nlcc = nlcc_
|
|
zp = Zval
|
|
etotps = 0.0
|
|
ecutrho=0.0
|
|
ecutwfc=0.0
|
|
if ( lmax_ == lloc) then
|
|
lmax = lmax_-1
|
|
else
|
|
lmax = lmax_
|
|
end if
|
|
nbeta= lmax_
|
|
mesh = comp(0)%nmesh
|
|
ntwfc= nwfs
|
|
allocate( elsw(ntwfc), ocw(ntwfc), lchiw(ntwfc) )
|
|
do i=1, nwfs
|
|
lchiw(i) = lchi(i)
|
|
ocw(i) = oc(i)
|
|
elsw(i) = els(i)
|
|
end do
|
|
call set_dft_from_name(dft)
|
|
iexch = get_iexch()
|
|
icorr = get_icorr()
|
|
igcx = get_igcx()
|
|
igcc = get_igcc()
|
|
|
|
allocate(rab(mesh))
|
|
allocate( r(mesh))
|
|
r = comp(0)%grid
|
|
rab = r * log( comp(0)%amesh )
|
|
|
|
if (nlcc) then
|
|
allocate (rho_atc(mesh))
|
|
rho_atc(:) = rho_atc_(:) / (4.d0*3.141592653589793d0)
|
|
end if
|
|
|
|
allocate (vloc0(mesh))
|
|
! the factor 2 converts from Hartree to Rydberg
|
|
vloc0 = 2.d0*comp(lloc)%pot
|
|
|
|
if (nbeta > 0) then
|
|
|
|
allocate(ikk2(nbeta), lll(nbeta))
|
|
kkbeta=mesh
|
|
do ir = 1,mesh
|
|
if ( r(ir) > rmax ) then
|
|
kkbeta=ir
|
|
exit
|
|
end if
|
|
end do
|
|
ikk2(:) = kkbeta
|
|
allocate(aux(kkbeta))
|
|
allocate(betar(mesh,nbeta))
|
|
allocate(qfunc(mesh,nbeta,nbeta))
|
|
allocate(dion(nbeta,nbeta))
|
|
allocate(qqq (nbeta,nbeta))
|
|
qfunc(:,:,:)=0.0d0
|
|
dion(:,:) =0.d0
|
|
qqq(:,:) =0.d0
|
|
iv=0
|
|
do i=1,nwfs
|
|
l=lchi(i)
|
|
if (l.ne.lloc) then
|
|
iv=iv+1
|
|
lll(iv)=l
|
|
do ir=1,kkbeta
|
|
! FHI potentials are in Hartree
|
|
betar(ir,iv) = 2.d0 * comp(l)%wfc(ir) * &
|
|
( comp(l)%pot(ir) - comp(lloc)%pot(ir) )
|
|
aux(ir) = comp(l)%wfc(ir) * betar(ir,iv)
|
|
end do
|
|
call simpson(kkbeta,aux,rab,vll)
|
|
dion(iv,iv) = 1.0d0/vll
|
|
end if
|
|
enddo
|
|
|
|
end if
|
|
|
|
allocate (rho_at(mesh))
|
|
rho_at = 0.d0
|
|
do i=1,nwfs
|
|
l=lchi(i)
|
|
rho_at = rho_at + ocw(i) * comp(l)%wfc ** 2
|
|
end do
|
|
|
|
allocate (chi(mesh,ntwfc))
|
|
do i=1,ntwfc
|
|
chi(:,i) = comp(i-1)%wfc(:)
|
|
end do
|
|
! ----------------------------------------------------------
|
|
write (6,'(a)') 'Pseudopotential successfully converted'
|
|
! ----------------------------------------------------------
|
|
return
|
|
end subroutine convert_fhi
|