libMesh::PetscNonlinearSolver< T > Class Template Reference

#include <petsc_nonlinear_solver.h>

Inheritance diagram for libMesh::PetscNonlinearSolver< T >:

Classes

class  ComputeLineSearchObject
 

Public Types

typedef NonlinearImplicitSystem sys_type
 

Public Member Functions

 PetscNonlinearSolver (sys_type &system)
 
 ~PetscNonlinearSolver ()
 
virtual void clear () override
 
virtual void init (const char *name=nullptr) override
 
SNES snes ()
 
virtual std::pair< unsigned int, Realsolve (SparseMatrix< T > &, NumericVector< T > &, NumericVector< T > &, const double, const unsigned int) override
 
virtual void print_converged_reason () override
 
SNESConvergedReason get_converged_reason ()
 
virtual int get_total_linear_iterations () override
 
virtual unsigned get_current_nonlinear_iteration_number () const override
 
void set_residual_zero_out (bool state)
 
void set_jacobian_zero_out (bool state)
 
void use_default_monitor (bool state)
 
void set_snesmf_reuse_base (bool state)
 
bool initialized () const
 
const sys_typesystem () const
 
sys_typesystem ()
 
void attach_preconditioner (Preconditioner< T > *preconditioner)
 
void set_solver_configuration (SolverConfiguration &solver_configuration)
 
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())
 
static std::string get_info ()
 
static void print_info (std::ostream &out=libMesh::out)
 
static unsigned int n_objects ()
 
static void enable_print_counter_info ()
 
static void disable_print_counter_info ()
 

Public Attributes

std::unique_ptr< ComputeLineSearchObjectlinesearch_object
 
void(* residual )(const NumericVector< Number > &X, NumericVector< Number > &R, sys_type &S)
 
NonlinearImplicitSystem::ComputeResidualresidual_object
 
NonlinearImplicitSystem::ComputeResidualfd_residual_object
 
NonlinearImplicitSystem::ComputeResidualmffd_residual_object
 
void(* jacobian )(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)
 
NonlinearImplicitSystem::ComputeJacobianjacobian_object
 
void(* matvec )(const NumericVector< Number > &X, NumericVector< Number > *R, SparseMatrix< Number > *J, sys_type &S)
 
NonlinearImplicitSystem::ComputeResidualandJacobianresidual_and_jacobian_object
 
void(* bounds )(NumericVector< Number > &XL, NumericVector< Number > &XU, sys_type &S)
 
NonlinearImplicitSystem::ComputeBoundsbounds_object
 
void(* nullspace )(std::vector< NumericVector< Number > *> &sp, sys_type &S)
 
NonlinearImplicitSystem::ComputeVectorSubspacenullspace_object
 
void(* transpose_nullspace )(std::vector< NumericVector< Number > *> &sp, sys_type &S)
 
NonlinearImplicitSystem::ComputeVectorSubspacetranspose_nullspace_object
 
void(* nearnullspace )(std::vector< NumericVector< Number > *> &sp, sys_type &S)
 
NonlinearImplicitSystem::ComputeVectorSubspacenearnullspace_object
 
void(* user_presolve )(sys_type &S)
 
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)
 
NonlinearImplicitSystem::ComputePostCheckpostcheck_object
 
unsigned int max_nonlinear_iterations
 
unsigned int max_function_evaluations
 
Real absolute_residual_tolerance
 
Real relative_residual_tolerance
 
Real absolute_step_tolerance
 
Real relative_step_tolerance
 
unsigned int max_linear_iterations
 
Real initial_linear_tolerance
 
Real minimum_linear_tolerance
 
bool converged
 

Protected Types

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

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)
 
void increment_destructor_count (const std::string &name)
 

Protected Attributes

SNES _snes
 
SNESConvergedReason _reason
 
PetscInt _n_linear_iterations
 
unsigned _current_nonlinear_iteration_number
 
bool _zero_out_residual
 
bool _zero_out_jacobian
 
bool _default_monitor
 
bool _snesmf_reuse_base
 
sys_type_system
 
bool _is_initialized
 
Preconditioner< T > * _preconditioner
 
SolverConfiguration_solver_configuration
 
const Parallel::Communicator_communicator
 

Static Protected Attributes

static Counts _counts
 
static Threads::atomic< unsigned int > _n_objects
 
static Threads::spin_mutex _mutex
 
static bool _enable_print_counter = true
 

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, MatStructure *msflag, void *ctx)
 
PetscErrorCode libmesh_petsc_snes_jacobian (SNES snes, Vec x, Mat jac, Mat pc, void *ctx)
 

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 93 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 117 of file reference_counter.h.

◆ sys_type

The type of system

Definition at line 99 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 547 of file petsc_nonlinear_solver.C.

547  :
548  NonlinearSolver<T>(system_in),
549  linesearch_object(nullptr),
550  _reason(SNES_CONVERGED_ITERATING/*==0*/), // Arbitrary initial value...
553  _zero_out_residual(true),
554  _zero_out_jacobian(true),
555  _default_monitor(true),
556  _snesmf_reuse_base(true)
557 {
558 }
std::unique_ptr< ComputeLineSearchObject > linesearch_object

◆ ~PetscNonlinearSolver()

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

Destructor.

Definition at line 563 of file petsc_nonlinear_solver.C.

564 {
565  this->clear ();
566 }

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 71 of file nonlinear_solver.C.

72 {
73  if (this->_is_initialized)
74  libmesh_error_msg("Preconditioner must be attached before the solver is initialized!");
75 
76  _preconditioner = preconditioner;
77 }
Preconditioner< T > * _preconditioner

◆ 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 37 of file nonlinear_solver.C.

38 {
39  // Build the appropriate solver
40  switch (solver_package)
41  {
42 
43 #ifdef LIBMESH_HAVE_PETSC
44  case PETSC_SOLVERS:
45  return libmesh_make_unique<PetscNonlinearSolver<T>>(s);
46 #endif // LIBMESH_HAVE_PETSC
47 
48 #if defined(LIBMESH_TRILINOS_HAVE_NOX) && defined(LIBMESH_TRILINOS_HAVE_EPETRA)
49  case TRILINOS_SOLVERS:
50  return libmesh_make_unique<NoxNonlinearSolver<T>>(s);
51 #endif
52 
53  default:
54  libmesh_error_msg("ERROR: Unrecognized solver package: " << solver_package);
55  }
56 }
TRILINOS_SOLVERS
Definition: libmesh.C:244

◆ 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 694 of file petsc_nonlinear_solver.C.

697 {
698  PetscErrorCode ierr;
699  std::vector<NumericVector<Number> *> sp;
700  if (computeSubspaceObject)
701  (*computeSubspaceObject)(sp, this->system());
702  else
703  (*computeSubspace)(sp, this->system());
704 
705  *msp = PETSC_NULL;
706  if (sp.size())
707  {
708  Vec * modes;
709  PetscScalar * dots;
710  PetscInt nmodes = cast_int<PetscInt>(sp.size());
711 
712 #if PETSC_RELEASE_LESS_THAN(3,5,0)
713  ierr = PetscMalloc2(nmodes,Vec,&modes,nmodes,PetscScalar,&dots);
714 #else
715  ierr = PetscMalloc2(nmodes,&modes,nmodes,&dots);
716 #endif
717  LIBMESH_CHKERR(ierr);
718 
719  for (PetscInt i=0; i<nmodes; ++i)
720  {
721  PetscVector<T> * pv = cast_ptr<PetscVector<T> *>(sp[i]);
722  Vec v = pv->vec();
723 
724  ierr = VecDuplicate(v, modes+i);
725  LIBMESH_CHKERR(ierr);
726 
727  ierr = VecCopy(v,modes[i]);
728  LIBMESH_CHKERR(ierr);
729  }
730 
731  // Normalize.
732  ierr = VecNormalize(modes[0],PETSC_NULL);
733  LIBMESH_CHKERR(ierr);
734 
735  for (PetscInt i=1; i<nmodes; i++)
736  {
737  // Orthonormalize vec[i] against vec[0:i-1]
738  ierr = VecMDot(modes[i],i,modes,dots);
739  LIBMESH_CHKERR(ierr);
740 
741  for (PetscInt j=0; j<i; j++)
742  dots[j] *= -1.;
743 
744  ierr = VecMAXPY(modes[i],i,dots,modes);
745  LIBMESH_CHKERR(ierr);
746 
747  ierr = VecNormalize(modes[i],PETSC_NULL);
748  LIBMESH_CHKERR(ierr);
749  }
750 
751  ierr = MatNullSpaceCreate(this->comm().get(), PETSC_FALSE, nmodes, modes, msp);
752  LIBMESH_CHKERR(ierr);
753 
754  for (PetscInt i=0; i<nmodes; ++i)
755  {
756  ierr = VecDestroy(modes+i);
757  LIBMESH_CHKERR(ierr);
758  }
759 
760  ierr = PetscFree2(modes,dots);
761  LIBMESH_CHKERR(ierr);
762  }
763 }
const Parallel::Communicator & comm() const
const sys_type & system() const
PetscErrorCode ierr

◆ 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 571 of file petsc_nonlinear_solver.C.

572 {
573  if (this->initialized())
574  {
575  this->_is_initialized = false;
576 
577  PetscErrorCode ierr=0;
578 
579  ierr = LibMeshSNESDestroy(&_snes);
580  LIBMESH_CHKERR(ierr);
581 
582  // Reset the nonlinear iteration counter. This information is only relevant
583  // *during* the solve(). After the solve is completed it should return to
584  // the default value of 0.
586  }
587 }
PetscErrorCode ierr

◆ comm()

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

