libMesh::NewmarkSolver Class Reference

#include <newmark_solver.h>

Inheritance diagram for libMesh::NewmarkSolver:

Public Types

typedef UnsteadySolver Parent
 
typedef DifferentiableSystem sys_type
 

Public Member Functions

 NewmarkSolver (sys_type &s)
 
virtual ~NewmarkSolver ()
 
virtual void advance_timestep () override
 
virtual void adjoint_advance_timestep () override
 
virtual void compute_initial_accel ()
 
void project_initial_accel (FunctionBase< Number > *f, FunctionBase< Gradient > *g=nullptr)
 
void set_initial_accel_avail (bool initial_accel_set)
 
virtual Real error_order () const override
 
virtual void solve () override
 
virtual bool element_residual (bool request_jacobian, DiffContext &) override
 
virtual bool side_residual (bool request_jacobian, DiffContext &) override
 
virtual bool nonlocal_residual (bool request_jacobian, DiffContext &) override
 
void set_beta (Real beta)
 
void set_gamma (Real gamma)
 
virtual unsigned int time_order () const override
 
virtual void init () override
 
virtual void init_data () override
 
virtual void reinit () override
 
virtual void retrieve_timestep () override
 
void project_initial_rate (FunctionBase< Number > *f, FunctionBase< Gradient > *g=nullptr)
 
Number old_solution_rate (const dof_id_type global_dof_number) const
 
Number old_solution_accel (const dof_id_type global_dof_number) const
 
Number old_nonlinear_solution (const dof_id_type global_dof_number) const
 
virtual Real du (const SystemNorm &norm) const override
 
virtual bool is_steady () const override
 
virtual void before_timestep ()
 
const sys_typesystem () const
 
sys_typesystem ()
 
virtual std::unique_ptr< DiffSolver > & diff_solver ()
 
virtual std::unique_ptr< LinearSolver< Number > > & linear_solver ()
 
void set_solution_history (const SolutionHistory &_solution_history)
 
bool is_adjoint () const
 
void set_is_adjoint (bool _is_adjoint_value)
 

Static Public Member Functions

static std::string get_info ()
 
static void print_info (std::ostream &out=libMesh::out)
 
static unsigned int n_objects ()
 
static void enable_print_counter_info ()
 
static void disable_print_counter_info ()
 

Public Attributes

std::unique_ptr< NumericVector< Number > > old_local_nonlinear_solution
 
bool quiet
 
unsigned int reduce_deltat_on_diffsolver_failure
 

Protected Types

typedef bool(DifferentiablePhysics::* ResFuncType) (bool, DiffContext &)
 
typedef void(DiffContext::* ReinitFuncType) (Real)
 
typedef std::map< std::string, std::pair< unsigned int, unsigned int > > Counts
 

Protected Member Functions

virtual bool _general_residual (bool request_jacobian, DiffContext &, ResFuncType mass, ResFuncType damping, ResFuncType time_deriv, ResFuncType constraint, ReinitFuncType reinit)
 
void increment_constructor_count (const std::string &name)
 
void increment_destructor_count (const std::string &name)
 

Protected Attributes

Real _beta
 
Real _gamma
 
bool _is_accel_solve
 
bool _initial_accel_set
 
std::unique_ptr< NumericVector< Number > > _old_local_solution_rate
 
std::unique_ptr< NumericVector< Number > > _old_local_solution_accel
 
bool first_solve
 
bool first_adjoint_step
 
std::unique_ptr< DiffSolver_diff_solver
 
std::unique_ptr< LinearSolver< Number > > _linear_solver
 
sys_type_system
 
std::unique_ptr< SolutionHistorysolution_history
 

Static Protected Attributes

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

Detailed Description

This class defines a Newmark time integrator for second order (in time) DifferentiableSystems. There are two parameters $\gamma$ and $\beta$ (defaulting to 0.5 and 0.25, respectively).

Note
Projectors are included for the initial velocity and acceleration; the initial displacement can be set by calling System::project_solution().

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
Paul T. Bauman
Date
2015

Definition at line 46 of file newmark_solver.h.

Member Typedef Documentation

◆ Counts

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

Data structure to log the information. The log is identified by the class name.

Definition at line 117 of file reference_counter.h.

◆ Parent

The parent class

Definition at line 52 of file newmark_solver.h.

◆ ReinitFuncType

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

Definition at line 273 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 271 of file time_solver.h.

◆ sys_type

The type of system

Definition at line 65 of file time_solver.h.

Constructor & Destructor Documentation

◆ NewmarkSolver()

libMesh::NewmarkSolver::NewmarkSolver ( sys_type s)
explicit

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

Definition at line 25 of file newmark_solver.C.

◆ ~NewmarkSolver()

libMesh::NewmarkSolver::~NewmarkSolver ( )
virtual

Destructor.

Definition at line 33 of file newmark_solver.C.

34 {}

Member Function Documentation

◆ _general_residual()

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

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

Definition at line 205 of file newmark_solver.C.

