TABLE OF CONTENTS


ABINIT/m_pawxmlps [ Modules ]

[ Top ] [ 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-2022 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

21 #include "libpaw.h"
22 
23 module m_pawxmlps
24 
25  USE_DEFS
26  USE_MSG_HANDLING
27  USE_MEMORY_PROFILING
28 
29 #if defined LIBPAW_HAVE_FOX
30  use fox_sax
31 #endif
32 
33  use m_pawrad     , only : pawrad_type, pawrad_init, pawrad_free, pawrad_ifromr, bound_deriv
34  use m_paw_numeric, only : paw_spline, paw_splint
35 
36  implicit none
37 
38  private
39 
40 !Procedures used for the Fortran reader
41  public  ::  rdpawpsxml
42  public  ::  rdpawpsxml_header
43  public  ::  rdpawpsxml_core
44 
45 !Procedures used for the FoX reader (called from xml_parser in response to events)
46 #if defined LIBPAW_HAVE_FOX
47  public :: paw_begin_element1
48  public :: paw_end_element1
49  public :: pawdata_chunk
50 #endif
51 
52 ! private procedures and global variables.
53  private ::  paw_rdfromline
54 
55 ! This type of real is used by the datatypes below
56  integer,parameter,private :: dpxml = selected_real_kind(14)
57 
58 ! The maximum length of a record in a file connected for sequential access.
59  integer,parameter,private :: XML_RECL = 50000
60  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

115 type, public              :: gaussian_expansion_t
116   logical             :: tread=.false.
117   integer             :: ngauss=0
118   character(len=6)    :: state=' '
119   real(dpxml), dimension(2, NGAUSSIAN_MAX) :: factors
120   real(dpxml), dimension(2, NGAUSSIAN_MAX) :: expos
121 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.

SOURCE

346 subroutine paw_begin_element1(namespaceURI,localName,name,attributes)
347 
348 character(len=*),intent(in)   :: namespaceURI,localName,name
349 type(dictionary_t),intent(in) :: attributes
350 
351 character(len=100)  :: msg,value
352 integer ::iaewf=0,iproj=0,ipswf=0,iprojfit=0,igauss=0
353 !Just to fool abirules
354  value=localName
355  value=namespaceURI
356 
357 select case(name)
358 
359       case ("paw_setup")
360         paw_setuploc%tread=.true.
361         igrid=0;ishpf=0
362         paw_setuploc%rpaw=-1.d0
363         LIBPAW_DATATYPE_ALLOCATE(grids,(10))
364         LIBPAW_DATATYPE_ALLOCATE(shpf,(7))
365         value = getValue(attributes,"version")
366         write(std_out,'(3a)') "Processing a PSEUDO version ",trim(value)," XML file"
367         paw_setuploc%version=trim(value)
368 
369 
370       case ("atom")
371          paw_setuploc%atom%tread=.true.
372          value = getValue(attributes,"symbol")
373          if (value == "" ) then
374            msg="Cannot determine atomic symbol"
375            LIBPAW_ERROR(msg)
376          end if
377          paw_setuploc%atom%symbol = trim(value)
378 
379          value = getValue(attributes,"Z") 
380          if (value == "" ) then
381            msg="Cannot determine znucl"
382            LIBPAW_ERROR(msg)
383          end if
384          read(unit=value,fmt=*) paw_setuploc%atom%znucl
385 
386          value = getValue(attributes,"core")
387          if (value == "" ) then
388            msg="Cannot determine zion"
389            LIBPAW_ERROR(msg)
390          end if
391          read(unit=value,fmt=*) paw_setuploc%atom%zion
392 
393          value = getValue(attributes,"valence")
394          if (value == "" ) then
395            msg="Cannot determine zval"
396            LIBPAW_ERROR(msg)
397          end if
398          read(unit=value,fmt=*) paw_setuploc%atom%zval
399 
400 
401       case ("xc_functional")
402          paw_setuploc%xc_functional%tread=.true.
403          value = getValue(attributes,"type")
404          if (value == "" ) then
405            msg="Cannot determine xc-functional-type"
406            LIBPAW_ERROR(msg)
407          end if
408          paw_setuploc%xc_functional%functionaltype = trim(value)
409 
410          value = getValue(attributes,"name")
411          if (value == "" ) then
412            msg="Cannot determine xc-functional-name "
413            LIBPAW_ERROR(msg)
414          end if
415          paw_setuploc%xc_functional%name= trim(value)
416 
417       case ("generator")
418          paw_setuploc%generator%tread=.true.
419          in_generator =.true.
420          value = getValue(attributes,"type")
421          if (value == "" ) value = "unknown"
422          paw_setuploc%generator%gen = trim(value)
423 
424          value = getValue(attributes,"name")
425          if (value == "" ) value = "unknown"
426          paw_setuploc%generator%name = trim(value)
427 
428       case ("PAW_radius")
429          value = getValue(attributes,"rpaw")
430          if (value == "" ) then
431            msg="Cannot determine rpaw"
432            LIBPAW_ERROR(msg)
433          end if
434          read(unit=value,fmt=*) paw_setuploc%rpaw
435 
436       case ("valence_states")
437          paw_setuploc%valence_states%tread=.true.
438          in_valenceStates=.true.
439          ival=0
440          lmax=0
441          LIBPAW_DATATYPE_ALLOCATE(valstate,(50))
442 
443       case ("state")
444          ival=ival+1
445 
446          value = getValue(attributes,"n")
447          if (value == "" ) then 
448            valstate(ival)%nn=-1    
449          else
450            read(unit=value,fmt=*) valstate(ival)%nn
451          end if
452  
453          value = getValue(attributes,"l")
454          if (value == "" ) then
455            msg="Cannot determine l"
456            LIBPAW_ERROR(msg)
457          end if
458          read(unit=value,fmt=*) valstate(ival)%ll
459          if(valstate(ival)%ll>lmax) lmax=valstate(ival)%ll
460 
461          value = getValue(attributes,"f")
462          if (value == "" ) then 
463            valstate(ival)%ff=-1.d0
464          else
465            read(unit=value,fmt=*) valstate(ival)%ff
466          end if
467 
468          value = getValue(attributes,"rc")
469          if (value == "" ) then
470            msg="Cannot determine rc"
471            LIBPAW_ERROR(msg)
472          end if
473          read(unit=value,fmt=*) valstate(ival)%rc
474 
475          value = getValue(attributes,"e")
476          if (value == "" ) then
477            msg="Cannot determine e"
478            LIBPAW_ERROR(msg)
479          end if
480          read(unit=value,fmt=*) valstate(ival)%ee
481 
482          value = getValue(attributes,"id")
483          if (value == "" ) value = "unknown"
484          valstate(ival)%id = trim(value)
485 
486       case ("radial_grid")
487          igrid=igrid+1
488          value = getValue(attributes,"eq")
489          if (value == "" ) value = "unknown"
490          grids(igrid)%eq = trim(value)
491 
492          value = getValue(attributes,"a")
493          if (value == "" ) then
494            grids(igrid)%aa=0.d0
495          else
496            read(unit=value,fmt=*) grids(igrid)%aa
497          end if
498 
499          value = getValue(attributes,"n")
500          if (value == "" ) then
501            grids(igrid)%nn=0
502          else
503            read(unit=value,fmt=*) grids(igrid)%nn
504          end if
505 
506          value = getValue(attributes,"d")
507          if (value == "" ) then
508            grids(igrid)%dd=0.d0
509          else
510            read(unit=value,fmt=*) grids(igrid)%dd
511          end if
512 
513          value = getValue(attributes,"b")
514          if (value == "" ) then
515            grids(igrid)%bb=0.d0
516          else
517            read(unit=value,fmt=*) grids(igrid)%bb
518          end if
519 
520          value = getValue(attributes,"istart")
521          if (value == "" ) then
522            msg="Cannot determine istart"
523            LIBPAW_ERROR(msg)
524          end if
525          read(unit=value,fmt=*) grids(igrid)%istart
526 
527          value = getValue(attributes,"iend")
528          if (value == "" ) then
529            msg="Cannot determine iend"
530            LIBPAW_ERROR(msg)
531          end if
532          read(unit=value,fmt=*) grids(igrid)%iend
533 
534          value = getValue(attributes,"id")
535          if (value == "" ) value = "unknown"
536          grids(igrid)%id = trim(value)
537 
538 end select
539 
540 select case(name)
541       case ("shape_function")
542          paw_setuploc%shape_function%tread=.true.
543          value = getValue(attributes,"type")
544          if (value == "" ) value = "unknown"
545          paw_setuploc%shape_function%gtype = trim(value)
546 
547          value = getValue(attributes,"grid")
548          paw_setuploc%shape_function%grid=trim(value)
549          if (value /= "" ) then
550            paw_setuploc%shape_function%gtype ="num" 
551            do ii=1,igrid
552              if(trim(paw_setuploc%shape_function%grid)==trim(grids(ii)%id)) then
553                mesh_size=grids(ii)%iend-grids(ii)%istart+1
554              end if
555            end do
556            ishpf=ishpf+1
557            LIBPAW_ALLOCATE(shpf(ishpf)%data,(mesh_size))
558            rp=>shpf(ishpf)
559            in_data=.true.
560            ndata = 0
561          end if
562 
563          value = getValue(attributes,"rc")
564          if (value == "" ) then
565            if(paw_setuploc%shape_function%gtype /="num") then
566               msg="Cannot determine rc"
567               LIBPAW_ERROR(msg)
568            end if
569          else
570            read(unit=value,fmt=*) paw_setuploc%shape_function%rc
571          end if
572 
573          value = getValue(attributes,"lamb")
574          if (value == "" ) then
575            paw_setuploc%shape_function%lamb=0
576          else
577            read(unit=value,fmt=*) paw_setuploc%shape_function%lamb
578          end if
579 
580       case ("pseudo_partial_wave")
581          ipswf=ipswf+1
582          paw_setuploc%pseudo_partial_wave(ipswf)%tread=.true.
583          value = getValue(attributes,"grid")
584          if (value == "" ) value = "unknown"
585          paw_setuploc%idgrid = trim(value)
586          paw_setuploc%pseudo_partial_wave(ipswf)%grid=trim(value)
587 
588          value = getValue(attributes,"state")
589          if (value == "" ) then
590            msg="Cannot determine pseudo_partial_wave state"
591            LIBPAW_ERROR(msg)
592          end if
593          paw_setuploc%pseudo_partial_wave(ipswf)%state=trim(value)
594 
595          do ii=1,igrid
596            if(trim(paw_setuploc%pseudo_partial_wave(ipswf)%grid)==trim(grids(ii)%id)) then
597              mesh_size=grids(ii)%iend-grids(ii)%istart+1
598            end if
599          end do
600 
601          LIBPAW_ALLOCATE(paw_setuploc%pseudo_partial_wave(ipswf)%data,(mesh_size))
602          rp=>paw_setuploc%pseudo_partial_wave(ipswf)
603          if(ipswf==paw_setuploc%valence_states%nval) ipswf=0
604          in_data=.true.
605          ndata = 0
606 
607       case ("ae_partial_wave")
608          iaewf=iaewf+1
609          paw_setuploc%ae_partial_wave(iaewf)%tread=.true.
610          value = getValue(attributes,"grid")
611          if (value == "" ) value = "unknown"
612          paw_setuploc%ae_partial_wave(iaewf)%grid=trim(value)
613 
614          value = getValue(attributes,"state")
615          if (value == "" ) then
616            LIBPAW_ERROR("Cannot determine ae_partial_wave state")
617          end if
618          paw_setuploc%ae_partial_wave(iaewf)%state=trim(value)
619 
620          do ii=1,igrid
621            if(trim(paw_setuploc%ae_partial_wave(iaewf)%grid)==trim(grids(ii)%id)) then
622              mesh_size=grids(ii)%iend-grids(ii)%istart+1
623            end if
624          end do
625 
626          LIBPAW_ALLOCATE(paw_setuploc%ae_partial_wave(iaewf)%data,(mesh_size))
627          rp=>paw_setuploc%ae_partial_wave(iaewf)
628          if(iaewf==paw_setuploc%valence_states%nval) iaewf=0
629          in_data=.true.
630          ndata = 0
631 
632       case ("projector_function")
633          iproj=iproj+1
634          paw_setuploc%projector_function(iproj)%tread=.true.
635          value = getValue(attributes,"grid")
636          if (value == "" ) value = "unknown"
637          paw_setuploc%projector_function(iproj)%grid=trim(value)
638 
639          value = getValue(attributes,"state")
640          if (value == "" ) then
641            msg="Cannot determine projector_function state"
642            LIBPAW_ERROR(msg)
643          end if
644          paw_setuploc%projector_function(iproj)%state=trim(value)
645 
646          do ii=1,igrid
647            if(trim(paw_setuploc%projector_function(iproj)%grid)==trim(grids(ii)%id)) then
648              mesh_size=grids(ii)%iend-grids(ii)%istart+1
649            end if
650          end do
651 
652          LIBPAW_ALLOCATE(paw_setuploc%projector_function(iproj)%data,(mesh_size))
653          rp=>paw_setuploc%projector_function(iproj)
654          if(iproj==paw_setuploc%valence_states%nval) iproj=0
655          in_data=.true.
656          ndata = 0
657 
658       case ("projector_fit")
659          if(.not.allocated(paw_setuploc%projector_fit)) then
660             LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%projector_fit,(paw_setuploc%valence_states%nval))
661          end if
662 
663          iprojfit=iprojfit+1
664          paw_setuploc%projector_fit(iprojfit)%tread=.true.
665          value = getValue(attributes,"state")
666          if (value == "" ) then
667            msg="Cannot determine projector_fit state"
668            LIBPAW_ERROR(msg)
669          end if
670          paw_setuploc%projector_fit(iprojfit)%state=trim(value)
671 
672          if(iprojfit==paw_setuploc%valence_states%nval) iprojfit=0
673          igauss = 0
674 
675       case ("gaussian")
676          igauss = igauss + 1
677          value = getValue(attributes,"factor")
678          if (value == "" ) then
679            msg="Cannot determine gaussian factor"
680            LIBPAW_ERROR(msg)
681          end if
682          read(value(2:100), *) paw_setuploc%projector_fit(iprojfit)%factors(1, igauss)
683          read(value(index(value, ',') + 1:100), *) paw_setuploc%projector_fit(iprojfit)%factors(2, igauss)
684          value = getValue(attributes,"exponent")
685          if (value == "" ) then
686            msg="Cannot determine gaussian exponent"
687            LIBPAW_ERROR(msg)
688          end if
689          read(value(2:100), *) paw_setuploc%projector_fit(iprojfit)%expos(1, igauss)
690          read(value(index(value, ',') + 1:100), *) paw_setuploc%projector_fit(iprojfit)%expos(2, igauss)
691 
692      case ("ae_core_density")
693          paw_setuploc%ae_core_density%tread=.true.
694           value = getValue(attributes,"grid")
695          if (value == "" ) value = "unknown"
696          paw_setuploc%ae_core_density%grid=trim(value)
697 
698          value = getValue(attributes,"state")
699          if (value == "" ) value = "unknown"
700          paw_setuploc%ae_core_density%state=trim(value)
701 
702          do ii=1,igrid
703            if(trim(paw_setuploc%ae_core_density%grid)==trim(grids(ii)%id)) then
704              mesh_size=grids(ii)%iend-grids(ii)%istart+1
705            end if
706          end do
707 
708          LIBPAW_ALLOCATE(paw_setuploc%ae_core_density%data,(mesh_size))
709          rp=>paw_setuploc%ae_core_density
710          in_data=.true.
711          ndata = 0
712 
713      case ("pseudo_core_density")
714          paw_setuploc%pseudo_core_density%tread=.true.
715           value = getValue(attributes,"grid")
716          if (value == "" ) value = "unknown"
717          paw_setuploc%pseudo_core_density%grid=trim(value)
718 
719          value = getValue(attributes,"state")
720          if (value == "" ) value = "unknown"
721          paw_setuploc%pseudo_core_density%state=trim(value)
722 
723          do ii=1,igrid
724            if(trim(paw_setuploc%pseudo_core_density%grid)==trim(grids(ii)%id)) then
725              mesh_size=grids(ii)%iend-grids(ii)%istart+1
726            end if
727          end do
728 
729          LIBPAW_ALLOCATE(paw_setuploc%pseudo_core_density%data,(mesh_size))
730          rp=>paw_setuploc%pseudo_core_density
731          in_data=.true.
732          ndata = 0
733 
734      case ("pseudo_valence_density")
735          paw_setuploc%pseudo_valence_density%tread=.true.
736           value = getValue(attributes,"grid")
737          if (value == "" ) value = "unknown"
738          paw_setuploc%pseudo_valence_density%grid=trim(value)
739 
740          value = getValue(attributes,"state")
741          if (value == "" ) value = "unknown"
742          paw_setuploc%pseudo_valence_density%state=trim(value)
743 
744          do ii=1,igrid
745            if(trim(paw_setuploc%pseudo_valence_density%grid)==trim(grids(ii)%id)) then
746              mesh_size=grids(ii)%iend-grids(ii)%istart+1
747            end if
748          end do
749 
750          LIBPAW_ALLOCATE(paw_setuploc%pseudo_valence_density%data,(mesh_size))
751          rp=>paw_setuploc%pseudo_valence_density
752          in_data=.true.
753          ndata = 0
754 
755      case ("zero_potential")
756          paw_setuploc%zero_potential%tread=.true.
757           value = getValue(attributes,"grid")
758          if (value == "" ) value = "unknown"
759          paw_setuploc%zero_potential%grid=trim(value)
760 
761          value = getValue(attributes,"state")
762          if (value == "" ) value = "unknown"
763          paw_setuploc%zero_potential%state=trim(value)
764 
765          do ii=1,igrid
766            if(trim(paw_setuploc%zero_potential%grid)==trim(grids(ii)%id)) then
767              mesh_size=grids(ii)%iend-grids(ii)%istart+1
768            end if
769          end do
770 
771          LIBPAW_ALLOCATE(paw_setuploc%zero_potential%data,(mesh_size))
772          rp=>paw_setuploc%zero_potential
773          in_data=.true.
774          ndata = 0
775 
776      case ("LDA_minus_half_potential")
777          paw_setuploc%LDA_minus_half_potential%tread=.true.
778           value = getValue(attributes,"grid")
779          if (value == "" ) value = "unknown"
780          paw_setuploc%LDA_minus_half_potential%grid=trim(value)
781 
782          value = getValue(attributes,"state")
783          if (value == "" ) value = "unknown"
784          paw_setuploc%LDA_minus_half_potential%state=trim(value)
785 
786          do ii=1,igrid
787            if(trim(paw_setuploc%LDA_minus_half_potential%grid)==trim(grids(ii)%id)) then
788              mesh_size=grids(ii)%iend-grids(ii)%istart+1
789            end if
790          end do
791 
792          LIBPAW_ALLOCATE(paw_setuploc%LDA_minus_half_potential%data,(mesh_size))
793          rp=>paw_setuploc%LDA_minus_half_potential
794          in_data=.true.
795          ndata = 0
796 
797      case ("ae_core_kinetic_energy_density")
798          paw_setuploc%ae_core_kinetic_energy_density%tread=.true.
799           value = getValue(attributes,"grid")
800          if (value == "" ) value = "unknown"
801          paw_setuploc%ae_core_kinetic_energy_density%grid=trim(value)
802 
803          value = getValue(attributes,"state")
804          if (value == "" ) value = "unknown"
805          paw_setuploc%ae_core_kinetic_energy_density%state=trim(value)
806 
807          do ii=1,igrid
808            if(trim(paw_setuploc%ae_core_kinetic_energy_density%grid)==trim(grids(ii)%id)) then
809              mesh_size=grids(ii)%iend-grids(ii)%istart+1
810            end if
811          end do
812 
813          LIBPAW_ALLOCATE(paw_setuploc%ae_core_kinetic_energy_density%data,(mesh_size))
814          rp=>paw_setuploc%ae_core_kinetic_energy_density
815          in_data=.true.
816          ndata = 0
817 
818      case ("pseudo_core_kinetic_energy_density")
819          paw_setuploc%pseudo_core_kinetic_energy_density%tread=.true.
820           value = getValue(attributes,"grid")
821          if (value == "" ) value = "unknown"
822          paw_setuploc%pseudo_core_kinetic_energy_density%grid=trim(value)
823 
824          value = getValue(attributes,"state")
825          if (value == "" ) value = "unknown"
826          paw_setuploc%pseudo_core_kinetic_energy_density%state=trim(value)
827 
828          do ii=1,igrid
829            if(trim(paw_setuploc%pseudo_core_kinetic_energy_density%grid)==trim(grids(ii)%id)) then
830              mesh_size=grids(ii)%iend-grids(ii)%istart+1
831            end if
832          end do
833 
834          LIBPAW_ALLOCATE(paw_setuploc%pseudo_core_kinetic_energy_density%data,(mesh_size))
835          rp=>paw_setuploc%pseudo_core_kinetic_energy_density
836          in_data=.true.
837          ndata = 0
838 
839      case ("kresse_joubert_local_ionic_potential")
840          paw_setuploc%kresse_joubert_local_ionic_potential%tread=.true.
841           value = getValue(attributes,"grid")
842          if (value == "" ) value = "unknown"
843          paw_setuploc%kresse_joubert_local_ionic_potential%grid=trim(value)
844 
845          value = getValue(attributes,"state")
846          if (value == "" ) value = "unknown"
847          paw_setuploc%kresse_joubert_local_ionic_potential%state=trim(value)
848 
849          do ii=1,igrid
850            if(trim(paw_setuploc%kresse_joubert_local_ionic_potential%grid)==trim(grids(ii)%id)) then
851              mesh_size=grids(ii)%iend-grids(ii)%istart+1
852            end if
853          end do
854 
855          LIBPAW_ALLOCATE(paw_setuploc%kresse_joubert_local_ionic_potential%data,(mesh_size))
856          rp=>paw_setuploc%kresse_joubert_local_ionic_potential
857          in_data=.true.
858          ndata = 0
859 
860      case ("blochl_local_ionic_potential")
861          paw_setuploc%blochl_local_ionic_potential%tread=.true.
862           value = getValue(attributes,"grid")
863          if (value == "" ) value = "unknown"
864          paw_setuploc%blochl_local_ionic_potential%grid=trim(value)
865 
866          value = getValue(attributes,"state")
867          if (value == "" ) value = "unknown"
868          paw_setuploc%blochl_local_ionic_potential%state=trim(value)
869 
870          do ii=1,igrid
871            if(trim(paw_setuploc%blochl_local_ionic_potential%grid)==trim(grids(ii)%id)) then
872              mesh_size=grids(ii)%iend-grids(ii)%istart+1
873            end if
874          end do
875 
876          LIBPAW_ALLOCATE(paw_setuploc%blochl_local_ionic_potential%data,(mesh_size))
877          rp=>paw_setuploc%blochl_local_ionic_potential
878          in_data=.true.
879          ndata = 0
880 
881     case ("kinetic_energy_differences")
882          paw_setuploc%kinetic_energy_differences%tread=.true.
883          mesh_size=paw_setuploc%valence_states%nval*paw_setuploc%valence_states%nval
884          LIBPAW_ALLOCATE(paw_setuploc%kinetic_energy_differences%data,(mesh_size))
885          rp=>paw_setuploc%kinetic_energy_differences
886          in_data=.true.
887          ndata = 0
888 
889 end select
890 
891 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.

