1 #ifndef DUNE_FEM_PRECONDITIONERWRAPPER_HH 2 #define DUNE_FEM_PRECONDITIONERWRAPPER_HH 5 #include <dune/common/shared_ptr.hh> 12 #include <dune/istl/operators.hh> 13 #include <dune/istl/preconditioners.hh> 15 #include <dune/istl/paamg/amg.hh> 16 #include <dune/istl/paamg/pinfo.hh> 26 template<
class M,
class X,
class Y,
int l=1>
27 class FemSeqILU0 :
public SeqILU0<M,X,Y,l>
29 typedef SeqILU0<M,X,Y,l> BaseType ;
31 typedef typename X::field_type field_type;
32 FemSeqILU0 (
const M& A,
int iter, field_type w)
35 FemSeqILU0 (
const M& A, field_type w)
41 template<
class MatrixObject,
class X,
class Y >
42 class FemDiagonalPreconditioner :
public Preconditioner<X,Y>
45 typedef typename MatrixObject :: ColumnDiscreteFunctionType DiscreteFunctionType ;
47 typedef DiagonalPreconditionerBase< DiscreteFunctionType, MatrixObject, true > PreconditionerType ;
50 PreconditionerType diagonalPrecon_;
53 typedef typename X::field_type field_type;
54 FemDiagonalPreconditioner(
const MatrixObject& mObj )
55 : diagonalPrecon_( mObj )
59 virtual void pre (X& x, Y& b) {}
62 virtual void apply (X& v,
const Y& d)
64 diagonalPrecon_.applyToISTLBlockVector( d, v );
68 virtual void post (X& x) {}
75 template<
class MatrixImp>
76 class PreconditionerWrapper
77 :
public Preconditioner<typename MatrixImp :: RowBlockVectorType,
78 typename MatrixImp :: ColBlockVectorType>
80 typedef MatrixImp MatrixType;
81 typedef typename MatrixType :: RowBlockVectorType X;
82 typedef typename MatrixType :: ColBlockVectorType Y;
85 typedef typename MatrixType :: BaseType ISTLMatrixType ;
87 typedef typename MatrixType :: CollectiveCommunictionType
88 CollectiveCommunictionType;
95 typedef MatrixAdapter< ISTLMatrixType, X, Y > OperatorType;
97 mutable std::shared_ptr< OperatorType > op_;
100 typedef Preconditioner<X,Y> PreconditionerInterfaceType;
101 mutable std::shared_ptr<PreconditionerInterfaceType> preconder_;
106 template <
class XImp,
class YImp>
109 inline static void apply(PreconditionerInterfaceType& preconder,
110 XImp& v,
const YImp& d)
114 inline static void copy(XImp& v,
const YImp& d)
119 template <
class XImp>
120 struct Apply<XImp,XImp>
123 inline static void apply(PreconditionerInterfaceType& preconder,
124 XImp& v,
const XImp& d)
126 preconder.apply( v ,d );
129 inline static void copy(X& v,
const Y& d)
136 typedef X domain_type;
138 typedef Y range_type;
140 typedef typename X::field_type field_type;
144 category=SolverCategory::sequential };
147 PreconditionerWrapper (
const PreconditionerWrapper& org)
149 , preconder_(org.preconder_)
155 PreconditionerWrapper()
162 template <
class PreconditionerType>
163 PreconditionerWrapper(MatrixType & matrix,
166 const PreconditionerType* p )
168 , preconder_( new PreconditionerType( matrix, iter, relax ) )
176 template <
class PreconditionerType>
177 PreconditionerWrapper(MatrixType & matrix,
178 PreconditionerType* p )
186 template <
class PreconditionerType>
187 PreconditionerWrapper(MatrixType & matrix,
190 const PreconditionerType* p ,
191 const CollectiveCommunictionType& comm )
195 : op_( new OperatorType( matrix ) )
197 , preconder_( createAMGPreconditioner(comm, iter, relax, p) )
203 virtual void pre (X& x, Y& b)
210 preconder_->pre(x,b);
216 virtual void apply (X& v,
const Y& d)
221 Apply<X,Y>::apply( *preconder_ , v, d );
225 Apply<X,Y>::copy( v, d );
230 virtual void post (X& x)
241 template <
class Smoother>
242 PreconditionerInterfaceType*
243 createAMGPreconditioner(
const CollectiveCommunictionType& comm,
244 int iter, field_type relax,
const Smoother* )
246 typedef typename Dune::FieldTraits< field_type>::real_type real_type;
247 typedef typename std::conditional< std::is_convertible<field_type,real_type>::value,
248 Dune::Amg::FirstDiagonal, Dune::Amg::RowSum >::type Norm;
249 typedef Dune::Amg::CoarsenCriterion<
250 Dune::Amg::UnSymmetricCriterion<ISTLMatrixType, Norm > > Criterion;
252 typedef typename Dune::Amg::SmootherTraits<Smoother>::Arguments SmootherArgs;
254 SmootherArgs smootherArgs;
256 smootherArgs.iterations = iter;
257 smootherArgs.relaxationFactor = relax ;
259 int coarsenTarget=1200;
260 Criterion criterion(15,coarsenTarget);
261 criterion.setDefaultValuesIsotropic(2);
262 criterion.setAlpha(.67);
263 criterion.setBeta(1.0e-8);
264 criterion.setMaxLevel(10);
278 typedef Dune::Amg::AMG<OperatorType, X, Smoother> AMG;
279 return new AMG(*op_, criterion, smootherArgs);
284 #endif // end HAVE_DUNE_ISTL 290 #endif // #ifndef DUNE_FEM_PRECONDITIONERWRAPPER_HH
Definition: coordinate.hh:4