References _beta, _gamma, _is_accel_solve, libMesh::TimeSolver::_system, libMesh::DenseVector< T >::add(), libMesh::DifferentiableSystem::deltat, libMesh::DiffContext::elem_solution_accel_derivative, libMesh::DiffContext::elem_solution_derivative, libMesh::DiffContext::elem_solution_rate_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_accel(), libMesh::DiffContext::get_elem_solution_rate(), libMesh::DifferentiableSystem::get_physics(), libMesh::UnsteadySolver::old_nonlinear_solution(), libMesh::SecondOrderUnsteadySolver::old_solution_accel(), libMesh::SecondOrderUnsteadySolver::old_solution_rate(), libMesh::Real, libMesh::DenseVector< T >::size(), libMesh::DenseMatrix< T >::swap(), and libMesh::System::use_fixed_solution.

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

212 {
213  unsigned int n_dofs = context.get_elem_solution().size();
214 
215  // We might need to save the old jacobian in case one of our physics
216  // terms later is unable to update it analytically.
217  DenseMatrix<Number> old_elem_jacobian(n_dofs, n_dofs);
218 
219  // Local velocity at old time step
220  DenseVector<Number> old_elem_solution_rate(n_dofs);
221  for (unsigned int i=0; i != n_dofs; ++i)
222  old_elem_solution_rate(i) =
223  old_solution_rate(context.get_dof_indices()[i]);
224 
225  // The user is computing the initial acceleration
226  // So upstream we've swapped _system.solution and _old_local_solution_accel
227  // So we need to give the context the correct entries since we're solving for
228  // acceleration here.
229  if (_is_accel_solve)
230  {
231  // System._solution is actually the acceleration right now so we need
232  // to reset the elem_solution to the right thing, which in this case
233  // is the initial guess for displacement, which got swapped into
234  // _old_solution_accel vector
235  DenseVector<Number> old_elem_solution(n_dofs);
236  for (unsigned int i=0; i != n_dofs; ++i)
237  old_elem_solution(i) =
238  old_solution_accel(context.get_dof_indices()[i]);
239 
240  context.elem_solution_derivative = 0.0;
241  context.elem_solution_rate_derivative = 0.0;
242  context.elem_solution_accel_derivative = 1.0;
243 
244  // Acceleration is currently the unknown so it's already sitting
245  // in elem_solution() thanks to FEMContext::pre_fe_reinit
246  context.get_elem_solution_accel() = context.get_elem_solution();
247 
248  // Now reset elem_solution() to what the user is expecting
249  context.get_elem_solution() = old_elem_solution;
250 
251  context.get_elem_solution_rate() = old_elem_solution_rate;
252 
253  // The user's Jacobians will be targeting derivatives w.r.t. u_{n+1}.
254  // Although the vast majority of cases will have the correct analytic
255  // Jacobians in this iteration, since we reset elem_solution_derivative*,
256  // if there are coupled/overlapping problems, there could be
257  // mismatches in the Jacobian. So we force finite differencing for
258  // the first iteration.
259  request_jacobian = false;
260  }
261  // Otherwise, the unknowns are the displacements and everything is straight
262  // forward and is what you think it is
263  else
264  {
265  if (request_jacobian)
266  old_elem_jacobian.swap(context.get_elem_jacobian());
267 
268  // Local displacement at old timestep
269  DenseVector<Number> old_elem_solution(n_dofs);
270  for (unsigned int i=0; i != n_dofs; ++i)
271  old_elem_solution(i) =
272  old_nonlinear_solution(context.get_dof_indices()[i]);
273 
274  // Local acceleration at old time step
275  DenseVector<Number> old_elem_solution_accel(n_dofs);
276  for (unsigned int i=0; i != n_dofs; ++i)
277  old_elem_solution_accel(i) =
278  old_solution_accel(context.get_dof_indices()[i]);
279 
280  // Convenience
282 
283  context.elem_solution_derivative = 1.0;
284 
285  // Local velocity at current time step
286  // v_{n+1} = gamma/(beta*Delta t)*(x_{n+1}-x_n)
287  // + (1-(gamma/beta))*v_n
288  // + (1-gamma/(2*beta))*(Delta t)*a_n
289  context.elem_solution_rate_derivative = (_gamma/(_beta*dt));
290 
291  context.get_elem_solution_rate() = context.get_elem_solution();
292  context.get_elem_solution_rate() -= old_elem_solution;
293  context.get_elem_solution_rate() *= context.elem_solution_rate_derivative;
294  context.get_elem_solution_rate().add( (1.0-_gamma/_beta), old_elem_solution_rate);
295  context.get_elem_solution_rate().add( (1.0-_gamma/(2.0*_beta))*dt, old_elem_solution_accel);
296 
297 
298 
299  // Local acceleration at current time step
300  // a_{n+1} = (1/(beta*(Delta t)^2))*(x_{n+1}-x_n)
301  // - 1/(beta*Delta t)*v_n
302  // - (1/(2*beta)-1)*a_n
303  context.elem_solution_accel_derivative = 1.0/(_beta*dt*dt);
304 
305  context.get_elem_solution_accel() = context.get_elem_solution();
306  context.get_elem_solution_accel() -= old_elem_solution;
307  context.get_elem_solution_accel() *= context.elem_solution_accel_derivative;
308  context.get_elem_solution_accel().add(-1.0/(_beta*dt), old_elem_solution_rate);
309  context.get_elem_solution_accel().add(-(1.0/(2.0*_beta)-1.0), old_elem_solution_accel);
310 
311  // Move the mesh into place first if necessary, set t = t_{n+1}
312  (context.*reinit_func)(1.);
313  }
314 
315  // If a fixed solution is requested, we'll use x_{n+1}
317  context.get_elem_fixed_solution() = context.get_elem_solution();
318 
319  // Get the time derivative at t_{n+1}, F(u_{n+1})
320  bool jacobian_computed = (_system.get_physics()->*time_deriv)(request_jacobian, context);
321 
322  // Damping at t_{n+1}, C(u_{n+1})
323  jacobian_computed = (_system.get_physics()->*damping)(jacobian_computed, context) &&
324  jacobian_computed;
325 
326  // Mass at t_{n+1}, M(u_{n+1})
327  jacobian_computed = (_system.get_physics()->*mass)(jacobian_computed, context) &&
328  jacobian_computed;
329 
330  // Add the constraint term
331  jacobian_computed = (_system.get_physics()->*constraint)(jacobian_computed, context) &&
332  jacobian_computed;
333 
334  // Add back (or restore) the old jacobian
335  if (request_jacobian)
336  {
337  if (jacobian_computed)
338  context.get_elem_jacobian() += old_elem_jacobian;
339  else
340  context.get_elem_jacobian().swap(old_elem_jacobian);
341  }
342 
343  return jacobian_computed;
344 }
Number old_solution_accel(const dof_id_type global_dof_number) const
Number old_nonlinear_solution(const dof_id_type global_dof_number) const
sys_type & _system
Definition: time_solver.h:258
bool use_fixed_solution
Definition: system.h:1493
Number old_solution_rate(const dof_id_type global_dof_number) const
const DifferentiablePhysics * get_physics() const
Definition: diff_system.h:182
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ adjoint_advance_timestep()

void libMesh::NewmarkSolver::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::UnsteadySolver.

Definition at line 101 of file newmark_solver.C.

102 {
103  libmesh_not_implemented();
104 }

◆ advance_timestep()

void libMesh::NewmarkSolver::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::UnsteadySolver.

Definition at line 43 of file newmark_solver.C.

References _beta, _gamma, libMesh::SecondOrderUnsteadySolver::_old_local_solution_accel, libMesh::SecondOrderUnsteadySolver::_old_local_solution_rate, libMesh::TimeSolver::_system, libMesh::NumericVector< T >::add(), libMesh::UnsteadySolver::advance_timestep(), libMesh::NumericVector< T >::clone(), libMesh::DifferentiableSystem::deltat, libMesh::UnsteadySolver::first_solve, libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::System::get_vector(), libMesh::SecondOrderUnsteadySolver::old_solution_accel(), libMesh::SecondOrderUnsteadySolver::old_solution_rate(), and libMesh::System::solution.

44 {
45  // We need to update velocity and acceleration before
46  // we update the nonlinear solution (displacement) and
47  // delta_t
48 
49  NumericVector<Number> & old_solution_rate =
50  _system.get_vector("_old_solution_rate");
51 
52  NumericVector<Number> & old_solution_accel =
53  _system.get_vector("_old_solution_accel");
54 
55  if (!first_solve)
56  {
57  NumericVector<Number> & old_nonlinear_soln =
58  _system.get_vector("_old_nonlinear_solution");
59 
60  NumericVector<Number> & nonlinear_solution =
61  *(_system.solution);
62 
63  // We need to cache the new solution_rate before updating the old_solution_rate
64  // so we can update acceleration with the proper old_solution_rate
65  // v_{n+1} = gamma/(beta*Delta t)*(x_{n+1}-x_n)
66  // - ((gamma/beta)-1)*v_n
67  // - (gamma/(2*beta)-1)*(Delta t)*a_n
68  std::unique_ptr<NumericVector<Number>> new_solution_rate = nonlinear_solution.clone();
69  (*new_solution_rate) -= old_nonlinear_soln;
70  (*new_solution_rate) *= (_gamma/(_beta*_system.deltat));
71  new_solution_rate->add( (1.0-_gamma/_beta), old_solution_rate );
72  new_solution_rate->add( (1.0-_gamma/(2.0*_beta))*_system.deltat, old_solution_accel );
73 
74  // a_{n+1} = (1/(beta*(Delta t)^2))*(x_{n+1}-x_n)
75  // - 1/(beta*Delta t)*v_n
76  // - (1-1/(2*beta))*a_n
77  std::unique_ptr<NumericVector<Number>> new_solution_accel = old_solution_accel.clone();
78  (*new_solution_accel) *= -(1.0/(2.0*_beta)-1.0);
79  new_solution_accel->add( -1.0/(_beta*_system.deltat), old_solution_rate );
80  new_solution_accel->add( 1.0/(_beta*_system.deltat*_system.deltat), nonlinear_solution );
81  new_solution_accel->add( -1.0/(_beta*_system.deltat*_system.deltat), old_nonlinear_soln );
82 
83  // Now update old_solution_rate
84  old_solution_rate = (*new_solution_rate);
85  old_solution_accel = (*new_solution_accel);
86  }
87 
88  // Localize updated vectors
89  old_solution_rate.localize
92 
93  old_solution_accel.localize
96 
97  // Now we can finish advancing the timestep
99 }
Number old_solution_accel(const dof_id_type global_dof_number) const
virtual std::unique_ptr< NumericVector< T > > clone() const =0
std::unique_ptr< NumericVector< Number > > _old_local_solution_accel
sys_type & _system
Definition: time_solver.h:258
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:774
std::unique_ptr< NumericVector< Number > > solution
Definition: system.h:1523
Number old_solution_rate(const dof_id_type global_dof_number) const
virtual void advance_timestep() override
std::unique_ptr< NumericVector< Number > > _old_local_solution_rate
const DofMap & get_dof_map() const
Definition: system.h:2049
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:450

