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.
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 pointp2
: 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 directionp1
: first pointp2
: 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 cavityfactor
: 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 cavityfactor
: 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 pointp2
: 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 directionp1
: first pointp2
: 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 cavityfactor
: 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 cavityfactor
: the scaling factor for the diagonal elements
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 derivativesProfilePolicy
: 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 p1p1
: first pointp2
: 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 p2p1
: first pointp2
: 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 pointp2
: 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 pointp2
: 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 pointprobe
: 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 pointp2
: second point
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 pointprobe
: 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 directionp1
: first pointp2
: 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 cavityfactor
: 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 cavityfactor
: 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 pointprobe
: 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 directionp1
: first pointp2
: 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 cavityfactor
: 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 cavityfactor
: 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 pointprobe
: 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 directionp1
: first pointp2
: 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 cavityfactor
: 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 cavityfactor
: 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 tensorseuler_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 pointprobe
: 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 directionp1
: first pointp2
: 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 cavityfactor
: 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 cavityfactor
: 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:
- a left-side dielectric constant;
- a right-side dielectric constant;
- an interface layer width;
- 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 momentumsp
: source pointpp
: probe pointCr12
: 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 pointpp
: 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 constante2
: right-side dielectric constantw
: width of the interface layerc
: center of the diffuse layero
: center of the spherel
: 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 chargeprobe
: 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 chargeprobe
: 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 chargeprobe
: 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 directionp1
: first pointp2
: 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 directionp1
: first pointp2
: 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 directionp1
: first pointp2
: second point
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 pointpp
: 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 directionp1
: first pointp2
: 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 cavityfactor
: 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 cavityfactor
: 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 constante2
: right-side dielectric constantw
: width of the interface layerc
: 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