Anna Urba
email: anna.urbanska@imag.fr

Towards an exact adaptive algorithm for the determinant of a rational matrix

In this paper we propose several strategies for the exact computation of the determinant of a rational matrix. First, we use the Chinese Remaindering Theorem and the rational reconstruction to recover the rational determinant from its modular images. Then we show a preconditioning for the determinant which allows us to skip the rational reconstruction process and reconstruct an integer result. We compare those approaches with matrix preconditioning which allow us to treat integer instead of rational matrices. This allows us to introduce integer determinant algorithms to the rational determinant problem. In particular, we discuss the applicability of the adaptive determinant algorithm of [9] and compare it with the integer Chinese Remaindering scheme. We present an analysis of the complexity of the strategies and evaluate their experimental performance on numerous examples. This experience allows us to develop an adaptive strategy which would choose the best solution at the run time, depending on matrix properties. All strategies have been implemented in LinBox linear algebra library.

INTRODUCTION

The determinant computation is one of the core problems in linear algebra. To our knowledge, the problem of the exact computation of the determinant of a rational matrix (i.e a matrix with rational entries) has not yet been widely studied. In general, exact algorithms can be used everywhere where large precision is required. For example, the determinant can be too close to 0 or ±∞ and thus cannot be computed by floating point precision algorithms. In the case of ill-conditioned matrices symbolic methods can be preferred as rounding errors can spoil the computation. It can also be interesting to compare the use of decimal and continued fractions approximations of the entries of realvalued matrices. Continued fractions are the best approximants with small denominators, see [START_REF] Gathen | Modern Computer Algebra[END_REF]Ch. 4]. In this paper, we will try to face the question of how efficient an exact determinant computation can be in both cases.

LinBox library [START_REF] Dumas | LinBox: A Generic Library for Exact Linear Algebra[END_REF] implements exact algorithms for the determinant computation in the case of modular and integer domains. By using fast modular routines [START_REF] Dumas | FFLAS: Finite field linear algebra subroutines[END_REF][START_REF] Dumas | FFPACK: finite field linear algebra package[END_REF] it can offer solutions an order of magnitude faster than other existing implementations [START_REF] Dumas | An introspective algorithm for the integer determinant[END_REF]. We apply these procedures to the computation of the determinant of a rational matrix.

Rational field arithmetics is implemented in GMP [START_REF]GNU Multiprecision Package[END_REF] and Givaro [START_REF]Givaro library[END_REF] libraries. In general, rational numbers are difficult to treat from the exact computation point of view. Mainly, the size of the numerator and denominator can increase very quickly with every addition and multiplication. When we add or multiply two fractions with numerators and denominators bounded by M , the numerator and denominator of the result are bounded by O(M 2 ). Moreover, one addition requires 3, and one multiplication requires 2 integer products, as well as a gcd computation. Therefore, the cost of an exact matrix-vector or matrix-matrix product can be prohibitive in practice. This prohibits the use of the rational field Q in most exact linear algebra algorithms which rely on matrix-matrix or matrixvector products.

However, the cost of computing a modular image of a rational number a b , where a, b are of moderate size, should be comparable with the cost of computing a modular image of a large integer number. This allows us to compute a modular image of a rational matrix at a reasonable cost and thus enables us to use modular procedures.

To compute the determinant of a rational matrix A = a ij b ij , bij > 0 the problem of matrix storage has to be considered. First, we can store the entries of A as rational numbers. Furthermore, one could store the common denominator D(A) of all entries of A and an integer matrix A ′ given by the formula A = 1 D(A) A ′ . This approach can be useful in the case when the entries of A are decimal fractions and D(A) can be set to a power of 10. But if we only assume that the values |aij |, bij are less than M , both D(A) and A ′ are bounded by O(M n 2 ). Still, we may store the common denominator for each row (column) separately. Then the integer vectors Ãi are given by the equation Ai = 1 D i Ãi, where Ai is the matrix row (column) and Di is the common denominator of its entries. Vectors Ãi form matrix Ã, the norm of which is bounded by O(M n ). The product πDi gives a more accurate approximation of the denominator D(det(A)) than D(A) n . The purpose of this paper is to propose the strategies to compute the denominator of a rational matrix. All approaches are based on modular computation. Depending on the matrix storage determinant and/or matrix preconditioning is proposed. The resulting algorithms can use the rational reconstruction [START_REF] Gathen | Modern Computer Algebra[END_REF]Ch.5] and/or existing integer determinant algorithms.

The rest of the paper is organized as follows. In section 2 we give a short description of the existing algorithms for the rational reconstruction and the integer determinant problem. In section 3 we present the main result i.e. two preconditioning strategies and four new algorithms to compute the rational determinant. The cost of the algorithms can be described in terms of the number of modular images of A and modular determinant computations needed. Depending on the strategy, the cost of the rational reconstruction or p-adic lifting is taken into account. In section 4 we discuss the cost of computing a modular image of a matrix and the overall cost of the algorithms. In section 5 we present the experimental results and discuss the best choice of the strategy in practice. We conclude the paper by proposing some mixed solutions in section 6.

EXISTING ALGORITHMS

The aim of this section is to introduce the algorithms that will be used later in section 3. In subsection 2.1 we give a short description of the rational reconstruction procedure. On the example of padic system solving [START_REF] Dixon | Exact Solution of Linear Equations Using P -Adic Expansions[END_REF], we present the application of this procedure to the computation of a rational solution. We show how to change the procedure in the case of early terminated reconstruction [START_REF] Wang | A p-adic Algorithm for Univariate Partial Fractions[END_REF] and give the complexity estimation in this case. Then in subsection 2.2 we present the classical CRA algorithm for the determinant and its modifications by [START_REF] Abbott | Fast deterministic computation of determinants of dense matrices[END_REF] and [START_REF] Dumas | An introspective algorithm for the integer determinant[END_REF].

Rational reconstruction and its application

A modular image of a rational number a b mod M can be computed by taking the modular images of a and b and applying the modular division. This fact can be written as

a b = u mod M ⇔ a = bu mod M.
It should be noticed is that the opposite procedure can also be performed. One can reconstruct the fraction a b where gcd(a, b) = 1, b > 0 from it modular image u. The solution is usually not unique but when we additionally require that |a| < N 2 , b ≤ M N , then there exists at most one solution, see [START_REF] Gathen | Modern Computer Algebra[END_REF]Ch.5].

The solution to the rational reconstruction problem can be computed by applying the extended Euclidean algorithm EEA which searches for the gcd of M and u. The procedure Ratrec(a,b,u,M, N, D) takes as the input modulus M , u ∈ Z and the bounds N and D, and returns a fraction a b = u mod M such that |a| < N, b < D or FAIL if no such solution exists. The worst case complexity of Ratrec is thus the same as for the EEA algorithm i.e. O log 2 (M ) for the classical algorithm and O(log(M ) log(log(M ))) for the fast Euclidean algorithm, see [START_REF] Gathen | Modern Computer Algebra[END_REF]Ch.11]. We will use the notation EEA(M ) for the complexity of the Extended Euclidean Algorithm with entries bounded by M . In many application, the cost of rational reconstruction is usually small compared with the cost of computing u and M . The general scheme is to recursively compute u k , M k , where

M k = p1p2 • • • p k or M k = p k until M k > 2N D
and then to apply the rational reconstruction. The complexity of the procedure depends on the number k of steps, which can be quite large Reducing the number of steps can be the easiest way to enhance the performance of the algorithm.

This can be seen on the example of the Dixon algorithm [START_REF] Dixon | Exact Solution of Linear Equations Using P -Adic Expansions[END_REF] to solve a linear system Ax = b of integer equations. Let N, D be the bound for the numerator and denominator of x. In the classical approach we compute the p-adic approximation in k > log(N ) + log(D) + 1 steps and then reconstruct the result, which gives the complexity O m 3 (log(m) + log( A )) 2 when we use the bound of Hadamard for D, N and assume b ∈ O(1). See [START_REF] Mulders | Diophantine Linear System Solving[END_REF] for a detailed complexity study. In fact, the number of entries in x which we need to reconstruct can often be reduced, see [START_REF] Dumas | LinBox: A Generic Library for Exact Linear Algebra[END_REF].

One should however notice, that the bounds N and D can be much bigger than the actual result. The idea is therefore to apply the rational reconstruction periodically and check the solution for correctness. If M k = p k is the modulus in the current step, the method of Wang [START_REF] Wang | A p-adic Algorithm for Univariate Partial Fractions[END_REF] prompts us to set

M k 2
as the current bound for numerators and the denominator in Ratrec. The algorithm is guaranteed to return the result if

M k > 2 max(N (x) 2 , D(x) 2 ),
where N (x), D(x) are the values of the numerator and the denominator. In the opposite case,

Ratrec(a, b, u, M k , M k 2 , M k
2 ) should fail with large probability. If we apply Wang's idea to the p-adic lifting we can reduce the number of steps to k = 2 log p (max(N (x), D(x)))+1 and the complexity becomes O m ω + m 2 k log(m A ) Current work on this field focus on further reducing the number of steps in the case when N (x) ≪ D(x) or D(x) ≪ N (x). A purely heuristic idea is to use the bounds

M k 2 N D , M k 2 D N instead of M k 2 .
For other approaches, see [START_REF] Khodadad | Fast rational function reconstruction[END_REF][START_REF] Olesh | The vector rational function reconstruction problem[END_REF].

