Helper classes and functions

../_images/utils.png

Input

class Input

A wrapper class for the Getkw Library C++ bindings.

An Input object is to be used as the unique point of access to user-provided input: input > parsed input (Python script) > Input object (contains all the input data) Definition of input parameters is to be done in the Python script and in this class. They must be specified as private data members with public accessor methods (get-ters). Most of the data members are anyway accessed through the input wrapping struct-s In general, no mutator methods (set-ters) should be needed, exceptions to this rule should be carefully considered.

Author
Roberto Di Remigio
Date
2013

Public Functions

Input()

Default constructor.

Input(const std::string &filename)

Constructor from human-readable input file name.

Input(const PCMInput &host_input)

Constructor from host input structs.

std::string units() const

Accessor methods.

Top-level section input

std::string cavityType() const

Cavity section input.

void molecule(const Molecule &m)

This method sets the molecule and the list of spheres.

Solvent solvent() const

Medium section input.

std::string greenInsideType() const

Green’s function section input.

std::string providedBy() const

Keeps track of who did the parsing: the API or the host program.

cavityData cavityParams()

Get-ters for input wrapping structs.

Private Functions

void reader(const PCMInput &host_input)

Read host data structures (host-side syntactic input parsing) into Input object. It provides access to a limited number of options only, basically the ones that can be filled into the cavityInput, solverInput and greenInput data structures. Lengths and areas are expected to be in Angstrom/Angstrom^2 and will hence be converted to au/au^2.

Note
Specification of the solvent by name overrides any input given through the greenInput data structure!
Warning
The cavity can only be built in the “Implicit” mode, i.e. by grabbing the coordinates for the sphere centers from the host program. Atomic coordinates are expected to be in au! The “Atoms” and “Explicit” methods are only available using the explicit parsing by our Python script of a separate input file.

void semanticCheck() const

Perform semantic input parsing aka sanity check

Private Members

std::string units_

Units of measure.

int CODATAyear_

Year of the CODATA set to be used.

std::string type_

The type of cavity.

std::string cavFilename_

Filename for the .npz cavity restart file.

std::string dyadicFilename_

Filename for the wavelet cavity dyadic file.

int patchLevel_

Wavelet cavity patch level.

double coarsity_

Wavelet cavity coarsity.

double area_

GePol and TsLess cavities average element area.

double minDistance_

TsLess cavity minimal distance between sampling points.

int derOrder_

TsLess cavity maximum derivative order of switch function.

bool scaling_

Whether the radii should be scaled by 1.2.

std::string radiiSet_

The set of radii to be used.

double minimalRadius_

Minimal radius of an added sphere.

std::string mode_

How the API should get the coordinates of the sphere centers.

std::vector<int> atoms_

List of selected atoms with custom radii.

std::vector<double> radii_

List of radii attached to the selected atoms.

std::vector<Sphere> spheres_

List of spheres for fully custom cavity generation.

Molecule molecule_

Molecule or atomic aggregate.

Solvent solvent_

The solvent for a vacuum/uniform dielectric run.

bool hasSolvent_

Whether the medium was initialized from a solvent object.

std::string solverType_

The solver type.

int equationType_

The integral equation type (wavelet solvers)

double correction_

Correction factor (C-PCM)

bool hermitivitize_

Whether the PCM matrix should be hermitivitized (collocation solvers)

bool isDynamic_

Whether the dynamic PCM matrix should be used.

double probeRadius_

Solvent probe radius.

int integratorType_

Type of integrator for the diagonal of the boundary integral operators.

std::string greenInsideType_

The Green’s function type inside the cavity.

std::string greenOutsideType_

The Green’s function type outside the cavity.

int derivativeInsideType_

How to calculate Green’s function derivatives inside the cavity.

int derivativeOutsideType_

How to calculate Green’s function derivatives outside the cavity.

double epsilonInside_

Permittivity inside the cavity.

double epsilonStaticOutside_

Static permittivity outside the cavity.

double epsilonDynamicOutside_

Dynamic permittivity outside the cavity.

double epsilonReal_

Real part of the metal NP permittivity.

double epsilonImaginary_

Imaginary part of the metal NP permittivity.

std::vector<double> spherePosition_

Center of the spherical metal NP.

double sphereRadius_

Radius of the spherical metal NP.

double epsilonStatic1_

Diffuse interface: static permittivity inside the interface.

