libMesh
Classes | Public Types | Public Member Functions | Static Public Member Functions | Public Attributes | Protected Types | Protected Member Functions | Protected Attributes | Static Protected Attributes | Friends | List of all members
libMesh::PetscNonlinearSolver< T > Class Template Reference

This class provides an interface to PETSc iterative solvers that is compatible with the libMesh NonlinearSolver<> More...

#include <petsc_nonlinear_solver.h>

Inheritance diagram for libMesh::PetscNonlinearSolver< T >:
[legend]

Classes

class  ComputeLineSearchObject
 Abstract base class to be used to implement a custom line-search algorithm. More...
 

Public Types

typedef NonlinearImplicitSystem sys_type
 The type of system. More...
 

Public Member Functions

 PetscNonlinearSolver (sys_type &system)
 Constructor. More...
 
 ~PetscNonlinearSolver ()
 Destructor. More...
 
virtual void clear () override
 Release all memory and clear data structures. More...
 
virtual void init (const char *name=nullptr) override
 Initialize data structures if not done so already. More...
 
SNES snes ()
 
virtual std::pair< unsigned int, Realsolve (SparseMatrix< T > &, NumericVector< T > &, NumericVector< T > &, const double, const unsigned int) override
 Call the Petsc solver. More...
 
virtual void print_converged_reason () override
 Prints a useful message about why the latest nonlinear solve con(di)verged. More...
 
SNESConvergedReason get_converged_reason ()
 
virtual int get_total_linear_iterations () override
 Get the total number of linear iterations done in the last solve. More...
 
virtual unsigned get_current_nonlinear_iteration_number () const override
 
void set_residual_zero_out (bool state)
 Set if the residual should be zeroed out in the callback. More...
 
void set_jacobian_zero_out (bool state)
 Set if the jacobian should be zeroed out in the callback. More...
 
void use_default_monitor (bool state)
 Set to true to use the libMesh's default monitor, set to false to use your own. More...
 
void set_snesmf_reuse_base (bool state)
 Set to true to let PETSc reuse the base vector. More...
 
bool snes_mf_reuse_base () const
 
void set_computing_base_vector (bool computing_base_vector)
 Set whether we are computing the base vector for matrix-free finite-differencing. More...
 
bool computing_base_vector () const
 
void setup_default_monitor ()
 Setup the default monitor if required. More...
 
virtual bool reuse_preconditioner () const override
 Getter for preconditioner reuse. More...
 
virtual unsigned int reuse_preconditioner_max_linear_its () const override
 Getter for the maximum iterations flag for preconditioner reuse. More...
 
virtual void force_new_preconditioner () override
 Immediately force a new preconditioner, even if reuse is set. More...
 
bool initialized () const
 
const sys_typesystem () const
 
sys_typesystem ()
 
void attach_preconditioner (Preconditioner< T > *preconditioner)
 Attaches a Preconditioner object to be used during the linear solves. More...
 
void set_solver_configuration (SolverConfiguration &solver_configuration)
 Set the solver configuration object. More...
 
virtual void set_reuse_preconditioner (bool reuse)
 Set the reuse preconditioner flag. More...
 
virtual void set_reuse_preconditioner_max_linear_its (unsigned int i)
 Set the reuse_preconditioner_max_linear_its parameter. More...
 
virtual void set_exact_constraint_enforcement (bool enable)
 Enable (or disable; it is true by default) exact enforcement of constraints at the solver level, correcting any constrained DoF coefficients in current_local_solution as well as applying nonlinear residual and Jacobian terms based on constraint equations. More...
 
bool exact_constraint_enforcement ()
 
const Parallel::Communicatorcomm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Static Public Member Functions

static std::unique_ptr< NonlinearSolver< T > > build (sys_type &s, const SolverPackage solver_package=libMesh::default_solver_package())
 Builds a NonlinearSolver using the nonlinear solver package specified by solver_package. More...
 
static std::string get_info ()
 Gets a string containing the reference information. More...
 
static void print_info (std::ostream &out_stream=libMesh::out)
 Prints the reference information, by default to libMesh::out. More...
 
static unsigned int n_objects ()
 Prints the number of outstanding (created, but not yet destroyed) objects. More...
 
static void enable_print_counter_info ()
 Methods to enable/disable the reference counter output from print_info() More...
 
static void disable_print_counter_info ()
 

Public Attributes

std::unique_ptr< ComputeLineSearchObjectlinesearch_object
 A callable object that can be used to specify a custom line-search. More...
 
void(* residual )(const NumericVector< Number > &X, NumericVector< Number > &R, sys_type &S)
 Function that computes the residual R(X) of the nonlinear system at the input iterate X. More...
 
NonlinearImplicitSystem::ComputeResidualresidual_object
 Object that computes the residual R(X) of the nonlinear system at the input iterate X. More...
 
NonlinearImplicitSystem::ComputeResidualfd_residual_object
 Object that computes the residual R(X) of the nonlinear system at the input iterate X for the purpose of forming a finite-differenced Jacobian. More...
 
NonlinearImplicitSystem::ComputeResidualmffd_residual_object
 Object that computes the residual R(X) of the nonlinear system at the input iterate X for the purpose of forming Jacobian-vector products via finite differencing. More...
 
void(* jacobian )(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)
 Function that computes the Jacobian J(X) of the nonlinear system at the input iterate X. More...
 
NonlinearImplicitSystem::ComputeJacobianjacobian_object
 Object that computes the Jacobian J(X) of the nonlinear system at the input iterate X. More...
 
void(* matvec )(const NumericVector< Number > &X, NumericVector< Number > *R, SparseMatrix< Number > *J, sys_type &S)
 Function that computes either the residual \( R(X) \) or the Jacobian \( J(X) \) of the nonlinear system at the input iterate \( X \). More...
 
NonlinearImplicitSystem::ComputeResidualandJacobianresidual_and_jacobian_object
 Object that computes either the residual \( R(X) \) or the Jacobian \( J(X) \) of the nonlinear system at the input iterate \( X \). More...
 
void(* bounds )(NumericVector< Number > &XL, NumericVector< Number > &XU, sys_type &S)
 Function that computes the lower and upper bounds XL and XU on the solution of the nonlinear system. More...
 
NonlinearImplicitSystem::ComputeBoundsbounds_object
 Object that computes the bounds vectors \( XL \) and \( XU \). More...
 
void(* nullspace )(std::vector< NumericVector< Number > *> &sp, sys_type &S)
 Function that computes a basis for the Jacobian's nullspace – the kernel or the "zero energy modes" – that can be used in solving a degenerate problem iteratively, if the solver supports it (e.g., PETSc's KSP). More...
 
NonlinearImplicitSystem::ComputeVectorSubspacenullspace_object
 A callable object that computes a basis for the Jacobian's nullspace – the kernel or the "zero energy modes" – that can be used in solving a degenerate problem iteratively, if the solver supports it (e.g., PETSc's KSP). More...
 
void(* transpose_nullspace )(std::vector< NumericVector< Number > *> &sp, sys_type &S)
 Function that computes a basis for the transpose Jacobian's nullspace – when solving a degenerate problem iteratively, if the solver supports it (e.g., PETSc's KSP), it is used to remove contributions outside of R(jac) More...
 
NonlinearImplicitSystem::ComputeVectorSubspacetranspose_nullspace_object
 A callable object that computes a basis for the transpose Jacobian's nullspace – when solving a degenerate problem iteratively, if the solver supports it (e.g., PETSc's KSP), it is used to remove contributions outside of R(jac) More...
 
void(* nearnullspace )(std::vector< NumericVector< Number > *> &sp, sys_type &S)
 Function that computes a basis for the Jacobian's near nullspace – the set of "low energy modes" – that can be used for AMG coarsening, if the solver supports it (e.g., ML, PETSc's GAMG). More...
 
NonlinearImplicitSystem::ComputeVectorSubspacenearnullspace_object
 A callable object that computes a basis for the Jacobian's near nullspace – the set of "low energy modes" – that can be used for AMG coarsening, if the solver supports it (e.g., ML, PETSc's GAMG). More...
 
void(* user_presolve )(sys_type &S)
 Customizable function pointer which users can attach to the solver. More...
 
void(* postcheck )(const NumericVector< Number > &old_soln, NumericVector< Number > &search_direction, NumericVector< Number > &new_soln, bool &changed_search_direction, bool &changed_new_soln, sys_type &S)
 Function that performs a "check" on the Newton search direction and solution after each nonlinear step. More...
 
NonlinearImplicitSystem::ComputePostCheckpostcheck_object
 A callable object that is executed after each nonlinear iteration. More...
 
NonlinearImplicitSystem::ComputePreCheckprecheck_object
 
unsigned int max_nonlinear_iterations
 Maximum number of non-linear iterations. More...
 
unsigned int max_function_evaluations
 Maximum number of function evaluations. More...
 
double absolute_residual_tolerance
 The NonlinearSolver should exit after the residual is reduced to either less than absolute_residual_tolerance or less than relative_residual_tolerance times the initial residual. More...
 
double relative_residual_tolerance
 
double divergence_tolerance
 The NonlinearSolver should exit if the residual becomes greater than the initial residual times the divergence_tolerance. More...
 
double absolute_step_tolerance
 The NonlinearSolver should exit after the full nonlinear step norm is reduced to either less than absolute_step_tolerance or less than relative_step_tolerance times the largest nonlinear solution which has been seen so far. More...
 
double relative_step_tolerance
 
unsigned int max_linear_iterations
 Each linear solver step should exit after max_linear_iterations is exceeded. More...
 
double initial_linear_tolerance
 Any required linear solves will at first be done with this tolerance; the NonlinearSolver may tighten the tolerance for later solves. More...
 
double minimum_linear_tolerance
 The tolerance for linear solves is kept above this minimum. More...
 
bool converged
 After a call to solve this will reflect whether or not the nonlinear solve was successful. More...
 

Protected Types

typedef std::map< std::string, std::pair< unsigned int, unsigned int > > Counts
 Data structure to log the information. More...
 

Protected Member Functions

void build_mat_null_space (NonlinearImplicitSystem::ComputeVectorSubspace *computeSubspaceObject, void(*)(std::vector< NumericVector< Number > *> &, sys_type &), MatNullSpace *)
 
void increment_constructor_count (const std::string &name) noexcept
 Increments the construction counter. More...
 
void increment_destructor_count (const std::string &name) noexcept
 Increments the destruction counter. More...
 

Protected Attributes

WrappedPetsc< SNES > _snes
 Nonlinear solver context. More...
 
SNESConvergedReason _reason
 Store the reason for SNES convergence/divergence for use even after the _snes has been cleared. More...
 
PetscInt _n_linear_iterations
 Stores the total number of linear iterations from the last solve. More...
 
unsigned _current_nonlinear_iteration_number
 Stores the current nonlinear iteration number. More...
 
bool _zero_out_residual
 true to zero out residual before going into application level call-back, otherwise false More...
 
bool _zero_out_jacobian
 true to zero out jacobian before going into application level call-back, otherwise false More...
 
bool _default_monitor
 true if we want the default monitor to be set, false for no monitor (i.e. More...
 
bool _snesmf_reuse_base
 True, If we want the base vector to be used for differencing even if the function provided to SNESSetFunction() is not the same as that provided to MatMFFDSetFunction(). More...
 
bool _computing_base_vector
 Whether we are computing the base vector for matrix-free finite differencing. More...
 
bool _setup_reuse
 Whether we've triggered the preconditioner reuse. More...
 
bool _reuse_preconditioner
 Whether we should reuse the linear preconditioner. More...
 
bool _exact_constraint_enforcement
 Whether we should enforce exact constraints globally during a solve. More...
 
unsigned int _reuse_preconditioner_max_linear_its
 Number of linear iterations to retain the preconditioner. More...
 
sys_type_system
 A reference to the system we are solving. More...
 
bool _is_initialized
 Flag indicating if the data structures have been initialized. More...
 
Preconditioner< T > * _preconditioner
 Holds the Preconditioner object to be used for the linear solves. More...
 
SolverConfiguration_solver_configuration
 Optionally store a SolverOptions object that can be used to set parameters like solver type, tolerances and iteration limits. More...
 
const Parallel::Communicator_communicator
 

Static Protected Attributes

static Counts _counts
 Actually holds the data. More...
 
static Threads::atomic< unsigned int_n_objects
 The number of objects. More...
 
static Threads::spin_mutex _mutex
 Mutual exclusion object to enable thread-safe reference counting. More...
 
static bool _enable_print_counter = true
 Flag to control whether reference count information is printed when print_info is called. More...
 

Friends

ResidualContext libmesh_petsc_snes_residual_helper (SNES snes, Vec x, void *ctx)
 
PetscErrorCode libmesh_petsc_snes_residual (SNES snes, Vec x, Vec r, void *ctx)
 
PetscErrorCode libmesh_petsc_snes_fd_residual (SNES snes, Vec x, Vec r, void *ctx)
 
PetscErrorCode libmesh_petsc_snes_mffd_residual (SNES snes, Vec x, Vec r, void *ctx)
 
PetscErrorCode libmesh_petsc_snes_jacobian (SNES snes, Vec x, Mat jac, Mat pc, void *ctx)
 
PetscErrorCode libmesh_petsc_snes_precheck (SNESLineSearch, Vec X, Vec Y, PetscBool *changed, void *context)
 

Detailed Description

template<typename T>
class libMesh::PetscNonlinearSolver< T >

This class provides an interface to PETSc iterative solvers that is compatible with the libMesh NonlinearSolver<>

Author
Benjamin Kirk
Date
2002-2007

Definition at line 80 of file petsc_nonlinear_solver.h.

Member Typedef Documentation

◆ Counts

typedef std::map<std::string, std::pair<unsigned int, unsigned int> > libMesh::ReferenceCounter::Counts
protectedinherited

Data structure to log the information.

The log is identified by the class name.

Definition at line 119 of file reference_counter.h.

◆ sys_type

The type of system.

Definition at line 86 of file petsc_nonlinear_solver.h.

Constructor & Destructor Documentation

◆ PetscNonlinearSolver()

template<typename T >
libMesh::PetscNonlinearSolver< T >::PetscNonlinearSolver ( sys_type system)
explicit

Constructor.

Initializes Petsc data structures

Definition at line 727 of file petsc_nonlinear_solver.C.

727  :
728  NonlinearSolver<T>(system_in),
729  linesearch_object(nullptr),
730  _reason(SNES_CONVERGED_ITERATING/*==0*/), // Arbitrary initial value...
733  _zero_out_residual(true),
734  _zero_out_jacobian(true),
735  _default_monitor(true),
736  _snesmf_reuse_base(true),
738  _setup_reuse(false)
739 {
740 }
unsigned _current_nonlinear_iteration_number
Stores the current nonlinear iteration number.
bool _zero_out_residual
true to zero out residual before going into application level call-back, otherwise false ...
bool _setup_reuse
Whether we&#39;ve triggered the preconditioner reuse.
bool _zero_out_jacobian
true to zero out jacobian before going into application level call-back, otherwise false ...
PetscInt _n_linear_iterations
Stores the total number of linear iterations from the last solve.
SNESConvergedReason _reason
Store the reason for SNES convergence/divergence for use even after the _snes has been cleared...
std::unique_ptr< ComputeLineSearchObject > linesearch_object
A callable object that can be used to specify a custom line-search.
bool _default_monitor
true if we want the default monitor to be set, false for no monitor (i.e.
bool _computing_base_vector
Whether we are computing the base vector for matrix-free finite differencing.
bool _snesmf_reuse_base
True, If we want the base vector to be used for differencing even if the function provided to SNESSet...

◆ ~PetscNonlinearSolver()

template<typename T >
libMesh::PetscNonlinearSolver< T >::~PetscNonlinearSolver ( )
default

Destructor.

Member Function Documentation

◆ attach_preconditioner()

template<typename T>
void libMesh::NonlinearSolver< T >::attach_preconditioner ( Preconditioner< T > *  preconditioner)
inherited

Attaches a Preconditioner object to be used during the linear solves.

Definition at line 74 of file nonlinear_solver.C.

75 {
76  libmesh_error_msg_if(this->_is_initialized, "Preconditioner must be attached before the solver is initialized!");
77 
78  _preconditioner = preconditioner;
79 }
Preconditioner< T > * _preconditioner
Holds the Preconditioner object to be used for the linear solves.
bool _is_initialized
Flag indicating if the data structures have been initialized.

◆ build()

template<typename T >
std::unique_ptr< NonlinearSolver< T > > libMesh::NonlinearSolver< T >::build ( sys_type s,
const SolverPackage  solver_package = libMesh::default_solver_package() 
)
staticinherited

Builds a NonlinearSolver using the nonlinear solver package specified by solver_package.

Definition at line 40 of file nonlinear_solver.C.

41 {
42  // Build the appropriate solver
43  switch (solver_package)
44  {
45 
46 #ifdef LIBMESH_HAVE_PETSC
47  case PETSC_SOLVERS:
48  return std::make_unique<PetscNonlinearSolver<T>>(s);
49 #endif // LIBMESH_HAVE_PETSC
50 
51 #if defined(LIBMESH_TRILINOS_HAVE_NOX) && defined(LIBMESH_TRILINOS_HAVE_EPETRA)
52  case TRILINOS_SOLVERS:
53  return std::make_unique<NoxNonlinearSolver<T>>(s);
54 #endif
55 
56  default:
57  libmesh_error_msg("ERROR: Unrecognized solver package: " << solver_package);
58  }
59 }

◆ build_mat_null_space()

template<typename T>
void libMesh::PetscNonlinearSolver< T >::build_mat_null_space ( NonlinearImplicitSystem::ComputeVectorSubspace computeSubspaceObject,
void(*)(std::vector< NumericVector< Number > * > &, sys_type &)  ,
MatNullSpace *   
)
protected

Definition at line 885 of file petsc_nonlinear_solver.C.

888 {
889  parallel_object_only();
890 
891  PetscErrorCode ierr;
892  std::vector<NumericVector<Number> *> sp;
893  if (computeSubspaceObject)
894  (*computeSubspaceObject)(sp, this->system());
895  else
896  (*computeSubspace)(sp, this->system());
897 
898  *msp = LIBMESH_PETSC_NULLPTR;
899  if (sp.size())
900  {
901  PetscInt nmodes = cast_int<PetscInt>(sp.size());
902 
903  std::vector<Vec> modes(nmodes);
904  std::vector<PetscScalar> dots(nmodes);
905 
906  for (PetscInt i=0; i<nmodes; ++i)
907  {
908  auto pv = cast_ptr<PetscVector<T> *>(sp[i]);
909 
910  ierr = VecDuplicate(pv->vec(), &modes[i]);
911  LIBMESH_CHKERR(ierr);
912 
913  ierr = VecCopy(pv->vec(), modes[i]);
914  LIBMESH_CHKERR(ierr);
915  }
916 
917  // Normalize.
918  ierr = VecNormalize(modes[0], LIBMESH_PETSC_NULLPTR);
919  LIBMESH_CHKERR(ierr);
920 
921  for (PetscInt i=1; i<nmodes; i++)
922  {
923  // Orthonormalize vec[i] against vec[0:i-1]
924  ierr = VecMDot(modes[i], i, modes.data(), dots.data());
925  LIBMESH_CHKERR(ierr);
926 
927  for (PetscInt j=0; j<i; j++)
928  dots[j] *= -1.;
929 
930  ierr = VecMAXPY(modes[i], i, dots.data(), modes.data());
931  LIBMESH_CHKERR(ierr);
932 
933  ierr = VecNormalize(modes[i], LIBMESH_PETSC_NULLPTR);
934  LIBMESH_CHKERR(ierr);
935  }
936 
937  ierr = MatNullSpaceCreate(this->comm().get(), PETSC_FALSE, nmodes, modes.data(), msp);
938  LIBMESH_CHKERR(ierr);
939 
940  for (PetscInt i=0; i<nmodes; ++i)
941  {
942  ierr = VecDestroy(&modes[i]);
943  LIBMESH_CHKERR(ierr);
944  }
945  }
946 }
const Parallel::Communicator & comm() const
const sys_type & system() const

◆ clear()

template<typename T >
void libMesh::PetscNonlinearSolver< T >::clear ( )
overridevirtual

Release all memory and clear data structures.

Reimplemented from libMesh::NonlinearSolver< T >.

Definition at line 750 of file petsc_nonlinear_solver.C.

751 {
752  if (this->initialized())
753  {
754  this->_is_initialized = false;
755 
756  // If we don't need the preconditioner next time
757  // retain the original behavior of clearing the data
758  // between solves.
759  if (!(reuse_preconditioner()))
760  {
761  // SNESReset really ought to work but replacing destroy() with
762  // SNESReset causes a very slight change in behavior that
763  // manifests as two failed MOOSE tests...
764  _snes.destroy();
765  }
766 
767  // Reset the nonlinear iteration counter. This information is only relevant
768  // *during* the solve(). After the solve is completed it should return to
769  // the default value of 0.
771  }
772 }
unsigned _current_nonlinear_iteration_number
Stores the current nonlinear iteration number.
WrappedPetsc< SNES > _snes
Nonlinear solver context.
bool _is_initialized
Flag indicating if the data structures have been initialized.
virtual bool reuse_preconditioner() const override
Getter for preconditioner reuse.
void destroy()
Must be specialized to call the appropriate XXXDestroy() routine in order for a WrappedPetsc<T> objec...

◆ comm()

const Parallel::Communicator& libMesh::ParallelObject::comm ( ) const
inlineinherited
Returns
A reference to the Parallel::Communicator object used by this mesh.

Definition at line 97 of file parallel_object.h.

References libMesh::ParallelObject::_communicator.

Referenced by libMesh::__libmesh_petsc_diff_solver_jacobian(), libMesh::__libmesh_petsc_diff_solver_monitor(), libMesh::__libmesh_petsc_diff_solver_residual(), libMesh::__libmesh_tao_equality_constraints(), libMesh::__libmesh_tao_equality_constraints_jacobian(), libMesh::__libmesh_tao_gradient(), libMesh::__libmesh_tao_hessian(), libMesh::__libmesh_tao_inequality_constraints(), libMesh::__libmesh_tao_inequality_constraints_jacobian(), libMesh::__libmesh_tao_objective(), libMesh::MeshRefinement::_coarsen_elements(), libMesh::ExactSolution::_compute_error(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::Partitioner::_find_global_index_by_pid_map(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::SlepcEigenSolver< libMesh::Number >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< Number >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< libMesh::Number >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< Number >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< Number >::_petsc_shell_matrix_mult_add(), libMesh::MeshRefinement::_refine_elements(), libMesh::MeshRefinement::_smooth_flags(), libMesh::DofMap::add_constraints_to_send_list(), add_cube_convex_hull_to_mesh(), libMesh::PetscDMWrapper::add_dofs_helper(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::TransientRBConstruction::add_IC_to_RB_space(), libMesh::EigenSystem::add_matrices(), libMesh::System::add_matrix(), libMesh::RBConstruction::add_scaled_matrix_and_vector(), libMesh::System::add_variable(), libMesh::System::add_variables(), libMesh::System::add_vector(), libMesh::MeshTools::Modification::all_tri(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::DofMap::allgather_recursive_constraints(), libMesh::TransientRBConstruction::allocate_data_structures(), libMesh::RBConstruction::allocate_data_structures(), libMesh::TransientRBConstruction::assemble_affine_expansion(), libMesh::FEMSystem::assemble_qoi(), libMesh::Nemesis_IO::assert_symmetric_cmaps(), libMesh::MeshCommunication::assign_global_indices(), libMesh::Partitioner::assign_partitioning(), libMesh::MeshTools::Generation::build_extrusion(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::EquationSystems::build_parallel_solution_vector(), libMesh::PetscDMWrapper::build_section(), libMesh::PetscDMWrapper::build_sf(), libMesh::MeshBase::cache_elem_data(), libMesh::System::calculate_norm(), libMesh::DofMap::check_dirichlet_bcid_consistency(), libMesh::RBConstruction::compute_Fq_representor_innerprods(), libMesh::RBConstruction::compute_max_error_bound(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::RBConstruction::compute_output_dual_innerprods(), libMesh::RBConstruction::compute_residual_dual_norm_slow(), libMesh::RBSCMConstruction::compute_SCM_bounds_on_training_set(), libMesh::DofMap::computed_sparsity_already(), libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::ContinuationSystem::ContinuationSystem(), libMesh::MeshBase::copy_constraint_rows(), libMesh::ExodusII_IO::copy_elemental_solution(), libMesh::ExodusII_IO::copy_nodal_solution(), libMesh::ExodusII_IO::copy_scalar_solution(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::MeshTools::create_bounding_box(), libMesh::DofMap::create_dof_constraints(), libMesh::MeshTools::create_nodal_bounding_box(), libMesh::MeshRefinement::create_parent_error_vector(), libMesh::MeshTools::create_processor_bounding_box(), libMesh::MeshTools::create_subdomain_bounding_box(), libMesh::PetscMatrix< libMesh::Number >::create_submatrix_nosort(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::RBEIMEvaluation::distribute_bfs(), libMesh::DofMap::distribute_dofs(), DMlibMeshFunction(), DMlibMeshJacobian(), DMlibMeshSetSystem_libMesh(), DMVariableBounds_libMesh(), libMesh::DTKSolutionTransfer::DTKSolutionTransfer(), libMesh::MeshRefinement::eliminate_unrefined_patches(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_interiors(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_nodes(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_sides(), libMesh::TransientRBConstruction::enrich_RB_space(), libMesh::EpetraVector< T >::EpetraVector(), AssembleOptimization::equality_constraints(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::MeshRefinement::flag_elements_by_elem_fraction(), libMesh::MeshRefinement::flag_elements_by_error_fraction(), libMesh::MeshRefinement::flag_elements_by_error_tolerance(), libMesh::MeshRefinement::flag_elements_by_mean_stddev(), libMesh::MeshRefinement::flag_elements_by_nelem_target(), libMesh::RBEIMEvaluation::gather_bfs(), libMesh::DofMap::gather_constraints(), libMesh::MeshfreeInterpolation::gather_remote_data(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::RBEIMEvaluation::get_eim_basis_function_node_value(), libMesh::RBEIMEvaluation::get_eim_basis_function_side_value(), libMesh::RBEIMEvaluation::get_eim_basis_function_value(), libMesh::MeshBase::get_info(), libMesh::System::get_info(), libMesh::DofMap::get_info(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::RBEIMConstruction::get_max_abs_value(), libMesh::RBEIMConstruction::get_node_max_abs_value(), libMesh::RBEIMEvaluation::get_parametrized_function_node_value(), libMesh::RBEIMEvaluation::get_parametrized_function_side_value(), libMesh::RBEIMEvaluation::get_parametrized_function_value(), libMesh::RBEIMConstruction::get_random_point(), AssembleOptimization::inequality_constraints(), AssembleOptimization::inequality_constraints_jacobian(), libMesh::LocationMap< T >::init(), libMesh::TimeSolver::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), libMesh::ExodusII_IO_Helper::initialize(), libMesh::OptimizationSystem::initialize_equality_constraints_storage(), libMesh::OptimizationSystem::initialize_inequality_constraints_storage(), libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set(), libMesh::RBEIMConstruction::inner_product(), integrate_function(), libMesh::MeshTools::libmesh_assert_consistent_distributed(), libMesh::MeshTools::libmesh_assert_consistent_distributed_nodes(), libMesh::MeshTools::libmesh_assert_contiguous_dof_ids(), libMesh::MeshTools::libmesh_assert_equal_connectivity(), libMesh::MeshTools::libmesh_assert_equal_points(), libMesh::MeshTools::libmesh_assert_parallel_consistent_new_node_procids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_topology_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_boundary_ids(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_flags(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_p_levels(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::MeshTools::libmesh_assert_valid_unique_ids(), libMesh::libmesh_petsc_linesearch_shellfunc(), libMesh::libmesh_petsc_preconditioner_apply(), libMesh::libmesh_petsc_recalculate_monitor(), libMesh::libmesh_petsc_snes_fd_residual(), libMesh::libmesh_petsc_snes_jacobian(), libMesh::libmesh_petsc_snes_mffd_interface(), libMesh::libmesh_petsc_snes_mffd_residual(), libMesh::libmesh_petsc_snes_postcheck(), libMesh::libmesh_petsc_snes_precheck(), libMesh::libmesh_petsc_snes_residual(), libMesh::libmesh_petsc_snes_residual_helper(), libMesh::MeshRefinement::limit_level_mismatch_at_edge(), libMesh::MeshRefinement::limit_level_mismatch_at_node(), libMesh::MeshRefinement::limit_overrefined_boundary(), libMesh::MeshRefinement::limit_underrefined_boundary(), libMesh::LinearImplicitSystem::LinearImplicitSystem(), main(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshCommunication::make_elems_parallel_consistent(), libMesh::MeshRefinement::make_flags_parallel_consistent(), libMesh::MeshCommunication::make_new_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_new_nodes_parallel_consistent(), libMesh::MeshCommunication::make_node_bcids_parallel_consistent(), libMesh::MeshCommunication::make_node_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_unique_ids_parallel_consistent(), libMesh::MeshCommunication::make_nodes_parallel_consistent(), libMesh::MeshCommunication::make_p_levels_parallel_consistent(), libMesh::MeshRefinement::make_refinement_compatible(), libMesh::TransientRBConstruction::mass_matrix_scaled_matvec(), libMesh::FEMSystem::mesh_position_set(), libMesh::TriangulatorInterface::MeshedHole::MeshedHole(), LinearElasticityWithContact::move_mesh(), libMesh::DistributedMesh::n_active_elem(), libMesh::MeshTools::n_active_levels(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::DofMap::n_constrained_dofs(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::CondensedEigenSystem::n_global_non_condensed_dofs(), libMesh::MeshTools::n_levels(), MixedOrderTest::n_neighbor_links(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::SparsityPattern::Build::n_nonzeros(), libMesh::MeshTools::n_p_levels(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::RBEIMEvaluation::node_distribute_bfs(), libMesh::RBEIMEvaluation::node_gather_bfs(), libMesh::RBEIMConstruction::node_inner_product(), libMesh::MeshBase::operator==(), libMesh::DistributedMesh::parallel_max_elem_id(), libMesh::DistributedMesh::parallel_max_node_id(), libMesh::ReplicatedMesh::parallel_max_unique_id(), libMesh::DistributedMesh::parallel_max_unique_id(), libMesh::DistributedMesh::parallel_n_elem(), libMesh::DistributedMesh::parallel_n_nodes(), libMesh::SparsityPattern::Build::parallel_sync(), libMesh::BoundaryInfo::parallel_sync_node_ids(), libMesh::BoundaryInfo::parallel_sync_side_ids(), libMesh::MeshTools::paranoid_n_levels(), libMesh::Partitioner::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::petsc_auto_fieldsplit(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshBase::prepare_for_use(), libMesh::DofMap::print_dof_constraints(), libMesh::DofMap::process_mesh_constraint_rows(), libMesh::Partitioner::processor_pairs_to_interface_nodes(), libMesh::InterMeshProjection::project_system_vectors(), FEMParameters::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::EquationSystems::read(), libMesh::ExodusII_IO::read_header(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::System::read_header(), libMesh::RBEIMEvaluation::read_in_interior_basis_functions(), libMesh::RBEIMEvaluation::read_in_node_basis_functions(), libMesh::RBEIMEvaluation::read_in_side_basis_functions(), libMesh::RBEvaluation::read_in_vectors_from_multiple_files(), libMesh::System::read_legacy_data(), libMesh::TransientRBConstruction::read_riesz_representors_from_files(), libMesh::RBConstruction::read_riesz_representors_from_files(), libMesh::System::read_SCALAR_dofs(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::System::read_serialized_vector(), libMesh::Nemesis_IO_Helper::read_var_names_impl(), libMesh::MeshBase::recalculate_n_partitions(), libMesh::MeshRefinement::refine_and_coarsen_elements(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::DistributedMesh::renumber_nodes_and_elements(), LinearElasticityWithContact::residual_and_jacobian(), OverlappingAlgebraicGhostingTest::run_ghosting_test(), OverlappingCouplingGhostingTest::run_sparsity_pattern_test(), scale_mesh_and_plot(), libMesh::DofMap::scatter_constraints(), libMesh::CheckpointIO::select_split_config(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::send_and_insert_dof_values(), libMesh::TransientRBConstruction::set_error_temporal_data(), libMesh::Partitioner::set_interface_node_processor_ids_BFS(), libMesh::Partitioner::set_interface_node_processor_ids_linear(), libMesh::Partitioner::set_interface_node_processor_ids_petscpartitioner(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Partitioner::set_parent_processor_ids(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::PetscDiffSolver::setup_petsc_data(), libMesh::RBEIMEvaluation::side_distribute_bfs(), libMesh::RBEIMEvaluation::side_gather_bfs(), libMesh::RBEIMConstruction::side_inner_product(), libMesh::Partitioner::single_partition(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::split_mesh(), libMesh::RBEIMConstruction::store_eim_solutions_for_training_set(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), ConstraintOperatorTest::test1DCoarseningNewNodes(), ConstraintOperatorTest::test1DCoarseningOperator(), libMesh::MeshRefinement::test_level_one(), MeshfunctionDFEM::test_mesh_function_dfem(), MeshfunctionDFEM::test_mesh_function_dfem_grad(), MeshFunctionTest::test_p_level(), libMesh::MeshRefinement::test_unflagged(), DofMapTest::testBadElemFECombo(), SystemsTest::testBlockRestrictedVarNDofs(), BoundaryInfoTest::testBoundaryOnChildrenErrors(), MeshInputTest::testExodusIGASidesets(), MeshTriangulationTest::testFoundCenters(), PointLocatorTest::testLocator(), BoundaryInfoTest::testMesh(), PointLocatorTest::testPlanar(), MeshTriangulationTest::testPoly2TriRefinementBase(), SystemsTest::testProjectCubeWithMeshFunction(), BoundaryInfoTest::testRenumber(), CheckpointIOTest::testSplitter(), MeshInputTest::testTetgenIO(), MeshTriangulationTest::testTriangulatorInterp(), MeshTriangulationTest::testTriangulatorMeshedHoles(), libMesh::MeshTools::total_weight(), libMesh::RBConstruction::train_reduced_basis_with_POD(), libMesh::MeshFunctionSolutionTransfer::transfer(), libMesh::MeshfreeSolutionTransfer::transfer(), libMesh::Poly2TriTriangulator::triangulate(), libMesh::TransientRBConstruction::truth_assembly(), libMesh::RBConstruction::truth_assembly(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::TransientRBConstruction::update_RB_initial_condition_all_N(), libMesh::TransientRBConstruction::update_RB_system_matrices(), libMesh::RBConstruction::update_RB_system_matrices(), libMesh::TransientRBConstruction::update_residual_terms(), libMesh::RBConstruction::update_residual_terms(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::VTKIO::write_nodal_data(), libMesh::RBEIMEvaluation::write_out_interior_basis_functions(), libMesh::RBEIMEvaluation::write_out_node_basis_functions(), libMesh::RBEIMEvaluation::write_out_side_basis_functions(), libMesh::RBEvaluation::write_out_vectors(), libMesh::TransientRBConstruction::write_riesz_representors_to_files(), libMesh::RBConstruction::write_riesz_representors_to_files(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::RBDataSerialization::RBEvaluationSerialization::write_to_file(), libMesh::RBDataSerialization::TransientRBEvaluationSerialization::write_to_file(), libMesh::RBDataSerialization::RBEIMEvaluationSerialization::write_to_file(), and libMesh::RBDataSerialization::RBSCMEvaluationSerialization::write_to_file().

98  { return _communicator; }
const Parallel::Communicator & _communicator

◆ computing_base_vector()

template<typename T>
bool libMesh::PetscNonlinearSolver< T >::computing_base_vector ( ) const
inline
Returns
whether we are computing the base vector for matrix-free finite-differencing

Definition at line 188 of file petsc_nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_mffd_interface(), and libMesh::PetscNonlinearSolver< Number >::set_computing_base_vector().

188 { return _computing_base_vector; }
bool _computing_base_vector
Whether we are computing the base vector for matrix-free finite differencing.

◆ disable_print_counter_info()

void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited

Definition at line 100 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

101 {
102  _enable_print_counter = false;
103  return;
104 }
static bool _enable_print_counter
Flag to control whether reference count information is printed when print_info is called...

◆ enable_print_counter_info()

void libMesh::ReferenceCounter::enable_print_counter_info ( )
staticinherited

Methods to enable/disable the reference counter output from print_info()

Definition at line 94 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

95 {
96  _enable_print_counter = true;
97  return;
98 }
static bool _enable_print_counter
Flag to control whether reference count information is printed when print_info is called...

◆ exact_constraint_enforcement()

template<typename T>
bool libMesh::NonlinearSolver< T >::exact_constraint_enforcement ( )
inlineinherited

Definition at line 403 of file nonlinear_solver.h.

404  {
406  }
bool _exact_constraint_enforcement
Whether we should enforce exact constraints globally during a solve.

◆ force_new_preconditioner()

template<typename T >
void libMesh::PetscNonlinearSolver< T >::force_new_preconditioner ( )
overridevirtual

Immediately force a new preconditioner, even if reuse is set.

Reimplemented from libMesh::NonlinearSolver< T >.

Definition at line 1242 of file petsc_nonlinear_solver.C.

1243 {
1244  // Easiest way is just to clear everything out
1245  this->_is_initialized = false;
1246  _snes.destroy();
1247  _setup_reuse = false;
1248 }
WrappedPetsc< SNES > _snes
Nonlinear solver context.
bool _setup_reuse
Whether we&#39;ve triggered the preconditioner reuse.
bool _is_initialized
Flag indicating if the data structures have been initialized.
void destroy()
Must be specialized to call the appropriate XXXDestroy() routine in order for a WrappedPetsc<T> objec...

◆ get_converged_reason()

template<typename T >
SNESConvergedReason libMesh::PetscNonlinearSolver< T >::get_converged_reason ( )
Returns
The currently-available (or most recently obtained, if the SNES object has been destroyed) convergence reason.

Refer to PETSc docs for the meaning of different SNESConvergedReasons.

Definition at line 1199 of file petsc_nonlinear_solver.C.

1200 {
1201  PetscErrorCode ierr=0;
1202 
1203  if (this->initialized())
1204  {
1205  ierr = SNESGetConvergedReason(_snes, &_reason);
1206  LIBMESH_CHKERR(ierr);
1207  }
1208 
1209  return _reason;
1210 }
WrappedPetsc< SNES > _snes
Nonlinear solver context.
SNESConvergedReason _reason
Store the reason for SNES convergence/divergence for use even after the _snes has been cleared...

◆ get_current_nonlinear_iteration_number()

template<typename T>
virtual unsigned libMesh::PetscNonlinearSolver< T >::get_current_nonlinear_iteration_number ( ) const
inlineoverridevirtual
Returns
The current nonlinear iteration number if called during the solve(), for example by the user-specified residual or Jacobian function.

Implements libMesh::NonlinearSolver< T >.

Definition at line 151 of file petsc_nonlinear_solver.h.

unsigned _current_nonlinear_iteration_number
Stores the current nonlinear iteration number.

◆ get_info()

std::string libMesh::ReferenceCounter::get_info ( )
staticinherited

Gets a string containing the reference information.

Definition at line 47 of file reference_counter.C.

References libMesh::ReferenceCounter::_counts, and libMesh::Quality::name().

Referenced by libMesh::ReferenceCounter::print_info().

48 {
49 #if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
50 
51  std::ostringstream oss;
52 
53  oss << '\n'
54  << " ---------------------------------------------------------------------------- \n"
55  << "| Reference count information |\n"
56  << " ---------------------------------------------------------------------------- \n";
57 
58  for (const auto & [name, cd] : _counts)
59  oss << "| " << name << " reference count information:\n"
60  << "| Creations: " << cd.first << '\n'
61  << "| Destructions: " << cd.second << '\n';
62 
63  oss << " ---------------------------------------------------------------------------- \n";
64 
65  return oss.str();
66 
67 #else
68 
69  return "";
70 
71 #endif
72 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
static Counts _counts
Actually holds the data.

◆ get_total_linear_iterations()

template<typename T >
int libMesh::PetscNonlinearSolver< T >::get_total_linear_iterations ( )
overridevirtual

Get the total number of linear iterations done in the last solve.

Implements libMesh::NonlinearSolver< T >.

Definition at line 1213 of file petsc_nonlinear_solver.C.

1214 {
1215  return _n_linear_iterations;
1216 }
PetscInt _n_linear_iterations
Stores the total number of linear iterations from the last solve.

◆ increment_constructor_count()

void libMesh::ReferenceCounter::increment_constructor_count ( const std::string &  name)
inlineprotectednoexceptinherited

Increments the construction counter.

Should be called in the constructor of any derived class that will be reference counted.

Definition at line 183 of file reference_counter.h.

References libMesh::err, libMesh::BasicOStreamProxy< charT, traits >::get(), libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCountedObject< RBParametrized >::ReferenceCountedObject().

184 {
185  libmesh_try
186  {
187  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
188  std::pair<unsigned int, unsigned int> & p = _counts[name];
189  p.first++;
190  }
191  libmesh_catch (...)
192  {
193  auto stream = libMesh::err.get();
194  stream->exceptions(stream->goodbit); // stream must not throw
195  libMesh::err << "Encountered unrecoverable error while calling "
196  << "ReferenceCounter::increment_constructor_count() "
197  << "for a(n) " << name << " object." << std::endl;
198  std::terminate();
199  }
200 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
OStreamProxy err
static Counts _counts
Actually holds the data.
streamT * get()
Rather than implement every ostream/ios/ios_base function, we&#39;ll be lazy and make esoteric uses go th...
spin_mutex spin_mtx
A convenient spin mutex object which can be used for obtaining locks.
Definition: threads.C:30

◆ increment_destructor_count()

void libMesh::ReferenceCounter::increment_destructor_count ( const std::string &  name)
inlineprotectednoexceptinherited

Increments the destruction counter.

Should be called in the destructor of any derived class that will be reference counted.

Definition at line 207 of file reference_counter.h.

References libMesh::err, libMesh::BasicOStreamProxy< charT, traits >::get(), libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCountedObject< RBParametrized >::~ReferenceCountedObject().

208 {
209  libmesh_try
210  {
211  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
212  std::pair<unsigned int, unsigned int> & p = _counts[name];
213  p.second++;
214  }
215  libmesh_catch (...)
216  {
217  auto stream = libMesh::err.get();
218  stream->exceptions(stream->goodbit); // stream must not throw
219  libMesh::err << "Encountered unrecoverable error while calling "
220  << "ReferenceCounter::increment_destructor_count() "
221  << "for a(n) " << name << " object." << std::endl;
222  std::terminate();
223  }
224 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
OStreamProxy err
static Counts _counts
Actually holds the data.
streamT * get()
Rather than implement every ostream/ios/ios_base function, we&#39;ll be lazy and make esoteric uses go th...
spin_mutex spin_mtx
A convenient spin mutex object which can be used for obtaining locks.
Definition: threads.C:30

◆ init()

template<typename T >
void libMesh::PetscNonlinearSolver< T >::init ( const char *  name = nullptr)
overridevirtual

Initialize data structures if not done so already.

May assign a name to the solver in some implementations

Implements libMesh::NonlinearSolver< T >.

Definition at line 775 of file petsc_nonlinear_solver.C.

776 {
777  parallel_object_only();
778 
779  // Initialize the data structures if not done so already.
780  if (!this->initialized())
781  {
782  this->_is_initialized = true;
783 
784  PetscErrorCode ierr=0;
785 
786  // Make only if we don't already have a retained snes
787  // hanging around from the last solve
788  if (!_snes) {
789  ierr = SNESCreate(this->comm().get(), _snes.get());
790  LIBMESH_CHKERR(ierr);
791  }
792 
793  // I believe all of the following can be safely repeated
794  // even on an old snes instance from the last solve
795 
796  if (name)
797  {
798  ierr = SNESSetOptionsPrefix(_snes, name);
799  LIBMESH_CHKERR(ierr);
800  }
801 
802  // Attaching a DM to SNES.
803  {
804  WrappedPetsc<DM> dm;
805  ierr = DMCreate(this->comm().get(), dm.get()); LIBMESH_CHKERR(ierr);
806  ierr = DMSetType(dm, DMLIBMESH); LIBMESH_CHKERR(ierr);
807  ierr = DMlibMeshSetSystem(dm, this->system()); LIBMESH_CHKERR(ierr);
808  if (name)
809  {
810  ierr = DMSetOptionsPrefix(dm, name); LIBMESH_CHKERR(ierr);
811  }
812  ierr = DMSetFromOptions(dm); LIBMESH_CHKERR(ierr);
813  ierr = DMSetUp(dm); LIBMESH_CHKERR(ierr);
814  ierr = SNESSetDM(_snes, dm); LIBMESH_CHKERR(ierr);
815  // SNES now owns the reference to dm.
816  }
817 
819 
820  // If the SolverConfiguration object is provided, use it to set
821  // options during solver initialization.
822  if (this->_solver_configuration)
823  {
825  }
826 
827  if (this->_preconditioner)
828  {
829  KSP ksp;
830  ierr = SNESGetKSP (_snes, &ksp);
831  LIBMESH_CHKERR(ierr);
832  PC pc;
833  ierr = KSPGetPC(ksp,&pc);
834  LIBMESH_CHKERR(ierr);
835 
836  this->_preconditioner->init();
837 
838  PCSetType(pc, PCSHELL);
839  PCShellSetContext(pc,(void *)this->_preconditioner);
840 
841  //Re-Use the shell functions from petsc_linear_solver
842  PCShellSetSetUp(pc,libmesh_petsc_preconditioner_setup);
843  PCShellSetApply(pc,libmesh_petsc_preconditioner_apply);
844  }
845  }
846 
847 
848  // Tell PETSc about our linesearch "post-check" function, but only
849  // if the user has provided one. There seem to be extra,
850  // unnecessary residual calculations if a postcheck function is
851  // attached for no reason.
852  if (this->postcheck || this->postcheck_object)
853  {
854  SNESLineSearch linesearch;
855  PetscErrorCode ierr = SNESGetLineSearch(_snes, &linesearch);
856  LIBMESH_CHKERR(ierr);
857 
858  ierr = SNESLineSearchSetPostCheck(linesearch, libmesh_petsc_snes_postcheck, this);
859  LIBMESH_CHKERR(ierr);
860  }
861 
862  if (this->precheck_object)
863  {
864  SNESLineSearch linesearch;
865  PetscErrorCode ierr = SNESGetLineSearch(_snes, &linesearch);
866  LIBMESH_CHKERR(ierr);
867 
868  ierr = SNESLineSearchSetPreCheck(linesearch, libmesh_petsc_snes_precheck, this);
869  LIBMESH_CHKERR(ierr);
870  }
871 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
PETSC_EXTERN PetscErrorCode DMlibMeshSetSystem(DM, libMesh::NonlinearImplicitSystem &)
Any functional implementation of the DMlibMesh API must compose the following functions with the DM o...
virtual void set_options_during_init()
Apply options during initialization of a solver.
WrappedPetsc< SNES > _snes
Nonlinear solver context.
Preconditioner< T > * _preconditioner
Holds the Preconditioner object to be used for the linear solves.
const Parallel::Communicator & comm() const
PetscErrorCode libmesh_petsc_snes_postcheck(SNESLineSearch, Vec x, Vec y, Vec w, PetscBool *changed_y, PetscBool *changed_w, void *context)
SolverConfiguration * _solver_configuration
Optionally store a SolverOptions object that can be used to set parameters like solver type...
const sys_type & system() const
bool _is_initialized
Flag indicating if the data structures have been initialized.
PetscErrorCode libmesh_petsc_preconditioner_apply(PC pc, Vec x, Vec y)
This function is called by PETSc to actually apply the preconditioner.
void setup_default_monitor()
Setup the default monitor if required.
communicator & get()
PetscErrorCode libmesh_petsc_preconditioner_setup(PC pc)
This function is called by PETSc to initialize the preconditioner.
NonlinearImplicitSystem::ComputePostCheck * postcheck_object
A callable object that is executed after each nonlinear iteration.
friend PetscErrorCode libmesh_petsc_snes_precheck(SNESLineSearch, Vec X, Vec Y, PetscBool *changed, void *context)
void(* postcheck)(const NumericVector< Number > &old_soln, NumericVector< Number > &search_direction, NumericVector< Number > &new_soln, bool &changed_search_direction, bool &changed_new_soln, sys_type &S)
Function that performs a "check" on the Newton search direction and solution after each nonlinear ste...
NonlinearImplicitSystem::ComputePreCheck * precheck_object

◆ initialized()

template<typename T>
bool libMesh::NonlinearSolver< T >::initialized ( ) const
inlineinherited
Returns
true if the data structures are initialized, false otherwise.

Definition at line 83 of file nonlinear_solver.h.

83 { return _is_initialized; }
bool _is_initialized
Flag indicating if the data structures have been initialized.

◆ n_objects()

static unsigned int libMesh::ReferenceCounter::n_objects ( )
inlinestaticinherited

Prints the number of outstanding (created, but not yet destroyed) objects.

Definition at line 85 of file reference_counter.h.

References libMesh::ReferenceCounter::_n_objects.

Referenced by libMesh::LibMeshInit::~LibMeshInit().

86  { return _n_objects; }
static Threads::atomic< unsigned int > _n_objects
The number of objects.

◆ n_processors()

processor_id_type libMesh::ParallelObject::n_processors ( ) const
inlineinherited
Returns
The number of processors in the group.

Definition at line 103 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, libMesh::libmesh_assert(), and TIMPI::Communicator::size().

Referenced by libMesh::Partitioner::_find_global_index_by_pid_map(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::DofMap::add_constraints_to_send_list(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::DistributedMesh::add_elem(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::DistributedMesh::add_node(), libMesh::System::add_vector(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::DofMap::allgather_recursive_constraints(), libMesh::FEMSystem::assembly(), libMesh::Nemesis_IO::assert_symmetric_cmaps(), libMesh::Partitioner::assign_partitioning(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::DistributedMesh::clear(), libMesh::DistributedMesh::clear_elems(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::ExodusII_IO::copy_scalar_solution(), libMesh::Nemesis_IO::copy_scalar_solution(), libMesh::UnstructuredMesh::create_pid_mesh(), libMesh::MeshTools::create_processor_bounding_box(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_scalar_dofs(), libMesh::DistributedMesh::DistributedMesh(), libMesh::EnsightIO::EnsightIO(), libMesh::RBEIMEvaluation::gather_bfs(), libMesh::MeshBase::get_info(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::ExodusII_IO_Helper::initialize(), libMesh::DistributedMesh::insert_elem(), libMesh::MeshTools::libmesh_assert_contiguous_dof_ids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_new_node_procids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_topology_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_boundary_ids(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::DofMap::local_variable_indices(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshBase::n_active_elem_on_proc(), libMesh::MeshBase::n_elem_on_proc(), libMesh::MeshBase::n_nodes_on_proc(), libMesh::RBEIMEvaluation::node_gather_bfs(), libMesh::Partitioner::partition(), libMesh::MeshBase::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::DofMap::prepare_send_list(), libMesh::DofMap::print_dof_constraints(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::CheckpointIO::read_connectivity(), libMesh::XdrIO::read_header(), libMesh::CheckpointIO::read_nodes(), libMesh::System::read_parallel_data(), libMesh::System::read_SCALAR_dofs(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::System::read_serialized_vector(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::Partitioner::repartition(), OverlappingFunctorTest::run_partitioner_test(), libMesh::DofMap::scatter_constraints(), libMesh::DistributedMesh::set_next_unique_id(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), WriteVecAndScalar::setupTests(), libMesh::RBEIMEvaluation::side_gather_bfs(), DistributedMeshTest::testRemoteElemError(), CheckpointIOTest::testSplitter(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO_Helper::write_nodal_coordinates(), libMesh::VTKIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::System::write_parallel_data(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::XdrIO::write_serialized_nodesets().

104  {
105  processor_id_type returnval =
106  cast_int<processor_id_type>(_communicator.size());
107  libmesh_assert(returnval); // We never have an empty comm
108  return returnval;
109  }
const Parallel::Communicator & _communicator
processor_id_type size() const
uint8_t processor_id_type
libmesh_assert(ctx)

◆ print_converged_reason()

template<typename T >
void libMesh::PetscNonlinearSolver< T >::print_converged_reason ( )
overridevirtual

Prints a useful message about why the latest nonlinear solve con(di)verged.

Reimplemented from libMesh::NonlinearSolver< T >.

Definition at line 1189 of file petsc_nonlinear_solver.C.

1190 {
1191 
1192  libMesh::out << "Nonlinear solver convergence/divergence reason: "
1193  << SNESConvergedReasons[this->get_converged_reason()] << std::endl;
1194 }
SNESConvergedReason get_converged_reason()
OStreamProxy out

◆ print_info()

void libMesh::ReferenceCounter::print_info ( std::ostream &  out_stream = libMesh::out)
staticinherited

Prints the reference information, by default to libMesh::out.

Definition at line 81 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter, and libMesh::ReferenceCounter::get_info().

Referenced by libMesh::LibMeshInit::~LibMeshInit().

82 {
84  out_stream << ReferenceCounter::get_info();
85 }
static std::string get_info()
Gets a string containing the reference information.
static bool _enable_print_counter
Flag to control whether reference count information is printed when print_info is called...

◆ processor_id()

processor_id_type libMesh::ParallelObject::processor_id ( ) const
inlineinherited
Returns
The rank of this processor in the group.

Definition at line 114 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, and TIMPI::Communicator::rank().

Referenced by libMesh::BoundaryInfo::_find_id_maps(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::DistributedMesh::add_elem(), libMesh::BoundaryInfo::add_elements(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::DistributedMesh::add_node(), libMesh::MeshTools::Modification::all_tri(), libMesh::DofMap::allgather_recursive_constraints(), libMesh::FEMSystem::assembly(), libMesh::Nemesis_IO::assert_symmetric_cmaps(), libMesh::Partitioner::assign_partitioning(), libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::Partitioner::build_graph(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::EquationSystems::build_parallel_solution_vector(), libMesh::DistributedMesh::clear(), libMesh::DistributedMesh::clear_elems(), libMesh::ExodusII_IO_Helper::close(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::compute_communication_map_parameters(), libMesh::Nemesis_IO_Helper::compute_internal_and_border_elems_and_internal_nodes(), libMesh::RBConstruction::compute_max_error_bound(), libMesh::Nemesis_IO_Helper::compute_node_communication_maps(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::ExodusII_IO::copy_elemental_solution(), libMesh::ExodusII_IO::copy_nodal_solution(), libMesh::ExodusII_IO::copy_scalar_solution(), libMesh::Nemesis_IO::copy_scalar_solution(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::ExodusII_IO_Helper::create(), libMesh::DistributedMesh::delete_elem(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DofMap::distribute_scalar_dofs(), libMesh::DistributedMesh::DistributedMesh(), libMesh::DofMap::end_dof(), libMesh::DofMap::end_old_dof(), libMesh::EnsightIO::EnsightIO(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::find_dofs_to_send(), libMesh::MeshFunction::find_element(), libMesh::MeshFunction::find_elements(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::DofMap::first_dof(), libMesh::DofMap::first_old_dof(), libMesh::RBEIMEvaluation::gather_bfs(), libMesh::Nemesis_IO_Helper::get_cmap_params(), libMesh::Nemesis_IO_Helper::get_eb_info_global(), libMesh::Nemesis_IO_Helper::get_elem_cmap(), libMesh::Nemesis_IO_Helper::get_elem_map(), libMesh::MeshBase::get_info(), libMesh::DofMap::get_info(), libMesh::Nemesis_IO_Helper::get_init_global(), libMesh::Nemesis_IO_Helper::get_init_info(), libMesh::Nemesis_IO_Helper::get_loadbal_param(), libMesh::DofMap::get_local_constraints(), libMesh::Nemesis_IO_Helper::get_node_cmap(), libMesh::Nemesis_IO_Helper::get_node_map(), libMesh::Nemesis_IO_Helper::get_ns_param_global(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), libMesh::SparsityPattern::Build::handle_vi_vj(), libMesh::LaplaceMeshSmoother::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), HeatSystem::init_data(), libMesh::ExodusII_IO_Helper::initialize(), libMesh::ExodusII_IO_Helper::initialize_element_variables(), libMesh::ExodusII_IO_Helper::initialize_global_variables(), libMesh::ExodusII_IO_Helper::initialize_nodal_variables(), libMesh::DistributedMesh::insert_elem(), libMesh::DofMap::is_evaluable(), libMesh::SparsityPattern::Build::join(), libMesh::TransientRBEvaluation::legacy_write_offline_data_to_files(), libMesh::RBSCMEvaluation::legacy_write_offline_data_to_files(), libMesh::RBEvaluation::legacy_write_offline_data_to_files(), libMesh::MeshTools::libmesh_assert_consistent_distributed(), libMesh::MeshTools::libmesh_assert_consistent_distributed_nodes(), libMesh::MeshTools::libmesh_assert_contiguous_dof_ids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DofMap::local_variable_indices(), main(), libMesh::MeshRefinement::make_coarsening_compatible(), AugmentSparsityOnInterface::mesh_reinit(), libMesh::TriangulatorInterface::MeshedHole::MeshedHole(), libMesh::MeshBase::n_active_local_elem(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::DofMap::n_local_dofs(), libMesh::System::n_local_dofs(), libMesh::MeshBase::n_local_elem(), libMesh::MeshBase::n_local_nodes(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::RBEIMEvaluation::node_gather_bfs(), libMesh::SparsityPattern::Build::operator()(), libMesh::DistributedMesh::own_node(), libMesh::BoundaryInfo::parallel_sync_node_ids(), libMesh::BoundaryInfo::parallel_sync_side_ids(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::DofMap::print_dof_constraints(), libMesh::DofMap::process_mesh_constraint_rows(), libMesh::Nemesis_IO_Helper::put_cmap_params(), libMesh::Nemesis_IO_Helper::put_elem_cmap(), libMesh::Nemesis_IO_Helper::put_elem_map(), libMesh::Nemesis_IO_Helper::put_loadbal_param(), libMesh::Nemesis_IO_Helper::put_node_cmap(), libMesh::Nemesis_IO_Helper::put_node_map(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read(), libMesh::EquationSystems::read(), libMesh::ExodusII_IO_Helper::read_elem_num_map(), libMesh::ExodusII_IO_Helper::read_global_values(), libMesh::ExodusII_IO::read_header(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), libMesh::DynaIO::read_mesh(), libMesh::ExodusII_IO_Helper::read_node_num_map(), libMesh::System::read_parallel_data(), libMesh::TransientRBConstruction::read_riesz_representors_from_files(), libMesh::RBConstruction::read_riesz_representors_from_files(), libMesh::System::read_SCALAR_dofs(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::System::read_serialized_data(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::System::read_serialized_vector(), libMesh::System::read_serialized_vectors(), libMesh::Nemesis_IO_Helper::read_var_names_impl(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::DistributedMesh::renumber_nodes_and_elements(), libMesh::DofMap::scatter_constraints(), libMesh::CheckpointIO::select_split_config(), libMesh::DistributedMesh::set_next_unique_id(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::RBEIMEvaluation::side_gather_bfs(), ExodusTest< elem_type >::test_read_gold(), ExodusTest< elem_type >::test_write(), MeshInputTest::testAbaqusRead(), MeshInputTest::testCopyElementSolutionImpl(), MeshInputTest::testCopyElementVectorImpl(), MeshInputTest::testCopyNodalSolutionImpl(), DefaultCouplingTest::testCoupling(), PointNeighborCouplingTest::testCoupling(), MeshInputTest::testDynaFileMappings(), MeshInputTest::testDynaNoSplines(), MeshInputTest::testDynaReadElem(), MeshInputTest::testDynaReadPatch(), MeshInputTest::testExodusFileMappings(), MeshInputTest::testExodusIGASidesets(), MeshInputTest::testExodusWriteElementDataFromDiscontinuousNodalData(), MeshInputTest::testLowOrderEdgeBlocks(), SystemsTest::testProjectMatrix1D(), SystemsTest::testProjectMatrix2D(), SystemsTest::testProjectMatrix3D(), BoundaryInfoTest::testShellFaceConstraints(), MeshInputTest::testSingleElementImpl(), WriteVecAndScalar::testSolution(), CheckpointIOTest::testSplitter(), MeshInputTest::testTetgenIO(), libMesh::MeshTools::total_weight(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::DTKAdapter::update_variable_values(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::CheckpointIO::write(), libMesh::EquationSystems::write(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), libMesh::ExodusII_IO_Helper::write_element_values(), libMesh::ExodusII_IO_Helper::write_element_values_element_major(), libMesh::ExodusII_IO_Helper::write_elements(), libMesh::ExodusII_IO_Helper::write_elemset_data(), libMesh::ExodusII_IO_Helper::write_elemsets(), libMesh::ExodusII_IO::write_global_data(), libMesh::ExodusII_IO_Helper::write_global_values(), libMesh::System::write_header(), libMesh::ExodusII_IO::write_information_records(), libMesh::ExodusII_IO_Helper::write_information_records(), libMesh::ExodusII_IO_Helper::write_nodal_coordinates(), libMesh::UCDIO::write_nodal_data(), libMesh::VTKIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_common(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::ExodusII_IO_Helper::write_nodal_values(), libMesh::ExodusII_IO_Helper::write_nodeset_data(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::ExodusII_IO_Helper::write_nodesets(), libMesh::RBEIMEvaluation::write_out_interior_basis_functions(), libMesh::RBEIMEvaluation::write_out_node_basis_functions(), libMesh::RBEIMEvaluation::write_out_side_basis_functions(), write_output_solvedata(), libMesh::System::write_parallel_data(), libMesh::RBConstruction::write_riesz_representors_to_files(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bc_names(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::System::write_serialized_data(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::XdrIO::write_serialized_subdomain_names(), libMesh::System::write_serialized_vector(), libMesh::System::write_serialized_vectors(), libMesh::ExodusII_IO_Helper::write_sideset_data(), libMesh::Nemesis_IO_Helper::write_sidesets(), libMesh::ExodusII_IO_Helper::write_sidesets(), libMesh::ExodusII_IO::write_timestep(), libMesh::ExodusII_IO_Helper::write_timestep(), and libMesh::ExodusII_IO::write_timestep_discontinuous().

115  { return cast_int<processor_id_type>(_communicator.rank()); }
processor_id_type rank() const
const Parallel::Communicator & _communicator

◆ reuse_preconditioner()

template<typename T >
bool libMesh::PetscNonlinearSolver< T >::reuse_preconditioner ( ) const
overridevirtual

Getter for preconditioner reuse.

Reimplemented from libMesh::NonlinearSolver< T >.

Definition at line 1230 of file petsc_nonlinear_solver.C.

1231 {
1232  return this->_reuse_preconditioner;
1233 }
bool _reuse_preconditioner
Whether we should reuse the linear preconditioner.

◆ reuse_preconditioner_max_linear_its()

template<typename T >
unsigned int libMesh::PetscNonlinearSolver< T >::reuse_preconditioner_max_linear_its ( ) const
overridevirtual

Getter for the maximum iterations flag for preconditioner reuse.

Reimplemented from libMesh::NonlinearSolver< T >.

Definition at line 1236 of file petsc_nonlinear_solver.C.

Referenced by libMesh::libmesh_petsc_recalculate_monitor().

1237 {
1239 }
unsigned int _reuse_preconditioner_max_linear_its
Number of linear iterations to retain the preconditioner.

◆ set_computing_base_vector()

template<typename T>
void libMesh::PetscNonlinearSolver< T >::set_computing_base_vector ( bool  computing_base_vector)
inline

Set whether we are computing the base vector for matrix-free finite-differencing.

Definition at line 183 of file petsc_nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_mffd_interface().

bool _computing_base_vector
Whether we are computing the base vector for matrix-free finite differencing.

◆ set_exact_constraint_enforcement()

template<typename T>
virtual void libMesh::NonlinearSolver< T >::set_exact_constraint_enforcement ( bool  enable)
inlinevirtualinherited

Enable (or disable; it is true by default) exact enforcement of constraints at the solver level, correcting any constrained DoF coefficients in current_local_solution as well as applying nonlinear residual and Jacobian terms based on constraint equations.

This is probably only safe to disable if user code is setting nonlinear residual and Jacobian terms based on constraint equations at an element-by-element level, by combining the asymmetric_constraint_rows option with the residual_constrain_element_vector processing option in DofMap.

Definition at line 398 of file nonlinear_solver.h.

399  {
401  }
bool _exact_constraint_enforcement
Whether we should enforce exact constraints globally during a solve.

◆ set_jacobian_zero_out()

template<typename T>
void libMesh::PetscNonlinearSolver< T >::set_jacobian_zero_out ( bool  state)
inline

Set if the jacobian should be zeroed out in the callback.

Definition at line 162 of file petsc_nonlinear_solver.h.

162 { _zero_out_jacobian = state; }
bool _zero_out_jacobian
true to zero out jacobian before going into application level call-back, otherwise false ...

◆ set_residual_zero_out()

template<typename T>
void libMesh::PetscNonlinearSolver< T >::set_residual_zero_out ( bool  state)
inline

Set if the residual should be zeroed out in the callback.

Definition at line 157 of file petsc_nonlinear_solver.h.

157 { _zero_out_residual = state; }
bool _zero_out_residual
true to zero out residual before going into application level call-back, otherwise false ...

◆ set_reuse_preconditioner()

template<typename T >
void libMesh::NonlinearSolver< T >::set_reuse_preconditioner ( bool  reuse)
virtualinherited

Set the reuse preconditioner flag.

Definition at line 94 of file nonlinear_solver.C.

95 {
96  _reuse_preconditioner = reuse;
97 }
bool _reuse_preconditioner
Whether we should reuse the linear preconditioner.

◆ set_reuse_preconditioner_max_linear_its()

template<typename T >
void libMesh::NonlinearSolver< T >::set_reuse_preconditioner_max_linear_its ( unsigned int  i)
virtualinherited

Set the reuse_preconditioner_max_linear_its parameter.

Definition at line 106 of file nonlinear_solver.C.

107 {
109 }
unsigned int _reuse_preconditioner_max_linear_its
Number of linear iterations to retain the preconditioner.

◆ set_snesmf_reuse_base()

template<typename T>
void libMesh::PetscNonlinearSolver< T >::set_snesmf_reuse_base ( bool  state)
inline

Set to true to let PETSc reuse the base vector.

Definition at line 172 of file petsc_nonlinear_solver.h.

172 { _snesmf_reuse_base = state; }
bool _snesmf_reuse_base
True, If we want the base vector to be used for differencing even if the function provided to SNESSet...

◆ set_solver_configuration()

template<typename T >
void libMesh::NonlinearSolver< T >::set_solver_configuration ( SolverConfiguration solver_configuration)
inherited

Set the solver configuration object.

Definition at line 82 of file nonlinear_solver.C.

83 {
84  _solver_configuration = &solver_configuration;
85 }
SolverConfiguration * _solver_configuration
Optionally store a SolverOptions object that can be used to set parameters like solver type...

◆ setup_default_monitor()

template<typename T >
void libMesh::PetscNonlinearSolver< T >::setup_default_monitor ( )

Setup the default monitor if required.

Definition at line 1219 of file petsc_nonlinear_solver.C.

1220 {
1221  if (_default_monitor)
1222  {
1223  PetscErrorCode ierr = SNESMonitorSet (_snes, libmesh_petsc_snes_monitor,
1224  this, LIBMESH_PETSC_NULLPTR);
1225  LIBMESH_CHKERR(ierr);
1226  }
1227 }
WrappedPetsc< SNES > _snes
Nonlinear solver context.
PetscErrorCode libmesh_petsc_snes_monitor(SNES, PetscInt its, PetscReal fnorm, void *)
bool _default_monitor
true if we want the default monitor to be set, false for no monitor (i.e.

◆ snes()

template<typename T >
SNES libMesh::PetscNonlinearSolver< T >::snes ( )
Returns
The raw PETSc snes context pointer.

Definition at line 875 of file petsc_nonlinear_solver.C.

Referenced by libMesh::libmesh_petsc_snes_mffd_interface().

876 {
877  this->init();
878  return _snes;
879 }
WrappedPetsc< SNES > _snes
Nonlinear solver context.
virtual void init(const char *name=nullptr) override
Initialize data structures if not done so already.

◆ snes_mf_reuse_base()

template<typename T>
bool libMesh::PetscNonlinearSolver< T >::snes_mf_reuse_base ( ) const
inline
Returns
Whether we are reusing the nonlinear function evaluation as the base for doing matrix-free approximation of the Jacobian action

Definition at line 178 of file petsc_nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_mffd_interface().

178 { return _snesmf_reuse_base; }
bool _snesmf_reuse_base
True, If we want the base vector to be used for differencing even if the function provided to SNESSet...

◆ solve()

template<typename T>
std::pair< unsigned int, Real > libMesh::PetscNonlinearSolver< T >::solve ( SparseMatrix< T > &  pre_in,
NumericVector< T > &  x_in,
NumericVector< T > &  r_in,
const double  ,
const unsigned  int 
)
overridevirtual

Call the Petsc solver.

It calls the method below, using the same matrix for the system and preconditioner matrices.

Implements libMesh::NonlinearSolver< T >.

Definition at line 950 of file petsc_nonlinear_solver.C.

955 {
956  parallel_object_only();
957 
958  LOG_SCOPE("solve()", "PetscNonlinearSolver");
959  this->init ();
960 
961  // Make sure the data passed in are really of Petsc types
962  PetscMatrix<T> * pre = cast_ptr<PetscMatrix<T> *>(&pre_in);
963  PetscVector<T> * x = cast_ptr<PetscVector<T> *>(&x_in);
964  PetscVector<T> * r = cast_ptr<PetscVector<T> *>(&r_in);
965 
966  PetscErrorCode ierr=0;
967  PetscInt n_iterations =0;
968  // Should actually be a PetscReal, but I don't know which version of PETSc first introduced PetscReal
969  Real final_residual_norm=0.;
970 
971  // We don't want to do this twice because it resets
972  // SNESSetLagPreconditioner
973  if ((reuse_preconditioner()) && (!_setup_reuse))
974  {
975  _setup_reuse = true;
976  ierr = SNESSetLagPreconditionerPersists(_snes, PETSC_TRUE);
977  LIBMESH_CHKERR(ierr);
978  // According to the PETSC 3.16.5 docs -2 is a magic number which
979  // means "recalculate the next time you need it and then not again"
980  ierr = SNESSetLagPreconditioner(_snes, -2);
981  LIBMESH_CHKERR(ierr);
982  // Add in our callback which will trigger recalculating
983  // the preconditioner when we hit reuse_preconditioner_max_linear_its
984  ierr = SNESMonitorSet(_snes, &libmesh_petsc_recalculate_monitor,
985  this,
986  NULL);
987  LIBMESH_CHKERR(ierr);
988  }
989  else if (!(reuse_preconditioner()))
990  // This covers the case where it was enabled but was then disabled
991  {
992  ierr = SNESSetLagPreconditionerPersists(_snes, PETSC_FALSE);
993  LIBMESH_CHKERR(ierr);
994  if (_setup_reuse)
995  {
996  _setup_reuse = false;
997  ierr = SNESMonitorCancel(_snes);
998  LIBMESH_CHKERR(ierr);
999  // Readd default monitor
1001  }
1002  }
1003 
1004  ierr = SNESSetFunction (_snes, r->vec(), libmesh_petsc_snes_residual, this);
1005  LIBMESH_CHKERR(ierr);
1006 
1007  // Only set the jacobian function if we've been provided with something to call.
1008  // This allows a user to set their own jacobian function if they want to
1009  if (this->jacobian || this->jacobian_object || this->residual_and_jacobian_object)
1010  {
1011  ierr = SNESSetJacobian (_snes, pre->mat(), pre->mat(), libmesh_petsc_snes_jacobian, this);
1012  LIBMESH_CHKERR(ierr);
1013  }
1014 
1015 
1016  // Only set the nullspace if we have a way of computing it and the result is non-empty.
1017  if (this->nullspace || this->nullspace_object)
1018  {
1019  WrappedPetsc<MatNullSpace> msp;
1020  this->build_mat_null_space(this->nullspace_object, this->nullspace, msp.get());
1021  if (msp)
1022  {
1023  ierr = MatSetNullSpace(pre->mat(), msp);
1024  LIBMESH_CHKERR(ierr);
1025  }
1026  }
1027 
1028  // Only set the transpose nullspace if we have a way of computing it and the result is non-empty.
1030  {
1031 #if PETSC_VERSION_LESS_THAN(3,6,0)
1032  libmesh_warning("MatSetTransposeNullSpace is only supported for PETSc >= 3.6, transpose nullspace will be ignored.");
1033 #else
1034  WrappedPetsc<MatNullSpace> msp;
1036  if (msp)
1037  {
1038  ierr = MatSetTransposeNullSpace(pre->mat(), msp);
1039  LIBMESH_CHKERR(ierr);
1040  }
1041 #endif
1042  }
1043 
1044  // Only set the nearnullspace if we have a way of computing it and the result is non-empty.
1045  if (this->nearnullspace || this->nearnullspace_object)
1046  {
1047  WrappedPetsc<MatNullSpace> msp;
1048  this->build_mat_null_space(this->nearnullspace_object, this->nearnullspace, msp.get());
1049 
1050  if (msp)
1051  {
1052  ierr = MatSetNearNullSpace(pre->mat(), msp);
1053  LIBMESH_CHKERR(ierr);
1054  }
1055  }
1056 
1057  // Have the Krylov subspace method use our good initial guess rather than 0
1058  KSP ksp;
1059  ierr = SNESGetKSP (_snes, &ksp);
1060  LIBMESH_CHKERR(ierr);
1061 
1062  // Set the tolerances for the iterative solver. Use the user-supplied
1063  // tolerance for the relative residual & leave the others at default values
1064  ierr = KSPSetTolerances (ksp, this->initial_linear_tolerance, PETSC_DEFAULT,
1065  PETSC_DEFAULT, this->max_linear_iterations);
1066  LIBMESH_CHKERR(ierr);
1067 
1068  // Set the tolerances for the non-linear solver.
1069  ierr = SNESSetTolerances(_snes, this->absolute_residual_tolerance, this->relative_residual_tolerance,
1071  LIBMESH_CHKERR(ierr);
1072 
1073  // Set the divergence tolerance for the non-linear solver
1074 #if !PETSC_VERSION_LESS_THAN(3,8,0)
1075  ierr = SNESSetDivergenceTolerance(_snes, this->divergence_tolerance);
1076  LIBMESH_CHKERR(ierr);
1077 #endif
1078 
1079  //Pull in command-line options
1080 #if PETSC_VERSION_LESS_THAN(3,7,0)
1081  KSPSetFromOptions(ksp);
1082 #endif
1083  SNESSetFromOptions(_snes);
1084 
1085  if (this->user_presolve)
1086  this->user_presolve(this->system());
1087 
1088  //Set the preconditioning matrix
1089  if (this->_preconditioner)
1090  {
1091  this->_preconditioner->set_matrix(pre_in);
1092  this->_preconditioner->init();
1093  }
1094 
1095  // If the SolverConfiguration object is provided, use it to override
1096  // solver options.
1097  if (this->_solver_configuration)
1098  {
1100  }
1101 
1102  // In PETSc versions before 3.5.0, it is not possible to call
1103  // SNESSetUp() before the solution and rhs vectors are initialized, as
1104  // this triggers the
1105  //
1106  // "Solution vector cannot be right hand side vector!"
1107  //
1108  // error message. It is also not possible to call SNESSetSolution()
1109  // in those versions of PETSc to work around the problem, since that
1110  // API was removed in 3.0.0 and only restored in 3.6.0. The
1111  // overzealous check was moved out of SNESSetUp in PETSc 3.5.0
1112  // (petsc/petsc@154060b), so this code block should be safe to use
1113  // in 3.5.0 and later.
1114 #if !PETSC_VERSION_LESS_THAN(3,6,0)
1115  ierr = SNESSetSolution(_snes, x->vec());
1116  LIBMESH_CHKERR(ierr);
1117 #endif
1118  ierr = SNESSetUp(_snes);
1119  LIBMESH_CHKERR(ierr);
1120 
1121  Mat J;
1122  ierr = SNESGetJacobian(_snes, &J, LIBMESH_PETSC_NULLPTR,
1123  LIBMESH_PETSC_NULLPTR,
1124  LIBMESH_PETSC_NULLPTR);
1125  LIBMESH_CHKERR(ierr);
1126  ierr = MatMFFDSetFunction(J, libmesh_petsc_snes_mffd_interface, this);
1127  LIBMESH_CHKERR(ierr);
1128 #if !PETSC_VERSION_LESS_THAN(3,8,4)
1129 #ifndef NDEBUG
1130  // If we're in debug mode, do not reuse the nonlinear function evaluation as the base for doing
1131  // matrix-free approximations of the Jacobian action. Instead if the user requested that we reuse
1132  // the base, we'll check the base function evaluation and compare it to the nonlinear residual
1133  // evaluation. If they are different, then we'll error and inform the user that it's unsafe to
1134  // reuse the base
1135  ierr = MatSNESMFSetReuseBase(J, PETSC_FALSE);
1136 #else
1137  // Resue the residual vector from SNES
1138  ierr = MatSNESMFSetReuseBase(J, static_cast<PetscBool>(_snesmf_reuse_base));
1139 #endif
1140  LIBMESH_CHKERR(ierr);
1141 #endif
1142 
1143  SNESLineSearch linesearch;
1144  if (linesearch_object)
1145  {
1146  ierr = SNESGetLineSearch(_snes, &linesearch);
1147  LIBMESH_CHKERR(ierr);
1148  ierr = SNESLineSearchSetType(linesearch, SNESLINESEARCHSHELL);
1149  LIBMESH_CHKERR(ierr);
1150  ierr = SNESLineSearchShellSetUserFunc(linesearch, libmesh_petsc_linesearch_shellfunc, this);
1151  LIBMESH_CHKERR(ierr);
1152  }
1153 
1154  ierr = SNESSolve (_snes, LIBMESH_PETSC_NULLPTR, x->vec());
1155  LIBMESH_CHKERR(ierr);
1156 
1157  ierr = SNESGetIterationNumber(_snes, &n_iterations);
1158  LIBMESH_CHKERR(ierr);
1159 
1160  ierr = SNESGetLinearSolveIterations(_snes, &_n_linear_iterations);
1161  LIBMESH_CHKERR(ierr);
1162 
1163  // SNESGetFunction has been around forever and should work on all
1164  // versions of PETSc. This is also now the recommended approach
1165  // according to the documentation for the PETSc 3.5.1 release:
1166  // http://www.mcs.anl.gov/petsc/documentation/changes/35.html
1167  Vec f;
1168  ierr = SNESGetFunction(_snes, &f, 0, 0);
1169  LIBMESH_CHKERR(ierr);
1170  ierr = VecNorm(f, NORM_2, pPR(&final_residual_norm));
1171  LIBMESH_CHKERR(ierr);
1172 
1173  // Get and store the reason for convergence
1174  SNESGetConvergedReason(_snes, &_reason);
1175 
1176  //Based on Petsc 2.3.3 documentation all diverged reasons are negative
1177  this->converged = (_reason >= 0);
1178 
1179  // Reset data structure
1180  this->clear();
1181 
1182  // return the # of its. and the final residual norm.
1183  return std::make_pair(n_iterations, final_residual_norm);
1184 }
friend PetscErrorCode libmesh_petsc_snes_residual(SNES snes, Vec x, Vec r, void *ctx)
PetscReal * pPR(T *ptr)
Definition: petsc_macro.h:184
NonlinearImplicitSystem::ComputeResidualandJacobian * residual_and_jacobian_object
Object that computes either the residual or the Jacobian of the nonlinear system at the input itera...
bool converged
After a call to solve this will reflect whether or not the nonlinear solve was successful.
WrappedPetsc< SNES > _snes
Nonlinear solver context.
unsigned int max_function_evaluations
Maximum number of function evaluations.
Preconditioner< T > * _preconditioner
Holds the Preconditioner object to be used for the linear solves.
void(* user_presolve)(sys_type &S)
Customizable function pointer which users can attach to the solver.
friend PetscErrorCode libmesh_petsc_snes_jacobian(SNES snes, Vec x, Mat jac, Mat pc, void *ctx)
virtual void clear() override
Release all memory and clear data structures.
virtual void configure_solver()=0
Apply solver options to a particular solver.
unsigned int max_linear_iterations
Each linear solver step should exit after max_linear_iterations is exceeded.
SolverConfiguration * _solver_configuration
Optionally store a SolverOptions object that can be used to set parameters like solver type...
void(* jacobian)(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)
Function that computes the Jacobian J(X) of the nonlinear system at the input iterate X...
bool _setup_reuse
Whether we&#39;ve triggered the preconditioner reuse.
const sys_type & system() const
NonlinearImplicitSystem::ComputeJacobian * jacobian_object
Object that computes the Jacobian J(X) of the nonlinear system at the input iterate X...
PetscInt _n_linear_iterations
Stores the total number of linear iterations from the last solve.
void(* transpose_nullspace)(std::vector< NumericVector< Number > *> &sp, sys_type &S)
Function that computes a basis for the transpose Jacobian&#39;s nullspace – when solving a degenerate pr...
void setup_default_monitor()
Setup the default monitor if required.
void(* nullspace)(std::vector< NumericVector< Number > *> &sp, sys_type &S)
Function that computes a basis for the Jacobian&#39;s nullspace – the kernel or the "zero energy modes" ...
void(* nearnullspace)(std::vector< NumericVector< Number > *> &sp, sys_type &S)
Function that computes a basis for the Jacobian&#39;s near nullspace – the set of "low energy modes" – ...
SNESConvergedReason _reason
Store the reason for SNES convergence/divergence for use even after the _snes has been cleared...
PetscErrorCode libmesh_petsc_snes_mffd_interface(void *ctx, Vec x, Vec r)
double divergence_tolerance
The NonlinearSolver should exit if the residual becomes greater than the initial residual times the d...
PetscErrorCode libmesh_petsc_linesearch_shellfunc(SNESLineSearch linesearch, void *ctx)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::unique_ptr< ComputeLineSearchObject > linesearch_object
A callable object that can be used to specify a custom line-search.
void build_mat_null_space(NonlinearImplicitSystem::ComputeVectorSubspace *computeSubspaceObject, void(*)(std::vector< NumericVector< Number > *> &, sys_type &), MatNullSpace *)
virtual void init(const char *name=nullptr) override
Initialize data structures if not done so already.
double absolute_residual_tolerance
The NonlinearSolver should exit after the residual is reduced to either less than absolute_residual_t...
virtual bool reuse_preconditioner() const override
Getter for preconditioner reuse.
NonlinearImplicitSystem::ComputeVectorSubspace * nearnullspace_object
A callable object that computes a basis for the Jacobian&#39;s near nullspace – the set of "low energy m...
double initial_linear_tolerance
Any required linear solves will at first be done with this tolerance; the NonlinearSolver may tighten...
unsigned int max_nonlinear_iterations
Maximum number of non-linear iterations.
NonlinearImplicitSystem::ComputeVectorSubspace * transpose_nullspace_object
A callable object that computes a basis for the transpose Jacobian&#39;s nullspace – when solving a dege...
bool _snesmf_reuse_base
True, If we want the base vector to be used for differencing even if the function provided to SNESSet...
PetscErrorCode libmesh_petsc_recalculate_monitor(SNES snes, PetscInt it, PetscReal norm, void *mctx)
NonlinearImplicitSystem::ComputeVectorSubspace * nullspace_object
A callable object that computes a basis for the Jacobian&#39;s nullspace – the kernel or the "zero energ...

◆ system() [1/2]

template<typename T>
const sys_type& libMesh::NonlinearSolver< T >::system ( ) const
inlineinherited

◆ system() [2/2]

template<typename T>
sys_type& libMesh::NonlinearSolver< T >::system ( )
inlineinherited
Returns
A writable reference to the system we are solving.

Definition at line 278 of file nonlinear_solver.h.

278 { return _system; }
sys_type & _system
A reference to the system we are solving.

◆ use_default_monitor()

template<typename T>
void libMesh::PetscNonlinearSolver< T >::use_default_monitor ( bool  state)
inline

Set to true to use the libMesh's default monitor, set to false to use your own.

Definition at line 167 of file petsc_nonlinear_solver.h.

167 { _default_monitor = state; }
bool _default_monitor
true if we want the default monitor to be set, false for no monitor (i.e.

Friends And Related Function Documentation

◆ libmesh_petsc_snes_fd_residual

template<typename T>
PetscErrorCode libmesh_petsc_snes_fd_residual ( SNES  snes,
Vec  x,
Vec  r,
void *  ctx 
)
friend

Definition at line 256 of file petsc_nonlinear_solver.C.

257  {
258  ResidualContext rc = libmesh_petsc_snes_residual_helper(snes, x, ctx);
259 
260  libmesh_parallel_only(rc.sys.comm());
261 
262  libmesh_assert(r);
263  PetscVector<Number> R(r, rc.sys.comm());
264 
265  if (rc.solver->_zero_out_residual)
266  R.zero();
267 
268  if (rc.solver->fd_residual_object != nullptr)
269  rc.solver->fd_residual_object->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
270 
271  else if (rc.solver->residual_object != nullptr)
272  rc.solver->residual_object->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
273 
274  else
275  libmesh_error_msg("Error! Unable to compute residual for forming finite difference Jacobian!");
276 
277  // Synchronize PETSc x to local solution since the local solution may be changed due to the constraints
278  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(rc.sys.solution.get());
279  PetscVector<Number> X_global(x, rc.sys.comm());
280 
281  X_global.swap(X_sys);
282  rc.sys.update();
283  X_global.swap(X_sys);
284 
285  R.close();
286 
287  if (rc.solver->_exact_constraint_enforcement)
288  {
289  rc.sys.get_dof_map().enforce_constraints_on_residual(rc.sys, &R, rc.sys.current_local_solution.get());
290  R.close();
291  }
292 
293  return rc.ierr;
294  }
template class LIBMESH_EXPORT PetscVector< Number >
friend ResidualContext libmesh_petsc_snes_residual_helper(SNES snes, Vec x, void *ctx)
libmesh_assert(ctx)
void * ctx

◆ libmesh_petsc_snes_jacobian

template<typename T>
PetscErrorCode libmesh_petsc_snes_jacobian ( SNES  snes,
Vec  x,
Mat  jac,
Mat  pc,
void *  ctx 
)
friend

Definition at line 447 of file petsc_nonlinear_solver.C.

448  {
449  LOG_SCOPE("jacobian()", "PetscNonlinearSolver");
450 
451  PetscErrorCode ierr=0;
452 
454 
455  // No way to safety-check this cast, since we got a void *...
457  static_cast<PetscNonlinearSolver<Number> *> (ctx);
458 
459  libmesh_parallel_only(solver->comm());
460 
461  // Get the current iteration number from the snes object,
462  // store it in the PetscNonlinearSolver object for possible use
463  // by the user's Jacobian function.
464  {
465  PetscInt n_iterations = 0;
466  ierr = SNESGetIterationNumber(snes, &n_iterations);
467  LIBMESH_CHKERR2(solver->comm(),ierr);
468  solver->_current_nonlinear_iteration_number = cast_int<unsigned>(n_iterations);
469  }
470 
471  //-----------------------------------------------------------------------------
472  // if the user has provided both function pointers and objects only the pointer
473  // will be used, so catch that as an error
474  libmesh_error_msg_if(solver->jacobian && solver->jacobian_object,
475  "ERROR: cannot specify both a function and object to compute the Jacobian!");
476 
477  libmesh_error_msg_if(solver->matvec && solver->residual_and_jacobian_object,
478  "ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
479 
480  NonlinearImplicitSystem & sys = solver->system();
481 
482  PetscMatrix<Number> PC(pc, sys.comm());
483  PetscMatrix<Number> Jac(jac, sys.comm());
484  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
485  PetscVector<Number> X_global(x, sys.comm());
486 
487  // We already computed the Jacobian during the residual evaluation
488  if (solver->residual_and_jacobian_object)
489  {
490  auto & sys_mat = static_cast<PetscMatrix<Number> &>(sys.get_system_matrix());
491 
492  // We could be doing matrix-free
493  if (jac && jac != sys_mat.mat())
494  Jac.close();
495  if (pc && pc != sys_mat.mat())
496  PC.close();
497 
498  return ierr;
499  }
500 
501  // Set the dof maps
502  PC.attach_dof_map(sys.get_dof_map());
503  Jac.attach_dof_map(sys.get_dof_map());
504 
505  // Use the systems update() to get a good local version of the parallel solution
506  X_global.swap(X_sys);
507  sys.update();
508  X_global.swap(X_sys);
509 
510  // Enforce constraints (if any) exactly on the
511  // current_local_solution. This is the solution vector that is
512  // actually used in the computation of the residual below, and is
513  // not locked by debug-enabled PETSc the way that "x" is.
514  if (solver->_exact_constraint_enforcement)
515  sys.get_dof_map().enforce_constraints_exactly(sys, sys.current_local_solution.get());
516 
517  if (solver->_zero_out_jacobian)
518  PC.zero();
519 
520 
521  if (solver->jacobian != nullptr)
522  solver->jacobian(*sys.current_local_solution.get(), PC, sys);
523 
524  else if (solver->jacobian_object != nullptr)
525  solver->jacobian_object->jacobian(*sys.current_local_solution.get(), PC, sys);
526 
527  else if (solver->matvec != nullptr)
528  solver->matvec(*sys.current_local_solution.get(), nullptr, &PC, sys);
529 
530  else
531  libmesh_error_msg("Error! Unable to compute residual and/or Jacobian!");
532 
533  PC.close();
534  if (solver->_exact_constraint_enforcement)
535  {
536  sys.get_dof_map().enforce_constraints_on_jacobian(sys, &PC);
537  PC.close();
538  }
539 
540  Jac.close();
541 
542  return ierr;
543  }
template class LIBMESH_EXPORT PetscMatrix< Number >
template class LIBMESH_EXPORT PetscNonlinearSolver< Number >
template class LIBMESH_EXPORT PetscVector< Number >
libmesh_assert(ctx)
void * ctx

◆ libmesh_petsc_snes_mffd_residual

template<typename T>
PetscErrorCode libmesh_petsc_snes_mffd_residual ( SNES  snes,
Vec  x,
Vec  r,
void *  ctx 
)
friend

Definition at line 309 of file petsc_nonlinear_solver.C.

310  {
311  ResidualContext rc = libmesh_petsc_snes_residual_helper(snes, x, ctx);
312 
313  libmesh_parallel_only(rc.sys.comm());
314 
315  libmesh_assert(r);
316  PetscVector<Number> R(r, rc.sys.comm());
317 
318  if (rc.solver->_zero_out_residual)
319  R.zero();
320 
321  if (rc.solver->mffd_residual_object != nullptr)
322  rc.solver->mffd_residual_object->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
323 
324  else if (rc.solver->residual_object != nullptr)
325  rc.solver->residual_object->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
326 
327  else
328  libmesh_error_msg("Error! Unable to compute residual for forming finite differenced"
329  "Jacobian-vector products!");
330 
331  // Synchronize PETSc x to local solution since the local solution may be changed due to the constraints
332  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(rc.sys.solution.get());
333  PetscVector<Number> X_global(x, rc.sys.comm());
334 
335  X_global.swap(X_sys);
336  rc.sys.update();
337  X_global.swap(X_sys);
338 
339  R.close();
340 
341  if (rc.solver->_exact_constraint_enforcement)
342  {
343  rc.sys.get_dof_map().enforce_constraints_on_residual(rc.sys, &R, rc.sys.current_local_solution.get());
344  R.close();
345  }
346 
347  return rc.ierr;
348  }
template class LIBMESH_EXPORT PetscVector< Number >
friend ResidualContext libmesh_petsc_snes_residual_helper(SNES snes, Vec x, void *ctx)
libmesh_assert(ctx)
void * ctx

◆ libmesh_petsc_snes_precheck

template<typename T>
PetscErrorCode libmesh_petsc_snes_precheck ( SNESLineSearch  ,
Vec  X,
Vec  Y,
PetscBool *  changed,
void *  context 
)
friend

Definition at line 662 of file petsc_nonlinear_solver.C.

663  {
664  LOG_SCOPE("precheck()", "PetscNonlinearSolver");
665 
666  PetscErrorCode ierr = 0;
667 
668  // PETSc almost certainly initializes these to false already, but
669  // it doesn't hurt to be explicit.
670  *changed = PETSC_FALSE;
671 
672  libmesh_assert(context);
673 
674  // Cast the context to a NonlinearSolver object.
676  static_cast<PetscNonlinearSolver<Number> *> (context);
677 
678  libmesh_parallel_only(solver->comm());
679 
680  // It's possible that we don't need to do anything at all, in
681  // that case return early...
682  if (!solver->precheck_object)
683  return ierr;
684 
685  // The user sets these flags in his/her postcheck function to
686  // indicate whether they changed something.
687  bool
688  petsc_changed = false;
689 
690  auto & sys = solver->system();
691  auto & x_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
692  PetscVector<Number> petsc_x(X, sys.comm());
693  PetscVector<Number> petsc_y(Y, sys.comm());
694 
695  // Use the systems update() to get a good local version of the parallel solution
696  petsc_x.swap(x_sys);
697  sys.update();
698  petsc_x.swap(x_sys);
699 
700  // Enforce constraints (if any) exactly on the
701  // current_local_solution. This is the solution vector that is
702  // actually used in the computation of residuals and Jacobians, and is
703  // not locked by debug-enabled PETSc the way that "x" is.
704  libmesh_assert(sys.current_local_solution.get());
705  auto & local_soln = *sys.current_local_solution.get();
706  if (solver->_exact_constraint_enforcement)
707  sys.get_dof_map().enforce_constraints_exactly(sys, &local_soln);
708 
709  solver->precheck_object->precheck(local_soln,
710  petsc_y,
711  petsc_changed,
712  sys);
713 
714  // Record whether the user changed the solution or the search direction.
715  if (petsc_changed)
716  *changed = PETSC_TRUE;
717 
718  return ierr;
719  }
template class LIBMESH_EXPORT PetscNonlinearSolver< Number >
template class LIBMESH_EXPORT PetscVector< Number >
libmesh_assert(ctx)

◆ libmesh_petsc_snes_residual

template<typename T>
PetscErrorCode libmesh_petsc_snes_residual ( SNES  snes,
Vec  x,
Vec  r,
void *  ctx 
)
friend

Definition at line 173 of file petsc_nonlinear_solver.C.

174  {
175  ResidualContext rc = libmesh_petsc_snes_residual_helper(snes, x, ctx);
176 
177  libmesh_parallel_only(rc.sys.comm());
178 
179  libmesh_assert(r);
180  PetscVector<Number> R(r, rc.sys.comm());
181 
182  if (rc.solver->_zero_out_residual)
183  R.zero();
184 
185  //-----------------------------------------------------------------------------
186  // if the user has provided both function pointers and objects only the pointer
187  // will be used, so catch that as an error
188  libmesh_error_msg_if(rc.solver->residual && rc.solver->residual_object,
189  "ERROR: cannot specify both a function and object to compute the Residual!");
190 
191  libmesh_error_msg_if(rc.solver->matvec && rc.solver->residual_and_jacobian_object,
192  "ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
193 
194  if (rc.solver->residual != nullptr)
195  rc.solver->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
196 
197  else if (rc.solver->residual_object != nullptr)
198  rc.solver->residual_object->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
199 
200  else if (rc.solver->matvec != nullptr)
201  rc.solver->matvec (*rc.sys.current_local_solution.get(), &R, nullptr, rc.sys);
202 
203  else if (rc.solver->residual_and_jacobian_object != nullptr)
204  {
205  auto & jac = rc.sys.get_system_matrix();
206 
207  if (rc.solver->_zero_out_jacobian)
208  jac.zero();
209 
210  rc.solver->residual_and_jacobian_object->residual_and_jacobian(
211  *rc.sys.current_local_solution.get(), &R, &jac, rc.sys);
212 
213  jac.close();
214  if (rc.solver->_exact_constraint_enforcement)
215  {
216  rc.sys.get_dof_map().enforce_constraints_on_jacobian(rc.sys, &jac);
217  jac.close();
218  }
219  }
220 
221  else
222  libmesh_error_msg("Error! Unable to compute residual and/or Jacobian!");
223 
224 
225  // Synchronize PETSc x to local solution since the local solution may be changed due to the constraints
226  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(rc.sys.solution.get());
227  PetscVector<Number> X_global(x, rc.sys.comm());
228 
229  X_global.swap(X_sys);
230  rc.sys.update();
231  X_global.swap(X_sys);
232 
233  R.close();
234 
235  if (rc.solver->_exact_constraint_enforcement)
236  {
237  rc.sys.get_dof_map().enforce_constraints_on_residual(rc.sys, &R, rc.sys.current_local_solution.get());
238  R.close();
239  }
240 
241  return rc.ierr;
242  }
template class LIBMESH_EXPORT PetscVector< Number >
friend ResidualContext libmesh_petsc_snes_residual_helper(SNES snes, Vec x, void *ctx)
libmesh_assert(ctx)
void * ctx

◆ libmesh_petsc_snes_residual_helper

template<typename T>
ResidualContext libmesh_petsc_snes_residual_helper ( SNES  snes,
Vec  x,
void *  ctx 
)
friend

Definition at line 54 of file petsc_nonlinear_solver.C.

55 {
56  LOG_SCOPE("residual()", "PetscNonlinearSolver");
57 
58  PetscErrorCode ierr = 0;
59 
60  libmesh_assert(x);
62 
63  // No way to safety-check this cast, since we got a void *...
65  static_cast<PetscNonlinearSolver<Number> *> (ctx);
66 
67  libmesh_parallel_only(solver->comm());
68 
69  // Get the current iteration number from the snes object,
70  // store it in the PetscNonlinearSolver object for possible use
71  // by the user's residual function.
72  {
73  PetscInt n_iterations = 0;
74  ierr = SNESGetIterationNumber(snes, &n_iterations);
75  LIBMESH_CHKERR2(solver->comm(),ierr);
76  solver->_current_nonlinear_iteration_number = cast_int<unsigned>(n_iterations);
77  }
78 
79  NonlinearImplicitSystem & sys = solver->system();
80 
81  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
82 
83  PetscVector<Number> X_global(x, sys.comm());
84 
85  // Use the system's update() to get a good local version of the
86  // parallel solution. This operation does not modify the incoming
87  // "x" vector, it only localizes information from "x" into
88  // sys.current_local_solution.
89  X_global.swap(X_sys);
90  sys.update();
91  X_global.swap(X_sys);
92 
93  // Enforce constraints (if any) exactly on the
94  // current_local_solution. This is the solution vector that is
95  // actually used in the computation of the residual below, and is
96  // not locked by debug-enabled PETSc the way that "x" is.
97  if (solver->_exact_constraint_enforcement)
98  sys.get_dof_map().enforce_constraints_exactly(sys, sys.current_local_solution.get());
99 
100  return ResidualContext(solver, sys, ierr);
101 }
template class LIBMESH_EXPORT PetscNonlinearSolver< Number >
template class LIBMESH_EXPORT PetscVector< Number >
libmesh_assert(ctx)
void * ctx

Member Data Documentation

◆ _communicator

const Parallel::Communicator& libMesh::ParallelObject::_communicator
protectedinherited

◆ _computing_base_vector

template<typename T>
bool libMesh::PetscNonlinearSolver< T >::_computing_base_vector
protected

Whether we are computing the base vector for matrix-free finite differencing.

Definition at line 284 of file petsc_nonlinear_solver.h.

Referenced by libMesh::PetscNonlinearSolver< Number >::computing_base_vector(), and libMesh::PetscNonlinearSolver< Number >::set_computing_base_vector().

◆ _counts

ReferenceCounter::Counts libMesh::ReferenceCounter::_counts
staticprotectedinherited

Actually holds the data.

Definition at line 124 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::get_info().

◆ _current_nonlinear_iteration_number

template<typename T>
unsigned libMesh::PetscNonlinearSolver< T >::_current_nonlinear_iteration_number
protected

◆ _default_monitor

template<typename T>
bool libMesh::PetscNonlinearSolver< T >::_default_monitor
protected

true if we want the default monitor to be set, false for no monitor (i.e.

user code can use their own)

Definition at line 268 of file petsc_nonlinear_solver.h.

Referenced by libMesh::PetscNonlinearSolver< Number >::use_default_monitor().

◆ _enable_print_counter

bool libMesh::ReferenceCounter::_enable_print_counter = true
staticprotectedinherited

Flag to control whether reference count information is printed when print_info is called.

Definition at line 143 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::disable_print_counter_info(), libMesh::ReferenceCounter::enable_print_counter_info(), and libMesh::ReferenceCounter::print_info().

◆ _exact_constraint_enforcement

template<typename T>
bool libMesh::NonlinearSolver< T >::_exact_constraint_enforcement
protectedinherited

◆ _is_initialized

template<typename T>
bool libMesh::NonlinearSolver< T >::_is_initialized
protectedinherited

Flag indicating if the data structures have been initialized.

Definition at line 433 of file nonlinear_solver.h.

Referenced by libMesh::NonlinearSolver< Number >::initialized().

◆ _mutex

Threads::spin_mutex libMesh::ReferenceCounter::_mutex
staticprotectedinherited

Mutual exclusion object to enable thread-safe reference counting.

Definition at line 137 of file reference_counter.h.

◆ _n_linear_iterations

template<typename T>
PetscInt libMesh::PetscNonlinearSolver< T >::_n_linear_iterations
protected

Stores the total number of linear iterations from the last solve.

Definition at line 248 of file petsc_nonlinear_solver.h.

◆ _n_objects

Threads::atomic< unsigned int > libMesh::ReferenceCounter::_n_objects
staticprotectedinherited

The number of objects.

Print the reference count information when the number returns to 0.

Definition at line 132 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::n_objects(), libMesh::ReferenceCounter::ReferenceCounter(), and libMesh::ReferenceCounter::~ReferenceCounter().

◆ _preconditioner

template<typename T>
Preconditioner<T>* libMesh::NonlinearSolver< T >::_preconditioner
protectedinherited

Holds the Preconditioner object to be used for the linear solves.

Definition at line 438 of file nonlinear_solver.h.

◆ _reason

template<typename T>
SNESConvergedReason libMesh::PetscNonlinearSolver< T >::_reason
protected

Store the reason for SNES convergence/divergence for use even after the _snes has been cleared.

Note
print_converged_reason() will always try to get the current reason with SNESGetConvergedReason(), but if the SNES object has already been cleared, it will fall back on this stored value. This value is therefore necessarily not cleared by the clear() function.

Definition at line 243 of file petsc_nonlinear_solver.h.

◆ _reuse_preconditioner

template<typename T>
bool libMesh::NonlinearSolver< T >::_reuse_preconditioner
protectedinherited

Whether we should reuse the linear preconditioner.

Definition at line 412 of file nonlinear_solver.h.

◆ _reuse_preconditioner_max_linear_its

template<typename T>
unsigned int libMesh::NonlinearSolver< T >::_reuse_preconditioner_max_linear_its
protectedinherited

Number of linear iterations to retain the preconditioner.

Definition at line 423 of file nonlinear_solver.h.

◆ _setup_reuse

template<typename T>
bool libMesh::PetscNonlinearSolver< T >::_setup_reuse
protected

Whether we've triggered the preconditioner reuse.

Definition at line 289 of file petsc_nonlinear_solver.h.

◆ _snes

template<typename T>
WrappedPetsc<SNES> libMesh::PetscNonlinearSolver< T >::_snes
protected

Nonlinear solver context.

Definition at line 231 of file petsc_nonlinear_solver.h.

◆ _snesmf_reuse_base

template<typename T>
bool libMesh::PetscNonlinearSolver< T >::_snesmf_reuse_base
protected

True, If we want the base vector to be used for differencing even if the function provided to SNESSetFunction() is not the same as that provided to MatMFFDSetFunction().

https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/SNES/MatSNESMFSetReuseBase.html

Definition at line 275 of file petsc_nonlinear_solver.h.

Referenced by libMesh::PetscNonlinearSolver< Number >::set_snesmf_reuse_base(), and libMesh::PetscNonlinearSolver< Number >::snes_mf_reuse_base().

◆ _solver_configuration

template<typename T>
SolverConfiguration* libMesh::NonlinearSolver< T >::_solver_configuration
protectedinherited

Optionally store a SolverOptions object that can be used to set parameters like solver type, tolerances and iteration limits.

Definition at line 444 of file nonlinear_solver.h.

◆ _system

template<typename T>
sys_type& libMesh::NonlinearSolver< T >::_system
protectedinherited

A reference to the system we are solving.

Definition at line 428 of file nonlinear_solver.h.

Referenced by libMesh::NonlinearSolver< Number >::system().

◆ _zero_out_jacobian

template<typename T>
bool libMesh::PetscNonlinearSolver< T >::_zero_out_jacobian
protected

true to zero out jacobian before going into application level call-back, otherwise false

Definition at line 263 of file petsc_nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_jacobian(), libMesh::libmesh_petsc_snes_residual(), and libMesh::PetscNonlinearSolver< Number >::set_jacobian_zero_out().

◆ _zero_out_residual

template<typename T>
bool libMesh::PetscNonlinearSolver< T >::_zero_out_residual
protected

true to zero out residual before going into application level call-back, otherwise false

Definition at line 258 of file petsc_nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_fd_residual(), libMesh::libmesh_petsc_snes_mffd_residual(), libMesh::libmesh_petsc_snes_residual(), and libMesh::PetscNonlinearSolver< Number >::set_residual_zero_out().

◆ absolute_residual_tolerance

template<typename T>
double libMesh::NonlinearSolver< T >::absolute_residual_tolerance
inherited

The NonlinearSolver should exit after the residual is reduced to either less than absolute_residual_tolerance or less than relative_residual_tolerance times the initial residual.

Users should increase any of these tolerances that they want to use for a stopping condition.

Definition at line 305 of file nonlinear_solver.h.

◆ absolute_step_tolerance

template<typename T>
double libMesh::NonlinearSolver< T >::absolute_step_tolerance
inherited

The NonlinearSolver should exit after the full nonlinear step norm is reduced to either less than absolute_step_tolerance or less than relative_step_tolerance times the largest nonlinear solution which has been seen so far.

Users should increase any of these tolerances that they want to use for a stopping condition.

Note
Not all NonlinearSolvers support relative_step_tolerance!

Definition at line 328 of file nonlinear_solver.h.

◆ bounds

template<typename T>
void(* libMesh::NonlinearSolver< T >::bounds) (NumericVector< Number > &XL, NumericVector< Number > &XU, sys_type &S)
inherited

Function that computes the lower and upper bounds XL and XU on the solution of the nonlinear system.

Definition at line 190 of file nonlinear_solver.h.

◆ bounds_object

template<typename T>
NonlinearImplicitSystem::ComputeBounds* libMesh::NonlinearSolver< T >::bounds_object
inherited

Object that computes the bounds vectors \( XL \) and \( XU \).

Definition at line 196 of file nonlinear_solver.h.

◆ converged

template<typename T>
bool libMesh::NonlinearSolver< T >::converged
inherited

After a call to solve this will reflect whether or not the nonlinear solve was successful.

Definition at line 352 of file nonlinear_solver.h.

◆ divergence_tolerance

template<typename T>
double libMesh::NonlinearSolver< T >::divergence_tolerance
inherited

The NonlinearSolver should exit if the residual becomes greater than the initial residual times the divergence_tolerance.

Users should adjust this tolerances to prevent divergence of the NonlinearSolver.

Definition at line 315 of file nonlinear_solver.h.

◆ fd_residual_object

template<typename T>
NonlinearImplicitSystem::ComputeResidual* libMesh::NonlinearSolver< T >::fd_residual_object
inherited

Object that computes the residual R(X) of the nonlinear system at the input iterate X for the purpose of forming a finite-differenced Jacobian.

Definition at line 143 of file nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_fd_residual().

◆ initial_linear_tolerance

template<typename T>
double libMesh::NonlinearSolver< T >::initial_linear_tolerance
inherited

Any required linear solves will at first be done with this tolerance; the NonlinearSolver may tighten the tolerance for later solves.

Definition at line 341 of file nonlinear_solver.h.

◆ jacobian

template<typename T>
void(* libMesh::NonlinearSolver< T >::jacobian) (const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)
inherited

Function that computes the Jacobian J(X) of the nonlinear system at the input iterate X.

Definition at line 156 of file nonlinear_solver.h.

Referenced by libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), and libMesh::libmesh_petsc_snes_jacobian().

◆ jacobian_object

template<typename T>
NonlinearImplicitSystem::ComputeJacobian* libMesh::NonlinearSolver< T >::jacobian_object
inherited

Object that computes the Jacobian J(X) of the nonlinear system at the input iterate X.

Definition at line 164 of file nonlinear_solver.h.

Referenced by libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), and libMesh::libmesh_petsc_snes_jacobian().

◆ linesearch_object

template<typename T>
std::unique_ptr<ComputeLineSearchObject> libMesh::PetscNonlinearSolver< T >::linesearch_object

A callable object that can be used to specify a custom line-search.

Definition at line 204 of file petsc_nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_linesearch_shellfunc().

◆ matvec

template<typename T>
void(* libMesh::NonlinearSolver< T >::matvec) (const NumericVector< Number > &X, NumericVector< Number > *R, SparseMatrix< Number > *J, sys_type &S)
inherited

Function that computes either the residual \( R(X) \) or the Jacobian \( J(X) \) of the nonlinear system at the input iterate \( X \).

Note
Either R or J could be nullptr.

Definition at line 173 of file nonlinear_solver.h.

Referenced by libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::libmesh_petsc_snes_jacobian(), and libMesh::libmesh_petsc_snes_residual().

◆ max_function_evaluations

template<typename T>
unsigned int libMesh::NonlinearSolver< T >::max_function_evaluations
inherited

Maximum number of function evaluations.

Definition at line 293 of file nonlinear_solver.h.

◆ max_linear_iterations

template<typename T>
unsigned int libMesh::NonlinearSolver< T >::max_linear_iterations
inherited

Each linear solver step should exit after max_linear_iterations is exceeded.

Definition at line 335 of file nonlinear_solver.h.

◆ max_nonlinear_iterations

template<typename T>
unsigned int libMesh::NonlinearSolver< T >::max_nonlinear_iterations
inherited

Maximum number of non-linear iterations.

Definition at line 288 of file nonlinear_solver.h.

◆ mffd_residual_object

template<typename T>
NonlinearImplicitSystem::ComputeResidual* libMesh::NonlinearSolver< T >::mffd_residual_object
inherited

Object that computes the residual R(X) of the nonlinear system at the input iterate X for the purpose of forming Jacobian-vector products via finite differencing.

Definition at line 150 of file nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_mffd_residual().

◆ minimum_linear_tolerance

template<typename T>
double libMesh::NonlinearSolver< T >::minimum_linear_tolerance
inherited

The tolerance for linear solves is kept above this minimum.

Definition at line 346 of file nonlinear_solver.h.

◆ nearnullspace

template<typename T>
void(* libMesh::NonlinearSolver< T >::nearnullspace) (std::vector< NumericVector< Number > * > &sp, sys_type &S)
inherited

Function that computes a basis for the Jacobian's near nullspace – the set of "low energy modes" – that can be used for AMG coarsening, if the solver supports it (e.g., ML, PETSc's GAMG).

Definition at line 233 of file nonlinear_solver.h.

◆ nearnullspace_object

template<typename T>
NonlinearImplicitSystem::ComputeVectorSubspace* libMesh::NonlinearSolver< T >::nearnullspace_object
inherited

A callable object that computes a basis for the Jacobian's near nullspace – the set of "low energy modes" – that can be used for AMG coarsening, if the solver supports it (e.g., ML, PETSc's GAMG).

Definition at line 240 of file nonlinear_solver.h.

◆ nullspace

template<typename T>
void(* libMesh::NonlinearSolver< T >::nullspace) (std::vector< NumericVector< Number > * > &sp, sys_type &S)
inherited

Function that computes a basis for the Jacobian's nullspace – the kernel or the "zero energy modes" – that can be used in solving a degenerate problem iteratively, if the solver supports it (e.g., PETSc's KSP).

Definition at line 204 of file nonlinear_solver.h.

◆ nullspace_object

template<typename T>
NonlinearImplicitSystem::ComputeVectorSubspace* libMesh::NonlinearSolver< T >::nullspace_object
inherited

A callable object that computes a basis for the Jacobian's nullspace – the kernel or the "zero energy modes" – that can be used in solving a degenerate problem iteratively, if the solver supports it (e.g., PETSc's KSP).

Definition at line 212 of file nonlinear_solver.h.

◆ postcheck

template<typename T>
void(* libMesh::NonlinearSolver< T >::postcheck) (const NumericVector< Number > &old_soln, NumericVector< Number > &search_direction, NumericVector< Number > &new_soln, bool &changed_search_direction, bool &changed_new_soln, sys_type &S)
inherited

Function that performs a "check" on the Newton search direction and solution after each nonlinear step.

See documentation for the NonlinearImplicitSystem::ComputePostCheck object for more information about the calling sequence.

Definition at line 254 of file nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_postcheck().

◆ postcheck_object

template<typename T>
NonlinearImplicitSystem::ComputePostCheck* libMesh::NonlinearSolver< T >::postcheck_object
inherited

A callable object that is executed after each nonlinear iteration.

Allows the user to modify both the search direction and the solution vector in an application-specific way.

Definition at line 266 of file nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_postcheck().

◆ precheck_object

template<typename T>
NonlinearImplicitSystem::ComputePreCheck* libMesh::NonlinearSolver< T >::precheck_object
inherited

Definition at line 268 of file nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_precheck().

◆ relative_residual_tolerance

template<typename T>
double libMesh::NonlinearSolver< T >::relative_residual_tolerance
inherited

Definition at line 306 of file nonlinear_solver.h.

◆ relative_step_tolerance

template<typename T>
double libMesh::NonlinearSolver< T >::relative_step_tolerance
inherited

Definition at line 329 of file nonlinear_solver.h.

◆ residual

template<typename T>
void(* libMesh::NonlinearSolver< T >::residual) (const NumericVector< Number > &X, NumericVector< Number > &R, sys_type &S)
inherited

Function that computes the residual R(X) of the nonlinear system at the input iterate X.

Definition at line 129 of file nonlinear_solver.h.

Referenced by libMesh::Problem_Interface::computeF(), and libMesh::libmesh_petsc_snes_residual().

◆ residual_and_jacobian_object

template<typename T>
NonlinearImplicitSystem::ComputeResidualandJacobian* libMesh::NonlinearSolver< T >::residual_and_jacobian_object
inherited

Object that computes either the residual \( R(X) \) or the Jacobian \( J(X) \) of the nonlinear system at the input iterate \( X \).

Note
Either R or J could be nullptr.

Definition at line 185 of file nonlinear_solver.h.

Referenced by libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::libmesh_petsc_snes_jacobian(), and libMesh::libmesh_petsc_snes_residual().

◆ residual_object

template<typename T>
NonlinearImplicitSystem::ComputeResidual* libMesh::NonlinearSolver< T >::residual_object
inherited

Object that computes the residual R(X) of the nonlinear system at the input iterate X.

Definition at line 137 of file nonlinear_solver.h.

Referenced by libMesh::Problem_Interface::computeF(), libMesh::libmesh_petsc_snes_fd_residual(), libMesh::libmesh_petsc_snes_mffd_residual(), and libMesh::libmesh_petsc_snes_residual().

◆ transpose_nullspace

template<typename T>
void(* libMesh::NonlinearSolver< T >::transpose_nullspace) (std::vector< NumericVector< Number > * > &sp, sys_type &S)
inherited

Function that computes a basis for the transpose Jacobian's nullspace – when solving a degenerate problem iteratively, if the solver supports it (e.g., PETSc's KSP), it is used to remove contributions outside of R(jac)

Definition at line 219 of file nonlinear_solver.h.

◆ transpose_nullspace_object

template<typename T>
NonlinearImplicitSystem::ComputeVectorSubspace* libMesh::NonlinearSolver< T >::transpose_nullspace_object
inherited

A callable object that computes a basis for the transpose Jacobian's nullspace – when solving a degenerate problem iteratively, if the solver supports it (e.g., PETSc's KSP), it is used to remove contributions outside of R(jac)

Definition at line 226 of file nonlinear_solver.h.

◆ user_presolve

template<typename T>
void(* libMesh::NonlinearSolver< T >::user_presolve) (sys_type &S)
inherited

Customizable function pointer which users can attach to the solver.

Gets called prior to every call to solve().

Definition at line 246 of file nonlinear_solver.h.


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