Integer Determinant Algorithms

For an integer matrix A one has several alternatives to compute the determinant. The classical approach is to use Chinese Remaindering Algorithm (CRA) to reconstruct the value from sufficiently many modular images. The modular determinant is computed by LU factorization in the time O(n ω ), where n is the matrix dimension.Each step of the algorithm consist of computation mod pi and a reconstruction of the determinant mod p1 • • • pi by the Chinese Remaindering Theorem. The computation is stopped when the early termination (ET) condition is fulfilled i.e. the reconstructed result rests the same for several iterations. The algorithm is Monte Carlo type, where the probability of success is controlled by the number of repetitions. See [START_REF] Dumas | On Efficient Sparse Integer Matrix Smith Normal Form Computations[END_REF][START_REF] Dumas | An introspective algorithm for the integer determinant[END_REF] for a detailed description.

A mixture of CRA loop and Dixon p-adic lifting is used to compute the integer determinant in [START_REF] Abbott | Fast deterministic computation of determinants of dense matrices[END_REF] and in the hybrid algorithm of [START_REF] Dumas | An introspective algorithm for the integer determinant[END_REF]. The principle is to reduce the value reconstructed by CRA algorithm by computing a large fraction of the determinant. By solving several linear systems we can compute some largest invariant factors sm, . . . sm-i. Their product π is potentially a large part of the determinant. An early terminated CRA loop which reconstructs det(A)/π mod p0p1 • • • pi usually requires only a few modular determinant computations. Informally, the algorithm can be described as follows.

1. For i = 0 to k do (a) Solve Axi = bi by Dixon p-adic lifting to find sm, . . . sm-i; Here, k should not exceed the expected number of invariant factors which is O( log(n)) see [START_REF] Dumas | An introspective algorithm for the integer determinant[END_REF]. The expected complexity of the hybrid determinant algorithm [START_REF] Dumas | An introspective algorithm for the integer determinant[END_REF] for random dense matrices is O n 3 log 2.5 (n A ) . In the worst case (step 2) we can choose between the CRA algorithm and the algorithms of [START_REF] Storjohann | The shifted number system for fast linear algebra on integer matrices[END_REF][START_REF] Eberly | On computing the determinant and Smith form of an integer matrix[END_REF][START_REF] Kaltofen | On the complexity of computing determinants[END_REF]. In fact, in the expected case we do not need to run this step. The experiment proved that thanks the adaptive solutions this algorithm performs better than other implementation for a larger group of matrices.

(b) π = sm • • • sm-i; (c)

RATIONAL DETERMINANT AL-GORITHMS

The algorithms to compute the rational determinant are based on the ideas described in section 2. We present four main strategies to compute the rational determinant. They all use CRA which allows us to compute the determinant of the matrix modulo a product p1 • • • p k of primes. Then the first variant uses the rational reconstruction to obtain the rational result. In order to make use of Early Termination condition we have to precondition the determinant to obtain its integer multiplication. Preconditioning of the matrix allows us to use the integers determinant algorithms. The application of two determinant algorithms is studied here. The common requirements for all algorithm are shown in 3.0. The algorithms are Monte Carlo type due to the early termination used. The effectiveness of our methods depends heavily on the number of modular determinants computed and thus on the bound N and D for the numerator and the denominator of the determinant. One can compute D as the product of lcm of all denominators in a row (or a column). Then N can be computed as D • H, where H is the Hadamard bound for matrix A. One should notice that the bounds can be largely overestimated. Thus, we proposed output-dependant approach which allows us to reduce the number of iteration.

The first idea is to employ the CRA scheme and compute the determinant for the modular images of a rational matrix. In the case when the determinant is rational, early termination condition never holds. Instead, we have to compute the bounds D and N for the denominator and numerator of the determinant. As soon as the product of primes M = p1 • • • p k overcomes 2N D we can apply rational reconstruction and reconstruct the determinant from the modular image. We can also use an output dependent rational reconstruction as described in section 2.1. This strategy is presented as algorithm RatLU. An early termination in the rational case would required applying the rational reconstruction from time to time with the bounds N = D = M k 2 and wait for the result to re-occur. This leads to solution when M > 2 max{n 2 , d 2 }, where n, d are the numerator and denominator of the determinant.

The second method can use the denominator bound D to make the CRA loop look for an integer value. Again, we compute the modular image of a rational matrix A but this time we call CRA to look for D × det(A) which is integer. Now the classic ET condition can be used and the result is obtained as soon as M > n D d . The effectiveness of this method Algorithm 3.1 RatLU 1: i = 0, k = 0, n = 0, d = 1, M = 1, u = 0; 2: repeat 3: ++i; Get pi from P ; 4:

Compute Ai = A mod pi; 5:

Compute ui = det(Ai); 6:

Reconstruct u = det(A) mod M pi using M, u, ui, pi, M = M pi; 7: Compute Ai = A mod pi;

if i = k 2 then 8: s = Ratrec(n, d, u, M, M 2 , M 2
5: Compute ui = D • det(Ai); 6: reconstruct u = D • det(A) mod M pi using M, u, ui, pi, M = M • pi 7: if ET holds then return u gcd(u,D) , D gcd(u,D) ; end if 8: until M > 2N D 9: return u gcd(u,D) , D gcd(u,D) ;
The last two strategies require an integer matrix à which can be obtained by preconditioning the rational matrix A. In order to obtain an integer matrix, the easiest way would be to take matrix A ′ = D(A)A, where D(A) is the common denominator of all entries. In the general case, where the entries of A are fractions a ij b ij with numerator and denominator bounded by A , this is not the best choice as the size of D(A) can be as large as O( A m 2 ). This causes log( A ′ ) to be O(m 2 ). Moreover, the denominator approximation is D(A) m in this case, which is O(m 3 ) in size. We have already defined a tighter bound for the denominator of det(A) by πDi, which is O(m 2 ) in size. Now, if we want to use the integer matrix à then we can precondition A by taking à = A diag(Di), where Di are the common denominators of the rows (or à = diag(Di)A, where Di are the common denomina-tors of the columns). For the preconditioned matrix à all integer determinant algorithms can be applied. In particular the hybrid determinant algorithm of [START_REF] Dumas | An introspective algorithm for the integer determinant[END_REF] can be used. The drawback of this approach is the size of the coefficients of à compared to A, see section 5 table 1. This forced us to use early terminated rational reconstruction for system solving in the Dixon p-adic lifting algorithm. The strategies that use the CRA loop or the hybrid algorithm are presented as algorithms PrecMatLU and PrecMatDixon respectively.

COMPLEXITY ANALYSIS

In this section we study the complexity of the algorithms presented in section 3. In subsection 4.1 we present the analysis of the general case, where we assume that the entries of the matrix are fractions with numerators and denominators bounded by A . Then, in subsection 4.2, we will focus on two special cases i.e. matrices of decimal fractions and Hilbert matrices.

The complexity of the strategies described in section 3 depends on the number of iterations required by the while loop of CRA. Then, depending on the strategy, we have to include the cost of computing the homomorphic image of the matrix, the cost of the rational reconstruction or the cost of p-adic lifting. If we use the early termination condition, the number of steps required for the computation of det(A) depends on the values: m -the size of the matrix, n, d -the real values of the numerator and denominator of det(A) and D -the bound for the denominator. The cost of homomorphic imaging depends on the maximum norm of the matrix i.e. A = max{ aij , bij } and à .

General case

We start this section by the analysis of the rational homomorphic imaging schemes. We have the following lemma.

Lemma 4.1. Let p be a word-size prime. Then the complexity of computing the modular image at p for a rational matrix A is O(m 2 (log( A ))+EEA(p)) word operations.

Proof. For a matrix without a pattern we compute an image for all m 2 entries. For a rational fraction the cost is O(log( A )) for the computation of the modular image of the numerator and denominator and EEA(p) = O(log(p) log(log(p))) for the modular inverse computation by fast extended Euclidean algorithm. Therefore for a word-size A the cost of computing the image is O(1) yet important, due to the constant for computing the inverse of an element mod p.

For the integer case, the cost is log( à )). We can notice that log( à ) can be O(m log( A )) in the worst case, so the complexity of homomorphic imaging in terms of m is O(m 2 ) for the rational and O(m 3 ) in the integer case. But if à < p the cost of imaging for one element is 1. Thus, if both à and A are less than p, the complexity of the homomorphic imaging becomes m 2 EEA(p) for the rational and m 2 for the integer case. In this case, it is better to use integer imaging. On the other hand, if matrix A is structured, for example it is Hankel-type, we have the complexity m EEA(p) for rational imaging. Due to the preconditioning, we loose the structure pattern for à and the complexity of integer imaging rests without change. Finally we notice, that for sparse matrices with Ω elements, we can take Ω instead of m 2 in the complexity formula.

Putting it together we have the following theorem. Special care should be taken if we consider the use of alg. PrecMatDixon. As à can potentially be O ∼ (m) in size and with a pessimistic bound on x, its worst case complexity can be O ∼ (log(m 4 )), which is worse than for the CRA computation. Nevertheless, the gain of computing sm can be important, as it is the case in the HybridDet algorithm, see [START_REF] Dumas | An introspective algorithm for the integer determinant[END_REF].

Complexity in the special cases

