mirror of
https://github.com/QuantumPackage/qp2.git
synced 2024-08-15 00:38:34 +02:00
working on complex determinants
This commit is contained in:
parent
702ba79af8
commit
6d12abf088
@ -248,6 +248,11 @@ BEGIN_PROVIDER [ double precision, one_e_spin_density_mo, (mo_num,mo_num) ]
|
|||||||
END_PROVIDER
|
END_PROVIDER
|
||||||
|
|
||||||
subroutine set_natural_mos
|
subroutine set_natural_mos
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
implicit none
|
implicit none
|
||||||
BEGIN_DOC
|
BEGIN_DOC
|
||||||
! Set natural orbitals, obtained by diagonalization of the one-body density matrix
|
! Set natural orbitals, obtained by diagonalization of the one-body density matrix
|
||||||
@ -274,6 +279,11 @@ subroutine set_natural_mos
|
|||||||
|
|
||||||
end
|
end
|
||||||
subroutine save_natural_mos
|
subroutine save_natural_mos
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
implicit none
|
implicit none
|
||||||
BEGIN_DOC
|
BEGIN_DOC
|
||||||
! Save natural orbitals, obtained by diagonalization of the one-body density matrix in
|
! Save natural orbitals, obtained by diagonalization of the one-body density matrix in
|
||||||
@ -292,11 +302,19 @@ BEGIN_PROVIDER [ double precision, c0_weight, (N_states) ]
|
|||||||
if (N_states > 1) then
|
if (N_states > 1) then
|
||||||
integer :: i
|
integer :: i
|
||||||
double precision :: c
|
double precision :: c
|
||||||
|
if (is_complex) then
|
||||||
|
do i=1,N_states
|
||||||
|
c0_weight(i) = 1.d-31
|
||||||
|
c = maxval(cdabs(psi_coef_complex(:,i) * psi_coef_complex(:,i)))
|
||||||
|
c0_weight(i) = 1.d0/(c+1.d-20)
|
||||||
|
enddo
|
||||||
|
else
|
||||||
do i=1,N_states
|
do i=1,N_states
|
||||||
c0_weight(i) = 1.d-31
|
c0_weight(i) = 1.d-31
|
||||||
c = maxval(psi_coef(:,i) * psi_coef(:,i))
|
c = maxval(psi_coef(:,i) * psi_coef(:,i))
|
||||||
c0_weight(i) = 1.d0/(c+1.d-20)
|
c0_weight(i) = 1.d0/(c+1.d-20)
|
||||||
enddo
|
enddo
|
||||||
|
endif
|
||||||
c = 1.d0/minval(c0_weight(:))
|
c = 1.d0/minval(c0_weight(:))
|
||||||
do i=1,N_states
|
do i=1,N_states
|
||||||
c0_weight(i) = c0_weight(i) * c
|
c0_weight(i) = c0_weight(i) * c
|
||||||
@ -398,8 +416,23 @@ subroutine get_occupation_from_dets(istate,occupation)
|
|||||||
ASSERT (istate <= N_states)
|
ASSERT (istate <= N_states)
|
||||||
|
|
||||||
occupation = 0.d0
|
occupation = 0.d0
|
||||||
double precision, external :: u_dot_u
|
|
||||||
|
|
||||||
|
if (is_complex) then
|
||||||
|
double precision, external :: u_dot_u_complex
|
||||||
|
norm_2 = 1.d0/u_dot_u_complex(psi_coef_complex(1,istate),N_det)
|
||||||
|
|
||||||
|
do i=1,N_det
|
||||||
|
c = cdabs(psi_coef_complex(i,istate)*psi_coef_complex(i,istate))*norm_2
|
||||||
|
call bitstring_to_list_ab(psi_det(1,1,i), list, n_elements, N_int)
|
||||||
|
do ispin=1,2
|
||||||
|
do j=1,n_elements(ispin)
|
||||||
|
ASSERT ( list(j,ispin) < mo_num )
|
||||||
|
occupation( list(j,ispin) ) += c
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
else
|
||||||
|
double precision, external :: u_dot_u
|
||||||
norm_2 = 1.d0/u_dot_u(psi_coef(1,istate),N_det)
|
norm_2 = 1.d0/u_dot_u(psi_coef(1,istate),N_det)
|
||||||
|
|
||||||
do i=1,N_det
|
do i=1,N_det
|
||||||
@ -412,5 +445,6 @@ subroutine get_occupation_from_dets(istate,occupation)
|
|||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
endif
|
||||||
end
|
end
|
||||||
|
|
||||||
|
311
src/determinants/density_matrix_complex.irp.f
Normal file
311
src/determinants/density_matrix_complex.irp.f
Normal file
@ -0,0 +1,311 @@
|
|||||||
|
BEGIN_PROVIDER [ complex*16, one_e_dm_mo_alpha_average_complex, (mo_num,mo_num) ]
|
||||||
|
&BEGIN_PROVIDER [ complex*16, one_e_dm_mo_beta_average_complex, (mo_num,mo_num) ]
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! $\alpha$ and $\beta$ one-body density matrix for each state
|
||||||
|
END_DOC
|
||||||
|
integer :: i
|
||||||
|
one_e_dm_mo_alpha_average_complex = (0.d0,0.d0)
|
||||||
|
one_e_dm_mo_beta_average_complex = (0.d0,0.d0)
|
||||||
|
do i = 1,N_states
|
||||||
|
one_e_dm_mo_alpha_average_complex(:,:) += one_e_dm_mo_alpha_complex(:,:,i) * state_average_weight(i)
|
||||||
|
one_e_dm_mo_beta_average_complex(:,:) += one_e_dm_mo_beta_complex(:,:,i) * state_average_weight(i)
|
||||||
|
enddo
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, one_e_dm_mo_diff_complex, (mo_num,mo_num,2:N_states) ]
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Difference of the one-body density matrix with respect to the ground state
|
||||||
|
END_DOC
|
||||||
|
integer :: i,j, istate
|
||||||
|
|
||||||
|
do istate=2,N_states
|
||||||
|
do j=1,mo_num
|
||||||
|
do i=1,mo_num
|
||||||
|
one_e_dm_mo_diff_complex(i,j,istate) = &
|
||||||
|
one_e_dm_mo_alpha_complex(i,j,istate) - one_e_dm_mo_alpha_complex(i,j,1) +&
|
||||||
|
one_e_dm_mo_beta_complex (i,j,istate) - one_e_dm_mo_beta_complex (i,j,1)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, one_e_dm_mo_spin_index_complex, (mo_num,mo_num,N_states,2) ]
|
||||||
|
implicit none
|
||||||
|
integer :: i,j,ispin,istate
|
||||||
|
ispin = 1
|
||||||
|
do istate = 1, N_states
|
||||||
|
do j = 1, mo_num
|
||||||
|
do i = 1, mo_num
|
||||||
|
one_e_dm_mo_spin_index_complex(i,j,istate,ispin) = one_e_dm_mo_alpha_complex(i,j,istate)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
ispin = 2
|
||||||
|
do istate = 1, N_states
|
||||||
|
do j = 1, mo_num
|
||||||
|
do i = 1, mo_num
|
||||||
|
one_e_dm_mo_spin_index_complex(i,j,istate,ispin) = one_e_dm_mo_beta_complex(i,j,istate)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, one_e_dm_dagger_mo_spin_index_complex, (mo_num,mo_num,N_states,2) ]
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
! implicit none
|
||||||
|
! integer :: i,j,ispin,istate
|
||||||
|
! ispin = 1
|
||||||
|
! do istate = 1, N_states
|
||||||
|
! do j = 1, mo_num
|
||||||
|
! one_e_dm_dagger_mo_spin_index(j,j,istate,ispin) = 1 - one_e_dm_mo_alpha(j,j,istate)
|
||||||
|
! do i = j+1, mo_num
|
||||||
|
! one_e_dm_dagger_mo_spin_index(i,j,istate,ispin) = -one_e_dm_mo_alpha(i,j,istate)
|
||||||
|
! one_e_dm_dagger_mo_spin_index(j,i,istate,ispin) = -one_e_dm_mo_alpha(i,j,istate)
|
||||||
|
! enddo
|
||||||
|
! enddo
|
||||||
|
! enddo
|
||||||
|
!
|
||||||
|
! ispin = 2
|
||||||
|
! do istate = 1, N_states
|
||||||
|
! do j = 1, mo_num
|
||||||
|
! one_e_dm_dagger_mo_spin_index(j,j,istate,ispin) = 1 - one_e_dm_mo_beta(j,j,istate)
|
||||||
|
! do i = j+1, mo_num
|
||||||
|
! one_e_dm_dagger_mo_spin_index(i,j,istate,ispin) = -one_e_dm_mo_beta(i,j,istate)
|
||||||
|
! one_e_dm_dagger_mo_spin_index(j,i,istate,ispin) = -one_e_dm_mo_beta(i,j,istate)
|
||||||
|
! enddo
|
||||||
|
! enddo
|
||||||
|
! enddo
|
||||||
|
!
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, one_e_dm_mo_alpha_complex, (mo_num,mo_num,N_states) ]
|
||||||
|
&BEGIN_PROVIDER [ complex*16, one_e_dm_mo_beta_complex, (mo_num,mo_num,N_states) ]
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! $\alpha$ and $\beta$ one-body density matrix for each state
|
||||||
|
! $\gamma_{\mu\nu} = \langle\Psi|a_{\nu}^{\dagger}a_{\mu}|\Psi\rangle$
|
||||||
|
! $\gamma_{\mu\nu} = \langle a_{\nu} \Psi|a_{\mu} \Psi\rangle$
|
||||||
|
! $\gamma_{\mu\nu} = \sum_{IJ} c^*_J c_I \langle a_{\nu} I|a_{\mu} J\rangle$
|
||||||
|
END_DOC
|
||||||
|
|
||||||
|
integer :: j,k,l,m,k_a,k_b
|
||||||
|
integer :: occ(N_int*bit_kind_size,2)
|
||||||
|
complex*16 :: ck, cl, ckl
|
||||||
|
double precision :: phase
|
||||||
|
integer :: h1,h2,p1,p2,s1,s2, degree
|
||||||
|
integer(bit_kind) :: tmp_det(N_int,2), tmp_det2(N_int)
|
||||||
|
integer :: exc(0:2,2),n_occ(2)
|
||||||
|
complex*16, allocatable :: tmp_a(:,:,:), tmp_b(:,:,:)
|
||||||
|
integer :: krow, kcol, lrow, lcol
|
||||||
|
|
||||||
|
PROVIDE psi_det
|
||||||
|
|
||||||
|
one_e_dm_mo_alpha_complex = (0.d0,0.d0)
|
||||||
|
one_e_dm_mo_beta_complex = (0.d0,0.d0)
|
||||||
|
!$OMP PARALLEL DEFAULT(NONE) &
|
||||||
|
!$OMP PRIVATE(j,k,k_a,k_b,l,m,occ,ck, cl, ckl,phase,h1,h2,p1,p2,s1,s2, degree,exc,&
|
||||||
|
!$OMP tmp_a, tmp_b, n_occ, krow, kcol, lrow, lcol, tmp_det, tmp_det2)&
|
||||||
|
!$OMP SHARED(psi_det,psi_coef_complex,N_int,N_states,elec_alpha_num, &
|
||||||
|
!$OMP elec_beta_num,one_e_dm_mo_alpha_complex,one_e_dm_mo_beta_complex,N_det,&
|
||||||
|
!$OMP mo_num,psi_bilinear_matrix_rows,psi_bilinear_matrix_columns,&
|
||||||
|
!$OMP psi_bilinear_matrix_transp_rows, psi_bilinear_matrix_transp_columns,&
|
||||||
|
!$OMP psi_bilinear_matrix_order_reverse, psi_det_alpha_unique, psi_det_beta_unique,&
|
||||||
|
!$OMP psi_bilinear_matrix_values_complex, psi_bilinear_matrix_transp_values_complex,&
|
||||||
|
!$OMP N_det_alpha_unique,N_det_beta_unique,irp_here)
|
||||||
|
allocate(tmp_a(mo_num,mo_num,N_states), tmp_b(mo_num,mo_num,N_states) )
|
||||||
|
tmp_a = (0.d0,0.d0)
|
||||||
|
!$OMP DO SCHEDULE(dynamic,64)
|
||||||
|
do k_a=1,N_det
|
||||||
|
krow = psi_bilinear_matrix_rows(k_a)
|
||||||
|
ASSERT (krow <= N_det_alpha_unique)
|
||||||
|
|
||||||
|
kcol = psi_bilinear_matrix_columns(k_a)
|
||||||
|
ASSERT (kcol <= N_det_beta_unique)
|
||||||
|
|
||||||
|
tmp_det(1:N_int,1) = psi_det_alpha_unique(1:N_int,krow)
|
||||||
|
tmp_det(1:N_int,2) = psi_det_beta_unique (1:N_int,kcol)
|
||||||
|
|
||||||
|
! Diagonal part
|
||||||
|
! -------------
|
||||||
|
|
||||||
|
call bitstring_to_list_ab(tmp_det, occ, n_occ, N_int)
|
||||||
|
do m=1,N_states
|
||||||
|
ck = cdabs(psi_bilinear_matrix_values_complex(k_a,m)*psi_bilinear_matrix_values_complex(k_a,m))
|
||||||
|
do l=1,elec_alpha_num
|
||||||
|
j = occ(l,1)
|
||||||
|
tmp_a(j,j,m) += ck
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
if (k_a == N_det) cycle
|
||||||
|
l = k_a+1
|
||||||
|
lrow = psi_bilinear_matrix_rows(l)
|
||||||
|
lcol = psi_bilinear_matrix_columns(l)
|
||||||
|
! Fix beta determinant, loop over alphas
|
||||||
|
do while ( lcol == kcol )
|
||||||
|
tmp_det2(:) = psi_det_alpha_unique(:, lrow)
|
||||||
|
call get_excitation_degree_spin(tmp_det(1,1),tmp_det2,degree,N_int)
|
||||||
|
if (degree == 1) then
|
||||||
|
exc = 0
|
||||||
|
call get_single_excitation_spin(tmp_det(1,1),tmp_det2,exc,phase,N_int)
|
||||||
|
call decode_exc_spin(exc,h1,p1,h2,p2)
|
||||||
|
! h1 occ in k
|
||||||
|
! p1 occ in l
|
||||||
|
do m=1,N_states
|
||||||
|
ckl = dconjg(psi_bilinear_matrix_values_complex(k_a,m))*psi_bilinear_matrix_values_complex(l,m) * phase
|
||||||
|
tmp_a(h1,p1,m) += dconjg(ckl)
|
||||||
|
tmp_a(p1,h1,m) += ckl
|
||||||
|
enddo
|
||||||
|
endif
|
||||||
|
l = l+1
|
||||||
|
if (l>N_det) exit
|
||||||
|
lrow = psi_bilinear_matrix_rows(l)
|
||||||
|
lcol = psi_bilinear_matrix_columns(l)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
enddo
|
||||||
|
!$OMP END DO NOWAIT
|
||||||
|
|
||||||
|
!$OMP CRITICAL
|
||||||
|
one_e_dm_mo_alpha_complex(:,:,:) = one_e_dm_mo_alpha_complex(:,:,:) + tmp_a(:,:,:)
|
||||||
|
!$OMP END CRITICAL
|
||||||
|
deallocate(tmp_a)
|
||||||
|
|
||||||
|
tmp_b = (0.d0,0.d0)
|
||||||
|
!$OMP DO SCHEDULE(dynamic,64)
|
||||||
|
do k_b=1,N_det
|
||||||
|
krow = psi_bilinear_matrix_transp_rows(k_b)
|
||||||
|
ASSERT (krow <= N_det_alpha_unique)
|
||||||
|
|
||||||
|
kcol = psi_bilinear_matrix_transp_columns(k_b)
|
||||||
|
ASSERT (kcol <= N_det_beta_unique)
|
||||||
|
|
||||||
|
tmp_det(1:N_int,1) = psi_det_alpha_unique(1:N_int,krow)
|
||||||
|
tmp_det(1:N_int,2) = psi_det_beta_unique (1:N_int,kcol)
|
||||||
|
|
||||||
|
! Diagonal part
|
||||||
|
! -------------
|
||||||
|
|
||||||
|
call bitstring_to_list_ab(tmp_det, occ, n_occ, N_int)
|
||||||
|
do m=1,N_states
|
||||||
|
ck = cdabs(psi_bilinear_matrix_transp_values_complex(k_b,m)*psi_bilinear_matrix_transp_values_complex(k_b,m))
|
||||||
|
do l=1,elec_beta_num
|
||||||
|
j = occ(l,2)
|
||||||
|
tmp_b(j,j,m) += ck
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
if (k_b == N_det) cycle
|
||||||
|
l = k_b+1
|
||||||
|
lrow = psi_bilinear_matrix_transp_rows(l)
|
||||||
|
lcol = psi_bilinear_matrix_transp_columns(l)
|
||||||
|
! Fix beta determinant, loop over alphas
|
||||||
|
do while ( lrow == krow )
|
||||||
|
tmp_det2(:) = psi_det_beta_unique(:, lcol)
|
||||||
|
call get_excitation_degree_spin(tmp_det(1,2),tmp_det2,degree,N_int)
|
||||||
|
if (degree == 1) then
|
||||||
|
exc = 0
|
||||||
|
call get_single_excitation_spin(tmp_det(1,2),tmp_det2,exc,phase,N_int)
|
||||||
|
call decode_exc_spin(exc,h1,p1,h2,p2)
|
||||||
|
do m=1,N_states
|
||||||
|
ckl = dconjg(psi_bilinear_matrix_transp_values_complex(k_b,m))*psi_bilinear_matrix_transp_values_complex(l,m) * phase
|
||||||
|
tmp_b(h1,p1,m) += dconjg(ckl)
|
||||||
|
tmp_b(p1,h1,m) += ckl
|
||||||
|
enddo
|
||||||
|
endif
|
||||||
|
l = l+1
|
||||||
|
if (l>N_det) exit
|
||||||
|
lrow = psi_bilinear_matrix_transp_rows(l)
|
||||||
|
lcol = psi_bilinear_matrix_transp_columns(l)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
enddo
|
||||||
|
!$OMP END DO NOWAIT
|
||||||
|
!$OMP CRITICAL
|
||||||
|
one_e_dm_mo_beta_complex(:,:,:) = one_e_dm_mo_beta_complex(:,:,:) + tmp_b(:,:,:)
|
||||||
|
!$OMP END CRITICAL
|
||||||
|
|
||||||
|
deallocate(tmp_b)
|
||||||
|
!$OMP END PARALLEL
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, one_e_dm_mo_complex, (mo_num,mo_num) ]
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! One-body density matrix
|
||||||
|
END_DOC
|
||||||
|
one_e_dm_mo_complex = one_e_dm_mo_alpha_average_complex + one_e_dm_mo_beta_average_complex
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, one_e_spin_density_mo_complex, (mo_num,mo_num) ]
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! $\rho(\alpha) - \rho(\beta)$
|
||||||
|
END_DOC
|
||||||
|
one_e_spin_density_mo_complex = one_e_dm_mo_alpha_average_complex - one_e_dm_mo_beta_average_complex
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, one_e_spin_density_ao_complex, (ao_num,ao_num) ]
|
||||||
|
BEGIN_DOC
|
||||||
|
! One body spin density matrix on the |AO| basis : $\rho_{AO}(\alpha) - \rho_{AO}(\beta)$
|
||||||
|
! todo: verify that this is correct for complex
|
||||||
|
! equivalent to using mo_to_ao_no_overlap?
|
||||||
|
END_DOC
|
||||||
|
implicit none
|
||||||
|
integer :: i,j,k,l
|
||||||
|
complex*16 :: dm_mo
|
||||||
|
|
||||||
|
one_e_spin_density_ao_complex = (0.d0,0.d0)
|
||||||
|
do k = 1, ao_num
|
||||||
|
do l = 1, ao_num
|
||||||
|
do i = 1, mo_num
|
||||||
|
do j = 1, mo_num
|
||||||
|
dm_mo = one_e_spin_density_mo_complex(j,i)
|
||||||
|
! if(dabs(dm_mo).le.1.d-10)cycle
|
||||||
|
one_e_spin_density_ao_complex(l,k) += dconjg(mo_coef_complex(k,i)) * mo_coef_complex(l,j) * dm_mo
|
||||||
|
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, one_e_dm_ao_alpha_complex, (ao_num,ao_num) ]
|
||||||
|
&BEGIN_PROVIDER [ complex*16, one_e_dm_ao_beta_complex, (ao_num,ao_num) ]
|
||||||
|
BEGIN_DOC
|
||||||
|
! One body density matrix on the |AO| basis : $\rho_{AO}(\alpha), \rho_{AO}(\beta)$.
|
||||||
|
END_DOC
|
||||||
|
implicit none
|
||||||
|
integer :: i,j,k,l
|
||||||
|
complex*16 :: mo_alpha,mo_beta
|
||||||
|
|
||||||
|
one_e_dm_ao_alpha = (0.d0,0.d0)
|
||||||
|
one_e_dm_ao_beta = (0.d0,0.d0)
|
||||||
|
do k = 1, ao_num
|
||||||
|
do l = 1, ao_num
|
||||||
|
do i = 1, mo_num
|
||||||
|
do j = 1, mo_num
|
||||||
|
mo_alpha = one_e_dm_mo_alpha_average_complex(j,i)
|
||||||
|
mo_beta = one_e_dm_mo_beta_average_complex(j,i)
|
||||||
|
! if(dabs(dm_mo).le.1.d-10)cycle
|
||||||
|
one_e_dm_ao_alpha_complex(l,k) += dconjg(mo_coef_complex(k,i)) * mo_coef_complex(l,j) * mo_alpha
|
||||||
|
one_e_dm_ao_beta_complex(l,k) += dconjg(mo_coef_complex(k,i)) * mo_coef_complex(l,j) * mo_beta
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
|
@ -1,4 +1,5 @@
|
|||||||
subroutine u_0_S2_u_0_complex(e_0,u_0,n,keys_tmp,Nint,N_st,sze_8)
|
subroutine u_0_S2_u_0_complex(e_0,u_0,n,keys_tmp,Nint,N_st,sze_8)
|
||||||
|
!todo: modify/implement for complex
|
||||||
print*,irp_here,' not implemented for complex'
|
print*,irp_here,' not implemented for complex'
|
||||||
stop -1
|
stop -1
|
||||||
! use bitmasks
|
! use bitmasks
|
||||||
@ -28,6 +29,7 @@ end
|
|||||||
|
|
||||||
|
|
||||||
subroutine S2_u_0_complex(v_0,u_0,n,keys_tmp,Nint)
|
subroutine S2_u_0_complex(v_0,u_0,n,keys_tmp,Nint)
|
||||||
|
!todo: modify/implement for complex
|
||||||
print*,irp_here,' not implemented for complex'
|
print*,irp_here,' not implemented for complex'
|
||||||
stop -1
|
stop -1
|
||||||
! use bitmasks
|
! use bitmasks
|
||||||
@ -46,6 +48,7 @@ subroutine S2_u_0_complex(v_0,u_0,n,keys_tmp,Nint)
|
|||||||
end
|
end
|
||||||
|
|
||||||
subroutine S2_u_0_nstates_complex(v_0,u_0,n,keys_tmp,Nint,N_st,sze_8)
|
subroutine S2_u_0_nstates_complex(v_0,u_0,n,keys_tmp,Nint,N_st,sze_8)
|
||||||
|
!todo: modify/implement for complex
|
||||||
print*,irp_here,' not implemented for complex'
|
print*,irp_here,' not implemented for complex'
|
||||||
stop -1
|
stop -1
|
||||||
! use bitmasks
|
! use bitmasks
|
||||||
@ -180,6 +183,7 @@ end
|
|||||||
|
|
||||||
|
|
||||||
subroutine get_uJ_s2_uI_complex(psi_keys_tmp,psi_coefs_tmp,n,nmax_coefs,nmax_keys,s2,nstates)
|
subroutine get_uJ_s2_uI_complex(psi_keys_tmp,psi_coefs_tmp,n,nmax_coefs,nmax_keys,s2,nstates)
|
||||||
|
!todo: modify/implement for complex
|
||||||
print*,irp_here,' not implemented for complex'
|
print*,irp_here,' not implemented for complex'
|
||||||
stop -1
|
stop -1
|
||||||
! implicit none
|
! implicit none
|
||||||
@ -232,6 +236,7 @@ end
|
|||||||
|
|
||||||
|
|
||||||
subroutine i_S2_psi_minilist_complex(key,keys,idx_key,N_minilist,coef,Nint,Ndet,Ndet_max,Nstate,i_S2_psi_array)
|
subroutine i_S2_psi_minilist_complex(key,keys,idx_key,N_minilist,coef,Nint,Ndet,Ndet_max,Nstate,i_S2_psi_array)
|
||||||
|
!todo: modify/implement for complex
|
||||||
print*,irp_here,' not implemented for complex'
|
print*,irp_here,' not implemented for complex'
|
||||||
stop -1
|
stop -1
|
||||||
! use bitmasks
|
! use bitmasks
|
||||||
|
@ -159,3 +159,148 @@ subroutine get_single_excitation_from_fock(det_1,det_2,h,p,spin,phase,hij)
|
|||||||
|
|
||||||
end
|
end
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [complex*16, fock_operator_closed_shell_ref_bitmask_complex, (mo_num, mo_num) ]
|
||||||
|
implicit none
|
||||||
|
integer :: i0,j0,i,j,k0,k
|
||||||
|
integer :: n_occ_ab(2)
|
||||||
|
integer :: occ(N_int*bit_kind_size,2)
|
||||||
|
integer :: n_occ_ab_virt(2)
|
||||||
|
integer :: occ_virt(N_int*bit_kind_size,2)
|
||||||
|
integer(bit_kind) :: key_test(N_int)
|
||||||
|
integer(bit_kind) :: key_virt(N_int,2)
|
||||||
|
complex*16 :: accu
|
||||||
|
|
||||||
|
call bitstring_to_list_ab(ref_closed_shell_bitmask, occ, n_occ_ab, N_int)
|
||||||
|
do i = 1, N_int
|
||||||
|
key_virt(i,1) = full_ijkl_bitmask(i)
|
||||||
|
key_virt(i,2) = full_ijkl_bitmask(i)
|
||||||
|
key_virt(i,1) = xor(key_virt(i,1),ref_closed_shell_bitmask(i,1))
|
||||||
|
key_virt(i,2) = xor(key_virt(i,2),ref_closed_shell_bitmask(i,2))
|
||||||
|
enddo
|
||||||
|
complex*16, allocatable :: array_coulomb(:),array_exchange(:)
|
||||||
|
allocate (array_coulomb(mo_num),array_exchange(mo_num))
|
||||||
|
call bitstring_to_list_ab(key_virt, occ_virt, n_occ_ab_virt, N_int)
|
||||||
|
! docc ---> virt single excitations
|
||||||
|
do i0 = 1, n_occ_ab(1)
|
||||||
|
i=occ(i0,1)
|
||||||
|
do j0 = 1, n_occ_ab_virt(1)
|
||||||
|
j = occ_virt(j0,1)
|
||||||
|
! <ia|ja>
|
||||||
|
call get_mo_two_e_integrals_coulomb_ii_complex(i,j,mo_num,array_coulomb,mo_integrals_map,mo_integrals_map_2)
|
||||||
|
! <ia|aj>
|
||||||
|
call get_mo_two_e_integrals_exch_ii_complex(i,j,mo_num,array_exchange,mo_integrals_map,mo_integrals_map_2)
|
||||||
|
accu = (0.d0,0.d0)
|
||||||
|
do k0 = 1, n_occ_ab(1)
|
||||||
|
k = occ(k0,1)
|
||||||
|
accu += 2.d0 * array_coulomb(k) - array_exchange(k)
|
||||||
|
enddo
|
||||||
|
fock_operator_closed_shell_ref_bitmask_complex(i,j) = accu + mo_one_e_integrals_complex(i,j)
|
||||||
|
!fock_operator_closed_shell_ref_bitmask_complex(j,i) = dconjg(accu) + mo_one_e_integrals_complex(j,i)
|
||||||
|
fock_operator_closed_shell_ref_bitmask_complex(j,i) = dconjg(fock_operator_closed_shell_ref_bitmask_complex(i,j))
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! virt ---> virt single excitations
|
||||||
|
do i0 = 1, n_occ_ab_virt(1)
|
||||||
|
i=occ_virt(i0,1)
|
||||||
|
do j0 = 1, n_occ_ab_virt(1)
|
||||||
|
j = occ_virt(j0,1)
|
||||||
|
call get_mo_two_e_integrals_coulomb_ii_complex(i,j,mo_num,array_coulomb,mo_integrals_map,mo_integrals_map_2)
|
||||||
|
call get_mo_two_e_integrals_exch_ii_complex(i,j,mo_num,array_exchange,mo_integrals_map,mo_integrals_map_2)
|
||||||
|
accu = (0.d0,0.d0)
|
||||||
|
do k0 = 1, n_occ_ab(1)
|
||||||
|
k = occ(k0,1)
|
||||||
|
accu += 2.d0 * array_coulomb(k) - array_exchange(k)
|
||||||
|
enddo
|
||||||
|
fock_operator_closed_shell_ref_bitmask_complex(i,j) = accu+ mo_one_e_integrals_complex(i,j)
|
||||||
|
fock_operator_closed_shell_ref_bitmask_complex(j,i) = dconjg(accu)+ mo_one_e_integrals_complex(j,i)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! docc ---> docc single excitations
|
||||||
|
do i0 = 1, n_occ_ab(1)
|
||||||
|
i=occ(i0,1)
|
||||||
|
do j0 = 1, n_occ_ab(1)
|
||||||
|
j = occ(j0,1)
|
||||||
|
call get_mo_two_e_integrals_coulomb_ii_complex(i,j,mo_num,array_coulomb,mo_integrals_map,mo_integrals_map_2)
|
||||||
|
call get_mo_two_e_integrals_exch_ii_complex(i,j,mo_num,array_exchange,mo_integrals_map,mo_integrals_map_2)
|
||||||
|
accu = (0.d0,0.d0)
|
||||||
|
do k0 = 1, n_occ_ab(1)
|
||||||
|
k = occ(k0,1)
|
||||||
|
accu += 2.d0 * array_coulomb(k) - array_exchange(k)
|
||||||
|
enddo
|
||||||
|
fock_operator_closed_shell_ref_bitmask_complex(i,j) = accu+ mo_one_e_integrals_complex(i,j)
|
||||||
|
fock_operator_closed_shell_ref_bitmask_complex(j,i) = dconjg(accu)+ mo_one_e_integrals_complex(j,i)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
deallocate(array_coulomb,array_exchange)
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
subroutine get_single_excitation_from_fock_complex(det_1,det_2,h,p,spin,phase,hij)
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
integer,intent(in) :: h,p,spin
|
||||||
|
double precision, intent(in) :: phase
|
||||||
|
integer(bit_kind), intent(in) :: det_1(N_int,2), det_2(N_int,2)
|
||||||
|
complex*16, intent(out) :: hij
|
||||||
|
integer(bit_kind) :: differences(N_int,2)
|
||||||
|
integer(bit_kind) :: hole(N_int,2)
|
||||||
|
integer(bit_kind) :: partcl(N_int,2)
|
||||||
|
integer :: occ_hole(N_int*bit_kind_size,2)
|
||||||
|
integer :: occ_partcl(N_int*bit_kind_size,2)
|
||||||
|
integer :: n_occ_ab_hole(2),n_occ_ab_partcl(2)
|
||||||
|
integer :: i0,i
|
||||||
|
complex*16 :: buffer_c(mo_num),buffer_x(mo_num)
|
||||||
|
do i=1, mo_num
|
||||||
|
buffer_c(i) = big_array_coulomb_integrals_complex(i,h,p)
|
||||||
|
buffer_x(i) = big_array_exchange_integrals_complex(i,h,p)
|
||||||
|
enddo
|
||||||
|
do i = 1, N_int
|
||||||
|
differences(i,1) = xor(det_1(i,1),ref_closed_shell_bitmask(i,1))
|
||||||
|
differences(i,2) = xor(det_1(i,2),ref_closed_shell_bitmask(i,2))
|
||||||
|
hole(i,1) = iand(differences(i,1),ref_closed_shell_bitmask(i,1))
|
||||||
|
hole(i,2) = iand(differences(i,2),ref_closed_shell_bitmask(i,2))
|
||||||
|
partcl(i,1) = iand(differences(i,1),det_1(i,1))
|
||||||
|
partcl(i,2) = iand(differences(i,2),det_1(i,2))
|
||||||
|
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)
|
||||||
|
hij = fock_operator_closed_shell_ref_bitmask_complex(h,p)
|
||||||
|
! holes :: direct terms
|
||||||
|
do i0 = 1, n_occ_ab_hole(1)
|
||||||
|
i = occ_hole(i0,1)
|
||||||
|
hij -= buffer_c(i)
|
||||||
|
enddo
|
||||||
|
do i0 = 1, n_occ_ab_hole(2)
|
||||||
|
i = occ_hole(i0,2)
|
||||||
|
hij -= buffer_c(i)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! holes :: exchange terms
|
||||||
|
do i0 = 1, n_occ_ab_hole(spin)
|
||||||
|
i = occ_hole(i0,spin)
|
||||||
|
hij += buffer_x(i)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! particles :: direct terms
|
||||||
|
do i0 = 1, n_occ_ab_partcl(1)
|
||||||
|
i = occ_partcl(i0,1)
|
||||||
|
hij += buffer_c(i)
|
||||||
|
enddo
|
||||||
|
do i0 = 1, n_occ_ab_partcl(2)
|
||||||
|
i = occ_partcl(i0,2)
|
||||||
|
hij += buffer_c(i)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
! particles :: exchange terms
|
||||||
|
do i0 = 1, n_occ_ab_partcl(spin)
|
||||||
|
i = occ_partcl(i0,spin)
|
||||||
|
hij -= buffer_x(i)
|
||||||
|
enddo
|
||||||
|
hij = hij * phase
|
||||||
|
|
||||||
|
end
|
||||||
|
|
||||||
|
@ -1581,9 +1581,12 @@ subroutine get_excitation_degree_vector(key1,key2,degree,Nint,sze,idx)
|
|||||||
end
|
end
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
double precision function diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
|
double precision function diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
use bitmasks
|
use bitmasks
|
||||||
implicit none
|
implicit none
|
||||||
BEGIN_DOC
|
BEGIN_DOC
|
||||||
@ -2292,3 +2295,592 @@ subroutine connected_to_hf(key_i,yes_no)
|
|||||||
yes_no = .True.
|
yes_no = .True.
|
||||||
endif
|
endif
|
||||||
end
|
end
|
||||||
|
|
||||||
|
|
||||||
|
!==============================================================================!
|
||||||
|
! !
|
||||||
|
! Complex !
|
||||||
|
! !
|
||||||
|
!==============================================================================!
|
||||||
|
|
||||||
|
|
||||||
|
subroutine i_H_j_s2_complex(key_i,key_j,Nint,hij,s2)
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Returns $\langle i|H|j \rangle$ and $\langle i|S^2|j \rangle$
|
||||||
|
! where $i$ and $j$ are determinants.
|
||||||
|
END_DOC
|
||||||
|
integer, intent(in) :: Nint
|
||||||
|
integer(bit_kind), intent(in) :: key_i(Nint,2), key_j(Nint,2)
|
||||||
|
double precision, intent(out) :: hij, s2
|
||||||
|
|
||||||
|
integer :: exc(0:2,2,2)
|
||||||
|
integer :: degree
|
||||||
|
double precision :: get_two_e_integral
|
||||||
|
integer :: m,n,p,q
|
||||||
|
integer :: i,j,k
|
||||||
|
integer :: occ(Nint*bit_kind_size,2)
|
||||||
|
double precision :: diag_H_mat_elem, phase
|
||||||
|
integer :: n_occ_ab(2)
|
||||||
|
PROVIDE mo_two_e_integrals_in_map mo_integrals_map big_array_exchange_integrals
|
||||||
|
|
||||||
|
ASSERT (Nint > 0)
|
||||||
|
ASSERT (Nint == N_int)
|
||||||
|
ASSERT (sum(popcnt(key_i(:,1))) == elec_alpha_num)
|
||||||
|
ASSERT (sum(popcnt(key_i(:,2))) == elec_beta_num)
|
||||||
|
ASSERT (sum(popcnt(key_j(:,1))) == elec_alpha_num)
|
||||||
|
ASSERT (sum(popcnt(key_j(:,2))) == elec_beta_num)
|
||||||
|
|
||||||
|
hij = 0.d0
|
||||||
|
s2 = 0d0
|
||||||
|
!DIR$ FORCEINLINE
|
||||||
|
call get_excitation_degree(key_i,key_j,degree,Nint)
|
||||||
|
integer :: spin
|
||||||
|
select case (degree)
|
||||||
|
case (2)
|
||||||
|
call get_double_excitation(key_i,key_j,exc,phase,Nint)
|
||||||
|
! Single alpha, single beta
|
||||||
|
if (exc(0,1,1) == 1) then
|
||||||
|
if ( (exc(1,1,1) == exc(1,2,2)).and.(exc(1,1,2) == exc(1,2,1)) ) then
|
||||||
|
s2 = -phase
|
||||||
|
endif
|
||||||
|
if(exc(1,1,1) == exc(1,2,2) )then
|
||||||
|
hij = phase * big_array_exchange_integrals(exc(1,1,1),exc(1,1,2),exc(1,2,1))
|
||||||
|
else if (exc(1,2,1) ==exc(1,1,2))then
|
||||||
|
hij = phase * big_array_exchange_integrals(exc(1,2,1),exc(1,1,1),exc(1,2,2))
|
||||||
|
else
|
||||||
|
hij = phase*get_two_e_integral( &
|
||||||
|
exc(1,1,1), &
|
||||||
|
exc(1,1,2), &
|
||||||
|
exc(1,2,1), &
|
||||||
|
exc(1,2,2) ,mo_integrals_map)
|
||||||
|
endif
|
||||||
|
! Double alpha
|
||||||
|
else if (exc(0,1,1) == 2) then
|
||||||
|
hij = phase*(get_two_e_integral( &
|
||||||
|
exc(1,1,1), &
|
||||||
|
exc(2,1,1), &
|
||||||
|
exc(1,2,1), &
|
||||||
|
exc(2,2,1) ,mo_integrals_map) - &
|
||||||
|
get_two_e_integral( &
|
||||||
|
exc(1,1,1), &
|
||||||
|
exc(2,1,1), &
|
||||||
|
exc(2,2,1), &
|
||||||
|
exc(1,2,1) ,mo_integrals_map) )
|
||||||
|
! Double beta
|
||||||
|
else if (exc(0,1,2) == 2) then
|
||||||
|
hij = phase*(get_two_e_integral( &
|
||||||
|
exc(1,1,2), &
|
||||||
|
exc(2,1,2), &
|
||||||
|
exc(1,2,2), &
|
||||||
|
exc(2,2,2) ,mo_integrals_map) - &
|
||||||
|
get_two_e_integral( &
|
||||||
|
exc(1,1,2), &
|
||||||
|
exc(2,1,2), &
|
||||||
|
exc(2,2,2), &
|
||||||
|
exc(1,2,2) ,mo_integrals_map) )
|
||||||
|
endif
|
||||||
|
case (1)
|
||||||
|
call get_single_excitation(key_i,key_j,exc,phase,Nint)
|
||||||
|
!DIR$ FORCEINLINE
|
||||||
|
call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint)
|
||||||
|
! Single alpha
|
||||||
|
if (exc(0,1,1) == 1) then
|
||||||
|
m = exc(1,1,1)
|
||||||
|
p = exc(1,2,1)
|
||||||
|
spin = 1
|
||||||
|
! Single beta
|
||||||
|
else
|
||||||
|
m = exc(1,1,2)
|
||||||
|
p = exc(1,2,2)
|
||||||
|
spin = 2
|
||||||
|
endif
|
||||||
|
call get_single_excitation_from_fock(key_i,key_j,p,m,spin,phase,hij)
|
||||||
|
|
||||||
|
case (0)
|
||||||
|
double precision, external :: diag_S_mat_elem
|
||||||
|
s2 = diag_S_mat_elem(key_i,Nint)
|
||||||
|
hij = diag_H_mat_elem(key_i,Nint)
|
||||||
|
end select
|
||||||
|
end
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
subroutine i_H_j_complex(key_i,key_j,Nint,hij)
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Returns $\langle i|H|j \rangle$ where $i$ and $j$ are determinants.
|
||||||
|
END_DOC
|
||||||
|
integer, intent(in) :: Nint
|
||||||
|
integer(bit_kind), intent(in) :: key_i(Nint,2), key_j(Nint,2)
|
||||||
|
double precision, intent(out) :: hij
|
||||||
|
|
||||||
|
integer :: exc(0:2,2,2)
|
||||||
|
integer :: degree
|
||||||
|
double precision :: get_two_e_integral
|
||||||
|
integer :: m,n,p,q
|
||||||
|
integer :: i,j,k
|
||||||
|
integer :: occ(Nint*bit_kind_size,2)
|
||||||
|
double precision :: diag_H_mat_elem, phase
|
||||||
|
integer :: n_occ_ab(2)
|
||||||
|
PROVIDE mo_two_e_integrals_in_map mo_integrals_map big_array_exchange_integrals
|
||||||
|
|
||||||
|
ASSERT (Nint > 0)
|
||||||
|
ASSERT (Nint == N_int)
|
||||||
|
ASSERT (sum(popcnt(key_i(:,1))) == elec_alpha_num)
|
||||||
|
ASSERT (sum(popcnt(key_i(:,2))) == elec_beta_num)
|
||||||
|
ASSERT (sum(popcnt(key_j(:,1))) == elec_alpha_num)
|
||||||
|
ASSERT (sum(popcnt(key_j(:,2))) == elec_beta_num)
|
||||||
|
|
||||||
|
|
||||||
|
hij = 0.d0
|
||||||
|
!DIR$ FORCEINLINE
|
||||||
|
call get_excitation_degree(key_i,key_j,degree,Nint)
|
||||||
|
integer :: spin
|
||||||
|
select case (degree)
|
||||||
|
case (2)
|
||||||
|
call get_double_excitation(key_i,key_j,exc,phase,Nint)
|
||||||
|
if (exc(0,1,1) == 1) then
|
||||||
|
! Single alpha, single beta
|
||||||
|
if(exc(1,1,1) == exc(1,2,2) )then
|
||||||
|
hij = phase * big_array_exchange_integrals(exc(1,1,1),exc(1,1,2),exc(1,2,1))
|
||||||
|
else if (exc(1,2,1) ==exc(1,1,2))then
|
||||||
|
hij = phase * big_array_exchange_integrals(exc(1,2,1),exc(1,1,1),exc(1,2,2))
|
||||||
|
else
|
||||||
|
hij = phase*get_two_e_integral( &
|
||||||
|
exc(1,1,1), &
|
||||||
|
exc(1,1,2), &
|
||||||
|
exc(1,2,1), &
|
||||||
|
exc(1,2,2) ,mo_integrals_map)
|
||||||
|
endif
|
||||||
|
else if (exc(0,1,1) == 2) then
|
||||||
|
! Double alpha
|
||||||
|
hij = phase*(get_two_e_integral( &
|
||||||
|
exc(1,1,1), &
|
||||||
|
exc(2,1,1), &
|
||||||
|
exc(1,2,1), &
|
||||||
|
exc(2,2,1) ,mo_integrals_map) - &
|
||||||
|
get_two_e_integral( &
|
||||||
|
exc(1,1,1), &
|
||||||
|
exc(2,1,1), &
|
||||||
|
exc(2,2,1), &
|
||||||
|
exc(1,2,1) ,mo_integrals_map) )
|
||||||
|
else if (exc(0,1,2) == 2) then
|
||||||
|
! Double beta
|
||||||
|
hij = phase*(get_two_e_integral( &
|
||||||
|
exc(1,1,2), &
|
||||||
|
exc(2,1,2), &
|
||||||
|
exc(1,2,2), &
|
||||||
|
exc(2,2,2) ,mo_integrals_map) - &
|
||||||
|
get_two_e_integral( &
|
||||||
|
exc(1,1,2), &
|
||||||
|
exc(2,1,2), &
|
||||||
|
exc(2,2,2), &
|
||||||
|
exc(1,2,2) ,mo_integrals_map) )
|
||||||
|
endif
|
||||||
|
case (1)
|
||||||
|
call get_single_excitation(key_i,key_j,exc,phase,Nint)
|
||||||
|
!DIR$ FORCEINLINE
|
||||||
|
call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint)
|
||||||
|
if (exc(0,1,1) == 1) then
|
||||||
|
! Single alpha
|
||||||
|
m = exc(1,1,1)
|
||||||
|
p = exc(1,2,1)
|
||||||
|
spin = 1
|
||||||
|
else
|
||||||
|
! Single beta
|
||||||
|
m = exc(1,1,2)
|
||||||
|
p = exc(1,2,2)
|
||||||
|
spin = 2
|
||||||
|
endif
|
||||||
|
call get_single_excitation_from_fock(key_i,key_j,p,m,spin,phase,hij)
|
||||||
|
|
||||||
|
case (0)
|
||||||
|
hij = diag_H_mat_elem(key_i,Nint)
|
||||||
|
end select
|
||||||
|
end
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
subroutine i_H_j_verbose_complex(key_i,key_j,Nint,hij,hmono,hdouble,phase)
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Returns $\langle i|H|j \rangle$ where $i$ and $j$ are determinants.
|
||||||
|
END_DOC
|
||||||
|
integer, intent(in) :: Nint
|
||||||
|
integer(bit_kind), intent(in) :: key_i(Nint,2), key_j(Nint,2)
|
||||||
|
double precision, intent(out) :: hij,hmono,hdouble,phase
|
||||||
|
|
||||||
|
integer :: exc(0:2,2,2)
|
||||||
|
integer :: degree
|
||||||
|
double precision :: get_two_e_integral
|
||||||
|
integer :: m,n,p,q
|
||||||
|
integer :: i,j,k
|
||||||
|
integer :: occ(Nint*bit_kind_size,2)
|
||||||
|
double precision :: diag_H_mat_elem
|
||||||
|
integer :: n_occ_ab(2)
|
||||||
|
logical :: has_mipi(Nint*bit_kind_size)
|
||||||
|
double precision :: mipi(Nint*bit_kind_size), miip(Nint*bit_kind_size)
|
||||||
|
PROVIDE mo_two_e_integrals_in_map mo_integrals_map
|
||||||
|
|
||||||
|
ASSERT (Nint > 0)
|
||||||
|
ASSERT (Nint == N_int)
|
||||||
|
ASSERT (sum(popcnt(key_i(:,1))) == elec_alpha_num)
|
||||||
|
ASSERT (sum(popcnt(key_i(:,2))) == elec_beta_num)
|
||||||
|
ASSERT (sum(popcnt(key_j(:,1))) == elec_alpha_num)
|
||||||
|
ASSERT (sum(popcnt(key_j(:,2))) == elec_beta_num)
|
||||||
|
|
||||||
|
hij = 0.d0
|
||||||
|
hmono = 0.d0
|
||||||
|
hdouble = 0.d0
|
||||||
|
!DIR$ FORCEINLINE
|
||||||
|
call get_excitation_degree(key_i,key_j,degree,Nint)
|
||||||
|
select case (degree)
|
||||||
|
case (2)
|
||||||
|
call get_double_excitation(key_i,key_j,exc,phase,Nint)
|
||||||
|
if (exc(0,1,1) == 1) then
|
||||||
|
! Single alpha, single beta
|
||||||
|
hij = phase*get_two_e_integral( &
|
||||||
|
exc(1,1,1), &
|
||||||
|
exc(1,1,2), &
|
||||||
|
exc(1,2,1), &
|
||||||
|
exc(1,2,2) ,mo_integrals_map)
|
||||||
|
else if (exc(0,1,1) == 2) then
|
||||||
|
! Double alpha
|
||||||
|
hij = phase*(get_two_e_integral( &
|
||||||
|
exc(1,1,1), &
|
||||||
|
exc(2,1,1), &
|
||||||
|
exc(1,2,1), &
|
||||||
|
exc(2,2,1) ,mo_integrals_map) - &
|
||||||
|
get_two_e_integral( &
|
||||||
|
exc(1,1,1), &
|
||||||
|
exc(2,1,1), &
|
||||||
|
exc(2,2,1), &
|
||||||
|
exc(1,2,1) ,mo_integrals_map) )
|
||||||
|
|
||||||
|
else if (exc(0,1,2) == 2) then
|
||||||
|
! Double beta
|
||||||
|
hij = phase*(get_two_e_integral( &
|
||||||
|
exc(1,1,2), &
|
||||||
|
exc(2,1,2), &
|
||||||
|
exc(1,2,2), &
|
||||||
|
exc(2,2,2) ,mo_integrals_map) - &
|
||||||
|
get_two_e_integral( &
|
||||||
|
exc(1,1,2), &
|
||||||
|
exc(2,1,2), &
|
||||||
|
exc(2,2,2), &
|
||||||
|
exc(1,2,2) ,mo_integrals_map) )
|
||||||
|
endif
|
||||||
|
case (1)
|
||||||
|
call get_single_excitation(key_i,key_j,exc,phase,Nint)
|
||||||
|
!DIR$ FORCEINLINE
|
||||||
|
call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint)
|
||||||
|
has_mipi = .False.
|
||||||
|
if (exc(0,1,1) == 1) then
|
||||||
|
! Single alpha
|
||||||
|
m = exc(1,1,1)
|
||||||
|
p = exc(1,2,1)
|
||||||
|
do k = 1, elec_alpha_num
|
||||||
|
i = occ(k,1)
|
||||||
|
if (.not.has_mipi(i)) then
|
||||||
|
mipi(i) = get_two_e_integral(m,i,p,i,mo_integrals_map)
|
||||||
|
miip(i) = get_two_e_integral(m,i,i,p,mo_integrals_map)
|
||||||
|
has_mipi(i) = .True.
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
do k = 1, elec_beta_num
|
||||||
|
i = occ(k,2)
|
||||||
|
if (.not.has_mipi(i)) then
|
||||||
|
mipi(i) = get_two_e_integral(m,i,p,i,mo_integrals_map)
|
||||||
|
has_mipi(i) = .True.
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
|
||||||
|
do k = 1, elec_alpha_num
|
||||||
|
hdouble = hdouble + mipi(occ(k,1)) - miip(occ(k,1))
|
||||||
|
enddo
|
||||||
|
do k = 1, elec_beta_num
|
||||||
|
hdouble = hdouble + mipi(occ(k,2))
|
||||||
|
enddo
|
||||||
|
|
||||||
|
else
|
||||||
|
! Single beta
|
||||||
|
m = exc(1,1,2)
|
||||||
|
p = exc(1,2,2)
|
||||||
|
do k = 1, elec_beta_num
|
||||||
|
i = occ(k,2)
|
||||||
|
if (.not.has_mipi(i)) then
|
||||||
|
mipi(i) = get_two_e_integral(m,i,p,i,mo_integrals_map)
|
||||||
|
miip(i) = get_two_e_integral(m,i,i,p,mo_integrals_map)
|
||||||
|
has_mipi(i) = .True.
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
do k = 1, elec_alpha_num
|
||||||
|
i = occ(k,1)
|
||||||
|
if (.not.has_mipi(i)) then
|
||||||
|
mipi(i) = get_two_e_integral(m,i,p,i,mo_integrals_map)
|
||||||
|
has_mipi(i) = .True.
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
|
||||||
|
do k = 1, elec_alpha_num
|
||||||
|
hdouble = hdouble + mipi(occ(k,1))
|
||||||
|
enddo
|
||||||
|
do k = 1, elec_beta_num
|
||||||
|
hdouble = hdouble + mipi(occ(k,2)) - miip(occ(k,2))
|
||||||
|
enddo
|
||||||
|
|
||||||
|
endif
|
||||||
|
hmono = mo_one_e_integrals(m,p)
|
||||||
|
hij = phase*(hdouble + hmono)
|
||||||
|
|
||||||
|
case (0)
|
||||||
|
phase = 1.d0
|
||||||
|
hij = diag_H_mat_elem(key_i,Nint)
|
||||||
|
end select
|
||||||
|
end
|
||||||
|
|
||||||
|
|
||||||
|
subroutine i_H_psi_complex(key,keys,coef,Nint,Ndet,Ndet_max,Nstate,i_H_psi_array)
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Computes $\langle i|H|Psi \rangle = \sum_J c_J \langle i | H | J \rangle$.
|
||||||
|
!
|
||||||
|
! Uses filter_connected_i_H_psi0 to get all the $|J \rangle$ to which $|i \rangle$
|
||||||
|
! is connected.
|
||||||
|
! The i_H_psi_minilist is much faster but requires to build the
|
||||||
|
! minilists.
|
||||||
|
END_DOC
|
||||||
|
integer, intent(in) :: Nint, Ndet,Ndet_max,Nstate
|
||||||
|
integer(bit_kind), intent(in) :: keys(Nint,2,Ndet)
|
||||||
|
integer(bit_kind), intent(in) :: key(Nint,2)
|
||||||
|
double precision, intent(in) :: coef(Ndet_max,Nstate)
|
||||||
|
double precision, intent(out) :: i_H_psi_array(Nstate)
|
||||||
|
|
||||||
|
integer :: i, ii,j
|
||||||
|
double precision :: phase
|
||||||
|
integer :: exc(0:2,2,2)
|
||||||
|
double precision :: hij
|
||||||
|
integer, allocatable :: idx(:)
|
||||||
|
|
||||||
|
ASSERT (Nint > 0)
|
||||||
|
ASSERT (N_int == Nint)
|
||||||
|
ASSERT (Nstate > 0)
|
||||||
|
ASSERT (Ndet > 0)
|
||||||
|
ASSERT (Ndet_max >= Ndet)
|
||||||
|
allocate(idx(0:Ndet))
|
||||||
|
|
||||||
|
i_H_psi_array = 0.d0
|
||||||
|
|
||||||
|
call filter_connected_i_H_psi0(keys,key,Nint,Ndet,idx)
|
||||||
|
if (Nstate == 1) then
|
||||||
|
|
||||||
|
do ii=1,idx(0)
|
||||||
|
i = idx(ii)
|
||||||
|
!DIR$ FORCEINLINE
|
||||||
|
call i_H_j(keys(1,1,i),key,Nint,hij)
|
||||||
|
i_H_psi_array(1) = i_H_psi_array(1) + coef(i,1)*hij
|
||||||
|
enddo
|
||||||
|
|
||||||
|
else
|
||||||
|
|
||||||
|
do ii=1,idx(0)
|
||||||
|
i = idx(ii)
|
||||||
|
!DIR$ FORCEINLINE
|
||||||
|
call i_H_j(keys(1,1,i),key,Nint,hij)
|
||||||
|
do j = 1, Nstate
|
||||||
|
i_H_psi_array(j) = i_H_psi_array(j) + coef(i,j)*hij
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
endif
|
||||||
|
|
||||||
|
end
|
||||||
|
|
||||||
|
|
||||||
|
subroutine i_H_psi_minilist_complex(key,keys,idx_key,N_minilist,coef,Nint,Ndet,Ndet_max,Nstate,i_H_psi_array)
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
integer, intent(in) :: Nint, Ndet,Ndet_max,Nstate,idx_key(Ndet), N_minilist
|
||||||
|
integer(bit_kind), intent(in) :: keys(Nint,2,Ndet)
|
||||||
|
integer(bit_kind), intent(in) :: key(Nint,2)
|
||||||
|
double precision, intent(in) :: coef(Ndet_max,Nstate)
|
||||||
|
double precision, intent(out) :: i_H_psi_array(Nstate)
|
||||||
|
|
||||||
|
integer :: i, ii,j, i_in_key, i_in_coef
|
||||||
|
double precision :: phase
|
||||||
|
integer :: exc(0:2,2,2)
|
||||||
|
double precision :: hij
|
||||||
|
integer, allocatable :: idx(:)
|
||||||
|
BEGIN_DOC
|
||||||
|
! Computes $\langle i|H|\Psi \rangle = \sum_J c_J \langle i|H|J\rangle$.
|
||||||
|
!
|
||||||
|
! Uses filter_connected_i_H_psi0 to get all the $|J \rangle$ to which $|i \rangle$
|
||||||
|
! is connected. The $|J\rangle$ are searched in short pre-computed lists.
|
||||||
|
END_DOC
|
||||||
|
|
||||||
|
ASSERT (Nint > 0)
|
||||||
|
ASSERT (N_int == Nint)
|
||||||
|
ASSERT (Nstate > 0)
|
||||||
|
ASSERT (Ndet > 0)
|
||||||
|
ASSERT (Ndet_max >= Ndet)
|
||||||
|
allocate(idx(0:Ndet))
|
||||||
|
i_H_psi_array = 0.d0
|
||||||
|
|
||||||
|
call filter_connected_i_H_psi0(keys,key,Nint,N_minilist,idx)
|
||||||
|
if (Nstate == 1) then
|
||||||
|
|
||||||
|
do ii=1,idx(0)
|
||||||
|
i_in_key = idx(ii)
|
||||||
|
i_in_coef = idx_key(idx(ii))
|
||||||
|
!DIR$ FORCEINLINE
|
||||||
|
call i_H_j(keys(1,1,i_in_key),key,Nint,hij)
|
||||||
|
! TODO : Cache misses
|
||||||
|
i_H_psi_array(1) = i_H_psi_array(1) + coef(i_in_coef,1)*hij
|
||||||
|
enddo
|
||||||
|
|
||||||
|
else
|
||||||
|
|
||||||
|
do ii=1,idx(0)
|
||||||
|
i_in_key = idx(ii)
|
||||||
|
i_in_coef = idx_key(idx(ii))
|
||||||
|
!DIR$ FORCEINLINE
|
||||||
|
call i_H_j(keys(1,1,i_in_key),key,Nint,hij)
|
||||||
|
do j = 1, Nstate
|
||||||
|
i_H_psi_array(j) = i_H_psi_array(j) + coef(i_in_coef,j)*hij
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
|
endif
|
||||||
|
|
||||||
|
end
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
subroutine i_H_j_single_spin_complex(key_i,key_j,Nint,spin,hij)
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Returns $\langle i|H|j \rangle$ where $i$ and $j$ are determinants differing by
|
||||||
|
! a single excitation.
|
||||||
|
END_DOC
|
||||||
|
integer, intent(in) :: Nint, spin
|
||||||
|
integer(bit_kind), intent(in) :: key_i(Nint,2), key_j(Nint,2)
|
||||||
|
double precision, intent(out) :: hij
|
||||||
|
|
||||||
|
integer :: exc(0:2,2)
|
||||||
|
double precision :: phase
|
||||||
|
|
||||||
|
PROVIDE big_array_exchange_integrals mo_two_e_integrals_in_map
|
||||||
|
|
||||||
|
call get_single_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint)
|
||||||
|
call get_single_excitation_from_fock(key_i,key_j,exc(1,1),exc(1,2),spin,phase,hij)
|
||||||
|
end
|
||||||
|
|
||||||
|
subroutine i_H_j_double_spin_complex(key_i,key_j,Nint,hij)
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Returns $\langle i|H|j \rangle$ where $i$ and $j$ are determinants differing by
|
||||||
|
! a same-spin double excitation.
|
||||||
|
END_DOC
|
||||||
|
integer, intent(in) :: Nint
|
||||||
|
integer(bit_kind), intent(in) :: key_i(Nint), key_j(Nint)
|
||||||
|
double precision, intent(out) :: hij
|
||||||
|
|
||||||
|
integer :: exc(0:2,2)
|
||||||
|
double precision :: phase
|
||||||
|
double precision, external :: get_two_e_integral
|
||||||
|
|
||||||
|
PROVIDE big_array_exchange_integrals mo_two_e_integrals_in_map
|
||||||
|
call get_double_excitation_spin(key_i,key_j,exc,phase,Nint)
|
||||||
|
hij = phase*(get_two_e_integral( &
|
||||||
|
exc(1,1), &
|
||||||
|
exc(2,1), &
|
||||||
|
exc(1,2), &
|
||||||
|
exc(2,2), mo_integrals_map) - &
|
||||||
|
get_two_e_integral( &
|
||||||
|
exc(1,1), &
|
||||||
|
exc(2,1), &
|
||||||
|
exc(2,2), &
|
||||||
|
exc(1,2), mo_integrals_map) )
|
||||||
|
end
|
||||||
|
|
||||||
|
subroutine i_H_j_double_alpha_beta_complex(key_i,key_j,Nint,hij)
|
||||||
|
!todo: modify/implement for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Returns $\langle i|H|j \rangle$ where $i$ and $j$ are determinants differing by
|
||||||
|
! an opposite-spin double excitation.
|
||||||
|
END_DOC
|
||||||
|
integer, intent(in) :: Nint
|
||||||
|
integer(bit_kind), intent(in) :: key_i(Nint,2), key_j(Nint,2)
|
||||||
|
double precision, intent(out) :: hij
|
||||||
|
|
||||||
|
integer :: exc(0:2,2,2)
|
||||||
|
double precision :: phase, phase2
|
||||||
|
double precision, external :: get_two_e_integral
|
||||||
|
|
||||||
|
PROVIDE big_array_exchange_integrals mo_two_e_integrals_in_map
|
||||||
|
|
||||||
|
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
|
||||||
|
hij = phase * big_array_exchange_integrals(exc(1,1,1),exc(1,1,2),exc(1,2,1))
|
||||||
|
else if (exc(1,2,1) == exc(1,1,2)) then
|
||||||
|
hij = phase * big_array_exchange_integrals(exc(1,2,1),exc(1,1,1),exc(1,2,2))
|
||||||
|
else
|
||||||
|
hij = phase*get_two_e_integral( &
|
||||||
|
exc(1,1,1), &
|
||||||
|
exc(1,1,2), &
|
||||||
|
exc(1,2,1), &
|
||||||
|
exc(1,2,2) ,mo_integrals_map)
|
||||||
|
endif
|
||||||
|
end
|
||||||
|
@ -10,6 +10,7 @@ spindeterminants
|
|||||||
psi_coef_matrix_rows integer (spindeterminants_n_det)
|
psi_coef_matrix_rows integer (spindeterminants_n_det)
|
||||||
psi_coef_matrix_columns integer (spindeterminants_n_det)
|
psi_coef_matrix_columns integer (spindeterminants_n_det)
|
||||||
psi_coef_matrix_values double precision (spindeterminants_n_det,spindeterminants_n_states)
|
psi_coef_matrix_values double precision (spindeterminants_n_det,spindeterminants_n_states)
|
||||||
|
psi_coef_matrix_values_complex double precision (2,spindeterminants_n_det,spindeterminants_n_states)
|
||||||
n_svd_coefs integer
|
n_svd_coefs integer
|
||||||
psi_svd_alpha double precision (spindeterminants_n_det_alpha,spindeterminants_n_svd_coefs,spindeterminants_n_states)
|
psi_svd_alpha double precision (spindeterminants_n_det_alpha,spindeterminants_n_svd_coefs,spindeterminants_n_states)
|
||||||
psi_svd_beta double precision (spindeterminants_n_det_beta,spindeterminants_n_svd_coefs,spindeterminants_n_states)
|
psi_svd_beta double precision (spindeterminants_n_det_beta,spindeterminants_n_svd_coefs,spindeterminants_n_states)
|
||||||
|
@ -307,8 +307,12 @@ integer function get_index_in_psi_det_beta_unique(key,Nint)
|
|||||||
|
|
||||||
end
|
end
|
||||||
|
|
||||||
|
|
||||||
subroutine write_spindeterminants
|
subroutine write_spindeterminants
|
||||||
|
!todo: modify for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
use bitmasks
|
use bitmasks
|
||||||
implicit none
|
implicit none
|
||||||
integer(8), allocatable :: tmpdet(:,:)
|
integer(8), allocatable :: tmpdet(:,:)
|
||||||
@ -350,7 +354,11 @@ subroutine write_spindeterminants
|
|||||||
call ezfio_set_spindeterminants_psi_det_beta(psi_det_beta_unique)
|
call ezfio_set_spindeterminants_psi_det_beta(psi_det_beta_unique)
|
||||||
deallocate(tmpdet)
|
deallocate(tmpdet)
|
||||||
|
|
||||||
|
if (is_complex) then
|
||||||
|
call ezfio_set_spindeterminants_psi_coef_matrix_values_complex(psi_bilinear_matrix_values_complex)
|
||||||
|
else
|
||||||
call ezfio_set_spindeterminants_psi_coef_matrix_values(psi_bilinear_matrix_values)
|
call ezfio_set_spindeterminants_psi_coef_matrix_values(psi_bilinear_matrix_values)
|
||||||
|
endif
|
||||||
call ezfio_set_spindeterminants_psi_coef_matrix_rows(psi_bilinear_matrix_rows)
|
call ezfio_set_spindeterminants_psi_coef_matrix_rows(psi_bilinear_matrix_rows)
|
||||||
call ezfio_set_spindeterminants_psi_coef_matrix_columns(psi_bilinear_matrix_columns)
|
call ezfio_set_spindeterminants_psi_coef_matrix_columns(psi_bilinear_matrix_columns)
|
||||||
|
|
||||||
@ -370,6 +378,18 @@ end
|
|||||||
|
|
||||||
det_alpha_norm = 0.d0
|
det_alpha_norm = 0.d0
|
||||||
det_beta_norm = 0.d0
|
det_beta_norm = 0.d0
|
||||||
|
if (is_complex) then
|
||||||
|
do k=1,N_det
|
||||||
|
i = psi_bilinear_matrix_rows(k)
|
||||||
|
j = psi_bilinear_matrix_columns(k)
|
||||||
|
f = 0.d0
|
||||||
|
do l=1,N_states
|
||||||
|
f += cdabs(psi_bilinear_matrix_values_complex(k,l)*psi_bilinear_matrix_values_complex(k,l)) * state_average_weight(l)
|
||||||
|
enddo
|
||||||
|
det_alpha_norm(i) += f
|
||||||
|
det_beta_norm(j) += f
|
||||||
|
enddo
|
||||||
|
else
|
||||||
do k=1,N_det
|
do k=1,N_det
|
||||||
i = psi_bilinear_matrix_rows(k)
|
i = psi_bilinear_matrix_rows(k)
|
||||||
j = psi_bilinear_matrix_columns(k)
|
j = psi_bilinear_matrix_columns(k)
|
||||||
@ -380,6 +400,7 @@ end
|
|||||||
det_alpha_norm(i) += f
|
det_alpha_norm(i) += f
|
||||||
det_beta_norm(j) += f
|
det_beta_norm(j) += f
|
||||||
enddo
|
enddo
|
||||||
|
endif
|
||||||
det_alpha_norm = det_alpha_norm
|
det_alpha_norm = det_alpha_norm
|
||||||
det_beta_norm = det_beta_norm
|
det_beta_norm = det_beta_norm
|
||||||
|
|
||||||
@ -392,8 +413,35 @@ END_PROVIDER
|
|||||||
! !
|
! !
|
||||||
!==============================================================================!
|
!==============================================================================!
|
||||||
|
|
||||||
BEGIN_PROVIDER [ double precision, psi_bilinear_matrix_values, (N_det,N_states) ]
|
BEGIN_PROVIDER [ double precision, psi_bilinear_matrix_values, (N_det,N_states) ]
|
||||||
&BEGIN_PROVIDER [ integer, psi_bilinear_matrix_rows , (N_det) ]
|
use bitmasks
|
||||||
|
PROVIDE psi_bilinear_matrix_rows
|
||||||
|
do k=1,N_det
|
||||||
|
do l=1,N_states
|
||||||
|
psi_bilinear_matrix_values(k,l) = psi_coef(k,l)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
do l=1,N_states
|
||||||
|
call dset_order(psi_bilinear_matrix_values(1,l),psi_bilinear_matrix_order,N_det)
|
||||||
|
enddo
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, psi_bilinear_matrix_values_complex, (N_det,N_states) ]
|
||||||
|
use bitmasks
|
||||||
|
PROVIDE psi_bilinear_matrix_rows
|
||||||
|
do k=1,N_det
|
||||||
|
do l=1,N_states
|
||||||
|
psi_bilinear_matrix_values_complex(k,l) = psi_coef_complex(k,l)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
do l=1,N_states
|
||||||
|
call cdset_order(psi_bilinear_matrix_values_complex(1,l),psi_bilinear_matrix_order,N_det)
|
||||||
|
enddo
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ integer, psi_bilinear_matrix_rows , (N_det) ]
|
||||||
&BEGIN_PROVIDER [ integer, psi_bilinear_matrix_columns, (N_det) ]
|
&BEGIN_PROVIDER [ integer, psi_bilinear_matrix_columns, (N_det) ]
|
||||||
&BEGIN_PROVIDER [ integer, psi_bilinear_matrix_order , (N_det) ]
|
&BEGIN_PROVIDER [ integer, psi_bilinear_matrix_order , (N_det) ]
|
||||||
use bitmasks
|
use bitmasks
|
||||||
@ -408,10 +456,13 @@ END_PROVIDER
|
|||||||
END_DOC
|
END_DOC
|
||||||
integer :: i,j,k, l
|
integer :: i,j,k, l
|
||||||
integer(bit_kind) :: tmp_det(N_int,2)
|
integer(bit_kind) :: tmp_det(N_int,2)
|
||||||
integer, external :: get_index_in_psi_det_sorted_bit
|
! integer, external :: get_index_in_psi_det_sorted_bit
|
||||||
|
|
||||||
|
if (is_complex) then
|
||||||
PROVIDE psi_coef_sorted_bit
|
PROVIDE psi_coef_sorted_bit_complex
|
||||||
|
else
|
||||||
|
PROVIDE psi_coef_sorted_bit
|
||||||
|
endif
|
||||||
|
|
||||||
integer*8, allocatable :: to_sort(:)
|
integer*8, allocatable :: to_sort(:)
|
||||||
integer, external :: get_index_in_psi_det_alpha_unique
|
integer, external :: get_index_in_psi_det_alpha_unique
|
||||||
@ -427,9 +478,6 @@ END_PROVIDER
|
|||||||
ASSERT (j>0)
|
ASSERT (j>0)
|
||||||
ASSERT (j<=N_det_beta_unique)
|
ASSERT (j<=N_det_beta_unique)
|
||||||
|
|
||||||
do l=1,N_states
|
|
||||||
psi_bilinear_matrix_values(k,l) = psi_coef(k,l)
|
|
||||||
enddo
|
|
||||||
psi_bilinear_matrix_rows(k) = i
|
psi_bilinear_matrix_rows(k) = i
|
||||||
psi_bilinear_matrix_columns(k) = j
|
psi_bilinear_matrix_columns(k) = j
|
||||||
to_sort(k) = int(N_det_alpha_unique,8) * int(j-1,8) + int(i,8)
|
to_sort(k) = int(N_det_alpha_unique,8) * int(j-1,8) + int(i,8)
|
||||||
@ -445,11 +493,6 @@ END_PROVIDER
|
|||||||
!$OMP SINGLE
|
!$OMP SINGLE
|
||||||
call iset_order(psi_bilinear_matrix_columns,psi_bilinear_matrix_order,N_det)
|
call iset_order(psi_bilinear_matrix_columns,psi_bilinear_matrix_order,N_det)
|
||||||
!$OMP END SINGLE
|
!$OMP END SINGLE
|
||||||
!$OMP DO
|
|
||||||
do l=1,N_states
|
|
||||||
call dset_order(psi_bilinear_matrix_values(1,l),psi_bilinear_matrix_order,N_det)
|
|
||||||
enddo
|
|
||||||
!$OMP END DO
|
|
||||||
!$OMP END PARALLEL
|
!$OMP END PARALLEL
|
||||||
deallocate(to_sort)
|
deallocate(to_sort)
|
||||||
ASSERT (minval(psi_bilinear_matrix_rows) == 1)
|
ASSERT (minval(psi_bilinear_matrix_rows) == 1)
|
||||||
@ -514,8 +557,71 @@ BEGIN_PROVIDER [ integer, psi_bilinear_matrix_columns_loc, (N_det_beta_unique+1)
|
|||||||
|
|
||||||
END_PROVIDER
|
END_PROVIDER
|
||||||
|
|
||||||
BEGIN_PROVIDER [ double precision, psi_bilinear_matrix_transp_values, (N_det,N_states) ]
|
BEGIN_PROVIDER [ double precision, psi_bilinear_matrix_transp_values, (N_det,N_states) ]
|
||||||
&BEGIN_PROVIDER [ integer, psi_bilinear_matrix_transp_rows , (N_det) ]
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Transpose of :c:data:`psi_bilinear_matrix`
|
||||||
|
!
|
||||||
|
! $D_\beta^\dagger.C^\dagger.D_\alpha$
|
||||||
|
!
|
||||||
|
! Rows are $\alpha$ determinants and columns are $\beta$, but the matrix is stored in row major
|
||||||
|
! format.
|
||||||
|
END_DOC
|
||||||
|
integer :: i,j,k,l
|
||||||
|
|
||||||
|
PROVIDE psi_bilinear_matrix_transp_rows
|
||||||
|
|
||||||
|
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(i,j,k,l)
|
||||||
|
do l=1,N_states
|
||||||
|
!$OMP DO
|
||||||
|
do k=1,N_det
|
||||||
|
psi_bilinear_matrix_transp_values (k,l) = psi_bilinear_matrix_values (k,l)
|
||||||
|
enddo
|
||||||
|
!$OMP ENDDO NOWAIT
|
||||||
|
enddo
|
||||||
|
!$OMP END PARALLEL
|
||||||
|
!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(l)
|
||||||
|
do l=1,N_states
|
||||||
|
call dset_order(psi_bilinear_matrix_transp_values(1,l),psi_bilinear_matrix_transp_order,N_det)
|
||||||
|
enddo
|
||||||
|
!$OMP END PARALLEL DO
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, psi_bilinear_matrix_transp_values_complex, (N_det,N_states) ]
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Transpose of :c:data:`psi_bilinear_matrix`
|
||||||
|
!
|
||||||
|
! $D_\beta^\dagger.C^\dagger.D_\alpha$
|
||||||
|
!
|
||||||
|
! Rows are $\alpha$ determinants and columns are $\beta$, but the matrix is stored in row major
|
||||||
|
! format.
|
||||||
|
END_DOC
|
||||||
|
integer :: i,j,k,l
|
||||||
|
|
||||||
|
PROVIDE psi_bilinear_matrix_transp_rows
|
||||||
|
|
||||||
|
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(i,j,k,l)
|
||||||
|
do l=1,N_states
|
||||||
|
!$OMP DO
|
||||||
|
do k=1,N_det
|
||||||
|
psi_bilinear_matrix_transp_values_complex (k,l) = psi_bilinear_matrix_values_complex (k,l)
|
||||||
|
enddo
|
||||||
|
!$OMP ENDDO NOWAIT
|
||||||
|
enddo
|
||||||
|
!$OMP END PARALLEL
|
||||||
|
!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(l)
|
||||||
|
do l=1,N_states
|
||||||
|
call cdset_order(psi_bilinear_matrix_transp_values_complex(1,l),psi_bilinear_matrix_transp_order,N_det)
|
||||||
|
enddo
|
||||||
|
!$OMP END PARALLEL DO
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ integer, psi_bilinear_matrix_transp_rows , (N_det) ]
|
||||||
&BEGIN_PROVIDER [ integer, psi_bilinear_matrix_transp_columns, (N_det) ]
|
&BEGIN_PROVIDER [ integer, psi_bilinear_matrix_transp_columns, (N_det) ]
|
||||||
&BEGIN_PROVIDER [ integer, psi_bilinear_matrix_transp_order , (N_det) ]
|
&BEGIN_PROVIDER [ integer, psi_bilinear_matrix_transp_order , (N_det) ]
|
||||||
use bitmasks
|
use bitmasks
|
||||||
@ -530,18 +636,15 @@ END_PROVIDER
|
|||||||
END_DOC
|
END_DOC
|
||||||
integer :: i,j,k,l
|
integer :: i,j,k,l
|
||||||
|
|
||||||
PROVIDE psi_coef_sorted_bit
|
if (is_complex) then
|
||||||
|
PROVIDE psi_coef_sorted_bit_complex
|
||||||
|
else
|
||||||
|
PROVIDE psi_coef_sorted_bit
|
||||||
|
endif
|
||||||
|
|
||||||
integer*8, allocatable :: to_sort(:)
|
integer*8, allocatable :: to_sort(:)
|
||||||
allocate(to_sort(N_det))
|
allocate(to_sort(N_det))
|
||||||
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(i,j,k,l)
|
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(i,j,k,l)
|
||||||
do l=1,N_states
|
|
||||||
!$OMP DO
|
|
||||||
do k=1,N_det
|
|
||||||
psi_bilinear_matrix_transp_values (k,l) = psi_bilinear_matrix_values (k,l)
|
|
||||||
enddo
|
|
||||||
!$OMP ENDDO NOWAIT
|
|
||||||
enddo
|
|
||||||
!$OMP DO
|
!$OMP DO
|
||||||
do k=1,N_det
|
do k=1,N_det
|
||||||
psi_bilinear_matrix_transp_columns(k) = psi_bilinear_matrix_columns(k)
|
psi_bilinear_matrix_transp_columns(k) = psi_bilinear_matrix_columns(k)
|
||||||
@ -563,11 +666,6 @@ END_PROVIDER
|
|||||||
call i8radix_sort(to_sort, psi_bilinear_matrix_transp_order, N_det,-1)
|
call i8radix_sort(to_sort, psi_bilinear_matrix_transp_order, N_det,-1)
|
||||||
call iset_order(psi_bilinear_matrix_transp_rows,psi_bilinear_matrix_transp_order,N_det)
|
call iset_order(psi_bilinear_matrix_transp_rows,psi_bilinear_matrix_transp_order,N_det)
|
||||||
call iset_order(psi_bilinear_matrix_transp_columns,psi_bilinear_matrix_transp_order,N_det)
|
call iset_order(psi_bilinear_matrix_transp_columns,psi_bilinear_matrix_transp_order,N_det)
|
||||||
!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(l)
|
|
||||||
do l=1,N_states
|
|
||||||
call dset_order(psi_bilinear_matrix_transp_values(1,l),psi_bilinear_matrix_transp_order,N_det)
|
|
||||||
enddo
|
|
||||||
!$OMP END PARALLEL DO
|
|
||||||
deallocate(to_sort)
|
deallocate(to_sort)
|
||||||
ASSERT (minval(psi_bilinear_matrix_transp_columns) == 1)
|
ASSERT (minval(psi_bilinear_matrix_transp_columns) == 1)
|
||||||
ASSERT (minval(psi_bilinear_matrix_transp_rows) == 1)
|
ASSERT (minval(psi_bilinear_matrix_transp_rows) == 1)
|
||||||
@ -641,7 +739,30 @@ BEGIN_PROVIDER [ double precision, psi_bilinear_matrix, (N_det_alpha_unique,N_de
|
|||||||
enddo
|
enddo
|
||||||
END_PROVIDER
|
END_PROVIDER
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ complex*16, psi_bilinear_matrix_complex, (N_det_alpha_unique,N_det_beta_unique,N_states) ]
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Coefficient matrix if the wave function is expressed in a bilinear form :
|
||||||
|
!
|
||||||
|
! $D_\alpha^\dagger.C.D_\beta$
|
||||||
|
END_DOC
|
||||||
|
integer :: i,j,k,istate
|
||||||
|
psi_bilinear_matrix_complex = (0.d0,0.d0)
|
||||||
|
do k=1,N_det
|
||||||
|
i = psi_bilinear_matrix_rows(k)
|
||||||
|
j = psi_bilinear_matrix_columns(k)
|
||||||
|
do istate=1,N_states
|
||||||
|
psi_bilinear_matrix_complex(i,j,istate) = psi_bilinear_matrix_values_complex(k,istate)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
subroutine create_wf_of_psi_bilinear_matrix(truncate)
|
subroutine create_wf_of_psi_bilinear_matrix(truncate)
|
||||||
|
!todo: modify for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
use bitmasks
|
use bitmasks
|
||||||
implicit none
|
implicit none
|
||||||
BEGIN_DOC
|
BEGIN_DOC
|
||||||
@ -713,6 +834,11 @@ subroutine create_wf_of_psi_bilinear_matrix(truncate)
|
|||||||
end
|
end
|
||||||
|
|
||||||
subroutine generate_all_alpha_beta_det_products
|
subroutine generate_all_alpha_beta_det_products
|
||||||
|
!todo: modify for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
implicit none
|
implicit none
|
||||||
BEGIN_DOC
|
BEGIN_DOC
|
||||||
! Creates a wave function from all possible $\alpha \times \beta$ determinants
|
! Creates a wave function from all possible $\alpha \times \beta$ determinants
|
||||||
@ -856,6 +982,11 @@ end
|
|||||||
|
|
||||||
|
|
||||||
subroutine copy_psi_bilinear_to_psi(psi, isize)
|
subroutine copy_psi_bilinear_to_psi(psi, isize)
|
||||||
|
!todo: modify for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
implicit none
|
implicit none
|
||||||
BEGIN_DOC
|
BEGIN_DOC
|
||||||
! Overwrites :c:data:`psi_det` and :c:data:`psi_coef` with the wave function
|
! Overwrites :c:data:`psi_det` and :c:data:`psi_coef` with the wave function
|
||||||
@ -1292,6 +1423,11 @@ END_TEMPLATE
|
|||||||
|
|
||||||
|
|
||||||
subroutine wf_of_psi_bilinear_matrix(truncate)
|
subroutine wf_of_psi_bilinear_matrix(truncate)
|
||||||
|
!todo: modify for complex
|
||||||
|
if (is_complex) then
|
||||||
|
print*,irp_here,' not implemented for complex'
|
||||||
|
stop -1
|
||||||
|
endif
|
||||||
use bitmasks
|
use bitmasks
|
||||||
implicit none
|
implicit none
|
||||||
BEGIN_DOC
|
BEGIN_DOC
|
||||||
|
@ -409,6 +409,7 @@ BEGIN_TEMPLATE
|
|||||||
SUBST [ X, type ]
|
SUBST [ X, type ]
|
||||||
; real ;;
|
; real ;;
|
||||||
d ; double precision ;;
|
d ; double precision ;;
|
||||||
|
cd; complex*16 ;;
|
||||||
i ; integer ;;
|
i ; integer ;;
|
||||||
i8; integer*8 ;;
|
i8; integer*8 ;;
|
||||||
i2; integer*2 ;;
|
i2; integer*2 ;;
|
||||||
|
@ -3,14 +3,14 @@
|
|||||||
current:
|
current:
|
||||||
|
|
||||||
general:
|
general:
|
||||||
i_h_j_complex
|
check for dependence on psi_det_sorted, clean up providers
|
||||||
diag_h_mat_elem if is_complex
|
|
||||||
|
|
||||||
|
|
||||||
determinants:
|
determinants:
|
||||||
(done) connected_to_ref.irp.f
|
(done) connected_to_ref.irp.f
|
||||||
(done) create_excitations.irp.f
|
(done) create_excitations.irp.f
|
||||||
(****) density_matrix.irp.f
|
(done?)density_matrix{,_complex}.irp.f
|
||||||
|
no one_e_dm_dagger_mo_spin_index_complex
|
||||||
|
need to test for complex
|
||||||
(done) determinants_bitmasks.irp.f
|
(done) determinants_bitmasks.irp.f
|
||||||
(****) determinants{_complex}.irp.f
|
(****) determinants{_complex}.irp.f
|
||||||
mostly done
|
mostly done
|
||||||
@ -41,15 +41,24 @@ determinants:
|
|||||||
(done) psi_cas{,_complex}.irp.f
|
(done) psi_cas{,_complex}.irp.f
|
||||||
might be able to combine some providers??
|
might be able to combine some providers??
|
||||||
(done) psi_energy_mono_elec.irp.f
|
(done) psi_energy_mono_elec.irp.f
|
||||||
(****) ref_bitmask.irp.f
|
(done) ref_bitmask.irp.f
|
||||||
(****) s2{,_complex}.irp.f
|
(****) s2{,_complex}.irp.f
|
||||||
(****) single_excitations.irp.f
|
made copies of needed functions for complex
|
||||||
|
still need to do implementation
|
||||||
|
(done) single_excitations.irp.f
|
||||||
(****) single_excitation_two_e.irp.f
|
(****) single_excitation_two_e.irp.f
|
||||||
(****) slater_rules.irp.f
|
(****) slater_rules.irp.f
|
||||||
|
made copies of needed functions for complex
|
||||||
|
still need to do implementation
|
||||||
(****) slater_rules_wee_mono.irp.f
|
(****) slater_rules_wee_mono.irp.f
|
||||||
(done) sort_dets_ab.irp.f
|
(done) sort_dets_ab.irp.f
|
||||||
spindeterminants.ezfio_config
|
spindeterminants.ezfio_config
|
||||||
(****) spindeterminants.irp.f
|
need svd complex?
|
||||||
|
(done?) spindeterminants.irp.f
|
||||||
|
separated psi_bilinear_matrix_values from psi_bilinear_matrix_{rows,columns,order}
|
||||||
|
new provider for psi_bilinear_matrix_values_complex
|
||||||
|
same for bilinear matrix transp (no conjugate)
|
||||||
|
done except for specific functions that are commented with todo
|
||||||
(****) two_e_density_matrix.irp.pouet
|
(****) two_e_density_matrix.irp.pouet
|
||||||
(done) utils.irp.f
|
(done) utils.irp.f
|
||||||
(****) zmq.irp.f
|
(****) zmq.irp.f
|
||||||
|
Loading…
Reference in New Issue
Block a user