Finite Element Mesh

Mesh management classes More...

Classes

class  Domain
 To store and treat finite element geometric information. More...
 
class  Edge
 To describe an edge. More...
 
class  Element
 To store and treat finite element geometric information. More...
 
class  Figure
 To store and treat a figure (or shape) information. More...
 
class  Rectangle
 To store and treat a rectangular figure. More...
 
class  Brick
 To store and treat a brick (parallelepiped) figure. More...
 
class  Circle
 To store and treat a circular figure. More...
 
class  Sphere
 To store and treat a sphere. More...
 
class  Ellipse
 To store and treat an ellipsoidal figure. More...
 
class  Triangle
 To store and treat a triangle. More...
 
class  Polygon
 To store and treat a polygonal figure. More...
 
class  Grid
 To manipulate structured grids. More...
 
class  Mesh
 To store and manipulate finite element meshes. More...
 
class  MeshAdapt
 To adapt mesh in function of given solution. More...
 
class  NodeList
 Class to construct a list of nodes having some common properties. More...
 
class  ElementList
 Class to construct a list of elements having some common properties. More...
 
class  SideList
 Class to construct a list of sides having some common properties. More...
 
class  EdgeList
 Class to construct a list of edges having some common properties. More...
 
class  Node
 To describe a node. More...
 
class  Partition
 To partition a finite element mesh into balanced submeshes. More...
 
class  Side
 To store and treat finite element sides (edges in 2-D or faces in 3-D) More...
 

Macros

#define GRAPH_MEMORY   1000000
 Memory necessary to store matrix graph.
 
#define MAX_NB_ELEMENTS   10000
 Maximal Number of elements.
 
#define MAX_NB_NODES   10000
 Maximal number of nodes.
 
#define MAX_NB_SIDES   30000
 Maximal number of sides in.
 
#define MAX_NB_EDGES   30000
 Maximal Number of edges.
 
#define MAX_NBDOF_NODE   6
 Maximum number of DOF supported by each node.
 
#define MAX_NBDOF_SIDE   6
 Maximum number of DOF supported by each side.
 
#define MAX_NBDOF_EDGE   2
 Maximum number of DOF supported by each edge.
 
#define MAX_NB_ELEMENT_NODES   20
 Maximum number of nodes by element.
 
#define MAX_NB_ELEMENT_EDGES   10
 Maximum number of edges by element.
 
#define MAX_NB_SIDE_NODES   9
 Maximum number of nodes by side.
 
#define MAX_NB_ELEMENT_SIDES   8
 Maximum number of sides by element.
 
#define MAX_NB_ELEMENT_DOF   27
 Maximum number of dof by element.
 
#define MAX_NB_SIDE_DOF   4
 Maximum number of dof by side.
 
#define MAX_NB_INT_PTS   20
 Maximum number of integration points in element.
 
#define MAX_NB_MATERIALS   10
 Maximum number of materials.
 
#define TheNode   (*theNode)
 
#define TheElement   (*theElement)
 
#define TheSide   (*theSide)
 
#define TheEdge   (*theEdge)
 
#define ElementLoop(m)   for (auto const& theElement: (m).theElements)
 
#define SideLoop(m)   for (auto const& theSide: (m).theSides)
 
#define EdgeLoop(m)   for (auto const& theEdge: (m).theEdges)
 
#define NodeLoop(m)   for (auto const& theNode: (m).theNodes)
 
#define BoundaryNodeLoop(m)   for (auto const& theNode: (m).theBoundaryNodes)
 
#define BoundarySideLoop(m)   for (auto const& theSide: (m).theBoundarySides)
 
#define ElementNodeLoop(el, nd)   for (auto const& nd: (el).theNodes)
 
#define ElementSideLoop(el, sd)   for (auto const& sd: (el).theSides)
 
#define theNodeLabel   theNode->n()
 
#define theSideLabel   theSide->n()
 A macro that returns side label in a loop using macro MeshSides
 
#define theSideNodeLabel(i)   theSide->getNodeLabel(i)
 A macro that returns label of i-th node of side using macro MeshSides
 
#define theElementLabel   theElement->n()
 A macro that returns element label in a loop using macro MeshElements
 
#define theElementNodeLabel(i)   theElement->getNodeLabel(i)
 A macro that returns label of i-th node of element using macro MeshElements
 

Functions

ostream & operator<< (ostream &s, const Edge &ed)
 Output edge data.
 
ostream & operator<< (ostream &s, const Element &el)
 Output element data.
 
Figure operator&& (const Figure &f1, const Figure &f2)
 Function to define a Figure instance as the intersection of two Figure instances.
 
Figure operator|| (const Figure &f1, const Figure &f2)
 Function to define a Figure instance as the union of two Figure instances.
 
Figure operator- (const Figure &f1, const Figure &f2)
 Function to define a Figure instance as the set subtraction of two Figure instances.
 
ostream & operator<< (ostream &s, const Material &m)
 Output material data.
 
ostream & operator<< (ostream &s, const Mesh &ms)
 Output mesh data.
 
ostream & operator<< (ostream &s, const MeshAdapt &a)
 Output MeshAdapt class data.
 
ostream & operator<< (ostream &s, const NodeList &nl)
 Output NodeList instance.
 
ostream & operator<< (ostream &s, const ElementList &el)
 Output ElementList instance.
 
ostream & operator<< (ostream &s, const SideList &sl)
 Output SideList instance.
 
ostream & operator<< (ostream &s, const EdgeList &el)
 Output EdgeList instance.
 
size_t Label (const Node &nd)
 Return label of a given node.
 
size_t Label (const Element &el)
 Return label of a given element.
 
size_t Label (const Side &sd)
 Return label of a given side.
 
size_t Label (const Edge &ed)
 Return label of a given edge.
 
size_t NodeLabel (const Element &el, size_t n)
 Return global label of node local label in element.
 
size_t NodeLabel (const Side &sd, size_t n)
 Return global label of node local label in side.
 
Point< real_t > Coord (const Node &nd)
 Return coordinates of a given node.
 
int Code (const Node &nd, size_t i=1)
 Return code of a given (degree of freedom of) node.
 
int Code (const Element &el)
 Return code of a given element.
 
int Code (const Side &sd, size_t i=1)
 Return code of a given (degree of freedom of) side.
 
bool operator== (const Element &el1, const Element &el2)
 Check equality between 2 elements.
 
bool operator== (const Side &sd1, const Side &sd2)
 Check equality between 2 sides.
 
void DeformMesh (Mesh &mesh, const Vect< real_t > &u, real_t rate=0.2)
 Calculate deformed mesh using a displacement field.
 
void MeshToMesh (Mesh &m1, Mesh &m2, const Vect< real_t > &u1, Vect< real_t > &u2, size_t nx, size_t ny=0, size_t nz=0, size_t dof=1)
 Function to redefine a vector defined on a mesh to a new mesh.
 
void MeshToMesh (const Vect< real_t > &u1, Vect< real_t > &u2, size_t nx, size_t ny=0, size_t nz=0, size_t dof=1)
 Function to redefine a vector defined on a mesh to a new mesh.
 
void MeshToMesh (Mesh &m1, Mesh &m2, const Vect< real_t > &u1, Vect< real_t > &u2, const Point< real_t > &xmin, const Point< real_t > &xmax, size_t nx, size_t ny, size_t nz, size_t dof=1)
 Function to redefine a vector defined on a mesh to a new mesh.
 
real_t getMaxSize (const Mesh &m)
 Return maximal size of element edges for given mesh.
 
real_t getMinSize (const Mesh &m)
 Return minimal size of element edges for given mesh.
 
real_t getMinElementMeasure (const Mesh &m)
 Return minimal measure (length, area or volume) of elements of given mesh.
 
real_t getMaxElementMeasure (const Mesh &m)
 Return maximal measure (length, area or volume) of elements of given mesh.
 
real_t getMinSideMeasure (const Mesh &m)
 Return minimal measure (length or area) of sides of given mesh.
 
real_t getMaxSideMeasure (const Mesh &m)
 Return maximal measure (length or area) of sides of given mesh.
 
real_t getMeanElementMeasure (const Mesh &m)
 Return average measure (length, area or volume) of elements of given mesh.
 
real_t getMeanSideMeasure (const Mesh &m)
 Return average measure (length or area) of sides of given mesh.
 
void setNodeCodes (Mesh &m, const string &exp, int code, size_t dof=1)
 Assign a given code to all nodes satisfying a boolean expression using node coordinates.
 
