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::EulerSolver Class Reference

This class defines a theta-method Euler (defaulting to Backward Euler with theta = 1.0) solver to handle time integration of DifferentiableSystems. More...

#include <euler_solver.h>

Inheritance diagram for libMesh::EulerSolver:
[legend]

Public Types

typedef FirstOrderUnsteadySolver Parent
 The parent class. More...
 
typedef DifferentiableSystem sys_type
 The type of system. More...
 

Public Member Functions

 EulerSolver (sys_type &s)
 Constructor. More...
 
virtual ~EulerSolver ()
 Destructor. More...
 
virtual Real error_order () const override
 Error convergence order: 2 for Crank-Nicolson, 1 otherwise. More...
 
virtual bool element_residual (bool request_jacobian, DiffContext &) override
 This method uses the DifferentiablePhysics' element_time_derivative() and element_constraint() to build a full residual on an element. More...
 
virtual bool side_residual (bool request_jacobian, DiffContext &) override
 This method uses the DifferentiablePhysics' side_time_derivative() and side_constraint() to build a full residual on an element's side. More...
 
virtual bool nonlocal_residual (bool request_jacobian, DiffContext &) override
 This method uses the DifferentiablePhysics' nonlocal_time_derivative() and nonlocal_constraint() to build a full residual for non-local terms. More...
 
virtual void integrate_qoi_timestep () override
 A method to integrate the system::QoI functionals. More...
 
virtual void integrate_adjoint_refinement_error_estimate (AdjointRefinementEstimator &adjoint_refinement_error_estimator, ErrorVector &QoI_elementwise_error) override
 A method to compute the adjoint refinement error estimate at the current timestep. More...
 
virtual unsigned int time_order () const override
 
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_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...
 
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 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

Real theta
 The value for the theta method to employ: 1.0 corresponds to backwards Euler, 0.0 corresponds to forwards Euler, 0.5 corresponds to Crank-Nicolson. More...
 
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

virtual bool _general_residual (bool request_jacobian, DiffContext &, ResFuncType mass, ResFuncType damping, ResFuncType time_deriv, ResFuncType constraint, ReinitFuncType reinit, bool compute_second_order_eqns)
 This method is the underlying implementation of the public residual methods. More...
 
void prepare_accel (DiffContext &context)
 If there are second order variables in the system, then we also prepare the accel for those variables so the user can treat them as such. More...
 
bool compute_second_order_eqns (bool compute_jacobian, DiffContext &c)
 If there are second order variables, then we need to compute their residual equations and corresponding Jacobian. More...
 
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 class defines a theta-method Euler (defaulting to Backward Euler with theta = 1.0) solver to handle time integration of DifferentiableSystems.

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

Definition at line 43 of file euler_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.

◆ Parent

The parent class.

Definition at line 49 of file euler_solver.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

◆ EulerSolver()

libMesh::EulerSolver::EulerSolver ( sys_type s)
explicit

Constructor.

Requires a reference to the system to be solved.

Definition at line 32 of file euler_solver.C.

34 {
35 }
Real theta
The value for the theta method to employ: 1.0 corresponds to backwards Euler, 0.0 corresponds to forw...
Definition: euler_solver.h:117
FirstOrderUnsteadySolver(sys_type &s)
Constructor.

◆ ~EulerSolver()

libMesh::EulerSolver::~EulerSolver ( )
virtualdefault

Destructor.

Member Function Documentation

◆ _general_residual()

bool libMesh::EulerSolver::_general_residual ( bool  request_jacobian,
DiffContext context,
ResFuncType  mass,
ResFuncType  damping,
ResFuncType  time_deriv,
ResFuncType  constraint,
ReinitFuncType  reinit,
bool  compute_second_order_eqns 
)
protectedvirtual

This method is the underlying implementation of the public residual methods.

Definition at line 102 of file euler_solver.C.

References libMesh::TimeSolver::_system, libMesh::DenseVector< T >::add(), libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), libMesh::DifferentiableSystem::deltat, libMesh::DiffContext::elem_solution_derivative, libMesh::DiffContext::elem_solution_rate_derivative, libMesh::DiffContext::fixed_solution_derivative, libMesh::DiffContext::get_dof_indices(), libMesh::DiffContext::get_elem_fixed_solution(), libMesh::DiffContext::get_elem_jacobian(), libMesh::DiffContext::get_elem_solution(), libMesh::DiffContext::get_elem_solution_rate(), libMesh::DifferentiableSystem::get_physics(), libMesh::UnsteadySolver::old_nonlinear_solution(), libMesh::FirstOrderUnsteadySolver::prepare_accel(), libMesh::DenseVector< T >::size(), libMesh::DenseVector< T >::swap(), libMesh::DenseMatrix< T >::swap(), theta, and libMesh::System::use_fixed_solution.

Referenced by element_residual(), nonlocal_residual(), and side_residual().

110 {
111  unsigned int n_dofs = context.get_elem_solution().size();
112 
113  // We might need to save the old jacobian in case one of our physics
114  // terms later is unable to update it analytically.
115  DenseMatrix<Number> old_elem_jacobian(n_dofs, n_dofs);
116  if (request_jacobian)
117  old_elem_jacobian.swap(context.get_elem_jacobian());
118 
119  // Local nonlinear solution at old timestep
120  DenseVector<Number> old_elem_solution(n_dofs);
121  for (unsigned int i=0; i != n_dofs; ++i)
122  old_elem_solution(i) =
123  old_nonlinear_solution(context.get_dof_indices()[i]);
124 
125  // Local time derivative of solution
126  context.get_elem_solution_rate() = context.get_elem_solution();
127  context.get_elem_solution_rate() -= old_elem_solution;
128  context.elem_solution_rate_derivative = 1 / _system.deltat;
129  context.get_elem_solution_rate() *=
130  context.elem_solution_rate_derivative;
131 
132  // If we are asked to compute residuals for second order variables,
133  // we also populate the acceleration part so the user can use that.
135  this->prepare_accel(context);
136 
137  // Local nonlinear solution at time t_theta
138  DenseVector<Number> theta_solution(context.get_elem_solution());
139  theta_solution *= theta;
140  theta_solution.add(1. - theta, old_elem_solution);
141 
142  context.elem_solution_derivative = theta;
143  context.fixed_solution_derivative = theta;
144 
145  // If a fixed solution is requested, we'll use theta_solution
147  context.get_elem_fixed_solution() = theta_solution;
148 
149  // Move theta_->elem_, elem_->theta_
150  context.get_elem_solution().swap(theta_solution);
151 
152  // Move the mesh into place first if necessary, set t = t_{\theta}
153  (context.*reinit_func)(theta);
154 
155  // Get the time derivative at t_theta
156  bool jacobian_computed =
157  (_system.get_physics()->*time_deriv)(request_jacobian, context);
158 
159  jacobian_computed = (_system.get_physics()->*mass)(jacobian_computed, context) &&
160  jacobian_computed;
161 
162  // If we have second-order variables, we need to get damping terms
163  // and the velocity equations
165  {
166  jacobian_computed = (_system.get_physics()->*damping)(jacobian_computed, context) &&
167  jacobian_computed;
168 
169  jacobian_computed = this->compute_second_order_eqns(jacobian_computed, context) &&
170  jacobian_computed;
171  }
172 
173  // Restore the elem position if necessary, set t = t_{n+1}
174  (context.*reinit_func)(1);
175 
176  // Move elem_->elem_, theta_->theta_
177  context.get_elem_solution().swap(theta_solution);
178  context.elem_solution_derivative = 1;
179 
180  // Add the constraint term
181  jacobian_computed = (_system.get_physics()->*constraint)(jacobian_computed, context) &&
182  jacobian_computed;
183 
184  // Add back (or restore) the old jacobian
185  if (request_jacobian)
186  {
187  if (jacobian_computed)
188  context.get_elem_jacobian() += old_elem_jacobian;
189  else
190  context.get_elem_jacobian().swap(old_elem_jacobian);
191  }
192 
193  return jacobian_computed;
194 }
bool compute_second_order_eqns(bool compute_jacobian, DiffContext &c)
If there are second order variables, then we need to compute their residual equations and correspondi...
void prepare_accel(DiffContext &context)
If there are second order variables in the system, then we also prepare the accel for those variables...
Real theta
The value for the theta method to employ: 1.0 corresponds to backwards Euler, 0.0 corresponds to forw...
Definition: euler_solver.h:117
Number old_nonlinear_solution(const dof_id_type global_dof_number) const
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
bool use_fixed_solution
A boolean to be set to true by systems using elem_fixed_solution, for optional use by e...
Definition: system.h:1543
const DifferentiablePhysics * get_physics() const
Definition: diff_system.h:181

◆ adjoint_advance_timestep()

void libMesh::UnsteadySolver::adjoint_advance_timestep ( )
overridevirtualinherited

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(), libMesh::UnsteadySolver::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)
overridevirtualinherited

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

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, libMesh::UnsteadySolver::first_solve, libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::System::get_vector(), libMesh::NumericVector< T >::localize(), libMesh::UnsteadySolver::old_local_nonlinear_solution, libMesh::System::solution, libMesh::TimeSolver::solution_history, and libMesh::System::time.

Referenced by libMesh::NewmarkSolver::advance_timestep(), and libMesh::UnsteadySolver::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 {}

◆ compute_second_order_eqns()

bool libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns ( bool  compute_jacobian,
DiffContext c 
)
protectedinherited

If there are second order variables, then we need to compute their residual equations and corresponding Jacobian.

