A multi-criteria approach to approximate solution of multiple-choice knapsack problem

We propose a method for finding approximate solutions to multiple-choice knapsack problems. To this aim we transform the multiple-choice knapsack problem into a bi-objective optimization problem whose solution set contains solutions of the original multiple-choice knapsack problem. The method relies on solving a series of suitably defined linearly scalarized bi-objective problems. The novelty which makes the method attractive from the computational point of view is that we are able to solve explicitly those linearly scalarized bi-objective problems with the help of the closed-form formulae. The method is computationally analyzed on a set of large-scale problem instances (test problems) of two categories: uncorrelated and weakly correlated. Computational results show that after solving, in average 10 scalarized bi-objective problems, the optimal value of the original knapsack problem is approximated with the accuracy comparable to the accuracies obtained by the greedy algorithm and an exact algorithm. More importantly, the respective approximate solution to the original knapsack problem (for which the approximate optimal value is attained) can be found without resorting to the dynamic programming. In the test problems, the number of multiple-choice constraints ranges up to hundreds with hundreds variables in each constraint.

In this paper, we adopt the convention that a vector x is a column vector, and hence the transpose of x, denoted by x T , is a row vector. Problem (M CKP ) is of the form x ij ∈ {0, 1} i = 1, ..., k, j = 1, ..., n i }.
By using the above notations, problem (M CKP ) can be equivalently rewritten in the vector form where p and c are vectors from R n , p := (p 11 , p 12 , ..., p 1n1 , p 21 , ..., p 2n2 , ...., p k1 , p k2 , ...p kn k ) T c := (c 11 , c 12 , ..., c 1n1 , c 21 , ..., c 2n2 , ...., c k1 , c k2 , ...c kn k ) T , and for any vectors u, v ∈ R n , the scalar product u T v is defined in the usual way as u T v := n i=1 u i v i . The feasible set F to problem (M CKP ) is defined by a single linear inequality constraint and the constraint x ∈ X, i.e., The optimal value of problem (M CKP ) is equal to max x∈F p T x and the solution set S * is given as Problem (M CKP ) is N P-hard. The approaches to solving (M CKP ) can be: heuristics [1,12], exact methods providing upper bounds for the optimal value of the profit together with the corresponding approximate solutions [26], exact methods providing solutions [18]. There are algorithms that efficiently solve (M CKP ) without sorting and reduction [8,28] or with sorting and reduction [4]. Solving (M CKP ) with a linear relaxation (by neglecting the constrains x ij ∈ {0, 1}, i = 1, ..., k, j = 1, ..., n i ) gives upper bounds on the value of optimal profit. Upper bounds can be also obtained with the help of the Lagrange relaxation. These facts and other features of (M CKP ) are described in details in monographs [13,19].
Exact branch-and-bound methods [6] (integer programming), even those using commercial optimization software (e.g., LINGO, CPLEX) can have troubles with solving large (M CKP ) problems. A branch-and-bound algorithm with a quick solution of the relaxation of reduced problems was proposed by Sinha and Zoltners [27]. Dudziński and Walukiewicz proposed an algorithm with pseudo-polynomial complexity [5].
Algorithms that use dynamic programming require integer values of data and for large-scale problems require large amount of memory for backtracking (finding solutions in set X), see also the monograph [19]. The algorithm we propose does not need the data to be integer numbers.
Heuristic algorithms, based on solving linear (or continuous) relaxation of (M CKP ) and dynamic programming [7,22,24] are reported to be fast, but have limitations typical for dynamic programming.
The most recent approach "reduce and solve" [2,10] is based on reducing the problem by proposed pseudo cuts and then solving the reduced problems by a Mixed Integer Programming (M IP ) solver.
In the present paper, we propose a new exact (not heuristic) method which provides approximate optimal profits together with the corresponding approximate solutions. The method is based on multi-objective optimization techniques. Namely, we start by formulating a linear bi-objective problem (BP ) related to the original problem (M CKP ). After investigating the relationships between (M CKP ) and (BP ) problems, we propose an algorithm for solving (M CKP ) via a series of scalarized linear bi-objective problems (BS(λ)).
The main advantage of the proposed method is that the scalarized linear bi-objective problems (BS(λ)) can be explicitly solved by exploiting the structure of the set X. Namely, these scalarized problems can be decomposed into k independent subproblems the solutions of which are given by simple closed-form formulas. This feature of our method is particularly suitable for parallelization. It allows to generate solutions of scalarized problems in an efficient and fast way.
The experiments show that the method we propose generates very quickly an outcomex ∈ F which is an approximate solution to (M CKP ). Moreover, lower bound (LB) and upper bound (UB) for the optimal profit are provided.
The obtained approximate solutionx ∈ F could serve as a good starting point for other, e.g., heuristic or exact algorithms for finding an optimal solution to the problem (M CKP ).
The organization of the paper is as follows. In Section 2, we provide preliminary facts on multi-objective optimization problems and we formulate a bi-objective optimization problem (BP ) associated with (M CKP ). In Section 3, we investigate the relationships between the problem (BP ) and the original problem (M CKP ). In Section 4, we formulate scalarized problems (BS(λ)) for bi-objective problem (BP ) and we provide closed-form formulae for solutions to problems (BS(λ)) by decomposing them into k independent subproblems (BS(λ)) i , i = 1, ..., k. In Section 5, we present our method (together with the pseudo-code) which provides a lower bound (LB) for the optimal profit together with the corresponding approximate feasible solutionx ∈ F to (M CKP ) for which the bound (LB) is attained. In Section 6, we report on the results of numerical experiments. The last section concludes.

