Marius Bozga
email: bozga@imag.fr

Radu Iosif
email: iosif@imag.fr

Filip Konecny
email: ikonecny@fit.vutbr.cz

Filip Konečný

Fast Acceleration of Ultimately Periodic Relations

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

The verification of safety properties of infinite-state systems (such as device drivers, communication protocols, control software, etc.) requires the computation of the set of reachable states, starting with an initial state from a given (possibly infinite) set. There are currently two ways of doing this: (i) compute a finite representation of an overapproximation of the set of reachable states, by applying a widening operator at each step, or (ii) attempt to compute precisely the transitive closure of the transition relation; the set of reachable states is the image of the set of initial states via the transitive closure. The first approach is guaranteed to terminate, but the abstraction usually introduces imprecision that may blur the verification result. On the other hand, the second approach, although precise, is not guaranteed to terminate -the problem of verifying safety properties being, in general, undecidable.

In practice, one usually tries to combine the two approaches and benefit from the advantages of both. To this end, it is important to know for which classes of transition relations it is possible to compute the transitive closure precisely and fast -the relations falling outside these classes being dealt with using suitable abstractions. To the best of our knowledge, the three main classes of integer relations for which transitive closures can be computed precisely in finite time are: [START_REF] Annichini | Symbolic techniques for parametric reasoning about counter and clock systems[END_REF] difference bounds constraints [START_REF] Comon | Multiple Counters Automata, Safety Analysis and Presburger Arithmetic[END_REF][START_REF] Bozga | Flat parametric counter automata[END_REF], (2) octagons [START_REF] Miné | The octagon abstract domain[END_REF][START_REF] Bozga | Iterating octagons[END_REF], and (3) finite monoid affine transformations [START_REF] Boigelot | Symbolic Methods for Exploring Infinite State Spaces[END_REF][START_REF] Finkel | How to compose presburger-accelerations: Applications to broadcast protocols[END_REF]. For these three classes, the transitive closures can be moreover defined in Presburger arithmetic.

The contributions of this paper are three-fold. First, we study the three non-trivial classes of relations mentioned in the previous and show that they are ultimately periodic, i.e. that each relation R in these classes can be mapped into an integer matrix M R such that the sequence {M R k } ∞ k=0 is periodic. The proof of the fact that a sequence of matrices is ultimately periodic relies on a result from tropical semiring theory [START_REF] Schutter | On the ultimate behavior of the sequence of consecutive powers of a matrix in the max-plus algebra[END_REF].

This provides shorter proofs to the fact that the transitive closures for these classes can be effectivelly computed, and that they are Presburger definable.

Second, from the practical side, the algorithm introduced in this paper computes the transitive closure of difference bounds and octagonal relations up to four orders of magnitude faster than our original methods from [START_REF] Bozga | Flat parametric counter automata[END_REF][START_REF] Bozga | Iterating octagons[END_REF], and also scales much better in the number of variables. The experimental comparison with the FAST tool [START_REF] Bardin | Fast extended release[END_REF] for difference bounds relations shows that large relations (> 50 variables), causing FAST to run out of memory, can now be handled by our implementation in less than 8 seconds, on average. We currently do not have a full implementation of the finite monoid affine transformation class, which is needed in order to compare our method with tools like FAST [START_REF] Bardin | Fast extended release[END_REF], LASH [START_REF] Wolper | Verifying systems with infinite but regular state spaces[END_REF], or TReX [START_REF] Annichini | Trex: A tool for reachability analysis of complex systems[END_REF], for this class of relations.

Related Work

Early attempts to apply Model Checking techniques to the verification of infinite-state systems consider the problem of accelerating transition relations by successive under-approximations, without any guarantee of termination. For systems with integer variables, the acceleration of affine relations has been considered primarily in the works of Annichini et. al [START_REF] Annichini | Symbolic techniques for parametric reasoning about counter and clock systems[END_REF], Boigelot [START_REF] Boigelot | Symbolic Methods for Exploring Infinite State Spaces[END_REF], and Finkel and Leroux [START_REF] Finkel | How to compose presburger-accelerations: Applications to broadcast protocols[END_REF]. Finite monoid affine relations have been first studied by Boigelot [START_REF] Boigelot | Symbolic Methods for Exploring Infinite State Spaces[END_REF], who shows that the finite monoid property is decidable, and that the transitive closure is Presburger definable in this case. On what concerns non-deterministic transition relations, difference bounds constraints appear in the context of timed automata verification. The transitive closure of a difference bounds constraint is shown to be Presburger definable first by Comon and Jurski [START_REF] Comon | Multiple Counters Automata, Safety Analysis and Presburger Arithmetic[END_REF]. Their proof was subsequently simplified and extended to parametric difference bounds constraints in [START_REF] Bozga | Flat parametric counter automata[END_REF]. We also showed that octagonal relations can be accelerated precisely, and that the transitive closure is also Presburger definable [START_REF] Bozga | Iterating octagons[END_REF]. The proofs of ultimate periodicity from this paper are based on some of our previous results [START_REF] Bozga | Flat parametric counter automata[END_REF][START_REF] Bozga | Iterating octagons[END_REF]. For difference bounds constraints, the proof from [START_REF] Bozga | Flat parametric counter automata[END_REF] was simplified using a result from tropical semiring theory [START_REF] Schutter | On the ultimate behavior of the sequence of consecutive powers of a matrix in the max-plus algebra[END_REF]. Roadmap The paper is organized as follows: Section 2 gives the definition of ultimately periodic relations, Section 3 describes the algorithm for computing transitive closures of ultimately periodic relations, in general, Section 4 describes three instances of the algorithm, Section 5 presents the experimental results, and Section 6 concludes. Missing proofs are deferred to [START_REF] Bozga | Fast Acceleration of Ultimately Periodic Relations[END_REF] due to reasons of space.

Preliminaries

We denote by Z, N and N + the sets of integers, positive (including zero) and strictly positive integers, respectivelly. The first order additive theory of integers is known as Presburger Arithmetic. The tropical semiring is defined as T = (Z ∞ , min, +, ∞, 0) [START_REF] Schutter | On the ultimate behavior of the sequence of consecutive powers of a matrix in the max-plus algebra[END_REF], where Z ∞ = Z ∪ {∞}, with the extended arithmetic operations x + ∞ = ∞, min(x, ∞) = x, for all x ∈ Z, where min(x, y) denotes the minimum between the values x and y. For two square matrices A, B ∈ S m×m , we define

(A + B) ij = A ij + B ij and (A × B) ij = min m k=1 (a ik + b kj ), for all 1 ≤ i, j ≤ m.
Let I ∈ T m×m be the identity matrix, i.e. I ii = 0 and

I ij = ∞, for all 1 ≤ i, j ≤ m, i = j. Definition 1. [13] An infinite sequence {s k } ∞ k=0 ∈ T is called ultimately periodic if: ∃K ∃c > 0 ∃λ 0 , λ 1 , . . . , λ c-1 ∈ T . s (k+1)c+i = λ i + s kc+i
for all k ≥ K and i = 0, 1, . . . , c -1. The smallest c and λ 0 , λ 1 , . . . , λ c-1 for which the above holds are called the period and rates of {s k } ∞ k=0 , respectivelly.

Example 1. The sequence

σ k = {3k+1 | k = 2l, l ≥ 2} ∪ {5k+3 | k = 2l+1, l ≥ 2}
is ultimately periodic, with K = 4, period c = 2 and rates λ 0 = 6, λ 1 = 10.

⊓ ⊔

A sequence of matrices {A k } ∞ k=0 ∈ T m×m is said to be ultimately periodic if, for all

1 ≤ i, j ≤ m, the sequence {(A k ) ij } ∞ k=0 is ultimately periodic. A matrix A ∈ T m×m is called ultimately periodic if the sequence {A k } ∞
k=1 is ultimately periodic, where A 0 = I and A k = A × A k-1 , for any k > 0. It is known that, every matrix is ultimately periodic in the tropical semiring [START_REF] Schutter | On the ultimate behavior of the sequence of consecutive powers of a matrix in the max-plus algebra[END_REF].

We have the following characterization of ultimately periodic sequences of matrices:

Lemma 1. A sequence of matrices {A k } ∞ k=1 ∈ T m×m is

ultimately periodic if and only if:

∃K ∃c > 0 ∃Λ 0 , Λ 1 , . . . , Λ c-1 ∈ T m×m . A (k+1)c+i = Λ i + A kc+i for all k ≥ K and i = 0, 1, . . . , c -1.
If A ∈ T m×m is a square matrix and n ∈ T, we define the matrix

(n•A) ij = n•A ij , for all 1 ≤ i, j ≤ m. If k is a parameter (typically interpreted over T), then T[k]
denotes the set of all terms where k may occur, built from the constants and operators of T. For instance

, if A, B ∈ T m×m , then k • A + B ∈ T[k] m×m denotes the matrix of terms (k • A + B) ij = k • A ij + B ij , for all 1 ≤ i, j ≤ m.

Ultimately Periodic Relations

Let x = {x 1 , x 2 , . . . , x N } be a set of variables, N > 0, and let x ′ = {x ′ 1 , x ′ 2 , . . . , x ′ N }. A relation is an arithmetic formula R(x, x ′ ) with free variables x ∪ x ′ . We say that two relations R and R ′ are equivalent, denoted R ⇔ R ′ if under all valuations of x and x ′ , R is true if and only if R ′ is true. A relation is called consistent if and only if there exist valuations of x and x ′ under which it holds. We denote a consistent relation R by writing R false, and an inconsistent relation by writing R ⇔ false. The composition of two relations is defined as R • R ′ ≡ ∃y . R(x, y) ∧ R ′ (y, x ′ ). Let I be the identity relation x∈x x ′ = x. We define R 0 ≡ I and

R n ≡ R n-1 • R, for any n > 0. With these notations, R * ≡ ∞ i=0 R i denotes the transitive closure of R. A relation R is called ω-consistent if R n is consistent for all n > 0.
For the rest of this section, let R be a class of relations 3 .

Definition 2. A relation R(x, x ′ ) ∈ R is called ultimately periodic if

and only if either:

1. there exists i 0 ≥ 0 such that R i0 is inconsistent, or 2. for all i ≥ 0, R i is consistent, and there exists two functions:

σ : R → T m×m ⊥ mapping each consistent relation in R into a m × m matrix of T, for some m > 0, and each inconsistent relation into ⊥.

ρ :

T m×m → R mapping each m × m matrix of T into a relation in R, such that ρ(σ(R)) ⇔ R, for each consistent relation R ∈ R such that the infinite sequence of matrices {σ(R i )} ∞ i=0 ∈ T m×m is ultimately pe- riodic.
Notice that the first condition of the definition implies that σ(R i ) = ⊥, for all i ≥ i 0 . If each relation R ∈ R is ultimately periodic, then R is called ultimately periodic as well.

The following lemma gives an alternative characterization of ω-consistent ultimately periodic relations.

Lemma 2. An ω-consistent relation R is ultimately periodic if and only if there exist K ≥ 0, b ≥ 0, c > 0 and Λ 0 , Λ 1 , . . . , Λ c-1 ∈ T m×m such that the following hold:

