Removed an old tentative to improve the noncollinear-gga case.

git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@6331 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
dalcorso 2010-01-29 16:52:32 +00:00
parent b4658b2b37
commit 28209081c6
4 changed files with 0 additions and 212 deletions

View File

@ -5,7 +5,6 @@
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
#define __OLD_NONCOLIN_GGA
!-----------------------------------------------------------------------
subroutine setup_dgc
!-----------------------------------------------------------------------
@ -70,7 +69,6 @@ subroutine setup_dgc
fac = 1.d0 / DBLE (nspin_gga)
IF (noncolin.and.domag) THEN
allocate(rhogout(ngm,nspin_mag))
#ifdef __OLD_NONCOLIN_GGA
call compute_rho(rho%of_r,rhoout,segni,nrxx)
DO is = 1, nspin_gga
!
@ -87,42 +85,6 @@ subroutine setup_dgc
rhogout(1,is), ngm, g, nl, grho(1,1,is) )
!
END DO
#else
call compute_rho_new(rho%of_r,rhoout,segni,nrxx,ux)
do is=1,nspin_mag
rhogout(:,is) = rho%of_g(:,is)
enddo
if (nlcc_any) then
rhogout(:,1) = rhog_core(:) + rho%of_g(:,1)
endif
do is = 1, nspin_mag
call gradrho (nrx1, nrx2, nrx3, nr1, nr2, nr3, nrxx, rhogout(1, is), &
ngm, g, nl, gmag (1, 1, is) )
enddo
DO is=1,nspin_gga
IF (is==1) seg0=0.5_dp
IF (is==2) seg0=-0.5_dp
DO ipol=1,3
grho(ipol,:,is) = 0.5_dp*gmag(ipol,:,1)
DO ir=1,nrxx
seg=seg0*segni(ir)
rhoout(ir,is) = fac*rho_core(ir) + 0.5_dp*rho%of_r(ir,1)
amag=sqrt(rho%of_r(ir,2)**2+rho%of_r(ir,3)**2+rho%of_r(ir,4)**2)
IF (amag>1.d-12) THEN
rhoout(ir,is)=rhoout(ir,is)+seg*amag
DO jpol=2,4
grho(ipol,ir,is)=grho(ipol,ir,is)+ seg*rho%of_r(ir,jpol)* &
gmag(ipol,ir,jpol)/amag
END DO
END IF
END DO
END DO
END DO
! write(6,*) 'setup dgc'
! do k=2,2
! write(6,'(3f20.5)') gmag(3,k,1), grho(3,k,1)
! enddo
#endif
DEALLOCATE(rhogout)
ELSE
do is = 1, nspin_gga

View File

@ -46,7 +46,6 @@ compute_becsum.o \
compute_deff.o \
compute_dip.o \
compute_rho.o \
compute_rho_new.o \
compute_fes_grads.o \
compute_scf.o \
compute_qdipol.o \

View File

@ -1,106 +0,0 @@
!
! Copyright (C) 2005 PWSCF group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
SUBROUTINE compute_rho_new(rho,rhoout,segni,nrxx,ux)
!
! This subroutine diagonalizes the spin density matrix and gives as output
! the spin up and spin down components of the charge
!
USE kinds, ONLY : dp
USE constants, ONLY: pi
USE io_global, ONLY : stdout
implicit none
REAL(DP) :: ux(3) ! fixed direction to calculate signs
INTEGER :: nrxx ! input: the dimension of the mesh
REAL(DP) :: rho(nrxx,4),segni(nrxx),rhoout(nrxx,2)
! input: the four components of the charge
! output: keep track of the spin direction
! output: the spin up and spin down charge
COMPLEX(DP) :: re(2,2) ! the rotated density matrix
REAL(DP) :: eps,amag,uxmod ! pi, a small number and the atan function
REAL(DP) :: m(3) ! magnetization
REAL(DP) :: ux0(3) ! magnetization
LOGICAL :: negative
INTEGER :: ir, i ! counter on mesh points
INTEGER :: it,counter(3)
eps=1.d-12
uxmod=sqrt(ux(1)**2+ux(2)**2+ux(3)**2)
if (uxmod<eps) then
do i=1,3
ux(i)=i
ux0(i)=i
enddo
counter(3)=0
do it=1,3
counter(it)=0
do ir=1,nrxx
amag=SQRT(rho(ir,2)**2+rho(ir,3)**2+rho(ir,4)**2)
if (amag > eps) then
do i=1,3
m(i)=rho(ir,i+1)/amag
enddo
if (it==2) then
if (ABS(ux(1)*m(1)+ux(2)*m(2)+ux(3)*m(3)).lt.1.d-3) then
do i=1,3
ux(i)=ux(i)+0.5d0*m(i)
enddo
uxmod=SQRT(ux(1)**2+ux(2)**2+ux(3)**2)
do i=1,3
ux(i)=ux(i)/uxmod
enddo
endif
else
if (ABS(ux(1)*m(1)+ux(2)*m(2)+ux(3)*m(3)).lt.1.d-3) &
counter(it)=counter(it)+1
endif
endif
enddo
if (counter(1)==0) goto 100
enddo
if (counter(1).lt.counter(3)) then
do i=1,3
ux(i)=ux0(i)
enddo
endif
! WRITE( stdout,'(5x,"it, counter:",3i5)') it,counter(1),counter(3)
WRITE( stdout,'(5x,"Fixed quantization axis: ", 3f12.6)') &
ux(1), ux(2), ux(3)
endif
100 continue
do ir=1,nrxx
amag=SQRT(rho(ir,2)**2+rho(ir,3)**2+rho(ir,4)**2)
re(1,1)=0.5_DP*(rho(ir,1)+amag)
re(2,2)=0.5_DP*(rho(ir,1)-amag)
if (amag.lt.eps) then
negative=.false.
else
do i=1,3
m(i)=rho(ir,i+1)/amag
enddo
negative=(m(1)*ux(1)+m(2)*ux(2)+m(3)*ux(3))>0.d0
endif
if (negative) then
rhoout(ir,2)= DBLE(re(1,1))
rhoout(ir,1)= DBLE(re(2,2))
segni(ir)=-1.d0
else
rhoout(ir,1)= DBLE(re(1,1))
rhoout(ir,2)= DBLE(re(2,2))
segni(ir)=1.d0
endif
enddo
return
end subroutine compute_rho_new

