12#ifndef MORTAR_EVALUATOR_HPP_
13#define MORTAR_EVALUATOR_HPP_
15#include <dune/istl/matrixmatrix.hh>
16#include <dune/istl/solvers.hh>
48 l2.resize(lambda.size());
64 l2.resize(lambda.size());
66 for (
size_t i=
A.N(); i < y.size(); ++i)
67 y[i] += alpha*l2[i-
A.N()];
70 Dune::SolverCategory::Category category()
const override
72 return Dune::SolverCategory::sequential;
This implements a operator evaluation for the schur mortar-block S = B^T*A^-1*B !
Definition mortar_evaluator.hpp:27
const Matrix & B
Reference to the mortar coupling matrix.
Definition mortar_evaluator.hpp:77
void apply(const Vector &x, Vector &y) const override
Apply the multiplier block.
Definition mortar_evaluator.hpp:41
void applyscaleadd(field_type alpha, const Vector &x, Vector &y) const override
Apply the multiplier block with an embedded axpy.
Definition mortar_evaluator.hpp:57
MortarEvaluator(const Matrix &A_, const Matrix &B_)
Constructor.
Definition mortar_evaluator.hpp:32
const Matrix & A
Reference to A matrix.
Definition mortar_evaluator.hpp:76
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
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