# Outer enclosures of nonlinear mapping with degenerate ellipsoids

Morgan Louedec, L Jaulin, Christophe Viel

► To cite this version:

Morgan Louedec, L Jaulin, Christophe Viel. Outer enclosures of nonlinear mapping with degenerate ellipsoids. IFAC ACNDC 2024, 2024. hal-04467673

HAL Id: hal-04467673

https://hal.science/hal-04467673v1

Submitted on 20 Feb 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

Public Domain
# Outer enclosures of nonlinear mapping with degenerate ellipsoids*

M. Louédec* L. Jaulin** C. Viel***

* Lab-STICC, ENSTA-Bretagne, 2 rue François Verny 29200 Brest,
France (e-mail: morgan.louedec@centraliens-nantes.org)

** Lab-STICC, ENSTA-Bretagne, 2 rue François Verny 29200 Brest,
France (e-mail: lucjaulin@gmail.com)

*** CNRS, Lab-STICC, F-29806, Brest, France (e-mail:
christophe.viel@ensta-bretagne.fr)

**Abstract:** The propagation of ellipsoids is a tool that can be used to study the stability and the convergence of non-linear systems with time. Recent works have presented a method to enclose the propagation of ellipsoids via nonlinear functions. While this method is effective in high dimensions, it requires the Jacobian matrix of the function to be invertible. In this paper, we propose to generalize this method to include non-invertible Jacobians, using degenerate ellipsoids. The generalized method is then applied to a nonlinear and/or synchronous hybrid system.

**Keywords:** N-dimensional systems, Stability of nonlinear systems, p-invariance, degenerate ellipsoids, guaranteed enclosures, set propagation, hybrid system.

## 1. INTRODUCTION

Non-degenerate ellipsoids have been used in a variety of approaches to study nonlinear systems in terms of reachability analysis Kousik et al. (2022), stability analysis Lian and Wu (2020); Rauh et al. (2022), and state estimation Althoff and Rath (2021); Valiauga et al. (2021). These approaches often involve the propagation of ellipsoids via nonlinear mappings. Since the result of this propagation is often no longer an ellipsoid, the result is often approximated by inner and outer ellipsoidal enclosures.

These approaches follow previous works such as Polyak et al. (2006); Messerer and Diehl (2021) which studied linear or linearized systems where the propagation of ellipsoids follow an affine arithmetic and where stability problems are often equivalent to LMI constraints.

Zonotope set type is an alternative to ellipsoids and is used in some papers such as Kochdumper and Althoff (2021); Rego et al. (2020). According to Althoff and Rath (2021), for linear time-invariant systems, the zonotopic method has better tightness than ellipsoidal methods but higher computational costs when using set propagation. Thus, ellipsoidal methods are often more relevant to study high-dimensional systems.

Several works such as Dompierre et al. (2017); Blanchini and Miani (2015); Poznyak et al. (2014); Calbert et al. (2022) have proposed various mathematical tools for operations on ellipsoids. In the recent paper Rauh and Jaulin (2021), a numerical guaranteed method has been proposed to compute the ellipsoidal enclosures of a nonlinear mapping. The polynomial complexity of this

method allows the study of high-dimensional mappings. This method can also study mapping whose analytical expression is unknown, provided that the Jacobian matrix of the mapping is enclosed in a known box. The main assumption of this method is that the Jacobian matrix is invertible. In Louedec et al. (2023), this method has been used to prove the existence of positive invariant ellipsoids with respect to n-dimensional time-independent nonlinear systems.

However, some mappings do not correspond to this assumption. The Jacobian matrix can be non-invertible when there are equality constraints. This singularity is common in systems with discrete measurement updates or shock effects. Some papers, such as Balandin et al. (2020); Becis-Aubry (2021), have proposed the use of degenerate ellipsoids in the singular case when the system is linear. Degenerate ellipsoids are harder to use in practice because they don't have the quadratic form of non-degenerate matrix because they don't have an invertible shape matrix. Linear methods work well in practice with degenerate ellipsoids because they don't need matrix inversion. However, for non-linear systems, the method of Rauh and Jaulin (2021) relies on matrix inversion.

