complex determinants

src/determinants/occ_pattern.irp.f

@@ -529,8 +529,12 @@ subroutine make_s2_eigenfunction
 
   if (update) then
     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
     TOUCH N_det psi_coef psi_det psi_occ_pattern N_occ_pattern
     endif
+  endif
   call write_time(6)
 
 end

src/determinants/prune_wf.irp.f

@@ -25,7 +25,11 @@ BEGIN_PROVIDER [ logical, pruned, (N_det) ]
  else
 
   ndet_new = max(1,int( dble(N_det) * (1.d0 - pruning) + 0.5d0 ))
+  if (is_complex) then
+    thr = psi_average_norm_contrib_sorted_complex(ndet_new)
+  else
   thr = psi_average_norm_contrib_sorted(ndet_new)
+  endif
   do i=1, N_det
     pruned(i) = psi_average_norm_contrib(i) < thr
   enddo

src/determinants/psi_cas.irp.f

@@ -150,7 +150,20 @@ END_PROVIDER
  double precision :: hij,norm,u_dot_v
   psi_cas_energy = 0.d0
 
-
+  if (is_complex) then
+    complex*16 :: hij_c
+    do k = 1, N_states
+      norm = 0.d0
+      do i = 1, N_det_cas_complex
+        norm += cdabs(psi_cas_coef_complex(i,k) * psi_cas_coef_complex(i,k))
+        do j = 1, N_det_cas_complex
+          !TODO: accum imag parts to ensure that sum is zero?
+          psi_cas_energy(k) += dble(dconjg(psi_cas_coef_complex(i,k)) * psi_cas_coef_complex(j,k) * H_matrix_cas_complex(i,j))
+        enddo
+      enddo
+      psi_cas_energy(k) = psi_cas_energy(k) /norm
+    enddo
+  else
   do k = 1, N_states
    norm = 0.d0
    do i = 1, N_det_cas
@@ -161,6 +174,7 @@ END_PROVIDER
    enddo
    psi_cas_energy(k) = psi_cas_energy(k) /norm
   enddo
+  endif
 
 END_PROVIDER
 

src/determinants/psi_cas_complex.irp.f

