MADNESS  0.10.1
Classes | Public Member Functions | Static Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | Private Types | Private Attributes | Friends | List of all members
madness::Nemo Class Reference

The Nemo class. More...

#include <nemo.h>

Inheritance diagram for madness::Nemo:
Inheritance graph
[legend]
Collaboration diagram for madness::Nemo:
Collaboration graph
[legend]

Classes

struct  NemoCalculationParameters
 class holding parameters for a nemo calculation beyond the standard dft parameters from moldft More...
 

Public Member Functions

 Nemo (World &world, const commandlineparser &parser)
 ctor More...
 
std::vector< vecfuncTcompute_all_cphf ()
 solve the CPHF equation for all displacements More...
 
vecfuncT compute_cphf_parallel_term (const size_t iatom, const int iaxis) const
 
tensorT compute_fock_matrix (const vecfuncT &nemo, const tensorT &occ) const
 compute the Fock matrix from scratch More...
 
Tensor< double > compute_IR_intensities (const Tensor< double > &normalmodes, const vecfuncT &dens_pt) const
 compute the IR intensities in the double harmonic approximation More...
 
bool do_ac () const
 
bool do_pcm () const
 
bool do_symmetry () const
 
AC< 3 > get_ac () const
 
std::shared_ptr< SCFget_calc () const
 
const NemoCalculationParametersget_param () const
 
PCM get_pcm () const
 
projector_irrep get_symmetry_projector () const
 return the symmetry_projector More...
 
Tensor< double > gradient (const Tensor< double > &x)
 compute the nuclear gradients More...
 
Tensor< double > hessian (const Tensor< double > &x)
 returns the molecular hessian matrix at structure x More...
 
bool is_dft () const
 
real_function_3d kinetic_energy_potential (const vecfuncT &nemo) const
 
real_function_3d make_ddensity (const real_function_3d &rhonemo, const int axis) const
 make the derivative of the density More...
 
real_function_3d make_density (const Tensor< double > &occ, const vecfuncT &nemo) const
 make the density (alpha or beta) More...
 
real_function_3d make_density (const tensorT &occ, const vecfuncT &bra, const vecfuncT &ket, const bool refine=false) const
 make the density using different bra and ket vectors More...
 
virtual std::shared_ptr< Fock< double, 3 > > make_fock_operator () const
 construct the fock operator based on the calculation parameters (K or XC?) More...
 
real_function_3d make_laplacian_density (const real_function_3d &rhonemo) const
 the Laplacian of the density More...
 
real_function_3d make_sigma (const real_function_3d &rho1, const real_function_3d &rho2) const
 compute the reduced densities sigma (gamma) for GGA functionals More...
 
Moleculemolecule ()
 return a reference to the molecule More...
 
Moleculemolecule () const
 return a reference to the molecule More...
 
std::string name () const
 
bool provides_gradient () const
 Override this to return true if the derivative is implemented. More...
 
Tensor< double > purify_hessian (const Tensor< double > &hessian) const
 purify and symmetrize the hessian More...
 
bool selftest ()
 
vecfuncT solve_cphf (const size_t iatom, const int iaxis, const Tensor< double > fock, const vecfuncT &guess, const vecfuncT &rhsconst, const Tensor< double > incomplete_hessian, const vecfuncT &parallel, const SCFProtocol &p, const std::string &xc_data) const
 solve the CPHF equations for the nuclear displacements More...
 
virtual double value ()
 
virtual double value (const Tensor< double > &x)
 Should return the value of the objective function. More...
 
- Public Member Functions inherited from madness::NemoBase
 NemoBase (World &w)
 
virtual ~NemoBase ()
 
bool check_convergence (const std::vector< double > energies, const std::vector< double > oldenergies, const double bsh_norm, const double delta_density, const CalculationParameters &param, const double econv, const double dconv) const
 
