Weichao Liang

Nina H Amini

Paolo Mason

On estimation and feedback control of spin-1 2 systems with unknown initial states

Keywords: Spin systems, Quantum filtering, Stochastic master equations, Fidelity

In this paper, we consider stochastic master equations describing the evolutions of quantum systems interacting with electromagnetic fields undergoing continuous-time measurements. In particular, we study feedback control of quantum spin-1 2 systems in the case of unawareness of initial states and in presence of measurement imperfections. We prove that the fidelity between the true quantum filter and its estimation (with arbitrary initial state) converges to one under appropriate assumption on the feedback controller. This shows the asymptotic convergence of such filters. In the more general case of spin-J systems, we discuss heuristically the asymptotic behavior of the true quantum filter and the associated estimation, and the possibility of exponentially stabilizing such systems towards an eigenvector of the measurement operator by an appropriate feedback.

INTRODUCTION

Classical filtering (see, e.g [START_REF] Kallianpur | Stochastic filtering theory[END_REF]; [START_REF] Xiong | An introduction to stochastic filtering theory[END_REF]) optimizes the estimation of the state of a classical system from noisy observations. The quantum analogue was developed in the 1960s by [START_REF] Davies | Quantum stochastic processes[END_REF][START_REF] Davies | Quantum theory of open systems[END_REF] and extended by Belavkin in the 1980s [START_REF] Belavkin | On the theory of controlling observable quantum systems[END_REF][START_REF] Belavkin | Nondemolition measurements, nonlinear filtering and dynamic programming of quantum stochastic processes[END_REF][START_REF] Belavkin | Quantum stochastic calculus and quantum nonlinear filtering[END_REF], relying on quantum probability theory and quantum stochastic calculus [START_REF] Hudson | Quantum Ito's formula and stochastic evolutions[END_REF]; [START_REF] Hudson | An introduction to quantum stochastic calculus and some of its applications[END_REF]; [START_REF] Meyer | Quantum probability for probabilists[END_REF]). The modern treatment of quantum filtering has been established in [START_REF] Bouten | A discrete invitation to quantum filtering and feedback control[END_REF]. Roughly speaking, quantum filtering theory gives a matrix-valued stochastic differential equation called stochastic master equation, to describe the time evolution of the state of an open quantum system interacting with an electromagnetic field under homodyne detection.

In real experiments, different types of imperfections, such as detection inefficiencies and unawareness of initial states, may be present (see e.g., [START_REF] Sayrin | Realtime quantum feedback prepares and stabilizes photon number states[END_REF]). In the case of unawareness of initial states and known physical parameters such as detection inefficiencies, one considers estimated quantum filters which are designed based on measurements. The main question is whether the estimated filter "forgets" its initial state and has the same asymptotic behaviour as if it was correctly initialized. This problem can be posed when a feedback depending on the estimated quantum filter is applied. Since the observation process depends on the true quantum state, we deal with coupled stochastic master equations whose asymptotic behavior is not at all a trivial problem.

For uncontrolled system, this problem has been investigated in some recent papers. In the series of papers van [START_REF] Van Handel | Filtering, stability, and robustness[END_REF]Handel ( , 2009a,b),b), a sufficient observability condition has been established so that the convergence of the estimation towards the true state is guaranteed. However, such condition is not easy to verify even if the system is finite dimensional. In [START_REF] Diósi | Coupled Ito equations of continuous quantum state measurement and estimation[END_REF], the problem of convergence is addressed assuming either the true quantum filter or its estimated state is always pure. The authors prove that the fidelity is a sub-martingale for this case. Then, in [START_REF] Rouchon | Fidelity is a sub-martingale for discrete-time quantum filters[END_REF], by applying Uhlmann's technique (Nielsen and Chuang, 2000, Theorem 9.4), the author shows that the fidelity between the state of the discrete-time quantum filter and its associated estimated state is a sub-martingale via a Kraus map. However, this sub-martingale property of the fidelity cannot ensure the convergence of the estimated filter state towards the true one. In [START_REF] Amini | On stability of continuous-time quantum filters[END_REF] the authors show that the fidelity is a sub-martingale for continuous-time quantum filters with perfect measurement for arbitrary mixed states. By the quantum repeated interaction approach, such result has been extended to the continuous-time jump-diffusion stochastic master equations with general measurement imperfections in [START_REF] Amini | Stability of continuous-time quantum filters with measurement imperfections[END_REF]. Then, in [START_REF] Benoist | Large time behavior and convergence rate for quantum filters under standard non demolition conditions[END_REF], for quantum filters described by jump-diffusion stochastic differential equations, the authors show that when the control input is turned off, under perfect Quantum Non-Demolition (QND) measurements and a non-degeneracy assumption, the convergence is ensured.

Concerning the feedback stabilization of discrete-time QND measures, in [START_REF] Mirrahimi | Feedback generation of quantum Fock states by discrete QND measures[END_REF]; [START_REF] Amini | Feedback stabilization of discrete-time quantum systems subject to non-demolition measurements with imperfections and delays[END_REF], the authors show the robustness of their feedback strategies with respect to initial states. Following [START_REF] Bouten | On the separation principle in quantum control[END_REF], this is referred to as quantum separation principle. For continuous-time case, stabilization results for the case of unknown initial states have not been addressed so far. On the other hand, for angular momentum systems and without initialization imperfections, in [START_REF] Mirrahimi | Stabilizing feedback controls for quantum systems[END_REF], the authors show the asymptotic convergence towards a chosen eigenvector of the measurement operator J z . In the same context, in Liang et al. (2019a), we provide general conditions on the feedback controller and a local Lyapunov type condition which ensure exponential convergence towards a chosen eigenvector of J z .

In this paper, we first analyze the dynamics of quantum spin-1 2 systems in presence of feedback control. We suppose imperfections in measurements and unawareness of the initial state. We show that the true filter and filter estimate have the same asymptotic behavior under appropriate assumption on the feedback. For spin-J systems, we discuss heuristically the asymptotic behaviors of the true quantum filter and its associated estimated filter and the possibility of exponentially stabilizing such systems towards an eigenvector of the measurement operator J z by a candidate feedback controller. We precise that our results don't prove the robustness properties of feedback strategies considered in Liang et al. (2019a[START_REF] Liang | On exponential stabilization of spin-1 2 systems[END_REF]. Numerical simulations are provided in order to illustrate our results and to support the efficiency of the proposed candidate feedback.

MODEL DESCRIPTION

Here, we consider quantum spin-1 2 systems. The stochastic master equations describing the evolution of the true system state and the corresponding estimated state are given as follows, dρ t = F ût (ρ t )dt+L(ρ t )dt+G(ρ t ) dY t -2 ηM Tr(σ z ρ t )dt , dρ t = F ût (ρ t )dt+L(ρ t )dt+G(ρ t ) dY t -2 ηM Tr(σ z ρt )dt , where

• the quantum state of the spin-1 2 system is denoted as ρ, and belongs to the space S 2 := {ρ ∈ C 2×2 | ρ = ρ * , Tr(ρ) = 1, ρ ≥ 0}. The estimated state is denoted as ρ ∈ S 2 , • the matrices σ x , σ y and σ z correspond to the Pauli matrices.

• F u (ρ) := -i[ωσ z + uσ y , ρ], L(ρ) := M (σ z ρσ z -ρ) and G(ρ) := √ ηM σ z ρ + ρσ z -2Tr(σ z ρ)ρ . • Y t
denotes the observation process of the quantum spin-1 2 system, which is a continuous semi-martingale whose quadratic variation is given by [Y, Y ] t = t. Its dynamics satisfies dY t = dW t + 2 √ ηM Tr(σ z ρ t )dt, where W t is a one-dimensional standard Wiener process,

• ût := u(ρ t ) denotes the feedback controller as a function of the estimated state ρt , • ω is the difference between the energies of the excited state and the ground state, η ∈ [0, 1] describes the efficiency of the detector, and M > 0 is the strength of the interaction between the system and the probe.

By replacing dY t = dW t + 2 √ ηM Tr(σ z ρ t )dt in the equation above, we obtain the following matrix-valued stochastic differential equations describing the time evolution of the pair (ρ

t , ρt ) ∈ S 2 × S 2 , dρ t = F ût (ρ t )dt+L(ρ t )dt+G(ρ t )dW t ,
(1)

dρ t = F ût (ρ t )dt+L(ρ t )dt+2 ηM G(ρ t )Tr σ z (ρ t -ρt ) dt +G(ρ t )dW t . (2)
If u ∈ C 1 (S 2 , R), the existence and uniqueness of the solution of ( 1)-( 2) can be shown by similar arguments as in (Mirrahimi and van Handel, 2007, Proposition 3.5). Recall that a density operator can be uniquely characterized by the Bloch sphere coordinates (x, y, z) as

ρ = 1 + xσ x + yσ y + zσ z 2 = 1 2 1 + z x -iy x + iy 1 -z .
The vector (x, y, z) belongs to the ball 1) expressed in the Bloch sphere coordinates takes the following form

