Commit bc9b9653 authored by Matthias Redies's avatar Matthias Redies

Merge branch 'develop' into RMA_olap

parents f4ef8549 715be85d
......@@ -4,6 +4,11 @@ if (EXISTS "${CMAKE_BINARY_DIR}/config.cmake")
include("${CMAKE_BINARY_DIR}/config.cmake")
endif()
if(${CMAKE_VERSION} VERSION_LESS "3.10.0")
message("Your cmake Version is rather old. Please consider to update cmake")
message("More modern cmake versions might be needed to make sure your compiler settings are detected correctly")
endif()
# sometimes cmake clears CMAKE_Fortran_FLAGS during project()
set(CMAKE_Fortran_FLAGS_backup ${CMAKE_Fortran_FLAGS})
set(CMAKE_Fortran_FLAGS "")
......
......@@ -156,7 +156,7 @@ SUBROUTINE cdnval(eig_id, fmpi,kpts,jspin,noco,nococonv,input,banddos,cell,atoms
CALL denCoeffs%init(atoms,sphhar,jsp_start,jsp_end)
! The last entry in denCoeffsOffdiag%init is l_fmpl. It is meant as a switch to a plot of the full magnet.
! density without the atomic sphere approximation for the magnet. density. It is not completely implemented (lo's missing).
CALL denCoeffsOffdiag%init(atoms,noco,banddos,sphhar,noco%l_mtnocopot.OR.noco%l_mperp)
CALL denCoeffsOffdiag%init(atoms,noco,sphhar,banddos%l_jDOS,noco%l_mtnocopot.OR.noco%l_mperp)
CALL force%init1(input,atoms)
CALL orb%init(atoms,noco,jsp_start,jsp_end)
......
This diff is collapsed.
......@@ -17,7 +17,7 @@ types_enparaXML.f90 types_forcetheo_data.f90
types_mpinp.f90 types_gfinp.F90 types_hub1inp.f90
types_xml.f90
calculator.f constants.f90 mpi_bc_tool.F90 ../include/inputSchema.h ../io/xml/dropInputSchema.c ../io/xml/inputSchema_old.h ../io/xml/xmlInterfaceWrapper.c
../math/d_wigner.F90 ../math/inv3.f90 ../math/grule.f90
../math/d_wigner.F90 ../math/inv3.f90 ../math/grule.f90 ../global/sort.f90
)
#Set module directories
include_directories("${CMAKE_CURRENT_BINARY_DIR}/modules/fleurinput")
......
......@@ -15,8 +15,8 @@ MODULE m_mpi_bc_tool
!have the same shape as the one on irank
INTERFACE mpi_bc
MODULE PROCEDURE mpi_bc_int,mpi_bc_int1,mpi_bc_int2,mpi_bc_int3,mpi_bc_int4,mpi_bc_int5
MODULE PROCEDURE mpi_bc_real_3,mpi_bc_real_fixed2,mpi_bc_real,mpi_bc_real1,mpi_bc_real2,mpi_bc_real3,mpi_bc_real4,mpi_bc_real5
MODULE PROCEDURE mpi_bc_complex,mpi_bc_complex1,mpi_bc_complex2,mpi_bc_complex3,mpi_bc_complex4,mpi_bc_complex5
MODULE PROCEDURE mpi_bc_real_3,mpi_bc_real_fixed2,mpi_bc_real,mpi_bc_real1,mpi_bc_real2,mpi_bc_real3,mpi_bc_real4,mpi_bc_real5,mpi_bc_real6,mpi_bc_real7
MODULE PROCEDURE mpi_bc_complex,mpi_bc_complex1,mpi_bc_complex2,mpi_bc_complex3,mpi_bc_complex4,mpi_bc_complex5,mpi_bc_complex6,mpi_bc_complex7
MODULE PROCEDURE mpi_bc_logical,mpi_bc_logical1,mpi_bc_logical2
MODULE PROCEDURE mpi_bc_character_fixed1, mpi_bc_character1, mpi_bc_character2
END INTERFACE mpi_bc
......@@ -447,6 +447,68 @@ CONTAINS
IF (ierr.NE.0) CALL judft_error("MPI_BCAST failed")
END SUBROUTINE mpi_bc_real5
SUBROUTINE mpi_bc_real6(r,irank,mpi_comm)
IMPLICIT NONE
REAL ,ALLOCATABLE,INTENT(INOUT) :: r(:,:,:,:,:,:)
INTEGER,INTENT(IN) :: irank,mpi_comm
INTEGER:: ierr=0,ilow(6),iup(6),myrank
#ifdef CPP_MPI
ilow=0;iup=0
CALL MPI_COMM_RANK(mpi_comm,myrank,ierr)
IF (ALLOCATED(r).AND.myrank==irank) THEN
ilow=LBOUND(r)
iup=UBOUND(r)
END IF
CALL MPI_BCAST(ilow,6,MPI_INTEGER,irank,mpi_comm,ierr)
CALL MPI_BCAST(iup,6,MPI_INTEGER,irank,mpi_comm,ierr)
IF (ALL(Ilow==0).AND.ALL(Iup==0)) THEN
RETURN
ENDIF
IF (myrank.NE.irank) THEN
IF (ALLOCATED(r)) DEALLOCATE(r)
ALLOCATE(r(ilow(1):iup(1),ilow(2):iup(2),ilow(3):iup(3),ilow(4):iup(4),ilow(5):iup(5),&
ilow(6):iup(6)))
ENDIF
CALL MPI_BCAST(r,SIZE(r),MPI_DOUBLE_PRECISION,irank,mpi_comm,ierr)
#endif
IF (ierr.NE.0) CALL judft_error("MPI_BCAST failed")
END SUBROUTINE mpi_bc_real6
SUBROUTINE mpi_bc_real7(r,irank,mpi_comm)
IMPLICIT NONE
REAL ,ALLOCATABLE,INTENT(INOUT) :: r(:,:,:,:,:,:,:)
INTEGER,INTENT(IN) :: irank,mpi_comm
INTEGER:: ierr=0,ilow(7),iup(7),myrank
#ifdef CPP_MPI
ilow=0;iup=0
CALL MPI_COMM_RANK(mpi_comm,myrank,ierr)
IF (ALLOCATED(r).AND.myrank==irank) THEN
ilow=LBOUND(r)
iup=UBOUND(r)
END IF
CALL MPI_BCAST(ilow,7,MPI_INTEGER,irank,mpi_comm,ierr)
CALL MPI_BCAST(iup,7,MPI_INTEGER,irank,mpi_comm,ierr)
IF (ALL(Ilow==0).AND.ALL(Iup==0)) THEN
RETURN
ENDIF
IF (myrank.NE.irank) THEN
IF (ALLOCATED(r)) DEALLOCATE(r)
ALLOCATE(r(ilow(1):iup(1),ilow(2):iup(2),ilow(3):iup(3),ilow(4):iup(4),ilow(5):iup(5),&
ilow(6):iup(6),ilow(7):iup(7)))
ENDIF
CALL MPI_BCAST(r,SIZE(r),MPI_DOUBLE_PRECISION,irank,mpi_comm,ierr)
#endif
IF (ierr.NE.0) CALL judft_error("MPI_BCAST failed")
END SUBROUTINE mpi_bc_real7
!
! And Complex!!
!
......@@ -616,6 +678,66 @@ CONTAINS
IF (ierr.NE.0) CALL judft_error("MPI_BCAST failed")
END SUBROUTINE mpi_bc_complex5
SUBROUTINE mpi_bc_complex6(c,irank,mpi_comm)
IMPLICIT NONE
COMPLEX,ALLOCATABLE,INTENT(INOUT) :: c(:,:,:,:,:,:)
INTEGER,INTENT(IN) :: irank,mpi_comm
INTEGER:: ierr=0,ilow(6),iup(6),myrank
#ifdef CPP_MPI
iup=0;ilow=0
CALL MPI_COMM_RANK(mpi_comm,myrank,ierr)
IF (ALLOCATED(c).AND.myrank==irank) THEN
ilow=LBOUND(c)
iup=UBOUND(c)
END IF
CALL MPI_BCAST(ilow,6,MPI_INTEGER,irank,mpi_comm,ierr)
CALL MPI_BCAST(iup,6,MPI_INTEGER,irank,mpi_comm,ierr)
IF (ALL(Ilow==0).AND.ALL(Iup==0)) THEN
RETURN
ENDIF
IF (myrank.NE.irank) THEN
IF (ALLOCATED(c)) DEALLOCATE(c)
ALLOCATE(c(ilow(1):iup(1),ilow(2):iup(2),ilow(3):iup(3),ilow(4):iup(4),ilow(5):iup(5),ilow(6):iup(6)))
ENDIF
CALL MPI_BCAST(c,SIZE(c),MPI_DOUBLE_COMPLEX,irank,mpi_comm,ierr)
#endif
IF (ierr.NE.0) CALL judft_error("MPI_BCAST failed")
END SUBROUTINE mpi_bc_complex6
SUBROUTINE mpi_bc_complex7(c,irank,mpi_comm)
IMPLICIT NONE
COMPLEX,ALLOCATABLE,INTENT(INOUT) :: c(:,:,:,:,:,:,:)
INTEGER,INTENT(IN) :: irank,mpi_comm
INTEGER:: ierr=0,ilow(7),iup(7),myrank
#ifdef CPP_MPI
iup=0;ilow=0
CALL MPI_COMM_RANK(mpi_comm,myrank,ierr)
IF (ALLOCATED(c).AND.myrank==irank) THEN
ilow=LBOUND(c)
iup=UBOUND(c)
END IF
CALL MPI_BCAST(ilow,7,MPI_INTEGER,irank,mpi_comm,ierr)
CALL MPI_BCAST(iup,7,MPI_INTEGER,irank,mpi_comm,ierr)
IF (ALL(Ilow==0).AND.ALL(Iup==0)) THEN
RETURN
ENDIF
IF (myrank.NE.irank) THEN
IF (ALLOCATED(c)) DEALLOCATE(c)
ALLOCATE(c(ilow(1):iup(1),ilow(2):iup(2),ilow(3):iup(3),ilow(4):iup(4),ilow(5):iup(5),&
ilow(6):iup(6),ilow(7):iup(7)))
ENDIF
CALL MPI_BCAST(c,SIZE(c),MPI_DOUBLE_COMPLEX,irank,mpi_comm,ierr)
#endif
IF (ierr.NE.0) CALL judft_error("MPI_BCAST failed")
END SUBROUTINE mpi_bc_complex7
SUBROUTINE mpi_bc_character_fixed1(irank,mpi_comm,c)
IMPLICIT NONE
......
This diff is collapsed.
......@@ -251,7 +251,7 @@ CONTAINS
ENDIF
ELSE
DO n = 1, kpts%numSpecialPoints
WRITE (fh, 209) TRIM(ADJUSTL(kpts%specialPointNames(n))), kpts%specialPoints(:, n)
WRITE (fh, 211) TRIM(ADJUSTL(kpts%specialPointNames(n))), kpts%specialPoints(:, n)
211 FORMAT(' <specialPoint name="', a, '">', f10.6, ' ', f10.6, ' ', f10.6, '</specialPoint>')
END DO
END IF
......
......@@ -19,8 +19,8 @@ MODULE m_types_noco
LOGICAL:: l_alignMT = .FALSE.
LOGICAL:: l_sourceFree = .FALSE.
LOGICAL:: l_scaleMag = .FALSE.
LOGICAL :: l_RelaxAlpha=.FALSE.
LOGICAL :: l_RelaxBeta=.FALSE.
LOGICAL:: l_RelaxAlpha=.FALSE.
LOGICAL:: l_RelaxBeta=.FALSE.
REAL :: mag_scale=1.0
REAL :: mix_b=0.0
REAL :: mix_RelaxWeightOffD=1.0
......
......@@ -30,10 +30,10 @@ MODULE m_greensfBZint
INTEGER :: i_gf,l,lp,atomType,atomTypep,indUnique
INTEGER :: natom,natomp,natomp_start,natomp_end,natom_start,natom_end
INTEGER :: i_elem
INTEGER :: i_elem,i_elemLO,nLO,imatSize
INTEGER :: spin1,spin2,ispin,spin_start,spin_end
COMPLEX :: phase
REAL :: atomFactor
REAL :: atomFactor,atomDiff(3)
LOGICAL :: l_sphavg
COMPLEX, ALLOCATABLE :: im(:,:,:,:,:)
......@@ -53,15 +53,26 @@ MODULE m_greensfBZint
atomType = gfinp%elem(i_gf)%atomType
atomTypep = gfinp%elem(i_gf)%atomTypep
l_sphavg = gfinp%elem(i_gf)%l_sphavg
atomDiff(:) = gfinp%elem(i_gf)%atomDiff(:)
atomFactor = MERGE(1.0,1.0/atoms%neq(atomType),l.NE.lp)
atomFactor = MERGE(1.0,atomFactor,atomType.NE.atomTypep)
i_elem = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg,indUnique=indUnique)
i_elem = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg,indUnique=indUnique)
i_elemLO = gfinp%uniqueElements(atoms,ind=i_gf,lo=.TRUE.,l_sphavg=l_sphavg,indUnique=indUnique)
IF(i_gf/=indUnique) CYCLE
nLO = 0
imatSize = 1
IF(.NOT.l_sphavg) THEN
imatSize = 4
nLO = gfinp%elem(i_gf)%countLOs(atoms)
IF(nLO/=0) THEN
imatSize = 4+4*nLO+nLO**2
ENDIF
ENDIF
ALLOCATE(im(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,nBands,&
MERGE(1,4,l_sphavg),spin_start:spin_end),source=cmplx_0)
imatSize,spin_start:spin_end),source=cmplx_0)
natom_start = SUM(atoms%neq(:atomType-1)) + 1
natom_end = MERGE(SUM(atoms%neq(:atomType-1)) + 1,SUM(atoms%neq(:atomType)),l.NE.lp)
......@@ -79,8 +90,8 @@ MODULE m_greensfBZint
DO natomp = natomp_start, natomp_end
!Phase factor for intersite elements (Does nothing atm)
IF(natom.NE.natomp) THEN
phase = exp(ImagUnit*dot_product(kpts%bk(:,ikpt),gfinp%elem(i_gf)%atomDiff))
IF(ANY(ABS(atomDiff).GT.1e-12)) THEN
phase = exp(ImagUnit*dot_product(kpts%bk(:,ikpt),atomDiff(:)))
ELSE
phase = cmplx_1
ENDIF
......@@ -108,8 +119,8 @@ MODULE m_greensfBZint
l_sphavg,atoms,denCoeffsOffdiag,eigVecCoeffs,im(:,:,:,:,ispin))
ENDIF
CALL greensfSym(ikpt_i,i_elem,natom,l,natom.EQ.natomp.AND.l.EQ.lp,l_sphavg,&
ispin,sym,atomFactor,phase,im(:,:,:,:,ispin),greensfBZintCoeffs)
CALL greensfSym(ikpt_i,i_elem,i_elemLO,nLO,natom,l,natom.EQ.natomp.AND.l.EQ.lp.AND.ALL(ABS(atomDiff).LT.1e-12),&
l_sphavg,ispin,sym,atomFactor,phase,im(:,:,:,:,ispin),greensfBZintCoeffs)
ENDDO
ENDDO !natomp
......
......@@ -26,9 +26,9 @@ MODULE m_greensfCalcImagPart
#include"cpp_double.h"
INTEGER :: ikpt_i,ikpt,nBands,jsp,i_gf
INTEGER :: ikpt_i,ikpt,nBands,jsp,i_gf,nLO,imatSize
INTEGER :: l,lp,m,mp,iBand,ie,j,eGrid_start,eGrid_end
INTEGER :: indUnique,i_elem
INTEGER :: indUnique,i_elem,imat,iLO,iLOp,i_elemLO
LOGICAL :: l_zero,l_sphavg
REAL :: del,eb,wtkpt
COMPLEX :: fac,weight
......@@ -89,16 +89,27 @@ MODULE m_greensfCalcImagPart
l_sphavg = gfinp%elem(i_gf)%l_sphavg
i_elem = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg,indUnique=indUnique)
i_elemLO = gfinp%uniqueElements(atoms,ind=i_gf,lo=.TRUE.,l_sphavg=l_sphavg,indUnique=indUnique)
IF(i_gf/=indUnique) CYCLE
nLO = 0
imatSize = 1
IF(.NOT.l_sphavg) THEN
imatSize = 4
nLO = gfinp%elem(i_gf)%countLOs(atoms)
IF(nLO/=0) THEN
imatSize = 4+4*nLO+nLO**2
ENDIF
ENDIF
!$OMP parallel default(none) &
!$OMP shared(input,greensfBZintCoeffs,greensfImagPart) &
!$OMP shared(i_elem,l,lp,ikpt_i,nBands,eMesh,l_sphavg)&
!$OMP shared(i_elem,i_elemLO,nLO,l,lp,ikpt_i,nBands,eMesh,l_sphavg,imatSize)&
!$OMP shared(del,eb,eig,weights,indBound,fac,wtkpt,spin_ind) &
!$OMP private(ie,m,mp,iBand,j,eGrid_start,eGrid_end,weight,imag,imagReal,l_zero)
ALLOCATE(imag(SIZE(eMesh),MERGE(1,4,l_sphavg)),source=cmplx_0)
ALLOCATE(imagReal(SIZE(eMesh),MERGE(1,4,l_sphavg)),source=0.0)
!$OMP private(ie,iLO,iLOp,imat,m,mp,iBand,j,eGrid_start,eGrid_end,weight,imag,imagReal,l_zero)
ALLOCATE(imag(SIZE(eMesh),imatSize),source=cmplx_0)
ALLOCATE(imagReal(SIZE(eMesh),imatSize),source=0.0)
!$OMP do collapse(2)
DO mp = -lp, lp
DO m = -l, l
......@@ -140,6 +151,23 @@ MODULE m_greensfCalcImagPart
imag(ie,2) = imag(ie,2) + weight * greensfBZintCoeffs%dd(iBand,m,mp,i_elem,ikpt_i,spin_ind)
imag(ie,3) = imag(ie,3) + weight * greensfBZintCoeffs%ud(iBand,m,mp,i_elem,ikpt_i,spin_ind)
imag(ie,4) = imag(ie,4) + weight * greensfBZintCoeffs%du(iBand,m,mp,i_elem,ikpt_i,spin_ind)
IF(nLO>0) THEN
imat = 0
DO iLO = 1, nLO
imat = imat + 4
imag(ie,imat+1) = imag(ie,imat+1) + weight * greensfBZintCoeffs%uulo(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind)
imag(ie,imat+2) = imag(ie,imat+2) + weight * greensfBZintCoeffs%ulou(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind)
imag(ie,imat+3) = imag(ie,imat+3) + weight * greensfBZintCoeffs%dulo(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind)
imag(ie,imat+4) = imag(ie,imat+4) + weight * greensfBZintCoeffs%ulod(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind)
ENDDO
imat = 0
DO iLO = 1, nLO
DO iLOp = 1, nLO
imat = imat + 1
imag(ie,4 + 4*nLO+imat) = imag(ie,4 + 4*nLO+imat) + weight * greensfBZintCoeffs%uloulop(iBand,m,mp,imat,i_elemLO,ikpt_i,spin_ind)
ENDDO
ENDDO
ENDIF
ENDIF
ELSE IF(eGrid_start==1 .AND. eGrid_end==SIZE(eMesh)) THEN!Here we always use the tetrahedron method
......@@ -156,6 +184,28 @@ MODULE m_greensfCalcImagPart
weights(:,iBand),1,imag(:,3),1)
CALL CPP_BLAS_caxpy(SIZE(eMesh),greensfBZintCoeffs%du(iBand,m,mp,i_elem,ikpt_i,spin_ind),&
weights(:,iBand),1,imag(:,4),1)
IF(nLO>0) THEN
imat = 0
DO iLO = 1, nLO
imat = imat + 4
CALL CPP_BLAS_caxpy(SIZE(eMesh),greensfBZintCoeffs%uulo(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind),&
weights(:,iBand),1,imag(:,imat+1),1)
CALL CPP_BLAS_caxpy(SIZE(eMesh),greensfBZintCoeffs%ulou(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind),&
weights(:,iBand),1,imag(:,imat+2),1)
CALL CPP_BLAS_caxpy(SIZE(eMesh),greensfBZintCoeffs%dulo(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind),&
weights(:,iBand),1,imag(:,imat+3),1)
CALL CPP_BLAS_caxpy(SIZE(eMesh),greensfBZintCoeffs%ulod(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind),&
weights(:,iBand),1,imag(:,imat+4),1)
ENDDO
imat = 0
DO iLO = 1, nLO
DO iLOp = 1, nLO
imat = imat + 1
CALL CPP_BLAS_caxpy(SIZE(eMesh),greensfBZintCoeffs%uloulop(iBand,m,mp,imat,i_elemLO,ikpt_i,spin_ind),&
weights(:,iBand),1,imag(:,4 + 4*nLO+imat),1)
ENDDO
ENDDO
ENDIF
ENDIF
ELSE
IF(l_sphavg) THEN
......@@ -170,6 +220,28 @@ MODULE m_greensfCalcImagPart
* greensfBZintCoeffs%ud(iBand,m,mp,i_elem,ikpt_i,spin_ind)
imag(eGrid_start:eGrid_end,4) = imag(eGrid_start:eGrid_end,4) + weights(eGrid_start:eGrid_end,iBand)&
* greensfBZintCoeffs%du(iBand,m,mp,i_elem,ikpt_i,spin_ind)
IF(nLO>0) THEN
imat = 0
DO iLO = 1, nLO
imat = imat + 4
imag(eGrid_start:eGrid_end,imat+1) = imag(eGrid_start:eGrid_end,imat+1) + weights(eGrid_start:eGrid_end,iBand) &
* greensfBZintCoeffs%uulo(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind)
imag(eGrid_start:eGrid_end,imat+2) = imag(eGrid_start:eGrid_end,imat+2) + weights(eGrid_start:eGrid_end,iBand) &
* greensfBZintCoeffs%ulou(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind)
imag(eGrid_start:eGrid_end,imat+3) = imag(eGrid_start:eGrid_end,imat+3) + weights(eGrid_start:eGrid_end,iBand) &
* greensfBZintCoeffs%dulo(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind)
imag(eGrid_start:eGrid_end,imat+4) = imag(eGrid_start:eGrid_end,imat+4) + weights(eGrid_start:eGrid_end,iBand) &
* greensfBZintCoeffs%ulod(iBand,m,mp,iLO,i_elemLO,ikpt_i,spin_ind)
ENDDO
imat = 0
DO iLO = 1, nLO
DO iLOp = 1, nLO
imat = imat + 1
imag(eGrid_start:eGrid_end,4 + 4*nLO+imat) = imag(eGrid_start:eGrid_end,4 + 4*nLO+imat) + weights(eGrid_start:eGrid_end,iBand) &
* greensfBZintCoeffs%uloulop(iBand,m,mp,imat,i_elemLO,ikpt_i,spin_ind)
ENDDO
ENDDO
ENDIF
ENDIF
ENDIF
......@@ -182,11 +254,30 @@ MODULE m_greensfCalcImagPart
CALL CPP_BLAS_saxpy(SIZE(eMesh),1.0,imagReal(:,2),1,greensfImagPart%dd(:,m,mp,i_elem,spin_ind),1)
CALL CPP_BLAS_saxpy(SIZE(eMesh),1.0,imagReal(:,3),1,greensfImagPart%ud(:,m,mp,i_elem,spin_ind),1)
CALL CPP_BLAS_saxpy(SIZE(eMesh),1.0,imagReal(:,4),1,greensfImagPart%du(:,m,mp,i_elem,spin_ind),1)
IF(nLO>0) THEN
imat = 0
DO iLO = 1, nLO
imat = imat + 4
CALL CPP_BLAS_saxpy(SIZE(eMesh),1.0,imagReal(:,imat+1),1,greensfImagPart%uulo(:,m,mp,iLO,i_elemLO,spin_ind),1)
CALL CPP_BLAS_saxpy(SIZE(eMesh),1.0,imagReal(:,imat+2),1,greensfImagPart%ulou(:,m,mp,iLO,i_elemLO,spin_ind),1)
CALL CPP_BLAS_saxpy(SIZE(eMesh),1.0,imagReal(:,imat+3),1,greensfImagPart%dulo(:,m,mp,iLO,i_elemLO,spin_ind),1)
CALL CPP_BLAS_saxpy(SIZE(eMesh),1.0,imagReal(:,imat+4),1,greensfImagPart%ulod(:,m,mp,iLO,i_elemLO,spin_ind),1)
ENDDO
imat = 0
DO iLO = 1, nLO
DO iLOp = 1, nLO
imat = imat + 1
CALL CPP_BLAS_saxpy(SIZE(eMesh),1.0,imagReal(:,4 + 4*nLO+imat),1,greensfImagPart%uloulop(:,m,mp,iLO,iLOp,i_elemLO,spin_ind),1)
ENDDO
ENDDO
ENDIF
ENDIF
ENDDO!m
ENDDO!mp
!$OMP end do
DEALLOCATE(imagReal,imag)
!$OMP end parallel
ENDDO!i_gf
CALL timestop("Green's Function: Imaginary Part")
......
......@@ -27,28 +27,27 @@ MODULE m_greensfCalcRealPart
CONTAINS
SUBROUTINE greensfCalcRealPart(atoms,gfinp,input,sym,noco,vTot,enpara,fmpi,hub1inp,ef,greensfImagPart,g)
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_gfinp), INTENT(IN) :: gfinp
TYPE(t_sym), INTENT(IN) :: sym
TYPE(t_noco), INTENT(IN) :: noco
TYPE(t_input), INTENT(IN) :: input
TYPE(t_potden), INTENT(IN) :: vTot
TYPE(t_enpara), INTENT(IN) :: enpara
TYPE(t_hub1inp), INTENT(IN) :: hub1inp
TYPE(t_mpi), INTENT(IN) :: fmpi
REAL, INTENT(IN) :: ef
TYPE(t_greensfImagPart),INTENT(INOUT) :: greensfImagPart
TYPE(t_greensf), INTENT(INOUT) :: g(:)
INTEGER :: i_gf,i_elem,ie,l,m,mp,nType,indUnique
SUBROUTINE greensfCalcRealPart(atoms,gfinp,input,sym,noco,usdus,denCoeffsOffDiag,fmpi,ef,greensfImagPart,g)
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_gfinp), INTENT(IN) :: gfinp
TYPE(t_sym), INTENT(IN) :: sym
TYPE(t_noco), INTENT(IN) :: noco
TYPE(t_usdus), INTENT(IN) :: usdus
TYPE(t_denCoeffsOffDiag), INTENT(IN) :: denCoeffsOffDiag
TYPE(t_input), INTENT(IN) :: input
TYPE(t_mpi), INTENT(IN) :: fmpi
REAL, INTENT(IN) :: ef
TYPE(t_greensfImagPart), INTENT(INOUT) :: greensfImagPart
TYPE(t_greensf), INTENT(INOUT) :: g(:)
INTEGER :: i_gf,i_elem,ie,l,m,mp,nType,indUnique,nLO,iLO,iLOp,i_elemLO
INTEGER :: jspin,nspins,ipm,kkcut,lp,nTypep,refCutoff
INTEGER :: spin_cut,contourShape
INTEGER :: i_gf_start,i_gf_end,spin_start,spin_end
INTEGER :: n_gf_task,extra
LOGICAL :: l_onsite,l_fixedCutoffset,l_sphavg
REAL :: fac,del,eb,et,fixedCutoff
REAL :: fac,del,eb,et,fixedCutoff,atomDiff(3)
REAL, ALLOCATABLE :: eMesh(:),imag(:)
!Get the information on the real axis energy mesh
......@@ -69,6 +68,7 @@ MODULE m_greensfCalcRealPart
l_fixedCutoffset = g(i_gf)%elem%l_fixedCutoffset
fixedCutoff = g(i_gf)%elem%fixedCutoff
refCutoff = g(i_gf)%elem%refCutoff
atomDiff(:) = g(i_gf)%elem%atomDiff(:)
i_elem = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg,indUnique=indUnique)
......@@ -78,8 +78,8 @@ MODULE m_greensfCalcRealPart
greensfImagPart%kkintgr_cutoff(i_gf,:,:) = greensfImagPart%kkintgr_cutoff(indUnique,:,:)
ELSE
!Is the current element suitable for automatic finding of the cutoff
l_onsite = nType.EQ.nTypep.AND.l.EQ.lp
IF(l_onsite.AND..NOT.l_fixedCutoffset.AND.refCutoff==-1 .AND.gfinp%countLOs(atoms,i_gf)==0) THEN
l_onsite = nType.EQ.nTypep.AND.l.EQ.lp.AND.ALL(ABS(atomDiff(:)).LT.1e-12)
IF(l_onsite.AND..NOT.l_fixedCutoffset.AND.refCutoff==-1 .AND. g(i_gf)%elem%countLOs(atoms)==0) THEN
!
!Check the integral over the fDOS to define a cutoff for the Kramer-Kronigs-Integration
! with LOs I just use a fixed cutoff or reference otherwise I would need to check whether
......@@ -91,7 +91,7 @@ MODULE m_greensfCalcRealPart
!Onsite element with radial dependence
CALL kk_cutoffRadial(greensfImagPart%uu(:,:,:,i_elem,:),greensfImagPart%ud(:,:,:,i_elem,:),&
greensfImagPart%du(:,:,:,i_elem,:),greensfImagPart%dd(:,:,:,i_elem,:),&
noco,atoms,vTot,enpara,fmpi,hub1inp,gfinp%l_mperp,l,nType,input,eMesh,&
noco,usdus,denCoeffsOffDiag,gfinp%l_mperp,l,nType,input,eMesh,&
greensfImagPart%kkintgr_cutoff(i_gf,:,:),greensfImagPart%scalingFactorRadial(i_elem,:))
ENDIF
ELSE IF (l_fixedCutoffset) THEN
......@@ -117,12 +117,14 @@ MODULE m_greensfCalcRealPart
refCutoff = g(i_gf)%elem%refCutoff
l_sphavg = g(i_gf)%elem%l_sphavg
i_elem = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg,indUnique=indUnique)
i_elemLO = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg,indUnique=indUnique,lo=.TRUE.)
nLO = g(i_gf)%elem%countLOs(atoms)
IF(refCutoff/=-1) THEN
!Overwrite cutoff with reference from other elements
greensfImagPart%kkintgr_cutoff(i_gf,:,:) = greensfImagPart%kkintgr_cutoff(refCutoff,:,:)
ENDIF
CALL greensfImagPart%scale(i_elem,l_sphavg)
CALL greensfImagPart%scale(i_elem,i_elemLO,l_sphavg,nLO)
ENDDO
CALL timestop("Green's Function: Integration Cutoff")
ENDIF
......@@ -184,8 +186,10 @@ MODULE m_greensfCalcRealPart
nTypep = g(i_gf)%elem%atomTypep
l_sphavg = g(i_gf)%elem%l_sphavg
contourShape = gfinp%contour(g(i_gf)%elem%iContour)%shape
nLO = g(i_gf)%elem%countLOs(atoms)
i_elem = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg)
i_elemLO = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg,lo=.TRUE.)
CALL timestart("Green's Function: Kramer-Kronigs-Integration")
DO jspin = spin_start, spin_end
......@@ -193,7 +197,7 @@ MODULE m_greensfCalcRealPart
DO m= -l,l
DO mp= -lp,lp
IF(greensfImagPart%checkEmpty(i_elem,m,mp,jspin,l_sphavg)) THEN
IF(greensfImagPart%checkEmpty(i_elem,i_elemLO,nLO,m,mp,jspin,l_sphavg)) THEN
CALL g(i_gf)%resetSingleElem(m,mp,jspin,ipm)
CYCLE
ENDIF
......@@ -217,6 +221,30 @@ MODULE m_greensfCalcRealPart
imag = greensfImagPart%applyCutoff(i_elem,i_gf,m,mp,jspin,l_sphavg,imat=4)
CALL kkintgr(imag,eMesh,g(i_gf)%contour%e,(ipm.EQ.2),&
g(i_gf)%du(:,m,mp,jspin,ipm),int_method(contourShape))
!KKT for LOs
IF(nLO>0) THEN
DO iLO = 1, nLO
imag = greensfImagPart%applyCutoff(i_elemLO,i_gf,m,mp,jspin,l_sphavg,imat=1,iLO=iLO)
CALL kkintgr(imag,eMesh,g(i_gf)%contour%e,(ipm.EQ.2),&
g(i_gf)%uulo(:,m,mp,iLO,jspin,ipm),int_method(contourShape))
imag = greensfImagPart%applyCutoff(i_elemLO,i_gf,m,mp,jspin,l_sphavg,imat=2,iLO=iLO)
CALL kkintgr(imag,eMesh,g(i_gf)%contour%e,(ipm.EQ.2),&
g(i_gf)%ulou(:,m,mp,iLO,jspin,ipm),int_method(contourShape))
imag = greensfImagPart%applyCutoff(i_elemLO,i_gf,m,mp,jspin,l_sphavg,imat=3,iLO=iLO)
CALL kkintgr(imag,eMesh,g(i_gf)%contour%e,(ipm.EQ.2),&
g(i_gf)%dulo(:,m,mp,iLO,jspin,ipm),int_method(contourShape))
imag = greensfImagPart%applyCutoff(i_elemLO,i_gf,m,mp,jspin,l_sphavg,imat=4,iLO=iLO)
CALL kkintgr(imag,eMesh,g(i_gf)%contour%e,(ipm.EQ.2),&
g(i_gf)%ulod(:,m,mp,iLO,jspin,ipm),int_method(contourShape))
DO iLOp = 1, nLO
imag = greensfImagPart%applyCutoff(i_elemLO,i_gf,m,mp,jspin,l_sphavg,iLO=iLO,iLOp=iLop)
CALL kkintgr(imag,eMesh,g(i_gf)%contour%e,(ipm.EQ.2),&
g(i_gf)%uloulop(:,m,mp,iLO,iLOp,jspin,ipm),int_method(contourShape))
ENDDO
ENDDO
ENDIF
ENDIF
ENDDO
ENDDO
......
This diff is collapsed.
......@@ -19,15 +19,19 @@ MODULE m_greensfSpinDiag
TYPE(t_eigVecCoeffs), INTENT(IN) :: eigVecCoeffs
COMPLEX, INTENT(INOUT) :: im(-lmaxU_const:,-lmaxU_const:,:,:)
INTEGER :: m,mp,lm,lmp,ilo,ilop
INTEGER :: m,mp,lm,lmp,ilo,ilop,nLO_ind,nLOp_ind,imat