Hermes2D
3.0
|
Represents an arbitrary function defined on an element. More...
#include <function.h>
Public Member Functions | |
Function () | |
Default constructor. | |
virtual | ~Function () |
unsigned char | get_num_components () const |
Returns the number of components of the function being represented by the class. | |
virtual const Scalar * | get_fn_values (int component=0) const |
Returns function values. More... | |
virtual const Scalar * | get_dx_values (int component=0) const |
Returns the x partial derivative. More... | |
virtual const Scalar * | get_dy_values (int component=0) const |
Returns the y partial derivative. More... | |
Scalar * | deep_copy_array (int component=0, int item=0) const |
Returns function values. More... | |
Quad2D * | get_quad_2d () const |
Returns the current quadrature points. | |
void | set_quad_order (unsigned short order, unsigned short mask=H2D_FN_DEFAULT) |
virtual const Scalar * | get_values (int component, int item) const |
virtual int | get_fn_order () const |
Returns the polynomial degree of the function being represented by the class. | |
virtual void | push_transform (int son) |
virtual void | pop_transform () |
virtual void | set_active_element (Element *e) |
Sets the active element. | |
virtual void | set_transform (uint64_t idx) |
Public Member Functions inherited from Hermes::Hermes2D::Transformable | |
Element * | get_active_element () const |
uint64_t | get_transform () const |
Protected Member Functions | |
virtual void | set_quad_2d (Quad2D *quad_2d) |
Selects the quadrature points in which the function will be evaluated. More... | |
virtual void | reset_transform () |
Empties the stack, loads identity transform. | |
virtual void | force_transform (uint64_t sub_idx, Trf *ctm) |
For internal use only. | |
virtual void | free ()=0 |
Frees all precalculated tables. | |
virtual int | get_edge_fn_order (unsigned char edge) const |
Returns the polynomial degree of the function at given edge. To be overridden in derived classes. More... | |
virtual void | precalculate (unsigned short order, unsigned short mask) |
precalculates the current function at the current integration points. | |
void | invalidate_values () |
Protected Member Functions inherited from Hermes::Hermes2D::Transformable | |
double | get_transform_jacobian () const |
Trf * | get_ctm () const |
unsigned int | get_depth () const |
Protected Attributes | |
Scalar | values [H2D_MAX_SOLUTION_COMPONENTS][H2D_NUM_FUNCTION_VALUES][H2D_MAX_INTEGRATION_POINTS_COUNT] |
The data. | |
bool | values_valid |
Flag that the data are not 'dirty'. | |
int | order |
Current function polynomial order. | |
unsigned char | num_components |
Number of vector components. | |
Quad2D * | quads [H2D_MAX_QUADRATURES] |
List of available quadratures. | |
int | cur_quad |
Active quadrature (index into 'quads') | |
Protected Attributes inherited from Hermes::Hermes2D::Transformable | |
Element * | element |
The active element. | |
Trf * | ctm |
Current sub-element transformation matrix. | |
uint64_t | sub_idx |
Sub-element transformation index. | |
Trf | stack [21] |
Transformation matrix stack. | |
unsigned int | top |
Stack top. | |
Static Protected Attributes | |
static int | idx2mask [H2D_NUM_FUNCTION_VALUES][2] |
Index to mask table. More... | |
Static Protected Attributes inherited from Hermes::Hermes2D::Transformable | |
static const uint64_t | H2D_MAX_IDX = (1ULL << 3 * H2D_MAX_TRN_LEVEL) - 1 |
Friends | |
template<typename T > | |
class | KellyTypeAdapt |
template<typename T > | |
class | RefinementSelectors::H1ProjBasedSelector |
template<typename T > | |
class | RefinementSelectors::L2ProjBasedSelector |
template<typename T > | |
class | RefinementSelectors::HcurlProjBasedSelector |
template<typename T > | |
class | Adapt |
template<typename T > | |
class | DiscontinuousFunc |
template<typename T > | |
class | DiscreteProblem |
class | CurvMap |
template<typename T > | |
class | Func |
template<typename T > | |
class | Filter |
template<typename T > | |
class | SimpleFilter |
template<typename T > | |
class | DXDYFilter |
class | ComplexFilter |
class | VonMisesFilter |
HERMES_API Func< double > * | init_fn (PrecalcShapeset *fu, RefMap *rm, const int order) |
Init the shape function for the evaluation of the volumetric/surface integral (transformation of values). | |
template<typename T > | |
HERMES_API Func< T > * | init_fn (MeshFunction< T > *fu, const int order) |
Additional Inherited Members | |
Static Public Attributes inherited from Hermes::Hermes2D::Transformable | |
static const unsigned int | H2D_MAX_TRN_LEVEL = 15 |
If this changes, NeighborSearch::H2D_MAX_NEIGHBORS must change too. | |
Static Protected Member Functions inherited from Hermes::Hermes2D::Transformable | |
static void | push_transforms (std::set< Transformable * > &transformables, int son) |
static void | pop_transforms (std::set< Transformable * > &transformables) |
Represents an arbitrary function defined on an element.
The Function class is an abstraction of a function defined in integration points on an element. You first specify what quadrature tables you want to use (set_quad_2d()) and select an element (Transformable::set_active_element()). Then you select concrete integration points (set_quad_order()) and obtain the function values by calling one of the functions get_fn_values(), get_dx_values(), etc.
This class is a template for RealFunction and ScalarFunction, depending of which type the function values are. For example, shape functions are always Real (see PrecalcShapeset), while the solution can be complex (see Solution).
The design goal for this class is to define a single common interface for functions used as integrands in the weak formulation. It should not matter whether you are integrating a shape function or, for example, a previous solution of the PDE in time-dependent problems. Ideally, you should also be able to apply the bilinear form not only to shape functions during assembling, but also to the solution when calculating energy norms etc. The last feature is unfortunately limited to Real code, because a PDE solution can be complex (hence Solution inherits from ScalarFunction), but shape functions are Real and for efficiency the bilinear form only takes RealFunction arguments.
Since this class inherits from Transformable, you can obtain function values in integration points transformed to sub-areas of the current element (see push_transform(), pop_transform()).
Definition at line 106 of file function.h.
|
virtual |
Default destructor. Added to delete Function::sub_tables.
Definition at line 42 of file function.cpp.
Scalar * Hermes::Hermes2D::Function< Scalar >::deep_copy_array | ( | int | component = 0 , |
int | item = 0 |
||
) | const |
Returns function values.
component[in] | The component of the function (0 or 1). |
Definition at line 212 of file function.cpp.
|
virtual |
Returns the x partial derivative.
component[in] | The component of the function (0 or 1). |
Reimplemented in Hermes::Hermes2D::PrecalcShapesetAssembling.
Definition at line 175 of file function.cpp.
Referenced by Hermes::Hermes2D::init_fn_preallocated(), Hermes::Hermes2D::RefinementSelectors::H1ProjBasedSelector< Scalar >::precalc_ref_solution(), Hermes::Hermes2D::RefinementSelectors::HcurlProjBasedSelector< Scalar >::precalc_ref_solution(), and Hermes::Hermes2D::LinearFilter< Scalar >::precalculate().
|
virtual |
Returns the y partial derivative.
component[in] | The component of the function (0 or 1). |
Reimplemented in Hermes::Hermes2D::PrecalcShapesetAssembling.
Definition at line 182 of file function.cpp.
Referenced by Hermes::Hermes2D::init_fn_preallocated(), Hermes::Hermes2D::RefinementSelectors::H1ProjBasedSelector< Scalar >::precalc_ref_solution(), Hermes::Hermes2D::RefinementSelectors::HcurlProjBasedSelector< Scalar >::precalc_ref_solution(), and Hermes::Hermes2D::LinearFilter< Scalar >::precalculate().
|
protectedvirtual |
Returns the polynomial degree of the function at given edge. To be overridden in derived classes.
edge[in] | Edge at which the order should be evaluated. (0-3) |
Definition at line 53 of file function.cpp.
|
virtual |
Returns function values.
component[in] | The component of the function (0 or 1). |
Reimplemented in Hermes::Hermes2D::PrecalcShapesetAssembling.
Definition at line 168 of file function.cpp.
Referenced by Hermes::Hermes2D::init_fn_preallocated(), Hermes::Hermes2D::RefinementSelectors::L2ProjBasedSelector< Scalar >::precalc_ref_solution(), Hermes::Hermes2D::RefinementSelectors::H1ProjBasedSelector< Scalar >::precalc_ref_solution(), Hermes::Hermes2D::RefinementSelectors::HcurlProjBasedSelector< Scalar >::precalc_ref_solution(), and Hermes::Hermes2D::LinearFilter< Scalar >::precalculate().
|
protected |
With changed sub-element mapping, or an element, or anything else there comes the need for a change of the current values. This invalidates the current values. See values_valid
Definition at line 135 of file function.cpp.
|
virtual |
See Transformable::pop_transform. Internal.
Reimplemented from Hermes::Hermes2D::Transformable.
Reimplemented in Hermes::Hermes2D::ComplexFilter, Hermes::Hermes2D::Filter< Scalar >, Hermes::Hermes2D::Filter< double >, and Hermes::Hermes2D::Filter< std::complex< double > >.
Definition at line 122 of file function.cpp.
Referenced by Hermes::Hermes2D::Filter< Scalar >::pop_transform(), and Hermes::Hermes2D::ComplexFilter::pop_transform().
|
virtual |
See Transformable::push_transform. Internal.
Reimplemented from Hermes::Hermes2D::Transformable.
Reimplemented in Hermes::Hermes2D::ComplexFilter, Hermes::Hermes2D::Filter< Scalar >, Hermes::Hermes2D::Filter< double >, and Hermes::Hermes2D::Filter< std::complex< double > >.
Definition at line 115 of file function.cpp.
Referenced by Hermes::Hermes2D::RefinementSelectors::ProjBasedSelector< Scalar >::calc_projection_errors(), Hermes::Hermes2D::Filter< Scalar >::push_transform(), and Hermes::Hermes2D::ComplexFilter::push_transform().
|
protectedvirtual |
Selects the quadrature points in which the function will be evaluated.
It is possible to switch back and forth between different quadrature points: no precalculated values are freed. The standard quadrature is always selected by default already.
quad_2d[in] | The quadrature points. |
Reimplemented in Hermes::Hermes2D::ComplexFilter, Hermes::Hermes2D::MeshFunction< Scalar >, Hermes::Hermes2D::MeshFunction< double >, Hermes::Hermes2D::MeshFunction< std::complex< double > >, Hermes::Hermes2D::PrecalcShapeset, Hermes::Hermes2D::Filter< Scalar >, Hermes::Hermes2D::Filter< double >, and Hermes::Hermes2D::Filter< std::complex< double > >.
Definition at line 79 of file function.cpp.
Referenced by Hermes::Hermes2D::PrecalcShapeset::set_quad_2d(), and Hermes::Hermes2D::UExtFunction< Scalar >::UExtFunction().
void Hermes::Hermes2D::Function< Scalar >::set_quad_order | ( | unsigned short | order, |
unsigned short | mask = H2D_FN_DEFAULT |
||
) |
Activates an integration rule of the specified order. Subsequent calls to get_values(), get_dx_values() etc. will be returning function values at these points.
order[in] | Integration rule order. |
mask[in] | A combination of one or more of the constants H2D_FN_VAL, H2D_FN_DX, H2D_FN_DY, ... |
Definition at line 59 of file function.cpp.
Referenced by Hermes::Hermes2D::RefinementSelectors::ProjBasedSelector< Scalar >::calc_projection_errors(), Hermes::Hermes2D::CurvMapStatic::calculate_bubble_projection_matrix(), and Hermes::Hermes2D::init_fn_preallocated().
|
virtual |
Sets the current transform at once as if it was created by multiple calls to push_transform().
idx[in] | The number of the sub-element, as returned by get_transform(). |
Reimplemented from Hermes::Hermes2D::Transformable.
Definition at line 155 of file function.cpp.
Referenced by Hermes::Hermes2D::NeighborSearch< Scalar >::init_ext_fn(), and Hermes::Hermes2D::CurvMap::update_refmap_coeffs().
|
staticprotected |
Index to mask table.
Definition at line 231 of file function.h.