LCPQ/quantum_package
10
0
Fork 0
mirror of https://github.com/LCPQ/quantum_package synced 2024-12-22 20:35:19 +01:00
Code Issues Releases Wiki Activity

MRCC is parallelized

Browse Source
This commit is contained in:
Anthony Scemama 2015-04-09 23:59:06 +02:00
parent 0f2bdedb90
commit b5b211aaea
5 changed files with 85 additions and 48 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

37
src/Dets/determinants.irp.f
View File

@ -467,33 +467,52 @@ END_PROVIDER
! to accelerate the search of a random determinant in the wave
! function.
END_DOC
call sort_dets_by_det_search_key(N_det, psi_det, psi_coef, &
psi_det_sorted_bit, psi_coef_sorted_bit)
END_PROVIDER
subroutine sort_dets_by_det_search_key(Ndet, det_in, coef_in, det_out, coef_out)
use bitmasks
implicit none
integer, intent(in) :: Ndet
integer(bit_kind), intent(in) :: det_in (N_int,2,psi_det_size)
double precision , intent(in) :: coef_in(psi_det_size,N_states)
integer(bit_kind), intent(out) :: det_out (N_int,2,psi_det_size)
double precision , intent(out) :: coef_out(psi_det_size,N_states)
BEGIN_DOC
! Determinants are sorted are sorted according to their det_search_key.
! Useful to accelerate the search of a random determinant in the wave
! function.
END_DOC
integer :: i,j,k
integer, allocatable :: iorder(:)
integer*8, allocatable :: bit_tmp(:)
integer*8, external :: det_search_key
allocate ( iorder(N_det), bit_tmp(N_det) )
allocate ( iorder(Ndet), bit_tmp(Ndet) )
do i=1,N_det
do i=1,Ndet
iorder(i) = i
!$DIR FORCEINLINE
bit_tmp(i) = det_search_key(psi_det(1,1,i),N_int)
bit_tmp(i) = det_search_key(det_in(1,1,i),N_int)
enddo
call i8sort(bit_tmp,iorder,N_det)
call i8sort(bit_tmp,iorder,Ndet)
!DIR$ IVDEP
do i=1,N_det
do i=1,Ndet
do j=1,N_int
psi_det_sorted_bit(j,1,i) = psi_det(j,1,iorder(i))
psi_det_sorted_bit(j,2,i) = psi_det(j,2,iorder(i))
det_out(j,1,i) = det_in(j,1,iorder(i))
det_out(j,2,i) = det_in(j,2,iorder(i))
enddo
do k=1,N_states
psi_coef_sorted_bit(i,k) = psi_coef(iorder(i),k)
coef_out(i,k) = coef_in(iorder(i),k)
enddo
enddo
deallocate(iorder, bit_tmp)
END_PROVIDER
end
subroutine int_of_3_highest_electrons( det_in, res, Nint )
implicit none

2
src/MRCC/H_apply.irp.f
View File

