HyTeG
Matrix-Free Finite Elements
Loading...
Searching...
No Matches
hyteg::dg::DGFunction< ValueType > Class Template Referencefinal

#include <DGFunction.hpp>

Inheritance diagram for hyteg::dg::DGFunction< ValueType >:
Collaboration diagram for hyteg::dg::DGFunction< ValueType >:

Public Types

template<typename VType >
using FunctionType = DGFunction< VType >
 
typedef FunctionTrait< DGFunction< ValueType > >::Tag Tag
 
typedef ValueType valueType
 

Public Member Functions

 DGFunction (const std::string &name, const std::shared_ptr< PrimitiveStorage > &storage, uint_t minLevel, uint_t maxLevel, const std::shared_ptr< DGBasisInfo > &basis, uint_t initialPolyDegree, BoundaryCondition boundaryCondition=BoundaryCondition::create0123BC())
 Instantiates a DGFunction.
 
void add (const std::vector< ValueType > scalars, const std::vector< std::reference_wrapper< const DGFunction< ValueType > > > &functions, uint_t level, DoFType flag=All) const
 
void add (const ValueType scalar, uint_t level, DoFType flag=All) const
 
void applyDirichletBoundaryConditions (const std::shared_ptr< DGForm > &form, uint_t level)
 Applies the Dirichlet boundary conditions to this function, treating it as the right-hand side of the linear system that corresponds to the form object.
 
void assign (const std::vector< ValueType > &scalars, const std::vector< std::reference_wrapper< const DGFunction< ValueType > > > &functions, uint_t level, DoFType flag=All)
 Assigns a linear combination of multiple VolumeDoFFunctions to this.
 
std::shared_ptr< DGBasisInfobasis () const
 Returns the associated DGBasisInfo object.
 
void communicate (uint_t level) const
 Updates ghost-layers.
 
template<typename OtherValueType >
void copyBoundaryConditionFromFunction (const DGFunction< OtherValueType > &other)
 
void copyFrom (const DGFunction< ValueType > &other, const uint_t &level, const std::map< PrimitiveID, uint_t > &localPrimitiveIDsToRank, const std::map< PrimitiveID, uint_t > &otherPrimitiveIDsToRank) const
 
ValueType dotGlobal (const DGFunction< ValueType > &rhs, uint_t level, DoFType flag=All) const
 Evaluates the (global) dot product. Involves communication and has to be called collectively.
 
ValueType dotLocal (const DGFunction< ValueType > &rhs, uint_t level, DoFType flag=All) const
 Evaluates the dot product on all local DoFs.
 
void enableTiming (const std::shared_ptr< walberla::WcTimingTree > &timingTree)
 
void enumerate (uint_t level) const
 Assigns unique values to all data points.
 
void enumerate (uint_t level, ValueType &offset) const
 Like enumerate() but starting with a value given by the offset parameter.
 
bool evaluate (const Point3D &physicalCoords, uint_t level, ValueType &value, real_t searchToleranceRadius=1e-05) const
 Evaluate finite element function at a specific coordinates.
 
void evaluateLinearFunctional (const std::function< real_t(const Point3D &) > &f, uint_t level)
 Evaluates the linear functional.
 
void evaluateOnMicroElement (const Point3D &coordinates, uint_t level, const PrimitiveID &cellID, hyteg::indexing::Index elementIndex, celldof::CellType cellType, ValueType &value) const
 Evaluate finite element function on a specific micro-cell.
 
void evaluateOnMicroElement (const Point3D &coordinates, uint_t level, const PrimitiveID &faceID, hyteg::indexing::Index elementIndex, facedof::FaceType faceType, ValueType &value) const
 Evaluate finite element function on a specific micro-face.
 
BoundaryCondition getBoundaryCondition () const
 Returns the boundary conditions of this function.
 
uint_t getDimension () const
 Query function object for the dimension of the field it represents.
 
const std::string & getFunctionName () const
 
ValueType getMax (uint_t level, bool mpiReduce=true) const
 Returns the max DoF value.
 
