Xavier Caruso

Random matrix over a DVR and LU factorization

Let R be a discrete valuation ring (DVR) and K be its fraction field. If M is a matrix over R admitting a LU decomposition, it could happen that the entries of the factors L and U do not lie in R, but just in K.

Having a good control on the valuations of these entries is very important for algorithmic applications. In the paper, we prove that in average these valuations are not too large and explain how one can apply this result to provide an efficient algorithm computing a basis of a coherent sheaf over A 1 K from the knowledge of its stalks.

Throughout the paper, we fix a ring R equipped with a discrete valuation v R : R → N ∪ {∞}. We assume that v R is normalized so that it takes the value 1 and we fix an element π ∈ R such that v R (π) = 1. We also assume that R is complete with respect to the distance defined by v R . The residue field of R and its fraction field are denoted by k and K respectively. The valuation v R extends uniquely to K and we continue to denote this extension by v R . We finally set q = Card k and assume that q is finite. Two typical examples of this are (1) R = Z p (the ring of p-adic integers) equipped with the usual p-adic valuation and (2) R = k[[x]] (the ring of power series) where k is a finite field.

If d is a positive integer, we denote by Ω the ring of square matrices of size d with coefficients in R. It is a compact additive group whose Haar measure is denoted by µ. We assume that µ is normalized so that (Ω, µ) is a probability space (i.e. µ(Ω) = 1). Thus, it makes sense to study some statistics on Ω. Surprinsingly, the literature around this subject seems to be very poor. Nevertheless related questions were already addressed by Abdel-Ghaffar in [START_REF] Abdel-Ghaffar | The determinant of random power series matrices over finite fields[END_REF] and Evans in [START_REF] Evans | Elementary divisors and determinants of random matrices over a local field[END_REF]: the main result of [START_REF] Abdel-Ghaffar | The determinant of random power series matrices over finite fields[END_REF] is the computation of the law of the random variable "valuation of the determinant" (in the case where R is a power series ring but his argument works for a more general discrete valuation ring) whereas, in [START_REF] Evans | Elementary divisors and determinants of random matrices over a local field[END_REF], Evans studies the random variable "valuation of the elementary divisors".

In this paper, we are mainly interested in the random variable V L : "valuation of the L-part in the LU decomposition". We give several estimations of its law, its expected value and its standard deviation. Roughly speaking, we prove that E[V L ] = log q d + O(1) and σ(V L ) = O(1) where the notation O(1) refers to a quantity bounded by a universal constant. We also bound from above the probability that V L deviates from its expectation. For more precise statements, we refer to Theorem 1.1, Theorem 1.2 and Corollary 1.3 in the introduction of §1.

In §2, we move to algorithmic applications. Firstable, we propose in §2.1 a stable algorithm to compute a LU decomposition of a matrix over R (unfortunately standard Gauss elimination is far for being stable) and analyze closely the losses of precision it generates in average (which turn out to be optimal in some sense). §2.2 is devoted to a study of the notion of "simultaneous PLU decomposition", which will play an important role for our next (and main) application presented in §2.3. This application is of geometric nature. We let X = A1 K be the affine line over K. Recall that a coherent subsheaf of F ⊂ O d X (where d is some integer) is determined by the data of all its stalks F x ⊂ O d X,x at all closed points x ∈ X. Furthermore, we know that all such sheaves F as above admit a global basis. In §2.3, we describe an algorithm that computes a basis of F knowing all its stalks and, once again, analyze its stability (which will turn out to be rather good).

Some statistics related to LU decomposition

If A is a (commutative) ring, we shall denote by M d (A) the ring of square d×d matrices. Recall that we have endowed Ω = M d (R) with its Haar measure. Choosing a matrix at random with respect to this measure is just choosing independently each entry at random with respect to the Haar measure on R. Furthermore, since R is complete, every element x ∈ R can be written uniquely as an infinite sum x = ∞ i=0 a i π i where the coefficients a i 's are taken in a fixed set R ⊂ R of representatives of elements of k (i.e. the restriction to R of the canonical projection R → k is bijective) and conversely, any sum ∞ i=0 a i π i as above converges and then defines an element in R. With this description, generating a random element (with respect to the Haar measure) of R is just choosing at random all coefficients a i 's in R independently and uniformly.

We shall say that a matrix M ∈ M d (K) admits a LU decomposition if it can be factorized as a product L(M ) • U (M ) where:

• L(M ) is a unit 1 lower triangular matrix with coefficients in K, and

• U (M ) is a upper triangular matrix with coefficients in K.

We underline that, even if M has coefficients in R, we do not require that L(M ) and U (M ) belong to M d (R). Here are some well known facts: [START_REF] Abdel-Ghaffar | The determinant of random power series matrices over finite fields[END_REF] an invertible matrix M admits a LU decomposition if and only if all its principal minors do not vanish and (2) when it exists, a LU decomposition is unique (i.e. the matrices L(M ) and U (M ) are uniquely determined). We will consider L and U as two partially defined functions on M d (K). For ω ∈ Ω such that L(ω) is defined, let us denote by V L (ω) the opposite of the smallest valuation of an entry of L(ω). The aim of this section is to study the random variable V L . Here are the main results we will prove.

Theorem 1.1. Setting E(q, d) = ∞ v=1 1 -(1 -q -v ) d
(1)

we have E(q, d) -1 q-1 < E[V L ] E(q, d). Furthermore, the distance between E(q, d) and log q d is bounded by 1 log(2) (and by 1 if q 3).

Theorem 1.2. For all positive real number ℓ, we have:

P |V L -log q d -1 2 | > ℓ + 1 2 q q -1 • q -ℓ • 2 + ℓ • log q .
Corollary 1.3. The standard deviation of V L is bounded by an explicit universal constant (which can be chosen equal to 6.5).

Some useful tools

This subsection gathers some preliminaries to the proof of Theorem 1.1, Theorem 1.2 and Corollary 1.3. We first recall some basic facts about LU decomposition, then introduce the random variables V i,j 's (which will play a crucial role in the sequel) and finally prove several important properties of them.

Cramer's rule for LU decomposition

Let M ∈ M d (K). A useful formula for our purpose is an analogue of Cramer's rule which gives a closed expression of the entries of L(M ) as a quotient of two determinants. This formula appears for instance in [START_REF] Householder | The Theory of Matrices in Numerical Analysis[END_REF], §1.4; let us recall it briefly. If I and J are two subsets of {1, . . . , d}, we denote by M I,J the submatrix of M obtained by deleting all columns and rows whose index are not in I and J respectively. The i-th principal minor is then the determinant of the matrix M I,I where I = {1, . . . , i}; we will denote it by δ i (M ). With these notations, we have:

if i > j : L(M ) i,j = det M I,J δ j (M )
where I = {1, . . . , j -1, i} and J = {1, . . . , j}

if i j : U (M ) i,j = det M I,J δ i-1 (M ) (2)
where I = {1, . . . , i} and J = {1, . . . , i -1, j}.

The proof of these formulas is not difficult. For Formula (2), note that L(M ) I,J •U (M ) J,J = M I,J provided that J has the particular shape J = {1, . . . , j}; then, passing to the determinant, we get det L(M ) I,J • det U (M ) J,J = det M I,J and the desired relation follows by combining these equalities for I = J and I = {1, . . . , j -1, i}. The proof of Formula (3) is similar.

1.1.2 The random variables V i,j
The aim of this paragraph is to define a collection of mututally independent random variables V i,j : Ω → N ∪ {∞} (1 i j d); they will be very useful in the sequel to study V L . We first construct a collection of random variables X i,j : Ω → R (1 i j d). The construction goes by induction on j. We start with a matrix ω in Ω. We first define X 1,1 (ω) to be the top left entry of ω. We then enter in the second round (i.e. j = 2). As before, we begin by letting X 1,2 (ω) denote the (1, 2)-th entry of ω but, before defining W 2,2 (ω) we do the two following modifications on the matrix ω:

• if the valuation of X 1,2 (ω) is less than the valuation of X 1,1 (ω), we swap the two first columns of ω and, then

• we clear the (1, 2)-th entry of ω by adding to its second column a suitable multiple of its first column (note that it is always possible because if the top left entry -which serves as pivot -vanishes, so does the (1, 2)-th entry).

Doing these operations, the coefficient of Ω in position (2, 2) may have changed and we define X 2,2 (ω) to be the new (2, 2)-th entry of ω. The general induction step works along the same ideas. Assume that, after the (j -1)-th step, we have ended up with a matrix ω such that ω i ′ ,j ′ = 0 when i ′ < j ′ < j. We define X i,j (ω) by induction on i by applying the following process successively for i = 1, 2, . . . , j -1:

• first, we set X i,j (ω) to the (i, j)-th entry of (the current) ω;

• second, if the valuation of X i,j (ω) is less than the valuation of the (i, i)-th entry of (the current) ω, we swap the first row of ω with its i-th one;

• third, we clear the (i, j)-th entry of ω by adding to its j-th column a suitable multiple of its first column.

We finally let X j,j (ω) denote the j-th diagonal entry of (the current) ω. For all (i, j) with 1 i j d, we also set V i,j = v R (X i,j ) and V i = V i,i . The V i,j 's take values in N ∪ {∞} and they are finite almost everywhere. Algorithm 1 summarizes the construction of the V i,j 's.

Proposition 1.4. The random variables X i,j (1 i j d) are uniformely distributed and mutually independent.

Proof. Set I = { (i, j) | 1 i j d }. Suppose we are given a family x = (x i,j ) (i,j)∈I of elements of R. We consider the following set:

Ω(x) = ω ∈ Ω | X i,j (ω) = x i,j , ∀(i, j) ∈ I .
Algorithm 1: The construction of the random variables V i,j 's Notation: : ⋆ ω i,j denotes the (i, j)-th entry of ω : ⋆ ω j denotes the j-th row of ω

1 for j from 1 to d do 2 for i from 1 to j -1 do 3 V i,j ← v R (ω i,j ); 4 if v R (ω i,j ) < v R (ω i,i
) then swap ω j and ω i ;

5 if ω i,i = 0 then ω j ← ω j - ωi,j ωi,i • ω i ; 6 V j,j ← v R (ω j,j ); Set v i,j = v R (x i,j
) and for all i, let j i denote the first index in {1, . . . , d

+ 1 -i} such that v i,i-1+ji is equal to min(v i,i , v i,i+1 , . . . , v i,d
). This sequence of integers (j i ) is the code of a certain permutation σ of {1, . . . , d} defined by the following rule. We write all the integers between 1 and d. We define σ(1) to be the j 1 -th written integer (that is j 1 ) and we erase it. We then define σ(2) to be the j 2 -th integer which remains written (that is j 2 is j 2 < j 1 and j 2 + 1 otherwise), we erase it and we continue. Let I x denote the subset of {1, . . . , d} 2 consisting of couples (i, j) such that i > σ -1 (j). One can check that it has cardinality d(d-1)

2

. Consider the function f x : Ω(x) → R Ix mapping ω to the family (ω i,j ) (i,j)∈Ix . Following the construction of the X i,j 's, one can check that f x is a bijection. Now, we globalize the previous construction. Let U be a subset of R I containing a distinguished element x and such that I y = I x for all y ∈ U . With this assumption, the collection of functions f y 's (y varying in U ) defines a bijection between

