Commit 5ef716dd authored by Matthias Redies's avatar Matthias Redies

TB09 converges, but lattice wrong

parent 53438780
......@@ -203,7 +203,7 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
perform_MetaGGA = ALLOCATED(EnergyDen%mt) &
.AND. (xcpot%exc_is_MetaGGA() .or. xcpot%vx_is_MetaGGA())
if(perform_MetaGGA) then
call set_kinED(mpi, sphhar, atoms, core_den, val_den, xcpot, &
call set_kinED(mpi, sphhar, atoms, sym, core_den, val_den, xcpot, &
input, noco, stars, cell, outDen, EnergyDen, vTot)
endif
#ifdef CPP_MPI
......
......@@ -280,12 +280,7 @@ CONTAINS
ALLOCATE(vxc_tmp(SIZE(vxc,2),SIZE(vxc,1)));vxc_tmp=0.0
ALLOCATE(vx_tmp(SIZE(vx,2),SIZE(vx,1)));vx_tmp=0.0
write (*,*) "vgen call"
if(present(kinED_KS)) then
write (*,*) "max(abs(kED)", maxval(abs(kinED_KS))
endif
IF (xcpot%vc_is_GGA()) THEN
write (*,*) "somehow in GGA branch"
IF (.NOT.PRESENT(grad)) CALL judft_error("Bug: You called get_vxc for a GGA potential without providing derivatives")
ALLOCATE(vsigma,mold=grad%vsigma)
!where(abs(grad%sigma)<1E-9) grad%sigma=1E-9
......@@ -298,8 +293,6 @@ CONTAINS
vxc_tmp=vx_tmp
ENDIF
ELSE !LDA potentials
write (*,*) "Is MGGA = ", xcpot%vx_is_MetaGGA()
write (*,*) "kinED_KS present = ", present(kinED_KS)
if(xcpot%vx_is_MetaGGA() .and. present(kinED_KS)) then
kinED_libXC = transpose(kinED_KS + 0.25 * grad%laplace)
......@@ -307,9 +300,6 @@ CONTAINS
allocate(tmp_vlapl, mold=vx_tmp)
allocate(tmp_vtau, mold=vx_tmp)
write (*,*) "using = ", trim(xc_f03_func_info_get_name(&
xc_f03_func_get_info(xcpot%vxc_func_x)))
call xc_f03_mgga_vxc(xcpot%vxc_func_x, size(rh,1), transpose(rh), &
grad%sigma, transpose(grad%laplace), kinED_libxc,&
vx_tmp, tmp_vsig, tmp_vlapl, tmp_vtau)
......@@ -317,15 +307,11 @@ CONTAINS
idx = find_first_normal(vx_tmp)+1
vx_tmp(:,:idx) = 0.0
CALL xc_f03_lda_vxc(initial_lda_func(jspins), idx, TRANSPOSE(rh(:idx,:)), vx_tmp(:,:idx))
write (*,*) "vx has nan = ", any(ieee_is_nan(vx_tmp))
else
write (*,*) "using = ", trim(xc_f03_func_info_get_name(&
xc_f03_func_get_info(initial_lda_func(jspins))))
CALL xc_f03_lda_vxc(initial_lda_func(jspins), SIZE(rh,1), TRANSPOSE(rh), vx_tmp)
endif
IF (xcpot%func_vxc_id_c>0) THEN
CALL xc_f03_lda_vxc(xcpot%vxc_func_c, SIZE(rh,1), TRANSPOSE(rh), vxc_tmp)
write (*,*) "added correlation"
vxc_tmp=vxc_tmp+vx_tmp
ENDIF
ENDIF
......@@ -344,6 +330,7 @@ CONTAINS
idx = size(vec, dim=2)
do while(all(.not. ieee_is_nan(vec(:,idx))))
idx = idx - 1
if(idx == 0) exit
enddo
END FUNCTION find_first_normal
......
This diff is collapsed.
......@@ -92,6 +92,5 @@ CONTAINS
ENDIF
CALL finish_pw_grid()
END SUBROUTINE vis_xc
END MODULE m_vis_xc
......@@ -51,8 +51,7 @@
TYPE(t_xcpot_inbuild) :: xcpot_tmp
TYPE(t_potden) :: vTot_tmp
TYPE(t_sphhar) :: tmp_sphhar
REAL, ALLOCATABLE :: ch(:,:), core_den_rs(:,:), val_den_rs(:,:), ED_rs(:,:), &
vTot_rs(:,:), vTot0_rs(:,:)
REAL, ALLOCATABLE :: ch(:,:)
INTEGER :: n,nsp,nt,jr, loc_n
INTEGER :: i, j, idx, cnt
REAL :: divi
......
......@@ -180,13 +180,14 @@ CONTAINS
res = matmul(transpose(cell%bmat), vec)
end function internal_to_rez
subroutine set_kinED(mpi, sphhar, atoms, core_den, val_den, xcpot, &
subroutine set_kinED(mpi, sphhar, atoms, sym, core_den, val_den, xcpot, &
input, noco, stars, cell, den, EnergyDen, vTot)
use m_types
implicit none
TYPE(t_mpi),INTENT(IN) :: mpi
TYPE(t_sphhar),INTENT(IN) :: sphhar
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_sym), INTENT(IN) :: sym
TYPE(t_potden),INTENT(IN) :: core_den, val_den
CLASS(t_xcpot),INTENT(INOUT) :: xcpot
TYPE(t_input),INTENT(IN) :: input
......@@ -195,12 +196,12 @@ CONTAINS
TYPE(t_cell),INTENT(IN) :: cell
TYPE(t_potden),INTENT(IN) :: den, EnergyDen, vTot
call set_kinED_is(xcpot, input, noco, stars, cell, den, EnergyDen, vTot)
call set_kinED_mt(mpi, sphhar, atoms, core_den, val_den, &
call set_kinED_is(xcpot, input, noco, stars, sym, cell, den, EnergyDen, vTot)
call set_kinED_mt(mpi, sphhar, atoms, sym, core_den, val_den, &
xcpot, EnergyDen, input, vTot)
end subroutine set_kinED
subroutine set_kinED_is(xcpot, input, noco, stars, cell, den, EnergyDen, vTot)
subroutine set_kinED_is(xcpot, input, noco, stars, sym, cell, den, EnergyDen, vTot)
use m_types
use m_pw_tofrom_grid
implicit none
......@@ -208,12 +209,15 @@ CONTAINS
TYPE(t_input),INTENT(IN) :: input
TYPE(t_noco),INTENT(IN) :: noco
TYPE(t_stars),INTENT(IN) :: stars
TYPE(t_sym), INTENT(IN) :: sym
TYPE(t_cell),INTENT(IN) :: cell
TYPE(t_potden),INTENT(IN) :: den, EnergyDen, vTot
!local arrays
REAL, ALLOCATABLE :: den_rs(:,:), ED_rs(:,:), vTot_rs(:,:)
TYPE(t_gradients) :: tmp_grad
CALL init_pw_grid(xcpot,stars,sym,cell)
CALL pw_to_grid(xcpot, input%jspins, noco%l_noco, stars, &
cell, EnergyDen%pw, tmp_grad, ED_rs)
......@@ -222,11 +226,13 @@ CONTAINS
CALL pw_to_grid(xcpot, input%jspins, noco%l_noco, stars, &
cell, den%pw, tmp_grad, den_rs)
CALL finish_pw_grid()
xcpot%kinED%is = ED_RS - vTot_RS * den_RS
xcpot%kinED%set = .True.
end subroutine set_kinED_is
subroutine set_kinED_mt(mpi, sphhar, atoms, core_den, val_den, &
subroutine set_kinED_mt(mpi, sphhar, atoms, sym, core_den, val_den, &
xcpot, EnergyDen, input, vTot)
use m_types
use m_mt_tofrom_grid
......@@ -234,6 +240,7 @@ CONTAINS
TYPE(t_mpi),INTENT(IN) :: mpi
TYPE(t_sphhar),INTENT(IN) :: sphhar
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_sym), INTENT(IN) :: sym
TYPE(t_potden),INTENT(IN) :: core_den, val_den, EnergyDen, vTot
CLASS(t_xcpot),INTENT(INOUT) :: xcpot
TYPE(t_input),INTENT(IN) :: input
......@@ -251,7 +258,14 @@ CONTAINS
n_start=1
n_stride=1
#endif
CALL init_mt_grid(input%jspins,atoms,sphhar,xcpot,sym)
loc_n = 0
allocate(ED_rs(atoms%nsp()*atoms%jmtd, input%jspins))
allocate(vTot_rs, mold=ED_rs)
allocate(vTot0_rs, mold=ED_rs)
allocate(core_den_rs, mold=ED_rs)
allocate(val_den_rs, mold=ED_rs)
call xcpot%kinED%alloc_mt(atoms%nsp()*atoms%jmtd, input%jspins, &
n_start, atoms%ntype, n_stride)
loc_n = 0
......@@ -262,10 +276,6 @@ CONTAINS
allocate(vTot_mt(lbound(vTot%mt, dim=1):ubound(vTot%mt, dim=1),&
lbound(vTot%mt, dim=2):ubound(vTot%mt, dim=2),&
lbound(vTot%mt, dim=4):ubound(vTot%mt, dim=4)))
write (*,*) "lbound vTot_mt = ", lbound(vTot_mt)
write (*,*) "ubound vTot_mt = ", ubound(vTot_mt)
write (*,*) "lbound vTot%mt = ", lbound(vTot%mt)
write (*,*) "ubound vTot%mt = ", ubound(vTot%mt)
endif
do jr=1,atoms%jri(n)
......@@ -289,5 +299,6 @@ CONTAINS
xcpot%kinED%mt(:,:,loc_n) = ED_RS - (vTot0_rs * core_den_rs + vTot_rs * val_den_rs)
enddo
xcpot%kinED%set = .True.
CALL finish_mt_grid()
end subroutine set_kinED_mt
END MODULE m_metagga
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