mpqc::lcao::gaussian Namespace Reference

Namespaces

 detail
 
 precomputed
 
 utility
 

Classes

class  AOFactory
 Gaussian AO integral factory using ab initio Hamiltonian. More...
 
class  AtomicBasis
 
class  Basis
 
class  BasisData
 
class  BasisPairData
 
struct  DFTaskGemm
 
class  DirectDFIntegralBuilder
 
class  DirectIntegralBuilder
 
struct  FMMShellData
 This holds shell-specific data relevant to FMM and related methods. More...
 
struct  FMMShellPairData
 This holds shell-pair-specific data relevant to FMM and related methods. More...
 
class  GaussianFactory
 
class  IntegralBuilder
 Builds integrals from an array of bases and an integral engine pool. More...
 
class  Q2Matrix
 Class which holds shell set information for screening. More...
 
class  QQRScreen
 Implements QQR screener for shell-sets of integrals. More...
 
class  SchwarzScreen
 Class for Schwarz based screening. More...
 
struct  ShellData
 aggregate of the Libint2 and SQVl shell data More...
 
struct  ShellPairData
 aggregate of the Libint2 and SQVl shell-pair data More...
 
class  SQVlScreen
 Implements SQVlScreen screener for shell-sets of (0|12) integrals. More...
 
struct  SQVlShellData
 This holds shell-specific data relevant to the SVQ(l) integral estimator described in DOI 10.1063/1.4917519. More...
 
struct  SQVlShellPairData
 This holds shell-specific data relevant to the SVQ(l) integral estimator described in DOI 10.1063/1.4917519. More...
 

Typedefs

using Shell = libint2::Shell
 
using ShellVec = std::vector< Shell >
 
using OrbitalBasisRegistry = OrbitalRegistry< std::shared_ptr< const Basis > >
 
using shellpair_data_accessor_t = std::function< const libint2::ShellPair &(std::size_t, std::size_t)>
 
using q_vector = std::vector< std::pair< double, std::array< double, 3 > >>
 
using Shell = Basis::Shell
 
using ShellVec = std::vector< Shell >
 
template<typename E >
using ShrPool = std::shared_ptr< mpqc::utility::TSPool< E > >
 
template<unsigned long N>
using BasisRefArray = std::array< std::reference_wrapper< const Basis >, N >
 
template<unsigned long N>
using BasisShrArray = std::array< std::shared_ptr< const Basis >, N >
 
using BasisRefVector = std::vector< std::reference_wrapper< const Basis > >
 
using BasisShrVector = std::vector< std::shared_ptr< const Basis > >
 
template<typename T >
using OrdTileVec = std::vector< std::pair< unsigned long, T > >
 
template<unsigned long N>
using ShellVecPtrArray = std::array< ShellVec const *, N >
 

Functions

int64_t max_am (ShellVec const &)
 Returns the maximum angular momement of any shell in the vector. More...
 
int64_t max_nprim (ShellVec const &)
 Returns the maximum number of primatives of any shell in the vector. More...
 
int64_t nfunctions (ShellVec const &)
 Returns the maximum number of primatives of any shell in the vector. More...
 
std::vector< std::vector< libint2::Shell > > reblock_basis (std::vector< libint2::Shell > shells, std::size_t blocksize)
 
template<typename Tile , typename Policy >
AOFactory< Tile, Policy > & to_ao_factory (lcao::AOFactory< Tile, Policy > &factory)
 
template<typename Tile , typename Policy >
std::shared_ptr< AOFactory< Tile, Policy > > to_ao_factory (const std::shared_ptr< lcao::AOFactory< Tile, Policy >> &factory)
 
 MPQC_EXTERN_TEMPLATE (class AOFactory<>)
 
template<typename OrbSpace >
void to_molden (std::string_view fname_prefix, const OrbSpace &orbspace, const WavefunctionWorld &wfn_world)
 prints an orbital space to a Molden file More...
 
libint2::Operator to_libint2_operator (Operator::Type mpqc_oper)
 
libint2::any to_libint2_operator_params (Operator::Type mpqc_oper, const Molecule &molecule, const std::map< Operator::Type, libint2::any > *oper_params)
 
libint2::scalar_type to_libint2_scale_factor (Operator::Type mpqc_oper, const std::map< Operator::Type, libint2::any > *oper_params)
 
std::shared_ptr< const Basisindex_to_basis (const OrbitalBasisRegistry &basis_registry, const OrbitalIndex &index)
 given OrbitalIndex, find the corresponding basis More...
 
std::vector< std::vector< size_t > > parallel_compute_shellpair_list (madness::World &world, const Basis &basis1, const Basis &basis2, double threshold=1e-12, double engine_precision=0.0)
 This computes non-negligible shell pair list; shells i and j form a non-negligible pair if they share a center or the Frobenius norm of their overlap is greater than threshold. More...
 
template<typename Bases >
std::shared_ptr< ThresholdedScreenermake_screener (madness::World &world, const ShrPool< libint2::Engine > &engine, Bases &&bases, const std::string &screen, double screen_threshold)
 factory for ThresholdedScreener objects More...
 
template<typename Policy , bool is_real = true>
TA::DistArray< TA::TensorD, Policy > tensorZ_to_tensorD (const TA::DistArray< TA::TensorZ, Policy > &complex_array)
 This takes real or imaginary part from a complex array. More...
 
template<typename Tile , typename Policy >
TA::DistArray< Tile, Policy > compute_shellblock_norm (const Basis &bs0, const Basis &bs1, const TA::DistArray< Tile, Policy > &D)
 
template<typename Tile , typename Policy >
TA::DistArray< Tile, Policy > compute_shellblock_norm (const Basis &bs0, const TA::DistArray< Tile, Policy > &D)
 
template<typename Tile = TA::TensorD, typename ComputeTile = Tile, typename Engine = libint2::Engine, typename Bases >
DirectArray< Tile, TA::SparsePolicy > direct_sparse_integrals (madness::World &world, ShrPool< Engine > shr_pool, Bases &&bases, std::shared_ptr< Screener > screen=std::make_shared< Screener >(Screener{}), std::function< Tile(ComputeTile &&)> op=TA::detail::Noop< Tile, ComputeTile, true >(), std::function< shellpair_data_accessor_t(const Basis *, const Basis *)> make_shellpair_data_accessor={}, std::shared_ptr< const math::PetiteList > plist=math::PetiteList::make_trivial())
 Construct direct integral tensors in parallel with screening. More...
 
