All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
mortar_schur_precond.hpp
1 //==============================================================================
11 //==============================================================================
12 #ifndef SCHUR_PRECOND_HPP_
13 #define SCHUR_PRECOND_HPP_
14 
15 #include <dune/istl/preconditioners.hh>
16 #include <dune/istl/matrixmatrix.hh>
20 
21 namespace Opm {
22 namespace Elasticity {
23 
36  template<class PrecondElasticityBlock>
37 class MortarSchurPre : public Dune::Preconditioner<Vector,Vector> {
38  public:
39  // define the category
40  enum {
42  category=Dune::SolverCategory::sequential
43  };
44 
50  MortarSchurPre(const Matrix& P_, const Matrix& B_,
51  PrecondElasticityBlock& Apre_, bool symmetric_=false) :
52  Apre(Apre_), B(B_), N(B.N()), M(B.M()),
53  Lpre(P_, false, false), symmetric(symmetric_)
54  {
55  }
56 
58  virtual ~MortarSchurPre()
59  {
60  }
61 
63  virtual void pre(Vector& x, Vector& b)
64  {
65  // extract displacement DOFs
66  Vector tempx, tempb;
67  MortarUtils::extractBlock(tempx, x, N);
68  MortarUtils::extractBlock(tempb, b, N);
69  Apre.pre(tempx, tempb);
70  MortarUtils::injectBlock(x, tempx, N);
71  MortarUtils::injectBlock(b, tempb, N);
72  }
73 
77  virtual void apply(Vector& v, const Vector& d)
78  {
79  // multiplier block (second row)
80  Vector temp;
81  MortarUtils::extractBlock(temp, d, M, N);
82  Dune::InverseOperatorResult r;
83  Vector temp2;
84  temp2.resize(temp.size(), false);
85  Lpre.apply(temp2, temp, r);
86  MortarUtils::injectBlock(v, temp2, M, N);
87 
88  // elasticity block (first row)
89  MortarUtils::extractBlock(temp, d, N);
90  if (!symmetric)
91  B.mmv(temp2, temp);
92 
93  temp2.resize(N, false);
94  temp2 = 0;
95  Apre.apply(temp2, temp);
96  MortarUtils::injectBlock(v, temp2, N);
97  }
98 
100  virtual void post(Vector& x)
101  {
102  Apre.post(x);
103  }
104  protected:
106  PrecondElasticityBlock& Apre;
107 
109  const Matrix& B;
110 
112  int N;
113 
115  int M;
116 
118  LUSolver Lpre;
119 
121  bool symmetric;
122 };
123 
124 }
125 }
126 
127 #endif
Helper class with some matrix operations.
Preconditioners for elasticity upscaling.
Mortar helper class.
bool symmetric
Whether or not to use a symmetric preconditioner.
Definition: mortar_schur_precond.hpp:121
int N
Number of displacement DOFs.
Definition: mortar_schur_precond.hpp:112
virtual ~MortarSchurPre()
Destructor.
Definition: mortar_schur_precond.hpp:58
The category the preconditioner is part of.
Definition: mortar_schur_precond.hpp:42
virtual void pre(Vector &x, Vector &b)
Preprocess preconditioner.
Definition: mortar_schur_precond.hpp:63
virtual void post(Vector &x)
Dummy post-process function.
Definition: mortar_schur_precond.hpp:100
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
This implements a Schur-decomposition based preconditioner for the mortar-elasticity system [A B] [B&#39;...
Definition: mortar_schur_precond.hpp:37
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
MortarSchurPre(const Matrix &P_, const Matrix &B_, PrecondElasticityBlock &Apre_, bool symmetric_=false)
Constructor.
Definition: mortar_schur_precond.hpp:50
PrecondElasticityBlock & Apre
The preconditioner for the elasticity operator.
Definition: mortar_schur_precond.hpp:106
const Matrix & B
The mortar coupling matrix.
Definition: mortar_schur_precond.hpp:109
LUSolver Lpre
Linear solver for the multiplier block.
Definition: mortar_schur_precond.hpp:118
int M
Number of multiplier DOFs.
Definition: mortar_schur_precond.hpp:115
virtual void apply(Vector &v, const Vector &d)
Applies the preconditioner.
Definition: mortar_schur_precond.hpp:77