By the precedent remarks it should be visible, that the analysis of the strategies should be divided into two main cases. One would consist of the matrices, whose entries are given by decimal fraction, or more generally, where the common denominator of all entries, the common denominator of the rows and the norm of A are of the same order i.e. D(A) = O(Di) = O( A ). In the other case matrix entries are given as fractions with different denominators. We will study the complexity of the algorithms on the example of Hilbert matrices.

In the case of matrices of decimal fractions let us further assume that A is O(1). This would be the case of numerous ill-conditioned matrices emerging from different applications in science and engineering. In order to better describe the differences between the algorithms, we include the cost of EEA when it is relevant. The theorem is a straightforward consequence of theorem 4.2. where k, x are as in theorem 4.2.

The analysis suggests that the algorithm PrecMatLU should be better than PrecDetLU (see 4.1 for the homomorphic image complexity). The performance analysis in section 5 confirms this observation. Furthermore, as long as the Smith form of à is simple, we encourage the use of strategy PrecMatDixon. In particular, we can establish an equivalence between matrices A of random decimal fractions with e decimal places taken randomly an uniformly from the interval [0, 1] and matrices Ã, à < 10 e . This allows us to use the expected complexity of the hybrid algorithms of [START_REF] Dumas | An introspective algorithm for the integer determinant[END_REF] as the expected complexity of the rational determinant computation by alg. PrecMatDixon. Also, the preconditioning should be used instead of strategy RatLU. For more details see section 5.

The other group consists of matrices with rational entries given by fractions with very different denominators. As a model case we can consider Hilbert matrices. Hilbert matrices are the matrices of the form Hm = [ 1 i+j-1 ]i,j=1..m. They are benchmarks examples for many numerical methods. The formula for the determinant of a Hilbert matrix is well known and is given by the equation

1 det(Hm) = Π m-1 k=1 (2k + 1) 2k k 2 .
Theorem 4.4. The complexities for rational determinant strategies in the case of Hilbert matrices are In the case of Hilbert matrices algorithm PrecDetLUhas the best time complexity and also performed best in the experiments, see section 5. Since the numerator is equal to 1, we only have to recover the size of the over-approximation. Experimental results show, that its size is equal to about 8% of the denominator size. Therefore, alg. PrecDetLU, PrecMatLU perform about 25 times less iterations than RatLU. As for the algorithm PrecMatDixon, the study of the Smith form of Hm has revealed that it is quite complex, with about 2 √ m nontrivial factors and the size log(sm( Hm)) equal O(m). Thus, it is not worth computing PrecMatDixon due to the high cost of the algorithm and poor gain.

PERFORMANCE COMPARISON

In this section we present the experimental results for four strategies from section 3. We have tested the performance of four strategies on three matrix sets: random, ill-conditioned and Hilbert matrices.

We generated the random matrices using Matlab procedure rand. The entries of the matrices are decimal fractions with 6 decimal places chosen randomly from the interval [0, 1]. The determinant of the resulting matrices is large in the absolute value. The result of the numerical procedure of Matlab is ±∞.

Ill-conditioned matrices have been chosen from the Matrix Market [START_REF]Matrix Market[END_REF] Harwell-Boeing collection. We chose three sets: Grenoble, Astroph and Bcsstruc3. Grenoble set represents the results of the simulation of computer systems. The sizes of the matrices varies from 115 to 1107 and the condition numbers range from 1.5 • 10 2 in the case of the smallest matrix to 9.7 • 10 7 for the biggest. The decimal precision of the entries depends on the matrix and ranges from 1 to 5 decimal places. The determinants are close to 0. For these matrices, Matlab procedure det computes the result correctly up to the 5th decimal place. Since matrix entries seem to be represented as rounded expansions of rational numbers, we computed the determinant of the matrices "as is" and then we took continued fractions approximants of the entries with the same precision as the decimal fractions.

Astroph set describes the process of nonlinear radiative transfer and statistical equilibrium in astrophysics. The condition number is 3.6 • 10 17 for the small 180 × 180 matrix and 1.7 • 10 14 for the 765 × 765 one. The result of Matlab computation is -∞. Bcsstruc3 gives dynamic analyses in structural engineering. All matrices are symmetric. The condition number is about 10 11 for matrices 19 and 20 and 10 5 for matrix 22. The result of Matlab computation is ∞.

We split the analysis of the performance of the algorithms in three phases. First, we will consider the cost of rational-modular vs. integer-modular imaging and compare it with the results for A and à . Then we will take a look on the numerator and denominator approximations D and N computed by our algorithms. Finally, we give the timings for all strategies and compare their performance.

As we can see in table 1, the time of computing an integer image can be several times shorter than for the rational image provided that the size of preconditioned matrix is still small. This is not the case for Hilbert matrices of dimension ≥ 250 , when the time of rational image computation is better. Furthermore, for structured matrices, like Hilbert, we can reduce the number of images computed. For a Hankel-type matrix, there are only 2n -1 images to compute, which makes the cost of imaging negligible.

The performance of the algorithms depends on the accuracy of denominator approximation used. For the bound D = πDi, the resulting size of the overapproximation is shown in table 2, column 4. In algorithm PrecMatDixon we additionally approximate the numerator by computing sm( Ã). In this case we are interested in the value App(N ) = sm( Ã) and D d n sm( Ã) which we compute instead of the numerator. As we can see in the table, the quality of the approximation of the denominator depends on the matrix and ranges from 1-2% in the case of sparse matrices in the Grenoble set, to 80% for Bccstk matrices. For Hilbert matrices the approximation is quite efficient, the over-approximation is always less than 10%. Table 3 The timings for all algorithms are shown in table 3. The results for Hilbert matrices agree with the complexity estimation in Thm. 4.4. Note that alg. PrecMatDixon is usually the best for the matrices from MatrixMarket collection.

For the Grenoble set, the approximation by continued fractions allowed quite well, in our opinion, to reconstruct the orginal rational matrix connected to the problem. Despite the difference in properties, the running times for the decimal and continued fractions variants were simmilar. However, although the matrices were close in the maximum norm, the determinants ratio reached as much as 2 in the case of grenoble1107.

In figure 5 we present the results of the determinant computation for Hilbert matrices. We compare the timings for algorithm RatLU, PrecDetLU, Prec-MatLU, PrecMatDixon, and the Maple LinearAlgebra::Determinant algorithm with method=rational. Hankel structure of the matrix.

CONCLUSIONS

It this paper we have presented four strategies for exact computation of the determinant of a rational matrix. We have evaluated the performance of these algorithms on several sets of matrices. The performance of the algorithms suggests that there exists a clear division between the matrices given as a rational approximation (by decimal fractions) of real valued matrices and the matrices with a great diversity of the denominators of the entries. For the first case, matrix preconditioning which leads to a integer matrix is proposed, which allows us to use integer determinant algorithms, see solution PrecMatDixon. For the second case, determinant preconditioning is preferred, which does not lead to matrix coefficient blow-up. In general, preconditioning proved more useful than rational reconstruction tools, although better early termination methods where the modulus M is linear in the size of the output n and d can bring a change, see [START_REF] Khodadad | Fast rational function reconstruction[END_REF][START_REF] Olesh | The vector rational function reconstruction problem[END_REF].

An adaptive solution should be able to choose the best storage method and homomorphic imaging scheme, and work independently of the determinant overapproximation.

We propose the following solution, which incorporates the elements of all algorithms This algorithm can be further developed to compute other invariant factors as in alg. PrecMat-Dixon if relevant. Notice, that the cost of introducing solution RatLU to the adaptive algorithm is virtually that of rational reconstruction.

Further work can include intertwining algorithms RatLU and PrecDetLU to include the use of less exact determinant preconditioners, which potentially are not a multiple of d. The aim would to reduce a factor of the denominator by preconditioning and reconstruct the remaining part by rational reconstruction. The strategy should be effective, if the over-approximation caused by preconditioning is reduced but a large fraction of the denominator is obtained at the same time. For example, D = πDi/gcd(Di) could be considered.

Further work can then focus on the implementation of the solution in the case of sparse matrices and on the parallelization of the algorithms.

In this paper we have considered the case of dense matrices in the analysis of the complexity of the strategies as well as in the implementation. However, sparse matrix counterparts of the algorithms can also be used. For the modular determinant computation one could used the algorithm of Wiedemann [START_REF] Wiedemann | Solving sparse linear equations over Finite Fields[END_REF] that computes the determinant by finding the characteristic polynomial of the matrix. In alg. PrecMatDixon the sparse solver of [START_REF] Eberly | Solving Sparse Integer Linear Systems[END_REF] can be used.

The strategies described in this paper contain elements that allow parallelization. This concerns in particular the CRA loop, where several iterations can be performed at the same time, see [START_REF] Dumas | Calcul parallele du polynome minimal entier en Athapascan-1 et Linbox[END_REF]. The question of an optimally distributed early termination in the case of integer Chinese reconstruction (alg. PrecDetLU, PrecMatLU, PrecMatDixon) as well as the rational reconstruction (alg. RatLU)

has not yet been addressed. For a parallel p-adic lifting for alg. PrecMatDixon, see [START_REF] Dumas | Exact Solution to Large Sparse Integer Linear Systems[END_REF].

