From f7a7ba2a3e2920062b0e64e24e376714a156cc36 Mon Sep 17 00:00:00 2001 From: Kevin Gasperich Date: Tue, 25 Feb 2020 09:11:16 -0600 Subject: [PATCH] started complex h_apply --- src/determinants/h_apply.irp.f | 140 ++++++++++++++++++++-- src/determinants/h_apply_nozmq.template.f | 5 +- src/determinants/occ_pattern.irp.f | 4 +- src/generators_full/generators.irp.f | 41 ++++++- src/utils_complex/qp2-pbc-diff.txt | 5 + 5 files changed, 180 insertions(+), 15 deletions(-) diff --git a/src/determinants/h_apply.irp.f b/src/determinants/h_apply.irp.f index 1c79bc75..cfbbf847 100644 --- a/src/determinants/h_apply.irp.f +++ b/src/determinants/h_apply.irp.f @@ -6,6 +6,7 @@ type H_apply_buffer_type integer :: sze integer(bit_kind), pointer :: det(:,:,:) double precision , pointer :: coef(:,:) + complex*16 , pointer :: coef_complex(:,:) double precision , pointer :: e2(:,:) end type H_apply_buffer_type @@ -32,11 +33,16 @@ type(H_apply_buffer_type), pointer :: H_apply_buffer(:) H_apply_buffer(iproc)%sze = sze allocate ( & H_apply_buffer(iproc)%det(N_int,2,sze), & - H_apply_buffer(iproc)%coef(sze,N_states), & H_apply_buffer(iproc)%e2(sze,N_states) & ) + if (is_complex) then + allocate(H_apply_buffer(iproc)%coef_complex(sze,N_states)) + H_apply_buffer(iproc)%coef_complex = (0.d0,0.d0) + else + allocate(H_apply_buffer(iproc)%coef(sze,N_states)) + H_apply_buffer(iproc)%coef = 0.d0 + endif H_apply_buffer(iproc)%det = 0_bit_kind - H_apply_buffer(iproc)%coef = 0.d0 H_apply_buffer(iproc)%e2 = 0.d0 call omp_init_lock(H_apply_buffer_lock(1,iproc)) !$OMP END PARALLEL @@ -59,6 +65,7 @@ subroutine resize_H_apply_buffer(new_size,iproc) integer, intent(in) :: new_size, iproc integer(bit_kind), pointer :: buffer_det(:,:,:) double precision, pointer :: buffer_coef(:,:) + complex*16, pointer :: buffer_coef_complex(:,:) double precision, pointer :: buffer_e2(:,:) integer :: i,j,k integer :: Ndet @@ -74,9 +81,14 @@ subroutine resize_H_apply_buffer(new_size,iproc) ASSERT (iproc < nproc) allocate ( buffer_det(N_int,2,new_size), & - buffer_coef(new_size,N_states), & buffer_e2(new_size,N_states) ) - buffer_coef = 0.d0 + if (is_complex) then + allocate(buffer_coef_complex(new_size,N_states)) + buffer_coef_complex = (0.d0,0.d0) + else + allocate(buffer_coef(new_size,N_states)) + buffer_coef = 0.d0 + endif buffer_e2 = 0.d0 do i=1,min(new_size,H_apply_buffer(iproc)%N_det) do k=1,N_int @@ -89,6 +101,15 @@ subroutine resize_H_apply_buffer(new_size,iproc) deallocate(H_apply_buffer(iproc)%det) H_apply_buffer(iproc)%det => buffer_det + if (is_complex) then + do k=1,N_states + do i=1,min(new_size,H_apply_buffer(iproc)%N_det) + buffer_coef_complex(i,k) = H_apply_buffer(iproc)%coef_complex(i,k) + enddo + enddo + deallocate(H_apply_buffer(iproc)%coef_complex) + H_apply_buffer(iproc)%coef_complex => buffer_coef_complex + else do k=1,N_states do i=1,min(new_size,H_apply_buffer(iproc)%N_det) buffer_coef(i,k) = H_apply_buffer(iproc)%coef(i,k) @@ -96,6 +117,7 @@ subroutine resize_H_apply_buffer(new_size,iproc) enddo deallocate(H_apply_buffer(iproc)%coef) H_apply_buffer(iproc)%coef => buffer_coef + endif do k=1,N_states do i=1,min(new_size,H_apply_buffer(iproc)%N_det) @@ -119,6 +141,7 @@ subroutine copy_H_apply_buffer_to_wf END_DOC integer(bit_kind), allocatable :: buffer_det(:,:,:) double precision, allocatable :: buffer_coef(:,:) + complex*16, allocatable :: buffer_coef_complex(:,:) integer :: i,j,k integer :: N_det_old @@ -128,7 +151,12 @@ subroutine copy_H_apply_buffer_to_wf ASSERT (N_int > 0) ASSERT (N_det > 0) - allocate ( buffer_det(N_int,2,N_det), buffer_coef(N_det,N_states) ) + allocate ( buffer_det(N_int,2,N_det)) + if (is_complex) then + allocate(buffer_coef_complex(N_det,N_states)) + else + allocate(buffer_coef(N_det,N_states)) + endif ! Backup determinants j=0 @@ -142,6 +170,17 @@ subroutine copy_H_apply_buffer_to_wf N_det_old = j ! Backup coefficients + if (is_complex) then + do k=1,N_states + j=0 + do i=1,N_det + if (pruned(i)) cycle ! Pruned determinants + j += 1 + buffer_coef_complex(j,k) = psi_coef_complex(i,k) + enddo + ASSERT ( j == N_det_old ) + enddo + else do k=1,N_states j=0 do i=1,N_det @@ -151,6 +190,7 @@ subroutine copy_H_apply_buffer_to_wf enddo ASSERT ( j == N_det_old ) enddo + endif ! Update N_det N_det = N_det_old @@ -170,14 +210,57 @@ subroutine copy_H_apply_buffer_to_wf ASSERT (sum(popcnt(psi_det(:,1,i))) == elec_alpha_num) ASSERT (sum(popcnt(psi_det(:,2,i))) == elec_beta_num ) enddo + if (is_complex) then + do k=1,N_states + do i=1,N_det_old + psi_coef_complex(i,k) = buffer_coef_complex(i,k) + enddo + enddo + else do k=1,N_states do i=1,N_det_old psi_coef(i,k) = buffer_coef(i,k) enddo enddo + endif ! Copy new buffers + if (is_complex) then + !$OMP PARALLEL DEFAULT(SHARED) & + !$OMP PRIVATE(j,k,i) FIRSTPRIVATE(N_det_old) & + !$OMP SHARED(N_int,H_apply_buffer,psi_det,psi_coef_complex,N_states,psi_det_size) + j=0 + !$ j=omp_get_thread_num() + do k=0,j-1 + N_det_old += H_apply_buffer(k)%N_det + enddo + do i=1,H_apply_buffer(j)%N_det + do k=1,N_int + psi_det(k,1,i+N_det_old) = H_apply_buffer(j)%det(k,1,i) + psi_det(k,2,i+N_det_old) = H_apply_buffer(j)%det(k,2,i) + enddo + ASSERT (sum(popcnt(psi_det(:,1,i+N_det_old))) == elec_alpha_num) + ASSERT (sum(popcnt(psi_det(:,2,i+N_det_old))) == elec_beta_num ) + enddo + do k=1,N_states + do i=1,H_apply_buffer(j)%N_det + psi_coef_complex(i+N_det_old,k) = H_apply_buffer(j)%coef_complex(i,k) + enddo + enddo + !$OMP BARRIER + H_apply_buffer(j)%N_det = 0 + !$OMP END PARALLEL + SOFT_TOUCH N_det psi_det psi_coef_complex + + logical :: found_duplicates + call remove_duplicates_in_psi_det(found_duplicates) + do k=1,N_states + call normalize(psi_coef_complex(1,k),N_det) + enddo + SOFT_TOUCH N_det psi_det psi_coef_complex + else + !$OMP PARALLEL DEFAULT(SHARED) & !$OMP PRIVATE(j,k,i) FIRSTPRIVATE(N_det_old) & !$OMP SHARED(N_int,H_apply_buffer,psi_det,psi_coef,N_states,psi_det_size) @@ -210,7 +293,8 @@ subroutine copy_H_apply_buffer_to_wf call normalize(psi_coef(1,k),N_det) enddo SOFT_TOUCH N_det psi_det psi_coef - + + endif end subroutine remove_duplicates_in_psi_det(found_duplicates) @@ -274,7 +358,30 @@ subroutine remove_duplicates_in_psi_det(found_duplicates) enddo !$OMP END DO !$OMP END PARALLEL - + + if (is_complex) then + if (found_duplicates) then + k=0 + do i=1,N_det + if (.not.duplicate(i)) then + k += 1 + psi_det(:,:,k) = psi_det_sorted_bit (:,:,i) + psi_coef_complex(k,:) = psi_coef_sorted_bit_complex(i,:) + else + if (sum(cdabs(psi_coef_sorted_bit_complex(i,:))) /= 0.d0 ) then + psi_coef_complex(k,:) = psi_coef_sorted_bit_complex(i,:) + endif + endif + enddo + N_det = k + psi_det_sorted_bit(:,:,1:N_det) = psi_det(:,:,1:N_det) + psi_coef_sorted_bit_complex(1:N_det,:) = psi_coef_complex(1:N_det,:) + TOUCH N_det psi_det psi_coef_complex psi_det_sorted_bit psi_coef_sorted_bit_complex c0_weight + endif + psi_det = psi_det_sorted + psi_coef_complex = psi_coef_sorted_complex + SOFT_TOUCH psi_det psi_coef_complex psi_det_sorted_bit psi_coef_sorted_bit_complex + else if (found_duplicates) then k=0 do i=1,N_det @@ -296,6 +403,7 @@ subroutine remove_duplicates_in_psi_det(found_duplicates) psi_det = psi_det_sorted psi_coef = psi_coef_sorted SOFT_TOUCH psi_det psi_coef psi_det_sorted_bit psi_coef_sorted_bit + endif deallocate (duplicate,bit_tmp) end @@ -329,11 +437,19 @@ subroutine fill_H_apply_buffer_no_selection(n_selected,det_buffer,Nint,iproc) ASSERT (sum(popcnt(H_apply_buffer(iproc)%det(:,1,i+H_apply_buffer(iproc)%N_det)) )== elec_alpha_num) ASSERT (sum(popcnt(H_apply_buffer(iproc)%det(:,2,i+H_apply_buffer(iproc)%N_det))) == elec_beta_num) enddo + if (is_complex) then + do j=1,N_states + do i=1,N_selected + H_apply_buffer(iproc)%coef_complex(i+H_apply_buffer(iproc)%N_det,j) = (0.d0,0.d0) + enddo + enddo + else do j=1,N_states do i=1,N_selected H_apply_buffer(iproc)%coef(i+H_apply_buffer(iproc)%N_det,j) = 0.d0 enddo enddo + endif H_apply_buffer(iproc)%N_det = new_size do i=1,H_apply_buffer(iproc)%N_det ASSERT (sum(popcnt(H_apply_buffer(iproc)%det(:,1,i)) )== elec_alpha_num) @@ -343,6 +459,11 @@ subroutine fill_H_apply_buffer_no_selection(n_selected,det_buffer,Nint,iproc) end subroutine push_pt2(zmq_socket_push,pt2,norm_pert,H_pert_diag,i_generator,N_st,task_id) + !todo: modify for complex + if (is_complex) then + print*,irp_here,' not implemented for complex' + stop -1 + endif use f77_zmq implicit none BEGIN_DOC @@ -404,6 +525,11 @@ IRP_ENDIF end subroutine pull_pt2(zmq_socket_pull,pt2,norm_pert,H_pert_diag,i_generator,N_st,n,task_id) + !todo: modify for complex + if (is_complex) then + print*,irp_here,' not implemented for complex' + stop -1 + endif use f77_zmq implicit none BEGIN_DOC diff --git a/src/determinants/h_apply_nozmq.template.f b/src/determinants/h_apply_nozmq.template.f index bd261bbe..6d769556 100644 --- a/src/determinants/h_apply_nozmq.template.f +++ b/src/determinants/h_apply_nozmq.template.f @@ -17,8 +17,11 @@ subroutine $subroutine($params_main) double precision, allocatable :: fock_diag_tmp(:,:) $initialization + if (is_complex) then + PROVIDE H_apply_buffer_allocated mo_two_e_integrals_in_map psi_det_generators psi_coef_generators_complex + else PROVIDE H_apply_buffer_allocated mo_two_e_integrals_in_map psi_det_generators psi_coef_generators - + endif call wall_time(wall_0) diff --git a/src/determinants/occ_pattern.irp.f b/src/determinants/occ_pattern.irp.f index 9e0ccd33..8596a3c1 100644 --- a/src/determinants/occ_pattern.irp.f +++ b/src/determinants/occ_pattern.irp.f @@ -513,7 +513,7 @@ subroutine make_s2_eigenfunction N_det_new += 1 det_buffer(:,:,N_det_new) = d(:,:,j) if (N_det_new == bufsze) then - call fill_H_apply_buffer_no_selection(bufsze,det_buffer,N_int,ithread) + call fill_h_apply_buffer_no_selection(bufsze,det_buffer,N_int,ithread) N_det_new = 0 endif enddo @@ -528,7 +528,7 @@ subroutine make_s2_eigenfunction !$OMP END PARALLEL if (update) then - call copy_H_apply_buffer_to_wf + call copy_h_apply_buffer_to_wf if (is_complex) then TOUCH N_det psi_coef_complex psi_det psi_occ_pattern N_occ_pattern else diff --git a/src/generators_full/generators.irp.f b/src/generators_full/generators.irp.f index 7f18947f..f6a42fad 100644 --- a/src/generators_full/generators.irp.f +++ b/src/generators_full/generators.irp.f @@ -22,20 +22,35 @@ BEGIN_PROVIDER [ integer, N_det_generators ] call write_int(6,N_det_generators,'Number of generators') END_PROVIDER - BEGIN_PROVIDER [ integer(bit_kind), psi_det_generators, (N_int,2,psi_det_size) ] -&BEGIN_PROVIDER [ double precision, psi_coef_generators, (psi_det_size,N_states) ] +BEGIN_PROVIDER [ integer(bit_kind), psi_det_generators, (N_int,2,psi_det_size) ] implicit none BEGIN_DOC ! For Single reference wave functions, the generator is the ! Hartree-Fock determinant END_DOC psi_det_generators(1:N_int,1:2,1:N_det) = psi_det_sorted(1:N_int,1:2,1:N_det) - psi_coef_generators(1:N_det,1:N_states) = psi_coef_sorted(1:N_det,1:N_states) END_PROVIDER +BEGIN_PROVIDER [ double precision, psi_coef_generators, (psi_det_size,N_states) ] + implicit none + BEGIN_DOC + ! For Single reference wave functions, the generator is the + ! Hartree-Fock determinant + END_DOC + psi_coef_generators(1:N_det,1:N_states) = psi_coef_sorted(1:N_det,1:N_states) +END_PROVIDER + +BEGIN_PROVIDER [ complex*16, psi_coef_generators_complex, (psi_det_size,N_states) ] + implicit none + BEGIN_DOC + ! For Single reference wave functions, the generator is the + ! Hartree-Fock determinant + END_DOC + psi_coef_generators_complex(1:N_det,1:N_states) = psi_coef_sorted_complex(1:N_det,1:N_states) +END_PROVIDER + BEGIN_PROVIDER [ integer(bit_kind), psi_det_sorted_gen, (N_int,2,psi_det_size) ] -&BEGIN_PROVIDER [ double precision, psi_coef_sorted_gen, (psi_det_size,N_states) ] &BEGIN_PROVIDER [ integer, psi_det_sorted_gen_order, (psi_det_size) ] implicit none @@ -44,10 +59,26 @@ END_PROVIDER ! Hartree-Fock determinant END_DOC psi_det_sorted_gen = psi_det_sorted - psi_coef_sorted_gen = psi_coef_sorted psi_det_sorted_gen_order = psi_det_sorted_order END_PROVIDER +BEGIN_PROVIDER [ double precision, psi_coef_sorted_gen, (psi_det_size,N_states) ] + implicit none + BEGIN_DOC + ! For Single reference wave functions, the generator is the + ! Hartree-Fock determinant + END_DOC + psi_coef_sorted_gen = psi_coef_sorted +END_PROVIDER + +BEGIN_PROVIDER [ complex*16, psi_coef_sorted_gen_complex, (psi_det_size,N_states) ] + implicit none + BEGIN_DOC + ! For Single reference wave functions, the generator is the + ! Hartree-Fock determinant + END_DOC + psi_coef_sorted_gen_complex = psi_coef_sorted_complex +END_PROVIDER BEGIN_PROVIDER [integer, degree_max_generators] implicit none diff --git a/src/utils_complex/qp2-pbc-diff.txt b/src/utils_complex/qp2-pbc-diff.txt index cbe87d05..cdc7a7f7 100644 --- a/src/utils_complex/qp2-pbc-diff.txt +++ b/src/utils_complex/qp2-pbc-diff.txt @@ -22,6 +22,11 @@ determinants: (done) filter_connected.irp.f (done) fock_diag.irp.f (****) h_apply.irp.f + added coef_complex to h_apply_buffer_type + either coef or coef_complex will remain unallocated + (if this causes problems (it shouldn't), maybe just allocate unused one with size 1?) + check {push,pull}_pt2 + pt2, norm_pert, h_pert_diag types? (should be real? documentation?) (****) h_apply_nozmq.template.f (****) h_apply.template.f (****) h_apply_zmq.template.f