template<typename Tile = TA::TensorD, typename ComputeTile = Tile, typename Engine , typename Idx , typename Bases >
DirectArray< Tile, TA::SparsePolicy > direct_sparse_integrals (madness::World &world, ShrPool< Engine > shr_pool, Bases &&bases, std::vector< std::pair< Idx, float >> const &user_provided_norms, bool user_provided_norms_perelem, std::shared_ptr< Screener > screen=std::make_shared< Screener >(Screener{}), std::function< Tile(ComputeTile &&)> op=TA::detail::Noop< Tile, ComputeTile, true >(), std::function< shellpair_data_accessor_t(const Basis *, const Basis *)> make_shellpair_data_accessor={}, std::shared_ptr< const math::PetiteList > plist=math::PetiteList::make_trivial())
 Construct direct integral tensors in parallel with screening. More...
 
template<typename Tile = TA::TensorD, typename ComputeTile = Tile, typename Engine , typename Bases >
DirectArray< Tile, TA::SparsePolicy, Engine > untruncated_direct_sparse_integrals (madness::World &world, ShrPool< Engine > shr_pool, Bases &&bases, std::shared_ptr< Screener > screen=std::make_shared< Screener >(Screener{}), std::function< Tile(ComputeTile &&)> op=TA::detail::Noop< Tile, ComputeTile, true >(), std::function< shellpair_data_accessor_t(const Basis *, const Basis *)> make_shellpair_data_accessor={}, std::shared_ptr< const math::PetiteList > plist=math::PetiteList::make_trivial())
 Construct direct integral tensors in parallel with screening. More...
 
template<typename Tile = TA::TensorD, typename ComputeTile = Tile, typename Engine , typename Bases >
DirectArray< Tile, TA::DensePolicy > direct_dense_integrals (madness::World &world, ShrPool< Engine > shr_pool, Bases &&bases, std::shared_ptr< Screener > screen=std::make_shared< Screener >(Screener{}), std::function< Tile(ComputeTile &&)> op=TA::detail::Noop< Tile, ComputeTile, true >(), std::function< shellpair_data_accessor_t(const Basis *, const Basis *)> make_shellpair_data_accessor={}, std::shared_ptr< const math::PetiteList > plist=math::PetiteList::make_trivial())
 Construct direct dense integral tensors in parallel with screening. More...
 
template<typename Tile >
DirectArray< Tile, TA::DensePolicy > df_direct_integrals (const TA::DistArray< Tile, TA::DensePolicy > &bra, const TA::DistArray< Tile, TA::DensePolicy > &ket, Formula::Notation notation=Formula::Notation::Chemical)
 
template<typename Tile >
DirectArray< Tile, TA::SparsePolicy > df_direct_integrals (const TA::DistArray< Tile, TA::SparsePolicy > &bra, const TA::DistArray< Tile, TA::SparsePolicy > &ket, Formula::Notation notation=Formula::Notation::Chemical)
 
template<typename Tile , typename ComputeTile = Tile, typename Engine = libint2::Engine>
std::shared_ptr< IntegralBuilder< Tile, ComputeTile, Engine > > make_integral_builder (ShrPool< Engine > shr_epool, const BasisShrVector &bases, std::shared_ptr< Screener > shr_screen, std::function< Tile(ComputeTile &&)> op, std::function< shellpair_data_accessor_t(const Basis *, const Basis *)> make_shellpair_data_accessor={}, std::shared_ptr< const math::PetiteList > plist=math::PetiteList::make_trivial())
 Function to make detection of template parameters easier, see IntegralBuilder for details. More...
 
template<typename Tile , typename ComputeTile = Tile, typename Engine = libint2::Engine>
std::shared_ptr< DirectIntegralBuilder< Tile, ComputeTile, Engine > > make_direct_integral_builder (madness::World &world, ShrPool< Engine > shr_epool, const BasisShrVector &bases, std::shared_ptr< Screener > shr_screen, std::function< Tile(ComputeTile &&)> op, std::function< shellpair_data_accessor_t(const Basis *, const Basis *)> make_shellpair_data_accessor={}, std::shared_ptr< const math::PetiteList > plist=math::PetiteList::make_trivial())
 Function to make detection of template parameters easier, see DirectIntegralBuilder for details. More...
 
template<typename Bases >
libint2::Engine make_engine (const libint2::Operator &oper, Bases &&bases, libint2::BraKet braket=libint2::BraKet::invalid, libint2::any oper_params=libint2::any(), double engine_precision=std::numeric_limits< double >::epsilon(), const libint2::scalar_type &scalar=libint2::scalar_type(1.0))
 
q_vector make_q (Molecule const &mol)
 
template<typename Bases >
ShrPool< libint2::Engine > make_engine_pool (const libint2::Operator &oper, Bases &&bases, libint2::BraKet braket=libint2::BraKet::invalid, libint2::any oper_params=libint2::any(), const libint2::scalar_type &scalar=libint2::scalar_type(1.0), double engine_precision=std::numeric_limits< double >::epsilon())
 Make an pool of Engine objects. More...
 
template<typename E , typename Bases >
QQRScreen create_qqr_screener (madness::World &world, const ShrPool< E > &engs, Bases &&bases, double threshold, SchwarzScreen::norm_op_type norm_func=detail::l2_norm)
 
std::unique_ptr< SchwarzScreenmake_schwarz_screener (std::tuple< const Basis & > bra, std::tuple< const Basis &, const Basis & > ket, double threshold)
 
bool operator== (const Q2Matrix &first, const Q2Matrix &second)
 
template<typename E , typename Bases >
SchwarzScreen create_schwarz_screener (madness::World &world, ShrPool< E > const &eng, Bases &&bs_array, double thresh, SchwarzScreen::norm_op_type norm_func=detail::l2_norm)
 Creates a Schwarz Screener. More...
 