◆ 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 167 of file time_solver.h.

167 {}

◆ compute_initial_accel()

void libMesh::NewmarkSolver::compute_initial_accel ( )
virtual

This method uses the specified initial displacement and velocity to compute the initial acceleration $a_0$.

Definition at line 106 of file newmark_solver.C.

References libMesh::TimeSolver::_diff_solver, _is_accel_solve, libMesh::TimeSolver::_system, libMesh::System::get_vector(), set_initial_accel_avail(), libMesh::System::solution, and libMesh::System::update().

107 {
108  // We need to compute the initial acceleration based off of
109  // the initial position and velocity and, thus, acceleration
110  // is the unknown in diff_solver and not the displacement. So,
111  // We swap solution and acceleration. NewmarkSolver::_general_residual
112  // will check _is_accel_solve and provide the correct
113  // values to the FEMContext assuming this swap was made.
114  this->_is_accel_solve = true;
115 
116  //solution->accel, accel->solution
117  _system.solution->swap(_system.get_vector("_old_solution_accel"));
118  _system.update();
119 
120  this->_diff_solver->solve();
121 
122  // solution->solution, accel->accel
123  _system.solution->swap(_system.get_vector("_old_solution_accel"));
124  _system.update();
125 
126  // We're done, so no longer doing an acceleration solve
127  this->_is_accel_solve = false;
128 
129  this->set_initial_accel_avail(true);
130 }
std::unique_ptr< DiffSolver > _diff_solver
Definition: time_solver.h:248
sys_type & _system
Definition: time_solver.h:258
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:774
std::unique_ptr< NumericVector< Number > > solution
Definition: system.h:1523
void set_initial_accel_avail(bool initial_accel_set)
virtual void update()
Definition: system.C:408

◆ 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 182 of file time_solver.h.

References libMesh::TimeSolver::_diff_solver.

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

182 { return _diff_solver; }
std::unique_ptr< DiffSolver > _diff_solver
Definition: time_solver.h:248

◆ disable_print_counter_info()

void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited

Definition at line 106 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

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

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

◆ 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 227 of file unsteady_solver.C.

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

228 {
229 
230  std::unique_ptr<NumericVector<Number>> solution_copy =
231  _system.solution->clone();
232 
233  solution_copy->add(-1., _system.get_vector("_old_nonlinear_solution"));
234 
235  solution_copy->close();
236 
237  return _system.calculate_norm(*solution_copy, norm);
238 }
sys_type & _system
Definition: time_solver.h:258
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:774
std::unique_ptr< NumericVector< Number > > solution
Definition: system.h:1523
Real calculate_norm(const NumericVector< Number > &v, unsigned int var, FEMNormType norm_type, std::set< unsigned int > *skip_dimensions=nullptr) const
Definition: system.C:1378

◆ element_residual()

bool libMesh::NewmarkSolver::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 163 of file newmark_solver.C.

References libMesh::DifferentiablePhysics::_eulerian_time_deriv(), _general_residual(), libMesh::DifferentiablePhysics::damping_residual(), libMesh::DiffContext::elem_reinit(), libMesh::DifferentiablePhysics::element_constraint(), and libMesh::DifferentiablePhysics::mass_residual().

