Skip to content

GitLab

Commit dfec2957 authored by Gregor Michalicek's avatar Gregor Michalicek
Browse files

Introduce potden type to some files + some cleanup

parent 9b164562
Expand all Hide whitespace changes
Inline Side-by-side
Showing with 888 additions and 955 deletions
global/types.F90
View file @ dfec2957
......@@ -986,11 +986,11 @@ CONTAINS
ALLOCATE(pd%vacxy(nmzxyd,n2d-1,2,jsp),stat=err(4))
ENDIF
IF (ANY(err>0)) CALL judft_error("Not enough memory allocating potential or density")
pd%pw=0.0
pd%pw=CMPLX(0.0,0.0)
pd%mt=0.0
IF (PRESENT(nmzd)) THEN
pd%vacz=0.0
pd%vacxy=0.0
pd%vacxy=CMPLX(0.0,0.0)
ENDIF
END SUBROUTINE init_potden_simple
......
optional/cdnsp.f90
View file @ dfec2957
......@@ -18,7 +18,7 @@
& stars,sym,oneD,cell,DIMENSION)
USE m_intgr, ONLY : intgr3
USE m_constants, ONLY : pi_const
USE m_constants
USE m_cdn_io
USE m_types
IMPLICIT NONE
......@@ -32,19 +32,19 @@
TYPE(t_oneD),INTENT(IN) :: oneD
TYPE(t_cell),INTENT(IN) :: cell
TYPE(t_dimension),INTENT(IN) :: DIMENSION
! ..
! local type instances
TYPE(t_potden) :: den
! .. Local Scalars ..
REAL dummy,p,pp,qtot1,qtot2,spmtot,qval,sfp,fermiEnergyTemp
INTEGER i,iter,ivac,j,k,lh,n,na,jsp_new
INTEGER i,ivac,j,k,lh,n,na,jsp_new
INTEGER ios
LOGICAL n_exist,l_qfix
! ..
! .. Local Arrays ..
REAL rhoc(atoms%jmtd,atoms%ntype,dimension%jspd)
REAL tec(atoms%ntype,dimension%jspd),qintc(atoms%ntype,dimension%jspd)
COMPLEX :: cdom(1),cdomvz(1,1),cdomvxy(1,1,1)
COMPLEX, ALLOCATABLE :: qpw(:,:),rhtxy(:,:,:,:)
REAL , ALLOCATABLE :: rho(:,:,:,:),rht(:,:,:)
CHARACTER(len=140), ALLOCATABLE :: clines(:)
CHARACTER(len=140) :: lineread
! ..
......@@ -54,13 +54,17 @@
IF (input%jspins/=2) CALL juDFT_error&
& ("cdnsp: set jspins = 2 and remove fl7para!",calledby&
& ="cdnsp")
ALLOCATE ( rho(atoms%jmtd,0:sphhar%nlhd,atoms%ntype,input%jspins),qpw(stars%ng3,input%jspins) )
ALLOCATE ( rhtxy(vacuum%nmzxy,stars%ng2-1,2,input%jspins),rht(vacuum%nmz,2,input%jspins) )
CALL den%init(stars,atoms,sphhar,vacuum,oneD,DIMENSION%jspd,.FALSE.,POTDEN_TYPE_DEN)
ALLOCATE (den%cdom(1),den%cdomvz(1,1),den%cdomvxy(1,1,1))
ALLOCATE (den%mmpMat(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,MAX(1,atoms%n_u),input%jspins))
den%mmpMat = CMPLX(0.0,0.0)
input%jspins=1
CALL readCoreDensity(input,atoms,dimension,rhoc,tec,qintc)
CALL readDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,CDN_ARCHIVE_TYPE_CDN1_const,&
CDN_INPUT_DEN_const,0,fermiEnergyTemp,l_qfix,iter,rho,qpw,rht,rhtxy,cdom,cdomvz,cdomvxy)
CDN_INPUT_DEN_const,0,fermiEnergyTemp,l_qfix,den%iter,den%mt,den%pw,den%vacz,&
den%vacxy,den%cdom,den%cdomvz,den%cdomvxy)
input%jspins=2
qval = 0.
......@@ -71,45 +75,46 @@
!
DO n = 1,atoms%ntype
DO j = 1,atoms%jri(n)
rho(j,0,n,1) = rho(j,0,n,1) - rhoc(j,n,1)/sfp
den%mt(j,0,n,1) = den%mt(j,0,n,1) - rhoc(j,n,1)/sfp
ENDDO
! WRITE (16,FMT='(8f10.4)') (rho(i,0,n,1),i=1,16)
CALL intgr3(rho(1,0,n,1),atoms%rmsh(1,n),atoms%dx(n),atoms%jri(n),qval)
! WRITE (16,FMT='(8f10.4)') (den%mt(i,0,n,1),i=1,16)
CALL intgr3(den%mt(1,0,n,1),atoms%rmsh(1,n),atoms%dx(n),atoms%jri(n),qval)
p = (atoms%bmu(n)+sfp*qval)/ (2.*sfp*qval)
pp = 1. - p
DO j = 1,atoms%jri(n)
rho(j,0,n,jsp_new) = pp*rho(j,0,n,1) + rhoc(j,n,1)/ (2.*sfp)
rho(j,0,n,1) = p*rho(j,0,n,1) + rhoc(j,n,1)/ (2.*sfp)
den%mt(j,0,n,jsp_new) = pp*den%mt(j,0,n,1) + rhoc(j,n,1)/ (2.*sfp)
den%mt(j,0,n,1) = p*den%mt(j,0,n,1) + rhoc(j,n,1)/ (2.*sfp)
ENDDO
DO lh = 1,sphhar%nlh(atoms%ntypsy(na))
DO j = 1,atoms%jri(n)
rho(j,lh,n,jsp_new) = pp*rho(j,lh,n,1)
rho(j,lh,n,1) = p*rho(j,lh,n,1)
den%mt(j,lh,n,jsp_new) = pp*den%mt(j,lh,n,1)
den%mt(j,lh,n,1) = p*den%mt(j,lh,n,1)
ENDDO
ENDDO
na = na + atoms%neq(n)
ENDDO
DO k = 1,stars%ng3
qpw(k,jsp_new) = 0.5 * qpw(k,1)
qpw(k,1) = qpw(k,jsp_new)
den%pw(k,jsp_new) = 0.5 * den%pw(k,1)
den%pw(k,1) = den%pw(k,jsp_new)
ENDDO
IF (input%film) THEN
DO ivac = 1,vacuum%nvac
DO j = 1, vacuum%nmz
rht(j,ivac,jsp_new) = 0.5 * rht(j,ivac,1)
rht(j,ivac,1) = rht(j,ivac,jsp_new)
den%vacz(j,ivac,jsp_new) = 0.5 * den%vacz(j,ivac,1)
den%vacz(j,ivac,1) = den%vacz(j,ivac,jsp_new)
ENDDO
DO k = 2, stars%ng2
DO j = 1,vacuum%nmzxy
rhtxy(j,k-1,ivac,jsp_new) = 0.5 * rhtxy(j,k-1,ivac,1)
rhtxy(j,k-1,ivac,1) = rhtxy(j,k-1,ivac,jsp_new)
den%vacxy(j,k-1,ivac,jsp_new) = 0.5 * den%vacxy(j,k-1,ivac,1)
den%vacxy(j,k-1,ivac,1) = den%vacxy(j,k-1,ivac,jsp_new)
ENDDO
ENDDO
ENDDO
ENDIF
! ----> write the spin-polarized density
CALL writeDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,CDN_ARCHIVE_TYPE_CDN1_const,&
CDN_INPUT_DEN_const,0,-1.0,0.0,.FALSE.,iter,rho,qpw,rht,rhtxy,cdom,cdomvz,cdomvxy)
CDN_INPUT_DEN_const,0,-1.0,0.0,.FALSE.,den%iter,den%mt,den%pw,den%vacz,den%vacxy,&
den%cdom,den%cdomvz,den%cdomvxy)
!
! -----> This part is only used for testing th e magnetic moment in
! -----> each sphere
......@@ -117,8 +122,8 @@
DO n = 1,atoms%ntype
qtot1=0.00
qtot2=0.00
CALL intgr3(rho(1,0,n,1),atoms%rmsh(1,n),atoms%dx(n),atoms%jri(n),qtot1)
CALL intgr3(rho(1,0,n,jsp_new),atoms%rmsh(1,n),atoms%dx(n),atoms%jri(n),qtot2)
CALL intgr3(den%mt(1,0,n,1),atoms%rmsh(1,n),atoms%dx(n),atoms%jri(n),qtot1)
CALL intgr3(den%mt(1,0,n,jsp_new),atoms%rmsh(1,n),atoms%dx(n),atoms%jri(n),qtot2)
spmtot=sfp*(qtot1-qtot2)
WRITE (6,'('' moment in sphere '',2x,'':'',f8.4)') spmtot
ENDDO
......@@ -150,7 +155,7 @@
WRITE (40,'(a)') TRIM(clines(i))
ENDDO
DEALLOCATE (clines,rho,qpw,rhtxy,rht)
DEALLOCATE (clines)
CLOSE(40)
ENDIF
!
......
optional/flipcdn.f90
View file @ dfec2957
......@@ -4,7 +4,7 @@
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------
MODULE m_flipcdn
MODULE m_flipcdn
! *******************************************************
! this subroutine reads the charge density and flips the
! magnetic moment within the m.t.sphere for each atom
......@@ -17,170 +17,164 @@
!
! Extension to multiple U per atom type by G.M. 2017
! *******************************************************
CONTAINS
SUBROUTINE flipcdn(&
& atoms,input,vacuum,sphhar,&
& stars,sym,oneD,cell,&
& l_noco)
USE m_cdn_io
USE m_types
IMPLICIT NONE
TYPE(t_stars),INTENT(IN) :: stars
TYPE(t_vacuum),INTENT(IN) :: vacuum
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_sphhar),INTENT(IN) :: sphhar
TYPE(t_input),INTENT(IN) :: input
TYPE(t_sym),INTENT(IN) :: sym
TYPE(t_oneD),INTENT(IN) :: oneD
TYPE(t_cell),INTENT(IN) :: cell
LOGICAL,INTENT(IN) :: l_noco
! .. Local Scalars ..
REAL rhodummy,rhodumms,fermiEnergyTemp
INTEGER i,iter,nt,j,lh,na,mp,ispin,urec,itype,m,i_u
INTEGER archiveType
CHARACTER(len=8) iop,dop
LOGICAL n_exist,l_qfix
! ..
! .. Local Arrays ..
COMPLEX, ALLOCATABLE :: n_mmp(:,:,:,:),qpw(:,:),rhtxy(:,:,:,:)
REAL , ALLOCATABLE :: rho(:,:,:,:),rht(:,:,:)
COMPLEX, ALLOCATABLE :: cdom(:),cdomvz(:,:),cdomvxy(:,:,:)
CHARACTER(len=80), ALLOCATABLE :: clines(:)
CHARACTER(len=8) name(10)
! ..
!atoms%jmtd = MAXVAL(atoms%jri(:))
!sphhar%nlhd = MAXVAL(sphhar%nlh(:))
ALLOCATE ( rho(atoms%jmtd,0:sphhar%nlhd,atoms%ntype,input%jspins),qpw(stars%ng3,input%jspins) )
ALLOCATE ( rhtxy(vacuum%nmzxy,stars%ng2-1,2,input%jspins),rht(vacuum%nmz,2,input%jspins) )
archiveType = CDN_ARCHIVE_TYPE_CDN1_const
IF (l_noco) THEN
ALLOCATE( cdom(stars%ng3) )
ALLOCATE( cdomvz(vacuum%nmz,2),cdomvxy(vacuum%nmzxy,stars%ng2-1,2) )
archiveType = CDN_ARCHIVE_TYPE_NOCO_const
ENDIF
! ---> read the charge density
CALL readDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,&
CDN_INPUT_DEN_const,0,fermiEnergyTemp,l_qfix,iter,rho,qpw,rht,rhtxy,cdom,cdomvz,cdomvxy)
! ---> flip cdn for each atom with nflip=-1
!
na = 1
DO itype = 1, atoms%ntype
IF (atoms%nflip(itype).EQ.-1) THEN
! ---> spherical and non-spherical m.t. charge density
DO lh = 0,sphhar%nlh(atoms%ntypsy(na))
DO j = 1,atoms%jri(itype)
rhodummy = rho(j,lh,itype,1)
rho(j,lh,itype,1) = rho(j,lh,itype,input%jspins)
rho(j,lh,itype,input%jspins) = rhodummy
END DO
END DO
ELSE IF (atoms%nflip(itype).EQ.-2) THEN
DO lh = 0,sphhar%nlh(atoms%ntypsy(na))
DO j = 1,atoms%jri(itype)
rhodummy = rho(j,lh,itype,1) + rho(j,lh,itype,input%jspins)
rhodumms = rho(j,lh,itype,1) - rho(j,lh,itype,input%jspins)
rho(j,lh,itype,1) = 0.5 * (rhodummy + atoms%bmu(itype)*rhodumms)
rho(j,lh,itype,input%jspins) = 0.5 * (rhodummy - atoms%bmu(itype)*rhodumms )
END DO
END DO
END IF
na = na + atoms%neq(itype)
ENDDO
! ----> write the spin-polarized density
CALL writeDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,CDN_INPUT_DEN_const,&
0,-1.0,0.0,.FALSE.,iter,rho,qpw,rht,rhtxy,cdom,cdomvz,cdomvxy)
!
! for lda+U: flip n-matrix
!
IF (atoms%n_u.GT.0) THEN
INQUIRE (file='n_mmp_mat',exist=n_exist)
IF (n_exist) THEN
OPEN (69,file='n_mmp_mat',status='old',form='formatted')
ALLOCATE (n_mmp(-3:3,-3:3,atoms%n_u,2))
READ (69,9000) n_mmp
! flip ...
DO i_u = 1, atoms%n_u
itype = atoms%lda_u(i_u)%atomType
IF (atoms%nflip(itype).EQ.-1) THEN
DO m = -3,3
DO mp = -3,3
rhodummy = n_mmp(m,mp,i_u,1)
n_mmp(m,mp,i_u,1) = n_mmp(m,mp,i_u,input%jspins)
n_mmp(m,mp,i_u,input%jspins) = rhodummy
END DO
END DO
ELSE IF (atoms%nflip(itype).EQ.-2) THEN
DO m = -3,3
DO mp = -3,3
rhodummy = n_mmp(m,mp,i_u,1) + n_mmp(m,mp,i_u,input%jspins)
rhodumms = n_mmp(m,mp,i_u,1) - n_mmp(m,mp,i_u,input%jspins)
n_mmp(m,mp,i_u,1) = 0.5 * (rhodummy + atoms%bmu(itype) * rhodumms)
n_mmp(m,mp,i_u,input%jspins) = 0.5 * (rhodummy - atoms%bmu(itype) * rhodumms)
END DO
END DO
END IF
END DO
! flip ...
REWIND (69)
WRITE (69,9000) n_mmp
9000 FORMAT(7f20.13)
!
DEALLOCATE ( n_mmp )
ENDIF
ENDIF
!-lda+U
!
!---> read enpara and flip lines
!
INQUIRE(file='enpara',exist=n_exist)
IF (n_exist) THEN
OPEN(40,file ='enpara',status='old',form='formatted')
j = 2
DO itype = 1, atoms%ntype
j = j + 1
IF (atoms%nlo(itype)>0) j = j + 2
ENDDO
IF (input%film) j = j + 1
ALLOCATE (clines(2*j))
DO i = 1, 2*j
READ (40,'(a)') clines(i)
ENDDO
REWIND 40
i = 0
DO ispin = 1,input%jspins
i = i + 2
WRITE (40,'(a)') TRIM(clines(i-1))
WRITE (40,'(a)') TRIM(clines(i))
DO itype = 1, atoms%ntype
i = i + 1
m = i
IF (atoms%nflip(itype)==-1) m = MOD(i+j,2*j)
IF (m==0) m = 2*j
WRITE (40,'(a)') TRIM(clines(m))
IF (atoms%nlo(itype)>0) THEN
WRITE (40,'(a)') TRIM(clines(m+1))
WRITE (40,'(a)') TRIM(clines(m+2))
i = i + 2
ENDIF
ENDDO
IF (input%film) THEN
i = i + 1
WRITE (40,'(a)') TRIM(clines(i))
ENDIF
ENDDO
DEALLOCATE (clines,rho,qpw,rhtxy,rht)
IF (l_noco) THEN
DEALLOCATE (cdom,cdomvz,cdomvxy)
ENDIF
CLOSE(40)
ENDIF
END SUBROUTINE flipcdn
END MODULE m_flipcdn
CONTAINS
SUBROUTINE flipcdn(atoms,input,vacuum,sphhar,&
stars,sym,oneD,cell,l_noco)
USE m_constants
USE m_cdn_io
USE m_types
IMPLICIT NONE
TYPE(t_stars),INTENT(IN) :: stars
TYPE(t_vacuum),INTENT(IN) :: vacuum
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_sphhar),INTENT(IN) :: sphhar
TYPE(t_input),INTENT(INOUT) :: input
TYPE(t_sym),INTENT(IN) :: sym
TYPE(t_oneD),INTENT(IN) :: oneD
TYPE(t_cell),INTENT(IN) :: cell
LOGICAL,INTENT(IN) :: l_noco
! Local type instance
TYPE(t_potden) :: den
! Local Scalars
REAL :: rhodummy,rhodumms,fermiEnergyTemp
INTEGER :: i,nt,j,lh,na,mp,ispin,urec,itype,m,i_u
INTEGER :: archiveType
LOGICAL :: n_exist,l_qfix,l_error
! Local Arrays
CHARACTER(len=80), ALLOCATABLE :: clines(:)
CALL den%init(stars,atoms,sphhar,vacuum,oneD,input%jspins,.FALSE.,POTDEN_TYPE_DEN)
IF(l_noco) THEN
ALLOCATE(den%cdom(stars%ng3))
ALLOCATE(den%cdomvz(vacuum%nmz,2),den%cdomvxy(vacuum%nmzxy,stars%ng2-1,2))
archiveType = CDN_ARCHIVE_TYPE_NOCO_const
ELSE
ALLOCATE (den%cdom(1),den%cdomvz(1,1),den%cdomvxy(1,1,1))
archiveType = CDN_ARCHIVE_TYPE_CDN1_const
END IF
ALLOCATE (den%mmpMat(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,MAX(1,atoms%n_u),input%jspins))
den%mmpMat = CMPLX(0.0,0.0)
! read the charge density
CALL readDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,&
CDN_INPUT_DEN_const,0,fermiEnergyTemp,l_qfix,den%iter,den%mt,den%pw,&
den%vacz,den%vacxy,den%cdom,den%cdomvz,den%cdomvxy)
! flip cdn for each atom with nflip=-1
na = 1
DO itype = 1, atoms%ntype
IF (atoms%nflip(itype).EQ.-1) THEN
! spherical and non-spherical m.t. charge density
DO lh = 0,sphhar%nlh(atoms%ntypsy(na))
DO j = 1,atoms%jri(itype)
rhodummy = den%mt(j,lh,itype,1)
den%mt(j,lh,itype,1) = den%mt(j,lh,itype,input%jspins)
den%mt(j,lh,itype,input%jspins) = rhodummy
END DO
END DO
ELSE IF (atoms%nflip(itype).EQ.-2) THEN
DO lh = 0,sphhar%nlh(atoms%ntypsy(na))
DO j = 1,atoms%jri(itype)
rhodummy = den%mt(j,lh,itype,1) + den%mt(j,lh,itype,input%jspins)
rhodumms = den%mt(j,lh,itype,1) - den%mt(j,lh,itype,input%jspins)
den%mt(j,lh,itype,1) = 0.5 * (rhodummy + atoms%bmu(itype)*rhodumms)
den%mt(j,lh,itype,input%jspins) = 0.5 * (rhodummy - atoms%bmu(itype)*rhodumms )
END DO
END DO
END IF
na = na + atoms%neq(itype)
END DO
! for LDA+U: flip density matrix
IF (isDensityMatrixPresent().AND.atoms%n_u>0) THEN
CALL readDensityMatrix(input,atoms,den%mmpMat,l_error)
DO i_u = 1, atoms%n_u
itype = atoms%lda_u(i_u)%atomType
IF (atoms%nflip(itype).EQ.-1) THEN
DO m = -3,3
DO mp = -3,3
rhodummy = den%mmpMat(m,mp,i_u,1)
den%mmpMat(m,mp,i_u,1) = den%mmpMat(m,mp,i_u,input%jspins)
den%mmpMat(m,mp,i_u,input%jspins) = rhodummy
END DO
END DO
ELSE IF (atoms%nflip(itype).EQ.-2) THEN
DO m = -3,3
DO mp = -3,3
rhodummy = den%mmpMat(m,mp,i_u,1) + den%mmpMat(m,mp,i_u,input%jspins)
rhodumms = den%mmpMat(m,mp,i_u,1) - den%mmpMat(m,mp,i_u,input%jspins)
den%mmpMat(m,mp,i_u,1) = 0.5 * (rhodummy + atoms%bmu(itype) * rhodumms)
den%mmpMat(m,mp,i_u,input%jspins) = 0.5 * (rhodummy - atoms%bmu(itype) * rhodumms)
END DO
END DO
END IF
END DO
END IF
! write the spin-polarized density
CALL writeDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,CDN_INPUT_DEN_const,&
0,-1.0,0.0,.FALSE.,den%iter,den%mt,den%pw,den%vacz,den%vacxy,den%cdom,den%cdomvz,den%cdomvxy)
! write density matrix for LDA+U
IF (isDensityMatrixPresent().AND.atoms%n_u>0) THEN
OPEN (69,file='n_mmp_mat',status='replace',form='formatted')
WRITE (69,'(7f20.13)') den%mmpMat(:,:,:,:)
IF (input%ldauLinMix) THEN
WRITE (69,'(2(a6,f5.3))') 'alpha=',input%ldauMixParam,'spinf=',input%ldauSpinf
END IF
CLOSE (69)
END IF
! read enpara and flip lines
INQUIRE(file='enpara',exist=n_exist)
IF (n_exist) THEN
OPEN(40,file ='enpara',status='old',form='formatted')
j = 2
DO itype = 1, atoms%ntype
j = j + 1
IF (atoms%nlo(itype)>0) j = j + 2
END DO
IF (input%film) j = j + 1
ALLOCATE (clines(2*j))
DO i = 1, 2*j
READ (40,'(a)') clines(i)
END DO
REWIND 40
i = 0
DO ispin = 1,input%jspins
i = i + 2
WRITE (40,'(a)') TRIM(clines(i-1))
WRITE (40,'(a)') TRIM(clines(i))
DO itype = 1, atoms%ntype
i = i + 1
m = i
IF (atoms%nflip(itype)==-1) m = MOD(i+j,2*j)
IF (m==0) m = 2*j
WRITE (40,'(a)') TRIM(clines(m))
IF (atoms%nlo(itype)>0) THEN
WRITE (40,'(a)') TRIM(clines(m+1))
WRITE (40,'(a)') TRIM(clines(m+2))
i = i + 2
END IF
END DO
IF (input%film) THEN
i = i + 1
WRITE (40,'(a)') TRIM(clines(i))
END IF
END DO
DEALLOCATE (clines)
CLOSE(40)
END IF
END SUBROUTINE flipcdn
END MODULE m_flipcdn
optional/pldngen.f90
View file @ dfec2957
This diff is collapsed.
optional/stden.f90
View file @ dfec2957
This diff is collapsed.
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment