Commit d42a5b4c authored by Henning Janssen's avatar Henning Janssen

Several fixes for spin-off diagonal elements

parent 836c675c
......@@ -235,7 +235,7 @@ SUBROUTINE cdnval(eig_id, mpi,kpts,jspin,noco,input,banddos,cell,atoms,enpara,st
IF (atoms%n_gf.GT.0.AND.(input%tria.OR.input%gfTet)) THEN
CALL timestart("TetrahedronWeights")
CALL tetrahedronInit(ikpt,kpts,input,SIZE(ev_list),results%eig(ev_list,:,jsp),&
greensfCoeffs,results%ef,resWeights,dosWeights,dosBound)
greensfCoeffs,results%ef,resWeights,dosWeights,dosBound)
CALL timestop("TetrahedronWeights")
ENDIF
DO ispin = jsp_start, jsp_end
......
......@@ -88,11 +88,11 @@ MODULE m_dosWeights
l_bloechl = .TRUE.
ALLOCATE( dos_weights(g%ne) )
dos_weights = 0.0
e_ind(:,1) = g%ne+1
e_ind(:,2) = 0
e_ind(:,1) = 0
e_ind(:,2) = g%ne+1
weights = 0.0
fac = count(kpts%bkp(:).EQ.ikpt)
fac = MERGE(1,count(kpts%bkp(:).EQ.ikpt),kpts%nkptf.EQ.0)
DO itet = 1, kpts%ntet
IF(ALL(kpts%ntetra(1:4,itet).NE.ikpt)) CYCLE !search for the tetrahedra containing ikpt
......
......@@ -52,6 +52,7 @@ MODULE m_greensfImag21
IF(.NOT.input%l_gfsphavg) CALL juDFT_error("NOCO-offdiagonal + Radial dependence of onsite-GF not implemented",calledby="onsite21")
l_tria = (input%tria.OR.input%gfTet).AND..NOT.input%l_hist
DO i_gf = 1, atoms%n_gf
nType = atoms%gfelem(i_gf)%atomType
l = atoms%gfelem(i_gf)%l
......@@ -95,20 +96,16 @@ MODULE m_greensfImag21
!
!Contribution from states
!
im(ie,m,mp,1) = im(ie,m,mp,1) + (weight *&
CONJG(eigVecCoeffs%acof(ib,lmp,natom,1)) * eigVecCoeffs%acof(ib,lm,natom,2) * denCoeffsOffdiag%uu21n(l,nType)&
+ CONJG(eigVecCoeffs%bcof(ib,lmp,natom,1)) * eigVecCoeffs%bcof(ib,lm,natom,2) * denCoeffsOffdiag%dd21n(l,nType)&
+ CONJG(eigVecCoeffs%acof(ib,lmp,natom,1)) * eigVecCoeffs%bcof(ib,lm,natom,2) * denCoeffsOffdiag%ud21n(l,nType)&
+ CONJG(eigVecCoeffs%bcof(ib,lmp,natom,1)) * eigVecCoeffs%acof(ib,lm,natom,2) * denCoeffsOffdiag%du21n(l,nType))
im(ie,m,mp,1) = im(ie,m,mp,1) + weight *&
(CONJG(eigVecCoeffs%acof(ib,lmp,natom,2)) * eigVecCoeffs%acof(ib,lm,natom,1) * denCoeffsOffdiag%uu21n(l,nType)&
+ CONJG(eigVecCoeffs%bcof(ib,lmp,natom,2)) * eigVecCoeffs%bcof(ib,lm,natom,1) * denCoeffsOffdiag%dd21n(l,nType)&
+ CONJG(eigVecCoeffs%acof(ib,lmp,natom,2)) * eigVecCoeffs%bcof(ib,lm,natom,1) * denCoeffsOffdiag%ud21n(l,nType)&
+ CONJG(eigVecCoeffs%bcof(ib,lmp,natom,2)) * eigVecCoeffs%acof(ib,lm,natom,1) * denCoeffsOffdiag%du21n(l,nType))
IF(.NOT.input%l_gfsphavg) THEN
im(ie,m,mp,2) = im(ie,m,mp,2) + (weight *&
conjg(eigVecCoeffs%acof(ib,lmp,natom,2))*eigVecCoeffs%acof(ib,lm,natom,1))
im(ie,m,mp,3) = im(ie,m,mp,3) + (weight *&
conjg(eigVecCoeffs%bcof(ib,lmp,natom,2))*eigVecCoeffs%bcof(ib,lm,natom,1))
im(ie,m,mp,4) = im(ie,m,mp,4) + (weight *&
conjg(eigVecCoeffs%acof(ib,lmp,natom,2))*eigVecCoeffs%bcof(ib,lm,natom,1))
im(ie,m,mp,5) = im(ie,m,mp,5) + (weight *&
conjg(eigVecCoeffs%bcof(ib,lmp,natom,2))*eigVecCoeffs%acof(ib,lm,natom,1))
im(ie,m,mp,2) = im(ie,m,mp,2) + weight * conjg(eigVecCoeffs%acof(ib,lmp,natom,2)) * eigVecCoeffs%acof(ib,lm,natom,1)
im(ie,m,mp,3) = im(ie,m,mp,3) + weight * conjg(eigVecCoeffs%bcof(ib,lmp,natom,2)) * eigVecCoeffs%bcof(ib,lm,natom,1)
im(ie,m,mp,4) = im(ie,m,mp,4) + weight * conjg(eigVecCoeffs%acof(ib,lmp,natom,2)) * eigVecCoeffs%bcof(ib,lm,natom,1)
im(ie,m,mp,5) = im(ie,m,mp,5) + weight * conjg(eigVecCoeffs%bcof(ib,lmp,natom,2)) * eigVecCoeffs%acof(ib,lm,natom,1)
ENDIF
!
!Contribution from local Orbitals
......@@ -116,7 +113,7 @@ MODULE m_greensfImag21
DO ilo = 1, atoms%nlo(nType)
IF(atoms%llo(ilo,nType).EQ.l) THEN
!TODO: Something is wrong here for noco%l_soc and noco%l_mperp with a local orbital on the same l
im(ie,m,mp,1) = im(ie,m,mp,1) + (weight *&
im(ie,m,mp,1) = im(ie,m,mp,1) + weight *&
(conjg(eigVecCoeffs%acof(ib,lmp,natom,2))*eigVecCoeffs%ccof(m,ib,ilo,natom,1) *&
denCoeffsOffDiag%uulo21n(ilo,nType) +&
conjg(eigVecCoeffs%ccof(mp,ib,ilo,natom,2))*eigVecCoeffs%acof(ib,lm,natom,1) *&
......@@ -124,7 +121,7 @@ MODULE m_greensfImag21
conjg(eigVecCoeffs%bcof(ib,lmp,natom,2))*eigVecCoeffs%ccof(m,ib,ilo,natom,1) *&
denCoeffsOffDiag%dulo21n(ilo,nType) +&
conjg(eigVecCoeffs%ccof(mp,ib,ilo,natom,2))*eigVecCoeffs%bcof(ib,lm,natom,1) *&
denCoeffsOffDiag%ulod21n(ilo,nType)))
denCoeffsOffDiag%ulod21n(ilo,nType))
DO ilop = 1, atoms%nlo(nType)
IF(atoms%llo(ilop,nType).EQ.l) THEN
im(ie,m,mp,1) = im(ie,m,mp,1) + (weight * denCoeffsOffDiag%uloulop21n(ilo,ilop,nType) *&
......@@ -143,7 +140,7 @@ MODULE m_greensfImag21
!Rotate the eqivalent atom into the irreducible brillouin zone
fac = 1.0/(sym%invarind(natom)*atoms%neq(nType))
IF(sym%invarind(natom).EQ.0) CALL juDFT_error("No symmetry operations",calledby="greensfImag")
DO imat = 1, 1
DO imat = 1, MERGE(1,5,input%l_gfsphavg)
DO ie = 1, greensfCoeffs%ne
DO it = 1, sym%invarind(natom)
is = sym%invarop(natom,it)
......
......@@ -40,7 +40,7 @@ MODULE m_kkintgr
IMPLICIT NONE
!Information about the integrand
REAL, INTENT(IN) :: im(ne) !part of the green's function on the real axis
REAL, INTENT(IN) :: im(ne) !Imaginary part of the green's function on the real axis
REAL, INTENT(IN) :: eb !Bottom energy cutoff
REAL, INTENT(IN) :: del !Energy step on the real axis
INTEGER, INTENT(IN) :: ne !Number of energy points on the real axis
......
......@@ -158,7 +158,7 @@ MODULE m_hubbard1_io
!Only write the exchange splitting here if its not zero to not conflict with possible additional args
IF(exc.NE.0.0) THEN
CALL comment(input_iounit,"Exchange splitting",1)
!Sign??
!The sign flip is just a convention between the solver and the DFT calculation
CALL writeValue(input_iounit,"Exc",-exc)
ENDIF
!---------------------------------------------------------
......
......@@ -321,7 +321,7 @@ MODULE m_hubbard1_setup
WRITE (6,*) results%e_ldau
ENDIF
ELSE IF(mpi%irank.NE.0) THEN
results%e_ldau = MERGE(results%lda_u,0.0,atoms%n_u>0)
results%e_ldau = MERGE(results%e_ldau,0.0,atoms%n_u>0)
pot%mmpMat(:,:,atoms%n_u+1:atoms%n_hia+atoms%n_u,:) = CMPLX(0.0,0.0)
!If we are on a different mpi%irank and no solver is linked we need to call juDFT_end here if the solver was not run
!kind of a weird workaround (replace with something better)
......@@ -340,7 +340,7 @@ MODULE m_hubbard1_setup
!There is nothing to be done yet just set the potential correction to 0
WRITE(*,*) "No density matrix and GF found -> skipping LDA+HIA"
pot%mmpMat(:,:,atoms%n_u+1:atoms%n_hia+atoms%n_u,:) = CMPLX(0.0,0.0)
results%e_ldau = MERGE(results%lda_u,0.0,atoms%n_u>0)
results%e_ldau = MERGE(results%e_ldau,0.0,atoms%n_u>0)
ENDIF
#ifdef CPP_MPI
!Broadcast both the potential and the density matrix here
......
......@@ -441,6 +441,23 @@ CONTAINS
IF(mpi%irank.EQ.0) CALL CPP_BLAS_scopy(n,r_b,1,greensfCoeffs%ud(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,jspin),1)
ENDIF
DEALLOCATE(r_b)
IF(input%l_gfmperp.AND.jspin.EQ.1) THEN
n = greensfCoeffs%ne*atoms%n_gf*(2*lmaxU_const+1)**2
ALLOCATE(r_b(n))
CALL MPI_REDUCE(greensfCoeffs%projdos(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,3),r_b,n,CPP_MPI_REAL,MPI_SUM,0,MPI_COMM_WORLD,ierr)
IF(mpi%irank.EQ.0) CALL CPP_BLAS_scopy(n,r_b,1,greensfCoeffs%projdos(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,3),1)
IF(.NOT.input%l_gfsphavg) THEN
CALL MPI_REDUCE(greensfCoeffs%uu(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,3),r_b,n,CPP_MPI_REAL,MPI_SUM,0,MPI_COMM_WORLD,ierr)
IF(mpi%irank.EQ.0) CALL CPP_BLAS_scopy(n,r_b,1,greensfCoeffs%uu(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,3),1)
CALL MPI_REDUCE(greensfCoeffs%du(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,3),r_b,n,CPP_MPI_REAL,MPI_SUM,0,MPI_COMM_WORLD,ierr)
IF(mpi%irank.EQ.0) CALL CPP_BLAS_scopy(n,r_b,1,greensfCoeffs%du(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,3),1)
CALL MPI_REDUCE(greensfCoeffs%dd(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,3),r_b,n,CPP_MPI_REAL,MPI_SUM,0,MPI_COMM_WORLD,ierr)
IF(mpi%irank.EQ.0) CALL CPP_BLAS_scopy(n,r_b,1,greensfCoeffs%dd(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,3),1)
CALL MPI_REDUCE(greensfCoeffs%ud(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,3),r_b,n,CPP_MPI_REAL,MPI_SUM,0,MPI_COMM_WORLD,ierr)
IF(mpi%irank.EQ.0) CALL CPP_BLAS_scopy(n,r_b,1,greensfCoeffs%ud(:,:,-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,3),1)
ENDIF
DEALLOCATE(r_b)
ENDIF
ENDIF
IF(PRESENT(greensf)) THEN
n = greensf%nz*atoms%n_gf*(2*lmaxU_const+1)**2*2
......
......@@ -269,7 +269,7 @@ MODULE m_types_greensf
END IF
IF(mpi%irank.EQ.0) THEN
!Write out the information about the energy contour to outfile
!Write out the information about the energy contour
WRITE(6,"(A)") "---------------------------------------------"
WRITE(6,"(A)") " Green's function energy contour"
WRITE(6,"(A)") "---------------------------------------------"
......@@ -292,7 +292,7 @@ MODULE m_types_greensf
WRITE(6,"(A)") "Equidistant Contour for DOS calculations: "
WRITE(6,1030) this%nz, input%gf_sigma
WRITE(6,"(A)") "Energy limits (rel. to fermi energy): "
WRITE(6,1040) eb,et
WRITE(6,1040) eb-ef,et-ef
CASE default
END SELECT
......@@ -307,7 +307,7 @@ MODULE m_types_greensf
1000 FORMAT("Using energy contour mode: ", I1)
1010 FORMAT("nz: ", I5.1,"; nmatsub: ", I5.1,"; n1: ", I5.1,"; n2: ", I5.1,"; n3: ", I5.1)
1020 FORMAT("nz: ", I5.1," alpha: ", f8.4, "(not doing anything atm)")
1020 FORMAT("nz: ", I5.1," alpha: ", f8.4)
1030 FORMAT("n: ", I5.1,"; sigma: ", f8.4)
1040 FORMAT("eb: ", f8.4,"; et: ",f8.4)
1050 FORMAT(2f8.4," weight: ",2e11.4)
......@@ -393,10 +393,10 @@ MODULE m_types_greensf
DO ispin = MERGE(1,spin,l_full), MERGE(ispin_end,spin,l_full)
!Find the right quadrant in gmat according to the spin index
spin1 = MERGE(ispin,1,ispin.NE.3)
spin1 = MERGE(spin1,2,ispin.NE.4)
spin2 = MERGE(ispin,2,ispin.NE.3)
spin2 = MERGE(spin2,1,ispin.NE.4)
spin1 = MERGE(ispin,2,ispin.NE.3)
spin1 = MERGE(spin1,1,ispin.NE.4)
spin2 = MERGE(ispin,1,ispin.NE.3)
spin2 = MERGE(spin2,2,ispin.NE.4)
spin_ind = MERGE(ispin,1,input%jspins.EQ.2)
spin_ind = MERGE(3,spin_ind,ispin.EQ.4)
......
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