Drake
IntegratorBase< T > Class Template Referenceabstract

An abstract class for an integrator for ODEs and DAEs as represented by a Drake System. More...

#include <drake/systems/analysis/integrator_base.h>

Inheritance diagram for IntegratorBase< T >:
[legend]

Public Types

enum  StepResult {
  kReachedPublishTime = 1, kReachedZeroCrossing = 2, kReachedUpdateTime = 3, kTimeHasAdvanced = 4,
  kReachedBoundaryTime = 5, kReachedStepLimit = 6
}
 Status returned by StepOnceAtMost(). More...
 

Public Member Functions

 IntegratorBase (const System< T > &system, Context< T > *context=nullptr)
 Maintains references to the system being integrated and the context used to specify the initial conditions for that system (if any). More...
 
virtual ~IntegratorBase ()=default
 Destructor. More...
 
virtual bool supports_error_estimation () const =0
 Indicates whether an integrator supports error estimation. More...
 
void set_fixed_step_mode (bool flag)
 Sets an integrator with error control to fixed step mode. More...
 
bool get_fixed_step_mode () const
 Gets whether an integrator is running in fixed step mode. More...
 
void set_target_accuracy (double accuracy)
 Request that the integrator attempt to achieve a particular accuracy for the continuous portions of the simulation. More...
 
double get_target_accuracy () const
 Gets the target accuracy. More...
 
double get_accuracy_in_use () const
 Gets the accuracy in use by the integrator. More...
 
void set_maximum_step_size (const T &max_step_size)
 Sets the maximum step size that may be taken by this integrator. More...
 
const T & get_maximum_step_size () const
 Gets the maximum step size that may be taken by this integrator. More...
 
void Reset ()
 Resets the integrator to initial values, i.e., default construction values. More...
 
void Initialize ()
 An integrator must be initialized before being used. More...
 
void request_initial_step_size_target (const T &step_size)
 Request that the first attempted integration step have a particular size. More...
 
const T & get_initial_step_size_target () const
 Gets the target size of the first integration step. More...
 
StepResult IntegrateAtMost (const T &publish_dt, const T &update_dt, const T &boundary_dt)
 Integrates the system forward in time by a single step with step size subject to integration error tolerances (assuming that the integrator supports error estimation). More...
 
double get_stretch_factor () const
 Gets the stretch factor (> 1), which is multiplied by the maximum (typically user-designated) integration step size to obtain the amount that the integrator is able to stretch the maximum time step toward hitting an upcoming publish or update event in IntegrateAtMost(). More...
 
void IntegrateWithMultipleSteps (const T &dt)
 Stepping function for integrators operating outside of Simulator that advances the continuous state exactly by dt. More...
 
void IntegrateWithSingleFixedStep (const T &dt)
 Stepping function for integrators operating outside of Simulator that advances the continuous state exactly by dt and using a single fixed step. More...
 
const T & get_ideal_next_step_size () const
 Return the step size the integrator would like to take next, based primarily on the integrator's accuracy prediction. More...
 
const Context< T > & get_context () const
 Returns a const reference to the internally-maintained Context holding the most recent state in the trajectory. More...
 
Context< T > * get_mutable_context ()
 Returns a mutable pointer to the internally-maintained Context holding the most recent state in the trajectory. More...
 
void reset_context (Context< T > *context)
 Replace the pointer to the internally-maintained Context with a different one. More...
 
const System< T > & get_system () const
 Gets a constant reference to the system that is being integrated (and was provided to the constructor of the integrator). More...
 
bool is_initialized () const
 Indicates whether the integrator has been initialized. More...
 
virtual int get_error_estimate_order () const =0
 Derived classes must override this function to return the order of the integrator's error estimate. More...
 
const T & get_previous_integration_step_size () const
 Gets the size of the last (previous) integration step. More...
 
const ContinuousState< T > * get_error_estimate () const
 Gets the error estimate (used only for integrators that support error estimation). More...
 
Does not allow copy, move, or assignment
 IntegratorBase (const IntegratorBase &)=delete
 
IntegratorBaseoperator= (const IntegratorBase &)=delete
 
 IntegratorBase (IntegratorBase &&)=delete
 
IntegratorBaseoperator= (IntegratorBase &&)=delete
 
Methods for minimum integration step size selection and behavior

Variable step integrators reduce their step sizes as needed to achieve requirements such as specified accuracy or step convergence. However, it is not possible to take an arbitrarily small step. Normally integrators choose an appropriate minimum step and throw an exception if the requirements can't be achieved without going below that. Methods in this section allow you to influence two aspects of this procedure:

  • you can increase the minimum step size, and
  • you can control whether an exception is thrown if a smaller step would have been needed to achieve the aforementioned integrator requirements.

By default, integrators allow a very small minimum step which can result in long run times. Setting a larger minimum can be helpful as a diagnostic to figure out what aspect of your simulation is requiring small steps. You can set the minimum to what should be a "reasonable" minimum based on what you know about the physical system. You will then get an std::runtime_error exception thrown at any point in time where your model behaves unexpectedly (due to, e.g., a discontinuity in the derivative evaluation function).

If you disable the exception (via set_throw_on_minimum_step_size_violation(false)), the integrator will simply proceed with a step of the minimum size: accuracy is guaranteed only when the minimum step size is not violated. Beware that there can be no guarantee about the magnitude of any errors introduced by violating the accuracy "requirements" in this manner, so disabling the exception should be done warily.

Details

Because time is maintained to finite precision, there is an absolute minimum step size h_floor required to avoid roundoff error. The integrator will never take a step smaller than h_floor. We calculate h_floor=max(ε,ε⋅t), where t is the current time and ε is a small multiple of machine precision, typically a number like 1e-14. Note that h_floor necessarily grows with time; if that is a concern you should limit how long your simulations are allowed to run without resetting time.

You may request a larger minimum step size h_min. Then at every time t, the integrator determines a "working" minimum h_work=max(h_min,h_floor). If the step size selection algorithm determines that a step smaller than h_work is needed to meet accuracy or other needs, then a std::runtime_error exception will be thrown and the simulation halted. On the other hand, if you have suppressed the exception (again, via set_throw_on_minimum_step_size_violation(false)), the integration will continue, taking a step of size h_work.

