Green’s Functions

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

../_images/green.png

IGreensFunction

class IGreensFunction

Interface for Green’s function classes.

We define as Green’s function a function:

\[ G(\mathbf{r}, \mathbf{r}^\prime) : \mathbb{R}^6 \rightarrow \mathbb{R} \]
Green’s functions and their directional derivatives appear as kernels of the \(\mathcal{S}\) and \(\mathcal{D}\) integral operators. Forming the matrix representation of these operators requires performing integrations over surface finite elements. Since these Green’s functions present a Coulombic divergence, the diagonal elements of the operators will diverge unless appropriately formulated. This is possible, but requires explicit access to the subtype of this abstract base object. This justifies the need for the singleLayer and doubleLayer functions. The code uses the Non-Virtual Interface (NVI) idiom.
Author
Luca Frediani and Roberto Di Remigio
Date
2012-2016

Subclassed by GreensFunction< DerivativeTraits, ProfilePolicy >, GreensFunction< DerivativeTraits, Anisotropic >, GreensFunction< DerivativeTraits, IntegratorPolicy, Sharp, SphericalSharp< DerivativeTraits, IntegratorPolicy > >, GreensFunction< DerivativeTraits, Uniform >, GreensFunction< DerivativeTraits, Yukawa >, GreensFunction< Numerical, ProfilePolicy >

Unnamed Group

