Commit d0536308 authored by Matthias Redies's avatar Matthias Redies

Merge branch 'develop' into RMA_olap

parents 4431b218 dba22843
......@@ -163,7 +163,7 @@ SUBROUTINE cdnval(eig_id, fmpi,kpts,jspin,noco,nococonv,input,banddos,cell,atoms
!Greens function always considers the empty states
IF(gfinp%n>0 .AND. PRESENT(greensfImagPart)) THEN
IF(greensfImagPart%l_calc) &
CALL greensfBZintCoeffs%init(gfinp,input,noco,jsp_start,jsp_end,&
CALL greensfBZintCoeffs%init(gfinp,atoms,input,noco,jsp_start,jsp_end,&
SIZE(cdnvalJob%k_list),SIZE(cdnvalJob%ev_list))
ENDIF
......
......@@ -51,11 +51,11 @@ SUBROUTINE rotateMagnetToSpinAxis(vacuum,sphhar,stars&
IF(noco%alph_inp(i).NE.0.0) nonZeroAngles=.TRUE.
IF(noco%beta_inp(i).NE.0.0) nonZeroAngles=.TRUE.
END DO
IF(l_firstIt.AND.nonZeroAngles) THEN
! Rotates cdn by given noco angles in first iteration. WARNING: If you want to continue/restart a calculation with MT relaxation set noco angles to 0!
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,zeros,nococonv%beta,den)
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,nococonv%alph,zeros,den)
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,zeros,noco%beta_inp,den)
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,noco%alph_inp,zeros,den)
!Setting angles to zero since we want our spinQ Axis to remain the same when writing out the density.
nococonv%alph=zeros
nococonv%beta=zeros
......@@ -79,17 +79,33 @@ SUBROUTINE rotateMagnetToSpinAxis(vacuum,sphhar,stars&
IF (abs(diffT(i)).LE.eps) diffT(i)=0.0
IF (abs(diffP(i)).LE.eps) diffP(i)=0.0
END DO
!Mixing of angles part
IF(noco%mix_RelaxAngle.NE.1.0) THEN
diffT=noco%mix_RelaxAngle*diffT
diffP=noco%mix_RelaxAngle*diffP
!Which angles to keep constant from input den
IF (.NOT.noco%l_RelaxAlpha) THEN
IF(l_firstIT) THEN
nococonv%alphRlx=diffP
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,-diffP,zeros,den)
END IF
END IF
! Rotate cdn by direction of magnetization so it alings with spin quantization axis.
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,-diffP,zeros,den)
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,zeros,-diffT,den)
IF (noco%l_RelaxAlpha.OR..NOT.l_firstIT) CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,-diffP,zeros,den)
!Which angles to keep constant from input den
IF (.NOT.noco%l_RelaxBeta) THEN
IF(l_firstIT) THEN
nococonv%betaRlx=diffT
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,zeros,-diffT,den)
END IF
END IF
IF (noco%l_RelaxBeta.OR..NOT.l_firstIT) CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,zeros,-diffT,den)
!Store angles in nococonv
nococonv%beta=nococonv%beta+diffT
nococonv%alph=nococonv%alph+diffP
IF (noco%l_RelaxBeta) nococonv%beta=nococonv%beta+diffT
IF (.NOT.noco%l_RelaxBeta) nococonv%beta=nococonv%betaRlx
IF (noco%l_RelaxAlpha) nococonv%alph=nococonv%alph+diffP
IF(.NOT.noco%l_RelaxAlpha) nococonv%alph=nococonv%alphRlx
! Set angles to zero if too low. (Prevent numerical rubbish to appear)
DO i=1, atoms%ntype
IF (abs(nococonv%beta(i)).LE.eps) nococonv%beta(i)=0.0
......@@ -126,12 +142,12 @@ SUBROUTINE rotateMagnetFromSpinAxis(noco,nococonv,vacuum,sphhar,stars&
CAlL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,zeros,nococonv%beta,inDen)
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,nococonv%alph,zeros,inDen)
IF (present(den)) THEN
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,zeros,nococonv%beta,den)
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,nococonv%alph,zeros,den)
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,zeros,nococonv%beta,den)
CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco,oneD,cell,nococonv%alph,zeros,den)
END IF
! Nococonv is zero now since rotation has been reverted.
nococonv%alph=zeros
nococonv%beta=zeros
nococonv%alph=zeros
nococonv%beta=zeros
END SUBROUTINE rotateMagnetFromSpinAxis
......
......@@ -2,9 +2,22 @@ if module list 2>&1 | grep -q mpt
then
export FC=mpif90
export CC=mpicc
export CXX=mpicxx
#FLEUR_LIBRARIES="-lxml2;-lmkl_scalapack_lp64;-lmkl_blacs_sgimpt_lp64"
FLEUR_LIBRARIES="-lxml2"
elif module list 2>&1 | grep -q openmpi
then
export FC=mpif90
export CC=mpicc
FLEUR_LIBRARIES="-lxml2;-lmkl_scalapack_lp64;-lmkl_blacs_openmpi_lp64"
else
#export FC=mpif90
#export CC=mpicc
export FC=mpiifort
export CC=mpiicc
FLEUR_LIBRARIES="-lxml2;-lmkl_scalapack_lp64;-lmkl_blacs_intelmpi_lp64"
fi
FLEUR_LIBRARIES="-lxml2;-lmkl_scalapack_lp64;-lmkl_blacs_sgimpt_lp64"
echo "Wannier is switched off-manually"
CLI_USE_WANNIER=FALSE
echo "ChASE is switched off-manually"
......
......@@ -77,7 +77,7 @@ if (DEFINED CLI_FLEUR_USE_EDSOLVER)
endif()
if (FLEUR_USE_ARPACK)
set(EDSOLVER_BUILD_FLEUR ON CACHE BOOL "Build EDsolver with juDFT")
set(EDSOLVER_BUILD_UTILS ON CACHE BOOL "Build auxiliary programs for the EDsolver")
set(EDSOLVER_BUILD_UTILS OFF CACHE BOOL "Build auxiliary programs for the EDsolver")
if (CMAKE_BUILD_TYPE MATCHES Debug)
set(EDSOLVER_DEBUG_BUILD ON CACHE BOOL "Build a debug version")
endif()
......
......@@ -268,19 +268,7 @@ CONTAINS
WRITE(attributes(1),'(f20.10)') results%ef
WRITE(attributes(2),'(a)') 'Htr'
IF (fmpi%irank.EQ.0) CALL writeXMLElement('FermiEnergy',(/'value','units'/),attributes(1:2))
ENDIF
!IF(.not.input%eig66(1))THEN
!Put w_iks into eig-file
DO jsp = 1,nspins
DO k = 1,kpts%nkpt
IF (fmpi%irank == 0) CALL write_eig(eig_id,k,jsp,w_iks=results%w_iks(:,k,jsp))
#ifdef CPP_MPI
CALL MPI_BARRIER(fmpi%mpi_comm,ierr)
#endif
ENDDO
ENDDO
!ENDIF
ENDIF
RETURN
8020 FORMAT (/,'FERMIE:',/,&
......
......@@ -78,7 +78,6 @@ MODULE m_fertetra
!Now that the fermi energy is guaranteed to be in the interval [lowBound,upperBound]
!We use the bisection method to find it
!-----------------------------------------------------------------------------------
DO WHILE(upperBound-lowBound.GT.1e-10)
ef = (lowBound+upperBound)/2.0
......@@ -105,6 +104,7 @@ MODULE m_fertetra
!---------------------------------------------------------------------
!Calculate final occupation and weights for individual kpts and bands
!---------------------------------------------------------------------
ef = (lowBound+upperBound)/2.0
dfermi = 0.0
DO jspin = 1, jspins
!-------------------------------------------------------
......
......@@ -312,6 +312,7 @@ SUBROUTINE read_xml_atoms(this,xml)
xpath=''
IF(xml%getNumberOfNodes(TRIM(ADJUSTL(xPaths))//'/force')==1) xpath=xpaths
IF(xml%getNumberOfNodes(TRIM(ADJUSTL(xPathg))//'/force')==1) xpath=xpathg
this%relax(:,n) = 0
IF (xpath.NE.'') THEN
this%l_geo(n) = evaluateFirstBoolOnly(xml%getAttributeValue(TRIM(ADJUSTL(xPath))//'/force/@calculate'))
valueString = xml%getAttributeValue(TRIM(ADJUSTL(xPath))//'/force/@relaxXYZ')
......@@ -321,7 +322,6 @@ SUBROUTINE read_xml_atoms(this,xml)
IF (relaxZ) this%relax(3,n) = 1
ELSE
this%l_geo(n) = .FALSE.
this%relax(:,n) = 0
END IF
!LO stuff
this%nlo(n) = 0
......
This diff is collapsed.
......@@ -19,9 +19,10 @@ MODULE m_types_noco
LOGICAL:: l_alignMT = .FALSE.
LOGICAL:: l_sourceFree = .FALSE.
LOGICAL:: l_scaleMag = .FALSE.
LOGICAL :: l_RelaxAlpha=.FALSE.
LOGICAL :: l_RelaxBeta=.FALSE.
REAL :: mag_scale=1.0
REAL :: mix_b=0.0
REAL :: mix_RelaxAngle=1.0
REAL :: mix_RelaxWeightOffD=1.0
LOGICAL:: l_spav= .FALSE.
REAL :: theta_inp=0.0
......@@ -65,12 +66,12 @@ MODULE m_types_noco
CALL mpi_bc(this%mag_scale ,rank,mpi_comm)
CALL mpi_bc(rank,mpi_comm,this%qss_inp)
CALL mpi_bc(this%mix_b,rank,mpi_comm)
CALL mpi_bc(this%mix_RelaxAngle,rank,mpi_comm)
CALL mpi_bc(this%l_RelaxAlpha,rank,mpi_comm)
CALL mpi_bc(this%l_RelaxBeta,rank,mpi_comm)
CALL mpi_bc(this%mix_RelaxWeightOffD,rank,mpi_comm)
CALL mpi_bc(this%l_spav,rank,mpi_comm)
CALL mpi_bc(this%theta_inp,rank,mpi_comm)
CALL mpi_bc(this%phi_inp,rank,mpi_comm)
CALL mpi_bc(this%l_relax,rank,mpi_comm)
CALL mpi_bc(this%alph_inp,rank,mpi_comm)
CALL mpi_bc(this%beta_inp,rank,mpi_comm)
......@@ -112,12 +113,13 @@ MODULE m_types_noco
this%l_mperp = evaluateFirstBoolOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@l_mperp'))
this%l_constr = evaluateFirstBoolOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@l_constr'))
this%l_mtNocoPot = evaluateFirstBoolOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@l_mtNocoPot'))
this%l_alignMT = evaluateFirstBoolOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@l_alignMT'))
this%l_sourceFree = evaluateFirstBoolOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@l_sourceFree'))
this%l_scaleMag = evaluateFirstBoolOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@l_scaleMag'))
this%mag_scale = evaluateFirstOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@mag_scale'))
this%mix_b = evaluateFirstOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@mix_b'))
this%mix_RelaxAngle = evaluateFirstOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@mix_RelaxAngle'))
this%l_alignMT = evaluateFirstBoolOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@l_RelaxMT'))
this%l_RelaxBeta = evaluateFirstBoolOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@l_RelaxBeta'))
this%l_RelaxAlpha = evaluateFirstBoolOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@l_RelaxAlpha'))
this%mix_RelaxWeightOffD = evaluateFirstOnly(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@mix_RelaxWeightOffD'))
valueString = TRIM(ADJUSTL(xml%GetAttributeValue(TRIM(ADJUSTL(xPathA))//'/qss')))
READ(valueString,*) this%qss_inp
......
......@@ -7,7 +7,7 @@ MODULE m_genMTBasis
CONTAINS
SUBROUTINE genMTBasis(atoms,enpara,vTot,fmpi,iType,jspin,usdus,f,g,flo,l_dftspinpol)
SUBROUTINE genMTBasis(atoms,enpara,vTot,fmpi,iType,jspin,usdus,f,g,flo,l_dftspinpol,l_writeArg)
USE m_types
USE m_constants
USE m_radfun
......@@ -30,6 +30,7 @@ CONTAINS
REAL, INTENT(INOUT) :: f(atoms%jmtd,2,0:atoms%lmaxd)
REAL, INTENT(INOUT) :: g(atoms%jmtd,2,0:atoms%lmaxd)
REAL, INTENT(INOUT) :: flo(atoms%jmtd,2,atoms%nlod)
LOGICAL, OPTIONAL, INTENT(IN) :: l_writeArg
INTEGER :: l,nodeu,noded
......@@ -39,7 +40,14 @@ CONTAINS
LOGICAL :: l_write,l_hia
REAL :: vrTmp(atoms%jmtd)
INTEGER :: i
l_write=fmpi%irank==0
IF(PRESENT(l_writeArg)) THEN
l_write = l_writeArg
ELSE
l_write = .TRUE.
ENDIF
l_write=l_write .and. fmpi%irank==0
!$ l_write = l_write .and. omp_get_num_threads()==1
......
......@@ -3,7 +3,7 @@ greensf/greensfSpinDiag.f90
greensf/greensfSpinOffDiag.f90
greensf/greensfBZint.f90
greensf/greensfSym.f90
greensf/greensfCalcImagPart.f90
greensf/greensfCalcImagPart.F90
greensf/greensfCalcRealPart.F90
greensf/greensfPostProcess.F90
greensf/kkintgr.f90
......
......@@ -9,19 +9,22 @@ MODULE m_excSplitting
CONTAINS
SUBROUTINE excSplitting(gfinp,input,greensfImagPart,ef)
SUBROUTINE excSplitting(gfinp,atoms,input,greensfImagPart,ef)
TYPE(t_gfinp), INTENT(IN) :: gfinp
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_input), INTENT(IN) :: input
TYPE(t_greensfImagPart), INTENT(IN) :: greensfImagPart
REAL, INTENT(IN) :: ef
INTEGER :: i_gf,i_elem,indUnique,ispin,kkcut,m,l,lp,atomType,atomTypep,ie
LOGICAL :: l_sphavg
REAL :: excSplit,del
REAL, ALLOCATABLE :: eMesh(:)
REAL, ALLOCATABLE :: eMesh(:), imag(:)
REAL, ALLOCATABLE :: intCOM(:,:), intNorm(:,:)
CHARACTER(LEN=20) :: attributes(4)
IF(.NOT.gfinp%checkSphavg()) RETURN
CALL gfinp%eMesh(ef,del=del,eMesh=eMesh)
......@@ -41,13 +44,14 @@ MODULE m_excSplitting
lp = gfinp%elem(i_gf)%lp
atomType = gfinp%elem(i_gf)%atomType
atomTypep = gfinp%elem(i_gf)%atomTypep
l_sphavg = gfinp%elem(i_gf)%l_sphavg
!Only onsite exchange splitting
IF(l /= lp) CYCLE
IF(atomType /= atomTypep) CYCLE
IF(.NOT.gfinp%l_sphavg) CYCLE
IF(.NOT.l_sphavg) CYCLE
i_elem = gfinp%uniqueElements(ind=i_gf,indUnique=indUnique)
i_elem = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg,indUnique=indUnique)
IF(i_gf /= indUnique) CYCLE
!-------------------------------------------------
......@@ -61,15 +65,13 @@ MODULE m_excSplitting
DO ispin = 1, input%jspins
intCOM = 0.0
intNorm = 0.0
kkcut = greensfImagPart%kkintgr_cutoff(i_gf,ispin,2)
DO m = -l, l
DO ie = 1, kkcut
intCOM(ie,ispin) = intCOM(ie,ispin) + eMesh(ie)*greensfImagPart%sphavg(ie,m,m,i_elem,ispin)
intNorm(ie,ispin) = intNorm(ie,ispin) + greensfImagPart%sphavg(ie,m,m,i_elem,ispin)
ENDDO
imag = greensfImagPart%applyCutoff(i_elem,i_gf,m,m,ispin,l_sphavg)
intCOM(:,ispin) = intCOM(:,ispin) + eMesh*imag
intNorm(:,ispin) = intNorm(:,ispin) + imag
ENDDO
excSplit = excSplit + (-1)**(ispin) * 1.0/(trapz(intNorm(:,ispin),del,kkcut)) &
* trapz(intCOM(:,ispin),del,kkcut)
excSplit = excSplit + (-1)**(ispin) * 1.0/(trapz(intNorm(:,ispin),del,SIZE(eMesh))) &
* trapz(intCOM(:,ispin),del,SIZE(eMesh))
ENDDO
WRITE(oUnit,'(A,I4,A,I4,A,f10.4,A)') ' atom: ', atomType, ' l: ', l,&
' DeltaExc: ',excSplit * hartree_to_ev_const, ' eV'
......@@ -79,7 +81,8 @@ MODULE m_excSplitting
WRITE(attributes(2),'(i0)') l
WRITE(attributes(3),'(f12.7)') excSplit * hartree_to_ev_const
WRITE(attributes(4),'(a2)') 'eV'
CALL writeXMLElementForm('excSplit',['atomType','l ','Delta ','unit '],attributes,reshape([8,1,5,4,6,1,12,2],[4,2]))
CALL writeXMLElementForm('excSplit',['atomType','l ','Delta ','unit '],&
attributes,reshape([8,1,5,4,6,1,12,2],[4,2]))
ENDDO
WRITE(oUnit,'(/)')
......
......@@ -34,6 +34,7 @@ MODULE m_greensfBZint
INTEGER :: spin1,spin2,ispin,spin_start,spin_end
COMPLEX :: phase
REAL :: atomFactor
LOGICAL :: l_sphavg
COMPLEX, ALLOCATABLE :: im(:,:,:,:,:)
spin_start = MERGE(1,jspin,gfinp%l_mperp)
......@@ -44,9 +45,6 @@ MODULE m_greensfBZint
CALL timestart("Green's Function: Brillouin-Zone-Integration")
ALLOCATE(im(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,nBands,&
MERGE(1,4,gfinp%l_sphavg),spin_start:spin_end),source=cmplx_0)
DO i_gf = 1, gfinp%n
!Get the information about the current element
......@@ -54,20 +52,29 @@ MODULE m_greensfBZint
lp = gfinp%elem(i_gf)%lp
atomType = gfinp%elem(i_gf)%atomType
atomTypep = gfinp%elem(i_gf)%atomTypep
l_sphavg = gfinp%elem(i_gf)%l_sphavg
atomFactor = MERGE(1.0,1.0/atoms%neq(atomType),l.NE.lp)
atomFactor = MERGE(1.0,atomFactor,atomType.NE.atomTypep)
i_elem = gfinp%uniqueElements(ind=i_gf,indUnique=indUnique)
i_elem = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg,indUnique=indUnique)
IF(i_gf/=indUnique) CYCLE
ALLOCATE(im(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,nBands,&
MERGE(1,4,l_sphavg),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)
!Loop over equivalent atoms
DO natom = natom_start , natom_end
!Only perform the second atom loop if we calculate intersite elements
!natomp_start = MERGE(natom,SUM(atoms%neq(:atomTypep-1)) + 1,atomType==atomTypep)
!natomp_end = MERGE(natom,SUM(atoms%neq(:atomTypep)) ,atomType==atomTypep)
!Deactivate this loop (notice natomp_end) (only calculate intersite between representative atoms)
natomp_start = MERGE(natom,SUM(atoms%neq(:atomTypep-1)) + 1,atomType==atomTypep)
natomp_end = MERGE(natom,SUM(atoms%neq(:atomTypep)) ,atomType==atomTypep)
natomp_end = MERGE(natom,SUM(atoms%neq(:atomTypep-1)) + 1,atomType==atomTypep)
DO natomp = natomp_start, natomp_end
......@@ -75,7 +82,7 @@ MODULE m_greensfBZint
IF(natom.NE.natomp) THEN
phase = exp(ImagUnit*dot_product(kpts%bk(:,ikpt),gfinp%elem(i_gf)%atomDiff))
ELSE
phase = 1.0
phase = cmplx_1
ENDIF
!l-offdiagonal phase
......@@ -94,21 +101,20 @@ MODULE m_greensfBZint
!Spin diagonal elements
IF(spin1==spin2) THEN
CALL greensfSpinDiag(nBands,l,lp,natom,natomp,atomType,atomTypep,spin1,&
gfinp%l_sphavg,atoms,usdus,eigVecCoeffs,im(:,:,:,:,ispin))
l_sphavg,atoms,usdus,eigVecCoeffs,im(:,:,:,:,ispin))
ELSE
!Spin offdiagonal elements
CALL greensfSpinOffDiag(nBands,l,lp,natom,natomp,atomType,atomTypep,spin1,spin2,&
gfinp%l_sphavg,atoms,denCoeffsOffdiag,eigVecCoeffs,im(:,:,:,:,ispin))
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)
ENDDO
CALL greensfSym(ikpt_i,i_elem,natom,l,natom.EQ.natomp.AND.l.EQ.lp,gfinp%l_sphavg,&
spin_start,spin_end,sym,atomFactor,phase,im,greensfBZintCoeffs)
ENDDO !natomp
ENDDO !natom
DEALLOCATE(im)
ENDDO !i_gf
CALL timestop("Green's Function: Brillouin-Zone-Integration")
......
This diff is collapsed.
......@@ -6,6 +6,7 @@ MODULE m_greensfPostProcess
USE m_greensfCalcRealPart
USE m_greensf_io
USE m_occmtx
USE m_greensfTorgue
USE m_excSplitting
USE m_crystalfield
USE m_genMTBasis
......@@ -16,7 +17,7 @@ MODULE m_greensfPostProcess
CONTAINS
SUBROUTINE greensfPostProcess(greensFunction,greensfImagPart,atoms,gfinp,input,sym,noco,mpi,&
nococonv,vTot,enpara,hub1inp,hub1data,results)
nococonv,vTot,enpara,hub1inp,sphhar,hub1data,results)
!contains all the modules for calculating properties from the greens function
......@@ -28,6 +29,7 @@ MODULE m_greensfPostProcess
TYPE(t_mpi), INTENT(IN) :: mpi
TYPE(t_nococonv), INTENT(IN) :: nococonv
TYPE(t_hub1inp), INTENT(IN) :: hub1inp
TYPE(t_sphhar), INTENT(IN) :: sphhar
TYPE(t_results), INTENT(IN) :: results
TYPE(t_potden), INTENT(IN) :: vTot
TYPE(t_enpara), INTENT(IN) :: enpara
......@@ -37,7 +39,9 @@ MODULE m_greensfPostProcess
INTEGER i_gf,nType,l,lp,atomType,atomTypep,i_elem,indUnique,jspin,ierr
COMPLEX mmpmat(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,gfinp%n,3)
LOGICAL l_sphavg,l_check
REAL :: torgue(3)
REAL, ALLOCATABLE :: u(:,:,:,:,:,:),udot(:,:,:,:,:,:)
REAL, ALLOCATABLE :: uun21(:,:),udn21(:,:),dun21(:,:),ddn21(:,:)
REAL, ALLOCATABLE :: f(:,:,:),g(:,:,:), flo(:,:,:)
......@@ -66,9 +70,9 @@ MODULE m_greensfPostProcess
CALL crystal_field(atoms,gfinp,hub1inp,input,nococonv,greensfImagPart,vTot,results%ef,hub1data)
ENDIF
CALL excSplitting(gfinp,input,greensfImagPart,results%ef)
CALL excSplitting(gfinp,atoms,input,greensfImagPart,results%ef)
IF(.NOT.gfinp%l_sphavg) THEN
IF(gfinp%checkRadial()) THEN
CALL timestart("Green's Function: Radial Functions")
ALLOCATE (f(atoms%jmtd,2,0:atoms%lmaxd),source=0.0)
ALLOCATE (g(atoms%jmtd,2,0:atoms%lmaxd),source=0.0)
......@@ -85,8 +89,10 @@ MODULE m_greensfPostProcess
lp = gfinp%elem(i_gf)%lp
atomType = gfinp%elem(i_gf)%atomType
atomTypep = gfinp%elem(i_gf)%atomTypep
l_sphavg = gfinp%elem(i_gf)%l_sphavg
IF(l_sphavg) CYCLE
i_elem = gfinp%uniqueElements(ind=i_gf,indUnique=indUnique)
i_elem = gfinp%uniqueElements(atoms,ind=i_gf,l_sphavg=l_sphavg,indUnique=indUnique)
IF(i_gf/=indUnique) THEN
u(:,:,:,:,:,i_gf) = u(:,:,:,:,:,indUnique)
......@@ -94,7 +100,7 @@ MODULE m_greensfPostProcess
ELSE
DO jspin = 1, input%jspins
CALL genMTBasis(atoms,enpara,vTot,mpi,atomType,jspin,usdus,&
f,g,flo,hub1inp%l_dftspinpol)
f,g,flo,hub1inp%l_dftspinpol,l_writeArg=.FALSE.)
u(:,:,1,1,jspin,i_gf) = f(:,:,l)
u(:,:,2,1,jspin,i_gf) = f(:,:,lp)
......@@ -103,7 +109,7 @@ MODULE m_greensfPostProcess
udot(:,:,2,1,jspin,i_gf) = g(:,:,lp)
CALL genMTBasis(atoms,enpara,vTot,mpi,atomTypep,jspin,usdus,&
f,g,flo,hub1inp%l_dftspinpol)
f,g,flo,hub1inp%l_dftspinpol,l_writeArg=.FALSE.)
u(:,:,1,2,jspin,i_gf) = f(:,:,l)
u(:,:,2,2,jspin,i_gf) = f(:,:,lp)
......@@ -133,22 +139,37 @@ MODULE m_greensfPostProcess
lp = greensFunction(i_gf)%elem%lp
atomType = greensFunction(i_gf)%elem%atomType
atomTypep = greensFunction(i_gf)%elem%atomTypep
l_sphavg = gfinp%elem(i_gf)%l_sphavg
IF(l.NE.lp) CYCLE
IF(atomType.NE.atomTypep) CYCLE
IF(gfinp%l_sphavg) THEN
CALL occmtx(greensFunction(i_gf),gfinp,input,mmpmat(:,:,i_gf,:),l_write=.TRUE.,check=.TRUE.)
l_check = gfinp%countLOs(atoms,i_gf)==0 !If there are SCLOs present the occupations can get bigger than 1
IF(l_sphavg) THEN
CALL occmtx(greensFunction(i_gf),gfinp,input,mmpmat(:,:,i_gf,:),l_write=.TRUE.,check=l_check)
ELSE IF(.NOT.gfinp%l_mperp) THEN
CALL occmtx(greensFunction(i_gf),gfinp,input,mmpmat(:,:,i_gf,:),&
ddn=usdus%ddn(l,atomType,:),l_write=.TRUE.,check=.TRUE.)
ddn=usdus%ddn(l,atomType,:),l_write=.TRUE.,check=l_check)
ELSE
CALL occmtx(greensFunction(i_gf),gfinp,input,mmpmat(:,:,i_gf,:),&
ddn=usdus%ddn(l,atomType,:),uun21=uun21(l,atomType),&
udn21=udn21(l,atomType),dun21=dun21(l,atomType),&
ddn21=ddn21(l,atomType),l_write=.TRUE.,check=.TRUE.)
ddn21=ddn21(l,atomType),l_write=.TRUE.,check=l_check)
ENDIF
ENDDO
CALL timestop("Green's Function: Occupation")
IF(ANY(gfinp%numTorgueElems>0)) THEN
CALL timestart("Green's Function: Torgue")
CALL openXMLElementNoAttributes('torgueCalculation')
WRITE(oUnit,'(/,A)') 'Torgue Calculation:'
WRITE(oUnit,'(/,A)') '------------------------'
DO atomType = 1, atoms%nType
IF(gfinp%numTorgueElems(atomType)==0) CYCLE
CALL greensfTorgue(greensFunction(gfinp%torgueElem(atomType,:gfinp%numTorgueElems(atomType))),vTot,&
sphhar,atoms,sym,noco,nococonv,input,enpara,hub1inp,mpi,atomType,torgue)
ENDDO
CALL closeXMLElement('torgueCalculation')
CALL timestop("Green's Function: Torgue")
ENDIF
#ifdef CPP_HDF
CALL timestart("Green's Function: IO/Write")
CALL openGreensFFile(greensf_fileID, input, gfinp, atoms)
......
......@@ -4,10 +4,12 @@ MODULE m_greensfSym
USE m_types
USE m_symMMPmat
IMPLICIT NONE
CONTAINS
SUBROUTINE greensfSym(ikpt_i,i_elem,natom,l,l_onsite,l_sphavg,spin_start,spin_end,&
sym,atomFactor,phase,im,greensfBZintCoeffs)
SUBROUTINE greensfSym(ikpt_i,i_elem,natom,l,l_onsite,l_sphavg,ispin,&
sym,atomFactor,addPhase,im,greensfBZintCoeffs)
INTEGER, INTENT(IN) :: ikpt_i
INTEGER, INTENT(IN) :: i_elem
......@@ -15,45 +17,44 @@ MODULE m_greensfSym
INTEGER, INTENT(IN) :: l
LOGICAL, INTENT(IN) :: l_onsite
LOGICAL, INTENT(IN) :: l_sphavg
INTEGER, INTENT(IN) :: spin_start
INTEGER, INTENT(IN) :: spin_end
INTEGER, INTENT(IN) :: ispin
TYPE(t_sym), INTENT(IN) :: sym
REAL, INTENT(IN) :: atomFactor
COMPLEX, INTENT(IN) :: phase
COMPLEX, INTENT(IN) :: im(-lmaxU_const:,-lmaxU_const:,:,:,:)
COMPLEX, INTENT(IN) :: addPhase
COMPLEX, INTENT(IN) :: im(-lmaxU_const:,-lmaxU_const:,:,:)
TYPE(t_greensfBZintCoeffs), INTENT(INOUT) :: greensfBZintCoeffs
INTEGER imat,iBand
COMPLEX, ALLOCATABLE :: imSym(:,:,:)
COMPLEX, ALLOCATABLE :: imSym(:,:)
!$OMP parallel default(none) &
!$OMP shared(ikpt_i,i_elem,natom,l,l_onsite,l_sphavg)&
!$OMP shared(spin_start,spin_end,sym,atomFactor,phase,im,greensfBZintCoeffs)&
!$OMP shared(ispin,sym,atomFactor,addPhase,im,greensfBZintCoeffs)&
!$OMP private(imat,iBand,imSym)
ALLOCATE(imSym(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,spin_start:spin_end),source=cmplx_0)
ALLOCATE(imSym(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const),source=cmplx_0)
!$OMP do collapse(2)
DO imat = 1, SIZE(im,4)
DO iBand = 1, SIZE(im,3)
IF(l_onsite) THEN !These rotations are only available for the onsite elements
imSym = symMMPmat(im(:,:,iBand,imat,:),sym,natom,l)