SOURCE

 912 subroutine paw_end_element1(namespaceURI,localName,name)
 913 
 914 character(len=*),intent(in) :: namespaceURI,localName,name
 915 character(len=100) :: msg,value
 916 
 917 !Just to fool abirules 
 918  value=localName
 919  value=namespaceURI
 920  
 921 select case(name)
 922 
 923       case ("generator")
 924          in_generator = .false.
 925 
 926       case ("valence_states")
 927         in_valenceStates = .false.
 928         if(ival>50) then
 929           msg="ival>50"
 930           LIBPAW_ERROR(msg)
 931         end if
 932         if(ival>0)then
 933           LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%valence_states%state,(ival))
 934           paw_setuploc%valence_states%state(ival)%tread=.true.
 935           paw_setuploc%valence_states%nval=ival
 936           do ii=1,ival
 937             paw_setuploc%valence_states%state(ii)=valstate(ii)
 938           end do
 939         end if
 940         LIBPAW_DATATYPE_DEALLOCATE(valstate)
 941         if(.not.allocated(paw_setuploc%ae_partial_wave)) then
 942           LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%ae_partial_wave,(paw_setuploc%valence_states%nval))
 943         end if
 944         if(.not.allocated(paw_setuploc%pseudo_partial_wave)) then
 945           LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%pseudo_partial_wave,(paw_setuploc%valence_states%nval))
 946         end if
 947         if(.not.allocated(paw_setuploc%projector_function)) then
 948           LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%projector_function,(paw_setuploc%valence_states%nval))
 949         end if
 950 
 951       case ("paw_setup")
 952         if(igrid>10) then
 953           msg="igrid>10"
 954           LIBPAW_ERROR(msg)
 955         end if
 956         LIBPAW_DATATYPE_ALLOCATE(paw_setuploc%radial_grid,(igrid))
 957         paw_setuploc%radial_grid(igrid)%tread=.true.
 958         paw_setuploc%ngrid=igrid
 959         do ii=1,igrid
 960           paw_setuploc%radial_grid(ii)=grids(ii)
 961         end do
 962         LIBPAW_DATATYPE_DEALLOCATE(grids)
 963         do ii=1,igrid
 964           if(trim(paw_setuploc%shape_function%grid)==trim(paw_setuploc%radial_grid(ii)%id)) then
 965             mesh_size=paw_setuploc%radial_grid(ii)%iend-paw_setuploc%radial_grid(ii)%istart+1
 966           end if
 967         end do
 968         if(ishpf>10) then
 969           msg="ishpf>7"
 970           LIBPAW_ERROR(msg)
 971         end if
 972         LIBPAW_ALLOCATE(paw_setuploc%shape_function%data,(mesh_size,ishpf))
 973         do ii=1,ishpf
 974           paw_setuploc%shape_function%data(:,ii)=shpf(ii)%data(:)
 975           LIBPAW_DEALLOCATE(shpf(ii)%data)
 976         end do
 977         LIBPAW_DATATYPE_DEALLOCATE(shpf)
 978 
 979       case ("shape_function")
 980         in_data=.false.
 981 
 982       case ("pseudo_partial_wave")
 983         in_data=.false.
 984 
 985       case ("ae_partial_wave")
 986         in_data=.false.
 987 
 988       case ("projector_function")
 989         in_data=.false.
 990 
 991       case ("ae_core_density")
 992         in_data=.false.
 993 
 994       case ("pseudo_core_density")
 995         in_data=.false.
 996 
 997       case ("pseudo_valence_density")
 998         in_data=.false.
 999 
1000       case ("zero_potential")
1001         in_data=.false.
1002 
1003       case ("LDA_minus_half_potential")
1004         in_data=.false.
1005 
1006       case ("ae_core_kinetic_energy_density")
1007         in_data=.false.
1008 
1009       case ("pseudo_core_kinetic_energy_density")
1010         in_data=.false.
1011 
1012       case ("kresse_joubert_local_ionic_potential")
1013         in_data=.false.
1014 
1015       case ("blochl_local_ionic_potential")
1016         in_data=.false.
1017 
1018       case ("kinetic_energy_differences")
1019         in_data=.false.
1020 
1021 end select
1022 
1023 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

SOURCE

1474  subroutine paw_rdfromline(keyword,line,output,ierr)
1475 
1476 !Arguments ---------------------------------------------
1477   character(len=*), intent(in) :: keyword,line
1478   character(len=*), intent(out) :: output
1479   integer, intent(out) :: ierr
1480 !Local variables ---------------------------------------
1481   character(len=len(line)) :: temp
1482   integer :: pos,pos2
1483 
1484 ! *********************************************************************
1485 
1486  ierr=1;output=""
1487  pos=index(line,trim(keyword))
1488  if (pos>0) then
1489    temp=line(pos+len_trim(keyword):len_trim(line))
1490    pos=index(temp,char(34))
1491    if (pos>0) then
1492      pos2=index(temp(pos+1:len_trim(temp)),char(34))
1493      if (pos2>0) then
1494        output=temp(pos+1:pos+pos2-1)
1495      end if
1496    end if
1497  end if
1498 
1499  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

SOURCE

