diff --git a/src/determinants/utils.irp.f b/src/determinants/utils.irp.f index 20d9e1e5..7e6ad39c 100644 --- a/src/determinants/utils.irp.f +++ b/src/determinants/utils.irp.f @@ -43,4 +43,57 @@ BEGIN_PROVIDER [ double precision, S2_matrix_all_dets,(N_det,N_det) ] !$OMP END PARALLEL DO END_PROVIDER +subroutine print_energy_components() + implicit none + BEGIN_DOC +! Prints the different components of the energy. + END_DOC + integer, save :: ifirst = 0 + double precision :: Vee, Ven, Vnn, Vecp, T, f + integer :: i,j,k + Vnn = nuclear_repulsion + + print *, 'Energy components' + print *, '=================' + print *, '' + do k=1,N_states + + Ven = 0.d0 + Vecp = 0.d0 + T = 0.d0 + + do j=1,mo_num + do i=1,mo_num + f = one_e_dm_mo_alpha(i,j,k) + one_e_dm_mo_beta(i,j,k) + Ven = Ven + f * mo_integrals_n_e(i,j) + Vecp = Vecp + f * mo_pseudo_integrals(i,j) + T = T + f * mo_kinetic_integrals(i,j) + enddo + enddo + Vee = psi_energy(k) - Ven - Vecp - T + + if (ifirst == 0) then + ifirst = 1 + print *, 'Vnn : Nucleus-Nucleus potential energy' + print *, 'Ven : Electron-Nucleus potential energy' + print *, 'Vee : Electron-Electron potential energy' + print *, 'Vecp : Potential energy of the pseudo-potentials' + print *, 'T : Electronic kinetic energy' + print *, '' + endif + + print *, 'State ', k + print *, '---------' + print *, '' + print *, 'Vnn = ', Vnn + print *, 'Ven = ', Ven + print *, 'Vee = ', Vee + print *, 'Vecp = ', Vecp + print *, 'T = ', T + print *, '' + enddo + + print *, '' + +end diff --git a/src/iterations/print_summary.irp.f b/src/iterations/print_summary.irp.f index a8037982..ad87bc8e 100644 --- a/src/iterations/print_summary.irp.f +++ b/src/iterations/print_summary.irp.f @@ -31,18 +31,19 @@ subroutine print_summary(e_,pt2_,error_,variance_,norm_,n_det_,n_occ_pattern_,n_ write(fmt,*) '(''# ============'',', N_states_p, '(1X,''=============================''))' write(*,fmt) - write(fmt,*) '(12X,', N_states_p, '(6X,A7,1X,I6,10X))' + write(fmt,*) '(13X,', N_states_p, '(6X,A7,1X,I6,10X))' write(*,fmt) ('State',k, k=1,N_states_p) write(fmt,*) '(''# ============'',', N_states_p, '(1X,''=============================''))' write(*,fmt) - write(fmt,*) '(A12,', N_states_p, '(1X,F14.8,15X))' + write(fmt,*) '(A13,', N_states_p, '(1X,F14.8,15X))' write(*,fmt) '# E ', e_(1:N_states_p) if (N_states_p > 1) then write(*,fmt) '# Excit. (au)', e_(1:N_states_p)-e_(1) write(*,fmt) '# Excit. (eV)', (e_(1:N_states_p)-e_(1))*27.211396641308d0 endif write(fmt,*) '(A13,', 2*N_states_p, '(1X,F14.8))' - write(*,fmt) '# PT2'//pt2_string, (pt2_(k), error_(k), k=1,N_states_p) + write(*,fmt) '# PT2 '//pt2_string, (pt2_(k), error_(k), k=1,N_states_p) + write(*,fmt) '# rPT2'//pt2_string, (pt2_(k)*f(k), error_(k)*f(k), k=1,N_states_p) write(*,'(A)') '#' write(*,fmt) '# E+PT2 ', (e_(k)+pt2_(k),error_(k), k=1,N_states_p) write(*,fmt) '# E+rPT2 ', (e_(k)+pt2_(k)*f(k),error_(k)*f(k), k=1,N_states_p) @@ -97,5 +98,7 @@ subroutine print_summary(e_,pt2_,error_,variance_,norm_,n_det_,n_occ_pattern_,n_ enddo endif + call print_energy_components() + end subroutine