template<typename E , typename Bases >
SQVlScreen create_sqvl_screener (madness::World &world, const ShrPool< E > &engs, Bases &&bases, double threshold, SchwarzScreen::norm_op_type norm_func=detail::l2_norm)
 
template<typename E , typename Tile = TA::TensorD, typename Bases >
TA::DistArrayVector< Tile, TA::SparsePolicy > sparse_xyz_integrals (madness::World &world, ShrPool< E > shr_pool, Bases &&bases, std::function< Tile(Tile &&)> op=TA::detail::Noop< Tile, Tile, true >{})
 Construct set of sparse integral tensors in parallel. More...
 
template<typename E , typename Bases , typename Tile = TA::TensorD>
TA::DistArrayVector< Tile, TA::DensePolicy > dense_xyz_integrals (madness::World &world, ShrPool< E > shr_pool, Bases &&bases, std::function< Tile(Tile &&)> op=TA::detail::Noop< Tile, Tile, true >{})
 Construct dense integral tensors from sets in parallel. More...
 
template<typename Tile , typename Policy , typename E , typename Bases >
TA::DistArrayVector< Tile, std::enable_if_t< std::is_same< Policy, TA::DensePolicy >::value, TA::DensePolicy > > xyz_integrals (madness::World &world, ShrPool< E > shr_pool, Bases &&bases, std::function< Tile(Tile &&)> op=TA::detail::Noop< Tile, Tile, true >{})
 
template<typename Tile , typename Policy , typename E , typename Bases >
TA::DistArrayVector< Tile, std::enable_if_t< std::is_same< Policy, TA::SparsePolicy >::value, TA::SparsePolicy > > xyz_integrals (madness::World &world, ShrPool< E > shr_pool, Bases &&bases, std::function< Tile(Tile &&)> op=TA::detail::Noop< Tile, Tile, true >{})
 
template<typename E , typename Tile = TA::TensorD, typename Bases >
TA::DistArray< Tile, TA::SparsePolicy > sparse_integrals (madness::World &world, ShrPool< E > shr_pool, Bases &&bases, std::shared_ptr< Screener > screen=std::make_shared< Screener >(Screener{}), std::function< Tile(Tile &&)> op=TA::detail::Noop< Tile, Tile, true >{}, std::function< shellpair_data_accessor_t(const Basis *, const Basis *)> make_shellpair_data_accessor={}, std::shared_ptr< const math::PetiteList > plist=math::PetiteList::make_trivial())
 Construct sparse integral tensors in parallel. More...
 
template<typename E , typename Tile = TA::TensorD, typename Bases >
TA::DistArray< Tile, TA::DensePolicy > dense_integrals (madness::World &world, ShrPool< E > shr_pool, Bases &&bases, std::shared_ptr< Screener > screen=std::make_shared< Screener >(Screener{}), std::function< Tile(Tile &&)> op=TA::detail::Noop< Tile, Tile, true >{}, std::function< shellpair_data_accessor_t(const Basis *, const Basis *)> make_shellpair_data_accessor={}, std::shared_ptr< const math::PetiteList > plist=math::PetiteList::make_trivial())
 Construct a dense integral tensor in parallel. More...
 
template<typename Tile , typename Policy , typename Bases >
TA::DistArray< Tile, typename std::enable_if< std::is_same< Policy, TA::SparsePolicy >::value, TA::SparsePolicy >::type > compute_integrals (madness::World &world, ShrPool< libint2::Engine > &engine, Bases &&bases, std::shared_ptr< Screener > p_screen=std::make_shared< Screener >(Screener{}), std::function< Tile(Tile &&)> op=TA::detail::Noop< Tile, Tile, true >{})
 
template<typename Tile , typename Policy , typename Bases >
TA::DistArray< Tile, typename std::enable_if< std::is_same< Policy, TA::DensePolicy >::value, TA::DensePolicy >::type > compute_integrals (madness::World &world, ShrPool< libint2::Engine > &engine, Bases &&bases, std::shared_ptr< Screener > p_screen=std::make_shared< Screener >(Screener{}), std::function< Tile(Tile &&)> op=TA::detail::Noop< Tile, Tile, true >{})
 
template<typename Bases >
BasisShrVector to_basis_shr_vector (Bases &&bases)
 
template<typename Bases >
BasisRefVector to_basis_ref_vector (Bases &&bases)
 
template<typename BS >
const Basisto_basis_ref (BS &&bs)
 
template<typename Bases >
TA::TiledRange make_trange (Bases &&bases)
 
template<typename Tile , typename Array >
void set_array (std::vector< std::pair< unsigned long, Tile >> &tiles, Array &a)
 
template<typename Tile , typename Policy >
void set_array (std::vector< std::pair< unsigned long, Tile >> &tiles, TiledArray::DistArrayVector< Tile, Policy > &a)
 
template<typename Engs , typename Tile , typename ComputeTile , typename Policy >
void soad_task (Engs eng_pool, const std::array< int64_t, 2 > &tile_idx, const std::array< std::vector< libint2::Shell > const *, 2 > &obs_row_and_col_ptrs, std::vector< libint2::Shell > const *min_bs, const RowMatrixXd *D, TA::DistArray< Tile, Policy > *F, std::function< Tile(ComputeTile &&)> op, const std::array< shellpair_data_accessor_t, 3 > *shellpair_data_accessors, const std::array< std::shared_ptr< ThresholdedScreener >, 2 > *jk_screeners)
 
template<typename Tile , typename ComputeTile , typename Policy >
TA::DistArrayVector< Tile, typename std::enable_if< std::is_same< Policy, TA::SparsePolicy >::value, TA::SparsePolicy >::type > fock_from_soad (madness::World &world, Molecule const &clustered_mol, Basis const &obs, Basis const &minbs, TA::DistArray< Tile, Policy > const &H, int const n_components, std::function< shellpair_data_accessor_t(const Basis *bs0, const Basis *bs1)> make_shellpair_data_accessor={}, std::function< Tile(ComputeTile &&)> op=TA::detail::Noop< Tile, ComputeTile, true >())
 
