TABLE OF CONTENTS
ABINIT/m_paw_efield [ Modules ]
NAME
m_paw_efield
FUNCTION
This module contains routines related to the treatment of electric field in the PAW approach.
COPYRIGHT
Copyright (C) 2018-2022 ABINIT group (FJ, PH) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 MODULE m_paw_efield 23 24 use defs_basis 25 use m_abicore 26 use m_errors 27 use m_time, only : timab 28 use m_xmpi, only : xmpi_sum 29 30 use m_pawtab, only : pawtab_type 31 use m_pawrhoij, only : pawrhoij_type 32 33 implicit none 34 35 private 36 37 !public procedures. 38 public :: pawpolev ! Compute the PAW on-site term for polarization 39 40 CONTAINS !========================================================================================
ABINIT/pawpolev [ Functions ]
NAME
pawpolev
FUNCTION
Compute the PAW term for polarization, named expected value term
COPYRIGHT
Copyright (C) 1998-2022 ABINIT group (FJ, PH) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt.
INPUTS
comm_atom=--optional-- MPI communicator over atoms my_natom=number of atoms treated by current processor natom=number of atoms in cell. ntypat = number of atom types pawrhoij(natom) <type(pawrhoij_type)>= paw rhoij occupancies and related data pawtab(ntypat) <type(pawtab_type)>=paw tabulated starting data
OUTPUT
pelev(3)= electronic polarisation. expectation value term (PAW only)
SOURCE
72 subroutine pawpolev(my_natom,natom,ntypat,pawrhoij,pawtab,pelev,& 73 & comm_atom) ! optional argument (parallelism) 74 75 implicit none 76 77 !Arguments --------------------------------------------- 78 !scalars 79 integer,intent(in) :: my_natom,natom,ntypat 80 integer,optional,intent(in) :: comm_atom 81 !arrays 82 real(dp),intent(out) :: pelev(3) 83 type(pawrhoij_type),intent(in) :: pawrhoij(my_natom) 84 type(pawtab_type),intent(in) :: pawtab(ntypat) 85 86 87 !Local variables --------------------------------------- 88 !scalars 89 integer :: iatom,idir,ierr,irhoij,ispden,itypat,jrhoij,klmn 90 logical :: paral_atom 91 real(dp) :: c1,ro_dlt 92 !arrays 93 integer,dimension(3) :: idirindx = (/4,2,3/) 94 real(dp) :: tsec(2) 95 96 ! ************************************************************************* 97 98 DBG_ENTER("COLL") 99 100 call timab(560,1,tsec) 101 102 if (my_natom>0) then 103 ABI_CHECK(pawrhoij(1)%qphase==1,'pawpolev not supposed to be called with qphase/=1!') 104 end if 105 106 !Check for parallelism over atoms 107 paral_atom=(present(comm_atom).and.(my_natom/=natom)) 108 109 !note that when vector r is expanded in real spherical harmonics, the factor 110 !sqrt(four_pi/three) appears, as in the following 111 !x = sqrt(four_pi/three)*r*S_{1,1} 112 !y = sqrt(four_pi/three)*r*S_{1,-1} 113 !z = sqrt(four_pi/three)*r*S_{1,0} 114 ! 115 !the moments pawtab()%qijl do not include such normalization factors 116 !see pawinit.F90 for their definition and computation 117 118 c1=sqrt(four_pi/three) 119 120 pelev=zero 121 do idir=1,3 122 do iatom=1,my_natom 123 itypat=pawrhoij(iatom)%itypat 124 do ispden=1,pawrhoij(iatom)%nspden 125 jrhoij=1 126 do irhoij=1,pawrhoij(iatom)%nrhoijsel 127 klmn=pawrhoij(iatom)%rhoijselect(irhoij) 128 ro_dlt=pawrhoij(iatom)%rhoijp(jrhoij,ispden)*pawtab(itypat)%dltij(klmn) 129 pelev(idir)=pelev(idir)+ro_dlt*c1*pawtab(itypat)%qijl(idirindx(idir),klmn) 130 jrhoij=jrhoij+pawrhoij(iatom)%cplex_rhoij 131 end do 132 end do 133 end do 134 end do 135 136 if (paral_atom) then 137 call xmpi_sum(pelev,comm_atom,ierr) 138 end if 139 140 call timab(560,2,tsec) 141 142 DBG_EXIT("COLL") 143 144 end subroutine pawpolev