In this paper we have developed and compared four strategies to compute the rational determinant of a matrix. We have proposed two preconditioning methods that allow us to transfer the problem from rational to integer domain. We believe that the approach described in this article can also be applied in other problems of exact computation in rational numbers such as rank computation or system solving.

Requirements 3. 0

 0 Require: A -an m × m rational matrix; Require: Di, i = 1 . . . m -the common denominator of the entries of the ith row (column); Require: N, D -the bounds for the numerator and the denominator of det(A), D = πDi; Require: A set P of random primes; Require: Ratrec(a, b, u, M, N, D) -a procedure which reconstructs a b = u mod M, a < N, b < D or returns FAIL. Ensure: det(A) -the determinant of the matrix.

Algorithm 3. 3

 3 PrecMatLU 1: i = 0; M = 1; u = 0; 2: Compute A = A diag(Di) (or diag(Di)A) 3: repeat 4: Get pi from P ; 5: Compute Ai = A mod pi; 6: Compute ui = det(Ai); 7: reconstruct u = det(A) mod M pi using M, u, ui, pi, M = M • pi 8: if ET holds then return u gcd(u,D) , D gcd(u,D) ; end if 9: until M > 2N D 10: return u gcd(u,D) , D gcd(u,D) ; Algorithm 3.4 PrecMatDixon 1: Compute A = A diag(Di) (or diag(Di)A); 2: Compute u = det(A) by HybridDet [

Theorem 4 . 2 . 1 2

 421 The worst case complexities of the strategies for computing the determinant of a rational matrix A of size m are1. O k(m 2 log( A ) + m ω ) +O ∼ k √ k for RatLU,where O ∼ hides some log(k) factors; 2. O log( D d n)(m 2 log( A ) + m ω ) for PrecDetLU; 3. O log( D d n)(m 2 log( à ) + m ω ) for PrecMatLU; 4. O ∼ (x(m 2 (log(m)+log( à ))+mx m 2 log( à )+m ω )) for PrecMatDixon, where sm = sm( Ã) and x ∈ m(log(m à b ) is the size of solution to Ãx = b. Here à is equal to A diag(Di) as in section 3; n, d are the numerator and denominator of det(A) and k = O(max(log(n), log(d))). Proof. The complexities can be obtained by a careful examination of the number of CRA steps. The result for alg. RatLU takes into account the cost of the rational reconstruction which is performed at most O( √ k) times. In alg. PrecMat-Dixon we introduce x to estimate the cost of early terminated p-adic lifting. The size of x can generally vary depending on the choice of b but is O(m log(m à b )) in the worst case. To further evaluate the worst case complexity of alg. PrecMat-Dixon we assumed that HybridDet continues to use CRA loop in the worst case. Thus the number of iterations O(log( D d n sm )) and the complexity.

Theorem 4 . 3 . 1 2

 431 The complexities of the strategies in the case when A = O( ˜ A ) = O(1) are:1. O ∼ k(m 2 EEA(p) + m ω + k √ k for alg. RatLU; 2. O ∼ log( D d n)(m 2 EEA(p) + m ω ) for alg. PrecDetLU; 3. O ∼ log( D d n)(m 2 + m ω ) for alg. PrecMatLU; 4. O ∼ x(m 2 log(m) + mx )) + log( D d n sm )(m 2 + m ω ) for alg. PrecMatDixon.

1 .

 1 O m 2 log(m)(m ω + m log(m)) for alg. RatLU; 2. O m ω+2 log(m) for alg. PrecDetLU; 3. O m 5 ) log(m) for alg. PrecMatLU; 4. O(smm 3 log 2 (m) + m 5 log(m)) for alg. Prec-MatDixon. Proof. One should notice that log( 1 det(Hm) ) is O(m 2 log(m)). The size of entries of Hm is log( Hm ) = O(log(m)) and log( Hm ) = O(m log(m)).

1 . 2 . 4 .

 124 Compute D = πDi, Ã; set N = 1; If log p ( Ã < C) compute N = sm( Ã ) -see alg. PrecMatDixon 3. Compute the modular image of the rational matrix A and integer matrix Ã, determine whether to use PrecDetLU or PrecMatLU based on the timings. Run the ET CRA loop for D N • det(A) using PrecDetLU or PrecMatLU.

5 .

 5 From time to time check by rational reconstruction the early termination condition on det(A) -see RatLU.

  Experimental results show that it is sufficient in practice, see sec. 5 table 2. This strategy is presented as algorithm PrecDetLU.

		);++k;
	9:	if s = FAIL then return n, d; end if
	10:	end if
	11: until M > 2N D
	12: status = Ratrec(n, d, u, M, N, D);
	13: if status = FAIL then return n, d; end if
	depends therefore on the exactness of denominator
	bound D. Algorithm 3.2 PrecDetLU
	1: i = 0; M = 1; u = 0;
	2: repeat
	3:	++i;Get pi from P ;
	4:

Table 1 :

 1 shows that despite Comparison of the times (in seconds) for homomorphic imaging are given in columns RatIm (for rational) and IntIm (for integer). The ratio of the timings is given in column 3. Last two columns give the size of entries for A and Ã. Matrix size is included in its name.the size of the over-approximation, preconditioning allow us to gain enough to beat the naive RatLU algorithm. If the size of à is small, as is the case for sparse matrices, we can compute sm( Ã) at a relatively low cost and efficiently approximate the numerator.

	A	RatIm IntIm IntIm/RatImlog( A )log( Ã )
	bccstk817 0.145870.03126	4.66696	60	66
	bccstk485 0.051890.01123	4.61980	65	69
	bccstk138 0.002800.00050	5.53681	42	42
	mmca180 0.008080.00120	6.74795	77	76
	mccf765	0.132220.03215	4.11219	70	68
	grenoble115 0.001620.00019	8.51887	19	19
	grenoble185 0.007460.00096	7.8125	19	19
	grenoble216a 0.010550.00145	7.25	2	1
	grenoble216b 0.0105 0.00106	9.90055	19	19
	grenoble343 0.0264 0.00507	5.21053	2	1
	grenoble512 0.0588 0.0126	4.66667	2	1
	grenoble1107 0.267620.05682	4.70958	16	16
	random200 0.037 0.003	11.692	19	19
	random500 0.330 0.028	11.831	19	19
	random800 0.599 0.071	8.436	19	19
	random1000 0.934 0.111	8.452	19	19
	hilbert100 0.004140.00255	1.62264	7	289
	hilbert200 0.021740.01984	1.09552	8	567
	hilbert250 0.034810.03629	0.95942	8	714
	hilbert300 0.050930.05967	0.85350	9	847
	hilbert400 0.093070.13343	0.69756	9	1134
	hilbert600 0.214850.41759	0.51450	10	1711
	hilbert800 0.388390.94920	0.40917	10	2294
	hilbert1000 0.614251.81285	0.33883	10	2866

Table 2 :

 2 The size of the numerator n and denominator d of det(A), the size of the denominator over-approximation D/d computed by PrecDetLU and PrecMatLU; the numerator approximation App(n) obtained as sm in PrecMatDixon, and the size of the part remaining to compute. sm depends on n and the over-approximation D/d.

	Matrix	RatLU PrecDetLU PrecMatLU PrecMatDixon
	bccstk817	*	789.02	553.624	318.62
	bccstk485 278.964 143.888	95.836	57.144
	bccstk138	4.12	1.868	1.324	0.764
	mmca180	14.404	5.896	3.644	1.604
	mccf765	*	585.724	416.352	128.24
	grenoble115 1.444	0.591813	0.456	0.288
	grenoble185	5.86	2.34	1.456	0.468
	grenoble216a 1.052	0.268	0.248	0.26
	grenoble216b 10.448	3.852	2.204	2.128
	grenoble343 4.292	0.924	0.832	0.732
	grenoble512 14.844	2.868	2.48	1.072
	grenoble1107	*	698.436	519.368	367.448
	random200 24.096	10.776	3.996	2.980
	random500 432.448 180.448	71.492	54.996
	random800 1715.316 789.154	331.008	205.188
	random1000	*	1572.024	662.956	403.232
	hilbert100	17.860	0.664	0.548	0.712
	hilbert200 330.280	11.104	10.52	11.312
	hilbert300	*	59.144	65.236	66.872
	hilbert400	*	200.844	252.676	265.276
	hilbert600	*	1072.754 1664.738	1735.574
	hilbert800	*	3476.188	6299.98	8830.372
	hilbert1000	*	8870.534 18466.348	19328.66

Table 3 :

 3 Timing comparison for 4 rational determinant strategies. All times in seconds. Best times in bold.

					Rational determinant computation -timings comparison
										RatLU
		800								PrecDetLU PrecDetLUsym
										PrecMatLU
		700								PrecMatDixon Maple
		600
	Time (s)	400 500
		300
		200
		100
		0
		50	100	150	200	250	300	350	400	450	500
						n
	Figure 1: Comparison of the timings for
	the exact computation of the rational de-
	terminant of Hilbert matrices. The results
	for algorithms RatLU, PrecDetLU, Prec-
	MatLU and PrecMatDixon implemented in
	LinBox and Maple Determinant procedure
	are shown. Algorithm PrecDetLU is used
	in the classic and symmetric variant, which
	takes into account the Hankel structure of
	the matrix. All times in seconds.


Laurent Gourvès

Jérôme Monnot

Lydia Tlilane

Approximate Tradeoffs on Matroids 1