void setBoundaryNodeCodes (Mesh &m, const string &exp, int code, size_t dof=1)
 Assign a given code to all nodes on boundary that satisfy a boolean expression using node coordinates.
 
int NodeInElement (const Node *nd, const Element *el)
 Say if a given node belongs to a given element.
 
int NodeInSide (const Node *nd, const Side *sd)
 Say if a given node belongs to a given side.
 
int SideInElement (const Side *sd, const Element *el)
 Say if a given side belongs to a given element.
 
ostream & operator<< (ostream &s, const Node &nd)
 Output node data.
 
ostream & operator<< (ostream &s, const Side &sd)
 Output side data.
 

Detailed Description

Mesh management classes

Macro Definition Documentation

◆ BoundaryNodeLoop

#define BoundaryNodeLoop (   m)    for (auto const& theNode: (m).theBoundaryNodes)

A macro to loop on mesh nodes m: Instance of Mesh

Note
: Each iteration updates the pointer theNode to current Node

◆ BoundarySideLoop

#define BoundarySideLoop (   m)    for (auto const& theSide: (m).theBoundarySides)

A macro to loop on mesh boundary sides m: Instance of Mesh

Note
: Each iteration updates the pointer theSide to current Node

◆ EdgeLoop

#define EdgeLoop (   m)    for (auto const& theEdge: (m).theEdges)

A macro to loop on mesh edges m : Instance of Mesh

Note
: Each iteration updates the pointer theEdge to current Edge

◆ ElementLoop

#define ElementLoop (   m)    for (auto const& theElement: (m).theElements)

A macro to loop on mesh elements m : Instance of Mesh

Note
: Each iteration updates the pointer theElement to current Element

◆ ElementNodeLoop

#define ElementNodeLoop (   el,
  nd 
)    for (auto const& nd: (el).theNodes)

A macro to loop on element nodes

Parameters
elInstance of Element
ndPointer to pointed node

◆ ElementSideLoop

#define ElementSideLoop (   el,
  sd 
)    for (auto const& sd: (el).theSides)

A macro to loop on element sides

Parameters
elInstance of Element
sdPointer to pointed side

◆ GRAPH_MEMORY

#define GRAPH_MEMORY   1000000

Memory necessary to store matrix graph.

This value is necessary only if nodes are to be renumbered.

◆ NodeLoop

#define NodeLoop (   m)    for (auto const& theNode: (m).theNodes)

A macro to loop on mesh nodes m : Instance of Mesh

Note
: Each iteration updates the pointer theNode to current Node

◆ SideLoop

#define SideLoop (   m)    for (auto const& theSide: (m).theSides)

A macro to loop on mesh sides m : Instance of Mesh

Note
: Each iteration updates the pointer theSide to current Element

◆ TheEdge

#define TheEdge   (*theEdge)

A macro that gives the instance pointed by theEdge

◆ TheElement

#define TheElement   (*theElement)

A macro that gives the instance pointed by theElement

◆ TheNode

#define TheNode   (*theNode)

A macro that gives the instance pointed by theNode

◆ theNodeLabel

#define theNodeLabel   theNode->n()

A macro that returns node label in a loop using macro MeshNodes

◆ TheSide

#define TheSide   (*theSide)

A macro that gives the instance pointed by theSide

Function Documentation

◆ Code() [1/3]

int Code ( const Element el)

Return code of a given element.

Parameters
[in]elReference to Element instance
Returns
Code of element
Author
Rachid Touzani

◆ Code() [2/3]

int Code ( const Node nd,
size_t  i = 1 
)

Return code of a given (degree of freedom of) node.

Parameters
[in]ndReference to Node instance
[in]iLabel of dof [Default: 1]
Returns
Code of dof of node
Author
Rachid Touzani

◆ Code() [3/3]

int Code ( const Side sd,
size_t  i = 1 
)

Return code of a given (degree of freedom of) side.

Parameters
[in]sdReference to Side instance
[in]iLabel of dof [Default: 1]
Returns
Code of dof of side
Author
Rachid Touzani

◆ Coord()

Point< real_t > Coord ( const Node nd)

Return coordinates of a given node.

Parameters
[in]ndReference to Node instance
Returns
Coordinates of node
Author
Rachid Touzani

