[skip-CI] More cleanup, GTH PP re-introduced, fixed a crash with PAW

Modules/read_pseudo.f90

@@ -50,11 +50,8 @@ SUBROUTINE readpp ( input_dft, printout, ecutwfc_pp, ecutrho_pp )
   USE xc_lib,       ONLY: xclib_get_id
   USE radial_grids, ONLY: deallocate_radial_grid, nullify_radial_grid
   USE clib_wrappers,     ONLY: md5_from_file
-  USE read_upf_v1_module,   ONLY: read_upf_v1
-  USE read_upf_new_module,  ONLY: read_upf_new
-  USE upf_auxtools, ONLY: upf_get_pp_format, upf_check_atwfc_norm
+  USE upf_auxtools, ONLY: upf_check_atwfc_norm
   USE upf_to_internal,  ONLY: add_upf_grid, set_upf_q
-  USE read_uspp_module, ONLY: readvan, readrrkj
   USE m_gth,            ONLY: readgth
   !
   IMPLICIT NONE
@@ -133,9 +130,14 @@ SUBROUTINE readpp ( input_dft, printout, ecutwfc_pp, ecutrho_pp )
      !
      CALL mp_bcast (ierr,ionode_id,intra_image_comm)
      !
-     IF ( ierr > 0 ) CALL errore('readpp', &
-          'file '//TRIM(file_pseudo)//' not readable',1)
-     !! Unrecoverable error
+     IF ( ierr == -7 ) THEN
+        !! FIXME: GTH PP files must be read from all processors
+        CALL  readgth( file_pseudo, nt, upf(nt), ierr )
+     END IF
+     IF ( ierr > 0 ) THEN
+        !! Unrecoverable error
+        CALL errore('readpp', 'file '//TRIM(file_pseudo)//' not readable',1)
+     END IF
      !
      CALL upf_bcast(upf(nt), ionode, ionode_id, intra_image_comm)
      !! Success: broadcast the pseudopotential to all processors
@@ -314,7 +316,8 @@ SUBROUTINE upf_bcast(upf, ionode, ionode_id, comm)
   CALL mp_bcast (upf%r,   ionode_id, comm )
   CALL mp_bcast (upf%rab, ionode_id, comm )
   !
-  IF ( upf%nlcc) THEN
+  IF ( upf%nlcc. OR. upf%tpawp ) THEN
+     !! FIXME: PAW uses the pseudo-core charge even when nlcc is not present
      IF ( .NOT. ionode) ALLOCATE( upf%rho_atc(upf%mesh) )
      CALL mp_bcast (upf%rho_atc, ionode_id, comm )
   END IF

upflib/gth.f90

@@ -518,7 +518,7 @@ subroutine deallocate_gth( lflag )
   !
 end subroutine deallocate_gth
 !-----------------------------------------------------------------------
-subroutine readgth (iunps, np, upf, ierr)
+subroutine readgth (psfile, np, upf, ierr)
   !-----------------------------------------------------------------------
   !
   USE upf_kinds,    ONLY: dp
@@ -532,13 +532,13 @@ subroutine readgth (iunps, np, upf, ierr)
   !
   ! I/O
   TYPE (pseudo_upf) :: upf
-  integer, intent(in) :: iunps
+  character(LEN=*), intent(in) :: psfile
   integer, intent(in) :: np
   integer, intent(out):: ierr
   !
   ! Local variables
-  integer  :: ios, pspdat, pspcod, pspxc, lmax, lloc, mmax, ii, jj, ll, nn, nnonloc, &
-              nprl, os, ns, iv, jv
+  integer  :: iunps, ios, pspdat, pspcod, pspxc, lmax, lloc, mmax, &
+              ii, jj, ll, nn, nnonloc, nprl, os, ns, iv, jv
   real(dp) :: rcore, qcore, rc2, prefact, znucl, r2well, rloc, rrl, cc(4)
   character(len=256)            :: info
   character(len=  1), parameter :: ch10=char(10)
@@ -584,6 +584,8 @@ subroutine readgth (iunps, np, upf, ierr)
   allocate(upf%lchi(upf%nwfc))
   upf%lchi(:) = 0
 
+  open(newunit=iunps, file=psfile, form='formatted', status='old', iostat = ios)
+  if ( ios .ne. 0 ) go to 400
   read (iunps, '(a)', end=400, err=400, iostat=ios) info
   read (iunps, *, err=400) znucl, upf%zp, pspdat
   if (upf%zp <= 0._dp .or. upf%zp > 100 ) then
@@ -677,6 +679,7 @@ subroutine readgth (iunps, np, upf, ierr)
       upf%rho_atc(ii) = prefact * exp(-0.5_dp * upf%r(ii)**2 / rc2)
     enddo
   end if
+  close (unit=iunps)
   !
   allocate(upf%lll(upf%nbeta), upf%els_beta(upf%nbeta), upf%dion(upf%nbeta,upf%nbeta))
   allocate(upf%rcut(upf%nbeta), upf%rcutus(upf%nbeta), upf%kbeta(upf%nbeta))

upflib/read_ps.f90

@@ -20,6 +20,7 @@ SUBROUTINE read_ps_new ( psfile, upf, ierr, printout )
   !!    ierr = -4  old Vanderbilt formatted USPP (deprecated)
   !!    ierr = -5  old RRKJ3 USPP format (deprecated)
   !!    ierr = -6  old PWscf NCPP format (deorecated)
+  !!    ierr = -7  Goedecker-Teter-Hutter NCPP
   !! Should be executed on a single processor
   !
   USE upf_io,             ONLY: stdout
@@ -74,6 +75,9 @@ SUBROUTINE read_ps_new ( psfile, upf, ierr, printout )
      ELSE IF (psfile (lm5:l) =='.RRKJ3') THEN
         CALL readrrkj (iunps, upf, ierr)
         IF ( ierr == 0 ) ierr = -5
+     ELSE IF (psfile (lm3:l) =='.gth' .OR. psfile(lm3:l) == '.GTH' ) THEN
+        !! FIXME: should be done in the same way as for the other cases
+        ierr = -7
      ELSE
         CALL read_ncpp (iunps, upf, ierr)
         IF ( ierr == 0 ) ierr = -6
@@ -108,6 +112,8 @@ SUBROUTINE read_ps_new ( psfile, upf, ierr, printout )
            WRITE( stdout, "('file format is RRKJ3')")
         CASE(-6)
            WRITE( stdout, "('file format is old PWscf NC format')")
+        CASE(-7)
+           WRITE( stdout, "('file format is GTH (Goedecker-Teter-Hutter)')")
         CASE DEFAULT
            WRITE( stdout, "('file format could not be determined')")
         END SELECT

upflib/read_upf_new.f90

@@ -81,6 +81,8 @@ CONTAINS
     CALL read_pp_mesh ( upf )
     !
     allocate ( upf%rho_atc(upf%mesh) )
+    !! FIXME: this is needed only if the nonlinear core correction is used,
+    !! FIXME: but with PAW the pseudo-core charge is used also if no nlcc
     IF(upf%nlcc) then
        CALL xmlr_readtag( capitalize_if_v2('pp_nlcc'), &
             upf%rho_atc(:) )

upflib/upf_auxtools.f90

@@ -14,45 +14,9 @@ module upf_auxtools
   implicit none
   private
 
-  public :: upf_get_pp_format
   public :: upf_check_atwfc_norm
 
 contains
-
-!
-!-----------------------------------------------------------------------
-FUNCTION upf_get_pp_format (psfile) result(pp_format)
-  !-----------------------------------------------------------------------
-  IMPLICIT NONE
-  INTEGER :: pp_format
-  CHARACTER (LEN=*) :: psfile
-  INTEGER :: l, lm3, lm4, lm5
-  !
-  l = LEN_TRIM (psfile)
-  lm3 = max(l-3,1)
-  lm4 = max(l-4,1)
-  lm5 = max(l-5,1)
-  IF (psfile (lm3:l) =='.xml' .OR. psfile (lm3:l) =='.XML') THEN
-     pp_format = 0
-  ELSE IF (psfile (lm3:l) =='.upf' .OR. psfile (lm3:l) =='.UPF') THEN
-     pp_format = 1
-  ELSE IF (psfile (lm3:l) =='.vdb' .OR. psfile (lm3:l) =='.van') THEN
-     pp_format = 2
-  ELSE IF (psfile (lm3:l) =='.gth') THEN
-     pp_format = 3
-  ELSE IF (psfile (lm3:l) =='.cpi' .OR. psfile (l-3:l) =='.fhi') THEN
-     pp_format = 6
-  ELSE IF (psfile (lm4:l) =='.cpmd') THEN
-     pp_format = 7
-  ELSE IF (psfile (lm4:l) =='.psml' ) THEN
-     pp_format = 8
-  ELSE IF (psfile (lm5:l) =='.RRKJ3') THEN
-     pp_format = 4
-  ELSE
-     pp_format = 5
-  END IF
-  !
-END FUNCTION upf_get_pp_format
 !
 !---------------------------------------------------------------
 SUBROUTINE upf_check_atwfc_norm(upf,psfile)