MADNESS  0.10.1
Classes | Public Types | Public Member Functions | Protected Attributes | Private Member Functions | Private Attributes | List of all members
madness::NuclearCorrelationFactor Class Referenceabstract

ABC for the nuclear correlation factors. More...

#include <correlationfactor.h>

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

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_3dU1vec () 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_3dU1_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 Moleculemolecule
 the molecule More...
 
double vtol
 the threshold for initial projection More...
 
Worldworld
 the world More...
 

Detailed Description

ABC for the nuclear correlation factors.

Member Typedef Documentation

◆ functorT

Member Enumeration Documentation

◆ corrfactype

Enumerator
None 
GradientalGaussSlater 
GaussSlater 
LinearSlater 
Polynomial 
Slater 
poly4erfc 
Two 
Adhoc 

Constructor & Destructor Documentation

◆ NuclearCorrelationFactor()

madness::NuclearCorrelationFactor::NuclearCorrelationFactor ( World world,
const Molecule mol 
)
inline

ctor

Parameters
[in]worldthe world
[in]molmolecule with the sites of the nuclei

◆ ~NuclearCorrelationFactor()

virtual madness::NuclearCorrelationFactor::~NuclearCorrelationFactor ( )
inlinevirtual

virtual destructor

Member Function Documentation

◆ apply_U()

virtual real_function_3d madness::NuclearCorrelationFactor::apply_U ( const real_function_3d rhs) const
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.

◆ dsmoothed_unitvec()

coord_3d madness::NuclearCorrelationFactor::dsmoothed_unitvec ( const coord_3d xyz,
const int  axis,
double  smoothing = 0.0 
) const
inline

derivative of smoothed unit vector wrt the electronic coordinate

note the sign change for exchanging nuclear and electronic coordinates

\[ \frac{\partial \vec n}{\partial x} = -\frac{\partial \vec n}{\partial X} \]

the derivative wrt x is given by

\[ \frac{\partial\vec n}{\partial x} = \left\{\frac{\left(r^2-x^2\right) \mathrm{erf}\left(\frac{r}{s}\right)}{r^3} +\frac{2 x^2 e^{-\frac{r^2}{s^2}}}{\sqrt{\pi } r^2 s}, x y \left(\frac{2 e^{-\frac{r^2}{s^2}}}{\sqrt{\pi } r^2 s} -\frac{\mathrm{erf}\left(\frac{r}{s}\right)}{r^3}\right), x z \left(\frac{2 e^{-\frac{r^2}{s^2}}}{\sqrt{\pi } r^2 s} -\frac{\mathrm{erf}\left(\frac{r}{s}\right)}{r^3}\right)\right\} \]

References axis, eprec, madness::Vector< T, N >::normf(), p(), madness::constants::pi, and madness::r2().

◆ function()

virtual real_function_3d madness::NuclearCorrelationFactor::function ( ) const
inlinevirtual

return the nuclear correlation factor

References vtol, and world.

◆ initialize()

void madness::NuclearCorrelationFactor::initialize ( const double  vtol1)
inline

◆ inverse()

virtual real_function_3d madness::NuclearCorrelationFactor::inverse ( ) const
inlinevirtual

return the inverse nuclear correlation factor

References vtol, and world.

◆ S()

virtual double madness::NuclearCorrelationFactor::S ( const double &  r,
const double &  Z 
) const
privatepure virtual

the correlation factor S wrt a given atom

Parameters
[in]rthe distance of the req'd coord to the nucleus
[in]Zthe nuclear charge
Returns
the nuclear correlation factor S_A(r_1A)

Implemented in madness::AdhocNuclearCorrelationFactor, madness::PseudoNuclearCorrelationFactor, madness::Polynomial< N >, madness::poly4erfc, madness::Slater, madness::LinearSlater, madness::GradientalGaussSlater, and madness::GaussSlater.

◆ smoothed_unitvec()

coord_3d madness::NuclearCorrelationFactor::smoothed_unitvec ( const coord_3d xyz,
double  smoothing = 0.0 
) const
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)

\[ \vec n = \frac{\partial \rho}{\partial x} = -\frac{\partial \rho}{\partial X} \]

\[ \vec n = \left\{\frac{x \mathrm{erf}\left(\frac{r}{s}\right)}{r}, \frac{y \mathrm{erf}\left(\frac{r}{s}\right)}{r}, \frac{z \mathrm{erf}\left(\frac{r}{s}\right)}{r}\right\} \]

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().

◆ Sp()

virtual coord_3d madness::NuclearCorrelationFactor::Sp ( const coord_3d vr1A,
const double &  Z 
) const
privatepure virtual

the partial derivative of correlation factor S' wrt the cartesian coordinates

Parameters
[in]vr1Athe vector of the req'd coord to the nucleus
[in]Zthe nuclear charge
Returns
the gradient of the nuclear correlation factor S'_A(r_1A)