1. σ(R (n+1)c+i ) = Λ i + σ(R nc+i ), for all n ≥ K. 2. R nc+b+i ⇔ ρ(n • Λ i + σ(R b+i )), for all n ≥ 0.
for all i = 0, 1, . . . , c -1, where σ and ρ are the functions from Def. 2.

Proof. By Lemma 1, if R is ω-consistent, then it is ultimately periodic if and only if

∃K ∃c > 0 ∃Λ 0 , Λ 1 , . . . , Λ c-1 ∈ T N ×N . σ(R (k+1)c+i ) = Λ i + σ(R kc+i )
for all k ≥ K and i = 0, 1, . . . , c -1. By induction on k ≥ K, one shows first that

R kc+i ⇔ ρ(Λ i k-K + σ(R Kc+i )), ∀k ≥ K Let b = Kc. By replacing k -K with k, we obtain R kc+b+i ⇔ ρ(Λ i k + σ(R b+i ))
, ∀k ≥ 0 ⊓ ⊔ For practical reasons related to the representation of R * , we are interested in finding the symbolic expression R k , where k is a parameter (because R * ≡ ∃k . R k ). Notice that the second point of lemma 2 can be used to compute the expression R k symbolically (as a formula over x, x ′ and k), assuming that we are given a function, call it π : T[k] m×m → R(k), where R(k) is the class of all parametric relations over x, x ′ and k. Intuitivelly, π is the parametric counterpart of the ρ function from Def. 2, mapping a matrix of terms over k into a parametric relation R(x, x ′ , k). Concrete definitions of π will be given in Section 4.

Computing Transitive Closures of Ultimately Periodic Relations

In this section we give a generic algorithm that computes the transitive closure of a given ultimately periodic relation. The algorithm needs to be instantiated for a specific class R of ultimately periodic relations by providing the mappings σ, ρ (Def. 2) and π (the parametric counterpart of ρ) as discussed in the previous. Next, in Section 4, we show how this algorithm can be used for accelerating three classes of relations: difference bounds, octagons, and finite monoid affine transformations.

Fig. 1 shows the generic framework for computing transitive closures. The input to the algorithm is a relation R, and the mappings σ : R → T m×m , ρ : T m×m → R, and π : T[k] m×m → R(k). The algorithm is guaranteed to terminate if R is ultimately periodic, as it will be explained in the following.

The main idea of the algorithm is to discover the prefix b and period c of the sequence {σ(R i )} ∞ i=0 -cf. the second point of lemma 2. If R is ultimately periodic, such values are guaranteed to exist. The dove-tailing enumeration on lines 1 and 2 is guaranteed to yield the smallest value c for which the sequence is shown to be periodic 4 .

Second, the algorithm attempts to compute the first rate of the sequence (line 6), by comparing the matrices σ(R b ), σ(R c+b ) and σ(R 2c+b ). If the difference Λ between σ(R c+b ) and σ(R b ) equals the difference between σ(R 2c+b ) and σ(R c+b ), then Λ is a valid candidate for the first rate of the progression (see lemma 2). Notice that the consistency check on line 4 is needed to ensure that we apply σ to consistent relationsotherwise, the relation is not ω-consistent, and the algorithm returns directly the transitive closure, i.e. the finite disjunction

kc+b-1 i=0 R i , 0 ≤ k ≤ 2 (line 4).
Once a candidate Λ for the initial rate was found, the test Q 1 on line 7 is used to check that R is ultimately periodic and ω-consistent. Notice that the characterization of ultimately periodic relations from lemma 2 cannot be applied here, since R n is not known in general, for arbitrary n ≥ 0. The condition used here is local, i.e. it needs only the relation R b , for a typically small constant b ≥ 0. The next lemma establishes the correctness of the criterion used by Q 1 :

Lemma 3. An ω-consistent relation R is ultimately periodic if and only if ∃b ∃c > 0 ∃Λ 0 , Λ 1 , . . . , Λ c-1 ∈ T m×m . ρ(n•Λ i +σ(R b+i ))•R c ⇔ ρ((n+1)•Λ i +σ(R b+i ))
for all n ≥ 0 and i = 0, 1, . . . , c -1, where σ and ρ are the functions from Def. 2. Moreover, Λ 0 , Λ 1 , . . . , Λ c-1 satisfy the equivalences of Lemma 2.

Proof. "⇒" If R is ω-consistent and ultimately periodic, by Lemma 2, there exist b ≥ 0, c > 0 and Λ 0 , Λ 1 , . . . , Λ c-1 ∈ T m×m such that R kc+b+i ⇔ ρ(Λ i k + σ(R b+i ))
for all k ≥ 0 and i = 0, 1, . . . , c -1. We have:

R (k+1)c+b+i ⇔ R kc+b+i • R c ρ(Λ i k+1 + σ(R b+i )) ⇔ ρ(Λ i k + σ(R b+i )) • R c
"⇐" We prove the equivalent condition of Lemma 2 by induction on k ≥ 0. The base case k = 0 is immediate. The induction step is as follows:

R (k+1)c+b+i ⇔ R kc+b+i • R c ⇔ ρ(Λ i k + σ(R b+i )) • R c , by the induction hypothesis ⇔ ρ(Λ i k+1 + σ(R b+i ))
⊓ ⊔ The universal query Q 1 on line 7 is in general handled by procedures that are specific to the class of relations R we work with. Notice furthermore that Q 1 can be handled symbolically by checking the validity of the first order formula:

∀k . π(k • Λ + σ(R b )) • R c ⇔ π((k + 1) • Λ + σ(R b ))
, where π is the parametric counterpart of ρ. Next, in Section 4, we detail two ways in which this test can be performed efficiently (for difference bounds and octagonal relations), without resorting to external proof engines, such as SMT or Presburger solvers. foreach

k := 0, 1, 2 do 4. if R kc+b ⇔ false then return R * ≡ kc+b-1 i=0 R i 5. endfor 6. if exists Λ ∈ T m×m : σ(R c+b ) = Λ + σ(R b ) and σ(R 2c+b ) = Λ + σ(R c+b ) then 7. if forall n ≥ 0 : ρ(n • Λ + σ(R b )) • R c ⇔ ρ((n + 1) • Λ + σ(R b )) false (Q1) then 8. return R * ≡ b-1 i=0 R i ∨ ∃k ≥ 0 . c-1 i=0 π(k • Λ + σ(R b )) • R i 9. else if exists n ≥ 0 : ρ(n • Λ + σ(R b )) • R c ⇔ false (Q2) then 10. let n0 = min{n | ρ(n • Λ + σ(R b )) • R c ⇔ false} 11. if forall n ∈ [0, n0 -1] : ρ(n • Λ + σ(R b )) • R c ⇔ ρ((n + 1) • Λ + σ(R b )) then 12. return R * ≡ b-1 i=0 R i ∨ n 0 -1 n=0 c-1 i=0 ρ(n • Λ + σ(R b )) • R i 13. endif 14.
endif 15. endfor 16. endfor Fig. 1: Transitive Closure Algorithm If the universal query on line 7 holds, the rate Λ can be used now to express the transitive closure (line 8) as a finite disjunction over the prefix ( b-1 i=0 R i ) followed by a formula defining an arbitrary number of iterations (∃k .

c-1 i=0 π(k • Λ + σ(R b )) • R i ).
Note that the formula on line 8 defines indeed the transitive closure of R, as a consequence of lemma 2. Moreover, this is a formula of Presburger arithmetic, provided that the classes of relations R and R(k) are Presburger definable.

Otherwise, if Q 1 does not hold, there are two possibilities: either (i) Λ is actually not the first rate of the sequence {σ(R i )} ∞ i=0 for given b ≥ 0 and c > 0, or (ii) the relation is not ω-consistent. In the first case, we need to reiterate with another prefix-period pair, which will give us another candidate Λ.

In the second case, R m becomes inconsistent, for some m > 0 -in this case the computation of its transitive closure is possible, in principle, by taking the disjunction of all powers of R up to m. However, in practice this may take a long time, if m is large. In order to speed up the computation, we check whether:

-ρ(n • Λ + σ(R b )) • R c is inconsistent (line 9)
; the existential query Q 2 (and namely finding the smallest value for which it holds) is dealt with in Section 4, specifically for the classes of difference bounds and octagonal relations. -R is periodic with first rate Λ between 0 and n 0 -1 (line 11), where n 0 is the smallest n satisfying the first point (line 10).

If both conditions above hold, then m = (n 0 + 1)c + b is the smallest value for which R m becomes inconsistent, and moreover, R is periodic with rate Λ between 0 and m. If this is the case, we compute the transitive closure using the period Λ and return (line 12). The following theorem can be proved along the lines of the discussion above:

Theorem 1. If R is an ultimately periodic relation, the algorithm in Fig. 1 eventually terminates and returns the transitive closure of R.

Some Ultimately Periodic Classes of Arithmetic Relations

This section is dedicated to the application of the transitive closure computation algorithm from the previous section (Fig. 1) to three classes of arithmetic relations, for which the transitive closure is Presburger-definable: difference bounds relations [START_REF] Bozga | Flat parametric counter automata[END_REF], octagonal relations [START_REF] Bozga | Iterating octagons[END_REF], and finite monoid affine transformations [START_REF] Boigelot | Symbolic Methods for Exploring Infinite State Spaces[END_REF].

In order to apply the transitive closure computation method, one needs to address two issues. First, the class of relations considered needs to be proved ultimately periodic (or else, our algorithm is not guaranteed to terminate). The proofs rely mostly on the fact that any matrix A is ultimately periodic in T [START_REF] Schutter | On the ultimate behavior of the sequence of consecutive powers of a matrix in the max-plus algebra[END_REF] (see Section 2 for the definition of ultimately periodic matrices).

Second, the queries Q 1 and Q 2 (Fig. 1) need to be answered efficiently, by avoiding excessive calls to external decision procedures. In theory, all these queries can be expressed in Presburger arithmetic, for the classes of difference constraints, octagons and affine transformations, yet in practice we would like to avoid as much as possible using Presburger solvers, due to reasons of high complexity. For the classes of difference bounds and octagons, we give direct decision methods for handling these queries. The class of affine transformations without guards can also be dealt with by simply checking equivalence between Diophantine systems, whereas the general case still needs to be handled by a Presburger solver.

Difference Constraints

Let x = {x 1 , x 2 , ..., x N } be a set of variables ranging over Z.

Definition 3. A formula φ(x) is a difference bounds constraint if it is equivalent to a finite conjunction of atomic propositions of the form x i -x j ≤ a ij , 1 ≤ i, j ≤ N, i = j, where a ij ∈ Z.
For example, x = y + 5 is a difference bounds constraint, as it is equivalent to x -y ≤ 5 ∧ y -x ≤ -5. Let R db denote the class of difference bound relations. Difference bounds constraints are alternatively represented as matrices or, equivalently, weighted graphs.

Given a difference bounds constraint φ, a difference bounds matrix (DBM) representing φ is a matrix