@ -18,7 +18,7 @@ s.data["decls_main"] += """
s.data["finalization"] = ""
s.data["copy_buffer"] = ""
s.data["generate_psi_guess"] = ""
s.data["size_max"] = "256"
s.data["size_max"] = "3072"
print s

6
src/MRCC/mrcc.irp.f
View File

@ -14,14 +14,14 @@ subroutine run
print *, 'CAS'
print *, '==='
do i=1,N_det_cas
print *, psi_cas_coefs(i,:)
print *, psi_cas_coef(i,:)
call debug_det(psi_cas(1,1,i),N_int)
enddo
! print *, 'SD'
! print *, '=='
! do i=1,N_det_non_cas
! print *, psi_non_cas_coefs(i,:)
! print *, psi_non_cas_coef(i,:)
! call debug_det(psi_non_cas(1,1,i),N_int)
! enddo
call write_double(6,ci_energy(1),"Initial CI energy")
@ -56,9 +56,9 @@ subroutine run_mrcc
print *, ''
E_old = sum(ci_energy_dressed)
call diagonalize_ci_dressed
call write_double(6,ci_energy_dressed(1),"MRCC energy")
E_new = sum(ci_energy_dressed)
delta_E = dabs(E_new - E_old)
call write_double(6,ci_energy_dressed(1),"MRCC energy")
enddo
call ezfio_set_mrcc_energy(ci_energy_dressed(1))
! call save_wavefunction

58
src/MRCC/mrcc_dress.irp.f
View File

@ -1,3 +1,17 @@
use omp_lib
BEGIN_PROVIDER [ integer(omp_lock_kind), psi_cas_lock, (psi_det_size) ]
implicit none
BEGIN_DOC
! Locks on CAS determinants to fill delta_ij
END_DOC
integer :: i
do i=1,psi_det_size
call omp_init_lock( psi_cas_lock(i) )
enddo
END_PROVIDER
subroutine mrcc_dress(delta_ij_,Ndet,i_generator,n_selected,det_buffer,Nint,iproc)
use bitmasks
implicit none
@ -7,11 +21,8 @@ subroutine mrcc_dress(delta_ij_,Ndet,i_generator,n_selected,det_buffer,Nint,ipro
double precision, intent(inout) :: delta_ij_(Ndet,Ndet,*)
integer(bit_kind), intent(in) :: det_buffer(Nint,2,n_selected)
integer :: i,j,k,l,m
logical :: is_in_wavefunction
integer :: i,j,k,l
integer :: degree_alpha(psi_det_size)
integer :: degree_I(psi_det_size)
integer :: idx_I(0:psi_det_size)
integer :: idx_alpha(0:psi_det_size)
logical :: good
@ -19,9 +30,8 @@ subroutine mrcc_dress(delta_ij_,Ndet,i_generator,n_selected,det_buffer,Nint,ipro
integer :: N_tq, c_ref ,degree
integer :: connected_to_ref
call find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq,N_tq)
double precision :: hIk, hIl, hla, dIk(N_states), dka(N_states), dIa(N_states)
double precision :: hIk, hla, hIl, dIk(N_states), dka(N_states), dIa(N_states)
double precision, allocatable :: dIa_hla(:,:)
double precision :: haj, phase, phase2
double precision :: f(N_states), ci_inv(N_states)
integer :: exc(0:2,2,2)
@ -30,7 +40,9 @@ subroutine mrcc_dress(delta_ij_,Ndet,i_generator,n_selected,det_buffer,Nint,ipro
integer :: iint, ipos
integer :: i_state, k_sd, l_sd, i_I, i_alpha
call find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq,N_tq)
allocate (dIa_hla(N_states,Ndet))
! |I>
@ -60,12 +72,15 @@ subroutine mrcc_dress(delta_ij_,Ndet,i_generator,n_selected,det_buffer,Nint,ipro
! <I|H|k>
call i_h_j(psi_cas(1,1,i_I),psi_non_cas(1,1,idx_alpha(k_sd)),Nint,hIk)
do i_state=1,N_states
dIk(i_state) = hIk * lambda_mrcc(idx_alpha(k_sd),i_state)
dIk(i_state) = hIk * lambda_mrcc(i_state,idx_alpha(k_sd))
enddo
! |l> = Exc(k -> alpha) |I>
call get_excitation(psi_non_cas(1,1,idx_alpha(k_sd)),tq(1,1,i_alpha),exc,degree,phase,Nint)
call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)
tmp_det(1:Nint,1:2) = psi_cas(1,1,i_I)
do k=1,N_int
tmp_det(k,1) = psi_cas(k,1,i_I)
tmp_det(k,2) = psi_cas(k,2,i_I)
enddo
! Hole (see list_to_bitstring)
iint = ishft(h1-1,-bit_kind_shift) + 1
ipos = h1-ishft((iint-1),bit_kind_shift)-1
@ -75,7 +90,7 @@ subroutine mrcc_dress(delta_ij_,Ndet,i_generator,n_selected,det_buffer,Nint,ipro
iint = ishft(p1-1,-bit_kind_shift) + 1
ipos = p1-ishft((iint-1),bit_kind_shift)-1
tmp_det(iint,s1) = ibset(tmp_det(iint,s1),ipos)
if (degree == 2) then
if (degree_alpha(k_sd) == 2) then
! Hole (see list_to_bitstring)
iint = ishft(h2-1,-bit_kind_shift) + 1
ipos = h2-ishft((iint-1),bit_kind_shift)-1
@ -97,7 +112,7 @@ subroutine mrcc_dress(delta_ij_,Ndet,i_generator,n_selected,det_buffer,Nint,ipro
call get_excitation(psi_cas(1,1,i_I),psi_non_cas(1,1,idx_alpha(l_sd)),exc,degree,phase2,Nint)
call i_h_j(psi_cas(1,1,i_I),psi_non_cas(1,1,idx_alpha(l_sd)),Nint,hIl)
do i_state=1,N_states
dka(i_state) = hIl * lambda_mrcc(idx_alpha(l_sd),i_state) * phase * phase2
dka(i_state) = hIl * lambda_mrcc(i_state,idx_alpha(l_sd)) * phase * phase2
enddo
exit
endif
@ -108,20 +123,28 @@ subroutine mrcc_dress(delta_ij_,Ndet,i_generator,n_selected,det_buffer,Nint,ipro
enddo
do i_state=1,N_states
ci_inv(i_state) = 1.d0/psi_cas_coefs(i_I,i_state)
ci_inv(i_state) = 1.d0/psi_cas_coef(i_I,i_state)
enddo
do l_sd=1,idx_alpha(0)
k_sd = idx_alpha(l_sd)
call i_h_j(tq(1,1,i_alpha),psi_non_cas(1,1,idx_alpha(l_sd)),Nint,hla)
do m=1,N_states
delta_ij_(idx_non_cas(k_sd),idx_cas(i_I),m) += dIa(m) * hla
delta_ij_(idx_cas(i_I),idx_non_cas(k_sd),m) += dIa(m) * hla
delta_ij_(idx_cas(i_I),idx_cas(i_I),m) -= dIa(m) * hla * ci_inv(m) * psi_non_cas_coefs(k_sd,m)
do i_state=1,N_states
dIa_hla(i_state,k_sd) = dIa(i_state) * hla
enddo
enddo
call omp_set_lock( psi_cas_lock(i_I) )
do l_sd=1,idx_alpha(0)
k_sd = idx_alpha(l_sd)
do i_state=1,N_states
delta_ij_(idx_non_cas(k_sd),idx_cas(i_I),i_state) += dIa_hla(i_state,k_sd)
delta_ij_(idx_cas(i_I),idx_non_cas(k_sd),i_state) += dIa_hla(i_state,k_sd)
delta_ij_(idx_cas(i_I),idx_cas(i_I),i_state) -= dIa_hla(i_state,k_sd) * ci_inv(i_state) * psi_non_cas_coef(k_sd,i_state)
enddo
enddo
call omp_unset_lock( psi_cas_lock(i_I) )
enddo
enddo
deallocate (dIa_hla)
end
@ -141,7 +164,6 @@ subroutine mrcc_dress_simple(delta_ij_non_cas_,Ndet_non_cas,i_generator,n_select
integer(bit_kind), intent(in) :: det_buffer(Nint,2,n_selected)
integer :: i,j,k,m
integer :: new_size
logical :: is_in_wavefunction
integer :: degree(psi_det_size)
integer :: idx(0:psi_det_size)
logical :: good

14
src/MRCC/mrcc_utils.irp.f
View File

@ -1,5 +1,5 @@
BEGIN_PROVIDER [ double precision, lambda_mrcc, (psi_det_size,n_states) ]
BEGIN_PROVIDER [ double precision, lambda_mrcc, (N_states,psi_det_size) ]
implicit none
BEGIN_DOC
! cm/<Psi_0|H|D_m>
@ -7,16 +7,12 @@ BEGIN_PROVIDER [ double precision, lambda_mrcc, (psi_det_size,n_states) ]
integer :: i,k
double precision :: ihpsi(N_states)
do i=1,N_det_non_cas
call i_h_psi(psi_non_cas(1,1,i), psi_cas, psi_cas_coefs, N_int, N_det_cas, &
size(psi_cas_coefs,1), n_states, ihpsi)
call i_h_psi(psi_non_cas(1,1,i), psi_cas, psi_cas_coef, N_int, N_det_cas, &
size(psi_cas_coef,1), n_states, ihpsi)
double precision :: hij
do k=1,N_states
if (dabs(ihpsi(k)) < 1.d-6) then
lambda_mrcc(i,k) = 0.d0
else
lambda_mrcc(i,k) = psi_non_cas_coefs(i,k)/ihpsi(k)
lambda_mrcc(i,k) = min( lambda_mrcc (i,k),0.d0 )
endif
lambda_mrcc(k,i) = psi_non_cas_coef(i,k)/ihpsi(k)
lambda_mrcc(k,i) = min( lambda_mrcc (k,i),0.d0 )
enddo
enddo
END_PROVIDER