src/AFQMC/Hamiltonians/HSPotential_Helpers.cpp
// In this case, c0/c1 refer to the ranges in the expanded list of CVs, e.g. 2*n/2*n+1
In this case, c0/c1 refer to the ranges in the expanded list of CVs, e.g.
In this case c0 / c1 refer to the ranges in the expanded list of CVs e . g .
void count_nnz(double cut, std::size_t& nik, std::size_t& nki, int c0, int c1, SpVType_shm_csr_matrix::reference const& Lik, SpVType_shm_csr_matrix::reference const& Lki)
{
  using ma::conj;
  using std::abs;
  using std::size_t;
  assert(c1>=c0);
  nik = nki = size_t(0);
  if(c0==c1) return;
  auto cik = std::lower_bound( to_address(Lik.non_zero_indices2_data()),
                               to_address((Lik.non_zero_indices2_data()
                                               +Lik.num_non_zero_elements())),
                               c0/2);
  auto cik_end = std::lower_bound( cik,
                                   to_address((Lik.non_zero_indices2_data()
                                               +Lik.num_non_zero_elements())),
                                   (c1+1)/2);
  auto vik = to_address(Lik.non_zero_values_data()) +
                    std::distance(to_address(Lik.non_zero_indices2_data()),cik);
  auto cki = std::lower_bound( to_address(Lki.non_zero_indices2_data()),
                               to_address((Lki.non_zero_indices2_data()
                                               +Lki.num_non_zero_elements())),
                               c0/2);
  auto cki_end = std::lower_bound( cki,
                                   to_address((Lki.non_zero_indices2_data()
                                               +Lki.num_non_zero_elements())),
                                   (c1+1)/2);
  auto vki = to_address(Lki.non_zero_values_data()) +
                    std::distance(to_address(Lki.non_zero_indices2_data()),cki);
  while( cik!=cik_end && cki!=cki_end ) {
    if( *cik == *cki ) { 
      if(abs(*vik + ma::conj(*vki)) > cut && 
          2*(*cik)>=c0 && 2*(*cik)<c1) { nik++; nki++; }
#if defined(QMC_COMPLEX)
      if(abs(*vik - ma::conj(*vki)) > cut && 
          2*(*cik)+1>=c0 && 2*(*cik)+1<c1) { nik++; nki++; }
#endif
      ++cik;
      ++vik;
      ++cki;
      ++vki;
    } else if( *cik < *cki )  { 
#if defined(QMC_COMPLEX)
      if(abs(*vik) > cut) {
        if(2*(*cik)>=c0 && 2*(*cik)<c1) { ++nik; ++nki; }
        if(2*(*cik)+1>=c0 && 2*(*cik)+1<c1) { ++nik; ++nki; }
      }
#else
      if(abs(*vik) > cut && 2*(*cik)>=c0 && 2*(*cik)<c1) { ++nik; ++nki; }
#endif
      ++cik;
      ++vik;
    } else {
#if defined(QMC_COMPLEX)
      if(abs(*vki) > cut) {
        if(2*(*cki)>=c0 && 2*(*cki)<c1) { ++nik; ++nki; }
        if(2*(*cki)+1>=c0 && 2*(*cki)+1<c1) { ++nik; ++nki; }
      }
#else
      if(abs(*vki) > cut && 2*(*cki)>=c0 && 2*(*cki)<c1) { ++nik; ++nki; }
#endif
      ++cki;
      ++vki;
    }
  }
  while( cik!=cik_end ) {
#if defined(QMC_COMPLEX)
    if(abs(*vik) > cut) {
      if(2*(*cik)>=c0 && 2*(*cik)<c1) { ++nik; ++nki; }
      if(2*(*cik)+1>=c0 && 2*(*cik)+1<c1) { ++nik; ++nki; }
    }
#else
    if(abs(*vik) > cut && 2*(*cik)>=c0 && 2*(*cik)<c1) { ++nik; ++nki; }
#endif
    ++cik;
    ++vik;
  }
  while( cki!=cki_end ) {
#if defined(QMC_COMPLEX)
    if(abs(*vki) > cut) {
      if(2*(*cki)>=c0 && 2*(*cki)<c1) { ++nik; ++nki; }
      if(2*(*cki)+1>=c0 && 2*(*cki)+1<c1) { ++nik; ++nki; }
    }
#else
    if(abs(*vki) > cut && 2*(*cki)>=c0 && 2*(*cki)<c1) { ++nik; ++nki; }
#endif
    ++cki;
    ++vki;
  }
}
void count_nnz ( double cut , std :: size_t & nik , std :: size_t & nki , int c0 , int c1 , SpVType_shm_csr_matrix :: reference const & Lik , SpVType_shm_csr_matrix :: reference const & Lki ) { using ma :: conj ; using std :: abs ; using std :: size_t ; assert ( c1 >= c0 ) ; nik = nki = size_t ( 0 ) ; if ( c0 == c1 ) return ; auto cik = std :: lower_bound ( to_address ( Lik . non_zero_indices2_data ( ) ) , to_address ( ( Lik . non_zero_indices2_data ( ) + Lik . num_non_zero_elements ( ) ) ) , c0 / 2 ) ; auto cik_end = std :: lower_bound ( cik , to_address ( ( Lik . non_zero_indices2_data ( ) + Lik . num_non_zero_elements ( ) ) ) , ( c1 + 1 ) / 2 ) ; auto vik = to_address ( Lik . non_zero_values_data ( ) ) + std :: distance ( to_address ( Lik . non_zero_indices2_data ( ) ) , cik ) ; auto cki = std :: lower_bound ( to_address ( Lki . non_zero_indices2_data ( ) ) , to_address ( ( Lki . non_zero_indices2_data ( ) + Lki . num_non_zero_elements ( ) ) ) , c0 / 2 ) ; auto cki_end = std :: lower_bound ( cki , to_address ( ( Lki . non_zero_indices2_data ( ) + Lki . num_non_zero_elements ( ) ) ) , ( c1 + 1 ) / 2 ) ; auto vki = to_address ( Lki . non_zero_values_data ( ) ) + std :: distance ( to_address ( Lki . non_zero_indices2_data ( ) ) , cki ) ; while ( cik != cik_end && cki != cki_end ) { if ( * cik == * cki ) { if ( abs ( * vik + ma :: conj ( * vki ) ) > cut && 2 * ( * cik ) >= c0 && 2 * ( * cik ) < c1 ) { nik ++ ; nki ++ ; } #if defined ( QMC_COMPLEX ) if ( abs ( * vik - ma :: conj ( * vki ) ) > cut && 2 * ( * cik ) + 1 >= c0 && 2 * ( * cik ) + 1 < c1 ) { nik ++ ; nki ++ ; } #endif ++ cik ; ++ vik ; ++ cki ; ++ vki ; } else if ( * cik < * cki ) { #if defined ( QMC_COMPLEX ) if ( abs ( * vik ) > cut ) { if ( 2 * ( * cik ) >= c0 && 2 * ( * cik ) < c1 ) { ++ nik ; ++ nki ; } if ( 2 * ( * cik ) + 1 >= c0 && 2 * ( * cik ) + 1 < c1 ) { ++ nik ; ++ nki ; } } #else if ( abs ( * vik ) > cut && 2 * ( * cik ) >= c0 && 2 * ( * cik ) < c1 ) { ++ nik ; ++ nki ; } #endif ++ cik ; ++ vik ; } else { #if defined ( QMC_COMPLEX ) if ( abs ( * vki ) > cut ) { if ( 2 * ( * cki ) >= c0 && 2 * ( * cki ) < c1 ) { ++ nik ; ++ nki ; } if ( 2 * ( * cki ) + 1 >= c0 && 2 * ( * cki ) + 1 < c1 ) { ++ nik ; ++ nki ; } } #else if ( abs ( * vki ) > cut && 2 * ( * cki ) >= c0 && 2 * ( * cki ) < c1 ) { ++ nik ; ++ nki ; } #endif ++ cki ; ++ vki ; } } while ( cik != cik_end ) { #if defined ( QMC_COMPLEX ) if ( abs ( * vik ) > cut ) { if ( 2 * ( * cik ) >= c0 && 2 * ( * cik ) < c1 ) { ++ nik ; ++ nki ; } if ( 2 * ( * cik ) + 1 >= c0 && 2 * ( * cik ) + 1 < c1 ) { ++ nik ; ++ nki ; } } #else if ( abs ( * vik ) > cut && 2 * ( * cik ) >= c0 && 2 * ( * cik ) < c1 ) { ++ nik ; ++ nki ; } #endif ++ cik ; ++ vik ; } while ( cki != cki_end ) { #if defined ( QMC_COMPLEX ) if ( abs ( * vki ) > cut ) { if ( 2 * ( * cki ) >= c0 && 2 * ( * cki ) < c1 ) { ++ nik ; ++ nki ; } if ( 2 * ( * cki ) + 1 >= c0 && 2 * ( * cki ) + 1 < c1 ) { ++ nik ; ++ nki ; } } #else if ( abs ( * vki ) > cut && 2 * ( * cki ) >= c0 && 2 * ( * cki ) < c1 ) { ++ nik ; ++ nki ; } #endif ++ cki ; ++ vki ; } }
In this case, c0/c1 refer to the ranges in the expanded list of CVs, e.g.
In this case c0 / c1 refer to the ranges in the expanded list of CVs e . g .
// ignoring factor of 0.5 to keep it consistent with the old code for now // both Lik and Lki have components on Chol Vec *cik==*cki // Lik + Lki* and Lki + Lik* // Lik - Lki* and Lki - Lik* // not on the same chol vector, only operate on the smallest