Public Member Functions | |
def | __init__ |
def | reset |
def | getStatus |
def | step |
def | solve |
def | getSolutionGroup |
def | getPreviousSolutionGroup |
def | getNumIterations |
def | getList |
def | __init__ |
def | reset |
def | getStatus |
def | step |
def | solve |
def | getSolutionGroup |
def | getPreviousSolutionGroup |
def | getNumIterations |
def | getList |
Public Attributes | |
this |
Nonlinear solver based on a rank-1 tensor method. Solves $F(x)=0$ using a rank-1 tensor method and a linesearch globalization. At the kth nonlinear iteration, the solver does the following: Computes the tensor direction $ d_T $ by finding the root or smallest magnitude minimizer of the local model \\[ M_T(x_k+d) = F_k + J_kd + a_k(s_k^Td)^2, \\] where \\[ a_k = 2(F_{k-1} - F_k - J_ks_k) / (s_k^Ts_k)^2 \\] and \\[ s_k = s_{k-1} - s_k. \\] Modifies the step according to a global strategy and updates $x$ as $x_{k+1} = x_k + d(\\lambda) $ via a linesearch method, where $ d(\\lambda) $ is some function of $ \\lambda $. For instance, the curvilinear step $ d_{\\lambda T} $ is a function of the linesearch parameter $ \\lambda $ and is a parametric step that spans the directions of the tensor step and the Newton step. At $ \\lambda=1 $, the curvilinear step equals the full tensor step, and as $ \\lambda $ nears 0, the curvilinear step approaches the Newton direction. This step provides a monotonic decrease in the norm of the local tensor model as $ \\lambda $ varies from 0 to 1. The solver iterates until the status tests (see NOX::StatusTest) determine either failure or convergence. Input Parameters To use this solver, set the "Nonlinear Solver" parameter to be "Tensor Based". Then, specify the following sublists with the appropriate parameters as indicated below. "Direction" - Sublist of the direction parameters, passed to the NOX::Direction::Factory constructor. Defaults to an empty list. "Method" - Name of the direction to be computed in this solver. "Tensor" and "Newton" are the only two valid choices. A sublist by this name specifies all of the parameters to be passed to the linear solver. See below under "Linear Solver". "Rescue Bad Newton Solve" (Boolean) - If the linear solve does not meet the tolerance specified by the forcing term, then use the step anyway. Defaults to true. "Linear Solver" - Sublist for the specific linear solver parameters that are passed to NOX::Abstract::Group::computeNewton() and NOX::Abstract::Group::applyJacobianInverse(). "Linear Solver" is itself a sublist of the list specified in "Method" above (i.e., "Tensor" or "Newton"). Below is a partial list of standard parameters usually available in common linear solvers. Check with the specific linear solver being used for other parameters. "Max Iterations" - Maximum number of Arnoldi iterations (also max Krylov space dimension) "Tolerance" - Relative tolerance for solving local model [default = 1e-4] "Output Frequency" - Print output at every number of iterations [default = 20] "Line Search" - Sublist of the line search parameters. Because the tensor step is not guaranteed to be a descent direction on the function, not all "basic" line search approaches would be appropriate. Thus, the LineSearch classes available to Newton's method (e.g., Polynomial, More-Thuente) are not used here. Instead, this solver class approriately handles technical considerations for tensor methods with its own set of global strategies. The following parameters specify the specific options for this line search: "Method" - Name of the line search available to tensor methods Valid choices are: "Curvilinear" - Backtrack along the "curvilinear" path that spans the tensor direction and the Newton direction and that maintains monotonicity on the tensor model. Recommended because it tends to be more robust and efficient than the other choices. [Default] "Standard" - Backtrack along tensor direction unless it is not a descent direction, in which case backtrack along Newton direction. "Dual" - Backtrack along both the Newton and tensor directions and choose the better of the two. "Full Step" - Only use the full step and do not backtrack along both the Newton and tensor directions and choose the better of the two. "Lambda selection" - Flag for how to calculate the next linesearch parameter lambda. Valid choices are "Quadratic" and "Halving" (default). Quadratic constructs a quadratic interpolating polynomial from the last trial point and uses the minimum of this function as the next trial lambda (bounded by 0.1). Halving divides the linesearch parameter by 2 before each trial, which is simpler but tends to generate longer steps than quadratic. "Default Step" - Starting value of the linesearch parameter (defaults to 1.0) "Minimum Step" - Minimum acceptable linesearch parameter before the linesearch terminates (defaults to 1.0e-12). If there are many linesearch failures, then lowering this value is one thing to try. "Recovery Step Type" - Determines the step size to take when the line search fails. Choices are: "Constant" [default] - Uses a constant value set in "Recovery Step". "Last Computed Step" - Uses the last value computed by the line search algorithm. "Recovery Step" - Step parameter to take when the line search fails (defaults to value for "Default Step") "Max Iters" - Maximum number of iterations (i.e., backtracks) "Solver Options" - Sublist of general solver options. "User Defined Pre/Post Operator" is supported. See NOX::Parameter::PrePostOperator for more details. Output Parameters Every time solve() is called, a sublist for output parameters called "Output" will be created and will contain the following parameters: "Nonlinear Iterations" - Number of nonlinear iterations "2-Norm of Residual" - L-2 norm of the final residual $ F(x_k) $. References B. W. Bader, Tensor-Krylov methods for solving large-scale systems of nonlinear equations, Ph.D. Thesis, 2003, University of Colorado, Boulder, Colorado. B. W. Bader, Tensor-Krylov methods for solving large-scale systems of nonlinear equations, submitted to SIAM J. Numer. Anal. B. W. Bader and R. B. Schnabel, Curvilinear linesearch for tensor methods, SISC, 25(2):604-622. R. B. Schnabel and P. D. Frank, Tensor methods for nonlinear equations, SIAM J. Numer. Anal., 21(5):815-843. Brett Bader (SNL 9233) C++ includes: NOX_Solver_TensorBased.H
def PyTrilinos::NOX::Solver::TensorBased::__init__ | ( | self, | ||
args | ||||
) |
__init__(self, Teuchos::RCP<(NOX::Abstract::Group)> grp, Teuchos::RCP<(NOX::StatusTest::Generic)> tests, Teuchos::RCP<(Teuchos::ParameterList)> params) -> TensorBased NOX::Solver::TensorBased::TensorBased(const Teuchos::RCP< NOX::Abstract::Group > &grp, const Teuchos::RCP< NOX::StatusTest::Generic > &tests, const Teuchos::RCP< Teuchos::ParameterList > ¶ms) Constructor. See reset() for description.
def PyTrilinos::NOX::Solver::TensorBased::__init__ | ( | self, | ||
args | ||||
) |
__init__(self, Teuchos::RCP<(NOX::Abstract::Group)> grp, Teuchos::RCP<(NOX::StatusTest::Generic)> tests, Teuchos::RCP<(Teuchos::ParameterList)> params) -> TensorBased NOX::Solver::TensorBased::TensorBased(const Teuchos::RCP< NOX::Abstract::Group > &grp, const Teuchos::RCP< NOX::StatusTest::Generic > &tests, const Teuchos::RCP< Teuchos::ParameterList > ¶ms) Constructor. See reset() for description.
def PyTrilinos::NOX::Solver::TensorBased::getList | ( | self, | ||
args | ||||
) |
getList(self) -> ParameterList const Teuchos::ParameterList & NOX::Solver::TensorBased::getList() const Return a refernece to the solver parameters.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::getList | ( | self, | ||
args | ||||
) |
getList(self) -> ParameterList const Teuchos::ParameterList & NOX::Solver::TensorBased::getList() const Return a refernece to the solver parameters.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::getNumIterations | ( | self, | ||
args | ||||
) |
getNumIterations(self) -> int int NOX::Solver::TensorBased::getNumIterations() const Get number of iterations.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::getNumIterations | ( | self, | ||
args | ||||
) |
getNumIterations(self) -> int int NOX::Solver::TensorBased::getNumIterations() const Get number of iterations.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::getPreviousSolutionGroup | ( | self, | ||
args | ||||
) |
getPreviousSolutionGroup(self) -> Group const NOX::Abstract::Group & NOX::Solver::TensorBased::getPreviousSolutionGroup() const Return a reference to the previous solution group.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::getPreviousSolutionGroup | ( | self, | ||
args | ||||
) |
getPreviousSolutionGroup(self) -> Group const NOX::Abstract::Group & NOX::Solver::TensorBased::getPreviousSolutionGroup() const Return a reference to the previous solution group.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::getSolutionGroup | ( | self, | ||
args | ||||
) |
getSolutionGroup(self) -> Group const NOX::Abstract::Group & NOX::Solver::TensorBased::getSolutionGroup() const Return a reference to the current solution group.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::getSolutionGroup | ( | self, | ||
args | ||||
) |
getSolutionGroup(self) -> Group const NOX::Abstract::Group & NOX::Solver::TensorBased::getSolutionGroup() const Return a reference to the current solution group.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::getStatus | ( | self, | ||
args | ||||
) |
getStatus(self) -> StatusType NOX::StatusTest::StatusType NOX::Solver::TensorBased::getStatus() Check current convergence and failure status.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::getStatus | ( | self, | ||
args | ||||
) |
getStatus(self) -> StatusType NOX::StatusTest::StatusType NOX::Solver::TensorBased::getStatus() Check current convergence and failure status.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::reset | ( | self, | ||
args | ||||
) |
reset(self, Vector initialGuess, Teuchos::RCP<(NOX::StatusTest::Generic)> tests) reset(self, Vector initialGuess) void NOX::Solver::TensorBased::reset(const NOX::Abstract::Vector &initialGuess) Resets the solver and sets a new initial guess.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::reset | ( | self, | ||
args | ||||
) |
reset(self, Vector initialGuess, Teuchos::RCP<(NOX::StatusTest::Generic)> tests) reset(self, Vector initialGuess) void NOX::Solver::TensorBased::reset(const NOX::Abstract::Vector &initialGuess) Resets the solver and sets a new initial guess.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::solve | ( | self, | ||
args | ||||
) |
solve(self) -> StatusType NOX::StatusTest::StatusType NOX::Solver::TensorBased::solve() Solve the nonlinear problem and return final status. By "solve", we call iterate() until the NOX::StatusTest value is either NOX::StatusTest::Converged or NOX::StatusTest::Failed.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::solve | ( | self, | ||
args | ||||
) |
solve(self) -> StatusType NOX::StatusTest::StatusType NOX::Solver::TensorBased::solve() Solve the nonlinear problem and return final status. By "solve", we call iterate() until the NOX::StatusTest value is either NOX::StatusTest::Converged or NOX::StatusTest::Failed.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::step | ( | self, | ||
args | ||||
) |
step(self) -> StatusType NOX::StatusTest::StatusType NOX::Solver::TensorBased::step() Do one nonlinear step in the iteration sequence and return status.
Reimplemented from PyTrilinos::NOX::Solver::Generic.
def PyTrilinos::NOX::Solver::TensorBased::step | ( | self, | ||
args | ||||
) |
step(self) -> StatusType NOX::StatusTest::StatusType NOX::Solver::TensorBased::step() Do one nonlinear step in the iteration sequence and return status.
Reimplemented from PyTrilinos::NOX::Solver::Generic.