LCPQ/quantum_package
10
0
Fork 0
mirror of https://github.com/LCPQ/quantum_package synced 2024-07-22 10:47:33 +02:00
Code Issues Releases Wiki Activity

Working on davidson

Browse Source
This commit is contained in:
Anthony Scemama 2016-09-23 16:16:48 +02:00
parent 18ff53e063
commit c8a5cf37cd
4 changed files with 58 additions and 47 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1
configure vendored
View File

@ -543,7 +543,6 @@ def recommendation():
print ""
print "Finally :"
print " ninja"
print " make -C ocaml"
print ""
print "You can install more plugin with the qp_module.py install command"
print "PS : For more info on compiling the code, read the README.md"

50
plugins/Full_CI_ZMQ/fci_zmq.irp.f
View File

@ -3,12 +3,10 @@ program fci_zmq
integer :: i,j,k
logical, external :: detEq
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer :: N_st, degree
integer(bit_kind) :: chk
double precision, allocatable :: pt2(:)
integer :: degree
N_st = N_states
allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st))
allocate (pt2(N_states))
pt2 = 1.d0
diag_algorithm = "Lapack"
@ -24,20 +22,22 @@ program fci_zmq
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
do k=1,N_states
print*,'State ',k
print *, 'PT2 = ', pt2(k)
print *, 'E = ', CI_energy(k)
print *, 'E+PT2 = ', CI_energy(k) + pt2(k)
print *, '-----'
enddo
endif
double precision :: i_H_psi_array(N_states),diag_H_mat_elem,h,i_O1_psi_array(N_states)
double precision :: E_CI_before(N_states)
integer :: n_det_before
print*,'Beginning the selection ...'
E_CI_before = CI_energy
E_CI_before(1:N_states) = CI_energy(1:N_states)
do while ( (N_det < N_det_max) .and. (maxval(abs(pt2(1:N_st))) > pt2_max) )
do while ( (N_det < N_det_max) .and. (maxval(abs(pt2(1:N_states))) > pt2_max) )
n_det_before = N_det
call ZMQ_selection(max(1024-N_det, N_det), pt2)
@ -56,14 +56,13 @@ program fci_zmq
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
do k = 1, N_states
print*,'State ',k
print *, 'PT2 = ', pt2(k)
print *, 'E = ', CI_energy(k)
print *, 'E(before)+PT2 = ', E_CI_before(k)+pt2(k)
do k=1, N_states
print*,'State ',k
print *, 'PT2 = ', pt2(k)
print *, 'E = ', CI_energy(k)
print *, 'E(before)+PT2 = ', E_CI_before(k)+pt2(k)
enddo
print *, '-----'
E_CI_before = CI_energy
if(N_states.gt.1)then
print*,'Variational Energy difference'
do i = 2, N_states
@ -76,7 +75,7 @@ program fci_zmq
print*,'Delta E = ',E_CI_before(i)+ pt2(i) - (E_CI_before(1) + pt2(1))
enddo
endif
E_CI_before = CI_energy
E_CI_before(1:N_states) = CI_energy(1:N_states)
call ezfio_set_full_ci_energy(CI_energy)
enddo
N_det = min(N_det_max,N_det)
@ -86,15 +85,18 @@ program fci_zmq
print*,'Last iteration only to compute the PT2'
threshold_selectors = 1.d0
threshold_generators = 0.9999d0
E_CI_before = CI_energy
E_CI_before(1:N_states) = CI_energy(1:N_states)
call ZMQ_selection(1, pt2)
print *, 'Final step'
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', E_CI_before
print *, 'E+PT2 = ', E_CI_before+pt2
print *, '-----'
do k=1,N_states
print *, 'State', k
print *, 'PT2 = ', pt2
print *, 'E = ', E_CI_before
print *, 'E+PT2 = ', E_CI_before+pt2
print *, '-----'
enddo
call ezfio_set_full_ci_energy_pt2(E_CI_before+pt2)
endif
call save_wavefunction

23
src/Davidson/diagonalization.irp.f
View File

