1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*- 2 // vi: set et ts=4 sw=2 sts=2: 3 #ifndef DUNE_ISTL_PRECONDITIONER_HH 4 #define DUNE_ISTL_PRECONDITIONER_HH 5 6 #include <dune/common/exceptions.hh> 7 8 #include "solvercategory.hh" 9 10 namespace Dune { 11 /** 12 * @addtogroup ISTL_Prec 13 * @{ 14 */ 15 //===================================================================== 16 /*! \brief Base class for matrix free definition of preconditioners. 17 18 Note that the operator, which is the basis for the preconditioning, 19 is supplied to the preconditioner from the outside in the 20 constructor or some other method. 21 22 This interface allows the encapsulation of all parallelization 23 aspects into the preconditioners. 24 25 \tparam X Type of the update 26 \tparam Y Type of the defect 27 */ 28 //===================================================================== 29 template<class X, class Y> 30 class Preconditioner { 31 public: 32 //! \brief The domain type of the preconditioner. 33 typedef X domain_type; 34 //! \brief The range type of the preconditioner. 35 typedef Y range_type; 36 //! \brief The field type of the preconditioner. 37 typedef typename X::field_type field_type; 38 39 /*! \brief Prepare the preconditioner. 40 41 A solver solves a linear operator equation A(x)=b by applying 42 one or several steps of the preconditioner. The method pre() 43 is called before the first apply operation. 44 b and x are right hand side and solution vector of the linear 45 system respectively. It may. e.g., scale the system, allocate memory or 46 compute a (I)LU decomposition. 47 Note: The ILU decomposition could also be computed in the constructor 48 or with a separate method of the derived method if several 49 linear systems with the same matrix are to be solved. 50 51 \note if a preconditioner is copied (e.g. for a second thread) 52 again the pre() method has to be called to ensure proper memory 53 mangement. 54 55 \code 56 X x(0.0); 57 Y b = ...; // rhs 58 Preconditioner<X,Y> prec(...); 59 prec.pre(x,b); // prepare the preconditioner 60 prec.apply(x,b); // can be called multiple times now... 61 prec.post(x); // cleanup internal state 62 \endcode 63 64 \param x The left hand side of the equation. 65 \param b The right hand side of the equation. 66 */ 67 virtual void pre (X& x, Y& b) = 0; 68 69 /*! \brief Apply one step of the preconditioner to the system A(v)=d. 70 71 On entry v=0 and d=b-A(x) (although this might not be 72 computed in that way. On exit v contains the update, i.e 73 one step computes \f$ v = M^{-1} d \f$ where \f$ M \f$ is the 74 approximate inverse of the operator \f$ A \f$ characterizing 75 the preconditioner. 76 \param[out] v The update to be computed 77 \param d The current defect. 78 */ 79 virtual void apply (X& v, const Y& d) = 0; 80 81 /*! \brief Clean up. 82 83 This method is called after the last apply call for the 84 linear system to be solved. Memory may be deallocated safely 85 here. x is the solution of the linear equation. 86 87 \param x The right hand side of the equation. 88 */ 89 virtual void post (X& x) = 0; 90 91 //! Category of the preconditioner (see SolverCategory::Category) category() const92 virtual SolverCategory::Category category() const 93 #if DUNE_ISTL_SUPPORT_OLD_CATEGORY_INTERFACE 94 { 95 DUNE_THROW(Dune::Exception,"It is necessary to implement the category method in a derived classes, in the future this method will pure virtual."); 96 } 97 #else 98 = 0; 99 #endif 100 101 //! every abstract base class has a virtual destructor ~Preconditioner()102 virtual ~Preconditioner () {} 103 104 }; 105 106 /** 107 * @} 108 */ 109 } 110 #endif 111