Reference documentation for deal.II version GIT relicensing-218-g1f873f3929 2024-03-28 15:00:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Classes | Public Member Functions | Public Attributes | Private Attributes | List of all members
TrilinosWrappers::NOXSolver< VectorType > Class Template Reference

#include <deal.II/trilinos/nox.h>

Inheritance diagram for TrilinosWrappers::NOXSolver< VectorType >:
Inheritance graph
[legend]

Classes

struct  AdditionalData
 

Public Member Functions

 NOXSolver (AdditionalData &additional_data, const Teuchos::RCP< Teuchos::ParameterList > &parameters=Teuchos::rcp(new Teuchos::ParameterList))
 
 ~NOXSolver ()
 
void clear ()
 
unsigned int solve (VectorType &solution)
 

Public Attributes

std::function< void(const VectorType &u, VectorType &F)> residual
 
std::function< void(const VectorType &current_u)> setup_jacobian
 
std::function< void(const VectorType &current_u)> setup_preconditioner
 
std::function< void(const VectorType &x, VectorType &y)> apply_jacobian
 
std::function< void(const VectorType &y, VectorType &x, const double tolerance)> solve_with_jacobian
 
std::function< int(const VectorType &y, VectorType &x, const double tolerance)> solve_with_jacobian_and_track_n_linear_iterations
 
std::function< SolverControl::State(const unsigned int i, const double norm_f, const VectorType &current_u, const VectorType &f)> check_iteration_status
 
std::function< bool()> update_preconditioner_predicate
 

Private Attributes

AdditionalData additional_data
 
const Teuchos::RCP< Teuchos::ParameterList > parameters
 
unsigned int n_residual_evaluations
 
unsigned int n_jacobian_applications
 
unsigned int n_nonlinear_iterations
 
unsigned int n_last_linear_iterations
 
std::exception_ptr pending_exception
 

Detailed Description

template<typename VectorType>
class TrilinosWrappers::NOXSolver< VectorType >

