Christophe Alias Prism
email: christophe.alias@prism.uvsq.fr

Denis Barthou Prism

Algorithm Recognition based on Demand-Driven Dataflow Analysis

Algorithm recognition is an important problem in program analysis, optimization, and more particularly program comprehension. Basically, one would like to submit a piece of code, and get an answer like "Lines 10 to 23 are an implementation of Gauss-Jordan elimination". Existing approaches often perform a bottom-up recognition, using a database describing many implementation variations of a given algorithm, in a format making them difficult to maintain.

In this paper, we present a new method to recognize algorithm templates in a program. We perform a top-down recognition, stopping at template variables and thus limiting the need for many program variations. The code to analyze is given in SSA form, and the -functions are computed on-demand, when computation is possible, along the candidate slices. Criteria Wills DiMartino Metzger Us Cost high high low high Adding patterns is... difficult difficult easy easy

Introduction

Algorithm recognition is an old problem in computer science. Basically, one would like to submit a piece of code to an analyzer, and get answers like "Lines 10 to 23 are an implementation of Gaussian elimination". Such a facility would enable many important techniques: ¡ Program comprehension and reverse engineering: we try to recognize in the program some codes that have a significant meaning.

¡

Program optimization: if we have the necessary items in our library, we may replace lines 10 to 23 by a hand optimized version, or by a sparse version, or a parallel version. If we are bold enough, we may even replace the relevant part of the code by a completely different implementation, as for instance an iterative solver.

¡

Program verification: if we know that the program specification asks for Gaussian elimination and the analyzer does not find it, we may suspect an error.

¡

Hardware-software codesign: if we recognize in the source program a piece of code for which we have a hardware implementation (e.g. as a coprocessor or an Intellectual Property component) we can remove the code and replace it by an activation of the hardware.

Simple cases of algorithm recognition have already been solved, mostly using pattern matching as the basic technique. An example is reduction recognition, which is included in many parallelizing compilers. A reduction is the application of an associative commutative operator to a data set. It can be detected by normalizing the input program, then matching it with a set of patterns which should include the most common associative operators (addition, multiplication, and, or, max, min ...). See [START_REF] Redon | Detection of scans in the polytope model[END_REF] and its references. One issue of pattern matching is that it does not handle well the variations in the control or data flow of the program. Several complex pattern matching techniques have been proposed (see the recent book by Metzger and Wen [START_REF] Metzger | Automatic Algorithm Recognition: A New Approach to Program Optimization[END_REF] and its references) in order to tackle this issue. An alternative approach, explored by [START_REF] Barthou | On the equivalence of two systems of affine recurrence equations[END_REF][START_REF] Barthou | On the equivalence of two systems of affine recurrence equations[END_REF] and [START_REF] Shashidhar | Geometric model checking: An automatic verification technique for loop and data reuse transformations[END_REF] is to consider system of affine recurrences as the starting point of the algorithm recognition. From this normal form the method described in [START_REF] Barthou | On the equivalence of two systems of affine recurrence equations[END_REF] is able to find the equivalence of two programs, modulo transformations such as variable hoisting, data expansion/shrinking, affine transformations of the iteration domain, or common sub-expression optimizations.

All these methods recognize only algorithms that have exactly the same semantics as the code they match. Many algorithms however are better described in generic terms, abstracting away the details of implementation. For instance, Gaussian elimination is one instance of the wellknown algebraic path problem (APP), as the Warshall's transitive closure algorithm and Floyd's shortest path algorithm are also instances of this same APP. The only difference is the underlying algebraic structure. The only way to handle them by the previous methods is to consider one different pattern for each instantiation. Such generic algorithms are called algorithm templates and many efficient implementations of templates have been proposed. See [START_REF] Siek | The matrix template library: A generic programming approach to high performance numerical linear algebra[END_REF] for matrix manipulations, [START_REF] Lee | The generic graph component library[END_REF] for graph algorithms or [START_REF] Tayoudjameni | Derivation of systolic algorithm path problem by recurrence transformations[END_REF] for the APP, to name a few.

The aim of this paper is to present a method for algorithm template recognition, based on the preliminary works described in [START_REF] Alias | On the recognition of algorithm templates[END_REF][START_REF] Barthou | On the equivalence of two systems of affine recurrence equations[END_REF], enlarging the scope of the technique in [START_REF] Alias | On the recognition of algorithm templates[END_REF] from static control programs [START_REF] Feautrier | Dataflow analysis of scalar and array references[END_REF] to any program. Given a program written in SSA form [START_REF] Cytron | Efficiently computing static single assignment form and the control dependence graph[END_REF] and an algorithm template, our goal is to find out, if it exists, an instantiation of the template corresponding to the program. The -functions of the SSA form capture some dataflow information that may be needed for the matching to succeed. The on-demand evaluation of the -functions during the matching procedure provides a convenient way to combine the power of the previous analysis to a larger range of applications.

As in most algorithm recognition methods, the first step is to normalize the given program as much as possible. One candidate for such a normalization is the conversion to a System of Recurrence Equations (SRE). Algorithm templates are given as static control programs, whereas the code to analyze is in SSA form. It has been shown that static control programs [START_REF] Feautrier | Dataflow analysis of scalar and array references[END_REF] can be automatically converted to Systems of Affine Recurrence Equations (SARE), and such a conversion already was the first step in [START_REF] Redon | Detection of scans in the polytope model[END_REF]. As for codes in SSA form, we propose a way to encode it into an SRE. The next step is to design an equivalence test between SREs and SARE templates. This is the main theme of this paper.

Section 2 introduces some essential definitions about SREs and the transformation procedure from SSA form to SRE. Section 3 defines the problem at hand and Sections 4 and 5 build the semi-algorithm performing the match between the template and the code. Comparison with related works is performed in Section 7.

Background

We present the definition of systems of recurrence equations which are used in the rest of the paper, as well as how to transform a program in SSA form into an SRE.

The code in Figure 1.a provides a toy example for the notions presented in this section. This is a computation of the bit vector representing the local maxima in the input array I.

Systems of Recurrence Equations

Systems of recurrence equations are a convenient way to represent algorithms: they already eliminate some syntactic aspects of the programs since they represent the computations with explicit dataflow information. Moreover, systems of affine recurrence equations, SAREs, can be obtained from static control programs by reaching definition analysis [START_REF] Feautrier | Dataflow analysis of scalar and array references[END_REF]. The basic reference on SAREs is [START_REF] Darte | Scheduling and automatic Parallelization[END_REF] ¡ £¢ ¤ ¦ ¨ ) " 6 0 %¢ 7 21 ) " 6( 3 %1 2 4¨" 6¢ 7 %7 54 6£ for all dependences in the SARE. If the domains are bounded, a schedule exists iff the given SARE has no dependence cycle. The scheduling problem for parametrically bounded SAREs is undecidable [START_REF] Saouter | Computability of recurrence equations[END_REF]. However, the existence of affine schedules for SAREs is decidable [START_REF] Feautrier | Some efficient solutions to the affine scheduling problem, II, multidimensional time[END_REF]. Note that in general, these schedules have a parametric latency. We only consider in this paper SAREs with a schedule.

¡ ¢ ¢ 4 ¤£ ¦¥ ¢ £ ¨ § " © A © A ¢ 7 ¢ ! #" ¦¥ ¢ £ ¨ § " 6 ¦¥ ¢ $ %£ & & '© A ¢ 7 ( G H ¥ " 3
A SARE template has the same definition as a SARE, except that in the definition of the clauses, We assume that " 71 8£ .

From SSA to SRE

Consider first a program in SSA form, where only scalars appear in the -functions. We assume in the following that lower and upper bounds of do loops are static (they depend on surrounding loop counters and on parameters). The procedure builds for each assignment in the program a clause of an SRE.

Perform scalar expansion. When a scalar, s, has an upward-exposed use in a loop, rewrite this use into

M0 := I[0] S0 := 0 do i := 1, n S1 := (S0,S4) M1 := (M0,M3) if (M1 < I[i]) then M2 := I[i] S2 := S1+1 endif S3 := (S1,S2) M3 := (M1,M2) S4 := S3*2 enddo OUTPUT := (S0,S4[n])

Figure 3. SSA form of the code s[prec(i)]

where i is the iteration vector of the use statement. prec is a function providing the previous value of this vector, in the execution order. The defining domains of the SRE built by this procedure are polyhedra, as in SAREs. However, the dependence functions may be non affine.

Perform an if-convertion on the program.

For each assignment S[i] := t in the program, consider the intersection of the surrounding loop domains with the predicates of the if. This set may be defined by non-affine constraints in the loop counters and parameters. Let

Fig. 3 describes the SSA form of the code in Fig. 1.a and Fig. 4 shows its represention as a system of recurrence equations, as given by the previous algorithm.

The SSA form does not take into account the indices of the arrays. In order to transform a statement of the form A3 := (A1,A2) with A1, A2 and A3 arrays, into a clause of an SRE, instead of applying the previous steps, we compute the source of all elements fe g B , for all indices g

, with an instance-wise dataflow analysis such as [START_REF] Feautrier | Dataflow analysis of scalar and array references[END_REF]. The code must then be a static control program and there is no more the advantages of an on-demand dataflow analysis.

¡

¢ £ ¢ ! #" @ ¢¡ ¢ 3 £¡ " @ £¤ A @ ¦¥ ¢ 5$ ¨£ 7 £ ¢ ! #" @ ¥ B ¢ 3 B ¢ § © ¦ © " &@ ¢¡ ¢ X4 6£ 7 £ ¢ ! #" @ e ¢ 3 ¥ " @ ¡ ¢ & @ ¥ ¢ 7 £ ¢ ! #" @ ¥ ¢ 3 @ e ¢ X ¤ ! £ £ ¢ ! #" ¡ ¢ £ e " ¤ e ¢ 5$ ¨£ 7 £ ¢ ! #" ¥ ¢ £ B ¢ § ¨© ¦ © " " ¢ 7 £ ¢ ! #" e ¢ £ ¥ " ¡ ¢ & ¥ ¢ 7 G $# &% B" @ ¤ @ ¥ " 3 7

Matching problem

We want to check whether a program ' is an instance of a template ( , where ( is a static control program, and ' is any program with scalar variables. A template is con- sidered as a program with some redefinable functions and inputs, as templates in C++ operating on parametric types.

For a program ' in SRE form, ' 0) denotes the term ob- tained by completely unfolding the recurrence of ' . This term is well-defined for a program ' with a schedule (no infinite loop). Likewise, for a template ( with a schedule, the term ( 1) exists. The matching problem is to find a substitution 2 of the template variables such that 2 " ( 3) 7 is syntactically equal to ' 4) . Such a substitution, if it exists, is called a unifier of the matching problem between ( and ' .

This problem depends clearly on the underlying algebra. It is clear, however, that equivalence in the Herbrand universe implies equivalence in all conforming algebras. We only consider in this paper equivalence in the initial algebra (no semantics). In particular -functions are not given any interpretation during the resolution of the matching problem.

The matching problem between two SAREs has been proved undecidable in [START_REF] Alias | On the recognition of algorithm templates[END_REF], therefore the matching problem between a SARE template and a SRE, which is at least as complex, is also undecidable.

Principle of the algorithm

In order to give the intuition of the method, we study the case of the template and code in Fig. 1.a when "

. The generalization will be described in the following section.

We can completely unfold the recurrence of the template:

5 76 98 A@ CB D@ EB DF HG I P RQ SF TG VU WP YX ¦ba cQ dF HG a P X &èa ).

