Namespaces

namespace pcm

Typedefs

typedef boost::container::flat_map<std::string, Eigen::VectorXd> SurfaceFunctionMap
typedef SurfaceFunctionMap::value_type SurfaceFunctionPair
typedef SurfaceFunctionMap::const_iterator SurfaceFunctionMapConstIter

Functions

void printer(const std::string &message)
void printer(const std::ostringstream &stream)
void initMolecule(const Input &inp, const Symmetry &pg, int nuclei, const Eigen::VectorXd &charges, const Eigen::Matrix3Xd &centers, Molecule &molecule)
void initSpheresAtoms(const Input &inp, const Eigen::Matrix3Xd &sphereCenter_, std::vector<Sphere> &spheres_)
unsigned int pcmsolver_get_version(void)
void print(const PCMInput &inp)
template <typename Source>
std::string to_string(const Source &arg)
class Meddle
#include <Meddle.hpp>

Contains functions exposing an interface to the module internals.

Public Functions

Meddle(pcmsolver_reader_t input_reading, int nr_nuclei, double charges[], double coordinates[], int symmetry_info[], const PCMInput &host_input)

Constructor.

The molecular point group information is passed as an array of 4 integers: number of generators, first, second and third generator respectively. Generators map to integers as in table :ref:symmetry-ops

Parameters
  • input_reading -

    input processing strategy

  • nr_nuclei -

    number of atoms in the molecule

  • charges -

    atomic charges

  • coordinates -

    atomic coordinates

  • symmetry_info -

    molecular point group information

  • host_input -

    input to the module, as read by the host

PCMSolverIndex getCavitySize() const

Getter for the number of finite elements composing the molecular cavity.

Return
the size of the cavity

PCMSolverIndex getIrreducibleCavitySize() const

Getter for the number of irreducible finite elements composing the molecular cavity.

Return
the number of irreducible finite elements

void getCenters(double centers[]) const

Getter for the centers of the finite elements composing the molecular cavity.

Parameters
  • centers -

    array holding the coordinates of the finite elements centers

void getCenter(int its, double center[]) const

Getter for the center of the i-th finite element.

Parameters
  • its -

    index of the finite element

  • center -

    array holding the coordinates of the finite element center

void getAreas(double areas[]) const

Getter for the areas/weights of the finite elements.

Parameters
  • context -

    the PCM context object

  • areas -

    array holding the weights/areas of the finite elements

void computeASC(const char *mep_name, const char *asc_name, int irrep) const

Computes ASC given a MEP and the desired irreducible representation.

Parameters
  • mep_name -

    label of the MEP surface function

  • asc_name -

    label of the ASC surface function

  • irrep -

    index of the desired irreducible representation The module uses the surface function concept to handle potentials and charges. Given labels for each, this function retrieves the MEP and computes the corresponding ASC.

void computeResponseASC(const char *mep_name, const char *asc_name, int irrep) const

Computes response ASC given a MEP and the desired irreducible representation.

Parameters
  • mep_name -

    label of the MEP surface function

  • asc_name -

    label of the ASC surface function

  • irrep -

    index of the desired irreducible representation If Nonequilibrium = True in the input, calculates a response ASC using the dynamic permittivity. Falls back to the solver with static permittivity otherwise.

double computePolarizationEnergy(const char *mep_name, const char *asc_name) const

Computes the polarization energy.

Return
the polarization energy This function calculates the dot product of the given MEP and ASC vectors.
Parameters
  • mep_name -

    label of the MEP surface function

  • asc_name -

    label of the ASC surface function

void getSurfaceFunction(PCMSolverIndex size, double values[], const char *name) const

Retrieves data wrapped in a given surface function.

Parameters
  • size -

    the size of the surface function

  • values -

    the values wrapped in the surface function

  • name -

    label of the surface function

void setSurfaceFunction(PCMSolverIndex size, double values[], const char *name) const

Sets a surface function given data and label.

Parameters
  • size -

    the size of the surface function

  • values -

    the values to be wrapped in the surface function

  • name -

    label of the surface function

void saveSurfaceFunctions() const

Dumps all currently saved surface functions to NumPy arrays in .npy files.

void saveSurfaceFunction(const char *name) const

Dumps a surface function to NumPy array in .npy file.

Parameters
  • name -

    label of the surface function

void loadSurfaceFunction(const char *name) const

Loads a surface function from a .npy file.

Parameters
  • name -

    label of the surface function

void printInfo() const

Prints citation and set up information.

void writeTimings() const

Writes timing results for the API.

Private Functions

void initInput(pcmsolver_reader_t input_reading, int nr_nuclei, double charges[], double coordinates[], int symmetry_info[], const PCMInput &host_input)

Initialize input_

void initCavity()

Initialize cavity_

void initStaticSolver()

Initialize static solver K_0_

void initDynamicSolver()

Initialize dynamic solver K_d_

void mediumInfo(IGreensFunction *gf_i, IGreensFunction *gf_o) const

Collect info on medium

Private Members

Input input_

Input object

Cavity *cavity_

Cavity

PCMSolver *K_0_

Solver with static permittivity

PCMSolver *K_d_

Solver with dynamic permittivity

std::ostringstream infoStream_

PCMSolver set up information

bool hasDynamic_

Whether K_d_ was initialized

SurfaceFunctionMap functions_

SurfaceFunction map

namespace interfaces

Typedefs

typedef StateType

state vector for the differential equation integrator

typedef ProfileEvaluator

sort of a function pointer to the dielectric profile evaluation function

struct IntegratorParameters
#include <InterfacesImpl.hpp>

holds parameters for the integrator

Public Members

double eps_abs_

Absolute tolerance level

double eps_rel_

Relative tolerance level

double factor_x_

Weight of the state

double factor_dxdt_

Weight of the state derivative

double r_0_

Lower bound of the integration interval

double r_infinity_

Upper bound of the integration interval

double observer_step_

Time step between observer calls

class LnTransformedRadial
#include <InterfacesImpl.hpp>

system of ln-transformed first-order radial differential equations

Provides a handle to the system of differential equations for the integrator. The dielectric profile comes in as a boost::function object.

Author
Roberto Di Remigio
Date
2015

Public Functions

LnTransformedRadial(const ProfileEvaluator &e, int lval)

Constructor from profile evaluator and angular momentum

void operator()(const StateType &rho, StateType &drhodr, const double r)

Provides a functor for the evaluation of the system of first-order ODEs needed by Boost.Odeint The second-order ODE and the system of first-order ODEs are reported in the manuscript.

Parameters
  • rho -

    state vector holding the function and its first derivative

  • drhodr -

    state vector holding the first and second derivative

  • r -

    position on the integration grid

Private Members

ProfileEvaluator eval_

Dielectric profile function and derivative evaluation

int l_

Angular momentum