B := {(x, y, z) ∈ R 3 | x 2 + y 2 + z 2 ≤ 1}. The stochastic differential equation (
dx t = -ω eg y t - M 2 x t + ût z t dt-ηM x t z t dW t , ( 3a
)
dy t = ω eg x t - M 2 y t dt -ηM y t z t dW t , (3b)
dz t = -û t x t dt + ηM (1 -z 2 t )dW t .
(3c) The stochastic differential equation ( 2) in the Bloch sphere coordinates is given by,

dx t = -ω eg ŷt - M 2 xt + ût ẑt + ηM xt ẑt (ẑ t -z t ) dt -ηM xt ẑt dW t , (4a)
dŷ t = ω eg xt - M 2 ŷt + ηM ŷt ẑt (ẑ t -z t ) dt -ηM ŷt ẑt dW t , (4b) dẑ t = -û t xt -ηM (1 -ẑ2 t )(ẑ t -z t ) dt + ηM (1 -ẑ2 t )dW t . (4c) 3. CONVERGENCE PROPERTY OF QUANTUM SPIN-1 2 SYSTEMS
We focus on the fidelity F(ρ, ρ) which defines a "distance" between the state ρ and its estimation ρ. In the two-level case, the fidelity can be written in the following form

F(ρ, ρ) = Tr(ρρ) + 2 det(ρ) det(ρ).
The fidelity in the Bloch sphere coordinates is given by

F(ρ, ρ) = F(v, v) = 1 2 1+v v+ (1 -v 2 )(1 -v 2 ) ,
where v := (x, y, z) denotes the true state and v := (x, ŷ, ẑ) denotes the estimated state in Bloch sphere coordinates.

We have two following special cases,

(1) if F(ρ, ρ) = 1, then v = v; (2) if F(ρ, ρ) = 0, then v + v = 0 and v 2 = v 2 = 1.
In order to apply the Itô formula on the fidelity F(ρ, ρ), we need to show the unattainability of the boundary for ρ and ρ. By straightforward calculations, we can show that

{ρ ∈ S 2 | det(ρ) = 0} = {ρ ∈ S 2 | Tr(ρ 2 ) = 1},
(5) which means that the boundary ∂S 2 is equal to the set of all pure states P. The following lemma states some invariance properties for Equations ( 1)-(2). Lemma 3.1. If ρ 0 > 0, then P(ρ t > 0, ∀t ≥ 0) = 1. Moreover, if η = 1, ∂S 2 × S 2 is a.s. invariant for Equations ( 1)-( 2). The same results hold true for S 2 ×∂S 2 .

Proof. The dynamics of the purification function S(ρ t ) := 1 -Tr(ρ 2 t ) is given by

dS(ρ t ) =M (1 -η)(1 -z 2 t ) 2 -(1 -ηz 2 t )S(ρ t ) dt
-2 ηM z t S(ρ t )dW t . Then, if η = 1, it is obvious that the set of all pure states P for Equation ( 1) is a.s. invariant.

Next, let us prove the first part of the lemma. Given ε > 0, consider any C 2 function on S such that

V (ρ) = 1 S(ρ) , if S(ρ) > ε.
We find

L V (ρ) = M 1 + 3ηz 2 -(1 -η) 1 -z 2 2S(ρ) V (ρ), ≤ 4M V (ρ) if S(ρ) > ε,
where L V denotes the infinitesimal generator of V (see e.g., [START_REF] Mao | Stochastic differential equations and applications[END_REF]). To conclude the proof, one applies standard arguments (see e.g., (Liang et al., 2019a, Lemma 4.1)). Roughly speaking, by setting f (ρ, t) = e -4M t V (ρ), one has L f ≤ 0 whenever S(ρ) > ε. From this fact one proves that the probability of S(ρ) becoming zero in a finite fixed time T is proportional to ε and, being the latter arbitrary, it must be 0. Due to the equality (5), P(ρ t > 0, ∀t ≥ 0) = 1 when ρ 0 > 0. The last part of the lemma can be proved in the same manner.

Next, we analyze the behavior of ρ t (resp. ρt ) whenever the corresponding initial datum ρ 0 (resp. ρ0 ) lies at the boundary of S 2 . Denote ρ g :=diag(1, 0) and ρ e :=diag(0, 1), which are the eigenvectors of σ z and pure states. Lemma 3.2. Assume η ∈ (0, 1) and u ∈ C 1 (S 2 , R). Suppose that ρ0 lies in ∂S 2 \ {ρ g , ρ e }, then ρt > 0 for all t > 0 almost surely. Moreover, if ρ0 ∈ {ρ g , ρ e } and u(ρ 0 ) = 0 then, ρt > 0 for all t > 0 almost surely. In addition, under the assumption ρ 0 ∈ ∂S 2 \ {ρ g , ρ e }, then ρ t > 0 for all t > 0 almost surely. Also, if u(ρ g )u(ρ e ) = 0 then, for all ρ 0 ∈ ∂S 2 , ρ t exits the boundary in finite time and stays in the interior of S 2 almost surely.

Proof. First, consider the purification function S(ρ) := 1 -Tr(ρ 2 ) for Equation (2), whose dynamics is given by

dS(ρ t ) =M (1 -η)(1 -ẑ2 t ) 2 -(1 -ηẑ 2 t )S(ρ t ) -4η(z t -ẑt )ẑ t S(ρ t ) dt -2 ηM ẑt S(ρ t )dW t . Now, assume ρ0 ∈ ∂S 2 \ {B (ρ g ) ∪ B (ρ e )}}. By compact- ness, there exists a ζ > 0 such that 1 2 M (1 -η)(1 -ẑ2 ) ≥ ζ. Define τ := inf{t > 0| ρt / ∈ ∂S 2 \ {B (ρ g ) ∪ B (ρ e )}
, for all ρ0 ∈ ∂S 2 \ {B (ρ g ) ∪ B (ρ e )} and t > 0, by Itô formula, we have

E(S(ρ t∧τ )) = E t∧τ 0 1 2 M (1 -η)(1 -ẑ2 s )ds ≥ ζE(t∧τ ).
By continuity and the definition of τ , S(ρ t∧τ ) = 0 almost surely. This implies that E(t ∧ τ ) = 0. Since we have E(t ∧ τ ) ≥ tP(τ ≥ t) we deduce that P(τ ≥ t) = 0 for all t > 0.

Due to the arbitrariness of , if ρ0 ∈ ∂S 2 \ {ρ g , ρ e } then ρt exits the boundary immediately. Combining with the strong Markov property and Lemma 3.1, we have ρt > 0 for all t > 0, almost surely. Moreover, if ρ0 ∈ {ρ g , ρ e } then by the condition u(ρ 0 ) = 0 we deduce the same result.

For the case ρ 0 ∈ ∂S 2 \ {ρ g , ρ e }, the above arguments can be repeated so that ρ t > 0 for all t > 0 almost surely. Moreover, if u(ρ g )u(ρ e ) = 0 and u ∈ C 1 (S 2 , R), then there exists a neighborhood of ρ e denoted by B re (ρ e ) and a neighborhood of ρ g denoted by B rg (ρ g ) such that, u(ρ) = 0 for all ρ ∈ B re (ρ e ) ∪ B rg (ρ g ). By applying the same arguments as in (Liang et al., 2019a Then for all (ρ 0 , ρ0 ) ∈ S 2 × S 2 , either ρt converges to {ρ g , ρ e } or F(ρ t , ρt ) converges to one, almost surely. In particular, if u(ρ e )u(ρ g ) = 0 then F(ρ t , ρt ) converges to one almost surely.

Proof. In order to study the asymptotic behavior of F(ρ, ρ), in the following we will apply Itô formula. For this purpose, we first need to show the C 2 regularity of F on appropriate invariant sets. If η = 1, then by Lemma 3.1, ∂S 2 ×∂S 2 , ∂S 2 ×int(S 2 ), int(S 2 )×∂S 2 and int(S 2 )×int(S 2 ) are invariant for the coupled system (1)-( 2) almost surely. Moreover, if either ρ or ρ belongs to the boundary of S 2 , then the fidelity takes the form F(ρ, ρ) = Tr(ρρ), which is a C 2 function. For the case η ∈ (0, 1), under the assumptions that u ∈ C 1 (S 2 , R) and u(ρ e )u(ρ g ) = 0, by Lemma 3.2, (ρ t , ρt ) exits the boundary in finite time and stays in int(S 2 ) × int(S 2 ) afterwards almost surely. Note that the fidelity function is C 2 in int(S 2 ) × int(S 2 ).

Consider the Lyapunov function V(ρ, ρ) := 1-F(ρ, ρ). Denote Ξ := (1 -v 2 )(1 -v 2 ). For any u ∈ C 1 (S 2 , R), the infinitesimal generator of F(ρ, ρ) is given by

L F (ρ, ρ) = M (1 -η) 4Ξ (1 -ẑ2 )(1 -v 2 ) + (1 -z 2 )(1 -v 2 ) + 2ẑ 2 (1 -v v -Ξ)Ξ -2(1 -z ẑ)Ξ + M 2 (1 -ẑ2 )(1 -v v -Ξ). (6)
In particular, for η = 0, we have

L F (ρ, ρ) = M 2 (1 -ẑ2 )(1 -v 2 ) + (1 -z 2 )(1 -v 2 ) 2Ξ + z ẑ -v v -Ξ ≥ M 2 (1 -ẑ2 )(1 -v 2 ) + (1 -z 2 )(1 -v 2 ) 2Ξ -( v 2 -z 2 )( v 2 -ẑ2 ) -Ξ = M 4Ξ ( v 2 -ẑ2 )(1 -v 2 ) -( v 2 -z 2 )(1 -v 2 ) 2 .
For η = 1, we have

L F (ρ, ρ) = M 2 (1 -ẑ2 )(1 -v v -Ξ) = M (1 -ẑ2 ) 1 -F (ρ, ρ) .
Therefore, for all η ∈ [0, 1] and (ρ, ρ) ∈ int(S 2 ) × int(S 2 ), we have L F(ρ, ρ) ≥ 0 which implies that L V(ρ, ρ) ≤ 0. By the stochastic LaSalle-type theorem in [START_REF] Mao | Stochastic versions of the LaSalle theorem[END_REF], we deduce that lim t→∞ L F(ρ t , ρt ) = 0 almost surely. Since L F(ρ, ρ) for any η ∈ (0, 1] can be written as a convex combination of the expressions of L F(ρ, ρ) for η = 0 and η = 1 given above, we have that either |ẑ| converges to one or F(ρ, ρ) converges to one almost surely. This concludes the proof of the first part of the proposition. The additional assumption u(ρ e )u(ρ g ) = 0 rules out the first possibility, completing the proof of the proposition.

Remark 3.4. The fact that the fidelity is a sub-martingale in int(S 2 ) × int(S 2 ), as shown in the above proof, suggests that the convergence of F (ρ t , ρt ) to one should not require the technical condition u(ρ e )u(ρ g ) = 0. This will be discussed more in detail in Section 4.

Simulations

In this section, we illustrate Proposition 3.3 through simulations of the system ( 1)-( 2) in the case u(ρ) ≡ 1 and with parameters ω = 0.3, η = 0.3 and M = 1. We set ρ e as the initial state of the quantum filter and ρ g as the initial state of the quantum filter estimate. In Fig. 1, we simulate the fidelity F(ρ t , ρt ). Fig. 2 represents the behavior of a sample trajectory ρ t and its corresponding estimation ρt in Bloch sphere coordinates. Fig. 2. Behavior of the trajectories ρt and ρt with the feedback law u(ρ) = ẑ2 , starting at (ρ e , ρ g ), when ω = 0.3, η = 0.3 and M = 1. The black curve corresponds to a sample trajectory ρt, the blue point represents its end point; the red curve corresponds to the estimated trajectory ρt, the magenta point represents its end point.

ASYMPTOTICS AND FEEDBACK CONTROL OF THE COUPLED SPIN-J SYSTEMS

In this section, we discuss the asymptotic behavior of the true quantum filter and its estimate for spin-J systems with unknown initial states. The stochastic master equations are given by

dρ t = F ût (ρ t )dt+L(ρ t )dt+G(ρ t ) dY t -2 ηM Tr(J z ρ t )dt , dρ t = F ût (ρ t )dt+L(ρ t )dt+G(ρ t ) dY t -2 ηM Tr(J z ρt )dt ,
where

• the quantum state of the spin system is denoted as ρ, and belongs to the space

S N := {ρ ∈ C N ×N | ρ = ρ * , Tr(ρ) = 1, ρ ≥ 0}. The estimated state is denoted as ρ ∈ S N , • F u (ρ) := -i[ωJ z +u J y , ρ], L(ρ) := M 2 (2J z ρJ z -J 2 z ρ- ρJ 2 z ) and G(ρ) := √ ηM J z ρ + ρJ z -2Tr(J z ρ)ρ , • Y t
denotes the observation process of the quantum spin system, which is a continuous semi-martingale whose quadratic variation is given by [Y, Y ] t = t. Its dynamics satisfies dY t = dW t + 2 √ ηM Tr(J z ρ t )dt, where W t is a one-dimensional standard Wiener process,

• ût := u(ρ t ) denotes the feedback controller as a function of the estimated state ρt , • J z is the (self-adjoint) angular momentum along the z axis, and it is defined by J z e n = (J -n)e n , n ∈ {0, . . . , 2J}, where J = N -1 2 represents the fixed angular momentum and {e 0 , . . . , e 2J } corresponds to an orthonormal basis of C N . With respect to this basis, the matrix form of J z is given by

J z =       J J -1 . . . -J + 1 -J       ,
We define the pure states ρ n := e n e * n for n ∈ {0, . . . , 2J} corresponding to the eigenvectors of J z .

• J y is the (self-adjoint) angular momentum along the y axis, and it is defined by

J y e n = -ic n e n-1 + ic n+1 e n+1 , n ∈ {0, . . . , 2J},
where

c m = 1 2 (2J + 1 -m)m.
The matrix form of J y is given by

J y =       0 -ic 1 ic 1 0 -ic 2 . . . . . . . . . ic 2J-1 0 -ic 2J ic 2J 0       ,
• η ∈ (0, 1] measures the efficiency of the detectors, M > 0 is the strength of the interaction between the system and the probe, and ω ≥ 0 is a parameter characterizing the free Hamiltonian.

By replacing dY t = dW t + 2 √ ηM Tr(J z ρ t )dt in the above equation, we obtain the following matrix-valued stochastic differential equations describing the time evolution of the pair (ρ t , ρt ) ∈ S N × S N ,

dρ t = F ût (ρ t )dt + L(ρ t )dt + G(ρ t )dW t , (7)
dρ t = F ût (ρ t )dt + L(ρ t )dt + 2 ηM G(ρ t )Tr J z (ρ t -ρt ) dt + G(ρ t )dW t . (8) If u ∈ C 1 (S N , R
), the existence and uniqueness of solutions of ( 7)-( 8) can be shown by similar arguments as in (Mirrahimi and van Handel, 2007, Proposition 3.5).

Note that, if we turn off the feedback controller, there are N 2 equilibria (ρ n , ρ m ) with n, m ∈ {0, . . . , 2J} for the coupled system. However, since the system (7) satisfies the non-demolition condition and the measurement operator J z satisfies the non-degeneracy condition [START_REF] Benoist | Large time behavior and convergence rate for quantum filters under standard non demolition conditions[END_REF], Definition 2 and Assumption (ND)), we may state the following result. Theorem 4.1. [START_REF] Benoist | Large time behavior and convergence rate for quantum filters under standard non demolition conditions[END_REF]). If u ≡ 0 and η = 1, (ρ t , ρt ) converges exponentially towards the set {(ρ 0 , ρ 0 ), . . . , (ρ 2J , ρ 2J )}.

