2007-10-03 00:54:13 +08:00
|
|
|
!
|
2009-08-01 22:26:46 +08:00
|
|
|
! Copyright (C) 2006-2007 Quantum ESPRESSO group
|
2006-03-07 21:26:52 +08:00
|
|
|
! 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 .
|
|
|
|
!
|
|
|
|
!---------------------------------------------------------------------
|
|
|
|
MODULE read_uspp_module
|
|
|
|
!---------------------------------------------------------------------
|
|
|
|
!
|
|
|
|
! routines reading ultrasoft pseudopotentials in older formats:
|
|
|
|
! Vanderbilt's code and Andrea's RRKJ3 format
|
|
|
|
!
|
|
|
|
USE kinds, ONLY: DP
|
2007-10-24 03:47:26 +08:00
|
|
|
USE parameters, ONLY: lmaxx, lqmax
|
2007-10-08 15:33:57 +08:00
|
|
|
USE io_global, ONLY: stdout
|
2006-03-07 21:26:52 +08:00
|
|
|
USE funct, ONLY: set_dft_from_name, dft_is_hybrid, dft_is_meta, &
|
|
|
|
set_dft_from_indices
|
|
|
|
!
|
|
|
|
! Variables above are not modified, variables below are
|
|
|
|
!
|
2007-10-08 15:33:57 +08:00
|
|
|
USE uspp_param, ONLY: oldvan
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
IMPLICIT NONE
|
|
|
|
SAVE
|
|
|
|
PRIVATE
|
|
|
|
PUBLIC :: readvan, readrrkj
|
|
|
|
!
|
|
|
|
CONTAINS
|
|
|
|
!---------------------------------------------------------------------
|
2007-10-03 00:54:13 +08:00
|
|
|
subroutine readvan( iunps, is, upf )
|
2006-03-07 21:26:52 +08:00
|
|
|
!---------------------------------------------------------------------
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
! Read Vanderbilt pseudopotential from unit "iunps"
|
|
|
|
! for species "is" into the structure "upf"
|
2006-03-07 21:26:52 +08:00
|
|
|
! info on DFT level in module "funct"
|
|
|
|
!
|
|
|
|
! ------------------------------------------------------
|
|
|
|
! Important:
|
|
|
|
! ------------------------------------------------------
|
|
|
|
! The order of all l-dependent objects is always s,p,d
|
|
|
|
! ------------------------------------------------------
|
|
|
|
! potentials, e.g. vloc_at, are really r*v(r)
|
|
|
|
! wave funcs, e.g. chi, are really proportional to r*psi(r)
|
|
|
|
! and are normalized so int (chi**2) dr = 1
|
|
|
|
! thus psi(r-vec)=(1/r)*chi(r)*y_lm(theta,phi)
|
|
|
|
! conventions carry over to beta, etc
|
|
|
|
! charge dens, e.g. rho_atc, really 4*pi*r**2*rho
|
|
|
|
!
|
|
|
|
! ------------------------------------------------------
|
|
|
|
! Notes on qfunc and qfcoef:
|
|
|
|
! ------------------------------------------------------
|
|
|
|
! Since Q_ij(r) is the product of two orbitals like
|
|
|
|
! psi_{l1,m1}^star * psi_{l2,m2}, it can be decomposed by
|
|
|
|
! total angular momentum L, where L runs over | l1-l2 | ,
|
|
|
|
! | l1-l2 | +2 , ... , l1+l2. (L=0 is the only component
|
|
|
|
! needed by the atomic program, which assumes spherical
|
|
|
|
! charge symmetry.)
|
|
|
|
!
|
|
|
|
! Recall qfunc(r) = y1(r) * y2(r) where y1 and y2 are the
|
|
|
|
! radial parts of the wave functions defined according to
|
|
|
|
!
|
|
|
|
! psi(r-vec) = (1/r) * y(r) * Y_lm(r-hat) .
|
|
|
|
!
|
|
|
|
! For each total angular momentum L, we pseudize qfunc(r)
|
|
|
|
! inside rc as:
|
|
|
|
!
|
|
|
|
! qfunc(r) = r**(L+2) * [ a_1 + a_2*r**2 + a_3*r**4 ]
|
|
|
|
!
|
|
|
|
! in such a way as to match qfunc and its 1'st derivative at
|
|
|
|
! rc, and to preserve
|
|
|
|
!
|
|
|
|
! integral dr r**L * qfunc(r) ,
|
|
|
|
!
|
|
|
|
! i.e., to preserve the L'th moment of the charge. The array
|
|
|
|
! qfunc has been set inside rc to correspond to this pseudized
|
|
|
|
! version using the minimal L, namely L = | l1-l2 | (e.g., L=0
|
|
|
|
! for diagonal elements). The coefficients a_i (i=1,2,3)
|
|
|
|
! are stored in the array qfcoef(i,L+1,j,k) for each L so that
|
|
|
|
! the correctly pseudized versions of qfunc can be reconstructed
|
|
|
|
! for each L. (Note that for given l1 and l2, only the values
|
|
|
|
! L = | l1-l2 | , | l1-l2 | +2 , ... , l1+l2 are ever used.)
|
|
|
|
! ------------------------------------------------------
|
|
|
|
!
|
2006-12-14 16:53:47 +08:00
|
|
|
USE constants, ONLY : fpi
|
2007-10-03 00:54:13 +08:00
|
|
|
USE pseudo_types
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
implicit none
|
|
|
|
!
|
|
|
|
! First the arguments passed to the subroutine
|
2007-10-03 00:54:13 +08:00
|
|
|
!
|
|
|
|
TYPE (pseudo_upf) :: upf
|
2006-03-07 21:26:52 +08:00
|
|
|
integer &
|
|
|
|
& is, &! The number of the pseudopotential
|
|
|
|
& iunps ! The unit of the pseudo file
|
|
|
|
!
|
|
|
|
! Local variables
|
|
|
|
|
|
|
|
real(DP) &
|
|
|
|
& exfact, &! index of the exchange and correlation used
|
|
|
|
& etotpseu, &! total pseudopotential energy
|
|
|
|
& eloc, &! energy of the local potential
|
|
|
|
& dummy, &! dummy real variable
|
2007-10-03 00:54:13 +08:00
|
|
|
& rinner1, &! rinner if only one is present
|
2006-03-07 21:26:52 +08:00
|
|
|
& rcloc ! the cut-off radius of the local potential
|
2007-10-06 16:23:39 +08:00
|
|
|
real(DP), allocatable:: &
|
2007-10-24 03:47:26 +08:00
|
|
|
& ee(:), &! the energy of the valence states
|
|
|
|
& rc(:), &! the cut-off radii of the pseudopotential
|
2007-10-06 16:23:39 +08:00
|
|
|
& eee(:), &! energies of the beta function
|
|
|
|
& ddd(:,:) ! the screened D_{\mu,\nu} parameters
|
|
|
|
integer, allocatable :: &
|
2007-10-24 03:47:26 +08:00
|
|
|
& nnlz(:), &! The nlm values of the valence states
|
2007-10-06 16:23:39 +08:00
|
|
|
& iptype(:) ! more recent parameters
|
2006-03-07 21:26:52 +08:00
|
|
|
integer &
|
2007-10-08 15:33:57 +08:00
|
|
|
& iver(3), &! contains the version of generating code
|
2006-03-07 21:26:52 +08:00
|
|
|
& idmy(3), &! contains the date of creation of the pseudo
|
|
|
|
& ifpcor, &! for core correction, 0 otherwise
|
|
|
|
& ios, &! integer variable for I/O control
|
|
|
|
& i, &! dummy counter
|
|
|
|
& keyps, &! the type of pseudopotential. Only US allowed
|
2007-10-06 16:23:39 +08:00
|
|
|
& irel, &! says if the pseudopotential is relativistic
|
2007-10-03 00:54:13 +08:00
|
|
|
& ifqopt, &! level of Q optimization
|
2006-03-07 21:26:52 +08:00
|
|
|
& npf, &! as above
|
|
|
|
& nang, &! number of angular momenta in pseudopotentials
|
|
|
|
& lloc, &! angular momentum of the local part of PPs
|
|
|
|
& lp, &! counter on Q angular momenta
|
|
|
|
& l, &! counter on angular momenta
|
|
|
|
& iv, jv, ijv, &! beta function counter
|
|
|
|
& ir ! mesh points counter
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
character(len=20) :: title
|
2006-03-07 21:26:52 +08:00
|
|
|
character(len=60) fmt
|
|
|
|
!
|
|
|
|
! We first check the input variables
|
|
|
|
!
|
2007-10-17 04:49:00 +08:00
|
|
|
if (is <= 0) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore('readvan','routine called with wrong 1st argument', 1)
|
|
|
|
if (iunps <= 0 .or. iunps >= 100000) &
|
|
|
|
call errore('readvan','routine called with wrong 2nd argument', 1)
|
|
|
|
!
|
2013-06-01 18:33:31 +08:00
|
|
|
read(iunps, *, err=100, iostat=ios ) &
|
2007-10-08 15:33:57 +08:00
|
|
|
(iver(i),i=1,3), (idmy(i),i=1,3)
|
|
|
|
write(upf%generated, &
|
|
|
|
"('Generated by Vanderbilt code, v. ',i1,'.',i1,'.',i1)") iver
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-08 15:33:57 +08:00
|
|
|
if ( iver(1) > 7 .or. iver(1) < 1 .or. &
|
|
|
|
iver(2) > 9 .or. iver(2) < 0 .or. &
|
|
|
|
iver(3) > 9 .or. iver(3) < 0 ) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore('readvan','wrong file version read',1)
|
|
|
|
!
|
|
|
|
read( iunps, '(a20,3f15.9)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
title, upf%zmesh, upf%zp, exfact
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%psd = title(1:2)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( upf%zmesh < 1 .or. upf%zmesh > 100.0_DP) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore( 'readvan','wrong zmesh read', is )
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( upf%zp <= 0.0_DP .or. upf%zp > 100.0_DP) &
|
2007-03-19 03:24:56 +08:00
|
|
|
call errore('readvan','wrong atomic charge read', is )
|
2006-03-07 21:26:52 +08:00
|
|
|
if ( exfact < -6 .or. exfact > 6) &
|
|
|
|
& call errore('readvan','Wrong xc in pseudopotential',1)
|
|
|
|
! convert from "our" conventions to Vanderbilt conventions
|
2007-10-03 00:54:13 +08:00
|
|
|
call dftname_cp (nint(exfact), upf%dft)
|
|
|
|
call set_dft_from_name( upf%dft )
|
2006-03-07 21:26:52 +08:00
|
|
|
IF ( dft_is_meta() ) &
|
|
|
|
CALL errore( 'readvan ', 'META-GGA not implemented', 1 )
|
|
|
|
!
|
|
|
|
read( iunps, '(2i5,1pe19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%nwfc, upf%mesh, etotpseu
|
2007-10-24 03:47:26 +08:00
|
|
|
if ( upf%nwfc < 0 ) &
|
|
|
|
call errore( 'readvan', 'wrong nchi read', upf%nwfc )
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( upf%mesh < 0 ) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore( 'readvan','wrong mesh', is )
|
|
|
|
!
|
|
|
|
! info on pseudo eigenstates - energies are not used
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%oc(upf%nwfc), upf%lchi(upf%nwfc) )
|
2007-10-24 03:47:26 +08:00
|
|
|
ALLOCATE ( nnlz(upf%nwfc), ee(upf%nwfc) )
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(i5,2f15.9)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
( nnlz(iv), upf%oc(iv), ee(iv), iv=1,upf%nwfc )
|
|
|
|
do iv = 1, upf%nwfc
|
2006-03-07 21:26:52 +08:00
|
|
|
i = nnlz(iv) / 100
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%lchi(iv) = nnlz(iv)/10 - i * 10
|
2006-03-07 21:26:52 +08:00
|
|
|
enddo
|
|
|
|
read( iunps, '(2i5,f15.9)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
keyps, ifpcor, rinner1
|
|
|
|
upf%nlcc = (ifpcor == 1)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! keyps= 0 --> standard hsc pseudopotential with exponent 4.0
|
|
|
|
! 1 --> standard hsc pseudopotential with exponent 3.5
|
|
|
|
! 2 --> vanderbilt modifications using defaults
|
|
|
|
! 3 --> new generalized eigenvalue pseudopotentials
|
|
|
|
! 4 --> frozen core all-electron case
|
|
|
|
if ( keyps < 0 .or. keyps > 4 ) then
|
|
|
|
call errore('readvan','wrong keyps',keyps)
|
|
|
|
else if (keyps == 4) then
|
|
|
|
call errore('readvan','keyps not implemented',keyps)
|
|
|
|
end if
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%tvanp = (keyps == 3)
|
2008-01-16 02:06:09 +08:00
|
|
|
upf%tpawp = .false.
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! Read information on the angular momenta, and on Q pseudization
|
|
|
|
! (version > 3.0)
|
|
|
|
!
|
2007-10-08 15:33:57 +08:00
|
|
|
if (iver(1) >= 3) then
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(2i5,f9.5,2i5,f9.5)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
nang, lloc, eloc, ifqopt, upf%nqf, dummy
|
2006-03-07 21:26:52 +08:00
|
|
|
!!! PWSCF: lmax(is)=nang, lloc(is)=lloc
|
|
|
|
!
|
|
|
|
! NB: In the Vanderbilt atomic code the angular momentum goes
|
|
|
|
! from 1 to nang
|
|
|
|
!
|
2007-10-24 03:47:26 +08:00
|
|
|
if ( nang < 0 ) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore(' readvan', 'Wrong nang read', nang)
|
|
|
|
if ( lloc == -1 ) lloc = nang+1
|
|
|
|
if ( lloc > nang+1 .or. lloc < 0 ) &
|
|
|
|
call errore( 'readvan', 'wrong lloc read', is )
|
2007-10-06 16:23:39 +08:00
|
|
|
if ( upf%nqf < 0 ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
call errore(' readvan', 'Wrong nqf read', upf%nqf)
|
|
|
|
if ( ifqopt < 0 ) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore( 'readvan', 'wrong ifqopt read', is )
|
2007-10-17 04:49:00 +08:00
|
|
|
else
|
|
|
|
! old format: no distinction between nang and nchi
|
|
|
|
nang = upf%nwfc
|
2006-03-07 21:26:52 +08:00
|
|
|
end if
|
|
|
|
!
|
|
|
|
! Read and test the values of rinner (version > 5.1)
|
|
|
|
! rinner = radius at which to cut off partial core or q_ij
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%rinner(2*nang-1) )
|
2007-10-08 15:33:57 +08:00
|
|
|
if (10*iver(1)+iver(2) >= 51) then
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
read( iunps, *, err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
(upf%rinner(lp), lp=1,2*nang-1 )
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
do lp = 1, 2*nang-1
|
2007-10-03 00:54:13 +08:00
|
|
|
if (upf%rinner(lp) < 0.0_DP) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore('readvan','Wrong rinner read', is )
|
|
|
|
enddo
|
2007-10-08 15:33:57 +08:00
|
|
|
else if (iver(1) > 3) then
|
2006-03-07 21:26:52 +08:00
|
|
|
do lp = 2, 2*nang-1
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%rinner(lp)=rinner1
|
2006-03-07 21:26:52 +08:00
|
|
|
end do
|
|
|
|
end if
|
|
|
|
!
|
2007-10-08 15:33:57 +08:00
|
|
|
if (iver(1) >= 4) &
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(i5)',err=100, iostat=ios ) irel
|
|
|
|
!
|
|
|
|
! set the number of angular momentum terms in q_ij to read in
|
|
|
|
!
|
2007-10-08 15:33:57 +08:00
|
|
|
if (iver(1) == 1) then
|
2006-03-07 21:26:52 +08:00
|
|
|
oldvan(is) = .TRUE.
|
|
|
|
! old format: no optimization of q_ij => 3-term taylor series
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%nqf=3
|
|
|
|
upf%nqlc=5
|
2007-10-08 15:33:57 +08:00
|
|
|
else if (iver(1) == 2) then
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%nqf=3
|
|
|
|
upf%nqlc = 2*nang - 1
|
2006-03-07 21:26:52 +08:00
|
|
|
else
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%nqlc = 2*nang - 1
|
2006-03-07 21:26:52 +08:00
|
|
|
end if
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( upf%nqlc > lqmax .or. upf%nqlc < 0 ) &
|
|
|
|
call errore(' readvan', 'Wrong nqlc read', upf%nqlc )
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-24 03:47:26 +08:00
|
|
|
ALLOCATE ( rc(nang) )
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
|
|
|
( rc(l), l=1,nang )
|
|
|
|
!
|
|
|
|
! reads the number of beta functions
|
|
|
|
!
|
|
|
|
read( iunps, '(2i5)', err=100, iostat=ios ) &
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
upf%nbeta, upf%kkbeta
|
2007-10-03 00:54:13 +08:00
|
|
|
!
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
ALLOCATE ( upf%kbeta(upf%nbeta) )
|
|
|
|
upf%kbeta(:) = upf%kkbeta
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-06 16:23:39 +08:00
|
|
|
if( upf%nbeta < 0 ) &
|
|
|
|
call errore( 'readvan','nbeta wrong', is )
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
if( upf%kkbeta > upf%mesh .or. upf%kkbeta < 0 ) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore( 'readvan','kkbeta wrong or too large', is )
|
|
|
|
!
|
|
|
|
! Now reads the main Vanderbilt parameters
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%lll(upf%nbeta) )
|
|
|
|
ALLOCATE ( upf%beta(upf%mesh,upf%nbeta) )
|
|
|
|
ALLOCATE ( upf%dion(upf%nbeta,upf%nbeta), upf%qqq(upf%nbeta,upf%nbeta) )
|
2007-10-06 16:23:39 +08:00
|
|
|
ALLOCATE ( upf%qfunc(upf%mesh,upf%nbeta*(upf%nbeta+1)/2) )
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%qfcoef(upf%nqf, upf%nqlc, upf%nbeta, upf%nbeta) )
|
2007-10-06 16:23:39 +08:00
|
|
|
ALLOCATE ( eee(upf%nbeta), ddd(upf%nbeta,upf%nbeta) )
|
2007-10-03 00:54:13 +08:00
|
|
|
do iv=1,upf%nbeta
|
|
|
|
read( iunps, '(i5)',err=100, iostat=ios ) upf%lll(iv)
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)',err=100, iostat=ios ) &
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
eee(iv), ( upf%beta(ir,iv), ir=1,upf%kkbeta )
|
|
|
|
do ir=upf%kkbeta+1,upf%mesh
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%beta(ir,iv)=0.0_DP
|
|
|
|
enddo
|
|
|
|
if ( upf%lll(iv) > lmaxx .or. upf%lll(iv) < 0 ) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore( 'readvan', 'lll wrong or too large ', is )
|
2007-10-03 00:54:13 +08:00
|
|
|
do jv=iv,upf%nbeta
|
2006-11-03 20:21:28 +08:00
|
|
|
!
|
|
|
|
! the symmetric matric Q_{nb,mb} is stored in packed form
|
|
|
|
! Q(iv,jv) => qfunc(ijv) as defined below (for jv >= iv)
|
|
|
|
!
|
|
|
|
ijv = jv * (jv-1) / 2 + iv
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%dion(iv,jv), ddd(iv,jv), upf%qqq(iv,jv), &
|
2007-10-06 16:23:39 +08:00
|
|
|
(upf%qfunc(ir,ijv),ir=1,upf%kkbeta), &
|
2007-10-03 00:54:13 +08:00
|
|
|
((upf%qfcoef(i,lp,iv,jv),i=1,upf%nqf),lp=1,upf%nqlc)
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
do ir=upf%kkbeta+1,upf%mesh
|
2007-10-06 16:23:39 +08:00
|
|
|
upf%qfunc(ir,ijv)=0.0_DP
|
2007-10-03 00:54:13 +08:00
|
|
|
enddo
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! Use the symmetry of the coefficients
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( iv /= jv ) then
|
|
|
|
upf%dion(jv,iv)=upf%dion(iv,jv)
|
|
|
|
upf%qqq(jv,iv) =upf%qqq(iv,jv)
|
|
|
|
upf%qfcoef(:,:,jv,iv)=upf%qfcoef(:,:,iv,jv)
|
|
|
|
end if
|
2006-03-07 21:26:52 +08:00
|
|
|
enddo
|
|
|
|
enddo
|
This is a quite complex check-in, but actually not very much is done. Changelog follows.
LP
UPF file format updated completely, UPFv2 introduced:
* ld1.x can still produce old format, with the switch upf_v1_format=.true. in inputp
this is disabled by default, but we can discuss if it should be the opposite.
* pw.x cp.x and all utilities should notice no difference
* some utilities in upftools still need to be updated, anyway conversion UPFv1 to UPFv2
is very easy, so this should be no big issue
* starting from now to produce an UPF file you need to fill the pseudo_upf derivedd type
and feed it to write_upf woutine in upf_module (Modules/upf.f90)
* extensive use of iotk
I have tried to make the new format as self contained as possible, e.g. there should be
minimal need for post-processing after the data is read, no more reconstruction of known
quantities, and no more odd syntax to save negligible quantity of space. Also the human
readable section is a bit richer, all the rest is more machine readable.
I hope this will not cause any throuble, and tried really hard to, all examples and all
tests works as fine as before and gives (what really looks like) the same results.
Other changes that I needed to make:
* radial grids are now allocatable, they management is a bit less of a hack too
* paw and uspp augmentation are stored in the same place
* paw print total all-electron energy if all atoms are paw, not very useful, but nice
* most of the pseudopotential-writing reading files have been renamed to some more
logical name, I spare you the list. E.g. read_oldpseudo -> read_pseudo_rrkj3
* paw_t derived type was only used in atomic, so I have put it there (as the pseudo_type
module take ages to recompile it was awkward to leave it there).
PAW tests inserted in test/ there are 6 of them, as a consequence I have also put 5 paw
pseudopotentials in the pseudo/ directory.
I will update the PAW scf examples soon, by deleting them (as running a pw with a PAW
pseudopotential requires no option at all). PAW generation examples should be updated.
A lot of small bugfixes here & there mostly uninitialized variables or unallocated
pointers used as subrotuine arguments.
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4769 c92efa57-630b-4861-b058-cf58834340f0
2008-04-03 23:50:43 +08:00
|
|
|
!
|
|
|
|
! Set additional, not present, variables to dummy values
|
|
|
|
ALLOCATE(upf%els(upf%nwfc))
|
|
|
|
upf%els(:) = 'nX'
|
|
|
|
ALLOCATE(upf%els_beta(upf%nbeta))
|
|
|
|
upf%els_beta(:) = 'nX'
|
|
|
|
ALLOCATE(upf%rcut(upf%nbeta), upf%rcutus(upf%nbeta))
|
|
|
|
upf%rcut(:) = 0._dp
|
|
|
|
upf%rcutus(:) = 0._dp
|
|
|
|
|
2007-10-06 16:23:39 +08:00
|
|
|
DEALLOCATE (ddd)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! for versions later than 7.2
|
|
|
|
!
|
2007-10-08 15:33:57 +08:00
|
|
|
if (10*iver(1)+iver(2) >= 72) then
|
2007-10-06 16:23:39 +08:00
|
|
|
ALLOCATE (iptype(upf%nbeta))
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(6i5)',err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
(iptype(iv), iv=1,upf%nbeta)
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(i5,f15.9)',err=100, iostat=ios ) &
|
|
|
|
npf, dummy
|
2007-10-06 16:23:39 +08:00
|
|
|
DEALLOCATE (iptype)
|
2006-03-07 21:26:52 +08:00
|
|
|
end if
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
! read the local potential
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%vloc(upf%mesh) )
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)',err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
rcloc, ( upf%vloc(ir), ir=1,upf%mesh )
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! If present reads the core charge rho_atc(r)=4*pi*r**2*rho_core(r)
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( upf%nlcc ) then
|
|
|
|
ALLOCATE ( upf%rho_atc(upf%mesh) )
|
2007-10-08 15:33:57 +08:00
|
|
|
if (iver(1) >= 7) &
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) dummy
|
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
( upf%rho_atc(ir), ir=1,upf%mesh )
|
2006-03-07 21:26:52 +08:00
|
|
|
endif
|
|
|
|
!
|
|
|
|
! Read the screened local potential (not used)
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%rho_at(upf%mesh) )
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
(upf%rho_at(ir), ir=1,upf%mesh)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! Read the valence atomic charge
|
|
|
|
!
|
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
(upf%rho_at(ir), ir=1,upf%mesh)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! Read the logarithmic mesh (if version > 1)
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%r(upf%mesh), upf%rab(upf%mesh) )
|
2007-10-08 15:33:57 +08:00
|
|
|
if (iver(1) >1) then
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)',err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
(upf%r(ir),ir=1,upf%mesh)
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)',err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
(upf%rab(ir),ir=1,upf%mesh)
|
2006-03-07 21:26:52 +08:00
|
|
|
else
|
|
|
|
!
|
|
|
|
! generate herman-skillman mesh (if version = 1)
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
call herman_skillman_grid &
|
|
|
|
( upf%mesh, upf%zmesh, upf%r, upf%rab )
|
2006-03-07 21:26:52 +08:00
|
|
|
end if
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
! convert vloc to the conventions used in the rest of the code
|
2006-03-07 21:26:52 +08:00
|
|
|
! (as read from Vanderbilt's format it is r*v_loc(r))
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
do ir = 2, upf%mesh
|
|
|
|
upf%vloc (ir) = upf%vloc (ir) / upf%r(ir)
|
2006-03-07 21:26:52 +08:00
|
|
|
enddo
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%vloc (1) = upf%vloc (2)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! set rho_atc(r)=rho_core(r) (without 4*pi*r^2 factor,
|
|
|
|
! for compatibility with rho_atc in the non-US case)
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
if (upf%nlcc) then
|
|
|
|
upf%rho_atc(1) = 0.0_DP
|
|
|
|
do ir=2,upf%mesh
|
|
|
|
upf%rho_atc(ir) = upf%rho_atc(ir)/fpi/upf%r(ir)**2
|
2006-03-07 21:26:52 +08:00
|
|
|
enddo
|
|
|
|
end if
|
|
|
|
!
|
|
|
|
! Read the wavefunctions of the atom
|
2007-10-03 00:54:13 +08:00
|
|
|
!
|
2007-10-08 15:33:57 +08:00
|
|
|
if (iver(1) >= 7) then
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, *, err=100, iostat=ios ) i
|
2007-10-03 00:54:13 +08:00
|
|
|
if (i /= upf%nwfc) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore('readvan','unexpected or unimplemented case',1)
|
|
|
|
end if
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%chi(upf%mesh, upf%nwfc) )
|
2007-10-08 15:33:57 +08:00
|
|
|
if (iver(1) >= 6) &
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, *, err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
( (upf%chi(ir,iv), ir=1,upf%mesh), iv=1,upf%nwfc )
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-08 15:33:57 +08:00
|
|
|
if (iver(1) == 1) then
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2006-03-15 19:23:03 +08:00
|
|
|
! old version: read the q_l(r) and fit them with the Vanderbilt's form
|
2007-10-03 00:54:13 +08:00
|
|
|
!
|
|
|
|
call fit_qrl ( )
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
end if
|
|
|
|
!
|
|
|
|
! Here we write on output information on the pseudopotential
|
|
|
|
!
|
|
|
|
WRITE( stdout,200) is
|
|
|
|
200 format (/4x,60('=')/4x,'| pseudopotential report', &
|
|
|
|
& ' for atomic species:',i3,11x,'|')
|
|
|
|
WRITE( stdout,300) 'pseudo potential version', &
|
2007-10-08 15:33:57 +08:00
|
|
|
iver(1), iver(2), iver(3)
|
2006-03-07 21:26:52 +08:00
|
|
|
300 format (4x,'| ',1a30,3i4,13x,' |' /4x,60('-'))
|
2007-10-03 00:54:13 +08:00
|
|
|
WRITE( stdout,400) title, upf%dft
|
2006-03-07 21:26:52 +08:00
|
|
|
400 format (4x,'| ',2a20,' exchange-corr |')
|
2007-10-03 00:54:13 +08:00
|
|
|
WRITE( stdout,500) upf%zmesh, is, upf%zp, exfact
|
2006-03-07 21:26:52 +08:00
|
|
|
500 format (4x,'| z =',f5.0,4x,'zv(',i2,') =',f5.0,4x,'exfact =', &
|
|
|
|
& f10.5, 9x,'|')
|
|
|
|
WRITE( stdout,600) ifpcor, etotpseu
|
|
|
|
600 format (4x,'| ifpcor = ',i2,10x,' atomic energy =',f10.5, &
|
|
|
|
& ' Ry',6x,'|')
|
|
|
|
WRITE( stdout,700)
|
|
|
|
700 format(4x,'| index orbital occupation energy',14x,'|')
|
2007-10-03 00:54:13 +08:00
|
|
|
WRITE( stdout,800) ( iv, nnlz(iv), upf%oc(iv), ee(iv), iv=1,upf%nwfc )
|
2007-10-24 03:47:26 +08:00
|
|
|
DEALLOCATE (ee, nnlz)
|
2006-03-07 21:26:52 +08:00
|
|
|
800 format(4x,'|',i5,i11,5x,f10.2,f12.2,15x,'|')
|
2007-10-08 15:33:57 +08:00
|
|
|
if (iver(1) >= 3 .and. nang > 0) then
|
2013-06-04 21:49:32 +08:00
|
|
|
IF (nang < 4) THEN
|
|
|
|
write(fmt,900) 2*nang-1, 40-8*(2*nang-2)
|
|
|
|
ELSE
|
|
|
|
write(fmt,900) 2*nang-1, 1
|
|
|
|
ENDIF
|
2010-12-08 00:17:49 +08:00
|
|
|
900 format('(4x,"| rinner =",',i1,'f8.4,',i2,'x,"|")')
|
2007-10-03 00:54:13 +08:00
|
|
|
WRITE( stdout,fmt) (upf%rinner(lp),lp=1,2*nang-1)
|
2006-03-07 21:26:52 +08:00
|
|
|
end if
|
|
|
|
WRITE( stdout,1000)
|
|
|
|
1000 format(4x,'| new generation scheme:',32x,'|')
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
WRITE( stdout,1100) upf%nbeta, upf%kkbeta, rcloc
|
2007-10-03 00:54:13 +08:00
|
|
|
1100 format(4x,'| nbeta = ',i2,5x,'kkbeta =',i5,5x,'rcloc =',f10.4,4x,&
|
|
|
|
& '|'/4x,'| ibeta l epsilon rcut',25x,'|')
|
|
|
|
do iv = 1, upf%nbeta
|
|
|
|
lp=upf%lll(iv)+1
|
|
|
|
WRITE( stdout,1200) iv,upf%lll(iv),eee(iv),rc(lp)
|
2006-03-07 21:26:52 +08:00
|
|
|
1200 format(4x,'|',5x,i2,6x,i2,4x,2f7.2,25x,'|')
|
|
|
|
enddo
|
|
|
|
WRITE( stdout,1300)
|
|
|
|
1300 format (4x,60('='))
|
|
|
|
!
|
2007-10-24 03:47:26 +08:00
|
|
|
DEALLOCATE (eee, rc)
|
2006-03-07 21:26:52 +08:00
|
|
|
return
|
|
|
|
100 call errore('readvan','error reading pseudo file', abs(ios) )
|
2007-10-03 00:54:13 +08:00
|
|
|
!
|
|
|
|
CONTAINS
|
2006-03-07 21:26:52 +08:00
|
|
|
!-----------------------------------------------------------------------
|
2007-10-03 00:54:13 +08:00
|
|
|
subroutine fit_qrl ( )
|
2006-03-07 21:26:52 +08:00
|
|
|
!-----------------------------------------------------------------------
|
|
|
|
!
|
|
|
|
! find coefficients qfcoef that fit the pseudized qrl in US PP
|
|
|
|
! these coefficients are written to file in newer versions of the
|
|
|
|
! Vanderbilt PP generation code but not in some ancient versions
|
|
|
|
!
|
|
|
|
implicit none
|
|
|
|
!
|
2006-03-15 19:23:03 +08:00
|
|
|
real (kind=DP), allocatable :: qrl(:,:), a(:,:), ainv(:,:), b(:), x(:)
|
2006-03-07 21:26:52 +08:00
|
|
|
real (kind=DP) :: deta
|
|
|
|
integer :: iv, jv, ijv, lmin, lmax, l, ir, irinner, i,j
|
|
|
|
!
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
allocate ( a(upf%nqf,upf%nqf), ainv(upf%nqf,upf%nqf) )
|
|
|
|
allocate ( b(upf%nqf), x(upf%nqf) )
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
ALLOCATE ( qrl(upf%kkbeta, upf%nqlc) )
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
do iv=1,upf%nbeta
|
|
|
|
do jv=iv,upf%nbeta
|
2006-03-15 19:23:03 +08:00
|
|
|
!
|
|
|
|
! original version, assuming lll(jv) >= lll(iv)
|
|
|
|
! lmin=lll(jv,is)-lll(iv,is)+1
|
|
|
|
! lmax=lmin+2*lll(iv,is)
|
|
|
|
! note that indices run from 1 to Lmax+1, not from 0 to Lmax
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
lmin = ABS( upf%lll(jv) - upf%lll(iv) ) + 1
|
|
|
|
lmax = upf%lll(jv) + upf%lll(iv) + 1
|
2006-03-15 19:23:03 +08:00
|
|
|
IF ( lmin < 1 .OR. lmax > SIZE(qrl,2)) &
|
|
|
|
CALL errore ('fit_qrl', 'bad 2rd dimension for array qrl', 1)
|
|
|
|
!
|
|
|
|
! read q_l(r) for all l
|
|
|
|
!
|
|
|
|
read(iunps,*, err=100) &
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
( (qrl(ir,l),ir=1,upf%kkbeta), l=lmin,lmax)
|
2007-10-03 00:54:13 +08:00
|
|
|
!
|
2007-10-06 16:23:39 +08:00
|
|
|
ijv = jv * (jv-1) / 2 + iv
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
do l=lmin,lmax
|
|
|
|
!
|
|
|
|
! reconstruct rinner
|
|
|
|
!
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
do ir=upf%kkbeta,1,-1
|
2007-10-06 16:23:39 +08:00
|
|
|
if ( abs(qrl(ir,l)-upf%qfunc(ir,ijv)) > 1.0d-6) go to 10
|
2006-03-07 21:26:52 +08:00
|
|
|
end do
|
|
|
|
10 irinner = ir+1
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%rinner(l) = upf%r(irinner)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! least square minimization: find
|
|
|
|
! qrl = sum_i c_i r^{l+1}r^{2i-2} for r < rinner
|
|
|
|
!
|
2007-06-12 01:13:15 +08:00
|
|
|
a(:,:) = 0.0_DP
|
|
|
|
b(:) = 0.0_DP
|
2007-10-03 00:54:13 +08:00
|
|
|
do i = 1, upf%nqf
|
2006-03-07 21:26:52 +08:00
|
|
|
do ir=1,irinner
|
2007-10-03 00:54:13 +08:00
|
|
|
b(i) = b(i) + upf%r(ir)**(2*i-2+l+1) * qrl(ir,l)
|
2006-03-07 21:26:52 +08:00
|
|
|
end do
|
2007-10-03 00:54:13 +08:00
|
|
|
do j = i, upf%nqf
|
2006-03-07 21:26:52 +08:00
|
|
|
do ir=1,irinner
|
2007-10-03 00:54:13 +08:00
|
|
|
a(i,j) = a(i,j) + upf%r(ir)**(2*i-2+l+1) * &
|
|
|
|
upf%r(ir)**(2*j-2+l+1)
|
2006-03-07 21:26:52 +08:00
|
|
|
end do
|
2006-03-15 19:23:03 +08:00
|
|
|
if (j > i) a(j,i) = a(i,j)
|
2006-03-07 21:26:52 +08:00
|
|
|
end do
|
|
|
|
end do
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
call invmat (upf%nqf, a, ainv, deta)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
do i = 1, upf%nqf
|
|
|
|
upf%qfcoef(i,l,iv,jv) = dot_product(ainv(i,:),b(:))
|
|
|
|
if (iv /= jv) upf%qfcoef(i,l,jv,iv) = upf%qfcoef(i,l,iv,jv)
|
2006-03-07 21:26:52 +08:00
|
|
|
end do
|
|
|
|
end do
|
|
|
|
end do
|
|
|
|
end do
|
|
|
|
!
|
2006-03-15 19:23:03 +08:00
|
|
|
deallocate ( qrl, x, b , ainv, a )
|
|
|
|
return
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2006-03-15 19:23:03 +08:00
|
|
|
100 call errore('readvan','error reading Q_L(r)', 1 )
|
2006-03-07 21:26:52 +08:00
|
|
|
end subroutine fit_qrl
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
end subroutine readvan
|
2006-03-07 21:26:52 +08:00
|
|
|
!-----------------------------------------------------------------------
|
2007-10-03 00:54:13 +08:00
|
|
|
SUBROUTINE herman_skillman_grid (mesh,z,r,rab)
|
2006-03-07 21:26:52 +08:00
|
|
|
!-----------------------------------------------------------------------
|
|
|
|
!
|
|
|
|
! generate Herman-Skillman radial grid (obsolescent)
|
2007-06-12 01:13:15 +08:00
|
|
|
! c - 0.88534138/z**(1/3)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
IMPLICIT NONE
|
|
|
|
!
|
|
|
|
INTEGER mesh
|
2007-06-12 01:13:15 +08:00
|
|
|
REAL(DP) :: z, r(mesh), rab(mesh)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-06-12 01:13:15 +08:00
|
|
|
REAL(DP) :: deltax,pi
|
|
|
|
INTEGER :: nblock,i,j,k
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-06-12 01:13:15 +08:00
|
|
|
pi=4.0_DP*ATAN(1.0_DP)
|
2006-03-07 21:26:52 +08:00
|
|
|
nblock = mesh/40
|
|
|
|
i=1
|
2007-06-12 01:13:15 +08:00
|
|
|
r(i)=0.0_DP
|
|
|
|
deltax=0.0025_DP*0.5_DP*(3.0_DP*pi/4.0_DP)**(2.0_DP/3.0_DP)/z**(1.0_DP/3.0_DP)
|
2006-03-07 21:26:52 +08:00
|
|
|
DO j=1,nblock
|
|
|
|
DO k=1,40
|
|
|
|
i=i+1
|
|
|
|
r(i)=r(i-1)+deltax
|
|
|
|
rab(i)=deltax
|
|
|
|
END DO
|
|
|
|
deltax=deltax+deltax
|
|
|
|
END DO
|
|
|
|
!
|
|
|
|
RETURN
|
|
|
|
END SUBROUTINE herman_skillman_grid
|
|
|
|
!
|
|
|
|
!---------------------------------------------------------------------
|
2007-10-03 00:54:13 +08:00
|
|
|
subroutine readrrkj( iunps, is, upf )
|
2006-03-07 21:26:52 +08:00
|
|
|
!---------------------------------------------------------------------
|
|
|
|
!
|
|
|
|
! This routine reads Vanderbilt pseudopotentials produced by the
|
|
|
|
! code of Andrea Dal Corso. Hard PPs are first generated
|
|
|
|
! according to the Rabe Rappe Kaxiras Johannopoulos recipe.
|
|
|
|
! Ultrasoft PP's are subsequently generated from the hard PP's.
|
|
|
|
!
|
|
|
|
! Output parameters in module "uspp_param"
|
|
|
|
! info on DFT level in module "dft"
|
|
|
|
!
|
2006-12-14 16:53:47 +08:00
|
|
|
USE constants, ONLY : fpi
|
2007-10-03 00:54:13 +08:00
|
|
|
USE pseudo_types
|
2006-12-14 16:53:47 +08:00
|
|
|
!
|
2006-03-07 21:26:52 +08:00
|
|
|
implicit none
|
|
|
|
!
|
|
|
|
! First the arguments passed to the subroutine
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
TYPE (pseudo_upf) :: upf
|
2006-03-07 21:26:52 +08:00
|
|
|
integer :: &
|
|
|
|
is, &! The index of the pseudopotential
|
|
|
|
iunps ! the unit from with pseudopotential is read
|
|
|
|
!
|
|
|
|
! Local variables
|
|
|
|
!
|
|
|
|
integer:: iexch, icorr, igcx, igcc
|
|
|
|
|
|
|
|
integer:: &
|
2006-11-03 20:21:28 +08:00
|
|
|
nb,mb, ijv,&! counters on beta functions
|
2006-03-07 21:26:52 +08:00
|
|
|
n, &! counter on mesh points
|
|
|
|
ir, &! counters on mesh points
|
|
|
|
pseudotype,&! the type of pseudopotential
|
|
|
|
ios, &! I/O control
|
|
|
|
ndum, &! dummy integer variable
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
l ! counter on angular momentum
|
2006-03-07 21:26:52 +08:00
|
|
|
real(DP):: &
|
|
|
|
x, &! auxiliary variable
|
|
|
|
etotps, &! total energy of the pseudoatom
|
|
|
|
rdum ! dummy real variable
|
|
|
|
!
|
|
|
|
logical :: &
|
2006-09-12 20:55:30 +08:00
|
|
|
rel ! if true the atomic calculation is relativistic
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
character(len=75) :: &
|
2006-09-12 20:55:30 +08:00
|
|
|
titleps ! the title of the pseudo
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
integer :: &
|
|
|
|
lmax ! max angular momentum
|
2006-09-12 20:55:30 +08:00
|
|
|
character(len=2) :: &
|
|
|
|
adum ! dummy character variable
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! We first check the input variables
|
|
|
|
!
|
2007-10-17 04:49:00 +08:00
|
|
|
if (is <= 0) &
|
2007-06-21 23:27:34 +08:00
|
|
|
call errore('readrrkj','routine called with wrong 1st argument', 1)
|
2006-03-07 21:26:52 +08:00
|
|
|
if (iunps <= 0 .or. iunps >= 100000) &
|
2007-06-21 23:27:34 +08:00
|
|
|
call errore('readrrkj','routine called with wrong 2nd argument', 1)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
read( iunps, '(a75)', err=100, iostat=ios ) &
|
|
|
|
titleps
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%psd = titleps(7:8)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
read( iunps, '(i5)',err=100, iostat=ios ) &
|
|
|
|
pseudotype
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%tvanp = (pseudotype == 3)
|
2008-01-16 02:06:09 +08:00
|
|
|
upf%tpawp = .false.
|
2007-10-08 15:33:57 +08:00
|
|
|
|
|
|
|
if ( upf%tvanp ) then
|
|
|
|
upf%generated = &
|
|
|
|
"RRKJ3 Ultrasoft PP, generated by Andrea Dal Corso code"
|
2006-03-07 21:26:52 +08:00
|
|
|
else
|
2007-10-08 15:33:57 +08:00
|
|
|
upf%generated = &
|
|
|
|
"RRKJ3 norm-conserving PP, generated by Andrea Dal Corso code"
|
2006-03-07 21:26:52 +08:00
|
|
|
endif
|
This is a quite complex check-in, but actually not very much is done. Changelog follows.
LP
UPF file format updated completely, UPFv2 introduced:
* ld1.x can still produce old format, with the switch upf_v1_format=.true. in inputp
this is disabled by default, but we can discuss if it should be the opposite.
* pw.x cp.x and all utilities should notice no difference
* some utilities in upftools still need to be updated, anyway conversion UPFv1 to UPFv2
is very easy, so this should be no big issue
* starting from now to produce an UPF file you need to fill the pseudo_upf derivedd type
and feed it to write_upf woutine in upf_module (Modules/upf.f90)
* extensive use of iotk
I have tried to make the new format as self contained as possible, e.g. there should be
minimal need for post-processing after the data is read, no more reconstruction of known
quantities, and no more odd syntax to save negligible quantity of space. Also the human
readable section is a bit richer, all the rest is more machine readable.
I hope this will not cause any throuble, and tried really hard to, all examples and all
tests works as fine as before and gives (what really looks like) the same results.
Other changes that I needed to make:
* radial grids are now allocatable, they management is a bit less of a hack too
* paw and uspp augmentation are stored in the same place
* paw print total all-electron energy if all atoms are paw, not very useful, but nice
* most of the pseudopotential-writing reading files have been renamed to some more
logical name, I spare you the list. E.g. read_oldpseudo -> read_pseudo_rrkj3
* paw_t derived type was only used in atomic, so I have put it there (as the pseudo_type
module take ages to recompile it was awkward to leave it there).
PAW tests inserted in test/ there are 6 of them, as a consequence I have also put 5 paw
pseudopotentials in the pseudo/ directory.
I will update the PAW scf examples soon, by deleting them (as running a pw with a PAW
pseudopotential requires no option at all). PAW generation examples should be updated.
A lot of small bugfixes here & there mostly uninitialized variables or unallocated
pointers used as subrotuine arguments.
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4769 c92efa57-630b-4861-b058-cf58834340f0
2008-04-03 23:50:43 +08:00
|
|
|
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(2l5)',err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
rel, upf%nlcc
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(4i5)',err=100, iostat=ios ) &
|
|
|
|
iexch, icorr, igcx, igcc
|
2007-10-03 00:54:13 +08:00
|
|
|
!
|
|
|
|
! workaround to keep track of which dft was read
|
|
|
|
! See also upf2internals
|
|
|
|
!
|
|
|
|
write( upf%dft, "('INDEX:',4i1)") iexch,icorr,igcx,igcc
|
2011-04-27 23:18:18 +08:00
|
|
|
call set_dft_from_indices(iexch,icorr,igcx,igcc, 0) ! Cannot read nonlocal in this format
|
2006-03-07 21:26:52 +08:00
|
|
|
|
|
|
|
read( iunps, '(2e17.11,i5)') &
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%zp, etotps, lmax
|
|
|
|
if ( upf%zp < 1 .or. upf%zp > 100 ) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore('readrrkj','wrong potential read',is)
|
|
|
|
!
|
|
|
|
read( iunps, '(4e17.11,i5)',err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%xmin, rdum, upf%zmesh, upf%dx, upf%mesh
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( upf%mesh < 0) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore('readrrkj', 'wrong mesh',is)
|
|
|
|
!
|
|
|
|
read( iunps, '(2i5)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%nwfc, upf%nbeta
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-06 16:23:39 +08:00
|
|
|
if ( upf%nbeta < 0) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore('readrrkj', 'wrong nbeta', is)
|
2007-10-24 03:47:26 +08:00
|
|
|
if ( upf%nwfc < 0 ) &
|
2006-03-07 21:26:52 +08:00
|
|
|
call errore('readrrkj', 'wrong nchi', is)
|
|
|
|
!
|
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
( rdum, nb=1,upf%nwfc )
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
( rdum, nb=1,upf%nwfc )
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%oc(upf%nwfc), upf%lchi(upf%nwfc), upf%lll(upf%nwfc) )
|
|
|
|
!
|
|
|
|
do nb=1,upf%nwfc
|
2006-03-07 21:26:52 +08:00
|
|
|
read(iunps,'(a2,2i3,f6.2)',err=100,iostat=ios) &
|
2007-10-03 00:54:13 +08:00
|
|
|
adum, ndum, upf%lchi(nb), upf%oc(nb)
|
|
|
|
upf%lll(nb)=upf%lchi(nb)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! oc < 0 distinguishes between bound states from unbound states
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( upf%oc(nb) <= 0.0_DP) upf%oc(nb) = -1.0_DP
|
2006-03-07 21:26:52 +08:00
|
|
|
enddo
|
|
|
|
!
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
ALLOCATE ( upf%kbeta(upf%nbeta) )
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%dion(upf%nbeta,upf%nbeta), upf%qqq(upf%nbeta,upf%nbeta) )
|
|
|
|
ALLOCATE ( upf%beta(upf%mesh,upf%nbeta) )
|
2007-10-06 16:23:39 +08:00
|
|
|
ALLOCATE ( upf%qfunc(upf%mesh,upf%nbeta*(upf%nbeta+1)/2) )
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
upf%kkbeta = 0
|
2007-10-03 00:54:13 +08:00
|
|
|
do nb=1,upf%nbeta
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
read ( iunps, '(i6)',err=100, iostat=ios ) upf%kbeta(nb)
|
|
|
|
upf%kkbeta = MAX ( upf%kkbeta, upf%kbeta(nb) )
|
2006-03-07 21:26:52 +08:00
|
|
|
read ( iunps, '(1p4e19.11)',err=100, iostat=ios ) &
|
The following pseudopotential-related variables in module uspp_param:
zp, psd, dion, betar, jjj, qqq, qfunc, qfcoef, vloc_at, rinner,
nbeta, kkbeta, nqf, nqlc, lll, tvanp
have been replaced by the corresponding variables in structure 'upf'.
There shouldn't be any side effects, but who knows. There is still a
copy of the above variables that will be removed sooner or later.
Basically : variable([i,j,k,..,]n) => upf(n)%variable [(i,j,k,..)]
Note that upf%qfunc has for the time being three indices instead of two,
and that upf%kkbeta is the analogous of kkbeta and not what it used to be.
The logic of this operation will be clearer when it will be completed
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4308 c92efa57-630b-4861-b058-cf58834340f0
2007-10-05 17:26:23 +08:00
|
|
|
( upf%beta(ir,nb), ir=1,upf%kbeta(nb))
|
|
|
|
do ir=upf%kbeta(nb)+1,upf%mesh
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%beta(ir,nb)=0.0_DP
|
2006-03-07 21:26:52 +08:00
|
|
|
enddo
|
|
|
|
do mb=1,nb
|
2006-11-03 20:21:28 +08:00
|
|
|
!
|
|
|
|
! the symmetric matric Q_{nb,mb} is stored in packed form
|
|
|
|
! Q(nb,mb) => qfunc(ijv) as defined below (for mb <= nb)
|
2007-10-06 16:23:39 +08:00
|
|
|
!
|
2006-11-03 20:21:28 +08:00
|
|
|
ijv = nb * (nb - 1) / 2 + mb
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%dion(nb,mb)
|
|
|
|
if (pseudotype == 3) then
|
2006-03-07 21:26:52 +08:00
|
|
|
read(iunps,'(1p4e19.11)',err=100,iostat=ios) &
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%qqq(nb,mb)
|
2006-03-07 21:26:52 +08:00
|
|
|
read(iunps,'(1p4e19.11)',err=100,iostat=ios) &
|
2007-10-06 16:23:39 +08:00
|
|
|
(upf%qfunc(n,ijv),n=1,upf%mesh)
|
2006-03-07 21:26:52 +08:00
|
|
|
else
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%qqq(nb,mb)=0.0_DP
|
2007-10-06 16:23:39 +08:00
|
|
|
upf%qfunc(:,ijv)=0.0_DP
|
2006-03-07 21:26:52 +08:00
|
|
|
endif
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( mb /= nb ) then
|
|
|
|
upf%dion(mb,nb)=upf%dion(nb,mb)
|
|
|
|
upf%qqq(mb,nb)=upf%qqq(nb,mb)
|
|
|
|
end if
|
2006-03-07 21:26:52 +08:00
|
|
|
enddo
|
|
|
|
enddo
|
|
|
|
!
|
|
|
|
! reads the local potential
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%vloc(upf%mesh) )
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)',err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
rdum, ( upf%vloc(ir), ir=1,upf%mesh )
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! reads the atomic charge
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%rho_at(upf%mesh) )
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
( upf%rho_at(ir), ir=1,upf%mesh )
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! if present reads the core charge
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( upf%nlcc ) then
|
|
|
|
ALLOCATE ( upf%rho_atc(upf%mesh) )
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
( upf%rho_atc(ir), ir=1,upf%mesh )
|
2006-03-07 21:26:52 +08:00
|
|
|
endif
|
|
|
|
!
|
|
|
|
! read the pseudo wavefunctions of the atom
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%chi(upf%mesh, upf%nwfc) )
|
2006-03-07 21:26:52 +08:00
|
|
|
read( iunps, '(1p4e19.11)', err=100, iostat=ios ) &
|
2007-10-03 00:54:13 +08:00
|
|
|
((upf%chi(ir,nb),ir=1,upf%mesh),nb=1,upf%nwfc)
|
2006-03-07 21:26:52 +08:00
|
|
|
!
|
|
|
|
! set several variables for compatibility with the rest of the code
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
upf%nqf=0
|
|
|
|
upf%nqlc=2*lmax+1
|
|
|
|
if ( upf%nqlc > lqmax .or. upf%nqlc < 0 ) &
|
|
|
|
call errore(' readrrkj', 'Wrong nqlc', upf%nqlc )
|
|
|
|
ALLOCATE ( upf%rinner(upf%nqlc) )
|
|
|
|
do l=1,upf%nqlc
|
|
|
|
upf%rinner(l)=0.0_DP
|
2006-03-07 21:26:52 +08:00
|
|
|
enddo
|
|
|
|
!
|
|
|
|
! compute the radial mesh
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
ALLOCATE ( upf%r(upf%mesh), upf%rab(upf%mesh) )
|
|
|
|
do ir = 1, upf%mesh
|
|
|
|
x = upf%xmin + DBLE(ir-1) * upf%dx
|
|
|
|
upf%r(ir) = EXP(x) / upf%zmesh
|
|
|
|
upf%rab(ir) = upf%dx * upf%r(ir)
|
2006-03-07 21:26:52 +08:00
|
|
|
end do
|
|
|
|
!
|
|
|
|
! set rho_atc(r)=rho_core(r) (without 4*pi*r^2 factor)
|
|
|
|
!
|
2007-10-03 00:54:13 +08:00
|
|
|
if ( upf%nlcc ) then
|
|
|
|
do ir=1,upf%mesh
|
|
|
|
upf%rho_atc(ir) = upf%rho_atc(ir)/fpi/upf%r(ir)**2
|
2006-03-07 21:26:52 +08:00
|
|
|
enddo
|
|
|
|
end if
|
|
|
|
!
|
This is a quite complex check-in, but actually not very much is done. Changelog follows.
LP
UPF file format updated completely, UPFv2 introduced:
* ld1.x can still produce old format, with the switch upf_v1_format=.true. in inputp
this is disabled by default, but we can discuss if it should be the opposite.
* pw.x cp.x and all utilities should notice no difference
* some utilities in upftools still need to be updated, anyway conversion UPFv1 to UPFv2
is very easy, so this should be no big issue
* starting from now to produce an UPF file you need to fill the pseudo_upf derivedd type
and feed it to write_upf woutine in upf_module (Modules/upf.f90)
* extensive use of iotk
I have tried to make the new format as self contained as possible, e.g. there should be
minimal need for post-processing after the data is read, no more reconstruction of known
quantities, and no more odd syntax to save negligible quantity of space. Also the human
readable section is a bit richer, all the rest is more machine readable.
I hope this will not cause any throuble, and tried really hard to, all examples and all
tests works as fine as before and gives (what really looks like) the same results.
Other changes that I needed to make:
* radial grids are now allocatable, they management is a bit less of a hack too
* paw and uspp augmentation are stored in the same place
* paw print total all-electron energy if all atoms are paw, not very useful, but nice
* most of the pseudopotential-writing reading files have been renamed to some more
logical name, I spare you the list. E.g. read_oldpseudo -> read_pseudo_rrkj3
* paw_t derived type was only used in atomic, so I have put it there (as the pseudo_type
module take ages to recompile it was awkward to leave it there).
PAW tests inserted in test/ there are 6 of them, as a consequence I have also put 5 paw
pseudopotentials in the pseudo/ directory.
I will update the PAW scf examples soon, by deleting them (as running a pw with a PAW
pseudopotential requires no option at all). PAW generation examples should be updated.
A lot of small bugfixes here & there mostly uninitialized variables or unallocated
pointers used as subrotuine arguments.
git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@4769 c92efa57-630b-4861-b058-cf58834340f0
2008-04-03 23:50:43 +08:00
|
|
|
! Set additional, not present, variables to dummy values
|
|
|
|
allocate(upf%els(upf%nwfc))
|
|
|
|
upf%els(:) = 'nX'
|
|
|
|
allocate(upf%els_beta(upf%nbeta))
|
|
|
|
upf%els_beta(:) = 'nX'
|
|
|
|
allocate(upf%rcut(upf%nbeta), upf%rcutus(upf%nbeta))
|
|
|
|
upf%rcut(:) = 0._dp
|
|
|
|
upf%rcutus(:) = 0._dp
|
|
|
|
!
|
2006-03-07 21:26:52 +08:00
|
|
|
return
|
|
|
|
100 call errore('readrrkj','Reading pseudo file',abs(ios))
|
|
|
|
end subroutine readrrkj
|
2007-10-03 00:54:13 +08:00
|
|
|
!
|
2006-03-07 21:26:52 +08:00
|
|
|
end module read_uspp_module
|