double epsilonDynamic1_

Diffuse interface: dynamic permittivity inside the interface.

double epsilonStatic2_

Diffuse interface: static permittivity outside the interface.

double epsilonDynamic2_

Diffuse interface: dynamic permittivity outside the interface.

double center_

Center of the diffuse interface.

double width_

Width of the diffuse interface.

int profileType_

Profile chosen for the diffuse interface.

int maxL_

Maximum angular momentum.

std::vector<double> origin_

Center of the dielectric sphere.

std::vector<double> geometry_

Molecular geometry.

std::string providedBy_

Who performed the syntactic input parsing.

cavityData cavData_

Input wrapping struct for the cavity.

greenData insideGreenData_

Input wrapping struct for the Green’s function inside.

greenData outsideStaticGreenData_

Input wrapping struct for the Green’s function outside (static permittivity)

greenData outsideDynamicGreenData_

Input wrapping struct for the Green’s function outside (dynamic permittivity)

solverData solverData_

Input wrapping struct for the solver.

Friends

std::ostream &operator<<(std::ostream &os, const Input &input)

Operators operator<<

Sphere

class Sphere

Class describing a sphere.

Author
Roberto Di Remigio
Date
2011

Public Functions

void scale(double scaling)

Scale sphere to other units.

Atom

class Atom

Class describing an atom.

This class contains all the radii sets available in the module. They can be obtained by the client through static functions initializing the sets via the lazy evaluation idiom.

Author
Roberto Di Remigio
Date
2011

Public Functions

double atomRadius() const

Returns the atomic radius in Angstrom.

Public Static Functions

std::vector<Atom> &initBondi()

Returns a reference to a vector<Atom> containing Bondi van der Waals radii.

The van der Waals radii are taken from: A. Bondi, J. Phys. Chem. 68, 441-451 (1964) complemented with the ones reported in: M. Mantina, A. C. Chamberlin, R. Valero, C. J. Cramer, D. G. Truhlar, J. Phys. Chem. A, 113, 5806-5812 (2009) We are here using Angstrom as in the papers.

std::vector<Atom> &initUFF()

Returns a reference to a vector<Atom> containing UFF radii.

The UFF set of radii is taken from: A. Rappé, C. J. Casewit, K. S. Colwell, W. A. Goddard, W. M. Skiff, J. Am. Chem. Soc., 114, 10024-10035 (1992) We are here using Angstrom as in the paper.

Molecule

class Molecule

Class representing a molecule or general aggregate of atoms.

This class is based on the similar class available in the Mints library of Psi4

Author
Roberto Di Remigio
Date
2014

Unnamed Group

Molecule &operator=(const Molecule &other)

Operators Assignment operator.

Public Functions

Molecule()

Default constructor Initialize a dummy molecule, e.g. as placeholder, see Cavity.cpp loadCavity method.

Molecule(int nat, const Eigen::VectorXd &chg, const Eigen::VectorXd &masses, const Eigen::Matrix3Xd &geo, const std::vector<Atom> &at, const std::vector<Sphere> &sph)

Constructor from full molecular data.

This initializes the molecule in C1 symmetry

Parameters
  • nat -

    number of atoms

  • chg -

    vector of atomic charges

  • masses -

    vector of atomic masses

  • geo -

    molecular geometry (format nat*3)

  • at -

    vector of Atom objects

  • sph -

    vector of Sphere objects

Molecule(int nat, const Eigen::VectorXd &chg, const Eigen::VectorXd &masses, const Eigen::Matrix3Xd &geo, const std::vector<Atom> &at, const std::vector<Sphere> &sph, int nr_gen, int gen[3])

Constructor from full molecular data, plus number of generators and generators.

This initializes the molecule in the symmetry prescribed by nr_gen and gen. See documentation of the Symmetry object for the conventions.

Parameters
  • nat -

    number of atoms

  • chg -

    vector of atomic charges

  • masses -

    vector of atomic masses

  • geo -

    molecular geometry (format nat*3)

  • at -

    vector of Atom objects

  • sph -

    vector of Sphere objects

  • nr_gen -

    number of molecular point group generators

  • gen -

    molecular point group generators

Molecule(int nat, const Eigen::VectorXd &chg, const Eigen::VectorXd &masses, const Eigen::Matrix3Xd &geo, const std::vector<Atom> &at, const std::vector<Sphere> &sph, const Symmetry &pg)

