To store and manipulate finite element meshes. More...

#include <Mesh.h>

Public Member Functions

 Mesh ()
 Default constructor (Empty mesh)
 
 Mesh (const string &file, bool bc=false, int opt=NODE_DOF, int nb_dof=1)
 Constructor using a mesh file.
 
 Mesh (real_t xmin, real_t xmax, size_t nb_el, bool bc=false, size_t p=1, size_t nb_dof=1, int c1=0, int c2=0)
 Constructor for a 1-D mesh. The domain is the interval [xmin,xmax].
 
 Mesh (const Grid &g, int opt=QUADRILATERAL)
 Constructor for a uniform finite difference grid given by and instance of class Grid.
 
 Mesh (const Grid &g, int shape, int opt)
 Constructor of dual mesh for a uniform finite difference grid given by and instance of class Grid.
 
 Mesh (real_t xmin, real_t xmax, size_t ne, int c1, int c2, int p=1, size_t nb_dof=1)
 Constructor for a uniform 1-D finite element mesh.
 
 Mesh (real_t xmin, real_t xmax, real_t ymin, real_t ymax, size_t nx, size_t ny, int cx0, int cxN, int cy0, int cyN, int opt=0, size_t nb_dof=1)
 Constructor for a uniform 2-D structured finite element mesh.
 
 Mesh (real_t xmin, real_t xmax, real_t ymin, real_t ymax, real_t zmin, real_t zmax, size_t nx, size_t ny, size_t nz, int cx0, int cxN, int cy0, int cyN, int cz0, int czN, int opt=0, size_t nb_dof=1)
 Constructor for a uniform 3-D structured finite element mesh.
 
 Mesh (const Mesh &m, const Point< real_t > &x_bl, const Point< real_t > &x_tr)
 Constructor that extracts the mesh of a rectangular region from an initial mesh.
 
 Mesh (const Mesh &mesh, int opt, size_t dof1, size_t dof2, bool bc=false)
 Constructor that copies the input mesh and selects given degrees of freedom.
 
 Mesh (const Mesh &ms)
 Copy Constructor.
 
 ~Mesh ()
 Destructor.
 
void setDim (size_t dim)
 Define space dimension. Normally, between 1 and 3.
 
void Add (Node *nd)
 Add a node to mesh.
 
void Add (Element *el)
 Add an element to mesh.
 
void Add (Side *sd)
 Add a side to mesh.
 
void Add (Edge *ed)
 Add an edge to mesh.
 
Meshoperator*= (real_t a)
 Operator *=
 
void get (const string &mesh_file)
 Read mesh data in file.
 
void get (const string &mesh_file, int ff, int nb_dof=1)
 Read mesh data in file with giving its format.
 
void setDOFSupport (int opt, int nb_nodes=1)
 Define supports of degrees of freedom.
 
void setNbDOFPerNode (size_t nb_dof=1)
 Define number of degrees of freedom for each node.
 
void setPointInDomain (Point< real_t > x)
 Define a point in the domain. This function makes sense only if boundary mesh is given without internal mesh (Case of Boundary Elements)
 
void removeImposedDOF ()
 Eliminate equations corresponding to imposed DOF.
 
size_t NumberEquations (size_t dof=0)
 Renumber Equations.
 
size_t NumberEquations (size_t dof, int c)
 Renumber Equations.
 
int getAllSides (int opt=0)
 Determine all mesh sides.
 
size_t getNbSideNodes () const
 Return the number of nodes on each side.
 
size_t getNbElementNodes () const
 Return the number of nodes in each element.
 
int getBoundarySides ()
 Determine all boundary sides.
 
int createBoundarySideList ()
 Create list of boundary sides.
 
int getBoundaryNodes ()
 Determine all boundary nodes.
 
int createInternalSideList ()
 Create list of internal sides (not on the boundary).
 
int getAllEdges ()
 Determine all edges.
 
void getNodeNeighborElements ()
 Create node neighboring elements.
 
void getElementNeighborElements ()
 Create element neighboring elements.
 
void setMaterial (int code, const string &mname)
 Associate material to code of element.
 
void Reorder (size_t m=GRAPH_MEMORY)
 Renumber mesh nodes according to reverse Cuthill Mc Kee algorithm.
 
void Add (size_t num, real_t *x)
 Add a node by giving its label and an array containing its coordinates.
 
