mirror of
https://github.com/LCPQ/quantum_package
synced 2024-08-25 05:51:46 +02:00
Cleaned MRCC
This commit is contained in:
parent
2d1f40cae7
commit
2a2e099bca
228
plugins/MRCC_Utils/amplitudes.irp.f
Normal file
228
plugins/MRCC_Utils/amplitudes.irp.f
Normal file
@ -0,0 +1,228 @@
|
|||||||
|
BEGIN_PROVIDER [ integer, n_exc_active ]
|
||||||
|
&BEGIN_PROVIDER [ integer, active_pp_idx, (hh_nex) ]
|
||||||
|
&BEGIN_PROVIDER [ integer, active_hh_idx, (hh_nex) ]
|
||||||
|
&BEGIN_PROVIDER [ logical, is_active_exc, (hh_nex) ]
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! is_active_exc : True if the excitation involves at least one active MO
|
||||||
|
!
|
||||||
|
! n_exc_active : Number of active excitations : Number of excitations without the inactive ones.
|
||||||
|
!
|
||||||
|
! active_hh_idx :
|
||||||
|
!
|
||||||
|
! active_pp_idx :
|
||||||
|
END_DOC
|
||||||
|
integer :: hh, pp, II
|
||||||
|
integer :: ind
|
||||||
|
logical :: ok
|
||||||
|
integer(bit_kind) :: myDet(N_int, 2), myMask(N_int, 2)
|
||||||
|
|
||||||
|
integer, allocatable :: pathTo(:)
|
||||||
|
integer, external :: searchDet
|
||||||
|
|
||||||
|
allocate(pathTo(N_det_non_ref))
|
||||||
|
|
||||||
|
pathTo(:) = 0
|
||||||
|
is_active_exc(:) = .false.
|
||||||
|
n_exc_active = 0
|
||||||
|
|
||||||
|
do hh = 1, hh_shortcut(0)
|
||||||
|
do pp = hh_shortcut(hh), hh_shortcut(hh+1)-1
|
||||||
|
do II = 1, N_det_ref
|
||||||
|
|
||||||
|
call apply_hole_local(psi_ref(1,1,II), hh_exists(1, hh), myMask, ok, N_int)
|
||||||
|
if(.not. ok) cycle
|
||||||
|
|
||||||
|
call apply_particle_local(myMask, pp_exists(1, pp), myDet, ok, N_int)
|
||||||
|
if(.not. ok) cycle
|
||||||
|
|
||||||
|
ind = searchDet(psi_non_ref_sorted(1,1,1), myDet(1,1), N_det_non_ref, N_int)
|
||||||
|
if(ind == -1) cycle
|
||||||
|
|
||||||
|
ind = psi_non_ref_sorted_idx(ind)
|
||||||
|
if(pathTo(ind) == 0) then
|
||||||
|
pathTo(ind) = pp
|
||||||
|
else
|
||||||
|
is_active_exc(pp) = .true.
|
||||||
|
is_active_exc(pathTo(ind)) = .true.
|
||||||
|
end if
|
||||||
|
end do
|
||||||
|
end do
|
||||||
|
end do
|
||||||
|
|
||||||
|
do hh = 1, hh_shortcut(0)
|
||||||
|
do pp = hh_shortcut(hh), hh_shortcut(hh+1)-1
|
||||||
|
if(is_active_exc(pp)) then
|
||||||
|
n_exc_active = n_exc_active + 1
|
||||||
|
active_hh_idx(n_exc_active) = hh
|
||||||
|
active_pp_idx(n_exc_active) = pp
|
||||||
|
end if
|
||||||
|
end do
|
||||||
|
end do
|
||||||
|
|
||||||
|
deallocate(pathTo)
|
||||||
|
|
||||||
|
print *, n_exc_active, "inactive excitations /", hh_nex
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ integer, active_excitation_to_determinants_idx, (0:N_det_ref+1, n_exc_active) ]
|
||||||
|
&BEGIN_PROVIDER [ double precision, active_excitation_to_determinants_val, (N_states,N_det_ref+1, n_exc_active) ]
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Sparse matrix A containing the matrix to transform the active excitations to
|
||||||
|
! determinants : A | \Psi_0 > = | \Psi_SD >
|
||||||
|
END_DOC
|
||||||
|
integer :: s, ppp, pp, hh, II, ind, wk, i
|
||||||
|
integer, allocatable :: lref(:)
|
||||||
|
integer(bit_kind) :: myDet(N_int,2), myMask(N_int,2)
|
||||||
|
double precision :: phase
|
||||||
|
logical :: ok
|
||||||
|
integer, external :: searchDet
|
||||||
|
|
||||||
|
|
||||||
|
!$OMP PARALLEL default(none) shared(psi_non_ref, hh_exists, pp_exists, N_int,&
|
||||||
|
!$OMP active_excitation_to_determinants_val, active_excitation_to_determinants_idx)&
|
||||||
|
!$OMP shared(hh_shortcut, psi_ref_coef, N_det_non_ref, psi_non_ref_sorted, &
|
||||||
|
!$OMP psi_non_ref_sorted_idx, psi_ref, N_det_ref, N_states)&
|
||||||
|
!$OMP shared(is_active_exc, active_hh_idx, active_pp_idx, n_exc_active)&
|
||||||
|
!$OMP private(lref, pp, II, ok, myMask, myDet, ind, phase, wk, ppp, hh, s)
|
||||||
|
allocate(lref(N_det_non_ref))
|
||||||
|
!$OMP DO schedule(static,10)
|
||||||
|
do ppp=1,n_exc_active
|
||||||
|
active_excitation_to_determinants_val(:,:,ppp) = 0d0
|
||||||
|
active_excitation_to_determinants_idx(:,ppp) = 0
|
||||||
|
pp = active_pp_idx(ppp)
|
||||||
|
hh = active_hh_idx(ppp)
|
||||||
|
lref = 0
|
||||||
|
do II = 1, N_det_ref
|
||||||
|
call apply_hole_local(psi_ref(1,1,II), hh_exists(1, hh), myMask, ok, N_int)
|
||||||
|
if(.not. ok) cycle
|
||||||
|
call apply_particle_local(myMask, pp_exists(1, pp), myDet, ok, N_int)
|
||||||
|
if(.not. ok) cycle
|
||||||
|
ind = searchDet(psi_non_ref_sorted(1,1,1), myDet(1,1), N_det_non_ref, N_int)
|
||||||
|
if(ind /= -1) then
|
||||||
|
call get_phase(myDet(1,1), psi_ref(1,1,II), phase, N_int)
|
||||||
|
if (phase > 0.d0) then
|
||||||
|
lref(psi_non_ref_sorted_idx(ind)) = II
|
||||||
|
else
|
||||||
|
lref(psi_non_ref_sorted_idx(ind)) = -II
|
||||||
|
endif
|
||||||
|
end if
|
||||||
|
end do
|
||||||
|
wk = 0
|
||||||
|
do i=1, N_det_non_ref
|
||||||
|
if(lref(i) > 0) then
|
||||||
|
wk += 1
|
||||||
|
do s=1,N_states
|
||||||
|
active_excitation_to_determinants_val(s,wk, ppp) = psi_ref_coef(lref(i), s)
|
||||||
|
enddo
|
||||||
|
active_excitation_to_determinants_idx(wk, ppp) = i
|
||||||
|
else if(lref(i) < 0) then
|
||||||
|
wk += 1
|
||||||
|
do s=1,N_states
|
||||||
|
active_excitation_to_determinants_val(s,wk, ppp) = -psi_ref_coef(-lref(i), s)
|
||||||
|
enddo
|
||||||
|
active_excitation_to_determinants_idx(wk, ppp) = i
|
||||||
|
end if
|
||||||
|
end do
|
||||||
|
active_excitation_to_determinants_idx(0,ppp) = wk
|
||||||
|
end do
|
||||||
|
!$OMP END DO
|
||||||
|
deallocate(lref)
|
||||||
|
!$OMP END PARALLEL
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ integer, mrcc_AtA_ind, (N_det_ref * n_exc_active) ]
|
||||||
|
&BEGIN_PROVIDER [ double precision, mrcc_AtA_val, (N_states, N_det_ref * n_exc_active) ]
|
||||||
|
&BEGIN_PROVIDER [ integer, mrcc_col_shortcut, (n_exc_active) ]
|
||||||
|
&BEGIN_PROVIDER [ integer, mrcc_N_col, (n_exc_active) ]
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! A is active_excitation_to_determinants in At.A
|
||||||
|
END_DOC
|
||||||
|
integer :: AtA_size, i,k
|
||||||
|
integer :: at_roww, at_row, wk, a_coll, a_col, r1, r2, s
|
||||||
|
double precision, allocatable :: t(:), A_val_mwen(:,:)
|
||||||
|
integer, allocatable :: A_ind_mwen(:)
|
||||||
|
|
||||||
|
mrcc_AtA_ind(:) = 0
|
||||||
|
mrcc_AtA_val(:,:) = 0.d0
|
||||||
|
mrcc_col_shortcut(:) = 0
|
||||||
|
mrcc_N_col(:) = 0
|
||||||
|
AtA_size = 0
|
||||||
|
|
||||||
|
|
||||||
|
!$OMP PARALLEL default(none) shared(k, active_excitation_to_determinants_idx,&
|
||||||
|
!$OMP active_excitation_to_determinants_val, hh_nex) &
|
||||||
|
!$OMP private(at_row, a_col, t, i, r1, r2, wk, A_ind_mwen, A_val_mwen, a_coll, at_roww)&
|
||||||
|
!$OMP shared(N_states,mrcc_col_shortcut, mrcc_N_col, AtA_size, mrcc_AtA_val, mrcc_AtA_ind, n_exc_active, active_pp_idx)
|
||||||
|
allocate(A_val_mwen(N_states,hh_nex), A_ind_mwen(hh_nex), t(N_states))
|
||||||
|
|
||||||
|
!$OMP DO schedule(dynamic, 100)
|
||||||
|
do at_roww = 1, n_exc_active ! hh_nex
|
||||||
|
at_row = active_pp_idx(at_roww)
|
||||||
|
wk = 0
|
||||||
|
if(mod(at_roww, 100) == 0) print *, "AtA", at_row, "/", hh_nex
|
||||||
|
|
||||||
|
do a_coll = 1, n_exc_active
|
||||||
|
a_col = active_pp_idx(a_coll)
|
||||||
|
t(:) = 0d0
|
||||||
|
r1 = 1
|
||||||
|
r2 = 1
|
||||||
|
do while ((active_excitation_to_determinants_idx(r1, at_roww) /= 0).and.(active_excitation_to_determinants_idx(r2, a_coll) /= 0))
|
||||||
|
if(active_excitation_to_determinants_idx(r1, at_roww) > active_excitation_to_determinants_idx(r2, a_coll)) then
|
||||||
|
r2 = r2+1
|
||||||
|
else if(active_excitation_to_determinants_idx(r1, at_roww) < active_excitation_to_determinants_idx(r2, a_coll)) then
|
||||||
|
r1 = r1+1
|
||||||
|
else
|
||||||
|
do s=1,N_states
|
||||||
|
t(s) = t(s) - active_excitation_to_determinants_val(s,r1, at_roww) * active_excitation_to_determinants_val(s,r2, a_coll)
|
||||||
|
enddo
|
||||||
|
r1 = r1+1
|
||||||
|
r2 = r2+1
|
||||||
|
end if
|
||||||
|
end do
|
||||||
|
|
||||||
|
if(a_col == at_row) then
|
||||||
|
do s=1,N_states
|
||||||
|
t(s) = t(s) + 1.d0
|
||||||
|
enddo
|
||||||
|
end if
|
||||||
|
if(sum(abs(t)) /= 0.d0) then
|
||||||
|
wk += 1
|
||||||
|
A_ind_mwen(wk) = a_col
|
||||||
|
do s=1,N_states
|
||||||
|
A_val_mwen(s,wk) = t(s)
|
||||||
|
enddo
|
||||||
|
end if
|
||||||
|
end do
|
||||||
|
|
||||||
|
if(wk /= 0) then
|
||||||
|
!$OMP CRITICAL
|
||||||
|
mrcc_col_shortcut(at_roww) = AtA_size+1
|
||||||
|
mrcc_N_col(at_roww) = wk
|
||||||
|
if (AtA_size+wk > size(mrcc_AtA_ind,1)) then
|
||||||
|
print *, AtA_size+wk , size(mrcc_AtA_ind,1)
|
||||||
|
stop 'too small'
|
||||||
|
endif
|
||||||
|
do i=1,wk
|
||||||
|
mrcc_AtA_ind(AtA_size+i) = A_ind_mwen(i)
|
||||||
|
do s=1,N_states
|
||||||
|
mrcc_AtA_val(s,AtA_size+i) = A_val_mwen(s,i)
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
AtA_size += wk
|
||||||
|
!$OMP END CRITICAL
|
||||||
|
end if
|
||||||
|
end do
|
||||||
|
!$OMP END DO NOWAIT
|
||||||
|
deallocate (A_ind_mwen, A_val_mwen, t)
|
||||||
|
!$OMP END PARALLEL
|
||||||
|
|
||||||
|
print *, "ATA SIZE", ata_size
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
@ -614,207 +614,60 @@ END_PROVIDER
|
|||||||
END_PROVIDER
|
END_PROVIDER
|
||||||
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ double precision, dIj_unique, (hh_shortcut(hh_shortcut(0)+1)-1, N_states) ]
|
BEGIN_PROVIDER [ double precision, dIj_unique, (hh_nex, N_states) ]
|
||||||
&BEGIN_PROVIDER [ double precision, rho_mrcc, (N_det_non_ref, N_states) ]
|
&BEGIN_PROVIDER [ double precision, rho_mrcc, (N_det_non_ref, N_states) ]
|
||||||
implicit none
|
implicit none
|
||||||
logical :: ok
|
logical :: ok
|
||||||
integer :: i, j, k, s, II, pp, ppp, hh, ind, wk, nex, a_col, at_row
|
integer :: i, j, k, s, II, pp, ppp, hh, ind, wk, a_col, at_row
|
||||||
integer, external :: searchDet, unsortedSearchDet
|
integer, external :: searchDet, unsortedSearchDet
|
||||||
integer(bit_kind) :: myDet(N_int, 2), myMask(N_int, 2)
|
integer(bit_kind) :: myDet(N_int, 2), myMask(N_int, 2)
|
||||||
integer :: N, INFO, AtA_size, r1, r2
|
integer :: N, INFO, r1, r2
|
||||||
double precision , allocatable :: AtB(:), AtA_val(:), A_val(:,:), x(:), x_new(:), A_val_mwen(:)
|
double precision , allocatable :: AtB(:), x(:), x_new(:), A_val_mwen(:,:), t(:)
|
||||||
double precision :: t, norm, cx, res
|
double precision :: norm, cx, res
|
||||||
integer, allocatable :: A_ind(:,:), lref(:), AtA_ind(:), A_ind_mwen(:), col_shortcut(:), N_col(:)
|
integer, allocatable :: lref(:), A_ind_mwen(:)
|
||||||
double precision :: phase
|
double precision :: phase
|
||||||
|
! double precision , allocatable :: mrcc_AtA_val(:,:)
|
||||||
|
! integer, allocatable :: mrcc_AtA_ind(:), col_shortcut(:), , mrcc_N_col(:)
|
||||||
|
|
||||||
|
|
||||||
integer, allocatable :: pathTo(:), active_hh_idx(:), active_pp_idx(:)
|
|
||||||
logical, allocatable :: active(:)
|
|
||||||
double precision, allocatable :: rho_mrcc_init(:,:)
|
double precision, allocatable :: rho_mrcc_init(:,:)
|
||||||
integer :: nactive
|
integer :: a_coll, at_roww
|
||||||
|
|
||||||
|
print *, "TI", hh_nex, N_det_non_ref
|
||||||
|
|
||||||
|
|
||||||
nex = hh_shortcut(hh_shortcut(0)+1)-1
|
|
||||||
print *, "TI", nex, N_det_non_ref
|
|
||||||
|
|
||||||
allocate(pathTo(N_det_non_ref), active(nex))
|
|
||||||
allocate(active_pp_idx(nex), active_hh_idx(nex))
|
|
||||||
allocate(rho_mrcc_init(N_det_non_ref, N_states))
|
allocate(rho_mrcc_init(N_det_non_ref, N_states))
|
||||||
|
|
||||||
pathTo = 0
|
allocate(x(hh_nex), AtB(hh_nex))
|
||||||
active = .false.
|
x = 0d0
|
||||||
nactive = 0
|
allocate(x_new(hh_nex))
|
||||||
|
|
||||||
|
|
||||||
do hh = 1, hh_shortcut(0)
|
|
||||||
do pp = hh_shortcut(hh), hh_shortcut(hh+1)-1
|
|
||||||
do II = 1, N_det_ref
|
|
||||||
call apply_hole_local(psi_ref(1,1,II), hh_exists(1, hh), myMask, ok, N_int)
|
|
||||||
if(.not. ok) cycle
|
|
||||||
call apply_particle_local(myMask, pp_exists(1, pp), myDet, ok, N_int)
|
|
||||||
if(.not. ok) cycle
|
|
||||||
ind = searchDet(psi_non_ref_sorted(1,1,1), myDet(1,1), N_det_non_ref, N_int)
|
|
||||||
if(ind == -1) cycle
|
|
||||||
ind = psi_non_ref_sorted_idx(ind)
|
|
||||||
if(pathTo(ind) == 0) then
|
|
||||||
pathTo(ind) = pp
|
|
||||||
else
|
|
||||||
active(pp) = .true.
|
|
||||||
active(pathTo(ind)) = .true.
|
|
||||||
end if
|
|
||||||
end do
|
|
||||||
end do
|
|
||||||
end do
|
|
||||||
|
|
||||||
do hh = 1, hh_shortcut(0)
|
|
||||||
do pp = hh_shortcut(hh), hh_shortcut(hh+1)-1
|
|
||||||
if(active(pp)) then
|
|
||||||
nactive = nactive + 1
|
|
||||||
active_hh_idx(nactive) = hh
|
|
||||||
active_pp_idx(nactive) = pp
|
|
||||||
end if
|
|
||||||
end do
|
|
||||||
end do
|
|
||||||
|
|
||||||
print *, nactive, "inact/", size(active)
|
|
||||||
|
|
||||||
allocate(A_ind(0:N_det_ref+1, nactive), A_val(N_det_ref+1, nactive))
|
|
||||||
allocate(AtA_ind(N_det_ref * nactive), AtA_val(N_det_ref * nactive))
|
|
||||||
allocate(x(nex), AtB(nex))
|
|
||||||
allocate(N_col(nactive), col_shortcut(nactive))
|
|
||||||
allocate(x_new(nex))
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
do s=1,N_states
|
do s=1,N_states
|
||||||
|
|
||||||
A_val = 0d0
|
AtB(:) = 0.d0
|
||||||
A_ind = 0
|
!$OMP PARALLEL default(none) shared(k, psi_non_ref_coef, active_excitation_to_determinants_idx,&
|
||||||
AtA_ind = 0
|
!$OMP active_excitation_to_determinants_val, x, N_det_ref, hh_nex, N_det_non_ref) &
|
||||||
AtB = 0d0
|
|
||||||
AtA_val = 0d0
|
|
||||||
x = 0d0
|
|
||||||
N_col = 0
|
|
||||||
col_shortcut = 0
|
|
||||||
|
|
||||||
!$OMP PARALLEL default(none) shared(psi_non_ref, hh_exists, pp_exists, N_int, A_val, A_ind)&
|
|
||||||
!$OMP shared(s, hh_shortcut, psi_ref_coef, N_det_non_ref, psi_non_ref_sorted, psi_non_ref_sorted_idx, psi_ref, N_det_ref)&
|
|
||||||
!$OMP shared(active, active_hh_idx, active_pp_idx, nactive) &
|
|
||||||
!$OMP private(lref, pp, II, ok, myMask, myDet, ind, phase, wk, ppp, hh)
|
|
||||||
allocate(lref(N_det_non_ref))
|
|
||||||
!$OMP DO schedule(static,10)
|
|
||||||
do ppp=1,nactive
|
|
||||||
pp = active_pp_idx(ppp)
|
|
||||||
hh = active_hh_idx(ppp)
|
|
||||||
lref = 0
|
|
||||||
do II = 1, N_det_ref
|
|
||||||
call apply_hole_local(psi_ref(1,1,II), hh_exists(1, hh), myMask, ok, N_int)
|
|
||||||
if(.not. ok) cycle
|
|
||||||
call apply_particle_local(myMask, pp_exists(1, pp), myDet, ok, N_int)
|
|
||||||
if(.not. ok) cycle
|
|
||||||
ind = searchDet(psi_non_ref_sorted(1,1,1), myDet(1,1), N_det_non_ref, N_int)
|
|
||||||
if(ind /= -1) then
|
|
||||||
call get_phase(myDet(1,1), psi_ref(1,1,II), phase, N_int)
|
|
||||||
if (phase > 0.d0) then
|
|
||||||
lref(psi_non_ref_sorted_idx(ind)) = II
|
|
||||||
else
|
|
||||||
lref(psi_non_ref_sorted_idx(ind)) = -II
|
|
||||||
endif
|
|
||||||
end if
|
|
||||||
end do
|
|
||||||
wk = 0
|
|
||||||
do i=1, N_det_non_ref
|
|
||||||
if(lref(i) > 0) then
|
|
||||||
wk += 1
|
|
||||||
A_val(wk, ppp) = psi_ref_coef(lref(i), s)
|
|
||||||
A_ind(wk, ppp) = i
|
|
||||||
else if(lref(i) < 0) then
|
|
||||||
wk += 1
|
|
||||||
A_val(wk, ppp) = -psi_ref_coef(-lref(i), s)
|
|
||||||
A_ind(wk, ppp) = i
|
|
||||||
end if
|
|
||||||
end do
|
|
||||||
A_ind(0,ppp) = wk
|
|
||||||
end do
|
|
||||||
!$OMP END DO
|
|
||||||
deallocate(lref)
|
|
||||||
!$OMP END PARALLEL
|
|
||||||
|
|
||||||
|
|
||||||
print *, 'Done building A_val, A_ind'
|
|
||||||
|
|
||||||
AtA_size = 0
|
|
||||||
col_shortcut = 0
|
|
||||||
N_col = 0
|
|
||||||
integer :: a_coll, at_roww
|
|
||||||
|
|
||||||
|
|
||||||
!$OMP PARALLEL default(none) shared(k, psi_non_ref_coef, A_ind, A_val, x, N_det_ref, nex, N_det_non_ref)&
|
|
||||||
!$OMP private(at_row, a_col, t, i, j, r1, r2, wk, A_ind_mwen, A_val_mwen, a_coll, at_roww)&
|
!$OMP private(at_row, a_col, t, i, j, r1, r2, wk, A_ind_mwen, A_val_mwen, a_coll, at_roww)&
|
||||||
!$OMP shared(col_shortcut, N_col, AtB, AtA_size, AtA_val, AtA_ind, s, nactive, active_pp_idx)
|
!$OMP shared(N_states,mrcc_col_shortcut, mrcc_N_col, AtB, mrcc_AtA_val, mrcc_AtA_ind, s, n_exc_active, active_pp_idx)
|
||||||
allocate(A_val_mwen(nex), A_ind_mwen(nex))
|
allocate(A_val_mwen(N_states,hh_nex), A_ind_mwen(hh_nex), t(N_states))
|
||||||
|
|
||||||
!$OMP DO schedule(dynamic, 100)
|
!$OMP DO schedule(dynamic, 100)
|
||||||
do at_roww = 1, nactive ! nex
|
do at_roww = 1, n_exc_active ! hh_nex
|
||||||
at_row = active_pp_idx(at_roww)
|
at_row = active_pp_idx(at_roww)
|
||||||
wk = 0
|
do i=1,active_excitation_to_determinants_idx(0,at_roww)
|
||||||
if(mod(at_roww, 100) == 0) print *, "AtA", at_row, "/", nex
|
AtB(at_row) = AtB(at_row) + psi_non_ref_coef(active_excitation_to_determinants_idx(i, at_roww), s) * active_excitation_to_determinants_val(s,i, at_roww)
|
||||||
do i=1,A_ind(0,at_roww)
|
|
||||||
j = active_pp_idx(i)
|
|
||||||
AtB(at_row) = AtB(at_row) + psi_non_ref_coef(A_ind(i, at_roww), s) * A_val(i, at_roww)
|
|
||||||
end do
|
end do
|
||||||
|
|
||||||
do a_coll = 1, nactive
|
|
||||||
a_col = active_pp_idx(a_coll)
|
|
||||||
t = 0d0
|
|
||||||
r1 = 1
|
|
||||||
r2 = 1
|
|
||||||
do while ((A_ind(r1, at_roww) /= 0).and.(A_ind(r2, a_coll) /= 0))
|
|
||||||
if(A_ind(r1, at_roww) > A_ind(r2, a_coll)) then
|
|
||||||
r2 = r2+1
|
|
||||||
else if(A_ind(r1, at_roww) < A_ind(r2, a_coll)) then
|
|
||||||
r1 = r1+1
|
|
||||||
else
|
|
||||||
t = t - A_val(r1, at_roww) * A_val(r2, a_coll)
|
|
||||||
r1 = r1+1
|
|
||||||
r2 = r2+1
|
|
||||||
end if
|
|
||||||
end do
|
|
||||||
|
|
||||||
if(a_col == at_row) then
|
|
||||||
t = t + 1.d0
|
|
||||||
end if
|
|
||||||
if(t /= 0.d0) then
|
|
||||||
wk += 1
|
|
||||||
A_ind_mwen(wk) = a_col
|
|
||||||
A_val_mwen(wk) = t
|
|
||||||
end if
|
|
||||||
end do
|
|
||||||
|
|
||||||
if(wk /= 0) then
|
|
||||||
!$OMP CRITICAL
|
|
||||||
col_shortcut(at_roww) = AtA_size+1
|
|
||||||
N_col(at_roww) = wk
|
|
||||||
if (AtA_size+wk > size(AtA_ind,1)) then
|
|
||||||
print *, AtA_size+wk , size(AtA_ind,1)
|
|
||||||
stop 'too small'
|
|
||||||
endif
|
|
||||||
do i=1,wk
|
|
||||||
AtA_ind(AtA_size+i) = A_ind_mwen(i)
|
|
||||||
AtA_val(AtA_size+i) = A_val_mwen(i)
|
|
||||||
enddo
|
|
||||||
AtA_size += wk
|
|
||||||
!$OMP END CRITICAL
|
|
||||||
end if
|
|
||||||
end do
|
end do
|
||||||
!$OMP END DO NOWAIT
|
!$OMP END DO NOWAIT
|
||||||
deallocate (A_ind_mwen, A_val_mwen)
|
deallocate (A_ind_mwen, A_val_mwen)
|
||||||
!$OMP END PARALLEL
|
!$OMP END PARALLEL
|
||||||
|
|
||||||
print *, "ATA SIZE", ata_size
|
|
||||||
x = 0d0
|
x = 0d0
|
||||||
|
|
||||||
|
|
||||||
do a_coll = 1, nactive
|
do a_coll = 1, n_exc_active
|
||||||
a_col = active_pp_idx(a_coll)
|
a_col = active_pp_idx(a_coll)
|
||||||
X(a_col) = AtB(a_col)
|
X(a_col) = AtB(a_col)
|
||||||
end do
|
end do
|
||||||
@ -827,7 +680,7 @@ END_PROVIDER
|
|||||||
!$OMP DO schedule(static, 1)
|
!$OMP DO schedule(static, 1)
|
||||||
do hh = 1, hh_shortcut(0)
|
do hh = 1, hh_shortcut(0)
|
||||||
do pp = hh_shortcut(hh), hh_shortcut(hh+1)-1
|
do pp = hh_shortcut(hh), hh_shortcut(hh+1)-1
|
||||||
if(active(pp)) cycle
|
if(is_active_exc(pp)) cycle
|
||||||
lref = 0
|
lref = 0
|
||||||
do II=1,N_det_ref
|
do II=1,N_det_ref
|
||||||
call apply_hole_local(psi_ref(1,1,II), hh_exists(1, hh), myMask, ok, N_int)
|
call apply_hole_local(psi_ref(1,1,II), hh_exists(1, hh), myMask, ok, N_int)
|
||||||
@ -872,11 +725,11 @@ END_PROVIDER
|
|||||||
!$OMP END DO
|
!$OMP END DO
|
||||||
|
|
||||||
!$OMP DO
|
!$OMP DO
|
||||||
do a_coll = 1, nactive !: nex
|
do a_coll = 1, n_exc_active !: hh_nex
|
||||||
a_col = active_pp_idx(a_coll)
|
a_col = active_pp_idx(a_coll)
|
||||||
cx = 0d0
|
cx = 0d0
|
||||||
do i=col_shortcut(a_coll), col_shortcut(a_coll) + N_col(a_coll) - 1
|
do i=mrcc_col_shortcut(a_coll), mrcc_col_shortcut(a_coll) + mrcc_N_col(a_coll) - 1
|
||||||
cx = cx + x(AtA_ind(i)) * AtA_val(i)
|
cx = cx + x(mrcc_AtA_ind(i)) * mrcc_AtA_val(s,i)
|
||||||
end do
|
end do
|
||||||
x_new(a_col) = AtB(a_col) + cx * factor
|
x_new(a_col) = AtB(a_col) + cx * factor
|
||||||
end do
|
end do
|
||||||
@ -887,12 +740,12 @@ END_PROVIDER
|
|||||||
|
|
||||||
|
|
||||||
res = 0.d0
|
res = 0.d0
|
||||||
do a_coll=1,nactive ! nex
|
do a_coll=1,n_exc_active ! hh_nex
|
||||||
a_col = active_pp_idx(a_coll)
|
a_col = active_pp_idx(a_coll)
|
||||||
do j=1,N_det_non_ref
|
do j=1,N_det_non_ref
|
||||||
i = A_ind(j,a_coll)
|
i = active_excitation_to_determinants_idx(j,a_coll)
|
||||||
if (i==0) exit
|
if (i==0) exit
|
||||||
rho_mrcc(i,s) = rho_mrcc(i,s) + A_val(j,a_coll) * X_new(a_col)
|
rho_mrcc(i,s) = rho_mrcc(i,s) + active_excitation_to_determinants_val(s,j,a_coll) * X_new(a_col)
|
||||||
enddo
|
enddo
|
||||||
res = res + (X_new(a_col) - X(a_col))*(X_new(a_col) - X(a_col))
|
res = res + (X_new(a_col) - X(a_col))*(X_new(a_col) - X(a_col))
|
||||||
X(a_col) = X_new(a_col)
|
X(a_col) = X_new(a_col)
|
||||||
@ -1051,6 +904,7 @@ END_PROVIDER
|
|||||||
! Avoid numerical instabilities
|
! Avoid numerical instabilities
|
||||||
f = min(f,2.d0)
|
f = min(f,2.d0)
|
||||||
f = max(f,-2.d0)
|
f = max(f,-2.d0)
|
||||||
|
f = 1.d0
|
||||||
|
|
||||||
norm = norm + f*f *rho_mrcc(i,s)*rho_mrcc(i,s)
|
norm = norm + f*f *rho_mrcc(i,s)*rho_mrcc(i,s)
|
||||||
rho_mrcc(i,s) = f
|
rho_mrcc(i,s) = f
|
||||||
@ -1180,9 +1034,21 @@ end function
|
|||||||
|
|
||||||
|
|
||||||
BEGIN_PROVIDER [ integer*2, hh_exists, (4, N_hh_exists) ]
|
BEGIN_PROVIDER [ integer*2, hh_exists, (4, N_hh_exists) ]
|
||||||
&BEGIN_PROVIDER [ integer, hh_shortcut, (0:N_hh_exists + 1) ]
|
|
||||||
&BEGIN_PROVIDER [ integer*2, pp_exists, (4, N_pp_exists) ]
|
&BEGIN_PROVIDER [ integer*2, pp_exists, (4, N_pp_exists) ]
|
||||||
|
&BEGIN_PROVIDER [ integer, hh_shortcut, (0:N_hh_exists + 1) ]
|
||||||
|
&BEGIN_PROVIDER [ integer, hh_nex ]
|
||||||
implicit none
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
!
|
||||||
|
! hh_exists :
|
||||||
|
!
|
||||||
|
! pp_exists :
|
||||||
|
!
|
||||||
|
! hh_shortcut :
|
||||||
|
!
|
||||||
|
! hh_nex : Total number of excitation operators
|
||||||
|
!
|
||||||
|
END_DOC
|
||||||
integer*2,allocatable :: num(:,:)
|
integer*2,allocatable :: num(:,:)
|
||||||
integer :: exc(0:2, 2, 2), degree, n, on, s, l, i
|
integer :: exc(0:2, 2, 2), degree, n, on, s, l, i
|
||||||
integer*2 :: h1, h2, p1, p2
|
integer*2 :: h1, h2, p1, p2
|
||||||
@ -1248,6 +1114,7 @@ end function
|
|||||||
end if
|
end if
|
||||||
end do
|
end do
|
||||||
end do
|
end do
|
||||||
|
hh_nex = hh_shortcut(hh_shortcut(0)+1)-1
|
||||||
END_PROVIDER
|
END_PROVIDER
|
||||||
|
|
||||||
|
|
||||||
|
@ -66,7 +66,7 @@ subroutine print_cas_coefs
|
|||||||
print *, 'CAS'
|
print *, 'CAS'
|
||||||
print *, '==='
|
print *, '==='
|
||||||
do i=1,N_det_cas
|
do i=1,N_det_cas
|
||||||
print *, psi_cas_coef(i,:)
|
print *, (psi_cas_coef(i,j), j=1,N_states)
|
||||||
call debug_det(psi_cas(1,1,i),N_int)
|
call debug_det(psi_cas(1,1,i),N_int)
|
||||||
enddo
|
enddo
|
||||||
call write_double(6,ci_energy(1),"Initial CI energy")
|
call write_double(6,ci_energy(1),"Initial CI energy")
|
||||||
|
@ -513,7 +513,7 @@ subroutine i_H_j(key_i,key_j,Nint,hij)
|
|||||||
integer :: occ(Nint*bit_kind_size,2)
|
integer :: occ(Nint*bit_kind_size,2)
|
||||||
double precision :: diag_H_mat_elem, phase,phase_2
|
double precision :: diag_H_mat_elem, phase,phase_2
|
||||||
integer :: n_occ_ab(2)
|
integer :: n_occ_ab(2)
|
||||||
PROVIDE mo_bielec_integrals_in_map mo_integrals_map
|
PROVIDE mo_bielec_integrals_in_map mo_integrals_map big_array_exchange_integrals
|
||||||
|
|
||||||
ASSERT (Nint > 0)
|
ASSERT (Nint > 0)
|
||||||
ASSERT (Nint == N_int)
|
ASSERT (Nint == N_int)
|
||||||
|
Loading…
Reference in New Issue
Block a user