NOX Development
|
Calculates a search direction using the Nonlinear Conjugate Gradient method. More...
#include <NOX_Direction_NonlinearCG.H>
Public Member Functions | |
NonlinearCG (const Teuchos::RCP< NOX::GlobalData > &gd, Teuchos::ParameterList ¶ms) | |
Constructor. | |
virtual | ~NonlinearCG () |
Destructor. | |
virtual bool | reset (const Teuchos::RCP< NOX::GlobalData > &gd, Teuchos::ParameterList &p) |
derived | |
virtual bool | compute (Abstract::Vector &dir, Abstract::Group &grp, const Solver::Generic &solver) |
derived | |
virtual bool | compute (NOX::Abstract::Vector &dir, NOX::Abstract::Group &grp, const NOX::Solver::LineSearchBased &solver) |
Same as compute(NOX::Abstract::Vector&,
NOX::Abstract::Group&, const NOX::Solver::Generic&) | |
![]() | |
Generic () | |
Constructor. | |
virtual | ~Generic () |
Destructor. | |
Calculates a search direction using the Nonlinear Conjugate Gradient method.
Calculates the direction
where
This method provides a generalization of Linear CG to nonlinear problems. It does this by computing a search direction using an expression analogous to that of Linear CG. The negative of the current residual vector,
To use this direction, specify that the "Method" is "NonlinearCG" in the "Direction" sublist of the parameters that are passed to the solver (see NOX::Direction::Manager for more information on choosing the search direction).
The following options may be specified in the "Nonlinear CG" sublist of the "Direction" sublist of the solver parameters.
"Orthogonalize" can be either of:
"Fletcher-Reeves" [default] -
"Polak-Ribiere" -
where
These comprise the two most popular choices for orthogonalization. Both reduce to the linear result for linear problems. "Polak-Ribiere" provides an implied restart capability by setting
"Precondition" can be either "On" or "Off" [default]: determines whether or not to compute and apply preconditioner
"Restart Frequency" - An integer specification of the number of nonlinear iterations between restarts [default = 10]. Restart corresponds to setting
References
information about both linear and nonlinear conjugate gradient methods can be found in Chapter 5 of:
Nocedal & Wright, "Numerical Optimization", Springer-Verlag, New York, 1999. Though presented within the context of nonlinear optimization, the connection to nonlinear systems of equations is made by the correspondence
Another useful useful reference is:
|
virtual |
derived
Implements NOX::Direction::Generic.
References NOX::Abstract::Group::applyRightPreconditioning(), NOX::Abstract::Vector::clone(), NOX::Abstract::Group::computeF(), NOX::Abstract::Group::computeJacobian(), NOX::Abstract::Group::getF(), NOX::Solver::Generic::getNumIterations(), NOX::Solver::Generic::getPreviousSolutionGroup(), NOX::Abstract::Group::getX(), NOX::Abstract::Vector::innerProduct(), NOX::Abstract::Group::isJacobian(), NOX::Abstract::Group::Ok, NOX::Utils::OuterIteration, NOX::Abstract::Vector::scale(), NOX::ShapeCopy, NOX::Abstract::Vector::update(), and NOX::Utils::Warning.
|
virtual |
Same as compute(NOX::Abstract::Vector&, NOX::Abstract::Group&, const NOX::Solver::Generic&)
Enables direct support for line search based solvers for the purpose of efficiency since the LineSearchBased object has a getStep() function that some directions require.
If it is not redefined in the derived class, it will just call the compute with the NOX::Solver::Generic argument.
Reimplemented from NOX::Direction::Generic.
References NOX::Direction::Generic::compute().
|
virtual |