dune-istl 2.9.1
Loading...
Searching...
No Matches
twolevelmethod.hh
Go to the documentation of this file.
1// SPDX-FileCopyrightText: Copyright (C) DUNE Project contributors, see file LICENSE.md in module root
2// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
3// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
4// vi: set et ts=4 sw=2 sts=2:
5#ifndef DUNE_ISTL_TWOLEVELMETHOD_HH
6#define DUNE_ISTL_TWOLEVELMETHOD_HH
7
8#include <tuple>
9
11#include"amg.hh"
12#include"galerkin.hh"
13#include<dune/istl/solver.hh>
14
22namespace Dune
23{
24namespace Amg
25{
26
36template<class FO, class CO>
38{
39public:
44 typedef FO FineOperatorType;
48 typedef typename FineOperatorType::range_type FineRangeType;
52 typedef typename FineOperatorType::domain_type FineDomainType;
61 typedef typename CoarseOperatorType::range_type CoarseRangeType;
65 typedef typename CoarseOperatorType::domain_type CoarseDomainType;
70 std::shared_ptr<CoarseOperatorType>& getCoarseLevelOperator()
71 {
72 return operator_;
73 }
79 {
80 return rhs_;
81 }
82
88 {
89 return lhs_;
90 }
100 virtual void moveToCoarseLevel(const FineRangeType& fineRhs)=0;
110 virtual void moveToFineLevel(FineDomainType& fineLhs)=0;
118 virtual void createCoarseLevelSystem(const FineOperatorType& fineOperator)=0;
119
121 virtual LevelTransferPolicy* clone() const =0;
122
125
126 protected:
132 std::shared_ptr<CoarseOperatorType> operator_;
133};
134
140template<class O, class C>
142 : public LevelTransferPolicy<O,O>
143{
145public:
147 typedef C Criterion;
149
151 : criterion_(crit)
152 {}
153
154 void createCoarseLevelSystem(const O& fineOperator)
155 {
156 prolongDamp_ = criterion_.getProlongationDampingFactor();
160 Dune::Amg::EdgeProperties,Dune::IdentityMap,Dune::IdentityMap> PropertiesGraph;
161 MatrixGraph mg(fineOperator.getmat());
162 PropertiesGraph pg(mg,Dune::IdentityMap(),Dune::IdentityMap());
163 typedef NegateSet<typename ParallelInformation::OwnerSet> OverlapFlags;
164
165 aggregatesMap_ = std::make_shared<AggregatesMap>(pg.maxVertex()+1);
166
167 int noAggregates, isoAggregates, oneAggregates, skippedAggregates;
168
169 std::tie(noAggregates, isoAggregates, oneAggregates, skippedAggregates) =
170 aggregatesMap_->buildAggregates(fineOperator.getmat(), pg, criterion_, true);
171 std::cout<<"no aggregates="<<noAggregates<<" iso="<<isoAggregates<<" one="<<oneAggregates<<" skipped="<<skippedAggregates<<std::endl;
172 // misuse coarsener to renumber aggregates
174 typedef std::vector<bool>::iterator Iterator;
175 typedef Dune::IteratorPropertyMap<Iterator, Dune::IdentityMap> VisitedMap;
176 std::vector<bool> excluded(fineOperator.getmat().N(), false);
177 VisitedMap vm(excluded.begin(), Dune::IdentityMap());
179 std::size_t aggregates = renumberer.coarsen(pinfo, pg, vm,
180 *aggregatesMap_, pinfo,
181 noAggregates);
182 std::vector<bool>& visited=excluded;
183
184 typedef std::vector<bool>::iterator Iterator;
185
186 for(Iterator iter= visited.begin(), end=visited.end();
187 iter != end; ++iter)
188 *iter=false;
189 matrix_.reset(productBuilder.build(mg, vm,
191 *aggregatesMap_,
192 aggregates,
193 OverlapFlags()));
194 productBuilder.calculate(fineOperator.getmat(), *aggregatesMap_, *matrix_, pinfo, OverlapFlags());
195 this->lhs_.resize(this->matrix_->M());
196 this->rhs_.resize(this->matrix_->N());
197 this->operator_ = std::make_shared<O>(*matrix_);
198 }
199
206
213
215 {
216 return new AggregationLevelTransferPolicy(*this);
217 }
218
219private:
220 typename O::matrix_type::field_type prolongDamp_;
221 std::shared_ptr<AggregatesMap> aggregatesMap_;
222 Criterion criterion_;
223 std::shared_ptr<typename O::matrix_type> matrix_;
224};
225
232template<class O, class S, class C>
234{
235public:
237 typedef O Operator;
239 typedef typename O::range_type X;
241 typedef C Criterion;
243 typedef S Smoother;
254 : smootherArgs_(args), criterion_(c)
255 {}
258 : coarseOperator_(other.coarseOperator_), smootherArgs_(other.smootherArgs_),
259 criterion_(other.criterion_)
260 {}
261private:
268 struct AMGInverseOperator : public InverseOperator<X,X>
269 {
270 AMGInverseOperator(const typename AMGType::Operator& op,
271 const Criterion& crit,
272 const typename AMGType::SmootherArgs& args)
273 : amg_(op, crit,args), first_(true)
274 {}
275
276 void apply(X& x, X& b, [[maybe_unused]] double reduction, [[maybe_unused]] InverseOperatorResult& res)
277 {
278 if(first_)
279 {
280 amg_.pre(x,b);
281 first_=false;
282 x_=x;
283 }
284 amg_.apply(x,b);
285 }
286
287 void apply(X& x, X& b, InverseOperatorResult& res)
288 {
289 return apply(x,b,1e-8,res);
290 }
291
293 virtual SolverCategory::Category category() const
294 {
295 return amg_.category();
296 }
297
298 ~AMGInverseOperator()
299 {
300 if(!first_)
301 amg_.post(x_);
302 }
303 AMGInverseOperator(const AMGInverseOperator& other)
304 : x_(other.x_), amg_(other.amg_), first_(other.first_)
305 {
306 }
307 private:
308 X x_;
309 AMGType amg_;
310 bool first_;
311 };
312
313public:
315 typedef AMGInverseOperator CoarseLevelSolver;
316
324 template<class P>
326 {
327 coarseOperator_=transferPolicy.getCoarseLevelOperator();
328 AMGInverseOperator* inv = new AMGInverseOperator(*coarseOperator_,
329 criterion_,
330 smootherArgs_);
331
332 return inv; //std::shared_ptr<InverseOperator<X,X> >(inv);
333
334 }
335
336private:
338 std::shared_ptr<Operator> coarseOperator_;
340 SmootherArgs smootherArgs_;
342 Criterion criterion_;
343};
344
350template<class FO, class CSP, class S>
352 public Preconditioner<typename FO::domain_type, typename FO::range_type>
353{
354public:
358 typedef typename CoarseLevelSolverPolicy::CoarseLevelSolver CoarseLevelSolver;
367 typedef typename FineOperatorType::range_type FineRangeType;
371 typedef typename FineOperatorType::domain_type FineDomainType;
376 typedef typename CSP::Operator CoarseOperatorType;
380 typedef typename CoarseOperatorType::range_type CoarseRangeType;
384 typedef typename CoarseOperatorType::domain_type CoarseDomainType;
388 typedef S SmootherType;
389
405 std::shared_ptr<SmootherType> smoother,
407 CoarseOperatorType>& policy,
408 CoarseLevelSolverPolicy& coarsePolicy,
409 std::size_t preSteps=1, std::size_t postSteps=1)
410 : operator_(&op), smoother_(smoother),
411 preSteps_(preSteps), postSteps_(postSteps)
412 {
413 policy_ = policy.clone();
414 policy_->createCoarseLevelSystem(*operator_);
415 coarseSolver_=coarsePolicy.createCoarseLevelSolver(*policy_);
416 }
417
419 : operator_(other.operator_), coarseSolver_(new CoarseLevelSolver(*other.coarseSolver_)),
420 smoother_(other.smoother_), policy_(other.policy_->clone()),
421 preSteps_(other.preSteps_), postSteps_(other.postSteps_)
422 {}
423
425 {
426 // Each instance has its own policy.
427 delete policy_;
428 delete coarseSolver_;
429 }
430
432 {
433 smoother_->pre(x,b);
434 }
435
436 void post([[maybe_unused]] FineDomainType& x)
437 {}
438
440 {
441 FineDomainType u(v);
442 FineRangeType rhs(d);
443 LevelContext context;
445 context.pinfo=&info;
446 context.lhs=&u;
447 context.update=&v;
448 context.smoother=smoother_;
449 context.rhs=&rhs;
450 context.matrix=operator_;
451 // Presmoothing
452 presmooth(context, preSteps_);
453 //Coarse grid correction
454 policy_->moveToCoarseLevel(*context.rhs);
456 coarseSolver_->apply(policy_->getCoarseLevelLhs(), policy_->getCoarseLevelRhs(), res);
457 *context.lhs=0;
458 policy_->moveToFineLevel(*context.lhs);
459 *context.update += *context.lhs;
460 // Postsmoothing
461 postsmooth(context, postSteps_);
462
463 }
464
467 {
469 }
470
471private:
475 struct LevelContext
476 {
478 typedef S SmootherType;
480 std::shared_ptr<SmootherType> smoother;
482 FineDomainType* lhs;
483 /*
484 * @brief The right hand side holding the current residual.
485 *
486 * This is passed to the smoother as the right hand side.
487 */
488 FineRangeType* rhs;
494 FineDomainType* update;
502 const FineOperatorType* matrix;
503 };
504 const FineOperatorType* operator_;
506 CoarseLevelSolver* coarseSolver_;
508 std::shared_ptr<S> smoother_;
510 LevelTransferPolicy<FO,typename CSP::Operator>* policy_;
512 std::size_t preSteps_;
514 std::size_t postSteps_;
515};
516}// end namespace Amg
517}// end namespace Dune
518
520#endif
Define general, extensible interface for inverse operators.
Provides a class for building the galerkin product based on a aggregation scheme.
The AMG preconditioner.
Define general, extensible interface for operators. The available implementation wraps a matrix.
void pre(Domain &x, Range &b)
Prepare the preconditioner.
Definition amg.hh:801
G::MutableMatrix * build(G &fineGraph, V &visitedMap, const ParallelInformation &pinfo, AggregatesMap< typename G::VertexDescriptor > &aggregates, const typename G::Matrix::size_type &size, const Set &copy)
Calculates the coarse matrix via a Galerkin product.
Definition galerkin.hh:563
SmootherTraits< Smoother >::Arguments SmootherArgs
The argument type for the construction of the smoother.
Definition amg.hh:100
Operator Operator
The matrix operator type.
Definition amg.hh:73
void presmooth(LevelContext &levelContext, size_t steps)
Apply pre smoothing on the current level.
Definition smoother.hh:406
void post(Domain &x)
Clean up.
Definition amg.hh:1139
void postsmooth(LevelContext &levelContext, size_t steps)
Apply post smoothing on the current level.
Definition smoother.hh:428
void apply(Domain &v, const Range &d)
Apply one step of the preconditioner to the system A(v)=d.
Definition amg.hh:888
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition amg.hh:194
void calculate(const M &fine, const AggregatesMap< V > &aggregates, M &coarse, const I &pinfo, const O &copy)
Calculate the galerkin product.
Definition allocator.hh:11
Class providing information about the mapping of the vertices onto aggregates.
Definition aggregates.hh:560
Parallel algebraic multigrid based on agglomeration.
Definition amg.hh:65
Class representing the properties of an ede in the matrix graph.
Definition dependency.hh:39
Class representing a node in the matrix graph.
Definition dependency.hh:126
Definition galerkin.hh:118
The (undirected) graph of a matrix.
Definition graph.hh:51
Attaches properties to the edges and vertices of a graph.
Definition graph.hh:978
VertexDescriptor maxVertex() const
Get the maximal vertex descriptor.
Definition indicescoarsener.hh:36
Definition pinfo.hh:28
The default class for the smoother arguments.
Definition smoother.hh:38
Definition transfer.hh:32
Abstract base class for transfer between levels and creation of the coarse level system.
Definition twolevelmethod.hh:38
CO CoarseOperatorType
The linear operator of the finel level system. Has to be derived from AssembledLinearOperator.
Definition twolevelmethod.hh:57
virtual void moveToCoarseLevel(const FineRangeType &fineRhs)=0
Transfers the data to the coarse level.
FineOperatorType::range_type FineRangeType
The type of the range of the fine level operator.
Definition twolevelmethod.hh:48
virtual void createCoarseLevelSystem(const FineOperatorType &fineOperator)=0
Algebraically creates the coarse level system.
CoarseOperatorType::range_type CoarseRangeType
The type of the range of the coarse level operator.
Definition twolevelmethod.hh:61
virtual ~LevelTransferPolicy()
Destructor.
Definition twolevelmethod.hh:124
CoarseDomainType lhs_
The coarse level lhs.
Definition twolevelmethod.hh:130
virtual LevelTransferPolicy * clone() const =0
Clone the current object.
CoarseDomainType & getCoarseLevelLhs()
Get the coarse level left hand side.
Definition twolevelmethod.hh:87
std::shared_ptr< CoarseOperatorType > operator_
the coarse level linear operator.
Definition twolevelmethod.hh:132
CoarseRangeType rhs_
The coarse level rhs.
Definition twolevelmethod.hh:128
virtual void moveToFineLevel(FineDomainType &fineLhs)=0
Updates the fine level linear system after the correction of the coarse levels system.
std::shared_ptr< CoarseOperatorType > & getCoarseLevelOperator()
Get the coarse level operator.
Definition twolevelmethod.hh:70
CoarseRangeType & getCoarseLevelRhs()
Get the coarse level right hand side.
Definition twolevelmethod.hh:78
FO FineOperatorType
The linear operator of the finel level system. Has to be derived from AssembledLinearOperator.
Definition twolevelmethod.hh:44
CoarseOperatorType::domain_type CoarseDomainType
The type of the domain of the coarse level operator.
Definition twolevelmethod.hh:65
FineOperatorType::domain_type FineDomainType
The type of the domain of the fine level operator.
Definition twolevelmethod.hh:52
A LeveTransferPolicy that used aggregation to construct the coarse level system.
Definition twolevelmethod.hh:143
C Criterion
Definition twolevelmethod.hh:147
AggregationLevelTransferPolicy(const Criterion &crit)
Definition twolevelmethod.hh:150
AggregationLevelTransferPolicy * clone() const
Clone the current object.
Definition twolevelmethod.hh:214
void moveToFineLevel(typename FatherType::FineDomainType &fineLhs)
Updates the fine level linear system after the correction of the coarse levels system.
Definition twolevelmethod.hh:207
void moveToCoarseLevel(const typename FatherType::FineRangeType &fineRhs)
Definition twolevelmethod.hh:200
SequentialInformation ParallelInformation
Definition twolevelmethod.hh:148
LevelTransferPolicy< O, O > FatherType
Definition twolevelmethod.hh:146
void createCoarseLevelSystem(const O &fineOperator)
Algebraically creates the coarse level system.
Definition twolevelmethod.hh:154
A policy class for solving the coarse level system using one step of AMG.
Definition twolevelmethod.hh:234
OneStepAMGCoarseSolverPolicy(const SmootherArgs &args, const Criterion &c)
Constructs the coarse solver policy.
Definition twolevelmethod.hh:253
AMGInverseOperator CoarseLevelSolver
The type of solver constructed for the coarse level.
Definition twolevelmethod.hh:315
OneStepAMGCoarseSolverPolicy(const OneStepAMGCoarseSolverPolicy &other)
Copy constructor.
Definition twolevelmethod.hh:257
O::range_type X
The type of the range and domain of the operator.
Definition twolevelmethod.hh:239
C Criterion
The type of the crition used for the aggregation within AMG.
Definition twolevelmethod.hh:241
Dune::Amg::SmootherTraits< S >::Arguments SmootherArgs
The type of the arguments used for constructing the smoother.
Definition twolevelmethod.hh:245
O Operator
The type of the linear operator used.
Definition twolevelmethod.hh:237
AMG< Operator, X, Smoother > AMGType
The type of the AMG construct on the coarse level.
Definition twolevelmethod.hh:247
CoarseLevelSolver * createCoarseLevelSolver(P &transferPolicy)
Constructs a coarse level solver.
Definition twolevelmethod.hh:325
S Smoother
The type of the smoother used in AMG.
Definition twolevelmethod.hh:243
Definition twolevelmethod.hh:353
CoarseOperatorType::range_type CoarseRangeType
The type of the range of the coarse level operator.
Definition twolevelmethod.hh:380
FineOperatorType::domain_type FineDomainType
The type of the domain of the fine level operator.
Definition twolevelmethod.hh:371
TwoLevelMethod(const TwoLevelMethod &other)
Definition twolevelmethod.hh:418
void pre(FineDomainType &x, FineRangeType &b)
Definition twolevelmethod.hh:431
FO FineOperatorType
The linear operator of the finel level system. Has to be derived from AssembledLinearOperator.
Definition twolevelmethod.hh:363
CoarseLevelSolverPolicy::CoarseLevelSolver CoarseLevelSolver
The type of the coarse level solver.
Definition twolevelmethod.hh:358
void apply(FineDomainType &v, const FineRangeType &d)
Definition twolevelmethod.hh:439
CSP CoarseLevelSolverPolicy
The type of the policy for constructing the coarse level solver.
Definition twolevelmethod.hh:356
CoarseOperatorType::domain_type CoarseDomainType
The type of the domain of the coarse level operator.
Definition twolevelmethod.hh:384
TwoLevelMethod(const FineOperatorType &op, std::shared_ptr< SmootherType > smoother, const LevelTransferPolicy< FineOperatorType, CoarseOperatorType > &policy, CoarseLevelSolverPolicy &coarsePolicy, std::size_t preSteps=1, std::size_t postSteps=1)
Constructs a two level method.
Definition twolevelmethod.hh:404
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition twolevelmethod.hh:466
FineOperatorType::range_type FineRangeType
The type of the range of the fine level operator.
Definition twolevelmethod.hh:367
~TwoLevelMethod()
Definition twolevelmethod.hh:424
CSP::Operator CoarseOperatorType
The linear operator of the finel level system. Has to be derived from AssembledLinearOperator.
Definition twolevelmethod.hh:376
void post(FineDomainType &x)
Definition twolevelmethod.hh:436
S SmootherType
The type of the fine level smoother.
Definition twolevelmethod.hh:388
Base class for matrix free definition of preconditioners.
Definition preconditioner.hh:32
Statistics about the application of an inverse operator.
Definition solver.hh:48
Abstract base class for all solvers.
Definition solver.hh:99
Category
Definition solvercategory.hh:23
@ sequential
Category for sequential solvers.
Definition solvercategory.hh:25