Constructor from full molecular data and point group.

This initializes the molecule in the symmetry prescribed by pg.

Parameters
  • nat -

    number of atoms

  • chg -

    vector of atomic charges

  • masses -

    vector of atomic masses

  • geo -

    molecular geometry (format nat*3)

  • at -

    vector of Atom objects

  • sph -

    vector of Sphere objects

  • pg -

    the molecular point group (a Symmetry object)

Molecule(const std::vector<Sphere> &sph)

Constructor from list of spheres.

Molecule is treated as an aggregate of spheres. We do not have information on the atomic species involved in the aggregate. Charges are set to 1.0; masses are set based on the radii; geometry is set from the list of spheres. All the atoms are dummy atoms. The point group is C1.

Warning
This constructor is to be used exclusively when initializing the Molecule in EXPLICIT mode, i.e. when the user specifies explicitly spheres centers and radii.
Parameters
  • sph -

    list of spheres

Molecule(const Molecule &other)

Copy constructor.

void translate(const Eigen::Vector3d &translationVector)

Given a vector, carries out translation of the molecule.

Parameters
  • translationVector -

    The translation vector.

void moveToCOM()

Performs translation to the Center of Mass Frame.

void rotate(const Eigen::Matrix3d &rotationMatrix)

Given a matrix, carries out rotation of the molecule.

Parameters
  • rotationMatrix -

    The matrix representing the rotation.

void moveToPAF()

Performs rotation to the Principal Axes Frame.

Private Members

size_t nAtoms_

The number of atoms in the molecule.

Eigen::VectorXd charges_

