MADNESS
0.10.1
|
ABC for the nuclear correlation factors. More...
#include <correlationfactor.h>
Classes | |
class | R_functor |
class | RX_functor |
compute the derivative of R wrt the displacement of atom A, coord axis More... | |
class | square_times_V_derivative_functor |
class | square_times_V_functor |
class | U1_atomic_functor |
U1 functor for a specific atom. More... | |
class | U1_dot_U1_functor |
functor for a local U1 dot U1 potential More... | |
class | U1_functor |
functor for the local part of the U1 potential – NOTE THE SIGN More... | |
class | U1X_functor |
compute the derivative of U1 wrt the displacement of atom A, coord axis More... | |
class | U2_atomic_functor |
U2 functor for a specific atom. More... | |
class | U2_functor |
class | U2X_functor |
compute the derivative of U2 wrt the displacement of atom A More... | |
class | U3_atomic_functor |
U3 functor for a specific atom. More... | |
class | U3_functor |
class | U3X_functor |
compute the derivative of U3 wrt the displacement of atom A, coord axis More... | |
Public Types | |
enum | corrfactype { None , GradientalGaussSlater , GaussSlater , LinearSlater , Polynomial , Slater , poly4erfc , Two , Adhoc } |
typedef std::shared_ptr< FunctionFunctorInterface< double, 3 > > | functorT |
Public Member Functions | |
NuclearCorrelationFactor (World &world, const Molecule &mol) | |
ctor More... | |
virtual | ~NuclearCorrelationFactor () |
virtual destructor More... | |
virtual real_function_3d | apply_U (const real_function_3d &rhs) const |
apply the regularized potential U_nuc on a given function rhs More... | |
coord_3d | dsmoothed_unitvec (const coord_3d &xyz, const int axis, double smoothing=0.0) const |
derivative of smoothed unit vector wrt the electronic coordinate More... | |
virtual real_function_3d | function () const |
return the nuclear correlation factor More... | |
void | initialize (const double vtol1) |
initialize the regularized potentials U1 and U2 More... | |
virtual real_function_3d | inverse () const |
return the inverse nuclear correlation factor More... | |
coord_3d | smoothed_unitvec (const coord_3d &xyz, double smoothing=0.0) const |
smoothed unit vector for the computation of the U1 potential More... | |
virtual real_function_3d | square () const |
return the square of the nuclear correlation factor More... | |
virtual real_function_3d | square_times_V_derivative (const int iatom, const int axis) const |
virtual double | Sr_div_S (const double &r, const double &Z) const =0 |
virtual double | Srr_div_S (const double &r, const double &Z) const =0 |
virtual double | Srrr_div_S (const double &r, const double &Z) const =0 |
virtual corrfactype | type () const =0 |
virtual const real_function_3d | U1 (const int axis) const |
return the U1 term of the correlation function More... | |
std::vector< real_function_3d > | U1vec () const |
return the U1 functions in a vector More... | |
virtual const real_function_3d | U2 () const |
return the U2 term of the correlation function More... | |
virtual double | U2X_spherical (const double &r, const double &Z, const double &rcut) const |
derivative of the U2 potential wrt nuclear coordinate X (spherical part) More... | |
Protected Attributes | |
std::vector< real_function_3d > | U1_function |
the three components of the U1 potential More... | |
real_function_3d | U2_function |
the purely local U2 potential, having absorbed the nuclear pot V_nuc More... | |
Private Member Functions | |
virtual double | S (const double &r, const double &Z) const =0 |
the correlation factor S wrt a given atom More... | |
virtual coord_3d | Sp (const coord_3d &vr1A, const double &Z) const =0 |
the partial derivative of correlation factor S' wrt the cartesian coordinates More... | |
virtual double | Spp_div_S (const double &r, const double &Z) const =0 |
the regularized potential wrt a given atom wrt the cartesian coordinate More... | |
Private Attributes | |
double | eprec |
smoothing of the potential/step function More... | |
const Molecule & | molecule |
the molecule More... | |
double | vtol |
the threshold for initial projection More... | |
World & | world |
the world More... | |
ABC for the nuclear correlation factors.
typedef std::shared_ptr< FunctionFunctorInterface<double,3> > madness::NuclearCorrelationFactor::functorT |
|
inline |
ctor
[in] | world | the world |
[in] | mol | molecule with the sites of the nuclei |
|
inlinevirtual |
virtual destructor
|
inlinevirtual |
apply the regularized potential U_nuc on a given function rhs
Reimplemented in madness::PseudoNuclearCorrelationFactor.
References axis, madness::Function< T, NDIM >::compress(), madness::Function< T, NDIM >::truncate(), madness::truncate(), U1(), U2(), and world.
|
inline |
derivative of smoothed unit vector wrt the electronic coordinate
note the sign change for exchanging nuclear and electronic coordinates
the derivative wrt x is given by
References axis, eprec, madness::Vector< T, N >::normf(), p(), madness::constants::pi, and madness::r2().
|
inlinevirtual |
|
inline |
initialize the regularized potentials U1 and U2
References axis, madness::Function< T, NDIM >::set_thresh(), thresh, madness::Function< T, NDIM >::truncate(), U1_function, U2_function, vtol, and world.
|
inlinevirtual |
|
privatepure virtual |
the correlation factor S wrt a given atom
[in] | r | the distance of the req'd coord to the nucleus |
[in] | Z | the nuclear charge |
Implemented in madness::AdhocNuclearCorrelationFactor, madness::PseudoNuclearCorrelationFactor, madness::Polynomial< N >, madness::poly4erfc, madness::Slater, madness::LinearSlater, madness::GradientalGaussSlater, and madness::GaussSlater.
|
inline |
smoothed unit vector for the computation of the U1 potential
note the identity for exchanging nuclear and electronic coordinates (there is a sign change, unlike for the smoothed potential)
References eprec, madness::Molecule::get_eprec(), molecule, madness::Vector< T, N >::normf(), madness::constants::pi, and xi.
Referenced by madness::GaussSlater::Sp(), madness::GradientalGaussSlater::Sp(), madness::LinearSlater::Sp(), madness::Slater::Sp(), madness::poly4erfc::Sp(), and madness::Polynomial< N >::Sp().
|
privatepure virtual |
the partial derivative of correlation factor S' wrt the cartesian coordinates
[in] | vr1A | the vector of the req'd coord to the nucleus |
[in] | Z | the nuclear charge |
Implemented in madness::AdhocNuclearCorrelationFactor, madness::PseudoNuclearCorrelationFactor, madness::Polynomial< N >, madness::poly4erfc, madness::Slater, madness::LinearSlater, madness::GradientalGaussSlater, and madness::GaussSlater.
|
privatepure virtual |
the regularized potential wrt a given atom wrt the cartesian coordinate
S" is the Cartesian Laplacian applied on the NCF. Note the difference to Srr_div_S, which is the second derivative wrt the distance rho. this is: -S"/S - Z/r
[in] | r | the distance of the req'd coord to the nucleus |
[in] | Z | the nuclear charge |
Implemented in madness::AdhocNuclearCorrelationFactor, madness::PseudoNuclearCorrelationFactor, madness::Polynomial< N >, madness::poly4erfc, madness::Slater, madness::LinearSlater, madness::GradientalGaussSlater, and madness::GaussSlater.
|
inlinevirtual |
|
inlinevirtual |
|
pure virtual |
first derivative of the NCF with respect to the relative distance rho
where the distance of the electron to the nucleus A is given by
Implemented in madness::AdhocNuclearCorrelationFactor, madness::PseudoNuclearCorrelationFactor, madness::Polynomial< N >, madness::poly4erfc, madness::Slater, madness::LinearSlater, madness::GradientalGaussSlater, and madness::GaussSlater.
Referenced by madness::NuclearCorrelationFactor::U1_dot_U1_functor::operator()(), madness::GaussSlater::U2X_spherical(), madness::GradientalGaussSlater::U2X_spherical(), madness::Slater::U2X_spherical(), and madness::Polynomial< N >::U2X_spherical().
|
pure virtual |
second derivative of the NCF with respect to the relative distance rho
where the distance of the electron to the nucleus A is given by
Implemented in madness::AdhocNuclearCorrelationFactor, madness::PseudoNuclearCorrelationFactor, madness::Polynomial< N >, madness::poly4erfc, madness::Slater, madness::LinearSlater, madness::GradientalGaussSlater, and madness::GaussSlater.
Referenced by madness::GaussSlater::U2X_spherical(), madness::GradientalGaussSlater::U2X_spherical(), madness::Slater::U2X_spherical(), and madness::Polynomial< N >::U2X_spherical().
|
pure virtual |
third derivative of the NCF with respect to the relative distance rho
where the distance of the electron to the nucleus A is given by
Implemented in madness::AdhocNuclearCorrelationFactor, madness::PseudoNuclearCorrelationFactor, madness::Polynomial< N >, madness::poly4erfc, madness::Slater, madness::LinearSlater, madness::GradientalGaussSlater, and madness::GaussSlater.
Referenced by madness::GaussSlater::U2X_spherical(), madness::GradientalGaussSlater::U2X_spherical(), madness::Slater::U2X_spherical(), and madness::Polynomial< N >::U2X_spherical().
|
pure virtual |
|
inlinevirtual |
return the U1 term of the correlation function
References axis, and U1_function.
Referenced by madness::AdhocNuclearCorrelationFactor::AdhocNuclearCorrelationFactor(), and apply_U().
|
inline |
return the U1 functions in a vector
References U1_function.
|
inlinevirtual |
return the U2 term of the correlation function
Reimplemented in madness::PseudoNuclearCorrelationFactor.
References U2_function.
Referenced by madness::AdhocNuclearCorrelationFactor::AdhocNuclearCorrelationFactor(), apply_U(), and madness::PseudoNuclearCorrelationFactor::apply_U().
|
inlinevirtual |
derivative of the U2 potential wrt nuclear coordinate X (spherical part)
need to reimplement this for all derived classes due to the range for r -> 0, where the singular terms cancel. With
returns the term in the parenthesis without the the derivative of rho
Reimplemented in madness::PseudoNuclearCorrelationFactor, madness::Polynomial< N >, madness::poly4erfc, madness::Slater, madness::GradientalGaussSlater, and madness::GaussSlater.
References MADNESS_EXCEPTION, madness::print(), madness::World::rank(), and world.
|
private |
smoothing of the potential/step function
Referenced by dsmoothed_unitvec(), smoothed_unitvec(), and madness::PseudoNuclearCorrelationFactor::Spp_div_S().
|
private |
the molecule
Referenced by smoothed_unitvec(), and square_times_V_derivative().
|
protected |
the three components of the U1 potential
Referenced by madness::AdhocNuclearCorrelationFactor::AdhocNuclearCorrelationFactor(), initialize(), U1(), and U1vec().
|
protected |
the purely local U2 potential, having absorbed the nuclear pot V_nuc
Referenced by madness::AdhocNuclearCorrelationFactor::AdhocNuclearCorrelationFactor(), initialize(), and U2().
|
private |
the threshold for initial projection
Referenced by function(), initialize(), inverse(), square(), and square_times_V_derivative().
|
private |
the world
Referenced by madness::AdhocNuclearCorrelationFactor::AdhocNuclearCorrelationFactor(), madness::GaussSlater::GaussSlater(), madness::GradientalGaussSlater::GradientalGaussSlater(), madness::LinearSlater::LinearSlater(), madness::poly4erfc::poly4erfc(), madness::Polynomial< N >::Polynomial(), madness::PseudoNuclearCorrelationFactor::PseudoNuclearCorrelationFactor(), madness::Slater::Slater(), apply_U(), function(), initialize(), inverse(), square(), square_times_V_derivative(), and U2X_spherical().