void DeleteNode (size_t label)
 Remove a node given by its label.
 
void DeleteElement (size_t label)
 Remove an element given by its label.
 
void DeleteSide (size_t label)
 Remove a side given by its label.
 
void Delete (Node *nd)
 Remove a node given by its pointer.
 
void Delete (Element *el)
 Remove a node given by its pointer.
 
void Delete (Side *sd)
 Remove a side given by its pointer.
 
void Delete (Edge *ed)
 Remove an edge given by its pointer.
 
void RenumberNode (size_t n1, size_t n2)
 Renumber a node.
 
void RenumberElement (size_t n1, size_t n2)
 Renumber an element.
 
void RenumberSide (size_t n1, size_t n2)
 Renumber a side.
 
void RenumberEdge (size_t n1, size_t n2)
 Renumber an edge.
 
void setList (const vector< Node * > &nl)
 Initialize list of mesh nodes using the input vector.
 
void setList (const vector< Element * > &el)
 Initialize list of mesh elements using the input vector.
 
void setList (const vector< Side * > &sl)
 Initialize list of mesh sides using the input vector.
 
void Rescale (real_t sx, real_t sy=0., real_t sz=0.)
 Rescale mesh by multiplying node coordinates by constants.
 
size_t getDim () const
 Return space dimension.
 
size_t getNbNodes () const
 Return number of nodes.
 
size_t getNbMarkedNodes () const
 Return number of marked nodes.
 
size_t getNbVertices () const
 Return number of vertices.
 
size_t getNbDOF () const
 Return total number of degrees of freedom (DOF)
 
size_t getNbEq () const
 Return number of equations.
 
size_t getNbEq (int i) const
 Return number of equations for the i-th set of degrees of freedom.
 
size_t getNbElements () const
 Return number of elements.
 
size_t getNbSides () const
 Return number of sides.
 
size_t getNbEdges () const
 Return number of sides.
 
size_t getNbBoundarySides () const
 Return number of boundary sides.
 
size_t getNbInternalSides () const
 Return number of internal sides.
 
size_t getNbMat () const
 Return number of materials.
 
void AddMidNodes (int g=0)
 Add mid-side nodes.
 
Point< real_t > getMaxCoord () const
 Return maximum coordinates of nodes.
 
Point< real_t > getMinCoord () const
 Return minimum coordinates of nodes.
 
void set (Node *nd)
 Replace node in the mesh.
 
void set (Element *el)
 Replace element in the mesh.
 
void set (Side *sd)
 Choose side in the mesh.
 
bool NodesAreDOF () const
 Return information about DOF type.
 
bool SidesAreDOF () const
 Return information about DOF type.
 
bool EdgesAreDOF () const
 Return information about DOF type.
 
bool ElementsAreDOF () const
 Return information about DOF type.
 
int getDOFSupport () const
 Return information on dof support Return an integer according to enumerated values: NODE_DOF, ELEMENT_DOF SIDE_DOF.
 
void Deform (const Vect< real_t > &u, real_t rate=0.2)
 Deform mesh according to a displacement vector.
 
void put (const string &mesh_file) const
 Write mesh data on file.
 
void save (const string &mesh_file) const
 Write mesh data on file in various formats.
 
bool withImposedDOF () const
 Return true if imposed DOF count in equations, false if not.
 
bool isStructured () const
 Return true is mesh is structured, false if not.
 
size_t getNodeNewLabel (size_t n) const
 Return new label of node of a renumbered node.
 
void getList (vector< Node * > &nl) const
 Fill vector nl with list of pointers to nodes.
 
void getList (vector< Element * > &el) const
 Fill vector el with list of pointers to elements.
 
void getList (vector< Side * > &sl) const
 Fill vector sl with list of pointers to sides.
 
NodegetPtrNode (size_t i) const
 Return pointer to node with label i.
 
NodegetNode (size_t i) const
 Return refenrece to node with label i
 
ElementgetPtrElement (size_t i) const
 Return pointer to element with label i
 
SidegetPtrSide (size_t i) const
 Return pointer to side with label i
 
SidegetSide (size_t i) const
 Return reference to side with label i
 
EdgegetPtrEdge (size_t i) const
 Return pointer to edge with label i
 
size_t getNodeLabel (size_t i) const
 Return label of i-th node.
 
size_t getElementLabel (size_t i) const
 Return label of i-th element.
 
