libMesh::DiffSolver Class Referenceabstract

#include <diff_solver.h>

Inheritance diagram for libMesh::DiffSolver:

Public Types

enum  SolveResult {
  INVALID_SOLVE_RESULT = 0, CONVERGED_NO_REASON = 1, CONVERGED_ABSOLUTE_RESIDUAL = 2, CONVERGED_RELATIVE_RESIDUAL = 4,
  CONVERGED_ABSOLUTE_STEP = 8, CONVERGED_RELATIVE_STEP = 16, DIVERGED_NO_REASON = 32, DIVERGED_MAX_NONLINEAR_ITERATIONS = 64,
  DIVERGED_BACKTRACKING_FAILURE = 128, DIVERGED_LINEAR_SOLVER_FAILURE = 256
}
 
typedef ImplicitSystem sys_type
 

Public Member Functions

 DiffSolver (sys_type &s)
 
virtual ~DiffSolver ()
 
virtual void init ()
 
virtual void reinit ()
 
virtual unsigned int solve ()=0
 
unsigned int total_outer_iterations ()
 
unsigned int total_inner_iterations ()
 
unsigned int solve_result ()
 
const sys_typesystem () const
 
sys_typesystem ()
 
const Parallel::Communicatorcomm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Static Public Member Functions

static std::unique_ptr< DiffSolverbuild (sys_type &s)
 
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

unsigned int max_linear_iterations
 
unsigned int max_nonlinear_iterations
 
bool quiet
 
bool verbose
 
bool continue_after_max_iterations
 
bool continue_after_backtrack_failure
 
Real absolute_residual_tolerance
 
Real relative_residual_tolerance
 
Real absolute_step_tolerance
 
Real relative_step_tolerance
 
Real initial_linear_tolerance
 
Real minimum_linear_tolerance
 
std::unique_ptr< LinearSolutionMonitorlinear_solution_monitor
 

Protected Types

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

Protected Member Functions

void increment_constructor_count (const std::string &name)
 
void increment_destructor_count (const std::string &name)
 

Protected Attributes

Real max_solution_norm
 
Real max_residual_norm
 
unsigned int _outer_iterations
 
unsigned int _inner_iterations
 
sys_type_system
 
unsigned int _solve_result
 
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
 

Detailed Description

This is a generic class that defines a solver to handle ImplicitSystem classes, including NonlinearImplicitSystem and DifferentiableSystem A user can define a solver by deriving from this class and implementing certain functions.

This class is part of the new DifferentiableSystem framework, which is still experimental. Users of this framework should beware of bugs and future API changes.

Author
Roy H. Stogner
Date
2006-2010

Definition at line 71 of file diff_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 78 of file diff_solver.h.

Member Enumeration Documentation

◆ SolveResult