M φ ∈ T N ×N such that (M φ ) ij = a ij , if x i -x j ≤ a ij is an atomic proposition in φ, and ∞, otherwise. Dually, if M ∈ T N ×N is a DBM, the corresponding difference bounds constraint is ∆ M ≡ Mij <∞ x i -x j ≤ M ij .
A DBM M is said to be consistent if and only if its corresponding constraint ϕ M is consistent. An elementary path in a DBM M is a sequence of indices 1 ≤ i 1 , i 2 , . . . , i k ≤ N , where i 1,...,k-1 are pairwise distinct, such that M ij ij+1 < ∞, for all 1 ≤ j < k. An elementary path is called an elementary cycle if moreover i 1 = i k . An elementary cycle is said to be strictly negative if k-1 j=1 M ij ij+1 < 0. A DBM M is inconsistent if and only if it has a strictly negative elementary cycle -a proof can be found in [START_REF] Miné | The octagon abstract domain[END_REF]. The next definition gives a canonical form for consistent DBMs.

Definition 4. A consistent DBM M ∈ T N ×N is said to be closed if and only if M ii = 0 and M ij ≤ M ik + M kj , for all 1 ≤ i, j, k ≤ N .
Given a consistent DBM M , we denote by M * the (unique) closed DBM such that ϕ M ⇔ ϕ M * . It is well-known that, if M is consistent, then M * is unique, and can be computed from M in time O(N 3 ), by the classical Floyd-Warshall algorithm. Moreover, if M is a consistent DBM, we have, for all 1 ≤ i, j ≤ N :

M * ij = min k-1 l=0 M i l i l+1 i = i 0 . . . i k-1 = j is an elementary path in M (1)
The closed form of DBMs is needed for the elimination of existentially quantified variables -if φ is a difference bounds constraint, then ∃x . φ is also a difference bounds constraint [START_REF] Miné | The octagon abstract domain[END_REF]. Consequently, we have that the class of difference bounds relations is closed under relational composition:

R 1 (x, x ′ )•R 2 (x, x ′ ) ≡ ∃y . R 1 (x, y)∧R 2 (y, x ′ ).
Difference Bounds Relations are Ultimately Periodic Given a consistent difference bounds relation With these definitions, the algorithm in Fig. 1 will return the transitive closure of a difference bounds relation R, provided that the sequence {σ(R i )} ∞ i=0 is ultimately periodic. If R is not ω-consistent then, by Def. 2, it is ultimately periodic. We consider henceforth that R is ω-consistent, i.e. σ(R i ) = M R i , for all i ≥ 0.

R(x, x ′ ) ∈ R db , let σ(R) = M R ∈ T 2N ×2N be
For a difference bounds relation R, we define the directed graph G R , whose set of vertices is the set x ∪ x ′ , and in which there is an edge from x i to x j labeled α ij if and only if the atomic proposition

x i -x j ≤ α ij occurs in R. Clearly, M R is the incidence matrix of G R .
Next, we define the concatenation of G R with itself as the disjoint union of two copies of G R , in which the x ′ vertices of the second copy overlap with the x vertices of the first copy. Then R m corresponds to the graph G m R , obtained by concatenating the graph of R to itself m > 0 times. Since R db is closed under relational composition, then R m ∈ R db , and moreover we have:

1≤i,j≤N x i -x j ≤ min{x 0 i - → x 0 j } ∧ x ′ i -x ′ j ≤ min{x m i - → x m j } ∧ x i -x ′ j ≤ min{x 0 i - → x m j } ∧ x ′ i -x j ≤ min{x m i - → x 0 j }
where min{x p i -→ x q j } is the minimal weight of all paths between the extremal vertices x p i and x q j in G m R , for p, q ∈ {0, m}. In other words, we have the equalities from Fig. 2 (a), for all 1 ≤ i, j ≤ N .

(MRm )i,j = min{x 0 i -→ x 0 j } (MRm )i+N,j+N = min{x m i -→ x m j } (MRm )i,j+N = min{x 0 i -→ x m j } (MRm )i+N,j = min{x m i -→ x 0 j } (a) min{x 0 i -→ x 0 j } = (M m R ) I ef (x i ),F ef (x j ) min{x m i -→ x m j } = (M m R ) I eb (x i ),F eb (x j ) min{x 0 i -→ x m j } = (M m R ) I of (x i ),F of (x j ) min{x m i -→ x 0 j } = (M m R ) I ob (x i ),F ob (x j ) (b)
Fig. 2 As proved in [START_REF] Bozga | Flat parametric counter automata[END_REF], the paths between x p i and x q j , for arbitrary 1 ≤ i, j ≤ N and p, q ∈ {0, m}, can be seen as words (over a finite alphabet of subgraphs of G m R ) recognized by a finite weighted automaton of size up to 5 N . For space reasons, the definition of this automaton is detailed in [START_REF] Bozga | Fast Acceleration of Ultimately Periodic Relations[END_REF].

Let M R be the incidence matrix of this automaton. By the construction of M R , for each variable x ∈ x, there are eight indices, denoted as5 I of (x), I ob (x), I ef (x), I eb (x), F of (x), F ob (x), F ef (x), F eb (x) ∈ {1, . . . , 5 N }, such that all relations from Fig. 2 (b) hold, for all 1 ≤ i, j ≤ N . Intuitivelly, all paths from x 0 i to x 0 j are recognized by the automaton with I ef (x i ) and F ef (x j ) as the initial and final states, respectivelly. The same holds for the other pairs of indices, from Fig. 2 (b). It is easy to see (as an immediate consequence of the interpretation of the matrix product in T) that, for any m > 0, the matrix M m R gives the minimal weight among all paths, of length m, between x p i and x q j , for any 1 ≤ i, j ≤ N and p, q ∈ {0, m}. But the sequence {M m R } ∞ m=0 is ultimately periodic, since every matrix is ultimately periodic in T [START_REF] Schutter | On the ultimate behavior of the sequence of consecutive powers of a matrix in the max-plus algebra[END_REF]. By equating the relations from Fig. 2 (a) with the ones from Fig. 2 (b), we obtain that the sequence

{σ(R m )} ∞ m=0 = {M R m } ∞ m=0
is ultimately periodic as well. In conclusion, the algorithm from Fig. 1 will terminate on difference bounds relations. Moreover, the result is a formula in Presburger arithmetic. This also simplifies the proof that transitive closures of difference bounds relations are Presburger definable, from [START_REF] Bozga | Flat parametric counter automata[END_REF], since the minimal paths of length m within the weighted automaton recognizing the paths of G m R correspond in fact to elements of the m-th power of M R (the incidence matrix of the automaton) in T.

Checking ω-Consistency and Inconsistency of Difference Bounds Relations For a difference bounds relation R(x, x ′ ) ∈ R db and a matrix Λ ∈ T 2N ×2N , we give methods to decide the queries Q 1 and Q 2 (lines 7 and 9 in Fig. 1) efficiently. To this end, we consider the class of parametric difference bounds relations. From now on, let k ∈ x be a variable interpreted over N + .

Definition 5. A formula φ(x, k) is a parametric difference bounds constraint if it is equivalent to a finite conjunction of atomic propositions of the form x

i -x j ≤ a ij • k + b ij , for some 1 ≤ i, j ≤ N , i = j, where a ij , b ij ∈ Z.
The class of parametric difference bounds relations with parameter k is denoted as R db (k). A parametric difference bounds constraint φ(k) can be represented by a matrix M φ [k] of linear terms, where

(M φ [k]) ij = a ij • k + b ij if x i -x j ≤ a ij • k + b ij occurs in φ,
M (k) ≡ M [k]ij =∞ x i -x j ≤ M [k] ij . With these considerations, we define π(M [k]) = ∆ M (k). Clearly, we have π(k • Λ + σ(R b )) ∈ R db (k), for R ∈ R db , b ≥ 0 and Λ ∈ T 2N ×2N .
Parametric DBMs do not have a closed form, since in general, the minimum of two linear terms in k (for all valuations of k) cannot be expressed again as a linear term. According to [START_REF] Annichini | Symbolic techniques for parametric reasoning about counter and clock systems[END_REF], one can define the closed form of a parametric DBM as a matrix of terms of the form min{a i •k+b i } m i=1 , for some a i , b i ∈ Z and m > 0. Then the query Q 1 can be written as a conjunction of formulae of the form ∀k > 0 .

min{a i • k + b i } m i=1 = a 0 • k + b 0 .
The following lemma gives a way to decide the validity of such formulae:

Lemma 4. Given ℓ, a 0 , a 1 , . . . , a m , b 0 , b 1 , . . . , b m ∈ Z, the following are equivalent: 1. ∀k ≥ ℓ . min{a i • k + b i } m i=1 = a 0 • k + b 0 2. m i=1 (a i = a 0 ∧ b i = b 0 ) ∧ m j=1 (a 0 ≤ a j ∧ a 0 • ℓ + b 0 ≤ a j • ℓ + b j )
In analogy to the non-parametric case, the inconsistency of a parametric difference bounds constraint φ(k) amounts to the existence of a strictly negative elementary cycle in M φ [k], for some valuation k ∈ N + . We are also interested in finding the smallest value for which such a cycle exists. The following lemma gives this value.

Lemma 5. Let φ(x, k) be a parametric difference bounds constraint and M φ [k] be its associated matrix. For some

a ij , b ij ∈ Z, let {a ij • k + b ij } mi j=1 , i = 1, .
. . , 2N be the set of terms denoting weights of elementary cycles going through i. Then φ is inconsistent for some ℓ ∈ N and k ≥ ℓ if and only if there exists

1 ≤ i ≤ 2N and 1 ≤ j ≤ m i such that either (i) a ij < 0 or (ii) a ij ≥ 0 ∧ a ij • ℓ + b ij < 0 holds. Moreover, the smallest value for which φ becomes inconsistent is min 2N i=1 {min mi j=1 γ ij }, where γ ij = max(ℓ, ⌊- bij aij ⌋ + 1), if a ij < 0, γ ij = ℓ, if a ij ≥ 0 ∧ a ij • ℓ + b ij < 0, and γ ij = ∞, otherwise.

Octagons

Let x = {x 1 , x 2 , ..., x N } be a set of variables ranging over Z.

Definition 6. A formula φ(x) is an octagonal constraint if it is equivalent to a finite conjunction of terms of the form ±x

i ± x j ≤ a ij , 2x i ≤ b i , or -2x i ≤ c i , where a ij , b i , c i ∈ Z and 1 ≤ i, j ≤ N, i = j.
The class of octagonal relations is denoted by R oct in the following. We represent octagons as difference bounds constraints over the set of variables y = {y 1 , y 2 , . . . , y 2N }, with the convention that y 2i-1 stands for x i and y 2i for -x i , respectively. For example, the octagonal constraint

x 1 + x 2 = 3 is represented as y 1 -y 4 ≤ 3 ∧ y 2 -y 3 ≤ -3.
To handle the y variables in the following, we define ī = i -1, if i is even, and ī = i + 1 if i is odd. Obviously, we have ī = i, for all i ∈ Z, i ≥ 0. We denote by φ the difference bounds formula φ[y 1 /x 1 , y 2 / -x 1 , . . . , y 2n-1 /x n , y 2n / -x n ] over y. The following equivalence relates φ and φ :

φ(x) ⇔ (∃y 2 , y 4 , . . . , y 2N . φ ∧ N i=1 y 2i-1 + y 2i = 0)[x 1 /y 1 , . . . , x n /y 2N -1 ] (2)
An octagonal constraint φ is equivalently represented by the DBM M φ ∈ T 2N ×2N , corresponding to φ. We say that a DBM M ∈ T 2N ×2N is coherent iff M ij = Mj ī for all 1 ≤ i, j ≤ 2N . This property is needed since any atomic proposition x i -x j ≤ a, in φ can be represented as both y 2i-1 -y 2j-1 ≤ a and y 2j -

y 2i ≤ a, 1 ≤ i, j ≤ N . Dually, a coherent DBM M ∈ T 2N ×2N corresponds to the octagonal constraint Ω M : 1≤i,j≤N (x i -x j ≤ M 2i-1,2j-1 ∧ x i + x j ≤ M 2i-1,2j ∧ -x i -x j ≤ M 2i,2j-1 ) (3)
A coherent DBM M is said to be octagonal-consistent if and only if Ω M is consistent. Definition 7. An octagonal-consistent coherent DBM M ∈ T 2N ×2N is said to be tightly closed if and only if the following hold:

1. M ii = 0, ∀1 ≤ i ≤ 2N 3. M ij ≤ M ik + M kj , ∀1 ≤ i, j, k ≤ 2N 2. M iī is even, ∀1 ≤ i ≤ 2N 4. M ij ≤ ⌊ Miī 2 ⌋ + ⌊ Mj j 2 ⌋, ∀1 ≤ i, j ≤ 2N
The following theorem from [START_REF] Bagnara | An improved tight closure algorithm for integer octagonal constraints[END_REF] provides an effective way of testing consistency and computing the tight closure of a coherent DBM. Moreover, it shows that the tight closure of a given DBM is unique and can also be computed in time O(N 3 ).

Theorem 2. [3] Let M ∈ T 2N ×2N be a coherent DBM. Then M is octagonal-consistent if and only if M is consistent and ⌊ Miī 2 ⌋ + ⌊ Mj j 2 ⌋ ≥ 0, for all 1 ≤ i, j ≤ 2N, i = j. Moreover, the tight closure of M is the DBM M t ∈ T 2N ×2N defined as M t ij = min M * ij , M * iī 2 + M * jj 2 , for all 1 ≤ i, j ≤ 2N , where M * ∈ T 2N ×2N is the closure of M .
The tight closure of an octagonal constraint is needed for existential quantifier elimination, and ultimately, for proving that the class of octagonal relations is closed under composition [START_REF] Bozga | Iterating octagons[END_REF].

Octagonal Relations are Ultimately Periodic Given a consistent octagonal relation

R(x, x ′ ) let σ(R) = M R . Dually, for any coherent DBM M ∈ T 4N ×4N , let ρ(M ) = Ω M . Clearly, ρ(σ(R)) ⇔ R, as required by Def. 2.
In order to prove that the class R oct of octagonal relations is ultimately periodic, we need to prove that the sequence {σ(R m )} ∞ m=0 is ultimately periodic, for an arbitrary relation R ∈ R oct . It is sufficient to consider only the case where R is ω-consistent, hence σ(R m ) = M R m , for all m ≥ 0. We rely in the following on the main result of [START_REF] Bozga | Iterating octagons[END_REF], which establishes a relation between M R m (the octagonal DBM corresponding to the m-th iteration of R) and M R m (the DBM corresponding to the m-th iteration of R ∈ R db ), for m > 0:

(M R m ) ij = min (M R m ) ij , (M R m )iī 2 + (M R m )j j 2 , for all 1 ≤ i, j ≤ 4N ( * )
This relation is in fact a generalization of the tight closure expression from theorem 2, from m = 1 to any m > 0.

In Section 4.1 it was shown that difference bounds relations are ultimately periodic. In particular, this means that the sequence {M R m } ∞ m=0 , corresponding to the iteration of the difference bounds relation R, is ultimately periodic. To prove that the sequence {M R m } ∞ m=0 is also ultimately periodic, it is sufficient to show that: the minimum and the sum of two ultimately periodic sequences are ultimately periodic, and also that the integer half of an ultimately periodic sequence is also ultimately periodic. Together with the above relation ( * ), lemma 6 proves that R oct is ultimately periodic.

Checking ω-Consistency and Inconsistency of Octagonal Relations This section describes an efficient method of deciding the queries Q 1 and Q 2 (lines 7 and 9 in Fig. 1) for the class of octagonal relations. In order to deal with these queries symbolically,we need to consider first the class R oct (k) of octagonal relations with parameter k. In the rest of this section, let k ∈ x be a variable ranging over N + . Definition 8. Then a formula φ(x, z) is a parametric octagonal constraint if it is equivalent to a finite conjunction of terms of the form

±x i ±x j ≤ a ij •k+b ij , 2x i ≤ c i •k+d i , or -2x i ≤ c ′ i • k + d ′ i , where a ij , b ij , c i , d i , c ′ i , d ′ i ∈ Z and 1 ≤ i, j ≤ N, i = j.
A parametric octagon φ(x, k) is represented by a matrix M φ [k]T[k] 2N ×2N of linear terms over k, and viceversa, a matrix M [k] ∈ T[k] 2N ×2N corresponds to a parametric octagon Ω M (k). We define π(M [k]) = Ω M (k). As in the case of difference bounds constraints, one notices that π(k

• Λ + σ(R b )) ∈ R oct (k), for R ∈ R oct , b ≥ 0 and Λ ∈ T 4N ×4N .
The composition of parametric octagonal relations (from e.g. Q 1 ) requires the computation of the tight closure in the presence of parameters. According to theorem 2, the parametric tight closure can be expressed as a matrix of elements of the form min{t i (k)} m i=1 , where t i (k) are either: (i) linear terms, i.e.

t i (k) = a i • k + b i , or (ii) sums of halved linear terms, i.e. t i (k) = ⌊ ai•k+bi 2 ⌋ + ⌊ ci•k+di 2 ⌋.
The main idea is to split a halved linear term of the form

⌊ a•k+b 2 ⌋, k > 0 into two linear terms a • k + ⌊ b 2 ⌋ and a • k + ⌊ b-a 2 ⌋
, corresponding to the cases of k > 0 being even or odd, respectivelly. This is justified by the following equivalence:

{⌊ a•k+b 2 ⌋ | k > 0} = {a • k + ⌊ b 2 ⌋ | k > 0} ∪ {a • k + ⌊ b-a 2 ⌋ | k > 0}
Hence, an expression of the form min{t i (k)} m i=1 yields two expressions min{t e i (k)} m i=1 , for even k, and min{t o i (k)} m i=1 , for odd k, where t e i and t o i , 1 ≤ i ≤ m, are effectively computable linear terms. With these considerations, Q 1 (for octagonal relations) is equivalent to a conjunction of equalities of the form ∀k > 0 . min{t

• i (k)} m i=1 = t • 0 (k),
• ∈ {e, o}. Now we can apply lemma 4 to the right-hand sides of the equivalences above, to give efficient equivalent conditions for deciding Q 1 .

The query Q 2 is, according to theorem 2, equivalent to finding either (i) a strictly negative cycle in a parametric octagonal DBM M [k], or (ii) a pair of indices

1 ≤ i, j ≤ 4N, i = j such that ⌊ M [k]iī 2 ⌋ + ⌊ M [k]j j 2 ⌋ < 0.
Considering that the set of terms corresponding to the two cases above is

T = {a i • k + b i } m i=1 ∪ {⌊ ci•k+di 2 ⌋ + ⌊ ei•k+fi

2

⌋} p i=1 , we split each term t ∈ T into two matching linear terms, and obtain, equivalently:

T e,o = {α e i • k + β e i } m+p i=1 ∪ {α o i • k + β o i } m+p i=1
Now we can apply lemma 5, and compute the minimal value for which a term t ∈ T e,o becomes negative, i.e. n 0 = min m+p i=1 min(2γ e i , 2γ o i -1), where γ

• i = max(1, ⌊- β • i α • i ⌋+ 1), if α • i < 0, 1 if α • i ≥ 0 ∧ α • i + β • i < 0, and
∞, otherwise, for • ∈ {e, o}.

Finite Monoid Affine Transformations

The class of affine transformations is one of the most general classes of deterministic transition relations involving integer variables. If x = x 1 , . . . , x N is a vector of variables ranging over Z, an affine transformation is a relation of the form:

T ≡ x ′ = A ⊗ x + b ∧ φ(x) (4)
where A ∈ Z N ×N , b ∈ Z N , φ is a Presburger formula, and ⊗ stands for the standard matrix multiplication in Z.

The affine transformation is said to have the finite monoid property [START_REF] Boigelot | Symbolic Methods for Exploring Infinite State Spaces[END_REF][START_REF] Finkel | How to compose presburger-accelerations: Applications to broadcast protocols[END_REF] if the monoid M A , ⊗ , where M A = {A ⊗ i | i ≥ 0} is finite. In this case, we also say that A is finite monoid. Here A ⊗ 0 = I N and A ⊗ i = A ⊗ A ⊗ i-1 , for i > 0. Intuitivelly, the finite monoid property is equivalent to the fact that A has finitely many powers (for the standard integer multiplication) that repeat periodically. It is easy to see that A is finite monoid if and only if there exists p ≥ 0 and l > 0 such that A ⊗ p = A ⊗ p+l , i.e. M A = {A ⊗ 0 , . . . , A ⊗ p , . . . , A ⊗ p+l-1 }.

If A is finite monoid, it can be shown that T * can be defined in Presburger arithmetic [START_REF] Boigelot | Symbolic Methods for Exploring Infinite State Spaces[END_REF][START_REF] Finkel | How to compose presburger-accelerations: Applications to broadcast protocols[END_REF]. We achieve the same result below, by showing that finite monoid affine transformations are ultimately periodic relations. As a byproduct, the transitive closure of such relations can also be computed by the algorithm in Fig. 1.

An affine tranformation T (4) can be equivalently written in the homogeneous form:

T ≡ x ′ h = A h ⊗ x h ∧ φ h (x h ) where A h ≡ A b 0 . . . 0 1 where x h = x 1 , . . . x N , x N +1 with x N +1 ∈ x being a fresh variable and φ h (x h ) ≡ φ(x) ∧ x N +1 = 1.
In general, the k-th iteration of an affine transformation can be expressed as:

T k ≡ x ′ h = A h ⊗ k ⊗ x h ∧ ∀0 ≤ ℓ < k . φ h (A h ⊗ ℓ ⊗ x h ) (5)
Notice that, if x

h denotes the initial values of x h , the values of x h at the ℓ-th iteration are x