size_t getSideLabel (size_t i) const
 Return label of i-th side.
 
size_t getEdgeLabel (size_t i) const
 Return label of i-th edge.
 
int getShape () const
 Determine shape of elements Return Shape index (see enum ElementShape) if all elements have the same shape, 0 if not.
 
Elementoperator() (size_t i) const
 Operator () : Return pointer to i-th element.
 
Nodeoperator[] (size_t i) const
 Operator [] : Return pointer to i-th node.
 
size_t operator() (size_t i, size_t n) const
 Operator () : Return pointer to i-th node of n-th element.
 
Meshoperator= (Mesh &ms)
 Operator = : Assign a Mesh instance.
 

Detailed Description

To store and manipulate finite element meshes.

Class Mesh enables defining as an object a finite element mesh. A finite element mesh is characterized by its nodes, elements and sides. Each of these types of data constitutes a class in the OFELI library.

The standard procedure to introduce the finite element mesh is to provide an input file containing its data. For this, we have defined our own mesh data file (following the XML syntax). Of course, a developer can write his own function to read his finite element mesh file using the methods in Mesh.

Author
Rachid Touzani

Constructor & Destructor Documentation

◆ Mesh() [1/10]

Mesh ( const string &  file,
bool  bc = false,
int  opt = NODE_DOF,
int  nb_dof = 1 
)

Constructor using a mesh file.

Parameters
[in]fileFile containing mesh data. The extension of the file yields the file format: The extension .m implies OFELI file format and .msh implies GMSH msh file.
[in]bcFlag to remove (true) or not (false) imposed Degrees of Freedom [default: false]
[in]optType of DOF support: To choose among enumerated values NODE_DOF, SIDE_DOF or ELEMENT_DOF.
Say if degrees of freedom (unknowns) are supported by nodes, sides or elements.
[in]nb_dofNumber of degrees of freedom per node [Default: 1].

◆ Mesh() [2/10]

Mesh ( real_t  xmin,
real_t  xmax,
size_t  nb_el,
bool  bc = false,
size_t  p = 1,
size_t  nb_dof = 1,
int  c1 = 0,
int  c2 = 0 
)

Constructor for a 1-D mesh. The domain is the interval [xmin,xmax].

Parameters
[in]xminValue of xmin
[in]xmaxValue of xmax
[in]nb_elNumber of elements to generate
[in]bcFlag to remove (true) or not (false) imposed Degrees of Freedom [default: false]
[in]pDegree of finite element polynomial [default: 1]
[in]nb_dofNumber of degrees of freedom for each node [default: 1]
[in]c1Code to assign to first node (at x=xmin) [default: 0]
[in]c2Code to assign to last node (at x=xmax) [default: 0]

◆ Mesh() [3/10]

Mesh ( const Grid g,
int  opt = QUADRILATERAL 
)

Constructor for a uniform finite difference grid given by and instance of class Grid.

Parameters
[in]gGrid instance
[in]optOptional value to say which type of elements to generate
  • TRIANGLE: Mesh elements are triangles
  • QUADRILATERAL: Mesh elements are quadrilaterals [default]

◆ Mesh() [4/10]

Mesh ( const Grid g,
int  shape,
int  opt 
)

Constructor of dual mesh for a uniform finite difference grid given by and instance of class Grid.

Parameters
[in]gGrid instance
[in]shapeValue to say which type of elements to generate
  • TRIANGLE: Mesh elements are triangles
  • QUADRILATERAL: Mesh elements are quadrilaterals [default]
[in]optThis argument can take any value. It is here only to distinguish from the other constructor using Grid instance.
Remarks
This constructor is to be used to obtain a dual mesh from a structured grid. It is mainly useful if a cell centered finite volume method is used.

◆ Mesh() [5/10]

Mesh ( real_t  xmin,
real_t  xmax,
size_t  ne,
int  c1,
int  c2,
int  p = 1,
size_t  nb_dof = 1 
)

Constructor for a uniform 1-D finite element mesh.

The domain is the line (xmin,xmax)

Parameters
[in]xminMinimal coordinate
[in]xmaxMaximal coordinate
[in]neNumber of elements
[in]c1Code for the first node (x=xmin)
[in]c2Code for the last node (x=xmax)
[in]pDegree of approximation polynomial [Default: 1].
[in]nb_dofNumber of degrees of freedom per node [Default: 1].
Remarks
The option p can be set to 1 if the user intends to use finite differences.

