odetoolbox package

odetoolbox._analysis(indict, disable_stiffness_check: bool = False, disable_analytic_solver: bool = False, preserve_expressions: Union[bool, Iterable[str]] = False, simplify_expression: Optional[str] = None, log_level: Union[str, int] = 30) → Tuple[List[Dict[KT, VT]], odetoolbox.system_of_shapes.SystemOfShapes, List[odetoolbox.shapes.Shape]]

Like analysis(), but additionally returns shape_sys and shapes.

For internal use only!

odetoolbox._find_analytically_solvable_equations(shape_sys, shapes, parameters=None)

Find which equations can be solved analytically (and, conversely, which cannot).

Perform dependency analysis and plot dependency graph.

odetoolbox._find_variable_definition(indict, name: str, order: int) → Optional[str]

Find the definition (as a string in the input dictionary) of variable named name with order order, and return it as a string. Return None if a definition by that name and order could not be found.

odetoolbox._from_json_to_shapes(indict, parameters=None) → Tuple[List[odetoolbox.shapes.Shape], Dict[sympy.core.symbol.Symbol, str]]

Process the input, construct Shape instances.

Parameters:indict – ODE-toolbox input dictionary.
odetoolbox._get_all_first_order_variables(indict) → Iterable[str]

Return a list of variable names, containing those variables that were defined as a first-order ordinary differential equation in the input.

odetoolbox._init_logging(log_level: Union[str, int] = 30)

Initialise message logging.

Parameters:log_level – Sets the logging threshold. Logging messages which are less severe than log_level will be ignored. Log levels can be provided as an integer or string, for example “INFO” (more messages) or “WARN” (fewer messages). For a list of valid logging levels, see https://docs.python.org/3/library/logging.html#logging-levels
odetoolbox._read_global_config(indict)

Process global configuration options.

odetoolbox.analysis(indict, disable_stiffness_check: bool = False, disable_analytic_solver: bool = False, preserve_expressions: Union[bool, Iterable[str]] = False, simplify_expression: Optional[str] = None, log_level: Union[str, int] = 30) → List[Dict[KT, VT]]

The main entry point of the ODE-toolbox API.

Parameters:
  • indict – Input dictionary for the analysis. For details, see https://ode-toolbox.readthedocs.io/en/master/#input
  • disable_stiffness_check – Whether to perform stiffness checking.
  • disable_analytic_solver – Set to True to return numerical solver recommendations, and no propagators, even for ODEs that are analytically tractable.
  • preserve_expressions – Set to True, or a list of strings corresponding to individual variable names, to disable internal rewriting of expressions, and return same output as input expression where possible. Only applies to variables specified as first-order differential equations.
  • simplify_expression – For all expressions expr that are rewritten internally: the contents of this parameter string are evaluated with eval() in Python to obtain the final output expression. Override for custom expression simplification steps. Example: "sympy.logcombine(sympy.powsimp(sympy.expand(expr)))".
  • log_level – Sets the logging threshold. Logging messages which are less severe than log_level will be ignored. Log levels can be provided as an integer or string, for example “INFO” (more messages) or “WARN” (fewer messages). For a list of valid logging levels, see https://docs.python.org/3/library/logging.html#logging-levels
Returns:

The result of the analysis. For details, see https://ode-toolbox.readthedocs.io/en/latest/index.html#output

Submodules

odetoolbox.analytic_integrator module

class odetoolbox.analytic_integrator.AnalyticIntegrator(solver_dict, spike_times: Optional[Dict[str, List[float]]] = None, enable_caching: bool = True)

Bases: odetoolbox.integrator.Integrator

Integrate a dynamical system by means of the propagators returned by ODE-toolbox.

__init__(solver_dict, spike_times: Optional[Dict[str, List[float]]] = None, enable_caching: bool = True)
Parameters:
  • solve_dict – The results dictionary returned by a call to odetoolbox.analysis().
  • spike_times – For each variable, used as a key, the list of times at which a spike occurs.
  • enable_caching – Allow caching of results between requested times.
