scnc
/
siesta
mirror of https://gitlab.com/siesta-project/siesta.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

New documentation regarding the density matrix with SO added

This commit is contained in:
Javier Junquera 2022-07-18 13:27:32 +02:00
parent a8bf5db931
commit 617d9a09e2
3 changed files with 66 additions and 3 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

19
Src/m_occ_proj_dftu.F90
View File

@ -687,6 +687,25 @@ subroutine occ_proj_dftu
! Transform the density matrix between the two atomic orbital
! in a complex (2x2) matrix
! When SO coupling is considered, the density matrix and the
! Hamiltonian must be globally Hermitian
! (see last sentence of Section 7 of the
! technical SIESTA paper in JPCM 14, 2745 (2002).
! However, in diag3k, the sign of the imaginary part of the
! (up,down) matrix elements is changed:
! Dnew(ind,1) = Dnew(ind,1) + real(D11,dp)
! Dnew(ind,2) = Dnew(ind,2) + real(D22,dp)
! Dnew(ind,3) = Dnew(ind,3) + real(D12,dp)
! Dnew(ind,4) = Dnew(ind,4) - aimag(D12)
! Dnew(ind,5) = Dnew(ind,5) + aimag(D11)
! Dnew(ind,6) = Dnew(ind,6) + aimag(D22)
! Dnew(ind,7) = Dnew(ind,7) + real(D21,dp)
! Dnew(ind,8) = Dnew(ind,8) + aimag(D21)
! In the subroutines to compute the corresponding DFT+U
! matrix elements:
! - We change locally the sign of this imaginary part, so the
! density matrix recovers all its properties
! - The occupations are computed with the "good" DM.
Dij(1) =cmplx( Di(jo,1), Di(jo,5), kind = dp )
Dij(2) =cmplx( Di(jo,2), Di(jo,6), kind = dp )
Dij(3) =cmplx( Di(jo,3), -1.0_dp * Di(jo,4), kind = dp )

23
Src/m_pot_dftu.F90
View File

@ -914,6 +914,24 @@ subroutine dftu_so_hamil_2( H_dftu_so, fal, stressl )
! Store in complex format the density matrix between
! orbitals mu and nu
! When SO coupling is considered, the density matrix and the Hamiltonian
! must be globally Hermitian (see last sentence of Section 7 of the
! technical SIESTA paper in JPCM 14, 2745 (2002).
! However, in diag3k, the sign of the imaginary part of the (up,down)
! matrix elements is changed:
! Dnew(ind,1) = Dnew(ind,1) + real(D11,dp)
! Dnew(ind,2) = Dnew(ind,2) + real(D22,dp)
! Dnew(ind,3) = Dnew(ind,3) + real(D12,dp)
! Dnew(ind,4) = Dnew(ind,4) - aimag(D12)
! Dnew(ind,5) = Dnew(ind,5) + aimag(D11)
! Dnew(ind,6) = Dnew(ind,6) + aimag(D22)
! Dnew(ind,7) = Dnew(ind,7) + real(D21,dp)
! Dnew(ind,8) = Dnew(ind,8) + aimag(D21)
! In the subroutines to compute the corresponding DFT+U matrix elements:
! - We change locally the sign of this imaginary part, so the
! density matrix recovers all its properties
! - The Potential matrix elements are computed with the "good" DM.
Dscf_cmplx_1 = cmplx(Di(jo,1),Di(jo,5), dp)
Dscf_cmplx_2 = cmplx(Di(jo,2),Di(jo,6), dp)
Dscf_cmplx_3 = cmplx(Di(jo,3),-Di(jo,4), dp)
@ -1207,12 +1225,17 @@ subroutine dftu_so_hamil_2( H_dftu_so, fal, stressl )
E_dftu_so = 0.0_dp
do ind = 1, maxnh
! Note the change in the imaginary part of the (up,down) component,
! as discussed above
Dscf_cmplx_1 = cmplx(Dscf(ind,1),Dscf(ind,5), dp)
Dscf_cmplx_2 = cmplx(Dscf(ind,2),Dscf(ind,6), dp)
Dscf_cmplx_3 = cmplx(Dscf(ind,3),-Dscf(ind,4), dp)
Dscf_cmplx_4 = cmplx(Dscf(ind,7),Dscf(ind,8), dp)
! Compute the energy according to the Equations developed in
! the doxygen documentation
! Here, we are computing the trace of the the potential times the
! density matrix, considering that both of them are (2x2) matrices
! and the global hermitian of the density matrix
E_dftu_so = E_dftu_so + &
& 0.5_dp * ( real( H_dftu_so_Hubbard(ind,1)*conjg(Dscf_cmplx_1), dp) + &
& real( H_dftu_so_Hubbard(ind,2)*conjg(Dscf_cmplx_2), dp) + &

27
Src/setup_hamiltonian.F
View File

@ -235,15 +235,36 @@
! Non-SCF part of total energy
call update_E0()
! Hubbard term for LDA+U: energy, forces, stress and matrix elements ....
! Hubbard term for DFT+U: energy, forces, stress and matrix elements ....
if( switch_dftu ) then
call re_alloc(fal, 1, 3, 1, na_u, 'fal', 'setup_hamiltonian')
if ( spin%NCol ) then
call die('LDA+U cannot be used with non-collinear spin.')
call die('DFT+U cannot be used with non-collinear spin.')
else if ( spin%SO ) then
H_dftu_so => val(H_dftu_so_2D)
call dftu_so_hamil( H_dftu_so, fal, stressl )
! When SO coupling is considered, the density matrix and the Hamiltonian
! must be globally Hermitian (see last sentence of Section 7 of the
! technical SIESTA paper in JPCM 14, 2745 (2002).
! However, in diag3k, the sign of the imaginary part of the (up,down)
! matrix elements is changed:
! Dnew(ind,1) = Dnew(ind,1) + real(D11,dp)
! Dnew(ind,2) = Dnew(ind,2) + real(D22,dp)
! Dnew(ind,3) = Dnew(ind,3) + real(D12,dp)
! Dnew(ind,4) = Dnew(ind,4) - aimag(D12)
! Dnew(ind,5) = Dnew(ind,5) + aimag(D11)
! Dnew(ind,6) = Dnew(ind,6) + aimag(D22)
! Dnew(ind,7) = Dnew(ind,7) + real(D21,dp)
! Dnew(ind,8) = Dnew(ind,8) + aimag(D21)
! In the subroutines to compute the corresponding DFT+U matrix elements:
! - We change locally the sign of this imaginary part, so the
! density matrix recovers all its properties
! - The Potential matrix elements are computed with the "good" DM.
! - Here, at the time of adding up the new Hamiltonian matrix elements
! to the potential, we change the sign of the imaginary part of the
! (up,down) component
!
!------- H(u,u)
H(:,1) = H(:,1) + real(H_dftu_so(:,1), dp)
H(:,5) = H(:,5) + dimag(H_dftu_so(:,1))
@ -343,7 +364,7 @@
! ind = listhptr(io) + j
! jo = listh(ind)
! write(6,'(a,6i7,8f12.5)')
! . ' Node, Nodes, io, jo, j, ind, H-LDAU = ' ,
! . ' Node, Nodes, io, jo, j, ind, H-DFTU = ' ,
! . Node, Nodes, io, jo, j, ind, H_dftu_so(ind,:)
! enddo
! enddo