Multi-objective optimization problems
Let f i : R n → R, i = 1, ..., k, be functions defined on R n and Ω ⊂ R n be a subset in R n .
The multi-objective optimization problem is defined as where the symbol V max means that solutions to problem (P ) are understood in the sense of Pareto efficiency defined in Definition 2.1.
The problem (P ) where all the functions f i , i = 1, ..., k are linear is called a linear multi-objective optimization problem.

Remark 2.1
The bi-objective problem in the sense that Pareto efficient solution sets (as subsets of the feasible set Ω) coincide and Pareto elements (the images in R 2 of Pareto efficient solutions) differ in sign in the second component.

A bi-objective optimization problem related to (M CKP )
In relation to the original multiple-choice knapsack problem (M CKP ), we consider the linear bi-objective binary optimization problem (BP 1) of the form In this problem the left-hand side of the linear inequality constraint c T x ≤ b of (M CKP ) becomes a second criterion and the constraint set reduces to the set X. There are two-fold motivations of considering the bi-objective problem (BP 1).
First motivation comes from the fact that in (M CKP ) the inequality is usually seen as a budget (in general: a resource) constraint with the lefthand-side to be preferably not greater than a given available budget b. In the bi-objective problem (BP 1), this requirement is represented through the minimization of ni j=1 c ij x ij . In Theorem 3.1 of Section 3, we show that under relatively mild conditions among solutions of the bi-objective problem (BP 1) (or the equivalent problem (BP )) there are solutions to problem (M CKP ).
Second motivation is important from the algorithmic point of view and is related to the fact that in the proposed algorithm we are able to exploit efficiently the specific structure of the constraint set X which contains k linear equality constraints (each one referring to a different group of variables) and the binary conditions only. More precisely, the set X can be represented as the Cartesian product of the sets i.e., Note that due to the presence of the budget inequality constraint the feasible set F of problem (M CKP ) cannot be represented in the form analogous to (3).
According to Remark 2.1, problem (BP 1) can be equivalently reformulated in the form