Under some circumstances the integrator may legitimately take a step of size h smaller than your specified h_min, although never smaller than h_floor. For example, occasionally the integrator may reach an event or time limit that occurs a very short time after the end of a previous step, necessitating that a tiny "sliver" of a step be taken to complete the interval. That does not indicate an error, and required accuracy and convergence goals are achieved. Larger steps can resume immediately afterwards. Another circumstance is when one of the integrator's stepping methods is called directly requesting a very small step, for example IntegrateWithMultipleSteps(h). No exception will be thrown in either of these cases.

void set_requested_minimum_step_size (const T &min_step_size)
 Sets the requested minimum step size h_min that may be taken by this integrator. More...
 
const T & get_requested_minimum_step_size () const
 Gets the requested minimum step size h_min for this integrator. More...
 
void set_throw_on_minimum_step_size_violation (bool throws)
 Sets whether the integrator should throw a std::runtime_error exception when the integrator's step size selection algorithm determines that it must take a step smaller than the minimum step size (for, e.g., purposes of error control). More...
 
bool get_throw_on_minimum_step_size_violation () const
 Reports the current setting of the throw_on_minimum_step_size_violation flag. More...
 
get_working_minimum_step_size () const
 Gets the current value of the working minimum step size h_work(t) for this integrator, which may vary with the current time t as stored in the integrator's context. More...
 
Integrator statistics methods.

These methods allow the caller to manipulate and query integrator statistics.

Generally speaking, the larger the integration step taken, the faster a simulation will run. These methods allow querying (and resetting) the integrator statistics as one means of determining how to make a simulation run faster.

void ResetStatistics ()
 Forget accumulated statistics. More...
 
int64_t get_num_substep_failures () const
 Gets the number of failed sub-steps (implying one or more step size reductions was required to permit solving the necessary nonlinear system of equations). More...
 
int64_t get_num_step_shrinkages_from_substep_failures () const
 Gets the number of step size shrinkages due to sub-step failures (e.g., integrator convergence failures) since the last call to ResetStatistics() or Initialize(). More...
 
int64_t get_num_step_shrinkages_from_error_control () const
 Gets the number of step size shrinkages due to failure to meet targeted error tolerances, since the last call to ResetStatistics or Initialize(). More...
 
int64_t get_num_derivative_evaluations () const
 Returns the number of ODE function evaluations (calls to CalcTimeDerivatives()) since the last call to ResetStatistics() or Initialize(). More...
 
const T & get_actual_initial_step_size_taken () const
 The actual size of the successful first step. More...
 
const T & get_smallest_adapted_step_size_taken () const
 The size of the smallest step taken as the result of a controlled integration step adjustment since the last Initialize() or ResetStatistics() call. More...
 
const T & get_largest_step_size_taken () const
 The size of the largest step taken since the last Initialize() or ResetStatistics() call. More...
 
int64_t get_num_steps_taken () const
 The number of integration steps taken since the last Initialize() or ResetStatistics() call. More...
 
Methods for weighting state variable errors

This group of methods describes how errors for state variables with heterogeneous units are weighted in the context of error-controlled integration.

This is an advanced topic and most users can simply specify desired accuracy and accept the default state variable weights.

A collection of state variables is generally defined in heterogenous units (e.g. length, angles, velocities, energy). Some of the state variables cannot even be expressed in meaningful units, like quaternions. Certain integrators provide an estimate of the absolute error made in each state variable during an integration step. These errors must be properly weighted to obtain an "accuracy" with respect to each particular variable. These per-variable accuracy determinations can be compared against the user's requirements and used to select an appropriate size for the next step [Sherman 2011]. The weights are normally determined automatically using the system's characteristic dimensions, so most users can stop reading now! Custom weighting is primarily useful for performance improvement; an optimal weighting would allow an error-controlled integrator to provide the desired level of accuracy across all state variables without wasting computation achieving superfluous accuracy for some of those variables.

Users interested in more precise control over state variable weighting may use the methods in this group to access and modify weighting factors for individual state variables. Changes to these weights can only be made prior to integrator initialization or as a result of an event being triggered and then followed by re-initialization.

Relative versus absolute accuracy

State variable integration error, as estimated by an integrator, is an absolute quantity with the same units as the variable. At each time step we therefore need to determine an absolute error that would be deemed "good enough", i.e. satisfies the user's accuracy requirement. If a variable is maintained to a relative accuracy then that "good enough" value is defined to be the required accuracy a (a fraction like 0.001) times the current value of the variable, as long as that value is far from zero. For variables maintained to an absolute accuracy, or relative variables that are at or near zero (where relative accuracy would be undefined or too strict, respectively), we need a different way to determine the "good enough" absolute error. The methods in this section control how that absolute error value is calculated.

How to choose weights

The weight wᵢ for a state variable xᵢ should be chosen so that the product wᵢ * dxᵢ is unitless, and in particular is 1 when dxᵢ represents a "unit effect" of state variable xᵢ; that is, the change in xᵢ that produces a unit change in some quantity of interest in the system being simulated. Why unity (1)? Aside from normalizing the values, unity "grounds" the weighted error to the user-specified accuracy. A weighting can be applied individually to each state variable, but typically it is done approximately by combining the known type of the variable (e.g. length, angle) with a "characteristic scale" for that quantity. For example, if a "characteristic length" for the system being simulated is 0.1 meters, and x₀ is a length variable measured in meters, then w₀ should be 10 so that w₀*dx₀=1 when dx₀=0.1. For angles representing pointing accuracy (say a camera direction) we typically assume a "characteristic angle" is one radian (about 60 degrees), so if x₁ is a pointing direction then w₁=1 is an appropriate weight. We can now scale an error vector e=[dx₀ dx₁] to a unitless fractional error vector f=[w₀*dx₀ w₁*dx₁]. Now to achieve a given accuracy a, say a=.0001, we need only check that |fᵢ|<=a for each element i of f. Further, this gives us a quantitative measure of "worst accuracy" that we can use to increase or reduce size of the next attempted step, so that we will just achieve the required accuracy but not much more. We'll be more precise about this below.

