Keywords: Hybrid systems, DAE, index, nonstandard analysis Systèmes hybrides, DAE, index, analyse non standard

Hybrid systems modelers exhibit a number of difficulties related to the mix of continuous and discrete dynamics and sensitivity to the discretization scheme. Modular modeling, where subsystems models can be simply assembled with no rework, calls for using Differential Algebraic Equations (DAE). In turn, DAE are strictly more difficult than ODE. ¶ In most modeling and simulation tools, before simulation can occur, sophisticated pre-processing is applied to DAE systems based on the notion of differentiation index. Graph based algorithms such as the one originally proposed by Pantelides are efficient at finding the index, structurally (i.e., outside some exceptional values for the system parameters). The differentiation index for DAE explicitly relies on everything being differentiable. Therefore, extensions to hybrid systems must be done with caution-to our knowledge, no such extension exists. We propose to rely on non-standard analysis for this. Non-standard analysis formalizes differential equations as discrete step transition systems with infinitesimal time basis. We can thus bring hybrid DAE systems to their non-standard form, where the notion of difference index can be firmly used-the difference index of a difference Algebraic Equation (dAE) is an easy transposition of the differentiation index, in which forward shift replaces differentiation. We prove that the differentiation index of a DAE is structurally equal to the difference index of its non-standard interpretation, which is a dAE. We can thus propose the difference index of the non-standard semantics of a hybrid DAE system, as a consistent extension of both the differentiation index of DAE and the difference index of dAE. It turns out that the index theory for (discrete time) dAE systems is interesting in itself and raises new issues. We discuss graph based algorithms à la Pantelides for computing the dAE index and discuss examples.

The booming of system design during the last decade has called for drastic changes in techniques and tools used for physical system modeling. So far the dominant technology for physical system modeling has been and still is Simulink, 1 2 in which both discrete and continuous time systems must be expressed in state space form-we give here the continuous time version of it:

ẋ = f (x, u) y = g(x, u) (1)
where u is the input vector, x the state vector, and y the output vector. Systems of the form (1) consist of Ordinary Differential Equations (ODE). They compose as input/output functions, which requires that no integrator-free loop is created as the result of the composition. Unfortunately, commonly encountered physical systems are naturally modeled by balance equations (Ohm and Kirchhoff laws for electrical networks, Lagrange equations for mechanical systems, thermal systems), which naturally leads to integrator-free loops. Manually moving from balance equations to state-space forms like (1) quickly becomes a significant burden. Since such a translation, from balance equations to state-space form, must be performed globally for the entire system, the reuse of partial models from libraries gets impaired [START_REF] Mattsson | Physical system modeling with Modelica[END_REF], [START_REF] Fritzson | Introduction to Modeling and Simulation of Technical and Physical Systems with Modelica[END_REF], [START_REF]Principles of Object-Oriented Modeling and Simulation with Modelica 3.3[END_REF].

The need for an adapted modeling framework was already recognized in the 70's by the engineering community with the proposal of bond graphs [START_REF] Thoma | Bond graphs: introduction and applications[END_REF], [START_REF] Gawthrop | Metamodelling: bond graphs and dynamic systems[END_REF], [START_REF] Ljung | Modeling of Dynamic Systems[END_REF]. Bond graphs guide the development of models from first physical principles by identifying generic concepts such as power, effort, and flow. The causality analysis of a bond graph allows deducing who in the model should be considered as an input, state, or output. Model reuse occurs at the level of bond graphs, whereas the causality analysis is performed automatically and globally for the system, prior to simulation. The causality analysis of bond graphs is a graphical problem.

The mathematical counterpart of models from physical principles is that of a Differential Algebraic Equation (DAE) [START_REF] Petzold | ODE methods for the solution of differentialalgebraic systems[END_REF], [START_REF] Petzold | Differential algebraic equations are not ODEs[END_REF], of which a simple instance is-again in continuous time:

0 = f ( ẋ, x, u) 0 ≤ g(x, u)
Here we have shown a combination of equality and inequality (also called bilateral and unilateral) constraints; x is the state vector collecting all variables subject to differentiation, whereas u is the algebraic vector collecting other variables. Observe that we do not mention inputs nor outputs. Unlike state-space systems, DAE systems compose with no restriction. Modeling languages based on DAEs were initially proposed in the early 80's by H. Elmqvist [START_REF] Elmqvist | A structured model language for large continuous systems[END_REF] from the automatic control group lead by K.J. Aström in Lund. The Modelica language [START_REF] Mattsson | Physical system modeling with Modelica[END_REF], [START_REF] Fritzson | Introduction to Modeling and Simulation of Technical and Physical Systems with Modelica[END_REF], [START_REF]Principles of Object-Oriented Modeling and Simulation with Modelica 3.3[END_REF] subsequently resulted from this effort, grounding the birth of the Modelica community. 3 1 www.mathworks.com/products/simulink/ 2 http://www.mathworks.com/help/simulink/ 3 https://www.modelica.org/ There is, however, no free lunch, and L.R. Petzold, one of the mathematicians having grounded the foundations of DAE theory, entitled one of his papers: "DAEs are not ODEs" [START_REF] Petzold | Differential algebraic equations are not ODEs[END_REF]. As an example, consider the following system, which we state in both continuous time (left handside) and discrete time (right handside) versions:

ẋ = f (x, u) 0 = g(x) x • = f (x, u) 0 = g(x) (2)
In (2) both x and u are scalar signals, taking their values in R. For the discrete time version, x • denotes the forward shifted version of signal x, i.e., such that x • k = x k+1 for every discrete instant k = 0, 1, 2, . . . Systems of the form (2) exhibit so-called latent constraints. For the discrete time version, additional constraint 0 = g(x • ) follows by shift invariance, which furthers constrains the triple (x • , x, u). Similarly, for the continuous time version, differentiating the second equation yields the additional constraint 0 = d dt g(x) = g (x) ẋ = g (x)f (x, u) where g denotes the derivative of g, which furthers constrains the triple ( ẋ, x, u). To summarize, [START_REF] Aubry | Pseudo-symplectic Runge-Kutta Methods[END_REF] should rather be reformulated by making the latent constraints explicit (shown in red):

   ẋ = f (x, u) 0 = g(x) 0 = g (x) ẋ    x • = f (x, u) 0 = g(x) 0 = g(x • ) (3)
Further shifting the last equation in the discrete time version of (3) would bring the fresh variable x •• , denoted by x •2 , thus resulting in no further constraint on the triple (x • , x, u). Similarly, further differentiating the last equation of the continuous time version of (3) does not bring anything useful. Focus once more on the discrete time version of (3), which we rewrite by substituting f (x, u) for x • in its last equation:

   x • = f (x, u) 0 = g(x) 0 = g(f (x, u)) (4)
Assume the current value of state x is consistent, meaning that constraint g(x) = 0 is satisfied. Then, the last equation determines u as an output of the system and then, the pair (x, u) determines the next value for the state x • by using the first equation, ensuring that x • remains consistent. This was an instance of the kind of causality analysis we can perform to infer the input/output status of the different algebraic variables. The same analysis can be performed for the continuous time version of (3), replacing forward shift by differentiation. One point should be noticed: it could happen that the last two equations of (4) are redundant. This would occur if f (x, u) ≡ x. In this case, our causality analysis can still be performed but leads to incorrect conclusions: u is not an output any more as it disappears from the equations of the system. The kind of causality analysis that we performed leads to valid conclusions but for exceptional instances of the functions f and g: properties inferred in this way are called structural.

The process of revealing latent constraints by shifting or differentiating is called index reduction and the number of times it has to be done is called the differentiation or difference index, for continuous time and discrete time systems, respectively. The notion of differentiation index was proposed in the late 80's, see [START_REF] Campbell | A computational method for general higher index nonlinear singular systems of differential equations[END_REF], [START_REF] Campbell | The index of general nonlinear DAEs[END_REF]. Index reduction goes along with causality analysis, performed by graph based algorithms proposed for the first time by Pantelides [START_REF] Pantelides | The consistent initialization of differential-algebraic systems[END_REF] and subsequently developed and modified by several authors [START_REF] Erik | Index reduction in Differential-Algebraic Equations using dummy derivatives[END_REF]. Some physical modeling tools use index reduction to prepare for the work of the solvers, which are designed to simulate DAE of low index only (1, 2, or at most 3). Some solvers address the issue of latent constraints (including invariants such as occurring in Hamiltonian or Symplectic systems) in an implicit way, by using appropriate discretization schemes [START_REF] Aubry | Pseudo-symplectic Runge-Kutta Methods[END_REF]. Nevertheless, index reduction was originally proposed to find consistent initial conditions for the DAE and still remains useful for that purpose. In addition, preprocessing required by some solvers for so-called non-smooth systems [START_REF] Vincent | Numerical Methods for Nonsmooth Dynamical Systems[END_REF] amount to computing the relative degree, which appears closely related to the index. Surprisingly enough, no notion of index exists for hybrid extensions of DAE. When equations are used in both the continuous dynamics (in each mode) and the reset of states at mode transitions, it becomes unclear how such resets must be performed, thus leading to latest physical systems modelers to refuse certain, otherwise mathematically well defined, models.

The intuition conjectures that the right notion of index should imply that, in each mode, the index should coincide with the known notion of index for DAE. This looks fine also when transitions between modes are isolated. How to handle cascades of mode transitions, however, seems out of reach of the intuition.

In this work, we show that the above intuition is indeed incorrect and we provide a formal definition of the index, for hybrid DAE systems. While the notion of index remains global and "mode independent", the dependence on modes actually appears while performing the graph based causality analysis supporting index evaluation. Our approach is as follows:

1) We observe that a notion of difference index can be defined for discrete time algebraic dynamical systems, we call them dAE, (i.e., involving constraints) by just borrowing, to discrete time dynamics, the notion of differentiation index. 2) We use the nonstandard analysis interpretation of DAE-a nonstandard interpretation of a DAE consists in interpreting the derivative as a difference operator using an "infinitesimal" step. This interpretation is exact in a certain sense, it is not an approximation unlike discretization schemes-it is not effective, however (there is no free lunch). We show that the two notions of index for a DAE (the classical differentiation index and the nonstandard difference index) coincide. 3) We take the nonstandard interpretation of a hybrid DAE system and consider its difference index. Thanks to the previous item, this notion subsumes the differentiation index of DAE and the difference index of dAE.

In fact, our study does not address the "true" index [START_REF] Campbell | The index of general nonlinear DAEs[END_REF] (see the forthcoming Definition 2), but rather the structural index. For a linear DAE system, the index equals the structural index almost everywhere when the non-zero coefficients vary over a neighborhood, and this extends to nonlinear DAE systems by considering the tangent linear DAE. Highly efficient graph based algorithms exist for computing the structural index, including the well-known first one proposed by Pantelides [START_REF] Pantelides | The consistent initialization of differential-algebraic systems[END_REF]. This graph-based analysis proposed by Pantelides is indeed a causality analysis, very much related to the causality analysis of bond graphs.

As our study concerns the structural notions of index, the term index will refer to the structural index in the sequel, whereas we use the term true index to refer to the original notion of differentiation index.

The paper is organized as follows. The background on differentiation index theory for DAE is recalled in Section II.

After recalling what the index is following the original definition given in Campbell and Gear [START_REF] Campbell | The index of general nonlinear DAEs[END_REF], we introduce the structural index and relate it to the notion of structurally nonsingular matrices. Structural properties of matrices are characterized by the existence of certain permutation matrices, whose quest relies on graph based algorithms. Section III introduces the index theory for difference Algebraic Systems (dAE), the discrete time counterpart of DAE; dAE can also be seen as discrete time dynamical systems subject to a transition relation (not function); we focus on dAE for which the transition relation is specified as a set of equations. While the definition of the (true) index is a straightforward translation from the continuous time case, new issues arise, due to the consideration of nonsmooth systems. We devote subsection III-E to the discussion of illustration examples, thus motivating the new techniques of guarded causality analysis we develop in Section IV-these techniques turn out to be closely reminiscent of those used in the compilation of synchronous languages [START_REF] Benveniste | The synchronous languages 12 years later[END_REF], [START_REF] Benveniste | Compositionality in dataflow synchronous languages: Specification and distributed code generation[END_REF]. In Section V we recall how DAE can be given a discrete time (dAE) interpretation by using nonstandard analysis with infinitesimals, which yields an alternative way of defining the index, based on this dAE interpretation. We show that, for smooth DAE systems, the two notions of structural index coincide, which allows us to rely on nonstandard semantics for the index theory of hybrid DAE systems. The latter is developed in Section VI. We propose a mini-formalism of guarded equations that is expressive enough to cover practical examples, including inequality constraints, also called "unilateral", and variations thereof. We give the nonstandard semantics of such systems, which provides as a byproduct a notion of structural index together with the graph based algorithm for computing it. Examples are revisited in Section VII and hints for effective modular algorithms for performing causality analyses and computing the index are sketched in Section VIII.

II. BACKGROUND ON DIFFERENTIATION INDEX THEORY

In this section we recall the background on differentiation index theory. We pay particular attention to the foundations of graph based algorithms associated to the quest for the structural index. To this end, we recall the background on structurally nonsingular matrices with some detail and we give proofs.

A. DAE differentiation index

The basic reference is the work of Campbell and Gear [START_REF] Campbell | The index of general nonlinear DAEs[END_REF], see also Mattsson and Söderlin [25]. Let R denote the set of reals, Z = {. . . , -1, 0, +1, . . . } the set of integers, and N = {0, +1, . . . } the set of non-negative integers. In this section we consider time invariant DAE problems of the following form:4

F (x, ẋ) = 0 (5)
where x takes its values in R n and F in R m . In the sequel, F x and F ẋ denote the partial derivatives of F with respect to the first and second variables of F , respectively. For the following definition, the reader is referred to [START_REF] Campbell | The index of general nonlinear DAEs[END_REF].

Definition 1: DAE (5) is solvable 5 in the connected open set Ω ⊂ R 2n if there are connected open sets Λ ⊂ R ρ and I ⊂ R and a function (t, λ) → Φ(t, λ) such that: 1) Θ(t, λ) = (t, Φ(t, λ)) is a diffeomorphism of I×Λ into R n+1 . 2) Φ(t, λ) is a solution of (5) for each value of λ. 3) (Φ(t, λ), d dt Φ(t, λ)) ∈ Ω for every λ ∈ Λ and t ∈ I. 4) If x(t) is a solution of (5) such that (x(t), ẋ(t)) ∈ Ω for some t ∈ I, then it holds that x(t) = Φ(t, λ) for some λ ∈ Λ. A pair (t, x) is called consistent if x = Φ(t, λ)
holds for some λ.

