1! 2! Copyright (C) 2020 Quantum ESPRESSO group 3! This file is distributed under the terms of the 4! GNU General Public License. See the file `License' 5! in the root directory of the present distribution, 6! or http://www.gnu.org/copyleft/gpl.txt . 7! 8!-------------------------------------------------------- 9MODULE read_upf_new_module 10 !----------------------------------------------------- 11 !! this module contains the simplified code for reading 12 !! pseudopotential files in either UPF v.2 or xml 13 ! 14 USE xmltools 15 USE upf_kinds, ONLY: dp 16 USE pseudo_types, ONLY: pseudo_upf, pseudo_config 17 ! 18 LOGICAL :: v2 19 !! true if UPF v.2 version, false if new UPF with xml schema 20 INTEGER :: iun 21 !! unit for reading data 22 ! 23 PUBLIC 24 ! 25CONTAINS 26 ! 27 !------------------------------------------------+ 28 SUBROUTINE read_upf_new (filename, upf, ierr) ! 29 !---------------------------------------------+ 30 !! Reads pseudopotential in UPF format (either v.2 or upf_schema). 31 !! Derived-type variable *upf* store in output the data read from file. 32 !! File *filename* is opened and closed inside the routine 33 ! 34 IMPLICIT NONE 35 CHARACTER(len=*), INTENT(IN) :: filename 36 !! i/o filename 37 TYPE(pseudo_upf),INTENT(OUT) :: upf 38 !! the derived type storing the pseudo data 39 INTEGER, INTENT(OUT) :: ierr 40 !! ierr=0 : xml schema, ierr=-2: UPF v.2 41 !! ierr=-81: error reading PP file 42 ! 43 iun = xml_openfile ( filename ) 44 IF ( iun == -1 ) CALL upf_error('read_upf', 'cannot open file',1) 45 call xmlr_opentag ( 'qe_pp:pseudo', IERR = ierr ) 46 if ( ierr == 0 ) then 47 v2 =.false. 48 else if ( ierr == -1 ) then 49 rewind (iun) 50 call xmlr_opentag ( 'UPF', IERR = ierr ) 51 if ( ierr == 0 ) then 52 v2 =.true. 53 ierr = -2 54 CALL get_attr ( 'version', upf%nv ) 55 end if 56 end if 57 if ( ierr /= 0 .and. ierr /= -2 ) then 58 call xml_closefile( ) 59 ierr = -81 60 return 61 end if 62 ! 63 ! The header sections differ a lot between UPF v.2 and UPF with schema 64 ! 65 IF ( v2 ) THEN 66 CALL read_pp_header_v2 ( upf ) 67 ELSE 68 CALL read_pp_header_schema ( upf ) 69 END IF 70 ! compatibility 71 upf%is_gth = .false. 72 upf%is_multiproj = .true. 73 ! 74 ! From here on the format of v2 and schema do not differ much: 75 ! the most frequent difference is capitalization of tags 76 ! (see function capitalize_if_v2) 77 ! 78 CALL read_pp_mesh ( upf ) 79 ! 80 allocate ( upf%rho_atc(upf%mesh) ) 81 IF(upf%nlcc) then 82 CALL xmlr_readtag( capitalize_if_v2('pp_nlcc'), & 83 upf%rho_atc(:) ) 84 else 85 upf%rho_atc(:) = 0.0_dp 86 end if 87 IF( .NOT. upf%tcoulombp) then 88 allocate ( upf%vloc(upf%mesh) ) 89 CALL xmlr_readtag( capitalize_if_v2('pp_local'), & 90 upf%vloc(:), ierr ) 91 ! 92 ! existing PP files may have pp_nlcc first, pp_local later, 93 ! but also the other way round - check that everything was right 94 ! 95 if ( ierr /= 0 ) then 96 ierr = -81 97 return 98 end if 99 end if 100 ! 101 CALL read_pp_semilocal ( upf ) 102 ! 103 CALL read_pp_nonlocal ( upf ) 104 ! 105 CALL read_pp_pswfc ( upf ) 106 ! 107 CALL read_pp_full_wfc ( upf ) 108 ! 109 allocate( upf%rho_at(1:upf%mesh) ) 110 CALL xmlr_readtag( capitalize_if_v2('pp_rhoatom'), & 111 upf%rho_at(1:upf%mesh) ) 112 ! 113 CALL read_pp_spinorb ( upf ) 114 ! 115 CALL read_pp_paw ( upf ) 116 ! 117 CALL read_pp_gipaw ( upf ) 118 ! 119 ! close initial tag, qe_pp:pseudo or UPF 120 ! 121 CALL xmlr_closetag ( ) 122 ! 123 CALL xml_closefile ( ) 124 ! 125 END SUBROUTINE read_upf_new 126 ! 127 FUNCTION capitalize_if_v2 ( strin ) RESULT ( strout ) 128 ! 129 ! returns a capitalized string for UPF v.2, the same string otherwise 130 ! (UPF v.2 uses capitalized tags, UPF with schema use lowercase) 131 ! 132 USE upf_utils, ONLY: capital 133 IMPLICIT NONE 134 CHARACTER(LEN=*) :: strin 135 ! 136 INTEGER :: n 137 CHARACTER(LEN=:), ALLOCATABLE :: strout 138 ! 139 IF ( v2 ) THEN 140 strout = '' 141 DO n = 1,LEN_TRIM(strin) 142 strout = strout // capital(strin(n:n)) 143 END DO 144 ELSE 145 strout = TRIM(strin) 146 END IF 147 ! 148 END FUNCTION capitalize_if_v2 149 !-------------------------------------------------------- 150 SUBROUTINE read_pp_header_schema ( upf ) 151 !-------------------------------------------------------- 152 ! 153 IMPLICIT NONE 154 TYPE(pseudo_upf), INTENT(INOUT) :: upf ! the pseudo data 155 ! 156 CALL xmlr_opentag( capitalize_if_v2('pp_header') ) 157 ! 158 CALL xmlr_readtag( 'element', upf%psd ) 159 CALL xmlr_readtag( 'z_valence', upf%zp ) 160 CALL xmlr_readtag( 'type', upf%typ ) 161 CALL xmlr_readtag( 'functional', upf%dft ) 162 CALL xmlr_readtag( 'relativistic', upf%rel ) 163 CALL xmlr_readtag( 'is_ultrasoft', upf%tvanp ) 164 CALL xmlr_readtag( 'is_paw', upf%tpawp ) 165 CALL xmlr_readtag( 'is_coulomb', upf%tcoulombp ) 166 CALL xmlr_readtag( 'has_so', upf%has_so ) 167 CALL xmlr_readtag( 'has_wfc', upf%has_wfc ) 168 CALL xmlr_readtag( 'has_gipaw', upf%has_gipaw ) 169 CALL xmlr_readtag( 'paw_as_gipaw', upf%paw_as_gipaw) 170 CALL xmlr_readtag( 'core_correction', upf%nlcc) 171 CALL xmlr_readtag( 'total_psenergy', upf%etotps ) 172 CALL xmlr_readtag( 'wfc_cutoff', upf%ecutwfc ) 173 CALL xmlr_readtag( 'rho_cutoff', upf%ecutrho ) 174 CALL xmlr_readtag( 'l_max', upf%lmax ) 175 CALL xmlr_readtag( 'l_max_rho', upf%lmax_rho ) 176 CALL xmlr_readtag( 'l_local', upf%lloc ) 177 CALL xmlr_readtag( 'mesh_size', upf%mesh ) 178 CALL xmlr_readtag( 'number_of_wfc', upf%nwfc ) 179 CALL xmlr_readtag( 'number_of_proj', upf%nbeta ) 180 ! 181 CALL xmlr_closetag( ) 182 ! 183 END SUBROUTINE read_pp_header_schema 184 ! 185 !-------------------------------------------------------- 186 SUBROUTINE read_pp_header_v2 ( upf ) 187 !-------------------------------------------------------- 188 ! 189 IMPLICIT NONE 190 TYPE(pseudo_upf), INTENT(INOUT) :: upf ! the pseudo data 191 ! 192 CHARACTER(LEN=1) :: dummy 193 ! 194 CALL xmlr_readtag ( capitalize_if_v2('pp_header'), dummy ) 195 CALL get_attr ('generated', upf%generated) 196 CALL get_attr ('author', upf%author) 197 CALL get_attr ('date', upf%date) 198 CALL get_attr ('comment', upf%comment) 199 CALL get_attr ('element', upf%psd) 200 CALL get_attr ('pseudo_type', upf%typ) 201 CALL get_attr ('relativistic', upf%rel) 202 CALL get_attr ('is_ultrasoft', upf%tvanp) 203 CALL get_attr ('is_paw', upf%tpawp) 204 CALL get_attr ('is_coulomb', upf%tcoulombp) 205 CALL get_attr ('has_so', upf%has_so) 206 CALL get_attr ('has_wfc', upf%has_wfc) 207 CALL get_attr ('has_gipaw', upf%has_gipaw) 208 CALL get_attr ('paw_as_gipaw', upf%paw_as_gipaw) 209 CALL get_attr ('core_correction', upf%nlcc) 210 CALL get_attr ('functional', upf%dft) 211 CALL get_attr ('z_valence', upf%zp) 212 CALL get_attr ('total_psenergy', upf%etotps) 213 CALL get_attr ('wfc_cutoff', upf%ecutwfc) 214 CALL get_attr ('rho_cutoff', upf%ecutrho) 215 CALL get_attr ('l_max', upf%lmax) 216 CALL get_attr ('l_max_rho', upf%lmax_rho) 217 CALL get_attr ('l_local', upf%lloc) 218 CALL get_attr ('mesh_size', upf%mesh) 219 CALL get_attr ('number_of_wfc', upf%nwfc) 220 CALL get_attr ('number_of_proj', upf%nbeta ) 221 ! 222 END SUBROUTINE read_pp_header_v2 223 ! 224 !-------------------------------------------------------- 225 SUBROUTINE read_pp_mesh ( upf ) 226 !-------------------------------------------------------- 227 ! 228 IMPLICIT NONE 229 TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data 230 integer :: mesh 231 ! 232 CALL xmlr_opentag( capitalize_if_v2('pp_mesh') ) 233 CALL get_attr ( 'mesh', mesh ) 234 if ( mesh == 0 ) THEN 235 call upf_error('read_pp_mesh',& 236 'mesh size missing, using the one in header',-1) 237 else if ( mesh /= upf%mesh ) THEN 238 call upf_error('read_pp_mesh',& 239 'mismatch in mesh size, discarding the one in header',-1) 240 upf%mesh = mesh 241 end if 242 CALL get_attr ( 'dx' , upf%dx ) 243 CALL get_attr ( 'xmin', upf%xmin ) 244 CALL get_attr ( 'rmax', upf%rmax ) 245 CALL get_attr ( 'zmesh', upf%zmesh ) 246 allocate ( upf%r(1:upf%mesh) ) 247 CALL xmlr_readtag( capitalize_if_v2('pp_r'), upf%r(1:upf%mesh) ) 248 allocate ( upf%rab(1:upf%mesh) ) 249 CALL xmlr_readtag( capitalize_if_v2('pp_rab'), upf%rab(1:upf%mesh) ) 250 ! 251 CALL xmlr_closetag( ) ! end pp_mesh 252 ! 253 END SUBROUTINE read_pp_mesh 254 ! 255 !-------------------------------------------------------- 256 SUBROUTINE read_pp_semilocal ( upf ) 257 !-------------------------------------------------------- 258 ! 259 IMPLICIT NONE 260 TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data 261 ! 262 INTEGER :: nb, ind, l, j, ierr 263 CHARACTER(LEN=8) :: tag 264 real(dp), allocatable :: vnl(:) 265 ! 266 IF ( upf%typ == "SL" ) THEN 267 ! 268 IF ( upf%has_so ) then 269 ALLOCATE(upf%vnl(upf%mesh,0:upf%lmax,2)) 270 else 271 ALLOCATE(upf%vnl(upf%mesh,0:upf%lmax,1)) 272 end if 273 allocate ( vnl(1:upf%mesh) ) 274 CALL xmlr_opentag( capitalize_if_v2('pp_semilocal') ) 275 ! 276 tag = 'vnl' 277 DO nb = 1,upf%nbeta 278 IF ( v2 ) THEN 279 ! NOTA BENE: v2 format follows available PP files, written 280 ! using original write_upf_v2; not FoX-based write_upf_v2 281 IF ( nb - 1 == upf%lloc ) CYCLE 282 tag = 'PP_VNL.'//i2c(nb-1) 283 END IF 284 CALL xmlr_readtag( tag, vnl, ierr ) 285 if ( ierr /= 0 ) & 286 call upf_error('read_pp_semilocal','error reading SL PPs',1) 287 CALL get_attr ( 'l', l) 288 ind = 1 289 IF ( upf%has_so ) then 290 CALL get_attr ( 'j', j) 291 IF ( l > 0 .AND. ABS(j-l-0.5_dp) < 0.001_dp ) ind = 2 292 ! FIXME: what about spin-orbit case for v.2 upf? 293 if ( v2 ) & 294 call upf_error('read_pp_semilocal','check spin-orbit',1) 295 END IF 296 upf%vnl(:,l,ind) = vnl(:) 297 END DO 298 deallocate ( vnl ) 299 ! 300 CALL xmlr_closetag( ) ! end pp_semilocal 301 ! 302 END IF 303 ! 304 END SUBROUTINE read_pp_semilocal 305 ! 306 !-------------------------------------------------------- 307 SUBROUTINE read_pp_nonlocal ( upf ) 308 !-------------------------------------------------------- 309 ! 310 IMPLICIT NONE 311 TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data 312 ! 313 LOGICAL :: isnull 314 INTEGER :: nb, ind, l, l_, ln, lm, mb, nmb 315 CHARACTER(LEN=15) :: tag 316 REAL(dp), ALLOCATABLE :: aux(:) 317 ! 318 nb = upf%nbeta 319 IF ( nb == 0 ) nb = 1 320 ALLOCATE (upf%beta(upf%mesh,nb) ) 321 ALLOCATE (upf%els_beta(nb), & 322 upf%lll(nb), & 323 upf%kbeta(nb), & 324 upf%rcut(nb), & 325 upf%rcutus(nb), & 326 upf%dion(nb,nb), & 327 upf%qqq(nb,nb) ) 328 ! 329 IF (upf%has_so) ALLOCATE( upf%jjj(upf%nbeta)) 330 ! 331 IF ( upf%nbeta == 0 ) THEN 332 upf%nqf = 0 333 upf%nqlc= 0 334 upf%kkbeta = 0 335 upf%qqq_eps=-1.0_dp 336 RETURN 337 END IF 338 ! 339 CALL xmlr_opentag( capitalize_if_v2('pp_nonlocal') ) 340 ! 341 DO nb = 1,upf%nbeta 342 ! 343 IF ( v2 ) THEN 344 tag = 'PP_BETA.'//i2c(nb) 345 ELSE 346 tag = 'pp_beta' 347 END IF 348 CALL xmlr_readtag( tag, upf%beta(1:upf%mesh,nb) ) 349 CALL get_attr('index', mb) 350 ! not-so-strict test: index is absent or incorrect in some UPF v.2 files 351 IF ( .NOT. v2 .AND. nb /= mb ) & 352 CALL upf_error('read_pp_nonlocal','mismatch',nb) 353 CALL get_attr('label', upf%els_beta(nb)) 354 CALL get_attr('angular_momentum', upf%lll(nb)) 355 IF ( .NOT. v2 .AND. upf%has_so ) & 356 CALL get_attr('tot_ang_mom', upf%jjj(nb)) 357 CALL get_attr('cutoff_radius_index', upf%kbeta(nb)) 358 CALL get_attr('cutoff_radius', upf%rcut(nb)) 359 CALL get_attr('ultrasoft_cutoff_radius', upf%rcutus(nb)) 360 ! 361 END DO 362 ! 363 ! pp_dij (D_lm matrix) 364 ! 365 CALL xmlr_opentag( capitalize_if_v2 ('pp_dij') ) 366 READ(iun,*) upf%dion(1:upf%nbeta,1:upf%nbeta) 367 CALL xmlr_closetag( ) 368 ! 369 ! pp_augmentation 370 ! 371 IF (upf%tvanp .or. upf%tpawp) THEN 372 CALL xmlr_opentag( capitalize_if_v2('pp_augmentation') ) 373 ! 374 IF ( v2 ) THEN 375 CALL get_attr ( 'q_with_l', upf%q_with_l ) 376 CALL get_attr ( 'nqf', upf%nqf ) 377 CALL get_attr ( 'nqlc', upf%nqlc ) 378 IF (upf%tpawp) THEN 379 CALL get_attr ( 'shape', upf%paw%augshape ) 380 CALL get_attr ( 'cutoff_r', upf%paw%raug ) 381 CALL get_attr ( 'cutoff_r_index', upf%paw%iraug ) 382 CALL get_attr ( 'augmentation_epsilon', upf%qqq_eps ) 383 CALL get_attr ( 'l_max_aug', upf%paw%lmax_aug ) 384 ENDIF 385 ELSE 386 CALL xmlr_readtag( 'q_with_l', upf%q_with_l ) 387 CALL xmlr_readtag( 'nqf', upf%nqf ) 388 CALL xmlr_readtag( 'nqlc', upf%nqlc ) 389 IF (upf%tpawp) THEN 390 CALL xmlr_readtag( 'shape', upf%paw%augshape ) 391 CALL xmlr_readtag( 'cutoff_r', upf%paw%raug ) 392 CALL xmlr_readtag( 'cutoff_r_index', upf%paw%iraug ) 393 CALL xmlr_readtag( 'augmentation_epsilon', upf%qqq_eps ) 394 CALL xmlr_readtag( 'l_max_aug', upf%paw%lmax_aug ) 395 ENDIF 396 ENDIF 397 ! 398 CALL xmlr_opentag( capitalize_if_v2('pp_q') ) 399 READ(iun,*) upf%qqq(1:upf%nbeta,1:upf%nbeta) 400 CALL xmlr_closetag( ) 401 ! 402 IF ( upf%tpawp ) THEN 403 CALL xmlr_opentag( capitalize_if_v2('pp_multipoles') ) 404 ALLOCATE ( upf%paw%augmom(1:upf%nbeta,1:upf%nbeta,0:2*upf%lmax) ) 405 READ(iun,*) upf%paw%augmom(1:upf%nbeta,1:upf%nbeta,0:2*upf%lmax) 406 CALL xmlr_closetag () 407 ENDIF 408 ! 409 ! read polinomial coefficients for Q_ij expansion at small radius 410 ! 411 IF ( upf%nqlc == 0 ) upf%nqlc = 2*upf%lmax+1 412 ALLOCATE( upf%rinner( upf%nqlc ) ) 413 IF ( v2 .AND. upf%nqf > 0) THEN 414 ALLOCATE ( upf%qfcoef(upf%nqf, upf%nqlc, upf%nbeta, upf%nbeta) ) 415 CALL xmlr_opentag('PP_QFCOEF') 416 READ(iun,*) upf%qfcoef 417 CALL xmlr_closetag () 418 CALL xmlr_readtag('PP_RINNER',upf%rinner) 419 ELSE IF ( upf%nqf == 0 ) THEN 420 ALLOCATE( upf%qfcoef(1,1,1,1) ) 421 upf%qfcoef =0.0_dp 422 ENDIF 423 ! 424 ! Read augmentation charge Q_ij 425 ! 426 IF( upf%q_with_l ) THEN 427 ALLOCATE( upf%qfuncl(upf%mesh,upf%nbeta*(upf%nbeta+1)/2,0:2*upf%lmax) ) 428 upf%qfuncl(:,:,:) = 0.0_dp 429 ! NOTE: it would be wiser to dimension qfuncl as (:,:,0:upf%lmax) 430 ! and store the q_l(r) with index l=L/2 (see loop_on_l below) 431 ! This would save some storage and avoid "holes" in the array 432 ! that may be a source of trouble if not initialized to zero 433 ELSE 434 ALLOCATE ( upf%qfunc(upf%mesh,upf%nbeta*(upf%nbeta+1)/2) ) 435 upf%qfunc (:,:) = 0.0_dp 436 END IF 437 ALLOCATE ( aux(upf%mesh) ) 438 loop_on_nb: DO nb = 1,upf%nbeta 439 ln = upf%lll(nb) 440 loop_on_mb: DO mb = nb,upf%nbeta 441 lm = upf%lll(mb) 442 IF( upf%q_with_l ) THEN 443 loop_on_l: DO l = abs(ln-lm),ln+lm,2 ! only even terms 444 isnull = .FALSE. 445 IF( upf%tpawp ) isnull = (abs(upf%paw%augmom(nb,mb,l)) < upf%qqq_eps) 446 IF(isnull) CYCLE loop_on_l 447 IF ( v2 ) THEN 448 tag = 'PP_QIJL.'//i2c(nb)//'.'//i2c(mb)//'.'//i2c(l) 449 ELSE 450 tag = 'pp_qijl' 451 END IF 452 CALL xmlr_readtag( tag, aux ) 453 CALL get_attr ('composite_index', nmb) 454 IF ( nmb /= mb*(mb-1)/2 + nb ) & 455 CALL upf_error ('read_pp_nonlocal','mismatch',1) 456 CALL get_attr ('angular_momentum', l_) 457 IF ( l /= l_ ) CALL upf_error ('read_pp_nonlocal','mismatch',2) 458 upf%qfuncl(:,nmb,l) = aux(:) 459 IF (upf%tpawp) upf%qfuncl(upf%paw%iraug+1:,nmb,l) = 0._DP 460 ENDDO loop_on_l 461 ELSE 462 isnull = .FALSE. 463 IF ( upf%tpawp ) isnull = ( abs(upf%qqq(nb,mb)) < upf%qqq_eps ) 464 IF (isnull) CYCLE loop_on_mb 465 IF ( v2 ) THEN 466 tag = 'PP_QIJ.'//i2c(nb)//'.'//i2c(mb) 467 ELSE 468 tag = 'pp_qij' 469 END IF 470 CALL xmlr_readtag( tag, aux ) 471 CALL get_attr ('composite_index', nmb) 472 IF ( nmb /= mb*(mb-1)/2 + nb ) & 473 CALL upf_error ('read_pp_nonlocal','mismatch',3) 474 upf%qfunc(:,nmb) = aux(:) 475 ! 476 ENDIF 477 ENDDO loop_on_mb 478 ENDDO loop_on_nb 479 ! 480 DEALLOCATE (aux) 481 CALL xmlr_closetag( ) ! end pp_augmentation 482 ! 483 END IF 484 CALL xmlr_closetag( ) ! end pp_nonlocal 485 ! 486 ! Maximum radius of beta projector: outer radius to integrate 487 upf%kkbeta = MAXVAL(upf%kbeta(1:upf%nbeta)) 488 ! For PAW, augmentation charge may extend a bit further: 489 IF(upf%tpawp) upf%kkbeta = MAX(upf%kkbeta, upf%paw%iraug) 490 ! 491 END SUBROUTINE read_pp_nonlocal 492 ! 493 !-------------------------------------------------------- 494 SUBROUTINE read_pp_pswfc ( upf ) 495 !-------------------------------------------------------- 496 ! 497 IMPLICIT NONE 498 TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data 499 ! 500 INTEGER :: nw, ind, l 501 CHARACTER(LEN=8) :: tag 502 ! 503 allocate ( upf%chi(1:upf%mesh,upf%nwfc) ) 504 allocate ( upf%els(upf%nwfc), & 505 upf%oc(upf%nwfc), & 506 upf%lchi(upf%nwfc), & 507 upf%nchi(upf%nwfc), & 508 upf%rcut_chi(upf%nwfc), & 509 upf%rcutus_chi(upf%nwfc), & 510 upf%epseu(upf%nwfc) ) 511 IF ( upf%has_so ) THEN 512 allocate ( upf%nn(upf%nwfc) ) 513 allocate ( upf%jchi(upf%nwfc) ) 514 END IF 515 ! 516 CALL xmlr_opentag( capitalize_if_v2('pp_pswfc') ) 517 DO nw=1,upf%nwfc 518 IF ( v2 ) THEN 519 tag = 'PP_CHI.'//i2c(nw) 520 ELSE 521 tag = 'pp_chi' 522 END IF 523 CALL xmlr_readtag( tag, upf%chi(1:upf%mesh,nw) ) 524 call get_attr('index', ind) 525 ! not-so-strict test: index is absent or incorrect in some UPF v.2 files 526 if ( .NOT. v2 .AND. ind /= nw ) & 527 call upf_error('read_pp_pswfc','mismatch reading PSWFC', nw) 528 call get_attr( 'label', upf%els(nw) ) 529 call get_attr( 'l', upf%lchi(nw) ) 530 IF ( .not. v2 .and. upf%has_so ) THEN 531 call get_attr( 'nn', upf%nn(nw) ) 532 call get_attr( 'jchi', upf%jchi(nw) ) 533 END IF 534 call get_attr( 'occupation', upf%oc(nw) ) 535 call get_attr( 'n', upf%nchi(nw) ) 536 call get_attr( 'pseudo_energy', upf%epseu(nw) ) 537 call get_attr( 'cutoff_radius', upf%rcut_chi(nw) ) 538 call get_attr( 'ultrasoft_cutoff_radius', upf%rcutus_chi(nw) ) 539 END DO 540 CALL xmlr_closetag( ) ! end pp_pswfc 541 ! 542 END SUBROUTINE read_pp_pswfc 543 ! 544 !-------------------------------------------------------- 545 SUBROUTINE read_pp_full_wfc ( upf ) 546 !-------------------------------------------------------- 547 ! 548 IMPLICIT NONE 549 TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data 550 ! 551 INTEGER :: nb, mb 552 CHARACTER(LEN=15) :: tag 553 ! 554 IF ( upf%has_wfc ) THEN 555 ! 556 ALLOCATE (upf%aewfc(1:upf%mesh,upf%nbeta) ) 557 CALL xmlr_opentag( capitalize_if_v2('pp_full_wfc') ) 558 ! 559 DO nb = 1, upf%nbeta 560 IF ( v2 ) THEN 561 tag = 'PP_AEWFC.'//i2c(nb) 562 ELSE 563 tag = 'pp_aewfc' 564 END IF 565 CALL xmlr_readtag( tag, upf%aewfc(1:upf%mesh,nb) ) 566 CALL get_attr ('index',mb) 567 ! not-so-strict test (and two more below): 568 ! index may be absent or incorrect in some UPF v.2 files 569 IF ( .NOT. v2 .AND. nb /= mb ) CALL upf_error('read_pp_full_wfc','mismatch',1) 570 END DO 571 ! 572 IF ( upf%has_so .AND. upf%tpawp ) THEN 573 ALLOCATE (upf%paw%aewfc_rel(1:upf%mesh,upf%nbeta) ) 574 DO nb = 1, upf%nbeta 575 IF ( v2 ) THEN 576 tag = 'PP_AEWFC_rel.'//i2c(nb) 577 ELSE 578 tag = 'pp_aewfc_rel' 579 END IF 580 CALL xmlr_readtag(tag, upf%paw%aewfc_rel(1:upf%mesh,nb) ) 581 CALL get_attr ('index',mb) 582 IF ( .NOT. v2 .AND. nb /= mb ) CALL upf_error('read_pp_full_wfc','mismatch',2) 583 END DO 584 END IF 585 ! 586 ALLOCATE (upf%pswfc(1:upf%mesh,upf%nbeta) ) 587 DO nb = 1, upf%nbeta 588 IF ( v2 ) THEN 589 tag = 'PP_PSWFC.'//i2c(nb) 590 ELSE 591 tag = 'pp_pswfc' 592 END IF 593 CALL xmlr_readtag(tag, upf%pswfc(1:upf%mesh,nb) ) 594 CALL get_attr ('index',mb) 595 IF ( .NOT. v2 .AND. nb /= mb ) CALL upf_error('read_pp_full_wfc','mismatch',3) 596 END DO 597 ! 598 CALL xmlr_closetag( ) 599 ! 600 END IF 601 ! 602 END SUBROUTINE read_pp_full_wfc 603 ! 604 !-------------------------------------------------------- 605 SUBROUTINE read_pp_spinorb ( upf ) 606 !-------------------------------------------------------- 607 ! 608 IMPLICIT NONE 609 TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data 610 INTEGER :: nw, nb, ierr 611 CHARACTER(LEN=1) :: dummy 612 ! 613 IF ( .NOT. v2 .OR. .NOT. upf%has_so ) RETURN 614 ! 615 CALL xmlr_opentag( 'PP_SPIN_ORB' ) 616 DO nw = 1,upf%nwfc 617 CALL xmlr_readtag( 'PP_RELWFC.'//i2c(nw), dummy ) 618 CALL get_attr( 'index' , nb ) 619 ! not-so-strict test: index absent or incorrect in some UPF v.2 files 620 IF ( .NOT. v2 .AND. nb /= nw ) CALL upf_error('read_pp_spinorb','mismatch',1) 621 CALL get_attr( 'els', upf%els(nw) ) 622 CALL get_attr( 'nn', upf%nn(nw) ) 623 CALL get_attr( 'lchi', upf%lchi(nw) ) 624 CALL get_attr( 'jchi', upf%jchi(nw) ) 625 CALL get_attr( 'oc', upf%oc(nw) ) 626 ENDDO 627 ! 628 DO nb = 1,upf%nbeta 629 CALL xmlr_readtag( 'PP_RELBETA.'//i2c(nb), dummy, ierr ) 630 ! 631 ! existing PP files may have pp_relbeta first, pp_relwfc later, 632 ! but also the other way round - check that everything was right 633 ! 634 if ( ierr /= 0 ) then 635 ierr = -81 636 return 637 end if 638 CALL get_attr( 'index' , nw ) 639 IF ( .NOT.v2 .AND. nb /= nw ) CALL upf_error('read_pp_spinorb','mismatch',2) 640 CALL get_attr( 'lll', upf%lll(nb) ) 641 CALL get_attr( 'jjj', upf%jjj(nb) ) 642 ENDDO 643 CALL xmlr_closetag () ! end pp_spin_orb 644 ! 645 END SUBROUTINE read_pp_spinorb 646 ! 647 !-------------------------------------------------------- 648 SUBROUTINE read_pp_paw ( upf ) 649 !-------------------------------------------------------- 650 ! 651 IMPLICIT NONE 652 TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data 653 INTEGER :: nb, mb 654 ! 655 IF ( .NOT. upf%tpawp ) RETURN 656 ! 657 CALL xmlr_opentag( capitalize_if_v2('pp_paw') ) 658 CALL get_attr ('paw_data_format', upf%paw_data_format) 659 CALL get_attr ('core_energy', upf%paw%core_energy) 660 ! Full occupation (not only > 0 ones) 661 ALLOCATE (upf%paw%oc(upf%nbeta) ) 662 ALLOCATE (upf%paw%ae_rho_atc(upf%mesh) ) 663 ALLOCATE (upf%paw%ae_vloc(upf%mesh) ) 664 CALL xmlr_readtag( capitalize_if_v2('pp_occupations'), & 665 upf%paw%oc(1:upf%nbeta) ) 666 ! All-electron core charge 667 CALL xmlr_readtag( capitalize_if_v2('pp_ae_nlcc'), & 668 upf%paw%ae_rho_atc(1:upf%mesh) ) 669 ! All-electron local potential 670 CALL xmlr_readtag( capitalize_if_v2('pp_ae_vloc'), & 671 upf%paw%ae_vloc(1:upf%mesh) ) 672 CALL xmlr_closetag () ! end pp_paw 673 ! 674 ALLOCATE(upf%paw%pfunc(upf%mesh, upf%nbeta,upf%nbeta) ) 675 upf%paw%pfunc(:,:,:) = 0._dp 676 IF (upf%has_so) THEN 677 ALLOCATE(upf%paw%pfunc_rel(upf%mesh, upf%nbeta,upf%nbeta) ) 678 upf%paw%pfunc_rel(:,:,:) = 0._dp 679 ENDIF 680 DO nb=1,upf%nbeta 681 DO mb=1,nb 682 upf%paw%pfunc (1:upf%mesh, nb, mb) = & 683 upf%aewfc(1:upf%mesh, nb) * upf%aewfc(1:upf%mesh, mb) 684 IF (upf%has_so) THEN 685 upf%paw%pfunc_rel (1:upf%paw%iraug, nb, mb) = & 686 upf%paw%aewfc_rel(1:upf%paw%iraug, nb) * & 687 upf%paw%aewfc_rel(1:upf%paw%iraug, mb) 688! 689! The small component is added to pfunc. pfunc_rel is useful only 690! to add a small magnetic contribution 691! 692 upf%paw%pfunc (1:upf%paw%iraug, nb, mb) = & 693 upf%paw%pfunc (1:upf%paw%iraug, nb, mb) + & 694 upf%paw%pfunc_rel (1:upf%paw%iraug, nb, mb) 695 ENDIF 696 upf%paw%pfunc(upf%paw%iraug+1:,nb,mb) = 0._dp 697 ! 698 upf%paw%pfunc (1:upf%mesh, mb, nb) = upf%paw%pfunc (1:upf%mesh, nb, mb) 699 IF (upf%has_so) upf%paw%pfunc_rel (1:upf%mesh, mb, nb) = & 700 upf%paw%pfunc_rel (1:upf%mesh, nb, mb) 701 ENDDO 702 ENDDO 703 ! 704 ! Pseudo wavefunctions (not only the ones for oc > 0) 705 ! All-electron wavefunctions 706 ALLOCATE(upf%paw%ptfunc(upf%mesh, upf%nbeta,upf%nbeta) ) 707 upf%paw%ptfunc(:,:,:) = 0._dp 708 DO nb=1,upf%nbeta 709 DO mb=1,upf%nbeta 710 upf%paw%ptfunc (1:upf%mesh, nb, mb) = & 711 upf%pswfc(1:upf%mesh, nb) * upf%pswfc(1:upf%mesh, mb) 712 upf%paw%ptfunc(upf%paw%iraug+1:,nb,mb) = 0._dp 713 ! 714 upf%paw%ptfunc (1:upf%mesh, mb, nb) = upf%paw%ptfunc (1:upf%mesh, nb, mb) 715 ENDDO 716 ENDDO 717 ! 718 END SUBROUTINE read_pp_paw 719 !-------------------------------------------------------- 720 SUBROUTINE read_pp_gipaw ( upf ) 721 !-------------------------------------------------------- 722 ! 723 IMPLICIT NONE 724 TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data 725 ! 726 INTEGER :: nb, mb 727 CHARACTER(LEN=24) :: tag 728 ! 729 IF (.NOT. upf%has_gipaw) RETURN 730 ! 731 CALL xmlr_opentag( capitalize_if_v2('pp_gipaw') ) 732 CALL get_attr ('gipaw_data_format', upf%gipaw_data_format ) 733 IF ( v2 ) THEN 734 CALL xmlr_opentag( 'PP_GIPAW_CORE_ORBITALS') 735 CALL get_attr ('number_of_core_orbitals', upf%gipaw_ncore_orbitals) 736 ELSE 737 CALL xmlr_readtag ('number_of_core_orbitals', upf%gipaw_ncore_orbitals) 738 IF ( .NOT. upf%paw_as_gipaw) & 739 CALL xmlr_readtag( 'number_of_valence_orbitals', upf%gipaw_wfs_nchannels) 740 END IF 741 ALLOCATE ( upf%gipaw_core_orbital(upf%mesh,upf%gipaw_ncore_orbitals) ) 742 ALLOCATE ( upf%gipaw_core_orbital_n(upf%gipaw_ncore_orbitals) ) 743 ALLOCATE ( upf%gipaw_core_orbital_el(upf%gipaw_ncore_orbitals) ) 744 ALLOCATE ( upf%gipaw_core_orbital_l(upf%gipaw_ncore_orbitals) ) 745 DO nb = 1,upf%gipaw_ncore_orbitals 746 IF ( v2 ) THEN 747 tag = "PP_GIPAW_CORE_ORBITAL."//i2c(nb) 748 ELSE 749 tag = 'pp_gipaw_core_orbital' 750 END IF 751 CALL xmlr_readtag( tag, upf%gipaw_core_orbital(1:upf%mesh,nb) ) 752 CALL get_attr ('index', mb) 753 IF ( nb /= mb ) CALL upf_error('read_pp_gipaw','mismatch',1) 754 CALL get_attr ('label', upf%gipaw_core_orbital_el(nb) ) 755 CALL get_attr ('n', upf%gipaw_core_orbital_n(nb) ) 756 CALL get_attr ('l', upf%gipaw_core_orbital_l(nb) ) 757 END DO 758 IF ( v2 ) CALL xmlr_closetag ( ) 759 ! 760 IF ( upf%paw_as_gipaw) THEN 761 ! 762 ! PAW as GIPAW case: all-electron and pseudo-orbitals not read here 763 ! 764 upf%gipaw_wfs_nchannels = upf%nbeta 765 ALLOCATE ( upf%gipaw_wfs_el(upf%gipaw_wfs_nchannels) ) 766 ALLOCATE ( upf%gipaw_wfs_ll(upf%gipaw_wfs_nchannels) ) 767 ALLOCATE ( upf%gipaw_wfs_rcut(upf%gipaw_wfs_nchannels) ) 768 ALLOCATE ( upf%gipaw_wfs_rcutus(upf%gipaw_wfs_nchannels) ) 769 ALLOCATE ( upf%gipaw_wfs_ae(upf%mesh,upf%gipaw_wfs_nchannels) ) 770 ALLOCATE ( upf%gipaw_wfs_ps(upf%mesh,upf%gipaw_wfs_nchannels) ) 771 DO nb = 1,upf%gipaw_wfs_nchannels 772 upf%gipaw_wfs_el(nb) = upf%els_beta(nb) 773 upf%gipaw_wfs_ll(nb) = upf%lll(nb) 774 upf%gipaw_wfs_ae(:,nb) = upf%aewfc(:,nb) 775 ENDDO 776 DO nb = 1,upf%gipaw_wfs_nchannels 777 upf%gipaw_wfs_ps(:,nb) = upf%pswfc(:,nb) 778 ENDDO 779 ALLOCATE ( upf%gipaw_vlocal_ae(upf%mesh) ) 780 ALLOCATE ( upf%gipaw_vlocal_ps(upf%mesh) ) 781 upf%gipaw_vlocal_ae(:)= upf%paw%ae_vloc(:) 782 upf%gipaw_vlocal_ps(:)= upf%vloc(:) 783 DO nb = 1,upf%gipaw_wfs_nchannels 784 upf%gipaw_wfs_rcut(nb)=upf%rcut(nb) 785 upf%gipaw_wfs_rcutus(nb)=upf%rcutus(nb) 786 ENDDO 787 ! 788 ELSE 789 ! 790 ! Read valence all-electron and pseudo orbitals 791 ! 792 IF ( v2 ) THEN 793 CALL xmlr_opentag( 'PP_GIPAW_ORBITALS' ) 794 CALL get_attr( 'number_of_valence_orbitals', & 795 upf%gipaw_wfs_nchannels ) 796 END IF 797 ALLOCATE ( upf%gipaw_wfs_el(upf%gipaw_wfs_nchannels) ) 798 ALLOCATE ( upf%gipaw_wfs_ll(upf%gipaw_wfs_nchannels) ) 799 ALLOCATE ( upf%gipaw_wfs_rcut(upf%gipaw_wfs_nchannels) ) 800 ALLOCATE ( upf%gipaw_wfs_rcutus(upf%gipaw_wfs_nchannels) ) 801 ALLOCATE ( upf%gipaw_wfs_ae(upf%mesh,upf%gipaw_wfs_nchannels) ) 802 ALLOCATE ( upf%gipaw_wfs_ps(upf%mesh,upf%gipaw_wfs_nchannels) ) 803 ALLOCATE ( upf%gipaw_vlocal_ae(upf%mesh) ) 804 ALLOCATE ( upf%gipaw_vlocal_ps(upf%mesh) ) 805 DO nb = 1,upf%gipaw_wfs_nchannels 806 IF ( v2 ) THEN 807 tag = "PP_GIPAW_ORBITAL."//i2c(nb) 808 ELSE 809 tag = 'pp_gipaw_orbital' 810 END IF 811 CALL xmlr_opentag( tag ) 812 CALL get_attr ('index', mb) 813 IF ( nb /= mb ) CALL upf_error('read_pp_gipaw','mismatch',2) 814 CALL get_attr ('label', upf%gipaw_wfs_el(nb) ) 815 CALL get_attr ('l', upf%gipaw_wfs_ll(nb) ) 816 CALL get_attr ('cutoff_radius', upf%gipaw_wfs_rcut(nb) ) 817 CALL get_attr ('ultrasoft_cutoff_radius', upf%gipaw_wfs_rcutus(nb) ) 818 CALL xmlr_readtag( capitalize_if_v2('pp_gipaw_wfs_ae'), & 819 upf%gipaw_wfs_ae(1:upf%mesh,nb) ) 820 CALL xmlr_readtag( capitalize_if_v2('pp_gipaw_wfs_ps'),& 821 upf%gipaw_wfs_ps(1:upf%mesh,nb) ) 822 CALL xmlr_closetag () 823 END DO 824 IF ( v2 ) CALL xmlr_closetag( ) 825 ! 826 ! Read all-electron and pseudo local potentials 827 ! 828 CALL xmlr_opentag( capitalize_if_v2('pp_gipaw_vlocal') ) 829 CALL xmlr_readtag( capitalize_if_v2('pp_gipaw_vlocal_ae'), & 830 upf%gipaw_vlocal_ae(1:upf%mesh) ) 831 CALL xmlr_readtag( capitalize_if_v2('pp_gipaw_vlocal_ps'), & 832 upf%gipaw_vlocal_ps(1:upf%mesh) ) 833 CALL xmlr_closetag () 834 END IF 835 CALL xmlr_closetag () ! end pp_gipaw 836 ! 837 END SUBROUTINE read_pp_gipaw 838 ! 839END MODULE read_upf_new_module 840