To unify two terms, we can apply Huet's procedure [START_REF] Huet | A unification algorithm for typed -calculus[END_REF]. Basically, we try to decompose the problem by applying the rule:

C" gf E %7 ih C" fE Sj 7 lk E ¡ mh E j ¡ ¢n $ &$ '$ n E So !h E j o until terms of the form § " Wf E %7 9h C" fE Sj 7
are obtained, where § is a template variable. Then we try to construct an unifier by trying:

¡ Projections: § " %$ &$ &$ G 3¨$ &$ &$ 97 ¥ 6G £¨; ¡ Imitation: § " f G 37 C" § ¡ f G $ '$ &$ § o f G 7 , with the § new function variables.
When decomposition is not possible because the two head operators are different, then this is a failure, the program and template are not equivalent. When one of these operators is function then this function may be computed in order to avoid the failure of the method. If one of the variable in the -function is never in dependence with the variable defined by the -function, then a dependence analysis can prove it. Otherwise, more complex dataflow analyses are necessary [START_REF] Feautrier | Dataflow analysis of scalar and array references[END_REF][START_REF] Pugh | The Omega test: a fast and practical integer programming algorithm for dependence analysis[END_REF][START_REF] Maydan | Array dataflow analysis and its use in array privatization[END_REF]. The -function cannot be evaluated if the code containing it is outside the scope of the dataflow analysis. In this case, the -function remains unevaluated and the unification fails, the template is then said non-equivalent to the program. This is a safe approximation. Likewise, B

-functions could be evaluated if there were a dataflow analysis able to handle in an exact way non-affine constraints.

On the example, we try to decompose the head operator in both expressions: for the code, this is £% , for the template this is . Thus ¦% is computed: as we have supposed that " 1 , a dependence analysis can prove that G $# does not depend on its first argument, £% B" G 5 qp 7 § rp . This is equivalent to prove that at least one iteration of the loop is executed. The head operator is now in both terms.

By decomposition, we have to match

§ " F § " $ &$ '$ 7 © A 7 3h ¥ " ¡ " $ '$ &$ 97 B ¥" $ &$ '$ 7 )7
. By applying the rule imitate, we define § " FG qp 7 as ¥ " § ¡ " G 5 qp 7 § ¥ " G 5 qp 7 )7

where § ¡ and § ¥ are new template variables and the next step is a decomposition for ¥ . Proceeding in the same manner, we find out the following solution: § " FG qp 7 ¥ ¥ " £¡ B" G d B7 B ¥" ¦¡ B" G 7 %7 54 6£ 7 and © A ! . This is one of the possible solutions which can be obtained by an exhaustive application of the proce-dure. This solution is not very helpful as it is. The semantics of the functions is:

¡ " FG 5 p 7 H P ¢ ¡ % Q TS VU XW mp U Xa bU G ¥ B" G 5 qp 7 H P ¡ ¢ ¢ ¢ Q TS VU W p #U Xa bU G $ Therefore, the function § of the template of Fig.1.b is, in a more natural way: § " G 5 qp 7 ¥ £ £ £ £ £ £ £ £ £ £ £ £ £ £ £ £ £ £ H P ¡ ¢ ¢¤ £ ¢ Q TS VU XW £ £ £ £ £ £ £ £ £ £ H P ¡ ! ¢ ¢¤ £ ¢ Q TS VU XW ¥ ¦ £ £ £ £ £ £ H P ¢ § ¨ Q bS U XW G d U Xa bU © 4 6£ U Xa bU 0c U Xa bU £ £ £ £ £ £ H P ¢ ¡ % Q TS VU W YG d U `a TU 0 Computation of ¡ , used by §
, can be added likewise. In fact, as templates and programs can depend on a parameter (here n), computed terms can have a parametric length, making Huet's algorithm not applicable as it is. Let us try nevertheless to build a unification tree. We start with state G H h G # , then we unfold the recurrence, applying on demand Huet's rules. We then obtain an unification tree -with a parametric number of states -which gives us the solution. This tree as many similar nodes, and as in [START_REF] Barthou | On the equivalence of two systems of affine recurrence equations[END_REF], the idea is to capture this unification tree into an automaton with a finite number of states and then analyze it to construct the final set of solutions.

Our algorithm

As said above, the idea of this semi-algorithm is to implement the procedure with an automaton and to analyze the automaton, without executing it, in order to construct the set of unifiers. The automaton, a Memory State Automaton (MSA) is described below.

Memory State Automaton

Definition

The state of an MSA has two parts: an element of a finite set and a vector of integers. The vector associated to state is denoted f and the full state is f . The dimension of tf is determined by and is noted " f . f . There is an edge from to in an MSA iff f (' 0)

. Let 1 ¤ f 32 be the initial state of the automaton. A state f is reachable iff there exists a finite sequence of tran- sitions from the initial state to f 4 : 5 76 w Q 98 98 @8 gQ 6 BA Q DC FE ¨Q 98 F8 98 gQ GC 7E 9H PI B 6 BA 8

RQ TS C 7E 9U WV ¨Q DC FE FU YX acb dE FU eV 9f E FU qX g8

The reachable set of , noted f , is the set of vectors f such that f 4 is reachable from the initial state.

Computing the Reachability Relation

One method for computing the reachability relation consists in characterizing all possible paths in the MSA, then computing the relation associated to each path and "summing" the results. This can be done by associating a letter from a new alphabet to each edge of the MSA. This results in a finite state automaton on the given alphabet. Familiar algorithms [START_REF] Autebert | Context-free languages and push-down automata[END_REF]

Construction of the matching MSA

Let us consider ( h ' a matching problem. We explain in this section how to build the corresponding MSA, by first describing the general form of a state, and then the different kind of transitions between two states. Next section gives an algorithm to analyze the matching MSA in order to compute the set of unifiers.

States

Each state of our matching MSA has two part: a clause 2 !E h E j with 2 a substitution, and a vector of integers f , which is the concatenation of f h and f h j

, where f h is the lhs iteration vector, and f h j is the rhs iteration vector. The initial state is

i G H f h !h G # f h j
, where G H and G # are corresponding outputs of ( and ' . Its reaching set is

p f h f h j q f h f h j &r
. The final states are either:

¡ © 0 f h h E j where ©
is an input variable of the template

( ; ¡ c
if there is no unifier.

Transitions

We describe thereafter the transitions:

¡ Decompose: From a state with label:

2

C" Wf E " f h 7 %7 h j " fE j " f h j 7 )7

starts a transition to each state:

2 E ¨" f h 7 h E j ¨" f h j 7
with the firing relation:

i p f h k f h f h j k f h j &r provided j
. All these transitions constitute an and-branching. If

2 C" f E " f h 7 )7 h 4" fE Sj " f h j 7 %7
If is not computable then start a transition to c