We consider problems where a solution is evaluated with a couple. Each coordinate of this couple represents an agent's utility. Due to the possible conflicts, it is unlikely that one feasible solution is optimal for both agents. Then, a natural aim is to find tradeoffs. We investigate tradeoff solutions with guarantees for the agents.The focus is on discrete problems having a matroid structure. We provide polynomial-time deterministic algorithms which achieve several guarantees and we prove that some guarantees are not possible to reach.

Introduction

This paper deals with the existence and computation of a solution s which is common to two agents. The interest of agent i ∈ {1, 2} over the set of possible solutions is captured by a utility function ui. When these functions are conflicting, it is unlikely that a feasible solution s, such that u1(s) and u2(s) are both nearly optimal, exists. So one has to make a tradeoff.

A natural way to cope with several functions is to aggregate them in a weighted sum. For example, which solution s maximizes the utilitarian function f λ (s) := λu1(s)+(1-λ)u2(s) for some λ ∈ [0, 1]? Unfortunately, this approach has two issues. The first issue is about computation: finding a solution which optimizes f λ may be puzzling when u1 and u2, though separately solvable, require completely different algorithms. The second issue is that an optimum to f λ may lead to unbalanced solutions. If s * i denotes the solution that maximizes ui(s * i ) then it is possible that a solution s, though optimal for f λ (s), satisfies u1(s)/u1(s * 1 ) ≈ 1 and u2(s)/u2(s * 2 ) ≈ 0 (or conversely u1(s)/u1(s * 1 ) ≈ 0 and u2(s)/u2(s * 2 ) ≈ 1). This pathological case indicates that s can be unfair, i.e. close to optimality for one agent, and very far from optimality for the other agent.

Then we address the following questions: For which lower bounds on u1(s)/u1(s * 1 ) and u2(s)/u2(s * 2 ) a solution s is guaranteed to exists? Which algorithm can cope with a possibly different nature of the agent's utility functions and such that non trivial a priori lower bounds on u1(s)/u1(s * 1 ) and u2(s)/u2(s * 2 ) can be derived? In this article, we seek for (α, β)-approximate algorithms, i.e. algorithms returning a solution s such that u1(s)/u1(s * 1 ) ≥ α and u2(s)/u2(s * 2 ) ≥ β for every instance. Of particular interest are the vectors (α, β) which are optimal in the sense of Pareto. There are several papers that deal with (α, β)-approximate algorithms, including [START_REF] Hansen | Bicriterion path problems[END_REF][START_REF] Warburton | Approximation of pareto optima in multiple-objective, shortest-path problems[END_REF][START_REF] Papadimitriou | On the approximability of trade-offs and optimal access of web sources[END_REF][START_REF] Stein | On the existence of schedules that are nearoptimal for both makespan and total weighted completion time[END_REF][START_REF] Goel | Simultaneous optimization via approximate majorization for concave profits or convex costs[END_REF][START_REF] Escoffier | Fair solutions for some multiagent optimization problems[END_REF].

The article is devoted to problems having a matroid structure (defined in Section 3). There is a rich literature on matroids [START_REF] Schrijver | Combinatorial Optimization: Polyhedra and Efficiency[END_REF][START_REF] Oxley | Matroid Theory[END_REF][START_REF] Korte | Combinatorial Optimization: Theory and Algorithms[END_REF]. 1 This research has been supported by the project ANR-09-BLAN-0361 "GUaranteed Efficiency for PAReto optimal solutions Determination (GUEPARD)" 2 CNRS, LAMSADE UMR 7243 3 PLS, Université Paris-Dauphine, email: firstname.lastname@dauphine.fr They model many practical situations (e.g. schedules, forests of a graph), possess many remarkable structural properties and admit elegant polynomial-time algorithms. These are motivations for studying tradeoffs on matroids.

Throughout the paper we assume that the utility of the first agent is additive while the utility of the second agent is more complex; it is a particular submodular function which captures more elaborate preferences.

A typical example of the instances covered in this article is the following. Each label ℓ has a non negative gain, denoted by g(ℓ), which captures a tourist's interest for any activity ai such that L(ai) = ℓ.

Example 1 A group of tourists is visiting

If A denotes the subset of activities actually scheduled then the travel agency would like to maximize its profit a∈A w(a). On the other hand, the tourists's viewpoint is different since they want ℓ∈{L(a):a∈A} g(ℓ) to be maximized. Here we assume that a tourist's utility increases when he does an activity of a new kind (no activity with the same label was done before).

The tourist/travel agency problem has a matroid structure usually called tranversal and explained below.

Related work and contribution

In this article, we study bicriteria approximation of labeled matroid where the utility of the two agents are in conflict. Hence, it is unlikely that the best strategy for the first agent is also the best one for the second agent, as illustrated by the tourist/travel-agency problem. One way of tackling this problem is to approximate the Pareto set ie., the set of non-dominated solutions (any improvement on one objective induces a deterioration on another objective). This point of view has been studied in the literature in several papers [START_REF] Stein | On the existence of schedules that are nearoptimal for both makespan and total weighted completion time[END_REF][START_REF] Ravi | The constrained minimum spanning tree problem (extended abstract)[END_REF][START_REF] Papadimitriou | On the approximability of trade-offs and optimal access of web sources[END_REF][START_REF] Berger | Budgeted matching and budgeted matroid intersection via the gasoline puzzle[END_REF]. When we produce a unique solution to approximate the problem, this approach is similar in the spirit to the notion of max-min fairness [START_REF] Asadpour | An approximation algorithm for max-min fair allocation of indivisible goods[END_REF][START_REF] Bertsimas | The price of fairness[END_REF][START_REF] Kalai | Other solutions to Nash's bargaining problem[END_REF][START_REF] Kumar | Fairness measures for resource allocation[END_REF][START_REF] Goel | Simultaneous optimization via approximate majorization for concave profits or convex costs[END_REF]. Fairness has been initially considered in economics and social choice theory where fairness notions based on some axiomatic characterizations such as proportional fairness, envy-freeness and max-min fairness [START_REF] Bouveret | Efficiency and envy-freeness in fair division of indivisible goods: Logical representation and complexity[END_REF][START_REF] Marsh | Equity measurement in facility location analysis: a review and framework[END_REF][START_REF] Lipton | On approximately fair allocations of indivisible goods[END_REF][START_REF] Kumar | Fairness measures for resource allocation[END_REF][START_REF] Kalai | Other solutions to Nash's bargaining problem[END_REF][START_REF] Nash | The bargaining problem[END_REF]. The goal of the max-min fairness criterion is to maximize the satisfaction of the least satisfied agent where the individual utilities of each agent is normalized in order to lie on the same scale [START_REF] Bertsimas | The price of fairness[END_REF][START_REF] Kalai | Other solutions to Nash's bargaining problem[END_REF]. An (α, β)approximation gives the satisfaction of the two agents so the satisfaction of the least satisfied agent is min{α, β}.

In this article we first recall definitions on matroids (Section 3). Section 4 presents the model studied in this paper and Section 5 gives the main definitions on multicriteria approximation.

Then we extend a result of [START_REF] Escoffier | Fair solutions for some multiagent optimization problems[END_REF] which deals with two agents having additive utility functions and willing to build a common spanning tree (a particular matroid problem). In Section 6, we propose a (1/2, 1/4)-approximation algorithm which simulates a natural process where two agents build a common solution. The general computational complexity is mentioned in Section 6 and followed by an algorithm finding a particular lexicographic optimum. In Section 8, we study a particular case, called uniform case, and show the following results: on the one hand, we produce within polynomial-time, a ( k-1 k , 1 k )-approximation in the uniform case, for any positive integer k given as the input, and on the other hand, we exhibit some instances without any (α, 1α)-approximation in the uniform case, for any α ∈ (0; 1) with α / ∈ { k-1 k : k is a positive integer }. Some open questions are indicated in Section 9.

Matroids

Matroids play an important role in combinatorial optimization and graph theory. We briefly mention some basic definitions, properties and algorithms on matroids and refer the reader to [START_REF] Schrijver | Combinatorial Optimization: Polyhedra and Efficiency[END_REF][START_REF] Korte | Combinatorial Optimization: Theory and Algorithms[END_REF][START_REF] Oxley | Matroid Theory[END_REF] for deeper expositions.

A matroid M = (X, F ) consists of a finite set of n elements X and a collection F of subsets of X such that:

(i) ∅ ∈ F , (ii) if F2 ⊆ F1 and F1 ∈ F then F2 ∈ F , (iii) for every couple F1, F2 ∈ F such that |F1| < |F2|, ∃ x ∈ F2\F1 such that F1 ∪ {x} ∈ F .
By induction (iii) is equivalent to

(iii ′ ) for every couple F1, F2 ∈ F where |F1| < |F2|, ∃A ⊆ F2\F1 with |A| = |F2| -|F1| such that F1 ∪ A ∈ F .
The elements of F are called independent, the element of 2 X \ F dependent. Inclusionwise minimal dependent sets are called circuits and inclusionwise maximal independent sets are called bases. All bases of a matroid M have the same cardinality r(M ), defined as the rank of M . Given a matroid M = (X, F ) and a subset

