newtonmethod.hh
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
27 #ifndef EWOMS_NEWTON_METHOD_HH
28 #define EWOMS_NEWTON_METHOD_HH
29 
30 #include "nullconvergencewriter.hh"
31 
34 #include <ewoms/common/timer.hh>
36 
37 #include <opm/material/densead/Math.hpp>
38 
39 #include <opm/common/Unused.hpp>
40 #include <opm/common/Exceptions.hpp>
41 #include <opm/common/ErrorMacros.hpp>
42 
43 #include <dune/istl/istlexception.hh>
44 #include <dune/common/classname.hh>
45 #include <dune/common/version.hh>
46 #include <dune/common/parallel/mpihelper.hh>
47 
48 #include <iostream>
49 #include <sstream>
50 
51 #include <unistd.h>
52 
53 namespace Ewoms {
54 // forward declaration of classes
55 template <class TypeTag>
57 }
58 
59 namespace Ewoms {
60 // forward declaration of property tags
61 namespace Properties {
65 
68 
70 NEW_PROP_TAG(Problem);
71 
73 NEW_PROP_TAG(Model);
74 
76 NEW_PROP_TAG(Scalar);
77 
80 
82 NEW_PROP_TAG(SolutionVector);
83 
85 NEW_PROP_TAG(PrimaryVariables);
86 
88 NEW_PROP_TAG(EnableConstraints);
89 
91 NEW_PROP_TAG(Constraints);
92 
94 NEW_PROP_TAG(GlobalEqVector);
95 
97 NEW_PROP_TAG(EqVector);
98 
100 NEW_PROP_TAG(Linearizer);
101 
103 NEW_PROP_TAG(JacobianMatrix);
104 
106 NEW_PROP_TAG(LinearSolverBackend);
107 
109 NEW_PROP_TAG(NewtonVerbose);
110 
112 NEW_PROP_TAG(NewtonConvergenceWriter);
113 
116 NEW_PROP_TAG(NewtonWriteConvergence);
117 
120 NEW_PROP_TAG(ConvergenceWriter);
121 
128 NEW_PROP_TAG(NewtonRawTolerance);
129 
132 NEW_PROP_TAG(NewtonMaxError);
133 
142 NEW_PROP_TAG(NewtonTargetIterations);
143 
145 NEW_PROP_TAG(NewtonMaxIterations);
146 
147 // set default values for the properties
150 SET_BOOL_PROP(NewtonMethod, NewtonWriteConvergence, false);
151 SET_BOOL_PROP(NewtonMethod, NewtonVerbose, true);
152 SET_SCALAR_PROP(NewtonMethod, NewtonRawTolerance, 1e-8);
153 // set the abortion tolerace to some very large value. if not
154 // overwritten at run-time this basically disables abortions
155 SET_SCALAR_PROP(NewtonMethod, NewtonMaxError, 1e100);
156 SET_INT_PROP(NewtonMethod, NewtonTargetIterations, 10);
157 SET_INT_PROP(NewtonMethod, NewtonMaxIterations, 18);
158 } // namespace Properties
159 } // namespace Ewoms
160 
161 namespace Ewoms {
169 template <class TypeTag>
170 class NewtonMethod
171 {
172  typedef typename GET_PROP_TYPE(TypeTag, NewtonMethod) Implementation;
173  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
174  typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
175  typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
176  typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
177 
178  typedef typename GET_PROP_TYPE(TypeTag, SolutionVector) SolutionVector;
179  typedef typename GET_PROP_TYPE(TypeTag, GlobalEqVector) GlobalEqVector;
180  typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
181  typedef typename GET_PROP_TYPE(TypeTag, Constraints) Constraints;
182  typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
183  typedef typename GET_PROP_TYPE(TypeTag, Linearizer) Linearizer;
184  typedef typename GET_PROP_TYPE(TypeTag, JacobianMatrix) JacobianMatrix;
185  typedef typename GET_PROP_TYPE(TypeTag, LinearSolverBackend) LinearSolverBackend;
186  typedef typename GET_PROP_TYPE(TypeTag, NewtonConvergenceWriter) ConvergenceWriter;
187 
188  typedef typename Dune::MPIHelper::MPICommunicator Communicator;
189  typedef Dune::CollectiveCommunication<Communicator> CollectiveCommunication;
190 
191 public:
192  NewtonMethod(Simulator& simulator)
193  : simulator_(simulator)
194  , endIterMsgStream_(std::ostringstream::out)
195  , linearSolver_(simulator)
196  , comm_(Dune::MPIHelper::getCommunicator())
197  , convergenceWriter_(asImp_())
198  {
199  lastError_ = 1e100;
200  error_ = 1e100;
201  tolerance_ = EWOMS_GET_PARAM(TypeTag, Scalar, NewtonRawTolerance);
202 
203  numIterations_ = 0;
204  }
205 
209  static void registerParameters()
210  {
211  LinearSolverBackend::registerParameters();
212 
213  EWOMS_REGISTER_PARAM(TypeTag, bool, NewtonVerbose,
214  "Specify whether the Newton method should inform "
215  "the user about its progress or not");
216  EWOMS_REGISTER_PARAM(TypeTag, bool, NewtonWriteConvergence,
217  "Write the convergence behaviour of the Newton "
218  "method to a VTK file");
219  EWOMS_REGISTER_PARAM(TypeTag, int, NewtonTargetIterations,
220  "The 'optimum' number of Newton iterations per "
221  "time step");
222  EWOMS_REGISTER_PARAM(TypeTag, int, NewtonMaxIterations,
223  "The maximum number of Newton iterations per time "
224  "step");
225  EWOMS_REGISTER_PARAM(TypeTag, Scalar, NewtonRawTolerance,
226  "The maximum raw error tolerated by the Newton"
227  "method for considering a solution to be "
228  "converged");
229  EWOMS_REGISTER_PARAM(TypeTag, Scalar, NewtonMaxError,
230  "The maximum error tolerated by the Newton "
231  "method to which does not cause an abort");
232  }
233 
238  bool converged() const
239  { return error_ <= tolerance(); }
240 
244  Problem& problem()
245  { return simulator_.problem(); }
246 
250  const Problem& problem() const
251  { return simulator_.problem(); }
252 
256  Model& model()
257  { return simulator_.model(); }
258 
262  const Model& model() const
263  { return simulator_.model(); }
264 
269  int numIterations() const
270  { return numIterations_; }
271 
279  void setIterationIndex(int value)
280  { numIterations_ = value; }
281 
286  Scalar tolerance() const
287  { return tolerance_; }
288 
293  void setTolerance(Scalar value)
294  { tolerance_ = value; }
295 
302  bool apply()
303  {
304  // Clear the current line using an ansi escape
305  // sequence. For an explanation see
306  // http://en.wikipedia.org/wiki/ANSI_escape_code
307  const char *clearRemainingLine = "\n";
308  if (isatty(fileno(stdout))) {
309  static const char blubb[] = { 0x1b, '[', 'K', '\r', 0 };
310  clearRemainingLine = blubb;
311  }
312 
313  // make sure all timers are prestine
314  prePostProcessTimer_.halt();
315  linearizeTimer_.halt();
316  solveTimer_.halt();
317  updateTimer_.halt();
318 
319  SolutionVector& nextSolution = model().solution(/*historyIdx=*/0);
320  SolutionVector currentSolution(nextSolution);
321  GlobalEqVector solutionUpdate(nextSolution.size());
322 
323  Linearizer& linearizer = model().linearizer();
324 
325  Ewoms::TimerGuard prePostProcessTimerGuard(prePostProcessTimer_);
326 
327  // tell the implementation that we begin solving
328  prePostProcessTimer_.start();
329  asImp_().begin_(nextSolution);
330  prePostProcessTimer_.stop();
331 
332  try {
333  Ewoms::TimerGuard innerPrePostProcessTimerGuard(prePostProcessTimer_);
334  Ewoms::TimerGuard linearizeTimerGuard(linearizeTimer_);
335  Ewoms::TimerGuard updateTimerGuard(updateTimer_);
336  Ewoms::TimerGuard solveTimerGuard(solveTimer_);
337 
338  // execute the method as long as the implementation thinks
339  // that we should do another iteration
340  while (asImp_().proceed_()) {
341  // linearize the problem at the current solution
342 
343  // notify the implementation that we're about to start
344  // a new iteration
345  prePostProcessTimer_.start();
346  asImp_().beginIteration_();
347  prePostProcessTimer_.stop();
348 
349  // make the current solution to the old one
350  currentSolution = nextSolution;
351 
352  if (asImp_().verbose_()) {
353  std::cout << "Linearize: r(x^k) = dS/dt + div F - q; M = grad r"
354  << clearRemainingLine
355  << std::flush;
356  }
357 
358  // do the actual linearization
359  linearizeTimer_.start();
360  asImp_().linearize_();
361  linearizeTimer_.stop();
362 
363  // notify the implementation of the successful linearization on order to
364  // give it the chance to update the error and thus to terminate the
365  // Newton method without the need of solving the last linearization.
366  updateTimer_.start();
367  auto& M = linearizer.matrix();
368  auto& b = linearizer.residual();
369  linearSolver_.prepareRhs(M, b);
370  asImp_().preSolve_(currentSolution, b);
371  updateTimer_.stop();
372 
373  if (!asImp_().proceed_()) {
374  if (asImp_().verbose_() && isatty(fileno(stdout)))
375  std::cout << clearRemainingLine
376  << std::flush;
377 
378  // tell the implementation that we're done with this iteration
379  prePostProcessTimer_.start();
380  asImp_().endIteration_(nextSolution, currentSolution);
381  prePostProcessTimer_.stop();
382 
383  break;
384  }
385 
386  // solve the resulting linear equation system
387  if (asImp_().verbose_()) {
388  std::cout << "Solve: M deltax^k = r"
389  << clearRemainingLine
390  << std::flush;
391  }
392 
393  solveTimer_.start();
394  solutionUpdate = 0;
395  linearSolver_.prepareMatrix(M);
396  bool converged = linearSolver_.solve(solutionUpdate);
397  solveTimer_.stop();
398 
399  if (!converged) {
400  solveTimer_.stop();
401  if (asImp_().verbose_())
402  std::cout << "Newton: Linear solver did not converge\n" << std::flush;
403 
404  prePostProcessTimer_.start();
405  asImp_().failed_();
406  prePostProcessTimer_.stop();
407 
408  return false;
409  }
410 
411  // update the solution
412  if (asImp_().verbose_()) {
413  std::cout << "Update: x^(k+1) = x^k - deltax^k"
414  << clearRemainingLine
415  << std::flush;
416  }
417 
418  // update the current solution (i.e. uOld) with the delta
419  // (i.e. u). The result is stored in u
420  updateTimer_.start();
421  asImp_().postSolve_(nextSolution,
422  currentSolution,
423  b,
424  solutionUpdate);
425  asImp_().update_(nextSolution, currentSolution, solutionUpdate, b);
426  updateTimer_.stop();
427 
428  if (asImp_().verbose_() && isatty(fileno(stdout)))
429  // make sure that the line currently holding the cursor is prestine
430  std::cout << clearRemainingLine
431  << std::flush;
432 
433  // tell the implementation that we're done with this iteration
434  prePostProcessTimer_.start();
435  asImp_().endIteration_(nextSolution, currentSolution);
436  prePostProcessTimer_.stop();
437  }
438  }
439  catch (const Dune::Exception& e)
440  {
441  if (asImp_().verbose_())
442  std::cout << "Newton method caught exception: \""
443  << e.what() << "\"\n" << std::flush;
444 
445  prePostProcessTimer_.start();
446  asImp_().failed_();
447  prePostProcessTimer_.stop();
448 
449  return false;
450  }
451  catch (const Opm::NumericalProblem& e)
452  {
453  if (asImp_().verbose_())
454  std::cout << "Newton method caught exception: \""
455  << e.what() << "\"\n" << std::flush;
456 
457  prePostProcessTimer_.start();
458  asImp_().failed_();
459  prePostProcessTimer_.stop();
460 
461  return false;
462  }
463 
464  // clear current line on terminal
465  if (asImp_().verbose_() && isatty(fileno(stdout)))
466  std::cout << clearRemainingLine
467  << std::flush;
468 
469  // tell the implementation that we're done
470  prePostProcessTimer_.start();
471  asImp_().end_();
472  prePostProcessTimer_.stop();
473 
474  // print the timing summary of the time step
475  if (asImp_().verbose_()) {
476  Scalar elapsedTot =
477  linearizeTimer_.realTimeElapsed()
478  + solveTimer_.realTimeElapsed()
479  + updateTimer_.realTimeElapsed();
480  std::cout << "Linearization/solve/update time: "
481  << linearizeTimer_.realTimeElapsed() << "("
482  << 100 * linearizeTimer_.realTimeElapsed()/elapsedTot << "%)/"
483  << solveTimer_.realTimeElapsed() << "("
484  << 100 * solveTimer_.realTimeElapsed()/elapsedTot << "%)/"
485  << updateTimer_.realTimeElapsed() << "("
486  << 100 * updateTimer_.realTimeElapsed()/elapsedTot << "%)"
487  << "\n" << std::flush;
488  }
489 
490 
491  // if we're not converged, tell the implementation that we've failed
492  if (!asImp_().converged()) {
493  prePostProcessTimer_.start();
494  asImp_().failed_();
495  prePostProcessTimer_.stop();
496  return false;
497  }
498 
499  // if we converged, tell the implementation that we've succeeded
500  prePostProcessTimer_.start();
501  asImp_().succeeded_();
502  prePostProcessTimer_.stop();
503 
504  return true;
505  }
506 
515  Scalar suggestTimeStepSize(Scalar oldTimeStep) const
516  {
517  // be aggressive reducing the time-step size but
518  // conservative when increasing it. the rationale is
519  // that we want to avoid failing in the next time
520  // integration which would be quite expensive
521  if (numIterations_ > targetIterations_()) {
522  Scalar percent = Scalar(numIterations_ - targetIterations_())
523  / targetIterations_();
524  return oldTimeStep / (1.0 + percent);
525  }
526 
527  Scalar percent = Scalar(targetIterations_() - numIterations_)
528  / targetIterations_();
529  return oldTimeStep * (1.0 + percent / 1.2);
530  }
531 
536  std::ostringstream& endIterMsg()
537  { return endIterMsgStream_; }
538 
543  void eraseMatrix()
544  { linearSolver_.eraseMatrix(); }
545 
549  LinearSolverBackend& linearSolver()
550  { return linearSolver_; }
551 
555  const LinearSolverBackend& linearSolver() const
556  { return linearSolver_; }
557 
558  const Ewoms::Timer& prePostProcessTimer() const
559  { return prePostProcessTimer_; }
560 
561  const Ewoms::Timer& linearizeTimer() const
562  { return linearizeTimer_; }
563 
564  const Ewoms::Timer& solveTimer() const
565  { return solveTimer_; }
566 
567  const Ewoms::Timer& updateTimer() const
568  { return updateTimer_; }
569 
570 protected:
574  bool verbose_() const
575  {
576  return EWOMS_GET_PARAM(TypeTag, bool, NewtonVerbose) && (comm_.rank() == 0);
577  }
578 
585  void begin_(const SolutionVector& u OPM_UNUSED)
586  {
587  numIterations_ = 0;
588 
589  if (EWOMS_GET_PARAM(TypeTag, bool, NewtonWriteConvergence))
590  convergenceWriter_.beginTimeStep();
591  }
592 
597  {
598  problem().beginIteration();
599  lastError_ = error_;
600  }
601 
605  void linearize_()
606  { model().linearizer().linearize(); }
607 
608  void preSolve_(const SolutionVector& currentSolution OPM_UNUSED,
609  const GlobalEqVector& currentResidual)
610  {
611  const auto& constraintsMap = model().linearizer().constraintsMap();
612  lastError_ = error_;
613 
614  // calculate the error as the maximum weighted tolerance of
615  // the solution's residual
616  error_ = 0;
617  for (unsigned dofIdx = 0; dofIdx < currentResidual.size(); ++dofIdx) {
618  // do not consider auxiliary DOFs for the error
619  if (dofIdx >= model().numGridDof() || model().dofTotalVolume(dofIdx) <= 0.0)
620  continue;
621 
622  // also do not consider DOFs which are constraint
623  if (enableConstraints_()) {
624  if (constraintsMap.count(dofIdx) > 0)
625  continue;
626  }
627 
628  const auto& r = currentResidual[dofIdx];
629  for (unsigned eqIdx = 0; eqIdx < r.size(); ++eqIdx)
630  error_ = Opm::max(std::abs(r[eqIdx] * model().eqWeight(dofIdx, eqIdx)), error_);
631  }
632 
633  // take the other processes into account
634  error_ = comm_.max(error_);
635 
636  // make sure that the error never grows beyond the maximum
637  // allowed one
638  if (error_ > EWOMS_GET_PARAM(TypeTag, Scalar, NewtonMaxError))
639  OPM_THROW(Opm::NumericalProblem,
640  "Newton: Error " << error_
641  << " is larger than maximum allowed error of "
642  << EWOMS_GET_PARAM(TypeTag, Scalar, NewtonMaxError));
643  }
644 
658  void postSolve_(const SolutionVector& nextSolution OPM_UNUSED,
659  const SolutionVector& currentSolution OPM_UNUSED,
660  const GlobalEqVector& currentResidual OPM_UNUSED,
661  const GlobalEqVector& solutionUpdate OPM_UNUSED)
662  { }
663 
678  void update_(SolutionVector& nextSolution,
679  const SolutionVector& currentSolution,
680  const GlobalEqVector& solutionUpdate,
681  const GlobalEqVector& currentResidual)
682  {
683  const auto& constraintsMap = model().linearizer().constraintsMap();
684 
685  // first, write out the current solution to make convergence
686  // analysis possible
687  asImp_().writeConvergence_(currentSolution, solutionUpdate);
688 
689  // make sure not to swallow non-finite values at this point
690  if (!std::isfinite(solutionUpdate.one_norm()))
691  OPM_THROW(Opm::NumericalProblem, "Non-finite update!");
692 
693  size_t numGridDof = model().numGridDof();
694  for (unsigned dofIdx = 0; dofIdx < numGridDof; ++dofIdx) {
695  if (enableConstraints_()) {
696  if (constraintsMap.count(dofIdx) > 0) {
697  const auto& constraints = constraintsMap.at(dofIdx);
698  asImp_().updateConstraintDof_(dofIdx,
699  nextSolution[dofIdx],
700  constraints);
701  }
702  else
703  asImp_().updatePrimaryVariables_(dofIdx,
704  nextSolution[dofIdx],
705  currentSolution[dofIdx],
706  solutionUpdate[dofIdx],
707  currentResidual[dofIdx]);
708  }
709  else
710  asImp_().updatePrimaryVariables_(dofIdx,
711  nextSolution[dofIdx],
712  currentSolution[dofIdx],
713  solutionUpdate[dofIdx],
714  currentResidual[dofIdx]);
715  }
716 
717  // update the DOFs of the auxiliary equations
718  size_t numDof = model().numTotalDof();
719  for (size_t dofIdx = numGridDof; dofIdx < numDof; ++dofIdx) {
720  nextSolution[dofIdx] = currentSolution[dofIdx];
721  nextSolution[dofIdx] -= solutionUpdate[dofIdx];
722  }
723  }
724 
728  void updateConstraintDof_(unsigned globalDofIdx OPM_UNUSED,
729  PrimaryVariables& nextValue,
730  const Constraints& constraints)
731  { nextValue = constraints; }
732 
736  void updatePrimaryVariables_(unsigned globalDofIdx OPM_UNUSED,
737  PrimaryVariables& nextValue,
738  const PrimaryVariables& currentValue,
739  const EqVector& update,
740  const EqVector& currentResidual OPM_UNUSED)
741  {
742  nextValue = currentValue;
743  nextValue -= update;
744  }
745 
752  void writeConvergence_(const SolutionVector& currentSolution,
753  const GlobalEqVector& solutionUpdate)
754  {
755  if (EWOMS_GET_PARAM(TypeTag, bool, NewtonWriteConvergence)) {
756  convergenceWriter_.beginIteration();
757  convergenceWriter_.writeFields(currentSolution, solutionUpdate);
758  convergenceWriter_.endIteration();
759  }
760  }
761 
768  void endIteration_(const SolutionVector& nextSolution OPM_UNUSED,
769  const SolutionVector& currentSolution OPM_UNUSED)
770  {
771  ++numIterations_;
772  problem().endIteration();
773 
774  if (asImp_().verbose_()) {
775  std::cout << "Newton iteration " << numIterations_ << ""
776  << " error: " << error_
777  << endIterMsg().str() << "\n" << std::flush;
778  }
779 
780  endIterMsgStream_.str("");
781  }
782 
786  bool proceed_() const
787  {
788  if (asImp_().numIterations() < 1)
789  return true; // we always do at least one full iteration
790  else if (asImp_().converged()) {
791  // we are below the specified tolerance, so we don't have to
792  // do more iterations
793  return false;
794  }
795  else if (asImp_().numIterations() >= asImp_().maxIterations_()) {
796  // we have exceeded the allowed number of steps. If the
797  // error was reduced by a factor of at least 4,
798  // in the last iterations we proceed even if we are above
799  // the maximum number of steps
800  return error_ * 4.0 < lastError_;
801  }
802 
803  return true;
804  }
805 
810  void end_()
811  {
812  if (EWOMS_GET_PARAM(TypeTag, bool, NewtonWriteConvergence))
813  convergenceWriter_.endTimeStep();
814  }
815 
821  void failed_()
822  { numIterations_ = targetIterations_() * 2; }
823 
829  void succeeded_()
830  {}
831 
832  // optimal number of iterations we want to achieve
833  int targetIterations_() const
834  { return EWOMS_GET_PARAM(TypeTag, int, NewtonTargetIterations); }
835  // maximum number of iterations we do before giving up
836  int maxIterations_() const
837  { return EWOMS_GET_PARAM(TypeTag, int, NewtonMaxIterations); }
838 
839  static bool enableConstraints_()
840  { return GET_PROP_VALUE(TypeTag, EnableConstraints); }
841 
842  Simulator& simulator_;
843 
844  Ewoms::Timer prePostProcessTimer_;
845  Ewoms::Timer linearizeTimer_;
846  Ewoms::Timer solveTimer_;
847  Ewoms::Timer updateTimer_;
848 
849  std::ostringstream endIterMsgStream_;
850 
851  Scalar error_;
852  Scalar lastError_;
853  Scalar tolerance_;
854 
855  // actual number of iterations done so far
856  int numIterations_;
857 
858  // the linear solver
859  LinearSolverBackend linearSolver_;
860 
861  // the collective communication used by the simulation (i.e. fake
862  // or MPI)
863  CollectiveCommunication comm_;
864 
865  // the object which writes the convergence behaviour of the Newton
866  // method to disk
867  ConvergenceWriter convergenceWriter_;
868 
869 private:
870  Implementation& asImp_()
871  { return *static_cast<Implementation *>(this); }
872  const Implementation& asImp_() const
873  { return *static_cast<const Implementation *>(this); }
874 };
875 
876 } // namespace Ewoms
877 
878 #endif
void beginIteration_()
Indicates the beginning of a Newton iteration.
Definition: newtonmethod.hh:596
A simple class which makes sure that a timer gets stopped if an exception is thrown.
Definition: timerguard.hh:40
void updatePrimaryVariables_(unsigned globalDofIdx OPM_UNUSED, PrimaryVariables &nextValue, const PrimaryVariables &currentValue, const EqVector &update, const EqVector &currentResidual OPM_UNUSED)
Update a single primary variables object.
Definition: newtonmethod.hh:736
bool proceed_() const
Returns true iff another Newton iteration should be done.
Definition: newtonmethod.hh:786
#define SET_BOOL_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant boolean value.
Definition: propertysystem.hh:361
LinearSolverBackend & linearSolver()
Returns the linear solver backend object for external use.
Definition: newtonmethod.hh:549
void start()
Start counting the time resources used by the simulation.
Definition: timer.hh:62
void setIterationIndex(int value)
Set the index of current iteration.
Definition: newtonmethod.hh:279
void begin_(const SolutionVector &u OPM_UNUSED)
Called before the Newton method is applied to an non-linear system of equations.
Definition: newtonmethod.hh:585
Definition: baseauxiliarymodule.hh:37
void succeeded_()
Called if the Newton method was successful.
Definition: newtonmethod.hh:829
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:189
std::ostringstream & endIterMsg()
Message that should be printed for the user after the end of an iteration.
Definition: newtonmethod.hh:536
Scalar tolerance() const
Return the current tolerance at which the Newton method considers itself to be converged.
Definition: newtonmethod.hh:286
void failed_()
Called if the Newton method broke down.
Definition: newtonmethod.hh:821
const Problem & problem() const
Returns a reference to the object describing the current physical problem.
Definition: newtonmethod.hh:250
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
void linearize_()
Linearize the global non-linear system of equations.
Definition: newtonmethod.hh:605
Problem & problem()
Return the object which specifies the pysical setup of the simulation.
Definition: simulator.hh:202
#define NEW_TYPE_TAG(...)
Define a new type tag.
Definition: propertysystem.hh:169
void halt()
Stop the measurement reset all timing values.
Definition: timer.hh:99
void writeConvergence_(const SolutionVector &currentSolution, const GlobalEqVector &solutionUpdate)
Write the convergence behaviour of the newton method to disk.
Definition: newtonmethod.hh:752
void end_()
Indicates that we&#39;re done solving the non-linear system of equations.
Definition: newtonmethod.hh:810
#define SET_INT_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant integer value.
Definition: propertysystem.hh:345
bool apply()
Run the Newton method.
Definition: newtonmethod.hh:302
void postSolve_(const SolutionVector &nextSolution OPM_UNUSED, const SolutionVector &currentSolution OPM_UNUSED, const GlobalEqVector &currentResidual OPM_UNUSED, const GlobalEqVector &solutionUpdate OPM_UNUSED)
Update the error of the solution given the previous iteration.
Definition: newtonmethod.hh:658
#define EWOMS_REGISTER_PARAM(TypeTag, ParamType, ParamName, Description)
Register a run-time parameter.
Definition: parametersystem.hh:68
This file provides the infrastructure to retrieve run-time parameters.
int numIterations() const
Returns the number of iterations done since the Newton method was invoked.
Definition: newtonmethod.hh:269
A convergence writer for the Newton method which does nothing.
Provides an encapsulation to measure the system time.
Definition: timer.hh:48
const LinearSolverBackend & linearSolver() const
Returns the linear solver backend object for external use.
Definition: newtonmethod.hh:555
bool converged() const
Returns true if the error of the solution is below the tolerance.
Definition: newtonmethod.hh:238
double realTimeElapsed() const
Return the real time [s] elapsed during the periods the timer was active since the last reset...
Definition: timer.hh:121
void updateConstraintDof_(unsigned globalDofIdx OPM_UNUSED, PrimaryVariables &nextValue, const Constraints &constraints)
Update the primary variables for a degree of freedom which is constraint.
Definition: newtonmethod.hh:728
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:99
The multi-dimensional Newton method.
Definition: newtonmethod.hh:56
bool verbose_() const
Returns true if the Newton method ought to be chatty.
Definition: newtonmethod.hh:574
static void registerParameters()
Register all run-time parameters for the Newton method.
Definition: newtonmethod.hh:209
Provides an encapsulation to measure the system time.
void endIteration_(const SolutionVector &nextSolution OPM_UNUSED, const SolutionVector &currentSolution OPM_UNUSED)
Indicates that one Newton iteration was finished.
Definition: newtonmethod.hh:768
Provides the magic behind the eWoms property system.
Scalar suggestTimeStepSize(Scalar oldTimeStep) const
Suggest a new time-step size based on the old time-step size.
Definition: newtonmethod.hh:515
A simple class which makes sure that a timer gets stopped if an exception is thrown.
void update_(SolutionVector &nextSolution, const SolutionVector &currentSolution, const GlobalEqVector &solutionUpdate, const GlobalEqVector &currentResidual)
Update the current solution with a delta vector.
Definition: newtonmethod.hh:678
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:75
#define NEW_PROP_TAG(PTagName)
Define a property tag.
Definition: propertysystem.hh:247
Problem & problem()
Returns a reference to the object describing the current physical problem.
Definition: newtonmethod.hh:244
const Model & model() const
Returns a reference to the numeric model.
Definition: newtonmethod.hh:262
#define SET_TYPE_PROP(EffTypeTagName, PropTagName,...)
Set a property which defines a type.
Definition: propertysystem.hh:377
Model & model()
Returns a reference to the numeric model.
Definition: newtonmethod.hh:256
void eraseMatrix()
Causes the solve() method to discared the structure of the linear system of equations the next time i...
Definition: newtonmethod.hh:543
A convergence writer for the Newton method which does nothing.
Definition: nullconvergencewriter.hh:51
void setTolerance(Scalar value)
Set the current tolerance at which the Newton method considers itself to be converged.
Definition: newtonmethod.hh:293
double stop()
Stop counting the time resources.
Definition: timer.hh:73
#define SET_SCALAR_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant scalar value.
Definition: propertysystem.hh:394