TABLE OF CONTENTS


ABINIT/mklocl [ Functions ]

[ Top ] [ Functions ]

NAME

 mklocl

FUNCTION

 This method is a wrapper for mklocl_recipspace and mklocl_realspace.
 It does some consistency checks before calling one of the two methods.

 Optionally compute :
  option=1 : local ionic potential throughout unit cell
  option=2 : contribution of local ionic potential to E gradient wrt xred
  option=3 : contribution of local ionic potential to
                stress tensor (only with reciprocal space computations)
  option=4 : contribution of local ionic potential to
                second derivative of E wrt xred  (only with reciprocal space computations)

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  if(option==3) eei=local pseudopotential part of total energy (hartree)
  gmet(3,3)=reciprocal space metric ($\textrm{Bohr}^{-2}$).
  gprimd(3,3)=reciprocal space dimensional primitive translations
  gsqcut=cutoff on $|G|^2$: see setup1 for definition (doubled sphere).
  mgfft=maximum size of 1D FFTs
  mpi_enreg=information about MPI parallelization
  natom=number of atoms in unit cell.
  nattyp(ntypat)=number of atoms of each type in cell.
  nfft=(effective) number of FFT grid points (for this processor)
  ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft
  nspden=number of spin-density components
  ntypat=number of types of atoms.
  option= (see above)
  pawtab(ntypat) <type(pawtab_type)>=paw tabulated starting data
  ph1d(2,3*(2*mgfft+1)*natom)=1-dim structure factor phase information.
  psps <type(pseudopotential_type)>=variables related to pseudopotentials
  qprtrb(3)= integer wavevector of possible perturbing potential
   in basis of reciprocal lattice translations
  rhog(2,nfft)=electron density rho(G) (electrons/$\textrm{Bohr}^3$)
    (needed if option==2 or if option==4)
  rhor(nfft,nspden)=electron density in electrons/bohr**3.
    (needed if option==2 or if option==4)
  rprimd(3,3)=dimensional primitive translations in real space (bohr)
  ucvol=unit cell volume ($\textrm{Bohr}^3$).
  vprtrb(2)=complex amplitude of possible perturbing potential; if nonzero,
   perturbing potential is added of the form
   $V(G)=(vprtrb(1)+I*vprtrb(2))/2$ at the values G=qprtrb and
   $(vprtrb(1)-I*vprtrb(2))/2$ at $G=-qprtrb$ (integers)
  xred(3,natom)=reduced dimensionless atomic coordinates

OUTPUT

  (if option==1) vpsp(nfft)=local crystal pseudopotential in real space.
  (if option==2) grtn(3,natom)=grads of Etot wrt tn.
  (if option==3) lpsstr(6)=components of local psp part of stress tensor
   (Cartesian coordinates, symmetric tensor) in hartree/$\textrm{bohr}^3$
   Store 6 unique components in order 11, 22, 33, 32, 31, 21
  (if option==4) dyfrlo(3,3,natom)=d2 Eei/d tn(i)/d tn(j).  (Hartrees)

SIDE EFFECTS

NOTES

 Note that the present routine is tightly connected to the dfpt_vlocal.f routine,
 that compute the derivative of the local ionic potential
 with respect to one atomic displacement. The argument list
 and the internal loops to be considered were sufficiently different
 as to make the two routine different.

PARENTS

      forces,prcref,prcref_PMA,respfn,setvtr

CHILDREN

      mklocl_realspace,mklocl_recipspace,mklocl_wavelets,wvl_rho_abi2big
      xred2xcart

SOURCE

 81 #if defined HAVE_CONFIG_H
 82 #include "config.h"
 83 #endif
 84 
 85 #include "abi_common.h"
 86 
 87 subroutine mklocl(dtset, dyfrlo,eei,gmet,gprimd,grtn,gsqcut,lpsstr,mgfft,&
 88 &  mpi_enreg,natom,nattyp,nfft,ngfft,nspden,ntypat,option,pawtab,ph1d,psps,qprtrb,&
 89 &  rhog,rhor,rprimd,ucvol,vprtrb,vpsp,wvl,wvl_den,xred)
 90 
 91  use defs_basis
 92  use defs_datatypes
 93  use defs_abitypes
 94  use defs_wvltypes
 95  use m_profiling_abi
 96  use m_errors
 97 
 98  use m_geometry,   only : xred2xcart
 99  use m_pawtab,     only : pawtab_type