(ℓ) h = A h ⊗ ℓ ⊗ x (0)
h . Moreover, we need to ensure that all guards up to (and including) the (k -1)-th step are satisfied, i.e. φ h (A h ⊗ ℓ ⊗ x h ), for all 0 ≤ ℓ < k. For the rest of the section we fix A and b, as in [START_REF] Bardin | Fast extended release[END_REF]. The encoding of a consistent affine transformation T is defined as σ(T ) = A h ∈ T (N +1)×(N +1) . Dually, for some M ∈ T[k] (N +1)×(N +1) , we define:

π(M ) : ∃x N +1 , x ′ N +1 . x ′ h = M ⊗ x h ∧ ∀0 ≤ ℓ < k . φ h (M [ℓ/k] ⊗ x h )
where M [ℓ/k] denotes the matrix M in which each occurrence of k is replaced by ℓ. In contrast with the previous cases (Section 4.1 and Section 4.2), only M is not sufficient here to recover the relation π(M ) -φ needs to be remembered as well 6 .

With these definitions, we have σ(T k ) = A ⊗ h k , for all k > 0 -as an immediate consequence of (5). The next lemma proves that the class of finite monoid affine relations is ultimately periodic.

Lemma 7. Given a finite monoid matrix

A ∈ Z N ×N and a vector b ∈ Z N , the se- quence {A ⊗ h k } ∞ k=0 is ultimately periodic.
The queries Q 1 and Q 2 (lines 7 and 9 in Fig. 1) in the case of finite monoid affine transformations, are expressible in Presburger arithmetic. These problems could be simplified in the case of affine transformations without guards, i.e T ≡ x ′ = Ax + b. The transformation is, in this case, ω-consistent. Consequently, Q 1 reduces to an equivalence between two homogeneous systems x ′ h = A 1h ⊗ x h and x ′ h = A 2h ⊗ x h . This is true if and only if A 1h = A 2h . The query Q 2 becomes trivially false in this case.

We have implemented the transitive closure algorithm from Fig. 1 within the FLATA toolset [11], a framework we develop for the analysis of counter systems. We compared the performance of this algorithm with our older transitive closure computation methods for difference bounds [START_REF] Bozga | Flat parametric counter automata[END_REF] and octagonal relations [START_REF] Bozga | Iterating octagons[END_REF]. We currently lack experimental data for finite monoid relations (namely, a comparison with existing tools such as FAST [START_REF] Bardin | Fast extended release[END_REF], LASH [START_REF] Wolper | Verifying systems with infinite but regular state spaces[END_REF] or TReX [START_REF] Annichini | Trex: A tool for reachability analysis of complex systems[END_REF] on this class), as our implementation of finite monoid affine transformation class is still underway.

Table 1 shows the results of the comparison between the older algorithms described in [START_REF] Bozga | Flat parametric counter automata[END_REF][START_REF] Bozga | Iterating octagons[END_REF] (denoted as old) and the algorithm in Fig. 1 for difference bounds relations d 1,...,6 and octagonal relations o 1,...,6 . The tests have been performed on both compact (minimum number of constraints) and canonical (i.e. closed, for difference bounds and tightly closed, for octagons) relations. The speedup column gives the ratio between the old and new execution times. The experiments were performed on a 2.53GHz machine with 2.9GB of memory. (x ≥ 0) ∧ (y ≥ 0) ∧ (x ′ ≥ 0) ∧ (y ′ ≥ 0) 2.5 > 10 6 ∞ > 10 6 ∞ ∧(x + y ≤ 1) ∧ (x ′ + y ′ ≤ 1) ∧ (x -1 ≤ x ′ ) ∧(x ′ ≤ x + 1) ∧ (y -1 ≤ y ′ ) ∧ (y ′ ≤ y + 1) Table 1: Comparison with older algorithms on difference bounds and octagons. Times are in milliseconds.

As shown in Table 1, the maximum observed speedup is almost 10 5 for difference bounds (d 4 in canonical form) and of the order of four for octagons. For the relations d 5 (canonical form), d 6 and o 6 the computation using older methods took longer than 10 6 msec. It is also worth noticing that the highest execution time with the new method was of 2.5 msec.

Table 2 compares FLATA with FAST [START_REF] Bardin | Fast extended release[END_REF] on counter systems with one self loop labeled with a randomly generated deterministic difference bound relation. We generated 50 such relations for each size N = 10, 15, 20, 25, 50, 100. Notice that FAST usually runs out of memory for more than 25 variables, whereas FLATA can handle 100 variables in reasonable time (less than 8 seconds on average).

1 .

 1 foreach b := 0, 1, 2, . . . do 2. foreach c := 1, 2, . . . , b do 3.

  the characteristic DBM of R, and for any M ∈ T 2N ×2N , let ρ(M ) = ∆ M ∈ R db be the difference bounds relation corresponding to R. Clearly, ρ(σ(R)) ⇔ R, as required by Def. 2.

  and ∞ otherwise. Dually, a matrix M [k] of linear terms corresponds to the formula ∆

Lemma 6 . 2 ∞

 62 Let {s m } ∞ m=0 and {t m } ∞ m=0 be two ultimately periodic sequences. Then the sequences {min(s m , t m )} ∞ m=0 , {s m + t m } ∞ m=0 and sm m=0 are ultimately periodic as well.

∞ d 6 (o 2 o 3

 623 Relationnew compact canonical old speedup old speedup d 0 (xx ′ = -1) ∧ (x = y ′ ) 0.18 0.7 3.89 38.77 215.39d 1 (xx ′ = -1) ∧ (x ′ = y ′ )0.18 18.18 101.0 38.77 215.39d 2 (xx ′ = -1) ∧ (x = y ′ ) ∧ (xz ′ ≤ 5) ∧ (z = z ′ ) 1.2 26.5 22.1 33431.2 27859.3 d 3 (xx ′ = -1) ∧ (x = y ′ ) ∧ (xz ≤ 5) ∧ (z = z ′ )0.6 32.7 54.5 33505.5 55841.7d 4 (xx ′ = -1) ∧ (x = y) ∧ (xz ≤ 5) ∧ (z = z ′ ) 0.5 702.3 1404.6 48913.8 97827.6 d 5 (a = c) ∧ (b = a ′ ) ∧ (b = b ′ ) ∧ (c = c ′ ) 1.8 5556.6 3087.0 > 10 6 ab ′ ≤ -1) ∧ (ae ′ ≤ -2) ∧ (ba ′ ≤ -2) 5.6 > 10 6 ∞ > 10 6 ∞ ∧(b -c ′ ≤ -1) ∧ (cb ′ ≤ -2) ∧ (cd ′ ≤ -1) ∧(d -c ′ ≤ -2) ∧ (de ′ ≤ -1 ∧ ea ′ ≤ -1) ∧(e -d ′ ≤ -2) ∧ (a ′ -b ≤ 4) ∧ (a ′ -c ≤ 3) ∧(b ′ -c ≤ 4 ∧ b ′ -d ≤ 3) ∧ (c ′ -d ≤ 4) ∧ (c ′ -e ≤ 3) ∧(d ′ -a ≤ 3 ∧ d ′ -e ≤ 4) ∧ (e ′ -a ≤ 4) ∧ (e ′ -b ≤ 3) o 1 (x + x ′ = 1) (x + y ′ ≤ -1) ∧ (-yx ′ ≤ -2) (x ≤ x ′ ) ∧ (x + y ′ ≤ -1) ∧ (-yx ′ ≤ -2) o 4 (x + y ≤ 5) ∧ (-x + x ′ ≤ -2) ∧ (-y + y ′ ≤ -3) 0.21 3.67 17.48 13.52 64.38 o 5 (x + y ≤ 1) ∧ (-x ≤ 0) ∧ (-y ≤ 0) 1.2 20050.9 16709.1 > 10 6 ∞ o 6

A class of relations is usually defined by syntactic conditions.

The nested loop from Fig.1will always yield a pair (b, c) such that b ≥ c. To ensure that b is also minimal, and thus cover up the case b < c, once the smallest period c has been detected at prefix b = c, we need to also try all prefixes b = c -1, c -2, . . . , 0.

Paths between x 0 and y m (x m and y 0 ) are called odd forward (backward) in[START_REF] Bozga | Flat parametric counter automata[END_REF], whereas paths between x 0 and y 0 (x m and y m ) are called even forward (backward). Hence the indices of , ob, ef and eb.

This incurs a slight modification of the algorithm presented in Fig.1.[START_REF] Schutter | On the ultimate behavior of the sequence of consecutive powers of a matrix in the max-plus algebra[END_REF]

⋆ This work was supported by the French national project ANR-09-SEGI-016 VERIDYC, by the Czech Science Foundation (projects P103/10/0306 and 102/09/H042), the Czech Ministry of Education (projects COST OC10009 and MSM 0021630528), and the internal FIT BUT grant FIT-10-1.

Conclusion

We presented a new, scalable algorithm for computing the transitive closure of ultimately periodic relations. We show that this algorithm is applicable to difference bounds, octagonal and finite monoid affine relations, as all three classes are shown to be ultimately periodic. Experimental results show great improvement in the time needed to compute transitive closures of difference bounds and octagonal relations.


P Bourgos
email: bourgos@imag.fr

A Basu
email: basu@imag.fr

M Bozga
email: bozga@imag.fr

S Bensalem
email: bensalem@imag.fr

J Sifakis
email: sifakis@imag.fr

K Huang
email: huangk@vlsi.zju.edu.cn

Rigorous System Level Modeling and Analysis of Mixed HW/SW Systems

A grand challenge in complex embedded systems design is developing methods and tools for modeling and analyzing the behavior of an application software running on multicore or distributed platforms. We propose a rigorous method and a tool chain that allows to obtain a faithful model representing the behavior of a mixed hardware/software system from a model of its application software and a model of its underlying hardware architecture. The system model can be simulated and analyzed for validation of both functional and extra-functional properties. The tool chain uses DOL (Distributed Operation Layer [1]) as the frontend for specifying the application software and hardware architecture, and BIP (Behavior Interaction Priority [2]) as the modeling and analysis framework. It is illustrated through the construction of system models of MJPEG and MPEG2 decoder applications running on MPARM, a multicore architecture.

I. INTRODUCTION

Performance of embedded applications strongly depends on features of the underlying hardware platform. In contrast to performance of application software running on a single core, getting the maximum throughput out of multicore processors demands application software to be designed taking parallelism into account from scratch. This is needed to catch up with the fast growth of computing capacity due to the foreseeable exponential increase of physical parallelism. But programming, testing and verifying parallel software with currently existing tools is notoriously hard, even for experts. There are no rigorous techniques for deriving global model of a given system from models of its application software and its execution platform.

Application software must be programmed for performance, in a platform independent way, exhibiting all potential parallelism. Its implementation must deal with mapping the specified application-level parallelism onto platformlevel (threads, cores, processors) on an as-needed/asavailable basis. Actually, this mapping would need to be adapted dynamically as applications must scale up or down

