File MatrixPartitioned.h#
Sparse matrix that is partitioned in:
unknown DOFs
prescribed DOFs
- 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 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
-
MatrixPartitioned() = default#
-
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].
Protected Attributes
-
Eigen::SparseMatrix<double> m_Auu#
The matrix.
-
Eigen::SparseMatrix<double> m_Aup#
The matrix.
-
Eigen::SparseMatrix<double> m_Apu#
The matrix.
-
Eigen::SparseMatrix<double> m_App#
The matrix.
-
array_type::tensor<size_t, 2> m_part#
Renumbered DOFs per node, such that:
making is much simpler to slice.iiu = arange(nnu) iip = nnu + arange(nnp)
-
array_type::tensor<size_t, 1> m_part1d#
Map real DOF to DOF in partitioned system.
The partitioned system is defined as:
Similar toiiu = arange(nnu) iip = nnu + arange(nnp)
m_part
but for a 1d sequential list of DOFs.
Private Functions
-
inline void reaction_p_impl(const array_type::tensor<double, 1> &x_u, const array_type::tensor<double, 1> &x_p, array_type::tensor<double, 1> &b_p) const#
-
inline Eigen::VectorXd AsDofs_u(const array_type::tensor<double, 1> &dofval) const#
-
inline Eigen::VectorXd AsDofs_u(const array_type::tensor<double, 2> &nodevec) const#
-
inline Eigen::VectorXd AsDofs_p(const array_type::tensor<double, 1> &dofval) const#
-
inline Eigen::VectorXd AsDofs_p(const array_type::tensor<double, 2> &nodevec) const#
Private Members
- friend MatrixBase< MatrixPartitioned >
- friend MatrixPartitionedBase< MatrixPartitioned >
Friends
- friend class MatrixPartitionedSolver
-
MatrixPartitioned() = default#
-
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
-
MatrixPartitionedSolver() = default#
-
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].
Private Functions
-
template<class T>
inline void solve_nodevec_impl(MatrixPartitioned &A, const T &b, T &x)#
-
template<class T>
inline void solve_dofval_impl(MatrixPartitioned &A, const T &b, T &x)#
-
template<class T>
inline void solve_u_impl(MatrixPartitioned &A, const T &b_u, const T &x_p, T &x_u)#
-
inline void factorize(MatrixPartitioned &A)#
compute inverse (evaluated by “solve”)
-
MatrixPartitionedSolver() = default#
-
class MatrixPartitioned : public GooseFEM::MatrixPartitionedBase<MatrixPartitioned>