◆ DeformMesh()

void DeformMesh ( Mesh mesh,
const Vect< real_t > &  u,
real_t  a = 0.2 
)

Calculate deformed mesh using a displacement field.

Parameters
[in,out]meshMesh instance. On output, node coordinates are modified to take into account the displacement
[in]uDisplacement field at nodes
[in]aMaximal deformation rate. [Default: 1]. A typical value is 0.2 (i.e. 20%).
Author
Rachid Touzani

◆ getMaxElementMeasure()

real_t getMaxElementMeasure ( const Mesh m)

Return maximal measure (length, area or volume) of elements of given mesh.

Parameters
[in]mReference to mesh instance
Author
Rachid Touzani

◆ getMaxSideMeasure()

real_t getMaxSideMeasure ( const Mesh m)

Return maximal measure (length or area) of sides of given mesh.

Parameters
[in]mReference to mesh instance
Note
Use this function only if sides are present in the mesh and for 2-D meshes
Author
Rachid Touzani

◆ getMaxSize()

real_t getMaxSize ( const Mesh m)

Return maximal size of element edges for given mesh.

Parameters
[in]mReference to mesh instance
Author
Rachid Touzani

◆ getMeanElementMeasure()

real_t getMeanElementMeasure ( const Mesh m)

Return average measure (length, area or volume) of elements of given mesh.

Parameters
[in]mReference to mesh instance
Author
Rachid Touzani

◆ getMeanSideMeasure()

real_t getMeanSideMeasure ( const Mesh m)

Return average measure (length or area) of sides of given mesh.

Parameters
[in]mReference to mesh instance
Note
Use this function only if sides are present in the mesh and for 2-D meshes
Author
Rachid Touzani

◆ getMinElementMeasure()

real_t getMinElementMeasure ( const Mesh m)

Return minimal measure (length, area or volume) of elements of given mesh.

Parameters
[in]mReference to mesh instance
Author
Rachid Touzani

◆ getMinSideMeasure()

real_t getMinSideMeasure ( const Mesh m)

Return minimal measure (length or area) of sides of given mesh.

Parameters
[in]mReference to mesh instance
Note
Use this function only if sides are present in the mesh and for 2-D meshes
Author
Rachid Touzani

◆ getMinSize()

real_t getMinSize ( const Mesh m)

Return minimal size of element edges for given mesh.

Parameters
[in]mReference to mesh instance
Author
Rachid Touzani

◆ Label() [1/4]

size_t Label ( const Edge ed)

Return label of a given edge.

Parameters
[in]edReference to Edge instance
Returns
Label of edge
Author
Rachid Touzani

◆ Label() [2/4]

size_t Label ( const Element el)

Return label of a given element.

Parameters
[in]elReference to Element instance
Returns
Label of element
Author
Rachid Touzani

◆ Label() [3/4]

size_t Label ( const Node nd)

Return label of a given node.

Parameters
[in]ndReference to Node instance
Returns
Label of node
Author
Rachid Touzani

◆ Label() [4/4]

size_t Label ( const Side sd)

Return label of a given side.

Parameters
[in]sdReference to Side instance
Returns
Label of side
Author
Rachid Touzani

◆ MeshToMesh() [1/3]

void MeshToMesh ( const Vect< real_t > &  u1,
Vect< real_t > &  u2,
size_t  nx,
size_t  ny = 0,
size_t  nz = 0,
size_t  dof = 1 
)

Function to redefine a vector defined on a mesh to a new mesh.

The program interpolates (piecewise linear) first the vector on a finer structured grid. Then the values on the new mesh nodes are computed.

Remarks
For efficiency the number of grid cells must be large enough so that interpolation provides efficient accuracy
Parameters
[in]u1Input vector of nodal values defined on first mesh. This vector instance must contain Mesh instance
[out]u2Output vector of nodal values defined on second mesh. This vector instance must contain Mesh instance
[in]nxNumber of cells in the x-direction in the fine structured grid
[in]nyNumber of cells in the y-direction in the fine structured grid The default value of ny is 0, i.e. a 1-D grid
[in]nzNumber of cells in the z-direction in the fine structured grid The default value of nz is 0, i.e. a 1-D or 2-D grid
[in]dofLabel of degree of freedom of vector u. Only this dof is considered. [Default: 1]
Note
The input vector u1 is a one degree of freedom per node vector, i.e. its size must be equal (or greater than) the total number of nodes of mesh m1. The size of vector u2 is deduced from the mesh m2
Author
Rachid Touzani

