1! 2! Copyright (C) 2008-2012 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 write_upf_schema_module 10 !---------------------------------------------------------------------------=! 11 ! this module handles the writing of pseudopotential data 12 ! ... declare modules 13#if defined (__use_fox) 14 USE upf_kinds, ONLY: DP 15 USE pseudo_types, ONLY: pseudo_upf, pseudo_config, deallocate_pseudo_config 16 USE Fox_wxml 17 ! 18 IMPLICIT NONE 19 ! 20 PRIVATE 21 PUBLIC :: write_upf_schema 22 23CONTAINS 24 25 !-------------------------------+ 26 SUBROUTINE write_upf_schema(xf, upf, conf, u_input) 27 !----------------------------+ 28 ! Write pseudopotential in UPF format version 2, using FoX 29 ! 30 IMPLICIT NONE 31 TYPE(xmlf_t),INTENT(INOUT) :: xf ! xmlfile for writing 32 TYPE(pseudo_upf),INTENT(IN) :: upf ! the pseudo data 33 ! optional: configuration used to generate the pseudopotential 34 TYPE(pseudo_config),OPTIONAL,INTENT(IN) :: conf 35 ! optional: unit pointing to input file containing generation data 36 INTEGER, OPTIONAL, INTENT(IN):: u_input 37 INTEGER :: irow 38 CHARACTER(LEN=*),PARAMETER :: QE_PP_URI="http://www.quantum-espresso.org/ns/qes/qe_pp-1.0",& 39 QE_PP_LOC=QE_PP_URI//" "//& 40 "http://www.quantum-espresso.org/ns/qes/qe_pp-1.0.xsd",& 41 XSD_VERSION = "QE_PP-1.0" 42 ! 43 ! Initialize the file 44 CALL xml_DeclareNameSpace(XF = xf, PREFIX = "xsi", nsUri="http://www.w3.org/2001/XMLSchema-instance") 45 CALL xml_DeclareNameSpace(XF = xf, PREFIX = "qe_pp", nsUri=QE_PP_URI) 46 CALL xml_NewElement(XF=xf, name="qe_pp:pseudo") 47 CALL xml_addAttribute(XF=xf, name="xsi:schemalocation", VALUE = QE_PP_LOC) 48 ! 49 CALL xml_NewElement(xf,"xsd_version") 50 CALL xml_AddCharacters(xf, XSD_VERSION) 51 CALL xml_endElement(xf, "xsd_version") 52 ! 53 ! 54 CALL write_info(xf, upf, conf, u_input) 55 ! Write machine-readable header 56 CALL write_header(xf, upf) 57 ! Write radial grid mesh 58 CALL write_mesh(xf, upf) 59 ! Write non-linear core correction charge 60 IF(upf%nlcc) THEN 61 CALL xml_newElement(xf, 'pp_nlcc') 62 CALL xml_addAttribute(xf, 'size', upf%mesh) 63 DO irow = 1, upf%mesh, 4 64 CALL xml_addNewLine(xf) 65 CALL xml_addCharacters(xf, upf%rho_atc(irow:min(irow-1+4,upf%mesh)) , fmt = 's16') 66 END DO 67 CALL xml_addNewLine(xf) 68 CALL xml_endElement(xf, 'pp_nlcc') 69 END IF 70 ! Write local potential 71 IF(.not. upf%tcoulombp) THEN 72 CALL xml_newElement(xf, 'pp_local') 73 CALL xml_addAttribute(xf,'size', upf%mesh) 74 DO irow = 1, upf%mesh, 4 75 CALL xml_addNewLine(xf) 76 CALL xml_addCharacters(xf, upf%vloc(irow:min(irow-1+4,upf%mesh)) , fmt = 's16') 77 END DO 78 CALL xml_addNewLine(xf) 79 CALL xml_endElement(xf, 'pp_local') 80 ENDIF 81 ! Write potentials in semilocal form (if existing) 82 IF ( upf%typ == "SL" ) CALL write_semilocal(xf, upf) 83 ! Write nonlocal components: projectors, augmentation, hamiltonian elements 84 CALL write_nonlocal(xf, upf) 85 ! Write initial pseudo wavefunctions 86 ! (usually only wfcs with occupancy > 0) 87 CALL write_pswfc(xf, upf) 88 ! If included, write all-electron and pseudo wavefunctions 89 CALL write_full_wfc(xf, upf) 90 ! Write valence atomic density (used for initial density) 91 CALL xml_NewElement(xf, 'pp_rhoatom') 92 CALL xml_addAttribute(xf,'size', upf%mesh) 93 DO irow = 1, upf%mesh, 4 94 CALL xml_addNewLine(xf) 95 CALL xml_addCharacters(xf, upf%rho_at(irow:min(irow-1+4,upf%mesh)) , fmt = 's16') 96 END DO 97 CALL xml_addNewLine(xf) 98 CALL xml_endElement(xf, 'pp_rhoatom') 99 ! Write additional data for PAW (All-electron charge, wavefunctions, vloc..) 100 CALL write_paw(xf, upf) 101 ! Write additional data for GIPAW reconstruction 102 CALL write_gipaw(xf, upf) 103 ! 104 ! Close the file (not the unit!) 105 CALL xml_endElement(xf,'qe_pp:pseudo') 106 CALL xml_Close(xf) 107 108 CONTAINS 109 ! 110 SUBROUTINE write_info(u, upf, conf, u_input) 111 ! Write human-readable header 112 IMPLICIT NONE 113#include "qe_version.h" 114 TYPE(xmlf_t),INTENT(INOUT) :: u ! write to xml file u 115 TYPE(pseudo_upf),INTENT(IN) :: upf ! the pseudo data 116 ! optional: configuration used to generate the pseudopotential 117 TYPE(pseudo_config),OPTIONAL,INTENT(IN) :: conf 118 INTEGER, OPTIONAL, INTENT(IN) :: u_input ! read input data from u_input 119 ! 120 INTEGER :: nb ! aux counter 121 INTEGER :: ierr ! /= 0 if something went wrong 122 CHARACTER(LEN=4096) :: char_buf 123 CHARACTER(len=256) :: line 124 LOGICAL :: opnd 125 ! 126 CALL xml_NewElement( u, "pp_info") 127 CALL xml_NewElement(u,"generated") 128 CALL xml_addCharacters(u, TRIM(upf%generated), parsed = .FALSE. ) 129 CALL xml_EndElement(u,"generated") 130 ! 131 CALL xml_NewElement(u,"creator") 132 CALL xml_addAttribute(u,name="NAME",value="QE Atomic Code") 133 CALL xml_addAttribute(u,name= "VERSION", value = version_number ) 134 CALL xml_addCharacters(u,TRIM(upf%author)) 135 CALL xml_EndElement(u, 'creator') 136 ! 137 CALL xml_NewElement(u,"created") 138 CALL xml_addAttribute(u, name="DATE", VALUE=TRIM(upf%date)) 139 CALL xml_endElement(u,'created') 140 IF ( PRESENT(u_input) ) THEN 141 ! 142 ! copy content of input file used in pseudopotential generation 143 ! 144 INQUIRE (unit=u_input, opened=opnd) 145 IF (opnd) THEN 146 CALL xml_NewElement(u,"input") 147 CALL xml_addAttribute(u, name = "program", value = "ld1.x") 148 149 REWIND (unit=u_input) 150 char_buf="" 151 read_write_loop: DO 152 READ (u_input, '(A)',end=20,err=25) line 153 char_buf = TRIM(char_buf) // new_line('a') // trim(line) 154 CYCLE read_write_loop 15525 CALL upf_error('write_upf_schema::write_inputfile', 'problem writing input data',-1) 15620 EXIT read_write_loop 157 END DO read_write_loop 158 char_buf = TRIM(char_buf)// new_line('a') 159 CALL xml_addCharacters(u, TRIM(char_buf), parsed = .FALSE.) 160 CALL xml_endElement(u, 'input') 161 ELSE 162 CALL upf_error('write_upf_v2::write_inputfile', 'input file not open',-1) 163 END IF 164 ! 165 END IF 166 167 168 CALL xml_NewElement(u,'type') 169 CALL xml_addCharacters(u,TRIM(upf%typ)) 170 CALL xml_EndElement(u,'type') 171 CALL xml_NewElement(u,'relativistic_effects') 172 IF (TRIM(upf%rel)=='full') THEN 173 CALL xml_addCharacters(u,'full') 174 ELSE IF (TRIM(upf%rel)=='scalar') THEN 175 CALL xml_addCharacters(u,'scalar') 176 ELSE 177 CALL xml_addCharacters(u,'none') 178 ENDIF 179 CALL xml_EndElement(u, 'relativistic_effects') 180 CALL xml_NewElement(u,'element') 181 CALL xml_addCharacters(u,TRIM(upf%psd)) 182 CALL xml_endElement(u,'element') 183 CALL xml_NewElement(u,'functional') 184 CALL xml_addCharacters(u,TRIM(upf%dft)) 185 CALL xml_endElement(u,'functional') 186 CALL xml_newElement(u,'suggested_basis') 187 CALL xml_addAttribute(u,name='ecutwfc',value = upf%ecutwfc) 188 IF (upf%tpawp .OR. upf%tvanp ) CALL xml_addAttribute(u,name='ecutrho',value=upf%ecutrho) 189 CALL xml_endElement(u,'suggested_basis') 190 DO nb =1, upf%nwfc 191 IF( upf%oc(nb) >= 0._dp) THEN 192 CALL xml_newElement(u, name = "valence_orbital") 193 CALL xml_addAttribute(u, name='nl', value = upf%els(nb)) 194 CALL xml_addAttribute(u, name = 'pn', value = upf%nchi(nb) ) 195 CALL xml_addAttribute(u, name = 'l', value = upf%lchi(nb) ) 196 ! 197 CALL xml_newElement (u, name = "occupation") 198 CALL xml_addCharacters(u, upf%oc(nb) ) 199 CALL xml_endElement(u, "occupation") 200 CALL xml_newElement(u, "Rcut") 201 CALL xml_addCharacters(u, upf%rcut_chi(nb)) 202 CALL xml_endElement(u, "Rcut") 203 IF (upf%rcutus_chi(nb) > 0.d0) THEN 204 CALL xml_newElement(u, "RcutUS") 205 CALL xml_addCharacters(u, upf%rcutus_chi(nb)) 206 CALL xml_endElement(u, "RcutUS") 207 END IF 208 CALL xml_newElement(u, "Epseu") 209 CALL xml_addCharacters(u, upf%epseu(nb)) 210 CALL xml_endElement(u,"Epseu") 211 CALL xml_endElement(u,"valence_orbital") 212 END IF 213 END DO 214 215 216 IF( present(conf) ) THEN 217 CALL xml_newElement(u, "generation_configuration") 218 char_buf = "" 219 DO nb = 1,conf%nwfs 220 WRITE(line, '(4x,a2,2i3,f6.2,2f11.3,1f13.6)') & 221 conf%els(nb), conf%nns(nb), & 222 conf%lls(nb), conf%ocs(nb), conf%rcut(nb), & 223 conf%rcutus(nb), conf%enls(nb) 224 char_buf = TRIM(char_buf) // new_line('a') // TRIM(line) 225 ENDDO 226 WRITE(line,'(4x,2a)') 'Pseudization used: ',TRIM(conf%pseud) 227 char_buf = TRIM(char_buf) // new_line('a') // TRIM(line) // new_line('a') 228 CALL xml_addCharacters(u, TRIM(char_buf), parsed = .FALSE.) 229 CALL xml_endElement(u,"generation_configuration") 230 ENDIF 231 232 IF(TRIM(upf%comment) /= ' ') THEN 233 CALL xml_addComment (u, TRIM(upf%comment), WS_SIGNIFICANT=.TRUE. ) 234 END IF 235 ! 236 ! 237 CALL xml_endElement(u,"pp_info") 238 ! 239 RETURN 240100 CALL upf_error('write_upf_schema::write_info', 'Writing pseudo file', 1) 241 ! 242 END SUBROUTINE write_info 243 ! 244 ! 245 SUBROUTINE write_header(u, upf) 246 IMPLICIT NONE 247 TYPE(xmlf_t), INTENT(INOUT) :: u ! i/o unit 248 TYPE(pseudo_upf),INTENT(IN) :: upf ! the pseudo data 249 INTEGER :: ierr ! /= 0 if something went wrong 250 251 ! 252 INTEGER :: nw 253 ! 254 ! Write HEADER section with some initialization data 255 ! 256 CALL xml_newElement(u, 'pp_header') 257 CALL xml_newElement(u,'element') 258 CALL xml_addCharacters(u, TRIM(upf%psd)) 259 CALL xml_endElement(u,'element') 260 ! 261 CALL xml_newElement(u,'z_valence') 262 CALL xml_addCharacters(u, upf%zp) 263 CALL xml_endElement(u,'z_valence') 264 ! 265 CALL xml_newElement(u, 'type') 266 CALL xml_addCharacters(u, TRIM(upf%typ)) 267 CALL xml_endElement(u,'type') 268 CALL xml_newElement(u,'functional') 269 CALL xml_addCharacters(u, TRIM(upf%dft)) 270 CALL xml_endElement(u,'functional') 271 ! 272 CALL xml_newElement(u, 'relativistic') 273 CALL xml_addCharacters(u, TRIM(upf%rel)) 274 CALL xml_endElement(u,'relativistic') 275 ! 276 CALL xml_newElement(u,'is_ultrasoft') 277 CALL xml_addCharacters(u, upf%tvanp) 278 CALL xml_endElement(u,'is_ultrasoft') 279 ! 280 CALL xml_newElement(u,'is_paw') 281 CALL xml_addCharacters(u, upf%tpawp) 282 CALL xml_endElement(u,'is_paw') 283 ! 284 CALL xml_newElement(u,'is_coulomb') 285 CALL xml_addCharacters(u, upf%tcoulombp) 286 CALL xml_endElement(u,'is_coulomb') 287 ! 288 CALL xml_newElement(u,'has_so') 289 CALL xml_addCharacters(u, upf%has_so) 290 CALL xml_endElement(u,'has_so') 291 ! 292 CALL xml_newElement(u,'has_wfc') 293 CALL xml_addCharacters(u, upf%has_wfc) 294 CALL xml_endElement(u,'has_wfc') 295 ! 296 CALL xml_newElement(u,'has_gipaw') 297 CALL xml_addCharacters(u, upf%has_gipaw) 298 CALL xml_endElement(u,'has_gipaw') 299 ! 300 !Emine 301 CALL xml_newElement(u,'paw_as_gipaw') 302 CALL xml_addCharacters(u, upf%paw_as_gipaw) 303 CALL xml_endElement(u,'paw_as_gipaw') 304 ! 305 CALL xml_newElement(u,'core_correction') 306 CALL xml_addCharacters(u, upf%nlcc) 307 CALL xml_endElement(u,'core_correction') 308 ! 309 CALL xml_newElement(u,'total_psenergy') 310 CALL xml_addCharacters(u, upf%etotps) 311 CALL xml_endElement(u,'total_psenergy') 312 ! 313 CALL xml_newElement(u,'wfc_cutoff') 314 CALL xml_addCharacters(u, upf%ecutwfc) 315 CALL xml_endElement(u,'wfc_cutoff') 316 ! 317 CALL xml_newElement(u,'rho_cutoff') 318 CALL xml_addCharacters(u, upf%ecutrho) 319 CALL xml_endElement(u,'rho_cutoff') 320 ! 321 CALL xml_newElement(u,'l_max') 322 CALL xml_addCharacters(u, upf%lmax) 323 CALL xml_endElement(u,'l_max') 324 ! 325 CALL xml_newElement(u,'l_max_rho') 326 CALL xml_addCharacters(u, upf%lmax_rho) 327 CALL xml_endElement(u,'l_max_rho') 328 ! 329 CALL xml_newElement(u,'l_local') 330 CALL xml_addCharacters(u, upf%lloc) 331 CALL xml_endElement(u,'l_local') 332 ! 333 CALL xml_newElement(u,'mesh_size') 334 CALL xml_addCharacters(u, upf%mesh) 335 CALL xml_endElement(u,'mesh_size') 336 ! 337 CALL xml_newElement(u,'number_of_wfc') 338 CALL xml_addCharacters(u, upf%nwfc) 339 CALL xml_endElement(u,'number_of_wfc') 340 ! 341 CALL xml_newElement(u,'number_of_proj') 342 CALL xml_addCharacters(u, upf%nbeta) 343 CALL xml_endElement(u,'number_of_proj') 344 ! 345 CALL xml_endElement(u, 'pp_header') 346 RETURN 347 END SUBROUTINE write_header 348 ! 349 SUBROUTINE write_mesh(u, upf) 350 IMPLICIT NONE 351 TYPE(xmlf_t), INTENT(INOUT) :: u ! i/o unit 352 TYPE(pseudo_upf),INTENT(IN) :: upf ! the pseudo data 353 INTEGER :: ierr ! /= 0 if something went wrong 354 355 ! 356 ! 357 CALL xml_NewElement(u, 'pp_mesh') 358 IF ( upf%dx .GT. 0.d0) CALL xml_addAttribute(u, name='dx', value = upf%dx , fmt = 's16') 359 IF ( upf%mesh .GT. 0 ) CALL xml_addAttribute(u, name='mesh', value = upf%mesh) 360 IF (upf%dx .GT. 0.d0) CALL xml_addAttribute(u, name='xmin', value = upf%xmin, fmt = 's16') 361 IF (upf%rmax .GT. 0.d0) CALL xml_addAttribute(u, name='rmax', value = upf%rmax, fmt = 's16') 362 IF (upf%zmesh .GT.0.d0 )CALL xml_addAttribute(u, name='zmesh',value = upf%zmesh, fmt = 's16') 363 364 CALL xml_NewElement(u, 'pp_r' ) 365 DO irow =1, upf%mesh, 4 366 CALL xml_addNewLine(u) 367 CALL xml_addCharacters(u, upf%r(irow:min(irow-1+4,upf%mesh) ) , fmt='s16') 368 END DO 369 CALL xml_addNewLine(xf) 370 CALL xml_endElement(u,'pp_r') 371 CALL xml_NewElement(u, 'pp_rab') 372 DO irow = 1, upf%mesh, 4 373 CALL xml_addNewLine(u) 374 CALL xml_addCharacters(u, upf%rab(irow:min(irow-1+4, upf%mesh)), fmt = 's16') 375 END DO 376 CALL xml_addNewLine(xf) 377 CALL xml_endElement(u, 'pp_rab') 378 ! 379 380 CALL xml_endElement(u, 'pp_mesh') 381 ! 382 RETURN 383 END SUBROUTINE write_mesh 384 ! 385 SUBROUTINE write_semilocal(u, upf) 386 IMPLICIT NONE 387 TYPE(xmlf_t),INTENT(INOUT) :: u ! i/o unit 388 TYPE(pseudo_upf),INTENT(IN) :: upf ! the pseudo data 389 INTEGER :: ierr ! /= 0 if something went wrong 390 ! 391 INTEGER :: nb, l, ind 392 ! 393 CALL xml_newElement(u, 'pp_semilocal') 394 ! 395 ! Write V_l(r) 396 DO nb = 1,upf%nbeta 397 l = upf%lll(nb) 398 ind = 1 399 CALL xml_newElement(u, 'vnl') 400 CALL xml_addAttribute(u, name = 'l', value = l) 401 IF ( upf%has_so ) THEN 402 CALL xml_addAttribute(u, 'j', upf%jjj(nb)) 403 IF ( l > 0 .AND. ABS (upf%jjj(nb)-l-0.5_dp) < 0.001_dp) ind = 2 404 ENDIF 405 DO irow = 1, upf%mesh, 4 406 CALL xml_addNewLine(u) 407 CALL xml_addCharacters(u, upf%vnl(irow:min(irow-1+4, upf%mesh),l, ind), fmt = 's16') 408 END DO 409 CALL xml_addNewLine(xf) 410 CALL xml_endElement(u, 'vnl') 411 END DO 412 ! 413 CALL xml_endElement(u, 'pp_semilocal') 414 ! 415 END SUBROUTINE write_semilocal 416 ! 417 SUBROUTINE write_nonlocal(u, upf) 418 IMPLICIT NONE 419 TYPE(xmlf_t),INTENT(INOUT) :: u ! xml file descriptor 420 TYPE(pseudo_upf),INTENT(IN) :: upf ! the pseudo data 421 INTEGER :: ierr ! /= 0 if something went wrong 422 ! 423 INTEGER :: nb,mb,ln,lm,l,nmb 424 LOGICAL :: isnull 425 REAL(DP),ALLOCATABLE :: tmp_dbuffer(:) 426 ! 427 IF (upf%tcoulombp) RETURN 428 ! 429 CALL xml_NewElement(u, 'pp_nonlocal') 430 ! 431 ! Write the projectors: 432 DO nb = 1,upf%nbeta 433 CALL xml_NewElement(u, 'pp_beta') 434 CALL xml_addAttribute(u, 'index', nb ) 435 CALL xml_addAttribute(u, 'label', upf%els_beta(nb)) 436 CALL xml_addAttribute(u, 'angular_momentum', upf%lll(nb)) 437 IF (upf%has_so) & 438 CALL xml_addAttribute(u, 'tot_ang_mom', upf%jjj(nb)) 439 CALL xml_addAttribute(u, 'cutoff_radius_index', upf%kbeta(nb)) 440 CALL xml_addAttribute(u, 'cutoff_radius', upf%rcut(nb), fmt = 's16') 441 CALL xml_addAttribute(u, 'ultrasoft_cutoff_radius',upf%rcutus(nb), fmt = 's16') 442 DO irow = 1, upf%mesh, 4 443 CALL xml_addNewLine(u) 444 CALL xml_addCharacters(u, upf%beta(irow:min(irow-1+4,upf%mesh), nb), fmt = 's16') 445 END DO 446 CALL xml_addNewLine(xf) 447 CALL xml_endElement(u, 'pp_beta') 448 ENDDO 449 ! 450 ! Write the hamiltonian terms D_ij 451 452 CALL xml_newElement(u, 'pp_dij') 453 CALL xml_addAttribute(u,'columns', upf%nbeta) 454 CALL xml_addAttribute(u, 'rows', upf%nbeta ) 455 nb = upf%nbeta*upf%nbeta 456 ALLOCATE(tmp_dbuffer(nb)) 457 tmp_dbuffer = reshape(upf%dion,[nb]) 458 DO irow = 1, nb, 4 459 CALL xml_addNewLine(u) 460 CALL xml_addCharacters(u, tmp_dbuffer(irow:min(irow-1+4,nb)), fmt ='s16') 461 END DO 462 CALL xml_addNewLine(xf) 463 DEALLOCATE(tmp_dbuffer) 464 CALL xml_endElement(u, 'pp_dij') 465 ! 466 ! Write the augmentation charge section 467 augmentation : & 468 IF(upf%tvanp .or. upf%tpawp) THEN 469 CALL xml_newElement(u, 'pp_augmentation') 470 CALL xml_newElement(u, 'nqf') 471 CALL xml_addCharacters(u, upf%nqf) 472 CALL xml_endElement(u, 'nqf') 473 CALL xml_newElement(u, 'q_with_l') 474 CALL xml_addCharacters(u, upf%q_with_l) 475 CALL xml_endElement(u, 'q_with_l') 476 CALL xml_newElement(u, 'nqlc') 477 CALL xml_addCharacters(u, upf%nqlc) 478 CALL xml_endElement(u, 'nqlc') 479 IF (upf%tpawp) THEN 480 CALL xml_newElement(u, 'shape') 481 CALL xml_addCharacters(u, upf%paw%augshape) 482 CALL xml_endElement(u, 'shape') 483 CALL xml_newElement(u, 'cutoff_r') 484 CALL xml_addCharacters(u, upf%paw%raug) 485 CALL xml_endElement(u, 'cutoff_r') 486 CALL xml_newElement(u, 'cutoff_r_index') 487 CALL xml_addCharacters(u, upf%paw%iraug) 488 CALL xml_endElement(u, 'cutoff_r_index') 489 CALL xml_newElement(u, 'l_max_aug') 490 CALL xml_addCharacters(u, upf%paw%lmax_aug) 491 CALL xml_endElement(u, 'l_max_aug') 492 CALL xml_newElement(u, 'augmentation_epsilon') 493 CALL xml_addCharacters(u, upf%qqq_eps) 494 CALL xml_endElement(u, 'augmentation_epsilon') 495 ENDIF 496 ! 497 ! Write the integrals of the Q functions 498 CALL xml_newElement(u, 'pp_q') 499 nb = upf%nbeta*upf%nbeta 500 CALL xml_addAttribute(u, 'size',nb) 501 ALLOCATE(tmp_dbuffer (nb)) 502 tmp_dbuffer = reshape(upf%qqq, [nb]) 503 DO irow =1, nb, 4 504 CALL xml_addNewLine(u) 505 CALL xml_addCharacters(u, tmp_dbuffer(irow:min(irow-1+4, nb)), fmt='s16') 506 END DO 507 CALL xml_addNewLine(xf) 508 DEALLOCATE(tmp_dbuffer) 509 CALL xml_endElement(u, 'pp_q') 510 ! 511 ! Write charge multipoles (only if PAW) 512 IF ( upf%tpawp ) THEN 513 CALL xml_addComment(u, 'augmentation charge multipoles ( only for PAW) ' //& 514 'multipole array dims = (nbeta,nbeta,2*lmax+1)') 515 CALL xml_newElement(u, 'pp_multipoles') 516 CALL xml_addAttribute(u, 'nbeta', upf%nbeta) 517 CALL xml_addAttribute(u, 'lmax', upf%lmax) 518 nb = upf%nbeta*upf%nbeta*(2*upf%lmax+1) 519 ALLOCATE (tmp_dbuffer(nb)) 520 tmp_dbuffer = reshape(upf%paw%augmom,[nb]) 521 DO irow = 1, nb, 4 522 CALL xml_addNewLine(u) 523 CALL xml_addCharacters(u, tmp_dbuffer(irow: min(irow-1+4, nb) ), fmt='s16') 524 END DO 525 CALL xml_addNewLine(u) 526 DEALLOCATE(tmp_dbuffer) 527 CALL xml_endElement(u, 'pp_multipoles') 528 ENDIF 529 ! 530 ! Write augmentation charge Q_ij 531 loop_on_nb: DO nb = 1,upf%nbeta 532 ln = upf%lll(nb) 533 loop_on_mb: DO mb = nb,upf%nbeta 534 lm = upf%lll(mb) 535 nmb = mb * (mb-1) /2 + nb 536 IF( upf%q_with_l ) THEN 537 loop_on_l: DO l = abs(ln-lm),ln+lm,2 ! only even terms 538 isnull = .FALSE. 539 IF( upf%tpawp ) isnull = (abs(upf%paw%augmom(nb,mb,l)) < upf%qqq_eps) 540 IF( isnull) CYCLE loop_on_l 541 CALL xml_NewElement(u, 'pp_qijl') 542 CALL xml_addAttribute(u, 'first_index', nb ) 543 CALL xml_addAttribute(u, 'second_index', mb) 544 CALL xml_addAttribute(u, 'composite_index', nmb) 545 CALL xml_addAttribute(u, 'angular_momentum', l) 546 CALL xml_addAttribute(u, 'size', upf%mesh ) 547 DO irow = 1, upf%mesh, 4 548 CALL xml_addNewLine(u) 549 CALL xml_addCharacters(u,upf%qfuncl(irow:min(irow-1+4,upf%mesh),nmb,l), fmt = 's16') 550 END DO 551 CALL xml_addNewLine(xf) 552 CALL xml_EndElement(u, 'pp_qijl') 553 ENDDO loop_on_l 554 ELSE 555 isnull = .FALSE. 556 IF ( upf%tpawp ) isnull = ( abs(upf%qqq(nb,mb)) < upf%qqq_eps ) 557 IF (isnull) CYCLE loop_on_mb 558 CALL xml_NewElement(u, 'pp_qij') 559 CALL xml_addAttribute(u, 'size', upf%mesh) 560 CALL xml_addAttribute(u, 'first_index', nb ) 561 CALL xml_addAttribute(u, 'second_index', mb) 562 CALL xml_addAttribute(u, 'composite_index', nmb) 563 DO irow = 1, upf%mesh, 4 564 CALL xml_addNewLine(u) 565 CALL xml_addCharacters(u, upf%qfunc(irow:min(irow-1+4, upf%mesh), nmb), fmt = 's16') 566 END DO 567 CALL xml_addNewLine(xf) 568 CALL xml_endElement(u, 'pp_qij') 569 ! 570 ENDIF 571 ENDDO loop_on_mb 572 ENDDO loop_on_nb 573 ! 574 CALL xml_endElement(u, 'pp_augmentation') 575 ! 576 ENDIF augmentation 577 ! 578 CALL xml_endElement(u, 'pp_nonlocal') 579 ! 580 RETURN 581 END SUBROUTINE write_nonlocal 582 ! 583 SUBROUTINE write_pswfc(u, upf) 584 IMPLICIT NONE 585 TYPE(xmlf_t),INTENT(INOUT) :: u ! i/o unit descriptor 586 TYPE(pseudo_upf),INTENT(IN) :: upf ! the pseudo data 587 INTEGER :: ierr ! /= 0 if something went wrong 588 589 ! 590 INTEGER :: nw 591 ! 592 CALL xml_newElement(u, 'pp_pswfc') 593 ! 594 DO nw = 1,upf%nwfc 595 CALL xml_newElement(u, 'pp_chi') 596 CALL xml_addAttribute(u, 'size', upf%mesh) 597 CALL xml_addAttribute(u, 'index', nw ) 598 CALL xml_addAttribute(u, 'label', upf%els(nw)) 599 CALL xml_addAttribute(u, 'l', upf%lchi(nw)) 600 IF ( upf%has_so) THEN 601 CALL xml_addAttribute(u, 'nn', upf%nn(nw)) 602 CALL xml_addAttribute (u, 'jchi', upf%jchi(nw)) 603 END IF 604 CALL xml_addAttribute(u, 'occupation', upf%oc(nw)) 605 IF ( upf%nchi(nw) .GT. upf%lchi(nw) ) CALL xml_addAttribute(u, 'n', upf%nchi(nw)) 606 IF ( upf%epseu(nw) .GT. 0.0_DP ) CALL xml_addAttribute(u, 'pseudo_energy', upf%epseu(nw)) 607 IF ( upf%rcut_chi(nw) .GT. 0.0_DP) CALL xml_addAttribute(u, 'cutoff_radius', upf%rcut_chi(nw)) 608 IF ( upf%rcutus_chi(nw) .GT. 0.0_DP) CALL xml_addAttribute(u, 'ultrasoft_cutoff_radius', upf%rcutus_chi(nw)) 609 DO irow =1, upf%mesh, 4 610 CALL xml_addNewLine(u) 611 CALL xml_addCharacters(u, upf%chi(irow:min(irow-1+4,upf%mesh),nw) , fmt = 's16') 612 END DO 613 CALL xml_addNewLine(xf) 614 CALL xml_endElement(u, 'pp_chi') 615 ENDDO 616 ! 617 CALL xml_endElement(u, 'pp_pswfc') 618 ! 619 RETURN 620 END SUBROUTINE write_pswfc 621 ! 622 SUBROUTINE write_full_wfc(u, upf) 623 IMPLICIT NONE 624 TYPE(xmlf_t),INTENT(INOUT) :: u ! i/o unit descriptor 625 TYPE(pseudo_upf),INTENT(IN) :: upf ! the pseudo data 626 INTEGER :: ierr ! /= 0 if something went wrong 627 ! 628 INTEGER :: nb 629 630 IF(.not. upf%has_wfc) RETURN 631 CALL xml_NewElement(u, 'pp_full_wfc') 632 CALL xml_addAttribute(u, 'number_of_wfc', upf%nbeta ) 633 ! All-electron wavefunctions corresponding to beta functions 634 DO nb = 1,upf%nbeta 635 CALL xml_NewElement(u, 'pp_aewfc') 636 CALL xml_AddAttribute(u, 'index', nb ) 637 CALL xml_AddAttribute(u, 'label', upf%els_beta(nb)) 638 CALL xml_AddAttribute(u, 'l', upf%lll(nb)) 639 DO irow = 1, upf%mesh, 4 640 CALL xml_addNewLine(u) 641 CALL xml_addCharacters(u, upf%aewfc(irow:min(irow-1+4, upf%mesh), nb), fmt = 's16') 642 END DO 643 CALL xml_addNewLine(xf) 644 CALL xml_endElement(u, 'pp_aewfc') 645 ENDDO 646 IF (upf%has_so.and.upf%tpawp) THEN 647 DO nb = 1,upf%nbeta 648 CALL xml_NewElement(u, 'pp_aewfc_rel') 649 CALL xml_AddAttribute(u, 'size', upf%mesh) 650 CALL xml_addAttribute(u, 'index', nb ) 651 CALL xml_addAttribute(u, 'label', upf%els_beta(nb)) 652 CALL xml_addAttribute(u, 'l', upf%lll(nb)) 653 CALL xml_addAttribute(u, 'j', upf%jjj(nb)) 654 DO irow = 1, upf%mesh, 4 655 CALL xml_addNewLine(u) 656 CALL xml_addCharacters(u, upf%paw%aewfc_rel(irow:min(irow-1+4, upf%mesh), nb) , fmt = 's16') 657 END DO 658 CALL xml_addNewLine(xf) 659 CALL xml_endElement(u, 'pp_aewfc_rel') 660 ENDDO 661 ENDIF 662 ! Pseudo wavefunctions 663 DO nb = 1,upf%nbeta 664 CALL xml_newElement(u, 'pp_pswfc') 665 CALL xml_addAttribute(u, 'size', upf%mesh) 666 CALL xml_addAttribute(u, 'index', nb ) 667 CALL xml_addAttribute(u, 'label', upf%els_beta(nb)) 668 CALL xml_addAttribute(u, 'l', upf%lll(nb)) 669 DO irow = 1, upf%mesh, 4 670 CALL xml_addNewLine(u) 671 CALL xml_AddCharacters(u, upf%pswfc(irow:min(irow-1+4, upf%mesh), nb) , fmt = 's16') 672 END DO 673 CALL xml_addNewLine(xf) 674 CALL xml_endElement(u, 'pp_pswfc') 675 ENDDO 676 ! Finalize 677 CALL xml_endELement(u, 'pp_full_wfc') 678 END SUBROUTINE write_full_wfc 679 680 SUBROUTINE write_paw(u, upf) 681 IMPLICIT NONE 682 TYPE(xmlf_t),INTENT(INOUT) :: u ! i/o unit descriptor 683 TYPE(pseudo_upf),INTENT(IN) :: upf ! the pseudo data 684 INTEGER :: ierr ! /= 0 if something went wrong 685 ! 686 INTEGER :: nb 687 688 IF (.not. upf%tpawp ) RETURN 689 CALL xml_newElement(u, 'pp_paw') 690 CALL xml_addAttribute(u, 'paw_data_format', upf%paw_data_format ) 691 CALL xml_addAttribute(u, 'core_energy', upf%paw%core_energy) 692 ! Full occupation (not only > 0 ones) 693 CALL xml_newElement(u, 'pp_occupations') 694 CALL xml_addAttribute(u, 'size', upf%nbeta) 695 DO irow = 1, upf%nbeta, 4 696 CALL xml_addNewLine(u) 697 CALL xml_addCharacters(u, upf%paw%oc(irow:min(irow-1+4,upf%nbeta)), fmt = 's16') 698 END DO 699 CALL xml_addNewLine(xf) 700 CALL xml_endElement(u, 'pp_occupations') 701 ! All-electron core charge 702 CALL xml_newElement(u, 'pp_ae_nlcc') 703 CALL xml_AddAttribute(u, 'size' , upf%mesh) 704 DO irow = 1, upf%mesh, 4 705 CALL xml_addNewLine(u) 706 CALL xml_addCharacters(u, upf%paw%ae_rho_atc(irow:min(irow-1+4, upf%mesh)), fmt ='s16') 707 END DO 708 CALL xml_addNewLine(xf) 709 CALL xml_endElement(u, 'pp_ae_nlcc') 710 ! All-electron local potential 711 CALL xml_newElement(u, 'pp_ae_vloc') 712 CALL xml_addAttribute(u, 'size', upf%mesh) 713 DO irow = 1, upf%mesh, 4 714 CALL xml_addNewLine(u) 715 CALL xml_addCharacters(u, upf%paw%ae_vloc(irow:min(irow-1+4, upf%mesh)), fmt ='s16') 716 END DO 717 CALL xml_addNewLine(xf) 718 CALL xml_endElement(u, 'pp_ae_vloc') 719 ! 720 CALL xml_endElement(u, 'pp_paw') 721 ! 722 RETURN 723 END SUBROUTINE write_paw 724 ! 725 SUBROUTINE write_gipaw(u, upf) 726 IMPLICIT NONE 727 TYPE(xmlf_t),INTENT(INOUT) :: u ! i/o unit descriptor 728 TYPE(pseudo_upf),INTENT(IN) :: upf ! the pseudo data 729 INTEGER :: ierr ! /= 0 if something went wrong 730 ! 731 INTEGER :: nb 732 IF (.not. upf%has_gipaw) RETURN 733 734 CALL xml_newElement(u, 'pp_gipaw') 735 CALL xml_addAttribute(u, 'gipaw_data_format', upf%gipaw_data_format ) 736 CALL xml_NewElement (u, 'number_of_core_orbitals') 737 CALL xml_AddCharacters(u, upf%gipaw_ncore_orbitals) 738 CALL xml_endElement(u, 'number_of_core_orbitals') 739 IF ( .NOT. upf%paw_as_gipaw ) THEN 740 CALL xml_newElement (u, 'number_of_valence_orbitals') 741 CALL xml_AddCharacters(u, upf%gipaw_wfs_nchannels) 742 CALL xml_endElement(u, 'number_of_valence_orbitals') 743 END IF 744 DO nb = 1,upf%gipaw_ncore_orbitals 745 CALL xml_NewElement(u, 'pp_gipaw_core_orbital') 746 CALL xml_addAttribute(u, 'size', upf%mesh) 747 CALL xml_addAttribute(u, 'index', nb ) 748 CALL xml_addAttribute(u, 'label', upf%gipaw_core_orbital_el(nb)) 749 CALL xml_addAttribute(u, 'n', upf%gipaw_core_orbital_n(nb)) 750 CALL xml_addAttribute(u, 'l', upf%gipaw_core_orbital_l(nb)) 751 DO irow = 1, upf%mesh, 4 752 CALL xml_addNewLine(u) 753 CALL xml_addCharacters(u, upf%gipaw_core_orbital(irow:min(irow-1+4,upf%mesh), nb), fmt = 's16') 754 END DO 755 CALL xml_addNewLine(xf) 756 CALL xml_EndElement(u, 'pp_gipaw_core_orbital') 757 ENDDO 758 ! 759 ! Only write core orbitals in the PAW as GIPAW case 760 IF (upf%paw_as_gipaw) THEN 761 CALL xml_EndElement(u, 'pp_gipaw') 762 RETURN 763 ENDIF 764 ! 765 ! Write valence all-electron and pseudo orbitals 766 ! 767 DO nb = 1,upf%gipaw_wfs_nchannels 768 CALL xml_NewElement(u, 'pp_gipaw_orbital') 769 CALL xml_addAttribute(u, 'index', nb ) 770 CALL xml_addAttribute(u, 'label', upf%gipaw_wfs_el(nb)) 771 CALL xml_addAttribute(u, 'l', upf%gipaw_wfs_ll(nb)) 772 CALL xml_addAttribute(u, 'cutoff_radius', upf%gipaw_wfs_rcut(nb), fmt = 's16') 773 CALL xml_addAttribute(u, 'ultrasoft_cutoff_radius', upf%gipaw_wfs_rcutus(nb), fmt = 's16') 774 ! 775 CALL xml_NewElement(u, 'pp_gipaw_wfs_ae') 776 CALL xml_addAttribute(u, 'size', upf%mesh) 777 DO irow = 1, upf%mesh, 4 778 CALL xml_addNewLine(u) 779 CALL xml_addCharacters(u, upf%gipaw_wfs_ae(irow:min(irow-1+4, upf%mesh), nb), fmt = 's16') 780 END DO 781 CALL xml_addNewLine(xf) 782 CALL xml_endElement (u, 'pp_gipaw_wfs_ae') 783 CALL xml_NewElement(u, 'pp_gipaw_wfs_ps') 784 CALL xml_addAttribute(u, 'size', upf%mesh) 785 DO irow = 1, upf%mesh, 4 786 CALL xml_addNewLine(u) 787 CALL xml_addCharacters(u, upf%gipaw_wfs_ps(irow:min(irow-1+4, upf%mesh), nb), fmt = 's16') 788 END DO 789 CALL xml_addNewLine(xf) 790 CALL xml_endElement (u, 'pp_gipaw_wfs_ps') 791 ! 792 CALL xml_endElement(u, 'pp_gipaw_orbital' ) 793 ENDDO 794 ! 795 ! Write all-electron and pseudo local potentials 796 CALL xml_NewElement(u, 'pp_gipaw_vlocal') 797 CALL xml_NewElement(u, 'pp_gipaw_vlocal_ae') 798 CALL xml_addAttribute(u, 'size', upf%mesh) 799 DO irow = 1, upf%mesh, 4 800 CALL xml_addNewLine(u) 801 CALL xml_addCharacters(u, upf%gipaw_vlocal_ae(irow:min(irow-1+4, upf%mesh)), fmt = 's16') 802 END DO 803 CALL xml_addNewLine(xf) 804 CALL xml_endElement(u, 'pp_gipaw_vlocal_ae') 805 CALL xml_NewElement(u, 'pp_gipaw_vlocal_ps') 806 CALL xml_addAttribute(u, 'size', upf%mesh) 807 DO irow = 1, upf%mesh, 4 808 CALL xml_addNewLine(u) 809 CALL xml_addCharacters(u, upf%gipaw_vlocal_ps(irow:min(irow-1+4, upf%mesh)), fmt = 's16') 810 END DO 811 CALL xml_addNewLine(xf) 812 CALL xml_endElement(u, 'pp_gipaw_vlocal_ps') 813 CALL xml_endElement(u, 'pp_gipaw_vlocal') 814 ! 815 CALL xml_endElement(u, 'pp_gipaw') 816 817 RETURN 818 END SUBROUTINE write_gipaw 819 ! 820END SUBROUTINE write_upf_schema 821#endif 822 823END MODULE write_upf_schema_module 824 825