![]() |
deal.II version GIT relicensing-3083-g7b89508ac7 2025-04-18 12:50:00+00:00
|
#include <deal.II/numerics/nonlinear.h>
Classes | |
class | AdditionalData |
Public Member Functions | |
NonlinearSolverSelector () | |
NonlinearSolverSelector (const AdditionalData &additional_data) | |
NonlinearSolverSelector (const AdditionalData &additional_data, const MPI_Comm &mpi_communicator) | |
void | select (const typename AdditionalData::SolverType &type) |
void | set_data (const AdditionalData &additional_data) |
void | set_data (const typename TrilinosWrappers::NOXSolver< VectorType >::AdditionalData &additional_data, const Teuchos::RCP< Teuchos::ParameterList > ¶meters=Teuchos::rcp(new Teuchos::ParameterList)) |
void | set_data (const typename SUNDIALS::KINSOL< VectorType >::AdditionalData &additional_data) |
void | set_data (const typename PETScWrappers::NonlinearSolverData &additional_data) |
void | solve (VectorType &initial_guess_and_solution) |
Public Attributes | |
std::function< void(VectorType &)> | reinit_vector |
std::function< void(const VectorType &src, VectorType &dst)> | residual |
std::function< void(const VectorType ¤t_u)> | setup_jacobian |
std::function< void(const VectorType &rhs, VectorType &dst, const double tolerance)> | solve_with_jacobian |
Private Member Functions | |
void | solve_with_petsc (VectorType &initial_guess_and_solution) |
void | solve_with_petsc (PETScWrappers::MPI::Vector &initial_guess_and_solution) |
Private Attributes | |
AdditionalData | additional_data |
MPI_Comm | mpi_communicator |
SUNDIALS::KINSOL< VectorType >::AdditionalData | additional_data_kinsol |
TrilinosWrappers::NOXSolver< VectorType >::AdditionalData | additional_data_nox |
Teuchos::RCP< Teuchos::ParameterList > | parameters_nox |
PETScWrappers::NonlinearSolverData | additional_data_petsc_snes |
This class offers a unified interface to several nonlinear solver implementations like KINSOL, NOX, and SNES.
KINSOL nonlinear solvers are part of the SUNDIALS package, while NOX is part of Trilinos and SNES is part of PETSc, respectively. If no solver is manually specified, this class will automaticlaly choose one of the available solvers depending on the enabled dependencies.
By calling the solve
function of this NonlinearSolverSelector
, it selects the solve
function of that Solver
that was specified in the constructor of this class, similar to the SolverSelector class.
An example of code one would run with this class is the following:
Definition at line 80 of file nonlinear.h.
NonlinearSolverSelector< VectorType >::NonlinearSolverSelector | ( | ) |
Constructor, filling in default values
Definition at line 495 of file nonlinear.h.
NonlinearSolverSelector< VectorType >::NonlinearSolverSelector | ( | const AdditionalData & | additional_data | ) |
Constructor, selecting the solver and other parameters specified in additional_data
.
Definition at line 502 of file nonlinear.h.
NonlinearSolverSelector< VectorType >::NonlinearSolverSelector | ( | const AdditionalData & | additional_data, |
const MPI_Comm & | mpi_communicator | ||
) |
Constructor.
additional_data | NonlinearSolverSelector configuration data |
mpi_communicator | MPI communicator used by the nonlinear solver. |
Definition at line 513 of file nonlinear.h.
void NonlinearSolverSelector< VectorType >::select | ( | const typename AdditionalData::SolverType & | type | ) |
Select a new nonlinear solver. All solver names used in this class are all lower case.
Definition at line 526 of file nonlinear.h.
void NonlinearSolverSelector< VectorType >::set_data | ( | const AdditionalData & | additional_data | ) |
Set the generic additional data for all nonlinear solvers.
Definition at line 412 of file nonlinear.h.
void NonlinearSolverSelector< VectorType >::set_data | ( | const typename TrilinosWrappers::NOXSolver< VectorType >::AdditionalData & | additional_data, |
const Teuchos::RCP< Teuchos::ParameterList > & | parameters = Teuchos::rcp(new Teuchos::ParameterList) |
||
) |
Set the additional data for NOX. See TrilinosWrappers::NOXSolver for more information.
Definition at line 557 of file nonlinear.h.
void NonlinearSolverSelector< VectorType >::set_data | ( | const typename SUNDIALS::KINSOL< VectorType >::AdditionalData & | additional_data | ) |
Set the additional data for KINSOL. See SUNDIALS::KINSOL for more information.
Definition at line 572 of file nonlinear.h.
void NonlinearSolverSelector< VectorType >::set_data | ( | const typename PETScWrappers::NonlinearSolverData & | additional_data | ) |
Set the additional data for SNES. See PETScWrappers::NonlinearSolverData for more information.
void NonlinearSolverSelector< VectorType >::solve | ( | VectorType & | initial_guess_and_solution | ) |
Solve the nonlinear system.
The content of initial_guess_and_solution
is used as an initial guess and the final solution is stored in the same vector.
Definition at line 622 of file nonlinear.h.
|
private |
Solve with PETSc SNES. Internal functions specialized for PETSc Vectors. This is necessary to only instantiate the SNES solvers with PETSc Vectors, as this is the only vector type currently supported.
Definition at line 583 of file nonlinear.h.
|
private |
Definition at line 595 of file nonlinear.h.
std::function<void(VectorType &)> NonlinearSolverSelector< VectorType >::reinit_vector |
A function object that users need to supply and that is intended to reinitize the given vector to its correct size, and block structure (if block vectors are used), along with any other properties necessary.
Definition at line 292 of file nonlinear.h.
std::function<void(const VectorType &src, VectorType &dst)> NonlinearSolverSelector< VectorType >::residual |
A function object that users should supply and that is intended to compute the residual dst = F(src)
.
Definition at line 307 of file nonlinear.h.
std::function<void(const VectorType ¤t_u)> NonlinearSolverSelector< VectorType >::setup_jacobian |
A function object that users may supply and that is intended to prepare the linear solver for subsequent calls to solve_with_jacobian).
The job of setup_jacobian() is to prepare the linear solver for subsequent calls to solve_with_jacobian(), in the solution of linear systems Ax = b. The exact nature of this system depends on the SolutionStrategy that has been selected.
In the cases strategy = SolutionStrategy::newton or SolutionStrategy::linesearch, A is the Jacobian J = \partial F/\partial u. If strategy = SolutionStrategy::picard, A is the approximate Jacobian matrix L.
current_u | Current value of u |
Definition at line 335 of file nonlinear.h.
std::function< void(const VectorType &rhs, VectorType &dst, const double tolerance)> NonlinearSolverSelector< VectorType >::solve_with_jacobian |
A function object that users may supply and that is intended to solve a linear system with the Jacobian matrix.
[in] | rhs | The system right hand side to solve for. |
[out] | dst | The solution of J^{-1} * \texttt{src}. |
[in] | tolerance | The tolerance with which to solve the linear system of equations. |
Definition at line 360 of file nonlinear.h.
|
private |
NonlinearSolverSelector configuration data.
Definition at line 366 of file nonlinear.h.
|
private |
The MPI communicator to be used by this solver, if any.
Definition at line 372 of file nonlinear.h.
|
private |
KINSOL configuration data
Definition at line 378 of file nonlinear.h.
|
private |
NOX configuration data
Definition at line 386 of file nonlinear.h.
|
private |
Definition at line 387 of file nonlinear.h.
|
private |
PETSc SNES configuration data
Definition at line 395 of file nonlinear.h.