TABLE OF CONTENTS


ABINIT/m_paw_efield [ Modules ]

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

17 #if defined HAVE_CONFIG_H
18 #include "config.h"
19 #endif
20 
21 #include "abi_common.h"
22 
23 MODULE m_paw_efield
24 
25  use defs_basis
26  use m_abicore
27  use m_errors
28  use m_time, only : timab
29  use m_xmpi, only : xmpi_sum
30 
31  use m_pawtab,   only : pawtab_type
32  use m_pawrhoij, only : pawrhoij_type
33 
34  implicit none
35 
36  private
37 
38 !public procedures.
39  public :: pawpolev ! Compute the PAW on-site term for polarization
40 
41 CONTAINS  !========================================================================================

ABINIT/pawpolev [ Functions ]

[ Top ] [ Functions ]

NAME

 pawpolev

FUNCTION

 Compute the PAW term for polarization, named expected value term

COPYRIGHT

 Copyright (C) 1998-2018 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)

PARENTS

      berryphase_new

CHILDREN

      timab,xmpi_sum

SOURCE

 79 subroutine pawpolev(my_natom,natom,ntypat,pawrhoij,pawtab,pelev,&
 80 &                   comm_atom) ! optional argument (parallelism)
 81 
 82 
 83 !This section has been created automatically by the script Abilint (TD).
 84 !Do not modify the following lines by hand.
 85 #undef ABI_FUNC
 86 #define ABI_FUNC 'pawpolev'
 87 !End of the abilint section
 88 
 89  implicit none
 90 
 91 !Arguments ---------------------------------------------
 92 !scalars
 93  integer,intent(in) :: my_natom,natom,ntypat
 94  integer,optional,intent(in) :: comm_atom
 95 !arrays
 96  real(dp),intent(out) :: pelev(3)
 97  type(pawrhoij_type),intent(in) :: pawrhoij(my_natom)
 98  type(pawtab_type),intent(in) :: pawtab(ntypat)
 99 
100 
101 !Local variables ---------------------------------------
102 !scalars
103  integer :: iatom,idir,ierr,irhoij,ispden,itypat,jrhoij,klmn
104  logical :: paral_atom
105  real(dp) :: c1,ro_dlt
106 !arrays
107  integer,dimension(3) :: idirindx = (/4,2,3/)
108  real(dp) :: tsec(2)
109 
110 ! *************************************************************************
111 
112  DBG_ENTER("COLL")
113 
114  call timab(560,1,tsec)
115 
116 !Check for parallelism over atoms
117  paral_atom=(present(comm_atom).and.(my_natom/=natom))
118 
119 !note that when vector r is expanded in real spherical harmonics, the factor
120 !sqrt(four_pi/three) appears, as in the following
121 !x = sqrt(four_pi/three)*r*S_{1,1}
122 !y = sqrt(four_pi/three)*r*S_{1,-1}
123 !z = sqrt(four_pi/three)*r*S_{1,0}
124 !
125 !the moments pawtab()%qijl do not include such normalization factors
126 !see pawinit.F90 for their definition and computation
127 
128  c1=sqrt(four_pi/three)
129 
130  pelev=zero
131  do idir=1,3
132    do iatom=1,my_natom
133      itypat=pawrhoij(iatom)%itypat
134      do ispden=1,pawrhoij(iatom)%nspden
135        jrhoij=1
136        do irhoij=1,pawrhoij(iatom)%nrhoijsel
137          klmn=pawrhoij(iatom)%rhoijselect(irhoij)
138          ro_dlt=pawrhoij(iatom)%rhoijp(jrhoij,ispden)*pawtab(itypat)%dltij(klmn)
139          pelev(idir)=pelev(idir)+ro_dlt*c1*pawtab(itypat)%qijl(idirindx(idir),klmn)
140          jrhoij=jrhoij+pawrhoij(iatom)%cplex
141        end do
142      end do
143    end do
144  end do
145 
146  if (paral_atom) then
147    call xmpi_sum(pelev,comm_atom,ierr)
148  end if
149 
150  call timab(560,2,tsec)
151 
152  DBG_EXIT("COLL")
153 
154 end subroutine pawpolev