qp2/plugins/local/fci_tc_bi/diagonalize_ci.irp.f

84 lines
2.9 KiB
Fortran

! ---
subroutine diagonalize_CI_tc_bi_ortho(ndet, E_tc, norm, pt2_data, print_pt2)
BEGIN_DOC
! Replace the coefficients of the CI states by the coefficients of the
! eigenstates of the CI matrix
END_DOC
use selection_types
implicit none
integer, intent(inout) :: ndet ! number of determinants from before
double precision, intent(inout) :: E_tc(N_states), norm(N_states) ! E and norm from previous wave function
type(pt2_type) , intent(in) :: pt2_data ! PT2 from previous wave function
logical, intent(in) :: print_pt2
integer :: i, j,k
double precision:: pt2_minus,pt2_plus,pt2_tot, pt2_abs,pt1_norm,rpt2_tot
double precision :: error_pt2_minus, error_pt2_plus, error_pt2_tot, error_pt2_abs
PROVIDE mo_l_coef mo_r_coef
! print*,'*****'
! print*,'New wave function information'
! print*,'N_det tc = ',N_det
! do k = 1, N_states
! print*,'************'
! print*,'State ',k
! pt2_plus = pt2_data % variance(k)
! pt2_minus = pt2_data % pt2(k)
! pt2_abs = pt2_plus - pt2_minus
! pt2_tot = pt2_plus + pt2_minus
!
! pt1_norm = pt2_data % overlap(k,k)
! rpt2_tot = pt2_tot / (1.d0 + pt1_norm)
!
!
! print*,'norm_ground_left_right_bi_orth = ',norm_ground_left_right_bi_orth(k)
! print*,'eigval_right_tc = ',eigval_right_tc_bi_orth(k)
! print*,'*****'
!
! if(print_pt2) then
! print*,'*****'
! print*,'previous wave function info'
! print*,'norm(before) = ',norm
! print*,'E(before) = ',E_tc
! print*,'PT1 norm = ',dsqrt(pt1_norm)
! print*,'PT2 = ',pt2_tot
! print*,'rPT2 = ',rpt2_tot
! print*,'|PT2| = ',pt2_abs
! print*,'Positive PT2 = ',pt2_plus
! print*,'Negative PT2 = ',pt2_minus
! print*,'E(before) + PT2 = ',E_tc + pt2_tot/norm
! print*,'E(before) +rPT2 = ',E_tc + rpt2_tot/norm
! write(*,'(A28,X,I10,X,100(F16.8,X))')'Ndet,E,E+PT2,E+RPT2,|PT2|=',ndet,E_tc ,E_tc + pt2_tot/norm,E_tc + rpt2_tot/norm,pt2_minus, pt2_plus
! print*,'*****'
! endif
! enddo
do k = 1, N_states
E_tc(k) = eigval_right_tc_bi_orth(k)
norm(k) = norm_ground_left_right_bi_orth(k)
enddo
psi_energy(1:N_states) = eigval_right_tc_bi_orth(1:N_states) - nuclear_repulsion
psi_s2(1:N_states) = s2_eigvec_tc_bi_orth(1:N_states)
ndet = N_det
do j = 1, N_states
do i = 1, N_det
psi_l_coef_bi_ortho(i,j) = leigvec_tc_bi_orth(i,j)
psi_r_coef_bi_ortho(i,j) = reigvec_tc_bi_orth(i,j)
psi_coef(i,j) = dabs(psi_l_coef_bi_ortho(i,j) * psi_r_coef_bi_ortho(i,j))
enddo
enddo
SOFT_TOUCH eigval_left_tc_bi_orth eigval_right_tc_bi_orth leigvec_tc_bi_orth reigvec_tc_bi_orth norm_ground_left_right_bi_orth
SOFT_TOUCH psi_l_coef_bi_ortho psi_r_coef_bi_ortho psi_coef psi_energy psi_s2
call save_tc_bi_ortho_wavefunction()
end
! ---