/usr/include/dune/istl/operators.hh is in libdune-istl-dev 2.2.1-2.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 | // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_ISTLOPERATORS_HH
#define DUNE_ISTLOPERATORS_HH
#include<cmath>
#include<complex>
#include<iostream>
#include<iomanip>
#include<string>
#include"solvercategory.hh"
namespace Dune {
/**
* @defgroup ISTL_Operators Operator concept
* @ingroup ISTL_Solvers
*
* The solvers in ISTL do not work on matrices directly. Instead we use
* an abstract operator concept. This allows for using matrix-free
* operators, i.e. operators that are not stored as matrices in any
* form. Thus our solver algorithms can easily be turned into matrix-free
* solvers just by plugging in matrix-free representations of linear
* operators and preconditioners.
*/
/** @addtogroup ISTL_Operators
@{
*/
/** \file
\brief Define general, extensible interface for operators.
The available implementation wraps a matrix.
*/
//=====================================================================
// Abstract operator interface
//=====================================================================
/*!
@brief A linear operator.
Abstract base class defining a linear operator \f$ A : X\to Y\f$,
i.e. \f$ A(\alpha x) = \alpha A(x) \f$ and
\f$ A(x+y) = A(x)+A(y)\f$ hold. The
simplest solvers just need the application \f$ A(x)\f$ of
the operator.
- enables on the fly computation through operator concept. If explicit
representation of the operator is required use AssembledLinearOperator
- Some inverters may need an explicit formation of the operator
as a matrix, e.g. BiCGStab, ILU, AMG, etc. In that case use the
derived class
*/
template<class X, class Y>
class LinearOperator {
public:
//! The type of the domain of the operator.
typedef X domain_type;
//! The type of the range of the operator.
typedef Y range_type;
//! The field type of the operator.
typedef typename X::field_type field_type;
/*! \brief apply operator to x: \f$ y = A(x) \f$
The input vector is consistent and the output must also be
consistent on the interior+border partition.
*/
virtual void apply (const X& x, Y& y) const = 0;
//! apply operator to x, scale and add: \f$ y = y + \alpha A(x) \f$
virtual void applyscaleadd (field_type alpha, const X& x, Y& y) const = 0;
//! every abstract base class has a virtual destructor
virtual ~LinearOperator () {}
};
/*!
\brief A linear operator exporting itself in matrix form.
Linear Operator that exports the operator in
matrix form. This is needed for certain solvers, such as
LU decomposition, ILU preconditioners or BiCG-Stab (because
of multiplication with A^T).
*/
template<class M, class X, class Y>
class AssembledLinearOperator : public LinearOperator<X,Y> {
public:
//! export types, usually they come from the derived class
typedef M matrix_type;
typedef X domain_type;
typedef Y range_type;
typedef typename X::field_type field_type;
//! get matrix via *
virtual const M& getmat () const = 0;
//! every abstract base class has a virtual destructor
virtual ~AssembledLinearOperator () {}
};
//=====================================================================
// Implementation for ISTL-matrix based operator
//=====================================================================
/*!
\brief Adapter to turn a matrix into a linear operator.
Adapts a matrix to the assembled linear operator interface
*/
template<class M, class X, class Y>
class MatrixAdapter : public AssembledLinearOperator<M,X,Y>
{
public:
//! export types
typedef M matrix_type;
typedef X domain_type;
typedef Y range_type;
typedef typename X::field_type field_type;
//! define the category
enum {category=SolverCategory::sequential};
//! constructor: just store a reference to a matrix
MatrixAdapter (const M& A) : _A_(A) {}
//! apply operator to x: \f$ y = A(x) \f$
virtual void apply (const X& x, Y& y) const
{
_A_.mv(x,y);
}
//! apply operator to x, scale and add: \f$ y = y + \alpha A(x) \f$
virtual void applyscaleadd (field_type alpha, const X& x, Y& y) const
{
_A_.usmv(alpha,x,y);
}
//! get matrix via *
virtual const M& getmat () const
{
return _A_;
}
private:
const M& _A_;
};
/** @} end documentation */
} // end namespace
#endif
|