Ω(U ) = { ω ∈ Ω | (X i,j (ω)) (i,j)∈I ∈ U } and U × R Ix .
It is morever easy to check that this bijection preserves the measure; in other words

P[(X i,j ) (i,j)∈I ∈ U ] = µ(U ) (4)
where µ denotes the Haar measure on R I . But, since the function v R is locally constant on R\{0}, any open subset U ⊂ (R\{0}) I can be written as a disjoint union of subsets U ′ on which the function y → I y is constant. Therefore the equality (4) holds for all these U . Since furthemore the complement of (R\{0}) I in R I is a measure-zero set, the equality [START_REF] Evans | Elementary divisors and determinants of random matrices over a local field[END_REF] holds for all open subset U of R I .

Corollary 1.5. The random variables V i,j (1 i j d) are mutually independent and they all follow a geometric law of parameter 1 -q -1 (i.e. they take value v with probability (1 -q)q v-1 ).

Proof. Clear after Proposition 1.4.

Another interest of the V i,j 's is that they are closely related to V L . The following Proposition precises this relationship.

Proposition 1.6. We have max(V 1 , V 2 , . . . , V d ) -v R (det) V L max(V 1 , . . . , V d-1 ) (recall that V i = V i,i by definition).
Proof. Let ω ∈ Ω. To avoid confusion, agree to call T j (ω) the matrix ω computed by Algorithm 1 (run with ω as input) after the j-th iteration of the main loop and reserv the notation ω for the matrix we have started with. It follows from the construction that T j (ω) has the following particular shape: if i ′ < j ′ j, then the (i ′ , j ′ )-th entry of ω j vanishes. Moreover, clearly, T j (ω) is obtained from ω by performing successive elementary operations on the first j columns. Therefore, if J = {1, . . . , j} and if I is a subset of {1, . . . , d} of cardinality J, we have det ω I,J = ± det(ω j ) I,J . In particular these two determinants have the same valuation. Fix a couple (i, j) such that 1 j i d and set I = {1, . . . , j -1, i}, J = {1, . . . , j}. From Formula ( 2) and what we have said before, we derive:

v R (L(ω) i,j ) = v R (det T j (ω) I,J ) -v R (det T j (ω) J,J ) = v R (T j (ω) i,j ) -v R (T j (ω) j,j ) = v R (T j (ω) i,j ) -V j (ω).
Since all coefficients of ω j lie in R, so does its determinant. It follows that v R (T j (ω) i,j ) 0 and consequently that v R (L(ω) i,j ) -V j (ω), which proves the second inequality. To establish the first one, note that ω and T j (ω) share the same determinant up to a sign. Thus there must exist an index i, necessarily not less than j, such that v R (T j (ω) i,j ) v R (det ω). For this particular i, we have v R (L(ω) i,j ) v R (det ω) -V j (ω) and then V L (ω) V j (ω) -v R (det ω). The conclusion follows.

Remark 1.7. In the same way, we can prove that the valuation of the i-th minor of ω ∈ Ω is equal to j i=1 min(V i,i (ω), V i,i+1 (ω), . . . , V i,j (ω)). Combining this with Corollary 1.5, one can easily recover Abdel-Ghaffar's formula

d i=1 1
q i -1 (see Theorem 3 of [START_REF] Abdel-Ghaffar | The determinant of random power series matrices over finite fields[END_REF]) giving the expected value of the random variable "valuation of the determinant".

Proof of the main results

Estimation of the expected value

This subsection is devoted to the proof of Theorem 1.1.

Estimation of the expected value of

V L Let V = max(V 1 , V 2 , . . . , V d ).
The event "V < v" occurs if and only if V i,i < v for all index i, and Corollary 1.5 shows that it happens with probability (1 -q -v ) d . The expected value of V is then equal to

∞ v=1 P[V v] = ∞ v=1 1 -(1 -q -v ) d that is exactly E(q, v). On the other hand, Proposition 1.6 implies that E[V ] -E[v R (det)] E[V L ] E[V ].
Moreover, by Abdel-Ghaffar's Theorem, we know that the expected value of v R (det) is given by

d i=1 1
q i -1 and hence is less that

d i=1 1
q i < ∞ i=1 1
q i = 1 q-1 .
The first part of Theorem 1.1 is proved.

Estimation of E(q, v) Consider the function f :

x → 1 -(1 -q -x ) d .
It is decreasing on the interval [0, ∞) and therefore one can write:

∞ 0 f (x)dx E(q, d) ∞ 1 f (x)dx -1 + ∞ 0 f (x)dx.
Doing the substitution y = 1 -q -x , we get:

∞ 0 f (x)dx = 1 log q • 1 0 1 -y d 1 -y dy = 1 log q • 1 0 (1 + y + y 2 + • • • + y d-1 )dy = H d log q
where

H d = 1 + 1 2 + • • • + 1 d
is the harmonic series. It is well known that γ + log d H d 1 + log d where γ is the Euler's constant. Therefore E(q, d) is almost equal to log q d, the error term being bounded by a universal constant. The second part of Theorem 1.1 follows.

Some additional remarks

We would like first to emphasize that the difference E(q, d) -log q d does not converge to 0 when q and/or d goes to infinity. Indeed the following Lemma shows that, when log q d is far from an integer and q is large, E(q, d) might be closer to the integral part of log q d than to log q d itself.

Lemma 1.8. For all q qnd d,

E(q, d) -[log q d] < q q -1 • q -dist(log q d, N)
where [log q d] and dist(log q d, N) denotes respectively the integral part and the distance to N of log q d.

Proof. We claim that the function f : x → 1 -(1 -q -x ) d satisfies:

1 - q x d f (x) d q x , for all x 0. (5)
Indeed, the second inequality directly comes from the standard inequality (1 + t) d 1 + td whereas the first one is a consequence of AM-GM inequality applied with the numbers dq -x and 1 -q -x , 1 -q -x , . . . , 1q -x (d times). If v 0 = [log q d], we then get v 0 + v0 v=1

q v d E(q, d) v 0 + ∞ v=v0+1 d
q v , which gives:

- q q -1 • q v0 d E(q, d) -v 0 q q -1 • d q v0+1 .
The Lemma follows from this.

Let us end this paragraph by a last remark: the sum E(q, d) can also be exactly computed. Indeed, we have:

E(q, d) = ∞ v=1 1 -(1 -q -v ) d = ∞ v=1 d k=1 (-1) k-1 d k q -vk = d k=1 (-1) k-1 d k • 1 q k -1 .
Nevertheless, this expression does not yield the order of magnitude of E(q, d); indeed, each term in the latter sum (the one over k) can individualy be very large whereas the sum itself grows rather slowly.

Estimation of the law of V L

We now start the proof of Theorem 1.2. The strategy is quite clear: we use Corollary 1.5 and Proposition 1.6 to bound from below and from above the distribution function of V L . First, let us investigate the consequences of the inequality V L V (where we recall that we have set

V = max(V 1 , . . . , V d ))
. For all (nonnegative) real number x, it implies that:

P[V L < x] P[V < x] = d i=1 P[V i < x] (1 -q -x ) d 1 -d • q -x . (6)
It is a bit more tricky to use the other inequality

V L V -v R (det) because v R (det)
and the V i 's are certainly not independent (cf Remark 1.7). Nevertheless, one can pick two nonnegative real numbers x and t and consider the event E x,t : "V > x + t and v R (det) t". It is clear that V -v R (det) is always greater than x when E x,t occurs. Thus we have:

P[V L x] P[V -v R (det) x] P[E x,t ] 1 -P[V x + t] -P[v R (det) > t]. (7)
Moreover we know that P[V x + t] (1 -q -x-t ) d and from Abdel-Ghaffar's result (see [START_REF] Abdel-Ghaffar | The determinant of random power series matrices over finite fields[END_REF]), we

derive P[V d > t] q -t+2
q-1 . Indeed, Abdel-Ghaffar Theorem states that for all integer v, the equality

P[v R (det) v] = (1 -q -v-1 )(1 -q -v-2 ) • • • (1 -q -v-d ) holds. In particular P[v R (det) v] 1 - d i=1 q -v-i 1 -q -v q-1 . Taking v = [t],
we get the claimed result. Putting these inputs in (7), we obtain:

P[V L x] 1 -(1 -q -x-t ) d - q -t+2 q -1 .
This estimation being true for all t, one can optimize it on t. For simplicity, let us define u = 1 -q -x-t ; the variable u now varies in [1 -q -x , 1], and for all u in this range, one have

P[V L > v] 1 -f (u) where f (u) = u d + dλ(1 -u), λ = q x+2
d(q-1) . Assume that λ < 1. A quick study of f shows that it is minimal when u = u 0 = λ 1/(d-1) . Moreover, one can check (using AG-MG inequality for instance) that u 0 always lies in the interval

[1 -q -x , 1]. It follows that P[V L x] 1 -f (u 0 ) = 1 -λ • (d -(d -1)u 0 ).
We can further simplify this formula and write a bound depending only on λ. For this, remark that λ

(1 + log λ d-1 ) d-1 . Raising to the power d -1, we find u 0 1 + log λ d-1 and then:

P[V L x] λ(1 -log λ) where λ = q x+2 d(q -1) . ( 8
)
We are now ready to prove Theorem 1.2. Let ℓ be a positive real number and define v 0 = log q d -1 2 . Applying Formulas ( 6) and ( 8)

with x = v 0 + (ℓ + 1 2 ) and x = v 0 -(ℓ + 1
2 ) respectively, we find:

P[V L v 0 + (ℓ + 1 2 )] q -ℓ and P[V L v 0 -(ℓ + 1 2 )] q q -1 • q -ℓ • 1 -log q q -1 + ℓ • log q q q -1 • q -ℓ • (1 + ℓ • log q).
Theorem 1.2 follows by adding these two inequalities. Corollary 1.3 can be now easily deduced. Indeed, note that the function v → E((V L -v) 2 ) is maximal when v is equal to the expected value of V L and the value taken at this optimal point is the variance of V L . It is then enough to bound the expected value of (V L -v 0 ) 2 , which can be done as follows:

E[(V L -v 0 ) 2 ] = ∞ 0 P (V L -v 0 ) 2 x • dx 1 4 + ∞ 0 P (V L -v 0 ) 2 (ℓ + 1 2 ) 2 • (2ℓ + 1) • dℓ 1 4 + q q -1 • ∞ 0 q -ℓ • (2 + ℓ • log q) • (2ℓ + 1) • dℓ = 1 4 + q q -1 • 3 log q + 8 log 2 q .
The standard deviation of V L is then always less than σ(q) = 1 4 + q q-1 • ( 3 log q + 8 log 2 q ). The function σ is decreasing on [2, ∞) and then bounded from above by its value at 2 (which is < 6.5). Note furthemore that when q goes to infinity, σ(q) = 1 2 + O( 1 log q ).

Generalization to block LU decomposition

Let d = (d 1 , . . . , d r ) be a tuple of positive integers such that d 1 + • • • + d r = d. By definition, a block LU decomposition of type d of a matrix M ∈ M d (K) is a factorization M = L d (M ) • U d (M )
where L d (M ) and U d (M ) are respectively block unit lower triangular and block upper triangular with respect to the partition d:

L d (M ) =       I d1 0 • • • 0 ⋆ I d2 . . . . . . . . . . . . . . . 0 ⋆ • • • ⋆ I dr       and U d (M ) =       ⋆ • • • • • • ⋆ 0 . . . . . . . . . . . . . . . . . . 0 • • • 0 ⋆      
where the s-th block has size d s and, for an integer n, I n denotes the identity matrix of size n. Of course, a block LU decomposition of type (1, 1, . . . , 1) is nothing but a standard LU decomposition and every matrix M ∈ M d (K) admits a block LU decomposition of type (d), which is simply M = I d •M . As in the standard case, a LU decomposition of type d is unique (when it exists) -which justifies the notations L d (M ) and U d (M ) -and, an invertible matrix M admits such a decomposition if and only if, for all i ∈ {1, . . . , s}, its d i -th principal minor does not vanish. For d as before and ω ∈ Ω, we let V L,d denote the opposite of the smallest valuation of an entry of L d (ω). This defines a random variable V L,d : Ω → N ∪ {∞} for each d. The aim of this subsection is to study them. Following the same strategy as in the standard case (i.e. d = (1, . . . , 1)), our first task is to establish a link between V L,d and the random variables V i,j defined in §1.1.2. To shorten notations, we set

I s = {d 1 + • • • + d s-1 + 1, . . . , d 1 + • • • + d s } and recall that if M ∈ M d (K)
and I, J ⊂ {1, . . . , d}, we denote by M I,J the submatrix of M consisting of entries whose row index and column index are in I and J respectively. For all s ∈ {1, . . . , t}, we further introduce:

V d,s = i∈Is min(V i,i , V i,i+1 , . . . , V i,d1+•••+ds ).
Corollary 1.5 learns us that the V d,s 's are mutually independant for s varying between 1 et r (and d remains fixed). The following Lemma shows that their laws are also precisely known.

Lemma 1.9. For all s ∈ {1, . . . , r} and all integer v, we have:

P[V d,s v] = (1 -q -v-1 )(1 -q -v-2 ) • • • (1 -q -v-ds ).
Proof. Throughout this proof, we set

a = d 1 + • • • + d s-1 , b = d 1 + . . . + d s and, for i ∈ {1, . . . , d s } and i ′ = b + 1 -i, W i = min(V i ′ ,i ′ , V i ′ ,i ′ +1 , . . . , V i ′ ,b
). It follows from Corollary 1.5 that W i follows a geometric law of parameter (1 -q -i ) and furthermore that the W i 's (1 i d s ) are mutually independant. For all ℓ ∈ {1, . . . , d s }, define moreover

S ℓ = W 1 + • • • + W ℓ . Clearly S ds = V d,s .
We will prove by induction on the couple (ℓ, v) (lexicographically ordered) that:

P[S ℓ v] = (1 -q -v-1 )(1 -q -v-2 ) • • • (1 -q -v-ℓ ).
For ℓ = 1, the statement is true. Assume now that it is true for all (ℓ ′ , v ′ ) with ℓ ′ < ℓ or ℓ ′ = ℓ and v ′ < v.

The strategy is to decompose the event "S ℓ v" in two parts according to the vanishing or the nonvanishing of W ds-ℓ+1 . Clearly, if W ds-ℓ+1 = 0, we have S ℓ = S ℓ-1 . On the other hand, if we know for sure that W ds-ℓ+1 does not vanish, one can subtract 1 to it and get this way a new random variable which still follows of a geometric law with the same parameter. Hence, one can write:

P[S ℓ v] = P[W ds-ℓ+1 = 0] • P[S ℓ-1 v] + P[W ds-ℓ+1 > 0] • P[S ℓ v -1] = (1 -q -ℓ ) • P[S ℓ-1 v] + q -ℓ • P[S ℓ v -1].
Replacing P[S ℓ-1 v] and P[S ℓ v -1] by their values (coming from the induction hypothesis), we get the desired result.

Remark 1.10. Alternatively, one can notice that V d,s follows the same law as the variable "determinant of a random matrix of size d s " and then conclude by Abdel-Ghaffar's Theorem. Actually the proof we have presented above is very inspired by Abdel-Ghaffar's one.

Proposition 1.11. We have max(V d,1 , . . . , V d,r ) -v R (det) V L,d max(V d,1 , . . . , V d,r ).
Proof. We follow the lines of the proof of Proposition 1.6. To avoid confusion, we begin by letting T j (ω) denote the matrix ω computed by Algorithm 1 (run with ω as input) after the j-th iteration of the main loop. Pick some s ∈ {1, . . . , r} and set j(s) = d 1 + . . . + d s . We are going to prove the two following statements from which the Proposition will follow directly:

• the determinant of the square

d s × d s matrix T j(s) (ω) Is,Is has valuation V d,s (ω);
• for all t ∈ {s, . . . , r}, we have the identity

L d (ω) It,Is • T j(s) (ω) Is,Is = T j(s) (ω) It,Is .
The first assertion is easily proved. Indeed, by construction, the submatrix T j(s) (ω) {1,...,j(s)},{1,...,j(s)} is lower triangular and that its i-th diagonal entry has valuation min(V i,i , V i,i+1 , . . . , V i,j(s) ). To prove the second assertion, we first remark that, up to replacing ω by ω + π N for a sufficiently large integer N , one may assume that ω is invertible. All the matrices T j(s) (ω) Is,Is are then also invertible. Consider the matrix L ∈ M d (K) whose i-th column is the i-th column of T j(s) (ω) where s is the unique such that i ∈ I s . It is apparently lower block triangular with respect to the partition d. Furthermore, noting that, for i ∈ I s , the i-th column of T j(s) (ω) is a linear combination of the first j(s) columns of ω, we see that L -1 • ω is upper block triangular. Hence, if D is the diagonal block matrix:

D =    T j(1) (ω) I1,I1 . . . T j(r) (ω) Ir ,Ir    the factorization ω = (LD -1 ) • (DL -1 ω) is the LU decomposition of type d of ω. Therefore L d (ω) = LD -1 .
Our claim follows directly from this.

From Lemma 1.9, we easily derive that q

-v P[V d,s v] q q-1 • q -v .
Arguing then as in §1.2, one can prove analogues of the results we have shown before concerning the random variable V L : the expected value of V L,d is equal to log q s + O(1), its standard deviation is a O(1) (where the notation O(1) stands for a quantity bounded by a universal constant which can be made explicit) and, actually, we even have a more precise (but also more technical) estimation of its law in the spirit of Theorem 1.2.

LU decomposition over a DVR: algorithmic issues

LU decomposition is a very basic and important tool when we are doing algorithmics involving matrices, and especially matrices over a complete DVR. But unfortunately, on some particular inputs, computing it may cause important numerical instability; it is the case for instance if the top left entry of the input matrix has a very large valuation (compared to the other entries). The first aim of this second section, is to study this phemonemon; more precisely, following the ideas of §1, we will design a new algorithm to compute LU decomposition (see Algorithm 2) and show that the set of unpleasant inputs for which it is numerically unstable is very small.

In particular, we may expect that if Algorithm 2 is called as a subroutine by an other probabilistic algorithm, it will not never generate important losses of precision. In § §2.2 and 2.3, we will illustrate this idea on a particular example: we will propose a probabilistic stable algorithm (based on LU decomposition) whose aim is to compute a basis of a coherent over A 1 K (where K is the fraction field of a complete DVR) from the knowledge of all its stalks. We keep the general notations of §1: let R be a discrete valuation ring whose valuation v R : R → N ∪ {∞} is assumed to be surjective. Let π be an element of R of valuation 1. Let k (resp. K) denote the residue field (resp. the fraction field) of R and set q = Card k. We recall that v R extends uniquely to K and that, in a slight abuse of notations, we continue to denote by v R this extended map. We recall also that we have set Ω = M d (R) and that this space is endowed with its Haar measure. For ω ∈ Ω, denote by W i (ω) the valuation of the i-th principal minor of ω and set W = max(W 1 , . . . , W d ). Thanks to Abdel-Ghaffar's Theorem (see [START_REF] Abdel-Ghaffar | The determinant of random power series matrices over finite fields[END_REF]), the law of the W i 's is known:

P[W i v] = (1 -q -v-1 )(1 -q -v-2 ) • • • (1 -q -v-i
) for all v 0 and i ∈ {1, . . . , d}. From this, we derive

P[W i > v] q -v
q-1 and then:

P[W > v] d • q -v q -1 (9)
for all nonnegative integer v. Adding all these probabilities, one finds

E[W ] = log q d + O(1)
where, as usual, the notation O(1) refers to a quantity bounded by a universal constant.

Loss of precision in LU decomposition

By Formula (2), we know that the entries of L(M ) can be all expressed as the quotient of one minor by one principal minor. Noting that if x and y are both known with precision O(π N ) and if y has valuation v, the quotient x y is known with precision at least O(π N -2v ), one may expect that a good algorithm computing the LU factorization of ω would shrink the initial precision by a factor π2•W (M) .

Unfortunately, a quick experiment shows that the naive algorithm based on usual Gauss elimination generates losses of precision much more important than that. For example, on a random input matrix M ∈ M 25 (Z 5 ) given with precision O(5 N ), it outputs a matrix L which is in average known up to precision O(5 N -c ) where c ≃ 10 whereas the mean value of 2 • W (M ) is only ≃ 2 • log q d = 4. For matrices of size d = 125, the deviation is amplified: we find c ≃ 50... to be compared to 2 • log q d = 6.

A first simple solution

Our starting remark is the following: it follows from Cramer like formulae [START_REF] Caruso | Semi-simplifiée modulo p des représentations semi-stables : une approche algorithmique[END_REF]

that if M are M ′ are two matrices in M d (R) congruent modulo π N (for some positive integer N ) such that W (M ) < N , that W (M ′ ) = W (M ) and L i,j (M ) ≡ L i,j (M ′ ) (mod π N -2•Wi(M) )
for all i, j ∈ {1, . . . , d} with i > j. In particular, under the previous assumptions, we have

L(M ) ≡ L(M ′ ) (mod π N -2•W (M)
). This result suggests the following method to compute L(M ) with a correct precision when M is a matrix known with precision O(π N ):

• we lift M to a matrix M ′ known with precision O(π N ′ ) for some N ′ > N ;

• we compute W (M ′ ) and L(M ′ ) with our favorite algorithm (e.g. Gauss elimination) 2 ;

• we answer

L(M ) = L(M ′ ) + O(π N -2•W (M ′ ) )
.

By what we have said before, our answer L(M ) is always correct. Furthemore, if N ′ is sufficiently large, then L(M ′ ) will be known with precision at least O(π N -2•W (M ′ ) ) and L(M ) itself will be known with precision O(π N -2•W (M ′ ) ).

It then remains to find a suitable value for N ′ . Of course, it will strongly depend on the algorithm we use to compute L(M ′ ). Let us study a bit the case of Gauss elimination. Since the successive pivots appearing during the elimination have valuations W 1 (M ′ ), W 2 (M ′ ), . . . , W d (M ′ ) and since we are only dividing by pivots, the maximal loss of precision is bounded from above by 2•(W

1 (M ′ )+• • •+W d (M ′ )).
In other terms, using Gauss elimination, one can certainly compute L(M ) with precision O(π

N -2•(W1(M ′ )+•••+W d (M ′ )) ).
As a consequence, it is enough to choose N ′ so that:

N ′ -N 2 • W 1 (M ) + • • • + W d (M ) -W (M ) .
However, at the very beginning, we have not computed the W i (M )'s yet. So we cannot figure out at this moment what is the best value of N ′ (i.e. the smallest one satisfying the above inequality). Nevertheless, we know that in average W i (M ′ ) ≃ 1 q and W (M ) ≃ log q d. To begin with, we can then try to take

