diff --git a/Makefile b/Makefile index d73e3a4eeae6873a2632eb2e2af59fdfe57ae525..9d482efec0516a41f36d9edd8acdfe69b2bfc794 100644 --- a/Makefile +++ b/Makefile @@ -127,6 +127,7 @@ CSRC = $(ALLCSRC) \ $(SRCDIR)/benchmarks/memory.c \ $(SRCDIR)/benchmarks/blind_abstraction.c \ $(SRCDIR)/benchmarks/observer_abstraction.c \ + $(SRCDIR)/benchmarks/luenberger_observer.c \ $(SRCDIR)/main.c \ $(SRCDIR)/random.c \ $(SRCDIR)/serial.c \ diff --git a/doc/references.bib b/doc/references.bib index 5d802b0723be9febb033010f237709b2d8e193c4..62b5c0d4f21eb9849be06931c3f8d3983308a129 100644 --- a/doc/references.bib +++ b/doc/references.bib @@ -7,3 +7,14 @@ pages = {t.b.d.}, notes = {Accepted}, } + +@article{CSS21-Dkhil, + author = {Dkhil, Monia and Dinh, Thach Ngoc and Wang, Zhenhua and Raïssi, Tarek and Amairi, Messaoud}, + journal = {IEEE Control Systems Letters}, + title = {{Interval Estimation for Discrete-Time Switched Linear Systems Based on $L_\infty$ Observer and Ellipsoid Analysis}}, + year = {2021}, + number = {1}, + pages = {13--18}, + volume = {5}, + doi = {10/gq5wdm}, +} diff --git a/scripts/benchmarks/luenberger_observer.py b/scripts/benchmarks/luenberger_observer.py new file mode 100644 index 0000000000000000000000000000000000000000..4a29bfca7b7b83349d0f8125588e7c9350267ae2 --- /dev/null +++ b/scripts/benchmarks/luenberger_observer.py @@ -0,0 +1,93 @@ +## This file is part of the execution-time 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/execution-time +## +## This program is free software: you can redistribute it and/or modify +## it under the terms of the GNU General Public License as published by +## the Free Software Foundation, either version 3 of the License, or +## (at your option) any later version. +## +## This program is distributed in the hope that it will be useful, +## but WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +## GNU General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with this program. If not, see <https://www.gnu.org/licenses/>. + +### +### @file benchmarks/luenberger_observer.py +### +### @brief Provides @p luenberger_observer.LuenbergerObserver, the @p benchmark.Benchmark +### subclass for the Luenberger observer benchmark +### +### @author Tim Rheinfels <tim.rheinfels@fau.de> +### + +import numpy as np +import json + +from benchmark import Benchmark +from data import decode + +### +### @brief Encapsulates the Luenberger observer execution time benchmark +### +class LuenbergerObserver(Benchmark): + + ### + ### @brief Validates the computations performed on the STM32 by running them + ### again in numpy and comparing the results. + ### + ### @see benchmark.Benchmark._validate for parameters + ### + def _validate(self, benchmark, key): + n_y = int(key) + + for run in benchmark['tests']: + ALC = decode(run['setup']['ALC']) + L = decode(run['setup']['L']) + M = decode(run['setup']['M']) + x = decode(run['setup']['x']) + X = decode(run['setup']['X']) + y = decode(run['setup']['y']) + + if n_y > 0: + assert(ALC.shape == (n_y, n_y)) + assert(L.shape == (n_y, n_y)) + assert(M.shape == (n_y, n_y)) + assert(np.allclose(M, M.T) and np.all(np.linalg.eigvalsh(M) > 0.0)) + assert(X.shape == (n_y, n_y)) + assert(np.allclose(X, X.T) and np.all(np.linalg.eigvalsh(X) > 0.0)) + else: + assert(ALC.size == 0) + assert(L.size == 0) + assert(M.size == 0) + assert(X.size == 0) + ALC = np.zeros((0, 0)) + L = np.zeros((0, 0)) + M = np.zeros((0, 0)) + X = np.zeros((0, 0)) + + assert(x.shape == (n_y,)) + assert(y.shape == (n_y,)) + + for i in range(benchmark['repetition_count']): + x = ALC @ x + L @ y + X = ALC @ X @ ALC.T + delta = np.sqrt( np.trace(X) / np.trace(M) ) + X = (1.0 + 1.0 / delta) * X + (1.0 + delta) * M + + x_p = decode(run['teardown']['x_p']) + X_p = decode(run['teardown']['X_p']) + + assert(np.allclose(x, x_p)) + assert(np.allclose(X, X_p)) + + ### + ### @brief Constructor + ### + ### @see benchmark.Benchmark.__init__ for parameters + ### + def __init__(self, data): + Benchmark.__init__(self, data, 'Luenberger Observer', r'luenberger_observer_([0-9]+)', lambda m: m.group(1)) diff --git a/src/benchmarks.c b/src/benchmarks.c index 9fa7d76dfcbd616ac08964408977fac0b90181cc..47f43ef94484d7a9c5c58a603b6e64572370b6a9 100644 --- a/src/benchmarks.c +++ b/src/benchmarks.c @@ -29,6 +29,7 @@ #include <benchmarks/blind_abstraction.h> #include <benchmarks/observer_abstraction.h> +#include <benchmarks/luenberger_observer.h> void benchmarks_run(void) { @@ -43,5 +44,6 @@ void benchmarks_run(void) // Comparative execution time benchmarks benchmarks_run_blind_abstraction(); benchmarks_run_observer_abstraction(); + benchmarks_run_luenberger_observer(); } diff --git a/src/benchmarks/luenberger_observer.c b/src/benchmarks/luenberger_observer.c new file mode 100644 index 0000000000000000000000000000000000000000..891b5d6fbf3038cb39379ce8677678e2a3a81d10 --- /dev/null +++ b/src/benchmarks/luenberger_observer.c @@ -0,0 +1,175 @@ +// This file is part of the execution-time 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/execution-time +// +// This program is free software: you can redistribute it and/or modify +// it under the terms of the GNU General Public License as published by +// the Free Software Foundation, either version 3 of the License, or +// (at your option) any later version. +// +// This program is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +// GNU General Public License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this program. If not, see <https://www.gnu.org/licenses/>. + +/// +/// @file benchmarks/luenberger_observer.c +/// +/// @brief Provides the implementation of the Luenberger observer execution time benchmark +/// +/// @author Tim Rheinfels <tim.rheinfels@fau.de> +/// + +#include "luenberger_observer.h" + +#include <stdlib.h> +#include <string.h> + +#include <linalg.h> +#include <random.h> +#include <serial.h> +#include <benchmarks/campaign.h> + +static float *ALC; ///< Observer dynamics matrix +static float *L; ///< Observer gain matrix +static float *M; ///< Disturbance ellipsoid's shape matrix +static float *x; ///< Ellipsoid's center vector +static float *x_p; ///< Memory for transforming the ellipsoid's center vector +static float *X; ///< Ellipsoid's shape matrix +static float *X_p; ///< Memory for transforming the ellipsoid's shape matrix +static float *y; ///< Measurement vector +static float *Ly; ///< Correction introduced by the measurement +static float trace_M; ///< Trace of the disturbance ellipsoid's shape matrix + +/// +/// @brief Configures the memory layout and randomly initializes the data +/// +static void setup_matrices(size_t n_y) +{ + ALC = benchmarks_memory; + L = ALC + n_y*n_y; + M = L + n_y*n_y; + x = M + n_y*n_y; + X = x + n_y; + y = X + n_y*n_y; + + x_p = y + n_y; + X_p = x_p + n_y; + Ly = X_p + n_y*n_y; + + random_matrix(ALC, n_y, n_y); + random_matrix(L, n_y, n_y); + random_matrix_spd(M, n_y); + random_matrix_spd(X, n_y); + random_matrix(x, n_y, 1u); + random_matrix(y, n_y, 1u); + memset(x_p, 0x00u, n_y * sizeof(float)); + memset(X_p, 0x00u, n_y * n_y * sizeof(float)); + memset(Ly, 0x00u, n_y * sizeof(float)); + + LINALG_TRACE(M, trace_M, n_y); + + serial_print_matrix("ALC", ALC, n_y, n_y); + serial_print_matrix("L", L, n_y, n_y); + serial_print_matrix("M", M, n_y, n_y); + serial_print_vector("x", x, n_y); + serial_print_matrix("X", X, n_y, n_y); + serial_print_vector("y", y, n_y); +} + +void benchmarks_run_luenberger_observer(void) +{ + #define SETUP(N_Y) { \ + setup_matrices(N_Y); \ + } + + #define TEARDOWN(N_Y) { \ + serial_print_vector("x_p", x, N_Y); \ + serial_print_matrix("X_p", X, N_Y, N_Y); \ + } + + #define CAMPAIGN(N_Y) { \ + BENCHMARKS_CAMPAIGN("luenberger_observer_" #N_Y, \ + { \ + LINALG_MV(ALC, x, x_p, 1.0f, 0.0f, N_Y, N_Y); \ + LINALG_MV(L, y, x_p, 1.0f, 1.0f, N_Y, N_Y); \ + \ + LINALG_MM(ALC, X, X_p, 1.0f, 0.0f, N_Y, N_Y, N_Y); \ + LINALG_MM_T(X_p, ALC, X, 1.0f, 0.0f, N_Y, N_Y, N_Y); \ + float trace_X; \ + LINALG_TRACE(X, trace_X, N_Y); \ + float delta = sqrtf(trace_X / trace_M); \ + LINALG_MA(X, M, X, 1.0f + 1.0f/delta, 1.0f + delta, 0.0f, N_Y, N_Y); \ + }, \ + SETUP(N_Y), \ + TEARDOWN(N_Y), \ + 10u, \ + 10u \ + ) \ + } + + CAMPAIGN(1); + CAMPAIGN(2); + CAMPAIGN(3); + CAMPAIGN(4); + CAMPAIGN(5); + CAMPAIGN(6); + CAMPAIGN(7); + CAMPAIGN(8); + CAMPAIGN(9); + CAMPAIGN(10); + CAMPAIGN(11); + CAMPAIGN(12); + CAMPAIGN(13); + CAMPAIGN(14); + CAMPAIGN(15); + CAMPAIGN(16); + CAMPAIGN(17); + CAMPAIGN(18); + CAMPAIGN(19); + CAMPAIGN(20); + CAMPAIGN(21); + CAMPAIGN(22); + CAMPAIGN(23); + CAMPAIGN(24); + CAMPAIGN(25); + CAMPAIGN(26); + CAMPAIGN(27); + CAMPAIGN(28); + CAMPAIGN(29); + CAMPAIGN(30); + CAMPAIGN(31); + CAMPAIGN(32); + CAMPAIGN(33); + CAMPAIGN(34); + CAMPAIGN(35); + CAMPAIGN(36); + CAMPAIGN(37); + CAMPAIGN(38); + CAMPAIGN(39); + CAMPAIGN(40); + CAMPAIGN(41); + CAMPAIGN(42); + CAMPAIGN(43); + CAMPAIGN(44); + CAMPAIGN(45); + CAMPAIGN(46); + CAMPAIGN(47); + CAMPAIGN(48); + CAMPAIGN(49); + CAMPAIGN(50); + CAMPAIGN(51); + CAMPAIGN(52); + CAMPAIGN(53); + CAMPAIGN(54); + CAMPAIGN(55); + CAMPAIGN(56); + CAMPAIGN(57); + CAMPAIGN(58); + CAMPAIGN(59); + CAMPAIGN(60); + +} diff --git a/src/benchmarks/luenberger_observer.h b/src/benchmarks/luenberger_observer.h new file mode 100644 index 0000000000000000000000000000000000000000..3def90d30074f978203c6b7937df872237aaff52 --- /dev/null +++ b/src/benchmarks/luenberger_observer.h @@ -0,0 +1,44 @@ +// This file is part of the execution-time 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/execution-time +// +// This program is free software: you can redistribute it and/or modify +// it under the terms of the GNU General Public License as published by +// the Free Software Foundation, either version 3 of the License, or +// (at your option) any later version. +// +// This program is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +// GNU General Public License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this program. If not, see <https://www.gnu.org/licenses/>. + +/// +/// @file benchmarks/luenberger_observer.h +/// +/// @brief Provides the interface for the Luenberger observer execution time benchmark +/// +/// @author Tim Rheinfels <tim.rheinfels@fau.de> +/// + +#ifndef BENCHMARKS_LUENBERGEROBSERVER_H +#define BENCHMARKS_LUENBERGEROBSERVER_H + +#include <benchmarks/memory.h> + +/// +/// @brief Memory required for running the Luenberger observer execution time benchmark +/// +#define BENCHMARKS_LUENBERGER_OBSERVER_WORK_SIZE (5u*N_Y_MAX*N_Y_MAX + 4u*N_Y_MAX) + +/// +/// @brief Runs the Luenberger observer execution time benchmark +/// +/// This function measures the execution time needed in order to update the +/// interval Luenberger observer proposed in @cite CSS21-Dkhil. +/// +void benchmarks_run_luenberger_observer(void); + +#endif // BENCHMARKS_LUENBERGEROBSERVER_H diff --git a/src/benchmarks/memory.c b/src/benchmarks/memory.c index 84a95f59df5d4f06436cb1fb0df4acbed88c702b..caaa789f5523b865c17971fe5d0bdfa3e40227a2 100644 --- a/src/benchmarks/memory.c +++ b/src/benchmarks/memory.c @@ -30,6 +30,7 @@ #include <benchmarks/blind_abstraction.h> #include <benchmarks/observer_abstraction.h> +#include <benchmarks/luenberger_observer.h> /// /// @brief Size of the padding in bytes used for checking memory violations @@ -51,7 +52,8 @@ /// static float work[2u * PADDING + MAX(BENCHMARKS_BLIND_ABSTRACTION_WORK_SIZE, - BENCHMARKS_OBSERVER_ABSTRACTION_WORK_SIZE) + MAX(BENCHMARKS_OBSERVER_ABSTRACTION_WORK_SIZE, + BENCHMARKS_LUENBERGER_OBSERVER_WORK_SIZE)) ]; float * const benchmarks_memory = work + PADDING;