1242 subroutine paw_setup_copy(paw_setupin,paw_setupout)
1243 
1244 !Arguments ------------------------------------
1245 !scalars
1246  type(paw_setup_t),intent(in) :: paw_setupin
1247  type(paw_setup_t),intent(out) :: paw_setupout
1248 
1249 !Local variables-------------------------------
1250 !scalars
1251  integer :: ii,sz1,sz2
1252 
1253 ! *********************************************************************
1254 
1255 !scalars
1256  paw_setupout%version=paw_setupin%version
1257  paw_setupout%tread=paw_setupin%tread
1258  paw_setupout%ngrid=paw_setupin%ngrid
1259  paw_setupout%idgrid=paw_setupin%idgrid
1260  paw_setupout%optortho=paw_setupin%optortho
1261  paw_setupout%rpaw=paw_setupin%rpaw
1262  paw_setupout%ex_cc=paw_setupin%ex_cc
1263  paw_setupout%lamb_shielding=paw_setupin%lamb_shielding
1264  paw_setupout%atom%tread=paw_setupin%atom%tread
1265  paw_setupout%atom%symbol=paw_setupin%atom%symbol
1266  paw_setupout%atom%znucl=paw_setupin%atom%znucl
1267  paw_setupout%atom%zion=paw_setupin%atom%zion
1268  paw_setupout%atom%zval=paw_setupin%atom%zval
1269  paw_setupout%xc_functional%tread=paw_setupin%xc_functional%tread
1270  paw_setupout%xc_functional%functionaltype=paw_setupin%xc_functional%functionaltype
1271  paw_setupout%xc_functional%name=paw_setupin%xc_functional%name
1272  paw_setupout%generator%tread=paw_setupin%generator%tread
1273  paw_setupout%generator%gen=paw_setupin%generator%gen
1274  paw_setupout%generator%name=paw_setupin%generator%name
1275  paw_setupout%valence_states%tread=paw_setupin%valence_states%tread
1276  paw_setupout%valence_states%nval=paw_setupin%valence_states%nval
1277  paw_setupout%shape_function%tread=paw_setupin%shape_function%tread
1278  paw_setupout%shape_function%gtype=paw_setupin%shape_function%gtype
1279  paw_setupout%shape_function%grid=paw_setupin%shape_function%grid
1280  paw_setupout%shape_function%rc=paw_setupin%shape_function%rc
1281  paw_setupout%shape_function%lamb=paw_setupin%shape_function%lamb
1282  paw_setupout%ae_core_density%tread=paw_setupin%ae_core_density%tread
1283  paw_setupout%ae_core_density%grid=paw_setupin%ae_core_density%grid
1284  paw_setupout%ae_core_density%state=paw_setupin%ae_core_density%state
1285  paw_setupout%pseudo_core_density%tread=paw_setupin%pseudo_core_density%tread
1286  paw_setupout%pseudo_core_density%grid=paw_setupin%pseudo_core_density%grid
1287  paw_setupout%pseudo_core_density%state=paw_setupin%pseudo_core_density%state
1288  paw_setupout%pseudo_valence_density%tread=paw_setupin%pseudo_valence_density%tread
1289  paw_setupout%pseudo_valence_density%grid=paw_setupin%pseudo_valence_density%grid
1290  paw_setupout%pseudo_valence_density%state=paw_setupin%pseudo_valence_density%state
1291  paw_setupout%zero_potential%tread=paw_setupin%zero_potential%tread
1292  paw_setupout%zero_potential%grid=paw_setupin%zero_potential%grid
1293  paw_setupout%zero_potential%state=paw_setupin%zero_potential%state
1294  paw_setupout%LDA_minus_half_potential%tread=paw_setupin%LDA_minus_half_potential%tread
1295  paw_setupout%LDA_minus_half_potential%grid=paw_setupin%LDA_minus_half_potential%grid
1296  paw_setupout%LDA_minus_half_potential%state=paw_setupin%LDA_minus_half_potential%state
1297  paw_setupout%ae_core_kinetic_energy_density%tread=&
1298 &     paw_setupin%ae_core_kinetic_energy_density%tread
1299  paw_setupout%ae_core_kinetic_energy_density%grid=&
1300 &     paw_setupin%ae_core_kinetic_energy_density%grid
1301  paw_setupout%ae_core_kinetic_energy_density%state=&
1302 &     paw_setupin%ae_core_kinetic_energy_density%state
1303  paw_setupout%pseudo_core_kinetic_energy_density%tread=&
1304 &     paw_setupin%pseudo_core_kinetic_energy_density%tread
1305  paw_setupout%pseudo_core_kinetic_energy_density%grid=&
1306 &     paw_setupin%pseudo_core_kinetic_energy_density%grid
1307  paw_setupout%pseudo_core_kinetic_energy_density%state=&
1308 &     paw_setupin%pseudo_core_kinetic_energy_density%state
1309  paw_setupout%kresse_joubert_local_ionic_potential%tread=&
1310 &    paw_setupin%kresse_joubert_local_ionic_potential%tread
1311  paw_setupout%kresse_joubert_local_ionic_potential%grid=&
1312 &    paw_setupin%kresse_joubert_local_ionic_potential%grid
1313  paw_setupout%kresse_joubert_local_ionic_potential%state=&
1314 &    paw_setupin%kresse_joubert_local_ionic_potential%state
1315  paw_setupout%blochl_local_ionic_potential%tread=&
1316 &    paw_setupin%blochl_local_ionic_potential%tread
1317  paw_setupout%blochl_local_ionic_potential%grid=&
1318 &    paw_setupin%blochl_local_ionic_potential%grid
1319  paw_setupout%blochl_local_ionic_potential%state=&
1320 &    paw_setupin%blochl_local_ionic_potential%state
1321  paw_setupout%kinetic_energy_differences%tread=paw_setupin%kinetic_energy_differences%tread
1322  paw_setupout%kinetic_energy_differences%grid=paw_setupin%kinetic_energy_differences%grid
1323  paw_setupout%kinetic_energy_differences%state=paw_setupin%kinetic_energy_differences%state
1324  paw_setupout%exact_exchange_matrix%tread=paw_setupin%exact_exchange_matrix%tread
1325 ! allocatable arrays
1326  if (allocated(paw_setupin%shape_function%data)) then
1327    sz1=size(paw_setupin%shape_function%data,1)
1328    sz2=size(paw_setupin%shape_function%data,2)
1329    LIBPAW_ALLOCATE(paw_setupout%shape_function%data,(sz1,sz2))
1330    paw_setupout%shape_function%data=paw_setupin%shape_function%data
1331  end if
1332  if (allocated(paw_setupin%ae_core_density%data)) then
1333    sz1=size(paw_setupin%ae_core_density%data,1)
1334    LIBPAW_ALLOCATE(paw_setupout%ae_core_density%data,(sz1))
1335    paw_setupout%ae_core_density%data=paw_setupin%ae_core_density%data
1336  end if
1337  if (allocated(paw_setupin%pseudo_core_density%data)) then
1338    sz1=size(paw_setupin%pseudo_core_density%data,1)
1339    LIBPAW_ALLOCATE(paw_setupout%pseudo_core_density%data,(sz1))
1340    paw_setupout%pseudo_core_density%data=paw_setupin%pseudo_core_density%data
1341  end if
1342  if (allocated(paw_setupin%pseudo_valence_density%data)) then
1343    sz1=size(paw_setupin%pseudo_valence_density%data,1)
1344    LIBPAW_ALLOCATE(paw_setupout%pseudo_valence_density%data,(sz1))
1345    paw_setupout%pseudo_valence_density%data=paw_setupin%pseudo_valence_density%data
1346  end if
1347  if (allocated(paw_setupin%zero_potential%data)) then
1348    sz1=size(paw_setupin%zero_potential%data,1)
1349    LIBPAW_ALLOCATE(paw_setupout%zero_potential%data,(sz1))
1350    paw_setupout%zero_potential%data=paw_setupin%zero_potential%data
1351  end if
1352  if (allocated(paw_setupin%LDA_minus_half_potential%data)) then
1353    sz1=size(paw_setupin%LDA_minus_half_potential%data,1)
1354    LIBPAW_ALLOCATE(paw_setupout%LDA_minus_half_potential%data,(sz1))
1355    paw_setupout%LDA_minus_half_potential%data=paw_setupin%LDA_minus_half_potential%data
1356  end if
1357  if (allocated(paw_setupin%ae_core_kinetic_energy_density%data)) then
1358    sz1=size(paw_setupin%ae_core_kinetic_energy_density%data,1)
1359    LIBPAW_ALLOCATE(paw_setupout%ae_core_kinetic_energy_density%data,(sz1))
1360    paw_setupout%ae_core_kinetic_energy_density%data=paw_setupin%ae_core_kinetic_energy_density%data
1361  end if
1362  if (allocated(paw_setupin%pseudo_core_kinetic_energy_density%data)) then
1363    sz1=size(paw_setupin%pseudo_core_kinetic_energy_density%data,1)
1364    LIBPAW_ALLOCATE(paw_setupout%pseudo_core_kinetic_energy_density%data,(sz1))
1365    paw_setupout%pseudo_core_kinetic_energy_density%data=paw_setupin%pseudo_core_kinetic_energy_density%data
1366  end if
1367  if (allocated(paw_setupin%kresse_joubert_local_ionic_potential%data)) then
1368    sz1=size(paw_setupin%kresse_joubert_local_ionic_potential%data,1)
1369    LIBPAW_ALLOCATE(paw_setupout%kresse_joubert_local_ionic_potential%data,(sz1))
1370    paw_setupout%kresse_joubert_local_ionic_potential%data=paw_setupin%kresse_joubert_local_ionic_potential%data
1371  end if
1372  if (allocated(paw_setupin%blochl_local_ionic_potential%data)) then
1373    sz1=size(paw_setupin%blochl_local_ionic_potential%data,1)
1374    LIBPAW_ALLOCATE(paw_setupout%blochl_local_ionic_potential%data,(sz1))
1375    paw_setupout%blochl_local_ionic_potential%data=paw_setupin%blochl_local_ionic_potential%data
1376  end if
1377  if (allocated(paw_setupin%exact_exchange_matrix%data)) then
1378    sz1=size(paw_setupin%exact_exchange_matrix%data,1)
1379    LIBPAW_ALLOCATE(paw_setupout%exact_exchange_matrix%data,(sz1))
1380    paw_setupout%exact_exchange_matrix%data=paw_setupin%exact_exchange_matrix%data
1381  end if
1382  if (allocated(paw_setupin%kinetic_energy_differences%data)) then
1383    sz1=size(paw_setupin%kinetic_energy_differences%data,1)
1384    LIBPAW_ALLOCATE(paw_setupout%kinetic_energy_differences%data,(sz1))
1385    paw_setupout%kinetic_energy_differences%data=paw_setupin%kinetic_energy_differences%data
1386  end if
1387  if(allocated( paw_setupin%radial_grid)) then
1388    sz1=size(paw_setupin%radial_grid,1)
1389    LIBPAW_DATATYPE_ALLOCATE(paw_setupout%radial_grid,(sz1))
1390    paw_setupout%radial_grid=paw_setupin%radial_grid
1391  end if
1392  if(allocated(paw_setupin%valence_states%state)) then
1393    sz1=size(paw_setupin%valence_states%state,1)
1394    LIBPAW_DATATYPE_ALLOCATE(paw_setupout%valence_states%state,(sz1))
1395    paw_setupout%valence_states%state=paw_setupin%valence_states%state
1396  end if
1397 
1398  if (allocated( paw_setupin%ae_partial_wave)) then
1399    sz1=size(paw_setupin%ae_partial_wave,1)
1400    LIBPAW_DATATYPE_ALLOCATE(paw_setupout%ae_partial_wave,(sz1))
1401    do ii=1,paw_setupin%valence_states%nval
1402      paw_setupout%ae_partial_wave(ii)%tread=paw_setupin%ae_partial_wave(ii)%tread
1403      paw_setupout%ae_partial_wave(ii)%grid=paw_setupin%ae_partial_wave(ii)%grid
1404      paw_setupout%ae_partial_wave(ii)%state=paw_setupin%ae_partial_wave(ii)%state
1405      if(allocated( paw_setupin%ae_partial_wave(ii)%data)) then
1406        sz1=size(paw_setupin%ae_partial_wave(ii)%data,1)
1407        LIBPAW_ALLOCATE(paw_setupout%ae_partial_wave(ii)%data,(sz1))
1408        paw_setupout%ae_partial_wave(ii)%data=paw_setupin%ae_partial_wave(ii)%data
1409      end if
1410    end do
1411  end if 
1412  if (allocated( paw_setupin%pseudo_partial_wave)) then
1413    sz1=size(paw_setupin%pseudo_partial_wave,1)
1414    LIBPAW_DATATYPE_ALLOCATE(paw_setupout%pseudo_partial_wave,(sz1))
1415    do ii=1,paw_setupin%valence_states%nval
1416      paw_setupout%pseudo_partial_wave(ii)%tread=paw_setupin%pseudo_partial_wave(ii)%tread
1417      paw_setupout%pseudo_partial_wave(ii)%grid=paw_setupin%pseudo_partial_wave(ii)%grid
1418      paw_setupout%pseudo_partial_wave(ii)%state=paw_setupin%pseudo_partial_wave(ii)%state
1419      if(allocated( paw_setupin%pseudo_partial_wave(ii)%data)) then
1420        sz1=size(paw_setupin%pseudo_partial_wave(ii)%data,1)
1421        LIBPAW_ALLOCATE(paw_setupout%pseudo_partial_wave(ii)%data,(sz1))
1422        paw_setupout%pseudo_partial_wave(ii)%data=paw_setupin%pseudo_partial_wave(ii)%data
1423      end if
1424    end do
1425  end if 
1426   if (allocated( paw_setupin%projector_function)) then
1427    sz1=size(paw_setupin%projector_function,1)
1428    LIBPAW_DATATYPE_ALLOCATE(paw_setupout%projector_function,(sz1))
1429    do ii=1,paw_setupin%valence_states%nval
1430      paw_setupout%projector_function(ii)%tread=paw_setupin%projector_function(ii)%tread
1431      paw_setupout%projector_function(ii)%grid=paw_setupin%projector_function(ii)%grid
1432      paw_setupout%projector_function(ii)%state=paw_setupin%projector_function(ii)%state
1433      if(allocated( paw_setupin%projector_function(ii)%data)) then
1434        sz1=size(paw_setupin%projector_function(ii)%data,1)
1435        LIBPAW_ALLOCATE(paw_setupout%projector_function(ii)%data,(sz1))
1436        paw_setupout%projector_function(ii)%data=paw_setupin%projector_function(ii)%data
1437      end if
1438    end do
1439  end if 
1440   if (allocated( paw_setupin%projector_fit)) then
1441    sz1=size(paw_setupin%projector_fit,1)
1442    LIBPAW_DATATYPE_ALLOCATE(paw_setupout%projector_fit,(sz1))
1443    do ii=1,paw_setupin%valence_states%nval
1444      paw_setupout%projector_fit(ii)%tread=paw_setupin%projector_fit(ii)%tread
1445      paw_setupout%projector_fit(ii)%ngauss=paw_setupin%projector_fit(ii)%ngauss
1446      paw_setupout%projector_fit(ii)%state=paw_setupin%projector_fit(ii)%state
1447      paw_setupout%projector_fit(ii)%factors=paw_setupin%projector_fit(ii)%factors
1448      paw_setupout%projector_fit(ii)%expos=paw_setupin%projector_fit(ii)%expos
1449    end do
1450  end if 
1451 
1452 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.

SOURCE

