12#ifndef SCHUR_PRECOND_HPP_
13#define SCHUR_PRECOND_HPP_
15#include <dune/istl/preconditioners.hh>
16#include <dune/istl/matrixmatrix.hh>
36 template<
class PrecondElasticityBlock>
45 PrecondElasticityBlock& Apre_,
bool symmetric_=
false) :
63 Apre.pre(tempx, tempb);
76 Dune::InverseOperatorResult r;
78 temp2.resize(temp.size());
79 Lpre.apply(temp2, temp, r);
89 Apre.apply(temp2, temp);
99 Dune::SolverCategory::Category category()
const override
101 return Dune::SolverCategory::sequential;
This implements a Schur-decomposition based preconditioner for the mortar-elasticity system [A B] [B'...
Definition mortar_schur_precond.hpp:37
void post(Vector &x) override
Dummy post-process function.
Definition mortar_schur_precond.hpp:94
virtual ~MortarSchurPre()
Destructor.
Definition mortar_schur_precond.hpp:52
PrecondElasticityBlock & Apre
The preconditioner for the elasticity operator.
Definition mortar_schur_precond.hpp:106
int M
Number of multiplier DOFs.
Definition mortar_schur_precond.hpp:115
void pre(Vector &x, Vector &b) override
Preprocess preconditioner.
Definition mortar_schur_precond.hpp:57
const Matrix & B
The mortar coupling matrix.
Definition mortar_schur_precond.hpp:109
bool symmetric
Whether or not to use a symmetric preconditioner.
Definition mortar_schur_precond.hpp:121
LUSolver Lpre
Linear solver for the multiplier block.
Definition mortar_schur_precond.hpp:118
MortarSchurPre(const Matrix &P_, const Matrix &B_, PrecondElasticityBlock &Apre_, bool symmetric_=false)
Constructor.
Definition mortar_schur_precond.hpp:44
int N
Number of displacement DOFs.
Definition mortar_schur_precond.hpp:112
void apply(Vector &v, const Vector &d) override
Applies the preconditioner.
Definition mortar_schur_precond.hpp:71
static void extractBlock(Vector &x, const Vector &y, int len, int start=0)
Extract a range of indices from a vector.
Definition mortar_utils.hpp:25
static void injectBlock(Vector &x, const Vector &y, int len, int start=0)
Inject a range of indices into a vector.
Definition mortar_utils.hpp:36
Preconditioners for elasticity upscaling.
Helper class with some matrix operations.
Dune::BCRSMatrix< Dune::FieldMatrix< double, 1, 1 > > Matrix
A sparse matrix holding our operator.
Definition matrixops.hpp:27
Dune::BlockVector< Dune::FieldVector< double, 1 > > Vector
A vector holding our RHS.
Definition matrixops.hpp:33
Inverting small matrices.
Definition ImplicitAssembly.hpp:43