165 {
166  return this->_general_residual(request_jacobian,
167  context,
173 }
virtual bool damping_residual(bool request_jacobian, DiffContext &)
Definition: diff_physics.h:374
virtual void elem_reinit(Real)
Definition: diff_context.h:76
bool _eulerian_time_deriv(bool request_jacobian, DiffContext &)
virtual bool mass_residual(bool request_jacobian, DiffContext &)
Definition: diff_physics.h:318
virtual bool _general_residual(bool request_jacobian, DiffContext &, ResFuncType mass, ResFuncType damping, ResFuncType time_deriv, ResFuncType constraint, ReinitFuncType reinit)
virtual bool element_constraint(bool request_jacobian, DiffContext &)
Definition: diff_physics.h:142

◆ enable_print_counter_info()

void libMesh::ReferenceCounter::enable_print_counter_info ( )
staticinherited

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

Definition at line 100 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

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

◆ error_order()

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

Error convergence order: 2 for $\gamma=0.5$, 1 otherwise

Implements libMesh::UnsteadySolver.

Definition at line 36 of file newmark_solver.C.

References _gamma.

37 {
38  if (_gamma == 0.5)
39  return 2.;
40  return 1.;
41 }

◆ get_info()

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

Gets a string containing the reference information.

Definition at line 47 of file reference_counter.C.

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

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

48 {
49 #if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
50 
51  std::ostringstream oss;
52 
53  oss << '\n'
54  << " ---------------------------------------------------------------------------- \n"
55  << "| Reference count information |\n"
56  << " ---------------------------------------------------------------------------- \n";
57 
58  for (const auto & pr : _counts)
59  {
60  const std::string name(pr.first);
61  const unsigned int creations = pr.second.first;
62  const unsigned int destructions = pr.second.second;
63 
64  oss << "| " << name << " reference count information:\n"
65  << "| Creations: " << creations << '\n'
66  << "| Destructions: " << destructions << '\n';
67  }
68 
69  oss << " ---------------------------------------------------------------------------- \n";
70 
71  return oss.str();
72 
73 #else
74 
75  return "";
76 
77 #endif
78 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:42

◆ increment_constructor_count()

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

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

Definition at line 181 of file reference_counter.h.

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

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

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

◆ increment_destructor_count()

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

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

Definition at line 194 of file reference_counter.h.

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

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

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

◆ init()

void libMesh::SecondOrderUnsteadySolver::init ( )
overridevirtualinherited

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

Reimplemented from libMesh::UnsteadySolver.

Definition at line 34 of file second_order_unsteady_solver.C.

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

35 {
37 
38  _system.add_vector("_old_solution_rate");
39  _system.add_vector("_old_solution_accel");
40 }
sys_type & _system
Definition: time_solver.h:258
NumericVector< Number > & add_vector(const std::string &vec_name, const bool projections=true, const ParallelType type=PARALLEL)
Definition: system.C:661
virtual void init() override

◆ init_data()

void libMesh::SecondOrderUnsteadySolver::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::UnsteadySolver.

Definition at line 42 of file second_order_unsteady_solver.C.

References libMesh::SecondOrderUnsteadySolver::_old_local_solution_accel, libMesh::SecondOrderUnsteadySolver::_old_local_solution_rate, libMesh::TimeSolver::_system, libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::GHOSTED, libMesh::UnsteadySolver::init_data(), libMesh::System::n_dofs(), libMesh::System::n_local_dofs(), and libMesh::SERIAL.

43 {
45 
46 #ifdef LIBMESH_ENABLE_GHOSTED
49  GHOSTED);
50 
53  GHOSTED);
54 #else
57 #endif
58 }
virtual void init_data() override
std::unique_ptr< NumericVector< Number > > _old_local_solution_accel
dof_id_type n_local_dofs() const
Definition: system.C:187
dof_id_type n_dofs() const
Definition: system.C:150
sys_type & _system
Definition: time_solver.h:258
std::unique_ptr< NumericVector< Number > > _old_local_solution_rate
const DofMap & get_dof_map() const
Definition: system.h:2049
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:450

◆ 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 233 of file time_solver.h.

References libMesh::TimeSolver::_is_adjoint.

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

234  { return _is_adjoint; }

◆ is_steady()

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

This is not a steady-state solver.

Implements libMesh::TimeSolver.

Definition at line 154 of file unsteady_solver.h.

154 { return false; }

◆ 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 187 of file time_solver.h.

References libMesh::TimeSolver::_linear_solver.

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

187 { return _linear_solver; }
std::unique_ptr< LinearSolver< Number > > _linear_solver
Definition: time_solver.h:253

◆ n_objects()

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

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

Definition at line 83 of file reference_counter.h.

References libMesh::ReferenceCounter::_n_objects.

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

◆ nonlocal_residual()

bool libMesh::NewmarkSolver::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 191 of file newmark_solver.C.

References _general_residual(), libMesh::DifferentiablePhysics::nonlocal_constraint(), libMesh::DifferentiablePhysics::nonlocal_damping_residual(), libMesh::DifferentiablePhysics::nonlocal_mass_residual(), libMesh::DiffContext::nonlocal_reinit(), and libMesh::DifferentiablePhysics::nonlocal_time_derivative().

