TABLE OF CONTENTS


ABINIT/m_pspheads [ Modules ]

[ Top ] [ Modules ]

NAME

 m_pspheads

FUNCTION

  Functions used to read the pseudopotential header of each psp file, in order to initialize pspheads(1:npsp).

COPYRIGHT

  Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR, FrD, AF, MT, FJ, MJV)
  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_pspheads
28 
29  use defs_basis
30  use defs_datatypes
31  use m_abicore
32  use m_errors
33  use m_hash_md5
34  use m_psxml2ab
35 #if defined HAVE_PSML
36  use m_psml
37 #endif
38 #if defined HAVE_BIGDFT
39   use BigDFT_API, only: atomic_info,psp_from_data
40 #endif
41  use m_atomdata
42  use pseudo_pwscf ! pwscf module with all data explicit!
43  use funct_pwscf  ! pwscf module for naming xc functionals
44  use m_xmpi
45 
46  use m_time,     only : timab
47  use m_io_tools, only : open_file
48  use m_fstrings, only : basename, lstrip, sjoin, startswith
49  use m_read_upf_pwscf,  only : read_pseudo
50  use m_pawpsp,   only : pawpsp_read_header_xml,pawpsp_read_pawheader
51  use m_pawxmlps, only : rdpawpsxml,rdpawpsxml_header, paw_setup_free,paw_setuploc
52 
53  implicit none
54 
55  private

ABINIT/pspheads_comm [ Functions ]

[ Top ] [ Functions ]

NAME

 pspheads_comm

FUNCTION

 Communicate pspheads to all processors

INPUTS

  npsp=number of pseudopotentials
  test_paw=0 if no PAW, 1 if PAW

SIDE EFFECTS

  pspheads(npsp)=<type pspheader_type>=all the important information from the
   pseudopotential file headers, as well as the psp file names. On one processor at input,
   on all processors at output

PARENTS

      m_ab7_invars_f90

CHILDREN

      timab,xmpi_bcast

SOURCE

