Commit fbd7a926 authored by Matthias Redies's avatar Matthias Redies

fixed jspin=2 initial guess for libxc

parent 4c3403d8
......@@ -219,9 +219,8 @@
ENDDO
ENDDO
IF (xcpot%is_gga()) THEN
CALL potl0(&
xcpot,DIMENSION%msh,DIMENSION%jspd,input%jspins,n,&
atoms%dx(ntyp),rad,rhoss, vxc)
CALL potl0(xcpot,DIMENSION%msh,DIMENSION%jspd,input%jspins,n,&
atoms%dx(ntyp), rad, rhoss, vxc)
ELSE
CALL xcpot%get_vxc(input%jspins,rhoss,vxc,vx)
ENDIF
......
......@@ -311,7 +311,7 @@ SUBROUTINE stden(mpi,sphhar,stars,atoms,sym,DIMENSION,vacuum,&
END DO
CALL qsf(vacuum%delz,sigm,vacpar(ivac),vacuum%nmz,0)
denz1 = den%vacz(1,ivac,ispin) ! get estimate for potential at vacuum boundary
CALL xcpot%get_vxc_start(1,denz1,vacpot,vacxpot)
CALL xcpot%get_vxc(1,denz1,vacpot,vacxpot)
! seems to be the best choice for 1D not to substract vacpar
IF (.NOT.oneD%odi%d1) THEN
vacpot = vacpot - fpi_const*vacpar(ivac)
......
<?xml version="1.0" encoding="UTF-8" standalone="no"?>
<fleurInput fleurInputVersion="0.27">
<comment>
Fe fcc 2-atom uc
</comment>
<calculationSetup>
<cutoffs Kmax="3.40000000" Gmax="10.20000000" GmaxXC="8.50000000" numbands="0"/>
<scfLoop itmax="20" maxIterBroyd="99" imix="Anderson" alpha=".05000000" spinf="2.00000000"/>
<coreElectrons ctail="F" frcor="F" kcrel="0"/>
<magnetism jspins="2" l_noco="T" l_J="F" swsp="F" lflip="F"/>
<soc theta=".00000000" phi=".00000000" l_soc="F" spav="F" off="F" soc66="T"/>
<nocoParams l_ss="F" l_mperp="F" l_constr="F" l_disp="F" sso_opt="FFF" mix_b="0.5" thetaJ="0.0" nsh="0">
<qss>0.0 0.0 0.0</qss>
</nocoParams>
<expertModes gw="0" pot8="F" eig66="F" lpr="0" isec1="99" secvar="F"/>
<geometryOptimization l_f="F" xa="2.00000000" thetad="330.00000000" epsdisp=".00001000" epsforce=".00001000"/>
<bzIntegration valenceElectrons="16.00000000" mode="hist" fermiSmearingEnergy=".00100000">
<kPointCount count="5" gamma="F"/>
</bzIntegration>
<energyParameterLimits ellow="-.80000000" elup="1.00000000"/>
</calculationSetup>
<cell>
<symmetryFile filename="sym.out"/>
<bulkLattice scale="1.000000000000" latnam="squ">
<a1>4.82246838</a1>
<c>6.82000000</c>
</bulkLattice>
</cell>
<xcFunctional name="libxc" relativisticCorrections="F">
<libXC exchange="101" correlation="130" />
</xcFunctional>
<atomSpecies>
<species name="Fe-1" element="Fe" atomicNumber="26" coreStates="7" magMom="2.20000000" flipSpin="T">
<mtSphere radius="2.35000000" gridPoints="565" logIncrement=".02100000"/>
<atomicCutoffs lmax="8" lnonsphr="6"/>
<energyParameters s="4" p="4" d="3" f="4"/>
</species>
</atomSpecies>
<atomGroups>
<atomGroup species="Fe-1">
<relPos>0.0 0.0 0.0</relPos>
<force calculate="T" relaxXYZ="TTT"/>
<nocoParams l_relax="F" l_magn="F" M="0.0" alpha="0.0" beta="0.0" b_cons_x="0.0" b_cons_y="0.0"/>
</atomGroup>
<atomGroup species="Fe-1">
<relPos>0.5 0.5 0.5</relPos>
<force calculate="T" relaxXYZ="TTT"/>
<nocoParams l_relax="F" l_magn="F" M="0.0" alpha="0.0" beta="3.1415926536" b_cons_x="0.0" b_cons_y="0.0"/>
</atomGroup>
</atomGroups>
<output dos="F" band="F" vacdos="F" slice="F">
<checks vchk="F" cdinf="F" disp="F"/>
<densityOfStates ndir="0" minEnergy="-.50000000" maxEnergy=".50000000" sigma=".01500000"/>
<vacuumDOS layers="0" integ="F" star="F" nstars="0" locx1=".00000000" locy1=".00000000" locx2=".00000000" locy2=".00000000" nstm="0" tworkf=".00000000"/>
<plotting iplot="F" score="F" plplot="F"/>
<chargeDensitySlicing numkpt="0" minEigenval=".00000000" maxEigenval=".00000000" nnne="0" pallst="F"/>
<specialOutput form66="F" eonly="F" bmt="F"/>
</output>
</fleurInput>
16 16 T ! nop,nop2,symor
! 1
1 0 0 0.00000
0 1 0 0.00000
0 0 1 0.00000
! 2
-1 0 0 0.00000
0 1 0 0.00000
0 0 1 0.00000
! 3
1 0 0 0.00000
0 -1 0 0.00000
0 0 1 0.00000
! 4
-1 0 0 0.00000
0 -1 0 0.00000
0 0 1 0.00000
! 5
0 -1 0 0.00000
-1 0 0 0.00000
0 0 1 0.00000
! 6
0 -1 0 0.00000
1 0 0 0.00000
0 0 1 0.00000
! 7
0 1 0 0.00000
-1 0 0 0.00000
0 0 1 0.00000
! 8
0 1 0 0.00000
1 0 0 0.00000
0 0 1 0.00000
! 9
1 0 0 0.00000
0 1 0 0.00000
0 0 -1 0.00000
! 10
-1 0 0 0.00000
0 1 0 0.00000
0 0 -1 0.00000
! 11
1 0 0 0.00000
0 -1 0 0.00000
0 0 -1 0.00000
! 12
-1 0 0 0.00000
0 -1 0 0.00000
0 0 -1 0.00000
! 13
0 -1 0 0.00000
-1 0 0 0.00000
0 0 -1 0.00000
! 14
0 -1 0 0.00000
1 0 0 0.00000
0 0 -1 0.00000
! 15
0 1 0 0.00000
-1 0 0 0.00000
0 0 -1 0.00000
! 16
0 1 0 0.00000
1 0 0 0.00000
0 0 -1 0.00000
$test_name="Fleur Fe bct non-collinear XML";
$test_code="Fleur";
%test_requirements=("SOC",0,"complex",1);
$test_stages=1;
$test_desc=<<EOF
Simple test of Fleur with two steps:
1.Generate a starting density
2.Run 20 iterations and compare convergence, fermi-energy & total energy
EOF
;
#juDFT Testscript
jt::copyfile("files/inp.xml",$workdir);
jt::copyfile("files/sym.out",$workdir);
jt::testrun("$executable -xmlInput",$workdir);
#now test output
$result=jt::test_fileexists("$workdir/out");
$result+=jt::test_fileexists("$workdir/cdn1");
$result+=jt::test_grepexists("$workdir/out","total charge");
$result+=jt::test_grepnumber("$workdir/out","qfix","qfix= *([^ ]*)",1.0,0.00001);
$result=jt::test_grepexists("$workdir/out","it= 20 is completed");
$result+=jt::test_grepnumber("$workdir/out","new fermi energy",".*: *([^ ]*)",0.326,0.005);
$result+=jt::test_grepnumber("$workdir/out","total energy=",".*= *([^ ]*)",-2545.607,0.005);
$result+=jt::test_grepnumber("$workdir/out","distance of charge densities for it= 20",": *([^ ]*)",0.0000,0.09);
$result+=jt::test_grepnumber("$workdir/out","mm 2",".*mm 2 *([^ ]*)",1.73,0.03);
jt::stageresult($workdir,$result,"1");
......@@ -32,8 +32,8 @@ endif()
#Tests for LibXC
if (${FLEUR_USE_LIBXC})
set(Testdirs ${Testdirs} CuBulkLibXC)
set(ParTestdirs ${ParTestdirs} CuBulkLibXC)
set(Testdirs ${Testdirs} CuBulkLibXC Fe_bct_LibXC)
set(ParTestdirs ${ParTestdirs} CuBulkLibXC Fe_bct_LibXC)
endif()
#The serial tests
if (${FLEUR_USE_SERIAL})
......
......@@ -25,7 +25,6 @@ MODULE m_types_xcpot
PROCEDURE :: is_hybrid=>xcpot_is_hybrid
PROCEDURE :: get_exchange_weight=>xcpot_get_exchange_weight
PROCEDURE :: get_vxc=>xcpot_get_vxc
PROCEDURE :: get_vxc_start=>xcpot_get_vxc_start
PROCEDURE :: get_exc=>xcpot_get_exc
PROCEDURE,NOPASS :: alloc_gradients=>xcpot_alloc_gradients
END TYPE t_xcpot
......@@ -93,17 +92,6 @@ CONTAINS
TYPE(t_gradients),OPTIONAL,INTENT(INOUT)::grad
END SUBROUTINE xcpot_get_vxc
SUBROUTINE xcpot_get_vxc_start(xcpot,jspins,rh,vxc,vx,grad)
CLASS(t_xcpot),INTENT(IN) :: xcpot
INTEGER, INTENT (IN) :: jspins
!--> charge density
REAL,INTENT (IN) :: rh(:,:)
!---> xc potential
REAL, INTENT (OUT) :: vxc (:,:),vx(:,:)
TYPE(t_gradients),OPTIONAL,INTENT(INOUT)::grad
END SUBROUTINE xcpot_get_vxc_start
SUBROUTINE xcpot_get_exc(xcpot,jspins,rh,exc,grad)
CLASS(t_xcpot),INTENT(IN) :: xcpot
INTEGER, INTENT (IN) :: jspins
......
......@@ -47,7 +47,6 @@ MODULE m_types_xcpot_inbuild
PROCEDURE :: is_hybrid=>xcpot_is_hybrid
PROCEDURE :: get_exchange_weight=>xcpot_get_exchange_weight
PROCEDURE :: get_vxc=>xcpot_get_vxc
PROCEDURE :: get_vxc_start=>xcpot_get_vxc_start
PROCEDURE :: get_exc=>xcpot_get_exc
!not overloaded
PROCEDURE :: get_name=>xcpot_get_name
......@@ -134,19 +133,6 @@ CONTAINS
END FUNCTION xcpot_get_exchange_weight
SUBROUTINE xcpot_get_vxc_start(xcpot,jspins,rh, vxc,vx, grad)
CLASS(t_xcpot_inbuild),INTENT(IN) :: xcpot
INTEGER, INTENT (IN) :: jspins
REAL,INTENT (IN) :: rh(:,:)
REAL, INTENT (OUT) :: vx (:,:)
REAL, INTENT (OUT) :: vxc(:,:)
TYPE(t_gradients), INTENT(INOUT), OPTIONAL :: grad
! there is no difference for inbuild case
call xcpot%get_vxc(jspins, rh, vxc, vx, grad)
END SUBROUTINE xcpot_get_vxc_start
SUBROUTINE xcpot_get_vxc(xcpot,jspins,rh, vxc,vx, grad)
!
USE m_xcxal, ONLY : vxcxal
......
......@@ -27,7 +27,6 @@ MODULE m_types_xcpot_libxc
PROCEDURE :: is_hybrid=>xcpot_is_hybrid
PROCEDURE :: get_exchange_weight=>xcpot_get_exchange_weight
PROCEDURE :: get_vxc=>xcpot_get_vxc
PROCEDURE :: get_vxc_start=>xcpot_get_vxc_start
PROCEDURE :: get_exc=>xcpot_get_exc
PROCEDURE,NOPASS :: alloc_gradients=>xcpot_alloc_gradients
!Not overloeaded...
......@@ -129,23 +128,6 @@ END FUNCTION xcpot_is_LDA
#endif
END FUNCTION xcpot_get_exchange_weight
SUBROUTINE xcpot_get_vxc_start(xcpot,jspins,rh, vxc,vx, grad)
IMPLICIT NONE
CLASS(t_xcpot_libxc),INTENT(IN) :: xcpot
INTEGER, INTENT (IN) :: jspins
REAL,INTENT (IN) :: rh(:,:) !Dimensions here
REAL, INTENT (OUT) :: vx (:,:) !points,spin
REAL, INTENT (OUT ) :: vxc(:,:) !
! optional arguments for GGA
TYPE(t_gradients),OPTIONAL,INTENT(INOUT)::grad
IF (xcpot%is_lda()) THEN
call xcpot%get_vxc(jspins,rh, vxc,vx, grad)
ELSE
ENDIF
END SUBROUTINE xcpot_get_vxc_start
!***********************************************************************
SUBROUTINE xcpot_get_vxc(xcpot,jspins,rh, vxc,vx, grad)
!***********************************************************************
......
set(fleur_F90 ${fleur_F90}
xc-pot/libxc_postprocess_gga.f90
xc-pot/potl0.f90
)
set(fleur_F77 ${fleur_F77}
......@@ -16,7 +17,6 @@ xc-pot/gaunt.f
xc-pot/grdchlh.f
xc-pot/mkgl0.f
xc-pot/pbecor2.f
xc-pot/potl0.f
xc-pot/relcor.f
xc-pot/vxcepbe.f
xc-pot/vxcl91.f
......
MODULE m_potl0
c ******************************************************************
c evaluate the xc-potential vxc for charge density and its
c gradients,dens,... only for nonmagnetic.
c ******************************************************************
CONTAINS
SUBROUTINE potl0(
> xcpot,mshd,jspd,jspins,msh,dx,rad,dens,
< vxc)
MODULE m_potl0
! ******************************************************************
! evaluate the xc-potential vxc for charge density and its
! gradients,dens,... only for nonmagnetic.
! ******************************************************************
CONTAINS
SUBROUTINE potl0(xcpot,mshd,jspd,jspins,msh,dx,rad,dens, &
vxc)
USE m_grdchlh
USE m_mkgl0
USE m_types
IMPLICIT NONE
c ..
CLASS(t_xcpot),intent(in)::xcpot
INTEGER, INTENT (IN) :: jspins,jspd,mshd,msh
REAL, INTENT (IN) :: dx
REAL, INTENT (IN) :: rad(msh),dens(mshd,jspd)
REAL, INTENT (OUT):: vxc(mshd,jspd)
c .. previously untyped names ..
! .. previously untyped names ..
INTEGER,PARAMETER :: ndvgrd=6
TYPE(t_gradients)::grad
INTEGER i,ispin
REAL, ALLOCATABLE :: drr(:,:),ddrr(:,:)
REAL :: vx(mshd,jspd)
ALLOCATE ( drr(mshd,jspd),ddrr(mshd,jspd))
ALLOCATE ( drr(mshd,jspd),ddrr(mshd,jspd))
!
!--> evaluate gradients of dens.
!
CALL xcpot%alloc_gradients(msh,1,grad)
CALL xcpot%alloc_gradients(msh,jspins,grad)
DO ispin = 1, jspins
CALL grdchlh(
> 1,1,msh,dx,rad,dens(1,ispin),ndvgrd,
< drr(1,ispin),ddrr(1,ispin))
CALL grdchlh(1,1,msh,dx,rad,dens(1,ispin),ndvgrd,&
drr(1,ispin),ddrr(1,ispin))
ENDDO
CALL mkgl0(
> mshd,msh,jspd,jspins,rad,dens,drr,ddrr,
< grad)
CALL mkgl0(mshd,msh,jspd,jspins,rad,dens,drr,ddrr,&
grad)
!
! --> calculate the potential.
!
call xcpot%get_vxc(jspins,dens,vxc,vx,grad)
call xcpot%get_vxc(jspins, dens, vxc, vx, grad)
DEALLOCATE ( drr,ddrr )
END SUBROUTINE potl0
END MODULE m_potl0
END SUBROUTINE potl0
END MODULE m_potl0
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