domain: math info

text:
Andreas Rauh
email: andreas.rauh@uni-oldenburg.de

Marit Lahme
email: marit.lahme@uni-oldenburg.de

Simon Rohou
email: simon.rohou@ensta-bretagne.fr

LUC JAULIN b , THACH Ngoc Dinh

R A Ïssi

AND Mohamed Fnadi
email: mohamed.fnadi@univ-littoral.fr

© A Rauh

L Jaulin
email: lucjaulin@gmail.com

T N Dinh

T Raïssi

OFFLINE AND ONLINE USE OF INTERVAL AND SET-BASED APPROACHES FOR CONTROL AND STATE ESTIMATION: A REVIEW OF METHODOLOGICAL APPROACHES AND THEIR APPLICATION

Keywords: Interval analysis, uncertain systems, parameter identification, state estimation, guaranteed stabilizing control, model-predictive control, interval observers

Control and state estimation procedures need to be robust against imprecisely known parameters, uncertainty in initial conditions, and external disturbances. Interval methods and other set-based techniques form the basis for the implementation of powerful approaches that can be used to identify parameters of dynamic system models in the presence of the aforementioned types of uncertainty. Moreover, they are applicable to a verified feasibility and stability analysis of controllers and state estimators. In addition to these offline approaches for analysis, interval and set-based methods have also been developed in recent years, which allow to solve the associated design tasks and to implement reliable techniques that are applicable online. The latter approaches include set-based model-predictive control, online parameter adaptation techniques for nonlinear variablestructure and backstepping controllers, interval observers, and fault diagnosis techniques. This paper provides an overview of the methodological background and reviews numerous practical applications for which interval and other set-valued approaches have been employed successfully.

Introduction

Interval analysis has become an active field of research over the past decades. Its fundamentals can especially be traced back to the works of R.E. Moore and his famous book published in 1966 [START_REF] Moore | Interval Arithmetic[END_REF].

In its widest sense, interval analysis provides an approach to perform computations on finitely large sets in such a way that the determined solutions contain all the solutions of the mathematically stated problem with certainty [START_REF] Jaulin | Applied Interval Analysis[END_REF][START_REF] Mayer | Interval Analysis and Automatic Result Verification[END_REF]. In such a way, it is not only possible to account for numerical errors due to the limited precision of the numbers representable in floating point arithmetic on a CPU, but also for discretization and method errors. The latter arise inevitably when solving, for example, ordinary differential equations with the help of temporal series expansions truncated after a finite order [START_REF] Bünger | A Taylor model toolbox for solving ODEs implemented in MATLAB/INTLAB[END_REF][START_REF] Kapela | CAPD::DynSys: A flexible C++ toolbox for rigorous numerical analysis of dynamical systems[END_REF][START_REF] Lohner | Enclosing the solutions of ordinary initial and boundary value problems[END_REF][START_REF] Nedialkov | Interval tools for ODEs and DAEs[END_REF][START_REF] Nedialkov | Implementing a rigorous ode solver through literate programming[END_REF]. In addition to these effects, which are small (at least) for a single evaluation step, but can cause severe damage if neglected, the same techniques also allow for computing with wider sets that result from imperfectly known system parameters and initial conditions in (technical) applications. This lack of knowledge is commonly referred to as epistemic uncertainty [START_REF] Der Kiureghian | Aleatory or epistemic? Does it matter?[END_REF]. As such, interval analysis provides a mathematically sound way to handle epistemic uncertainty that can furthermore also be traced back to other kinds of a lack of knowledge of the actual system behavior, for example, friction phenomena in rigid body-body dynamics of fluidic applications.

In addition to the application domains mentioned above, further high-level functionalities have been developed on the basis of interval analysis. These contain techniques for solving sets of nonlinear algebraic equations with the help of Newton-type methods [START_REF] Krawczyk | Newton-Algorithmen zur Bestimmung von Nullstellen mit Fehlerschranken[END_REF][START_REF] Neumaier | Interval Methods for Systems of Equations[END_REF] or the identification of feasible solution sets that satisfy (nonlinear) constraints in a guaranteed manner [START_REF] Dit Sandretto | Certified calibration of a cabledriven robot using interval contractor programming[END_REF]. In the latter class of problems, the SIVIA algorithm [START_REF] Jaulin | Applied Interval Analysis[END_REF][START_REF] Walter | Guaranteed characterization of stability domains via set inversion[END_REF] (set inversion via interval analysis) as well as contractor approaches form the fundamental background, commonly being combined with Newton-type methods. Moreover, extensions exist that allow for handling outliers in a suitable manner by relaxed set intersections [START_REF] Jaulin | Robust set-membership state estimation; application to underwater robotics[END_REF].

Besides the advantageous features mentioned above that interval techniques always provide guaranteed enclosures of a mathematically formulated problem, they also come with the disadvantage of a certain degree of pessimism. This pessimism is caused by the following two fundamental issues [START_REF] Cornelius | Computing the range of values of real functions with accuracy higher than second order[END_REF][START_REF] Jaulin | Applied Interval Analysis[END_REF][START_REF] Mayer | Interval Analysis and Automatic Result Verification[END_REF]:

• The dependency effect: The difference of two intervals with numerically identical bounds generally does not simplify to the value zero as shown in

[0 ; 1] -[0 ; 1] = [-1 ; 1] . (1.1)
This is due to the fact that the difference in the equation above does not carry any information on whether both arguments are identical quantities or not. Therefore, variables are generally treated as independent in any arithmetic operation according to

x • y|x ∈ [x] := [x ; x] , y ∈ [y] := y ; y , where • ∈ {+, -, •, /} , (1.2)
unless simplifications can be carried out on the basis of symbolic formula manipulation, which allows to exploit knowledge that the two arguments x and y represent exactly the same (physical) quantity. Therefore, in classical interval analysis, there do not exist any set-based inversions for the operations of set-based addition and multiplication. • The wrapping effect: Transformations of sets, for example, rotations and scalings even in linear matrix-vector products as in the recursive evaluation of discrete-time dynamic systems

[x](t k+1 ) = A • [x](t k ) (1.
3) with the exemplarily chosen initial conditions and dynamics matrix

x(t 0 ) ∈ [x](t 0 ) = [-1 ; 1] [-1 ; 1] , A = 1 2
may lead to an exponential increase in the volume of the solution enclosures

x(t k ) ∈ 2 k 2 • [-1 ; 1] [-1 ; 1] (1.5)
if axis-aligned interval box enclosures were determined in each step (1.3) despite the fact that the true volume would stay constant in the example above because the matrix A in (1.4) represents a volume-preserving rotation of the initial box by 45 • in each step.

In this example, changing the evaluation order according to

[x](t k+1 ) = A k+1 • [x](t 0 ) , (1.6)
where the integer-order powers of the point-valued matrix A are firstly determined, allows to solve this problem [START_REF] Krasnochtanova | Interval methods as a simulation tool for the dynamics of biological wastewater treatment processes with parameter uncertainties[END_REF]. However, there is no general solution available, especially in cases where the transformations are nonlinear or the matrices A may themselves include intervals in each element. Then, a large variety of other set-representations may turn out to be more efficient, for example, the use of ellipsoids [START_REF] Becis-Aubry | Ellipsoidal constrained state estimation in presence of bounded disturbances[END_REF][START_REF] Kurzhanskii | Ellipsoidal Calculus for Estimation and Control[END_REF] (being mapped exactly to an ellipsoid in pointvalued linear transformations), (constrained) zonotopes [START_REF] Ifqir | Zonotopic set-membership state estimation for switched systems[END_REF][START_REF] Rego | Joint state and parameter estimation based on constrained zonotopes[END_REF], zonotope bundles [START_REF] Althoff | Zonotope bundles for the efficient computation of reachable sets[END_REF], or polytopes [START_REF] Rauh | An ellipsoidal predictor-corrector state estimation scheme for linear continuous-time systems with bounded parameters and bounded measurement errors[END_REF][START_REF] Rohou | Exact bounded-error continuous-time linear state estimator[END_REF]. In nonlinear settings, the computationally inexpensive thick ellipsoid approach recently published by two of the authors provides the possibility to detect pessimism in the guaranteed outer enclosure by performing a comparison with an inner ellipsoidal enclosure belonging to the solution set with certainty [START_REF] Rauh | Union and intersection operators for thick ellipsoid state enclosures: Application to bounded-error discrete-time state observer design[END_REF][START_REF] Rauh | A computationally inexpensive algorithm for determining outer and inner enclosures of nonlinear mappings of ellipsoidal domains[END_REF]. Moreover, computationally more demanding techniques such as Taylor model arithmetic [START_REF] Bünger | Shrink wrapping for Taylor models revisited[END_REF][START_REF] Bünger | A Taylor model toolbox for solving ODEs implemented in MATLAB/INTLAB[END_REF][START_REF] Bünger | Preconditioning of Taylor models, implementation and test cases[END_REF][START_REF] Hoefkens | Rigorous Numerical Analysis with High-Order Taylor Models[END_REF] provide further means to reduce the wrapping effect that is omnipresent in any dynamic system simulation of solution finding to algebraic systems of equations or during the integration of sets of differential equations.

Finally, coordinate transformations, as derived in [START_REF] Kersten | State-space transformations of uncertain systems with purely real and conjugate-complex eigenvalues into a cooperative form[END_REF][START_REF] Kersten | Transformation of uncertain linear fractional order differential equations into a cooperative form[END_REF][START_REF] Mazenc | Asymptotically stable interval observers for planar systems with complex poles[END_REF] for systems with either real or complex eigenvalues after their linearization, as well as the exploitation of specific monotonicity properties, allow for reducing both the pessimism discussed above and the computational effort required for its reduction. An example for the exploitation of monotonicity properties is the design of interval observers [START_REF] Efimov | Interval state observer for nonlinear time varying systems[END_REF][START_REF] Raïssi | Interval state estimation for a class of nonlinear systems[END_REF] that allows for computing elementwise lower and upper bounds for the state trajectories according to the following scheme [START_REF] Smith | Monotone Dynamical Systems: An Introduction to the Theory of Competitive and Cooperative Systems[END_REF].

An autonomous set of ordinary differential equations

ẋ(t) = f (x(t)) , x ∈ R n (1.7)
is cooperative according to [START_REF] Angeli | Monotone control systems[END_REF][START_REF] Hirsch | On the nonchaotic nature of monotone dynamical systems[END_REF][START_REF] Smith | Monotone Dynamical Systems: An Introduction to the Theory of Competitive and Cooperative Systems[END_REF] if the property

x ⟨1⟩ i (t) ≥ x ⟨2⟩ i (t) (1.8)
holds for all elements i ∈ {1, . . . , n} of the two vectors with x ⟨1⟩ (0) and x ⟨2⟩ (0) being the initial conditions which satisfy again the inequalities

x ⟨1⟩ i (t 0 ) ≥ x ⟨2⟩ i (t 0
) for all i ∈ {1, . . . , n} .

(1.9)

For continuous-time systems, cooperativity can be checked by the sufficient sign conditions

J i,j (x) ≥ 0 , i, j ∈ {1, . . . , n} , i ̸ = j (1.10)
for all off-diagonal elements of the Jacobian

J (x) = ∂f ∂x (x) (1.11)
of the right-hand side of the state equations evaluated for all reachable states x = x(t).

Matrices satisfying this non-negativity property for the off-diagonal elements are also denoted as Metzler matrices in the literature [START_REF] Raïssi | Some recent results on the design and implementation of interval observers for uncertain systems[END_REF][START_REF] Raïssi | Interval state estimation for a class of nonlinear systems[END_REF][START_REF] Rauh | Interval and linear matrix inequality techniques for reliable control of linear continuous-time cooperative systems with applications to heat transfer[END_REF][START_REF] Smith | Monotone Dynamical Systems: An Introduction to the Theory of Competitive and Cooperative Systems[END_REF].

In the case of linear uncertain systems, this cooperativity property simplifies to an elementwise sign condition for the off-diagonal elements of the lower interval bounds of the uncertain dynamics matrix A ∈ A ; A in the dynamic system model

ẋ(t) = A • x(t) .
(1.12)

Guaranteed bounds for the state trajectories x(t) can be obtained in the form [x] (t) = [v(t) ; w(t)] so that all reachable states are included by an evaluation of the following coupled set of state equations (which are a direct consequence of [START_REF] Müller | Über die Eindeutigkeit der Integrale eines Systems gewöhnlicher Differenzialgleichungen und die Konvergenz einer Gattung von Verfahren zur Approximation dieser Integrale[END_REF])

vi (t) = inf ([A i,i ] • v i (t)) + inf    n j=1 j̸ =i [A i,j ] • [z j ] (t)    ẇi (t) = sup ([A i,i ] • w i (t)) + sup    n j=1 j̸ =i [A i,j ] • [z j ] (t)    (1.13)
with the resulting state bounds

[z i ] (t) = [v i (t) ; w i (t)] . (1.14)
In the equations above, the operators inf (•) and sup (•) denote the lower and upper interval bounds of the corresponding arguments. Note that couplings between the vectors v(t) and w(t) can be ignored in (1.13) if the system is positive, i.e., v i (t) ≥ 0 holds for all t ≥ 0, see for example [START_REF] Efimov | Interval state observer for nonlinear time varying systems[END_REF][START_REF] Raïssi | Interval state estimation for a class of nonlinear systems[END_REF]. There, this property is exploited for the design of a Luenberger-type interval observer, which allows for guaranteed state estimation in the case of bounded measurement uncertainty.

Moreover, recent developments of the so-called mixed monotonicity approach and remainder form decompositions [START_REF] Khajenejad | Interval observer synthesis for locally lipschitz nonlinear dynamical systems via mixed-monotone decompositions[END_REF][START_REF] Khajenejad | Tight remainder-form decomposition functions with applications to constrained reachability and guaranteed state estimation[END_REF] can be exploited as a relaxation of the aforementioned evaluation technique based on cooperativity. In this case, a cooperative subsystem is extracted from the dynamic system model, while the non-cooperative part is shifted in to an additive bounded error term. In such a way, the possibility for extracting bounding systems is preserved and the need for computing with general sets is still avoided.

For practical implementation of set-valued calculus, numerous software libraries exist. The following short list is given without any claim for completeness: IntLab [START_REF] Rump | IntLab -INTerval LABoratory[END_REF] for Matlab and Octave, the Octave interval package [START_REF] Heimlich | Octave forge interval[END_REF], Versoft [START_REF] Rohn | VERSOFT[END_REF], Juliaintervals [START_REF] Sanders | IntervalArithmetic.jl[END_REF], C-XSC [START_REF] Krämer | XSC Languages (C-XSC, PASCAL-XSC) -Scientific computing with validation, arithmetic requirements, hardware solution and language support[END_REF], DynIbex [START_REF] Dit Sandretto | DynIBEX: A differential constraint library for studying dynamical systems[END_REF], the CORA library [START_REF] Althoff | An introduction to CORA[END_REF], and codac [START_REF] Rohou | The Codac library: A catalog of domains and contractors[END_REF].

This paper is structured as follows. Sec. 2 gives an overview of the fundamentals of identifying parameters for dynamic systems with the help of interval and other setbased approaches. These identification results are the basis for the variable-structure control approaches reviewed in Sec. 3 as well as for the guaranteed model-predictive control approaches in Sec. [START_REF] Althoff | An introduction to CORA[END_REF]. Thereafter, we present approaches for set-based state estimation in Sec. 5, applicable to the identification of the open-circuit voltage characteristic of Lithiumion batteries, and an ellipsoidal technique in Sec. 6 that allows for treating set-valued and stochastic uncertainty and disturbances in a unified manner when solving state estimation tasks. Finally, conclusions and a brief outlook on future work are given in Sec. 7.

System Modeling and Verified Parameter Identification

The estimation of the parameters of both finite-dimensional sets of discrete-time difference equations as well as continuous-time ordinary differential equations is crucial for the design and implementation of all model-based control and state estimation procedures. Identification experiments are typically performed (for asymptotically stable systems) with predefined, time-dependent actuation signals so that measurements of selected outputs can be gathered at certain discrete points of time t k , k ∈ N 0 . Here, t 0 = 0 is assumed without loss of generality. Due to the collection of measurements at discrete points of time, it is necessary that parameter identification routines allow for a propagation of state information between two subsequent measurement instants t k-1 and t k , k ≥ 1. The following description of the set-based parameter identification approach is based on the more detailed discussions that can be found in [START_REF] Rauh | Interval methods and contractor-based branch-and-bound procedures for verified parameter identification of quasi-linear cooperative system models[END_REF].

2.1. Predictor-Corrector Framework for the Verified Identification of Parameters of Dynamic Systems. As summarized in Fig. 1, a first class of verified parameter identification routines is characterized by an observer-based structure. This structure is closely linked to the evaluation steps that are performed in classical state estimators such as Luenberger-type observers or (Extended) Kalman Filters [START_REF] Julier | A new approach for the nonlinear transformation of means and covariances in filters and estimators[END_REF][START_REF] Kalman | A new approach to linear filtering and prediction problems[END_REF]. Hence, prediction and correction steps are separated from each other and evaluated in an alternating manner.

