mpqc::lcao::gaussian::AOFactory< Tile, Policy > Class Template Reference

Documentation

template<typename Tile, typename Policy>
class mpqc::lcao::gaussian::AOFactory< Tile, Policy >

Gaussian AO integral factory using ab initio Hamiltonian.

This class computes Gaussian AO integrals exactly and using density fitting

Public Types

using AOFactoryBase = lcao::AOFactory< Tile, Policy >
 
using TArray = typename AOFactoryBase::TArray
 
using DirectTArray = typename AOFactoryBase::DirectTArray
 
using gtg_params_t = std::vector< std::pair< double, double > >
 
using Op = std::function< Tile(Tile &&)>
 

Public Member Functions

 AOFactory ()=default
 
 AOFactory (const KeyVal &kv)
 The KeyVal constructor. More...
 
 AOFactory (AOFactory &&)=default
 
AOFactoryoperator= (AOFactory &&)=default
 
virtual ~AOFactory () noexcept=default
 
TArray compute (const Formula &formula) override
 
DirectTArray compute_direct (const Formula &formula) override
 
const std::string & screen () const
 
double screen_threshold () const
 
double screen_threshold_push (double thr)
 
double screen_threshold_pop ()
 
std::shared_ptr< Screenerscreener (const std::wstring &formula) const
 
template<typename U , typename Bases >
TA::DistArray< Tile, std::enable_if_t< std::is_same< U, TA::SparsePolicy >::value, TA::SparsePolicy > > compute_integrals (madness::World &world, ShrPool< libint2::Engine > &engine, Bases &bases, std::shared_ptr< Screener > p_screen)
 
template<typename U , typename Bases >
TA::DistArray< Tile, std::enable_if_t< std::is_same< U, TA::DensePolicy >::value, TA::DensePolicy > > compute_integrals (madness::World &world, ShrPool< libint2::Engine > &engine, Bases &&bases, std::shared_ptr< Screener > p_screen)
 
template<typename U , typename Bases >
DirectArray< Tile, std::enable_if_t< std::is_same< U, TA::SparsePolicy >::value, TA::SparsePolicy > > compute_direct_integrals (madness::World &world, ShrPool< libint2::Engine > &engine, Bases &&bases, std::shared_ptr< Screener > p_screen, std::shared_ptr< const math::PetiteList > plist)
 
template<typename U , typename Bases >
DirectArray< Tile, std::enable_if_t< std::is_same< U, TA::DensePolicy >::value, TA::DensePolicy > > compute_direct_integrals (madness::World &world, ShrPool< libint2::Engine > &engine, Bases &bases, std::shared_ptr< Screener > p_screen, std::shared_ptr< const math::PetiteList > plist)
 

Member Typedef Documentation

◆ AOFactoryBase

template<typename Tile , typename Policy >
using mpqc::lcao::gaussian::AOFactory< Tile, Policy >::AOFactoryBase = lcao::AOFactory<Tile, Policy>

◆ DirectTArray

template<typename Tile , typename Policy >
using mpqc::lcao::gaussian::AOFactory< Tile, Policy >::DirectTArray = typename AOFactoryBase::DirectTArray

◆ gtg_params_t

template<typename Tile , typename Policy >
using mpqc::lcao::gaussian::AOFactory< Tile, Policy >::gtg_params_t = std::vector<std::pair<double, double> >

◆ Op

template<typename Tile , typename Policy >
using mpqc::lcao::gaussian::AOFactory< Tile, Policy >::Op = std::function<Tile(Tile&&)>

Op is a function pointer that convert Tile to Tile using Op = Tile (*) (Tile&&);

◆ TArray

template<typename Tile , typename Policy >
using mpqc::lcao::gaussian::AOFactory< Tile, Policy >::TArray = typename AOFactoryBase::TArray

Constructor & Destructor Documentation

◆ AOFactory() [1/3]

template<typename Tile , typename Policy >
mpqc::lcao::gaussian::AOFactory< Tile, Policy >::AOFactory ( )
default

◆ AOFactory() [2/3]

template<typename Tile , typename Policy >
mpqc::lcao::gaussian::AOFactory< Tile, Policy >::AOFactory ( const KeyVal kv)

The KeyVal constructor.

Parameters
kvthe KeyVal object, it will be queried for all keywords of the Factory ctor as well as the following additional keywords:
Keyword Type Default Description
screen string dist method of screening; valid values are schwarz (2-norm based), schwarz_inf (infinity-norm based), dist (2-norm based distance-dependent screeners: QQR for 4-center integrals and SQVl for 3-center integrals), dist_inf (same as dist but infinity-norm based), or empty string to turn screening off
screen_threshold real 1e-12 screening threshold for integrals (integral tiles with per-element norm smaller than this value wil be neglected)
precision real std::numeric_limits<double>::epsilon() integral precision
iterative_inv_sqrt bool false use iterative inverse square root
f12_param string ttg-6g[1] Te_n-no-type F12 correlation factor (defined if aux_basis exists in OrbitalBasisRegistry); valid format is ttg-Ng[A] where N and A are nonzero integer and positive real parameters defining the Slater-type correlation factor as a linear combination of N Gaussian geminals: $ - \exp(- A r_{12})/A \approx \sum\limits_{i=1}^N c_i \exp(- \alpha_i r_{12}^2) $.
See also
mpqc::lcao::f12::ttg_ng_fit
Warning
These keywords should be prefixed by "wfn_world:" (not recommended) or defined directly in this KeyVal object