Condition 2) expresses that λ acts as a daemon solving the possible nondeterminism. Condition 4) states that λ captures all the nondeterminism. Indeed, λ parameterizes consistent initial conditions, which, in turn, determine solutions of (5). Systems with exogeneous inputs, of the form, e.g.:

F (x, ẋ, u) = 0 (6)
are used in control and when composing subsystems to form larger systems. Systems of the form (6) are a specialization of (5) by putting y = (x, u) and reformulating it as a DAE with larger state y. Systems of the form (6) leave generally some freedom on exogeneous u (subject to the constraints) when selecting solution Φ(t, λ).

The true differentiation index for DAE [START_REF] Benveniste | The synchronous languages 12 years later[END_REF] is defined as follows. The kth derivative array associated to (5) is:

     F (x, ẋ) d dt F (x, ẋ) . . .x, ẋ) d k dt k F (x, ẋ)      = def      F (0) (x, ẋ, w) F (1) (x, ẋ, w) . . .x, ẋ, w) F (k) (x, ẋ, w)      = def F k (x, v, w) (7)
where

v = def ẋ (8)
and w = def (x (2) , . . . , x (k+1) )

where we recall that d dt F (x, ẋ) = F x (x, ẋ) ẋ + F ẋ(x, ẋ)ẍ, and so on for higher degree derivatives. In [START_REF] Diener | Analyse non standard[END_REF] x (1) = ẋ, x (2) = ẍ, x (3) = . . . denote the successive derivatives of x.

The reason for considering the kth derivative array equations is the following. The additional equations d j dt j F (x, ẋ) = 0 added when forming the array are implied by the original DAE system, due to its time-invariance. These additional equations, however, add new equations and new variables, namely some components of x (2) involved in these new equations. Some of the new equations may not bring fresh variables, but only reuse previous variables, which they further constrain. In this case, latent constraints get revealed. Such latent constraints were not "visible" in the original formulation of the DAE system but are a consequence of its time-invariance.

Following again [START_REF] Campbell | The index of general nonlinear DAEs[END_REF], a value x is called consistent for (7) if there exists (v, w) such that

F k (x, v, w) = 0 (10)
seen as an algebraic equation. Given a consistent value x for [START_REF] Campbell | A computational method for general higher index nonlinear singular systems of differential equations[END_REF], algebraic equation ( 10) will generally have a set of solutions for (v, w). Definition 2: Assume that DAE (5) is solvable. The true differentiation index of this DAE, denoted by ν D , is the smallest index k such that v is uniquely determined by the algebraic equation [START_REF]Nonstandard analysis and its applications[END_REF] for any consistent value x for [START_REF] Campbell | A computational method for general higher index nonlinear singular systems of differential equations[END_REF].

That is, for k ≥ ν D , the map

x → ∃w.F k (x, v, w) = 0 (11)
defines v as a deterministic function of x. Since ( 7) is equivalent to the original DAE, (11) determines ẋ and is, therefore, sort of an ODE that can be solved for x. At this point, it is worth detailing how quantification over time applies to the map defined by [START_REF] Elmqvist | A structured model language for large continuous systems[END_REF]:

∀t ∈ R : x t → ∃w t .F k (x t , v t , w t ) = 0 (12)
By the implicit function theorem, [START_REF] Elmqvist | A structured model language for large continuous systems[END_REF] is, locally around a triple (v o , w o , x o ), equivalent to the following condition:

the map x → ∃w.Av + Cw + Ex = 0 defines v as a deterministic function of x (13)
where A, C, E are the Jacobians

A = (F k ) v , C = (F k ) w , E = (F k ) x , at the considered triple (v o , w o , x o ).
As an illustration example consider the well known pendulum example in Cartesian coordinates:

ẍ = T x ÿ = T y -g L 2 = x 2 + y 2 (14)
In this example, x, ẋ, y, ẏ are the states and T is an algebraic variable. The form (5) for DAE system ( 14) is:

ẋ = u u = T x ẏ = v v = T y -g L 2 = x 2 + y 2 (15) RR n°8630
This is not of index 0 since the Jacobian with respect to ẋ, u, ẏ, v, T is singular, leaving T undefined since the last row is identically zero:

      1 0 0 0 0 0 1 0 0 -x 0 0 1 0 0 0 0 0 1 -y 0 0 0 0 0      
So we must differentiate the system. Differentiating, in [START_REF] Lee | Leveraging synchronous language principles for heterogeneous modeling and design of embedded systems[END_REF], the last equation twice yields:

ẋ = u (i1) u = T x (i2) ẏ = v (ii1) v = T y -g (ii2) L 2 = x 2 + y 2 (iii) 0 = ẋx + ẏy (iv) 0 = ux + ẋ2 + ẏ2 + vy (v) (16)
Unknowns of highest derivative order are ẋ, u, ẏ, v, T . Rewriting all equations (i-v) in the form 0 = . . . yields the following Jacobian for the equations involving ẋ, u, ẏ, v, T :

             ∂(i1) ∂ ẋ ∂(i1) ∂ u ∂(i1) ∂ ẏ ∂(i1) ∂ v ∂(i1) ∂T ∂(i2) ∂ ẋ ∂(i2) ∂ u ∂(i2) ∂ ẏ ∂(i2) ∂ v ∂(i2) ∂T ∂(ii1) ∂ ẋ ∂(ii1) ∂ u ∂(ii1) ∂ ẏ ∂(ii1) ∂ v ∂(ii1) ∂T ∂(ii2) ∂ ẋ ∂(ii2) ∂ u ∂(ii2) ∂ ẏ ∂(ii2) ∂ v ∂(ii2) ∂T ∂(iv) ∂ ẋ ∂(iv) ∂ u ∂(iv) ∂ ẏ ∂(iv) ∂ v ∂(iv) ∂T ∂(v) ∂ ẋ ∂(v) ∂ u ∂(v) ∂ ẏ ∂(v) ∂ v ∂(v) ∂T              =         1 0 0 0 0 0 1 0 0 -x 0 0 1 0 0 0 0 0 1 -y x 0 y 0 0 2x x 2y y 0        
which, by reordering the rows, yields the Jacobian:

        x 0 y 0 0 1 0 0 0 0 0 1 0 0 -x 0 0 1 0 0 2x x 2y y 0 0 0 0 1 -y         (17)
Removing the first (red) row yields, under the condition y =0, a structurally nonsingular Jacobian, meaning that it is nonsingular but for exceptional values for the pair (x, y). Hence, ẋ, u, ẏ, v, T is determined as a function of other variables. When y is zero or close to zero, then, due to constraint (iii), x is not small and we can exchange the roles of x and y. Hence, the index was found equal to 2.

In principle, Definition 2 requires that we not only differentiate (iii) twice, but also (i1-ii2). This would, however, introduce fresh variables x (2) , x (3) , u (2) , u (3) , y (2) , y (3) , v (2) , v (3) , which enter the w of [START_REF] Elmqvist | A structured model language for large continuous systems[END_REF]; eliminating this w is simply achieved by ignoring the differentiation of (i1-ii2).

Following [START_REF]Nonstandard analysis and its applications[END_REF], a consistent value for the tuple x, u, y, v must satisfy the following equations, obtained by substituting ẋ using (i1) and ẏ using (ii1) in (iv):

L 2 = x 2 + y 2 (iii) 0 = ux + vy (iv) (18)
The remaining equations form an ODE with highest order derivatives ẋ, u, ẏ, v, T , since the Jacobian is structurally nonsingular (outside a neighborhood of y = 0):

ẋ = u (i1) u = T x (i2) ẏ = v (ii1) v = T y -g (ii2) 0 = ux + u 2 + v 2 + vy (v) (19)
The reasoning regarding the Jacobian obtained by erasing the first row of matrix [START_REF] Lindstrøm | An invitation to nonstandard analysis[END_REF] refers to so-called structural properties of matrices, which we recall now.

B. Structurally nonsingular matrices

Structurally nonsingular matrices play a central role in finding the differentiation degree using graph based algorithms. We thus recall some basic material in this section.

Say that a property P (x 1 , . . . , x k ) involving the real variables x 1 , . . . , x k holds almost everywhere if it holds for every x 1 , . . . , x k outside a subset of R k having empty interior. Matrix A is called structurally onto if it remains almost everywhere onto (surjective) when its non-zero entries vary over some neighborhood. Square n×n-matrix P is a permutation matrix if and only if there exists a permutation σ of the set {1, . . . , n} such that p ij = 1 if j = σ(i) and p ij = 0 otherwise. Pre-and post-multiplication of a matrix A by a permutation matrix results in permuting the rows and columns of this matrix. The following result holds [?]:

Lemma 1: A is structurally onto if and only if there exist two permutation matrices P and Q, of appropriate dimensions, such that PAQ = B 1 B 2 , where B 1 is square with nonzero diagonal.

Proof: See Appendix A1.

Lemma 1 specializes to the following simpler result:

Corollary 1: A is structurally nonsingular if and only if PA has a nonzero diagonal (all entries of the diagonal are nonzero) for some permutation matrix P .

Matrix A is structurally nonsingular if and only if the linear equation Av=y, where v is the unknown, has a unique solution for any y, for almost all values for the non-zero entries of matrix A. We will, however, be interested in the more general equation, to be solved for v:

∃w : Av + Cw + x = [ A C ] v w + x = 0 (20)
where x has dimension n, A is a n×p-matrix, C is a n×qmatrix, and v, w are of appropriate dimensions. In other words, when considering equation Av+Cw+x=0, we see v as the unknown (or output), x as the input, and w as a don't carewe are not interested in the particular value for w, we just want it to exist. We formalize this next. Definition 3: Say that the pair (A, C) is structurally nonsingular if equation [START_REF] Fritzson | Introduction to Modeling and Simulation of Technical and Physical Systems with Modelica[END_REF] uniquely defines v as a function of x when x is consistent, almost everywhere when the non zero entries of (A, C) vary over some neighborhood. The following result generalizes Corollary 1:

Lemma 2: Pair (A, C) is structurally nonsingular if and only if there exists a permutation matrix P , of dimension n, such that

P A =   A 1 A 2 A 3   and P C =   C 1 0 0   ( 21
)
where A 2 is square with nonzero diagonal and C 1 is structurally onto. Proof: See Appendix A2.

C. Graph based algorithms

In this section, we recall the well-known graph based algorithm, originally due to Pantelides [START_REF] Pantelides | The consistent initialization of differential-algebraic systems[END_REF], which implements the quest for pivoting, i.e., for the two permutation matrices A and C of Lemma 2-we specify the algorithm but otherwise pay no attention to its efficient implementation, unlike the original work of Pantelides.

We assume three disjoint sets X , V, W of underlying realvalued variables-we also indicate the generic notation for their elements: x∈X , v∈V, w∈W, and we set X = def X V W where denotes the disjoint union. For E a scalar equation, i.e., of the form F (. . . ) = 0 where F is a real-valued linear function of variables belonging to X, write z ∈ E to indicate that variable z occurs in equation E with a coefficient that is structurally non-zero. Let S = i∈I E i be a system consisting of the conjunction of a set {E i | i ∈ I} of scalar equations over a finite subset of variables belonging to X. For convenience, we shall identify S with its set

{E i | i ∈ I} of constitutive equations. Define
x ∈ S if and only if x ∈ E i holds for some i ∈ I and denote by ∈ S the set of variables x ∈ X such that x ∈ S. By decomposing

∈ S = ∈ S ∩ X X S ∈ S ∩ V V S ∈ S ∩ W W S
and collecting all variables z ∈ X S into a vector x, all variables z ∈ V S into a vector v, and all variables z ∈ W S into a vector w, system S uniquely defines an equation of the form [START_REF] Fritzson | Introduction to Modeling and Simulation of Technical and Physical Systems with Modelica[END_REF]

Call consistent causality of S any directed bipartite graph P S such that: 1) Its set of vertices is Z S, for some Z such that X S V S ⊆ Z ⊆ ∈ S ; 2) P S covers P S : for every z ∈ Z, z--E ∈ P S if and only if either z→E ∈ P S or E→z ∈ P S ; 3) P C has the single assignment property, meaning that, for every z, there exists at most one equation E such that E→z ∈ P C ; 4) Graph P S is circuitfree-hence its transitive closure is a partial order S on Z S, we call it the causality order induced by P S ; 5) X S ⊇ min( S ), the set of minimal elements of S ; 6) For every pair (z, w) ∈ (X S V S )×(W S ∩ Z), it never holds that w z. What does the existence of a consistent causality ensure? By condition 1), P S involves as vertices all variables belonging to X S (the "candidate inputs"), all variables belonging to V S (the "outputs"), plus some additional variables belonging to Z. By conditions 2) and 4), P S amounts to selecting a consistent orientation for P S , making it circuitfree: this defines the pivoting order between variables and equations, since, by Condition 3, at most one equation defines each variable. Condition 5) states that inputs for the system are found within X S : since they are minimal for the order S , they get evaluated first (when reading the values of the inputs). Finally, condition 6) expresses that no variable belonging to X S V S requires the prior evaluation of a variable belonging to W S , for its own evaluation. Thus eliminating the w variables simply proceeds by discarding certain equations.

The above discussion is formalized by the following theorem, where we associate with S the equation of the form [START_REF] Fritzson | Introduction to Modeling and Simulation of Technical and Physical Systems with Modelica[END_REF] it defines, together with its pair (A, C) of matrices:

Theorem 1: System S possesses a consistent causality if and only if the pair (A, C) is structurally nonsingular.

Proof: See Appendix A3.

The above results apply to the nonlinear equation ∃w.F (x, v, w) = 0 by considering the Jacobian ∇F at a solution (x o , v o , w o ) of F (x, v, w) = 0.

III. INDEX OF DIFFERENCE ALGEBRAIC EQUATIONS

In this section we translate the theory of differentiation index to discrete time systems consisting of difference Algebraic Equations (dAE), the discrete time counterpart of DAEs. The translation simply consists in replacing the derivative operator by the forward shift operator. While the definition of the true index is a straightforward translation, getting graph based algorithms is much more involved, as we shall see.

A. Difference Algebraic Equations (dAE)

A difference Algebraic Equation (dAE) consists of a constraint relating tuples of variables x and x • :

(x, x • ) ∈ C, where C ⊆ D×D (24)
meaning that D is the domain of both x and x • .

Definition 5: A solution of dAE (24) is any infinite sequence {x k | k ∈ Z} satisfying ∀k ∈ Z : (x k , x k+1 ) ∈ C
and ( 24) is solvable if solutions for it exist.

Definition 5 makes dAE [START_REF] Campbell | The index of general nonlinear DAEs[END_REF] shift-invariant and expresses that

x • is the forward shifted version of x: x • k = x k+1 (25)
For convenience, we write in the sequel

C(x, x • ) (26)
instead of [START_REF] Campbell | The index of general nonlinear DAEs[END_REF].