The relationships between (BP ) and (M CP K)
Starting from the multiple-choice knapsack problem (M CKP ) of the form in the present section we analyse relationships between problems (M CKP ) and (BP ).
We start with a basic observation.
On the other hand, if b ≥ max x∈X c T x, (M CKP ) is trivially solvable. Thus, in the sequel we assume that Let P max := max x∈X p T x, i.e., P max is the maximal value of the function p T x on the set X. The following observations are essential for further considerations.
1. First, among the elements of X which realize the maximal value P max , there exists at least one which is feasible for (M CKP ), i.e., there exists Then, clearly, x p solves (M CKP ). 2. Second, none of elements which realize the maximal value P max is feasible for (M CKP ), i.e., for every x p ∈ X, p T x p = P max we have c T x p > b, i.e., any x p realizing the maximal value P max is infeasible for (M CKP ), i.e.
In the sequel, we concentrate on Case 2, characterized by (7). This case is related to problem (BP ). To see this let us introduce some additional notations. Let x cmin ∈ X and x pmax ∈ X be defined as Let S bo be the set of all Pareto solutions to the bi-objective problem (BP ), (c.f. Definition 2.1). The following lemma holds.
Lemma 3.1 Assume that we are in Case 2, i.e., condition (7) holds. There exists a Pareto solution to the bi-objective optimization problem (BP ),x ∈ S bo which is feasible to problem (M CKP ), i.e. c Tx ≤ b which amounts tox ∈ F .
Proof. According to Definition 2.1, both x pmax ∈ X and x cmin ∈ X are Pareto efficient solutions to (BP ), i.e., there is no Moreover, by (7), In view of (8), which means thatx is feasible to problem (M CKP ), which concludes the proof. Now we are ready to formulate the result establishing the relationship between solutions of (M CKP ) and Pareto efficient solutions of (BP ) in the case where the condition (7) holds.
Proof. Observe first that, by Lemma 3.1, there exist x ∈ S bo satisfying the constraint c T x ≤ b, i.e., condition ( * ) is not dummy. By contradiction, suppose that a feasible element x * ∈ F , i.e., x * ∈ X, c T x * ≤ b, is not a solution to (M CKP ), i.e., there exists an x 1 ∈ X, such that We show that x * cannot satisfy condition ( * ). If c T x 1 ≤ c T x * , then x * is not a Pareto solution to (BP ), i.e., x * ∈ S bo , and x * does not satisfy condition ( * ).
If x 1 ∈ S bo , then because x * ∈ S bo , then, according to (9), x * cannot satisfy condition ( * ). If and x * does not satisfy condition ( * ), a contradiction which completes the proof. Theorem 3.1 says that under condition (7) any solution to (BP ) satisfying condition ( * ) solves problem (M CKP ). General relations between constrained optimization and multi-objective programming were investigated in [15].
Basing ourselves on Theorem 3.1, in Section 5 we provide an algorithm for finding x ∈ S bo , a Pareto solution to (BP ), which is feasible to problem (M CKP ) and for which the condition ( * ) is either satisfied or is, in some sense, as close as possible to be satisfied. In this latter case, the algorithm provides upper and lower bounds for the optimal value of (M CKP ).

Decomposition of the scalarized bi-objective problem (BP )
In the present section, we consider problem (BS(λ 1 , λ 2 )) defined by (10) which is a linear scalarization of problem (BP ). In our algorithm BISSA, presented in Section 5, we obtain an approximate feasible solution to (M CKP ) by solving a (usually very small) number of problems of the form (BS(λ 1 , λ 2 )). The main advantage of basing our algorithm on problems (BS(λ 1 , λ 2 )) is that they are explicitly solvable by simple closed-form expressions (17) .
For problem (BP ) the following classical scalarization result holds.
Theorem 4.1 [9,20] If there exist λ > 0, = 1, 2, such that x * ∈ X is a solution to the scalarized problem then x * is a Pareto efficient solution to problem (BP ).
Without loosing generality we can assume that 2 l=1 λ = 1. In the sequel, we consider, for 0 < λ < 1, scalarized problems of the form Remark 4.1 According to Theorem 4.1, solutions to problems need not be Pareto efficient because the weights are not both positive. However, there exist Pareto efficient solutions to (BP ) among solutions to these problems. Namely, there exist ε 1 > 0 and ε 2 > 0 such that solutions to problems are Pareto efficient solutions to problems (12), respectively. Suitable ε 1 and ε 2 will be determined in the next subsection.

Decomposition
Due to the highly structured form of the set X and the possibility of representing X in the form (3), we can provide explicit formulae for solving problems (BS(λ)). To this aim we decompose problems (BS(λ)) as follows.
Recall that by using the notation (4) we can put any x ∈ X in the form By solving problems (BS(λ)) i , i = 1, ..., k, we find their solutionsx i . We shall show thatx := (x 1 , ...,x k ) T solves (BS(λ)). Thus, problem (11) is decomposed into k subproblems (13), the solutions of which form solutions to (11). Note that similar decomposed problems with feasible sets X i and another objective functions have already been considered in [3] in relation to multidimensional multiple-choice knapsack problems. Now we give a closed-form formulae for solutions of (BS(λ)) i . For i = 1, .., , k, let and let 1 ≤ j * i ≤ n i be the index number for which the value V i is attained, i.e., We show thatx i := (0, .., 1 is a solution to (BS(λ)) i and is a solution to (BS(λ)). The optimal value of (BS(λ)) is Namely, the following proposition holds.
Proof. Clearly,x i are feasible for (BS(λ)) i , i = 1, ..., k, becausex i is of the form (16) and hence belongs to the set X i which is the constraint set of (BS(λ)) i . Consequently,x * defined by (17) is feasible for (BS(λ)) because all the components are binary and the linear equality constraints ni j=1 x i j = 1, i = 1, 2, ..., k are satisfied. To see thatx i are also optimal for (BS(λ)) i , i = 1, ..., k, suppose by the contrary, that there exists 1 ≤ i ≤ k and an element y ∈ R ni which is feasible for (BS(λ)) i with the value of the objective function strictly greater than the value atx i , i.e., This, however, would mean that there exists an index 1 ≤ j ≤ n i such that contrary to the definition of j * .
To see thatx * is optimal for (BS(λ)), suppose by the contrary, that there exists an element y ∈ R n which is feasible for (BS(λ)) and the value of the objective function at y is strictly greater than the value of the objective function atx * , i.e., In the same way as previously, we get the contradiction with the definition of the components ofx * given by (17).
Let us observe that each optimization problem (BS(λ)) i can be solved in time O(n i ), hence problem (BS(λ)) can be solved in time O(n), where n = k i=1 n i . Clearly, one can have more than one solution to (BS(λ)) i , i = 1, ..., k. In the next section, according to Theorem 3.1, from among all the solutions of (BS(λ)) we choose the one for which the value of the second criterion is greater than and as close as possible to −b.
Note that by using Proposition 4.1, one can easily solve problems (P 1) and (P 2) defined in Remark 4.1, i.e., by applying (18) we immediately get the optimal values of (P 1) and (P 2) and by (17), we find their solutionsx 1 andx 2 , respectively. Proposition 4.1 and formula (17) allows to find ε 1 > 0 and ε 2 > 0 as defined in Remark 4.1. By (17), it is easy to find elementsx 1 ,x 2 ∈ X such that For any 1 ≤ i ≤ k, the submaximal values of a linear function (d i ) T x i on X i can be found by: ordering first the coefficients of the function (d i ) decreasingly, and next observing that the submaximal (i.e., smaller than maximal but as close as possible to the maximal) value of (d i ) on X i is attained for We have the following fact.

Bi-objective Approximate Solution Search Algorithm (BISSA) for solving (M CKP )
In this section, we propose the bi-objective approximate solution search algorithm BISSA, for finding an elementx ∈ F which is an approximate solution to (M CKP ). The algorithm relies on solving a series of problems (BS(λ)) defined by (11)   According to Theorem 4.1, each solution to (BS(λ)) solves the linear biobjective optimization problem (BP ), According to Theorem 3.1, any Pareto efficient solution x * to problem (BP ) which is feasible to (M CKP ), i.e., (−c) T x * ≥ −b, and satisfies condition ( * ), i.e., solves problem (M CKP ). Since problems (BS(λ)) are defined with the help of linear scalarization, we are not able, in general, to enumerate all x ∈ S bo such that (−c) T x + b ≥ 0 in order to find an x * which satisfy condition ( * ).
On the other hand, by using linear scalarization, we are able to decompose and easily solve problems (BS(λ)). The BISSA algorithm aims at finding a Pareto efficient solutionx ∈ X to (BP ) which is feasible to (M CKP ), i.e., c Tx ≤ b for which the value of b − c Tx is as small as possible (but not necessarily minimal) and approaches condition ( * ) of Theorem 3.1 as close as possible.
Here, we give a description of the BISSA algorithm. The first step of the algorithm (lines 1-5) is to find solutions to problems (P 1) and (P 2) as well as their outcomes. The solutions are the extreme Pareto solutions to problem (BP ). Those points named (a 1 , b 1 ) 0 and (a 2 , b 2 ) 0 are presented in Fig. 3. Then (lines [6][7][8][9], in order to assert whether a solution to problem (M CKP ) exists or not, a basic checking is made against value −b. If the algorithm reaches line 10, no solution has been found yet, and we can begin the exploration of the search space.
We calculate λ according to line 13. The value of λ is the slope of the straight line joining (a 1 , b 1 ) and (a 2 , b 2 ). At the same time it is the scalarization parameter defining the problem (BS(λ)) (formula (11)). The outcome of the solution to problem (BS(λ)) cannot lie below the straight line determined by points (a 1 , b 1 ) and (a 2 , b 2 ). It must lie on or above this line, as it is the Pareto efficient solution to problem (BP ). Then, problem (BS(λ)) is solved (line 14) by using formulae (16) and (17). Next, in lines 15-27 of the repeat-until loop a scanning of the search space is conducted to find solutions to problem (BP ) which are feasible to problem (M CKP ). If there exist solutions with outcomes lying above the straight line determined by λ (the condition in line 15 is true), either the narrowing of the search space is made (by determining new points  (a 1 , b 1 ) and (a 2 , b 2 ), see Fig.3, and points with upper index equal to 1), and the loop continues, or the solution to problem (M CKP ) is derived. If not, the solution x from set S which outcome lies above the line determined by −b (the feasible solution to problem (M CKP )) and for which value f 2 (x) + b is minimal in this set, is an approximate solution (x) to problem (M CKP ), and the loop terminates. Finally (line 28), the upper bound f 1 (x) + u on the profit value of exact solution to problem (M CKP ) is calculated.
Algorithm 1 BISSA -Approximate solution search to (M CKP ) 1: Calculate ε 1 , ε 2 according to (20) 2: Assume that f 1 (x) = p T x and f 2 (x) = (−c) T x 3: Solve (P 1) according to (18) and (17) x 1 a solution to (P 1) 4: Solve (P 2) according to (17) and (18) x 2 a solution to (P 2) 5: λ := Solve (BS(λ)) according to (16) and (17) x a solution, opt the optimal value, S the solution set to (BS(λ)) 15: if opt > α then 16: if f 2 (x) > −b then 17: The BISSA algorithm finds either an exact solution to problem (M CKP ), or (after reaching line 27) a lower bound (LB) with its solutionx and an upper bound (UB) (see Fig.3). A solution found by the algorithm is, in general, only an approximate solution to problem (M CKP ) because a triangle (called fur- ther the triangle of uncertainty) determined by points (f 1 (x), f 2 (x)), (f 1 (x) + u, −b), (f 1 (x), −b) may contain other Pareto outcomes (candidates for outcomes of exact solutions to problem (M CKP )) which the proposed algorithm is not able to derive. The reason is that we use a scalarization technique based on weighted sums of criteria functions to obtain Pareto solutions to problem (BP ).
Let us recall that each instance of the optimization problem (BS(λ)) can be solved in time O(n), but the number of these instances solved by the proposed algorithm depends on the size of the problem (values k and n i ) and the data.

Computational experiments
Most publicly available test instances refer not to the (M CKP ) problem (let us recall, that there is only one inequality or budget constraint in the problem we consider) but to multi-dimensional knapsack problems. Due to this fact we generate new random instances (available from the authors on request). However, to compare solutions obtained by the BISSA algorithm to the exact solutions we used the minimal algorithm for the multiple-choice knapsack problem [22] which we call EXACT and its implementation in C [23]. The EXACT algorithm gives the profit value of the optimal solution as well as the solution obtained by the greedy algorithm for the (M CKP ) problem, so the quality of the BISSA algorithm approximate solutions can be assessed in terms of the difference or relative difference between profit values of approximate solutions and exact ones.
Since the difficulty of knapsack problems (see, e.g., the monograph [19]) depends on the correlation between profits and weights of items, we conducted two computational experiments: Experiment 1 with uncorrelated data instances (easy to solve) and Experiment 2 with weakly correlated data instances (more difficult to solve) (c.f. [11]). The explanation why weakly correlated problems are more difficult to solve by the BISSA algorithm than uncorrelated ones we give later.
To prepare test problems (data instances) we used a method proposed in [22] and our own procedure for calculating total cost values. The BISSA algorithm has been implemented in C. The implementation of BISSA algorithm was run on off-the-shelf laptop (2GHz AMD processor, Windows 10), and the implementation of EXACT algorithm was run on PC machine (4x3.2GHz Intel processor, Linux). The running time for BISSA and EXACT algorithms for each of the test problems was below one second.
The contents of the tables columns containing experiment results is as follows.
1 -problem no. 2 -profit of the exact solution found by the EXACT algorithm. 3 -profit of the approximate solution found by the BISSA algorithm. 4 -difference between 2 and 3. 5 -relative (%) difference between 2 and 3. 6 -upper bound for (M CKP ) found by the BISSA algorithm. 7 -the difference between the upper bound and profit of the approximate solution. 8 -the relative difference between the upper bound and profit of the approximate solution. 9 -upper bound for (M CKP ) found by the greedy algorithm. 10 -number of (BS(λ)) problems solved by the BISSA algorithm.

Experiment 2 -weakly correlated (wco) data instances
We generated 10 test problems assuming that k = 20 and n i = 20, i = 1, ..., k (problem set (wco, 20, 20)). For each test problem costs (c ij ) of items in set N i were randomly distributed (according to the uniform distribution) in [1, R], R = 10000, and profits of items (p ij ) in this set were randomly distributed in [c ij − 10, c ij + 10], such that p ij ≥ 1. Profits and costs of items were integers. For each test problem the total cost b was calculated as for Experiment 1.
The results for problem set (wco, 20, 20) are given in Table 4.
In the case of uncorrelated data instances, the BISSA algorithm was able to find approximate solutions (and profit values) to problems with 10000 binary variables in reasonable time. The relative difference between profit values of exact and approximate solutions are small for each of the test problems. Upper bounds found by the BISSA algorithm are almost the same as upper bounds found by the greedy algorithm for (M CKP ). Even for the problem set (unc, 1000, 10) number of (BS(λ)) problems solved by the BISSA algorithm is relatively small in regards to number of decision variables.
In the case of weakly correlated data instances, the BISSA algorithm solved problems with 400 binary variables in reasonable time. The relative difference between profit values of exact and approximate solutions is, in average, greater than for uncorrelated test problems. As one can see in Table 4, upper bounds found by the BISSA algorithm are almost the same as upper bounds found by the greedy algorithm for (M CKP ). The reason why the BISSA algorithm solves weakly correlated instances with a significantly smaller number of variables than for uncorrelated ones in reasonable time is as follows. In line 24 of the BISSA algorithm, in order to find an elementx, we have to go through the solution set S to the problem (BS(λ)) (the complete scan of set S according to values of the second objective function of problem (BP )). For weakly correlated data instances the cardinality of the set S may be large even for problems belonging to class (wco, 30, 30). We conducted experiments for problem class (wco, 30, 30). For the most difficult test problem in this class, the cardinality of solution set S to the problem (BS(λ)) was 199,065,600. For greater weakly correlated problems that number may be even larger.

Conclusions and future works
A new approximate method of solving multiple-choice knapsack problems by replacing the budget constraint with the second objective function has been presented. Such a relaxation of the original problem allows to the smart scanning of the decision space by quick solving of the binary linear optimization problem (it is possible by the decomposition of this problem to independently solved easy subproblems). Let us note that our method can also be used for finding an upper bound for the multi-dimensional multiple-choice knapsack   Table 3 Obtained results for Experiment 1, problem set (unc, 1000, 10). problem (M M CKP ) via the relaxation obtained by summing up all the linear inequality constraints [1]. The method can be compared to greedy algorithm for multiple-choice knapsack problems which also finds, in general, an approximate solution and an upper bound.
Two preliminary computational experiments have been conducted to check how the proposed algorithm behaves for simple to solve (uncorrelated) instances and hard to solve (weakly correlated) instances. The results have