uint_t getMaxLevel ()
 Query function object for maximal level on which it defined.
 
ValueType getMaxMagnitude (uint_t level, bool mpiReduce=true) const
 Returns the max absolute DoF value.
 
ValueType getMin (uint_t level, bool mpiReduce=true) const
 Returns the min DoF value.
 
uint_t getMinLevel ()
 Query function object for minimal level on which it defined.
 
uint_t getNumberOfGlobalDoFs (uint_t level, const MPI_Comm &communicator=walberla::mpi::MPIManager::instance() ->comm(), const bool &onRootOnly=false) const
 Returns the number of DoFs.
 
uint_t getNumberOfLocalDoFs (uint_t level) const
 Returns the number of DoFs that are allocated on this process.
 
std::shared_ptr< PrimitiveStoragegetStorage () const
 
void interpolate (const std::function< ValueType(const hyteg::Point3D &) > &expr, uint_t level, DoFType flag=All) const
 
void interpolate (const std::vector< std::function< ValueType(const hyteg::Point3D &) > > &expressions, uint_t level, DoFType flag=All) const
 
void interpolate (ValueType constant, uint_t level, DoFType flag=All) const
 
void multElementwise (const std::vector< std::reference_wrapper< const DGFunction< ValueType > > > &functions, uint_t level, DoFType flag=All) const
 
int polynomialDegree (const PrimitiveID &primitiveID) const
 Returns the polynomial degree of the function on the passed primitive.
 
void setBoundaryCondition (BoundaryCondition bc)
 Sets the boundary conditions of this function.
 
ValueType sumGlobal (uint_t level, DoFType flag=All) const
 Evaluates the sum on all DoFs.
 
void swap (const DGFunction< ValueType > &other, const uint_t &level, const DoFType &flag=All) const
 
std::shared_ptr< volumedofspace::VolumeDoFFunction< ValueType > > volumeDoFFunction () const
 Returns the internally stored VolumeDoFFunction.
 
void toVector (const DGFunction< idx_t > &numerator, const std::shared_ptr< VectorProxy > &vec, uint_t level, DoFType flag) const
 conversion to/from linear algebra representation
 
void fromVector (const DGFunction< idx_t > &numerator, const std::shared_ptr< VectorProxy > &vec, uint_t level, DoFType flag) const
 conversion to/from linear algebra representation
 

Static Public Member Functions

static std::vector< std::string > getFunctionNames ()
 
static std::map< uint_t, uint_t > getLevelWiseFunctionCounter ()
 
static uint_t getNumFunctions ()
 

Protected Member Functions

void startTiming (const std::string &timerString) const
 
void stopTiming (const std::string &timerString) const
 

Protected Attributes

std::map< uint_t, std::shared_ptr< communication::BufferedCommunicator > > additiveCommunicators_
 
std::map< uint_t, std::shared_ptr< communication::BufferedCommunicator > > communicators_
 
std::string functionName_
 
std::shared_ptr< walberla::WcTimingTreetimingTree_
 

Private Attributes

std::shared_ptr< DGBasisInfobasis_
 
BoundaryCondition boundaryCondition_
 
uint_t maxLevel_
 
uint_t minLevel_
 
std::string name_
 
std::map< PrimitiveID, uint_t > polyDegreesPerPrimitive_
 
std::shared_ptr< PrimitiveStoragestorage_
 
std::shared_ptr< volumedofspace::VolumeDoFFunction< ValueType > > volumeDoFFunction_
 

Static Private Attributes

static std::vector< std::string > functionNames_
 
static std::map< uint_t, uint_t > levelWiseFunctionCounter_
 

Member Typedef Documentation

◆ FunctionType

template<typename ValueType >
template<typename VType >
using hyteg::dg::DGFunction< ValueType >::FunctionType = DGFunction< VType >

◆ Tag

typedef FunctionTrait<DGFunction< ValueType > >::Tag hyteg::Function< DGFunction< ValueType > >::Tag
inherited

◆ valueType