193 {
194  return this->_general_residual(request_jacobian,
195  context,
201 }
virtual bool nonlocal_time_derivative(bool request_jacobian, DiffContext &)
Definition: diff_physics.h:209
virtual bool nonlocal_constraint(bool request_jacobian, DiffContext &)
Definition: diff_physics.h:227
virtual bool _general_residual(bool request_jacobian, DiffContext &, ResFuncType mass, ResFuncType damping, ResFuncType time_deriv, ResFuncType constraint, ReinitFuncType reinit)
virtual bool nonlocal_damping_residual(bool request_jacobian, DiffContext &)
Definition: diff_physics.h:406
virtual bool nonlocal_mass_residual(bool request_jacobian, DiffContext &c)
virtual void nonlocal_reinit(Real)
Definition: diff_context.h:94

◆ 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 216 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 libMesh::EulerSolver::_general_residual(), libMesh::Euler2Solver::_general_residual(), and _general_residual().

218 {
219  libmesh_assert_less (global_dof_number, _system.get_dof_map().n_dofs());
220  libmesh_assert_less (global_dof_number, old_local_nonlinear_solution->size());
221 
222  return (*old_local_nonlinear_solution)(global_dof_number);
223 }
std::unique_ptr< NumericVector< Number > > old_local_nonlinear_solution
sys_type & _system
Definition: time_solver.h:258
dof_id_type n_dofs() const
Definition: dof_map.h:574
const DofMap & get_dof_map() const
Definition: system.h:2049

◆ old_solution_accel()

Number libMesh::SecondOrderUnsteadySolver::old_solution_accel ( const dof_id_type  global_dof_number) const
inherited
Returns
The solution acceleration at the previous time step, $\ddot{u}_n$, for the specified global DOF.

Definition at line 116 of file second_order_unsteady_solver.C.

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

Referenced by _general_residual(), advance_timestep(), project_initial_accel(), and libMesh::SecondOrderUnsteadySolver::reinit().

118 {
119  libmesh_assert_less (global_dof_number, _system.get_dof_map().n_dofs());
120  libmesh_assert_less (global_dof_number, _old_local_solution_accel->size());
121 
122  return (*_old_local_solution_accel)(global_dof_number);
123 }
std::unique_ptr< NumericVector< Number > > _old_local_solution_accel
sys_type & _system
Definition: time_solver.h:258
dof_id_type n_dofs() const
Definition: dof_map.h:574
const DofMap & get_dof_map() const
Definition: system.h:2049

◆ old_solution_rate()

Number libMesh::SecondOrderUnsteadySolver::old_solution_rate ( const dof_id_type  global_dof_number) const
inherited
Returns
The solution rate at the previous time step, $\dot{u}_n$, for the specified global DOF.

Definition at line 107 of file second_order_unsteady_solver.C.

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

Referenced by _general_residual(), advance_timestep(), libMesh::SecondOrderUnsteadySolver::project_initial_rate(), and libMesh::SecondOrderUnsteadySolver::reinit().

109 {
110  libmesh_assert_less (global_dof_number, _system.get_dof_map().n_dofs());
111  libmesh_assert_less (global_dof_number, _old_local_solution_rate->size());
112 
113  return (*_old_local_solution_rate)(global_dof_number);
114 }
sys_type & _system
Definition: time_solver.h:258
dof_id_type n_dofs() const
Definition: dof_map.h:574
std::unique_ptr< NumericVector< Number > > _old_local_solution_rate
const DofMap & get_dof_map() const
Definition: system.h:2049

◆ print_info()

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

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

Definition at line 87 of file reference_counter.C.

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

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

◆ project_initial_accel()

void libMesh::NewmarkSolver::project_initial_accel ( FunctionBase< Number > *  f,
FunctionBase< Gradient > *  g = nullptr 
)

Specify non-zero initial acceleration. Should be called before solve(). This is an alternative to compute_initial_acceleration() if the initial acceleration is actually known. The function value f and its gradient g are user-provided cloneable functors. A gradient g is only required/used for projecting onto finite element spaces with continuous derivatives.

Definition at line 132 of file newmark_solver.C.

References libMesh::TimeSolver::_system, libMesh::System::get_vector(), libMesh::SecondOrderUnsteadySolver::old_solution_accel(), libMesh::System::project_vector(), and set_initial_accel_avail().

134 {
135  NumericVector<Number> & old_solution_accel =
136  _system.get_vector("_old_solution_accel");
137 
139 
140  this->set_initial_accel_avail(true);
141 }
Number old_solution_accel(const dof_id_type global_dof_number) const
sys_type & _system
Definition: time_solver.h:258
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:774
void set_initial_accel_avail(bool initial_accel_set)
void project_vector(NumericVector< Number > &new_vector, FunctionBase< Number > *f, FunctionBase< Gradient > *g=nullptr, int is_adjoint=-1) const

◆ project_initial_rate()

void libMesh::SecondOrderUnsteadySolver::project_initial_rate ( FunctionBase< Number > *  f,
FunctionBase< Gradient > *  g = nullptr 
)
inherited

Specify non-zero initial velocity. Should be called before solve(). The function value f and its gradient g are user-provided cloneable functors. A gradient g is only required/used for projecting onto finite element spaces with continuous derivatives.

Definition at line 98 of file second_order_unsteady_solver.C.

References libMesh::TimeSolver::_system, libMesh::System::get_vector(), libMesh::SecondOrderUnsteadySolver::old_solution_rate(), and libMesh::System::project_vector().

100 {
101  NumericVector<Number> & old_solution_rate =
102  _system.get_vector("_old_solution_rate");
103 
105 }
sys_type & _system
Definition: time_solver.h:258
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:774
Number old_solution_rate(const dof_id_type global_dof_number) const
void project_vector(NumericVector< Number > &new_vector, FunctionBase< Number > *f, FunctionBase< Gradient > *g=nullptr, int is_adjoint=-1) const

◆ reinit()

void libMesh::SecondOrderUnsteadySolver::reinit ( )
overridevirtualinherited

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

Reimplemented from libMesh::UnsteadySolver.

Definition at line 60 of file second_order_unsteady_solver.C.

References libMesh::SecondOrderUnsteadySolver::_old_local_solution_accel, libMesh::SecondOrderUnsteadySolver::_old_local_solution_rate, libMesh::TimeSolver::_system, libMesh::System::get_dof_map(), libMesh::DofMap::get_send_list(), libMesh::System::get_vector(), libMesh::GHOSTED, libMesh::System::n_dofs(), libMesh::System::n_local_dofs(), libMesh::SecondOrderUnsteadySolver::old_solution_accel(), libMesh::SecondOrderUnsteadySolver::old_solution_rate(), libMesh::UnsteadySolver::reinit(), and libMesh::SERIAL.

61 {
63 
64 #ifdef LIBMESH_ENABLE_GHOSTED
67  GHOSTED);
68 
71  GHOSTED);
72 #else
75 #endif
76 
77  // localize the old solutions
78  NumericVector<Number> & old_solution_rate =
79  _system.get_vector("_old_solution_rate");
80 
81  NumericVector<Number> & old_solution_accel =
82  _system.get_vector("_old_solution_accel");
83 
84  old_solution_rate.localize
87 
88  old_solution_accel.localize
91 }
virtual void reinit() override
Number old_solution_accel(const dof_id_type global_dof_number) const
std::unique_ptr< NumericVector< Number > > _old_local_solution_accel
dof_id_type n_local_dofs() const
Definition: system.C:187
dof_id_type n_dofs() const
Definition: system.C:150
sys_type & _system
Definition: time_solver.h:258
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:774
Number old_solution_rate(const dof_id_type global_dof_number) const
std::unique_ptr< NumericVector< Number > > _old_local_solution_rate
const DofMap & get_dof_map() const
Definition: system.h:2049
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:450