N ′ = N + ⌈ 2d q ⌉
and see what happens: we do the computation with this particular N ′ , we determine the W i (M )'s, if the above inequality is fulfilled, we are done, otherwise, we determine the right N ′ and redo the computation. Actually, it could happen -but it is very rare -that the first precision O(π N ′ ) does not allow us to determine some of the W i (M )'s; in that case, we just guess a new larger N ′ , try with it and continue like this until it works.

Let us finally analyze the complexity of this method in the favorable case where N ′ = N + 2d q is enough. In order to fix notations, we assume moreover that doing basic operations (i.e. additions, substractions, multiplications and divisions) in R with precision π N requires O((N log q) α ) bit operations where α is some constant 3 , necessarily greater than or equal to 1. Since the complexity of Gauss elimination is O(d 3 ) operations in the base ring, our method needs:

O d 3 • (N + d q ) α • log α q bit operations. If d ≪ qN
, it is quite nice. However, if the opposite situation when d ≫ qN , the dominant term in the above complexity is d 3+α , which is very large and actually not really acceptable for many practical applications.

A stable algorithm to compute LU decomposition

In this subsection, we propose and study a different method to compute LU decomposition which has the advantage of not requiring to increase the precision at any time and whose complexity is comparable to Gauss elimination. Our algorithm is strongly inspired by the constructions of §1 and especially those of §1.1.2. Here is it: : ⋆ A i,j denotes the (i, j)-th entry of a matrix A : ⋆ A j denotes the j-th row of A

1 ω ← M ; 2 L ← identity matrix of size d × d; 3 for j from 1 to d do 4 for i from 1 to j -1 do 5 if v R (ω i,j ) < v R (ω i,i
) then swap ω j and ω i ;

6 if ω i,i = 0 then s ← ωi,j ωi,i lifted to precision O(π N ); ω j ← ω j -s • ω i ; 7 v ← j k=1 v R (ω k,k ); 8 for i from j + 1 to d do L i,j ←
ωi,j ωj,j + O(π N -v-max(0,vR(ωj,j )-vR(ωi,j )) ); 9 return L;

A first important remark related to Algorithm 2 is the following: at each step, all entries of ω are known with precision O(π N ). Indeed, ω itself is updated only on line 6 and the corresponding computation does not affect the precision (because s has been lifted modulo π N previously).

Correctness of Algorithm 2

We fix an integer j ∈ {1, . . . , d} and focus on the matrix ω computed by the algorithm after the j-th iteration of the main loop. It is clear that it is obtained from M by performing a sequence of elementary operations on its j first columns. Thus, for all i > j, we have

L(M ) i,j = det ω (j) I,J det ω (j) J,J
where I = {1, . . . , j -1, i} and J = {1, . . . , j}. On the other hand, by construction, ω (j) I,J and ω (j) J,J are two upper triangular matrices modulo π N . Their determinants are then congruent to the product of their diagonal entries modulo π N . Therefore:

L(M ) i,j = ω (j) 1,1 • • • ω (j) j-1,j-1 • ω (j) i,j + O(π N ) ω (j) 1,1 • • • ω (j) j-1,j-1 • ω (j) j,j + O(π N ) .
Of course, the value of this quotient is

ω (j) i,j ω (j) j,j
up to some precision. To compute this precision, it is easier to work with relative precision (i.e. the difference between the absolute precision and the valuation); indeed, we know that the relative precision of a quotient is equal to the minimum between the relative precisions of the numerator and the numerator. In our case, if we set v = v R (ω

(j) 1,1 ) + • • • + v R (ω (j) j,j ) and w = v R (ω (j) i,j ) -v R (ω (j) j,j
), the relative precision of the numerator (resp. the denominator) is N -(v + w) (resp. N -v). Thus, the relative precision of the quotient is N -v -max(0, w) and its absolute precision is then N -v + min(0, w) (since its valuation is w). The value L i,j computed by Algorithm 2, together with its precision, are then correct.

Precision issues Keeping the previous notations, one certainly have w -v and then N -v+min(0, w) N -2v. In other words, the (i, j)-th entry of the matrix L returned by the Algorithm 2 is known with precision at least O(p N -2Vj (M) ) (recall that V j (M ) denotes the valuation of ω j,j at the end of the j-th loop, i.e. our previous v). The maximal loss of precision is then bounded above by 2

• max(V 1 (M ), . . . , V d (M )).
By the results of §1, we know that the mean of this upper bound is close to 2 • log q d, that is the value we expected.

Algorithm 2 and Hermite normal form

Let us denote by H ′ (M ) the matrix ω computed at the end of the execution of Algorithm 2. It worths remarking that H ′ (M ) has a lot of things to do with the Hermite normal form of M . Let us first agree on the definition of the Hermite normal form of M : throughout this paper, it will refer to the unique lower triangular matrix whose diagonal entries are powers of π and which is right-equivalent to M (it means that H(M ) is obtained from M by multiplying on a right by a unimodular matrix). We will denote it by H(M ).

Proposition 2.1. Let M ∈ M d (R) known with precision π N . We assume that all diagonal entries of H ′ (M ) are not congruent to 0 modulo π N and, for all j ∈ {1, . . . , d}, we write H ′ j,j (M ) = p vj u j where v j is a nonnegative integer and u j is a unit. For all i, j ∈ {1, . . . , d}, we then have:

if i < j : H i,j (M ) = 0 if i = j : H i,j (M ) = π vj if i > j : H i,j (M ) ≡ u -1 j • H ′ i,j (M ) (mod π N -vj ).
Remark 2.2. Keeping the notations of the Proposition, it is clear that u j is only known modulo π N -vj . The congruence of the Proposition is then, by nature, the best one can expect.

Proof. Let W ′ (M ) be the matrix W ′ computed by Algorithm 2. One can easily check that W ′ (M ) is unimodular and moreover that H ′ (M ) = M • W ′ (M ). Consequently the Hermite normal form of M is equal to the Hermite normal form of H ′ (M ).

On the other hand, we know that H ′ (M ) has a very particular shape: firstly, it is lower triangular modulo π N and secondly, by assumption, its diagonal entries are not divisible by π N . Thus, H(M ) is obtained from H ′ (M ) by clearing one by one its entries lying above the diagonal and by dividing its j-th column by u j . But, if H i,j (M ) = π N v i,j (for some pair (i, j) with i < j), one clears the (i, j)-th entry of H ′ (M ) by doing the following elementary operation on columns:

H ′ j (M ) ← H ′ j (M ) -π N -v u -1 i v i,j H ′ i (M )
. Hence clearings do not affect the value of H ′ i,j (M ) modulo π N -vj . The Proposition follows easily from this observation.

The notion of L'V' decomposition

The L-part of the LU decomposition has of course very nice abstract properties but unfortunately does not behave very well regarding to precision. Indeed, as we have seen before, if a matrix M is known modulo π N , it is not true that L(M ) is known with the same precision. But, beyond that, the precision data attached to L(M ) is not uniform in the sense that all entries of L(M ) are not known with the same precision. In order to tackle this problem, we introduce the following definition.

Definition 2.3. Let M ∈ M d (R) and N be a positive integer. A L'V' decomposition of M modulo π N is a couple of d × d matrices (L ′ , V ′ ) such that L ′ ≡ M V ′ (mod π N
) and L ′ and V ′ are lower triangular modulo π N and upper triangular modulo π N respectively.

If there exists a diagonal entry of L ′ which is congruent to 0 modulo π N , (L ′ , V ′ ) is said to be degenerate. Otherwise, it is nondegenerate. Remark 2.4. It is easy to see that if (L ′ , V ′ ) is nondegenerate, then all diagonal entries of V ′ are not congruent to 0 modulo π N as well.

It is not difficult to modify Algorithm 2 so that it computes a L'V' decomposition modulo π N ; we end up this way with Algorithm 3. On the other hand, it is worth noting that L'V' decomposition is closely Algorithm 3: An algorithm to compute a L'V' decomposition

Input: : A matrix M of size d × d known with precision O(π N ) Output: : A L'V' decomposition of M modulo π N 1 ω ← M ; 2 L ′ , V ′ ← two new matrices of size d × d; 3 W ′ ← identity matrix of size d × d; 4 for j from 1 to d do 5 for i from 1 to j -1 do 6 if v R (ω i,j ) < v R (ω i,i
) then swap ω j and ω i ; swap W ′ j and W ′ i ;

7 if ω i,i = 0 then 8 s ← ωi,j
ωi,i lifted to precision O(π N );

9 ω j ← ω j -s • ω i ; W ′ j ← W ′ j -s • W ′ i ; 10 L ′ j ← ω j ; V ′ j ← W ′ j ; 11 return L ′ , V ′ ;
related to LU decomposition. The following proposition makes this statement precise. Proposition 2.5. Let M ∈ M d (R), N be a positive integer and (L ′ , V ′ ) be a nondegenerate L'V' decomposition of M modulo π N . Then M admits a LU decomposition and for all (i, j) with 1 i < j d, one have:

L i,j (M ) ≡ L ′ i,j L ′ j,j (mod π N -vj -max(0,vR(L ′ j,j )-vR(L ′ i,j )) ) with v j = j-1 k=1 v R (L ′ k,k ) -v R (V ′ k,k ).
Proof. Left to the reader (the arguments are very similar to those detailed in §2.1.2).

Of course, executing first Algorithm 3 and then applying the result of Proposition 2.5 is almost the same than running directly Algorithm 2. Nevertheless splitting Algorithm 2 in two parts can be very useful for some applications (we will see an example of this in §2.1.5) because, as we have already said before, the pair (L ′ , V ′ ) is generally easier to manipulate than L(M ) since it carries a flat precision (and, in addition, it consists of two integral matrices if M is itself integral).

Complexity and Hafner-McCauley's algorithm

It is easily seen that the asymptotic complexity of Algorithm 2 is O(d 3 ) (operations in the base ring R) where d denotes the size of the input matrix. It is then similar to the complexity of usual Gauss elimination whereas it is true that our Algorithm 2 runs a little bit more slowly because it basically makes more swaps and copies.

When precision is not an issue (e.g. when we are working over an exact ring), Hafner and McCauley showed in [START_REF] Hafner | Asymptotically fast triangularization of matrices over rings[END_REF] how to reduce the computation of the LU decomposition to matrix multiplication and got this way a nice recursive algorithm that computes the LU decomposition of a matrix in only O(d ω ) operations where ω is the exponent for matrix multiplication 4 . The aim of this subsection is to extend Hafner-McCauley's algorithm in our setting where we want to take care of precision.

A preliminary result about Algorithm 3 Roughly speaking, Algorithm 3 clears the entries of ω lying above the diagonal in the colexigographic order. We would like to study what happens if we decide to clear these entries in a different order. Definition 2.6. Let a < b be two positive integers and set I a,b = {(i, j) ∈ N 2 | a i < j b}. A total order on I a,b is nice if:

• for 1 i i ′ < j d, one always have (i, j) (i ′ , j), and

• for all (i, j) and (i ′ , j ′ ) ∈ I a,b such that j i ′ , one have (i, j) (i ′ , j ′ ).

Remark 2.7. It is easy to check that the colexicographic order on I a,b is nice. However, it is not the only one: the lexicographic order, for instance, is nice as well. One can also build recursively nice orders on I a,b as follows. Fix an integer c between a and b and pick 1 and 2 two nice orders defined on I a,c and I c+1,b respectively. Consider also a third order 3 defined on the cartesian product {a, . . . , c} × {c + 1, . . . , b} and satisfying the first condition of Definition 2.6. Now define a new order on I a,b by agreeing that I a,c {a, . . . , c} × {c + 1, . . . , b} I c+1,b5 and furthermore that agrees with 1 , 2 and 3 on I a,c , I c+1,b and {a, . . . , c} × {c + 1, . . . , b} respectively. A quick check then shows that is nice as well.

If is a nice order on

I 1,d = {(i, j) ∈ N 2 | 1 i < j
d}, let us agree to use the expression "to execute Algorithm 3 with respect to " to mean that we execute this algorithm but, instead of running through all (i, j) ∈ I 1,d according to the colexicographic order, we run through these pairs according to and execute line 10 when i = j -1.

Proposition 2.8. When they are called on the same input, Algorithm 3 and Algorithm 3 executed with respect to a nice order return the same answer.

Proof. Easy check.

Description of the algorithm Suppose that we are given a matrix M ∈ M d (R) known with precision O(π N ). The basic idea (which comes from Hafner and McCauley) is to obtain a recursive algorithm to compute the LU decomposition and, doing this, to replace as much as possible elementary operations on rows by matrix multiplication. Moreover, in order to avoid many problems related to precision, it would be really better to work with L'V' decomposition instead of LU decomposition. Actually, for the purpose of the recursion, we will not just need the matrices L ′ and V ′ but also H ′ (which is the matrix ω at the end of the execution; see §2.1.3) and W ′ . The prototype of the algorithm we want to design is then:

LV : M → (L ′ , V ′ , H ′ , W ′ ).
Proposition 2.8, together with the recursive construction of a nice order detailed in Remark 2.7 suggests the following strategy for a recursive implementation of LV: It turns out that the first step can be computed in a recursive way. Precisely, we decompose M as a block matrix

M = M 1 M 2 M 3 M 4 (where M 1 has size d ′ × d ′ )
we call recursively the routine LV on the input M 1 and then recover the matrices L ′ , V ′ and ω (as they have to be just after the first step) using the following formulas:

V ′ = V ′ 1 0 0 I ; L ′ = M V ′ = L ′ 1 M 2 M 3 • V ′ 1 M 4 ; ω = H ′ 1 M 2 M 3 • W ′ 1 M 4
where the quadruple

(L ′ 1 , V ′ 1 , H ′ 1 , W ′ 1 )
is the output of the recursive call of LV. Last but not least: remark furthermore that, proceeding this way, we are replacing elementary operators (on the columns of M 3 ) by matrix multiplication (by V ′ 1 and W ′ 1 ). It is exactly the benefit we were looking for! Let us now focus on step 2. With the notations above, it consists in clearing all the entries of M 2 (using eventually the diagonal entries of H ′ 1 as pivots). Of course, this can be done just by running the corresponding part of Algorithm 3. Nevertheless, we do not want to proceed exactly along these lines but we would like instead to use a recursive version of this algorithm in order to take advantage again of the complexity of the matrix multiplication. Writing such a recursive version is actually very similar to what we have done before. In order to have more coherent notations, let us rename H ′ 1 and M 2 to X and Y respectively and write:

X = X 1 0 X 3 X 4 ; Y = Y 1 Y 2 Y 3 Y 4 .
Note that X 1 and X 4 are then lower triangular modulo π N . We can then proceed recursively along the following lines:

1. we clear Y 1 using X 1 as pivot;

2. we clear Y 2 using the new X 1 as pivot;

3. we clear Y 3 using X 4 as pivot;

4. we clear Y 4 using the new X 4 as pivot.

Each of these steps can be done recursively. As in the previous case, we just need to be careful and let our recursive routine return not only the new matrix X gotten after clearing Y but also the transformation matrix T such that ( X Y ) • T ≡ ( X ′ 0 ) where X ′ is the new X mentionned previously. Indeed, this matrix is needed to update X and Y after each step.

A brief study of complexity Let us denote T ′ (d) the complexity of the clearing algorithm we have just described (i.e. the number of elementary operations on R performed by this algorithm when the size of the input matrices is d) and by T (d) the complexity of our complete recursive algorithm computing a LV decomposition. From the description of these algorithms, we find:

T (d) = 2 • T ( d 2 ) + T ′ ( d 2 ) + O(d ω ) (10) T ′ (d) = 4 • T ′ ( d 2 ) + O(d ω ) (11)
where we recall that ω is the exponant of the complexity of matrix multiplication. Since a d × d matrix have d 2 entries, one certainly have ω 2. For simplicity, we assume that ω > 2 (we recall that the fastest asymptotic algorithm known today corresponds to ω ≃ 2.376). It is then a classical exercise to deduce from the recursion formula (11 The aim of this paragraph is to prove that this algorithm works quite well in the following sense: not only it finds very quickly a matrix P that does the job, but it furthermore finds quickly a matrix P for which all entries of P , P -1 and the L i 's are kwown with a good precision and do not have a too small valuation. Our precise result can be stated as follows.

) that T ′ (d) = O(d ω ).
Theorem 2.9. Let n be some positive integer. Suppose that for all m ∈ {1, . . . , n} we are given a matrix M m ∈ M d (K) together with a finite sequence d m = (d m,1 , . . . , d m,rm ) of positive integers whose sum equals d. Let ε be a real number in (0, 1) and take v an integer log q ( r1+•••+rn q-1

) -log q ε. Then, a random matrix ω ∈ Ω satisfies the following conditions with probability at least 1 -ε:

• ω is invertible in M d (K) and ω -1 ∈ π -v M d (R);
• for all m ∈ {1, . . . , n}, the matrix ωM m admit a block LU decomposition of type d m and (where V d,s is the random variable defined in §1.3) and set:

L d m (ωM m ) ∈ π -v M d (R); Moreover if the M m 's all lie in M d (R),
W = max(v R (det), W d 1 ,1,1 , . . . , W d 1 ,1,r1-1 , W d 2 ,2,1 , . . . , W d 2 ,2,r2-1 , . . . , W d n ,n,1 , . . . , W d n ,n,rn-1 ).
Since M m is invertible in Ω, Lemma 1.9 implies that

P[W d m ,m,s v] = (1 -q -v-1 ) • • • (1 -q -v-ds ) 1 - q -v q -1
for all m, s and v. Furthermore by Adbel-Khaffar's Theorem (see Theorem 1 of [START_REF] Abdel-Ghaffar | The determinant of random power series matrices over finite fields[END_REF]) we also know the law of the random variable v R (det); we have

P[v R (det) v] = (1 + q -v-1 )(1 + q -v-2 ) • (1 + q -v-d ).
Let us simplify this Formula and just remember that

P[v R (det) v] 1 -(q -v-1 + q -v-2 + • • • + q -v-d ) 1 -q -v q-1 .
We can now estimate the law of W as follows:

P[W > v] P[v R (det) > v] + 1 m n 1 s<rm P[W m,i > v] q -v • r 1 + • • • + r n -n q -1 + 1 q -1 ε. Proposition 1.11 shows that L d m (ωM m ) ∈ π -W (ω) M d (R)
for all ω ∈ Ω and all m ∈ {1, . . . , n} and, on the other hand, it is clear that ω

-1 ∈ π -W (ω) M d (R) because v R (det ω) W (ω). It is enough to conclude the proof.
The second assertion (concerning precision) is now clear.

Remark 2.11. One may wonder if the bound log q ( r1+•••+rn q-1 ) -log q ε is sharp. Actually, it cannot be for any data of (M 1 , . . . , M n ). Indeed, an integer v satisfies the required conditions of Theorem 2.9 for the families (M, . . . , M ) and (d, . . . , d) (n times) if and only if it satisfies the same conditions for the family reduced to the unique matrix M . So if M 1 = • • • = M n and d 1 = . . . = d n , one can certainly improve the bound log q ( rn q-1 ) -log q ε by removing the facteur n in the first log. Nevertheless, by using similar methods as those of §1, one can prove, first, that the result of Theorem 2.9 fails if v ≪ log q ( max(r1,...,rn) q-1 ) -log q ε and, second, that if M 1 , . . . , M n are themselves chosen randomly, it even fails for v ≪ log q ( r1+•••+rn q-1 ) -log q ε (i.e. the given bound is sharp). Lemma 2.14. Let f : Rd → Rr be a surjective map. Then

P f (x) ∈ π v Rr | x ∈ R d q -rv
for all nonnegative integer v.

Proof. Let k denote the residue field of R; it is a finite extension of k. Since f is surjective, it induces a surjective k-linear map f : kd → kr over the residue field. Moreover, the image of f is generated over k by

f (k d ). Thus dim k f (k d ) = dim k f ( kd ) = r.
This fact implies the existence of a k-linear map ḡ : kr → k r such that the composite ḡ • f : k d → k r is surjective. Let g : Rr → R r be any R-linear lifting of ḡ.

The R-linear morphism h = g • f : R d → R r induces a surjection over the residue field and thus is itself surjective. Furthermore, it is clear that h(x) is divisible by

π v if f (x) is. Hence, for x staying in R d , we have P[f (x) ∈ π v Rr ] P[h(x) ∈ π v R r
] and we are reduced to prove the Lemma with f replaced by h. (In other words, we may assume that K = K.) By the structure Theorem for finitely generated modules over a principal domain (recall that R is a principal domain), there exists a basis (e 1 , . . . , e d ) of R d such that the first (d -r) vectors e 1 , . . . , e d-r form a basis of ker h. Now, using that h is surjective, we easily see that a vector x = d i=0 x i e i ∈ R d satisfies h(x) ∈ π v R r if and only if x i is divisible by π v for all i > d -r. But, the probability that such an event occurs is q -rv and we are done.

Corollary 2.15. For all integers v 1 , . . . , v jm(s) , we have:

P[W i v i , ∀i] jm(s) i=1 q -(r+1-i)vi .
Proof. For ω ∈ Ω and i ∈ {1, . . . , j m (s)}, we denote by ω i the i-th row of ω (and consider it as a vector of R d ) and by M (i) (ω) the matrix defined above. Let F (1) be the submatrix of M m consisting of its first j m (s) columns and let f (1) : Rd → Rjm(s) be the R-linear map whose matrix is t F (1) . The fact that M m is invertible implies that f (1) is surjective. Lemma 2.14 applied to f (1) yields:

P[W 1 v 1 ] q -jm(s)v1 . (12)
Now remember that M (2) (ω) is obtained from M (1) (ω) by performing a sequence of elementary operations on columns. It then exists a matrix P (1) (ω) such that M (2) (ω) = M (1) (ω) • P (1) (ω). Clearly P (1) (ω) depends only on ω 1 and we will denote it P (1) (ω 1 ) in the sequel. Set F (2) (ω 1 ) = F (1) • P (1) (ω 1 ) and let f (2) (ω 1 ) : Rd → Rjm(s) denote the map whose matrix is t F (2) (ω 1 ). It is surjective and one can then apply Lemma 2.14 to the composite pr jm(s)-1 • f (2) (ω 1 ) where pr jm(s)-1 : Rjm(s) → Rjm(s)-1 is the projection on the first coordinates. It gives

P[W 2 v 2 | ω 1 = x 1 ]
q -(jm(s)-1)v2 for all x 1 ∈ R d . Integrating now over x 1 and using (12), we get

P[W 1 v 1 and W 2 v 2 ] q -jm(s)v1 • q -(jm(s)-1)v2 .
The Corollary follows by repeating j m (s) times the previous argument.

If we denote by δ s,m (ω) the valuation of the determinant of the j m (s)-th principal minor of ωM m , Corollary 2.15 allows us to do the following computation:

P[δ s,m > v] v1,...,v jm(s) 0 v1+•••+v jm (s) =v+1 P[W i v i , ∀i] v1,...,v jm (s) 0 v1+•••+v jm(s) =v+1 q -(v1+2v2+•••+jm(s)v jm (s) ) = q -v-1 v2,...,v jm (s) 0 v2+•••+v jm (s) v+1 q -(v2+2v3+•••+(jm(s)-1)v jm (s) ) q -v-1 v2,...,v jm(s) 0 q -(v2+2v3+•••+(jm(s)-1)v jm (s) ) = q -v-1 • ∞ v2=0 q -v2 • ∞ v3=0 q -2v3 • • • ∞ v jm(s) =0 q -(jm(s)-1)v jm (s) = q -v-1 • (1 -q -1 ) -1 (1 -q -2 ) -1 • • • (1 -q -jm(s)+1 ) -1 q -v • Π(q) -1 .
It is time now to free s and m: summing the above estimation over all possible s and m, we find that δ(ω) = max s,m δ s,m (ω) is greater than v -which implies that L d (ωM m ) does not lie in π -v M d (R ′ )with probability at most (r

1 + • • • + r n -n) • q -v • Π(q) -1
and consequently that ω does not satisfy the conditions of the first statement of Theorem 2.12 with probability at most:

q -v • 1 q -1 + r 1 + • • • + r n -n Π(q) q -v • r 1 + • • • + r n Π(q) . Hence if v is chosen log q ( r1+•••+rn Π(q)
) -log q (ε), this probability is less than ε: the first part of Theorem 2.12 is proved.

The second part now follows easily: indeed, we know that, on the input ωM m , the Algorithm 2 decreases the precision by a factor that cannot exceed π 2 maxi Vi(ωMm) and so, a fortiori, by a factor that cannot exceed π 2δ(M) where δ is the random variable defined above. The conclusion follows from this.

Remark 2.16. The bound of Theorem 2.12 is sharp if M 1 , . . . , M n are chosen randomly among all square d × d matrices with coefficients in R. However, it is not true in general and it is even not true if M 1 , . . . , M n are chosen randomly among all matrices over R. Indeed, in that case, using results of §1, one can prove that, in average, the better possible bound for v is given by:

[ K : K] -1 • log q ( r1+...+rn Π(q) ) -log q (ε) + O(1)
with an extra factor [ K : K] -1 , which can be very small.

Modules over K[X] and sheaves over A 1

K

Let X denote an affine curve over K and A = K[X] be the ring of regular functions over X. It is well known that the category of coherent sheaves over X is equivalent to that of finitely generated modules over A. In particular, the data of a submodule M ⊂ A d (for some fixed integer d) is equivalent to the data of a coherent subsheaf M ⊂ O d X . Nevertheless, these two objects are of different nature and we would like to represent them in two different ways:

• a submodule M ⊂ A d by a matrix of generators

• a subsheaf M ⊂ O d
X by the data of the stalk M x ⊂ O d X,x for each closed point x ∈ X (note that this inclusion is not trivial for only a finite number of points x).

Since these objects are supposed to be equivalent, it is natural to ask if one can find an efficient way to go from one representation to the other. Actually going from the global description to the local one is quite easy: it suffices to localize at each point x. Contrariwise, going in the opposite direction is not so obvious and will be discuss now.

From now on, we assume for simplicity that X is the affine line A 1 K (and leave to the reader the exercise to extend our constructions to a more general setting). With this extra assumption, the ring A is nothing but the ring of univariate polynomials with coefficients in K.

Rephrasing our problem in concrete terms

For all irreducible polynomials P ∈ K[X], let A P denote the completion of A for the P -adic topology, that is A P = lim ← -r A/P r A. Concretely A P can be identified with a ring of power series with coefficients in the residue field K P = A/P A in one indeterminate X P . This variable X P should be thought as "X -a P " where a P is a (fixed) root of P in K P . Under the identification A P ≃ K P [[X P ]], the natural embedding A → A P is just the Taylor expansion at a P :

F (X) → ∞ i=0 F (i) (a P ) i! • X i P .
Let P 1 , . . . , P n be the minimal polynomials of a 1 . . . , a n respectively and, for simplicity, set A m = A Pm . The question we have addressed earlier is then equivalent to the following: given, for all m ∈ {1, . . . , n}, a submodule M m ⊂ A d m free of maximal rank, how can one find explicitely a A-module M ⊂ A d such that A m ⊗ M = M m (as a submodule of A d m ) for all m and A P ⊗ M = A d P for all other P ? One can actually rephrase again this question in very concrete terms by taking basis everywhere. Indeed, if B is A or one of the A m 's, any free submodule of B d of rank d can certainly be represented by a square d × d matrix with coefficients in B: the module is recovered from the matrix by taking the span of its column vectors. Note furthermore that two matrices G and H defines the same module if and only if there exists an invertible matrix P with coefficients in B such that G = HP ; if this property holds, we shall say that G and H are right-equivalent. Since all our base rings are principal domains, we know that any matrix G ∈ M d (B) admits a factorization G = M DN where M and N are two invertible matrices, D is diagonal and each diagonal entry of D divides the next one. Up to replacing G by a right-equivalent matrix, one can furthermore assume that N is the identity matrix, i.e. that G has the particular form G = M D. Moreover, if B is one of the A m 's, it is safe to assume that the diagonal entries of D are all some powers of the variable X m since all nonvanishing element of A m can be written as a product of an invertible element with a power of X m . In that case the data of D is then reduced to that of a nondecreasing sequence of integers n 1 . . . n d .

With all these remarks, our question becomes: ) ; iv) for all irreducible polynomial P ∈ K[X] which is not one of the P m 's, the matrix M D is rightequivalent to the identity matrix over A P .

The answer

We consider, for all m ∈ {1, . . . , n}, an invertible matrix M m ∈ M d (A m ) together with a nondecreasing sequence of d integers e m,1 • • • e m,d . Our aim is to construct a couple (M, D) satisfying the Conditions i), ii), iii) and iv) above. Firstable, we define the matrix D as follows:

D = Diag(P e1,1 1 • • • P en,1 n , . . . , P e 1,d 1 • • • P e n,d n
).

It clearly satisfies Condition ii).

Lemma 2.18. Let M ∈ M d (A). a) Assume that, for all m ∈ {1, . . . , n}, the matrix M considered as an element of

M d (A m ) (via the natural embedding A → A m ) is congruent to M m modulo X e m,d +1 m . Then, the couple (M, D) satisfies Condition iii). b) Assume moreover that M is invertible in M d (A).
Then the couple (M, D) satisfies Conditions i), ii), iii) and iv).

Proof. Note that in the ring A m , the polynomial P m is equal to the product of X m by a unit whereas all other P m ′ 's (for m ′ = m) are invertible. We deduce from this that D is right-equivalent to D m over M m . Hence, our first hypothesis implies that M D is right-equivalent to a matrix congruent to M m D m modulo X e m,d +1 m

. In other words, there exists a matrix We assume now that M is invertible. Then, clearly, Condition i) holds. Moreover, we have already seen that Conditions ii) are iii) are fulfilled. It is then enough to prove Condition iv). Let P ∈ K[X] be an irreducible polynomial different from all the P m 's. All P m 's are then invertible in A P and, consequently, so is the matrix D. Since M is itself invertible, the product M D belongs to GL d (A P ) and is then rightequivalent to the identity matrix.

Q ∈ GL d (M m ) such that M D is right equivalent to: M m D m + X e m,d +1 m Q = M m D m • I d + X m • Diag(X
It is actually not difficult to produce a matrix M satisfying the assumption of the Lemma 2.18.a). Indeed, the identification

A m /X e m,d +1 m A m ≃ A/P e m,d +1 m A shows that the congruence M ≡ M m (mod X e m,d +1 m ) is equivalent to M ≡ M ′ m (mod P e m,d +1 m
) for a certain matrix M ′ m ∈ M d (A). Hence, finding a convenient M is just a direct application of the Chinese Remainder Theorem (recall that all P m 's are irreducible and pairwise distinct polynomials).

Producing a matrix M satisfying also the second assumption of Lemma 2.18 is a bit more tricky but can be achieved using block LU decomposition. For m ∈ {1, . . . , n}, let r m be the numbers of differents values taken by the sequence (e m,1 , . . . , e m,d ) and let d m,s (1 s r m ) denote the number of times this sequence takes its i-th smallest value. We then have: ). Since M m is invertible, so is U m . Let V m be the matrix obtained from U m by multiplying its (i, j)-th entry by X em,j -em,i m (note that the exponent is always nonnegative when the (i, j)-th entry of U m does not vanish). Obviously, V m is again upper triangular and its diagonal entries are equal to those of U m . Thus U m and V m share the same determinant and V m is invertible. Moreover, we check that for all m. Such a matrix is apparently invertible and also satisfies the assumption in part a) of Lemma 2.18. We can then simply take M = L. Now let us go back to the general case where some M m might not have a block LU decomposition of type d m . In that case, we denote by M m (0) ∈ M d (K m ) the image of M m under the canonical projection

e m,1 = • • • = e m,
M m D m = L m D m V m ,
A m → A m /P m A m = K m (or, equivalently, A m ≃ K m [[X m ]] → K m ).
The coefficients of M m (0) then all lie in K m , which is a finite extension of K. We can therefore apply Theorem 2.12 which implies in particular the existence of a matrix ω ∈ GL ) for all m. The matrix M = ω -1 L then satisfies the two assumptions of Lemma 2.18. Hence, it satisfies also the conclusions of this Lemma and we have solved our problem. Algorithm 5 summarizes the different steps of the proposed solution. Of course, if ω is not invertible or one of the ωM m 's does not admit a block LU decomposition of the required type, Algorithm 5 fails. If it happens, we simply rerun the algorithm again and again until it works: it follows from Theorem 2.12 that we will get the desired answer quite fast.

Let us analyze quickly how much precision is loss in average by this method. In order to fix ideas, let us assume that the entries of the matrix M m are explicitely given as polynomials in K m [X m ] (eventually modulo X e m,d +1 m

) and that all these polynomials are known with precision O(π N ) for some integer N . For Algorithm 5: A solution to Problem 2.17 ) for all m; 7 return (ω -1 L, D); simplicity, we assume moreover that M m (0) has coefficients in the ring of integers R m of K m and that it is invertible in M d (R m ) 6 . Set as before j m (s) = d m,1 + . . . + d m,s and, for all admissible pair (m, s), let denote by D m,s the determinant of the j m (s)-th principal minor of ωM m . Define also δ(ω) to the maximum of all v R (D m,s (0)) when m and s run over all the possibilities. By the proof of Theorem 2.12, we know that δ(ω) is less than v = log q ( 2 Π(q) ) + log q (r 1 + • • • + r n )