565 subroutine pspheads_comm(npsp,pspheads,test_paw)
566 
567 
568 !This section has been created automatically by the script Abilint (TD).
569 !Do not modify the following lines by hand.
570 #undef ABI_FUNC
571 #define ABI_FUNC 'pspheads_comm'
572 !End of the abilint section
573 
574  implicit none
575 
576 !Arguments ------------------------------------
577  integer,intent(in) :: npsp
578  integer,intent(inout) :: test_paw
579  type(pspheader_type),intent(inout) :: pspheads(npsp)
580 
581 !Local variables-------------------------------
582 #if defined HAVE_MPI
583 !scalars
584  integer,parameter :: master=0
585  integer :: ierr,comm
586 !arrays
587  integer,allocatable :: list_int(:)
588  real(dp) :: tsec(2)
589  real(dp),allocatable :: list_dpr(:)
590  character(len=fnlen),allocatable :: list_char(:)
591 #endif
592 
593 !*************************************************************************
594 
595 #if defined HAVE_MPI
596  call timab(48,1,tsec)
597 
598  comm = xmpi_world
599 
600 !Broadcast the characters (file names and titles)
601  ABI_ALLOCATE(list_char,(3*npsp))
602  list_char(1:npsp)=pspheads(1:npsp)%filpsp
603  list_char(npsp+1:2*npsp)=pspheads(1:npsp)%title
604  list_char(2*npsp+1:3*npsp)=pspheads(1:npsp)%md5_checksum
605 
606  call xmpi_bcast(list_char,master,comm,ierr)
607 
608  pspheads(1:npsp)%filpsp=list_char(1:npsp)
609  pspheads(1:npsp)%title=list_char(npsp+1:2*npsp)
610  pspheads(1:npsp)%md5_checksum=list_char(2*npsp+1:3*npsp)(1:md5_slen)
611  ABI_DEALLOCATE(list_char)
612 
613 !Brodcast the integers
614  ABI_ALLOCATE(list_int,(1+13*npsp))
615  list_int(1        :   npsp) = pspheads(1:npsp)%nproj(0)
616  list_int(1+   npsp: 2*npsp) = pspheads(1:npsp)%nproj(1)
617  list_int(1+ 2*npsp: 3*npsp) = pspheads(1:npsp)%nproj(2)
618  list_int(1+ 3*npsp: 4*npsp) = pspheads(1:npsp)%nproj(3)
619  list_int(1+ 4*npsp: 5*npsp) = pspheads(1:npsp)%lmax
620  list_int(1+ 5*npsp: 6*npsp) = pspheads(1:npsp)%xccc
621  list_int(1+ 6*npsp: 7*npsp) = pspheads(1:npsp)%pspxc
622  list_int(1+ 7*npsp: 8*npsp) = pspheads(1:npsp)%pspdat
623  list_int(1+ 8*npsp: 9*npsp) = pspheads(1:npsp)%pspcod
624  list_int(1+ 9*npsp:10*npsp) = pspheads(1:npsp)%pspso
625  list_int(1+10*npsp:11*npsp) = pspheads(1:npsp)%nprojso(1)
626  list_int(1+11*npsp:12*npsp) = pspheads(1:npsp)%nprojso(2)
627  list_int(1+12*npsp:13*npsp) = pspheads(1:npsp)%nprojso(3)
628  list_int(1+13*npsp)         = test_paw
629 
630  call xmpi_bcast(list_int,master,comm,ierr)
631 
632  pspheads(1:npsp)%nproj(0)   = list_int(1        :   npsp)
633  pspheads(1:npsp)%nproj(1)   = list_int(1+   npsp: 2*npsp)
634  pspheads(1:npsp)%nproj(2)   = list_int(1+ 2*npsp: 3*npsp)
635  pspheads(1:npsp)%nproj(3)   = list_int(1+ 3*npsp: 4*npsp)
636  pspheads(1:npsp)%lmax       = list_int(1+ 4*npsp: 5*npsp)
637  pspheads(1:npsp)%xccc       = list_int(1+ 5*npsp: 6*npsp)
638  pspheads(1:npsp)%pspxc      = list_int(1+ 6*npsp: 7*npsp)
639  pspheads(1:npsp)%pspdat     = list_int(1+ 7*npsp: 8*npsp)
640  pspheads(1:npsp)%pspcod     = list_int(1+ 8*npsp: 9*npsp)
641  pspheads(1:npsp)%pspso      = list_int(1+ 9*npsp:10*npsp)
642  pspheads(1:npsp)%nprojso(1) = list_int(1+10*npsp:11*npsp)
643  pspheads(1:npsp)%nprojso(2) = list_int(1+11*npsp:12*npsp)
644  pspheads(1:npsp)%nprojso(3) = list_int(1+12*npsp:13*npsp)
645  test_paw                    = list_int(1+13*npsp)
646  ABI_DEALLOCATE(list_int)
647 
648 !Unbeliveable, this cannot be sent with the others, for woopy
649  ABI_ALLOCATE(list_int,(npsp))
650  list_int(1:npsp) = pspheads(1:npsp)%usewvl
651  call xmpi_bcast(list_int,master,comm,ierr)
652  pspheads(1:npsp)%usewvl     = list_int(1:npsp)
653  ABI_DEALLOCATE(list_int)
654 
655 
656 !Broadcast zionpsp and znuclpsp
657  ABI_ALLOCATE(list_dpr,(7*npsp))
658  list_dpr(1       :  npsp) = pspheads(1:npsp)%zionpsp
659  list_dpr(1+  npsp:2*npsp) = pspheads(1:npsp)%znuclpsp
660  list_dpr(1+2*npsp:3*npsp) = pspheads(1:npsp)%GTHradii(0)
661  list_dpr(1+3*npsp:4*npsp) = pspheads(1:npsp)%GTHradii(1)
662  list_dpr(1+4*npsp:5*npsp) = pspheads(1:npsp)%GTHradii(2)
663  list_dpr(1+5*npsp:6*npsp) = pspheads(1:npsp)%GTHradii(3)
664  list_dpr(1+6*npsp:7*npsp) = pspheads(1:npsp)%GTHradii(4)
665 
666  call xmpi_bcast(list_dpr,master,comm,ierr)
667 
668  pspheads(1:npsp)%zionpsp     = list_dpr(1       :  npsp)
669  pspheads(1:npsp)%znuclpsp    = list_dpr(1+  npsp:2*npsp)
670  pspheads(1:npsp)%GTHradii(0) = list_dpr(1+2*npsp:3*npsp)
671  pspheads(1:npsp)%GTHradii(1) = list_dpr(1+3*npsp:4*npsp)
672  pspheads(1:npsp)%GTHradii(2) = list_dpr(1+4*npsp:5*npsp)
673  pspheads(1:npsp)%GTHradii(3) = list_dpr(1+5*npsp:6*npsp)
674  pspheads(1:npsp)%GTHradii(4) = list_dpr(1+6*npsp:7*npsp)
675  ABI_DEALLOCATE(list_dpr)
676 
677 !Broadcast additional integers for PAW psps (testpaw was spread, previously)
678  if (test_paw==1) then
679    ABI_ALLOCATE(list_int,(6*npsp))
680    list_int(1       :  npsp)=pspheads(1:npsp)%pawheader%basis_size
681    list_int(1+  npsp:2*npsp)=pspheads(1:npsp)%pawheader%l_size
682    list_int(1+2*npsp:3*npsp)=pspheads(1:npsp)%pawheader%lmn_size
683    list_int(1+3*npsp:4*npsp)=pspheads(1:npsp)%pawheader%mesh_size
684    list_int(1+4*npsp:5*npsp)=pspheads(1:npsp)%pawheader%pawver
685    list_int(1+5*npsp:6*npsp)=pspheads(1:npsp)%pawheader%shape_type
686 
687    call xmpi_bcast(list_int,master,comm,ierr)
688 
689    pspheads(1:npsp)%pawheader%basis_size=list_int(1       :  npsp)
690    pspheads(1:npsp)%pawheader%l_size    =list_int(1+  npsp:2*npsp)
691    pspheads(1:npsp)%pawheader%lmn_size  =list_int(1+2*npsp:3*npsp)
692    pspheads(1:npsp)%pawheader%mesh_size =list_int(1+3*npsp:4*npsp)
693    pspheads(1:npsp)%pawheader%pawver    =list_int(1+4*npsp:5*npsp)
694    pspheads(1:npsp)%pawheader%shape_type=list_int(1+5*npsp:6*npsp)
695    ABI_DEALLOCATE(list_int)
696 
697 !  broadcast rpaw values
698    ABI_ALLOCATE(list_dpr,(2*npsp))
699 
700    list_dpr(1       :  npsp) = pspheads(1:npsp)%pawheader%rpaw
701    list_dpr(1+1*npsp:2*npsp) = pspheads(1:npsp)%pawheader%rshp
702 
703    call xmpi_bcast(list_dpr,master,comm,ierr)
704 
705    pspheads(1:npsp)%pawheader%rpaw = list_dpr(1       :  npsp)
706    pspheads(1:npsp)%pawheader%rshp = list_dpr(1+  npsp:2*npsp)
707 
708    ABI_DEALLOCATE(list_dpr)
709  end if
710 
711  call timab(48,2,tsec)
712 
713 #else
714 !Code to use unused dummy arguments
715  if(pspheads(1)%lmax == -10) pspheads(1)%lmax=-10
716  if(test_paw == -1) test_paw = -1
717 #endif
718 
719 end subroutine pspheads_comm