Some subtleties for second-order dynamic systems

Systems governed by 2nd-order differential equations are typically split into second order (configuration) variables q, and rate (velocity) variables v, where the time derivatives qdot of q are linearly related to v by the kinematic differential equation qdot = dq/dt = N(q)*v. Velocity variables are chosen to be physically significant, but configuration variables may be chosen for convenience and do not necessarily have direct physical interpretation. For examples, quaternions are chosen as a numerically stable orientation representation. This is problematic for choosing weights which must be done by physical reasoning as sketched above. We resolve this by introducing the notion of "quasi-coordinates" ꝗ (pronounced "qbar") which are defined by the equation ꝗdot = dꝗ/dt = v. Other than time scaling, quasi-coordinates have the same units as their corresponding velocity variables. That is, for weighting we need to think of the configuration coordinates in the same physical space as the velocity variables; weight those by their physical significance; and then map back to an instantaneous weighting on the actual configuration variables q. This mapping is performed automatically; you need only to be concerned about physical weightings.

Note that generalized quasi-coordinates can only be defined locally for a particular configuration q. There is in general no meaningful set of n generalized coordinates which can be differentiated with respect to time to yield v. For example, the Hairy Ball Theorem implies that it is not possible for three orientation variables to represent all 3D rotations without singularities, yet three velocity variables can represent angular velocity in 3D without singularities.

To summarize, separate weights can be provided for each of

  • n generalized quasi-coordinates (configuration variables in the velocity variable space), and
  • nz miscellaneous continuous state variables z.

Weights on the generalized velocity variables v (= dꝗ/dt) are derived directly from the weights on , weighted by a characteristic time. Weights on the actual nq generalized coordinates can be calculated efficiently from weights on the quasi-coordinates (details below).

How the weights are used

The errors in the and z variables are weighted by the diagonal elements of diagonal weighting matrices Wꝗ and Wz, respectively. (The block-diagonal weighting matrix Wq on the original generalized coordinates q is calculated from N and Wꝗ; see below.) In the absence of other information, the default for all weighting values is one, so Wꝗ and Wz are n × n and nz × nz identity matrices. The weighting matrix Wv for the velocity variables is just Wv = τ*Wꝗ where τ is a "characteristic time" for the system, that is, a quantity in time units that represents a significant evolution of the trajectory. This serves to control the accuracy with which velocity is determined relative to configuration. Note that larger values of τ are more conservative since they increase the velocity weights. Typically we use τ=1.0 or 0.1 seconds for human-scale mechanical systems.

