MPQC  3.0.0-alpha
sc::DensityFitting Class Reference

Decomposition by density fitting with respect to some kernel. More...

#include <chemistry/qc/lcao/df.h>

Inheritance diagram for sc::DensityFitting:
sc::SavableState sc::DescribedClass sc::RefCount sc::PermutedDensityFitting sc::TransformedDensityFitting

Public Types

enum  SolveMethod {
  SolveMethod_InverseCholesky = 0, SolveMethod_Cholesky = 1, SolveMethod_RefinedCholesky = 2, SolveMethod_InverseBunchKaufman = 3,
  SolveMethod_BunchKaufman = 4, SolveMethod_RefinedBunchKaufman = 5
}
 
typedef TwoBodyMOIntsRuntimeUnion23 MOIntsRuntime
 

Public Member Functions

 DensityFitting (const Ref< MOIntsRuntime > &rtime, const std::string &kernel_key, SolveMethod solver, const Ref< OrbitalSpace > &space1, const Ref< OrbitalSpace > &space2, const Ref< GaussianBasisSet > &fitting_basis, bool ri=false)
 Determines density fitting of product (space1 space2) in terms of fitting basis. More...
 
 DensityFitting (StateIn &)
 
void save_data_state (StateOut &)
 Save the base classes (with save_data_state) and the members in the same order that the StateIn CTOR initializes them. More...
 
const Ref< MOIntsRuntime > & runtime () const
 
const Ref< Integral > & integral () const
 
const Ref< OrbitalSpace > & space1 () const
 
const Ref< OrbitalSpace > & space2 () const
 
const Ref< GaussianBasisSet > & fbasis () const
 
const std::string & kernel_key () const
 
SolveMethod solver () const
 
const Ref< DistArray4 > & C () const
 returns the fitting coeffcients
 
RefSCDimension product_dimension () const
 produces RefSCDimension that corresponds to the product space
 
virtual void compute ()
 
- Public Member Functions inherited from sc::SavableState
SavableStateoperator= (const SavableState &)
 
void save_state (StateOut &)
 Save the state of the object as specified by the StateOut object. More...
 
void save_object_state (StateOut &)
 This can be used for saving state when the exact type of the object is known for both the save and the restore. More...
 
virtual void save_vbase_state (StateOut &)
 Save the virtual bases for the object. More...
 
- Public Member Functions inherited from sc::DescribedClass
 DescribedClass (const DescribedClass &)
 
DescribedClassoperator= (const DescribedClass &)
 
ClassDescclass_desc () const MPQC__NOEXCEPT
 This returns the unique pointer to the ClassDesc corresponding to the given type_info object. More...
 
const char * class_name () const
 Return the name of the object's exact type.
 
int class_version () const
 Return the version of the class.
 
virtual void print (std::ostream &=ExEnv::out0()) const
 Print the object.
 
Ref< DescribedClassref ()
 Return this object wrapped up in a Ref smart pointer. More...
 
- Public Member Functions inherited from sc::RefCount
size_t identifier () const
 Return the unique identifier for this object that can be compared for different objects of different types. More...
 
int lock_ptr () const
 Lock this object.
 
int unlock_ptr () const
 Unlock this object.
 
void use_locks (bool inVal)
 start and stop using locks on this object
 
refcount_t nreference () const
 Return the reference count.
 
refcount_t reference ()
 Increment the reference count and return the new count.
 
refcount_t dereference ()
 Decrement the reference count and return the new count.
 
int managed () const
 
void unmanage ()
 Turn off the reference counting mechanism for this object. More...
 

Static Public Member Functions

static int definite_kernel (TwoBodyOper::type kernel_type, const Ref< IntParams > &kernel_params)
 Test whether a density-fitting kernel has definite sign. More...
 
static int definite_kernel (TwoBodyOper::type kernel_type, const std::string &kernel_params_key)
 
- Static Public Member Functions inherited from sc::SavableState
static void save_state (SavableState *s, StateOut &)
 
static SavableStaterestore_state (StateIn &si)
 Restores objects saved with save_state. More...
 
static SavableStatekey_restore_state (StateIn &si, const char *keyword)
 Like restore_state, but keyword is used to override values while restoring.
 
