Namespace GooseFEM#
-
namespace GooseFEM
Toolbox to perform finite element computations.
Functions
-
template<class T, class R>
inline void asTensor(const T &arg, R &ret) “Broadcast” a scalar stored in an array (e.g.
[r, s]
) to the same scalar of all tensor components of a tensor of certain rank (e.g. for rank 2:[r, s, i, j]
).- Parameters
arg – An array with scalars.
ret – Corresponding array with tensors.
-
template<class T, class S>
inline auto AsTensor(const T &arg, const S &shape) “Broadcast” a scalar stored in an array (e.g.
[r, s]
) to the same scalar of all tensor components of a tensor of certain rank (e.g. for rank 2:[r, s, i, j]
).- Parameters
arg – An array with scalars.
shape – The shape of the added tensor dimensions (e.g.:
[i, j]
).
- Returns
Corresponding array with tensors.
-
template<class T, class I, size_t L>
inline auto AsTensor(const T &arg, const I (&shape)[L]) “Broadcast” a scalar stored in an array (e.g.
[r, s]
) to the same scalar of all tensor components of a tensor of certain rank (e.g. for rank 2:[r, s, i, j]
).- Parameters
arg – An array with scalars.
shape – The shape of the added tensor dimensions (e.g.:
[i, j]
).
- Returns
Corresponding array with tensors.
-
template<size_t rank, class T>
inline auto AsTensor(const T &arg, size_t n) “Broadcast” a scalar stored in an array (e.g.
[r, s]
) to the same scalar of all tensor components of a tensor of certain rank (e.g. for rank 2:[r, s, n, n]
).- Template Parameters
rank – Number of tensor dimensions (number of dimensions to add to the input).
- Parameters
arg – An array with scalars.
n – The shape along each of the added dimensions.
- Returns
Corresponding array with tensors.
-
template<class T>
inline auto AsTensor(size_t rank, const T &arg, size_t n) “Broadcast” a scalar stored in an array (e.g.
[r, s]
) to the same scalar of all tensor components of a tensor of certain rank (e.g. for rank 2:[r, s, n, n]
).- Parameters
rank – Number of tensor dimensions (number of dimensions to add to the input).
arg – An array with scalars.
n – The shape along each of the added dimensions.
- Returns
Corresponding array with tensors.
-
template<class T>
inline T as3d(const T &arg) Zero-pad columns to a matrix until is that shape
[m, 3]
.- Parameters
arg – A “nodevec” (
arg.shape(1) <= 3
).- Returns
Corresponding “nodevec” in 3-d (
ret.shape(1) == 3
)
-
template<class T>
inline bool is_unique(const T &arg) Returns true is a list is unique (has not duplicate items).
- Parameters
arg – Array-like.
- Returns
true
ifarg
is unique.
-
inline std::string version()
Return version string, e.g.
"0.8.0"
- Returns
String.
-
class Matrix : public GooseFEM::MatrixBase<Matrix>
- #include <GooseFEM/Matrix.h>
Sparse matrix.
See GooseFEM::Vector() for bookkeeping definitions.
Public Functions
-
inline Matrix(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<size_t, 2> &dofs)
Constructor.
-
inline const Eigen::SparseMatrix<double> &data() const
Pointer to data.
-
inline void set(const array_type::tensor<size_t, 1> &rows, const array_type::tensor<size_t, 1> &cols, const array_type::tensor<double, 2> &matrix)
Overwrite matrix.
- Parameters
rows – Row numbers [m].
cols – Column numbers [n].
matrix – Data entries
matrix(i, j)
forrows(i), cols(j)
[m, n].
-
inline void add(const array_type::tensor<size_t, 1> &rows, const array_type::tensor<size_t, 1> &cols, const array_type::tensor<double, 2> &matrix)
Add matrix.
- Parameters
rows – Row numbers [m].
cols – Column numbers [n].
matrix – Data entries
matrix(i, j)
forrows(i), cols(j)
[m, n].
-
inline Matrix(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<size_t, 2> &dofs)
-
template<class D>
class MatrixBase - #include <GooseFEM/Matrix.h>
CRTP base class for a matrix.
Subclassed by GooseFEM::MatrixPartitionedBase< MatrixDiagonalPartitioned >, GooseFEM::MatrixPartitionedBase< MatrixPartitioned >, GooseFEM::MatrixPartitionedBase< MatrixPartitionedTyings >, GooseFEM::MatrixPartitionedBase< D >
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
inline size_t nelem() const
Number of elements.
- Returns
Unsigned integer.
-
inline size_t nne() const
Number of nodes per element.
- Returns
Unsigned integer.
-
inline size_t nnode() const
Number of nodes.
- Returns
Unsigned integer.
-
inline size_t ndim() const
Number of dimensions.
- Returns
Unsigned integer.
-
inline size_t ndof() const
Number of DOFs.
- Returns
Unsigned integer.
-
inline const array_type::tensor<size_t, 2> &dofs() const
DOFs per node.
-
inline const array_type::tensor<size_t, 2> &conn() const
Connectivity.
-
inline std::array<size_t, 3> shape_elemmat() const
Shape of “elemmat”.
-
template<class T>
inline void assemble(const T &elemmat) Assemble from “elemmat”.
-
inline array_type::tensor<double, 2> Todense() const
Copy as dense matrix.
-
template<class T>
inline void todense(T &ret) const Copy to dense matrix.
-
inline array_type::tensor<double, 2> Dot(const array_type::tensor<double, 2> &x) const
Dot-product \( b_i = A_{ij} x_j \).
-
inline array_type::tensor<double, 1> Dot(const array_type::tensor<double, 1> &x) const
Dot-product \( b_i = A_{ij} x_j \).
-
inline void dot(const array_type::tensor<double, 2> &x, array_type::tensor<double, 2> &b) const
Dot-product \( b_i = A_{ij} x_j \).
-
inline void dot(const array_type::tensor<double, 1> &x, array_type::tensor<double, 1> &b) const
Dot-product \( b_i = A_{ij} x_j \).
-
using derived_type = D
-
class MatrixDiagonal : public GooseFEM::MatrixBase<MatrixDiagonal>, public GooseFEM::MatrixDiagonalBase<MatrixDiagonal>
- #include <GooseFEM/MatrixDiagonal.h>
Diagonal matrix.
Warning: assemble() ignores all off-diagonal terms.
See Vector() for bookkeeping definitions.
Public Functions
-
template<class C, class D>
inline MatrixDiagonal(const C &conn, const D &dofs) Constructor.
- Template Parameters
C – e.g.
array_type::tensor<size_t, 2>
D – e.g.
array_type::tensor<size_t, 2>
- Parameters
-
inline void set(const array_type::tensor<double, 1> &A)
Set all (diagonal) matrix components.
- Parameters
A – The matrix [ndof].
-
inline const array_type::tensor<double, 1> &Todiagonal() const
Copy as diagonal matrix.
- Returns
[ndof].
-
inline const array_type::tensor<double, 1> &data() const
Underlying matrix.
- Returns
[ndof].
-
template<class C, class D>
-
template<class D>
class MatrixDiagonalBase - #include <GooseFEM/MatrixDiagonal.h>
CRTP base class for a partitioned matrix with tying.
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
inline array_type::tensor<double, 2> Solve(const array_type::tensor<double, 2> &b)
Solve \( x = A^{-1} b \).
Note that this does not involve a conversion to DOFs.
In case of GooseFEM::MatrixDiagonalPartitioned under the hood, schematically: \( x_u = A_{uu}^{-1} (b_u - A_{up} * x_p) \equiv A_{uu}^{-1} b_u \) (again, no conversion to DOFs is needed). Use GooseFEM::MatrixDiagonalPartitioned::Reaction() to get reaction forces.
- Parameters
b – nodevec [nelem, ndim].
- Returns
x nodevec [nelem, ndim].
-
inline array_type::tensor<double, 1> Solve(const array_type::tensor<double, 1> &b)
Solve \( x = A^{-1} b \).
For GooseFEM::MatrixDiagonalPartitioned under the hood solved \( x_u = A_{uu}^{-1} (b_u - A_{up} * x_p) \equiv A_{uu}^{-1} b_u \). Use GooseFEM::MatrixDiagonalPartitioned::Reaction() to get reaction forces.
- Parameters
b – dofval [ndof].
- Returns
x dofval [ndof].
-
inline void solve(const array_type::tensor<double, 2> &b, array_type::tensor<double, 2> &x)
Solve \( x = A^{-1} b \).
For GooseFEM::MatrixDiagonalPartitioned under the hood solved \( x_u = A_{uu}^{-1} (b_u - A_{up} * x_p) \equiv A_{uu}^{-1} b_u \). Use GooseFEM::MatrixDiagonalPartitioned::Reaction() to get reaction forces.
- Parameters
b – nodevec [nelem, ndim].
x – (overwritten) nodevec [nelem, ndim].
-
inline void solve(const array_type::tensor<double, 1> &b, array_type::tensor<double, 1> &x)
Solve \( x = A^{-1} b \).
For GooseFEM::MatrixDiagonalPartitioned under the hood solved \( x_u = A_{uu}^{-1} (b_u - A_{up} * x_p) \equiv A_{uu}^{-1} b_u \). Use GooseFEM::MatrixDiagonalPartitioned::Reaction() to get reaction forces.
- Parameters
b – nodevec [nelem, ndim].
x – (overwritten) nodevec [nelem, ndim].
-
using derived_type = D
-
class MatrixDiagonalPartitioned : public GooseFEM::MatrixPartitionedBase<MatrixDiagonalPartitioned>, public GooseFEM::MatrixDiagonalBase<MatrixDiagonalPartitioned>
- #include <GooseFEM/MatrixDiagonalPartitioned.h>
Diagonal and partitioned matrix.
See Vector() for bookkeeping definitions.
Public Functions
-
inline MatrixDiagonalPartitioned(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<size_t, 2> &dofs, const array_type::tensor<size_t, 1> &iip)
Constructor.
-
inline void set(const array_type::tensor<double, 1> &A)
Set all (diagonal) matrix components.
- Parameters
A – The matrix [ndof].
-
inline array_type::tensor<double, 1> data() const
Assemble to diagonal matrix (involves copies).
- Returns
[ndof].
-
inline const array_type::tensor<double, 1> &data_uu() const
Pointer to data.
- Returns
[nnu].
-
inline const array_type::tensor<double, 1> &data_pp() const
Pointer to data.
- Returns
[nnu].
-
inline array_type::tensor<double, 1> Todiagonal() const
Pointer to data.
- Returns
[nnu].
-
inline array_type::tensor<double, 1> Dot_u(const array_type::tensor<double, 1> &x_u, const array_type::tensor<double, 1> &x_p) const
- Todo:
Decide if this function should be kept.
-
inline void dot_u(const array_type::tensor<double, 1> &x_u, const array_type::tensor<double, 1> &x_p, array_type::tensor<double, 1> &b_u) const
- Todo:
Decide if this function should be kept.
-
inline array_type::tensor<double, 1> Dot_p(const array_type::tensor<double, 1> &x_u, const array_type::tensor<double, 1> &x_p) const
- Todo:
Decide if this function should be kept.
-
inline void dot_p(const array_type::tensor<double, 1> &x_u, const array_type::tensor<double, 1> &x_p, array_type::tensor<double, 1> &b_p) const
- Todo:
Decide if this function should be kept.
-
inline array_type::tensor<double, 1> Solve_u(const array_type::tensor<double, 1> &b_u, const array_type::tensor<double, 1> &x_p)
-
inline void solve_u(const array_type::tensor<double, 1> &b_u, const array_type::tensor<double, 1> &x_p, array_type::tensor<double, 1> &x_u)
-
inline MatrixDiagonalPartitioned(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<size_t, 2> &dofs, const array_type::tensor<size_t, 1> &iip)
-
class MatrixPartitioned : public GooseFEM::MatrixPartitionedBase<MatrixPartitioned>
- #include <GooseFEM/MatrixPartitioned.h>
Sparse matrix partitioned in an unknown and a prescribed part.
In particular: \( \begin{bmatrix} A_{uu} & A_{up} \\ A_{pu} & A_{pp} \end{bmatrix} \)
See VectorPartitioned() for bookkeeping definitions.
Public Functions
-
inline MatrixPartitioned(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<size_t, 2> &dofs, const array_type::tensor<size_t, 1> &iip)
Constructor.
-
inline const Eigen::SparseMatrix<double> &data_uu() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_up() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_pu() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_pp() const
Pointer to data.
-
inline void set(const array_type::tensor<size_t, 1> &rows, const array_type::tensor<size_t, 1> &cols, const array_type::tensor<double, 2> &matrix)
Overwrite matrix.
- Parameters
rows – Row numbers [m].
cols – Column numbers [n].
matrix – Data entries
matrix(i, j)
forrows(i), cols(j)
[m, n].
-
inline void add(const array_type::tensor<size_t, 1> &rows, const array_type::tensor<size_t, 1> &cols, const array_type::tensor<double, 2> &matrix)
Add matrix.
- Parameters
rows – Row numbers [m].
cols – Column numbers [n].
matrix – Data entries
matrix(i, j)
forrows(i), cols(j)
[m, n].
-
inline MatrixPartitioned(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<size_t, 2> &dofs, const array_type::tensor<size_t, 1> &iip)
-
template<class D>
class MatrixPartitionedBase : public GooseFEM::MatrixBase<D> - #include <GooseFEM/Matrix.h>
CRTP base class for a partitioned matrix.
Subclassed by GooseFEM::MatrixPartitionedTyingsBase< MatrixPartitionedTyings >, GooseFEM::MatrixPartitionedTyingsBase< D >
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
inline size_t nnu() const
Number of unknown DOFs.
- Returns
Unsigned integer.
-
inline size_t nnp() const
Number of prescribed DOFs.
- Returns
Unsigned integer.
-
inline const array_type::tensor<size_t, 1> &iiu() const
Unknown DOFs.
- Returns
[nnu].
-
inline const array_type::tensor<size_t, 1> &iip() const
Prescribed DOFs.
- Returns
[nnp].
-
inline array_type::tensor<double, 2> Reaction(const array_type::tensor<double, 2> &x, const array_type::tensor<double, 2> &b) const
Right-hand-size for corresponding to the prescribed DOFs:
\( b_p = A_{pu} * x_u + A_{pp} * x_p \)
and assemble them to the appropriate places in “nodevec”.
-
inline array_type::tensor<double, 1> Reaction(const array_type::tensor<double, 1> &x, const array_type::tensor<double, 1> &b) const
Same as Reaction(const array_type::tensor<double, 2>&, const array_type::tensor<double, 2>&), but of “dofval” input and output.
-
inline void reaction(const array_type::tensor<double, 2> &x, array_type::tensor<double, 2> &b) const
Same as Reaction(const array_type::tensor<double, 2>&, const array_type::tensor<double, 2>&), but inserting in-place.
-
inline void reaction(const array_type::tensor<double, 1> &x, array_type::tensor<double, 1> &b) const
Same as Reaction(const array_type::tensor<double, 1>&, const array_type::tensor<double, 1>&), but inserting in-place.
-
inline array_type::tensor<double, 1> Reaction_p(const array_type::tensor<double, 1> &x_u, const array_type::tensor<double, 1> &x_p) const
Same as Reaction(const array_type::tensor<double, 1>&, const array_type::tensor<double, 1>&), but with partitioned input and output.
-
inline void reaction_p(const array_type::tensor<double, 1> &x_u, const array_type::tensor<double, 1> &x_p, array_type::tensor<double, 1> &b_p) const
Same as Reaction_p(const array_type::tensor<double, 1>&, const array_type::tensor<double,
1>&), but writing to preallocated output.
-
using derived_type = D
-
template<class Solver = Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>
class MatrixPartitionedSolver : public MatrixSolverBase<MatrixPartitionedSolver<Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>>, public MatrixSolverPartitionedBase<MatrixPartitionedSolver<Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>> - #include <GooseFEM/MatrixPartitioned.h>
Solve \( x_u = A_{uu}^{-1} (b_u - A_{up} * x_p) \) for
A
of the MatrixPartitioned() class.You can solve for multiple right-hand-sides using one factorisation.
For “nodevec” input
x
is used to read \( x_p \), while \( x_u \) is written. See MatrixPartitioned::Reaction() to get \( b_p \).Public Functions
-
template<class M>
inline array_type::tensor<double, 1> Solve_u(M &A, const array_type::tensor<double, 1> &b_u, const array_type::tensor<double, 1> &x_p) Solve \( x = A^{-1} b \).
- Parameters
A – GooseFEM (sparse) matrix, see e.g. GooseFEM::MatrixPartitioned().
b_u – unknown dofval [nnu].
x_p – prescribed dofval [nnp]
- Returns
x_u unknown dofval [nnu].
-
template<class M>
inline void solve_u(M &A, const array_type::tensor<double, 1> &b_u, const array_type::tensor<double, 1> &x_p, array_type::tensor<double, 1> &x_u) Same as Solve \( x = A^{-1} b \).
- Parameters
A – GooseFEM (sparse) matrix, see e.g. GooseFEM::MatrixPartitioned().
b_u – unknown dofval [nnu].
x_p – prescribed dofval [nnp]
x_u – (overwritten) unknown dofval [nnu].
-
template<class M>
-
class MatrixPartitionedTyings : public GooseFEM::MatrixPartitionedTyingsBase<MatrixPartitionedTyings>
- #include <GooseFEM/MatrixPartitionedTyings.h>
Sparse matrix from with dependent DOFs are eliminated, and the remaining (small) independent system is partitioned in an unknown and a prescribed part.
In particular:
\( A_{ii} = \begin{bmatrix} A_{uu} & A_{up} \\ A_{pu} & A_{pp} \end{bmatrix} \)
See VectorPartitionedTyings() for bookkeeping definitions.
Public Functions
-
inline MatrixPartitionedTyings(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<size_t, 2> &dofs, const Eigen::SparseMatrix<double> &Cdu, const Eigen::SparseMatrix<double> &Cdp)
Constructor.
- Parameters
Cdu – See Tyings::Periodic::Cdu().
Cdp – See Tyings::Periodic::Cdp().
-
inline const Eigen::SparseMatrix<double> &data_uu() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_up() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_pu() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_pp() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_ud() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_pd() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_du() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_dp() const
Pointer to data.
-
inline const Eigen::SparseMatrix<double> &data_dd() const
Pointer to data.
-
inline MatrixPartitionedTyings(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<size_t, 2> &dofs, const Eigen::SparseMatrix<double> &Cdu, const Eigen::SparseMatrix<double> &Cdp)
-
template<class D>
class MatrixPartitionedTyingsBase : public GooseFEM::MatrixPartitionedBase<D> - #include <GooseFEM/Matrix.h>
CRTP base class for a partitioned matrix with tying.
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
inline size_t nni() const
Number of independent DOFs.
- Returns
Unsigned integer.
-
inline size_t nnd() const
Number of dependent DOFs.
- Returns
Unsigned integer.
-
inline const array_type::tensor<size_t, 1> &iii() const
Independent DOFs.
- Returns
[nnu].
-
inline const array_type::tensor<size_t, 1> &iid() const
Dependent DOFs.
- Returns
[nnp].
-
using derived_type = D
-
template<class Solver = Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>
class MatrixPartitionedTyingsSolver : public MatrixSolverBase<MatrixPartitionedTyingsSolver<Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>>, public MatrixSolverPartitionedBase<MatrixPartitionedTyingsSolver<Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>> - #include <GooseFEM/MatrixPartitionedTyings.h>
Solver for MatrixPartitionedTyings().
This solver class can be used to solve for multiple right-hand-sides using one factorisation.
Solving proceeds as follows:
\( A' = A_{ii} + A_{id} * C_{di} + C_{di}^T * A_{di} + C_{di}^T * A_{dd} * C_{di} \)
\( b' = b_i + C_{di}^T * b_d \)
\( x_u = A'_{uu} \ ( b'_u - A'_{up} * x_p ) \)
\( x_i = \begin{bmatrix} x_u \\ x_p \end{bmatrix} \)
\( x_d = C_{di} * x_i \)
Public Functions
-
inline array_type::tensor<double, 1> Solve_u(MatrixPartitionedTyings &A, const array_type::tensor<double, 1> &b_u, const array_type::tensor<double, 1> &b_d, const array_type::tensor<double, 1> &x_p)
Same as Solve(MatrixPartitionedTyings&, const array_type::tensor<double, 2>&, const
array_type::tensor<double, 2>&), but with partitioned input and output.
- Parameters
A – sparse matrix, see MatrixPartitionedTyings().
b_u – unknown dofval [nnu].
b_d – dependent dofval [nnd].
x_p – prescribed dofval [nnp]
- Returns
x_u unknown dofval [nnu].
-
inline void solve_u(MatrixPartitionedTyings &A, const array_type::tensor<double, 1> &b_u, const array_type::tensor<double, 1> &b_d, const array_type::tensor<double, 1> &x_p, array_type::tensor<double, 1> &x_u)
Same as Solve_u(MatrixPartitionedTyings&, const array_type::tensor<double, 1>&, constarray_type::tensor<double, 1>&, const array_type::tensor<double, 1>&), but writing to pre-allocated output.
- Parameters
A – sparse matrix, see MatrixPartitionedTyings().
b_u – unknown dofval [nnu].
b_d – dependent dofval [nnd].
x_p – prescribed dofval [nnp]
x_u – (overwritten) unknown dofval [nnu].
-
inline array_type::tensor<double, 1> Solve_u(MatrixPartitionedTyings &A, const array_type::tensor<double, 1> &b_u, const array_type::tensor<double, 1> &b_d, const array_type::tensor<double, 1> &x_p)
-
template<class Solver = Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>
class MatrixSolver : public MatrixSolverBase<MatrixSolver<Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>>, public MatrixSolverSingleBase<MatrixSolver<Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>> - #include <GooseFEM/Matrix.h>
Solve \( x = A^{-1} b \), for
A
of the GooseFEM::Matrix() class.You can solve for multiple right-hand-sides using one factorisation.
-
template<class D>
class MatrixSolverBase - #include <GooseFEM/Matrix.h>
CRTP base class for a solver class.
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
template<class M>
inline void solve(M &A, const array_type::tensor<double, 2> &b, array_type::tensor<double, 2> &x) Solve \( x = A^{-1} b \).
- Parameters
A – GooseFEM (sparse) matrix, see e.g. GooseFEM::Matrix().
b – nodevec [nelem, ndim].
x – (overwritten) nodevec [nelem, ndim].
-
template<class M>
inline void solve(M &A, const array_type::tensor<double, 1> &b, array_type::tensor<double, 1> &x) Solve \( x = A^{-1} b \).
- Parameters
A – GooseFEM (sparse) matrix, see e.g. GooseFEM::Matrix().
b – dofval [ndof].
x – (overwritten) dofval [ndof].
-
using derived_type = D
-
template<class D>
class MatrixSolverPartitionedBase - #include <GooseFEM/Matrix.h>
CRTP base class for a extra functions for a partitioned solver class.
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
template<class M>
inline array_type::tensor<double, 2> Solve(M &A, const array_type::tensor<double, 2> &b, const array_type::tensor<double, 2> &x) Solve \( x = A^{-1} b \).
- Parameters
A – GooseFEM (sparse) matrix, see e.g. GooseFEM::Matrix().
b – nodevec [nelem, ndim].
x – nodevec [nelem, ndim].
- Returns
x nodevec [nelem, ndim].
-
template<class M>
inline array_type::tensor<double, 1> Solve(M &A, const array_type::tensor<double, 1> &b, const array_type::tensor<double, 1> &x) Solve \( x = A^{-1} b \).
- Parameters
A – GooseFEM (sparse) matrix, see e.g. GooseFEM::Matrix().
b – dofval [ndof].
x – nodevec [nelem, ndim].
- Returns
x dofval [ndof].
-
using derived_type = D
-
template<class D>
class MatrixSolverSingleBase - #include <GooseFEM/Matrix.h>
CRTP base class for a solver class.
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
template<class M>
inline array_type::tensor<double, 2> Solve(M &A, const array_type::tensor<double, 2> &b) Solve \( x = A^{-1} b \).
- Parameters
A – GooseFEM (sparse) matrix, see e.g. GooseFEM::Matrix().
b – nodevec [nelem, ndim].
- Returns
x nodevec [nelem, ndim].
-
template<class M>
inline array_type::tensor<double, 1> Solve(M &A, const array_type::tensor<double, 1> &b) Solve \( x = A^{-1} b \).
- Parameters
A – GooseFEM (sparse) matrix, see e.g. GooseFEM::Matrix().
b – dofval [ndof].
- Returns
x dofval [ndof].
-
using derived_type = D
-
class Vector
- #include <GooseFEM/Vector.h>
Class to switch between storage types.
In particular:
Subclassed by GooseFEM::VectorPartitioned, GooseFEM::VectorPartitionedTyings
Public Functions
-
inline size_t nelem() const
- Returns
Number of elements.
-
inline size_t nne() const
- Returns
Number of nodes per element.
-
inline size_t nnode() const
- Returns
Number of nodes.
-
inline size_t ndim() const
- Returns
Number of dimensions.
-
inline size_t ndof() const
- Returns
Number of DOFs.
-
inline const array_type::tensor<size_t, 2> &conn() const
-
inline const array_type::tensor<size_t, 2> &dofs() const
-
template<class T>
inline T Copy(const T &nodevec_src, const T &nodevec_dest) const Copy “nodevec” to another “nodevec”.
-
template<class T>
inline void copy(const T &nodevec_src, T &nodevec_dest) const Copy “nodevec” to another “nodevec”.
-
template<class T>
inline array_type::tensor<double, 1> AsDofs(const T &arg) const Convert “nodevec” or “elemvec” to “dofval” (overwrite entries that occur more than once).
-
template<class T, class R>
inline void asDofs(const T &arg, R &ret) const Convert “nodevec” or “elemvec” to “dofval” (overwrite entries that occur more than once).
-
template<class T>
inline array_type::tensor<double, 2> AsNode(const T &arg) const Convert “dofval” or “elemvec” to “nodevec” (overwrite entries that occur more than once).
-
template<class T, class R>
inline void asNode(const T &arg, R &ret) const Convert “dofval” or “elemvec” to “nodevec” (overwrite entries that occur more than once).
-
template<class T>
inline array_type::tensor<double, 3> AsElement(const T &arg) const Convert “dofval” or “nodevec” to “elemvec” (overwrite entries that occur more than once).
-
template<class T, class R>
inline void asElement(const T &arg, R &ret) const Convert “dofval” or “nodevec” to “elemvec” (overwrite entries that occur more than once).
-
template<class T>
inline array_type::tensor<double, 1> AssembleDofs(const T &arg) const Assemble “nodevec” or “elemvec” to “dofval” (adds entries that occur more that once).
-
template<class T, class R>
inline void assembleDofs(const T &arg, R &ret) const Assemble “nodevec” or “elemvec” to “dofval” (adds entries that occur more that once).
-
template<class T>
inline array_type::tensor<double, 2> AssembleNode(const T &arg) const Assemble “elemvec” to “nodevec” (adds entries that occur more that once.
-
template<class T, class R>
inline void assembleNode(const T &arg, R &ret) const Assemble “elemvec” to “nodevec” (adds entries that occur more that once.
-
inline std::array<size_t, 3> shape_elemmat() const
Shape of “elemmat”.
-
inline array_type::tensor<double, 1> allocate_dofval() const
Allocated “dofval”.
- Returns
[ndof]
-
inline array_type::tensor<double, 1> allocate_dofval(double val) const
Allocated and initialised “dofval”.
- Parameters
val – value to which to initialise.
- Returns
[ndof]
-
inline array_type::tensor<double, 2> allocate_nodevec() const
Allocated “nodevec”.
-
inline array_type::tensor<double, 2> allocate_nodevec(double val) const
Allocated and initialised “nodevec”.
-
inline array_type::tensor<double, 3> allocate_elemvec() const
Allocated “elemvec”.
-
inline array_type::tensor<double, 3> allocate_elemvec(double val) const
Allocated and initialised “elemvec”.
-
inline array_type::tensor<double, 3> allocate_elemmat() const
Allocated “elemmat”.
-
inline size_t nelem() const
-
class VectorPartitioned : public GooseFEM::Vector
- #include <GooseFEM/VectorPartitioned.h>
Class to switch between storage types, based on a mesh and DOFs that are partitioned in:
To this end some internal re-ordering of the DOFs has to be done, as follows:
which is relevant only if you interact using partitioned DOF-lists (“dofval_u” or “dofval_p”).iiu() -> arange(nnu()) iip() -> nnu() + arange(nnp())
The “dofval”, “nodevec”, and “elemvec” are all stored in the ‘normal’ order.
For reference:
Public Functions
-
inline VectorPartitioned(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<size_t, 2> &dofs, const array_type::tensor<size_t, 1> &iip)
Constructor.
-
inline size_t nnu() const
- Returns
Number of unknown DOFs.
-
inline size_t nnp() const
- Returns
Number of prescribed DOFs.
-
inline const array_type::tensor<size_t, 1> &iiu() const
- Returns
Unknown DOFs [nnu].
-
inline const array_type::tensor<size_t, 1> &iip() const
- Returns
Prescribed DOFs [nnp].
-
inline array_type::tensor<bool, 2> dofs_is_u() const
Per DOF (see Vector::dofs()) list if unknown (“u”).
- Returns
Boolean “nodevec”.
-
inline array_type::tensor<bool, 2> dofs_is_p() const
Per DOF (see Vector::dofs()) list if prescribed (“p”).
- Returns
Boolean “nodevec”.
-
inline array_type::tensor<double, 2> Copy_u(const array_type::tensor<double, 2> &nodevec_src, const array_type::tensor<double, 2> &nodevec_dest) const
Copy unknown DOFs from “nodevec” to another “nodvec”:
the other DOFs are taken fromnodevec_dest[vector.dofs_is_u()] = nodevec_src
nodevec_dest
:nodevec_dest[vector.dofs_is_p()] = nodevec_dest
-
inline void copy_u(const array_type::tensor<double, 2> &nodevec_src, array_type::tensor<double, 2> &nodevec_dest) const
Copy unknown DOFs from “nodevec” to another “nodvec”:
the other DOFs are taken fromnodevec_dest[vector.dofs_is_u()] = nodevec_src
nodevec_dest
:nodevec_dest[vector.dofs_is_p()] = nodevec_dest
-
inline array_type::tensor<double, 2> Copy_p(const array_type::tensor<double, 2> &nodevec_src, const array_type::tensor<double, 2> &nodevec_dest) const
Copy prescribed DOFs from “nodevec” to another “nodvec”:
the other DOFs are taken fromnodevec_dest[vector.dofs_is_p()] = nodevec_src
nodevec_dest
:nodevec_dest[vector.dofs_is_u()] = nodevec_dest
-
inline void copy_p(const array_type::tensor<double, 2> &nodevec_src, array_type::tensor<double, 2> &nodevec_dest) const
Copy prescribed DOFs from “nodevec” to another “nodvec”:
the other DOFs are taken fromnodevec_dest[vector.dofs_is_p()] = nodevec_src
nodevec_dest
:nodevec_dest[vector.dofs_is_u()] = nodevec_dest
-
inline array_type::tensor<double, 1> DofsFromParitioned(const array_type::tensor<double, 1> &dofval_u, const array_type::tensor<double, 1> &dofval_p) const
Combine unknown and prescribed “dofval” into a single “dofval” list.
-
inline void dofsFromParitioned(const array_type::tensor<double, 1> &dofval_u, const array_type::tensor<double, 1> &dofval_p, array_type::tensor<double, 1> &dofval) const
Combine unknown and prescribed “dofval” into a single “dofval” list.
-
inline array_type::tensor<double, 2> NodeFromPartitioned(const array_type::tensor<double, 1> &dofval_u, const array_type::tensor<double, 1> &dofval_p) const
Combine unknown and prescribed “dofval” into a single “dofval” list and directly convert to “nodeval” without a temporary (overwrite entries that occur more than once).
-
inline void nodeFromPartitioned(const array_type::tensor<double, 1> &dofval_u, const array_type::tensor<double, 1> &dofval_p, array_type::tensor<double, 2> &nodevec) const
Combine unknown and prescribed “dofval” into a single “dofval” list and directly convert to “nodeval” without a temporary (overwrite entries that occur more than once).
-
inline array_type::tensor<double, 3> ElementFromPartitioned(const array_type::tensor<double, 1> &dofval_u, const array_type::tensor<double, 1> &dofval_p) const
Combine unknown and prescribed “dofval” into a single “dofval” list and directly convert to “elemvec” without a temporary (overwrite entries that occur more than once).
-
inline void elementFromPartitioned(const array_type::tensor<double, 1> &dofval_u, const array_type::tensor<double, 1> &dofval_p, array_type::tensor<double, 3> &elemvec) const
Combine unknown and prescribed “dofval” into a single “dofval” list and directly convert to “elemvec” without a temporary (overwrite entries that occur more than once).
-
inline array_type::tensor<double, 1> AsDofs_u(const array_type::tensor<double, 1> &dofval) const
Extract the unknown “dofval”:
dofval[iiu()]
-
inline void asDofs_u(const array_type::tensor<double, 1> &dofval, array_type::tensor<double, 1> &dofval_u) const
Extract the unknown “dofval”:
dofval[iiu()]
-
inline array_type::tensor<double, 1> AsDofs_u(const array_type::tensor<double, 2> &nodevec) const
Convert “nodevec” to “dofval” (overwrite entries that occur more than once) and extract the unknown “dofval” without a temporary.
-
inline void asDofs_u(const array_type::tensor<double, 2> &nodevec, array_type::tensor<double, 1> &dofval_u) const
Convert “nodevec” to “dofval” (overwrite entries that occur more than once) and extract the unknown “dofval” without a temporary.
-
inline array_type::tensor<double, 1> AsDofs_u(const array_type::tensor<double, 3> &elemvec) const
Convert “elemvec” to “dofval” (overwrite entries that occur more than once) and extract the unknown “dofval” without a temporary.
-
inline void asDofs_u(const array_type::tensor<double, 3> &elemvec, array_type::tensor<double, 1> &dofval_u) const
Convert “elemvec” to “dofval” (overwrite entries that occur more than once) and extract the unknown “dofval” without a temporary.
-
inline array_type::tensor<double, 1> AsDofs_p(const array_type::tensor<double, 1> &dofval) const
Extract the prescribed “dofval”:
dofval[iip()]
-
inline void asDofs_p(const array_type::tensor<double, 1> &dofval, array_type::tensor<double, 1> &dofval_p) const
Extract the prescribed “dofval”:
dofval[iip()]
-
inline array_type::tensor<double, 1> AsDofs_p(const array_type::tensor<double, 2> &nodevec) const
Convert “nodevec” to “dofval” (overwrite entries that occur more than once) and extract the prescribed “dofval” without a temporary.
-
inline void asDofs_p(const array_type::tensor<double, 2> &nodevec, array_type::tensor<double, 1> &dofval_p) const
Convert “nodevec” to “dofval” (overwrite entries that occur more than once) and extract the prescribed “dofval” without a temporary.
-
inline array_type::tensor<double, 1> AsDofs_p(const array_type::tensor<double, 3> &elemvec) const
Convert “elemvec” to “dofval” (overwrite entries that occur more than once) and extract the prescribed “dofval” without a temporary.
-
inline void asDofs_p(const array_type::tensor<double, 3> &elemvec, array_type::tensor<double, 1> &dofval_p) const
Convert “elemvec” to “dofval” (overwrite entries that occur more than once) and extract the prescribed “dofval” without a temporary.
-
inline VectorPartitioned(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<size_t, 2> &dofs, const array_type::tensor<size_t, 1> &iip)
-
class VectorPartitionedTyings : public GooseFEM::Vector
- #include <GooseFEM/VectorPartitionedTyings.h>
Class to switch between storage types.
In particular:
”elemvec”: nodal vectors stored per element [nelem, nne, ndim].
”dofval”: DOF values [ndof].
”dofval_u”: DOF values (Unknown),
== dofval[iiu]
, [nnu].”dofval_p”: DOF values (Prescribed),
== dofval[iiu]
, [nnp].
Public Functions
-
template<class E, class M>
inline VectorPartitionedTyings(const E &conn, const E &dofs, const M &Cdu, const M &Cdp, const M &Cdi) Constructor.
- Template Parameters
E – e.g.
array_type::tensor<size_t, 2>
M – e.g.
Eigen::SparseMatrix<double>
- Parameters
Cdu – See Tyings::Periodic::Cdu().
Cdp – See Tyings::Periodic::Cdp().
Cdi – See Tyings::Periodic::Cdi().
-
inline size_t nnd() const
- Returns
Number of dependent DOFs.
-
inline size_t nni() const
- Returns
Number of independent DOFs.
-
inline size_t nnu() const
- Returns
Number of independent unknown DOFs.
-
inline size_t nnp() const
- Returns
Number of independent prescribed DOFs.
-
inline const array_type::tensor<size_t, 1> &iid() const
Dependent DOFs (list of DOF numbers).
- Returns
Pointer.
-
inline const array_type::tensor<size_t, 1> &iii() const
Independent DOFs (list of DOF numbers).
- Returns
Copy.
-
inline const array_type::tensor<size_t, 1> &iiu() const
Independent unknown DOFs (list of DOF numbers).
- Returns
Pointer.
-
inline const array_type::tensor<size_t, 1> &iip() const
Independent prescribed DOFs (list of DOF numbers).
- Returns
Pointer.
-
template<class T>
inline void copy_p(const T &dofval_src, T &dofval_dest) const Copy (part of) “dofval” to another “dofval”: dofval_dest[iip()] = dofval_src[iip()].
-
template<class T>
inline array_type::tensor<double, 1> AsDofs_i(const T &nodevec) const Convert to “dofval” (overwrite entries that occur more than once).
Only the independent DOFs are retained.
-
namespace array_type
Container type.
By default
array_type::tensor
is used. Otherwise:#define GOOSEFEM_USE_XTENSOR_PYTHON
to usext::pytensor
-
namespace Element
Element quadrature and interpolation.
Functions
-
inline array_type::tensor<double, 3> asElementVector(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<double, 2> &nodevec)
Convert nodal vector with (“nodevec”, shape:
[nnode, ndim]
) to nodal vector stored per element (“elemvec”, shape:[nelem, nne, ndim]
).- Parameters
conn – Connectivity.
nodevec – “nodevec”.
- Returns
“elemvec”.
-
inline array_type::tensor<double, 2> assembleNodeVector(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<double, 3> &elemvec)
Assemble nodal vector stored per element (“elemvec”, shape
[nelem, nne, ndim]
) to nodal vector (“nodevec”, shape[nnode, ndim]
).- Parameters
conn – Connectivity.
elemvec – “elemvec”.
- Returns
“nodevec”.
-
template<class E>
inline bool isSequential(const E &dofs) Check that DOFs leave no holes.
- Parameters
dofs – DOFs (“nodevec”)
- Returns
true
if there are no holds.
-
bool isDiagonal(const array_type::tensor<double, 3> &elemmat)
Check that all of the matrices stored per elemmat (shape:
[nelem, nne * ndim, nne * ndim]
) are diagonal.- Parameters
elemmat – Element-vectors (“elemmat”)
- Returns
true
if all element matrices are diagonal.
-
template<class D>
class QuadratureBase - #include <GooseFEM/Element.h>
CRTP base class for quadrature.
Subclassed by GooseFEM::Element::QuadratureBaseCartesian< Quadrature >, GooseFEM::Element::QuadratureBaseCartesian< QuadratureAxisymmetric >, GooseFEM::Element::QuadratureBaseCartesian< QuadraturePlanar >, GooseFEM::Element::QuadratureBaseCartesian< D >
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
inline auto nelem() const
Number of elements.
- Returns
Scalar.
-
inline auto nne() const
Number of nodes per element.
- Returns
Scalar.
-
inline auto ndim() const
Number of dimensions for node vectors.
- Returns
Scalar.
-
inline auto tdim() const
Number of dimensions for integration point tensors.
- Returns
Scalar.
-
inline auto nip() const
Number of integration points.
- Returns
Scalar.
-
template<class T, class R>
inline void asTensor(const T &arg, R &ret) const Convert “qscalar” to “qtensor” of certain rank.
Fully allocated output passed as reference, use AsTensor to allocate and return data.
- Parameters
arg – A “qscalar”.
ret – A “qtensor”.
-
template<size_t rank, class T>
inline auto AsTensor(const T &arg) const Convert “qscalar” to “qtensor” of certain rank.
- Parameters
arg – A “qscalar”.
- Returns
“qtensor”.
-
template<class T>
inline auto AsTensor(size_t rank, const T &arg) const Convert “qscalar” to “qtensor” of certain rank.
- Parameters
rank – Tensor rank.
arg – A “qscalar”.
- Returns
“qtensor”.
-
inline auto shape_elemvec() const -> std::array<size_t, 3>
Get the shape of an “elemvec”.
-
inline auto shape_elemvec(size_t arg) const -> std::array<size_t, 3>
Get the shape of an “elemvec”.
-
inline auto shape_elemmat() const -> std::array<size_t, 3>
Get the shape of an “elemmat”.
-
template<size_t rank = 0>
inline auto shape_qtensor() const -> std::array<size_t, 2 + rank> Get the shape of a “qtensor” of a certain rank (0 = scalar, 1, vector, 2 = 2nd-order tensor, etc.).
-
inline auto shape_qtensor(size_t rank) const -> std::vector<size_t>
Get the shape of a “qtensor” of a certain rank (0 = scalar, 1, vector, 2 = 2nd-order tensor, etc.).
-
template<size_t trank>
inline auto shape_qtensor(size_t rank, size_t arg) const -> std::array<size_t, 2 + trank> Get the shape of a “qtensor” of a certain rank (0 = scalar, 1, vector, 2 = 2nd-order tensor, etc.).
-
inline auto shape_qtensor(size_t rank, size_t arg) const -> std::vector<size_t>
Get the shape of a “qtensor” of a certain rank (0 = scalar, 1, vector, 2 = 2nd-order tensor, etc.).
-
inline auto shape_qscalar() const -> std::array<size_t, 2>
Get the shape of a “qscalar” (a “qtensor” of rank 0)
-
inline auto shape_qvector() const -> std::array<size_t, 3>
Get the shape of a “qvector” (a “qtensor” of rank 1)
-
inline auto shape_qvector(size_t arg) const -> std::array<size_t, 3>
Get the shape of a “qvector” (a “qtensor” of rank 1)
-
template<class R>
inline auto allocate_elemvec() const Get an allocated
array_type::tensor
to store a “elemvec”.Note: the container is not (zero-)initialised.
- Template Parameters
R – value-type of the array, e.g.
double
.- Returns
xt::xarray
container of the correct shape.
-
template<class R>
inline auto allocate_elemvec(R val) const Get an allocated and initialised
xt::xarray
to store a “elemvec”.- Template Parameters
R – value-type of the array, e.g.
double
.- Parameters
val – The value to which to initialise all items.
- Returns
array_type::tensor
container of the correct shape.
-
template<class R>
inline auto allocate_elemmat() const Get an allocated
array_type::tensor
to store a “elemmat”.Note: the container is not (zero-)initialised.
- Template Parameters
R – value-type of the array, e.g.
double
.- Returns
xt::xarray
container of the correct shape.
-
template<class R>
inline auto allocate_elemmat(R val) const Get an allocated and initialised
xt::xarray
to store a “elemmat”.- Template Parameters
R – value-type of the array, e.g.
double
.- Parameters
val – The value to which to initialise all items.
- Returns
array_type::tensor
container of the correct shape.
-
template<size_t rank = 0, class R>
inline auto allocate_qtensor() const Get an allocated
array_type::tensor
to store a “qtensor” of a certain rank (0 = scalar, 1, vector, 2 = 2nd-order tensor, etc.).Default: rank = 0, a.k.a. scalar. Note: the container is not (zero-)initialised.
-
template<size_t rank = 0, class R>
inline auto allocate_qtensor(R val) const Get an allocated and initialised
array_type::tensor
to store a “qtensor” of a certain rank (0 = scalar, 1, vector, 2 = 2nd-order tensor, etc.).Default: rank = 0, a.k.a. scalar.
- Template Parameters
R – value-type of the array, e.g.
double
.- Parameters
val – The value to which to initialise all items.
- Returns
array_type::tensor
container of the correct shape (and rank).
-
template<class R>
inline auto allocate_qtensor(size_t rank) const Get an allocated
xt::xarray
to store a “qtensor” of a certain rank (0 = scalar, 1, vector, 2 = 2nd-order tensor, etc.).Note: the container is not (zero-)initialised.
- Template Parameters
R – value-type of the array, e.g.
double
.- Parameters
rank – The tensor rank.
- Returns
xt::xarray
container of the correct shape.
-
template<class R>
inline auto allocate_qtensor(size_t rank, R val) const Get an allocated and initialised
xt::xarray
to store a “qtensor” of a certain rank (0 = scalar, 1, vector, 2 = 2nd-order tensor, etc.).- Template Parameters
R – value-type of the array, e.g.
double
.- Parameters
rank – The tensor rank.
val – The value to which to initialise all items.
- Returns
array_type::tensor
container of the correct shape (and rank).
-
template<class R>
inline auto allocate_qscalar() const Get an allocated
array_type::tensor
to store a “qscalar” (a “qtensor” of rank 0).Note: the container is not (zero-)initialised.
- Template Parameters
R – value-type of the array, e.g.
double
.- Returns
xt::xarray
container of the correct shape.
-
template<class R>
inline auto allocate_qscalar(R val) const Get an allocated and initialised
xt::xarray
to store a “qscalar” (a “qtensor” of rank 0).- Template Parameters
R – value-type of the array, e.g.
double
.- Parameters
val – The value to which to initialise all items.
- Returns
array_type::tensor
container of the correct shape (and rank).
-
using derived_type = D
-
template<class D>
class QuadratureBaseCartesian : public GooseFEM::Element::QuadratureBase<D> - #include <GooseFEM/Element.h>
CRTP base class for interpolation and quadrature for a generic element in Cartesian coordinates.
Naming convention:
elemmat
: matrices stored per element, [nelem, nne * ndim, nne * ndim]elemvec
: nodal vectors stored per element, [nelem, nne, ndim]
Public Types
-
using derived_type = D
Underlying type.
Public Functions
-
template<class T>
inline void update_x(const T &x) Update the nodal positions.
This recomputes:
the shape functions,
the shape function gradients (in local and global) coordinates,
the integration points volumes. Under the small deformations assumption this function should not be called.
- Parameters
x – nodal coordinates (
elemvec
). Shape should match the earlier definition.
-
inline auto GradN() const -> const array_type::tensor<double, 4>&
Shape function gradients (in global coordinates).
-
inline auto dV() const -> const array_type::tensor<double, 2>&
Integration volume.
-
template<class T>
inline auto InterpQuad_vector(const T &elemvec) const -> array_type::tensor<double, 3> Interpolate element vector and evaluate at each quadrature point.
\( \vec{u}(\vec{x}_q) = N_i^e(\vec{x}) \vec{u}_i^e \)
-
template<class T, class R>
inline void interpQuad_vector(const T &elemvec, R &qvector) const Same as InterpQuad_vector(), but writing to preallocated return.
-
template<class T>
inline auto GradN_vector(const T &elemvec) const -> array_type::tensor<double, 4> Element-by-element: dyadic product of the shape function gradients and a nodal vector.
Typical input: nodal displacements. Typical output: quadrature point strains. Within one element::
Note that the functions and their gradients are precomputed upon construction, or updated when calling update_x().for e in range(nelem): for q in range(nip): for m in range(nne): qtensor(e, q, i, j) += dNdx(e, q, m, i) * elemvec(e, m, j)
-
template<class T, class R>
inline void gradN_vector(const T &elemvec, R &qtensor) const Same as GradN_vector(), but writing to preallocated return.
-
template<class T>
inline auto GradN_vector_T(const T &elemvec) const -> array_type::tensor<double, 4> The transposed output of GradN_vector().
Within one element::
for e in range(nelem): for q in range(nip): for m in range(nne): qtensor(e, q, j, i) += dNdx(e, q, m, i) * elemvec(e, m, j)
-
template<class T, class R>
inline void gradN_vector_T(const T &elemvec, R &qtensor) const Same as GradN_vector_T(), but writing to preallocated return.
-
template<class T>
inline auto SymGradN_vector(const T &elemvec) const -> array_type::tensor<double, 4> The symmetric output of GradN_vector().
Without one element::
for e in range(nelem): for q in range(nip): for m in range(nne): qtensor(e, q, i, j) += 0.5 * dNdx(e, q, m, i) * elemvec(e, m, j) qtensor(e, q, j, i) += 0.5 * dNdx(e, q, m, i) * elemvec(e, m, j)
-
template<class T, class R>
inline void symGradN_vector(const T &elemvec, R &qtensor) const Same as SymGradN_vector(), but writing to preallocated return.
-
template<class T>
inline auto Int_N_vector_dV(const T &qvector) const -> array_type::tensor<double, 3> Element-by-element: integral of a continuous vector-field.
\( \vec{f}_i^e = \int N_i^e(\vec{x}) \vec{f}(\vec{x}) d\Omega_e \)
which is integration numerically as follows
\( \vec{f}_i^e = \sum\limits_q N_i^e(\vec{x}_q) \vec{f}(\vec{x}_q) \)
-
template<class T, class R>
inline void int_N_vector_dV(const T &qvector, R &elemvec) const Same as Int_N_vector_dV(), but writing to preallocated return.
-
template<class T>
inline auto Int_N_scalar_NT_dV(const T &qscalar) const -> array_type::tensor<double, 3> Element-by-element: integral of the scalar product of the shape function with a scalar.
Within one one element::
withfor e in range(nelem): for q in range(nip): for m in range(nne): for n in range(nne): elemmat(e, m * ndim + i, n * ndim + i) += N(e, q, m) * qscalar(e, q) * N(e, q, n) * dV(e, q)
i
a tensor dimension. Note that the functions and their gradients are precomputed upon construction, or updated when calling update_x().
-
template<class T, class R>
inline void int_N_scalar_NT_dV(const T &qscalar, R &elemmat) const Same as Int_N_scalar_NT_dV(), but writing to preallocated return.
-
template<class T>
inline auto Int_gradN_dot_tensor2_dV(const T &qtensor) const -> array_type::tensor<double, 3> Element-by-element: integral of the dot product of the shape function gradients with a second order tensor.
Typical input: stress. Typical output: nodal force. Within one one element::
withfor e in range(nelem): for q in range(nip): for m in range(nne): elemvec(e, m, j) += dNdx(e, q, m, i) * qtensor(e, q, i, j) * dV(e, q)
i
andj
tensor dimensions. Note that the functions and their gradients are precomputed upon construction, or updated when calling update_x().
-
template<class T, class R>
inline void int_gradN_dot_tensor2_dV(const T &qtensor, R &elemvec) const Same as Int_gradN_dot_tensor2_dV(), but writing to preallocated return.
-
template<class T>
inline auto Int_gradN_dot_tensor4_dot_gradNT_dV(const T &qtensor) const -> array_type::tensor<double, 3> Element-by-element: integral of the dot products of the shape function gradients with a fourth order tensor.
Typical input: stiffness tensor. Typical output: stiffness matrix. Within one one element::
withfor e in range(nelem): for q in range(nip): for m in range(nne): for n in range(nne): elemmat(e, m * ndim + j, n * ndim + k) += dNdx(e,q,m,i) * qtensor(e,q,i,j,k,l) * dNdx(e,q,n,l) * dV(e,q)
i
,j
,k
, andl
tensor dimensions. Note that the functions and their gradients are precomputed upon construction, or updated when calling update_x().
-
template<class T, class R>
inline void int_gradN_dot_tensor4_dot_gradNT_dV(const T &qtensor, R &elemmat) const Same as Int_gradN_dot_tensor4_dot_gradNT_dV(), but writing to preallocated return.
-
namespace Hex8
8-noded hexahedral element in 3d (GooseFEM::Mesh::ElementType::Hex8).
-
class Quadrature : public GooseFEM::Element::QuadratureBaseCartesian<Quadrature>
- #include <GooseFEM/ElementHex8.h>
Interpolation and quadrature.
Fixed dimensions:
ndim = 3
: number of dimensions.nne = 8
: number of nodes per element.
Naming convention:
elemmat
: matrices stored per element, [nelem, nne * ndim, nne * ndim]elemvec
: nodal vectors stored per element, [nelem, nne, ndim]
Public Functions
-
template<class T>
inline Quadrature(const T &x) Constructor: use default Gauss integration.
The following is pre-computed during construction:
the shape functions,
the shape function gradients (in local and global) coordinates,
the integration points volumes. They can be reused without any cost. They only have to be recomputed when the nodal position changes (note that they are assumed to be constant under a small-strain assumption). In that case use update_x() to update the nodal positions and to recompute the above listed quantities.
- Parameters
x – nodal coordinates (
elemvec
).
-
template<class T, class X, class W>
inline Quadrature(const T &x, const X &xi, const W &w) Constructor with custom integration.
The following is pre-computed during construction:
the shape functions,
the shape function gradients (in local and global) coordinates,
the integration points volumes. They can be reused without any cost. They only have to be recomputed when the nodal position changes (note that they are assumed to be constant under a small-strain assumption). In that case use update_x() to update the nodal positions and to recompute the above listed quantities.
-
namespace Gauss
gauss quadrature: quadrature points such that integration is exact for these bi-linear elements::
Functions
-
inline size_t nip()
Number of integration points:
nip = nne = 8
- Returns
unsigned int
-
inline array_type::tensor<double, 2> xi()
Integration point coordinates (local coordinates).
- Returns
Coordinates [nip, ndim], with
ndim = 3
.
-
inline array_type::tensor<double, 1> w()
Integration point weights.
- Returns
Coordinates [nip].
-
inline size_t nip()
-
namespace Nodal
nodal quadrature: quadrature points coincide with the nodes.
The order is the same as in the connectivity.
Functions
-
inline size_t nip()
Number of integration points:
nip = nne = 8
- Returns
unsigned int
-
inline array_type::tensor<double, 2> xi()
Integration point coordinates (local coordinates).
- Returns
Coordinates [nip,
ndim
], withndim = 3
.
-
inline array_type::tensor<double, 1> w()
Integration point weights.
- Returns
Coordinates [nip].
-
inline size_t nip()
-
class Quadrature : public GooseFEM::Element::QuadratureBaseCartesian<Quadrature>
-
namespace Quad4
4-noded quadrilateral element in 2d (GooseFEM::Mesh::ElementType::Quad4).
-
class Quadrature : public GooseFEM::Element::QuadratureBaseCartesian<Quadrature>
- #include <GooseFEM/ElementQuad4.h>
Interpolation and quadrature.
Fixed dimensions:
ndim = 2
: number of dimensions.nne = 4
: number of nodes per element.
Naming convention:
elemmat
: matrices stored per element, [nelem, nne * ndim, nne * ndim]elemvec
: nodal vectors stored per element, [nelem, nne, ndim]
Public Functions
-
template<class T>
inline Quadrature(const T &x) Constructor: use default Gauss integration.
The following is pre-computed during construction:
the shape functions,
the shape function gradients (in local and global) coordinates,
the integration points volumes. They can be reused without any cost. They only have to be recomputed when the nodal position changes (note that they are assumed to be constant under a small-strain assumption). In that case use update_x() to update the nodal positions and to recompute the above listed quantities.
- Parameters
x – nodal coordinates (
elemvec
).
-
template<class T, class X, class W>
inline Quadrature(const T &x, const X &xi, const W &w) Constructor with custom integration.
The following is pre-computed during construction:
the shape functions,
the shape function gradients (in local and global) coordinates,
the integration points volumes. They can be reused without any cost. They only have to be recomputed when the nodal position changes (note that they are assumed to be constant under a small-strain assumption). In that case use update_x() to update the nodal positions and to recompute the above listed quantities.
-
class QuadratureAxisymmetric : public GooseFEM::Element::QuadratureBaseCartesian<QuadratureAxisymmetric>
- #include <GooseFEM/ElementQuad4Axisymmetric.h>
Interpolation and quadrature.
Fixed dimensions:
ndim = 2
: number of dimensions.tdim = 3
: number of dimensions or tensor.nne = 4
: number of nodes per element.
Naming convention:
elemmat
: matrices stored per element, [nelem, nne * ndim, nne * ndim]elemvec
: nodal vectors stored per element, [nelem, nne, ndim]
Public Functions
-
template<class T>
inline QuadratureAxisymmetric(const T &x) Constructor: use default Gauss integration.
The following is pre-computed during construction:
the shape functions,
the shape function gradients (in local and global) coordinates,
the integration points volumes. They can be reused without any cost. They only have to be recomputed when the nodal position changes (note that they are assumed to be constant under a small-strain assumption). In that case use update_x() to update the nodal positions and to recompute the above listed quantities.
- Parameters
x – nodal coordinates (
elemvec
).
-
template<class T, class X, class W>
inline QuadratureAxisymmetric(const T &x, const X &xi, const W &w) Constructor with custom integration.
The following is pre-computed during construction:
the shape functions,
the shape function gradients (in local and global) coordinates,
the integration points volumes. They can be reused without any cost. They only have to be recomputed when the nodal position changes (note that they are assumed to be constant under a small-strain assumption). In that case use update_x() to update the nodal positions and to recompute the above listed quantities.
-
inline const array_type::tensor<double, 6> &B() const
Get the B-matrix (shape function gradients) (in global coordinates).
Note that the functions and their gradients are precomputed upon construction, or updated when calling update_x().
-
class QuadraturePlanar : public GooseFEM::Element::QuadratureBaseCartesian<QuadraturePlanar>
- #include <GooseFEM/ElementQuad4Planar.h>
Interpolation and quadrature.
Similar to Element::Quad4::Quadrature with the only different that quadrature point tensors are 3d (“plane strain”) while the mesh is 2d.
Fixed dimensions:
ndim = 2
: number of dimensions.tdim = 3
: number of dimensions or tensor.nne = 4
: number of nodes per element.
Naming convention:
elemmat
: matrices stored per element, [nelem, nne * ndim, nne * ndim]elemvec
: nodal vectors stored per element, [nelem, nne, ndim]
Public Functions
-
template<class T>
inline QuadraturePlanar(const T &x, double thick = 1.0) Constructor: use default Gauss integration.
The following is pre-computed during construction:
the shape functions,
the shape function gradients (in local and global) coordinates,
the integration points volumes. They can be reused without any cost. They only have to be recomputed when the nodal position changes (note that they are assumed to be constant under a small-strain assumption). In that case use update_x() to update the nodal positions and to recompute the above listed quantities.
- Parameters
x – nodal coordinates (
elemvec
).thick – out-of-plane thickness (incorporated in the element volume).
-
template<class T, class X, class W>
inline QuadraturePlanar(const T &x, const X &xi, const W &w, double thick = 1.0) Constructor with custom integration.
The following is pre-computed during construction:
the shape functions,
the shape function gradients (in local and global) coordinates,
the integration points volumes. They can be reused without any cost. They only have to be recomputed when the nodal position changes (note that they are assumed to be constant under a small-strain assumption). In that case use update_x() to update the nodal positions and to recompute the above listed quantities.
-
namespace Gauss
Gauss quadrature: quadrature points such that integration is exact for this bi-linear element:
+ ----------- + | | | 3 2 | | | | 0 1 | | | + ----------- +
Functions
-
inline size_t nip()
Number of integration points:
nip = nne = 4
- Returns
unsigned int
-
inline array_type::tensor<double, 2> xi()
Integration point coordinates (local coordinates).
- Returns
Coordinates [nip,
ndim
], withndim = 2
.
-
inline array_type::tensor<double, 1> w()
Integration point weights.
- Returns
Coordinates [nip].
-
inline size_t nip()
-
namespace MidPoint
midpoint quadrature: quadrature points in the middle of the element::
+ ------- + | | | 0 | | | + ------- +
Functions
-
inline size_t nip()
Number of integration points::
nip = 1
- Returns
unsigned int
-
inline array_type::tensor<double, 2> xi()
Integration point coordinates (local coordinates).
- Returns
Coordinates [nip,
ndim
], withndim = 2
.
-
inline array_type::tensor<double, 1> w()
Integration point weights.
- Returns
Coordinates [nip].
-
inline size_t nip()
-
namespace Nodal
nodal quadrature: quadrature points coincide with the nodes.
The order is the same as in the connectivity:
3 -- 2 | | 0 -- 1
Functions
-
inline size_t nip()
Number of integration points::
nip = nne = 4
- Returns
unsigned int
-
inline array_type::tensor<double, 2> xi()
Integration point coordinates (local coordinates).
- Returns
Coordinates [nip,
ndim
], withndim = 2
.
-
inline array_type::tensor<double, 1> w()
Integration point weights.
- Returns
Coordinates [nip].
-
inline size_t nip()
-
class Quadrature : public GooseFEM::Element::QuadratureBaseCartesian<Quadrature>
-
inline array_type::tensor<double, 3> asElementVector(const array_type::tensor<size_t, 2> &conn, const array_type::tensor<double, 2> &nodevec)
-
namespace Iterate
Support function for iterations in end-user programs.
-
class StopList
- #include <GooseFEM/Iterate.h>
Class to perform a residual check based on the last “n” iterations.
A typical usage is in dynamic simulations where equilibrium is checked based on a force residual. Fluctuations could however be responsible for this criterion to be triggered too early. By checking several time-steps such case can be avoided.
Public Functions
-
inline StopList(size_t n = 1)
Constructor.
- Parameters
n – Number of consecutive iterations to consider.
-
inline void reset()
Reset all residuals to infinity.
-
inline void reset(size_t n)
Reset all residuals to infinity, and change the number of residuals to check.
- Parameters
n – Number of consecutive iterations to consider.
-
inline void roll_insert(double res)
Roll the list with the residuals, and add a new residual to the end.
In Python code this function corresponds to::
I.e. the residual ofresiduals = residuals[1:] + [new_residual]
n
iterations ago will be forgotten.- Parameters
res – New residual to add to the list of residuals.
-
inline bool descending() const
Check of the sequence of
n
residuals is in descending order.- Returns
true
if then
residuals are in descending order.
-
inline bool all_less(double tol) const
Check of the sequence of
n
residuals are all below a tolerance.- Parameters
tol – Tolerance.
- Returns
true
if alln
residuals are less than the tolerance.
-
inline const std::vector<double> &data() const
Get the historic residuals.
-
inline const std::vector<double> &get() const
Get the historic residuals.
-
inline StopList(size_t n = 1)
-
class StopList
-
namespace Mesh
Generic mesh operations, and simple mesh definitions.
Enums
-
enum class ElementType
Enumerator for element-types.
Values:
-
enumerator Unknown
Unknown element-type.
-
enumerator Quad4
Quadrilateral: 4-noded element in 2-d.
-
enumerator Hex8
Hexahedron: 8-noded element in 3-d.
-
enumerator Tri3
Triangle: 3-noded element in 2-d.
-
enumerator Unknown
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
-
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.
-
template<class CA, class EA, class NA, class CB, class EB, class NB>
-
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.
-
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.
-
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.
-
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
-
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.
-
template<class T>
-
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
-
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.
-
template<class T>
-
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.
-
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].
-
inline Vstack(bool check_overlap = true, double rtol = 1e-5, double atol = 1e-8)
-
namespace Hex8
Simple meshes of 8-noded hexahedral elements in 3d (ElementType::Hex8).
-
class FineLayer : public GooseFEM::Mesh::RegularBase3d<FineLayer>
- #include <GooseFEM/MeshHex8.h>
Mesh with fine middle layer, and coarser elements towards the top and bottom.
Public Functions
-
inline FineLayer(size_t nelx, size_t nely, size_t nelz, double h = 1.0, size_t nfine = 1)
Constructor.
- Parameters
nelx – Number of elements (along the middle layer) in horizontal (x) direction.
nely – Approximate equivalent number of elements in vertical (y) direction.
nelz – Number of elements (along the middle layer) in depth (z) direction.
h – Edge size (width == height == depth) of elements along the weak layer.
nfine – Extra number of fine layers around the middle layer. By default the element size is kept smaller than the distance to the middle layer.
-
inline array_type::tensor<size_t, 1> elementsMiddleLayer() const
Elements in the middle (fine) layer.
- Returns
List of element numbers (copy, involves computation).
-
inline FineLayer(size_t nelx, size_t nely, size_t nelz, double h = 1.0, size_t nfine = 1)
-
class Regular : public GooseFEM::Mesh::RegularBase3d<Regular>
- #include <GooseFEM/MeshHex8.h>
Regular mesh: equi-sized elements.
Public Functions
-
inline Regular(size_t nelx, size_t nely, size_t nelz, double h = 1.0)
Constructor.
- Parameters
nelx – Number of elements in horizontal (x) direction.
nely – Number of elements in vertical (y) direction.
nelz – Number of elements in vertical (z) direction.
h – Edge size (width == height == depth).
-
inline Regular(size_t nelx, size_t nely, size_t nelz, double h = 1.0)
-
class FineLayer : public GooseFEM::Mesh::RegularBase3d<FineLayer>
-
namespace Quad4
Simple meshes of 4-noded quadrilateral elements in 2d (ElementType::Quad4).
-
class FineLayer : public GooseFEM::Mesh::RegularBase2d<FineLayer>
- #include <GooseFEM/MeshQuad4.h>
Mesh with fine middle layer, and coarser elements towards the top and bottom.
Public Functions
-
inline FineLayer(size_t nelx, size_t nely, double h = 1.0, size_t nfine = 1)
Constructor.
- Parameters
nelx – Number of elements (along the middle layer) in horizontal (x) direction.
nely – Approximate equivalent number of elements in vertical (y) direction.
h – Edge size (width == height) of elements along the weak layer.
nfine – Extra number of fine layers around the middle layer. By default the element size is kept smaller than the distance to the middle layer.
-
template<class C, class E, std::enable_if_t<xt::is_xexpression<C>::value, bool> = true>
inline FineLayer(const C &coor, const E &conn) Reconstruct class for given coordinates / connectivity.
- 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]
withndim == 2
.conn – Connectivity
[nne, nne]
withnne == 4
.
-
inline const array_type::tensor<size_t, 1> &elemrow_nhx() const
Edge size in x-direction of a block, in units of h, per row of blocks.
Note that a block is equal to an element except in refinement layers where it contains three elements.
- Returns
List of size equal to the number of rows of blocks.
-
inline const array_type::tensor<size_t, 1> &elemrow_nhy() const
Edge size in y-direction of a block, in units of h, per row of blocks.
Note that a block is equal to an element except in refinement layers where it contains three elements.
- Returns
List of size equal to the number of rows of blocks.
-
inline const array_type::tensor<int, 1> &elemrow_type() const
Per row of blocks:
-1
: normal layer0
: transition layer to match coarse and finer element on the previous/next row.- Returns
List of size equal to the number of rows of blocks.
-
inline const array_type::tensor<size_t, 1> &elemrow_nelem() const
Number of elements per row of blocks.
Note that a block is equal to an element except in refinement layers where it contains three elements.
- Returns
List of size equal to the number of rows of blocks.
-
inline array_type::tensor<size_t, 1> elementsMiddleLayer() const
Elements in the middle (fine) layer.
- Returns
List of element numbers.
-
inline array_type::tensor<size_t, 1> elementsLayer(size_t layer) const
Elements along a layer.
- Returns
List of element numbers.
-
inline array_type::tensor<size_t, 1> elementgrid_ravel(std::vector<size_t> start_stop_rows, std::vector<size_t> start_stop_cols) const
Select region of elements from ‘matrix’ of element numbers.
- Returns
List of element numbers.
-
inline array_type::tensor<size_t, 1> elementgrid_around_ravel(size_t e, size_t size, bool periodic = true)
Select region of elements from ‘matrix’ of element numbers around an element: square box with edge-size
(2 * size + 1) * h
, aroundelement
.- Parameters
e – The element around which to select elements.
size – Edge size of the square box encapsulating the selected element.
periodic – Assume the mesh periodic.
- Returns
List of elements.
-
inline array_type::tensor<size_t, 1> elementgrid_leftright(size_t e, size_t left, size_t right, bool periodic = true)
Select region of elements from ‘matrix’ of element numbers around an element: left/right from
element
(on the same layer).- Parameters
e – The element around which to select elements.
left – Number of elements to select to the left.
right – Number of elements to select to the right.
periodic – Assume the mesh periodic.
- Returns
List of elements.
-
inline array_type::tensor<size_t, 1> roll(size_t n)
Mapping to ‘roll’ periodically in the x-direction,.
- Returns
element mapping, such that: new_elemvar = elemvar[elem_map]
-
inline FineLayer(size_t nelx, size_t nely, double h = 1.0, size_t nfine = 1)
-
class Regular : public GooseFEM::Mesh::RegularBase2d<Regular>
- #include <GooseFEM/MeshQuad4.h>
Regular mesh: equi-sized elements.
Public Functions
-
inline Regular(size_t nelx, size_t nely, double h = 1.0)
Constructor.
- Parameters
nelx – Number of elements in horizontal (x) direction.
nely – Number of elements in vertical (y) direction.
h – Edge size (width == height).
-
inline array_type::tensor<size_t, 2> elementgrid() const
Element numbers as ‘matrix’.
-
inline Regular(size_t nelx, size_t nely, double h = 1.0)
-
namespace Map
Mesh mappings.
-
class FineLayer2Regular
- #include <GooseFEM/MeshQuad4.h>
Map a FineLayer mesh to a Regular mesh.
The element size of the Regular corresponds to the smallest elements of the FineLayer mesh (along the middle layer).
Public Functions
-
inline FineLayer2Regular(const GooseFEM::Mesh::Quad4::FineLayer &mesh)
Constructors.
- Parameters
mesh – The FineLayer mesh.
-
inline GooseFEM::Mesh::Quad4::FineLayer fineLayerMesh() const
Obtain the FineLayer mesh (copy of the mesh passed to the constructor).
- Returns
mesh.
-
inline std::vector<std::vector<size_t>> map() const
Get element-mapping: elements of the Regular mesh per element of the FineLayer mesh.
The number of Regular elements varies between elements of the FineLayer mesh.
- Returns
[nelem_finelayer, ?]
-
inline std::vector<std::vector<double>> mapFraction() const
To overlap fraction for each item in the mapping in map().
- Returns
[nelem_finelayer, ?]
-
inline GooseFEM::Mesh::Quad4::FineLayer getFineLayerMesh() const
Obtain the FineLayer mesh (copy of the mesh passed to the constructor).
- Returns
mesh.
-
inline std::vector<std::vector<size_t>> getMap() const
Get element-mapping: elements of the Regular mesh per element of the FineLayer mesh.
The number of Regular elements varies between elements of the FineLayer mesh.
- Returns
[nelem_finelayer, ?]
-
inline std::vector<std::vector<double>> getMapFraction() const
To overlap fraction for each item in the mapping in getMap().
- Returns
[nelem_finelayer, ?]
-
template<class T, size_t rank>
inline array_type::tensor<T, rank> mapToRegular(const array_type::tensor<T, rank> &data) const Map element quantities to Regular.
The mapping is a bit simplistic: no interpolation is involved, the function just accounts the fraction of overlap between the FineLayer element and the Regular element. The mapping is such that::
ret[e_regular, ...] <- arg[e_finelayer, ...]
- Template Parameters
T – type of the data (e.g.
double
).rank – rank of the data.
- Parameters
data – data.
- Returns
mapped data.
-
inline FineLayer2Regular(const GooseFEM::Mesh::Quad4::FineLayer &mesh)
-
class RefineRegular
- #include <GooseFEM/MeshQuad4.h>
Refine a Regular mesh: subdivide elements in several smaller elements.
Public Functions
-
inline RefineRegular(const GooseFEM::Mesh::Quad4::Regular &mesh, size_t nx, size_t ny)
Constructor.
- Parameters
mesh – the coarse mesh.
nx – for each coarse element: number of fine elements in x-direction.
ny – for each coarse element: number of fine elements in y-direction.
-
inline size_t nx() const
For each coarse element: number of fine elements in x-direction.
- Returns
unsigned int (same as used in constructor)
-
inline size_t ny() const
For each coarse element: number of fine elements in y-direction.
- Returns
unsigned int (same as used in constructor)
-
inline GooseFEM::Mesh::Quad4::Regular coarseMesh() const
Obtain the coarse mesh (copy of the mesh passed to the constructor).
- Returns
mesh
-
inline const array_type::tensor<size_t, 2> &map() const
Get element-mapping: elements of the fine mesh per element of the coarse mesh.
-
inline GooseFEM::Mesh::Quad4::Regular getCoarseMesh() const
Obtain the coarse mesh (copy of the mesh passed to the constructor).
- Returns
mesh
-
inline const array_type::tensor<size_t, 2> &getMap() const
Get element-mapping: elements of the fine mesh per element of the coarse mesh.
-
template<class T, size_t rank>
inline array_type::tensor<T, rank> meanToCoarse(const array_type::tensor<T, rank> &data) const Compute the mean of the quantity define on the fine mesh when mapped on the coarse mesh.
- Template Parameters
T – type of the data (e.g.
double
).rank – rank of the data.
- Parameters
data – the data [nelem_fine, …]
- Returns
the average data of the coarse mesh [nelem_coarse, …]
-
template<class T, size_t rank, class S>
inline array_type::tensor<T, rank> averageToCoarse(const array_type::tensor<T, rank> &data, const array_type::tensor<S, rank> &weights) const Compute the average of the quantity define on the fine mesh when mapped on the coarse mesh.
- Template Parameters
T – type of the data (e.g.
double
).rank – rank of the data.
S – type of the weights (e.g.
double
).
- Parameters
data – the data [nelem_fine, …]
weights – the weights [nelem_fine, …]
- Returns
the average data of the coarse mesh [nelem_coarse, …]
-
template<class T, size_t rank>
inline array_type::tensor<T, rank> mapToFine(const array_type::tensor<T, rank> &data) const Map element quantities to the fine mesh.
The mapping is a bit simplistic: no interpolation is involved. The mapping is such that::
ret[e_fine, ...] <- data[e_coarse, ...]
- Template Parameters
T – type of the data (e.g.
double
).rank – rank of the data.
- Parameters
data – the data.
- Returns
mapped data.
-
inline RefineRegular(const GooseFEM::Mesh::Quad4::Regular &mesh, size_t nx, size_t ny)
-
class FineLayer2Regular
-
class FineLayer : public GooseFEM::Mesh::RegularBase2d<FineLayer>
-
namespace Tri3
Simple meshes of and mesh operations for triangular elements of type ElementType::Tri3.
Functions
-
inline array_type::tensor<int, 1> getOrientation(const array_type::tensor<double, 2> &coor, const array_type::tensor<size_t, 2> &conn)
Read the orientation of a mesh of triangular elements of type ElementType::Tri3.
- Parameters
coor – Nodal coordinates [nnode, ndim].
conn – Connectivity [nelem, nne].
- Returns
Orientation (-1 or +1) per element [nelem].
-
inline array_type::tensor<size_t, 2> setOrientation(const array_type::tensor<double, 2> &coor, const array_type::tensor<size_t, 2> &conn, const array_type::tensor<int, 1> &val, int orientation = -1)
Set the orientation of a mesh of triangular elements of type ElementType::Tri3.
For efficiency this function reuses the output of getOrientation().
- Parameters
coor – Nodal coordinates [nnode, ndim].
conn – Connectivity [nelem, nne].
val – Current orientation per element (output of getOrientation()) [nelem].
orientation – Target orientation (applied to all elements).
- Returns
Connectivity (order of nodes-per-element may have changed) [nelem, nne].
-
inline array_type::tensor<size_t, 2> setOrientation(const array_type::tensor<double, 2> &coor, const array_type::tensor<size_t, 2> &conn, int orientation = -1)
Set the orientation of a mesh of triangular elements of type ElementType::Tri3.
- Parameters
coor – Nodal coordinates [nnode, ndim].
conn – Connectivity [nelem, nne].
orientation – Target orientation (applied to all elements).
- Returns
Connectivity (order of nodes-per-element may have changed) [nelem, nne].
-
class Regular : public GooseFEM::Mesh::RegularBase2d<Regular>
- #include <GooseFEM/MeshTri3.h>
Regular grid of squares, with each square cut into two triangular elements.
Public Functions
-
inline Regular(size_t nelx, size_t nely, double h = 1.0)
Constructor.
- Parameters
nelx – Number of elements in x-direction.
nely – Number of elements in y-direction.
h – Edge-size (of the squares, and thus of two of three element-edges).
-
inline Regular(size_t nelx, size_t nely, double h = 1.0)
-
inline array_type::tensor<int, 1> getOrientation(const array_type::tensor<double, 2> &coor, const array_type::tensor<size_t, 2> &conn)
-
enum class ElementType
-
namespace Tyings
Tools to store and apply nodal/DOF tyings.
-
class Control
- #include <GooseFEM/TyingsPeriodic.h>
Add control nodes to an existing system.
Public Functions
-
template<class C, class N>
inline Control(const C &coor, const N &dofs) Constructor.
- Template Parameters
- Parameters
coor – Nodal coordinates [nnode, ndim].
dofs – DOF-numbers per node [nnode, ndim].
-
inline const array_type::tensor<double, 2> &coor() const
Nodal coordinates, for the system with control nodes added to it.
- Returns
[nnode + ndim, ndim], with nnode the number of nodes of the original system.
-
inline const array_type::tensor<size_t, 2> &dofs() const
DOF-numbers per node, for the system with control nodes added to it.
- Returns
[nnode + ndim, ndim], with nnode the number of nodes of the original system.
-
inline const array_type::tensor<size_t, 2> &controlDofs() const
DOF-numbers of each control node.
- Returns
[ndim, ndim].
-
inline const array_type::tensor<size_t, 1> &controlNodes() const
Node-numbers of the control nodes.
- Returns
[ndim].
-
template<class C, class N>
-
class Periodic
- #include <GooseFEM/TyingsPeriodic.h>
Nodal tyings per periodic boundary conditions.
The idea is that the displacement of all DOFs of a node are tied to another node and to the average displacement gradient. The latter is applied/measured using ‘virtual’ control nodes.
Consider the DOF list \( u \) renumbered such that it is split up in independent and dependent DOFs as follows
\( u = \begin{bmatrix} u_i \\ u_d \end{bmatrix}\)
whereby the independent DOFs are furthermore split up in unknown and prescribed nodes as follows
\( u_i = \begin{bmatrix} u_u \\ u_p \end{bmatrix}\)
such that
\( u = \begin{bmatrix} u_u \\ u_p \\ u_d \end{bmatrix}\)
- Todo:
Document how the DOFs are tied to the control nodes, and what the has to do with the mean.
Public Functions
-
template<class C, class D, class S, class T>
inline Periodic(const C &coor, const D &dofs, const S &control_dofs, const T &nodal_tyings) Constructor.
- Template Parameters
- Parameters
coor – Nodal coordinates [nnode, ndim].
dofs – DOF-numbers per node [nnode, ndim].
control_dofs – DOF-numbers per control node [ndim, ndim].
nodal_tyings – List of nodal tyings, see nodal_tyings(). [ntyings, 2].
-
template<class C, class D, class S, class T, class U>
inline Periodic(const C &coor, const D &dofs, const S &control_dofs, const T &nodal_tyings, const U &iip) Constructor.
- Template Parameters
- Parameters
coor – Nodal coordinates [nnode, ndim].
dofs – DOF-numbers per node [nnode, ndim].
control_dofs – DOF-numbers per control node [ndim, ndim].
nodal_tyings – List of nodal tyings, see nodal_tyings(). [ntyings, 2].
iip – List of prescribed DOF-numbers.
-
inline size_t nnd() const
- Returns
Number of dependent DOFs.
-
inline size_t nni() const
- Returns
Number of independent DOFs.
-
inline size_t nnu() const
- Returns
Number of independent unknown DOFs.
-
inline size_t nnp() const
- Returns
Number of independent prescribed DOFs.
-
inline const array_type::tensor<size_t, 2> &dofs() const
- Returns
DOF-numbers per node, as used internally (after renumbering), [nnode, ndim].
-
inline const array_type::tensor<size_t, 2> &control() const
- Returns
DOF-numbers for each control node, as used internally (after renumbering), [ndim, ndim].
-
inline const array_type::tensor<size_t, 2> &nodal_tyings() const
Return the applied nodal tyings.
Per tying (row) two node numbers are specified, according to the convention (independent, dependent).
- Returns
[ntyings, 2].
-
inline array_type::tensor<size_t, 1> iid() const
Dependent DOFs.
- Returns
List of DOF numbers.
-
inline array_type::tensor<size_t, 1> iii() const
Independent DOFs.
- Returns
List of DOF numbers.
-
inline array_type::tensor<size_t, 1> iiu() const
Independent unknown DOFs.
- Returns
List of DOF numbers.
-
inline array_type::tensor<size_t, 1> iip() const
Independent prescribed DOFs.
- Returns
List of DOF numbers.
-
inline Eigen::SparseMatrix<double> Cdi() const
Return tying matrix such as to get the dependent DOFs \( u_d \) from the independent DOFs \( u_i \) as follows.
\( u_d = C_{di} u_i \)
Note that this can be further partitioned in
\( u_d = C_{du} u_u + C_{dp} u_p \)
- Returns
Sparse matrix.
-
inline Eigen::SparseMatrix<double> Cdu() const
Unknown part of the partitioned tying matrix, see Cdi().
- Returns
Sparse matrix.
-
inline Eigen::SparseMatrix<double> Cdp() const
Prescribed part of the partitioned tying matrix, see Cdi().
- Returns
Sparse matrix.
-
class Control
-
template<class T, class R>