B. The true index of a dAE

Referring to (9), we define the kth difference array equations associated to (26), which collects forward shifted versions of constraint C(x, x • ):

     C(x, x • ) C • (x, x • ) . . ., x • ) C •k (x, x • )      = def      C (0) (x, x • , w) C (1) (x, x • , w) . . ., x • , w) C (k) (x, x • , w)      = def C k (x, v, w) (27)
where

v = def x • (28) w = def (x •2 , . . . , x •k+1 ) (29) and x •k+1 = def (x •k ) • A value x is called consistent for (27) if there exists (v, w) such that C k (x, v, w) holds, (30)
seen as an algebraic equation. Given a consistent value x for (27), algebraic equation (30) will generally have a set of solutions for (v, w). Writing w k instead of w in (27), the chain of sets

V k = def {v | ∃w k : C k (x, v, w k )} is decreasing for set inclusion.
Having finite index for the considered dAE means that this chain becomes a singleton for some finite value of k, and then remains so. The true difference index of a dAE system is defined with reference to the true differentiation index for DAE systems:

Definition 6: Assume that dAE (24) is solvable. The true difference index of this dAE, denoted by ν d , is the smallest index k such that v is uniquely determined by the algebraic equation (30) for any consistent value x for (27).

That is, for k ≥ ν d , the map

x → ∃w.C k (x, v, w) (31)
defines v as a deterministic function of x. Since ( 27) is equivalent to the original dAE, (31) determines x • and is, therefore, an OdE (Ordinary difference Equation), i.e., a transition system that can be directly executed. In (31), quantification over time applies as in [START_REF] Gawthrop | Metamodelling: bond graphs and dynamic systems[END_REF], namely, (31) expands as

∀n ∈ N : x n → ∃w n .C k (x n , v n , w n ) (32)
Again, we are rather interested in the structural difference index, simply referred to as difference index in the sequel, which is the almost everywhere value of the true difference index when the equations have their parameters modified while preserving the involvement/non-involvement of the variables in the equations.

C. Constructive semantics

Suppose that C has > 0 index, i.e., is not an OdE (a transition system in the usual sense). Then C k involves in its second block-row the constraint C • and thus, condition (31) ensures that x • is also consistent. That is, if the system has a finite index k, considering the array C k ensures that, at each reaction, current values for the variables and next values for the states are selected in a way that no deadlock will be caused in the future. For discrete time systems, the index appears as the needed look ahead horizon to be taken into account when generating the successives transitions of the system. The execution scheme of C is then specified by the following constructive semantics:

Constructive Semantics 1: 1) Initial condition: find a consistent x 0 , i.e., such that ∃w, ∃v : C k (x 0 , v, w) holds.

2) Non terminating while loop: for every n ≥ 0, a) assuming x n consistent, find the unique v n solution for v of ∃w :

C k (x n , v, w) holds; b) set x n+1 = v n , which is consistent by construction; c) repeat.
2 Constructive Semantics 1 specifies how runs of C are effectively constructed while time progresses, assuming a (static) constraint solver at hand. So far Execution scheme 1 is global, so that the entire burden is on the shoulders of the constraint solvers. In the next sections we show how calls for constraint solvers can be made "local" for dAE specified as systems of equations, by developing causality analyses related to the notion of structural difference index.

D. Index of smooth dAE

By the implicit function theorem, for the special case in which constraint

C k (x, v, w) has the form F k (x, v, w) = 0 for F k smooth, (31) is, locally around a triple (v o , w o , x o ), equivalent to the following condition: the map x → ∃w.Av + Cw + Ex = 0 defines v as a deterministic function of x (33) RR n°8630
where A, C, E are the Jacobians

A = (F k ) v , C = (F k ) w , E = (F k ) x , at the considered triple (v o , w o , x o ).
Hence Lemma 2 can be invoked for a structural checking of (33) based on causality analyses, from which the (structural) difference index and then the constructive semantics follow. We do not develop details for this since our interest is rather in non-smooth dAEs.

E. Index of non-smooth dAE: examples

For smooth systems, index definition and analysis reduces to the structural analysis of Jacobians. Graph based algorithms rely on the so-called "incidence graph" consisting of branches linking each equation to the variables involved in it. The algorithms search for a consistent direction for each branch, making the graph directed and circuitfree. Theorem 1 provided the justification of the graph based algorithm of Section II-C for index evaluation.

In our study of hybrid systems, we will need to consider dynamical systems with modes. A mode is characterized by a predicate over the system variables. A hybrid system has a smooth dynamics in each of its different modes-the reader is referred to Section VI for a detailed description of this. The consequence is that we cannot restrict ourselves to smooth systems. The following simple examples illustrate the new difficulties arising with the consideration of non-smooth systems.

1) Guarded equations: Consider the following system of static equations:

E 0 : b = [x > 0] E 3 : if b then F 1 (u, v) = 0 else F 2 (u, v) = 0 (34)
and call E 1 and E 2 the equations F 1 (u, v)=0 and F 2 (u, v)=0, respectively. Observe that E 0 is a function whose output is b.

Then, a possible consistent causality for dAE system (34) in the sense of Definition 4 consists of the following directed branches:

x→E 0 →b b→E 3 u→E 3 →v (35)
(Exchanging u and v in the last causality would do as well.)

The last two directed branches shown in the third line of (35) are legitimate according to Definition 4. Now, assume that equations E 1 and E 2 have indeed the form E 1 : v=f 1 (u) and E 2 : u=f 2 (v), where f 1 and f 2 are functions. Since equations E 1 and E 2 are indeed functions whose respective inputs are u and v, it is not true that E 3 will uniquely determinate v when the values of b and u are given. This is only valid if guard b is true. Thus the directed graph (35) does not help finding the index.

What was wrong? The point is that, for non-smooth equations, it is not true that any variable involved in it can be assigned the status of output, even structurally. How can we fix this? The solution is obvious for this case: replace (35) by the following guarded consistent causality, in which predicate b = [x>0] acts as a guard:

x→E 0 →b b→E 3 if b then u→E 3 →v else v→E 3 →u (36)
The actual consistent causality depends on the value of guard b. Observe that, in this guarded consistent causality, b is evaluated prior to being used as a guard. So, our first observation is that causality analyses must be guard dependent.

2) Unilateral constraint: A second interesting example is that of a unilateral constraint:

E 1 : x • = f (x, u) E 2 : 0 ≤ g(x) (37)
No consistent causality can be found for (37). The first difference array associated with (37) consists in adding the latent constraint E • 2 , obtained by shifting E 2 forward:

E 1 : x • = f (x, u) E 2 : 0 ≤ g(x) E • 2 : 0 ≤ g(x • ) (38) Substituting x • by f (x, u) in E • 2 , which we rename E 3 , yields by setting h(x, u) = def g(f (x, u)): E 1 : x • = f (x, u) E 2 : 0 ≤ g(x) E 3 : 0 ≤ h(x, u) (39)
The consistent causality associated with (39) is thus the following, where equations have been reordered to match causality:

x consistent for E 2 (40) x→E • 2 →u (41) (x, u)→E 1 →x • (42)
where "x consistent for E 2 " means that 0 ≤ g(x) and (x, u)→E 1 means x→E 1 , u→E 1 ; (41) abstracts the unilateral constraint E • 2 acting on u and requires a solver. A valid execution scheme at each time step for (39) would be the following one, formalized as a constructive semantics:

Constructive Semantics 2: 1) (40): assume x consistent for E 2 ; 2) (41): use a suitable static constraint solver proposing a u such that h(x, u) ≥ 0; 3) (42): having (x, u) compute x • , which is guaranteed consistent for the next time step.

3) The need for atomic sets of equations: Here we continue the discussion of the example of unilateral constraint in Section III-E2 by developing a special but important point.

So far unilateral constraints are not part of our syntax of guarded equations. We could consider including them. There are, however, many more candidate primitives to be includedfor instance the so-called complementarity constraints we discuss in the sequel. To avoid expanding our syntax with too many primitives, can we instead expand unilateral constraints in terms of guarded equations? The answer is: yes! g(x)≥0 expands as

   E 21 : b = [g(x)≤0] E 22 : if b then g(x)=0 F22 (43)
Unfortunately, the causality constraints induced by the set {E 21 , E 22 } of equations has a circuit:

x→E 21 →b→E 22 , if b then F 22 →x (
{E 21 , E 22 }→(x, b) (45)
We now apply this technique to system (39), which expands as follows:

E 1 : x • = f (x, u) E 21 : b = [0 ≥ g(x)] E 22 : if b then 0 = g(x) E 31 : c = [0 ≥ h(x, u)] E 32 : if c then 0 = h(x, u) F32(x,u) (46)
An attempt for the consistent causality of (46) is:

(x, u)→E 1 →x • (47) x consistent for {E 21 , E 22 }, E 21 →b (48) (x, u)→E 31 →c→E 32 , if c then (x, u)→F 32 →u (49)
Causality constraints collected in (49) were derived by abstracting each individual equations as its causality constraints. Subset of equations {E 31 , E 32 } exhibits a causality circuit, so we consider it as atomic. Applying to (49) the same transformation as we did from (44) to (45), yields:

x→{E 31 , E 32 }→(u, c) (50)
thus resulting in the following consistent causality for (46):

(x, u)→E 1 →x • (x, b) consistent for {E 21 , E 22 } x→{E 31 , E 32 }→(u, c)
which coincides with (40)-(42) if we abstract b and c away. In reinforcing (49) to (50), we decided to have x as the source and u as a target. We could have equally well done the opposite:

u→{E 31 , E 32 }→(x, c) (51)
In other words, the atom {E 31 , E 32 } is assigned two possible causality constraints, namely

x→{E 31 , E 32 }→(u, c) or u→{E 31 , E 32 }→(x, c). ( 52
)
Which one is ultimately kept depends on the other causality constraints, when attempting to avoid circuits. To summarize this discussion, when dealing with a minimal circuit in the set of causality constraints, the following must be performed: 1) propose for it a set of candidate circuitfree causality constraints; 2) taking the causality constraints induced by the whole set of equations into account, select a compatible candidate (causing no global circuit), if any. In this work, we are not proposing a systematic procedure for implementing Task 1). We can see 1) as manual task to be performed for a given library of primitive constraints (e.g., the unilateral constraint and the complementarity constraint discussed below). Grouping to an atom the equations involved in a circuit and applying the above procedure allows invoking "local solvers" as part of the execution of dAE.

4) Complementarity condition:

The following is a simple case of non-smooth dAE system, namely a dynamical system subject to a so-called complementarity condition, usually written as 0≤U (x)⊥V (y)≥0:

U (x)≥0 and V (y)≥0 and U (x)V (y) = 0 F (x, y) = 0 (53)
where x and y have dimensions possibly >1, U and V are Rvalued functions, and the constraint F (x, y)=0 is assumed to make the overall system (53) nonsingular. Such systems are frequently encountered in circuits (e.g., perfect diodes) and mechanics (e.g., dry friction, or contact). The complementary condition is a multi-mode constraint. We reformulate example (53) using guarded equations as follows:

E 1 : 0 = F (x, y) E 21 : b U = [U (x)>0] E 22 : b V = [V (y)>0] E 23 : if b U then 0 = V (y) else 0 = U (x) E 24 : if b V then 0 = U (x) else 0 = V (y) (54)
Regard E 2 = def 1≤i≤4 E 2i as an atomic system. Its consistent causality is derived similarly as for system (46), where E U and E V denote the equations 0=U (x) and 0=V (y), respectively:

x→E 21 →b U y→E 22 →b V b U →E 23 , if b U then E V →y else E U →x b V →E 24 , if b V then E U →x else E V →y (55)
In (55) we are missing the important side information that b U ∧ b V = F. We add it as an assertion:

assert b U ∧ b V = F x→E 21 →b U y→E 22 →b V b U →E 23 , if b U then E V →y else E U →x b V →E 24 , if b V then E U →x else E V →y (56)
We are now ready to include E 1 in our causality analysis.

The entire system (54) is now considered atomic, with the RR n°8630

following set of causality constraints:

assert b U ∧ b V = F x→E 21 →b U y→E 22 →b V b U →E 23 , if b U then E V →y→E 1 →x b U →E 23 , if b U else E U →x→E 1 →y b V →E 24 , if b V then E U →x→E 1 →y b U →E 23 , if b U else E V →y→E 1 →x (57)
Summary of this informal discussion:

• We will systematically expand unilateral constraints into guarded equations; • While doing so, we may need to "glue together" some equations and see them as an atomic system; • Atomic systems are given at once a consistent causality.

This consistent causality may involve assertions relating guards, which may somtimes be abstracted as a nondeterministic choice between alternatives. In the next section we formalize the above approach.

IV. INDEX OF NON-SMOOTH DAE: THEORY

In a first stage, we ignore guards. So the causality analyses we develop will be guardless. This will raise difficulties as we shall see. Thus, in a second stage, we account for guards.

A. Guardless causality analysis

We assume an underlying set of variables, together with the generic notation for its elements:

x ∈ X For x ∈ X , we denote by D(x) the domain of x and we set D(X ) = x∈X D(x). Then, for k any positive integer, we consider a fresh copy

X •k of X and we write X • instead of X •1 ; the elements of X •k are x •k . In particular, D(x • ) = D(x) and D(X • ) = D(X ).
Next, we assume an underlying set

E of primitive constraints E ⊆ D(X )×D(X •
), which we call equations in the sequel. Some subsets of this set of equations can be specified as being atomic. Atomic subsets are assumed to be pairwise disjoint. Atoms consisting of a singleton are identified with the single primitive equation they contain. Thus, for our theoretical development, E will be an underlying set of atoms generically denoted by the symbol E. For convenience, the term "equation" will often be used instead of "atom" in the sequel.

The dAE systems we consider are finite conjunctions of equations C = i∈I E i . For convenience, we shall often identify system C with its set {E i |i∈I} of constitutive equations.

Write

x ∈ E to indicate that x occurs in E and denote by

∈ E
the set of all variables x such that x ∈ E. If E is an equation of the form F = 0 where F is smooth, x ∈ E is equivalent to saying that the partial derivative ∂F/∂x is structurally nonzero. For E ∈ E and x ∈ E, we consider the predicate

x ∈ out E (58)
to specify that the value of x must be known prior to solving equation E (intuitively, x cannot be seen as an output of equation E).

Comment 1: Predicate (58) is not a semantic property of equation E; it is rather an additional specification that must be provided together with the semantics of E. For example, for E : x+y=0, neither x nor y satisfies (58). For E : b=[x > 0], we specify x ∈ out E. For the guarded equation

E : if b then F 1 (u, v)=0 else F 2 (u, v)=0
we specify b ∈ out E. For the atom E = def {E 31 , E 32 } of system (46), we specify c ∈ out E, which is equivalent to the set (52) of two candidate causalities.

2

For C a dAE system, define

x ∈ C if and only if x ∈ E holds for some E∈C and denote by ∈ C the set of variables z ∈ X X • such that z ∈ C. Define, for k a non-negative integer:

X C = def ∈ C ∩ X V C = def ∈ C ∩ X • W C = def x •l | 1≤l≤k, x • ∈ V C , x ∈ X C ∪ x •l | 2≤l≤k, x ∈ V C X C = def X C V C W C (59)
Variables belonging to V C are the state variables of C and variables belonging to X C are its other variables. Variables belonging to W C are the additional variables resulting from shifting forward the former ones and integer k is the given bound for those shifts. The following definition is a variation of Definition 4 to account for specified causality constraints:

Definition 7: Let C be a dAE system and C = i∈I E i be its k-th difference array following (27). We associate to array C its Pantelides graph P C , which is a bipartite graph defined by the following rules, where we identify C with its set of constitutive equations {E i | i∈I}:

• Its set of vertices is X C C; • For every z ∈ X C , using notation [START_REF] Petzold | Differential algebraic equations are not ODEs[END_REF]:

E ∈ C and z ∈ E z--E ∈ P C Call consistent causality of C any directed bipartite graph P C such that 1) Its set of vertices is Z C, for some Z ⊆ X C such that X C V C ⊆ Z; RR n°8630
( C ) ⊆ X C ; 7) For every pair (z, w) ∈ (X C V C )×(W C ∩ Z), it never
holds that w C z. The index of dAE system C is re-defined as the minimal integer k such that its k-th array C possesses a consistent causality.

The novelty with respect to Definition 4 lies in condition 3), which accounts for the fact that, in searching for a consistent causality, some constraints are subject to causality specifications, see Comment 1. Definition 4 is justified by the following weakening of Theorem 1, which follows from a straightforward pivoting argument:

Theorem 2: Consider the following condition for the considered consistent causality: for every branch E→z ∈ P C , z is uniquely determined as an output of E, for any consistent choice of values for the other variables involved in E.

(60) Assuming (60), if C possesses a consistent causality, then property (31) holds, and thus C possesses a constructive semantics.

If the system is smooth, condition (60) holds structurally. The example (34)-( 36) in section III-E illustrated how condition (60) may be violated by non-smooth systems. See the discussion sitting between (35) and (36), where the use of guards was motivated. In the next section, we formalize this by refining our development with the consideration of guards. In passing, we describe the constructive semantics, which yields the execution scheme.

B. Guarded causality analysis

We assume a subset B⊂X of guards and we augment E with two new kinds of equations:

Assertions over guards: In the sequel, dAE systems are pairs consisting of a conjunction of possibly guarded equations, and a conjunction of assertions over the involved guards:

C = i∈I E i , j∈J assert P j (64)
To properly capture the intuition behind (62), we assume the following about this operation: Axiom 1: The following holds:

∀E ∈ E =⇒ [ if T then E ] ≡ E (65) ∀E ∈ E =⇒ [ if F then E ] ≡ (66) ∀E ∈ E =⇒ b(E) ∈ out E (67)
(65) states that, if b is true, equation "if b then E" reduces to E. (66) states that, if the guard b is false, then "if b then E" collapses to the trivial equation setting no constraint at all, i.e., having no variable involved in it: ∈ = ∅. Focus finally on condition (67) by making the form of E explicit: E: "if b then E". Then, condition (67) expresses that guard b must be evaluated prior to the equation it guards, see (58) for the definition of ∈ out . Observe that this is a nontrivial restriction, since it prevents if b then E from being a fixpoint in both b and the variables of E, jointly.

We now formalize the notion of guarded Pantelides graph, which provide the graph abstraction of guarded equations. For (b, x, E) ∈ B×X ×E, we consider the following guarded bipartite branch and guarded directed bipartite branch:

if b then x--E and if b then x→E or if b then E→x (68)
(We make no difference between non directed branches x--E and E--x.) We consider the counterpart of Axiom 1: Axiom 2: The following holds, for any (possibly directed) bipartite branch π:

[ if T then π ] ≡ π (69) [ if F then π ] ≡ (70)
Due to Axiom 2, we feel free to identify the non guarded branch π with the guarded branch "if T then π". For π = if b then x--E (respectively "if b then x→E") a branch, define its base:

π = def x--E respectively π = def x→E
and, for the directed branch, its reversal:

← π = def if b then E→x
Due to Axiom 2, we feel free to identify the non guarded branch π with the guarded branch "if T then π". Accordingly, for a non guarded branch π, its base π coincides with π itself.

A guarded (directed) bipartite graph P is a pair consisting of a finite set of guarded (directed) bipartite branches, and a finite conjunction of assertions over its involved guards:

P = {π i | i∈I} , j∈J assert P j (71)
In the sequel, we simply write π ∈ P to indicate that π is a branch of the graph part of P.

If P is directed, say that π 1 , . . . , π k is a path (circuit) of P if π 1 , . . . , π k is a path (circuit). For π 1 , . . . , π k a path in P such that π i = if b i then π i , we define its guard b(π 1 , . . . , π k ) by b(π 1 , . . . , π k ) = def 1≤i≤k b i (72)
and say that P is guarded circuitfree if all its circuits have a false guard-taking assertions of P into account. If P is circuitfree, then the transitive closure P of the graph part of P is a partial order, for any configuration of the guards. Say that P has the guarded single assignment property if the assertions of P imply that, for every z, ∈L b = F holds, where {b | ∈L} is the set of guards of all the branches "if b then E →z" of P ending at z. Using these prerequisites, Definition 7 is reformulated: Definition 8: Let C be a dAE system and C = i∈I E i , j∈J assert P j be its k-th difference array following (27), together with its set of assertions. We associate to array C its guarded Pantelides graph P C , which is a guarded bipartite graph keeping the assertions of C and otherwise obtained by applying the following rules, where we identify C with its set of constitutive equations {E i | i∈I}:

E ∈ C [ b(E)--E ] ∈ P C [ if b then E ] ∈ C and z ∈ E [ if b then z--E ] ∈ P C (73)
Call guarded consistent causality of C any guarded directed bipartite graph P C keeping the assertions of C and otherwise such that, with reference to notations (59):

1) Its set of vertices is Z C, for some Z ⊆ X C such that X C V C ⊆ Z;
2) P C covers P C : for every π, π ∈ P C either π ∈ P C or ← π∈ P C

3) Causality specifications are respected:

[if b then z--E] ∈ P C and z ∈ out E [if b then z→E] ∈ P C
4) Causality distributes over guards:

[ z→[ if b then E ] ] ∈ P C [ if b then z→E ] ∈ P C 5
(z, w) ∈ (X C V C )×(W C ∩ Z), it never
holds that w C z. The index of dAE system C is re-defined as the minimal integer k such that the k-th difference array C possesses a consistent causality. Observe that, due to condition (67) of Axiom 1, entailment 3), applied with b being the constant true and E replaced by E = "if b then E", yields:

E ∈ C [ b(E)→E ] ∈ P C (74)
(Compare with the first entailment rule in (73).)

C. Constructive semantics

The following result justifies the consideration of Definition 8:

Theorem 3: Assuming (60), if C possesses a consistent causality according to Definition 8, then property (31) holds. Said differently, Theorem 2 still holds with Definition 8 substituted for Definition 7. The important point about refining Definition 7 into Definition 8 is that condition (60), which ensures the correctness of the graphical abstraction, is much easier satisfying, see the paragraph sitting just before Section IV-B.

Proof: By condition 7) and entailment (74), vertices belonging to min( C ) can only be of two categories:

• either non guarded equations having no inputs and only outputs; such equations belong to min( C ) whatever configuration of the guards is; • or variables belonging to X C that are output of no equation; such variables belong to min( C ) whatever configuration of the guards is. In particular, the set min( C ) is independent of the configuration of the guards and, thus, condition 7) is indeed meaningful.

The evaluation of all variables of C at a considered instant proceeds according to the following constructive semantics:

Constructive Semantics 3: P is a running subgraph of P C and a running subrelation of C :

1) Initialization: P := P C and := C ; 2) While P has non-empty set of vertices, do: a) evaluate variables or equations belonging to min( ); the variables are evaluated by reading their values; thanks to Condition 5) of Definition 8 and condition (60), evaluating the equations fixes the values of their immediate sucessors in ; b) for those evalued variables that are guards, replace them by their value T or F and then apply simplifying rules (65) and (66) of Axiom 1; Condition 4 of Definition 8 ensures that no causality is lost by removing the trivial guards T; c) erase, from P, the vertices belonging to min( ) and adjacent edges and redefine accordingly. 2

The key observation is that, at each round of Algorithm 3, min( ) is independent of the configuration of the guards that remain to be evaluated. By Condition 2) of Definition 8, every equation and variable of array C gets evaluated. Finally, due to Condition 8) of Definition 8, unnecessary variables w are eliminated by discarding the equation defining them. This finishes the proof of Theorem 3.

Constructive Semantics 3 can be seen as an interpreter of system C.

Comment 3 (atoms):

Recall that some of the E's involved in the considered dAE system are atomic systems of equations representing constraints for which a solver is needed. The principle is that each such constraint E comes with its guarded Pantelides graph P E , possibly enhanced with a set of assertions on the guards of E and a set of specifications of the form x ∈ out E, see (58). By using the context of E, if a consistent causality can be found, a directed graph P E results, which in turn specifies the variables for which the solver must return a solution of E. See Section III-E3 for a detailed example.

Comment 4 (nested guards): Capturing multi-mode systems with arbitrarily nested modes would require that E is further extended to become closed under the mapping (b, E) → [if b then E], meaning that nested guards should be considered-something forbidden by our current extension of E, see (62). Extending Definition 8 and Theorem 3 to cover this is feasible but more involved. This is left for further work.

Comment 5 (multiple clocked dAE systems): So the dAE as defined in Section III-A are "single-clocked" in that all variables possess a value at every instant. It is well known from synchronous languages that multiple-clocked discrete time systems can be made multiple-clocked by adding, to every data type, a distinguished symbol "⊥" denoting the absence of the considered variable at the considered instant [START_REF] Benveniste | The synchronous languages 12 years later[END_REF]. In a given run, the subset of instants at which a given variables is present, i.e., =⊥, is called its clock. Clocks can be seen as a type system and these types can either be verified or synthesized [START_REF] Benveniste | The synchronous languages 12 years later[END_REF]. Of course, the same techniques apply to dAE systems. We did not consider this issue here since it would bring yet another layer of technicality.

D. Back to dAE examples

In this section we revisit the non-smooth dAE examples of Section III-E and we add the most basic example, namely the change of mode upon detection of a zero-crossing.

1) Zero-crossing:

The following example of continuoustime ODE with event-based reset is the simplest case of hybrid system, i.e., continuous time system with mode change:

if b then x = h(x -) else ẋ = f (x), (75)
where boolean signal b selects the events of zero-crossings of some smooth function g(x), i.e., at any instant when g crosses zero from below. The default dynamics is the ODE ẋ = f (x) and x = h(x -) yields the reset value for x at the instants of zero-crossing, where x -denotes the left-limit of trajectory x, i.e., the value x had just before the reset.

We consider here example (75), albeit in its dAE form. Three variations can be considered, depending on how shifts get substituted for the left-limit and derivative operators:

if b then x = h( • x) else x • = f (x) b = Q(x, x • ) (76) if b then x = h( • x) else x • = f (x) b = Q( • x, x) (77) if b then x • = h(x) else x • = f (x) b = Q( • x, x) (78)
where Q is the predicate specifying the zero-crossings of g:

Q(z, x) = def [g(z)≤0] ∧ [g(x) > 0] (79)
We successively analyze the above three variants.

Analyzing (76): Since (76) involves two successive instants, we first put it in state space form:

   z • = x if b then x = h(z) else x • = f (x) b = Q(x, x • ) (80)
which possesses no consistent causality and is, therefore, not of index 0. Denoting by E(b, z, x, x • ) the guarded equation of (80), the next difference array for consideration is the following, where the "----" indicate the separation between the original and shifted dynamics:

          0 = z • -x E(b, z, x, x • ) 0 = b -Q(x, x • ) --------- 0 = z •2 -x • E(b • , z • , x • , x •2 ) 0 = b • -Q(x • , x •2 )          
According to the definition of the 2-nd array, variables b • , x •2 , z •2 have to be eliminated. Unfortunately, considering this array and even further increasing it does not help finding a consistent causality for the pair (b, x • ). The index is thus found to be infinite. Let us try to fix this problem by specifying the predicate evaluation to be directed, see (58) and Comment 1 for the definition of ∈ out :

   0 = z • -x E(b, z, x, x • ) (x, x • ) ∈ out in b = Q(x, x • )
As a consequence, we cannot select x • as an output of the first equation, since this would yield a causality circuit. The conclusion is that (77) is a bad model for a dAE version of zero-crossing.

Analyzing (77): Repeating the same for (77) and using previous notation E(b, z, x, x • ) yields:

   E 1 : 0 = z • -x E 2 : E(b, z, x, x • ) E 3 : (z, x) ∈ out in b = Q(z, x) (81)
At this point we make an attempt to use guardless causality analysis, based on Definition 7. According to this definition, the following is a correct consistent causality for (81):

   x→E 1 →z • (b, z, x)→E 2 →x • (z, x)→E 3 →b • (82)
Consistent causality (82) is optimistic and erroneous, however, in that it does not satisfy the important condition (60). While it yields a seemingly correct order, it makes the assumption that x • can be determined as an output of equation E 2 . This is, however, wrong at zero-crossing events, i.e., when b = T. Thus, to guarantee that Theorem 3 applies, Definition 8 must be invoked and, thus, the following guarded Pantelides graph must be used instead (note the non-directed branches in the third and fourth lines):

           x→E 1 →z • b→E 2 if b then b→E 2,T --(z, x) if b else b→E 2,F --(x, x • ) (z, x)→E 3 →b (83)
where E 2,T and E 2,F denote the if-branch and else-branch of guarded equation E 2 . Unfortunately, no consistent causality can be found for (83) since every attempt is not circuitfree by exhibiting the circuit x→b→x having guard b, which is not false in general. Thus, (77) is not appropriate either.

Analyzing (78): A cascade of two successive events occurs, namely: the observation of b = T at some instant, followed by a reset of x at the next instant. Putting (78) in state space form yields:

   z • = x if b then x • = h(x) else x • = f (x) b = Q(z, x) (84)
for which the following consistent causality is found, showing that the index is 0:

       x→E 1 →z • b→E 2 if b then (b, x)→E 2,T →x • else (b, x)→E 2,F →x • (z, x)→E 3 →b
We can simplify the guarded branch by merging its two alternatives, which finally yields:

   x→E 1 →z • (b, x)→E 2 →x • (z, x)→E 3 →b
