autokoopman.benchmark package#

Submodules#

autokoopman.benchmark.bio2 module#

class autokoopman.benchmark.bio2.Bio2#

Bases: SymbolicContinuousSystem

Spring pendulum

A nine-dimensional continuous model which is adapted from a biological system given in [1].

Initial set:

x1 in [0.99,1.01] x2 in [0.99,1.01] x3 in [0.99,1.01] x4 in [0.99,1.01] x5 in [0.99,1.01] x6 in [0.99,1.01] x7 in [0.99,1.01] x8 in [0.99,1.01] x9 in [0.99,1.01]

References
  1. Klipp, R. Herwig, A. Kowald, C. Wierling, H. Lehrach. Systems Biology in Practice: Concepts, Implementation and Application. Wiley-Blackwell, 2005.

autokoopman.benchmark.fhn module#

class autokoopman.benchmark.fhn.FitzHughNagumo(i_ext=0.0, r=1.0, a=0.3, b=0.3, tau=3.0)#

Bases: SymbolicContinuousSystem

FitzHugh-Nagumo System (FHN)

The FHN Model is an example of a relaxation oscillator.

\[\begin{split}\left\{\begin{array}{l} \dot{x}_{1}=x_0 - \frac{x_0^3}{3} - x_1 + R I_{ext} \\ \dot{x}_{2}=\frac{x_0 + a - b x_1}{\tau} \\ \end{array}\right.\end{split}\]
Parameters:
  • i_ext – external stimulus

  • r – resistance

  • a – model parameter \(a\)

  • b – model parameter \(b\)

  • tau – model parameter \(\tau\)

Setting \(a=b=0\) creates the Van der Pol Oscillator.

Example:
import autokoopman.benchmark.fhn as fhn

sys = fhn.FitzHughNagumo()
traj = sys.solve_ivp(
    initial_state=[2.0, 1.5],
    tspan=(0.0, 20.0),
    sampling_period=0.01
)
traj.states
# array([[ 2.        ,  1.5       ],
#        [ 1.97862133,  1.50612782],
#        [ 1.95779684,  1.51217922],
#        ...,
#        [-1.32849614, -0.81944446],
#        [-1.32576847, -0.82204749],
#        [-1.32303553, -0.82463882]])
References

FitzHugh, R. (1961). Impulses and physiological states in theoretical models of nerve membrane. Biophysical journal, 1(6), 445-466.

autokoopman.benchmark.lalo20 module#

class autokoopman.benchmark.lalo20.LaubLoomis#

Bases: SymbolicContinuousSystem

Laub-Loomis Benchmark (LALO20)

They are boxes centered at \(x_1(0) = 1.2\), \(x_2(0) = 1.05\), \(x_3(0) = 1.5\), \(x_4(0) = 2.4\), \(x_5(0) = 1\), \(x_6(0) = 0.1\), \(x_7(0) = 0.45\). The range of the box in the \(i\)-th dimension is defined by the interval \([x_i(0) − W, x_i(0) + W ]\).

We consider \(W = 0.01\), \(W = 0.05\), and \(W =0.1\). For \(W =0.01\) and \(W =0.05\) we consider the unsafe region defined by \(x4 \ge 4.5\), while for \(W = 0.1\), the unsafe set is defined by \(x4 \ge 5\). The time horizon for all cases is \([0, 20]\).

The dynamics are defined as

\[\begin{split}\left\{\begin{array}{l} \dot{x}_{1}=1.4 x_{3}-0.9 x_{1} \\ \dot{x}_{2}=2.5 x_{5}-1.5 x_{2} \\ \dot{x}_{3}=0.6 x_{7}-0.8 x_{2} x_{3} \\ \dot{x}_{4}=2-1.3 x_{3} x_{4} \\ \dot{x}_{5}=0.7 x_{1}-x_{4} x_{5} \\ \dot{x}_{6}=0.3 x_{1}-3.1 x_{6} \\ \dot{x}_{7}=1.8 x_{6}-1.5 x_{2} x_{7} \end{array}\right.\end{split}\]

The system is asymptotically stable and the equilibrium is the origin.

Example:
import autokoopman.benchmark.lalo20 as lalo

sys = lalo.LaubLoomis()
traj = sys.solve_ivp(
    initial_state=[1.2, 1.05, 1.5, 2.4, 1.0, 0.1, 0.45],
    tspan=(0.0, 20.0),
    sampling_period=0.1
)
traj.states
# array([[1.2       , 1.05      , 1.5       , ..., 1.        , 0.1       ,
#    0.45      ],
#   [1.29071503, 1.11992937, 1.39927254, ..., 0.87447922, 0.10557573,
#    0.39928611],
#   [1.3600053 , 1.15715223, 1.29895902, ..., 0.79376949, 0.11169569,
#    0.35464381],
#   ...,
#   [0.89473053, 0.36803281, 0.58535322, ..., 0.22894195, 0.08595487,
#    0.28538999],
#   [0.89608814, 0.37001616, 0.58517889, ..., 0.22991276, 0.08609767,
#    0.28511588],
#   [0.89729548, 0.37195084, 0.58490685, ..., 0.23082109, 0.0862578 ,
#    0.28478635]])
References

Geretti, L., Sandretto, J. A. D., Althoff, M., Benet, L., Chapoutot, A., Chen, X., … & Schilling, C. (2020). ARCH-COMP20 category report: Continuous and hybrid systems with nonlinear dynamics. EPiC Series in Computing, 74, 49-75. pp 59

autokoopman.benchmark.pendulum module#

class autokoopman.benchmark.pendulum.PendulumWithInput(g=9.81, l=1.0, beta=0.0)#

Bases: SymbolicContinuousSystem

Simple Pendulum with Constant Torque Input

We model a pendulum with the equation:

\[l^2 \ddot{\theta} + \beta \theta + g l \sin \theta = \tau\]

This leads to the state space form:

\[\begin{split}\left\{\begin{array}{l} \dot{x}_{1} = x_2 \\ \dot{x}_{2} = - g / l \sin x_1 - 2 \beta x_2 + u_1 \\ \end{array}\right.\end{split}\]

Note that \(\beta=b / m\) kn this formulation: http://underactuated.mit.edu/pend.html .

autokoopman.benchmark.prde20 module#

class autokoopman.benchmark.prde20.ProdDestr(alpha=0.3)#

Bases: SymbolicContinuousSystem

Production Destruction System (PRDE20)

In this model, \(x_0\) is the nutrients, \(x_1\) the phytoplankton and \(x_2\) the detritus. The constraints are \(x_0(t), x_1(t), x_2(t)\) are positive and \(x_0(t) + x_1(t) + x_2(t) = 10\) for all t.

The initial condition \(x_0(0) = 9.98\), \(x_1(0) = 0.01\) and \(x_2(0) = 0.01\)

\(I: x_0(0) \in [9.5, 10.0]\), i.e., uncertainty on the initial condition;

From Michaelis-Menten Theory, the evolution function is

\[\begin{split}\begin{bmatrix} \dot x_0 \\ \dot x_1 \\ \dot x_2 \end{bmatrix} = \begin{bmatrix} \frac{-x_0 x_1}{1 + x_0} \\ \frac{x_0 x_1}{1 + x_0} - 0.3 x_1 \\ 0.3 x_1 \end{bmatrix}\end{split}\]
Parameters:

alpha – model parameter

Example:
import autokoopman.benchmark.prde20 as prde

sys = prde.ProdDestr()
traj = sys.solve_ivp(
    initial_state=[9.98, 0.01, 0.01],
    tspan=(0.0, 20.0),
    sampling_period=0.01
)
traj.states
# array([[9.98000000e+00, 1.00000000e-02, 1.00000000e-02],
#        [9.97990883e+00, 1.00610783e-02, 1.00300915e-02],
#        [9.97981710e+00, 1.01225295e-02, 1.00603668e-02],
#        ...,
#        [3.34608394e-08, 4.41976926e-01, 9.55802304e+00],
#        [3.33135031e-08, 4.40652987e-01, 9.55934698e+00],
#        [3.31672496e-08, 4.39333011e-01, 9.56066696e+00]])
References

Geretti, L., Sandretto, J. A. D., Althoff, M., Benet, L., Chapoutot, A., Chen, X., … & Schilling, C. (2020). ARCH-COMP20 category report: Continuous and hybrid systems with nonlinear dynamics. EPiC Series in Computing, 74, 49-75. pp 53

autokoopman.benchmark.robe21 module#

class autokoopman.benchmark.robe21.RobBench(alpha=0.4, beta=100, gamma=1000)#

Bases: SymbolicContinuousSystem

Robertson chemical reaction benchmark (ROBE21)

As proposed by Robertson [31], this chemical reaction system models the kinetics of an auto-catalytic reaction. x, y and z are the (positive) concentrations of the species, with the assumption that x+y+z = 1. Here alpha is a small constant, while beta and gamma take on large values. The initial condition is always \(x(0) = 1\), \(y(0) = 0\) and \(z(0) = 0\).

\(I: x_0(0) \in [9.5, 10.0]\), i.e., uncertainty on the initial condition;

In this benchmark we fix alpha = 0.4 and analyze the system under three different pairs of values for beta and gamma: 1. beta = 10^2 , gamma = 10^3 2. beta=10^3, gamma = 10^5 3. beta = 10^3, gamma = 10^7

We are interested in computing the reachable tube until t = 40, to see how the integration scheme holds under the stiff behavior. No verification objective is enforced.

For each of the three setups, the following three measures are required: 1. the execution time for evolution; 2. the number of integration steps taken; 3. the width of the sum of the concentrations s = x + y + z at the final time. Additionally, a figure with s (in the [0.999, 1.001] range) w.r.t. time overlaid for the three setups should be shown.

As proposed by Robertson, this chemical reaction system models the kinetics of an auto- catalytic reaction.

References

Geretti, L., dit Sandretto, J. A., Althoff, M., Benet, L., Chapoutot, A., Collins, P., … & Wetzlinger, M. (2021). ARCH-COMP21 Category Report: Continuous and Hybrid Systems with Nonlinear Dynamics. EPiC Series in Computing, 80, 32-54. H. H. Robertson. The solution of a set of reaction rate equations. In ”Numerical analysis: an introduction”, page 178–182. Academic Press, 1966.

autokoopman.benchmark.spring module#

class autokoopman.benchmark.spring.Spring(g=9.81)#

Bases: SymbolicContinuousSystem

Spring pendulum

We study the behavior of the planar spring-pendulum described in [1]. It consists of a solid ball of mass \(m\) and a spring of natural length \(L\). The spring constant is \(k\).

We study the evolutions of the length r of the spring and the angle \(\theta\) between the spring and the vertical. They are modeled by the following differential equations.

The constants are set as \(k=2\), \(L=1\), and \(g=9.8\).

\[\begin{split}\begin{bmatrix} \dot r \\ \dot \theta \\ \dot d_r \\ \dot d_{\theta} \end{bmatrix} = \begin{bmatrix} d_r \\ d_{\theta} \\ r d_{\theta}^2 + g * \cos (\theta) -2 (r - 1) \\ \left( -2 d_r d_{\theta} - g \sin (\theta) \right) / 2 \end{bmatrix}\end{split}\]
Initial set:

r in [1.19,1.21] theta in [0.49,0.51] dr in [0,0] dtheta in [0,0]

References
    1. Meiss. Differential Dynamical Systems (Monographs on Mathematical Modeling and Computation), Book 14, SIAM publishers, 2007.

Module contents#