libMesh
Public Types | Public Member Functions | Static Public Member Functions | Public Attributes | Protected Types | Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
libMesh::UnsteadySolver Class Referenceabstract

This is a generic class that defines a solver to handle time integration of DifferentiableSystems. More...

#include <unsteady_solver.h>

Inheritance diagram for libMesh::UnsteadySolver:
[legend]

Public Types

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

Public Member Functions

 UnsteadySolver (sys_type &s)
 Constructor. More...
 
virtual ~UnsteadySolver ()
 Destructor. More...
 
virtual void init () override
 The initialization function. More...
 
virtual void init_adjoints () override
 Add adjoint vectors and old_adjoint_vectors as per the indices of QoISet. More...
 
virtual void init_data () override
 The data initialization function. More...
 
virtual void reinit () override
 The reinitialization function. More...
 
virtual void solve () override
 This method solves for the solution at the next timestep. More...
 
virtual void advance_timestep () override
 This method advances the solution to the next timestep, after a solve() has been performed. More...
 
void update ()
 
virtual std::pair< unsigned int, Realadjoint_solve (const QoISet &qoi_indices) override
 This method solves for the adjoint solution at the next adjoint timestep (or a steady state adjoint solve) More...
 
virtual void adjoint_advance_timestep () override
 This method advances the adjoint solution to the previous timestep, after an adjoint_solve() has been performed. More...
 
virtual void retrieve_timestep () override
 This method retrieves all the stored solutions at the current system.time. More...
 
virtual void integrate_qoi_timestep () override
 A method to integrate the system::QoI functionals. More...
 
virtual void integrate_adjoint_sensitivity (const QoISet &qois, const ParameterVector &parameter_vector, SensitivityData &sensitivities) override
 A method to integrate the adjoint sensitivity w.r.t a given parameter vector. More...
 
virtual void integrate_adjoint_refinement_error_estimate (AdjointRefinementEstimator &, ErrorVector &) override
 A method to compute the adjoint refinement error estimate at the current timestep. More...
 
virtual Real error_order () const =0
 This method should return the expected convergence order of the (non-local) error of the time discretization scheme - e.g. More...
 
virtual unsigned int time_order () const =0
 
Number old_nonlinear_solution (const dof_id_type global_dof_number) const
 
virtual Real du (const SystemNorm &norm) const override
 Computes the size of ||u^{n+1} - u^{n}|| in some norm. More...
 
virtual bool is_steady () const override
 This is not a steady-state solver. More...
 
void set_first_adjoint_step (bool first_adjoint_step_setting)
 A setter for the first_adjoint_step boolean. More...
 
void set_first_solve (bool first_solve_setting)
 
virtual bool element_residual (bool request_jacobian, DiffContext &)=0
 This method uses the DifferentiablePhysics element_time_derivative(), element_constraint(), and mass_residual() to build a full residual on an element. More...
 
virtual bool side_residual (bool request_jacobian, DiffContext &)=0
 This method uses the DifferentiablePhysics side_time_derivative(), side_constraint(), and side_mass_residual() to build a full residual on an element's side. More...
 
virtual bool nonlocal_residual (bool request_jacobian, DiffContext &)=0
 This method uses the DifferentiablePhysics nonlocal_time_derivative(), nonlocal_constraint(), and nonlocal_mass_residual() to build a full residual of non-local terms. More...
 
virtual void before_timestep ()
 This method is for subclasses or users to override to do arbitrary processing between timesteps. More...
 
const sys_typesystem () const
 
sys_typesystem ()
 
virtual std::unique_ptr< DiffSolver > & diff_solver ()
 An implicit linear or nonlinear solver to use at each timestep. More...
 
virtual std::unique_ptr< LinearSolver< Number > > & linear_solver ()
 An implicit linear solver to use for adjoint and sensitivity problems. More...
 
void set_solution_history (const SolutionHistory &_solution_history)
 A setter function users will employ if they need to do something other than save no solution history. More...
 
SolutionHistoryget_solution_history ()
 A getter function that returns a reference to the solution history object owned by TimeSolver. More...
 
bool is_adjoint () const
 Accessor for querying whether we need to do a primal or adjoint solve. More...
 
void set_is_adjoint (bool _is_adjoint_value)
 Accessor for setting whether we need to do a primal or adjoint solve. More...
 
virtual Real last_completed_timestep_size ()
 Returns system.deltat if fixed timestep solver is used, the complete timestep size (sum of all substeps) if the adaptive time solver is used. More...
 

Static Public Member Functions

static std::string get_info ()
 Gets a string containing the reference information. More...
 
static void print_info (std::ostream &out_stream=libMesh::out)
 Prints the reference information, by default to libMesh::out. More...
 
static unsigned int n_objects ()
 Prints the number of outstanding (created, but not yet destroyed) objects. More...
 
static void enable_print_counter_info ()
 Methods to enable/disable the reference counter output from print_info() More...
 
static void disable_print_counter_info ()
 

Public Attributes

std::shared_ptr< NumericVector< Number > > old_local_nonlinear_solution
 Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be shared between different derived class instances, as in e.g. More...
 
bool quiet
 Print extra debugging information if quiet == false. More...
 
unsigned int reduce_deltat_on_diffsolver_failure
 This value (which defaults to zero) is the number of times the TimeSolver is allowed to halve deltat and let the DiffSolver repeat the latest failed solve with a reduced timestep. More...
 

Protected Types

typedef bool(DifferentiablePhysics::* ResFuncType) (bool, DiffContext &)
 Definitions of argument types for use in refactoring subclasses. More...
 
typedef void(DiffContext::* ReinitFuncType) (Real)
 
typedef std::map< std::string, std::pair< unsigned int, unsigned int > > Counts
 Data structure to log the information. More...
 

Protected Member Functions

void increment_constructor_count (const std::string &name) noexcept
 Increments the construction counter. More...
 
void increment_destructor_count (const std::string &name) noexcept
 Increments the destruction counter. More...
 

Protected Attributes

bool first_solve
 A bool that will be true the first time solve() is called, and false thereafter. More...
 
bool first_adjoint_step
 A bool that will be true the first time adjoint_advance_timestep() is called, (when the primal solution is to be used to set adjoint boundary conditions) and false thereafter. More...
 
std::vector< std::unique_ptr< NumericVector< Number > > > old_adjoints
 A vector of pointers to vectors holding the adjoint solution at the last time step. More...
 
Real last_step_deltat
 We will need to move the system.time around to ensure that residuals are built with the right deltat and the right time. More...
 
Real next_step_deltat
 
std::unique_ptr< DiffSolver_diff_solver
 An implicit linear or nonlinear solver to use at each timestep. More...
 
std::unique_ptr< LinearSolver< Number > > _linear_solver
 An implicit linear solver to use for adjoint problems. More...
 
sys_type_system
 A reference to the system we are solving. More...
 
std::unique_ptr< SolutionHistorysolution_history
 A std::unique_ptr to a SolutionHistory object. More...
 
Real last_deltat
 The deltat for the last completed timestep before the current one. More...
 

Static Protected Attributes

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

Detailed Description

This is a generic class that defines a solver to handle time integration of DifferentiableSystems.

A user can define a solver for unsteady problems 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
2008

Definition at line 50 of file unsteady_solver.h.

Member Typedef Documentation

◆ Counts

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

Data structure to log the information.

The log is identified by the class name.

Definition at line 119 of file reference_counter.h.

◆ ReinitFuncType

typedef void(DiffContext::* libMesh::TimeSolver::ReinitFuncType) (Real)
protectedinherited

Definition at line 327 of file time_solver.h.

◆ ResFuncType

typedef bool(DifferentiablePhysics::* libMesh::TimeSolver::ResFuncType) (bool, DiffContext &)
protectedinherited

Definitions of argument types for use in refactoring subclasses.

Definition at line 325 of file time_solver.h.

◆ sys_type

The type of system.

Definition at line 69 of file time_solver.h.

Constructor & Destructor Documentation

◆ UnsteadySolver()

libMesh::UnsteadySolver::UnsteadySolver ( sys_type s)
explicit

Constructor.

Requires a reference to the system to be solved.

Definition at line 37 of file unsteady_solver.C.

References libMesh::make_range(), libMesh::System::n_qois(), and old_adjoints.

38  : TimeSolver(s),
40  first_solve (true),
41  first_adjoint_step (true)
42 {
43  old_adjoints.resize(s.n_qois());
44 
45  // Set the old adjoint pointers to nullptrs
46  // We will use this nullness to skip the initial time instant,
47  // when there is no older adjoint.
48  for(auto j : make_range(s.n_qois()))
49  {
50  old_adjoints[j] = nullptr;
51  }
52 }
std::vector< std::unique_ptr< NumericVector< Number > > > old_adjoints
A vector of pointers to vectors holding the adjoint solution at the last time step.
TimeSolver(sys_type &s)
Constructor.
Definition: time_solver.C:36
std::shared_ptr< NumericVector< Number > > old_local_nonlinear_solution
Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be...
static std::unique_ptr< NumericVector< Number > > build(const Parallel::Communicator &comm, const SolverPackage solver_package=libMesh::default_solver_package())
Builds a NumericVector on the processors in communicator comm using the linear solver package specifi...
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
bool first_adjoint_step
A bool that will be true the first time adjoint_advance_timestep() is called, (when the primal soluti...
bool first_solve
A bool that will be true the first time solve() is called, and false thereafter.

◆ ~UnsteadySolver()

libMesh::UnsteadySolver::~UnsteadySolver ( )
virtualdefault

Destructor.

Member Function Documentation

◆ adjoint_advance_timestep()

void libMesh::UnsteadySolver::adjoint_advance_timestep ( )
overridevirtual

This method advances the adjoint solution to the previous timestep, after an adjoint_solve() has been performed.

This will be done before every UnsteadySolver::adjoint_solve().

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::AdaptiveTimeSolver, and libMesh::NewmarkSolver.

Definition at line 239 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::DifferentiableSystem::deltat, libMesh::System::get_adjoint_solution(), libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::System::get_vector(), libMesh::NumericVector< T >::localize(), libMesh::make_range(), libMesh::System::n_qois(), old_local_nonlinear_solution, libMesh::TimeSolver::solution_history, and libMesh::System::time.

240 {
241  // Call the store function to store the adjoint we have computed (or
242  // for first_adjoint_step, the adjoint initial condition) in this
243  // time step for the time instance.
244  solution_history->store(true, _system.time);
245 
246  // Before moving to the next time instant, copy over the current adjoint solutions into _old_adjoint_solutions
247  for(auto i : make_range(_system.n_qois()))
248  {
249  std::string old_adjoint_solution_name = "_old_adjoint_solution";
250  old_adjoint_solution_name+= std::to_string(i);
251  NumericVector<Number> & old_adjoint_solution_i = _system.get_vector(old_adjoint_solution_name);
252  NumericVector<Number> & adjoint_solution_i = _system.get_adjoint_solution(i);
253  old_adjoint_solution_i = adjoint_solution_i;
254  }
255 
257 
258  // Retrieve the primal solution vectors at this new (or for
259  // first_adjoint_step, initial) time instance. These provide the
260  // data to solve the adjoint problem for the next time instance.
261  solution_history->retrieve(true, _system.time);
262 
263  // Dont forget to localize the old_nonlinear_solution !
264  _system.get_vector("_old_nonlinear_solution").localize
267 }
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1595
unsigned int n_qois() const
Number of currently active quantities of interest.
Definition: system.h:2516
std::unique_ptr< SolutionHistory > solution_history
A std::unique_ptr to a SolutionHistory object.
Definition: time_solver.h:319
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
std::shared_ptr< NumericVector< Number > > old_local_nonlinear_solution
Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be...
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
Definition: diff_system.h:278
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
NumericVector< Number > & get_adjoint_solution(unsigned int i=0)
Definition: system.C:1193
const DofMap & get_dof_map() const
Definition: system.h:2293
template class LIBMESH_EXPORT NumericVector< Number >
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:511
const NumericVector< Number > & get_vector(std::string_view vec_name) const
Definition: system.C:918
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.

◆ adjoint_solve()

std::pair< unsigned int, Real > libMesh::UnsteadySolver::adjoint_solve ( const QoISet qoi_indices)
overridevirtual

This method solves for the adjoint solution at the next adjoint timestep (or a steady state adjoint solve)

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::AdaptiveTimeSolver, and libMesh::TwostepTimeSolver.

Definition at line 228 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::DifferentiableSystem::deltat, and libMesh::TimeSolver::last_deltat.

229 {
230  std::pair<unsigned int, Real> adjoint_output = _system.ImplicitSystem::adjoint_solve(qoi_indices);
231 
232  // Record the deltat we used for this adjoint timestep. This was determined completely
233  // by SolutionHistory::retrieve methods. The adjoint_solve methods should never change deltat.
235 
236  return adjoint_output;
237 }
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
Definition: diff_system.h:278
Real last_deltat
The deltat for the last completed timestep before the current one.
Definition: time_solver.h:332

◆ advance_timestep()

void libMesh::UnsteadySolver::advance_timestep ( )
overridevirtual

This method advances the solution to the next timestep, after a solve() has been performed.

Often this will be done after every UnsteadySolver::solve(), but adaptive mesh refinement and/or adaptive time step selection may require some solve() steps to be repeated.

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::AdaptiveTimeSolver, and libMesh::NewmarkSolver.

Definition at line 192 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::DifferentiableSystem::deltat, first_solve, libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::System::get_vector(), libMesh::NumericVector< T >::localize(), old_local_nonlinear_solution, libMesh::System::solution, libMesh::TimeSolver::solution_history, and libMesh::System::time.

Referenced by libMesh::NewmarkSolver::advance_timestep(), and solve().

193 {
194  // The first access of advance_timestep happens via solve, not user code
195  // It is used here to store any initial conditions data
196  if (!first_solve)
197  {
198  // We call advance_timestep in user code after solve, so any solutions
199  // we will be storing will be for the next time instance
201  }
202  else
203  {
204  // We are here because of a call to advance_timestep that happens
205  // via solve, the very first solve. All we are doing here is storing
206  // the initial condition. The actual solution computed via this solve
207  // will be stored when we call advance_timestep in the user's timestep loop
208  first_solve = false;
209  }
210 
211  // If the user has attached a memory or file solution history object
212  // to the solver, this will store the current solution indexed with
213  // the current time
214  solution_history->store(false, _system.time);
215 
216  NumericVector<Number> & old_nonlinear_soln =
217  _system.get_vector("_old_nonlinear_solution");
218  NumericVector<Number> & nonlinear_solution =
219  *(_system.solution);
220 
221  old_nonlinear_soln = nonlinear_solution;
222 
223  old_nonlinear_soln.localize
226 }
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1595
std::unique_ptr< SolutionHistory > solution_history
A std::unique_ptr to a SolutionHistory object.
Definition: time_solver.h:319
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
std::shared_ptr< NumericVector< Number > > old_local_nonlinear_solution
Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be...
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
Definition: diff_system.h:278
std::unique_ptr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1573
const DofMap & get_dof_map() const
Definition: system.h:2293
template class LIBMESH_EXPORT NumericVector< Number >
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:511
const NumericVector< Number > & get_vector(std::string_view vec_name) const
Definition: system.C:918
bool first_solve
A bool that will be true the first time solve() is called, and false thereafter.

◆ before_timestep()

virtual void libMesh::TimeSolver::before_timestep ( )
inlinevirtualinherited

This method is for subclasses or users to override to do arbitrary processing between timesteps.

Definition at line 205 of file time_solver.h.

205 {}

◆ diff_solver()

virtual std::unique_ptr<DiffSolver>& libMesh::TimeSolver::diff_solver ( )
inlinevirtualinherited

An implicit linear or nonlinear solver to use at each timestep.

Reimplemented in libMesh::AdaptiveTimeSolver.

Definition at line 220 of file time_solver.h.

References libMesh::TimeSolver::_diff_solver.

Referenced by libMesh::TimeSolver::adjoint_solve(), adjust_linear_solvers(), libMesh::TimeSolver::init(), libMesh::TimeSolver::init_data(), libMesh::TimeSolver::reinit(), and libMesh::TimeSolver::solve().

220 { return _diff_solver; }
std::unique_ptr< DiffSolver > _diff_solver
An implicit linear or nonlinear solver to use at each timestep.
Definition: time_solver.h:302

◆ disable_print_counter_info()

void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited

Definition at line 100 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

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

◆ du()

Real libMesh::UnsteadySolver::du ( const SystemNorm norm) const
overridevirtual

Computes the size of ||u^{n+1} - u^{n}|| in some norm.

Note
While you can always call this function, its result may or may not be very meaningful. For example, if you call this function right after calling advance_timestep() then you'll get a result of zero since old_nonlinear_solution is set equal to nonlinear_solution in this function.

Implements libMesh::TimeSolver.

Definition at line 348 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::System::calculate_norm(), libMesh::System::get_vector(), libMesh::TensorTools::norm(), and libMesh::System::solution.

349 {
350 
351  std::unique_ptr<NumericVector<Number>> solution_copy =
352  _system.solution->clone();
353 
354  solution_copy->add(-1., _system.get_vector("_old_nonlinear_solution"));
355 
356  solution_copy->close();
357 
358  return _system.calculate_norm(*solution_copy, norm);
359 }
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
std::unique_ptr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1573
Real calculate_norm(const NumericVector< Number > &v, unsigned int var, FEMNormType norm_type, std::set< unsigned int > *skip_dimensions=nullptr) const
Definition: system.C:1672
auto norm(const T &a) -> decltype(std::abs(a))
Definition: tensor_tools.h:74
const NumericVector< Number > & get_vector(std::string_view vec_name) const
Definition: system.C:918

◆ element_residual()

virtual bool libMesh::TimeSolver::element_residual ( bool  request_jacobian,
DiffContext  
)
pure virtualinherited

This method uses the DifferentiablePhysics element_time_derivative(), element_constraint(), and mass_residual() to build a full residual on an element.

What combination

it uses will depend on the type of solver. See the subclasses for more details.

Implemented in libMesh::EigenTimeSolver, libMesh::NewmarkSolver, libMesh::AdaptiveTimeSolver, libMesh::Euler2Solver, libMesh::SteadySolver, and libMesh::EulerSolver.

Referenced by libMesh::FEMSystem::numerical_elem_jacobian().

◆ enable_print_counter_info()

void libMesh::ReferenceCounter::enable_print_counter_info ( )
staticinherited

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

Definition at line 94 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

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

◆ error_order()

virtual Real libMesh::UnsteadySolver::error_order ( ) const
pure virtual

This method should return the expected convergence order of the (non-local) error of the time discretization scheme - e.g.

2 for the O(deltat^2) Crank-Nicholson, or 1 for the O(deltat) Backward Euler.

Useful for adaptive timestepping schemes.

Implemented in libMesh::NewmarkSolver, libMesh::AdaptiveTimeSolver, libMesh::Euler2Solver, and libMesh::EulerSolver.

◆ get_info()

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

Gets a string containing the reference information.

Definition at line 47 of file reference_counter.C.

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

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

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

◆ get_solution_history()

SolutionHistory & libMesh::TimeSolver::get_solution_history ( )
inherited

A getter function that returns a reference to the solution history object owned by TimeSolver.

Definition at line 124 of file time_solver.C.

References libMesh::TimeSolver::solution_history.

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

125 {
126  return *solution_history;
127 }
std::unique_ptr< SolutionHistory > solution_history
A std::unique_ptr to a SolutionHistory object.
Definition: time_solver.h:319

◆ increment_constructor_count()

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

Increments the construction counter.

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

Definition at line 183 of file reference_counter.h.

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

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

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

◆ increment_destructor_count()

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

Increments the destruction counter.

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

Definition at line 207 of file reference_counter.h.

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

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

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

◆ init()

void libMesh::UnsteadySolver::init ( )
overridevirtual

The initialization function.

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

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::AdaptiveTimeSolver, and libMesh::SecondOrderUnsteadySolver.

Definition at line 60 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::System::add_vector(), and libMesh::TimeSolver::init().

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

61 {
63 
64  _system.add_vector("_old_nonlinear_solution");
65 }
NumericVector< Number > & add_vector(std::string_view vec_name, const bool projections=true, const ParallelType type=PARALLEL)
Adds the additional vector vec_name to this system.
Definition: system.C:751
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
virtual void init()
The initialization function.
Definition: time_solver.C:72

◆ init_adjoints()

void libMesh::UnsteadySolver::init_adjoints ( )
overridevirtual

Add adjoint vectors and old_adjoint_vectors as per the indices of QoISet.

Reimplemented from libMesh::TimeSolver.

Definition at line 68 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::System::add_vector(), libMesh::GHOSTED, libMesh::TimeSolver::init_adjoints(), libMesh::make_range(), and libMesh::System::n_qois().

69 {
71 
72  // Add old adjoint solutions
73  // To keep the number of vectors consistent between the primal and adjoint
74  // time loops, we will also add the adjoint rhs vector during initialization
75  for(auto i : make_range(_system.n_qois()))
76  {
77  std::string old_adjoint_solution_name = "_old_adjoint_solution";
78  old_adjoint_solution_name+= std::to_string(i);
79  _system.add_vector(old_adjoint_solution_name, false, GHOSTED);
80 
81  std::string adjoint_rhs_name = "adjoint_rhs";
82  adjoint_rhs_name+= std::to_string(i);
83  _system.add_vector(adjoint_rhs_name, false, GHOSTED);
84  }
85 
86 }
unsigned int n_qois() const
Number of currently active quantities of interest.
Definition: system.h:2516
NumericVector< Number > & add_vector(std::string_view vec_name, const bool projections=true, const ParallelType type=PARALLEL)
Adds the additional vector vec_name to this system.
Definition: system.C:751
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
virtual void init_adjoints()
Initialize any adjoint related data structures, based on the number of qois.
Definition: time_solver.C:83

◆ init_data()

void libMesh::UnsteadySolver::init_data ( )
overridevirtual

The data initialization function.

This method is used to initialize internal data structures after the underlying System has been initialized

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::SecondOrderUnsteadySolver.

Definition at line 89 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::GHOSTED, libMesh::TimeSolver::init_data(), libMesh::System::n_dofs(), libMesh::System::n_local_dofs(), old_local_nonlinear_solution, and libMesh::SERIAL.

Referenced by libMesh::SecondOrderUnsteadySolver::init_data().

90 {
92 
93 #ifdef LIBMESH_ENABLE_GHOSTED
96  GHOSTED);
97 #else
99 #endif
100 }
virtual void init_data()
The data initialization function.
Definition: time_solver.C:97
dof_id_type n_local_dofs() const
Definition: system.C:150
dof_id_type n_dofs() const
Definition: system.C:113
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
std::shared_ptr< NumericVector< Number > > old_local_nonlinear_solution
Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be...
const DofMap & get_dof_map() const
Definition: system.h:2293
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:511

◆ integrate_adjoint_refinement_error_estimate()

void libMesh::UnsteadySolver::integrate_adjoint_refinement_error_estimate ( AdjointRefinementEstimator ,
ErrorVector  
)
overridevirtual

A method to compute the adjoint refinement error estimate at the current timestep.

int_{tstep_start}^{tstep_end} R(u^h,z) dt The user provides an initialized ARefEE object. Fills in an ErrorVector that contains the weighted sum of errors from all the QoIs and can be used to guide AMR. CURRENTLY ONLY SUPPORTED for Backward Euler.

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::EulerSolver, libMesh::FirstOrderUnsteadySolver, libMesh::TwostepTimeSolver, libMesh::AdaptiveTimeSolver, and libMesh::Euler2Solver.

Definition at line 331 of file unsteady_solver.C.

332 {
333  libmesh_not_implemented();
334 }

◆ integrate_adjoint_sensitivity()

void libMesh::UnsteadySolver::integrate_adjoint_sensitivity ( const QoISet qois,
const ParameterVector parameter_vector,
SensitivityData sensitivities 
)
overridevirtual

A method to integrate the adjoint sensitivity w.r.t a given parameter vector.

int_{tstep_start}^{tstep_end} dQ/dp dt = int_{tstep_start}^{tstep_end} ( / p) - ( R (u,z) / p ) dt The trapezoidal rule is used to numerically integrate the timestep.

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::AdaptiveTimeSolver, and libMesh::TwostepTimeSolver.

Definition at line 280 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::ImplicitSystem::adjoint_qoi_parameter_sensitivity(), libMesh::DifferentiableSystem::deltat, libMesh::make_range(), libMesh::Real, libMesh::System::remove_vector(), retrieve_timestep(), libMesh::ParameterVector::size(), libMesh::QoISet::size(), and libMesh::System::time.

