1
0
mirror of https://gitlab.com/scemama/qp_plugins_scemama.git synced 2024-06-02 03:15:25 +02:00

Complete FCI

This commit is contained in:
Anthony Scemama 2021-06-07 15:15:01 +02:00
parent 6333429d20
commit 53f956e204
6 changed files with 98 additions and 12 deletions

View File

@ -0,0 +1,6 @@
[energy]
type: double precision
doc: Calculated Selected |FCI| energy
interface: ezfio
size: (determinants.n_states)

3
devel/fci_complete/NEED Normal file
View File

@ -0,0 +1,3 @@
davidson_undressed
hartree_fock
determinants

View File

@ -0,0 +1,4 @@
===
mpn
===

View File

@ -0,0 +1,8 @@
program fci_complete
call generate_fci_space
call diagonalize_ci
end

View File

@ -0,0 +1,75 @@
subroutine generate_fci_space
use bitmasks
implicit none
BEGIN_DOC
! Generates the complete FCI space
END_DOC
integer :: i, sze, ncore
integer(bit_kind) :: o(N_int,2)
integer(bit_kind) :: u, coremask
if (mo_num > 64) then
stop 'No more than 64 MOs'
endif
ncore = 0
coremask = 0_bit_kind
do i=1,mo_num
if (trim(mo_class(i)) == 'Core') then
ncore += 1
coremask = ibset(coremask,i-1)
endif
enddo
o(1,1) = iand(full_ijkl_bitmask(1),not(coremask))
o(1,2) = 0_bit_kind
call configuration_to_dets_size(o,n_det_alpha_unique,elec_alpha_num-ncore,N_int)
TOUCH n_det_alpha_unique
integer :: k,n,m, t, t1, t2
k=0
n = elec_alpha_num
m = mo_num - n
n = n
u = shiftl(1_bit_kind,n) -1
do while (u < shiftl(1_bit_kind,n+m))
if (iand(coremask, u) == coremask) then
k = k+1
psi_det_alpha_unique(1,k) = u
endif
t = ior(u,u-1)
t1 = t+1
IRP_IF WITHOUT_TRAILZ
t2 = shiftr((iand(not(t),t1)-1), popcnt(ieor(u,u-1)))
IRP_ELSE
t2 = shiftr((iand(not(t),t1)-1), trailz(u)+1)
IRP_ENDIF
u = ior(t1,t2)
enddo
call configuration_to_dets_size(o,n_det_beta_unique,elec_beta_num-ncore,N_int)
TOUCH n_det_beta_unique
k=0
n = elec_beta_num
m = mo_num - n
u = shiftl(1_bit_kind,n) -1
do while (u < shiftl(1_bit_kind,n+m))
if (iand(coremask, u) == coremask) then
k = k+1
psi_det_beta_unique(1,k) = u
endif
t = ior(u,u-1)
t1 = t+1
t2 = shiftr((iand(not(t),t1)-1), trailz(u)+1)
u = ior(t1,t2)
enddo
call generate_all_alpha_beta_det_products
print *, 'Ndet = ', N_det
end

View File

@ -227,18 +227,6 @@ subroutine run
Ept2 += ctmp*ctmp / (E0 - (H_diag(i,j)+nuclear_repulsion) )
enddo
enddo
! do j=n_selected+1,n_det_beta_unique
! do i=n_selected+1,n_det_alpha_unique
! ctmp = 0.d0
! do l=1,n_selected
! do k=1,n_selected
! ctmp = ctmp + H(k,l,i,j) * psi_postsvd(k,l)
! enddo
! enddo
! psi_postsvd(i,j) = ctmp / (E0 - (H(i,j,i,j)+nuclear_repulsion) )
! Ept2 += ctmp*ctmp / (E0 - (H(i,j,i,j)+nuclear_repulsion) )
! enddo
! enddo
tol_energy = dabs(E_prev - E0)
print '(I5, 3(3X, F20.10))', it_svd, E0, E0 + Ept2, tol_energy
@ -403,6 +391,7 @@ subroutine const_H_uv(Uref, Vref, H, H_diag, n_selected)
if ( degree > 2) cycle
call i_H_j(det1, det2, N_int, h12)
if (h12 == 0.d0) cycle
do m=1,nb
tmp3(m,j,l) = tmp3(m,j,l) + h12 * tmp1(m,i) * tmp1(m,k)
@ -410,6 +399,7 @@ subroutine const_H_uv(Uref, Vref, H, H_diag, n_selected)
do n=1,n_selected
c_tmp = h12 * Vref(j,n)
if (c_tmp == 0.d0) cycle
do m=1,n_selected
H0(k,l,m,n) = H0(k,l,m,n) + c_tmp * tmp1(m,i)
enddo