Dune Core Modules (2.5.2)

preconditioners.hh
Go to the documentation of this file.
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_PRECONDITIONERS_HH
4#define DUNE_ISTL_PRECONDITIONERS_HH
5
6#include <cmath>
7#include <complex>
8#include <iostream>
9#include <iomanip>
10#include <string>
11
12#include <dune/common/unused.hh>
13
14#include "preconditioner.hh"
15#include "solver.hh"
16#include "solvercategory.hh"
17#include "istlexception.hh"
18#include "matrixutils.hh"
19#include "gsetc.hh"
20#include "ilu.hh"
21
22
23namespace Dune {
66 template<class O, int c>
68 public Preconditioner<typename O::domain_type, typename O::range_type>
69 {
70 public:
72 typedef typename O::domain_type domain_type;
74 typedef typename O::range_type range_type;
76 typedef typename range_type::field_type field_type;
77 typedef O InverseOperator;
78
79 // define the category
80 enum {
82 category=c
83 };
84
89 InverseOperator2Preconditioner(InverseOperator& inverse_operator)
90 : inverse_operator_(inverse_operator)
91 {}
92
94 {}
95
96 void apply(domain_type& v, const range_type& d)
97 {
99 range_type copy(d);
100 inverse_operator_.apply(v, copy, res);
101 }
102
104 {}
105
106 private:
107 InverseOperator& inverse_operator_;
108 };
109
110 //=====================================================================
111 // Implementation of this interface for sequential ISTL-preconditioners
112 //=====================================================================
113
114
126 template<class M, class X, class Y, int l=1>
127 class SeqSSOR : public Preconditioner<X,Y> {
128 public:
130 typedef M matrix_type;
132 typedef X domain_type;
134 typedef Y range_type;
136 typedef typename X::field_type field_type;
137
138 // define the category
139 enum {
142 };
143
151 SeqSSOR (const M& A, int n, field_type w)
152 : _A_(A), _n(n), _w(w)
153 {
155 }
156
162 virtual void pre (X& x, Y& b)
163 {
166
167 }
168
174 virtual void apply (X& v, const Y& d)
175 {
176 for (int i=0; i<_n; i++) {
177 bsorf(_A_,v,d,_w,BL<l>());
178 bsorb(_A_,v,d,_w,BL<l>());
179 }
180 }
181
187 virtual void post (X& x)
188 {
190 }
191
192 private:
194 const M& _A_;
196 int _n;
198 field_type _w;
199 };
200
201
202
214 template<class M, class X, class Y, int l=1>
215 class SeqSOR : public Preconditioner<X,Y> {
216 public:
218 typedef M matrix_type;
220 typedef X domain_type;
222 typedef Y range_type;
224 typedef typename X::field_type field_type;
225
226 // define the category
227 enum {
230 };
231
239 SeqSOR (const M& A, int n, field_type w)
240 : _A_(A), _n(n), _w(w)
241 {
243 }
244
250 virtual void pre (X& x, Y& b)
251 {
254 }
255
261 virtual void apply (X& v, const Y& d)
262 {
263 this->template apply<true>(v,d);
264 }
265
274 template<bool forward>
275 void apply(X& v, const Y& d)
276 {
277 if(forward)
278 for (int i=0; i<_n; i++) {
279 bsorf(_A_,v,d,_w,BL<l>());
280 }
281 else
282 for (int i=0; i<_n; i++) {
283 bsorb(_A_,v,d,_w,BL<l>());
284 }
285 }
286
292 virtual void post (X& x)
293 {
295 }
296
297 private:
299 const M& _A_;
301 int _n;
303 field_type _w;
304 };
305
306
317 template<class M, class X, class Y, int l=1>
318 class SeqGS : public Preconditioner<X,Y> {
319 public:
321 typedef M matrix_type;
323 typedef X domain_type;
325 typedef Y range_type;
327 typedef typename X::field_type field_type;
328
329 // define the category
330 enum {
333 };
334
342 SeqGS (const M& A, int n, field_type w)
343 : _A_(A), _n(n), _w(w)
344 {
346 }
347
353 virtual void pre (X& x, Y& b)
354 {
357 }
358
364 virtual void apply (X& v, const Y& d)
365 {
366 for (int i=0; i<_n; i++) {
367 dbgs(_A_,v,d,_w,BL<l>());
368 }
369 }
370
376 virtual void post (X& x)
377 {
379 }
380
381 private:
383 const M& _A_;
385 int _n;
387 field_type _w;
388 };
389
390
401 template<class M, class X, class Y, int l=1>
403 public:
405 typedef M matrix_type;
407 typedef X domain_type;
409 typedef Y range_type;
411 typedef typename X::field_type field_type;
412
413 // define the category
414 enum {
417 };
418
426 SeqJac (const M& A, int n, field_type w)
427 : _A_(A), _n(n), _w(w)
428 {
430 }
431
437 virtual void pre (X& x, Y& b)
438 {
441 }
442
448 virtual void apply (X& v, const Y& d)
449 {
450 for (int i=0; i<_n; i++) {
451 dbjac(_A_,v,d,_w,BL<l>());
452 }
453 }
454
460 virtual void post (X& x)
461 {
463 }
464
465 private:
467 const M& _A_;
469 int _n;
471 field_type _w;
472 };
473
474
475
487 template<class M, class X, class Y, int l=1>
489 public:
491 typedef typename std::remove_const<M>::type matrix_type;
493 typedef X domain_type;
495 typedef Y range_type;
497 typedef typename X::field_type field_type;
498
499 // define the category
500 enum {
503 };
504
511 SeqILU0 (const M& A, field_type w)
512 : ILU(A) // copy A
513 {
514 _w =w;
516 }
517
523 virtual void pre (X& x, Y& b)
524 {
527 }
528
534 virtual void apply (X& v, const Y& d)
535 {
536 bilu_backsolve(ILU,v,d);
537 v *= _w;
538 }
539
545 virtual void post (X& x)
546 {
548 }
549
550 private:
552 field_type _w;
554 matrix_type ILU;
555 };
556
557
571 template<class M, class X, class Y, int l=1>
573 public:
575 typedef typename std::remove_const<M>::type matrix_type;
577 typedef X domain_type;
579 typedef Y range_type;
581 typedef typename X::field_type field_type;
582
583 // define the category
584 enum {
587 };
588
596 SeqILUn (const M& A, int n, field_type w)
597 : ILU(A.N(),A.M(),M::row_wise)
598 {
599 _n = n;
600 _w = w;
601 bilu_decomposition(A,n,ILU);
602 }
603
609 virtual void pre (X& x, Y& b)
610 {
613 }
614
620 virtual void apply (X& v, const Y& d)
621 {
622 bilu_backsolve(ILU,v,d);
623 v *= _w;
624 }
625
631 virtual void post (X& x)
632 {
634 }
635
636 private:
638 matrix_type ILU;
640 int _n;
642 field_type _w;
643 };
644
645
646
655 template<class X, class Y>
657 public:
659 typedef X domain_type;
661 typedef Y range_type;
663 typedef typename X::field_type field_type;
664
665 // define the category
666 enum {
669 };
670
677 {
678 _w = w;
679 }
680
686 virtual void pre (X& x, Y& b)
687 {
690 }
691
697 virtual void apply (X& v, const Y& d)
698 {
699 v = d;
700 v *= _w;
701 }
702
708 virtual void post (X& x)
709 {
711 }
712
713 private:
715 field_type _w;
716 };
717
718
719
722} // end namespace
723
724#endif
Turns an InverseOperator into a Preconditioner.
Definition: preconditioners.hh:69
O::range_type range_type
The range type of the preconditioner.
Definition: preconditioners.hh:74
void apply(domain_type &v, const range_type &d)
Apply one step of the preconditioner to the system A(v)=d.
Definition: preconditioners.hh:96
O::domain_type domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:72
range_type::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:76
@ category
The category the preconditioner is part of.
Definition: preconditioners.hh:82
InverseOperator2Preconditioner(InverseOperator &inverse_operator)
Construct the preconditioner from the solver.
Definition: preconditioners.hh:89
void pre(domain_type &, range_type &)
Prepare the preconditioner.
Definition: preconditioners.hh:93
void post(domain_type &)
Clean up.
Definition: preconditioners.hh:103
Abstract base class for all solvers.
Definition: solver.hh:79
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:26
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioner.hh:33
Richardson preconditioner.
Definition: preconditioners.hh:656
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:663
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:661
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:686
Richardson(field_type w=1.0)
Constructor.
Definition: preconditioners.hh:676
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:708
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:659
virtual void apply(X &v, const Y &d)
Apply the precondioner.
Definition: preconditioners.hh:697
Sequential Gauss Seidel preconditioner.
Definition: preconditioners.hh:318
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:376
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:353
SeqGS(const M &A, int n, field_type w)
Constructor.
Definition: preconditioners.hh:342
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:325
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:327
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition: preconditioners.hh:364
M matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:321
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:323
Sequential ILU0 preconditioner.
Definition: preconditioners.hh:488
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:545
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:523
virtual void apply(X &v, const Y &d)
Apply the preconditoner.
Definition: preconditioners.hh:534
SeqILU0(const M &A, field_type w)
Constructor.
Definition: preconditioners.hh:511
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:495
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:493
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:497
std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:491
Sequential ILU(n) preconditioner.
Definition: preconditioners.hh:572
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:609
std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:575
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:631
SeqILUn(const M &A, int n, field_type w)
Constructor.
Definition: preconditioners.hh:596
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:577
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:579
virtual void apply(X &v, const Y &d)
Apply the precondioner.
Definition: preconditioners.hh:620
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:581
The sequential jacobian preconditioner.
Definition: preconditioners.hh:402
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:460
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition: preconditioners.hh:448
M matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:405
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:411
SeqJac(const M &A, int n, field_type w)
Constructor.
Definition: preconditioners.hh:426
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:437
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:407
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:409
Sequential SOR preconditioner.
Definition: preconditioners.hh:215
M matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:218
void apply(X &v, const Y &d)
Apply the preconditioner in a special direction.
Definition: preconditioners.hh:275
SeqSOR(const M &A, int n, field_type w)
Constructor.
Definition: preconditioners.hh:239
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:220
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:292
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:250
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition: preconditioners.hh:261
@ category
The category the preconditioner is part of.
Definition: preconditioners.hh:229
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:222
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:224
Sequential SSOR preconditioner.
Definition: preconditioners.hh:127
virtual void post(X &x)
Clean up.
Definition: preconditioners.hh:187
SeqSSOR(const M &A, int n, field_type w)
Constructor.
Definition: preconditioners.hh:151
X::field_type field_type
The field type of the preconditioner.
Definition: preconditioners.hh:136
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:132
M matrix_type
The matrix type the preconditioner is for.
Definition: preconditioners.hh:130
@ category
The category the preconditioner is part of.
Definition: preconditioners.hh:141
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition: preconditioners.hh:174
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: preconditioners.hh:162
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:134
void bsorb(const M &A, X &x, const Y &b, const K &w)
SSOR step.
Definition: gsetc.hh:591
void dbjac(const M &A, X &x, const Y &b, const K &w)
Jacobi step.
Definition: gsetc.hh:603
void dbgs(const M &A, X &x, const Y &b, const K &w)
GS step.
Definition: gsetc.hh:567
void bilu_backsolve(const M &A, X &v, const Y &d)
LU backsolve with stored inverse.
Definition: ilu.hh:96
void bilu_decomposition(const M &A, int n, M &ILU)
Definition: ilu.hh:156
void bilu0_decomposition(M &A)
compute ILU decomposition of A. A is overwritten by its decomposition
Definition: ilu.hh:35
void bsorf(const M &A, X &x, const Y &b, const K &w)
SOR step.
Definition: gsetc.hh:579
Simple iterative methods like Jacobi, Gauss-Seidel, SOR, SSOR, etc. in a generic way.
Some handy generic functions for ISTL matrices.
Dune namespace.
Definition: alignment.hh:11
Define general, extensible interface for inverse operators.
compile-time parameter for block recursion depth
Definition: gsetc.hh:40
static void check(const Matrix &mat)
Check whether the a matrix has diagonal values on blocklevel recursion levels.
Definition: matrixutils.hh:88
Statistics about the application of an inverse operator.
Definition: solver.hh:32
@ sequential
Category for sequential solvers.
Definition: solvercategory.hh:21
Definition of the DUNE_UNUSED macro for the case that config.h is not available.
#define DUNE_UNUSED_PARAMETER(parm)
A macro to mark intentionally unused function parameters with.
Definition: unused.hh:18
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 13, 23:29, 2024)