_update_step(delta_t, initial_values)

Apply propagator to update the state, starting from initial_values, by timestep delta_t.

Parameters:
  • delta_t – Timestep to take.
  • initial_values – A dictionary mapping variable names (as strings) to initial value expressions.
disable_cache_update()

Disallow caching of results between requested times.

enable_cache_update()

Allow caching of results between requested times.

get_all_variable_symbols()
get_value(t)

Get numerical solution of the dynamical system at time t.

Parameters:t – The time to compute the solution for.
reset()

Reset time to zero and state to initial values.

set_initial_values(vals)

Set initial values, i.e. the state of the system at \(t = 0\). This will additionally cause the system state to be reset to \(t = 0\) and the new initial conditions.

Parameters:vals – New initial values.

odetoolbox.config module

class odetoolbox.config.Config

Bases: object

Static class to store global configuration options.

Options are stored in the static dict config. Access using either Config().key or Config()[key] or Config.config[key].

config = {'differential_order_symbol': '__d', 'expression_simplification_threshold': 1000, 'input_time_symbol': 't', 'integration_accuracy_abs': 1e-06, 'integration_accuracy_rel': 1e-06, 'max_step_size': 999.0, 'output_timestep_symbol': '__h', 'sim_time': 0.1, 'simplify_expression': 'sympy.simplify(expr)'}
keys()

odetoolbox.dependency_graph_plotter module

class odetoolbox.dependency_graph_plotter.DependencyGraphPlotter

Bases: object

Use graphviz to plot a dependency graph between state variables.

classmethod plot_graph(shapes, dependency_edges, node_is_lin, fn=None)

Plot graph and write to file.

Parameters:
  • shapes – List of Shape instances.
  • dependency_edges – List of edges returned from dependency analysis.
  • node_is_lin – List of Booleans returned from dependency analysis.
  • fn – Filename to write PNG image as.

odetoolbox.integrator module

class odetoolbox.integrator.Integrator

Bases: object

Integrate a dynamical system by means of the propagators returned by ODE-toolbox (base class).

all_variable_symbols = []
get_sorted_spike_times()

Returns a global, sorted list of spike times.

Return all_spike_times:
 A sorted list of all spike times for all variables.
Return all_spike_times_sym:
 For the spike at time all_spike_times[i], the variables to which that spike applies are listed in all_spike_times_sym[i].
get_spike_times()

Get spike times.

Return spike_times:
 For each variable, used as a key, the list of spike times associated with it.
set_spike_times(spike_times: Optional[Dict[str, List[float]]])

Internally converts to a global, sorted list of spike times.

Parameters:spike_times – For each variable, used as a key, the list of spike times associated with it.

odetoolbox.mixed_integrator module

class odetoolbox.mixed_integrator.MixedIntegrator(numeric_integrator, system_of_shapes, shapes, analytic_solver_dict=None, parameters=None, spike_times=None, random_seed=123, max_step_size=inf, integration_accuracy_abs=1e-06, integration_accuracy_rel=1e-06, sim_time=1.0, alias_spikes=False, debug_plot_dir: Optional[str] = None)

Bases: odetoolbox.integrator.Integrator

Mixed numeric+analytic integrator. Supply with a result from ODE-toolbox analysis; calculates numeric approximation of the solution.

__init__(numeric_integrator, system_of_shapes, shapes, analytic_solver_dict=None, parameters=None, spike_times=None, random_seed=123, max_step_size=inf, integration_accuracy_abs=1e-06, integration_accuracy_rel=1e-06, sim_time=1.0, alias_spikes=False, debug_plot_dir: Optional[str] = None)
Parameters:
  • numeric_integrator – A method from the GSL library for evolving ODEs, e.g. odeiv.step_rk4
  • system_of_shapes – Dynamical system to solve.
  • shapes – List of shapes in the dynamical system.
  • analytic_solver_dict – Analytic solver dictionary from ODE-toolbox analysis result.
  • parameters – Dictionary mapping parameter name (as string) to value expression.
  • spike_times – For each variable, used as a key, the list of times at which a spike occurs.
  • random_seed – Random number generator seed.
  • max_step_size – The maximum step size taken by the integrator.
  • integration_accuracy_abs – Absolute integration accuracy.
  • integration_accuracy_rel – Relative integration accuracy.
  • sim_time – How long to simulate for.
  • alias_spikes – Whether to alias spike times to the numerical integration grid. False means that precise integration will be used for spike times whenever possible. True means that after taking a timestep \(dt\) and arriving at \(t\), spikes from \(\langle t - dt, t]\) will only be processed at time \(t\).
  • debug_plot_dir – If given, enable debug plotting to this directory. If enabled, matplotlib is imported and used for plotting.
integrate_ode(initial_values=None, h_min_lower_bound=5e-09, raise_errors=True, debug=False)

This function computes the average step size and the minimal step size that a given integration method from GSL uses to evolve a certain system of ODEs during a certain simulation time, integration method from GSL and spike train.

Parameters:
  • initial_values – A dictionary mapping variable names (as strings) to initial value expressions.
  • h_min_lower_bound – The minimum acceptable step size.
  • raise_errors – Stop and raise exception when error occurs, or try to continue.
  • debug – Return extra values useful for debugging.
Returns:

Average and minimal step size, and elapsed wall clock time.
integrator_debug_plot(t_log, h_log, y_log, dir)
numerical_jacobian(t, y, params)

Compute the numerical values of the Jacobian matrix at the current time t and state y.

If the dynamics of variables \(x_1, \ldots, x_N\) is defined as \(x_i' = f_i\), then row \(i\) of the Jacobian matrix \(\mathbf{J}_i = \left[\begin{matrix}\frac{\partial f_i}{\partial x_0} & \cdots & \frac{\partial f_i}{\partial x_N}\end{matrix}\right]\).

Parameters:
  • t – Current time.
  • y – Current state vector of the dynamical system.
  • params – GSL parameters (not used here).
Returns:

Jacobian matrix \(\mathbf{J}\).
step(t, y, params)

“Stepping function”: compute the (numerical) value of the derivative of y over time, at the current time t and state y.

Parameters:
  • t – Current time.
  • y – Current state vector of the dynamical system.
  • params – GSL parameters (not used here).
Returns:

Updated state vector
exception odetoolbox.mixed_integrator.ParametersIncompleteException

Bases: Exception

Thrown in case not all parameters are assigned a numerical value before integration was attempted.

odetoolbox.plot_helper module

odetoolbox.plot_helper.import_matplotlib()

Try to import and configure matplotlib. Returns the “mpl” and “plt” packages if the import was successful, or return (None, None) if unsuccessful.

odetoolbox.shapes module

exception odetoolbox.shapes.MalformedInputException

Bases: Exception

Thrown in case an error occurred while processing malformed input.

class odetoolbox.shapes.Shape(symbol, order, initial_values, derivative_factors, inhom_term=0.0, nonlin_term=0.0, lower_bound=None, upper_bound=None)

Bases: object

This class provides a canonical representation of a shape function independently of the way in which the user specified the shape. It assumes a differential equation of the general form (where bracketed superscript \({}^{(n)}\) indicates the \(n\)-th derivative with respect to time):

\[x^{(n)} = N + \sum_{i=0}^{n-1} c_i x^{(i)}\]

Any constant or nonlinear part is here contained in the term N.

In the input and output, derivatives are indicated by adding one prime (single quotation mark) for each derivative order. For example, in the expression

x''' = c0 + c1*x + c2*x' + c3*x'' + x*y + x**2

the symbol of the ODE would be x (i.e. without any qualifiers), order would be 3, derivative_factors would contain the linear part in the form of the list [c1, c2, c3], the inhomogeneous term would be c0 and the nonlinear part x*y + x**2 is stored in nonlin_term.

__init__(symbol, order, initial_values, derivative_factors, inhom_term=0.0, nonlin_term=0.0, lower_bound=None, upper_bound=None)

Perform type and consistency checks and assign arguments to member variables.

Parameters:
  • symbol – Symbolic name of the shape without additional qualifiers like prime symbols or similar.
  • order – Order of the ODE representing the shape.
  • initial_values – Initial values of the ODE representing the shape. The dict contains order many key-value pairs: one for each derivative that occurs in the ODE. The keys are strings created by concatenating the variable symbol with as many single quotation marks (’) as the derivation order. The values are SymPy expressions.
  • derivative_factors – The factors for the derivatives that occur in the ODE. This list has to contain :path:`order` many values, i.e. one for each derivative that occurs in the ODE. The values have to be in ascending order, i.e. [c1, c2, c3] for the given example.
  • inhom_term – Inhomogeneous part of the ODE representing the shape, i.e. c0 for the given example.
  • nonlin_term – Nonlinear part of the ODE representing the shape, i.e. x*y + x**2 for the given example.
__str__()

Return str(self).

classmethod _parse_defining_expression(s: str) → Tuple[str, int, str]

Parse a defining expression, for example, if the ODE-toolbox JSON input file contains the snippet:

::
{
“expression”: “h’ = -g / tau**2 - 2 * h / tau”

}

then the corresponding defining expression is "h' = -g / tau**2 - 2 * h / tau".

This function parses that string and returns the variable name (h), the derivative order (1) and the right-hand side expression (-g / tau**2 - 2 * h / tau").

_sympy_autowrap_helpers = [('Min', -Abs(x - y)/2 + Abs(x + y)/2, [x, y]), ('Max', Abs(x - y)/2 + Abs(x + y)/2, [x, y]), ('Heaviside', (x + Abs(x))/(2*Abs(x) + 1.0e-300), [x])]
_sympy_globals = {'DiracDelta': DiracDelta, 'E': E, 'Float': <class 'sympy.core.numbers.Float'>, 'Function': Function, 'Heaviside': Heaviside, 'Integer': <class 'sympy.core.numbers.Integer'>, 'Pow': <class 'sympy.core.power.Pow'>, 'Symbol': <class 'sympy.core.symbol.Symbol'>, 'acos': acos, 'acosh': acosh, 'asin': asin, 'asinh': asinh, 'atanh': atanh, 'cos': cos, 'cosh': cosh, 'e': E, 'exp': exp, 'log': log, 'max': Max, 'min': Min, 'power': <class 'sympy.core.power.Pow'>, 'sin': sin, 'sinh': sinh, 't': t, 'tan': tan, 'tanh': tanh}
classmethod from_function(symbol: str, definition, max_t=100, max_order=4, all_variable_symbols=None, debug=False) → odetoolbox.shapes.Shape

Create a Shape object given a function of time.

Only functions of time that have a homogeneous ODE equivalent are supported (inhomogeneous ODE functions are not supported).

For a complete description of the algorithm, see https://ode-toolbox.readthedocs.io/en/latest/index.html#converting-direct-functions-of-time Note that this function uses sympy.simplify() rather than the custom simplification expression, because the latter risks that we fail to recognize certain shapes.

Parameters:
  • symbol – The variable name of the shape (e.g. “alpha”, “I”)
  • definition – The definition of the shape (e.g. “(e/tau_syn_in) * t *  exp(-t/tau_syn_in)”)
  • all_variable_symbols – All known variable symbols. None or list of string.
Returns:

The canonical representation of the postsynaptic shape
Example: 
>>> Shape.from_function("I_in", "(e/tau) * t * exp(-t/tau)")
classmethod from_json(indict, all_variable_symbols=None, parameters=None, _debug=False)

Create a Shape instance from an input dictionary.

Parameters:
  • indict – Input dictionary, i.e. one element of the “dynamics” list supplied in the ODE-toolbox input dictionary.
  • all_variable_symbols – All known variable symbols. None or list of string.
  • parameters – An optional dictionary of parameters to their defining expressions.
classmethod from_ode(symbol: str, definition: str, initial_values: dict, all_variable_symbols=None, lower_bound=None, upper_bound=None, parameters=None, debug=False, **kwargs) → odetoolbox.shapes.Shape

Create a Shape object given an ODE and initial values.

Note that shapes are only aware of their own state variables: if an equation for \(x\) depends on another state variable \(y\) of another shape, then \(y\) will be assumed to be a parameter and the term will appear in the inhomogeneous component of \(x\).

Parameters:
  • symbol – The symbol (variable name) of the ODE.
  • definition – The definition of the ODE.
  • initial_values – A dictionary mapping initial values to expressions.
  • all_variable_symbols – All known variable symbols. None or list of string.
Example: 
>>> Shape.from_ode("alpha",
                  "-1/tau**2 * shape_alpha -2/tau * shape_alpha'",
                  {"alpha" : "0", "alpha'" : "e/tau", "0"]})
get_all_variable_symbols(shapes=None, derivative_symbol="'") → List[sympy.core.symbol.Symbol]

Get all variable symbols for this shape and all other shapes in shapes, without duplicates, in no particular order.

Returns:all_symbols
get_initial_value(sym: str)

Get the initial value corresponding to the variable symbol.

Parameters:sym – String representation of a sympy symbol, e.g. “V_m’”
get_state_variables(derivative_symbol="'") → List[sympy.core.symbol.Symbol]

Get all variable symbols for this shape, ordered according to derivative order, up to the shape’s order \(N\): [sym, dsym/dt, d^2sym/dt^2, …, d^(N-1)sym/dt^(N-1)]

is_homogeneous() → bool
Returns:False if and only if the shape has a nonzero right-hand side.
is_lin_const_coeff() → bool
Returns:True if and only if the shape is linear and constant coefficient.
is_lin_const_coeff_in(symbols, parameters=None)
Returns:True if and only if the shape is linear and constant coefficient in those variables passed in symbols.
reconstitute_expr() → sympy.core.expr.Expr

Recreate right-hand side expression from internal representation (linear coefficients, inhomogeneous, and nonlinear parts).

static split_lin_inhom_nonlin(expr, x, parameters=None)

Split an expression into a linear, inhomogeneous and nonlinear part.

For example, in the expression

x''' = c0 + c1*x + c2*x' + c3*x'' + x*y + x**2

lin_factors would contain the linear part in the form of the list [c1, c2, c3], the inhomogeneous term would be c0 and the nonlinear part x*y + x**2 is returned as nonlin_term.

All parameters in parameters are assumed to be constants.

odetoolbox.shapes.is_constant_term(term, parameters: Optional[Mapping[sympy.core.symbol.Symbol, str]] = None)
Returns:True if and only if this term contains only numerical values and parameters; False otherwise.

odetoolbox.singularity_detection module

class odetoolbox.singularity_detection.SingularityDetection

Bases: object

Singularity detection for generated propagator matrix.

Some ordinary differential equations (ODEs) can be solved analytically: an expression for the solution can be readily derived by algebraic manipulation. This allows us to formulate an “exact integrator”, that yields the next state of the system given the current state and the timestep Δt, to floating point (machine) precision [1].

In some cases, an ODE is analytically tractable, but vulnerable to an edge case condition in the generated propagator matrices. Consider the following example: Let the system of ODEs be given by

\[y' = A \cdot y\]

Then the propagator matrix for a timestep \(\Delta t\) is

\[P = \exp(A \cdot \Delta t)\]

which we can use to advance the system

\[y(t + \Delta t) = P \cdot y(t)\]

If \(A\) is of the form:

\[\begin{split}\begin{bmatrix} -a & 0 & 0\\ 1 & -a & 0\\ 0 & 1 & -b \end{bmatrix}\end{split}\]

Then the generated propagator matrix contains denominators that include the factor \(a - b\). When the parameters are chosen such that \(a = b\), a singularity (division by zero fault) occurs. However, the singularity is readily avoided if we assume that \(a = b\) before generating the propagator, i.e. we start out with the matrix

\[\begin{split}\begin{bmatrix} -a & 0 & 0\\ 1 & -a & 0\\ 0 & 1 & -a \end{bmatrix}\end{split}\]

The resulting propagator contains no singularities.

This class detects the potential occurrence of such singularities (potential division by zero) in the generated propagator matrix, which occur under certain choices of parameter values. These choices are reported as “conditions” by the find_singularities() function.

References

[1]Stefan Rotter, Markus Diesmann. Exact digital simulation of time-invariant linear systems with applications to neuronal modeling. Neurobiologie und Biophysik, Institut für Biologie III, Universität Freiburg, Freiburg, Germany Biol. Cybern. 81, 381-402 (1999)
static _filter_valid_conditions(cond, A: sympy.matrices.dense.MutableDenseMatrix)
static _flatten_conditions(cond)

Return a list with conditions in the form of dictionaries

static _generate_singularity_conditions(A: sympy.matrices.dense.MutableDenseMatrix)

The function solve returns a list where each element is a dictionary. And each dictionary entry (condition: expression) corresponds to a condition at which that expression goes to zero. If the expression is quadratic, like let’s say “x**2-1” then the function ‘solve() returns two dictionaries in a list. each dictionary corresponds to one solution. We are then collecting these lists in our own list called ‘condition’.

static _is_matrix_defined_under_substitution(A: sympy.matrices.dense.MutableDenseMatrix, cond: Mapping[KT, VT_co]) → bool

Function to check if a matrix is defined (i.e. does not contain NaN or infinity) after we perform a given set of subsitutions.

Parameters:
  • A (sympy.Matrix) – input matrix
  • cond (Mapping) – mapping from expression that is to be subsituted, to expression to put in its place
static find_singularities(P: sympy.matrices.dense.MutableDenseMatrix, A: sympy.matrices.dense.MutableDenseMatrix)

Find singularities in the propagator matrix \(P\) given the system matrix \(A\).

Parameters:
  • P (sympy.Matrix) – propagator matrix to check for singularities
  • A (sympy.Matrix) – system matrix

odetoolbox.spike_generator module

class odetoolbox.spike_generator.SpikeGenerator

Bases: object

classmethod _generate_homogeneous_poisson_spikes(T: float, rate: float, min_isi: float = 1e-06)

Generate spike trains for the given simulation length. Uses a Poisson distribution to create biologically realistic characteristics of the spike-trains.

Parameters:
  • T – Spikes are generated in the window \([0, T]\).
  • min_isi – Minimum time between two consecutive spikes.
Returns:

spike_times: For each symbol: a list with spike times
Return type:

Dict(str -> List(float))
classmethod _generate_regular_spikes(T: float, rate: float)

Generate spike trains for the given simulation length.

Parameters:T – Spikes are generated in the window \(\langle 0, T]\).
Returns:For each symbol: a list with spike times
Return type:Dict(str -> List(float))
classmethod spike_times_from_json(stimuli, sim_time) → Mapping[str, List[float]]

Read or generate spike times according to a JSON specification.

Returns:spike_times: Each variable symbol is a key in the dictionary, and the list of spike times is the corresponding value. Symbol names use derivative_symbol to indicate differential order.

odetoolbox.stiffness module

class odetoolbox.stiffness.StiffnessTester(system_of_shapes, shapes, analytic_solver_dict=None, parameters=None, stimuli=None, random_seed=123, max_step_size=inf, integration_accuracy_abs=1e-06, integration_accuracy_rel=1e-06, sim_time=100.0, alias_spikes=False)

Bases: object

__init__(system_of_shapes, shapes, analytic_solver_dict=None, parameters=None, stimuli=None, random_seed=123, max_step_size=inf, integration_accuracy_abs=1e-06, integration_accuracy_rel=1e-06, sim_time=100.0, alias_spikes=False)
Parameters:
  • system_of_shapes – Dynamical system to solve.
  • shapes – List of shapes in the dynamical system.
  • analytic_solver_dict – Analytic solver dictionary from ODE-toolbox analysis result.
  • parameters – Dictionary mapping parameter name (as string) to value expression.
  • stimuli – Dictionary containing spiking stimuli.
  • random_seed – Random number generator seed.
  • max_step_size – The maximum step size taken by the integrator.
  • integration_accuracy_abs – Absolute integration accuracy.
  • integration_accuracy_rel – Relative integration accuracy.
  • sim_time – How long to simulate for.
  • alias_spikes – Whether to alias spike times to the numerical integration grid. False means that precise integration will be used for spike times whenever possible. True means that after taking a timestep \(dt\) and arriving at \(t\), spikes from \(\langle t - dt, t]\) will only be processed at time \(t\).
_draw_decision(step_min_imp, step_min_exp, step_average_imp, step_average_exp, machine_precision_dist_ratio=10, avg_step_size_ratio=6)

Decide which is the best integrator to use.

For details, see https://ode-toolbox.readthedocs.io/en/latest/index.html#numeric-solver-selection-criteria

Parameters:
  • step_min_imp – Minimum recommended step size returned for implicit solver.
  • step_min_exp – Minimum recommended step size returned for explicit solver.
  • step_average_imp – Average recommended step size returned for implicit solver.
  • step_average_exp – Average recommended step size returned for explicit solver.
_evaluate_integrator(integrator, h_min_lower_bound=1e-12, raise_errors=True, debug=True)

This function computes the average step size and the minimal step size that a given integration method from GSL uses to evolve a certain system of ODEs during a certain simulation time, integration method from GSL and spike train for a given maximal stepsize.

This function will reset the numpy random seed.

Parameters:
  • integrator – A method from the GSL library for evolving ODEs, e.g. odeiv.step_rk4.
  • h_min_lower_bound – The minimum acceptable step size. Integration will terminate with an error if this step size is reached.
  • raise_errors – Stop and raise exception when error occurs, or try to continue.
Return h_min:

Minimum recommended step size.
Return h_avg:

Average recommended step size.
Return runtime:

Wall clock runtime.
check_stiffness(raise_errors=False)

Perform stiffness testing: use implicit and explicit solvers to simulate the dynamical system, then decide which is the better solver to use.

For details, see https://ode-toolbox.readthedocs.io/en/latest/index.html#numeric-solver-selection-criteria

Returns:Either “implicit”, “explicit” or “warning”.
Return type:str
random_seed

odetoolbox.sympy_helpers module

class odetoolbox.sympy_helpers.SympyPrinter(settings=None)

Bases: sympy.printing.str.StrPrinter

_print_Exp1(expr)
_print_Function(expr)

Overrides base class method to print min() and max() functions in lowercase.

odetoolbox.sympy_helpers._custom_simplify_expr(expr: str)

Custom expression simplification

odetoolbox.sympy_helpers._find_in_matrix(A, el)
odetoolbox.sympy_helpers._is_sympy_type(var)
odetoolbox.sympy_helpers._is_zero(x)

Check if a sympy expression is equal to zero.

In the ideal case, we would like to use sympy.simplify() to do simplification of an expression before comparing it to zero. However, for expressions of moderate size (e.g. a few dozen terms involving exp() functions), it becomes unbearably slow. We therefore use this internal function, so that the simplification function can be easily switched over.

Tests by expand_mul only, suitable for polynomials and rational functions.

Ref.: https://github.com/sympy/sympy PR #13877 by @normalhuman et al. merged on Jan 27, 2018

odetoolbox.system_of_shapes module

exception odetoolbox.system_of_shapes.PropagatorGenerationException

Bases: Exception

Thrown in case an error occurs while generating propagators.

class odetoolbox.system_of_shapes.SystemOfShapes(x: sympy.matrices.dense.MutableDenseMatrix, A: sympy.matrices.dense.MutableDenseMatrix, b: sympy.matrices.dense.MutableDenseMatrix, c: sympy.matrices.dense.MutableDenseMatrix, shapes: List[odetoolbox.shapes.Shape])

Bases: object

Represent a dynamical system in the canonical form \(\mathbf{x}' = \mathbf{Ax} + \mathbf{b} + \mathbf{c}\).

__init__(x: sympy.matrices.dense.MutableDenseMatrix, A: sympy.matrices.dense.MutableDenseMatrix, b: sympy.matrices.dense.MutableDenseMatrix, c: sympy.matrices.dense.MutableDenseMatrix, shapes: List[odetoolbox.shapes.Shape])

Initialize a dynamical system in the canonical form \(\mathbf{x}' = \mathbf{Ax} + \mathbf{b} + \mathbf{c}\).

Parameters:
  • x – Vector containing variable symbols.
  • A – Matrix containing linear part.
  • b – Vector containing inhomogeneous part (constant term).
  • c – Vector containing nonlinear part.
classmethod from_shapes(shapes: List[odetoolbox.shapes.Shape], parameters=None)

Construct the global system matrix \(\mathbf{A}\) and inhomogeneous part (constant term) \(\mathbf{b}\) and nonlinear part \(\mathbf{c}\) on the basis of all shapes in shapes, and such that

\[\mathbf{x}' = \mathbf{Ax} + \mathbf{b} + \mathbf{c}\]
generate_numeric_solver(state_variables=None)

Generate the symbolic expressions for numeric integration state updates; return as JSON.

generate_propagator_solver()

Generate the propagator matrix and symbolic expressions for propagator-based updates; return as JSON.

get_connected_symbols(idx: int) → List[sympy.core.symbol.Symbol]

Extract all symbols belonging to a shape with symbol self.x_[idx] from the system matrix.

For example, if symbol i is x, and symbol j is y, and the system is:

\[\begin{split}\frac{dx}{dt} &= y\\ \frac{dy}{dt} &= y' = -\frac{1}{\tau^2} x - \frac{2}{\tau} y\end{split}\]

Then get_connected_symbols() for symbol x would return [x, y], and get_connected_symbols() for y would return the same.

get_dependency_edges()
get_initial_value(sym)
get_jacobian_matrix()

Get the Jacobian matrix as symbolic expressions. Entries in the matrix are sympy expressions.

If the dynamics of variables \(x_1, \ldots, x_N\) is defined as \(x_i' = f_i\), then row \(i\) of the Jacobian matrix \(\mathbf{J}_i = \left[\begin{matrix}\frac{\partial f_i}{\partial x_0} & \cdots & \frac{\partial f_i}{\partial x_N}\end{matrix}\right]\).

get_lin_cc_symbols(E, parameters=None)

Retrieve the variable symbols of those shapes that are linear and constant coefficient. In the case of a higher-order shape, will return all the variable symbols with "__d" suffixes up to the order of the shape.

get_shape_by_symbol(sym: str) → Optional[odetoolbox.shapes.Shape]
get_sub_system(symbols)

Return a new SystemOfShapes instance which discards all symbols and equations except for those in symbols. This is probably only sensible when the elements in symbols do not dependend on any of the other symbols that will be thrown away.

propagate_lin_cc_judgements(node_is_lin, E)

Propagate: if a node depends on a node that is not linear and constant coefficient, it cannot be linear and constant coefficient.

Parameters:
  • node_is_lin – Initial assumption about whether node is linear and constant coefficient.
  • E – List of edges returned from dependency analysis.
reconstitute_expr(state_variables=None)

Reconstitute a sympy expression from a system of shapes (which is internally encoded in the form \(\mathbf{x}' = \mathbf{Ax} + \mathbf{b} + \mathbf{c}\)).

Before returning, the expression is simplified using a custom series of steps, passed via the simplify_expression argument (see the ODE-toolbox documentation for more details).

shape_order_from_system_matrix(idx: int) → int

Determine shape differential order from system matrix of symbol self.x_[idx]

odetoolbox.system_of_shapes.off_diagonal_is_zero(row: int, A) → bool