Commit 5370997f authored by Daniel Wortmann's avatar Daniel Wortmann

Merge branch 'develop' of iffgit.fz-juelich.de:fleur/fleur into develop

parents fdac2fb9 ac09716d
...@@ -25,4 +25,5 @@ cdn/qpw_to_nmt.f90 ...@@ -25,4 +25,5 @@ cdn/qpw_to_nmt.f90
cdn/slab_dim.f90 cdn/slab_dim.f90
cdn/slabgeom.f90 cdn/slabgeom.f90
cdn/vacden.F90 cdn/vacden.F90
cdn/genNewNocoInp.f90
) )
This diff is collapsed.
!--------------------------------------------------------------------------------
! 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_genNewNocoInp
CONTAINS
SUBROUTINE genNewNocoInp(input,atoms,jij,noco,noco_new)
USE m_juDFT
USE m_types
USE m_constants
USE m_rwnoco
IMPLICIT NONE
TYPE(t_input),INTENT(IN) :: input
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_jij),INTENT(IN) :: jij
TYPE(t_noco),INTENT(IN) :: noco
TYPE(t_noco),INTENT(INOUT) :: noco_new
INTEGER :: iAtom, iType
REAL :: alphdiff
IF (.NOT.noco%l_mperp) THEN
CALL juDFT_error ("genNewNocoInp without noco%l_mperp" ,calledby ="genNewNocoInp")
END IF
iAtom = 1
DO iType = 1, atoms%ntype
IF (noco%l_ss) THEN
alphdiff = 2.0*pi_const*(noco%qss(1)*atoms%taual(1,iAtom) + &
noco%qss(2)*atoms%taual(2,iAtom) + &
noco%qss(3)*atoms%taual(3,iAtom) )
noco_new%alph(iType) = noco%alph(iType) - alphdiff
DO WHILE (noco_new%alph(iType) > +pi_const)
noco_new%alph(iType)= noco_new%alph(iType) - 2.0*pi_const
END DO
DO WHILE (noco_new%alph(iType) < -pi_const)
noco_new%alph(iType)= noco_new%alph(iType) + 2.0*pi_const
END DO
ELSE
noco_new%alph(iType) = noco%alph(iType)
END IF
iatom= iatom + atoms%neq(iType)
END DO
OPEN (24,file='nocoinp',form='formatted', status='old')
REWIND (24)
CALL rw_noco_write(atoms,jij,noco_new, input)
CLOSE (24)
END SUBROUTINE genNewNocoInp
END MODULE m_genNewNocoInp
...@@ -9,6 +9,7 @@ cdn_mt/abcof.F90 ...@@ -9,6 +9,7 @@ cdn_mt/abcof.F90
cdn_mt/abcof3.F90 cdn_mt/abcof3.F90
cdn_mt/abcrot2.f90 cdn_mt/abcrot2.f90
cdn_mt/cdnmt.f90 cdn_mt/cdnmt.f90
cdn_mt/cdncore.F90
cdn_mt/magMoms.f90 cdn_mt/magMoms.f90
cdn_mt/orbMagMoms.f90 cdn_mt/orbMagMoms.f90
cdn_mt/orb_comp2.f90 cdn_mt/orb_comp2.f90
......
!--------------------------------------------------------------------------------
! 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_cdncore
CONTAINS
SUBROUTINE cdncore(results,mpi,dimension,oneD,sliceplot,input,vacuum,noco,sym,&
stars,cell,sphhar,atoms,vTot,outDen,stdn,svdn)
USE m_constants
USE m_cdn_io
USE m_cdnovlp
USE m_cored
USE m_coredr
USE m_types
USE m_xmlOutput
USE m_magMoms
USE m_orbMagMoms
#ifdef CPP_MPI
USE m_mpi_bc_coreden
#endif
IMPLICIT NONE
TYPE(t_results),INTENT(INOUT) :: results
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
TYPE(t_sym),INTENT(IN) :: sym
TYPE(t_stars),INTENT(IN) :: stars
TYPE(t_cell),INTENT(IN) :: cell
TYPE(t_sphhar),INTENT(IN) :: sphhar
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_potden),INTENT(IN) :: vTot
TYPE(t_potden),INTENT(INOUT) :: outDen
REAL, INTENT(INOUT) :: stdn(atoms%ntype,dimension%jspd)
REAL, INTENT(INOUT) :: svdn(atoms%ntype,dimension%jspd)
INTEGER :: jspin, n, iType
REAL :: seig, rhoint, momint
LOGICAL, PARAMETER :: l_st=.FALSE.
REAL :: rh(dimension%msh,atoms%ntype,dimension%jspd)
REAL :: qint(atoms%ntype,dimension%jspd)
REAL :: tec(atoms%ntype,DIMENSION%jspd)
REAL :: rhTemp(dimension%msh,atoms%ntype,dimension%jspd)
results%seigc = 0.0
IF (mpi%irank.EQ.0) THEN
DO jspin = 1,input%jspins
DO n = 1,atoms%ntype
svdn(n,jspin) = outDen%mt(1,0,n,jspin) / (sfp_const*atoms%rmsh(1,n)*atoms%rmsh(1,n))
END DO
END DO
END IF
IF (input%kcrel.EQ.0) THEN
! Generate input file ecore for subsequent GW calculation
! 11.2.2004 Arno Schindlmayr
IF ((input%gw.eq.1 .or. input%gw.eq.3).AND.(mpi%irank.EQ.0)) THEN
OPEN (15,file='ecore',status='unknown', action='write',form='unformatted')
END IF
rh = 0.0
tec = 0.0
qint = 0.0
IF (input%frcor) THEN
IF (mpi%irank.EQ.0) THEN
CALL readCoreDensity(input,atoms,dimension,rh,tec,qint)
END IF
#ifdef CPP_MPI
CALL mpi_bc_coreDen(mpi,atoms,input,dimension,rh,tec,qint)
#endif
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)
results%seigc = results%seigc + seig
END IF
END IF
DO jspin = 1,input%jspins
IF (mpi%irank.EQ.0) THEN
DO n = 1,atoms%ntype
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
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 IF
END DO
END IF
IF (input%kcrel.EQ.0) THEN
IF (mpi%irank.EQ.0) THEN
CALL writeCoreDensity(input,atoms,dimension,rhTemp,tec,qint)
END IF
IF ((input%gw.eq.1 .or. input%gw.eq.3).AND.(mpi%irank.EQ.0)) CLOSE(15)
END IF
END SUBROUTINE cdncore
END MODULE m_cdncore
...@@ -20,6 +20,7 @@ global/savewigner.f ...@@ -20,6 +20,7 @@ global/savewigner.f
set(fleur_F90 ${fleur_F90} set(fleur_F90 ${fleur_F90}
global/constants.f90 global/constants.f90
global/checkdop.F90 global/checkdop.F90
global/checkdopall.f90
global/chkmt.f90 global/chkmt.f90
global/convn.f90 global/convn.f90
global/enpara.f90 global/enpara.f90
......
MODULE m_checkdop MODULE m_checkdop
CONTAINS CONTAINS
SUBROUTINE checkdop(& SUBROUTINE checkdop(&
& p,np,n,na,ivac,iflag,jsp,cdn,& & p,np,n,na,ivac,iflag,jsp,&
& DIMENSION,atoms,sphhar,stars,sym,& & DIMENSION,atoms,sphhar,stars,sym,&
& vacuum,cell,oneD,& & vacuum,cell,oneD,potden)
& fpw,fr,fxy,fz)
! ************************************************************ ! ************************************************************
! subroutines checks the continuity of coulomb * ! subroutines checks the continuity of coulomb *
! potential or valence charge density * ! potential or valence charge density *
...@@ -14,6 +13,7 @@ ...@@ -14,6 +13,7 @@
! YM: this routine doesn't really work in the vacuum in 1D case yet ! YM: this routine doesn't really work in the vacuum in 1D case yet
! ************************************************************ ! ************************************************************
USE m_juDFT
USE m_starf, ONLY : starf2,starf3 USE m_starf, ONLY : starf2,starf3
USE m_angle USE m_angle
USE m_ylm USE m_ylm
...@@ -23,32 +23,37 @@ ...@@ -23,32 +23,37 @@
IMPLICIT NONE IMPLICIT NONE
! .. ! ..
! .. Scalar Arguments .. ! .. Scalar Arguments ..
TYPE(t_dimension),INTENT(IN):: dimension TYPE(t_dimension),INTENT(IN) :: dimension
type(t_sphhar),intent(in):: sphhar type(t_sphhar),intent(in) :: sphhar
TYPE(t_stars),INTENT(IN) :: stars TYPE(t_stars),INTENT(IN) :: stars
TYPE(t_atoms),INTENT(IN) :: atoms TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_sym),INTENT(IN) :: sym TYPE(t_sym),INTENT(IN) :: sym
TYPE(t_vacuum),INTENT(IN):: vacuum TYPE(t_vacuum),INTENT(IN) :: vacuum
TYPE(t_oneD),INTENT(IN) :: oneD TYPE(t_oneD),INTENT(IN) :: oneD
TYPE(t_cell),INTENT(IN) :: cell TYPE(t_cell),INTENT(IN) :: cell
TYPE(t_potden),INTENT(IN) :: potden
INTEGER, INTENT (IN) :: iflag,ivac,n,na,np,jsp INTEGER, INTENT (IN) :: iflag,ivac,n,na,np,jsp
LOGICAL, INTENT (IN) :: cdn
!-odim !-odim
!+odim !+odim
! .. Array Arguments .. ! .. Array Arguments ..
REAL, INTENT (IN) :: p(3,DIMENSION%nspd) REAL, INTENT (IN) :: p(3,DIMENSION%nspd)
REAL, INTENT (IN) :: fr(atoms%jmtd,0:sphhar%nlhd,atoms%ntype,DIMENSION%jspd),fz(vacuum%nmzd,2,DIMENSION%jspd)
COMPLEX, INTENT (IN) :: fpw(stars%ng3,DIMENSION%jspd),fxy(vacuum%nmzxyd,oneD%odi%n2d-1,2,DIMENSION%jspd)
! .. ! ..
! .. Local Scalars .. ! .. Local Scalars ..
REAL av,dms,rms,s,ir2,help,phi REAL av,dms,rms,s,ir2,help,phi
INTEGER i,j,k,lh,mem,nd,lm,ll1,nopa ,gz,m INTEGER i,j,k,lh,mem,nd,lm,ll1,nopa ,gz,m
COMPLEX ic COMPLEX ic
LOGICAL l_cdn
! .. ! ..
! .. Local Arrays .. ! .. Local Arrays ..
COMPLEX sf2(stars%ng2),sf3(stars%ng3),ylm( (atoms%lmaxd+1)**2 ) COMPLEX sf2(stars%ng2),sf3(stars%ng3),ylm( (atoms%lmaxd+1)**2 )
REAL rcc(3),v1(DIMENSION%nspd),v2(DIMENSION%nspd),x(3),ri(3) REAL rcc(3),v1(DIMENSION%nspd),v2(DIMENSION%nspd),x(3),ri(3)
l_cdn = .FALSE. ! By default we assume that the input is a potential.
IF (potden%potdenType.LE.0) CALL juDFT_error('unknown potden type', calledby='checkdop')
IF (potden%potdenType.GT.1000) l_cdn = .TRUE. ! potdenTypes > 1000 are reserved for densities
! .. ! ..
! .. ! ..
#ifdef __TOS_BGQ__ #ifdef __TOS_BGQ__
...@@ -75,11 +80,11 @@ ...@@ -75,11 +80,11 @@
ENDIF ENDIF
v1(j) = 0.0 v1(j) = 0.0
DO k = 1,stars%ng3 DO k = 1,stars%ng3
v1(j) = v1(j) + REAL(fpw(k,jsp)*sf3(k))*stars%nstr(k) v1(j) = v1(j) + REAL(potden%pw(k,jsp)*sf3(k))*stars%nstr(k)
ENDDO ENDDO
ENDDO ENDDO
! ---> vacuum part ! ---> vacuum part
IF (cdn) THEN IF (l_cdn) THEN
WRITE (6,FMT=9000) ivac WRITE (6,FMT=9000) ivac
WRITE (16,FMT=9000) ivac WRITE (16,FMT=9000) ivac
ELSE ELSE
...@@ -91,18 +96,18 @@ ...@@ -91,18 +96,18 @@
CALL starf2(& CALL starf2(&
& sym%nop2,stars%ng2,stars%kv2,sym%mrot,sym%symor,sym%tau,p(1,j),sym%invtab,& & sym%nop2,stars%ng2,stars%kv2,sym%mrot,sym%symor,sym%tau,p(1,j),sym%invtab,&
& sf2)!keep & sf2)!keep
v2(j) = fz(1,ivac,jsp) v2(j) = potden%vacz(1,ivac,jsp)
DO k = 2,stars%ng2 DO k = 2,stars%ng2
v2(j) = v2(j) + REAL(fxy(1,k-1,ivac,jsp)*sf2(k))*stars%nstr2(k) v2(j) = v2(j) + REAL(potden%vacxy(1,k-1,ivac,jsp)*sf2(k))*stars%nstr2(k)
ENDDO ENDDO
ELSE ELSE
!-odim !-odim
v2(j) = fz(1,ivac,jsp) v2(j) = potden%vacz(1,ivac,jsp)
phi = angle(p(1,j),p(2,j)) phi = angle(p(1,j),p(2,j))
DO k = 2,oneD%odi%nq2 DO k = 2,oneD%odi%nq2
m = oneD%odi%kv(2,k) m = oneD%odi%kv(2,k)
gz = oneD%odi%kv(1,k) gz = oneD%odi%kv(1,k)
v2(j) = v2(j) + REAL(fxy(1,k-1,ivac,jsp)*& v2(j) = v2(j) + REAL(potden%vacxy(1,k-1,ivac,jsp)*&
& EXP(ic*m*phi)*EXP(ic*cell%bmat(3,3)*gz*p(3,j)))*oneD%odi%nst2(k) & EXP(ic*m*phi)*EXP(ic*cell%bmat(3,3)*gz*p(3,j)))*oneD%odi%nst2(k)
ENDDO ENDDO
!+odim !+odim
...@@ -136,11 +141,11 @@ ...@@ -136,11 +141,11 @@
! !
v1(j) = 0.0 v1(j) = 0.0
DO k = 1,stars%ng3 DO k = 1,stars%ng3
v1(j) = v1(j) + REAL(fpw(k,jsp)*sf3(k))*stars%nstr(k) v1(j) = v1(j) + REAL(potden%pw(k,jsp)*sf3(k))*stars%nstr(k)
ENDDO ENDDO
ENDDO ENDDO
! ----> m.t. part ! ----> m.t. part
IF (cdn) THEN IF (l_cdn) THEN
WRITE (6,FMT=9010) n WRITE (6,FMT=9010) n
WRITE (16,FMT=9010) n WRITE (16,FMT=9010) n
ELSE ELSE
...@@ -148,7 +153,7 @@ ...@@ -148,7 +153,7 @@
WRITE (16,FMT=8010) n WRITE (16,FMT=8010) n
ENDIF ENDIF
ir2 = 1.0 ir2 = 1.0
IF (cdn) ir2 = 1.0 / ( atoms%rmt(n)*atoms%rmt(n) ) IF (l_cdn) ir2 = 1.0 / ( atoms%rmt(n)*atoms%rmt(n) )
nd = atoms%ntypsy(na) nd = atoms%ntypsy(na)
nopa = atoms%ngopr(na) nopa = atoms%ngopr(na)
IF (oneD%odi%d1) THEN IF (oneD%odi%d1) THEN
...@@ -190,7 +195,7 @@ ...@@ -190,7 +195,7 @@
lm = ll1 + sphhar%mlh(mem,lh,nd) lm = ll1 + sphhar%mlh(mem,lh,nd)
s = s + REAL( sphhar%clnu(mem,lh,nd)* ylm(lm) ) s = s + REAL( sphhar%clnu(mem,lh,nd)* ylm(lm) )
ENDDO ENDDO
help = help + fr(atoms%jri(n),lh,n,jsp) * s help = help + potden%mt(atoms%jri(n),lh,n,jsp) * s
ENDDO ENDDO
v2(j) = help * ir2 v2(j) = help * ir2
IF (j.LE.8) THEN IF (j.LE.8) THEN
......
!--------------------------------------------------------------------------------
! 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_checkdopall
CONTAINS
SUBROUTINE checkDOPAll(input,dimension,sphhar,stars,atoms,sym,vacuum,oneD,&
cell,potden,ispin)
USE m_sphpts
USE m_checkdop
USE m_types
USE m_cylpts
USE m_points
IMPLICIT NONE
TYPE(t_input),INTENT(IN) :: input
TYPE(t_dimension),INTENT(IN) :: dimension
TYPE(t_sphhar),intent(in) :: sphhar
TYPE(t_stars),INTENT(IN) :: stars
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_sym),INTENT(IN) :: sym
TYPE(t_vacuum),INTENT(IN) :: vacuum
TYPE(t_oneD),INTENT(IN) :: oneD
TYPE(t_cell),INTENT(IN) :: cell
TYPE(t_potden),INTENT(IN) :: potden
INTEGER, INTENT(IN) :: ispin
INTEGER :: npd, nat, n, ivac
REAL :: signum
REAL :: xp(3,dimension%nspd)
IF ((input%film).AND.(.NOT.oneD%odi%d1)) THEN
!---> vacuum boundaries
npd = min(dimension%nspd,25)
CALL points(xp,npd)
DO ivac = 1,vacuum%nvac
signum = 3.0 - 2.0*ivac
xp(3,:npd) = signum*cell%z1/cell%amat(3,3)
CALL checkdop(xp,npd,0,0,ivac,1,ispin,dimension,atoms,&
sphhar,stars,sym,vacuum,cell,oneD,potden)
END DO
ELSE IF (oneD%odi%d1) THEN
!-odim
npd = min(dimension%nspd,25)
CALL cylpts(xp,npd,cell%z1)
CALL checkdop(xp,npd,0,0,ivac,1,ispin,dimension,atoms,&
sphhar,stars,sym,vacuum,cell,oneD,potden)
!+odim
END IF
!---> m.t. boundaries
nat = 1
DO n = 1, atoms%ntype
CALL sphpts(xp,dimension%nspd,atoms%rmt(n),atoms%pos(1,atoms%nat))
CALL checkdop(xp,dimension%nspd,n,nat,0,-1,ispin,&
dimension,atoms,sphhar,stars,sym,vacuum,cell,oneD,potden)
nat = nat + atoms%neq(n)
END DO
END SUBROUTINE checkDOPAll
END MODULE m_checkdopall
...@@ -20,12 +20,11 @@ MODULE m_constants ...@@ -20,12 +20,11 @@ MODULE m_constants
REAL, PARAMETER :: eVac0Default_const = -0.25 REAL, PARAMETER :: eVac0Default_const = -0.25
CHARACTER(len=9), PARAMETER :: version_const = 'fleur 27' CHARACTER(len=9), PARAMETER :: version_const = 'fleur 27'
INTEGER, PARAMETER :: POTDEN_TYPE_OTHER = 0 INTEGER, PARAMETER :: POTDEN_TYPE_OTHER = 0 ! POTDEN_TYPE <= 0 ==> undefined
INTEGER, PARAMETER :: POTDEN_TYPE_POTTOT = 1 INTEGER, PARAMETER :: POTDEN_TYPE_POTTOT = 1 ! 0 < POTDEN_TYPE <= 1000 ==> potential
INTEGER, PARAMETER :: POTDEN_TYPE_POTCOUL = 2 INTEGER, PARAMETER :: POTDEN_TYPE_POTCOUL = 2
INTEGER, PARAMETER :: POTDEN_TYPE_POTX = 3 INTEGER, PARAMETER :: POTDEN_TYPE_POTX = 3
INTEGER, PARAMETER :: POTDEN_TYPE_DEN = 1001 INTEGER, PARAMETER :: POTDEN_TYPE_DEN = 1001 ! 1000 < POTDEN_TYPE ==> density
CHARACTER(2),DIMENSION(0:103),PARAMETER :: namat_const=(/& CHARACTER(2),DIMENSION(0:103),PARAMETER :: namat_const=(/&
'va',' H','He','Li','Be',' B',' C',' N',' O',' F','Ne',& 'va',' H','He','Li','Be',' B',' C',' N',' O',' F','Ne',&
......
...@@ -12,31 +12,26 @@ MODULE m_qfix ...@@ -12,31 +12,26 @@ MODULE m_qfix
CONTAINS CONTAINS
SUBROUTINE qfix( stars,atoms,sym,vacuum,& SUBROUTINE qfix( stars,atoms,sym,vacuum,&
sphhar,input,cell,oneD,qpw,rhtxy,rho,rht,l_printData,force_fix,fix) sphhar,input,cell,oneD,den,l_printData,force_fix,fix)
USE m_types USE m_types
USE m_cdntot USE m_cdntot
USE m_xmlOutput USE m_xmlOutput
IMPLICIT NONE IMPLICIT NONE
! ..
! .. Scalar Arguments .. ! .. Scalar Arguments ..
TYPE(t_stars),INTENT(IN) :: stars TYPE(t_stars),INTENT(IN) :: stars
TYPE(t_atoms),INTENT(IN) :: atoms TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_sym),INTENT(IN) :: sym TYPE(t_sym),INTENT(IN) :: sym
TYPE(t_vacuum),INTENT(IN):: vacuum TYPE(t_vacuum),INTENT(IN) :: vacuum
TYPE(t_sphhar),INTENT(IN):: sphhar TYPE(t_sphhar),INTENT(IN) :: sphhar
TYPE(t_input),INTENT(IN) :: input TYPE(t_input),INTENT(IN) :: input
TYPE(t_oneD),INTENT(IN) :: oneD TYPE(t_oneD),INTENT(IN) :: oneD
TYPE(t_cell),INTENT(IN) :: cell TYPE(t_cell),INTENT(IN) :: cell
LOGICAL,INTENT(IN) :: l_printData,force_fix TYPE(t_potden),INTENT(INOUT) :: den
REAL, INTENT (OUT) :: fix LOGICAL,INTENT(IN) :: l_printData,force_fix
! .. REAL, INTENT (OUT) :: fix
! .. Array Arguments ..
COMPLEX,INTENT (INOUT) :: qpw(stars%ng3,input%jspins)
COMPLEX,INTENT (INOUT) :: rhtxy(vacuum%nmzxyd,oneD%odi%n2d-1,2,input%jspins)
REAL, INTENT (INOUT) :: rho(atoms%jmtd,0:sphhar%nlhd,atoms%ntype,input%jspins)
REAL, INTENT (INOUT) :: rht(vacuum%nmzd,2,input%jspins)
! ..
! .. Local Scalars .. ! .. Local Scalars ..
LOGICAL :: l_qfixfile,fixtotal LOGICAL :: l_qfixfile,fixtotal
LOGICAL :: l_firstcall=.true. LOGICAL :: l_firstcall=.true.
...@@ -57,7 +52,7 @@ CONTAINS ...@@ -57,7 +52,7 @@ CONTAINS
! In this case do nothing except when forced to fix! ! In this case do nothing except when forced to fix!
CALL cdntot( stars,atoms,sym, vacuum,input,cell,oneD,& CALL cdntot( stars,atoms,sym, vacuum,input,cell,oneD,&
qpw,rho,rht,.TRUE., qtot,qis) den%pw,den%mt,den%vacz,.TRUE., qtot,qis)
!The total nucleii charge !The total nucleii charge
zc=SUM(atoms%neq(:)*atoms%zatom(:)) zc=SUM(atoms%neq(:)*atoms%zatom(:))
...@@ -71,30 +66,30 @@ CONTAINS ...@@ -71,30 +66,30 @@ CONTAINS
DO n = 1,atoms%ntype DO n = 1,atoms%ntype
lh = sphhar%nlh(atoms%ntypsy(na)) lh = sphhar%nlh(atoms%ntypsy(na))
jm = atoms%jri(n) jm = atoms%jri(n)
rho(:jm,0:lh,n,:) = fix*rho(:jm,0:lh,n,:) den%mt(:jm,0:lh,n,:) = fix*den%mt(:jm,0:lh,n,:)
na = na + atoms%neq(n) na = na + atoms%neq(n)
ENDDO ENDDO
qpw(:stars%ng3,:) = fix*qpw(:stars%ng3,:) den%pw(:stars%ng3,:) = fix*den%pw(:stars%ng3,:)
IF (input%film) THEN IF (input%film) THEN
rht(:vacuum%nmz,:vacuum%nvac,:) = fix*rht(:vacuum%nmz,:vacuum%nvac,:) den%vacz(:vacuum%nmz,:vacuum%nvac,:) = fix*den%vacz(:vacuum%nmz,:vacuum%nvac,:)
rhtxy(:vacuum%nmzxy,:oneD%odi%nq2-1,:vacuum%nvac,:) = fix*& den%vacxy(:vacuum%nmzxy,:oneD%odi%nq2-1,:vacuum%nvac,:) = fix*&
rhtxy(:vacuum%nmzxy,:oneD%odi%nq2-1,:vacuum%nvac,:) den%vacxy(:vacuum%nmzxy,:oneD%odi%nq2-1,:vacuum%nvac,:)
END IF END IF
WRITE (6,FMT=8000) zc,fix WRITE (6,FMT=8000) zc,fix
IF (ABS(fix-1.0)<1.E-6) RETURN !no second calculation of cdntot as nothing was fixed IF (ABS(fix-1.0)<1.E-6) RETURN !no second calculation of cdntot as nothing was fixed
CALL openXMLElementNoAttributes('fixedCharges') CALL openXMLElementNoAttributes('fixedCharges')
CALL cdntot( stars,atoms,sym, vacuum,input,cell,oneD,& CALL cdntot( stars,atoms,sym, vacuum,input,cell,oneD,&
qpw,rho,rht,l_printData, qtot,qis) den%pw,den%mt,den%vacz,l_printData, qtot,qis)
CALL closeXMLElement('fixedCharges')