Furthermore, the following holds true. Lemma 4.2. Consider the coupled system ( 7)-( 8) with η ∈ (0, 1] and u ∈ C 1 (S N , R). Then the fidelity F(ρ, ρ) is a sub-martingale. Indeed, note that the value L F(ρ, ρ) only depends on the instantaneous value of the Hamiltonian ωJ z +u(ρ)J y . Then by (Amini et al., 2014, Theorem 5), one has L F(ρ, ρ) ≥ 0 for the constant Hamiltonian instantaneously equal to ωJ z + u(ρ)J y , which proves the result. Remark 4.3. In Proposition 3.3, by explicit computations, we have shown that for spin-1 2 systems L F(ρ, ρ) ≥ 0, which is consistent with Lemma 4.2. This implied that L F(ρ, ρ) should converge to zero almost surely. Moreover, we were able to analyze the zeros of L F(ρ, ρ) and deduce some information about the asymptotic behavior of F(ρ, ρ). Similar results for spin-J systems cannot be obtained easily because of the complicated form of L F(ρ, ρ).

The above results give some intuitions about the asymptotic behaviors of ρ and ρ at least for the case η = 1 or for spin-1 2 systems. In these cases, without loss of generality, one may assume that the fidelity F(ρ 0 , ρ0 ) is close to one. If η = 1, this can be obtained for example by turning off the feedback controller for a large enough time in view of Theorem 4.1. For spin-1 2 systems, one can exploit Proposition 3.3 and apply a feedback controller satisfying u(ρ g )u(ρ e ) = 0 for a large enough time. Now, by Lemma 4.2, we have E(F(ρ t , ρt )) ≥ F(ρ 0 , ρ0 ) which leads to the conclusion that E(F(ρ t , ρt )) is close to one, independently of the chosen feedback law.

In the following, we discuss the possibility of designing a feedback controller u(ρ) which stabilizes exponentially almost surely the coupled system (7)-( 8) towards a given target state (ρ n, ρ n) with n ∈ {0, . . . , 2J}. We note that the two subsystems ( 7)-( 8) share the same feedback controller u(ρ), only depending on the estimated state ρ. Hence, if we suppose that the feedback u satisfies the assumption u(ρ n) = 0 and u(ρ k ) = 0 for all k = n, then the coupled system ( 7)-( 8) possesses exactly N equilibria, given by (ρ k , ρ n) for k ∈ {0, . . . , 2J}.

Concerning feedback exponential stabilization of the coupled spin-J systems ( 7)-( 8) towards (ρ n, ρ n), we propose a conjecture inspired by (Liang et al., 2019a, Theorem 6.4 and Theorem 6.5). These results were developed for the case ρ 0 = ρ0 . Conjecture 4.4. Consider the coupled system ( 7)-( 8) with (ρ 0 , ρ0 ) ∈ S N × S N \ (ρ n, ρ n) and assume η ∈ (0, 1). Then, the feedback controller u n( ρ) = α 1 -Tr(ρρ n) β , α > 0, β ≥ 1, (9) exponentially stabilizes (ρ t , ρt ) to (ρ n, ρ n) almost surely for the special case n ∈ {0, 2J} with sample Lyapunov exponent less or equal than -ηM . Moreover, the feedback controller u n( ρ) = α J -n -Tr(J z ρ) β , α > 0, β ≥ 1, (10) exponentially stabilizes (ρ t , ρt ) to (ρ n, ρ n) almost surely for the general case n ∈ {0, . . . , 2J} with sample Lyapunov exponent less or equal than -ηM/2 for n ∈ {1, . . . , 2J -1} and -ηM for n ∈ {0, 2J}.

Simulations

Here, we illustrate the above conjecture through simulations for a three-level quantum spin system. In Fig. 3, we consider the target state (ρ 0 , ρ 0 ), and the Lyapunov function V 0 (ρ, ρ) = 1 -Tr(ρρ 0 )Tr(ρρ 0 ). We define the distance for the coupled system via the Bures distance d B (•, •) in S 3 [START_REF] Bengtsson | Geometry of quantum states: an introduction to quantum entanglement[END_REF]

) as follows d B ((ρ, ρ), ρ 0 ) := d B (ρ, ρ 0 ) + d B (ρ, ρ 0 ). We have √ 2 4 V 0 (ρ, ρ) ≤ d B ((ρ, ρ), ρ 0 ) ≤ √ 2V 0 (ρ, ρ). In

CONCLUSION

In this paper, for spin-1 2 systems, we have shown that under appropriate assumption on the feedback controller, we can guarantee the same asymptotics for the quantum filter and its estimate. For the general case of spin-J systems, we provided an heuristic approach regarding the asymptotic behavior of ρ and ρ and feedback exponential stabilization of the coupled system (7)-( 8). Future goals will be to Also, exponential stabilization of entangled states such as GHZ states, with unknown initial states is included in our research lines (see [START_REF] Liang | On exponential stabilization of two-qubit systems[END_REF] for exponential stabilization of two-qubit systems with perfect initialization).

Fig. 1 .

 1 Fig. 1. Convergence of the fidelity F (ρt, ρt) towards one with the feedback law u(ρ) = ẑ2 starting at (ρ 0 , ρ0 ) = (ρ e , ρ g ), when ω = 0.3, η = 0.3 and M = 1: the black curve represents the mean value of 100 arbitrary samples.

  Fig 4, we choose (ρ 1 , ρ 1 ) as the target state and V 1 (ρ, ρ)= k =1 Tr(ρρ k ) + Tr(ρρ k ) as the Lyapunov function.We can show that√ 2 2 V 1 (ρ, ρ) ≤ d B ((ρ, ρ), ρ 1 ) ≤ √ 2V 1 (ρ, ρ).

Fig. 3 .

 3 Fig. 3. Exponential stabilization of a three-level quantum spin system towards (ρ 0 , ρ 0 ) with the feedback law (9) starting at (ρ 0 , ρ0 ) = (ρ 2 , ρ 1 ) with ω = 0.3, η = 0.3, M = 1, α = 5 and β = 2. The black curve represents the mean value of 100 arbitrary sample trajectories, the red and blue curves represent the exponential references with exponents -ηM/2 and -ηM respectively. The figures at the bottom are the semi-log versions of the ones at the top.

Fig. 4 .

 4 Fig. 4. Exponential stabilization of a three-level quantum spin system towards (ρ 1 , ρ 1 ) with the feedback law (10) starting at (ρ 0 , ρ0 ) = diag(0.2, 0.2, 0.6), diag(0.8, 0.1, 0.1) with ω = 0.3, η = 0.3, M = 1, α = 2 and β = 2. The black curve represents the mean value of 100 arbitrary sample trajectories, the red curve represents the exponential reference with exponent -ηM/2. The figures at the bottom are the semi-log versions of the ones at the top.

  provide a rigorous proof of Conjecture 4.4 and present systematic methods ensuring exponential stabilization of the coupled spin systems. Other research lines include generalization to more general forms of filter equations which are driven by Wiener and/or Poisson processes.

  , Lemma 6.1), we can show that if ρ 0 ∈ ∂S 2 , ρt enters in B re (ρ e )∪B rg (ρ g ) in finite time almost surely, which means that u(ρ t ) becomes non-zero in finite time almost surely. As a consequence, ρ t exits the boundary and stays in the interior of S almost surely. The proof is then complete. Proposition 3.3. Assume η ∈ (0, 1] and let u ∈ C 1 (S 2 , R).

The authors thank the support from Agence Nationale de la Recherche projects Q-COAST ANR-19-CE48-0003 and QUACO ANR-17-CE40-0007.


Inderjeet Singh

Manarshhjot Singh

Ismail Bensekrane

Othman Lakhal

Rochdi Merzouki
email: rochdi.merzouki@polytech-lille.fr

Curve-based Approach for Optimal Trajectory Planning with Optimal Energy Consumption: application to Wheeled Mobile Robots

Wheeled transportation systems have been used since the dawn of civilization. The great amount of experience we have with wheeled transportation systems gives the Wheeled Mobile Robots (WMRs) an edge over other types of robots. This has allowed WMRs to be the prime candidates for many applications including, but not limited to, military and space applications. Sometimes these vehicles have to work in complex environments with complicated constraints to finish their assigned tasks. This paper presents a methodology to generate a smooth feasible trajectory for WMR using Pythagorean Hodograph (PH) curves. Generated trajectory not only avoids the static obstacles but also conforms to minimum energy consumption and travel time possible within the kinematic constraints of the WMR. Simulations are presented to check the performance of the proposed approach. The proposed method is experimentally validated using a WMR 'Robotino'.

INTRODUCTION

Accurate trajectory planning is a crucial task for WMRs to navigate safely and efficiently from the start point to the finish point autonomously. The path undertaken by them should be smooth, short, and collision-free. The selected velocity profile on the generated path should assure the minimum energy consumption to increase the availability of the WMR.

Related Work