. Otherwise compute its value (for instance with the methods in [START_REF] Feautrier | Dataflow analysis of scalar and array references[END_REF][START_REF] Maydan | Array dataflow analysis and its use in array privatization[END_REF], which is of the form:

¡ f G 4$ G ¨¡ f h j 4¤ ¦
p f h k f h f h j k f h j f h j ¤ ¦ ¨r
All these transitions constitute an ¢ " i -branching.

¡

Generalize: From a state with label:

2 § 1 C" f h 7 h
6E j create a transition to state:

2 § f h h 6E j
with the firing relation:

p f h k 1 4" f h 7 f h j k f h j Wr
There is a similar rule for the rhs. This rule is important because it normalizes the form of the states, ensuring a finite number of states [START_REF] Alias | On the recognition of algorithm templates[END_REF].

¡ Compute: If § f h E ¨¡ f h ¤ ¦ ¨,
p f h k f h f h j k f h j f h ¤ ¦ ¨r
All these transitions constitute an ¢ " i -branching. There is a similar rule for the rhs. One can remark that the -functions are only computed when they are needed by the Decompose rule. Therefore, this limitates the cost due to dataflow analysis. Moreover, a -function can be imitated as it is. Hence, it allows us to handle programs with non-computable -functions.

Analysis of the matching MSA

We have now to analyze the matching MSA in order to decide whether the program is an instantiation of the template, and to find out the set of unifiers. This can be done by the following algorithm:

Algorithm 2 Match SRE IN : A matching problem ( h ' .

OUT : A set of unifiers p @ ¤£ ¦¥ "¡ $ &$ '$ '@ ¤£ ¦¥ Do r

E . -add § k E #¡ C" f h f h j 7 ¤ to @ ¤£ ¦¥ ¡
If there is no functional incoherence in @ ¤£ ¦¥ Then emit( @ ¤£ ¦¥

);

Until they is no more £ -branching to enable.

Step 2 allows us to eliminate unaccessible nodes due to systematic application of Compute and Compute rules.

Step 3 avoids unnecessary work by eliminating the nodes which do not produce unifiers. Step 4 studies all possible combinations of unifiers in the remaining automaton. For each definition of a -function it collects all couples " f h f h j 7 to which the definition is applied. The final solution is discarded if a variable § have two distinct definitions for a value of f h (functional incoherence).

Rule Compute does not change the soundness result proved in [START_REF] Alias | On the recognition of algorithm templates[END_REF]. However, our algorithm is unable to find functions with a recursive definition, so it is not complete. Another restriction is the construction of the transitive closure of a relation which occurs when a reaching set is computed. This is not an effective procedure [START_REF] Kelly | Transitive closure of infinite graphs and its applications[END_REF] and the algorithm works only when transitive closures are computable.

An example

Let us apply our algorithm to the matching problem between the template of Figure 1.b and the program of Figure 1.a. Applying on demand the rules described in Section 5.2, we obtain an MSA of 72 states. For presentation reasons, we will just describe a small, but significant part of the MSA. Figure 5 give the first part of the MSA. For sake of clarity, arcs are labeled by shortened no- " 3 3h % " @ ¤ '@ ¥ " 3 7

¢¡ ¤£ ¥ ¢¦ D © o £ ¥ ¨ § © 2 © o D © £ ¥ § © 2 ! © o " # D © %$ ¡ '& © " )( ¥ £ ¥ ¨ § © 2 © o " # D © 0$ ¡ 1& © " )( ¥ £ ¥ ¨ 2 . . . # D © )$ ¡ '& © " )( ¥ £ ¥ ¨ 2 © o # D © 0$ ¡ 1& © 3( ¥ £ ¥ ¨ © 4 # D © )$ ¡ 5& © " )( ¥ £ ¥ ¨ 26 © 74 8( ¥ # D © %$ ¡ '& © " £ ¥ ¨ 26 © 74 . . . U 79 H U A@ 0B H 2C 2 H @ Ü 4 9 H U 4 @ DFE 'G IH 'P ! H 1Q 'E
, label

¢ "

indicates a generalization represented by the relation p " %$ '$ 7 ik " 6¢ &$ 97 q ¢ " r .

¢

is computed with respect to its definition in the SARE template. We will only consider ¢ 01 , which represents relation p " ¢ %¢ j 7 k " ¢ %¢ j 7 q ¢ d1 r

. Because § " ¢ $ C£ & & '© A ¢ 7 h % " @ ¦¤ '@ £¥ " 3 7 is a ¢ 'R ¢ Di - ¢ 5R 5¢ Di
pair (the head symbol of lhs is the unary function

G k G 3
), % must be computed. This explains the two next transitions. @ £¥ " 3 is generalized in @ £¥ ¢ j

, then computed, and gives @ e ¢ j !

. The two terms are decomposed, then @ e ¢ j is computed and gives the

¥ 5S bG ¢ UT ¥ IS - ¢ 5R 5¢ Di pair § " ¢ $ 6£ & & '© A ¢ 7 h ¥ " @ ¡ ¢ j & @ ¥ ¢ j 7
. We can now try to build a value for § by applying Huet's rules. Figure 6 gives the part of the MSA obtained while try- p " F" " 7 9r so the following substitution is added to @ £ ¦¥ : c ed f hg 2i 2p 0q p sr xq 7t s It Vu D wv yx q s ' '

ing § k G p 3$ G . Since ¥ is not computable, we ob- D © 0$ ¡ £ ¥ ¨ § B V Ẅ© 4 B © 74 D © £ ¥ ¨ § B g© 4 B © 4 " # & © ¤ '& © ¡ " )( ¥ £ ¥ § B Ẅ© 4 B © 4 # D © )$ ¡ '& © " )( ¥ £ ¥ § B g© 4 B © 4 " W W U VX Y U eV Ü79 ¨U `@ 0B
In the same way, the reaching set of the final state © ¢ 9h @ e ¢ j is p " " " 7 9r so the following substitution is added to @ ¤£ ¦¥ : g t Y 0d f y t u V 8r xq t s It u wv x q 7 s 5 5

Finally, the following solution is emitted: c d f hg 2i 2p 0q p r xq 7t s It u v x q s ' ' g t Y %d f y p t u r xq 7t s It u v x q s ' '

The solution emitted when the imitation arc is enabled is: c d f hg 2i 2p %q h vu q 7i l £ ts xq 7i £ 9 r xq 7t s It u v x q " s q )q Uq Uq q 7 s I 5 g t i d 0d f r W r xq 7t s It u v x q " s q 5

Which corresponds to the solution given in section 4.

Related works

It is a well-known fact that programmers spend most of their time to understand programs. Hence, tools which facilitate program comprehension are important. They can be divided in two complementary categories: The first category includes software visualization tools (see for example [START_REF] Telea | An open visualization toolkit for reverse architecting[END_REF] or [START_REF] Quigley | Experience with FADE for the visualization and abstraction of software views[END_REF]) which allows the user to understand software architecture by providing a graphical presentation of modules dependencies and searching facilities. Pintzger's revealer [START_REF] Pinzger | Revealer: A lexical pattern matcher for architecture recovery[END_REF] extracts parts of the source code which corresponds to a given pattern by performing a lexical analysis. The user can then quickly find parts of the code containing characteristic statements, like for example the declaration of a Java server socket.

We are interested in automatic comprehension of program. But what does it really mean to "understand a program"? Biggerstaff et al. [START_REF] Biggerstaff | The concept assignment problem in program understanding[END_REF] defined the concept assignment problem as "[...] a process of recognizing concepts within a computer program and building up an 'understanding' of the program by relating recognized concepts to portions of the program, its operational context and to one another.". Since semantic equivalence is not decidable, the concept assignment problem which is at least as difficult, is not decidable. There exists several approaches to recognize a concept in a program. We will evaluate these with respect to the following criteria: The last one is the most important. Indeed, the main difficulties come from various codes implementing the same algorithm. The different kinds of variation that can be found in a program are mentioned in [START_REF] Metzger | Automatic Algorithm Recognition: A New Approach to Program Optimization[END_REF]. To sum them up, they are: ¡ Organization variation: any permutation of independent statements and introduction of temporary variables ; ¡ Data structure variation: the same computation with a different data structure ; ¡ Control variation: any loop transformation, e.g. fusion, splitting, skewing, . . . or control transformation e.g. if-conversion or dead-code suppression; ¡ Semantic variation: pieces of code equivalent modulo a theory defining properties on operators, e.g. associativity/commutativity of 4 . Wills [START_REF] Wills | Automated Program Recognition by Graph Parsing[END_REF] represents programs by a particular kind of dependence graph called flow-graphs, and patterns by flowgraph grammar rules. The recognition is performed by parsing the program's graph according to the grammar rules. We finally obtain a parsing tree which represents a hierarchical description of a plausible project of the program. This approach is a pure bottom-up code-driven analysis based on exact graph matching. Wills has shown that the flow-graph recognition is NP-complete, and argues that even if the cost of her algorithm is exponential in the worst case, it is feasible to apply it to practical partial program recognition. Patterns are represented by grammars rules, encoding a hierarchy among them, but making the pattern base difficult to maintain. Organization variation is partially supported. Temporary variables can be handled by adding specific rules. All others variations can be handled only if they are explicitly described in the pattern base.

Di Martino et al. [START_REF] Martino | PAP recognizer: A tool for automatic recognition of parallelizable patterns[END_REF] propose an approach similar to Wills'. Flow graphs and patterns are coded by prolog clauses, then the recognition is performed by SLDresolution. The differences are essentially in the abstract program representation (flow graph).

Metzger [START_REF] Metzger | Automatic Algorithm Recognition: A New Approach to Program Optimization[END_REF] normalizes the program's AST by applying various heuristics, and then compares it to the pattern AST. Of course this approach is scalable and low cost. Moreover, the patterns are given in the source form. This cannot handle much program variations though, in particular data structure variations. Organization variations are not handled by the algorithm itself, but by pre-treatments applied to the program. In the same way, the control variations supported are bounded to pre-treatments. One can remark that semantic variations are partly taken into account, with commutativity of usual operators.

Our approach is able to handle all variations described, except for semantic variation. Organization and data structure variations are normalized during the conversion to the SRE. Since our algorithm checks on-demand the order in which operations appear in the template and the program, it copes with any transformation preserving the order of the operations, in particular the control variations. Templates are given as source code, so new patterns can be added to the base without difficulties. However, our algorithm is expensive because of the reaching set computations and dataflow analyses. The evaluation of different approaches is summarized in Figure 7.

Conclusion and Future work

Algorithm templates represent programming models that convey genericity, portability, that can be easily customized by the programmer to suit its need and at the same time have efficient implementations. Algorithm template recognition thus appears as a promising tool for code comprehension, validation and optimization. In this paper, we have presented an approach that provides such recognition for templates described by systems of affine recurrent equations that can be applied on any code. As a consequence, our analysis is able to recognize algorithms obtained by composition of other algorithms, since templates can be composed with other templates. While other analyses [START_REF] Wills | Automated Program Recognition by Graph Parsing[END_REF] could recognize an algorithm made of several known algorithms, ours works also for unknown algorithms.

The method presented relies on a precise knowledge of the dataflow in the program. The on-demand dataflow analysis approach reduces the piece of code that must be in the scope of an exact dataflow analysis to the portion only needed by the template. Non static control programs can thus be analyzed with this method, provided the difficult parts of the program correspond to the template variables.

In future work, we will investigate the feasibility of the approach on benchmark applications, with respect to the assumptions that have been made and by extending the existing prototype developed for the equivalence of SAREs. In order to make algorithm template recognition more efficient, several further developpements can be explored: faster heuristics, use of slicing techniques in order to target only good candidates for template recognition, and recognition of templates organized in a hierachical way: the template described in this paper is an instance of a template of reduction. This could drastically reduce the cost of the recognition with respect to the number of templates considered.

Fig. 1 .Figure 1 .

 11 Figure 1. (a). Code to analyze (b). Template

  The example of Fig.2illustrates the transformation from the program of Fig.1.b to a SARE. The output, G IH is set to the

Figure 2 .

 2 Figure 2. SARE template

  defined during the matching procedure. The SARE template representing the template of Fig.1.b is given in Fig.2.1.

Algorithm 1

 1 SSA to SRE IN : A program in SSA form. OUT : A SRE 1. Give a different number to each -function 2. For each while loop, add a loop counter going from £ to an upper bound denoted 1 C" ¢ 7 where ¢ is the iteration vector of the surrounding loops.

¦

  denote the non-affine part of this set and 9 the affine part. We define the clause ¢ ¤ 9 A@ ¢ 3 CB (D " FE %7where B 5D is a new function symbol. Its semantics is B (D " FG 37 IH P ¢ ¤ §¦ RQ TS VU XW YG 6U Xa bU dc where c denotes an undefined value.

Figure 4 .

 4 Figure 4. SRE form of the code

A

  transition in an MSA has three elements: a start state, , an arrival state , and a firing relation ¢f in o "! $# o % . A transition from f to & can occur only if & f ¤

  Consider a state with label:

¡Proposition 1

 1 Huet's rules produce an or-branching between each Project and Imitate. The firing relation is i since they do not modify the index variables. ¡ Substitute: Simply replace an already defined template variable by its value. Firing relation is i . Rule Compute does not change the following result proved in[START_REF] Alias | On the recognition of algorithm templates[END_REF]: Let ( h ' be a matching problem, and be its corresponding MSA. Then the number of states of is finite.

3 .--

 3 Propagate c by applying rules c n G G n c c and c ¨ § G G © § c 6G Initialize the set with the reaching set of the definition node; Visit the following nodes, then add to the reaching sets of the nodes which precede the substitution of § by

Figure 5 .

 5 Figure 5. First part of the example MSA

Figure 6 .

 6 Figure 6. Part of MSA obtained while trying § k G p 3$ G tain failures states ( c ) which are final states. While trying § k G p 3$ p, we obtain © ¢ h @ e ¢ j, which is a final state. To achieve the construction of the MSA, we should try imitation of rhs head symbol (here ¥) by try-

¡

  Tools which help the user to understand a program ; ¡ Tools which try to understand automatically parts of a program by performing subcomponent extraction.

¡

  Maintainability of pattern base. Can we easily add our own patterns ?¡ Program variations handled.

Figure 7 .

 7 Figure 7. Evaluation of different approaches

index variables in the equations are called pa- rameters of the SRE; Domains can be finite sets, paramet- rically bounded (the domains are finite but their sizes de- pend on unbounded parameters), or infinite;

  . Let us recall that a static control program is a program which manipulate arrays, iteration domains are bounded -polyhedra and index functions of arrays are affine w.r.t. index variables. Figure 1.b gives an example of static control program.

	Definition 1 A System of Recurrence Equations is a set of
	equations called clauses, of the form:
		¡ £¢ ¥¤ §¦ © ¢ !#" %$ &$ '$ )( 0 1 32 45" 6¢ 87 8 $ '$ &$ 97 @$	(1)
	where and	¢ is an index vector, ( are arrays, B¨i s a function symbol, and ¦ A¨a domain of integer vectors, 1 2 C¨i s
	a function of indices. We introduce the following defini-
	tions: free ¦ 0D denotes the union of all the sets ¦ ¨, for all E , defining the clauses of ; Functions 1 2 C¨a re called dependence functions. Arrays
	that do not appear in the left-hand side (lhs) of any clause
	are called the inputs of the SRE. The outputs are special ar-
	rays defined in a lhs of some clauses. Note that there can
	be several output variables in a SRE. When the domains are union of -polyhedra and functions 1 F2 4¨a re affine w.r.t.
	index variables, then the SRE is a SARE.
	Moreover a SRE must satisfy the single assignment prop-
	erty, i.e. each value of	is defined uniquely, and we assume
	that all values of arrays which are not inputs are defined in
	the SRE.

  Note that in the preceding expressions, @ e ! such that the two expressions are the same, syntactically. Hence, the problem of matching is reduced to a semi-unification problem, between the terms of the template and of the program, where unknowns are template variables

		its value will be denoted c	is not defined,
		(here § and © ), and closed variables are program inputs (here , £ £ &
	For the program, G gf is h pi &q Dr ts "h vu ¦q Yh xw xq Dr ts y d ps q h xw xq Yr ts dy t d t £ v E " p £
	with	y p t d Ch pu &q Dh w xq r ts dy p t r W t £ s eq Yh xw xq Dr ts y r W v 9

. The problem is to find out the possible values of § and © and £

  allow one to associate to each state a regular expression representing all paths from the initial state to the current state. To obtain the reachability relation from such a regular expression, replace each letter by the corresponding firing relation, concatenation by relation composition, alternation by union and Kleene star by transitive closure. The reachable set is obtained by composing the result with the reachable set of the initial state.


Nadjet Bourdache
email: nadjet.bourdache@gmail.com

Patrice Perny
email: patrice.perny@lip6.fr

Anytime Algorithms for Adaptive Robust Optimization with OWA and WOWA

Keywords: Robust optimization, Preference Elicitation, OWA, WOWA, Ranking Algorithms, Assignment Problem

We consider optimization problems in graphs where the utilities of solutions depend on different scenarios. In this context, we study incremental approaches for the determination of robust solutions, i.e. solutions yielding good outcomes in all scenarios. Our approach consists in interleaving adaptive preference elicitation methods aiming to assess the attitude of the Decision Maker towards robustness or risk with combinatorial optimization algorithms aiming to determine a robust solution. Our work focuses on the use of ordered weighted average (OWA) and weighted ordered weighted average (WOWA) to respectively model preferences under uncertainty and risk while accounting for the idea of robustness. These models are parameterized by weighting coefficients or weighting functions that must be fitted to the value system of the Decision Maker. We introduce and justify anytime algorithms for the adaptive elicitation of these parameters until a robust solution can be determined. We also test these algorithms on the robust assignment problem.

Introduction

The practice of decision support in complex environments has shown the importance of developing new models and algorithms for optimization under partial information. Uncertainty is pervasive in discrete optimization and various contributions concern robust optimization problems [START_REF] Kouvelis | Robust discrete optimization and its applications[END_REF] in which multiple instances of the same problem must be considered simultaneously. These instances correspond to possible scenarios and the goal is to determine a feasible solution that remains as good as possible in all scenarios.

In graph optimization, several problems have been revisited under this perspective. For example, assuming that uncertainty only impacts the valuation of the graph (and not on its structure), several contributions address the discrete scenario case, and propose various reformulations of shortest-path problems, assignment problems, minimum spanning tree problems, as a min-max or min-max regret optimization problems, see [START_REF] Kouvelis | Robust discrete optimization and its applications[END_REF][START_REF] Murthy | Solving min-max shortest-path problems on a network[END_REF][START_REF] Yu | On the robust shortest path problem[END_REF][START_REF] Yu | Min-max optimization of several classical discrete optimization problems[END_REF][START_REF] Aissi | Complexity of the minmax and minmax regret assignment problems[END_REF][START_REF] Vladimir | On the robust assignment problem under a fixed number of cost scenarios[END_REF]. For the same problems, a similar work has been carried out in the case of graphs valued with interval data, which correponds to an infinity of possible scenarios [START_REF] Kouvelis | Robust discrete optimization and its applications[END_REF][START_REF] Yaman | The robust spanning tree problem with interval data[END_REF][START_REF] Ionut | On the complexity of the robust spanning tree problem with interval data[END_REF][START_REF] Montemanni | An exact algorithm for the robust shortest path problem with interval data[END_REF][START_REF] Aissi | Complexity of the minmax and minmax regret assignment problems[END_REF][START_REF] Wu | Algorithms for the min-max regret generalized assignment problem with interval data[END_REF]. In this paper we consider the discrete scenario case and propose new (interactive) approaches for robust optimization, with an implementation on the assignment problem. Let us introduce an instance of the robust assignment problem, for illustrative purpose:

Example 1 We consider a problem with 4 items that must be assigned to 4 agents (one item per agent and one agent per item). The utility of every item for every agent depends on the context which remains unknown. Three scenarios {s 1 , s 2 , s 3 } are considered leading to the following three utility matrices: where U s is the utility matrix giving in row i and column j the utility u s ij of item j for agent i in scenario s, for i = 1, . . . , 4, j = 1, . . . , 4, s = 1, . . . , 3. We want to find an assignment that remains as good as possible in all possible scenarios.

U1 =    
Any solution of an n × n assignment problem is characterized by a one-to-one mapping α defined from the set of agents to the set of items and associating item α(i) to agent i for i = 1, . . . n. Equivalently assignment α will be represented by the set of arcs {(i, α(i)), i = 1, . . . , n} in the bi-partite assignment graph. Given an assignment α, its total utility in scenario s is denoted as u s (α) and defined by u s (α) = n i=1 u s iα(i) , s = 1, . . . , n. If q distinct scenarios are considered, then any assignment α is characterized by the utility vector u(α) = (u 1 (α), . . . , u q (α)) representing the possible utilities in the q scenarios. For example, in the problem described in Example 1, the assignment α = {(1, 3), (2, 1), [START_REF] Benabbou | Incremental weight elicitation for multiobjective state space search[END_REF][START_REF] Ionut | On the complexity of the robust spanning tree problem with interval data[END_REF], [START_REF] Boutilier | Constraintbased optimization and utility elicitation using the minimax decision criterion[END_REF][START_REF] Boutilier | Constraintbased optimization and utility elicitation using the minimax decision criterion[END_REF]} leads to the utility vector u(α) = [START_REF] Yu | On the robust shortest path problem[END_REF][START_REF] Heinen | Fairness and rank-weighted utilitarianism in resource allocation[END_REF][START_REF] Yu | Min-max optimization of several classical discrete optimization problems[END_REF] whereas α = {(1, 3), [START_REF] Ionut | On the complexity of the robust spanning tree problem with interval data[END_REF][START_REF] Boutilier | Constraintbased optimization and utility elicitation using the minimax decision criterion[END_REF], [START_REF] Benabbou | Incremental weight elicitation for multiobjective state space search[END_REF][START_REF] Aissi | Complexity of the minmax and minmax regret assignment problems[END_REF], (4, 2)} leads to the utility vector u(α ) = [START_REF] Perny | A decision-theoretic approach to robust optimization in multivalued graphs[END_REF][START_REF] Ramsay | Monotone regression spline in action[END_REF][START_REF] Menahem | The dual theory of choice under risk[END_REF]. Hence the comparison of α and α amounts to comparing utility vectors u(α) and u(α ).

In order to be able to compare the utility vectors attached to feasible solutions, we need to assess the attitude of the Decision Maker (DM) toward uncertainty or risk (the latter situation occurring when the probabilities of the scenarios are known), and to decide how the outcomes attached to the q scenarios must be aggregated to define the overall value of a solution. One standard criterion for decision making under uncertainty is the Laplace criterion that consists in maximizing the average utility taken over all scenarios. In Example 1, the optimal assignment w.r.t. the Laplace criterion is α with an average of 65/3.

A more cautious approach consists in choosing a solution maximizing the utility measured in its worst scenario. This corresponds to max-min optimization, a standard approach in robust optimization (equivalent to min-max optimization when costs are considered). In Example 1, the max-min optimal assignment is α introduced above. It can be seen as more robust than α since it guarantees a utility greater or equal to 18 whereas we could obtain 8 with α in scenario s 2 .

However the max-min criterion is often considered as overpessimistic and poorly discriminating. In particular, we may obtain solutions having a much lower average utility than with the Laplace criterion. On the other hand, by compensating high and low utilities, the Laplace criterion does not provide any guarantee on the robustness of solutions. In order to have a better flexibility in modelling the DM's attitude toward uncertainty, and to find compromise attitudes between overpessimism and full compensation, we will use an ordered weighted average (OWA) to aggregate the utilities obtained in the different scenarios. OWA operators have been introduced by Yager [START_REF] Yager | On ordered weighted averaging aggregation operators in multicriteria decision making[END_REF] and axiomatically justified in the context of robust discrete optimization in [START_REF] Perny | A decision-theoretic approach to robust optimization in multivalued graphs[END_REF]. They are also widely used in fair multiagent optimization for their ability to generate Pareto-optimal solutions with well-balanced profiles [START_REF] Ogryczak | Multicriteria models for fair resource allocation[END_REF][START_REF] Lesca | Lp solvable models for multiagent fair allocation problems[END_REF][START_REF] Heinen | Fairness and rank-weighted utilitarianism in resource allocation[END_REF].

As we shall see later more formally, the OWA of a given vector is a kind of weighted sum where the weights are not attached to positions of components in the vector but to their ranks. This allows the importance attached to good or bad outcomes to be controlled and provides a continuum of attitudes in the aggregation, ranging from the Laplace criterion (modelled by the average) to pure pessimism (modelled by the min). These various attitudes are defined by the weighting vector parameterizing the OWA model.

When the probabilities of the scenarios are known, the scenarios do not play symmetric roles. In this case we will use the WOWA model [START_REF] Torra | The weighted owa operator[END_REF] also known as Yaari's model [START_REF] Menahem | The dual theory of choice under risk[END_REF] which is a non-symmetric extension of OWA allowing weights to be attached to components (here scenarios). In the latter model, the DM's attitude toward risk is controlled by a probability weighting function.

Whether OWA or WOWA is used, assessing the weighting parameter of the model is a critical issue. The aim of this paper is to propose an incremental elicitation method for facilitating the parametrization of these models in order to determine solutions that are well fitted to the value system of the DM. Our aim is not to determine precise values of these parameters for the DM, prior to the optimization stage. Instead, we propose to interleave preference queries with the exploration of solutions in order to progressively reduce the uncertainty attached to these weighting parameters until a robust solution can be found. These preference queries consist in asking the DM to compare between pairs of solutions, his preferences can then be easily translated into linear constraints that we can add to our model. By integrating the elicitation to the exploration, we aim to save a large part of the elicitation burden.

The paper is organized as follows: we introduce in Section 2 some background on OWA and WOWA. Then we present in Section 3 a ranking algorithm for the determination of possibly optimal solutions when the weighting parameters of OWA or WOWA are imprecise. In Section 4 we present an incremental elicitation algorithm to determine or approximate an OWA-optimal or WOWA-optimal solution. In Section 5 we implement the proposed approach on the robust assignment problem and present numerical tests, as well as some preliminary results on the robust shortest path problem.

Models for Robust Optimization

In this section, we consider a general robust optimization problem involving q distinct scenarios and we discuss the evaluation of solutions represented by vectors of type (x 1 , . . . , x q ) where x i represents the utility of solution x in scenario i. In particular we recall some background on OWA and WOWA operators.

OWA Optimization. In robust optimization problems with discrete scenarios, one basically looks for Pareto-optimal solutions having a well-balanced utility profile. This vision of robustness in the context of uncertainty can be related to the notion of fairness in social choice (the scenarios acting as different agents providing different views on solutions). In particular, the preference of the DM is expected to satisfy monotonicity w.r.t. -transfers, a standard axiom used in inequality measurement that reads as follows: for all j, k j = k, for all such that 0 < x k -x j < , (x 1 , . . . , x j + , . . . , x k -, . . . , x q ) (x 1 , . . . , x j , . . . , x k , . . . , x q ).

Moving from a solution to another one using such -transfers contributes to reducing the utility gap between pairs of scenarios and thus makes the DM better off. It is well known that the minimal preference relation satisfying this condition is the Lorenz dominance relation (L-dominance) [START_REF] Albert W Marshall | Inequalities: theory of majorization and its applications[END_REF] defined by x L y if and only if L(x) P L(y) where P is the Pareto dominance and L(x) is the Lorenz vector the i th component of which is defined by

L i (x) = i k=1 x σ(k)
, where σ is the permutation of (1, . . . , q) that reorders the components of x by increasing order (x σ(i) ≤ x σ(i+1) , i = 1, . . . , q -1). However, L-dominance is a partial order and many solutions remain incomparable. A natural way to extend this partial order, axiomatically justified in the context of robust optimization [START_REF] Perny | A decision-theoretic approach to robust optimization in multivalued graphs[END_REF], is to resort to Ordered Weighted Averages (OWA for short).

OWA is an aggregation function that weights the components of a vector in function of their rank. Let w ∈ R q + be a weighting vector. The OWA defined by w reads: f (x, w) = q i=1 w i x σ(i) . This function is symmetric because the weights are not attached to the components of x but to the components of the reordered vector (x σ(1) , . . . , x σ(q) ). The OWA family of aggregation functions f (x, w), w ∈ R q + includes the minimum, the maximmum, the median and all order statistics as particular cases.

OWA is widely used in fair optimization because it enables a linear extension of the Lorenz Dominance order. Remark indeed that x σ(i) = L i (x) -L i-1 (x) for all i > 1. Hence we have: f (x, w) = q-1 i=1 (w i -w i+1 )L i (x)+w q L q (x). Thus, function f is nothing else but a linear combination of the components of the Lorenz vector. Then, if the weights are decreasing (i.e. w i > w i+1 for all i = 1, . . . , q -1), then the coefficients w i -w i+1 are strictly positive. Therefore, an f -optimal solution is necessary L-non-dominated. Thus, OWA used with strictly decreasing weights w i leads to Pareto-optimal solutions that cannot be improved in terms of -transfers reducing inequalities. This favours solutions having a well-balanced utility vector. Considering the properties of OWA recalled above, robust optimization under uncertainty can soundly be reformulated as the problem of maximizing f (x, w) for some w with decreasing components.

Example 1 (continued) Let us compare solutions α and α such that u(α) = [START_REF] Yu | On the robust shortest path problem[END_REF][START_REF] Heinen | Fairness and rank-weighted utilitarianism in resource allocation[END_REF][START_REF] Yu | Min-max optimization of several classical discrete optimization problems[END_REF] and u(α ) = (18, 20, 25) using f (., w) with w = (1/2, 1/3, 1/6) we obtain f (u(α), w) = 8/2+28/3+29/6 = 18.17 and similarily f (u(α ), w) = 119/6 = 19.83. Here, we observe that α is preferred to α because the average loss in utility incurred when passing form u(α ) to u(α) is compensated by the improvement of the minimum (18 against 8). Remark that, although α is the max-min optimal solution in Example 1, it is not OWA-optimal. The OWA-optimal assignment for the chosen weights is α = {(1, 1), [START_REF] Ionut | On the complexity of the robust spanning tree problem with interval data[END_REF][START_REF] Boutilier | Constraintbased optimization and utility elicitation using the minimax decision criterion[END_REF], [START_REF] Benabbou | Incremental weight elicitation for multiobjective state space search[END_REF][START_REF] Ionut | On the complexity of the robust spanning tree problem with interval data[END_REF], (4, 3)} with u(α ) = [START_REF] Wang | Incremental utility elicitation with the minimax regret decision criterion[END_REF][START_REF] Wu | Algorithms for the min-max regret generalized assignment problem with interval data[END_REF][START_REF] Ogryczak | Multicriteria models for fair resource allocation[END_REF] and f (u(α ), w) = 121/6 = 20.17. This is a compromise between α and α that improves the average utility of α but slightly downgrades the worst case utility.

WOWA Optimization. As explained before, one typical property of OWA is to be a symmetric aggregator. This property seems natural when the same attention or importance is attached to every scenario. This is no longer the case when the probabilities (p 1 , . . . , p q ) of the scenarios are known. In such cases we may naturally consider a weighted extension of OWA defined as follows:

g(x, ϕ) = q i=1 [x σ(i) -x σ(i-1) ]ϕ( q k=i p σ(k) ) (1) = q i=1 [ϕ( q k=i p σ(k) ) -ϕ( q k=i+1 p σ(k) )]x σ(i) (2)
where x σ(0) = 0; function ϕ is strictly increasing on the unit interval. A solution x is as least as good as a solution y when g(x, ϕ) ≥ g(y, ϕ). This formulation is known as the Yaari's model in the literature on decision under risk because it has been introduced and axiomatically justified by Yaari [START_REF] Menahem | The dual theory of choice under risk[END_REF] in this context. Since ϕ is increasing on the unit interval, the preferences induced by the Yaari's model are monotonic with respect to first-order stochastic dominance (FSD). This means that if two solutions x and y are such that the probability G x (t) = P [u(x) > t] is greater or equal than the probability G y (t) = P [u(y) ≥ t] for all t then g(x, ϕ) ≥ g(y, ϕ). Moreover, a necessary and sufficient condition for these preferences to be monotonic with respect to second order stochastic dominance (SSD) is the convexity of ϕ. This means that x SSD y implies g(x, ϕ) ≥ g(y, ϕ) whenever ϕ is convex, where x SSD y means that

u -∞ G x (t)dt ≥ u -∞ G y (t)
dt for all u. This property established in [START_REF] Chew | Risk aversion in the theory of expected utility with rank dependent probabilities[END_REF] has a major importance in the context of robust optimization since monotonicity with respect to SSD is the standard model of strong risk aversion. This is due to the fact that SSD is strongly related to the existence of mean-preserving spread of the utility distribution increasing the risk attached to a solution (just as Lorenz dominance relates to the existence of -transfers reducing inequalities); for more details see [START_REF] Joseph | Increasing risk. i. a definition[END_REF]. For this reason, we shall use convex functions ϕ to account for risk aversion in robust optimization problem with probabilities over scenarios.

Example 2 We come back to Example 1 and assume now that the probabilities of the 3 scenarios are given by p = (1/2, 1/6, 1/3). Let us compare solutions α and α such that u(α) = (29, 8, 28) and u(α ) = (18, 20, 25) using Yaari's model with ϕ(x) = x 2 . We have: g(α, ϕ) = 8 + (28 -8)ϕ(5/6) + (29 -28)ϕ(1/2) = 22.1 and g(α , ϕ) = 18 + (20 -18)ϕ(1/2) + (25 -20)ϕ(1/3) = 19.1. Therefore α is preferred to α . The fact that scenario s 1 has a high probability compared to the others gives the advantage to α even if it is more risky than α . Note that the conclusion would be different if the requirement of risk-aversion were strengthened, using a function ϕ of higher convexity. For example, if ϕ(x) = x 4 we have g(α, ϕ) = 17.7 whereas g(α , ϕ) = 18.2. In this case α is preferred to α.

Note that when scenarios are equiprobable then p i = 1/q and the coefficient of x σ(i) in Equation( 2) becomes:

w i = [ϕ( q-i+1 q )-ϕ( q-i q )
] which is now constant (independent of σ). In this case the Yaari's model reduces to a standard OWA. Hence Yaari's model can be seen as a weighted generalization of OWA which explains the name WOWA due to Torra [START_REF] Torra | The weighted owa operator[END_REF] and used in this paper.

Considering the properties of WOWA recalled above, robust optimization under multiple scenarios of known probabilities can soundly be reformulated as the problem of maximizing g(x, ϕ) for some proper convex weighting function ϕ. The determination of an f -optimal or g-optimal solution in graph optimization is generally a challenging problem. Standard constructive algorithms based on dynamic programming or greedy search do not apply directly because neither f -optimality nor g-optimality satisfy the Bellman principle. An f -optimal solution for a given vector w can include sub-optimal subsolutions due to the nonlinearity of f with respect to outcomes. For example, when w = {(1, 0, . . . , 0)}, f -optimization is nothing else but -min maximization ( -max minimization for costs) which is known to be NP-hard for assignment problems but also for shortest-path problems [START_REF] Kouvelis | Robust discrete optimization and its applications[END_REF][START_REF] Aissi | Complexity of the minmax and minmax regret assignment problems[END_REF]. The same remarks holds for g-optimality. Yet some solution methods are available both for OWA and WOWA optimization, e.g., [START_REF] Ogryczak | Multicriteria models for fair resource allocation[END_REF][START_REF] Lesca | Lp solvable models for multiagent fair allocation problems[END_REF].

Optimization with Imprecise Parameters

Possibly OWA-optimal Solutions. It is often difficult to precisely determine the weighting vector w to be used in the OWA model. Indeed, the only prior information we have is that the weighting vector is positive and strictly decreasing, providing an exponential number of possible weights, and then an exponential number of preference queries to be assessed. Yet, in most cases, some preference information is available, putting some constraints on the set of admissible weighting vectors. Note that, any judgment of type "i prefer x to y" (where x and y are two feasible utility vectors) translates into the following inequality: f (x, w) ≥ f (y, w), which is a linear constraint in w bounding the set of admissible weighting vectors W . Thus, we define the uncertainty set W as a convex polyhedron including all weighting vectors compatible with the preference information collected so far. Given the uncertainty set W and a set X of utility vectors attached to feasible solutions, we define P O W (X) as the set of possibly f -optimal solutions in X, i.e., the elements of X which are f -optimal for some weighting vector w in W . More formally: ∀X ⊆ R n , P O W (X) = w∈W arg max x∈X f (x, w).

The uncertainty set W allows a W -dominance relation over utility vectors to be defined as follows: x W y ⇐⇒ [∀w ∈ W, f (x, w) > f (y, w)]. This relation can be extended to set-wise dominance: a solution x ∈ X is said to be dominated by a set Y ⊆ X (denoted Y W x) if there exists y ∈ Y such that y W x. If Y is explicitly defined, we can decide whether a solution x is dominated by Y in polynomial time by testing whether min w∈W max y∈Y f (y, w)-f (x, w) > 0 which is done by solving a linear program. Hence, when X is given explicitly, P O W (X) can be computed in polynomial time. It is sufficient to iteratively remove from X any solution x i dominated by X, for i = 1, . . . , n. This was shown in [START_REF] Benabbou | Incremental weight elicitation for multiobjective state space search[END_REF] for decision models based on the minimization of a weighted sum. We kept the same idea while adapting the algorithm for the OWA maximization problem.

In robust optimization problems, the set X is only implicitly defined and computing P O W (X) is more difficult. To overcome the problem we now introduce a ranking approach to compute P O W (X). This approach consists of three steps: 1) linear scalarization: a scalar valued instance of the problem is constructed by replacing utility vectors attached to the edges of the graph by their average over all scenarios. 2) ranking: we perform an enumeration of solutions by decreasing order of utilities; several algorithms are available in the literature to rank the solutions of an optimization problem by decreasing order of preferences, e.g., Murty algorithm [START_REF] Katta | An algorithm for ranking all the assignments in order of increasing costs[END_REF] for assignment problems or Eppstein algorithm [START_REF] Eppstein | Finding the k shortest paths[END_REF] for shortest path problems. 3) stopping condition: we stop the enumeration when we can prove that all possibly optimal solutions have been enumerated. The stopping condition used in step 3 is justified by the following propositions: Proposition 1 For any weighting vector w ∈ R q + with decreasing weights such that q i=1 w i = 1, we have, for all x ∈ R q , f (x, w) ≤ 1 q q i=1 x i .

Proposition 1, that directly derives from the result presented in proposition 3 of [START_REF] Galand | Exact algorithms for owa-optimization in multiobjective spanning tree problems[END_REF], allows us to establish the following result:

Proposition 2 Let X be the set of all feasible vectors, and let X k = {x 1 , . . . , x k } be the list of the k best elements of X ordered by decreasing average. We have:

max x∈P O W (X k ) min w∈W f (x, w) > 1 q q i=1 x k i ⇒ P O W (X) ⊆ P O W (X k ).
Proof. We show that when the "if" condition holds at step k, then any element that does not belong to X k cannot be optimal. Let us consider

x ∈ X \ X k . Then if max x∈P O W (X k ) min w∈W f (x, w) > 1 q q i=1
x k i then there exists y ∈ X k such that min w∈W f (y, w) > 1 q q i=1 x k i . Moreover since x comes after x k in the ordered enumeration, we have 1 q q i=1 x k i ≥ 1 q q i=1 x i ≥ f (x, w) for any w ∈ W . Hence, for any w ∈ W we have f (y, w) > f (x, w) and therefore x / ∈ P O W (X). This shows that P O W (X) ⊆ X k . Moreover an element in P O W (X) cannot be

W -dominated in X k since X k ⊆ X. Hence P O W (X) ⊆ P O W (X k ).
Proposition 2 provides a stopping condition for a ranking algorithm based on the mean of utilities to determine the set P O W (X). If at step k of the ranking algorithm the condition is fulfilled, then all solutions that would be enumerated after step k are W -dominated by a solution in P O W (X k ). The enumeration can be stopped since the elements we are looking for are all included in P O W (X k ) which can easily be computed (as previously explained) since X k is explicitly known at step k. To summarize, when X is implicitly defined, the set P O W (X) can still be determined by ranking its elements by increasing average utilities. This ranking is stopped as soon as the condition described above is activated. Remark that, for the robust version of assignment problem [START_REF] Kouvelis | Robust discrete optimization and its applications[END_REF][START_REF] Aissi | Complexity of the minmax and minmax regret assignment problems[END_REF], we cannot expect the stopping condition to be activated within a polynomial number of steps since the problem of determining P O W (X) is already NP-hard when W is reduced to the single vector {(1, . . . , 0)}. Nevertheless, we shall see in the section dedicated to numerical tests that the stopping condition is activated after a reasonable number of iterations on average. Moreover, the stopping condition can be relaxed into

1 q q i=1 x k i -max x∈P O W (X k ) min w∈W f (z, w) ≤ δ
where δ is a positive threshold representing the maximum admissible error. This possibly saves multiple iterations while providing good approximations of P O W (X).

Possibly WOWA-optimal Solutions. The ranking approach described above can also be adapted to compute P O Φ (X) for g-optimization under strong risk aversion. In this case ϕ is assumed to be convex. Consistently with the approach proposed for OWA, we consider now an uncertainty set Φ of admissible convex functions (functions compatible with the preferences observed so far). We will show later that, under mild hypothesis, Φ can also be represented by a convex polyhedron. Given the uncertainty set Φ, we define P O Φ (X) as the set of possibly g-optimal solutions in X, i.e., the elements of X which are g-optimal for at least one function ϕ in Φ. More formally: ∀X ⊆ R n , P O Φ (X) = ϕ∈Φ arg max x∈X g(x, ϕ). The associated dominance relation is defined by: x Φ y ⇔ [∀ϕ ∈ Φ, g(x, ϕ) > g(y, ϕ)]. This relation obviously extends to set-wise dominance. Deciding whether a solution y is dominated by a set X amounts to testing whether min ϕ∈Φ max x∈X g(X, ϕ) -g(y, ϕ) > 0 which may be done by linear programming provided that Φ is represented by a convex polyhedron. Moreover, when X is explicitly defined, the set P O Φ (X) can easily be computed by iteratively eliminating dominated elements, as done for P O W (X). Now, in order to install a ranking procedure to determine P O Φ (X) when X is implicitly defined, we need to establish a counterpart of Propositions 1 and 2 for g-optimization. This is exactly the role of the two following propositions: Proposition 3 For any convex function ϕ : [0, 1] → [0, 1] such that ϕ(0) = 0 and ϕ(1) = 1 we have, for all x ∈ R q , g(x, ϕ)

≤ q i=1 p i x i . Proof. Since ϕ is convex we have ϕ(ta + (1 -t)b) ≤ tϕ(a) + (1 -t)ϕ(b) for all a, b, t ∈ [0, 1]. Setting a = 1 and b = 0 we obtain: ϕ(t) ≤ t for all t ∈ [0, 1]. Hence ϕ( n k=i p σ(k) ) ≤ n k=i p σ(k) for i = 1, . . . , q. Therefore we have: g(x, ϕ) ≤ q i=1 x σ(i) -x σ(i-1) q k=i p σ(k) by Equation (1) since x σ(i) -x σ(i-1) ≥ 0 for i = 1, . . . , q. Hence we obtain : g(x, ϕ) ≤ q i=1 q k=i p σ(k) - q k=i+1 p σ(k) x σ(i) = q i=1 p σ(i) x σ(i) = q i=1 p i x i .
Proposition 4 Let X be the set of all feasible vectors and X k = {x 1 , . . . , x k } be the list of the k best elements of X ordered by decreasing expected utility. Let p 1 , . . . , p q be the probabilities of the q scenarios. The following property holds:

max x∈P O W (X k ) min ϕ∈Φ g(z, ϕ) > q i=1 p i x k i ⇒ P O Φ (X) ⊆ P O Φ (X k ).
The proof is very similar to the one of Proposition 2 and is deliberately omitted. To complete the parallel with the approach proposed for computing P O W (X), when X implicitly defined, the set P O Φ (X) can be determined by ranking the solutions of the instance obtained by replacing the utility vectors attached to the edges of the graph by their expected utility. This approach applies to problems for which a ranking algorithm is available, in particular, assignment problems, shortest path and minimum spanning tree problems.

Interleaving Elicitation and Ranking

Incremental Elicitation. As shown in Example 1, preferences induced by f (•, w) or g(•, ϕ) models may be sensitive to variations of their respective parameters, w and ϕ. It is therefore necessary to design elicitation procedures aiming to reduce the uncertainty set W (resp. Φ) introduced in the previous section. The elicitation of these parameters may require numerous preference queries if it is performed independently on the problem instance to be solved. For this reason, it is preferable to interleave elicitation and search. We suggest inserting preference queries in the ranking algorithm presented above in order to progressively enrich the set of preference statements and the list of constraints defining W (resp. Φ). This will iteratively reduce the W (resp. Φ), and therefore the set P O W (X) (resp. P O Φ (X)) until the obtention of a necessarily optimal solution, i.e., a solution that is f -optimal (resp. g-optimal) for all remaining parameter values in the uncertainty set. This incremental elicitation process should save a large part of the elicitation burden since an optimal solution can be identified although the parameters of the models remain largely imprecise.

Regret Minimization. We want to design an anytime algorithm that can return a valid solution to the problem even if it is interrupted before it ends. To make such a recommendation upon request at any step of the algorithm we use a standard regret based elicitation approach [START_REF] Wang | Incremental utility elicitation with the minimax regret decision criterion[END_REF] based on the following definitions:

PMR(x, y, W ) = max w∈W {f (y, w) -f (x, w)} (3) MR(x, X, W ) = max y∈X PMR(x, y) (4) MMR(X, W ) = min x∈X MR(x, X, W ) (5)
The pairwise regret PMR(x, y, W ) is the maximum regret of choosing x instead of y, defined as the maximum gap of OWA values. When X is explicitly known, this regret can easily be computed by linear programming since f is linear in w. The maximum regret (MR) attached to a solution x is the maximum regret of choosing x instead of any other solution. Finally, the minimax regret (MMR)

is the minimal MR regret over X. If the algorithm is stopped at a given step with a set X of possibly optimal solutions, then we shall recommend an element x in X that achieves the MMR. This solution will be named the MMR solution hereafter. As suggested in [START_REF] Wang | Incremental utility elicitation with the minimax regret decision criterion[END_REF], these regrets can also be used to select informative preference queries during an incremental elicitation process to iteratively reduce the MMR to zero. An efficient strategy introduced in [START_REF] Boutilier | Constraintbased optimization and utility elicitation using the minimax decision criterion[END_REF] under the name the Current Solution Strategy (CSS for short) consists in asking the DM to compare the current MMR solution x * with its strongest challenger defined by y * =arg max y∈X PMR(x * , y, W ). Whatever the answer to this query, a new constraint will be derived, further restricting the set W .

A similar approach could be implemented for WOWA optimization, using regrets PMR(x, y, Φ), MR(x, X, Φ) and MMR(X, Φ) that simply derive from (3-5) by substituting f (•, w) by g(•, ϕ). However the optimization of such regrets might be challenging because we have to optimize over a continuous set of weighting functions. To overcome this problem, we use a spline representation of function ϕ. Spline functions are piecewise polynomials whose elements connect with a high degree of smoothness. They are widely used in data interpolation due to their ability to approximate complex shapes [START_REF] Ramsay | Monotone regression spline in action[END_REF]. Interestingly enough, spline functions can be generated by linear combinations of basis spline functions. This allows to reduce the elicitation of a spline function to the determination of its weights in the spline basis. The use of spline representations for function ϕ in WOWA model has been recently introduced in [START_REF] Perny | Incremental preference elicitation for decision making under risk with the rank-dependent utility model[END_REF]. It enables an efficient incremental elicitation of the model to describe DM's preferences over probabilistic distributions. The proposed construction relies on the definition of ϕ as a convex combination of m basis spline functions of degree 3, increasing from 0 to 1 on the unit interval, and known as I-spline functions [START_REF] Ramsay | Monotone regression spline in action[END_REF]. More precisely we have: ϕ(x) = m j=1 b j I j (x) where I j (x), j = 1, . . . , m are the basic spline functions (see [START_REF] Perny | Incremental preference elicitation for decision making under risk with the rank-dependent utility model[END_REF][START_REF] Ramsay | Monotone regression spline in action[END_REF] for a formal definition of I j ).

Note however that this construction does not completely fit to our context because we have the additional constraint that ϕ must be convex, in order to enforce strong risk aversion, as explained in the previous section. To overcome this problem, we use another spline basis to generate spline functions that are both increasing and convex on the unit interval. To this end, ϕ(x) is defined by ϕ(x) = m j=1 φ j C j (x), j = 1, . . . , m where C j are C-spline functions defined as the normalized integrals of the I-spline functions. More precisely: C j (x) = x 0 I j (t)dt/ 1 0 I j (t)dt. As the integrals of positive and increasing functions, C j functions are increasing and convex. Moreover we have C j (0) = 0 and C j (1) = 1 for all j. Therefore, ϕ(x) will also be increasing and convex since coefficients φ j will be constrained to be non-negative. This is the model we use hereafter because g defined in this way is a linear function of coefficients φ j , j = 1, . . . , m, a key property for regret optimization. Hence, any preference constraint of type g(x, ϕ) ≥ g(y, ϕ) translates into a linear equation in coefficients φ j ; thus Φ can be soundly defined as the convex polyhedron of vectors (φ 1 , . . . , φ m ) compatible with the preference statements collected so far. We give below on Figure 1 the I-spline basis and the associated C-spline basis used to generate ϕ. Adaptative Elicitation in a Ranking Algorithm. Our adaptive elicitation procedure of OWA's weights uses progressive reductions of MMR values to discriminate between the solutions generated by the ranking algorithm. The baseline of our algorithm is a ranking procedure enumerating the solutions by decreasing average utilities (or increasing average costs). During the enumeration, if the current solution is not dominated by the previous ones, it is inserted in a bag of solution named X. When |X| = sb, preferences queries are generated to discriminate between the element of X and W is reduced accordingly until MMR(X, W ) = 0. Then all solutions with a strictly positive MR value are removed from X. The process is iterated until the stopping condition introduced in Proposition 2 holds. Due to this proposition, we know that, at this point, X includes all possibly optimal solutions. Moreover, since there exists x * ∈ X such that MR(x * , X, W ) = 0 (since MMR(X, W ) = 0), x * dominates all other solutions in X for all w ∈ W and is therefore necessarily optimal. This is the solution returned by the algorithm. The pseudo-code is given in Algorithm 1.

Function next(G, k) generates the k th best solution x k in the ranking process and returns its outcome vector (one component per scenario) or the empty vector if there are no solutions left. Function MMR(X, W ) returns the best MR value in X for the uncertainty set W . We keep asking queries, applying CSS, until the MMR becomes zero. The stopping condition of the ranking process holds when current-distance = 0 where current-distance is defined in line 20. So δ should be set to 0. However, in practice, the stopping condition can be relaxed by using a positive tolerance threshold δ in line 4 in order to save many iterations. More generally the algorithm can be stopped at any step k of the ranking (anytime property). The current MMR solution x * will be returned as the current best solution. It is necessarily optimal within the set {x 1 , . . . , x k } of solutions enumerated so far. Moreover, all solutions coming after x k in the ranking have an f -value lower than the average utility of x k . Hence current-distance provides an upper bound on the gap to optimality in the case of an early interruption of the ranking process. Consequently, when current-distance is less than δ, the gap to optimality for x * is at most δ.

Algorithm 1: MMR based elicitation in ranking for OWA optimization

Input: G = (V, E) undirected graph; ui(e) 0 edge utility in scenario i; q: number of scenarios, sb: size of bags; δ: tolerance threshold Output: x * : a near-optimal solution 1 X ← {}, W ← {w ∈ R q + : q i=1 wi = 1}, We use bags of size 5 (sb = 5) and simulate the answers of the DM using an OWA with w = (1/2, 1/3, 1/6). The best solution according to the mean value is x 1 = {(1, 3), (2, 1), (3, 2), (4, 4)} with utility vector u(x 1 ) = [START_REF] Yu | On the robust shortest path problem[END_REF][START_REF] Heinen | Fairness and rank-weighted utilitarianism in resource allocation[END_REF][START_REF] Yu | Min-max optimization of several classical discrete optimization problems[END_REF] and a mean at 21.67. Then, a lower bound of the OWA value of x 1 is obtained by minimizing f (u(x 1 ), w) = 8w 1 + 28w 2 + 29w 3 over all possible weighting vectors. We obtain 8, hence current-distance = 21.67 -8 > 0 and another iteration is necessary. The ordered enumeration continues until step 5 where the algorithm computes the MMR value which is strictlty positive. So a preference query is asked to the DM: the MMR solution x 2 = {(1, 1), (2, 4), (3, 2), (4, 3)} such that u(x 2 ) = [START_REF] Wang | Incremental utility elicitation with the minimax regret decision criterion[END_REF][START_REF] Wu | Algorithms for the min-max regret generalized assignment problem with interval data[END_REF][START_REF] Ogryczak | Multicriteria models for fair resource allocation[END_REF] must be compared to its best challenger x 4 chosen as explained before and such that u(x 4 ) = [START_REF] Perny | A decision-theoretic approach to robust optimization in multivalued graphs[END_REF][START_REF] Ramsay | Monotone regression spline in action[END_REF][START_REF] Menahem | The dual theory of choice under risk[END_REF]. The DM prefers x 2 therefore x 4 is removed from X and the constraint f (x 2 , w) ≥ f (y 4 , w) is added to the definition of polyhedron W i.e. 17w 1 + 23w 2 + 24w 3 ≥ 18w 1 + 20w 2 + 25w 3 or equivalently -w 1 + 3w 2 -w 3 ≥ 0. At this time, the MMR is still positive so the algorithm asks another query, this time between x 2 and x 1 (chosen in the same way as x 4 ), the DM still prefers x 2 , the algorithm proceeds similarly, removing x 1 and adding a constraint according to this preference. After this second question, the MMR is equal to 0, therefore all solutions but x 2 are removed from X. Then the ranking algorithm continues until step 15 without inserting new solutions in X because all of them are dominated by x 2 . At step 15, the stopping condition is activated and the algorithm returns the optimal solution x 2 .

x 1 ← next(G, 0), k ← 2, 2 current-distance ← 1 q q i=1 x k i -minw∈W f (x 1 , w) 3 enumeration-completed ← false 4 while (current-distance > δ and enumeration-completed = false) do 5 x k ← next(G, k -1) 6 if x k = emptyV ector then enumeration-completed ← true 7 else 8 if not(X W {x k }) then 9 X ← X ∪ {x k } 10 if |X| = sb and MMR(X, W ) > 0 then
A variant of Algorithm 1 can be used for WOWA optimization. It is sufficient to replace f by g and to modify the definition of regrets accordingly. The definition of the current distance must also be adapted not only by substituting f by g but also by replacing the weighted average by q i=1 p i x k i . The correctness of this variant derives from Proposition 4 instead of Proposition 2.

Numerical Tests

We have implemented Algorithm 1 on the robust assignment problem using both OWA and WOWA models. For these tests, function next(G, k) (line 5) was implemented with Murty's algorithm [START_REF] Katta | An algorithm for ranking all the assignments in order of increasing costs[END_REF]. The complexity of Murthy's algorithm to rank assignments by increasing utility is O(Kn 4 ) for n agents and K enumerations. We used the Gurobi library of Python to solve the linear programs required for dominance tests and regret minimization. During the elicitation steps, the DM's answers to preference queries are simulated using a hidden OWA or WOWA model. We evaluate the performance of the algorithm in terms of computation time, number of preference queries and number of ranking steps, we performed tests on multiple instances of different size and number of scenarios. For every case, performances are averaged over 20 runs. The tests are performed on a Intel Core i7-4770 CPU with 11GB of RAM. Table 1 shows the results obtained for OWA and WOWA elicitation and optimization. The performance are obtained for an error threshold δ set to 10% of the initial error (obtained for solution x 1 ), with a bag size of 10 and a timeout of 20 minutes (1200s). Time is given in seconds in the tables and the gap is the maximal error attached to the returned solution, expressed as a percentage of the range of the valuation scale. To test the generality of the approach, we also made some preliminary tests on the robust shortest path problem. In this problem, the arcs of the graph are valued by cost vectors with q components corresponding to the cost of the arc in the different scenarios. The definition of OWA and WOWA aggregators as well as Algorithm 1 have been modified to fit to minimization problems. For ranking the paths by increasing average costs, we used a lazy version of Eppstein's algorithm introduced in [START_REF] Víctor | A lazy version of Eppstein's k shortest paths algorithm[END_REF] with a complexity of O(m + n log n) in the worst case for a graph with n nodes and m arcs. To give an idea of the performance of Algorithm 1 on robust shortest path problems, we solve instances including 500 nodes and 5 scenarios in 548 seconds and 5 preference queries with a gap of 0.1 on average (after 28460 steps on average). In the case of 10 scenarios, the algorithm needs 681 seconds on average on similar graphs, 6.3 queries with a gap of 0.1 on average (after 41796 steps on average).

Conclusion

We have introduced a new adaptive elicitation approach for OWA and WOWA optimization and tested its practical efficiency on robust assignment problems and on robust shortest path problems in the discrete scenarios case. Our approach is quite general and applies to any other optimization problem for which an efficient ranking algorithm is known. An interesting extension of this work would be to design a similar approach for the incremental elicitation of Choquet integrals under the constraint of convex capacity (this is a more general model to account for robustness in optimization under uncertainty, including OWA and WOWA as special cases). The implementation of the ranking approach for Choquet integrals and the definition of a valid stopping condition are challenging questions because the capacity is imprecisely known.

Fig. 1 .

 1 Fig. 1. I-spline and C-spline cubic functions Ij(x) and Cj(x) for m = 5

16 X

 16 query to the DM (selected by CSS) 13 restrict W according to the answer 14 update regrets PMR, MR and MMR 15 until MMR(X, W ) = 0 ← {x ∈ X : MR(x, W ) = 0} return x * = arg minx∈X M R(x, X, W ), current-distance Example 1 (continued) Let us briefly illustrate the behavior of Algorithm 1 on Example 1.

Table 1 .

 1 Tests for OWA and WOWA optimization on the robust assignment problem

	Number of	q = 3	q = 5	q = 10
		agents time steps queries gap time steps queries gap time steps queries gap
	OWA	20 60 100	274 13607 1.3 1200 7236 4 1200 2075 5.1	0 719 34009 2.3 0.5 940 37739 1.1 1200 6945 6.7 5 1200 6382 11.6 4.7 2 1.4 2 1200 2550 11.3 6 1200 2216 17.5 7.5
	WOWA	20 60 100	270 5338 10.3 1200 5446 8.95 8.1 1200 3192 5.3 7.6 1200 4381 5.4 6.5 1 475 4657 9.45 1.1 865 6743 9.2 1.7 1200 2978 6.15 10.3 1200 2302 5 15.9 1200 1901 4.9 11.3