281 {
282  // CURRENTLY using the trapezoidal rule to integrate each timestep
283  // (f(t_j) + f(t_j+1))/2 (t_j+1 - t_j)
284  // Fix me: This function needs to be moved to the EulerSolver classes like the
285  // other integrate_timestep functions, and use an integration rule consistent with
286  // the theta method used for the time integration.
287 
288  // Get t_j
289  Real time_left = _system.time;
290 
291  // Left side sensitivities to hold f(t_j)
292  SensitivityData sensitivities_left(qois, _system, parameter_vector);
293 
294  // Get f(t_j)
295  _system.adjoint_qoi_parameter_sensitivity(qois, parameter_vector, sensitivities_left);
296 
297  // Advance to t_j+1
299 
300  // Get t_j+1
301  Real time_right = _system.time;
302 
303  // Right side sensitivities f(t_j+1)
304  SensitivityData sensitivities_right(qois, _system, parameter_vector);
305 
306  // Remove the sensitivity rhs vector from system since we did not write it to file and it cannot be retrieved
307  _system.remove_vector("sensitivity_rhs0");
308 
309  // Retrieve the primal and adjoint solutions at the current timestep
311 
312  // Get f(t_j+1)
313  _system.adjoint_qoi_parameter_sensitivity(qois, parameter_vector, sensitivities_right);
314 
315  // Remove the sensitivity rhs vector from system since we did not write it to file and it cannot be retrieved
316  _system.remove_vector("sensitivity_rhs0");
317 
318  // Get the contributions for each sensitivity from this timestep
319  const auto pv_size = parameter_vector.size();
320  for (auto i : make_range(qois.size(_system)))
321  for (auto j : make_range(pv_size))
322  sensitivities[i][j] = ( (sensitivities_left[i][j] + sensitivities_right[i][j])/2. )*(time_right - time_left);
323 }
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1595
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
void remove_vector(std::string_view vec_name)
Removes the additional vector vec_name from this system.
Definition: system.C:846
virtual void retrieve_timestep() override
This method retrieves all the stored solutions at the current system.time.
virtual void adjoint_qoi_parameter_sensitivity(const QoISet &qoi_indices, const ParameterVector &parameters, SensitivityData &sensitivities) override
Solves for the derivative of each of the system&#39;s quantities of interest q in qoi[qoi_indices] with r...
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
Definition: diff_system.h:278
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134

◆ integrate_qoi_timestep()

void libMesh::UnsteadySolver::integrate_qoi_timestep ( )
overridevirtual

A method to integrate the system::QoI functionals.

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::EulerSolver, libMesh::FirstOrderUnsteadySolver, libMesh::AdaptiveTimeSolver, libMesh::TwostepTimeSolver, and libMesh::Euler2Solver.

Definition at line 325 of file unsteady_solver.C.

326 {
327  libmesh_not_implemented();
328 }

◆ is_adjoint()

bool libMesh::TimeSolver::is_adjoint ( ) const
inlineinherited

Accessor for querying whether we need to do a primal or adjoint solve.

Definition at line 277 of file time_solver.h.

References libMesh::TimeSolver::_is_adjoint.

Referenced by libMesh::FEMSystem::build_context().

278  { return _is_adjoint; }
bool _is_adjoint
This boolean tells the TimeSolver whether we are solving a primal or adjoint problem.
Definition: time_solver.h:340

◆ is_steady()

virtual bool libMesh::UnsteadySolver::is_steady ( ) const
inlineoverridevirtual

This is not a steady-state solver.

Implements libMesh::TimeSolver.

Definition at line 194 of file unsteady_solver.h.

194 { return false; }

◆ last_completed_timestep_size()

Real libMesh::TimeSolver::last_completed_timestep_size ( )
virtualinherited

Returns system.deltat if fixed timestep solver is used, the complete timestep size (sum of all substeps) if the adaptive time solver is used.

Returns the change in system.time, deltat, for the last timestep which was successfully completed. This only returns the outermost step size in the case of nested time solvers. If no time step has yet been successfully completed, then returns system.deltat.

Reimplemented in libMesh::AdaptiveTimeSolver.

Definition at line 160 of file time_solver.C.

References libMesh::TimeSolver::last_deltat.

161 {
162  return last_deltat;
163 }
Real last_deltat
The deltat for the last completed timestep before the current one.
Definition: time_solver.h:332

◆ linear_solver()

virtual std::unique_ptr<LinearSolver<Number> >& libMesh::TimeSolver::linear_solver ( )
inlinevirtualinherited

An implicit linear solver to use for adjoint and sensitivity problems.

Reimplemented in libMesh::AdaptiveTimeSolver.

Definition at line 225 of file time_solver.h.

References libMesh::TimeSolver::_linear_solver.

Referenced by libMesh::TimeSolver::init(), libMesh::TimeSolver::init_data(), and libMesh::TimeSolver::reinit().

225 { return _linear_solver; }
std::unique_ptr< LinearSolver< Number > > _linear_solver
An implicit linear solver to use for adjoint problems.
Definition: time_solver.h:307

◆ n_objects()

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

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

Definition at line 85 of file reference_counter.h.

References libMesh::ReferenceCounter::_n_objects.

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

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

◆ nonlocal_residual()

virtual bool libMesh::TimeSolver::nonlocal_residual ( bool  request_jacobian,
DiffContext  
)
pure virtualinherited

This method uses the DifferentiablePhysics nonlocal_time_derivative(), nonlocal_constraint(), and nonlocal_mass_residual() to build a full residual of non-local terms.

What combination it uses will depend on the type of solver. See the subclasses for more details.

Implemented in libMesh::NewmarkSolver, libMesh::EigenTimeSolver, libMesh::Euler2Solver, libMesh::AdaptiveTimeSolver, libMesh::SteadySolver, and libMesh::EulerSolver.

Referenced by libMesh::FEMSystem::numerical_nonlocal_jacobian().

◆ old_nonlinear_solution()

Number libMesh::UnsteadySolver::old_nonlinear_solution ( const dof_id_type  global_dof_number) const
Returns
The old nonlinear solution for the specified global DOF.

Definition at line 337 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::System::get_dof_map(), libMesh::DofMap::n_dofs(), and old_local_nonlinear_solution.

Referenced by libMesh::EulerSolver::_general_residual(), libMesh::Euler2Solver::_general_residual(), libMesh::NewmarkSolver::_general_residual(), and libMesh::FEMPhysics::eulerian_residual().

339 {
340  libmesh_assert_less (global_dof_number, _system.get_dof_map().n_dofs());
341  libmesh_assert_less (global_dof_number, old_local_nonlinear_solution->size());
342 
343  return (*old_local_nonlinear_solution)(global_dof_number);
344 }
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
dof_id_type n_dofs() const
Definition: dof_map.h:659
std::shared_ptr< NumericVector< Number > > old_local_nonlinear_solution
Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be...
const DofMap & get_dof_map() const
Definition: system.h:2293

