Commit 3db65716 authored by Henning Janssen's avatar Henning Janssen

Moved phase factors for spin-offdiagonal elements of lda+u/gf to t_sym

parent b5433f1a
......@@ -12,7 +12,7 @@ CONTAINS
SUBROUTINE cdnval(eig_id, mpi,kpts,jspin,noco,input,banddos,cell,atoms,enpara,stars,&
vacuum,sphhar,sym,vTot,oneD,cdnvalJob,den,regCharges,dos,results,&
moments,gfinp,hub1inp,hub1data,coreSpecInput,mcd,slab,orbcomp,greensfCoeffs,angle)
moments,gfinp,hub1inp,hub1data,coreSpecInput,mcd,slab,orbcomp,greensfCoeffs)
!************************************************************************************
! This is the FLEUR valence density generator
......@@ -87,7 +87,6 @@ SUBROUTINE cdnval(eig_id, mpi,kpts,jspin,noco,input,banddos,cell,atoms,enpara,st
TYPE(t_slab), OPTIONAL, INTENT(INOUT) :: slab
TYPE(t_orbcomp), OPTIONAL, INTENT(INOUT) :: orbcomp
TYPE(t_greensfCoeffs), OPTIONAL, INTENT(INOUT) :: greensfCoeffs
REAL, OPTIONAL, INTENT(IN) :: angle(:)
! Scalar Arguments
INTEGER, INTENT(IN) :: eig_id, jspin
......@@ -247,11 +246,12 @@ SUBROUTINE cdnval(eig_id, mpi,kpts,jspin,noco,input,banddos,cell,atoms,enpara,st
CALL abcof(input,atoms,sym,cell,lapw,noccbd,usdus,noco,ispin,oneD,&
eigVecCoeffs%acof(:,0:,:,ispin),eigVecCoeffs%bcof(:,0:,:,ispin),&
eigVecCoeffs%ccof(-atoms%llod:,:,:,:,ispin),zMat,eig,force)
IF (atoms%n_u.GT.0) CALL n_mat(atoms,sym,noccbd,usdus,ispin,we,eigVecCoeffs,den%mmpMat(:,:,:,ispin))
IF (atoms%n_u.GT.0.AND.noco%l_mperp.AND.(ispin==jsp_end)) CALL n_mat21(atoms,sym,angle,noccbd,we,denCoeffsOffdiag,&
eigVecCoeffs,den%mmpMat(:,:,:,3))
IF (atoms%n_u.GT.0.AND.noco%l_mperp.AND.(ispin==jsp_end)) &
CALL n_mat21(atoms,sym,noccbd,we,denCoeffsOffdiag,eigVecCoeffs,den%mmpMat(:,:,:,3))
IF (gfinp%n.GT.0) CALL bzIntegrationGF(atoms,gfinp,sym,input,angle,ispin,noccbd,dosWeights,resWeights,dosBound,kpts%wtkpt(ikpt),&
IF (gfinp%n.GT.0) CALL bzIntegrationGF(atoms,gfinp,sym,input,ispin,noccbd,dosWeights,resWeights,dosBound,kpts%wtkpt(ikpt),&
results%ef,eig,denCoeffsOffdiag,usdus,eigVecCoeffs,greensfCoeffs,ispin==jsp_end)
! perform Brillouin zone integration and summation over the
......
MODULE m_nmat21
! ************************************************************
! This subroutine calculates the density matrix n^{s}_{m,m'}
! for a given atom 'n' and l-quantum number 'l'. The l's for
! all atoms are stored in lda_u(), if lda_u()<0, no +U is used.
! For details see Eq.(12) of Shick et al. PRB 60, 10765 (1999)
! Part of the LDA+U package G.B., Oct. 2000
! ************************************************************
CONTAINS
SUBROUTINE n_mat21(atoms,sym,angle,ne,we,denCoeffsOffdiag,eigVecCoeffs,n_mmp)
!------------------------------------------------------------
!This subroutine calculates the density matrix n^{s}_{m,m'}
!for a given atom 'n' and l-quantum number 'l'. The l's for
!all atoms are stored in lda_u(), if lda_u()<0, no +U is used.
!For details see Eq.(12) of Shick et al. PRB 60, 10765 (1999)
!Part of the LDA+U package G.B., Oct. 2000
!------------------------------------------------------------
USE m_types
USE m_constants
USE m_types
USE m_constants
IMPLICIT NONE
CONTAINS
IMPLICIT NONE
SUBROUTINE n_mat21(atoms,sym,ne,we,denCoeffsOffdiag,eigVecCoeffs,n_mmp)
TYPE(t_sym), INTENT(IN) :: sym
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_eigVecCoeffs), INTENT(IN) :: eigVecCoeffs
TYPE(t_denCoeffsOffDiag), INTENT(IN) :: denCoeffsOffdiag
REAL, INTENT(IN) :: angle(:)
INTEGER, INTENT(IN) :: ne
REAL, INTENT(IN) :: we(:)!(input%neig)
COMPLEX, INTENT(INOUT) :: n_mmp(-lmaxU_const:,-lmaxU_const:,:)
......@@ -113,7 +114,7 @@ MODULE m_nmat21
ENDDO
nr_tmp = matmul( transpose( conjg(d_tmp) ) , n_tmp)
n1_tmp = matmul( nr_tmp, d_tmp )
phase = exp(ImagUnit*angle(isi))
phase = exp(ImagUnit*sym%phase(isi))
DO m = -l,l
DO mp = -l,l
n_mmp(m,mp,i_u) = n_mmp(m,mp,i_u) + conjg(n1_tmp(m,mp)) * fac * phase
......
......@@ -33,7 +33,9 @@ MODULE m_types_sym
!No of 2D-sym ops
INTEGER ::nop2
!Wigner matrix for lda+u
COMPLEX, ALLOCATABLE:: d_wgn(:, :, :, :)
COMPLEX, ALLOCATABLE :: d_wgn(:, :, :, :)
!Phase factors for offidagonal lda+u
REAL, ALLOCATABLE :: phase(:)
!
! Atom sepecific stuff
!
......@@ -85,6 +87,7 @@ CONTAINS
call mpi_bc(this%multab,rank,mpi_comm)
call mpi_bc(this%nop2,rank,mpi_comm)
call mpi_bc(this%d_wgn,rank,mpi_comm)
call mpi_bc(this%phase,rank,mpi_comm)
call mpi_bc(this%invsatnr,rank,mpi_comm)
call mpi_bc(this%invarop,rank,mpi_comm)
call mpi_bc(this%invarind,rank,mpi_comm)
......
......@@ -4,7 +4,6 @@ set(fleur_F77 ${fleur_F77}
set(fleur_F90 ${fleur_F90}
greensf/greensfImag.f90
greensf/greensfImag21.f90
greensf/angles.f90
greensf/onsite.f90
greensf/kkintgr.f90
greensf/kk_cutoff.f90
......
......@@ -11,14 +11,13 @@ MODULE m_gfcalc
USE m_hybridization
USE m_crystalfield
IMPLICIT NONE
CONTAINS
SUBROUTINE bzIntegrationGF(atoms,gfinp,sym,input,angle,ispin,nbands,dosWeights,resWeights,&
indBound,wtkpt,ef,eig,denCoeffsOffdiag,usdus,eigVecCoeffs,greensfCoeffs,l21)
SUBROUTINE bzIntegrationGF(atoms,gfinp,sym,input,ispin,nbands,dosWeights,resWeights,indBound,&
wtkpt,ef,eig,denCoeffsOffdiag,usdus,eigVecCoeffs,greensfCoeffs,l21)
USE m_greensfImag
USE m_greensfImag21
......@@ -40,13 +39,12 @@ MODULE m_gfcalc
REAL, INTENT(IN) :: dosWeights(:,:) !Precalculated tetrahedron weights for the current k-point
INTEGER, INTENT(IN) :: indBound(:,:) !Gives the range where the tetrahedron weights are non-zero
REAL, INTENT(IN) :: eig(:) !Eigenvalues for the current k-point
REAL, INTENT(IN) :: angle(:) !Phases for spin-offdiagonal part
CALL timestart("Greens Function: Imaginary Part")
CALL greensfImag(atoms,gfinp,sym,input,ispin,nbands,dosWeights,resWeights,indBound,&
wtkpt,ef,eig,usdus,eigVecCoeffs,greensfCoeffs)
IF(gfinp%l_mperp.AND.l21) THEN
CALL greensfImag21(atoms,gfinp,sym,angle,input,nbands,dosWeights,resWeights,indBound,&
CALL greensfImag21(atoms,gfinp,sym,input,nbands,dosWeights,resWeights,indBound,&
wtkpt,ef,eig,denCoeffsOffdiag,eigVecCoeffs,greensfCoeffs)
ENDIF
CALL timestop("Greens Function: Imaginary Part")
......@@ -55,7 +53,7 @@ MODULE m_gfcalc
SUBROUTINE postProcessGF(greensf,greensfCoeffs,atoms,gfinp,input,sym,noco,vTot,hub1inp,hub1data,results,angle)
SUBROUTINE postProcessGF(greensf,greensfCoeffs,atoms,gfinp,input,sym,noco,vTot,hub1inp,hub1data,results)
!contains all the modules for calculating properties from the greens function
USE m_onsite
......@@ -69,7 +67,6 @@ MODULE m_gfcalc
TYPE(t_hub1inp), INTENT(IN) :: hub1inp
TYPE(t_results), INTENT(IN) :: results
TYPE(t_potden), INTENT(IN) :: vTot
REAL, INTENT(IN) :: angle(:)
TYPE(t_hub1data), INTENT(INOUT) :: hub1data
TYPE(t_greensfCoeffs), INTENT(INOUT) :: greensfCoeffs
TYPE(t_greensf), INTENT(INOUT) :: greensf
......@@ -80,9 +77,9 @@ MODULE m_gfcalc
CALL timestart("Green's Function: Postprocess")
CALL rot_projDOS(sym,atoms,gfinp,input,angle,greensfCoeffs)
CALL rot_projDOS(sym,atoms,gfinp,input,greensfCoeffs)
!Perform the Kramer-Kronigs-Integration if we only have the imaginary part at this point
CALL calc_onsite(atoms,gfinp,input,sym,noco,results%ef,angle,greensfCoeffs,greensf)
CALL calc_onsite(atoms,gfinp,input,sym,noco,results%ef,greensfCoeffs,greensf)
!-------------------------------------------------------------
! Calculate various properties from the greens function
!-------------------------------------------------------------
......
......@@ -24,12 +24,11 @@ MODULE m_greensfImag21
CONTAINS
SUBROUTINE greensfImag21(atoms,gfinp,sym,angle,input,nbands,dosWeights,resWeights,ind,wtkpt,ef,eig,denCoeffsOffDiag,eigVecCoeffs,greensfCoeffs)
SUBROUTINE greensfImag21(atoms,gfinp,sym,input,nbands,dosWeights,resWeights,ind,wtkpt,ef,eig,denCoeffsOffDiag,eigVecCoeffs,greensfCoeffs)
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_gfinp), INTENT(IN) :: gfinp
TYPE(t_sym), INTENT(IN) :: sym
REAL, INTENT(IN) :: angle(:)
TYPE(t_input), INTENT(IN) :: input
TYPE(t_eigVecCoeffs), INTENT(IN) :: eigVecCoeffs
TYPE(t_denCoeffsOffDiag), INTENT(IN) :: denCoeffsOffDiag
......
......@@ -27,7 +27,7 @@ MODULE m_onsite
CONTAINS
SUBROUTINE calc_onsite(atoms,gfinp,input,sym,noco,ef,angle,greensfCoeffs,g)
SUBROUTINE calc_onsite(atoms,gfinp,input,sym,noco,ef,greensfCoeffs,g)
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_gfinp), INTENT(IN) :: gfinp
......@@ -37,9 +37,7 @@ MODULE m_onsite
TYPE(t_noco), INTENT(IN) :: noco
TYPE(t_input), INTENT(IN) :: input
REAL, INTENT(IN) :: ef
REAL, INTENT(IN) :: angle(:)
!-Local Scalars
INTEGER i_gf,ie,l,m,mp,nType,jspin,ipm,kkcut,lp,nTypep,spin_cut,nn,natom
REAL fac,del,eb,et
INTEGER it,is,isi
......@@ -123,7 +121,7 @@ MODULE m_onsite
calc_mat = matmul( transpose( conjg(d_mat) ) , &
g21(ie,:,:))
calc_mat = matmul( calc_mat, d_mat )
phase = exp(ImagUnit*angle(isi))
phase = exp(ImagUnit*sym%phase(isi))
DO m = -l,l
DO mp = -l,l
g%gmmpMat(ie,m,mp,3,ipm,i_gf) =&
......
......@@ -12,14 +12,13 @@ MODULE m_rot_gf
CONTAINS
SUBROUTINE rot_projDOS(sym,atoms,gfinp,input,angle,greensfCoeffs)
SUBROUTINE rot_projDOS(sym,atoms,gfinp,input,greensfCoeffs)
TYPE(t_sym), INTENT(IN) :: sym
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_gfinp), INTENT(IN) :: gfinp
TYPE(t_input), INTENT(IN) :: input
TYPE(t_greensfCoeffs), INTENT(INOUT) :: greensfCoeffs
REAL, INTENT(IN) :: angle(:)
COMPLEX, ALLOCATABLE :: curr_dos(:,:,:),calc_mat(:,:,:)
COMPLEX d_mat(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const)
......@@ -69,7 +68,7 @@ MODULE m_rot_gf
d_mat(m,mp) = sym%d_wgn(m,mp,l,isi)
ENDDO
ENDDO
phase = MERGE(exp(ImagUnit*angle(isi)),CMPLX(1.0,0.0),ispin.EQ.3)
phase = MERGE(exp(ImagUnit*sym%phase(isi)),CMPLX(1.0,0.0),ispin.EQ.3)
DO ie = 1, gfinp%ne
calc_mat(ie,:,:) = matmul( transpose( conjg(d_mat) ) , curr_dos(ie,:,:))
calc_mat(ie,:,:) = matmul( calc_mat(ie,:,:), d_mat )
......
......@@ -32,4 +32,5 @@ init/prp_xcfft.f90
init/stepf.F90
init/strgn.f90
init/lapw_dim.F90
init/angles.f90
)
MODULE m_angles
CONTAINS
USE m_types
USE m_constants
!Calculate the correct angles for the rotation in spin-space (needed for LDA+U and Greens functions with noco)
IMPLICIT NONE
SUBROUTINE angles(sym,angle)
CONTAINS
USE m_types
USE m_constants
!Calculate the correct phase shifts for the rotation in spin-space
!(needed for LDA+U and Green's functions with noco%l_mperp)
IMPLICIT NONE
SUBROUTINE angles(sym)
TYPE(t_sym), INTENT(IN) :: sym
REAL, INTENT(OUT) :: angle(:)
TYPE(t_sym), INTENT(INOUT) :: sym
INTEGER iop,d,t
angle = 0.0
sym%phase = 0.0
DO iop = 1, sym%nop
d = det(sym%mrot(:,:,iop))
t = (sym%mrot(1,1,iop)+sym%mrot(2,2,iop)+sym%mrot(3,3,iop)) * d
IF(t.EQ.-1) angle(iop) = 1.0
IF(t.EQ.0) angle(iop) = 2.0/3.0
IF(t.EQ.1) angle(iop) = 1.0/2.0
IF(t.EQ.2) angle(iop) = 1.0/3.0
IF(t.EQ.3) angle(iop) = 0.0
angle(iop) = d*angle(iop)*pi_const
IF(t.EQ.-1) sym%phase(iop) = 1.0
IF(t.EQ.0) sym%phase(iop) = 2.0/3.0
IF(t.EQ.1) sym%phase(iop) = 1.0/2.0
IF(t.EQ.2) sym%phase(iop) = 1.0/3.0
IF(t.EQ.3) sym%phase(iop) = 0.0
sym%phase(iop) = d*sym%phase(iop)*pi_const
ENDDO
END SUBROUTINE
......
......@@ -13,6 +13,8 @@ CONTAINS
!Generates missing symmetry info.
!tau,mrot and nop have to be specified alread
USE m_dwigner
USE m_angles !Phase factors for spin-offdiagonal lda+u
USE m_constants
USE m_mapatom
USE m_od_mapatom
use m_ptsym
......@@ -44,10 +46,16 @@ CONTAINS
END IF
!Generated wigner symbols for LDA+U (includes DFT+HubbardI)
IF (ALLOCATED(sym%d_wgn)) DEALLOCATE(sym%d_wgn)
ALLOCATE(sym%d_wgn(-3:3,-3:3,3,sym%nop))
IF (atoms%n_u+gfinp%n.GT.0) THEN
CALL d_wigner(sym%nop,sym%mrot,cell%bmat,3,sym%d_wgn)
IF(ALLOCATED(sym%d_wgn)) DEALLOCATE(sym%d_wgn)
ALLOCATE(sym%d_wgn(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,lmaxU_const,sym%nop))
IF(atoms%n_u+gfinp%n.GT.0) THEN !replace with atoms%n_u+gfinp%n
CALL d_wigner(sym%nop,sym%mrot,cell%bmat,lmaxU_const,sym%d_wgn)
!For spin-offdiagonal parts, we need additional phase factors
IF(noco%l_mperp) THEN
IF(ALLOCATED(sym%phase)) DEALLOCATE(sym%phase)
ALLOCATE(sym%phase(sym%nop),source=0.0)
CALL angles(sym)
ENDIF
END IF
!Atom specific symmetries
......
......@@ -103,7 +103,6 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
INTEGER(HID_T) :: banddosFile_id
#endif
LOGICAL :: l_error, perform_MetaGGA
REAL :: angle(sym%nop)
CALL regCharges%init(input,atoms)
CALL dos%init(input,atoms,kpts,vacuum)
......@@ -116,12 +115,10 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
!Only calculate the greens function when needed
CALL greensfCoeffs%init(gfinp,input,atoms,noco)
CALL gfinp%eContour(results%ef,mpi%irank,gOnsite%nz,gOnsite%e,gOnsite%de)
gOnsite%gmmpMat = 0.0
CALL gOnsite%reset(gfinp)
IF(atoms%n_hia.GT.0.AND.mpi%irank==0.AND.PRESENT(hub1data)) hub1data%mag_mom = 0.0
ENDIF
IF(gfinp%n+atoms%n_u.GT.0.AND.noco%l_mperp) CALL angles(sym,angle)
CALL outDen%init(stars, atoms, sphhar, vacuum, noco, input%jspins, POTDEN_TYPE_DEN)
CALL EnergyDen%init(stars, atoms, sphhar, vacuum, noco, input%jspins, POTDEN_TYPE_EnergyDen)
......@@ -139,12 +136,13 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
IF (sliceplot%slice) CALL cdnvalJob%select_slice(sliceplot,results,input,kpts,noco,jspin)
CALL cdnval(eig_id,mpi,kpts,jspin,noco,input,banddos,cell,atoms,enpara,stars,vacuum,&
sphhar,sym,vTot,oneD,cdnvalJob,outDen,regCharges,dos,results,moments,gfinp,&
hub1inp,hub1data,coreSpecInput,mcd,slab,orbcomp,greensfCoeffs,angle)
hub1inp,hub1data,coreSpecInput,mcd,slab,orbcomp,greensfCoeffs)
END DO
IF(PRESENT(gOnsite).AND.mpi%irank.EQ.0) THEN
IF(gfinp%n.GT.0) THEN
CALL postProcessGF(gOnsite,greensfCoeffs,atoms,gfinp,input,sym,noco,vTot,hub1inp,hub1data,results,angle)
CALL postProcessGF(gOnsite,greensfCoeffs,atoms,gfinp,input,sym,noco,vTot,&
hub1inp,hub1data,results)
ENDIF
ENDIF
......
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment