system

The System class manages the simulation (integration) of a system whose equations are given by KanesMethod.

Many of the attributes are also properties, and can be directly modified.

Here is the procedure for using this class.

  1. specify your options either via the constructor or via the attributes.

  2. optionally, call generate_ode_function() if you want to customize how the ODE function is generated.

  3. call integrate() to simulate your system.

The simplest usage of this class is as follows. First, we need a KanesMethod object on which we have already invoked kanes_equations():

km = KanesMethod(...)
km.kanes_equations(force_list, body_list)
times = np.linspace(0, 5, 100)
sys = System(km, times=times)
sys.integrate()

In this case, we use defaults for the numerical values of the constants, specified quantities, initial conditions, etc. You probably won’t like these defaults. You can also specify such values via constructor keyword arguments or via the attributes:

sys = System(km,
             initial_conditions={dynamicsymbol('q1'): 0.5},
             times=times)
sys.constants = {symbol('m'): 5.0}
sys.integrate()

To double-check the constants, specifieds, states and times in your problem, look at these properties:

sys.constants_symbols
sys.specifieds_symbols
sys.states
sys.times

In this case, the System generates the numerical ode function for you behind the scenes. If you want to customize how this function is generated, you must call generate_ode_function on your own:

sys = System(KM)
sys.generate_ode_function(generator='cython')
sys.integrate()
class pydy.system.System(eom_method, constants=None, specifieds=None, ode_solver=None, initial_conditions=None, times=None)[source]

See the class’s attributes for a description of the arguments to this constructor.

The parameters to this constructor are all attributes of the System. Actually, they are properties. With the exception of eom_method, these attributes can be modified directly at any future point.

Parameters
eom_methodsympy.physics.mechanics.KanesMethod

You must have called KanesMethod.kanes_equations() before constructing this System.

constantsdict, optional (default: all 1.0)

This dictionary maps SymPy Symbol objects to floats.

specifiedsdict, optional (default: all 0)

This dictionary maps SymPy Functions of time objects, or tuples of them, to floats, NumPy arrays, or functions of the state and time.

ode_solverfunction, optional (default: scipy.integrate.odeint)

This function computes the derivatives of the states.

initial_conditionsdict, optional (default: all zero)

This dictionary maps SymPy Functions of time objects to floats.

timesarray_like, shape(n,), optional

An array_like object, which contains time values over which equations are integrated. It has to be supplied before System.integrate() can be called.

__init__(eom_method, constants=None, specifieds=None, ode_solver=None, initial_conditions=None, times=None)[source]
property constants

A dict that provides the numerical values for the constants in the problem (all non-dynamics symbols). Keys are the symbols for the constants, and values are floats. Constants that are not specified in this dict are given a default value of 1.0.

property constants_symbols

A set of the symbolic constants (not functions of time) in the system.

property coordinates

Returns a list of the symbolic functions of time representing the system’s generalized coordinates.

property eom_method

This is a sympy.physics.mechanics.KanesMethod. The method used to generate the equations of motion. Read-only.

property evaluate_ode_function

A function generated by generate_ode_function that computes the state derivatives:

x’ = evaluate_ode_function(x, t, args=(…))

This function is used by the ode_solver.

generate_ode_function(**kwargs)[source]

Calls pydy.codegen.ode_function_generators.generate_ode_function with the appropriate arguments, and sets the evaluate_ode_function attribute to the resulting function.

Parameters
kwargs

All other kwargs are passed onto pydy.codegen.ode_function_generators.generate_ode_function(). Don’t specify the specifieds keyword argument though; the System class takes care of those.

Returns
evaluate_ode_functionfunction

A function which evaluates the derivaties of the states.

property initial_conditions

Initial conditions for all states (coordinates and speeds). Keys are the symbols for the coordinates and speeds, and values are floats. Coordinates or speeds that are not specified in this dict are given a default value of 0.0.

integrate(**solver_kwargs)[source]

Integrates the equations evaluate_ode_function() using ode_solver.

It is necessary to have first generated an ode function. If you have not done so, we do so automatically by invoking generate_ode_function(). However, if you want to customize how this function is generated (e.g., change the generator to cython), you can call generate_ode_function() on your own (before calling integrate()).

Parameters
**solver_kwargs

Optional arguments that are passed on to the ode_solver.

Returns
x_historynp.array, shape(num_integrator_time_steps, 2)

The trajectory of states (coordinates and speeds) through the requested time interval. num_integrator_time_steps is either len(times) if len(times) > 2, or is determined by the ode_solver.

property ode_solver

A function that performs forward integration. It must have the same signature as scipy.integrate.odeint, which is:

x_history = ode_solver(f, x0, t, args=(args,))

where f is a function f(x, t, args), x0 are the initial conditions, x_history is the state time history, x is the state, t is the time, and args is a keyword argument takes arguments that are then passed to f. The default solver is odeint.

property specifieds

A dict that provides numerical values for the specified quantities in the problem (all dynamicsymbols that are not defined by the equations of motion). There are two possible formats. (1) is more flexible, but (2) is more efficient (by a factor of 3).

(1) Keys are the symbols for the specified quantities, or a tuple of symbols, and values are the floats, arrays of floats, or functions that generate the values. If a dictionary value is a function, it must have the same signature as f(x, t), the ode right-hand-side function (see the documentation for the ode_solver attribute). You needn’t provide values for all specified symbols. Those for which you do not give a value will default to 0.0.

(2) There are two keys: ‘symbols’ and ‘values’. The value for ‘symbols’ is an iterable of all the specified quantities in the order that you have provided them in ‘values’. Values is an ndarray, whose length is len(sys.specifieds_symbols), or a function of x and t that returns an ndarray (also of length len(sys.specifieds_symbols)). NOTE: You must provide values for all specified symbols. In this case, we do not provide default values.

NOTE: If you switch formats with the same instance of System, you must call generate_ode_function() before calling integrate() again.

Examples

Here are examples for (1). Keys can be individual symbols, or a tuple of symbols. Length of a value must match the length of the corresponding key. Values can be functions that return iterables:

   sys = System(km)
   sys.specifieds = {(a, b, c): np.ones(3), d: lambda x, t: -3 * x[0]}
   sys.specifieds = {(a, b, c): lambda x, t: np.ones(3)}

Here are examples for (2):

   sys.specifieds = {'symbols': (a, b, c, d),
                     'values': np.ones(4)}
   sys.specifieds = {'symbols': (a, b, c, d),
                     'values': lambda x, t: np.ones(4)}
property specifieds_symbols

A set of the dynamicsymbols you must specify.

property speeds

Returns a list of the symbolic functions of time representing the system’s generalized speeds.

property states

Returns a list of the symbolic functions of time representing the system’s states, i.e. generalized coordinates plus the generalized speeds. These are in the same order as used in integration (as passed into evaluate_ode_function) and match the order of the mass matrix and forcing vector.

property times

An array-like object, containing time values over which the equations of motion are integrated, numerically.

The object should be in a format which the integration module to be used can accept.