A vector of dimension (# atoms) containing the charges.

Eigen::VectorXd masses_

A vector of dimension (# atoms) containing the masses.

Eigen::Matrix3Xd geometry_

Molecular geometry, in cartesian coordinates. The dimensions are (# atoms * 3) Units are Bohr.

std::vector<Atom> atoms_

A container for all the atoms composing the molecule.

std::vector<Sphere> spheres_

A container for the spheres composing the molecule.

rotorType rotor_

The molecular rotor type.

Symmetry pointGroup_

The molecular point group.

Solvent

class Solvent

Class describing a solvent.

A Solvent object contains all the solvent-related experimental data needed to set up the Green’s functions and the non-electrostatic terms calculations.

Author
Roberto Di Remigio
Date
2011

Public Types

typedef std::map<std::string, Solvent> SolventMap

typedef for the map between solvent name and Solvent object.

Public Static Functions

Solvent::SolventMap &initSolventMap()

Returns the map between solvent names and Solvent objects.

This map contains solvent data taken from the DALTON2011 internal implementation of the Polarizable Continuum Model.

SurfaceFunction

class SurfaceFunction

A basic surface function class.

This class is basically a wrapper around vectors containing electrostatic potentials and apparent charges. Use judiciously, i.e. DO NOT use it directly in the core classes (cavities, solvers) to avoid high coupling.

Author
Luca Frediani, Roberto Di Remigio
Date
2012, 2013

Public Functions

SurfaceFunction(const SurfaceFunction &other)

Copy constructor.

SurfaceFunction &operator=(SurfaceFunction other)

Assignment operator.

double operator*(const SurfaceFunction &other) const

Multiplication operator: product of two SurfaceFunctions version (scalar product of the values vectors).

SurfaceFunction &operator+=(const SurfaceFunction &other)

Addition-assignment operator.

SurfaceFunction &operator-=(const SurfaceFunction &other)

Subtraction-assignment operator.

SurfaceFunction &operator*=(double scaling)

Multiplication-assignment operator. Defined only for the uniform scaling case.

SurfaceFunction &operator/=(double scaling)

Division-assignment operator. Defined only for the uniform scaling case.

Symmetry

class Symmetry

Contains very basic info about symmetry (only Abelian groups)

Just a wrapper around a vector containing the generators of the group

Author
Roberto Di Remigio
Date
2014

Public Functions

Symmetry()

Default constructor sets up C1 point group.

Private Members

int nrGenerators_

Number of generators

int generators_[3]

Generators

int nrIrrep_

Number of irreps

Mathematical utilities

Functions

template <size_t nBits>
int parity(std::bitset<nBits> bitrep)

Calculate the parity of the bitset as defined by: bitrep[0] XOR bitrep[1] XOR ... XOR bitrep[nBits-1]

Parameters
  • bitrep -

    a bitset

Template Parameters
  • nBits -

    lenght of the input bitset

double parity(unsigned int i)

Returns parity of input integer. The parity is defined as the result of using XOR on the bitrep of the given integer. For example: 2 -> 010 -> 0^1^0 = 1 -> -1.0 6 -> 110 -> 1^1^0 = 0 -> 1.0

Parameters
  • i -

    an integer, usually an index for an irrep or a symmetry operation

It can also be interpreted as the action of a given operation on the Cartesian axes: zyx Parity 0 000 E 1.0 1 001 Oyz -1.0 2 010 Oxz -1.0 3 011 C2z 1.0 4 100 Oxy -1.0 5 101 C2y 1.0 6 110 C2x 1.0 7 111 i -1.0

bool isZero(double value, double threshold)

Returns true if value is less or equal to threshold

Parameters
  • value -

    the value to be checked

  • threshold -

    the threshold

bool numericalZero(double value)

Returns true if value is less than 1.0e-14

Parameters
  • value -

    the value to be checked

void symmetryBlocking(Eigen::MatrixXd &matrix, size_t cavitySize, int ntsirr, int nr_irrep)
void symmetryPacking(std::vector<Eigen::MatrixXd> &blockedMatrix, const Eigen::MatrixXd &fullMatrix, int dimBlock, int nrBlocks)

Parameters
  • blockedMatrix -

    the result of packing fullMatrix

  • fullMatrix -

    the matrix to be packed

  • dimBlock -

    the dimension of the square blocks

  • nrBlocks -

    the number of square blocks

template <typename Derived>
void hermitivitize(Eigen::MatrixBase<Derived> &matrix_)

Given matrix_ returns 0.5 * (matrix_ + matrix_^dagger)

Parameters
  • matrix_ -

    the matrix to be hermitivitized

Template Parameters
  • Derived -

    the numeric type of matrix_ elements

template <typename T>
Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic> getFromRawBuffer(size_t _rows, size_t _columns, void *_inData)

Returns an Eigen matrix of type T, with dimensions _rows*_columns.

Warning! This function assumes that the raw buffer is in column-major order as in Fortran.

Parameters
  • _rows -

    the number of rows.

  • _columns -

    the number of columns.

  • _inData -

    the raw data buffer.

Template Parameters
  • T -

    the data type of the matrix to be returned.

void eulerRotation(Eigen::Matrix3d &R_, const Eigen::Vector3d &eulerAngles_)

Build rotation matrix between two reference frames given the Euler angles.

We assume the convention \( R = Z_3 X_2 Z_1 \) for the ordering of the extrinsic elemental rotations (see http://en.wikipedia.org/wiki/Euler_angles) The Euler angles are given in the order \( \phi, \theta, \psi \). If we write \( c_i, s_i \,\, i = 1, 3 \) for their cosines and sines the rotation matrix will be:

\[\begin{split} R = \begin{pmatrix} c_1c_3 - s_1c_2s_3 & -s_1c_3 - c_1c_2s_3 & s_2s_3 \\ c_1s_3 + s_1c_2c_3 & -s_1s_3 + c_1c_2c_3 & -s_2c_3 \\ s_1s_2 & c_1s_2 & c_2 \end{pmatrix} \end{split}\]
Eigen’s geometry module is used to calculate the rotation matrix
Parameters
  • R_ -

    the rotation matrix

  • eulerAngles_ -

    the Euler angles, in degrees, describing the rotation

double linearInterpolation(const double point, const std::vector<double> &grid, const std::vector<double> &function)

Return value of function defined on grid at an arbitrary point.

This function finds the nearest values for the given point and performs a linear interpolation.

Warning
This function assumes that grid has already been sorted!
Parameters
  • point -

    where the function has to be evaluated

  • grid -

    holds points on grid where function is known

  • function -

    holds known function values

double splineInterpolation(const double point, const std::vector<double> &grid, const std::vector<double> &function)

Return value of function defined on grid at an arbitrary point.

This function finds the nearest values for the given point and performs a cubic spline interpolation.

Warning
This function assumes that grid has already been sorted!
Parameters
  • point -

    where the function has to be evaluated

  • grid -

    holds points on grid where function is known

  • function -

    holds known function values