Hermes2D  3.0
Hermes::Hermes2D::Mesh Class Reference

Represents a finite element mesh. Typical usage: MeshSharedPtr mesh; Hermes::Hermes2D::MeshReaderH2DXML mloader; try {  mloader.load("mesh.xml", &mesh); } catch(Exceptions::MeshLoadFailureException& e) {  e.print_msg();  return -1; }. More...

#include <mesh.h>

+ Inheritance diagram for Hermes::Hermes2D::Mesh:

Classes

class  BoundaryMarkersConversion
 
class  ElementMarkersConversion
 
class  MarkersConversion
 
class  ReferenceMeshCreator
 Class for creating reference mesh. More...
 

Public Member Functions

void init (int size=H2D_DEFAULT_HASH_SIZE)
 
virtual bool isOkay () const
 State querying helpers.
 
std::string getClassName () const
 
bool rescale (double x_ref, double y_ref)
 Rescales the mesh.
 
void copy (MeshSharedPtr mesh)
 Creates a copy of another mesh.
 
void copy_base (MeshSharedPtr mesh)
 Copies the coarsest elements of another mesh.
 
void copy_converted (MeshSharedPtr mesh)
 Copies the active elements of a converted mesh.
 
void create (int nv, double2 *verts, int nt, int3 *tris, std::string *tri_markers, int nq, int4 *quads, std::string *quad_markers, int nm, int2 *mark, std::string *boundary_markers)
 Creates a mesh from given vertex, triangle, quad, and marker arrays.
 
Elementelement_on_physical_coordinates (double x, double y)
 
double get_marker_area (int marker)
 
Elementget_element (int id) const
 Retrieves an element by its id number.
 
int get_num_elements () const
 Returns the total number of elements stored.
 
int get_num_base_elements () const
 Returns the number of base mesh elements.
 
int get_num_used_base_elements () const
 Returns the number of used base mesh elements.
 
int get_num_active_elements () const
 Returns the current number of active elements in the mesh.
 
int get_max_element_id () const
 Returns the maximum node id number plus one.
 
int get_num_vertex_nodes () const
 Returns the number of vertex nodes.
 
int get_num_edge_nodes () const
 Returns the number of edge nodes.
 
void get_bounding_box (double &bottom_left_x, double &bottom_left_y, double &top_right_x, double &top_right_y)
 
Elementget_element_fast (int id) const
 For internal use.
 
unsigned get_seq () const
 For internal use.
 
void refine_element_id (int id, int refinement=0)
 
void refine_all_elements (int refinement=0, bool mark_as_initial=false)
 
void refine_by_criterion (int(*criterion)(Element *e), int depth=1, bool mark_as_initial=false)
 
void refine_towards_vertex (int vertex_id, int depth=1, bool mark_as_initial=false)
 
void refine_towards_boundary (std::vector< std::string > markers, int depth=1, bool aniso=true, bool mark_as_initial=false)
 
void refine_towards_boundary (std::string marker, int depth=1, bool aniso=true, bool mark_as_initial=false)
 
void refine_in_area (std::string marker, int depth=1, int refinement=0, bool mark_as_initial=false)
 Refines all element sharing the marker passed.
 
void refine_in_areas (std::vector< std::string > markers, int depth=1, int refinement=0, bool mark_as_initial=false)
 Refines all element sharing the markers passed.
 
int * regularize (int n)
 
void unrefine_element_id (int id)
 
void unrefine_all_elements (bool keep_initial_refinements=true)
 
void free ()
 Frees all data associated with the mesh.
 
void set_seq (unsigned seq)
 For internal use.
 
void initial_single_check ()
 For internal use.
 
const ElementMarkersConversionget_element_markers_conversion () const
 
const BoundaryMarkersConversionget_boundary_markers_conversion () const
 
ElementMarkersConversionget_element_markers_conversion ()
 
BoundaryMarkersConversionget_boundary_markers_conversion ()
 
Elementcreate_quad (int marker, Node *v0, Node *v1, Node *v2, Node *v3, CurvMap *cm, int id=-1)
 
Elementcreate_triangle (int marker, Node *v0, Node *v1, Node *v2, CurvMap *cm, int id=-1)
 
void refine_element (Element *e, int refinement)
 
void refine_quad (Element *e, int refinement, Element **sons_out=nullptr)
 
void refine_triangle_to_triangles (Element *e, Element **sons=nullptr)
 
- Public Member Functions inherited from Hermes::Hermes2D::HashTable
Nodeget_node (int id) const
 Retrieves a node by its id number.
 
int get_max_node_id () const
 Returns the maximum node id number plus one.
 
Nodepeek_vertex_node (int p1, int p2) const
 Returns a vertex node with parent id's p1 and p2 if it exists, nullptr otherwise.
 
Nodepeek_edge_node (int p1, int p2) const
 Returns an edge node with parent id's p1 and p2 if it exists, nullptr otherwise.
 
Array< Node > & get_nodes ()
 

Static Public Member Functions

static double vector_length (double a_1, double a_2)
 Computing vector length.
 
static bool same_line (double p_1, double p_2, double q_1, double q_2, double r_1, double r_2)
 Checking whether the points p, q, r lie on the same line.
 
static bool is_convex (double a_1, double a_2, double b_1, double b_2)
 Checking whether the angle of vectors 'a' and 'b' is between zero and Pi.
 
static void check_triangle (int i, Node *&v0, Node *&v1, Node *&v2)
 
static void check_quad (int i, Node *&v0, Node *&v1, Node *&v2, Node *&v3)
 

Public Attributes

MeshHashGridmeshHashGrid
 
std::map< int, MarkerArea * > marker_areas
 
int nbase
 
int ntopvert
 
int ninitial
 
int nactive
 
ElementMarkersConversion element_markers_conversion
 
BoundaryMarkersConversion boundary_markers_conversion
 
Array< Elementelements
 
unsigned seq
 
std::vector< std::pair
< unsigned int, int > > 
refinements
 

Friends

class MeshHashGrid
 
class MeshReaderH2D
 
class MeshReaderH2DBSON
 
class MeshReaderH2DXML
 
class MeshReaderH1DXML
 
class MeshReaderExodusII
 
class DiscreteProblem< double >
 
class DiscreteProblem< std::complex< double > >
 
template<typename Scalar >
class DiscreteProblemDGAssembler
 
class WeakForm< double >
 
class WeakForm< std::complex< double > >
 
template<typename Scalar >
class Adapt
 
class KellyTypeAdapt< double >
 
class KellyTypeAdapt< std::complex< double > >
 
template<typename Scalar >
class Solution
 
template<typename Scalar >
class Filter
 
template<typename Scalar >
class MeshFunction
 
class RefMap
 
class Traverse
 
class Transformable
 
class Curved
 
class Views::MeshView
 
template<typename Scalar >
class RefinementSelectors::Selector
 
template<typename Scalar >
class RefinementSelectors::POnlySelector
 
template<typename Scalar >
class RefinementSelectors::HOnlySelector
 
template<typename Scalar >
class RefinementSelectors::OptimumSelector
 
template<typename Scalar >
class RefinementSelectors::ProjBasedSelector
 
template<typename Scalar >
class RefinementSelectors::H1ProjBasedSelector
 
template<typename Scalar >
class RefinementSelectors::L2ProjBasedSelector
 
template<typename Scalar >
class RefinementSelectors::HcurlProjBasedSelector
 
class PrecalcShapeset
 
template<typename Scalar >
class Space
 
template<typename Scalar >
class H1Space
 
template<typename Scalar >
class HcurlSpace
 
template<typename Scalar >
class HdivSpace
 
template<typename Scalar >
class L2Space
 
class Views::ScalarView
 
class Views::Orderizer
 
class EggShell
 

Additional Inherited Members

- Static Public Attributes inherited from Hermes::Hermes2D::HashTable
static const int H2D_DEFAULT_HASH_SIZE = 0x8000
 32K entries
 
- Protected Member Functions inherited from Hermes::Hermes2D::HashTable
int get_num_nodes () const
 Returns the total number of nodes stored.
 
Nodeget_vertex_node (int p1, int p2)
 
Nodeget_edge_node (int p1, int p2)
 
void init (int size=H2D_DEFAULT_HASH_SIZE)
 
void copy (const HashTable *ht)
 Copies another hash table contents.
 
void rebuild ()
 Reconstructs the hashtable, after, e.g., the nodes have been loaded from a file.
 
void free ()
 Frees all memory used by the instance.
 
void remove_vertex_node (int id)
 Removes a vertex node with parent id's p1 and p2.
 
void remove_edge_node (int id)
 Removes an edge node with parent id's p1 and p2.
 
- Protected Attributes inherited from Hermes::Hermes2D::HashTable
Array< Nodenodes
 Array storing all nodes.
 

Detailed Description

Represents a finite element mesh. Typical usage: MeshSharedPtr mesh; Hermes::Hermes2D::MeshReaderH2DXML mloader; try {  mloader.load("mesh.xml", &mesh); } catch(Exceptions::MeshLoadFailureException& e) {  e.print_msg();  return -1; }.

Definition at line 61 of file mesh.h.

Member Function Documentation

Element * Hermes::Hermes2D::Mesh::element_on_physical_coordinates ( double  x,
double  y 
)

Returns the element pointer located at physical coordinates x, y.

Parameters
[in]xPhysical x-coordinate.
[in]yPhysical y-coordinate.
[in]x_referenceOptional parameter, in which the x-coordinate of x in the reference domain will be returned.
[in]y_referenceOptional parameter, in which the y-coordinate of y in the reference domain will be returned.

Definition at line 1279 of file mesh.cpp.