◆ MeshToMesh() [2/3]

void MeshToMesh ( Mesh m1,
Mesh m2,
const Vect< real_t > &  u1,
Vect< real_t > &  u2,
const Point< real_t > &  xmin,
const Point< real_t > &  xmax,
size_t  nx,
size_t  ny,
size_t  nz,
size_t  dof = 1 
)

Function to redefine a vector defined on a mesh to a new mesh.

The program interpolates (piecewise linear) first the vector on a finer structured grid. Then the values on the new mesh nodes are computed. In this function the grid rectangle is defined so that this one can cover only a submesh of m1.

Remarks
For efficiency the number of grid cells must be large enough so that interpolation provides efficient accuracy
Parameters
[in]m1Reference to the first mesh instance
[out]m2Reference to the second mesh instance
[in]u1Input vector of nodal values defined on first mesh
[out]u2Output vector of nodal values defined on second mesh
[in]xminPoint instance containing minimal coordinates of the rectangle that defines the grid
[in]xmaxPoint instance containing maximal coordinates of the rectangle that defines the grid
[in]nxNumber of cells in the x-direction in the fine structured grid
[in]nyNumber of cells in the y-direction in the fine structured grid The default value of ny is 0, i.e. a 1-D grid
[in]nzNumber of cells in the z-direction in the fine structured grid The default value of nz is 0, i.e. a 1-D or 2-D grid
[in]dofLabel of degree of freedom of vector u. Only this dof is considered. [Default: 1]
Note
The input vector u1 is a one degree of freedom per node vector, i.e. its size must be equal (or greater than) the total number of nodes of mesh m1. The size of vector u2 is deduced from the mesh m2
Author
Rachid Touzani

◆ MeshToMesh() [3/3]

void MeshToMesh ( Mesh m1,
Mesh m2,
const Vect< real_t > &  u1,
Vect< real_t > &  u2,
size_t  nx,
size_t  ny = 0,
size_t  nz = 0,
size_t  dof = 1 
)

Function to redefine a vector defined on a mesh to a new mesh.

The program interpolates (piecewise linear) first the vector on a finer structured grid. Then the values on the new mesh nodes are computed.

Remarks
For efficiency the number of grid cells must be large enough so that interpolation provides efficient accuracy
Parameters
[in]m1Reference to the first mesh instance
[out]m2Reference to the second mesh instance
[in]u1Input vector of nodal values defined on first mesh
[out]u2Output vector of nodal values defined on second mesh
[in]nxNumber of cells in the x-direction in the fine structured grid
[in]nyNumber of cells in the y-direction in the fine structured grid The default value of ny is 0, i.e. a 1-D grid
[in]nzNumber of cells in the z-direction in the fine structured grid The default value of nz is 0, i.e. a 1-D or 2-D grid
[in]dofLabel of degree of freedom of vector u. Only this dof is considered. [Default: 1]
Note
The input vector u1 is a one degree of freedom per node vector, i.e. its size must be equal (or greater than) the total number of nodes of mesh m1. The size of vector u2 is deduced from the mesh m2
Author
Rachid Touzani

◆ NodeInElement()

int NodeInElement ( const Node nd,
const Element el 
)

Say if a given node belongs to a given element.

Parameters
[in]ndPointer to Node
[in]elPointer to Element
Returns
Local label of the node if this one is found, 0 if not.
Author
Rachid Touzani

◆ NodeInSide()

int NodeInSide ( const Node nd,
const Side sd 
)

Say if a given node belongs to a given side.

Parameters
[in]ndPointer to Node
[in]sdPointer to Side
Returns
Local label of the node if this one is found, 0 if not.
Author
Rachid Touzani

◆ NodeLabel() [1/2]

size_t NodeLabel ( const Element el,
size_t  n 
)

Return global label of node local label in element.

Parameters
[in]elReference to Element instance
[in]nLocal label of node in element
Returns
Global label of node
Author
Rachid Touzani

◆ NodeLabel() [2/2]

size_t NodeLabel ( const Side sd,
size_t  n 
)

Return global label of node local label in side.