◆ retrieve_timestep()

void libMesh::SecondOrderUnsteadySolver::retrieve_timestep ( )
overridevirtualinherited

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

Reimplemented from libMesh::UnsteadySolver.

Definition at line 93 of file second_order_unsteady_solver.C.

94 {
95  libmesh_not_implemented();
96 }

◆ set_beta()

void libMesh::NewmarkSolver::set_beta ( Real  beta)
inline

Setter for $ \beta $

Definition at line 149 of file newmark_solver.h.

References _beta.

150  { _beta = beta; }

◆ set_gamma()

void libMesh::NewmarkSolver::set_gamma ( Real  gamma)
inline

Setter for $ \gamma $.

Note
The method is second order only for $ \gamma = 1/2 $.

Definition at line 157 of file newmark_solver.h.

References _gamma.

158  { _gamma = gamma; }

◆ set_initial_accel_avail()

void libMesh::NewmarkSolver::set_initial_accel_avail ( bool  initial_accel_set)

Allow the user to (re)set whether the initial acceleration is available. This is not needed if either compute_initial_accel() or project_initial_accel() are called. This is useful is the user is restarting a calculation and the acceleration is available from the restart.

Definition at line 143 of file newmark_solver.C.

References _initial_accel_set.

Referenced by compute_initial_accel(), and project_initial_accel().

144 {
145  _initial_accel_set = initial_accel_set;
146 }

◆ 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 240 of file time_solver.h.

References libMesh::TimeSolver::_is_adjoint.

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

241  { _is_adjoint = _is_adjoint_value; }

◆ 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 97 of file time_solver.C.

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

98 {
99  solution_history = _solution_history.clone();
100 }
std::unique_ptr< SolutionHistory > solution_history
Definition: time_solver.h:265

◆ side_residual()

bool libMesh::NewmarkSolver::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 177 of file newmark_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().