1122 subroutine paw_setup_free(paw_setupin)
1123 
1124 !Arguments ------------------------------------
1125 !scalars
1126  type(paw_setup_t),intent(inout) :: paw_setupin
1127 
1128 !Local variables-------------------------------
1129  integer :: ii
1130 
1131 ! *********************************************************************
1132  
1133  paw_setupin%tread=.false.
1134  paw_setupin%atom%tread=.false.
1135  paw_setupin%xc_functional%tread=.false.
1136  paw_setupin%generator%tread=.false.
1137  paw_setupin%valence_states%tread=.false.
1138  paw_setupin%shape_function%tread=.false.
1139  paw_setupin%ae_core_density%tread=.false.
1140  paw_setupin%pseudo_core_density%tread=.false.
1141  paw_setupin%pseudo_valence_density%tread=.false.
1142  paw_setupin%zero_potential%tread=.false.
1143  paw_setupin%LDA_minus_half_potential%tread=.false.
1144  paw_setupin%ae_core_kinetic_energy_density%tread=.false.
1145  paw_setupin%pseudo_core_kinetic_energy_density%tread=.false.
1146  paw_setupin%kresse_joubert_local_ionic_potential%tread=.false.
1147  paw_setupin%blochl_local_ionic_potential%tread=.false.
1148  paw_setupin%kinetic_energy_differences%tread=.false.
1149  paw_setupin%exact_exchange_matrix%tread=.false.
1150 
1151  if(allocated( paw_setupin%shape_function%data)) then
1152    LIBPAW_DEALLOCATE(paw_setupin%shape_function%data)
1153  end if
1154  if(allocated( paw_setupin%ae_core_density%data)) then
1155    LIBPAW_DEALLOCATE(paw_setupin%ae_core_density%data)
1156  end if
1157  if(allocated( paw_setupin%pseudo_core_density%data)) then
1158    LIBPAW_DEALLOCATE(paw_setupin%pseudo_core_density%data)
1159  end if
1160  if(allocated( paw_setupin%pseudo_valence_density%data)) then
1161    LIBPAW_DEALLOCATE(paw_setupin%pseudo_valence_density%data)
1162  end if
1163  if(allocated( paw_setupin%zero_potential%data)) then
1164    LIBPAW_DEALLOCATE(paw_setupin%zero_potential%data)
1165  end if
1166  if(allocated( paw_setupin%LDA_minus_half_potential%data)) then
1167    LIBPAW_DEALLOCATE(paw_setupin%LDA_minus_half_potential%data)
1168  end if
1169  if(allocated( paw_setupin%ae_core_kinetic_energy_density%data)) then
1170    LIBPAW_DEALLOCATE(paw_setupin%ae_core_kinetic_energy_density%data)
1171  end if
1172  if(allocated( paw_setupin%pseudo_core_kinetic_energy_density%data)) then
1173    LIBPAW_DEALLOCATE(paw_setupin%pseudo_core_kinetic_energy_density%data)
1174  end if
1175  if(allocated( paw_setupin%kresse_joubert_local_ionic_potential%data)) then
1176    LIBPAW_DEALLOCATE(paw_setupin%kresse_joubert_local_ionic_potential%data)
1177  end if
1178  if(allocated( paw_setupin%blochl_local_ionic_potential%data)) then
1179    LIBPAW_DEALLOCATE(paw_setupin%blochl_local_ionic_potential%data)
1180  end if
1181  if(allocated( paw_setupin%kinetic_energy_differences%data)) then
1182    LIBPAW_DEALLOCATE(paw_setupin%kinetic_energy_differences%data)
1183  end if
1184  if(allocated( paw_setupin%exact_exchange_matrix%data)) then
1185    LIBPAW_DEALLOCATE(paw_setupin%exact_exchange_matrix%data)
1186  end if
1187  if (allocated( paw_setupin%ae_partial_wave)) then
1188    do ii=1,paw_setupin%valence_states%nval
1189      if(allocated( paw_setupin%ae_partial_wave(ii)%data)) then
1190        LIBPAW_DEALLOCATE(paw_setupin%ae_partial_wave(ii)%data)
1191      end if
1192    end do
1193    LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%ae_partial_wave)
1194  end if
1195  if (allocated( paw_setupin%pseudo_partial_wave)) then
1196    do ii=1,paw_setupin%valence_states%nval
1197      if(allocated( paw_setupin%pseudo_partial_wave(ii)%data)) then
1198        LIBPAW_DEALLOCATE(paw_setupin%pseudo_partial_wave(ii)%data)
1199      end if
1200    end do
1201    LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%pseudo_partial_wave)
1202  end if
1203  if (allocated( paw_setupin%projector_function)) then
1204    do ii=1,paw_setupin%valence_states%nval
1205      if(allocated( paw_setupin%projector_function(ii)%data)) then
1206        LIBPAW_DEALLOCATE(paw_setupin%projector_function(ii)%data)
1207      end if
1208    end do
1209    LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%projector_function)
1210  end if
1211  if (allocated( paw_setupin%projector_fit)) then
1212    LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%projector_fit)
1213  end if
1214  if(allocated(paw_setupin%valence_states%state)) then
1215    LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%valence_states%state)
1216  end if
1217  if(allocated( paw_setupin%radial_grid)) then
1218    LIBPAW_DATATYPE_DEALLOCATE(paw_setupin%radial_grid)
1219  end if
1220 
1221 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   real(dpxml)                  :: lamb_shielding=0.0D0
260   character(len=4)             :: idgrid
261   character(len=12)            :: optortho
262   type(atom_t)                 :: atom
263   type(xc_functional_t)        :: xc_functional
264   type(generator_t)            :: generator
265   type(valence_states_t)       :: valence_states
266   type(radial_grid_t), allocatable :: radial_grid(:)
267   type(shape_function_t)       :: shape_function
268   type(radialfunc_t)           :: ae_core_density
269   type(radialfunc_t)           :: pseudo_core_density
270   type(radialfunc_t)           :: pseudo_valence_density
271   type(radialfunc_t)           :: zero_potential
272   type(radialfunc_t)           :: LDA_minus_half_potential
273   type(radialfunc_t)           :: ae_core_kinetic_energy_density
274   type(radialfunc_t)           :: pseudo_core_kinetic_energy_density
275   type(radialfunc_t),allocatable :: ae_partial_wave(:)
276   type(radialfunc_t),allocatable :: pseudo_partial_wave(:)
277   type(radialfunc_t),allocatable :: projector_function(:)
278   type(gaussian_expansion_t),allocatable :: projector_fit(:)
279   type(radialfunc_t)           :: kresse_joubert_local_ionic_potential
280   type(radialfunc_t)           :: blochl_local_ionic_potential
281   type(radialfunc_t)           :: kinetic_energy_differences
282   type(radialfunc_t)           :: exact_exchange_matrix
283 end type paw_setup_t
284 
285  public :: paw_setup_free  ! Free memory
286  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)

SOURCE

1044 subroutine pawdata_chunk(chunk)
1045 
1046 character(len=*),intent(in) :: chunk
1047 
1048 integer :: ii,ntokens,status,last_pos
1049 logical :: in_token
1050 character(len=len(chunk))  :: str
1051 character(len=50)  :: msg
1052 character(len=1)  :: cc
1053 real(dpxml),pointer :: x(:)
1054 
1055 
1056 if (len_trim(chunk) == 0) RETURN     ! skip empty chunk
1057 
1058 if (in_data) then
1059   str = chunk ; x => rp%data
1060 
1061 ! Check the contents of the string and find the number of tokens it contains
1062 ! The standard separator is generalized whitespace (space, tab, CR, or LF)
1063   in_token=.false.;ntokens=0;last_pos=0
1064   do ii=1,len_trim(str)
1065     cc=str(ii:ii)
1066     if (in_token) then
1067       if (cc==char(9).or.cc==char(10).or.cc==char(13).or.cc==char(32)) then
1068         in_token = .false.
1069         if (cc==char(10).or.cc==char(13)) str(ii:ii) = " "
1070       else
1071         last_pos=ii
1072       end if
1073     else
1074       if (cc==char(9).or.cc==char(10).or.cc==char(13).or.cc==char(32)) then
1075         if (cc==char(10).or.cc==char(13)) str(ii:ii) = " "
1076       else
1077         in_token=.true.
1078         last_pos=ii
1079         ntokens=ntokens + 1
1080       end if
1081     end if
1082   end do
1083 
1084   if ((ndata+ntokens)>size(x)) then 
1085     msg="data array full"
1086     LIBPAW_ERROR(msg)
1087   end if
1088 
1089 ! Take the string and turn it into useful reals
1090   read(unit=str(1:last_pos),fmt=*,iostat=status) x(ndata+1:ndata+ntokens)
1091   if (status/=0) then
1092     msg="real conversion error"
1093     LIBPAW_ERROR(msg)
1094   end if
1095   ndata=ndata+ntokens
1096 
1097 end if
1098 
1099 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

 96 type, public              :: radialfunc_t
 97   logical             :: tread=.false.
 98   character(len=6)    :: grid=' '  !vz_z
 99   character(len=6)    :: state=' ' !vz_z
100   real(dpxml),allocatable :: data(:)
101 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

72 type, public :: radial_grid_t
73   logical           :: tread=.false.
74   character(len=20) :: eq
75   character(len=4)  :: id
76   real(dpxml)       :: aa
77   real(dpxml)       :: bb
78   real(dpxml)       :: dd
79   integer           :: nn
80   integer           :: istart
81   integer           :: iend
82 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

SOURCE