In this framework, the detection of feasible parameter domains is performed by means of the following three-stage procedure [START_REF] Dötschel | Thermal behavior of high-temperature fuel cells: Reliable parameter identification and interval-based sliding mode control[END_REF][START_REF] Rauh | Experimental parameter identification for a control-oriented model of the thermal behavior of high-temperature fuel cells[END_REF]: (1) A verified evaluation of the set of state equations is performed between two subsequent measurement points. In the case of discrete-time state equations, the corresponding expressions are evaluated with the help of interval arithmetic (or one of the alternative set-valued approaches mentioned in the introduction of this article), while verified solvers for initial value problems for ordinary differential equations [START_REF] Dit Sandretto | DynIBEX: A differential constraint library for studying dynamical systems[END_REF][START_REF] Dit Sandretto | Validated explicit and implicit Runge-Kutta methods[END_REF][START_REF] Auer | Validated modeling of mechanical systems with SmartMOBILE: Improvement of performance by ValEncIA-IVP[END_REF][START_REF] Kapela | CAPD::DynSys: A flexible C++ toolbox for rigorous numerical analysis of dynamical systems[END_REF][START_REF] Kühn | Rigorous error bounds for the initial value problem based on defect estimation[END_REF][START_REF] Lin | Validated solutions of initial value problems for parametric ODEs[END_REF][START_REF] Nedialkov | Computing Rigorous Bounds on the Solution of an Initial Value Problem for an Ordinary Differential Equation[END_REF][START_REF] Nedialkov | Interval tools for ODEs and DAEs[END_REF][START_REF] Nedialkov | Implementing a rigorous ode solver through literate programming[END_REF] are required if continuous-time processes are taken into consideration. In any case, problem-specific approaches for the reduction of overestimation [START_REF] Rauh | Quantification of overestimation in interval simulations of uncertain systems[END_REF] -such as interval subdivision, preconditioning of state equations, symbolic simplifications, or enclosures allowing to trace the shape of the reachable state domains in a tight manner -shall be employed to minimize the pessimism in the computed state and output boundaries that is caused by both the dependency and wrapping effects [START_REF] Jaulin | Applied Interval Analysis[END_REF][START_REF] Lohner | Enclosing the solutions of ordinary initial and boundary value problems[END_REF]. (2) Information resulting from the prediction of state intervals up to the point at which new measurements are available are combined with information obtained from sensors. The applicable techniques are:

• A direct intersection of the enclosures of predicted state variables with directly measured state intervals (in terms of an additive superposition of a point-valued measurement with the respective tolerance interval); • The use of contractors and set inversion approaches to determine an inverse mapping from the measured outputs to the internal system states and parameters [START_REF] Ifqir | Zonotopic set-membership state estimation for switched systems[END_REF]. Example are forward-backward contractors [START_REF] Jaulin | Applied Interval Analysis[END_REF] initialized with the predicted state intervals that are intersected with the associated measurement outputs or interval Newton techniques in the case of nonlinear output equations or for output equations involving more than one state variable. (3) Parameter subdomains are eliminated either by the contractors for which details are presented in the following section or after a subdivision of parameter intervals with associated multiple evaluations of the state equations. In the latter case, a parameter subinterval is guaranteed to be inconsistent with the information provided by the specified system model and the available measurements if the result of the state prediction substituted into the measurement model and the interval-valued measured data do not overlap. For a visualization of this exclusion text, see Fig. 2, where a branch-and-bound procedure on the basis of multiple simulations of the state equations over the complete horizon of gathered measured data is depicted. The corresponding consistency tests correspond to those also used in the following subsection. 2.2. Simulation-Based Parameter Identification Procedures. Simulation-based parameter identification procedures, making use of multiple evaluations of the set of state equations for subintervals of the complete initial domain of possible system parameters, are a second option for parameter identification. The underlying branch-and-bound procedure1 relies on the following assumptions:

• measured data are available at discrete points of time, • worst-case bounds for measurement tolerances are known a-priori (modeling of the sensor uncertainty), • the model structure to be parameterized is assumed to be structurally correct,

• outer enclosures for the domains of possibly uncertain initial states are given, and • interval bounds on the uncertain parameters are known in terms of conservative overapproximations.

To exclude parts of the domains of uncertain initial states and parts of the a-priori given parameter intervals, a subdivision procedure of the respective domains is executed [START_REF] Csendes | Subdivision direction selection in interval methods for global optimization[END_REF]. To detect intervals for the parameterization of the system model which are guaranteed to be inconsistent, directly measured and simulated state intervals are intersected at identical points of time. If state variables are not directly measured, the contractor procedures mentioned in the previous subsection are equally applicable here. Parameter interval subdivisions are performed according to the following distinction of cases, where the directions in which the determined parameter intervals are subdivided with the help of a sensitivity analysis [START_REF] Dötschel | Thermal behavior of high-temperature fuel cells: Reliable parameter identification and interval-based sliding mode control[END_REF]. For accelerations by both additional (physics-inspired) contractors and GPUbased parallelizations, the readers are referred to [START_REF] Auer | Experiments-based parameter identification on the GPU for cooperative systems[END_REF][START_REF] Rauh | Interval methods and contractor-based branch-and-bound procedures for verified parameter identification of quasi-linear cooperative system models[END_REF]:

time t t 0 t 1 t 2 t 3 ... measured data y m (t k ) [ y m ](t 0 ) [ y m ](t 1 ) [ y m ](t 3 ) [ y m ](t 2 ) y m (t k )∈[ y m ](t k )
(a) Measurement information of the system output ym(t), sampled at the equidistant instants t0, t1, t2, . . ..

time t t 0 t 1 t 2 t 3 ... measured data y m (t k ) simulated output enclosure [ y m ](t 0 ) [ y m ](t 1 ) [ y m ](t 3 ) [ y m ](t 2 )
(b) Output trajectory including interval uncertainty for a possibly feasible (undecided) system parameterization (Case 1).

time t t 0 t 1 t 2 t 3 ... measured data y m (t k ) [ y m ](t 0 ) [ y m ](t 1 ) [ y m ](t 3 ) [ y m ](t 2 )
(c) Output trajectory including interval uncertainty for a guaranteed feasible system parameterization (Case 2).

time t t 0 t 1 t 2 t 3 ... measured data y m (t k ) [ y m ](t 0 ) [ y m ](t 1 ) [ y m ](t 3 ) [ y m ](t 2 )
(d) Output trajectory including interval uncertainty for a guaranteed infeasible system parameterization (Case 3).

Figure 2.

Simulation-based distinction between possibly feasible, guaranteed feasible, and infeasible system parameterizations [START_REF] Dötschel | Thermal behavior of high-temperature fuel cells: Reliable parameter identification and interval-based sliding mode control[END_REF][START_REF] Rauh | Interval methods and contractor-based branch-and-bound procedures for verified parameter identification of quasi-linear cooperative system models[END_REF].

Case 1: Parameter subintervals are yet undecided if their corresponding predicted output intervals overlap at least partially for all available sensors at each point of time at which measurements are available and are never outside the range of measured data. Case 2: A parameter subinterval is guaranteed to be consistent if the corresponding output intervals are true subintervals of the interval-valued measured data for each of the available sensors at each point of time. This subinterval is no longer investigated but stored in a list of guaranteed admissible boxes. Case 3: A parameter subinterval is guaranteed to be inconsistent if the corresponding output interval lies outside the range of the interval-valued measured data. This has to hold true for at least one of the available sensors and for at least one point in time. The resulting parameter subinterval is excluded from further evaluations and the underlying simulation can be aborted as soon as the inconsistency is detected.

Interval-Based Variable-Structure Control Implementation

As mentioned in the previous section, the verified parameter identification is the prerequisite for the implementation of guaranteed stabilizing control procedures. These include linear matrix inequality methods for systems with polytopic and norm-bounded parameter uncertainty [START_REF] Boyd | Linear Matrix Inequalities in System and Control Theory[END_REF][START_REF] Rauh | Interval and linear matrix inequality techniques for reliable control of linear continuous-time cooperative systems with applications to heat transfer[END_REF], interval-based gain scheduling techniques [START_REF] Kersten | Interval methods for robust gain scheduling controllers: An lmi-based approach[END_REF], and variable structure as well as backstepping-type controllers [START_REF] Rauh | Interval methods for robust sliding mode control synthesis of high-temperature fuel cells with state and input constraints[END_REF][START_REF] Rauh | Interval-based sliding mode control design for solid oxide fuel cells with state and actuator constraints[END_REF][START_REF] Rauh | Reliable control of high-temperature fuel cell systems using interval-based sliding mode techniques[END_REF][START_REF] Senkel | Experimental and numerical validation of a reliable sliding mode control strategy considering uncertainty with interval arithmetic[END_REF] without and with hard state constraints.

In this section, we review a robust control procedure that extends the sliding mode design methodology by means of the online use of interval analysis to prevent the violation of state constraints with certainty [START_REF] Rauh | Interval methods for robust sliding mode control synthesis of high-temperature fuel cells with state and input constraints[END_REF][START_REF] Rauh | Interval methods for variable-structure control of dynamic systems with state constraints[END_REF][START_REF] Senkel | Experimental and numerical validation of a reliable sliding mode control strategy considering uncertainty with interval arithmetic[END_REF]. For a similar approach, applicable to the backstepping design methodology, the reader is referred to [START_REF] Rauh | Interval-based implementation of robust variable-structure and backstepping controllers of single-input single-output systems[END_REF].

3.1. Fundamental Sliding Mode Control Laws. Both the treatment of strict inequality constraints and bounded interval uncertainty can be combined with first-and second-order sliding mode techniques [START_REF] Eker | Second-order sliding mode control with experimental application[END_REF][START_REF] Fridman | Higher order sliding modes[END_REF].

3.1.1. First-Order Sliding Mode Control. For a summary of the fundamental sliding mode control design of single-input single-output systems, consider the n-th order set of ordinary differential equations

ẋ(t) = ẋ1 (t) . . . ẋn-1 (t) ẋn (t) T = x 2 (t) . . . x n (t) a (x(t), p) + b (x(t), p)•v(t) T (3.1)
in nonlinear controller canonical form with the state vector x(t) ∈ R n and the control input

(b (x(t), p) ̸ = 0) v(t) = -a (x(t), p) + u(t) b (x(t), p) ∈ R . (3.2)
For this system model, the output variable is given by

y(t) = x 1 (t) . (3.3)
Due to the assumed structure of the system (3.1), the output variable (3.3) has the relative degree n [START_REF] Marquez | Nonlinear Control Systems[END_REF], corresponding to the property that the n-th time derivative x (n) 1 (t) of the system output is the lowest-order derivative explicitly depending on either of the control inputs v(t) and u(t). Therefore, the output y(t) corresponds to a (trivial) flat system output [START_REF] Fliess | Flatness and defect of nonlinear systems: Introductory theory and examples[END_REF] and an exact tracking of a sufficiently smooth desired trajectory x 1,d (t) can be achieved in the case of exactly known parameters.

On the basis of this desired trajectory, the corresponding tracking error (r = 0) and its r-th time derivative 2 are given by ξ(r)

1 = x (r) 1 -x (r) 1,d , r ∈ {0, 1, . . . , n} , (3.4)
with the sliding surface

s := s(t) = n-1 r=0 α r ξ(r) 1 (3.5)
in which the highest-order coefficient is normalized according to α n-1 = 1. To guarantee asymptotic stability for states exactly located on this sliding surface, the parameters α r have to fulfill the necessary and sufficient stability conditions for a Hurwitz polynomial of the order n -1.

A classical first-order sliding mode control can be derived with the help of the quadratic radially unbounded candidate for a Lyapunov function

V ⟨I⟩ = 1 2 s 2 > 0 for s ̸ = 0 . (3.6)
2 Time arguments are omitted, whenever the meaning is non-ambiguous.

For s ̸ = 0, (global) asymptotic stability of the dynamic system corresponds to the (global) negative definiteness of the corresponding time derivative

V ⟨I⟩ = s • ṡ = n-1 r=0 α r ξ(r) 1 • n-1 r=0 α r ξ(r+1) 1 < 0 . (3.7)
Following the derivation of the variable-structure control law according to [START_REF] Rauh | Interval-based sliding mode control design for solid oxide fuel cells with state and actuator constraints[END_REF][START_REF] Rauh | Reliable control of high-temperature fuel cell systems using interval-based sliding mode techniques[END_REF], the right-hand side of the inequality (3.7) is replaced by the more conservative formulation

n-1 r=0 α r ξ(r) 1 • n-1 r=0 α r ξ(r+1) 1 < -η • |s| = -η • n-1 r=0 α r ξ(r) 1 • sign (s) (3.8)
which guarantees global asymptotic stability for η > 0. Here, the actual choice of η significantly influences the dynamics and the maximum absolute values of the control signal in the reaching phase characterized by s ̸ = 0. As soon as s = 0 has been reached after a finite time duration, the control amplitudes depend (in the unperturbed case) on the actual choice of the reference trajectory x 1,d and on the coefficients α r . These latter values also have a major influence on the control amplitudes if non-modeled errors and disturbances act on the system dynamics and if the error signals ξ(r) 1 are corrupted by non-negligible measurement noise or state reconstruction errors, preventing the perfect tracking of the desired trajectory x 1,d (t) and making it impossible that s ≡ 0 holds perfectly after the end of the reaching phase.

The derivation of the control law is completed by enforcing that the second factor in (3.8) becomes proportional to the sign of the actual value of s according to

n-2 r=0 α r ξ(r+1) 1 + u -x (n) 1,d + η • sign (s) = -β • sign (s) with β > 0 . (3.9)
Using η := η + β > 0, the control law results in

u(t) = u ⟨I⟩ = x (n) 1,d - n-2 r=0 α r ξ(r+1) 1 -η • sign (s) . (3.10) 3.1.2.
Second-Order Sliding Mode Control. For a second-order sliding mode, both s = 0 and ṡ = 0 have to be ensured by the feedback controller [START_REF] Eker | Second-order sliding mode control with experimental application[END_REF][START_REF] Fridman | Higher order sliding modes[END_REF]. This can be achieved by adding a first-order lag dynamics on the left-hand side of

γ 1 ṡ + γ 0 s = n-1 r=-1 α r ξ(r) 1 with ξ(-1) 1 (t) := t 0 ξ1 (τ )dτ . (3.11)
For the sake of asymptotic stability, the coefficients γ 0 and γ 1 need to be strictly positive.

For α -1 = 0, this sliding surface has a proportional and differentiating characteristic, while α -1 ̸ = 0 means that the integral ξ(-1)

1
of the output error is additionally fed back in the definition of the sliding variable s.

To ensure both s = 0 and ṡ = 0 for the closed-loop control system, the Lyapunov function candidate, which was previously chosen as a function solely depending on the value of s, is redefined in the form

V ⟨II⟩ = 1 2 • s 2 + λ ṡ2
with the scaling factor λ > 0 .

(3.12)

The computation of the time derivative of (3.12) results in

V ⟨II⟩ = s • ṡ + ṡ • - λγ 0 γ 1 ṡ + λ γ 1 n r=0 α r-1 ξ(r) 1 < 0 , (3.13)
where λ = γ 1 > 0 is used without loss of generality [START_REF] Rauh | Interval methods for robust sliding mode control synthesis of high-temperature fuel cells with state and input constraints[END_REF]. Also according to [START_REF] Rauh | Interval methods for robust sliding mode control synthesis of high-temperature fuel cells with state and input constraints[END_REF], the stabilization of the closed-loop system towards s = 0 can be achieved by setting

V ⟨II⟩ < -η 1 • | ṡ| -η 2 • |s| • | ṡ| = -ṡ • sign ( ṡ) • (η 1 + η 2 • |s|) . (3.14)
This finally leads to the nonlinear feedback controller [21, generalized form of Eqs. ( 22), ( 23)]

u = u ⟨II⟩ = x (n) 1,d + 1 α n-1 • γ 0 ṡ -s - n-1 r=0 α r-1 ξ(r) 1 -sign ( ṡ) • (η 1 + η2 • |s|) (3.15)
with ηi ≥ η i > 0 for both i ∈ {1, 2}.

3.2.

Extension by One-and Two-Sided Barrier Functions. Both control laws u ⟨I⟩ and u ⟨II⟩ can be extended by one-and two-sided barrier functions [START_REF] Tee | Barrier Lyapunov functions for the control of output-constrained nonlinear systems[END_REF].

3.2.1. One-Sided State Constraints. For the case of a one-sided barrier, it is necessary that the (time-varying) state (respectively output) constraint

x 1 < x1,max := x 1,d + ∆x 1,max with ∆x 1,max > 0 (3.16)
is not violated for any point of time t > 0 if the initial conditions for the state vector x(t) at t = 0 are also compatible with this constraint. Moreover, it is necessary that the sliding surface s = 0 (equivalent to x 1 = x 1,d ) lies within the admissible operating range defined in (3.16).

Then, the extended Lyapunov function ansatz

V ⟨j,A⟩ = V ⟨j⟩ + V ⟨A⟩ > 0 for s ̸ = 0 with (3.17)

V ⟨A⟩ = ρ V • ln σ V •x 1,max x1,max -x 1 and x 1 < x1,max (3.18)
is introduced for both alternatives j ∈ {I, II}. In (3.18), ρ V > 0 is selected so that the singularity x1,max -x 1 = 0 represents a repelling potential and that control constraints are not violated for usual operating conditions. Then, the term V ⟨j⟩ has dominating influence in the neighborhood of s = 0, while σ V > 0 can be used to adapt the steepness of the barrier function near its singularity.

In the case of the first-order sliding mode, the time derivative of (3.17) can be computed as

V ⟨j,A⟩ = V ⟨j⟩ + V ⟨A⟩ < 0 with V ⟨A⟩ = ρ V x1,max • -x 1 • ẋ1,max + ẋ1 •x 1,max x1,max -x 1 , (3.19)
where V ⟨A⟩ does not explicitly depend on the system input u.

In analogy to the fundamental first-order sliding mode control law u ⟨I⟩ (t) derived from (3.9), the inequality

s • n-2 r=0 α r ξ(r+1) 1 + u -x (n) 1,d + η • sign (s) + 1 s • V ⟨A⟩ -β•sign(s) < 0 (3.20)
is obtained. Following the same steps as in the derivation of u ⟨I⟩ (t) yields the control law