The residual equation will simply be \( \dot{u} - v = 0 \), where \( u \) is the second order variable add by the user and \( v \) is the variable added by the time-solver as the "velocity" variable.

Definition at line 32 of file first_order_unsteady_solver.C.

References libMesh::TimeSolver::_system, compute_jacobian(), libMesh::DiffContext::get_dof_indices(), libMesh::DiffContext::get_elem_jacobian(), libMesh::DiffContext::get_elem_residual(), libMesh::DiffContext::get_elem_solution_derivative(), libMesh::DiffContext::get_elem_solution_rate_derivative(), libMesh::FEMContext::get_element_fe(), libMesh::FEMContext::get_element_qrule(), libMesh::DifferentiableSystem::get_second_order_dot_var(), libMesh::DiffContext::get_system(), libMesh::FEMContext::interior_rate(), libMesh::FEMContext::interior_value(), libMesh::DifferentiablePhysics::is_second_order_var(), libMesh::libmesh_assert(), libMesh::make_range(), libMesh::QBase::n_points(), libMesh::DiffContext::n_vars(), libMesh::Variable::type(), and libMesh::System::variable().

Referenced by _general_residual(), libMesh::Euler2Solver::_general_residual(), element_residual(), libMesh::Euler2Solver::element_residual(), nonlocal_residual(), and libMesh::Euler2Solver::nonlocal_residual().

33 {
34  FEMContext & context = cast_ref<FEMContext &>(c);
35 
36  unsigned int n_qpoints = context.get_element_qrule().n_points();
37 
38  for (auto var : make_range(context.n_vars()))
39  {
40  if (!this->_system.is_second_order_var(var))
41  continue;
42 
43  unsigned int dot_var = this->_system.get_second_order_dot_var(var);
44 
45  // We're assuming that the FE space for var and dot_var are the same
46  libmesh_assert( context.get_system().variable(var).type() ==
47  context.get_system().variable(dot_var).type() );
48 
49  FEBase * elem_fe = nullptr;
50  context.get_element_fe( var, elem_fe );
51 
52  const std::vector<Real> & JxW = elem_fe->get_JxW();
53 
54  const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
55 
56  const unsigned int n_dofs = cast_int<unsigned int>
57  (context.get_dof_indices(dot_var).size());
58 
59  DenseSubVector<Number> & Fu = context.get_elem_residual(var);
60  DenseSubMatrix<Number> & Kuu = context.get_elem_jacobian( var, var );
61  DenseSubMatrix<Number> & Kuv = context.get_elem_jacobian( var, dot_var );
62 
63  for (unsigned int qp = 0; qp != n_qpoints; ++qp)
64  {
65  Number udot, v;
66  context.interior_rate(var, qp, udot);
67  context.interior_value(dot_var, qp, v);
68 
69  for (unsigned int i = 0; i < n_dofs; i++)
70  {
71  Fu(i) += JxW[qp]*(udot-v)*phi[i][qp];
72 
73  if (compute_jacobian)
74  {
75  Number rate_factor = JxW[qp]*context.get_elem_solution_rate_derivative()*phi[i][qp];
76  Number soln_factor = JxW[qp]*context.get_elem_solution_derivative()*phi[i][qp];
77 
78  Kuu(i,i) += rate_factor*phi[i][qp];
79  Kuv(i,i) -= soln_factor*phi[i][qp];
80 
81  for (unsigned int j = i+1; j < n_dofs; j++)
82  {
83  Kuu(i,j) += rate_factor*phi[j][qp];
84  Kuu(j,i) += rate_factor*phi[j][qp];
85 
86  Kuv(i,j) -= soln_factor*phi[j][qp];
87  Kuv(j,i) -= soln_factor*phi[j][qp];
88  }
89  }
90  }
91  }
92  }
93 
94  return compute_jacobian;
95 }
bool is_second_order_var(unsigned int var) const
Definition: diff_physics.h:522
void compute_jacobian(const NumericVector< Number > &, SparseMatrix< Number > &J, NonlinearImplicitSystem &system)
Definition: assembly.C:315
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
unsigned int get_second_order_dot_var(unsigned int var) const
For a given second order (in time) variable var, this method will return the index to the correspondi...
Definition: diff_system.C:306
libmesh_assert(ctx)
FEGenericBase< Real > FEBase
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

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

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

bool libMesh::EulerSolver::element_residual ( bool  request_jacobian,
DiffContext context 
)
overridevirtual

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

What combination it uses will depend on theta.

Implements libMesh::TimeSolver.

Definition at line 53 of file euler_solver.C.

References libMesh::DifferentiablePhysics::_eulerian_time_deriv(), _general_residual(), libMesh::TimeSolver::_system, libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), libMesh::DifferentiablePhysics::damping_residual(), libMesh::DiffContext::elem_reinit(), libMesh::DifferentiablePhysics::element_constraint(), libMesh::DifferentiableSystem::get_physics(), libMesh::DifferentiablePhysics::have_second_order_vars(), and libMesh::DifferentiablePhysics::mass_residual().

55 {
57 
58  return this->_general_residual(request_jacobian,
59  context,
66 }
bool compute_second_order_eqns(bool compute_jacobian, DiffContext &c)
If there are second order variables, then we need to compute their residual equations and correspondi...
virtual bool damping_residual(bool request_jacobian, DiffContext &)
Subtracts a damping vector contribution on elem from elem_residual.
Definition: diff_physics.h:360
virtual bool _general_residual(bool request_jacobian, DiffContext &, ResFuncType mass, ResFuncType damping, ResFuncType time_deriv, ResFuncType constraint, ReinitFuncType reinit, bool compute_second_order_eqns)
This method is the underlying implementation of the public residual methods.
Definition: euler_solver.C:102
virtual void elem_reinit(Real)
Gives derived classes the opportunity to reinitialize data (FE objects in FEMSystem, for example) needed for an interior integration at a new point within a timestep.
Definition: diff_context.h:77
bool _eulerian_time_deriv(bool request_jacobian, DiffContext &)
This method simply combines element_time_derivative() and eulerian_residual(), which makes its addres...
Definition: diff_physics.C:96
virtual bool mass_residual(bool request_jacobian, DiffContext &)
Subtracts a mass vector contribution on elem from elem_residual.
Definition: diff_physics.h:302
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
const DifferentiablePhysics * get_physics() const
Definition: diff_system.h:181
virtual bool element_constraint(bool request_jacobian, DiffContext &)
Adds the constraint contribution on elem to elem_residual.
Definition: diff_physics.h:144

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

Real libMesh::EulerSolver::error_order ( ) const
overridevirtual

Error convergence order: 2 for Crank-Nicolson, 1 otherwise.

Implements libMesh::UnsteadySolver.

Definition at line 43 of file euler_solver.C.

References theta.

44 {
45  if (theta == 0.5)
46  return 2.;
47  return 1.;
48 }
Real theta
The value for the theta method to employ: 1.0 corresponds to backwards Euler, 0.0 corresponds to forw...
Definition: euler_solver.h:117

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

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

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

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(), libMesh::UnsteadySolver::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::EulerSolver::integrate_adjoint_refinement_error_estimate ( AdjointRefinementEstimator adjoint_refinement_error_estimator,
ErrorVector QoI_elementwise_error 
)
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. The integration scheme used is consistent with the theta used for the unsteady solution. CURRENTLY ONLY SUPPORTED for Backward Euler.

Implements libMesh::FirstOrderUnsteadySolver.

Definition at line 254 of file euler_solver.C.

References libMesh::TimeSolver::_system, std::abs(), libMesh::NumericVector< T >::clone(), libMesh::DifferentiableSystem::deltat, libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::System::get_adjoint_solution(), libMesh::AdjointRefinementEstimator::get_global_QoI_error_estimate(), libMesh::QoISet::has_index(), libMesh::index_range(), libMesh::UnsteadySolver::last_step_deltat, libMesh::make_range(), libMesh::System::n_qois(), libMesh::UnsteadySolver::next_step_deltat, libMesh::UnsteadySolver::old_adjoints, libMesh::AdjointRefinementEstimator::qoi_set(), libMesh::Real, libMesh::UnsteadySolver::retrieve_timestep(), libMesh::System::set_qoi_error_estimate(), std::sqrt(), libMesh::NumericVector< T >::swap(), theta, libMesh::System::time, libMesh::TOLERANCE, and libMesh::System::update().

255 {
256  // Currently, we only support this functionality when Backward-Euler time integration is used.
257  if (theta != 1.0)
258  libmesh_not_implemented();
259 
260  // There are two possibilities regarding the integration rule we need to use for time integration.
261  // If we have a instantaneous QoI, then we need to use a left sided Riemann sum, otherwise the trapezoidal rule for temporally smooth QoIs.
262 
263  // Create left and right error estimate vectors of the right size
264  std::vector<Number> qoi_error_estimates_left(_system.n_qois());
265  std::vector<Number> qoi_error_estimates_right(_system.n_qois());
266 
267  // Get t_j
268  Real time_left = _system.time;
269 
270  // Get f(t_j)
271  ErrorVector QoI_elementwise_error_left;
272 
273  // If we are at the very initial step, the error contribution is zero,
274  // otherwise the old ajoint vector has been filled and we are the left end
275  // of a subsequent timestep or sub-timestep
276  if(old_adjoints[0] != nullptr)
277  {
278  // For evaluating the residual, we need to use the deltat that was used
279  // to get us to this solution, so we save the current deltat as next_step_deltat
280  // and set _system.deltat to the last completed deltat.
283 
284  // The adjoint error estimate expression for a backwards facing step
285  // scheme needs the adjoint for the last time instant, so save the current adjoint for future use
286  for (auto j : make_range(_system.n_qois()))
287  {
288  // Swap for residual weighting
290  }
291 
292  _system.update();
293 
294  // The residual has to be evaluated at the last time
296 
297  adjoint_refinement_error_estimator.estimate_error(_system, QoI_elementwise_error_left);
298 
299  // Shift the time back
301 
302  // Swap back the current and old adjoints
303  for (auto j : make_range(_system.n_qois()))
304  {
306  }
307 
308  // Set the system deltat back to what it should be to march to the next time
310 
311  }
312  else
313  {
314  for (auto i : index_range(QoI_elementwise_error))
315  QoI_elementwise_error_left[i] = 0.0;
316  }
317 
318  // Also get the left side contributions for the spatially integrated errors for all the QoIs in the QoI set
319  for (auto j : make_range(_system.n_qois()))
320  {
321  // Skip this QoI if not in the QoI Set
322  if (adjoint_refinement_error_estimator.qoi_set().has_index(j))
323  {
324  // If we are at the initial time, the error contribution is zero
326  {
327  qoi_error_estimates_left[j] = adjoint_refinement_error_estimator.get_global_QoI_error_estimate(j);
328  }
329  else
330  {
331  qoi_error_estimates_left[j] = 0.0;
332  }
333  }
334  }
335 
336  // Advance to t_j+1
338 
339  // Get t_j+1
340  Real time_right = _system.time;
341 
342  // We will need to use the last step deltat for the weighted residual evaluation
344 
345  // The adjoint error estimate expression for a backwards facing step
346  // scheme needs the adjoint for the last time instant, so save the current adjoint for future use
347  for (auto j : make_range(_system.n_qois()))
348  {
350  }
351 
352  // Retrieve the state and adjoint vectors for the next time instant
354 
355  // Swap for residual weighting
356  for (auto j : make_range(_system.n_qois()))
357  {
359  }
360 
361  // Swap out the deltats as we did for the left side
364 
365  // Get f(t_j+1)
366  ErrorVector QoI_elementwise_error_right;
367 
368  _system.update();
369 
370  // The residual has to be evaluated at the last time
372 
373  adjoint_refinement_error_estimator.estimate_error(_system, QoI_elementwise_error_right);
374 
375  // Shift the time back
377 
378  // Set the system deltat back to what it needs to be able to march to the next time
380 
381  // Swap back now that the residual weighting is done
382  for (auto j : make_range(_system.n_qois()))
383  {
385  }
386 
387  // Also get the right side contributions for the spatially integrated errors for all the QoIs in the QoI set
388  for (auto j : make_range(_system.n_qois()))
389  {
390  // Skip this QoI if not in the QoI Set
391  if (adjoint_refinement_error_estimator.qoi_set().has_index(j))
392  {
393  qoi_error_estimates_right[j] = adjoint_refinement_error_estimator.get_global_QoI_error_estimate(j);
394  }
395  }
396 
397  // Error contribution from this timestep
398  for (auto i : index_range(QoI_elementwise_error))
399  QoI_elementwise_error[i] = float(((QoI_elementwise_error_right[i] + QoI_elementwise_error_left[i])/2)
400  * (time_right - time_left));
401 
402  // QoI set spatially integrated errors contribution from this timestep
403  for (auto j : make_range(_system.n_qois()))
404  {
405  // Skip this QoI if not in the QoI Set
406  if (adjoint_refinement_error_estimator.qoi_set().has_index(j))
407  {
408  _system.set_qoi_error_estimate(j, ( (1.0 - theta)*qoi_error_estimates_left[j] + theta*qoi_error_estimates_right[j] )*last_step_deltat);
409  }
410  }
411 
412 }
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1595
static constexpr Real TOLERANCE
unsigned int n_qois() const
Number of currently active quantities of interest.
Definition: system.h:2516
Real theta
The value for the theta method to employ: 1.0 corresponds to backwards Euler, 0.0 corresponds to forw...
Definition: euler_solver.h:117
virtual std::unique_ptr< NumericVector< T > > clone() const =0
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:53
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:57
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
std::vector< std::unique_ptr< NumericVector< Number > > > old_adjoints
A vector of pointers to vectors holding the adjoint solution at the last time step.
virtual void retrieve_timestep() override
This method retrieves all the stored solutions at the current system.time.
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
Definition: diff_system.h:278
virtual void update()
Update the local values to reflect the solution on neighboring processors.
Definition: system.C:493
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void set_qoi_error_estimate(unsigned int qoi_index, Number qoi_error_estimate)
Definition: system.C:2354
virtual void swap(NumericVector< T > &v)
Swaps the contents of this with v.
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
Real last_step_deltat
We will need to move the system.time around to ensure that residuals are built with the right deltat ...
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

◆ integrate_adjoint_sensitivity()

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

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(), libMesh::UnsteadySolver::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::EulerSolver::integrate_qoi_timestep ( )
overridevirtual

A method to integrate the system::QoI functionals.

Implements libMesh::FirstOrderUnsteadySolver.

Definition at line 196 of file euler_solver.C.

References libMesh::TimeSolver::_system, libMesh::ExplicitSystem::assemble_qoi(), libMesh::DifferentiableSystem::deltat, libMesh::System::get_qoi_value(), libMesh::make_range(), libMesh::System::n_qois(), libMesh::UnsteadySolver::retrieve_timestep(), libMesh::System::set_qoi(), theta, and libMesh::System::time.

197 {
198  // We use a numerical integration scheme consistent with the theta used for the timesolver.
199 
200  // Zero out the system.qoi vector
201  for (auto j : make_range(_system.n_qois()))
202  {
203  _system.set_qoi(j, 0.0);
204  }
205 
206  // Left and right side contributions
207  std::vector<Number> left_contribution(_system.n_qois(), 0.0);
208  Number time_left = 0.0;
209  std::vector<Number> right_contribution(_system.n_qois(), 0.0);
210  Number time_right = 0.0;
211 
212  time_left = _system.time;
213 
214  // Base class assumes a direct steady evaluation
215  this->_system.assemble_qoi();
216 
217  // Also get the spatially integrated errors for all the QoIs in the QoI set
218  for (auto j : make_range(_system.n_qois()))
219  {
220  left_contribution[j] = _system.get_qoi_value(j);
221  }
222 
223  // Advance to t_j+1
225 
226  time_right = _system.time;
227 
228  // Load the solution at the next timestep
230 
231  // Zero out the system.qoi vector
232  for (auto j : make_range(_system.n_qois()))
233  {
234  _system.set_qoi(j, 0.0);
235  }
236 
237  // Base class assumes a direct steady evaluation
238  this->_system.assemble_qoi();
239 
240  for(auto j : make_range(_system.n_qois()))
241  {
242  right_contribution[j] = _system.get_qoi_value(j);
243  }
244 
245  // Combine the left and right side contributions as per the specified theta
246  // theta = 0.5 (Crank-Nicholson) gives the trapezoidal rule.
247  for (auto j : make_range(_system.n_qois()))
248  {
249  _system.set_qoi(j, ( ((1.0 - theta)*left_contribution[j]) + (theta*right_contribution[j]) )*(time_right - time_left));
250  }
251 }
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1595
Number get_qoi_value(unsigned int qoi_index) const
Definition: system.C:2334
unsigned int n_qois() const
Number of currently active quantities of interest.
Definition: system.h:2516
Real theta
The value for the theta method to employ: 1.0 corresponds to backwards Euler, 0.0 corresponds to forw...
Definition: euler_solver.h:117
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
virtual void retrieve_timestep() override
This method retrieves all the stored solutions at the current system.time.
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
Definition: diff_system.h:278
virtual void assemble_qoi(const QoISet &qoi_indices=QoISet()) override
Prepares qoi for quantity of interest assembly, then calls user qoi function.
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
void set_qoi(unsigned int qoi_index, Number qoi_value)
Definition: system.C:2326

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

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

bool libMesh::EulerSolver::nonlocal_residual ( bool  request_jacobian,
DiffContext context 
)
overridevirtual

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

What combination it uses will depend on theta.

Implements libMesh::TimeSolver.

Definition at line 85 of file euler_solver.C.

References _general_residual(), libMesh::TimeSolver::_system, libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), libMesh::DifferentiableSystem::have_second_order_scalar_vars(), libMesh::DifferentiablePhysics::nonlocal_constraint(), libMesh::DifferentiablePhysics::nonlocal_damping_residual(), libMesh::DifferentiablePhysics::nonlocal_mass_residual(), libMesh::DiffContext::nonlocal_reinit(), and libMesh::DifferentiablePhysics::nonlocal_time_derivative().