Definition at line 89 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::BoundaryInfo::_find_id_maps(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult_add(), libMesh::EquationSystems::_read_impl(), libMesh::MeshRefinement::_refine_elements(), libMesh::MeshRefinement::_smooth_flags(), libMesh::PetscDMWrapper::add_dofs_helper(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::ImplicitSystem::add_matrix(), libMesh::System::add_vector(), libMesh::UnstructuredMesh::all_second_order(), libMesh::MeshTools::Modification::all_tri(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assemble_qoi(), libMesh::MeshCommunication::assign_global_indices(), libMesh::DofMap::attach_matrix(), 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_dims(), libMesh::System::calculate_norm(), libMesh::DofMap::check_dirichlet_bcid_consistency(), libMesh::PetscDMWrapper::check_section_n_dofs(), 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::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::ExodusII_IO::copy_elemental_solution(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::MeshTools::create_bounding_box(), 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::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), DMlibMeshFunction(), DMlibMeshJacobian(), DMlibMeshSetSystem_libMesh(), DMVariableBounds_libMesh(), libMesh::MeshRefinement::eliminate_unrefined_patches(), libMesh::EpetraVector< T >::EpetraVector(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::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_nelem_target(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::DofMap::get_info(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::LocationMap< T >::init(), libMesh::TimeSolver::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), libMesh::EigenSystem::init_data(), libMesh::EigenSystem::init_matrices(), libMesh::OptimizationSystem::initialize_equality_constraints_storage(), libMesh::OptimizationSystem::initialize_inequality_constraints_storage(), 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_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_snes_fd_residual(), libMesh::libmesh_petsc_snes_jacobian(), libMesh::libmesh_petsc_snes_mffd_residual(), libMesh::libmesh_petsc_snes_postcheck(), 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::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_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::FEMSystem::mesh_position_set(), libMesh::DistributedMesh::n_active_elem(), libMesh::MeshTools::n_active_levels(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::CondensedEigenSystem::n_global_non_condensed_dofs(), libMesh::MeshTools::n_levels(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::MeshTools::n_p_levels(), libMesh::BoundaryInfo::n_shellface_conds(), 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::MeshTools::paranoid_n_levels(), libMesh::petsc_auto_fieldsplit(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshBase::prepare_for_use(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), 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::MeshBase::recalculate_n_partitions(), libMesh::MeshRefinement::refine_and_coarsen_elements(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::CheckpointIO::select_split_config(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::PetscDiffSolver::setup_petsc_data(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::split_mesh(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), libMesh::MeshRefinement::test_level_one(), libMesh::MeshRefinement::test_unflagged(), libMesh::MeshTools::total_weight(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), 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().

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

◆ disable_print_counter_info()

void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited

Definition at line 106 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

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

107 {
108  _enable_print_counter = false;
109  return;
110 }

◆ 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 100 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

101 {
102  _enable_print_counter = true;
103  return;
104 }

◆ 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 985 of file petsc_nonlinear_solver.C.

986 {
987  PetscErrorCode ierr=0;
988 
989  if (this->initialized())
990  {
991  ierr = SNESGetConvergedReason(_snes, &_reason);
992  LIBMESH_CHKERR(ierr);
993  }
994 
995  return _reason;
996 }
PetscErrorCode ierr

◆ 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 164 of file petsc_nonlinear_solver.h.

◆ 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 & pr : _counts)
59  {
60  const std::string name(pr.first);
61  const unsigned int creations = pr.second.first;
62  const unsigned int destructions = pr.second.second;
63 
64  oss << "| " << name << " reference count information:\n"
65  << "| Creations: " << creations << '\n'
66  << "| Destructions: " << destructions << '\n';
67  }
68 
69  oss << " ---------------------------------------------------------------------------- \n";
70 
71  return oss.str();
72 
73 #else
74 
75  return "";
76 
77 #endif
78 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:42

◆ 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 999 of file petsc_nonlinear_solver.C.

1000 {
1001  return _n_linear_iterations;
1002 }

◆ increment_constructor_count()

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

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

Definition at line 181 of file reference_counter.h.

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

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

182 {
183  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
184  std::pair<unsigned int, unsigned int> & p = _counts[name];
185 
186  p.first++;
187 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:42
spin_mutex spin_mtx
Definition: threads.C:29

◆ increment_destructor_count()

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

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

Definition at line 194 of file reference_counter.h.

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

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

195 {
196  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
197  std::pair<unsigned int, unsigned int> & p = _counts[name];
198 
199  p.second++;
200 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:42
spin_mutex spin_mtx
Definition: threads.C:29

◆ 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 592 of file petsc_nonlinear_solver.C.

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

593 {
594  // Initialize the data structures if not done so already.
595  if (!this->initialized())
596  {
597  this->_is_initialized = true;
598 
599  PetscErrorCode ierr=0;
600 
601  ierr = SNESCreate(this->comm().get(),&_snes);
602  LIBMESH_CHKERR(ierr);
603 
604  if (name)
605  {
606  ierr = SNESSetOptionsPrefix(_snes, name);
607  LIBMESH_CHKERR(ierr);
608  }
609 
610 #if !PETSC_RELEASE_LESS_THAN(3,3,0)
611  // Attaching a DM to SNES.
612  DM dm;
613  ierr = DMCreate(this->comm().get(), &dm);LIBMESH_CHKERR(ierr);
614  ierr = DMSetType(dm,DMLIBMESH);LIBMESH_CHKERR(ierr);
615  ierr = DMlibMeshSetSystem(dm,this->system());LIBMESH_CHKERR(ierr);
616  if (name)
617  {
618  ierr = DMSetOptionsPrefix(dm,name); LIBMESH_CHKERR(ierr);
619  }
620  ierr = DMSetFromOptions(dm); LIBMESH_CHKERR(ierr);
621  ierr = DMSetUp(dm); LIBMESH_CHKERR(ierr);
622  ierr = SNESSetDM(this->_snes, dm); LIBMESH_CHKERR(ierr);
623  // SNES now owns the reference to dm.
624  ierr = DMDestroy(&dm); LIBMESH_CHKERR(ierr);
625 
626 #endif
627 
628  if (_default_monitor)
629  {
630  ierr = SNESMonitorSet (_snes, libmesh_petsc_snes_monitor,
631  this, PETSC_NULL);
632  LIBMESH_CHKERR(ierr);
633  }
634 
635  // If the SolverConfiguration object is provided, use it to set
636  // options during solver initialization.
637  if (this->_solver_configuration)
638  {
640  }
641 
642 #if PETSC_VERSION_LESS_THAN(3,1,0)
643  // Cannot call SNESSetOptions before SNESSetFunction when using
644  // any matrix free options with PETSc 3.1.0+
645  ierr = SNESSetFromOptions(_snes);
646  LIBMESH_CHKERR(ierr);
647 #endif
648 
649  if (this->_preconditioner)
650  {
651  KSP ksp;
652  ierr = SNESGetKSP (_snes, &ksp);
653  LIBMESH_CHKERR(ierr);
654  PC pc;
655  ierr = KSPGetPC(ksp,&pc);
656  LIBMESH_CHKERR(ierr);
657 
658  this->_preconditioner->init();
659 
660  PCSetType(pc, PCSHELL);
661  PCShellSetContext(pc,(void *)this->_preconditioner);
662 
663  //Re-Use the shell functions from petsc_linear_solver
664  PCShellSetSetUp(pc,libmesh_petsc_preconditioner_setup);
665  PCShellSetApply(pc,libmesh_petsc_preconditioner_apply);
666  }
667  }
668 
669 
670  // Tell PETSc about our linesearch "post-check" function, but only
671  // if the user has provided one. There seem to be extra,
672  // unnecessary residual calculations if a postcheck function is
673  // attached for no reason.
674  if (this->postcheck || this->postcheck_object)
675  {
676 #if PETSC_VERSION_LESS_THAN(3,3,0)
677  PetscErrorCode ierr = SNESLineSearchSetPostCheck(_snes, libmesh_petsc_snes_postcheck, this);
678  LIBMESH_CHKERR(ierr);
679 
680 #else
681  SNESLineSearch linesearch;
682  PetscErrorCode ierr = SNESGETLINESEARCH(_snes, &linesearch);
683  LIBMESH_CHKERR(ierr);
684 
685  ierr = SNESLineSearchSetPostCheck(linesearch, libmesh_petsc_snes_postcheck, this);
686  LIBMESH_CHKERR(ierr);
687 #endif
688  }
689 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:42
PetscErrorCode libmesh_petsc_snes_postcheck(#if PETSC_VERSION_LESS_THAN(3, 3, 0) SNES, Vec x, Vec y, Vec w, void *context, PetscBool *changed_y, PetscBool *changed_w #else SNESLineSearch, Vec x, Vec y, Vec w, PetscBool *changed_y, PetscBool *changed_w, void *context #endif)
Preconditioner< T > * _preconditioner
const Parallel::Communicator & comm() const
SolverConfiguration * _solver_configuration
PetscErrorCode libmesh_petsc_snes_monitor(SNES, PetscInt its, PetscReal fnorm, void *)
const sys_type & system() const
PetscErrorCode libmesh_petsc_preconditioner_setup(void *ctx)
PetscErrorCode DMlibMeshSetSystem(DM dm, libMesh::NonlinearImplicitSystem &sys)
NonlinearImplicitSystem::ComputePostCheck * postcheck_object
PetscErrorCode ierr
PetscErrorCode libmesh_petsc_preconditioner_apply(void *ctx, Vec x, Vec y)
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)

◆ initialized()

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

Definition at line 92 of file nonlinear_solver.h.

◆ n_objects()

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

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

Definition at line 83 of file reference_counter.h.

References libMesh::ReferenceCounter::_n_objects.

84  { return _n_objects; }
static Threads::atomic< unsigned int > _n_objects

◆ n_processors()

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

Definition at line 95 of file parallel_object.h.

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

Referenced by libMesh::BoundaryInfo::_find_id_maps(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::DistributedMesh::add_elem(), libMesh::DistributedMesh::add_node(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assembly(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::DistributedMesh::clear(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::UnstructuredMesh::create_pid_mesh(), libMesh::MeshTools::create_processor_bounding_box(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshBase::get_info(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), libMesh::Nemesis_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::MeshBase::partition(), libMesh::PetscLinearSolver< T >::PetscLinearSolver(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), 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::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), 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().

96  { return cast_int<processor_id_type>(_communicator.size()); }
processor_id_type size() const
Definition: communicator.h:175
const Parallel::Communicator & _communicator

◆ 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 975 of file petsc_nonlinear_solver.C.

976 {
977 
978  libMesh::out << "Nonlinear solver convergence/divergence reason: "
979  << SNESConvergedReasons[this->get_converged_reason()] << std::endl;
980 }
SNESConvergedReason get_converged_reason()
OStreamProxy out(std::cout)

◆ print_info()

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

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

Definition at line 87 of file reference_counter.C.

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

88 {
90  out_stream << ReferenceCounter::get_info();
91 }
static std::string get_info()

◆ processor_id()

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

Definition at line 101 of file parallel_object.h.

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

Referenced by libMesh::BoundaryInfo::_find_id_maps(), libMesh::EquationSystems::_read_impl(), 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::UnstructuredMesh::all_second_order(), libMesh::MeshTools::Modification::all_tri(), libMesh::FEMSystem::assembly(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::DistributedMesh::clear(), 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::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::MeshTools::correct_node_proc_ids(), libMesh::ExodusII_IO_Helper::create(), libMesh::DistributedMesh::delete_elem(), libMesh::DistributedMesh::delete_node(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DistributedMesh::DistributedMesh(), libMesh::DofMap::end_dof(), libMesh::DofMap::end_old_dof(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshFunction::find_element(), libMesh::MeshFunction::find_elements(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::DofMap::first_dof(), libMesh::DofMap::first_old_dof(), 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::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::SystemSubsetBySubdomain::init(), 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::DofMap::last_dof(), 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(), libMesh::MeshRefinement::make_coarsening_compatible(), 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::SparsityPattern::Build::operator()(), libMesh::DistributedMesh::own_node(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), 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::ExodusII_IO_Helper::read_elem_num_map(), libMesh::ExodusII_IO_Helper::read_global_values(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), libMesh::ExodusII_IO_Helper::read_node_num_map(), libMesh::System::read_parallel_data(), 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::DistributedMesh::renumber_dof_objects(), libMesh::CheckpointIO::select_split_config(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::MeshTools::total_weight(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::Parallel::Packing< T >::unpack(), libMesh::DistributedMesh::update_parallel_id_counts(), 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_elements(), 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::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::ExodusII_IO_Helper::write_nodal_values(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::ExodusII_IO_Helper::write_nodesets(), libMesh::System::write_parallel_data(), 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::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().

102  { return cast_int<processor_id_type>(_communicator.rank()); }
const Parallel::Communicator & _communicator
processor_id_type rank() const
Definition: communicator.h:173

◆ 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 175 of file petsc_nonlinear_solver.h.

◆ 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 170 of file petsc_nonlinear_solver.h.

◆ 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 185 of file petsc_nonlinear_solver.h.

◆ 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 80 of file nonlinear_solver.C.

81 {
82  _solver_configuration = &solver_configuration;
83 }
SolverConfiguration * _solver_configuration

◆ snes()

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

Definition at line 126 of file petsc_nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_mffd_interface().

126 { this->init(); return _snes; }
virtual void init(const char *name=nullptr) override

◆ 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 768 of file petsc_nonlinear_solver.C.

773 {
774  LOG_SCOPE("solve()", "PetscNonlinearSolver");
775  this->init ();
776 
777  // Make sure the data passed in are really of Petsc types
778  PetscMatrix<T> * pre = cast_ptr<PetscMatrix<T> *>(&pre_in);
779  PetscVector<T> * x = cast_ptr<PetscVector<T> *>(&x_in);
780  PetscVector<T> * r = cast_ptr<PetscVector<T> *>(&r_in);
781 
782  PetscErrorCode ierr=0;
783  PetscInt n_iterations =0;
784  // Should actually be a PetscReal, but I don't know which version of PETSc first introduced PetscReal
785  Real final_residual_norm=0.;
786 
787  ierr = SNESSetFunction (_snes, r->vec(), libmesh_petsc_snes_residual, this);
788  LIBMESH_CHKERR(ierr);
789 
790  // Only set the jacobian function if we've been provided with something to call.
791  // This allows a user to set their own jacobian function if they want to
792  if (this->jacobian || this->jacobian_object || this->residual_and_jacobian_object)
793  {
794  ierr = SNESSetJacobian (_snes, pre->mat(), pre->mat(), libmesh_petsc_snes_jacobian, this);
795  LIBMESH_CHKERR(ierr);
796  }
797 
798 
799 #if !PETSC_VERSION_LESS_THAN(3,3,0)
800  // Only set the nullspace if we have a way of computing it and the result is non-empty.
801  if (this->nullspace || this->nullspace_object)
802  {
803  MatNullSpace msp;
804  this->build_mat_null_space(this->nullspace_object, this->nullspace, &msp);
805  if (msp)
806  {
807  ierr = MatSetNullSpace(pre->mat(), msp);
808  LIBMESH_CHKERR(ierr);
809  ierr = MatNullSpaceDestroy(&msp);
810  LIBMESH_CHKERR(ierr);
811  }
812  }
813 
814  // Only set the transpose nullspace if we have a way of computing it and the result is non-empty.
816  {
817 #if PETSC_VERSION_LESS_THAN(3,6,0)
818  libmesh_warning("MatSetTransposeNullSpace is only supported for PETSc >= 3.6, transpose nullspace will be ignored.");
819 #else
820  MatNullSpace msp = PETSC_NULL;
822  if (msp)
823  {
824  ierr = MatSetTransposeNullSpace(pre->mat(), msp);
825  LIBMESH_CHKERR(ierr);
826  ierr = MatNullSpaceDestroy(&msp);
827  LIBMESH_CHKERR(ierr);
828  }
829 #endif
830  }
831 
832  // Only set the nearnullspace if we have a way of computing it and the result is non-empty.
833  if (this->nearnullspace || this->nearnullspace_object)
834  {
835  MatNullSpace msp = PETSC_NULL;
836  this->build_mat_null_space(this->nearnullspace_object, this->nearnullspace, &msp);
837 
838  if (msp)
839  {
840  ierr = MatSetNearNullSpace(pre->mat(), msp);
841  LIBMESH_CHKERR(ierr);
842  ierr = MatNullSpaceDestroy(&msp);
843  LIBMESH_CHKERR(ierr);
844  }
845  }
846 #endif
847 
848  // Have the Krylov subspace method use our good initial guess rather than 0
849  KSP ksp;
850  ierr = SNESGetKSP (_snes, &ksp);
851  LIBMESH_CHKERR(ierr);
852 
853  // Set the tolerances for the iterative solver. Use the user-supplied
854  // tolerance for the relative residual & leave the others at default values
855  ierr = KSPSetTolerances (ksp, this->initial_linear_tolerance, PETSC_DEFAULT,
856  PETSC_DEFAULT, this->max_linear_iterations);
857  LIBMESH_CHKERR(ierr);
858 
859  // Set the tolerances for the non-linear solver.
860  ierr = SNESSetTolerances(_snes, this->absolute_residual_tolerance, this->relative_residual_tolerance,
862  LIBMESH_CHKERR(ierr);
863 
864  //Pull in command-line options
865 #if PETSC_VERSION_LESS_THAN(3,7,0)
866  KSPSetFromOptions(ksp);
867 #endif
868  SNESSetFromOptions(_snes);
869 
870  if (this->user_presolve)
871  this->user_presolve(this->system());
872 
873  //Set the preconditioning matrix
874  if (this->_preconditioner)
875  {
876  this->_preconditioner->set_matrix(pre_in);
877  this->_preconditioner->init();
878  }
879 
880  // If the SolverConfiguration object is provided, use it to override
881  // solver options.
882  if (this->_solver_configuration)
883  {
885  }
886 
887  // In PETSc versions before 3.5.0, it is not possible to call
888  // SNESSetUp() before the solution and rhs vectors are initialized, as
889  // this triggers the
890  //
891  // "Solution vector cannot be right hand side vector!"
892  //
893  // error message. It is also not possible to call SNESSetSolution()
894  // in those versions of PETSc to work around the problem, since that
895  // API was removed in 3.0.0 and only restored in 3.6.0. The
896  // overzealous check was moved out of SNESSetUp in PETSc 3.5.0
897  // (petsc/petsc@154060b), so this code block should be safe to use
898  // in 3.5.0 and later.
899 #if !PETSC_VERSION_LESS_THAN(3,5,0)
900 #if !PETSC_VERSION_LESS_THAN(3,6,0)
901  ierr = SNESSetSolution(_snes, x->vec());
902  LIBMESH_CHKERR(ierr);
903 #endif
904  ierr = SNESSetUp(_snes);
905  LIBMESH_CHKERR(ierr);
906 
907  Mat J;
908  ierr = SNESGetJacobian(_snes,&J,PETSC_NULL,PETSC_NULL,PETSC_NULL);
909  LIBMESH_CHKERR(ierr);
910  ierr = MatMFFDSetFunction(J, libmesh_petsc_snes_mffd_interface, this);
911  LIBMESH_CHKERR(ierr);
912 #if !PETSC_VERSION_LESS_THAN(3, 8, 4)
913  // Resue the residual vector from SNES
914  ierr = MatSNESMFSetReuseBase(J, static_cast<PetscBool>(_snesmf_reuse_base));
915  LIBMESH_CHKERR(ierr);
916 #endif
917 #endif
918 
919 #if PETSC_VERSION_LESS_THAN(3, 3, 0)
920  if (linesearch_object)
921  libmesh_error_msg("Line search setter interface introduced in petsc version 3.3!");
922 #else
923  SNESLineSearch linesearch;
924  if (linesearch_object)
925  {
926  ierr = SNESGETLINESEARCH(_snes, &linesearch);
927  LIBMESH_CHKERR(ierr);
928  ierr = SNESLineSearchSetType(linesearch, SNESLINESEARCHSHELL);
929  LIBMESH_CHKERR(ierr);
930  ierr = SNESLineSearchShellSetUserFunc(linesearch, libmesh_petsc_linesearch_shellfunc, this);
931  LIBMESH_CHKERR(ierr);
932  }
933 #endif
934 
935  ierr = SNESSolve (_snes, PETSC_NULL, x->vec());
936  LIBMESH_CHKERR(ierr);
937 
938  ierr = SNESGetIterationNumber(_snes,&n_iterations);
939  LIBMESH_CHKERR(ierr);
940 
941  ierr = SNESGetLinearSolveIterations(_snes, &_n_linear_iterations);
942  LIBMESH_CHKERR(ierr);
943 
944  // Enforce constraints exactly now that the solve is done. We have
945  // been enforcing them on the current_local_solution during the
946  // solve, but now need to be sure they are enforced on the parallel
947  // solution vector as well.
949 
950  // SNESGetFunction has been around forever and should work on all
951  // versions of PETSc. This is also now the recommended approach
952  // according to the documentation for the PETSc 3.5.1 release:
953  // http://www.mcs.anl.gov/petsc/documentation/changes/35.html
954  Vec f;
955  ierr = SNESGetFunction(_snes, &f, 0, 0);
956  LIBMESH_CHKERR(ierr);
957  ierr = VecNorm(f, NORM_2, &final_residual_norm);
958  LIBMESH_CHKERR(ierr);
959 
960  // Get and store the reason for convergence
961  SNESGetConvergedReason(_snes, &_reason);
962 
963  //Based on Petsc 2.3.3 documentation all diverged reasons are negative
964  this->converged = (_reason >= 0);
965 
966  this->clear();
967 
968  // return the # of its. and the final residual norm.
969  return std::make_pair(n_iterations, final_residual_norm);
970 }
friend PetscErrorCode libmesh_petsc_snes_residual(SNES snes, Vec x, Vec r, void *ctx)
NonlinearImplicitSystem::ComputeResidualandJacobian * residual_and_jacobian_object
friend PetscErrorCode libmesh_petsc_snes_jacobian(SNES snes, Vec x, Mat *jac, Mat *pc, MatStructure *msflag, void *ctx)
unsigned int max_function_evaluations
Preconditioner< T > * _preconditioner
void(* user_presolve)(sys_type &S)
virtual void configure_solver()=0
unsigned int max_linear_iterations
SolverConfiguration * _solver_configuration
void(* jacobian)(const NumericVector< Number > &X, SparseMatrix< Number > &J, sys_type &S)
const sys_type & system() const
NonlinearImplicitSystem::ComputeJacobian * jacobian_object
void(* transpose_nullspace)(std::vector< NumericVector< Number > *> &sp, sys_type &S)
void(* nullspace)(std::vector< NumericVector< Number > *> &sp, sys_type &S)
void(* nearnullspace)(std::vector< NumericVector< Number > *> &sp, sys_type &S)
PetscErrorCode libmesh_petsc_snes_mffd_interface(void *ctx, Vec x, Vec r)
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
PetscErrorCode ierr
void build_mat_null_space(NonlinearImplicitSystem::ComputeVectorSubspace *computeSubspaceObject, void(*)(std::vector< NumericVector< Number > *> &, sys_type &), MatNullSpace *)
virtual void init(const char *name=nullptr) override
NonlinearImplicitSystem::ComputeVectorSubspace * nearnullspace_object
const DofMap & get_dof_map() const
Definition: system.h:2049
unsigned int max_nonlinear_iterations
NonlinearImplicitSystem::ComputeVectorSubspace * transpose_nullspace_object
void enforce_constraints_exactly(const System &system, NumericVector< Number > *v=nullptr, bool homogeneous=false) const
NonlinearImplicitSystem::ComputeVectorSubspace * nullspace_object

◆ system() [1/2]

◆ 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 285 of file nonlinear_solver.h.

285 { return _system; }

◆ 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 180 of file petsc_nonlinear_solver.h.

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 201 of file petsc_nonlinear_solver.C.

202  {
203  ResidualContext rc = libmesh_petsc_snes_residual_helper(snes, x, ctx);
204 
205  libmesh_assert(r);
206  PetscVector<Number> R(r, rc.sys.comm());
207 
208  if (rc.solver->_zero_out_residual)
209  R.zero();
210 
211  if (rc.solver->fd_residual_object != nullptr)
212  rc.solver->fd_residual_object->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
213 
214  else if (rc.solver->residual_object != nullptr)
215  rc.solver->residual_object->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
216 
217  else
218  libmesh_error_msg("Error! Unable to compute residual for forming finite difference Jacobian!");
219 
220  R.close();
221 
222  return rc.ierr;
223  }
friend ResidualContext libmesh_petsc_snes_residual_helper(SNES snes, Vec x, void *ctx)

◆ libmesh_petsc_snes_jacobian [1/2]

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

◆ libmesh_petsc_snes_jacobian [2/2]

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

◆ 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 238 of file petsc_nonlinear_solver.C.

239  {
240  ResidualContext rc = libmesh_petsc_snes_residual_helper(snes, x, ctx);
241 
242  libmesh_assert(r);
243  PetscVector<Number> R(r, rc.sys.comm());
244 
245  if (rc.solver->_zero_out_residual)
246  R.zero();
247 
248  if (rc.solver->mffd_residual_object != nullptr)
249  rc.solver->mffd_residual_object->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
250 
251  else if (rc.solver->residual_object != nullptr)
252  rc.solver->residual_object->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
253 
254  else
255  libmesh_error_msg("Error! Unable to compute residual for forming finite differenced"
256  "Jacobian-vector products!");
257 
258  R.close();
259 
260  return rc.ierr;
261  }
friend ResidualContext libmesh_petsc_snes_residual_helper(SNES snes, Vec x, void *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 150 of file petsc_nonlinear_solver.C.

151  {
152  ResidualContext rc = libmesh_petsc_snes_residual_helper(snes, x, ctx);
153 
154  libmesh_assert(r);
155  PetscVector<Number> R(r, rc.sys.comm());
156 
157  if (rc.solver->_zero_out_residual)
158  R.zero();
159 
160  //-----------------------------------------------------------------------------
161  // if the user has provided both function pointers and objects only the pointer
162  // will be used, so catch that as an error
163  if (rc.solver->residual && rc.solver->residual_object)
164  libmesh_error_msg("ERROR: cannot specify both a function and object to compute the Residual!");
165 
166  if (rc.solver->matvec && rc.solver->residual_and_jacobian_object)
167  libmesh_error_msg("ERROR: cannot specify both a function and object to compute the combined Residual & Jacobian!");
168 
169  if (rc.solver->residual != nullptr)
170  rc.solver->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
171 
172  else if (rc.solver->residual_object != nullptr)
173  rc.solver->residual_object->residual(*rc.sys.current_local_solution.get(), R, rc.sys);
174 
175  else if (rc.solver->matvec != nullptr)
176  rc.solver->matvec (*rc.sys.current_local_solution.get(), &R, nullptr, rc.sys);
177 
178  else if (rc.solver->residual_and_jacobian_object != nullptr)
179  rc.solver->residual_and_jacobian_object->residual_and_jacobian (*rc.sys.current_local_solution.get(), &R, nullptr, rc.sys);
180 
181  else
182  libmesh_error_msg("Error! Unable to compute residual and/or Jacobian!");
183 
184  R.close();
185 
186  return rc.ierr;
187  }
friend ResidualContext libmesh_petsc_snes_residual_helper(SNES snes, Vec x, 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 64 of file petsc_nonlinear_solver.C.

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

Member Data Documentation

◆ _communicator

◆ _counts

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

◆ _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 245 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 141 of file reference_counter.h.

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

◆ _is_initialized

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

Flag indicating if the data structures have been initialized.

Definition at line 366 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 135 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 225 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 130 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 371 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 220 of file petsc_nonlinear_solver.h.

◆ _snes

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

Nonlinear solver context

Definition at line 208 of file petsc_nonlinear_solver.h.

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

◆ _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 252 of file petsc_nonlinear_solver.h.

Referenced by libMesh::PetscNonlinearSolver< Number >::set_snesmf_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 377 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 361 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 240 of file petsc_nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_jacobian(), 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 235 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>
Real 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 312 of file nonlinear_solver.h.

◆ absolute_step_tolerance

template<typename T>
Real 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 326 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 199 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 205 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 350 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 152 of file nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_fd_residual().

◆ initial_linear_tolerance

template<typename T>
Real 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 339 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 165 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 173 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 201 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 182 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 300 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 333 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 295 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 159 of file nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_mffd_residual().

◆ minimum_linear_tolerance

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

The tolerance for linear solves is kept above this minimum

Definition at line 344 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 242 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 249 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 213 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 221 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 263 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 275 of file nonlinear_solver.h.

Referenced by libMesh::libmesh_petsc_snes_postcheck().

◆ relative_residual_tolerance

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

Definition at line 313 of file nonlinear_solver.h.

◆ relative_step_tolerance

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

Definition at line 327 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 138 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 194 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 146 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 228 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 235 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 255 of file nonlinear_solver.h.


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