MPQC  3.0.0-alpha

A shell of Gaussian functions. More...

#include <chemistry/qc/basis/gaussshell.h>

Inheritance diagram for sc::GaussianShell:
sc::DescribedClass sc::RefCount sc::GaussianBasisSet::Shell

Public Types

enum  PrimitiveType { Normalized, Unnormalized }
 
enum  GaussianType { Cartesian, Pure }
 

Public Member Functions

 GaussianShell (const std::vector< unsigned int > &am, const std::vector< bool > &puream, const std::vector< double > &exps, const std::vector< double > &contr_coefs, PrimitiveType pt=GaussianShell::Normalized, bool normalize_shell=true)
 Constructs GaussianShell programmatically. More...
 
 GaussianShell (const Ref< KeyVal > &kv, int pure=-1)
 Construct a GaussianShell object from KeyVal input. More...
 
 GaussianShell (const GaussianShell &other)
 Copy constructor (deep :-)
 
DEPRECATED GaussianShell (int ncn, int nprm, double *e, int *am, int *pure, double **c, PrimitiveType pt=GaussianShell::Normalized, bool do_normalize_shell=true)
 
DEPRECATED GaussianShell (int ncn, int nprm, double *e, int *am, GaussianType pure, double **c, PrimitiveType pt=GaussianShell::Normalized)
 
unsigned int nprimitive () const
 The number of primitive Gaussian shells.
 
unsigned int ncontraction () const
 The number of contractions formed from the primitives.
 
unsigned int nfunction () const
 The number of basis functions.
 
int max_angular_momentum () const
 The maximum angular momentum in the shell.
 
int min_angular_momentum () const
 The minimum angular momentum in the shell.
 
int max_cartesian () const
 The maximum number of Cartesian functions in any contraction.
 
const std::vector< unsigned int > & am () const
 The angular momenta of contractions.
 
unsigned int am (int con) const
 The angular momentum of the given contraction.
 
unsigned int max_am () const
 The maximum angular momentum of any contraction.
 
unsigned int min_am () const
 The minimum angular momentum of any contraction.
 
char amchar (int con) const
 The character symbol for the angular momentum of the given contraction.
 
unsigned int nfunction (int con) const
 The number of basis functions coming from the given contraction.
 
unsigned int ncartesian () const
 The total number of functions if this shell was Cartesian.
 
unsigned int ncartesian_with_aminc (int aminc) const
 The total number of Cartesian functions if this shift is applied to all of the angular momentums.
 
unsigned int ncartesian (int con) const
 The number of Cartesian functions for the given contraction.
 
bool is_cartesian (int con) const
 Returns nonzero if contraction con is Cartesian.
 
bool has_cartesian () const
 Returns nonzero if any contraction is Cartesian.
 
bool is_pure (int con) const
 Returns true if contraction con is solid harmonics.
 
const std::vector< bool > & is_pure () const
 Vector of booleans that indicate whether each contraction is solid harmonics.
 
bool has_pure () const
 Returns true if any contraction is solid harmonics.
 
int contraction_to_function (int c) const
 Returns the number of the first function in the given contraction.
 
int function_to_contraction (int f) const
 Returns the contraction to which this function belongs.
 
double coefficient_unnorm (int con, int prim) const
 Returns the contraction coef for unnormalized primitives.
 
double coefficient_norm (int con, int prim) const
 Returns the contraction coef for normalized primitives.
 
const std::vector< double > & coefficient_unnorm_block () const
 returns coefficients for unnormalization primitives, in block form
 
double exponent (int iprim) const
 Returns the exponents of the given primitive.
 
const std::vector< double > & exponents () const
 Returns the exponents.
 
int values (CartesianIter **, SphericalTransformIter **, const SCVector3 &r, double *basis_values)
 Compute the values for this shell at position r. More...
 
int grad_values (CartesianIter **, SphericalTransformIter **, const SCVector3 &R, double *g_values, double *basis_values=0) const
 Like values(...), but computes gradients of the basis function values, too.
 
int hessian_values (CartesianIter **, SphericalTransformIter **, const SCVector3 &R, double *h_values, double *g_values=0, double *basis_values=0) const
 Like values(...), but computes first and second derivatives of the basis function values, too.
 
double relative_overlap (const Ref< Integral > &, int con, int func1, int func2) const
 Returns the intra-generalized-contraction overlap matrix element <con func1|con func2> within an arbitrary constant for the shell.
 
double relative_overlap (int con, int a1, int b1, int c1, int a2, int b2, int c2) const
 Returns the intra-generalized-contraction overlap matrix element <con func1|con func2> within an arbitrary constant for the shell. More...
 
bool equiv (const GaussianShell &s) const
 Returns true if this and the argument are equivalent.
 
double extent (double threshold) const
 Returns a radius. More...
 
double monobound (double r) const
 Returns a bound for the basis function. More...
 
void print (std::ostream &=ExEnv::out0()) const
 Print the object.
 
- 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.
 
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 double epsilon ()
 
static GaussianShell unit ()
 

Static Public Attributes

static const char * amtypes
 
static const char * AMTYPES
 

Friends

void ToStateOut (const GaussianShell &s, StateOut &so, int &count)
 writes GaussianShell to sc::StateOut
 
