LIBINT 2.9.0
Classes | Typedefs | Functions
libint2::detail Namespace Reference

Most basic type – TwoPRep_11_11 – has one bfs for each particle in bra and ket. More...

Classes

struct  __initializer
 
class  compressed_pair
 
struct  CoreEvalScratch
 some evaluators need thread-local scratch, but most don't More...
 
struct  CoreEvalScratch< GaussianGmEval< Real, -1 > >
 GaussianGmEval<Real,-1> needs extra scratch data. More...
 
class  ext_stack_allocator
 allocator that uses an externally-managed stack-allocated array for allocations up to max_size, for larger allocations uses heap. More...
 
struct  has_static_size
 
struct  has_static_size< std::array< T, N > >
 
struct  has_static_size< T[N]>
 
struct  IsSharedPtr
 
struct  IsSharedPtrHelper
 Can be used to determine whether a type is a std::shared_ptr. More...
 
struct  IsSharedPtrHelper< std::shared_ptr< T > >
 
class  managed_singleton
 
struct  scale
 
struct  scale< Real, 2 >
 
struct  scale< Real, 4 >
 

Typedefs

template<typename Base , typename T >
using disable_if_same_or_derived
 

Functions

template<typename Bra >
bool is_nonderiv_ss_product (Bra &&bra)
 
constexpr double pow10 (long exp)
 
template<typename Real >
std::vector< Real > make_df_of_i_minux_1 (int imax)
 
template<typename Real >
std::vector< std::vector< Real > > make_cart_coeffs (int lmax)
 
int notxyz (int a, int b)
 
std::pair< int, int > notxyz (int a)
 
template<typename Real >
void subshell_occvec (Real *&occvec, size_t size, size_t &ne)
 computes orbital occupation numbers for a subshell of size size created by smearing no more than ne electrons (corresponds to spherical averaging)
 
std::vector< Shelluncontract (const std::vector< Shell > &shells)
 Uncontracts a set of shells.
 
double gamma_function (const double x)
 returns $ \Gamma(x) $
 
double alpha_eff (const Shell &shell1, const Shell &shell2, const int L)
 Computes an effective exponent for a product of two primitive gaussians as as described in J.
 
std::vector< Shellproduct_functions (const std::vector< Shell > &primitive_shells)
 Creates a set of product functions from a set of primitive shells.
 
std::vector< size_t > map_shell_to_basis_function (const std::vector< libint2::Shell > &shells)
 returns a hash map of shell indices to basis function indices
 
Eigen::MatrixXd compute_coulomb_matrix (const std::vector< Shell > &shells)
 Computes the Coulomb matrix $ (\mu|rij^{-1}|\nu) $ for a set of shells.
 
std::vector< std::vector< Shell > > split_by_L (const std::vector< Shell > &shells)
 Splits a set of shells by angular momentum.
 
std::vector< Shellpivoted_cholesky_in_L (const std::vector< Shell > &shells, const double cholesky_threshold)
 Computes the reduced set of product functions via pivoted Cholesky decomposition as described in J.
 
__libint2_engine_inline std::vector< Engine::compute2_ptr_type > init_compute2_ptrs ()
 
std::atomic< bool > & verbose_accessor ()
 
std::ostream *& verbose_stream_accessor ()
 
std::atomic< SHGShellOrdering > & solid_harmonics_ordering_accessor ()
 
ScreeningMethoddefault_screening_method_accessor ()
 
template<class T , std::size_t N>
bool operator== (const ext_stack_allocator< T, N > &x, const ext_stack_allocator< T, N > &y) noexcept
 
template<class T , std::size_t N>
bool operator!= (const ext_stack_allocator< T, N > &x, const ext_stack_allocator< T, N > &y) noexcept
 
template<class T1 , class T2 >
compressed_pair< T1, T2 > make_compressed_pair (const T1 &x, const T2 &y)
 
template<class C >
constexpr auto size (const C &c) -> decltype(c.size())
 
template<class T , std::size_t N>
constexpr std::size_t size (const T(&array)[N]) noexcept
 
template<typename T , typename A >
void resize (std::vector< T, A > &x, std::size_t n)
 

Detailed Description

Most basic type – TwoPRep_11_11 – has one bfs for each particle in bra and ket.

Note that GenIntegralSet is initialized with an abstract type libint2::BFSet, from which BFS derives.

Typedef Documentation

◆ disable_if_same_or_derived

