BEM++  2.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Groups Pages
Public Types | Public Member Functions | List of all members
Bempp::EvaluationOptions Class Reference

Options controlling evaluation of potentials. More...

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

Public Types

enum  { AUTO = -1 }
 

Public Member Functions

 EvaluationOptions ()
 Constructor.
 
void setMaxThreadCount (int maxThreadCount)
 Set the maximum number of threads used during evaluation of potentials. More...
 
BEMPP_DEPRECATED void switchToTbb (int maxThreadCount=AUTO)
 Set the maximum number of threads used during evaluation of potentials. More...
 
const ParallelizationOptionsparallelizationOptions () const
 Return current parallelization options.
 
void setVerbosityLevel (VerbosityLevel::Level level)
 Set the verbosity level. More...
 
VerbosityLevel::Level verbosityLevel () const
 Return the verbosity level.
 

Evaluation mode

enum  Mode { DENSE, ACA }
 Possible evaluation modes. More...
 
void switchToDenseMode ()
 Use dense-matrix representations of elementary potential operators. More...
 
void switchToAcaMode (const AcaOptions &acaOptions)
 Use adaptive cross approximation (ACA) to obtain hierarchical-matrix representations of potential operators. More...
 
Mode evaluationMode () const
 Return current evaluation mode. More...
 
const AcaOptionsacaOptions () const
 Return the current adaptive cross approximation (ACA) settings. More...
 

Detailed Description

Options controlling evaluation of potentials.

Member Enumeration Documentation

Possible evaluation modes.

Enumerator
DENSE 

Assemble dense matrices.

ACA 

Assemble hierarchical matrices using adaptive cross approximation (ACA).

Member Function Documentation

const AcaOptions & Bempp::EvaluationOptions::acaOptions ( ) const

Return the current adaptive cross approximation (ACA) settings.

Note
These settings are only used in the ACA evaluation mode, i.e. when evaluationMode() returns ACA.

Referenced by switchToAcaMode().

EvaluationOptions::Mode Bempp::EvaluationOptions::evaluationMode ( ) const

Return current evaluation mode.

The evaluation mode can be changed by calling switchToDenseMode() or switchToAcaMode().

Referenced by Bempp::ElementaryPotentialOperator< BasisFunctionType_, KernelType_, ResultType_ >::evaluateAtPoints().

void Bempp::EvaluationOptions::setMaxThreadCount ( int  maxThreadCount)

Set the maximum number of threads used during evaluation of potentials.

maxThreadCount must be a positive number or AUTO. In the latter case the number of threads is determined automatically.

Referenced by switchToTbb().

void Bempp::EvaluationOptions::setVerbosityLevel ( VerbosityLevel::Level  level)

Set the verbosity level.

This setting determines the amount of information printed out by functions from BEM++.

void Bempp::EvaluationOptions::switchToAcaMode ( const AcaOptions acaOptions)

Use adaptive cross approximation (ACA) to obtain hierarchical-matrix representations of potential operators.

Parameters
[in]acaOptionsParameters influencing the ACA algorithm.

In this mode, evaluation of potentials always entails the construction of a hierarchical-matrix representation of the potential operator, as defined in the documentation of switchToDenseMode(). The potential given by a particular charge distribution is then calculated by left-multiplying the vector of its expansion coefficients by the above hierarchical matrix.

Note
Our (BEM++ developers') numerical experiments indicate that to maximize the efficiency of potential evaluation, the acaOptions.eta parameter should be chosen much smaller than the default 1.2 (which is adequate for the ACA of discrete weak forms). For potential evaluation, eta = 0.4 seems to work well.

References ACA, and acaOptions().

void Bempp::EvaluationOptions::switchToDenseMode ( )

Use dense-matrix representations of elementary potential operators.

This is the default evaluation mode. If it is active, potentials are evaluated in the following way.

(1) When the potential needs to be evaluated only for a single charge distribution $\psi(y)$, as is the case in the PotentialOperator::evaluateAtPoints() and evaluateOnGrid() functions, the defining formula of the potential

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

is approximated with a numerical quadrature rule

\[ k(x) \approx \sum_{j=1}^N w_j F[x, \psi(y_j)], \]

where $y_j \in \Gamma$ and $w_j$ are the quadrature points and weights. The values of $\psi(y)$ at the quadrature points $y_j$ are precalculated and stored, and then the potential is sequentially evaluated at each requested point $x_i$ using the formula above. The values of any kernels involved in $F$ are evaluated once for each pair $(x_i, y_j)$ and then discarded.

(2) If the potentials due to multiple charge distributions need to be evaluated at a fixed set of points, the PotentialOperator::assemble() function can be used to generate an AssembledPotentialOperator object storing a (dense) matrix representation of the potential operator, with the element (c * i, j) containing the value of cth component of the potential produced at ith point by the charge distribution equal to jth basis function of a certain function space. The potential generated by a specific charge distribution expanded in this space is then calculated by left-multiplying the vector of its expansion coefficients by the matrix representation of the potential operator.

References DENSE.

void Bempp::EvaluationOptions::switchToTbb ( int  maxThreadCount = AUTO)

Set the maximum number of threads used during evaluation of potentials.

Deprecated:
Use setMaxThreadCount() instead.

References setMaxThreadCount().


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