Enumeration return type for the solve() function. Multiple SolveResults may be combined (OR'd) in the single return. To test which ones are present, just AND the return value with any of the SolveResult flags defined below.

Enumerator
INVALID_SOLVE_RESULT 

A default or invalid solve result. This usually means no solve has occurred yet.

CONVERGED_NO_REASON 

The solver converged but no particular reason is specified.

CONVERGED_ABSOLUTE_RESIDUAL 

The DiffSolver achieved the desired absolute residual tolerance.

CONVERGED_RELATIVE_RESIDUAL 

The DiffSolver achieved the desired relative residual tolerance.

CONVERGED_ABSOLUTE_STEP 

The DiffSolver achieved the desired absolute step size tolerance.

CONVERGED_RELATIVE_STEP 

The DiffSolver achieved the desired relative step size tolerance.

DIVERGED_NO_REASON 

The DiffSolver diverged but no particular reason is specified.

DIVERGED_MAX_NONLINEAR_ITERATIONS 

The DiffSolver reached the maximum allowed number of nonlinear iterations before satisfying any convergence tests.

DIVERGED_BACKTRACKING_FAILURE 

The DiffSolver failed to find a descent direction by backtracking (See newton_solver.C)

DIVERGED_LINEAR_SOLVER_FAILURE 

The linear solver used by the DiffSolver failed to find a solution.

Definition at line 223 of file diff_solver.h.

Constructor & Destructor Documentation

◆ DiffSolver()

libMesh::DiffSolver::DiffSolver ( sys_type s)

Constructor. Requires a reference to the system to be solved.

Definition at line 28 of file diff_solver.C.

28  :
29  ParallelObject(s),
32  quiet(true),
33  verbose(false),
46  _system (s),
48 {
49 }
bool continue_after_max_iterations
Definition: diff_solver.h:175
ParallelObject(const Parallel::Communicator &comm_in)
unsigned int max_nonlinear_iterations
Definition: diff_solver.h:157
Real absolute_residual_tolerance
Definition: diff_solver.h:192
static const Real TOLERANCE
unsigned int _solve_result
Definition: diff_solver.h:327
unsigned int _inner_iterations
Definition: diff_solver.h:314
sys_type & _system
Definition: diff_solver.h:319
unsigned int max_linear_iterations
Definition: diff_solver.h:149
bool continue_after_backtrack_failure
Definition: diff_solver.h:181
unsigned int _outer_iterations
Definition: diff_solver.h:309
Real relative_residual_tolerance
Definition: diff_solver.h:193

◆ ~DiffSolver()

virtual libMesh::DiffSolver::~DiffSolver ( )
inlinevirtual

Destructor.

Definition at line 97 of file diff_solver.h.

97 {}

Member Function Documentation

◆ build()

std::unique_ptr< DiffSolver > libMesh::DiffSolver::build ( sys_type s)
static

Factory method. Requires a reference to the system to be solved.

Returns
A NewtonSolver by default.

Definition at line 53 of file diff_solver.C.

Referenced by libMesh::TimeSolver::init().

54 {
55  return libmesh_make_unique<NewtonSolver>(s);
56 }

◆ 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_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

◆ 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()

void libMesh::DiffSolver::init ( )
virtual

The initialization function. This method is used to initialize internal data structures before a simulation begins.

Reimplemented in libMesh::PetscDiffSolver, and libMesh::NewtonSolver.

Definition at line 69 of file diff_solver.C.

References max_residual_norm, and max_solution_norm.

Referenced by libMesh::NewtonSolver::init(), and libMesh::PetscDiffSolver::init().

70 {
71  // Reset the max_step_size and max_residual_norm for a new problem
72  max_solution_norm = 0.;
73  max_residual_norm = 0.;
74 }

◆ 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_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

◆ reinit()

void libMesh::DiffSolver::reinit ( )
virtual

The reinitialization function. This method is used after changes in the mesh.

Reimplemented in libMesh::NewtonSolver, and libMesh::PetscDiffSolver.

Definition at line 60 of file diff_solver.C.

References max_residual_norm, and max_solution_norm.

Referenced by libMesh::PetscDiffSolver::reinit(), and libMesh::NewtonSolver::reinit().

61 {
62  // Reset the max_step_size and max_residual_norm for a new mesh
63  max_solution_norm = 0.;
64  max_residual_norm = 0.;
65 }

◆ solve()

virtual unsigned int libMesh::DiffSolver::solve ( )
pure virtual

This method performs a solve. What occurs in this method will depend on the type of solver. See the subclasses for more details.

Implemented in libMesh::PetscDiffSolver, and libMesh::NewtonSolver.

◆ solve_result()

unsigned int libMesh::DiffSolver::solve_result ( )
inline
Returns
The value of the SolveResult from the last solve.

Definition at line 133 of file diff_solver.h.

References _solve_result.

133 { return _solve_result; }
unsigned int _solve_result
Definition: diff_solver.h:327

◆ system() [1/2]

const sys_type& libMesh::DiffSolver::system ( ) const
inline
Returns
A constant reference to the system we are solving.

Definition at line 138 of file diff_solver.h.

References _system.

Referenced by libMesh::__libmesh_petsc_diff_solver_jacobian(), and libMesh::__libmesh_petsc_diff_solver_residual().

138 { return _system; }
sys_type & _system
Definition: diff_solver.h:319

◆ system() [2/2]

sys_type& libMesh::DiffSolver::system ( )
inline
Returns
A writable reference to the system we are solving.

Definition at line 143 of file diff_solver.h.

References _system.

143 { return _system; }
sys_type & _system
Definition: diff_solver.h:319

◆ total_inner_iterations()

unsigned int libMesh::DiffSolver::total_inner_iterations ( )
inline
Returns
The number of "inner" (e.g. Krylov) iterations required by the last solve.

Definition at line 128 of file diff_solver.h.

References _inner_iterations.

128 { return _inner_iterations; }
unsigned int _inner_iterations
Definition: diff_solver.h:314

◆ total_outer_iterations()

unsigned int libMesh::DiffSolver::total_outer_iterations ( )
inline
Returns
The number of "outer" (e.g. quasi-Newton) iterations required by the last solve.

Definition at line 122 of file diff_solver.h.

References _outer_iterations.

122 { return _outer_iterations; }
unsigned int _outer_iterations
Definition: diff_solver.h:309

Member Data Documentation

◆ _communicator

◆ _counts

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

◆ _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().

◆ _inner_iterations

unsigned int libMesh::DiffSolver::_inner_iterations
protected

The number of inner iterations used by the last solve

Definition at line 314 of file diff_solver.h.

Referenced by libMesh::NewtonSolver::solve(), and total_inner_iterations().

◆ _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_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().

◆ _outer_iterations

unsigned int libMesh::DiffSolver::_outer_iterations
protected

The number of outer iterations used by the last solve

Definition at line 309 of file diff_solver.h.

Referenced by libMesh::NewtonSolver::line_search(), libMesh::NewtonSolver::solve(), and total_outer_iterations().

◆ _solve_result

unsigned int libMesh::DiffSolver::_solve_result
protected

Initialized to zero. solve_result is typically set internally in the solve() function before it returns. When non-zero, solve_result tells the result of the latest solve. See enum definition for description.

Definition at line 327 of file diff_solver.h.

Referenced by libMesh::NewtonSolver::line_search(), libMesh::NewtonSolver::solve(), solve_result(), and libMesh::NewtonSolver::test_convergence().

◆ _system

◆ absolute_residual_tolerance

Real libMesh::DiffSolver::absolute_residual_tolerance

The DiffSolver 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 192 of file diff_solver.h.

Referenced by libMesh::NewtonSolver::print_convergence(), libMesh::NewtonSolver::solve(), and libMesh::NewtonSolver::test_convergence().

◆ absolute_step_tolerance

Real libMesh::DiffSolver::absolute_step_tolerance

The DiffSolver 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.

Definition at line 204 of file diff_solver.h.

Referenced by libMesh::NewtonSolver::print_convergence(), libMesh::NewtonSolver::solve(), and libMesh::NewtonSolver::test_convergence().

◆ continue_after_backtrack_failure

bool libMesh::DiffSolver::continue_after_backtrack_failure

Defaults to false, telling the DiffSolver to throw an error when the backtracking scheme fails to find a descent direction.

Definition at line 181 of file diff_solver.h.

Referenced by libMesh::NewtonSolver::line_search().

◆ continue_after_max_iterations

bool libMesh::DiffSolver::continue_after_max_iterations

Defaults to true, telling the DiffSolver to continue rather than exit when a solve has reached its maximum number of nonlinear iterations.

Definition at line 175 of file diff_solver.h.

Referenced by libMesh::NewtonSolver::solve().

◆ initial_linear_tolerance

Real libMesh::DiffSolver::initial_linear_tolerance

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

Definition at line 211 of file diff_solver.h.

Referenced by libMesh::NewtonSolver::solve().

◆ linear_solution_monitor

std::unique_ptr<LinearSolutionMonitor> libMesh::DiffSolver::linear_solution_monitor

Pointer to functor which is called right after each linear solve

Definition at line 289 of file diff_solver.h.

Referenced by libMesh::__libmesh_petsc_diff_solver_monitor(), and libMesh::NewtonSolver::solve().

◆ max_linear_iterations

unsigned int libMesh::DiffSolver::max_linear_iterations

Each linear solver step should exit after max_linear_iterations is exceeded.

Definition at line 149 of file diff_solver.h.

Referenced by libMesh::ContinuationSystem::continuation_solve(), libMesh::NewtonSolver::solve(), and libMesh::ContinuationSystem::solve_tangent().

◆ max_nonlinear_iterations

unsigned int libMesh::DiffSolver::max_nonlinear_iterations

The DiffSolver should exit in failure if max_nonlinear_iterations is exceeded and continue_after_max_iterations is false, or should end the nonlinear solve if max_nonlinear_iterations is exceeded and continue_after_max_iterations is true.

Definition at line 157 of file diff_solver.h.

Referenced by libMesh::ContinuationSystem::continuation_solve(), libMesh::NewtonSolver::solve(), and libMesh::ContinuationSystem::update_solution().

◆ max_residual_norm

Real libMesh::DiffSolver::max_residual_norm
protected

The largest nonlinear residual which the DiffSolver has yet seen will be stored here, to be used for stopping criteria based on relative_residual_tolerance

Definition at line 304 of file diff_solver.h.

Referenced by init(), libMesh::NewtonSolver::print_convergence(), reinit(), libMesh::NewtonSolver::solve(), and libMesh::NewtonSolver::test_convergence().

◆ max_solution_norm

Real libMesh::DiffSolver::max_solution_norm
protected

The largest solution norm which the DiffSolver has yet seen will be stored here, to be used for stopping criteria based on relative_step_tolerance

Definition at line 297 of file diff_solver.h.

Referenced by init(), libMesh::NewtonSolver::print_convergence(), reinit(), libMesh::NewtonSolver::solve(), and libMesh::NewtonSolver::test_convergence().

◆ minimum_linear_tolerance

Real libMesh::DiffSolver::minimum_linear_tolerance

The tolerance for linear solves is kept above this minimum

Definition at line 216 of file diff_solver.h.

Referenced by libMesh::NewtonSolver::solve().

◆ quiet

bool libMesh::DiffSolver::quiet

The DiffSolver should not print anything to libMesh::out unless quiet is set to false; default is true.

Definition at line 163 of file diff_solver.h.

Referenced by libMesh::NewtonSolver::line_search(), and libMesh::NewtonSolver::solve().

◆ relative_residual_tolerance

Real libMesh::DiffSolver::relative_residual_tolerance

◆ relative_step_tolerance

Real libMesh::DiffSolver::relative_step_tolerance

◆ verbose


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