ITL_Preconditioner.h 5.1 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
33
34
35
36
37
38
39
40
41
42
43
44
45
#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>

namespace AMDiS {

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

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

52
      virtual Vec member_solve(const Vec& vin) const =0 ;
53

54
      virtual Vec member_adjoint_solve(const Vec& vin) const =0 ;
55
56
  };

57
58
59
60
61
62
  template< typename Vec >
  inline Vec
  solve(const ITL_BasePreconditioner< Vec >& P, const Vec& vin)
  {
    return P.member_solve(vin);
  }
63

64
65
66
67
68
69
  template< typename Vec >
  inline Vec
  adjoint_solve(const ITL_BasePreconditioner< Vec >& P, const Vec& vin)
  {
    return P.member_adjoint_solve(vin);
  }
70
71
72
73
74
75

  /**
   * \ingroup Solver
   * 
   * \brief Wrapper for using ITL preconditioners in AMDiS.
   */
76
77
  template <typename Preconditioner, typename MTLVector, typename MTLMatrix >
  class ITL_Preconditioner : public ITL_BasePreconditioner< MTLVector >
78
79
  {
  public:
80
    ITL_Preconditioner(const MTLMatrix& A) : precon(A) {}
81
   
82
83
    MTLVector
    member_solve(const MTLVector& vin) const
84
85
86
87
    {
      return solve(precon, vin);
    }

88
89
    MTLVector
    member_adjoint_solve(const MTLVector& vin) const
90
91
92
93
94
    {
      return adjoint_solve(precon, vin);
    }

    /// Creator class used in the OEMSolverMap.
95
    class Creator : public CreatorInterface< ITL_BasePreconditioner< MTLVector > >
96
97
    {
    public:
Thomas Witkowski's avatar
Thomas Witkowski committed
98
      virtual ~Creator() {}
99
100
101
102
      
      /** \brief
       * Creates an ITL preconditioner
       */
103
104
      ITL_BasePreconditioner< MTLVector >* create(const MTLMatrix& A) { 
	return new ITL_Preconditioner<Preconditioner, MTLVector, MTLMatrix>(A);
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
      }
    };

  private:
      Preconditioner precon;
  };


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

  /**
   * \ingroup Solver
   * 
   * \brief
   * Diagonal preconditioner. 
   */
  class DiagonalPreconditioner 
124
    : public ITL_Preconditioner<itl::pc::diagonal<DOFMatrix::base_matrix_type>, MTLTypes::MTLVector, MTLTypes::MTLMatrix >
125
  {
126
    typedef ITL_Preconditioner<itl::pc::diagonal<DOFMatrix::base_matrix_type>, MTLTypes::MTLVector, MTLTypes::MTLMatrix > base;
127

Thomas Witkowski's avatar
Thomas Witkowski committed
128
  public:
129
    DiagonalPreconditioner(const MTLTypes::MTLMatrix& A) : base(A) {}
130
131
132
133
134
135
136
137
138
  };

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

143
  public:
144
      IdentityPreconditioner(const MTLTypes::MTLMatrix& A) : base(A) {}
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
  };


  // ============================================================================
  // ===== 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 
163
      : public ITL_Preconditioner< itl::pc::ilu_0<DOFMatrix::base_matrix_type>, MTLTypes::MTLVector, MTLTypes::MTLMatrix >
164
165
166
  {};

  class ICPreconditioner 
167
      : public ITL_Preconditioner< itl::pc::ic_0<DOFMatrix::base_matrix_type>, MTLTypes::MTLVector, MTLTypes::MTLMatrix >
168
169
170
171
172
173
  {};


} // namespace AMDiS

#endif // AMDIS_ITL_PRECONDITIONER_H