Antoine Jouglet
email: antoine.jouglet@utc.fr

Dritan Nace
email: dritan.nace@utc.fr

Christophe Outteryck
email: christophe.outteryck@a-sis.com

Timetabling of sorting slots in a logistic warehouse

Keywords: warehouse management, sorting slots, scheduling, timetabling

come

Introduction

Logistic warehouses are becoming more and more complex installations: flows are continually increasing, systems are becoming more automated and logistic processes globalized. In an increasingly competitive environment, warehouse management software (comprising Warehouse Management Systems (WMS) and Ware-house Control Systems (WCS)) is becoming a vital part of operational decision making.

Since the introduction of automation systems in production and distribution environments a large number of research studies into operational decision systems have been published (see for example [3-5, 8, 9] for surveys in this area). In this paper we study a sub-problem which occurs at the end of a logistic stream in a warehouse and concerns the timetabling of sorting slots. To our knowledge, this problem has so far not been studied in the literature.

We consider a warehouse composed of several order preparation shops and a set of orders to be prepared over a given time horizon. Each order is composed of several order lines. An order line corresponds to a quantity of work to be processed in one of the preparation shops. Each order is associated with the particular transportation truck that is to transport it to its final customer. This makes it necessary for the subset of all orders for a given truck to be totally processed before the truck's scheduled departure. Once an order has been prepared it is dispatched to a sorting slot, that is to say an accumulation area where the order waits to be loaded onto the truck. The difficulty of the sorting slot timetabling problem is due to a specific resource constraint: for a given truck, from the time the first order destined for the truck begins to be prepared, and up until the truck's scheduled departure, a precise number of sorting slots must be allocated. Since several orders corresponding to different trucks are prepared simultaneously, and since the number of sorting slots is limited, timetabling the use of sorting slots is essential. The WCS dynamically makes decisions concerning the order lines to be processed in each preparation shop, and it needs to know which of these order lines can be started at a given time without adversely affecting the processing of other order lines by using up available sorting slots (see Figure 1). In other words, sorting slots need to be timetabled to ensure that all orders can be processed over the considered time horizon.

For instance, consider a simple example with only one preparation shop, 3 trucks and a capacity of 5 sorting slots. Truck 1 needs 3 sorting slots to process its orders and must leave the warehouse at time 100. Truck 2 needs 2 sorting slots to process its orders and must leave the warehouse at time 80. Truck 3 needs 2 sorting slots to process its orders and must leave the warehouse at time 50. Suppose that preparation shop ol 1 ol 2 . . . Fig. 2 A non-optimal use of sorting slots, illustrating why a timetable for sorting slots needs to be established.

the WCS has no timetable for the use of sorting slots and that it begins (at time 0) to start an order line ol 1 destined for truck 1 and, immediately afterwards (at time 1), an order line ol 2 associated with truck 2. This means that 5 sorting slots are entirely allocated to trucks 1 and 2 at least from time 1 to time 80 (the departure time of truck 2) (see Figure 2). Unfortunately, this makes it impossible to process order lines for truck 3, scheduled to depart earlier.

preparation shop ol 3 ol 2 . . . ol 1 . . . Suppose now that a timetable for using sorting slots has been established, stipulating that 2 sorting slots are to be used by truck 2 from time 0 to time 80, 2 are to be used by truck 3 from time 0 to time 50, and 3 are to be used by truck 1 from time 50 to time 100. This timetable must make it possible for all order lines to be processed. The WCS knows which order lines it is authorized to dynamically start at a given time t, and which order lines belong to the trucks using sorting slots over an interval of time including t. In this example, order line ol 1 cannot be started before time 50 (see Figure 3).

Note that it is not our concern to schedule the order lines in the preparation shops. That is the role of the WCS, which must make this kind of decision according to the current state (not precisely predictable) of the whole system. Most of the time, a preparation shop is a very complex system involving several subsystems. One of the functions of the WCS is to start the appropriate order lines at the right time in order to balance the work between all its subsystems. This is why we shall not be assuming below that we have complete information about the set of orders. Instead, we have access to macro-information corresponding to the amount of work related to each truck in each preparation shop. Thus, as we will see in the next section, the orders corresponding to a particular truck can be described by a set of operations (one for each preparation shop) whose processing times correspond to these amounts of work. However, it should be remembered that the schedule of these operations is not really of any importance. We simply need to be certain that for an established timetable there exists a possible scheduling of these operations. In other words, we have to find a timetable for the use of sorting slots, given an operation scheduling feasibility constraint.

The remainder of the paper is organized as follows. We first define the problem in Section 2. Two operational cases are considered: the non-preemptive case in which order lines for the same preparation shop and for the same transportation truck must be processed in a single operation, and the preemptive case in which this constraint is released. We then provide a study of the complexity of the problem in Section 3: we state that some particular cases of the problem are polynomially solvable, while the general problem is NP-complete [START_REF] Garey | Computers and intractability: a guide to the theory of NP-completeness[END_REF]. In the remaining sections, we propose optimization methods to solve it. Thus, in Section 4, we describe a dominance rule allowing integer linear programming formulations and constraintbased scheduling methods [START_REF] Ph | Constraint-based scheduling, applying constraint programming to scheduling problems[END_REF] to be built to solve the problem in practice. These methods are described in Section 5 for the non-preemptive case and in Section 6 for the preemptive case. Some experimental results are provided in Section 7 and we conclude in Section 8.

Problem definition

The problem in hand can be formalized as follows. Let n be the number of trucks and let m be the number of preparation shops in the warehouse. From a schedulingbased modeling point of view, a truck is assimilated to a job and a preparation shop is assimilated to a machine. Thus, let N = {1, . . . , n} be the set of jobs (trucks) and let M = {1, . . . , m} be the set of machines (preparation shops).

Each job i is composed of m operations {o i1 , o i2 , . . . , o im }. The processing time p ij of operation o ij corresponds to the duration of the whole set of order lines which have to be processed by the machine (preparation shop) j that corresponds to truck i, i.e. to the total amount of work which has to be done by shop j for truck i. Each of the m machines therefore has to process n operations (one per job). Let d i be the known departure date of truck i. In a feasible solution of the problem, each job i is then considered to have finished by time d i , while the starttime s i of i is to be determined such that all operations of i can be performed. The set of sorting slots is assimilated as a cumulative resource of capacity E, E being the number of sorting slots available in the warehouse. Thus, job i needs a given number n i of sorting slots from among E, from time s i up until time d i . Note that the period over which the resource constraint applies depends only on the start-time of the job, in contrast to classical scheduling problems where it depends on periods during which the job is actually being processed.

Let s ij and c ij be the variables representing the start time and the completion time respectively of operation o ij in a particular solution. In a feasible solution, each operation o ij can start after the start of job i and has to be completed before the departure time of the associated truck: s ij ≥ s i and c ij ≤ d i . Each machine is disjunctive and can process only one operation at a time. Note that unlike in a scheduling shop problem (see e.g. [START_REF] Brucker | Scheduling Algorithms[END_REF]), operations for the same job can be processed simultaneously (and, in fact, very often will be). Two cases are considered:

-The non-preemptive case where, once it has begun to be processed by a machine, an operation o ij must be continued until its completion (

s ij + p ij = c ij ).
This implies that all order lines in a preparation shop associated with a truck will be dealt with sequentially. -The preemptive case, where interruptions to the processing of operations are possible. Such interruptions are useful, since they can allow part of another more urgent operation associated with another truck to be executed. This may be desirable when more flexibility is preferred in preparation shops. In this case s ij is taken as the first time an element of operation o ij is processed by machine j, while c ij is the date at which the operation is completed. Thus we have

s ij + p ij ≤ c ij .
The problem is to determine the values of variables {s 1 , . . . , s n } such that no more than E sorting slots can be used at a time, while the scheduling of operations is feasible. If variables {s 1 , s 2 , . . . , s n } can be fixed such that all the previous constraints hold, it follows that n i sorting slots from among E will be used from s i to d i for truck i. Thus, the values of variables {s 1 , s 2 , . . . , s n } and data {d 1 , . . . , d n } provide the timetable for the use of the sorting slots. From now on we shall refer to this problem as "the search problem".

Of course, if several solutions exist, we should like to find the best solution from an operational point of view. In this context, the longer the span of time over which a truck can use sorting slots, the more flexibility the WCS has for processing the orders for the truck in question. Greater overall versatility is thus ensured by using the sorting slots for as long as possible, i.e. by leaving as few idle times as possible in the use of sorting slots. Let H = max i∈{1,...,n} d i be the time horizon of an instance of the problem. One interesting objective is therefore to minimize H × E -i∈N n i (d i -s i ), which is equivalent to minimizing i∈N n i s i .

From now on we shall refer to the search for a solution that minimizes i∈N n i s i as "the optimization problem".

We shall be using the following additional notation. Let T = {d i |i ∈ N }∪{0} = {t 0 = 0, . . . , t w = H} be the set of departure dates, their values t i being indexed in non-decreasing order, i.e. t 0 = 0 < t 1 < . . . < t w = H. Moreover, for all i ∈ N , let [i] be the position of d i in the ordered set T , i.e. [i] = k such that k ∈ {1, . . . , w} and d i = t k .

Some complexity results

In this section we propose several complexity results for the general problem and for some of its special cases.

The special preemptive case where we have n i=1 n i ≤ E is polynomially solvable. The resource constraint according to the number of sorting slots disappears in this case, and it only remains to determine whether the scheduling of the operations is feasible, which will enable us to solve m single-machine scheduling problems with deadlines (1| di |_) which are polynomially solvable [START_REF] Garey | Computers and intractability: a guide to the theory of NP-completeness[END_REF].

The other special case in which any p ij = 1 is also polynomially solvable. Here the operations can be scheduled identically on all the machines. The problem is thus reduced to a problem with m = 1. We then have to determine whether the operations can be scheduled without violating the resource constraint on sorting slots. To this end we iteratively schedule operations from H to 0 in decreasing order of departure times. When several operations have the same departure time, they are scheduled in decreasing order of the number of sorting slots. The process stops when all operations are scheduled or when the resource constraint on the sorting slots is no longer satisfied. If the latter is the case, then the instance is not feasible.

Proposition 1 The preemptive case with m = 1 is N P -complete.

Proof To prove this result, we state that the KNAPSACK problem [START_REF] Garey | Computers and intractability: a guide to the theory of NP-completeness[END_REF] can be reduced to the preemptive case with m = 1. The KNAPSACK problem can be stated as follows: INSTANCE: A positive integer B. Two sets {a i |i ∈ I = {1, . . . , u}} and {b i |i ∈ I = {1, . . . , u}} of u strictly positive integers. Two positive integers A and B. QUESTION: Is there a subset

I ′ ⊆ I such that i∈I ′ a i ≥ A and i∈I ′ b i ≤ B ?
The decision version of the preemptive sorting slot problem with m = 1 can be stated as follows: INSTANCE: A set N = {1, . . . , n} of n jobs and a positive number E of sorting slots. Each job i ∈ N has a completion time d i , requires n i sorting slots in order to be processed and has one operation whose processing time is p i . QUESTION: Is there a preemptive single-machine schedule of the operations such that no more than E sorting slots are used at any one time, and such that each operation is completed before the completion time of its job?

Consider an instance of the KNAPSACK problem. We build an instance of the preemptive sorting slot problem with n = u + 1 jobs, E = i∈I a i and we set H = 1 + i∈I b i . For i ∈ {1, . . . , u} we set d i = H, p i = b i and n i = a i . We set d n = H -B, p n = 1 and n n = A. Note that this reduction is polynomial.