with probability at least 1 2 . In many concrete situations, it is not easy to compute exactly the r m 's but it will nevertheless in general quite simple to estimate them. Indeed, going back to the definition, it is clear that r m is less than both d and e m,d and these latter quantites are natural parameters on which we will in general have a good control (cf [START_REF] Caruso | Semi-simplifiée modulo p des représentations semi-stables : une approche algorithmique[END_REF], §3.2 for a concrete example). From now on, we assume that all matrices M m computed on line 5 satisfy this estimation. If this property does not hold, we simply agree to rerun Algorithm 5 until the desired property holds.

The next step is to measure the size of the denominators appearing in the following nonconstant coefficients. In order to do this, we introduce a new parameter w by requiring that all matrices M m have coefficients in the ring R m,w defined as the image of R m [ Xm for all m ′ = m.

In order to bound the loss of precision as we would like to do, we assume for simplicity that all P m 's are entirely known. We introduce again two new parameters. The first one is an integer v 1 for which we require that the image of R[X]/P e d,m +1 m under the isomorphism (13) contains

π v1 • • • R m [X m ]/X e d,m
+1 m for all m. The second parameter is the integer v 2 defined as the opposite of the smallest valuation of a coefficient of the unique polynomial C m of degree < n m=1 (e d,m + 1) deg P m satisfying the above condition. Now, remember that we have proved that L m are known up to an element of π N -2v • R m,v+w . It is then a fortiori known up to an element of π N -2v-e(v+w) where e = max(e 1,d , . . . , e n,d ). Inverting the isomorphism (13), we find that L ′ m is certainly known modulo π N -2v-e m,d (v+w)-v1 . Finally the formula L = C 1 L ′ 1 + . . . + C n L ′ n shows that L is known with precision O(π N -2v-e(v+w)-v1-v2 ) (recall that we have assumed that the P m 's -and consequently the C m 's -are known with infinite precision). The total loss of precision of Algorithm 5 is then bounded by 2v + e(v + w) + v 1 + v 2 .

Remark 2.20. The parameters v 1 and v 2 are not easy to estimate in general. One can nevertheless keep in mind the following: v 1 measures the ramification of the roots a m of the P m 's and v 2 measures the distance between these roots. For instance, to be more precise, one can easily prove that if all a m lie in the ring of integers of an unramified extension K ′ of K then one can just take v 1 = 0. If in addition the a m 's are pairwise distinct in the residue field of K ′ (which is the same as to be distinct in the residue field of K since K ′ /K is unramified), one can also take v 2 = 0. In that very particular case, the computation of line 6 does not generate any loss of precision. We refer to [START_REF] Caruso | Semi-simplifiée modulo p des représentations semi-stables : une approche algorithmique[END_REF] for a quite different example where the constants v 1 and v 2 do not vanish but stay nevertheless under control.

Here is a final important remark. Algorithm 5 still works if, instead of computing (the L-part) of the block LU decomposition of ωM m , we compute a unit lower triangular (and not block unit lower triangular) L m such that there exists a block upper triangular (with respect to d) matrix U m with the property that M m = L m U m . Indeed, the knowledge of these L m 's is enough to compute L (which need to be only unit lower triangular) and then to conclude using Lemma 2.18. This remark is important because Algorithm 4 spends some time in line 10 in clearing entries in order to make the computed matrix L block unit lower triangular instead of simply unit lower triangular. In other words, commenting the line 10 in Algorithm 4 speeds up the execution of Algorithm 5 but do not have any influence on its correctness.

Algorithm 2 :

 2 A stable algorithm to compute the L-part of the LU decomposition Input: : A matrix M of size d × d known with precision O(π N ) Output: : The L-part of the LU decomposition of M Notations: : ⋆ d is the dimension of the matrix M

1 .

 1 we start the computation of a L'V' decomposition of M but stop it after d ′ columns for some d ′ < d (e.g. d ′ = [ d 2 ]); 2. we clear all the entries in the d ′ × (d -d ′ ) top right corner of the matrix ω we have ended up after the first step; 3. we finally compute a L'V' decomposition of the (d -d ′ ) × (d -d ′ ) bottom right corner of ω.

  Knowing this, equation (10) becomes T (d) = 2 • T ( d 2 ) + O(d ω ) and then yields T (d) = O(d ω ) as expected.

  are invertible in this ring and are only known with precision O(π N ), one can furthermore require (without changing the probability) that, on each input ωM m , Algorithm 4 outputs L d m (ωM m ) with precision at least O(π N -2v ). Remark 2.10. It is also possible to bound the loss of precision if we drop the hypothesis of inversibility of the M m 's. The valuations of their determinants then enter into the scene. (The exercise is left to the reader.) Proof. Let us begin by proving the first assertion. Replacing eventually M m by M m +π N for a large integer N , one may assume that all M m 's are invertible in M d (K). Furthermore, since multiplying M m (1 m n) on the right by an upper triangular matrix with coefficients in R does not change the matrix L d m (ωM m ) for any ω ∈ M d (K), one can even safely assume that all M m 's are invertible in M d (R) = Ω. For all m ∈ {1, . . . , n} and all s ∈ {1, . . . , r m }, define W d,s,m : Ω → N ∪ {∞} by W d m ,s,m (ω) = V d,s (ωM m )

Problem 2 .

 2 17. Given for all m, an invertible matrix M m ∈ M d (A m ) and a nondecreasing sequence of d integers e m,1 • • • e m,d , how can one construct explicitely a couple (M, D) of matrices over A such that: i) the matrix M is invertible in M d (A) ; ii) the matrix D is diagonal and each of its diagonal entry divides the next one ; iii) for all m ∈ {1, . . . , n}, the matrix M D is right-equivalent to M m D m over A m where D m =

  e m,d -em,1 m , . . . , X e m,d -e m,d-1 m , 1) • Q where I d is of course the identity matrix. The last factor (the one between brackets) is a matrix over A m congruent to identity modulo X m . It is therefore invertible. It follows that M D is right-equivalent to M m D m , and part a) of the Lemma is proved.

  dm,1 < e m,dm,1+1 = e m,dm,1+2 = • • • < e m,dm,1+dm,2 < e m,dm,1+dm,2+1 = • • • Assume now for a moment that all M m 's admit a block LU factorization M m = L m U m of type d m = (d m,1 , . . . , d m,rm

  from what we derive that M m D m is right-equivalent to L m D m . Since all L m 's are unit lower triangular, there certainly exists a unit lower triangular matrix L ∈ M d (A) which is congruent to L m modulo X e m,d +1 m

  d (K) such that ω • M m (0) has a block LU decomposition of type d m for all m. Lemma 2.19 below shows that this decomposition lifts to a LU decomposition of type d m of ωM m . Lemma 2.19. Let L be a finite extension of K. Pick M ∈ M d (L[[Y ]]) and denote by M (0) its image in M d (L) under the projection L[[Y ]] → L, Y → 0. Assume that M (0) is invertible and admits a block LU decomposition of type d for a certain partition d of d. Then M also does. Proof. Write d = (d 1 , . . . , d r ) and set as usual j(s) = d 1 + • • • + d s for all s ∈ {1, . . . , r}. It is enough to check that, for all s, the j(s)-th principal minor of M , say δ s (M ), is invertible in L[[Y ]]. But, L[[Y ]] being a local ring, δ s (M ) is invertible if and only if its image δ s (M )(0) is. Now remark that this image is nothing but the corresponding minor of M (0): in other words δ i (M )(0) = δ i (M (0)). The invertibility of M (0) together with the fact that it has a block LU decomposition of type d shows that δ i (M (0)) is invertible in L and we are done. We are now in position to argue as above. For all m, write ωM m = L m U m the block LU decomposition of M m of type d. By the Chinese Remainder Theorem, there exists a unit lower triangular matrix L with coefficients in A such that L ≡ L m (mod X e m,d +1 m

2 ω 5 L

 25 ← a random matrix in M d (R); 3 for m from 1 to n do 4 compute d m = (d m,1 , . . . , d m,rm ); m ← L d m (ωM m ) (computed by Algorithm 4); 6 L ← a unit lower triangular matrix in M d (A) such that L ≡ L m (mod X e m,d +1 m

F•

  π w ] in the quotient ring K m [X m ]/X e m,d +1 m .Clearly D m,s belongs to R m,w for all (m, s) and, by we have said before, it has a representant whose constant coefficient has a valuation less than v. Its inverse D -1 m,s then belongs to π -v • R m,v+w and is known up to an element of π N -2v • R m,v+w . All entries of L m will consequently be known with this precision.It remains to analyze the line 6 of Algorithm 5. Note that the matrix L we want to compute can be expressed in terms of the L m 's by the formulaL = C 1 L ′ 1 + • • • + C n L ′ n where:• L ′ m is a matrix with coefficients in K[X] whose reduction modulo P e d,m +1 m corresponds to L m via the natural isomorphism:K[X]/P e d,m +1 m ∼ -→ K m [X m ]/X e d,m +1 m , F (X) → e d,mi=0(i) (a m ) C m is a polynomialcongruent to 1 modulo P e d,m +1 m and divisible by P e d,m ′ +1 m ′

It means that all diagonal entries are equal to 1.

Generally, these two computations can be done simultaneously. It happens in particular if one uses Gauss elimination.

In usual situations, one can take α = 1 + ε for all positive real number ε.

Nowadays, the best known value for ω is 2.376 but, unfortunalety, the corresponding algorithm due to Coppersmith and Winograd (see[START_REF] Coppersmith | Matrix Multiplication via Arithmetic Progressions[END_REF]) is not efficient in practice (even for very large d) because the constant hidden in the O is quite large. A good compromise is to use classical Strassen's algorithm whose asymptotic complexity is a little bit worse -exactly O(d log 2 7 ) -but which is easy to implement and works very well in practice.

By this inequality, we mean that elements in Ia,c are all less than those in {a, . . . , c} × {c + 1, . . . , b} and, in the same way, that the latter elements are less than any pair in I c+1,b .

Otherwise, we would need to take in account the valuation of det Mm(0) as in Remarks 2.10 and 2.13.

for i from 1 to j -1 do 6 if v R (ω i,j ) < v R (ω i,i ) then swap ω j and ω i ; 7 if ω i,i = 0 then s ← ωi,j ωi,i lifted to precision O(π N ); ω j ← ω j -s • ω i ; 8 if j = j 0 + d s then 9 for j ′ from j 0 + 1 to j 0 + d s do L j ′ ← 1 ω j ′ ,j ′ • ω j ′ ; 10 for j ′ from j 0 + 1 to

14 return L;

Block LU decomposition

The results of §2.1.2 extend to block LU decomposition using §1.

). In average, the loss of precision is then bounded by