template<typename ValueType >
typedef ValueType hyteg::dg::DGFunction< ValueType >::valueType

Constructor & Destructor Documentation

◆ DGFunction()

template<typename ValueType >
hyteg::dg::DGFunction< ValueType >::DGFunction ( const std::string & name,
const std::shared_ptr< PrimitiveStorage > & storage,
uint_t minLevel,
uint_t maxLevel,
const std::shared_ptr< DGBasisInfo > & basis,
uint_t initialPolyDegree,
BoundaryCondition boundaryCondition = BoundaryCondition::create0123BC() )

Instantiates a DGFunction.

Parameters
namename of the function for output / debugging purposes
storagePrimitiveStorage that this lives on
minLevelmin refinement level to allocate
maxLevelmax refinement level to allocate
basispolynomial basis information through corresponding DGBasisInfo object
initialPolyDegreepolynomial degree to initialize with globally

Member Function Documentation

◆ add() [1/2]

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::add ( const std::vector< ValueType > scalars,
const std::vector< std::reference_wrapper< const DGFunction< ValueType > > > & functions,
uint_t level,
DoFType flag = All ) const
inline

◆ add() [2/2]

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::add ( const ValueType scalar,
uint_t level,
DoFType flag = All ) const
inline

◆ applyDirichletBoundaryConditions()

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::applyDirichletBoundaryConditions ( const std::shared_ptr< DGForm > & form,
uint_t level )

Applies the Dirichlet boundary conditions to this function, treating it as the right-hand side of the linear system that corresponds to the form object.

◆ assign()

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::assign ( const std::vector< ValueType > & scalars,
const std::vector< std::reference_wrapper< const DGFunction< ValueType > > > & functions,
uint_t level,
DoFType flag = All )
inline

Assigns a linear combination of multiple VolumeDoFFunctions to this.

◆ basis()

template<typename ValueType >
std::shared_ptr< DGBasisInfo > hyteg::dg::DGFunction< ValueType >::basis ( ) const
inline

Returns the associated DGBasisInfo object.

◆ communicate()

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::communicate ( uint_t level) const
inline

Updates ghost-layers.

◆ copyBoundaryConditionFromFunction()

template<typename ValueType >
template<typename OtherValueType >
void hyteg::dg::DGFunction< ValueType >::copyBoundaryConditionFromFunction ( const DGFunction< OtherValueType > & other)
inline

◆ copyFrom()

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::copyFrom ( const DGFunction< ValueType > & other,
const uint_t & level,
const std::map< PrimitiveID, uint_t > & localPrimitiveIDsToRank,
const std::map< PrimitiveID, uint_t > & otherPrimitiveIDsToRank ) const
inline

◆ dotGlobal()

template<typename ValueType >
ValueType hyteg::dg::DGFunction< ValueType >::dotGlobal ( const DGFunction< ValueType > & rhs,
uint_t level,
DoFType flag = All ) const
inline

Evaluates the (global) dot product. Involves communication and has to be called collectively.

◆ dotLocal()

template<typename ValueType >
ValueType hyteg::dg::DGFunction< ValueType >::dotLocal ( const DGFunction< ValueType > & rhs,
uint_t level,
DoFType flag = All ) const
inline

Evaluates the dot product on all local DoFs.

No communication is involved and the results may be different on each process.

◆ enableTiming()

void hyteg::Function< DGFunction< ValueType > >::enableTiming ( const std::shared_ptr< walberla::WcTimingTree > & timingTree)
inlineinherited

◆ enumerate() [1/2]

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::enumerate ( uint_t level) const

Assigns unique values to all data points.

Increments by 1 per DoF.

Mainly used for sparse matrix assembly to get global integer identifiers for all DoFs.

Parameters
levelrefinement level to enumerate

◆ enumerate() [2/2]

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::enumerate ( uint_t level,
ValueType & offset ) const

Like enumerate() but starting with a value given by the offset parameter.

Parameters
levelrefinement level to enumerate
offsetvalue to start from, the next "free" value is returned

◆ evaluate()