179 {
180  return this->_general_residual(request_jacobian,
181  context,
187 }
virtual bool side_constraint(bool request_jacobian, DiffContext &)
Definition: diff_physics.h:191
virtual bool side_damping_residual(bool request_jacobian, DiffContext &)
Definition: diff_physics.h:391
virtual bool _general_residual(bool request_jacobian, DiffContext &, ResFuncType mass, ResFuncType damping, ResFuncType time_deriv, ResFuncType constraint, ReinitFuncType reinit)
virtual void elem_side_reinit(Real)
Definition: diff_context.h:82
virtual bool side_mass_residual(bool request_jacobian, DiffContext &)
Definition: diff_physics.h:335
virtual bool side_time_derivative(bool request_jacobian, DiffContext &)
Definition: diff_physics.h:171

◆ solve()

void libMesh::NewmarkSolver::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::UnsteadySolver.

Definition at line 148 of file newmark_solver.C.

References _initial_accel_set, and libMesh::UnsteadySolver::solve().

149 {
150  // First, check that the initial accel was set one way or another
151  if (!_initial_accel_set)
152  {
153  std::string error = "ERROR: Must first set initial acceleration using one of:\n";
154  error += "NewmarkSolver::compute_initial_accel()\n";
155  error += "NewmarkSolver::project_initial_accel()\n";
156  libmesh_error_msg(error);
157  }
158 
159  // That satisfied, now we can solve
161 }
virtual void solve() override

◆ system() [1/2]

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

Definition at line 172 of file time_solver.h.

References libMesh::TimeSolver::_system.

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

172 { return _system; }
sys_type & _system
Definition: time_solver.h:258

◆ system() [2/2]

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

Definition at line 177 of file time_solver.h.

References libMesh::TimeSolver::_system.

177 { return _system; }
sys_type & _system
Definition: time_solver.h:258

◆ time_order()

virtual unsigned int libMesh::SecondOrderUnsteadySolver::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 53 of file second_order_unsteady_solver.h.

54  { return 2; }

Member Data Documentation

◆ _beta

Real libMesh::NewmarkSolver::_beta
protected

The value for the $\beta$ to employ. Method is unconditionally stable for $ \beta \ge \frac{1}{4} \left( \gamma + \frac{1}{2}\right)^2 $

Definition at line 167 of file newmark_solver.h.

Referenced by _general_residual(), advance_timestep(), and set_beta().

◆ _counts

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

◆ _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 248 of file time_solver.h.

Referenced by 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 141 of file reference_counter.h.

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

◆ _gamma

Real libMesh::NewmarkSolver::_gamma
protected

The value for $\gamma$ to employ. Newmark is 2nd order iff $\gamma=0.5$.

Definition at line 173 of file newmark_solver.h.

Referenced by _general_residual(), advance_timestep(), error_order(), and set_gamma().

◆ _initial_accel_set

bool libMesh::NewmarkSolver::_initial_accel_set
protected

This method requires an initial acceleration. So, we force the user to call either compute_initial_accel or project_initial_accel to set the initial acceleration.

Definition at line 187 of file newmark_solver.h.

Referenced by set_initial_accel_avail(), and solve().

◆ _is_accel_solve

bool libMesh::NewmarkSolver::_is_accel_solve
protected

Need to be able to indicate to _general_residual if we are doing an acceleration solve or not.

Definition at line 179 of file newmark_solver.h.

Referenced by _general_residual(), and compute_initial_accel().

◆ _linear_solver

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

An implicit linear solver to use for adjoint problems.

Definition at line 253 of file time_solver.h.

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

◆ _mutex

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

Mutual exclusion object to enable thread-safe reference counting.

Definition at line 135 of file reference_counter.h.

◆ _n_objects

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

The number of objects. Print the reference count information when the number returns to 0.

Definition at line 130 of file reference_counter.h.

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

◆ _old_local_solution_accel

std::unique_ptr<NumericVector<Number> > libMesh::SecondOrderUnsteadySolver::_old_local_solution_accel
protectedinherited

◆ _old_local_solution_rate

std::unique_ptr<NumericVector<Number> > libMesh::SecondOrderUnsteadySolver::_old_local_solution_rate
protectedinherited

◆ _system

sys_type& libMesh::TimeSolver::_system
protectedinherited

A reference to the system we are solving.

Definition at line 258 of file time_solver.h.

Referenced by libMesh::EulerSolver::_general_residual(), libMesh::Euler2Solver::_general_residual(), libMesh::SteadySolver::_general_residual(), _general_residual(), libMesh::UnsteadySolver::adjoint_advance_timestep(), advance_timestep(), libMesh::AdaptiveTimeSolver::advance_timestep(), libMesh::UnsteadySolver::advance_timestep(), compute_initial_accel(), libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), libMesh::UnsteadySolver::du(), libMesh::EulerSolver::element_residual(), libMesh::Euler2Solver::element_residual(), libMesh::EigenTimeSolver::element_residual(), libMesh::SecondOrderUnsteadySolver::init(), libMesh::UnsteadySolver::init(), libMesh::TimeSolver::init(), libMesh::EigenTimeSolver::init(), libMesh::SecondOrderUnsteadySolver::init_data(), libMesh::UnsteadySolver::init_data(), libMesh::TimeSolver::init_data(), libMesh::EulerSolver::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(), 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(), and libMesh::TimeSolver::system().

◆ 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 168 of file unsteady_solver.h.

Referenced by libMesh::UnsteadySolver::adjoint_advance_timestep().

◆ first_solve

bool libMesh::UnsteadySolver::first_solve
protectedinherited

◆ old_local_nonlinear_solution

◆ quiet

bool libMesh::TimeSolver::quiet
inherited

Print extra debugging information if quiet == false.

Definition at line 192 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 221 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

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 265 of file time_solver.h.

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


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