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