template<typename T , std::size_t NDIM>
Function< typename Tensor< T >::scalar_type, NDIMcompute_density (const std::vector< Function< T, NDIM > > nemo) const
 
Tensor< double > compute_gradient (const real_function_3d &rhonemo, const Molecule &molecule) const
 compute the nuclear gradients More...
 
template<typename T , std::size_t NDIM>
double compute_kinetic_energy (const std::vector< Function< T, NDIM > > &nemo) const
 compute kinetic energy as square of the "analytical" expectation value More...
 
template<typename T , std::size_t NDIM>
double compute_kinetic_energy1 (const std::vector< Function< T, NDIM > > &nemo) const
 compute kinetic energy as square of the "analytical" derivative of the orbitals More...
 
template<typename T , std::size_t NDIM>
double compute_kinetic_energy1a (const std::vector< Function< T, NDIM > > &nemo) const
 compute kinetic energy as square of the "analytical" derivative of the orbitals More...
 
template<typename T , std::size_t NDIM>
double compute_kinetic_energy2 (const std::vector< Function< T, NDIM > > &nemo) const
 compute kinetic energy as direct derivative of the orbitals (probably imprecise) More...
 
void construct_nuclear_correlation_factor (const Molecule &molecule, const std::shared_ptr< PotentialManager > pm, const std::pair< std::string, double > ncf_parameter)
 
std::shared_ptr< NuclearCorrelationFactorget_ncf_ptr () const
 create an instance of the derived object based on the input parameters More...
 
virtual void invalidate_factors_and_potentials ()
 
virtual bool need_recompute_factors_and_potentials (const double thresh) const
 
template<typename T , std::size_t NDIM>
void orthonormalize (std::vector< Function< T, NDIM > > &nemo, const Function< double, NDIM > metric=Function< double, NDIM >(), const double trantol=FunctionDefaults< NDIM >::get_thresh() *0.01) const
 orthonormalize the vectors More...
 
- Public Member Functions inherited from madness::OptimizationTargetInterface
virtual ~OptimizationTargetInterface ()
 
double test_gradient (Tensor< double > &x, double value_precision, bool doprint=true)
 Numerical test of the derivative ... optionally prints to stdout, returns max abs error. More...
 
virtual void value_and_gradient (const Tensor< double > &x, double &value, Tensor< double > &gradient)
 Reimplement if more efficient to evaluate both value and gradient in one call. More...
 

Static Public Member Functions

static void help ()
 
static void print_parameters ()
 
static void smoothen (real_function_3d &f)
 smooth a function by projecting it onto k-1 and then average with k More...
 
- Static Public Member Functions inherited from madness::NemoBase
template<typename T , std::size_t NDIM>
static void normalize (std::vector< Function< T, NDIM > > &nemo, const Function< double, NDIM > metric=Function< double, NDIM >())
 normalize the nemos More...
 
template<typename T >
static Tensor< TQ2 (const Tensor< T > &s)
 

Public Attributes

NemoCalculationParameters param
 
- Public Attributes inherited from madness::NemoBase
std::shared_ptr< NuclearCorrelationFactorncf
 the nuclear correlation factor More...
 
real_function_3d R
 the nuclear correlation factor More...
 
real_function_3d R_square
 the square of the nuclear correlation factor More...
 
Worldworld
 

Protected Member Functions

std::vector< double > compute_energy_regularized (const vecfuncT &nemo, const vecfuncT &Jnemo, const vecfuncT &Knemo, const vecfuncT &Unemo) const
 given nemos, compute the HF energy using the regularized expressions for T and V More...
 
void compute_nemo_potentials (const vecfuncT &nemo, vecfuncT &Jnemo, vecfuncT &Knemo, vecfuncT &xcnemo, vecfuncT &pcmnemo, vecfuncT &Unemo) const
 compute the reconstructed orbitals, and all potentials applied on nemo More...
 
real_function_3d get_coulomb_potential (const vecfuncT &psi) const
 return the Coulomb potential More...
 
template<typename T , size_t NDIM>
void load_function (std::vector< Function< T, NDIM > > &f, const std::string name) const
 load a function More...
 
vecfuncT localize (const vecfuncT &nemo, const double dconv, const bool randomize) const
 localize the nemo orbitals More...
 
vecfuncT make_cphf_constant_term (const size_t iatom, const int iaxis, const vecfuncT &R2nemo, const real_function_3d &rhonemo) const
 compute the constant term for the CPHF equations More...
 
Tensor< double > make_incomplete_hessian () const
 compute the incomplete hessian More...
 
Tensor< double > make_incomplete_hessian_response_part (const std::vector< vecfuncT > &xi) const
 compute the complementary incomplete hessian More...
 
void make_plots (const real_function_3d &f, const std::string &name="function") const
 
template<typename solverT >
void rotate_subspace (World &world, const tensorT &U, solverT &solver, int lo, int nfunc) const
 rotate the KAIN subspace (cf. SCF.cc) More...
 
template<typename T , size_t NDIM>
void save_function (const std::vector< Function< T, NDIM > > &f, const std::string name) const
 save a function More...
 
void set_protocol (const double thresh)
 adapt the thresholds consistently to a common value More...
 
double solve (const SCFProtocol &proto)
 solve the HF equations More...
 
double trantol () const
 return the threshold for vanishing elements in orbital rotations More...
 

Protected Attributes

AC< 3 > ac
 asymptotic correction for DFT More...
 
std::shared_ptr< SCFcalc
 
std::shared_ptr< real_convolution_3dpoisson
 a poisson solver More...
 
projector_irrep symmetry_projector
 

Private Types

typedef std::shared_ptr< real_convolution_3dpoperatorT
 

Private Attributes

double coords_sum
 sum of square of coords at last solved geometry More...
 
PCM pcm
 polarizable continuum model More...
 

Friends

class PNO
 
class TDHF
 

Detailed Description

The Nemo class.

Member Typedef Documentation

◆ poperatorT

typedef std::shared_ptr<real_convolution_3d> madness::Nemo::poperatorT
private

Constructor & Destructor Documentation

◆ Nemo()

madness::Nemo::Nemo ( World world,
const commandlineparser parser 
)

Member Function Documentation

◆ compute_all_cphf()

std::vector< vecfuncT > madness::Nemo::compute_all_cphf ( )

◆ compute_cphf_parallel_term()

vecfuncT madness::Nemo::compute_cphf_parallel_term ( const size_t  iatom,
const int  iaxis 
) const

this function computes that part of the orbital response that is parallel to the occupied space.

\[ F^X = F^\perp + F^\parallel \]

If no NCF's are used F^\parallel vanishes. If NCF's are used this term does not vanish because the derivatives of the NCF does not vanish, and it is given by

\[ F_i^\parallel = -\frac{1}{2}\sum_k|F_k ><F_k | (R^2)^X | F_i> \]

References calc, madness::matrix_inner(), madness::mul(), madness::transform(), madness::truncate(), and madness::NemoBase::world.

Referenced by compute_all_cphf().

◆ compute_energy_regularized()

std::vector< double > madness::Nemo::compute_energy_regularized ( const vecfuncT nemo,
const vecfuncT Jnemo,
const vecfuncT Knemo,
const vecfuncT Unemo 
) const
protected

◆ compute_fock_matrix()

tensorT madness::Nemo::compute_fock_matrix ( const vecfuncT nemo,
const tensorT occ 
) const

◆ compute_IR_intensities()

Tensor< double > madness::Nemo::compute_IR_intensities ( const Tensor< double > &  normalmodes,
const vecfuncT dens_pt 
) const

◆ compute_nemo_potentials()

void madness::Nemo::compute_nemo_potentials ( const vecfuncT nemo,
vecfuncT Jnemo,
vecfuncT Knemo,
vecfuncT xcnemo,
vecfuncT pcmnemo,
vecfuncT Unemo 
) const
protected