This is an example of the manipulations we can perform on guarded branches to avoid enumerating modes while computing the guarded incidence matrix of Pantelides' graph.

2) Unilateral constraint: The unilateral constraint in discrete time was studied in Sections III-E2 and III-E3. The developments performed there, exactly follow the theory of Definition 8, Theorem 3, and Constructive Semantics 3.

3) Complementarity condition : Complementarity condition (53) of Section III-E4 is reformulated as the following dAE system with guards and assertions:

E 1 : 0 = F (x, y) A 2 : assert b U ∧ b V = F E 21 : b U = [U (x)>0] E 22 : b V = [V (x)>0] E 23 : if b U then 0 = V (y) else 0 = U (x) E 24 : if b V then 0 = U (x) F24 else 0 = V (y) F23 (85)
With reference to (53), we have added the assertion that guards b U and b V are never true simultaneously. Clearly, this can be deduced from the equations defining the system. Still, we state it explicitly so it can be used when constructing the causality analysis.

No consistent causality can be found for it other than by viewing the entire system (85) but its first equation as an atom. Doing so is indeed consistent with (85) being itself, in the litterature of so-called non-smooth systems, seen as a primitive for which dedicated solvers are developed [START_REF] Vincent | Numerical Methods for Nonsmooth Dynamical Systems[END_REF]. Depending on the particular form for the first equation, we may need to extend (85) into a larger array, by shifting forward the atom.

The following set of causality constraints is found for the atom (A 2 , E 2 ) = def (A 2 , {E 21 , E 22 , E 23 , E 24 }):

A 2 : assert b U ∧ b V = F E 2 :            x→E 2 →(b U , b V ) b U →E 23 if b U then b U →F 23 →y else b U →F 24 →x b V →E 24 if b V then b V →F 24 →x else b V →F 23 →y
which boils down to

A 2 : assert b U ∧ b V = F E 2 :    x→E 2 →(b U , b V ) b U →E 23 , if b U else b U →F 24 →x b V →E 24 , if b V else b V →F 23 →y (86)
Combining (86) with equation E 1 yields the following set of causality constraints for the overall atomic system (85):

A 2 : assert b U ∧ b V = F {E 1 , E 2 } :    x→E 2 →(b U , b V ) b U →E 23 , if b U else b U →F 24 →x→E 1 →y b V →E 24 , if b V else b V →F 23 →y→E 1 →x
Since assertion A 2 implies ¬b U ∨ ¬b V = T, atomic system {E 1 , E 2 } entirely determines the pair (x, y) as its output. A correct abstraction for (86) is thus

{E 1 , E 2 }→(x, y)
Consider next the following variation of (85), in which E 1 has been modified so that x is now a state:

E 1 : x • = f (x, y) A 2 : assert b U ∧ b V = F E 21 : b U = [U (x) > 0] E 22 : b V = [V (y) > 0] E 23 : if b U then 0 = V (y) else 0 = U (x) E 24 : if b V then 0 = U (x) else 0 = V (y) (87)
This is no longer of index 0, so we must shift the complementarity condition. Shifting its causality analysis yields the following, where we have eliminated the unnecessary y • and taken into account that x is assumed consistent (and thus is no longer occurring in the causality analysis):

A 2 : assert b U ∧ b V = F A • 2 : assert b • U ∧ b • V = F E 2 :    x→E 2 →(b U , b V ) b U →E 23 , if b U then b U →F 23 →y b V →E 24 , if b V else b V →F 23 →y E • 2 :    x • →E • 2 →(b • U , b • V ) b • U →E • 23 , if b • U else b • U →F • 24 →x • b • V →E • 24 , if b • V then b • V →F • 24 →x •
Taking E 1 into account yields:

A 2 : assert b U ∧ b V = F A • 2 : assert b • U ∧ b • V = F (88) atom {E 1 , E 2 , E • 2 } :                x→E 2 →(b U , b V ) x • →E • 2 →(b • U , b • V ) b U →E 23 , if b U then b U →F 23 →y→E 1 →x • b V →E 24 , if b V else b V →F 23 →y→E 1 →x • b • U →E • 23 , if b • U else b • U →F • 24 →x • →E 1 →y b • V →E • 24 , if b • V then b • V →F • 24 →x
• →E 1 →y which determinates (y, x • ) as its ouput if the following condition holds:

b U ∨ ¬b V ∨ ¬b • U ∨ b • V = T (89)
Unfortunately, assuming assertion A 2 ∧ A • 2 still does not make (89) a tautology. Extending the array (88) will not solve the problem. This suggests that our assertion A 2 may not be strong enough. Let us strenghten A 2 as follows:

A 2 : assert b U ∧ b V = F assert [b U = b • U ] ∨ [b V = b • V ] = T ( 90
)
which amounts to assuming that the two unilateral constraints do not change their status saturated/unsaturated simultaneously. Then, assuming (90) makes (89) a tautology, showing that the atom (A 2 , {E 1 , E 2 }) has index 1 under (90).

V. INDEX REDUCTION AND NONSTANDARD SEMANTICS

So far we have developed a theory of index for dAE with associated graph based algorithms. A key step in our agenda of developing a theory of index for hybrid DAE systems is to cast smooth DAE systems into the world of dAEs and see what our dAE index theory gives. Our aim is to show that the two ways of defining an index for this class of systems actually coincide. Casting smooth DAE systems into the world of dAEs is performed by interpreting DAE using nonstandard analysis [START_REF] Robinson | Nonstandard Analysis[END_REF], [START_REF] Diener | Analyse non standard[END_REF], [START_REF] Lindstrøm | An invitation to nonstandard analysis[END_REF], [START_REF]Nonstandard analysis and its applications[END_REF], [START_REF] Bliudze | Un cadre formel pour l'étude des systèmes industriels complexes: un exemple basé sur l'infrastructure de l'UMTS[END_REF].

A. Nonstandard analysis for the engineer

Here we collect the minimum needed to follow the subsequent developments (excluding the proofs). Corresponding details are found in the Appendix B. We will use N and R, which are the nonstandard extensions of N and R, respectively. N is a totally ordered set obtained by extending N with infinite integers, greater than any standard finite integer. R is obtained by extending R in two ways: 1) by adding infinitesimals, which are nonzero elements closer to zero than any standard (ordinary) nonzero real, and 2) infinite reals, greater than any standard finite real; to every finite real x we then add its nonstandard neighbours x+ where is any infinitesimal. R is totally ordered, and usual operations (+, ×, etc.) and relations on reals extend to it. For any finite positive real t and any finite nonstandard positive real ∂, there exists n ∈ N such that n×∂ ≤ t < (n + 1)×∂. This applies in particular when ∂ is infinitesimal, in which case n must be an infinite nonstandard integer.

This suggests to consider the following set as a time domain:

T = {k∂ | k ∈ Z} RR n°8630
where ∂ is an infinitesimal time step and Z is the set of nonstandard positive or negative integers. Elements of T will be denoted by the symbol τ . Time domain T is totally ordered and ranges from minus infinity to plus infinity. It is both "continuous": any finite real standard number t has an element τ ∈ T such that |t -τ | is infinitesimal; and "discrete": each τ = n×∂ ∈ T has a unique immediate predecessor • τ = def (n -1)×∂ and immediate successor τ • = def (n + 1)×∂.

B. Nonstandard semantics of DAEs

The nonstandard semantics of a DAE is obtained by applying the following substitution rules:

ẋ ← 1 ∂ (x • -x) ẍ ← 1 ∂ 2 (x •2 -2 x • + x) x (3) ← 1 ∂ 3 (x •3 -3 x •2 + 3 x • -x) x (4) ← . . . (91)
where x • , x •2 , x •3 , . . . are defined as follows (we will need both the pre operator • x and post operator x • in the sequel):

• τ = def τ -∂ , τ • = def τ + ∂ • x τ = def x• τ , x • τ = def x τ • •m x = def • ( •m-1 x) , x •m = def (x •m-1 ) • (92)
Applying (91) to the pendulum example (15) yields:

x • = x + ∂ × u u • = u + ∂ × T x y • = y + ∂ × v v • = v + ∂ × (T y -g) L 2 = x 2 + y 2 (93)

C. The two notions of index coincide

We now establish a link between index reduction for DAE and index reduction for dAE. Return to the pendulum example. Highest degree shifted variables are x • , u • , y • , v • , T . Corresponding Jacobian is singular, thus the difference degree is strictly positive. Forward shifting the last equation two times yields:

x • = x + ∂ × u (i1) u • = u + ∂ × T x (i2) y • = y + ∂ × v (ii1) v • = v + ∂ × (T y -g) (ii2) L 2 = x 2 + y 2 (iii) L 2 = (x • ) 2 + (y • ) 2 (iv) L 2 = (x •2 ) 2 + (y •2 ) 2 (v) (94)
Substituting, in (iv, v), x • and y • by using (i1) and (ii1), and reorganizing the result yields:

x • = x + ∂ × u (i1) u • = u + ∂ × T x (i2) y • = y + ∂ × v (ii1) v • = v + ∂ × (T y -g) (ii2) L 2 = (x • + ∂ × u • ) 2 + (y • + ∂ × v • ) 2 (v) ---------- L 2 = x 2 + y 2 (iii) L 2 = (x + ∂ × u) 2 + (y + ∂ × v) 2 (iv) (95)
The first group of equations has structurally nonsingular Jacobian with respect to x • , u • , y • , v • , T , and is thus a dAE of index 0. dAE system (93) has thus index 2. All of this suggests that we have the following "pseudo-equation":

index(NS(DAE)) structurally = index(DAE)
where NS(DAE) denotes the nonstandard interpretation of DAE, seen as a dAE. Let us formalize this. In the following, n denotes the dimension of variable v.

Definition 9: Say that the following property holds:

x → ∃w.F (x, v, w) = 0 structurally defines v (96)
if the pair (F v (x, v, w), F w (x, v, w)) of Jacobians is structurally nonsingular, see Definition 3 and Lemma 2.

Informally, if we perturb F while respecting the independence of a given equation with respect to a given variable, the map x → ∃w.F (x, v, w) = 0 defines v almost everywhere. Next, consider

• the differential array F k (x, v, w) associated with DAE F (x, ẋ) = 0 following [START_REF] Campbell | A computational method for general higher index nonlinear singular systems of differential equations[END_REF], and • the difference array F NS k (x, v, w) associated with F NS (x, x • ) = 0 following (27), where F NS (x, x • ) is the nonstandard semantics of F (x, ẋ). Let us relate these two arrays. In doing so we use different notations for the variables involved in F , written x, v, w, and the variables involved in F NS , written x, v, w-we do not need to distinguish the x-variable for the two cases. Coordinates of w are w 2 . . .w k and similarly for w. We have, with reference to [START_REF] Campbell | A computational method for general higher index nonlinear singular systems of differential equations[END_REF] and ( 27):

F NS k (x, v, w) =       F (0) NS (x, v, w) F (1) NS (x, v, w) . . . F (k) NS (x, v, w)       F k (x, v, w) =      F (0) (x, v, w) F (1) (x, v, w) . . . F (k) (x, v, w)      (97)
Using (91), we get that (x, v, w) and (x, v, w) are related by the following formulas, where

i m = i! m!(i -m)! denote the binomial coefficients: v = 1 ∂ (v -x) w 2 = 1 ∂ 2 ( w2 -2 v + x) w 3 = 1 ∂ 3 ( w3 -3 w2 + 3 v -x) . . . w i = 1 ∂ i i j=0 (-1) i-j i j wj (98)
with the convention that w1 = v and w0 = x. Then,

F (0) (x, v, w) = def F (x, v) = F (x, 1 ∂ (v -x)) = def F (0) NS (x, v, w) . (99)
F (i) (x, v, w) is obtained by differentiating i-times F (x, ẋ)
with respect to time and then replacing the successive derivatives of x by v, w 2 , w 3 , . . . , w i . Similarly,

F (i) NS (x, v, w) is ob- tained by shifting forward i-times F (x, 1 ∂ (x • -x))
and then replacing the succesive forward shifts of x by v, w2 , w3 , . . . , wi . Hence, F NS k (x, v, w) and F k (x, v, w) are related by

F (0) (x, v, w) = F (0) NS (x, v, w) F (1) (x, v, w) = 1 ∂ F (1)
NS (x, v, w) -F

NS (x, v, w)

F (2) (x, v, w) = 1 ∂ 2 F (2)
NS (x, v, w) -2F

NS (x, v, w)

1 ∂ 2 ( + F (0)
NS (x, v, w) . . .

F (i) (x, v, w) = 1 ∂ i i j=0 (-1) i-j i j F (j) NS (x, v, w) (100)
Accordingly, if we define the square lower triangular matrix matrix L(p) by for j≤i :

L(p) ij = (-1) i-j ∂ i i j
and set

L = L(k) ⊗ I n (101)
where ⊗ denotes the Kronecker product and I n is the n×nidentity matrix, 6 the following holds:

Lemma 3: The matrix L defined in (101) is block lower tri- 6 The Kronecker product A⊗B is the block matrix whose ijth block is

a ij B.
angular (BLT), invertible and with BLT inverse. Furthermore,

       v w 2 w 3 . . . w k        = L        v w2 w3 . . . wk        +        -1/∂ 1/∂ 2 -1/∂ 3 . . .        x (102)           F (1) NS F (2) NS F (3) NS
. . .

F (k) NS           = L           F (1) NS F (2) NS F (3) NS
. . .

F (k) NS           +           -1/∂ 1/∂ 2 -1/∂ 3 . . . . . .           F (0) NS (103)
On the other hand, the map

(x, v, w) → F k (x, v, w)
factorizes as the following composition of maps:

(x, v, w) (102) → (x, v, w) → F NS k (x, v, w) (103)
→ F k (x, v, w)
Using the formula giving the Jacobian of the composition of maps, we get

∇ x,v,w F k = I n 0 0 L ∇ x,v, wF NS k I n 0 0 L -1
where the notation ∇ x,v,w f stands for the Jacobian of f with respect to the triple (x, v, w).

Theorem 4: The following two properties are equivalent:

x → ∃w.F k (x, v, w) = 0 structurally defines v (104) x → ∃w.F NS k (x, v, w) = 0 structurally defines v (105)
Proof: By Lemma 2, property (104) holds if and only if 21), for some two permutation matrices P and