Path planning of wheeled mobile robots has been done using various algorithms like Dijkstra's algorithm, A* algorithm, Rapidly exploring random tree (RRT) algorithm, Probabilistic road map (PRM) algorithm, Potential fields algorithm etc.. [START_REF] Yang | Survey of robot 3d path planning algorithms[END_REF]. Most of the methods generate a path consisting of straight lines with sharp turns. Therefore, the path is unnatural for wheeled robots to follow as they can not take sudden sharp turns. To avoid this problem, curve based approaches are used to make a smooth path. [START_REF] Ravankar | Path smoothing techniques in robot navigation: State-of-the-art, current and future challenges[END_REF] gives a review of different curve based approaches: Bezier curves, cubic splines, B-spline, NURBS, etc. In 1990, Pythagorean Hodograph (PH) curves were introduced by Farouki and Sakkalis [START_REF] Farouki | Pythagorean hodographs[END_REF] to overcome the various drawbacks of the above-mentioned curves. Therefore, in [START_REF] Bruyninckx | Path planning for mobile and hyper-redundant robots using pythagorean hodograph curves[END_REF], PH curves are used for path planning of mobile robots. Most of the PH based path planning has been done by varying the curvature of the PH curve to avoid the obstacles. However, the curvature variation leads to a significant increase in the path length. Therefore, [START_REF] Shah | Obstacles avoidance in the frame work of pythagorean hodograph based path planning[END_REF] proposed a methodology to avoid obstacles using PH curves by taking a different combination of direction vectors at the start and endpoint. The method considered only some combinations of direction vectors from four different quadrants. This paper overcomes this limitation and proposes an approach to optimize the direction vectors at start and destination (end) points to create an obstaclefree path. Further, motion planning of WMR is needed on the planned path so that it should consume minimum energy to finish its task in minimum travel time. For the energy consumption of autonomous vehicles, most of the works use the A* algorithm to generate a path on which minimum energy consumption of the vehicle is considered [START_REF] Datouo | Optimal motion planning for minimizing energy consumption of wheeled mobile robots[END_REF]. In [START_REF] Ondrúška | Scheduled perception for energyefficient path following[END_REF], the minimum energy consumption technique is presented for a known path. In this work, parameters of the PH curve are optimized to generate a trajectory with minimum energy consumption in minimum travel time. In the present work, we propose an optimization approach:

• to optimize the direction vectors at the start and endpoints of the PH curve to avoid known obstacles. • to generate a trajectory using PH curves by consideration of minimum energy consumption and travel time of WMR.

Outline

The paper is organized as follows: Section 2 discusses the formulation of a quintic PH curve for path planning. The energy consumption of wheeled robots is discussed in section 3. Optimization methodology to generate a smooth and collision-free path has been discussed in section 4. Section 5 presents the discretization of a path into more than one PH curve. Section 6 provides simulation as well as experimental results for the proposed methodology. Finally, section 7 provides the conclusions and future perspectives of the work.

PYTHAGOREAN HODOGRAPH CURVES FOR PATH PLANNING

Quintic PH curves are used for path planning because they are a more appropriate compromise between complexity and smoothness. Formulation of the planar PH-curve used for path planning, originating from an initial position to the final (destination) position is discussed using normal polynomial parametric curves as (Singh et al., 2018b,a): Let r(h) be the polynomial curve, r(h) = (x(h), y(h)); 0 ≤ h ≤ 1 (1) where h is the curvilinear coordinate of the curve. It means start and final points of the curve are:

P s = r(0) = (x(0), y(0)) at h = 0 P f = r(1) = (x(1), y(1)) at h = 1 (2)
The curve joining derivatives of the polynomial curve is defined as the 'hodograph' of the curve. This is parallel to the tangent to the curve and is given as follows: r (h) = (x (h), y (h))

(3) Let L(h) represents the length of the path constructed using the curve r(h):

L(h) = 1 0 |r (h)| dh = 1 0 x (h) 2 + y (h) 2 dh (4)
From eq. ( 4), the square root sign can be eliminated, if:

x (h) 2 + y (h) 2 = σ(h) 2
(5) Eq. ( 5) is known as Pythagorean law or condition. Therefore, if the first derivative (hodograph) of parametric polynomials, satisfies the pythagorean condition at all the points, the curve is known as Pythagorean Hodographs. Now eq. ( 4) can be written as:

L(h) = 1 0 |σ(h)| dh (6)
This also allows the closed form solution of the length of the path.

For the construction of a PH curve r with minimum bending energy, four boundary conditions are required [START_REF] Singh | Modeling of continuum manipulators using pythagorean hodograph curves[END_REF]):

• P s (x s , y s ), start point of the curve.

• d s (θ s ), unit direction vector at the start point P s .

• P f (x f , y f ), end point of the curve.

• d f (θ f ), unit direction vector at the end point P f .

Using these four boundary conditions, six control points of a quintic PH curve are computed to generate a curve with minimum bending and is explained in [START_REF] Singh | Modeling of continuum manipulators using pythagorean hodograph curves[END_REF]. A basic quintic PH is shown in Fig. 1(a) with its four boundary conditions and six control points from P 0 to P 5 . Four boundary conditions (P s , P f , d s , d f ) fix the four control points of a quintic PH curve and there are two free control points (Fig. 1(a)). Due to these two free control points there are multiple solutions to generate quintic PH curve for four fixed boundary conditions. Formulation of a quintic PH uses [START_REF] Singh | Modeling of continuum manipulators using pythagorean hodograph curves[END_REF][START_REF] Farouki | Pythagorean hodographs[END_REF]) a condition to generate a unique curve with minimum bending energy for the above mentioned four boundary conditions. Hence, a direct solution for minimum energy curve is possible instead of an iterative solution.

The general equation of quintic PH curve r(h) with control points is as follows:

r(h) =        1 h h 2 h 3 h 4 h 5        T        1 0 0 0 0 0 -5 5 0 0 0 0 10 -20 10 0 0 0 -10 30 -30 10 0 0 5 -20 30 -20 5 0 -1 5 -10 10 -5 1               P 0 P 1 P 2 P 3 P 4 P 5        (7)
The curvature of curve r(h) is given as:

κ = |r × r | |r | 3 (8)
Start (P s ) and end (P f ) points of a path are fixed but direction vectors (d s (θ s ) and d f (θ f )) at these points can vary. Figures 1(b),1(c) shows the effect of direction vectors on the path. Also, Fig. 1(d) illustrates that by variation of direction vectors, different possible paths can be achieved to avoid obstacles. Section 4 explains the method to select optimal direction vectors according to our requirements.

ENERGY CONSUMPTION OF WHEELED ROBOTS

The total energy consumption of a WMR is a sum of the total mechanical (kinetic + potential) energy of the WMR and the various losses. For the current formulation, the path is assumed as planar, hence the potential energy is unchanged. The losses associated with WMR include the frictional losses at the wheel and electrical losses. V and V are the linear velocity and acceleration of the WMR. ω and ω represent the angular velocity and angular acceleration of the WMR while steering. P s represents the static power consumption of the embedded computer and different sensors. g is the standard acceleration due to gravity.

OPTIMIZATION STRATEGY FOR WMR

The optimization problem in this paper is recognized as a bounded, constrained, multi-variable, multi-objective, feasibility study. Hence, a multi-objective Genetic Algorithm is used to solve the problem.

Objective Functions and Constraints

The objective functions to be minimized for the trajectory generation are generally expressions containing significant physical parameters related to the WMR's behavior and to the efficiency of the system. Hereafter, we particularly consider the following functions:

F 1 = E T otal , F 2 = T (
10) where, T is the travel time, E T otal represents the total energy consumption. As discussed earlier PH based path generation depends only on start and end direction vectors: r = f 1 (θ s , θ f ) (11) From eq. 9, energy consumption of a WMR depends on the linear and angular velocity:

E T otal = f 2 (V, ω) (12
) but for a selected path, angular velocity of the WMR is the function of linear velocity and the path.

ω = f 3 (r, V ) (13
) Therefore, from eqs. 11-13, the energy consumption depends on start and end direction vectors (θ s , θ f ) of the path (r) and the velocity profile (V ) as:

E T otal = f 4 (V, θ s , θ f ) (14
) Therefore, to minimize energy consumption, both the optimal direction vectors (θ s , θ f ) and optimal velocity profile (V ) is required for the generated path. Physical limitations: The search space for the variables (start and end direction vectors) is not continuous because some combinations of variables give a path that passes through the obstacles. Therefore these solutions are considered infeasible. The solution from the optimization Although due to the cyclic nature of angles there are no bounds on the direction vectors, but, the search space for direction vectors is bounded to one cycle.

0 • ≤ θ s ≤ 360 • 0 • ≤ θ f ≤ 360 • (16)
The velocity profile and acceleration evaluated at the center of gravity are bounded by the dynamic limitations of the robot and also the path:

V min ≤ V ≤ V max Vmin ≤ V ≤ Vmax (17)
In this work, Genetic Algorithm (GA) is used to find the global optimal solution.

Methodology

The method includes the calculation of an optimal path and an optimal velocity profile for the path. This is achieved using the GA algorithm, to find multiple Paretooptimal solutions in a multi-objective optimization problem. From the generated Pareto front, a solution (objective point) is selected (Fig. 4) according to our requirement regarding energy consumption and travel time.

Pareto front generation:

The procedure for generating optimal Pareto front for WMR path and velocity profiles is described in Fig. 2. The algorithm aims to find the solution of the path and velocity profiles ensuring the compromise between the travel time T and the energy consumption E T otal . A PH path is created using the information of direction vectors as discussed in section 2. Further, a velocity profile is generated using a cubic spline interpolating a set of intermediate control points uniformly distributed along the time scale [0, T ], as shown in Fig. 3. The positions of these control points are investigated within the limits imposed by kinematic constraints in eq. 17.

Compromise point selection: The optimization process provides Pareto front. A compromise point is selected from Pareto front in regards of the ideal point whose coordinates are minimum energy consumption and minimum travel time [min(T ), min(E)]. The selected objective point corresponds to the nearest point to the ideal point as shown in Fig. 4. As the range of energy consumption and travel time is different and can not be compared, therefore both Figure 2 represents the block diagram of optimization algorithm. This depicts that, using GA, direction vectors (θ s , θ f ) and velocity profile (V ) are optimized to generate a trajectory considering minimum energy consumption with minimum travel time by WMR.

! V(t) t V max V 1 V 2 V 3 V 4 V 5 V 6 V 7 V m 0 t

DISCRETIZATION OF PATH

The PH generation method discussed in Section 2 gives a PH curve of order 5. This is a major limitation because, for a particular set of obstacles, a fifth-order PH might not be flexible enough to generate a feasible path. In such a condition, it is proposed that the curve is disjointed to form multiple curves. To achieve this, intermediate points are created. An intermediate point is initiated at the endpoint of the path and shifted using the optimization technique. The need can be explained using Fig. 5. In the figure, Path 1 is generated using a single PH curve but Path 2 which requires a sharp turn is a combination of two PH curves joined at the Intermediate point.

Hence, the problem becomes a bi-level optimization problem. In this problem, a search algorithm for optimal path and the velocity profile is embedded inside a search algorithm for the optimum location of an intermediate point.

The structure of the optimization problem is discussed in Fig. 6.

RESULTS AND DISCUSSIONS

The presented methodology is applied on mobile robot Robotino. The robot is powered by two 12V batteries capable for powering the robot for upto two hours. The mechanical specifications of the robot are given in table 1. An environment is created with obstacles, shown in Fig. 7 to plan a planar trajectory of the Robotino from start to the destination point. The algorithm provides the trajectory of the center of mass of the Robotino. This trajectory must always maintain a specific distance from any obstacle. In current case, this distance is equal to the radius of Robotino. Hence, a clearance equal to this radius is provided around the obstacles represented by dotted lines. The available work environment after consideration of clearance is provided in Fig. 7. The bounds for Robotino are as follows:

0.01 ≤ V ≤ 1.325(ms -1 ) 0 • ≤ θ s ≤ 360 • 0 • ≤ θ f ≤ 360 • (18)
Also, the linear acceleration constraint for robotino is given as: -0.05 ≤ V ≤ 1.02(ms -2 ) (19) A trajectory is planned from start to the destination point using the proposed approach (Fig. 7). This simulated trajectory ensures the minimum energy consumption of the Robotino in the minimum travel time. Here, trajectory is made up of two PH curves. Also, as per discritization algorithm, continuity between the two PH curves is ensured. The final values of the optimization variables are: θ s = 2.

