From 2a5cb23fd9f4029fefac377c6a8da26fe81e1e21 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Tue, 24 Nov 2015 11:16:14 +0100 Subject: [PATCH] Added template for pt2_equations --- plugins/Perturbation/Moller_plesset.irp.f | 50 ---------- ...stein_nesbet.irp.f => pt2_equations.irp.f} | 94 +++++++++++++++---- 2 files changed, 76 insertions(+), 68 deletions(-) delete mode 100644 plugins/Perturbation/Moller_plesset.irp.f rename plugins/Perturbation/{epstein_nesbet.irp.f => pt2_equations.irp.f} (55%) diff --git a/plugins/Perturbation/Moller_plesset.irp.f b/plugins/Perturbation/Moller_plesset.irp.f deleted file mode 100644 index 38fe6610..00000000 --- a/plugins/Perturbation/Moller_plesset.irp.f +++ /dev/null @@ -1,50 +0,0 @@ -subroutine pt2_moller_plesset(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_st) - use bitmasks - implicit none - integer, intent(in) :: Nint,ndet,n_st - integer(bit_kind), intent(in) :: det_pert(Nint,2) - double precision , intent(out) :: c_pert(n_st),e_2_pert(n_st),H_pert_diag(N_st) - double precision :: i_H_psi_array(N_st) - - BEGIN_DOC - ! compute the standard Moller-Plesset perturbative first order coefficient and second order energetic contribution - ! - ! for the various n_st states. - ! - ! c_pert(i) = /(difference of orbital energies) - ! - ! e_2_pert(i) = ^2/(difference of orbital energies) - ! - END_DOC - - integer :: i,j - double precision :: diag_H_mat_elem - integer :: exc(0:2,2,2) - integer :: degree - double precision :: phase,delta_e,h - integer :: h1,h2,p1,p2,s1,s2 - ASSERT (Nint == N_int) - ASSERT (Nint > 0) - call get_excitation(ref_bitmask,det_pert,exc,degree,phase,Nint) - if (degree == 2) then - call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2) - delta_e = Fock_matrix_diag_mo(h1) + Fock_matrix_diag_mo(h2) - & - (Fock_matrix_diag_mo(p1) + Fock_matrix_diag_mo(p2)) - delta_e = 1.d0/delta_e - else if (degree == 1) then - call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2) - delta_e = Fock_matrix_diag_mo(h1) - Fock_matrix_diag_mo(p1) - delta_e = 1.d0/delta_e - else - delta_e = 0.d0 - endif - - call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det,psi_selectors_size,n_st,i_H_psi_array) - h = diag_H_mat_elem(det_pert,Nint) - 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 diff --git a/plugins/Perturbation/epstein_nesbet.irp.f b/plugins/Perturbation/pt2_equations.irp.f similarity index 55% rename from plugins/Perturbation/epstein_nesbet.irp.f rename to plugins/Perturbation/pt2_equations.irp.f index ede75df9..79bba2be 100644 --- a/plugins/Perturbation/epstein_nesbet.irp.f +++ b/plugins/Perturbation/pt2_equations.irp.f @@ -1,14 +1,9 @@ -subroutine pt2_epstein_nesbet(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist) +BEGIN_TEMPLATE + +subroutine pt2_epstein_nesbet ($arguments) use bitmasks implicit none - integer, intent(in) :: Nint,ndet,N_st - integer(bit_kind), intent(in) :: det_pert(Nint,2) - double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st) - double precision :: i_H_psi_array(N_st) - - integer, intent(in) :: N_minilist - integer, intent(in) :: idx_minilist(0:N_det_selectors) - integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) + $declarations BEGIN_DOC ! compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution @@ -23,6 +18,7 @@ subroutine pt2_epstein_nesbet(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_s integer :: i,j double precision :: diag_H_mat_elem, h + double precision :: i_H_psi_array(N_st) PROVIDE selection_criterion ASSERT (Nint == N_int) @@ -49,17 +45,10 @@ subroutine pt2_epstein_nesbet(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_s end -subroutine pt2_epstein_nesbet_2x2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist) +subroutine pt2_epstein_nesbet_2x2 ($arguments) use bitmasks implicit none - integer, intent(in) :: Nint,ndet,N_st - integer(bit_kind), intent(in) :: det_pert(Nint,2) - double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st) - double precision :: i_H_psi_array(N_st) - - integer, intent(in) :: N_minilist - integer, intent(in) :: idx_minilist(0:N_det_selectors) - integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) + $declarations BEGIN_DOC ! compute the Epstein-Nesbet 2x2 diagonalization coefficient and energetic contribution @@ -74,6 +63,7 @@ subroutine pt2_epstein_nesbet_2x2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet integer :: i,j double precision :: diag_H_mat_elem,delta_e, h + double precision :: i_H_psi_array(N_st) ASSERT (Nint == N_int) ASSERT (Nint > 0) PROVIDE CI_electronic_energy @@ -104,3 +94,71 @@ subroutine pt2_epstein_nesbet_2x2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet enddo end + +subroutine pt2_moller_plesset ($arguments) + use bitmasks + implicit none + $declarations + + BEGIN_DOC + ! compute the standard Moller-Plesset perturbative first order coefficient and second order energetic contribution + ! + ! for the various n_st states. + ! + ! c_pert(i) = /(difference of orbital energies) + ! + ! e_2_pert(i) = ^2/(difference of orbital energies) + ! + END_DOC + + integer :: i,j + double precision :: diag_H_mat_elem + 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(ref_bitmask,det_pert,exc,degree,phase,Nint) + if (degree == 2) then + call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2) + delta_e = Fock_matrix_diag_mo(h1) + Fock_matrix_diag_mo(h2) - & + (Fock_matrix_diag_mo(p1) + Fock_matrix_diag_mo(p2)) + delta_e = 1.d0/delta_e + else if (degree == 1) then + call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2) + delta_e = Fock_matrix_diag_mo(h1) - Fock_matrix_diag_mo(p1) + delta_e = 1.d0/delta_e + else + delta_e = 0.d0 + endif + + call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det,psi_selectors_size,n_st,i_H_psi_array) + h = diag_H_mat_elem(det_pert,Nint) + 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 +SUBST [ arguments, declarations ] + +det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist ; + + integer, intent(in) :: Nint + integer, intent(in) :: ndet + integer, intent(in) :: N_st + integer, intent(in) :: N_minilist + integer(bit_kind), intent(in) :: det_pert(Nint,2) + double precision , intent(out) :: c_pert(N_st) + double precision , intent(out) :: e_2_pert(N_st) + double precision, intent(out) :: H_pert_diag(N_st) + integer, intent(in) :: idx_minilist(0:N_det_selectors) + integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors) +;; + + +END_TEMPLATE +