Hermes2D
3.0
|
Should be exactly the same as is the count of enum ShapesetType. More...
#include <shapeset.h>
Public Types | |
typedef double(* | shape_fn_t )(double, double) |
Shape-function function type. Internal. | |
Public Member Functions | |
unsigned short | get_order (int index, ElementMode2D mode) const |
virtual Shapeset * | clone ()=0 |
unsigned char | get_num_components () const |
Returns 2 if this is a vector shapeset, 1 otherwise. | |
unsigned short | get_max_order () const |
Returns the maximum poly degree for all shape functions. | |
unsigned short | get_min_order () const |
virtual unsigned short | get_max_index (ElementMode2D mode) const =0 |
Returns the highest shape function index. | |
short | get_vertex_index (int vertex, ElementMode2D mode) const |
Returns the index of a vertex shape function associated with the specified vertex. | |
short | get_edge_index (unsigned char edge, unsigned short ori, unsigned short order, ElementMode2D mode) const |
virtual SpaceType | get_space_type () const =0 |
virtual unsigned char | get_id () const =0 |
double | get_value (int n, int index, double x, double y, unsigned short component, ElementMode2D mode) |
double | get_fn_value (int index, double x, double y, unsigned short component, ElementMode2D mode) |
double | get_dx_value (int index, double x, double y, unsigned short component, ElementMode2D mode) |
double | get_dy_value (int index, double x, double y, unsigned short component, ElementMode2D mode) |
double | get_fn_value_0_tri (int index, double x, double y) |
The most used calls are distinguished for optimization. | |
double | get_dx_value_0_tri (int index, double x, double y) |
The most used calls are distinguished for optimization. | |
double | get_dy_value_0_tri (int index, double x, double y) |
The most used calls are distinguished for optimization. | |
double | get_fn_value_0_quad (int index, double x, double y) |
The most used calls are distinguished for optimization. | |
double | get_dx_value_0_quad (int index, double x, double y) |
The most used calls are distinguished for optimization. | |
double | get_dy_value_0_quad (int index, double x, double y) |
The most used calls are distinguished for optimization. | |
double | get_dxx_value (int index, double x, double y, unsigned short component, ElementMode2D mode) |
double | get_dyy_value (int index, double x, double y, unsigned short component, ElementMode2D mode) |
double | get_dxy_value (int index, double x, double y, unsigned short component, ElementMode2D mode) |
virtual unsigned short | get_num_bubbles (unsigned short order, ElementMode2D mode) const |
Returns the number of bubble functions for an element of the given order. | |
Protected Member Functions | |
virtual short * | get_bubble_indices (unsigned short order, ElementMode2D mode) const |
Returns a complete set of indices of bubble functions for an element of the given order. | |
int | get_constrained_edge_index (unsigned char edge, unsigned short order, unsigned short ori, unsigned short part, ElementMode2D mode) const |
double2 * | get_ref_vertex (int vertex, ElementMode2D mode) |
Returns the coordinates of the reference domain vertices. | |
double * | calculate_constrained_edge_combination (unsigned short order, unsigned short part, unsigned short ori, ElementMode2D mode) |
double * | get_constrained_edge_combination (unsigned short order, unsigned short part, unsigned short ori, unsigned short &nitems, ElementMode2D mode) |
void | free_constrained_edge_combinations () |
Releases all cached coefficients. | |
double | get_constrained_value (int n, int index, double x, double y, unsigned short component, ElementMode2D mode) |
Protected Attributes | |
shape_fn_t *** | shape_table [6] |
short ** | vertex_indices |
short *** | edge_indices |
short *** | bubble_indices |
unsigned short ** | bubble_count |
unsigned short ** | index_to_order |
double2 | ref_vert [H2D_MAX_SOLUTION_COMPONENTS][H2D_MAX_NUMBER_VERTICES] |
unsigned char | max_order |
unsigned char | min_order |
unsigned char | num_components |
unsigned short | ebias |
2 for H1 shapesets, 0 for H(curl) shapesets. It is the order of the More... | |
double ** | comb_table |
unsigned short | table_size |
Should be exactly the same as is the count of enum ShapesetType.
Defines a set of shape functions.
This class stores mainly the definitions of the polynomials for all shape functions, but also their polynomial degrees, their types (vertex, edge, bubble), and contains mechanisms for the calculation and storage of constrained shape functions.
The class returns shape function values for both triangles and quads, depending on what mode is requested.
Each shape function is assigned a unique number - 'index'. For standard shape functions, index is positive and not greater than the value returned by get_max_index(). Negative index values are reserved for constrained shape functions. These are special edge functions designed to fit standard edge functions along a portion of an edge, used in the construction of FE spaces on meshes with hanging nodes.
The usage of Shapeset is simple: you first obtain the index of the shape function you are interested in, be it a vertex function for a given vertex, or an edge function of a given order, etc. Then you call one of the functions get_fn_value(), get_dx_value(), etc. All shape functions are defined on the reference domain. For triangles, this is the standard triangle (-1,-1), (1,-1), (-1,1), and for quads this is the square (-1,1)^2.
The polynomial degree (or 'order') is an integer typically in the range[1-10] for H1 shapesets and[0-10] for H(curl) shapesets. Quadrilaterals are allowed to have different orders in the x and y directions (of the reference domain). The 'order' for quads thus has to be formed with the macro H2D_MAKE_QUAD_ORDER(), see h2d_common.h.
Vertex shape functions in H1 shapesets are also regarded as edge functions of orders 0 and 1. This simplifies constraint calculations and BC projections.
Shape functions are always Real-valued.
Very important: Should one create a custom shapeset, or enlarge the existing ones, one must increase the value of H2D_MAX_LOCAL_BASIS_SIZE.
Definition at line 95 of file shapeset.h.
|
protected |
numbering of edge intervals: (the variable 'part') -+- -+- -+- | | 13 | | 5 | -+- | | 12 | 1 | -+- -+- finer interval: | | 11 | | 4 | -+- p = (p + 1) * 2 (+1) | | 10 | -+- -+- -+- ... etc. | | 9 | | 3 | -+- | | 8 | 0 | -+- -+- | | 7 | | 2 | -+- | | 6 | -+- -+- -+- Constrained edge functions are constructed by subtracting the linear part (ie., two vertex functions) from the constraining edge function and expressing the rest as a linear combination of standard edge functions. This function determines the coefficients of such linear combination by forming and solving a simple linear system.
Definition at line 27 of file shapeset.cpp.
|
protected |
Returns the coefficients for the linear combination forming a constrained edge function. This function performs the storage (caching) of these coefficients, so that they can be calculated only once.
Definition at line 102 of file shapeset.cpp.
|
protected |
Returns the index of a constrained edge function. 'part' is 0 or 1 for edge halves, 2, 3, 4, 5 for edge quarters, etc. See shapeset.cpp.
Definition at line 242 of file shapeset.cpp.
|
protected |
Constructs the linear combination of edge functions, forming a constrained edge function.
Definition at line 163 of file shapeset.cpp.
short Hermes::Hermes2D::Shapeset::get_edge_index | ( | unsigned char | edge, |
unsigned short | ori, | ||
unsigned short | order, | ||
ElementMode2D | mode | ||
) | const |
Returns the index of an edge function associated with the specified edge and of the requested order. 'ori' can be 0 or 1 and determines edge orientation (this is for shapesets with non-symmetric edge functions).
Definition at line 208 of file shapeset.cpp.
Referenced by Hermes::Hermes2D::RefMap::set_active_element(), and Hermes::Hermes2D::RefMap::untransform().
|
pure virtual |
Returns shapeset identifier. Internal.
Implemented in Hermes::Hermes2D::L2ShapesetTaylor, Hermes::Hermes2D::H1ShapesetJacobi, Hermes::Hermes2D::L2ShapesetLegendre, Hermes::Hermes2D::HdivShapesetLegendre, and Hermes::Hermes2D::HcurlShapesetGradLeg.
Referenced by Hermes::Hermes2D::PrecalcShapesetAssembling::PrecalcShapesetAssembling(), Hermes::Hermes2D::HcurlSpace< Scalar >::set_shapeset(), Hermes::Hermes2D::HdivSpace< Scalar >::set_shapeset(), Hermes::Hermes2D::H1Space< Scalar >::set_shapeset(), and Hermes::Hermes2D::L2Space< Scalar >::set_shapeset().
unsigned short Hermes::Hermes2D::Shapeset::get_order | ( | int | index, |
ElementMode2D | mode | ||
) | const |
Returns the polynomial degree of the specified shape function. If on quads, it returns encoded orders. The orders has to be decoded through macros H2D_GET_H_ORDER and H2D_GET_V_ORDER.
Definition at line 256 of file shapeset.cpp.
Referenced by Hermes::Hermes2D::PrecalcShapeset::get_edge_fn_order(), and Hermes::Hermes2D::PrecalcShapeset::set_active_shape().
|
pure virtual |
Returns space type. Internal.
Implemented in Hermes::Hermes2D::L2ShapesetTaylor, Hermes::Hermes2D::HcurlShapesetGradLeg, Hermes::Hermes2D::L2ShapesetLegendre, Hermes::Hermes2D::H1ShapesetJacobi, and Hermes::Hermes2D::HdivShapesetLegendre.
Referenced by Hermes::Hermes2D::PrecalcShapeset::get_space_type(), and Hermes::Hermes2D::Space< Scalar >::init_empty_space().
double Hermes::Hermes2D::Shapeset::get_value | ( | int | n, |
int | index, | ||
double | x, | ||
double | y, | ||
unsigned short | component, | ||
ElementMode2D | mode | ||
) |
Obtains the value of the given shape function. (x,y) is a coordinate in the reference domain, component is 0 for Scalar shapesets and 0 or 1 for vector shapesets.
Definition at line 264 of file shapeset.cpp.
Referenced by Hermes::Hermes2D::PrecalcShapeset::precalculate().
|
protected |
2 for H1 shapesets, 0 for H(curl) shapesets. It is the order of the
first edge function.
Definition at line 189 of file shapeset.h.