@@ -0,0 +1,145 @@
+use bitmasks
+
+ BEGIN_PROVIDER [ integer(bit_kind), psi_cas_complex, (N_int,2,psi_det_size) ]
+&BEGIN_PROVIDER [ complex*16, psi_cas_coef_complex,  (psi_det_size,n_states) ]
+&BEGIN_PROVIDER [ integer, idx_cas_complex, (psi_det_size) ]
+&BEGIN_PROVIDER [ integer, N_det_cas_complex ]
+  implicit none
+  BEGIN_DOC
+  ! |CAS| wave function, defined from the application of the |CAS| bitmask on the
+  ! determinants. idx_cas gives the indice of the |CAS| determinant in psi_det.
+  END_DOC
+  integer                        :: i, k, l
+  logical                        :: good
+  n_det_cas_complex = 0
+  do i=1,N_det
+    do l = 1, N_states
+     psi_cas_coef_complex(i,l) = (0.d0,0.d0)
+    enddo
+    good = .True.
+    do k=1,N_int
+      good = good .and. (                                          &
+          iand(not(act_bitmask(k,1)), psi_det(k,1,i)) ==         &
+          iand(not(act_bitmask(k,1)), hf_bitmask(k,1)) ) .and. (  &
+          iand(not(act_bitmask(k,2)), psi_det(k,2,i)) ==         &
+          iand(not(act_bitmask(k,2)), hf_bitmask(k,2)) )
+    enddo
+    if (good) then
+      exit
+    endif
+    if (good) then
+      n_det_cas_complex = n_det_cas_complex+1
+      do k=1,N_int
+        psi_cas_complex(k,1,n_det_cas_complex) = psi_det(k,1,i)
+        psi_cas_complex(k,2,n_det_cas_complex) = psi_det(k,2,i)
+      enddo
+      idx_cas(n_det_cas_complex) = i
+      do k=1,N_states
+        psi_cas_coef_complex(n_det_cas_complex,k) = psi_coef_complex(i,k)
+      enddo
+    endif
+  enddo
+  call write_int(6,n_det_cas_complex, 'Number of determinants in the CAS')
+
+END_PROVIDER
+
+
+ BEGIN_PROVIDER [ integer(bit_kind), psi_cas_sorted_bit_complex, (N_int,2,psi_det_size) ]
+&BEGIN_PROVIDER [ complex*16, psi_cas_coef_sorted_bit_complex, (psi_det_size,N_states) ]
+ implicit none
+ BEGIN_DOC
+ ! |CAS| determinants sorted to accelerate the search of a random determinant in the wave
+ ! function.
+ END_DOC
+ call sort_dets_by_det_search_key_complex(n_det_cas_complex, psi_cas_complex, psi_cas_coef_complex, size(psi_cas_coef_complex,1), &
+     psi_cas_sorted_bit_complex, psi_cas_coef_sorted_bit_complex, N_states)
+
+END_PROVIDER
+
+
+
+ BEGIN_PROVIDER [ integer(bit_kind), psi_non_cas_complex,  (N_int,2,psi_det_size) ]
+&BEGIN_PROVIDER [ complex*16, psi_non_cas_coef,_complex (psi_det_size,n_states) ]
+&BEGIN_PROVIDER [ integer, idx_non_cas_complex,  (psi_det_size) ]
+&BEGIN_PROVIDER [ integer, N_det_non_cas_complex ]
+ implicit none
+ BEGIN_DOC
+  ! Set of determinants which are not part of the |CAS|, defined from the application
+  ! of the |CAS| bitmask on the determinants.
+  ! idx_non_cas gives the indice of the determinant in psi_det.
+ END_DOC
+ integer                        :: i_non_cas,j,k
+ integer                        :: degree
+ logical                        :: in_cas
+ i_non_cas =0
+ do k=1,N_det
+   in_cas = .False.
+   do j=1,N_det_cas_complex
+     call get_excitation_degree(psi_cas_complex(1,1,j), psi_det(1,1,k), degree, N_int)
+     if (degree == 0) then
+       in_cas = .True.
+       exit
+     endif
+   enddo
+   if (.not.in_cas) then
+     double precision :: hij
+     i_non_cas += 1
+     do j=1,N_int
+       psi_non_cas_complex(j,1,i_non_cas) = psi_det(j,1,k)
+       psi_non_cas_complex(j,2,i_non_cas) = psi_det(j,2,k)
+     enddo
+     do j=1,N_states
+       psi_non_cas_coef_complex(i_non_cas,j) = psi_coef_complex(k,j)
+     enddo
+     idx_non_cas_complex(i_non_cas) = k
+   endif
+ enddo
+ N_det_non_cas_complex = i_non_cas
+END_PROVIDER
+
+ BEGIN_PROVIDER [ integer(bit_kind), psi_non_cas_sorted_bit_complex, (N_int,2,psi_det_size) ]
+&BEGIN_PROVIDER [ complex*16, psi_non_cas_coef_sorted_bit_complex, (psi_det_size,N_states) ]
+ implicit none
+ BEGIN_DOC
+ ! |CAS| determinants sorted to accelerate the search of a random determinant in the wave
+ ! function.
+ END_DOC
+ !TODO: should this be n_det_non_cas_complex?
+ call sort_dets_by_det_search_key_complex(N_det_cas_complex, psi_non_cas_complex, psi_non_cas_coef_complex, size(psi_non_cas_coef_complex,1), &
+     psi_non_cas_sorted_bit_complex, psi_non_cas_coef_sorted_bit_complex, N_states)
+
+END_PROVIDER
+
+
+BEGIN_PROVIDER [complex*16, H_matrix_cas_complex, (N_det_cas_complex,N_det_cas_complex)]
+ implicit none
+ integer :: i,j
+ complex*16 :: hij
+  do i = 1, N_det_cas_complex
+   do j = 1, N_det_cas_complex
+    call i_h_j_complex(psi_cas_complex(1,1,i),psi_cas_complex(1,1,j),N_int,hij)
+    H_matrix_cas_complex(i,j) = hij
+   enddo
+  enddo
+END_PROVIDER
+
+ BEGIN_PROVIDER [complex*16, psi_coef_cas_diagonalized_complex, (N_det_cas_complex,N_states)]
+&BEGIN_PROVIDER [double precision, psi_cas_energy_diagonalized_complex, (N_states)]
+ implicit none
+ integer :: i,j
+  double precision, allocatable  :: eigenvalues(:)
+  complex*16, allocatable  :: eigenvectors(:,:)
+  allocate (eigenvectors(size(H_matrix_cas,1),N_det_cas))
+  allocate (eigenvalues(N_det_cas))
+  call lapack_diag_complex(eigenvalues,eigenvectors,                       &
+      H_matrix_cas_complex,size(H_matrix_cas_complex,1),N_det_cas_complex)
+  do i = 1, N_states
+   psi_cas_energy_diagonalized_complex(i) = eigenvalues(i)
+   do j = 1, N_det_cas_complex
+    psi_coef_cas_diagonalized_complex(j,i) = eigenvectors(j,i)
+   enddo
+  enddo
+
+
+ END_PROVIDER
+

src/utils_complex/qp2-pbc-diff.txt

@@ -13,12 +13,14 @@ determinants:
       ezfio_set_determinants_psi_coef_complex_qp_edit? (need ocaml?)
       psi_coef_{max,min}?
       save_wavefunction_specified{,_complex} qp_edit save?
-
+    psi_energy_mono_elec
     diag_h_mat_elem for complex
+    ...
 
   DONE
     create_excitations
       build_singly_excited_wavefunction{_complex}
+    ...
 
 
 -------------------------------------------------------------------------------------