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
|
! Copyright (C) 2008 F. Bultmark, F. Cricchio and L. Nordstrom.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.
!BOP
! !ROUTINE: writetm3
! !INTERFACE:
subroutine writetm3
! !USES:
use modmain
use moddftu
use modtest
use modvars
! !DESCRIPTION:
! Decompose the density matrix into 3-index tensor moments and write to
! {\tt TENSMOM.OUT}. See {\it Phys. Rev. B} {\bf 80}, 035121 (2009) and
! {\it J. Phys.: Condens. Matter} {\bf 7} 9947 (1995). See also the routines
! {\tt tm2todm} and {\tt tm3todm}.
!
! !REVISION HISTORY:
! Created April 2008 (F. Cricchio and L. Nordstrom)
! Updated, December 2021 (JKD)
!EOP
!BOC
implicit none
! local variables
integer is,ia,ias,idu
integer l,p,k,r,t
real(8) ehx,sm,t0,t1
! automatic arrays
real(8) wkpr(-lmmaxdm:lmmaxdm)
complex(8) zkpr(-lmmaxdm:lmmaxdm)
complex(8) gamma(lmmaxdm,2,lmmaxdm,2)
! external functions
real(8), external :: trzhmm
wkpr(:)=0.d0
! open TENSMOM.OUT file
open(50,file='TENSMOM'//trim(filext),form='FORMATTED',action='WRITE')
write(50,'("Density matrix decomposition in coupled tensor moments")')
write(50,'("Components are in the spherical basis")')
write(50,'(" (see Phys. Rev. B. 80, 035121 (2009))")')
! loop over DFT+U entries
do idu=1,ndftu
is=isldu(1,idu)
l=isldu(2,idu)
! scale factor for conventional normalisation
t0=sqrt(dble((2*l+1)*2))
do ia=1,natoms(is)
ias=idxas(ia,is)
write(50,*)
write(50,'("Species : ",I4," (",A,"), atom : ",I4)') is,trim(spsymb(is)),ia
write(50,'(" l = ",I1)') l
do k=0,2*l
do p=0,1
do r=abs(k-p),k+p
! decompose density matrix in 3-index tensor moment components
call dmtotm3(l,k,p,r,lmmaxdm,dmatmt(:,:,:,:,ias),wkpr)
! determine the contribution to muffin-tin Hartree + exchange energy
call tm3todm(l,k,p,r,lmmaxdm,wkpr,gamma)
ehx=0.5d0*trzhmm(lmmaxdm*2,gamma,vmatmt(:,:,:,:,ias))
! write out tensor moment components and vector magnitude
write(50,*)
write(50,'(" k = ",I1,", p = ",I1,", r = ",I1)') k,p,r
if (tm3vdl) then
! G. van der Laan complex convention
call tm3rtoz(l,k,p,r,lmmaxdm,wkpr,zkpr)
do t=-r,r
write(50,'(" t = ",I2," : ",2F14.8)') t,zkpr(t)
end do
else
! new real convention
sm=0.d0
do t=-r,r
t1=t0*wkpr(t)
write(50,'(" t = ",I2," : ",F14.8)') t,t1
sm=sm+t1**2
end do
sm=sqrt(sm)
write(50,'(" magnitude : ",F14.8)') sm
end if
write(50,'(" Hartree + exchange energy : ",F14.8)') ehx
! write to VARIABLES.OUT if required, but only on the last iteration
if (wrtvars.and.tlast) then
call writevars('wkpr',n1=is,n2=ia,n3=l,n4=k,n5=p,n6=r,nv=2*r+1, &
rva=wkpr(-r:r))
end if
end do
end do
end do
! end loop over atoms and species
end do
end do
close(50)
! write last entry of tensor moment components to test file if required
if (test) then
call writetest(820,'Coupled tensor moments',nv=size(wkpr),tol=2.d-2,rva=wkpr)
end if
end subroutine
!EOC
|
