TABLE OF CONTENTS
- ABINIT/m_pspheads
- ABINIT/pspheads_comm
- m_pspheads/inpspheads
- m_pspheads/pawpsxml2ab
- m_pspheads/updft_to_ixc
- m_pspheads/upf1_to_psphead
- m_pspheads/upf2_jl2srso
- m_pspheads/upf2_to_psphead
- m_pspheads/upfxc2abi
ABINIT/m_pspheads [ 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-2024 ABINIT group (DCA, XG, GMR, FrD, AF, MT, FJ, MJV, MG, DRH) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
17 #if defined HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "abi_common.h" 22 23 MODULE m_pspheads 24 25 use defs_basis 26 use m_abicore 27 use m_errors 28 use m_xmpi 29 use m_atomdata 30 use m_hash_md5 31 use m_psxml2ab 32 #if defined HAVE_LIBPSML 33 use m_psml 34 #endif 35 #if defined HAVE_BIGDFT 36 use BigDFT_API, only: atomic_info, psp_from_data 37 #endif 38 39 use defs_datatypes, only : pspheader_type 40 use m_time, only : timab 41 use m_io_tools, only : open_file 42 use m_numeric_tools,only : simpson 43 use m_fstrings, only : basename, lstrip, sjoin, startswith, atoi, itoa, ftoa, toupper, next_token 44 use m_pawpsp, only : pawpsp_read_header_xml,pawpsp_read_pawheader 45 use m_pawxmlps, only : rdpawpsxml,rdpawpsxml_header, paw_setup_free,paw_setuploc 46 use pseudo_types, only : pseudo_upf, deallocate_pseudo_upf !, pseudo_config 47 use read_upf_new_module, only : read_upf_new 48 49 implicit none 50 51 private
ABINIT/pspheads_comm [ 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
SOURCE
550 subroutine pspheads_comm(npsp,pspheads,test_paw) 551 552 !Arguments ------------------------------------ 553 integer,intent(in) :: npsp 554 integer,intent(inout) :: test_paw 555 type(pspheader_type),intent(inout) :: pspheads(npsp) 556 557 !Local variables------------------------------- 558 #if defined HAVE_MPI 559 !scalars 560 integer,parameter :: master=0 561 integer :: ierr,comm 562 !arrays 563 integer,allocatable :: list_int(:) 564 real(dp) :: tsec(2) 565 real(dp),allocatable :: list_dpr(:) 566 character(len=fnlen),allocatable :: list_char(:) 567 #endif 568 569 !************************************************************************* 570 571 #if defined HAVE_MPI 572 call timab(48,1,tsec) 573 574 comm = xmpi_world 575 576 ! Broadcast the characters (file names and titles) 577 ABI_MALLOC(list_char,(3*npsp)) 578 list_char(1:npsp)=pspheads(1:npsp)%filpsp 579 list_char(npsp+1:2*npsp)=pspheads(1:npsp)%title 580 list_char(2*npsp+1:3*npsp)=pspheads(1:npsp)%md5_checksum 581 582 call xmpi_bcast(list_char,master,comm,ierr) 583 584 pspheads(1:npsp)%filpsp=list_char(1:npsp) 585 pspheads(1:npsp)%title=list_char(npsp+1:2*npsp) 586 pspheads(1:npsp)%md5_checksum=list_char(2*npsp+1:3*npsp)(1:md5_slen) 587 ABI_FREE(list_char) 588 589 ! Brodcast the integers 590 ABI_MALLOC(list_int,(1+13*npsp)) 591 list_int(1 : npsp) = pspheads(1:npsp)%nproj(0) 592 list_int(1+ npsp: 2*npsp) = pspheads(1:npsp)%nproj(1) 593 list_int(1+ 2*npsp: 3*npsp) = pspheads(1:npsp)%nproj(2) 594 list_int(1+ 3*npsp: 4*npsp) = pspheads(1:npsp)%nproj(3) 595 list_int(1+ 4*npsp: 5*npsp) = pspheads(1:npsp)%lmax 596 list_int(1+ 5*npsp: 6*npsp) = pspheads(1:npsp)%xccc 597 list_int(1+ 6*npsp: 7*npsp) = pspheads(1:npsp)%pspxc 598 list_int(1+ 7*npsp: 8*npsp) = pspheads(1:npsp)%pspdat 599 list_int(1+ 8*npsp: 9*npsp) = pspheads(1:npsp)%pspcod 600 list_int(1+ 9*npsp:10*npsp) = pspheads(1:npsp)%pspso 601 list_int(1+10*npsp:11*npsp) = pspheads(1:npsp)%nprojso(1) 602 list_int(1+11*npsp:12*npsp) = pspheads(1:npsp)%nprojso(2) 603 list_int(1+12*npsp:13*npsp) = pspheads(1:npsp)%nprojso(3) 604 list_int(1+13*npsp) = test_paw 605 606 call xmpi_bcast(list_int,master,comm,ierr) 607 608 pspheads(1:npsp)%nproj(0) = list_int(1 : npsp) 609 pspheads(1:npsp)%nproj(1) = list_int(1+ npsp: 2*npsp) 610 pspheads(1:npsp)%nproj(2) = list_int(1+ 2*npsp: 3*npsp) 611 pspheads(1:npsp)%nproj(3) = list_int(1+ 3*npsp: 4*npsp) 612 pspheads(1:npsp)%lmax = list_int(1+ 4*npsp: 5*npsp) 613 pspheads(1:npsp)%xccc = list_int(1+ 5*npsp: 6*npsp) 614 pspheads(1:npsp)%pspxc = list_int(1+ 6*npsp: 7*npsp) 615 pspheads(1:npsp)%pspdat = list_int(1+ 7*npsp: 8*npsp) 616 pspheads(1:npsp)%pspcod = list_int(1+ 8*npsp: 9*npsp) 617 pspheads(1:npsp)%pspso = list_int(1+ 9*npsp:10*npsp) 618 pspheads(1:npsp)%nprojso(1) = list_int(1+10*npsp:11*npsp) 619 pspheads(1:npsp)%nprojso(2) = list_int(1+11*npsp:12*npsp) 620 pspheads(1:npsp)%nprojso(3) = list_int(1+12*npsp:13*npsp) 621 test_paw = list_int(1+13*npsp) 622 ABI_FREE(list_int) 623 624 ! Unbeliveable, this cannot be sent with the others, for woopy 625 ABI_MALLOC(list_int,(npsp)) 626 list_int(1:npsp) = pspheads(1:npsp)%usewvl 627 call xmpi_bcast(list_int,master,comm,ierr) 628 pspheads(1:npsp)%usewvl = list_int(1:npsp) 629 ABI_FREE(list_int) 630 631 ! Broadcast zionpsp and znuclpsp 632 ABI_MALLOC(list_dpr,(7*npsp)) 633 list_dpr(1 : npsp) = pspheads(1:npsp)%zionpsp 634 list_dpr(1+ npsp:2*npsp) = pspheads(1:npsp)%znuclpsp 635 list_dpr(1+2*npsp:3*npsp) = pspheads(1:npsp)%GTHradii(0) 636 list_dpr(1+3*npsp:4*npsp) = pspheads(1:npsp)%GTHradii(1) 637 list_dpr(1+4*npsp:5*npsp) = pspheads(1:npsp)%GTHradii(2) 638 list_dpr(1+5*npsp:6*npsp) = pspheads(1:npsp)%GTHradii(3) 639 list_dpr(1+6*npsp:7*npsp) = pspheads(1:npsp)%GTHradii(4) 640 641 call xmpi_bcast(list_dpr,master,comm,ierr) 642 643 pspheads(1:npsp)%zionpsp = list_dpr(1 : npsp) 644 pspheads(1:npsp)%znuclpsp = list_dpr(1+ npsp:2*npsp) 645 pspheads(1:npsp)%GTHradii(0) = list_dpr(1+2*npsp:3*npsp) 646 pspheads(1:npsp)%GTHradii(1) = list_dpr(1+3*npsp:4*npsp) 647 pspheads(1:npsp)%GTHradii(2) = list_dpr(1+4*npsp:5*npsp) 648 pspheads(1:npsp)%GTHradii(3) = list_dpr(1+5*npsp:6*npsp) 649 pspheads(1:npsp)%GTHradii(4) = list_dpr(1+6*npsp:7*npsp) 650 ABI_FREE(list_dpr) 651 652 ! Broadcast additional integers for PAW psps (testpaw was sent, previously) 653 if (test_paw==1) then 654 ABI_MALLOC(list_int,(6*npsp)) 655 list_int(1 : npsp)=pspheads(1:npsp)%pawheader%basis_size 656 list_int(1+ npsp:2*npsp)=pspheads(1:npsp)%pawheader%l_size 657 list_int(1+2*npsp:3*npsp)=pspheads(1:npsp)%pawheader%lmn_size 658 list_int(1+3*npsp:4*npsp)=pspheads(1:npsp)%pawheader%mesh_size 659 list_int(1+4*npsp:5*npsp)=pspheads(1:npsp)%pawheader%pawver 660 list_int(1+5*npsp:6*npsp)=pspheads(1:npsp)%pawheader%shape_type 661 662 call xmpi_bcast(list_int,master,comm,ierr) 663 664 pspheads(1:npsp)%pawheader%basis_size=list_int(1 : npsp) 665 pspheads(1:npsp)%pawheader%l_size =list_int(1+ npsp:2*npsp) 666 pspheads(1:npsp)%pawheader%lmn_size =list_int(1+2*npsp:3*npsp) 667 pspheads(1:npsp)%pawheader%mesh_size =list_int(1+3*npsp:4*npsp) 668 pspheads(1:npsp)%pawheader%pawver =list_int(1+4*npsp:5*npsp) 669 pspheads(1:npsp)%pawheader%shape_type=list_int(1+5*npsp:6*npsp) 670 ABI_FREE(list_int) 671 672 ! broadcast rpaw values 673 ABI_MALLOC(list_dpr,(2*npsp)) 674 675 list_dpr(1 : npsp) = pspheads(1:npsp)%pawheader%rpaw 676 list_dpr(1+1*npsp:2*npsp) = pspheads(1:npsp)%pawheader%rshp 677 678 call xmpi_bcast(list_dpr,master,comm,ierr) 679 680 pspheads(1:npsp)%pawheader%rpaw = list_dpr(1 : npsp) 681 pspheads(1:npsp)%pawheader%rshp = list_dpr(1+ npsp:2*npsp) 682 683 ABI_FREE(list_dpr) 684 end if 685 686 call timab(48,2,tsec) 687 688 #else 689 ! Code to use unused dummy arguments 690 if(pspheads(1)%lmax == -10) pspheads(1)%lmax=-10 691 if(test_paw == -1) test_paw = -1 692 #endif 693 694 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
SOURCE
82 subroutine inpspheads(filnam, npsp, pspheads, ecut_tmp) 83 84 !Arguments ------------------------------------ 85 !scalars 86 integer,intent(in) :: npsp 87 !arrays 88 real(dp),intent(inout) :: ecut_tmp(3,2,10) 89 character(len=fnlen), intent(in) :: filnam(npsp) 90 type(pspheader_type),intent(inout) :: pspheads(npsp) !vz_i 91 92 !Local variables------------------------------- 93 !In case a xc core correction is to be taken into account, 94 !the n1xccc value will be given by n1xccc_default. Otherwise it is set to 0. 95 !scalars 96 integer,parameter :: n1xccc_default=2501 97 integer :: extension_switch 98 integer :: idum,ii,ilmax,ipsp,lang,lmax,mmax,mpsang,n1xccc,nmesh 99 integer :: pspcod,pspso,test_paw,usexml,unt,useupf 100 real(dp) :: al,e990,e999,fchrg,qchrg,r1,rchrg,rr 101 character(len=3) :: testxc 102 character(len=500) :: msg,errmsg 103 character(len=70) :: testxml 104 character(len=80) :: pspline 105 !arrays 106 integer :: nproj(0:3),nprojso(1:3) 107 integer,allocatable :: orb(:) 108 real(dp) :: hdum(3) 109 #if defined HAVE_BIGDFT 110 !new variables for wvl+paw 111 character(len=2) :: symbol 112 integer :: iasctype,nzatom, nelpsp, npspcode_,ixc_ 113 real(dp) :: rcov,ehomo 114 real(dp) :: psppar(0:4,0:6) 115 logical :: exists 116 #endif 117 #if defined HAVE_LIBPSML 118 character(len=3) :: atmsymb 119 character(len=30) :: creator 120 #endif 121 122 !************************************************************************* 123 124 test_paw=0 125 126 do ipsp=1,npsp 127 128 pspheads(ipsp)%filpsp=trim(filnam(ipsp)) 129 130 ! Check if the file is written in XML 131 usexml = 0 132 if (open_file(filnam(ipsp), msg, newunit=unt, form="formatted", status="old") /= 0) then 133 ABI_ERROR(msg) 134 end if 135 136 rewind(unit=unt, err=10, iomsg=errmsg) 137 read(unt, "(a)", err=10, iomsg=errmsg) testxml 138 139 if(testxml(1:5)=='<?xml')then 140 usexml = 1 141 read(unt,*, err=10, iomsg=errmsg) testxml 142 if(testxml(1:4)=='<paw')then 143 test_paw = 1 144 else 145 test_paw = 0 146 end if 147 else 148 usexml = 0 149 end if 150 151 ! Check if pseudopotential file is a QE UPF2 file 152 ! "<UPF version="2.0.1"> 153 useupf = 0 154 if (testxml(1:4) == '<UPF') then 155 ii = index(testxml, '"') 156 if (ii /= 0) then 157 useupf = atoi(testxml(ii+1:ii+1)) 158 else 159 ABI_ERROR(sjoin("Cannot find version attribute in UPF2 file:", filnam(ipsp))) 160 end if 161 end if 162 163 close(unit=unt, err=10, iomsg=errmsg) 164 165 ! Check if pseudopotential file is a QE UPF1 file 166 if (useupf == 0) then 167 if (open_file(filnam(ipsp), msg, newunit=unt, form="formatted", status="old") /= 0) then 168 ABI_ERROR(msg) 169 end if 170 171 rewind(unit=unt, err=10, iomsg=errmsg) 172 read(unt,*,err=10,iomsg=errmsg) testxml ! just a string, no relation to xml. 173 if (testxml(1:9)=='<PP_INFO>') then 174 useupf = 1 175 else 176 useupf = 0 177 end if 178 close(unit=unt,err=10,iomsg=errmsg) 179 end if 180 181 ! Read the header of the pseudopotential file 182 if (usexml /= 1 .and. useupf == 0) then 183 ! Open the psp file and read a normal abinit style header 184 if (open_file(filnam(ipsp), msg, newunit=unt, form='formatted', status='old') /= 0) then 185 ABI_ERROR(msg) 186 end if 187 rewind (unit=unt, err=10, iomsg=errmsg) 188 189 ! Read the three first lines 190 read(unt, '(a)', err=10, iomsg=errmsg) pspheads(ipsp)%title 191 read(unt,*, err=10, iomsg=errmsg)pspheads(ipsp)%znuclpsp,pspheads(ipsp)%zionpsp,pspheads(ipsp)%pspdat 192 read(unt,*, err=10, iomsg=errmsg)pspheads(ipsp)%pspcod,pspheads(ipsp)%pspxc,pspheads(ipsp)%lmax,idum,mmax 193 194 pspcod=pspheads(ipsp)%pspcod 195 lmax=pspheads(ipsp)%lmax 196 write(msg,'(a,f5.1,a,i4,a,i4)')' read the values zionpsp=',pspheads(ipsp)%zionpsp,' , pspcod=',pspcod,' , lmax=',lmax 197 call wrtout(std_out,msg,'PERS') 198 199 nproj(0:3)=0 ; nprojso(1:3)=0 200 201 pspheads(ipsp)%xccc=0 202 pspheads(ipsp)%pspso=0 203 204 else if (usexml==1 .and. test_paw==0) then 205 #if defined HAVE_LIBPSML 206 write(msg,'(4a)') & 207 '- inpspheads : Reading pseudopotential header in XML form from ',ch10,& 208 & '- ',trim(filnam(ipsp)) 209 call wrtout([std_out, ab_out], msg) 210 211 ! could pass pspheads(ipsp) directly and fill all of it in psxml2ab 212 call psxml2abheader( filnam(ipsp), pspheads(ipsp), atmsymb, creator, 1) 213 214 ! save some stuff locally for this ipsp 215 pspcod = pspheads(ipsp)%pspcod 216 lmax = pspheads(ipsp)%lmax 217 nproj = pspheads(ipsp)%nproj 218 nprojso = pspheads(ipsp)%nprojso 219 220 #else 221 write(msg, '(2a)') "XML norm-conserving pseudopotential has been input,", & 222 " but abinit is not compiled with libPSML support. Reconfigure and recompile." 223 ABI_ERROR(msg) 224 #endif 225 226 else if(usexml==1.and.test_paw==1)then 227 228 write(msg,'(4a)') & 229 '- inpspheads : Reading pseudopotential header in XML form from ',ch10,& 230 '- ',trim(filnam(ipsp)) 231 call wrtout([std_out, ab_out], msg) 232 233 call pawpsxml2ab(filnam(ipsp),ecut_tmp(:,:,ipsp), pspheads(ipsp),1) 234 pspcod=17; pspheads(ipsp)%pspcod=pspcod 235 236 else if (useupf > 0) then 237 pspheads(ipsp)%xccc = n1xccc_default ! will be set to 0 if no nlcc 238 239 if (useupf == 1) then 240 pspheads(ipsp)%pspcod = 11 241 call upf1_to_psphead(filnam(ipsp), pspheads(ipsp)%znuclpsp, pspheads(ipsp)%zionpsp, pspheads(ipsp)%pspxc, & 242 pspheads(ipsp)%lmax, pspheads(ipsp)%xccc, nproj, nprojso) 243 244 ! FIXME: generalize for SO pseudos 245 pspheads(ipsp)%pspso = 0 246 247 else 248 ! UPF2 format 249 pspheads(ipsp)%pspcod = 12 250 call upf2_to_psphead(filnam(ipsp), pspheads(ipsp)%znuclpsp, pspheads(ipsp)%zionpsp, pspheads(ipsp)%pspxc, & 251 pspheads(ipsp)%lmax, pspheads(ipsp)%xccc, nproj, nprojso) 252 253 pspheads(ipsp)%pspso = merge(2, 0, any(nprojso > 0)) 254 end if 255 256 pspcod = pspheads(ipsp)%pspcod 257 lmax = pspheads(ipsp)%lmax 258 end if 259 260 !write(std_out,*) "pspheads(ipsp)%znuclpsp", pspheads(ipsp)%znuclpsp 261 !write(std_out,*) "pspheads(ipsp)%zionpsp", pspheads(ipsp)%zionpsp 262 !write(std_out,*) "pspheads(ipsp)%pspcod", pspheads(ipsp)%pspcod 263 !write(std_out,*) "pspheads(ipsp)%pspxc", pspheads(ipsp)%pspxc 264 !write(std_out,*) "pspheads(ipsp)%lmax", pspheads(ipsp)%lmax 265 266 ! Initialize nproj, nprojso, pspso, as well as xccc, for each type of psp 267 pspheads(ipsp)%GTHradii = zero 268 269 if(pspcod==1 .or. pspcod==4)then 270 271 ! Teter format 272 do ilmax=0,lmax 273 read(unt,*, err=10, iomsg=errmsg) lang,e990,e999,nproj(ilmax) 274 read(unt,*, err=10, iomsg=errmsg) 275 end do 276 read(unt,*, err=10, iomsg=errmsg) rchrg,fchrg,qchrg 277 if (fchrg>1.d-15) pspheads(ipsp)%xccc=n1xccc_default 278 279 else if(pspcod==2)then 280 281 ! GTH pseudopotentials 282 read(unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%GTHradii(0) !rloc 283 read(unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%GTHradii(1),hdum(1),hdum(2) 284 if(abs(hdum(1))>1.d-9) nproj(0)=1 285 if(abs(hdum(2))>1.d-9) nproj(0)=2 286 read(unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%GTHradii(2),hdum(3) 287 if(abs(hdum(3))>1.d-9) nproj(1)=1 288 289 else if(pspcod==3)then 290 291 ! HGH pseudopotentials 292 read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%GTHradii(0) !rloc 293 do ilmax=0,lmax 294 read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%GTHradii(ilmax + 1),hdum(1),hdum(2),hdum(3) 295 if (abs(hdum(1))>1.d-9)nproj(ilmax)=1 296 if (abs(hdum(2))>1.d-9)nproj(ilmax)=2 297 if (abs(hdum(3))>1.d-9)nproj(ilmax)=3 298 if (ilmax>0.and.ilmax<3) then 299 read (unt,*, err=10, iomsg=errmsg) hdum(1),hdum(2),hdum(3) 300 if (abs(hdum(1))>1.d-9)nprojso(ilmax)=1 301 if (abs(hdum(2))>1.d-9)nprojso(ilmax)=2 302 if (abs(hdum(3))>1.d-9)nprojso(ilmax)=3 303 if(nprojso(ilmax)>0)pspheads(ipsp)%pspso=2 304 end if 305 if (ilmax==3) then 306 read (unt,*, err=10, iomsg=errmsg) hdum(1) 307 if (abs(hdum(1))>1.d-9)nprojso(3)=1 308 if(nprojso(3)>0)pspheads(ipsp)%pspso=2 309 end if 310 end do 311 312 else if(pspcod==5)then 313 314 ! PHONEY pseudopotentials 315 ! read parameter for Hamman grid 316 pspso=1 317 read (unt,fmt=*,err=50,end=50) r1,al,pspso 318 50 continue 319 do ilmax=0,lmax 320 read (unt,*, err=10, iomsg=errmsg) lang,e990,e999,nproj(ilmax) 321 read (unt,*, err=10, iomsg=errmsg) 322 if (ilmax>0.and.pspso/=1) then 323 read (unt,*, err=10, iomsg=errmsg) lang,e990,e999,nprojso(ilmax) 324 read (unt,*, err=10, iomsg=errmsg) 325 pspheads(ipsp)%pspso=pspso 326 ! Meaning of pspso internally to ABINIT has been changed in v5.4 327 ! So file must contain pspso 1, but ABINIT will have pspso 0. 328 if(pspso==1)pspheads(ipsp)%pspso=0 329 end if 330 end do 331 read (unt,*, err=10, iomsg=errmsg) rchrg,fchrg,qchrg 332 if (fchrg>1.d-15) pspheads(ipsp)%xccc=n1xccc_default 333 334 else if(pspcod==6)then 335 336 ! FHI pseudopotentials 337 read (unt, '(a3)') testxc 338 ! Note: prior to version 2.2, this 4th line started with 4-- , 339 ! and no core-correction was available. 340 if(testxc/='4--')then 341 backspace(unt, err=10, iomsg=errmsg) 342 read (unt,*, err=10, iomsg=errmsg) rchrg,fchrg,qchrg 343 else 344 fchrg=0.0_dp 345 end if 346 if (fchrg>1.d-15) pspheads(ipsp)%xccc=n1xccc_default 347 ! XG020728 : Should take lloc into account ?? 348 do ilmax=0,lmax 349 nproj(ilmax)=1 350 end do 351 352 else if(pspcod==7)then 353 354 ! PAW pseudopotentials 355 test_paw=1;pspheads(ipsp)%pawheader%pawver=1 356 read (unt,'(a80)', err=10, iomsg=errmsg) pspline 357 pspline=adjustl(pspline) 358 if (pspline(1:3)=="paw".or.pspline(1:3)=="PAW") & 359 & read(unit=pspline(4:80),fmt=*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%pawver 360 if (pspheads(ipsp)%pawheader%pawver==1) then ! Compatibility with Abinit v4.2.x 361 read (unit=pspline,fmt=*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%basis_size,& 362 pspheads(ipsp)%pawheader%lmn_size 363 ABI_MALLOC(orb,(pspheads(ipsp)%pawheader%basis_size)) 364 orb(:)=0 365 read (unt,*, err=10, iomsg=errmsg) (orb(ii), ii=1,pspheads(ipsp)%pawheader%basis_size) 366 read (unt,*, err=10, iomsg=errmsg) 367 read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%rpaw 368 pspheads(ipsp)%pawheader%rshp=pspheads(ipsp)%pawheader%rpaw 369 read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%mesh_size 370 read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%shape_type 371 if (pspheads(ipsp)%pawheader%shape_type==3) pspheads(ipsp)%pawheader%shape_type=-1 372 else 373 read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%basis_size,pspheads(ipsp)%pawheader%lmn_size 374 ABI_MALLOC(orb,(pspheads(ipsp)%pawheader%basis_size)) 375 orb(:)=0 376 read (unt,*, err=10, iomsg=errmsg) (orb(ii), ii=1,pspheads(ipsp)%pawheader%basis_size) 377 pspheads(ipsp)%pawheader%mesh_size=mmax 378 read (unt,*, err=10, iomsg=errmsg) nmesh 379 do ii=1,nmesh 380 read(unt,*, err=10, iomsg=errmsg) 381 end do 382 read (unt,*, err=10, iomsg=errmsg) pspheads(ipsp)%pawheader%rpaw 383 pspheads(ipsp)%pawheader%rshp=pspheads(ipsp)%pawheader%rpaw 384 read (unt,'(a80)', err=10, iomsg=errmsg) pspline 385 pspline=adjustl(pspline); write(std_out,*) pspline 386 read(unit=pspline,fmt=*) pspheads(ipsp)%pawheader%shape_type 387 if (pspheads(ipsp)%pawheader%pawver==2.and.& 388 pspheads(ipsp)%pawheader%shape_type==3) pspheads(ipsp)%pawheader%shape_type=-1 389 if (pspheads(ipsp)%pawheader%pawver>=3.and.pspheads(ipsp)%pawheader%shape_type==-1) then 390 rr=zero;read(unit=pspline,fmt=*,err=20,end=20) ii,rr 391 20 continue 392 if (rr>=tol8) pspheads(ipsp)%pawheader%rshp=rr 393 end if 394 end if 395 do ilmax=0,lmax 396 do ii=1,pspheads(ipsp)%pawheader%basis_size 397 if(orb(ii)==ilmax) nproj(ilmax)=nproj(ilmax)+1 398 end do 399 end do 400 pspheads(ipsp)%pawheader%l_size=2*maxval(orb)+1 401 pspheads(ipsp)%xccc=1 ! We suppose apriori that cc is used (but n1xccc is not used in PAW) 402 ABI_FREE(orb) 403 404 #if defined HAVE_BIGDFT 405 ! WVL+PAW case, need to define GTHradii 406 if(pspheads(ipsp)%usewvl==1) then 407 ! Obtain the HGH parameters by default from BigDFT 408 409 call atomic_info(int(pspheads(ipsp)%znuclpsp), int(pspheads(ipsp)%zionpsp), & 410 symbol = symbol, ehomo = ehomo, rcov = rcov, nsccode = iasctype) 411 412 ! I use the XC: Perdew, Burke & Ernzerhof as default, since 413 ! other XC potentials may not be in the BigDFT table. 414 ixc_=1 415 call psp_from_data(symbol, nzatom, nelpsp, npspcode_, ixc_, psppar, exists) 416 if(.not. exists) then 417 write(msg,'(4a)')ch10,& 418 "Chemical element not found in BigDFT table",ch10,& 419 "Action: upgrade BigDFT table" 420 ABI_BUG(msg) 421 end if 422 ! 423 ! pspheads(ipsp)%pawheader%rpaw/4.0d0 424 pspheads(ipsp)%GTHradii(0)=psppar(0,0) !rloc 425 pspheads(ipsp)%GTHradii(1)=psppar(1,0) !rrs 426 pspheads(ipsp)%GTHradii(2)=psppar(2,0) !rrp 427 ! pspheads(ipsp)%GTHradii(1) = one / sqrt(abs(two * ehomo)) 428 ! write(*,*)pspheads(ipsp)%GTHradii(:) 429 end if 430 #endif 431 432 else if(pspcod==8)then 433 434 ! DRH pseudopotentials 435 read(unt,*, err=10, iomsg=errmsg) rchrg,fchrg,qchrg 436 if (fchrg>1.d-15) pspheads(ipsp)%xccc=n1xccc_default 437 read(unt,*, err=10, iomsg=errmsg) nproj(0:lmax) 438 read(unt,*, err=10, iomsg=errmsg) extension_switch 439 if(any(extension_switch == [2, 3])) then 440 pspso=2 441 read(unt,*,err=10,iomsg=errmsg) nprojso(1:lmax) 442 else 443 pspso=0 444 end if 445 pspheads(ipsp)%pspso=pspso 446 447 else if(pspcod==9)then 448 ! placeholder: nothing to do everything is read above 449 450 else if(pspcod==10)then 451 452 ! HGH pseudopotentials, full h/k matrices 453 read (unt,*,err=10,iomsg=errmsg) pspheads(ipsp)%GTHradii(0) !rloc 454 read (unt,*,err=10,iomsg=errmsg) idum 455 if(idum-1/=lmax) then 456 ABI_ERROR("in inpspheads: nnonloc-1 /= lmax") 457 end if 458 do ilmax=0,lmax 459 read (unt,*,err=10,iomsg=errmsg) & 460 pspheads(ipsp)%GTHradii(ilmax + 1),nproj(ilmax),(hdum(idum),idum=1,nproj(ilmax)) 461 do idum=2,nproj(ilmax) !skip the rest of h_ij 462 read (unt,*, err=10, iomsg=errmsg) 463 end do 464 if (ilmax==0) cycle 465 nprojso(ilmax)=nproj(ilmax) 466 if(nprojso(ilmax)>0)then 467 pspheads(ipsp)%pspso=2 468 do idum=1,nprojso(ilmax) !skip the rest of k_ij 469 read (unt,*, err=10, iomsg=errmsg) 470 end do 471 end if 472 end do 473 474 else if (any(pspcod == [11, 12, 17])) then 475 ! already done above 476 477 else 478 write(msg, '(a,i0,4a)' )& 479 'The pseudopotential code (pspcod) read from file is ',pspcod,ch10,& 480 'This value is not allowed.',ch10,& 481 'Action: use a correct pseudopotential file.' 482 ABI_ERROR(msg) 483 end if ! pspcod 484 485 ! Store in pspheads 486 if (pspcod /= 17) then 487 pspheads(ipsp)%nproj(0:3)=nproj(0:3) 488 pspheads(ipsp)%nprojso(1:3)=nprojso(1:3) 489 end if 490 !write(std_out,'(a,*(i0,1x))') 'nproj = ', pspheads(ipsp)%nproj(:) 491 !write(std_out,'(a,*(i0,1x))') 'nprojso = ', pspheads(ipsp)%nprojso(:) 492 493 close(unt) 494 495 ! Compute md5 checksum 496 pspheads(ipsp)%md5_checksum = md5_sum_from_file(filnam(ipsp)) 497 end do ! ipsp=1,npsp 498 499 ! Note that mpsang is the max of 1+lmax, with minimal value 1 (even for local psps, at present) 500 mpsang=1 501 n1xccc=pspheads(1)%xccc 502 do ii=1,npsp 503 mpsang=max(pspheads(ii)%lmax+1,mpsang) 504 n1xccc=max(pspheads(ii)%xccc,n1xccc) 505 end do 506 507 write(msg,'(2a,i0,a,i0,a)')ch10,' inpspheads: deduce mpsang = ',mpsang,', n1xccc = ',n1xccc,'.' 508 call wrtout(std_out,msg,'PERS') 509 510 ! Test: if one psp is PAW, all must be 511 if (test_paw==1) then 512 do ipsp=1,npsp 513 if (all(pspheads(ipsp)%pspcod /= [7, 17])) then 514 write(msg, '(5a)' )& 515 'One pseudopotential is PAW (pspcod=7 or 17) !',ch10,& 516 'All pseudopotentials must be PAW (this is not the case here) !',ch10,& 517 'Action: use only PAW pseudopotential files.' 518 ABI_ERROR(msg) 519 end if 520 end do 521 end if 522 523 return 524 525 ! Handle IO error 526 10 continue 527 ABI_ERROR(errmsg) 528 529 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
SOURCE
715 subroutine pawpsxml2ab(filnam, ecut_tmp, pspheads, option) 716 717 !Arguments ------------------------------------ 718 !scalars 719 integer, intent(in) :: option 720 character(len=fnlen), intent(in) :: filnam 721 type(pspheader_type),intent(inout) :: pspheads !vz_i 722 !arrays 723 real(dp),intent(inout) :: ecut_tmp(3,2) 724 725 !Local variables------------------------------- 726 integer :: ii,il,lloc,lmax,pspcod,pspxc 727 real(dp) :: r2well,zionpsp,znuclpsp 728 ! character(len=100) :: xclibxc, msg 729 730 ! ********************************************************************* 731 732 if (option==1) then 733 call rdpawpsxml_header(ecut_tmp,filnam,paw_setuploc) 734 paw_setuploc%idgrid= paw_setuploc%radial_grid(1)%id 735 else 736 call rdpawpsxml(filnam,paw_setuploc) 737 end if 738 739 call pawpsp_read_header_xml(lloc,lmax,pspcod, pspxc,paw_setuploc,r2well,zionpsp,znuclpsp) 740 741 pspheads%lmax=lmax 742 pspheads%pspxc=pspxc 743 pspheads%zionpsp=zionpsp 744 pspheads%znuclpsp=znuclpsp 745 746 call pawpsp_read_pawheader(pspheads%pawheader%basis_size,& 747 pspheads%lmax,pspheads%pawheader%lmn_size,& 748 pspheads%pawheader%l_size,pspheads%pawheader%mesh_size,& 749 pspheads%pawheader%pawver,paw_setuploc,pspheads%pawheader%rpaw,& 750 pspheads%pawheader%rshp,pspheads%pawheader%shape_type) 751 752 pspheads%nproj=0 753 do il=0,pspheads%lmax 754 do ii=1,pspheads%pawheader%basis_size 755 if(paw_setuploc%valence_states%state(ii)%ll==il) pspheads%nproj(il)=pspheads%nproj(il)+1 756 end do 757 end do 758 759 pspheads%nprojso=0 760 pspheads%pspdat=27061961 761 pspheads%pspso=1 762 pspheads%xccc=1 763 pspheads%title=paw_setuploc%atom%symbol 764 765 if (option==1) call paw_setup_free(paw_setuploc) 766 767 end subroutine pawpsxml2ab
m_pspheads/updft_to_ixc [ Functions ]
[ Top ] [ m_pspheads ] [ Functions ]
NAME
updft_to_ixc
FUNCTION
Returns the abinit internal `ixc` from `dft` string with XC functional in QE format.
SOURCE
1182 integer function upfdft_to_ixc(dft, ixc, msg) result(ierr) 1183 1184 !Arguments ------------------------------------ 1185 character(len=*),intent(in) :: dft 1186 character(len=*),intent(out) :: msg 1187 integer,intent(out) :: ixc 1188 1189 !Local variables------------------------------- 1190 integer :: start !, ii 1191 character(len=500) :: x_name, c_name, gcx_name, gcc_name 1192 1193 !************************************************************************* 1194 1195 ! This list taken from oncvpsp/src/upfout.f90 1196 ! It should be OK as long as the UPF2 NC pseudos are generated with oncvpsp 1197 ! but it does not cover all QE possibilities. 1198 ierr = 0; msg = "" 1199 ixc = 0 1200 select case (dft) 1201 case ("PZ") 1202 ixc = -001009 1203 case ("PBE") 1204 ixc = -101130 !; ixc = 11 1205 case ("PW91") 1206 ixc = -109134 1207 case ("PBESOL") 1208 ixc = -116133 1209 case ("REVPBE") 1210 ixc = -102130 1211 case ("BP") 1212 ixc = -106132 1213 case ("BLYP") 1214 ixc = -106131 1215 case ("WC") 1216 ixc = -118130 1217 case ('SLA PW NOGX NOGC') ! string produced by oncvpsp3 1218 ixc = -1012 1219 case default 1220 ierr = 1 1221 end select 1222 1223 ! Extract substrings with 1224 ! 1) exchange 1225 ! 2) correlation 1226 ! 3) gradient correction, exchange 1227 ! 4) gradient correction, correlation 1228 if (ierr == 1) then 1229 ierr = 0; start = 1 1230 ABI_CHECK(next_token(dft, start, x_name) == 0 , "Error reading x_name") 1231 ABI_CHECK(next_token(dft, start, c_name) == 0 , "Error reading c_name") 1232 ABI_CHECK(next_token(dft, start, gcx_name) == 0 , "Error reading gcx_name") 1233 ABI_CHECK(next_token(dft, start, gcc_name) == 0 , "Error reading gcc_name") 1234 !call remove_non_ascii(gcc_name) 1235 !print *, "dft: `", trim(dft), "`" 1236 !print *, "x_name: `", trim(x_name), "`, c_name: `", trim(c_name), & 1237 ! "`, gcx_name: `", trim(gcx_name), "`, gcc_name: `", trim(gcc_name), "`" 1238 1239 if (x_name == "SLA" .and. c_name == "PW") then 1240 !print *, "in first if", gcx_name == "NOGX", trim(gcc_name) == "NOGC" 1241 !print *, "len_trim(gcc_name)", len_trim(gcc_name) 1242 if (gcx_name == "NOGX" .and. gcc_name == "NOGC") then 1243 ixc = -1012 1244 else 1245 ierr = 1 1246 !print *, "in second ierr" 1247 end if 1248 else 1249 ierr = 1 1250 end if 1251 end if 1252 1253 if (ierr == 1) then 1254 write(msg, "(5a)") & 1255 "Cannot find ABINIT ixc value corresponding to QE dft string: `", trim(dft), "`", ch10, & 1256 "Please update mapping in m_pspheads/upfdft_to_ixc." 1257 end if 1258 1259 end function upfdft_to_ixc
m_pspheads/upf1_to_psphead [ Functions ]
[ Top ] [ m_pspheads ] [ Functions ]
NAME
upf1_to_psphead
FUNCTION
This routine wraps a call to a PWSCF module, which reads in a UPF1 (PWSCF / Espresso) format pseudopotential, then transfers data for the HEADER of abinit psps only!
INPUTS
filpsp = name of file with UPF1 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
SOURCE
793 subroutine upf1_to_psphead(filpsp, znucl, zion, pspxc, lmax_, n1xccc, nproj_l, nprojso_l) 794 795 use m_read_upf_pwscf, only : read_pseudo 796 use pseudo_pwscf ! pwscf module with all data explicit! 797 798 !Arguments ------------------------------- 799 character(len=fnlen), intent(in) :: filpsp 800 integer,intent(inout) :: n1xccc 801 integer,intent(out) :: pspxc, lmax_ 802 real(dp),intent(out) :: znucl, zion 803 !arrays 804 integer,intent(out) :: nproj_l(0:3) 805 integer,intent(out) :: nprojso_l(1:3) 806 807 !Local variables ------------------------- 808 integer :: iproj, ll, iunit 809 character(len=500) :: msg 810 type(atomdata_t) :: atom 811 812 ! ********************************************************************* 813 814 if (open_file(filpsp, msg, newunit=iunit, status='old',form='formatted') /= 0) then 815 ABI_ERROR(msg) 816 end if 817 818 ! read in psp data to static data in pseudo module, for ipsx == 1 819 call read_pseudo(1,iunit) 820 close (iunit) 821 822 ! copy over to abinit internal arrays and vars 823 ! FIXME: The API is broken. It does not recognize PBEsol 824 ! should use upfdft_to_ixc 825 call upfxc2abi(dft(1), pspxc) 826 lmax_ = lmax(1) 827 call atomdata_from_symbol(atom,psd(1)) 828 znucl = atom%znucl 829 zion = zp(1) 830 831 nproj_l = 0 832 do iproj = 1, nbeta(1) 833 ll = lll(iproj,1) 834 nproj_l(ll) = nproj_l(ll) + 1 835 end do 836 837 nprojso_l = 0 !FIXME deal with so 838 !do iproj = 1, nbeta(1) 839 !nprojso_l(ll+1) = nprojso_l(ll+1) + 1 840 !end do 841 842 if (.not. nlcc(1)) n1xccc = 0 843 844 end subroutine upf1_to_psphead
m_pspheads/upf2_jl2srso [ Functions ]
[ Top ] [ m_pspheads ] [ Functions ]
NAME
upf2_jl2srso
FUNCTION
INPUTS
OUTPUT
SOURCE
958 subroutine upf2_jl2srso(upf, nproj_l, nprojso_l, vsr, esr, vso, eso) 959 960 !Arguments ------------------------------- 961 type(pseudo_upf),intent(in) :: upf 962 !arrays 963 integer,intent(out) :: nproj_l(0:3), nprojso_l(1:3) 964 real(dp),allocatable,intent(out) :: vsr(:,:,:), esr(:,:), vso(:,:,:), eso(:,:) 965 966 !Local variables ------------------------- 967 integer :: iprj, ii, ll, l1, il, ik, lmax, mmax, mxprj 968 real(dp) :: jtot, eprmin !eps_srso, 969 !character(len=500) :: msg 970 ! arrays 971 integer :: irc6(6),nproj6(6), done_ilk(6,2) 972 real(dp),allocatable :: vkb(:,:,:,:), evkb(:,:,:) 973 974 ! ********************************************************************* 975 976 lmax = upf%lmax; mmax = upf%mesh 977 nproj_l = 0; nprojso_l = 0 978 979 ! Pseudo in j = l + s representation. 980 irc6 = zero; nproj6 = zero 981 do iprj=1,upf%nbeta 982 ll = upf%lll(iprj) 983 nproj6(ll+1) = nproj6(ll+1) + 1 984 !irc6(ll+1) = max(upf%kbeta(iprj), irc6(ll+1)) 985 irc6(ll+1) = mmax 986 end do 987 988 ! Divide by two for l > 0 as this is sr_so_r expects. 989 nproj6(2:) = nproj6(2:) / 2 990 mxprj = maxval(nproj6) 991 992 ABI_MALLOC(vkb, (mmax,mxprj,4,2)) 993 ABI_MALLOC(evkb, (mxprj,4,2)) 994 ABI_MALLOC(vsr, (mmax,2*mxprj,4)) 995 ABI_MALLOC(esr, (2*mxprj,4)) 996 ABI_MALLOC(vso, (mmax,2*mxprj,4)) 997 ABI_MALLOC(eso, (2*mxprj,4)) 998 999 done_ilk = 0 1000 do iprj=1,upf%nbeta 1001 jtot = upf%jjj(iprj) 1002 ll = upf%lll(iprj) 1003 il = ll + 1 1004 if (ll == 0) then 1005 ik = 1 1006 else 1007 ! l+1/2 --> ik 1, l-1/2 --> ik 2 1008 if (abs(jtot - (ll + half)) < tol6) then 1009 ik = 1 1010 else if (abs(jtot - (ll - half)) < tol6) then 1011 ik = 2 1012 else 1013 ABI_ERROR(sjoin("Cannot detect ik index from jtot:", ftoa(jtot))) 1014 end if 1015 end if 1016 1017 done_ilk(il, ik) = done_ilk(il, ik) + 1 1018 ii = done_ilk(il, ik) 1019 evkb(ii,il,ik) = upf%dion(iprj,iprj) * half ! convert from Rydberg to Ha 1020 vkb(:,ii,il,ik) = upf%beta(:,iprj) 1021 end do 1022 1023 call sr_so_r(lmax, irc6, nproj6, upf%r, mmax, mxprj, evkb, vkb, vsr, esr, vso, eso) 1024 1025 ! MG: This is done in oncvpsp 3.3 but not in oncvpsp4 1026 ! drop sr, so orthonormal projectors with neglibible coefficients 1027 ! modify cutoff if desired 1028 1029 eprmin=2.0d-5 1030 write(std_out,'(/a,1p,e10.2,a)') 'Orthonormal projectors with coefficients <', & 1031 eprmin,' Ha will be dropped' 1032 1033 do l1=1,lmax+1 1034 if(abs(esr(3,l1))<eprmin) esr(3,l1)=0.0d0 1035 if(abs(esr(4,l1))<eprmin) esr(4,l1)=0.0d0 1036 if(abs(eso(3,l1))<eprmin) eso(3,l1)=0.0d0 1037 if(abs(eso(4,l1))<eprmin) eso(4,l1)=0.0d0 1038 end do 1039 1040 #if 0 1041 ! MG: This is done in oncvpsp 4 but not in oncvpsp 3.3 1042 ! set smallest components to zero (following the approach used in oncvpsp) 1043 eps_srso=1.0d-3 1044 do l1=1,lmax+1 1045 if (nproj6(l1) > 0) then 1046 do iprj=2,2*nproj6(l1) 1047 if (abs(esr(iprj,l1)) < eps_srso*abs(esr(1,l1))) esr(iprj,l1) = 0.0d0 1048 if (l1 == 1) cycle 1049 if (abs(eso(iprj,l1)) < eps_srso*abs(eso(1,l1))) eso(iprj,l1) = 0.0d0 1050 end do 1051 end if 1052 end do 1053 #endif 1054 1055 ! set up projector number for sr_so calculations based on non-zero coefficients 1056 ! note that energies and projectors have been sorted sr_so_r 1057 ! so the relevant projectors are packed in the first positions. 1058 do l1=1,lmax+1 1059 ll = l1 - 1 1060 do ii=1,2*nproj6(l1) 1061 if (abs(esr(ii,l1)) > 0.0d0) nproj_l(ll) = nproj_l(ll) + 1 1062 if (abs(eso(ii,l1)) > 0.0d0) nprojso_l(ll) = nprojso_l(ll) + 1 1063 end do 1064 !write(std_out, '(a,3(i0,1x))')' ll, nproj_l, nprojso_l',ll, nproj_l(ll), nprojso_l(ll) 1065 end do 1066 1067 ABI_FREE(vkb) 1068 ABI_FREE(evkb) 1069 1070 end subroutine upf2_jl2srso
m_pspheads/upf2_to_psphead [ Functions ]
[ Top ] [ m_pspheads ] [ Functions ]
NAME
upf2_to_psphead
FUNCTION
This routine wraps a call to a PWSCF module, which reads in a UPF2 (PWSCF / Espresso) format pseudopotential, then transfers data for the HEADER of abinit psps only!
INPUTS
filpsp = name of file with UPF1 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
SOURCE
870 subroutine upf2_to_psphead(filpsp, znucl, zion, pspxc, lmax, n1xccc, nproj_l, nprojso_l) 871 872 !Arguments ------------------------------- 873 character(len=fnlen), intent(in) :: filpsp 874 integer,intent(inout) :: n1xccc 875 integer,intent(out) :: pspxc, lmax 876 real(dp),intent(out) :: znucl, zion 877 !arrays 878 integer,intent(out) :: nproj_l(0:3), nprojso_l(1:3) 879 880 !Local variables ------------------------- 881 integer :: ierr , iprj, ll, mmax, irad 882 real(dp) :: amesh, damesh 883 character(len=500) :: msg 884 logical :: linear_mesh 885 type(pseudo_upf) :: upf 886 type(atomdata_t) :: atom 887 ! arrays 888 real(dp),allocatable :: vsr(:,:,:), esr(:,:), vso(:,:,:), eso(:,:) 889 890 ! ********************************************************************* 891 892 ! See also https://github.com/QEF/qeschemas/blob/master/UPF/qe_pp-0.99.xsd 893 call read_upf_new(filpsp, upf, ierr) 894 ABI_CHECK(ierr == 0, sjoin("read_upf_new returned ierr:", itoa(ierr))) 895 896 call atomdata_from_symbol(atom, upf%psd) 897 znucl = atom%znucl 898 zion = upf%zp 899 lmax = upf%lmax 900 mmax = upf%mesh 901 902 ! Consistency check 903 ABI_CHECK(upf%typ == "NC", sjoin("Only NC pseudos in UPF2 format are supported while type is:", upf%typ)) 904 ABI_CHECK(upfdft_to_ixc(upf%dft, pspxc, msg) == 0, msg) 905 if (.not. upf%nlcc) n1xccc = 0 906 907 ! Check that rad grid is linear starting at zero 908 linear_mesh = .True. 909 amesh = upf%r(2) - upf%r(1); damesh = zero 910 do irad=2,mmax-1 911 damesh = max(damesh, abs(upf%r(irad)+amesh-upf%r(irad+1))) 912 end do 913 linear_mesh = damesh < tol8 914 915 if (.not. linear_mesh .or. abs(upf%r(1)) > tol16) then 916 write(msg,'(3a)')& 917 'Assuming pseudized valence charge given on linear radial mesh starting at zero.',ch10,& 918 'Action: check your pseudopotential file.' 919 ABI_ERROR(msg) 920 end if 921 922 nproj_l = 0; nprojso_l = 0 923 924 if (.not. upf%has_so) then 925 ! Scalar case 926 do iprj=1,upf%nbeta 927 ll = upf%lll(iprj) 928 nproj_l(ll) = nproj_l(ll) + 1 929 end do 930 931 else 932 ! Pseudo in j = l + s representation. 933 call upf2_jl2srso(upf, nproj_l, nprojso_l, vsr, esr, vso, eso) 934 935 ABI_FREE(vsr) 936 ABI_FREE(esr) 937 ABI_FREE(vso) 938 ABI_FREE(eso) 939 end if 940 941 call deallocate_pseudo_upf(upf) 942 943 end subroutine upf2_to_psphead
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
SOURCE
1095 subroutine upfxc2abi(dft, pspxc) 1096 1097 use funct_pwscf ! pwscf module for naming xc functionals 1098 1099 !Arguments ------------------------------- 1100 character(len=*), intent(in) :: dft 1101 integer, intent(out) :: pspxc 1102 1103 !Local variables ------------------------- 1104 integer :: iexch,icorr,igcx,igcc 1105 integer :: totalindex, offset 1106 1107 ! ********************************************************************* 1108 !extract from char*20 :: dft(:) 1109 !### The following has been copied from pwscf src/Modules/upf_to_internal.f90: 1110 !workaround for rrkj format - it contains the indices, not the name 1111 if ( dft(1:6)=='INDEX:') then 1112 read( dft(7:10), '(4i1)') iexch,icorr,igcx,igcc 1113 call set_dft_from_indices(iexch,icorr,igcx,igcc) 1114 else 1115 call set_dft_from_name( dft ) 1116 iexch = get_iexch() 1117 icorr = get_icorr() 1118 igcx = get_igcx() 1119 igcc = get_igcc() 1120 end if 1121 1122 !reset dft string to avoid stray spaces 1123 call set_dft_from_indices(iexch,icorr,igcx,igcc) 1124 write(std_out,'(a)') ' upf2abinit: XC string from pseudopotential is :' 1125 write(std_out,'(3a)') '>', dft, '<' 1126 ABI_WARNING("upfxc2abi is not guaranteed to return the right ixc from QE XC string e.g. PBEsol. Please crosscheck!") 1127 1128 offset = 100 1129 totalindex = offset*offset*offset*iexch + offset*offset*icorr + offset*igcx + igcc 1130 select case (totalindex) 1131 case (00000000) !(" NOX NOC NOGX NOGC") ! no xc 1132 pspxc = 0 1133 case (01010000) !(" SLA PZ NOGX NOGC") ! slater exchange + Perdew Zunger 1134 pspxc = 2 1135 case (01050000) !(" SLA WIG NOGX NOGC") ! slater exchange + Wigner corr 1136 pspxc = 4 1137 case (01060000) !(" SLA HL NOGX NOGC") ! Hedin + Lundqvist 1138 pspxc = 5 1139 case (02000000) !(" SL1 NOC NOGX NOGC") ! full slater exchange 1140 pspxc = 6 1141 case (01040000) !(" SLA PW NOGX NOGC") ! slater exchange + Perdew Wang 1142 pspxc = 7 1143 case (01000000) !(" SLA NOC NOGX NOGC") ! Perdew Wang + no corr 1144 pspxc = 8 1145 case (01040304) !(" SLA PW PBX PBC") ! LDA + PBE GGA 1146 pspxc = 11 ! PBE 1147 case (01000300) !(" SLA NOC PBX NOGC") ! exchange part of PBE GGA 1148 pspxc = 12 1149 case (01040404) !(" SLA PW RPB PBC") ! rev PBE 1150 pspxc = 14 1151 case (00000505) !(" NOX NOC HTCH HTCH") ! HTCH 120 1152 pspxc = 17 1153 case (01030103) !(" SLA LYP B88 BLYP") ! BLYP 1154 pspxc = -106131 1155 case (01040101) !(" SLA PW B88 P86") ! BP86 1156 pspxc = -106132 1157 case (00030603) !(" NOX LYP OPTX BLYP") ! OLYP 1158 pspxc = -110131 1159 ! FIXME: important cases left to be patched with libxc: 1160 ! vosko wilkins nusair 1161 ! ortiz ballone 1162 ! pbe0 1163 ! Gunnarson-Lunqvist 1164 ! make general approach: check gradient parts first, then lda. 1165 ! event. check if they are consistent. 1166 case default 1167 ABI_ERROR('upf2abinit: XC functional not recognized') 1168 end select 1169 1170 end subroutine upfxc2abi