void Hermes::Hermes2D::Mesh::get_bounding_box ( double &  bottom_left_x,
double &  bottom_left_y,
double &  top_right_x,
double &  top_right_y 
)

Get the mesh bounding box.

Parameters
[out]bottom_left_xBottom left corner - x coordinate.
[out]bottom_left_yBottom left corner - y coordinate.
[out]top_right_xTop right corner - x coordinate.
[out]top_right_yTop right corner - y coordinate.

Definition at line 293 of file mesh.cpp.

void Hermes::Hermes2D::Mesh::init ( int  size = H2D_DEFAULT_HASH_SIZE)

Initializes the mesh.

Parameters
size[in]Hash table size; must be a power of two.

Definition at line 1188 of file mesh.cpp.

Referenced by copy_base(), and create().

void Hermes::Hermes2D::Mesh::refine_all_elements ( int  refinement = 0,
bool  mark_as_initial = false 
)

Refines all elements.

Parameters
refinement[in]Same meaning as in refine_element_id().

Definition at line 752 of file mesh.cpp.

Referenced by Hermes::Hermes2D::Mesh::ReferenceMeshCreator::create_ref_mesh().

void Hermes::Hermes2D::Mesh::refine_by_criterion ( int(*)(Element *e)  criterion,
int  depth = 1,
bool  mark_as_initial = false 
)

Selects elements to refine according to a given criterion and performs 'depth' levels of refinements. The criterion function receives a pointer to an element to be considered. It must return -1 if the element is not to be refined, 0 if it should be refined uniformly, 1 if it is a quad and should be split horizontally, 2 if it is a quad and should be split vertically, and 3 if it is a triangle and should be split into three quads.

Definition at line 776 of file mesh.cpp.

Referenced by refine_towards_boundary(), and refine_towards_vertex().

void Hermes::Hermes2D::Mesh::refine_element_id ( int  id,
int  refinement = 0 
)

Refines an element.

Parameters
id[in]Element id number.
refinement[in]Ignored for triangles. If the element is a quad, 0 means refine in both directions, 1 means refine horizontally (with respect to the reference domain), 2 means refine vertically.

Definition at line 740 of file mesh.cpp.

Referenced by refine_by_criterion(), and regularize().

void Hermes::Hermes2D::Mesh::refine_quad ( Element e,
int  refinement,
Element **  sons_out = nullptr 
)

Refines a quad element into four quads, or two quads (horizontally or vertically. If mesh != nullptr, the new_ elements are incorporated into the mesh. The option mesh == nullptr is used to perform adaptive numerical quadrature. If sons_out != nullptr, pointers to the new_ elements will be saved there.

Definition at line 505 of file mesh.cpp.

void Hermes::Hermes2D::Mesh::refine_towards_boundary ( std::vector< std::string >  markers,
int  depth = 1,
bool  aniso = true,
bool  mark_as_initial = false 
)

Performs repeated refinements of elements touching a part of the boundary marked by 'marker'. Elements touching both by an edge or by a vertex are refined. 'aniso' allows or disables anisotropic splits of quads.

Definition at line 845 of file mesh.cpp.

void Hermes::Hermes2D::Mesh::refine_towards_vertex ( int  vertex_id,
int  depth = 1,
bool  mark_as_initial = false 
)

Performs repeated refinements of elements containing the given vertex. A mesh graded towards the vertex is created.

Definition at line 804 of file mesh.cpp.

int * Hermes::Hermes2D::Mesh::regularize ( int  n)

Regularizes the mesh by refining elements with hanging nodes of degree more than 'n'. As a result, n-irregular mesh is obtained. If n = 0, completely regular mesh is created. In this case, however, due to incompatible refinements, the element refinement hierarchy is removed and all elements become top-level elements. Also, total regularization does not work on curved elements. Returns an array of new_ element parents which can be passed to Space::distribute_orders(). The array must be deallocated with free_with_check().

Definition at line 2682 of file mesh.cpp.

void Hermes::Hermes2D::Mesh::unrefine_all_elements ( bool  keep_initial_refinements = true)

Unrefines all elements with immediate active sons. In effect, this shaves off one layer of refinements from the mesh. If done immediately after refine_all_elements(), this function reverts the mesh to its original state. However, it is not exactly an inverse to refine_all_elements().

Definition at line 995 of file mesh.cpp.

void Hermes::Hermes2D::Mesh::unrefine_element_id ( int  id)

Recursively removes all son elements of the given element and makes it active.

Definition at line 981 of file mesh.cpp.

Referenced by unrefine_all_elements().

Member Data Documentation

std::vector<std::pair<unsigned int, int> > Hermes::Hermes2D::Mesh::refinements

Vector for storing refinements in order to be able to save/load meshes with identical element IDs. Refinement "-1" stands for unrefinement.

Definition at line 435 of file mesh.h.

Referenced by copy(), and free().


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