◆ Mesh() [6/10]

Mesh ( real_t  xmin,
real_t  xmax,
real_t  ymin,
real_t  ymax,
size_t  nx,
size_t  ny,
int  cx0,
int  cxN,
int  cy0,
int  cyN,
int  opt = 0,
size_t  nb_dof = 1 
)

Constructor for a uniform 2-D structured finite element mesh.

The domain is the rectangle (xmin,xmax)x(ymin,ymax)

Parameters
[in]xminMinimal x-coordinate
[in]xmaxMaximal x-coordinate
[in]yminMinimal y-coordinate
[in]ymaxMaximal y-coordinate
[in]nxNumber of subintervals on the x-axis
[in]nyNumber of subintervals on the y-axis
[in]cx0Code for nodes generated on the line x=x0 if >0, for sides on this line if <0
[in]cxNCode for nodes generated on the line x=xN if >0, for sides on this line if <0
[in]cy0Code for nodes generated on the line y=y0 if >0, for sides on this line if <0
[in]cyNCode for nodes generated on the line y=yN if >0, for sides on this line if <0
[in]optFlag to generate elements as well (if not zero) [Default: 0]. If the flag is not 0, it can take one of the enumerated values: TRIANGLE or QUADRILATERAL, with obvious meaning.
[in]nb_dofNumber of degrees of freedom per node [Default: 1].
Remarks
The option opt can be set to 0 if the user intends to use finite differences.

◆ Mesh() [7/10]

Mesh ( real_t  xmin,
real_t  xmax,
real_t  ymin,
real_t  ymax,
real_t  zmin,
real_t  zmax,
size_t  nx,
size_t  ny,
size_t  nz,
int  cx0,
int  cxN,
int  cy0,
int  cyN,
int  cz0,
int  czN,
int  opt = 0,
size_t  nb_dof = 1 
)

Constructor for a uniform 3-D structured finite element mesh.

The domain is the parallepiped (xmin,xmax)x(ymin,ymax)x(zmin,zmax)

Parameters
[in]xminMinimal x-coordinate
[in]xmaxMaximal x-coordinate
[in]yminMinimal y-coordinate
[in]ymaxMaximal y-coordinate
[in]zminMinimal z-coordinate
[in]zmaxMaximal z-coordinate
[in]nxNumber of subintervals on the x-axis
[in]nyNumber of subintervals on the y-axis
[in]nzNumber of subintervals on the z-axis
[in]cx0Code for nodes generated on the line x=xmin if >0, for sides on this line if <0
[in]cxNCode for nodes generated on the line x=xmax if >0, for sides on this line if <0
[in]cy0Code for nodes generated on the line y=ymin if >0, for sides on this line if <0
[in]cyNCode for nodes generated on the line y=ymax if >0, for sides on this line if <0
[in]cz0Code for nodes generated on the line z=zmin if >0, for sides on this line if <0
[in]czNCode for nodes generated on the line z=zmax if >0, for sides on this line if <0
[in]optFlag to generate elements as well (if not zero) [Default: 0]. If the flag is not 0, it can take one of the enumerated values: HEXAHEDRON or TETRAHEDRON, with obvious meaning.
[in]nb_dofNumber of degrees of freedom per node [Default: 1].
Remarks
The option opt can be set to 0 if the user intends to use finite differences.

◆ Mesh() [8/10]

Mesh ( const Mesh m,
const Point< real_t > &  x_bl,
const Point< real_t > &  x_tr 
)

Constructor that extracts the mesh of a rectangular region from an initial mesh.

This constructor is useful for zooming purposes for instance.

Parameters
[in]mInitial mesh from which the submesh is extracted
[in]x_blCoordinate of bottom left vertex of the rectangle
[in]x_trCoordinate of top right vertex of the rectangle

◆ Mesh() [9/10]

Mesh ( const Mesh mesh,
int  opt,
size_t  dof1,
size_t  dof2,
bool  bc = false 
)

Constructor that copies the input mesh and selects given degrees of freedom.

This constructor is to be used for coupled problems where each subproblem uses a choice of degrees of freedom.