Implemented in madness::AdhocNuclearCorrelationFactor, madness::PseudoNuclearCorrelationFactor, madness::Polynomial< N >, madness::poly4erfc, madness::Slater, madness::LinearSlater, madness::GradientalGaussSlater, and madness::GaussSlater.

◆ Spp_div_S()

virtual double madness::NuclearCorrelationFactor::Spp_div_S ( const double &  r,
const double &  Z 
) const
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

Parameters
[in]rthe distance of the req'd coord to the nucleus
[in]Zthe nuclear charge
Returns
the Laplacian of the nuclear correlation factor divided by the correlation factor minus the nuclear potential

Implemented in madness::AdhocNuclearCorrelationFactor, madness::PseudoNuclearCorrelationFactor, madness::Polynomial< N >, madness::poly4erfc, madness::Slater, madness::LinearSlater, madness::GradientalGaussSlater, and madness::GaussSlater.

◆ square()

virtual real_function_3d madness::NuclearCorrelationFactor::square ( ) const
inlinevirtual

return the square of the nuclear correlation factor

References R2, vtol, and world.

◆ square_times_V_derivative()

virtual real_function_3d madness::NuclearCorrelationFactor::square_times_V_derivative ( const int  iatom,
const int  axis 
) const
inlinevirtual

return the square of the nuclear correlation factor multiplied with the derivative of the nuclear potential for the specified atom

Returns
R^2 * \frac{\partial Z_A/r_{1A}}{\partial X_A}

References axis, madness::func(), molecule, R2, vtol, and world.

◆ Sr_div_S()

virtual double madness::NuclearCorrelationFactor::Sr_div_S ( const double &  r,
const double &  Z 
) const
pure virtual

◆ Srr_div_S()

virtual double madness::NuclearCorrelationFactor::Srr_div_S ( const double &  r,
const double &  Z 
) const
pure virtual

◆ Srrr_div_S()

virtual double madness::NuclearCorrelationFactor::Srrr_div_S ( const double &  r,
const double &  Z 
) const
pure virtual

◆ type()

virtual corrfactype madness::NuclearCorrelationFactor::type ( ) const
pure virtual

◆ U1()

virtual const real_function_3d madness::NuclearCorrelationFactor::U1 ( const int  axis) const
inlinevirtual

return the U1 term of the correlation function

References axis, and U1_function.

Referenced by madness::AdhocNuclearCorrelationFactor::AdhocNuclearCorrelationFactor(), and apply_U().

◆ U1vec()

std::vector<real_function_3d> madness::NuclearCorrelationFactor::U1vec ( ) const
inline

return the U1 functions in a vector

References U1_function.

◆ U2()

virtual const real_function_3d madness::NuclearCorrelationFactor::U2 ( ) const
inlinevirtual

◆ U2X_spherical()

virtual double madness::NuclearCorrelationFactor::U2X_spherical ( const double &  r,
const double &  Z,
const double &  rcut 
) const
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

\[ \rho = \left| \vec r- \vec R_A \right| \]

returns the term in the parenthesis without the the derivative of rho

\[ \frac{\partial U_2}{\partial X_A} = \frac{\partial \rho}{\partial X} \left(-\frac{1}{2}\frac{S''' S - S'' S'}{S^2} + \frac{1}{\rho^2}\frac{S'}{S} - \frac{1}{\rho} \frac{S''S - S'^2}{S^2} + \frac{Z_A}{\rho^2}\right) \]

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.

Member Data Documentation

◆ eprec

double madness::NuclearCorrelationFactor::eprec
private

smoothing of the potential/step function

Referenced by dsmoothed_unitvec(), smoothed_unitvec(), and madness::PseudoNuclearCorrelationFactor::Spp_div_S().

◆ molecule

const Molecule& madness::NuclearCorrelationFactor::molecule
private

the molecule

Referenced by smoothed_unitvec(), and square_times_V_derivative().

◆ U1_function

std::vector<real_function_3d> madness::NuclearCorrelationFactor::U1_function
protected

the three components of the U1 potential

Referenced by madness::AdhocNuclearCorrelationFactor::AdhocNuclearCorrelationFactor(), initialize(), U1(), and U1vec().

◆ U2_function

real_function_3d madness::NuclearCorrelationFactor::U2_function
protected

the purely local U2 potential, having absorbed the nuclear pot V_nuc

Referenced by madness::AdhocNuclearCorrelationFactor::AdhocNuclearCorrelationFactor(), initialize(), and U2().

◆ vtol

double madness::NuclearCorrelationFactor::vtol
private

the threshold for initial projection

Referenced by function(), initialize(), inverse(), square(), and square_times_V_derivative().

◆ world

World& madness::NuclearCorrelationFactor::world
private

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