Cavities

We will here describe the inheritance hierarchy for generating cavities, in order to use and extend it properly. The runtime creation of cavity objects relies on the Factory Method pattern [GHJV94][Ale01], implemented through the generic Factory class.

../_images/cavity.png

Cavity

class Cavity

Abstract Base Class for cavities.

This class represents a cavity made of spheres, its surface being discretized in terms of finite elements.

Author
Krzysztof Mozgawa
Date
2011

Subclassed by GePolCavity, RestartCavity

Public Functions

Cavity()

Default constructor.

PCMSolver, an API for the Polarizable Continuum Model Copyright (C) 2016 Roberto Di Remigio, Luca Frediani and collaborators.

This file is part of PCMSolver.

PCMSolver is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

PCMSolver is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License along with PCMSolver. If not, see http://www.gnu.org/licenses/.

For information on the complete list of contributors to the PCMSolver API, see: http://pcmsolver.readthedocs.io/

Cavity(const Sphere &sph)

Constructor from a single sphere.

Only used when we have to deal with a single sphere, i.e. in the unit tests

Parameters
  • sph: the sphere

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

Constructor from list of spheres.

Only used when we have to deal with a single sphere, i.e. in the unit tests

Parameters
  • sph: the list of spheres

Cavity(const Molecule &molec)

Constructor from Molecule.

Parameters
  • molec: the molecular aggregate

void saveCavity(const std::string &fname = "cavity.npz")

Save cavity specification to file.

The cavity specification contains: 0. the number of finite elements, nElements;

  1. the weight of the finite elements, elementArea;
  2. the radius of the finite elements, elementRadius;
  3. the centers of the finite elements, elementCenter;
  4. the normal vectors relative to the centers, elementNormal. Each of these objects is saved in a separate .npy binary file and compressed into one .npz file. Notice that this is just the minimal set of data needed to restart an energy calculation.

void loadCavity(const std::string &fname = "cavity.npz")

Load cavity specification from file.

Protected Attributes

std::vector<Sphere> spheres_

List of spheres.

Molecule molecule_

The molecule to be wrapped by the cavity.

PCMSolverIndex nElements_

Number of finite elements generated.

PCMSolverIndex nIrrElements_

Number of irreducible finite elements.

bool built

Whether the cavity has been built.

Eigen::Matrix3Xd elementCenter_

Coordinates of elements centers.

Eigen::Matrix3Xd elementNormal_

Outward-pointing normal vectors to the elements centers.

Eigen::VectorXd elementArea_

Elements areas.

int nSpheres_

Number of spheres.

Eigen::Matrix3Xd elementSphereCenter_

Centers of the sphere the elements belong to.

Eigen::VectorXd elementRadius_

Radii of the sphere the elements belong to.

Eigen::Matrix3Xd sphereCenter_

Spheres centers.

Eigen::VectorXd sphereRadius_

Spheres radii.

std::vector<Element> elements_

List of finite elements.

Symmetry pointGroup_

Molecular point group.

Private Functions

virtual void makeCavity() = 0

Creates the cavity and discretizes its surface.

Has to be implemented by classes lower down in the inheritance hierarchy

GePolCavity

class GePolCavity

A class for GePol cavity.

This class is an interface to the Fortran code PEDRA for the generation of cavities according to the GePol algorithm.

Author
Krzysztof Mozgawa, Roberto Di Remigio
Date
2011, 2016

PCMSolver, an API for the Polarizable Continuum Model Copyright (C) 2016 Roberto Di Remigio, Luca Frediani and collaborators.

This file is part of PCMSolver.

PCMSolver is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

PCMSolver is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License along with PCMSolver. If not, see http://www.gnu.org/licenses/.

For information on the complete list of contributors to the PCMSolver API, see: http://pcmsolver.readthedocs.io/

Inherits from Cavity

Public Functions

GePolCavity(const Molecule &molec, double a, double pr, double minR, const std::string &suffix = "")

PCMSolver, an API for the Polarizable Continuum Model Copyright (C) 2016 Roberto Di Remigio, Luca Frediani and collaborators.

This file is part of PCMSolver.

PCMSolver is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

PCMSolver is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License along with PCMSolver. If not, see http://www.gnu.org/licenses/.

For information on the complete list of contributors to the PCMSolver API, see: http://pcmsolver.readthedocs.io/

Private Functions

virtual void makeCavity()

Creates the cavity and discretizes its surface.

Has to be implemented by classes lower down in the inheritance hierarchy

void build(const std::string &suffix, int maxts, int maxsp, int maxvert)

Driver for PEDRA Fortran module.

Parameters
  • suffix: for the cavity.off and PEDRA.OUT files, the PID will also be added
  • maxts: maximum number of tesserae
  • maxsp: maximum number of spheres (original + added)
  • maxvert: maximum number of vertices

void writeOFF(const std::string &suffix)

Writes the cavity.off file for visualizing the cavity.

Parameters
  • suffix: for the cavity.off The full name of the visualization file will be cavity.off_suffix_PID

RestartCavity

class RestartCavity

A class for Restart cavity.

Author
Roberto Di Remigio
Date
2014

Inherits from Cavity

Public Functions

virtual void makeCavity()

Creates the cavity and discretizes its surface.

Has to be implemented by classes lower down in the inheritance hierarchy

Private Functions

std::ostream &printCavity(std::ostream &os)

PCMSolver, an API for the Polarizable Continuum Model Copyright (C) 2016 Roberto Di Remigio, Luca Frediani and collaborators.

This file is part of PCMSolver.

PCMSolver is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

PCMSolver is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License along with PCMSolver. If not, see http://www.gnu.org/licenses/.

For information on the complete list of contributors to the PCMSolver API, see: http://pcmsolver.readthedocs.io/