Parameters
[in]meshInitial mesh from which the submesh is extracted
[in]optType of DOF support: To choose among enumerated values NODE_DOF, SIDE_DOF or ELEMENT_DOF.
[in]dof1Label of first degree of freedom to select to the output mesh
[in]dof2Label of last degree of freedom to select to the output mesh
[in]bcFlag to remove (true) or not (false) imposed Degrees of Freedom [Default: false]

◆ Mesh() [10/10]

Mesh ( const Mesh ms)

Copy Constructor.

Parameters
[in]msMesh instance to copy

Member Function Documentation

◆ Add() [1/5]

void Add ( Edge ed)

Add an edge to mesh.

Parameters
[in]edPointer to Edge to add

◆ Add() [2/5]

void Add ( Element el)

Add an element to mesh.

Parameters
[in]elPointer to Element to add

◆ Add() [3/5]

void Add ( Node nd)

Add a node to mesh.

Parameters
[in]ndPointer to Node to add

◆ Add() [4/5]

void Add ( Side sd)

Add a side to mesh.

Parameters
[in]sdPointer to Side to add

◆ Add() [5/5]

void Add ( size_t  num,
real_t *  x 
)

Add a node by giving its label and an array containing its coordinates.

Parameters
[in]numLabel of node to add
[in]xC-array of node coordinates

◆ AddMidNodes()

void AddMidNodes ( int  g = 0)

Add mid-side nodes.

This is function is valid for triangles only

Parameters
[in]gOption to say of barycentre node is to be added (>0) or not (=0)

◆ createBoundarySideList()

int createBoundarySideList ( )

Create list of boundary sides.

This function is useful to loop over boundary sides without testing Once this one is called, the function getNbBoundarySides() is available. Moreover, looping over boundary sides is available via the member functions topBoundarySide() and getBoundarySide()

Returns
Number of boundary sides.

◆ createInternalSideList()

int createInternalSideList ( )

Create list of internal sides (not on the boundary).

This function is useful to loop over internal sides without testing Once this one is called, the function getNbInternalSides() is available. Moreover, looping over internal sides is available via the member functions topInternalSide() and getInternalSide()

Returns
n Number of internal sides.

◆ Deform()

void Deform ( const Vect< real_t > &  u,
real_t  rate = 0.2 
)

Deform mesh according to a displacement vector.

This function modifies node coordinates according to given displacement vector and given rate

Parameters
[in]uDisplacement vector
[in]rateMaximal rate of deformation of resulting mesh. Its default value is 0.2, i.e. The resulting mesh has a maximum of deformation rate of 20%

◆ Delete() [1/4]

void Delete ( Edge ed)

Remove an edge given by its pointer.

This function does not release the space previously occupied

Parameters
[in]edPointer to edge to delete

◆ Delete() [2/4]

void Delete ( Element el)

Remove a node given by its pointer.

This function does not release the space previously occupied

Parameters
[in]elPointer to element to delete

◆ Delete() [3/4]

void Delete ( Node nd)

Remove a node given by its pointer.

This function does not release the space previously occupied

Parameters
[in]ndPointer to node to delete

◆ Delete() [4/4]

void Delete ( Side sd)

Remove a side given by its pointer.

This function does not release the space previously occupied

Parameters
[in]sdPointer to side to delete

◆ DeleteElement()

void DeleteElement ( size_t  label)

Remove an element given by its label.

This function does not release the space previously occupied

Parameters
[in]labelLabel of element to delete

◆ DeleteNode()

void DeleteNode ( size_t  label)

Remove a node given by its label.

This function does not release the space previously occupied

Parameters
[in]labelLabel of node to delete

◆ DeleteSide()

void DeleteSide ( size_t  label)

Remove a side given by its label.

This function does not release the space previously occupied

Parameters
[in]labelLabel of side to delete

◆ EdgesAreDOF()

bool EdgesAreDOF ( ) const

Return information about DOF type.

Returns
true if DOF are supported by edges, false otherwise

◆ ElementsAreDOF()

bool ElementsAreDOF ( ) const

Return information about DOF type.

Returns
true if DOF are supported by elements, false otherwise

◆ get() [1/2]

void get ( const string &  mesh_file)

Read mesh data in file.

Mesh file must be in OFELI format. See "File Formats" page

Parameters
[in]mesh_fileMesh file name

◆ get() [2/2]

void get ( const string &  mesh_file,
int  ff,
int  nb_dof = 1 
)

Read mesh data in file with giving its format.

