- Email: [email protected]
The inverse-parametric knapsack problem
EUROPEAN JOURNAL OF OPERATIONAL RESEARCH ELSEVIER
European Journal of Operational Research 83 (1995) 376-393
The inverse-parametric knapsack problem Rainer E. Burkard *, Ulrich Pferschy Institut fiir Mathematik B, TU Graz, Kopernikusgasse 24, A-8010 Graz, Austria
Abstract
This paper deals with parametric knapsack problems where the costs resp. weights are replaced by linear functions depending on a parameter t. The aim is to find the smallest parameter t * such that the optimal solution value of the knapsack problem is equal to a prespecified solution value. For this inverse-parametric problem pseudopolynomial algorithms are developed and search methods making use of the special properties of the parametric value function are constructed. Using computational experiments the behaviour of these algorithms is investigated and the favourable practical performance of different search methods is exhibited. Keywords: Parametric programming; Inverse-parametric; Knapsack problem
I. Introduction
Parametrizations of integer programs have been the subject of active research for many years. Usually the cost coefficients are substituted by (linear) cost functions depending on a parameter t. The corresponding optimal value function assigns to every parameter value t o the optimal solution value of the corresponding problem with fixed costs given by evaluating all cost functions at t o. A different way to look at this setup is the inverse-parametric problem where we do not want to know the complete value function but have a fixed 'target' value given in advance and ask at which parameter t * the optimal value function attains at least this value. In this case the parameter is frequently interpreted as a time variable, and we look for the earliest date a certain threshold value is attained. An instance of an explicitly treated inverse-parametric problem is the quickest path problem (Chen and Chin [4], cf. also Lee and Papadopoulou [9]). In this problem a certain amount of units are given which are to be sent in the fastest possible way along a single path from the source to the sink in a network. The transmission time along a path is given by the sum of transmission times over all arcs of the path. Burkard, Dlaska and Klinz [2] solved the problem of finding the quickest flow using the whole network by computing the earliest time a certain flow value can be sent through the network. The resulting
This research was supported by the Fonds zur F6rderung der wissenschaftlichen Forschung, Project P8971-PHY. * Corresponding author. 0377-2217/95/$09.50 © 1995 Elsevier Science B.V. All rights reserved SSDI 0377-2217(95)00014-3
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
377
optimal value function giving for each time t the maximum possible flow during t time periods is piecewise linear, convex and monotone. Determining the point t * where this function reaches the desired flow value is equivalent to an inverse-parametric problem. Closely related is the quickest disjoint flow problem (Burkard, Dlaska and KeUerer [1]) and the multiterminal version of the quickest flow problem (Hoppe and Tardos [6]) where units of flow are sent between different sources and sinks. (Inverse-)parametric approaches are closely related to fractional programs. Detailed elaborations of this approach are given in Ibaraki [7] and Megiddo [12]. Sometimes approximation algorithms can be applied to the resulting parametric problem. Hashizume et al. [5] give a framework where the accuracy of the parametric solution carries over to the quality of the solution of the original problem. In this paper we investigate inverse-parametric problems where the underlying combinatorial structure is given by a knapsack problem. The related fractional knapsack problem was subject of a study by Ishii et al. [8] and has also been discussed in [5]. An explicit inverse-parametric instance can arise, e.g. from real estate marketing where plots correspond to items and area to weight.
Example. Consider an estate agent who wants to get a costumer interested in a large investment. Due to legal provisions against extensive land holdings, the prospective costumer cannot own more than a certain amount of land. The estate agent has a large number of different plots available, for each of which he can give a good estimation about the future development of its value (because of his experience and knowledge about planned projects and developments). The costumer originally thought of investing into fixed k year government bonds. To persuade the costumer, the agent wants to know in what time (possibly less than k years) and by which selection of plots the same earnings guaranteed by the bond investment can be achieved in real estate. One objective of this research is the development of pseudopolynomial algorithms for solving inverse-parametric knapsack problems using elements of Megiddo's method for fractional programming. Based on properties of the underlying value function, a second approach leads to the construction of search methods which performed convincingly in computational experiments. In the following Section 2 we recall the knapsack problem and its parametrizations. In Section 3 the inverse-parametric knapsack problem is introduced and two solution methods are presented. Extensive computational experiments for both methods are reported in Section 4, followed by final remarks and conclusions.
2. The parametric knapsack problem We briefly describe different parametrizations of the general knapsack problem. The knapsack problem is defined as
(K) max
s.t.
n E CiX i i=l ~ WiX i ~_~b, i=1
xi e {0, 1},
i = 1,...,n,
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995)376-393
378
where c i E ~+ is the value of item i, wi is the positive weight of item i and b the capacity of the knapsack. Recall that one classical method to solve a knapsack problem is the application of dynamic programming (see, e.g. [15]). The knapsack function
z(k,d)=max
cix i ~ w i x i < d , i=1
xi~{O, 1 } , i = l ..... n
i=1
is computed for k = 1. . . . . n and d = 0 . . . . , b using the well known property
z ( k , d) = m a x ( z ( k -
1, d), z ( k , d - wk) + ck).
(1)
The optimal solution of (K) is given by z(n, b) and can be found in pseudopolynomial time O(nb). Parametrization of the right hand side, i.e. replacement of the knapsack capacity b by an arbitrary function preserves the structure of the standard knapsack problem and does not pose a new problem.
Parametric item values If the item costs are parametrized by linear functions we get the following parametric cost knapsack problem: (PKC) max
~.~ (cit + di)xi i=1
s.t.
~ WiXi ~_~b, i=1
X i E { O , 1},
i = 1 . . . . . n,
where c i and d i are the real coefficients of the linear cost function of item i such that cit + d i > 0 for t > 0 (i.e. ci>O, di>O). For (PKC) the optimal value function v(t)=max
~ ( c i t + d i ) x i s.t. i=1
wixi
is piecewise linear and convex. For general parametric programs it is known that the number of breakpoints, i.e. values of t where the slope of v(t) changes, is exponential in the number of variables. Examples of such 'bad' problem instances are given, e.g. by Zadeh [17] and Murty [14]. One may hope that a parametric integer program with only one single constraint shows a more favourable behaviour. But this is not the case as shown by Carstensen. Theorem 1 (Carstensen [3]). For any n, there exists a parametric cost knapsack problem with ½(9n 2 - 7n) variables such that the corresponding value function v(t ) has 2 n - 1 breakpoints in the interval ( - 2 n, 2~). For the dynamic programming scheme the following interpretation of this result can be given: The function z(1, d) is 0 if d < w 1 and a linear function (namely Clt + d 1) otherwise. To compute z(2, d) we have to take the maximum of two linear functions in every step. This generates a piecewise linear function with one breakpoint. In each subsequent step the maximum of two piecewise linear functions is computed, which may in the worst case double the number of breakpoints.
R.E. Burkard, U. Pferschy / European Journal of Operational Research 83 (1995) 376-393
379
Parametric item weights Parametrizing the weights we get the parametric weight knapsack problem which has the form (PKW) max
•
CiX i
i=1
s.t.
~
(wit+vi)xi<_b,
i=1
x i ~ { 0 , 1 },
i = 1 . . . . . n,
where w i and v i are the real coefficients of a linear function representing the weight of item i such that wit+c~>Ofor t>0(i.e.w i>0,v i>0). In this case each knapsack function z ( k , d) is a monotone decreasing staircase function consisting of intervals, where the function is constant. It is shown by the next theorem that the number of these intervals can be exponential.
Theorem 2 (Woeginger [16]). For any n, there exists a parametric weight knapsack problem with all coefficients in the interval (0, 2 n) which has 2 ~ different values in the optimal value function for t ~ [1, 2"]. Proof. Consider the following construction: Set c i = wi = 2 i, v i = 0 , i = 0 . . . . . n - 1, and b = 2 " - 1 . Because the costs are directly proportional to the weights, for any t an optimal knapsack solution is attained by putting as much weight as possible into the knapsack. It is a basic fact that every integer in the interval [1, 2 n - 1] can be represented in a unique way as a sum of distinct values from {1, 2, 4 . . . . . 2"-1}. Hence, for a given integer solution value v ~ [1, 2 n - l ] there exists a unique set I such that ]~i~12 i = v. We consider the p a r a m e t e r t = (2" - 1)/v and get
Y'~ (wit + vi) = Y', 2i(2" - 1 ) / v = 2" - 1 = b. i~l
i~l
As the knapsack is completely filled its value is optimal. According to the above construction 2 n - 1 different solution values can be attained. Setting t = 2 n the optimal value of the knapsack remains 0. [] Based on these parametrizations we will concentrate in the next section on inverse-parametric problems.
3. The inverse-parametric knapsack problem Let us consider a knapsack problem with parametric costs and a value function v(t). For a given target value u * the inverse-parametric knapsack problem asks for the smallest value of the p a r a m e t e r t (if it exists) such that the optimal value of (PKC) is equal to v *" (IPK)
m i n { t l u ( t ) = v*}. In the following two different approaches for solving problem (IPK) are described. At first we develop pseudopolynomial methods where the time complexities depend on the n u m b e r of items n and on the size of the knapsack b.
380
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
The second subsection deals with various search methods based on specific properties of the value function. The latter approach proved to be highly successful in computational experiments. An inverse-parametric problem can be also defined for knapsack problems with parametric weights. In this case the optimal solution identifies the largest parameter such that the corresponding staircase function value is greater or equal to the desired solution value. To avoid the computation of the complete value function of (PKW), methods similar to those presented in Section 3.1 can be used leading to pseudopolynomial algorithms, although we will not elaborate on the details. The use of the same search methods as described in Section 3.2 is not possible because we lack the crucial property of convexity. Other heuristic methods, however, can be developed which are based on properties of the value function like its monotonicity. 3.1. Pseudopolynomial algorithms The algorithms in this subsection are based on the method of Megiddo for solving fractional problems, which was first described in [12] and further elaborated in [13]. It was used in various applications and also for inverse-parametric problems (see, e.g. [2]). It is widely believed, however, to be only of theoretical significance. Its relatively large overhead prevented implementations in many applications. Therefore, no experience exists apparently about the practical behaviour of Megiddo's method. A discussion about our computational experiences is given in Section 4.2. After a short outline of the basic idea of Megiddo's method we describe in more detail how this approach can be used for the inverse-parametric knapsack problem. Megiddo 's method Let A be an algorithm which solves a combinatorial optimization problem P just by additions, subtractions and comparisons of its cost coefficients. Let us assume that the corresponding parametric problem, where the cost coefficients are replaced by linear functions, has a monotone increasing value function. This parametric problem can be solved by an extension of algorithm A: Additions and subtractions are performed with linear functions. A comparison of two linear functions is made by use of their intersection point ~ with a different ordering on either side. To know whether t * will be greater or smaller than ~, the objective value ~7 of problem P with fixed costs given by parameter }" is computed using algorithm A. If ~ > v the result of the comparison is the order of the two linear functions on the left side of t, otherwise the contrary holds. After termination of algorithm A, the optimal solution value is given as a linear function and t * can be computed immediately. The time complexity of this scheme is given by the number of comparisons in A times the running time of A, and is therefore bounded by O(T2(A)). To derive a pseudopolynomial time bound for the (IPK), we use the dynamic programming scheme to adapt Megiddo's method for the inverse-parametric knapsack problem. Megiddo's method requires in general a monotone value function. This can be achieved for the (PKC) by restricting the slopes of the linear costs c i to be positive. Nevertheless, the general case can be treated in a likewise manner by checking at every comparison whether the test parameter ? lies in the decreasing or increasing part of the convex value function. However, to simplify presentation we stick to the monotone case. To solve (IPK) we start with an interval P = [tmi . . . . . . tmax] containing t * and guarantee in every step that the computation performed so far leads to a correct solution of the corresponding (PKC) for the current interval. As mentioned in the interpretation of Theorem 1, the computation of the knapsack functions for (IPK) starts with linear functions for z(1, d). To compute z(2, d) two instances of z(1, *) are compared; the linear cost function of item 2 is added to one of them. This is carried out by determining the
R.E. Burkard, U. Pferschy/ European Journal of OperationalResearch83 (1995) 376-393
381
intersection point i of the two functions. Then a knapsack problem ( / ( ) with fixed costs 6 i := cff + d i is solved by a subroutine. Its optimal value is denoted by t3. If /3 > v the monotonicity of the value function implies that the test parameter i is too large. Therefore the parameter interval P is restricted to [tmin,... , t'], and the precedence which is valid left of i is reported to the main computation. In the opposite case (b < v), we will react conversely setting P = [? . . . . , tmax] and reporting the ordering which is present right of ?. If b = v we have solved (IPK) by chance and stop with t * = i. Sometimes it is admissible to stop even if t3 is within an e-range of v. Obviously, the computation of t3 is useful only if the intersection point ? lies in P. Otherwise i is too large or too small, and the result of the comparison follows without effort. This case will occur very frequently in practice and saves many knapsack computations. Throughout this procedure all z ( k , d) are linear functions. At the end we have attained an interval P necessarily containing i and a linear function v(t) = ott+/3 for z(n, b) which is the optimal solution of (PKC) for all t ~ P. Hence, t * can be computed from at * +/3 = v* The worst case running time of this straight forward algorithm is bounded by O(nZb 2) with O(b) space requirement.
Improvement 1 An improvement of Megiddo's method described already in the original paper [12] consists of postponing the decision on a comparison until it is really necessary to know its outcome. Up to this point the generated intersection points can be 'collected'. When we finally have to use the result of a previous comparison, binary search is applied to the set of collected points yielding a smaller interval on which the comparisons can be uniquely performed. It depends on the specific problem whether this strategy can be applied in a way to improve not only the practical behaviour but also the theoretical complexity bounds. The case of the (IPK) provides a natural way to combine such intersection points. While computing z ( k , * ) only data from z ( k - 1, *) is used. Therefore the b comparisons of stage k do not have to be carried out at once but only after completion of the whole stage ( z ( k , 1)) up to z ( k , b). At the end of each such iteration we have an interval P containing (in the worst case) b new intersection points. To find out the smallest possible interval for t *, binary search will perform at most O(log b) knapsack computations. Hence, the time bound is improved to O(nZb log b) without increasing the asymptotic space requirement. Improvement 2 Obviously, the above idea cannot be extended to stage k + 1 because the computation of z ( k + 1, * ) uses the values of z ( k , *). However, the intersection points of several stages can be collected if the knapsack functions are not required to be linear. This means that computing z ( k , * ) we do not carry out the comparisons required to find the maximum according to (1) but let z ( k , *) consist of two linear pieces. In stage k + 1 the maximum of two functions of the preceding stage is determined, which may in the worst case double the number of linear pieces. Continuing this process up to stage k + s at most O(b 2 s) intersection points are collected. In general, we can partition the n stages of the computation into [n/s] phases containing s stages each. (The last phase may contain fewer stages.) At the end of each phase binary search is applied to the set of generated intersection points to find the smallest interval containing t *. The knapsack functions of the last stage of each phase are updated by reducing the piecewise linear functions to one piece identified by the binary search computations. The total running time is bounded by [b2 s + nb log(b2S)]n/s.
382
R.E. Burkard, U. Pferschy/ European Journal of OperationalResearch 83 (1995)376-393
The first term denotes the time to compute all intersection points of one phase, the second term the effort spent on the knapsack computations during binary search. To minimize this expression, the size of the phases is set to s := log(n log b) which yields Theorem 3. The time complexity of problem ( I P K ) / s bounded by O(n2b log b / l o g ( n log b)).
An obvious drawback of this extended collection of intersection points is the increase of required space to O(nb log b). Improvement 3 Instead of solving the knapsack problems (I() to optimality, approximation algorithms may be used to generate upper and lower bounds t7u and t3~. For the knapsack problem, bounds of high quality can be computed in reasonable time (for details see Section 4.1). After performing the approximation, we test whether the bounds can be used: If /~1 ~> V (and hence t3 > v) or ~u < v ( ~ / ~ < v), the corresponding operations can be performed as described above. Otherwise tT~< v < t3u and the exact solution t3 has to be computed. This method can also be applied for the knapsack problems occurring during the binary search steps. In Section 4.2 we will describe how these different variants behave in practical problems. 3.2. Search methods The nature of the optimal value function of (PKC) (piecewise linear and convex) inspires the development of search methods based on these properties (cf. also [2]). These methods will intuitively and practically turn out to be very efficient. They can solve problems of considerable size amazingly fast. A detailed performance analysis of the different algorithms described in this subsection is given in Section 4.3. As in the previous subsection, the description will be restricted to monotone increasing value functions by dealing with positive slopes c i only. The general case can be treated in a way similar to the approach described below, but with more complicated technical details. The general idea in this subsection is to generate a sequence of smaller and smaller parameter intervals P each containing t *. The number of iterations is defined as the number of evaluations of knapsack problems for fixed parameter values. Obviously, the total effort of such an approach depends almost exclusively on the number of iterations. Other operations performed to compute a new interval are comparatively insignificant. At first we lay out a general framework which can be used for different versions of algorithms. Some of them will be described below. General Search Procedure We denote the optimal value function of (PKC) by v(t). The aim of this procedure is to compute t * such that v(t * ) = v. Let tes denote the linear function implied by the segment of v(t) valid for the parameter s. L is a constant function defined by L ( t ) = v. We keep the steps of the iteration as general as possible and will discuss possible realization afterwards. Initialization step: Step il. Set P = [train,... , tmax] , {tmin, tmax are chosen such that t * ~ P.} Step i2. improve (tmin, tmax); Iteration step: Step 1. Choose t' ~ P.
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
383
v(t)
et,n~
L(t)
train
imroved train
t*
[ L ( t ) O etm~x [
n e....
tmax
L(t) n e..,.
Fig. 1. The value function v(t) of (IPK) and improvements of the bounds.
Step 2. /3 := k n a p s a c k (/~); {solution of a knapsack problem with costs ci? + di.}
Step 3. if /3 > v
then
tma x :---- t';
if /3 < v then tmi n : = t; i f /3 = v t h e n t * := ~'; STOP; Step 4. improve (tmin, tmax); Step 5. if tmax - tmin < ~ t h e n STOP; else goto Step 1; Most of the steps in the above algorithm require some comment or elaboration. Note that Steps 4 and i2 are optional. For a better understanding of the bound improvement procedure improve and the role of ft see Fig. 1. The more practical aspects are discussed in Section 4.3. Initialization: ad Step il. If nothing is known about a possible interval P, a knapsack problem can be solved for t = 0. Comparing the solution value attained at this point with v yields an u p p e r or lower bound of P. For the opposite end of P a bound can be found, e.g. by applying exponential search: Move
384
R.E. Burkard, U. Pferschy/ European Journal of OperationalResearch 83 (1995)376-393
from t = 0 in direction of the missing endpoint with a (large) step width A. If the knapsack problem for t = A yields no bound for P this procedure is iterated and A multiplied in each step by a factor (e.g. greater than two). This will generate bounds for P in a small number of iterations depending on the chosen factor. ad Step i2. The improvement of the bounds is described in Step 4 below. Iteration: ad Step 1. Some of the various possibilities of 'guessing' a parameter ~' are the following: 1 • A-partitioning: ~ : = A t m i n q- (1 - A ) t m a x , A ~ ( 0 , 1) (e.g. A = 7 =, binary search). • Newton type search: ~ := intersection of /~i with L(t). • Sub-approximation: ~':= intersection of tQtm.n and t~tm~x. ad Step 2. The knapsack computation can ~e performed by any available algorithm. The most efficient codes available are Branch and Bound methods (see, e.g. Martello and Toth [10]). Considering the conditions in Step 3, approximation algorithms may be used to get lower and upper bounds for t3. As described in the last subsection, (I() must be solved to optimality only if v lies between these bounds. This will occur usually only in the last few iterations of the algorithm when ~" is already close to t * ad Step 3. As mentioned in the last subsection we are sometimes allowed to stop with a solution satisfying I ~ - v I < e. ad Step 4. Upper bounds for t * are given by the intersection of t~t ~ and get with L(t). By adding up the linear cost functions of the items included in the optimal ~solution o~ (I(), g1 can be computed easily. It is therefore an immediate byproduct of Step 2. However, in this case it is necessary to compute in Step 2 not only the optimal solution value t3, but also the optimal packing scheme. Major difficulties arise if an approximation algorithm is used in Step 2. Constructing tet with an approximate value function leads to considerably less attractive bound improvements. The lower bound of P can be improved by the intersection of L ( t ) with the linear function implied by (train, V(tmin)) and (tmax, V(tmax)) (see Fig. 1). ad Step 5. A better stopping criterion is attained by comparing the slopes of gt and £t . If they are identical, we have identified a linear piece of the value function and get t * by intersecting . max . .~it with L(t). Various combinations of the options mentioned above, especially different choices of ~ and the influence of the bounds, are investigated in the next section. It is intuitively clear that in all cases described above, P will converge towards an interval containing only one linear piece of the value function which guarantees termination. This is obvious for the A-partition, and one can easily see that Newton type search and Sub-approximation eliminate at least one breakpoint from P in every iteration. As the number of breakpoints of v(t) is finite these procedures will necessarily terminate. We conclude that there exists a variety of search methods which are related to each other but may be subject to modification and extension. A comparison of their behaviour is given in Section 4.3. For special applications, additional knowledge about the nature of v(t) may be taken into account for the selection of an algorithm.
4. Numerical experiments In this section we report on implementations of the developed algorithms. We describe their typical performance and interpret their practical behaviour. First the technical details of the experiments are given, then the pseudopolynomial and search algorithms are discussed separately. It should be clear that these two approaches can hardly be compared as they are based on completely different concepts and with different aims. Not surprisingly the latter outperforms the former by a huge margin.
R.E. Burkard, U. Pferschy / European Journal of Operational Research 83 (1995) 376-393
385
The results of the experiments depend strongly on the knapsack problem and cannot be compared to other inverse-parametric problems. Moreover, no other implementations of closely related problems are known to the authors.
4.1. Details of the computational experiments Generation of test data To avoid the treatment of additional subcases we restricted the problems to the first quadrant, i.e. the feasible range for t * is R + and hence v(t) is positive. The number of items n and the knapsack capacity b were fixed in advance for each example. The weights w i were chosen uniformly from [1 . . . . . 500]. To generate random linear functions we chose d i uniformly from (0, MAXVAL) where M A X V A L was set to 500. This upper bound of the distances on the vertical axis makes the trivial case t * = 0 less likely to appear. The slopes c i were chosen such that the angles between the graph of cit + d~ and the horizontal axis are distributed uniformly between 0 and 1 ~-rr. Restricting these angles to smaller ranges naturally increased the value of t * but did not change the characteristic properties of the value function. The choice of v influences the behaviour of the computation. We tested different situations with v fixed as M A X V A L • n • R A T I O and different values of R A T I O between 0.1 and 2. The behaviour of the pseudopolynomial algorithms did not seem to be affected by the magnitude of RATIO, whereas a dependence was observed for the search method which will be explained later. To get less regular and more structured examples which are more difficult to solve, it is common for knapsack problems to generate items with correlated costs and weights. This can not be done for parametric knapsack problems, as the cost vs. weight ratio changes with the parameter. Evaluating knapsack problems Solving the arising knapsack problems is the most time consuming part of the algorithms based on search methods, but only a minor factor of the pseudopolynomial programs. As required by most knapsack routines the items were sorted at first with respect to the cost/weight ratio by a quicksort routine. In the approximation algorithm, lower bounds were determined by a greedy algorithm, which packs all items consecutively into the knapsack if they fit in, without ever removing an inserted item again. In the same procedure the following upper bound given by Martello and Toth [11] was computed: U=
~.~ cj + m a x j=l
c s[
Ws+l J
-
, Ws-1J
where b = b - E j= s-1l Wj and s denotes the 'critical' item (i.e. the first item not fitting into the knapsack when using the greedy algorithm). This bound can be computed with small effort after sorting and gives good estimates for larger problems. To compute optimal knapsack solutions, a code taken from Martello and Toth [10] (there called MT1) was used. Although this is not the fastest known code for large problems, it is still reasonably efficient if the items are presorted, which is necessary anyway for the approximation procedures. The code was translated from Fortran to C and simplified for computing only the optimal knapsack value. For some variants of the search methods, the packing scheme also had to be reported. The behaviour of the knapsack procedure for different test data was consistent with the results reported in [10]. Obviously, the generation of items within a smaller weight range makes the problem more difficult (and takes more time) because a larger number of items is included in an optimal knapsack. For problems with 'large' items the set of candidates for an optimal packing is naturally smaller and an optimal solution can be found faster.
386
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
A similar argument shows that the solution of knapsack problems corresponding to small parameter values takes more time because the cost coefficients are closer to each other, whereas problems with high parameter values have larger gaps between the costs making a faster solution possible. As these properties are not specific for the inverse-parametric knapsack problem we fixed the range of the weights as described above and will not discuss these aspects any further.
Test environment The experiments were performed on a DIGITAL DECstation 5000/240. All algorithms were implemented in the standard C programming language. For all sets of parameters in Tables 1-5, 20 different random examples were generated and the average performance data reported. Tests with the same parameters in different tables originate from the same example. Time is always taken as CPU-time in seconds.
4.2. Experiments with pseudopolynomial algorithms Dimensions of the problems The number of items n was chosen from 1000 up to 10000. The knapsack capacity b was set between 2 000 and 15 000. Larger problems caused in some cases excessive running times and were solved only with search methods. As the magnitude of v had no substantial influence on the behaviour of algorithms, we fixed RATIO = 0.8 in this subsection. Experiments were performed on three types of algorithms as described in Section 3.1. Let us recall their characteristics: 1. Algorithm PPI: Carry out each comparison immediately: Either ? ~ P and the precedence is obvious or a knapsack problem is solved, which determines the outcome of the comparison and yields a new interval P. 2. Algorithm PP2: Collect all intersection points from one stage. Perform binary search to get a smaller interval P on which the comparisons can be determined. 3. Algorithm PP3: Collect the intersection points of maxstage stages. The knapsack functions computed during this collection phase are piecewise linear and may consist of up to 2 to the power of maxstage segments. Perform binary search as above. Reduce each function of the last stage to the segment valid in the new interval. Each algorithm was performed by trying at first approximate solutions for the arising knapsack problems. If the approximation did not settle the outcome of the comparison, the exact knapsack solution was computed. As most of the running time of the approximation algorithm is spent on the sorting of the items, which is also necessary for the exact solution, the amount of 'wasted' time used for trying the approximation in vain was very small. However, the gain of time achieved by using the approximation instead of the exact solution justified this approach. The necessity to compute exact solutions arose mainly during the last few breakpoint calculations. Once a complete knapsack computation had to be performed, the subsequent breakpoints tended to make further complete solutions necessary more frequently, due to the small range of P at this point of the computation. We tried to exploit this behaviour by switching off the approximation module after one or several 'failures' and continue with exact solutions. However, this variant performed hardly differently, sometimes slightly worse, than the program always trying the approximation. Table 1 shows the running time data. It turned out that the amount of time spent on solving or approximating the knapsack problems was very small compared to the running time of the main program.
387
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
Table 1 Performance of PP1 and PP2: For both algorithms total running time, time used for the approximate resp. exact solution of the arising knapsack problems and the number of generated breakpoints are given n b PP1 PP2
1000 1000 1000 2000 2000 2000 4000 4000 4000 6000 6000 6000 6000 8000 8000 8000 10000 10000
2000 4000 6000 2000 4000 6000 2000 4000 6000 6000 8000 10000 12000 8000 12000 16000 10000 15000
Total time (sec.)
Knapsack approx, time
Knapsack solution time
Breakpoints
Total time (sec.)
Knapsack approx, time
Knapsack solution time
Breakpoints
6.6 13.4 21.9 13.3 28.3 43.9 26.5 56.5 88.9 132.1 179.5 224.5 269.4 238.4 366.2 494.2 374.3 575.1
0.22 0.30 0.28 0.56 0.60 0.63 1.17 1.44 1.44 2.49 2.71 2.45 1.97 3.51 3.37 3.15 5.65 5.98
0.04 0.05 0.07 0.08 0.16 0.16 0.19 0.17 0.22 0.40 0.46 0.46 0.42 0.66 0.68 0.62 0.79 0.98
20.6 28.2 28.5 25.1 28.9 30.0 23.9 27.6 28.2 31.1 33.9 31.3 25.5 31.8 30.4 28.7 38.2 41.4
8.2 16,5 27.2 16,5 34.7 54.2 33.0 70.3 108.3 163.9 221.1 279.2 335.2 293.4 446.2 605.8 465.3 706.2
0.23 0.29 0.28 0.59 0.60 0.68 1.25 1.46 1.47 2.71 2.55 2.45 2.11 3.68 3.57 3.29 6.15 6.46
0.04 0.05 0.07 0.09 0.16 0.16 0.19 0.18 0.20 0.41 0.45 0.50 0.48 0.70 0.67 0.68 0.71 0.79
21.7 28.6 28.7 26.1 28.7 31.4 25.0 28.2 28.4 33.1 31.6 30.8 27.3 33.2 31.4 29.8 23.4 43.1
A n i n t e r e s t i n g a s p e c t o f t h e e x p e r i m e n t s is t h e n u m b e r o f g e n e r a t e d i n t e r s e c t i o n p o i n t s ~ a n d w h e r e t h e y w e r e d e t e c t e d by t h e a l g o r i t h m . It t u r n e d o u t t h a t t h e t o t a l n u m b e r of e n c o u n t e r e d b r e a k p o i n t s was surprisingly small (see T a b l e 1). A vast m a j o r i t y o f this n u m b e r was g e n e r a t e d in t h e first few stages w h e n t h e large s t a r t i n g interval was quickly r e d u c e d to a m u c h s m a l l e r one. A f t e r only 3 % o f t h e stages t h e r a n g e of P was always r e d u c e d to t * + 5 % a n d a f t e r less t h a n 10% o f t h e s t a g e s t * was a p p r o x i m a t e d by P w i t h i n a r a n g e o f 1%. H a r d l y any b r e a k p o i n t s o c c u r r e d in t h e s e c o n d half o f t h e stages. This m e a n s t h a t a small n u m b e r o f a r b i t r a r y i t e m s gives a close d e s c r i p t i o n o f t *. T h e r e f o r e , t h e a l g o r i t h m s c o u l d b e s t o p p e d e a r l y a n d w o u l d p r o d u c e a r e a s o n a b l e a p p r o x i m a t i o n to t h e o p t i m a l p a r a m e t e r . D a t a a b o u t this a s p e c t is given in T a b l e 2. N o t e t h a t P P 2 exhibits a b e h a v i o u r similar to PP1. T h e fact t h a t t h e n u m b e r o f e n c o u n t e r e d b r e a k p o i n t s was small in t o t a l a n d very small for a single s t a g e (the l a r g e s t n u m b e r o f d i f f e r e n t b r e a k p o i n t s f o u n d within o n e stage was 12) m a d e t h e t h e o r e t i c a l i m p r o v e m e n t of a l g o r i t h m P P 2 useless, as t h e c o l l e c t i o n o f all b r e a k p o i n t s o f o n e stage usually p r o d u c e d only 5 to 8 points. Thus, the a d d i t i o n a l effort o f collecting t h e s e b r e a k p o i n t s a n d u p d a t i n g t h e d a t a s t r u c t u r e at t h e e n d o f e a c h b r e a k p o i n t e v a l u a t i o n was l a r g e r t h a n t h e gain r e a c h e d by saving a few b r e a k p o i n t c o m p u t a t i o n s . F u r t h e r m o r e , s o m e o f t h e b r e a k p o i n t s in P P 2 w o u l d n o t have b e e n d e t e c t e d at all if P w e r e u p d a t e d i m m e d i a t e l y a f t e r a b r e a k p o i n t was f o u n d (as d o n e by PP1). In m a n y cases t h e s a m e b r e a k p o i n t was c o m p u t e d for several consecutive e n t r i e s o f a stage. This is i n d e e d p l a u s i b l e c o n s i d e r i n g t h a t t h e d e c i s i o n o f t e n has to b e d r a w n b e t w e e n t h e s a m e two p a c k a g e s c h e m e s if t h e k n a p s a c k c a p a c i t y is c h a n g e d by one. T h e s e m u l t i p l e b r e a k p o i n t s w e r e i m m e d i a t e l y e l i m i n a t e d by PP1 b u t h a d to b e s t o r e d in P P 2 for u p d a t i n g at t h e e n d o f t h e stage. T o c h e c k w h e t h e r t h e s e a d v a n t a g e s o f PP1 also h o l d in c o m p a r e n c e with PP3, w e t e s t e d this m o r e
ILE. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
388
Table 2 Stage where the last breakpoint is generated and first stage where P is contained in t * +5% resp. t * + 1%
n
b
PP1
PP2
Last breakpoint
Breakpoint t * + 5%
Breakpoint t * + 1%
Last breakpoint
Breakpoint t * + 5%
Breakpoint t * + 1%
1000 1000 1000 2 000 2 000 2000 4000 4000 4000
2 000 4 000 6000 2000 4 000 6000 2000 4000 6000
666 647 733 1066 1359 1456 3065 2265 1717
29 22 21 55 26 23 47 28 22
111 54 34 156 66 41 105 56 55
666 647 733 1066 1359 1456 3065 2265 1717
28 22 21 50 25 23 47 26 20
111 51 33 147 64 41 110 56 55
6000 6000 6000 6000 8000 8 000 8 000 10 000 10 000
6000 8000 10000 12000 8000 12 000 16 000 10 000 15 000
3198 3455 3747 2218 4521 4 057 3122 4 865 4 680
32 27 19 18 29 24 22 24 23
60 48 35 43 43 42 33 45 49
3198 4510 3747 2218 4521 4 057 3122 4 865 4 680
32 24 18 18 29 24 22 24 23
60 45 35 40 43 47 33 45 49
s o p h i s t i c a t e d a l g o r i t h m with d i f f e r e n t c r i t e r i a to d e t e r m i n e how m a n y s t a g e s s h o u l d b e c o l l e c t e d for e a c h p h a s e . T h r e e p a r a m e t e r s w e r e used: 1. maxsegrnent: M a x i m a l n u m b e r o f d i f f e r e n t line s e g m e n t s o f which o n e k n a p s a c k f u n c t i o n z ( k , d) consists. 2. maxbreak: M a x i m a l n u m b e r o f d i f f e r e n t b r e a k p o i n t s g e n e r a t e d in o n e p h a s e . 3. maxstage: M a x i m a l n u m b e r o f s t a g e s in o n e p h a s e . T h e c o l l e c t i o n p r o c e s s was s t o p p e d ( a n d t h e p h a s e e s t a b l i s h e d t h e r e b y e v a l u a t e d ) as soon as o n e o f t h e s e limits was r e a c h e d . A s i n d i c a t e d a b o v e a l a r g e n u m b e r o f b r e a k p o i n t s a r e f o u n d in t h e first few stages. T h e e x p e r i m e n t s s h o w e d t h a t it is m o r e efficient n o t to collect m a n y b r e a k p o i n t s at t h e b e g i n n i n g b u t to d e c r e a s e t h e interval P very e a r l y by s t a r t i n g with p h a s e s consisting o f only a few stages by c h o o s i n g e i t h e r maxsegment very small (e.g. < 4) o r limiting maxbreak (e.g. < b / 1 0 0 ) . B e c a u s e t h e g e n e r a t e d e x a m p l e s h a d a fairly r e g u l a r s t r u c t u r e , k n a p s a c k f u n c t i o n s consisting o f a n e x c e p t i o n a l l y large n u m b e r o f line s e g m e n t s d i d n o t occur, a n d b o u n d i n g any o f t h e two p a r a m e t e r s y i e l d e d r o u g h l y t h e s a m e results. H o w e v e r , by varying maxbreak a m o r e d e t a i l e d analysis o f t h e b e h a v i o u r o f t h e a l g o r i t h m c o u l d b e p e r f o r m e d . T h e r e f o r e , n o extensive e x a m i n a t i o n c o n c e r n i n g t h e i n f l u e n c e o f maxsegment was p e r formed. By t h e a b o v e s t r a t e g y P will b e very small a f t e r at m o s t t h r e e o r f o u r p h a s e s , which d r a m a t i c a l l y r e d u c e s t h e n u m b e r o f b r e a k p o i n t s f o u n d in t h e c o m p u t a t i o n t h e r e a f t e r . T h e b r e a k p o i n t s c o m p u t e d a f t e r t h e s e initial p h a s e s l e d to k n a p s a c k f u n c t i o n s consisting only o f a small n u m b e r o f line s e g m e n t s ( h a r d l y ever r e a c h i n g maxsegment). D u r i n g this m a i n p e r i o d o f t h e a l g o r i t h m it p r o v e d to b e m o s t efficient to c o m b i n e all stages into only o n e o r two p h a s e s . N o t e t h a t t h e last p h a s e d o e s n o t have to b e e v a l u a t e d at all b e c a u s e t * c a n b e f o u n d d i r e c t l y by i n t e r s e c t i n g z ( k , d) with L(t), which m a k e s a very l a r g e last p h a s e e s p e c i a l l y attractive. Although the different schemes of stage combination produced quite different patterns of behaviour
R.E. Burkard, U. Pferschy / European Journal of Operational Research 83 (1995) 376-393
389
Table 3 Behaviour of PP3 with different stage collection methods: Depending on maxbreak and maxstage the total running time, the time spent on the arising Knapsack subproblems, the number of generated breakpoints and the number of phases is given n
1000
b
4000
max-
max-
break
stage
Total time
40
200 400 1333 400 1333 400 1333
200 400 2 000
6 000
60
300 600 4000
8 000
80
400 800
Knapsack time
Number of breakpoints
Number of phases
22.3 22.3 22.2 23.1 23.0 23.8 24.1
0.29 0.26 0.25 0.26 0.20 0.24 0.22
23.4 22.1 21.2 21.8 16.9 20.2 18.4
7.7 5.6 4.5 4.8 3.1 4.3 3.1
300 600 2 000 600 2 000 600 2 000
72.8 72.9 72.6 73.8 73.6 74.7 74.9
0.68 0.61 0.48 0.59 0.46 0.57 0.39
25.2 22.7 18.5 22.0 17.8 21.9 15.6
9.4 6.3 3.9 5.4 3.0 5.2 2.6
400 800 2666 800 2666 800 2666
195.4 195.3 196.5 197.8 199.2 202.7 202.3
1.84 1.61 1.49 1.57 1.47 1.64 1.53
30.2 26.1 24.5 25.4 23.8 26.3 24.7
12.4 7.9 5.1 7.2 4.3 7.3 4.0
with notable changes in the number and structure of generated breakpoints, they surprisingly had only minor differences in the total running time (see Table 3). As the running-time of PP3 was clearly inferior to PP1, we only give the data of some instances with different stage collecting parameters in Table 3. The performance of PP2 was somewhere between the other two variants but closer towards PP1, which was overall the significantly fastest pseudopolynomial approach. Summarizing, we conclude that the more complicated but asymptotically faster algorithms which save a number of knapsack computations, are clearly inferior to the straight-forward approach. This may not correlate with the behaviour of other applications of Megiddo-based algorithms (of which no practical experiments are known to the authors), but is due to the fact that on one hand the subproblems in our case were efficiently solvable, and on the other hand, only a moderate number of these subproblems had to be solved. Hence, the core of the algorithm (solving the problem with parametric costs within a given parameter range) consumed most of the running time. Therefore, strategies to decrease the number of subproblems did not improve the performance of the algorithm but due to the more complicated data structures and bookkeepng procedures actually increased the total running time. An improvement of the main part of the algorithm does not seem to be at hand when following the Megiddo method.
4.3. Experiments with search methods Dimensions of the problems Much larger problems could be solved in reasonable time by search methods than by pseudopolynomial algorithms. The number of items n was chosen from 5 000 up to 40 000 and the knapsack capacity b set between 5 000 and 80 000.
390
R.E. Burkard, U. Pferschy/ European Journal of Operat•nal Research 83 (1995)376-393
Table 4 Performance of search methods with different values of RATIO n
b
RATIO
20 000
40 000
0.1 1 2
40000
80000
0.1 1 2
Method (BSS)
(BSB)
(NSL)
(NSU)
(SAS)
time iterations time iterations time iterations
4.29 10.2 2.23 5.8 1.85 5.2
2.87 7.2 1.86 4.9 1.83 5.0
1.58 4.5 1.37 4.l 1.38 4.2
1.38 4.1 0.92 3.2 0.91 3.0
3.42 8.4 2.32 6.0 2.28 6.0
time iterations time iterations time iterations
9.71 10.8 4.93 6.1 3.95 5.0
5.93 7.2 3.95 4.9 3.94 4.9
3.55 4.7 2.96 4.0 2.95 4.0
2.97 4.1 1.98 3.1 1.96 3.0
6.92 8.0 4.90 5.9 4.96 6.0
As mentioned above the solution of knapsack problems takes more time for smaller parameter values, whereas problems corresponding to large parameters are usually easier to solve. Therefore, the knapsack computations for problems with smaller values of RATIO tended to take longer because they naturally imply smaller values for t *. Furthermore, also the number of iterations was slightly higher for low RATIOs, especially for the binary search based methods (BSS) and (BSB), which is due to the fact that the search methods work better for strongly ascending parts of the value functions, which are present for higher parameter values (see Table 4 for selected examples). The initial interval P was set to [0. . . . . T] with a large upper bound T fixed according to previous experiences (e.g. T = 2 000). Because of the rapid improvement of P in the first iterations the choice of T hardly effected the performance of the algorithms and is not critical in practice. To choose a new search value ? inside P = [tmin. . . . . tmax], we tested three variants as described in Section 3.2. The A-partitioning was restricted to the classical binary search (A = ½) because other versions were either inferior (for A close to 0 or 1) or performed similarly to the bisection. More important than the choice of ? is how to take advantage of the properties of the search function, especially how to decrease the length of P between the knapsack computations by the bounds given in Step 4 of the general search procedure. P was updated by setting
tmax:=min{~tmmnL(t),
/,maxn L(t)}
and tmin :~
[linear function through
(tmin, U(tmin) ) a n d (tmax, U ( t m a x ) ) ] n L ( t ) .
As indicated in Step 5, any search procedure can be stopped as soon as the slopes of b~tm and Ct~" are equal (or differ only by a small factor). In this case intersecting t~tm, with L(t) yields t~e optimal parameter t * The improvement of P by the above bounding techniques is especially useful in the first few iterations when P is still quite large. Thereby, a considerable number of knapsack evaluations can be saved, although part of this effect is diminished by the necessity to determine ft, which requires the exact knapsack solution at the given point and makes the use of approximation methods for the knapsack computations impossible. After a couple of iterations the interval P was usually rather small which made the use of the bounds less effective. As also the approximation algorithms were not precise enough in this phase of the
R.E. Burkard, U. Pferschy/ EuropeanJournalof OperationalResearch83 (1995)376-393
391
a l g o r i t h m s , t h e exact k n a p s a c k s o l u t i o n h a d to b e c o m p u t e d anyway in m o s t cases. H e n c e , t h e b o u n d c o m p u t a t i o n s d i d n o t c a u s e m u c h a d d i t i o n a l effort a n d w e r e c a r r i e d o u t until t h e e n d o f t h e search. I n c l u d i n g t h e c a l c u l a t i o n o f gt a n d t h e b o u n d s t h e m a i n p r o g r a m t o o k always less t h a n 5 % o f t h e t o t a l r u n n i n g t i m e with t h e r e m a i n d e r b e i n g s p e n t o n t h e k n a p s a c k c o m p u t a t i o n s . H e n c e , we r e p o r t only t h e t o t a l r u n n i n g times. C h a r a c t e r i z e d by t h e choice o f t" six d i f f e r e n t a l g o r i t h m s w e r e finally t e s t e d with t h e following properties and behaviour: (BS): pure binary search. This v a r i a n t r e a c h e d t e r m i n a t i o n only if tmax - tmin < e. O f c o u r s e t h e n u m b e r o f i t e r a t i o n s d e p e n d e d m a i n l y o n t h e size o f e ( r e a s o n a b l e choices s e e m to be, e.g. e = 0.001) a n d was a l m o s t i d e n t i c a l for d i f f e r e n t p r o b l e m s . T h e possibility to use a p p r o x i m a t i o n a l g o r i t h m s for t h e arising k n a p s a c k p r o b l e m s gave s o m e a d v a n t a g e s to (BS), b u t still t h e p e r f o r m a n c e o f t h e other algorithms could not be matched. (BSS): binary search with slope computation. (BSS) is d e r i v e d f r o m (BS) by a d d i n g as s t o p p i n g c r i t e r i o n t h e e - e q u i v a l e n c e o f t h e gt s l o p e s as d e s c r i b e d above. This saves a n u m b e r o f i t e r a t i o n s at t h e e n d o f t h e a l g o r i t h m a n d l e a d s to a s h o r t e r r u n n i n g time t h a n (BS), a l t h o u g h the c o m p l e t e s o l u t i o n o f t h e k n a p s a c k p r o b l e m s is n e c e s s a r y in every iteration. (BSB): binary search with slope computation and bounding. By a d d i n g to (BSS) t h e b o u n d i n g p r o c e d u r e s d e s c r i b e d above, a very efficient a l g o r i t h m s is a t t a i n e d . T h e t o t a l n u m b e r o f i t e r a t i o n s was always less t h a n 10. N o t e t h a t t°t is only a v a i l a b l e for a b o u n d o f P which was g e n e r a t e d as s e a r c h v a l u e f in a p r e v i o u s i t e r a t i o n , b e c a u s e k n a p s a c k solutions a r e not k n o w n for interval b o u n d s /impr d e r i v e d by b o u n d i m p r o v e m e n t s . H e n c e , in p r a c t i c e we can use only r e l a x e d v e r s i o n s o f t h e full i m p r o v e m e n t by using, e.g. for timpr t h e l i n e a r f u n c t i o n m a t c h i n g v ( t ) at timpr with t h e slope t a k e n f r o m gton w h e r e torig is t h e o r i g i n a l b o u n d . (NSL): Newton type search starting from O. T h e N e w t o n m e t h o d c a n n o t m a k e use of any b o u n d s , b e c a u s e •
g
Table 5 Total running time and number of iterations for different search methods (RATIO = 0.8) n
b 5 000
5 000
5 000
10000
10 000
15 000
10000
30000
20 000
20 000
20000
40 000
20 000
60 000
40 000
40 000
40000
60 000
40 000
80 000
Method time iterations time iterations time iterations time iterations time iterations time iterations time iterations time iterations time iterations time iterations
(BS)
(BSS)
(BSB)
(NSL)
(NSU)
(SAS)
0.83 11.8 0.84 11.9 1.90 12.0 1.91 12.0 4.26 12.0 4.28 12.0 4.25 12.0 9.27 12.0 9.31 12.1 9.32 12.2
0.41 5.6 0.53 6.7 1.08 6.2 1.24 6.8 2.12 5.6 2.40 6.3 2.59 6.5 4.95 6.1 4.98 6.0 4.96 6.1
0.37 4.9 0.38 4.9 0.82 5.0 0.83 5.0 1.82 4.9 1.84 4.9 1.83 5.0 3.91 4.9 3.97 5.0 3.92 5.0
0.28 4.1 0.29 4.2 0.63 4.2 0.64 4.2 1.37 4.1 1.39 4.2 1.39 4.3 2.94 4.0 2.95 4.2 2.99 4.3
0.19 3.2 0.20 3.2 0.48 3.3 0.42 3.1 0.92 3.2 0.94 3.2 0.92 3.0 1.98 3.0 1.98 3.1 1.99 3.0
0.48 6.1 0.45 6.0 1.03 6.0 1.05 5.9 2.33 6.0 2.37 6.2 2.31 6.1 4.79 6.1 4.59 5.9 4.63 6.0
392
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
only the slope at the last ~' is used and no interval improvement makes sense. Still it performed very well and competed successfully with (BSB). (NSU): Newton type search starting f r o m T. In most cases the Newton method took one or two iterations less if started on the right hand side of t * rather than at 0. This approach was the best method in all experiments. There was a tendency observed that (NSU) took less iterations than (NSL) if the slopes c i were chosen between 0 and infinity, which causes a very steep value function, whereas (NSL) was superior if the slopes were restricted to a smaller interval, e.g. (0 . . . . . tan('rr/16)). (SAS): Sub-approximation search. The performance of this approach was varying more than the other methods. Although it sometimes produced solutions surprisingly fast, in some cases it was the worst of all variants. However, also for this method the number of iterations was very stable. As can be seen on Table 5 the algorithms performed regularly with a modest sublinear increase of the running time with respect to the problem size. We summarize that the inverse-parametric knapsack problem can be solved very efficiently in a small number of iterations by search methods exploiting the properties of the optimal value function. In contrast to the pseudopolynomial methods, in this case the more sophisticated (though still not too complicate) algorithms produced the best results. Newton methods usually performed best, followed by (BSB). This might be explained by the regular structure of random test examples. Their smooth value function made the approximation of the function by ft as employed by the Newton method highly successful. We suspect that in examples with more irregular data, the advantages of (BSB) will gain more importance. However, sets of linear functions generally lead to more regular shapes of value functions than constant cost coefficients. Degenerate examples usually lead to difficult problems only for a special range of v and are easy to solve for all other target values.
5. Conclusions In this paper we treated the inverse-parametric problem where the costs (resp. the weights) of a knapsack problem are replaced by linear functions in a parameter t and the smallest possible parameter t * has to be found such that the optimal solution of the knapsack problem taken at t * is equal to a given input value. Based on a general method of Megiddo, we constructed a pseudopolynomial algorithm for the parametric cost case with running time O(n2b log b / l o g ( n log b)). Exploiting the specific properties of the parametric value function, we developed various geometrically motivated search methods. Computational studies showed that the more complicated and theoretically improved pseudopolynomial methods performed worse than a straight-forward approach which is due to the fact that the solution of the arising subproblems (knapsack problems) can either be avoided in many cases or computed efficiently. Comparative analysis of the search methods showed that even large problems could be solved very fast, especially by Newton-based approaches.
References [1] Burkard, R.E., Dlaska, K., and Kellerer, H., "The quickest disjoint flow problem", Report 189, Institute of Mathematics, Technical UniversityGraz, Austria, 1991.
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
393
[2] Burkard, R.E., Dlaska, K., and Klinz, B., "The quickest flow problem", Zeitschrift fiir Operations Research 37 (1993) 31-58. [3] Carstensen, P.J., "Complexity of some parametric integer and network programming problems", Mathematical Programming 26 (1983) 64-75. [4] Chen, Y.L., and Chin, Y.H., "The quickest path problem", Computers & Operations Research 17 (1990) 153-161. [5] Hashizume, S., Fukushima, M., Katoh, N., and Ibaraki, T., "Approximation algorithms for combinatorial fractional programming problems", Mathematical Programming 37 (1987) 255-267. [6] Hoppe, B., and Tardos, E., "The quickest transshipment problem", to appear in: Proc. of the Sixth AnnualACM-SIAM Symp. on Discrete Algorithms, 1995. [7] Ibaraki, T., "Parametric approaches to fractional programs", Mathematical Programming 26 (1983) 345-362. [8] Ishii, H., Ibaraki, T., and Mine, H., "Fractional knapsack problems", Mathematical Programming 13 (1977) 255-271. [9] Lee, D.T., and Papadopoulou, E., "The all-pairs quickest path problem", Information Processing Letters 45 (1993) 261-267. [10] Martello, S., and Toth, P., Knapsack Problems, Wiley, New York, 1990. [11] Martello, S., and Toth, P., "An upper bound for the zero-one knapsack problem and a branch and bound algorithm", European Journal of Operational Research 1 (1977) 169-175. [12] Megiddo, N., "Combinatorial optimization with rational objective functions", Mathematics of Operations Research 4 (1979) 414-424. [13] Megiddo, N., "Applying parallel computation algorithms in the design of serial algorithms", Journal of the ACM 30 (1983) 852-865. [14] Murty, K.G., "Computational complexity of parametric linear programming", Mathematical Programming 19 (1980) 213-219. [15] Nemhauser, G.L., and Wolsey, L.A., Integer and Combinatorial Optimization, Wiley, New York, 1988. [16] Woeginger, G.J., Private communication, 1994. [17] Zadeh, N., "A bad network problem for the simplex method and other minimum cost flow algorithms", Mathematical Programming 5 (1973) 255-266.
European Journal of Operational Research 83 (1995) 376-393
The inverse-parametric knapsack problem Rainer E. Burkard *, Ulrich Pferschy Institut fiir Mathematik B, TU Graz, Kopernikusgasse 24, A-8010 Graz, Austria
Abstract
This paper deals with parametric knapsack problems where the costs resp. weights are replaced by linear functions depending on a parameter t. The aim is to find the smallest parameter t * such that the optimal solution value of the knapsack problem is equal to a prespecified solution value. For this inverse-parametric problem pseudopolynomial algorithms are developed and search methods making use of the special properties of the parametric value function are constructed. Using computational experiments the behaviour of these algorithms is investigated and the favourable practical performance of different search methods is exhibited. Keywords: Parametric programming; Inverse-parametric; Knapsack problem
I. Introduction
Parametrizations of integer programs have been the subject of active research for many years. Usually the cost coefficients are substituted by (linear) cost functions depending on a parameter t. The corresponding optimal value function assigns to every parameter value t o the optimal solution value of the corresponding problem with fixed costs given by evaluating all cost functions at t o. A different way to look at this setup is the inverse-parametric problem where we do not want to know the complete value function but have a fixed 'target' value given in advance and ask at which parameter t * the optimal value function attains at least this value. In this case the parameter is frequently interpreted as a time variable, and we look for the earliest date a certain threshold value is attained. An instance of an explicitly treated inverse-parametric problem is the quickest path problem (Chen and Chin [4], cf. also Lee and Papadopoulou [9]). In this problem a certain amount of units are given which are to be sent in the fastest possible way along a single path from the source to the sink in a network. The transmission time along a path is given by the sum of transmission times over all arcs of the path. Burkard, Dlaska and Klinz [2] solved the problem of finding the quickest flow using the whole network by computing the earliest time a certain flow value can be sent through the network. The resulting
This research was supported by the Fonds zur F6rderung der wissenschaftlichen Forschung, Project P8971-PHY. * Corresponding author. 0377-2217/95/$09.50 © 1995 Elsevier Science B.V. All rights reserved SSDI 0377-2217(95)00014-3
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
377
optimal value function giving for each time t the maximum possible flow during t time periods is piecewise linear, convex and monotone. Determining the point t * where this function reaches the desired flow value is equivalent to an inverse-parametric problem. Closely related is the quickest disjoint flow problem (Burkard, Dlaska and KeUerer [1]) and the multiterminal version of the quickest flow problem (Hoppe and Tardos [6]) where units of flow are sent between different sources and sinks. (Inverse-)parametric approaches are closely related to fractional programs. Detailed elaborations of this approach are given in Ibaraki [7] and Megiddo [12]. Sometimes approximation algorithms can be applied to the resulting parametric problem. Hashizume et al. [5] give a framework where the accuracy of the parametric solution carries over to the quality of the solution of the original problem. In this paper we investigate inverse-parametric problems where the underlying combinatorial structure is given by a knapsack problem. The related fractional knapsack problem was subject of a study by Ishii et al. [8] and has also been discussed in [5]. An explicit inverse-parametric instance can arise, e.g. from real estate marketing where plots correspond to items and area to weight.
Example. Consider an estate agent who wants to get a costumer interested in a large investment. Due to legal provisions against extensive land holdings, the prospective costumer cannot own more than a certain amount of land. The estate agent has a large number of different plots available, for each of which he can give a good estimation about the future development of its value (because of his experience and knowledge about planned projects and developments). The costumer originally thought of investing into fixed k year government bonds. To persuade the costumer, the agent wants to know in what time (possibly less than k years) and by which selection of plots the same earnings guaranteed by the bond investment can be achieved in real estate. One objective of this research is the development of pseudopolynomial algorithms for solving inverse-parametric knapsack problems using elements of Megiddo's method for fractional programming. Based on properties of the underlying value function, a second approach leads to the construction of search methods which performed convincingly in computational experiments. In the following Section 2 we recall the knapsack problem and its parametrizations. In Section 3 the inverse-parametric knapsack problem is introduced and two solution methods are presented. Extensive computational experiments for both methods are reported in Section 4, followed by final remarks and conclusions.
2. The parametric knapsack problem We briefly describe different parametrizations of the general knapsack problem. The knapsack problem is defined as
(K) max
s.t.
n E CiX i i=l ~ WiX i ~_~b, i=1
xi e {0, 1},
i = 1,...,n,
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995)376-393
378
where c i E ~+ is the value of item i, wi is the positive weight of item i and b the capacity of the knapsack. Recall that one classical method to solve a knapsack problem is the application of dynamic programming (see, e.g. [15]). The knapsack function
z(k,d)=max
cix i ~ w i x i < d , i=1
xi~{O, 1 } , i = l ..... n
i=1
is computed for k = 1. . . . . n and d = 0 . . . . , b using the well known property
z ( k , d) = m a x ( z ( k -
1, d), z ( k , d - wk) + ck).
(1)
The optimal solution of (K) is given by z(n, b) and can be found in pseudopolynomial time O(nb). Parametrization of the right hand side, i.e. replacement of the knapsack capacity b by an arbitrary function preserves the structure of the standard knapsack problem and does not pose a new problem.
Parametric item values If the item costs are parametrized by linear functions we get the following parametric cost knapsack problem: (PKC) max
~.~ (cit + di)xi i=1
s.t.
~ WiXi ~_~b, i=1
X i E { O , 1},
i = 1 . . . . . n,
where c i and d i are the real coefficients of the linear cost function of item i such that cit + d i > 0 for t > 0 (i.e. ci>O, di>O). For (PKC) the optimal value function v(t)=max
~ ( c i t + d i ) x i s.t. i=1
wixi
is piecewise linear and convex. For general parametric programs it is known that the number of breakpoints, i.e. values of t where the slope of v(t) changes, is exponential in the number of variables. Examples of such 'bad' problem instances are given, e.g. by Zadeh [17] and Murty [14]. One may hope that a parametric integer program with only one single constraint shows a more favourable behaviour. But this is not the case as shown by Carstensen. Theorem 1 (Carstensen [3]). For any n, there exists a parametric cost knapsack problem with ½(9n 2 - 7n) variables such that the corresponding value function v(t ) has 2 n - 1 breakpoints in the interval ( - 2 n, 2~). For the dynamic programming scheme the following interpretation of this result can be given: The function z(1, d) is 0 if d < w 1 and a linear function (namely Clt + d 1) otherwise. To compute z(2, d) we have to take the maximum of two linear functions in every step. This generates a piecewise linear function with one breakpoint. In each subsequent step the maximum of two piecewise linear functions is computed, which may in the worst case double the number of breakpoints.
R.E. Burkard, U. Pferschy / European Journal of Operational Research 83 (1995) 376-393
379
Parametric item weights Parametrizing the weights we get the parametric weight knapsack problem which has the form (PKW) max
•
CiX i
i=1
s.t.
~
(wit+vi)xi<_b,
i=1
x i ~ { 0 , 1 },
i = 1 . . . . . n,
where w i and v i are the real coefficients of a linear function representing the weight of item i such that wit+c~>Ofor t>0(i.e.w i>0,v i>0). In this case each knapsack function z ( k , d) is a monotone decreasing staircase function consisting of intervals, where the function is constant. It is shown by the next theorem that the number of these intervals can be exponential.
Theorem 2 (Woeginger [16]). For any n, there exists a parametric weight knapsack problem with all coefficients in the interval (0, 2 n) which has 2 ~ different values in the optimal value function for t ~ [1, 2"]. Proof. Consider the following construction: Set c i = wi = 2 i, v i = 0 , i = 0 . . . . . n - 1, and b = 2 " - 1 . Because the costs are directly proportional to the weights, for any t an optimal knapsack solution is attained by putting as much weight as possible into the knapsack. It is a basic fact that every integer in the interval [1, 2 n - 1] can be represented in a unique way as a sum of distinct values from {1, 2, 4 . . . . . 2"-1}. Hence, for a given integer solution value v ~ [1, 2 n - l ] there exists a unique set I such that ]~i~12 i = v. We consider the p a r a m e t e r t = (2" - 1)/v and get
Y'~ (wit + vi) = Y', 2i(2" - 1 ) / v = 2" - 1 = b. i~l
i~l
As the knapsack is completely filled its value is optimal. According to the above construction 2 n - 1 different solution values can be attained. Setting t = 2 n the optimal value of the knapsack remains 0. [] Based on these parametrizations we will concentrate in the next section on inverse-parametric problems.
3. The inverse-parametric knapsack problem Let us consider a knapsack problem with parametric costs and a value function v(t). For a given target value u * the inverse-parametric knapsack problem asks for the smallest value of the p a r a m e t e r t (if it exists) such that the optimal value of (PKC) is equal to v *" (IPK)
m i n { t l u ( t ) = v*}. In the following two different approaches for solving problem (IPK) are described. At first we develop pseudopolynomial methods where the time complexities depend on the n u m b e r of items n and on the size of the knapsack b.
380
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
The second subsection deals with various search methods based on specific properties of the value function. The latter approach proved to be highly successful in computational experiments. An inverse-parametric problem can be also defined for knapsack problems with parametric weights. In this case the optimal solution identifies the largest parameter such that the corresponding staircase function value is greater or equal to the desired solution value. To avoid the computation of the complete value function of (PKW), methods similar to those presented in Section 3.1 can be used leading to pseudopolynomial algorithms, although we will not elaborate on the details. The use of the same search methods as described in Section 3.2 is not possible because we lack the crucial property of convexity. Other heuristic methods, however, can be developed which are based on properties of the value function like its monotonicity. 3.1. Pseudopolynomial algorithms The algorithms in this subsection are based on the method of Megiddo for solving fractional problems, which was first described in [12] and further elaborated in [13]. It was used in various applications and also for inverse-parametric problems (see, e.g. [2]). It is widely believed, however, to be only of theoretical significance. Its relatively large overhead prevented implementations in many applications. Therefore, no experience exists apparently about the practical behaviour of Megiddo's method. A discussion about our computational experiences is given in Section 4.2. After a short outline of the basic idea of Megiddo's method we describe in more detail how this approach can be used for the inverse-parametric knapsack problem. Megiddo 's method Let A be an algorithm which solves a combinatorial optimization problem P just by additions, subtractions and comparisons of its cost coefficients. Let us assume that the corresponding parametric problem, where the cost coefficients are replaced by linear functions, has a monotone increasing value function. This parametric problem can be solved by an extension of algorithm A: Additions and subtractions are performed with linear functions. A comparison of two linear functions is made by use of their intersection point ~ with a different ordering on either side. To know whether t * will be greater or smaller than ~, the objective value ~7 of problem P with fixed costs given by parameter }" is computed using algorithm A. If ~ > v the result of the comparison is the order of the two linear functions on the left side of t, otherwise the contrary holds. After termination of algorithm A, the optimal solution value is given as a linear function and t * can be computed immediately. The time complexity of this scheme is given by the number of comparisons in A times the running time of A, and is therefore bounded by O(T2(A)). To derive a pseudopolynomial time bound for the (IPK), we use the dynamic programming scheme to adapt Megiddo's method for the inverse-parametric knapsack problem. Megiddo's method requires in general a monotone value function. This can be achieved for the (PKC) by restricting the slopes of the linear costs c i to be positive. Nevertheless, the general case can be treated in a likewise manner by checking at every comparison whether the test parameter ? lies in the decreasing or increasing part of the convex value function. However, to simplify presentation we stick to the monotone case. To solve (IPK) we start with an interval P = [tmi . . . . . . tmax] containing t * and guarantee in every step that the computation performed so far leads to a correct solution of the corresponding (PKC) for the current interval. As mentioned in the interpretation of Theorem 1, the computation of the knapsack functions for (IPK) starts with linear functions for z(1, d). To compute z(2, d) two instances of z(1, *) are compared; the linear cost function of item 2 is added to one of them. This is carried out by determining the
R.E. Burkard, U. Pferschy/ European Journal of OperationalResearch83 (1995) 376-393
381
intersection point i of the two functions. Then a knapsack problem ( / ( ) with fixed costs 6 i := cff + d i is solved by a subroutine. Its optimal value is denoted by t3. If /3 > v the monotonicity of the value function implies that the test parameter i is too large. Therefore the parameter interval P is restricted to [tmin,... , t'], and the precedence which is valid left of i is reported to the main computation. In the opposite case (b < v), we will react conversely setting P = [? . . . . , tmax] and reporting the ordering which is present right of ?. If b = v we have solved (IPK) by chance and stop with t * = i. Sometimes it is admissible to stop even if t3 is within an e-range of v. Obviously, the computation of t3 is useful only if the intersection point ? lies in P. Otherwise i is too large or too small, and the result of the comparison follows without effort. This case will occur very frequently in practice and saves many knapsack computations. Throughout this procedure all z ( k , d) are linear functions. At the end we have attained an interval P necessarily containing i and a linear function v(t) = ott+/3 for z(n, b) which is the optimal solution of (PKC) for all t ~ P. Hence, t * can be computed from at * +/3 = v* The worst case running time of this straight forward algorithm is bounded by O(nZb 2) with O(b) space requirement.
Improvement 1 An improvement of Megiddo's method described already in the original paper [12] consists of postponing the decision on a comparison until it is really necessary to know its outcome. Up to this point the generated intersection points can be 'collected'. When we finally have to use the result of a previous comparison, binary search is applied to the set of collected points yielding a smaller interval on which the comparisons can be uniquely performed. It depends on the specific problem whether this strategy can be applied in a way to improve not only the practical behaviour but also the theoretical complexity bounds. The case of the (IPK) provides a natural way to combine such intersection points. While computing z ( k , * ) only data from z ( k - 1, *) is used. Therefore the b comparisons of stage k do not have to be carried out at once but only after completion of the whole stage ( z ( k , 1)) up to z ( k , b). At the end of each such iteration we have an interval P containing (in the worst case) b new intersection points. To find out the smallest possible interval for t *, binary search will perform at most O(log b) knapsack computations. Hence, the time bound is improved to O(nZb log b) without increasing the asymptotic space requirement. Improvement 2 Obviously, the above idea cannot be extended to stage k + 1 because the computation of z ( k + 1, * ) uses the values of z ( k , *). However, the intersection points of several stages can be collected if the knapsack functions are not required to be linear. This means that computing z ( k , * ) we do not carry out the comparisons required to find the maximum according to (1) but let z ( k , *) consist of two linear pieces. In stage k + 1 the maximum of two functions of the preceding stage is determined, which may in the worst case double the number of linear pieces. Continuing this process up to stage k + s at most O(b 2 s) intersection points are collected. In general, we can partition the n stages of the computation into [n/s] phases containing s stages each. (The last phase may contain fewer stages.) At the end of each phase binary search is applied to the set of generated intersection points to find the smallest interval containing t *. The knapsack functions of the last stage of each phase are updated by reducing the piecewise linear functions to one piece identified by the binary search computations. The total running time is bounded by [b2 s + nb log(b2S)]n/s.
382
R.E. Burkard, U. Pferschy/ European Journal of OperationalResearch 83 (1995)376-393
The first term denotes the time to compute all intersection points of one phase, the second term the effort spent on the knapsack computations during binary search. To minimize this expression, the size of the phases is set to s := log(n log b) which yields Theorem 3. The time complexity of problem ( I P K ) / s bounded by O(n2b log b / l o g ( n log b)).
An obvious drawback of this extended collection of intersection points is the increase of required space to O(nb log b). Improvement 3 Instead of solving the knapsack problems (I() to optimality, approximation algorithms may be used to generate upper and lower bounds t7u and t3~. For the knapsack problem, bounds of high quality can be computed in reasonable time (for details see Section 4.1). After performing the approximation, we test whether the bounds can be used: If /~1 ~> V (and hence t3 > v) or ~u < v ( ~ / ~ < v), the corresponding operations can be performed as described above. Otherwise tT~< v < t3u and the exact solution t3 has to be computed. This method can also be applied for the knapsack problems occurring during the binary search steps. In Section 4.2 we will describe how these different variants behave in practical problems. 3.2. Search methods The nature of the optimal value function of (PKC) (piecewise linear and convex) inspires the development of search methods based on these properties (cf. also [2]). These methods will intuitively and practically turn out to be very efficient. They can solve problems of considerable size amazingly fast. A detailed performance analysis of the different algorithms described in this subsection is given in Section 4.3. As in the previous subsection, the description will be restricted to monotone increasing value functions by dealing with positive slopes c i only. The general case can be treated in a way similar to the approach described below, but with more complicated technical details. The general idea in this subsection is to generate a sequence of smaller and smaller parameter intervals P each containing t *. The number of iterations is defined as the number of evaluations of knapsack problems for fixed parameter values. Obviously, the total effort of such an approach depends almost exclusively on the number of iterations. Other operations performed to compute a new interval are comparatively insignificant. At first we lay out a general framework which can be used for different versions of algorithms. Some of them will be described below. General Search Procedure We denote the optimal value function of (PKC) by v(t). The aim of this procedure is to compute t * such that v(t * ) = v. Let tes denote the linear function implied by the segment of v(t) valid for the parameter s. L is a constant function defined by L ( t ) = v. We keep the steps of the iteration as general as possible and will discuss possible realization afterwards. Initialization step: Step il. Set P = [train,... , tmax] , {tmin, tmax are chosen such that t * ~ P.} Step i2. improve (tmin, tmax); Iteration step: Step 1. Choose t' ~ P.
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
383
v(t)
et,n~
L(t)
train
imroved train
t*
[ L ( t ) O etm~x [
n e....
tmax
L(t) n e..,.
Fig. 1. The value function v(t) of (IPK) and improvements of the bounds.
Step 2. /3 := k n a p s a c k (/~); {solution of a knapsack problem with costs ci? + di.}
Step 3. if /3 > v
then
tma x :---- t';
if /3 < v then tmi n : = t; i f /3 = v t h e n t * := ~'; STOP; Step 4. improve (tmin, tmax); Step 5. if tmax - tmin < ~ t h e n STOP; else goto Step 1; Most of the steps in the above algorithm require some comment or elaboration. Note that Steps 4 and i2 are optional. For a better understanding of the bound improvement procedure improve and the role of ft see Fig. 1. The more practical aspects are discussed in Section 4.3. Initialization: ad Step il. If nothing is known about a possible interval P, a knapsack problem can be solved for t = 0. Comparing the solution value attained at this point with v yields an u p p e r or lower bound of P. For the opposite end of P a bound can be found, e.g. by applying exponential search: Move
384
R.E. Burkard, U. Pferschy/ European Journal of OperationalResearch 83 (1995)376-393
from t = 0 in direction of the missing endpoint with a (large) step width A. If the knapsack problem for t = A yields no bound for P this procedure is iterated and A multiplied in each step by a factor (e.g. greater than two). This will generate bounds for P in a small number of iterations depending on the chosen factor. ad Step i2. The improvement of the bounds is described in Step 4 below. Iteration: ad Step 1. Some of the various possibilities of 'guessing' a parameter ~' are the following: 1 • A-partitioning: ~ : = A t m i n q- (1 - A ) t m a x , A ~ ( 0 , 1) (e.g. A = 7 =, binary search). • Newton type search: ~ := intersection of /~i with L(t). • Sub-approximation: ~':= intersection of tQtm.n and t~tm~x. ad Step 2. The knapsack computation can ~e performed by any available algorithm. The most efficient codes available are Branch and Bound methods (see, e.g. Martello and Toth [10]). Considering the conditions in Step 3, approximation algorithms may be used to get lower and upper bounds for t3. As described in the last subsection, (I() must be solved to optimality only if v lies between these bounds. This will occur usually only in the last few iterations of the algorithm when ~" is already close to t * ad Step 3. As mentioned in the last subsection we are sometimes allowed to stop with a solution satisfying I ~ - v I < e. ad Step 4. Upper bounds for t * are given by the intersection of t~t ~ and get with L(t). By adding up the linear cost functions of the items included in the optimal ~solution o~ (I(), g1 can be computed easily. It is therefore an immediate byproduct of Step 2. However, in this case it is necessary to compute in Step 2 not only the optimal solution value t3, but also the optimal packing scheme. Major difficulties arise if an approximation algorithm is used in Step 2. Constructing tet with an approximate value function leads to considerably less attractive bound improvements. The lower bound of P can be improved by the intersection of L ( t ) with the linear function implied by (train, V(tmin)) and (tmax, V(tmax)) (see Fig. 1). ad Step 5. A better stopping criterion is attained by comparing the slopes of gt and £t . If they are identical, we have identified a linear piece of the value function and get t * by intersecting . max . .~it with L(t). Various combinations of the options mentioned above, especially different choices of ~ and the influence of the bounds, are investigated in the next section. It is intuitively clear that in all cases described above, P will converge towards an interval containing only one linear piece of the value function which guarantees termination. This is obvious for the A-partition, and one can easily see that Newton type search and Sub-approximation eliminate at least one breakpoint from P in every iteration. As the number of breakpoints of v(t) is finite these procedures will necessarily terminate. We conclude that there exists a variety of search methods which are related to each other but may be subject to modification and extension. A comparison of their behaviour is given in Section 4.3. For special applications, additional knowledge about the nature of v(t) may be taken into account for the selection of an algorithm.
4. Numerical experiments In this section we report on implementations of the developed algorithms. We describe their typical performance and interpret their practical behaviour. First the technical details of the experiments are given, then the pseudopolynomial and search algorithms are discussed separately. It should be clear that these two approaches can hardly be compared as they are based on completely different concepts and with different aims. Not surprisingly the latter outperforms the former by a huge margin.
R.E. Burkard, U. Pferschy / European Journal of Operational Research 83 (1995) 376-393
385
The results of the experiments depend strongly on the knapsack problem and cannot be compared to other inverse-parametric problems. Moreover, no other implementations of closely related problems are known to the authors.
4.1. Details of the computational experiments Generation of test data To avoid the treatment of additional subcases we restricted the problems to the first quadrant, i.e. the feasible range for t * is R + and hence v(t) is positive. The number of items n and the knapsack capacity b were fixed in advance for each example. The weights w i were chosen uniformly from [1 . . . . . 500]. To generate random linear functions we chose d i uniformly from (0, MAXVAL) where M A X V A L was set to 500. This upper bound of the distances on the vertical axis makes the trivial case t * = 0 less likely to appear. The slopes c i were chosen such that the angles between the graph of cit + d~ and the horizontal axis are distributed uniformly between 0 and 1 ~-rr. Restricting these angles to smaller ranges naturally increased the value of t * but did not change the characteristic properties of the value function. The choice of v influences the behaviour of the computation. We tested different situations with v fixed as M A X V A L • n • R A T I O and different values of R A T I O between 0.1 and 2. The behaviour of the pseudopolynomial algorithms did not seem to be affected by the magnitude of RATIO, whereas a dependence was observed for the search method which will be explained later. To get less regular and more structured examples which are more difficult to solve, it is common for knapsack problems to generate items with correlated costs and weights. This can not be done for parametric knapsack problems, as the cost vs. weight ratio changes with the parameter. Evaluating knapsack problems Solving the arising knapsack problems is the most time consuming part of the algorithms based on search methods, but only a minor factor of the pseudopolynomial programs. As required by most knapsack routines the items were sorted at first with respect to the cost/weight ratio by a quicksort routine. In the approximation algorithm, lower bounds were determined by a greedy algorithm, which packs all items consecutively into the knapsack if they fit in, without ever removing an inserted item again. In the same procedure the following upper bound given by Martello and Toth [11] was computed: U=
~.~ cj + m a x j=l
c s[
Ws+l J
-
, Ws-1J
where b = b - E j= s-1l Wj and s denotes the 'critical' item (i.e. the first item not fitting into the knapsack when using the greedy algorithm). This bound can be computed with small effort after sorting and gives good estimates for larger problems. To compute optimal knapsack solutions, a code taken from Martello and Toth [10] (there called MT1) was used. Although this is not the fastest known code for large problems, it is still reasonably efficient if the items are presorted, which is necessary anyway for the approximation procedures. The code was translated from Fortran to C and simplified for computing only the optimal knapsack value. For some variants of the search methods, the packing scheme also had to be reported. The behaviour of the knapsack procedure for different test data was consistent with the results reported in [10]. Obviously, the generation of items within a smaller weight range makes the problem more difficult (and takes more time) because a larger number of items is included in an optimal knapsack. For problems with 'large' items the set of candidates for an optimal packing is naturally smaller and an optimal solution can be found faster.
386
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
A similar argument shows that the solution of knapsack problems corresponding to small parameter values takes more time because the cost coefficients are closer to each other, whereas problems with high parameter values have larger gaps between the costs making a faster solution possible. As these properties are not specific for the inverse-parametric knapsack problem we fixed the range of the weights as described above and will not discuss these aspects any further.
Test environment The experiments were performed on a DIGITAL DECstation 5000/240. All algorithms were implemented in the standard C programming language. For all sets of parameters in Tables 1-5, 20 different random examples were generated and the average performance data reported. Tests with the same parameters in different tables originate from the same example. Time is always taken as CPU-time in seconds.
4.2. Experiments with pseudopolynomial algorithms Dimensions of the problems The number of items n was chosen from 1000 up to 10000. The knapsack capacity b was set between 2 000 and 15 000. Larger problems caused in some cases excessive running times and were solved only with search methods. As the magnitude of v had no substantial influence on the behaviour of algorithms, we fixed RATIO = 0.8 in this subsection. Experiments were performed on three types of algorithms as described in Section 3.1. Let us recall their characteristics: 1. Algorithm PPI: Carry out each comparison immediately: Either ? ~ P and the precedence is obvious or a knapsack problem is solved, which determines the outcome of the comparison and yields a new interval P. 2. Algorithm PP2: Collect all intersection points from one stage. Perform binary search to get a smaller interval P on which the comparisons can be determined. 3. Algorithm PP3: Collect the intersection points of maxstage stages. The knapsack functions computed during this collection phase are piecewise linear and may consist of up to 2 to the power of maxstage segments. Perform binary search as above. Reduce each function of the last stage to the segment valid in the new interval. Each algorithm was performed by trying at first approximate solutions for the arising knapsack problems. If the approximation did not settle the outcome of the comparison, the exact knapsack solution was computed. As most of the running time of the approximation algorithm is spent on the sorting of the items, which is also necessary for the exact solution, the amount of 'wasted' time used for trying the approximation in vain was very small. However, the gain of time achieved by using the approximation instead of the exact solution justified this approach. The necessity to compute exact solutions arose mainly during the last few breakpoint calculations. Once a complete knapsack computation had to be performed, the subsequent breakpoints tended to make further complete solutions necessary more frequently, due to the small range of P at this point of the computation. We tried to exploit this behaviour by switching off the approximation module after one or several 'failures' and continue with exact solutions. However, this variant performed hardly differently, sometimes slightly worse, than the program always trying the approximation. Table 1 shows the running time data. It turned out that the amount of time spent on solving or approximating the knapsack problems was very small compared to the running time of the main program.
387
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
Table 1 Performance of PP1 and PP2: For both algorithms total running time, time used for the approximate resp. exact solution of the arising knapsack problems and the number of generated breakpoints are given n b PP1 PP2
1000 1000 1000 2000 2000 2000 4000 4000 4000 6000 6000 6000 6000 8000 8000 8000 10000 10000
2000 4000 6000 2000 4000 6000 2000 4000 6000 6000 8000 10000 12000 8000 12000 16000 10000 15000
Total time (sec.)
Knapsack approx, time
Knapsack solution time
Breakpoints
Total time (sec.)
Knapsack approx, time
Knapsack solution time
Breakpoints
6.6 13.4 21.9 13.3 28.3 43.9 26.5 56.5 88.9 132.1 179.5 224.5 269.4 238.4 366.2 494.2 374.3 575.1
0.22 0.30 0.28 0.56 0.60 0.63 1.17 1.44 1.44 2.49 2.71 2.45 1.97 3.51 3.37 3.15 5.65 5.98
0.04 0.05 0.07 0.08 0.16 0.16 0.19 0.17 0.22 0.40 0.46 0.46 0.42 0.66 0.68 0.62 0.79 0.98
20.6 28.2 28.5 25.1 28.9 30.0 23.9 27.6 28.2 31.1 33.9 31.3 25.5 31.8 30.4 28.7 38.2 41.4
8.2 16,5 27.2 16,5 34.7 54.2 33.0 70.3 108.3 163.9 221.1 279.2 335.2 293.4 446.2 605.8 465.3 706.2
0.23 0.29 0.28 0.59 0.60 0.68 1.25 1.46 1.47 2.71 2.55 2.45 2.11 3.68 3.57 3.29 6.15 6.46
0.04 0.05 0.07 0.09 0.16 0.16 0.19 0.18 0.20 0.41 0.45 0.50 0.48 0.70 0.67 0.68 0.71 0.79
21.7 28.6 28.7 26.1 28.7 31.4 25.0 28.2 28.4 33.1 31.6 30.8 27.3 33.2 31.4 29.8 23.4 43.1
A n i n t e r e s t i n g a s p e c t o f t h e e x p e r i m e n t s is t h e n u m b e r o f g e n e r a t e d i n t e r s e c t i o n p o i n t s ~ a n d w h e r e t h e y w e r e d e t e c t e d by t h e a l g o r i t h m . It t u r n e d o u t t h a t t h e t o t a l n u m b e r of e n c o u n t e r e d b r e a k p o i n t s was surprisingly small (see T a b l e 1). A vast m a j o r i t y o f this n u m b e r was g e n e r a t e d in t h e first few stages w h e n t h e large s t a r t i n g interval was quickly r e d u c e d to a m u c h s m a l l e r one. A f t e r only 3 % o f t h e stages t h e r a n g e of P was always r e d u c e d to t * + 5 % a n d a f t e r less t h a n 10% o f t h e s t a g e s t * was a p p r o x i m a t e d by P w i t h i n a r a n g e o f 1%. H a r d l y any b r e a k p o i n t s o c c u r r e d in t h e s e c o n d half o f t h e stages. This m e a n s t h a t a small n u m b e r o f a r b i t r a r y i t e m s gives a close d e s c r i p t i o n o f t *. T h e r e f o r e , t h e a l g o r i t h m s c o u l d b e s t o p p e d e a r l y a n d w o u l d p r o d u c e a r e a s o n a b l e a p p r o x i m a t i o n to t h e o p t i m a l p a r a m e t e r . D a t a a b o u t this a s p e c t is given in T a b l e 2. N o t e t h a t P P 2 exhibits a b e h a v i o u r similar to PP1. T h e fact t h a t t h e n u m b e r o f e n c o u n t e r e d b r e a k p o i n t s was small in t o t a l a n d very small for a single s t a g e (the l a r g e s t n u m b e r o f d i f f e r e n t b r e a k p o i n t s f o u n d within o n e stage was 12) m a d e t h e t h e o r e t i c a l i m p r o v e m e n t of a l g o r i t h m P P 2 useless, as t h e c o l l e c t i o n o f all b r e a k p o i n t s o f o n e stage usually p r o d u c e d only 5 to 8 points. Thus, the a d d i t i o n a l effort o f collecting t h e s e b r e a k p o i n t s a n d u p d a t i n g t h e d a t a s t r u c t u r e at t h e e n d o f e a c h b r e a k p o i n t e v a l u a t i o n was l a r g e r t h a n t h e gain r e a c h e d by saving a few b r e a k p o i n t c o m p u t a t i o n s . F u r t h e r m o r e , s o m e o f t h e b r e a k p o i n t s in P P 2 w o u l d n o t have b e e n d e t e c t e d at all if P w e r e u p d a t e d i m m e d i a t e l y a f t e r a b r e a k p o i n t was f o u n d (as d o n e by PP1). In m a n y cases t h e s a m e b r e a k p o i n t was c o m p u t e d for several consecutive e n t r i e s o f a stage. This is i n d e e d p l a u s i b l e c o n s i d e r i n g t h a t t h e d e c i s i o n o f t e n has to b e d r a w n b e t w e e n t h e s a m e two p a c k a g e s c h e m e s if t h e k n a p s a c k c a p a c i t y is c h a n g e d by one. T h e s e m u l t i p l e b r e a k p o i n t s w e r e i m m e d i a t e l y e l i m i n a t e d by PP1 b u t h a d to b e s t o r e d in P P 2 for u p d a t i n g at t h e e n d o f t h e stage. T o c h e c k w h e t h e r t h e s e a d v a n t a g e s o f PP1 also h o l d in c o m p a r e n c e with PP3, w e t e s t e d this m o r e
ILE. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
388
Table 2 Stage where the last breakpoint is generated and first stage where P is contained in t * +5% resp. t * + 1%
n
b
PP1
PP2
Last breakpoint
Breakpoint t * + 5%
Breakpoint t * + 1%
Last breakpoint
Breakpoint t * + 5%
Breakpoint t * + 1%
1000 1000 1000 2 000 2 000 2000 4000 4000 4000
2 000 4 000 6000 2000 4 000 6000 2000 4000 6000
666 647 733 1066 1359 1456 3065 2265 1717
29 22 21 55 26 23 47 28 22
111 54 34 156 66 41 105 56 55
666 647 733 1066 1359 1456 3065 2265 1717
28 22 21 50 25 23 47 26 20
111 51 33 147 64 41 110 56 55
6000 6000 6000 6000 8000 8 000 8 000 10 000 10 000
6000 8000 10000 12000 8000 12 000 16 000 10 000 15 000
3198 3455 3747 2218 4521 4 057 3122 4 865 4 680
32 27 19 18 29 24 22 24 23
60 48 35 43 43 42 33 45 49
3198 4510 3747 2218 4521 4 057 3122 4 865 4 680
32 24 18 18 29 24 22 24 23
60 45 35 40 43 47 33 45 49
s o p h i s t i c a t e d a l g o r i t h m with d i f f e r e n t c r i t e r i a to d e t e r m i n e how m a n y s t a g e s s h o u l d b e c o l l e c t e d for e a c h p h a s e . T h r e e p a r a m e t e r s w e r e used: 1. maxsegrnent: M a x i m a l n u m b e r o f d i f f e r e n t line s e g m e n t s o f which o n e k n a p s a c k f u n c t i o n z ( k , d) consists. 2. maxbreak: M a x i m a l n u m b e r o f d i f f e r e n t b r e a k p o i n t s g e n e r a t e d in o n e p h a s e . 3. maxstage: M a x i m a l n u m b e r o f s t a g e s in o n e p h a s e . T h e c o l l e c t i o n p r o c e s s was s t o p p e d ( a n d t h e p h a s e e s t a b l i s h e d t h e r e b y e v a l u a t e d ) as soon as o n e o f t h e s e limits was r e a c h e d . A s i n d i c a t e d a b o v e a l a r g e n u m b e r o f b r e a k p o i n t s a r e f o u n d in t h e first few stages. T h e e x p e r i m e n t s s h o w e d t h a t it is m o r e efficient n o t to collect m a n y b r e a k p o i n t s at t h e b e g i n n i n g b u t to d e c r e a s e t h e interval P very e a r l y by s t a r t i n g with p h a s e s consisting o f only a few stages by c h o o s i n g e i t h e r maxsegment very small (e.g. < 4) o r limiting maxbreak (e.g. < b / 1 0 0 ) . B e c a u s e t h e g e n e r a t e d e x a m p l e s h a d a fairly r e g u l a r s t r u c t u r e , k n a p s a c k f u n c t i o n s consisting o f a n e x c e p t i o n a l l y large n u m b e r o f line s e g m e n t s d i d n o t occur, a n d b o u n d i n g any o f t h e two p a r a m e t e r s y i e l d e d r o u g h l y t h e s a m e results. H o w e v e r , by varying maxbreak a m o r e d e t a i l e d analysis o f t h e b e h a v i o u r o f t h e a l g o r i t h m c o u l d b e p e r f o r m e d . T h e r e f o r e , n o extensive e x a m i n a t i o n c o n c e r n i n g t h e i n f l u e n c e o f maxsegment was p e r formed. By t h e a b o v e s t r a t e g y P will b e very small a f t e r at m o s t t h r e e o r f o u r p h a s e s , which d r a m a t i c a l l y r e d u c e s t h e n u m b e r o f b r e a k p o i n t s f o u n d in t h e c o m p u t a t i o n t h e r e a f t e r . T h e b r e a k p o i n t s c o m p u t e d a f t e r t h e s e initial p h a s e s l e d to k n a p s a c k f u n c t i o n s consisting only o f a small n u m b e r o f line s e g m e n t s ( h a r d l y ever r e a c h i n g maxsegment). D u r i n g this m a i n p e r i o d o f t h e a l g o r i t h m it p r o v e d to b e m o s t efficient to c o m b i n e all stages into only o n e o r two p h a s e s . N o t e t h a t t h e last p h a s e d o e s n o t have to b e e v a l u a t e d at all b e c a u s e t * c a n b e f o u n d d i r e c t l y by i n t e r s e c t i n g z ( k , d) with L(t), which m a k e s a very l a r g e last p h a s e e s p e c i a l l y attractive. Although the different schemes of stage combination produced quite different patterns of behaviour
R.E. Burkard, U. Pferschy / European Journal of Operational Research 83 (1995) 376-393
389
Table 3 Behaviour of PP3 with different stage collection methods: Depending on maxbreak and maxstage the total running time, the time spent on the arising Knapsack subproblems, the number of generated breakpoints and the number of phases is given n
1000
b
4000
max-
max-
break
stage
Total time
40
200 400 1333 400 1333 400 1333
200 400 2 000
6 000
60
300 600 4000
8 000
80
400 800
Knapsack time
Number of breakpoints
Number of phases
22.3 22.3 22.2 23.1 23.0 23.8 24.1
0.29 0.26 0.25 0.26 0.20 0.24 0.22
23.4 22.1 21.2 21.8 16.9 20.2 18.4
7.7 5.6 4.5 4.8 3.1 4.3 3.1
300 600 2 000 600 2 000 600 2 000
72.8 72.9 72.6 73.8 73.6 74.7 74.9
0.68 0.61 0.48 0.59 0.46 0.57 0.39
25.2 22.7 18.5 22.0 17.8 21.9 15.6
9.4 6.3 3.9 5.4 3.0 5.2 2.6
400 800 2666 800 2666 800 2666
195.4 195.3 196.5 197.8 199.2 202.7 202.3
1.84 1.61 1.49 1.57 1.47 1.64 1.53
30.2 26.1 24.5 25.4 23.8 26.3 24.7
12.4 7.9 5.1 7.2 4.3 7.3 4.0
with notable changes in the number and structure of generated breakpoints, they surprisingly had only minor differences in the total running time (see Table 3). As the running-time of PP3 was clearly inferior to PP1, we only give the data of some instances with different stage collecting parameters in Table 3. The performance of PP2 was somewhere between the other two variants but closer towards PP1, which was overall the significantly fastest pseudopolynomial approach. Summarizing, we conclude that the more complicated but asymptotically faster algorithms which save a number of knapsack computations, are clearly inferior to the straight-forward approach. This may not correlate with the behaviour of other applications of Megiddo-based algorithms (of which no practical experiments are known to the authors), but is due to the fact that on one hand the subproblems in our case were efficiently solvable, and on the other hand, only a moderate number of these subproblems had to be solved. Hence, the core of the algorithm (solving the problem with parametric costs within a given parameter range) consumed most of the running time. Therefore, strategies to decrease the number of subproblems did not improve the performance of the algorithm but due to the more complicated data structures and bookkeepng procedures actually increased the total running time. An improvement of the main part of the algorithm does not seem to be at hand when following the Megiddo method.
4.3. Experiments with search methods Dimensions of the problems Much larger problems could be solved in reasonable time by search methods than by pseudopolynomial algorithms. The number of items n was chosen from 5 000 up to 40 000 and the knapsack capacity b set between 5 000 and 80 000.
390
R.E. Burkard, U. Pferschy/ European Journal of Operat•nal Research 83 (1995)376-393
Table 4 Performance of search methods with different values of RATIO n
b
RATIO
20 000
40 000
0.1 1 2
40000
80000
0.1 1 2
Method (BSS)
(BSB)
(NSL)
(NSU)
(SAS)
time iterations time iterations time iterations
4.29 10.2 2.23 5.8 1.85 5.2
2.87 7.2 1.86 4.9 1.83 5.0
1.58 4.5 1.37 4.l 1.38 4.2
1.38 4.1 0.92 3.2 0.91 3.0
3.42 8.4 2.32 6.0 2.28 6.0
time iterations time iterations time iterations
9.71 10.8 4.93 6.1 3.95 5.0
5.93 7.2 3.95 4.9 3.94 4.9
3.55 4.7 2.96 4.0 2.95 4.0
2.97 4.1 1.98 3.1 1.96 3.0
6.92 8.0 4.90 5.9 4.96 6.0
As mentioned above the solution of knapsack problems takes more time for smaller parameter values, whereas problems corresponding to large parameters are usually easier to solve. Therefore, the knapsack computations for problems with smaller values of RATIO tended to take longer because they naturally imply smaller values for t *. Furthermore, also the number of iterations was slightly higher for low RATIOs, especially for the binary search based methods (BSS) and (BSB), which is due to the fact that the search methods work better for strongly ascending parts of the value functions, which are present for higher parameter values (see Table 4 for selected examples). The initial interval P was set to [0. . . . . T] with a large upper bound T fixed according to previous experiences (e.g. T = 2 000). Because of the rapid improvement of P in the first iterations the choice of T hardly effected the performance of the algorithms and is not critical in practice. To choose a new search value ? inside P = [tmin. . . . . tmax], we tested three variants as described in Section 3.2. The A-partitioning was restricted to the classical binary search (A = ½) because other versions were either inferior (for A close to 0 or 1) or performed similarly to the bisection. More important than the choice of ? is how to take advantage of the properties of the search function, especially how to decrease the length of P between the knapsack computations by the bounds given in Step 4 of the general search procedure. P was updated by setting
tmax:=min{~tmmnL(t),
/,maxn L(t)}
and tmin :~
[linear function through
(tmin, U(tmin) ) a n d (tmax, U ( t m a x ) ) ] n L ( t ) .
As indicated in Step 5, any search procedure can be stopped as soon as the slopes of b~tm and Ct~" are equal (or differ only by a small factor). In this case intersecting t~tm, with L(t) yields t~e optimal parameter t * The improvement of P by the above bounding techniques is especially useful in the first few iterations when P is still quite large. Thereby, a considerable number of knapsack evaluations can be saved, although part of this effect is diminished by the necessity to determine ft, which requires the exact knapsack solution at the given point and makes the use of approximation methods for the knapsack computations impossible. After a couple of iterations the interval P was usually rather small which made the use of the bounds less effective. As also the approximation algorithms were not precise enough in this phase of the
R.E. Burkard, U. Pferschy/ EuropeanJournalof OperationalResearch83 (1995)376-393
391
a l g o r i t h m s , t h e exact k n a p s a c k s o l u t i o n h a d to b e c o m p u t e d anyway in m o s t cases. H e n c e , t h e b o u n d c o m p u t a t i o n s d i d n o t c a u s e m u c h a d d i t i o n a l effort a n d w e r e c a r r i e d o u t until t h e e n d o f t h e search. I n c l u d i n g t h e c a l c u l a t i o n o f gt a n d t h e b o u n d s t h e m a i n p r o g r a m t o o k always less t h a n 5 % o f t h e t o t a l r u n n i n g t i m e with t h e r e m a i n d e r b e i n g s p e n t o n t h e k n a p s a c k c o m p u t a t i o n s . H e n c e , we r e p o r t only t h e t o t a l r u n n i n g times. C h a r a c t e r i z e d by t h e choice o f t" six d i f f e r e n t a l g o r i t h m s w e r e finally t e s t e d with t h e following properties and behaviour: (BS): pure binary search. This v a r i a n t r e a c h e d t e r m i n a t i o n only if tmax - tmin < e. O f c o u r s e t h e n u m b e r o f i t e r a t i o n s d e p e n d e d m a i n l y o n t h e size o f e ( r e a s o n a b l e choices s e e m to be, e.g. e = 0.001) a n d was a l m o s t i d e n t i c a l for d i f f e r e n t p r o b l e m s . T h e possibility to use a p p r o x i m a t i o n a l g o r i t h m s for t h e arising k n a p s a c k p r o b l e m s gave s o m e a d v a n t a g e s to (BS), b u t still t h e p e r f o r m a n c e o f t h e other algorithms could not be matched. (BSS): binary search with slope computation. (BSS) is d e r i v e d f r o m (BS) by a d d i n g as s t o p p i n g c r i t e r i o n t h e e - e q u i v a l e n c e o f t h e gt s l o p e s as d e s c r i b e d above. This saves a n u m b e r o f i t e r a t i o n s at t h e e n d o f t h e a l g o r i t h m a n d l e a d s to a s h o r t e r r u n n i n g time t h a n (BS), a l t h o u g h the c o m p l e t e s o l u t i o n o f t h e k n a p s a c k p r o b l e m s is n e c e s s a r y in every iteration. (BSB): binary search with slope computation and bounding. By a d d i n g to (BSS) t h e b o u n d i n g p r o c e d u r e s d e s c r i b e d above, a very efficient a l g o r i t h m s is a t t a i n e d . T h e t o t a l n u m b e r o f i t e r a t i o n s was always less t h a n 10. N o t e t h a t t°t is only a v a i l a b l e for a b o u n d o f P which was g e n e r a t e d as s e a r c h v a l u e f in a p r e v i o u s i t e r a t i o n , b e c a u s e k n a p s a c k solutions a r e not k n o w n for interval b o u n d s /impr d e r i v e d by b o u n d i m p r o v e m e n t s . H e n c e , in p r a c t i c e we can use only r e l a x e d v e r s i o n s o f t h e full i m p r o v e m e n t by using, e.g. for timpr t h e l i n e a r f u n c t i o n m a t c h i n g v ( t ) at timpr with t h e slope t a k e n f r o m gton w h e r e torig is t h e o r i g i n a l b o u n d . (NSL): Newton type search starting from O. T h e N e w t o n m e t h o d c a n n o t m a k e use of any b o u n d s , b e c a u s e •
g
Table 5 Total running time and number of iterations for different search methods (RATIO = 0.8) n
b 5 000
5 000
5 000
10000
10 000
15 000
10000
30000
20 000
20 000
20000
40 000
20 000
60 000
40 000
40 000
40000
60 000
40 000
80 000
Method time iterations time iterations time iterations time iterations time iterations time iterations time iterations time iterations time iterations time iterations
(BS)
(BSS)
(BSB)
(NSL)
(NSU)
(SAS)
0.83 11.8 0.84 11.9 1.90 12.0 1.91 12.0 4.26 12.0 4.28 12.0 4.25 12.0 9.27 12.0 9.31 12.1 9.32 12.2
0.41 5.6 0.53 6.7 1.08 6.2 1.24 6.8 2.12 5.6 2.40 6.3 2.59 6.5 4.95 6.1 4.98 6.0 4.96 6.1
0.37 4.9 0.38 4.9 0.82 5.0 0.83 5.0 1.82 4.9 1.84 4.9 1.83 5.0 3.91 4.9 3.97 5.0 3.92 5.0
0.28 4.1 0.29 4.2 0.63 4.2 0.64 4.2 1.37 4.1 1.39 4.2 1.39 4.3 2.94 4.0 2.95 4.2 2.99 4.3
0.19 3.2 0.20 3.2 0.48 3.3 0.42 3.1 0.92 3.2 0.94 3.2 0.92 3.0 1.98 3.0 1.98 3.1 1.99 3.0
0.48 6.1 0.45 6.0 1.03 6.0 1.05 5.9 2.33 6.0 2.37 6.2 2.31 6.1 4.79 6.1 4.59 5.9 4.63 6.0
392
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
only the slope at the last ~' is used and no interval improvement makes sense. Still it performed very well and competed successfully with (BSB). (NSU): Newton type search starting f r o m T. In most cases the Newton method took one or two iterations less if started on the right hand side of t * rather than at 0. This approach was the best method in all experiments. There was a tendency observed that (NSU) took less iterations than (NSL) if the slopes c i were chosen between 0 and infinity, which causes a very steep value function, whereas (NSL) was superior if the slopes were restricted to a smaller interval, e.g. (0 . . . . . tan('rr/16)). (SAS): Sub-approximation search. The performance of this approach was varying more than the other methods. Although it sometimes produced solutions surprisingly fast, in some cases it was the worst of all variants. However, also for this method the number of iterations was very stable. As can be seen on Table 5 the algorithms performed regularly with a modest sublinear increase of the running time with respect to the problem size. We summarize that the inverse-parametric knapsack problem can be solved very efficiently in a small number of iterations by search methods exploiting the properties of the optimal value function. In contrast to the pseudopolynomial methods, in this case the more sophisticated (though still not too complicate) algorithms produced the best results. Newton methods usually performed best, followed by (BSB). This might be explained by the regular structure of random test examples. Their smooth value function made the approximation of the function by ft as employed by the Newton method highly successful. We suspect that in examples with more irregular data, the advantages of (BSB) will gain more importance. However, sets of linear functions generally lead to more regular shapes of value functions than constant cost coefficients. Degenerate examples usually lead to difficult problems only for a special range of v and are easy to solve for all other target values.
5. Conclusions In this paper we treated the inverse-parametric problem where the costs (resp. the weights) of a knapsack problem are replaced by linear functions in a parameter t and the smallest possible parameter t * has to be found such that the optimal solution of the knapsack problem taken at t * is equal to a given input value. Based on a general method of Megiddo, we constructed a pseudopolynomial algorithm for the parametric cost case with running time O(n2b log b / l o g ( n log b)). Exploiting the specific properties of the parametric value function, we developed various geometrically motivated search methods. Computational studies showed that the more complicated and theoretically improved pseudopolynomial methods performed worse than a straight-forward approach which is due to the fact that the solution of the arising subproblems (knapsack problems) can either be avoided in many cases or computed efficiently. Comparative analysis of the search methods showed that even large problems could be solved very fast, especially by Newton-based approaches.
References [1] Burkard, R.E., Dlaska, K., and Kellerer, H., "The quickest disjoint flow problem", Report 189, Institute of Mathematics, Technical UniversityGraz, Austria, 1991.
R.E. Burkard, U. Pferschy/ European Journal of Operational Research 83 (1995) 376-393
393
[2] Burkard, R.E., Dlaska, K., and Klinz, B., "The quickest flow problem", Zeitschrift fiir Operations Research 37 (1993) 31-58. [3] Carstensen, P.J., "Complexity of some parametric integer and network programming problems", Mathematical Programming 26 (1983) 64-75. [4] Chen, Y.L., and Chin, Y.H., "The quickest path problem", Computers & Operations Research 17 (1990) 153-161. [5] Hashizume, S., Fukushima, M., Katoh, N., and Ibaraki, T., "Approximation algorithms for combinatorial fractional programming problems", Mathematical Programming 37 (1987) 255-267. [6] Hoppe, B., and Tardos, E., "The quickest transshipment problem", to appear in: Proc. of the Sixth AnnualACM-SIAM Symp. on Discrete Algorithms, 1995. [7] Ibaraki, T., "Parametric approaches to fractional programs", Mathematical Programming 26 (1983) 345-362. [8] Ishii, H., Ibaraki, T., and Mine, H., "Fractional knapsack problems", Mathematical Programming 13 (1977) 255-271. [9] Lee, D.T., and Papadopoulou, E., "The all-pairs quickest path problem", Information Processing Letters 45 (1993) 261-267. [10] Martello, S., and Toth, P., Knapsack Problems, Wiley, New York, 1990. [11] Martello, S., and Toth, P., "An upper bound for the zero-one knapsack problem and a branch and bound algorithm", European Journal of Operational Research 1 (1977) 169-175. [12] Megiddo, N., "Combinatorial optimization with rational objective functions", Mathematics of Operations Research 4 (1979) 414-424. [13] Megiddo, N., "Applying parallel computation algorithms in the design of serial algorithms", Journal of the ACM 30 (1983) 852-865. [14] Murty, K.G., "Computational complexity of parametric linear programming", Mathematical Programming 19 (1980) 213-219. [15] Nemhauser, G.L., and Wolsey, L.A., Integer and Combinatorial Optimization, Wiley, New York, 1988. [16] Woeginger, G.J., Private communication, 1994. [17] Zadeh, N., "A bad network problem for the simplex method and other minimum cost flow algorithms", Mathematical Programming 5 (1973) 255-266.