m_pspheads/inpspheads [ Functions ]

[ Top ] [ m_pspheads ] [ Functions ]

NAME

 inpspheads

FUNCTION

 Read the pseudopotential header of each psp file, in order to initialize pspheads(1:npsp).

INPUTS

  npsp=number of pseudopotentials

OUTPUT

  pspheads(npsp)=<type pspheader_type>=all the important information from the
  pseudopotential file headers, as well as the psp file names
  ecut_tmp(3,2,npsp)= possible ecut values as read in psp files

PARENTS

      m_ab7_invars_f90

CHILDREN

      atomic_info,pawpsxml2ab
      psp_from_data,psxml2abheader,upfheader2abi,wrtout

SOURCE

 90 subroutine inpspheads(filnam,npsp,pspheads,ecut_tmp)
 91 
 92 
 93 !This section has been created automatically by the script Abilint (TD).
 94 !Do not modify the following lines by hand.
 95 #undef ABI_FUNC
 96 #define ABI_FUNC 'inpspheads'
 97 !End of the abilint section
 98 
 99  implicit none
100 
101 !Arguments ------------------------------------
102 !scalars
103  integer,intent(in) :: npsp
104 !arrays
105  real(dp),intent(inout) :: ecut_tmp(3,2,10)
106  character(len=fnlen), intent(in) :: filnam(npsp)
107  type(pspheader_type),intent(inout) :: pspheads(npsp) !vz_i
108 
109 !Local variables-------------------------------
110 !In case a xc core correction is to be taken into account,
111 !the n1xccc value will be given by n1xccc_default. Otherwise it is set to 0.
112 !scalars
113  integer,parameter :: n1xccc_default=2501
114  integer :: extension_switch
115  integer :: idum,ii,ilmax,ipsp,lang,lmax,mmax,mpsang,n1xccc,nmesh
116  integer :: pspcod,pspso,test_paw,usexml,unt,useupf !,pspxc
117  real(dp) :: al,e990,e999,fchrg,qchrg,r1,rchrg,rr ! ,rp,rs
118  character(len=3) :: testxc
119  character(len=500) :: message,errmsg
120  character(len=70) :: testxml
121  character(len=80) :: pspline
122 !arrays
123  integer :: nproj(0:3),nprojso(1:3)
124  integer,allocatable :: orb(:)
125  real(dp) :: hdum(3)
126 #if defined HAVE_BIGDFT
127  !new variables for wvl+paw
128  character(len=2) :: symbol
129  integer :: iasctype,nzatom, nelpsp, npspcode_,ixc_
130  real(dp) :: rcov,ehomo
131  real(dp) :: psppar(0:4,0:6)
132  logical :: exists
133 #endif
134 #if defined HAVE_PSML
135  character(len=3) :: atmsymb
136  character(len=30) :: creator
137 #endif
138 
139 !*************************************************************************
140 
141  test_paw=0
142 
143  do ipsp=1,npsp
144 
145 !  Check if the file is written in XML
146    pspheads(ipsp)%filpsp=trim(filnam(ipsp))
147 
148    usexml = 0
149    if (open_file(filnam(ipsp),message,newunit=unt,form="formatted",status="old") /= 0) then
150      MSG_ERROR(message)
151    end if
152 
153    rewind(unit=unt, err=10, iomsg=errmsg)
154    read(unt,*, err=10, iomsg=errmsg) testxml
155    if(testxml(1:5)=='<?xml')then
156      usexml = 1
157      read(unt,*, err=10, iomsg=errmsg) testxml
158      if(testxml(1:4)=='<paw')then
159        test_paw = 1
160      else
161        test_paw = 0
162      end if
163    else
164      usexml = 0
165    end if
166 
167    close(unit=unt, err=10, iomsg=errmsg)
168 
169 !  Check if pseudopotential file is a Q-espresso UPF file
170    useupf = 0
171    if (open_file(filnam(ipsp),message,newunit=unt,form="formatted",status="old") /= 0) then
172      MSG_ERROR(message)
173    end if
174 
175    rewind(unit=unt,err=10,iomsg=errmsg)
176    read(unt,*,err=10,iomsg=errmsg) testxml ! just a string, no relation to xml.
177    if (testxml(1:9)=='<PP_INFO>') then
178      useupf = 1
179    else
180      useupf = 0
181    end if
182    close(unit=unt,err=10,iomsg=errmsg)
183 
184 !  Read the header of the pseudopotential file
185    if (usexml /= 1 .and. useupf /= 1) then
186 !    Open the psp file and read a normal abinit style header
187      if (open_file(filnam(ipsp), message, newunit=unt, form='formatted', status='old') /= 0) then
188        MSG_ERROR(message)
189      end if
190 
191      rewind (unit=unt, err=10, iomsg=errmsg)
192 
193 !    Read the three first lines
194      read(unt, '(a)', err=10, iomsg=errmsg) pspheads(ipsp)%title
195      read(unt,*, err=10, iomsg=errmsg)pspheads(ipsp)%znuclpsp,pspheads(ipsp)%zionpsp,pspheads(ipsp)%pspdat
196      read(unt,*, err=10, iomsg=errmsg)pspheads(ipsp)%pspcod,pspheads(ipsp)%pspxc,pspheads(ipsp)%lmax,idum,mmax
197 
198      pspcod=pspheads(ipsp)%pspcod
199      lmax=pspheads(ipsp)%lmax
200      write(message,'(a,f5.1,a,i4,a,i4)')'  read the values zionpsp=',pspheads(ipsp)%zionpsp,' , pspcod=',pspcod,' , lmax=',lmax
201      call wrtout(std_out,message,'PERS')
202 
203      nproj(0:3)=0 ; nprojso(1:3)=0
204 
205      pspheads(ipsp)%xccc=0
206      pspheads(ipsp)%pspso=0
207 
208    else if (usexml==1 .and. test_paw==0) then
209 #if defined HAVE_PSML
210      write(message,'(a,a)')  &
211 &     '- inpspheads : Reading pseudopotential header in XML form from ', trim(filnam(ipsp))
212      call wrtout(ab_out,message,'COLL')
213      call wrtout(std_out,  message,'COLL')
214 
215 ! could pass pspheads(ipsp) directly and fill all of it in psxml2ab
216      call psxml2abheader( filnam(ipsp), pspheads(ipsp), atmsymb, creator, 1)
217 
218 ! save some stuff locally for this ipsp
219      pspcod = pspheads(ipsp)%pspcod
220      lmax   = pspheads(ipsp)%lmax
221      nproj = pspheads(ipsp)%nproj
222      nprojso = pspheads(ipsp)%nprojso
223 
224 #else
225      write(message, '(a,a)') "XML norm-conserving pseudopotential has been input,", &
226 &     " but abinit is not compiled with libPSML support. Reconfigure and recompile."
227      MSG_ERROR(message)
228 #endif
229 
230    else if(usexml==1.and.test_paw==1)then
231 
232      write(message,'(a,a)')  &
233 &     '- inpspheads : Reading pseudopotential header in XML form from ', trim(filnam(ipsp))
234      call wrtout(ab_out,message,'COLL')
235      call wrtout(std_out,  message,'COLL')
236      call pawpsxml2ab(filnam(ipsp),ecut_tmp(:,:,ipsp), pspheads(ipsp),1)
237      pspcod=17; pspheads(ipsp)%pspcod=pspcod
238 
239    else if (useupf == 1) then
240      pspheads(ipsp)%pspcod = 11
241 
242      pspheads(ipsp)%xccc  = n1xccc_default ! will be set to 0 if no nlcc
243 !    call upfoctheader2abi (filnam(ipsp),  &
244      call upfheader2abi (filnam(ipsp),  &
245 &     pspheads(ipsp)%znuclpsp, &
246 &     pspheads(ipsp)%zionpsp,  &
247 &     pspheads(ipsp)%pspxc,    &
248 &     pspheads(ipsp)%lmax,     &
249 &     pspheads(ipsp)%xccc,     &
250 &     nproj, nprojso)
251 
252      pspcod = pspheads(ipsp)%pspcod
253      lmax   = pspheads(ipsp)%lmax
254 !    FIXME : generalize for SO pseudos
255      pspheads(ipsp)%pspso = 0
256    end if
257 
258 !  DEBUG
259 !  write(std_out,*) pspheads(ipsp)%znuclpsp
260 !  write(std_out,*) pspheads(ipsp)%zionpsp
261 !  write(std_out,*) pspheads(ipsp)%pspcod
262 !  write(std_out,*) pspheads(ipsp)%pspxc
263 !  write(std_out,*) pspheads(ipsp)%lmax
264 !  stop
265 !  ENDDEBUG
266 
267 !  Initialize nproj, nprojso, pspso, as well as xccc, for each type of psp
268    pspheads(ipsp)%GTHradii = zero
269 
270    if(pspcod==1 .or. pspcod==4)then
271 
272 !    Teter format
273      do ilmax=0,lmax
274        read(unt,*, err=10, iomsg=errmsg) lang,e990,e999,nproj(ilmax)
275        read(unt,*, err=10, iomsg=errmsg)
276      end do
277      read(unt,*, err=10, iomsg=errmsg) rchrg,fchrg,qchrg
278      if (fchrg>1.d-15) pspheads(ipsp)%xccc=n1xccc_default
279 
280    else if(pspcod==2)then
281 
282 !    GTH pseudopotentials
283      read(unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%GTHradii(0) !rloc
284      read(unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%GTHradii(1),hdum(1),hdum(2)
285      if(abs(hdum(1))>1.d-9) nproj(0)=1
286      if(abs(hdum(2))>1.d-9) nproj(0)=2
287      read(unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%GTHradii(2),hdum(3)
288      if(abs(hdum(3))>1.d-9) nproj(1)=1
289 
290    else if(pspcod==3)then
291 
292 !    HGH pseudopotentials
293      read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%GTHradii(0) !rloc
294      do ilmax=0,lmax
295        read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%GTHradii(ilmax + 1),hdum(1),hdum(2),hdum(3)
296        if (abs(hdum(1))>1.d-9)nproj(ilmax)=1
297        if (abs(hdum(2))>1.d-9)nproj(ilmax)=2
298        if (abs(hdum(3))>1.d-9)nproj(ilmax)=3
299        if (ilmax>0.and.ilmax<3) then
300          read (unt,*, err=10, iomsg=errmsg) hdum(1),hdum(2),hdum(3)
301          if (abs(hdum(1))>1.d-9)nprojso(ilmax)=1
302          if (abs(hdum(2))>1.d-9)nprojso(ilmax)=2
303          if (abs(hdum(3))>1.d-9)nprojso(ilmax)=3
304          if(nprojso(ilmax)>0)pspheads(ipsp)%pspso=2
305        end if
306        if (ilmax==3) then
307          read (unt,*, err=10, iomsg=errmsg) hdum(1)
308          if (abs(hdum(1))>1.d-9)nprojso(3)=1
309          if(nprojso(3)>0)pspheads(ipsp)%pspso=2
310        end if
311      end do
312 
313    else if(pspcod==5)then
314 
315 !    PHONEY pseudopotentials
316 !    read parameter for Hamman grid
317      pspso=1
318      read (unt,fmt=*,err=50,end=50) r1,al,pspso
319      50 continue
320      do ilmax=0,lmax
321        read (unt,*, err=10, iomsg=errmsg) lang,e990,e999,nproj(ilmax)
322        read (unt,*, err=10, iomsg=errmsg)
323        if (ilmax>0.and.pspso/=1) then
324          read (unt,*, err=10, iomsg=errmsg) lang,e990,e999,nprojso(ilmax)
325          read (unt,*, err=10, iomsg=errmsg)
326          pspheads(ipsp)%pspso=pspso
327 !        Meaning of pspso internally to ABINIT has been changed in v5.4
328 !        So : file must contain pspso 1 , but ABINIT will have pspso 0 .
329          if(pspso==1)pspheads(ipsp)%pspso=0
330        end if
331      end do
332      read (unt,*, err=10, iomsg=errmsg) rchrg,fchrg,qchrg
333      if (fchrg>1.d-15) pspheads(ipsp)%xccc=n1xccc_default
334 
335    else if(pspcod==6)then
336 
337 !    FHI pseudopotentials
338      read (unt, '(a3)') testxc
339 !    Note : prior to version 2.2, this 4th line started with  4--  ,
340 !    and no core-correction was available.
341      if(testxc/='4--')then
342        backspace(unt, err=10, iomsg=errmsg)
343        read (unt,*, err=10, iomsg=errmsg) rchrg,fchrg,qchrg
344      else
345        fchrg=0.0_dp
346      end if
347      if (fchrg>1.d-15) pspheads(ipsp)%xccc=n1xccc_default
348 !    XG020728 : Should take lloc into account ??
349      do ilmax=0,lmax
350        nproj(ilmax)=1
351      end do
352 
353    else if(pspcod==7)then
354 
355 !    PAW pseudopotentials
356      test_paw=1;pspheads(ipsp)%pawheader%pawver=1
357      read (unt,'(a80)', err=10, iomsg=errmsg) pspline
358      pspline=adjustl(pspline)
359      if (pspline(1:3)=="paw".or.pspline(1:3)=="PAW") &
360 &     read(unit=pspline(4:80),fmt=*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%pawver
361      if (pspheads(ipsp)%pawheader%pawver==1) then   ! Compatibility with Abinit v4.2.x
362        read (unit=pspline,fmt=*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%basis_size,&
363 &       pspheads(ipsp)%pawheader%lmn_size
364        ABI_ALLOCATE(orb,(pspheads(ipsp)%pawheader%basis_size))
365        orb(:)=0
366        read (unt,*, err=10, iomsg=errmsg) (orb(ii), ii=1,pspheads(ipsp)%pawheader%basis_size)
367        read (unt,*, err=10, iomsg=errmsg)
368        read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%rpaw
369        pspheads(ipsp)%pawheader%rshp=pspheads(ipsp)%pawheader%rpaw
370        read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%mesh_size
371        read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%shape_type
372        if (pspheads(ipsp)%pawheader%shape_type==3) pspheads(ipsp)%pawheader%shape_type=-1
373      else
374        read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%basis_size,&
375 &       pspheads(ipsp)%pawheader%lmn_size
376        ABI_ALLOCATE(orb,(pspheads(ipsp)%pawheader%basis_size))
377        orb(:)=0
378        read (unt,*, err=10, iomsg=errmsg) (orb(ii), ii=1,pspheads(ipsp)%pawheader%basis_size)
379        pspheads(ipsp)%pawheader%mesh_size=mmax
380        read (unt,*, err=10, iomsg=errmsg) nmesh
381        do ii=1,nmesh
382          read(unt,*, err=10, iomsg=errmsg)
383        end do
384        read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%rpaw
385        pspheads(ipsp)%pawheader%rshp=pspheads(ipsp)%pawheader%rpaw
386        read (unt,'(a80)', err=10, iomsg=errmsg) pspline
387        pspline=adjustl(pspline); write(std_out,*) pspline
388        read(unit=pspline,fmt=*) pspheads(ipsp)%pawheader%shape_type
389        if (pspheads(ipsp)%pawheader%pawver==2.and.&
390 &       pspheads(ipsp)%pawheader%shape_type==3) pspheads(ipsp)%pawheader%shape_type=-1
391        if (pspheads(ipsp)%pawheader%pawver>=3.and.pspheads(ipsp)%pawheader%shape_type==-1) then
392          rr=zero;read(unit=pspline,fmt=*,err=20,end=20) ii,rr
393          20       continue
394          if (rr>=tol8) pspheads(ipsp)%pawheader%rshp=rr
395        end if
396      end if
397      do ilmax=0,lmax
398        do ii=1,pspheads(ipsp)%pawheader%basis_size
399          if(orb(ii)==ilmax) nproj(ilmax)=nproj(ilmax)+1
400        end do
401      end do
402      pspheads(ipsp)%pawheader%l_size=2*maxval(orb)+1
403      pspheads(ipsp)%xccc=1  ! We suppose apriori that cc is used (but n1xccc is not used in PAW)
404      ABI_DEALLOCATE(orb)
405 
406 !    WVL+PAW case, need to define GTHradii
407 #if defined HAVE_BIGDFT
408      if(pspheads(ipsp)%usewvl==1) then
409 !      Obtain the HGH parameters by default from BigDFT
410 
411        call atomic_info(int(pspheads(ipsp)%znuclpsp), int(pspheads(ipsp)%zionpsp), &
412 &       symbol = symbol, ehomo = ehomo, rcov = rcov, nsccode = iasctype)
413 
414 !      I use the XC: Perdew, Burke & Ernzerhof  as default, since
415 !      other XC potentials may not be in the BigDFT table.
416        ixc_=1
417        call psp_from_data(symbol, nzatom, nelpsp, npspcode_, ixc_, psppar, exists)
418        if(.not. exists) then
419          write(message,'(4a)')ch10,&
420 &         "Chemical element not found in BigDFT table",ch10,&
421 &         "Action: upgrade BigDFT table"
422          MSG_BUG(message)
423        end if
424 !
425 !      pspheads(ipsp)%pawheader%rpaw/4.0d0
426        pspheads(ipsp)%GTHradii(0)=psppar(0,0) !rloc
427        pspheads(ipsp)%GTHradii(1)=psppar(1,0) !rrs
428        pspheads(ipsp)%GTHradii(2)=psppar(2,0) !rrp
429 !      pspheads(ipsp)%GTHradii(1) = one / sqrt(abs(two * ehomo))
430 !      write(*,*)pspheads(ipsp)%GTHradii(:)
431      end if
432 #endif
433 
434    else if(pspcod==8)then
435 
436 !    DRH pseudopotentials
437      read(unt,*, err=10, iomsg=errmsg) rchrg,fchrg,qchrg
438      if (fchrg>1.d-15) pspheads(ipsp)%xccc=n1xccc_default
439      read(unt,*, err=10, iomsg=errmsg) nproj(0:lmax)
440      read(unt,*, err=10, iomsg=errmsg) extension_switch
441      if(any(extension_switch == [2, 3])) then
442        pspso=2
443        read(unt,*,err=10,iomsg=errmsg) nprojso(1:lmax)
444      else
445        pspso=0
446      end if
447      pspheads(ipsp)%pspso=pspso
448 
449    else if(pspcod==9)then
450 ! placeholder: nothing to do everything is read above
451 
452    else if(pspcod==10)then
453 
454 !    HGH pseudopotentials, full h/k matrices
455      read (unt,*,err=10,iomsg=errmsg) pspheads(ipsp)%GTHradii(0) !rloc
456      read (unt,*,err=10,iomsg=errmsg) idum
457      if(idum-1/=lmax) then
458        MSG_ERROR("in inpspheads: nnonloc-1 /= lmax")
459      end if
460      do ilmax=0,lmax
461        read (unt,*,err=10,iomsg=errmsg) &
462        pspheads(ipsp)%GTHradii(ilmax + 1),nproj(ilmax),(hdum(idum),idum=1,nproj(ilmax))
463        do idum=2,nproj(ilmax) !skip the rest of h_ij
464          read (unt,*, err=10, iomsg=errmsg)
465        end do
466        if (ilmax==0) cycle
467        nprojso(ilmax)=nproj(ilmax)
468        if(nprojso(ilmax)>0)then
469          pspheads(ipsp)%pspso=2
470          do idum=1,nprojso(ilmax) !skip the rest of k_ij
471            read (unt,*, err=10, iomsg=errmsg)
472          end do
473        end if
474      end do
475 
476    else if (pspcod == 11.or.pspcod == 17) then
477 !    already done above
478    else
479 
480      write(message, '(a,i0,4a)' )&
481 &     'The pseudopotential code (pspcod) read from file is ',pspcod,ch10,&
482 &     'This value is not allowed (should be between 1 and 10). ',ch10,&
483 &     'Action: use a correct pseudopotential file.'
484      MSG_ERROR(message)
485    end if ! pspcod=...
486 
487 !  Store in pspheads
488    if (pspcod /= 17) then
489      pspheads(ipsp)%nproj(0:3)=nproj(0:3)
490      pspheads(ipsp)%nprojso(1:3)=nprojso(1:3)
491    end if
492 !  DEBUG
493 !   write(message,'(a,10i6)') 'nproj = ', pspheads(ipsp)%nproj(:)
494 !   call wrtout(std_out,message,'PERS')
495 !   write(message,'(a,10i6)') 'nprojso = ', pspheads(ipsp)%nprojso(:)
496 !   call wrtout(std_out,message,'PERS')
497 !  ENDDEBUG
498 
499    close(unt)
500 
501    ! Compute md5 checksum
502    pspheads(ipsp)%md5_checksum = md5_sum_from_file(filnam(ipsp))
503  end do ! ipsp=1,npsp
504 
505 !Note that mpsang is the max of 1+lmax, with minimal value 1 (even for local psps, at present)
506 !mpsang=max(maxval(pspheads(1:npsp)%lmax)+1,1) ! Likely troubles with HP compiler
507 !n1xccc=maxval(pspheads(1:npsp)%xccc)
508  mpsang=1
509  n1xccc=pspheads(1)%xccc
510  do ii=1,npsp
511    mpsang=max(pspheads(ii)%lmax+1,mpsang)
512    n1xccc=max(pspheads(ii)%xccc,n1xccc)
513  end do
514 
515  write(message,'(2a,i0,a,i0,a)')ch10,&
516 & ' inpspheads: deduce mpsang = ',mpsang,', n1xccc = ',n1xccc,'.'
517  call wrtout(std_out,message,'PERS')
518 
519 !Test: if one psp is PAW, all must be
520  if (test_paw==1) then
521    do ipsp=1,npsp
522      if (all(pspheads(ipsp)%pspcod /= [7, 17])) then
523        write(message, '(5a)' )&
524 &       'One pseudopotential is PAW (pspcod=7 or 17) !',ch10,&
525 &       'All pseudopotentials must be PAW (this is not the case here) !',ch10,&
526 &       'Action: use only PAW pseudopotential files.'
527        MSG_ERROR(message)
528      end if
529    end do
530  end if
531 
532  return
533 
534  ! Handle IO error
535  10 continue
536  MSG_ERROR(errmsg)
537 
538 end subroutine inpspheads

m_pspheads/pawpsxml2ab [ Functions ]

[ Top ] [ m_pspheads ] [ Functions ]

NAME

 pawpsxml2ab

FUNCTION

  From a XML format pseudopotential file which has already been read in,
  convert to abinit internal datastructures.

INPUTS

  ecut_tmp(3,2)= possible ecut values as read in psp files
  filenam= input file name (atomicdata XML)
  option= 1 if header only is read; 0 if the whole data are read

OUTPUT

 pspheads data structure is filled

PARENTS

      inpspheads

CHILDREN

      pawpsp_read_header_xml,pawpsp_read_pawheader

SOURCE

746 subroutine pawpsxml2ab( filnam,ecut_tmp, pspheads,option)
747 
748 
749 !This section has been created automatically by the script Abilint (TD).
750 !Do not modify the following lines by hand.
751 #undef ABI_FUNC
752 #define ABI_FUNC 'pawpsxml2ab'
753 !End of the abilint section
754 
755  implicit none
756 
757 !Arguments ------------------------------------
758 !scalars
759  integer, intent(in) :: option
760  character(len=fnlen), intent(in) :: filnam
761  type(pspheader_type),intent(inout) :: pspheads !vz_i
762 !arrays
763  real(dp),intent(inout) :: ecut_tmp(3,2)
764 
765 !Local variables-------------------------------
766 !scalars
767  integer :: ii,il,lloc,lmax,pspcod,pspxc,unt
768  real(dp) :: r2well,zionpsp,znuclpsp
769 ! character(len=100) :: xclibxc
770 ! character(len=500) :: message
771 !arrays
772 
773 ! *********************************************************************
774  
775  if (option==1) then
776    call rdpawpsxml_header(ecut_tmp,filnam,paw_setuploc)
777    paw_setuploc%idgrid= paw_setuploc%radial_grid(1)%id
778  else
779    call rdpawpsxml(filnam,paw_setuploc)
780  end if
781 
782  call pawpsp_read_header_xml(lloc,lmax,pspcod,&
783 &   pspxc,paw_setuploc,r2well,zionpsp,znuclpsp)
784 
785  pspheads%lmax=lmax
786  pspheads%pspxc=pspxc
787  pspheads%zionpsp=zionpsp
788  pspheads%znuclpsp=znuclpsp
789 
790  call pawpsp_read_pawheader(pspheads%pawheader%basis_size,&
791 &   pspheads%lmax,pspheads%pawheader%lmn_size,&
792 &   pspheads%pawheader%l_size,pspheads%pawheader%mesh_size,&
793 &   pspheads%pawheader%pawver,paw_setuploc,pspheads%pawheader%rpaw,&
794 &   pspheads%pawheader%rshp,pspheads%pawheader%shape_type)
795 
796  pspheads%nproj=0
797  do il=0,pspheads%lmax
798    do ii=1,pspheads%pawheader%basis_size
799      if(paw_setuploc%valence_states%state(ii)%ll==il) pspheads%nproj(il)=pspheads%nproj(il)+1
800    end do
801  end do
802  pspheads%nprojso=0
803  pspheads%pspdat=27061961
804  pspheads%pspso=1
805  pspheads%xccc=1
806  pspheads%title=paw_setuploc%atom%symbol
807 
808  if (option==1) then
809    call paw_setup_free(paw_setuploc)
810  end if
811 
812 end subroutine pawpsxml2ab

m_pspheads/upfheader2abi [ Functions ]

[ Top ] [ m_pspheads ] [ Functions ]

NAME

 upfheader2abi

FUNCTION

  This routine wraps a call to a PWSCF module, which reads in
  a UPF (PWSCF / Espresso) format pseudopotential, then transfers
  data for the HEADER of abinit psps only!

INPUTS

  filpsp = name of file with UPF data

OUTPUT

  pspxc = index of xc functional for this pseudo
  lmax_ = maximal angular momentum
  znucl = charge of species nucleus
  zion = valence charge
  n1xccc = default number of points. Set to 0 if no nlcc is present
  nproj_l= number of projectors for each channel
  nprojso_l= number of projectors for each channel for SO correction projectors

PARENTS

      inpspheads

CHILDREN

      set_dft_from_indices,set_dft_from_name

SOURCE

844 subroutine upfheader2abi (filpsp, znucl, zion, pspxc, lmax_, n1xccc, nproj_l, nprojso_l)
845 
846 
847 !This section has been created automatically by the script Abilint (TD).
848 !Do not modify the following lines by hand.
849 #undef ABI_FUNC
850 #define ABI_FUNC 'upfheader2abi'
851 !End of the abilint section
852 
853   implicit none
854 
855 !Arguments -------------------------------
856   character(len=fnlen), intent(in) :: filpsp
857   integer, intent(inout) :: n1xccc
858   integer, intent(out) :: pspxc, lmax_
859   real(dp), intent(out) :: znucl, zion
860   !arrays
861   integer, intent(out) :: nproj_l(0:3)
862   integer, intent(out) :: nprojso_l(1:3)
863 
864 !Local variables -------------------------
865   integer :: iproj, ll, iunit
866   character(len=500) :: msg
867   type(atomdata_t) :: atom
868 
869 ! *********************************************************************
870 
871 !call pwscf routine for reading in UPF
872  if (open_file(filpsp, msg, newunit=iunit, status='old',form='formatted') /= 0) then
873    MSG_ERROR(msg)
874  end if
875 
876 !read in psp data to static data in pseudo module, for ipsx == 1
877  call read_pseudo(1,iunit)
878  close (iunit)
879 
880 !copy over to abinit internal arrays and vars
881  call upfxc2abi(dft(1), pspxc)
882  lmax_ = lmax(1)
883  call atomdata_from_symbol(atom,psd(1))
884  znucl = atom%znucl
885  zion = zp(1)
886 
887  nproj_l = 0
888  do iproj = 1, nbeta(1)
889    ll = lll(iproj,1)
890    nproj_l(ll) = nproj_l(ll) + 1
891  end do
892 
893  nprojso_l = 0 !FIXME deal with so
894 !do iproj = 1, nbeta(1)
895 !nprojso_l(ll+1) = nprojso_l(ll+1) + 1
896 !end do
897 
898  if (.not. nlcc(1)) n1xccc = 0
899 
900 end subroutine upfheader2abi

m_pspheads/upfxc2abi [ Functions ]

[ Top ] [ m_pspheads ] [ Functions ]

NAME

 upfxc2abi

FUNCTION

  This routine wraps a call to an OCTOPUS module, which reformats
  a UPF (PWSCF / Espresso) string describing XC functionals,
  and returns the abinit internal code pspxc

INPUTS

  dft = string with x/c functionals from PWSCF format

OUTPUT

  pspxc = index of xc functional for this pseudo

NOTES

   FIXME: extend to more functionals with libxc
   Could be included in separate module, eg read_upf_pwscf or funct_pwscf
   Left without defs_basis or calls to abinit routines ON PURPOSE

PARENTS

      upf2abinit,upfheader2abi

CHILDREN

      set_dft_from_indices,set_dft_from_name

SOURCE

 931 subroutine upfxc2abi(dft, pspxc)
 932 
 933 
 934 !This section has been created automatically by the script Abilint (TD).
 935 !Do not modify the following lines by hand.
 936 #undef ABI_FUNC
 937 #define ABI_FUNC 'upfxc2abi'
 938 !End of the abilint section
 939 
 940   implicit none
 941 
 942 !Arguments -------------------------------
 943   character(len=20), intent(in) :: dft
 944   integer, intent(out) :: pspxc
 945 
 946 !Local variables -------------------------
 947   integer :: iexch,icorr,igcx,igcc
 948   integer :: totalindex, offset
 949 
 950 ! *********************************************************************
 951 !extract from char*20 :: dft(:)
 952 !###  The following has been copied from pwscf src/Modules/upf_to_internal.f90:
 953 !workaround for rrkj format - it contains the indices, not the name
 954  if ( dft(1:6)=='INDEX:') then
 955    read( dft(7:10), '(4i1)') iexch,icorr,igcx,igcc
 956    call set_dft_from_indices(iexch,icorr,igcx,igcc)
 957  else
 958    call set_dft_from_name( dft )
 959    iexch = get_iexch()
 960    icorr = get_icorr()
 961    igcx = get_igcx()
 962    igcc = get_igcc()
 963  end if
 964 !reset dft string to avoid stray spaces
 965  call set_dft_from_indices(iexch,icorr,igcx,igcc)
 966  write(std_out,'(a)') ' upf2abinit: XC string from pseudopotential is :'
 967  write(std_out,'(3a)') '>', dft, '<'
 968 
 969  offset = 100
 970  totalindex = offset*offset*offset*iexch + offset*offset*icorr + offset*igcx + igcc
 971  select case (totalindex)
 972  case (00000000)  !(" NOX  NOC NOGX NOGC") ! no xc
 973    pspxc = 0
 974  case (01010000)  !(" SLA   PZ NOGX NOGC") ! slater exchange + Perdew Zunger
 975    pspxc = 2
 976  case (01050000)  !(" SLA  WIG NOGX NOGC") ! slater exchange + Wigner corr
 977    pspxc = 4
 978  case (01060000)  !(" SLA   HL NOGX NOGC") ! Hedin + Lundqvist
 979    pspxc = 5
 980  case (02000000)  !(" SL1  NOC NOGX NOGC") ! full slater exchange
 981    pspxc = 6
 982  case (01040000)  !(" SLA   PW NOGX NOGC") ! slater exchange + Perdew Wang
 983    pspxc = 7
 984  case (01000000)  !(" SLA  NOC NOGX NOGC") ! Perdew Wang + no corr
 985    pspxc = 8
 986  case (01040304)  !(" SLA   PW  PBX  PBC") ! LDA + PBE GGA
 987    pspxc = 11 ! PBE
 988  case (01000300)  !(" SLA  NOC  PBX NOGC") ! exchange part of PBE GGA
 989    pspxc = 12
 990  case (01040404)  !(" SLA   PW  RPB  PBC") ! rev PBE
 991    pspxc = 14
 992  case (00000505)  !(" NOX  NOC HTCH HTCH") ! HTCH 120
 993    pspxc = 17
 994  case (01030103)  !(" SLA  LYP  B88 BLYP") ! BLYP
 995    pspxc = -106131
 996  case (01040101)  !(" SLA   PW  B88  P86") ! BP86
 997    pspxc = -106132
 998  case (00030603)  !(" NOX  LYP OPTX BLYP") ! OLYP
 999    pspxc = -110131
1000 !    FIXME: important cases left to be patched with libxc:
1001 !    vosko wilkins nusair
1002 !    ortiz ballone
1003 !    pbe0
1004 !    Gunnarson-Lunqvist
1005 !    make general approach: check gradient parts first, then lda.
1006 !    event. check if they are consistent.
1007  case default
1008    MSG_ERROR('upf2abinit: XC functional not recognized')
1009  end select
1010 
1011 end subroutine upfxc2abi