TABLE OF CONTENTS
- ABINIT/m_pawxmlps
- m_pawxmlps/atom_t
- m_pawxmlps/gaussian_expansion_t
- m_pawxmlps/generator_t
- m_pawxmlps/paw_begin_element1
- m_pawxmlps/paw_end_element1
- m_pawxmlps/paw_rdfromline
- m_pawxmlps/paw_setup_copy
- m_pawxmlps/paw_setup_free
- m_pawxmlps/paw_setup_t
- m_pawxmlps/pawdata_chunk
- m_pawxmlps/radial_func_t
- m_pawxmlps/radial_grid_t
- m_pawxmlps/rdpawpsxml
- m_pawxmlps/rdpawpsxml_core
- m_pawxmlps/rdpawpsxml_header
- m_pawxmlps/shape_function_t
- m_pawxmlps/state_t
- m_pawxmlps/valence_states_t
- m_pawxmlps/xc_functional_t
ABINIT/m_pawxmlps [ Modules ]
NAME
m_pawxmlps
FUNCTION
This module reads a PAW pseudopotential file written in XML. Can use either FoX or pure Fortran routines.
COPYRIGHT
Copyright (C) 2005-2018 ABINIT group (MT, FJ) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
NOTES
FOR DEVELOPPERS: in order to preserve the portability of libPAW library, please consult ~abinit/src/??_libpaw/libpaw-coding-rules.txt
SOURCE
22 #include "libpaw.h" 23 24 module m_pawxmlps 25 26 USE_DEFS 27 USE_MSG_HANDLING 28 USE_MEMORY_PROFILING 29 30 #if defined LIBPAW_HAVE_FOX 31 use fox_sax 32 #endif 33 34 use m_pawrad , only : pawrad_type, pawrad_init, pawrad_free, bound_deriv 35 use m_paw_numeric, only : paw_spline, paw_splint 36 37 implicit none 38 39 private 40 41 !Procedures used for the Fortran reader 42 public :: rdpawpsxml 43 public :: rdpawpsxml_header 44 public :: rdpawpsxml_core 45 46 !Procedures used for the FoX reader (called from xml_parser in response to events) 47 #if defined LIBPAW_HAVE_FOX 48 public :: paw_begin_element1 49 public :: paw_end_element1 50 public :: pawdata_chunk 51 #endif 52 53 ! private procedures and global variables. 54 private :: paw_rdfromline 55 56 ! This type of real is used by the datatypes below 57 integer,parameter,private :: dpxml = selected_real_kind(14) 58 59 ! The maximum length of a record in a file connected for sequential access. 60 integer,parameter,private :: XML_RECL = 50000 61 integer,parameter,private :: NGAUSSIAN_MAX = 100
m_pawxmlps/atom_t [ Types ]
[ Top ] [ m_pawxmlps ] [ Types ]
NAME
atom_t
FUNCTION
Atom type for the FoX XML reader
SOURCE
233 type, public :: atom_t 234 logical :: tread=.false. 235 character(len=2) :: symbol 236 real(dpxml) :: znucl 237 real(dpxml) :: zion 238 real(dpxml) :: zval 239 end type atom_t
m_pawxmlps/gaussian_expansion_t [ Types ]
[ Top ] [ m_pawxmlps ] [ Types ]
NAME
radial_func_t
FUNCTION
Radial function type for the FoX XML reader
SOURCE
116 type, public :: gaussian_expansion_t 117 logical :: tread=.false. 118 integer :: ngauss=0 119 character(len=6) :: state=' ' 120 real(dpxml), dimension(2, NGAUSSIAN_MAX) :: factors 121 real(dpxml), dimension(2, NGAUSSIAN_MAX) :: expos 122 end type gaussian_expansion_t
m_pawxmlps/generator_t [ Types ]
[ Top ] [ m_pawxmlps ] [ Types ]
NAME
generator_t
FUNCTION
Generator type for the FoX XML reader
SOURCE
197 type, public :: generator_t 198 logical :: tread=.false. 199 character(len=20) :: gen 200 character(len=20) :: name 201 end type generator_t
m_pawxmlps/paw_begin_element1 [ Functions ]
[ Top ] [ m_pawxmlps ] [ Functions ]
NAME
begin_element
FUNCTION
Read an XML tag with a given name. Fills the present module private data.
INPUTS
namespaceURI = universal resource indicator for XML namespace??? Not used. localName = local equivalent of tag name?? Not used. name = name of XML tag which has been read in attributes = attributes of XML tag
OUTPUT
Fills private data in present module.
PARENTS
CHILDREN
SOURCE
348 subroutine paw_begin_element1(namespaceURI,localName,name,attributes) 349 350 351 !This section has been created automatically by the script Abilint (TD). 352 !Do not modify the following lines by hand. 353 #undef ABI_FUNC 354 #define ABI_FUNC 'paw_begin_element1' 355 !End of the abilint section 356 357 character(len=*),intent(in) :: namespaceURI,localName,name 358 type(dictionary_t),intent(in) :: attributes 359 360 character(len=100) :: msg,value 361 integer ::iaewf=0,iproj=0,ipswf=0,iprojfit=0,igauss=0 362 !Just to fool abirules 363 value=localName 364 value=namespaceURI 365 366 select case(name) 367 368 case ("paw_setup") 369 paw_setuploc%tread=.true. 370 igrid=0;ishpf=0 371 paw_setuploc%rpaw=-1.d0 372 LIBPAW_DATATYPE_ALLOCATE(grids,(10)) 373 LIBPAW_DATATYPE_ALLOCATE(shpf,(7)) 374 value = getValue(attributes,"version") 375 write(std_out,'(3a)') "Processing a PSEUDO version ",trim(value)," XML file" 376 paw_setuploc%version=trim(value) 377 378 379 case ("atom") 380 paw_setuploc%atom%tread=.true. 381 value = getValue(attributes,"symbol") 382 if (value == "" ) then 383 msg="Cannot determine atomic symbol" 384 MSG_ERROR(msg) 385 end if 386 paw_setuploc%atom%symbol = trim(value) 387 388 value = getValue(attributes,"Z") 389 if (value == "" ) then 390 msg="Cannot determine znucl" 391 MSG_ERROR(msg) 392 end if 393 read(unit=value,fmt=*) paw_setuploc%atom%znucl 394 395 value = getValue(attributes,"core") 396 if (value == "" ) then 397 msg="Cannot determine zion" 398 MSG_ERROR(msg) 399 end if 400 read(unit=value,fmt=*) paw_setuploc%atom%zion 401 402 value = getValue(attributes,"valence") 403 if (value == "" ) then 404 msg="Cannot determine zval" 405 MSG_ERROR(msg) 406 end if 407 read(unit=value,fmt=*) paw_setuploc%atom%zval 408 409 410 case ("xc_functional") 411 paw_setuploc%xc_functional%tread=.true. 412 value = getValue(attributes,"type") 413 if (value == "" ) then 414 msg="Cannot determine xc-functional-type" 415 MSG_ERROR(msg) 416 end if 417 paw_setuploc%xc_functional%functionaltype = trim(value) 418 419 value = getValue(attributes,"name") 420 if (value == "" ) then 421 msg="Cannot determine xc-functional-name " 422 MSG_ERROR(msg) 423 end if 424 paw_setuploc%xc_functional%name= trim(value) 425 426 case ("generator") 427 paw_setuploc%generator%tread=.true. 428 in_generator =.true. 429 value = getValue(attributes,"type") 430 if (value == "" ) value = "unknown" 431 paw_setuploc%generator%gen = trim(value) 432 433 value = getValue(attributes,"name") 434 if (value == "" ) value = "unknown" 435 paw_setuploc%generator%name = trim(value) 436 437 case ("PAW_radius") 438 value = getValue(attributes,"rpaw") 439 if (value == "" ) then 440 msg="Cannot determine rpaw" 441 MSG_ERROR(msg) 442 end if 443 read(unit=value,fmt=*) paw_setuploc%rpaw 444 445 case ("valence_states") 446 paw_setuploc%valence_states%tread=.true. 447 in_valenceStates=.true. 448 ival=0 449 lmax=0 450 LIBPAW_DATATYPE_ALLOCATE(valstate,(50)) 451 452 case ("state") 453 ival=ival+1 454 455 value = getValue(attributes,"n") 456 if (value == "" ) then 457 valstate(ival)%nn=-1 458 else 459 read(unit=value,fmt=*) valstate(ival)%nn 460 end if 461 462 value = getValue(attributes,"l") 463 if (value == "" ) then 464 msg="Cannot determine l" 465 MSG_ERROR(msg) 466 end if 467 read(unit=value,fmt=*) valstate(ival)%ll 468 if(valstate(ival)%ll>lmax) lmax=valstate(ival)%ll 469 470 value = getValue(attributes,"f") 471 if (value == "" ) then 472 valstate(ival)%ff=-1.d0 473 else 474 read(unit=value,fmt=*) valstate(ival)%ff 475 end if 476 477 value = getValue(attributes,"rc") 478 if (value == "" ) then 479 msg="Cannot determine rc" 480 MSG_ERROR(msg) 481 end if 482 read(unit=value,fmt=*) valstate(ival)%rc 483 484 value = getValue(attributes,"e") 485 if (value == "" ) then 486 msg="Cannot determine e" 487 MSG_ERROR(msg) 488 end if 489 read(unit=value,fmt=*) valstate(ival)%ee 490 491 value = getValue(attributes,"id") 492 if (value == "" ) value = "unknown" 493 valstate(ival)%id = trim(value) 494 495 case ("radial_grid") 496 igrid=igrid+1 497 value = getValue(attributes,"eq") 498 if (value == "" ) value = "unknown" 499 grids(igrid)%eq = trim(value) 500 501 value = getValue(attributes,"a") 502 if (value == "" ) then 503 grids(igrid)%aa=0.d0 504 else 505 read(unit=value,fmt=*) grids(igrid)%aa 506 end if 507 508 value = getValue(attributes,"n") 509 if (value == "" ) then 510 grids(igrid)%nn=0 511 else 512 read(unit=value,fmt=*) grids(igrid)%nn 513 end if 514 515 value = getValue(attributes,"d") 516 if (value == "" ) then 517 grids(igrid)%dd=0.d0 518 else 519 read(unit=value,fmt=*) grids(igrid)%dd 520 end if 521 522 value = getValue(attributes,"b") 523 if (value == "" ) then 524 grids(igrid)%bb=0.d0 525 else 526 read(unit=value,fmt=*) grids(igrid)%bb 527 end if 528 529 value = getValue(attributes,"istart") 530 if (value == "" ) then 531 msg="Cannot determine istart" 532 MSG_ERROR(msg) 533 end if 534 read(unit=value,fmt=*) grids(igrid)%istart 535 536 value = getValue(attributes,"iend") 537 if (value == "" ) then 538 msg="Cannot determine iend" 539 MSG_ERROR(msg) 540 end if 541 read(unit=value,fmt=*) grids(igrid)%iend 542 543 value = getValue(attributes,"id") 544 if (value == "" ) value = "unknown" 545 grids(igrid)%id = trim(value) 546 547 end select 548 549 select case(name) 550 case ("shape_function") 551 paw_setuploc%shape_function%tread=.true. 552 value = getValue(attributes,"type") 553 if (value == "" ) value = "unknown" 554 paw_setuploc%shape_function%gtype = trim(value) 555 556 value = getValue(attributes,"grid") 557 paw_setuploc%shape_function%grid=trim(value) 558 if (value /= "" ) then 559 paw_setuploc%shape_function%gtype ="num" 560 do ii=1,igrid 561 if(trim(paw_setuploc%shape_function%grid)==trim(grids(ii)%id)) then 562 mesh_size=grids(ii)%iend-grids(ii)%istart+1 563 end if 564 end do 565 ishpf=ishpf+1 566 LIBPAW_ALLOCATE(shpf(ishpf)%data,(mesh_size)) 567 rp=>shpf(ishpf) 568 in_data=.true. 569 ndata = 0 570 end if 571 572 value = getValue(attributes,"rc") 573 if (value == "" ) then 574 if(paw_setuploc%shape_function%gtype /="num") then 575 msg="Cannot determine rc" 576 MSG_ERROR(msg) 577 end if 578 else 579 read(unit=value,fmt=*) paw_setuploc%shape_function%rc 580 end if 581 582 value = getValue(attributes,"lamb") 583 if (value == "" ) then 584 paw_setuploc%shape_function%lamb=0 585 else 586 read(unit=value,fmt=*) paw_setuploc%shape_function%lamb 587 end if 588 589 case ("pseudo_partial_wave") 590 ipswf=ipswf+1 591 paw_setuploc%pseudo_partial_wave(ipswf)%tread=.true. 592 value = getValue(attributes,"grid") 593 if (value == "" ) value = "unknown" 594 paw_setuploc%idgrid = trim(value) 595 paw_setuploc%pseudo_partial_wave(ipswf)%grid=trim(value) 596 597 value = getValue(attributes,"state") 598 if (value == "" ) then 599 msg="Cannot determine pseudo_partial_wave state" 600 MSG_ERROR(msg) 601 end if 602 paw_setuploc%pseudo_partial_wave(ipswf)%state=trim(value) 603 604 do ii=1,igrid 605 if(trim(paw_setuploc%pseudo_partial_wave(ipswf)%grid)==trim(grids(ii)%id)) then 606 mesh_size=grids(ii)%iend-grids(ii)%istart+1 607 end if 608 end do 609 610 LIBPAW_ALLOCATE(paw_setuploc%pseudo_partial_wave(ipswf)%data,(mesh_size)) 611 rp=>paw_setuploc%pseudo_partial_wave(ipswf) 612 if(ipswf==paw_setuploc%valence_states%nval) ipswf=0 613 in_data=.true. 614 ndata = 0 615 616 case ("ae_partial_wave") 617 iaewf=iaewf+1 618 paw_setuploc%ae_partial_wave(iaewf)%tread=.true. 619 value = getValue(attributes,"grid") 620 if (value == "" ) value = "unknown" 621 paw_setuploc%ae_partial_wave(iaewf)%grid=trim(value) 622 623 value = getValue(attributes,"state") 624 if (value == "" ) then 625 msg="Cannot determine ae_partial_wave state" 626 MSG_ERROR(msg) 627 end if 628 paw_setuploc%ae_partial_wave(iaewf)%state=trim(value) 629 630 do ii=1,igrid 631 if(trim(paw_setuploc%ae_partial_wave(iaewf)%grid)==trim(grids(ii)%id)) then 632 mesh_size=grids(ii)%iend-grids(ii)%istart+1 633 end if 634 end do 635 636 LIBPAW_ALLOCATE(paw_setuploc%ae_partial_wave(iaewf)%data,(mesh_size)) 637 rp=>paw_setuploc%ae_partial_wave(iaewf) 638 if(iaewf==paw_setuploc%valence_states%nval) iaewf=0 639 in_data=.true. 640 ndata = 0 641 642 case ("projector_function") 643 iproj=iproj+1 644 paw_setuploc%projector_function(iproj)%tread=.true. 645 value = getValue(attributes,"grid") 646 if (value == "" ) value = "unknown" 647 paw_setuploc%projector_function(iproj)%grid=trim(value) 648 649 value = getValue(attributes,"state") 650 if (value == "" ) then 651 msg="Cannot determine projector_function state" 652 MSG_ERROR(msg) 653 end if 654 paw_setuploc%projector_function(iproj)%state=trim(value) 655 656 do ii=1,igrid 657 if(trim(paw_setuploc%projector_function(iproj)%grid)==trim(grids(ii)%id)) then 658 mesh_size=grids(ii)%iend-grids(ii)%istart+1 659 end if 660 end do 661 662 LIBPAW_ALLOCATE(paw_setuploc%projector_function(iproj)%data,(mesh_size)) 663 rp=>paw_setuploc%projector_function(iproj) 664 if(iproj==paw_setuploc%valence_states%nval) iproj=0 665 in_data=.true. 666 ndata = 0 667 668 case ("projector_fit") 669 if(.not.allocated(paw_setuploc%projector_fit)) then 670 LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%projector_fit,(paw_setuploc%valence_states%nval)) 671 end if 672 673 iprojfit=iprojfit+1 674 paw_setuploc%projector_fit(iprojfit)%tread=.true. 675 value = getValue(attributes,"state") 676 if (value == "" ) then 677 msg="Cannot determine projector_fit state" 678 MSG_ERROR(msg) 679 end if 680 paw_setuploc%projector_fit(iprojfit)%state=trim(value) 681 682 if(iprojfit==paw_setuploc%valence_states%nval) iprojfit=0 683 igauss = 0 684 685 case ("gaussian") 686 igauss = igauss + 1 687 value = getValue(attributes,"factor") 688 if (value == "" ) then 689 msg="Cannot determine gaussian factor" 690 MSG_ERROR(msg) 691 end if 692 read(value(2:100), *) paw_setuploc%projector_fit(iprojfit)%factors(1, igauss) 693 read(value(index(value, ',') + 1:100), *) paw_setuploc%projector_fit(iprojfit)%factors(2, igauss) 694 value = getValue(attributes,"exponent") 695 if (value == "" ) then 696 msg="Cannot determine gaussian exponent" 697 MSG_ERROR(msg) 698 end if 699 read(value(2:100), *) paw_setuploc%projector_fit(iprojfit)%expos(1, igauss) 700 read(value(index(value, ',') + 1:100), *) paw_setuploc%projector_fit(iprojfit)%expos(2, igauss) 701 702 case ("ae_core_density") 703 paw_setuploc%ae_core_density%tread=.true. 704 value = getValue(attributes,"grid") 705 if (value == "" ) value = "unknown" 706 paw_setuploc%ae_core_density%grid=trim(value) 707 708 value = getValue(attributes,"state") 709 if (value == "" ) value = "unknown" 710 paw_setuploc%ae_core_density%state=trim(value) 711 712 do ii=1,igrid 713 if(trim(paw_setuploc%ae_core_density%grid)==trim(grids(ii)%id)) then 714 mesh_size=grids(ii)%iend-grids(ii)%istart+1 715 end if 716 end do 717 718 LIBPAW_ALLOCATE(paw_setuploc%ae_core_density%data,(mesh_size)) 719 rp=>paw_setuploc%ae_core_density 720 in_data=.true. 721 ndata = 0 722 723 case ("pseudo_core_density") 724 paw_setuploc%pseudo_core_density%tread=.true. 725 value = getValue(attributes,"grid") 726 if (value == "" ) value = "unknown" 727 paw_setuploc%pseudo_core_density%grid=trim(value) 728 729 value = getValue(attributes,"state") 730 if (value == "" ) value = "unknown" 731 paw_setuploc%pseudo_core_density%state=trim(value) 732 733 do ii=1,igrid 734 if(trim(paw_setuploc%pseudo_core_density%grid)==trim(grids(ii)%id)) then 735 mesh_size=grids(ii)%iend-grids(ii)%istart+1 736 end if 737 end do 738 739 LIBPAW_ALLOCATE(paw_setuploc%pseudo_core_density%data,(mesh_size)) 740 rp=>paw_setuploc%pseudo_core_density 741 in_data=.true. 742 ndata = 0 743 744 case ("pseudo_valence_density") 745 paw_setuploc%pseudo_valence_density%tread=.true. 746 value = getValue(attributes,"grid") 747 if (value == "" ) value = "unknown" 748 paw_setuploc%pseudo_valence_density%grid=trim(value) 749 750 value = getValue(attributes,"state") 751 if (value == "" ) value = "unknown" 752 paw_setuploc%pseudo_valence_density%state=trim(value) 753 754 do ii=1,igrid 755 if(trim(paw_setuploc%pseudo_valence_density%grid)==trim(grids(ii)%id)) then 756 mesh_size=grids(ii)%iend-grids(ii)%istart+1 757 end if 758 end do 759 760 LIBPAW_ALLOCATE(paw_setuploc%pseudo_valence_density%data,(mesh_size)) 761 rp=>paw_setuploc%pseudo_valence_density 762 in_data=.true. 763 ndata = 0 764 765 case ("zero_potential") 766 paw_setuploc%zero_potential%tread=.true. 767 value = getValue(attributes,"grid") 768 if (value == "" ) value = "unknown" 769 paw_setuploc%zero_potential%grid=trim(value) 770 771 value = getValue(attributes,"state") 772 if (value == "" ) value = "unknown" 773 paw_setuploc%zero_potential%state=trim(value) 774 775 do ii=1,igrid 776 if(trim(paw_setuploc%zero_potential%grid)==trim(grids(ii)%id)) then 777 mesh_size=grids(ii)%iend-grids(ii)%istart+1 778 end if 779 end do 780 781 LIBPAW_ALLOCATE(paw_setuploc%zero_potential%data,(mesh_size)) 782 rp=>paw_setuploc%zero_potential 783 in_data=.true. 784 ndata = 0 785 786 case ("LDA_minus_half_potential") 787 paw_setuploc%LDA_minus_half_potential%tread=.true. 788 value = getValue(attributes,"grid") 789 if (value == "" ) value = "unknown" 790 paw_setuploc%LDA_minus_half_potential%grid=trim(value) 791 792 value = getValue(attributes,"state") 793 if (value == "" ) value = "unknown" 794 paw_setuploc%LDA_minus_half_potential%state=trim(value) 795 796 do ii=1,igrid 797 if(trim(paw_setuploc%LDA_minus_half_potential%grid)==trim(grids(ii)%id)) then 798 mesh_size=grids(ii)%iend-grids(ii)%istart+1 799 end if 800 end do 801 802 LIBPAW_ALLOCATE(paw_setuploc%LDA_minus_half_potential%data,(mesh_size)) 803 rp=>paw_setuploc%LDA_minus_half_potential 804 in_data=.true. 805 ndata = 0 806 807 case ("ae_core_kinetic_energy_density") 808 paw_setuploc%ae_core_kinetic_energy_density%tread=.true. 809 value = getValue(attributes,"grid") 810 if (value == "" ) value = "unknown" 811 paw_setuploc%ae_core_kinetic_energy_density%grid=trim(value) 812 813 value = getValue(attributes,"state") 814 if (value == "" ) value = "unknown" 815 paw_setuploc%ae_core_kinetic_energy_density%state=trim(value) 816 817 do ii=1,igrid 818 if(trim(paw_setuploc%ae_core_kinetic_energy_density%grid)==trim(grids(ii)%id)) then 819 mesh_size=grids(ii)%iend-grids(ii)%istart+1 820 end if 821 end do 822 823 LIBPAW_ALLOCATE(paw_setuploc%ae_core_kinetic_energy_density%data,(mesh_size)) 824 rp=>paw_setuploc%ae_core_kinetic_energy_density 825 in_data=.true. 826 ndata = 0 827 828 case ("pseudo_core_kinetic_energy_density") 829 paw_setuploc%pseudo_core_kinetic_energy_density%tread=.true. 830 value = getValue(attributes,"grid") 831 if (value == "" ) value = "unknown" 832 paw_setuploc%pseudo_core_kinetic_energy_density%grid=trim(value) 833 834 value = getValue(attributes,"state") 835 if (value == "" ) value = "unknown" 836 paw_setuploc%pseudo_core_kinetic_energy_density%state=trim(value) 837 838 do ii=1,igrid 839 if(trim(paw_setuploc%pseudo_core_kinetic_energy_density%grid)==trim(grids(ii)%id)) then 840 mesh_size=grids(ii)%iend-grids(ii)%istart+1 841 end if 842 end do 843 844 LIBPAW_ALLOCATE(paw_setuploc%pseudo_core_kinetic_energy_density%data,(mesh_size)) 845 rp=>paw_setuploc%pseudo_core_kinetic_energy_density 846 in_data=.true. 847 ndata = 0 848 849 case ("kresse_joubert_local_ionic_potential") 850 paw_setuploc%kresse_joubert_local_ionic_potential%tread=.true. 851 value = getValue(attributes,"grid") 852 if (value == "" ) value = "unknown" 853 paw_setuploc%kresse_joubert_local_ionic_potential%grid=trim(value) 854 855 value = getValue(attributes,"state") 856 if (value == "" ) value = "unknown" 857 paw_setuploc%kresse_joubert_local_ionic_potential%state=trim(value) 858 859 do ii=1,igrid 860 if(trim(paw_setuploc%kresse_joubert_local_ionic_potential%grid)==trim(grids(ii)%id)) then 861 mesh_size=grids(ii)%iend-grids(ii)%istart+1 862 end if 863 end do 864 865 LIBPAW_ALLOCATE(paw_setuploc%kresse_joubert_local_ionic_potential%data,(mesh_size)) 866 rp=>paw_setuploc%kresse_joubert_local_ionic_potential 867 in_data=.true. 868 ndata = 0 869 870 case ("blochl_local_ionic_potential") 871 paw_setuploc%blochl_local_ionic_potential%tread=.true. 872 value = getValue(attributes,"grid") 873 if (value == "" ) value = "unknown" 874 paw_setuploc%blochl_local_ionic_potential%grid=trim(value) 875 876 value = getValue(attributes,"state") 877 if (value == "" ) value = "unknown" 878 paw_setuploc%blochl_local_ionic_potential%state=trim(value) 879 880 do ii=1,igrid 881 if(trim(paw_setuploc%blochl_local_ionic_potential%grid)==trim(grids(ii)%id)) then 882 mesh_size=grids(ii)%iend-grids(ii)%istart+1 883 end if 884 end do 885 886 LIBPAW_ALLOCATE(paw_setuploc%blochl_local_ionic_potential%data,(mesh_size)) 887 rp=>paw_setuploc%blochl_local_ionic_potential 888 in_data=.true. 889 ndata = 0 890 891 case ("kinetic_energy_differences") 892 paw_setuploc%kinetic_energy_differences%tread=.true. 893 mesh_size=paw_setuploc%valence_states%nval*paw_setuploc%valence_states%nval 894 LIBPAW_ALLOCATE(paw_setuploc%kinetic_energy_differences%data,(mesh_size)) 895 rp=>paw_setuploc%kinetic_energy_differences 896 in_data=.true. 897 ndata = 0 898 899 end select 900 901 end subroutine paw_begin_element1
m_pawxmlps/paw_end_element1 [ Functions ]
[ Top ] [ m_pawxmlps ] [ Functions ]
NAME
end_element
FUNCTION
End XML tag effect: switches flags in private data of this module
INPUTS
namespaceURI = universal resource indicator for XML namespace??? Not used. localName = local equivalent of tag name?? Not used. name = name of XML tag which has been read in
OUTPUT
side effect: private data flags in present module are turned to .false.
PARENTS
CHILDREN
SOURCE
926 subroutine paw_end_element1(namespaceURI,localName,name) 927 928 929 !This section has been created automatically by the script Abilint (TD). 930 !Do not modify the following lines by hand. 931 #undef ABI_FUNC 932 #define ABI_FUNC 'paw_end_element1' 933 !End of the abilint section 934 935 character(len=*),intent(in) :: namespaceURI,localName,name 936 character(len=100) :: msg,value 937 938 !Just to fool abirules 939 value=localName 940 value=namespaceURI 941 942 select case(name) 943 944 case ("generator") 945 in_generator = .false. 946 947 case ("valence_states") 948 in_valenceStates = .false. 949 if(ival>50) then 950 msg="ival>50" 951 MSG_ERROR(msg) 952 end if 953 if(ival>0)then 954 LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%valence_states%state,(ival)) 955 paw_setuploc%valence_states%state(ival)%tread=.true. 956 paw_setuploc%valence_states%nval=ival 957 do ii=1,ival 958 paw_setuploc%valence_states%state(ii)=valstate(ii) 959 end do 960 end if 961 LIBPAW_DATATYPE_DEALLOCATE(valstate) 962 if(.not.allocated(paw_setuploc%ae_partial_wave)) then 963 LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%ae_partial_wave,(paw_setuploc%valence_states%nval)) 964 end if 965 if(.not.allocated(paw_setuploc%pseudo_partial_wave)) then 966 LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%pseudo_partial_wave,(paw_setuploc%valence_states%nval)) 967 end if 968 if(.not.allocated(paw_setuploc%projector_function)) then 969 LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%projector_function,(paw_setuploc%valence_states%nval)) 970 end if 971 972 case ("paw_setup") 973 if(igrid>10) then 974 msg="igrid>10" 975 MSG_ERROR(msg) 976 end if 977 LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%radial_grid,(igrid)) 978 paw_setuploc%radial_grid(igrid)%tread=.true. 979 paw_setuploc%ngrid=igrid 980 do ii=1,igrid 981 paw_setuploc%radial_grid(ii)=grids(ii) 982 end do 983 LIBPAW_DATATYPE_DEALLOCATE(grids) 984 do ii=1,igrid 985 if(trim(paw_setuploc%shape_function%grid)==trim(paw_setuploc%radial_grid(ii)%id)) then 986 mesh_size=paw_setuploc%radial_grid(ii)%iend-paw_setuploc%radial_grid(ii)%istart+1 987 end if 988 end do 989 if(ishpf>10) then 990 msg="ishpf>7" 991 MSG_ERROR(msg) 992 end if 993 LIBPAW_ALLOCATE(paw_setuploc%shape_function%data,(mesh_size,ishpf)) 994 do ii=1,ishpf 995 paw_setuploc%shape_function%data(:,ii)=shpf(ii)%data(:) 996 LIBPAW_DEALLOCATE(shpf(ii)%data) 997 end do 998 LIBPAW_DATATYPE_DEALLOCATE(shpf) 999 1000 case ("shape_function") 1001 in_data=.false. 1002 1003 case ("pseudo_partial_wave") 1004 in_data=.false. 1005 1006 case ("ae_partial_wave") 1007 in_data=.false. 1008 1009 case ("projector_function") 1010 in_data=.false. 1011 1012 case ("ae_core_density") 1013 in_data=.false. 1014 1015 case ("pseudo_core_density") 1016 in_data=.false. 1017 1018 case ("pseudo_valence_density") 1019 in_data=.false. 1020 1021 case ("zero_potential") 1022 in_data=.false. 1023 1024 case ("LDA_minus_half_potential") 1025 in_data=.false. 1026 1027 case ("ae_core_kinetic_energy_density") 1028 in_data=.false. 1029 1030 case ("pseudo_core_kinetic_energy_density") 1031 in_data=.false. 1032 1033 case ("kresse_joubert_local_ionic_potential") 1034 in_data=.false. 1035 1036 case ("blochl_local_ionic_potential") 1037 in_data=.false. 1038 1039 case ("kinetic_energy_differences") 1040 in_data=.false. 1041 1042 end select 1043 1044 end subroutine paw_end_element1
m_pawxmlps/paw_rdfromline [ Functions ]
[ Top ] [ m_pawxmlps ] [ Functions ]
NAME
paw_rdfromline
FUNCTION
Read the value of a keyword from a XML line
INPUTS
keyword= keyword which value has to be read line= string from which the data are read (line from a XML)
OUTPUT
ierr= error code output= (string) value of the keyword
PARENTS
m_pawxmlps
CHILDREN
paw_rdfromline
SOURCE
1541 subroutine paw_rdfromline(keyword,line,output,ierr) 1542 1543 1544 !This section has been created automatically by the script Abilint (TD). 1545 !Do not modify the following lines by hand. 1546 #undef ABI_FUNC 1547 #define ABI_FUNC 'paw_rdfromline' 1548 !End of the abilint section 1549 1550 implicit none 1551 1552 !Arguments --------------------------------------------- 1553 character(len=*), intent(in) :: keyword,line 1554 character(len=*), intent(out) :: output 1555 integer, intent(out) :: ierr 1556 !Local variables --------------------------------------- 1557 character(len=len(line)) :: temp 1558 integer :: pos,pos2 1559 1560 ! ********************************************************************* 1561 1562 ierr=1;output="" 1563 pos=index(line,trim(keyword)) 1564 if (pos>0) then 1565 temp=line(pos+len_trim(keyword):len_trim(line)) 1566 pos=index(temp,char(34)) 1567 if (pos>0) then 1568 pos2=index(temp(pos+1:len_trim(temp)),char(34)) 1569 if (pos2>0) then 1570 output=temp(pos+1:pos+pos2-1) 1571 end if 1572 end if 1573 end if 1574 1575 end subroutine paw_rdfromline
m_pawxmlps/paw_setup_copy [ Functions ]
[ Top ] [ m_pawxmlps ] [ Functions ]
NAME
paw_setup_copy
FUNCTION
Copy a paw_setup datastructure into another
INPUTS
paw_setupin<paw_setup_type>=input paw_setup datastructure
OUTPUT
paw_setupout<paw_setup_type>=output paw_setup datastructure
PARENTS
inpspheads
CHILDREN
paw_rdfromline
SOURCE
1296 subroutine paw_setup_copy(paw_setupin,paw_setupout) 1297 1298 1299 !This section has been created automatically by the script Abilint (TD). 1300 !Do not modify the following lines by hand. 1301 #undef ABI_FUNC 1302 #define ABI_FUNC 'paw_setup_copy' 1303 !End of the abilint section 1304 1305 implicit none 1306 1307 !Arguments ------------------------------------ 1308 !scalars 1309 type(paw_setup_t),intent(in) :: paw_setupin 1310 type(paw_setup_t),intent(out) :: paw_setupout 1311 1312 !Local variables------------------------------- 1313 !scalars 1314 integer :: ii,sz1,sz2 1315 1316 ! ********************************************************************* 1317 1318 !scalars 1319 paw_setupout%version=paw_setupin%version 1320 paw_setupout%tread=paw_setupin%tread 1321 paw_setupout%ngrid=paw_setupin%ngrid 1322 paw_setupout%idgrid=paw_setupin%idgrid 1323 paw_setupout%rpaw=paw_setupin%rpaw 1324 paw_setupout%ex_cc=paw_setupin%ex_cc 1325 paw_setupout%atom%tread=paw_setupin%atom%tread 1326 paw_setupout%atom%symbol=paw_setupin%atom%symbol 1327 paw_setupout%atom%znucl=paw_setupin%atom%znucl 1328 paw_setupout%atom%zion=paw_setupin%atom%zion 1329 paw_setupout%atom%zval=paw_setupin%atom%zval 1330 paw_setupout%xc_functional%tread=paw_setupin%xc_functional%tread 1331 paw_setupout%xc_functional%functionaltype=paw_setupin%xc_functional%functionaltype 1332 paw_setupout%xc_functional%name=paw_setupin%xc_functional%name 1333 paw_setupout%generator%tread=paw_setupin%generator%tread 1334 paw_setupout%generator%gen=paw_setupin%generator%gen 1335 paw_setupout%generator%name=paw_setupin%generator%name 1336 paw_setupout%valence_states%tread=paw_setupin%valence_states%tread 1337 paw_setupout%valence_states%nval=paw_setupin%valence_states%nval 1338 paw_setupout%shape_function%tread=paw_setupin%shape_function%tread 1339 paw_setupout%shape_function%gtype=paw_setupin%shape_function%gtype 1340 paw_setupout%shape_function%grid=paw_setupin%shape_function%grid 1341 paw_setupout%shape_function%rc=paw_setupin%shape_function%rc 1342 paw_setupout%shape_function%lamb=paw_setupin%shape_function%lamb 1343 paw_setupout%ae_core_density%tread=paw_setupin%ae_core_density%tread 1344 paw_setupout%ae_core_density%grid=paw_setupin%ae_core_density%grid 1345 paw_setupout%ae_core_density%state=paw_setupin%ae_core_density%state 1346 paw_setupout%pseudo_core_density%tread=paw_setupin%pseudo_core_density%tread 1347 paw_setupout%pseudo_core_density%grid=paw_setupin%pseudo_core_density%grid 1348 paw_setupout%pseudo_core_density%state=paw_setupin%pseudo_core_density%state 1349 paw_setupout%pseudo_valence_density%tread=paw_setupin%pseudo_valence_density%tread 1350 paw_setupout%pseudo_valence_density%grid=paw_setupin%pseudo_valence_density%grid 1351 paw_setupout%pseudo_valence_density%state=paw_setupin%pseudo_valence_density%state 1352 paw_setupout%zero_potential%tread=paw_setupin%zero_potential%tread 1353 paw_setupout%zero_potential%grid=paw_setupin%zero_potential%grid 1354 paw_setupout%zero_potential%state=paw_setupin%zero_potential%state 1355 paw_setupout%LDA_minus_half_potential%tread=paw_setupin%LDA_minus_half_potential%tread 1356 paw_setupout%LDA_minus_half_potential%grid=paw_setupin%LDA_minus_half_potential%grid 1357 paw_setupout%LDA_minus_half_potential%state=paw_setupin%LDA_minus_half_potential%state 1358 paw_setupout%ae_core_kinetic_energy_density%tread=& 1359 & paw_setupin%ae_core_kinetic_energy_density%tread 1360 paw_setupout%ae_core_kinetic_energy_density%grid=& 1361 & paw_setupin%ae_core_kinetic_energy_density%grid 1362 paw_setupout%ae_core_kinetic_energy_density%state=& 1363 & paw_setupin%ae_core_kinetic_energy_density%state 1364 paw_setupout%pseudo_core_kinetic_energy_density%tread=& 1365 & paw_setupin%pseudo_core_kinetic_energy_density%tread 1366 paw_setupout%pseudo_core_kinetic_energy_density%grid=& 1367 & paw_setupin%pseudo_core_kinetic_energy_density%grid 1368 paw_setupout%pseudo_core_kinetic_energy_density%state=& 1369 & paw_setupin%pseudo_core_kinetic_energy_density%state 1370 paw_setupout%kresse_joubert_local_ionic_potential%tread=& 1371 & paw_setupin%kresse_joubert_local_ionic_potential%tread 1372 paw_setupout%kresse_joubert_local_ionic_potential%grid=& 1373 & paw_setupin%kresse_joubert_local_ionic_potential%grid 1374 paw_setupout%kresse_joubert_local_ionic_potential%state=& 1375 & paw_setupin%kresse_joubert_local_ionic_potential%state 1376 paw_setupout%blochl_local_ionic_potential%tread=& 1377 & paw_setupin%blochl_local_ionic_potential%tread 1378 paw_setupout%blochl_local_ionic_potential%grid=& 1379 & paw_setupin%blochl_local_ionic_potential%grid 1380 paw_setupout%blochl_local_ionic_potential%state=& 1381 & paw_setupin%blochl_local_ionic_potential%state 1382 paw_setupout%kinetic_energy_differences%tread=paw_setupin%kinetic_energy_differences%tread 1383 paw_setupout%kinetic_energy_differences%grid=paw_setupin%kinetic_energy_differences%grid 1384 paw_setupout%kinetic_energy_differences%state=paw_setupin%kinetic_energy_differences%state 1385 paw_setupout%exact_exchange_matrix%tread=paw_setupin%exact_exchange_matrix%tread 1386 ! allocatable arrays 1387 if (allocated(paw_setupin%shape_function%data)) then 1388 sz1=size(paw_setupin%shape_function%data,1) 1389 sz2=size(paw_setupin%shape_function%data,2) 1390 LIBPAW_ALLOCATE(paw_setupout%shape_function%data,(sz1,sz2)) 1391 paw_setupout%shape_function%data=paw_setupin%shape_function%data 1392 end if 1393 if (allocated(paw_setupin%ae_core_density%data)) then 1394 sz1=size(paw_setupin%ae_core_density%data,1) 1395 LIBPAW_ALLOCATE(paw_setupout%ae_core_density%data,(sz1)) 1396 paw_setupout%ae_core_density%data=paw_setupin%ae_core_density%data 1397 end if 1398 if (allocated(paw_setupin%pseudo_core_density%data)) then 1399 sz1=size(paw_setupin%pseudo_core_density%data,1) 1400 LIBPAW_ALLOCATE(paw_setupout%pseudo_core_density%data,(sz1)) 1401 paw_setupout%pseudo_core_density%data=paw_setupin%pseudo_core_density%data 1402 end if 1403 if (allocated(paw_setupin%pseudo_valence_density%data)) then 1404 sz1=size(paw_setupin%pseudo_valence_density%data,1) 1405 LIBPAW_ALLOCATE(paw_setupout%pseudo_valence_density%data,(sz1)) 1406 paw_setupout%pseudo_valence_density%data=paw_setupin%pseudo_valence_density%data 1407 end if 1408 if (allocated(paw_setupin%zero_potential%data)) then 1409 sz1=size(paw_setupin%zero_potential%data,1) 1410 LIBPAW_ALLOCATE(paw_setupout%zero_potential%data,(sz1)) 1411 paw_setupout%zero_potential%data=paw_setupin%zero_potential%data 1412 end if 1413 if (allocated(paw_setupin%LDA_minus_half_potential%data)) then 1414 sz1=size(paw_setupin%LDA_minus_half_potential%data,1) 1415 LIBPAW_ALLOCATE(paw_setupout%LDA_minus_half_potential%data,(sz1)) 1416 paw_setupout%LDA_minus_half_potential%data=paw_setupin%LDA_minus_half_potential%data 1417 end if 1418 if (allocated(paw_setupin%ae_core_kinetic_energy_density%data)) then 1419 sz1=size(paw_setupin%ae_core_kinetic_energy_density%data,1) 1420 LIBPAW_ALLOCATE(paw_setupout%ae_core_kinetic_energy_density%data,(sz1)) 1421 paw_setupout%ae_core_kinetic_energy_density%data=paw_setupin%ae_core_kinetic_energy_density%data 1422 end if 1423 if (allocated(paw_setupin%pseudo_core_kinetic_energy_density%data)) then 1424 sz1=size(paw_setupin%pseudo_core_kinetic_energy_density%data,1) 1425 LIBPAW_ALLOCATE(paw_setupout%pseudo_core_kinetic_energy_density%data,(sz1)) 1426 paw_setupout%pseudo_core_kinetic_energy_density%data=paw_setupin%pseudo_core_kinetic_energy_density%data 1427 end if 1428 if (allocated(paw_setupin%kresse_joubert_local_ionic_potential%data)) then 1429 sz1=size(paw_setupin%kresse_joubert_local_ionic_potential%data,1) 1430 LIBPAW_ALLOCATE(paw_setupout%kresse_joubert_local_ionic_potential%data,(sz1)) 1431 paw_setupout%kresse_joubert_local_ionic_potential%data=paw_setupin%kresse_joubert_local_ionic_potential%data 1432 end if 1433 if (allocated(paw_setupin%blochl_local_ionic_potential%data)) then 1434 sz1=size(paw_setupin%blochl_local_ionic_potential%data,1) 1435 LIBPAW_ALLOCATE(paw_setupout%blochl_local_ionic_potential%data,(sz1)) 1436 paw_setupout%blochl_local_ionic_potential%data=paw_setupin%blochl_local_ionic_potential%data 1437 end if 1438 if (allocated(paw_setupin%exact_exchange_matrix%data)) then 1439 sz1=size(paw_setupin%exact_exchange_matrix%data,1) 1440 LIBPAW_ALLOCATE(paw_setupout%exact_exchange_matrix%data,(sz1)) 1441 paw_setupout%exact_exchange_matrix%data=paw_setupin%exact_exchange_matrix%data 1442 end if 1443 if (allocated(paw_setupin%kinetic_energy_differences%data)) then 1444 sz1=size(paw_setupin%kinetic_energy_differences%data,1) 1445 LIBPAW_ALLOCATE(paw_setupout%kinetic_energy_differences%data,(sz1)) 1446 paw_setupout%kinetic_energy_differences%data=paw_setupin%kinetic_energy_differences%data 1447 end if 1448 if(allocated( paw_setupin%radial_grid)) then 1449 sz1=size(paw_setupin%radial_grid,1) 1450 LIBPAW_DATATYPE_ALLOCATE(paw_setupout%radial_grid,(sz1)) 1451 paw_setupout%radial_grid=paw_setupin%radial_grid 1452 end if 1453 if(allocated(paw_setupin%valence_states%state)) then 1454 sz1=size(paw_setupin%valence_states%state,1) 1455 LIBPAW_DATATYPE_ALLOCATE(paw_setupout%valence_states%state,(sz1)) 1456 paw_setupout%valence_states%state=paw_setupin%valence_states%state 1457 end if 1458 1459 if (allocated( paw_setupin%ae_partial_wave)) then 1460 sz1=size(paw_setupin%ae_partial_wave,1) 1461 LIBPAW_DATATYPE_ALLOCATE(paw_setupout%ae_partial_wave,(sz1)) 1462 do ii=1,paw_setupin%valence_states%nval 1463 paw_setupout%ae_partial_wave(ii)%tread=paw_setupin%ae_partial_wave(ii)%tread 1464 paw_setupout%ae_partial_wave(ii)%grid=paw_setupin%ae_partial_wave(ii)%grid 1465 paw_setupout%ae_partial_wave(ii)%state=paw_setupin%ae_partial_wave(ii)%state 1466 if(allocated( paw_setupin%ae_partial_wave(ii)%data)) then 1467 sz1=size(paw_setupin%ae_partial_wave(ii)%data,1) 1468 LIBPAW_ALLOCATE(paw_setupout%ae_partial_wave(ii)%data,(sz1)) 1469 paw_setupout%ae_partial_wave(ii)%data=paw_setupin%ae_partial_wave(ii)%data 1470 end if 1471 end do 1472 end if 1473 if (allocated( paw_setupin%pseudo_partial_wave)) then 1474 sz1=size(paw_setupin%pseudo_partial_wave,1) 1475 LIBPAW_DATATYPE_ALLOCATE(paw_setupout%pseudo_partial_wave,(sz1)) 1476 do ii=1,paw_setupin%valence_states%nval 1477 paw_setupout%pseudo_partial_wave(ii)%tread=paw_setupin%pseudo_partial_wave(ii)%tread 1478 paw_setupout%pseudo_partial_wave(ii)%grid=paw_setupin%pseudo_partial_wave(ii)%grid 1479 paw_setupout%pseudo_partial_wave(ii)%state=paw_setupin%pseudo_partial_wave(ii)%state 1480 if(allocated( paw_setupin%pseudo_partial_wave(ii)%data)) then 1481 sz1=size(paw_setupin%pseudo_partial_wave(ii)%data,1) 1482 LIBPAW_ALLOCATE(paw_setupout%pseudo_partial_wave(ii)%data,(sz1)) 1483 paw_setupout%pseudo_partial_wave(ii)%data=paw_setupin%pseudo_partial_wave(ii)%data 1484 end if 1485 end do 1486 end if 1487 if (allocated( paw_setupin%projector_function)) then 1488 sz1=size(paw_setupin%projector_function,1) 1489 LIBPAW_DATATYPE_ALLOCATE(paw_setupout%projector_function,(sz1)) 1490 do ii=1,paw_setupin%valence_states%nval 1491 paw_setupout%projector_function(ii)%tread=paw_setupin%projector_function(ii)%tread 1492 paw_setupout%projector_function(ii)%grid=paw_setupin%projector_function(ii)%grid 1493 paw_setupout%projector_function(ii)%state=paw_setupin%projector_function(ii)%state 1494 if(allocated( paw_setupin%projector_function(ii)%data)) then 1495 sz1=size(paw_setupin%projector_function(ii)%data,1) 1496 LIBPAW_ALLOCATE(paw_setupout%projector_function(ii)%data,(sz1)) 1497 paw_setupout%projector_function(ii)%data=paw_setupin%projector_function(ii)%data 1498 end if 1499 end do 1500 end if 1501 if (allocated( paw_setupin%projector_fit)) then 1502 sz1=size(paw_setupin%projector_fit,1) 1503 LIBPAW_DATATYPE_ALLOCATE(paw_setupout%projector_fit,(sz1)) 1504 do ii=1,paw_setupin%valence_states%nval 1505 paw_setupout%projector_fit(ii)%tread=paw_setupin%projector_fit(ii)%tread 1506 paw_setupout%projector_fit(ii)%ngauss=paw_setupin%projector_fit(ii)%ngauss 1507 paw_setupout%projector_fit(ii)%state=paw_setupin%projector_fit(ii)%state 1508 paw_setupout%projector_fit(ii)%factors=paw_setupin%projector_fit(ii)%factors 1509 paw_setupout%projector_fit(ii)%expos=paw_setupin%projector_fit(ii)%expos 1510 end do 1511 end if 1512 1513 end subroutine paw_setup_copy
m_pawxmlps/paw_setup_free [ Functions ]
[ Top ] [ m_pawxmlps ] [ Functions ]
NAME
paw_setup_free
FUNCTION
Destroy a paw_setup datastructure
SIDE EFFECTS
paw_setup<paw_setup_type>=Datatype gathering information on XML paw setup.
PARENTS
abinit,inpspheads
CHILDREN
paw_rdfromline
SOURCE
1161 subroutine paw_setup_free(paw_setupin) 1162 1163 1164 !This section has been created automatically by the script Abilint (TD). 1165 !Do not modify the following lines by hand. 1166 #undef ABI_FUNC 1167 #define ABI_FUNC 'paw_setup_free' 1168 !End of the abilint section 1169 1170 implicit none 1171 1172 !Arguments ------------------------------------ 1173 !scalars 1174 type(paw_setup_t),intent(inout) :: paw_setupin 1175 1176 !Local variables------------------------------- 1177 integer :: ii 1178 1179 ! ********************************************************************* 1180 1181 paw_setupin%tread=.false. 1182 paw_setupin%atom%tread=.false. 1183 paw_setupin%xc_functional%tread=.false. 1184 paw_setupin%generator%tread=.false. 1185 paw_setupin%valence_states%tread=.false. 1186 paw_setupin%shape_function%tread=.false. 1187 paw_setupin%ae_core_density%tread=.false. 1188 paw_setupin%pseudo_core_density%tread=.false. 1189 paw_setupin%pseudo_valence_density%tread=.false. 1190 paw_setupin%zero_potential%tread=.false. 1191 paw_setupin%LDA_minus_half_potential%tread=.false. 1192 paw_setupin%ae_core_kinetic_energy_density%tread=.false. 1193 paw_setupin%pseudo_core_kinetic_energy_density%tread=.false. 1194 paw_setupin%kresse_joubert_local_ionic_potential%tread=.false. 1195 paw_setupin%blochl_local_ionic_potential%tread=.false. 1196 paw_setupin%kinetic_energy_differences%tread=.false. 1197 paw_setupin%exact_exchange_matrix%tread=.false. 1198 1199 if(allocated( paw_setupin%shape_function%data)) then 1200 LIBPAW_DEALLOCATE(paw_setupin%shape_function%data) 1201 end if 1202 if(allocated( paw_setupin%ae_core_density%data)) then 1203 LIBPAW_DEALLOCATE(paw_setupin%ae_core_density%data) 1204 end if 1205 if(allocated( paw_setupin%pseudo_core_density%data)) then 1206 LIBPAW_DEALLOCATE(paw_setupin%pseudo_core_density%data) 1207 end if 1208 if(allocated( paw_setupin%pseudo_valence_density%data)) then 1209 LIBPAW_DEALLOCATE(paw_setupin%pseudo_valence_density%data) 1210 end if 1211 if(allocated( paw_setupin%zero_potential%data)) then 1212 LIBPAW_DEALLOCATE(paw_setupin%zero_potential%data) 1213 end if 1214 if(allocated( paw_setupin%LDA_minus_half_potential%data)) then 1215 LIBPAW_DEALLOCATE(paw_setupin%LDA_minus_half_potential%data) 1216 end if 1217 if(allocated( paw_setupin%ae_core_kinetic_energy_density%data)) then 1218 LIBPAW_DEALLOCATE(paw_setupin%ae_core_kinetic_energy_density%data) 1219 end if 1220 if(allocated( paw_setupin%pseudo_core_kinetic_energy_density%data)) then 1221 LIBPAW_DEALLOCATE(paw_setupin%pseudo_core_kinetic_energy_density%data) 1222 end if 1223 if(allocated( paw_setupin%kresse_joubert_local_ionic_potential%data)) then 1224 LIBPAW_DEALLOCATE(paw_setupin%kresse_joubert_local_ionic_potential%data) 1225 end if 1226 if(allocated( paw_setupin%blochl_local_ionic_potential%data)) then 1227 LIBPAW_DEALLOCATE(paw_setupin%blochl_local_ionic_potential%data) 1228 end if 1229 if(allocated( paw_setupin%kinetic_energy_differences%data)) then 1230 LIBPAW_DEALLOCATE(paw_setupin%kinetic_energy_differences%data) 1231 end if 1232 if(allocated( paw_setupin%exact_exchange_matrix%data)) then 1233 LIBPAW_DEALLOCATE(paw_setupin%exact_exchange_matrix%data) 1234 end if 1235 if (allocated( paw_setupin%ae_partial_wave)) then 1236 do ii=1,paw_setupin%valence_states%nval 1237 if(allocated( paw_setupin%ae_partial_wave(ii)%data)) then 1238 LIBPAW_DEALLOCATE(paw_setupin%ae_partial_wave(ii)%data) 1239 end if 1240 end do 1241 LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%ae_partial_wave) 1242 end if 1243 if (allocated( paw_setupin%pseudo_partial_wave)) then 1244 do ii=1,paw_setupin%valence_states%nval 1245 if(allocated( paw_setupin%pseudo_partial_wave(ii)%data)) then 1246 LIBPAW_DEALLOCATE(paw_setupin%pseudo_partial_wave(ii)%data) 1247 end if 1248 end do 1249 LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%pseudo_partial_wave) 1250 end if 1251 if (allocated( paw_setupin%projector_function)) then 1252 do ii=1,paw_setupin%valence_states%nval 1253 if(allocated( paw_setupin%projector_function(ii)%data)) then 1254 LIBPAW_DEALLOCATE(paw_setupin%projector_function(ii)%data) 1255 end if 1256 end do 1257 LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%projector_function) 1258 end if 1259 if (allocated( paw_setupin%projector_fit)) then 1260 LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%projector_fit) 1261 end if 1262 if(allocated(paw_setupin%valence_states%state)) then 1263 LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%valence_states%state) 1264 end if 1265 if(allocated( paw_setupin%radial_grid)) then 1266 LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%radial_grid) 1267 end if 1268 1269 end subroutine paw_setup_free
m_pawxmlps/paw_setup_t [ Types ]
[ Top ] [ m_pawxmlps ] [ Types ]
NAME
paw_setup_t
FUNCTION
PAW setup type (contain all the data for a PAW setup)
SOURCE
253 type, public :: paw_setup_t 254 character(len=3) :: version 255 logical :: tread=.false. 256 integer :: ngrid 257 real(dpxml) :: rpaw 258 real(dpxml) :: ex_cc 259 character(len=4) :: idgrid 260 type(atom_t) :: atom 261 type(xc_functional_t) :: xc_functional 262 type(generator_t) :: generator 263 type(valence_states_t) :: valence_states 264 type(radial_grid_t), allocatable :: radial_grid(:) 265 type(shape_function_t) :: shape_function 266 type(radialfunc_t) :: ae_core_density 267 type(radialfunc_t) :: pseudo_core_density 268 type(radialfunc_t) :: pseudo_valence_density 269 type(radialfunc_t) :: zero_potential 270 type(radialfunc_t) :: LDA_minus_half_potential 271 type(radialfunc_t) :: ae_core_kinetic_energy_density 272 type(radialfunc_t) :: pseudo_core_kinetic_energy_density 273 type(radialfunc_t),allocatable :: ae_partial_wave(:) 274 type(radialfunc_t),allocatable :: pseudo_partial_wave(:) 275 type(radialfunc_t),allocatable :: projector_function(:) 276 type(gaussian_expansion_t),allocatable :: projector_fit(:) 277 type(radialfunc_t) :: kresse_joubert_local_ionic_potential 278 type(radialfunc_t) :: blochl_local_ionic_potential 279 type(radialfunc_t) :: kinetic_energy_differences 280 type(radialfunc_t) :: exact_exchange_matrix 281 end type paw_setup_t 282 283 public :: paw_setup_free ! Free memory 284 public :: paw_setup_copy ! Copy object
m_pawxmlps/pawdata_chunk [ Functions ]
[ Top ] [ m_pawxmlps ] [ Functions ]
NAME
pawdata_chunk
FUNCTION
Take a string and turn it into useful data structure (reals)
INPUTS
chunk=raw data for chunk of XML data
OUTPUT
SIDE EFFECTS
Copied and translated into module data (side effect)
PARENTS
CHILDREN
paw_rdfromline
SOURCE
1070 subroutine pawdata_chunk(chunk) 1071 1072 1073 !This section has been created automatically by the script Abilint (TD). 1074 !Do not modify the following lines by hand. 1075 #undef ABI_FUNC 1076 #define ABI_FUNC 'pawdata_chunk' 1077 !End of the abilint section 1078 1079 character(len=*),intent(in) :: chunk 1080 1081 integer :: ii,ntokens,status,last_pos 1082 logical :: in_token 1083 character(len=len(chunk)) :: str 1084 character(len=50) :: msg 1085 character(len=1) :: cc 1086 real(dpxml),pointer :: x(:) 1087 1088 1089 if (len_trim(chunk) == 0) RETURN ! skip empty chunk 1090 1091 if (in_data) then 1092 str = chunk ; x => rp%data 1093 1094 ! Check the contents of the string and find the number of tokens it contains 1095 ! The standard separator is generalized whitespace (space, tab, CR, or LF) 1096 in_token=.false.;ntokens=0;last_pos=0 1097 do ii=1,len_trim(str) 1098 cc=str(ii:ii) 1099 if (in_token) then 1100 if (cc==char(9).or.cc==char(10).or.cc==char(13).or.cc==char(32)) then 1101 in_token = .false. 1102 if (cc==char(10).or.cc==char(13)) str(ii:ii) = " " 1103 else 1104 last_pos=ii 1105 end if 1106 else 1107 if (cc==char(9).or.cc==char(10).or.cc==char(13).or.cc==char(32)) then 1108 if (cc==char(10).or.cc==char(13)) str(ii:ii) = " " 1109 else 1110 in_token=.true. 1111 last_pos=ii 1112 ntokens=ntokens + 1 1113 end if 1114 end if 1115 end do 1116 1117 if ((ndata+ntokens)>size(x)) then 1118 msg="data array full" 1119 MSG_ERROR(msg) 1120 end if 1121 1122 ! Take the string and turn it into useful reals 1123 read(unit=str(1:last_pos),fmt=*,iostat=status) x(ndata+1:ndata+ntokens) 1124 if (status/=0) then 1125 msg="real conversion error" 1126 MSG_ERROR(msg) 1127 end if 1128 ndata=ndata+ntokens 1129 1130 end if 1131 1132 end subroutine pawdata_chunk
m_pawxmlps/radial_func_t [ Types ]
[ Top ] [ m_pawxmlps ] [ Types ]
NAME
radial_func_t
FUNCTION
Radial function type for the FoX XML reader
SOURCE
97 type, public :: radialfunc_t 98 logical :: tread=.false. 99 character(len=6) :: grid=' ' !vz_z 100 character(len=6) :: state=' ' !vz_z 101 real(dpxml),allocatable :: data(:) 102 end type radialfunc_t
m_pawxmlps/radial_grid_t [ Types ]
[ Top ] [ m_pawxmlps ] [ Types ]
NAME
radial_grid_t
FUNCTION
Radial grid type for the FoX XML reader
SOURCE
73 type, public :: radial_grid_t 74 logical :: tread=.false. 75 character(len=20) :: eq 76 character(len=4) :: id 77 real(dpxml) :: aa 78 real(dpxml) :: bb 79 real(dpxml) :: dd 80 integer :: nn 81 integer :: istart 82 integer :: iend 83 end type radial_grid_t
m_pawxmlps/rdpawpsxml [ Functions ]
[ Top ] [ m_pawxmlps ] [ Functions ]
NAME
rdpawpsxml
FUNCTION
Read the PAW pseudopotential XML file generated by AtomPAW
INPUTS
filename= input file name (atomicdata XML)
OUTPUT
paw_setup=pseudopotential data structure
PARENTS
pawpsxml2ab
CHILDREN
paw_rdfromline
SOURCE
1996 subroutine rdpawpsxml(filename,paw_setup) 1997 1998 1999 !This section has been created automatically by the script Abilint (TD). 2000 !Do not modify the following lines by hand. 2001 #undef ABI_FUNC 2002 #define ABI_FUNC 'rdpawpsxml' 2003 !End of the abilint section 2004 2005 implicit none 2006 2007 !Arguments --------------------------------------------- 2008 character (len=fnlen),intent(in) :: filename 2009 type(paw_setup_t),intent(inout) :: paw_setup 2010 !Local variables --------------------------------------- 2011 integer :: funit, iaewf,ii,ipswf,iproj,ir,igrid,ival,ierr,ishpf,lmax,mesh_size,igauss,iprojfit 2012 logical :: endfile,found,endgauss 2013 character(len=100) :: msg 2014 character (len=XML_RECL) :: line,readline 2015 character (len=XML_RECL) :: strg 2016 character (len=30) :: strg1 2017 real(dp) :: rc(6) 2018 real(dp), allocatable :: shpf(:,:) 2019 type(state_t), pointer :: valstate (:) 2020 type(radial_grid_t), pointer :: grids (:) 2021 2022 ! ************************************************************************* 2023 2024 !Open the atomicdata XML file for reading 2025 open(newunit=funit,file=filename,form='formatted',status='old', recl=XML_RECL) 2026 2027 !Start a reading loop 2028 endfile=.false. 2029 found=.false. 2030 paw_setup%rpaw=-1.d0 2031 rc=-1.d0 2032 2033 do while ((.not.endfile).and.(.not.found)) 2034 read(funit,'(a)',err=10,end=10) readline 2035 line=adjustl(readline);goto 20 2036 10 line="";endfile=.true. 2037 20 continue 2038 2039 ! --Read VERSION 2040 if ((line(1:10)=='<paw_setup').or.(line(1:12)=='<paw_dataset')) then 2041 paw_setup%tread=.true. 2042 igrid=0;ishpf=0 2043 LIBPAW_DATATYPE_ALLOCATE(grids,(10)) 2044 2045 call paw_rdfromline(" version",line,strg,ierr) 2046 paw_setup%version=trim(strg) 2047 cycle 2048 end if 2049 2050 ! --Read TITLE, ATOMIC CHARGE AND CORE CHARGE 2051 if (line(1:6)=='<atom ') then 2052 paw_setup%atom%tread=.true. 2053 call paw_rdfromline(" symbol",line,strg,ierr) 2054 paw_setup%atom%symbol=trim(strg) 2055 call paw_rdfromline(" Z",line,strg,ierr) 2056 if (len(trim(strg))<=30) then 2057 strg1=trim(strg) 2058 read(unit=strg1,fmt=*) paw_setup%atom%znucl 2059 else 2060 read(unit=strg,fmt=*) paw_setup%atom%znucl 2061 end if 2062 call paw_rdfromline(" core",line,strg,ierr) 2063 if (len(trim(strg))<=30) then 2064 strg1=trim(strg) 2065 read(unit=strg1,fmt=*) paw_setup%atom%zion 2066 else 2067 read(unit=strg,fmt=*) paw_setup%atom%zion 2068 end if 2069 call paw_rdfromline(" valence",line,strg,ierr) 2070 if (len(trim(strg))<=30) then 2071 strg1=trim(strg) 2072 read(unit=strg1,fmt=*) paw_setup%atom%zval 2073 else 2074 read(unit=strg,fmt=*) paw_setup%atom%zval 2075 end if 2076 cycle 2077 end if 2078 2079 ! --Read EXCHANGE-CORRELATION TYPE 2080 if (line(1:14)=='<xc_functional') then 2081 paw_setup%xc_functional%tread=.true. 2082 call paw_rdfromline(" type",line,strg,ierr) 2083 paw_setup%xc_functional%functionaltype = trim(strg) 2084 call paw_rdfromline(" name",line,strg,ierr) 2085 paw_setup%xc_functional%name= trim(strg) 2086 cycle 2087 end if 2088 2089 ! --Read GENERATOR 2090 if (line(1:10)=='<generator') then 2091 paw_setup%generator%tread=.true. 2092 call paw_rdfromline(" type",line,strg,ierr) 2093 paw_setup%generator%gen = trim(strg) 2094 call paw_rdfromline(" name",line,strg,ierr) 2095 paw_setup%generator%name= trim(strg) 2096 cycle 2097 end if 2098 2099 ! --Read PAW RADIUS 2100 if (line(1:11)=='<PAW_radius') then 2101 call paw_rdfromline(" rpaw",line,strg,ierr) 2102 if (len(trim(strg))<=30) then 2103 strg1=trim(strg) 2104 read(unit=strg1,fmt=*) paw_setup%rpaw 2105 else 2106 read(unit=strg,fmt=*) paw_setup%rpaw 2107 end if 2108 cycle 2109 end if 2110 if (line(1:11)=='<paw_radius') then 2111 call paw_rdfromline(" rc",line,strg,ierr) 2112 if (len(trim(strg))<=30) then 2113 strg1=trim(strg) 2114 read(unit=strg1,fmt=*) paw_setup%rpaw 2115 else 2116 read(unit=strg,fmt=*) paw_setup%rpaw 2117 end if 2118 cycle 2119 end if 2120 2121 ! --Read BASIS SIZE, ORBITALS, RC AND OCCUPATIONS/STATE IDs 2122 if (line(1:16)=='<valence_states>') then 2123 paw_setup%valence_states%tread=.true. 2124 LIBPAW_DATATYPE_ALLOCATE(valstate,(50)) 2125 ival=0 2126 lmax=0 2127 do while (line(1:17)/='</valence_states>') 2128 read(funit,'(a)') readline;line=adjustl(readline) 2129 if (line(1:6)=='<state') then 2130 ival=ival+1 2131 if (ival>50) then 2132 close(funit) 2133 msg="Error in rdpawps1xml: basis size too large (>50)!" 2134 MSG_ERROR(msg) 2135 end if 2136 call paw_rdfromline(" n",line,strg,ierr) 2137 if (strg == "" ) then 2138 valstate(ival)%nn=-1 2139 else 2140 if (len(trim(strg))<=30) then 2141 strg1=trim(strg) 2142 read(unit=strg1,fmt=*) valstate(ival)%nn 2143 else 2144 read(unit=strg,fmt=*) valstate(ival)%nn 2145 end if 2146 end if 2147 call paw_rdfromline(" l",line,strg,ierr) 2148 if (len(trim(strg))<=30) then 2149 strg1=trim(strg) 2150 read(unit=strg1,fmt=*) valstate(ival)%ll 2151 else 2152 read(unit=strg,fmt=*) valstate(ival)%ll 2153 end if 2154 if(valstate(ival)%ll>lmax) lmax=valstate(ival)%ll 2155 call paw_rdfromline(" f",line,strg,ierr) 2156 if (strg == "" ) then 2157 valstate(ival)%ff=-1.d0 2158 else 2159 if (len(trim(strg))<=30) then 2160 strg1=trim(strg) 2161 read(unit=strg1,fmt=*) valstate(ival)%ff 2162 else 2163 read(unit=strg,fmt=*) valstate(ival)%ff 2164 end if 2165 end if 2166 call paw_rdfromline(" rc",line,strg,ierr) 2167 if (len(trim(strg))<=30) then 2168 strg1=trim(strg) 2169 read(unit=strg1,fmt=*) valstate(ival)%rc 2170 else 2171 read(unit=strg,fmt=*) valstate(ival)%rc 2172 end if 2173 call paw_rdfromline(" e",line,strg,ierr) 2174 if (len(trim(strg))<=30) then 2175 strg1=trim(strg) 2176 read(unit=strg1,fmt=*) valstate(ival)%ee 2177 else 2178 read(unit=strg,fmt=*) valstate(ival)%ee 2179 end if 2180 call paw_rdfromline(" id",line,strg,ierr) 2181 valstate(ival)%id = trim(strg) 2182 end if 2183 end do 2184 cycle 2185 end if 2186 2187 ! --Read MESH_STEP AND NUMBER OF POINTS 2188 if (line(1:12)=='<radial_grid')then 2189 igrid=igrid+1 2190 call paw_rdfromline(" eq",line,strg,ierr) 2191 grids(igrid)%eq = trim(strg) 2192 call paw_rdfromline(" a",line,strg,ierr) 2193 if (strg == "" ) then 2194 grids(igrid)%aa=0.d0 2195 else 2196 if (len(trim(strg))<=30) then 2197 strg1=trim(strg) 2198 read(unit=strg1,fmt=*) grids(igrid)%aa 2199 else 2200 read(unit=strg,fmt=*) grids(igrid)%aa 2201 end if 2202 end if 2203 call paw_rdfromline(" n",line,strg,ierr) 2204 if (strg == "" ) then 2205 grids(igrid)%nn=0 2206 else 2207 if (len(trim(strg))<=30) then 2208 strg1=trim(strg) 2209 read(unit=strg1,fmt=*) grids(igrid)%nn 2210 else 2211 read(unit=strg,fmt=*) grids(igrid)%nn 2212 end if 2213 end if 2214 call paw_rdfromline(" d",line,strg,ierr) 2215 if (strg == "" ) then 2216 grids(igrid)%dd=0.d0 2217 else 2218 if (len(trim(strg))<=30) then 2219 strg1=trim(strg) 2220 read(unit=strg1,fmt=*) grids(igrid)%dd 2221 else 2222 read(unit=strg,fmt=*) grids(igrid)%dd 2223 end if 2224 end if 2225 call paw_rdfromline(" b",line,strg,ierr) 2226 if (strg == "" ) then 2227 grids(igrid)%bb=0.d0 2228 else 2229 if (len(trim(strg))<=30) then 2230 strg1=trim(strg) 2231 read(unit=strg1,fmt=*) grids(igrid)%bb 2232 else 2233 read(unit=strg,fmt=*) grids(igrid)%bb 2234 end if 2235 end if 2236 call paw_rdfromline("istart",line,strg,ierr) 2237 if (len(trim(strg))<=30) then 2238 strg1=trim(strg) 2239 read(unit=strg1,fmt=*) grids(igrid)%istart 2240 else 2241 read(unit=strg,fmt=*) grids(igrid)%istart 2242 end if 2243 call paw_rdfromline("iend",line,strg,ierr) 2244 if (len(trim(strg))<=30) then 2245 strg1=trim(strg) 2246 read(unit=strg1,fmt=*) grids(igrid)%iend 2247 else 2248 read(unit=strg,fmt=*) grids(igrid)%iend 2249 end if 2250 call paw_rdfromline(" id",line,strg,ierr) 2251 grids(igrid)%id = trim(strg) 2252 if(igrid>10)then 2253 close(funit) 2254 msg="igrid>10" 2255 MSG_ERROR(msg) 2256 end if 2257 cycle 2258 end if 2259 2260 ! --Read SHAPE TYPE 2261 if (line(1:15)=='<shape_function') then 2262 paw_setup%shape_function%tread=.true. 2263 call paw_rdfromline(" type",line,strg,ierr) 2264 paw_setup%shape_function%gtype = trim(strg) 2265 call paw_rdfromline(" rc",line,strg,ierr) 2266 if (strg /= "" ) then 2267 if (len(trim(strg))<=30) then 2268 strg1=trim(strg) 2269 read(unit=strg1,fmt=*) paw_setup%shape_function%rc 2270 else 2271 read(unit=strg,fmt=*) paw_setup%shape_function%rc 2272 end if 2273 end if 2274 call paw_rdfromline(" lamb",line,strg,ierr) 2275 if (strg == "" ) then 2276 paw_setup%shape_function%lamb=0 2277 else 2278 if (len(trim(strg))<=30) then 2279 strg1=trim(strg) 2280 read(unit=strg1,fmt=*) paw_setup%shape_function%lamb 2281 else 2282 read(unit=strg,fmt=*) paw_setup%shape_function%lamb 2283 end if 2284 end if 2285 found=paw_setup%shape_function%tread 2286 call paw_rdfromline("grid",line,strg,ierr) 2287 paw_setup%shape_function%grid=trim(strg) 2288 if (strg /= "" ) then 2289 paw_setup%shape_function%gtype ="num" 2290 do ii=1,igrid 2291 if(trim(paw_setup%shape_function%grid)==trim(grids(ii)%id)) then 2292 mesh_size=grids(ii)%iend-grids(ii)%istart+1 2293 exit 2294 end if 2295 end do 2296 if(.not.allocated(shpf)) then 2297 LIBPAW_ALLOCATE(shpf,(mesh_size,7)) 2298 end if 2299 ishpf=ishpf+1 2300 read(funit,*) (shpf(ir,ishpf),ir=1,mesh_size) 2301 call paw_rdfromline(" l",line,strg,ierr) 2302 if (strg /= "" ) then 2303 found=.false. 2304 if(paw_setup%valence_states%tread) then 2305 if(ishpf==2*lmax+1) found=.true. 2306 else 2307 write(msg,'(a,a,a)')"the grids and the states must be read before the shapefunction",ch10,& 2308 & "Action: Modify your XML PAW data file" 2309 MSG_ERROR(msg) 2310 end if 2311 end if 2312 end if 2313 cycle 2314 end if 2315 2316 ! End of reading loop 2317 end do 2318 2319 if(igrid==0.or.ival==0) then 2320 write(msg,'(a,a,a)')"the grids and the states must be read before the shapefunction",ch10,& 2321 & "Action: Modify your XML PAW data file" 2322 MSG_ERROR(msg) 2323 end if 2324 if(ishpf>0)then 2325 LIBPAW_ALLOCATE(paw_setup%shape_function%data,(mesh_size,ishpf)) 2326 do ii=1,ishpf 2327 paw_setup%shape_function%data(:,ii)=shpf(:,ii) 2328 end do 2329 LIBPAW_DEALLOCATE(shpf) 2330 end if 2331 2332 if(ival>0)then 2333 LIBPAW_DATATYPE_ALLOCATE(paw_setup%valence_states%state,(ival)) 2334 paw_setup%valence_states%state(ival)%tread=.true. 2335 paw_setup%valence_states%nval=ival 2336 do ii=1,ival 2337 paw_setup%valence_states%state(ii)=valstate(ii) 2338 end do 2339 end if 2340 LIBPAW_DATATYPE_DEALLOCATE(valstate) 2341 if(.not.allocated(paw_setup%ae_partial_wave)) then 2342 LIBPAW_DATATYPE_ALLOCATE(paw_setup%ae_partial_wave,(paw_setup%valence_states%nval)) 2343 end if 2344 if(.not.allocated(paw_setup%pseudo_partial_wave)) then 2345 LIBPAW_DATATYPE_ALLOCATE(paw_setup%pseudo_partial_wave,(paw_setup%valence_states%nval)) 2346 end if 2347 if(.not.allocated(paw_setup%projector_function)) then 2348 LIBPAW_DATATYPE_ALLOCATE(paw_setup%projector_function,(paw_setup%valence_states%nval)) 2349 end if 2350 2351 LIBPAW_DATATYPE_ALLOCATE(paw_setup%radial_grid,(igrid)) 2352 paw_setup%radial_grid(igrid)%tread=.true. 2353 paw_setup%ngrid=igrid 2354 do ii=1,igrid 2355 paw_setup%radial_grid(ii)=grids(ii) 2356 end do 2357 LIBPAW_DATATYPE_DEALLOCATE(grids) 2358 2359 !Start a reading loop 2360 ipswf=0;iaewf=0;iproj=0;iprojfit=0 2361 endfile=.false. 2362 do while (.not.endfile) 2363 read(funit,'(a)',err=11,end=11) readline 2364 line=adjustl(readline);goto 21 2365 11 line="";endfile=.true. 2366 21 continue 2367 2368 ! --Read core density CORE_DENSITY 2369 if (line(1:16)=='<ae_core_density') then 2370 paw_setup%ae_core_density%tread=.true. 2371 call paw_rdfromline(" grid",line,strg,ierr) 2372 if (strg == "" ) strg = "unknown" 2373 paw_setup%ae_core_density%grid=trim(strg) 2374 do ii=1,paw_setup%ngrid 2375 if(trim(paw_setup%ae_core_density%grid)==trim(paw_setup%radial_grid(ii)%id)) then 2376 mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1 2377 exit 2378 end if 2379 end do 2380 call paw_rdfromline(" rc",line,strg,ierr) 2381 if (strg /= "" ) then 2382 if (len(trim(strg))<=30) then 2383 strg1=trim(strg) 2384 read(unit=strg1,fmt=*) rc(1) 2385 else 2386 read(unit=strg,fmt=*) rc(1) 2387 end if 2388 end if 2389 LIBPAW_ALLOCATE(paw_setup%ae_core_density%data,(mesh_size)) 2390 !MGNAG v7[62] 2391 ! Runtime Error: m_pawxmlps_cpp.f90, line 1657: 2392 ! Record too long for input bufferProgram terminated by I/O error on unit 9 2393 ! (File="/home/buildbot/ABINIT_OD/petrus_nag/gmatteo_7.7.1-training/tests/Psps_for_tests/Al.LDA",Formatted,Sequential) 2394 read(funit,*) (paw_setup%ae_core_density%data(ir),ir=1,mesh_size) 2395 cycle 2396 end if 2397 2398 ! --Read pseudized core density CORETAIL_DENSITY 2399 if (line(1:20)=='<pseudo_core_density') then 2400 paw_setup%pseudo_core_density%tread=.true. 2401 call paw_rdfromline(" grid",line,strg,ierr) 2402 if (strg == "" ) strg = "unknown" 2403 paw_setup%pseudo_core_density%grid=trim(strg) 2404 do ii=1,paw_setup%ngrid 2405 if(trim(paw_setup%pseudo_core_density%grid)==trim(paw_setup%radial_grid(ii)%id)) then 2406 mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1 2407 exit 2408 end if 2409 end do 2410 call paw_rdfromline(" rc",line,strg,ierr) 2411 if (strg /= "" ) then 2412 if (len(trim(strg))<=30) then 2413 strg1=trim(strg) 2414 read(unit=strg1,fmt=*) rc(2) 2415 else 2416 read(unit=strg,fmt=*) rc(2) 2417 end if 2418 end if 2419 LIBPAW_ALLOCATE(paw_setup%pseudo_core_density%data,(mesh_size)) 2420 read(funit,*) (paw_setup%pseudo_core_density%data(ir),ir=1,mesh_size) 2421 cycle 2422 end if 2423 2424 ! --Read pseudized valence density PSEUDO_VALENCE_DENSITY 2425 if (line(1:23)=='<pseudo_valence_density') then 2426 paw_setup%pseudo_valence_density%tread=.true. 2427 call paw_rdfromline(" grid",line,strg,ierr) 2428 if (strg == "" ) strg = "unknown" 2429 paw_setup%pseudo_valence_density%grid=trim(strg) 2430 do ii=1,paw_setup%ngrid 2431 if(trim(paw_setup%pseudo_valence_density%grid)==trim(paw_setup%radial_grid(ii)%id)) then 2432 mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1 2433 exit 2434 end if 2435 end do 2436 call paw_rdfromline(" rc",line,strg,ierr) 2437 if (strg /= "" ) then 2438 if (len(trim(strg))<=30) then 2439 strg1=trim(strg) 2440 read(unit=strg1,fmt=*) rc(3) 2441 else 2442 read(unit=strg,fmt=*) rc(3) 2443 end if 2444 end if 2445 LIBPAW_ALLOCATE(paw_setup%pseudo_valence_density%data,(mesh_size)) 2446 read(funit,*) (paw_setup%pseudo_valence_density%data(ir),ir=1,mesh_size) 2447 cycle 2448 end if 2449 2450 ! --Read Vbare potential VLOCFUN 2451 if (line(1:15)=='<zero_potential') then 2452 paw_setup%zero_potential%tread=.true. 2453 call paw_rdfromline(" grid",line,strg,ierr) 2454 if (strg == "" ) strg = "unknown" 2455 paw_setup%zero_potential%grid=trim(strg) 2456 do ii=1,paw_setup%ngrid 2457 if(trim(paw_setup%zero_potential%grid)==trim(paw_setup%radial_grid(ii)%id)) then 2458 mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1 2459 exit 2460 end if 2461 end do 2462 call paw_rdfromline(" rc",line,strg,ierr) 2463 if (strg /= "" ) then 2464 if (len(trim(strg))<=30) then 2465 strg1=trim(strg) 2466 read(unit=strg1,fmt=*) rc(4) 2467 else 2468 read(unit=strg,fmt=*) rc(4) 2469 end if 2470 end if 2471 LIBPAW_ALLOCATE(paw_setup%zero_potential%data,(mesh_size)) 2472 read(funit,*) (paw_setup%zero_potential%data(ir),ir=1,mesh_size) 2473 cycle 2474 end if 2475 ! --Read external potential 2476 if (line(1:25)=='<LDA_minus_half_potential') then 2477 paw_setup%LDA_minus_half_potential%tread=.true. 2478 call paw_rdfromline(" grid",line,strg,ierr) 2479 if (strg == "" ) strg = "unknown" 2480 paw_setup%LDA_minus_half_potential%grid=trim(strg) 2481 do ii=1,paw_setup%ngrid 2482 if(trim(paw_setup%LDA_minus_half_potential%grid)==trim(paw_setup%radial_grid(ii)%id)) then 2483 mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1 2484 exit 2485 end if 2486 end do 2487 call paw_rdfromline(" rc",line,strg,ierr) 2488 if (strg /= "" ) then 2489 if (len(trim(strg))<=30) then 2490 strg1=trim(strg) 2491 read(unit=strg1,fmt=*) rc(4) 2492 else 2493 read(unit=strg,fmt=*) rc(4) 2494 end if 2495 end if 2496 LIBPAW_ALLOCATE(paw_setup%LDA_minus_half_potential%data,(mesh_size)) 2497 read(funit,*) (paw_setup%LDA_minus_half_potential%data(ir),ir=1,mesh_size) 2498 cycle 2499 end if 2500 ! --Read Vloc for Abinit potential VLOC_ION 2501 if (line(1:37)=='<kresse_joubert_local_ionic_potential') then 2502 paw_setup%kresse_joubert_local_ionic_potential%tread=.true. 2503 call paw_rdfromline(" grid",line,strg,ierr) 2504 if (strg == "" ) strg = "unknown" 2505 paw_setup%kresse_joubert_local_ionic_potential%grid=trim(strg) 2506 do ii=1,paw_setup%ngrid 2507 if(trim(paw_setup%kresse_joubert_local_ionic_potential%grid)==trim(paw_setup%radial_grid(ii)%id)) then 2508 mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1 2509 exit 2510 end if 2511 end do 2512 call paw_rdfromline(" rc",line,strg,ierr) 2513 if (strg /= "" ) then 2514 if (len(trim(strg))<=30) then 2515 strg1=trim(strg) 2516 read(unit=strg1,fmt=*) rc(5) 2517 else 2518 read(unit=strg,fmt=*) rc(5) 2519 end if 2520 end if 2521 LIBPAW_ALLOCATE(paw_setup%kresse_joubert_local_ionic_potential%data,(mesh_size)) 2522 read(funit,*) (paw_setup%kresse_joubert_local_ionic_potential%data(ir),ir=1,mesh_size) 2523 cycle 2524 end if 2525 if (line(1:29)=='<blochl_local_ionic_potential') then 2526 paw_setup%blochl_local_ionic_potential%tread=.true. 2527 call paw_rdfromline(" grid",line,strg,ierr) 2528 if (strg == "" ) strg = "unknown" 2529 paw_setup%blochl_local_ionic_potential%grid=trim(strg) 2530 do ii=1,paw_setup%ngrid 2531 if(trim(paw_setup%blochl_local_ionic_potential%grid)==trim(paw_setup%radial_grid(ii)%id)) then 2532 mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1 2533 exit 2534 end if 2535 end do 2536 call paw_rdfromline(" rc",line,strg,ierr) 2537 if (strg /= "" ) then 2538 if (len(trim(strg))<=30) then 2539 strg1=trim(strg) 2540 read(unit=strg1,fmt=*) rc(6) 2541 else 2542 read(unit=strg,fmt=*) rc(6) 2543 end if 2544 end if 2545 LIBPAW_ALLOCATE(paw_setup%blochl_local_ionic_potential%data,(mesh_size)) 2546 read(funit,*) (paw_setup%blochl_local_ionic_potential%data(ir),ir=1,mesh_size) 2547 cycle 2548 end if 2549 2550 ! --Read WAVE FUNCTIONS PHI 2551 if (line(1:16)=='<ae_partial_wave') then 2552 iaewf=iaewf+1 2553 paw_setup%ae_partial_wave(iaewf)%tread=.true. 2554 call paw_rdfromline(" grid",line,strg,ierr) 2555 if (strg == "" ) strg = "unknown" 2556 paw_setup%ae_partial_wave(iaewf)%grid=trim(strg) 2557 call paw_rdfromline(" state",line,strg,ierr) 2558 paw_setup%ae_partial_wave(iaewf)%state=trim(strg) 2559 do ii=1,paw_setup%ngrid 2560 if(trim(paw_setup%ae_partial_wave(iaewf)%grid)==trim(paw_setup%radial_grid(ii)%id)) then 2561 mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1 2562 exit 2563 end if 2564 end do 2565 LIBPAW_ALLOCATE(paw_setup%ae_partial_wave(iaewf)%data,(mesh_size)) 2566 read(funit,*) (paw_setup%ae_partial_wave(iaewf)%data(ir),ir=1,mesh_size) 2567 cycle 2568 end if 2569 2570 ! --Read PSEUDO WAVE FUNCTIONS TPHI 2571 if (line(1:20)=='<pseudo_partial_wave') then 2572 ipswf=ipswf+1 2573 paw_setup%pseudo_partial_wave(ipswf)%tread=.true. 2574 call paw_rdfromline(" grid",line,strg,ierr) 2575 if (strg == "" ) strg = "unknown" 2576 paw_setup%idgrid = trim(strg) 2577 paw_setup%pseudo_partial_wave(ipswf)%grid=trim(strg) 2578 call paw_rdfromline(" state",line,strg,ierr) 2579 paw_setup%pseudo_partial_wave(ipswf)%state=trim(strg) 2580 do ii=1,paw_setup%ngrid 2581 if(trim(paw_setup%pseudo_partial_wave(ipswf)%grid)==trim(paw_setup%radial_grid(ii)%id)) then 2582 mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1 2583 exit 2584 end if 2585 end do 2586 LIBPAW_ALLOCATE(paw_setup%pseudo_partial_wave(ipswf)%data,(mesh_size)) 2587 read(funit,*) (paw_setup%pseudo_partial_wave(ipswf)%data(ir),ir=1,mesh_size) 2588 cycle 2589 end if 2590 2591 ! --Read PROJECTORS TPROJ 2592 if (line(1:19)=='<projector_function') then 2593 iproj=iproj+1 2594 paw_setup%projector_function(iproj)%tread=.true. 2595 call paw_rdfromline(" grid",line,strg,ierr) 2596 if (strg == "" ) strg = "unknown" 2597 paw_setup%projector_function(iproj)%grid=trim(strg) 2598 call paw_rdfromline(" state",line,strg,ierr) 2599 paw_setup%projector_function(iproj)%state=trim(strg) 2600 do ii=1,paw_setup%ngrid 2601 if(trim(paw_setup%projector_function(iproj)%grid)==trim(paw_setup%radial_grid(ii)%id)) then 2602 mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1 2603 exit 2604 end if 2605 end do 2606 LIBPAW_ALLOCATE(paw_setup%projector_function(iproj)%data,(mesh_size)) 2607 read(funit,*) (paw_setup%projector_function(iproj)%data(ir),ir=1,mesh_size) 2608 cycle 2609 end if 2610 2611 ! --Read PROJECTORS TPROJ as gaussian representations 2612 if (line(1:14)=='<projector_fit') then 2613 if(.not.allocated(paw_setup%projector_fit)) then 2614 LIBPAW_DATATYPE_ALLOCATE(paw_setup%projector_fit,(paw_setup%valence_states%nval)) 2615 end if 2616 iprojfit=iprojfit+1 2617 paw_setup%projector_fit(iprojfit)%tread=.true. 2618 call paw_rdfromline(" state",line,strg,ierr) 2619 paw_setup%projector_fit(iprojfit)%state=trim(strg) 2620 igauss = 0 2621 endgauss = .false. 2622 do while(.not. endgauss) 2623 read(funit,'(a)',err=12,end=12) readline 2624 line=adjustl(readline);goto 22 2625 12 line="";endgauss=.true. 2626 22 continue 2627 endgauss = (line(1:15)=='</projector_fit') 2628 if (line(1:9)=='<gaussian') then 2629 igauss = igauss + 1 2630 call paw_rdfromline(" factor",line,strg,ierr) 2631 read(strg(index(strg, '{') + 1:index(strg, ',') - 1), *) & 2632 & paw_setup%projector_fit(iprojfit)%factors(1, igauss) 2633 read(strg(index(strg, ',') + 1:index(strg, '}') - 1), *) & 2634 & paw_setup%projector_fit(iprojfit)%factors(2, igauss) 2635 call paw_rdfromline(" exponent",line,strg,ierr) 2636 read(strg(index(strg, '{') + 1:index(strg, ',') - 1), *) & 2637 & paw_setup%projector_fit(iprojfit)%expos(1, igauss) 2638 read(strg(index(strg, ',') + 1:index(strg, '}') - 1), *) & 2639 & paw_setup%projector_fit(iprojfit)%expos(2, igauss) 2640 end if 2641 end do 2642 paw_setup%projector_fit(iprojfit)%ngauss = igauss 2643 cycle 2644 end if 2645 2646 ! --Read Kinetic term KINETIC_ENERGY_MATRIX 2647 if (line(1:28)=='<kinetic_energy_differences>') then 2648 paw_setup%kinetic_energy_differences%tread=.true. 2649 mesh_size=paw_setup%valence_states%nval*paw_setup%valence_states%nval 2650 LIBPAW_ALLOCATE(paw_setup%kinetic_energy_differences%data,(mesh_size)) 2651 read(funit,*) (paw_setup%kinetic_energy_differences%data(ir),ir=1,mesh_size) 2652 cycle 2653 end if 2654 2655 ! --Read Exact exchange term EXACT_EXCHANGE_X_MATRIX 2656 if (line(1:25)=='<exact_exchange_X_matrix>') then 2657 paw_setup%exact_exchange_matrix%tread=.true. 2658 mesh_size=paw_setup%valence_states%nval*paw_setup%valence_states%nval 2659 LIBPAW_ALLOCATE(paw_setup%exact_exchange_matrix%data,(mesh_size)) 2660 read(funit,*) (paw_setup%exact_exchange_matrix%data(ir),ir=1,mesh_size) 2661 cycle 2662 end if 2663 2664 ! --Read Exact exchange core-core energy 2665 if (line(1:25)=='<exact_exchange core-core') then 2666 call paw_rdfromline(" core-core",line,strg,ierr) 2667 if (len(trim(strg))<=30) then 2668 strg1=trim(strg) 2669 read(unit=strg1,fmt=*) paw_setup%ex_cc 2670 else 2671 read(unit=strg,fmt=*) paw_setup%ex_cc 2672 end if 2673 cycle 2674 end if 2675 2676 2677 ! --Read the Atompaw input file 2678 ir=0 2679 if ((line(1:13)=='<!-- Program:').and.(ir==1)) then 2680 msg=" " 2681 do while ((msg(1:9)/=' Program:').and.(msg(1:8)/='Program:')) 2682 read(funit,'(a)') msg 2683 write(ab_out,'(a)') trim(msg) 2684 end do 2685 cycle 2686 end if 2687 2688 ! End of reading loop 2689 end do 2690 if(paw_setup%rpaw<0.d0) paw_setup%rpaw=maxval(rc) 2691 !Close the XML atomicdata file 2692 close(funit) 2693 2694 !Test flags: is anything OK ? 2695 found=paw_setup%atom%tread.and.paw_setup%valence_states%tread.and.& 2696 & paw_setup%xc_functional%tread.and.paw_setup%shape_function%tread 2697 2698 if (.not.paw_setup%atom%tread) then 2699 msg="ATOM SYMBOL not found !" 2700 MSG_WARNING(msg) 2701 end if 2702 if (.not.paw_setup%valence_states%tread) then 2703 msg="VALENCE STATES not found!" 2704 MSG_WARNING(msg) 2705 end if 2706 if (.not.paw_setup%xc_functional%tread) then 2707 msg="EXCHANGE/CORRELATION not found !" 2708 MSG_WARNING(msg) 2709 end if 2710 if (.not.paw_setup%shape_function%tread) then 2711 msg="SHAPE FUNCTION TYPE not found !" 2712 MSG_WARNING(msg) 2713 end if 2714 2715 if (.not.found) then 2716 msg="Aborting now" 2717 MSG_ERROR(msg) 2718 end if 2719 2720 end subroutine rdpawpsxml
m_pawxmlps/rdpawpsxml_core [ Functions ]
[ Top ] [ m_pawxmlps ] [ Functions ]
NAME
rdpawpsxml_core
FUNCTION
Read the core wavefunctions in the XML file generated by AtomPAW
INPUTS
filename= input file name (atomicdata XML) funit= input unit number
OUTPUT
paw_setup=pseudopotential data structure
PARENTS
m_pawpsp
CHILDREN
paw_rdfromline
SOURCE
2746 subroutine rdpawpsxml_core(energy_cor,filename,lcor,ncor,nphicor,pawrad,phi_cor) 2747 2748 2749 !This section has been created automatically by the script Abilint (TD). 2750 !Do not modify the following lines by hand. 2751 #undef ABI_FUNC 2752 #define ABI_FUNC 'rdpawpsxml_core' 2753 !End of the abilint section 2754 2755 implicit none 2756 2757 !Arguments --------------------------------------------- 2758 character (len=fnlen),intent(in) :: filename 2759 integer,intent(out) :: nphicor 2760 !arrays 2761 integer,allocatable,intent(inout) :: lcor(:),ncor(:) 2762 real(dp),allocatable,intent(inout) :: phi_cor(:,:),energy_cor(:) 2763 type(pawrad_type),intent(in) :: pawrad 2764 2765 2766 !Local variables --------------------------------------- 2767 integer :: funit,iaewf,ii,imeshae,imsh,ir,igrid,icor,ierr,maxmeshz,mesh_size,nmesh,shft 2768 logical :: endfile,found, tread 2769 type(pawrad_type) :: tmpmesh 2770 character(len=100) :: msg,version 2771 character (len=XML_RECL) :: line,readline 2772 character (len=XML_RECL) :: strg 2773 character (len=30) :: strg1 2774 integer,allocatable :: mesh_shift(:) 2775 real(dp),allocatable :: work(:),phitmp(:,:) 2776 character (len=20), allocatable :: gridwf(:),statewf(:) 2777 type(state_t),allocatable :: corestate (:) 2778 type(radial_grid_t),allocatable :: grids (:) 2779 type(pawrad_type),allocatable :: radmesh(:) 2780 2781 ! ************************************************************************* 2782 2783 !Open the atomicdata XML file for reading 2784 funit=100 2785 open(unit=funit,file=filename,form='formatted',status='old', recl=XML_RECL) 2786 2787 !Start a reading loop 2788 endfile=.false. 2789 found=.false. 2790 2791 2792 do while ((.not.endfile).and.(.not.found)) 2793 read(funit,'(a)',err=10,end=10) readline 2794 line=adjustl(readline);goto 20 2795 10 line="";endfile=.true. 2796 20 continue 2797 2798 ! --Read VERSION 2799 if (line(1:10)=='<paw_setup') then 2800 tread=.true. 2801 igrid=0 2802 LIBPAW_DATATYPE_ALLOCATE(grids,(10)) 2803 2804 call paw_rdfromline(" version",line,strg,ierr) 2805 version=trim(strg) 2806 cycle 2807 end if 2808 2809 ! --Read BASIS SIZE, ORBITALS, RC AND OCCUPATIONS/STATE IDs 2810 if (line(1:13)=='<core_states>') then 2811 tread=.true. 2812 LIBPAW_DATATYPE_ALLOCATE(corestate,(50)) 2813 icor=0 2814 do while (line(1:14)/='</core_states>') 2815 read(funit,'(a)') readline;line=adjustl(readline) 2816 if (line(1:6)=='<state') then 2817 icor=icor+1 2818 if (icor>50) then 2819 close(funit) 2820 msg="basis size too large (>50)!" 2821 MSG_ERROR(msg) 2822 end if 2823 call paw_rdfromline(" n",line,strg,ierr) 2824 if (strg == "" ) then 2825 corestate(icor)%nn=-1 2826 else 2827 if (len(trim(strg))<=30) then 2828 strg1=trim(strg) 2829 read(unit=strg1,fmt=*) corestate(icor)%nn 2830 else 2831 read(unit=strg,fmt=*) corestate(icor)%nn 2832 end if 2833 end if 2834 call paw_rdfromline(" l",line,strg,ierr) 2835 if (len(trim(strg))<=30) then 2836 strg1=trim(strg) 2837 read(unit=strg1,fmt=*) corestate(icor)%ll 2838 else 2839 read(unit=strg,fmt=*) corestate(icor)%ll 2840 end if 2841 call paw_rdfromline(" f",line,strg,ierr) 2842 if (strg == "" ) then 2843 corestate(icor)%ff=-1.d0 2844 else 2845 if (len(trim(strg))<=30) then 2846 strg1=trim(strg) 2847 read(unit=strg1,fmt=*) corestate(icor)%ff 2848 else 2849 read(unit=strg,fmt=*) corestate(icor)%ff 2850 end if 2851 end if 2852 call paw_rdfromline(" rc",line,strg,ierr) 2853 if (strg == "" ) then 2854 corestate(icor)%rc=-1 2855 else 2856 if (len(trim(strg))<=30) then 2857 strg1=trim(strg) 2858 read(unit=strg1,fmt=*) corestate(icor)%rc 2859 else 2860 read(unit=strg,fmt=*) corestate(icor)%rc 2861 end if 2862 end if 2863 call paw_rdfromline(" e",line,strg,ierr) 2864 if (len(trim(strg))<=30) then 2865 strg1=trim(strg) 2866 read(unit=strg1,fmt=*) corestate(icor)%ee 2867 else 2868 read(unit=strg,fmt=*) corestate(icor)%ee 2869 end if 2870 call paw_rdfromline(" id",line,strg,ierr) 2871 corestate(icor)%id = trim(strg) 2872 end if 2873 end do 2874 cycle 2875 end if 2876 2877 ! --Read MESH_STEP AND NUMBER OF POINTS 2878 if (line(1:12)=='<radial_grid')then 2879 igrid=igrid+1 2880 call paw_rdfromline(" eq",line,strg,ierr) 2881 grids(igrid)%eq = trim(strg) 2882 call paw_rdfromline(" a",line,strg,ierr) 2883 if (strg == "" ) then 2884 grids(igrid)%aa=0.d0 2885 else 2886 if (len(trim(strg))<=30) then 2887 strg1=trim(strg) 2888 read(unit=strg1,fmt=*) grids(igrid)%aa 2889 else 2890 read(unit=strg,fmt=*) grids(igrid)%aa 2891 end if 2892 end if 2893 call paw_rdfromline(" n",line,strg,ierr) 2894 if (strg == "" ) then 2895 grids(igrid)%nn=0 2896 else 2897 if (len(trim(strg))<=30) then 2898 strg1=trim(strg) 2899 read(unit=strg1,fmt=*) grids(igrid)%nn 2900 else 2901 read(unit=strg,fmt=*) grids(igrid)%nn 2902 end if 2903 end if 2904 call paw_rdfromline(" d",line,strg,ierr) 2905 if (strg == "" ) then 2906 grids(igrid)%dd=0.d0 2907 else 2908 if (len(trim(strg))<=30) then 2909 strg1=trim(strg) 2910 read(unit=strg1,fmt=*) grids(igrid)%dd 2911 else 2912 read(unit=strg,fmt=*) grids(igrid)%dd 2913 end if 2914 end if 2915 call paw_rdfromline(" b",line,strg,ierr) 2916 if (strg == "" ) then 2917 grids(igrid)%bb=0.d0 2918 else 2919 if (len(trim(strg))<=30) then 2920 strg1=trim(strg) 2921 read(unit=strg1,fmt=*) grids(igrid)%bb 2922 else 2923 read(unit=strg,fmt=*) grids(igrid)%bb 2924 end if 2925 end if 2926 call paw_rdfromline("istart",line,strg,ierr) 2927 if (len(trim(strg))<=30) then 2928 strg1=trim(strg) 2929 read(unit=strg1,fmt=*) grids(igrid)%istart 2930 else 2931 read(unit=strg,fmt=*) grids(igrid)%istart 2932 end if 2933 call paw_rdfromline("iend",line,strg,ierr) 2934 if (len(trim(strg))<=30) then 2935 strg1=trim(strg) 2936 read(unit=strg1,fmt=*) grids(igrid)%iend 2937 else 2938 read(unit=strg,fmt=*) grids(igrid)%iend 2939 end if 2940 call paw_rdfromline(" id",line,strg,ierr) 2941 grids(igrid)%id = trim(strg) 2942 if(igrid>10)then 2943 close(funit) 2944 msg="igrid>10" 2945 MSG_ERROR(msg) 2946 end if 2947 found=.true. 2948 cycle 2949 end if 2950 2951 ! End of reading loop 2952 end do 2953 2954 nphicor=icor 2955 nmesh=igrid 2956 if(nmesh>0)then 2957 LIBPAW_DATATYPE_ALLOCATE(radmesh,(nmesh)) 2958 LIBPAW_ALLOCATE(mesh_shift,(nmesh)) 2959 do imsh=1,nmesh 2960 radmesh(imsh)%mesh_type=-1 2961 radmesh(imsh)%rstep=zero 2962 radmesh(imsh)%lstep=zero 2963 mesh_shift(imsh)=0 2964 select case(trim(grids(imsh)%eq)) 2965 case("r=a*exp(d*i)") 2966 mesh_shift(imsh)=1 2967 radmesh(imsh)%mesh_type=3 2968 radmesh(imsh)%mesh_size=grids(imsh)%iend-grids(imsh)%istart+1+mesh_shift(imsh) 2969 radmesh(imsh)%rstep=grids(imsh)%aa 2970 radmesh(imsh)%lstep=grids(imsh)%dd 2971 case("r=a*i/(1-b*i)") 2972 write(msg, '(3a)' )& 2973 & ' the grid r=a*i/(1-b*i) is not implemented in ABINIT !',ch10,& 2974 & ' Action: check your psp file.' 2975 MSG_ERROR(msg) 2976 case("r=a*i/(n-i)") 2977 mesh_shift(imsh)=0 2978 radmesh(imsh)%mesh_type=5 2979 radmesh(imsh)%mesh_size=grids(imsh)%iend-grids(imsh)%istart+1+mesh_shift(imsh) 2980 radmesh(imsh)%rstep=grids(imsh)%aa 2981 radmesh(imsh)%lstep=dble(grids(imsh)%nn) 2982 case("r=a*(exp(d*i)-1)") 2983 mesh_shift(imsh)=0 2984 radmesh(imsh)%mesh_type=2 2985 radmesh(imsh)%mesh_size=grids(imsh)%iend-grids(imsh)%istart+1+mesh_shift(imsh) 2986 if(grids(imsh)%istart==1)radmesh(imsh)%mesh_size=radmesh(imsh)%mesh_size+1 2987 radmesh(imsh)%rstep=grids(imsh)%aa 2988 radmesh(imsh)%lstep=grids(imsh)%dd 2989 case("r=d*i") 2990 mesh_shift(imsh)=0 2991 radmesh(imsh)%mesh_type=1 2992 radmesh(imsh)%mesh_size=grids(imsh)%iend-grids(imsh)%istart+1+mesh_shift(imsh) 2993 if(grids(imsh)%istart==1)radmesh(imsh)%mesh_size=radmesh(imsh)%mesh_size+1 2994 radmesh(imsh)%rstep=grids(imsh)%dd 2995 case("r=(i/n+a)^5/a-a^4") 2996 write(msg, '(3a)' )& 2997 & ' the grid r=(i/n+a)^5/a-a^4 is not implemented in ABINIT !',ch10,& 2998 & ' Action: check your psp file.' 2999 MSG_ERROR(msg) 3000 end select 3001 end do 3002 end if 3003 3004 !Start a reading loop 3005 iaewf=0 3006 endfile=.false. 3007 LIBPAW_DATATYPE_ALLOCATE(gridwf,(nphicor)) 3008 LIBPAW_DATATYPE_ALLOCATE(statewf,(nphicor)) 3009 maxmeshz=maxval(radmesh(:)%mesh_size) 3010 LIBPAW_ALLOCATE(phitmp,(nphicor,maxmeshz)) 3011 3012 do while (.not.endfile) 3013 read(funit,'(a)',err=11,end=11) readline 3014 line=adjustl(readline);goto 21 3015 11 line="";endfile=.true. 3016 21 continue 3017 3018 ! --Read CORE WAVE FUNCTIONS PHI 3019 if (line(1:21)=='<ae_core_wavefunction') then 3020 iaewf=iaewf+1 3021 tread=.true. 3022 call paw_rdfromline(" grid",line,strg,ierr) 3023 if (strg == "" ) strg = "unknown" 3024 gridwf(iaewf)=trim(strg) 3025 call paw_rdfromline(" state",line,strg,ierr) 3026 statewf(iaewf)=trim(strg) 3027 do ii=1,nmesh 3028 if(trim(gridwf(iaewf))==trim(grids(ii)%id)) then 3029 mesh_size=grids(ii)%iend-grids(ii)%istart+1 3030 exit 3031 end if 3032 end do 3033 read(funit,*) (phitmp(iaewf,ir),ir=1,mesh_size) 3034 cycle 3035 end if 3036 ! End of reading loop 3037 end do 3038 3039 if(nphicor>0)then 3040 LIBPAW_ALLOCATE(ncor,(nphicor)) 3041 LIBPAW_ALLOCATE(lcor,(nphicor)) 3042 LIBPAW_ALLOCATE(energy_cor,(nphicor)) 3043 LIBPAW_ALLOCATE(phi_cor,(maxmeshz,nphicor)) 3044 do ii=1,nphicor 3045 ncor(ii)=corestate(ii)%nn 3046 lcor(ii)=corestate(ii)%ll 3047 energy_cor(ii)=corestate(ii)%ee 3048 do imsh=1,nmesh 3049 if(trim(gridwf(ii))==trim(grids(imsh)%id)) imeshae=imsh 3050 end do 3051 if ((pawrad%mesh_type/=radmesh(imeshae)%mesh_type) & 3052 & .or.(pawrad%rstep/=radmesh(imeshae)%rstep) & 3053 & .or.(pawrad%lstep/=radmesh(imeshae)%lstep)) then 3054 if(maxmeshz>pawrad%mesh_size) then 3055 write(msg, '(3a)' )& 3056 & ' rdpawpsxml_core:maxmeshz>pawrad%mesh_size',ch10,& 3057 & ' change pseudopotential' 3058 MSG_ERROR(msg) 3059 end if 3060 call pawrad_init(tmpmesh,mesh_size=maxmeshz,mesh_type=radmesh(imeshae)%mesh_type,& 3061 & rstep=radmesh(imeshae)%rstep,lstep=radmesh(imeshae)%lstep) 3062 LIBPAW_ALLOCATE(work,(maxmeshz)) 3063 call bound_deriv(phitmp(ii,:),tmpmesh,maxmeshz,radmesh(imeshae)%rstep,radmesh(imeshae)%lstep) 3064 call paw_spline(tmpmesh%rad,phitmp(ii,:),maxmeshz,radmesh(imeshae)%rstep,radmesh(imeshae)%lstep,work) 3065 call paw_splint(maxmeshz,tmpmesh%rad,phitmp(ii,:),work,maxmeshz,tmpmesh%rad(1:maxmeshz),phi_cor(1:maxmeshz,ii)) 3066 phi_cor(1:maxmeshz,ii)=phi_cor(1:maxmeshz,ii)*tmpmesh%rad(1:maxmeshz) 3067 LIBPAW_DEALLOCATE(work) 3068 call pawrad_free(tmpmesh) 3069 else 3070 shft=mesh_shift(imeshae) 3071 phi_cor(1+shft:radmesh(imeshae)%mesh_size,ii)=phitmp(ii,1:radmesh(imeshae)%mesh_size-shft) 3072 phi_cor(1+shft:maxmeshz,ii)=phi_cor(1+shft:maxmeshz,ii)*radmesh(imeshae)%rad(1+shft:maxmeshz) 3073 if (radmesh(imeshae)%mesh_size<maxmeshz) phi_cor(radmesh(imeshae)%mesh_size+1:maxmeshz,ii)=zero 3074 if (shft==1) phi_cor(1,ii)=zero 3075 end if 3076 end do 3077 end if 3078 3079 LIBPAW_DATATYPE_DEALLOCATE(radmesh) 3080 LIBPAW_DATATYPE_DEALLOCATE(mesh_shift) 3081 3082 LIBPAW_DATATYPE_DEALLOCATE(grids) 3083 LIBPAW_DATATYPE_DEALLOCATE(corestate) 3084 3085 LIBPAW_DATATYPE_DEALLOCATE(gridwf) 3086 LIBPAW_DATATYPE_DEALLOCATE(statewf) 3087 LIBPAW_DEALLOCATE(phitmp) 3088 3089 !Close the XML atomicdata file 3090 close(funit) 3091 3092 3093 end subroutine rdpawpsxml_core
m_pawxmlps/rdpawpsxml_header [ Functions ]
[ Top ] [ m_pawxmlps ] [ Functions ]
NAME
rdpawpsxml_header
FUNCTION
Read the header of a PAW pseudopotential XML file generated by AtomPAW
INPUTS
filename= input file name (atomicdata XML)
OUTPUT
paw_setup=pseudopotential data structure
PARENTS
pawpsxml2ab
CHILDREN
paw_rdfromline
SOURCE
1600 subroutine rdpawpsxml_header(ecut_tmp,filename,paw_setup) 1601 1602 1603 !This section has been created automatically by the script Abilint (TD). 1604 !Do not modify the following lines by hand. 1605 #undef ABI_FUNC 1606 #define ABI_FUNC 'rdpawpsxml_header' 1607 !End of the abilint section 1608 1609 implicit none 1610 1611 !Arguments --------------------------------------------- 1612 1613 character (len=fnlen),intent(in) :: filename 1614 real(dp), intent(inout) :: ecut_tmp(3,2) 1615 type(paw_setup_t),intent(inout) :: paw_setup 1616 !Local variables --------------------------------------- 1617 integer :: funit,iaewf,ii,ipswf,iproj,ir,igrid,ival,ierr,ishpf,lmax,mesh_size,igauss,iprojfit 1618 logical :: endfile,found,endgauss 1619 character(len=100) :: msg 1620 character (len=XML_RECL) :: line,readline 1621 character (len=XML_RECL) :: strg 1622 character (len=30) :: strg1 1623 real(dp) :: rc(6) 1624 real(dp), allocatable :: shpf(:,:) 1625 type(state_t), pointer :: valstate (:) 1626 type(radial_grid_t), pointer :: grids (:) 1627 1628 ! ************************************************************************* 1629 1630 !Open the atomicdata XML file for reading 1631 open(newunit=funit,file=filename,form='formatted',status='old', recl=XML_RECL) 1632 1633 !Start a reading loop 1634 endfile=.false. 1635 found=.false. 1636 paw_setup%rpaw=-1.d0 1637 rc=-1.d0 1638 1639 do while ((.not.endfile).and.(.not.found)) 1640 read(funit,'(a)',err=10,end=10) readline 1641 line=adjustl(readline);goto 20 1642 10 line="";endfile=.true. 1643 20 continue 1644 1645 ! --Read VERSION 1646 if ((line(1:10)=='<paw_setup').or.(line(1:12)=='<paw_dataset')) then 1647 paw_setup%tread=.true. 1648 igrid=0;ishpf=0 1649 LIBPAW_DATATYPE_ALLOCATE(grids,(10)) 1650 1651 call paw_rdfromline(" version",line,strg,ierr) 1652 paw_setup%version=trim(strg) 1653 cycle 1654 end if 1655 1656 ! --Read TITLE, ATOMIC CHARGE AND CORE CHARGE 1657 if (line(1:6)=='<atom ') then 1658 paw_setup%atom%tread=.true. 1659 call paw_rdfromline(" symbol",line,strg,ierr) 1660 paw_setup%atom%symbol=trim(strg) 1661 call paw_rdfromline(" Z",line,strg,ierr) 1662 if (len(trim(strg))<=30) then 1663 strg1=trim(strg) 1664 read(unit=strg1,fmt=*) paw_setup%atom%znucl 1665 else 1666 read(unit=strg,fmt=*) paw_setup%atom%znucl 1667 end if 1668 call paw_rdfromline(" core",line,strg,ierr) 1669 if (len(trim(strg))<=30) then 1670 strg1=trim(strg) 1671 read(unit=strg1,fmt=*) paw_setup%atom%zion 1672 else 1673 read(unit=strg,fmt=*) paw_setup%atom%zion 1674 end if 1675 call paw_rdfromline(" valence",line,strg,ierr) 1676 if (len(trim(strg))<=30) then 1677 strg1=trim(strg) 1678 read(unit=strg1,fmt=*) paw_setup%atom%zval 1679 else 1680 read(unit=strg,fmt=*) paw_setup%atom%zval 1681 end if 1682 cycle 1683 end if 1684 1685 ! --Read Ecut and Ecutdg 1686 if (line(1:8)=='<pw_ecut') then 1687 call paw_rdfromline(" low",line,strg,ierr) 1688 if (len(trim(strg))<=30) then 1689 strg1=trim(strg) 1690 read(unit=strg1,fmt=*) ecut_tmp(1,1) 1691 else 1692 read(unit=strg,fmt=*) ecut_tmp(1,1) 1693 end if 1694 call paw_rdfromline(" medium",line,strg,ierr) 1695 if (len(trim(strg))<=30) then 1696 strg1=trim(strg) 1697 read(unit=strg1,fmt=*) ecut_tmp(2,1) 1698 else 1699 read(unit=strg,fmt=*) ecut_tmp(2,1) 1700 end if 1701 call paw_rdfromline(" high",line,strg,ierr) 1702 if (len(trim(strg))<=30) then 1703 strg1=trim(strg) 1704 read(unit=strg1,fmt=*) ecut_tmp(3,1) 1705 else 1706 read(unit=strg,fmt=*) ecut_tmp(3,1) 1707 end if 1708 cycle 1709 end if 1710 1711 ! --Read EXCHANGE-CORRELATION TYPE 1712 if (line(1:14)=='<xc_functional') then 1713 paw_setup%xc_functional%tread=.true. 1714 call paw_rdfromline(" type",line,strg,ierr) 1715 paw_setup%xc_functional%functionaltype = trim(strg) 1716 call paw_rdfromline(" name",line,strg,ierr) 1717 paw_setup%xc_functional%name= trim(strg) 1718 cycle 1719 end if 1720 1721 ! --Read GENERATOR 1722 if (line(1:10)=='<generator') then 1723 paw_setup%generator%tread=.true. 1724 call paw_rdfromline(" type",line,strg,ierr) 1725 paw_setup%generator%gen = trim(strg) 1726 call paw_rdfromline(" name",line,strg,ierr) 1727 paw_setup%generator%name= trim(strg) 1728 cycle 1729 end if 1730 1731 ! --Read PAW RADIUS 1732 if (line(1:11)=='<PAW_radius') then 1733 call paw_rdfromline(" rpaw",line,strg,ierr) 1734 if (len(trim(strg))<=30) then 1735 strg1=trim(strg) 1736 read(unit=strg1,fmt=*) paw_setup%rpaw 1737 else 1738 read(unit=strg,fmt=*) paw_setup%rpaw 1739 end if 1740 cycle 1741 end if 1742 if (line(1:11)=='<paw_radius') then 1743 call paw_rdfromline(" rc",line,strg,ierr) 1744 if (len(trim(strg))<=30) then 1745 strg1=trim(strg) 1746 read(unit=strg1,fmt=*) paw_setup%rpaw 1747 else 1748 read(unit=strg,fmt=*) paw_setup%rpaw 1749 end if 1750 cycle 1751 end if 1752 1753 ! --Read BASIS SIZE, ORBITALS, RC AND OCCUPATIONS/STATE IDs 1754 if (line(1:16)=='<valence_states>') then 1755 paw_setup%valence_states%tread=.true. 1756 LIBPAW_DATATYPE_ALLOCATE(valstate,(50)) 1757 ival=0 1758 lmax=0 1759 do while (line(1:17)/='</valence_states>') 1760 read(funit,'(a)') readline;line=adjustl(readline) 1761 if (line(1:6)=='<state') then 1762 ival=ival+1 1763 if (ival>50) then 1764 close(funit) 1765 msg="Error in rdpawps1xml: basis size too large (>50)!" 1766 MSG_ERROR(msg) 1767 end if 1768 call paw_rdfromline(" n",line,strg,ierr) 1769 if (strg == "" ) then 1770 valstate(ival)%nn=-1 1771 else 1772 if (len(trim(strg))<=30) then 1773 strg1=trim(strg) 1774 read(unit=strg1,fmt=*) valstate(ival)%nn 1775 else 1776 read(unit=strg,fmt=*) valstate(ival)%nn 1777 end if 1778 end if 1779 call paw_rdfromline(" l",line,strg,ierr) 1780 if (len(trim(strg))<=30) then 1781 strg1=trim(strg) 1782 read(unit=strg1,fmt=*) valstate(ival)%ll 1783 else 1784 read(unit=strg,fmt=*) valstate(ival)%ll 1785 end if 1786 if(valstate(ival)%ll>lmax) lmax=valstate(ival)%ll 1787 call paw_rdfromline(" f",line,strg,ierr) 1788 if (strg == "" ) then 1789 valstate(ival)%ff=-1.d0 1790 else 1791 if (len(trim(strg))<=30) then 1792 strg1=trim(strg) 1793 read(unit=strg1,fmt=*) valstate(ival)%ff 1794 else 1795 read(unit=strg,fmt=*) valstate(ival)%ff 1796 end if 1797 end if 1798 call paw_rdfromline(" rc",line,strg,ierr) 1799 if (len(trim(strg))<=30) then 1800 strg1=trim(strg) 1801 read(unit=strg1,fmt=*) valstate(ival)%rc 1802 else 1803 read(unit=strg,fmt=*) valstate(ival)%rc 1804 end if 1805 call paw_rdfromline(" e",line,strg,ierr) 1806 if (len(trim(strg))<=30) then 1807 strg1=trim(strg) 1808 read(unit=strg1,fmt=*) valstate(ival)%ee 1809 else 1810 read(unit=strg,fmt=*) valstate(ival)%ee 1811 end if 1812 call paw_rdfromline(" id",line,strg,ierr) 1813 valstate(ival)%id = trim(strg) 1814 end if 1815 end do 1816 cycle 1817 end if 1818 1819 ! --Read MESH_STEP AND NUMBER OF POINTS 1820 if (line(1:12)=='<radial_grid')then 1821 igrid=igrid+1 1822 call paw_rdfromline(" eq",line,strg,ierr) 1823 grids(igrid)%eq = trim(strg) 1824 call paw_rdfromline(" a",line,strg,ierr) 1825 if (strg == "" ) then 1826 grids(igrid)%aa=0.d0 1827 else 1828 if (len(trim(strg))<=30) then 1829 strg1=trim(strg) 1830 read(unit=strg1,fmt=*) grids(igrid)%aa 1831 else 1832 read(unit=strg,fmt=*) grids(igrid)%aa 1833 end if 1834 end if 1835 call paw_rdfromline(" n",line,strg,ierr) 1836 if (strg == "" ) then 1837 grids(igrid)%nn=0 1838 else 1839 if (len(trim(strg))<=30) then 1840 strg1=trim(strg) 1841 read(unit=strg1,fmt=*) grids(igrid)%nn 1842 else 1843 read(unit=strg,fmt=*) grids(igrid)%nn 1844 end if 1845 end if 1846 call paw_rdfromline(" d",line,strg,ierr) 1847 if (strg == "" ) then 1848 grids(igrid)%dd=0.d0 1849 else 1850 if (len(trim(strg))<=30) then 1851 strg1=trim(strg) 1852 read(unit=strg1,fmt=*) grids(igrid)%dd 1853 else 1854 read(unit=strg,fmt=*) grids(igrid)%dd 1855 end if 1856 end if 1857 call paw_rdfromline(" b",line,strg,ierr) 1858 if (strg == "" ) then 1859 grids(igrid)%bb=0.d0 1860 else 1861 if (len(trim(strg))<=30) then 1862 strg1=trim(strg) 1863 read(unit=strg1,fmt=*) grids(igrid)%bb 1864 else 1865 read(unit=strg,fmt=*) grids(igrid)%bb 1866 end if 1867 end if 1868 call paw_rdfromline("istart",line,strg,ierr) 1869 if (len(trim(strg))<=30) then 1870 strg1=trim(strg) 1871 read(unit=strg1,fmt=*) grids(igrid)%istart 1872 else 1873 read(unit=strg,fmt=*) grids(igrid)%istart 1874 end if 1875 call paw_rdfromline("iend",line,strg,ierr) 1876 if (len(trim(strg))<=30) then 1877 strg1=trim(strg) 1878 read(unit=strg1,fmt=*) grids(igrid)%iend 1879 else 1880 read(unit=strg,fmt=*) grids(igrid)%iend 1881 end if 1882 call paw_rdfromline(" id",line,strg,ierr) 1883 grids(igrid)%id = trim(strg) 1884 if(igrid>10)then 1885 close(funit) 1886 msg="igrid>10" 1887 MSG_ERROR(msg) 1888 end if 1889 cycle 1890 end if 1891 1892 ! --Read SHAPE TYPE 1893 if (line(1:15)=='<shape_function') then 1894 paw_setup%shape_function%tread=.true. 1895 call paw_rdfromline(" type",line,strg,ierr) 1896 paw_setup%shape_function%gtype = trim(strg) 1897 call paw_rdfromline(" rc",line,strg,ierr) 1898 if (strg /= "" ) then 1899 if (len(trim(strg))<=30) then 1900 strg1=trim(strg) 1901 read(unit=strg1,fmt=*) paw_setup%shape_function%rc 1902 else 1903 read(unit=strg,fmt=*) paw_setup%shape_function%rc 1904 end if 1905 end if 1906 call paw_rdfromline(" lamb",line,strg,ierr) 1907 if (strg == "" ) then 1908 paw_setup%shape_function%lamb=0 1909 else 1910 if (len(trim(strg))<=30) then 1911 strg1=trim(strg) 1912 read(unit=strg1,fmt=*) paw_setup%shape_function%lamb 1913 else 1914 read(unit=strg,fmt=*) paw_setup%shape_function%lamb 1915 end if 1916 end if 1917 found=paw_setup%shape_function%tread 1918 call paw_rdfromline("grid",line,strg,ierr) 1919 paw_setup%shape_function%grid=trim(strg) 1920 if (strg /= "" ) then 1921 paw_setup%shape_function%gtype ="num" 1922 do ii=1,igrid 1923 if(trim(paw_setup%shape_function%grid)==trim(grids(ii)%id)) then 1924 mesh_size=grids(ii)%iend-grids(ii)%istart+1 1925 exit 1926 end if 1927 end do 1928 if(.not.allocated(shpf)) then 1929 LIBPAW_ALLOCATE(shpf,(mesh_size,7)) 1930 end if 1931 ishpf=ishpf+1 1932 read(funit,*) (shpf(ir,ishpf),ir=1,mesh_size) 1933 call paw_rdfromline(" l",line,strg,ierr) 1934 if (strg /= "" ) then 1935 found=.false. 1936 if(paw_setup%valence_states%tread) then 1937 if(ishpf==2*lmax+1) found=.true. 1938 else 1939 write(msg,'(a,a,a)')"the grids and the states must be read before the shapefunction",ch10,& 1940 & "Action: Modify your XML PAW data file" 1941 MSG_ERROR(msg) 1942 end if 1943 end if 1944 end if 1945 cycle 1946 end if 1947 1948 ! End of reading loop 1949 end do 1950 if(ival>0)then 1951 LIBPAW_DATATYPE_ALLOCATE(paw_setup%valence_states%state,(ival)) 1952 paw_setup%valence_states%state(ival)%tread=.true. 1953 paw_setup%valence_states%nval=ival 1954 do ii=1,ival 1955 paw_setup%valence_states%state(ii)=valstate(ii) 1956 end do 1957 end if 1958 LIBPAW_DATATYPE_DEALLOCATE(valstate) 1959 LIBPAW_DATATYPE_ALLOCATE(paw_setup%radial_grid,(igrid)) 1960 paw_setup%radial_grid(igrid)%tread=.true. 1961 paw_setup%ngrid=igrid 1962 do ii=1,igrid 1963 paw_setup%radial_grid(ii)=grids(ii) 1964 end do 1965 LIBPAW_DATATYPE_DEALLOCATE(grids) 1966 if(allocated(shpf)) then 1967 LIBPAW_DEALLOCATE(shpf) 1968 end if 1969 close(funit) 1970 1971 end subroutine rdpawpsxml_header
m_pawxmlps/shape_function_t [ Types ]
[ Top ] [ m_pawxmlps ] [ Types ]
NAME
shape_function_t
FUNCTION
Shape function type for the FoX XML reader
SOURCE
136 type, public :: shape_function_t 137 logical :: tread=.false. 138 character(len=20) :: gtype 139 real(dpxml) :: rc=0 140 character(len=6) :: grid 141 integer :: lamb 142 real(dpxml), allocatable :: data(:,:) 143 end type shape_function_t
m_pawxmlps/state_t [ Types ]
[ Top ] [ m_pawxmlps ] [ Types ]
NAME
state_t
FUNCTION
State type for the FoX XML reader
SOURCE
157 type, public :: state_t 158 logical :: tread=.false. 159 character(len=6) :: id 160 real(dpxml) :: ff 161 real(dpxml) :: rc 162 real(dpxml) :: ee 163 integer :: nn 164 integer :: ll 165 end type state_t
m_pawxmlps/valence_states_t [ Types ]
[ Top ] [ m_pawxmlps ] [ Types ]
NAME
valence_states_t
FUNCTION
Valence state type for the FoX XML reader
SOURCE
179 type, public :: valence_states_t 180 logical :: tread=.false. 181 integer :: nval 182 type(state_t),allocatable :: state(:) 183 end type valence_states_t
m_pawxmlps/xc_functional_t [ Types ]
[ Top ] [ m_pawxmlps ] [ Types ]
NAME
xc_functional_t function_t
FUNCTION
XC functional type for the FoX XML reader
SOURCE
215 type, public :: xc_functional_t 216 logical :: tread=.false. 217 character(len=12) :: functionaltype 218 character(len=100) :: name 219 end type xc_functional_t