87 {
89 
90  return this->_general_residual(request_jacobian,
91  context,
98 }
bool compute_second_order_eqns(bool compute_jacobian, DiffContext &c)
If there are second order variables, then we need to compute their residual equations and correspondi...
virtual bool nonlocal_time_derivative(bool request_jacobian, DiffContext &)
Adds any nonlocal time derivative contributions (e.g.
Definition: diff_physics.h:214
virtual bool _general_residual(bool request_jacobian, DiffContext &, ResFuncType mass, ResFuncType damping, ResFuncType time_deriv, ResFuncType constraint, ReinitFuncType reinit, bool compute_second_order_eqns)
This method is the underlying implementation of the public residual methods.
Definition: euler_solver.C:102
virtual bool nonlocal_constraint(bool request_jacobian, DiffContext &)
Adds any nonlocal constraint contributions (e.g.
Definition: diff_physics.h:233
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:312
bool have_second_order_scalar_vars() const
Check for any second order vars that are also belong to FEFamily::SCALAR.
Definition: diff_system.C:335
virtual bool nonlocal_damping_residual(bool request_jacobian, DiffContext &)
Subtracts any nonlocal damping vector contributions (e.g.
Definition: diff_physics.h:394
virtual bool nonlocal_mass_residual(bool request_jacobian, DiffContext &c)
Subtracts any nonlocal mass vector contributions (e.g.
Definition: diff_physics.C:57
virtual void nonlocal_reinit(Real)
Gives derived classes the opportunity to reinitialize data needed for nonlocal calculations at a new ...
Definition: diff_context.h:95

◆ old_nonlinear_solution()

Number libMesh::UnsteadySolver::old_nonlinear_solution ( const dof_id_type  global_dof_number) const
inherited
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 libMesh::UnsteadySolver::old_local_nonlinear_solution.

Referenced by _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

◆ prepare_accel()

void libMesh::FirstOrderUnsteadySolver::prepare_accel ( DiffContext context)
protectedinherited

If there are second order variables in the system, then we also prepare the accel for those variables so the user can treat them as such.

Definition at line 25 of file first_order_unsteady_solver.C.

References libMesh::DiffContext::elem_solution_accel_derivative, libMesh::DiffContext::get_elem_solution_accel(), libMesh::DiffContext::get_elem_solution_rate(), and libMesh::DiffContext::get_elem_solution_rate_derivative().

Referenced by _general_residual(), and libMesh::Euler2Solver::_general_residual().

26 {
27  context.get_elem_solution_accel() = context.get_elem_solution_rate();
28 
29  context.elem_solution_accel_derivative = context.get_elem_solution_rate_derivative();
30 }

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

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(), libMesh::UnsteadySolver::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 ( )
overridevirtualinherited

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(), libMesh::UnsteadySolver::old_local_nonlinear_solution, libMesh::TimeSolver::solution_history, and libMesh::System::time.

Referenced by libMesh::Euler2Solver::integrate_adjoint_refinement_error_estimate(), integrate_adjoint_refinement_error_estimate(), libMesh::UnsteadySolver::integrate_adjoint_sensitivity(), libMesh::Euler2Solver::integrate_qoi_timestep(), and 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)
inlineinherited

A setter for the first_adjoint_step boolean.

Needed for nested time solvers.

Definition at line 199 of file unsteady_solver.h.

References libMesh::UnsteadySolver::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)
inlineinherited

Definition at line 209 of file unsteady_solver.h.

References libMesh::UnsteadySolver::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()

bool libMesh::EulerSolver::side_residual ( bool  request_jacobian,
DiffContext context 
)
overridevirtual

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

What combination it uses will depend on theta.

Implements libMesh::TimeSolver.

Definition at line 70 of file euler_solver.C.

References _general_residual(), libMesh::DiffContext::elem_side_reinit(), libMesh::DifferentiablePhysics::side_constraint(), libMesh::DifferentiablePhysics::side_damping_residual(), libMesh::DifferentiablePhysics::side_mass_residual(), and libMesh::DifferentiablePhysics::side_time_derivative().

72 {
73  return this->_general_residual(request_jacobian,
74  context,
80  false);
81 }
virtual bool side_constraint(bool request_jacobian, DiffContext &)
Adds the constraint contribution on side of elem to elem_residual.
Definition: diff_physics.h:195
virtual bool side_damping_residual(bool request_jacobian, DiffContext &)
Subtracts a damping vector contribution on side of elem from elem_residual.
Definition: diff_physics.h:378
virtual bool _general_residual(bool request_jacobian, DiffContext &, ResFuncType mass, ResFuncType damping, ResFuncType time_deriv, ResFuncType constraint, ReinitFuncType reinit, bool compute_second_order_eqns)
This method is the underlying implementation of the public residual methods.
Definition: euler_solver.C:102
virtual void elem_side_reinit(Real)
Gives derived classes the opportunity to reinitialize data needed for a side integration at a new poi...
Definition: diff_context.h:83
virtual bool side_mass_residual(bool request_jacobian, DiffContext &)
Subtracts a mass vector contribution on side of elem from elem_residual.
Definition: diff_physics.h:320
virtual bool side_time_derivative(bool request_jacobian, DiffContext &)
Adds the time derivative contribution on side of elem to elem_residual.
Definition: diff_physics.h:174

◆ solve()

void libMesh::UnsteadySolver::solve ( )
overridevirtualinherited

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, libMesh::UnsteadySolver::advance_timestep(), libMesh::DifferentiableSystem::deltat, libMesh::DiffSolver::DIVERGED_BACKTRACKING_FAILURE, libMesh::DiffSolver::DIVERGED_MAX_NONLINEAR_ITERATIONS, libMesh::UnsteadySolver::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::FirstOrderUnsteadySolver::time_order ( ) const
inlineoverridevirtualinherited
Returns
The maximum order of time derivatives for which the UnsteadySolver subclass is capable of handling.

For example, EulerSolver will have time_order() = 1 and NewmarkSolver will have time_order() = 2.

Implements libMesh::UnsteadySolver.

Definition at line 90 of file first_order_unsteady_solver.h.

91  { return 1; }

◆ update()

void libMesh::UnsteadySolver::update ( )
inherited

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 libMesh::UnsteadySolver::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 libMesh::UnsteadySolver::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 _general_residual(), libMesh::Euler2Solver::_general_residual(), libMesh::SteadySolver::_general_residual(), libMesh::NewmarkSolver::_general_residual(), libMesh::AdaptiveTimeSolver::adjoint_advance_timestep(), libMesh::UnsteadySolver::adjoint_advance_timestep(), libMesh::TwostepTimeSolver::adjoint_solve(), libMesh::UnsteadySolver::adjoint_solve(), libMesh::TimeSolver::adjoint_solve(), libMesh::NewmarkSolver::advance_timestep(), libMesh::AdaptiveTimeSolver::advance_timestep(), libMesh::UnsteadySolver::advance_timestep(), libMesh::NewmarkSolver::compute_initial_accel(), libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), libMesh::UnsteadySolver::du(), element_residual(), libMesh::Euler2Solver::element_residual(), libMesh::EigenTimeSolver::element_residual(), libMesh::SecondOrderUnsteadySolver::init(), libMesh::UnsteadySolver::init(), libMesh::TimeSolver::init(), libMesh::EigenTimeSolver::init(), libMesh::UnsteadySolver::init_adjoints(), libMesh::TimeSolver::init_adjoints(), libMesh::SecondOrderUnsteadySolver::init_data(), libMesh::UnsteadySolver::init_data(), libMesh::TimeSolver::init_data(), libMesh::Euler2Solver::integrate_adjoint_refinement_error_estimate(), libMesh::TwostepTimeSolver::integrate_adjoint_refinement_error_estimate(), integrate_adjoint_refinement_error_estimate(), libMesh::TwostepTimeSolver::integrate_adjoint_sensitivity(), libMesh::SteadySolver::integrate_adjoint_sensitivity(), libMesh::UnsteadySolver::integrate_adjoint_sensitivity(), libMesh::Euler2Solver::integrate_qoi_timestep(), libMesh::TwostepTimeSolver::integrate_qoi_timestep(), integrate_qoi_timestep(), libMesh::SteadySolver::integrate_qoi_timestep(), nonlocal_residual(), libMesh::Euler2Solver::nonlocal_residual(), libMesh::EigenTimeSolver::nonlocal_residual(), libMesh::UnsteadySolver::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(), libMesh::UnsteadySolver::reinit(), libMesh::TimeSolver::reinit(), libMesh::UnsteadySolver::retrieve_timestep(), libMesh::EigenTimeSolver::side_residual(), libMesh::TwostepTimeSolver::solve(), libMesh::UnsteadySolver::solve(), libMesh::EigenTimeSolver::solve(), libMesh::TimeSolver::system(), and libMesh::UnsteadySolver::update().

◆ first_adjoint_step

bool libMesh::UnsteadySolver::first_adjoint_step
protectedinherited

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 libMesh::UnsteadySolver::set_first_adjoint_step().

◆ first_solve

bool libMesh::UnsteadySolver::first_solve
protectedinherited

◆ last_deltat

Real libMesh::TimeSolver::last_deltat
protectedinherited

◆ last_step_deltat

Real libMesh::UnsteadySolver::last_step_deltat
protectedinherited

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 integrate_adjoint_refinement_error_estimate().

◆ next_step_deltat

Real libMesh::UnsteadySolver::next_step_deltat
protectedinherited

◆ old_adjoints

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

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(), integrate_adjoint_refinement_error_estimate(), and libMesh::UnsteadySolver::UnsteadySolver().

◆ old_local_nonlinear_solution

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

◆ 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(), libMesh::UnsteadySolver::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 libMesh::UnsteadySolver::solve().

◆ solution_history

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

◆ theta

Real libMesh::EulerSolver::theta

The value for the theta method to employ: 1.0 corresponds to backwards Euler, 0.0 corresponds to forwards Euler, 0.5 corresponds to Crank-Nicolson.

Definition at line 117 of file euler_solver.h.

Referenced by _general_residual(), error_order(), integrate_adjoint_refinement_error_estimate(), and integrate_qoi_timestep().


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