X ′ ⊂ X, the restriction of M to X ′ , denoted by M|X ′ , is the structure (X ′ , F ′ ) where F ′ = {F ∈ F : F ⊆ X ′ }. If X ′ ∈ F , the contraction 2 of M by X ′ , denoted M/X ′ , is the structure (X \ X ′ , F ′ ) where F ′ = {F ⊆ X \ X ′ : F ∪ X ′ ∈ F :}. It is well known that M|X ′ and M/X ′ are matroids.
Typical examples of matroids are the following:

• The forests (set of edges which do not admit a cycle) of a multigraph G form a matroid usually called the graphic matroid of G.

A base in this matroid is a spanning tree.

• Given k disjoint sets E1, . . . , E k which forms a ground set E = ∪ k i=1 Ei and k non negative integers bi (i = 1..k), the sets F ⊆ E satisfying |F ∩ Ei| ≤ bi form a matroid usually called the partition matroid.

• Given k (not necessarily disjoint) sets E1, . . . , E k , subsets of a ground set E, a partial transversal is a set T ⊆ E that an injective map Φ : T → [1..k] satisfying t ∈ E Φ(t) exists.
Then (E, T ) where T = {T ∈ 2 E : T is a partial transversal of E} is a matroid usually called the transversal matroid.

Returning to Example 1, let E = {a1, a2, ..., am} a set of activities and Ei ⊆ E the subset of activities available on day i, i = 1..n. A set of activities T ⊆ E is feasible if it exists an injective mapping Φ : T → [1.

.n] combining at most one activity of Ei per day i, i = 1..n.

A matroid is said simple if no single element, or pair of elements, is a circuit [START_REF] Oxley | Matroid Theory[END_REF]. For example, the forests of a simple graph define a simple matroid.

When every element e ∈ X has a weight w(e) ∈ R + , a typioptimization problem consists of computing a base B ∈ F that maximizes e∈B w(e). This problem is solved by the following algorithm:

Algorithm 1 GREEDY Require: M = (X, F ), w : X → R + 1: Sort X = {e1, • • • , en} such that w(ei) ≥ w(ei+1), i = 1..n -1 2: Set F = ∅ 3: for i = 1 to n do 4: if F ∪ {ei} ∈ F then 5: F ← F ∪ {ei} 6:
end if 7: end for 8: return F Note that the execution of GREEDY on a forest matroid coincides with Kruskal's algorithm for maximum weight spanning trees.

We always assume that an independence oracle can decide within polynomial time whether a set F is independent or dependent. Given two matroids (X, F1) and (X, F2) defined over the same set of elements X, there are algorithms (more elaborate than GREEDY) to solve the following problems in polynomial time [START_REF] Korte | Combinatorial Optimization: Theory and Algorithms[END_REF][START_REF] Schrijver | Combinatorial Optimization: Polyhedra and Efficiency[END_REF]:

• find an independent F ∈ F1 ∩ F2 of maximum cardinality.

• when every element e ∈ X has a weight w(e) ∈ R + , find an independent F ∈ F1 ∩ F2 that maximizes e∈F w(e).

The model: two agents on a matroid

Let M = (X, F ) be a matroid and consider the following functions:

• w : X → R + where w(x) is called the weight of x ∈ X;
• L : X → {ℓ1, . . . , ℓp} where {ℓ1, . . . , ℓp} is the set of labels and L(e) is called the label of e; • g : L → R + where g(ℓ) is the gain of label ℓ.

Note that w is additive while g has a more general form (it is a particular submodular function).

For the ease of presentation we often write g(x) instead of g(L(x)) for x ∈ X. The labels of a set X ′ ⊆ X is a set denoted by L(X ′ ) and defined as x∈X ′ {L(x)}. In the tourist/travel-agency problem, the labels {ℓ1, . . . , ℓp} are the activity types.

We study a model where F is the set of feasible solutions. The first and second agent's utilities are w and g, respectively. A set F ∈ F is then evaluated with two objective functions: w(F ) = x∈F w(x) and g(F ) = l∈L(F ) g(l) which should be maximized. An instance is then a tuple M, w, L, g .

As previously mentioned, GREEDY finds a base F1 that maximizes w(F1). Finding F2 that maximizes g(F2) can be done in polynomial time via the search for a maximum cost intersection of two matroids. The first matroid is M = (X, F ). The second matroid M ′ = (X, F ′ ) is a partition matroid defined by the labels of X: F ′ := {S ⊂ X : |L(S)| = |S|}. Every element x ∈ X has a cost c(x) defined as g(L(x)). Use any appropriate algorithm to compute a maximum cost independent in F ∩F ′ (see [START_REF] Korte | Combinatorial Optimization: Theory and Algorithms[END_REF][START_REF] Schrijver | Combinatorial Optimization: Polyhedra and Efficiency[END_REF]) and complete it into a base F2 if necessary. It is not difficult to show that the resulting base maximizes g.

Non trivial approximation

When focusing on a particular agent i, we say that an algorithm A is ρ-approximate if, for every instance, A returns a solution s satisfying ui(s)/ui(s * i ) ≥ ρ. Here s * i is a solution which maximizes ui(s * i ) and ρ ∈ [0, 1] is called the approximation ratio or performance guarantee [START_REF] Vazirani | Approximation algorithms[END_REF].

When dealing with k ≥ 2 agents, we say that an algorithm A is (ρ1, . . . , ρ k )-approximate if, for every instance, A returns a solution s satisfying ui(s)/ui(s * i ) ≥ ρi for i = 1..k. Actually the vector (u1(s * 1 ), . . . , u k (s * k )) is often called the ideal point [START_REF] Ehrgott | Multicriteria Optimization[END_REF] since it is the of an unlikely feasible solution where optimality is reached for all agents. In this paper we propose algorithms which approximate this point.

Obviously, returning a solution that maximizes w (resp. g) gives a (1, 0)-approximation (resp. 1)-approximation) but we expect α and β to be positive so that the solution constitutes a non trivial tradeoff. Next examples show that, for the general model considered in this article (described in Section 4), there is no hope for a non trivial (α, β)-approximate tradeoff (such that α > 0 and β > 0) if we consider the whole class of matroids. In both examples, for every independent F , min{w(F ), g(F )} = 0 while maxF ∈F {w(F )} = maxF ∈F {g(F )} = 1. This means that either α = 0 or β = 0 for every feasible solution. However one can overcome this issue by considering a notion which generalizes the notion of simple matroid.

Definition 1 A matroid M = (X, F ) is said labeled-simple if, for every circuit C = {x1, x2} of size two, we have L(x1) = L(x2).
Note that the matroids of Examples 2 and 3 are not labeled-simple. More generally, a simple matroid is labeled-simple and when every label appears once, these two notions coincide.

As shown in the next section, excluding matroids which are not labeled-simple makes the existence of non trivial approximation possible.

A general greedy algorithm

We propose to analyze a simple extension of GREEDY which builds a tradeoff solution. This extension, called ALT-GREEDY, simulates a simple and natural process for the construction of a tradeoff.

At the beginning F = ∅ and the agents alternatively add an element e to F such that F + e ∈ F until F becomes a base. If it is the first agent's turn then we assume that he selects e that maximizes w(F + e). In a symmetric way, the second agent chooses e that maximizes g(F + e) during his turn.

We suppose that the first agent (the one who tries to maximize w(F )) plays first.

Theorem 1 ALT-GREEDY is (1/2, 1/4)-approximate for labeledsimple matroids.

Proof. Let B = {e1, . . . , er} be the base returned by ALT-GREEDY. Let us first focus on the ratio 1/2 for the first agent. Let B * = {e * 1 , . . . , e * r } be a base with maximum weight satisfying w(e * 1 ) ≥ . . . ≥ w(e * r ). Each element of B with odd index is inserted by the first agent who wants to maximize the total weight. Take i odd. ). Since the first agent selects the element ei that maximizes w(F + ei) because this agent is additive, we deduce that w(ei) ≥ w(ê) ≥ w(e * i ), for every odd i. We get that w(

When i = 1 the current solution F is ∅ (just before ei is inserted). Otherwise F = {e1, . . . , ei-1}. Set F ′ = {e * 1 } if i = 1, oth- erwise F ′ =
r i=1 ei) ≥ w( r i=1, odd ei) ≥ w( r i=1, odd e * i ) (1)
where the non negativity of an element's weight is used. Now observe that w(e * i ) ≥ Thus, ALT-GREEDY is (1/2, •)-approximate. Now, consider the second agent. Let B g be a base with maximum gain. Let E be the elements of B g whose label does not appear in the solution returned by ALT-GREEDY: E := {e ∈ B g : L(e) / ∈ L(B)}. We suppose that the first ν = |L(E)| elements of E have distinct labels and they are sorted by non increasing gain: g(e g 1 ) ≥ g(e g 2 ) ≥ . . . ≥ g(e g ν ) where E = {e g 1 , . . . , e g ν , . . . , e g |E| }. E is in F because E ⊆ B g and B g ∈ F . Note that g(E) = g({e g 1 , . . . , e g ν }) = ν i=1 g(e g i ). Suppose that during an even step of ALT-GREEDY, the second agent could not add an element with a new label to the current solution F . We know that F ∈ F because F ⊆ B. If |F | < |E| then, by property (iii) of matroids, one could add an element with a new label to F , contradiction. We deduce that |F | ≥ |E|. Hence F has at least ν elements. Let F ′ = {e1, . . . , eν } be the first ν elements of F (following the order by which they are inserted during the algorithm). Every element with an even index, within F ′ , was inserted by the second agent. By property (iii) of a matroid, we know that g(ei) ≥ g(e g i ) holds for every even i between 2 and ν. It follows that ν i=2, even g(ei) ≥ ν i=2, even g(e g i ). The first ν elements of E being sorted by non increasing gain, we know that 2 ν i=2, even g(e g i ) ≥ ν i=1 g(e g i )g(e g 1 ) where ν i=1 g(e g i ) =