◆ AOFactory() [3/3]

template<typename Tile , typename Policy >
mpqc::lcao::gaussian::AOFactory< Tile, Policy >::AOFactory ( AOFactory< Tile, Policy > &&  )
default

◆ ~AOFactory()

template<typename Tile , typename Policy >
virtual mpqc::lcao::gaussian::AOFactory< Tile, Policy >::~AOFactory ( )
virtualdefaultnoexcept

Member Function Documentation

◆ compute()

template<typename Tile , typename Policy >
AOFactory< Tile, Policy >::TArray mpqc::lcao::gaussian::AOFactory< Tile, Policy >::compute ( const Formula formula)
override

◆ compute_direct()

template<typename Tile , typename Policy >
AOFactory< Tile, Policy >::DirectTArray mpqc::lcao::gaussian::AOFactory< Tile, Policy >::compute_direct ( const Formula formula)
override

Functions to compute direct integral

◆ compute_direct_integrals() [1/2]

template<typename Tile , typename Policy >
template<typename U , typename Bases >
DirectArray<Tile, std::enable_if_t<std::is_same<U, TA::SparsePolicy>::value, TA::SparsePolicy> > mpqc::lcao::gaussian::AOFactory< Tile, Policy >::compute_direct_integrals ( madness::World &  world,
ShrPool< libint2::Engine > &  engine,
Bases &&  bases,
std::shared_ptr< Screener p_screen,
std::shared_ptr< const math::PetiteList >  plist 
)

◆ compute_direct_integrals() [2/2]

template<typename Tile , typename Policy >
template<typename U , typename Bases >
DirectArray<Tile, std::enable_if_t<std::is_same<U, TA::DensePolicy>::value, TA::DensePolicy> > mpqc::lcao::gaussian::AOFactory< Tile, Policy >::compute_direct_integrals ( madness::World &  world,
ShrPool< libint2::Engine > &  engine,
Bases &  bases,
std::shared_ptr< Screener p_screen,
std::shared_ptr< const math::PetiteList >  plist 
)

◆ compute_integrals() [1/2]

template<typename Tile , typename Policy >
template<typename U , typename Bases >
TA::DistArray<Tile, std::enable_if_t<std::is_same<U, TA::DensePolicy>::value, TA::DensePolicy> > mpqc::lcao::gaussian::AOFactory< Tile, Policy >::compute_integrals ( madness::World &  world,
ShrPool< libint2::Engine > &  engine,
Bases &&  bases,
std::shared_ptr< Screener p_screen 
)

◆ compute_integrals() [2/2]

template<typename Tile , typename Policy >
template<typename U , typename Bases >
TA::DistArray<Tile, std::enable_if_t<std::is_same<U, TA::SparsePolicy>::value, TA::SparsePolicy> > mpqc::lcao::gaussian::AOFactory< Tile, Policy >::compute_integrals ( madness::World &  world,
ShrPool< libint2::Engine > &  engine,
Bases &  bases,
std::shared_ptr< Screener p_screen 
)

◆ operator=()

template<typename Tile , typename Policy >
AOFactory& mpqc::lcao::gaussian::AOFactory< Tile, Policy >::operator= ( AOFactory< Tile, Policy > &&  )
default

◆ screen()

template<typename Tile , typename Policy >
const std::string& mpqc::lcao::gaussian::AOFactory< Tile, Policy >::screen ( ) const
inline

◆ screen_threshold()

template<typename Tile , typename Policy >
double mpqc::lcao::gaussian::AOFactory< Tile, Policy >::screen_threshold ( ) const
inline
Returns
the current value of screening threshold

◆ screen_threshold_pop()

template<typename Tile , typename Policy >
double mpqc::lcao::gaussian::AOFactory< Tile, Policy >::screen_threshold_pop

restores the previous value of the screening threshold

Returns
current value of screening threshold

◆ screen_threshold_push()

template<typename Tile , typename Policy >
double mpqc::lcao::gaussian::AOFactory< Tile, Policy >::screen_threshold_push ( double  thr)

updates the screening threshold value

Parameters
thrnew value of screening threshold
Returns
previous value of screening threshold

◆ screener()

template<typename Tile , typename Policy >
std::shared_ptr< Screener > mpqc::lcao::gaussian::AOFactory< Tile, Policy >::screener ( const std::wstring &  formula) const

The documentation for this class was generated from the following files: