Pascal Giorgi

Romain Lebreton

Online order basis algorithm and its impact on the block Wiedemann algorithm *

Order bases are a fundamental tool for linear algebra with polynomial coefficients. In particular, block Wiedemann methods are nowadays able to tackle large sparse matrix problems because they benefit from fast order basis algorithms. However, such fast algorithms suffer from two practical drawbacks: they are not designed for early termination and often require more knowledge on the input than necessary. In this paper, we propose an online algorithm for order basis which allows for both early termination and minimal input requirement while keeping quasi-optimal complexity in the order. Using this algorithm inside block Wiedemann methods leads to an improvement of their practical performance by a constant factor.

Introduction

Order bases (also called sigma bases) are the cornerstone of efficient algorithms with polynomial matrices over finite fields. Their use in algorithms for the determinant, column reduction [START_REF] Giorgi | On the complexity of polynomial matrix computations[END_REF] or minimal nullspace basis [START_REF] Zhou | Computing minimal nullspace bases[END_REF] has reduced these problems to the simpler computation of polynomial matrix multiplication. A complete panel of these reductions can be found in [START_REF] Jeannerod | Asymptotically fast polynomial matrix algorithms for multivariable systems[END_REF][START_REF] Zhou | Fast Order Basis and Kernel Basis Computation and Related Problems[END_REF].

Since their introduction in 1994 [START_REF] Beckermann | A uniform approach for the fast computation of matrix-type Padé approximants[END_REF], order basis computations have a quasioptimal algebraic cost in the order of approximation. Fast linear algebra has been included later on for square matrices [START_REF] Giorgi | On the complexity of polynomial matrix computations[END_REF], then for any matrix dimensions [START_REF] Zhou | Efficient algorithms for order basis computation[END_REF]. Historically, order basis were introduced to generalize Hermite-Padé approximants to the matrix case. It is nowadays well known that a minimal generator of a linearly generated scalar sequence is related to a specific Padé approximant. As a non trivial transposition to the matrix case, it has been shown in [START_REF] Turner | Black Box Linear Algebra with the LinBox Library[END_REF][START_REF] Villard | A study of coppersmith's block Wiedemann algorithm using matrix polynomials[END_REF] that order bases provide a fast solution to the problem of computing a minimal generator of a linearly generated matrix sequence. Note that many other approaches solve this problem, e.g. Toeplitz/Hankel solver, matrix Berlekamp-Massey or matrix extended Euclidean, but only few of them provide fast polynomial arithmetic together with fast matrix multiplication. We refer to [START_REF] Kaltofen | On the matrix Berlekamp-Massey algorithm[END_REF] for a classification of all possible algorithms with their corresponding reference and complexity.

Modern sparse linear algebra over a finite field relies on Coppersmith's block Wiedemann method [START_REF] Coppersmith | Solving Homogeneous Linear Equation Over GF(2) via Block Wiedemann Algorithm[END_REF], which in turn computes a minimal matrix generating polynomial. This is the case for linear systems solving [START_REF] Coppersmith | Solving Homogeneous Linear Equation Over GF(2) via Block Wiedemann Algorithm[END_REF] or for rank computations [START_REF] Turner | A block Wiedemann rank algorithm[END_REF]. For generic matrices with low rank deficiency, as in integer factorization, the general a priori bound on the degree of this polynomial is tight. However, a significant rank deficiency makes this bound loose, which introduces a computational overhead. This situation has been observed in [START_REF] Dumas | Parallel computation of the rank of large sparse matrices from algebraic k-theory[END_REF] when computing the rank of matrices from Algebraic K-theory using the block Wiedemann method of [START_REF] Turner | A block Wiedemann rank algorithm[END_REF]. In such a case, early termination techniques [START_REF] Kaltofen | Early termination in sparse interpolation algorithms[END_REF] are an important tool to probabilistically detect the degree of the minimal matrix generating polynomial. However no fast block Wiedemann method enables early termination in the order basis computation to the best of our knowledge.

Hence we focus in this paper on the application of order basis to the block Wiedemann framework in the context of early termination. In particular, existing fast order basis algorithms such as Algorithm PM-Basis of [START_REF] Giorgi | On the complexity of polynomial matrix computations[END_REF] does not allow for early termination.

Main results and outline of the paper. We contribute by giving two new algorithms, called iPM-Basis and oPM-Basis, which both have a quasi-linear time complexity in the degree of the approximation's order. After some background on order basis in Sections 2 and 3, we introduce in Section 4 Algorithm iPM-Basis, which is an iterative variant of PM-Basis more suited to early termination. This enables us to incorporate in block Wiedemann both early termination and quasi-linear time algorithm for the computation of the minimal matrix generating polynomial.

However, this is not optimal yet because we may have to compute more coefficients of the input of the order basis than necessary. Besides, the cost of computing this input is dominant in block Wiedemann, and so it is of great interest to minimize the number of required coefficients. Since online algorithms require minimal knowledge on the input, we propose a fast online order basis oPM-Basis in Section 5.

Then we study our order basis algorithms inside block Wiedemann method in the context of early termination. We analyze the complexity of our approach in Section 6 and make explicit the potential gains in Proposition 7. Finally Section 7 gives implementations and benchmarks and shows the positive impact of oPM-Basis in practice (see Figure 4).

Background

Let K be a field.

If A = i A i x i is in K[x] m×n or K[[x]
] m×n , and i, j are integers with i j, then we write

A i...j = A i + A i+1 x + • • • + A j-1 x j-i-1 , so that A i...j has degree less than j -i. Let A ∈ K[x] m×m and B ∈ K[x] m×n such that deg(A) < d, deg(B) < h and d h. Then, the middle product MP(A, B, d, h) of A and B is defined as the part D d-1...h of the product D := A • B, so that deg(MP(A, B, d, h)) < h -d + 1.
The balanced case corresponds to h = 2d -1; we denote this MP(A, B, d).

Row reduced matrix

The notion of shifted degree of a polynomial matrix is an extension of the classic polynomial degree where some scalings are introduced either by row or by column. In particular, the shifted s-row degree of a row vector v ∈ K[x] 1×n with s = [s 1 , . . . , s n ] ∈ Z n corresponds to max

1 i n (s i + deg v[i]) ∈ Z.
This naturally extends to matrices A ∈ K[x] m×n by considering the vector in Z m of s-row degrees of the rows of A.

A matrix R ∈ K[x] m×n is said to be s-row reduced if the s-row degree of R is less than or equal to the s-row degree of U R for any unimodular matrix U . Comparisons of s-row degrees use the lexicographic order on their sorted vectors. The existence of a minimal s-row degree matrix in the set of U R for U unimodular is ensured by [START_REF] Zhou | Fast Order Basis and Kernel Basis Computation and Related Problems[END_REF]Lemma 2.11]. One can ensure the unicity of row reduction by using the well known Popov form [START_REF] Beckermann | Shifted normal forms of polynomial matrices[END_REF].

The following lemma gives a criterion to certify that a matrix is s-row reduced. Let us denote by x s the diagonal matrix Diag(x s1 , . . . , x sn ) ∈ K[x] n×n . Lemma 1 ( [START_REF] Beckermann | Shifted normal forms of polynomial matrices[END_REF][START_REF] Zhou | Fast Order Basis and Kernel Basis Computation and Related Problems[END_REF]). Let u be the s-row degree of a matrix R ∈ K[x] m×n and T ∈ K m×n such that x -u Rx s = T + O(x -1 ) x→∞ . Then R is s-row reduced if and only if T has full rank.

The following lemma states the compatibility of s-row degree and row reducedness with matrix multiplication.

Lemma 2 ([22, Lemmas 2.14, 2.18]). Let A, B belong to K[x] m×m with B invertible, u-row reduced and of u-row degree v. Then the u-row degree of AB is equal to the v-row degree of A and AB

∈ K[x] m×m is u-row reduced if and only if A ∈ K[x] m×m is v-row reduced.

Order basis algorithms

Let F = i 0 F i x i ∈ K[[x]
] m×n be a matrix of power series with m n, σ a positive integer and (F, σ) be the K[x]-module v ∈ K[x] 1×m s.t. vF ≡ 0 mod x σ . A polynomial matrix P is a (left) order basis of F of order σ ∈ N * and shift s ∈ Z m if the rows of P form a basis of (F, σ) and P is s-row reduced (see [START_REF] Zhou | Fast Order Basis and Kernel Basis Computation and Related Problems[END_REF][START_REF] Zhou | Efficient algorithms for order basis computation[END_REF] for more details). Without loss of generality, we only consider in this paper the case n = Θ(m) as in [START_REF] Giorgi | On the complexity of polynomial matrix computations[END_REF]. Indeed the techniques of [START_REF] Zhou | Efficient algorithms for order basis computation[END_REF] allow to reduce the general situation to our particular case.

Nowadays there are mainly two algorithms for computing order bases. In this section, we propose a new presentation of these two algorithms of [START_REF] Giorgi | On the complexity of polynomial matrix computations[END_REF] using the framework of order basis of [START_REF] Zhou | Fast Order Basis and Kernel Basis Computation and Related Problems[END_REF].

Throughout the paper, the order basis algorithms, except Basis, will take as input a power series matrix F ∈ K[[x]] m×n , a shift vector s ∈ Z m and an approximation order σ ∈ N * , and return a (σ, s) order basis P of F and its s-row degree u.

Base case (order 1)

In order to simplify the presentation of all following order basis algorithms, we treat separately the case σ = 1. This case corresponds to finding a polynomial matrix M of minimal s-row degree such that M F ≡ 0 mod x. This problem roughly reduces to the computation of a row echelon form of the kernel of πF ∈ K m×n , where π is a permutation matrix such that the sequence π s is increasing. Hence, one can use the PLE decomposition [START_REF] Jeannerod | Rank-profile revealing gaussian elimination and the CUP matrix decomposition[END_REF] over K for this matter. In this decomposition, P is a permutation matrix, L is a lower unit triangular matrix and E is the row echelon form. Note that among all possible PLE decompositions, we use one revealing the row rank profile.

Let ω be the exponent in the complexity of matrix multiplication over fields, which satisfies ω < 2.3727 [START_REF] Williams | Multiplying matrices faster than Coppersmith-Winograd[END_REF].

Algorithm 1: Basis(F 0 , s)

Input: F 0 ∈ K m×n , s ∈ Z m Output: M ∈ K[x] m×m , u ∈ Z m 1:
Find a permutation π s.t. π s is sorted increasingly

M = xI r -GL -1 r I m-r τ -1 π 5: u = τ -1 π s + 1 r 0 n-r T 6:
return M, u Proposition 1. Algorithm Basis outputs a (1, s) order basis M of F 0 and its s-row degree u in O(mnr ω-2 ) operations in K, where r corresponds to the rank of F 0 .

Proof. For the correctness of the algorithm, we need to prove that M F 0 ≡ 0 mod x, det M = 0 and M is s-row reduced. By definition of the row echelon form, it is clear that the m-r last rows of E are zero. Consequently the last m-r rows of M F 0 are also zero. Since the first r rows of M are [xI r 0], it is clear that M F 0 ≡ 0 mod x. Then det M = ±x r = 0. Finally we assume, without loss of generality, that the permutations π and τ are equal to the identity. This assumption implies that s is already sorted increasingly and that u = [s 1 +1, . . . , s r +1, s r+1 , . . . , s m ] is the s-row degree of M . Therefore, we have x -u M x s = I m + O(x -1 ) x→∞ , which implies by Lemma 1 that M is s-row reduced and then concludes the proof of correctness of algorithm Basis. The complexity of Basis algorithm is dominated by the cost O(mnr ω-2 ) of PLE decomposition [START_REF] Jeannerod | Rank-profile revealing gaussian elimination and the CUP matrix decomposition[END_REF].

Quadratic algorithm

To compute a (σ, s)-order basis, Algorithm M-Basis presented in [START_REF] Giorgi | On the complexity of polynomial matrix computations[END_REF] works iteratively by induction on the approximation order.

Proposition 2. Algorithm M-Basis is correct and runs in O(m ω σ 2 ) operations in K.

Proof. Let us prove the correctness of Algorithm M-Basis by induction on the order σ. The case σ = 1 is already proven as it corresponds to Algorithm Basis.

Algorithm 2: M-Basis(F, σ, s)

Input: F ∈ K[[x]] m×n , s ∈ Z m and σ ∈ N * Output: P ∈ K[x] m×m and u ∈ Z m 1: P 0 , u 0 = Basis(F mod x, s ) 2:
for k = 1 to σ -1 do 3:

F = (x -k P k-1 F ) mod x 4: M k , u k = Basis(F , u k-1 )
5:

P k = M k P k-1 6: return P σ-1 , u σ-1
Assume P k-1 is a (k, s)-order basis of F , and let us prove that

M k P k-1 is a (k + 1, s)-order basis. By definition, F ∈ K m×n is the first non-zero term of P k-1 F . Since M k F ≡ 0 mod x, P k = M k P k-1 is a basis for (F, k + 1)
. One needs to show that P k is s-row reduced to finish the proof. By definition of order basis, P k-1 is s-row reduced and its s-row degree is u k-1 . Similarly, M k is u k-1 -row reduced and its u k-1 -row degree is u k . So it follows directly from Lemma 2 that P k is s-row reduced and its s-row degree is u k .

The complexity of M-Basis is dominated by step 5 which computes at most k + 1 products of m × m matrices. Hence, the cost of each loop is bounded by O(m ω k). Incorporating the latter cost into the 'for' loop gives the announced complexity for M-Basis.

Quasi-linear algorithm

A divide-and-conquer approach has been presented in [START_REF] Giorgi | On the complexity of polynomial matrix computations[END_REF], namely the PM-Basis algorithm. The idea is similar to the ones used in [START_REF] Beckermann | A uniform approach for the fast computation of matrix-type Padé approximants[END_REF][START_REF] Thomé | Subquadratic computation of vector generating polynomials and improvement of the block Wiedemann algorithm[END_REF]: it reduces the computation of one order basis of order σ to two order bases of order about σ/2. Algorithm 3: PM-Basis(F, σ, s)

Input: F ∈ K[[x]] m×n , s ∈ Z m and σ ∈ N * Output: P ∈ K[x] m×m and u ∈ Z m 1: if σ = 1 then 2:
return Basis(F mod x, s )

3: else 4: P l , u l = PM-Basis(F, σ/2 , s ) 5: F = MP(P l , F, σ/2 + 1, σ) 6: P h , u h = PM-Basis(F , σ/2 , u l ) 7: return P h • P l , u h
As explained in [START_REF] Giorgi | On the complexity of polynomial matrix computations[END_REF], this strategy allows to reduce the arithmetic complexity to O(m ω M(σ) log(σ)), where M denotes the arithmetic complexity of polynomial multiplication. Note that the proof of correctness of PM-Basis is a direct implication of Lemma 2. Contrary to M-Basis, the recursive structure of PM-Basis makes it difficult to enable early termination.

Fast iterative order basis

In this section, we give a variant of PM-Basis more suited to early termination. We present an original iterative variant iPM-Basis of the recursive algorithm

PM-Basis. Let k = r i=1 2 ni be the binary decomposition of k with n 1 < • • • < n r . The valuation ν 2 (k) in 2 of positive k corresponds to n 1 . By convention, ν 2 (0) = +∞. Algorithm 4: iPM-Basis(F, σ, s) Input: F ∈ K[[x]] m×n , s ∈ Z m and σ ∈ N * Output: P ∈ K[x] m×m and u ∈ Z m 1: F (∞) = F 2: M 0 , u 0 = Basis(F mod x, s ) 3: for k = 1 to σ -1 do 4: v = ν 2 (k), v = ν 2 (k -2 v ) 5: M (v) = ((M k-1 • M (0) ) • M (1) ) • • • M (v-1)
6:

F (v) = MP(M (v) , F (v ) , 2 v + 1, 2 v+1 ) 7: M k , u k = Basis(F (v) mod x, u k-1 ) 8:
Let r i=1 2 ni be the binary decomposition of σ 9:

M (n1) = ((M σ-1 • M (0) ) • M (1) ) • • • M (n1-1) 10: return (M (n1) • M (n2) ) • • • M (nr) , u σ-1
It is easy to add a stopping criteria to iPM-Basis in the 'for' loop after line 7. When σ is a power of two, Algorithms PM-Basis and iPM-Basis perform the same ordered sequence of polynomial matrix multiplication and call to Basis. Note that these algorithms differ in terms of memory management. Besides, our algorithm should store all the M k in one variable M and all the u k in u but these indexed variables will be useful in the proof.

Before proving that iPM-Basis is the iterative counterpart of PM-Basis, we give some hindsight on iPM-Basis. The polynomial matrices M k (of degree at most 1) coincide with those of M-Basis. Our order basis output consists in their product M σ-1 • • • M 0 , which is computed using the divide-and-conquer scheme of PM-Basis. For this matter, we have the variables

M (v) which store at step k the product M r-1 • • • M r-2 v where r = k -(k mod 2 v
). We perform the same products than a binary multiplication tree on M 0 , . . . , M σ-1 would do when σ is a power of two. But our approach decomposes these computations in steps and therefore yields an iterative variant of the binary multiplication tree recursive algorithm.

Proposition

3. For any F, s and t ∈ N, Algorithms PM-Basis and iPM-Basis have the same output and perform the same computations on the input F, 2 t , s.

Proof. We proceed by induction on t. When t = 0, both algorithms perform only a call to Basis (F mod x, s ). Now recursively, suppose the result true for inputs F, 2 r , s with r < t and let us prove it for the input F, 2 t , s. The program PM t := PM-Basis(F, 2 t , s) is made of 4 instructions, and we claim that they are mapped to subparts of It t := iPM-Basis(F, 2 t , s) as follows:

1. line 4 of PM t maps to the code It t,beg of It t from the beginning until line 5 for k = 2 t-1 ; 2. line 5 of PM t corresponds to line 6, k = 2 t-1 of It t ; 3. lines 6 and 7 of PM t are mapped to It t,end defined by the code from line 7, k = 2 t-1 to the end of It t . Our proposition follows from this former correspondence, which we now prove. Note that since σ is a power of two, iPM-Basis outputs M (σ) .

Proof of 1. First, the call to PM t-1 on line 4 of PM t corresponds to It t-1 by the inductive hypothesis. The latter code and It t share the same code until line 7, k = 2 t-1 -1. Then It t-1 outputs the product M (t-1) of line 9, which corresponds to the product of It t , line 5 for k = 2 t-1 .

Proof of 2. Let us prove that the inputs of both middle products are the same. First, M (t-1) and P l are equal because they are the respective outputs of It t-1 and PM t-1 which coincide by the inductive hypothesis. Concerning the second input, ) has not been overwritten since its definition F (∞) = F at line 1. So the middle products coincide and their outputs F and F (t-1) correspond.

F (v ) = F (∞) is equal to F because F (∞
Proof of 3. By recurrence hypothesis, line 6 of PM t corresponds to It t-1 := iPM-Basis(F (t-1) , 2 t-1 , u 2 t-1 -1 ). Thus, we will prove the correspondence of It t-1 followed by line 7 of PM t , with It t,end . We split both of these programs in 4 parts: 1) we will see that the initialization

F (∞) = F (t-1) of line 1 of It t-1 is compensated by direct access to F (t-1) instead of F (∞) in It t,end ; 2) the initialization of 2 of It t-1 corresponds to line 7, k = 2 t-1 of It t-1 ;
3) we will show later that respective 'for' loops of the two programs correspond; 4) the products of line 9 of It t-1 and It t,end differ only by one term, which corresponds to line 7 of PM t .

Let us prove that both respective 'for' loops of It t-1 and It t,end read the same variables. Let 1 k 2 t-1 -1 be the indices of the first 'for' loop and

k = k + 2 t-1 the corresponding indices inside It t,end . Notice that v = ν 2 (k) = ν 2 (k ). Similarly v = ν 2 (k -2 v ) = ν 2 (k -2 v ) except if k = 2 r in which case ν 2 (k -2 v ) = ∞ and ν 2 (k -2 v ) = t -1.
Therefore the variables (M (i) ) 0 i v and F (v ) that both program read coincide. Note that the variables F (v ) also coincide when k = 2 r because they both point to F (t-1) .

Corollary 1. Algorithm iPM-Basis outputs a (σ, s) order basis P of F and its s-row degree u in O(m ω M(σ) log(σ)) operations in K.

Proof. Proposition 3 proves both the complexity and correctness statements of iPM-Basis, but only for orders σ that are powers of two. We deduce our complexity estimate for general σ from the case σ = 2 r . Concerning correctness, we will prove that we match the output of M-Basis. The variables M k and u k of iPM-Basis coincide with those of M-Basis because PM-Basis and iPM-Basis perform the same computations when σ = 2 r and PM-Basis performs the same calls to Basis as M-Basis.

It remains to prove that we output

M σ-1 • • • M 0 like M-Basis to conclude. Let σ = r i=1 2 ni be the binary decomposition of σ with ν 2 (σ) = n 1 < • • • < n r .
For 0 j r, we define ϕ(k, j) = r i=r-j+1 2 ni as the sum of the jth highest bits of k. Then at the end of the algorithm, the variables

M (ni) equal to M ϕ(k,r-i+1)-1 • • • M ϕ(k,r-i) when 1 i r and their product is M σ-1 • • • M 0 as claimed.

Online order basis algorithm

The (σ, s)-order basis M σ-1 • • • M 0 of F is a function of F mod x σ ; this can be seen easily on Algorithm M-Basis. Moreover, during the intermediate steps of M-Basis(F, σ, s), only F mod x k+1 is read to compute M k . Therefore M-Basis allows for savings both on the input and on the computations in case of early termination. This is not the case for the fast algorithms of the family of PM-Basis. The speed of these algorithms is attained by performing computations in advance, which in turn require more knowledge on F . Consider the computation PM-Basis(F, 2 r , s) and let us focus on its intermediate steps which are made explicit in iPM-Basis. Then for any t r and any intermediate step k such that

2 t-1 < k 2 t , Algorithm iPM-Basis reads F mod x 2 t to compute M k (see line 6 for k = 2 t-1 where F (v ) = F ).
In some applications, such as the block Wiedemann method, the cost of computing the input F is dominant (see Section 6) and minimizing its required precision can have a practical impact. The purpose of this section is to give an algorithm that requires minimal knowledge on F at each intermediate step while keeping a quasi-optimal complexity in the order σ. Note that Algorithm M-Basis satisfies the first condition but not the second, and that PM-Basis (and iPM-Basis) is in the opposite situation.

Online algorithms

Let A be an algorithm, i be one of its inputs and o its only output. Assume that both i = [i 0 , . . . , i n ] and o = [o 0 , . . . , o n ] are given as sequences. We say that A is an online algorithm with respect to its input i if it reads at most i 0 , . . . , i j when computing o j for any 0 j n. This notion, first defined by [START_REF] Hennie | On-line Turing machine computations[END_REF], was popularized in the domain of Computer Algebra in [START_REF] Hoeven | Relax, but don't be too lazy[END_REF] under the name of relaxed algorithms.

In our context of order basis algorithms, we consider the input F as the sequence of its power series coefficients [F i ] i<σ and the output to be the sequence [M i ] i<σ . To put it differently, we say that an order basis algorithm is online if it requires minimal knowledge on F , that is if it reads at most F 0 , . . . , F k when computing M k .

We have seen that M-Basis is an online order basis algorithm. Our objective in this section is to derive a fast online order basis algorithm from iPM-Basis. We have noticed that the middle product of step 6 of iPM-Basis may read more entries of F than necessary. Therefore, we need to perform this step differently.

Shifted online middle product

Throughout this section, we will denote by R a ring and M a left R-module. Let A ∈ R[x] and B ∈ M [x] with A of degree less than d. This theoretical framework is meant to allow for the multiplication of A by B in our special case R = K m×m and M = K m×n .

We are interested in computing the middle product C = MP(A, B, d, h) under the following constraint: we cannot use the polynomial coefficients A i+d or B i+d before we have computed the coefficients C 0 , . . . , C i of C. Compared to online algorithms, there is a shift in the indices of A and B due to the fact that the ith coefficient of MP(A, B, d, h) is the (i + d -1)-th coefficient of A • B. We will call such an algorithm a shifted online algorithm for the middle product.

Note that for our problem, the reading constraint affects only B since A i+d = 0 by definition. Therefore, we are only interested in algorithms for the middle product that are shifted online with respect to B. Algorithms that are (shifted) online with respect to only one input out of two are called (shifted) half-line algorithms.

Let us focus on the balanced case h = 2d -1 for the rest of this section, where deg(MP(A, B, d)) < d. As in the offline (not online) case, we have a naive algorithm for the middle product which consists in computing (A•B) mod x 2d-1 and removing the unnecessary lower terms. We denote by H(d) the number of R-module operations on M necessary to multiply two power series at precision d (i.e. modulo x d ) by a half-line algorithm. Similarly to offline multiplication, we suppose that H is super-linear, i.e. that rH(s) H(rs) for any r, s ∈ N. The cost of the naive shifted online middle product algorithm is therefore SMP Naive (d) := H(2d -1).

We now propose a dedicated shifted online algorithm sMiddleProduct for the middle product that will save asymptotically a factor 2 in the complexity compared to the naive algorithm.

Algorithm 5: sMiddleProduct(A, B, d)

Input:

A ∈ R[x] of degree < d, B ∈ M [x] Output: MP(A, B, d) ∈ M [x] 1: C = (A 0...d • B 0...d ) div x d-1 2: for k = 1 to d -1 do 3: m = ν 2 (k) 4: C = C + MP(A 0...2 m+1 -1 , B d+k-2 m ...d+k , 2 m ) x k 5:
return C

The scheme of computation of our algorithm is given in Figure 1 In Figure 1 and Algorithm sMiddleProduct (see lines 3-4), we use the half-line multiplication of [START_REF] Hoeven | Relaxed multiplication using the middle product[END_REF] due to its good performances.

It turns out that the middle product of line 6 of iPM-Basis is unbalanced and does not exactly fit the settings of this section. However, the same ideas apply; we can cut the middle product into two parts, and then apply an offline multiplication on the lowest part and an half-line algorithm on the highest part.

Online order basis

For the purpose of a fast online order basis algorithm, we will use the lazy data structure to encode our polynomial matrices as in [START_REF] Berthomieu | Relaxed algorithms for padic numbers[END_REF][START_REF] Hoeven | Relax, but don't be too lazy[END_REF]. In this data structure, polynomials are given as a promise, and the evaluation of its coefficients is delayed.

We assume that we dispose of a class O sMiddleProduct that performs the shifted online middle product using the lazy data structure. In this setting, step k of Algorithm sMiddleProduct is performed when the kth coefficient of the middle product is asked for. By computing the middle product in this fashion, we obtain an online order basis algorithm oPM-Basis.

Algorithm 6: oPM-Basis(F, σ, s)

Input: F ∈ K[[x]] m×n , s ∈ Z m and σ ∈ N * Output: P ∈ K[x] m×m and u ∈ Z m 1: F (∞) = F 2:
M, u 0 = Basis(F mod x, s )

3: for k = 1 to σ -1 do 4: v = ν 2 (k), v = ν 2 (k -2 v ) 5: M (v) = ((M • M (0) ) • M (1) ) • • • M (v-1)
6:

F (v) = O sMiddleProduct (M (v) , F (v ) , 2 v + 1, 2 v+1 ) 7: M, u k = Basis(F (v) mod x, u k-1 ) 8:
Let r i=1 2 ni be the binary decomposition of σ 9:

M (n1) = ((M • M (0) ) • M (1) ) • • • M (n1-1) 10: return (M (n1) • M (n2) ) • • • M (nr) , u σ-1 Proposition 4. Algorithm oPM-Basis is correct and online. It runs in O(m ω M(σ) log 2 (σ)) operations in K. Proof. Noticing that F (v) = (M k-1 • • • M 0 F ) k...k+2 v ,
we deduce that the call to (F (v) mod x) on line 7 will read at most the entries F 0 , . . . , F k of F at step k.

The cost analysis is similar to the one of PM-Basis, except that middle products now cost

H(d) = O(M(d) log(d)) instead of M(d) [10]. Thus, Algorithm oPM-Basis requires O(m ω H(σ) log(σ)) = O(m ω M(σ) log 2 (σ)) operations in K, which is quasi-linear in σ.

Application to block Wiedemann

Let A ∈ K N ×N be a sparse matrix, which means that it has O(N log(N ) C ) non-zero entries for a constant C > 0. Block Wiedemann methods are useful for handling such matrices, and especially for linear systems solving [START_REF] Coppersmith | Solving Homogeneous Linear Equation Over GF(2) via Block Wiedemann Algorithm[END_REF] or for computing their rank [START_REF] Turner | A block Wiedemann rank algorithm[END_REF]. The main ingredient of block Wiedemann algorithms is the computation of the left minimal matrix generating polynomial Π AV U ∈ K[x] m×m of the matrix sequence S = (S i ) i∈N where S i = (U A i V ) ∈ K m×n and m n. Here U ∈ K m×N and V ∈ K N ×n are chosen randomly. Note that A can be preconditioned to ensure structural properties (e.g. generic rank profile, no nilpotent blocks of size greater than one in its Jordan form) [START_REF] Chen | Efficient matrix preconditioners for black box linear algebra[END_REF] and then allow a good probability of success for block Wiedemann method.

Definition 1. The matrix Π = d i=0 Π i x i ∈ K[x] m×m is a left matrix gener- ating polynomial of S if d i=0 Π i S i+j = 0 m×n , ∀j 0,
and det(Π) = 0. It is said to be minimal if it is row reduced and deg det(Π) is minimal. Furthermore, it is unique if it is in Popov form [START_REF] Kaltofen | On the matrix Berlekamp-Massey algorithm[END_REF].

As described in [START_REF] Turner | Black Box Linear Algebra with the LinBox Library[END_REF][START_REF] Villard | A study of coppersmith's block Wiedemann algorithm using matrix polynomials[END_REF], one can obtain Π AV U from a (σ, s) order basis P of the matrix series F = [

σ i=0 S T i x i | I m ] T ∈ K[[x]] 2m×n where s = [0 m , 1 m ].
It is well known since Coppersmith's work [START_REF] Coppersmith | Solving Homogeneous Linear Equation Over GF(2) via Block Wiedemann Algorithm[END_REF] that an order σ C = N/m + N/n + O(1) is sufficient in most cases. However, this bound may be loose for rank deficient matrices, as experimented in [START_REF] Dumas | Parallel computation of the rank of large sparse matrices from algebraic k-theory[END_REF].

Theorem 2.12 of [START_REF] Kaltofen | On the complexity of computing determinants[END_REF] gives a tighter bound by looking at the invariant factors of the characteristic matrix xI N -A. Let µ be the sum of the degrees of the n largest invariant factors of xI N -A. Then a minimal matrix generating polynomial of S can be derived from an order basis at precision σ KV = µ/m + µ/n +O(1) [START_REF] Kaltofen | On the complexity of computing determinants[END_REF][START_REF] Turner | A block Wiedemann rank algorithm[END_REF]. Note that µ gives a bound on the degree of the determinant of Π AV U .

When the values of µ and σ KV are not known, and assuming µ N , early termination is a good strategy to stop the course of the order basis algorithm. We use a heuristic method to detect that we have reached order σ KV . For example, we can test if the m smallest values of the s-row degree of P remain unchanged for a few consecutive orders. This condition is similar to getting consecutive zero discrepancies in the Berlekamp-Massey algorithm [START_REF] Kaltofen | Early termination in sparse interpolation algorithms[END_REF].

We will from now on set ourselves in the context of early termination in the block Wiedemann algorithm. We have three candidate algorithms for order basis in block Wiedemann that allow for early termination : the two online algorithms M-Basis and oPM-Basis, and the offline algorithm iPM-Basis. In the following, we only consider the sequential case as the analysis is equivalent in the parallel setting. Indeed, using t processors with t divides m will divide the cost for computing the sequence elements by t. The gain is similar for parallel computation of order basis.

We demonstrate below that the online algorithms M-Basis and oPM-Basis can improve offline block Wiedemann complexity by a constant factor, up to 2, under some conditions. Whereas this gain can only be observed with M-Basis when the ratio µ/N is small, we show that this is always the case with oPM-Basis whenever N is big. For the sake of simplicity we now assume that m = n. Computing S has dominant cost. Let δ be the minimal value such that a (δ, s)-order basis of F contains Π AV U . We know that 2µ/m < δ σ KV σ C [START_REF] Kaltofen | On the complexity of computing determinants[END_REF]. We will assume that our heuristic stops the order basis algorithm at order δ.

We now compare the complexity of the first two steps of the block Wiedemann algorithm : 1) computing the sequence S and 2) computing the (δ, s)-order basis, depending on the choice of order basis algorithm (M-Basis, iPM-Basis or oPM-Basis). We denote these complexities respectively by C M-Basis , C iPM-Basis and C oPM-Basis .

M-Basis and oPM-Basis are both online algorithms that only require the first δ coefficients of F . On the contrary, iPM-Basis requires the first 2 t coefficients of F when 2 t-1 < δ 2 t . Let X be the cost of computing one coefficient of the sequence F , i.e.

U A i V = U (A(A i-1 V )), that is X = Θ(mN (m + log(N ) C )). Then C M-Basis = δX + O(m ω δ 2 ), C iPM-Basis = 2 t X + O(m ω δ log(δ)), C oPM-Basis = δX + O(m ω δ log 2 (δ)).
Proposition 5. The ratio C iPM-Basis /C M-Basis tends to 2 t /δ when µ/N tends to zero.

Proof. First, one has C iPM-Basis C M-Basis = 2 t /δ+O(ε1) 1+O(ε2) where ε 1 := m ω log(δ)/X and ε 2 := m ω δ/X . Since ε 1 = O(ε 2 ), it is enough to prove that ε 2 tends to 0 to imply that C iPM-Basis /C M-Basis tends to 2 t /δ when µ/N tends to 0. Using X = Ω(m 2 N ) and δ = O(µ/m), we get that ε 2 = O(m ω-3 µ/N ) which tends to 0 as claimed. Proposition 6. The ratio C iPM-Basis /C oPM-Basis tends to 2 t /δ when N tends to infinity.

Proof. The proof is similar but with ε 2 := m ω log 2 (δ)/X . From ε 2 = O(m ω-2 log 2 (δ)/N ) and using δ = O(N ) and m = O(log(N )), we deduce that ε 2 = O(log 3 (N )/N ) which tends to 0 when N tends to infinity.

Analysis of the gain due to online algorithms. We will now assume that either we are using M-Basis and µ/N is small, or we are using oPM-Basis and N is large. Therefore, Propositions 5 and 6 state that we gain a constant factor K δ = 2 t /δ compared to offline block Wiedemann. This constant depends on δ and satisfies 1 K δ < 2. The upper bound is tight as it is easy to find δ such that K δ is arbitrarily close to 2. The following proposition states that the average gain is larger than 1.

Online middle product

Our fast shifted online middle product relies on the half-line product of [START_REF] Hoeven | Relaxed multiplication using the middle product[END_REF]. Such a choice is driven by its complexity, roughly 1 4 M(σ) log(σ), and its simplicity of implementation requiring only an efficient offline middle product. Therefore, we readily benefit from the optimizations of our polynomial matrix multiplication implementation.

Figure 2 illustrates the ratio of timings of the offline and shifted online middle product algorithms with respect to polynomial matrix multiplication. As expected, the offline middle product behaves exactly as polynomial matrix multiplication and the shifted online ratio behaves like 1/4 • log(σ).

Online order basis

As described in Sections 3 to 5, fast order basis algorithms reduce to polynomial matrix multiplication. Our implementations of iPM-Basis and oPM-Basis follow this reduction and therefore mainly use the online middle products of Section 7.2 and the polynomial matrix multiplication of Section 7.1. In practice, for small values of σ, it appears not surprisingly that using a fast order basis algorithm is not relevant, and that switching to a quadratic one as M-Basis improves the performance. From our experiments, we observe that this is always the case for σ 32. In order to speed up order basis implementations, we therefore use an x 32 -adic version of our algorithms. For iPM-Basis this corresponds to stopping the recursion at σ = 32, while for oPM-Basis this comes down to incrementing k by 32 in the 'for' loop. For both of them, the calls to Basis are then replaced with calls to M-Basis on the input series mod x 32 . In Figure 3, we report the relative performance of online order basis, i.e. M-Basis and oPM-Basis, against the fast offline order basis iPM-Basis. From this figure, one can first see that a quadratic approach, such as M-Basis , is not competitive w.r.t. its fast counterpart as soon as σ > 128. The high efficiency of our polynomial matrix operations (multiplication and middle product) makes it possible to use the fast variant very early. Furthermore, the larger the matrix dimensions the earlier fast variants outperform M-Basis.

As shown in Figure 3, oPM-Basis almost keeps up with the performance of iPM-Basis. The ratio of timings stays below 2 in the given range of orders, despite its theoretical logarithmic behavior. This makes our oPM-Basis algorithm an interesting approach when online computation makes sense, for example when computing the matrix power series coefficients is very expensive.

Block Wiedemann algorithms

As demonstrated in Section 6, using early termination strategy in block Wiedemann should benefit in practice from online order basis calculation. In order to exhibit such a behavior in practice, we compare the timings of computing both the required sequence elements U A i V and the corresponding order basis for any approximation order up to a given bound. We characterize which algorithm is best suited depending on how early the algorithm terminates.

We use a random sparse matrix A ∈ M 2 17 (F 6946817 ) with 20 elements per row and two random dense matrices U, V T ∈ M 16×2 17 (F 6946817 ). To give a fair benchmark, we perform the calculation of the sequence elements U A i V in parallel (on the columns of V ) while order basis calculation is done sequentially. Order basis could also benefited from parallel computations but this only affects the result by a constant factor and does not change our observation.

Figure 4 illustrates the measured times for the first two steps of block Wiedemann with early termination using the offline iPM-Basis and online M-Basis, oPM-Basis algorithms. Note that σ = 2 14 corresponds to the classic bound 2N/m.

One sees on this figure that whenever early termination leads us to stop at a value δ smaller than the a priori bound σ C , online algorithms can improve performance. Of course, this is only the case when δ is not close to a power of two, where iPM-Basis still remains the fastest. In other cases, M-Basis improves the performance only for very small approximant orders, which is a consequence of Proposition 5. As soon as δ is getting larger, the quadratic complexity of M-Basis is not competitive with the quasi-linearity of iPM-Basis, even with its extra U A i V 's.

In the case of oPM-Basis, most of the possible values for δ allow to be faster than using the fast offline approach. Indeed, it is almost always true that the extra log(δ) of oPM-Basis is negligible compared to few extra U A i V 's needed by iPM-Basis. Except for very few cases around powers of two, this makes oPM-Basis always the fastest algorithm, as predicted by Proposition 6.

2 : 3 :

 23 Compute the PLE decomposition τ LE of πF 0 Let r = rank(E) and write L = L r G I m-r 4:

when d = 8 .

 8 Decompose C = MP(A, B, d) in two parts C = L + xU where L = (A 0...d • B 0...d ) div x d-1 and U = (A 0...d-1 • B d...2d-1 ) mod x d-1 . At the first step of our algorithm, we are allowed to use A 0...d and B 0...d in order to compute C 0 = L 0 . We chose to anticipate some computations and compute the whole L via a classic offline multiplication algorithm at cost M(d). Let B := B d...2d-1 . For the computation of U , the reading constraint translates into the following: the computation of U i must read at most the coefficients B i of B . So we can use a half-line multiplication algorithm for this task, at cost H(d -1). Finally the addition C = L + xU costs O(d). Altogether, Algorithm requires SMP GL (d) R-module operations on M where SMP GL (d) := M(d) + H(d -1) + O(d) . Since d M(d) H(d) and 2H(d -1) H(2d -1), we deduce SMP GL (d) (1/2 + o(1)) SMP Naive (d).

Figure 1 :

 1 Figure 1: Shifted online middle product MP(A, B, 8)

Figure 2 :

 2 Figure 2: Relative performance of two middle products against polynomial matrix multiplication over Fp[x] 16×16 for p a 23-bit prime.

Figure 3 :

 3 Figure 3: Relative performance of two online order basis algorithms against polynomial matrix multiplication. The input matrix power series are in Fp[x] 32×16 for p a 23-bit prime.

Figure 4 :

 4 Figure 4: Timings of block Wiedemann with early termination for different order basis algorithms. A ∈ M 2 17 (Fp) with 20 elts/row, m = 16 and p a 23-bit prime.

Table 1 :

 1 Computation times of polynomial matrix multiplication in Fp[x] m×m .Here, p is a 23-bit FFT prime and matrices have degree d.

	m, log 2 (d) Mathemagix FLINT 2.4 our code
	16, 10	0.02s	0.26s	0.03s
	16, 11	0.04s	0.70s	0.06s
	16, 12	0.09s	1.68s	0.13s
	16, 13	0.20s	4.52s	0.28s
	128, 9	1.00s	26.21s	0.82s
	256, 8	4.00s	36.71s	1.75s
	512, 9	69.19s	465.66s	19.64s
	1024, 6	71.36s	115.52s	13.95s
	2048, 5	298.27s	263.88s	48.90s

This work has been supported by the French ANR under the grants HPAC (ANR-11-BS02-013), CATREL (ANR-12-BS02-001).

Proposition 7. Let δ be a random variable uniformly distributed between 1 and σ C = 2 r . Then the expected value of K δ tends to 2 ln 2 1.39 when r tends to infinity.

Proof. If δ were a random variable uniformly distributed between 2 t-1 + 1 and 2 t , the mean value of K δ would be

The actual mean value of K δ is

Practical performance

In this section, we discuss the implementation of our algorithms. The code is developed and distributed in the LinBox library (www.linalg.org). All the following benchmarks have been done on an Intel Core i7-4960HQ @ 2.6GHZ with 16Gb RAM and 4 cores. Only the computations of the block Wiedemann sequences have been done in parallel as it is the primary motivation for this approach.

Polynomial matrix multiplication

Our fast implementation of polynomial matrix multiplication relies on several optimizations. First, we use the FFLAS package, available through the LinBox library, as a primitive for matrix multiplication. This allows to reach the processor's peak performance for matrix multiplication over word size prime fields [START_REF] Dumas | Dense Linear Algebra over Word-Size Prime Fields: The FFLAS and FFPACK Packages[END_REF]. Second, we provide cache friendly DFT transforms for polynomials, which use SIMD vectorization, i.e. SSE4 instructions. As both matrices and polynomials must be stored contiguously to be cache friendly, we therefore implement two data structures, namely matrix of polynomials as well as polynomial of matrices, and efficient conversions between them. Table 1 gives a quick preview of the performance of our code. For comparison purposes, we provide the times of two standard libraries: FLINT (www. flintlib.org) and Mathemagix (www.mathemagix.org). Our code provides similar or better performance than these reference libraries. In particular, our optimizations on DFT are close to the ones of Mathemagix, which are known to be efficient. On the other side, for large matrices, our code performs even better thanks to the use of optimized BLAS libraries, which rely on the recent AVX2 instructions set of our processor. Note that Mathemagix does not yet provide such a support in their matrix_naive type. This emphasizes the benefit of the FFLAS approach, which always makes use of the most recent technology. Regarding FLINT's performance, we explain the gap in performance by our use of both cache friendly storage and vectorization.


A Ferdjoukh

A Baert

E Bourreau

A Chateau

R Coletta

C Nebut

Instantiation of meta-models constrained with OCL: A CSP approach

Keywords: Automated Model Generation, Constraint Satisfaction Problem (CSP), Object Constraint Language (OCL)

The automated generation of models that conform to a given meta-model is an important challenge in Model Driven Engineering, as well for model transformation testing, as for designing and exploring new metamodels. Amongst the main issues, we are mainly concerned by scalability, flexibility and a reasonable computing time. This paper presents an approach for model generation, which relies on Constraint Programming. After the translation of a meta-model into a CSP, our software generates models that conform to this metamodel, using a Constraint Solver. Our model also includes the most frequent types of OCL constraints. Since we are concerned by the relevance of the produced models, we describe a first attempt to improve them. We outperform the existing approaches from the mentioned point of view, and propose a configurable, easy-to-use and free-access tool, together with an on-line demonstrator.

Introduction

Model transformations are crucial in model driven development, and they must be accurately and carefully tested. Testing a model transformation is very similar to testing a classical program, except that the test data are more complex, since they are models. One issue in testing model transformations is then the automatic generation of test data, i.e. models. This task is complex first because the data to generate are themselves complex, and second because the models are constrained by constraints on the structure of the meta-model, and possibly by additional constraints coming from the transformation to test. A model generation mechanism is also needed to validate a metamodel [START_REF] Baudry | Barriers to Systematic Model Transformation Testing[END_REF]. Indeed, a meta-model is supposed to capture the essence of a domain, and the objective of its design is to define models. Therefore, it is of prime importance that the meta-model is adequately designed.

In our opinion, the main properties that are required for a proper model generation are: (i) scalability, which implies that the generation of many large models should take a small computational time, (ii) validity, meaning that the OCL constraints [START_REF] Omg | Object Constraint Language Specification, Version 2.4. Official Specification[END_REF] have to be taken into account, in order to generate only valid models, (iii) flexibility, implying that it can easily be parametrized, depending on the purpose of the automated generation, (iv) diversity of the solutions: ideally, models should be uniformly distributed on the solution space defined by the meta-model and the constraints. Model generation have been studied with several conceptual tools, but none of them fulfils the previously mentioned main properties. The contribution of this paper is a model generation mechanism based on CSP, that is closer to fitful all these properties. We propose a modelling of a meta-model as a CSP, which is as simple as possible, to minimize the number of variables and constraints in the CSP, and that thus provides a very quick generation process. This modelling allows OCL constraints to be introduced, thus produces only valid instances, i.e. models that conform to the input metamodels and that respect the OCL constraints of these meta-models. The approach is flexible, parametrizable, and allows additional constraints added by the user. Experiments of model generation from several meta-models taken from the literature validate our modelling. We also experiment the assistance to a meta-model designer in a real use case.

The rest of the paper is organized as follows. In section 2, we detail the main principles of Constraint Programming. In section 3, we analyse the existing work on model generation to identify the improvable aspects and extract the main issues. In section 4, we present our original CSP modelling of a meta-model. To complete this modelling, we detail the treatment of OCL constraints in section 5. After the presentation of an illustrative example in section 6, we detail our experiments performed on a benchmark of several meta-models, in section 7.

Background: CSP

This section is devoted to a very short introduction to the concepts and vocabulary of Constraint Satisfaction Problems (CSP). For a more complete view of this huge area, please refer to (Rossi et al., 2006). CSP are widely used to model many Artificial intelligence issues and combinatorial problems. Many real life and artificial intelligence problems can be modelized in CSP, for example: scheduling a sport competition, solving a board game like a sudoku or planning flights in an airport. A CSP is defined by Mackworth [START_REF] Mackworth | Consistency in Networks of Relations[END_REF] as follows: " We are given a set of variables, a domain of possible values for each variable, and a conjunction of constraints. Each constraint is a relation defined over a subset of the variables, limiting the combination of values that the variables in this subset can take. The goal is to find a consistent assignment of values from the domains to the variables so that all the constraints are satisfied simultaneously.". More formally, variables take their values in a finite domain. These domains are mapped on the set Z of integers.

Definition 1. A Constraint Network is composed of • a set of variables X = {x 1 , x 2 , . . . , x n }, • a domain on X , that is, a set D = {D(x 1 ), D(x 2 ), . . . , D(x n )}, where D(x i ) ⊂ Z
is a finite set of values that variable x i can take (its domain),

• and a set of constraints C = {C 1 ,C 2 , . . . ,C e }.

Constraints specify combinations of values that given subsets of variables must take or not take.

Definition 2. A constraint C i ∈ C is a boolean function involving a sequence of variables X (C i ) = {x i 1 , x i 2 , . . . , x i m } called its scope. The function is de- fined on Z m . A combination of values (or tuple) t ∈ Z m satisfies C i iff C i (t) = 1. Definition 3. Given a constraint network (X , D,C ), • an instantiation I on a set Y = {x 1 , . . . x k } ⊆ X is an assignment of values v 1 , . . . , v k to variables x 1 , . . . x k with v i ∈ D(x i ). • an instantiation I on Y violates C i iff X(C i ) ⊆ Y and I[X(C i )] ∈ C i .
• a solution to (X , D,C ) is an instantiation I on X which does not violate any constraint.

Example 1. Let us define a CSP (X ,D,C ), where :

X = {x 1 , x 2 , x 3 }, D = {{1, 2}, {1, 3}, {1, 2, 3}}.
We specify three constraints : C 1 :

x 1 = x 2 , C 2 : x 1 = x 3 , C 3 : x 2 = x 3 . The instantiation (x 1 , x 2 , x 3 ) = (2, 1, 3
) is a solution of (X ,D,C ).

Global Constraints

Global Constraints are one of the most important features in Constraint Satisfaction Problems. They capture a relation between a non-fixed number of variables. A Global Constraint provides shorthand for frequently recurring patterns and facilitates the work of the solver. Indeed, in a tree-search process, it allows dead-ends paths to be quickly detected and then speeds up the resolution.

allDifferent

The constraint allDifferent(Vars) is a global constraint that enforces all the variables of the set Vars to be different. This constraint is very useful, it occurs in many practical problems. Finding a solution for a CSP is performed by a program called a CSP solver. There are many solvers in the literature, such as Abscon, Choco and ECL i PS e solver. They essentially differ by their implementation language (Prolog for ECL i PS e , Java for the two others), the expressivity of the supported constraints and the complexity of the propagation algorithms implemented in the solver. To perform our experiments, we choose a solver allowing the use of a generic language for defining CSP such as XCSP [START_REF] Lecoutre | XML Representation of Constraint Networks: Format XCSP 2.1[END_REF], thus, we use Abscon solver [START_REF] Merchez | AbsCon: A prototype to solve CSPs with abstraction[END_REF].

Related work

Several approaches can be found in the literature for model generation, relying on several paradigms: Random graph, Graph grammars, Alloy constraints, Satisfiability Modulo Theory (SMT) formulas and CSP (Constraint Satisfaction Problems). In this section, we analyze those approaches considering the four quality criteria given in section 1.

In [START_REF] Mougenot | Uniform Random Generation of Huge Metamodel Instances[END_REF], Mougenot et al. use random graphs to generate models. First, they compute a spanning tree over the meta-model. Then, it is translated into a context-free grammar which is used in the generation process (ad-hoc algorithms may be used). This approach has a linear complexity, which makes it scaling very well. On the other side, it generates skeletons of models (Containment relationships only) and OCL constraints of meta-models can not be treated in such an approach.

An approach based on graph grammars is described by [START_REF] Ehrig | Generating Instance Models from Meta models[END_REF]. The idea of the approach is to translate all the components of a meta-model into grammar rules. Using these rules, the method generates nodes (representing class instances) or sub-graphs containing nodes related between them (representing references connecting two or more instances), leading to the construction of a graph representing the model. These rules are applied successively and randomly until the desired model is obtained. Unfortunately, this approach does not take into account OCL constraints. In addition, it is not flexible enough, since the parameterization of the process to add constraints is impossible.

In [START_REF] Sen | Automatic Model Generation Strategies for Model Transformation Testing[END_REF], Sen et al. present an approach based on Alloy constraints. The meta-model and its OCL constraints are translated into alloy language, then the Alloy solver generates solutions. The translating operation from a meta-model to Alloy constraints is quite simple. However, the approach suffers from a lack of scalability, since Alloy is dedicated to model checking and does not allow solutions to be quickly generated. In addition, the OCL constraints of a given meta-model are manually translated. This is not suitable when there are a lot of constraints to treat.

Cabot et al. use an approach based on CSP in [START_REF] González Pérez | EMFtoCSP: A Tool for the Lightweight Verification of EMF Models[END_REF], [START_REF] Cabot | Verification of UML/OCL Class Diagrams using Constraint Programming[END_REF]. It consists in the modelling of a meta-model and its OCL constraints into a CSP instance. Then, a CSP solver (namely ECL i PS e ) generates solutions. The CSP paradigm is very flexible and allows a large part of OCL constraints as well as additional constraints to be processed to improve the relevance of the generated models. Unfortunately, the proposed work suffers for a lack of performance. We think that it is due to a non-optimal modelling of the meta-model into CSP, for instance without using global constraints, and the naive transformation of OCL constraints. Indeed, all OCL constraints are translated in the same manner. In Constraints Programming, improve modelling can have a large influence on the performance. Using the appropriate constraint, such as, global constraint, reduce the number of variables are ideas to get this improvement. Wu et al. present in (Wu et al., 2013), an approach based on Satisfiability Modulo Theory. An SMT is a constraint problem which includes propositional satisfiability (SAT) and supports a richer language, for instance arithmetic operations, thus allowing more expressiveness for the constraints than with only SAT formulas. They translate a meta-model and its OCL constraints into SMT formulas. An SMT solver investigates for the satisfiability of these formulas. This approach automatically processes the OCL constraints of a meta-model. We were not able to evaluate the performances and the scalability of this solution because they generate only 1 or 2 instances per class in the paper.

In our previous work [START_REF] Ferdjoukh | A CSP Approach for Metamodel Instantiation[END_REF], we presented an approach to generate models from metamodels using CSP, translating a meta-model into a CSP instance and using the Abscon CSP solver to get solutions. In this paper, we present some central contributions to this preliminary work. The first one is to reduce the number of variables representing instances of meta-model classes by half compared to [START_REF] González Pérez | EMFtoCSP: A Tool for the Lightweight Verification of EMF Models[END_REF]. The second contribution concerns the treatment of references. The previous solution of Cabot et al. considered references as relationships between two different classes. We consider them as pointers. The third contribution exploits one of the most important features of CSP, the global constraints. The perspectives of [START_REF] Ferdjoukh | A CSP Approach for Metamodel Instantiation[END_REF] were mainly to improve the performances of our approach by optimizing the modelling into CSP, one step further, and to apply the tool on a wider range of meta-models. We were also concerned by the translation of the OCL constraints into CSP and the relevance of the generated models, together with the userfriendliness of the tool. For instance, it is important to visualize the produced solutions. Since the automated generation has to be integrated to more complex processes, it is very useful to allow an expertise on the solutions, to qualify them as relevant or not, or to score them according to some measure, in a given context. For the moment, automated expertise is not available, thus a human eye may be needed.

Meta-model to CSP

In this section we propose a modelling, as a CSP instance, of a meta-model that conform to Ecore. To figure out the inheritance between classes, we copy the references and the features from the super-classes to their sub-classes. We assume that each meta-model has a root class from which all the classes of the metamodel can be accessed through composition links. Definition 7. For each class c i , the interval D c i of instances is defined as follows: D c i = {M i-1 + 1, . . . , M i }, where M i = ∑ i j=1 size(c j ). Note 3. The root class of a meta-model has one and only one instance because of its particular role in the generated model.

Definition 4. A meta-model M is defined by {C l,F ,R }. C l = {c 1 , c 2 , . . . , c n } is a set of n classes, F = { f 1 , f 2 , . . . , f m } is a set of m features and R =
Example 3. To illustrate the modelling of the classes of a meta-model, we apply the process of the transformation on the meta-model of Figure 1

. C l = {root, A, B}, Minsize = 5, Maxsize = 10, size(root) = 1, size(A) = 6, size(B) = 8. D root = {1}, D A = {2, . . . , 7}, D B = {8, . . . , 15}. (1)
This means that there will be 6 different instances of the class A and 8 of the class B in the generated model. In addition, the set of values D B will be used to process the reference r of the meta-model.

Feature Modelling

To model a feature of simple type into CSP, we create a variable whose type is the domain. A variable is created for each instance of each class having a feature.

∀c ∈ C l, ∀ f ∈ F feature of c, ∀i ∈ {1, . . . , size(c)} : create a variable F c,i, f .
Notation 1. The data type of a feature f is denoted Type( f ).

The domain of the variables F c,i, f is given by: D(F c,i, f ) = Type( f ). The most frequently used data types are: Integer, Char and String.

Definition 8. The data type Enumeration is a set of values {v 1 , v 2 , . . . , v l } of the same type. A value v i is called a literal of the enumeration.

Enumerations are modelled as follows. Let f be a feature of type enum, where enum is an enumeration of l literals:

Type( f ) = enum ⇒ D(F c,i, f ) = {1, 2, . . . , l}.

Reference modelling

Considering references as pointers is quite simple and natural. While the modelling proposed by most of related works, see for example For each i ∈ {1, . . . , size(c)} and each reference r ∈ c.AllRe f erences, we create j variables Re f c,r i, j , where j ∈ {1, . . . ,U pperBound(r)}.

Note 4. The variable Re f c,r i, j is the value of j th instance of the reference r for i th instance of class c. Definition 10. The set of all types of the reference r, denoted S dst (r), is defined by : S dst (r) = r.EReferenceType∪ r.EReferenceType.getSubTypes(),

where: getSubTypes() designates the set of subtypes of a class and ERe f erenceType returns the reference destination class. It is needed to treat the references pointing abstract classes with a set of subclasses.

Definition 11. The upper bound of class domains, denoted by UDom, is defined as: Max c∈C l (Max(D c )).

We define the maximum number of optional reference instances by Max r∈R (U pperBound(r) -LowerBound(r)). We note it Max(R ).

To represent the non-allocated values (values of the variables which are not taken into account in the generated model), we define a set of joker values.

Definition 12. The set of joker values is defined by

jokers = {UDom + 1, . . . ,UDom + Max(R )}.
The domain of the variable Re f c,r i, j , denoted D(Re f c,r i, j ), is given as follows:

∪ c∈S dst (r) (D c ), if j < LowerBound(r), ∪ c∈S dst (r) (D c ) ∪ jokers, otherwise.
It means that the LowerBound(r) first variables must be allocated, thus there is no joker in their domain.

The other variables are optional, therefore their domain includes the set jokers.

Example 4. To model the reference r linking the two classes A and B in the meta-model of Figure 1 we create some variables. For each instance of class a, we create 3 variables {Re f A,r i,1 , Re f A,r i,2 , Re f A,r }. The domains of variables are D(Re f A,r i,1 ) = D(Re f A,r i,2 ) = D b = {8, . . . , 15} and D(Re f A,r i,3 ) = {8, . . . , 15} ∪ jokers, where jokers = {16, 17}. The third variable is optional, then its domain should include the set jokers for the case where this variable is not used. The three variables Re f A,r i, j , j ∈ {1, 2, 3} can be seen as pointers from instances of class A to instances of class B. An example of assignment of these variables is shown on Figure 2. When a variable Re f A,r i, j is equal to a value v, we consider that the instance of class A references the instance of class B number v.

OCL constraints' modelling

In this section we propose a modelling for a set of constructions of the OCL language. The modelled constructions are, in our opinion, some of the most common in the OCL constraints we met. They gather important operations such as navigation of references, loop operations, typing operations and constraints on the features of classes. The goal of this process is the automatic generation of models that respect the OCL constraints of a given meta-model. The modelling we propose is a kind-specific treatment of the OCL constraints. Indeed, each OCL construction has a specific CSP modelling. This allows us to reduce the number of CSP variables in some cases and to use global constraints in others. Our modelling is able to process a four valued OCL language since the CSP constraints we create guarantee the generation of valid models which consider the OCL constraints of a meta-model.

A constraint about one feature

This kind of OCL constraints applies a boolean expression on a class feature of the meta-model. These constraints have the following shape:

Context c inv: expr(t). c is a class of the meta-model and t a feature of c. To model such a constraint in CSP, we will apply, as a pre-treatment, the boolean expression expr(t) on the domain of t to obtain a new domain: D(F c,t ) = {e ∈ Type(t)|expr(e)}. Example 5. We consider the following OCL constraint on the feature marking of the class Place of the Petri nets meta-model in Figure 8:

Context Place inv: marking ≥ 0.
To model this constraint into CSP, we apply the boolean expression marking ≥ 0 on the domain D(marking) = [-20, 20] and we get the new domain: D new (marking) = [0, 20].

A constraint on more features

These constraints apply a boolean expression on two or more features of the same class or of different classes. It is not possible to process them by modifying the domain as in the previous case because each domain belongs to one feature and all the domains are independent. To model an OCL constraint of this kind, we proceed as follows:

Let Context c inv: expr(t 1 ,t 2 ) be an OCL constraint on two features of a class c.

First, we create a CSP predicate pred(expr) satisfying the boolean expression expr. Then, ∀i ∈ D c , ∀ the variables F c,i,t 1 and F c,i,t 1 , we create a constraint C, where: X (C) = {F c,i,t 1 , F c,i,t 2 } are its variables and pred(expr)) is its function.

Navigation of n references

The navigation of references is an important operation in the OCL language. It allows to jump from a class to another using a reference between these two classes.

Let Context C 0 inv: r 1 . . . . .r n .expr( f ) be an OCL constraint containing the navigation of n classes and a boolean expression on a feature of the last class as shown in Figure 3. To model this issue into CSP, we consider each reference r i as a pointer from the source class C i-1 to the target class C i . The expression expr( f ) is applied only on the instances of C n that are pointed by the instances of C n-1 and so on up to C 0 . We create the following CSP constraints:

∧ U pper(r i ) m i =1 (Re f c i-1 ,r i k i-1 ,m i = k i ) ∧ expr(F c n ,k n , f ), ∀k i ∈ D c i ,
where i ∈ {1, . . . , n}. Example 6. Let Context C 0 inv: r 1 .r 2 .t < 10 be an OCL constraint. A part of meta-model with the navigation of two references is given in Figure 4. In this example, we give the appropriate CSP constraints that we have to create in order to treat the previous OCL constraint:

The domains of the classes are:

D C 0 = {1, 2}, D C 1 = {3}, D C 2 = {4}.
We create the following constraints in CSP:

(Re f c 0 ,r 1 1,1 = 3) ∧ (Re f c 1 ,r 2 3,1 = 4) ∧ (F c 2 ,4,t < 10) (Re f c 0 ,r 1 2,1 = 3) ∧ (Re f c 1 ,r 2 3,1 = 4) ∧ (F c 2 ,4,t < 10)
Figure 4: A part of a meta-model with a navigation of two references.

Parallel branches' Navigation

We consider the OCL constraint Context C 0 inv: Expr(r 1 . . . . .r n . f , r 1 . . . . .r m . f ) containing the navigation of two parallel branches of n and m references and applying a boolean expression on two features f , f where f (resp. f ) belongs to the last class of the first (resp. second) branch of classes as shown in Figure 5.

To model this OCL constraint in CSP, we have to create the following constraints, ∀k i ∈ D c i of the first branch of classes and ∀k j ∈ D c j of the second branch of classes:

U pper(r i ) m i =1 (Re f c i-1 ,r i k i-1 ,m i = k i ) ∧ U pper(r j ) m j =1 (Re f c j-1 ,r j k j-1 ,m j = k j ) ∧ Expr(F c n ,k n , f , F c m ,k m , f ),
where i ∈ {1, . . . , n}, j ∈ {1, . . . , m}. Note 5. The modelling of Expr(F c n ,k n , f , F c m ,k m , f ) is the same as the treatment of an OCL constraint containing two features of different classes as explained in Section 5.2.

A loop operation: Forall

In the case of OCL constraints containing a loop operation, for example, forall, the idea is to identify the CSP variables on which we will apply the boolean expressions of these constraints. Then, we have to define the appropriate CSP constraints in order to process these boolean expressions.

We consider the part of the meta-model shown in Figure 6. It contains two classes C 0 , C 1 and one reference r 1 between them. An OCL constraint containing a forall loop operation is defined on this meta-model: The previous OCL constraint applies a boolean expression on the variables representing the feature f of the instances of class C 1 which are referenced by the instances of C 0 . The treatment of this constraint is divided into two operations: a navigation of one reference and the loop operation itself. The CSP constraints that should be created are the following:

Context C 0 inv: r 1 → f orall(expr(C 1 . f )).
∀i ∈ D C 0 , j ∈ [1,U pper(r 1 )], k ∈ D C 1 , (Re f C 0 ,r 1 i, j = k) ∧ Expr(F C 1 ,k, f ).
Our approach is able to process many other loops and operations over collections of variables, such as: exists, size, sum and count. It is possible to treat these operators because of the diversity of global constraints in CSP. For example, the global constraint sum(VARS,Res) computes the sum of the collection of variables VARS. Thus, we can match an operation in OCL language with a global constraint into CSP. Table 1 gives examples of these connections.

OCL Operation

Global Table 1: Some examples of correspondence between OCL and Global Constraints.

A typing operation: OCLType

The meta-model in Figure 7 contains n references r 1 , r 2 , . . . , r n from a class B to an abstract class C. This latter has n sub-classes A 1 , A 2 , . . . , A n . Context B inv: AllDiff(r 1 .oclType,. . . , r n .oclType) is an OCL constraint defined to make the types (A 1 , A 2 , . . . , A n ) of the references r 1 , r 2 , . . . , r n different.

To model such an OCL constraint, we do the following manipulations ∀ j ∈ [1,U pper(r)] and i ∈ {1, . . . , n}:

• Create n variables U i , where: U i ∈ D A i . • Create n variables V i , where: V i = i. • Create n Element constraints: Element(Re f B,C i, j , [U 1 , . . . ,U n ],V i ). • Create one allDifferent constraint: allDiffer- ent(V 1 , . . . ,V n ).
6 Illustrative example In this section, we use a meta-model describing Petri nets (Figure 8) to illustrate our approach. We explain step by step the modelling and the resolution process to get conform models, thus correct Petri nets. We developed a tool called Grimm for which a demonstration web version is available at: http: //info-demo.lirmm.fr/grimm/. To generate models conform to Petri nets, it proceeds as follows: For example:

• marking ≥ 0 ⇒ D(marking) = [0, 20]. • weight > 0 ⇒ D(weight) = [1, 20]. • ∀c 1 , c 2 instance of Place, c 1 .name = c 2 .name ⇒ allDifferent(name 1 , name 2 ).
6. Solve the resulting CSP instance with the abscon CSP solver and build the corresponding conform model (see Figure 9).

Experiments

In this section, we show the results of experiments we performed on a benchmark of meta-models. We focus on the following issues:

• Study the effect of the size and the structure of the meta-model on our generation process. Answer some questions: Are we able to generate models for larger meta-models with the same efficiency as for smaller ones ? What is the effect of the number of classes, or the number of references of the meta-model on the generation process ?

• Determine if our new modelling leads to a significant improvement of its execution time, compared to the one in [START_REF] Ferdjoukh | A CSP Approach for Metamodel Instantiation[END_REF]. In other words, does the removal of the variables representing instances of classes improve time performance ?

• Perform experiments on meta-models constrained with the recently supported OCL constraints. The goal is to compare the resolution time with and without the OCL constraints.

• Study how our tool can help designers during meta-model designing process. Indeed, during the conception of a meta-model, one needs to generate models to verify the correctness of the metamodel and the OCL constraints he designs. Is it possible to use our tool to achieve this goal ? Are results relevant and can they be applied in practical cases ?

In our experiments, we used 13 different metamodels extracted from various domains. These metamodels have different sizes (smaller and bigger ones). Some of them are from research and others are from industry. Table 2 gives their sizes (number of classes).

To perform the experiments we proceed on the metamodels as follows: First, vary the number of instances per class of the meta-model in [10,120] and deduce the total number of instances. The number of reference instances is equivalent to the number of class in-stances. Then, generate the CSP instance corresponding to the parameters in XCSP format [START_REF] Lecoutre | XML Representation of Constraint Networks: Format XCSP 2.1[END_REF]. Finally, solve it using the Abscon1 CSP solver.

Experiments on meta-models

We present the results of the experiments performed with our tool Grimm on a benchmark of metamodels2 . Table 2 shows the model generation time on 11 different meta-models by varying the number of instances contained in the resulting model. We observe that our approach is able to generate models that conform to meta-models of different sizes (number of classes and references) and domains. We notice that the generation is generally quicker for bigger metamodels. This is due to the number of classes in the meta-model. When this number is large, the number of instances for each one is smaller than in the case of meta-models with few classes. The other reason is related to CSP. The arity of our constraints (number of variables appearing in the constraint) grows more quickly in meta-models of small sizes. This slows down the resolution process.

In Table 3, we compare the efficiency of our approach with our previous work [START_REF] Ferdjoukh | A CSP Approach for Metamodel Instantiation[END_REF]. We observe an important improvement in performance. Resolution time gain is between 13% and 80% depending on the meta-model. These good results are the consequence of:

• Removal of an important number of variables, mainly those representing class instances. The number of variables is reduced by n rB×r in average, where n (resp. r) is the number of classes (resp. references) of the meta-model and rB is the upper bound of non-bounded references. For example, if rB = 10, the number of variables is reduced by 7.5% for the Petri net meta-model.

• Removal of a constraint of large arity on previously mentioned variables.

We compared our previous solution to Cabot et al. solution in [START_REF] Ferdjoukh | A CSP Approach for Metamodel Instantiation[END_REF]. We showed that our solution is more efficient than the solution using also CSP paradigm of [START_REF] Ferdjoukh | A CSP Approach for Metamodel Instantiation[END_REF]. This ratio represents the average resolution time with the previous version divided by the resolution time with the current version.

Meta-models constrained with OCL

In this section, we perform experiments on metamodels constrained with OCL constraints. We compare the resolution time for some meta-models in the two cases of taking into account the OCL constraints or not. Table 4 shows the results of these experiments.

We observe that the difference between the generation time with OCL constraints and the one without them is really small. In some cases, generation time is even smaller with OCL constraints. There are several reasons that may explain these results:

(1) The reduction of domains of some CSP variables reduces the search tree then the resolution time. (2) The use of global constraints (allDifferent and element) makes the solution efficient because these constraints propagate easily in the search tree.

(3) Our modelling of OCL does not create a lot of new variables, so the problem does not increase a lot.

Discussion: Assistant for designers

Here we describe how our tool can be used as an assistant for designers during meta-modelling process.

In order to validate a meta-model during its design, we proceed as follows:

(1) Generate a model that conform to the input meta-model using our tool. (2) If the expert detects an anomaly with the generated model, he has to correct the meta-model by: Correcting a reference cardinality, Adding classes or features or Adding OCL constraints. (3) Then, return to

Step (1) or finish. For more precision, we discuss the case of a meta-model corrected with the use of our tool and the eye of expert. The Feature meta-model is a meta-model designed by a PhD student. Using the Grimm model generation tool, this meta-model was corrected and published in a defended PhD thesis (AL-Msie'Deen, 2014). After an interview with the designer in order to understand the meta-model, we generated models that conform to the first version. The designer examined these models, detected an error, the correction leading to a new version of the meta-model. The process was repeated 4 times until getting a correct meta-model in the opinion of the user. The following types of errors were detected and corrected during this process: References with wrong cardinalities, Use of a class inheritance structure instead of enumeration features and The need of some OCL constraints to complete the meta-model.

Conclusion

Model generation is an important tool to test model transformations as well as to design meta-models. In this paper, we proposed a model generation mechanism, based on CSP. This paper develops three main contributions to automated model generation. The first contribution is to improve and refine our previous modelling, by notably reducing the amount of variables in the CSP. As shown by our experiments, it induces a very affordable computation time, allowing the perspective of a largescale generation. The second contribution addresses the treatment of OCL constraints. We do not handle the whole OCL, however we selected a subset of the OCL language and proposed efficient ways to transform them into additional constraints for the CSP modelling of the meta-model. The experiments show that adding those constraints do not significantly increase the resolution time, while it improves the relevance of the produced models. This is a particularly encouraging result, which illustrates the power of the CSP modeling in such a situation. Finally, we implemented a complete, easy-to-use and available tool, which can be used to generate models, and visualize them, from a meta-model provided in Ecore format. This tool allowed hidden constraints to be detected on the PetriNets meta-model, and helped to design a new meta-model for Features. It therefore proves its immediate usefulness and usability. A web demonstrator of our tool is available at: http: //info-demo.lirmm.fr/grimm/.

As perspectives, we are intending to continue the optimization of the solving time, especially in the generation of consecutive instances. During this sequence of generations we would like to introduce some additional constraints to separate the produced model by a given distance, in order to improve the diversity aspect of our tool. We may also cover a wider range of OCL constraints, and design an expert-driven add-on, to facilitate the adjustment of additional constraints leading to more realistic models.

Figure 1 :

 1 Figure 1: A meta-model containing three classes.

  [START_REF] González Pérez | EMFtoCSP: A Tool for the Lightweight Verification of EMF Models[END_REF], considers a reference instance as a pair of variables (a, b) where a class instance a references a class instance b, we propose to model a reference by only one variable associated to a class instance a. It will take the value assigned to the class instance b referenced by a. Let c be a class. The set of all references of c is c.AllRe f erences ⊂ R . Let r ∈ c.AllRe f erences be a reference of c. We denote by LowerBound(r) (resp. U pperBound(r)) the lower (resp. upper) bound of r. Definition 9. We define the reference bound of a meta-model M containing R a set of references, denoted Re f Bound, as the upper bound of unbounded references of M , in other words, the references which have * as upper bound.

Figure 2 :

 2 Figure 2: An example of value assignment for variables modelling the references. The instance of class A references three instances of class B. The dashed line indicates that this variable is optional why its domain includes jokers.

Figure 3 :

 3 Figure 3: A meta-model with a navigation of n References.

Figure 5 :

 5 Figure 5: A meta-model with a navigation of two branches of n and m references.

Figure 6 :

 6 Figure 6: A part of a meta-model containing two classes and one reference.

  Constraint → sum() = Res sum(Vars, Res) → count(Value) = Res count(Value,Vars, Res) → size() = Res Gcc(Vars,Vals, Res) → exists(Value) AtLeast(1,Vars,Value)

Figure 7 :

 7 Figure 7: Meta-model in which we want to type n classes.

Figure 8 :

 8 Figure 8: Petri nets meta-model.

  1. Put in three parameters to determine the size of generated model. The first and second are the couple (lB, uB) which are respectively the minimum and maximum number of instances for each class of the meta-model. The third one is the reference upper bound, i.e. the maximum number of instances for the non-bounded references. For example: lB = 2, uB = 2 and rB = 2. 2. For each concrete class c of the Petri nets meta-model, generate a random number lB ≤ size(c) ≤ uB. It represents the effective number of instances of class c. Then create the domains of each class. The size of the root class must be equal to 1. The domains are given by: D(PetriNet) = {1}, D(Transition) = {2, 3}, D(Place) = {4, 5}, D(Arc) = {6, 7}. After this step we can build the instances of all classes of the meta-model. 3. For each class instance, create the variables and the domains of all features including inherited ones. For example, D(marking) = [-20, 20]. 4. Create variables and domains of references. e.g. D(source) = D(Transition) ∪ D(Place). 5. Translate the OCL constraints of the meta-model.

Figure 9 :

 9 Figure 9: Generated Petri net written as an instance of Ecore (top) and with Petri nets syntax (bottom).

  element(Index,Table,Value), where Value is equal to the Index th item of Table, Value = Table[Index].

Note 1. The clique of O(n 2 ) binary constraints of difference can be replaced by a single allDifferent Global Constraint. Example 2. In example 1, the constraints C 1 , C 2 and C 3 can be replaced by the global constraint allDifferent(x 1 , x 2 , x 3 ) resulting in speed-up of the resolution.

2.1.2 element

The constraint element is a global constraint used to get a value from a table depending of an index. A global constraint of type element has this shape:

Table 2 :

 2 Experiments on meta-models of different sizes. The resolution time is given in CPU time (seconds). Detailed results: http://www2.lirmm.fr/ ˜ferdjoukh/ details.pdf.

	Models #class instances

Table 3 :

 3 Cabot et al.. Concerning the other existing approaches, we did not find any available tool for comparison. Resolution time ratio for the current modelling compared to our previous work

	Meta-model	#Classes % improvement
	Graph Colour. Petri nets ER DiaGraph Jess UML class Feature Ecore Business process BIBT E X Bmethod Royal&loyal Sad3	2 5 5 5 7 7 8 17 25 28 33 34 39	20.3 16.0 30.2 48.2 63.1 13.1 34.0 19.0 31.9 69.2 80.4 29.0 32.3

Table 4 :

 4 Comparing the resolution time for some metamodels constrained with OCL constraints. Time is given in seconds for average of 15 executions. #Inst.: total number of class instances, Diff: Difference.

	Meta-model	Petri nets	ER	Feature	Graph Colour.
	#Inst.	OCL ? Yes No Yes No Yes No Yes No OCL ? OCL ? OCL ?
	50 0.39 0.38 0.38 0.36 0.38 0.38 0.36 0.33
	100 0.55 0.54 0.54 0.47 0.55 0.54 0.43 0.39
	300 1.65 1.64 1.42 1.26 1.65 1.64 0.62 0.47
	600 14.3 14.4 7.79 7.22 14.3 14.5 0.85 0.60
	Diff.	0.02	0.2	0.02	0.12

Abscon solver : http://www.cril.univ-artois. fr/ ˜lecoutre/software.html.

Used meta-models: http://www2.lirmm.fr/ ˜ferdjoukh/english/research.html

Acknowledgment

The authors would like to thank Félix Vonthron and Olivier Perrier for the work during their internship.