LCPQ
/
quantum_package
mirror of https://github.com/LCPQ/quantum_package
10
0
Fork 0
Code Issues Releases Wiki Activity

qdpt test

This commit is contained in:
Anthony Scemama 2017-11-28 17:11:44 +01:00
parent d680c85b8c
commit ff8fcfa42b
2 changed files with 108 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
plugins/Full_CI/H_apply.irp.f
View File

@ -4,11 +4,13 @@ from generate_h_apply import *
s = H_apply("FCI")
s.set_selection_pt2("epstein_nesbet_2x2")
#s.set_selection_pt2("qdpt")
s.unset_skip()
print s
s = H_apply("FCI_PT2")
s.set_perturbation("epstein_nesbet_2x2")
#s.set_perturbation("qdpt")
s.unset_skip()
s.unset_openmp()
print s

106
plugins/Perturbation/pt2_equations.irp.f
View File

@ -45,6 +45,58 @@ subroutine pt2_epstein_nesbet ($arguments)
end
subroutine pt2_qdpt ($arguments)
use bitmasks
implicit none
$declarations
BEGIN_DOC
! compute the QDPT first order coefficient and second order energetic contribution
!
! for the various N_st states.
!
! c_pert(i) = <psi(i)|H|det_pert>/( <psi(i)|H|psi(i)> - <det_pert|H|det_pert> )
!
END_DOC
integer :: i,j
double precision :: diag_H_mat_elem_fock, h, E, diag_H_mat_elem, hij
double precision :: i_H_psi_array(N_st)
integer :: degree
double precision :: delta_E
PROVIDE selection_criterion
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
!call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
call i_H_psi_minilist(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array)
h = diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
c_pert = 0.d0
do j=1,N_det_selectors
call get_excitation_degree(det_ref, psi_selectors(1,1,j), degree, Nint)
if (degree > 2) then
E = diag_H_mat_elem(psi_selectors(1,1,j),Nint)
else
E = diag_H_mat_elem_fock(det_ref,det_ref,fock_diag_tmp,Nint)
endif
delta_E = E-h
! delta_E = electronic_energy(1) - h
call i_H_j(psi_selectors(1,1,j),det_pert,Nint,hij)
if (dabs(delta_e) > 1.d-3) then
do i =1,N_st
c_pert(i) += psi_selectors_coef(j,i) * hij / delta_e
enddo
endif
enddo
do i =1,N_st
e_2_pert(i) = c_pert(i)*i_H_psi_array(i)
H_pert_diag(i) = h*c_pert(i)*c_pert(i)
enddo
end
subroutine pt2_decontracted ($arguments)
use bitmasks
@ -249,6 +301,60 @@ subroutine pt2_moller_plesset ($arguments)
end
subroutine pt2_moller_plesset_general ($arguments)
use bitmasks
implicit none
$declarations
BEGIN_DOC
! compute the general Moller-Plesset perturbative first order coefficient and second order energetic contribution
!
! for the various n_st states.
!
! c_pert(i) = <psi(i)|H|det_pert>/(difference of orbital energies)
!
! e_2_pert(i) = <psi(i)|H|det_pert>^2/(difference of orbital energies)
!
END_DOC
integer :: i,j
double precision :: diag_H_mat_elem_fock
integer :: exc(0:2,2,2)
integer :: degree
double precision :: phase,delta_e,h
double precision :: i_H_psi_array(N_st)
integer :: h1,h2,p1,p2,s1,s2
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
call get_excitation(det_ref,det_pert,exc,degree,phase,Nint)
if (degree == 2) then
call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)
delta_e = (mo_energy_expval(1,h1,s1,1) - mo_energy_expval(1,p1,s1,2)) + &
(mo_energy_expval(1,h2,s2,1) - mo_energy_expval(1,p2,s2,2))
else if (degree == 1) then
call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)
delta_e = mo_energy_expval(1,h1,s1,1) - mo_energy_expval(1,p1,s1,2)
else
delta_e = 0.d0
endif
if (dabs(delta_e) > 1.d-10) then
delta_e = 1.d0/delta_e
call i_H_psi_minilist(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array)
h = diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
else
i_H_psi_array(:) = 0.d0
h = 0.d0
endif
do i =1,N_st
H_pert_diag(i) = h
c_pert(i) = i_H_psi_array(i) *delta_e
e_2_pert(i) = c_pert(i) * i_H_psi_array(i)
enddo
end
subroutine pt2_epstein_nesbet_SC2_projected ($arguments)
use bitmasks