dune-pdelab 2.7-git
Loading...
Searching...
No Matches
Public Types | Public Member Functions | List of all members
Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ > Class Template Reference

Newton solver for solving non-linear problems. More...

#include <dune/pdelab/solver/newton.hh>

Public Types

using GridOperator = GridOperator_
 Type of the grid operator.
 
using LinearSolver = LinearSolver_
 Type of the linear solver.
 
using Domain = typename GridOperator::Traits::Domain
 Type of the domain (solution)
 
using Range = typename GridOperator::Traits::Range
 Type of the range (residual)
 
using Jacobian = typename GridOperator::Traits::Jacobian
 Type of the Jacobian matrix.
 
using Real = typename Dune::FieldTraits< typename Domain::ElementType >::real_type
 Number type.
 
using Result = PDESolverResult< Real >
 Type of results.
 

Public Member Functions

const Resultresult () const
 Return results.
 
virtual void prepareStep (Domain &solution)
 
virtual void linearSolve ()
 
virtual void apply (Domain &solution)
 Solve the nonlinear problem using solution as initial guess and for storing the result.
 
virtual void updateDefect (Domain &solution)
 Update _residual and defect in _result.
 
void setVerbosityLevel (unsigned int verbosity)
 Set how much output you get.
 
unsigned int getVerbosityLevel () const
 Get verbosity level.
 
void setReduction (Real reduction)
 Set reduction Newton needs to achieve.
 
Real getReduction () const
 Get reduction.
 
void setAbsoluteLimit (Real absoluteLimit)
 Set absolute convergence limit.
 
Real getAbsoluteLimit () const
 
void setKeepMatrix (bool b)
 Set whether the jacobian matrix should be kept across calls to apply().
 
void setUseMaxNorm (bool b)
 Set whether to use the maximum norm for stopping criteria.
 
void setHangingNodeModifications (bool b)
 Does the problem have hanging nodes.
 
bool keepMatrix () const
 Return whether the jacobian matrix is kept across calls to apply().
 
void discardMatrix ()
 Discard the stored Jacobian matrix.
 
void setMinLinearReduction (Real minLinearReduction)
 Set the minimal reduction in the linear solver.
 
void setFixedLinearReduction (bool fixedLinearReduction)
 Set wether to use a fixed reduction in the linear solver.
 
void setReassembleThreshold (Real reassembleThreshold)
 Set a threshold, when the linear operator is reassembled.
 
void setParameters (const ParameterTree &parameterTree)
 Interpret a parameter tree as a set of options for the newton solver.
 
void setTerminate (std::shared_ptr< TerminateInterface > terminate)
 Set the termination criterion.
 
void setLineSearch (std::shared_ptr< LineSearchInterface< Domain > > lineSearch)
 Set the line search.
 
 NewtonMethod (const GridOperator &gridOperator, LinearSolver &linearSolver, const std::string &lineSearchStrategy="hackbusch_reusken")
 Construct Newton using default parameters.
 
 NewtonMethod (const GridOperator &gridOperator, LinearSolver &linearSolver, const ParameterTree &parameterTree, const std::string &lineSearchStrategy="hackbusch_reusken")
 Construct Newton passing a parameter tree.
 

Detailed Description

template<typename GridOperator_, typename LinearSolver_>
class Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >

Newton solver for solving non-linear problems.

Template Parameters
GridOperator_Grid operator for evaluation of resdidual and Jacobian
LinearSolver_Solver backend for solving linear system of equations
Examples
recipe-operator-splitting.cc.

Member Typedef Documentation

◆ Domain

template<typename GridOperator_ , typename LinearSolver_ >
using Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::Domain = typename GridOperator::Traits::Domain

Type of the domain (solution)

◆ GridOperator

template<typename GridOperator_ , typename LinearSolver_ >
using Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::GridOperator = GridOperator_

Type of the grid operator.

◆ Jacobian

template<typename GridOperator_ , typename LinearSolver_ >
using Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::Jacobian = typename GridOperator::Traits::Jacobian

Type of the Jacobian matrix.

◆ LinearSolver

template<typename GridOperator_ , typename LinearSolver_ >
using Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::LinearSolver = LinearSolver_

Type of the linear solver.

◆ Range

template<typename GridOperator_ , typename LinearSolver_ >
using Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::Range = typename GridOperator::Traits::Range

Type of the range (residual)

◆ Real

template<typename GridOperator_ , typename LinearSolver_ >
using Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::Real = typename Dune::FieldTraits<typename Domain::ElementType>::real_type

Number type.

◆ Result

template<typename GridOperator_ , typename LinearSolver_ >
using Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::Result = PDESolverResult<Real>

Type of results.

Constructor & Destructor Documentation

◆ NewtonMethod() [1/2]

template<typename GridOperator_ , typename LinearSolver_ >
Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::NewtonMethod ( const GridOperator gridOperator,
LinearSolver linearSolver,
const std::string &  lineSearchStrategy = "hackbusch_reusken" 
)
inline

Construct Newton using default parameters.

◆ NewtonMethod() [2/2]

template<typename GridOperator_ , typename LinearSolver_ >
Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::NewtonMethod ( const GridOperator gridOperator,
LinearSolver linearSolver,
const ParameterTree &  parameterTree,
const std::string &  lineSearchStrategy = "hackbusch_reusken" 
)
inline

Construct Newton passing a parameter tree.

Member Function Documentation

◆ apply()

template<typename GridOperator_ , typename LinearSolver_ >
virtual void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::apply ( Domain solution)
inlinevirtual

Solve the nonlinear problem using solution as initial guess and for storing the result.

◆ discardMatrix()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::discardMatrix ( )
inline

Discard the stored Jacobian matrix.

◆ getAbsoluteLimit()

template<typename GridOperator_ , typename LinearSolver_ >
Real Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::getAbsoluteLimit ( ) const
inline

◆ getReduction()

template<typename GridOperator_ , typename LinearSolver_ >
Real Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::getReduction ( ) const
inline

Get reduction.

◆ getVerbosityLevel()

template<typename GridOperator_ , typename LinearSolver_ >
unsigned int Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::getVerbosityLevel ( ) const
inline

Get verbosity level.

◆ keepMatrix()

template<typename GridOperator_ , typename LinearSolver_ >
bool Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::keepMatrix ( ) const
inline

Return whether the jacobian matrix is kept across calls to apply().

◆ linearSolve()

template<typename GridOperator_ , typename LinearSolver_ >
virtual void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::linearSolve ( )
inlinevirtual

◆ prepareStep()

template<typename GridOperator_ , typename LinearSolver_ >
virtual void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::prepareStep ( Domain solution)
inlinevirtual

◆ result()

template<typename GridOperator_ , typename LinearSolver_ >
const Result & Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::result ( ) const
inline

Return results.

◆ setAbsoluteLimit()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setAbsoluteLimit ( Real  absoluteLimit)
inline

Set absolute convergence limit.

◆ setFixedLinearReduction()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setFixedLinearReduction ( bool  fixedLinearReduction)
inline

Set wether to use a fixed reduction in the linear solver.

Note
If fixedLinearReduction is true, the linear reduction rate will always be fixed to minLinearReduction.

◆ setHangingNodeModifications()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setHangingNodeModifications ( bool  b)
inline

Does the problem have hanging nodes.

◆ setKeepMatrix()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setKeepMatrix ( bool  b)
inline

Set whether the jacobian matrix should be kept across calls to apply().

◆ setLineSearch()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setLineSearch ( std::shared_ptr< LineSearchInterface< Domain > >  lineSearch)
inline

Set the line search.

See getLineSearch() for already implemented line searches

◆ setMinLinearReduction()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setMinLinearReduction ( Real  minLinearReduction)
inline

Set the minimal reduction in the linear solver.

Note
with minLinearReduction > 0, the linear reduction will be determined as mininum of the minLinearReduction and the linear reduction needed to achieve second order Newton convergence. (As long as you are not using a fixed linear reduction)

◆ setParameters()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setParameters ( const ParameterTree &  parameterTree)
inline

Interpret a parameter tree as a set of options for the newton solver.

Possible parameters:

example configuration:

[newton_parameters]
reassemble_threshold = 0.1
absolute_limit = 1e-6
reduction = 1e-4
min_linear_reduction = 1e-3
[newton_parameters.terminate]
max_iterations = 15
[newton_parameters.line_search]
line_search_damping_factor = 0.7

and invocation in the code:

newton.setParameters(param.sub("NewtonParameters"));

This can also be used to set single parameters like this

Dune::ParameterTree ptree;
ptree["verbosity"] = "4";
newton.setParameters(ptree);

◆ setReassembleThreshold()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setReassembleThreshold ( Real  reassembleThreshold)
inline

Set a threshold, when the linear operator is reassembled.

We allow to keep the linear operator over several newton iterations. If the reduction in the newton drops below a given threshold the linear operator is reassembled to ensure convergence.

◆ setReduction()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setReduction ( Real  reduction)
inline

Set reduction Newton needs to achieve.

◆ setTerminate()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setTerminate ( std::shared_ptr< TerminateInterface terminate)
inline

Set the termination criterion.

◆ setUseMaxNorm()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setUseMaxNorm ( bool  b)
inline

Set whether to use the maximum norm for stopping criteria.

◆ setVerbosityLevel()

template<typename GridOperator_ , typename LinearSolver_ >
void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::setVerbosityLevel ( unsigned int  verbosity)
inline

Set how much output you get.

◆ updateDefect()

template<typename GridOperator_ , typename LinearSolver_ >
virtual void Dune::PDELab::NewtonMethod< GridOperator_, LinearSolver_ >::updateDefect ( Domain solution)
inlinevirtual

Update _residual and defect in _result.


The documentation for this class was generated from the following file: