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.
The top-level header, _i.e._ to be included in client code, is Green.hpp
.
The common interface to all Green’s function classes is specified by the IGreensFunction
class,
this is non-templated.
All other classes are templated.
The Green’s functions are registered to the factory based on a label encoding: type, derivative, and dielectric profile.
The only allowed labels must be listed in src/green/Green.hpp
. If they are not, they can not be selected at run time.
IGreensFunction¶
-
class
pcm
::
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 pcm::green::GreensFunction< DerivativeTraits, dielectric_profile::Anisotropic >, pcm::green::GreensFunction< DerivativeTraits, dielectric_profile::Sharp >, pcm::green::GreensFunction< DerivativeTraits, dielectric_profile::Uniform >, pcm::green::GreensFunction< DerivativeTraits, dielectric_profile::Yukawa >, pcm::green::GreensFunction< DerivativeTraits, ProfilePolicy >, pcm::green::GreensFunction< Stencil, 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
[in] p1
: first point[in] 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
[in] direction
: the direction[in] p1
: first point[in] 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
[in] e
: finite element on the cavity[in] 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
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
Unnamed Group
-
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
[in] p1
: first point[in] p2
: second point
-
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
[in] direction
: the direction[in] p1
: first point[in] p2
: second point
Unnamed Group
-
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
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
-
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
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
GreensFunction¶
-
template<typename
DerivativeTraits
, typenameProfilePolicy
>
classpcm::green
::
GreensFunction
: public pcm::IGreensFunction¶ 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
Unnamed Group
-
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
[in] normal_p1
: the normal vector to p1[in] p1
: first point[in] p2
: second point
-
double
derivativeProbe
(const Eigen::Vector3d &normal_p2, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const final¶ 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
[in] normal_p2
: the normal vector to p2[in] p1
: first point[in] 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
[in] p1
: first point[in] 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
[in] p1
: first point[in] p2
: second point
Public Functions
-
bool
uniform
() const final override¶ Whether the Green’s function describes a uniform environment
Protected Functions
-
DerivativeTraits
operator()
(DerivativeTraits *source, DerivativeTraits *probe) const = 0¶ Evaluates the Green’s function given a pair of points
- Parameters
[in] source
: the source point[in] probe
: the probe point
-
double
kernelS_impl
(const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const final override¶ 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
[in] p1
: first point[in] p2
: second point
Protected Attributes
-
double
delta_
¶ Step for numerical differentiation.
-
ProfilePolicy
profile_
¶ Permittivity profile.
Vacuum¶
-
template<typename
DerivativeTraits
= AD_directional>
classpcm::green
::
Vacuum
: public pcm::green::GreensFunction<DerivativeTraits, dielectric_profile::Uniform>¶ 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
Public Functions
-
double
permittivity
() const final override¶ Returns a dielectric permittivity
Private Functions
-
DerivativeTraits
operator()
(DerivativeTraits *sp, DerivativeTraits *pp) const override¶ Evaluates the Green’s function given a pair of points
- Parameters
[in] source
: the source point[in] probe
: the probe point
-
double
kernelD_impl
(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const override¶ 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
[in] direction
: the direction[in] p1
: first point[in] p2
: second point
-
double
singleLayer_impl
(const Element &e, double factor) const override¶ 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
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
-
double
doubleLayer_impl
(const Element &e, double factor) const override¶ Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.
- Parameters
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
UniformDielectric¶
-
template<typename
DerivativeTraits
= AD_directional>
classpcm::green
::
UniformDielectric
: public pcm::green::GreensFunction<DerivativeTraits, dielectric_profile::Uniform>¶ 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
Public Functions
-
double
permittivity
() const final override¶ Returns a dielectric permittivity
Private Functions
-
DerivativeTraits
operator()
(DerivativeTraits *sp, DerivativeTraits *pp) const override¶ Evaluates the Green’s function given a pair of points
- Parameters
[in] source
: the source point[in] probe
: the probe point
-
double
kernelD_impl
(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const override¶ 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
[in] direction
: the direction[in] p1
: first point[in] p2
: second point
-
double
singleLayer_impl
(const Element &e, double factor) const override¶ 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
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
-
double
doubleLayer_impl
(const Element &e, double factor) const override¶ Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.
- Parameters
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
IonicLiquid¶
-
template<typename
DerivativeTraits
= AD_directional>
classpcm::green
::
IonicLiquid
: public pcm::green::GreensFunction<DerivativeTraits, dielectric_profile::Yukawa>¶ 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
Public Functions
-
double
permittivity
() const final override¶ Returns a dielectric permittivity
Private Functions
-
DerivativeTraits
operator()
(DerivativeTraits *sp, DerivativeTraits *pp) const override¶ Evaluates the Green’s function given a pair of points
- Parameters
[in] source
: the source point[in] probe
: the probe point
-
double
kernelD_impl
(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const override¶ 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
[in] direction
: the direction[in] p1
: first point[in] p2
: second point
-
double
singleLayer_impl
(const Element&, double) const override¶ 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
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
-
double
doubleLayer_impl
(const Element&, double) const override¶ Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.
- Parameters
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
AnisotropicLiquid¶
-
template<typename
DerivativeTraits
= AD_directional>
classpcm::green
::
AnisotropicLiquid
: public pcm::green::GreensFunction<DerivativeTraits, dielectric_profile::Anisotropic>¶ 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
Public Functions
-
AnisotropicLiquid
(const Eigen::Vector3d &eigen_eps, const Eigen::Vector3d &euler_ang)¶ - Parameters
[in] eigen_eps
: eigenvalues of the permittivity tensors[in] euler_ang
: Euler angles in degrees
-
double
permittivity
() const final override¶ Returns a dielectric permittivity
Private Functions
-
DerivativeTraits
operator()
(DerivativeTraits *sp, DerivativeTraits *pp) const override¶ Evaluates the Green’s function given a pair of points
- Parameters
[in] source
: the source point[in] probe
: the probe point
-
double
kernelD_impl
(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const override¶ 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
[in] direction
: the direction[in] p1
: first point[in] p2
: second point
-
double
singleLayer_impl
(const Element&, double) const override¶ 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
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
-
double
doubleLayer_impl
(const Element&, double) const override¶ Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.
- Parameters
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
SphericalDiffuse¶
-
template<typename
ProfilePolicy
= dielectric_profile::OneLayerLog>
classpcm::green
::
SphericalDiffuse
: public pcm::green::GreensFunction<Stencil, ProfilePolicy>¶ Green’s function for a 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
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<detail::StateType, detail::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<detail::StateType, detail::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
[in] L
: angular momentum[in] sp
: source point[in] pp
: probe point[in] 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<detail::StateType, detail::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<detail::StateType, detail::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
[in] sp
: first point[in] 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
- Parameters
[in] e1
: left-side dielectric constant[in] e2
: right-side dielectric constant[in] w
: width of the interface layer[in] c
: center of the diffuse layer[in] o
: center of the sphere[in] l
: maximum value of angular momentum
-
double
permittivity
() const final override¶ Returns a dielectric permittivity
-
double
coefficientCoulomb
(const Eigen::Vector3d &source, const Eigen::Vector3d &probe) const¶ Returns Coulomb singularity separation coefficient.
- Parameters
[in] source
: location of the source charge[in] 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
[in] source
: location of the source charge[in] 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
[in] source
: location of the source charge[in] 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
[in] direction
: the direction[in] p1
: first point[in] 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
[in] direction
: the direction[in] p1
: first point[in] 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
[in] direction
: the direction[in] p1
: first point[in] p2
: second point
-
std::tuple<double, double>
epsilon
(const Eigen::Vector3d &point) const¶ Handle to the dielectric profile evaluation
Private Functions
-
Stencil
operator()
(Stencil *sp, Stencil *pp) const override¶ Evaluates the Green’s function given a pair of points
- Note
This takes care of the origin shift
- Parameters
[in] sp
: the source point[in] pp
: the probe point
-
double
kernelD_impl
(const Eigen::Vector3d &direction, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2) const override¶ 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
[in] direction
: the direction[in] p1
: first point[in] p2
: second point
-
double
singleLayer_impl
(const Element &e, double factor) const override¶ 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
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
-
double
doubleLayer_impl
(const Element &e, double factor) const override¶ Calculates an element of the diagonal of the matrix representation of the D operator using an approximate collocation formula.
- Parameters
[in] e
: finite element on the cavity[in] factor
: the scaling factor for the diagonal elements
-
void
initSphericalDiffuse
()¶ This calculates all the components needed to evaluate the Green’s function
Private Members
-
Eigen::Vector3d
origin_
¶ Center of the dielectric sphere