The same environment is created as shown in Fig. 8 for experimental validation. In this environment, a trajectory is planned using the proposed PH based approach as well as the widely used A-star approach. For the Astar approach, the optimal velocity profile is calculated using the same method as used for PH based approach discussed in Section 4.2. Figure 9 shows both of the planned trajectories in the same environment. Robotino is used to follow both of these trajectories. The energy consumption of Robotino is computed using its battery current and voltage. The comparison of data is given in table 2. Figure 10 shows the comparison of the optimal velocity profile as well as power consumption for the planned trajectory for the presented scenario (Fig. 7).

For the PH based path the high energy consumption is observed during high velocities, whereas for A-star based path, high energy consumption is observed during low velocity. It can be interpreted that for the PH curve most of the energy is consumed for the propulsion of the robot along the path but most of the energy for A-star trajectory is consumed when the robot turns to shift from one portion of the curve to the other. This shows that the vehicle is better able to utilize the available energy while following the PH based path. Therefore, for straight line and simple paths with negligible steering requirements Astar trajectory is expected to be energy efficient but for a complex path PH based trajectory is expected to be more energy-efficient. From the results, it is observed that the energy consumption and the travel time is significantly less when robot moves on the path generated using the proposed approach in-spite of a slight increase in path length.

As discussed in section 4.2, one solution is selected from pareto front considering an equal compromise between minimum energy and minimum travel time. However, according to the need, preference can be given to a point for either a lower energy consumption or a lower travel time.

CONCLUSIONS

This paper presents a curve-based approach to plan a smooth and power-efficient trajectory for a WMR. On account of the possibility of generating a higher-order curve, using lesser boundary conditions and availability of direct condition for minimum bending energy, PH curves are chosen to generate a path. The energy consumption and travel time of the robot are minimized by optimizing the PH based path and velocity profile. The paper also discusses the limitations of PH based path in case of using a single PH and a solution is proposed to overcome this limitation by including a serial of PH curves. Experiments are performed using Robotino to check and validate the performance of the proposed trajectory planning methodology. Our future work includes adapting this curve-based approach to the case of a dynamic environment.

  Fig. 1. Quintic PH curve (a) Basic quintic PH, (b) Effect of direction vectors on length of PH, (c) Effect of direction vectors on curvature of PH and (d) Role of direction vectors in obstacle avoidance

Fig. 2 .

 2 Fig. 2. Generation of Pareto Front in GA problem must be feasible. This is modeled as a constraint given below: {r ∩ O} = {φ} (15) where r is the PH based path and O is the set of obstacles.

FigFig. 4 .

 4 Fig. 3. Generation of optimal velocity profile

  Initially, a new intermediate point is initialized. PH curves are generated, first between the start and intermediate point and second between the intermediate and endpoint. The individual optimization function is evaluated for each curve as discussed previously in section 4.1. These objective functions are used to find optimal direction vectors and a velocity profile. The objective function to find the optimum location of the intermediate point is evaluated by adding individual objective functions of each curve. A constraint is added to assure the smoothness of the path. This is done by equating the direction vectors at the intermediate point for both the PH. The algorithm tries to find the best combination of direction vectors i.e. at the start, intermediate and endpoint. If a feasible combination is not found then the search algorithm for intermediate point continues. If the search for optimum intermediate point location exhausts without a feasible solution, another intermediate point must be added, thereby increasing the number of curves. The solution achieved before the introduction of a new curve is considered as the best possible solution using lesser number of curves. Therefore, this acts as an initial solution for the search for a new intermediate point problem. The new intermediate point must be placed to split the worst curve i.e. with the highest individual value of the objective function. The optimum position of the new intermediate point and the associated direction vectors are found. The process is continued and new intermediate points are added till a collision-free path is created.

Fig. 7 .

 7 Fig. 5. Need for discretization of path

FigFig. 9 .

 9 Fig. 8. Experimental set-up

Fig. 6 .

 6 Fig. 6. Block diagram of algorithm of discretization of path

Fig. 10 .

 10 Fig. 10. Comparison between PH and A-Star based trajectories: (a) Optimal velocity profile, (b) Power consumption

Table 1 .

 1 Robotino Parameters Moment of inertia 0.16245 Kgm -2

	Parameter	Value
	Robot radius	175 mm
	Height	210 mm
	Battery voltage	24 V
	Robot mass	11 Kg
	Max velocity	1.325 ms -1

Table 2 .

 2 Comparison of proposed approach with simulation and A-star algorithm

	Trajectory	Path Length (mm) Energy Consumption (J) Travel Time (sec)
	PH based (Experiment)	3612.44	136.81	23.89
	PH based (Simulation)	3592.10	133.01	22.40
	A-star based (Experiment)	3578.10	161.67	32.65


L Etienne
email: lucien.etienne@imt-lille-douai.fr

K.M.D Motchon

E Duviella
email: eric.duviella@imt-lille-douai.fr

K Guelton
email: kevin.guelton@univ-reims.fr

Impulsive observer design for switched linear systems with time varying sampling and (a)synchronous switching rules

Keywords: Continuous-discrete observers, linear switched systems, sampled-data systems, LMIs

Observer synthesis for switched linear systems with time-varying sampled-data measurements is addressed in this paper. Between sampling times, the behaviour of the observer under consideration has the particularity to be that of the estimated active mode while, at the sampling times a correction (or an update) is made using the estimated active mode and the sampled output signal. Two cases are discussed: i) the situation where no delay occurs in the estimation of the active mode and, ii) the more practical case where the estimated switching signal that defines the proposed observer is a delayed replica of the plant's switching signal. In these two situations, appropriate time-varying piecewise quadratic Lyapunov functions are used to establish convergence conditions of the proposed observer in the linear matrix inequalities (LMIs) framework.

INTRODUCTION