template<typename Tile , typename ComputeTile , typename Policy >
TA::DistArrayVector< Tile, typename std::enable_if< std::is_same< Policy, TA::DensePolicy >::value, TA::DensePolicy >::type > fock_from_soad (madness::World &world, Molecule const &clustered_mol, Basis const &obs, Basis const &minbs, TA::DistArray< Tile, Policy > const &H, int const n_components, std::function< shellpair_data_accessor_t(const Basis *bs0, const Basis *bs1)> make_shellpair_data_accessor={}, std::function< Tile(ComputeTile &&)> op=TA::detail::Noop< Tile, ComputeTile, true >())
 
std::ostream & operator<< (std::ostream &os, Basis::Factory const &f)
 
std::ostream & operator<< (std::ostream &os, Basis const &b)
 
Eigen::RowVectorXi sub_basis_map (const Basis &basis, const Basis &sub_basis)
 
std::shared_ptr< const Basismerge (const Basis &basis1, const Basis &basis2)
 
Basis parallel_make_basis (madness::World &world, const Basis::Factory &factory, const mpqc::Molecule &mol)
 
std::ostream & operator<< (std::ostream &os, AtomicBasis const &b)
 
bool operator== (const Basis &one, const Basis &two)
 
template<typename Op , typename... Args>
Basis reblock (Basis const &basis, Op op, Args... args)
 reblock allows for reblocking a basis More...
 

Typedef Documentation

◆ OrbitalBasisRegistry

Typedef of OrbitalBasisRegistry A Registry that maps OrbitalIndex to a Gaussian Basis

◆ q_vector

using mpqc::lcao::gaussian::q_vector = typedef std::vector<std::pair<double, std::array<double, 3> >>

◆ Shell [1/2]

typedef libint2::Shell mpqc::lcao::gaussian::Shell

◆ shellpair_data_accessor_t

using mpqc::lcao::gaussian::shellpair_data_accessor_t = typedef std::function<const libint2::ShellPair &(std::size_t, std::size_t)>

◆ ShellVec [1/2]

typedef std::vector< Shell > mpqc::lcao::gaussian::ShellVec

Function Documentation

◆ compute_integrals() [1/2]

template<typename Tile , typename Policy , typename Bases >
TA::DistArray< Tile, typename std::enable_if<std::is_same<Policy, TA::SparsePolicy>::value, TA::SparsePolicy>::type> mpqc::lcao::gaussian::compute_integrals ( madness::World &  world,
ShrPool< libint2::Engine > &  engine,
Bases &&  bases,
std::shared_ptr< Screener p_screen = std::make_shared<Screener>(Screener{}),
std::function< Tile(Tile &&)>  op = TA::detail::Noop<Tile, Tile, true>{} 
)

◆ compute_integrals() [2/2]

template<typename Tile , typename Policy , typename Bases >
TA::DistArray< Tile, typename std::enable_if<std::is_same<Policy, TA::DensePolicy>::value, TA::DensePolicy>::type> mpqc::lcao::gaussian::compute_integrals ( madness::World &  world,
ShrPool< libint2::Engine > &  engine,
Bases &&  bases,
std::shared_ptr< Screener p_screen = std::make_shared<Screener>(Screener{}),
std::function< Tile(Tile &&)>  op = TA::detail::Noop<Tile, Tile, true>{} 
)

◆ compute_shellblock_norm() [1/2]

template<typename Tile , typename Policy >
TA::DistArray<Tile, Policy> mpqc::lcao::gaussian::compute_shellblock_norm ( const Basis bs0,
const Basis bs1,
const TA::DistArray< Tile, Policy > &  D 
)

◆ compute_shellblock_norm() [2/2]

template<typename Tile , typename Policy >
TA::DistArray<Tile, Policy> mpqc::lcao::gaussian::compute_shellblock_norm ( const Basis bs0,
const TA::DistArray< Tile, Policy > &  D 
)

◆ create_qqr_screener()

template<typename E , typename Bases >
QQRScreen mpqc::lcao::gaussian::create_qqr_screener ( madness::World &  world,
const ShrPool< E > &  engs,
Bases &&  bases,
double  threshold,
SchwarzScreen::norm_op_type  norm_func = detail::l2_norm 
)

◆ create_schwarz_screener()

template<typename E , typename Bases >
SchwarzScreen mpqc::lcao::gaussian::create_schwarz_screener ( madness::World &  world,
ShrPool< E > const &  eng,
Bases &&  bs_array,
double  thresh,
SchwarzScreen::norm_op_type  norm_func = detail::l2_norm 
)

Creates a Schwarz Screener.

Template Parameters
Ean engine type
Basesa directly-addressable sequence of basis sets (or reference proxies, e.g. std::reference_wrapper)
Parameters
worldis a reference to the world in which the distributed screener data will be computed and stores
engis a reference to a ShrPool<E>
bs_arrayis a reference to a directly-addressable sequence of basis sets, if the length is 3 then DF integrals are assumed and the first basis is assumed to be the auxiliary basis if the length is 4 then it is assumed that four center screening is desired. There is no requirement that any basis sets be the same.
threshis the Schwarz Screening threshold
norm_funcis the function pointer that returns a norm, it should have a signature double (Eigen::Map<const RowMatrixXd> const &)

◆ create_sqvl_screener()

template<typename E , typename Bases >
SQVlScreen mpqc::lcao::gaussian::create_sqvl_screener ( madness::World &  world,
const ShrPool< E > &  engs,
Bases &&  bases,
double  threshold,
SchwarzScreen::norm_op_type  norm_func = detail::l2_norm 
)

◆ dense_integrals()

template<typename E , typename Tile = TA::TensorD, typename Bases >
TA::DistArray<Tile, TA::DensePolicy> mpqc::lcao::gaussian::dense_integrals ( madness::World &  world,
ShrPool< E >  shr_pool,
Bases &&  bases,
std::shared_ptr< Screener screen = std::make_shared<Screener>(Screener{}),
std::function< Tile(Tile &&)>  op = TA::detail::Noop<Tile, Tile, true>{},
std::function< shellpair_data_accessor_t(const Basis *, const Basis *)>  make_shellpair_data_accessor = {},
std::shared_ptr< const math::PetiteList >  plist = math::PetiteList::make_trivial() 
)