template<typename ValueType >
bool hyteg::dg::DGFunction< ValueType >::evaluate ( const Point3D & physicalCoords,
uint_t level,
ValueType & value,
real_t searchToleranceRadius = 1e-05 ) const

Evaluate finite element function at a specific coordinates.

In a parallel setting, the specified coordinate might not lie in the local subdomain.

Evaluation is performed in two steps:

  1. For all volume primitives of the local subdomain: If a point-tet (point-triangle in 2D) inclusion test succeeds, the function returns true and the finite-element function is evaluated.
  2. Skipped, if radius is negative. For all volume primitives of the local subdomain: A sphere-tet (circle-triangle) intersection test is performed, if successful returns true, and the finite-element function is extrapolated to the specified coordinate and evaluated, depending on the radius, this might introduce (large) errors.

If both tests fail, this function returns false, and no evaluation is performed (i.e. the returned, evaluated value is not set to anything meaningful).

Note that two parallel processes that return true, may return different values.

No communication is performed in this function. -> Does not need to be called collectively. -> Different values are returned on each process.

Parameters
physicalCoordscoordinates in physical domain where the function is to be evaluated
levelrefinement level
valuefunction value at the coordinate if search was successful
searchToleranceRadiusradius of the sphere (circle) for the second search phase, skipped if negative
Returns
true if the function was evaluated successfully, false otherwise

◆ evaluateLinearFunctional()

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::evaluateLinearFunctional ( const std::function< real_t(const Point3D &) > & f,
uint_t level )

Evaluates the linear functional.

\[ l( v ) = \int_\Omega f \cdot v \]

by integration over the local basis functions and writes the result into the vector, i.e.

\[ u_i \leftarrow \int_T f \cdot \phi_i \]

where \(\phi_i\) is the basis function associated with the DoF \(u_i\) and \(f\) a given analytical function.

◆ evaluateOnMicroElement() [1/2]

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::evaluateOnMicroElement ( const Point3D & coordinates,
uint_t level,
const PrimitiveID & cellID,
hyteg::indexing::Index elementIndex,
celldof::CellType cellType,
ValueType & value ) const

Evaluate finite element function on a specific micro-cell.

The standard evaluate() function does not take the discontinuity of the function into account. If the function is evaluated at a discontinuity with evaluate(), it is kind of "random" which of the neighboring elements is used to evaluate the polynom.

This is for example an issue for accurate visualization of discontinuities.

The present method solves this issue by taking a specific (micro-)element as input argument. This way the function can be evaluated on different elements at the same coordinate. Note that regardless of whether the specified coordinate is located on the element or not - the local polynomial is extrapolated.

Parameters
coordinateswhere the function shall be evaluated
levelrefinement level
cellIDthe macro-cell where the (micro-)element is located on
elementIndexthe logical index of the micro-element
cellTypethe type of the local micro-element
valuethe evaluation

◆ evaluateOnMicroElement() [2/2]

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::evaluateOnMicroElement ( const Point3D & coordinates,
uint_t level,
const PrimitiveID & faceID,
hyteg::indexing::Index elementIndex,
facedof::FaceType faceType,
ValueType & value ) const

Evaluate finite element function on a specific micro-face.

The standard evaluate() function does not take the discontinuity of the function into account. If the function is evaluated at a discontinuity with evaluate(), it is kind of "random" which of the neighboring elements is used to evaluate the polynom.

This is for example an issue for accurate visualization of discontinuities.

The present method solves this issue by taking a specific (micro-)element as input argument. This way the function can be evaluated on different elements at the same coordinate. Note that regardless of whether the specified coordinate is located on the element or not - the local polynomial is extrapolated.

Parameters
coordinateswhere the function shall be evaluated
levelrefinement level
faceIDthe macro-face where the (micro-)element is located on
elementIndexthe logical index of the micro-element
faceTypethe type of the local micro-element
valuethe evaluation

◆ fromVector()

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::fromVector ( const DGFunction< idx_t > & numerator,
const std::shared_ptr< VectorProxy > & vec,
uint_t level,
DoFType flag ) const

