Accelerated Davidson with connection map

src/CISD/cisd.irp.f

@@ -5,6 +5,9 @@ program cisd
   print *,   'HF      = ', HF_energy
   print *,  'N_states = ', N_states
   call H_apply_cisd
+! do i=1,N_det
+!   print '(100(X,O32))',  det_connections(:,i)
+! enddo
   print *,  'N_det = ', N_det
   do i = 1,N_states
    print *,  'energy  = ',CI_energy(i) 

src/Dets/README.rst

@@ -51,9 +51,10 @@ Documentation
 .. NEEDED_MODULES file.
 
 `copy_h_apply_buffer_to_wf <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/H_apply.irp.f#L113>`_
-  Undocumented
+  Copies the H_apply buffer to psi_coef. You need to touch psi_det, psi_coef and N_det
+  after calling this function.
 
-`fill_h_apply_buffer_no_selection <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/H_apply.irp.f#L197>`_
+`fill_h_apply_buffer_no_selection <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/H_apply.irp.f#L199>`_
   Fill the H_apply buffer with determiants for CISD
 
 `h_apply_buffer_allocated <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/H_apply.irp.f#L14>`_
@@ -69,18 +70,18 @@ Documentation
 `connected_to_ref <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/connected_to_ref.irp.f#L1>`_
   Undocumented
 
-`det_is_not_or_may_be_in_ref <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/connected_to_ref.irp.f#L196>`_
+`det_is_not_or_may_be_in_ref <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/connected_to_ref.irp.f#L188>`_
   If true, det is not in ref
   If false, det may be in ref
 
-`key_pattern_not_in_ref <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/connected_to_ref.irp.f#L229>`_
+`key_pattern_not_in_ref <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/connected_to_ref.irp.f#L222>`_
   Min and max values of the integers of the keys of the reference
 
-`davidson_converged <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/davidson.irp.f#L375>`_
+`davidson_converged <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/davidson.irp.f#L383>`_
   True if the Davidson algorithm is converged
 
-`davidson_criterion <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/davidson.irp.f#L365>`_
+`davidson_criterion <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/davidson.irp.f#L373>`_
-  Can be : [  energy  | residual | both ]
+  Can be : [  energy  | residual | both | wall_time | cpu_time | iterations ]
 
 `davidson_diag <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/davidson.irp.f#L18>`_
   Davidson diagonalization.
@@ -126,24 +127,41 @@ Documentation
 `davidson_sze_max <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/davidson.irp.f#L9>`_
   Max number of Davidson sizes
 
-`davidson_threshold <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/davidson.irp.f#L366>`_
+`davidson_threshold <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/davidson.irp.f#L374>`_
-  Can be : [  energy  | residual | both ]
+  Can be : [  energy  | residual | both | wall_time | cpu_time | iterations ]
 
-`n_det <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L11>`_
+`n_det <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L20>`_
   Number of determinants in the wave function
 
+`n_det_reference <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L75>`_
+  Number of determinants in the reference wave function
+
 `n_states <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L3>`_
   Number of states to consider
 
-`psi_coef <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L28>`_
+`psi_average_norm_contrib <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L84>`_
-  The wave function. Initialized with Hartree-Fock
+  Contribution of determinants to the state-averaged density
 
-`psi_det <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L27>`_
+`psi_average_norm_contrib_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L105>`_
-  The wave function. Initialized with Hartree-Fock
+  Wave function sorted by determinants (state-averaged)
 
-`psi_det_size <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L19>`_
+`psi_coef <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L47>`_
+  The wave function. Initialized with Hartree-Fock if the EZFIO file
+  is empty
+
+`psi_coef_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L104>`_
+  Wave function sorted by determinants (state-averaged)
+
+`psi_det <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L46>`_
+  The wave function. Initialized with Hartree-Fock if the EZFIO file
+  is empty
+
+`psi_det_size <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L38>`_
   Size of the psi_det/psi_coef arrays
 
+`psi_det_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants.irp.f#L103>`_
+  Wave function sorted by determinants (state-averaged)
+
 `double_exc_bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/determinants_bitmasks.irp.f#L40>`_
   double_exc_bitmask(:,1,i) is the bitmask for holes of excitation 1
   double_exc_bitmask(:,2,i) is the bitmask for particles of excitation 1
@@ -162,6 +180,22 @@ Documentation
   single_exc_bitmask(:,2,i) is the bitmask for particles
   for a given couple of hole/particle excitations i.
 
+`ci_eigenvectors <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/diagonalize_CI.irp.f#L36>`_
+  Eigenvectors/values of the CI matrix
+
+`ci_electronic_energy <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/diagonalize_CI.irp.f#L35>`_
+  Eigenvectors/values of the CI matrix
+
+`ci_energy <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/diagonalize_CI.irp.f#L18>`_
+  N_states lowest eigenvalues of the CI matrix
+
+`diag_algorithm <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/diagonalize_CI.irp.f#L1>`_
+  Diagonalization algorithm (Davidson or Lapack)
+
+`diagonalize_ci <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/diagonalize_CI.irp.f#L73>`_
+  Replace the coefficients of the CI states by the coefficients of the
+  eigenstates of the CI matrix
+
 `filter_connected <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/filter_connected.irp.f#L2>`_
   Filters out the determinants that are not connected by H
   .br
@@ -173,7 +207,9 @@ Documentation
   .br
   idx(0) is the number of determinants that interact with key1
 
-`filter_connected_i_h_psi0 <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/filter_connected.irp.f#L94>`_
+`filter_connected_davidson <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/filter_connected.irp.f#L95>`_
+  Filters out the determinants that are not connected by H
+  .br
   returns the array idx which contains the index of the
   .br
   determinants in the array key1 that interact
@@ -182,7 +218,16 @@ Documentation
   .br
   idx(0) is the number of determinants that interact with key1
 
-`filter_connected_i_h_psi0_sc2 <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/filter_connected.irp.f#L193>`_
+`filter_connected_i_h_psi0 <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/filter_connected.irp.f#L211>`_
+  returns the array idx which contains the index of the
+  .br
+  determinants in the array key1 that interact
+  .br
+  via the H operator with key2.
+  .br
+  idx(0) is the number of determinants that interact with key1
+
+`filter_connected_i_h_psi0_sc2 <http://github.com/LCPQ/quantum_package/tree/master/src/Dets/filter_connected.irp.f#L310>`_
   standard filter_connected_i_H_psi but returns in addition
   .br
   the array of the index of the non connected determinants to key1

src/Dets/connections.irp.f

@@ -0,0 +1,41 @@
+use bitmasks
+BEGIN_PROVIDER [ integer, N_con_int ]
+  implicit none
+  BEGIN_DOC
+  ! Number of integers to represent the connections between determinants
+  END_DOC
+  N_con_int = 1 + ishft(N_det-1,-13)
+END_PROVIDER
+
+BEGIN_PROVIDER [ integer*8, det_connections, (N_con_int,N_det) ]
+  implicit none
+  BEGIN_DOC
+  !
+  END_DOC
+  integer                        :: i,j
+  integer                        :: degree
+  integer                        :: j_int, j_k, j_l
+  !$OMP PARALLEL DEFAULT (NONE) &
+  !$OMP SHARED(N_det, N_con_int, psi_det,N_int, det_connections) &
+  !$OMP PRIVATE(i,j_int,j_k,j_l,j,degree)
+  !$OMP DO SCHEDULE(guided)
+  do i=1,N_det
+    do j_int=1,N_con_int
+      det_connections(j_int,i) = 0_8
+      j_k = ishft(j_int-1,13)
+      do j_l = j_k,min(j_k+8191,N_det), 128
+        do j = j_l+1,min(j_l+128,i)
+          call get_excitation_degree(psi_det(1,1,i),psi_det(1,1,j),degree,N_int)
+          if (degree < 3) then
+            det_connections(j_int,i) = ibset( det_connections(j_int,i), iand(63,ishft(j_l,-7)) )
+            exit
+          endif
+        enddo
+      enddo
+    enddo
+  enddo
+  !$OMP ENDDO
+  !$OMP ENDPARALLEL
+
+END_PROVIDER
+

src/Dets/davidson.irp.f