The research leading to these results has received funding from the European Community's Seventh Framework Programme [FP7/ [2007][2008][2009][2010][2011][2012][2013] under grant agreement no 248776 (PRO3D) and from ARTEMIS JU grant agreement ARTEMIS-2009-1-100230 (SMECY) according to the available resources of the execution platform. Moreover, efficiency and correctness are not the only concerns. Programmer productivity, that is, the programmer's ability to design correct software that gathers the maximum performance out of an arbitrary multicore platform with ease should not be neglected [START_REF] Asanovic | The landscape of parallel computing research: A view from berkeley[END_REF].

Achieving these goals requires a design flow based on a single semantic model. The design flow must be able to generate rigorous models of mixed hardware/software systems, suitable for analysis, design space exploration and automatic code generation. The main contribution of this paper is deriving a rigorous system model combining the application software and the architecture, which can be the basis for multiple objectives, such as functional verification, performance evaluation and code generation for target architectures.

We propose a system construction method that is both rigorous and allows a fine analysis of system dynamics. It is rigorous because it is based on formal models, have precise semantics and thus can be analyzed by using formal techniques. A system model is derived by progressively integrating constraints induced on an application software model by the underlying hardware architecture model. Both models are described in BIP [START_REF] Basu | Modeling Heterogeneous Real-time Components in BIP[END_REF], which is a formal component based modeling framework. In contrast to ad hoc modeling approaches, the system model is obtained from a BIP model of the application software and a description of the hardware architecture, by application of source-tosource transformations that are correct-by-construction [START_REF] Bonakdarpour | From high-level component-based models to distributed implementations[END_REF]. The final generated model is a mixed software-hardware model which provides the capability using a single model to simulate and apply formal verification techniques on it using the BIP framework.

Metro II [START_REF] Abhijit | A next-generation design framework for platform-based design[END_REF] is a platform-based design framework and provides a simulation backend based on SystemC. Octopus [START_REF] Twan | Model-driven design-space exploration for embedded systems: The octopus toolset[END_REF] allows design space exploration by stochastic simulation of task graphs. Both have connections to formal verification tools based on model checking. Most of the frameworks for mixed HW/SW systems are based on SystemC [START_REF] Grtker | System Design with SystemC[END_REF] as a language for modeling at various levels of abstractions. Various tools and associated design methodologies emerged e.g., SystemCoDesigner [START_REF] Haubelt | Systemcodesigner: automatic design space exploration and rapid prototyping from behavioral models[END_REF], Spade [START_REF] Lieverse | System level design with SPADE: an M-JPEG case study[END_REF], Sesame [START_REF] Erbas | A framework for system-level modeling and simulation of embedded systems architectures[END_REF] to cite only a few. All these focus and facilitate the construction of executable simulation models which, while being claimed cycle-accurate, do not rely on a formal foundation. For instance, such models cannot be used to check formally the correctness of the constructed system. There have been attempts on providing formal semantics to System-C models using tools like LusSy [START_REF] Moy | Lussy: A toolbox for the analysis of systems-on-a-chip at the transactional level[END_REF], however, they remain difficult to use mainly because of the limited expressiveness of the target formalism compared with a general purpose language.

One of the main needs for rigorous system model is performance evaluation. Simulation based methods use adhoc executable system models such as [START_REF] Kienhuis | An approach for quantitative analysis of applicationspecific dataflow architectures[END_REF] or models in SystemC [START_REF] Grtker | System Design with SystemC[END_REF], [START_REF] Moussa | Exploring sw performance using soc transaction-level modeling[END_REF]. The latter provide cycle-accurate results, but are not adequate for thorough exploration of hardware architecture dynamics and its effects on software execution. Furthermore, long simulation time is a major drawback. Trace-based co-simulation is used in Spade [START_REF] Lieverse | System level design with SPADE: an M-JPEG case study[END_REF], Sesame [START_REF] Erbas | A framework for system-level modeling and simulation of embedded systems architectures[END_REF]. There exist much faster techniques that work on abstract system models e.g., Real Time Calculus [START_REF] Thiele | Real-time calculus for scheduling hard real-time systems[END_REF] and SymTA/S [START_REF] Henia | System-level performance analysis -the SymTA/S approach[END_REF]. They use formal analytical models representing a system as a network of nodes exchanging streams. The dynamics of the execution platform is characterized by execution times. Nonetheless, these techniques allow only estimation of pessimistic worst-case measures (delays, buffer sizes, etc) and moreover, they require an abstract model of the application software. Building these abstract models represents a significant modeling effort and, if done through a manual process, the results are not guaranteed to be accurate. Similar drawbacks exists for performance analysis techniques based on Timed-Automata [START_REF] Salah | Compositional timing analysis[END_REF], [START_REF] Abdeddaim | Scheduling with timed automata[END_REF]. These can be used for modeling and solving scheduling problems. An approach combining simulation and analytic models is presented in [START_REF] Künzli | Combining simulation and formal methods for system-level performance analysis[END_REF], where simulation results can be propagated to analytic models and vice versa through well defined interfaces.

The paper is structured as follows. Section II presents the method and the main steps in the design flow, with a brief overview of the BIP framework and associated toolbox. The generation of the system model follows in section III. Section IV describes the performance estimation technique applied on the system model. Finally, experimental results are provided in section V. In section VI we conclude and discuss future work directions.

II. DESIGN FLOW

The flow of our method is illustrated in Figure 1. The method takes three inputs: (i) the application software, (ii) the hardware architecture and (iii) the mapping. We consider application software defined using the Kahn process network model [START_REF] Kahn | The semantics of a simple language for parallel programming[END_REF]. They consists of a set of deterministic processes communicating through FIFO channels by executing atomic read/write operations. The behavior of each process is a sequential program. We consider hardware architectures described as interconnections of computational and communication devices such as processors, buses and memories. Finally, we consider mappings that associate application software elements to hardware architecture, that is, processes to processors and FIFO channels to memories.

In this paper, we will focus on the generation of the system model. We will also describe one of its utilities, i.e., performance evaluation. The first stage of the method is the construction of the system model in BIP. The system model represents the application mapped on the hardware architecture. The system model is obtained by the three following steps:

1) the construction of a BIP model by automatic translation from the application software, 2) the construction of a BIP model by automatic translation from the hardware architecture, 3) the construction of the system model by source-tosource transformation of the previous two models and their composition according to the mapping. The second stage of the method is performance evaluation realized on the system model. We provide a simulation-based technique allowing the accurate estimation of real-time characteristics (response times, delays, latencies, throughputs, etc.) and particular indicators about the use of resources (bus conflicts, memory conflicts, etc.).

The performance evaluation method combines native (BIP) simulation of the system model with online code profiling on the target hardware architecture. That is, the (simulated) processing time required by the application code is computed during simulation, on demand, using the application object code for the target architecture and the processor weight table. The later provides the raw execution times for elementary (assembler) instructions.

The method is completely automated and has been implemented in a tool. The tool uses as inputs Distributed Operation Layer (DOL) [START_REF] Thiele | Mapping applications to tiled multiprocessor embedded systems[END_REF] specifications, that is, the application software, the hardware architecture and the mapping are described using the concrete formalisms available in the DOL framework. The method is realized using the BIP framework [START_REF] Basu | Modeling Heterogeneous Real-time Components in BIP[END_REF], [START_REF] Bliudze | A Notion of Glue Expressiveness for Component-Based Systems[END_REF], [START_REF] Basu | Rigorous component-based design using the BIP framework[END_REF] and the associated toolbox 1 . The BIP language is a notation which allows complex systems to be built by coordinating the behavior of a set of atomic components. The behavior is described as automata or Petri nets extended with data and functions described in C/C++. Transitions are labelled with ports (action names), guards (enabling conditions on the state of a component) as well as functions (computations on local data). The description of coordination between components is layered. It consists of interactions and priorities that characterizes the overall architecture of a component. Their combination confers BIP strong expressiveness that cannot be matched by other languages [START_REF] Bliudze | A Notion of Glue Expressiveness for Component-Based Systems[END_REF]. BIP has clean operational semantics that

III. DERIVING SYSTEM MODEL

The construction of the system model in BIP from the input DOL specification [START_REF] Thiele | Mapping applications to tiled multiprocessor embedded systems[END_REF] is done in three steps, as described in the following subsections.

A. Construction of Application Software Model in BIP

An application software in DOL [START_REF] Thiele | Mapping applications to tiled multiprocessor embedded systems[END_REF] is a process network that consists of three basic entities: SW-Process, SW-Channel, and SW-Connection, organized as described by the following abstract grammar: A write connection between output port j of a process P and a channel C is a pair (P.OutPort j , C.RecvPort). A read connection between input port i of process P and a channel C is a pair (C.SendPort, P.InPort i ). We assume that ports of channels are uniquely associated with ports of processes and vice versa.

Process behavior is described using C programs with a particular structure (see figure 3 for a concrete example). In general, the behavior of a process P is defined by an initial call of the P init() function followed by an endless loop calling the P fire() function. Communication is realized by using two particular primitives, namely write and read for respectively sending and receiving data to software channels. A read operation reads data from an input port, and a write operation writes data to an output port. The code may also call another special primitive, namely detach, in order to terminate the execution of the process. Example 1: An example process network is shown in figure 2. It has three SW-processes (generator, square and consumer), connected through two SW-channels (C1 and C2). The generator produces an integer and sends it to square, which squares it and send it to the consumer which prints the result. The description of square process is shown in figure 3. It defines the data structure for the process state, the function square init() to initialize the process state and the function square fire() to define the cyclic behavior of the process. The square process uses integer variables index and len. The function square fire defines a floating variable i, which holds the value read from the port IN. On every call of square fire, it reads a value for i, squares it, writes it to the port OUT and increments the counter index. The process terminates when index reaches len. The construction of the application software model in BIP is structural: every process and every channel are independently translated to atomic components in BIP and then connected according to their connections in the process network.

1) Translation of Software Processes into BIP:

The translation converts every software process to an atomic component in BIP. Each port is defined as a port in the atomic component. Data structures defined in the C functions are used as data in the atomic component. Control locations correspond to invocation of read/write primitives for which synchronization is required. Transitions are labeled by the port name associated with the primitives. Computation statements are added as actions of the transitions.

The translation requires the extraction of a control-flow graph from the C code. It starts by parsing the process code into an intermediate, annotated abstract syntax tree (AST). The translation to BIP is then completed in two steps. In the first step, the interaction points in the AST are identified, that is, each call to a read/write primitive is registered as an interaction point. The second step involves the construction of an explicit control flow graph and its representation as a finite state automaton extended with data in BIP. For every interaction point, a control location is created. An outgoing transition is added from this location, labeled by the port used in the read/write call. The transition models the primitive call and requires synchronization with a software channel.

Statements other than read/write calls are added as actions to the existing transitions. Let us notice that any functions that contain read/write calls (either directly or through nested calls) are inlined in the BIP automaton. Consequently, our translation is restricted to programs without communication calls occurring within recursive functions. Additional restrictions are, namely: no use of global variable; and no goto statement.

3) Translation of Connections into BIP:

Every connection in the application software is translated into a BIP connector which strongly synchronizes the corresponding ports. Connectors provide the transfer of data implementing the read and write operations. A connector implementing write transfers data from a process to a channel, whereas the one implementing read transfers data from a channel to a process.

B. Construction of Hardware Architecture Model in BIP

A hardware architecture consists of computational resources interconnected according to communication paths. Resources are used for computation (processors, memories) or for communication (buses). Communication paths define the connections between computational resources. More formally, we consider the family of hardware architectures that can be represented in DOL [START_REF] Thiele | Mapping applications to tiled multiprocessor embedded systems[END_REF] and are abstracted by the following grammar: The local memory of each tile is also connected to the shared bus. We consider the following three communication paths, ordered (write, read) as follows:

WP1 = ARM1.LB1.LM1 RP1 = LM1.LB1.ARM1 WP2 = ARM1.LB1.SB.SM RP2 = SM.SB.LB2.ARM2 WP3 = ARM2.LB2.LM2 RP3 = LM2.LB2.ARM2
The BIP model constructed from the hardware architecture represents explicitly, in an operational manner, the interconnect between the different resources as defined by the communication paths. This model is organized as a collection of bus, processor and memory components. Nonetheless, let us notice that, the processor and memory components are just empty, placeholder components. We introduce them in the BIP model of the hardware architecture only for the sake of clarity. They will be filled during the next step, that is, the construction of the system model.

Every bus component is concretely defined as a scheduled collection of communication path fragments. That is, for every read/write path going on a bus, we consider the path fragment defined by three atomic components, respectively:

• the MasterInterface (MI) component, which controls the access of the communication path on the bus and initiates the read/write operation. Depending on its position on the path, the master component receives data either from some software processes executing inside the processor or from the previous path segment. • the VirtualLink (VL) component, which models effectively the transfer of data over the bus, from the master once it gets access to the bus, towards the slave. • the SlaveInterface (SI) component, which acts like a buffer and is needed to connect further either to the next path fragment or to some FIFO buffers on the memory, depending on the position of the bus on the path. All the paths segments going over the same bus must share its transport capabilities according to some predefined bus policy. The scheduling can be of one of fixed-priority, round-robin or TDMA. We model it explicitly by using a HW-Bus-Scheduler component, which interacts with all the master interface components and ensures exclusive access for transmission of data, according to the policy selected. The HW-Bus-Scheduler acts as an arbiter to resolve the bus access conflicts.

All these components are predefined and belong to the BIP hardware library. They have identical interfaces for the transport of data, respectively ports RR/WR (Read/Write-Request), RA/WA (Read/Write-Acknowledge) to connect with upper components, and RB/WB (Read/Write-Begin), RE/WE (Read/Write-End) to connect with lower components on the path. In addition, the MI components use ports ACQ (Acquire) and REL (Release) to interact with the bus scheduler.

Finally, let us also notice that all these components are timed BIP components [START_REF] Basu | Modeling Heterogeneous Real-time Components in BIP[END_REF]. The VirtualLink components model the latency of the buffer. The Master/SlaveInterface components observe the time progress and can be used for observation purposes, as explained later in section IV.

Example Every connection is realized using BIP connectors which strongly synchronize the corresponding ports. The behavior of the connector implements the transfer of data, its address and size between the successive components, corresponding to the write and read operations.

Example 6: Figure 13 shows the BIP hardware model of the 2-Tile ARM architecture of example 4. Communication paths between the processors and the memories are implemented using the previously defined set of bus components.

C. Construction of the System Model in BIP

Given the BIP models of respectively the application software and hardware architecture, the construction of the BIP system model is completed in two steps:

1) transformation of components in the BIP application model, namely decomposing the SW-Channels into data buffers and read/write FIFO access routines, and consequently breaking the atomicity of the read/write operations in SW-Processes. 2) allocation of the transformed processes and FIFO routines on hardware processors and respectively data buffers on hardware memories according to the mapping, and eventually filling up the processor and memory placeholder components. Formally, the BIP system model conforms to the following abstract grammar: In order to deploy the application software on the architecture, we need a low level implementation model for the SW-Channels where the control and the data are dissociated and moreover, the read/write operations are no longer atomic.

Splitting software channels: Every SW-Channel in the application software is replaced by a composition of FIFO-Write, FIFO-Read and a FIFO-Buffer atomic components (figure 9). The two former components represent the control part of the software channel, that is, the hardware dependent software routines implementing the read/write operations. The latter component simply represents the buffer of data. All the three components FIFO-Read, FIFO-Write, FIFO-Buffer are predefined BIP components and belong to the BIP hardware dependent software library. The FIFO-Read component, illustrated in figure 10, implements the read operation on channels. It has the ports RR (Read-Request), RA (Read-Acknowledge) for its interaction with a software process read operation, and ports RB (Read-Begin), RE (Read-End) for its interaction with the buffer. The FIFO-Write component implements the write action in a similar manner. Let us notice that the two routines, FIFO-Write and FIFO-Read, require extra synchronization with each other in order to maintain a coherent value for the used space within the buffer. This is realized by using strong synchronization between two control ports, SIGSEM and UPDSEM. Moreover, they also use the ports REL and ACQ for interaction with the processor scheduler. These ports are used to release (resp. acquire) the processor whenever the read/write operation is suspended (resp. resumed) due to lack (resp. presence) of available data (or available space) in the buffer.

The FIFO-Buffer represents a passive component modeling the data storage. It has ports WB, WE and RB, RE for writing and reading respectively. The ports for writing (resp. reading) synchronizes with the FIFO-Write (resp. FIFO-Read) component.

We can prove that the proposed model is a correct implementation of the SW-Channel. That is, the composition is a refined model of the SW-Channel which fully preserves the input/output behavior of the software channel.

Transformation of software processes: The splitting of SW-Channels as described before requires the transformation of the software processes as well.

The first transformation consists in breaking atomicity of write and read operations. Every transition involving an input/output port X is split into two transitions, labeled by fresh ports, respectively XB (i.e., X-begin) and XE (i.e., Xend). This is obtained by adding new control locations for each read/write operations in the behavior of the process.

The second transformation, completely orthogonal to the first one, consists in adding interactions with the processor scheduler. This transformation is needed since several processes, together with their associated FIFO access routines, are potentially mapped on the same hardware processor and must use it in mutual exclusion. The ports ACQ and REL are added for interaction with the processor scheduler. The port ACQ is used for acquiring and REL is for releasing the processor. A process acquires the processor at the start of its behavior. It releases the processor on its termination. Let us mention that, the transformed model is a correct implementation of the initial model constructed from the application software. That is, it can be formally proven that the input/output behavior of every process is fully preserved by the transformation above.

2) Allocation according to mapping: Given an Application-Software and a Hardware-Architecture, a mapping Map associates software processes to hardware processors and software channels to memories, formally:

Mapping ::= Mapping-Item +

Mapping-Item

::= SW-Process → HW-Processor |

SW-Channel → HW-Memory

A mapping must be consistent. That is, for every writeconnection from process P to channel C in the application software, if the mapping associates P on processor H and C on memory M, there must exist a write-path of the form H Bus 1 . . . Bus n M in the hardware architecture. Similarly, for every read-connection from channel C to process P , there must exist a read-path of the form M Bus ′ 1 . . . Bus ′ m H. Example 8: For our example, we consider the following consistent mapping: generator

→ ARM1 C1 → LM1 square → ARM1 C2 → SM consumer → ARM2
The construction of the system model is completed as follows. For every hardware processor, we consider the composition of all transformed software processes mapped on it, together with all the FIFO routines required to access the FIFO buffers. These components are connected as defined by the transformed software model. Additionally, the composition includes a HW-CPU-Scheduler component which ensures mutual exclusion for execution on the processor.

Example 9: The structure of the ARM1 processor is shown in figure 12. It contains the generator and square processes together with their associated FIFO routines respectively, the FIFO-Write for writing on C1, the FIFO-Read for reading from C1 and the FIFO-Write for writing on C2. Moreover, for every memory component, we consider the union of all the FIFO buffers mapped onto it according to the mapping. Let us remark that no scheduling is done here: all the operations requiring access to memory are controlled by the processor and the bus, the memories being simple passive components, with no behavior.

Finally, the direct connections between the FIFO routines and the FIFO buffers which exist in the transformed software model are replaced by connections over the associated hardware communication paths. For example, the request/acknowledge connectors between a FIFO routine and the FIFO buffer (FB) are replaced by (i) request/acknowledge connectors from the FIFO routine to the master interface of the first bus of the associated hardware path and (ii) request/acknowledge connectors from the slave interface of the last bus of the path to the FIFO buffer.

We assume a high cache hit rate for the local variables of the processes mapped on a processor, and hence we do not model explicitly the allocation of process data in the memory. The memory is used only to model inter process data communications through the software FIFOs.

The system model can be seen as a refined implementation of the transformed BIP model of the application software according to hardware constraints. In fact, direct communication between components within the application software model have been replaced by multi-hop communication using hardware communication paths, along different buses. Moreover, mutual exclusion constraints are enforced between components running on the same hardware processors. These transformations do not impact the input/output behavior of the application. This can be proved by establishing a trace equivalence between the input and the transformed model. Nevertheless, the transformations reveal all the non-functional constraints the hardware architecture put on the execution due to contention for bus and memory access, bus access and transfer latencies, contention for processor, etc. These constraints are mandatory for an accurate performance evaluation of the application mapped on the hardware architecture. Example 10: Figure 13 shows the complete system model obtained for the mapping of the software application given in figure 6 to the hardware architecture of example 4 according to the mapping from example 8.

IV. PERFORMANCE ESTIMATION ON SYSTEM MODEL

We provide an infrastructure for performance estimation of the system model based on native BIP simulation. The process is dynamic and based on fine granular analysis of code generated for the target platform, using weight table profiling, as shown in figure 1.

A. Instrumenting the System Model

The system model is instrumented with the profiling API, embedded in the behavior of the SW-Processes. Every block of code, except the read/write calls, is instrumented by inserting profiling function calls at its start and at its end. These calls invoke the profiler which provides accurate execution times.

The instrumented BIP system model is used as such by the BIP tool-chain for compilation and execution using BIP native simulator. On execution, the profiler is invoked , which dynamically estimates the computation time of the current block of code of the SW-Processes. The estimated execution time is recorded by dedicated observers for delay measurements.

The observers added in the system model are timed BIP components and monitor both the computation and the communication delays. The communication latencies of the buses and memories are also recorded by separate sets of observers, considering the conflicts arising in the use of the buses and the memories.

B. Weight Table Profiling

We use standard tools for cross-compilation and coverage profiling of the source code for SW-Processes, generated from the system model using the BIP tool-chain. The source code is cross-compiled to generate the object code (assembly) for the target processor. The source code is also instrumented for coverage analysis. The profiler is parameterized by a weight-table, which characterizes the time of executing each elementary instruction on the target HW-Processor. The object code, instrumented sources and weight-table are used by the profiler dynamically during the simulation to estimate the execution time of transitions within processes.

V. EXPERIMENTS

The method described in section III has been implemented in a tool2 . It consists of two parts, the frontend that transforms the input specification into a system model, and the backend for performance estimation on the system model. The frontend uses an open source C parser called codegen 3 to parse C files that describe the behavior of the DOL processes into an intermediate model. This, along with the description of the hardware architecture and mapping information (XML description) is transformed into the system model. The backend uses gcov as a profiling tool for code coverage, and arm-rtems-g++ cross compiler for assembly code generation for ARM processors. The weighttable conforms to the ARM7 data sheet 4 .

We experimented the method on two applications: MJPEG [START_REF] Huang | Coupling MPARM with DOL[END_REF] and MPEG-2 [START_REF] Thiele | Mapping applications to tiled multiprocessor embedded systems[END_REF], [START_REF] Huang | Coupling MPARM with DOL[END_REF] described in subsections V-A and V-B respectively. We used the multiprocessor ARM (MPARM 5 ) with five tiles as the target architecture (a two tile MPARM is illustrated in figure 7). For the hardware model in BIP, we assumed all the local memories as SRAM with an access time of 2 cycles. The shared memory is a DRAM with an access time of 6 cycles. All CPU frequencies are assumed to be 200MHz. Communication paths are defined between all five processors using shared and local memories.

A. MJPEG Decoder

The MJPEG decoder application software reads a sequence of MJPEG frames and displays the decompressed video frames. The process network of the application is illustrated in figure 14 We experimented with eight different mappings to analyze their effect on the total computation and communication time for decoding a frame. The process and the sw channel mappings are illustrated on table I.

For the mappings described above, a system model contains about 50 BIP atomic components and 220 BIP connectors, and consists of approximately 6K lines of BIP code, generating around 19.5K lines of C code for simulation. The total computation and communication delays for decoding a frame for different mappings are shown in figure 15. Mapping (1) produces the worst computation time as all processes are mapped to a single processor. Mapping (2) uses two processors, still the performance does not improve much. But (3) gives much better performance as the computation load is balanced. The other mappings can not produce better performance as the load can not be further distributed, even if more processors are used. The communication overhead is reduced if we map more channels to the local memories of the processors. The bus and memory access conflicts are shown in figure 15. As more channels are mapped to the local memory, the shared bus contention is reduced. However, this might increase the local memory contention, as shown for [START_REF] Haubelt | Systemcodesigner: automatic design space exploration and rapid prototyping from behavioral models[END_REF].

B. MPEG2 Decoder

The MPEG2 decoder application decodes a set of moving pictures and associated audio information. We used an application case study where there are seven processes Dis-patchGops (DG), DispatchMb (DM), DispatchBlocks (DB), TransformBlock (TB), CollectBlocks (CB), CollectMb (CM) and CollectGops (CG) and six software channels C1, . . . , C6. The process and the sw channel mappings are illustrated on table II.

For the MPEG-2 case study the BIP System Model contains about 90 BIP atomic components, 340 BIP connectors and 30K lines of BIP code generating approximately 100K 17. The MPEG-2 process network is characterized as computationally intensive. The more we distribute the computational load to different CPUs, the smaller is the computational delay. Since the SW-channels are few, there is small difference in the communication delays between the different mappings, except for mapping [START_REF] Thiele | Mapping applications to tiled multiprocessor embedded systems[END_REF] where all processes and SW-channels are mapped on a single tile. However, as we distribute the processes into more tiles, the communication delay increases and more bus conflicts occur. The best throughput is achieved in Mapping [START_REF] Grtker | System Design with SystemC[END_REF] due to the usage of five CPUs and their local memories.

VI. CONCLUSION

The presented method allows generation of a correct-byconstruction model of a mixed hardware/software system from application software, a description of the hardware architecture and a mapping. The method is completely automated and supported by BIP tools. The system model is obtained by refining the application software model and composing it with the hardware architecture model. The composition is defined by the mapping. BIP instruments the incremental construction of the models. Its expressiveness allows the integration of architecture constraints into the application model without suffering complexity explosion.

The method clearly separates software and hardware design issues. It is also parameterized by design choices related to resource management such as scheduling policies, memory size and execution times. This allows mastering the complexity and appreciation of the impact of each parameter on system behavior.

When the generated system model is adequately instrumented with execution times, it can be used for performance analysis and design space exploration. Experimental results show the feasibility of the system model for fine granular analysis of the effects of architecture and mapping constraints on the system behavior. The method is tractable and allows design space exploration to determine optimal solutions.

Future work includes extension to other programming models for the application software and richer hardware architecture models that includes DMA (Direct Memory Access) Controller, Bus Bridge and Network on Chip communication. Moreover, we plan to include statistical model checking on system models consisting of multiple applications running on complex multicore architectures for performance analysis, as in [START_REF] Basu | Statistical abstraction and model-checking of large heterogeneous systems[END_REF].

Figure 1 .

 1 Figure 1. System Model Construction and Performance Evaluation

  SW-Process + . SW-Channel + . SW-Conn + SW-Process ::= SW-InPort * . SW-OutPort * . SW-Behavior SW-Channel ::= SW-RecvPort . SW-SendPort . SW-Channel-Behav SW-Conn ::= SW-Read-Conn | SW-Write-Conn SW-Write-Conn ::= SW-OutPort . SW-RecvPort SW-Read-Conn ::= SW-SendPort . SW-InPort SW-Behavior ::= a-C-program SW-Channel-Behav ::= FIFO-Param + Each software process P has input ports P.InPort i , output ports P.OutPort j and behavior P.Behavior. Each channel C has a single input port C.RecvPort and a single output port C.SendPort.

Figure 2 .

 2 Figure 2. An application software

Figure 3 .

 3 Figure 3. C code fragment of the square process

Figure 4 .

 4 Figure 4. The model of the square process as an atomic BIP component Example 2: Figure 4 shows the translation of the square process into an atomic component in BIP. The generated BIP component has ports IN, OUT, control locations L1, . . . L5 and variables index, len and i. Additional variables size and address are associated as parameters of the ports. Transitions are labeled by IN, OUT and τ , denoting an internal transition. At L2, it awaits synchronization through IN corresponding to the read primitive call. At L4 it awaits synchronization through OUT corresponding the write primitive call. At L1, internal transitions with guard model the conditional (if) statement. Exit of the process on a detach is modeled by the deadlocked location L5.

Figure 5 .

 5 Figure 5. SW-channel (FIFO) in BIP

Figure 6 .

 6 Figure 6. Application software model in BIP Example 3: The figure 6 provides the complete BIP model obtained from the application example given in figure 2. It consists of the BIP component generator sending data to square and consumer by using channels C1 and C2 respectively.

Figure 7 .

 7 Figure 7. A multi-core hardware architecture with two ARM tiles

  HW-Processor + . HW-Memory + . HW-Bus + HW-Processor ::= SW-Process (t) + . HdS + . HW-Cpu-Scheduler . SW-Conn + HdS ::= FIFO-Read | FIFO-Write SW-Conn ::= SW-Process-HdS | SW-Process-HW-Cpu-Scheduler | HdS-HW-Cpu-Scheduler HW-Memory ::= FIFO-Buffer + 1) Transformation of the BIP Application Model:

Figure 9 .

 9 Figure 9. Low-level implementation BIP model for software channels

Figure 10 .

 10 Figure 10. FIFO-Read component

Figure 11 .

 11 Figure 11. The transformed BIP model for the square process Example 7: The transformed behavior of the square process from figure 4 is provided in figure 11.Let us mention that, the transformed model is a correct implementation of the initial model constructed from the application software. That is, it can be formally proven that the input/output behavior of every process is fully preserved by the transformation above.2) Allocation according to mapping: Given an Application-Software and a Hardware-Architecture, a mapping Map associates software processes

Figure 12 .

 12 Figure 12. The BIP Model of the HW Processor ARM1

Figure 13 .

 13 Figure 13. The BIP system model of generator-square-consumer application software mapped into 2-tile ARM hardware architecture

Figure 15 .

 15 Figure 15. Mjpeg Performance Analysis Results

Figure 17 .

 17 Figure 17. Mpeg-2 Performance Analysis Results

  Every software channel is translated into a predefined BIP atomic component, as shown in figure5. It has ports recvPort and sendPort, and a single control location L1. It contains an array of data buff parametrized by size N . The variable x associated with recvPort gets the received value which is inserted into buff. The variable y associated with sendPort contains the value to be read next. The FIFO policy is implemented by using two indices i and j, for respectively insertion/deletion into/from the (circular) buffer buff.

	2) Translation of Software Channels into BIP:
	recvPort	i=0; j=0; count=0;	sendPort
	x			y
	[count<N] recvPort	L1	[count>0] sendPort
	buff[i]=x; count++; i=(i+1)%N	y=buff[j]; count--; j=(j+1)%N
	y=buff[j];	var: x, y, i, j, count, buff[N]
	The generated
	BIP component has ports IN, OUT, control locations L1,
	. . . L5 and variables index, len and i. Additional variables
	size and address are associated as parameters of the ports.
	Transitions are labeled by IN, OUT and τ , denoting an
	internal transition. At L2, it awaits synchronization through
	IN corresponding to the read primitive call. At L4 it awaits
	synchronization through OUT corresponding the write prim-
	itive call. At L1, internal transitions with guard model the
	conditional (if) statement. Exit of the process on a detach is
	modeled by the deadlocked location L5.

  . It contains five processes SplitStream (SS), SplitFrame (SF), IqzigzagIDCT (IDCT), MergeFrame (MF) and MergeStream (MS), and nine communication sw channels C1, . . . , C9.

	ARM1	ARM2	ARM3			ARM4	ARM5
	SplitStream	SplitFrame	IqzigzagIDCT	MergeFrame	MergeStream
		C1		C3	C7		C8
		C2		C4			C9
		Shared		C5
				C6
		Figure 14. MJPEG Decoder application and a mapping
		ARM 1		ARM 2	ARM 3		ARM 4	ARM 5
	1	all
	2	SS, SF , IQ	M F , M S
	3	SS, SF		IQ, M F , M S
	4	SS, SF		IQ	M F , M S
	5	SS, M S		SF	IQ		M F
	6	SS		SF	IQ		M F	M S
	7	SS, SF		IQ	M F , M S
	8	SS		SF	IQ		M F	M S
		Shared		LM 1	LM 2		LM 3	LM 4
	1			all
	2	C6, C7		C1, C2, C3, C4, C5	C8, C9
	3	C3, C4, C5, C6	C1, C2	C7, C8, C9
	4	C3, C4, C5, C6, C7	C1, C2			C8, C9
	5	all
	6	all
	7	C6, C7		C1, C2, C3, C4, C5			C8, C9
	8			C1, C2	C3, C4, C5, C6	C7	C8, C9

Table I MAPPING

 I DESCRIPTION OF THE PROCESSES AND THE SW CHANNELS

Table II MAPPING

 II DESCRIPTION OF THE PROCESSES AND THE SW CHANNELS lines of C code. The total computation and communication delays for decoding 5 frames for different mappings are shown in figure

http://www-verimag.imag.fr/BIP-System-Designer.html

http://think.ow2.org

http://www.datasheetarchive.com/ARM7-datasheet.html

http://www-micrel.deis.unibo.it/sitonew/research/mparm.html