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
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
Meiss. Differential Dynamical Systems (Monographs on Mathematical Modeling and Computation), Book 14, SIAM publishers, 2007.