File format can be chosen among a variety of choices. See "File Formats" page

Parameters
[in]mesh_fileMesh file name
[in]ffFile format: Integer to chose among enumerated values: OFELI_FF, GMSH, MATLAB, EASYMESH, GAMBIT, BAMG, NETGEN, TRIANGLE_FF
[in]nb_dofNumber of degrees of freedom per node (Default value: 1)

◆ getAllEdges()

int getAllEdges ( )

Determine all edges.

Returns
Number of all edges.

◆ getAllSides()

int getAllSides ( int  opt = 0)

Determine all mesh sides.

Returns
Number of all sides.

◆ getBoundaryNodes()

int getBoundaryNodes ( )

Determine all boundary nodes.

Returns
n Number of boundary nodes.

◆ getBoundarySides()

int getBoundarySides ( )

Determine all boundary sides.

Returns
Number of boundary sides.

◆ getEdgeLabel()

size_t getEdgeLabel ( size_t  i) const

Return label of i-th edge.

Parameters
[in]iEdge index

◆ getElementLabel()

size_t getElementLabel ( size_t  i) const

Return label of i-th element.

Parameters
[in]iElement index

◆ getElementNeighborElements()

void getElementNeighborElements ( )

Create element neighboring elements.

This function creates for each element the list of elements that share a side with it. Once this function is invoked, one can retrieve the list of neighboring elements of any element (Element::getNeigborElement)

◆ getList() [1/3]

void getList ( vector< Element * > &  el) const

Fill vector el with list of pointers to elements.

Parameters
[out]elInstance of class vector that contain on output the list

◆ getList() [2/3]

void getList ( vector< Node * > &  nl) const

Fill vector nl with list of pointers to nodes.

Parameters
[out]nlInstance of class vector that contain on output the list

◆ getList() [3/3]

void getList ( vector< Side * > &  sl) const

Fill vector sl with list of pointers to sides.

Parameters
[out]slInstance of class vector that contain on output the list

◆ getNbBoundarySides()

size_t getNbBoundarySides ( ) const

Return number of boundary sides.

This function is valid if member function getAllSides or getBoundarySides has been invoked before

◆ getNbInternalSides()

size_t getNbInternalSides ( ) const

Return number of internal sides.

This function is valid if member functions getAllSides and createInternalSideList have been invoked before

◆ getNodeLabel()

size_t getNodeLabel ( size_t  i) const

Return label of i-th node.

Parameters
[in]iNode index

◆ getNodeNeighborElements()

void getNodeNeighborElements ( )

Create node neighboring elements.

This function is generally useful when, for a numerical method, one looks for a given node to the list of elements that share this node. Once this function is invoked, one can retrieve the list of neighboring elements of any node (Node::getNeigEl)

◆ getSideLabel()

size_t getSideLabel ( size_t  i) const

Return label of i-th side.

Parameters
[in]iSide index

◆ NodesAreDOF()

bool NodesAreDOF ( ) const

Return information about DOF type.

Returns
true if DOF are supported by nodes, false otherwise

◆ NumberEquations() [1/2]

size_t NumberEquations ( size_t  dof,
int  c 
)

Renumber Equations.

Parameters
[in]dofLabel of degree of freedom for which numbering is performed.
[in]ccode for which degrees of freedom are enforced.

◆ NumberEquations() [2/2]

size_t NumberEquations ( size_t  dof = 0)

Renumber Equations.

Parameters
[in]dofLabel of degree of freedom for which numbering is performed. Default value (0) means that all degrees of freedom are taken into account

◆ operator*=()

Mesh & operator*= ( real_t  a)

Operator *=

Rescale mesh coordinates by myltiplying by a factor

Parameters
[in]aValue to multiply by

◆ put()

void put ( const string &  mesh_file) const

Write mesh data on file.

Parameters
[in]mesh_fileMesh file name

◆ RenumberEdge()

void RenumberEdge ( size_t  n1,
size_t  n2 
)

Renumber an edge.

Parameters
[in]n1Old label
[in]n2New label

◆ RenumberElement()

void RenumberElement ( size_t  n1,
size_t  n2 
)

Renumber an element.

Parameters
[in]n1Old label
[in]n2New label

◆ RenumberNode()

void RenumberNode ( size_t  n1,
size_t  n2 
)

Renumber a node.

