PTEROSOR/QuantumPackage
10
0
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-06-24 06:02:26 +02:00
Code Releases Activity

updated green for qp2

Browse Source
This commit is contained in:
Kevin Gasperich 2020-06-03 16:13:16 -05:00
parent 25d0cbaa75
commit 0fd6eb3897
10 changed files with 82 additions and 80 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
src/green/green.main.irp.f
View File

@ -15,7 +15,7 @@ end
subroutine psicoefprinttest
implicit none
integer :: i
TOUCH psi_coef
TOUCH psi_coef_complex
print *, 'printing ndet', N_det
end
subroutine print_lanczos_eigvals

54
src/green/hu0_hp.irp.f
View File

@ -595,22 +595,22 @@ subroutine i_h_j_double_spin_hp(key_i,key_j,Nint,ispin,hij_hp,N_hp,spin_hp,sign_
double precision :: phase_hp(N_hp)
integer :: exc(0:2,2)
double precision :: phase
complex*16, external :: get_mo_bielec_integral
complex*16, external :: mo_two_e_integral_complex
integer :: i1,i2,i3,i4,j2,j3,ii
PROVIDE big_array_exchange_integrals mo_bielec_integrals_in_map
PROVIDE big_array_exchange_integrals_complex mo_two_e_integrals_in_map
call get_double_excitation_spin(key_i,key_j,exc,phase,Nint)
hij0 = phase*(get_mo_bielec_integral( &
hij0 = phase*(mo_two_e_integral_complex( &
exc(1,1), &
exc(2,1), &
exc(1,2), &
exc(2,2), mo_integrals_map) - &
get_mo_bielec_integral( &
exc(2,2)) - &
mo_two_e_integral_complex( &
exc(1,1), &
exc(2,1), &
exc(2,2), &
exc(1,2), mo_integrals_map) )
exc(1,2)) )
ASSERT (exc(1,1) < exc(2,1))
ASSERT (exc(1,2) < exc(2,2))
@ -661,14 +661,14 @@ subroutine i_h_j_mono_spin_hp(key_i,key_j,Nint,spin,hij_hp,N_hp,spin_hp,sign_hp,
integer :: exc(0:2,2)
double precision :: phase
PROVIDE big_array_exchange_integrals mo_bielec_integrals_in_map
PROVIDE big_array_exchange_integrals_complex mo_two_e_integrals_in_map
call get_mono_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint)
call get_single_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint)
call get_mono_excitation_from_fock_hp(key_i,key_j,exc(1,1),exc(1,2),spin,phase,N_hp,hij_hp,spin_hp,sign_hp,idx_hp,allowed_hp)
call get_single_excitation_from_fock_hp(key_i,key_j,exc(1,1),exc(1,2),spin,phase,N_hp,hij_hp,spin_hp,sign_hp,idx_hp,allowed_hp)
end
subroutine get_mono_excitation_from_fock_hp(det_1,det_2,h,p,spin,phase,N_hp,hij_hp,spin_hp,sign_hp,idx_hp,allowed_hp)
subroutine get_single_excitation_from_fock_hp(det_1,det_2,h,p,spin,phase,N_hp,hij_hp,spin_hp,sign_hp,idx_hp,allowed_hp)
use bitmasks
implicit none
integer,intent(in) :: h,p,spin,N_hp
@ -699,37 +699,37 @@ subroutine get_mono_excitation_from_fock_hp(det_1,det_2,h,p,spin,phase,N_hp,hij_
enddo
call bitstring_to_list_ab(hole, occ_hole, n_occ_ab_hole, N_int)
call bitstring_to_list_ab(partcl, occ_partcl, n_occ_ab_partcl, N_int)
hij0 = fock_operator_closed_shell_ref_bitmask(h,p)
hij0 = fock_op_cshell_ref_bitmask_cplx(h,p)
! holes :: direct terms
do i0 = 1, n_occ_ab_hole(1)
i = occ_hole(i0,1)
hij0 -= big_array_coulomb_integrals(i,h,p) ! get_mo_bielec_integral_schwartz(h,i,p,i,mo_integrals_map)
hij0 -= big_array_coulomb_integrals_complex(i,h,p) ! get_mo_bielec_integral_schwartz(h,i,p,i,mo_integrals_map)
enddo
do i0 = 1, n_occ_ab_hole(2)
i = occ_hole(i0,2)
hij0 -= big_array_coulomb_integrals(i,h,p) !get_mo_bielec_integral_schwartz(h,i,p,i,mo_integrals_map)
hij0 -= big_array_coulomb_integrals_complex(i,h,p) !get_mo_bielec_integral_schwartz(h,i,p,i,mo_integrals_map)
enddo
! holes :: exchange terms
do i0 = 1, n_occ_ab_hole(spin)
i = occ_hole(i0,spin)
hij0 += big_array_exchange_integrals(i,h,p) ! get_mo_bielec_integral_schwartz(h,i,i,p,mo_integrals_map)
hij0 += big_array_exchange_integrals_complex(i,h,p) ! get_mo_bielec_integral_schwartz(h,i,i,p,mo_integrals_map)
enddo
! particles :: direct terms
do i0 = 1, n_occ_ab_partcl(1)
i = occ_partcl(i0,1)
hij0 += big_array_coulomb_integrals(i,h,p)!get_mo_bielec_integral_schwartz(h,i,p,i,mo_integrals_map)
hij0 += big_array_coulomb_integrals_complex(i,h,p)!get_mo_bielec_integral_schwartz(h,i,p,i,mo_integrals_map)
enddo
do i0 = 1, n_occ_ab_partcl(2)
i = occ_partcl(i0,2)
hij0 += big_array_coulomb_integrals(i,h,p) !get_mo_bielec_integral_schwartz(h,i,p,i,mo_integrals_map)
hij0 += big_array_coulomb_integrals_complex(i,h,p) !get_mo_bielec_integral_schwartz(h,i,p,i,mo_integrals_map)
enddo
! particles :: exchange terms
do i0 = 1, n_occ_ab_partcl(spin)
i = occ_partcl(i0,spin)
hij0 -= big_array_exchange_integrals(i,h,p)!get_mo_bielec_integral_schwartz(h,i,i,p,mo_integrals_map)
hij0 -= big_array_exchange_integrals_complex(i,h,p)!get_mo_bielec_integral_schwartz(h,i,i,p,mo_integrals_map)
enddo
low=min(h,p)
@ -771,7 +771,7 @@ subroutine get_mono_excitation_from_fock_hp(det_1,det_2,h,p,spin,phase,N_hp,hij_
hij_hp(ii) = 0.d0
cycle
else if (spin.eq.spin_hp(ii)) then
hij_hp(ii) = hij0 + sign_hp(ii) *(big_array_coulomb_integrals(idx_hp(ii),h,p) - big_array_exchange_integrals(idx_hp(ii),h,p))
hij_hp(ii) = hij0 + sign_hp(ii) *(big_array_coulomb_integrals_complex(idx_hp(ii),h,p) - big_array_exchange_integrals_complex(idx_hp(ii),h,p))
if ((low.lt.idx_hp(ii)).and.(high.gt.idx_hp(ii))) then
phase_hp(ii) = -1.d0
else
@ -779,7 +779,7 @@ subroutine get_mono_excitation_from_fock_hp(det_1,det_2,h,p,spin,phase,N_hp,hij_
endif
else
phase_hp(ii) = 1.d0
hij_hp(ii) = hij0 + sign_hp(ii) * big_array_coulomb_integrals(idx_hp(ii),h,p)
hij_hp(ii) = hij0 + sign_hp(ii) * big_array_coulomb_integrals_complex(idx_hp(ii),h,p)
endif
hij_hp(ii) = hij_hp(ii) * phase * phase_hp(ii)
enddo
@ -806,24 +806,24 @@ subroutine i_H_j_double_alpha_beta_hp(key_i,key_j,Nint,hij_hp,N_hp,spin_hp,sign_
integer :: lowhigh(2,2)
integer :: exc(0:2,2,2)
double precision :: phase, phase2
complex*16, external :: get_mo_bielec_integral
complex*16, external :: mo_two_e_integral_complex
PROVIDE big_array_exchange_integrals mo_bielec_integrals_in_map
PROVIDE big_array_exchange_integrals_complex mo_two_e_integrals_in_map
call get_mono_excitation_spin(key_i(1,1),key_j(1,1),exc(0,1,1),phase,Nint)
call get_mono_excitation_spin(key_i(1,2),key_j(1,2),exc(0,1,2),phase2,Nint)
call get_single_excitation_spin(key_i(1,1),key_j(1,1),exc(0,1,1),phase,Nint)
call get_single_excitation_spin(key_i(1,2),key_j(1,2),exc(0,1,2),phase2,Nint)
phase = phase*phase2
if (exc(1,1,1) == exc(1,2,2)) then
hij0 = big_array_exchange_integrals(exc(1,1,1),exc(1,1,2),exc(1,2,1))
hij0 = big_array_exchange_integrals_complex(exc(1,1,1),exc(1,1,2),exc(1,2,1))
else if (exc(1,2,1) == exc(1,1,2)) then
hij0 = big_array_exchange_integrals(exc(1,2,1),exc(1,1,1),exc(1,2,2))
hij0 = big_array_exchange_integrals_complex(exc(1,2,1),exc(1,1,1),exc(1,2,2))
else
hij0 = get_mo_bielec_integral( &
hij0 = mo_two_e_integral_complex( &
exc(1,1,1), &
exc(1,1,2), &
exc(1,2,1), &
exc(1,2,2) ,mo_integrals_map)
exc(1,2,2))
endif
!todo: clean this up

10
src/green/hu0_lanczos.irp.f
View File

@ -194,7 +194,7 @@ subroutine H_u_0_openmp_work_$N_int(v_t,u_t,sze,istart,iend,ishift,istep)
ASSERT (lrow <= N_det_alpha_unique)
tmp_det2(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, lrow)
call i_h_j_double_alpha_beta(tmp_det,tmp_det2,$N_int,hij)
call i_h_j_double_alpha_beta_complex(tmp_det,tmp_det2,$N_int,hij)
v_t(k_a) = v_t(k_a) + hij * u_t(l_a)
enddo
enddo
@ -264,7 +264,7 @@ subroutine H_u_0_openmp_work_$N_int(v_t,u_t,sze,istart,iend,ishift,istep)
ASSERT (lrow <= N_det_alpha_unique)
tmp_det2(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, lrow)
call i_H_j_mono_spin( tmp_det, tmp_det2, $N_int, 1, hij)
call i_h_j_single_spin_complex( tmp_det, tmp_det2, $N_int, 1, hij)
v_t(k_a) = v_t(k_a) + hij * u_t(l_a)
enddo
@ -280,7 +280,7 @@ subroutine H_u_0_openmp_work_$N_int(v_t,u_t,sze,istart,iend,ishift,istep)
lrow = psi_bilinear_matrix_rows(l_a)
ASSERT (lrow <= N_det_alpha_unique)
call i_H_j_double_spin( tmp_det(1,1), psi_det_alpha_unique(1, lrow), $N_int, hij)
call i_h_j_double_spin_complex( tmp_det(1,1), psi_det_alpha_unique(1, lrow), $N_int, hij)
v_t(k_a) = v_t(k_a) + hij * u_t(l_a)
enddo
@ -341,7 +341,7 @@ subroutine H_u_0_openmp_work_$N_int(v_t,u_t,sze,istart,iend,ishift,istep)
ASSERT (lcol <= N_det_beta_unique)
tmp_det2(1:$N_int,2) = psi_det_beta_unique (1:$N_int, lcol)
call i_H_j_mono_spin( tmp_det, tmp_det2, $N_int, 2, hij)
call i_h_j_single_spin_complex( tmp_det, tmp_det2, $N_int, 2, hij)
l_a = psi_bilinear_matrix_transp_order(l_b)
ASSERT (l_a <= N_det)
v_t(k_a) = v_t(k_a) + hij * u_t(l_a)
@ -357,7 +357,7 @@ subroutine H_u_0_openmp_work_$N_int(v_t,u_t,sze,istart,iend,ishift,istep)
lcol = psi_bilinear_matrix_transp_columns(l_b)
ASSERT (lcol <= N_det_beta_unique)
call i_H_j_double_spin( tmp_det(1,2), psi_det_beta_unique(1, lcol), $N_int, hij)
call i_h_j_double_spin_complex( tmp_det(1,2), psi_det_beta_unique(1, lcol), $N_int, hij)
l_a = psi_bilinear_matrix_transp_order(l_b)
ASSERT (l_a <= N_det)

13
src/green/lanczos.irp.f
View File

@ -28,7 +28,7 @@ END_PROVIDER
integer :: idx_homo_lumo(2), spin_homo_lumo(2)
logical :: has_idx,has_spin,has_sign,has_lanc
integer :: nlanc
! needs psi_det, mo_tot_num, N_int, mo_bielec_integral_jj, mo_mono_elec_integral_diag
! needs psi_det, mo_num, N_int, mo_bielec_integral_jj, mo_mono_elec_integral_diag
call ezfio_has_green_green_idx(has_idx)
call ezfio_has_green_green_spin(has_spin)
call ezfio_has_green_green_sign(has_sign)
@ -43,7 +43,7 @@ END_PROVIDER
stop 'problem with lanczos restart; need idx, spin, sign'
else
print*,'new lanczos calculation, finding homo/lumo'
call get_homo_lumo(psi_det(1:N_int,1:2,1),N_int,mo_tot_num,idx_homo_lumo,spin_homo_lumo)
call get_homo_lumo(psi_det(1:N_int,1:2,1),N_int,mo_num,idx_homo_lumo,spin_homo_lumo)
! homo
green_idx(1)=idx_homo_lumo(1)
@ -75,7 +75,7 @@ END_PROVIDER
! endif
! else
! print*,'new lanczos calculation, finding homo/lumo'
! call get_homo_lumo(psi_det(1:N_int,1:2,1),N_int,mo_tot_num,idx_homo_lumo,spin_homo_lumo)
! call get_homo_lumo(psi_det(1:N_int,1:2,1),N_int,mo_num,idx_homo_lumo,spin_homo_lumo)
!
! ! homo
! green_idx(1)=idx_homo_lumo(1)
@ -279,9 +279,9 @@ BEGIN_PROVIDER [ double precision, spectral_lanczos, (n_omega,n_green_vec) ]
logical :: has_ci_energy
double precision :: ref_energy_0
PROVIDE delta_omega alpha_lanczos beta_lanczos omega_list
call ezfio_has_full_ci_zmq_energy(has_ci_energy)
call ezfio_has_fci_energy(has_ci_energy)
if (has_ci_energy) then
call ezfio_get_full_ci_zmq_energy(ref_energy_0)
call ezfio_get_fci_energy(ref_energy_0)
else
print*,'no reference energy from full_ci_zmq, exiting'
stop
@ -469,7 +469,8 @@ subroutine lanczos_h_step_hp(uu,vv,work,sze,alpha_i,beta_i,ng,spin_hp,sign_hp,id
complex*16, intent(inout) :: uu(sze,ng),vv(sze,ng)
complex*16, intent(out) :: work(sze,ng)
double precision, intent(out) :: alpha_i(ng), beta_i(ng)
integer, intent(in) :: spin_hp(ng), sign_hp(ng), idx_hp(ng)
integer, intent(in) :: spin_hp(ng), idx_hp(ng)
double precision, intent(in) :: sign_hp(ng)
double precision, external :: dznrm2
complex*16, external :: u_dot_v_complex

2
src/green/plot-spec-dens.py
View File

@ -23,7 +23,7 @@ def printspec(ezdir,wmin,wmax,nw,eps):
gdir=ezdir+'/green/'
with open(gdir+'n_green_vec') as infile:
ngvec=int(infile.readline().strip())
with open(ezdir+'/full_ci_zmq/energy') as infile:
with open(ezdir+'/fci/energy') as infile:
e0=float(infile.readline().strip())
with open(gdir+'n_lanczos_complete') as infile:
nlanc=int(infile.readline().strip())

2
src/green/print_e_mo_debug.irp.f
View File

@ -11,5 +11,5 @@ subroutine routine
implicit none
integer :: i
read*,i
call print_mo_energies(psi_det(:,:,i),N_int,mo_tot_num)
call print_mo_energies(psi_det(:,:,i),N_int,mo_num)
end

6
src/green/print_h_debug.irp.f
View File

@ -53,7 +53,7 @@ subroutine routine
print*,'H matrix '
double precision :: s2
complex*16 :: ref_h_matrix
ref_h_matrix = h_matrix_all_dets(1,1)
ref_h_matrix = h_matrix_all_dets_complex(1,1)
print*,'HF like determinant energy = ',ref_bitmask_energy+nuclear_repulsion
print*,'Ref element of H_matrix = ',ref_h_matrix+nuclear_repulsion
print*,'Printing the H matrix ...'
@ -64,7 +64,7 @@ subroutine routine
!enddo
do i = 1, N_det
H_matrix_all_dets(i,i) += nuclear_repulsion
H_matrix_all_dets_complex(i,i) += nuclear_repulsion
enddo
!do i = 5,N_det
@ -79,7 +79,7 @@ subroutine routine
! TODO: change for complex
do i = 1, N_det
write(*,'(I3,X,A3,2000(E24.15))')i,' | ',H_matrix_all_dets(i,:)
write(*,'(I3,X,A3,2000(E24.15))')i,' | ',H_matrix_all_dets_complex(i,:)
enddo
! print*,''

2
src/green/print_h_omp_debug.irp.f
View File

@ -30,7 +30,7 @@ subroutine routine_omp
u_tmp(i,i)=(1.d0,0.d0)
enddo
call h_s2_u_0_nstates_openmp(v_tmp,s_tmp,u_tmp,n_st,h_size)
call h_s2_u_0_nstates_openmp_complex(v_tmp,s_tmp,u_tmp,n_st,h_size)
do i = 1, n_st
v_tmp(i,i) += nuclear_repulsion
enddo

70
src/green/utils_hp.irp.f
View File

@ -26,8 +26,8 @@ subroutine print_mo_energies(key_ref,nint,nmo)
enddo
call bitstring_to_list_ab(key_virt,virt,n_virt,nint)
e_mo(1:nmo,1)=mo_mono_elec_integral_diag(1:nmo)
e_mo(1:nmo,2)=mo_mono_elec_integral_diag(1:nmo)
e_mo(1:nmo,1)=mo_one_e_integrals_diag(1:nmo)
e_mo(1:nmo,2)=mo_one_e_integrals_diag(1:nmo)
do ispin=1,2
jspin=int_spin2(ispin)
@ -36,23 +36,23 @@ subroutine print_mo_energies(key_ref,nint,nmo)
is_occ(i,ispin)=1
do j0=i0+1,n_occ(ispin)
j=occ(j0,ispin)
e_mo(i,ispin) = e_mo(i,ispin) + mo_bielec_integral_jj_anti(i,j)
e_mo(j,ispin) = e_mo(j,ispin) + mo_bielec_integral_jj_anti(i,j)
e_mo(i,ispin) = e_mo(i,ispin) + mo_two_e_integrals_jj_anti(i,j)
e_mo(j,ispin) = e_mo(j,ispin) + mo_two_e_integrals_jj_anti(i,j)
enddo
do k=2,ispin
do j0=1,n_occ(jspin)
j=occ(j0,jspin)
e_mo(i,ispin) = e_mo(i,ispin) + mo_bielec_integral_jj(i,j)
e_mo(j,jspin) = e_mo(j,jspin) + mo_bielec_integral_jj(i,j) !can delete this and remove k level of loop
e_mo(i,ispin) = e_mo(i,ispin) + mo_two_e_integrals_jj(i,j)
e_mo(j,jspin) = e_mo(j,jspin) + mo_two_e_integrals_jj(i,j) !can delete this and remove k level of loop
enddo
enddo
do j0=1,n_virt(ispin)
j=virt(j0,ispin)
e_mo(j,ispin) = e_mo(j,ispin) + mo_bielec_integral_jj_anti(i,j)
e_mo(j,ispin) = e_mo(j,ispin) + mo_two_e_integrals_jj_anti(i,j)
enddo
do j0=1,n_virt(jspin)
j=virt(j0,jspin)
e_mo(j,jspin) = e_mo(j,jspin) + mo_bielec_integral_jj(i,j)
e_mo(j,jspin) = e_mo(j,jspin) + mo_two_e_integrals_jj(i,j)
enddo
enddo
enddo
@ -89,8 +89,8 @@ subroutine get_mo_energies(key_ref,nint,nmo,e_mo)
enddo
call bitstring_to_list_ab(key_virt,virt,n_virt,nint)
e_mo(1:nmo,1)=mo_mono_elec_integral_diag(1:nmo)
e_mo(1:nmo,2)=mo_mono_elec_integral_diag(1:nmo)
e_mo(1:nmo,1)=mo_one_e_integrals_diag(1:nmo)
e_mo(1:nmo,2)=mo_one_e_integrals_diag(1:nmo)
do ispin=1,2
jspin=int_spin2(ispin)
@ -98,23 +98,23 @@ subroutine get_mo_energies(key_ref,nint,nmo,e_mo)
i=occ(i0,ispin)
do j0=i0+1,n_occ(ispin)
j=occ(j0,ispin)
e_mo(i,ispin) = e_mo(i,ispin) + mo_bielec_integral_jj_anti(i,j)
e_mo(j,ispin) = e_mo(j,ispin) + mo_bielec_integral_jj_anti(i,j)
e_mo(i,ispin) = e_mo(i,ispin) + mo_two_e_integrals_jj_anti(i,j)
e_mo(j,ispin) = e_mo(j,ispin) + mo_two_e_integrals_jj_anti(i,j)
enddo
do k=2,ispin
do j0=1,n_occ(jspin)
j=occ(j0,jspin)
e_mo(i,ispin) = e_mo(i,ispin) + mo_bielec_integral_jj(i,j)
e_mo(j,jspin) = e_mo(j,jspin) + mo_bielec_integral_jj(i,j) !can delete this and remove k level of loop
e_mo(i,ispin) = e_mo(i,ispin) + mo_two_e_integrals_jj(i,j)
e_mo(j,jspin) = e_mo(j,jspin) + mo_two_e_integrals_jj(i,j) !can delete this and remove k level of loop
enddo
enddo
do j0=1,n_virt(ispin)
j=virt(j0,ispin)
e_mo(j,ispin) = e_mo(j,ispin) + mo_bielec_integral_jj_anti(i,j)
e_mo(j,ispin) = e_mo(j,ispin) + mo_two_e_integrals_jj_anti(i,j)
enddo
do j0=1,n_virt(jspin)
j=virt(j0,jspin)
e_mo(j,jspin) = e_mo(j,jspin) + mo_bielec_integral_jj(i,j)
e_mo(j,jspin) = e_mo(j,jspin) + mo_two_e_integrals_jj(i,j)
enddo
enddo
enddo
@ -524,17 +524,17 @@ subroutine ac_operator_phase(key_new,key_ref,iorb,ispin,Nint,phase)
if (ispin==1) then
parity_filled=0_bit_kind
else
parity_filled=iand(elec_alpha_num,1_bit_kind)
parity_filled=iand(int(elec_alpha_num,bit_kind),1_bit_kind)
endif
! get parity due to orbs in other ints (with lower indices)
do i=1,k-1
parity_filled = iand(popcnt(key_ref(i,ispin)),parity_filled)
parity_filled = iand(int(popcnt(key_ref(i,ispin)),bit_kind),parity_filled)
enddo
! get parity due to orbs in same int as iorb
! ishft(1_bit_kind,l)-1 has its l rightmost bits set to 1, other bits set to 0
parity_filled = iand(popcnt(iand(ishft(1_bit_kind,l)-1,key_ref(k,ispin))),parity_filled)
parity_filled = iand(int(popcnt(iand(ishft(1_bit_kind,l)-1,key_ref(k,ispin))),bit_kind),parity_filled)
phase = p(iand(1_bit_kind,parity_filled))
end
@ -585,30 +585,30 @@ subroutine a_operator_phase(key_new,key_ref,iorb,ispin,Nint,phase)
if (ispin==1) then
parity_filled=0_bit_kind
else
parity_filled=iand(elec_alpha_num,1_bit_kind)
parity_filled=iand(int(elec_alpha_num,bit_kind),1_bit_kind)
endif
! get parity due to orbs in other ints (with lower indices)
do i=1,k-1
parity_filled = iand(popcnt(key_ref(i,ispin)),parity_filled)
parity_filled = iand(int(popcnt(key_ref(i,ispin)),bit_kind),parity_filled)
enddo
! get parity due to orbs in same int as iorb
! ishft(1_bit_kind,l)-1 has its l rightmost bits set to 1, other bits set to 0
parity_filled = iand(popcnt(iand(ishft(1_bit_kind,l)-1,key_ref(k,ispin))),parity_filled)
parity_filled = iand(int(popcnt(iand(ishft(1_bit_kind,l)-1,key_ref(k,ispin))),bit_kind),parity_filled)
phase = p(iand(1_bit_kind,parity_filled))
end
BEGIN_PROVIDER [ double precision, mo_mono_elec_integral_diag,(mo_tot_num)]
implicit none
integer :: i
BEGIN_DOC
! diagonal elements of mo_mono_elec_integral array
END_DOC
print*,'Providing the mono electronic integrals (diagonal)'
do i = 1, mo_tot_num
mo_mono_elec_integral_diag(i) = real(mo_mono_elec_integral(i,i))
enddo
END_PROVIDER
!BEGIN_PROVIDER [ double precision, mo_mono_elec_integral_diag,(mo_num)]
! implicit none
! integer :: i
! BEGIN_DOC
! ! diagonal elements of mo_mono_elec_integral array
! END_DOC
! print*,'Providing the mono electronic integrals (diagonal)'
!
! do i = 1, mo_num
! mo_mono_elec_integral_diag(i) = real(mo_mono_elec_integral(i,i))
! enddo
!
!END_PROVIDER

1
src/utils/constants.include.F
View File

@ -7,6 +7,7 @@ double precision, parameter :: sqpi = dsqrt(dacos(-1.d0))
double precision, parameter :: pi_5_2 = 34.9868366552d0
double precision, parameter :: dfour_pi = 4.d0*dacos(-1.d0)
double precision, parameter :: dtwo_pi = 2.d0*dacos(-1.d0)
double precision, parameter :: inv_pi = 1.d0/dacos(-1.d0)
double precision, parameter :: inv_sq_pi = 1.d0/dsqrt(dacos(-1.d0))
double precision, parameter :: inv_sq_pi_2 = 0.5d0/dsqrt(dacos(-1.d0))
double precision, parameter :: thresh = 1.d-15