Wrapper around the nonlinear solver from the NOX package (https://docs.trilinos.org/dev/packages/nox/doc/html/index.html), targeting deal.II data structures.

The following code shows the steps how to use this class:

// set configuration
// Define ParameterList object for more options
// These specifications are the default but we include them for
// clarification
const Teuchos::RCP<Teuchos::ParameterList> parameters =
Teuchos::rcp(new Teuchos::ParameterList);
// Specify nonlinear solver type
parameters->set("Nonlinear Solver","Line Search Based");
// Specify method of line search
parameters->sublist("Line Search").set("Method","Full Step");
// Specify direction
parameters->sublist("Direction").set("Method","Newton")
// create nonlinear solver
// Set user functions to compute residual, to set up the Jacobian, and to
// apply the inverse of the Jacobian.
// Note that there are more functions that can be set.
solver.residual = [](const auto &u, auto &F) {...};
solver.setup_jacobian = [](const auto &u) {...};
solver.solve_with_jacobian =
[](const auto &u, auto &F, const auto) {...};
// solver nonlinear system with solution containing the initial guess and
// the final solution
solver.solve(solution);
AdditionalData additional_data
Definition nox.h:373
const Teuchos::RCP< Teuchos::ParameterList > parameters
Definition nox.h:380

The functions used in NOX are nearly identical to the functions in SUNDIALS::KINSOL with a few exceptions (for example, SUNDIALS::KINSOL requires a reinit() function where NOX does not). So check the SUNDIALS::KINSOL documentation for more precise details on how these functions are implemented.

Definition at line 87 of file nox.h.

Constructor & Destructor Documentation

◆ NOXSolver()

template<typename VectorType >
TrilinosWrappers::NOXSolver< VectorType >::NOXSolver ( AdditionalData additional_data,
const Teuchos::RCP< Teuchos::ParameterList > &  parameters = Teuchos::rcp(new Teuchos::ParameterList) 
)

Constructor, with class parameters set by the AdditionalData object.

Parameters
additional_dataNOX configuration data.
parametersMore specific NOX solver configuration.

If parameters is not filled, a Newton solver with a full step is used. An overview of possible parameters is given at https://docs.trilinos.org/dev/packages/nox/doc/html/parameters.html.

◆ ~NOXSolver()

template<typename VectorType >
TrilinosWrappers::NOXSolver< VectorType >::~NOXSolver ( )

Destructor.

Member Function Documentation

◆ clear()

template<typename VectorType >
void TrilinosWrappers::NOXSolver< VectorType >::clear ( )

Clear the internal state.

◆ solve()

template<typename VectorType >
unsigned int TrilinosWrappers::NOXSolver< VectorType >::solve ( VectorType &  solution)

Solve the nonlinear problem and return the number of nonlinear iterations executed.

Member Data Documentation

◆ residual

template<typename VectorType >
std::function<void(const VectorType &u, VectorType &F)> TrilinosWrappers::NOXSolver< VectorType >::residual

A function object that users should supply and that is intended to compute the residual \(F(u)\).

Note
This variable represents a user provided callback. See there for a description of how to deal with errors and other requirements and conventions. NOX can not deal with "recoverable" errors for this callback, so if it throws an exception of type RecoverableUserCallbackError, then this exception is treated like any other exception.

Definition at line 204 of file nox.h.

◆ setup_jacobian

template<typename VectorType >
std::function<void(const VectorType &current_u)> TrilinosWrappers::NOXSolver< VectorType >::setup_jacobian

A user function that sets up the Jacobian, based on the current solution current_u.

Note
This variable represents a user provided callback. See there for a description of how to deal with errors and other requirements and conventions. NOX can not deal with "recoverable" errors for this callback, so if it throws an exception of type RecoverableUserCallbackError, then this exception is treated like any other exception.

Definition at line 218 of file nox.h.

◆ setup_preconditioner

template<typename VectorType >
std::function<void(const VectorType &current_u)> TrilinosWrappers::NOXSolver< VectorType >::setup_preconditioner

A user function that sets up the preconditioner for inverting the Jacobian, based on the current solution current_u.

Note
This function is optional and is used when setup_jacobian is called and the preconditioner needs to be updated (see update_preconditioner_predicate and AdditionalData::threshold_nonlinear_iterations).
This variable represents a user provided callback. See there for a description of how to deal with errors and other requirements and conventions. NOX can not deal with "recoverable" errors for this callback, so if it throws an exception of type RecoverableUserCallbackError, then this exception is treated like any other exception.

Definition at line 237 of file nox.h.

◆ apply_jacobian

template<typename VectorType >
std::function<void(const VectorType &x, VectorType &y)> TrilinosWrappers::NOXSolver< VectorType >::apply_jacobian

A user function that applies the Jacobian \(\nabla_u F(u)\) to x and writes the result in y. The Jacobian to be used (i.e., more precisely: the linearization point \(u\) above) is the one computed when the setup_jacobian function was last called.

Note
This function is optional and is used in the case of certain configurations. For instance, this function is required if the polynomial line search (NOX::LineSearch::Polynomial) is chosen, whereas for the full step case (NOX::LineSearch::FullStep) it won't be called.
This variable represents a user provided callback. See there for a description of how to deal with errors and other requirements and conventions. NOX can not deal with "recoverable" errors for this callback, so if it throws an exception of type RecoverableUserCallbackError, then this exception is treated like any other exception.

Definition at line 259 of file nox.h.

◆ solve_with_jacobian

template<typename VectorType >
std::function< void(const VectorType &y, VectorType &x, const double tolerance)> TrilinosWrappers::NOXSolver< VectorType >::solve_with_jacobian

A user function that applies the inverse of the Jacobian \([\nabla_u F(u)]^{-1}\) to y and writes the result in x. The parameter tolerance specifies the error reduction if an iterative solver is used in applying the inverse matrix. The Jacobian to be used (i.e., more precisely: the linearization point \(u\) above) is the one computed when the setup_jacobian function was last called.

Note
This function is optional and is used in the case of certain configurations.
This variable represents a user provided callback. See there for a description of how to deal with errors and other requirements and conventions. NOX can deal with "recoverable" errors for this callback, if the NOX parameter "Newton/Rescue Bad Newton Solve" is set to true (which is, in fact, its default value). If this parameters is set to true, then exceptions of type RecoverableUserCallbackError are eaten for this callback and NOX can safely proceed with a recovery step. Exceptions of other types are still treated as "irrecoverable".

Definition at line 286 of file nox.h.

◆ solve_with_jacobian_and_track_n_linear_iterations

template<typename VectorType >
std::function< int(const VectorType &y, VectorType &x, const double tolerance)> TrilinosWrappers::NOXSolver< VectorType >::solve_with_jacobian_and_track_n_linear_iterations

A user function that applies the inverse of the Jacobian \([\nabla_u F(u)]^{-1}\) to y, writes the result in x and returns the number of linear iterations the linear solver needed. The parameter tolerance species the error reduction if an iterative solver is used. The Jacobian to be used (i.e., more precisely: the linearization point \(u\) above) is the one computed when the setup_jacobian function was last called.

Note
This function is used if solve_with_jacobian is not provided. Its return value is compared again AdditionalFlags::threshold_n_linear_iterations; if it is larger, the preconditioner will be built before the next linear system is solved. The use of this approach is predicated on the idea that one can keep using a preconditioner built earlier as long as it is a good preconditioner for the matrix currently in use – where "good" is defined as leading to a number of iterations to solve linear systems less than the threshold given by the current variable.
This variable represents a user provided callback. See there for a description of how to deal with errors and other requirements and conventions. NOX can deal with "recoverable" errors for this callback, if the NOX parameter "Newton/Rescue Bad Newton Solve" is set to true (which is, in fact, its default value). If this parameters is set to true, then exceptions of type RecoverableUserCallbackError are eaten for this callback and NOX can safely proceed with a recovery step. Exceptions of other types are still treated as "irrecoverable".

Definition at line 322 of file nox.h.

◆ check_iteration_status

template<typename VectorType >
std::function<SolverControl::State(const unsigned int i, const double norm_f, const VectorType &current_u, const VectorType &f)> TrilinosWrappers::NOXSolver< VectorType >::check_iteration_status

A user function that allows to check convergence in addition to ones checking the l2-norm and the number of iterations (see AdditionalData). It is run after each nonlinear iteration.

The input are the current iteration number i, the l2-norm norm_f of the residual vector, the current solution current_u, and the current residual vector f.

Note
This function is optional.
This variable represents a user provided callback. See there for a description of how to deal with errors and other requirements and conventions. NOX can not deal with "recoverable" errors for this callback, so if it throws an exception of type RecoverableUserCallbackError, then this exception is treated like any other exception.

Definition at line 347 of file nox.h.

◆ update_preconditioner_predicate

template<typename VectorType >
std::function<bool()> TrilinosWrappers::NOXSolver< VectorType >::update_preconditioner_predicate

A user function that, in addition to AdditionalData::threshold_nonlinear_iterations, allows to force to update the preconditioner. A reason for wanting to update the preconditioner is when the expected number of linear iterations is exceeded.

Note
This function is optional. If no function is attached, this means implicitly a return value of false.
This variable represents a user provided callback. See there for a description of how to deal with errors and other requirements and conventions. NOX can not deal with "recoverable" errors for this callback, so if it throws an exception of type RecoverableUserCallbackError, then this exception is treated like any other exception.

Definition at line 367 of file nox.h.

◆ additional_data

template<typename VectorType >
AdditionalData TrilinosWrappers::NOXSolver< VectorType >::additional_data
private

Additional data with basic settings.

Definition at line 373 of file nox.h.

◆ parameters

template<typename VectorType >
const Teuchos::RCP<Teuchos::ParameterList> TrilinosWrappers::NOXSolver< VectorType >::parameters
private

Additional data with advanced settings. An overview of possible parameters is given at https://docs.trilinos.org/dev/packages/nox/doc/html/parameters.html.

Definition at line 380 of file nox.h.

◆ n_residual_evaluations

template<typename VectorType >
unsigned int TrilinosWrappers::NOXSolver< VectorType >::n_residual_evaluations
private

A counter for the number of (accumulated) residual evaluations.

Definition at line 385 of file nox.h.

◆ n_jacobian_applications

template<typename VectorType >
unsigned int TrilinosWrappers::NOXSolver< VectorType >::n_jacobian_applications
private

A counter for the number of (accumulated) Jacobi applications.

Definition at line 390 of file nox.h.

◆ n_nonlinear_iterations

template<typename VectorType >
unsigned int TrilinosWrappers::NOXSolver< VectorType >::n_nonlinear_iterations
private

A counter for the number of (accumulated) nonlinear iterations.

Definition at line 395 of file nox.h.

◆ n_last_linear_iterations

template<typename VectorType >
unsigned int TrilinosWrappers::NOXSolver< VectorType >::n_last_linear_iterations
private

The number of linear iterations of the last Jacobian solve.

Definition at line 400 of file nox.h.

◆ pending_exception

template<typename VectorType >
std::exception_ptr TrilinosWrappers::NOXSolver< VectorType >::pending_exception
mutableprivate

A pointer to any exception that may have been thrown in user-defined call-backs and that we have to deal after the KINSOL function we call has returned.

Definition at line 407 of file nox.h.


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