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