Construct a dense integral tensor in parallel.

See also
sparse_integrals

◆ dense_xyz_integrals()

template<typename E , typename Bases , typename Tile = TA::TensorD>
TA::DistArrayVector<Tile, TA::DensePolicy> mpqc::lcao::gaussian::dense_xyz_integrals ( madness::World &  world,
ShrPool< E >  shr_pool,
Bases &&  bases,
std::function< Tile(Tile &&)>  op = TA::detail::Noop<Tile, Tile, true>{} 
)

Construct dense integral tensors from sets in parallel.

This is needed for integrals such as the dipole integrals that come as a set.

Parameters
shr_poolshould be a std::shared_ptr to an IntegralEnginePool
basesshould be a std::array of Basis, which will be copied.
opneeds to be a function or functor that takes a TA::TensorD && and returns any valid tile type. Op is copied so it can be moved.
auto t = [](TA::TensorD &&ten){return std::move(ten);};

◆ df_direct_integrals() [1/2]

template<typename Tile >
DirectArray<Tile, TA::DensePolicy> mpqc::lcao::gaussian::df_direct_integrals ( const TA::DistArray< Tile, TA::DensePolicy > &  bra,
const TA::DistArray< Tile, TA::DensePolicy > &  ket,
Formula::Notation  notation = Formula::Notation::Chemical 
)

construct direct dense LCAO integral from density fitting

◆ df_direct_integrals() [2/2]

template<typename Tile >
DirectArray<Tile, TA::SparsePolicy> mpqc::lcao::gaussian::df_direct_integrals ( const TA::DistArray< Tile, TA::SparsePolicy > &  bra,
const TA::DistArray< Tile, TA::SparsePolicy > &  ket,
Formula::Notation  notation = Formula::Notation::Chemical 
)

construct direct sparse LCAO integral from density fitting

◆ direct_dense_integrals()

template<typename Tile = TA::TensorD, typename ComputeTile = Tile, typename Engine , typename Bases >
DirectArray<Tile, TA::DensePolicy> mpqc::lcao::gaussian::direct_dense_integrals ( madness::World &  world,
ShrPool< Engine >  shr_pool,
Bases &&  bases,
std::shared_ptr< Screener screen = std::make_shared<Screener>(Screener{}),
std::function< Tile(ComputeTile &&)>  op = TA::detail::Noop<Tile, ComputeTile, true>(),
std::function< shellpair_data_accessor_t(const Basis *, const Basis *)>  make_shellpair_data_accessor = {},
std::shared_ptr< const math::PetiteList >  plist = math::PetiteList::make_trivial() 
)

Construct direct dense integral tensors in parallel with screening.

Same requirements on Op as those in Integral Builder

◆ direct_sparse_integrals() [1/2]

template<typename Tile = TA::TensorD, typename ComputeTile = Tile, typename Engine = libint2::Engine, typename Bases >
DirectArray<Tile, TA::SparsePolicy> mpqc::lcao::gaussian::direct_sparse_integrals ( madness::World &  world,
ShrPool< Engine >  shr_pool,
Bases &&  bases,
std::shared_ptr< Screener screen = std::make_shared<Screener>(Screener{}),
std::function< Tile(ComputeTile &&)>  op = TA::detail::Noop<Tile, ComputeTile, true>(),
std::function< shellpair_data_accessor_t(const Basis *, const Basis *)>  make_shellpair_data_accessor = {},
std::shared_ptr< const math::PetiteList >  plist = math::PetiteList::make_trivial() 
)

Construct direct integral tensors in parallel with screening.

Same requirements on Op as those in Integral Builder

◆ direct_sparse_integrals() [2/2]

template<typename Tile = TA::TensorD, typename ComputeTile = Tile, typename Engine , typename Idx , typename Bases >
DirectArray<Tile, TA::SparsePolicy> mpqc::lcao::gaussian::direct_sparse_integrals ( madness::World &  world,
ShrPool< Engine >  shr_pool,
Bases &&  bases,
std::vector< std::pair< Idx, float >> const &  user_provided_norms,
bool  user_provided_norms_perelem,
std::shared_ptr< Screener screen = std::make_shared<Screener>(Screener{}),
std::function< Tile(ComputeTile &&)>  op = TA::detail::Noop<Tile, ComputeTile, true>(),
std::function< shellpair_data_accessor_t(const Basis *, const Basis *)>  make_shellpair_data_accessor = {},
std::shared_ptr< const math::PetiteList >  plist = math::PetiteList::make_trivial() 
)

Construct direct integral tensors in parallel with screening.

Parameters
user_provided_normsis a user supplied replicated tensor with the norm estimates for the array. It is of type vector<std::pair<Idx, float>> where Idx is a type that provides random access via operator[] and the values are the indices of the tile. The float is then the norm estimate for that tile.
user_provided_norms_perelemwhether user_provided_norms are already scaled (per-element) or not

◆ fock_from_soad() [1/2]

template<typename Tile , typename ComputeTile , typename Policy >
TA::DistArrayVector< Tile, typename std::enable_if<std::is_same<Policy, TA::SparsePolicy>::value, TA::SparsePolicy>::type> mpqc::lcao::gaussian::fock_from_soad ( madness::World &  world,
Molecule const &  clustered_mol,
Basis const &  obs,
Basis const &  minbs,
TA::DistArray< Tile, Policy > const &  H,
int const  n_components,
std::function< shellpair_data_accessor_t(const Basis *bs0, const Basis *bs1)>  make_shellpair_data_accessor = {},
std::function< Tile(ComputeTile &&)>  op = TA::detail::Noop<Tile, ComputeTile, true>() 
)

fock matrix computed from soad for SparsePolicy

Template Parameters
Tilea Tile type for the result
ComputeTilethe storage tile for integrals
PolicyTA::SparsePolicy
Parameters
worldworld object
clustered_molmolecule class
obsbasis object
minbsminimal basis object
engsengine
HH matrix
n_componentsnumber of components (size) of returned DistArrayVector
opoperator to convert ComputeTile to Tile
Returns
Fock matrix from soad

◆ fock_from_soad() [2/2]