Parameters
[in]n1Old label
[in]n2New label

◆ RenumberSide()

void RenumberSide ( size_t  n1,
size_t  n2 
)

Renumber a side.

Parameters
[in]n1Old label
[in]n2New label

◆ Reorder()

void Reorder ( size_t  m = GRAPH_MEMORY)

Renumber mesh nodes according to reverse Cuthill Mc Kee algorithm.

Parameters
[in]mMemory size needed for matrix graph (default value is GRAPH_MEMORY, see OFELI_Config.h)

◆ Rescale()

void Rescale ( real_t  sx,
real_t  sy = 0.,
real_t  sz = 0. 
)

Rescale mesh by multiplying node coordinates by constants.

This function can be used e.g. for changing coordinate units

Parameters
[in]sxFactor to multiply by x coordinates
[in]syFactor to multiply by y coordinates [Default: sx]
[in]szFactor to multiply by z coordinates [Default: sx]

◆ save()

void save ( const string &  mesh_file) const

Write mesh data on file in various formats.

File format depends on the extension in file name

Parameters
[in]mesh_fileMesh file name If the extension is '.m', the output file is an OFELI file If the extension is '.gpl', the output file is a Gnuplot file If the extension is '.msh' or '.geo', the output file is a Gmsh file If the extension is '.vtk', the output file is a VTK file

◆ set() [1/3]

void set ( Element el)

Replace element in the mesh.

If the element label exists already, the existing element pointer will be replaced by the current one. If not, an error message is displayed.

Parameters
[in]elPointer to element

◆ set() [2/3]

void set ( Node nd)

Replace node in the mesh.

If the node label exists already, the existing node pointer will be replaced by the current one. If not, an error message is displayed.

Parameters
[in]ndPointer to node

◆ set() [3/3]

void set ( Side sd)

Choose side in the mesh.

If the side label exists already, the existing side pointer will be replaced by the current one. If not, an error message is displayed.

Parameters
[in]sdPointer to side

◆ setDim()

void setDim ( size_t  dim)

Define space dimension. Normally, between 1 and 3.

Parameters
[in]dimSpace dimension to set (must be between 1 and 3)

◆ setDOFSupport()

void setDOFSupport ( int  opt,
int  nb_nodes = 1 
)

Define supports of degrees of freedom.

Parameters
[in]optDOF type:
  • NODE_DOF: Degrees of freedom are supported by nodes
  • SIDE_DOF: Degrees of freedom are supported by sides
  • EDGE_DOF: Degrees of freedom are supported by edges
  • ELEMENT_DOF: Degrees of freedom are supported by elements
[in]nb_nodesNumber of nodes on sides or elements (default=1). This parameter is useful only if dofs are supported by sides or elements
Note
This member function creates all mesh sides if the option ELEMENT_DOF or SIDE_DOF is selected. So it not necessary to call getAllSides() after

◆ setList() [1/3]

void setList ( const vector< Element * > &  el)

Initialize list of mesh elements using the input vector.

Parameters
[in]elvector instance that contains the list of pointers to elements

◆ setList() [2/3]

void setList ( const vector< Node * > &  nl)

Initialize list of mesh nodes using the input vector.

Parameters
[in]nlvector instance that contains the list of pointers to nodes

◆ setList() [3/3]

void setList ( const vector< Side * > &  sl)

Initialize list of mesh sides using the input vector.

Parameters
[in]slvector instance that contains the list of pointers to sides

◆ setMaterial()

void setMaterial ( int  code,
const string &  mname 
)

Associate material to code of element.

Parameters
[in]codeElement code for which material is assigned
[in]mnameName of material

◆ setNbDOFPerNode()

void setNbDOFPerNode ( size_t  nb_dof = 1)

Define number of degrees of freedom for each node.

Parameters
[in]nb_dofNumber of degrees of freedom (unknowns) for each mesh node (Default value is 1)
Note
This function first declares nodes as unknown supports, sets the number of degrees of freedom and renumbers equations

◆ setPointInDomain()

void setPointInDomain ( Point< real_t >  x)

Define a point in the domain. This function makes sense only if boundary mesh is given without internal mesh (Case of Boundary Elements)

Parameters
[in]xCoordinates of point to define

◆ SidesAreDOF()

bool SidesAreDOF ( ) const

Return information about DOF type.

Returns
true if DOF are supported by sides, false otherwise