DUNE PDELab (git)

linesearch.hh
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_LINESEARCH_HH
4#define DUNE_PDELAB_SOLVER_LINESEARCH_HH
5
6#include <dune/pdelab/solver/newtonerrors.hh>
7
8
9namespace Dune::PDELab
10{
11
13 template <typename Domain>
15 {
16 public:
19
21 virtual void lineSearch(Domain&, const Domain&) = 0;
22
24 virtual void setParameters(const ParameterTree&) = 0;
25
27 virtual void printParameters() const
28 {
29 std::cout << "NewtonMethod::_lineSearch->printParameters() is not implemented." << std::endl;
30 }
31 };
32
33
35 template <typename Solver>
36 class LineSearchNone : public LineSearchInterface<typename Solver::Domain>
37 {
38 public:
39 using Domain = typename Solver::Domain;
40 using Real = typename Solver::Real;
41
42 LineSearchNone(Solver& solver) : _solver(solver) {}
43
45 virtual void lineSearch(Domain& solution, const Domain& correction) override
46 {
47 solution.axpy(-1.0, correction);
48 _solver.updateDefect(solution);
49 }
50
51 virtual void setParameters(const ParameterTree&) override {}
52
53 // print line search type
54 virtual void printParameters() const override
55 {
56 std::cout << "LineSearch.Type........... None" << std::endl;
57 }
58
59 private:
60 Solver& _solver;
61 };
62
63
70 template <typename Solver>
71 class LineSearchHackbuschReusken : public LineSearchInterface<typename Solver::Domain>
72 {
73 public:
74 using Domain = typename Solver::Domain;
75 using Real = typename Solver::Real;
76
77 LineSearchHackbuschReusken(Solver& solver, bool forceAcceptBest = false) :
78 _solver(solver), _forceAcceptBest(forceAcceptBest) {}
79
81 virtual void lineSearch(Domain& solution, const Domain& correction) override
82 {
83 if ((_solver.result().defect < _solver.getAbsoluteLimit())){
84 solution.axpy(-1.0, correction);
85 _solver.updateDefect(solution);
86 return;
87 }
88
89 auto verbosity = _solver.getVerbosityLevel();
90
91 if (verbosity >= 4)
92 std::cout << " Performing line search..." << std::endl;
93 Real lambda = 1.0;
94 Real bestLambda = 0.0;
95 Real bestDefect = _solver.result().defect;
96 Real previousDefect = _solver.result().defect;
97 bool converged = false;
98
99 if (not _previousSolution)
100 _previousSolution = std::make_shared<Domain>(solution);
101 else
102 *_previousSolution = solution;
103
104 for (unsigned int iteration = 0; iteration < _lineSearchMaxIterations; ++iteration){
105 if (verbosity >= 4)
106 std::cout << " trying line search damping factor: "
107 << std::setw(12) << std::setprecision(4) << std::scientific
108 << lambda
109 << std::endl;
110
111 solution.axpy(-lambda, correction);
112 _solver.updateDefect(solution);
113 if (verbosity >= 4){
114 if (not std::isfinite(_solver.result().defect))
115 std::cout << " NaNs detected" << std::endl;
116 } // ignore NaNs and try again with lower lambda
117
118 if (_solver.result().defect <= (1.0 - lambda/4) * previousDefect){
119 if (verbosity >= 4)
120 std::cout << " line search converged" << std::endl;
121 converged = true;
122 break;
123 }
124
125 if (_solver.result().defect < bestDefect){
126 bestDefect = _solver.result().defect;
127 bestLambda = lambda;
128 }
129
130 lambda *= _lineSearchDampingFactor;
131 solution = *_previousSolution;
132 }
133
134 if (not converged){
135 if (verbosity >= 4)
136 std::cout << " max line search iterations exceeded" << std::endl;
137
138 if (not (_acceptBest or _forceAcceptBest)){
139 solution = *_previousSolution;
140 _solver.updateDefect(solution);
141 DUNE_THROW( LineSearchError,
142 "LineSearch::lineSearch(): line search failed, "
143 "max iteration count reached, "
144 "defect did not improve enough");
145 }
146 else{
147 if (bestLambda == 0.0){
148 solution = *_previousSolution;
149 _solver.updateDefect(solution);
150 DUNE_THROW(LineSearchError,
151 "LineSearch::lineSearch(): line search failed, "
152 "max iteration count reached, "
153 "defect did not improve in any of the iterations");
154 }
155 if (bestLambda != lambda){
156 solution = *_previousSolution;
157 solution.axpy(-bestLambda, correction);
158 _solver.updateDefect(solution);
159 converged = true;
160 }
161 }
162 }
163
164 if (converged)
165 if (verbosity >= 4)
166 std::cout << " line search damping factor: "
167 << std::setw(12) << std::setprecision(4) << std::scientific
168 << lambda << std::endl;
169 }
170
171
172 /* \brief Set parameters
173 *
174 * Possible parameters are:
175 *
176 * - MaxIterations: Maximum number of line search iterations.
177 *
178 * - DampingFactor: Multiplier to line search parameter after each iteration.
179 *
180 * - AcceptBest: Accept the best line search parameter if
181 * there was any improvement, even if the convergence criterion was not
182 * reached.
183 */
184 virtual void setParameters(const ParameterTree& parameterTree) override
185 {
186 _lineSearchMaxIterations = parameterTree.get<unsigned int>("MaxIterations",
187 _lineSearchMaxIterations);
188 _lineSearchDampingFactor = parameterTree.get<Real>("DampingFactor",
189 _lineSearchDampingFactor);
190 _acceptBest = parameterTree.get<bool>("AcceptBest", _acceptBest);
191 }
192
193 virtual void printParameters() const override
194 {
195 std::cout << "LineSearch.Type........... Hackbusch-Reusken" << std::endl;
196 std::cout << "LineSearch.MaxIterations.. " << _lineSearchMaxIterations << std::endl;
197 std::cout << "LineSearch.DampingFactor.. " << _lineSearchDampingFactor << std::endl;
198 std::cout << "LineSearch.AcceptBest..... " << (_acceptBest or _forceAcceptBest) << std::endl;
199 }
200
201 private:
202 Solver& _solver;
203 std::shared_ptr<Domain> _previousSolution;
204
205 // Line search parameters
206 unsigned int _lineSearchMaxIterations = 10;
207 Real _lineSearchDampingFactor = 0.5;
208 bool _acceptBest = false;
209 bool _forceAcceptBest;
210 };
211
213 enum class LineSearchStrategy
214 {
215 noLineSearch,
216 hackbuschReusken,
217 hackbuschReuskenAcceptBest
218 };
219
220 // we put this into an emty namespace, so that we don't violate the one-definition-rule
221 namespace {
228 inline
229 LineSearchStrategy lineSearchStrategyFromString (const std::string& name)
230 {
231 if (name == "noLineSearch")
232 return LineSearchStrategy::noLineSearch;
233 if (name == "hackbuschReusken")
234 return LineSearchStrategy::hackbuschReusken;
235 if (name == "hackbuschReuskenAcceptBest")
236 return LineSearchStrategy::hackbuschReuskenAcceptBest;
237 DUNE_THROW(Exception,"Unkown line search strategy: " << name);
238 }
239 }
240
241
252 template <typename Solver>
253 std::shared_ptr<LineSearchInterface<typename Solver::Domain>>
254 createLineSearch(Solver& solver, LineSearchStrategy strategy)
255 {
256 if (strategy == LineSearchStrategy::noLineSearch){
257 auto lineSearch = std::make_shared<LineSearchNone<Solver>> (solver);
258 return lineSearch;
259 }
260 if (strategy == LineSearchStrategy::hackbuschReusken){
261 auto lineSearch = std::make_shared<LineSearchHackbuschReusken<Solver>> (solver);
262 return lineSearch;
263 }
264 if (strategy == LineSearchStrategy::hackbuschReuskenAcceptBest){
265 auto lineSearch = std::make_shared<LineSearchHackbuschReusken<Solver>> (solver, true);
266 std::cout << "Warning: linesearch hackbuschReuskenAcceptBest is deprecated and will be removed after PDELab 2.7.\n"
267 << " Please use 'hackbuschReusken' and add the parameter 'LineSearchAcceptBest : true'";
268 return lineSearch;
269 }
270 DUNE_THROW(Exception,"Unkown line search strategy");
271 }
272
273} // namespace Dune::PDELab
274
275#endif
Hackbusch-Reusken line search.
Definition: linesearch.hh:72
virtual void setParameters(const ParameterTree &parameterTree) override
Set parameters.
Definition: linesearch.hh:184
virtual void printParameters() const override
Print paramters.
Definition: linesearch.hh:193
virtual void lineSearch(Domain &solution, const Domain &correction) override
Do line search.
Definition: linesearch.hh:81
Abstract base class describing the line search interface.
Definition: linesearch.hh:15
virtual ~LineSearchInterface()
Every abstract base class should have a virtual destructor.
Definition: linesearch.hh:18
virtual void printParameters() const
Print paramters.
Definition: linesearch.hh:27
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: linesearch.hh:37
virtual void printParameters() const override
Print paramters.
Definition: linesearch.hh:54
virtual void lineSearch(Domain &solution, const Domain &correction) override
Do line search (in this case just update the solution)
Definition: linesearch.hh:45
virtual void setParameters(const ParameterTree &) override
Set parameters.
Definition: linesearch.hh:51
Hierarchical structure of string parameters.
Definition: parametertree.hh:37
std::string get(const std::string &key, const std::string &defaultValue) const
get value as string
Definition: parametertree.cc:188
#define DUNE_THROW(E, m)
Definition: exceptions.hh:218
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 23, 23:29, 2024)