File: rdmft.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 (120 lines) | stat: -rw-r--r-- 3,288 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
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120

! Copyright (C) 2007-2008 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.

!BOP
! !ROUTINE: rdmft
! !INTERFACE:
subroutine rdmft
! !USES:
use modmain
use modrdm
use modmpi
! !DESCRIPTION:
!  Main routine for one-body reduced density matrix functional theory (RDMFT).
!
! !REVISION HISTORY:
!   Created 2008 (Sharma)
!EOP
!BOC
implicit none
! initialise global variables
call init0
call init1
! generate q-point set and gclq array
call init2
! read density and potentials from file
call readstate
! Fourier transform Kohn-Sham potential to G-space
call genvsig
! generate the core wavefunctions and densities
call gencore
! read Fermi energy from file
call readfermi
! find the new linearisation energies
call linengy
! generate the APW and local-orbital radial functions and integrals
call genapwlofr
! generate the spin-orbit coupling radial functions
call gensocfr
! compute the kinetic energy of the core
call energykncr
! generate the first- and second-variational eigenvectors and eigenvalues
call genevfsv
! find the occupation numbers
call occupy
! generate the kinetic matrix elements in the first-variational basis
call genkmat(.true.,.false.)
! open information files (MPI master process only)
if (mp_mpi) then
  open(60,file='RDM_INFO.OUT',form='FORMATTED')
  open(61,file='RDMN_ENERGY.OUT',form='FORMATTED')
  open(62,file='RDMC_ENERGY.OUT',form='FORMATTED')
  if (spinpol) then
    open(63,file='RDMN_MOMENT.OUT',form='FORMATTED')
    open(64,file='RDMC_MOMENT.OUT',form='FORMATTED')
  end if
  open(65,file='RDM_ENERGY.OUT',form='FORMATTED')
! write out general information to RDM_INFO.OUT
  call writeinfo(60)
  write(60,*)
  write(60,'("+------------------------------+")')
  write(60,'("| Self-consistent loop started |")')
  write(60,'("+------------------------------+")')
end if
! begin main RDMFT self-consistent loop
do iscl=1,rdmmaxscl
  if (mp_mpi) then
    write(60,*)
    write(60,'("+--------------------+")')
    write(60,'("| Loop number : ",I4," |")') iscl
    write(60,'("+--------------------+")')
    flush(60)
    write(*,*)
    write(*,'("Info(rdmft): self-consistent loop number : ",I4)') iscl
  end if
! synchronise MPI processes
  call mpi_barrier(mpicom,ierror)
! minimisation over natural orbitals
  call rdmminc
! minimisation over occupation number
  call rdmminn
! compute the RDMFT 'eigenvalues'
  call rdmeval
  if (mp_mpi) then
    call rdmwriteengy(60)
    call writechg(60)
    call writeeval
! write out the total energy
    write(65,'(G18.10)') engytot
    flush(65)
  end if
end do
if (mp_mpi) then
  write(60,*)
  write(60,'("+------------------------------+")')
  write(60,'("| Self-consistent loop stopped |")')
  write(60,'("+------------------------------+")')
! write density to STATE.OUT
  call writestate
  write(60,*)
  write(60,'("Wrote STATE.OUT")')
  write(60,*)
  write(60,'("+----------------------------+")')
  write(60,'("| Elk version ",I1.1,".",I1.1,".",I2.2," stopped |")') version
  write(60,'("+----------------------------+")')
! close information files
  close(60)
  close(61)
  close(62)
  if (spinpol) then
    close(63)
    close(64)
  end if
  close(65)
end if
return
end subroutine
!EOC