compute the reconstructed orbitals, and all potentials applied on nemo

to use these potentials in the fock matrix computation they must be multiplied by the nuclear correlation factor

Parameters
[in]nemothe nemo orbitals
[out]JnemoCoulomb operator applied on the nemos
[out]Knemoexchange operator applied on the nemos
[out]pcmnemoPCM (solvent) potential applied on the nemos
[out]Unemoregularized nuclear potential applied on the nemos

to use these potentials in the fock matrix computation they must be multiplied by the nuclear correlation factor

Parameters
[in]nemothe nemo orbitals
[out]psithe reconstructed, full orbitals
[out]JnemoCoulomb operator applied on the nemos
[out]Knemoexchange operator applied on the nemos
[out]Vnemonuclear potential applied on the nemos
[out]Unemoregularized nuclear potential applied on the nemos

References ac, madness::AC< NDIM >::apply(), calc, madness::CalculationParameters::charge(), madness::PCM::compute_pcm_potential(), do_ac(), do_pcm(), madness::WorldGopInterface::fence(), madness::get_size(), madness::World::gop, K(), madness::XCOperator< T, NDIM >::make_xc_potential(), molecule(), param, pcm, madness::Coulomb< T, NDIM >::potential(), madness::CalculationParameters::print_level(), madness::scale(), madness::Exchange< T, NDIM >::set_symmetric(), madness::Exchange< T, NDIM >::set_taskq(), madness::Coulomb< T, NDIM >::set_taskq(), madness::World::size(), madness::stringify(), madness::timer::tag(), madness::truncate(), and madness::NemoBase::world.

Referenced by compute_fock_matrix(), and solve().

◆ do_ac()

bool madness::Nemo::do_ac ( ) const
inline

◆ do_pcm()

bool madness::Nemo::do_pcm ( ) const
inline

◆ do_symmetry()

bool madness::Nemo::do_symmetry ( ) const
inline

◆ get_ac()

AC<3> madness::Nemo::get_ac ( ) const
inline

References ac.

◆ get_calc()

std::shared_ptr<SCF> madness::Nemo::get_calc ( ) const
inline

◆ get_coulomb_potential()

real_function_3d madness::Nemo::get_coulomb_potential ( const vecfuncT psi) const
protected

◆ get_param()

const NemoCalculationParameters& madness::Nemo::get_param ( ) const
inline

References param.

◆ get_pcm()

PCM madness::Nemo::get_pcm ( ) const
inline

References pcm.

◆ get_symmetry_projector()

projector_irrep madness::Nemo::get_symmetry_projector ( ) const
inline

return the symmetry_projector

References symmetry_projector.

◆ gradient()

Tensor< double > madness::Nemo::gradient ( const Tensor< double > &  x)
virtual

◆ help()

static void madness::Nemo::help ( )
inlinestatic

References madness::print(), and madness::print_header2().

Referenced by main().

◆ hessian()

Tensor< double > madness::Nemo::hessian ( const Tensor< double > &  x)

◆ is_dft()

bool madness::Nemo::is_dft ( ) const
inline

◆ kinetic_energy_potential()

real_function_3d madness::Nemo::kinetic_energy_potential ( const vecfuncT nemo) const

compute the kinetic energy potential using Eq. (16) of R. A. King and N. C. Handy, “Kinetic energy functionals from the Kohn–Sham potential,” Phys. Chem. Chem. Phys., vol. 2, no. 22, pp. 5049–5056, 2000.

References madness::apply(), madness::binary_op(), calc, diff(), madness::dot(), madness::gaxpy(), get_calc(), madness::mul(), madness::norm2(), madness::constants::pi, poisson, madness::print(), madness::NemoBase::R, madness::NemoBase::R_square, madness::save(), madness::scale(), madness::SmoothingOperator3D(), madness::sub(), madness::sum(), and madness::NemoBase::world.