@ -505,24 +505,21 @@ subroutine davidson_diag_hjj(dets_in,u_in,H_jj,energies,dim_in,sze,N_st,N_st_dia
! Gram-Schmidt
! ------------
double precision :: c
double precision :: c(N_st_diag)
do k=1,N_st_diag
do iter2=1,iter
do l=1,N_st_diag
c = u_dot_v(U(1,k,iter+1),U(1,l,iter2),sze)
do i=1,sze
U(i,k,iter+1) = U(i,k,iter+1) - c * U(i,l,iter2)
enddo
enddo
enddo
do l=1,k-1
c = u_dot_v(U(1,k,iter+1),U(1,l,iter+1),sze)
do i=1,sze
U(i,k,iter+1) = U(i,k,iter+1) - c * U(i,l,iter+1)
enddo
call dgemv('T',sze,N_st_diag,1.d0,U(1,1,iter2),size(U,1), &
U(1,k,iter+1),1,0.d0,c,1)
call dgemv('N',sze,N_st_diag,-1.d0,U(1,1,iter2),size(U,1), &
c,1,1.d0,U(1,k,iter+1),1)
enddo
call dgemv('T',sze,N_st_diag,1.d0,U(1,1,iter+1),size(U,1), &
U(1,k,iter+1),1,0.d0,c,1)
call dgemv('N',sze,k-1,-1.d0,U(1,1,iter+1),size(U,1), &
c,1,1.d0,U(1,k,iter+1),1)
call normalize( U(1,k,iter+1), sze )
enddo
!DEBUG : CHECK OVERLAP
!print *, '==='

31
src/Davidson/u0Hu0.irp.f
View File

@ -43,7 +43,7 @@ subroutine H_u_0_nstates(v_0,u_0,H_jj,n,keys_tmp,Nint,N_st,sze_8)
double precision, intent(in) :: H_jj(n)
integer(bit_kind),intent(in) :: keys_tmp(Nint,2,n)
double precision :: hij
double precision, allocatable :: vt(:,:)
double precision, allocatable :: vt(:,:), ut(:,:)
integer :: i,j,k,l, jj,ii
integer :: i0, j0
@ -52,7 +52,12 @@ subroutine H_u_0_nstates(v_0,u_0,H_jj,n,keys_tmp,Nint,N_st,sze_8)
integer(bit_kind) :: sorted_i(Nint)
integer :: sh, sh2, ni, exa, ext, org_i, org_j, endi, istate
integer :: N_st_8
integer, external :: align_double
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: vt, ut
N_st_8 = align_double(N_st)
ASSERT (Nint > 0)
ASSERT (Nint == N_int)
@ -60,15 +65,23 @@ subroutine H_u_0_nstates(v_0,u_0,H_jj,n,keys_tmp,Nint,N_st,sze_8)
PROVIDE ref_bitmask_energy
allocate (shortcut(0:n+1,2), sort_idx(n,2), sorted(Nint,n,2), version(Nint,n,2))
allocate(ut(N_st_8,n))
v_0 = 0.d0
do i=1,n
do istate=1,N_st
ut(istate,i) = u_0(i,istate)
enddo
enddo
call sort_dets_ab_v(keys_tmp, sorted(1,1,1), sort_idx(1,1), shortcut(0,1), version(1,1,1), n, Nint)
call sort_dets_ba_v(keys_tmp, sorted(1,1,2), sort_idx(1,2), shortcut(0,2), version(1,1,2), n, Nint)
!$OMP PARALLEL DEFAULT(NONE) &
!$OMP PRIVATE(i,hij,j,k,jj,vt,ii,sh,sh2,ni,exa,ext,org_i,org_j,endi,sorted_i,istate)&
!$OMP SHARED(n,H_jj,u_0,keys_tmp,Nint,v_0,sorted,shortcut,sort_idx,version,N_st,sze_8)
allocate(vt(sze_8,N_st))
!$OMP SHARED(n,H_jj,keys_tmp,ut,Nint,v_0,sorted,shortcut,sort_idx,version,N_st,N_st_8)
allocate(vt(N_st_8,n))
Vt = 0.d0
!$OMP DO SCHEDULE(dynamic)
@ -102,8 +115,8 @@ subroutine H_u_0_nstates(v_0,u_0,H_jj,n,keys_tmp,Nint,N_st,sze_8)
if(ext <= 4) then
call i_H_j(keys_tmp(1,1,org_j),keys_tmp(1,1,org_i),Nint,hij)
do istate=1,N_st
vt (org_i,istate) = vt (org_i,istate) + hij*u_0(org_j,istate)
vt (org_j,istate) = vt (org_j,istate) + hij*u_0(org_i,istate)
vt (istate,org_i) = vt (istate,org_i) + hij*ut(istate,org_j)
vt (istate,org_j) = vt (istate,org_j) + hij*ut(istate,org_i)
enddo
endif
enddo
@ -125,8 +138,8 @@ subroutine H_u_0_nstates(v_0,u_0,H_jj,n,keys_tmp,Nint,N_st,sze_8)
if(ext == 4) then
call i_H_j(keys_tmp(1,1,org_j),keys_tmp(1,1,org_i),Nint,hij)
do istate=1,N_st
vt (org_i,istate) = vt (org_i,istate) + hij*u_0(org_j,istate)
vt (org_j,istate) = vt (org_j,istate) + hij*u_0(org_i,istate)
vt (istate,org_i) = vt (istate,org_i) + hij*ut(istate,org_j)
vt (istate,org_j) = vt (istate,org_j) + hij*ut(istate,org_i)
enddo
end if
end do
@ -137,7 +150,7 @@ subroutine H_u_0_nstates(v_0,u_0,H_jj,n,keys_tmp,Nint,N_st,sze_8)
!$OMP CRITICAL
do istate=1,N_st
do i=n,1,-1
v_0(i,istate) = v_0(i,istate) + vt(i,istate)
v_0(i,istate) = v_0(i,istate) + vt(istate,i)
enddo
enddo
!$OMP END CRITICAL
@ -150,7 +163,7 @@ subroutine H_u_0_nstates(v_0,u_0,H_jj,n,keys_tmp,Nint,N_st,sze_8)
v_0(i,istate) += H_jj(i) * u_0(i,istate)
enddo
enddo
deallocate (shortcut, sort_idx, sorted, version)
deallocate (shortcut, sort_idx, sorted, version, ut)
end
BEGIN_PROVIDER [ double precision, psi_energy, (N_states) ]