File Mesh.h#
Generic mesh operations.
- Copyright
Copyright 2017. Tom de Geus. All rights reserved.
- License: This project is released under the GNU Public License (GPLv3).
-
namespace GooseFEM
Toolbox to perform finite element computations.
-
namespace Mesh
Generic mesh operations, and simple mesh definitions.
Enums
Functions
-
template<class D>
inline std::vector<std::vector<size_t>> nodaltyings(const D &dofs)# List nodal tyings based on DOF-numbers per node.
- Parameters
dofs – DOFs per node [nnode, ndim].
- Returns
Nodes to which the nodes is connected (sorted) [nnode, …]
-
template<class S, class T>
inline ElementType defaultElementType(const S &coor, const T &conn)# Extract the element type based on the connectivity.
- Parameters
coor – Nodal coordinates [nnode, ndim].
conn – Connectivity [nelem, nne].
- Returns
-
inline array_type::tensor<size_t, 2> dofs(size_t nnode, size_t ndim)#
List with DOF-numbers in sequential order.
The output is a sequential list of DOF-numbers for each vector-component of each node. For example for 3 nodes in 2 dimensions the output is
\( \begin{bmatrix} 0 & 1 \\ 2 & 3 \\ 4 & 5 \end{bmatrix} \)
- Parameters
nnode – Number of nodes.
ndim – Number of dimensions.
- Returns
DOF-numbers.
-
template<class T>
inline T renumber(const T &dofs)# Renumber to lowest possible index (see GooseFEM::Mesh::Renumber).
- Parameters
dofs – DOF-numbers [nnode, ndim].
- Returns
Renumbered DOF-numbers.
-
template<class S, class T>
inline array_type::tensor<size_t, 2> overlapping(const S &coor_a, const T &coor_b, double rtol = 1e-5, double atol = 1e-8)# Find overlapping nodes.
The output has the following structure:
[ [nodes_from_mesh_a], [nodes_from_mesh_b], ]
- Parameters
coor_a – Nodal coordinates of mesh “a” [nnode, ndim].
coor_b – Nodal coordinates of mesh “b” [nnode, ndim].
rtol – Relative tolerance for position match.
atol – Absolute tolerance for position match.
- Returns
Overlapping nodes.
-
template<class E>
inline array_type::tensor<size_t, 1> coordination(const E &conn)# Number of elements connected to each node.
- Parameters
conn – Connectivity [nelem, nne].
- Returns
Coordination per node.
-
template<class E>
inline std::vector<std::vector<size_t>> elem2node(const E &conn, bool sorted = true)# Elements connected to each node.
- Parameters
conn – Connectivity [nelem, nne].
sorted – If
true
the list of elements for each node is sorted.
- Returns
Elements per node [nnode, …].
-
template<class E, class D>
inline std::vector<std::vector<size_t>> elem2node(const E &conn, const D &dofs, bool sorted = true)# Elements connected to each node.
- Parameters
conn – Connectivity [nelem, nne].
sorted – If
true
the list of elements for each node is sorted.dofs – DOFs per node, allowing accounting for periodicity [nnode, ndim].
- Returns
Elements per node [nnode, …].
-
template<class D>
inline std::vector<std::vector<size_t>> node2dof(const D &dofs, bool sorted = true)# Nodes connected to each DOF.
- Parameters
dofs – DOFs per node [nnode, ndim].
sorted – If
true
the list of nodes for each DOF is sorted.
- Returns
Nodes per DOF [ndof, …].
-
template<class C, class E>
inline array_type::tensor<double, 2> edgesize(const C &coor, const E &conn, ElementType type)# Return size of each element edge.
- Parameters
coor – Nodal coordinates.
conn – Connectivity.
type – ElementType.
- Returns
Edge-sizes per element.
-
template<class C, class E>
inline array_type::tensor<double, 2> edgesize(const C &coor, const E &conn)# Return size of each element edge.
The element-type is automatically determined, see defaultElementType().
- Parameters
coor – Nodal coordinates.
conn – Connectivity.
- Returns
Edge-sizes per element.
-
template<class C, class E>
inline array_type::tensor<double, 2> centers(const C &coor, const E &conn, ElementType type)# Coordinates of the center of each element.
- Parameters
coor – Nodal coordinates.
conn – Connectivity.
type – ElementType.
- Returns
Center of each element.
-
template<class C, class E>
inline array_type::tensor<double, 2> centers(const C &coor, const E &conn)# Coordinates of the center of each element.
The element-type is automatically determined, see defaultElementType().
- Parameters
coor – Nodal coordinates.
conn – Connectivity.
- Returns
Center of each element.
-
template<class T, class C, class E>
inline array_type::tensor<size_t, 1> elemmap2nodemap(const T &elem_map, const C &coor, const E &conn, ElementType type)# Convert an element-map to a node-map.
- Parameters
elem_map – Element-map such that
new_elvar = elvar[elem_map]
.coor – Nodal coordinates.
conn – Connectivity.
type – ElementType.
- Returns
Node-map such that
new_nodevar = nodevar[node_map]
-
template<class T, class C, class E>
inline array_type::tensor<size_t, 1> elemmap2nodemap(const T &elem_map, const C &coor, const E &conn)# Convert an element-map to a node-map.
The element-type is automatically determined, see defaultElementType().
- Parameters
elem_map – Element-map such that
new_elvar = elvar[elem_map]
.coor – Nodal coordinates.
conn – Connectivity.
- Returns
Node-map such that
new_nodevar = nodevar[node_map]
-
template<class C, class E>
inline array_type::tensor<double, 2> nodal_mass(const C &coor, const E &conn, ElementType type)# Compute the mass of each node in the mesh.
If nodes are not part of the connectivity the mass is set to zero, such that the center of gravity is simply::
average(coor, GooseFEM.Mesh.nodal_mass(coor, conn, type), axis=0);
- Template Parameters
C – e.g.
array_type::tensor<double, 2>
E – e.g.
array_type::tensor<size_t, 2>
- Parameters
coor – Nodal coordinates
[nnode, ndim]
.conn – Connectivity
[nelem, nne]
.type – ElementType.
- Returns
Center of gravity
[ndim]
.
-
template<class C, class E>
inline array_type::tensor<double, 2> nodal_mass(const C &coor, const E &conn)# Compute the mass of each node in the mesh.
If nodes are not part of the connectivity the mass is set to zero, such that the center of gravity is simply::
average(coor, GooseFEM.Mesh.nodal_mass(coor, conn), axis=0);
- Template Parameters
C – e.g.
array_type::tensor<double, 2>
E – e.g.
array_type::tensor<size_t, 2>
- Parameters
coor – Nodal coordinates
[nnode, ndim]
.conn – Connectivity
[nelem, nne]
.
- Returns
Center of gravity
[ndim]
.
-
template<class C, class E>
inline array_type::tensor<double, 1> center_of_gravity(const C &coor, const E &conn, ElementType type)# Compute the center of gravity of a mesh.
- Template Parameters
C – e.g.
array_type::tensor<double, 2>
E – e.g.
array_type::tensor<size_t, 2>
- Parameters
coor – Nodal coordinates
[nnode, ndim]
.conn – Connectivity
[nelem, nne]
.type – ElementType.
- Returns
Center of gravity
[ndim]
.
-
template<class C, class E>
inline array_type::tensor<double, 1> center_of_gravity(const C &coor, const E &conn)# Compute the center of gravity of a mesh.
- Template Parameters
C – e.g.
array_type::tensor<double, 2>
E – e.g.
array_type::tensor<size_t, 2>
- Parameters
coor – Nodal coordinates
[nnode, ndim]
.conn – Connectivity
[nelem, nne]
.
- Returns
Center of gravity
[ndim]
.
-
class ManualStitch
- #include <GooseFEM/Mesh.h>
Stitch two mesh objects, specifying overlapping nodes by hand.
Public Functions
-
ManualStitch() = default#
-
template<class CA, class EA, class NA, class CB, class EB, class NB>
inline ManualStitch(const CA &coor_a, const EA &conn_a, const NA &overlapping_nodes_a, const CB &coor_b, const EB &conn_b, const NB &overlapping_nodes_b, bool check_position = true, double rtol = 1e-5, double atol = 1e-8) - Parameters
coor_a – Nodal coordinates of mesh “a” [nnode, ndim].
conn_a – Connectivity of mesh “a” [nelem, nne].
overlapping_nodes_a – Node-numbers of mesh “a” that overlap with mesh “b” [n].
coor_b – Nodal coordinates of mesh “b” [nnode, ndim].
conn_b – Connectivity of mesh “b” [nelem, nne].
overlapping_nodes_b – Node-numbers of mesh “b” that overlap with mesh “a” [n].
check_position – If
true
the nodes are checked for position overlap.rtol – Relative tolerance for check on position overlap.
atol – Absolute tolerance for check on position overlap.
-
inline size_t nmesh() const
Number of sub meshes == 2.
- Returns
unsigned int
-
inline size_t nelem() const
Number of elements.
- Returns
unsigned int
-
inline size_t nnode() const
Number of nodes.
- Returns
unsigned int
-
inline size_t nne() const
Number of nodes-per-element.
- Returns
unsigned int
-
inline size_t ndim() const
Number of dimensions.
- Returns
unsigned int
-
inline const array_type::tensor<double, 2> &coor() const
Nodal coordinates [nnode, ndim].
- Returns
coordinates per node
-
inline const array_type::tensor<size_t, 2> &conn() const
-
- Returns
nodes per element
-
inline array_type::tensor<size_t, 2> dofs() const
DOF numbers for each node (numbered sequentially) [nnode, ndim].
- Returns
DOFs per node
-
inline std::vector<array_type::tensor<size_t, 1>> nodemap() const
Node-map per sub-mesh.
- Returns
nodes per mesh
-
inline std::vector<array_type::tensor<size_t, 1>> elemmap() const
Element-map per sub-mesh.
- Returns
elements per mesh
-
inline array_type::tensor<size_t, 1> nodemap(size_t mesh_index) const
- Parameters
mesh_index – Index of the mesh (“a” = 1, “b” = 1).
- Returns
Node-map for a given mesh.
-
inline array_type::tensor<size_t, 1> elemmap(size_t mesh_index) const
- Parameters
mesh_index – Index of the mesh (“a” = 1, “b” = 1).
- Returns
Element-map for a given mesh.
Private Members
-
array_type::tensor<double, 2> m_coor#
-
array_type::tensor<size_t, 2> m_conn#
-
array_type::tensor<size_t, 1> m_map_b#
-
size_t m_nnd_a#
-
size_t m_nel_a#
-
size_t m_nel_b#
-
ManualStitch() = default#
-
template<class D>
class RegularBase - #include <GooseFEM/Mesh.h>
CRTP base class for regular meshes.
Subclassed by GooseFEM::Mesh::RegularBase2d< FineLayer >, GooseFEM::Mesh::RegularBase2d< Regular >, GooseFEM::Mesh::RegularBase3d< FineLayer >, GooseFEM::Mesh::RegularBase3d< Regular >, GooseFEM::Mesh::RegularBase2d< D >, GooseFEM::Mesh::RegularBase3d< D >
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
inline auto nelem() const
Number of elements.
- Returns
unsigned int
-
inline auto nnode() const
Number of nodes.
- Returns
unsigned int
-
inline auto nne() const
Number of nodes-per-element == 4.
- Returns
unsigned int
-
inline auto ndim() const
Number of dimensions == 2.
- Returns
unsigned int
-
inline auto nelx() const
Number of elements in x-direction == width of the mesh in units of h.
- Returns
unsigned int
-
inline auto nely() const
Number of elements in y-direction == height of the mesh, in units of h,.
- Returns
unsigned int
-
inline auto h() const
Linear edge size of one ‘block’.
- Returns
double
-
inline auto elementType() const
The ElementType().
- Returns
element type
-
inline auto dofs() const
DOF numbers for each node (numbered sequentially) [nnode, ndim].
- Returns
DOFs per node
-
inline auto dofsPeriodic() const
DOF-numbers for the case that the periodicity if fully eliminated.
-
inline auto nodesPeriodic() const
Periodic node pairs, in two columns: (independent, dependent).
- Returns
[ntyings, ndim].
-
inline auto nodesOrigin() const
Reference node to use for periodicity, because all corners are tied to it.
- Returns
Node number.
Private Functions
-
inline auto derived_cast() -> derived_type&#
-
inline auto derived_cast() const -> const derived_type&#
-
using derived_type = D
-
template<class D>
class RegularBase2d : public GooseFEM::Mesh::RegularBase<D> - #include <GooseFEM/Mesh.h>
CRTP base class for regular meshes in 2d.
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
inline auto nodesBottomEdge() const
Nodes along the bottom edge (y = 0), in order of increasing x.
- Returns
List of node numbers.
-
inline auto nodesTopEdge() const
Nodes along the top edge (y = nely * h), in order of increasing x.
- Returns
List of node numbers.
-
inline auto nodesLeftEdge() const
Nodes along the left edge (x = 0), in order of increasing y.
- Returns
List of node numbers.
-
inline auto nodesRightEdge() const
Nodes along the right edge (x = nelx * h), in order of increasing y.
- Returns
List of node numbers.
-
inline auto nodesBottomOpenEdge() const
Nodes along the bottom edge (y = 0), without the corners (at x = 0 and x = nelx * h).
Same as: nodesBottomEdge()[1: -1].
- Returns
List of node numbers.
-
inline auto nodesTopOpenEdge() const
Nodes along the top edge (y = nely * h), without the corners (at x = 0 and x = nelx * h).
Same as: nodesTopEdge()[1: -1].
- Returns
List of node numbers.
-
inline auto nodesLeftOpenEdge() const
Nodes along the left edge (x = 0), without the corners (at y = 0 and y = nely * h).
Same as: nodesLeftEdge()[1: -1].
- Returns
List of node numbers.
-
inline auto nodesRightOpenEdge() const
Nodes along the right edge (x = nelx * h), without the corners (at y = 0 and y = nely * h).
Same as: nodesRightEdge()[1: -1].
- Returns
List of node numbers.
-
inline auto nodesBottomLeftCorner() const
The bottom-left corner node (at x = 0, y = 0).
Same as nodesBottomEdge()[0] and nodesLeftEdge()[0].
- Returns
Node number.
-
inline auto nodesBottomRightCorner() const
The bottom-right corner node (at x = nelx * h, y = 0).
Same as nodesBottomEdge()[-1] and nodesRightEdge()[0].
- Returns
Node number.
-
inline auto nodesTopLeftCorner() const
The top-left corner node (at x = 0, y = nely * h).
Same as nodesTopEdge()[0] and nodesRightEdge()[-1].
- Returns
Node number.
-
inline auto nodesTopRightCorner() const
The top-right corner node (at x = nelx * h, y = nely * h).
Same as nodesTopEdge()[-1] and nodesRightEdge()[-1].
- Returns
Node number.
-
inline auto nodesLeftBottomCorner() const
Alias of nodesBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesLeftTopCorner() const
Alias of nodesTopLeftCorner().
- Returns
Node number.
-
inline auto nodesRightBottomCorner() const
Alias of nodesBottomRightCorner().
- Returns
Node number.
-
inline auto nodesRightTopCorner() const
Alias of nodesTopRightCorner().
- Returns
Node number.
Private Functions
-
inline auto derived_cast() -> derived_type&#
-
inline auto derived_cast() const -> const derived_type&#
-
inline array_type::tensor<size_t, 2> nodesPeriodic_impl() const#
-
inline auto nodesOrigin_impl() const#
Friends
- friend class RegularBase< D >
-
using derived_type = D
-
template<class D>
class RegularBase3d : public GooseFEM::Mesh::RegularBase<D> - #include <GooseFEM/Mesh.h>
CRTP base class for regular meshes in 3d.
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
inline auto nelz() const
Number of elements in y-direction == height of the mesh, in units of h,.
- Returns
unsigned int
-
inline auto nodesBottom() const
Nodes along the bottom face (y = 0).
- Returns
List of node numbers.
-
inline auto nodesLeft() const
Nodes along the left face (x = 0).
- Returns
List of node numbers.
-
inline auto nodesRight() const
Nodes along the right face (x = nelx * h).
- Returns
List of node numbers.
-
inline auto nodesFront() const
Nodes along the front face (z = 0).
- Returns
List of node numbers.
-
inline auto nodesBack() const
Nodes along the back face (z = nelz * h).
- Returns
List of node numbers.
-
inline auto nodesFrontBottomEdge() const
Nodes along the edge at the intersection of the front and bottom faces (z = 0 and y = 0).
- Returns
List of node numbers.
-
inline auto nodesFrontTopEdge() const
Nodes along the edge at the intersection of the front and top faces (z = 0 and y = nely * h).
- Returns
List of node numbers.
-
inline auto nodesFrontLeftEdge() const
Nodes along the edge at the intersection of the front and left faces (z = 0 and x = 0).
- Returns
List of node numbers.
-
inline auto nodesFrontRightEdge() const
Nodes along the edge at the intersection of the front and right faces (z = 0 and x = nelx * h).
- Returns
List of node numbers.
-
inline auto nodesBackBottomEdge() const
Nodes along the edge at the intersection of the back and bottom faces (z = nelz * h and y = nely * h).
- Returns
List of node numbers.
-
inline auto nodesBackTopEdge() const
Nodes along the edge at the intersection of the back and top faces (z = nelz * h and x = 0).
- Returns
List of node numbers.
-
inline auto nodesBackLeftEdge() const
Nodes along the edge at the intersection of the back and left faces (z = nelz * h and x = nelx * h).
- Returns
List of node numbers.
-
inline auto nodesBackRightEdge() const
Nodes along the edge at the intersection of the back and right faces (? = nelz * h and ? = ?).
- Returns
List of node numbers.
-
inline auto nodesBottomLeftEdge() const
Nodes along the edge at the intersection of the bottom and left faces (y = 0 and x = 0).
- Returns
List of node numbers.
-
inline auto nodesBottomRightEdge() const
Nodes along the edge at the intersection of the bottom and right faces (y = 0 and x = nelx * h).
- Returns
List of node numbers.
-
inline auto nodesTopLeftEdge() const
Nodes along the edge at the intersection of the top and left faces (y = 0 and x = nelx * h).
- Returns
List of node numbers.
-
inline auto nodesTopRightEdge() const
Nodes along the edge at the intersection of the top and right faces (y = nely * h and x = nelx * h).
- Returns
List of node numbers.
-
inline auto nodesBottomFrontEdge() const
Alias of nodesFrontBottomEdge()
- Returns
List of node numbers.
-
inline auto nodesBottomBackEdge() const
Alias of nodesBackBottomEdge()
- Returns
List of node numbers.
-
inline auto nodesTopFrontEdge() const
Alias of nodesFrontTopEdge()
- Returns
List of node numbers.
-
inline auto nodesTopBackEdge() const
Alias of nodesBackTopEdge()
- Returns
List of node numbers.
-
inline auto nodesLeftBottomEdge() const
Alias of nodesBottomLeftEdge()
- Returns
List of node numbers.
-
inline auto nodesLeftFrontEdge() const
Alias of nodesFrontLeftEdge()
- Returns
List of node numbers.
-
inline auto nodesLeftBackEdge() const
Alias of nodesBackLeftEdge()
- Returns
List of node numbers.
-
inline auto nodesLeftTopEdge() const
Alias of nodesTopLeftEdge()
- Returns
List of node numbers.
-
inline auto nodesRightBottomEdge() const
Alias of nodesBottomRightEdge()
- Returns
List of node numbers.
-
inline auto nodesRightTopEdge() const
Alias of nodesTopRightEdge()
- Returns
List of node numbers.
-
inline auto nodesRightFrontEdge() const
Alias of nodesFrontRightEdge()
- Returns
List of node numbers.
-
inline auto nodesRightBackEdge() const
Alias of nodesBackRightEdge()
- Returns
List of node numbers.
-
inline auto nodesFrontFace() const
Nodes along the front face excluding edges.
Same as different between nodesFront() and [nodesFrontBottomEdge(), nodesFrontTopEdge(), nodesFrontLeftEdge(), nodesFrontRightEdge()]
- Returns
list of node numbers.
-
inline auto nodesBackFace() const
Nodes along the back face excluding edges.
Same as different between nodesBack() and [nodesBackBottomEdge(), nodesBackTopEdge(), nodesBackLeftEdge(), nodesBackRightEdge()]
- Returns
list of node numbers.
-
inline auto nodesLeftFace() const
Nodes along the left face excluding edges.
Same as different between nodesLeft() and [nodesFrontLeftEdge(), nodesBackLeftEdge(), nodesBottomLeftEdge(), nodesTopLeftEdge()]
- Returns
list of node numbers.
-
inline auto nodesRightFace() const
Nodes along the right face excluding edges.
Same as different between nodesRight() and [nodesFrontRightEdge(), nodesBackRightEdge(), nodesBottomRightEdge(), nodesTopRightEdge()]
- Returns
list of node numbers.
-
inline auto nodesBottomFace() const
Nodes along the bottom face excluding edges.
Same as different between nodesBottom() and [nodesBackBottomEdge(), nodesBackTopEdge(), nodesBackLeftEdge(), nodesBackRightEdge()]
- Returns
list of node numbers.
-
inline auto nodesTopFace() const
Nodes along the top face excluding edges.
Same as different between nodesTop() and [nodesFrontBottomEdge(), nodesFrontTopEdge(), nodesFrontLeftEdge(), nodesFrontRightEdge()]
- Returns
list of node numbers.
-
inline auto nodesFrontBottomOpenEdge() const
Same as nodesFrontBottomEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesFrontTopOpenEdge() const
Same as nodesFrontTopEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesFrontLeftOpenEdge() const
Same as nodesFrontLeftEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesFrontRightOpenEdge() const
Same as nodesFrontRightEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesBackBottomOpenEdge() const
Same as nodesBackBottomEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesBackTopOpenEdge() const
Same as nodesBackTopEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesBackLeftOpenEdge() const
Same as nodesBackLeftEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesBackRightOpenEdge() const
Same as nodesBackRightEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesBottomLeftOpenEdge() const
Same as nodesBottomLeftEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesBottomRightOpenEdge() const
Same as nodesBottomRightEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesTopLeftOpenEdge() const
Same as nodesTopLeftEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesTopRightOpenEdge() const
Same as nodesTopRightEdge() but without corners.
- Returns
List of node numbers.
-
inline auto nodesBottomFrontOpenEdge() const
Alias of nodesFrontBottomOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesBottomBackOpenEdge() const
Alias of nodesBackBottomOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesTopFrontOpenEdge() const
Alias of nodesFrontTopOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesTopBackOpenEdge() const
Alias of nodesBackTopOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesLeftBottomOpenEdge() const
Alias of nodesBottomLeftOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesLeftFrontOpenEdge() const
Alias of nodesFrontLeftOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesLeftBackOpenEdge() const
Alias of nodesBackLeftOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesLeftTopOpenEdge() const
Alias of nodesTopLeftOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesRightBottomOpenEdge() const
Alias of nodesBottomRightOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesRightTopOpenEdge() const
Alias of nodesTopRightOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesRightFrontOpenEdge() const
Alias of nodesFrontRightOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesRightBackOpenEdge() const
Alias of nodesBackRightOpenEdge().
- Returns
List of node numbers.
-
inline auto nodesFrontBottomLeftCorner() const
Front-Bottom-Left corner node.
- Returns
Node number.
-
inline auto nodesFrontBottomRightCorner() const
Front-Bottom-Right corner node.
- Returns
Node number.
-
inline auto nodesFrontTopLeftCorner() const
Front-Top-Left corner node.
- Returns
Node number.
-
inline auto nodesFrontTopRightCorner() const
Front-Top-Right corner node.
- Returns
Node number.
-
inline auto nodesBackBottomLeftCorner() const
Back-Bottom-Left corner node.
- Returns
Node number.
-
inline auto nodesBackBottomRightCorner() const
Back-Bottom-Right corner node.
- Returns
Node number.
-
inline auto nodesBackTopLeftCorner() const
Back-Top-Left corner node.
- Returns
Node number.
-
inline auto nodesBackTopRightCorner() const
Back-Top-Right corner node.
- Returns
Node number.
-
inline auto nodesFrontLeftBottomCorner() const
Alias of nodesFrontBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesBottomFrontLeftCorner() const
Alias of nodesFrontBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesBottomLeftFrontCorner() const
Alias of nodesFrontBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesLeftFrontBottomCorner() const
Alias of nodesFrontBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesLeftBottomFrontCorner() const
Alias of nodesFrontBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesFrontRightBottomCorner() const
Alias of nodesFrontBottomRightCorner().
- Returns
Node number.
-
inline auto nodesBottomFrontRightCorner() const
Alias of nodesFrontBottomRightCorner().
- Returns
Node number.
-
inline auto nodesBottomRightFrontCorner() const
Alias of nodesFrontBottomRightCorner().
- Returns
Node number.
-
inline auto nodesRightFrontBottomCorner() const
Alias of nodesFrontBottomRightCorner().
- Returns
Node number.
-
inline auto nodesRightBottomFrontCorner() const
Alias of nodesFrontBottomRightCorner().
- Returns
Node number.
-
inline auto nodesFrontLeftTopCorner() const
Alias of nodesFrontTopLeftCorner().
- Returns
Node number.
-
inline auto nodesTopFrontLeftCorner() const
Alias of nodesFrontTopLeftCorner().
- Returns
Node number.
-
inline auto nodesTopLeftFrontCorner() const
Alias of nodesFrontTopLeftCorner().
- Returns
Node number.
-
inline auto nodesLeftFrontTopCorner() const
Alias of nodesFrontTopLeftCorner().
- Returns
Node number.
-
inline auto nodesLeftTopFrontCorner() const
Alias of nodesFrontTopLeftCorner().
- Returns
Node number.
-
inline auto nodesFrontRightTopCorner() const
Alias of nodesFrontTopRightCorner().
- Returns
Node number.
-
inline auto nodesTopFrontRightCorner() const
Alias of nodesFrontTopRightCorner().
- Returns
Node number.
-
inline auto nodesTopRightFrontCorner() const
Alias of nodesFrontTopRightCorner().
- Returns
Node number.
-
inline auto nodesRightFrontTopCorner() const
Alias of nodesFrontTopRightCorner().
- Returns
Node number.
-
inline auto nodesRightTopFrontCorner() const
Alias of nodesFrontTopRightCorner().
- Returns
Node number.
-
inline auto nodesBackLeftBottomCorner() const
Alias of nodesBackBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesBottomBackLeftCorner() const
Alias of nodesBackBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesBottomLeftBackCorner() const
Alias of nodesBackBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesLeftBackBottomCorner() const
Alias of nodesBackBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesLeftBottomBackCorner() const
Alias of nodesBackBottomLeftCorner().
- Returns
Node number.
-
inline auto nodesBackRightBottomCorner() const
Alias of nodesBackBottomRightCorner().
- Returns
Node number.
-
inline auto nodesBottomBackRightCorner() const
Alias of nodesBackBottomRightCorner().
- Returns
Node number.
-
inline auto nodesBottomRightBackCorner() const
Alias of nodesBackBottomRightCorner().
- Returns
Node number.
-
inline auto nodesRightBackBottomCorner() const
Alias of nodesBackBottomRightCorner().
- Returns
Node number.
-
inline auto nodesRightBottomBackCorner() const
Alias of nodesBackBottomRightCorner().
- Returns
Node number.
-
inline auto nodesBackLeftTopCorner() const
Alias of nodesBackTopLeftCorner().
- Returns
Node number.
-
inline auto nodesTopBackLeftCorner() const
Alias of nodesBackTopLeftCorner().
- Returns
Node number.
-
inline auto nodesTopLeftBackCorner() const
Alias of nodesBackTopLeftCorner().
- Returns
Node number.
-
inline auto nodesLeftBackTopCorner() const
Alias of nodesBackTopLeftCorner().
- Returns
Node number.
-
inline auto nodesLeftTopBackCorner() const
Alias of nodesBackTopLeftCorner().
- Returns
Node number.
-
inline auto nodesBackRightTopCorner() const
Alias of nodesBackTopRightCorner().
- Returns
Node number.
-
inline auto nodesTopBackRightCorner() const
Alias of nodesBackTopRightCorner().
- Returns
Node number.
-
inline auto nodesTopRightBackCorner() const
Alias of nodesBackTopRightCorner().
- Returns
Node number.
-
inline auto nodesRightBackTopCorner() const
Alias of nodesBackTopRightCorner().
- Returns
Node number.
-
inline auto nodesRightTopBackCorner() const
Alias of nodesBackTopRightCorner().
- Returns
Node number.
Private Functions
-
inline auto derived_cast() -> derived_type&#
-
inline auto derived_cast() const -> const derived_type&#
-
inline array_type::tensor<size_t, 2> nodesPeriodic_impl() const#
-
inline auto nodesOrigin_impl() const#
Friends
- friend class RegularBase< D >
-
using derived_type = D
-
class Renumber
- #include <GooseFEM/Mesh.h>
Renumber indices to lowest possible index.
For example:
\( \begin{bmatrix} 0 & 1 \\ 5 & 4 \end{bmatrix} \)
is renumbered to
\( \begin{bmatrix} 0 & 1 \\ 3 & 2 \end{bmatrix} \)
Or, in pseudo-code, the result of this function is that:
dofs = renumber(dofs) sort(unique(dofs[:])) == range(max(dofs+1))
Note
One can use the wrapper function renumber(). This class gives more advanced features.
Public Functions
-
Renumber() = default#
-
template<class T>
inline Renumber(const T &dofs) - Parameters
dofs – DOF-numbers.
-
template<class T>
inline T apply(const T &list) const Apply renumbering to other set.
- Parameters
list – List of (DOF-)numbers.
- Returns
Renumbered list of (DOF-)numbers.
-
inline const array_type::tensor<size_t, 1> &index() const
Get the list needed to renumber, e.g.:
dofs_renumbered(i, j) = index(dofs(i, j))
- Returns
Renumber-index.
Private Members
-
array_type::tensor<size_t, 1> m_renum#
-
Renumber() = default#
-
class Reorder
- #include <GooseFEM/Mesh.h>
Reorder to lowest possible index, in specific order.
For example for
Reorder({iiu, iip})
after reordering:iiu = xt::range<size_t>(nnu); iip = xt::range<size_t>(nnp) + nnu;
Public Functions
-
Reorder() = default#
-
template<class T>
inline Reorder(const std::initializer_list<T> args) - Parameters
args – List of (DOF-)numbers.
-
template<class T>
inline T apply(const T &list) const Apply reordering to other set.
- Parameters
list – List of (DOF-)numbers.
- Returns
Reordered list of (DOF-)numbers.
-
inline const array_type::tensor<size_t, 1> &index() const
Get the list needed to reorder, e.g.:
dofs_reordered(i, j) = index(dofs(i, j))
- Returns
Reorder-index.
Private Members
-
array_type::tensor<size_t, 1> m_renum#
-
Reorder() = default#
-
class Stitch
- #include <GooseFEM/Mesh.h>
Stitch mesh objects, automatically searching for overlapping nodes.
Subclassed by GooseFEM::Mesh::Vstack
Public Functions
-
inline Stitch(double rtol = 1e-5, double atol = 1e-8)
- Parameters
rtol – Relative tolerance for position match.
atol – Absolute tolerance for position match.
-
template<class C, class E>
inline void push_back(const C &coor, const E &conn) Add mesh to be stitched.
- Parameters
coor – Nodal coordinates [nnode, ndim].
conn – Connectivity [nelem, nne].
-
inline size_t nmesh() const
Number of sub meshes.
- Returns
unsigned int
-
inline size_t nelem() const
Number of elements.
- Returns
unsigned int
-
inline size_t nnode() const
Number of nodes.
- Returns
unsigned int
-
inline size_t nne() const
Number of nodes-per-element.
- Returns
unsigned int
-
inline size_t ndim() const
Number of dimensions.
- Returns
unsigned int
-
inline const array_type::tensor<double, 2> &coor() const
Nodal coordinates [nnode, ndim].
- Returns
coordinates per node
-
inline const array_type::tensor<size_t, 2> &conn() const
-
- Returns
nodes per element
-
inline array_type::tensor<size_t, 2> dofs() const
DOF numbers for each node (numbered sequentially) [nnode, ndim].
- Returns
DOFs per node
-
inline std::vector<array_type::tensor<size_t, 1>> nodemap() const
Node-map per sub-mesh.
- Returns
nodes per mesh
-
inline std::vector<array_type::tensor<size_t, 1>> elemmap() const
Element-map per sub-mesh.
- Returns
elements per mesh
-
inline array_type::tensor<size_t, 1> nodemap(size_t mesh_index) const
The node numbers in the stitched mesh that are coming from a specific sub-mesh.
- Parameters
mesh_index – Index of the sub-mesh.
- Returns
List of node numbers.
-
inline array_type::tensor<size_t, 1> elemmap(size_t mesh_index) const
The element numbers in the stitched mesh that are coming from a specific sub-mesh.
- Parameters
mesh_index – Index of the sub-mesh.
- Returns
List of element numbers.
-
template<class T>
inline T nodeset(const T &set, size_t mesh_index) const Convert set of node-numbers for a sub-mesh to the stitched mesh.
- Parameters
set – List of node numbers.
mesh_index – Index of the sub-mesh.
- Returns
List of node numbers for the stitched mesh.
-
template<class T>
inline T elemset(const T &set, size_t mesh_index) const Convert set of element-numbers for a sub-mesh to the stitched mesh.
- Parameters
set – List of element numbers.
mesh_index – Index of the sub-mesh.
- Returns
List of element numbers for the stitched mesh.
-
template<class T>
inline T nodeset(const std::vector<T> &set) const Combine set of node numbers for an original to the final mesh (removes duplicates).
- Parameters
set – List of node numbers per mesh.
- Returns
List of node numbers for the stitched mesh.
-
template<class T>
inline T nodeset(std::initializer_list<T> set) const Combine set of node numbers for an original to the final mesh (removes duplicates).
- Parameters
set – List of node numbers per mesh.
- Returns
List of node numbers for the stitched mesh.
Protected Attributes
-
array_type::tensor<double, 2> m_coor#
-
array_type::tensor<size_t, 2> m_conn#
-
std::vector<array_type::tensor<size_t, 1>> m_map#
See nodemap(size_t)
-
double m_rtol#
Relative tolerance to find overlapping nodes.
-
double m_atol#
Absolute tolerance to find overlapping nodes.
-
inline Stitch(double rtol = 1e-5, double atol = 1e-8)
-
class Vstack : public GooseFEM::Mesh::Stitch
- #include <GooseFEM/Mesh.h>
Vertically stack meshes.
Public Functions
-
inline Vstack(bool check_overlap = true, double rtol = 1e-5, double atol = 1e-8)
- Parameters
check_overlap – Check if nodes are overlapping when adding a mesh.
rtol – Relative tolerance for position match.
atol – Absolute tolerance for position match.
-
template<class C, class E, class N>
inline void push_back(const C &coor, const E &conn, const N &nodes_bot, const N &nodes_top) Add a mesh to the top of the current stack.
Each time the current
nodes_bot
are stitched with the then highestnodes_top
.- Parameters
coor – Nodal coordinates [nnode, ndim].
conn – Connectivity [nelem, nne].
nodes_bot – Nodes along the bottom edge [n].
nodes_top – Nodes along the top edge [n].
Private Members
-
std::vector<array_type::tensor<size_t, 1>> m_nodes_bot#
Bottom nodes of each mesh (renumbered).
-
std::vector<array_type::tensor<size_t, 1>> m_nodes_top#
Top nodes of each mesh (renumbered).
-
bool m_check_overlap#
Check if nodes are overlapping when adding a mesh.
-
inline Vstack(bool check_overlap = true, double rtol = 1e-5, double atol = 1e-8)
-
template<class D>
-
namespace Mesh