◆ load_function()

template<typename T , size_t NDIM>
void madness::Nemo::load_function ( std::vector< Function< T, NDIM > > &  f,
const std::string  name 
) const
protected

◆ localize()

vecfuncT madness::Nemo::localize ( const vecfuncT nemo,
const double  dconv,
const bool  randomize 
) const
protected

◆ make_cphf_constant_term()

vecfuncT madness::Nemo::make_cphf_constant_term ( const size_t  iatom,
const int  iaxis,
const vecfuncT R2nemo,
const real_function_3d rhonemo 
) const
protected

◆ make_ddensity()

real_function_3d madness::Nemo::make_ddensity ( const real_function_3d rhonemo,
const int  axis 
) const

make the derivative of the density

$ \nabla\rho = 2R^X R \rho_R + R^2\nabla \rho_R $

Parameters
[in]rhonemothe regularized density
[in]axisthe component of the nabla operator
Returns
the gradient of the reconstructed density

References axis, madness::copy(), madness::NemoBase::R_square, and madness::NemoBase::world.

Referenced by hessian().

◆ make_density() [1/2]

real_function_3d madness::Nemo::make_density ( const Tensor< double > &  occ,
const vecfuncT nemo 
) const

◆ make_density() [2/2]

real_function_3d madness::Nemo::make_density ( const tensorT occ,
const vecfuncT bra,
const vecfuncT ket,
const bool  refine = false 
) const

make the density using different bra and ket vectors

e.g. for computing the perturbed density \sum_i \phi_i \phi_i^X or when using nemos: \sum_i R2nemo_i nemo_i

References calc, madness::compress(), madness::WorldGopInterface::fence(), madness::Function< T, NDIM >::gaxpy(), madness::World::gop, madness::mul(), madness::refine(), and madness::NemoBase::world.

◆ make_fock_operator()

std::shared_ptr< Fock< double, 3 > > madness::Nemo::make_fock_operator ( ) const
virtual

◆ make_incomplete_hessian()

Tensor< double > madness::Nemo::make_incomplete_hessian ( ) const
protected

compute the incomplete hessian

incomplete hessian is the nuclear-nuclear contribution, and the contribution from the second derivative of the nuclear potential, and also the derivative of the nuclear correlation factor. i.e. all contributions that do not contain the regularized perturbed density, but it will contain parts of the perturbed density

References get_calc(), madness::inner(), make_density(), molecule(), madness::Molecule::natom(), madness::Molecule::nuclear_repulsion_hessian(), madness::NemoBase::R_square, madness::refine(), and madness::NemoBase::world.

Referenced by compute_all_cphf().

◆ make_incomplete_hessian_response_part()

Tensor< double > madness::Nemo::make_incomplete_hessian_response_part ( const std::vector< vecfuncT > &  xi) const
protected

compute the complementary incomplete hessian

Parameters
[in]xithe response functions including the parallel part

References madness::_(), calc, madness::dot(), h(), madness::inner(), molecule(), madness::Molecule::natom(), madness::NemoBase::R_square, madness::NemoBase::world, and xi.

Referenced by compute_all_cphf().

◆ make_laplacian_density()

real_function_3d madness::Nemo::make_laplacian_density ( const real_function_3d rhonemo) const

the Laplacian of the density

The Laplacian should currently only be used for subsequent convolution with a Green's function (which is reasonably stable), but not on its own!

The Laplacian of the cuspy density is numerically fairly unstable:

  • a singular term may be rewritten using the nuclear potential (see below)
  • the Laplacian of the regularized density is still very noisy

It may be computed as

\[ \Delta \rho = \Delta (R^2 \rho_R) = \Delta (R^2) \rho_R + 2\nabla R \nabla \rho_R + R^2 \Delta \rho_R = 2 R^2 U1^2 \rho_R -4 R^2 ( U-V ) \rho_R + R^2 \Delta\rho_R \]

