1#include "f_defs.h" 2 3module pseudopot_m 4 5 use global_m 6 use kpoint_pool_m, only: kpoint_pool_t 7 use bgw_mpi_m, only: bgw_bcast 8 9 implicit none 10 11 private 12 13 type, abstract :: pseudopot_t 14 integer :: nat 15 integer :: nsp 16 integer :: nkb 17 integer :: nhm 18 integer :: ns 19 integer, allocatable :: ityp(:) 20 integer, allocatable :: nh(:) 21 real(DP), allocatable :: deeq(:,:,:,:) 22 complex(DPC), allocatable :: vkb(:,:,:) 23 contains 24 ! High-level routines 25 procedure(pp_init), deferred :: init 26 procedure(pp_free), deferred :: free 27 procedure(pp_prepare_kpoints), deferred :: prepare_kpoints 28 procedure(pp_get_nh_for_atom), deferred :: get_nh_for_atom 29 procedure(pp_get_vkb_for_atom), deferred :: get_vkb_for_atom 30 ! Aux. routines 31 procedure, nopass :: newunit => pp_newunit 32 endtype pseudopot_t 33 34 35 abstract interface 36 subroutine pp_init(this, fname, kpp, mf) 37 import 38 class(pseudopot_t), intent(out) :: this 39 character(len=*), intent(in) :: fname 40 type(kpoint_pool_t), intent(in) :: kpp 41 type(mf_header_t), intent(in) :: mf 42 end subroutine pp_init 43 end interface 44 45 46 abstract interface 47 subroutine pp_prepare_kpoints(this, ik, kpp, mf, gind_k2m) 48 import 49 class(pseudopot_t), intent(inout) :: this 50 integer, intent(in) :: ik 51 type(kpoint_pool_t), intent(in) :: kpp 52 type(mf_header_t), intent(in) :: mf 53 integer, intent(in) :: gind_k2m(:) 54 end subroutine pp_prepare_kpoints 55 end interface 56 57 58 abstract interface 59 subroutine pp_free(this) 60 import 61 class(pseudopot_t), intent(inout) :: this 62 end subroutine pp_free 63 end interface 64 65 66 abstract interface 67 integer function pp_get_nh_for_atom(this, iatom) 68 import 69 class(pseudopot_t), intent(in) :: this 70 integer, intent(in) :: iatom 71 end function pp_get_nh_for_atom 72 end interface 73 74 75 abstract interface 76 subroutine pp_get_vkb_for_atom(this, iatom, is, vkb_atom, kpp) 77 import 78 class(pseudopot_t), intent(in) :: this 79 integer, intent(in) :: iatom 80 integer, intent(in) :: is 81 complex(DPC), intent(out) :: vkb_atom(:,:) !<(ngkmax,nh) 82 type(kpoint_pool_t), intent(in) :: kpp 83 end subroutine pp_get_vkb_for_atom 84 end interface 85 86 87 !> Abstract type for pseudopotentials that explicitly have the KB in memory 88 !! for each k-point iteration. Right now, all PPs inherit from this type. 89 type, abstract, extends(pseudopot_t) :: pseudopot_explicit_t 90 !> Number of atoms that I own 91 integer :: nat_own 92 !> Number of KB projectors that I own 93 integer :: nkb_own 94 integer, allocatable :: atom_offset(:) 95 contains 96 procedure :: get_nh_for_atom => pp2_get_nh_for_atom 97 procedure :: get_vkb_for_atom => pp2_get_vkb_for_atom 98 endtype pseudopot_explicit_t 99 100 101 public :: pseudopot_t, pseudopot_explicit_t 102 103 104contains 105 106!------------------------------------------------------------------------------- 107! pseudopot_t 108 109integer function pp_newunit(unit) 110 integer, intent(out), optional :: unit 111 112 integer, parameter :: LUN_MIN=201, LUN_MAX=1000 113 logical :: opened 114 integer :: lun 115 116 pp_newunit = -1 117 do lun = LUN_MIN, LUN_MAX 118 inquire(unit=lun, opened=opened) 119 if (.not.opened) then 120 pp_newunit = lun 121 exit 122 endif 123 enddo 124 if (present(unit)) unit = pp_newunit 125 126end function pp_newunit 127 128 129!------------------------------------------------------------------------------- 130! pseudopot_explicit_t 131 132 133integer function pp2_get_nh_for_atom(this, iatom) 134 class(pseudopot_explicit_t), intent(in) :: this 135 integer, intent(in) :: iatom 136 137 PUSH_SUB(pp2_get_nh_for_atom) 138 139 pp2_get_nh_for_atom = this%nh(this%ityp(iatom)) 140 141 POP_SUB(pp2_get_nh_for_atom) 142 143end function pp2_get_nh_for_atom 144 145 146subroutine pp2_get_vkb_for_atom(this, iatom, is, vkb_atom, kpp) 147 class(pseudopot_explicit_t), intent(in) :: this 148 integer, intent(in) :: iatom 149 integer, intent(in) :: is 150 complex(DPC), intent(out) :: vkb_atom(:,:) !<(ngkmax,nh) 151 type(kpoint_pool_t), intent(in) :: kpp 152 153 integer :: nh, ipe, iatom_loc, atom_offset, ngk 154 155 PUSH_SUB(pp2_get_vkb_for_atom) 156 157 ! FHJ: KB projectors are distributed over the atom index round robin 158 ! among ranks in the k-point group. 159 nh = this%get_nh_for_atom(iatom) 160 ipe = INDXG2P(iatom, 1, kpp%inode, 0, kpp%npes) 161 if (ipe==kpp%inode) then 162 iatom_loc = INDXG2L(iatom, 1, kpp%inode, 0, kpp%npes) 163 atom_offset = this%atom_offset(iatom_loc) 164 ngk = size(this%vkb, 1) 165 vkb_atom(:ngk,:nh) = this%vkb(:,atom_offset+1:atom_offset+nh,is) 166 endif 167 call timacc(44,1) 168 call bgw_bcast(vkb_atom, ipe, kpp%comm) 169 call timacc(44,2) 170 171 POP_SUB(pp2_get_vkb_for_atom) 172 173end subroutine pp2_get_vkb_for_atom 174 175 176end module pseudopot_m 177