src/rdmdexcdn.f90

   1
   2 ! Copyright (C) 2007-2008 J. K. Dewhurst, S. Sharma and E. K. U. Gross.
   3 ! This file is distributed under the terms of the GNU General Public License.
   4 ! See the file COPYING for license details.
   5
   6 subroutine rdmdexcdn(dedn)
   7 ! calculates derivative of exchange-correlation energy w.r.t. occupation numbers
   8 use modmain
   9 implicit none
  10 ! arguments
  11 real(8), intent(inout) :: dedn(nstsv,nkpt)
  12 ! local variables
  13 integer ik1,ik2,ik3
  14 integer ist1,ist2,iv(3)
  15 ! parameter for calculating the functional derivatives
  16 real(8), parameter :: eps=1.d-12
  17 real(8) t1,t2,t3,t4
  18 ! external functions
  19 real(8) r3dist
  20 external r3dist
  21 if (rdmxctype.eq.0) return
  22
  23 ! calculate the pre-factor
  24 if (rdmxctype.eq.1) then
  25   t1=1.d0/occmax
  26 else if (rdmxctype.eq.2) then
  27   if (spinpol) then
  28     t1=rdmalpha
  29   else
  30     t1=2.d0*rdmalpha*(0.25d0)**rdmalpha
  31   end if
  32 else
  33   write(*,*)
  34   write(*,'("Error(rdmdexcdn): rdmxctype not defined : ",I8)') rdmxctype
  35   write(*,*)
  36 end if
  37
  38
  39 do ik1=1,nkpt
  40   do ist1=1,nstsv
  41     do ik2=1,nkptnr
  42 ! find the equivalent reduced k-point
  43       iv(:)=ivknr(:,ik2)
  44       ik3=ikmap(iv(1),iv(2),iv(3))
  45       do ist2=1,nstsv
  46 ! Hartree-Fock functional
  47         if (rdmxctype.eq.1) then
  48           t2=t1*occsv(ist2,ik3)
  49 ! SDLG functional
  50         else if (rdmxctype.eq.2) then
  51           t3=max(occsv(ist1,ik1),eps)
  52           t4=max(occsv(ist2,ik3),eps)
  53           t2=t1*(t4**rdmalpha)/(t3**(1.d0-rdmalpha))
  54         end if
  55         dedn(ist1,ik1)=dedn(ist1,ik1)+t2*vnlrdm(ist1,ik1,ist2,ik2)
  56       end do
  57     end do
  58   end do
  59 end do
  60 return
  61 end subroutine