The weighting matrices Wq, Wv, and Wz are used to compute a weighted infinity norm as follows. Although Wv and Wz are constant, the actual weightings may be state dependent for relative-error calculations. Define block diagonal error weighting matrix E=diag(Eq,Ev,Ez) as follows:

  Eq = Wq
  Ev: Ev(i,i) = { min(Wv(i,i), 1/|vᵢ|)     if vᵢ is relative
                { Wv(i,i)                  if vᵢ is absolute
  Ez: Ez(i,i) = { min(Wz(i,i), 1/|zᵢ|)     if zᵢ is relative
                { Wz(i,i)                  if zᵢ is absolute

(Ev and Ez are diagonal.) A v or z will be maintained to relative accuracy unless (a) it is "close" to zero (less than 1), or (b) the variable has been defined as requiring absolute accuracy. Position variables q are always maintained to absolute accuracy (see [Sherman 2011] for rationale).

Now given an error estimate vector e=[eq ev ez], the vector f=E*e can be considered to provide a unitless fractional error for each of the state variables. To achieve a given user-specified accuracy a, we require that norm_inf(f) <= a. That is, no element of f can have absolute value larger than a. We also use f to determine an ideal next step size using an appropriate integrator-specific computation.

Determining weights for q

The kinematic differential equations qdot=N(q)*v employ an nq × n matrix N. By construction, this relationship is invertible using N's left pseudo-inverse N⁺ so that v=N⁺ qdot and N⁺ N = I (the identity matrix); however, N N⁺ != I, as N has more rows than columns generally. [Nikravesh 1988] shows how such a matrix N can be determined and provides more information. Given this relationship between N and N⁺, we can relate weighted errors in configuration coordinates q to weighted errors in generalized quasi-coordinates , as the following derivation shows:

           v = N⁺ qdot         Inverse kinematic differential equation
       dꝗ/dt = N⁺ dq/dt        Use synonyms for v and qdot
          dꝗ = N⁺ dq           Change time derivatives to differentials
       Wꝗ dꝗ = Wꝗ N⁺ dq        Pre-multiply both sides by Wꝗ
     N Wꝗ dꝗ = N Wꝗ N⁺ dq      Pre-multiply both sides by N
     N Wꝗ dꝗ = Wq dq           Define Wq := N Wꝗ N⁺
      N Wꝗ v = Wq qdot         Back to time derivatives.

The last two equations show that Wq as defined above provides the expected relationship between the weighted or v variables in velocity space and the weighted q or qdot (resp.) variables in configuration space.

Finally, note that a diagonal entry of one of the weighting matrices can be set to zero to disable error estimation for that state variable (i.e., auxiliary variable or configuration/velocity variable pair), but that setting an entry to a negative value will cause an exception to be thrown when the integrator is initialized.

const Eigen::VectorXd & get_generalized_state_weight_vector () const
 Gets the weighting vector (equivalent to a diagonal matrix) applied to weighting both generalized coordinate and velocity state variable errors, as described in the group documentation. More...
 
Eigen::VectorBlock< Eigen::VectorXd > get_mutable_generalized_state_weight_vector ()
 Gets a mutable weighting vector (equivalent to a diagonal matrix) applied to weighting both generalized coordinate and velocity state variable errors, as described in the group documentation. More...
 
const Eigen::VectorXd & get_misc_state_weight_vector () const
 Gets the weighting vector (equivalent to a diagonal matrix) for weighting errors in miscellaneous continuous state variables z. More...
 
Eigen::VectorBlock< Eigen::VectorXd > get_mutable_misc_state_weight_vector ()
 Gets a mutable weighting vector (equivalent to a diagonal matrix) for weighting errors in miscellaneous continuous state variables z. More...
 

Protected Member Functions

virtual void DoResetStatistics ()
 Resets any statistics particular to a specific integrator. More...
 
void CalcTimeDerivatives (const Context< T > &context, ContinuousState< T > *dxdt)
 Evaluates the derivative function (and updates call statistics). More...
 
template<typename U >
void CalcTimeDerivatives (const System< U > &system, const Context< U > &context, ContinuousState< U > *dxdt)
 Evaluates the derivative function (and updates call statistics). More...
 
void set_accuracy_in_use (double accuracy)
 Sets the working ("in use") accuracy for this integrator. More...
 
void InitializeAccuracy (double default_accuracy, double loosest_accuracy, double max_step_fraction)
 Generic code for validating (and resetting, if need be) the integrator working accuracy for error controlled integrators. More...
 
bool StepOnceErrorControlledAtMost (const T &dt_max)
 Default code for advancing the continuous state of the system by a single step of dt_max (or smaller, depending on error control). More...
 
CalcStateChangeNorm (const ContinuousState< T > &dx_state) const
 Computes the infinity norm of a change in continuous state. More...
 
std::pair< bool, T > CalcAdjustedStepSize (const T &err, const T &attempted_step_size, bool *at_minimum_step_size) const
 Calculates adjusted integrator step sizes toward keeping state variables within error bounds on the next integration step. More...
 
virtual void DoInitialize ()
 Derived classes can override this method to perform special initialization. More...
 
virtual void DoReset ()
 Derived classes can override this method to perform routines when Reset() is called. More...
 
virtual bool DoStep (const T &dt)=0
 Derived classes must implement this method to (1) integrate the continuous portion of this system forward by a single step of size dt and (2) set the error estimate (via get_mutable_error_estimate()). More...
 
ContinuousState< T > * get_mutable_error_estimate ()
 Gets an error estimate of the state variables recorded by the last call to StepOnceFixedSize(). More...
 
void set_actual_initial_step_size_taken (const T &dt)
 
void set_smallest_adapted_step_size_taken (const T &dt)
 Sets the size of the smallest-step-taken statistic as the result of a controlled integration step adjustment. More...
 
void set_largest_step_size_taken (const T &dt)
 
void set_ideal_next_step_size (const T &dt)
 

Detailed Description

template<class T>
class drake::systems::IntegratorBase< T >

An abstract class for an integrator for ODEs and DAEs as represented by a Drake System.

Integrators solve initial value problems of the form:

ẋ(t) = f(t, x(t)) with f : ℝ × ℝⁿ → ℝⁿ

(i.e., f() is an ordinary differential equation) given initial conditions (t₀, x₀). Thus, integrators advance the continuous state of a dynamical system forward in time.

Apart from solving initial value problems, for which the integrator is a key component of a simulator, integrators can also be used to solve boundary value problems (via numerical methods like the Multiple Shooting Method) and trajectory optimization problems (via numerical methods like direct transcription). This class and its derivatives were developed primarily toward the former application (through IntegrateAtMost() and the Simulator class). However, the IntegratorBase architecture was developed to support these ancillary applications as well using the IntegrateWithMultipleSteps() and IntegrateWithSingleFixedStep() methods; the latter permits the caller to advance time using fixed steps in applications where variable stepping would be deleterious (e.g., direct transcription).

A natural question for a user to ask of an integrator is: Which scheme (method) should be applied to a particular problem? The answer is whichever one most quickly computes the solution to the desired accuracy! Selecting an integration scheme for a particular problem is presently an artform. As examples of some selection criteria: multistep methods generally work poorly when events (that require state reinitializations) are common, symplectic methods generally work well at maintaining stability for large integration steps, and stiff integrators are often best for computationally stiff systems. If ignorant as to the characteristics of a particular problem, it is often best to start with an explicit, Runge-Kutta type method. Statistics collected by the integrator can help diagnose performance issues and possibly point to the use of a different integration scheme.

Some systems are known to exhibit "computational stiffness", by which it is meant that (excessively) small integration steps are necessary for purposes of stability: in other words, steps must be taken smaller than that required to achieve a desired accuracy over a particular interval. Thus, the nature of computationally stiff problems is that the solution to the ODE is smooth in the interval of stiffness (in contrast, some problems possess such high frequency dynamics that very small steps are simply necessary to capture the solution accurately). Implicit integrators are the go-to approach for solving computationally stiff problems, but careful consideration is warranted. Implicit integrators typically require much more computation than non-implicit (explicit) integrators, stiffness might be an issue on only a very small time interval, and some problems might be only "moderately stiff". Put another way, applying an implicit integrator to a potentially stiff problem might not yield faster computation. The first chapter of [Hairer, 1996] illustrates the issues broached in this paragraph using various examples.

Established methods for integrating ordinary differential equations invariably make provisions for estimating the "local error" (i.e., the error over a small time interval) of a solution. Although the relationship between local error and global error (i.e., the accumulated error over multiple time steps) can be tenuous, such error estimates can allow integrators to work adaptively, subdividing time intervals as necessary (if, e.g., the system is particularly dynamic or stationary in an interval). Even for applications that do not recommend such adaptive integration- like direct transcription methods for trajectory optimization- error estimation allows the user to assess the accuracy of the solution.

IntegratorBase provides numerous settings and flags that can leverage problem-specific information to speed integration and/or improve integration accuracy. As an example, set_maximum_step_size() allows the user to prevent overly large integration steps (that integration error control alone might be insufficient to detect). As noted previously, IntegratorBase also collects a plethora of statistics that can be used to diagnose poor integration performance. For example, a large number of shrinkages due to error control could indicate that a system is computationally stiff.

  • [Hairer, 1996] E. Hairer and G. Wanner. Solving Ordinary Differential Equations II (Stiff and Differential-Algebraic Problems). Springer, 1996.
Template Parameters
TThe vector element type, which must be a valid Eigen scalar.

Member Enumeration Documentation

enum StepResult

Status returned by StepOnceAtMost().

When a step is successful, it will return an indication of what caused it to stop where it did. When unsuccessful it will throw an exception so you won't see any return value. When return of control is due ONLY to reaching a publish time, (status is kReachedPublishTime) the context may return an interpolated value at an earlier time.

Note: the simulation step must always end at an update time but can end after a publish time.

Enumerator
kReachedPublishTime 

Indicates a publish time has been reached but not an update time.

kReachedZeroCrossing 

Localized an event; this is the before state (interpolated).

kReachedUpdateTime 

Indicates that integration terminated at an update time.

kTimeHasAdvanced 

User requested control whenever an internal step is successful.

kReachedBoundaryTime 

Reached the desired integration time without reaching an update time.

kReachedStepLimit 

Took maximum number of steps without finishing integrating over the interval.

Constructor & Destructor Documentation

IntegratorBase ( const IntegratorBase< T > &  )
delete
IntegratorBase ( IntegratorBase< T > &&  )
delete
IntegratorBase ( const System< T > &  system,
Context< T > *  context = nullptr 
)
inlineexplicit

Maintains references to the system being integrated and the context used to specify the initial conditions for that system (if any).

Parameters
systemA reference to the system to be integrated; the integrator will maintain a reference to the system in perpetuity, so the integrator must not outlive the system.
contextA pointer to a writeable context (nullptr is ok, but a non-null pointer must be set before Initialize() is called). The integrator will advance the system state using the pointer to this context. The pointer to the context will be maintained internally. The integrator must not outlive the context.
virtual ~IntegratorBase ( )
virtualdefault

Destructor.

Here is the caller graph for this function:

Member Function Documentation

std::pair< bool, T > CalcAdjustedStepSize ( const T &  err,
const T &  attempted_step_size,
bool at_minimum_step_size 
) const
protected

Calculates adjusted integrator step sizes toward keeping state variables within error bounds on the next integration step.

Note that it is not guaranteed that the (possibly) reduced step size will keep state variables within error bounds; however, the process of (1) taking a trial integration step, (2) calculating the error, and (3) adjusting the step size can be repeated until convergence.

Parameters
errThe norm of the integrator error that was computed using attempted_step_size.
attempted_step_sizeThe step size that was attempted.
[in,out]at_minimum_step_sizeIf true on entry, the error control mechanism is not allowed to shrink the step because the integrator is stepping at the minimum step size (note that this condition will only occur if get_throw_on_minimum_step_size_violation() == false- an exception would be thrown otherwise). If true on entry and false on exit, the error control mechanism has managed to increase the step size above the working minimum; if true on entry and true on exit, error control would like to shrink the step size but cannot. If false on entry and true on exit, error control shrank the step to the working minimum step size.
Returns
a pair of types bool and T; the bool will be set to true if the integration step was to be considered successful and false otherwise. The value of the T type will be set to the recommended next step size.

Here is the call graph for this function:

Here is the caller graph for this function:

T CalcStateChangeNorm ( const ContinuousState< T > &  dx_state) const
protected

Computes the infinity norm of a change in continuous state.

We use the infinity norm to capture the idea that, by providing accuracy requirements, the user can indirectly specify error tolerances that act to limit the largest error in any state vector component.

Returns
the norm (a non-negative value)

Here is the call graph for this function:

Here is the caller graph for this function:

void CalcTimeDerivatives ( const Context< T > &  context,
ContinuousState< T > *  dxdt 
)
inlineprotected

Evaluates the derivative function (and updates call statistics).

Subclasses should call this function rather than calling system.CalcTimeDerivatives() directly.

Here is the caller graph for this function:

void CalcTimeDerivatives ( const System< U > &  system,
const Context< U > &  context,
ContinuousState< U > *  dxdt 
)
inlineprotected

Evaluates the derivative function (and updates call statistics).

Subclasses should call this function rather than calling system.CalcTimeDerivatives() directly. This version of this function exists to allow integrators to count AutoDiff'd systems in derivative function evaluations.

virtual void DoInitialize ( )
inlineprotectedvirtual

Derived classes can override this method to perform special initialization.

This method is called during the Initialize() method. This default method does nothing.

Here is the caller graph for this function:

virtual void DoReset ( )
inlineprotectedvirtual

Derived classes can override this method to perform routines when Reset() is called.

This default method does nothing.

Here is the caller graph for this function:

virtual void DoResetStatistics ( )
inlineprotectedvirtual

Resets any statistics particular to a specific integrator.

The default implementation of this function does nothing. If your integrator collects its own statistics, you should re-implement this method and reset them there.

Here is the caller graph for this function:

virtual bool DoStep ( const T &  dt)
protectedpure virtual

Derived classes must implement this method to (1) integrate the continuous portion of this system forward by a single step of size dt and (2) set the error estimate (via get_mutable_error_estimate()).

This method is called during the default Step() method.

Parameters
dtThe integration step to take.
Returns
true if successful, false if the integrator was unable to take a single step of size dt (due to, e.g., an integrator convergence failure).
Postcondition
If the time on entry is denoted t, the time and state will be advanced to t+dt if the method returns true; otherwise, the time and state should be reset to those at t.
Warning
It is expected that DoStep() will return true for some, albeit possibly very small, positive value of dt. The derived integrator's stepping algorithm can make this guarantee, for example, by switching to an algorithm not subject to convergence failures (e.g., explicit Euler) for very small step sizes.

Here is the caller graph for this function:

double get_accuracy_in_use ( ) const
inline

Gets the accuracy in use by the integrator.

This number may differ from the target accuracy if, for example, the user has requested an accuracy not attainable or not recommended for the particular integrator.

Here is the caller graph for this function:

const T& get_actual_initial_step_size_taken ( ) const
inline

The actual size of the successful first step.

Here is the caller graph for this function:

const Context<T>& get_context ( ) const
inline

Returns a const reference to the internally-maintained Context holding the most recent state in the trajectory.

This is suitable for publishing or extracting information about this trajectory step.

Here is the caller graph for this function:

const ContinuousState<T>* get_error_estimate ( ) const
inline

Gets the error estimate (used only for integrators that support error estimation).

If the integrator does not support error estimation, nullptr is returned.

Here is the caller graph for this function:

virtual int get_error_estimate_order ( ) const
pure virtual

Derived classes must override this function to return the order of the integrator's error estimate.

If the integrator does not provide an error estimate, the derived class implementation should return 0.

Implemented in ImplicitEulerIntegrator< T >, SemiExplicitEulerIntegrator< T >, RungeKutta3Integrator< T >, ExplicitEulerIntegrator< T >, and RungeKutta2Integrator< T >.

Here is the caller graph for this function:

bool get_fixed_step_mode ( ) const
inline

Gets whether an integrator is running in fixed step mode.

If the integrator does not support error estimation, this function will always return true. If the integrator runs in fixed step mode, it will always take the maximum step size directed (which may be that determined by get_maximum_step_size() or may be smaller, as directed by, e.g., Simulator for event handling purposes).

See also
set_fixed_step_mode()

Here is the caller graph for this function:

const Eigen::VectorXd& get_generalized_state_weight_vector ( ) const
inline

Gets the weighting vector (equivalent to a diagonal matrix) applied to weighting both generalized coordinate and velocity state variable errors, as described in the group documentation.

Only used for integrators that support error estimation.

Here is the caller graph for this function:

const T& get_ideal_next_step_size ( ) const
inline

Return the step size the integrator would like to take next, based primarily on the integrator's accuracy prediction.

This value will not be computed for integrators that do not support error estimation and NaN will be returned.

Here is the caller graph for this function:

const T& get_initial_step_size_target ( ) const
inline

Gets the target size of the first integration step.

You can find out what step size was actually used for the first integration step with get_actual_initial_step_size_taken().

See also
request_initial_step_size_target()

Here is the caller graph for this function:

const T& get_largest_step_size_taken ( ) const
inline

The size of the largest step taken since the last Initialize() or ResetStatistics() call.

Here is the caller graph for this function:

const T& get_maximum_step_size ( ) const
inline

Gets the maximum step size that may be taken by this integrator.

This is a soft maximum: the integrator may stretch it by as much as 1% to hit a discrete event.

See also
set_requested_minimum_step_size()

Here is the caller graph for this function:

const Eigen::VectorXd& get_misc_state_weight_vector ( ) const
inline

Gets the weighting vector (equivalent to a diagonal matrix) for weighting errors in miscellaneous continuous state variables z.

Only used for integrators that support error estimation.

Here is the caller graph for this function:

Context<T>* get_mutable_context ( )
inline

Returns a mutable pointer to the internally-maintained Context holding the most recent state in the trajectory.

Here is the caller graph for this function:

ContinuousState<T>* get_mutable_error_estimate ( )
inlineprotected

Gets an error estimate of the state variables recorded by the last call to StepOnceFixedSize().

If the integrator does not support error estimation, this function will return nullptr.

Eigen::VectorBlock<Eigen::VectorXd> get_mutable_generalized_state_weight_vector ( )
inline

Gets a mutable weighting vector (equivalent to a diagonal matrix) applied to weighting both generalized coordinate and velocity state variable errors, as described in the group documentation.

Only used for integrators that support error estimation. Returns a VectorBlock to make the values mutable without permitting changing the size of the vector. Requires re-initializing the integrator after calling this method; if Initialize() is not called afterward, an exception will be thrown when attempting to call StepOnceAtMost(). If the caller sets one of the entries to a negative value, an exception will be thrown when the integrator is initialized.

Eigen::VectorBlock<Eigen::VectorXd> get_mutable_misc_state_weight_vector ( )
inline

Gets a mutable weighting vector (equivalent to a diagonal matrix) for weighting errors in miscellaneous continuous state variables z.

Only used for integrators that support error estimation. Returns a VectorBlock to make the values mutable without permitting changing the size of the vector. Requires re-initializing the integrator after calling this method. If Initialize() is not called afterward, an exception will be thrown when attempting to call StepOnceAtMost(). If the caller sets one of the entries to a negative value, an exception will be thrown when the integrator is initialized.

int64_t get_num_derivative_evaluations ( ) const
inline

Returns the number of ODE function evaluations (calls to CalcTimeDerivatives()) since the last call to ResetStatistics() or Initialize().

This count includes all such calls including (1) those necessary to compute Jacobian matrices; (2) those used in rejected integrated steps (for, e.g., purposes of error control); (3) those used strictly for integrator error estimation; and (4) calls that exhibit little cost (due to results being cached).

int64_t get_num_step_shrinkages_from_error_control ( ) const
inline

Gets the number of step size shrinkages due to failure to meet targeted error tolerances, since the last call to ResetStatistics or Initialize().

int64_t get_num_step_shrinkages_from_substep_failures ( ) const
inline

Gets the number of step size shrinkages due to sub-step failures (e.g., integrator convergence failures) since the last call to ResetStatistics() or Initialize().

int64_t get_num_steps_taken ( ) const
inline

The number of integration steps taken since the last Initialize() or ResetStatistics() call.

int64_t get_num_substep_failures ( ) const
inline

Gets the number of failed sub-steps (implying one or more step size reductions was required to permit solving the necessary nonlinear system of equations).

const T& get_previous_integration_step_size ( ) const
inline

Gets the size of the last (previous) integration step.

If no integration steps have been taken, value will be NaN.

const T& get_requested_minimum_step_size ( ) const
inline

Gets the requested minimum step size h_min for this integrator.

See also
set_requested_minimum_step_size()
get_working_minimum_step_size(T)
const T& get_smallest_adapted_step_size_taken ( ) const
inline

The size of the smallest step taken as the result of a controlled integration step adjustment since the last Initialize() or ResetStatistics() call.

This value will be NaN for integrators without error estimation.

Here is the caller graph for this function:

double get_stretch_factor ( ) const
inline

Gets the stretch factor (> 1), which is multiplied by the maximum (typically user-designated) integration step size to obtain the amount that the integrator is able to stretch the maximum time step toward hitting an upcoming publish or update event in IntegrateAtMost().

See also
IntegrateAtMost()

Here is the caller graph for this function:

const System<T>& get_system ( ) const
inline

Gets a constant reference to the system that is being integrated (and was provided to the constructor of the integrator).

Here is the caller graph for this function:

double get_target_accuracy ( ) const
inline

Gets the target accuracy.

See also
get_accuracy_in_use()

Here is the caller graph for this function:

bool get_throw_on_minimum_step_size_violation ( ) const
inline

Reports the current setting of the throw_on_minimum_step_size_violation flag.

See also
set_throw_on_minimum_step_size_violation().
T get_working_minimum_step_size ( ) const
inline

Gets the current value of the working minimum step size h_work(t) for this integrator, which may vary with the current time t as stored in the integrator's context.

See this section for more detail.

Here is the caller graph for this function:

void Initialize ( )
inline

An integrator must be initialized before being used.

The pointer to the context must be set before Initialize() is called (or an std::logic_error will be thrown). If Initialize() is not called, an exception will be thrown when attempting to call StepOnceAtMost(). To reinitialize the integrator, Reset() should be called followed by Initialize().

Exceptions
std::logic_errorIf the context has not been set or a user-set parameter has been set illogically (i.e., one of the weighting matrix coefficients is set to a negative value- this check is only performed for integrators that support error estimation; the maximum step size is smaller than the minimum step size; the requested initial step size is outside of the interval [minimum step size, maximum step size]).
See also
Reset()

Here is the caller graph for this function:

void InitializeAccuracy ( double  default_accuracy,
double  loosest_accuracy,
double  max_step_fraction 
)
inlineprotected

Generic code for validating (and resetting, if need be) the integrator working accuracy for error controlled integrators.

This method is intended to be called from an integrator's DoInitialize() method.

Parameters
default_accuracya reasonable default accuracy setting for this integrator.
loosest_accuracythe loosest accuracy that this integrator should support.
max_step_fractiona fraction of the maximum step size to use when setting the integrator accuracy and the user has not specified accuracy directly.
Exceptions
std::logic_errorif neither the initial step size target nor the maximum step size has been set.
IntegratorBase< T >::StepResult IntegrateAtMost ( const T &  publish_dt,
const T &  update_dt,
const T &  boundary_dt 
)

Integrates the system forward in time by a single step with step size subject to integration error tolerances (assuming that the integrator supports error estimation).

The integrator must already have been initialized or an exception will be thrown. The context will be integrated forward by an amount that will never exceed the minimum of publish_dt, update_dt, and 1.01 * get_maximum_step_size().

Parameters
publish_dtThe step size, >= 0.0 (exception will be thrown if this is not the case) at which the next publish will occur.
update_dtThe step size, > 0.0 (exception will be thrown if this is not the case) at which the next update will occur.
boundary_dtThe step size, >= 0.0 (exception will be thrown if this is not the case) marking the end of the user-designated simulated interval.
Exceptions
std::logic_errorIf the integrator has not been initialized or one of publish_dt, update_dt, or boundary_dt is negative.
Returns
The reason for the integration step ending.
Warning
Users should generally not call this function directly; within simulation circumstances, users will typically call Simulator::StepTo(). In other circumstances, users will typically call IntegratorBase::IntegrateWithMultipleSteps().

This method at a glance:

  • For integrating ODEs/DAEs via Simulator
  • Supports fixed step and variable step integration schemes
  • Takes only a single step forward.

Here is the call graph for this function:

Here is the caller graph for this function:

void IntegrateWithMultipleSteps ( const T &  dt)
inline

Stepping function for integrators operating outside of Simulator that advances the continuous state exactly by dt.

This method is designed for integrator users that do not wish to consider publishing or discontinuous, mid-interval updates. This method will step the integrator multiple times, as necessary, to attain requested error tolerances and to ensure the integrator converges.

Warning
Users should simulate systems using Simulator::StepTo() in place of this function (which was created for off-simulation purposes), generally.
Parameters
dtThe non-negative integration step to take.
Exceptions
std::logic_errorIf the integrator has not been initialized or dt is negative.
See also
IntegrateAtMost(), which is designed to be operated by Simulator and accounts for publishing and state reinitialization.
IntegrateWithSingleStep(), which is also designed to be operated outside of Simulator, but throws an exception if the integrator cannot advance time by dt in a single step.

This method at a glance:

  • For integrating ODEs/DAEs not using Simulator
  • Supports fixed step and variable step integration schemes
  • Takes as many steps as necessary until time has advanced by dt
void IntegrateWithSingleFixedStep ( const T &  dt)
inline

Stepping function for integrators operating outside of Simulator that advances the continuous state exactly by dt and using a single fixed step.

This method is designed for integrator users that do not wish to consider publishing or discontinuous, mid-interval updates. One such example application is that of direct transcription for trajectory optimization, for which the integration process should be consistent: it should execute the same sequence of arithmetic operations for all values of the nonlinear programming variables. In keeping with the naming semantics of this function, error controlled integration is not supported (though error estimates will be computed for integrators that support that feature), which is a minimal requirement for "consistency".

Warning
Users should simulate systems using Simulator::StepTo() in place of this function (which was created for off-simulation purposes), generally.
Parameters
dtThe non-negative integration step to take.
Exceptions
std::logic_errorIf the integrator has not been initialized or dt is negative or if the integrator is not operating in fixed step mode.
std::runtime_errorIf the integrator was unable to take a step of the requested size.
See also
IntegrateAtMost(), which is designed to be operated by Simulator and accounts for publishing and state reinitialization.
IntegrateWithMultipleSteps(), which is also designed to be operated outside of Simulator, but will take as many integration steps as necessary until time has been stepped forward by dt.

This method at a glance:

  • For integrating ODEs/DAEs not using Simulator
  • Fixed step integration (no step size reductions for error control or integrator convergence)
  • Takes only a single step forward.
bool is_initialized ( ) const
inline

Indicates whether the integrator has been initialized.

IntegratorBase& operator= ( const IntegratorBase< T > &  )
delete
IntegratorBase& operator= ( IntegratorBase< T > &&  )
delete
void request_initial_step_size_target ( const T &  step_size)
inline

Request that the first attempted integration step have a particular size.

If no request is made, the integrator will estimate a suitable size for the initial step attempt. If the integrator does not support error control, this method will throw a std::logic_error (call supports_error_estimation() to verify before calling this method). For variable-step integration, the initial target will be treated as a maximum step size subject to accuracy requirements and event occurrences. You can find out what size actually worked with get_actual_initial_step_size_taken().

Exceptions
std::logic_errorIf the integrator does not support error estimation.

Here is the caller graph for this function:

void Reset ( )
inline

Resets the integrator to initial values, i.e., default construction values.

void reset_context ( Context< T > *  context)
inline

Replace the pointer to the internally-maintained Context with a different one.

This is useful for supplying a new set of initial conditions or wiping out the current context (by passing in a null pointer). You should invoke Initialize() after replacing the Context unless the context is null.

Parameters
contextThe pointer to the new context or nullptr to wipe out the current context without replacing it with another.
void ResetStatistics ( )
inline

Forget accumulated statistics.

These are reset to the values they have post construction or immediately after Initialize().

Here is the caller graph for this function:

void set_accuracy_in_use ( double  accuracy)
inlineprotected

Sets the working ("in use") accuracy for this integrator.

The working accuracy may not be equivalent to the target accuracy when the latter is too loose or tight for an integrator's capabilities.

See also
get_accuracy_in_use()
get_target_accuracy()

Here is the caller graph for this function:

void set_actual_initial_step_size_taken ( const T &  dt)
inlineprotected

Here is the caller graph for this function:

void set_fixed_step_mode ( bool  flag)
inline

Sets an integrator with error control to fixed step mode.

If the integrator runs in fixed step mode, it will always take the maximum step size directed (which may be that determined by get_maximum_step_size(), or may be smaller, as directed by, e.g., Simulator for event handling purposes).

Exceptions
std::logic_errorif integrator does not support error estimation and flag is set to false.
void set_ideal_next_step_size ( const T &  dt)
inlineprotected
void set_largest_step_size_taken ( const T &  dt)
inlineprotected

Here is the caller graph for this function:

void set_maximum_step_size ( const T &  max_step_size)
inline

Sets the maximum step size that may be taken by this integrator.

The integrator may stretch the maximum step size by as much as 1% to reach a discrete event. For fixed step integrators, all steps will be taken at the maximum step size unless an event would be missed.

Here is the caller graph for this function:

void set_requested_minimum_step_size ( const T &  min_step_size)
inline

Sets the requested minimum step size h_min that may be taken by this integrator.

No step smaller than this will be taken except under circumstances as described above. This setting will be ignored if it is smaller than the absolute minimum h_floor also described above. Default value is zero.

Parameters
min_step_sizea non-negative value. Setting this value to zero will cause the integrator to use a reasonable value instead (see get_working_minimum_step_size()).
See also
get_requested_minimum_step_size()
get_working_minimum_step_size()
void set_smallest_adapted_step_size_taken ( const T &  dt)
inlineprotected

Sets the size of the smallest-step-taken statistic as the result of a controlled integration step adjustment.

Here is the caller graph for this function:

void set_target_accuracy ( double  accuracy)
inline

Request that the integrator attempt to achieve a particular accuracy for the continuous portions of the simulation.

Otherwise a default accuracy is chosen for you. This may be ignored for fixed-step integration since accuracy control requires variable step sizes. You should call supports_error_estimation() to ensure that the integrator supports this capability before calling this function; if the integrator does not support it, this method will throw an exception.

Integrators vary in the range of accuracy (loosest to tightest) that they can support. If you request accuracy outside the supported range for the chosen integrator it will be quietly adjusted to be in range. You can find out the accuracy setting actually being used using get_accuracy_in_use().

The precise meaning of accuracy is a complicated discussion, but translates roughly to the number of significant digits you want in the results. By convention it is supplied as 10^-digits, meaning that an accuracy of 1e-3 provides about three significant digits. For more information, see [Sherman 2011].

Here is the caller graph for this function:

void set_throw_on_minimum_step_size_violation ( bool  throws)
inline

Sets whether the integrator should throw a std::runtime_error exception when the integrator's step size selection algorithm determines that it must take a step smaller than the minimum step size (for, e.g., purposes of error control).

Default is true. If false, the integrator will advance time and state using the minimum specified step size in such situations. See this section for more detail.

bool StepOnceErrorControlledAtMost ( const T &  dt_max)
protected

Default code for advancing the continuous state of the system by a single step of dt_max (or smaller, depending on error control).

This particular function is designed to be called directly by an error estimating integrator's DoStep() method to effect error-controlled integration. The integrator can effect error controlled integration without calling this method, if the implementer so chooses, but this default method is expected to function well in most circumstances.

Parameters
[in]dt_maxThe maximum step size to be taken. The integrator may take a smaller step than specified to satisfy accuracy requirements, to resolve integrator convergence problems, or to respect the integrator's maximum step size.
Exceptions
std::logic_errorif integrator does not support error estimation.
Note
This function will shrink the integration step as necessary whenever the integrator's DoStep() fails to take the requested step e.g., due to integrator convergence failure.
Returns
true if the full step of size dt_max is taken and false otherwise (i.e., a smaller step than dt_max was taken).

Here is the call graph for this function:

Here is the caller graph for this function:

virtual bool supports_error_estimation ( ) const
pure virtual

Indicates whether an integrator supports error estimation.

Without error estimation, target accuracy will be unused.

Implemented in ImplicitEulerIntegrator< T >, SemiExplicitEulerIntegrator< T >, RungeKutta3Integrator< T >, ExplicitEulerIntegrator< T >, and RungeKutta2Integrator< T >.

Here is the caller graph for this function:


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