A C = def ∇ v F k ∇ w F k satisfies (
Q. It turns out that, if L is block-lower triangular, [ A C ] satisfies (21) if and only if so does L[ A C ]L -1 .
Using an implicit Euler scheme: So far the substitution rules (91) correspond to applying an explicit Euler scheme with an infinitesimal step ∂, see for example (93). What happens if we use an implicit Euler scheme instead, meaning that the nonstandard semantics of a DAE becomes

F ( ẋ, x) = 0 ← F ( 1 ∂ (x • -x), x • ) = 0 (106)
For example, ODE ẋ = f (x) yields x • = x+∂×f (x • ). Using (106) instead of (91) amounts to replacing, in (100), x by x • while leaving the rest unchanged. The same substitution holds in (102) while (103) remains unchanged. The substitution x ← x • in (102) does not affect the block-triangular matrix L. The bottom line is that Lemma 3 and thus also Theorem 4 still hold if an implicit Euler scheme is used.

VI. HYBRID DAE In this section we apply our previous results fo Hybrid DAE systems. We first describe the class of hybrid DAE systems we consider. Then, we describe a set of primitives that is sufficient to specify hybrid DAE systems. Finally, we study their index.

A. Trajectories

We collect in a preamble the needed notations regarding trajectories of hybrid systems. The reader is referred to Figure 1. This figure shows the trajectory of a variable x of type real. We will be interested in a restricted class of hybrid systems, in which:

dom(x) = R real time t (t, 1) (t, 0) (t, 2) (v, 0) (v, 1) (v, 3) (u, 0) v (v, 2) u
• Variables of type real possess trajectories such as in Figure 1, meaning that instants of discontinuity are isolated (t and v in the figure); at t the dashed trajectory hits the red bullet labeled with (t, 0), and then jumps to the red bullet labeled with (t, 1) and then (t, 2); similarly, at time v, the trajectory has its second discontinuity, by successively hitting (v, 0), (v, 1), (v, 2), and (v, 3). the trajectory is continuous in the closed interval [t, v], with a beginning depicted by the red bullet labeled with (t, 2) and an end depicted by the red bullet labeled with (v, 0). During the open interval (t, v), points on the trajectory are labeled with time index (u, 0), indicating that the associated date of occurrence is u (thus, the second field 0 of the date serves nothing).

• The same holds for all variables having a domain equipped with a topology. In particular, for discrete domains (such as Booleans or Integers), the same holds if the domains are equipped with the discrete topology. Formally, we use as time index the set S = R × N, equipped with the lexicographic order defined by: (t, k) < (t , k ) if and only if, either t < t , or t = t and k < k . Elements of S are denoted by the symbol s, or explicitly as pairs s = (t, k) whenever needed. By convention:

We identify R and the subset {s = (t, 0) | t ∈ R} ⊂ S. (107) Time set S defines the so-called super-dense time, see [START_REF] Lee | Operational semantics of hybrid systems[END_REF], [START_REF] Lee | Leveraging synchronous language principles for heterogeneous modeling and design of embedded systems[END_REF]. It allows defining solutions for (110) in which finite (but possibly unbounded) cascades of mode changes can occur.

To every trajectory, we associate a function R t → n t ∈ N such that n t = 0, except for t belonging to some increasing sequence T = {t k | k∈Z} of instants of R with lim k→±∞ t k =±∞. For the trajectory shown in Figure 1, we have n t = 2, n v = 3, and otherwise n takes the value zero.

For s ∈ S, define x - s and x + s by

x - s = def lim t t x (t ,0) if s = (t, 0) x (t,k-1) if s = (t, k), k > 0 (108) x + s = def lim t t x (t ,0) if s = (t, n t ) x (t,k+1) if s = (t, k), k < n t (109)

B. Mode dependent dynamics

In this section we use the guarded equations introduced in Section IV-B and the reader is referred to this section for the notations. The basic form for a hybrid DAE system is the following, where I is some finite set of modes and i ranges over I:

S : ∀i ∈ I b i = B i (x -, x, ẋ) if b i then 0 = F i (x -, x, ẋ) (110)
where I is some finite index set, x denotes a n-tuple of real variables, denote the left-and right-limit of x at instant t, the F i 's are real-valued and smooth, and the B i 's are smooth boolean predicates over the listed variables. The first equation specifies that, in mode i, DAE F i = 0 must hold; observe that this is a fixpoint equation. Thus, b i is the guard of mode i and F i = 0 its dynamics. If two or more guards overlap, then the conjunction of the corresponding dynamics must hold. If the "or" of all guards is not the constant "true", then S is incompletely specified. Although we do not require that the different modes satisfy i∈I b i = T and b i ∧b j = F for any two i =j, this will generally be the case in practice.

The generic form (110) does not allow for nested guards, e.g., of the form if b then [if b then E]. Flattening can be applied to bring nested guards to the form (110). Flattening, however, is known not to be desirable from the point of view of computational complexity, as it may typically give raise to combinatorial explosion. Nested guards are thus desirable to include but are left for future work.

Examples:

We discuss here a few examples, showing the flexibility of generic form (110):

A first example is obtained by considering a DAE system with unilateral constraint:

0 ≤ F (x, ẋ) (111)
where x and F are as above. Equation ( 111)

rewrites b = [ 0≥F (x, ẋ) ] if b then 0=F (x, ẋ)
This system has two modes. In the first mode, b = T expresses that unilateral constraint (111) is active. In the second mode corresponding to the else alternative, the constraint is trivial, RR n°8630

expressing that (x, ẋ) is unconstrained when the inequality is strict in (111). So-called DAE systems with a complementarity conditions are a second example: U (x)≥0 and V (y)≥0 and U (x)V (y) = 0 F (x, y) = 0 (112)

where we assume that the second equation ensures that the system (112) is nonsingular. The first equation of (112)the complementarity condition-is often written in a compact form as

0 ≤ U (x) ⊥ V (y) ≥ 0 (113)
Such systems are encountered, e.g., in electric circuits with perfect diodes. System (112) can be given the generic form (110) in the following way:

0 = F (x, y) b U = [U (x)>0] b V = [V (x)>0] if b U else [0=U (x)] if b V else [0=V (y)] if b U then [0=V (y)] if b V then [0=U (x)] (114)
A third example is the zero-crossing example (75) of Section IV-D1, where boolean guard b selects the events of zero-crossings of some smooth function g(x).

2

We now define what a solution of (110) is. Let B denote the Boolean domain.

Definition 10: Hybrid DAE (110) is solvable if there exists a pair of functions (s, λ) → (Φ(s, λ), β(s, λ)), from S×Λ into (R n ×B I )∪{ }, where is the undefined value and Λ is some nonempty open set of R ρ , satisfying the following conditions, where β i , i∈I denote the components of β:

1) For every λ, there exists a function R t → n t (λ) ∈ N, such that (Φ((t, k), λ), β((t, k), λ)) = if and only if k > n t (λ). Furthermore, n t (λ) = 0, except for t belonging to some increasing sequence T (λ) = {t k (λ) | k∈Z} of instants of R with lim k→±∞ t k (λ)=±∞. We write t k instead of t k (λ) when no confusion results. 2) Regarding the modes:

• The function t→β((t, 0), λ) is constant over every open interval (t k , t k+1 ), and, letting b ∈ B I be the corresponding value, we have β((t k , n t k ), λ) = b and β((t k+1 , 0), λ) = b. 3) Regarding the state x, and using (108), ( 109) and (107):

• For every λ ∈ Λ and every interval (t k , t k+1 ), t→Φ(t, λ) is a diffeomorphism from (t k , t k+1 ) into R n , and, if β i (t, λ) = T, then F i Φ -(t, λ), Φ(t, λ), d dt Φ(t, λ) = 0 holds for every t ∈ (t k , t k+1 ).

• For every λ ∈ Λ and every t∈T (λ), then

Φ((t, 0), λ) = Φ -((t, 0), λ) Φ((t, n t ), λ) = Φ + ((t, 0), λ)
and, for every k>0 and i such that β i ((t, k), λ)=T, then x i = def Φ((t, k), λ) is a consistent value for F i (x -, x, ẋ) = 0. 4) If s→(x s , b s ) satisfies conditions 1) -3) above, then it holds that (x s , b s ) = (Φ(s, λ), β(s, λ)) for some λ.

Some comments are in order regarding Definition 10, with reference to its successive conditions:

1) The function n t (λ) speficies, for instant t, whether the system traverses a cascade of mode changes (n t (λ) > 0 being the length of this cascade), or not (n t (λ) = 0).

Cascades are assumed to be finite (but not necessarily uniformly bounded) and isolated from each other. 2) A mode begins at the last instant of the last seen cascade (where it is reset), and ends at the first instant of the next cascade. 3) The system dynamics at a persistent mode i (meaning that b i = T for some positive duration) is F i . 4) Parameter λ serves to parameterize the solutions by fixing the consistent resets, for each mode. Our definition of a solution for a hybrid DAE system generalizes the classical definition for (ODE based) hybrid systems. Observe that conditions for existence and/or uniqueness of solutions are delicate, particularly so because we have ruled out Zeno behaviors.

The theory of DAE differentiation index recalled in Section II-A deeply relies on differentiability, so it does not apply as such to (110). In contrast, the notion of difference index for dAE does not require differentiability. To circumvent the lack of differentiability of (110), we move to its nonstandard semantics. As a matter of fact, the nonstandard semantics of a DAE hybrid system is simple and clean defining.

C. Nonstandard hybrid DAE index

Using (91) and (92), the nonstandard semantics of (110) is the following constraint:

H NS ( • x, x, x ) = b i = B i ( • x, x, x ) if b i then 0 = F i ( • x, x, x ) (115)
where

x τ = x τ +∂ -x τ ∂ i.e., x = x • -x ∂ (116)
That is, we substitute, in S, the left limit x -by the previous value • x defined in (92) and the derivative ẋ by its nonstandard version defined in (91). In the following reasoning, the concepts and notations of Appendix B are used.

Theorem 5: Assume that hybrid DAE (110) is solvable. Then, every solution of (110) is the standardisation of some solution of (115).

Proof: See Appendix A4. The proof relies on the material of Appendix B.

We now study the difference index of dAE (115). Accordingly, using (116) we regard F as a function of the tuple ( • x, x, x • ). To simplify the notations, when no confusion can result, we write F for short instead of F ( • x, x, x • ) and similarly for B. With this convention, we have

H •k NS ( • x, x, v, w) = b •k i = B •k i ( • x, x, x ) if b •k i then 0 = F •k i ( • x, x, x ) (117)
where

v = x • and w = (x •2 , . . . , x •k+1 ) (118)
and x is defined by (116). Using (117), the difference array of dAE ( 115) is

H NS k ( • x, x, v, w) = def      H NS ( • x, x, v, w) H • NS ( • x, x, v, w) . . .( • x, x, v, w) H •k NS ( • x, x, v, w)      (119)
Array ( 119) is amenable to Definition 8 and Theorem 3, which yields a guarded Pantelides algorithm for computing the causality and performing index reduction. As announced in the introduction, the causality is dynamically defined, that is, it depends on the sequence of modes traversed by the system at the considered successive instants t, . . . , t •k .

VII. ANALYZING SOME EXAMPLES

In this section, we analyze some examples. Most of them were already investigated in Section IV-D under their dAE variation.

A. Zero-crossing

We reconsider the example (75) of Section IV-D:

if b then x = h(x -) else ẋ = f (x) b = zero-crossing of g(x) (120)
where boolean guard b selects the events of zero-crossings of some smooth function g(x). Our first task is to make the semantics of (120) precise, particularly at and around the zerocrossing. We do this now by defining different versions of the the nonstandard semantics of (120). The different versions for the nonstandard semantics are:

if b then x = h( • x) else x = f (x) b = Q(x, x • ) (121) if b then x = h( • x) else x = f (x) b = Q( • x, x) (122) if b then x • = h(x) else x = f (x) b = Q( • x, x) (123)
where

x = def x • -x ∂ and Q(z, x) = def [g(z)≤0] ∧ [g(x) > 0] (124)
Semantics ( 121) is closest to (120). Semantics (123) introduces micro-delays to ease causality-to this end, the predicate defining guard b has been shifted backward and reset has been shifted forward. Semantics (122) lies in between the former two. Semantics (121), (122), and (123) have the form (76), (77), and (78), respectively. The analysis developed in Section IV-D1 shows that the only suitable semantics for zero-crossing ( 120) is (123). The discussion at the end of Section IV-D1 carries over to semantics (123).

Conclusions regarding this example:

A first generic conclusion is that we must perform mode-dependent causality analyses by not abstracting predicate evaluation. This was developed in Section IV-B.

The second, specific, conclusion is that we should be careful in defining the semantics of events (such as zero-crossings). Not putting appropriate micro-delays in the nonstandard semantics results in semantics with causality circuits that cannot get repaired by using higher order difference arrays, evidencing that those semantics have infinite index.

A generic form of zero-crossing: Equation (75) specifies a system in which the state is reset each time some zero-crossing event occurs. The generic case is when the zero-crossing event triggers the move from one arbitrary mode to another arbitrary mode. A possible continuous-time form for this system is the following, where the two modes 0 and 1 are encoded by the value of the piecewise constant variable ξ (its dynamics is ξ = 0), which takes first the value 0 and then the value 1:

           init ξ = 0 if ξ = 0 then b = zero-crossing of g(x) if b then ξ = 1 else ξ = 0 if ξ = 0 then F 0 (x, ẋ) = 0 if ξ = 0 else if ξ = 1 then F 1 (x, ẋ) = 0 (125)
System (125) can be put in the generic form (110):

             init ξ = 0 if b then ξ = 1 else ξ = 0 if ξ=0 then F 0 (x, ẋ) = 0 b = zero-crossing of g(x) if ξ=1 then F 1 (x, ẋ) = 0 (126)
A similar analysis as for the simpler zero-crossing example shows that the appropriate nonstandard semantics for (126) is the following, where Q(z, x) is as in (124):

           init ξ = 0 if b then ξ • = 1 else ξ = 0 if [ξ=0] then F 0 (x, x )=0 b=Q( • x, x) if [ξ=1] then F 1 (x, x )=0 (127)
Depending on the particular form for the two dynamics F 0 and F 1 , a consistent causality may be found for (127) exactly as for the simple zero-crossing case. The key point in obtaining this is the additional delay applied to the reset of ξ in the equation if b then [ξ • =1], which ensures that the relation between x, b, and ξ is not fixpoint but can rather be properly scheduled for its evaluation.

B. Unilateral constraint

In this section we reconsider the example of unilateral constraint introduced in Section IV-D, albeit in continuous time:

ẋ = f (x, u) 0 ≤ g(x)
Here is a word for word translation to nonstandard semantics:

   x = f (x, u) b = [g(x)≤0] if b then g(x) = 0 (128)
where we recall that x is the nonstandard translation of ẋ given in (124). This formulation has the form studied in Section III-E3.

C. Complementarity condition

Replacing, in example (114), the derivatives ẋ and ẏ by their nonstandard form yields a system of the form studied in Section IV-D3.

D. Circuit breaker

Figure 2 shows a simple circuit breaker.

ON i OFF i ⊥ U -|u i |≥0 |j i |≥0 |j i |<J i 0=u i |j i |≥J i j 2 u 1 j 1 j u 2 V v 1 v 2 w L F 1 R F 2 L 2 Figure 2.
A simple circuit breaker. Top: the circuit. Bottom: the mode automaton for each fuse i = 1, 2. For the ON mode, the current must stay below a threshold J i , while in the OFF mode, the complementarity condition shown holds. We assume that initial conditions are: fuse F 1 is ON, fuse F 2 is OFF with j 2 = 0 (no current), and J 1 is small enough so that fuse F 1 eventually switches to OFF.

Mode-independent equations for this circuit are:

0 = j 1 + j 2 -j 0 = u 1 + v 1 + w -V 0 = u 1 + v 1 -u 2 -v 2 0 = v 1 -Rj 1 0 = Lj -w 0 = L 2 j 2 -v 2 (129)
whereas the equations of the two fuses are mode-dependent:

while in ON i : ⊥ |j i |<J i 0=u i event ON2OFF i : ⊥ |j i |≥J i while in OFF i : ⊥ U -|u i |≥0 |j i |≥0 (130)
We begin with an informal review of this model, by referring to the classical notion of DAE differentiation index, which we apply in each mode separately-we know this does not give the right notion of index for the whole system. Focus on the complementarity condition in OFF i mode. It involves u i , which is not a state, and j i , which is a state. Hence, the two saturated constraints yield dynamics of index 0 and 1 in the two complementary modes of OFF 1 , respectively. As a result, the system moves from index 1 to index 0-or-1 dynamics for the fuse F 2 as a consequence of the switching of F 1 to OFF mode.

Let us detail the nonstandard semantics, by reusing our approach (125)-(127) to handle the transition ON i →OFF i , which is by a zero-crossing:

           init ξ = 0 if b then ξ • = 1 else ξ = 0 if [ξ=0] then E ON1 b=Q( • x, x) if [ξ=1] then E OFF1 (131)
where E ON1 and E OFF1 are the nonstandard dynamics of the circuit in modes ON 1 and OFF 1 , respectively:

E ON1 :                        0 = j 1 + j 2 -j 0 = u 1 + v 1 + w -V 0 = u 1 + v 1 -u 2 -v 2 0 = v 1 -Rj 1 0 = Lj -w 0 = L 2 j 2 -v 2 0 = u 1 0 ≤ U -|u 2 | ⊥ |j 2 | ≥ 0 E OFF1 :                        0 = j 1 + j 2 -j 0 = u 1 + v 1 + w -V 0 = u 1 + v 1 -u 2 -v 2 0 = v 1 -Rj 1 0 = Lj -w 0 = L 2 j 2 -v 2 0 ≤ U -|u 1 | ⊥ |j 1 | ≥ 0 0 ≤ U -|u 2 | ⊥ |j 2 | ≥ 0
Only the last two equations differ (they are shown in red). The other equations are mode-independent. The two modes involve complementarity conditions. Hence, our first task is to identify the atomic subsets of equations, for the two modes ON 1 and OFF 1 . We first consider mode ON 1 , which we rewrite in three blocks as follows:

E ON1 :                                E 1 : 0 = [j 1 ] out + j 2 -j E 2 : 0 = [v 1 ] out -Rj 1 ---------- E 31 : 0 = L 2 j 2 -v 2 E 32 : 0 = u 1 + v 1 + w -V E 33 : 0 = u 1 + v 1 -u 2 -v 2 E 34 : 0 = u 1 E 35 : 0 ≤ U -|u 2 | ⊥ |j 2 | ≥ 0 ---------- E 4 : 0 = Lj -w
The states are j and j 2 . Accordingly, their current value is known and the first block inherits its consistent causality as shown. The second block collects the equation of the inductance, the two Kirchhoff voltage laws, the constraint on u 1 , and the complementarity condition: we make it an atom with state j 2 . We handle this atom as the second case of Section IV-D3, showing that it has outputs v 2 , u 1 , u 2 , w, j • 2 as desired. Accordingly, the consistent causality for this mode is the following, where E 3 = def {E 31 , . . . , E 35 }:

(j 2 , j)→E 1 →j 1 j 1 →E 2 →v 1 (j 2 , v 1 )→E 3 →(v 2 , u 1 , u 2 , w, j • 2 ) w→E 4 →j •
Performing the same for E OFF1 yields:

E OFF1 :                                E 1 : 0 = [j 1 ] out + j 2 -j E 2 : 0 = [v 1 ] out -Rj 1 ---------- E 31 : 0 = L 2 j 2 -v 2 E 32 : 0 = u 1 + v 1 + w -V E 33 : 0 = u 1 + v 1 -u 2 -v 2 E 34 : 0 ≤ U -|u 1 | ⊥ |j 1 | ≥ 0 E 35 : 0 ≤ U -|u 2 | ⊥ |j 2 | ≥ 0 ---------- E 4 : 0 = Lj -w
Focus on E 34 . When 0 < U -|u 1 | holds, this complementarity condition enforces j 1 =0, which, together with E 1 , makes this system overconstrained, hence singular. This mode is therefore non reachable and E 34 boils down to E 34 : U = u 1 and, thus, E OFF1 is handled the same way as E ON1 . The resulting scenario for this circuit is then: in a first phase, fuse F 1 is on and fuse F 2 is OFF and open (i.e., such that j 2 =0); the circuit remains so until j 1 crosses threshold J 1 ; then, fuse F 1 switches to mode OFF and leaks immediately (j 1 >0, u 1 =U ), while the status of fuse F 2 remains open.

Conclusions regarding this example:

• The system has index 1;

• Its causality analysis is mode-dependent.

E. Discussion

The analysis of the above examples shows that our notion of index for DAE hybrid systems is the adequate generalization of the notion of consistent causality for ODE hybrid systems. Our notion of index is particularly useful at handling nonintuitive situations where the index changes at events or cascades thereof. We have shown that it is important not to abstract modes while performing this causality analysis: the right notion of causality must be mode-dependent. Not doing so may result in overly optimistic causality analyses, see the zero-crossing example. In turn, the index itself remains globally defined for the overall system and is a unique integer (not mode-dependent).

VIII. ALGORITHMS

B, X , X, E, are as in Section IV-B. Recall that E denotes equations or atoms. In this section, we provide an encoding of guarded Pantelides graphs as systems of Boolean equations, which we call Pantelides systems. Boolean values F, T and encoded as the integers 0, 1. Then, we consider the following sets of variables with domain {0, 1}:

x(E) δ(E, x) where x ∈ X, E ∈ E (132)
Non-directed and directed branches of the Pantelides graph are encoded as follows:

x

(E) = 1 encodes branch E--x x(E)δ(E, x) = 0 encodes branch E→x (133)
For the first equation, x(E) takes the value

(E)) b = 0 encodes branch if b then E--x. (1 -x
Stating that x is an output of E when x is involved in E and guard b = T is written:

b x(E) δ(E, x) = 0
In particular, due to condition 3) of Definition 8, we have, for every guard b, the equation B defining it, and every other variable x involved in B,

       0 = 1 -x(B) 0 = 1 -b(B) 0 = 1 -δ(B, x) 0 = δ(B, b) (134)
which is equivalent to the single equation In addition to the above encoding of the guarded Pantelides graph using Pantelides systems, additional constraints are for consideration.

Axioms:

The following axioms must hold, in order for a Pantelides system to represent a consistent causality: for every x ∈ X , we have Axiom (136) states that, in any mode (specified by a given guard), any variable cannot be determined by more than one equation. Axiom (135) follows from the fact that a previous variable cannot be an output. 2

Goal: The following is the goal of index reduction. For every x ∈ X :

1 = E∈E,b∈B b x(E) (1 -δ(E, x • )) (137)
It expresses that shifted variables must be entirely, not partially, determined. 2

Observe that Axioms (135) and (136), as well as goal (137) do not involve guards. Hence (135)-(137) must hold for any reachable configuration of the guards.

IX. CONCLUSION

To our knowledge, no proper notion of index existed for hybrid DAE systems (multi-mode). By relying on nonstandard analysis we were able to propose such an extension.

Nonstandard analysis formalizes differential equations as discrete step transition systems with infinitesimal time basis. The nonstandard semantics of a hybrid DAE system yields a (discrete time) difference Algebraic Equation (dAE), for which the notion of difference index can be defined as well-the difference index of a difference Algebraic Equation (dAE) is an easy transposition of the differentiation index, in which forward shift replaces derivation. We proved that the differentiation index of a DAE is equal to the difference index of its nonstandard semantics, structurally. "Structurally" means that we are seeking for a "generic equality", i.e., an equality that is valid if the actual parameters defining the system dynamics remain outside some exceptional set. Thanks to this result, we can propose, as a definition of the index of a hybrid DAE system, the index of its nonstandard semantics (which is a dAE system) and this definition is a conservative extension of both the DAE index and the dAE index.

For both DAE and dAE systems, the structural index is computed by using graph based algorithms of Pantelides type. Due the the need for considering modes and their guards, computing the structural index of a dAE system or a hybrid DAE system requires a new guarded causality analysis, in which causalities are derived in a mode-dependent way. Whereas causality analyses must become mode dependent, the index in itself remains a global notion. We illustrated all of this on some examples.

The algorithms we proposed avoid listing the different modes explicitly. In our guarded causality analysis, branches of the graph are labeled by a predicate characterizing the set of all states in which the considered branch occurs in the Pantelides graph. Doing this is expected to cope with the combinatorial explosion resulting from composing a large number of multi-mode systems. Also, we encode the search for a guarded causality as an integer programming problem, thus making it possible to solve it in a modular way-our algorithms are modular in that they provide results that can be reused in subsequent contexts, not known in advance. This is not the case for the purely graph based algorithms, which require knowing the entire system before analyzing it.

ACKNOWLEDGEMENT: Philippe Chartier is indebted for fruitful discussions. Equation (149) characterizes the consistent x, whereas equation (148) is used to determine v, as a function of x. Now, since A 2 was obtained by pre-and post-multiplying A by permutation matrices and then taking a submatrix, the zero coefficients of A are preserved in A 2 (modulo permutation and extraction) when the nonzero coefficients of A vary over some neighborhood. Hence (A, C) structurally nonsingular implies that A 2 is structurally nonsingular too. Hence, by Corollary 1, there exists a permutation matrix P 2 such that P 2 A 2 has nonzero diagonal. Using (142), (146), and (147), and setting P = I 0 0 P 2 P 1 , where P 2 = P 2 0 0 I P 2 yields the desired permutation matrix P . This proves the only if part and the lemma is proved.

3) Proof of Theorem 1: We first prove the only if part. Assume that S possesses a causality analysis. Chose a total order on Z S that is an extension of order S . This total order induces a renumbering of the equations E i , thus defining a permutation matrix P , of dimension n, the number of equations in S. We move to the if part. Assume that (A, C) is structurally nonsingular. The reader is referred to the proof of Lemma 2, from which the following argument is reused. Using [START_REF]Principles of Object-Oriented Modeling and Simulation with Modelica 3.3[END_REF], equation Av+Cw+x=0 can be given the form (139,140,141). Eqn. (140) uniquely determines v as a partial function of x by using the same pivoting argument as for the proof of the if part of Lemma 1. Since v is fixed by (140), Eqn. (141) characterizes the consistent x. Eliminating w is performed by discarding eqn. (139), since C 1 is structurally onto. By using Corollary 1, we get a causality analysis for S. 4) Proof of Theorem 5: Consider the standard hybrid DAE system (110) and let (x s , b s ), s = (t, k) ∈ S, be a solution for it in the sense of Definition 10. Then, let λ be the parameter representing its initial condition. We show that this solution is the standardisation of some solution of hybrid dAE (115). Two cases can occur, see condition 1) of Definition 10:

Case 1: t ∈ T (λ). Then, by condition 1) of Definition 10, we can assume that s = (t, 0), which we identify with t. Let i ∈ I be the active mode, such that b i,t = 1 holds. We then have F i (x - t , x t , ẋt ) = 0. Pick T τ = t. We can thus assume τ = [t n ] for some sequence t n of reals converging to t. Hence,

∃N ∈ N : n≥N =⇒ t n ∈ T (λ), (151)
thus, b i,(tn,0) = 1. Therefore, setting x τ = [x tn ] and ẋτ = [ ẋtn ], we get b i (x - τ , x τ , ẋτ ) = [b i,tn ] = 1. Similarly, since F i is smooth and the solution of F i (x -, x, ẋ) = 0 is assumed to be infinitely differentiable, we get F i (x - tn , x tn , ẋtn ) = 0 for n≥N where N is as in (151), whence F i (x - τ , x τ , x τ ) ≈ 0 follows as well. For this case, we thus proved the existence of a solution (x τ , b τ ) for (115) such that (x t , b t ) = st(x τ , b τ ).

Case 2: t ∈ T (λ), meaning that one or more successive mode changes occur at t, so that the super-dense instants for consideration are (t, 0), (t, 1), . . . , (t, n t (λ)), on which the active mode changes at every instant. Suppose that b i,(t,0) = 1, whence x (t,0) is a consistent value for F i = 0. By condition 2 of Definition 10, we also have b i,(tn,0) = 1 for t n any increasing sequence converging to t and n large enough. Consider T τ = [t n ], we have st(τ ) = t. Setting x τ = def [x (tn,0) ] yields a consistent value for F i (x -, x, x ) = 0 and complementing it with b i (x - τ , x τ , x τ ) = 1 extends the solution of (115) at the considered τ . We further extend this solution for the subsequent non-standard instants τ + ∂, . . . , τ + n t ∂ as follows. First, observe that τ + ∂ ≈ • • • ≈ τ + n t ∂ ≈ t. Then, we simply extend the solution of (115) by setting b i,τ +k∂ = b i,(t,k) and x τ +k∂ = x (t,k) , where (x (t,k) , b i,(t,k) ) is the given (standard) solution of (110). This proves the theorem.

  assert P (b 1 , . . . , b k ) where b 1 , . . . , b k ∈ B (61) where P (b 1 , . . . , b k ) is a boolean expression over the listed guards; this assertion states that predicate P (b 1 , . . . , b k ) holds true in the dAE system where this assertion occurs.Guarded equations: if b then E (62) where (b, E) ∈ B×E; this guarded equation specifies that equation E is in force whenever guard b is true. Convention and associated notation: We use the generic notation E (63) to denote a possibly guarded equation, i.e., a guarded equation (61) or simply an element of E. For E a guarded equation of the form "if b then E", we denote by b(E) = def b its guard. For a non guarded equation E and considering Axiom 1 below, we set by convention b(E) = T.

Comment 2 :

 2 Entailment rule (74) forbids the consideration of systems in which guards are themselves solutions of fixpoint equations, such as in the following example:b = [ x > 0 ] if b then [ 0 = f (x) ]which means informally: if b then f (x)=0 where b=[x>0]. Such fixpoint equations cannot be solved by graph based abstractions. Numerical methods must be used-a natural candidate being relaxation, by which a sequence b 0 , x 1 , b 1 , x 2 , . . . of values is computed by using the two equations in a loop until (hopefully) fixpoint occurs. Graph based algorithms à la Pantelides do require condition (67) of Axiom 1.

Figure 1 .

 1 Figure 1. Trajectory of a hybrid system.

[ 1 -

 1 x(B)] + [1-b(B)] + [1-δ(B, x)] + δ(B, b) = 0

∀E ∈ E : 1 =) 1 ≥

 11 x(E)δ(E, • x) (135E∈E,b∈B b x(E) (1 -δ(E, x)) (136)

Write P 1 2

 12 Av + Cw + x = 0 rewrites 0 = A 1 v + C 1 w + P 11 x (144) 0 = A 2 v + P 12 x(145) Eliminating w in equation Av + Cw + x = 0 amounts to discarding (144) and keeping only (145). Now, let the nonzero coefficients of A and C vary over some neighborhood. The same reasoning used in the proof of the if part of Lemma 1 shows that a same matrix P 1 can work almost everywhere when the nonzero coefficients of A and C vary over some neighborhood. (146) Thus we can now eliminate w and focus on equation (145) in which P 11 and P 12 are fixed. Apply to (145) the same pivoting argument we used to eliminate w in equation Av+Cw+x = 0 yields a permutation matrix P 2 such that P 2 v + P 21 P 12 x (148) 0 = P 22 P 12 x (149)

+ P 1 x by 5 )Setting C 1 =

 51 Using conditions 5) and 6) of Definition 4, we can rewrite equation Av+Cw+x=0 as follows: by 1) and 6) :0 = A 1 v + [ C 11 C 12 ] w w : 0 = A 2 v + P 2 x by 5) : 0 = P 3 x (150)where C 11 is structurally nonsingular, w collects the variables of W S ∩ Z, possibly reordered, and permutation matrix P decomposes as C 11 C 12 shows that (150) is a particular form of (21), hence the pair (A, C) is structurally nonsingular by Lemma 2.

  Definition 4: To each system S={E i | i∈I}, we associate its Pantelides graph P S , which is a nondirected bipartite graph defined as follows, using notation (22):

			• Its set of vertices is ∈ S S;
			• For every z ∈ ∈ S :
			E ∈ S and z ∈ E
			z--E ∈ P S
			.
	In the sequel,
	Prop 1 Prop 2	stands for Prop 1 entails Prop 2	(22)

  44)This should not come as a surprise since we know that (43) is a constraint. This leads to considering the set of equations {E 21 , E 22 } as atomic, meaning that it must be evaluated at once. To reflect this on the causality analysis, simply perform the following transformations on (44): 1) Reinforce (44) to x→E 21 →b→E 22 →x; 2) Turn the so obtained circuit into the following consistent causality involving the atom {E 21 , E 22 }(x, b):

  2) P C covers P C : for every z ∈ Z, z--E ∈ P C either z→E ∈ P C or E→z ∈ P C C has the single assignment property, meaning that, for every z, there exists at most one equation E such that E→z ∈ P C ; 5) P C is circuitfree; its transitive closure C is thus a partial order on Z C, called the causality order induced by P C ; 6) Denoting by min( C ) the set of the minimal elements of C , we have min

	3) Causality specifications are respected:
	z--E ∈ P C and z ∈ out E
	z→E ∈ P C
	4) P

  ) P C has the guarded single assignment property; 6) P C is guarded circuitfree, hence we can consider the transitive closure C of the graph part of P C , which we call the guarded causality order induced by P C ; 7) Denoting by min( C ) the set of the vertices that are minimal elements of C , then min( C ) ⊆ X C C holds; 8) For every pair

RR n°8630

We may also consider F (t, x, ẋ) = 0, but dependence on time t can always be removed by making t an additional variable obeying ṫ = 1.

The term used in[START_REF] Campbell | The index of general nonlinear DAEs[END_REF] is geometrically solvable.

Observe that, as a consequence, µ cannot be sigma-additive (in contrast to probability measures or Radon measures) in that it is not true that µ( n An) = n µ(An) holds for an infinite denumerable sequence An of pairwise disjoint subsets of N.RR n°8630

APPENDIX

A. Collecting proofs 1) Proof of Lemma 1: We consider the linear equation Av=y, where y∈R n and v is the R m -valued unknown with coordinates v 1 , . . . , v m . We assume m≥n otherwise the lemma is trivial.

We first prove the if part. We want to prove that, with the condition of the lemma, for every y, there exists a v such that Av=y holds, almost everywhere when the non zero coefficients of A vary over some neighborhood. With a suitable renumbering of the coordinates of v, we can assume that Q is the identity matrix. Let P be a permutation matrix such that A = def P A = B 1 B 2 where B 1 has a nonzero diagonal, and let σ be the permutation of {1, . . . , n} associated with P . Since a nn =0 we can express v n in terms of v 1 . . . v n-1 , y σ(n) and then remove the last equation. This yields a reduced equation A v =y , where A is (n -1)×(m -1), y ∈R n-1 , and v is the R m-1 -valued unknown collecting v 1 . . . v n-1 , v n+1 , . . . , v m . Matrix A has diagonal entries equal to a kk =a kk -a kn /a nn . Since a kk and a nn are non zero, the set of entries a ij of matrix A causing a kk =0 for some k is exceptional in R n×m since it requires the condition a kn = a kk a nn to hold. The corresponding set of entries a ij of matrix A = P -1 A is exceptional as well, we denote by Ξ n this exceptional subset of R n×m . Thus we assume that the entries a ij of matrix A do not belong to exceptional set Ξ n . Since the diagonal entries of matrix A are all non zero, we can express v n-1 in terms of v 1 . . . v n-2 , y σ(n-1) and we are left with a further reduced equation A v = y for which the same argument as before applies, namely: if the entries a ij of matrix A do not belong to exceptional set Ξ n-1 , then the diagonal entries of A are all non zero. And so on, by reducing the number of equations and unknowns. We have thus proved the existence of exceptional sets

uniquely determined as a function of the coordinates of y, by pivoting. This concludes the if part, since set

For the only-if part, assume that A is structurally onto, meaning that A is almost everywhere onto when its non-zero entries vary over some neighborhood. Call V (A) the resulting set of neighbor matrices of A and let Ξ n be the exceptional set for the entries ãij of à outside which à is onto, for à ranging over V (A). For each onto Ã, there exists a permutation matrix

where B 1 is square invertible. Neighborhood V (A) decomposes as the union V (A) = σ V (A, σ), where V (A, σ) collects the matrices à for which decomposition (138) holds with the matrix Q σ representing permutation σ. Since the number of permutations is finite, at least one of the V (A, σ) has non-empty interior. Hence, without loss of generality we can assume that decomposition (138) holds with the same matrix Q. The map V (A) à → det(B 1 ) takes the value 0 on an exceptional set only. Now, we have det(B 1 ) = n i=1 bi1 det(B i1 ), where B ij is the submatrix of B 1 obtained by erasing row i and column j. Since B 1 is structurally nonsingular, there must be some i ∈ {1, . . . , n} such that b i1 =0 and B i1 is structurally nonsingular-otherwise we would have det(B 1 ) = 0 for any Ã∈W (A), where W (A) is some neighborhood contained in V (A). Let P 1 be the permutation matrix exchanging rows 1 and i, so we replace

and it is enough to prove the lemma for B 1 1 . Now, since we have b 1 11 = b 1 i1 =0 and B 1 11 is structurally nonsingular, we apply the same reasoning to B 1 11 , which yields a second permutation matrix P 2 and so on. The wanted left permutation matrix is P = P m . . .P 2 P 1 .

2) Proof of Lemma 2: Since we reuse the same techniques as for the proof of Lemma 1, some details will be omitted.

We begin with the if part. Using [START_REF]Principles of Object-Oriented Modeling and Simulation with Modelica 3.3[END_REF], equation Av+Cw+x=0 becomes

where For the only if part, consider equation Av + Cw + x = 0 and use pivoting to eliminate w. To this end, pick an equation, say the i th one, and a component of w, say the j th one, such that w j has a nonzero coefficient in the i th equation. Let P 1 be the permutation matrix exchanging the equations 1 and i, and let Q 1 be the permutation matrix exchanging the rows 1 and j in w. Then, discard the first equation and the first coordinate of the permuted w and repeat the same argument on the remaining equations and components of the permuted w. Performing pivoting as above yields a second pair (P 2 , Q 2 ), leaving invariant the first equation and the first component of w. The procedure repeats until step r≤q, where no more component of w remains that has a nonzero coefficient in some equation. Setting

yields

where r×q-matrix C 1 has rank r and allows using the r first permuted equations to eliminate w entirely. Accordingly,

RR n°8630

B. A primer on non-standard analysis

The background material of this section is used in proofs, so the reader can skip it for a first reading. The text is borrowed verbatim from [START_REF] Benveniste | Nonstandard semantics of hybrid systems modelers[END_REF].

1) Intuitive introduction: We begin with an intuitive introduction to the construction of the non-standard reals. The goal is to augment R ∪ {±∞} by adding, to each x in the set, a set of elements that are "infinitesimally close" to it. We will call the resulting set R. Another requirement is that all operations and relations defined on R should extend to R.

A first idea is to represent such additional numbers as convergent sequences of reals. For example, elements infinitesimally close to the real number zero are the sequences u n = 1/n, v n = 1/ √ n and w n = 1/n 2 . Observe that the above three sequences can be ordered: v n > u n > w n > 0 where 0 denotes the constant zero sequence. Of course, infinitely large elements (close to +∞) can also be considered, e.g., sequences x u = n, y n = √ n, and z n = n 2 . Unfortunately, this way of defining R does not yield a total order since two sequences converging to zero cannot always be compared: if u n and u n are two such sequences, the three sets {n | u n > u n }, {n | u n = u n }, and {n | u n < u n } may even all be infinitely large. The beautiful idea of Lindstrøm is to enforce that exactly one of the above sets is important and the other two can be neglected. This is achieved by fixing once and for all a finitely additive positive measure µ over the set N of integers with the following properties: 7 1) µ : 2 N → {0, 1};

2) µ(X) = 0 whenever X is finite; 3) µ(N) = 1. Now, once µ is fixed, one can compare any two sequences: for the above case, exactly one of the three sets must have µ-measure 1 and the others must have µ-measure 0. Thus, say that u > u , u = u , or u < u , if µ({n

respectively. Indeed, the same trick works for many other relations and operations on non-standard real numbers, as we shall see. We now proceed with a more formal presentation.

2) Non-standard domains: For I an arbitrary set, a filter F over I is a family of subsets of I such that:

1) the empty set does not belong to F, 2) P, Q ∈ F implies P ∩ Q ∈ F, and 3) P ∈ F and P ⊂ Q ⊆ I implies Q ∈ F. Consequently, F cannot contain both a set P and its complement P c . A filter that contains one of the two for any subset P ⊆ I is called an ultra-filter. At this point we recall Zorn's lemma, known to be equivalent to the axiom of choice:

Lemma 4 (Zorn's lemma): Any partially ordered set (X, ≤) such that any chain in X possesses an upper bound has a maximal element. A filter F over I is an ultra-filter if and only if it is maximal with respect to set inclusion. By Zorn's lemma, any filter F over I can be extended to an ultra-filter over I. Now, if I is 7 The existence of such a measure is non trivial and is explained later. infinite, the family of sets F = {P ⊆ I | P c is finite} is a free filter, meaning it contains no finite set. It can thus be extended to a free ultra-filter over I:

Lemma 5: Any infinite set has a free ultra-filter. Every free ultra-filter F over I uniquely defines, by setting µ(P ) = 1 if P ∈ F and otherwise 0, a finitely additive measure 8 µ : 2 I → {0, 1}, which satisfies µ(I) = 1 and, if P is finite, then µ(P ) = 0. Now, fix an infinite set I and a finitely additive measure µ over I as above. Let X be a set and consider the Cartesian product

Relation ≈ is an equivalence relation whose equivalence classes are denoted by [x i ] and we define

X is naturally embedded into X by mapping every x ∈ X to the constant tuple such that x i = x for every i ∈ I.

Any algebraic structure over X (group, ring, field) carries over to X by almost point-wise extension. In particular, if

The existence of an inverse for any non-zero element of a ring is indeed stated by the formula: ∀x (x = 0 ∨ ∃y (xy = 1)). More generally:

Lemma 6 (Transfer Principle): Every first order formula is true over X if and only if it is true over X.

3) Non-standard reals and integers: The above general construction can simply be applied to X = R and I = N. The result is denoted R; it is a field according to the transfer principle. By the same principle, R is totally ordered by

To prove this, let It is also of interest to apply the general construction (152) to X = I = N, which results in the set N of non-standard natural numbers. The non-standard set N differs from N by the addition of infinite natural numbers, which are equivalence classes of sequences of integers whose essential limit is +∞. 4) Integrals and ODE: Any sequence (g n ) of functions g n : R → R point-wise defines a function [g n ] : R → R by setting

A function R → R so obtained is called internal. Properties of and operations on ordinary functions extend point-wise to internal functions of R → R. The non-standard version of g : R → R is the internal function g = [g, g, g, . . . ]. The same notions apply to sets. An internal set A

Now, consider an infinite number N ∈ N and the set

By definition, if

The sum of g over A can be defined:

If t is as above, and f : R → R is a standard function, we obtain

Under the same assumptions, for any t ∈ [0, 1],

Now, consider the following ODE:

Assume (159) possesses a solution [0, 1] t → x(t) such that the function t → f (x(t), t) is continuous. Rewriting (159) in its equivalent integral form x(t) = x 0 + t 0 f (x(u), u)du and using (158) yields

The substitution in (160) of ∂ = 1/N , which is positive and infinitesimal, yields T = {t n = n∂ | n = 0, . . . , N }. The expression in parentheses on the right hand side of (160) is the piecewise-constant right-continuous function x(t), t ∈ [0, 1] such that, for n = 1, . . . , N :

x(t n ) = x(t n-1 ) + ∂ × f ( x(t n-1 ), t n-1 )

x(t 0 ) = x 0 (161) By (160), the solutions x, of ODE (159), and x, as computed by algorithm (161), are related by x = st( x). Formula (161) can be seen as a non-standard semantics for ODE (159); one which depends on the choice of infinitesimal step parameter ∂. Property (160), though, expresses the idea that all these nonstandard semantics are equivalent from the standard viewpoint regardless of the choice made for ∂. This fact is referred to as the standardization principle.

RR n°8630