ITL_Preconditioner.h 5.87 KB
Newer Older
1
2
3
4
// ============================================================================
// ==                                                                        ==
// == AMDiS - Adaptive multidimensional simulations                          ==
// ==                                                                        ==
5
// ==  http://www.amdis-fem.org                                              ==
6
7
// ==                                                                        ==
// ============================================================================
8
9
10
11
12
13
14
15
16
17
18
19
//
// Software License for AMDiS
//
// Copyright (c) 2010 Dresden University of Technology 
// All rights reserved.
// Authors: Simon Vey, Thomas Witkowski et al.
//
// This file is part of AMDiS
//
// See also license.opensource.txt in the distribution.


20
21
22
23
24
25

/** \file ITL_Preconditioner.h */

#ifndef AMDIS_ITL_PRECONDITIONER_H
#define AMDIS_ITL_PRECONDITIONER_H

26
#include "MTL4Types.h"
27
28
29
30
31
32
#include "DOFMatrix.h"
#include "CreatorInterface.h"

#include <boost/numeric/itl/itl.hpp>
#include <boost/numeric/itl/pc/ilu_0.hpp>
#include <boost/numeric/itl/pc/ic_0.hpp>
Praetorius, Simon's avatar
Praetorius, Simon committed
33
#include <boost/numeric/mtl/vector/assigner.hpp>
34
35
36
37
38
39
40
41
42
43
44
45
46

namespace AMDiS {

  // ============================================================================
  // ===== class ITL_Preconditioner ==============================================
  // ============================================================================

 
  /**
   * \ingroup Solver
   * 
   * \brief Common base class for wrappers to use ITL preconditioners in AMDiS.
   */
47
  template< class Vec >
48
49
50
51
52
  class ITL_BasePreconditioner 
  {
  public:
      virtual ~ITL_BasePreconditioner() {}

Praetorius, Simon's avatar
Praetorius, Simon committed
53
54
55
      // Old style AMDiS
      // virtual Vec member_solve(const Vec& vin) const =0 ;
      // virtual Vec member_adjoint_solve(const Vec& vin) const =0 ;
56

Praetorius, Simon's avatar
Praetorius, Simon committed
57
58
59
     // New style without copying the result vector
     virtual void solve(const Vec& x, Vec& y) const= 0;
     virtual void adjoint_solve(const Vec& x, Vec& y) const= 0;
60
61
  };

Praetorius, Simon's avatar
Praetorius, Simon committed
62
63


64
  template< typename Vec >
Praetorius, Simon's avatar
Praetorius, Simon committed
65
66
  inline //Vec
  itl::pc::solver<ITL_BasePreconditioner< Vec >, Vec, false>
67
68
  solve(const ITL_BasePreconditioner< Vec >& P, const Vec& vin)
  {
Praetorius, Simon's avatar
Praetorius, Simon committed
69
70
    return itl::pc::solver<ITL_BasePreconditioner< Vec >, Vec, false>(P, vin);
    //return P.member_solve(vin);
71
  }
72

73
  template< typename Vec >
Praetorius, Simon's avatar
Praetorius, Simon committed
74
75
  inline //Vec
  itl::pc::solver<ITL_BasePreconditioner< Vec >, Vec, true>
76
77
  adjoint_solve(const ITL_BasePreconditioner< Vec >& P, const Vec& vin)
  {
Praetorius, Simon's avatar
Praetorius, Simon committed
78
79
    return itl::pc::solver<ITL_BasePreconditioner< Vec >, Vec, true>(P, vin);
    //return P.member_adjoint_solve(vin);
80
  }
81
82
83
84
85
86

