TABLE OF CONTENTS


ABINIT/xcmult [ Functions ]

[ Top ] [ Functions ]

NAME

 xcmult

FUNCTION

 In the case of GGA, multiply the different gradient of spin-density
 by the derivative of the XC functional with respect
 to the norm of the gradient, then divide it by the
 norm of the gradient

COPYRIGHT

 Copyright (C) 2001-2018 ABINIT group (XG)
 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

  depsxc(nfft,nspgrad)=derivative of Exc with respect to the (spin-)density,
    or to the norm of the gradient of the (spin-)density,
    further divided by the norm of the gradient of the (spin-)density
   The different components of depsxc will be
   for nspden=1,         depsxc(:,1)=d(rho.exc)/d(rho)
         and if ngrad=2, depsxc(:,2)=1/2*1/|grad rho_up|*d(rho.exc)/d(|grad rho_up|)
                                      +   1/|grad rho|*d(rho.exc)/d(|grad rho|)
         (do not forget : |grad rho| /= |grad rho_up| + |grad rho_down|
   for nspden=2,         depsxc(:,1)=d(rho.exc)/d(rho_up)
                         depsxc(:,2)=d(rho.exc)/d(rho_down)
         and if ngrad=2, depsxc(:,3)=1/|grad rho_up|*d(rho.exc)/d(|grad rho_up|)
                         depsxc(:,4)=1/|grad rho_down|*d(rho.exc)/d(|grad rho_down|)
                         depsxc(:,5)=1/|grad rho|*d(rho.exc)/d(|grad rho|)
  nfft=(effective) number of FFT grid points (for this processor)
  ngrad = must be 2
  nspden=number of spin-density components
  nspgrad=number of spin-density and spin-density-gradient components

OUTPUT

  (see side effects)

SIDE EFFECTS

  rhonow(nfft,nspden,ngrad*ngrad)=
   at input :
    electron (spin)-density in real space and its gradient,
    either on the unshifted grid (if ishift==0,
      then equal to rhor), or on the shifted grid
     rhonow(:,:,1)=electron density in electrons/bohr**3
     rhonow(:,:,2:4)=gradient of electron density in el./bohr**4
   at output :
    rhonow(:,:,2:4) has been multiplied by the proper factor,
    described above.

PARENTS

      m_pawxc,rhotoxc

CHILDREN

SOURCE

 60 #if defined HAVE_CONFIG_H
 61 #include "config.h"
 62 #endif
 63 
 64 #include "abi_common.h"
 65 
 66 
 67 subroutine xcmult (depsxc,nfft,ngrad,nspden,nspgrad,rhonow)
 68 
 69  use defs_basis
 70  use m_profiling_abi
 71 
 72 !This section has been created automatically by the script Abilint (TD).
 73 !Do not modify the following lines by hand.
 74 #undef ABI_FUNC
 75 #define ABI_FUNC 'xcmult'
 76 !End of the abilint section
 77 
 78  implicit none
 79 
 80 !Arguments ------------------------------------
 81 !scalars
 82  integer,intent(in) :: nfft,ngrad,nspden,nspgrad
 83 !arrays
 84  real(dp),intent(in) :: depsxc(nfft,nspgrad)
 85  real(dp),intent(inout) :: rhonow(nfft,nspden,ngrad*ngrad)
 86 
 87 !Local variables-------------------------------
 88 !scalars
 89  integer :: idir,ifft
 90  real(dp) :: rho_tot,rho_up
 91 
 92 ! *************************************************************************
 93 
 94  do idir=1,3
 95 
 96    if(nspden==1)then
 97 !$OMP PARALLEL DO PRIVATE(ifft) SHARED(depsxc,idir,nfft,rhonow)
 98      do ifft=1,nfft
 99        rhonow(ifft,1,1+idir)=rhonow(ifft,1,1+idir)*depsxc(ifft,2)
100      end do
101 
102    else
103 
104 !    In the spin-polarized case, there are more factors to take into account
105 !$OMP PARALLEL DO PRIVATE(ifft,rho_tot,rho_up) SHARED(depsxc,idir,nfft,rhonow)
106      do ifft=1,nfft
107        rho_tot=rhonow(ifft,1,1+idir)
108        rho_up =rhonow(ifft,2,1+idir)
109        rhonow(ifft,1,1+idir)=rho_up *depsxc(ifft,3)         + rho_tot*depsxc(ifft,5)
110        rhonow(ifft,2,1+idir)=(rho_tot-rho_up)*depsxc(ifft,4)+ rho_tot*depsxc(ifft,5)
111      end do
112 
113    end if ! nspden==1
114 
115  end do ! End loop on directions
116 
117 end subroutine xcmult