G+Smo  24.08.0
Geometry + Simulation Modules
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
gsMinResQLP.h
Go to the documentation of this file.
1 
18 #pragma once
19 
21 
22 //TODOs: Create accessors for values of interest, like, Arnorm, xnorm etc.
23 //TODOs: Add the same flaging system the matlab code uses
24 //TODOs: Merge the minres and minresQLP code like petsc and matlab has
25 namespace gismo
26 {
27 
32 template<class T=real_t>
33 class gsMinResQLP : public gsIterativeSolver<T>
34 {
35 
36 public:
37  typedef gsIterativeSolver<T> Base;
38 
39  typedef gsMatrix<T> VectorType;
40 
41  typedef typename Base::LinOpPtr LinOpPtr;
42 
43  typedef memory::shared_ptr<gsMinResQLP> Ptr;
44  typedef memory::unique_ptr<gsMinResQLP> uPtr;
45 
50  template< typename OperatorType >
51  explicit gsMinResQLP( const OperatorType& mat,
52  const LinOpPtr& precond = LinOpPtr())
53  : Base(mat, precond), m_inexact_residual(false) {}
54 
59  template< typename OperatorType >
60  static uPtr make( const OperatorType& mat, const LinOpPtr& precond = LinOpPtr() )
61  { return uPtr( new gsMinResQLP(mat, precond) ); }
62 
63  bool initIteration( const VectorType& rhs, VectorType& x );
64  void finalizeIteration( VectorType& x );
65 
66  bool step( VectorType& x );
67 
70  {
72  opt.addSwitch( "InexactResidual",
73  "If true, the residual is estimated, not accurately computed.",
74  false );
75  return opt;
76  }
77 
80  {
81  Base::setOptions(opt);
82  m_inexact_residual = opt.askSwitch("InexactResidual", m_inexact_residual);
83  return *this;
84  }
85 
87  void setInexactResidual( bool flag ) { m_inexact_residual = flag; }
88 
90  std::ostream &print(std::ostream &os) const
91  {
92  os << "gsMinResQLP\n";
93  return os;
94  }
95 
96 private:
97  using Base::m_mat;
98  using Base::m_precond;
99  using Base::m_max_iters;
100  using Base::m_tol;
101  using Base::m_num_iter;
102  using Base::m_rhs_norm;
103  using Base::m_error;
104 
105  gsMatrix<T> m_rhs, r1, r2, r3, resvec, Aresvec,
106  v, xl2, negResidual,
107  wl, w, wl2;
108 
109 
110  T beta1, betal, beta, tau, taul, taul2, phi, betan, gmin, cs, sn,cr1,sr1,
111  cr2,sr2, deltan,gamma,gammal, gammal2, gammal3, eta,etal, etal2,
112  vepln, veplnl, veplnl2, u, ul, ul2, ul3, rnorm, rnorml, xnorm, xl2norm,
113  Axnorm, Anorm, Acond, relres, alpha, pnorm, dbar, delta,
114  epsilon, epsilonn, gbar, maxxnorm, gammal_QLP, gamma_QLP, vepln_QLP,
115  abs_gamma, gminl, gminl2, Arnorml, Arnorm, relAresl, relAres,
116  rootl, ul_QLP, u_QLP;;
117 
118 
119  index_t QLPiter;
120  bool m_inexact_residual;
121 };
122 
123 } // namespace gismo
124 
125 #ifndef GISMO_BUILD_LIB
126 #include GISMO_HPP_HEADER(gsMinResQLP.hpp)
127 #endif
const LinOpPtr m_mat
The matrix/operator to be solved for.
Definition: gsIterativeSolver.h:252
bool step(VectorType &x)
Perform one step, requires initIteration.
Definition: gsMinResQLP.hpp:109
virtual gsIterativeSolver & setOptions(const gsOptionList &opt)
Set the options based on a gsOptionList.
Definition: gsIterativeSolver.h:102
gsMinResQLP & setOptions(const gsOptionList &opt)
Set the options based on a gsOptionList.
Definition: gsMinResQLP.h:79
std::ostream & print(std::ostream &os) const
Prints the object as a string.
Definition: gsMinResQLP.h:90
gsOptionList::Real m_tol
The tolerance for m_error to be reached.
Definition: gsIterativeSolver.h:255
void finalizeIteration(VectorType &x)
Some post-processing might be required.
Definition: gsMinResQLP.hpp:306
#define index_t
Definition: gsConfig.h:32
T m_rhs_norm
The norm of the right-hand-side.
Definition: gsIterativeSolver.h:257
gsMinResQLP(const OperatorType &mat, const LinOpPtr &precond=LinOpPtr())
Constructor using a matrix (operator) and optionally a preconditionner.
Definition: gsMinResQLP.h:51
static gsOptionList defaultOptions()
Returns a list of default options.
Definition: gsMinResQLP.h:69
static gsOptionList defaultOptions()
Returns a list of default options.
Definition: gsIterativeSolver.h:92
index_t m_max_iters
The upper bound for the number of iterations.
Definition: gsIterativeSolver.h:254
bool initIteration(const VectorType &rhs, VectorType &x)
Init the iteration.
Definition: gsMinResQLP.hpp:57
static uPtr make(const OperatorType &mat, const LinOpPtr &precond=LinOpPtr())
Make function using a matrix (operator) and optionally a preconditionner.
Definition: gsMinResQLP.h:60
T m_error
The relative error as absolute_error/m_rhs_norm.
Definition: gsIterativeSolver.h:258
Abstract class for iterative solvers.
void setInexactResidual(bool flag)
If true, the residual is estimated, not accurately computed.
Definition: gsMinResQLP.h:87
bool askSwitch(const std::string &label, const bool &value=false) const
Reads value for option label from options.
Definition: gsOptionList.cpp:128
Abstract class for iterative solvers.
Definition: gsIterativeSolver.h:26
Class which holds a list of parameters/options, and provides easy access to them. ...
Definition: gsOptionList.h:32
void addSwitch(const std::string &label, const std::string &desc, const bool &value)
Adds a option named label, with description desc and value value.
Definition: gsOptionList.cpp:235
index_t m_num_iter
The number of iterations performed.
Definition: gsIterativeSolver.h:256
LinOpPtr m_precond
The preconditioner.
Definition: gsIterativeSolver.h:253
The minimal residual (MinRes-QLP) method.
Definition: gsMinResQLP.h:33