static SavableStatedir_restore_state (StateIn &si, const char *objectname, const char *keyword=0)
 

Protected Member Functions

const RefSymmSCMatrixkernel () const
 returns the kernel matrix in the fitting basis, i.e. $ (A|\hat{W}|B) $ .
 
const Ref< DistArray4 > & conjugateC () const
 returns the conjugate expansion matrix defined as $ \tilde{C}^{rs}_A \equiv (A|\hat{W}|rs) = \sum_B C^{rs}_B (A|\hat{W}|B) $ . More...
 
- Protected Member Functions inherited from sc::SavableState
 SavableState (const SavableState &)
 
 SavableState (StateIn &)
 Each derived class StateIn CTOR handles the restore corresponding to calling save_object_state, save_vbase_state, and save_data_state listed above. More...
 
- Protected Member Functions inherited from sc::RefCount
 RefCount (const RefCount &)
 
RefCountoperator= (const RefCount &)
 

Protected Attributes

Ref< DistArray4C_
 
RefSymmSCMatrix kernel_
 
Ref< DistArray4cC_
 

Detailed Description

Decomposition by density fitting with respect to some kernel.

The definition is as follows:

\[ |rs) = \sum_A C^{rs}_A |A) \quad , \]

where $ |A) $ belong to a fitting basis and $ C^{rs}_A $ are determined by solving the linear system

\[ (A|\hat{W}|rs) = \sum_B C^{rs}_B (A|\hat{W}|B) \quad . \]

The kernel $\hat{W} $ is any (positive-definite) operator.

Parameters
runtimeis the MOIntsRuntime object used to compute MO/AO integrals
kernel_keydenotes the kernel operator to be used for fitting. Currently only "1/r_{12}" is supported.
space1is space r in |rs) product space
space2is space s in |rs) product space
fitting_basisis the basis set in which to perform the fitting

Constructor & Destructor Documentation

◆ DensityFitting()

sc::DensityFitting::DensityFitting ( const Ref< MOIntsRuntime > &  rtime,
const std::string &  kernel_key,
SolveMethod  solver,
const Ref< OrbitalSpace > &  space1,
const Ref< OrbitalSpace > &  space2,
const Ref< GaussianBasisSet > &  fitting_basis,
bool  ri = false 
)

Determines density fitting of product (space1 space2) in terms of fitting basis.

Parameters
rtime
kernel_key
solver
space1
space2
fitting_basis
rifalse -> compute DF coeficients, true -> compute DF coefficients multiplied by the inverse square root of the fitting matrix (or of the matrix scaled by -1, for negative-definite kernels); only supported if using inverse-based Cholesky solver.
Returns

Member Function Documentation

◆ conjugateC()

const Ref<DistArray4>& sc::DensityFitting::conjugateC ( ) const
inlineprotected

returns the conjugate expansion matrix defined as $ \tilde{C}^{rs}_A \equiv (A|\hat{W}|rs) = \sum_B C^{rs}_B (A|\hat{W}|B) $ .

$ \tilde{C} $ is conjugate to $ C $ in the sense that their product returns reconstructed integrals: $ (pq|rs) \approx \sum_A \tilde{C}^A_{pq} C^{rs}_A $

◆ definite_kernel()

static int sc::DensityFitting::definite_kernel ( TwoBodyOper::type  kernel_type,
const Ref< IntParams > &  kernel_params 
)
static

Test whether a density-fitting kernel has definite sign.

Note
this is a preliminary implementation that makes common assumptions about the kernel_params argument
Parameters
kernel_typekernel type; only some operators are accepted as density fitting kernel
kernel_paramskernel parameters; must be in ParamsRegistry
Returns
+1 if the kernel is positive-definite, -1 if the kernel is negative-definite, 0 otherwise

◆ save_data_state()

void sc::DensityFitting::save_data_state ( StateOut )
virtual

Save the base classes (with save_data_state) and the members in the same order that the StateIn CTOR initializes them.

This must be implemented by the derived class if the class has data.

Reimplemented from sc::SavableState.

Reimplemented in sc::PermutedDensityFitting, and sc::TransformedDensityFitting.


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

Generated at Sun Jan 26 2020 23:24:04 for MPQC 3.0.0-alpha using the documentation package Doxygen 1.8.16.