where we can use the identity

\[ U=V + R^{-1}[T,R] -2 R (U-V) = \Delta R + 2\nabla R\dot \nabla \]

first term comes from the definition of the U potential as the commutator over the kinetic energy (aka the Laplacian)

Parameters
[in]rhonemothe regularized density \rho_R
Returns
the laplacian of the reconstructed density \Delta (R^2\rho_R)

References axis, madness::copy(), get_calc(), madness::constants::pi, poisson, madness::NemoBase::R_square, madness::save(), smoothen(), madness::truncate(), and madness::NemoBase::world.

◆ make_plots()

void madness::Nemo::make_plots ( const real_function_3d f,
const std::string &  name = "function" 
) const
inlineprotected

◆ make_sigma()

real_function_3d madness::Nemo::make_sigma ( const real_function_3d rho1,
const real_function_3d rho2 
) const

compute the reduced densities sigma (gamma) for GGA functionals

the reduced density is given by

\[ \sigma = \nabla\rho_1 \nabla\rho_2 = (\nabla R^2 \rho_{R,1} + R^2 \nabla \rho_{R,1}) \cdot (\nabla R^2 \rho_{R,2} + R^2 \nabla \rho_{R,2}) = 4R^4 U_1^2 \rho_{R,1}\rho_{R,2} + 2R^4 \vec U_1 \left(\rho_{R,1}\nabla \rho_{R,2} + \nabla\rho_{R,1} \rho_{R,2}\right) + R^4 \nabla \rho_{R,1}\cdot \nabla\rho_{R,2} \]

References madness::dot(), madness::FunctionDefaults< NDIM >::get_thresh(), madness::grad(), madness::NemoBase::R_square, madness::FunctionDefaults< NDIM >::set_thresh(), thresh, and madness::NemoBase::world.

◆ molecule() [1/2]

Molecule& madness::Nemo::molecule ( )
inlinevirtual

◆ molecule() [2/2]

Molecule& madness::Nemo::molecule ( ) const
inline

return a reference to the molecule

References calc.

◆ name()

std::string madness::Nemo::name ( ) const
inlinevirtual

Implements madness::QCPropertyInterface.

Reimplemented in madness::OEP.

Referenced by load_function(), make_plots(), and save_function().

◆ print_parameters()

static void madness::Nemo::print_parameters ( )
inlinestatic

◆ provides_gradient()

bool madness::Nemo::provides_gradient ( ) const
inlinevirtual

Override this to return true if the derivative is implemented.

Reimplemented from madness::OptimizationTargetInterface.

◆ purify_hessian()

Tensor< double > madness::Nemo::purify_hessian ( const Tensor< double > &  hessian) const

purify and symmetrize the hessian

The hessian should be symmetric, but it is not, because

\[ \langle i^{Y_B}|H^{X_A}|i\rangle \neq \langle i|H^{X_A}|i^{Y_B}\rangle \]

does holds analytically, but not numerically. If the two numbers differ, pick the more trustworthy, which is the one with a heavy atom causing the perturbed density and the light atom being the nuclear singularity.

Parameters
[in]hessianthe raw hessian
Returns
a symmetrized hessian

References calc, madness::copy(), diff(), hessian(), max, and madness::print().

Referenced by hessian().

◆ rotate_subspace()

template<typename solverT >
void madness::Nemo::rotate_subspace ( World world,
const tensorT U,
solverT &  solver,
int  lo,
int  nfunc 
) const
protected

◆ save_function()

template<typename T , size_t NDIM>
void madness::Nemo::save_function ( const std::vector< Function< T, NDIM > > &  f,
const std::string  name 
) const
protected

◆ selftest()

bool madness::Nemo::selftest ( )
inlinevirtual

Implements madness::QCPropertyInterface.

Reimplemented in madness::OEP.

◆ set_protocol()

void madness::Nemo::set_protocol ( const double  thresh)
inlineprotected