conversion to/from linear algebra representation

◆ getBoundaryCondition()

template<typename ValueType >
BoundaryCondition hyteg::dg::DGFunction< ValueType >::getBoundaryCondition ( ) const
inline

Returns the boundary conditions of this function.

◆ getDimension()

template<typename ValueType >
uint_t hyteg::dg::DGFunction< ValueType >::getDimension ( ) const
inlinevirtual

Query function object for the dimension of the field it represents.

Implements hyteg::Function< DGFunction< ValueType > >.

◆ getFunctionName()

const std::string & hyteg::Function< DGFunction< ValueType > >::getFunctionName ( ) const
inlineinherited

◆ getFunctionNames()

static std::vector< std::string > hyteg::Function< DGFunction< ValueType > >::getFunctionNames ( )
inlinestaticinherited

◆ getLevelWiseFunctionCounter()

static std::map< uint_t, uint_t > hyteg::Function< DGFunction< ValueType > >::getLevelWiseFunctionCounter ( )
inlinestaticinherited

◆ getMax()

template<typename ValueType >
ValueType hyteg::dg::DGFunction< ValueType >::getMax ( uint_t level,
bool mpiReduce = true ) const

Returns the max DoF value.

Parameters
levelrefinement level
mpiReduceif true, reduces over all processes (global max), if false returns the process local value

◆ getMaxLevel()

uint_t hyteg::Function< DGFunction< ValueType > >::getMaxLevel ( )
inlineinherited

Query function object for maximal level on which it defined.

◆ getMaxMagnitude()

template<typename ValueType >
ValueType hyteg::dg::DGFunction< ValueType >::getMaxMagnitude ( uint_t level,
bool mpiReduce = true ) const

Returns the max absolute DoF value.

Parameters
levelrefinement level
mpiReduceif true, reduces over all processes (global max magnitude), if false returns the process local value

◆ getMin()

template<typename ValueType >
ValueType hyteg::dg::DGFunction< ValueType >::getMin ( uint_t level,
bool mpiReduce = true ) const

Returns the min DoF value.

Parameters
levelrefinement level
mpiReduceif true, reduces over all processes (global min), if false returns the process local value

◆ getMinLevel()

uint_t hyteg::Function< DGFunction< ValueType > >::getMinLevel ( )
inlineinherited

Query function object for minimal level on which it defined.

◆ getNumberOfGlobalDoFs()

template<typename ValueType >
uint_t hyteg::dg::DGFunction< ValueType >::getNumberOfGlobalDoFs ( uint_t level,
const MPI_Comm & communicator = walberla::mpi::MPIManager::instance()->comm(),
const bool & onRootOnly = false ) const

Returns the number of DoFs.

Performs global reduction, must be called collectively.

Parameters
levelrefinement level
communicatorif required, a custom communicator can be passed
onRootOnlyif true, the result is only returned on the root process
Returns

◆ getNumberOfLocalDoFs()

template<typename ValueType >
uint_t hyteg::dg::DGFunction< ValueType >::getNumberOfLocalDoFs ( uint_t level) const

Returns the number of DoFs that are allocated on this process.

◆ getNumFunctions()

static uint_t hyteg::Function< DGFunction< ValueType > >::getNumFunctions ( )
inlinestaticinherited

◆ getStorage()

std::shared_ptr< PrimitiveStorage > hyteg::Function< DGFunction< ValueType > >::getStorage ( ) const
inlineinherited

◆ interpolate() [1/3]

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::interpolate ( const std::function< ValueType(const hyteg::Point3D &) > & expr,
uint_t level,
DoFType flag = All ) const
inline

◆ interpolate() [2/3]

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::interpolate ( const std::vector< std::function< ValueType(const hyteg::Point3D &) > > & expressions,
uint_t level,
DoFType flag = All ) const
inline

◆ interpolate() [3/3]

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::interpolate ( ValueType constant,
uint_t level,
DoFType flag = All ) const