void FromStateIn (GaussianShell &s, StateIn &si, int &count)
 reads GaussianShell from sc::StateIn
 

Additional Inherited Members

- Protected Member Functions inherited from sc::RefCount
 RefCount (const RefCount &)
 
RefCountoperator= (const RefCount &)
 

Detailed Description

A shell of Gaussian functions.

A shell is a set of functions with same quantum numbers, contraction coefficients, and exponents, and located on the common origin. GaussianShell does include the origin information.

See also
GaussianBasisSet::Shell

Constructor & Destructor Documentation

◆ GaussianShell() [1/4]

sc::GaussianShell::GaussianShell ( const std::vector< unsigned int > &  am,
const std::vector< bool > &  puream,
const std::vector< double > &  exps,
const std::vector< double > &  contr_coefs,
PrimitiveType  pt = GaussianShell::Normalized,
bool  normalize_shell = true 
)

Constructs GaussianShell programmatically.

The constructed shells has ncontr contractions composed of nprim primitive Gaussians.

Parameters
amangular momenta of each contracted shell, std::vector of ncontr nonnegative integers
pureamspherical/cartesian flag, std::vector of ncontr boolean values
expsexponents of each primitive, std::vector of nprim nonnegative double-precision numbers
contr_coefscontraction coefficients for each contraction, std::vector of ncontr times nprim double-precision numbers. contr_coefs[c*ncontr+p] is the contraction coefficient of primitive p in contraction c
ptindicates whether the input contraction coefficients are in terms of normalized primitive functions (i.e. do not include the normalization factors in them).
normalize_shellif true, the coefficients will be scaled to normalize each contraction to unity.

◆ GaussianShell() [2/4]

sc::GaussianShell::GaussianShell ( const Ref< KeyVal > &  kv,
int  pure = -1 
)

Construct a GaussianShell object from KeyVal input.

Parameters
kvthe KeyVal object
purethis is an optional parameter that can be used to override programmatically the "pure" keyword value in kv. If pure=0 Cartesian GaussianShell will be constructed; if pure=1 solid(spherical) harmonics GaussianShell will be constructed; any other value (or omitting this value) will use the kv object.

◆ GaussianShell() [3/4]

DEPRECATED sc::GaussianShell::GaussianShell ( int  ncn,
int  nprm,
double *  e,
int *  am,
int *  pure,
double **  c,
PrimitiveType  pt = GaussianShell::Normalized,
bool  do_normalize_shell = true 
)
Deprecated:
A GaussianShell constructor.

Users of GaussianShell must pass pointers to newed memory that is kept by GaussianShell and deleted by the destructor. The arguments for the following ctor are:

  • ncn is the number of contracted functions (1 except for SP and gen. con.)
  • nprm is the number of primitives
  • e gives the exponents (length nprm)
  • am gives the angular momentum (length ncn)
  • pure is 1 for pure am and 0 for cartesian (length ncn)
  • c are the contraction coefficients (C-style array of arrays!! not single block ncn by nprm)
  • pt describes whether the primitive functions are to be considered normalized or unnormalized. This effects whether or not c is manipulated to give the correct normalization.
  • If do_normalize_shell is true (the default), then the shell normalization constants will be folded into the coefficients.

◆ GaussianShell() [4/4]

DEPRECATED sc::GaussianShell::GaussianShell ( int  ncn,
int  nprm,
double *  e,
int *  am,
GaussianType  pure,
double **  c,
PrimitiveType  pt = GaussianShell::Normalized 
)
Deprecated:
A GaussianShell constructor.

In this ctor pure is either GaussianShell::Cartesian or Gaussian::Pure and all of the contracted functions are treated in that way. (The user doesn\'t need to compute generate a int*pure vector in this case.)

Member Function Documentation

◆ epsilon()

static double sc::GaussianShell::epsilon ( )
inlinestatic
Returns
numerical epsilon value; exponents/coefficients that differ by less than this are the same, similarly coefficients less than that in absolute magntide are zero (zero exponents are OK!!!)

◆ extent()

double sc::GaussianShell::extent ( double  threshold) const

Returns a radius.

All functions in the shell are below threshold outside this radius.

◆ monobound()

double sc::GaussianShell::monobound ( double  r) const

Returns a bound for the basis function.

This bound is defined so that it is positive and monotonically decreasing as a function of r.

◆ relative_overlap()

double sc::GaussianShell::relative_overlap ( int  con,
int  a1,
int  b1,
int  c1,
int  a2,
int  b2,
int  c2 
) const

Returns the intra-generalized-contraction overlap matrix element <con func1|con func2> within an arbitrary constant for the shell.

func1 and func2 are determined according to the axis exponents, a1, b1, c1, a2, b2, and c2.

◆ unit()

static GaussianShell sc::GaussianShell::unit ( )
static
Returns
the "unit" Gaussian shell, zero angular momentum, and zero exponents and unit contraction coefficient

◆ values()

int sc::GaussianShell::values ( CartesianIter **  ,
SphericalTransformIter **  ,
const SCVector3 r,
double *  basis_values 
)

Compute the values for this shell at position r.

The basis_values argument must be vector of length nfunction().


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

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