◆ print_info()

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

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

Definition at line 81 of file reference_counter.C.

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

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

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

◆ reinit()

void libMesh::UnsteadySolver::reinit ( )
overridevirtual

The reinitialization function.

This method is used to resize internal data vectors after a mesh change.

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::SecondOrderUnsteadySolver, and libMesh::AdaptiveTimeSolver.

Definition at line 104 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::System::get_vector(), libMesh::GHOSTED, libMesh::NumericVector< T >::localize(), libMesh::System::n_dofs(), libMesh::System::n_local_dofs(), old_local_nonlinear_solution, libMesh::TimeSolver::reinit(), and libMesh::SERIAL.

Referenced by libMesh::SecondOrderUnsteadySolver::reinit().

105 {
107 
108 #ifdef LIBMESH_ENABLE_GHOSTED
110  _system.get_dof_map().get_send_list(), false,
111  GHOSTED);
112 #else
114 #endif
115 
116  // localize the old solution
117  NumericVector<Number> & old_nonlinear_soln =
118  _system.get_vector("_old_nonlinear_solution");
119 
120  old_nonlinear_soln.localize
123 }
virtual void reinit()
The reinitialization function.
Definition: time_solver.C:54
dof_id_type n_local_dofs() const
Definition: system.C:150
dof_id_type n_dofs() const
Definition: system.C:113
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
std::shared_ptr< NumericVector< Number > > old_local_nonlinear_solution
Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be...
const DofMap & get_dof_map() const
Definition: system.h:2293
template class LIBMESH_EXPORT NumericVector< Number >
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:511
const NumericVector< Number > & get_vector(std::string_view vec_name) const
Definition: system.C:918
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.

◆ retrieve_timestep()

void libMesh::UnsteadySolver::retrieve_timestep ( )
overridevirtual

This method retrieves all the stored solutions at the current system.time.

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::SecondOrderUnsteadySolver, and libMesh::AdaptiveTimeSolver.

Definition at line 269 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::System::get_vector(), libMesh::NumericVector< T >::localize(), old_local_nonlinear_solution, libMesh::TimeSolver::solution_history, and libMesh::System::time.

Referenced by libMesh::Euler2Solver::integrate_adjoint_refinement_error_estimate(), libMesh::EulerSolver::integrate_adjoint_refinement_error_estimate(), integrate_adjoint_sensitivity(), libMesh::Euler2Solver::integrate_qoi_timestep(), and libMesh::EulerSolver::integrate_qoi_timestep().

270 {
271  // Retrieve all the stored vectors at the current time
272  solution_history->retrieve(false, _system.time);
273 
274  // Dont forget to localize the old_nonlinear_solution !
275  _system.get_vector("_old_nonlinear_solution").localize
278 }
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1595
std::unique_ptr< SolutionHistory > solution_history
A std::unique_ptr to a SolutionHistory object.
Definition: time_solver.h:319
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
std::shared_ptr< NumericVector< Number > > old_local_nonlinear_solution
Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be...
const DofMap & get_dof_map() const
Definition: system.h:2293
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:511
const NumericVector< Number > & get_vector(std::string_view vec_name) const
Definition: system.C:918
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.

◆ set_first_adjoint_step()

void libMesh::UnsteadySolver::set_first_adjoint_step ( bool  first_adjoint_step_setting)
inline

A setter for the first_adjoint_step boolean.

Needed for nested time solvers.

Definition at line 199 of file unsteady_solver.h.

References first_adjoint_step.

