Class GooseFEM::MatrixPartitionedTyingsSolver#
-
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>>>># 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)#