g(E). We deduce that

2g(F ′ ) ≥ 2 ν i=2, even g(ei) ≥ ν i=1 g(e g i ) -g(e g 1 ) = g(E) -g(e g 1 ) (2)
If L(e1) = L(e g 1 ) then g(B) ≥ g(e g 1 ), otherwise {e1, e g 1 } ∈ F because the matroid is labeled-simple (see Definition 1) and g(e2) ≥ g(e g 1 ) by property (iii) 3 In any case we get that

g(B) ≥ g(e g 1 ) (3)
Use F ′ ⊆ B and Inequalities ( 2) and ( 3) to derive

3g(B) ≥ 2g(F ′ ) + g(B) ≥ g(E) (4)
By definition, L(B g )\L(E) are the labels appearing in B. It follows that g(B) ≥ g(B g )g(E) that we add to Inequality (4) to get that 4g(B) ≥ g(B g ).

Thus, ALT-GREEDY is (•, 1/4)-approximate if, during an even step of ALT-GREEDY, the second agent could not add an element with a new label to the current solution F . Now suppose that for every turn of the second agent, it was possible to add an element with a new label. Rename the ξ = |L(B g )| first elements of B g such that L({e g 1 , . . . , e g ξ }) = L(B g ) and g(e 1 ) ≥ g(e g 2 ) ≥ . . . ≥ g(e g ξ ). Note that ξ i=1 g(e g i ) = g(B g ). property (iii) we have g(ei) ≥ g(e g i ) for i ≤ ξ and even. We deduce that g(B) ≥ ξ i=2, even g(ei) ≥ ξ i=2, even g(e g i ). Since the first ξ elements of B g are sorted by non increasing gain, we know that 2 ξ i=2, even g(e g i )

≥ ξ i=1 g(e g i ) -g(e 1 g ). Hence 2g(B) ≥ ξ i=1 g(e g i ) -g(e 1 g ) = g(B g ) -g(e 1
g ). Since Inequality (3) holds, we deduce that 3g(B) ≥ g(B g ) which is better than 4g(B) ≥ g(B g ). In conclusion, ALT-GREEDY is (1/2, 1/4)-approximate.

✷ Note that if we apply ALT-GREEDY by inverting the role of the agents (even steps are for the weight and odd steps are for the labels), then we can prove that ALT-GREEDY is (1/3, 1/3)-approximate for labeled-simple matroids.

Next example shows that the analysis of ALT-GREEDY is tight.

Example 4 Consider the graphic matroid of the following graph.

v1 v2 v3 v4 v5 (0 , ℓ 5 ) (0, ℓ1) (0, ℓ1) (1, ℓ1) ( 1 , ℓ2 ) ( 1 , ℓ 1 ) ( 0 , ℓ 4 ) (1, ℓ3)
Each edge e has a pair (w(e), L(e)). We suppose that g(ℓ2) = 0 and g(ℓi) = 1 for any i = 2. The tree {(v1, v5), (v2, v5), (v3, v4), (v4, v5)} has weight 4 while the tree {(v1, v4), (v2, v3), (v2, v4), (v2, v5)} uses the labels {ℓ1, ℓ3, ℓ4, ℓ5} so its gain is 4. ALT-GREEDY output a tree containing (v4, v5), (v1, v2), (v1, v5) and (v2, v3) (edges are listed by the order they enter the solution). This tree has weight 2 and its gain is 1.

General computational complexity and a particular solution

Theorem 1 is a constructive proof that every instance of the model admits a (1/2, 1/4)-approximate solution. More generally, given an instance (X, F ), w, L, g and two bounds kw and kg, what is the computational complexity of the following decision problem?

Π: Is there any F ∈ F such that w(F ) ≥ kw and g(F ) ≥ kg?

Π generalizes the minimal spanning tree problem with a side constraint which was shown NP-complete [START_REF] Aggarwal | Minimal spanning tree subject to a side constraint[END_REF], see also [START_REF] Ehrgott | On matroids with multiple objectives[END_REF] for matroids. Actually, by considering the graphic matroid where each edge e has a distinct label e, the problem dealt with in the paper is exactly the minimal spanning tree problem with a side constraint.

Though NP-complete in general, next results states that a particular Pareto optimal solution can be computed in polynomial time. Within the set of optimal solutions for the weight, let F * w,g be the one which maximizes g. F * w,g is a lexicographic Pareto optimal solution.

Theorem 2 F * w,g can be computed in polynomial time.

Before giving a proof of Theorem 2, let us give an intermediate result.

Lemma 1 Let M = (X, F ) be a matroid, B its set of bases and w : X → R + a weight function. Let Bw ⊆ B be the set of all bases which are optimal for w. Then Mw = (X, Fw) where Fw := {F ⊆ B : B ∈ Bw} is a matroid.

Due to space limitation, the proof of Lemma 1 is skipped. The proof relies on the next property which follows from results of [START_REF] Brualdi | Comments on bases in different structures[END_REF][START_REF] Schrijver | Combinatorial Optimization: Polyhedra and Efficiency[END_REF]. Since Mw and M ′ are two matroids (cf Lemma 1), this can be done in polynomial time [START_REF] Korte | Combinatorial Optimization: Theory and Algorithms[END_REF][START_REF] Schrijver | Combinatorial Optimization: Polyhedra and Efficiency[END_REF]. Complete F into a base B in a manner with elements of X sorted by non increasing weight (like in GREEDY). We claim that B = F * w,g . F ∈ Fw its completion leads to a base of maximum weight: B ∈ Bw. Observe that g( F ) = g( B) since otherwise F does not maximize g( F ). Now suppose, by contradiction, that a base B ∈ Bw is such that g( B) > g( F ). Retain exactly one element per label of B to get an independent F . Since F ⊆ B we have F ∈ Fw. Moreover g( F ) = g( B) because every label of B appears in F . It follows that g( F ) > g( F ), contradicting the optimality of F . ✷ Note that the complexity of the algorithm that finds F * w,g depends on Edmonds' Matroid Intersection Algorithm which is O(|X| 4 +Y * |X| 3 ) where Y is the complexity of the independence oracle [START_REF] Korte | Combinatorial Optimization: Theory and Algorithms[END_REF]. Y is not given explicitly, it depends on the matroid under consideration. In our study, we suppose that Y is a polynomial (see Section 3).

The uniform subcase

this section we consider a particular case, called uniform, where g(l) = 1 for all l ∈ L. In fact g(F ) = |L(F )| holds in this case so the second agent's goal is to maximize the number of distinct labels. As previously mentioned, finding a base B L that maximizes |L(B L | can be done in polynomial time [START_REF] Korte | Combinatorial Optimization: Theory and Algorithms[END_REF][START_REF] Schrijver | Combinatorial Optimization: Polyhedra and Efficiency[END_REF]. In the following L denotes |L(B L )|. In this section, we do not assume that the matroid is labeled-simple.

One can analyze ALT-GREEDY and show that it is (1/2, 1/3 + 1/r(X))-approximate for any matroid M = (X, M) in the uniform subcase. It is noteworthy that we do not need to restrict ourselves to labeled-simple matroids anymore. The proof is skipped because it follows the line of Theorem 1's proof and we are able to propose another algorithm, called 3-PHASES, with better guarantees.

Indeed, we will prove that 3-PHASES is ( k-1 k , 1 k )-approximate for every positive integer k taken as the input. So k = 2 gives a (1/2, 1/2)-approximation which improves the guarantees of ALT-GREEDY.

3-PHASES starts with an empty solution F and agent 1 adds ⌊ k-1 k L⌋ elements in a greedy manner (trying to maximize the weight). Afterwards agent 2 adds to F a set of at most ⌊ L k ⌋ new elements so that the number of new labels (ie. labels not present in the first phase) is maximized. During the last phase, agent 1 completes F in a greedy way by adding elements according to their weight. In the following, Bj denotes the elements inserted during phase j. Hence 3-PHASES returns the independent set B1 ∪ B2 ∪ B3.

The first and third phase of the algorithm are both greedy and clearly polynomial. Only the second phase the algorithm is not greedy. It consists of finding an independent of limited cardinality (at most ⌊L/k⌋ elements) in the intersection of two matroids.

first matroid is

M ′ = (X ′ , F ′ ) where X ′ := \ {x ∈ X : L(x) / ∈ L(B1)} and F ′ := {F X ′ : B1 ∪ F ∈ F }.
Actually, M ′ can be obtained by applying first the contraction of M by B1, M/B1 and then, the restriction of M/B1 to X ′ = X \ {x ∈ X : L(x) / ∈ L(B1)} (note that X ′ ⊆ X \ B1). the contraction and the restriction of matroids are matroids, M ′ is a matroid. The second matroid M2 = (X ′ , F2) is the partition matroid of X ′ induced by the labels, ie., assume that L(X) = {ℓ1, . . . , ℓp} and if Xi = {x ∈ X : L(x) = ℓi} for i = 1..p denotes the elements with label ℓi, then F2 = {F ⊆ X ′ : ∀i = 1..p, |Xi ∩ F | ≤ 1}.