200  {
201  first_adjoint_step = first_adjoint_step_setting;
202  }
bool first_adjoint_step
A bool that will be true the first time adjoint_advance_timestep() is called, (when the primal soluti...

◆ set_first_solve()

void libMesh::UnsteadySolver::set_first_solve ( bool  first_solve_setting)
inline

Definition at line 209 of file unsteady_solver.h.

References first_solve.

210  {
211  first_solve = first_solve_setting;
212  }
bool first_solve
A bool that will be true the first time solve() is called, and false thereafter.

◆ set_is_adjoint()

void libMesh::TimeSolver::set_is_adjoint ( bool  _is_adjoint_value)
inlineinherited

Accessor for setting whether we need to do a primal or adjoint solve.

Definition at line 284 of file time_solver.h.

References libMesh::TimeSolver::_is_adjoint.

Referenced by libMesh::DifferentiableSystem::adjoint_solve(), libMesh::FEMSystem::postprocess(), and libMesh::DifferentiableSystem::solve().

285  { _is_adjoint = _is_adjoint_value; }
bool _is_adjoint
This boolean tells the TimeSolver whether we are solving a primal or adjoint problem.
Definition: time_solver.h:340

◆ set_solution_history()

void libMesh::TimeSolver::set_solution_history ( const SolutionHistory _solution_history)
inherited

A setter function users will employ if they need to do something other than save no solution history.

Definition at line 119 of file time_solver.C.

References libMesh::SolutionHistory::clone(), and libMesh::TimeSolver::solution_history.

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

120 {
121  solution_history = _solution_history.clone();
122 }
std::unique_ptr< SolutionHistory > solution_history
A std::unique_ptr to a SolutionHistory object.
Definition: time_solver.h:319

◆ side_residual()

virtual bool libMesh::TimeSolver::side_residual ( bool  request_jacobian,
DiffContext  
)
pure virtualinherited

This method uses the DifferentiablePhysics side_time_derivative(), side_constraint(), and side_mass_residual() to build a full residual on an element's side.

What combination it uses will depend on the type of solver. See the subclasses for more details.

Implemented in libMesh::NewmarkSolver, libMesh::EigenTimeSolver, libMesh::AdaptiveTimeSolver, libMesh::Euler2Solver, libMesh::SteadySolver, and libMesh::EulerSolver.

Referenced by libMesh::FEMSystem::numerical_side_jacobian().

◆ solve()

void libMesh::UnsteadySolver::solve ( )
overridevirtual

This method solves for the solution at the next timestep.

Usually we will only need to solve one (non)linear system per timestep, but more complex subclasses may override this.

Reimplemented from libMesh::TimeSolver.

Reimplemented in libMesh::NewmarkSolver, libMesh::AdaptiveTimeSolver, and libMesh::TwostepTimeSolver.

Definition at line 127 of file unsteady_solver.C.

References libMesh::TimeSolver::_diff_solver, libMesh::TimeSolver::_system, advance_timestep(), libMesh::DifferentiableSystem::deltat, libMesh::DiffSolver::DIVERGED_BACKTRACKING_FAILURE, libMesh::DiffSolver::DIVERGED_MAX_NONLINEAR_ITERATIONS, first_solve, libMesh::TimeSolver::last_deltat, libMesh::out, libMesh::TimeSolver::quiet, and libMesh::TimeSolver::reduce_deltat_on_diffsolver_failure.

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

128 {
129  if (first_solve)
130  {
132  first_solve = false;
133  }
134 
135  unsigned int solve_result = _diff_solver->solve();
136 
137  // If we requested the UnsteadySolver to attempt reducing dt after a
138  // failed DiffSolver solve, check the results of the solve now.
140  {
141  bool backtracking_failed =
143 
144  bool max_iterations =
146 
147  if (backtracking_failed || max_iterations)
148  {
149  // Cut timestep in half
150  for (unsigned int nr=0; nr<reduce_deltat_on_diffsolver_failure; ++nr)
151  {
152  _system.deltat *= 0.5;
153  libMesh::out << "Newton backtracking failed. Trying with smaller timestep, dt="
154  << _system.deltat << std::endl;
155 
156  solve_result = _diff_solver->solve();
157 
158  // Check solve results with reduced timestep
159  bool backtracking_still_failed =
161 
162  bool backtracking_max_iterations =
164 
165  if (!backtracking_still_failed && !backtracking_max_iterations)
166  {
167  // Set the successful deltat as the last deltat
169 
170  if (!quiet)
171  libMesh::out << "Reduced dt solve succeeded." << std::endl;
172  return;
173  }
174  }
175 
176  // If we made it here, we still couldn't converge the solve after
177  // reducing deltat
178  libMesh::out << "DiffSolver::solve() did not succeed after "
180  << " attempts." << std::endl;
181  libmesh_convergence_failure();
182 
183  } // end if (backtracking_failed || max_iterations)
184  } // end if (reduce_deltat_on_diffsolver_failure)
185 
186  // Set the successful deltat as the last deltat
188 }
bool quiet
Print extra debugging information if quiet == false.
Definition: time_solver.h:230
std::unique_ptr< DiffSolver > _diff_solver
An implicit linear or nonlinear solver to use at each timestep.
Definition: time_solver.h:302
The DiffSolver reached the maximum allowed number of nonlinear iterations before satisfying any conve...
Definition: diff_solver.h:268
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
Definition: diff_system.h:278
unsigned int reduce_deltat_on_diffsolver_failure
This value (which defaults to zero) is the number of times the TimeSolver is allowed to halve deltat ...
Definition: time_solver.h:259
virtual void advance_timestep() override
This method advances the solution to the next timestep, after a solve() has been performed.
OStreamProxy out
The DiffSolver failed to find a descent direction by backtracking (See newton_solver.C)
Definition: diff_solver.h:274
Real last_deltat
The deltat for the last completed timestep before the current one.
Definition: time_solver.h:332
bool first_solve
A bool that will be true the first time solve() is called, and false thereafter.

◆ system() [1/2]

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

Definition at line 210 of file time_solver.h.

References libMesh::TimeSolver::_system.

Referenced by libMesh::TimeSolver::adjoint_solve(), libMesh::TimeSolver::reinit(), and libMesh::TimeSolver::solve().

210 { return _system; }
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312

◆ system() [2/2]

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

Definition at line 215 of file time_solver.h.

References libMesh::TimeSolver::_system.

215 { return _system; }
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312

◆ time_order()

virtual unsigned int libMesh::UnsteadySolver::time_order ( ) const
pure virtual

◆ update()

void libMesh::UnsteadySolver::update ( )

Definition at line 361 of file unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::System::get_vector(), libMesh::NumericVector< T >::localize(), and old_local_nonlinear_solution.

362 {
363  // Dont forget to localize the old_nonlinear_solution !
364  _system.get_vector("_old_nonlinear_solution").localize
367 }
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
std::shared_ptr< NumericVector< Number > > old_local_nonlinear_solution
Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be...
const DofMap & get_dof_map() const
Definition: system.h:2293
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:511
const NumericVector< Number > & get_vector(std::string_view vec_name) const
Definition: system.C:918
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.

Member Data Documentation

◆ _counts

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

Actually holds the data.

Definition at line 124 of file reference_counter.h.

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

◆ _diff_solver

std::unique_ptr<DiffSolver> libMesh::TimeSolver::_diff_solver
protectedinherited

An implicit linear or nonlinear solver to use at each timestep.

Definition at line 302 of file time_solver.h.

Referenced by libMesh::NewmarkSolver::compute_initial_accel(), libMesh::TimeSolver::diff_solver(), and solve().

◆ _enable_print_counter

bool libMesh::ReferenceCounter::_enable_print_counter = true
staticprotectedinherited

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

Definition at line 143 of file reference_counter.h.

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

◆ _linear_solver

std::unique_ptr<LinearSolver<Number> > libMesh::TimeSolver::_linear_solver
protectedinherited

An implicit linear solver to use for adjoint problems.

Definition at line 307 of file time_solver.h.

Referenced by libMesh::TimeSolver::linear_solver(), and libMesh::TimeSolver::reinit().

◆ _mutex

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

Mutual exclusion object to enable thread-safe reference counting.

Definition at line 137 of file reference_counter.h.

◆ _n_objects

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

The number of objects.

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

Definition at line 132 of file reference_counter.h.

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

◆ _system

sys_type& libMesh::TimeSolver::_system
protectedinherited

A reference to the system we are solving.

Definition at line 312 of file time_solver.h.

Referenced by libMesh::EulerSolver::_general_residual(), libMesh::Euler2Solver::_general_residual(), libMesh::SteadySolver::_general_residual(), libMesh::NewmarkSolver::_general_residual(), libMesh::AdaptiveTimeSolver::adjoint_advance_timestep(), adjoint_advance_timestep(), libMesh::TwostepTimeSolver::adjoint_solve(), adjoint_solve(), libMesh::TimeSolver::adjoint_solve(), libMesh::NewmarkSolver::advance_timestep(), libMesh::AdaptiveTimeSolver::advance_timestep(), advance_timestep(), libMesh::NewmarkSolver::compute_initial_accel(), libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), du(), libMesh::EulerSolver::element_residual(), libMesh::Euler2Solver::element_residual(), libMesh::EigenTimeSolver::element_residual(), libMesh::SecondOrderUnsteadySolver::init(), init(), libMesh::TimeSolver::init(), libMesh::EigenTimeSolver::init(), init_adjoints(), libMesh::TimeSolver::init_adjoints(), libMesh::SecondOrderUnsteadySolver::init_data(), init_data(), libMesh::TimeSolver::init_data(), libMesh::Euler2Solver::integrate_adjoint_refinement_error_estimate(), libMesh::TwostepTimeSolver::integrate_adjoint_refinement_error_estimate(), libMesh::EulerSolver::integrate_adjoint_refinement_error_estimate(), libMesh::TwostepTimeSolver::integrate_adjoint_sensitivity(), libMesh::SteadySolver::integrate_adjoint_sensitivity(), integrate_adjoint_sensitivity(), libMesh::Euler2Solver::integrate_qoi_timestep(), libMesh::TwostepTimeSolver::integrate_qoi_timestep(), libMesh::EulerSolver::integrate_qoi_timestep(), libMesh::SteadySolver::integrate_qoi_timestep(), libMesh::EulerSolver::nonlocal_residual(), libMesh::Euler2Solver::nonlocal_residual(), libMesh::EigenTimeSolver::nonlocal_residual(), old_nonlinear_solution(), libMesh::SecondOrderUnsteadySolver::old_solution_accel(), libMesh::SecondOrderUnsteadySolver::old_solution_rate(), libMesh::NewmarkSolver::project_initial_accel(), libMesh::SecondOrderUnsteadySolver::project_initial_rate(), libMesh::SecondOrderUnsteadySolver::reinit(), reinit(), libMesh::TimeSolver::reinit(), retrieve_timestep(), libMesh::EigenTimeSolver::side_residual(), libMesh::TwostepTimeSolver::solve(), solve(), libMesh::EigenTimeSolver::solve(), libMesh::TimeSolver::system(), and update().

