dune-pdelab  2.0.0
newton.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_NEWTON_HH
4 #define DUNE_PDELAB_NEWTON_HH
5 
6 #include <iostream>
7 #include <iomanip>
8 #include <cmath>
9 
10 #include <math.h>
11 
12 #include <dune/common/exceptions.hh>
13 #include <dune/common/ios_state.hh>
14 #include <dune/common/timer.hh>
15 #include <dune/common/parametertree.hh>
16 
18 
19 namespace Dune
20 {
21  namespace PDELab
22  {
23  // Exception classes used in NewtonSolver
24  class NewtonError : public Exception {};
25  class NewtonDefectError : public NewtonError {};
28  class NewtonNotConverged : public NewtonError {};
29 
30  // Status information of Newton's method
31  template<class RFType>
33  {
34  RFType first_defect; // the first defect
35  RFType defect; // the final defect
36  double assembler_time; // Cumulative time for matrix assembly
37  double linear_solver_time; // Cumulative time for linear sovler
38  int linear_solver_iterations; // Total number of linear iterations
39 
43  };
44 
45  template<class GOS, class TrlV, class TstV>
46  class NewtonBase
47  {
48  typedef GOS GridOperator;
49  typedef TrlV TrialVector;
50  typedef TstV TestVector;
51 
52  typedef typename TestVector::ElementType RFType;
53  typedef typename GOS::Traits::Jacobian Matrix;
54 
55 
56  public:
57  // export result type
59 
60  void setVerbosityLevel(unsigned int verbosity_level_)
61  {
62  if (gridoperator.trialGridFunctionSpace().gridView().comm().rank()>0)
63  verbosity_level = 0;
64  else
65  verbosity_level = verbosity_level_;
66  }
67 
68  protected:
69  GridOperator& gridoperator;
70  TrialVector *u;
72  unsigned int verbosity_level;
73  RFType prev_defect;
76  RFType reduction;
77  RFType abs_limit;
78 
79  NewtonBase(GridOperator& go, TrialVector& u_)
80  : gridoperator(go)
81  , u(&u_)
82  , verbosity_level(1)
83  {
84  if (gridoperator.trialGridFunctionSpace().gridView().comm().rank()>0)
85  verbosity_level = 0;
86  }
87 
88  NewtonBase(GridOperator& go)
89  : gridoperator(go)
90  , u(0)
91  , verbosity_level(1)
92  {
93  if (gridoperator.trialGridFunctionSpace().gridView().comm().rank()>0)
94  verbosity_level = 0;
95  }
96 
97  virtual ~NewtonBase() { }
98 
99  virtual bool terminate() = 0;
100  virtual void prepare_step(Matrix& A, TestVector& r) = 0;
101  virtual void line_search(TrialVector& z, TestVector& r) = 0;
102  virtual void defect(TestVector& r) = 0;
103  };
104 
105  template<class GOS, class S, class TrlV, class TstV>
106  class NewtonSolver : public virtual NewtonBase<GOS,TrlV,TstV>
107  {
108  typedef S Solver;
109  typedef GOS GridOperator;
110  typedef TrlV TrialVector;
111  typedef TstV TestVector;
112 
113  typedef typename TestVector::ElementType RFType;
114  typedef typename GOS::Traits::Jacobian Matrix;
115 
116  public:
118 
119  NewtonSolver(GridOperator& go, TrialVector& u_, Solver& solver_)
120  : NewtonBase<GOS,TrlV,TstV>(go,u_)
121  , solver(solver_)
122  , result_valid(false)
123  {}
124 
125  NewtonSolver(GridOperator& go, Solver& solver_)
126  : NewtonBase<GOS,TrlV,TstV>(go)
127  , solver(solver_)
128  , result_valid(false)
129  {}
130 
131  void apply();
132 
133  void apply(TrialVector& u_);
134 
135  const Result& result() const
136  {
137  if (!result_valid)
138  DUNE_THROW(NewtonError,
139  "NewtonSolver::result() called before NewtonSolver::solve()");
140  return this->res;
141  }
142 
143  protected:
144  virtual void defect(TestVector& r)
145  {
146  r = 0.0; // TODO: vector interface
147  this->gridoperator.residual(*this->u, r);
148  this->res.defect = this->solver.norm(r); // TODO: solver interface
149  if (!std::isfinite(this->res.defect))
150  DUNE_THROW(NewtonDefectError,
151  "NewtonSolver::defect(): Non-linear defect is NaN or Inf");
152  }
153 
154 
155  private:
156  void linearSolve(Matrix& A, TrialVector& z, TestVector& r) const
157  {
158  if (this->verbosity_level >= 4)
159  std::cout << " Solving linear system..." << std::endl;
160  z = 0.0; // TODO: vector interface
161  this->solver.apply(A, z, r, this->linear_reduction); // TODO: solver interface
162 
163  ios_base_all_saver restorer(std::cout); // store old ios flags
164 
165  if (!this->solver.result().converged) // TODO: solver interface
166  DUNE_THROW(NewtonLinearSolverError,
167  "NewtonSolver::linearSolve(): Linear solver did not converge "
168  "in " << this->solver.result().iterations << " iterations");
169  if (this->verbosity_level >= 4)
170  std::cout << " linear solver iterations: "
171  << std::setw(12) << solver.result().iterations << std::endl
172  << " linear defect reduction: "
173  << std::setw(12) << std::setprecision(4) << std::scientific
174  << solver.result().reduction << std::endl;
175  }
176 
177  Solver& solver;
178  bool result_valid;
179  };
180 
181  template<class GOS, class S, class TrlV, class TstV>
183  {
184  this->u = &u_;
185  apply();
186  }
187 
188  template<class GOS, class S, class TrlV, class TstV>
190  {
191  this->res.iterations = 0;
192  this->res.converged = false;
193  this->res.reduction = 1.0;
194  this->res.conv_rate = 1.0;
195  this->res.elapsed = 0.0;
196  this->res.assembler_time = 0.0;
197  this->res.linear_solver_time = 0.0;
198  this->res.linear_solver_iterations = 0;
199  result_valid = true;
200  Timer timer;
201 
202  try
203  {
204  TestVector r(this->gridoperator.testGridFunctionSpace());
205  this->defect(r);
206  this->res.first_defect = this->res.defect;
207  this->prev_defect = this->res.defect;
208 
209  if (this->verbosity_level >= 2)
210  {
211  // store old ios flags
212  ios_base_all_saver restorer(std::cout);
213  std::cout << " Initial defect: "
214  << std::setw(12) << std::setprecision(4) << std::scientific
215  << this->res.defect << std::endl;
216  }
217 
218  Matrix A(this->gridoperator);
219  TrialVector z(this->gridoperator.trialGridFunctionSpace());
220 
221  while (!this->terminate())
222  {
223  if (this->verbosity_level >= 3)
224  std::cout << " Newton iteration " << this->res.iterations
225  << " --------------------------------" << std::endl;
226 
227  Timer assembler_timer;
228  try
229  {
230  this->prepare_step(A,r);
231  }
232  catch (...)
233  {
234  this->res.assembler_time += assembler_timer.elapsed();
235  throw;
236  }
237  double assembler_time = assembler_timer.elapsed();
238  this->res.assembler_time += assembler_time;
239  if (this->verbosity_level >= 3)
240  std::cout << " matrix assembly time: "
241  << std::setw(12) << std::setprecision(4) << std::scientific
242  << assembler_time << std::endl;
243 
244  Timer linear_solver_timer;
245  try
246  {
247  this->linearSolve(A, z, r);
248  }
249  catch (...)
250  {
251  this->res.linear_solver_time += linear_solver_timer.elapsed();
252  this->res.linear_solver_iterations += this->solver.result().iterations;
253  throw;
254  }
255  double linear_solver_time = linear_solver_timer.elapsed();
256  this->res.linear_solver_time += linear_solver_time;
257  this->res.linear_solver_iterations += this->solver.result().iterations;
258 
259  try
260  {
261  this->line_search(z, r);
262  }
263  catch (NewtonLineSearchError)
264  {
265  if (this->reassembled)
266  throw;
267  if (this->verbosity_level >= 3)
268  std::cout << " line search failed - trying again with reassembled matrix" << std::endl;
269  continue;
270  }
271 
272  this->res.reduction = this->res.defect/this->res.first_defect;
273  this->res.iterations++;
274  this->res.conv_rate = std::pow(this->res.reduction, 1.0/this->res.iterations);
275 
276  // store old ios flags
277  ios_base_all_saver restorer(std::cout);
278 
279  if (this->verbosity_level >= 3)
280  std::cout << " linear solver time: "
281  << std::setw(12) << std::setprecision(4) << std::scientific
282  << linear_solver_time << std::endl
283  << " defect reduction (this iteration):"
284  << std::setw(12) << std::setprecision(4) << std::scientific
285  << this->res.defect/this->prev_defect << std::endl
286  << " defect reduction (total): "
287  << std::setw(12) << std::setprecision(4) << std::scientific
288  << this->res.reduction << std::endl
289  << " new defect: "
290  << std::setw(12) << std::setprecision(4) << std::scientific
291  << this->res.defect << std::endl;
292  if (this->verbosity_level == 2)
293  std::cout << " Newton iteration " << std::setw(2) << this->res.iterations
294  << ". New defect: "
295  << std::setw(12) << std::setprecision(4) << std::scientific
296  << this->res.defect
297  << ". Reduction (this): "
298  << std::setw(12) << std::setprecision(4) << std::scientific
299  << this->res.defect/this->prev_defect
300  << ". Reduction (total): "
301  << std::setw(12) << std::setprecision(4) << std::scientific
302  << this->res.reduction << std::endl;
303  }
304  }
305  catch(...)
306  {
307  this->res.elapsed = timer.elapsed();
308  throw;
309  }
310  this->res.elapsed = timer.elapsed();
311 
312  ios_base_all_saver restorer(std::cout); // store old ios flags
313 
314  if (this->verbosity_level == 1)
315  std::cout << " Newton converged after " << std::setw(2) << this->res.iterations
316  << " iterations. Reduction: "
317  << std::setw(12) << std::setprecision(4) << std::scientific
318  << this->res.reduction
319  << " (" << std::setprecision(4) << this->res.elapsed << "s)"
320  << std::endl;
321  }
322 
323  template<class GOS, class TrlV, class TstV>
324  class NewtonTerminate : public virtual NewtonBase<GOS,TrlV,TstV>
325  {
326  typedef GOS GridOperator;
327  typedef TrlV TrialVector;
328 
329  typedef typename TstV::ElementType RFType;
330 
331  public:
332  NewtonTerminate(GridOperator& go, TrialVector& u_)
333  : NewtonBase<GOS,TrlV,TstV>(go,u_)
334  , maxit(40)
335  , force_iteration(false)
336  {
337  this->reduction = 1e-8;
338  this->abs_limit = 1e-12;
339  }
340 
341  NewtonTerminate(GridOperator& go)
342  : NewtonBase<GOS,TrlV,TstV>(go)
343  , maxit(40)
344  , force_iteration(false)
345  {
346  this->reduction = 1e-8;
347  this->abs_limit = 1e-12;
348  }
349 
350  void setReduction(RFType reduction_)
351  {
352  this->reduction = reduction_;
353  }
354 
355  void setMaxIterations(unsigned int maxit_)
356  {
357  maxit = maxit_;
358  }
359 
360  void setForceIteration(bool force_iteration_)
361  {
362  force_iteration = force_iteration_;
363  }
364 
365  void setAbsoluteLimit(RFType abs_limit_)
366  {
367  this->abs_limit = abs_limit_;
368  }
369 
370  virtual bool terminate()
371  {
372  if (force_iteration && this->res.iterations == 0)
373  return false;
374  this->res.converged = this->res.defect < this->abs_limit
375  || this->res.defect < this->res.first_defect * this->reduction;
376  if (this->res.iterations >= maxit && !this->res.converged)
377  DUNE_THROW(NewtonNotConverged,
378  "NewtonTerminate::terminate(): Maximum iteration count reached");
379  return this->res.converged;
380  }
381 
382  private:
383  unsigned int maxit;
384  bool force_iteration;
385  };
386 
387  template<class GOS, class TrlV, class TstV>
388  class NewtonPrepareStep : public virtual NewtonBase<GOS,TrlV,TstV>
389  {
390  typedef GOS GridOperator;
391  typedef TrlV TrialVector;
392 
393  typedef typename TstV::ElementType RFType;
394  typedef typename GOS::Traits::Jacobian Matrix;
395 
396  public:
397  NewtonPrepareStep(GridOperator& go, TrialVector& u_)
398  : NewtonBase<GOS,TrlV,TstV>(go,u_)
399  , min_linear_reduction(1e-3)
400  , fixed_linear_reduction(0.0)
401  , reassemble_threshold(0.0)
402  {}
403 
404  NewtonPrepareStep(GridOperator& go)
405  : NewtonBase<GOS,TrlV,TstV>(go)
406  , min_linear_reduction(1e-3)
407  , fixed_linear_reduction(0.0)
408  , reassemble_threshold(0.0)
409  {}
410 
411  /* with min_linear_reduction > 0, the linear reduction will be
412  determined as mininum of the min_linear_reduction and the
413  linear_reduction needed to achieve second order
414  Newton convergence. */
415  void setMinLinearReduction(RFType min_linear_reduction_)
416  {
417  min_linear_reduction = min_linear_reduction_;
418  }
419 
420  /* with fixed_linear_reduction > 0, the linear reduction
421  rate will always be fixed to min_linear_reduction. */
422  void setFixedLinearReduction(bool fixed_linear_reduction_)
423  {
424  fixed_linear_reduction = fixed_linear_reduction_;
425  }
426 
427  void setReassembleThreshold(RFType reassemble_threshold_)
428  {
429  reassemble_threshold = reassemble_threshold_;
430  }
431 
432  virtual void prepare_step(Matrix& A, TstV& )
433  {
434  this->reassembled = false;
435  if (this->res.defect/this->prev_defect > reassemble_threshold)
436  {
437  if (this->verbosity_level >= 3)
438  std::cout << " Reassembling matrix..." << std::endl;
439  A = 0.0; // TODO: Matrix interface
440  this->gridoperator.jacobian(*this->u, A);
441  this->reassembled = true;
442  }
443 
444  if (fixed_linear_reduction == true)
445  this->linear_reduction = min_linear_reduction;
446  else {
447  // determine maximum defect, where Newton is converged.
448  RFType stop_defect =
449  std::max(this->res.first_defect * this->reduction,
450  this->abs_limit);
451 
452  /*
453  To achieve second order convergence of newton
454  we need a linear reduction of at least
455  current_defect^2/prev_defect^2.
456  For the last newton step a linear reduction of
457  1/10*end_defect/current_defect
458  is sufficient for convergence.
459  */
460  if ( stop_defect/(10*this->res.defect) >
461  this->res.defect*this->res.defect/(this->prev_defect*this->prev_defect) )
462  this->linear_reduction =
463  stop_defect/(10*this->res.defect);
464  else
465  this->linear_reduction =
466  std::min(min_linear_reduction,this->res.defect*this->res.defect/(this->prev_defect*this->prev_defect));
467  }
468 
469  this->prev_defect = this->res.defect;
470 
471  ios_base_all_saver restorer(std::cout); // store old ios flags
472 
473  if (this->verbosity_level >= 3)
474  std::cout << " requested linear reduction: "
475  << std::setw(12) << std::setprecision(4) << std::scientific
476  << this->linear_reduction << std::endl;
477  }
478 
479  private:
480  RFType min_linear_reduction;
481  bool fixed_linear_reduction;
482  RFType reassemble_threshold;
483  };
484 
485  template<class GOS, class TrlV, class TstV>
486  class NewtonLineSearch : public virtual NewtonBase<GOS,TrlV,TstV>
487  {
488  typedef GOS GridOperator;
489  typedef TrlV TrialVector;
490  typedef TstV TestVector;
491 
492  typedef typename TestVector::ElementType RFType;
493 
494  public:
498 
499  NewtonLineSearch(GridOperator& go, TrialVector& u_)
500  : NewtonBase<GOS,TrlV,TstV>(go,u_)
501  , strategy(hackbuschReusken)
502  , maxit(10)
503  , damping_factor(0.5)
504  {}
505 
506  NewtonLineSearch(GridOperator& go)
507  : NewtonBase<GOS,TrlV,TstV>(go)
508  , strategy(hackbuschReusken)
509  , maxit(10)
510  , damping_factor(0.5)
511  {}
512 
514  {
515  strategy = strategy_;
516  }
517 
518  void setLineSearchStrategy(std::string strategy_)
519  {
520  strategy = strategyFromName(strategy_);
521  }
522 
523  void setLineSearchMaxIterations(unsigned int maxit_)
524  {
525  maxit = maxit_;
526  }
527 
528  void setLineSearchDampingFactor(RFType damping_factor_)
529  {
530  damping_factor = damping_factor_;
531  }
532 
533  virtual void line_search(TrialVector& z, TestVector& r)
534  {
535  if (strategy == noLineSearch)
536  {
537  this->u->axpy(-1.0, z); // TODO: vector interface
538  this->defect(r);
539  return;
540  }
541 
542  if (this->verbosity_level >= 4)
543  std::cout << " Performing line search..." << std::endl;
544  RFType lambda = 1.0;
545  RFType best_lambda = 0.0;
546  RFType best_defect = this->res.defect;
547  TrialVector prev_u(*this->u); // TODO: vector interface
548  unsigned int i = 0;
549  ios_base_all_saver restorer(std::cout); // store old ios flags
550 
551  while (1)
552  {
553  if (this->verbosity_level >= 4)
554  std::cout << " trying line search damping factor: "
555  << std::setw(12) << std::setprecision(4) << std::scientific
556  << lambda
557  << std::endl;
558 
559  this->u->axpy(-lambda, z); // TODO: vector interface
560  try {
561  this->defect(r);
562  }
563  catch (NewtonDefectError)
564  {
565  if (this->verbosity_level >= 4)
566  std::cout << " Nans detected" << std::endl;
567  } // ignore NaNs and try again with lower lambda
568 
569  if (this->res.defect <= (1.0 - lambda/4) * this->prev_defect)
570  {
571  if (this->verbosity_level >= 4)
572  std::cout << " line search converged" << std::endl;
573  break;
574  }
575 
576  if (this->res.defect < best_defect)
577  {
578  best_defect = this->res.defect;
579  best_lambda = lambda;
580  }
581 
582  if (++i >= maxit)
583  {
584  if (this->verbosity_level >= 4)
585  std::cout << " max line search iterations exceeded" << std::endl;
586  switch (strategy)
587  {
588  case hackbuschReusken:
589  *this->u = prev_u;
590  this->defect(r);
591  DUNE_THROW(NewtonLineSearchError,
592  "NewtonLineSearch::line_search(): line search failed, "
593  "max iteration count reached, "
594  "defect did not improve enough");
596  if (best_lambda == 0.0)
597  {
598  *this->u = prev_u;
599  this->defect(r);
600  DUNE_THROW(NewtonLineSearchError,
601  "NewtonLineSearch::line_search(): line search failed, "
602  "max iteration count reached, "
603  "defect did not improve in any of the iterations");
604  }
605  if (best_lambda != lambda)
606  {
607  *this->u = prev_u;
608  this->u->axpy(-best_lambda, z);
609  this->defect(r);
610  }
611  break;
612  case noLineSearch:
613  break;
614  }
615  break;
616  }
617 
618  lambda *= damping_factor;
619  *this->u = prev_u; // TODO: vector interface
620  }
621  if (this->verbosity_level >= 4)
622  std::cout << " line search damping factor: "
623  << std::setw(12) << std::setprecision(4) << std::scientific
624  << lambda << std::endl;
625  }
626 
627  protected:
629  Strategy strategyFromName(const std::string & s) {
630  if (s == "noLineSearch")
631  return noLineSearch;
632  if (s == "hackbuschReusken")
633  return hackbuschReusken;
634  if (s == "hackbuschReuskenAcceptBest")
636  }
637 
638  private:
639  Strategy strategy;
640  unsigned int maxit;
641  RFType damping_factor;
642  };
643 
644  template<class GOS, class S, class TrlV, class TstV = TrlV>
645  class Newton : public NewtonSolver<GOS,S,TrlV,TstV>
646  , public NewtonTerminate<GOS,TrlV,TstV>
647  , public NewtonLineSearch<GOS,TrlV,TstV>
648  , public NewtonPrepareStep<GOS,TrlV,TstV>
649  {
650  typedef GOS GridOperator;
651  typedef S Solver;
652  typedef TrlV TrialVector;
653 
654  public:
655  Newton(GridOperator& go, TrialVector& u_, Solver& solver_)
656  : NewtonBase<GOS,TrlV,TstV>(go,u_)
657  , NewtonSolver<GOS,S,TrlV,TstV>(go,u_,solver_)
658  , NewtonTerminate<GOS,TrlV,TstV>(go,u_)
659  , NewtonLineSearch<GOS,TrlV,TstV>(go,u_)
660  , NewtonPrepareStep<GOS,TrlV,TstV>(go,u_)
661  {}
662  Newton(GridOperator& go, Solver& solver_)
663  : NewtonBase<GOS,TrlV,TstV>(go)
664  , NewtonSolver<GOS,S,TrlV,TstV>(go,solver_)
665  , NewtonTerminate<GOS,TrlV,TstV>(go)
666  , NewtonLineSearch<GOS,TrlV,TstV>(go)
667  , NewtonPrepareStep<GOS,TrlV,TstV>(go)
668  {}
670 
691  void setParameters(Dune::ParameterTree & param)
692  {
693  typedef typename TstV::ElementType RFType;
694  if (param.hasKey("VerbosityLevel"))
695  this->setVerbosityLevel(
696  param.get<unsigned int>("VerbosityLevel"));
697  if (param.hasKey("VerbosityLevel"))
698  this->setVerbosityLevel(
699  param.get<unsigned int>("VerbosityLevel"));
700  if (param.hasKey("Reduction"))
701  this->setReduction(
702  param.get<RFType>("Reduction"));
703  if (param.hasKey("MaxIterations"))
704  this->setMaxIterations(
705  param.get<unsigned int>("MaxIterations"));
706  if (param.hasKey("ForceIteration"))
707  this->setForceIteration(
708  param.get<bool>("ForceIteration"));
709  if (param.hasKey("AbsoluteLimit"))
710  this->setAbsoluteLimit(
711  param.get<RFType>("AbsoluteLimit"));
712  if (param.hasKey("MinLinearReduction"))
713  this->setMinLinearReduction(
714  param.get<RFType>("MinLinearReduction"));
715  if (param.hasKey("FixedLinearReduction"))
717  param.get<bool>("FixedLinearReduction"));
718  if (param.hasKey("ReassembleThreshold"))
720  param.get<RFType>("ReassembleThreshold"));
721  if (param.hasKey("LineSearchStrategy"))
722  this->setLineSearchStrategy(
723  param.get<std::string>("LineSearchStrategy"));
724  if (param.hasKey("LineSearchMaxIterations"))
726  param.get<unsigned int>("LineSearchMaxIterations"));
727  if (param.hasKey("LineSearchDampingFactor"))
729  param.get<RFType>("LineSearchDampingFactor"));
730  }
731  };
732  }
733 }
734 
735 #endif // DUNE_PDELAB_NEWTON_HH
Definition: newton.hh:106
RFType defect
Definition: newton.hh:35
void setMinLinearReduction(RFType min_linear_reduction_)
Definition: newton.hh:415
void setAbsoluteLimit(RFType abs_limit_)
Definition: newton.hh:365
Definition: newton.hh:46
virtual void prepare_step(Matrix &A, TestVector &r)=0
Strategy strategyFromName(const std::string &s)
Definition: newton.hh:629
NewtonSolver(GridOperator &go, Solver &solver_)
Definition: newton.hh:125
bool converged
Definition: solver.hh:32
RFType first_defect
Definition: newton.hh:34
NewtonResult()
Definition: newton.hh:40
void setLineSearchMaxIterations(unsigned int maxit_)
Definition: newton.hh:523
virtual ~NewtonBase()
Definition: newton.hh:97
virtual void line_search(TrialVector &z, TestVector &r)=0
Definition: newton.hh:32
void setLineSearchStrategy(std::string strategy_)
Definition: newton.hh:518
Newton(GridOperator &go, Solver &solver_)
Definition: newton.hh:662
RFType reduction
Definition: newton.hh:76
void setReassembleThreshold(RFType reassemble_threshold_)
Definition: newton.hh:427
void setLineSearchStrategy(Strategy strategy_)
Definition: newton.hh:513
virtual bool terminate()=0
NewtonTerminate(GridOperator &go, TrialVector &u_)
Definition: newton.hh:332
NewtonLineSearch(GridOperator &go, TrialVector &u_)
Definition: newton.hh:499
NewtonBase(GridOperator &go)
Definition: newton.hh:88
TrialVector * u
Definition: newton.hh:70
RFType abs_limit
Definition: newton.hh:77
Definition: solver.hh:30
virtual void defect(TestVector &r)
Definition: newton.hh:144
NewtonResult< RFType > Result
Definition: newton.hh:58
void setVerbosityLevel(unsigned int verbosity_level_)
Definition: newton.hh:60
bool reassembled
Definition: newton.hh:75
Newton(GridOperator &go, TrialVector &u_, Solver &solver_)
Definition: newton.hh:655
NewtonPrepareStep(GridOperator &go)
Definition: newton.hh:404
virtual void prepare_step(Matrix &A, TstV &)
Definition: newton.hh:432
Definition: newton.hh:645
void setReduction(RFType reduction_)
Definition: newton.hh:350
virtual void defect(TestVector &r)=0
NewtonResult< RFType > Result
Definition: newton.hh:117
NewtonLineSearch(GridOperator &go)
Definition: newton.hh:506
Definition: newton.hh:24
double assembler_time
Definition: newton.hh:36
Definition: newton.hh:388
int linear_solver_iterations
Definition: newton.hh:38
Definition: newton.hh:27
double linear_solver_time
Definition: newton.hh:37
GridOperator & gridoperator
Definition: newton.hh:69
Definition: newton.hh:324
Result res
Definition: newton.hh:71
Strategy
Definition: newton.hh:495
NewtonPrepareStep(GridOperator &go, TrialVector &u_)
Definition: newton.hh:397
NewtonTerminate(GridOperator &go)
Definition: newton.hh:341
virtual bool terminate()
Definition: newton.hh:370
Definition: newton.hh:28
void setParameters(Dune::ParameterTree &param)
interpret a parameter tree as a set of options for the newton solver
Definition: newton.hh:691
Definition: newton.hh:486
NewtonBase(GridOperator &go, TrialVector &u_)
Definition: newton.hh:79
Base class for all PDELab exceptions.
Definition: exceptions.hh:17
void setLineSearchDampingFactor(RFType damping_factor_)
Definition: newton.hh:528
Definition: newton.hh:25
const E & e
Definition: interpolate.hh:172
void setForceIteration(bool force_iteration_)
Definition: newton.hh:360
void setFixedLinearReduction(bool fixed_linear_reduction_)
Definition: newton.hh:422
virtual void line_search(TrialVector &z, TestVector &r)
Definition: newton.hh:533
NewtonSolver(GridOperator &go, TrialVector &u_, Solver &solver_)
Definition: newton.hh:119
unsigned int verbosity_level
Definition: newton.hh:72
void apply()
Definition: newton.hh:189
void setMaxIterations(unsigned int maxit_)
Definition: newton.hh:355
RFType prev_defect
Definition: newton.hh:73
const std::string s
Definition: function.hh:1103
const Result & result() const
Definition: newton.hh:135
RFType linear_reduction
Definition: newton.hh:74
unsigned int iterations
Definition: solver.hh:33