Switched systems has been the subject of a growing amount of interests in the past decade and are by now a well studied topic [START_REF] Daafouz | Stability analysis and control synthesis for switched systems: a switched lyapunov function approach[END_REF][START_REF] Liberzon | Basic problems in stability and design of switched systems[END_REF][START_REF] Liberzon | Switching in systems and control[END_REF]. Switched systems are a subclass of hybrid systems where different dynamics describe the different operating modes. The transition from one mode to another is given by a switching rule. The motivation for considering such a class of systems is vast and encompasses robust stabilization (modeling and control) of systems subject to sampling uncertainties and singularly perturbed systems. An important problem arising when considering switched systems is the stability problem. Indeed given a switching rule, stability (resp. instability) of all modes does not guaranty the stability (resp. instability) of the overall switched system [START_REF] Liberzon | Switching in systems and control[END_REF]. Several results exist on stability of switched systems in different settings (for a survey on the topic, see for instance [START_REF] Briat | Convex conditions for robust stability analysis and stabilization of linear aperiodic impulsive and sampled-data systems under dwell-time constraints[END_REF][START_REF] Allerhand | Robust stability and stabilization of linear switched systems with dwell time[END_REF][START_REF] Zhang | Asynchronously switched control of switched linear systems with average dwell time[END_REF][START_REF] Xiang | Necessary and sufficient condition for stability of switched uncertain linear systems under dwell-time constraint[END_REF][START_REF] Lin | Stability and stabilizability of switched linear systems: a survey of recent results[END_REF].

While full state feedback is often assumed to be measured, in many practical applications due to technical or economic constraints, such an hypothesis is not justified. Therefore observation for linear and nonlinear systems constitutes a large area of study. For this reason observer synthesis for switched systems has received some attention in the literature [START_REF] Alessandri | Design of observers for switched discrete-time linear systems[END_REF][START_REF] Pettersson | Designing switched observers for switched systems using multiple lyapunov functions and dwell-time switching[END_REF][START_REF] Barbot | Nonlinear observer for autonomous switching systems with jumps[END_REF][START_REF] Bejarano | Switched observers for switched linear systems with unknown inputs[END_REF]. However when the output is transmitted over a network the assumption of continuous transmission is no longer warranted and close attention has to be paid to the aperiodic and sampled nature of the output. This consideration has, in the last decade, led to numerous works in the field of observation with aperiodic sampling for both continuous time linear and nonlinear systems. For such sampled data systems, the so called impulsive observer (sometimes referred as continuous-discrete observer) offer interesting properties in terms of both performance and simplicity of analysis (see for instance [START_REF] Raff | Observers with impulsive dynamical behavior for linear and nonlinear continuoustime systems[END_REF][START_REF] Dinh | Continuous-discrete time observer design for lipschitz systems with sampled measurements[END_REF][START_REF] Etienne | Observer synthesis under time-varying sampling for lipschitz nonlinear systems[END_REF]. However, to the best of the author's knowledge, observation of continuous time switched systems subject to time varying sampling has not received attention. While some results have recently been proposed for discrete time systems [START_REF] Han | Asynchronous observer design for switched linear systems: a tube-based approach[END_REF], such missing analysis motivates the present study.

In this work we will address the problem of observer synthesis for switched linear systems with aperiodic sampled measurements. Between the sampling times, the observer is just a copy of the dynamics given by the estimated active mode. At sampling times, an instantaneous correction (or an update) is made using the sampled output signal as well as the estimated active mode. For the design of this observer, two cases are considered:

• The first case assumes that the switching sequence is known in real time. • For the second case, since in practice the identification of the active mode may takes some times thus introducing a delay [START_REF] Tian | Algebraic switching time identification for a class of linear hybrid systems[END_REF][START_REF] Yang | Active mode identification and continuous state estimation for switched linear systems with unknown inputs and slow switching signal[END_REF], we will also consider the asynchronous situation where the switching times are not instan-taneously estimated. In this context, we will assume that there can be a bounded delay between the real active mode and its estimation provided by the observer. Indeed if one considers sampled data outputs, then it is natural to see the switching signal (which can be seen as an output) as not continuously measured and available in real-time .

In these two situations, using appropriate time-varying piecewise quadratic Lyapunov functions, we propose convergence conditions for the impulsive observers in the linear matrix inequalities framework.

The rest of the paper is organized as follows. The considered class of switched linear systems and the impulsive observers under consideration are described in Section 2. LMIs-based conditions proposed for the observer synthesis, obtained from piecewise quadratic Lyapunov functions, are presented in Section 3. To show the relevance of the proposed observer synthesis methods, numerical examples are also provided in Section 4. Concluding remarks are finally given in Section 5.

Notation: For a vector or a matrix v, v ⊤ denotes its transpose. We define for a matrix A, He(A) := A + A ⊤ . R ≥0 corresponds to the set of non-negative real numbers. S n denotes the set of n × n symmetric matrices while S + n denotes the set of n × n positive definite symmetric matrices. For a matrix M ∈ S n , λ min (M ) and λ max (M ) denote its smallest and largest eigenvalues, respectively. For two positive definite (resp. positive semidefinite) matrices P and Q, we write P ≻ Q (resp. P ≽ Q) if P -Q is positive definite (resp. positive semidefinite).

DEFINITIONS AND PROBLEM STATEMENT

Let's consider the following class of switched linear systems:

   ẋ(t) = A σ(t) x(t) + B σ(t) u(t), t ̸ = s k , t ̸ = t k , x(s + k ) = x(s k ), y(t k ) = C σ(t k ) x(t k ), (1)
where ∀t ∈ [0 ; ∞), x (t) ∈ R n , u (t) ∈ R ℓ and y (t) ∈ R m are the state vector, the control input and the output vector. The signal σ : [0 ; ∞) -→ S ⊂ N is a piecewise constant function, continuous from the right, representing the switching signal of the system. {s k } k∈N is the strictly increasing sequence of the switching instants and {t k } k∈N is the sequence of sampling times. By convention, we set s 0 = t 0 = 0 and for technical simplicity we assume that ∀k ̸ = 0, {t k } k∈N ∩ {s k } k∈N = ∅. For every i ∈ S , A i , B i and C i are known matrices with appropriate dimensions. Solution of the dynamical system are defined iteratively.

The switching signal is assumed to admit an average dwelltime (ADT), i.e.: Assumption 1. There exists s ∈ R ⋆ + such that for every s, t ∈ [0 ; ∞) with t ≥ s, the number N s (s, t) of switching of σ in the interval [s ; t[ verifies N s (s, t) ≤ N 0 + t-s s for some N 0 ∈ N; s is the ADT and N 0 is called the chatter bound.

Furthermore, sampling sequences with ranged dwell are considered, that is: Assumption 2. There exists constants τ ≥ τ > 0 such that ∀k ∈ N, t k+1 -t k ∈ [τ ; τ ]. Remark 3. Note that for sampling sequences {t k } k∈N satisfying the ranged dwell hypothesis in Assumption 2, for every s, t ∈ [0 ; ∞) with t ≥ s, the number N c (s, t) of sampling t k in [s ; t) verifies N c (s, t) ≥ t-s τ -1. For the class of switched linear systems of the form (4), the concept of global uniform exponential stability (GUES) needed in our subsequent analysis is recalled. Definition 4. (GUES). System (1) is said to be globally uniformly exponentially stable with convergence rate β and overshoot

M if, ∀t ∈ [0 ; ∞) and ∀x 0 = x(0) ∈ R n , ∥x(t)∥ ≤ M exp (-β t) ∥x 0 ∥.
The impulsive observer under consideration is described by:

   ẋ(t) = A σ(t) x(t) + B σ(t) u(t), t ̸ = s k , t ̸ = t k , x(s + k ) = x(s k ), x(t + k ) = x(t k ) + G σ(t k ) ( C σ(t k ) x(t k ) -C σ(t k ) x(t k ) ) , (2)
where ∀t ∈ [0 ; ∞), x (t) ∈ R n is the observer state vector (i.e. the estimated state vector) and σ : [0 ; ∞) -→ S is an estimate of the real switching signal σ. The switching times of σ are denoted by the strictly increasing sequence {ŝ k } k∈N . It is assumed next that the estimated switching times are obtained with a maximum time delay δ ≥ 0, that is:

Assumption 5. ∀k ∈ N, ŝk -s k ≤ δ.
Let us stress that the identification method of the real active mode is outside the scope of this paper, but our work can be used for any method that can reconstruct the real state in an time upper bounded by δ.

Note that, according to (2), between sampling, the observer is a copy of the system that run the estimated active mode σ. Moreover, when a sampling occurs, an instantaneous correction is made using the sampled output and the estimated active mode σ (t).

A simpler specific case occurs when no delay in the active mode detection is present (i.e. δ = 0). In this case, the observer (2) yields to:

   ẋ(t) = A σ(t) x(t) + B σ(t) u(t), t ̸ = s k , t ̸ = t k x(s + k ) = x(s k ), x(t + k ) = x(t k ) + G σ(t k ) C σ(t k ) ( x(t k ) -x(t k )
) .

(3)

In this delay-free situation, we are interested by the following observer synthesis problem (OSP): Problem 6. Solving the observer synthesis problem (OSP) is to find a set of impulsive observation gains G i , i ∈ S for observer (3) such that the estimated error e = x -x is GUES.

When the estimated switching signal is the delayed replica of the system's switching signal, the following observer based stability analysis problem (OBSAP) will be examined: Problem 7. Solving the observer based stability analysis problem (OBSAP) is to find conditions such that (x ⊤ , x ⊤ -x⊤ ) ⊤ is GUES when both matrices K i and G i , i ∈ S are given.

MAIN RESULT

Case of synchronous switching signals

In this part, we consider the case where δ = 0 (i.e. σ and σ are synchronous). Considering x given by (1) and x defined by (3), the observation error has the following dynamics:

     ė(t) = A σ(t) e(t), t ̸ = s k , e(s + k ) = e(s k ), e(t + k ) = (I -G σ(t k ) C σ(t k ) )e(t - k ). (4a)
∀i ∈ S , ∀τ ∈ {0, τ }, P i (τ ) ≻ 0, (5) ∀i ∈ S , ∀τ ∈ {0, τ }, He(P i (τ ) A i ) + P i 1 + αP i (τ ) ≺ 0, (6) ∀i ∈ S , ∀τ ∈ {τ , τ }, [ µ c P i (τ ) P i 0 -L i C j * P i 0 ] ≻ 0, (7) ∀i, j ∈ S , ∀τ ∈ {0, τ }, P i (τ ) ≺ µ s P j (τ ) , (8) α > ln µ s s + ln µ c τ (9)
where ∀i ∈ S , ∀τ ∈ R, P i (τ ) = P i 0 + τ P i 1 , (10) then the system (4) with impulsive observation gains G i = (P i 0 ) -1 L i , i ∈ S is GUES (i.e. the OSP is solved) with the convergence rate:

β = 1 2 ( α - ln µ s s - ln µ c τ ) . ( 11
)
Proof 1. Let's consider the piecewise continuous Lyapunov function candidate

V (t) = e ⊤ (t) P σ(t) (t -τ (t)) e (t), t ∈ [0 ; ∞) where τ (t) = max{t k : t ≥ t k }
is a piece wise constant function that denotes at time t the last sampling that occurred. P σ(t) (t -τ (t)) is introduced in (10) i.e.:

P σ(t) (t -τ (t)) = P σ(t) 0 + (t -τ (t)) P σ(t) 1 . For any t ∈ [0 ; ∞), by construction of τ (t), one has t -τ (t) ∈ [0 ; τ ]
and thus by convexity argument, one deduce from condition (5) that V is a positive definite function.

Denote by { tk } k∈N the strictly increasing sequence associated to

{t k } k∈N ∪ {s k } k∈N . Let t ∈ ] tk ; tk+1 [ with k ∈ N.
By a simple differentiation of V and using the dynamic equation (4a) of the observation error, one can verify that V (t) = e ⊤ (t) Γ σ(t) e (t), where

Γ σ(t) = [ He ( P σ(t) (t -τ (t)) A σ(t) ) + P σ(t) 1
] .

Since t -τ (t) ∈ [0 ; τ ], using condition (6) and convexity argument, one gets that Γ σ(t) ≺ -α P σ(t) (t -τ (t)) and consequently, V (τ ) ≤ -αV (τ ). Integrating this differential inequality on ] tk ; tk+1 [ yields:

V (t) ≤ exp (-α (t -tk )) V ( t+ k ) , ∀τ ∈ ] tk ; tk+1 [. (12) Furthermore, since e ( s + k ) = e ( s - k )
, using (4b) one obtains

V ( s + k ) = e ⊤ ( s - k ) P σ(s + k ) ( s k -τ ( s - k )) e ( s - k ) , ∀k ∈ N. As s k -τ ( s - k ) ∈ [0 ; τ ]
, by a convexity argument, one can deduce from condition (8) that

P σ(s + k ) ( s k -τ ( s - k )) ≺ µ s P σ(s - k ) ( s k -τ ( s - k ))
and consequently, the relation

V ( s + k ) ≤ µ s V ( s - k ) , ∀k ∈ N, (13)
holds. We will now establish a similar relation for the sampling times. Since t + k -τ

( t + k ) = 0 and e ( t + k ) = e ( t - k )
, using (4c), one can see that

V ( t + k ) = e ⊤ ( t - k ) Ω ⊤ σ(t k ) ( P σ(t k ) 0 ) -1 Ω σ(t k ) e ( t - k )
where

Ω σ(t k ) = P σ(t k ) 0 -L σ(t k ) C σ(t k ) . This implies that V ( t + k ) -µ c V ( t - k ) = -e ( t - k ) Γ σ(t k ) e (t k ) with Γ σ(t k ) = µ c P σ(t - k ) ( t k -τ ( t - k )) -Ω ⊤ σ(t k ) ( P σ(t k ) 0 ) -1 Ω σ(t k ) . As t k - τ ( t - k ) = t k -t k-1 ∈ [τ ; τ ] (see Assumption 2) and σ(t - k ) = σ(t k ),
then by also applying a convex argument and the Schur complement to condition (7), one gets that

Γ σ(t k ) ≻ 0. Therefore, ∀k ∈ N, the relation V ( t + k ) ≤ µ c V ( t - k )
holds for the sampling times. This relation and that established in (13) imply that

V ( t+ k ) ≤ μk V ( t- k )
, ∀k ∈ N (14) where constants μk , k ∈ N, are defined by:

μk = { µ s if tk ∈ {s k } k∈N , µ c if tk ∈ {t k } k∈N . ( 15
)
Now passing to the limit in ( 12) when τ tends to tk+1 from below and using ( 14) yields:

V ( t- k+1 ) ≤ μk exp (-α ( tk+1 -tk )) V ( t- k )
, ∀k ∈ N. (16) Now, let t ∈ [0 ; ∞) and denote by N (0, t), the number of sampling and switching times in the interval [0 ; t), that is N (0, t) = N s (0, t) + N c (0, t). We recall that N s (0, t) and N c (0, t) are the number of switching times and sampling times in [0 ; t), respectively. By iterating the relation ( 16) for k = 0, 1, . . . N (0, t) -1, one gets:

V ( t- N (o,t) ) ≤   N (0,t)-1 ∏ k=0 μk   exp ( -α tN(0,t) ) V (0) . (17)
Furthermore, as t ∈ ] tN(0,t) ; tN(0,t)+1 [, by setting τ = t and k = tN(0,t) in (12) and using ( 14), one obtains

V (t) ≤ μN(0,t) exp ( -α ( t -tN(0,t) )) V ( t- N (0,t)
) , which implies according to (17) that:

V (t) ≤   N (0,t) ∏ k=0 μk   exp (-α t) V (0) . ( 18
)
It follows from the definition (15) of constants μk and the definition of N (0, t) that ∏ N (0,t)

k=0 μk = µ Ns(0,t) s µ Nc(0,t) c . From Assumption 1, N s (0, t) ≤ N 0 + t s and from Re- mark 3, N c (0, t) ≥ t τ -1. So, one has ∏ N (0,t) k=0 μk ≤ µ N0 s µ c exp {( ln µ s s + ln µ c τ ) t } . Therefore: V (t) ≤ µ N0 s µ c exp (-2 β t) V (0) . ( 19
)
where β is defined by ( 11). Moreover, one can see that: This concludes the proof. Remark 9. Note that (9) gives a trade-off between the convergence rate in continuous time α, the discrete time convergence rate µ c introduced by G σ , and the parameters µ s that accounts for instability introduced by the switching. As it is the case for systems without switching [START_REF] Raff | Observers with impulsive dynamical behavior for linear and nonlinear continuoustime systems[END_REF][START_REF] Etienne | Observer synthesis under time-varying sampling for lipschitz nonlinear systems[END_REF], no subsystem A i , i ∈ S , has to be stable to ensure the converge of the observer. However in order to verify the LMIs conditions, every pair (A i , C i ), i ∈ S , must be detectable.

V (0) ≤ θ max ∥e (0)∥ 2 , V (t) ≥ θ min ∥e (t)∥ 2 , ( 20

Case of asynchronous switching signals

We assume now that the active mode is not detected instantaneously but with a delay. Such an assumption is coherent with the fact that the output y(t) is not continuously monitored. Therefore looking at σ as an output it is clear that some delay in the active mode reconstruction/measurement is reasonable. Considering x given by (1) and x defined by (2). We define in this case the extended state z (t) = [ x ⊤ (t) e ⊤ (t) ] ⊤ . We also define the observer based feedback u(t) = K σ(t) x(t). Then the following extended dynamical system is defined:

     ż(t) = A σ σ z(t), t ̸ = s k , z(s + k ) = z(s k ), z(t + k ) = J σ(t k )σ(t k ) z(t k ), (21a) (21b) (21c)
where:

A σ σ = [ A σ + B σ K σ -B σ K σ A σ σ + B σσ K σ A σ -B σ σ K σ ] , J σσ = [ I 0 G σ(t k ) C σ σ I -G σ C σ ] , (22a) (22b)
with:

A σ σ = A σ -A σ , B σσ = B σ -B σ , C σ σ = C σ -C σ . ( 23
∀i ∈ S , ∀τ ∈ {0, τ }, Q i (τ ) ≻ 0, (24) ∀i, j ∈ S , ∀τ ∈ {0, τ }, He ( Q i (τ ) A ji ) + Q i 1 + α ij Q i (τ ) ≺ 0, ( 25
) ∀i, j ∈ S , ∀τ ∈ {τ , τ }, [ µ ij Q i (τ ) Q i 0 J ji * Q i 0 ] ≻ 0, (26) ∀i, j ∈ S , ∀τ ∈ {0, τ }, Q i (τ ) ≺ µ s Q j (τ ) , ( 27
) -α[1 - δ s ] + α δ s + lnµ s s + lnµ c τ [1 - δ s ] + δ lnµ c τ s < 0. ( 28
)
where:

Q i (τ ) = Q i 0 + τ Q i 1 , ( 29
) and

(α ij , µ ij ) = { (α, µ c ) if i = j, (-α, µ c ) if i ̸ = j, ( 30
)
then the system (21) is GUES (i.e. the OBSAP is solved) with the convergence rate:

β = 1 2 ( α[1 - δ s ] -α δ s - lnµ s s - lnµ c τ [1 - δ s ] - δ lnµ c τ s ) (31)
Proof 2. Steps of the proof are similar to that followed in the proof of Theorem 8. The main difference lies in the treatment of existing delays in the active mode detection. In this case, we consider the piecewise continuous Lyapunov function candidate

V (t) = z ⊤ (t) Q σ(t) (t -τ (t)) z(t), t ∈ [0 ; ∞) where τ (t) = max{t k : t ≥ t k } and Q σ(t) (t -τ (t)) = Q σ(t) 0 + (t - τ (t)) Q σ(t) 1
. Note that condition (24) implies that V is positive definite. Let's also consider { tk } k∈N , the strictly increasing sequence associated to the sequence

{t k } k∈N ∪ {s k } k∈N . Let t ∈ ] tk ; tk+1 [ with k ∈ N. One has V (t) = z ⊤ (t)Υ σ(t) z(t) with Υ σ(t) the matrix defined by Υ σ(t) = He ( Q σ(t) (t -τ (t)) A σ σ ) + Q σ(t) 1
. Using (25), one can verify that Υ σ(t) ≺ -α σ(t) σ(t) Q σ(t) (tτ (t)) and consequently, V (t) ≤ -α σ(t) σ(t) V (t). Integrating this differential inequality and using the definition (30) of constant α σ(t) σ(t) , one gets:

V (t) ≤ exp{-α M atch ( tk , τ ) + α M iss atch ( tk , τ )} V ( t+ k ) where M atch (τ 1 , τ 2 ) = ∫ τ2 τ1 1 σ(s)=σ(s) ds and M iss atch (τ 1 , τ 2 ) = ∫ τ2 τ1 1 σ(s)̸ =σ(s) ds, ∀τ 1 , τ 2 ∈ [0 ; ∞) with 1, the indicator function. Since M atch (τ 1 , τ 2 ) = τ 2 -τ 1 -M iss atch (τ 1 , τ 2 ), then V (τ ) ≤ exp{γ ( tk , τ )} V ( t+ k ), ∀τ ∈ ] tk ; tk+1 [ (32)
where the notation:

γ (τ 1 , τ 2 ) = -α (τ 2 -τ 1 ) + (α + α) M iss atch (τ 1 , τ 2 ) (33) is introduced for simplicity in the presentation. Note that M iss atch (τ 1 , τ 2 ) satisfies: M iss atch (τ 1 , τ 2 ) ≤ δ N s (τ 1 , τ 2 ) ≤ δ ( N 0 + τ 2 -τ 1 s ) . ( 34
)
Now, we will also determine in this case an upper bound of V ( t+ k ) similar to that obtained in (14). For the switching times, one can see that

V (s + k ) = z ⊤ (s - k )Q σ(s + k ) (s k - τ (s - k ))z(s - k ). Moreover, it follows from (27) that Q σ(s + k ) (s k - τ (s - k )) ≺ µ s Q σ(s - k ) (s k -τ (s - k )). Consequently: V (s + k ) ≤ µ s V (s - k ), ∀k ∈ N.
(35) For the sampling times, using (21c), one can verify that

V (t + k ) -µ σ(t - k )σ(t - k ) V (t - k ) = -z ⊤ (t - k ) Υ σ(t k ) z(t - k ) with Υ σ(t k ) defined by Υ σ(t k ) = µ σ(t - k )σ(t - k ) Q σ(t - k ) (t k -t k-1 ) - J ⊤ σ(t - k ) Q σ(t - k ) (0)J σ(t - k )
. Applying the Schur complement to (26), one can show that Υ σ(t k ) is positive definite and, consequently,

V (t + k ) ≤ µ σ(t - k )σ(t - k ) V (t - k )
, ∀k ∈ N. This inequality and (35) imply that

V ( t+ k ) ≤ μk V ( t- k ), ∀k ∈ N, (36)
where, in this case, μk is the scalar defined by:

μk =    µ s if tk ∈ {s k } k∈N , µ c if tk ∈ {t k } k∈N and σ( tk ) = σ( tk ), µ c if tk ∈ {t k } k∈N and σ( tk ) ̸ = σ( tk ).
(37) By passing to the limit in (32) when τ tends to tk+1 from below and using (36), one obtains:

V ( t- k+1 ) ≤ μk exp{γ ( tk , tk+1 )}V ( t- k ) , ∀k ∈ N. (38)
Now, let t ∈ [0 ; ∞) and denote by N (0, t) the number of sampling times and switching times in [0 ; t), i.e. N (0, t) = N s (0, t) + N c (0, t). Similarly to the proof of ( 18), by iterating relation ( 38) for k = 0, 1, . . . , N (0, t), and using after ( 32) and ( 36) rewritten for τ = t and k = tN(0,t) , one gets:

V (t) ≤   N (0,t) ∏ k=0 μk   exp{γ (0, t)}V (0) . ( 39
)
Let N match c (0, t) be the number of sampling times t k in [0 ; t) with no error in the estimation of the switching signal (i.e. σ(t k ) = σ(t k )) and let N match c (0, t) = N c (0, t) -N match c (0, t). It follows from the definition (37) of μk , that ∏ N (0,t)

k=0 μk = µ Ns(0,t) s µ N match c (0,t) c µ N match c (0,t) c . Since N s (0, t) ≤ N 0 + t s , N match c (0, t) ≤ M iss atch (0,t) τ ≤ δ τ ( N 0 + t s ) and N match c (0, t) ≥ t-M iss atch (0,t) τ ≥ 1 τ [ t ( 1 -δ s ) -δ N 0 ] , then: N (0,t) ∏ k=0 μk ≤ µ N0 s µ N 0 δ τ c µ -δ N 0 τ c exp{ β t}, ( 40
)
where β =

lnµ c τ [1 -δ s ] + δ lnµ c τ s
+ lnµs s . Furthermore, from (34) and (33),

γ (0, t) ≤ (α + α) δ N 0 + t [ - ( 1 - δ s ) α + α δ s ] .
Then combining (39) and ( 40) gives the inequality

V (t) ≤ M exp{-2 β t} V (0) , with β the constant in- troduced in (31) and M the constant defined by M = exp{(α + α) δ N 0 } µ N0 s µ N 0 δ τ c µ -δ N 0 τ c . Finally, using the re- lations V (0) ≤ θmax ∥z (0)∥ 2 and V (t) ≥ θmin ∥z (t)∥ 2 with constants θmax = max i∈S λ max ( Q i 0 ) and θmin = min i∈S min τ ∈[0;τ ] λ min ( Q i (τ ) ) , one obtains ∥z (t)∥ ≤ √ M θmax θmin exp{-β t} ∥z (0)∥ ,
which concludes the proof. Remark 11. The parameter (28) gives a trade off between admissible degradation of convergence of the Lyapunov functions due to both δ and µ s and the convergence ensured by α and µ c . Interestingly, when δ = 0 one recovers (9). Note that, while the existence of parameters α and α is classical for switched systems with delays in the active mode detection, the presence of µ c is not. This is due to the fact that this parameter accounts for impulsive observation when the system's and the observer's active modes mismatch.

Remark 12. The conditions given in both Theorems 8 and 10 amount to solve a set of bilinear matrix inequalities. This is in general a hard problem. However, by choosing a priori a small amount of parameters (namely (α, α, µ c , µ c , µ s ) for Theorem 10 when τ , τ and s are given), the problem at hand becomes LMIs that can be solved for large scale systems. Furthermore, note that this number of parameters is independent on the dimension of the system. Remark 13. Note that here the result obtained do not assume a link between t k , s k and ŝk if one assume that the switching signal is transmitted along the output, simpler condition might be derived, furthermore if one assume that C σ = C, then no attention has to be paid to potential error when triggering the impulsive observation and less conservative results shall be obtained as well.

4. SIMULATION RESULT Consider the following academic example (also considered in [START_REF] Zhang | Asynchronously switched control of switched linear systems with average dwell time[END_REF][START_REF] Zhang | Mode-identifying time estimation and switching-delay tolerant control for switched systems: An elementary time unit approach[END_REF][START_REF] Etienne | Stability analysis for switched linear systems with dwell time and delay in the active mode detection[END_REF] given by:

A 1 = [ 0.2 -0.5 0.5 -0.3 ] , B 1 = [ -0.4 1.8 ] , C 1 = [1, 0], A 2 = [ 0.2 0.3 -1 0.4 ] , B 2 = [ 0.1 1.5 ] , C 2 = [1, 0],
with the parameters τ = 0.3, τ = 0.4 and s = 7.

Syntheses

In this part, with the parameters α = 0.2 and µ s = 5 and µ c = 0.98, we synthesis an impulsive observer without delay in the active mode detection. As is clear in Theorem 8 it is not necessary for the system to be stable (stabilized Fig. 1.a) for the observer to converge (Fig. 1.b). Furthermore, the switching of the original system does not lead to instability of the observer. In Fig. 1.c one can see that t k+1 -t k is randomly drown in the interval [0.3., 0.4].

Analyse

Considering a known (emulation based) stabilizing controller K 1 = [2.1124, -0.9336], K 2 = [-1.3111, -1.3158].

Taking the same parameters as before with furthermore δ = 0.15, µ c = 5 and α = 3, we found that the condition of Theorem 10 are verified. We can see in Fig. 2 that with the proposed controller and observer gains the overall extended system is asymptotically stable event with delay in the active mode detection.

CONCLUSION

In this work we have studied the impulsive observer synthesis problem for switched system with synchronous and Then, an emulation based approach has been used to establish asymptotic stability of an observer based controller in the asynchronous situation. Both results are provided in terms of LMIs and allow for effective computation of the controller action or the convergence rate.

)

  Theorem 10. Consider the hybrid system (21) with switching times {s k } k∈N and sampling sequences {t k } k∈N satisfying Assumptions 1 and 2. If there exist symmetric matrices Q i m ∈ S 2n , m ∈ {0, 1}, i ∈ S , and scalars α, α ∈ R and µ s , µ c , µ c ∈ R + with µ s , µ c ≥ 1 and µ c ≤ 1 such that the following conditions are verified:

  Fig. 1. a)Trajectory of the system. b) Trajectory of the observation error (plain line measured states dotted lined unmeasured state). c)Inter sampling time.

  Consider the switched linear system (4) with switching times {s k } k∈N and sampling sequences {t k } k∈N satisfying Assumptions 1 and 2. If there exist symmetric matrices P i m ∈ S n , m ∈ {0, 1}, i ∈ S and the scalars α ∈ R, µ s ≥ 1 and µ c ∈ ]0 ; 1] such that the following statements hold:

	(4b)
	(4c)
	LMIs-based conditions are proposed in Theorem 8 to solve
	the OSP problem:
	Theorem 8.

ACKNOWLEDGMENTS

This work was partially funded by a CPER Data project, which is co-financed by European Union with the financial support of European Regional Development Fund (ERDF), French State and the French Region of Hautsde-France.


Andréa Le Beller

Mouhamad Nasser

Julie Traclet

Anne-Sophie Blanchet-Legens²

Caroline Kannengiesser

Jean-François Mornex

Vincent Cottin
email: vincent.cottin@chu-lyon.fr

Combined pulmonary fibrosis and emphysema in alpha-1-antitrypsin deficiency

Keywords: emphysema, pulmonary fibrosis : alpha-1-antitrypsin, genetics

2

Combined pulmonary fibrosis and emphysema (CPFE) is characterized by severe dyspnea at exertion, distinct physiology with preserved static and dynamic lung volumes contrasting with significantly altered gas exchange capacity, and a combination of upper lobe emphysema and features of pulmonary fibrosis predominating in the lower lobes on chest CT [START_REF] Cottin | Combined pulmonary fibrosis and emphysema: a distinct underrecognised entity[END_REF]. A variety of patterns have been described both at chest imaging and histopathology, including usual interstitial pneumonia (UIP), non-specific idiopathic pneumonia, desquamative interstitial pneumonia, and more recently smoking related interstitial fibrosis [START_REF] Chae | Differentiating Smoking-Related Interstitial Fibrosis (SRIF) from Usual Interstitial Pneumonia (UIP) with Emphysema Using CT Features Based on Pathologically Proven Cases[END_REF]. Prognosis is dismal in progressive forms, with a median survival of 3 to 5 years.

Mortality is mostly driven by the development of precapillary pulmonary hypertension, acute exacerbation of fibrosis, and lung cancer [START_REF] Cottin | Combined pulmonary fibrosis and emphysema: a distinct underrecognised entity[END_REF].

Very few genetic variants have been associated with CPFE. Herein we report on a case of CPFE in the context of alpha-1-antitrypsin (A1AT) deficiency (AATD) with progressive course leading to death.

Observation

A 62-year-old female, who had smoked less than 5 pack-years until the age of 30 years, was diagnosed with liver cirrhosis, which was then related to AATD, with a low serum level of A1AT (0.23 g/L), and a PiZZ protease inhibitor phenotype. She had progressive shortness of breath over the past few months as well. Pulmonary investigations at that time demonstrated predominantly paraseptal emphysema in the upper lobes at chest CT.

Pulmonary function tests were normal. Therefore, augmentation therapy was not introduced and the patient followed a watchful waiting strategy. Six years later, the patient (now aged 68) was admitted for progressive chronic dyspnea limiting daily activities. Lung examination found basal "Velcro" crackles. The chest CT now demonstrated bilateral fibrotic changes, with basal and subpleural predominance, consisting of reticulations, traction bronchiectasis, and honeycombing. Paraseptal emphysema in the upper lobes was unchanged. Forced vital capacity (FVC) was 2.52L (111% of predicted value), forced expiratory volume in one second was 2.11L (112% of predicted), and the diffusion capacity of the lungs for carbon monoxide was 42% of predicted. Autoimmune studies were negative. With these findings, the diagnosis of CPFE was established, in this case corresponding to idiopathic pulmonary fibrosis (IPF) combined with predominantly paraseptal emphysema. The stage of cirrhosis was Child-Pugh A, allowing initiation of antifibrotic therapy with nintedanib, later switched to pirfenidone due to disease progression.

Nevertheless, lung fibrosis inexorably progressed with an average decline of FVC of 270 mL/year. The patient died from chronic respiratory failure 5 years following CPFE diagnosis.

Because she had a sister with pulmonary fibrosis, genetic tests were performed. No mutation was found in the telomerase-related genes (TERT, TERC, RTEL1, DKC1, PARN, TINF2) and a panel of other genes associated with pulmonary fibrosis. Homozygous mutation p.E366K was found in the SERPINA1 gene coding for A1AT.

Discussion

The current case is the first reported association of CPFE syndrome and AATD.

Despite the wealth of reports on genetic predisposition to pulmonary fibrosis, large cohorts of patients with AATD, and series of patients with CPFE, the association between AATD and pulmonary fibrosis seems very rare. In the seminal description of CPFE as a syndrome, the level of A1AT was low (0.19 g/L) in only one of 65 patients with CPFE, who was PiZZ homozygous, and another patient with a PiMZ phenotype had a normal serum level of A1AT (1.6 g/L). Among homozygous twin sisters with a PiZZ phenotype, one had liver cirrhosis and the other developed IPF [START_REF] Calabro | Idiopathic Pulmonary Fibrosis Associated with Alpha1-Antitrypsin Deficiency: Concomitant Finding or Real Association?[END_REF]. A cluster of both familial pulmonary fibrosis and AATD (with PiSZ or PiMZ phenotypes in several family members yet normal serum levels of A1AT) was reported in one family, including one individual who had IPF, and one who had both IPF and liver cirrhosis [START_REF] Kim | alpha 1-antitrypsin deficiency and idiopathic pulmonary fibrosis in a family[END_REF]. Another patient had liver cirrhose due to AATD and IPF but no emphysema [START_REF] Afonso | A rare mutation on alpha-1 antitrypsin deficit and lung fibrosis: case report[END_REF]. A highly significant increase in the frequency of PiMZ phenotype was also found among patients with pulmonary fibrosis, both with and without rheumatoid arthritis [START_REF] Geddes | alpha 1-antitrypsin phenotypes in fibrosing alveolitis and rheumatoid arthritis[END_REF].

Our case presented a radiologic UIP pattern, combined with upper lobe paraseptal emphysema and without characteristic panlobular emphysema as commonly observed in AATD. In addition, large cysts with thickened walls were present in lung bases, suggesting association with smoking-related interstitial fibrosis [START_REF] Chae | Differentiating Smoking-Related Interstitial Fibrosis (SRIF) from Usual Interstitial Pneumonia (UIP) with Emphysema Using CT Features Based on Pathologically Proven Cases[END_REF]. The lung function profile at baseline was characteristic of the CPFE syndrome and supported the diagnosis [START_REF] Cottin | Combined pulmonary fibrosis and emphysema: a distinct underrecognised entity[END_REF]. However, the outcome was characterised by progressive restrictive lung physiology and eventually death despite treatment, in a timeframe consistent with IPF. Antifibrotic therapy could be prescribed as the stage of liver cirrhosis was Chid-Pugh A.

The pathophysiology of CPFE has not been fully elucidated. Tobacco smoking is by far the most common etiology. Other risk factors such as increasing age, male sex, and environmental or occupational exposures are shared with IPF. CPFE occurs also in nonsmoker subjects with connective tissue diseases and hypersensitivity pneumonitis. CPFE has been reported in association with gene variants implicated in pulmonary fibrosis, including mutations in telomerase-related genes [START_REF] Nunes | Is telomeropathy the explanation for combined pulmonary fibrosis and emphysema syndrome?: report of a family with TERT mutation[END_REF] and genes coding for surfactant proteins [START_REF] Cottin | Combined pulmonary fibrosis and emphysema syndrome associated with familial SFTPC mutation[END_REF]. However, less attention has been paid to potential genetic components generally associated with an increased risk of emphysema. AATD is a genetic disorder with more than 100 genetic variants described, and an estimated prevalence of 23 per 100,000 persons [START_REF] Greulich | The prevalence of diagnosed α1-antitrypsin deficiency and its comorbidities: results from a large population-based database[END_REF]. The estimated prevalence of IPF in France is 8.2 per 100,000 persons [START_REF] Duchemann | Prevalence and incidence of interstitial lung diseases in a multi-ethnic county of Greater Paris[END_REF]. It is noteworthy that mostly familial forms of IPF [START_REF] Calabro | Idiopathic Pulmonary Fibrosis Associated with Alpha1-Antitrypsin Deficiency: Concomitant Finding or Real Association?[END_REF][START_REF] Kim | alpha 1-antitrypsin deficiency and idiopathic pulmonary fibrosis in a family[END_REF] or CPFE (the present case) have been reported in combination with AATD. Indeed, our patient had a family history of pulmonary fibrosis, although no genetic variant could be identified in the proband.

What caused the coexistence of AATD and CPFE in our patient is unknown. The two conditions may well be associated by chance, with two genetic disorders clustered in the same family, as previously reported for AATD and IPF [START_REF] Kim | alpha 1-antitrypsin deficiency and idiopathic pulmonary fibrosis in a family[END_REF]. The patient was a smoker, which is a risk factor shared by emphysema and pulmonary fibrosis. It remains to be elucidated whether AATD might predispose to, or protect from pulmonary fibrosis. Although emphysema and fibrosis have distinct pathophysiology [START_REF] Hanaoka | Comparison of gene expression profiling between lung fibrotic and emphysematous tissues sampled from patients with combined pulmonary fibrosis and emphysema[END_REF], experimental data in pallid mice that are deficient in serum alpha-1-proteinase inhibitor, suggested that both may coexist after exposure to bleomycin and may be potentiated by AATD, potentially through neutrophil elastase [START_REF] Cavarra | Genetic deficiency of alpha1-PI in mice influences lung responses to bleomycin[END_REF]. Indeed, AATD increases sensitivity to oxidative stress and inflammatory mediators within the lung, and might conceivably predispose to a pro-inflammatory and pro-fibrotic response to lung injury. If true, however, one would expect to encounter a sizable number of subjects with AATD-related emphysema also harboring pulmonary fibrosis. Conversely, administration of alpha-1-proteinase inhibitor (A1AT) to bleomycin-treated hamsters reduced experimental pulmonary fibrosis [START_REF] Nagai | Administration of alpha 1-proteinase inhibitor ameliorates bleomycininduced pulmonary fibrosis in hamsters[END_REF]. A1AT supplementation therapy was not prescribed in our patient because IPF predominated over emphysema; however, whether it might have had an impact on the course of IPF, if prescribed, is unknown.

Altogether, the current data mostly support the coincidental association of two genetic conditions predisposing to emphysema (AATD) and to pulmonary fibrosis. Whether AATD may alter the risk of pulmonary fibrosis in either way remains to be determined.

Figure 1 .

 1 Figure 1. Axial CT-scan of the lung in a patient with CPFE. A: limited paraseptal emphysema in the left upper lobe, associated with mild thickening of the cyst walls. B: usual interstitial pneumonia pattern in the lung bases, with honeycombing, subpleural reticular opacities and traction bronchiectasis in the lung bases. C: worsening of fibrotic features two years later than C; D: large cyst with thick walls in the right lower lobe, suggesting a pattern of smoking-related interstitial fibrosis.

Acknowledgements. We thank Dr Corinne Villaudy (Lyon, France) for the clinical management in hepatology.