! ! Copyright (C) 2024 INTW group ! ! This file is part of INTW. ! ! INTW is free software: you can redistribute it and/or modify ! it under the terms of the GNU General Public License as published by ! the Free Software Foundation, either version 3 of the License, or ! (at your option) any later version. ! ! INTW is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! You should have received a copy of the GNU General Public License ! along with this program. If not, see <https://www.gnu.org/licenses/>. ! module intw_pseudo !! display: none !! !! This module contains the pseudo-potentials in intwpseudo type, !! and the subroutine read_all_pseudo for reading them from the PP files. !! use kinds, only: dp implicit none save ! variables public :: nt_max, intwpseudo, upf ! subroutines public :: read_all_pseudo private integer, parameter :: nt_max = 99 ! Max number of different atomic types ! NOTE: Haritz 14/03/2024: ! In QE the maximum number of atomic types is 10, however, I didn't find any limit in SIESTA. ! In intw we just set a limit of 2 digits due to the PP filename format (see read_all_pseudo below). type intwpseudo character(len=2) :: psd = ' ' ! Element label character(len=6) :: rel = ' ' ! relativistic: {no|scalar|full} logical :: nlcc ! Non linear core corrections real(kind=dp) :: zp ! z valence logical :: has_so ! If .true. includes spin-orbit ! integer :: mesh ! Number of points in the radial mesh real(kind=dp), dimension(:), allocatable :: r ! r(mesh) radial grid real(kind=dp), dimension(:), allocatable :: rab ! rab(mesh) dr(i)/di ! real(kind=dp) :: rcloc ! vloc = v_ae for r > rcloc integer :: lloc ! l of channel used to generate local potential real(kind=dp), dimension(:), allocatable :: vloc ! vloc(mesh) local atomic potential ! integer :: nbeta ! Number of beta projectors integer :: lmax ! Max l component in beta integer, dimension(:), allocatable :: kbeta ! kbeta(nbeta): number of grid points used for betas ! this defines the cutoff radius for each of them. integer :: kkbeta ! Max number of grid points used for betas integer, dimension(:), allocatable :: lll ! lll(nbeta) l of each projector real(kind=dp), dimension(:), allocatable :: jjj ! jjj(nbeta) j=l+1/2 or l-1/2 of beta real(kind=dp), dimension(:,:), allocatable :: beta ! beta(mesh,nbeta) projectors real(kind=dp), dimension(:,:), allocatable :: dion ! dion(nbeta,nbeta) atomic D_{mu,nu} end type intwpseudo type(intwpseudo), dimension(:), allocatable :: upf contains !--------------------------------------------------------------------- subroutine read_all_pseudo() ! use intw_utility, only: find_free_unit use intw_reading, only: ntyp use intw_input_parameters, only: outdir, prefix ! implicit none ! ! Local variables ! integer :: ios, nr, nt, nb, ir, nb1 integer :: io_unit, ierr character(256) :: file_pseudo character(256) :: dum character(1) :: tag1 character(2) :: tag2 allocate(upf(ntyp)) ! ierr = 0 do nt = 1, ntyp if ( nt>=1 .and. nt<=9 ) then write(tag1,"(i1)") nt write(*,"(a)") "| "//tag1//"-KBPP.txt"//" .. |" file_pseudo=trim(outdir)//trim(prefix)//".save.intw/"//tag1//"-KBPP.txt" else if ( nt>=10 .and. nt<=nt_max ) then write(tag2,"(i2)") nt write(*,"(a)") "| "//tag2//"-KBPP.txt"//" .. |" file_pseudo=trim(outdir)//trim(prefix)//".save.intw/"//tag2//"-KBPP.txt" else stop "ERROR: read_all_pseudo: ntyp > nt_max" end if io_unit = find_free_unit() open(unit=io_unit ,file=file_pseudo, status='old', form='formatted', iostat=ios) read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%psd write(*,"(a)") "| .. for the specie "//trim(upf(nt)%psd)//" |" read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%rel read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%has_so read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%nlcc read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%zp read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%lloc read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%lmax read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%nbeta allocate(upf(nt)%kbeta(upf(nt)%nbeta)) read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%kbeta upf(nt)%kkbeta = maxval( upf(nt)%kbeta(:) ) allocate(upf(nt)%lll(upf(nt)%nbeta)) read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%lll allocate(upf(nt)%jjj(upf(nt)%nbeta)) read(unit=io_unit,fmt="(a)",iostat=ierr) dum if (upf(nt)%has_so) then read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%jjj else read(unit=io_unit,fmt=*,iostat=ierr) dum end if read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%mesh read(unit=io_unit,fmt="(a)",iostat=ierr) dum read(unit=io_unit,fmt=*,iostat=ierr) upf(nt)%rcloc allocate(upf(nt)%dion( upf(nt)%nbeta,upf(nt)%nbeta)) read(unit=io_unit,fmt="(a)",iostat=ierr) dum do nb = 1, upf(nt)%nbeta read(unit=io_unit,fmt=*,iostat=ierr) (upf(nt)%dion(nb,nb1),nb1=1,upf(nt)%nbeta) end do nr=upf(nt)%mesh nb=upf(nt)%nbeta allocate(upf(nt)%r(nr), upf(nt)%rab(nr), upf(nt)%vloc(nr)) allocate(upf(nt)%beta(nr,nb)) read(unit=io_unit,fmt="(a)",iostat=ierr)dum do ir = 1, nr read(unit=io_unit,fmt=*,iostat=ierr) & upf(nt)%r(ir), upf(nt)%rab(ir), upf(nt)%vloc(ir), (upf(nt)%beta(ir,nb), nb=1,upf(nt)%nbeta) end do if (ierr /= 0) then write(unit=*,fmt=*)"ERROR reading PP, ", file_pseudo stop end if ! close(io_unit) end do !nt end subroutine read_all_pseudo end module intw_pseudo