Find a set S of maximum cardinality in F ∩ F ′ (this can be done in polynomial time [START_REF] Korte | Combinatorial Optimization: Theory and Algorithms[END_REF][START_REF] Schrijver | Combinatorial Optimization: Polyhedra and Efficiency[END_REF]) and if |S| > ⌊L/k⌋, retain a subset of only ⌊L/k⌋ elements. By property (ii) of a matroid, the resulting set is independent. Thus, the second phase of 3-PHASES runs in time.

Theorem 3 Let k a positive integer taken as the input, 3-PHASES is ( k-1 k , 1 k )-approximate in the uniform case.

Proof. The algorithm is clearly (0, 1)-approximate when k = 1 because no element is picked during the first phase and during the second phase one can insert ⌊L⌋ = L elements with distinct labels.

From now on, we suppose that k ≥ 2. Let B = {e1, . . . , er} be the base returned by 3-PHASES. The elements of B are numbered according to order by which they are inserted in the solution. Let us first focus on the ratio k- w(e * j ) for p = 1..k -1.

Summing up these inequalities gives

w(B1) = k-1 p=1
w(B p 1 ) ≥ (k -1)

|B 1 |+|B 2 | j=|B 1 |+1
w(e * j ).

Thus 1 k -1 w(B1) ≥ |B 1 |+|B 2 | j=|B 1 |+1
w(e * j )

Using again [START_REF] Bouveret | Efficiency and envy-freeness in fair division of indivisible goods: Logical representation and complexity[END_REF] gives

w(B1) ≥ |B 1 | j=1 w(e * j ) (7)
Adding ( 6) and ( 7), we get that

k k -1 w(B1) ≥ |B 1 |+|B 2 | j=+1 w(e * j ) (8)
Concerning the third phase, Property (iii) implies w(ei) ≥ w(e * i ) for r ≥ i ≥ |B1| + |B2| + 1. Summing up these inequalities gives w(B3) =

  Montpellier for n days and they ask a local travel agency to arrange the stay with (at most) one activity per day. There are m activities available: a1, • • • , am. Each ai has • a label L(ai) which characterizes the activity (e.g. visiting a monument, attending a show, eating at a gastronomic restaurant, etc); there are typically several activities with the same label. The labels constitute a partition of the activity set. • a list of days T (ai) during which ai can be scheduled; • a non negative weight w(ai) indicating what the travel agency earns when the activity is scheduled.

Example 2 Example 3

 23 Consider the matroid (X = {a, b}, F = {∅, {a}, {b}} where w(a) = 1, w(b) = 0, ℓ(a) = l, ℓ(b) = l ′ , g(l) = 0 and g(l ′ ) = 1. The rank is 1. Consider the matroid (X = {a, b, c}, F = {∅, {a}, {b}, {c}, {a, c}, {b, c}} where w(a) = 1, w(b) = w(c) = 0, ℓ(a) = ℓ(c) = l, ℓ(b) = l ′ , g(l) = 0 and g(l ′ ) = 1.

1 2 (

 2 w(e * i ) + w(e * i+1 )) because w(e * i ) ≥ w(e * i+1 ). It follows that w(B) = w( w(B * )/2

Property 1

 1 Let B1, B2 ∈ Bw with B1 = B2. Then, ∀e1 ∈ B1 \ B2, ∃e2 ∈ B2 \ B1 such that (B1 \ {e1}) ∪ {e2} and (B2 \ {e2}) ∪ {e1} are in Bw. Now, we are ready to give a proof of Theorem 2. Proof. Start with M = (X, F ) and define Mw = (X, Fw) as in Lemma 1. Let M ′ = (X, F ′ ) be a partition matroid where F ∈ F ′ iff F ∈ Fw and |L(F )| = |F |. Find F ∈ Fw ∩ F ′ with maximizes g( F ).

1 k 1 1

 11 for the first agent. Let B * = {e * 1 , . . . , e * r } be a base with maximum weight satisfying w(e * 1 ) ≥ . . . ≥ w(e * r ). By the third property of a matroid, we know that w(ei) ≥ w(e * i ) for i = 1..|B1| and |B1| = ⌊ k-1 k L⌋. The elements of B * being sorted by non increasing weight, we also havew(ei) ≥ w(e * j )(5)for any pair (i, j) ∈ {1, . . . , |B1|} × {|B1| + 1, . . . , |B1| + |B2|}. Note that B2 is the set of elements added during the second phase by the second agent and |B2| ≤ ⌊L/k⌋.Since |B1| = ⌊ k-1 k L⌋ ≥ (k -1)⌊ L k ⌋ ≥ (k -1)|B2|, one can partition B1 in k -1 disjoints sets B 1 1 , . . . , B k-so that every one contains at least |B2| elements. Using Inequality (5) we get that w(B p 1 ) ≥ |B 1 |+|B 2 | j=|B 1 |+1

r i=|B 1 |B 1

 11 |+|B 2 |+1 w(ei) ≥ r i=|B 1 |+|B 2 |+1 w(e *i ) that we plug in the previous inequality to get thatk k -1 w(B1) + w(B3) ≥ ) ≥ w(B * )

  {e * 1 , . . . , e * i }. By property (iii) of a matroid, and because |F ′ | > |F |, there exists an element ê ∈ F ′ \ F such that F + ê ∈ F . Using w(e

* 1 ) ≥ . . . ≥ w(e * r ), we know that w(ê) ≥ w(e * i

Actually, the contraction is defined in the literature for any X ′ ⊂ X, and when X ′ ∈ F the definition is similar to the one given in this paper.

Note that the rank is at least 2 by assumption.

Then 3-PHASES is ( k-1 k , •)-approximate. Now consider the second agent and let B L be a base of M with L labels (the maximum number). Let E L be a subset of B L containing one element per label in B L , ie. k L⌋ = ⌈L/k⌉-p elements with new labels are added during the second phase, since on the one hand, by Property (iii ′ ) ∃A ⊆ E \B1 with |A| = |E|-|B1| such that A∪B1 ∈ F (ie., A is an independent set of M ′ = (X ′ , F ′ ) where we recall that X ′ = X \ {x ∈ X : L(x) / ∈ L(B1)} and F ′ := {F ⊂ X ′ : B1 ∪ F ∈ F }) and on the other hand, A is an independent set of the partition matroid

Thus, A is a feasible solution of the intersection of matroids M ′ and M2. So,

⌈L/k⌉ ≥ L/k labels appears at the end of the second phase. The elements added during the last phase can only increase the number of labels so we know that in every case 3-PHASES is (•, 1/k)approximate.

In conclusion, 3-PHASES is ( k-1 k , 1 k )-approximate. ✷ Next example shows that the guarantees of 3-PHASES are essentially Pareto optimal.

Lemma 2 Let α ∈ (0; 1) and β ∈ (0; 1). There are instances without any (α, β)-approximation in the uniform case in the following cases:

Proof. We prove these results for a particular matroid, the graphic matroid. Let α > 0 and β > 0.

For (a). We prove that the result holds for α + β = 1, and then a fortiori for α+β ≥ 1. So, assume α+β = 1 and α = k-1 k for every positive integer k. Thus, there is a unique integer L ≥ 2 such

and (1α)L > 1. Consider the following multigraph GL, instance of the uniform case for the graphic matroid.

The ideal point of GL is x = (L, L) by considering trees given by the bottom edges for the weight and the top edges for the labels.

By contradiction, assume that T is a spanning tree of GL with weight w(T ≥ αL and number of labels |L(T )| ≥ (1α)L. If T at least L -1 edges of weight L/(L -1), then T contains one label. Hence, 1 = |L(T )| ≥ (1 α)L > 1, contradiction. Thus, T contains at most L -2 edges of weight L/(L -1) and then L(L -2)/(L -1) ≥ w(T ) ≥ αL > L(L -2)/(L -1), contradiction.

For (b). Consider the multigraph G2, ie., L = 2 described previously. If β > 1/2, then |L(T )| ≥ 2β > 1. So, |L(T )| = 2 and T is given by the top edges of G2. In this case, w(T ) = 0. ✷ Item (b) of Lemma 2 shows that β = 1/2 is a tradeoff on the number of labels that we can reach if some positive guarantee on the weight is achieved. Hence, Theorem 3 gives the best results that we can hope.

Conclusion and future directions

For the general case, we have proposed a polynomial-time (1/2, 1/4)-approximation for labeled-simple matroids (a (1/3, 1/3)-approximation exists also by inverting the agents' role). An important question consists of improving these approximations. We believe that a (1/2, 1/3)-approximation exists and that it offers Pareto optimal guarantees, ie., there are small instances without any (α, β)-approximation with α ≥ 1/2 and β ≥ 1/3.

For the uniform subcase we proposed polynomial-time deterministic algorithms which achieve several tradeoffs and we prove that some tradeoffs are not possible. Note that the complexity result of Section 7 does not hold for the uniform case so the exact complexity of the uniform case is open. A first step in this direction would be to consider the open problem of computing F * g,w (a Pareto optimal solution with maximum weight among solutions of maximum gain). Another open question is the case where both agents' utility functions are of type "g" (not only the second agent's utility function).