In this paper, we propose to generalize the method of Rauh and Jaulin (2021) by using degenerate ellipsoids to solve the singular case where the Jacobian matrix is not invertible. This generalization is made using the singular value decomposition of the shape matrix and the addition of new singular values when necessary. The performances of the generalized method in terms of computational complexity and pessimism are similar to those of the initial method.

The paper is structured in the following way. The key concepts and the original method are presented in Section 2.

<footnote>* This work has been supported by the French Defence Innovation Agency (AID) and the Brittany Region</footnote>
The generalized method is presented in Section 3. Section 4 presents an algorithmic implementation of the method using interval analysis. Finally, Section 5 discusses the performance of the method in applications with degenerate ellipsoids.

## 2. NOTATION AND HYPOTHESIS

The theory of the propagation method presented in this paper is a modified version of the method presented in Section 2.1. It relies on degenerate ellipsoids presented in Section 2.2, used to solve the problem presented in Section 2.3.

### 2.1 Initial method

The method provided in Rauh and Jaulin (2021) considers the multi-dimensional state space of state vector $x \in \mathbb{R}^n$, the nonlinear mapping

$$y = f(x), f: \mathbb{R}^n \to \mathbb{R}^n, \quad (1)$$

and the non-degenerate ellipsoids $\mathcal{E}(\mu, \Gamma)$, given by the form

$$\mathcal{E}(\mu, \Gamma) = \{x \in \mathbb{R}^n | (x - \mu)^T \Gamma^{-1} \Gamma^{-1} (x - \mu) \le 1\}, \quad (2)$$

with the midpoint $\mu \in \mathbb{R}^n$ and the positive definite shape matrix $\Gamma\Gamma^T > 0$.

Given the initial domain $x \in \mathcal{E}(\mu_x, \Gamma_x)$, the method can compute an enclosing ellipsoid $\mathcal{E}(\mu_y, \Gamma_y)$ such that $y \in \mathcal{E}(\mu_y, \Gamma_y)$ using Theorem 1.

**Theorem 1.** (Rauh and Jaulin, 2021, Theorem 1) Consider the nondegenerate ellipsoid $\mathcal{E}(\mu_x, \Gamma_x)$. For a differential function (1) with

$$A = \frac{\partial f}{\partial x}(\mu_x) \quad (3)$$

invertible, $\mathcal{E}(\mu_y, \Gamma_y)$ is an outer enclosure of the solution set $f(\mathcal{E}(\mu_x, \Gamma_x))$ with

$$\mu_y = f(\mu_x), \quad (4)$$

$$\Gamma_y = (1 + \rho) \cdot A \cdot \Gamma_x, \quad (5)$$

where

$$\rho = \max_{\|\tilde{x}\| \le 1} \|\tilde{b}(\tilde{x})\|, \quad (6)$$

and

$$\tilde{b}: \mathbb{R}^n \to \mathbb{R}^n \\ \tilde{x} \to \Gamma_x^{-1} \cdot A^{-1} \cdot (f(\Gamma_x \tilde{x} + \mu_x) - f(\mu_x)) - \tilde{x}. \quad (7)$$

### 2.2 Degenerate ellipsoids

While the ellipsoids are often described by the quadratic form (2), they can also be presented as an affine transformation of the unit sphere, given by

$$\mathcal{E}(\mu, \Gamma) = \{x \in \mathbb{R}^n | \exists \tilde{x} \in \mathbb{R}^n, x = \mu + \Gamma \cdot \tilde{x}, \|\tilde{x}\| \le 1\}, \\ = \mu + \Gamma \cdot \mathcal{E}(0, I_n). \quad (8)$$

A degenerate ellipsoid is an ellipsoid whose shape matrix $\Gamma$ is not invertible. While degenerate ellipsoids can't be described by the quadratic form (2), they can be defined as a singular affine transformation with (8). In the problem

presented in Section 2.1, degenerate ellipsoid can appear if the matrix $A$ or $\Gamma_x$ are not invertible.

The general method proposed in this paper uses the singular value decomposition (SVD) of the matrix $A\Gamma_x$, given by

$$A\Gamma_x = U \cdot \Sigma \cdot V^T, \quad (9)$$

with $U \in \mathbb{R}^{n \times n}$ and $V \in \mathbb{R}^{n \times n}$ orthonormal and $\Sigma \in \mathbb{R}^{n \times n}$ diagonal. The diagonal elements of $\Sigma$ are the singular values $\sigma_i \ge 0$.

For a matrix $M$, let us define $\text{pinv}(M)$ the pseudo-inverse of $M$.

### 2.3 Problem formulation

In this paper, we propose to generalize Theorem 1 by considering the case of degenerate ellipsoids, i.e. singular case where the Jacobian matrix is not invertible. Consider the mapping (1), the initial domain $x \in \mathcal{E}(\mu_x, \Gamma_x)$ given by (8), and the matrix $A$ given by (3). Assuming $\Gamma_x A$ is noninvertible, the problem is to find an outer enclosure $\mathcal{E}(\mu_y, \Gamma_y)$ such that $y \in \mathcal{E}(\mu_y, \Gamma_y)$.

## 3. GENERAL METHOD

The modifications of Theorem 1 to solve the problem of Section 2.3 are described in the Theorem 2 and it illustrated in Figure 1.

**Theorem 2.** Consider the non-degenerate ellipsoid

$\mathcal{E}(\mu_x, \Gamma_x)$. For a differential function (1) with

$$A = \frac{\partial f}{\partial x}(\mu_x), \quad (10)$$

$\mathcal{E}(\mu_y, \Gamma_y)$ is an outer enclosure of the solution set $f(\mathcal{E}(\mu_x, \Gamma_x))$ with

$$\mu_y = f(\mu_x), \quad (11)$$

$$\Gamma_y = (1 + \rho) \cdot M, \quad (12)$$

where

$$\rho = \max_{\|\tilde{x}\| \le 1} \|\tilde{b}(\tilde{x})\|, \quad (13)$$

and

$$\tilde{b}: \mathbb{R}^n \to \mathbb{R}^n \\ \tilde{x} \to W \cdot (f(\Gamma_x \tilde{x} + \mu_x) - \mu_y) - Z\tilde{x}, \quad (14)$$

and where the matrices $M \in \mathbb{R}^{n \times n}$, $W \in \mathbb{R}^{n \times n}$ and $Z \in \mathbb{R}^{n \times n}$ are given by

* case 1 (General): if $A\Gamma$ is invertible, then

$$M = A\Gamma_x, \quad (15)$$

$$W = M^{-1}, \quad (16)$$

$$Z = I_n. \quad (17)$$

* case 2 (Singularity): if $A\Gamma$ is non invertible, then

$$M = US, \quad (18)$$

$$W = \text{pinv}(S) U^T, \quad (19)$$

$$Z = W A \Gamma_x, \quad (20)$$

$$S = \text{diag}(s_1, s_2, \dots, s_n), \quad (21)$$
where

$$s_i = \begin{cases} \sigma_i & \text{if } \sigma_i > 0, \\ q_i \cdot \max_{x \in \mathcal{E}(\mu, \Gamma_x)} \left| e_i^T U^T (f(x) - \mu_y) \right| & \text{else,} \end{cases} \quad (22)$$

where $U$, $V$ and $\sigma_i$ are defined from the SVD decomposition (9), the unitary vector $e_i$ from the Cartesian base of $\mathbb{R}^n$ written $\{e_1, \dots, e_n\}$ and the design parameters $q_i > 0$.

**Proof.** Case 1 has been studied in Rauh and Jaulin (2021). Let us study case 2.

Let $x \in \mathcal{E}(\mu_x, \Gamma_x)$. From (8), there is $\tilde{x} \in \mathbb{R}^n$ such that $\|\tilde{x}\| \le 1$ and

$$x = \mu_x + \Gamma_x \cdot \tilde{x}. \quad (23)$$

Then, in order to have a vector $\tilde{y}$ which plays the role of $\tilde{x}$ for the point $y \in \mathcal{E}(\mu_y, \Gamma_y)$, we propose

$$\tilde{y} = Z\tilde{x} + b(\tilde{x}). \quad (24)$$

Let us show (24) is a relevant solution. From (14), one has

$$\tilde{y} = W \cdot (f(\Gamma_x \tilde{x} + \mu_x) - \mu_y), \quad (25)$$