u = u ⟨I,A⟩ = u ⟨I⟩ -s s 2 +ε • V ⟨A⟩ , (3.21)
in which 1 s in (3.20) has been approximated by s s 2 +ε with the small constant ε > 0 to ensure regularity of the control law u ⟨I,A⟩ (t) on the sliding surface s = 0. Moreover, this modification guarantees that the barrier function becomes inactive as soon as the control goal has been reached.

In analogy to (3.21), the second-order sliding mode control (3.15) can be extended by the barrier function (3.18). Following the same steps as in Eqs. (3.19)-(3.21) leads to

u = u ⟨II,A⟩ = u ⟨II⟩ -1 α n-1 • ṡ ṡ2 +ε • V ⟨A⟩ . (3.22)
3.2.2. Two-Sided State Constraints. For the case of two-sided state constraints, the Lyapunov functions V ⟨j⟩ , j ∈ {I, II}, are extended according to

V ⟨j,B⟩ = V ⟨j⟩ + V ⟨B⟩ > 0 for s ̸ = 0 , (3.23)
where the additive term V ⟨B⟩ is chosen in this paper so that state deviations |x 1 -x 1,d | ≥ χ are prevented. For the alternative option |s| ≥ χ, the reader is referred to [START_REF] Rauh | Interval methods for robust sliding mode control synthesis of high-temperature fuel cells with state and input constraints[END_REF]. Symmetric barriers can be enforced by specifying a barrier function in the form

V ⟨B⟩ = ρ V • ln χ2l χ2l -(x 1 -x 1,d) 2l with l ∈ N . (3.24)
Large values for l typically lead to the fact that resulting state trajectories come closer to the edges of the admissible operating range. As before, the time derivative of the additive term (3.24) is computed, which yields

V ⟨B⟩ = ρ V • 2l•(x 1 -x 1,d) 2l-1 •( ẋ1 -ẋ1,d) χ2l -(x 1 -x 1,d) 2l . (3.25)
In full analogy to Eqs. (3.19)-(3.22), the requirement V ⟨j,B⟩ < 0 for s ̸ = 0 (and ṡ ̸ = 0, resp.) leads to the extension

u ⟨I,B⟩ = u ⟨I⟩ -s -1 • V ⟨B⟩ (3.26)
of the first-order sliding mode controller and to

u ⟨II,B⟩ = u ⟨II⟩ -(α n-1 • ṡ) -1 • V ⟨B⟩ (3.27)
for the second-order case, where the same regularization strategies for the rational terms 1 s and 1 ṡ become necessary as before. The term ṡ (if required) is typically estimated by a low-pass filtered differentiation or by means of an observer.

Interval Extensions to Handle Bounded Parameter Uncertainty and State

Estimation Errors. To guarantee asymptotic stability despite bounded uncertainty in both the measured (resp., estimated) states x and system parameters p included in the control law, interval techniques can be applied at run-time during the execution of the previously derived control laws if they are extended according to [START_REF] Rauh | Interval methods for robust sliding mode control synthesis of high-temperature fuel cells with state and input constraints[END_REF][START_REF] Rauh | Interval-based sliding mode control design for solid oxide fuel cells with state and actuator constraints[END_REF]. x ∈ [x] = [x ; x] of the measured or estimated state vector. In addition, it is assumed that the system model is given as an n-th order set of ODEs (3.1) in nonlinear controller canonical form. If this is not the case, extended control approaches according to [START_REF] Senkel | Experimental and numerical validation of a reliable sliding mode control strategy considering uncertainty with interval arithmetic[END_REF] are applicable.

For the sake of controllability (and, hence, also for the existence of the following interval-based variable-structure controllers), it has to be guaranteed that

0 ̸ ∈ b ([x] , [p]) (3.28)
holds for all x ∈ [x] and p ∈ [p], describing, respectively, the reachable domain in the state-space and the guaranteed enclosure of all uncertain parameters. To handle the setvalued state and parameter uncertainty, the output tracking error and its r-th derivative are enclosed by the intervals ξ(r)

1 ∈ ξ(r) 1 = x (r) 1 -x (r) 1,d (3.29)
for each r ∈ {0, 1, . . . , n}. These tracking error intervals can be used to generalize the first-order variable-structure controller (without and with state constraints) according to

v ⟨I⟩ = -a([x],[p])+x (n) 1,d - n-2 r=0 αr• ξ(r+1) 1 -η•sign([s]) b([x],[p]) , (3.30)
v ⟨I,A⟩ = v ⟨I⟩ - 1 b([x],[p]) • [s] [s] 2 +ε • V ⟨A⟩ , (3.31)
and

v ⟨I,B⟩ = v ⟨I⟩ - 1 b([x],[p]) • [s] [s] 2 +ε • V ⟨B⟩ . (3.32)
Similarly, interval-based generalizations can be defined for all before-mentioned second-order formulations. For a detailed discussion of restrictions in the case that accounts for integrator extensions in the definition of the sliding surface s = 0, the reader is referred to [START_REF] Rauh | Interval methods for robust sliding mode control synthesis of high-temperature fuel cells with state and input constraints[END_REF]. In all expressions above, a

([x] , [p]), b ([x] , [p]
), V ⟨A⟩ , and [s] denote the intervaldependent evaluations of the corresponding entries of the state equations, the time derivatives of the barrier function, and the sliding surface, respectively. For the actual control implementation, all interval expressions are evaluated by means of a suitable interval library.

In previous work such as [START_REF] Rauh | Interval-based implementation of robust variable-structure and backstepping controllers of single-input single-output systems[END_REF][START_REF] Rauh | Interval methods for robust sliding mode control synthesis of high-temperature fuel cells with state and input constraints[END_REF][START_REF] Rauh | Interval-based sliding mode control design for solid oxide fuel cells with state and actuator constraints[END_REF][START_REF] Rauh | Interval methods for variable-structure control of dynamic systems with state constraints[END_REF], where this approach was applied successfully to the control of the thermal behavior of high-temperature solid oxide fuel cells as well as for the control of inverted pendulum systems, the C++ toolbox C-XSC [START_REF] Krämer | XSC Languages (C-XSC, PASCAL-XSC) -Scientific computing with validation, arithmetic requirements, hardware solution and language support[END_REF] was used. Required derivatives, necessary to transform general nonlinear state equations into the system representation assumed in this section, can easily be obtained with the help of algorithmic differentiation. The template-based library FADBAD++ has shown its efficiency for this purpose, because it cannot only be applied to expressions involving classical point-valued data types but also interval variables.

It has to be pointed out that the actual implementation of the control laws (3.30)-(3.32) is performed in such a way that a point-valued control signal v is chosen from the computed intervals so that it guarantees asymptotic stability regardless of the sign of b(x, p). According to [START_REF] Rauh | Interval methods for robust sliding mode control synthesis of high-temperature fuel cells with state and input constraints[END_REF], this is done by testing the negative definiteness of the Lyapunov function candidate for the infima and suprema v := inf{[v]} and v := sup{[v]}. To account for roundoff and representation errors, these values are inflated by a small constant ϵ > 0 to obtain the final set of candidates V := {v -ϵ, v + ϵ, v -ϵ, v + ϵ}, from which the signal with minimum absolute value is chosen that ensures V < 0 (or more generally V ⟨j,A⟩ < 0, V ⟨j,B⟩ < 0, j ∈ {I, II}) despite the considered interval uncertainty.

Guaranteed Model Predictive Control

Predictive control techniques are powerful approaches for controlling uncertain systems, where the corresponding dynamics are formulated initial value problems for finite-dimensional sets of nonlinear ordinary differential equations. Besides stabilizing the dynamics towards a desired reference trajectory, it is possible to formulate the involved optimization problem in such a way that constraints for the admissible state trajectories and system inputs are handled. Classical approaches for nonlinear model predictive control commonly do not account for interval bounds in initial states and system parameters. This issue is resolved in this section, by presenting a control procedure for interval-valued initial value problems in the form [START_REF] Fnadi | Exponential state enclosure techniques for the implementation of validated model predictive control[END_REF] 

      ẋ(t) = f (t, x(t), u, p) x 0 ∈ [x 0 ] ⊆ IR n u ∈ [u] ⊆ IR m p ∈ [p] ⊆ IR p , (4.1)
where the state vector is denoted by x(t), the vector of dynamic parameters by p, and the control vector by u. The sets

[x 0 ] = [x 10 ] . . . [x n0 ] T , [u] = [u 1 ] . . . [u m ]
T , and

[p] = [p 1 ] . . . [p p ]
T , expressed as interval boxes, are respectively the initial condition of the state vector, the interval-bounded input, and the set of feasible dynamic parameters. We assume that this interval-based initial value problem has a unique solution x(t)(t, x 0 , u, p) at t > 0 since f : R × R n × R m × R p → R n is supposed to be continuous in t and Lipschitz in x(t). For the use of the proposed solution algorithm, we further require that f is sufficiently smooth, i.e., of class C k . For an approach applicable to linear time-invariant and linear parameter-varying systems, which is similar to the one described subsequently, the reader is referred to [START_REF] Reis De Souza | Robust output feedback model predictive control of time-delayed systems using interval observers[END_REF][START_REF] Reis De Souza | Robust output feedback mpc for lpv systems using interval observers[END_REF][START_REF] Reis De Souza | Robust output feedback model predictive control for constrained linear systems via interval observers[END_REF].

