Commit 969a7a51 authored by Daniel Wortmann's avatar Daniel Wortmann

Merge branch 'develop' into 'xc-pot-refactor'

# Conflicts:
#   types/types.F90
parents ffaa9200 917abc43
This diff is collapsed.
......@@ -24,7 +24,7 @@ MODULE m_eparas
!
CONTAINS
SUBROUTINE eparas(jsp,atoms,noccbd, mpi,ikpt,ne,we,eig,skip_t,l_evp,eigVecCoeffs,&
usdus,regCharges,mcd,l_mcd)
usdus,regCharges,dos,mcd,l_mcd)
USE m_types
IMPLICIT NONE
TYPE(t_usdus), INTENT(IN) :: usdus
......@@ -32,6 +32,7 @@ CONTAINS
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_eigVecCoeffs), INTENT(IN) :: eigVecCoeffs
TYPE(t_regionCharges), INTENT(INOUT) :: regCharges
TYPE(t_dos), INTENT(INOUT) :: dos
TYPE(t_mcd), INTENT(INOUT) :: mcd
! ..
! .. Scalar Arguments ..
......@@ -64,9 +65,9 @@ CONTAINS
ENDIF
regCharges%ener(:,:,jsp) = 0.0
regCharges%sqal(:,:,jsp) = 0.0
regCharges%qal(:,:,:,jsp) = 0.0
regCharges%enerlo(:,:,jsp) = 0.0
regCharges%sqlo(:,:,jsp) = 0.0
dos%qal(:,:,:,ikpt,jsp) = 0.0
END IF
!
!---> l-decomposed density for each occupied state
......@@ -109,7 +110,7 @@ CONTAINS
ENDDO
ENDIF ! end MCD
ENDDO
regCharges%qal(l,n,i,jsp) = (suma+sumb*usdus%ddn(l,n,jsp))/atoms%neq(n)
dos%qal(l,n,i,ikpt,jsp) = (suma+sumb*usdus%ddn(l,n,jsp))/atoms%neq(n)
ENDDO
nt1 = nt1 + atoms%neq(n)
ENDDO
......@@ -122,8 +123,8 @@ CONTAINS
DO l = 0,3
DO n = 1,atoms%ntype
DO i = (skip_t+1),noccbd
regCharges%ener(l,n,jsp) = regCharges%ener(l,n,jsp) + regCharges%qal(l,n,i,jsp)*we(i)*eig(i)
regCharges%sqal(l,n,jsp) = regCharges%sqal(l,n,jsp) + regCharges%qal(l,n,i,jsp)*we(i)
regCharges%ener(l,n,jsp) = regCharges%ener(l,n,jsp) + dos%qal(l,n,i,ikpt,jsp)*we(i)*eig(i)
regCharges%sqal(l,n,jsp) = regCharges%sqal(l,n,jsp) + dos%qal(l,n,i,ikpt,jsp)*we(i)
ENDDO
ENDDO
ENDDO
......@@ -176,7 +177,7 @@ CONTAINS
! llo > 3 used for unoccupied states only
IF( l .GT. 3 ) CYCLE
DO i = 1,ne
regCharges%qal(l,ntyp,i,jsp)= regCharges%qal(l,ntyp,i,jsp) + ( 1.0/atoms%neq(ntyp) )* (&
dos%qal(l,ntyp,i,ikpt,jsp)= dos%qal(l,ntyp,i,ikpt,jsp) + ( 1.0/atoms%neq(ntyp) )* (&
qaclo(i,lo,ntyp)*usdus%uulon(lo,ntyp,jsp)+qbclo(i,lo,ntyp)*usdus%dulon(lo,ntyp,jsp) )
END DO
DO lop = 1,atoms%nlo(ntyp)
......@@ -184,7 +185,7 @@ CONTAINS
DO i = 1,ne
regCharges%enerlo(lo,ntyp,jsp) = regCharges%enerlo(lo,ntyp,jsp) +qlo(i,lop,lo,ntyp)*we(i)*eig(i)
regCharges%sqlo(lo,ntyp,jsp) = regCharges%sqlo(lo,ntyp,jsp) + qlo(i,lop,lo,ntyp)*we(i)
regCharges%qal(l,ntyp,i,jsp)= regCharges%qal(l,ntyp,i,jsp) + ( 1.0/atoms%neq(ntyp) ) *&
dos%qal(l,ntyp,i,ikpt,jsp)= dos%qal(l,ntyp,i,ikpt,jsp) + ( 1.0/atoms%neq(ntyp) ) *&
qlo(i,lop,lo,ntyp)*usdus%uloulopn(lop,lo,ntyp,jsp)
ENDDO
ENDIF
......
......@@ -9,25 +9,13 @@ MODULE m_int21
!-----------------------------------------------------------
CONTAINS
SUBROUTINE int_21(f,g,atoms,ityp,l,denCoeffsOffdiag)
USE m_types
IMPLICIT NONE
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_denCoeffsOffdiag), INTENT(INOUT) :: denCoeffsOffdiag
INTEGER, INTENT (IN) :: l,ityp
REAL, INTENT (IN) :: f(:,:,0:,:)!(atoms%jmtd,2,0:atoms%lmaxd,dimension%jspd)
REAL, INTENT (IN) :: g(:,:,0:,:)!(atoms%jmtd,2,0:atoms%lmaxd,dimension%jspd)
CALL int_21_arrays(f,g,atoms,ityp,l,denCoeffsOffdiag%uu21n,denCoeffsOffdiag%ud21n,&
denCoeffsOffdiag%du21n,denCoeffsOffdiag%dd21n)
END SUBROUTINE int_21
SUBROUTINE int_21_arrays(f,g,atoms,ityp,l,uu21n,ud21n,du21n,dd21n)
SUBROUTINE int_21(f,g,atoms,ityp,l,uu21n,ud21n,du21n,dd21n)
USE m_intgr, ONLY : intgr3
USE m_types
USE m_types_setup
IMPLICIT NONE
TYPE(t_atoms), INTENT(IN) :: atoms
INTEGER, INTENT (IN) :: l,ityp
......@@ -55,5 +43,6 @@ CONTAINS
+ g(:atoms%jri(ityp),2,l,2)*g(:atoms%jri(ityp),2,l,1)
CALL intgr3(uu_tmp,atoms%rmsh(:,ityp),atoms%dx(ityp),atoms%jri(ityp),dd21n(l,ityp))
END SUBROUTINE int_21_arrays
END SUBROUTINE int_21
END MODULE m_int21
......@@ -10,14 +10,17 @@ MODULE m_int21lo
!
!-----------------------------------------------------------
CONTAINS
SUBROUTINE int_21lo(f,g,atoms,n, flo,ilo,denCoeffsOffdiag)
SUBROUTINE int_21lo(f,g,atoms,n,flo,ilo,uulo21n,ulou21n,dulo21n,ulod21n,uloulop21n)
USE m_intgr, ONLY : intgr3
USE m_types
USE m_types
USE m_types_setup
IMPLICIT NONE
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_denCoeffsOffdiag), INTENT(INOUT) :: denCoeffsOffdiag
TYPE(t_atoms), INTENT(IN) :: atoms
REAL, INTENT(INOUT) :: uulo21n(atoms%nlod,atoms%ntype)
REAL, INTENT(INOUT) :: ulou21n(atoms%nlod,atoms%ntype)
REAL, INTENT(INOUT) :: dulo21n(atoms%nlod,atoms%ntype)
REAL, INTENT(INOUT) :: ulod21n(atoms%nlod,atoms%ntype)
REAL, INTENT(INOUT) :: uloulop21n(atoms%nlod,atoms%nlod,atoms%ntype)
! ..
! .. Scalar Arguments ..
INTEGER, INTENT (IN) :: ilo,n
......@@ -38,22 +41,22 @@ CONTAINS
DO iri = 1, atoms%jri(n)
uu_tmp(iri) = f(iri,1,l,2)*flo(iri,1,ilo,1)+ f(iri,2,l,2)*flo(iri,2,ilo,1)
ENDDO
CALL intgr3(uu_tmp,atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),denCoeffsOffdiag%uulo21n(ilo,n))
CALL intgr3(uu_tmp,atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),uulo21n(ilo,n))
DO iri = 1, atoms%jri(n)
uu_tmp(iri) = f(iri,1,l,1)*flo(iri,1,ilo,2)+ f(iri,2,l,1)*flo(iri,2,ilo,2)
ENDDO
CALL intgr3(uu_tmp,atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),denCoeffsOffdiag%ulou21n(ilo,n))
CALL intgr3(uu_tmp,atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),ulou21n(ilo,n))
!
! --> norm of product of du and ulo:
!
DO iri = 1, atoms%jri(n)
uu_tmp(iri) = g(iri,1,l,2)*flo(iri,1,ilo,1) + g(iri,2,l,2)*flo(iri,2,ilo,1)
ENDDO
CALL intgr3(uu_tmp,atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),denCoeffsOffdiag%dulo21n(ilo,n))
CALL intgr3(uu_tmp,atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),dulo21n(ilo,n))
DO iri = 1, atoms%jri(n)
uu_tmp(iri) = g(iri,1,l,1)*flo(iri,1,ilo,2) + g(iri,2,l,1)*flo(iri,2,ilo,2)
ENDDO
CALL intgr3(uu_tmp,atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),denCoeffsOffdiag%ulod21n(ilo,n))
CALL intgr3(uu_tmp,atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),ulod21n(ilo,n))
!
! --> norm of product of ulo and ulo':
!
......@@ -63,10 +66,10 @@ CONTAINS
DO iri = 1, atoms%jri(n)
uu_tmp(iri) = flo(iri,1,ilo,2)*flo(iri,1,ilop,1) + flo(iri,2,ilo,2)*flo(iri,2,ilop,1)
ENDDO
CALL intgr3(uu_tmp,atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),denCoeffsOffdiag%uloulop21n(ilo,ilop,n))
CALL intgr3(uu_tmp,atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),uloulop21n(ilo,ilop,n))
ELSE
denCoeffsOffdiag%uloulop21n(ilo,ilop,n) = 0.0
uloulop21n(ilo,ilop,n) = 0.0
ENDIF
ENDDO
......
......@@ -176,8 +176,8 @@ CONTAINS
! pgfft : contains the phases of the g-vectors of sph.
! isn : isn = +1, FFT transform for g-space to r-space
! isn = -1, vice versa
!
CALL timestart("pwden")
ALLOCATE(cwk(stars%ng3),ecwk(stars%ng3))
......@@ -715,5 +715,7 @@ CONTAINS
IF (input%l_f) DEALLOCATE ( kpsir,kpsii,ekin)
ENDIF
CALL timestop("pwden")
END SUBROUTINE pwden
END MODULE m_pwden
......@@ -5,7 +5,7 @@ MODULE m_qal21
!***********************************************************************
!
CONTAINS
SUBROUTINE qal_21(dimension,atoms,input,noccbd,noco,eigVecCoeffs,denCoeffsOffdiag,regCharges)
SUBROUTINE qal_21(dimension,atoms,input,noccbd,noco,eigVecCoeffs,denCoeffsOffdiag,ikpt,dos)
USE m_rotdenmat
USE m_types
......@@ -16,10 +16,10 @@ CONTAINS
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_eigVecCoeffs), INTENT(IN) :: eigVecCoeffs
TYPE(t_denCoeffsOffdiag), INTENT(IN) :: denCoeffsOffdiag
TYPE(t_regionCharges), INTENT(INOUT) :: regCharges
TYPE(t_dos), INTENT(INOUT) :: dos
! .. Scalar Arguments ..
INTEGER, INTENT (IN) :: noccbd
INTEGER, INTENT (IN) :: noccbd,ikpt
! .. Local Scalars ..
INTEGER i,l,lo,lop ,natom,nn,ntyp
......@@ -149,11 +149,10 @@ CONTAINS
state : DO i = 1, noccbd
lls : DO l = 0,3
CALL rot_den_mat(noco%alph(n),noco%beta(n),&
regCharges%qal(l,n,i,1),regCharges%qal(l,n,i,2),qal21(l,n,i))
dos%qal(l,n,i,ikpt,1),dos%qal(l,n,i,ikpt,2),qal21(l,n,i))
IF (.FALSE.) THEN
IF (n==1) WRITE(*,'(3i3,4f10.5)') l,n,i,qal21(l,n,i),&
regCharges%qal(l,n,i,:)
q_loc(1,1) = regCharges%qal(l,n,i,1); q_loc(2,2) = regCharges%qal(l,n,i,2)
IF (n==1) WRITE(*,'(3i3,4f10.5)') l,n,i,qal21(l,n,i),dos%qal(l,n,i,ikpt,:)
q_loc(1,1) = dos%qal(l,n,i,ikpt,1); q_loc(2,2) = dos%qal(l,n,i,ikpt,2)
q_loc(1,2) = qal21(l,n,i); q_loc(2,1) = CONJG(q_loc(1,2))
q_hlp = MATMUL( TRANSPOSE( CONJG(chi) ) ,q_loc)
q_loc = MATMUL(q_hlp,chi)
......
......@@ -5,14 +5,9 @@ MODULE m_vacden
! vacuum charge density. speed up by r. wu 1992
! *************************************************************
CONTAINS
SUBROUTINE vacden(&
vacuum,DIMENSION,stars,oneD,&
kpts,input,cell,atoms,noco,banddos,&
gvac1,gvac2,&
we,ikpt,jspin,vz,&
ne,lapw,&
evac,eig,den,qvac,qvlay,&
stcoeff,zMat)
SUBROUTINE vacden(vacuum,DIMENSION,stars,oneD,kpts,input,sym,cell,atoms,noco,banddos,&
gVacMap,we,ikpt,jspin,vz,ne,lapw,evac,eig,den,qvac,qvlay,&
stcoeff,zMat)
!***********************************************************************
! ****** change vacden(....,q) for vacuum density of states shz Jan.96
......@@ -61,10 +56,12 @@ CONTAINS
TYPE(t_vacuum),INTENT(IN) :: vacuum
TYPE(t_noco),INTENT(IN) :: noco
TYPE(t_stars),INTENT(IN) :: stars
TYPE(t_sym),INTENT(IN) :: sym
TYPE(t_cell),INTENT(IN) :: cell
TYPE(t_kpts),INTENT(IN) :: kpts
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_zMat),INTENT(IN) :: zMat
TYPE(t_gVacMap),INTENT(IN) :: gVacMap
TYPE(t_potden),INTENT(INOUT) :: den
! .. Scalar Arguments ..
INTEGER, INTENT (IN) :: jspin
......@@ -73,15 +70,14 @@ CONTAINS
INTEGER,PARAMETER :: n2max=13
REAL,PARAMETER :: emax=2.0/hartree_to_ev_const
! .. Array Arguments ..
REAL, INTENT(IN) :: evac(2,DIMENSION%jspd)
REAL, INTENT(OUT) :: qvlay(DIMENSION%neigd,vacuum%layerd,2,kpts%nkpt,DIMENSION%jspd)
REAL, INTENT(INOUT) :: qvac(DIMENSION%neigd,2,kpts%nkpt,DIMENSION%jspd)
REAL, INTENT(IN) :: we(DIMENSION%neigd)
REAL :: vz(vacuum%nmzd,2) ! Note this breaks the INTENT(IN) from cdnval. It may be read from a file in this subroutine.
REAL, INTENT(IN) :: evac(2,DIMENSION%jspd)
REAL, INTENT(OUT) :: qvlay(DIMENSION%neigd,vacuum%layerd,2,kpts%nkpt,DIMENSION%jspd)
REAL, INTENT(INOUT) :: qvac(DIMENSION%neigd,2,kpts%nkpt,DIMENSION%jspd)
REAL, INTENT(IN) :: we(DIMENSION%neigd)
REAL :: vz(vacuum%nmzd,2) ! Note this breaks the INTENT(IN) from cdnval. It may be read from a file in this subroutine.
! STM-Arguments
REAL, INTENT (IN) :: eig(DIMENSION%neigd)
INTEGER, INTENT (IN) :: gvac1(DIMENSION%nv2d),gvac2(DIMENSION%nv2d)
COMPLEX, INTENT (OUT):: stcoeff(vacuum%nstars,DIMENSION%neigd,vacuum%layerd,2)
REAL, INTENT (IN) :: eig(DIMENSION%neigd)
COMPLEX, INTENT (INOUT) :: stcoeff(vacuum%nstars,DIMENSION%neigd,vacuum%layerd,2,kpts%nkpt,input%jspins)
!
! local STM variables
INTEGER nv2(DIMENSION%jspd)
......@@ -143,7 +139,9 @@ CONTAINS
! \vec{a}_1,2 are the 2D lattice vectors
!
! **************************************************************************************************
!
CALL timestart("vacden")
ALLOCATE ( ac(DIMENSION%nv2d,DIMENSION%neigd,DIMENSION%jspd),bc(DIMENSION%nv2d,DIMENSION%neigd,DIMENSION%jspd),dt(DIMENSION%nv2d),&
& dte(DIMENSION%nv2d),du(vacuum%nmzd),ddu(vacuum%nmzd,DIMENSION%nv2d),due(vacuum%nmzd),&
& ddue(vacuum%nmzd,DIMENSION%nv2d),t(DIMENSION%nv2d),te(DIMENSION%nv2d),&
......@@ -172,8 +170,6 @@ CONTAINS
! ------------------
! WRITE (16,'(a,i2)') 'nstars=',nstars
stcoeff(:,:,:,:) = CMPLX(0.0,0.0)
! -----> set up mapping arrays
IF (noco%l_ss) THEN
jsp_start = 1
......@@ -244,17 +240,14 @@ CONTAINS
DO j=1, n2max
mapg2k(j)=j
DO i=1, nv2(jspin)
IF ( kvac1(i,jspin).EQ.gvac1(j) .AND.&
& kvac2(i,jspin).EQ.gvac2(j) ) &
& mapg2k(j)=i
IF (kvac1(i,jspin).EQ.gVacMap%gvac1d(j).AND.kvac2(i,jspin).EQ.gVacMap%gvac2d(j)) mapg2k(j)=i
END DO
END DO
END IF
!
!-dw
IF (noco%l_noco) THEN
OPEN (25,FILE='potmat',FORM='unformatted',&
& STATUS='old')
OPEN (25,FILE='potmat',FORM='unformatted',STATUS='old')
!---> skip the four components of the interstitial potential matrix
DO ipot = 1,3
READ (25)
......@@ -1203,7 +1196,7 @@ CONTAINS
!=============================================================
!
! calculate 1. to nstars. starcoefficient for each k and energy eigenvalue
! to stcoeff(ne,layer,ivac) if starcoeff=T (the star coefficient values are written to vacdos)
! to stcoeff(ne,layer,ivac,ikpt) if starcoeff=T (the star coefficient values are written to vacdos)
!
IF (vacuum%starcoeff .AND. banddos%vacdos) THEN
DO n=1,ne
......@@ -1234,9 +1227,9 @@ CONTAINS
uej = ue(vacuum%izlay(jj,1),l1,jspin)
t1 = aa*ui*uj + bb*uei*uej +ba*ui*uej + ab*uei*uj
IF (ind2.GE.2.AND.ind2.LE.vacuum%nstars) &
stcoeff(ind2-1,n,jj,ivac) = stcoeff(ind2-1,n,jj,ivac)+ t1*phs/stars%nstr2(ind2)
stcoeff(ind2-1,n,jj,ivac,ikpt,jspin) = stcoeff(ind2-1,n,jj,ivac,ikpt,jspin)+ t1*phs/stars%nstr2(ind2)
IF (ind2p.GE.2.AND.ind2p.LE.vacuum%nstars) &
stcoeff(ind2p-1,n,jj,ivac) = stcoeff(ind2p-1,n,jj,ivac) +CONJG(t1)*phs/stars%nstr2(ind2p)
stcoeff(ind2p-1,n,jj,ivac,ikpt,jspin) = stcoeff(ind2p-1,n,jj,ivac,ikpt,jspin) +CONJG(t1)*phs/stars%nstr2(ind2p)
END DO
END IF
ENDDO
......@@ -1252,5 +1245,16 @@ CONTAINS
DEALLOCATE (t_1,te_1,tei_1,u_1,ue_1)
END IF ! oneD%odi%d1
IF(vacuum%nvac.EQ.1) THEN
den%vacz(:,2,:) = den%vacz(:,1,:)
IF (sym%invs) THEN
den%vacxy(:,:,2,:) = CONJG(den%vacxy(:,:,1,:))
ELSE
den%vacxy(:,:,2,:) = den%vacxy(:,:,1,:)
END IF
END IF
CALL timestop("vacden")
END SUBROUTINE vacden
END MODULE m_vacden
......@@ -21,4 +21,6 @@ cdn_mt/rhonmt.f90
cdn_mt/rhonmt21.f90
cdn_mt/rhonmtlo.f90
cdn_mt/rhosphnlo.f90
cdn_mt/calcDenCoeffs.f90
)
......@@ -15,6 +15,7 @@ CONTAINS
USE m_abclocdn
USE m_ylm
USE m_types
USE m_juDFT
IMPLICIT NONE
TYPE(t_input),INTENT(IN) :: input
TYPE(t_usdus),INTENT(IN) :: usdus
......@@ -55,6 +56,8 @@ CONTAINS
REAL, ALLOCATABLE :: work_r(:)
COMPLEX, ALLOCATABLE :: work_c(:)
CALL timestart("abcof")
IF (zmat%l_real) THEN
IF (noco%l_soc.AND.sym%invs) CALL judft_error("BUG in abcof, SOC&INVS but real?")
IF (noco%l_noco) CALL judft_error("BUG in abcof, l_noco but real?")
......@@ -362,5 +365,7 @@ CONTAINS
ENDDO
ENDIF
CALL timestop("abcof")
END SUBROUTINE abcof
END MODULE m_abcof
!--------------------------------------------------------------------------------
! Copyright (c) 2018 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------
MODULE m_calcDenCoeffs
CONTAINS
SUBROUTINE calcDenCoeffs(atoms,sphhar,sym,we,noccbd,eigVecCoeffs,ispin,denCoeffs)
USE m_juDFT
USE m_types
USE m_rhomt
USE m_rhonmt
USE m_rhomtlo
USE m_rhonmtlo
IMPLICIT NONE
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_sphhar), INTENT(IN) :: sphhar
TYPE(t_sym), INTENT(IN) :: sym
TYPE(t_eigVecCoeffs), INTENT(IN) :: eigVecCoeffs
TYPE(t_denCoeffs), INTENT(INOUT) :: denCoeffs
REAL, INTENT(IN) :: we(noccbd)
INTEGER, INTENT(IN) :: noccbd
INTEGER, INTENT(IN) :: ispin
!---> set up coefficients for the spherical and
CALL timestart("cdnval: rhomt")
CALL rhomt(atoms,we,noccbd,eigVecCoeffs,denCoeffs,ispin)
CALL timestop("cdnval: rhomt")
!---> non-spherical m.t. density
CALL timestart("cdnval: rhonmt")
CALL rhonmt(atoms,sphhar,we,noccbd,sym,eigVecCoeffs,denCoeffs,ispin)
CALL timestop("cdnval: rhonmt")
!---> set up coefficients of the local orbitals and the
!---> flapw - lo cross terms for the spherical and
!---> non-spherical mt density
CALL timestart("cdnval: rho(n)mtlo")
CALL rhomtlo(atoms,noccbd,we,eigVecCoeffs,denCoeffs,ispin)
CALL rhonmtlo(atoms,sphhar,noccbd,we,eigVecCoeffs,denCoeffs,ispin)
CALL timestop("cdnval: rho(n)mtlo")
END SUBROUTINE calcDenCoeffs
END MODULE m_calcDenCoeffs
......@@ -8,7 +8,7 @@ MODULE m_cdncore
CONTAINS
SUBROUTINE cdncore(results,mpi,dimension,oneD,sliceplot,input,vacuum,noco,sym,&
SUBROUTINE cdncore(results,mpi,dimension,oneD,input,vacuum,noco,sym,&
stars,cell,sphhar,atoms,vTot,outDen,moments)
USE m_constants
......@@ -30,7 +30,6 @@ SUBROUTINE cdncore(results,mpi,dimension,oneD,sliceplot,input,vacuum,noco,sym,&
TYPE(t_mpi),INTENT(IN) :: mpi
TYPE(t_dimension),INTENT(IN) :: dimension
TYPE(t_oneD),INTENT(IN) :: oneD
TYPE(t_sliceplot),INTENT(IN) :: sliceplot
TYPE(t_input),INTENT(IN) :: input
TYPE(t_vacuum),INTENT(IN) :: vacuum
TYPE(t_noco),INTENT(IN) :: noco
......@@ -81,60 +80,58 @@ SUBROUTINE cdncore(results,mpi,dimension,oneD,sliceplot,input,vacuum,noco,sym,&
END IF
END IF
IF (.NOT.sliceplot%slice) THEN
!add in core density
IF (mpi%irank.EQ.0) THEN
IF (input%kcrel.EQ.0) THEN
DO jspin = 1,input%jspins
CALL cored(input,jspin,atoms,outDen%mt,dimension,sphhar,vTot%mt(:,0,:,jspin), qint,rh,tec,seig)
rhTemp(:,:,jspin) = rh(:,:,jspin)
results%seigc = results%seigc + seig
END DO
ELSE
CALL coredr(input,atoms,seig, outDen%mt,dimension,sphhar,vTot%mt(:,0,:,:),qint,rh)
!add in core density
IF (mpi%irank.EQ.0) THEN
IF (input%kcrel.EQ.0) THEN
DO jspin = 1,input%jspins
CALL cored(input,jspin,atoms,outDen%mt,dimension,sphhar,vTot%mt(:,0,:,jspin), qint,rh,tec,seig)
rhTemp(:,:,jspin) = rh(:,:,jspin)
results%seigc = results%seigc + seig
END IF
END DO
ELSE
CALL coredr(input,atoms,seig, outDen%mt,dimension,sphhar,vTot%mt(:,0,:,:),qint,rh)
results%seigc = results%seigc + seig
END IF
DO jspin = 1,input%jspins
IF (mpi%irank.EQ.0) THEN
DO n = 1,atoms%ntype
moments%stdn(n,jspin) = outDen%mt(1,0,n,jspin) / (sfp_const*atoms%rmsh(1,n)*atoms%rmsh(1,n))
END IF
DO jspin = 1,input%jspins
IF (mpi%irank.EQ.0) THEN
DO n = 1,atoms%ntype
moments%stdn(n,jspin) = outDen%mt(1,0,n,jspin) / (sfp_const*atoms%rmsh(1,n)*atoms%rmsh(1,n))
END DO
END IF
IF ((noco%l_noco).AND.(mpi%irank.EQ.0)) THEN
IF (jspin.EQ.2) THEN
!pk non-collinear (start)
!add the coretail-charge to the constant interstitial
!charge (star 0), taking into account the direction of
!magnetisation of this atom
DO iType = 1,atoms%ntype
rhoint = (qint(iType,1) + qint(iType,2)) /cell%volint/input%jspins/2.0
momint = (qint(iType,1) - qint(iType,2)) /cell%volint/input%jspins/2.0
!rho_11
outDen%pw(1,1) = outDen%pw(1,1) + rhoint + momint*cos(noco%beta(iType))
!rho_22
outDen%pw(1,2) = outDen%pw(1,2) + rhoint - momint*cos(noco%beta(iType))
!real part rho_21
outDen%pw(1,3) = outDen%pw(1,3) + cmplx(0.5*momint *cos(noco%alph(iType))*sin(noco%beta(iType)),0.0)
!imaginary part rho_21
outDen%pw(1,3) = outDen%pw(1,3) + cmplx(0.0,-0.5*momint *sin(noco%alph(iType))*sin(noco%beta(iType)))
END DO
!pk non-collinear (end)
END IF
IF ((noco%l_noco).AND.(mpi%irank.EQ.0)) THEN
IF (jspin.EQ.2) THEN
!pk non-collinear (start)
!add the coretail-charge to the constant interstitial
!charge (star 0), taking into account the direction of
!magnetisation of this atom
DO iType = 1,atoms%ntype
rhoint = (qint(iType,1) + qint(iType,2)) /cell%volint/input%jspins/2.0
momint = (qint(iType,1) - qint(iType,2)) /cell%volint/input%jspins/2.0
!rho_11
outDen%pw(1,1) = outDen%pw(1,1) + rhoint + momint*cos(noco%beta(iType))
!rho_22
outDen%pw(1,2) = outDen%pw(1,2) + rhoint - momint*cos(noco%beta(iType))
!real part rho_21
outDen%pw(1,3) = outDen%pw(1,3) + cmplx(0.5*momint *cos(noco%alph(iType))*sin(noco%beta(iType)),0.0)
!imaginary part rho_21
outDen%pw(1,3) = outDen%pw(1,3) + cmplx(0.0,-0.5*momint *sin(noco%alph(iType))*sin(noco%beta(iType)))
END DO
!pk non-collinear (end)
END IF
ELSE
IF (input%ctail) THEN
!+gu hope this works as well
CALL cdnovlp(mpi,sphhar,stars,atoms,sym,dimension,vacuum,&
cell,input,oneD,l_st,jspin,rh(:,:,jspin),&
outDen%pw,outDen%vacxy,outDen%mt,outDen%vacz)
ELSE IF (mpi%irank.EQ.0) THEN
DO iType = 1,atoms%ntype
outDen%pw(1,jspin) = outDen%pw(1,jspin) + qint(iType,jspin)/input%jspins/cell%volint
END DO
END IF
ELSE
IF (input%ctail) THEN
!+gu hope this works as well
CALL cdnovlp(mpi,sphhar,stars,atoms,sym,dimension,vacuum,&
cell,input,oneD,l_st,jspin,rh(:,:,jspin),&
outDen%pw,outDen%vacxy,outDen%mt,outDen%vacz)
ELSE IF (mpi%irank.EQ.0) THEN
DO iType = 1,atoms%ntype
outDen%pw(1,jspin) = outDen%pw(1,jspin) + qint(iType,jspin)/input%jspins/cell%volint
END DO
END IF
END DO
END IF
END IF
END DO
IF (input%kcrel.EQ.0) THEN
IF (mpi%irank.EQ.0) THEN
......
......@@ -10,8 +10,8 @@ MODULE m_cdnmt
! Philipp Kurz 2000-02-03
!***********************************************************************
CONTAINS
SUBROUTINE cdnmt(jspd,atoms,sphhar,llpd, noco,l_fmpl,jsp_start,jsp_end, epar,&
ello,vr,denCoeffs,usdus,orb,denCoeffsOffdiag,moments,rho)
SUBROUTINE cdnmt(jspd,atoms,sphhar,noco,jsp_start,jsp_end,enpara,&
vr,denCoeffs,usdus,orb,denCoeffsOffdiag,moments,rho)
use m_constants,only: sfp_const
USE m_rhosphnlo
USE m_radfun
......@@ -23,17 +23,14 @@ CONTAINS
TYPE(t_noco), INTENT(IN) :: noco
TYPE(t_sphhar), INTENT(IN) :: sphhar
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_enpara), INTENT(IN) :: enpara
TYPE(t_moments), INTENT(INOUT) :: moments
! .. Scalar Arguments ..
INTEGER, INTENT (IN) :: llpd
INTEGER, INTENT (IN) :: jsp_start,jsp_end,jspd
LOGICAL, INTENT (IN) :: l_fmpl
! ..
! .. Array Arguments ..
REAL, INTENT (IN) :: epar(0:atoms%lmaxd,atoms%ntype,jspd)
REAL, INTENT (IN) :: vr(atoms%jmtd,atoms%ntype,jspd)
REAL, INTENT (IN) :: ello(atoms%nlod,atoms%ntype,jspd)
REAL, INTENT (INOUT) :: rho(:,0:,:,:)!(toms%jmtd,0:sphhar%nlhd,atoms%ntype,jspd)
TYPE (t_orb), INTENT(IN) :: orb
TYPE (t_denCoeffs), INTENT(IN) :: denCoeffs
......@@ -58,7 +55,7 @@ CONTAINS
CALL timestart("cdnmt")
IF (noco%l_mperp) THEN
IF (l_fmpl) THEN
IF (denCoeffsOffdiag%l_fmpl) THEN
ALLOCATE ( rho21(atoms%jmtd,0:sphhar%nlhd,atoms%ntype) )
rho21(:,:,:) = cmplx(0.0,0.0)
ENDIF
......@@ -67,8 +64,8 @@ CONTAINS
!$OMP PARALLEL DEFAULT(none) &
!$OMP SHARED(usdus,rho,moments,rho21,qmtl) &
!$OMP SHARED(atoms,jsp_start,jsp_end,epar,vr,denCoeffs,sphhar,ello)&
!$OMP SHARED(orb,noco,l_fmpl,denCoeffsOffdiag,jspd)&
!$OMP SHARED(atoms,jsp_start,jsp_end,enpara,vr,denCoeffs,sphhar)&
!$OMP SHARED(orb,noco,denCoeffsOffdiag,jspd)&
!$OMP PRIVATE(itype,na,ispin,l,f,g,nodeu,noded,wronk,i,j,s,qmtllo,qmtt,nd,lh,lp,llp,cs)
IF (noco%l_mperp) THEN
ALLOCATE ( f(atoms%jmtd,2,0:atoms%lmaxd,jspd),g(atoms%jmtd,2,0:atoms%lmaxd,jspd) )
......@@ -88,7 +85,7 @@ CONTAINS
!---> spherical component
DO ispin = jsp_start,jsp_end
DO l = 0,atoms%lmax(itype)
CALL radfun(l,itype,ispin,epar(l,itype,ispin),vr(1,itype,ispin),atoms,&
CALL radfun(l,itype,ispin,enpara%el0(l,itype,ispin),vr(1,itype,ispin),atoms,&
f(1,1,l,ispin),g(1,1,l,ispin),usdus, nodeu,noded,wronk)
DO j = 1,atoms%jri(itype)
s = denCoeffs%uu(l,itype,ispin)*( f(j,1,l,ispin)*f(j,1,l,ispin)+f(j,2,l,ispin)*f(j,2,l,ispin) )&
......@@ -107,7 +104,7 @@ CONTAINS
CALL rhosphnlo(itype,atoms,sphhar,&
usdus%uloulopn(1,1,itype,ispin),usdus%dulon(1,itype,ispin),&
usdus%uulon(1,itype,ispin),ello(1,itype,ispin),&
usdus%uulon(1,itype,ispin),enpara%ello0(1,itype,ispin),&
vr(1,itype,ispin),denCoeffs%aclo(1,itype,ispin),denCoeffs%bclo(1,itype,ispin),&
denCoeffs%cclo(1,1,itype,ispin),denCoeffs%acnmt(0,1,1,itype,ispin),&
denCoeffs%bcnmt(0,1,1,itype,ispin),denCoeffs%ccnmt(1,1,1,itype,ispin),&
......@@ -178,7 +175,7 @@ CONTAINS
ENDDO
ENDDO
IF (l_fmpl) THEN
IF (denCoeffsOffdiag%l_fmpl) THEN
!---> the following part can be used to calculate the full magnet.
!---> density without the atomic sphere approximation for the
!---> magnet. density, e.g. for plotting.
......@@ -212,7 +209,7 @@ CONTAINS
ENDDO
ENDDO
ENDIF ! l_fmpl
ENDIF ! denCoeffsOffdiag%l_fmpl
ENDIF ! noco%l_mperp
ENDDO ! end of loop over atom types
......@@ -246,7 +243,7 @@ CONTAINS
!---> for testing: to plot the offdiag. part of the density matrix it
!---> is written to the file rhomt21. This file can read in pldngen.
IF (l_fmpl) THEN
IF (denCoeffsOffdiag%l_fmpl) THEN
OPEN (26,file='rhomt21',form='unformatted',status='unknown')
WRITE (26) rho21
CLOSE</