dune-pdelab  2.7-git
newtonlinesearch.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_PDELAB_SOLVER_NEWTONLINESEARCH_HH
4 #define DUNE_PDELAB_SOLVER_NEWTONLINESEARCH_HH
5 
7 
8 namespace Dune::PDELab
9 {
11  template <typename Domain>
13  {
14  public:
16  virtual ~LineSearchInterface () {}
17 
19  virtual void lineSearch(Domain&, const Domain&) = 0;
20 
22  virtual void setParameters(const ParameterTree&) = 0;
23  };
24 
25 
27  template <typename Newton>
28  class LineSearchNone : public LineSearchInterface<typename Newton::Domain>
29  {
30  public:
31  using Domain = typename Newton::Domain;
32  using Real = typename Newton::Real;
33 
34  LineSearchNone(Newton& newton) : _newton(newton) {}
35 
37  virtual void lineSearch(Domain& solution, const Domain& correction) override
38  {
39  solution.axpy(-1.0, correction);
40  _newton.updateDefect(solution);
41  }
42 
43  virtual void setParameters(const ParameterTree&) override {}
44 
45  private:
46  Newton& _newton;
47  };
48 
49 
56  template <typename Newton>
57  class LineSearchHackbuschReusken : public LineSearchInterface<typename Newton::Domain>
58  {
59  public:
60  using Domain = typename Newton::Domain;
61  using Real = typename Newton::Real;
62 
63  LineSearchHackbuschReusken(Newton& newton) : _newton(newton) {}
64 
66  virtual void lineSearch(Domain& solution, const Domain& correction) override
67  {
68  if ((_newton.result().defect < _newton.getAbsoluteLimit())){
69  solution.axpy(-1.0, correction);
70  _newton.updateDefect(solution);
71  return;
72  }
73 
74  auto verbosity = _newton.getVerbosityLevel();
75 
76  if (verbosity >= 4)
77  std::cout << " Performing line search..." << std::endl;
78  Real lambda = 1.0;
79  Real bestLambda = 0.0;
80  Real bestDefect = _newton.result().defect;
81  Real previousDefect = _newton.result().defect;
82  bool converged = false;
83 
84  if (not _previousSolution)
85  _previousSolution = std::make_shared<Domain>(solution);
86  else
87  *_previousSolution = solution;
88 
89  for (unsigned int iteration = 0; iteration < _lineSearchMaxIterations; ++iteration){
90  if (verbosity >= 4)
91  std::cout << " trying line search damping factor: "
92  << std::setw(12) << std::setprecision(4) << std::scientific
93  << lambda
94  << std::endl;
95 
96  solution.axpy(-lambda, correction);
97  _newton.updateDefect(solution);
98  if (verbosity >= 4){
99  if (not std::isfinite(_newton.result().defect))
100  std::cout << " NaNs detected" << std::endl;
101  } // ignore NaNs and try again with lower lambda
102 
103  if (_newton.result().defect <= (1.0 - lambda/4) * previousDefect){
104  if (verbosity >= 4)
105  std::cout << " line search converged" << std::endl;
106  converged = true;
107  break;
108  }
109 
110  if (_newton.result().defect < bestDefect){
111  bestDefect = _newton.result().defect;
112  bestLambda = lambda;
113  }
114 
115  lambda *= _lineSearchDampingFactor;
116  solution = *_previousSolution;
117  }
118 
119  if (not converged){
120  if (verbosity >= 4)
121  std::cout << " max line search iterations exceeded" << std::endl;
122 
123  if (not _acceptBest){
124  solution = *_previousSolution;
125  _newton.updateDefect(solution);
126  DUNE_THROW(NewtonLineSearchError,
127  "NewtonLineSearch::line_search(): line search failed, "
128  "max iteration count reached, "
129  "defect did not improve enough");
130  }
131  else{
132  if (bestLambda == 0.0){
133  solution = *_previousSolution;
134  _newton.updateDefect(solution);
135  DUNE_THROW(NewtonLineSearchError,
136  "NewtonLineSearch::line_search(): line search failed, "
137  "max iteration count reached, "
138  "defect did not improve in any of the iterations");
139  }
140  if (bestLambda != lambda){
141  solution = *_previousSolution;
142  solution.axpy(-bestLambda, correction);
143  _newton.updateDefect(solution);
144  converged = true;
145  }
146  }
147  }
148 
149  if (converged)
150  if (verbosity >= 4)
151  std::cout << " line search damping factor: "
152  << std::setw(12) << std::setprecision(4) << std::scientific
153  << lambda << std::endl;
154  }
155 
156 
157  /* \brief Set parameters
158  *
159  * Possible parameters are:
160  *
161  * - line_search_max_iterations: Maximum number of line search iterations.
162  *
163  * - line_search_damping_factor: Multiplier to line search parameter after each iteration.
164  *
165  * - line_search_accept_best: Accept the best line search parameter if
166  * there was any improvement, even if the convergence criterion was not
167  * reached.
168  */
169  virtual void setParameters(const ParameterTree& parameterTree) override
170  {
171  _lineSearchMaxIterations = parameterTree.get<unsigned int>("line_search_max_iterations",
172  _lineSearchMaxIterations);
173  _lineSearchDampingFactor = parameterTree.get<Real>("line_search_damping_factor",
174  _lineSearchDampingFactor);
175  _acceptBest = parameterTree.get<bool>("line_search_accept_best",
176  _acceptBest);
177  }
178 
179  private:
180  Newton& _newton;
181  std::shared_ptr<Domain> _previousSolution;
182 
183  // Line search parameters
184  unsigned int _lineSearchMaxIterations = 10;
185  Real _lineSearchDampingFactor = 0.5;
186  bool _acceptBest = false;
187  };
188 
191  {
194  };
195 
196 
204  {
205  if (name == "noLineSearch")
207  if (name == "hackbusch_reusken")
209  DUNE_THROW(Exception,"Unkown line search strategy: " << name);
210  }
211 
212 
223  template <typename Newton>
224  std::shared_ptr<LineSearchInterface<typename Newton::Domain>>
225  getLineSearch(Newton& newton, const std::string& name)
226  {
227  auto strategy = lineSearchStrategyFromString(name);
228  if (strategy == LineSearchStrategy::noLineSearch){
229  auto lineSearch = std::make_shared<LineSearchNone<Newton>> (newton);
230  return lineSearch;
231  }
232  if (strategy == LineSearchStrategy::hackbuschReusken){
233  auto lineSearch = std::make_shared<LineSearchHackbuschReusken<Newton>> (newton);
234  return lineSearch;
235  }
236  DUNE_THROW(Exception,"Unkown line search strategy");
237  }
238 
239 } // namespace Dune::PDELab
240 
241 #endif
Definition: adaptivity.hh:29
LineSearchStrategy lineSearchStrategyFromString(const std::string &name)
Get a LineSearchStrategy from a string identifier.
Definition: newtonlinesearch.hh:203
std::shared_ptr< LineSearchInterface< typename Newton::Domain > > getLineSearch(Newton &newton, const std::string &name)
Get a pointer to a line search.
Definition: newtonlinesearch.hh:225
LineSearchStrategy
Flags for different line search strategies.
Definition: newtonlinesearch.hh:191
@ noLineSearch
Definition: newtonlinesearch.hh:192
@ hackbuschReusken
Definition: newtonlinesearch.hh:193
Base class for all PDELab exceptions.
Definition: exceptions.hh:19
Abstract base class describing the line search interface.
Definition: newtonlinesearch.hh:13
virtual ~LineSearchInterface()
Every abstract base class should have a virtual destructor.
Definition: newtonlinesearch.hh:16
virtual void setParameters(const ParameterTree &)=0
Set parameters.
virtual void lineSearch(Domain &, const Domain &)=0
Do line search.
Class for simply updating the solution without line search.
Definition: newtonlinesearch.hh:29
typename Newton::Domain Domain
Definition: newtonlinesearch.hh:31
virtual void setParameters(const ParameterTree &) override
Set parameters.
Definition: newtonlinesearch.hh:43
typename Newton::Real Real
Definition: newtonlinesearch.hh:32
virtual void lineSearch(Domain &solution, const Domain &correction) override
Do line search (in this case just update the solution)
Definition: newtonlinesearch.hh:37
LineSearchNone(Newton &newton)
Definition: newtonlinesearch.hh:34
Hackbusch-Reusken line search.
Definition: newtonlinesearch.hh:58
typename Newton::Real Real
Definition: newtonlinesearch.hh:61
virtual void setParameters(const ParameterTree &parameterTree) override
Set parameters.
Definition: newtonlinesearch.hh:169
typename Newton::Domain Domain
Definition: newtonlinesearch.hh:60
LineSearchHackbuschReusken(Newton &newton)
Definition: newtonlinesearch.hh:63
virtual void lineSearch(Domain &solution, const Domain &correction) override
Do line search.
Definition: newtonlinesearch.hh:66