$$\stackrel{(23)}{=} W \cdot (f(x) - \mu_y), \quad (26)$$

$$\stackrel{(1)}{=} W \cdot (y - \mu_y), \quad (27)$$

$$\stackrel{(19)}{=} \text{pinv}(S) U^T \cdot (y - \mu_y). \quad (28)$$

Then, from (18), it follows that

$$M \cdot \tilde{y} = US\text{pinv}(S)U^T(y - \mu_y). \quad (29)$$

From Appendix (7.1), it can be shown that

$$S\text{pinv}(S)U^T(y - \mu_y) = U^T(y - \mu_y). \quad (30)$$

Using (30) and since $U$ is orthonormal, (29) can be rewritten as

$$M \cdot \tilde{y} = UU^T(y - \mu_y), \quad (31)$$

$$M \cdot \tilde{y} = (y - \mu_y), \quad (32)$$

$$y = \mu_y + M \cdot \tilde{y}, \quad (33)$$

$$y \stackrel{(12)}{=} \mu_y + \Gamma_y \cdot \left( \frac{1}{1+\rho} \cdot \tilde{y} \right) \quad (34)$$

Then, from Appendix (7.2), one has

$$\|Z\| \le 1. \quad (35)$$

Therefore, because $\|\tilde{x}\| \le 1$, one gets $\|Z\tilde{x}\| \le 1$. Then, from (24), (13) and the triangle inequality

$$\left\| \frac{1}{1+\rho} \cdot \tilde{y} \right\| \le \frac{1}{1+\rho} \cdot \left( \|Z\tilde{x}\| + \|\tilde{b}(\tilde{x})\| \right), \quad (36)$$

$$\le \frac{1+\rho}{1+\rho}, \quad (37)$$

$$\le 1. \quad (38)$$

Finally, from (8), (34) and (38), one has $y \in \mathcal{E}(\mu_y, \Gamma_y)$ which implies

$$f(\mathcal{E}(\mu_x, \Gamma_x)) \subseteq \mathcal{E}(\mu_y, \Gamma_y). \quad (39)$$

**Remark 3.** The parameter $q_i$ can be tuned to reduce the value of $\rho$ to reduce pessimism. One should consider $q_i = 1$ by default.

## 4. IMPLEMENTATION

Theorem 2 can be implemented in an algorithm which computes $\mathcal{E}(\mu_y, \Gamma_y)$, as in Rauh and Jaulin (2021). In

practice, the values of $\rho$ and $S$ are computed using tools from interval analysis such as those in Jaulin et al. (2001).

To solve (13), the parameter $\rho$ is overestimated by the upper interval bound of the box $[\tilde{b}]$ as

$$\rho \le \sup \left\{ \| [\tilde{b}] \| \right\}, \quad (40)$$

with

$$\tilde{b}(\tilde{x}) \in [\tilde{b}] := (W \cdot [J_f] \cdot \Gamma_x - Z) \cdot [\tilde{x}], \quad (41)$$

$$\tilde{x} \in [\tilde{x}] := [-1; 1]^n, \quad (42)$$

where $[J_f]$ is the interval extension of $\frac{\partial f}{\partial x}$ over the tightest axis-aligned box $[x]$ which contains the ellipsoid $\mathcal{E}(\mu_x, \Gamma_x)$. (41) is deduced from the center form representation of $[\tilde{b}]$ as presented in (Rauh and Jaulin, 2021, Equation (35)).

Then, to solve (22) when $\sigma_i = 0$, the parameter $s_i$ is overestimated by the upper interval bound of the box $[s_i]$ as

$$0 < s_i \le \sup \{[s_i]\}, \quad (43)$$

with

$$[s_i] := q_i \cdot e_i^T U^T \cdot [J_f] \cdot \Gamma_x \cdot [\tilde{x}]. \quad (44)$$

In practice, if $s_i = 0$ is known, selecting $q_i = 0$ limit pessimism and is recommended. Note that overestimation of $\rho$ and $s_i$ don't invalidate the proof of Theorem 2.

These new parameter estimations can then be implemented in Algorithm 1, a modified version of algorithm (Rauh and Jaulin, 2021, Algorithm 1). The parameter $\gamma > 0$ is used to detect singularity and should have a small value.