All other results proved previously for classical LU decomposition (relation with Hermite normal form, notion of L'V' decomposition, Hafner-McCauley's improvement) also extend almost verbatim to block LU decomposition. We will not explain it in details here (but let the exercise to the reader).

Simultaneous PLU decompositions

As we have already said before, a LU decomposition may fail to exist for some particular matrices. Nevertheless, it is well known that all matrices over a DVR admit a PLU decomposition (recall that a PLU decomposition of a matrix M is a factorization M = P LU where P is a permutation matrix and L and U are as before) and, in general, that several matrices P are possible.

Assume now that we do not pick just one matrix, but a (finite) family of matrices (M 1 , . . . , M n ). The question we would like to address is the following: does there exist a "simultaneous PLU decomposition" of the M i 's, that is PLU decomposition of each M i with the same matrix P . If we want as before P to be a permutation matrix, the answer is negative in general. However, if we relax this condition and require only that P is invertible (which is enough for certain applications, see §2.3 for a concrete example), the answer is positive (at least if the base field is infinite).

The aim of this section is to study this notion of "simultaneous (block) PLU decomposition" over a base field which is the fraction field of a discrete valuation ring.

The basic result

Let (M 1 , . . . , M n ) be a family of square d × d matrices over K and fix a partition d = (d 1 , . . . , d s ) of d. There exists an obvious probabilistic algorithm to compute a simultaneous block PLU decomposition (of type d) of the M i 's: we choose P at random and compute the block LU decomposition of the P -1 M i 's.

A slight generalization

For the application we have in mind (see §2.3), we will need a slight generalization of Theorem 2.9 where the matrices M 1 , . . . , M n on the one hand and the matrix ω on the other hand are not defined over the same field. Let K be a finite extension of K. A classical result asserts that the valuation v R extends uniquely to K. Let R be the ring of integers of K, that is the subset of K consisting of elements with nonnegative valuation. Set finally:

It is easy to check that q -1 -1 q-1 < Π(q) < q -1. Theorem 2.12. Let n be some positive integer. Suppose that for all m ∈ {1, . . . , n} we are given a matrix M m ∈ M d ( K) together with a finite sequence d m = (d m,1 , . . . , d m,rm ) of positive integers whose sum equals d. Let ε be a real number in (0, 1) and take v an integer log q ( r1+•••+rn

) -log q ε. Then, a random matrix ω ∈ Ω satisfies the following conditions with probability at least 1 -ε:

• for all m ∈ {1, . . . , n}, the matrix ωM m admit a block LU decomposition of type d and We now start the proof of Theorem 2.12; it will occupy the rest of this subsection. As in the proof of Theorem 2.9, we start with the first assertion and assume that M m is invertible in M d ( R) for all m. However, in our new settings, this fact no longer implies that ωM m runs over M d (R) when ω runs over M d (R). Thus, we can no longer work with the random variables ω → V d m ,i (ωM m ) and we need to modify a bit our strategy. Actually, since we just want to bound from above -and not from below -the valuation of the matrices L d m (ωM m ), we can argue first assuming that m is fixed and then add probabilities. Moreover, by the proof of Proposition 1.11 (see also Formula (2) when d m = (1, . . . , 1)), bounding the valuation of L d m (ωM m ) reduces to bounding the valuation of the (d m,1 + . . . + d m,s )-th minor of ωM m for all s ∈ {1, . . . , r m -1}. Thus we first fix m ∈ {1, . . . , n} and s ∈ {1, . . . , r m -1} and look for an upper bound for the valuation of the j m (s)-th minor of ωM m where, by definition, j m (s) = d m,1 + • • • + d m,s . For 1 i j m (s), we are going to define a random variable W i : Ω → 1 e N ∪ {∞} where e is the ramification index of K/K (i.e. the index of v R (K ⋆ ) is v R ( K⋆ )). The construction of the W i 's is achieved by applying the classical algorithm of LU decomposition. We pick ω ∈ Ω and first set M (1) = (ωM m ) {1,...,rm},{1,...,rm} ∈ M r ( R). Let j be the first index for which v R (M (1) 1,j ) is minimum among the valuations of all entries of the first line of M (1) . Let M (2) be the matrix obtained from M (1) by swapping the j-th column with the first one and by clearing all the entries of the first row (expect the first one) by pivoting, i.e. adding to each column (expect the first one) a suitable multiple of the first one. The matrix M (2)

.

We now continue this process: we select the first index j 2 for which v R (M

2,j ) is minimal, we obtain M (3) by putting the j-th column in the second position and clearing all the other entries on the second row. Repeating this again and again, we obtain a finite sequence M (1) , . . . , M (jm(s)) of matrices and the last one is lower triangular. For i ∈ {1, . . . , j m (s)}, we define W i (ω) as the valuation of the i-th diagonal entry of M (jm(s)) (or equivalently of M (j) for some j i). It is clear that the determinant of the j m (s)-th principal minor of ωM m is equal to W 1 (ω) + W 2 (ω) + • • • + W jm(s) (ω). We need to determine the law and the correlations between the W i 's. We begin by a Lemma.


Imen Tounsi
email: imen.tounsi@redcad.org

Mohamed Hadj Kacem
email: mohamed.hadjkacem@redcad.org

Ahmed Hadj Kacem
email: ahmed.hadjkacem@fsegs.rnu.tn

Khalil Drira
email: khalil@lass.fr

Towards an Approach for Modeling and Formalizing SOA Design Patterns with Event-B

Keywords: Design patterns, SoaML modeling, Event-B method

This paper introduces a formal architecture-centric approach, which allows first to model message-oriented SOA design patterns with the SoaML standard language, and second to formally specify these patterns at a high level of abstraction using the Event-B method. These two steps are performed before undertaking the effective coding of a design pattern providing correct by construction pattern-based software architectures. We implement our approach under the Rodin platform which we use to prove model consistency.

Introduction

Service-oriented architectures (SOA) is a technology that offers a model and an opportunity to solve problems related to the communication and the integration between heterogeneous applications [Erl, 2009]. Nevertheless these architectures are subject to some quality attribute failures (e.g., reliability, availability, and performance problems). Design patterns, as tested solutions to common design problems within a context, have been widely used to solve these weaknesses.

Most design patterns are presented in an informal way that can raise ambiguity and may lead to their incorrect usage. Patterns, proposed by the SOA design pattern community, are described with informal visual notations [Erl, 2009]. The intent of our approach is to model and formalize messageoriented SOA design patterns. These two steps are performed before undertaking the effective coding of a design pattern, so that the pattern in question will be correct by construction. Our approach allows to reuse correct SOA design patterns, hence we can save effort on proving pattern correctness.

In this paper, we introduce a formal architecture-centric approach. The key idea is to model SOA design patterns with the semi-formal Service oriented architecture Modeling Language (SoaML) in order to attribute a standard notation to SOA design patterns, then to formally specify them with the Event-B method. We have tested our approach with pattern examples. From these patterns we quote the Asynchronous Queuing pattern proposed by the SOA design pattern community. We implement the specifications under the Rodin platform which we use to prove model consistency. We provide both structural and behavioral features of SOA design patterns in modeling and formalizing steps. Structural features of a design pattern are generally specified by assertions on the existence of entity types in the pattern. The configuration of the entities is also described, in terms of the static relationships between them. Behavioral features are defined by assertions on the temporal orders of the messages exchanged between the entities [START_REF] Zhu | Laws of pattern composition[END_REF].

The rest of this paper is organized as follows. Section 2 gives an overview of our proposed approach. Section 3 discusses related work. Section 4 concludes and gives future work directions.

Approach overview

The proposed approach mainly consists of two steps. In the first step, SOA design patterns are modeled graphically with the SoaML language [START_REF] Omg | Service oriented architecture Modeling Language (SoaML) Specification[END_REF]. In the second step, the obtained graphical models are formalized with the Event-B method [Abrial, 2010] (Fig. 1).

Pattern Modeling

We provide a modeling solution for describing SOA design patterns using a visual notation based on the graphical SoaML standard language. To specify structural features, we use Participant diagram that allows modeling entities that make up the pattern's architecture, their types and their dependencies (connections). We also use ServiceInterface and MessageType diagrams to respectively model interfaces of entities and exchanged messages. To specify behavioral features, we use UML2.0 sequence diagram that provides a graphical notation to describe dynamic aspects of design patterns.

Entities, that make up the architecture of an SOA design pattern, can be either P articipants or Agents. It is a P articipant when it provides and/or consumes services and it is an Agent when it can adapt to and interact with its environment. Agents are also participants, providing and using services. Entities can have ports that constitute interaction points with their environment. The port type can be either service or request. The communication path between Services and Requests within an architecture is called ServiceChannel.

Pattern Formalization

For the formalization of SOA design patterns, we use the Event-B method.We use the Rodin Platform [START_REF] Abrial | Rodin: An Open Toolset for Modelling and Reasoning in Event-B[END_REF] in order to prove the correctness of the pattern specification.

A design pattern is described with structural features and behavioral features. Structural features are specified with one or several contexts P Ci and behavioral features are specified with one or several machines P M i. To reduce the complexity of pattern formalizations, we define specification levels. In the first level, we create a very abstract model (a context P C0 and a machine P M 0). In the next levels, we use the refinement techniques to gradually introduce detail and complexity into our model until obtaining the final pattern specification. When we move from Level(i) to Level(i+1), we add a new entity and its connections to the model. In Level(i+1), the context P Ci is extended with the context P C(i + 1) and the machine P M i is refined with the machine P M (i + 1). The refined machine sees the extended context.

Related work

Research connected to design patterns in the field of software architecture are mainly classified into three branches of work according to their architectural style. The first is about design patterns for Object-Oriented Architectures [START_REF] Gamma | Design Patterns: Elements of Reusable Object-Oriented Software[END_REF], the second is about design patterns for Enterprise Application Integration (EAI) [START_REF] Gregor | Enterprise Integration Patterns -Designing, Building, and Deploying Messaging Solutions[END_REF], and the third is about design patterns for Service Oriented Architectures (SOA) [Erl, 2009]. Most proposed patterns are presented in an informal way, using a combination of textual description and a graphical notations in order to make them easy to read and understand [START_REF] Gamma | Design Patterns: Elements of Reusable Object-Oriented Software[END_REF],Gregor and Bobby, 2003,Erl, 2009]. However, using these descriptions makes patterns ambiguous and may lack details. Therefore, several researches have proposed the formalization of these patterns. Since the most famous ones are those proposed by Gamma [START_REF] Gamma | Design Patterns: Elements of Reusable Object-Oriented Software[END_REF], most researches refer to these patterns [START_REF] Zhu | Laws of pattern composition[END_REF] [ [START_REF] Taibi | Formal specification of design pattern combination using BPSL[END_REF],Dong et al., 2007,Kim and Carrington, 2009]. These researches use several formal techniques to define pattern specifications like the Balanced Pattern Specification Language (BPSL) [START_REF] Taibi | Formal specification of design pattern combination using BPSL[END_REF], TLA [START_REF] Dong | Composing pattern-based components and verifying correctness[END_REF] and Object-Z [START_REF] Kim | [END_REF]Carrington, 2009] [Dong et al., 2007]. All these research address object-oriented design patterns, however in our research work we are interested in SOA design patterns and we are based on patterns defined by Erl [Erl, 2009]. Erl presents his patterns with an informal proprietary notation. So, in our work, we propose to model SOA design patterns with the SoaML standard language and we focus on their structural and behavioral features.

Conclusions

In this paper, we introduce a formal architecture-centric design approach supporting modeling and formalizing message-oriented SOA design patterns. The modeling phase allows to represent SOA design patterns with a graphical standard notation using the SoaML language. The formalization phase allows to formally specify both structural and behavioral features of these patterns at a high level of abstraction using the Event-B method. So far, we have implemented the elaborated specifications under the Rodin platform. We also illustrated our approach through a pattern example within the "Service messaging patterns" category. Currently, the transition from the SoaML modeling to the formal specification is achieved manually, we are working on automating this phase by implementing transformation rules.

Acknowledgments

This paper is done with the support of the Ministry of Higher Education and Scientific Research of Tunisia within the Tunisian-French scientific cooperation (DGRS/CNRS).

  Fig. 1. Approach overview