@@ -114,6 +114,8 @@ subroutine davidson_diag_hjj(dets_in,u_in,H_jj,energies,dim_in,sze,N_st,Nint,iun
   double precision               :: to_print(2,N_st)
   double precision               :: cpu, wall
   
+  PROVIDE N_con_int det_connections
+
   call write_time(iunit)
   call wall_time(wall)
   call cpu_time(cpu)

src/Dets/filter_connected.irp.f

@@ -91,6 +91,123 @@ subroutine filter_connected(key1,key2,Nint,sze,idx)
   idx(0) = l-1
 end
 
+
+subroutine filter_connected_davidson(key1,key2,Nint,sze,idx)
+  use bitmasks
+  implicit none
+  BEGIN_DOC
+  ! Filters out the determinants that are not connected by H
+  !
+  ! returns the array idx which contains the index of the 
+  !
+  ! determinants in the array key1 that interact 
+  !
+  ! via the H operator with key2.
+  !
+  ! idx(0) is the number of determinants that interact with key1
+  END_DOC
+  integer, intent(in)            :: Nint, sze
+  integer(bit_kind), intent(in)  :: key1(Nint,2,sze)
+  integer(bit_kind), intent(in)  :: key2(Nint,2)
+  integer, intent(out)           :: idx(0:sze)
+  
+  integer                        :: i,j,l
+  integer                        :: degree_x2
+  integer :: j_int, j_start
+  integer*8 :: itmp
+  
+  PROVIDE N_con_int det_connections
+  ASSERT (Nint > 0)
+  ASSERT (sze >= 0)
+
+  l=1
+  
+  if (Nint==1) then
+
+    i = idx(0)
+    do j_int=1,N_con_int
+      itmp = det_connections(j_int,i)
+      do while (itmp /= 0_8)
+        j_start = ishft(j_int-1,13) + ishft(trailz(itmp),7)
+        do j = j_start+1, min(j_start+128,i-1)
+          degree_x2 = popcnt(xor( key1(1,1,j), key2(1,1))) +             &
+              popcnt(xor( key1(1,2,j), key2(1,2)))
+          if (degree_x2 < 5) then
+            idx(l) = j
+            l = l+1
+          endif
+        enddo
+        itmp = iand(itmp-1_8,itmp)
+      enddo
+    enddo
+
+ !  l=1
+ !  !DIR$ LOOP COUNT (1000)
+ !  do i=1,sze
+ !    degree_x2 =  popcnt(xor( key1(1,1,i), key2(1,1))) +            &
+ !        popcnt(xor( key1(2,1,i), key2(2,1))) 
+ !    if (degree_x2 < 5) then
+ !      if (idx(l) /= i) then
+ !      print *,  l, idx(l), i
+ !      endif
+ !      idx(l) = i
+ !      l = l+1
+ !    endif
+ !  enddo
+    
+  else if (Nint==2) then
+    
+    !DIR$ LOOP COUNT (1000)
+    do i=1,sze
+      degree_x2 =  popcnt(xor( key1(1,1,i), key2(1,1))) +            &
+          popcnt(xor( key1(2,1,i), key2(2,1))) +                     &
+          popcnt(xor( key1(1,2,i), key2(1,2))) +                     &
+          popcnt(xor( key1(2,2,i), key2(2,2)))
+      if (degree_x2 < 5) then
+        idx(l) = i
+        l = l+1
+      endif
+    enddo
+    
+  else if (Nint==3) then
+    
+    !DIR$ LOOP COUNT (1000)
+    do i=1,sze
+      degree_x2 = popcnt(xor( key1(1,1,i), key2(1,1))) +             &
+          popcnt(xor( key1(1,2,i), key2(1,2))) +                     &
+          popcnt(xor( key1(2,1,i), key2(2,1))) +                     &
+          popcnt(xor( key1(2,2,i), key2(2,2))) +                     &
+          popcnt(xor( key1(3,1,i), key2(3,1))) +                     &
+          popcnt(xor( key1(3,2,i), key2(3,2)))
+      if (degree_x2 < 5) then
+        idx(l) = i
+        l = l+1
+      endif
+    enddo
+    
+  else
+    
+    !DIR$ LOOP COUNT (1000)
+    do i=1,sze
+      degree_x2 = 0
+      !DEC$ LOOP COUNT MIN(4)
+      do j=1,Nint
+        degree_x2 = degree_x2+ popcnt(xor( key1(j,1,i), key2(j,1))) +&
+            popcnt(xor( key1(j,2,i), key2(j,2)))
+        if (degree_x2 > 4) then
+          exit
+        endif
+      enddo
+      if (degree_x2 <= 5) then
+        idx(l) = i
+        l = l+1
+      endif
+    enddo
+    
+  endif
+  idx(0) = l-1
+end
+
 subroutine filter_connected_i_H_psi0(key1,key2,Nint,sze,idx)
   use bitmasks
   BEGIN_DOC

src/Dets/slater_rules.irp.f

@@ -864,7 +864,8 @@ subroutine H_u_0(v_0,u_0,H_jj,n,keys_tmp,Nint)
   Vt = 0.d0
   !$OMP DO SCHEDULE(guided)
   do i=1,n
-        call filter_connected(keys_tmp,keys_tmp(1,1,i),Nint,i-1,idx)
+        idx(0) = i
+        call filter_connected_davidson(keys_tmp,keys_tmp(1,1,i),Nint,i-1,idx)
         do jj=1,idx(0)
           j = idx(jj)
           call i_H_j(keys_tmp(1,1,j),keys_tmp(1,1,i),Nint,hij)

src/NEEDED_MODULES

@@ -1 +1 @@
-AOs Bitmask Electrons Ezfio_files MOs Nuclei Output Utils Hartree_Fock BiInts MonoInts MOGuess Dets DensityMatrix CISD Perturbation SingleRefMethod CISD_selected Selectors_full MP2 Generators_full Full_CI
+AOs Bitmask Electrons Ezfio_files MOs Nuclei Output Utils Hartree_Fock BiInts MonoInts MOGuess Dets DensityMatrix CISD Perturbation SingleRefMethod CISD_selected Selectors_full MP2 Generators_full Full_CI TestFock

src/TestFock/README.rst

@@ -0,0 +1,31 @@
+===============
+TestFock Module
+===============
+
+* This is a test
+Documentation
+=============
+
+.. Do not edit this section. It was auto-generated from the
+.. NEEDED_MODULES file.
+
+
+
+Needed Modules
+==============
+
+.. Do not edit this section. It was auto-generated from the
+.. NEEDED_MODULES file.
+
+* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
+* `BiInts <http://github.com/LCPQ/quantum_package/tree/master/src/BiInts>`_
+* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
+* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
+* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
+* `Hartree_Fock <http://github.com/LCPQ/quantum_package/tree/master/src/Hartree_Fock>`_
+* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
+* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
+* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
+* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
+* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
+