**Algorithm 1.** General Outer ellipsoidal enclosure

Require: $\{\mu_x, \Gamma_x, A, [J_f], \gamma\}$

1: **if** $\det(A\Gamma_x) > \gamma$ **then**

2: do as in (Rauh and Jaulin, 2021, Algorithm 1)

3: **else**

4: $U, \Sigma, V = \text{full\_svd\_decomposition}(A\Gamma_x)$

5: **for** $i \in [1: n]$ **do**

6: **if** $\sigma_i < \gamma$ **then**

7: $[s_i] = q_i \cdot e_i^T U^T \cdot [J_f] \cdot \Gamma_x \cdot [\tilde{x}]$

8: $s_i = \sup \{[s_i]\}$

9: **else**

10: $s_i = \sigma_i$

11: **end if**

12: **end for**

13: $M = US$

14: $W = S^+U^T$

15: $Z = S^+\Sigma V^T$

16: $[\tilde{b}] = (W \cdot [J_f] \cdot \Gamma_x - Z) \cdot u$

17: $\rho = \sup \left\{ \| [\tilde{b}] \| \right\}$

18: $\mu_y = f(\mu_x)$

19: $\Gamma_y = (1 + \rho) \cdot M$

20: **end if**

return $\{\mu_y, \Gamma_y\}$

## 5. APPLICATION ON NONLINEAR SYSTEM

A common case of singularity appears when one considers hybrid systems like robots, whose dynamics are continuous
<figure><img src="image_1.png" /><figcaption>Fig. 1. Illustration of Theorem 2 in the singular case.</figcaption></figure>

physical variables but their memory store discrete estimations.

Robots are usually controlled by their state estimation which is often corrected when the robots make measurements, and predicted between measurements.

When there is no measurement noise, the estimation correction can be perfect, making the estimation vector equal to the state vector. This correction is a common cause of singularity as it can create degenerate ellipsoids. Thus, to illustrate our method, we propose to apply the general ellipsoidal propagation to this type of system. The system is defined in Sections 5.1. Section 5.2 gives an example of ellipsoidal propagation. Section 5.3 presents another example where the ellipsoidal propagation is used to prove the p-invariance of the initial ellipsoid. Section 5.4 comments on the obtained results.

## 5.1 System description

Consider the dynamical system described by