  /**
   * \ingroup Solver
   * 
   * \brief Wrapper for using ITL preconditioners in AMDiS.
   */
87
88
  template <typename Preconditioner, typename MTLVector, typename MTLMatrix >
  class ITL_Preconditioner : public ITL_BasePreconditioner< MTLVector >
89
90
  {
  public:
91
    ITL_Preconditioner(const MTLMatrix& A) : precon(A) {}
Praetorius, Simon's avatar
Praetorius, Simon committed
92
#if 0   
93
94
    MTLVector
    member_solve(const MTLVector& vin) const
95
96
97
98
    {
      return solve(precon, vin);
    }

99
100
    MTLVector
    member_adjoint_solve(const MTLVector& vin) const
101
102
103
    {
      return adjoint_solve(precon, vin);
    }
Praetorius, Simon's avatar
Praetorius, Simon committed
104
105
106
107
108
109
110
111
112
113
114
115
#else
      
    void solve(const MTLVector& vin, MTLVector& vout) const
    {
      precon.solve(vin, vout);
   }

    void adjoint_solve(const MTLVector& vin, MTLVector& vout) const
    {
      precon.adjoint_solve(vin, vout);
   }
#endif
116
117

    /// Creator class used in the OEMSolverMap.
118
    class Creator : public CreatorInterface< ITL_BasePreconditioner< MTLVector > >
119
120
    {
    public:
Thomas Witkowski's avatar
Thomas Witkowski committed
121
      virtual ~Creator() {}
122
123
124
125
      
      /** \brief
       * Creates an ITL preconditioner
       */
126
127
      ITL_BasePreconditioner< MTLVector >* create(const MTLMatrix& A) { 
	return new ITL_Preconditioner<Preconditioner, MTLVector, MTLMatrix>(A);
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
      }
    };

  private:
      Preconditioner precon;
  };


  // ============================================================================
  // ===== class DiagonalPreconditioner =========================================
  // ============================================================================

  /**
   * \ingroup Solver
   * 
   * \brief
   * Diagonal preconditioner. 
   */
  class DiagonalPreconditioner 
147
    : public ITL_Preconditioner<itl::pc::diagonal<DOFMatrix::base_matrix_type>, MTLTypes::MTLVector, MTLTypes::MTLMatrix >
148
  {
149
    typedef ITL_Preconditioner<itl::pc::diagonal<DOFMatrix::base_matrix_type>, MTLTypes::MTLVector, MTLTypes::MTLMatrix > base;
150

Thomas Witkowski's avatar
Thomas Witkowski committed
151
  public:
152
    DiagonalPreconditioner(const MTLTypes::MTLMatrix& A) : base(A) {}
153
154
155
156
157
158
159
160
161
  };

  /**
   * \ingroup Solver
   * 
   * \brief
   * Identity preconditioner. Behaves like no preconditioning.
   */
  class IdentityPreconditioner 
162
    : public ITL_Preconditioner<itl::pc::identity<DOFMatrix::base_matrix_type>, MTLTypes::MTLVector, MTLTypes::MTLMatrix >
163
  {
164
    typedef ITL_Preconditioner<itl::pc::identity<DOFMatrix::base_matrix_type>, MTLTypes::MTLVector, MTLTypes::MTLMatrix > base;
Thomas Witkowski's avatar
Thomas Witkowski committed
165

166
  public:
167
      IdentityPreconditioner(const MTLTypes::MTLMatrix& A) : base(A) {}
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
  };


  // ============================================================================
  // ===== class ILUPreconditioner ==============================================
  // ============================================================================

  /**
   * \ingroup Solver
   * 
   * \brief
   * ILU (Incomplete LU factorization) preconditioner.
   *
   * The preconditioner is used from ITL. It corresponds for instance to "Iterative Methods for
   * Sparce Linear Systems", second edition, Yousef Saad. The preconditioner is
   * described in chapter 10.3 (algorithm 10.4).
   */
  class ILUPreconditioner 
186
      : public ITL_Preconditioner< itl::pc::ilu_0<DOFMatrix::base_matrix_type>, MTLTypes::MTLVector, MTLTypes::MTLMatrix >
187
188
189
  {};

  class ICPreconditioner 
190
      : public ITL_Preconditioner< itl::pc::ic_0<DOFMatrix::base_matrix_type>, MTLTypes::MTLVector, MTLTypes::MTLMatrix >
191
192
193
194
195
196
  {};


} // namespace AMDiS

#endif // AMDIS_ITL_PRECONDITIONER_H