1901  subroutine rdpawpsxml(filename,paw_setup)
1902 
1903 !Arguments ---------------------------------------------
1904  character (len=fnlen),intent(in) :: filename
1905  type(paw_setup_t),intent(inout) :: paw_setup
1906 !Local variables ---------------------------------------
1907  integer :: funit, iaewf,ii,ipswf,iproj,ir,igrid,ival,ierr,ishpf,lmax,mesh_size,igauss,iprojfit
1908  logical :: endfile,found,endgauss
1909  character(len=100) :: msg
1910  character (len=XML_RECL) :: line,readline
1911  character (len=XML_RECL) :: strg
1912  character (len=30) :: strg1
1913  real(dp) :: rc(6)
1914  real(dp), allocatable :: shpf(:,:)
1915  type(state_t), pointer :: valstate (:)
1916  type(radial_grid_t), pointer :: grids (:)
1917 
1918 ! *************************************************************************
1919 
1920 !Open the atomicdata XML file for reading
1921  open(newunit=funit,file=filename,form='formatted',status='old', recl=XML_RECL)
1922 
1923 !Start a reading loop
1924  endfile=.false.
1925  found=.false.
1926  paw_setup%rpaw=-1.d0
1927  rc=-1.d0
1928 
1929  do while ((.not.endfile).and.(.not.found))
1930    read(funit,'(a)',err=10,end=10) readline
1931    line=adjustl(readline);goto 20
1932    10 line="";endfile=.true.
1933    20 continue
1934 
1935 !  --Read VERSION
1936    if ((line(1:10)=='<paw_setup').or.(line(1:12)=='<paw_dataset')) then
1937      paw_setup%tread=.true.
1938      igrid=0;ishpf=0
1939      LIBPAW_DATATYPE_ALLOCATE(grids,(10))
1940 
1941      call paw_rdfromline(" version",line,strg,ierr)
1942      paw_setup%version=trim(strg)
1943      cycle
1944    end if
1945 
1946 !  --Read TITLE, ATOMIC CHARGE AND CORE CHARGE
1947    if (line(1:6)=='<atom ') then
1948      paw_setup%atom%tread=.true.
1949      call paw_rdfromline(" symbol",line,strg,ierr)
1950      paw_setup%atom%symbol=trim(strg)
1951      call paw_rdfromline(" Z",line,strg,ierr)
1952      if (len(trim(strg))<=30) then
1953        strg1=trim(strg)
1954        read(unit=strg1,fmt=*) paw_setup%atom%znucl
1955      else
1956        read(unit=strg,fmt=*) paw_setup%atom%znucl
1957      end if
1958      call paw_rdfromline(" core",line,strg,ierr)
1959      if (len(trim(strg))<=30) then
1960        strg1=trim(strg)
1961        read(unit=strg1,fmt=*) paw_setup%atom%zion
1962      else
1963        read(unit=strg,fmt=*) paw_setup%atom%zion
1964      end if
1965      call paw_rdfromline(" valence",line,strg,ierr)
1966      if (len(trim(strg))<=30) then
1967        strg1=trim(strg)
1968        read(unit=strg1,fmt=*) paw_setup%atom%zval
1969      else
1970        read(unit=strg,fmt=*) paw_setup%atom%zval
1971      end if
1972      cycle
1973    end if
1974 
1975 !  --Read EXCHANGE-CORRELATION TYPE
1976    if (line(1:14)=='<xc_functional') then
1977      paw_setup%xc_functional%tread=.true.
1978      call paw_rdfromline(" type",line,strg,ierr)
1979      paw_setup%xc_functional%functionaltype = trim(strg)
1980      call paw_rdfromline(" name",line,strg,ierr)
1981      paw_setup%xc_functional%name= trim(strg)
1982      cycle
1983    end if
1984 
1985 !  --Read GENERATOR
1986    if (line(1:10)=='<generator') then
1987      paw_setup%generator%tread=.true.
1988      call paw_rdfromline(" type",line,strg,ierr)
1989      paw_setup%generator%gen  = trim(strg)
1990      call paw_rdfromline(" name",line,strg,ierr)
1991      paw_setup%generator%name= trim(strg)
1992      cycle
1993    end if
1994 
1995 !  --Read PAW RADIUS
1996    if (line(1:11)=='<PAW_radius') then
1997      call paw_rdfromline(" rpaw",line,strg,ierr)
1998      if (len(trim(strg))<=30) then
1999        strg1=trim(strg)
2000        read(unit=strg1,fmt=*) paw_setup%rpaw
2001      else
2002        read(unit=strg,fmt=*) paw_setup%rpaw
2003      end if
2004      cycle
2005    end if
2006    if (line(1:11)=='<paw_radius') then
2007      call paw_rdfromline(" rc",line,strg,ierr)
2008      if (len(trim(strg))<=30) then
2009        strg1=trim(strg)
2010        read(unit=strg1,fmt=*) paw_setup%rpaw
2011      else
2012        read(unit=strg,fmt=*) paw_setup%rpaw
2013      end if
2014      cycle
2015    end if
2016 
2017 !  --Read BASIS SIZE, ORBITALS, RC AND OCCUPATIONS/STATE IDs
2018    if (line(1:16)=='<valence_states>') then
2019      paw_setup%valence_states%tread=.true.
2020      LIBPAW_DATATYPE_ALLOCATE(valstate,(50))
2021      ival=0
2022      lmax=0
2023      do while (line(1:17)/='</valence_states>')
2024        read(funit,'(a)') readline;line=adjustl(readline)
2025        if (line(1:6)=='<state') then
2026          ival=ival+1
2027          if (ival>50) then
2028            close(funit)
2029            msg="Error in rdpawps1xml: basis size too large (>50)!"
2030            LIBPAW_ERROR(msg)
2031          end if
2032          call paw_rdfromline(" n",line,strg,ierr)
2033          if (strg == "" ) then 
2034            valstate(ival)%nn=-1    
2035          else
2036            if (len(trim(strg))<=30) then
2037              strg1=trim(strg)
2038              read(unit=strg1,fmt=*) valstate(ival)%nn
2039            else
2040              read(unit=strg,fmt=*) valstate(ival)%nn
2041            end if
2042          end if
2043          call paw_rdfromline(" l",line,strg,ierr)
2044          if (len(trim(strg))<=30) then
2045            strg1=trim(strg)
2046            read(unit=strg1,fmt=*) valstate(ival)%ll
2047          else
2048            read(unit=strg,fmt=*) valstate(ival)%ll
2049          end if
2050          if(valstate(ival)%ll>lmax) lmax=valstate(ival)%ll
2051          call paw_rdfromline(" f",line,strg,ierr)
2052          if (strg == "" ) then 
2053            valstate(ival)%ff=-1.d0
2054          else
2055            if (len(trim(strg))<=30) then
2056              strg1=trim(strg)
2057              read(unit=strg1,fmt=*) valstate(ival)%ff
2058            else
2059              read(unit=strg,fmt=*) valstate(ival)%ff
2060            end if
2061          end if
2062          call paw_rdfromline(" rc",line,strg,ierr)
2063          if (len(trim(strg))<=30) then
2064            strg1=trim(strg)
2065            read(unit=strg1,fmt=*) valstate(ival)%rc
2066          else
2067            read(unit=strg,fmt=*) valstate(ival)%rc
2068          end if
2069          call paw_rdfromline(" e",line,strg,ierr)
2070          if (len(trim(strg))<=30) then
2071            strg1=trim(strg)
2072            read(unit=strg1,fmt=*) valstate(ival)%ee
2073          else
2074            read(unit=strg,fmt=*) valstate(ival)%ee
2075          end if
2076          call paw_rdfromline(" id",line,strg,ierr)
2077          valstate(ival)%id = trim(strg)
2078        end if
2079      end do
2080      cycle
2081    end if
2082 
2083 !  --Read MESH_STEP AND NUMBER OF POINTS
2084    if (line(1:12)=='<radial_grid')then
2085      igrid=igrid+1
2086      call paw_rdfromline(" eq",line,strg,ierr)
2087      grids(igrid)%eq = trim(strg)
2088      call paw_rdfromline(" a",line,strg,ierr)
2089      if (strg == "" ) then
2090        grids(igrid)%aa=0.d0
2091      else
2092        if (len(trim(strg))<=30) then
2093          strg1=trim(strg)
2094          read(unit=strg1,fmt=*) grids(igrid)%aa
2095        else
2096          read(unit=strg,fmt=*) grids(igrid)%aa
2097        end if
2098      end if
2099      call paw_rdfromline(" n",line,strg,ierr)
2100      if (strg == "" ) then
2101        grids(igrid)%nn=0
2102      else
2103        if (len(trim(strg))<=30) then
2104          strg1=trim(strg)
2105          read(unit=strg1,fmt=*) grids(igrid)%nn
2106        else
2107          read(unit=strg,fmt=*) grids(igrid)%nn
2108        end if
2109      end if
2110      call paw_rdfromline(" d",line,strg,ierr)
2111      if (strg == "" ) then
2112        grids(igrid)%dd=0.d0
2113      else
2114        if (len(trim(strg))<=30) then
2115          strg1=trim(strg)
2116          read(unit=strg1,fmt=*) grids(igrid)%dd
2117        else
2118          read(unit=strg,fmt=*) grids(igrid)%dd
2119        end if
2120      end if
2121      call paw_rdfromline(" b",line,strg,ierr)
2122      if (strg == "" ) then
2123        grids(igrid)%bb=0.d0
2124      else
2125        if (len(trim(strg))<=30) then
2126          strg1=trim(strg)
2127          read(unit=strg1,fmt=*) grids(igrid)%bb
2128        else
2129          read(unit=strg,fmt=*) grids(igrid)%bb
2130        end if
2131      end if
2132      call paw_rdfromline("istart",line,strg,ierr)
2133      if (len(trim(strg))<=30) then
2134        strg1=trim(strg)
2135        read(unit=strg1,fmt=*) grids(igrid)%istart
2136      else
2137        read(unit=strg,fmt=*) grids(igrid)%istart
2138      end if
2139      call paw_rdfromline("iend",line,strg,ierr)
2140      if (len(trim(strg))<=30) then
2141        strg1=trim(strg)
2142        read(unit=strg1,fmt=*) grids(igrid)%iend
2143      else
2144        read(unit=strg,fmt=*) grids(igrid)%iend
2145      end if
2146      call paw_rdfromline(" id",line,strg,ierr)
2147      grids(igrid)%id = trim(strg)
2148      if(igrid>10)then
2149        close(funit)
2150        msg="igrid>10"
2151        LIBPAW_ERROR(msg)
2152      end if
2153      cycle
2154    end if
2155 
2156 !  --Read SHAPE TYPE
2157    if (line(1:15)=='<shape_function') then
2158      paw_setup%shape_function%tread=.true.
2159      call paw_rdfromline(" type",line,strg,ierr)
2160      paw_setup%shape_function%gtype = trim(strg)
2161      call paw_rdfromline(" rc",line,strg,ierr)
2162      if (strg /= "" ) then
2163        if (len(trim(strg))<=30) then
2164          strg1=trim(strg)
2165          read(unit=strg1,fmt=*) paw_setup%shape_function%rc
2166        else
2167          read(unit=strg,fmt=*) paw_setup%shape_function%rc
2168        end if
2169      end if
2170      call paw_rdfromline(" lamb",line,strg,ierr)
2171      if (strg == "" ) then
2172        paw_setup%shape_function%lamb=0
2173      else
2174        if (len(trim(strg))<=30) then
2175          strg1=trim(strg)
2176          read(unit=strg1,fmt=*) paw_setup%shape_function%lamb
2177        else
2178          read(unit=strg,fmt=*) paw_setup%shape_function%lamb
2179        end if
2180      end if
2181      found=paw_setup%shape_function%tread
2182      call paw_rdfromline("grid",line,strg,ierr)
2183      paw_setup%shape_function%grid=trim(strg)
2184      if (strg /= "" ) then
2185        paw_setup%shape_function%gtype ="num" 
2186        do ii=1,igrid
2187          if(trim(paw_setup%shape_function%grid)==trim(grids(ii)%id)) then
2188            mesh_size=grids(ii)%iend-grids(ii)%istart+1
2189            exit
2190          end if
2191        end do
2192        if(.not.allocated(shpf)) then
2193          LIBPAW_ALLOCATE(shpf,(mesh_size,7))
2194        end if
2195        ishpf=ishpf+1
2196        read(funit,*) (shpf(ir,ishpf),ir=1,mesh_size)
2197        call paw_rdfromline(" l",line,strg,ierr)
2198        if (strg /= "" ) then
2199          found=.false.
2200          if(paw_setup%valence_states%tread) then
2201            if(ishpf==2*lmax+1) found=.true.
2202          else
2203            write(msg,'(a,a,a)')"the grids and the states must be read before the shapefunction",ch10,&
2204 &           "Action: Modify your XML PAW data file"
2205            LIBPAW_ERROR(msg)
2206          end if
2207        end if
2208      end if
2209      cycle
2210    end if
2211 
2212 !  End of reading loop
2213  end do
2214 
2215  if(igrid==0.or.ival==0) then
2216    write(msg,'(a,a,a)')"the grids and the states must be read before the shapefunction",ch10,&
2217 &   "Action: Modify your XML PAW data file"
2218    LIBPAW_ERROR(msg)
2219  end if
2220  if(ishpf>0)then
2221    LIBPAW_ALLOCATE(paw_setup%shape_function%data,(mesh_size,ishpf))
2222    do ii=1,ishpf
2223      paw_setup%shape_function%data(:,ii)=shpf(:,ii)
2224    end do
2225    LIBPAW_DEALLOCATE(shpf)
2226  end if
2227 
2228  if(ival>0)then
2229    LIBPAW_DATATYPE_ALLOCATE(paw_setup%valence_states%state,(ival))
2230    paw_setup%valence_states%state(ival)%tread=.true.
2231    paw_setup%valence_states%nval=ival
2232    do ii=1,ival
2233      paw_setup%valence_states%state(ii)=valstate(ii)
2234    end do
2235  end if
2236  LIBPAW_DATATYPE_DEALLOCATE(valstate)
2237  if(.not.allocated(paw_setup%ae_partial_wave)) then
2238    LIBPAW_DATATYPE_ALLOCATE(paw_setup%ae_partial_wave,(paw_setup%valence_states%nval))
2239  end if
2240  if(.not.allocated(paw_setup%pseudo_partial_wave)) then
2241    LIBPAW_DATATYPE_ALLOCATE(paw_setup%pseudo_partial_wave,(paw_setup%valence_states%nval))
2242  end if
2243  if(.not.allocated(paw_setup%projector_function)) then
2244    LIBPAW_DATATYPE_ALLOCATE(paw_setup%projector_function,(paw_setup%valence_states%nval))
2245  end if
2246 
2247  LIBPAW_DATATYPE_ALLOCATE(paw_setup%radial_grid,(igrid))
2248  paw_setup%radial_grid(igrid)%tread=.true.
2249  paw_setup%ngrid=igrid
2250  do ii=1,igrid
2251    paw_setup%radial_grid(ii)=grids(ii)
2252  end do
2253  LIBPAW_DATATYPE_DEALLOCATE(grids)
2254 
2255 !Start a reading loop
2256  ipswf=0;iaewf=0;iproj=0;iprojfit=0
2257  endfile=.false.
2258  do while (.not.endfile)
2259    read(funit,'(a)',err=11,end=11) readline
2260    line=adjustl(readline);goto 21
2261    11 line="";endfile=.true.
2262    21 continue
2263 
2264 !  --Read core density CORE_DENSITY
2265    if (line(1:16)=='<ae_core_density') then
2266      paw_setup%ae_core_density%tread=.true.
2267      call paw_rdfromline(" grid",line,strg,ierr)
2268      if (strg == "" ) strg = "unknown"
2269      paw_setup%ae_core_density%grid=trim(strg)
2270      do ii=1,paw_setup%ngrid
2271        if(trim(paw_setup%ae_core_density%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2272          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2273          exit
2274        end if
2275      end do
2276      call paw_rdfromline(" rc",line,strg,ierr)
2277      if (strg /= "" ) then
2278        if (len(trim(strg))<=30) then
2279          strg1=trim(strg)
2280          read(unit=strg1,fmt=*) rc(1)
2281        else
2282          read(unit=strg,fmt=*) rc(1)
2283        end if
2284      end if
2285      LIBPAW_ALLOCATE(paw_setup%ae_core_density%data,(mesh_size))
2286      !MGNAG v7[62]
2287      ! Runtime Error: m_pawxmlps_cpp.f90, line 1657: 
2288      ! Record too long for input bufferProgram terminated by I/O error on unit 9 
2289      ! (File="/home/buildbot/ABINIT_OD/petrus_nag/gmatteo_7.7.1-training/tests/Psps_for_tests/Al.LDA",Formatted,Sequential)
2290      read(funit,*) (paw_setup%ae_core_density%data(ir),ir=1,mesh_size)
2291      cycle
2292    end if
2293 
2294 !  --Read pseudized core density CORETAIL_DENSITY
2295    if (line(1:20)=='<pseudo_core_density') then
2296      paw_setup%pseudo_core_density%tread=.true.
2297      call paw_rdfromline(" grid",line,strg,ierr)
2298      if (strg == "" ) strg = "unknown"
2299      paw_setup%pseudo_core_density%grid=trim(strg)
2300      do ii=1,paw_setup%ngrid
2301        if(trim(paw_setup%pseudo_core_density%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2302          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2303          exit
2304        end if
2305      end do
2306      call paw_rdfromline(" rc",line,strg,ierr)
2307      if (strg /= "" ) then
2308        if (len(trim(strg))<=30) then
2309          strg1=trim(strg)
2310          read(unit=strg1,fmt=*) rc(2)
2311        else
2312          read(unit=strg,fmt=*) rc(2)
2313        end if
2314      end if
2315      LIBPAW_ALLOCATE(paw_setup%pseudo_core_density%data,(mesh_size))
2316      read(funit,*) (paw_setup%pseudo_core_density%data(ir),ir=1,mesh_size)
2317      cycle
2318    end if
2319 
2320 !  --Read core density CORE_DENSITY
2321    if (line(1:31)=='<ae_core_kinetic_energy_density') then
2322      paw_setup%ae_core_kinetic_energy_density%tread=.true.
2323      call paw_rdfromline(" grid",line,strg,ierr)
2324      if (strg == "" ) strg = "unknown"
2325      paw_setup%ae_core_kinetic_energy_density%grid=trim(strg)
2326      do ii=1,paw_setup%ngrid
2327        if(trim(paw_setup%ae_core_kinetic_energy_density%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2328          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2329          exit
2330        end if
2331      end do
2332      call paw_rdfromline(" rc",line,strg,ierr)
2333      if (strg /= "" ) then
2334        if (len(trim(strg))<=30) then
2335          strg1=trim(strg)
2336          read(unit=strg1,fmt=*) rc(1)
2337        else
2338          read(unit=strg,fmt=*) rc(1)
2339        end if
2340      end if
2341      LIBPAW_ALLOCATE(paw_setup%ae_core_kinetic_energy_density%data,(mesh_size))
2342      !MGNAG v7[62]
2343      ! Runtime Error: m_pawxmlps_cpp.f90, line 1657: 
2344      ! Record too long for input bufferProgram terminated by I/O error on unit 9 
2345      ! (File="/home/buildbot/ABINIT_OD/petrus_nag/gmatteo_7.7.1-training/tests/Psps_for_tests/Al.LDA",Formatted,Sequential)
2346      read(funit,*) (paw_setup%ae_core_kinetic_energy_density%data(ir),ir=1,mesh_size)
2347      cycle
2348    end if
2349 
2350 !  --Read pseudized core density CORETAIL_DENSITY
2351    if (line(1:35)=='<pseudo_core_kinetic_energy_density') then
2352      paw_setup%pseudo_core_kinetic_energy_density%tread=.true.
2353      call paw_rdfromline(" grid",line,strg,ierr)
2354      if (strg == "" ) strg = "unknown"
2355      paw_setup%pseudo_core_kinetic_energy_density%grid=trim(strg)
2356      do ii=1,paw_setup%ngrid
2357        if(trim(paw_setup%pseudo_core_kinetic_energy_density%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2358          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2359          exit
2360        end if
2361      end do
2362      call paw_rdfromline(" rc",line,strg,ierr)
2363      if (strg /= "" ) then
2364        if (len(trim(strg))<=30) then
2365          strg1=trim(strg)
2366          read(unit=strg1,fmt=*) rc(2)
2367        else
2368          read(unit=strg,fmt=*) rc(2)
2369        end if
2370      end if
2371      LIBPAW_ALLOCATE(paw_setup%pseudo_core_kinetic_energy_density%data,(mesh_size))
2372      read(funit,*) (paw_setup%pseudo_core_kinetic_energy_density%data(ir),ir=1,mesh_size)
2373      cycle
2374    end if
2375 
2376 !  --Read pseudized valence density PSEUDO_VALENCE_DENSITY
2377    if (line(1:23)=='<pseudo_valence_density') then
2378      paw_setup%pseudo_valence_density%tread=.true.
2379      call paw_rdfromline(" grid",line,strg,ierr)
2380      if (strg == "" ) strg = "unknown"
2381      paw_setup%pseudo_valence_density%grid=trim(strg)
2382      do ii=1,paw_setup%ngrid
2383        if(trim(paw_setup%pseudo_valence_density%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2384          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2385          exit
2386        end if
2387      end do
2388      call paw_rdfromline(" rc",line,strg,ierr)
2389      if (strg /= "" ) then
2390        if (len(trim(strg))<=30) then
2391          strg1=trim(strg)
2392          read(unit=strg1,fmt=*) rc(3)
2393        else
2394          read(unit=strg,fmt=*) rc(3)
2395        end if
2396      end if
2397      LIBPAW_ALLOCATE(paw_setup%pseudo_valence_density%data,(mesh_size))
2398      read(funit,*) (paw_setup%pseudo_valence_density%data(ir),ir=1,mesh_size)
2399      cycle
2400    end if
2401 
2402 !  --Read Vbare potential VLOCFUN
2403    if (line(1:15)=='<zero_potential') then
2404      paw_setup%zero_potential%tread=.true.
2405      call paw_rdfromline(" grid",line,strg,ierr)
2406      if (strg == "" ) strg = "unknown"
2407      paw_setup%zero_potential%grid=trim(strg)
2408      do ii=1,paw_setup%ngrid
2409        if(trim(paw_setup%zero_potential%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2410          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2411          exit
2412        end if
2413      end do
2414      call paw_rdfromline(" rc",line,strg,ierr)
2415      if (strg /= "" ) then
2416        if (len(trim(strg))<=30) then
2417          strg1=trim(strg)
2418          read(unit=strg1,fmt=*) rc(4)
2419        else
2420          read(unit=strg,fmt=*) rc(4)
2421        end if
2422      end if
2423      LIBPAW_ALLOCATE(paw_setup%zero_potential%data,(mesh_size))
2424      read(funit,*) (paw_setup%zero_potential%data(ir),ir=1,mesh_size)
2425      cycle
2426    end if
2427 
2428 !  --Read external potential
2429    if (line(1:25)=='<LDA_minus_half_potential') then
2430      paw_setup%LDA_minus_half_potential%tread=.true.
2431      call paw_rdfromline(" grid",line,strg,ierr)
2432      if (strg == "" ) strg = "unknown"
2433      paw_setup%LDA_minus_half_potential%grid=trim(strg)
2434      do ii=1,paw_setup%ngrid
2435        if(trim(paw_setup%LDA_minus_half_potential%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2436          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2437          exit
2438        end if
2439      end do
2440      call paw_rdfromline(" rc",line,strg,ierr)
2441      if (strg /= "" ) then
2442        if (len(trim(strg))<=30) then
2443          strg1=trim(strg)
2444          read(unit=strg1,fmt=*) rc(4)
2445        else
2446          read(unit=strg,fmt=*) rc(4)
2447        end if
2448      end if
2449      LIBPAW_ALLOCATE(paw_setup%LDA_minus_half_potential%data,(mesh_size))
2450      read(funit,*) (paw_setup%LDA_minus_half_potential%data(ir),ir=1,mesh_size)
2451      cycle
2452    end if
2453 
2454 !  --Read Vloc for Abinit potential VLOC_ION
2455    if (line(1:37)=='<kresse_joubert_local_ionic_potential') then
2456      paw_setup%kresse_joubert_local_ionic_potential%tread=.true.
2457      call paw_rdfromline(" grid",line,strg,ierr)
2458      if (strg == "" ) strg = "unknown"
2459      paw_setup%kresse_joubert_local_ionic_potential%grid=trim(strg)
2460      do ii=1,paw_setup%ngrid
2461        if(trim(paw_setup%kresse_joubert_local_ionic_potential%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2462          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2463          exit
2464        end if
2465      end do
2466      call paw_rdfromline(" rc",line,strg,ierr)
2467      if (strg /= "" ) then
2468        if (len(trim(strg))<=30) then
2469          strg1=trim(strg)
2470          read(unit=strg1,fmt=*) rc(5)
2471        else
2472          read(unit=strg,fmt=*) rc(5)
2473        end if
2474      end if
2475      LIBPAW_ALLOCATE(paw_setup%kresse_joubert_local_ionic_potential%data,(mesh_size))
2476      read(funit,*) (paw_setup%kresse_joubert_local_ionic_potential%data(ir),ir=1,mesh_size)
2477      cycle
2478    end if
2479    if (line(1:29)=='<blochl_local_ionic_potential') then
2480      paw_setup%blochl_local_ionic_potential%tread=.true.
2481      call paw_rdfromline(" grid",line,strg,ierr)
2482      if (strg == "" ) strg = "unknown"
2483      paw_setup%blochl_local_ionic_potential%grid=trim(strg)
2484      do ii=1,paw_setup%ngrid
2485        if(trim(paw_setup%blochl_local_ionic_potential%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2486          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2487          exit
2488        end if
2489      end do
2490      call paw_rdfromline(" rc",line,strg,ierr)
2491      if (strg /= "" ) then
2492        if (len(trim(strg))<=30) then
2493          strg1=trim(strg)
2494          read(unit=strg1,fmt=*) rc(6)
2495        else
2496          read(unit=strg,fmt=*) rc(6)
2497        end if
2498      end if
2499      LIBPAW_ALLOCATE(paw_setup%blochl_local_ionic_potential%data,(mesh_size))
2500      read(funit,*) (paw_setup%blochl_local_ionic_potential%data(ir),ir=1,mesh_size)
2501      cycle
2502    end if
2503 
2504 !  --Read WAVE FUNCTIONS PHI
2505    if (line(1:16)=='<ae_partial_wave') then
2506      iaewf=iaewf+1
2507      paw_setup%ae_partial_wave(iaewf)%tread=.true.
2508      call paw_rdfromline(" grid",line,strg,ierr)
2509      if (strg == "" ) strg = "unknown"
2510      paw_setup%ae_partial_wave(iaewf)%grid=trim(strg)
2511      call paw_rdfromline(" state",line,strg,ierr)
2512      paw_setup%ae_partial_wave(iaewf)%state=trim(strg)
2513      do ii=1,paw_setup%ngrid
2514        if(trim(paw_setup%ae_partial_wave(iaewf)%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2515          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2516          exit
2517        end if
2518      end do
2519      LIBPAW_ALLOCATE(paw_setup%ae_partial_wave(iaewf)%data,(mesh_size))
2520      read(funit,*) (paw_setup%ae_partial_wave(iaewf)%data(ir),ir=1,mesh_size)
2521      cycle
2522    end if
2523 
2524 !  --Read PSEUDO WAVE FUNCTIONS TPHI
2525    if (line(1:20)=='<pseudo_partial_wave') then
2526      ipswf=ipswf+1
2527      paw_setup%pseudo_partial_wave(ipswf)%tread=.true.
2528      call paw_rdfromline(" grid",line,strg,ierr)
2529      if (strg == "" ) strg = "unknown"
2530      paw_setup%idgrid = trim(strg)
2531      paw_setup%pseudo_partial_wave(ipswf)%grid=trim(strg)
2532      call paw_rdfromline(" state",line,strg,ierr)
2533      paw_setup%pseudo_partial_wave(ipswf)%state=trim(strg)
2534      do ii=1,paw_setup%ngrid
2535        if(trim(paw_setup%pseudo_partial_wave(ipswf)%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2536          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2537          exit
2538        end if
2539      end do
2540      LIBPAW_ALLOCATE(paw_setup%pseudo_partial_wave(ipswf)%data,(mesh_size))
2541      read(funit,*) (paw_setup%pseudo_partial_wave(ipswf)%data(ir),ir=1,mesh_size)
2542      cycle
2543    end if
2544 
2545 !  --Read PROJECTORS TPROJ
2546    if (line(1:19)=='<projector_function') then
2547      iproj=iproj+1
2548      paw_setup%projector_function(iproj)%tread=.true.
2549      call paw_rdfromline(" grid",line,strg,ierr)
2550      if (strg == "" ) strg = "unknown"
2551      paw_setup%projector_function(iproj)%grid=trim(strg)
2552      call paw_rdfromline(" state",line,strg,ierr)
2553      paw_setup%projector_function(iproj)%state=trim(strg)
2554      do ii=1,paw_setup%ngrid
2555        if(trim(paw_setup%projector_function(iproj)%grid)==trim(paw_setup%radial_grid(ii)%id)) then
2556          mesh_size=paw_setup%radial_grid(ii)%iend-paw_setup%radial_grid(ii)%istart+1
2557          exit
2558        end if
2559      end do
2560      LIBPAW_ALLOCATE(paw_setup%projector_function(iproj)%data,(mesh_size))
2561      read(funit,*) (paw_setup%projector_function(iproj)%data(ir),ir=1,mesh_size)
2562      cycle
2563    end if
2564 
2565 !  --Read PROJECTORS TPROJ as gaussian representations
2566    if (line(1:14)=='<projector_fit') then
2567      if(.not.allocated(paw_setup%projector_fit)) then
2568         LIBPAW_DATATYPE_ALLOCATE(paw_setup%projector_fit,(paw_setup%valence_states%nval))
2569      end if
2570      iprojfit=iprojfit+1
2571      paw_setup%projector_fit(iprojfit)%tread=.true.
2572      call paw_rdfromline(" state",line,strg,ierr)
2573      paw_setup%projector_fit(iprojfit)%state=trim(strg)
2574      igauss = 0
2575      endgauss = .false.
2576      do while(.not. endgauss)
2577         read(funit,'(a)',err=12,end=12) readline
2578         line=adjustl(readline);goto 22
2579 12      line="";endgauss=.true.
2580 22      continue
2581         endgauss = (line(1:15)=='</projector_fit')
2582         if (line(1:9)=='<gaussian') then
2583            igauss = igauss + 1
2584            call paw_rdfromline(" factor",line,strg,ierr)
2585            read(strg(index(strg, '{') + 1:index(strg, ',') - 1), *) &
2586                 & paw_setup%projector_fit(iprojfit)%factors(1, igauss)
2587            read(strg(index(strg, ',') + 1:index(strg, '}') - 1), *) &
2588                 & paw_setup%projector_fit(iprojfit)%factors(2, igauss)
2589            call paw_rdfromline(" exponent",line,strg,ierr)
2590            read(strg(index(strg, '{') + 1:index(strg, ',') - 1), *) &
2591                 & paw_setup%projector_fit(iprojfit)%expos(1, igauss)
2592            read(strg(index(strg, ',') + 1:index(strg, '}') - 1), *) &
2593                 & paw_setup%projector_fit(iprojfit)%expos(2, igauss)
2594         end if
2595      end do
2596      paw_setup%projector_fit(iprojfit)%ngauss = igauss
2597      cycle
2598    end if
2599 
2600 !  --Read Kinetic term KINETIC_ENERGY_MATRIX
2601    if (line(1:28)=='<kinetic_energy_differences>') then
2602      paw_setup%kinetic_energy_differences%tread=.true.
2603      mesh_size=paw_setup%valence_states%nval*paw_setup%valence_states%nval
2604      LIBPAW_ALLOCATE(paw_setup%kinetic_energy_differences%data,(mesh_size))
2605      read(funit,*) (paw_setup%kinetic_energy_differences%data(ir),ir=1,mesh_size)
2606      cycle
2607    end if
2608 
2609 !  --Read Exact exchange term EXACT_EXCHANGE_X_MATRIX
2610    if (line(1:25)=='<exact_exchange_X_matrix>') then
2611      paw_setup%exact_exchange_matrix%tread=.true.
2612      mesh_size=paw_setup%valence_states%nval*paw_setup%valence_states%nval
2613      LIBPAW_ALLOCATE(paw_setup%exact_exchange_matrix%data,(mesh_size))
2614      read(funit,*) (paw_setup%exact_exchange_matrix%data(ir),ir=1,mesh_size)
2615      cycle
2616    end if
2617 
2618 !  --Read Exact exchange core-core energy
2619    if (line(1:25)=='<exact_exchange core-core') then
2620      call paw_rdfromline(" core-core",line,strg,ierr)
2621      if (len(trim(strg))<=30) then
2622        strg1=trim(strg)
2623        read(unit=strg1,fmt=*) paw_setup%ex_cc
2624      else
2625        read(unit=strg,fmt=*) paw_setup%ex_cc
2626      end if
2627      cycle
2628    end if
2629 
2630    !  --Read Lamb shielding
2631    if (line(1:25)=='<lamb_shielding shielding') then
2632      call paw_rdfromline(" shielding",line,strg,ierr)
2633      if (len(trim(strg))<=30) then
2634        strg1=trim(strg)
2635        read(unit=strg1,fmt=*) paw_setup%lamb_shielding
2636      else
2637        read(unit=strg,fmt=*) paw_setup%lamb_shielding
2638      end if
2639      cycle
2640    end if
2641 
2642 
2643 !  --Read orthogonalisation scheme
2644    if (line(1:18)=='<orthogonalisation') then
2645      call paw_rdfromline(" scheme",line,strg,ierr)
2646      if (len(trim(strg))<=30) then
2647        strg1=trim(strg)
2648        read(unit=strg1,fmt=*) paw_setup%optortho
2649      else
2650        read(unit=strg,fmt=*) paw_setup%optortho
2651      end if
2652      cycle
2653    end if
2654 
2655 !  --Read the Atompaw input file
2656    ir=0
2657    if ((line(1:13)=='<!-- Program:').and.(ir==1)) then
2658      msg=" "
2659      do while ((msg(1:9)/=' Program:').and.(msg(1:8)/='Program:'))
2660        read(funit,'(a)') msg
2661        write(ab_out,'(a)') trim(msg)
2662      end do   
2663      cycle
2664    end if 
2665 
2666 !  End of reading loop
2667  end do
2668  if(paw_setup%rpaw<0.d0) paw_setup%rpaw=maxval(rc)
2669 !Close the XML atomicdata file
2670  close(funit)
2671 
2672 !Test flags: is anything OK ?
2673  found=paw_setup%atom%tread.and.paw_setup%valence_states%tread.and.&
2674 & paw_setup%xc_functional%tread.and.paw_setup%shape_function%tread
2675 
2676  if (.not.paw_setup%atom%tread) then
2677    msg="ATOM SYMBOL not found !"
2678    LIBPAW_WARNING(msg)
2679  end if
2680  if (.not.paw_setup%valence_states%tread) then
2681    msg="VALENCE STATES not found!"
2682    LIBPAW_WARNING(msg)
2683  end if
2684  if (.not.paw_setup%xc_functional%tread) then
2685    msg="EXCHANGE/CORRELATION not found !"
2686    LIBPAW_WARNING(msg)
2687  end if
2688  if (.not.paw_setup%shape_function%tread) then
2689    msg="SHAPE FUNCTION TYPE not found !"
2690    LIBPAW_WARNING(msg)
2691  end if
2692 
2693  if (.not.found) then
2694    msg="Aborting now"
2695    LIBPAW_ERROR(msg)
2696  end if
2697  
2698  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

SOURCE

2719  subroutine rdpawpsxml_core(energy_cor,filename, &
2720 &           lcor,ncor,nphicor,pawrad,phi_cor,kappacor)
2721 
2722 !Arguments ---------------------------------------------
2723  character (len=*),intent(in) :: filename
2724  integer,intent(out) :: nphicor
2725 !arrays
2726  integer,allocatable,intent(inout) :: lcor(:),ncor(:)
2727  integer,allocatable,intent(inout),optional :: kappacor(:)
2728  real(dp),allocatable,intent(inout) :: phi_cor(:,:),energy_cor(:)
2729  type(pawrad_type),intent(in) :: pawrad
2730 
2731 
2732 !Local variables ---------------------------------------
2733  integer :: funit,iaewf,ii,imeshae,imsh,ir,igrid,icor,ierr,maxmeshz,mesh_size,nmesh,shft
2734  logical :: endfile,found,tread,diracrel
2735  real(dp) :: yp1,ypn
2736  character(len=100) :: msg,version
2737  character (len=XML_RECL) :: line,readline
2738  character (len=XML_RECL) :: strg
2739  character (len=30) :: strg1
2740  integer,allocatable :: mesh_shift(:)
2741  real(dp),allocatable :: work(:),phitmp(:,:)
2742  character (len=20), allocatable :: gridwf(:),statewf(:)
2743  type(state_t),allocatable   :: corestate (:)
2744  type(radial_grid_t),allocatable   :: grids (:)
2745  type(pawrad_type),allocatable :: radmesh(:)
2746 
2747 ! *************************************************************************
2748 
2749 !Open the atomicdata XML file for reading
2750  funit=100
2751  open(unit=funit,file=filename,form='formatted',status='old', recl=XML_RECL)
2752 
2753 !Start a reading loop
2754  endfile=.false.
2755  found=.false.
2756  diracrel=.false.
2757 
2758  do while ((.not.endfile).and.(.not.found))
2759    read(funit,'(a)',err=10,end=10) readline
2760    line=adjustl(readline);goto 20
2761    10 line="";endfile=.true.
2762    20 continue
2763 
2764 !  --Read VERSION
2765    if (line(1:10)=='<paw_setup') then
2766      tread=.true.
2767      igrid=0
2768      LIBPAW_DATATYPE_ALLOCATE(grids,(10))
2769 
2770      call paw_rdfromline(" version",line,strg,ierr)
2771      version=trim(strg)
2772      cycle
2773    end if
2774 
2775 !  --Read GENERATOR
2776    if (line(1:10)=='<generator') then
2777      tread=.true.
2778      call paw_rdfromline(" type",line,strg,ierr)
2779      if(strg=="dirac-relativistic") then
2780        diracrel=.true.
2781      else
2782        if(present(kappacor)) then
2783          write(msg,'(a)') 'Error in pawpsp_read_core: To read kappa a diracrelativistic corewf file has to be provided!'
2784          LIBPAW_ERROR(msg)
2785        endif
2786      endif
2787      cycle
2788    end if
2789 
2790 !  --Read BASIS SIZE, ORBITALS, RC AND OCCUPATIONS/STATE IDs
2791    if (line(1:13)=='<core_states>') then
2792      tread=.true.
2793      LIBPAW_DATATYPE_ALLOCATE(corestate,(50))
2794      icor=0
2795      do while (line(1:14)/='</core_states>')
2796        read(funit,'(a)') readline;line=adjustl(readline)
2797        if (line(1:6)=='<state') then
2798          icor=icor+1
2799          if (icor>50) then
2800            close(funit)
2801            msg="basis size too large (>50)!"
2802            LIBPAW_ERROR(msg)
2803          end if
2804          call paw_rdfromline(" n",line,strg,ierr)
2805          if (strg == "" ) then
2806            corestate(icor)%nn=-1
2807          else
2808            if (len(trim(strg))<=30) then
2809              strg1=trim(strg)
2810              read(unit=strg1,fmt=*) corestate(icor)%nn
2811            else
2812              read(unit=strg,fmt=*) corestate(icor)%nn
2813            end if
2814          end if
2815          call paw_rdfromline(" l",line,strg,ierr)
2816          if (len(trim(strg))<=30) then
2817            strg1=trim(strg)
2818            read(unit=strg1,fmt=*) corestate(icor)%ll
2819          else
2820            read(unit=strg,fmt=*) corestate(icor)%ll
2821          end if
2822          if(diracrel) then!does not work if xml file is in the wrong order, which is bad xml, alternatives?
2823            call paw_rdfromline(" kappa",line,strg,ierr)
2824            if(present(kappacor)) then
2825              if (len(trim(strg))<=30) then
2826                strg1=trim(strg)
2827                read(unit=strg1,fmt=*) corestate(icor)%kk
2828              else
2829                read(unit=strg,fmt=*) corestate(icor)%kk
2830              end if
2831            endif
2832          endif
2833          call paw_rdfromline(" f",line,strg,ierr)
2834          if (strg == "" ) then
2835            corestate(icor)%ff=-1.d0
2836          else
2837            if (len(trim(strg))<=30) then
2838              strg1=trim(strg)
2839              read(unit=strg1,fmt=*) corestate(icor)%ff
2840            else
2841              read(unit=strg,fmt=*) corestate(icor)%ff
2842            end if
2843          end if
2844          call paw_rdfromline(" rc",line,strg,ierr)
2845          if (strg == "" ) then
2846            corestate(icor)%rc=-1
2847          else
2848            if (len(trim(strg))<=30) then
2849              strg1=trim(strg)
2850              read(unit=strg1,fmt=*) corestate(icor)%rc
2851            else
2852              read(unit=strg,fmt=*) corestate(icor)%rc
2853            end if
2854          end if
2855          call paw_rdfromline(" e",line,strg,ierr)
2856          if (len(trim(strg))<=30) then
2857            strg1=trim(strg)
2858            read(unit=strg1,fmt=*) corestate(icor)%ee
2859          else
2860            read(unit=strg,fmt=*) corestate(icor)%ee
2861          end if
2862          call paw_rdfromline(" id",line,strg,ierr)
2863          corestate(icor)%id = trim(strg)
2864        end if
2865      end do
2866      cycle
2867    end if
2868 
2869 !  --Read MESH_STEP AND NUMBER OF POINTS
2870    if (line(1:12)=='<radial_grid')then
2871      igrid=igrid+1
2872      call paw_rdfromline(" eq",line,strg,ierr)
2873      grids(igrid)%eq = trim(strg)
2874      call paw_rdfromline(" a",line,strg,ierr)
2875      if (strg == "" ) then
2876        grids(igrid)%aa=0.d0
2877      else
2878        if (len(trim(strg))<=30) then
2879          strg1=trim(strg)
2880          read(unit=strg1,fmt=*) grids(igrid)%aa
2881        else
2882          read(unit=strg,fmt=*) grids(igrid)%aa
2883        end if
2884      end if
2885      call paw_rdfromline(" n",line,strg,ierr)
2886      if (strg == "" ) then
2887        grids(igrid)%nn=0
2888      else
2889        if (len(trim(strg))<=30) then
2890          strg1=trim(strg)
2891          read(unit=strg1,fmt=*) grids(igrid)%nn
2892        else
2893          read(unit=strg,fmt=*) grids(igrid)%nn
2894        end if
2895      end if
2896      call paw_rdfromline(" d",line,strg,ierr)
2897      if (strg == "" ) then
2898        grids(igrid)%dd=0.d0
2899      else
2900        if (len(trim(strg))<=30) then
2901          strg1=trim(strg)
2902          read(unit=strg1,fmt=*) grids(igrid)%dd
2903        else
2904          read(unit=strg,fmt=*) grids(igrid)%dd
2905        end if
2906      end if
2907      call paw_rdfromline(" b",line,strg,ierr)
2908      if (strg == "" ) then
2909        grids(igrid)%bb=0.d0
2910      else
2911        if (len(trim(strg))<=30) then
2912          strg1=trim(strg)
2913          read(unit=strg1,fmt=*) grids(igrid)%bb
2914        else
2915          read(unit=strg,fmt=*) grids(igrid)%bb
2916        end if
2917      end if
2918      call paw_rdfromline("istart",line,strg,ierr)
2919      if (len(trim(strg))<=30) then
2920        strg1=trim(strg)
2921        read(unit=strg1,fmt=*) grids(igrid)%istart
2922      else
2923        read(unit=strg,fmt=*) grids(igrid)%istart
2924      end if
2925      call paw_rdfromline("iend",line,strg,ierr)
2926      if (len(trim(strg))<=30) then
2927        strg1=trim(strg)
2928        read(unit=strg1,fmt=*) grids(igrid)%iend
2929      else
2930        read(unit=strg,fmt=*) grids(igrid)%iend
2931      end if
2932      call paw_rdfromline(" id",line,strg,ierr)
2933      grids(igrid)%id = trim(strg)
2934      if(igrid>10)then
2935        close(funit)
2936        msg="igrid>10"
2937        LIBPAW_ERROR(msg)
2938      end if
2939      found=.true.
2940      cycle
2941    end if
2942 
2943 !  End of reading loop
2944  end do
2945 
2946  nphicor=icor
2947  nmesh=igrid
2948  if(nmesh>0)then
2949    LIBPAW_DATATYPE_ALLOCATE(radmesh,(nmesh))
2950    LIBPAW_ALLOCATE(mesh_shift,(nmesh))
2951    do imsh=1,nmesh
2952      radmesh(imsh)%mesh_type=-1
2953      radmesh(imsh)%rstep=zero
2954      radmesh(imsh)%lstep=zero
2955      mesh_shift(imsh)=0
2956      select case(trim(grids(imsh)%eq))
2957        case("r=a*exp(d*i)")
2958          mesh_shift(imsh)=1
2959          radmesh(imsh)%mesh_type=3
2960          radmesh(imsh)%mesh_size=grids(imsh)%iend-grids(imsh)%istart+1+mesh_shift(imsh)
2961          radmesh(imsh)%rstep=grids(imsh)%aa
2962          radmesh(imsh)%lstep=grids(imsh)%dd
2963        case("r=a*i/(1-b*i)")
2964          write(msg, '(3a)' )&
2965 &         '  the grid r=a*i/(1-b*i) is not implemented in ABINIT !',ch10,&
2966 &         '  Action: check your psp file.'
2967          LIBPAW_ERROR(msg)
2968        case("r=a*i/(n-i)")
2969          mesh_shift(imsh)=0
2970          radmesh(imsh)%mesh_type=5
2971          radmesh(imsh)%mesh_size=grids(imsh)%iend-grids(imsh)%istart+1+mesh_shift(imsh)
2972          radmesh(imsh)%rstep=grids(imsh)%aa
2973          radmesh(imsh)%lstep=dble(grids(imsh)%nn)
2974        case("r=a*(exp(d*i)-1)")
2975          mesh_shift(imsh)=0
2976          radmesh(imsh)%mesh_type=2
2977          radmesh(imsh)%mesh_size=grids(imsh)%iend-grids(imsh)%istart+1+mesh_shift(imsh)
2978          if(grids(imsh)%istart==1)radmesh(imsh)%mesh_size=radmesh(imsh)%mesh_size+1
2979          radmesh(imsh)%rstep=grids(imsh)%aa
2980          radmesh(imsh)%lstep=grids(imsh)%dd
2981        case("r=d*i")
2982          mesh_shift(imsh)=0
2983          radmesh(imsh)%mesh_type=1
2984          radmesh(imsh)%mesh_size=grids(imsh)%iend-grids(imsh)%istart+1+mesh_shift(imsh)
2985          if(grids(imsh)%istart==1)radmesh(imsh)%mesh_size=radmesh(imsh)%mesh_size+1
2986          radmesh(imsh)%rstep=grids(imsh)%dd
2987        case("r=(i/n+a)^5/a-a^4")
2988          write(msg, '(3a)' )&
2989 &       '  the grid r=(i/n+a)^5/a-a^4 is not implemented in ABINIT !',ch10,&
2990 &       '  Action: check your psp file.'
2991          LIBPAW_ERROR(msg)
2992      end select
2993    end do
2994  end if
2995 
2996 !Initialize radial meshes
2997  do imsh=1,nmesh
2998    call pawrad_init(radmesh(imsh))
2999  end do
3000 
3001  maxmeshz=maxval(radmesh(:)%mesh_size)
3002  LIBPAW_DATATYPE_ALLOCATE(gridwf,(nphicor))
3003  LIBPAW_DATATYPE_ALLOCATE(statewf,(nphicor))
3004  LIBPAW_ALLOCATE(phitmp,(maxmeshz,nphicor))
3005  phitmp(:,:)=zero
3006 
3007 !Start of reading loop
3008  iaewf=0 ; endfile=.false.
3009  do while (.not.endfile)
3010    read(funit,'(a)',err=11,end=11) readline
3011    line=adjustl(readline);goto 21
3012    11 line="";endfile=.true.
3013    21 continue
3014 
3015 !  --Read CORE WAVE FUNCTIONS PHI
3016    if (line(1:21)=='<ae_core_wavefunction') then
3017      iaewf=iaewf+1
3018      tread=.true.
3019      call paw_rdfromline(" grid",line,strg,ierr)
3020      if (strg == "" ) strg = "unknown"
3021      gridwf(iaewf)=trim(strg)
3022      call paw_rdfromline(" state",line,strg,ierr)
3023      statewf(iaewf)=trim(strg)
3024      do ii=1,nmesh
3025        if(trim(gridwf(iaewf))==trim(grids(ii)%id)) then
3026          mesh_size=grids(ii)%iend-grids(ii)%istart+1
3027          exit
3028        end if
3029      end do
3030      read(funit,*) (phitmp(ir,iaewf),ir=1,mesh_size)
3031      cycle
3032    end if
3033 !  End of reading loop
3034  end do
3035 
3036  if(nphicor>0)then
3037    LIBPAW_ALLOCATE(ncor,(nphicor))
3038    LIBPAW_ALLOCATE(lcor,(nphicor))
3039    LIBPAW_ALLOCATE(energy_cor,(nphicor))
3040    LIBPAW_ALLOCATE(phi_cor,(pawrad%mesh_size,nphicor))
3041    if (present(kappacor)) then
3042      LIBPAW_ALLOCATE(kappacor,(nphicor))
3043    endif
3044    phi_cor(:,:)=zero
3045    do ii=1,nphicor
3046      ncor(ii)=corestate(ii)%nn
3047      lcor(ii)=corestate(ii)%ll
3048      if(present(kappacor)) kappacor(ii)=corestate(ii)%kk
3049      energy_cor(ii)=corestate(ii)%ee
3050      do imsh=1,nmesh
3051        if(trim(gridwf(ii))==trim(grids(imsh)%id)) imeshae=imsh
3052      end do
3053      if ((pawrad%mesh_type/=radmesh(imeshae)%mesh_type) &
3054 &    .or.(pawrad%rstep/=radmesh(imeshae)%rstep) &
3055 &    .or.(pawrad%lstep/=radmesh(imeshae)%lstep)) then
3056        mesh_size=radmesh(imeshae)%mesh_size
3057        LIBPAW_ALLOCATE(work,(mesh_size))
3058        call bound_deriv(phitmp(1:mesh_size,ii),radmesh(imeshae),mesh_size,yp1,ypn)
3059        call paw_spline(radmesh(imeshae)%rad(1:mesh_size),phitmp(1:mesh_size,ii),mesh_size,yp1,ypn,work(1:mesh_size))
3060        ir=pawrad%mesh_size
3061        if (radmesh(imeshae)%rmax<pawrad%rmax+tol8) ir=pawrad_ifromr(pawrad,radmesh(imeshae)%rmax)-1
3062        call paw_splint(mesh_size,radmesh(imeshae)%rad(1:mesh_size),phitmp(1:mesh_size,ii),work(1:mesh_size),&
3063 &                      ir,pawrad%rad(1:ir),phi_cor(1:ir,ii))
3064        phi_cor(1:ir,ii)=phi_cor(1:ir,ii)*pawrad%rad(1:ir)
3065        LIBPAW_DEALLOCATE(work)
3066      else
3067        shft=mesh_shift(imeshae)
3068        mesh_size=min(radmesh(imeshae)%mesh_size,pawrad%mesh_size)
3069        phi_cor(1+shft:mesh_size,ii)=phitmp(1:mesh_size-shft,ii)*radmesh(imeshae)%rad(1:mesh_size-shft)
3070        if (shft==1) phi_cor(1,ii)=zero
3071      end if
3072    end do
3073  end if
3074 
3075  if (allocated(radmesh)) then
3076    call pawrad_free(radmesh)
3077    LIBPAW_DATATYPE_DEALLOCATE(radmesh)
3078  end if
3079  if (allocated(mesh_shift)) then
3080    LIBPAW_DATATYPE_DEALLOCATE(mesh_shift)
3081  end if
3082 
3083  if (allocated(grids)) then
3084    LIBPAW_DATATYPE_DEALLOCATE(grids)
3085  end if
3086  if (allocated(corestate)) then
3087    LIBPAW_DATATYPE_DEALLOCATE(corestate)
3088  end if
3089 
3090  if (allocated(gridwf)) then
3091    LIBPAW_DATATYPE_DEALLOCATE(gridwf)
3092  end if
3093  if (allocated(statewf)) then
3094    LIBPAW_DATATYPE_DEALLOCATE(statewf)
3095  end if
3096  if (allocated(phitmp)) then
3097    LIBPAW_DEALLOCATE(phitmp)
3098  end if
3099 
3100 !Close the XML atomicdata file
3101  close(funit)
3102 
3103 
3104  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

SOURCE

1519  subroutine rdpawpsxml_header(ecut_tmp,filename,paw_setup)
1520 
1521 !Arguments ---------------------------------------------
1522  
1523  character (len=fnlen),intent(in) :: filename
1524  real(dp), intent(inout) :: ecut_tmp(3,2)
1525  type(paw_setup_t),intent(inout) :: paw_setup
1526 !Local variables ---------------------------------------
1527  integer :: funit,ii,ir,igrid,ival,ierr,ishpf,lmax,mesh_size
1528  logical :: endfile,found
1529  character(len=100) :: msg
1530  character (len=XML_RECL) :: line,readline
1531  character (len=XML_RECL) :: strg
1532  character (len=30) :: strg1
1533  real(dp) :: rc(6)
1534  real(dp), allocatable :: shpf(:,:)
1535  type(state_t), pointer :: valstate (:)
1536  type(radial_grid_t), pointer :: grids (:)
1537 
1538 ! *************************************************************************
1539 
1540 !Open the atomicdata XML file for reading
1541  open(newunit=funit,file=filename,form='formatted',status='old', recl=XML_RECL)
1542 
1543 !Start a reading loop
1544  endfile=.false.
1545  found=.false.
1546  paw_setup%rpaw=-1.d0
1547  rc=-1.d0
1548 
1549  do while ((.not.endfile).and.(.not.found))
1550    read(funit,'(a)',err=10,end=10) readline
1551    line=adjustl(readline);goto 20
1552    10 line="";endfile=.true.
1553    20 continue
1554 
1555 !  --Read VERSION
1556    if ((line(1:10)=='<paw_setup').or.(line(1:12)=='<paw_dataset')) then
1557      paw_setup%tread=.true.
1558      igrid=0;ishpf=0
1559      LIBPAW_DATATYPE_ALLOCATE(grids,(10))
1560 
1561      call paw_rdfromline(" version",line,strg,ierr)
1562      paw_setup%version=trim(strg)
1563      cycle
1564    end if
1565 
1566 !  --Read TITLE, ATOMIC CHARGE AND CORE CHARGE
1567    if (line(1:6)=='<atom ') then
1568      paw_setup%atom%tread=.true.
1569      call paw_rdfromline(" symbol",line,strg,ierr)
1570      paw_setup%atom%symbol=trim(strg)
1571      call paw_rdfromline(" Z",line,strg,ierr)
1572      if (len(trim(strg))<=30) then
1573        strg1=trim(strg)
1574        read(unit=strg1,fmt=*) paw_setup%atom%znucl
1575      else
1576        read(unit=strg,fmt=*) paw_setup%atom%znucl
1577      end if
1578      call paw_rdfromline(" core",line,strg,ierr)
1579      if (len(trim(strg))<=30) then
1580        strg1=trim(strg)
1581        read(unit=strg1,fmt=*) paw_setup%atom%zion
1582      else
1583        read(unit=strg,fmt=*) paw_setup%atom%zion
1584      end if
1585      call paw_rdfromline(" valence",line,strg,ierr)
1586      if (len(trim(strg))<=30) then
1587        strg1=trim(strg)
1588        read(unit=strg1,fmt=*) paw_setup%atom%zval
1589      else
1590        read(unit=strg,fmt=*) paw_setup%atom%zval
1591      end if
1592      cycle
1593    end if
1594 
1595 !  --Read Ecut and Ecutdg
1596    if (line(1:8)=='<pw_ecut') then
1597      call paw_rdfromline(" low",line,strg,ierr)
1598      if (len(trim(strg))<=30) then
1599        strg1=trim(strg)
1600        read(unit=strg1,fmt=*) ecut_tmp(1,1)
1601      else
1602        read(unit=strg,fmt=*) ecut_tmp(1,1)
1603      end if
1604      call paw_rdfromline(" medium",line,strg,ierr)
1605      if (len(trim(strg))<=30) then
1606        strg1=trim(strg)
1607        read(unit=strg1,fmt=*) ecut_tmp(2,1)
1608      else
1609        read(unit=strg,fmt=*) ecut_tmp(2,1)
1610      end if
1611      call paw_rdfromline(" high",line,strg,ierr)
1612      if (len(trim(strg))<=30) then
1613        strg1=trim(strg)
1614        read(unit=strg1,fmt=*) ecut_tmp(3,1)
1615      else
1616        read(unit=strg,fmt=*) ecut_tmp(3,1)
1617      end if
1618      cycle
1619    end if
1620 
1621 !  --Read EXCHANGE-CORRELATION TYPE
1622    if (line(1:14)=='<xc_functional') then
1623      paw_setup%xc_functional%tread=.true.
1624      call paw_rdfromline(" type",line,strg,ierr)
1625      paw_setup%xc_functional%functionaltype = trim(strg)
1626      call paw_rdfromline(" name",line,strg,ierr)
1627      paw_setup%xc_functional%name= trim(strg)
1628      cycle
1629    end if
1630 
1631 !  --Read GENERATOR
1632    if (line(1:10)=='<generator') then
1633      paw_setup%generator%tread=.true.
1634      call paw_rdfromline(" type",line,strg,ierr)
1635      paw_setup%generator%gen  = trim(strg)
1636      call paw_rdfromline(" name",line,strg,ierr)
1637      paw_setup%generator%name= trim(strg)
1638      cycle
1639    end if
1640 
1641 !  --Read PAW RADIUS
1642    if (line(1:11)=='<PAW_radius') then
1643      call paw_rdfromline(" rpaw",line,strg,ierr)
1644      if (len(trim(strg))<=30) then
1645        strg1=trim(strg)
1646        read(unit=strg1,fmt=*) paw_setup%rpaw
1647      else
1648        read(unit=strg,fmt=*) paw_setup%rpaw
1649      end if
1650      cycle
1651    end if
1652    if (line(1:11)=='<paw_radius') then
1653      call paw_rdfromline(" rc",line,strg,ierr)
1654      if (len(trim(strg))<=30) then
1655        strg1=trim(strg)
1656        read(unit=strg1,fmt=*) paw_setup%rpaw
1657      else
1658        read(unit=strg,fmt=*) paw_setup%rpaw
1659      end if
1660      cycle
1661    end if
1662 
1663 !  --Read BASIS SIZE, ORBITALS, RC AND OCCUPATIONS/STATE IDs
1664    if (line(1:16)=='<valence_states>') then
1665      paw_setup%valence_states%tread=.true.
1666      LIBPAW_DATATYPE_ALLOCATE(valstate,(50))
1667      ival=0
1668      lmax=0
1669      do while (line(1:17)/='</valence_states>')
1670        read(funit,'(a)') readline;line=adjustl(readline)
1671        if (line(1:6)=='<state') then
1672          ival=ival+1
1673          if (ival>50) then
1674            close(funit)
1675            msg="Error in rdpawps1xml: basis size too large (>50)!"
1676            LIBPAW_ERROR(msg)
1677          end if
1678          call paw_rdfromline(" n",line,strg,ierr)
1679          if (strg == "" ) then 
1680            valstate(ival)%nn=-1    
1681          else
1682            if (len(trim(strg))<=30) then
1683              strg1=trim(strg)
1684              read(unit=strg1,fmt=*) valstate(ival)%nn
1685            else
1686              read(unit=strg,fmt=*) valstate(ival)%nn
1687            end if
1688          end if
1689          call paw_rdfromline(" l",line,strg,ierr)
1690          if (len(trim(strg))<=30) then
1691            strg1=trim(strg)
1692            read(unit=strg1,fmt=*) valstate(ival)%ll
1693          else
1694            read(unit=strg,fmt=*) valstate(ival)%ll
1695          end if
1696          if(valstate(ival)%ll>lmax) lmax=valstate(ival)%ll
1697          call paw_rdfromline(" f",line,strg,ierr)
1698          if (strg == "" ) then 
1699            valstate(ival)%ff=-1.d0
1700          else
1701            if (len(trim(strg))<=30) then
1702              strg1=trim(strg)
1703              read(unit=strg1,fmt=*) valstate(ival)%ff
1704            else
1705              read(unit=strg,fmt=*) valstate(ival)%ff
1706            end if
1707          end if
1708          call paw_rdfromline(" rc",line,strg,ierr)
1709          if (len(trim(strg))<=30) then
1710            strg1=trim(strg)
1711            read(unit=strg1,fmt=*) valstate(ival)%rc
1712          else
1713            read(unit=strg,fmt=*) valstate(ival)%rc
1714          end if
1715          call paw_rdfromline(" e",line,strg,ierr)
1716          if (len(trim(strg))<=30) then
1717            strg1=trim(strg)
1718            read(unit=strg1,fmt=*) valstate(ival)%ee
1719          else
1720            read(unit=strg,fmt=*) valstate(ival)%ee
1721          end if
1722          call paw_rdfromline(" id",line,strg,ierr)
1723          valstate(ival)%id = trim(strg)
1724        end if
1725      end do
1726      cycle
1727    end if
1728 
1729 !  --Read MESH_STEP AND NUMBER OF POINTS
1730    if (line(1:12)=='<radial_grid')then
1731      igrid=igrid+1
1732      call paw_rdfromline(" eq",line,strg,ierr)
1733      grids(igrid)%eq = trim(strg)
1734      call paw_rdfromline(" a",line,strg,ierr)
1735      if (strg == "" ) then
1736        grids(igrid)%aa=0.d0
1737      else
1738        if (len(trim(strg))<=30) then
1739          strg1=trim(strg)
1740          read(unit=strg1,fmt=*) grids(igrid)%aa
1741        else
1742          read(unit=strg,fmt=*) grids(igrid)%aa
1743        end if
1744      end if
1745      call paw_rdfromline(" n",line,strg,ierr)
1746      if (strg == "" ) then
1747        grids(igrid)%nn=0
1748      else
1749        if (len(trim(strg))<=30) then
1750          strg1=trim(strg)
1751          read(unit=strg1,fmt=*) grids(igrid)%nn
1752        else
1753          read(unit=strg,fmt=*) grids(igrid)%nn
1754        end if
1755      end if
1756      call paw_rdfromline(" d",line,strg,ierr)
1757      if (strg == "" ) then
1758        grids(igrid)%dd=0.d0
1759      else
1760        if (len(trim(strg))<=30) then
1761          strg1=trim(strg)
1762          read(unit=strg1,fmt=*) grids(igrid)%dd
1763        else
1764          read(unit=strg,fmt=*) grids(igrid)%dd
1765        end if
1766      end if
1767      call paw_rdfromline(" b",line,strg,ierr)
1768      if (strg == "" ) then
1769        grids(igrid)%bb=0.d0
1770      else
1771        if (len(trim(strg))<=30) then
1772          strg1=trim(strg)
1773          read(unit=strg1,fmt=*) grids(igrid)%bb
1774        else
1775          read(unit=strg,fmt=*) grids(igrid)%bb
1776        end if
1777      end if
1778      call paw_rdfromline("istart",line,strg,ierr)
1779      if (len(trim(strg))<=30) then
1780        strg1=trim(strg)
1781        read(unit=strg1,fmt=*) grids(igrid)%istart
1782      else
1783        read(unit=strg,fmt=*) grids(igrid)%istart
1784      end if
1785      call paw_rdfromline("iend",line,strg,ierr)
1786      if (len(trim(strg))<=30) then
1787        strg1=trim(strg)
1788        read(unit=strg1,fmt=*) grids(igrid)%iend
1789      else
1790        read(unit=strg,fmt=*) grids(igrid)%iend
1791      end if
1792      call paw_rdfromline(" id",line,strg,ierr)
1793      grids(igrid)%id = trim(strg)
1794      if(igrid>10)then
1795        close(funit)
1796        msg="igrid>10"
1797        LIBPAW_ERROR(msg)
1798      end if
1799      cycle
1800    end if
1801 
1802 !  --Read SHAPE TYPE
1803    if (line(1:15)=='<shape_function') then
1804      paw_setup%shape_function%tread=.true.
1805      call paw_rdfromline(" type",line,strg,ierr)
1806      paw_setup%shape_function%gtype = trim(strg)
1807      call paw_rdfromline(" rc",line,strg,ierr)
1808      if (strg /= "" ) then
1809        if (len(trim(strg))<=30) then
1810          strg1=trim(strg)
1811          read(unit=strg1,fmt=*) paw_setup%shape_function%rc
1812        else
1813          read(unit=strg,fmt=*) paw_setup%shape_function%rc
1814        end if
1815      end if
1816      call paw_rdfromline(" lamb",line,strg,ierr)
1817      if (strg == "" ) then
1818        paw_setup%shape_function%lamb=0
1819      else
1820        if (len(trim(strg))<=30) then
1821          strg1=trim(strg)
1822          read(unit=strg1,fmt=*) paw_setup%shape_function%lamb
1823        else
1824          read(unit=strg,fmt=*) paw_setup%shape_function%lamb
1825        end if
1826      end if
1827      found=paw_setup%shape_function%tread
1828      call paw_rdfromline("grid",line,strg,ierr)
1829      paw_setup%shape_function%grid=trim(strg)
1830      if (strg /= "" ) then
1831        paw_setup%shape_function%gtype ="num" 
1832        do ii=1,igrid
1833          if(trim(paw_setup%shape_function%grid)==trim(grids(ii)%id)) then
1834            mesh_size=grids(ii)%iend-grids(ii)%istart+1
1835            exit
1836          end if
1837        end do
1838        if(.not.allocated(shpf)) then
1839          LIBPAW_ALLOCATE(shpf,(mesh_size,7))
1840        end if
1841        ishpf=ishpf+1
1842        read(funit,*) (shpf(ir,ishpf),ir=1,mesh_size)
1843        call paw_rdfromline(" l",line,strg,ierr)
1844        if (strg /= "" ) then
1845          found=.false.
1846          if(paw_setup%valence_states%tread) then
1847            if(ishpf==2*lmax+1) found=.true.
1848          else
1849            write(msg,'(a,a,a)')"the grids and the states must be read before the shapefunction",ch10,&
1850 &           "Action: Modify your XML PAW data file"
1851            LIBPAW_ERROR(msg)
1852          end if
1853        end if
1854      end if
1855      cycle
1856    end if
1857 
1858 !  End of reading loop
1859  end do
1860  if(ival>0)then
1861    LIBPAW_DATATYPE_ALLOCATE(paw_setup%valence_states%state,(ival))
1862    paw_setup%valence_states%state(ival)%tread=.true.
1863    paw_setup%valence_states%nval=ival
1864    do ii=1,ival
1865      paw_setup%valence_states%state(ii)=valstate(ii)
1866    end do
1867  end if
1868  LIBPAW_DATATYPE_DEALLOCATE(valstate)
1869  LIBPAW_DATATYPE_ALLOCATE(paw_setup%radial_grid,(igrid))
1870  paw_setup%radial_grid(igrid)%tread=.true.
1871  paw_setup%ngrid=igrid
1872  do ii=1,igrid
1873    paw_setup%radial_grid(ii)=grids(ii)
1874  end do
1875  LIBPAW_DATATYPE_DEALLOCATE(grids)
1876  if(allocated(shpf)) then
1877    LIBPAW_DEALLOCATE(shpf)
1878  end if
1879  close(funit)
1880 
1881  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

135 type, public :: shape_function_t
136   logical              :: tread=.false.
137   character(len=20)    :: gtype
138   real(dpxml)          :: rc=0
139   character(len=6)     :: grid
140   integer              :: lamb
141   real(dpxml), allocatable :: data(:,:)
142 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

156 type, public :: state_t
157   logical          :: tread=.false.
158   character(len=6) :: id
159   real(dpxml)      :: ff
160   real(dpxml)      :: rc
161   real(dpxml)      :: ee
162   integer          :: nn
163   integer          :: ll
164   integer          :: kk
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