LCPQ/qp2
9
1
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-12-22 03:23:29 +01:00
Code Releases Activity

Merge branch 'dev-stable' of github.com:QuantumPackage/qp2 into dev-stable

Browse Source
This commit is contained in:
Anthony Scemama 2024-04-03 15:34:16 +02:00
commit 4190aee606
2 changed files with 36 additions and 18 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

24
src/ccsd/ccsd_t_space_orb_stoch.irp.f
View File

@ -181,8 +181,8 @@ subroutine ccsd_par_t_space_stoch(nO,nV,t1,t2,f_o,f_v,v_vvvo,v_vvoo,v_vooo,energ
integer :: nbuckets integer :: nbuckets
nbuckets = 100 nbuckets = 100
double precision, allocatable :: ED(:)
double precision, allocatable :: wsum(:) double precision, allocatable :: wsum(:)
allocate(wsum(nbuckets))
converged = .False. converged = .False.
Ncomputed = 0_8 Ncomputed = 0_8
@ -197,7 +197,8 @@ subroutine ccsd_par_t_space_stoch(nO,nV,t1,t2,f_o,f_v,v_vvvo,v_vvoo,v_vooo,energ
iright = Nabc iright = Nabc
integer*8, allocatable :: bounds(:,:) integer*8, allocatable :: bounds(:,:)
allocate (bounds(2,nbuckets)) allocate(wsum(nbuckets), ED(nbuckets), bounds(2,nbuckets))
ED(:) = 0.d0
do isample=1,nbuckets do isample=1,nbuckets
eta = 1.d0/dble(nbuckets) * dble(isample) eta = 1.d0/dble(nbuckets) * dble(isample)
ieta = binary_search(waccu,eta,Nabc) ieta = binary_search(waccu,eta,Nabc)
@ -233,7 +234,7 @@ subroutine ccsd_par_t_space_stoch(nO,nV,t1,t2,f_o,f_v,v_vvvo,v_vvoo,v_vooo,energ
enddo enddo
! Deterministic part ! Deterministic part
if (imin < Nabc) then if (imin <= Nabc) then
ieta=imin ieta=imin
sampled(ieta) = 0_8 sampled(ieta) = 0_8
a = abc(1,ieta) a = abc(1,ieta)
@ -254,7 +255,7 @@ subroutine ccsd_par_t_space_stoch(nO,nV,t1,t2,f_o,f_v,v_vvvo,v_vvoo,v_vooo,energ
! Stochastic part ! Stochastic part
call random_number(eta) call random_number(eta)
do isample=1,nbuckets do isample=1,nbuckets
if (imin >= bounds(2,isample)) then if (imin > bounds(2,isample)) then
cycle cycle
endif endif
ieta = binary_search(waccu,(eta + dble(isample-1))/dble(nbuckets),Nabc)+1 ieta = binary_search(waccu,(eta + dble(isample-1))/dble(nbuckets),Nabc)+1
@ -280,7 +281,7 @@ subroutine ccsd_par_t_space_stoch(nO,nV,t1,t2,f_o,f_v,v_vvvo,v_vvoo,v_vooo,energ
enddo enddo
call wall_time(t01) call wall_time(t01)
if ((t01-t00 > 1.0d0).or.(imin >= Nabc)) then if ((t01-t00 > 1.0d0).or.(imin > Nabc)) then
!$OMP TASKWAIT !$OMP TASKWAIT
call wall_time(t01) call wall_time(t01)
@ -300,8 +301,11 @@ subroutine ccsd_par_t_space_stoch(nO,nV,t1,t2,f_o,f_v,v_vvvo,v_vvoo,v_vooo,energ
do isample=1,nbuckets do isample=1,nbuckets
if (imin >= bounds(2,isample)) then if (imin > bounds(2,isample)) then
energy_det = energy_det + sum(memo(bounds(1,isample):bounds(2,isample))) if (ED(isample) == 0.d0) then
ED(isample) = sum(memo(bounds(1,isample):bounds(2,isample)))
endif
energy_det = energy_det + ED(isample)
scale = scale - wsum(isample) scale = scale - wsum(isample)
else else
exit exit
@ -310,12 +314,14 @@ subroutine ccsd_par_t_space_stoch(nO,nV,t1,t2,f_o,f_v,v_vvvo,v_vvoo,v_vooo,energ
isample = min(isample,nbuckets) isample = min(isample,nbuckets)
do ieta=bounds(1,isample), Nabc do ieta=bounds(1,isample), Nabc
w = dble(max(sampled(ieta),0_8)) if (sampled(ieta) < 0_8) cycle
w = dble(sampled(ieta))
tmp = w * memo(ieta) * Pabc(ieta) tmp = w * memo(ieta) * Pabc(ieta)
ET = ET + tmp ET = ET + tmp
ET2 = ET2 + tmp * memo(ieta) * Pabc(ieta) ET2 = ET2 + tmp * memo(ieta) * Pabc(ieta)
norm = norm + w norm = norm + w
enddo enddo
norm = norm/scale norm = norm/scale
if (norm > 0.d0) then if (norm > 0.d0) then
energy_stoch = ET / norm energy_stoch = ET / norm
@ -327,7 +333,7 @@ subroutine ccsd_par_t_space_stoch(nO,nV,t1,t2,f_o,f_v,v_vvvo,v_vvoo,v_vooo,energ
print '('' '',F20.8, '' '', ES12.4,'' '', F8.2,'' '')', eccsd+energy, dsqrt(variance/(norm-1.d0)), 100.*real(Ncomputed)/real(Nabc) print '('' '',F20.8, '' '', ES12.4,'' '', F8.2,'' '')', eccsd+energy, dsqrt(variance/(norm-1.d0)), 100.*real(Ncomputed)/real(Nabc)
endif endif
!$OMP END MASTER !$OMP END MASTER
if (imin >= Nabc) exit if (imin > Nabc) exit
enddo enddo
!$OMP END PARALLEL !$OMP END PARALLEL

22
src/davidson/diagonalization_hs2_dressed.irp.f
View File

@ -139,7 +139,7 @@ subroutine davidson_diag_hjj_sjj(dets_in,u_in,H_jj,s2_out,energies,dim_in,sze,N_
integer :: iter2, itertot integer :: iter2, itertot
double precision, allocatable :: y(:,:), h(:,:), h_p(:,:), lambda(:), s2(:) double precision, allocatable :: y(:,:), h(:,:), h_p(:,:), lambda(:), s2(:)
real, allocatable :: y_s(:,:) real, allocatable :: y_s(:,:)
double precision, allocatable :: s_(:,:), s_tmp(:,:) double precision, allocatable :: s_(:,:), s_tmp(:,:), prev_y(:,:)
double precision :: diag_h_mat_elem double precision :: diag_h_mat_elem
double precision, allocatable :: residual_norm(:) double precision, allocatable :: residual_norm(:)
character*(16384) :: write_buffer character*(16384) :: write_buffer
@ -288,6 +288,7 @@ subroutine davidson_diag_hjj_sjj(dets_in,u_in,H_jj,s2_out,energies,dim_in,sze,N_
h(N_st_diag*itermax,N_st_diag*itermax), & h(N_st_diag*itermax,N_st_diag*itermax), &
! h_p(N_st_diag*itermax,N_st_diag*itermax), & ! h_p(N_st_diag*itermax,N_st_diag*itermax), &
y(N_st_diag*itermax,N_st_diag*itermax), & y(N_st_diag*itermax,N_st_diag*itermax), &
prev_y(N_st_diag*itermax,N_st_diag*itermax), &
s_(N_st_diag*itermax,N_st_diag*itermax), & s_(N_st_diag*itermax,N_st_diag*itermax), &
s_tmp(N_st_diag*itermax,N_st_diag*itermax), & s_tmp(N_st_diag*itermax,N_st_diag*itermax), &
residual_norm(N_st_diag), & residual_norm(N_st_diag), &
@ -301,6 +302,10 @@ subroutine davidson_diag_hjj_sjj(dets_in,u_in,H_jj,s2_out,energies,dim_in,sze,N_
s_ = 0.d0 s_ = 0.d0
s_tmp = 0.d0 s_tmp = 0.d0
prev_y = 0.d0
do i = 1, N_st_diag*itermax
prev_y(i,i) = 1d0
enddo
ASSERT (N_st > 0) ASSERT (N_st > 0)
ASSERT (N_st_diag >= N_st) ASSERT (N_st_diag >= N_st)
@ -479,6 +484,10 @@ subroutine davidson_diag_hjj_sjj(dets_in,u_in,H_jj,s2_out,energies,dim_in,sze,N_
if (info > 0) then if (info > 0) then
! Numerical errors propagate. We need to reduce the number of iterations ! Numerical errors propagate. We need to reduce the number of iterations
itermax = iter-1 itermax = iter-1
! eigenvectors of the previous iteration
y = prev_y
shift2 = shift2 - N_st_diag
exit exit
endif endif
@ -553,7 +562,7 @@ subroutine davidson_diag_hjj_sjj(dets_in,u_in,H_jj,s2_out,energies,dim_in,sze,N_
do l = 1, N_st do l = 1, N_st
do k = 1, N_st_diag do k = 1, N_st_diag
do i = 1, sze do i = 1, sze
overlp(k+j-1,l) += U(i,l) * U(i,shift2+k) overlp(k+j-1,l) += u_in(i,l) * U(i,shift2+k)
enddo enddo
enddo enddo
enddo enddo
@ -576,7 +585,7 @@ subroutine davidson_diag_hjj_sjj(dets_in,u_in,H_jj,s2_out,energies,dim_in,sze,N_
enddo enddo
! Maximum overlap ! Maximum overlap
if (dabs(overlp(k,l)) > omax .and. .not. used .and. state_ok(k)) then if ((dabs(overlp(k,l)) > omax) .and. (.not. used) .and. state_ok(k)) then
omax = dabs(overlp(k,l)) omax = dabs(overlp(k,l))
idx = k idx = k
endif endif
@ -615,6 +624,9 @@ subroutine davidson_diag_hjj_sjj(dets_in,u_in,H_jj,s2_out,energies,dim_in,sze,N_
endif endif
enddo enddo
! Swapped eigenvectors
prev_y = y
! if (state_following) then ! if (state_following) then
! !
! overlap = -1.d0 ! overlap = -1.d0
@ -677,7 +689,7 @@ subroutine davidson_diag_hjj_sjj(dets_in,u_in,H_jj,s2_out,energies,dim_in,sze,N_
do i=1,sze do i=1,sze
U(i,shift2+k) = & U(i,shift2+k) = &
(lambda(k) * U(i,shift2+k) - W(i,shift2+k) ) & (lambda(k) * U(i,shift2+k) - W(i,shift2+k) ) &
/max(H_jj(i) - lambda (k),1.d-2) /max(dabs(H_jj(i) - lambda (k)),1.d-2) * dsign(1d0,H_jj(i) - lambda (k))
enddo enddo
if (k <= N_st) then if (k <= N_st) then
@ -792,7 +804,7 @@ subroutine davidson_diag_hjj_sjj(dets_in,u_in,H_jj,s2_out,energies,dim_in,sze,N_
residual_norm, & residual_norm, &
U, overlap, & U, overlap, &
h, y_s, S_d, & h, y_s, S_d, &
y, s_, s_tmp, & y, s_, s_tmp, prev_y, &
lambda & lambda &
) )
FREE nthreads_davidson FREE nthreads_davidson