$$ 
\begin{aligned}
\dot{\mathbf{x}} &= \mathbf{f}_x(\mathbf{x}, \hat{\mathbf{x}}, \boldsymbol{\alpha}), \\
\dot{\hat{\mathbf{x}}} &= \mathbf{f}_x(\hat{\mathbf{x}}, \hat{\mathbf{x}}, 0), \text{ for } t \in ]t_k; t_{k+1}[, \\
\hat{\mathbf{x}}(t_k^+) &= \mathbf{x}(t_k)
\end{aligned}
\quad (45) $$

where $\mathbf{x} \in \mathbb{R}^n$ is the state vector, $\hat{\mathbf{x}} \in \mathbb{R}^n$ is the state estimation, $\boldsymbol{\alpha} \in \mathbb{R}^n$ is a constant disturbance, and where $\mathbf{f}_x$ is differentiable and the time $t_k^+$ represents the time $t_k$ after the measurement. At every measurement time $t_k = k \cdot T$ with $k \in \mathbb{N}$ and the period $T > 0$, the piece-wise continuous estimation $\hat{\mathbf{x}}$ is corrected on the value of the state $\mathbf{x}(t_k)$.

The evolution of the system can be described by the augmented state vector $\mathbf{z} = [\mathbf{x}^T \hat{\mathbf{x}}^T \boldsymbol{\alpha}^T]^T \in \mathbb{R}^{3n}$ whose

evolution is described by

$$ \dot{\mathbf{z}} = \mathbf{f}_z(\mathbf{z}), \text{ for } t \in ]t_k; t_{k+1}[, \quad (46) $$

$$ \mathbf{z}(t_k^+) = \mathbf{P}_z \cdot \mathbf{z}(t_k), \quad (47) $$

where

$$ \mathbf{f}_z(\mathbf{z}) = \begin{bmatrix} \mathbf{f}_x(\mathbf{x}, \hat{\mathbf{x}}, \boldsymbol{\alpha}) \\ \mathbf{f}_x(\hat{\mathbf{x}}, \hat{\mathbf{x}}, 0) \\ 0 \end{bmatrix}, \quad (48) $$

and

$$ \mathbf{P}_z = \begin{bmatrix} 1 & 0 & 0 \\ 1 & 0 & 0 \\ 0 & 0 & 1 \end{bmatrix} \otimes \mathbf{I}_n. \quad (49) $$

Thus, the evolution of the system from one measurement time $t_k$ to the next one $t_{k+1}$ can be described by the nonlinear mapping

$$ \mathbf{z}_{k+1} = \mathbf{f}(\mathbf{z}_k) = \mathbf{z}_k + \phi(T, \mathbf{z}_k), \quad (50) $$

where $\phi$ is the flow function of (46).

In order to illustrate the system in low dimensions, let us consider a system with $n = 1$ and

$$ \mathbf{f}_x(\mathbf{x}, \hat{\mathbf{x}}, \boldsymbol{\alpha}) = 0.5 \cdot \sin(\mathbf{x}) - \hat{\mathbf{x}} + \boldsymbol{\alpha}. \quad (51) $$

Note that our method can easily treat high-dimension problems with polynomial calculation time, see, Louedec et al. (2023).

## 5.2 Ellipsoidal propagation

The first test consider the initial ellipsoid $\mathcal{E}(\boldsymbol{\mu}_0^+, \boldsymbol{\Gamma}_0^+)$ a time $t_0^+$ with

$$ \boldsymbol{\mu}_0^+ = [1 \ 1 \ 0]^T, \quad (52) $$

$$ \boldsymbol{\Gamma}_0^+ = 0.05 \cdot \begin{bmatrix} 1 & 0 & 0 \\ 1 & 0 & 0 \\ 0 & 0 & 1 \end{bmatrix}. \quad (53) $$
<figure><img src="image_2.png" /><figcaption>Fig. 2. Ellipsoidal propagation for the system (46) with $z(t_k) \in \mathcal{E}(\mu_k, \Gamma_k)$ and $z(t_k^+) \in \mathcal{E}(\mu_k^+, \Gamma_k^+)$ for $k \in \{0, 1, 2, 3\}$. Ellipsoids are projected on $\{\mathbf{x}, \hat{\mathbf{x}}\}$.</figcaption></figure>

The ellipsoids $\mathcal{E}(\mu_k, \Gamma_k)$ and $\mathcal{E}(\mu_k^+, \Gamma_k^+)$ are computed for respectively $t_k$ and $t_k^+$ with $k \in \{1, 2, 3\}$ using the Algorithm 1. The parameters $\mathbf{q} = [q_1 \ q_2 \ q_3]^T$ and $\gamma$ are empirically set to $\mathbf{q} = [0 \ 0 \ 10]^T$ and $\gamma = 10^{-6}$. The result is illustrated in Figure 2, where the ellipsoids are projected on the two first dimensions.

## 5.3 p-invariant neighborhood

The second test consider the initial ellipsoid $\mathcal{E}(\mu_0^+, \Gamma_0^+)$ with

$$ \mu_0^+ = [0 \ 0 \ 0]^T, \qquad (54) $$

$$ \Gamma_0^+ = 0.5 \cdot \begin{bmatrix} 1 & 0 & 0 \\ 1 & 0.5 & 0 \\ 0 & 0 & 0.1 \end{bmatrix}. \qquad (55) $$

The ellipsoids $\mathcal{E}(\mu_1, \Gamma_1)$ and $\mathcal{E}(\mu_1^+, \Gamma_1^+)$ are computed. The parameters $\mathbf{q}$ and $\gamma$ are unchanged.

Then, the inclusion $\mathcal{E}(\mu_1^+, \Gamma_1^+) \subseteq \mathcal{E}(\mu_0^+, \Gamma_0^+)$ is tested to verify that $\mathcal{E}(\mu_0^+, \Gamma_0^+)$ is p-invariant for the time $t_k^+$ according to (Soyer et al., 2020, Theorem IV.5). The ellipsoids are illustrated in Figure 3.

## 5.4 Results and discussion

As expected, the ellipsoids are degenerate at the times $t_k^+$ with $\hat{x}(t_k) = \mathbf{x}(t_k)$ but become non-degenerate at time $t_{k+1}$ because of the bias $\alpha$. Moreover, the nonlinearity of the system creates some pessimism in the ellipsoidal propagation. A smaller initial set leads to smaller pessimism. Pessimism has also been reduced by selecting the value of $\mathbf{q}$. For now, one needs to find it empirically, but a method to find it automatically will be the subject of future works. Note that pessimism is incremental and thus increase at every iteration.

<figure><img src="image_3.png" /><figcaption>Fig. 3. Propagation of the p-invariant ellipsoid $\mathcal{E}(\mu_0^+, \Gamma_0^+)$ until the next measurement of the system 46. Ellipsoids are projected on $\{\mathbf{x}, \hat{\mathbf{x}}\}$.</figcaption></figure>

# 6. CONCLUSION

In this paper, the ellipsoid propagation method from Rauh and Jaulin (2021) has been generalized for non-degenerate ellipsoids with similar computational time and precision. The method propagates degenerate ellipsoids through an n-dimensional nonlinear mapping. This method makes it possible to analyze high-dimensional synchronous nonlinear hybrid systems, which is a challenging topic in the field of mobile robotics.

Future work will study the stability of such a system using the proposed method, to prove the existence of p-invariant ellipsoids.

# REFERENCES

* Althoff, M. and Rath, J.J. (2021). Comparison of guaranteed state estimators for linear time-invariant systems. Automatica, 130, 109662.

* Balandin, D.V., Biryukov, R.S., and Kogan, M.M. (2020). Ellipsoidal reachable sets of linear time-varying continuous and discrete systems in control and estimation problems. Automatica, 116.

* Becis-Aubry, Y. (2021). Ellipsoidal constrained state estimation in presence of bounded disturbances. In 2021 European Control Conference (ECC), 555–560.

* Blanchini, F. and Miani, S. (2015). Set-Theoretic Methods in Control. Systems & Control: Foundations & Applications. Springer International Publishing, Cham.

* Calbert, J., Egidio, L.N., and Jungers, R.M. (2022). An Efficient Method to Verify the Inclusion of Ellipsoids. arXiv. 2211.06237 [math.OA].

* Dompierre, J., Mokwinski, Y., Vallet, M.G., and Guibault, F. (2017). On ellipse intersection and union with application to anisotropic mesh adaptation. Engineering with Computers, 33(4), 745–766.

* Jaulin, L., Kieffer, M., Didrit, O., and Walter, E. (2001). Applied Interval Analysis. Springer, London.
Kochdumper, N. and Althoff, M. (2021). Sparse Polyno- mial Zonotopes: A Novel Set Representation for Reach- ability Analysis. *IEEE Transactions on Automatic Con- trol*, 66(9), 4043–4058.

Kousik, S., Dai, A., and Gao, G.X. (2022). Ellipsotopes: Uniting Ellipsoids and Zonotopes for Reachability Analysis and Fault Detection. *IEEE Transactions on Automatic Control*, 1–13.

Lian, J. and Wu, F. (2020). Stabilization of Switched Lin- ear Systems Subject to Actuator Saturation via Invari- ant Semiellipsoids. *IEEE Transactions on Automatic Control*, 65(10), 4332–4339. Conference Name: IEEE Transactions on Automatic Control.

Louedec, M., Jaulin, L., and Viel, C. (2023). Computa- tional tractable guaranteed numerical method to study the stability of n-dimensional time-independent nonlin- ear systems with bounded perturbation. *Automatica*, 153, 110981.

Messerer, F. and Diehl, M. (2021). An Efficient Algorithm for Tube-based Robust Nonlinear Optimal Control with Optimal Linear Feedback. In *2021 60th IEEE Confer- ence on Decision and Control (CDC)*, 6714–6721. ISSN: 2576-2370.

Polyak, B.T., Nazin, A.V., Topunov, M.V., and Nazin, S.A. (2006). Rejection of Bounded Disturbances via Invariant Ellipsoids Technique. In *Proceedings of the 45th IEEE Conference on Decision and Control*, 1429– 1434. ISSN: 0191-2216.

Poznyak, A., Polyakov, A., and Azhmyakov, V. (2014). Attractive Ellipsoids in Robust Control. Systems & Con- trol: Foundations & Applications. Springer International Publishing, Cham.

Rauh, A., Bourgois, A., and Jaulin, L. (2022). Verifying Provable Stability Domains for Discrete-Time Systems Using Ellipsoidal State Enclosures. *Acta Cybernetica*.

Rauh, A. and Jaulin, L. (2021). A computationally in- expensive algorithm for determining outer and inner enclosures of nonlinear mappings of ellipsoidal domains. *International Journal of Applied Mathematics and Com- puter Science*, 31(3), 399–415.

Rego, B.S., Raffo, G.V., Scott, J.K., and Raimondo, D.M. (2020). Guaranteed methods based on constrained zonotopes for set-valued state estimation of nonlinear discrete-time systems. *Automatica*, 111, 108614.

Soyer, M., Olaru, S., and Fang, Z. (2020). From constraint satisfactions to periodic positive invariance for discrete- time systems. In *2020 59th IEEE Conference on Deci- sion and Control (CDC)*, 4547–4552. ISSN: 2576-2370.

Valiauga, P., Feng, X., Villanueva, M.E., Paulen, R., and Houska, B. (2021). Set-membership Estimation using Ellipsoidal Ensembles. *IFAC-PapersOnLine*, 54(3), 596– 601.

7. APPENDIX

7.1 Proof of equation (30)

Consider **S** from (21), **e**<sub>i</sub> from (22), **U** from (20), **y** from (2) and **μ**<sub>y</sub> from (11). the matrix **S** is diagonal, so its pseudo inverse pinv (**S**) is also diagonal and its diagonal elements are

$$
p_i = \begin{cases} \frac{1}{s_i} & \text{if } s_i > 0, \\ 0 & \text{else.} \end{cases} \tag{56}
$$

Let **c** = **SS**<sup>+</sup>**U**<sup>T</sup> (**y** − **μ**<sub>y</sub>) and **d** = **U**<sup>T</sup> (**y** − **μ**<sub>y</sub>). One has

$$
\begin{align}
e_i^T c &= e_i^T S S^+ U^T (\mathbf{y} - \boldsymbol{\mu}_y), \nonumber \\
&= s_i s_i^+ e_i^T U^T (\mathbf{y} - \boldsymbol{\mu}_y), \nonumber \\
&= s_i s_i^+ e_i^T d. \tag{57}
\end{align}
$$

From (22), if $e_i^T d = 0$, then $s_i s_i^+ = 0$, so $e_i^T c = e_i^T d$.
Moreover, if $e_i^T d \neq 0$, then $s_i s_i^+ = 1$, so $e_i^T c = e_i^T d$.
Therefore $c = d$ and (30) is verified.

7.2 Proof of equation (35)

One has

$$
\begin{align}
\| \mathbf{Z} \|^{(20)} &\stackrel{(19),(9)}{=} \| \mathbf{W} \mathbf{A} \mathbf{\Gamma}_x \|, \nonumber \\
&\stackrel{(19),(9)}{=} \| \operatorname{pinv}(\mathbf{S}) \mathbf{U}^T \mathbf{U} \mathbf{\Sigma} \mathbf{V}^T \|, \nonumber \\
&\stackrel{(19),(9)}{=} \| \operatorname{pinv}(\mathbf{S}) \mathbf{\Sigma} \mathbf{V}^T \| . \tag{58}
\end{align}
$$

Then, since **V** is orthonormal,

$$
\begin{align}
\| \mathbf{Z} \| &= \| \operatorname{pinv}(\mathbf{S}) \mathbf{\Sigma} \|, \nonumber \\
&\stackrel{(56)}{\le} \max_i (p_i \sigma_i). \tag{59}
\end{align}
$$

Moreover, from (21), if $\sigma_i \neq 0$, then $p_i = \frac{1}{\sigma_i}$. Thus
$$ \| \mathbf{Z} \| \le 1. $$