View File

@ -6,7 +6,6 @@
! or http://www.gnu.org/copyleft/gpl.txt .
!
!
#define __OLD_NONCOLIN_GGA
!----------------------------------------------------------------------------
SUBROUTINE gradcorr( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
!----------------------------------------------------------------------------
@ -67,10 +66,6 @@ SUBROUTINE gradcorr( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
ALLOCATE( vgg( nrxx, nspin0 ) )
ALLOCATE( vsave( nrxx, nspin ) )
ALLOCATE( segni( nrxx ) )
#ifdef __OLD_NONCOLIN_GGA
#else
ALLOCATE( gmag( 3, nrxx, nspin) )
#endif
vsave=v
v=0.d0
ENDIF
@ -79,7 +74,6 @@ SUBROUTINE gradcorr( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
!
! ... calculate the gradient of rho + rho_core in real space
!
#ifdef __OLD_NONCOLIN_GGA
IF ( nspin == 4 .AND. domag ) THEN
!
CALL compute_rho(rho,rhoout,segni,nrxx)
@ -110,53 +104,6 @@ SUBROUTINE gradcorr( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
rhogsum(1,is), ngm, g, nl, grho(1,1,is) )
!
END DO
#else
IF ( nspin == 4 .AND. domag ) THEN
!
CALL compute_rho_new(rho,rhoout,segni,nrxx,ux)
!
rhogsum(:,1) =rhog_core(:) + rhog(:,1)
!
CALL gradrho( nrx1, nrx2, nrx3, nr1, nr2, nr3, nrxx, &
rhogsum, ngm, g, nl, gmag )
DO is = 2, nspin
rhogsum(:,1) = rhog(:,is)
!
CALL gradrho( nrx1, nrx2, nrx3, nr1, nr2, nr3, nrxx, &
rhogsum, ngm, g, nl, gmag(1,1,is) )
END DO
DO is=1,nspin0
IF (is==1) seg=0.5_dp
IF (is==2) seg=-0.5_dp
DO ipol=1,3
grho(ipol,:,is) = 0.5_dp*gmag(ipol,:,1)
ENDDO
DO ir=1,nrxx
amag=sqrt(rho(ir,2)**2+rho(ir,3)**2+rho(ir,4)**2)
rhoout(ir,is) = fac*rho_core(ir) + 0.5_dp*rho(ir,1)
IF (amag>1.d-12) THEN
rhoout(ir,is)= rhoout(ir,is) + segni(ir)*seg*amag
DO ipol=1,3
DO jpol=2,4
grho(ipol,ir,is)=grho(ipol,ir,is)+ segni(ir)*seg*rho(ir,jpol)* &
gmag(ipol,ir,jpol)/amag
END DO
END DO
END IF
END DO
END DO
ELSE
DO is = 1, nspin0
!
rhoout(:,is) = fac * rho_core(:) + rho(:,is)
rhogsum(:,is) = fac * rhog_core(:) + rhog(:,is)
!
CALL gradrho( nrx1, nrx2, nrx3, nr1, nr2, nr3, nrxx, &
rhogsum(1,is), ngm, g, nl, grho(1,1,is) )
!
END DO
END IF
#endif
!
DEALLOCATE( rhogsum )
!
@ -241,20 +188,10 @@ SUBROUTINE gradcorr( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
!
zeta = ( rhoout(k,1) - rhoout(k,2) ) / rh
if (nspin.eq.4.and.domag) zeta=abs(zeta)*segni(k)
#ifdef __OLD_NONCOLIN_GGA
!
grh2 = ( grho(1,k,1) + grho(1,k,2) )**2 + &
( grho(2,k,1) + grho(2,k,2) )**2 + &
( grho(3,k,1) + grho(3,k,2) )**2
#else
if (nspin==4) then
grh2= gmag(1,k,1)**2+ gmag(2,k,1)**2+gmag(3,k,1)**2
else
grh2 = ( grho(1,k,1) + grho(1,k,2) )**2 + &
( grho(2,k,1) + grho(2,k,2) )**2 + &
( grho(3,k,1) + grho(3,k,2) )**2
endif
#endif
!
CALL gcc_spin( rh, zeta, grh2, sc, v1cup, v1cdw, v2c )
!
@ -353,10 +290,6 @@ SUBROUTINE gradcorr( rho, rhog, rho_core, rhog_core, etxc, vtxc, v )
DEALLOCATE( vgg )
DEALLOCATE( vsave )
DEALLOCATE( segni )
#ifdef __OLD_NONCOLIN_GGA
#else
DEALLOCATE( gmag )
#endif
ENDIF
!
RETURN