◆ first_adjoint_step

bool libMesh::UnsteadySolver::first_adjoint_step
protected

A bool that will be true the first time adjoint_advance_timestep() is called, (when the primal solution is to be used to set adjoint boundary conditions) and false thereafter.

Definition at line 226 of file unsteady_solver.h.

Referenced by libMesh::AdaptiveTimeSolver::adjoint_advance_timestep(), and set_first_adjoint_step().

◆ first_solve

bool libMesh::UnsteadySolver::first_solve
protected

A bool that will be true the first time solve() is called, and false thereafter.

Definition at line 220 of file unsteady_solver.h.

Referenced by libMesh::NewmarkSolver::advance_timestep(), libMesh::AdaptiveTimeSolver::advance_timestep(), advance_timestep(), set_first_solve(), libMesh::TwostepTimeSolver::solve(), and solve().

◆ last_deltat

Real libMesh::TimeSolver::last_deltat
protectedinherited

The deltat for the last completed timestep before the current one.

Definition at line 332 of file time_solver.h.

Referenced by libMesh::TwostepTimeSolver::adjoint_solve(), adjoint_solve(), libMesh::TimeSolver::last_completed_timestep_size(), libMesh::TwostepTimeSolver::solve(), and solve().

◆ last_step_deltat

Real libMesh::UnsteadySolver::last_step_deltat
protected

We will need to move the system.time around to ensure that residuals are built with the right deltat and the right time.

Definition at line 237 of file unsteady_solver.h.

Referenced by libMesh::Euler2Solver::integrate_adjoint_refinement_error_estimate(), and libMesh::EulerSolver::integrate_adjoint_refinement_error_estimate().

◆ next_step_deltat

Real libMesh::UnsteadySolver::next_step_deltat
protected

◆ old_adjoints

std::vector< std::unique_ptr<NumericVector<Number> > > libMesh::UnsteadySolver::old_adjoints
protected

A vector of pointers to vectors holding the adjoint solution at the last time step.

Definition at line 231 of file unsteady_solver.h.

Referenced by libMesh::Euler2Solver::integrate_adjoint_refinement_error_estimate(), libMesh::EulerSolver::integrate_adjoint_refinement_error_estimate(), and UnsteadySolver().

◆ old_local_nonlinear_solution

std::shared_ptr<NumericVector<Number> > libMesh::UnsteadySolver::old_local_nonlinear_solution

Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be shared between different derived class instances, as in e.g.

AdaptiveTimeSolver.

Definition at line 178 of file unsteady_solver.h.

Referenced by libMesh::AdaptiveTimeSolver::adjoint_advance_timestep(), adjoint_advance_timestep(), libMesh::AdaptiveTimeSolver::advance_timestep(), advance_timestep(), libMesh::AdaptiveTimeSolver::init(), init_data(), old_nonlinear_solution(), reinit(), retrieve_timestep(), and update().

◆ quiet

bool libMesh::TimeSolver::quiet
inherited

Print extra debugging information if quiet == false.

Definition at line 230 of file time_solver.h.

Referenced by libMesh::TwostepTimeSolver::solve(), solve(), and libMesh::EigenTimeSolver::solve().

◆ reduce_deltat_on_diffsolver_failure

unsigned int libMesh::TimeSolver::reduce_deltat_on_diffsolver_failure
inherited

This value (which defaults to zero) is the number of times the TimeSolver is allowed to halve deltat and let the DiffSolver repeat the latest failed solve with a reduced timestep.

Note
This has no effect for SteadySolvers.
You must set at least one of the DiffSolver flags "continue_after_max_iterations" or "continue_after_backtrack_failure" to allow the TimeSolver to retry the solve.

Definition at line 259 of file time_solver.h.

Referenced by libMesh::TwostepTimeSolver::solve(), and solve().

◆ solution_history

std::unique_ptr<SolutionHistory> libMesh::TimeSolver::solution_history
protectedinherited

A std::unique_ptr to a SolutionHistory object.

Default is NoSolutionHistory, which the user can override by declaring a different kind of SolutionHistory in the application

Definition at line 319 of file time_solver.h.

Referenced by adjoint_advance_timestep(), advance_timestep(), libMesh::TimeSolver::get_solution_history(), retrieve_timestep(), and libMesh::TimeSolver::set_solution_history().


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