100 
101 #if defined HAVE_BIGDFT
102  use BigDFT_API, only : ELECTRONIC_DENSITY
103  use m_abi2big, only : wvl_rho_abi2big
104 #endif
105 
106 !This section has been created automatically by the script Abilint (TD).
107 !Do not modify the following lines by hand.
108 #undef ABI_FUNC
109 #define ABI_FUNC 'mklocl'
110  use interfaces_67_common, except_this_one => mklocl
111 !End of the abilint section
112 
113  implicit none
114 
115 !Arguments ------------------------------------
116 !scalars
117  integer,intent(in) :: mgfft,natom,nfft,nspden,ntypat,option
118  real(dp),intent(in) :: eei,gsqcut,ucvol
119  type(MPI_type),intent(in) :: mpi_enreg
120  type(dataset_type),intent(in) :: dtset
121  type(pseudopotential_type),intent(in) :: psps
122  type(wvl_internal_type), intent(in) :: wvl
123  type(wvl_denspot_type), intent(inout) :: wvl_den
124  type(pawtab_type),intent(in)  :: pawtab(ntypat*psps%usepaw)
125 !arrays
126  integer,intent(in) :: nattyp(ntypat),ngfft(18),qprtrb(3)
127  real(dp),intent(in) :: gmet(3,3),gprimd(3,3),ph1d(2,3*(2*mgfft+1)*natom)
128  real(dp),intent(in) :: rhog(2,nfft),rprimd(3,3)
129  real(dp),intent(in) :: vprtrb(2),xred(3,natom)
130  real(dp),intent(in),target :: rhor(nfft,nspden)
131  real(dp),intent(out) :: dyfrlo(3,3,natom),grtn(3,natom),lpsstr(6)
132  real(dp),intent(inout) :: vpsp(nfft)
133 
134 !Local variables-------------------------------
135 !scalars
136  character(len=500) :: message
137 !arrays
138  real(dp),allocatable :: xcart(:,:)
139 #if defined HAVE_BIGDFT
140  real(dp),pointer :: rhor_ptr(:,:)
141 #endif
142 
143 ! *************************************************************************
144 
145  if (option < 1 .or. option > 4) then
146    write(message,'(a,i0,a,a)')&
147 &   'From the calling routine, option=',option,ch10,&
148 &   'The only allowed values are between 1 and 4.'
149    MSG_ERROR(message)
150  end if
151  if (option > 2 .and. .not.psps%vlspl_recipSpace) then
152    write(message,'(a,i0,a,a,a,a)')&
153 &   'From the calling routine, option=',option,ch10,&
154 &   'but the local part of the pseudo-potential is in real space.',ch10,&
155 &   'Action: set icoulomb = 0 to turn-off real space computations.'
156    MSG_ERROR(message)
157  end if
158  if (option > 2 .and. dtset%usewvl == 1) then
159    write(message,'(a,i0,a,a)')&
160 &   'From the calling routine, option=',option,ch10,&
161 &   'but this is not implemented yet from wavelets.'
162    MSG_ERROR(message)
163  end if
164 
165  if (dtset%usewvl == 0) then
166 !  Plane wave case
167    if (psps%vlspl_recipSpace) then
168      call mklocl_recipspace(dyfrlo,eei,gmet,gprimd,grtn,gsqcut,lpsstr,mgfft, &
169 &     mpi_enreg,psps%mqgrid_vl,natom,nattyp,nfft,ngfft, &
170 &     ntypat,option,dtset%paral_kgb,ph1d,psps%qgrid_vl,qprtrb,rhog,ucvol, &
171 &     psps%vlspl,vprtrb,vpsp)
172    else
173      call mklocl_realspace(grtn,dtset%icoulomb,mpi_enreg,natom,nattyp,nfft, &
174 &     ngfft,dtset%nscforder,nspden,ntypat,option,pawtab,psps,rhog,rhor, &
175 &     rprimd,dtset%typat,ucvol,dtset%usewvl,vpsp,xred)
176    end if
177  else
178 !  Store xcart for each atom
179    ABI_ALLOCATE(xcart,(3, dtset%natom))
180    call xred2xcart(dtset%natom, rprimd, xcart, xred)
181 !  Eventually retrieve density
182 #if defined HAVE_BIGDFT
183    if (option>1.and.wvl_den%denspot%rhov_is/=ELECTRONIC_DENSITY) then
184      rhor_ptr => rhor ! Just to bypass intent(inout)
185      call wvl_rho_abi2big(1,rhor_ptr,wvl_den)
186    end if
187 #endif
188 !  Wavelets case
189    call mklocl_wavelets(dtset%efield, grtn, mpi_enreg, dtset%natom, &
190 &   nfft, nspden, option, rprimd, vpsp, &
191 &   wvl_den, wvl, xcart)
192    ABI_DEALLOCATE(xcart)
193  end if
194 
195 end subroutine mklocl