KrylovRunner.hpp 3.93 KB
Newer Older
1
2
3
4
/** \file ITL_Runner.h */

#pragma once

5
6
#include <string>

7
8
9
10
11
// MTL4 headers
#include <boost/numeric/itl/itl.hpp>
#include <boost/numeric/mtl/mtl.hpp>

// AMDiS headers
12
13
14
#include <amdis/linear_algebra/PreconditionerInterface.hpp>
#include <amdis/linear_algebra/RunnerInterface.hpp>
#include <amdis/linear_algebra/SolverInfo.hpp>
15

16
17
#include <amdis/linear_algebra/mtl/ITL_Preconditioner.hpp>
#include <amdis/linear_algebra/mtl/BITL_Preconditioner.hpp>
18
19
20
21
22
23
24
25
26
27
28
29

namespace AMDiS
{
  /** \ingroup Solver
   *
   * \brief
   * Wrapper class for different MTL4 itl-solvers. These solvers
   * are parametrized by Matrix- and Vector.
   **/
  template <class ITLSolver, class Matrix, class Vector>
  class KrylovRunner : public RunnerInterface<Matrix, Vector>
  {
30
31
    using Super = RunnerInterface<Matrix, Vector>;
    using PreconBase = typename Super::PreconBase;
32

33
34
  public:
    /// Constructor.
35
    explicit KrylovRunner(std::string prefix)
36
37
38
      : itlSolver(prefix)
    {
      initPrecon(prefix);
39

40
41
      Parameters::get(prefix + "->max iteration", maxIter);
      Parameters::get(prefix + "->print cycle", printCycle);
42
43
    }

44
    /// Implements \ref RunnerInterface::getLeftPrecon(), Returns \ref L
45
    virtual std::shared_ptr<PreconBase> getLeftPrecon() override
46
47
    {
      return L;
48
49
    }

50
    /// Implements \ref RunnerInterface::getRightPrecon(), Returns \ref R
51
    virtual std::shared_ptr<PreconBase> getRightPrecon() override
52
53
54
    {
      return R;
    }
55

56
    /// Set a new left preconditioner \ref L
57
    void setLeftPrecon(std::shared_ptr<PreconBase> const& precon)
58
59
60
    {
      L = precon;
    }
61

62
    /// Set a new right preconditioner \ref R
63
    void setRightPrecon(std::shared_ptr<PreconBase> const& precon)
64
65
66
    {
      R = precon;
    }
67
68
69
70
71
72
73

    /// Implementation of \ref RunnerInterface::init()
    virtual void init(Matrix const& A) override
    {
      L->init(A);
      R->init(A);
    }
74

75
76
77
78
79
80
81
82
    /// Implementation of \ref RunnerInterface::exit()
    virtual void exit() override
    {
      L->exit();
      R->exit();
    }

    /// Implementation of \ref RunnerInterface::solve()
83
    virtual int solve(Matrix const& A, Vector& x, Vector const& b,
84
                      SolverInfo& solverInfo) override
85
    {
86
87
      test_exit(bool(L), "There is no left preconditioner");
      test_exit(bool(R), "There is no right preconditioner");
88

89
      auto r = Vector(b);
90
91
      if (two_norm(x) != 0)
        r -= A * x;
92

93
      // print information about the solution process
94
      itl::cyclic_iteration<double> iter(r, maxIter, solverInfo.getRelTolerance(),
95
                                                     solverInfo.getAbsTolerance(),
96
                                                     printCycle);
97
      iter.set_quite(solverInfo.getInfo() == 0);
98
      iter.suppress_resume(solverInfo.getInfo() == 0);
99
100

      int error = itlSolver(A, x, b, *L, *R, iter);
101
102
      solverInfo.setAbsResidual(iter.resid());
      solverInfo.setRelResidual(iter.relresid());
103
104
105
106

      return error;
    }

107
  private:
108
109
    /// Create left/right preconditioners from parameters given in the init-file
    void initPrecon(std::string prefix)
110
    {
111
112
113
114
      // Creator for the left preconditioner
      std::string leftPrecon = "", rightPrecon = "";
      Parameters::get(prefix + "->left precon", leftPrecon);
      Parameters::get(prefix + "->right precon", rightPrecon);
115

116
117
118
119
      auto leftCreator
        = CreatorMap<PreconBase>::getCreator(leftPrecon, prefix + "->left precon");
      auto rightCreator
        = CreatorMap<PreconBase>::getCreator(rightPrecon, prefix + "->right precon");
120

121
122
      L = leftCreator->create();
      R = rightCreator->create();
123
124
    }

125
  private:
126
127
    /// The itl solver object that performes the actual solve
    ITLSolver itlSolver;
128

129
    /// Pointer to the left and right preconditioner
130
    std::shared_ptr<PreconBase> L, R;
131

132
    /// The maximal number of iterations
133
    std::size_t maxIter = 1000;
134

135
    /// Print solver residuals every \ref printCycle 's iteration
136
    std::size_t printCycle = 100;
137

138
139
140
  };

} // end namespace AMDiS