Suppose that there exists a subset I ′ ⊆ I such that i∈I ′ a i ≥ A and i∈I ′ b i ≤ B. We process the operations as follows. Operations of jobs i ∈ I ′ are processed from time H -i∈I ′ p i to time H in any order. Operations of jobs i ∈ I ∪ {n} \ I ′ are processed from time 0 to time H -i∈I ′ p i . These jobs need i∈I\I ′ a i + A ≤ i∈I a i sorting slots to be processed. Since jobs in I ′ are totally processed after time H -i∈I ′ p i ≥ H -B (since i∈I ′ p i = i∈I ′ b i ≤ B), job n has released A sorting slots. Hence, there remain at most i∈I\I ′ a i occupied sorting slots and there are still i∈I ′ a i sorting slots available to be used from jobs in I ′ . Now suppose that there exits a valid preemptive schedule for the built instance of the preemptive sorting sort problem. Since job n requires at least A sorting slots from time H -B -1 to time H -B, then a subset I ′ ⊆ I is totally processed after time H -B letting at least A available sorting slots. This subset I ′ of jobs is then such that i∈I ′ a i ≥ A and i∈I ′ b i ≤ B.

Proposition 2

The non-preemptive case with m = 1 is strongly N P -complete.

Proof To prove this result, we first establish that the 3-PARTITION problem [START_REF] Garey | Computers and intractability: a guide to the theory of NP-completeness[END_REF] can be reduced to the general non-preemptive case with m = 1. The 3-PARTITION problem can be stated as follows: INSTANCE: A positive integer B. A set A = {a i |i ∈ I = {1, . . . , 3q}} of 3q integers with i∈I a i = qB and ∀i ∈ I, B/4 < a i < B/2. QUESTION: Is there a partition of I into q subsets {I 1 , . . . , I q } of cardinal 3 such that ∀k ∈ {1, . . . , q}, j∈I k a j = B ? The decision version of the non-preemptive sorting slot problem with m = 1 can be stated as follows: INSTANCE: A set N = {1, . . . , n} of n jobs and a positive number E of sorting slots. Each job i ∈ N has a completion time d i , needs n i sorting slots to be processed and has an operation whose processing time is p i . QUESTION: Is there a non-preemptive single-machine schedule of the operations such that no more than E sorting slots are used at a time and that each job is completed before the completion time of its job?

Consider an instance of the 3-PARTITION problem. We build an instance of the sorting slot problem with n = 5q jobs, E = 3q and we set H = qB + 2q. For i ∈ {1, . . . , 3q} we set d i = H, p i = a i and n i = 1. For k ∈ {0, . . . , q -1}, we set

d 3q+2k+1 = (B + 2)k + 1, d 3q+2k+2 = (B + 2)k + 2, p 3q+2k+1 = p 3q+2k+2 = 1 and n 3q+2k+1 = n 3q+2k+2 = 3(q -k)
, and so on. This reduction is polynomial.

Suppose that there exists a partition of I into q subsets {I 0 , . . . , I q-1 } of cardinal 3 such that ∀k ∈ {0, . . . , q -1} we have j∈I k a j = B. We build a feasible schedule as follows.