template<typename Tile , typename ComputeTile , typename Policy >
TA::DistArrayVector< Tile, typename std::enable_if<std::is_same<Policy, TA::DensePolicy>::value, TA::DensePolicy>::type> mpqc::lcao::gaussian::fock_from_soad ( madness::World &  world,
Molecule const &  clustered_mol,
Basis const &  obs,
Basis const &  minbs,
TA::DistArray< Tile, Policy > const &  H,
int const  n_components,
std::function< shellpair_data_accessor_t(const Basis *bs0, const Basis *bs1)>  make_shellpair_data_accessor = {},
std::function< Tile(ComputeTile &&)>  op = TA::detail::Noop<Tile, ComputeTile, true>() 
)

fock matrix computed from soad for DensePolicy

Template Parameters
ShrPool
Tilea Tile type for the result
ComputeTilethe storage tile for integrals
PolicyTA::DensePolicy
Parameters
worldworld object
clustered_molmolecule class
obsbasis object
minbsminimal basis object
engsengine
HH matrix
n_componentsnumber of components (size) of returned DistArrayVector
opoperator to convert ComputeTile to Tile
Returns
Fock matrix from soad

◆ index_to_basis()

std::shared_ptr< const Basis > mpqc::lcao::gaussian::index_to_basis ( const OrbitalBasisRegistry basis_registry,
const OrbitalIndex index 
)

given OrbitalIndex, find the corresponding basis

◆ make_direct_integral_builder()

template<typename Tile , typename ComputeTile = Tile, typename Engine = libint2::Engine>
std::shared_ptr<DirectIntegralBuilder<Tile, ComputeTile, Engine> > mpqc::lcao::gaussian::make_direct_integral_builder ( madness::World &  world,
ShrPool< Engine >  shr_epool,
const BasisShrVector bases,
std::shared_ptr< Screener shr_screen,
std::function< Tile(ComputeTile &&)>  op,
std::function< shellpair_data_accessor_t(const Basis *, const Basis *)>  make_shellpair_data_accessor = {},
std::shared_ptr< const math::PetiteList >  plist = math::PetiteList::make_trivial() 
)

Function to make detection of template parameters easier, see DirectIntegralBuilder for details.

◆ make_engine()

template<typename Bases >
libint2::Engine mpqc::lcao::gaussian::make_engine ( const libint2::Operator &  oper,
Bases &&  bases,
libint2::BraKet  braket = libint2::BraKet::invalid,
libint2::any  oper_params = libint2::any(),
double  engine_precision = std::numeric_limits<double>::epsilon(),
const libint2::scalar_type &  scalar = libint2::scalar_type(1.0) 
)

makes an engine for computing integrals of operator oper over bases bases

◆ make_engine_pool()

template<typename Bases >
ShrPool<libint2::Engine> mpqc::lcao::gaussian::make_engine_pool ( const libint2::Operator &  oper,
Bases &&  bases,
libint2::BraKet  braket = libint2::BraKet::invalid,
libint2::any  oper_params = libint2::any(),
const libint2::scalar_type &  scalar = libint2::scalar_type(1.0),
double  engine_precision = std::numeric_limits<double>::epsilon() 
)
inline

Make an pool of Engine objects.

Template Parameters
Basesa contiguous range of Basis objects or std::reference_wrapper's to (const) Basis objects
Parameters
operthe operator type
basesan object of Bases type
braketthe braket type; the default is to use the value returned by libint2::default_braket(oper)
oper_paramsthe operator parameters; the default is to use libint2::default_params(oper)
scalarintegrals will be scaled by this factor; the default value is 1
engine_precisionthe (absolute) precision target for the Engine
Returns
a std::shared_ptr to the Engine pool

◆ make_integral_builder()

template<typename Tile , typename ComputeTile = Tile, typename Engine = libint2::Engine>
std::shared_ptr<IntegralBuilder<Tile, ComputeTile, Engine> > mpqc::lcao::gaussian::make_integral_builder ( ShrPool< Engine >  shr_epool,
const BasisShrVector bases,
std::shared_ptr< Screener shr_screen,
std::function< Tile(ComputeTile &&)>  op,
std::function< shellpair_data_accessor_t(const Basis *, const Basis *)>  make_shellpair_data_accessor = {},
std::shared_ptr< const math::PetiteList >  plist = math::PetiteList::make_trivial() 
)

Function to make detection of template parameters easier, see IntegralBuilder for details.

◆ make_q()

q_vector mpqc::lcao::gaussian::make_q ( Molecule const &  mol)
inline

◆ make_schwarz_screener()

std::unique_ptr< SchwarzScreen > mpqc::lcao::gaussian::make_schwarz_screener ( std::tuple< const Basis & >  bra,
std::tuple< const Basis &, const Basis & >  ket,
double  threshold 
)

constructs an l2-based SchwarzScreen screener for (bra[0]|ket[0] ket[1]) integrals

Parameters
[in]bratuple of references to bra bases
[in]kettuple of references to ket bases
[in]thresholdscreening threshold; see SchwarzScreen::SchwarzScreen()

◆ make_screener()

template<typename Bases >
std::shared_ptr<ThresholdedScreener> mpqc::lcao::gaussian::make_screener ( madness::World &  world,
const ShrPool< libint2::Engine > &  engine,
Bases &&  bases,
const std::string &  screen,
double  screen_threshold 
)

factory for ThresholdedScreener objects

Template Parameters
Basesa directly-addressable sequence of basis sets (or reference proxies, e.g. std::reference_wrapper)
Parameters
worldis a reference to the world in which the distributed screener data will be computed and stores
engis a reference to a ShrPool<E>
bs_arrayis a reference to a directly-addressable sequence of basis sets, if the length is 3 then DF integrals are assumed and the first basis is assumed to be the auxiliary basis if the length is 4 then it is assumed that four center screening is desired. There is no requirement that any basis sets be the same.
screenstring type of the screener requested; valid values are schwarz (Schwarz screener using the L2 norm), schwarz_inf (Schwarz screener using the infinity norm), qqr (QQR screener using the L2 norm; 4-center only), qqr_inf (QQR screener using the infinity norm; 4-center only) sqvl (SQVl screener using the L2 norm; 3-center only), sqvl_inf (SQVl screener using the infinity norm; 3-center only), default (QQR screener for 4-center ints, SQVl screener for 3-center ints, both L2-based)
screen_thresholdis the screening threshold
Returns
the screener object