◆ multElementwise()

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::multElementwise ( const std::vector< std::reference_wrapper< const DGFunction< ValueType > > > & functions,
uint_t level,
DoFType flag = All ) const
inline

◆ polynomialDegree()

template<typename ValueType >
int hyteg::dg::DGFunction< ValueType >::polynomialDegree ( const PrimitiveID & primitiveID) const
inline

Returns the polynomial degree of the function on the passed primitive.

◆ setBoundaryCondition()

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::setBoundaryCondition ( BoundaryCondition bc)
inline

Sets the boundary conditions of this function.

◆ startTiming()

void hyteg::Function< DGFunction< ValueType > >::startTiming ( const std::string & timerString) const
inlineprotectedinherited

◆ stopTiming()

void hyteg::Function< DGFunction< ValueType > >::stopTiming ( const std::string & timerString) const
inlineprotectedinherited

◆ sumGlobal()

template<typename ValueType >
ValueType hyteg::dg::DGFunction< ValueType >::sumGlobal ( uint_t level,
DoFType flag = All ) const
inline

Evaluates the sum on all DoFs.

◆ swap()

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::swap ( const DGFunction< ValueType > & other,
const uint_t & level,
const DoFType & flag = All ) const
inline

◆ toVector()

template<typename ValueType >
void hyteg::dg::DGFunction< ValueType >::toVector ( const DGFunction< idx_t > & numerator,
const std::shared_ptr< VectorProxy > & vec,
uint_t level,
DoFType flag ) const

conversion to/from linear algebra representation

◆ volumeDoFFunction()

template<typename ValueType >
std::shared_ptr< volumedofspace::VolumeDoFFunction< ValueType > > hyteg::dg::DGFunction< ValueType >::volumeDoFFunction ( ) const
inline

Returns the internally stored VolumeDoFFunction.

Member Data Documentation

◆ additiveCommunicators_

std::map< uint_t, std::shared_ptr< communication::BufferedCommunicator > > hyteg::Function< DGFunction< ValueType > >::additiveCommunicators_
protectedinherited

◆ basis_

template<typename ValueType >
std::shared_ptr< DGBasisInfo > hyteg::dg::DGFunction< ValueType >::basis_
private

◆ boundaryCondition_

template<typename ValueType >
BoundaryCondition hyteg::dg::DGFunction< ValueType >::boundaryCondition_
private

◆ communicators_

std::map< uint_t, std::shared_ptr< communication::BufferedCommunicator > > hyteg::Function< DGFunction< ValueType > >::communicators_
protectedinherited

◆ functionName_

std::string hyteg::Function< DGFunction< ValueType > >::functionName_
protectedinherited

◆ functionNames_

std::vector< std::string > hyteg::Function< DGFunction< ValueType > >::functionNames_
staticprivateinherited

◆ levelWiseFunctionCounter_

std::map< uint_t, uint_t > hyteg::Function< DGFunction< ValueType > >::levelWiseFunctionCounter_
staticprivateinherited

◆ maxLevel_

template<typename ValueType >
uint_t hyteg::dg::DGFunction< ValueType >::maxLevel_
private

◆ minLevel_

template<typename ValueType >
uint_t hyteg::dg::DGFunction< ValueType >::minLevel_
private

◆ name_

template<typename ValueType >
std::string hyteg::dg::DGFunction< ValueType >::name_
private

◆ polyDegreesPerPrimitive_

template<typename ValueType >
std::map< PrimitiveID, uint_t > hyteg::dg::DGFunction< ValueType >::polyDegreesPerPrimitive_
private

◆ storage_

template<typename ValueType >
std::shared_ptr< PrimitiveStorage > hyteg::dg::DGFunction< ValueType >::storage_
private

◆ timingTree_

std::shared_ptr< walberla::WcTimingTree > hyteg::Function< DGFunction< ValueType > >::timingTree_
protectedinherited

◆ volumeDoFFunction_

template<typename ValueType >
std::shared_ptr< volumedofspace::VolumeDoFFunction< ValueType > > hyteg::dg::DGFunction< ValueType >::volumeDoFFunction_
private

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