Skip to content

Commit

Permalink
GAPW triplet excitation energies
Browse files Browse the repository at this point in the history
  • Loading branch information
@oschuett
Anna Hehn authored and oschuett committed Jun 26, 2023
1 parent 71da3ea commit bd8365e
Show file tree
Hide file tree
Showing 6 changed files with 162 additions and 9 deletions.
5 changes: 3 additions & 2 deletions src/qs_vxc_atom.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -546,6 +546,7 @@ SUBROUTINE calculate_xc_2nd_deriv_atom(rho_atom_set, rho1_atom_set, qs_env, xc_s
CALL section_vals_val_get(input, "DFT%TDDFPT%RES_ETYPE", &
i_val=res_etype)
END IF

xc_fun_section => section_vals_get_subs_vals(xc_section, &
"XC_FUNCTIONAL")
IF (lsd) THEN
Expand Down Expand Up @@ -710,11 +711,11 @@ SUBROUTINE calculate_xc_2nd_deriv_atom(rho_atom_set, rho1_atom_set, qs_env, xc_s
CALL xc_2nd_deriv_of_r(xc_section=xc_section, &
rho_set=rho_set_h, rho1_set=rho1_set_h, &
deriv_set=deriv_set, &
w=weight, vxc=vxc_h, vxg=vxg_h)
w=weight, vxc=vxc_h, vxg=vxg_h, do_triplet=do_triplet)
CALL xc_2nd_deriv_of_r(xc_section=xc_section, &
rho_set=rho_set_s, rho1_set=rho1_set_s, &
deriv_set=deriv_set, &
w=weight, vxc=vxc_s, vxg=vxg_s)
w=weight, vxc=vxc_s, vxg=vxg_s, do_triplet=do_triplet)

CALL get_rho_atom(rho_atom=rho1_atom, ga_Vlocal_gb_h=int_hh, ga_Vlocal_gb_s=int_ss)
IF (gradient_functional) THEN
Expand Down
10 changes: 6 additions & 4 deletions src/xc/xc.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2623,7 +2623,7 @@ SUBROUTINE xc_calc_2nd_deriv_analytical(v_xc, v_xc_tau, deriv_set, rho_set, rho1
INTEGER, DIMENSION(2, 3) :: bo
LOGICAL :: gradient_f, lsd, my_compute_virial, &
my_gapw, tau_f, laplace_f, rho_f
REAL(KIND=dp) :: fac, gradient_cut, tmp
REAL(KIND=dp) :: fac, gradient_cut, tmp, factor2
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :) :: dr1dr, dra1dra, drb1drb
REAL(kind=dp), DIMENSION(:, :, :), POINTER :: deriv_data, e_drhoa, e_drhob, &
e_drho, norm_drho, norm_drhoa, &
Expand Down Expand Up @@ -2661,7 +2661,9 @@ SUBROUTINE xc_calc_2nd_deriv_analytical(v_xc, v_xc_tau, deriv_set, rho_set, rho1
nspins = SIZE(v_xc)
lsd = ASSOCIATED(rho_set%rhoa)
fac = 0.0_dp
factor2 = 1.0_dp
IF (PRESENT(tddfpt_fac)) fac = tddfpt_fac
IF (PRESENT(tddfpt_fac)) factor2 = tddfpt_fac
bo = rho_set%local_bounds
Expand Down Expand Up @@ -2938,12 +2940,12 @@ SUBROUTINE xc_calc_2nd_deriv_analytical(v_xc, v_xc_tau, deriv_set, rho_set, rho1
DO idir = 1, 3
!$OMP PARALLEL DO PRIVATE(ia,ir) DEFAULT(NONE) &
!$OMP SHARED(bo,vxg,drho,v_drho,e_drho,drho1,idir) COLLAPSE(2)
!$OMP SHARED(bo,vxg,drho,v_drho,e_drho,drho1,idir,factor2) COLLAPSE(2)
DO ia = bo(1, 1), bo(2, 1)
DO ir = bo(1, 2), bo(2, 2)
vxg(idir, ia, ir, 1) = -drho(idir)%array(ia, ir, 1)*v_drho(1)%cr3d(ia, ir, 1)
IF (ASSOCIATED(e_drho)) THEN
vxg(idir, ia, ir, 1) = vxg(idir, ia, ir, 1) + drho1(idir)%array(ia, ir, 1)*e_drho(ia, ir, 1)
vxg(idir, ia, ir, 1) = vxg(idir, ia, ir, 1) + factor2*drho1(idir)%array(ia, ir, 1)*e_drho(ia, ir, 1)
END IF
END DO
END DO
Expand Down Expand Up @@ -3401,7 +3403,7 @@ END SUBROUTINE prepare_dr1dr
SUBROUTINE prepare_dr1dr_ab(dr1dr, drhoa, drhob, drho1a, drho1b, fac)
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :), &
INTENT(OUT) :: dr1dr
TYPE(cp_3d_r_cp_type), DIMENSION(3), INTENT(IN) :: drhoa, drhob, drho1a, drho1b
TYPE(cp_3d_r_cp_type), DIMENSION(3), INTENT(IN) :: drhoa, drhob, drho1a, drho1b
REAL(KIND=dp), INTENT(IN) :: fac
CHARACTER(len=*), PARAMETER :: routineN = 'prepare_dr1dr_ab'
Expand Down
11 changes: 8 additions & 3 deletions src/xc/xc_atom.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -337,9 +337,10 @@ END SUBROUTINE vxc_of_r_new
!> \param w ...
!> \param vxc ...
!> \param vxg ...
!> \param do_triplet ...
! **************************************************************************************************
SUBROUTINE xc_2nd_deriv_of_r(rho_set, rho1_set, xc_section, &
deriv_set, w, vxc, vxg)
deriv_set, w, vxc, vxg, do_triplet)
! As input of this routine one gets rho and drho on a one dimensional grid.
! The grid is the angular grid corresponding to a given point ir on the radial grid.
Expand All @@ -354,12 +355,13 @@ SUBROUTINE xc_2nd_deriv_of_r(rho_set, rho1_set, xc_section, &
REAL(dp), DIMENSION(:, :), POINTER :: w
REAL(dp), CONTIGUOUS, DIMENSION(:, :, :), POINTER :: vxc
REAL(dp), DIMENSION(:, :, :, :), POINTER :: vxg
LOGICAL, INTENT(IN), OPTIONAL :: do_triplet
CHARACTER(LEN=*), PARAMETER :: routineN = 'xc_2nd_deriv_of_r'
INTEGER :: handle, ispin, nspins
LOGICAL :: lsd
REAL(dp) :: drho_cutoff
REAL(dp) :: drho_cutoff, my_fac_triplet
TYPE(cp_sll_xc_deriv_type), POINTER :: pos
TYPE(pw_pool_type), POINTER :: pw_pool
TYPE(pw_type), DIMENSION(:), POINTER :: vxc_pw, vxc_tau_pw
Expand All @@ -373,6 +375,9 @@ SUBROUTINE xc_2nd_deriv_of_r(rho_set, rho1_set, xc_section, &
IF (ASSOCIATED(rho_set%rhoa)) THEN
lsd = .TRUE.
END IF
my_fac_triplet = 1.0_dp
IF ((PRESENT(do_triplet)) .AND. (do_triplet)) my_fac_triplet = -1.0_dp
CALL xc_rho_set_get(rho_set, drho_cutoff=drho_cutoff)
xc_fun_section => section_vals_get_subs_vals(xc_section, &
"XC_FUNCTIONAL")
Expand Down Expand Up @@ -402,7 +407,7 @@ SUBROUTINE xc_2nd_deriv_of_r(rho_set, rho1_set, xc_section, &
NULLIFY (vxc_tau_pw)
CALL xc_calc_2nd_deriv_analytical(vxc_pw, vxc_tau_pw, deriv_set, rho_set, rho1_set, pw_pool, &
xc_section, gapw=.TRUE., vxg=vxg)
xc_section, gapw=.TRUE., vxg=vxg, tddfpt_fac=my_fac_triplet)
DEALLOCATE (vxc_pw)
Expand Down
67 changes: 67 additions & 0 deletions tests/QS/regtest-tddfpt/H2O_GAPW_1_triplet.inp
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,67 @@
&GLOBAL
PROJECT H2O_GAPW
RUN_TYPE ENERGY
PRINT_LEVEL LOW
&END GLOBAL
&FORCE_EVAL
METHOD Quickstep
&PROPERTIES
&TDDFPT
NSTATES 3
MAX_ITER 10
MAX_KV 10
CONVERGENCE 1.0e-5
&XC
&XC_FUNCTIONAL PBE
&END XC_FUNCTIONAL
&END XC
RKS_TRIPLETS T
&END TDDFPT
&END PROPERTIES

&DFT
BASIS_SET_FILE_NAME EMSL_BASIS_SETS
POTENTIAL_FILE_NAME POTENTIAL
&MGRID
CUTOFF 200
&END MGRID
&QS
METHOD GAPW
&END QS
&SCF
MAX_SCF 40
SCF_GUESS ATOMIC
&END SCF
&XC
&XC_FUNCTIONAL PBE
&END XC_FUNCTIONAL
&END XC
&POISSON
PERIODIC NONE
POISSON_SOLVER WAVELET
&END
&END DFT
&SUBSYS
&CELL
ABC 6.0 6.0 6.0
PERIODIC NONE
&END CELL
&COORD
O 0.000000 0.000000 -0.065587 H2O
H 0.000000 -0.757136 0.520545 H2O
H 0.000000 0.757136 0.520545 H2O
&END COORD
&TOPOLOGY
&CENTER_COORDINATES
&END
&END
&KIND H
BASIS_SET ORB 6-311Gxx
POTENTIAL ALL
&END KIND
&KIND O
BASIS_SET ORB 6-311Gxx
POTENTIAL ALL
&END KIND
&END SUBSYS
&END FORCE_EVAL
76 changes: 76 additions & 0 deletions tests/QS/regtest-tddfpt/Ne_GAPW_triplet.inp
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,76 @@
&GLOBAL
PROJECT Ne
PRINT_LEVEL LOW
RUN_TYPE ENERGY
&END GLOBAL
&FORCE_EVAL
&PROPERTIES
&TDDFPT
&DIPOLE_MOMENTS
DIPOLE_FORM LENGTH
&END DIPOLE_MOMENTS
NSTATES 1
MAX_ITER 10
MAX_KV 60
CONVERGENCE 1.0e-5
&XC
&XC_FUNCTIONAL PBE
&END XC_FUNCTIONAL
2ND_DERIV_ANALYTICAL T
&END XC
RKS_TRIPLETS T
&END TDDFPT
&END PROPERTIES
METHOD Quickstep
&DFT
BASIS_SET_FILE_NAME EMSL_BASIS_SETS
BASIS_SET_FILE_NAME basis
POTENTIAL_FILE_NAME POTENTIAL
&QS
METHOD GAPW
EPS_DEFAULT 1.0E-17
EPS_PGF_ORB 1.0E-20
&END QS
&MGRID
CUTOFF 200
REL_CUTOFF 80
&END MGRID
&SCF
SCF_GUESS RESTART
&OT
PRECONDITIONER FULL_SINGLE_INVERSE
MINIMIZER DIIS
STEPSIZE 0.1
&END
&OUTER_SCF
MAX_SCF 20
EPS_SCF 1.0E-7
&END
MAX_SCF 20
EPS_SCF 1.0E-7
&END SCF
&POISSON
PERIODIC NONE
PSOLVER WAVELET
&END
&XC
&XC_FUNCTIONAL PBE
&END XC_FUNCTIONAL
2ND_DERIV_ANALYTICAL T
&END XC
&END DFT
&SUBSYS
&CELL
ABC 6.0 6.0 6.0
PERIODIC NONE
&END CELL
&COORD
Ne 2.5 2.5 2.5
&END COORD
&KIND Ne
BASIS_SET Ahlrichs-def2-SVP
POTENTIAL ALL
&END KIND
&END SUBSYS
&END FORCE_EVAL

2 changes: 2 additions & 0 deletions tests/QS/regtest-tddfpt/TEST_FILES
View file
Original file line number Diff line number Diff line change
Expand Up @@ -24,4 +24,6 @@ H2O_GAPW_4.inp 37 4.0E-06
H2O_GAPW_XC_1.inp 37 4.0E-06 0.619451E+00
H2O_GAPW_XC_2.inp 37 4.0E-06 0.812642E+00
H2O_GAPW_XC_3.inp 37 4.0E-06 0.836577E+00
H2O_GAPW_1_triplet.inp 37 4.0E-06 0.505991E+00
Ne_GAPW_triplet.inp 37 4.0E-06 0.168075E+01
#EOF

0 comments on commit bd8365e

Please sign in to comment.