◆ max_am()

int64_t mpqc::lcao::gaussian::max_am ( ShellVec const &  shell_vec)

Returns the maximum angular momement of any shell in the vector.

◆ max_nprim()

int64_t mpqc::lcao::gaussian::max_nprim ( ShellVec const &  shell_vec)

Returns the maximum number of primatives of any shell in the vector.

◆ merge()

std::shared_ptr< const Basis > mpqc::lcao::gaussian::merge ( const Basis basis1,
const Basis basis2 
)

merges two Basis objects by concatenating their shell cluster sequences.

Returns
a Basis object in which shells of basis1 and basis2

◆ MPQC_EXTERN_TEMPLATE()

mpqc::lcao::gaussian::MPQC_EXTERN_TEMPLATE ( class AOFactory<>  )

◆ nfunctions()

int64_t mpqc::lcao::gaussian::nfunctions ( ShellVec const &  shell_vec)

Returns the maximum number of primatives of any shell in the vector.

◆ operator<<() [1/3]

std::ostream & mpqc::lcao::gaussian::operator<< ( std::ostream &  os,
AtomicBasis const &  b 
)

construct a map that maps column of basis to column of sub_basis, sub_basis has to be a subset of basis

Warning
the value in index starts with 1, value 0 in index indicates this column is missing in sub_basis This approach uses N^2 algorithm //TODO need unit test for this
Parameters
basisBasis object
sub_basisBasis object, subset of basis
Returns
a vector of column id of sub_basis

◆ operator<<() [2/3]

std::ostream & mpqc::lcao::gaussian::operator<< ( std::ostream &  os,
Basis const &  b 
)

construct a map that maps column of basis to column of sub_basis, sub_basis has to be a subset of basis

Warning
the value in index starts with 1, value 0 in index indicates this column is missing in sub_basis This approach uses N^2 algorithm //TODO need unit test for this
Parameters
basisBasis object
sub_basisBasis object, subset of basis
Returns
a vector of column id of sub_basis

◆ operator<<() [3/3]

std::ostream & mpqc::lcao::gaussian::operator<< ( std::ostream &  os,
Basis::Factory const &  f 
)

construct a map that maps column of basis to column of sub_basis, sub_basis has to be a subset of basis

Warning
the value in index starts with 1, value 0 in index indicates this column is missing in sub_basis This approach uses N^2 algorithm //TODO need unit test for this
Parameters
basisBasis object
sub_basisBasis object, subset of basis
Returns
a vector of column id of sub_basis

◆ operator==() [1/2]

bool mpqc::lcao::gaussian::operator== ( const Basis one,
const Basis two 
)
inline

merges two Basis objects by concatenating their shell cluster sequences.

Returns
a Basis object in which shells of basis1 and basis2

◆ operator==() [2/2]

bool mpqc::lcao::gaussian::operator== ( const Q2Matrix first,
const Q2Matrix second 
)
inline

Compares two Q2Matrix objects. This does not handle transposes, i.e. it returns false if first is a transpose of second .

Parameters
[in]firsta Q2Matrix object
[in]seconda Q2Matrix object
Returns
true if first refers to the same bases and second and uses same operator type and norm op

◆ parallel_compute_shellpair_list()

std::vector< std::vector< size_t > > mpqc::lcao::gaussian::parallel_compute_shellpair_list ( madness::World &  world,
const Basis basis1,
const Basis basis2,
double  threshold = 1e-12,
double  engine_precision = 0.0 
)

This computes non-negligible shell pair list; shells i and j form a non-negligible pair if they share a center or the Frobenius norm of their overlap is greater than threshold.

Parameters
basis1a basis
basis2a basis
threshold
Returns
a list of pairs with key: shell index mapped value: a vector of shell indices

◆ parallel_make_basis()

Basis mpqc::lcao::gaussian::parallel_make_basis ( madness::World &  world,
const Basis::Factory factory,
const mpqc::Molecule mol 
)

construct Basis from a factory and a Molecule on process 0 and broadcast to the entire world

Parameters
worldthe madness::World
factorythe Basis::Factory object
molthe Molecule object
Returns
the Basis object
Exceptions
Uncomputableif Libint does not include this basis or it is not supported for one or more atoms.

◆ reblock()

template<typename Op , typename... Args>
Basis mpqc::lcao::gaussian::reblock ( Basis const &  basis,
Op  op,
Args...  args 
)

reblock allows for reblocking a basis

Warning
If reblocking a basis with the intent to use it with tensors computed with the old basis you must be careful not to reorder the shells.
Parameters
opshould be a function that takes a std::vector<Shell> and returns a std::vector<std::vector<Shell>> for use in initializing a Basis.

◆ reblock_basis()

std::vector<std::vector<libint2::Shell> > mpqc::lcao::gaussian::reblock_basis ( std::vector< libint2::Shell >  shells,
std::size_t  blocksize 
)

◆ soad_task()

template<typename Engs , typename Tile , typename ComputeTile , typename Policy >
void mpqc::lcao::gaussian::soad_task ( Engs  eng_pool,
const std::array< int64_t, 2 > &  tile_idx,
const std::array< std::vector< libint2::Shell > const *, 2 > &  obs_row_and_col_ptrs,
std::vector< libint2::Shell > const *  min_bs,
const RowMatrixXd D,
TA::DistArray< Tile, Policy > *  F,
std::function< Tile(ComputeTile &&)>  op,
const std::array< shellpair_data_accessor_t, 3 > *  shellpair_data_accessors,
const std::array< std::shared_ptr< ThresholdedScreener >, 2 > *  jk_screeners 
)

◆ sparse_integrals()

template<typename E , typename Tile = TA::TensorD, typename Bases >
TA::DistArray<Tile, TA::SparsePolicy> mpqc::lcao::gaussian::sparse_integrals ( madness::World &  world,
ShrPool< E >  shr_pool,
Bases &&  bases,
std::shared_ptr< Screener screen = std::make_shared<Screener>(Screener{}),
std::function< Tile(Tile &&)>  op = TA::detail::Noop<Tile, Tile, true>{},
std::function< shellpair_data_accessor_t(const Basis *, const Basis *)>  make_shellpair_data_accessor = {},
std::shared_ptr< const math::PetiteList >  plist = math::PetiteList::make_trivial() 
)

Construct sparse integral tensors in parallel.

Parameters
shr_poolshould be a std::shared_ptr to an IntegralTSPool
basesshould be a std::array of Basis, which will be copied.
opneeds to be a function or functor that takes Tile&& and returns an type. Op is copied so it can be moved.
auto t = [](TA::TensorD &&ten){return std::move(ten);};
screena std::shared_ptr to a Screener.
make_shellpair_data_accessora functor that makes shellpair data accessors; such accessors are to accept basis function indices, rather than shell indices

◆ sparse_xyz_integrals()

template<typename E , typename Tile = TA::TensorD, typename Bases >
TA::DistArrayVector<Tile, TA::SparsePolicy> mpqc::lcao::gaussian::sparse_xyz_integrals ( madness::World &  world,
ShrPool< E >  shr_pool,
Bases &&  bases,
std::function< Tile(Tile &&)>  op = TA::detail::Noop<Tile, Tile, true>{} 
)

Construct set of sparse integral tensors in parallel.

This is needed for integrals such as the dipole integrals that come as a set.

Template Parameters
Tilea (contiguous) tensor type to store the result
Parameters
shr_poolshould be a std::shared_ptr to an IntegralTSPool
basesshould be a std::array of Basis, which will be copied.
opneeds to be a function or functor that takes a TA::TensorD && and returns any valid tile type. Op is copied so it can be moved.
auto t = [](TA::TensorD &&ten){return std::move(ten);};

◆ sub_basis_map()

Eigen::RowVectorXi mpqc::lcao::gaussian::sub_basis_map ( const Basis basis,
const Basis sub_basis 
)

construct a map that maps column of basis to column of sub_basis, sub_basis has to be a subset of basis

Warning
the value in index starts with 1, value 0 in index indicates this column is missing in sub_basis This approach uses N^2 algorithm //TODO need unit test for this
Parameters
basisBasis object
sub_basisBasis object, subset of basis
Returns
a vector of column id of sub_basis

◆ tensorZ_to_tensorD()

template<typename Policy , bool is_real = true>
TA::DistArray<TA::TensorD, Policy> mpqc::lcao::gaussian::tensorZ_to_tensorD ( const TA::DistArray< TA::TensorZ, Policy > &  complex_array)

This takes real or imaginary part from a complex array.

Template Parameters
Policycan be TA::SparsePolicy or TA::DensePolicy
Is_realtrue if user requests real part, false if imaginary
Parameters
aTensorZ array
Returns
a TensorD array

◆ to_ao_factory() [1/2]

template<typename Tile , typename Policy >
std::shared_ptr<AOFactory<Tile, Policy> > mpqc::lcao::gaussian::to_ao_factory ( const std::shared_ptr< lcao::AOFactory< Tile, Policy >> &  factory)

◆ to_ao_factory() [2/2]

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

◆ to_libint2_operator()

libint2::Operator mpqc::lcao::gaussian::to_libint2_operator ( Operator::Type  mpqc_oper)

◆ to_libint2_operator_params()

libint2::any mpqc::lcao::gaussian::to_libint2_operator_params ( Operator::Type  mpqc_oper,
const Molecule molecule,
const std::map< Operator::Type, libint2::any > *  oper_params 
)

◆ to_libint2_scale_factor()

libint2::scalar_type mpqc::lcao::gaussian::to_libint2_scale_factor ( Operator::Type  mpqc_oper,
const std::map< Operator::Type, libint2::any > *  oper_params 
)

◆ to_molden()

template<typename OrbSpace >
void mpqc::lcao::gaussian::to_molden ( std::string_view  fname_prefix,
const OrbSpace &  orbspace,
const WavefunctionWorld wfn_world 
)

prints an orbital space to a Molden file

Template Parameters
OrbSpace
Parameters
orbspace
wfn_world

◆ untruncated_direct_sparse_integrals()

template<typename Tile = TA::TensorD, typename ComputeTile = Tile, typename Engine , typename Bases >
DirectArray<Tile, TA::SparsePolicy, Engine> mpqc::lcao::gaussian::untruncated_direct_sparse_integrals ( madness::World &  world,
ShrPool< Engine >  shr_pool,
Bases &&  bases,
std::shared_ptr< Screener screen = std::make_shared<Screener>(Screener{}),
std::function< Tile(ComputeTile &&)>  op = TA::detail::Noop<Tile, ComputeTile, true>(),
std::function< shellpair_data_accessor_t(const Basis *, const Basis *)>  make_shellpair_data_accessor = {},
std::shared_ptr< const math::PetiteList >  plist = math::PetiteList::make_trivial() 
)

Construct direct integral tensors in parallel with screening.

Same requirements on Op as those in Integral Builder.

I only plan to use this for CADF, no point in truncating tiles.

◆ xyz_integrals() [1/2]

template<typename Tile , typename Policy , typename E , typename Bases >
TA::DistArrayVector< Tile, std::enable_if_t<std::is_same<Policy, TA::DensePolicy>::value, TA::DensePolicy> > mpqc::lcao::gaussian::xyz_integrals ( madness::World &  world,
ShrPool< E >  shr_pool,
Bases &&  bases,
std::function< Tile(Tile &&)>  op = TA::detail::Noop<Tile, Tile, true>{} 
)

◆ xyz_integrals() [2/2]

template<typename Tile , typename Policy , typename E , typename Bases >
TA::DistArrayVector< Tile, std::enable_if_t<std::is_same<Policy, TA::SparsePolicy>::value, TA::SparsePolicy> > mpqc::lcao::gaussian::xyz_integrals ( madness::World &  world,
ShrPool< E >  shr_pool,
Bases &&  bases,
std::function< Tile(Tile &&)>  op = TA::detail::Noop<Tile, Tile, true>{} 
)