File MatrixDiagonal.h#
Diagonal matrix.
- Copyright
Copyright 2017. Tom de Geus. All rights reserved.
- License: This project is released under the GNU Public License (GPLv3).
-
namespace GooseFEM
Toolbox to perform finite element computations.
-
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
-
MatrixDiagonal() = default#
-
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].
Private Functions
-
inline void factorize()#
Inverse of the matrix.
Compute inverse (automatically evaluated by “solve”).
Private Members
- friend MatrixBase< MatrixDiagonal >
- friend MatrixDiagonalBase< MatrixDiagonal >
-
array_type::tensor<double, 1> m_A#
The matrix.
-
array_type::tensor<double, 1> m_inv#
-
MatrixDiagonal() = default#
-
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].
Private Functions
-
inline auto derived_cast() -> derived_type&#
-
inline auto derived_cast() const -> const derived_type&#
-
using derived_type = D
-
class MatrixDiagonal : public GooseFEM::MatrixBase<MatrixDiagonal>, public GooseFEM::MatrixDiagonalBase<MatrixDiagonal>