template<typename Base , typename T >
using libint2::detail::disable_if_same_or_derived
Initial value:
typename std::enable_if<
!std::is_base_of<Base, typename std::decay<T>::type>::value>::type

Function Documentation

◆ alpha_eff()

double libint2::detail::alpha_eff ( const Shell & shell1,
const Shell & shell2,
const int L )
inline

Computes an effective exponent for a product of two primitive gaussians as as described in J.

Chem. Theory Comput. 2021, 17, 6886−6900. $ \alpha_{eff} = \left[ \frac{\Gamma(L+2)\Gamma(l_\mu +l_\nu +
\frac{3}{2})}{\Gamma(L+ \frac{3}{2})\Gamma(l_\mu +l_\nu + 2)} \right]^2
(\alpha_\mu + \alpha_\nu) $

Parameters
shell1first primitive shell
shell2second primitive shell
Ltotal angular momentum of product function
Returns
effective exponent of product function

References libint2::Shell::alpha, libint2::Shell::contr, and gamma_function().

Referenced by product_functions().

◆ compute_coulomb_matrix()

Eigen::MatrixXd libint2::detail::compute_coulomb_matrix ( const std::vector< Shell > & shells)
inline

Computes the Coulomb matrix $ (\mu|rij^{-1}|\nu) $ for a set of shells.

Parameters
shellsset of shells
Returns
Coulomb matrix

References libint2::Conservative, libint2::default_params(), and map_shell_to_basis_function().

Referenced by pivoted_cholesky_in_L().

◆ pivoted_cholesky_in_L()

std::vector< Shell > libint2::detail::pivoted_cholesky_in_L ( const std::vector< Shell > & shells,
const double cholesky_threshold )
inline

Computes the reduced set of product functions via pivoted Cholesky decomposition as described in J.

Chem. Theory Comput. 2021, 17, 6886−6900. Cholesky decomposition is performed on the Coulomb matrix of the product functions and the pivot indices are sorted in ascending order of column sums of the Coulomb matrix. The reduced set of product functions is then constructed by selecting the product functions corresponding to the pivot indices.

Parameters
shellsset of shells
cholesky_thresholdthreshold for choosing a product function via pivoted Cholesky decomposition
Returns
reduced set of product functions

References compute_coulomb_matrix(), and libint2::pivoted_cholesky().

Referenced by libint2::DFBasisSetGenerator::reduced_shells().

◆ product_functions()

std::vector< Shell > libint2::detail::product_functions ( const std::vector< Shell > & primitive_shells)
inline

Creates a set of product functions from a set of primitive shells.

Each pair of primitive shells produces a set of product functions with an effective exponent ( $ \alpha_{eff} $ as described above) and angular momentum ranging from $ |l1-l2| \leq L \leq l1+l2 $ as described in J. Chem. Theory Comput. 2021, 17, 6886−6900.

Parameters
primitive_shellsset of primitive shells
Returns
set of product functions

References alpha_eff(), and product_functions().

Referenced by libint2::DFBasisSetGenerator::DFBasisSetGenerator(), libint2::DFBasisSetGenerator::DFBasisSetGenerator(), and product_functions().

◆ split_by_L()

std::vector< std::vector< Shell > > libint2::detail::split_by_L ( const std::vector< Shell > & shells)
inline

Splits a set of shells by angular momentum.

Parameters
shellsset of shells
Returns
vector of vectors of shells split by angular momentum L. The i-th vector contains all shells with angular momentum i

Referenced by libint2::DFBasisSetGenerator::reduced_shells().

◆ subshell_occvec()

template<typename Real >
void libint2::detail::subshell_occvec ( Real *& occvec,
size_t size,
size_t & ne )

computes orbital occupation numbers for a subshell of size size created by smearing no more than ne electrons (corresponds to spherical averaging)

Parameters
[in,out]occvecoccupation vector, increments by size on return
[in]sizethe size of the subshell
[in,out]nethe number of electrons, on return contains the number of "remaining" electrons

Referenced by libint2::sto3g_ao_occupation_vector().

◆ uncontract()

std::vector< Shell > libint2::detail::uncontract ( const std::vector< Shell > & shells)
inline

Uncontracts a set of shells.

Parameters
[in]shellsa vector of shells to be uncontracted
Returns
a vector of uncontracted shells

Referenced by libint2::DFBasisSetGenerator::DFBasisSetGenerator(), and libint2::DFBasisSetGenerator::DFBasisSetGenerator().