NOX  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Public Member Functions | Protected Types | Protected Member Functions | Protected Attributes | List of all members
NOX::Solver::TrustRegionBased Class Reference

Newton-like solver using a trust region. More...

#include <NOX_Solver_TrustRegionBased.H>

Inheritance diagram for NOX::Solver::TrustRegionBased:
Inheritance graph
[legend]
Collaboration diagram for NOX::Solver::TrustRegionBased:
Collaboration graph
[legend]

Public Member Functions

 TrustRegionBased (const Teuchos::RCP< NOX::Abstract::Group > &grp, const Teuchos::RCP< NOX::StatusTest::Generic > &tests, const Teuchos::RCP< Teuchos::ParameterList > &params)
 Constructor. More...
 
virtual ~TrustRegionBased ()
 Destructor.
 
virtual void reset (const NOX::Abstract::Vector &initialGuess, const Teuchos::RCP< NOX::StatusTest::Generic > &tests)
 Resets the solver, sets a new status test, and sets a new initial guess.
 
virtual void reset (const NOX::Abstract::Vector &initialGuess)
 Resets the solver and sets a new initial guess.
 
virtual void reset ()
 Resets the solver for another solve. This resets the counters and status only. Uses the final solution from the last solve as the initial guess for the next solve. More...
 
virtual NOX::StatusTest::StatusType step ()
 Do one nonlinear step in the iteration sequence and return status.
 
virtual NOX::StatusTest::StatusType solve ()
 Solve the nonlinear problem and return final status. More...
 
virtual const
NOX::Abstract::Group
getSolutionGroup () const
 Return a reference to the current solution group.
 
virtual const
NOX::Abstract::Group
getPreviousSolutionGroup () const
 Return a reference to the previous solution group.
 
virtual NOX::StatusTest::StatusType getStatus () const
 Returns the current status of the solver.
 
virtual int getNumIterations () const
 Get number of iterations.
 
virtual const
Teuchos::ParameterList
getList () const
 Return a reference to the solver parameters.
 
virtual Teuchos::RCP< const
NOX::Abstract::Group
getSolutionGroupPtr () const
 Return a RCP to the solution group.
 
virtual Teuchos::RCP< const
NOX::Abstract::Group
getPreviousSolutionGroupPtr () const
 Return a RCP to the previous solution group.
 
virtual Teuchos::RCP< const
Teuchos::ParameterList
getListPtr () const
 Return a RCP to the solver parameters.
 
virtual Teuchos::RCP< const
NOX::SolverStats
getSolverStatistics () const
 Return a RCP to the solver statistics.
 
- Public Member Functions inherited from NOX::Solver::Generic
 Generic ()
 Constructor (does nothing)
 
virtual ~Generic ()
 Destructor (does nothing)
 

Protected Types

enum  StepType { Newton, Cauchy, Dogleg }
 Enumerated list for each direction that may be required in the Trust region computation. More...
 

Protected Member Functions

virtual void init ()
 Print out initialization information and calcuation the RHS.
 
virtual void invalid (const std::string &param, double value) const
 Print and error message and throw and error.
 
virtual void printUpdate ()
 Prints the current iteration information.
 

Protected Attributes

Teuchos::RCP< NOX::GlobalDataglobalDataPtr
 Pointer to the global data object.
 
Teuchos::RCP< NOX::UtilsutilsPtr
 Printing Utils.
 
Teuchos::RCP
< NOX::Abstract::Group
solnPtr
 Current solution.
 
Teuchos::RCP
< NOX::Abstract::Group
oldSolnPtr
 Previous solution pointer. More...
 
Teuchos::RCP
< NOX::Abstract::Vector
newtonVecPtr
 Current search direction.pointer. More...
 
Teuchos::RCP
< NOX::Abstract::Vector
cauchyVecPtr
 Current search direction.pointer. More...
 
Teuchos::RCP
< NOX::Abstract::Vector
aVecPtr
 Extra vector used in computations. More...
 
Teuchos::RCP
< NOX::Abstract::Vector
bVecPtr
 Extra vector used in computations. More...
 
Teuchos::RCP
< NOX::StatusTest::Generic
testPtr
 Stopping test.
 
NOX::StatusTest::CheckType checkType
 Type of check to use for status tests. See NOX::StatusTest for more details.
 
Teuchos::RCP
< Teuchos::ParameterList
paramsPtr
 Input parameters.
 
Teuchos::RCP
< NOX::Direction::Generic
newtonPtr
 Newton Search Direction.
 
Teuchos::RCP
< NOX::Direction::Generic
cauchyPtr
 Cauchy Search Direction.
 
double radius
 Radius of the trust region.
 
double minRatio
 Minimum improvement ratio to accept step.
 
double initRadius
 Initial trust region radius.
 
double minRadius
 Minimum trust region radius.
 
double maxRadius
 Maximum trust region radius.
 
double contractTriggerRatio
 ratio < alpha triggers contraction
 
double expandTriggerRatio
 ratio > beta triggers expansion
 
double expandFactor
 Expansion factor.
 
double contractFactor
 Constraction factor.
 
double recoveryStep
 
double newF
 Value of $ f $ at current solution.
 
double oldF
 Value of $ f $ at previous solution.
 
double dx
 norm(xnew - xold)
 
int nIter
 Number of nonlinear iterations.
 
NOX::StatusTest::StatusType status
 Status of nonlinear solver.
 
StepType stepType
 Type of step to be taken.
 
Teuchos::RCP
< NOX::MeritFunction::Generic
meritFuncPtr
 Stores a user supplied merit function if supplied in the parameter list.
 
bool useAredPredRatio
 If set to true, the minimum improvement ratio condition uses an Ared/Pred approach.
 
Teuchos::RCP< NOX::Observerobserver
 Pointer to a user defined NOX::Observer object.
 

Detailed Description

Newton-like solver using a trust region.

Our goal is to solve: $ F(x) = 0, $ where $ F:\Re^n \rightarrow \Re^n $. Alternatively, we might say that we wish to solve

   $ \min f(x) \equiv \frac{1}{2} \|F(x)\|^2_2. $

The trust region subproblem (TRSP) at iteration $ k$ is given by

   $ \min \; m_k(s) \equiv f_k + g_k^T d + \frac{1}{2} d^T B_k d, \mbox{ s.t. } \|d\| \leq \Delta_k \quad \mbox{(TRSP)} $

where

The "improvement ratio" for a given step $ s $ is defined as

   $ \rho = \displaystyle\frac{ f(x_k) - f(x_k + d) } { m_k(0) - m_k(d) } $

An iteration consists of the following steps.

"Minimum Trust Region Radius" ( $\Delta_{\min}$) - Minimum allowable trust region radius. Defaults to 1.0e-6.

Output Paramters

A sublist for output parameters called "Output" will be created and contain the following parameters:

Member Enumeration Documentation

Enumerated list for each direction that may be required in the Trust region computation.

Enumerator
Newton 

Use the Newton direction.

Cauchy 

Use the Cauchy direction.

Dogleg 

Use the doglog direction.

Constructor & Destructor Documentation

TrustRegionBased::TrustRegionBased ( const Teuchos::RCP< NOX::Abstract::Group > &  grp,
const Teuchos::RCP< NOX::StatusTest::Generic > &  tests,
const Teuchos::RCP< Teuchos::ParameterList > &  params 
)

Member Function Documentation

void TrustRegionBased::reset ( )
virtual

Resets the solver for another solve. This resets the counters and status only. Uses the final solution from the last solve as the initial guess for the next solve.

NOTE: All NOX solvers will call reset() automatically at teh beginning of the solve() method. We add the reset() method to the solver interface for the application to call in case the application needs to reset counters and status manually before the next call to solve() is made.

Implements NOX::Solver::Generic.

References dx, nIter, status, and NOX::StatusTest::Unconverged.

Referenced by solve().

NOX::StatusTest::StatusType TrustRegionBased::solve ( )
virtual

Member Data Documentation

Teuchos::RCP<NOX::Abstract::Vector> NOX::Solver::TrustRegionBased::aVecPtr
protected

Extra vector used in computations.

We have both a pointer and a reference because we need to create a DERIVED object and then want to have a reference to it.

Referenced by step().

Teuchos::RCP<NOX::Abstract::Vector> NOX::Solver::TrustRegionBased::bVecPtr
protected

Extra vector used in computations.

We have both a pointer and a reference because we need to create a DERIVED object and then want to have a reference to it.

Referenced by step().

Teuchos::RCP<NOX::Abstract::Vector> NOX::Solver::TrustRegionBased::cauchyVecPtr
protected

Current search direction.pointer.

We have both a pointer and a reference because we need to create a DERIVED object and then want to have a reference to it.

Referenced by step().

Teuchos::RCP<NOX::Abstract::Vector> NOX::Solver::TrustRegionBased::newtonVecPtr
protected

Current search direction.pointer.

We have both a pointer and a reference because we need to create a DERIVED object and then want to have a reference to it.

Referenced by step().

Teuchos::RCP<NOX::Abstract::Group> NOX::Solver::TrustRegionBased::oldSolnPtr
protected

Previous solution pointer.

We have both a pointer and a reference because we need to create a DERIVED object and then want to have a reference to it.

Referenced by getPreviousSolutionGroup(), getPreviousSolutionGroupPtr(), and step().

double NOX::Solver::TrustRegionBased::recoveryStep
protected

Take a step of this length in the Newton direction if the trust-region search fails

Referenced by init(), and step().


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