BEM++  2.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Groups Pages
Public Types | Public Member Functions | Private Types | Private Member Functions | List of all members
Bempp::ElementaryPotentialOperator< BasisFunctionType_, KernelType_, ResultType_ > Class Template Referenceabstract

Elementary potential operator. More...

#include </home/wojtek/Projects/BEM/bempp-sven/bempp/lib/assembly/elementary_potential_operator.hpp>

Inheritance diagram for Bempp::ElementaryPotentialOperator< BasisFunctionType_, KernelType_, ResultType_ >:
Bempp::PotentialOperator< BasisFunctionType_, ResultType_ >

Public Types

typedef Base::BasisFunctionType BasisFunctionType
 Type of the values of the (components of the) basis functions into which functions acted upon by the operator are expanded. More...
 
typedef KernelType_ KernelType
 Type of the values of the (components of the) kernel functions.
 
typedef Base::ResultType ResultType
 Type of the values of the (components of the) potential. More...
 
typedef Base::CoordinateType CoordinateType
 Type used to represent coordinates. More...
 
typedef Base::QuadratureStrategy QuadratureStrategy
 Type of the appropriate instantiation of Fiber::QuadratureStrategy. More...
 
typedef
Fiber::EvaluatorForIntegralOperators
< ResultType
Evaluator
 Type of the appropriate instantiation of Fiber::EvaluatorForIntegralOperators.
 
typedef
Fiber::LocalAssemblerForPotentialOperators
< ResultType
LocalAssembler
 Type of the appropriate instantiation of Fiber::LocalAssemblerForPotentialOperators.
 
typedef
Fiber::CollectionOfShapesetTransformations
< CoordinateType
CollectionOfShapesetTransformations
 Type of the appropriate instantiation of Fiber::CollectionOfShapesetTransformations.
 
typedef
Fiber::CollectionOfShapesetTransformations
< CoordinateType
CollectionOfBasisTransformations
 Type of the appropriate instantiation of Fiber::CollectionOfBasisTransformations. More...
 
typedef
Fiber::CollectionOfKernels
< KernelType
CollectionOfKernels
 Type of the appropriate instantiation of Fiber::CollectionOfKernels.
 
typedef
Fiber::KernelTrialIntegral
< BasisFunctionType,
KernelType, ResultType
KernelTrialIntegral
 Type of the appropriate instantiation of Fiber::KernelTrialIntegral.
 
- Public Types inherited from Bempp::PotentialOperator< BasisFunctionType_, ResultType_ >
typedef BasisFunctionType_ BasisFunctionType
 Type of the values of the (components of the) basis functions into which functions acted upon by the operator are expanded.
 
typedef ResultType_ ResultType
 Type of the values of the (components of the) potential.
 
typedef ScalarTraits
< ResultType >::RealType 
CoordinateType
 Type used to represent coordinates.
 
typedef
Fiber::QuadratureStrategy
< BasisFunctionType,
ResultType, GeometryFactory
QuadratureStrategy
 Type of the appropriate instantiation of Fiber::QuadratureStrategy.
 

Public Member Functions

virtual std::auto_ptr
< InterpolatedFunction
< ResultType_ > > 
evaluateOnGrid (const GridFunction< BasisFunctionType, ResultType > &argument, const Grid &evaluationGrid, const QuadratureStrategy &quadStrategy, const EvaluationOptions &options) const
 Evaluate the potential of a given charge distribution on a prescribed grid. More...
 
virtual arma::Mat< ResultType_ > evaluateAtPoints (const GridFunction< BasisFunctionType, ResultType > &argument, const arma::Mat< CoordinateType > &evaluationPoints, const QuadratureStrategy &quadStrategy, const EvaluationOptions &options) const
 Evaluate the potential of a given charge distribution at prescribed points. More...
 
virtual
AssembledPotentialOperator
< BasisFunctionType_,
ResultType_ > 
assemble (const shared_ptr< const Space< BasisFunctionType > > &space, const shared_ptr< const arma::Mat< CoordinateType > > &evaluationPoints, const QuadratureStrategy &quadStrategy, const EvaluationOptions &options) const
 Create and return an AssembledPotentialOperator object. More...
 
virtual int componentCount () const
 Number of components of the values of the potential. More...
 
- Public Member Functions inherited from Bempp::PotentialOperator< BasisFunctionType_, ResultType_ >
virtual ~PotentialOperator ()
 Destructor.
 

Private Types

typedef PotentialOperator
< BasisFunctionType_,
ResultType_ > 
Base
 

Private Member Functions

virtual const CollectionOfKernelskernels () const =0
 Return the collection of kernel functions occurring in the integrand of this operator.
 
virtual const
CollectionOfBasisTransformations
trialTransformations () const =0
 Return the collection of transformations of the charge distribution that occur in the weak form of this operator.
 
virtual const KernelTrialIntegralintegral () const =0
 Return an object representing the integral used to evaluate the potential of a charge distribution. More...
 

Detailed Description

template<typename BasisFunctionType_, typename KernelType_, typename ResultType_>
class Bempp::ElementaryPotentialOperator< BasisFunctionType_, KernelType_, ResultType_ >

Elementary potential operator.

This class provides the interface for evaluation of a potential $k(x)$ defined by the formula

\[ k(x) = \int_\Gamma F[x, \psi(y)] \, \mathrm{d}\Gamma, \]

where the integration goes over a surface $\Gamma$ and the integrand $F$ depends on the coordinates of a point $x$ lying outside $\Gamma$ and a surface charge distribution $\psi(y)$ does not lie on $\Gamma$. The function $F[x, \psi(y)]$ is assumed to be linear in $psi(y)$. In the simplest and most common case, $F[x, \psi(y)]$ is just

\[ F[x, \psi(y)] = K(x, y) \psi(y), \]

where $K(x, y)$ is a kernel function. For more complex operators, $F$ might involve some transformations of the charge distribution (e.g. its surface divergence or curl), the kernel function might be a tensor, or $F$ might consist of several terms. The form of $F$ for a particular potential operator is determined by the implementation of the virtual functions integral(), kernels() and trialTransformations().

Template Parameters
BasisFunctionType_Type of the values of the (components of the) basis functions into which functions acted upon by the operator are expanded.
KernelType_Type of the values of the (components of the) kernel functions occurring in the integrand of the operator.
ResultType_Type of the values of the (components of the) potential.

All three template parameters can take the following values: float, double, std::complex<float> and std::complex<double>. All types must have the same precision: for instance, mixing float with std::complex<double> is not allowed. If either BasisFunctionType_ or KernelType_ is a complex type, then ResultType_ must be set to the same type.

Member Typedef Documentation

template<typename BasisFunctionType_, typename KernelType_, typename ResultType_>
typedef Base::BasisFunctionType Bempp::ElementaryPotentialOperator< BasisFunctionType_, KernelType_, ResultType_ >::BasisFunctionType

Type of the values of the (components of the) basis functions into which functions acted upon by the operator are expanded.

template<typename BasisFunctionType_, typename KernelType_, typename ResultType_>
typedef Fiber::CollectionOfShapesetTransformations<CoordinateType> Bempp::ElementaryPotentialOperator< BasisFunctionType_, KernelType_, ResultType_ >::CollectionOfBasisTransformations

Type of the appropriate instantiation of Fiber::CollectionOfBasisTransformations.

Deprecated:
This type is deprecated; use CollectionOfShapesetTransformations instead.
template<typename BasisFunctionType_, typename KernelType_, typename ResultType_>
typedef Base::CoordinateType Bempp::ElementaryPotentialOperator< BasisFunctionType_, KernelType_, ResultType_ >::CoordinateType

Type used to represent coordinates.

template<typename BasisFunctionType_, typename KernelType_, typename ResultType_>
typedef Base::QuadratureStrategy Bempp::ElementaryPotentialOperator< BasisFunctionType_, KernelType_, ResultType_ >::QuadratureStrategy

Type of the appropriate instantiation of Fiber::QuadratureStrategy.

template<typename BasisFunctionType_, typename KernelType_, typename ResultType_>
typedef Base::ResultType Bempp::ElementaryPotentialOperator< BasisFunctionType_, KernelType_, ResultType_ >::ResultType

Type of the values of the (components of the) potential.

Member Function Documentation

template<typename BasisFunctionType , typename KernelType , typename ResultType >
AssembledPotentialOperator< BasisFunctionType, ResultType > Bempp::ElementaryPotentialOperator< BasisFunctionType, KernelType, ResultType >::assemble ( const shared_ptr< const Space< BasisFunctionType > > &  space,
const shared_ptr< const arma::Mat< CoordinateType > > &  evaluationPoints,
const QuadratureStrategy quadStrategy,
const EvaluationOptions options 
) const
virtual

Create and return an AssembledPotentialOperator object.

The returned AssembledPotentialOperator object stores the values of the potentials generated at the points listed in the array evaluationPoints by the charge distributions equal to the individual basis functions of the space space. The object can afterwards be used to evaluate efficiently the potentials generated by multiple GridFunctions expanded in the space space.

Parameters
[in]spaceThe space whose basis functions will be taken as the charge distributions inducing the potentials to be evaluated.
[in]evaluationPoints2D array whose (i, j)th element is the ith coordinate of the jth point at which the potential should be evaluated. The first dimension of this array should be equal to space.grid().dimWorld().
[in]quadStrategyA QuadratureStrategy object controlling how the integrals will be evaluated.
[in]optionsEvaluation options. This parameter controls, notably, the format used to store the matrix of precalculated potential values: a dense matrix or an H-matrix.

Implements Bempp::PotentialOperator< BasisFunctionType_, ResultType_ >.

template<typename BasisFunctionType , typename KernelType , typename ResultType >
int Bempp::ElementaryPotentialOperator< BasisFunctionType, KernelType, ResultType >::componentCount ( ) const
virtual

Number of components of the values of the potential.

E.g. 1 for a scalar-valued potential, 3 for a vector-valued potential.

Implements Bempp::PotentialOperator< BasisFunctionType_, ResultType_ >.

template<typename BasisFunctionType , typename KernelType , typename ResultType >
arma::Mat< ResultType > Bempp::ElementaryPotentialOperator< BasisFunctionType, KernelType, ResultType >::evaluateAtPoints ( const GridFunction< BasisFunctionType, ResultType > &  argument,
const arma::Mat< CoordinateType > &  evaluationPoints,
const QuadratureStrategy quadStrategy,
const EvaluationOptions options 
) const
virtual

Evaluate the potential of a given charge distribution at prescribed points.

Parameters
[in]argumentArgument of the potential operator ( $\psi(y)$ in the notation above), represented by a grid function.
[in]evaluationPoints2D array whose (i, j)th element is the ith coordinate of the jth point at which the potential should be evaluated. The first dimension of this array should be equal to argument.grid().dimWorld().
[in]quadStrategyA QuadratureStrategy object controlling how the integrals will be evaluated.
[in]optionsEvaluation options.
Returns
A 2D array whose (i, j)th element is the ith component of the potential at the jth point.
Note
This function is not designed to yield accurate values of the potential on the surface $\Gamma$ containing the charge distribution, i.e. argument.grid(), even if the potential has a unique extension from $\Omega \setminus \Gamma$ to $\Gamma$. Hence values of the potential at any points belonging to $\Gamma$ can be badly wrong.

The current implementation does not yet take special measures to prevent loss of accuracy near $\Gamma$, either. Users are advised to increase the quadrature accuracy for points lying in the vicinity of $\Gamma$.

Implements Bempp::PotentialOperator< BasisFunctionType_, ResultType_ >.

References Bempp::EvaluationOptions::ACA, Bempp::AssembledPotentialOperator< BasisFunctionType, ResultType >::apply(), Bempp::EvaluationOptions::DENSE, Bempp::EvaluationOptions::evaluationMode(), Bempp::GridFunction< BasisFunctionType, ResultType >::grid(), and Bempp::GridFunction< BasisFunctionType, ResultType >::space().

template<typename BasisFunctionType , typename KernelType , typename ResultType >
std::auto_ptr< InterpolatedFunction< ResultType > > Bempp::ElementaryPotentialOperator< BasisFunctionType, KernelType, ResultType >::evaluateOnGrid ( const GridFunction< BasisFunctionType, ResultType > &  argument,
const Grid evaluationGrid,
const QuadratureStrategy quadStrategy,
const EvaluationOptions options 
) const
virtual

Evaluate the potential of a given charge distribution on a prescribed grid.

Parameters
[in]argumentArgument of the potential operator ( $\psi(y)$ in the notation above), represented by a grid function.
[in]evaluationGridGrid at whose vertices the potential will be evaluated. The grid may have arbitrary dimension, but must be embedded in a world of the same dimension as argument.grid().
[in]quadStrategyA QuadratureStrategy object controlling how the integrals will be evaluated.
[in]optionsEvaluation options.
Returns
The potential represented by a function interpolated on the vertices of evaluationGrid.
Note
This function is not designed to yield accurate values of the potential on the surface $\Gamma$ containing the charge distribution, i.e. argument.grid(), even if the potential has a unique extension from $\Omega \setminus \Gamma$ to $\Gamma$. Hence values of the potential at any vertices of evaluationGrid that coincide with $\Gamma$ can be badly wrong.

The current implementation does not yet take special measures to prevent loss of accuracy near $\Gamma$, either. If in doubt, increase the quadrature accuracy.

Implements Bempp::PotentialOperator< BasisFunctionType_, ResultType_ >.

References Bempp::Grid::dim(), Bempp::Grid::dimWorld(), Bempp::IndexSet::entityIndex(), Bempp::Entity< codim >::geometry(), Bempp::Geometry::getCenter(), Bempp::GridFunction< BasisFunctionType, ResultType >::grid(), and Bempp::Grid::leafView().

template<typename BasisFunctionType_, typename KernelType_, typename ResultType_>
virtual const KernelTrialIntegral& Bempp::ElementaryPotentialOperator< BasisFunctionType_, KernelType_, ResultType_ >::integral ( ) const
privatepure virtual

Return an object representing the integral used to evaluate the potential of a charge distribution.

Subclasses of KernelTrialIntegral implement functions that evaluate the integral using the data provided by a CollectionOfKernels representing the kernel functions occurring in the integrand and a CollectionOfBasisTransformations representing the charge-distribution transformations occurring in the integrand.

Implemented in Bempp::Laplace3dPotentialOperatorBase< Impl, BasisFunctionType_, ResultType_ >, Bempp::Laplace3dPotentialOperatorBase< Laplace3dDoubleLayerPotentialOperatorImpl< BasisFunctionType_, ResultType_ >, BasisFunctionType_, ResultType_ >, Bempp::Laplace3dPotentialOperatorBase< Laplace3dSingleLayerPotentialOperatorImpl< BasisFunctionType_, ResultType_ >, BasisFunctionType_, ResultType_ >, Bempp::ModifiedHelmholtz3dPotentialOperatorBase< Impl, BasisFunctionType_ >, Bempp::ModifiedHelmholtz3dPotentialOperatorBase< ModifiedHelmholtz3dDoubleLayerPotentialOperatorImpl< BasisFunctionType_ >, BasisFunctionType_ >, Bempp::ModifiedHelmholtz3dPotentialOperatorBase< ModifiedHelmholtz3dSingleLayerPotentialOperatorImpl< BasisFunctionType_ >, BasisFunctionType_ >, Bempp::Helmholtz3dPotentialOperatorBase< Impl, BasisFunctionType_ >, Bempp::Helmholtz3dPotentialOperatorBase< Helmholtz3dFarFieldDoubleLayerPotentialOperatorImpl< BasisFunctionType_ >, BasisFunctionType_ >, Bempp::Helmholtz3dPotentialOperatorBase< Helmholtz3dFarFieldSingleLayerPotentialOperatorImpl< BasisFunctionType_ >, BasisFunctionType_ >, Bempp::Helmholtz3dPotentialOperatorBase< Maxwell3dSingleLayerPotentialOperatorImpl< BasisFunctionType_ >, BasisFunctionType_ >, Bempp::Helmholtz3dPotentialOperatorBase< Maxwell3dDoubleLayerPotentialOperatorImpl< BasisFunctionType_ >, BasisFunctionType_ >, Bempp::Helmholtz3dPotentialOperatorBase< Helmholtz3dSingleLayerPotentialOperatorImpl< BasisFunctionType_ >, BasisFunctionType_ >, Bempp::Helmholtz3dPotentialOperatorBase< Helmholtz3dDoubleLayerPotentialOperatorImpl< BasisFunctionType_ >, BasisFunctionType_ >, Bempp::Helmholtz3dPotentialOperatorBase< Maxwell3dFarFieldSingleLayerPotentialOperatorImpl< BasisFunctionType_ >, BasisFunctionType_ >, and Bempp::Helmholtz3dPotentialOperatorBase< Maxwell3dFarFieldDoubleLayerPotentialOperatorImpl< BasisFunctionType_ >, BasisFunctionType_ >.


The documentation for this class was generated from the following files: