Commit 8ee7d48f authored by Gregor Michalicek's avatar Gregor Michalicek

Introduce t_regionCharges type to cdn/cdnval.F90

parent e89fd760
......@@ -3,7 +3,7 @@ MODULE m_cdnval
CONTAINS
SUBROUTINE cdnval(eig_id, mpi,kpts,jspin,sliceplot,noco, input,banddos,cell,atoms,enpara,stars,&
vacuum,dimension,sphhar,sym,obsolete,vTot,oneD,coreSpecInput,den,results,&
qvac,qvlay,qa21, chmom,clmom)
qvlay,qa21, chmom,clmom)
!
! ***********************************************************
! this subroutin is a modified version of cdnval.F.
......@@ -108,52 +108,38 @@ CONTAINS
TYPE(t_potden),INTENT(INOUT) :: den
! .. Scalar Arguments ..
INTEGER, INTENT (IN) :: eig_id,jspin
INTEGER, INTENT(IN) :: eig_id,jspin
! .. Array Arguments ..
COMPLEX, INTENT(INOUT) :: qa21(atoms%ntype)
REAL, INTENT (OUT) :: chmom(atoms%ntype,dimension%jspd),clmom(3,atoms%ntype,dimension%jspd)
REAL, INTENT (INOUT) :: qvac(dimension%neigd,2,kpts%nkpt,dimension%jspd)
REAL, INTENT (INOUT) :: qvlay(dimension%neigd,vacuum%layerd,2,kpts%nkpt,dimension%jspd)
COMPLEX, INTENT(INOUT) :: qa21(atoms%ntype)
REAL, INTENT(OUT) :: chmom(atoms%ntype,dimension%jspd)
REAL, INTENT(OUT) :: clmom(3,atoms%ntype,dimension%jspd)
REAL, INTENT(INOUT) :: qvlay(dimension%neigd,vacuum%layerd,2,kpts%nkpt,dimension%jspd)
#ifdef CPP_MPI
INCLUDE 'mpif.h'
logical :: mpi_flag, mpi_status
LOGICAL :: mpi_flag, mpi_status
#endif
! .. Local Scalars ..
TYPE(t_lapw):: lapw
INTEGER :: llpd
INTEGER i,ie,iv,ivac,j,k,l,n,ilo,isp,&
nbands,noccbd,nslibd,na,&
ikpt,jsp_start,jsp_end,ispin
INTEGER skip_t,skip_tt
INTEGER n_size,i_rec,n_rank ,ncored,n_start,n_end,noccbd_l,nbasfcn
LOGICAL l_fmpl,l_evp,l_orbcomprot,l_real, l_write
INTEGER :: llpd,ikpt,jsp_start,jsp_end,ispin
INTEGER :: i,ie,iv,ivac,j,k,l,n,ilo,isp,nbands,noccbd,nslibd,na
INTEGER :: skip_t,skip_tt, nkpt_extended
INTEGER :: n_size,i_rec,n_rank ,ncored,n_start,n_end,noccbd_l,nbasfcn
LOGICAL :: l_fmpl,l_evp,l_orbcomprot,l_real, l_write
! ...Local Arrays ..
INTEGER n_bands(0:dimension%neigd)
REAL eig(dimension%neigd)
REAL vz0(2)
INTEGER :: n_bands(0:dimension%neigd)
REAL :: eig(dimension%neigd)
REAL :: vz0(2)
!orbcomp
REAL, ALLOCATABLE :: orbcomp(:,:,:),qmtp(:,:)
REAL, ALLOCATABLE :: qis(:,:,:)
!-dw
REAL, ALLOCATABLE :: orbcomp(:,:,:),qmtp(:,:) ! orbcomp
INTEGER, ALLOCATABLE :: gvac1d(:),gvac2d(:)
INTEGER, ALLOCATABLE :: jsym(:),ksym(:)
REAL, ALLOCATABLE :: we(:)
! radial functions
REAL, ALLOCATABLE :: f(:,:,:,:),g(:,:,:,:),flo(:,:,:,:)
REAL, ALLOCATABLE :: sqlo(:,:,:)
REAL, ALLOCATABLE :: qal(:,:,:,:),sqal(:,:,:),ener(:,:,:)
REAL, ALLOCATABLE :: svac(:,:),pvac(:,:)
REAL, ALLOCATABLE :: enerlo(:,:,:),qmat(:,:,:,:)
REAL, ALLOCATABLE :: f(:,:,:,:),g(:,:,:,:),flo(:,:,:,:) ! radial functions
REAL, ALLOCATABLE :: qmat(:,:,:,:)
COMPLEX, ALLOCATABLE :: qstars(:,:,:,:)
TYPE (t_lapw) :: lapw
TYPE (t_orb) :: orb
TYPE (t_denCoeffs) :: denCoeffs
TYPE (t_denCoeffsOffdiag) :: denCoeffsOffdiag
......@@ -161,10 +147,9 @@ CONTAINS
TYPE (t_slab) :: slab
TYPE (t_eigVecCoeffs) :: eigVecCoeffs
TYPE (t_mcd) :: mcd
TYPE (t_usdus) :: usdus
TYPE (t_zMat) :: zMat
INTEGER :: nkpt_extended
TYPE (t_regionCharges) :: regCharges
TYPE (t_usdus) :: usdus
TYPE (t_zMat) :: zMat
l_real = sym%invs.AND.(.NOT.noco%l_soc).AND.(.NOT.noco%l_noco)
......@@ -190,15 +175,9 @@ CONTAINS
!---> at a time, otherwise for both spins:
ALLOCATE ( f(atoms%jmtd,2,0:atoms%lmaxd,jsp_start:jsp_end) ) ! Deallocation before mpi_col_den
ALLOCATE ( g(atoms%jmtd,2,0:atoms%lmaxd,jsp_start:jsp_end) )
ALLOCATE ( flo(atoms%jmtd,2,atoms%nlod,dimension%jspd) )
ALLOCATE ( jsym(dimension%neigd),ksym(dimension%neigd) )
ALLOCATE ( gvac1d(dimension%nv2d),gvac2d(dimension%nv2d) )
ALLOCATE ( qal(0:3,atoms%ntype,dimension%neigd,jsp_start:jsp_end) )
ALLOCATE ( sqal(0:3,atoms%ntype,jsp_start:jsp_end) )
ALLOCATE ( ener(0:3,atoms%ntype,jsp_start:jsp_end) )
ALLOCATE ( sqlo(atoms%nlod,atoms%ntype,jsp_start:jsp_end) )
ALLOCATE ( enerlo(atoms%nlod,atoms%ntype,jsp_start:jsp_end) )
ALLOCATE ( svac(2,jsp_start:jsp_end) )
ALLOCATE ( pvac(2,jsp_start:jsp_end) )
ALLOCATE ( qstars(vacuum%nstars,dimension%neigd,vacuum%layerd,2) )
! --> Initializations
......@@ -207,22 +186,17 @@ CONTAINS
CALL denCoeffs%init(atoms,sphhar,jsp_start,jsp_end)
CALL denCoeffsOffdiag%init(atoms,noco,sphhar,l_fmpl)
CALL force%init1(input,atoms)
IF ((l_fmpl).AND.(.not.noco%l_mperp)) CALL juDFT_error("for fmpl set noco%l_mperp = T!" ,calledby ="cdnval")
svac(:,:) = 0.0 ; pvac(:,:) = 0.0
sqal(:,:,:) = 0.0 ; ener(:,:,:) = 0.0
CALL regCharges%init(atoms,dimension,kpts,jsp_start,jsp_end)
CALL orb%init(atoms,noco,jsp_start,jsp_end)
CALL mcd%init1(banddos,dimension,input,atoms)
IF ((l_fmpl).AND.(.not.noco%l_mperp)) CALL juDFT_error("for fmpl set noco%l_mperp = T!" ,calledby ="cdnval")
! calculation of core spectra (EELS) initializations -start-
CALL corespec_init(input,atoms,coreSpecInput)
IF(l_cs.AND.(mpi%isize.NE.1)) CALL juDFT_error('EELS + MPI not implemented', calledby = 'cdnval')
IF(l_cs.AND.jspin.EQ.1) CALL corespec_gaunt()
! calculation of core spectra (EELS) initializations -end-
IF (mpi%irank==0) THEN
WRITE (6,FMT=8000) jspin
......@@ -231,25 +205,20 @@ CONTAINS
END IF
8000 FORMAT (/,/,10x,'valence density: spin=',i2)
CALL cdn_read0(eig_id,mpi%irank,mpi%isize,jspin,dimension%jspd,&
noco%l_noco,n_bands,n_size)
CALL cdn_read0(eig_id,mpi%irank,mpi%isize,jspin,dimension%jspd,noco%l_noco,n_bands,n_size)
#ifdef CPP_MPI
! Sinchronizes the RMA operations
CALL MPI_BARRIER(mpi%mpi_comm,ie)
#endif
ALLOCATE ( qis(dimension%neigd,kpts%nkpt,dimension%jspd))
skip_tt = dot_product(enpara%skiplo(:atoms%ntype,jspin),atoms%neq(:atoms%ntype))
IF (noco%l_soc.OR.noco%l_noco) skip_tt = 2 * skip_tt
!-lo
!---> set up l-dependent m.t. wavefunctions
na = 1
ncored = 0
l_write = input%cdinf.AND.mpi%irank==0
ALLOCATE ( flo(atoms%jmtd,2,atoms%nlod,dimension%jspd) )
DO n = 1,atoms%ntype
DO ispin = jsp_start, jsp_end
......@@ -276,7 +245,7 @@ CONTAINS
na = na + atoms%neq(n)
END DO
DEALLOCATE (flo)
DEALLOCATE (f,g,flo)
IF (input%film) vz0(:) = vTot%vacz(vacuum%nmz,:,jspin)
......@@ -287,7 +256,6 @@ CONTAINS
ALLOCATE ( orbcomp(dimension%neigd,23,atoms%nat) )
ALLOCATE ( qmtp(dimension%neigd,atoms%nat) )
IF (.NOT.input%film) qvac(:,:,:,jspin) = 0.0
ELSE
ALLOCATE(orbcomp(1,1,1),qmtp(1,1))
END IF
......@@ -295,11 +263,6 @@ CONTAINS
!--> loop over k-points: each can be a separate task
IF (kpts%nkpt < mpi%isize) THEN
l_evp = .true.
ener(:,:,:) = 0.0
sqal(:,:,:) = 0.0
qal(:,:,:,:) = 0.0
enerlo(:,:,:) = 0.0
sqlo(:,:,:) = 0.0
ELSE
l_evp = .false.
END IF
......@@ -475,7 +438,7 @@ CONTAINS
IF (.NOT.((jspin.EQ.2) .AND. noco%l_noco)) THEN
CALL timestart("cdnval: pwden")
CALL pwden(stars,kpts,banddos,oneD,input,mpi,noco,cell,atoms,sym,ikpt,&
jspin,lapw,noccbd,we,eig,den,qis,results,force%f_b8,zMat)
jspin,lapw,noccbd,we,eig,den,regCharges%qis,results,force%f_b8,zMat)
CALL timestop("cdnval: pwden")
END IF
!+new
......@@ -494,15 +457,15 @@ CONTAINS
CALL timestart("cdnval: vacden")
CALL vacden(vacuum,dimension,stars,oneD, kpts,input, cell,atoms,noco,banddos,&
gvac1d,gvac2d, we,ikpt,jspin,vTot%vacz(:,:,jspin),vz0, noccbd,lapw, enpara%evac0,eig,&
den,qvac,qvlay, qstars,zMat)
den,regCharges%qvac,qvlay, qstars,zMat)
CALL timestop("cdnval: vacden")
END IF
!---> perform Brillouin zone integration and summation over the
!---> bands in order to determine the vacuum energy parameters.
DO ispin = jsp_start,jsp_end
DO ivac = 1,vacuum%nvac
pvac(ivac,ispin)=pvac(ivac,ispin)+dot_product(eig(:noccbd)*qvac(:noccbd,ivac,ikpt,ispin),we(:noccbd))
svac(ivac,ispin)=svac(ivac,ispin)+dot_product(qvac(:noccbd,ivac,ikpt,ispin),we(:noccbd))
regCharges%pvac(ivac,ispin)=regCharges%pvac(ivac,ispin)+dot_product(eig(:noccbd)*regCharges%qvac(:noccbd,ivac,ikpt,ispin),we(:noccbd))
regCharges%svac(ivac,ispin)=regCharges%svac(ivac,ispin)+dot_product(regCharges%qvac(:noccbd,ivac,ikpt,ispin),we(:noccbd))
END DO
END DO
END IF
......@@ -536,13 +499,10 @@ CONTAINS
!---> atom and angular momentum
IF (.not.sliceplot%slice) THEN
CALL eparas(ispin,atoms,noccbd,mpi,ikpt,noccbd,we,eig,&
skip_t,l_evp,eigVecCoeffs,usdus,mcd,&
banddos%l_mcd,enerlo(1,1,ispin),sqlo(1,1,ispin),&
ener(0,1,ispin),sqal(0,1,ispin),&
qal(0:,:,:,ispin))
skip_t,l_evp,eigVecCoeffs,usdus,regCharges,mcd,banddos%l_mcd)
IF (noco%l_mperp.AND.(ispin == jsp_end)) THEN
CALL qal_21(atoms,input,noccbd,we,noco,eigVecCoeffs,denCoeffsOffdiag,qal,qmat)
CALL qal_21(atoms,input,noccbd,we,noco,eigVecCoeffs,denCoeffsOffdiag,regCharges,qmat)
END IF
END IF
!
......@@ -600,9 +560,7 @@ CONTAINS
IF (noco%l_mperp) THEN
CALL rhomt21(atoms,we,noccbd,eigVecCoeffs,denCoeffsOffdiag)
IF (l_fmpl) THEN
CALL rhonmt21(atoms,sphhar,we,noccbd,sym,eigVecCoeffs,denCoeffsOffdiag)
END IF
IF (l_fmpl) CALL rhonmt21(atoms,sphhar,we,noccbd,sym,eigVecCoeffs,denCoeffsOffdiag)
END IF
199 CONTINUE
......@@ -622,7 +580,7 @@ CONTAINS
END IF
!--dw now write k-point data to tmp_dos
CALL write_dos(eig_id,ikpt,jspin,qal(:,:,:,jspin),qvac(:,:,ikpt,jspin),qis(:,ikpt,jspin),&
CALL write_dos(eig_id,ikpt,jspin,regCharges%qal(:,:,:,jspin),regCharges%qvac(:,:,ikpt,jspin),regCharges%qis(:,ikpt,jspin),&
qvlay(:,:,:,ikpt,jspin),qstars,ksym,jsym,mcd%mcd,slab%qintsl,&
slab%qmtsl(:,:),qmtp(:,:),orbcomp)
......@@ -639,17 +597,15 @@ CONTAINS
END IF
END IF ! --> end "IF ((mod(i_rec-1,mpi%isize).EQ.mpi%irank).OR.l_evp) THEN"
END DO !---> end of k-point loop
DEALLOCATE (we,f,g)
DEALLOCATE (we)
!+t3e
#ifdef CPP_MPI
CALL timestart("cdnval: mpi_col_den")
DO ispin = jsp_start,jsp_end
CALL mpi_col_den(mpi,sphhar,atoms,oneD,stars,vacuum,&
input,noco,l_fmpl,ispin,llpd, den%vacxy(1,1,1,ispin),&
den%vacz(1,1,ispin),den%pw(1,ispin), ener(0,1,ispin),sqal(0,1,ispin),&
results,svac(1,ispin),pvac(1,ispin),denCoeffs,&
sqlo(1,1,ispin),enerlo(1,1,ispin),orb,&
denCoeffsOffdiag,den,den%mmpMat(:,:,:,jspin))
den%vacz(1,1,ispin),den%pw(1,ispin),regCharges,&
results,denCoeffs,orb,denCoeffsOffdiag,den,den%mmpMat(:,:,:,jspin))
END DO
CALL timestop("cdnval: mpi_col_den")
#endif
......@@ -675,16 +631,15 @@ CONTAINS
CALL cdnmt(dimension%jspd,atoms,sphhar,llpd,&
noco,l_fmpl,jsp_start,jsp_end,&
enpara%el0,enpara%ello0,vTot%mt(:,0,:,:),denCoeffs,&
usdus,orb,denCoeffsOffdiag,&
chmom,clmom,qa21,den%mt)
usdus,orb,denCoeffsOffdiag,chmom,clmom,qa21,den%mt)
DO ispin = jsp_start,jsp_end
IF (.NOT.sliceplot%slice) THEN
DO n=1,atoms%ntype
enpara%el1(0:3,n,ispin)=ener(0:3,n,ispin)/sqal(0:3,n,ispin)
IF (atoms%nlo(n)>0) enpara%ello1(:atoms%nlo(n),n,ispin)=enerlo(:atoms%nlo(n),n,ispin)/sqlo(:atoms%nlo(n),n,ispin)
enpara%el1(0:3,n,ispin)=regCharges%ener(0:3,n,ispin)/regCharges%sqal(0:3,n,ispin)
IF (atoms%nlo(n)>0) enpara%ello1(:atoms%nlo(n),n,ispin)=regCharges%enerlo(:atoms%nlo(n),n,ispin)/regCharges%sqlo(:atoms%nlo(n),n,ispin)
END DO
IF (input%film) enpara%evac1(:vacuum%nvac,ispin)=pvac(:vacuum%nvac,ispin)/svac(:vacuum%nvac,ispin)
IF (input%film) enpara%evac1(:vacuum%nvac,ispin)=regCharges%pvac(:vacuum%nvac,ispin)/regCharges%svac(:vacuum%nvac,ispin)
END IF
!---> check continuity of charge density
......@@ -717,8 +672,6 @@ CONTAINS
END IF
END IF ! end of (mpi%irank==0)
!+t3e
!Note: no deallocation anymore, we rely on Fortran08 :-)
IF ((jsp_end.EQ.input%jspins)) THEN
IF ((banddos%dos.OR.banddos%vacdos).AND.(banddos%ndir/=-2)) CALL juDFT_end("DOS OK",mpi%irank)
......
......@@ -24,14 +24,15 @@ MODULE m_eparas
!
CONTAINS
SUBROUTINE eparas(jsp,atoms,noccbd, mpi,ikpt,ne,we,eig,skip_t,l_evp,eigVecCoeffs,&
usdus,mcd,l_mcd,enerlo,sqlo,ener,sqal,qal)
usdus,regCharges,mcd,l_mcd)
USE m_types
IMPLICIT NONE
TYPE(t_usdus),INTENT(IN) :: usdus
TYPE(t_mpi),INTENT(IN) :: mpi
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_eigVecCoeffs),INTENT(IN) :: eigVecCoeffs
TYPE(t_mcd),INTENT(INOUT) :: mcd
TYPE(t_usdus), INTENT(IN) :: usdus
TYPE(t_mpi), INTENT(IN) :: mpi
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_eigVecCoeffs), INTENT(IN) :: eigVecCoeffs
TYPE(t_regionCharges), INTENT(INOUT) :: regCharges
TYPE(t_mcd), INTENT(INOUT) :: mcd
! ..
! .. Scalar Arguments ..
INTEGER, INTENT (IN) :: noccbd,jsp
......@@ -41,9 +42,6 @@ CONTAINS
! .. Array Arguments ..
REAL, INTENT (IN) :: eig(:)!(dimension%neigd),
REAL, INTENT (IN) :: we(noccbd)
REAL, INTENT (INOUT) :: enerlo(atoms%nlod,atoms%ntype),sqlo(atoms%nlod,atoms%ntype)
REAL, INTENT (INOUT) :: ener(0:3,atoms%ntype),sqal(0:3,atoms%ntype)
REAL, INTENT (INOUT) :: qal(0:,:,:)!(0:3,atoms%ntype,dimension%neigd)
! ..
! .. Local Scalars ..
......@@ -64,11 +62,11 @@ CONTAINS
IF (l_mcd) THEN
mcd%mcd(:,:,:) = 0.0
ENDIF
ener(:,:) = 0.0
sqal(:,:) = 0.0
qal(:,:,:) = 0.0
enerlo(:,:) = 0.0
sqlo(:,:) = 0.0
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
END IF
!
!---> l-decomposed density for each occupied state
......@@ -111,7 +109,7 @@ CONTAINS
ENDDO
ENDIF ! end MCD
ENDDO
qal(l,n,i) = (suma+sumb*usdus%ddn(l,n,jsp))/atoms%neq(n)
regCharges%qal(l,n,i,jsp) = (suma+sumb*usdus%ddn(l,n,jsp))/atoms%neq(n)
ENDDO
nt1 = nt1 + atoms%neq(n)
ENDDO
......@@ -124,8 +122,8 @@ CONTAINS
DO l = 0,3
DO n = 1,atoms%ntype
DO i = (skip_t+1),noccbd
ener(l,n) = ener(l,n) + qal(l,n,i)*we(i)*eig(i)
sqal(l,n) = sqal(l,n) + qal(l,n,i)*we(i)
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)
ENDDO
ENDDO
ENDDO
......@@ -178,15 +176,15 @@ CONTAINS
! llo > 3 used for unoccupied states only
IF( l .GT. 3 ) CYCLE
DO i = 1,ne
qal(l,ntyp,i)= qal(l,ntyp,i) + ( 1.0/atoms%neq(ntyp) )* (&
regCharges%qal(l,ntyp,i,jsp)= regCharges%qal(l,ntyp,i,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)
IF (atoms%llo(lop,ntyp).EQ.l) THEN
DO i = 1,ne
enerlo(lo,ntyp) = enerlo(lo,ntyp) +qlo(i,lop,lo,ntyp)*we(i)*eig(i)
sqlo(lo,ntyp) = sqlo(lo,ntyp) + qlo(i,lop,lo,ntyp)*we(i)
qal(l,ntyp,i)= qal(l,ntyp,i) + ( 1.0/atoms%neq(ntyp) ) *&
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) ) *&
qlo(i,lop,lo,ntyp)*usdus%uloulopn(lop,lo,ntyp,jsp)
ENDDO
ENDIF
......
......@@ -5,23 +5,24 @@ MODULE m_qal21
!***********************************************************************
!
CONTAINS
SUBROUTINE qal_21(atoms,input,noccbd,we,noco,eigVecCoeffs,denCoeffsOffdiag,qal,qmat)
SUBROUTINE qal_21(atoms,input,noccbd,we,noco,eigVecCoeffs,denCoeffsOffdiag,regCharges,qmat)
USE m_rotdenmat
USE m_types
IMPLICIT NONE
TYPE(t_input),INTENT(IN) :: input
TYPE(t_noco),INTENT(IN) :: noco
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_eigVecCoeffs),INTENT(IN) :: eigVecCoeffs
TYPE(t_input), INTENT(IN) :: input
TYPE(t_noco), INTENT(IN) :: noco
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_eigVecCoeffs), INTENT(IN) :: eigVecCoeffs
TYPE (t_denCoeffsOffdiag), INTENT(IN) :: denCoeffsOffdiag
TYPE(t_regionCharges), INTENT(INOUT) :: regCharges
! ..
! .. Scalar Arguments ..
INTEGER, INTENT (IN) :: noccbd
! ..
! .. Array Arguments ..
REAL, INTENT (INout) :: we(noccbd),qal(0:,:,:,:)!(0:3,atoms%ntype,DIMENSION%neigd,input%jspins)
REAL, INTENT (INout) :: we(noccbd)
REAL, INTENT (OUT) :: qmat(0:,:,:,:)!(0:3,atoms%ntype,DIMENSION%neigd,4)
TYPE (t_denCoeffsOffdiag), INTENT (IN) :: denCoeffsOffdiag
! ..
! .. Local Scalars ..
......@@ -153,11 +154,11 @@ CONTAINS
state : DO i = 1, noccbd
lls : DO l = 0,3
CALL rot_den_mat(noco%alph(n),noco%beta(n),&
qal(l,n,i,1),qal(l,n,i,2),qal21(l,n,i))
regCharges%qal(l,n,i,1),regCharges%qal(l,n,i,2),qal21(l,n,i))
IF (.FALSE.) THEN
IF (n==1) WRITE(*,'(3i3,4f10.5)') l,n,i,qal21(l,n,i),&
qal(l,n,i,:)
q_loc(1,1) = qal(l,n,i,1); q_loc(2,2) = qal(l,n,i,2)
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)
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)
......
......@@ -76,7 +76,7 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
!Local Arrays
REAL stdn(atoms%ntype,dimension%jspd),svdn(atoms%ntype,dimension%jspd)
REAL chmom(atoms%ntype,dimension%jspd),clmom(3,atoms%ntype,dimension%jspd)
REAL ,ALLOCATABLE :: qvac(:,:,:,:),qvlay(:,:,:,:,:)
REAL ,ALLOCATABLE :: qvlay(:,:,:,:,:)
COMPLEX,ALLOCATABLE :: qa21(:)
IF (mpi%irank.EQ.0) THEN
......@@ -84,12 +84,10 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
IF (l_enpara) OPEN (40,file ='enpara',form = 'formatted',status ='unknown')
ENDIF
ALLOCATE (qa21(atoms%ntype))
ALLOCATE (qvac(dimension%neigd,2,kpts%nkpt,dimension%jspd))
ALLOCATE (qvlay(dimension%neigd,vacuum%layerd,2,kpts%nkpt,dimension%jspd))
!initialize density arrays with zero
qa21(:) = cmplx(0.0,0.0)
qvac(:,:,:,:) = 0.0
qvlay(:,:,:,:,:) = 0.0
IF (mpi%irank.EQ.0) CALL openXMLElementNoAttributes('valenceDensity')
......@@ -105,7 +103,7 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
CALL cdnval(eig_id,&
mpi,kpts,jspin,sliceplot,noco,input,banddos,cell,atoms,enpara,stars,vacuum,dimension,&
sphhar,sym,obsolete,vTot,oneD,coreSpecInput,&
outDen,results,qvac,qvlay,qa21,chmom,clmom)
outDen,results,qvlay,qa21,chmom,clmom)
CALL timestop("cdngen: cdnval")
END DO
......@@ -151,7 +149,7 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
CALL mpi_bc_potden(mpi,stars,sphhar,atoms,input,vacuum,oneD,noco,outDen)
#endif
DEALLOCATE (qvac,qvlay,qa21)
DEALLOCATE (qvlay,qa21)
END SUBROUTINE cdngen
......
......@@ -10,11 +10,9 @@ MODULE m_mpi_col_den
!
CONTAINS
SUBROUTINE mpi_col_den(mpi,sphhar,atoms,oneD,stars,vacuum,&
input, noco,l_fmpl,jspin,llpd,rhtxy,rht,qpw,ener,&
sqal,results,svac,pvac,denCoeffs,sqlo,&
enerlo,orb,&
denCoeffsOffdiag,den,n_mmp)
!
input, noco,l_fmpl,jspin,llpd,rhtxy,rht,qpw,regCharges,&
results,denCoeffs,orb,denCoeffsOffdiag,den,n_mmp)
#include"cpp_double.h"
USE m_types
USE m_constants
......@@ -39,14 +37,12 @@ CONTAINS
! .. Array Arguments ..
COMPLEX, INTENT (INOUT) :: qpw(stars%ng3)
COMPLEX, INTENT (INOUT) :: rhtxy(vacuum%nmzxyd,oneD%odi%n2d-1,2)
REAL, INTENT (INOUT) :: rht(vacuum%nmzd,2)
REAL, INTENT (INOUT) :: ener(0:3,atoms%ntype),sqal(0:3,atoms%ntype)
REAL, INTENT (INOUT) :: svac(2),pvac(2)
REAL, INTENT (INOUT) :: sqlo(atoms%nlod,atoms%ntype),enerlo(atoms%nlod,atoms%ntype)
REAL, INTENT (INOUT) :: rht(vacuum%nmzd,2)
COMPLEX,INTENT(INOUT) :: n_mmp(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u)
TYPE (t_orb), INTENT(INOUT) :: orb
TYPE (t_denCoeffs), INTENT(INOUT) :: denCoeffs
TYPE (t_denCoeffsOffdiag), INTENT(INOUT) :: denCoeffsOffdiag
TYPE (t_regionCharges), INTENT(INOUT) :: regCharges
! ..
! .. Local Scalars ..
INTEGER :: n
......@@ -138,13 +134,13 @@ CONTAINS
!
n=4*atoms%ntype
ALLOCATE(r_b(n))
CALL MPI_REDUCE(ener,r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
CALL MPI_REDUCE(regCharges%ener(0:,:,jspin),r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
IF (mpi%irank.EQ.0) THEN
CALL CPP_BLAS_scopy(n, r_b, 1, ener, 1)
CALL CPP_BLAS_scopy(n, r_b, 1, regCharges%ener(0:,:,jspin), 1)
ENDIF
CALL MPI_REDUCE(sqal,r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
CALL MPI_REDUCE(regCharges%sqal(0:,:,jspin),r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
IF (mpi%irank.EQ.0) THEN
CALL CPP_BLAS_scopy(n, r_b, 1, sqal, 1)
CALL CPP_BLAS_scopy(n, r_b, 1, regCharges%sqal(0:,:,jspin), 1)
ENDIF
DEALLOCATE (r_b)
!
......@@ -154,13 +150,13 @@ CONTAINS
n=2
ALLOCATE(r_b(n))
CALL MPI_REDUCE(svac,r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
CALL MPI_REDUCE(regCharges%svac(:,jspin),r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
IF (mpi%irank.EQ.0) THEN
CALL CPP_BLAS_scopy(n, r_b, 1, svac, 1)
CALL CPP_BLAS_scopy(n, r_b, 1, regCharges%svac(:,jspin), 1)
ENDIF
CALL MPI_REDUCE(pvac,r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
CALL MPI_REDUCE(regCharges%pvac(:,jspin),r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
IF (mpi%irank.EQ.0) THEN
CALL CPP_BLAS_scopy(n, r_b, 1, pvac, 1)
CALL CPP_BLAS_scopy(n, r_b, 1, regCharges%pvac(:,jspin), 1)
ENDIF
DEALLOCATE (r_b)
......@@ -194,13 +190,13 @@ CONTAINS
IF (mpi%irank.EQ.0) THEN
CALL CPP_BLAS_scopy(n, r_b, 1, denCoeffs%bclo(:,:,jspin), 1)
ENDIF
CALL MPI_REDUCE(enerlo,r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
CALL MPI_REDUCE(regCharges%enerlo(:,:,jspin),r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
IF (mpi%irank.EQ.0) THEN
CALL CPP_BLAS_scopy(n, r_b, 1, enerlo, 1)
CALL CPP_BLAS_scopy(n, r_b, 1, regCharges%enerlo(:,:,jspin), 1)
ENDIF
CALL MPI_REDUCE(sqlo,r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
CALL MPI_REDUCE(regCharges%sqlo(:,:,jspin),r_b,n,CPP_MPI_REAL,MPI_SUM,0, MPI_COMM_WORLD,ierr)
IF (mpi%irank.EQ.0) THEN
CALL CPP_BLAS_scopy(n, r_b, 1, sqlo, 1)
CALL CPP_BLAS_scopy(n, r_b, 1, regCharges%sqlo(:,:,jspin), 1)
ENDIF
DEALLOCATE (r_b)
......
......@@ -154,7 +154,27 @@ PRIVATE
PROCEDURE,PASS :: init1 => mcd_init1
END TYPE t_mcd
PUBLIC t_orb, t_denCoeffs, t_denCoeffsOffdiag, t_force, t_slab, t_eigVecCoeffs, t_mcd
TYPE t_regionCharges
REAL, ALLOCATABLE :: qis(:,:,:)
REAL, ALLOCATABLE :: qal(:,:,:,:)
REAL, ALLOCATABLE :: sqal(:,:,:)
REAL, ALLOCATABLE :: ener(:,:,:)
REAL, ALLOCATABLE :: sqlo(:,:,:)
REAL, ALLOCATABLE :: enerlo(:,:,:)
REAL, ALLOCATABLE :: qvac(:,:,:,:)
REAL, ALLOCATABLE :: svac(:,:)
REAL, ALLOCATABLE :: pvac(:,:)
CONTAINS
PROCEDURE,PASS :: init => regionCharges_init
END TYPE t_regionCharges
PUBLIC t_orb, t_denCoeffs, t_denCoeffsOffdiag, t_force, t_slab, t_eigVecCoeffs, t_mcd, t_regionCharges
CONTAINS
......@@ -578,5 +598,46 @@ SUBROUTINE mcd_init1(thisMCD,banddos,dimension,input,atoms)
END SUBROUTINE mcd_init1
SUBROUTINE regionCharges_init(thisRegCharges,atoms,dimension,kpts,jsp_start,jsp_end)
USE m_types_setup
USE m_types_kpts
IMPLICIT NONE
CLASS(t_regionCharges), INTENT(INOUT) :: thisRegCharges
TYPE(t_atoms), INTENT(IN) :: atoms
TYPE(t_dimension), INTENT(IN) :: dimension
TYPE(t_kpts), INTENT(IN) :: kpts
INTEGER, INTENT(IN) :: jsp_start
INTEGER, INTENT(IN) :: jsp_end
ALLOCATE(thisRegCharges%qis(dimension%neigd,kpts%nkpt,dimension%jspd))
ALLOCATE(thisRegCharges%qal(0:3,atoms%ntype,dimension%neigd,jsp_start:jsp_end))
ALLOCATE(thisRegCharges%sqal(0:3,atoms%ntype,jsp_start:jsp_end))
ALLOCATE(thisRegCharges%ener(0:3,atoms%ntype,jsp_start:jsp_end))
ALLOCATE(thisRegCharges%sqlo(atoms%nlod,atoms%ntype,jsp_start:jsp_end))
ALLOCATE(thisRegCharges%enerlo(atoms%nlod,atoms%ntype,jsp_start:jsp_end))