For k ∈ {0, . . . , q -1}, operation of job 3q + 2k + 1 is scheduled from time (B + 2)k to time (B + 2)k + 1 and operation of job 3q + 2k + 2 is scheduled from time (B + 2)k + 1 to time (B + 2)k + 2. Note that the departure time of each of these jobs is then respected. Moreover, the constraint on sorting slots holds, since these jobs are totally processed on disjoint intervals and never use more than 3q = E sorting slots. Next, let k iteratively take values in {0, . . . , q -1}. Operations of jobs of indices in I k are scheduled in any order within the interval [(B+2)k+2, (B+2)(k+1)[. Note that the machine is never idle during this interval, since j∈I k p j = j∈I k a j = B and (B + 2)(k + 1) -(B + 2)k -2 = B. Moreover, at least 3(q -k + 1) sorting slots are available from time (B + 2)k + 2, since job 3q + 2k + 2 has just finished. 3 of these sorting slots are used by jobs in I k up until time H, leaving the 3(q -k) other slots for the processing of jobs 3q + 2(k + 1) + 1 and 3q + 2(k + 1) + 2 within the interval [(B + 2)(k + 1), (B + 2)(k + 1) + 2[. Now suppose that there exists a feasible non-preemptive schedule for the jobs. This implies that the machine is never idle on interval [0, H[ since i=5q i=1 p j = i=3q i=1 a j + 2q = qB + 2q = H. Recall that during interval [H -1, H[ all sorting slots are necessarily used by jobs {1, . . . , 3q}. Job 5q requires at least n 5q = n 3q+2(q-1)+2 = 3(q -(q -1)) = 3 sorting slots in order to be processed during the interval [H -B -1, H -B[. Thus, at least 3 operations from jobs among {1, . . . , 3q} must start during the interval [H -B, H[ (see Figure 4). Moreover, no more than 3 operations can be processed over this interval of size B, because ∀i ∈ {1, . . . , 3q}, p i = a i > B/4 and preempting is not possible. As the machine is never idle, there is no other choice for job 5q -1 but to use the same 3 sorting slots as job 5q. Since d 5q-1 = H -B -1 and d 5q = H -B, the operation of job 5q is necessarily scheduled from time H -B -1 to time H -B. Therefore exactly 3 operations from jobs x, y and z such that p x + p y + p z = B are necessarily processed by the machine in interval [H -B, H[. We then set I q-1 = {x, y, z}.

. . . . . . 5q -2 . . . 5q -1 . . . 5q I q-1 = {x, y, z} H -2B -2 H -B -1 H -B H B + 1 B 3q sorting slots 3 3 . . . . . . . . . . . .

Fig. 4

Reducing the 3-PARTITION problem to the sorting slot problem (the part above the time-axis represents the schedule of jobs on the machine, while the part below represents the use of the sorting slots over the time horizon). Case of I q-1 Note that unlike job 5q, job 5q -1 can use these [START_REF] Cormier | A review of warehouse models[END_REF] . . .

(B + 2)k + 1 (B + 2)k + 2 (B + 2)(k + 1) + 1 (B + 2)(k + 1) + 2 B +
3(q -k) 3(q -k -1)
Fig. 5 Reducing the 3-PARTITION problem to the sorting slot problem (the part above the time-axis represents the schedule of jobs on the machine while the below part represents the use of the sorting slots over the time horizon). General case.

Consider k ∈ {0, . . . , q-2} and let us assume that, for all z ∈ {k+1, . . . , q-1}, job {3q + 2z + 2} starts at (B + 2)z + 1 and is completed at (B + 2)z + 2, that 3(q -z) operations of jobs in {1, . . . , 3q} are totally processed during the interval [(B + 2)z + 2, H[ and that job 3q + 2z + 1 uses 3(q -z) sorting slots from its start to its departure time (B + 2)z + 1. Note that this assumption is true for k = q -1, as previously shown.

Job 3q + 2k + 2 needs at least 3(q -k) sorting slots over the interval [(B + 2)k + 1, (B + 2)k + 2[ for its execution to be possible. Exactly 3(q -k -1) jobs in {1, . . . , 3q} are totally processed after time (B + 2)(k + 1) + 2, leaving 3(q -k -1) sorting slots which are all used by job 3q+2(k+1)+2 = 3q+2k+4 over the interval [(B + 2)(k + 1) + 1, (B + 2)(k + 1) + 2[ and by job 3q + 2(k + 1) + 1 = 3q + 2k + 3 from its start time to its departure time (B + 2)(k + 1) + 1 (see Figure 5). At most 4 operations belonging to jobs from among {1, . . . , 3q} may be executed over the interval [(B + 2)k + 2, (B + 2)(k + 1) + 1[ of size B + 1 (the total processing time for 5 operations belonging to these jobs is necessarily strictly greater than B + 1). This means that at most 4 sorting slots may be freed for job 3q + 2k + 2, which is not sufficient, since 3(q-k) >= 6. The only way that 3(q-k) sorting slots might become available for job 3q + 2k + 2 before time (B + 2)k + 2 would be for 3 jobs to start in {1, . . . , 3q} and for job 3q+2(k+1)+1 to start during [(B+2)k+2, (B+2)(k+1)+1[, which would leave exactly the 3 + 3(q -k -1) = 3(q -k) sorting slots needed by job 3q + 2k + 2. However, there is no alternative to processing the operations belonging to job 3q +2k +2 from time (B +2)k +1 to time (B +2)k +2, because job 3q +2k +1 needs the same available sorting slots before its departure time (B + 2)k + 1. Since job 3q + 2k + 3 needs one unit of time to process its operation on the machine during [(B + 2)k + 2, (B + 2)(k + 1) + 1[, we must have 3 operations of jobs x, y and z in {1, . . . , 3q} such that p x +p y +p z = B being processed on the B remaining unit times of this interval, thus giving the subset I k = {x, y, z}.

This shows that there exists a partition of I into q subsets {I 0 , . . . , I q-1 } of cardinal 3 such that for all k ∈ {0, . . . , q -1} we have j∈I k a j = B.

Dominance rules

Proposition 3 Let S be a feasible solution of an instance of the problem in which a job i is executed from s i to d i . It is possible to transform S into a feasible solution S ′ in which job i uses n i sorting slots from

s ′ i = max{t k |t k ∈ T ∧ t k ≤ s i }
and where the processing of operations is not modified.

Proof In a feasible solution S, note that the times at which some sorting slots become free after use can occur only at points {t 1 , . . . , t k } of the considered time horizon. If the cumulative constraint on sorting slots for a job i is satisfied from

preparation shops o i2 o i1 t 0 = 0 t 1 s ′ i = tx s i / ∈ T d i = ty H = tw 5 sorting slots
Fig. 6 Starting times for sorting slots correspond to the different values in set T s i to d i , then the corresponding sorting slot n i can start to be used at time s ′ i = max{t k |t k ∈ T ∧t k ≤ s i } without violating the cumulative constraint (see Figure 6).

We can deduce the following dominance rules.

Proposition 4 There exists at least one solution of the search problem in which the starting times for sorting slots for the different jobs correspond to the values in set T .

Proof This proposition is a direct consequence of Proposition 3.

Proposition 5 All optimal solutions of the optimization problem are such that the starting times for sorting slots for the different jobs correspond to the values in set T .

Proof Consider an optimal solution of the optimization problem in which the start time s i of job i is not in T . As shown in Proposition 3, if the cumulative constraint on sorting slots is respected from s i to d i for a job i, then the corresponding sorting slot n i can start to be used at time

s ′ i = max{t k |t k ∈ T ∧ t k ≤ s i } without violating the cumulative constraint. Since s ′
i < s i , we obtain a new solution with a smaller value of the objective function ( i∈N n i s i ), contradicting the hypothesis that the initial solution was optimal.

The previous two propositions can be used in several ways. From a practical point of view, they offer a way of simplifying the timetabling of the sorting slots, since timetables can be established in terms of periods whose time points correspond to the departure times of the different trucks. Of course, the additional use (in principle, futile) of sorting slots by a job can be seen as additional versatility in the preparation process. From the point of view of the problem resolution, these dominance rules may be used to establish efficient models, as shown in the following sections.

5 Models and approaches for the non-preemptive case

A mixed integer linear programming approach

In this section, we propose a first mixed integer linear programming model which allows the problem to be solved.

We could have proposed a classical time-indexed mixed-integer formulation for modeling the use of sorting slots by the jobs, where, for example, some binary variables x i,t are used to determine at any time t whether or not job i is using sorting slots. However, thanks to the dominance rule presented in the previous section, we know that if a job i uses sorting slots at time t, then it can use the same sorting slots at time t -1, if t -1 is not the departure time of a job. It is consequently possible to avoid taking the constraint on sorting slots into account, other than at each interval [t k , t k+1 [, with k ∈ {0, . . . w -1}.

This can be achieved through the use of binary variables x i,k for each job i and each interval [t k , t k+1 [, which are equal to 1 if job i uses n i sorting slots on [t k-1 ..t k [, and 0 otherwise. The constraints linked to these variables can be formulated as follows.

x i,k = 0, ∀i ∈ N, ∀k ∈ {[i] + 1, . . . , w -1} (1)
x i,[i] = 1, ∀i ∈ N (2)
x i,k ≤ x i,k+1 , ∀i ∈ N, ∀k ∈ {0, . . . , [i] -1} (3)
i∈N x i,k n i ≤ E, ∀k ∈ {0, . . . , w -1} (4)
Equation [START_REF] Ph | Constraint-based scheduling, applying constraint programming to scheduling problems[END_REF] ensures that job i cannot use sorting slot after the departure d i . Equation ( 2) means that job i should use sorting slots in the period [t

[i]-1 , t [i] [ (in
the interval just before departure of truck i). Equation (3) ensures that there is no interruption in the use of sorting slots from the start of job i until the departure time d i . Equation (4) ensures that in each interval [t k , t k+1 [, the total number of used sorting slots is at most E. Note that these equations are also valid in the preemptive case and will be used in the next section.

To compute possible values for the starting times s i of the jobs we can introduce variables {s i |i ∈ N } into the MIP by adding the following constraints, implying that the starting time of a job i that does not use sorting slots in interval [t k , t k+1 [ with t k+1 < d i is greater than or equal to t k+1 :

s i ≥ (1 -x i,k )t k+1 , ∀i ∈ N, ∀k ∈ {0, . . . , [i]} (5)
We note that an operation o ij can be processed only in the periods in which job i is active, i.e. in intervals [t k ..t k+1 [ such that x i,k = 1. We then need to show that the scheduling of operations is possible according to the values of the job start times. For this we use a classical scheduling model, with binary variables b j il ∈ {0, 1}, j ∈ M, i, l ∈ N, i < l such that b j il = 1 if operation o ij starts before operation o lj on machine j, otherwise b j il = 0. The search version of the problem can be mathematically formulated as follows:

s ij + p ij ≤ d i , ∀i ∈ N, ∀j ∈ M (6)
s i ≤ s ij , ∀i ∈ N, ∀j ∈ M (7)
s ij + p ij ≤ s lj + Ω × (1 -b j il ), ∀j ∈ M, ∀i, l ∈ N, i < l (8)
s lj + p lj ≤ s ij + Ω × b j il , ∀j ∈ M, ∀i, l ∈ N, i < l (9) b j il ∈ {0, 1}, ∀j ∈ M, ∀i, l ∈ N, i < l (10)
Ω is an arbitrary big number which can be H. Equations ( 6) and [START_REF] Moore | An n job, one machine sequencing algorithm for minimizing the number of late jobs[END_REF] ensure that an operation o ij starts after s i and ends before d i . Equations ( 8) and ( 9) deal with the disjunctive constraints of operations on each machine and ensure that operations can be non-preemptively scheduled. If b j il = 1, then operation o ij starts before operation o lj on machine j, s ij + p ij ≤ s lj and s lj + p lj ≤ s ij + Ω. If, on the other hand, b j il = 0, then operation o lj starts before operation o ij on machine j, s ij + p ij ≤ s lj + Ω and s lj + p lj ≤ s ij .

The optimization version of the problem is obtained by setting the minimization of i∈N n i s i as objective function.

A constraint-based scheduling approach

Constraint programming has been widely applied in the area of scheduling [START_REF] Ph | Constraint-based scheduling, applying constraint programming to scheduling problems[END_REF]. In our constraint-based scheduling model, an activity of duration p ij is associated with each operation o ij and two variables s ij and c ij represent respectively its start time and its completion time. The constraint s ij + p ij = c ij is maintained (since we are in the non-preemptive case). Let Dom(V ) be the set of values which can be taken by a variable V . The minimum and maximum values of the domain of s ij are respectively denoted as:

-est ij = min v∈Dom(s ij ) v: the earliest starting time of operation o ij , -lst ij = max v∈Dom(s ij ) v: the latest starting time of operation o ij .
Similarly, the minimum and maximum values of the domain of c ij are respectively denoted as:

-eet ij = min v∈Dom(c ij ) v = est ij + p ij : the earliest completion time of operation o ij , -let ij = max v∈Dom(c ij ) v = lst ij + p ij : the latest completion time of operation o ij .
Thus, each operation o ij has to be executed within the time window [est ij , let ij ) whose bounds can be initialized with est ij = 0 and let ij = d i . Each machine j is modeled with a disjunctive unary resource, each operation o ij requiring the resource from its start time to its completion time. Edge-finding propagation techniques [START_REF] Ph | Constraint-based scheduling, applying constraint programming to scheduling problems[END_REF] which are able to adjust the time-windows of operations according to the disjunctive constraints between operations and their possible starting times are used. Each job i is also represented by an activity whose start time is denoted by the variable s i (its completion time being already fixed at d i ). In contrast to the activities corresponding to the operations, the durations of the activities corresponding to the jobs are not fixed, and these are denoted by the variable p i . The constraint s i + p i = d i is maintained. The sorting slots are represented by a discrete resource of capacity E. Each job i requires n i units of this resource from its start time to its completion time. Note that the domains of variables s i , p i , s ij and c ij will be reduced through the use of propagation algorithms as the search proceeds.

In the optimization problem an additional variable Φ represents the objective function Φ = n i=1 n i s i to be optimized. To find an optimal solution we use a branch and bound algorithm solving successive variants of the decision problem. At each iteration we try to improve the best known solution and to this end we add an additional constraint stipulating that Φ is lower than or equal to the best solution minus 1. At each new iteration the search resumes at the last-visited node of the previous iteration.

At each node of the enumeration procedure, the operation which can be scheduled at the latest time is scheduled at its latest start time. When the created branch fails, the operation is backward-postponed pending modification of its latest start time following the propagation of decisions. When all operations have been scheduled, the job which can start the earliest among the others, is scheduled at its earliest start time. When the created branch fails, the job is forward-postponed pending modification of its earliest start time following the propagation of decisions.

To improve the behavior of the enumeration procedure, we try to limit the enumeration of (partial) equivalent solutions. Consider a machine j, and suppose that 3 operations o xj , o yj , o zj are consecutively scheduled on machine j (i.e. c xj = s yj and c yj = s zj ). Suppose also that over the interval [s xj , c zj [ jobs x, y and z use sorting slots (i.e. s x ≤ s xj , s y ≤ s xj , s z ≤ s xj and d x ≥ c zj , d y ≥ c zj , d z ≥ c zj . Then clearly operations o xj , o yj and o zj can be permuted (in any order) and scheduled from the initial value of s xj without changing the validity (or the optimality, in the case of the optimization problem) of the solution (see Figure 7).

The following dominance rule can therefore be established:

Proposition 6
There is at least one (optimal) solution in which ∀j ∈ M and ∀x, y ∈ N such that c xj = s yj , at least one of the 3 conditions holds: Moreover, the start time of jobs x and y does not change. Therefore, if the initial solution was optimal, the obtained solution is also optimal.

This dominance rule is implemented as a propagation algorithm to ensure that only solutions belonging to the dominant subset of solutions are enumerated.

Models and appoaches for the preemptive case

In the previous section the start times of operations were computed simply in order to establish these operations could be scheduled. The values of the start times obtained were not really required. In the preemptive case this computation is unnecessary. Ensuring that no more than t k+1 -t k units of jobs are processed on each machine j ∈ M and within each interval [t k , t k+1 [ is sufficient to ensure that such a schedule exists. We now propose two models based on this property.

A flow-based mixed integer linear programming approach

In this model we propose a flow-based mixed integer linear programming. Just like in Section 5.1, this is implemented using binary variables for each job i and for each interval

[t k , t k+1 [. The variable x i,k is equal to 1 if job i uses n i sorting slots on [t k-1 ..t k [, otherwise it is 0.
The same constraints described by equations (1-4) are used on these variables to manage the constraint regarding the number of available sorting slots. In addition we use the set of variables {s i |i ∈ N } to represent the start time of jobs, along with the set of constraints described by Equation ( 5). Now, suppose that the values of variables x ik are known.The problem of the feasibility of operations on each machine j can then be seen as a kind of max flow problem. From this point of view, the network flow corresponding to each machine j is a connected, directed graph G(X, U ), where X is the set of vertices and U the set of edges. Set X is composed of a vertex Source, a vertex Sink, n vertices o ij with i ∈ {1, . . . , n} associated with the operations to be processed by machine j, and w vertices I k with k ∈ {0, . . . , w-1} corresponding to the consecutive intervals Figure 8).

I 0 = [t 0 , t 1 [, I 1 = [t 1 , t 2 [, . . . , I w-1 = [t w-1 , t w = H[ that partition interval [0, H[ (see
t 0 = 0 t k t k+1 tw = H I 0 I k = [t k , t k+1 [ I w-1
I w-1 p ij ∞ t k+1 -t k
Fig. 9 A network flow representation of the preemptive problem (each edge is labeled with its capacity).

Vertex Source is connected to each vertex o ij by an edge of capacity p ij . Each vertex I k is connected to vertex Sink by an edge of capacity t k+1 -t k . Finally, each vertex o ij is connected to vertex I k by an edge of infinite capacity if t k+1 ≤ d i (see Figure 9). It is easy to see that the preemptive scheduling of the operations on machine j is feasible if and only if the maximal flow which can be sent through the network is equal to i∈N p ij . This flow-max subproblem is modeled using variables f ijk , ∀i ∈ N, ∀j ∈ M, ∀k ∈ {0, . . . , w -1} corresponding to the quantity of operation o ij which is processed by machine j during interval [t k , t k+1 [. This subproblem is mathematically formulated by adding the following constraints:

f ijk ≤ x ik (t k+1 -t k ), ∀i ∈ N, ∀j ∈ M, ∀k ∈ {0, . . . , w -1} (11) k∈{0,...,w-1} f ijk = p ij , ∀i ∈ N, ∀j ∈ M (12) i∈N f ijk ≤ t k+1 -t k , ∀j ∈ M, ∀k ∈ {0, . . . w -1} (
y ik = 1, ∀i ∈ N ( 14
)
i∈N |d i >t k ∧d i ≤t l p ij l z=k y iz ≤ t l -t k , ∀j ∈ M, ∀k ∈ {0, . . . , w -1}, ∀l ∈ {k + 1, . . . , w} (15)
i∈N |d i >t k n i k l=0 y il ≤ E ∀k ∈ {0, . . . , w -1} (16)
Equation ( 14) ensures that job i begins to use sorting slots in only one interval. Equation (15) ensures that within each interval [t k , t l [ such that (t k , t l ) ∈ T 2 (with t k < t l ), the amount of operation which is processed by each machine j is less than or equal to t l -t k . Note that on machine j, the operations which are totally processed within interval [t k , t l [ belong to the jobs whose start times are in [t k , t l [ (i.e. those where l z=k y iz = 1) and whose completion times are in ]t k , t l ] (i.e. jobs i such that d i > t k ∧ d i ≤ t l ). Equation (16) ensures that the jobs which start at a date less than or equal to t k (i.e. the jobs i such that k l=0 y il = 1) and which end after t k (i.e. jobs i such that d i > t k ) use at most E sorting slots on [t k , t k+1 [. In conventional multiperiod models some term introduced in some constraint (period) will appear in all subsequent constraints (periods). In our model, the multiperiod aspect is to be found essentially in constraint ( 16), where each constraint indexed at k involves the previous one (indexed at k -1), subject to some condition on d i .

Recall that y ik = 1 means that job i will begin to use sorting slots at time t k . The optimization version of the problem is therefore obtained by adding the minimization of i∈N n i ( k∈{0,...,w-1} t k y ik ).

7 Experimental results

Generating benchmarks

All experimental results described in this paper were computed on a Dell Precision M2400 PC with an Intel(R) Core(TM) 2 Duo T9900 3.06Ghz, running Windows 7. The constraint programming model has been implemented on the top of IBM Ilog Scheduler 6.7. The MIP algorithm of IBM Ilog CPLEX 12.1 has been used to solve the mixed integer linear programming approaches. In aim to fairly compare the CP and MIP approaches and to obtain reproducible results, the parallelization option of Cplex has been disabled (Scheduler is not parallelized and parallelizing Cplex leads to important computational time variations when running on the same instances).

The methods have been successfully tested on a small number of industrial instances. However, it is well known that practical cases are often less of a challenge than benchmarks that have been generated specifically to test the limits of solvers. In this section we propose two schemes for generating instances. As shown in the following sections, these test instances would appear to be much harder to solve than instances encountered in practice.

Generating instances using the proof of complexity

In the first generation scheme, we propose drawing inspiration from the proof of the NP-completeness of the problem (see Section 3) with the 3-PARTITION problem for the non-preemptive version of the problem.

In this case, we have only one machine (m = 1). For given values of q and B, an instance has n = 5q jobs with 1 operation, and the number of sorting slots is equal to E = 3q. The time horizon is fixed at H = qB + 2q. For i ∈ {1, . . . , 3q} we set d i = H and n i = 1 and the processing time of the sole operation o i1 of job i is randomly generated in ]B/4, B/2[. For k ∈ {0, . . . , q -1}, we set

d 3q+2k+1 = (B + 2)k + 1, d 3q+2k+2 = (B + 2)k + 2, p 3q+2k+1 = p 3q+2k+2 = 1 and n 3q+2k+1 = n 3q+2k+2 = 3(q -k).
We generated two types of instance (feasible and unfeasible) using this scheme. In a first generation scheme (instances of this type will be referred to below as "3PART_F"), we make sure that there is a solution to the problem by randomly generating buckets of 3 operations, the sum of whose processing times is equal to B. Let lb and ub respectively be the smallest and largest integers in ]B/4, B/2[. For each bucket (a, b, c) of processing times, a is generated from a uniform distribution on interval [lb, min(B -2lb, ub)], then b is generated from a uniform distribution on interval [lb, min(B -a -lb, ub)], and finally c is set to B -a -b.

In a second generation scheme (instances of this type will be referred to below as "3PART_U"), the processing times are randomly generated from a uniform distribution on interval ]B/4, B/2[ (without applying a control as for instances of type 3PART_F). Then we keep only instances which are proved unfeasible by a solver.

For each type (3PART_F and 3PART_U), for each q ∈ {2, 4, 6} (i.e. n ∈ {10, 20, 30}) and for each B ∈ {28, 52, 100, 140} we generate 5 instances.

Generating more general instances

In this section we propose a more general generation scheme used both for the non-preemptive and the preemptive cases, and which can be either feasible or unfeasible (instances of this type will be referred to below as "GEN_F&U").

For each n ∈ {10, 20, 30, 40, 50}, m ∈ {1, 2, 5, 10}, we generated 1000 instances as follows.

For each instance, we randomly generate a number E of available sorting slots from the uniform distribution on interval [10, 20]. A maximum processing time P is randomly generated from the uniform distribution on interval [10, 50] and the processing times of operations o ij with i ∈ {1, . . . , n} and j ∈ {1, . . . , m} are then randomly generated from a uniform distribution on interval [0, P ]. We fix the time horizon H to (1 + γ) max j∈M { i∈N p ij }, where γ is randomly generated from a continuous distribution on interval [0, 0.25] and used to tighten or loosen the scheduling of operations on the bottleneck machine (i.e. the most loaded).

The completion time of the first job is fixed to H and its required number of sorting slots is generated from a decreasing triangular distribution on [1, 1 2 E], i.e. the required number of sorting slots can take a value in [1, 1 2 E], but the probability that any integer in this interval will be generated linearly decreases according to the value of that integer.

We then iteratively generate the completion time and the number of sorting slots required by the other jobs i as follows. Considering that jobs in {1, . . . i -1} (jobs with already computed departure times) are re-indexed in non-decreasing order of departure time, let k be the first index such that t k = max j∈{1,...m} p ij + k l=1 p lj ≤ d k . Thus t k corresponds to the earliest possible completion time that does not obviously lead to an unfeasible instance with respect to the departure times already generated (see for example the algorithm for solving 1|| U i in [START_REF] Moore | An n job, one machine sequencing algorithm for minimizing the number of late jobs[END_REF]). A completion time d i is then randomly generated from the increasing triangular distribution [t k , H] (higher values in the interval have a higher probability of being generated). Then, knowing that each job l ∈ {1, . . . , i -1} requires sorting slots at least from time d l -max j∈{1,...,m} p lj to time d l , we look at an upper bound of the maximum number

ni = E - max t∈[d i -max j∈{1,...,m} p ij ,d i ] {l|t∈[d l -max j∈{1,...,m} p lj ,d l ]} n l
of sorting slots which could be still available from time d i -max j∈{1,...,m} p ij to time d i (corresponding to the smallest interval of time over which job i will require its sorting slots). If ni = 0, another completion time d i is randomly generated as previously described. Otherwise, the number of sorting slots n i is randomly generated from a decreasing triangular distribution on interval [1, max(1, 1 2 ni )] (lower values in the interval have a higher probability of being generated).

Results for the non-preemptive case

In Table 1 results are shown for 3PART_F and 3PART_U instances tested with different methods proposed for the non-preemptive case. The results obtained with the mixed integer linear programming approach described in Section 5.1 are shown in column "MIP". To show the usefulness of the dominance rule described in Section 4, in column "ti_MIP" we also provide the results obtained by the same method but in a time-indexed version (that is without using the above dominance rule) where a variable x it is introduced for each t ∈ {0, 1, . . . , d i } meaning that job i uses sorting slots at time t if x it = 1. The results obtained with the constraintbased scheduling approach described in Section 5.2 are provided in column "CP". To show the usefulness of dominance rule described in Proposition 6, we also give (in column "CP_no_DR") the results obtained by the same method but without using the dominance rule. For 3PART_F and 3PART_U instances, we give results for the search problem. For 3PART_F instances, we also provide results for the optimization problem.

Note that when generating instances 3PART_U with n = 30, either the instances were feasible (and then not retained for the set of instances 3PART_U), or our approaches were not able to prove the feasibility (or the infeasibility) of the instances even after several hours of computation. This indicates that there is a clear limit to the size of unfeasible instances that can be solved within the time limit, and it explains why we provide results only up to n = 20 for instances 3PART_U.

The average number of generated nodes ("nodes") and the average computation time in milliseconds ("cpu") over the 5 generated instances are shown for each type of instance and each solution context (for 3PART_F instances, for each method, and for each couple (q, B)). A time limit of 100 seconds was used. The number of instances over 5 which were not solved within the time limit is given in column "uns". For unsolved instances, the number of nodes and the computation times were not taken into account in the computation of the averages. Then "_" in a column means that no instance over the 5 ones was solved within the time limit.

As expected, dominance rules are effective in reducing computation time and enabling a greater number of instances to be solved within the time limit. We can see that the CP method appears to be more efficient than the MIP method for these instances. We also remark that unfeasible instances are much more difficult to solve than the feasible ones, due to the proof of infeasibility to be established for these instances. Thus, from n = 20 and m = 1 onwards, the MIP method was unable to solve any of the unfeasible instances. Note that from n = 30 onwards we also begin to find unsolved instances for method CP.

In Table 2, we provide results obtained using the GEN_UF instances for the MIP and CP methods, for both the search and the optimization problems. For each method, for each resolution context, and for each couple (n, m), the average number of generated nodes ("nodes") and the average computation time in milliseconds ("cpu") over the generated instances are shown (1000 in the case of the search problem). Of course, whatever the method, the test on the optimization problem is only performed on the instances which have been proved to be feasible. Thus, column "nb" indicates the number of instances concerned. Once again, a time limit of 100 seconds was used. The percentage of instances which were not solved within the time limit is given in column "%uns". For unsolved instances, the number of nodes and the computational times were not taken into account in the computation of the averages. It will be remarked that the difficulty of the problem increases as the number n of jobs and the number m of machines increase. It is for this reason that we did not perform the tests for higher parameters m or n when more than 50% of the instances could not be solved within the time limit. In this case "_" in a column means that no test was performed for these parameters.

These results confirm that for the search problem the CP method outperforms the MIP method. However, this is not the case for the optimization problem, where MIP performs much better, as shown by the results for n = 10 and m = 5. Note that some instances with n = 10 and m = 5 remain unsolved for the search problem for both methods.

In Table 3 we compare the performance of the MIP and CP methods. For each method, we indicate the percentage of times where the method outperforms its rival (column "%best"), i.e. the percentage of times where the method found a solution either with a shorter computation time or where the other method failed to find a solution at all. We also indicate this proportion with respect to unfeasible instances ("%unf") and to feasible instances ("%feas"). We report the average difference in computation time (in milliseconds) when one method is better than the other ("cpu_dist"). Where the other method failed to find a solution a difference of 100s is taken. Column "eq" shows the percentage of instances for which both methods are equivalent. Finally "mix" shows what the results obtained would be, in terms of the percentage of unsolved instances ("%uns") and of the average computation time in milliseconds ("cpu"), if for each instance we were able to choose the best method. These results show that while CP is on average better than MIP, MIP dominates CP on a non-negligible set of instances, especially when the number of machines m is increased. For example, for n = 10 and m = 10, although neither method is able to solve more than 10% of the instances (see previous table), MIP is better than CP for 12% of the instances, taking on average 81 seconds less computational time than method CP. Then, taking the best method for each of these instances, only 0.4% of the instances still remain unsolved and the computation time is drastically reduced. This illustrates a certain complementarity between the two methods. From the results it is not possible to draw any conclusions regarding the best method and the feasibility or infeasibility of the instances.

Results for the preemptive case

For the preemptive case, Tables 4 and5 give the same measures for the flowbased mixed interger linear programming (FLOW) and the multiperiod (MULTIP) approaches methods as given above for the non-preemptive case (see tables 2 and 3). For the preemptive case, it would appear that the MULTIP method is better than the FLOW method. This is because this method uses fewer constraints with respect to the feasibility of the associated scheduling problem. We note that from n = 50 for the search problem and from n = 30 for the optimization problem, the problem becomes very hard to solve within the time limit. Table 5 (with the same colums as Table 3) reveals a complementarity between the two methods. When trying to determine the best method for each instance a certain balance can be observed. Note, however, that in cases where FLOW is better, the computation time for MULTIP is only marginally greater. The difference tends to be more acute for higher values of m and n. On the other hand, when MULTIP is better, the difference in computation time between the two methods tends to be large. Consequently, using the best method for each instance yields only slightly better results than using MULTIP exclusively. This enables us to solve only a small number of additional instances within the time limit, and reduces the computation time in only a limited number of instances and by a ratio of not more than 2. Note also that, for higher values of m and n, the FLOW method is of interest only for unfeasible instances.

Conclusion and extensions

In this paper we have formalized and studied a new combinatorial problem. While the study of this problem has a real application in warehouse management, its combinatorial structure is interesting from a theoretical point of view. We have shown the NP-completeness of the problem and we have proposed several methods for solving instances of the problem in both the preemptive and non-preemptive cases. Not only we can always determine which of the methods is better on average for the case in question, but we have also shown a certain complementarity between them. Moreover, although these methods are effective for instances encountered in practice, we have shown in our experimental study that the solvers can be challenged by theoretical instances of small dimensions. The technical part of the paper describes some initial approaches to solving the problems. As far as future research is concerned, these methods might be improved by working on the branching schemes of the methods and by attempting to find efficient lower bounds.

There are also several ways in which this research work might be extended. First, some other objective functions could be examined. For example, seeking a solution which balances between the different trucks the extra time available for using sorting slots might improve the robustness of solutions. Secondly, important questions arise when the instance in hand is not feasible. How can decision makers be assisted in making the right choices? Several possibilities might be investigated. For example, minimizing the (weighted) number of late jobs might help to determine which trucks might be canceled (ideally, as few as possible). Similarly, minimizing the total (weighted) tardiness can indicate which orders need to be canceled (again, as few as possible) when all truck departures are maintained.

While the approaches have been described in the context of solving a particular warehouse management problem, they can also be used to solve other optimization problems in which any particular resource is required during the processing of a set of operations on other resources, even if this processing is intermittent. For instance, a special case of multi-skill project management can fall in this category. Each project (job) is composed of tasks (operations) necessiting separated skills and resources (preparation shops). In top of this, one should affect some special resources (for instance, project manager, commercial agent, hotline number) from the start of the project until its delivery (a given deadline).

Fig. 1

 1 Fig. 1 Preparation of order lines in a warehouse

Fig. 3

 3 Fig. 3 Correct dynamic decisions by the WCS when the use of sorting slots is timetabled

  machine

Fig. 8

 8 Fig. 8 Intervals I 0 , I 1 , . . . , I w-1

  sorting slots starting at a date earlier than H -B -2 if no other job requires them.

	. . . . . . 3q + 2k +1	. . . 3q + 2k +2	part of {x, y, z}	3q + 2k +3	part of {x,y,z}	3q + 2k +4	machine

  Equivalent solutions when permuting operations o xj , o yj and o zj . Proof Consider a solution in which there exist two operations o xj and o yj scheduled on machine j such that y > x, s y ≤ s xj and d x ≥ c yj . It is possible to build an equivalent solution in which operation o yj is scheduled at the initial time s xj and where operation o xj is scheduled immediately after o yj and is completed at the initial value of c yj . Since in the new schedule operations o xj and o yj are scheduled in the interval during which their jobs are using sorting slots, and since the schedule of the other operations does not change, the obtained schedule remains valid.

	s xj			c zj
	. . . o xj	o yj	o zj	. . .	machine j
					sorting
					slots
	Fig. 7 1. x < y;
	2. s y > s xj ;
	3. d x < c yj .

  13) Equation (11) ensures that an operation o ij cannot be processed during an interval [t k , t k+1 [ in which job i does not use sorting slots. Equation (12) ensures that operation o ij is totally processed in [0, H[. Equation (13) ensures that machine j processes at most t k+1 -t k unit of operation on interval [t k , t k+1 [.The optimization version of the problem is obtained by adding the minimization of i∈N n i s i .6.2 A multiperiod approachIn this model we propose a kind of multiperiod approach. In the flow-based model described above, although we do not know exactly how the component parts of operations are scheduled in each interval [t k , t k+1 [, we nevertheless know what quantity of each operation o ij is scheduled. Actually, we can go further in the non-knowing thing by verifying only that the amount of operations processed by each machine j on each interval [t k , t l [ (with t k < t l ) is less than or equal to t l -t k .To implement this idea, we use binary variables y i,k , ∀i ∈ N, ∀k ∈ {0, . . . w -1} such that y ik is equal to 1 if job i begins to use sorting slots in interval [t k , t k+1 [ and to 0 in the other case. Note that y ik = 1 means that job i will use n i sorting slots in interval [t k , t [i] [. The search version of the problem can be mathematically formulated as follows:

	k∈{0,...,w-1}

Table 1

 1 Test of methods and dominance rules in the non-preemptive case on 3PART_U and 3PART_F instances search and optimization problems.

					3PART_U, search problem
			ti_MIP			MIP		CP_no_DR		CP
	q n B uns nodes	cpu uns nodes	cpu uns	nodes	cpu uns	nodes	cpu
	2 10 28	0	335	118	0	362	31	0	103	3	0	35	3
	52	0	293	165	0	302	28	0	101	9	0	34	3
	100	0	225	293	0	265	25	0	96	3	0	33	0
	140	0	198	389	0	247	25	0	95	3	0	33	3
	4 20 28	5	_	_	5	_	_	5	_	_	0	23349 1154
	52	5	_	_	5	_	_	0 2260940 44117	0	4767	315
	100	5	_	_	5	_	_	0 1174739 22326	0	2742	184
	140	5	_	_	5	_	_	0 419399 11681	0	1409	56
					3PART_F, search problem
			ti_MIP			MIP		CP_no_DR		CP
	q n B uns nodes	cpu uns nodes	cpu uns	nodes	cpu uns	nodes	cpu
	2 10 28	0	24	46	0	17	12	0	1	6	0	1	0
	52	0	14	81	0	26	9	0	3	3	0	3	3
	100	0	28	137	0	15	6	0	5	0	0	4	0
	140	0	33	196	0	28	9	0	2	0	0	1	0
	4 20 28	0 1861 20298	0 7967 2957	0	37592	727	0	302	12
	52	3 2342 22425	0 63116 22189	0	537	12	0	45	0
	100	3 1202 15319	1 16050 8213	0	16892	321	0	189	3
	140	2	411 38022	2 6438 13525	0	38673	742	0	437	53
	6 30 28	4	529 10498	4	203	234	1 190480 5815	1	933	43
	52	5	_	_	5	_	_	1	319	43	0 147328 4617
	100	5	_	_	5	_	_	1	65	19	0	29584	954
	140	5	_	_	5	_	_	3	5200	640	0	21674 1351
				3PART_F, optimization problem
			ti_MIP			MIP		CP_no_DR		CP
	q n B uns nodes	cpu uns nodes	cpu uns	nodes	cpu uns	nodes	cpu
	2 10 28	0	46	99	0	109	15	0	121	3	0	41	0
	52	0	42	143	0	10	9	0	113	0	0	39	0
	100	0	84	549	0	69	18	0	112	0	0	39	3
	140	0	112 1198	0	68	15	0	105	3	0	36	0
	4 20 28	5	_	_	0 16988 11606	5	_	_	0	67163 3981
	52	5	_	_	0	836 2003	5	_	_	0	18435	798
	100	5	_	_	2 58034 42765	3 2607247 48852	0	5572	393
	140	5	_	_	1 16983 22889	1 2641537 59198	0	5768	365
	6 30 28	5	_	_	4 25742 52260	5	_	_	5	_	_
	52	5	_	_	5	_	_	5	_	_	5	_	_
	100	5	_	_	5	_	_	5	_	_	3 1575292 50489
	140	5	_	_	5	_	_	5	_	_	3 1130425 61612

Table 2

 2 Test of methods MIP and CP (non-preemptive case) on GEN_F&U instances, search and optimization problems.

				MIP						CP
			search		optimization			search		optimization
	n m %uns nodes cpu nb %uns nodes cpu %uns nodes cpu nb %uns nodes cpu
	10 1	0	1085	120 707	0	1200	197	0	0.95	0 707	0	2940	32
	2	0.1	3430	477 698	0.72	4233	856	0	22	1 699	0 199104 2630
	5	3 10350 2611 656	9	8636 2709	2 19040 453 641	77 194517 4959
	10	11	9301 4251 646	16	7646 3335	10	9320 378 _	_	_
	20 1	41 20221 8518 348	51 14186 23616	0	126	4 455	63 1616711 25090
	2	50 11560 10006 _	_	_	_	0.3 19545 603 _	_	_	_
	5	54	3222 7079 _	_	_	_	22 47785 2412 _	_	_	_
	10	_	_	_ _	_	_	_	47 31122 2721 _	_	_	_

Table 3

 3 Comparing MIP and CP methods (non-preemptive case) on GEN_F&U instances, search problem.

			MIP				CP				mix
	n m %best %unf %feas cpu dist %best %unf %feas cpu dist %eq %uns cpu
	10 1	1	0	1	15	91	89	92	131	8	0 0.18
	2	0	0	0	70	95	89	97	606	5	0	1
	5	6	2	8	38376	90	85	92	5835	5	0	19
	10	12	7	15	81202	84	81	85	17250	4	0.4 235
	20 1	0	0	0	437	98	97	100	46917	2	0	3
	2	2	3	1	13739	92	88	99	58691	5	0.2 400
	5	16	22	14	47467	62	65	86	66877	8	14 2778

Table 4

 4 Test of methods FLOW and MULTIP (preemptive case) on GEN_F&U instances, search and optimization problems.

				FLOW						MULTIP
			search		optimization			search		optimization
	n m %uns nodes cpu nb %uns nodes cpu %uns nodes cpu nb %uns nodes cpu
	10 1	0	0.12	3 750	0	28	21	0	0.00	3 750	0	2	8
	2	0	0.11	4 730	0	21	23	0	0.00	5 730	0	1	9
	5	0	0.19	9 683	0	19	40	0	0.00	9 683	0	1	16
	10	0	0.21	19 716	0	19	76	0	0.04	17 716	0	1	29
	20 1	0	20	84 558	0	1172 6614	0	0.40	48 558	0	225	339
	2	0	56	419 474	1	1389 12317	0	0.40	75 474	0	366	759
	5	0.3	74 1790 402	10	1256 24843	0	3	179 402	0	502 1533
	10	2	77 4725 326	38	978 37469	0.1	4	424 326	0.31	557 3596
	30 1	3	235 4411 327	83	2181 44586	0	2	232 327	2	2174 6998
	2	8	229 9287 _	_	_	_	0	4	399 195	8	2488 11519
	5	26	80 12253 _	_	_	_	0	7 1147 117	29	2149 22029
	10	35	16 8146 _	_	_	_	0.2	3 2637 86	42	948 33311
	40 1	21	268 13270 _	_	_	_	0	3	602 131	21	3163 24774
	2	28	82 10916 _	_	_	_	0.1	3 1046 68	54	2109 32370
	5	41	7 5227 _	_	_	_	0	2 3577 _	_	_	_
	10	41	1 6516 _	_	_	_	0.6	1 11449 _	_	_	_
	50 1	40	115 12395 _	_	_	_	0	2 1270 42	81	1630 35564
	2	45	16 6623 _	_	_	_	0	3 2853 _	_	_	_
	5	42	0.30 5566 _	_	_	_	0.1	0.38 16007 _	_	_	_
	10	42	0.00 16380 _	_	_	_	23	0.00 40863 _	_	_	_

Table 5

 5 Comparing the FLOW and MULTIP methods (preemptive case) on GEN_F&U instances, search problem.

			FLOW			MULTIP			mix
	n m %best %unf %feas cpu dist %best %unf %feas cpu dist %eq %uns cpu
	10 1	20	14	22	15	20	16	21	15	60	0 0.02
	2	30	20	34	15	24	16	27	15	46	0 0.18
	5	31	38	28	15	31	27	33	16	38	0	4
	10	29	45	23	18	33	26	35	20	38	0	12
	20 1	59	64	55	32	28	24	31	196	13	0	30


Sonia Houssany-Pissot

Julien Rosencher

Philippe Allouch

Christophe Bensouda

Remy Pillière

Léa Cacoub

Christophe Caussin

Simon El-Hadad

Serge Makowski

Fabrizio Beverelli

Patrice Cacoub
email: patrice.cacoub@aphp.fr

Patrice Pr

MD, MSc Cacoub

Screening coronary artery disease with computed tomography angiogram should limit normal invasive coronary angiogram, regardless of pre-test probability

Keywords: invasive coronary angiogram, functional testing, computed tomography coronary angiogram, coronary artery disease. ICA: Invasive coronary angiography CAD: coronary artery disease CCTA: computed tomography angiogram PTP: pre-test probability CI: confidence interval

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

Screening coronary artery disease with computed tomography angiogram should limit normal invasive coronary angiogram, regardless of pre-test probability.

Introduction

Detection of stable coronary artery disease (CAD) remains a major public health problem. Most guidelines recommend a stepwise approach for decision making in patients with suspected stable CAD. The process begins (Step 1) with a clinical assessment of the probability that CAD is present in a particular patient, i.e. determination of pre-test probability (PTP). This approach uses clinicians' PTP of CAD along with the results of diagnostic tests to generate individualized post-test CAD probabilities for a given patient. The PTP is influenced by the prevalence of the disease in the population studied, as well as clinical features (including the presence of cardiovascular risk factors). Major determinants of PTP are age, gender and the nature of symptoms (1). Step 1 is followed by non-invasive testing to establish the diagnosis of CAD (Step 2). Usually, based on available non-invasive test results, patients who may benefit from invasive coronary angiogram (ICA) are selected. Recent guidelines recommend that patients with presumed stable CAD and intermediate PTP should be screened with non-invasive functional testing before undergoing ICA [START_REF] Fihn | ACC/AHA/AATS/PCNA/SCAI/STS focused update of the guideline for the diagnosis and management of patients with stable ischemic heart disease: a report of the American College of Cardiology/American Heart Association Task Force on Practice Guidelines, and the American Association for Thoracic Surgery, Preventive Cardiovascular Nurses Association, Society for Cardiovascular Angiography and Interventions, and Society of Thoracic Surgeons[END_REF][START_REF] Patel | Low diagnostic yield of elective coronary angiography[END_REF]. One of the guideline's objectives is to limit the number of unnecessary ICA, i.e. the number of patients who would undergo an ICA that would not lead to the diagnosis of CAD. Such strategy appears suboptimal as more than half of the patients undergo ICA without previous functional testing and more than one-third of ICA performed are normal [START_REF] Lin | Frequency of stress testing to document ischemia prior to elective percutaneous coronary intervention[END_REF]. In the very recently published ESC 2019 guidelines for CAD, the PTP of CAD based on age, gender and nature of symptoms have undergone major revisions, notably as they introduced a new phrase "Clinical likelihood of CAD that utilizes also various risk factors of CAD as PTP modifiers" [START_REF] Knuuti | ESC Scientific Document Group. 2019 ESC Guidelines for the diagnosis and management of chronic coronary syndromes[END_REF]. This is precisely what we have done in the present study, i.e. PTP based on clinical likelihood.

Coronary computed tomography angiogram (CCTA) has emerged as a non-invasive anatomical imaging to detect CAD. Several studies have reported high levels of sensitivity and specificity of CCTA compared to the gold standard [START_REF] Knuuti | ESC Scientific Document Group. 2019 ESC Guidelines for the diagnosis and management of chronic coronary syndromes[END_REF][START_REF] Miller | Diagnostic performance of coronary angiography by 64row CT[END_REF][START_REF] Budoff | Diagnostic performance of 64-multidetector row coronary computed tomographic angiography for evaluation of coronary artery stenosis in individuals without known coronary artery disease: results from the prospective multicenter ACCURACY (Assessment by Coronary Computed Tomographic Angiography of Individuals Undergoing Invasive Coronary Angiography) trial[END_REF]. The ability of CCTA to diagnose CAD in patients with suspected stable angina compared to stress-test has been demonstrated [START_REF] Investigators | CT coronary angiography in patients with suspected angina due to coronary heart disease (SCOT-HEART): an open-label, parallel-group, multicentre trial[END_REF]. Many studies have reported excellent operating characteristics of CCTA for CAD diagnosis and a positive impact on referrals for ICA, although multicentre studies have yielded mixed results [START_REF] Meijboom | Diagnostic accuracy of 64-slice computed tomography coronary angiography: a prospective, multicenter, multivendor study[END_REF][START_REF] Hoffmann | Prognostic value of noninvasive cardiovascular testing in patients with stable chest pain: insights from the PROMISE trial (Prospective Multicenter Imaging Study for Evaluation of Chest Pain)[END_REF]11,[START_REF] Haase | Diagnosis of obstructive coronary artery disease using computed tomography angiography in patients with stable chest pain depending on clinical probability and in clinically important subgroups: meta-analysis of individual patient data[END_REF]. However, in current guidelines, CCTA has been given limited place as an alternative to traditional functional testing, when ischemia inducing test is contraindicated or non-conclusive (1,2). To better estimate the diagnostic accuracy of CCTA as a first step-test to screen patients suspected to have a stable CAD in real-life practice, we performed a multicentre study comparing the prevalence rate of CAD on ICA in patients with either positive CCTA or positive functional tests, according to the pre-test probability value.

Methods

Patients

We retrospectively analyzed consecutive patients suspected of stable CAD, admitted for an ICA in three major interventional cardiology centres in Paris area, France between January 1 st , 2014 and June 1 st , 2018 (Ambroise Paré Clinic, Neuilly Sur Seine; Lagny Marne-la-Vallée Hospital, Lagny sur Marne; and Montsouris Mutualist Institute, Paris, FRANCE).

For the present analysis, patients should have had a positive functional testing or a positive CCTA indicating an ICA. Positive functional testing was defined by ischemia findings during stress or recovery, like patient chest pain, ECG modifications, left ventricle ejection fraction decrease, abnormal cinetic wall motion, and abnormal myocardial perfusion. CCTA findings were deemed positive if coronary artery stenosis ≥ 50% was reported, if the stenosis calcification was classified as severe, or if the coronary artery calcium score considering the Agatston method was too high (i.e. above 400) [START_REF] Mcclelland | Distribution of coronary artery calcium by race, gender, and age: results from the Multi-Ethnic Study of Atherosclerosis (MESA)[END_REF]. Data were prospectively collected on patient's demographic characteristics, cardiovascular risk factors, history of chronic renal failure, stroke and peripheral vascular disease, symptoms and non-invasive tests motivating ICA. For all included patients, we calculated the PTP according to the Consortium Basic Score [START_REF] Bittencourt | European Society of Cardiology-Recommended Coronary Artery Disease Consortium Pretest Probability Scores More Accurately Predict Obstructive Coronary Disease and Cardiovascular Events Than the Diamond and Forrester Score: The Partners Registry[END_REF]. Based on guidelines recommendations, patients were categorized in one of the 5 PTP following groups: [START_REF] Montalescot | ESC guidelines on the management of stable coronary artery disease: the Task Force on the management of stable coronary artery disease of the European Society of Cardiology[END_REF] (1) obstructive CAD (oCAD, refers to diameter stenosis on ICA), as an estimated stenosis of 50% or more in any major epicardial vessels, including side branches of at least 2 mm in diameter; (2) non-obstructive CAD (non-oCAD), as an estimated stenosis below 50% in any major epicardial vessels, including side branches of at least 2mm in diameter; and (3) normal ICA (no CAD), as estimated stenosis of less than 20% in any vessel. Most guidelines use a threshold of 50% to define a significant coronary artery stenosis [START_REF] Raff | Society of Cardiovascular Computed Tomography. SCCT guidelines for the interpretation and reporting of coronary computed tomographic angiography[END_REF], and the tighter the stenosis found the higher is the risk of ischemia [START_REF] Leipsic | SCCT guidelines for the interpretation and reporting of coronary CT angiography: a report of the Society of Cardiovascular Computed Tomography Guidelines Committee[END_REF]. All data were recorded prospectively before and just after ICA in the CardioReport™ software database (MediReport Ltd, Paris, France).

Statistical analysis

Continuous variables are presented with the median and interquartile range (IQR); categorical variables are presented with counts and proportions. We compared baseline demographic characteristics, cardiovascular risk factors, history of chronic renal failure, stroke and peripheral vascular disease and results of ICA of patients who previously had either a positive functional testing or a positive CCTA. We compared the main results of ICA (normal ICA vs. oCAD) according to the PTP value of CAD. Categorical variables were compared using the chi-square test. P values of less than 0.05 were considered to indicate statistical significance. Statistical analyses were performed using R software version 3.2.4 (The R Foundation for Statistical Computing; Vienna, Austria).

No extramural funding was used to support this work. The authors are solely responsible for the design and conduct of this study, all study analyses, the drafting and editing of the paper and its final contents.

Results

From 7465 consecutive patients suspected of stable CAD who had an ICA in our institutions during the study period, we excluded 2513 patients who have had neither functional testing nor CCTA prior to ICA. Therefore, we analyzed data of 4952 patients who underwent ICA following either a positive functional testing (3276, 66.2%) or a positive CCTA (1676, 33.8%). Prior positive functional tests included single photon emission computed tomography (SPECT) (1335 patients, 27%), exercise stress test (1179, 23.8%), stress echocardiography (603, 12.2%), and stress magnetic resonance imaging (MRI) (159, 3.2%). To estimate the pre-ICA probability of having a CAD, patients were categorized in the previously defined 5 PTP groups : [START_REF] Montalescot | ESC guidelines on the management of stable coronary artery disease: the Task Force on the management of stable coronary artery disease of the European Society of Cardiology[END_REF] 29.3%, p<0.0001) or a chronic kidney disease (0.1% vs. 1.5%, p<0.0001). The population of the CCTA group had a slightly lower median pretest probability of CAD compared to the functional test group (25.4% vs 26.5%, p=0.03). Among the symptoms motivating the non-invasive tests for searching a CAD, patients in the CCTA group compared to the functional test group had more frequently atypical angina (22.3% vs. 17.2%, p<0.0001) and less frequently non-specific symptoms (46.4% vs. 51.5%, p =0.0006).

Main results of ICA showed, considering the entire population, a normal coronary angiogram (819 patients, 16.5%), non-oCAD (1193 patients, 24.1%) or oCAD (2940 patients, 59.4%) (Table 2). Our primary objective was to analyze the rates of normal ICA in the entire population according to the results of functional testing vs. CCTA, then according to the PTP values. When we analyzed patients without taking into account the PTP group, CCTA always showed significantly better results than functional tests, i.e. less frequently normal ICA (7% vs. 16.5%), and more frequently non-oCAD (27.9% vs. 22.2%) and oCAD (65.1% vs. 56.4%)(all p<0.0001)(Table 2, Figure 1). When we considered the PTP values in the entire population, the rate of normal ICA ranged from 14.7% to 27.3% with lower rates in the groups with higher PTP values (Table 2). When we considered the different PTP values, CCTA always showed lower rates of normal ICA than the functional testing. Subjects with a very low risk of CAD (PTP < 15%) and a positive functional test prior to ICA had 4.7 times higher rate of normal ICA compared to those with a positive CCTA (95% CI=3.2-7.0, p<0.001)(Figure 2, suppl. Table 2). Subjects with lower intermediate-risk (PTP 15-35%) or higher intermediate-risk (PTP 35-50%) who had a positive functional test had 4.0 and 2.4 times higher rate of normal ICA compared to subjects with a positive CCTA, respectively (95% CI=2.7-5.8 and 95% CI=1.5-3.9, respectively; p<0.01). Subjects at high-risk (PTP 50-65%) and very high-risk of CAD (PTP > 65%) who had a positive functional test prior to ICA had 1.9 and 2.3 times the rate of normal ICA compared to subjects with a positive CCTA, respectively (95% CI=1.2-2.9 and 95% CI=1.4-3.6, respectively; p<0.01). Because of the known lowest performance of exercise tests in the diabetic population, we did the same analyses after excluding exercise ECG which showed very similar results (Supplementary Figures 1 and2). For the screening of oCAD, the ratio of positive functional testing/positive CCTA was significantly inferior to 1 only for the PTP groups <15% and 15-35% (p<0.001) (Figure 2, suppl. Table 2). CCTA, compared to functional tests significantly better predicted the presence of CAD (non-oCAD and oCAD) on ICA in patients with very low and low PTP values (i.e. between 0% and 35%). For patients with higher PTP values (>35%), there was a trend for better prediction efficacy for CCTA versus functional tests although not reaching statistical significance (Table 2).

Discussion

This study shows that, in a real-life setting, patients suspected to have a stable CAD and a positive non-invasive test had lower rates of normal ICA in the CCTA group compared to the functional testing group, regardless of the PTP value. Having a positive CCTA prior to ICA led less frequently to a normal ICA compared to having positive functional testing (from 4.0 times less frequently when PTP<15% to 2.3 times when PTP>65%). The rate of obstructive CAD was significantly higher in patients who had a positive CCTA compared to those with positive functional testing only for low-risk and lower intermediate risk patients (PTP<35%).

The major goal of performing a non-invasive test before ICA is to avoid an unnecessary invasive test. Actual guidelines recommend functional testing in the first intention in patients with intermediate PTP (15-65%) [START_REF] Montalescot | ESC guidelines on the management of stable coronary artery disease: the Task Force on the management of stable coronary artery disease of the European Society of Cardiology[END_REF][START_REF] Fihn | ACC/AHA/AATS/PCNA/SCAI/STS focused update of the guideline for the diagnosis and management of patients with stable ischemic heart disease: a report of the American College of Cardiology/American Heart Association Task Force on Practice Guidelines, and the American Association for Thoracic Surgery, Preventive Cardiovascular Nurses Association, Society for Cardiovascular Angiography and Interventions, and Society of Thoracic Surgeons[END_REF]. Because of its excellent negative predictive value, CCTA appeared as a good alternative in patients with a lower intermediate-risk PTP [START_REF] Budoff | Diagnostic performance of 64-multidetector row coronary computed tomographic angiography for evaluation of coronary artery stenosis in individuals without known coronary artery disease: results from the prospective multicenter ACCURACY (Assessment by Coronary Computed Tomographic Angiography of Individuals Undergoing Invasive Coronary Angiography) trial[END_REF][START_REF] Miller | Diagnostic performance of coronary angiography by 64row CT[END_REF]. In a recent meta-analysis, Knuuti and al. reported that CCTA was superior to rule out the presence of a CAD when PTP value was intermediate [START_REF] Miller | Diagnostic performance of coronary angiography by 64row CT[END_REF]. In this study, CCTA was superior to functional tests to avoid unnecessary ICA, regardless the PTP value, and this was mainly due to a lower rate of normal ICA (7.0% vs. 21.4%, respectively p<0.001). In a retrospective study evaluating symptomatic and asymptomatic patients without known CAD undergoing ICA, Patel et al.

showed an overall rate of 39.2% of normal ICA (3). Non-obstructive coronary angiogram or normal ICA was present in 40.6% in our population (in 34.9% of positive CCTA group vs. 43.6% of positive functional testing group, p <0.001). The rate of non-obstructive CAD was not statistically different between these groups (CCTA 27.9% vs. functional testing 22.2%, p=0.12). In the Prospective Multicenter Imaging Study for Evaluation of Chest Pain (PROMISE) Study, evaluating the outcomes of anatomical study with CCTA vs. functional study for patients with symptoms suggestive of CAD, the rate of non-oCAD was 28% and 52% in the CCTA and functional testing group, respectively [START_REF] Douglas | Outcomes of anatomical versus functional testing for coronary artery disease[END_REF].

The recommended strategy for CAD screening is based on the determination of PTP before choosing further investigations. Despite technical improvements and the development of radial access, ICA remains an invasive method with a rate of complications estimated between 0.6% and 4.2% [START_REF] Pristipino | Major improvement of percutaneous cardiovascular procedure outcomes with radial artery catheterisation: results from the PREVAIL study[END_REF][START_REF] Nakazato | Age-related risk of major adverse cardiac event risk and coronary artery disease extent and severity by coronary CT angiography: results from 15 187 patients from the International Multisite CONFIRM Study[END_REF]. However, retrospective data showed that 17.1% of low-risk patients and 15.9% of intermediate risk patients are directly referred to ICA without undergoing a non-invasive test [START_REF] Patel | Low diagnostic yield of elective coronary angiography[END_REF]. Up to 45% of patients are directly referred for percutaneous coronary intervention without prior functional testing (4). European guidelines recommend using the CAD consortium score while American guidelines recommend using the Diamond and Forrester score [START_REF] Montalescot | ESC guidelines on the management of stable coronary artery disease: the Task Force on the management of stable coronary artery disease of the European Society of Cardiology[END_REF][START_REF] Fihn | ACC/AHA/AATS/PCNA/SCAI/STS focused update of the guideline for the diagnosis and management of patients with stable ischemic heart disease: a report of the American College of Cardiology/American Heart Association Task Force on Practice Guidelines, and the American Association for Thoracic Surgery, Preventive Cardiovascular Nurses Association, Society for Cardiovascular Angiography and Interventions, and Society of Thoracic Surgeons[END_REF]. The latter score seems to overestimate the risk of obstructive CAD while CAD consortium score seems to improve the PTP estimation [START_REF] Bittencourt | European Society of Cardiology-Recommended Coronary Artery Disease Consortium Pretest Probability Scores More Accurately Predict Obstructive Coronary Disease and Cardiovascular Events Than the Diamond and Forrester Score: The Partners Registry[END_REF]. In the present study, 50.7% of patients in the very low-risk PTP group (58.9% in CCTA group vs. 45.6% in the functional testing group) had a non-oCAD or an oCAD. These data may question the relevance of these scores since guidelines consider that when PTP is under 15%, no further investigations are recommended. In the non-invasive strategy of CAD detection, using CCTA at step 1 appears an interesting option as it should permit to have two times less normal ICA without missing abnormal ICA (including oCAD and non-oCAD)[Figure 1]. Interestingly, the latest ESC guidance continues to suggest CCTA should be for lowintermediate risk and functional testing for intermediate-higher risk [START_REF] Knuuti | ESC Scientific Document Group. 2019 ESC Guidelines for the diagnosis and management of chronic coronary syndromes[END_REF]. The current results clearly discounts such approach and reinforce the United Kingdom National Institute for Clinical Excellence guidance which does not recommend the assessment of pretest probabilities but that CCTA should be the first line test of choice for patients with possible angina [START_REF] Kelion | The rationale for the primacy of coronary CT angiography in the National Institute for Health and Care Excellence (NICE) guideline (CG95) for the investigation of chest pain of recent onset[END_REF].

Data regarding whether CCTA leads more frequently to ICA are conflicting. Some studies showed that in patients with suspected stable CAD, compared to a SPECT or cardiac magnetic resonance, a positive anatomical test leads more frequently to ICA [START_REF] Investigators | CT coronary angiography in patients with suspected angina due to coronary heart disease (SCOT-HEART): an open-label, parallel-group, multicentre trial[END_REF][START_REF] Douglas | Outcomes of anatomical versus functional testing for coronary artery disease[END_REF][START_REF] Hachamovitch | Patient management after noninvasive cardiac imaging results from SPARC (Study of myocardial perfusion and coronary anatomy imaging roles in coronary artery disease)[END_REF][START_REF] Pontone | The STRATEGY Study (Stress Cardiac Magnetic Resonance Versus Computed Tomography Coronary Angiography for the Management of Symptomatic Revascularized Patients): Resources and Outcomes Impact[END_REF]. Similar results were observed in patients presenting with a low risk of acute coronary syndrome [START_REF] Linde | Cardiac computed tomography guided treatment strategy in patients with recent acute-onset chest pain: results from the randomised, controlled trial: CArdiac cT in the treatment of acute CHest pain (CATCH)[END_REF]. After a CCTA, the rate of revascularization was higher whereas the rate of myocardial infarction was lower as compared to functional testing [START_REF] Jørgensen | Functional Testing or Coronary Computed Tomography Angiography in Patients With Stable Coronary Artery Disease[END_REF]. In the SCOT-HEART study, the prognosis of patients screened for stable angina was better after CCTA compared to functional testing without resulting in increasing referral for ICA or coronary revascularization [START_REF] Investigators | Coronary CT Angiography and 5-Year Risk of Myocardial Infarction[END_REF]. However, CCTA did show an early increase of coronary revascularisation, but beyond the first year there was significantly less revascularisation. This confirms that, by using at a first step CCTA the right patients were identified earlier with an increased number of revascularisation because more oCAD was found, as in the present study. This should prevent patients progressing to needing downstream revascularisation, and hence less late revascularisation.

Recently, the Coronary Computed Tomographic Angiography for Selective Cardiac Catheterization (CONSERVE) trial showed that patients referred to ICA after a recommendation-based screening strategy compared to a systematic CCTA screening, had similar outcomes at one year, without increasing the number of ICA or revascularization [START_REF] Chang | Selective Referral Using CCTA Versus Direct Referral for Individuals Referred to Invasive Coronary Angiography for Suspected CAD: A Randomized, Controlled, Open-Label Trial[END_REF]. The authors suggested a revision of the current CAD management guidelines, consistently with the United Kingdom National Institute for Clinical Excellence guidelines [START_REF] Moss | The Updated NICE Guidelines: Cardiac CT as the First-Line Test for Coronary Artery Disease[END_REF].

Our study has several limitations. First, it was a retrospective study. So the comparison between anatomical and functional testing was not based on randomized inclusion. Second, as the study was not randomized, the proportions of each non-invasive functional testing were not equal. However, because the functional testing group were higher risk, it is even more surprising that CTCA managed to have lower rates of normal angiograms. Many exercise test and SPECT were realized in the CCTA and functional testing group whereas the number of stress tests (echocardiography and MRI) was low. This may have underestimated the interest of stress-testing, although it also reflects the small number of stress-MRI realized in a real-life setting. Patients in the functional testing group had slightly higher median PTP value and a more severe profile as they were more frequently male, suffered more frequently from diabetes, chronic kidney disease and had a higher BMI. In this study, we classified patients into 3 groups (oCAD, non-oCAD and normal ICA) instead of 2 groups (oCAD vs. non-oCAD).

Previous studies showed that in multivariable analyses, an increased risk of major adverse cardiovascular events at 2 years was observed for both non-oCAD and oCAD when compared with patients with normal CCTA [START_REF] Bittencourt | European Society of Cardiology-Recommended Coronary Artery Disease Consortium Pretest Probability Scores More Accurately Predict Obstructive Coronary Disease and Cardiovascular Events Than the Diamond and Forrester Score: The Partners Registry[END_REF]. Finally, exercise testing and as SPECT provide other functional information not given by CCTA such as exercise capacity, arrhythmias evaluation or blood pressure follow-up.

Conclusion

In a real-life condition, this study suggests that CCTA is a better alternative than functional testing to limit unnecessary ICA regardless of PTP value, as it should permit to have two times less normal ICA without missing abnormal ICA. This finding should weigh in the decision making to perform CCTA for screening patients suspected of stable CAD.

Funding: none.

Figure 1 .

 1 Figure 1. Results of invasive coronary angiogram according to the positivity of previous non-invasive tests for coronary artery

Figure 2 .

 2 Figure 2. Normal coronary angiogram or obstructive coronary artery disease rate ratio in patients with prior positive functional

Figure 1 .

 1 Figure 1. Results of invasive coronary angiogram according to the positivity of

Figure 2 .

 2 Figure 2. Normal coronary angiogram or obstructive coronary artery disease

Table 1 .

 1 Compared to the functional test group, patients in the CCTA group were found to be older (66.6 vs. 65.8 years, p=0.01), less frequently male (69.7 % vs. 76.5%, p<0.0001), and had less frequently a diabetes (18.4% vs.

	low risk [PTP < 15%; n=968,19.5%], (2) lower intermediate risk [PTP 15 to 35%;

Table 1 . Baseline characteristics of patients suspected of coronary artery disease after a positive coronary computed tomography angiogram or a positive functional test for cardiac ischemia.

 1

		Total	CCTA	Functional test
	Variables	N = 4952	N = 1676	N = 3276
					p-value
		(100%)	(33.8%)	(66.2%)
	Patient characteristics
	Age (years), mean (SD)	66.1 (10.3)	66.6 (10.0)	65.8 (10.4)	0.01
	Male gender, n (%)	3674 (74.2%)	1168 (69.7%)	2506 (76.5%)	< 0.001
	BMI kg/m 2 , median, [IQR]	26.2 [24.0,	26.1 [23.8, 29.4] 26.3 [24.1, 29.4]	0.02
		29.4]
	Cardiovascular risk factors
	Diabetes mellitus, n (%)	1270 (25.6%)	309 (18.4%)	961 (29.3%)	< 0.001
	Dyslipidemia, n (%)	2363 (47.7%)	791 (47.2%)	1572 (48.0%)	0.62
	Hypertension, n (%)	2194 (44.3%)	759 (45.3%)	1435 (43.8%)	0.34
	Current smoking, n (%)	2025 (40.9%)	698 (41.6%)	1327 (40.5%)	0.46
	Chronic kidney disease, n (%)	50 (1.0%)	2 (0.1%)	48 (1.5%)	< 0.001
	Peripheral vascular disease, n	312 (6.3%)	91 (5.4%)	221 (6.7%)	0.08
	(%)
	Previous stroke, n (%)	111 (2.2%)	35 (2.1%)	76 (2.3%)	0.67

PTP, pre-test probability; BMI, body mass index; SD, standard deviation; IQR, interquartile range; SPECT, single photon emission computed tomography; MRI, magnetic resonance imaging; ICA, invasive coronary angiogram; CAD, coronary artery disease; CCTA, computed tomography angiogram; NA,

Table 2 . Results of invasive coronary angiogram in patients with positive functional testing or positive coronary computed tomography angiogram, according to the pretest probability of coronary artery disease.

 2

		Normal coronary angiogram		Non-obstructive coronary artery disease	Obstructive coronary artery disease
				Functiona
			CCTA		P		CCTA	Functional	P		CCTA	Functional	P
		Total		l		Total				Total
								testing				testing
				testing
		16.5%	7.0%	16.5%	<0.00	24.1%	467	726	<0.00	59.4%	65.1%		<0.00
	All PTP	(819)	(118)	(819)	1	(1193)	(27.9%)	(22.2%)	1	(2940)	(1091)	56.4% (1849)	1
	PTP	27,3%	7.0%	33.1%	<0.00	26.2%	127		<0.00		58.9%		<0.00
	0-15%	(223)	(26)	(197)	1	(254)	(34.1%)	127 (21.3%)	1 22,5% (710)	(219)	45.6% (491)	1
	PTP	26,3%	5.6 %	22.3 %	<0.00	22.0%	122				69.8 %		<0.00
	15-35%	(215)	(28)	(187)	1	(294)	(24.5%)	172 (20.7%)	0.097 26,2% (827)	(347)	57.2 % (480)	1
	PTP	14,7%	7.5 %	18.3 %	<0.00	23.6%	68				65.6 %
	35-50%	(120)	(19)	(101)	1	(190)	(26.9%)	122 (22.1%)	0.160 15,7% (496)	(166)	59.7 % (330)	0.13

CCTA, computed tomography angiogram; PTP, pre-test probability

Table 1 : Baseline characteristics of patients suspected of coronary artery disease after a positive coronary computed tomography angiogram or a positive functional test for cardiac ischemia, according to each center.

 1

		Center 1	Center 2	Center 3	Total	p-
		(n=2184)	(n=1466)	(n=1302)	(n=4952)	valu
						e
	Patient
	characteristics
	Age (years), mean	66.5 (10.1)	65.3 (10.0)	66.2 (10.7)	66.1 (10.3)	<0.00
	(SD)					1
	Male gender, n (%)	1623 (74.3%)	1116 (76.1%)	935 (71.9%)	3674 (74.2%)	0.04
	BMI kg/m 2 , median,	25.9 [23.8,	27.1 [24.5,	26.0 [23.7,	26.2 [24.0,	<0.00
	[IQR]	28.7]	30.8]	28.8]	29.4]	1
	Diabetes mellitus, n	443 (20.3%)	479 (32.7%)	348 (26.7%)	1270 (25.6%) <0.00
	(%)					1
	Dyslipidemia, n (%)	955 (43.7%)	748 (51.0%)	660 (50.7%)	2363 (47.7%) <0.00
						1
	Hypertension, n (%)	965 (44.2%)	659 (45.0%)	570 (43.8%)	2194 (44.3%)	0.82
	Current smoking, n	879 (40.2%)	677 (46.2%)	469 (36.0%)	2025 (40.9%) <0.00
	(%)					1
	CKD, n (%)	13 (0.6%)	28 (1.9%)	9 (0.7%)	50 (1.0%)	<0.00
						1
	Vascular disease, n	112 (5.1%)	128 (8.7%)	72 (5.5%)	312 (6.3%)	<0.00
	(%)					1
	Previous stroke, n	40 (1.8%)	50 (3.4%)	21 (1.6%)	111 (2.2%)	<0.00
	(%)					1
	PTP, median [IQR]	38.0 [18.5,	38.9 [19.2,	37.0 [18.5,	38.0 [18.9,	0.37
		60.7]	61.2]	58.1]	60.4]
	Positive pre-ICA					<0.00
	tests					1
	SPECT	402 (18.4%)	693 (47.3%)	240 (18.4%)	1335 (27.0%)
	Exercise stress	525 (24.0%)	210 (14.3%)	444 (34.1%)	1179 (23.8%)
	test
	S t r e s s	269 (12.3%)	97 (6.6%)	237 (18.2%)	603 (12.2%)
	echocardiography

Table 2 .

 2 Rate ratio of positive functional testing/positive computed tomography angiogram in patients suspected of stable coronary artery disease, according to the pre-test probability of coronary artery disease and to the results of coronary

	angiogram.

Table 1 . Baseline characteristics of patients suspected of coronary artery disease after a positive coronary computed tomography angiogram (CCTA) or a positive functional test for cardiac ischemia.

 1

	Variables	Total N = 4952 (100%)	CCTA N = 1676 (33.8%)	Functional test (66.2%) N = 3276	p-value
	Patient characteristics
	Age (years), mean (SD)	66.1 (10.3)	66.6 (10.0)	65.8 (10.4)	0.01
	Male gender, n (%)	3674 (74.2%)	1168 (69.7%)	2506 (76.5%)	< 0.001
	BMI kg/m 2 , median, [IQR]	26.2 [24.0,	26.1 [23.8,	26.3 [24.1,	0.02
		29.4]	29.4]	29.4]
	Cardiovascular risk factors
	Diabetes mellitus, n (%)	1270 (25.6%)	309 (18.4%)	961 (29.3%)	< 0.001
	Dyslipidemia, n (%)	2363 (47.7%)	791 (47.2%)	1572 (48.0%)	0.62
	Hypertension, n (%)	2194 (44.3%)	759 (45.3%)	1435 (43.8%)	0.34
	Current smoking, n (%)	2025 (40.9%)	698 (41.6%)	1327 (40.5%)	0.46
	Chronic kidney disease, n	50 (1.0%)	2 (0.1%)	48 (1.5%)	< 0.001
	(%)
	Peripheral vascular disease,	312 (6.3%)	91 (5.4%)	221 (6.7%)	0.08
	n (%)
	Previous stroke, n (%)	111 (2.2%)	35 (2.1%)	76 (2.3%)	0.67
	PTP, median [IQR]	26.2 [14.0,	25.4 [13.5,	26.5 [14.4,	0.03
		46.2]	44.3]	47.1]
	Symptoms
	Typical Angina	1549 (31.3%)	525 (31.3%)	1024(31.3)	0.99
	Atypical angina	938 (18.9%)	374 (22.3%)	564 (17.2%)	< 0.001
	Non-specific symptoms	2465 (49.8%)	777 (46.4%)	1688 (51.5%)	< 0.001
	Prior positive function test for CAD
	SPECT		NA	1335 (40.7%)	-
	Exercise stress test		NA	1179 (36.0%)	-
	Stress echocardiography		NA	603 (18.4%)	-
	Stress MRI		NA	159 (4.9%)	-

Table 2 . Results of invasive coronary angiogram in patients with positive functional testing or positive coronary computed tomography angiogram (CCTA), according to the pretest probability (PTP) coronary artery disease. Normal coronary angiogram Non-obstructive coronary artery disease Obstructive coronary artery disease

 2

		Tota l	CCT A	Funct ional testin g	P Total	CCTA	Functi onal g testin	P	Total	CCTA	testing Functio nal	P
		16.5			<0								<0
		%	7.0%		.		467	726					.
	All	(819)	(118	16.5%	00	24.1%	(27.9	(22.2%	<0.	59.4%	65.1%	56.4%	00
	PTP		)	(819)	1	(1193)	%)	)	001	(2940)	(1091)	(1849)	1
		27,3			<0								<0
	PTP	%			.		127	127					.
	0-15	(223) 7.0%	33.1%	00	26.2%	(34.1	(21.3%	<0.	22,5%	58.9%	45.6%	00
	%		(26)	(197)	1	(254)	%)	)	001	(710)	(219)	(491)	1
		26,3			<0
	PTP	%	5.6		.		122	172
	15-35	(215)	%	22.3 %	00	22.0%	(24.5	(20.7%	0.0	26,2%	69.8 %
	%		(28)	(187)	1	(294)	%)	)	97	(827)	(347)

Supplementary Table 1: Baseline characteristics of patients suspected of coronary artery disease after a positive coronary computed tomography angiogram or a positive functional test for cardiac ischemia, according to each center.

  ICA, invasive coronary angiography; CAD, coronary artery disease; oCAD, obstructive coronary artery disease; CCTA, computed tomography angiogram; PTP, pre-test probability; GFR, glomerular filtration rate; BMI, body mass index; SD, standard deviation; IQR, interquartile range; SPECT, single photon emission computed tomography; MRI, magnetic resonance imaging; CKD, chronic kidney disease.

		Center 1	Center 2	Center 3	Total	p-
		(n=2184)	(n=1466)	(n=1302)	(n=4952)	valu
						e
	Patient
	characteristics
	Age (years), mean	66.5 (10.1)	65.3 (10.0)	66.2 (10.7)	66.1 (10.3)	<0.00
	(SD)					1
	Male gender, n (%)	1623 (74.3%)	1116 (76.1%)	935 (71.9%)	3674 (74.2%)	0.04
	BMI kg/m 2 , median,	25.9 [23.8,	27.1 [24.5,	26.0 [23.7,	26.2 [24.0,	<0.00
	[IQR]	28.7]	30.8]	28.8]	29.4]	1
	Diabetes mellitus, n	443 (20.3%)	479 (32.7%)	348 (26.7%)	1270 (25.6%) <0.00
	(%)					1
	Dyslipidemia, n (%)	955 (43.7%)	748 (51.0%)	660 (50.7%)	2363 (47.7%) <0.00
						1
	Hypertension, n (%)	965 (44.2%)	659 (45.0%)	570 (43.8%)	2194 (44.3%)	0.82
	Current smoking, n	879 (40.2%)	677 (46.2%)	469 (36.0%)	2025 (40.9%) <0.00
	(%)					1
	CKD, n (%)	13 (0.6%)	28 (1.9%)	9 (0.7%)	50 (1.0%)	<0.00
						1
	Vascular disease, n	112 (5.1%)	128 (8.7%)	72 (5.5%)	312 (6.3%)	<0.00
	(%)					1
	Previous stroke, n	40 (1.8%)	50 (3.4%)	21 (1.6%)	111 (2.2%)	<0.00
	(%)					1
	PTP, median [IQR]	38.0 [18.5,	38.9 [19.2,	37.0 [18.5,	38.0 [18.9,	0.37
		60.7]	61.2]	58.1]	60.4]
	Positive pre-ICA					<0.00
	tests					1
	SPECT	402 (18.4%)	693 (47.3%)	240 (18.4%)	1335 (27.0%)
	Exercise stress	525 (24.0%)	210 (14.3%)	444 (34.1%)	1179 (23.8%)
	test
	S t r e s s	269 (12.3%)	97 (6.6%)	237 (18.2%)	603 (12.2%)
	echocardiography
	Stress MRI	59 (2.7%)	59 (4.0%)	41 (3.1%)	159 (3.2%)
	CCTA	929 (42.5%)	407 (27.8%)	340 (26.1%)	1676 (33.8%)
	Results of ICA					<0.00
						1
	Normal ICA	279 (12.8%)	303 (20.7%)	237 (18.2%)	819 (16.5%)
	Non-oCAD	510 (23.4%)	343 (23.4%)	340 (26.1%)	1193 (24.1%)
	oCAD	1395 (63.9%)	820 (55.9%)	725 (55.7%)	2940 (59.4%)

Conflict of Interest statement: all authors declare to have nothing to declare.

Declaration of Helsinki:

All authors state that this study complies with the Declaration of Helsinki.

Supplementary Figure 1. Results of invasive coronary angiogram according to the positivity of previous non-invasive tests for coronary artery disease, after excluding exercise ECG.

Supplementary Figure 2. Normal coronary angiogram or obstructive coronary artery disease rate ratio in patients with