diff --git a/src/system_model.py b/src/system_model.py
new file mode 100644
index 0000000000000000000000000000000000000000..5fac65a9722da2b4c5c8ccb29cc7c7d7914795e4
--- /dev/null
+++ b/src/system_model.py
@@ -0,0 +1,400 @@
+## This file is part of the simulative evaluation for the qronos observer abstractions.
+## Copyright (C) 2022-2023 Tim Rheinfels <tim.rheinfels@fau.de>
+## See https://gitlab.cs.fau.de/qronos-state-abstractions/simulation
+##
+## Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
+##
+## 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
+##
+## 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
+##
+## 3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.
+##
+## THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS “AS IS” AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+
+###
+### @file system_model.py
+###
+### @brief Provides a class for the system model described in @cite ECRTS23-Rheinfels
+###
+### @author Tim Rheinfels <tim.rheinfels@fau.de>
+###
+
+import control
+import json
+import numpy as np
+import scipy as sp
+import scipy.linalg
+
+from ellipsoid import Ellipsoid
+from experiment.system import System
+from experiment.signal import Signal
+from util import Serializable
+
+###
+### @brief Encapsulates the switched linear system model described in @cite ECRTS23-Rheinfels, section 2.2
+###
+class SystemModel(System, Serializable):
+
+ ###
+ ### @brief Helper function for creating an @f$ n\ times m @f$ matrix from different input data
+ ###
+ ### @param A Input data. Either a matrix @f$ A \in \R^{n \times m} @f$, a scalar, or None
+ ### @param n Number of rows
+ ### @param m Number of columns
+ ###
+ ### @returns Either
+ ### - The matrix @f$ A \in \R^{n \times m} @f$ if @p A is a matrix,
+ ### - a scalar @f$ A \in \R @f$ encapsulated as an @f$ n \times m @f$ matrix if @p A is a sclar, or
+ ### - the zero matrix @f$ 0 \in \R^{n \times m} @f$ if @p A is None
+ ###
+ @staticmethod
+ def _promote_to_matrix(A, n, m):
+ if A is None:
+ return np.zeros((n, m))
+ elif isinstance(A, float) or isinstance(A, int):
+ return np.array([[float(A)]])
+ return A
+
+ ###
+ ### @brief Encapsualtes on of the switched system's modes
+ ###
+ class Mode(Serializable):
+
+ ###
+ ### @brief Constructor
+ ###
+ ### @param name Name to assign to the mode
+ ### @param A Dynamics matrix @f$ A \in \R^{n_x \times n_x} @f$ or None to omit
+ ### @param G Process disturbance input matrix @f$ G \in \R^{n_x \times n_d} @f$ or None to omit
+ ### @param E_D Process disturbance bounding @ref ellipsoid.Ellipsoid
+ ### @param C Measurement matrix @f$ C \in \R^{n_y \times n_x} @f$ or None to omit
+ ### @param H Measurement disturbance input matrix @f$ H \in \R^{n_y \times n_z} @f$ or None to omit
+ ### @param E_Z Measurement disturbance bounding @ref ellipsoid.Ellipsoid
+ ###
+ def __init__(self, name, A, G, E_D, C, H, E_Z):
+ assert(isinstance(name, str))
+
+ assert(A is not None)
+ A = SystemModel._promote_to_matrix(A, 0, 0)
+ n_x = A.shape[0]
+
+ G = SystemModel._promote_to_matrix(G, n_x, 0)
+ n_d = G.shape[1]
+
+ if E_D is None:
+ E_D = Ellipsoid()
+ assert(isinstance(E_D, Ellipsoid))
+ assert(E_D.n == n_d)
+
+ C = SystemModel._promote_to_matrix(C, 0, n_x)
+ assert(C.shape[1] == n_x)
+ n_y = C.shape[0]
+
+ H = SystemModel._promote_to_matrix(H, n_y, 0)
+ assert(H.shape[0] == n_y)
+ n_z = H.shape[1]
+
+ if E_Z is None:
+ E_Z = Ellipsoid()
+ assert(isinstance(E_Z, Ellipsoid))
+ assert(E_Z.n == n_z)
+
+ self._name = name
+ self._A = A
+ self._G = G
+ self._E_D = E_D
+ self._C = C
+ self._H = H
+ self._E_Z = E_Z
+ self._n_x = n_x
+ self._n_d = n_d
+ self._n_y = n_y
+ self._n_z = n_z
+
+ ###
+ ### @brief Getter for the mode's name
+ ###
+ ### @returns Mode's name
+ ###
+ @property
+ def name(self):
+ return self._name
+
+ ###
+ ### @brief Getter for the mode's dynamic matrix @f$ A \in \R^{n_x \times n_x} @f$
+ ###
+ ### @returns Mode's dynamic matrix @f$ A \in \R^{n_x \times n_x} @f$
+ ###
+ @property
+ def A(self):
+ return self._A
+
+ ###
+ ### @brief Getter for the mode's process disturbance input matrix @f$ G \in \R^{n_x \times n_d} @f$
+ ###
+ ### @returns Mode's process disturbance input matrix @f$ G \in \R^{n_x \times n_d} @f$
+ ###
+ @property
+ def G(self):
+ return self._G
+
+ ###
+ ### @brief Getter for the mode's process disturbance bounding @ref ellipsoid.Ellipsoid
+ ###
+ ### @returns Mode's process disturbance bounding @ref ellipsoid.Ellipsoid
+ ###
+ @property
+ def E_D(self):
+ return self._E_D
+
+ ###
+ ### @brief Getter for the mode's measurement matrix @f$ C \in \R^{n_y \times n_x} @f$
+ ###
+ ### @returns Mode's measurement matrix @f$ C \in \R^{n_y \times n_x} @f$
+ ###
+ @property
+ def C(self):
+ return self._C
+
+ ###
+ ### @brief Getter for the mode's measurement disturbance input matrix @f$ H \in \R^{n_y \times n_z} @f$
+ ###
+ ### @returns Mode's measurement disturbance input matrix @f$ H \in \R^{n_y \times n_z} @f$
+ ###
+ @property
+ def H(self):
+ return self._H
+
+ ###
+ ### @brief Getter for the mode's measurement disturbance bounding @ref ellipsoid.Ellipsoid
+ ###
+ ### @returns Mode's measurement disturbance bounding @ref ellipsoid.Ellipsoid
+ ###
+ @property
+ def E_Z(self):
+ return self._E_Z
+
+ ###
+ ### @brief Getter for the mode's state dimension @f$ n_x @f$
+ ###
+ ### @returns Mode's state dimension @f$ n_x @f$
+ ###
+ @property
+ def n_x(self):
+ return self._n_x
+
+ ###
+ ### @brief Getter for the mode's process disturbance dimension @f$ n_d @f$
+ ###
+ ### @returns Mode's process disturbance dimension @f$ n_d @f$
+ ###
+ @property
+ def n_d(self):
+ return self._n_d
+
+ ###
+ ### @brief Getter for the mode's measurement dimension @f$ n_y @f$
+ ###
+ ### @returns Mode's measurement dimension @f$ n_y @f$
+ ###
+ @property
+ def n_y(self):
+ return self._n_y
+
+ ###
+ ### @brief Getter for the mode's measurement disturbance dimension @f$ n_z @f$
+ ###
+ ### @returns Mode's measurement disturbance dimension @f$ n_z @f$
+ ###
+ @property
+ def n_z(self):
+ return self._n_z
+
+ ###
+ ### @brief Converts the mode's data into a dictionary
+ ###
+ ### @returns Dictionary characterizing the mode
+ ###
+ ### @see util.Serializable
+ ###
+ def to_dict(self):
+ return {
+ 'name': self.name,
+ 'A': self.A,
+ 'G': self.G,
+ 'E_D': self.E_D.to_dict(),
+ 'C': self.C,
+ 'H': self.H,
+ 'E_Z': self.E_Z.to_dict(),
+ 'n_x': self.n_x,
+ 'n_d': self.n_d,
+ 'n_y': self.n_y,
+ 'n_z': self.n_z,
+ }
+
+ ###
+ ### @brief Constructor
+ ###
+ ### @param name Name to assign to the system model
+ ### @param modes List of @ref system_model.SystemModel.Mode instances
+ ### @param E_X_0 An @ref ellipsoid.Ellipsoid bounding the initial state
+ ### @param C_s Safety output matrix @f$ C_s \in \R^{n_s \times n_x} @f$ or None to use identity
+ ### @param E_S Safety specification @ref ellipsoid.Ellipsoid
+ ###
+ def __init__(self, name, modes, E_X_0, C_s, E_S):
+ assert(isinstance(name, str))
+ assert(len(modes) > 0)
+ n_x = modes[0].n_x
+ n_d = 0
+ n_y = 0
+ n_z = 0
+ for mode in modes:
+ assert(isinstance(mode, SystemModel.Mode))
+ assert(mode.n_x == n_x)
+ n_d = max((n_d, mode.n_d))
+ n_y = max((n_y, mode.n_y))
+ n_z = max((n_z, mode.n_z))
+
+ assert(isinstance(E_X_0, Ellipsoid))
+ assert(E_X_0.n == n_x)
+ assert(isinstance(E_S, Ellipsoid))
+ assert(not E_S.degenerate)
+ n_s = E_S.n
+ if C_s is None:
+ C_s = np.eye(n_s)
+ elif isinstance(C_s, float) or isinstance(C_s, int):
+ C_s = np.array([[C_s]])
+ assert(isinstance(C_s, np.ndarray))
+ assert((C_s.ndim == 2) and (C_s.shape[0] == n_s) and (C_s.shape[1] == n_x))
+
+ self._name = name
+ self._modes = modes
+ self._E_X_0 = E_X_0
+ self._C_s = C_s
+ self._E_S = E_S
+ self._n_x = n_x
+ self._n_d = n_d
+ self._n_y = n_y
+ self._n_z = n_z
+ self._n_s = n_s
+
+ ###
+ ### @brief Getter for the system model's name
+ ###
+ ### @returns System model's name
+ ###
+ @property
+ def name(self):
+ return self._name
+
+ ###
+ ### @brief Getter for the system model's modes
+ ###
+ ### @returns System model's modes
+ ###
+ @property
+ def modes(self):
+ return self._modes
+
+ ###
+ ### @brief Getter for the system model's initial state @ref ellipsoid.Ellipsoid
+ ###
+ ### @returns System model's initial state @ref ellipsoid.Ellipsoid
+ ###
+ @property
+ def E_X_0(self):
+ return self._E_X_0
+
+ ###
+ ### @brief Getter for the system model's safety output matrix @f$ C_s \in \R^{n_s \times n_x} @f$
+ ###
+ ### @returns System model's safety output matrix @f$ C_s \in \R^{n_s \times n_x} @f$
+ ###
+ @property
+ def C_s(self):
+ return self._C_s
+
+ ###
+ ### @brief Getter for the system model's safety specification @ref ellipsoid.Ellipsoid
+ ###
+ ### @returns System model's safety specification @ref ellipsoid.Ellipsoid
+ ###
+ @property
+ def E_S(self):
+ return self._E_S
+
+ ###
+ ### @brief Getter for the system model's state dimension @f$ n_x @f$
+ ###
+ ### @returns System model's state dimension @f$ n_x @f$
+ ###
+ @property
+ def n_x(self):
+ return self._n_x
+
+ ###
+ ### @brief Getter for the system model's maximum process disturbance dimension @f$ n_d @f$
+ ###
+ ### @returns System model's maximum process disturbance dimension @f$ n_d @f$
+ ###
+ @property
+ def n_d(self):
+ return self._n_d
+
+ ###
+ ### @brief Getter for the system model's maximum measurement dimension @f$ n_y @f$
+ ###
+ ### @returns System model's maximum measurement dimension @f$ n_y @f$
+ ###
+ @property
+ def n_y(self):
+ return self._n_y
+
+ ###
+ ### @brief Getter for the system model's maximum measurement disturbance dimension @f$ n_z @f$
+ ###
+ ### @returns System model's maximum measurement disturbance dimension @f$ n_z @f$
+ ###
+ @property
+ def n_z(self):
+ return self._n_z
+
+ ###
+ ### @brief Getter for the system model's safety output dimension @f$ n_sw @f$
+ ###
+ ### @returns System model's safety output dimension @f$ n_sw @f$
+ ###
+ @property
+ def n_s(self):
+ return self._n_s
+
+ ###
+ ### @brief Getter for the number of modes @f$ n_\Sigma @f$
+ ###
+ ### @returns Number of modes @f$ n_\Sigma @f$
+ ###
+ @property
+ def n_Sigma(self):
+ return len(self._modes)
+
+ ###
+ ### @brief Converts the system model's data into a dictionary
+ ###
+ ### @returns Dictionary characterizing the system model
+ ###
+ ### @see util.Serializable
+ ###
+ def to_dict(self):
+ return {
+ 'name': self.name,
+ 'modes': [mode.to_dict() for mode in self.modes],
+ 'E_X_0': self.E_X_0.to_dict(),
+ 'C_s': self.C_s,
+ 'E_S': self.E_S.to_dict(),
+ 'n_x': self.n_x,
+ 'n_d': self.n_d,
+ 'n_y': self.n_y,
+ 'n_z': self.n_z,
+ 'n_s': self.n_s,
+ }
diff --git a/src/test/system_model/__init__.py b/src/test/system_model/__init__.py
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/src/test/system_model/constructor.py b/src/test/system_model/constructor.py
new file mode 100644
index 0000000000000000000000000000000000000000..991c93e157629b36b9798d362921749dffa2a048
--- /dev/null
+++ b/src/test/system_model/constructor.py
@@ -0,0 +1,211 @@
+## This file is part of the simulative evaluation for the qronos observer abstractions.
+## Copyright (C) 2022-2023 Tim Rheinfels <tim.rheinfels@fau.de>
+## See https://gitlab.cs.fau.de/qronos-state-abstractions/simulation
+##
+## Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
+##
+## 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
+##
+## 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
+##
+## 3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.
+##
+## THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS “AS IS” AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+
+###
+### @file test/system_model/constructor.py
+###
+### @brief Provides unit tests for @ref system_model.SystemModel.__init__
+###
+### @author Tim Rheinfels <tim.rheinfels@fau.de>
+###
+
+import numpy as np
+import unittest
+
+from ellipsoid import Ellipsoid
+from system_model import SystemModel
+
+###
+### @brief Encapsulates the unit tests for @ref system_model.SystemModel.__init__
+###
+class Test(unittest.TestCase):
+
+ ###
+ ### @brief Sets up the same rng for all tests
+ ###
+ def setUp(self):
+ self.rng = np.random.default_rng(0)
+
+ ###
+ ### @brief Draws a 1000 valid system models and checks if the constructor works properly for all of them
+ ###
+ def test_random(self):
+ for i in range(1000):
+ n_x = self.rng.integers(1, 11)
+ n_s = self.rng.integers(1, 11)
+ n_Sigma = self.rng.integers(1, 11)
+
+ def draw_random_matrix(rows, cols):
+ M = self.rng.uniform(-1, +1, size=(rows, cols))
+ if (rows == 1) and (cols == 1) and (self.rng.integers(0, 101) >= 50):
+ if self.rng.integers(0, 101) >= 50:
+ M = float(M)
+ else:
+ M = int(M)
+ return M
+
+ def draw_random_spd_matrix(n):
+ Pi = self.rng.uniform(-1, +1, size=(n, n))
+ Pi = (Pi + Pi.T) / 2.0 + n * np.eye(n)
+ return Pi
+
+ def check_matrix(M, M_ref):
+ self.assertTrue(isinstance(M, np.ndarray))
+ if M_ref is None:
+ self.assertEqual(M.size, 0)
+ elif isinstance(M_ref, float) or isinstance(M_ref, int):
+ self.assertEqual(M.size, 1)
+ self.assertEqual(M[0], M_ref)
+ else:
+ self.assertEqual(M.shape, M_ref.shape)
+ self.assertTrue(np.all(M == M_ref))
+
+ def check_ellipsoid(E, E_ref):
+ self.assertTrue(isinstance(E, Ellipsoid))
+ self.assertTrue(E.centered)
+ if E_ref is None:
+ self.assertEqual(E.n, 0)
+ else:
+ self.assertTrue(np.all(E.P == E_ref.P))
+
+ modes = []
+ n_d_max = 0
+ n_z_max = 0
+ n_y_max = 0
+ for j in range(n_Sigma):
+ n_d = self.rng.integers(0, 11)
+ n_z = self.rng.integers(0, 11)
+ n_y = self.rng.integers(0, 11)
+
+ E_D = Ellipsoid( draw_random_spd_matrix(n_d) )
+ E_Z = Ellipsoid( draw_random_spd_matrix(n_z) )
+
+ A = draw_random_matrix(n_x, n_x)
+
+ G = draw_random_matrix(n_x, n_d)
+ if self.rng.integers(0, 101) >= 80:
+ G = None
+ E_D = None
+ n_d = 0
+
+ C = draw_random_matrix(n_y, n_x)
+ if self.rng.integers(0, 101) >= 80:
+ C = None
+ n_y = 0
+
+ H = draw_random_matrix(n_y, n_z)
+ if self.rng.integers(0, 101) >= 80:
+ H = None
+ E_Z = None
+ n_z = 0
+
+ if n_d > n_d_max:
+ n_d_max = n_d
+
+ if n_z > n_z_max:
+ n_z_max = n_z
+
+ if n_y > n_y_max:
+ n_y_max = n_y
+
+ modes.append(
+ SystemModel.Mode(
+ 'test_mode_%d' % j,
+ A, G, E_D,
+ C, H, E_Z,
+ )
+ )
+
+ self.assertEqual(modes[-1].name, 'test_mode_%d' % j)
+
+ check_matrix(modes[-1].A, A)
+ check_matrix(modes[-1].G, G)
+ check_ellipsoid(modes[-1].E_D, E_D)
+
+ check_matrix(modes[-1].C, C)
+ check_matrix(modes[-1].H, H)
+ check_ellipsoid(modes[-1].E_Z, E_Z)
+
+ d = modes[-1].to_dict()
+
+ if E_D is None:
+ E_D = Ellipsoid()
+
+ if E_Z is None:
+ E_Z = Ellipsoid()
+
+ self.assertEqual(d['name'], 'test_mode_%d' % j)
+ self.assertTrue(np.all(d['A'] == A))
+ self.assertTrue(np.all(d['G'] == G))
+ self.assertTrue(np.all(d['E_D']['P'] == E_D.P))
+ self.assertTrue(np.all(d['E_D']['Pi'] == E_D.Pi))
+ self.assertTrue(np.all(d['E_D']['c'] == E_D.c))
+ self.assertTrue(np.all(d['C'] == C))
+ self.assertTrue(np.all(d['H'] == H))
+ self.assertTrue(np.all(d['E_Z']['P'] == E_Z.P))
+ self.assertTrue(np.all(d['E_Z']['Pi'] == E_Z.Pi))
+ self.assertTrue(np.all(d['E_Z']['c'] == E_Z.c))
+ self.assertEqual(d['n_x'], n_x)
+ self.assertEqual(d['n_d'], n_d)
+ self.assertEqual(d['n_y'], n_y)
+ self.assertEqual(d['n_z'], n_z)
+
+ C_s = draw_random_matrix(n_s, n_x)
+ if (n_s == n_x) and (self.rng.integers(0, 101) >= 80):
+ C_s = None
+
+ E_X_0 = Ellipsoid( draw_random_spd_matrix(n_x) )
+ E_S = Ellipsoid( draw_random_spd_matrix(n_s) )
+
+ system_model = SystemModel(
+ 'test_model',
+ modes,
+ E_X_0,
+ C_s,
+ E_S,
+ )
+
+ self.assertEqual(system_model.name, 'test_model')
+ self.assertEqual(system_model.modes, modes)
+ check_ellipsoid(system_model.E_X_0, E_X_0)
+ if C_s is not None:
+ check_matrix(system_model.C_s, C_s)
+ else:
+ self.assertTrue(np.all(system_model.C_s == np.eye(n_x)))
+ check_ellipsoid(system_model.E_S, E_S)
+ self.assertEqual(system_model.n_x, n_x)
+ self.assertEqual(system_model.n_d, n_d_max)
+ self.assertEqual(system_model.n_z, n_z_max)
+ self.assertEqual(system_model.n_y, n_y_max)
+ self.assertEqual(system_model.n_s, n_s)
+ self.assertEqual(system_model.n_Sigma, n_Sigma)
+
+ d = system_model.to_dict()
+
+ if C_s is None:
+ C_s = np.eye(n_x)
+
+ self.assertEqual(d['name'], 'test_model')
+ self.assertTrue(isinstance(d['modes'], list))
+ self.assertTrue(np.all(d['E_X_0']['P'] == E_X_0.P))
+ self.assertTrue(np.all(d['E_X_0']['Pi'] == E_X_0.Pi))
+ self.assertTrue(np.all(d['E_X_0']['c'] == E_X_0.c))
+ self.assertTrue(np.all(d['C_s'] == C_s))
+ self.assertTrue(np.all(d['E_S']['P'] == E_S.P))
+ self.assertTrue(np.all(d['E_S']['Pi'] == E_S.Pi))
+ self.assertTrue(np.all(d['E_S']['c'] == E_S.c))
+ self.assertEqual(d['n_x'], n_x)
+ self.assertEqual(d['n_d'], n_d_max)
+ self.assertEqual(d['n_y'], n_y_max)
+ self.assertEqual(d['n_s'], n_s)