4.1. Interval-Based Nonlinear Model Predictive Control. The interval-based nonlinear model predictive control approach presented in [START_REF] Fnadi | Experimental validation of a guaranteed nonlinear model predictive control[END_REF][START_REF] Fnadi | Exponential state enclosure techniques for the implementation of validated model predictive control[END_REF] is a generalization of the optimization approach in [START_REF] Rauh | Interval methods for optimal control[END_REF] and is based on the computation of interval bounds for the control sequence over a receding horizon which takes into account bounded uncertainties in the parameters of the dynamic system model and in the measured data. The control intervals are calculated such that the convergence to the set-point interval is ensured (i.e., x(t j ) → [x r ], ∀j, and all the state and input constraints are satisfied (i.e., x(t i ) ∈ [x(t i )] and u(t j ) ∈ [u(t j )], ∀i, j). The corresponding interval algorithm encompasses two stages [START_REF] Fnadi | Experimental validation of a guaranteed nonlinear model predictive control[END_REF], namely, • Filtering and branching: This first step provides a sequence of guaranteed input interval boxes at each time-step t k over the prediction horizon N p , denoted as [U(t

k )] = [u(t k )] × [u(t k+1 )] × . . . × [u(t k+Np-1 )
]. Branching and filtering procedures allow the computation of safe input intervals along the receding time horizon that satisfy the state constraints (i.e., ∀j, [x(t j )] ⊆ [x min (t j ), x max (t j )], where x min (t j ) and x max (t j ) are the bounds for the admissible domain for each state variable) and ensure convergence to the reference interval (i.e., [x k ] → [x r ]). • Interval optimization: Since safe inputs are computed over a finite time horizon, the optimization algorithm is launched to compute the optimal inputs [U(t k )] ⋆ by minimizing as much as possible a newly formulated interval objective function to reduce both the norm of the input intervals and the error between the predicted as well as the reference outputs.

Formulation and Minimization of Interval Cost Functions.

For the feasible intervals of input values that satisfy the constraints on the state trajectories and system inputs, an optimization procedure is required to find the optimal control box and -on this basis -a point-valued input u(t k ) that can be applied to the actuator of the considered system [START_REF] Fnadi | Experimental validation of a guaranteed nonlinear model predictive control[END_REF][START_REF] Fnadi | Exponential state enclosure techniques for the implementation of validated model predictive control[END_REF]. The continuous cost function to be minimized can be expressed over the prediction horizon T p = T c × N p , with the control update step size T c , as

J(x(t), u(t)) = t+Tp t F (x(τ ), u(τ ))dτ . (4.2)
Here, F is commonly chosen as a quadratic function that minimizes the norm of the inputs and the error between the predicted outputs x and the reference x r according to

F (x(τ ), u(τ )) = (x(τ ) -x r ) T Q (x(τ ) -x r ) + u(τ ) T Ru(τ ) , (4.3)
where Q and R are both positive (semi-)definite weighting matrices. Under the assumption of piecewise constant control inputs, the continuous objective function (4.2) takes the form

J(x(t), u(t)) = t+Tp t (x(τ ) -x r ) T Q (x(τ ) -x r ) + u(τ ) T Ru(τ ) dτ = Np k=1 t+kTc t+(k-1)Tc (x(τ ) -x r ) T Q (x(τ ) -x r ) dτ + Np k=1 t+kTc t+(k-1)Tc u(τ ) T Ru(τ )dτ . (4.4)
As shown in [START_REF] Fnadi | Experimental validation of a guaranteed nonlinear model predictive control[END_REF][START_REF] Fnadi | Exponential state enclosure techniques for the implementation of validated model predictive control[END_REF], the value of cost function (4.4) can be enclosed with the help of techniques from interval analysis. This is achieved by applying a validated integration method for the initial value problem under consideration. This method provides a list of tight enclosures [x 0 ], [x 1 ], . . . , [x K ] for the solution of the model (4.1). These enclosures are obtained starting from the initial conditions [x 0 ] with bounded, piecewise constant control sequences [u j ]. Here, each solution enclosure [x t ] is defined as a box so that x(t) ∈ [x t ] holds for all t ∈ [t + (k -1)T c ; t + kT c ]. By using these piecewise constant interval enclosures, bounds for the cost function are obtained according to

J(x t , u) ∈ T c k+Np-1 j=k [x [j:j+1] ] -x r T Q [x [j:j+1] ] -x r + [u j ] T R[u j ] ≤ sup    T c k+Np-1 j=k [x [j:j+1] ] -x r T Q [x [j:j+1] ] -x r + [u j ] T R[u j ]    , (4.5)
where [x [j:j+1] ] is a guaranteed enclosure of all reachable states over the temporal slice t ∈ [t j ; t j + T c ] and sup{•} signifies the upper bound of the corresponding interval.

Figure 3 gives an overview of the extension of the interval-based nonlinear modelpredictive control procedure by a further pre-stabilization of open-loop unstable plants with the help of linear matrix inequality (LMI) techniques. In such cases, the interval-based optimization procedure does not directly provide the overall control sequence but rather a feedforward control sequence [u ff ] that allows for optimizing the tracking behavior within the operating domain in which the underlying feedback controller possesses provable asymptotic stability properties [START_REF] Fnadi | Exponential state enclosure techniques for the implementation of validated model predictive control[END_REF]. This paper also demonstrates the successful implementation of the control procedure to the swing-up control of a rotary inverse pendulum.

Reference set-point

[𝐱 r ]

Guaranteed nonlinear model-predictive control Validated Simulation: Guaranteed State Enclosures

IVP-ODE [𝐮 1 ] IVP-ODE [𝐮 2 ] IVP-ODE [𝐮 𝑁 p -1 ] [𝐔 𝑘 ] = 𝐮 1 × ⋯ × [𝐮 𝑁 p -1 ]

Interval Optimization

An Interval Observer-Based Approach for the Identification of the Open-Circuit Voltage Characteristic of Battery Cells

In this section, a continuous-time interval observer design is combined with a set-based identification of state-dependent output nonlinearities for the online-identification of the model of a Lithium-ion battery. According to [START_REF] Lahme | Set-valued approach for the online identification of the open-circuit voltage of Lithium-ion batteries[END_REF], this approach makes use of an equivalent circuit model of the battery as described in [START_REF] Erdinç | A dynamic Lithium-ion battery model considering the effects of temperature and capacity fading[END_REF]. In a first stage, the state variables of the system model

ẋ(t) = A (σ(t)) • x(t) + b (σ(t)) • u(t) (5.1) =    0 0 0 0 -1 C TS (σ(t))•R TS (σ(t)) 0 0 0 -1 C TL (σ(t))•R TL (σ(t))    • x(t) +    -1 C Bat 1 C TS (σ(t)) 1 C TL (σ(t))    • u(t)
with the state vector

x(t) = σ(t) v TS (t) v TL (t) T (5.
2) are estimated. These state variables denote the state of charge of the battery as well as two voltages across resistor-capacitor sub-networks included in the equivalent circuit to represent delay phenomena between state-like changes of the terminal current and measurable variations of the terminal voltage. According to [START_REF] Lahme | Set-valued approach for the online identification of the open-circuit voltage of Lithium-ion batteries[END_REF], interval enclosures of the state variables are estimated after reformulation of both the state equations (5.1) and the output equation

ṽOC (σ(t)) = C (σ(t)) • x(t) (5.3) = η OC (σ(t)) • σ(t) = v OC (σ(t)) -v 0 -v 2 (5.4) = v 0 • e v 1 •σ(t) -1 σ(t) + v 3 + v 4 • σ(t) + v 5 • σ 2 (t) • σ(t) (5.5)
in a quasi-linear form in which state of charge dependent matrices and vectors A (σ(t)), b (σ(t)), and C (σ(t)), respectively, are separated from the state vector x(t) and the input variable u(t), which is represented by the terminal current i T .

In a second stage, an interval-based set intersection approach is subsequently used to identify the unknown characteristic of the open circuit voltage v OC (σ(t)). This overall approach is summarized in Fig. 4.

For the observer design, cooperativity of the system model (as discussed in the introduction of this article) as well as stability of the estimation error dynamics are ensured by choosing the observer gain matrix H according to [START_REF] Hildebrandt | Robust interval observer design for fractionalorder models with applications to state estimation of batteries[END_REF], in the form

H = h 1 0 0 T (5.6) with h 1 > 0.
In such a way, interval bounds x ∈ [x] := x ; x can be determined for the true, non-measurable state vector x by the decoupled lower and upper bounding systems [START_REF] Hildebrandt | Robust interval observer design for fractionalorder models with applications to state estimation of batteries[END_REF][START_REF] Raïssi | Some recent results on the design and implementation of interval observers for uncertain systems[END_REF], following Müller's theorem [START_REF] Müller | Über die Eindeutigkeit der Integrale eines Systems gewöhnlicher Differenzialgleichungen und die Konvergenz einer Gattung von Verfahren zur Approximation dieser Integrale[END_REF],

A O x + bu + Hy m ≤ ẋ ≤ A O x + bu + Hy m , x ∈ [x] , (5.7)
where Successive intersections of these interval boxes -at different points in time -are utilized to improve the approximation of the true v OC (σ(t)) characteristic according to Figure 5.

A O = A -HC ∈ A O ; A O (5.
For strategies, allowing for a reduction of the computational complexity due to a large number of interval boxes to be intersected, the interval merging routine according to [START_REF] Krasnochtanova | Interval methods as a simulation tool for the dynamics of biological wastewater treatment processes with parameter uncertainties[END_REF] possible by implementing a forward-backward contractor, as discussed from a methodological point of view in [START_REF] Jaulin | Applied Interval Analysis[END_REF]. These refinements have been demonstrated recently in [START_REF] Lahme | Online identification of the open-circuit voltage characteristic of Lithium-ion batteries with a contractor-based procedure[END_REF].

6. Combination of Set-Based and Stochastic Uncertainty Representations: An Application to Iterative-Learning Observer Design

In the final part of this article, a combination of stochastic and set-based (in this case, ellipsoidal) uncertainty representations is considered. This approach allows, on the one hand, for a rigorous quantification of predefined confidence levels in stochastic state estimation procedures. On the other hand, it allows for handling nonlinearities in such a way that the previously mentioned tolerance bounds are definitely not determined in an overly optimistic manner [START_REF] Rauh | Robust iterative learning observers based on a combination of stochastic estimation schemes and ellipsoidal calculus[END_REF][START_REF] Rauh | A stochastic design approach for iterative learning observers for the estimation of periodically recurring trajectories and disturbances[END_REF][START_REF] Rauh | Reliability assessment of an unscented kalman filter by using ellipsoidal enclosure techniques[END_REF].

In addition to the properties mentioned above, we present a formulation of the state estimation procedure that allows to iteratively enhance the estimation accuracy for the case of periodically recurring trajectories and disturbance profiles (e.g., due to errors caused by a deterministic model mismatch due to purposefully performed model simplifications). For that purpose, assume the discrete time system model

z k+1 = Φ (z k , p) • z k (6.1)
with the bounded parameter vector p ∈ p ; p . As discussed in the following, the vector z k ∈ R ñ consists of both system states and noise variables. Furthermore, assume that the ellipsoidal domain

E k µ k , Γ ′ k , r := z k ∈ R n (z k -µ k ) T Γ ′ k -T Γ ′ k -1 (z k -µ k ) ≤ r 2 (6.2)
specifies the confidence bound of a given percentage if the vector z k is normally distributed with the positive definite covariance matrix

Q ′ k = Γ ′ k Γ ′ k
T ≻ 0 and the ellipsoid midpoint

(expected value) µ k ∈ R ñ.
Here, the parameter r describes a magnification factor according to [START_REF] Wang | Confidence analysis of standard deviational ellipse and its extension into higher dimensional euclidean space[END_REF] associated with this confidence level.

6.1. Ellipsoidal Calculus for the Guaranteed Preduction of Covariance Ellipsoids. For the implementation of an ellipsoidal covariance prediction step (based on the references [START_REF] Rauh | Union and intersection operators for thick ellipsoid state enclosures: Application to bounded-error discrete-time state observer design[END_REF][START_REF] Rauh | Ellipsoidal enclosure techniques for a verified simulation of initial value problems for ordinary differential equations[END_REF][START_REF] Rauh | A computationally inexpensive algorithm for determining outer and inner enclosures of nonlinear mappings of ellipsoidal domains[END_REF]), we use the definitions

Γ k := r • Γ ′ k and Q k := r 2 • Q ′ k . (6.3)
For the compactness of notation, reformulate the system model (6.1) into the form

z k+1 = Φ (z k , p) • žk + Φ • µ k + Φ (z k , p) -Φ • µ k (6.4)
with z k = žk + µ k , where

z k ∈ E k = E k µ k , Γ ′ k , r . (6.5)
Here, E k denotes the uncertainty on the non-origin centered states z k ,

žk ∈ Ěk = Ěk 0, Γ ′ k , r (6.6)
the uncertainty on žk after shifting the ellipsoid to the origin, and

Φ = Φ (µ k , mid ([p])) (6.7)
is the midpoint approximation of the quasi-linear system matrix with

mid ([p]) = 1 2 • p + p . (6.8)
Let also □E k denote an axis-aligned enclosure of E k in the form of an ñ-dimensional interval box. An alternative to the definition (6.7) is given by

Φ = mid (Φ (□E k , [p])) , (6.9)
which is preferable if (6.7) and (6.9) strongly differ from each other in the case of large uncertainty.

Then, a confidence bound of the predicted states z k+1 of magnification r is given by the ellipsoid E k+1 µ k+1 , Γ ′ k+1 , r with the covariance

Q ′ k+1 = Γ ′ k+1 Γ ′ k+1
T for which the parameters are computed in the following steps.

P1 : Apply žk+1 = Φ (z k , p) • žk (6.10) to the ellipsoid Ěk in (6.6). The outer ellipsoid enclosure of the image set is described by an ellipsoid with the shape matrix

Qk+1 = α 2 k+1 • Γ k+1 • Γ T k+1 , (6.11)
where α k+1 ≥ 0 is the smallest value for which the LMI

M k+1 := Λ -Q -1 k Φ T (z k , p) • Φ-T Φ-1 • Φ (z k , p) -α 2 k+1 R k Λ ⪯ 0 (6.12)
is satisfied for all z k ∈ □E k and p ∈ [p] with

R k := Γ k • Γ T k . (6.13)
In (6.12), the symbol ⪯ denotes the negative semi-definiteness of the corresponding matrix expression and Λ is a preconditioning matrix chosen according to [START_REF] Rauh | Ellipsoidal enclosure techniques for a verified simulation of initial value problems for ordinary differential equations[END_REF]. P2 : Compute interval bounds for the term

b k = Φ (p) -Φ • µ k ∈ [b k ] (6.14)
which accounts for a non-zero ellipsoid midpoint with z k , Φ, and p defined according to (6.5), (6.7), and (6.8). Inflate the ellipsoid bound described by the shape matrix (6.11) according to [START_REF] Rauh | A computationally inexpensive algorithm for determining outer and inner enclosures of nonlinear mappings of ellipsoidal domains[END_REF] and its factorized shape matrix

Q k+1 = (1 + ρ O,k+1 ) 2 • Qk+1 , (6.15) ρ O,k+1 = sup α -1 k+1 • Γ -1 k • [b k ] . ( 6
Γ ′ k+1 = α k+1 • (1 + ρ O,k+1 ) • Φ • Γ ′ k . ( 6
x k+1 = A (x k , p) • x k + E (x k , p) • w k y k = C (x k , p) • x k + v k (6.19)
with the state vector x k ∈ R n , the measured output vector y k ∈ R m (m ≤ n), as well as the uncorrelated process and measurement noise vectors w k ∈ R nw and v k ∈ R nv , respectively. Let both noise vectors be normally distributed with the covariances C w,k and C v,k and vanishing mean. The iterative-learning observer is designed in this section to estimate the state vector x k as well as its uncertainty (expressed by its covariance) by a Kalman Filter-like procedure. This estimator does not only operate along the time domain k but also enhances the estimates successively from the trial i to the trial i + 1. Moreover, to allow for the identification of a systematic model mismatch, a lumped correction term δ k is added to the state equations (6.19) in the form

x k+1 = A (x k , p) • x k + E (x k , p) • w k + δ k . (6.20)
For the derivation of the learning-type framework, consider the trials ξ = i and ξ = i + 1 with the actually measured data y ξ m,k = C ξ k • x ξ k + v ξ k corresponding to the realizations of the general outputs y k in (6.19) for these two trials.

For the recursive formulation of the iterative-learning observer, we assume that both process and measurement noise are uncorrelated and normally distributed with zero mean. The superscript p denotes the result of the prediction step, while the superscript e refers to the estimation result as the outcome of the measurement-based innovation step.

Prediction

Step for State-and Parameter-Dependent Disturbance Inputs E (x k , p) [START_REF] Rauh | Robust iterative learning observers based on a combination of stochastic estimation schemes and ellipsoidal calculus[END_REF]. Define an ellipsoid for the augmented state vector

z k = x i k T x i+1 k T w i k T w i+1 k T T (6.21)
corresponding to the result of the preceding innovation step, augmented by the influence of the uncorrelated process noise in both trials i and i + 1 in the form

E e,i|i+1 k         µ e,i k µ e,i+1 k 0 0     , Γ e k , r     (6.22)
with as well as of the iteration-independent noise covariance C w,k , and the magnification factor r ≥ 1 as a user-defined degree of freedom. Then, the application of the system model

Γ e k =    Γ e,i|i+1 k
0 2n×2nw 0 2nw×2n C w,k 0 0 C w,k 1
Φ (z k , p) =     A x i k , p 0 n×n E x i k , p 0 n×nw 0 n×n A x i+1 k , p 0 n×nw E x i+1 k , p 0 nw×n 0 nw×n I nw×nw 0 nw×nw 0 nw×n 0 nw×n 0 nw×nw I nw×nw     (6.25)
to the ellipsoid (6.22) according to Sec. 6.1 yields the ellipsoid

E p,i|i+1 k+1         µ p,i k+1 µ p,i+1 k+1 0 0     , Γ p k+1 , r     (
i k = y i m,k -C i k • µ p,i k (6.27)
and

∆y i+1 k = y i+1 m,k -C i+1 k • µ p,i+1 k (6.28)
between the measured data in the trials i and i + 1, respectively, and the corresponding output forecasts based on the prediction step of the previous subsection. Using these output deviations, and assuming the standard detectability requirements of the Kalman Filter design to be satisfied, the expected values are updated according to [START_REF] Rauh | A stochastic design approach for iterative learning observers for the estimation of periodically recurring trajectories and disturbances[END_REF] by As shown in [START_REF] Rauh | A stochastic design approach for iterative learning observers for the estimation of periodically recurring trajectories and disturbances[END_REF], the estimation error covariance is minimized by the filter gain matrices H i+1 1,k and H i+1 2,k that are derived in the following, where the predicted covariance is partitioned in a blockwise manner according to

µ e,i k µ e,i+1 k = µ p,i k µ p,i+1 k + Hk • y i m,k y i+1 m,k -Hk Ck • µ p,i k µ p,i+1 k , ( 6
C p,i|i+1 k = C p A,k C p B,k ⋆ C p C,k , (6.36)
and ⋆ denotes blocks that can be inferred from the symmetry of the result. Moreover, the residual covariance is defined as is computed in Eq. (6.37), followed by an extraction of the upper left (2m × 2m) block that corresponds to the actually measured output quantities. This subblock extraction is performed by the multiplication with the projection matrix P in the first summand of Eq. (6.37).

P • Γe k • Γe k T • P T + Cv,k = C A,k C B,k ⋆ C C,k (6
The optimal iterative-learning observer gains, in the sense of a minimization of the estimation error covariance, jointly considering the trials i and i + 1 are given by

H i+1 1,k H i+1 2,k = C i k C p A,k + C i+1 k C p C,k T C i k C p B,k -C i+1 k C p C,k T • C A,k + C C,k ⋆ C B,k -C C,k C A,k -C B,k + C T B,k + C C,k -1 .
(6.42)

For a proof of this expression for the filter gain matrices, see [START_REF] Rauh | A stochastic design approach for iterative learning observers for the estimation of periodically recurring trajectories and disturbances[END_REF]. The reference [START_REF] Rauh | Robust iterative learning observers based on a combination of stochastic estimation schemes and ellipsoidal calculus[END_REF] demonstrates the successful application of this iterative-learning observer for the state and disturbance estimation of a Lithium-ion battery with periodically recurring input current profiles. Compared to the work [START_REF] Rauh | A stochastic design approach for iterative learning observers for the estimation of periodically recurring trajectories and disturbances[END_REF], where an Extended Kalman Filter methodology was used in the form of a point-valued linearization of the state and output equations during the covariance computation in the prediction and innovation stages, the combination with the ellipsoidal calculus allows for a significant improvement of the state estimation accuracy, also yielding significantly tighter uncertainty bounds.

Conclusions and Future Work

In this paper, a comprehensive overview of recent advances concerning the offline and online use of interval and other set-based approaches for control and state estimation has been given. We have equally addressed methodological approaches and selected fields of practical applications in our review, where we have put a focus on the tasks of parameter estimation, robust and model-predictive control, as well as state estimation procedures. Future work will focus especially on widening the field of applications by an investigation of control and state estimation tasks for distributed and large-scale, interconnected systems with parameter uncertainty and disturbances that are, among others, omnipresent in control scenarios in the Industry 4.0 as well as in the domain of energy systems. To cope with such systems, we will focus especially on the development of decentralized techniques for the parameter identification of large-scale system models, which will again form the basis for control and state estimator design. Moreover, we will further develop approaches that allow for removing restrictive monotonicity assumptions as already shown in works such as [START_REF] Khajenejad | Interval observer synthesis for locally lipschitz nonlinear dynamical systems via mixed-monotone decompositions[END_REF][START_REF] Khajenejad | Tight remainder-form decomposition functions with applications to constrained reachability and guaranteed state estimation[END_REF][START_REF] Wang | Interval observer design for uncertain discrete-time linear systems[END_REF].

Figure 1 .

 1 Figure1. Predictor-corrector scheme for combined state estimation and parameter identification[START_REF] Rauh | Interval methods and contractor-based branch-and-bound procedures for verified parameter identification of quasi-linear cooperative system models[END_REF].

Figure 3 .

 3 Figure 3. Overall structure of the interval-based nonlinear model predictive control approach, extended by an underlying pre-stabilization of the plant dynamics.

Figure 4 .

 4 Figure 4. Estimation of the open-circuit voltage and the state of charge.

Figure 5 .

 5 Figure 5. Identification of nonlinear dependencies using interval methods.

. 16 )

 16 P3 : Compute the updated ellipsoid midpoint µ k+1 = Φ • µ k (6.17

  )

  from pointwise evaluations of the quasi-linear system's output matrix according toC i k := C k µ p,i k , mid ([

. 37 )

 37 with the projection matrix P = I m×m 0 m×(n-m) and the trial-independent measurement noise feature introduced in[START_REF] Rauh | Robust iterative learning observers based on a combination of stochastic estimation schemes and ellipsoidal calculus[END_REF], aiming at a robustification of the innovation stage against nonlinearities, the matrix Γe k is obtained by the consideration of the quasi-linearity of the output equation with the help of the ellipsoid E is propagated through a quasi-linear system model in the form (6.1) with the associated system matrix (y i m,k , y i+1 m,k ∈ R m ) Φ (z k , p) -m)×m I (n-m)×(n-m) 0 (n-m)×n 0 (n-m)×n 0 (n-m)×m I (n-m)×(n-m) Φ (z k , p) ∈ R 2n×2n .The evaluation of this quasi-linear model yields an ellipsoid

Uncertain Environment Plant, including LMI-based pre-stabilization

	𝐮 min ; 𝐮 max	𝐱 min ; 𝐱 max
	𝐱 0
	𝐩
	𝑁 p
	𝑇 c
		Feedforward
		control 𝐮 ff
		Optimal box 𝐔 𝑘 ⋆
		Sensors
		[𝐱 𝑘 ]

  [START_REF] Auer | Experiments-based parameter identification on the GPU for cooperative systems[END_REF]) is a Metzler matrix[START_REF] Gennat | Computing guaranteed bounds for uncertain cooperative and monotone nonlinear systems[END_REF][START_REF] Kaczorek | Positive 1D and 2D Systems[END_REF][START_REF] Kersten | Transformation of uncertain linear fractional order differential equations into a cooperative form[END_REF] with non-negative off-diagonal elements and [y m ] := y m ; y m = y m + [-∆y m ; ∆y m ]

		(5.9)
	is the measured system output with bounded uncertainty.
	The estimated state of charge, enclosed by the interval [σ], together with the estimated
	open-circuit voltage
	ṽOC (t) = y * m (t) + vTS (t) + vTL (t) ,	(5.10)
	with vTS (t) ∈ [v TS (t)] and vTL (t) ∈ [v TL (t)], allow for constructing the interval box [Γ] =
	[σ(t)] × [ṽ OC (t)] with σ(t) ∈ [σ(t) ; σ(t)] and ṽOC (t) ∈ ṽOC (t) ; ṽOC (t) .

  .18) 6.2. Stochastic Iterative-Learning Observer Design for Quasi-Linear State Equations. Consider the quasi-linear discrete-time state-space representation

  .2.2. Innovation Step. The measurement-based innovation step makes use of the deviations ∆y

						6.26)
	in which the first 2n components of the midpoint represent the predicted expected value
	vector µ p,i k+1	T	µ p,i+1 k+1	T	T	p,i|i+1 and the covariance C k+1

is obtained by extracting the upper left (2n × 2n) block of the matrix product Γ p k+1 • Γ p k+1 T . For the simplified case, in which the disturbance input matrix E is constant, the reader is referred to

[START_REF] Rauh | Robust iterative learning observers based on a combination of stochastic estimation schemes and ellipsoidal calculus[END_REF]

.

6

Branching of the parameter domain by an interval subdivision approach with a subsequent simulationbased bounding of the output trajectories.