◆ smoothen()

static void madness::Nemo::smoothen ( real_function_3d f)
inlinestatic

smooth a function by projecting it onto k-1 and then average with k

kept it here for further testing

References madness::f, k, and madness::project().

Referenced by make_laplacian_density().

◆ solve()

double madness::Nemo::solve ( const SCFProtocol proto)
protected

◆ solve_cphf()

vecfuncT madness::Nemo::solve_cphf ( const size_t  iatom,
const int  iaxis,
const Tensor< double >  fock,
const vecfuncT guess,
const vecfuncT rhsconst,
const Tensor< double >  incomplete_hessian,
const vecfuncT parallel,
const SCFProtocol proto,
const std::string &  xc_data 
) const

solve the CPHF equations for the nuclear displacements

solve the CPHF equation for the nuclear displacements

this function computes that part of the orbital response that is orthogonal to the occupied space. If no NCF's are used this corresponds to the normal response. If NCF's are used the part parallel to the occupied space must be added!

\[ F^X = F^\perp + F^\parallel \]

cf parallel_CPHF()

Parameters
[in]iatomthe atom A to be moved
[in]iaxisthe coordinate X of iatom to be moved
Returns
\ket{i^X} or \ket{F^\perp}
Parameters
[in]iatomthe atom A to be moved
[in]iaxisthe coordinate X of iatom to be moved
Returns
\frac{\partial}{\partial X_A} \varphi

References madness::_(), madness::Tensor< T >::absmax(), madness::apply(), b, calc, madness::Coulomb< T, NDIM >::compute_potential(), madness::copy(), madness::SCFProtocol::dconv, madness::CalculationParameters::do_localize(), madness::dot(), madness::gaxpy(), get_calc(), guess(), h(), madness::inner(), is_dft(), K(), madness::CalculationParameters::lo(), make_density(), molecule(), madness::mul(), madness::Molecule::natom(), madness::detail::norm(), madness::norm2(), madness::norm2s(), param, madness::Coulomb< T, NDIM >::potential(), madness::print(), Q(), madness::NemoBase::R_square, madness::World::rank(), residual(), madness::scale(), madness::Exchange< T, NDIM >::set_bra_and_ket(), madness::Exchange< T, NDIM >::set_symmetric(), madness::World::size(), madness::CalculationParameters::spin_restricted(), madness::stringify(), madness::transform(), trantol(), madness::truncate(), V(), madness::NemoBase::world, xc, and xi.

Referenced by compute_all_cphf().

◆ trantol()

double madness::Nemo::trantol ( ) const
inlineprotected

return the threshold for vanishing elements in orbital rotations

References calc, and get_calc().

Referenced by rotate_subspace(), solve(), and solve_cphf().

◆ value() [1/2]

virtual double madness::Nemo::value ( )
inlinevirtual

Reimplemented in madness::OEP.

References calc, and value().

Referenced by main(), test_fock(), test_nemo(), and value().

◆ value() [2/2]

double madness::Nemo::value ( const Tensor< double > &  x)
virtual

Friends And Related Function Documentation

◆ PNO

friend class PNO
friend

◆ TDHF

friend class TDHF
friend

Member Data Documentation

◆ ac

AC<3> madness::Nemo::ac
protected

asymptotic correction for DFT

Referenced by compute_nemo_potentials(), get_ac(), and make_fock_operator().

◆ calc

std::shared_ptr<SCF> madness::Nemo::calc
protected

◆ coords_sum

double madness::Nemo::coords_sum
mutableprivate

sum of square of coords at last solved geometry

Referenced by value().

◆ param

NemoCalculationParameters madness::Nemo::param

◆ pcm

PCM madness::Nemo::pcm
private

◆ poisson

std::shared_ptr<real_convolution_3d> madness::Nemo::poisson
protected

◆ symmetry_projector

projector_irrep madness::Nemo::symmetry_projector
protected

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