scemama/qp_plugins_scemama
1
0
Fork 0
mirror of https://gitlab.com/scemama/qp_plugins_scemama.git synced 2024-11-07 06:33:40 +01:00
Code Issues Releases Wiki Activity

Rewrote qmc_e_curve

Browse Source
This commit is contained in:
Anthony Scemama 2023-05-24 15:23:39 +02:00
parent edc592e012
commit 8efcf7a649
3 changed files with 86 additions and 71 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

8
devel/dmc_dress/dressing_vector.irp.f
View File

@ -29,6 +29,14 @@
END_PROVIDER
BEGIN_PROVIDER [ double precision, dressing_delta, (1,1) ]
implicit none
BEGIN_DOC
! TODO
END_DOC
stop 'dressing_delta not implemented'
END_PROVIDER
BEGIN_PROVIDER [ double precision, dmc_delta_htc , (n_det) ]

147
devel/qmcchem/qmc_e_curve.irp.f
View File

@ -1,33 +1,12 @@
program e_curve
use bitmasks
implicit none
integer :: i,j,k, kk, nab, m, l
double precision :: norm, E, hij, num, ci, cj
integer :: i,j,k, nab, m, l
double precision :: norm
integer, allocatable :: iorder(:)
double precision , allocatable :: norm_sort(:)
double precision :: e_0(N_states)
PROVIDE mo_two_e_integrals_in_map mo_one_e_integrals
nab = n_det_alpha_unique+n_det_beta_unique
allocate ( norm_sort(0:nab), iorder(0:nab) )
double precision, allocatable :: u_t(:,:), v_t(:,:), s_t(:,:)
double precision, allocatable :: u_0(:,:), v_0(:,:)
norm_sort(0) = 0.d0
iorder(0) = 0
do i=1,n_det_alpha_unique
norm_sort(i) = det_alpha_norm(i)
iorder(i) = i
enddo
do i=1,n_det_beta_unique
norm_sort(i+n_det_alpha_unique) = det_beta_norm(i)
iorder(i+n_det_alpha_unique) = -i
enddo
call dsort(norm_sort(1),iorder(1),nab)
if (.not.read_wf) then
stop 'Please set read_wf to true'
@ -35,40 +14,67 @@ program e_curve
PROVIDE psi_bilinear_matrix_values nuclear_repulsion
if (.True.) then
print *, ''
print *, 'Energy: ', psi_energy(1)+nuclear_repulsion
endif
print *, ''
print *, '=============================='
print *, 'Energies at different cut-offs'
print *, '=============================='
print *, ''
print *, '=========================================================='
print '(A8,2X,A8,2X,A12,2X,A10,2X,A12)', 'Thresh.', 'Ndet', 'Cost', 'Norm', 'E'
print *, '=========================================================='
print *, '=============================================================================='
print '(A8,2X,A8,A8,A8,2X,A12,2X,A10,2X,A12)', 'Thresh.', 'Ndet', 'Na', 'Nb', 'Cost', 'Norm', 'E'
print *, '=============================================================================='
double precision :: thresh
integer(bit_kind), allocatable :: det_i(:,:), det_j(:,:)
integer :: na, nb
thresh = 1.d-10
na = n_det_alpha_unique
nb = n_det_beta_unique
do j=0,nab
i = iorder(j)
if (i<0) then
nb -= 1
do k=1,n_det
if (psi_bilinear_matrix_columns(k) == -i) then
psi_bilinear_matrix_values(k,1) = 0.d0
endif
enddo
else
na -= 1
do k=1,n_det
if (psi_bilinear_matrix_rows(k) == i) then
psi_bilinear_matrix_values(k,1) = 0.d0
endif
enddo
endif
if (thresh > norm_sort(j)) then
cycle
endif
nab = n_det_alpha_unique+n_det_beta_unique
allocate ( norm_sort(0:nab), iorder(0:nab) )
thresh = 1.d-8
do while (N_det > 1)
nab = n_det_alpha_unique+n_det_beta_unique
norm_sort(0) = 0.d0
iorder(0) = 0
do i=1,n_det_alpha_unique
norm_sort(i) = det_alpha_norm(i)
iorder(i) = i
enddo
do i=1,n_det_beta_unique
norm_sort(i+n_det_alpha_unique) = det_beta_norm(i)
iorder(i+n_det_alpha_unique) = -i
enddo
call dsort(norm_sort(1),iorder(1),nab)
na = n_det_alpha_unique
nb = n_det_beta_unique
do j=1,nab
i = iorder(j)
if ((i<0).and.(nb>1)) then
nb -= 1
do k=1,n_det
if (psi_bilinear_matrix_columns(k) == -i) then
psi_bilinear_matrix_values(k,1) = 0.d0
endif
enddo
else if ((i>0).and.(na>1)) then
na -= 1
do k=1,n_det
if (psi_bilinear_matrix_rows(k) == i) then
psi_bilinear_matrix_values(k,1) = 0.d0
endif
enddo
endif
if (thresh < norm_sort(j)) then
exit
endif
enddo
do k=1,N_states
psi_coef(1:N_det,k) = psi_bilinear_matrix_values(1:N_det,k)
@ -76,30 +82,31 @@ program e_curve
enddo
TOUCH psi_det psi_coef
m = 0
do k=1,n_det
if (psi_bilinear_matrix_values(k,1) /= 0.d0) then
m = m+1
endif
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
do m=1,n_det
if (psi_coef_sorted(m,1) == 0.d0) exit
enddo
if (m == 0) then
exit
N_det = m-1
TOUCH psi_det psi_coef N_det
! Avoid providing psi_energy
if (.True.) then
double precision :: cost0, cost
cost0 = elec_alpha_num**3 + elec_beta_num**3
cost = (na-1) * elec_alpha_num**2 + &
(nb-1) * elec_beta_num**2 + &
elec_alpha_num**3 + elec_beta_num**3
cost = cost/cost0
double precision :: u_dot_u
norm = dsqrt(u_dot_u(psi_coef(1,1),N_det))
print '(E9.1,2X,I8,I8,I8,2X,F10.2,2X,F10.8,2X,F15.10)', thresh, N_det, &
na, nb, cost, norm, psi_energy(1) + nuclear_repulsion
thresh = thresh * dsqrt(10.d0)
endif
E = E_0(1) + nuclear_repulsion
double precision :: cost0, cost
cost0 = elec_alpha_num**3 + elec_beta_num**3
cost = (na-1) * elec_alpha_num**2 + &
(nb-1) * elec_beta_num**2 + &
elec_alpha_num**3 + elec_beta_num**3
cost = cost/cost0
double precision :: u_dot_u
norm = dsqrt(u_dot_u(psi_coef(1,1),N_det))
print '(E9.1,2X,I8,2X,F10.2,2X,F10.8,2X,F15.10)', thresh, m, &
cost, norm, psi_energy(1)
thresh = thresh * dsqrt(10.d0)
enddo
print *, '=========================================================='

2
devel/trexio/uninstall
View File

@ -14,7 +14,7 @@ scripts_list="qp_import_trexio.py"
# ${QP_ROOT}/scripts/ directory.
for i in $scripts_list
do
find ${QP_ROOT}/scripts/$i -type l -delete
find ${QP_ROOT}/scripts/$i -type l -delete || :
done