Parameters
[in]sdReference to Side instance
[in]nLocal label of node in side
Returns
Global label of node
Author
Rachid Touzani

◆ operator&&()

Figure operator&& ( const Figure f1,
const Figure f2 
)

Function to define a Figure instance as the intersection of two Figure instances.

Parameters
[in]f1First Figure instance
[in]f2Second Figure instance
Returns
Updated resulting Figure instance

◆ operator-()

Figure operator- ( const Figure f1,
const Figure f2 
)

Function to define a Figure instance as the set subtraction of two Figure instances.

Parameters
[in]f1First Figure instance to subtract from
[in]f2Second Figure instance to subtract
Returns
Updated resulting Figure instance

◆ operator<<() [1/7]

ostream & operator<< ( ostream &  s,
const EdgeList el 
)

Output EdgeList instance.

Author
Rachid Touzani

◆ operator<<() [2/7]

ostream & operator<< ( ostream &  s,
const ElementList el 
)

Output ElementList instance.

Author
Rachid Touzani

◆ operator<<() [3/7]

ostream & operator<< ( ostream &  s,
const Mesh ms 
)

Output mesh data.

Level of mesh output depends on the global variable Verbosity

  • If Verbosity=0 or Verbosity=1, this function outputs only principal mesh parameters: number of nodes, number of elements, ...
  • If Verbosity>1, this function outputs in addition the list of 10 first nodes, elements and sides
  • If Verbosity>2, this function outputs in addition the list of 50 first nodes, elements and sides
  • If Verbosity>3, this function outputs all mesh data

◆ operator<<() [4/7]

ostream & operator<< ( ostream &  s,
const Node nd 
)

Output node data.

Author
Rachid Touzani

◆ operator<<() [5/7]

ostream & operator<< ( ostream &  s,
const NodeList nl 
)

Output NodeList instance.

Author
Rachid Touzani

◆ operator<<() [6/7]

ostream & operator<< ( ostream &  s,
const Side sd 
)

Output side data.

Author
Rachid Touzani

◆ operator<<() [7/7]

ostream & operator<< ( ostream &  s,
const SideList sl 
)

Output SideList instance.

Author
Rachid Touzani

◆ operator==() [1/2]

operator== ( const Element el1,
const Element el2 
)

Check equality between 2 elements.

Parameters
[in]el1Reference to first Side instance
[in]el2Reference to second Side instance
Returns
true is elements are equal, i.e. if they have the same nodes, false if not.
Author
Rachid Touzani

◆ operator==() [2/2]

bool operator== ( const Side sd1,
const Side sd2 
)

Check equality between 2 sides.

Parameters
[in]sd1Reference to first Side instance
[in]sd2Reference to second Side instance
Returns
true is sides are equal, i.e. if they have the same nodes, false if not.
Author
Rachid Touzani

◆ operator||()

Figure operator|| ( const Figure f1,
const Figure f2 
)

Function to define a Figure instance as the union of two Figure instances.

Parameters
[in]f1First Figure instance
[in]f2Second Figure instance
Returns
Updated resulting Figure instance

◆ setBoundaryNodeCodes()

void setBoundaryNodeCodes ( Mesh m,
const string &  exp,
int  code,
size_t  dof = 1 
)

Assign a given code to all nodes on boundary that satisfy a boolean expression using node coordinates.

Parameters
[in]mReference to mesh instance
[in]expRegular expression using x, y, and z coordinates of nodes, according to exprtk parser
[in]codeCode to assign
[in]dofDegree of freedom for which code is assigned [Default: 1]
Author
Rachid Touzani

◆ setNodeCodes()

void setNodeCodes ( Mesh m,
const string &  exp,
int  code,
size_t  dof = 1 
)

Assign a given code to all nodes satisfying a boolean expression using node coordinates.

Parameters
[in]mReference to mesh instance
[in]expRegular expression using x, y, and z coordinates of nodes, according to exprtk parser
[in]codeCode to assign
[in]dofDegree of freedom for which code is assigned [Default: 1]

◆ SideInElement()

int SideInElement ( const Side sd,
const Element el 
)

Say if a given side belongs to a given element.

Parameters
[in]sdPointer to Side
[in]elPointer to Element
Returns
Local label of the side if this one is found, 0 if not.
Author
Rachid Touzani