double kernelS(const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Methods to sample the Green’s function and its probe point directional derivative

Returns value of the kernel of the \(\mathcal{S}\) integral operator, i.e. the value of the Greens’s function for the pair of points p1, p2: \( G(\mathbf{p}_1, \mathbf{p}_2)\)

Note
This is the Non-Virtual Interface (NVI)
Parameters
  • p1: first point
  • p2: second point

double kernelD(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the kernel of the \(\mathcal{D}\) integral operator for the pair of points p1, p2: \( [\boldsymbol{\varepsilon}\nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)]\cdot \mathbf{n}_{\mathbf{p}_2}\) To obtain the kernel of the \(\mathcal{D}^\dagger\) operator call this methods with \(\mathbf{p}_1\) and \(\mathbf{p}_2\) exchanged and with \(\mathbf{n}_{\mathbf{p}_2} = \mathbf{n}_{\mathbf{p}_1}\)

Note
This is the Non-Virtual Interface (NVI)
Parameters
  • direction: the direction
  • p1: first point
  • p2: second point

Unnamed Group

double singleLayer(const Element &e, double factor) const

Methods to compute the diagonal of the matrix representation of the S and D operators by approximate collocation.

Calculates an element on the diagonal of the matrix representation of the S operator using an approximate collocation formula.

Note
This is the Non-Virtual Interface (NVI)
Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

double doubleLayer(const Element &e, double factor) const

Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.

Note
This is the Non-Virtual Interface (NVI)
Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

Unnamed Group

virtual double kernelS_impl(const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const = 0

Methods to sample the Green’s function and its probe point directional derivative

Returns value of the kernel of the \(\mathcal{S}\) integral operator, i.e. the value of the Greens’s function for the pair of points p1, p2: \( G(\mathbf{p}_1, \mathbf{p}_2)\)

Parameters
  • p1: first point
  • p2: second point

virtual double kernelD_impl(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const = 0

Returns value of the kernel of the \(\mathcal{D}\) integral operator for the pair of points p1, p2: \( [\boldsymbol{\varepsilon}\nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)]\cdot \mathbf{n}_{\mathbf{p}_2}\) To obtain the kernel of the \(\mathcal{D}^\dagger\) operator call this methods with \(\mathbf{p}_1\) and \(\mathbf{p}_2\) exchanged and with \(\mathbf{n}_{\mathbf{p}_2} = \mathbf{n}_{\mathbf{p}_1}\)

Parameters
  • direction: the direction
  • p1: first point
  • p2: second point

Unnamed Group

virtual double singleLayer_impl(const Element &e, double factor) const = 0

Methods to compute the diagonal of the matrix representation of the S and D operators by approximate collocation.

Calculates an element on the diagonal of the matrix representation of the S operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

virtual double doubleLayer_impl(const Element &e, double factor) const = 0

Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

Public Functions

virtual bool uniform() const = 0

Whether the Green’s function describes a uniform environment

virtual Permittivity permittivity() const = 0

Returns a dielectric permittivity profile

GreensFunction

template <typename DerivativeTraits, typename ProfilePolicy>
class GreensFunction

Templated interface for Green’s functions.

Author
Luca Frediani and Roberto Di Remigio
Date
2012-2016
Template Parameters
  • DerivativeTraits: evaluation strategy for the function and its derivatives
  • ProfilePolicy: dielectric profile type

Inherits from IGreensFunction

Unnamed Group

virtual double derivativeSource(const Eigen::Vector3d &normal_p1, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Methods to sample the Green’s function directional derivatives

Returns value of the directional derivative of the Greens’s function for the pair of points p1, p2: \( \nabla_{\mathbf{p_1}}G(\mathbf{p}_1, \mathbf{p}_2)\cdot \mathbf{n}_{\mathbf{p}_1}\) Notice that this method returns the directional derivative with respect to the source point.

Parameters
  • normal_p1: the normal vector to p1
  • p1: first point
  • p2: second point

virtual double derivativeProbe(const Eigen::Vector3d &normal_p2, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the directional derivative of the Greens’s function for the pair of points p1, p2: \( \nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)\cdot \mathbf{n}_{\mathbf{p}_2}\) Notice that this method returns the directional derivative with respect to the probe point.

Parameters
  • normal_p2: the normal vector to p2
  • p1: first point
  • p2: second point

Unnamed Group

Eigen::Vector3d gradientSource(const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Methods to sample the Green’s function gradients

Returns full gradient of Greens’s function for the pair of points p1, p2: \( \nabla_{\mathbf{p_1}}G(\mathbf{p}_1, \mathbf{p}_2)\) Notice that this method returns the gradient with respect to the source point.

Parameters
  • p1: first point
  • p2: second point

Eigen::Vector3d gradientProbe(const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns full gradient of Greens’s function for the pair of points p1, p2: \( \nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)\) Notice that this method returns the gradient with respect to the probe point.

Parameters
  • p1: first point
  • p2: second point

Public Functions

virtual bool uniform() const

Whether the Green’s function describes a uniform environment

virtual Permittivity permittivity() const

Returns a dielectric permittivity profile

Protected Functions

virtual DerivativeTraits operator()(DerivativeTraits *source, DerivativeTraits *probe) const = 0

Evaluates the Green’s function given a pair of points

Parameters
  • source: the source point
  • probe: the probe point

virtual double kernelS_impl(const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the kernel of the \(\mathcal{S}\) integral operator, i.e. the value of the Greens’s function for the pair of points p1, p2: \( G(\mathbf{p}_1, \mathbf{p}_2)\)

Note
Relies on the implementation of operator() in the subclasses and that is all subclasses need to implement. Thus this method is marked __final.
Parameters
  • p1: first point
  • p2: second point

Protected Attributes

double delta_

Step for numerical differentiation.

ProfilePolicy profile_

Permittivity profile.

Vacuum

template <typename DerivativeTraits = AD_directional>
class Vacuum

Green’s function for vacuum.

Author
Luca Frediani and Roberto Di Remigio
Date
2012-2016
Template Parameters
  • DerivativeTraits: evaluation strategy for the function and its derivatives

Inherits from GreensFunction< DerivativeTraits, Uniform >

Public Functions

Vacuum()

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

DerivativeTraits operator()(DerivativeTraits *source, DerivativeTraits *probe) const

Evaluates the Green’s function given a pair of points

Parameters
  • source: the source point
  • probe: the probe point

double kernelD_impl(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the kernel of the \(\mathcal{D}\) integral operator for the pair of points p1, p2: \( [\boldsymbol{\varepsilon}\nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)]\cdot \mathbf{n}_{\mathbf{p}_2}\) To obtain the kernel of the \(\mathcal{D}^\dagger\) operator call this methods with \(\mathbf{p}_1\) and \(\mathbf{p}_2\) exchanged and with \(\mathbf{n}_{\mathbf{p}_2} = \mathbf{n}_{\mathbf{p}_1}\)

Parameters
  • direction: the direction
  • p1: first point
  • p2: second point

double singleLayer_impl(const Element &e, double factor) const

Methods to compute the diagonal of the matrix representation of the S and D operators by approximate collocation.

Calculates an element on the diagonal of the matrix representation of the S operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

double doubleLayer_impl(const Element &e, double factor) const

Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

UniformDielectric

template <typename DerivativeTraits = AD_directional>
class UniformDielectric

Green’s function for uniform dielectric.

Author
Luca Frediani and Roberto Di Remigio
Date
2012-2016
Template Parameters
  • DerivativeTraits: evaluation strategy for the function and its derivatives

Inherits from GreensFunction< DerivativeTraits, Uniform >

Public Functions

UniformDielectric(double eps)

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

DerivativeTraits operator()(DerivativeTraits *source, DerivativeTraits *probe) const

Evaluates the Green’s function given a pair of points

Parameters
  • source: the source point
  • probe: the probe point

double kernelD_impl(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the kernel of the \(\mathcal{D}\) integral operator for the pair of points p1, p2: \( [\boldsymbol{\varepsilon}\nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)]\cdot \mathbf{n}_{\mathbf{p}_2}\) To obtain the kernel of the \(\mathcal{D}^\dagger\) operator call this methods with \(\mathbf{p}_1\) and \(\mathbf{p}_2\) exchanged and with \(\mathbf{n}_{\mathbf{p}_2} = \mathbf{n}_{\mathbf{p}_1}\)

Parameters
  • direction: the direction
  • p1: first point
  • p2: second point

double singleLayer_impl(const Element &e, double factor) const

Methods to compute the diagonal of the matrix representation of the S and D operators by approximate collocation.

Calculates an element on the diagonal of the matrix representation of the S operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

double doubleLayer_impl(const Element &e, double factor) const

Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

IonicLiquid

template <typename DerivativeTraits = AD_directional>
class IonicLiquid

Green’s functions for ionic liquid, described by the linearized Poisson-Boltzmann equation.

Author
Luca Frediani, Roberto Di Remigio
Date
2013-2016
Template Parameters
  • DerivativeTraits: evaluation strategy for the function and its derivatives

Inherits from GreensFunction< DerivativeTraits, Yukawa >

Public Functions

IonicLiquid(double eps, double k)

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

DerivativeTraits operator()(DerivativeTraits *source, DerivativeTraits *probe) const

Evaluates the Green’s function given a pair of points

Parameters
  • source: the source point
  • probe: the probe point

double kernelD_impl(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the kernel of the \(\mathcal{D}\) integral operator for the pair of points p1, p2: \( [\boldsymbol{\varepsilon}\nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)]\cdot \mathbf{n}_{\mathbf{p}_2}\) To obtain the kernel of the \(\mathcal{D}^\dagger\) operator call this methods with \(\mathbf{p}_1\) and \(\mathbf{p}_2\) exchanged and with \(\mathbf{n}_{\mathbf{p}_2} = \mathbf{n}_{\mathbf{p}_1}\)

Parameters
  • direction: the direction
  • p1: first point
  • p2: second point

double singleLayer_impl(const Element &e, double factor) const

Methods to compute the diagonal of the matrix representation of the S and D operators by approximate collocation.

Calculates an element on the diagonal of the matrix representation of the S operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

double doubleLayer_impl(const Element &e, double factor) const

Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

AnisotropicLiquid

template <typename DerivativeTraits = AD_directional>
class AnisotropicLiquid

Green’s functions for anisotropic liquid, described by a tensorial permittivity.

Author
Roberto Di Remigio
Date
2016
Template Parameters
  • DerivativeTraits: evaluation strategy for the function and its derivatives

Inherits from GreensFunction< DerivativeTraits, Anisotropic >

Public Functions

AnisotropicLiquid(const Eigen::Vector3d &eigen_eps, const Eigen::Vector3d &euler_ang)

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

Parameters
  • eigen_eps: eigenvalues of the permittivity tensors
  • euler_ang: Euler angles in degrees

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

DerivativeTraits operator()(DerivativeTraits *source, DerivativeTraits *probe) const

Evaluates the Green’s function given a pair of points

Parameters
  • source: the source point
  • probe: the probe point

double kernelD_impl(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the kernel of the \(\mathcal{D}\) integral operator for the pair of points p1, p2: \( [\boldsymbol{\varepsilon}\nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)]\cdot \mathbf{n}_{\mathbf{p}_2}\) To obtain the kernel of the \(\mathcal{D}^\dagger\) operator call this methods with \(\mathbf{p}_1\) and \(\mathbf{p}_2\) exchanged and with \(\mathbf{n}_{\mathbf{p}_2} = \mathbf{n}_{\mathbf{p}_1}\)

Parameters
  • direction: the direction
  • p1: first point
  • p2: second point

double singleLayer_impl(const Element &e, double factor) const

Methods to compute the diagonal of the matrix representation of the S and D operators by approximate collocation.

Calculates an element on the diagonal of the matrix representation of the S operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

double doubleLayer_impl(const Element &e, double factor) const

Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

SphericalDiffuse

template <typename ProfilePolicy = OneLayerTanh>
class SphericalDiffuse

Green’s function for a diffuse interface with spherical symmetry.

This class is general, in the sense that no specific dielectric profile has been set in its definition. In principle any profile that can be described by:

  1. a left-side dielectric constant;
  2. a right-side dielectric constant;
  3. an interface layer width;
  4. an interface layer center can be used to define a new diffuse interface with spherical symmetry. The origin of the dielectric sphere can be changed by means of the constructor. The solution of the differential equation defining the Green’s function is always** performed assuming that the dielectric sphere is centered in the origin of the coordinate system. Whenever the public methods are invoked to “sample” the Green’s function at a pair of points, a translation of the sampling points is performed first.
Author
Hui Cao, Ville Weijo, Luca Frediani and Roberto Di Remigio
Date
2010-2015
Template Parameters
  • ProfilePolicy: functional form of the diffuse layer

Inherits from GreensFunction< Numerical, ProfilePolicy >

Unnamed Group

int maxLGreen_

Parameters and functions for the calculation of the Green’s function, including Coulomb singularity

Maximum angular momentum in the __final summation over Legendre polynomials to obtain G

std::vector<RadialFunction<interfaces::StateType, interfaces::LnTransformedRadial, Zeta>> zeta_

First independent radial solution, used to build Green’s function.

Note
The vector has dimension maxLGreen_ and has r^l behavior

std::vector<RadialFunction<interfaces::StateType, interfaces::LnTransformedRadial, Omega>> omega_

Second independent radial solution, used to build Green’s function.

Note
The vector has dimension maxLGreen_ and has r^(-l-1) behavior

double imagePotentialComponent_impl(int L, const Eigen::Vector3d &sp, const Eigen::Vector3d &pp, double Cr12) const

Returns L-th component of the radial part of the Green’s function.

Note
This function shifts the given source and probe points by the location of the dielectric sphere.
Parameters
  • L: angular momentum
  • sp: source point
  • pp: probe point
  • Cr12: Coulomb singularity separation coefficient

Unnamed Group

int maxLC_

Parameters and functions for the calculation of the Coulomb singularity separation coefficient

Maximum angular momentum to obtain C(r, r’), needed to separate the Coulomb singularity

RadialFunction<interfaces::StateType, interfaces::LnTransformedRadial, Zeta> zetaC_

First independent radial solution, used to build coefficient.

Note
This is needed to separate the Coulomb singularity and has r^l behavior

RadialFunction<interfaces::StateType, interfaces::LnTransformedRadial, Omega> omegaC_

Second independent radial solution, used to build coefficient.

Note
This is needed to separate the Coulomb singularity and has r^(-l-1) behavior

double coefficient_impl(const Eigen::Vector3d &sp, const Eigen::Vector3d &pp) const

Returns coefficient for the separation of the Coulomb singularity.

Note
This function shifts the given source and probe points by the location of the dielectric sphere.
Parameters
  • sp: first point
  • pp: second point

Public Functions

SphericalDiffuse(double e1, double e2, double w, double c, const Eigen::Vector3d &o, int l)

Constructor for a one-layer interface PCMSolver, an API for the Polarizable Continuum Model Copyright (C) 2016 Roberto Di Remigio, Luca Frediani and collaborators.

Parameters
  • e1: left-side dielectric constant
  • e2: right-side dielectric constant
  • w: width of the interface layer
  • c: center of the diffuse layer
  • o: center of the sphere
  • l: maximum value of angular momentum

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/

double coefficientCoulomb(const Eigen::Vector3d &source, const Eigen::Vector3d &probe) const

Returns Coulomb singularity separation coefficient.

Parameters
  • source: location of the source charge
  • probe: location of the probe charge

double Coulomb(const Eigen::Vector3d &source, const Eigen::Vector3d &probe) const

Returns singular part of the Green’s function.

Parameters
  • source: location of the source charge
  • probe: location of the probe charge

double imagePotential(const Eigen::Vector3d &source, const Eigen::Vector3d &probe) const

Returns non-singular part of the Green’s function (image potential)

Parameters
  • source: location of the source charge
  • probe: location of the probe charge

double coefficientCoulombDerivative(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the directional derivative of the Coulomb singularity separation coefficient for the pair of points p1, p2: \( \nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)\cdot *\mathbf{n}_{\mathbf{p}_2}\) Notice that this method returns the directional derivative with respect to the probe point, thus assuming that the direction is relative to that point.

Parameters
  • direction: the direction
  • p1: first point
  • p2: second point

double CoulombDerivative(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the directional derivative of the singular part of the Greens’s function for the pair of points p1, p2: \( \nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)\cdot *\mathbf{n}_{\mathbf{p}_2}\) Notice that this method returns the directional derivative with respect to the probe point, thus assuming that the direction is relative to that point.

Parameters
  • direction: the direction
  • p1: first point
  • p2: second point

double imagePotentialDerivative(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the directional derivative of the non-singular part (image potential) of the Greens’s function for the pair of points p1, p2: \( \nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)\cdot *\mathbf{n}_{\mathbf{p}_2}\) Notice that this method returns the directional derivative with respect to the probe point, thus assuming that the direction is relative to that point.

Parameters
  • direction: the direction
  • p1: first point
  • p2: second point

pcm::tuple<double, double> epsilon(const Eigen::Vector3d &point) const

Handle to the dielectric profile evaluation

Private Functions

Numerical operator()(Numerical *sp, Numerical *pp) const

Evaluates the Green’s function given a pair of points

Note
This takes care of the origin shift
Parameters
  • sp: the source point
  • pp: the probe point

double kernelD_impl(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const

Returns value of the kernel of the \(\mathcal{D}\) integral operator for the pair of points p1, p2: \( [\boldsymbol{\varepsilon}\nabla_{\mathbf{p_2}}G(\mathbf{p}_1, \mathbf{p}_2)]\cdot \mathbf{n}_{\mathbf{p}_2}\) To obtain the kernel of the \(\mathcal{D}^\dagger\) operator call this methods with \(\mathbf{p}_1\) and \(\mathbf{p}_2\) exchanged and with \(\mathbf{n}_{\mathbf{p}_2} = \mathbf{n}_{\mathbf{p}_1}\)

Parameters
  • direction: the direction
  • p1: first point
  • p2: second point

double singleLayer_impl(const Element &e, double factor) const

Methods to compute the diagonal of the matrix representation of the S and D operators by approximate collocation.

Calculates an element on the diagonal of the matrix representation of the S operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

double doubleLayer_impl(const Element &e, double factor) const

Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.

Parameters
  • e: finite element on the cavity
  • factor: the scaling factor for the diagonal elements

void initProfilePolicy(double e1, double e2, double w, double c)

Initializes a one-layer profile

Parameters
  • e1: left-side dielectric constant
  • e2: right-side dielectric constant
  • w: width of the interface layer
  • c: center of the diffuse layer

void initSphericalDiffuse()

This calculates all the components needed to evaluate the Green’s function

Private Members

Eigen::Vector3d origin_

Center of the dielectric sphere