File: ksinvert.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster, sid
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (88 lines) | stat: -rw-r--r-- 2,357 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88

! Copyright (C) 2018 P. Elliott, J. K. Dewhurst, S. Sharma and E. K. U. Gross.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

subroutine ksinvert
use modmain
use modmpi
implicit none
! local variables
integer idm,it
! allocatable arrays
real(8), allocatable :: rhomt_(:,:),rhoir_(:)
real(8), allocatable :: magmt_(:,:,:),magir_(:,:)
real(8), allocatable :: rfmt(:,:),rfir(:)
! copy the existing density and magnetisation
allocate(rhomt_(npmtmax,natmtot),rhoir_(ngtot))
call rfcopy(rhomt,rhoir,rhomt_,rhoir_)
allocate(rfmt(npmtmax,natmtot),rfir(ngtot))
if (spinpol) then
  allocate(magmt_(npmtmax,natmtot,ndmag),magir_(ngtot,ndmag))
  do idm=1,ndmag
    call rfcopy(magmt(:,:,idm),magir(:,idm),magmt_(:,:,idm),magir_(:,idm))
  end do
end if
if (mp_mpi) then
  write(*,*)
  write(*,'("Info(ksinvert): inverting the Kohn-Sham equations")')
end if
do it=1,maxitksi
  if (mp_mpi.and.(mod(it,10).eq.0)) then
    write(*,'("Info(ksinvert): done ",I4," iterations of ",I4)') it,maxitksi
  end if
  call hmlrad
  call genevfsv
  call occupy
  call rhomag
! determine the residual and add it to the exchange-correlation potential
  call residual
! set the constant part of V_xc equal to that of LDA/GGA
  call rfint0(vxc0,vxcmt,vxcir)
! add the external and Hartree potentials
  call potks(.false.)
  call genvsig
end do
deallocate(rhomt_,rhoir_,rfmt,rfir)
if (spinpol) deallocate(magmt_,magir_)
! synchronise MPI processes
call mpi_barrier(mpicom,ierror)
return

contains

subroutine residual
implicit none
! local variables
integer is,ias,np
! external functions
real(8) rfinp
external rfinp
do ias=1,natmtot
  is=idxis(ias)
  np=npmt(is)
  rfmt(1:np,ias)=rhomt(1:np,ias)-rhomt_(1:np,ias)
  call daxpy(np,tauksi,rfmt(:,ias),1,vxcmt(:,ias),1)
end do
rfir(:)=rhoir(:)-rhoir_(:)
call daxpy(ngtot,tauksi,rfir,1,vxcir,1)
resoep=rfinp(rfmt,rfir,rfmt,rfir)
if (spinpol) then
  do idm=1,ndmag
    do ias=1,natmtot
      is=idxis(ias)
      np=npmt(is)
      rfmt(1:np,ias)=magmt(1:np,ias,idm)-magmt_(1:np,ias,idm)
      call daxpy(np,tauksi,rfmt(:,ias),1,bxcmt(:,ias,idm),1)
    end do
    rfir(:)=magir(:,idm)-magir_(:,idm)
    call daxpy(ngtot,tauksi,rfir,1,bxcir(:,idm),1)
    resoep=resoep+rfinp(rfmt,rfir,rfmt,rfir)
  end do
end if
resoep=resoep/omega
return
end subroutine

end subroutine