On the rectangular knapsack problem: approximation of a specific quadratic knapsack problem

In this article, we introduce the rectangular knapsack problem as a special case of the quadratic knapsack problem consisting in the maximization of the product of two separate knapsack profits subject to a cardinality constraint. We propose a polynomial time algorithm for this problem that provides a constant approximation ratio of 4.5. Our experimental results on a large number of artificially generated problem instances show that the average ratio is far from theoretical guarantee. In addition, we suggest refined versions of this approximation algorithm with the same time complexity and approximation ratio that lead to even better experimental results.


Introduction
The classical (linear) knapsack problem (KP) is a combinatorial optimization problem. Given a finite set {1, . . . , n} of items i with assigned profit and weight values p i and w i , respectively, and a finite capacity W , KP decides which items to include to maximize the total profit while satisfying a capacity constraint. The capacity and profit and weight values are all assumed to be positive integer and each item can be included at most once. KP is an N P-complete problem but is solvable in pseudo-polynomial time by dynamic programming. Several interesting and challenging variants of KP have been introduced in the past (see (Kellerer et al. 2004) for an overview). One of those is the quadratic knapsack problem (QKP).
In contrast to KP, the profit of a collection of items in QKP is not only determined by the sum of individual profits, but also by profits generated by pairwise combinations of items. This can be used to model the fact that two items may complement each other such that their profit is increased if both of them are selected. The quadratic objective still allows to model that the profit of two items is independent of each other by setting the combined profit to 0. In this case, including both items does not increase the profit over the sum of the individual profits. Furthermore, a negative combined profit value can model the fact that both items together are less profitable than the sum of individual profits. This might be the case, if both items are substitutes for each other and including both items is as profitable as including one.
The formulation of quadratic knapsack problems (QKP) is very general and, therefore, its range of applications is quite wide. For example, Johnson et al. (1993) present a problem in the context of compiler construction that can be formulated as a QKP. Moreover, QKPs have been discussed in the context of the location of airports, freight handling terminals, railway stations, and satellite stations (Rhys 1970;Witzgall 1975).
We present a variant of QKP which we call the rectangular knapsack problem (RKP). The profit matrix is built by the product of two vectors and the constraint is a cardinality constraint.
The main motivation for problem RKP arises when solving the cardinality constrained bi-objective knapsack problem (2oKP) where a, b ∈ N n , with a, b = 0 n = (0, 0, . . . , 0) ∈ N n , and k ∈ N, k < n. Instead of computing the set of efficient solutions for this bi-objective optimization problem, we want to find one (or several) representative nondominated point(s). Originally proposed by Zitzler and Thiele (1998) in the context of evolutionary algorithms, the hypervolume indicator is often used as a versatile quality measure of representation of the efficient set in multiobjetive optimization (c.f. Kuhn et al. 2016). The problem of finding one solution of 2oKP that maximizes the hypervolume, considering (0, 0) as reference point, is equivalent to RKP. The structure of RKP allows to formulate a polynomial time 4.5-approximation algorithm. For maximization problems, an algorithm is called a polynomial time ρ-approximation algorithm, if it computes a feasible solution in run time being polynomial in the encoding length of the input such that ρ ≥ OPT ALG .
Here, OPT denotes the optimal objective function value of the maximization problem and ALG the objective function value of the solution which is the output of the algorithm (Cormen et al. 2001). The remainder of this article is organized as follows: In Sect. 2, we give an introduction to quadratic knapsack problems. We introduce the rectangular knapsack problem in Sect. 3 and present upper and lower bounds. These bounds motivate an approximation algorithm that is formulated in Sect. 4, for which a constant approximation ratio ρ is proven. Furthermore, we also introduce improved implementations of this approximation method. In Sect. 5 we present a computational study of these algorithms and compare the realized approximation ratios to the theoretical bound of 4.5. Section 6 concludes this article.
2 Quadratic knapsack problems Gallo et al. (1980) first introduced the binary quadratic knapsack problem (QKP). It is a variant of the classical knapsack problem and can be concisely stated as follows: Given n items, the profit for including item i is given by the coefficient p ii . Additionally, a profit p i j + p ji is generated if both items i and j are selected. The values p i j (which are often assumed to be non-negative integers) can be compactly written in a profit matrix The profits p i j and p ji are either both realized, i.e., if items i and j are selected, or both not realized, i.e., if item i or item j is not selected. Hence, p i j and p ji can be assumed to be equally valued, which results in a symmetric matrix P.
As for the classical knapsack problem, each item i has a positive integral weight w i and the goal is to select a subset of items that maximizes the overall profit while the sum of weights does not exceed the given capacity W . As usual, the binary decision variable x i indicates if item i is selected, x i = 1, or not, x i = 0. Thus, QKP can be defined as follows: It is well-known that the quadratic knapsack problem (with knapsack constraint but also with cardinality constraint) is N P-complete in the strong sense, which can be shown by a polynomial reduction from the Clique-problem (Garey and Johnson 1979;Pisinger 2007).
The quadratic knapsack problem has been widely studied in the literature, see Pisinger (2007) for a comprehensive survey. Exact solution algorithms are mainly based on branch-and-bound (B&B) schemes. Besides the model of QKP, Gallo et al. (1980) also presented the first B&B algorithm for this optimization problem. Since then, several techniques have been presented to compute good upper bounds for B&B algorithms [e.g., Calmels 1996, Rodrigues et al. (2012): linearization of the quadratic problem; Caprara et al. (1999): upper planes; Billionnet et al. (1999): Lagrangian decomposition;Helmberg et al. (2000): semidefinite relaxation techniques] and to fix decision variables at their optimal value before applying the final optimization (Billionnet and Soutif 2004;Pisinger et al. 2007).
However, few results are known about the approximation of QKP. Since the problem is strongly N P-hard, a fully polynomial time approximation scheme (FPTAS) cannot be expected unless P = N P. Furthermore, it is unknown whether there exists an approximation with a constant approximation ratio for QKP. Taylor (2016) present an approximation algorithm based on an approach for the densest k-subgraph problem. They show that for ε > 0, QKP can be approximated with an approximation ratio in O(n 2 5 +ε ) and a run time of O(n 9 ε ). Rader and Woeginger (2002) prove that for a variant of QKP, where positive as well as negative profit coefficients p i j are considered, there does not exist any polynomial time approximation algorithm with finite worst case guarantee unless P = N P. Other approximation results concentrate on special cases of QKP (Pferschy and Schauer 2016): FPTAS for QKP on graphs of bounded tree width and polynomial time approximation scheme (PTAS) for graphs that do not contain any fixed graph H as a minor; Kellerer and Strusevich (2010) and Xu (2012): FPTAS for the symmetric quadratic knapsack problem; Pferschy and Schauer (2016): QKP on 3-book embeddable graphs is strongly N P-hard; Rader and Woeginger (2002): QKP on vertex series-parallel graphs is strongly N P-hard).

Rectangular knapsack problems
The (cardinality constrained) rectangular knapsack problem (RKP) is a variant of QKP which can be written as follows: where a, b ∈ N n , with a, b = 0 n , and k ∈ N, k < n. Note, that P = a b , i.e., rank(P) = 1, with p i j = a i b j and p ji = a j b i , i.e., in general, P is not symmetric. We assume that k ≥ 2. Otherwise, i.e., if k = 1, the problem reduces to finding the largest coefficient a i b i , for i ∈ {1, . . . , n}.
The rectangular objective function is formulated in analogy to the so-called Koopmans-Beckmann form of the quadratic assignment problem, see Burkard et al. (1998), which is also a particular case of the more general Lawler formulation. In both cases, the two respective four dimensional arrays of profit/cost coefficients are given as a product of two lower dimensional parameter matrices or vectors, respectively.
The complexity of RKP is an open question. However, we can prove that two different simple extensions of this problem are NP-hard: Again, the answer to Partition is "yes" if and only if the optimal solution has value , i.e., if the bound is tight. The answer is "no" in two cases: Either if there exists a partition but with |I | = n 2 . In this case inequality ( * ) is strict. Or if there does not exist a partition. In this case inequality ( * * ) is strict.

Illustrative interpretation
The denotation rectangular knapsack problem is motivated by the special structure of P given by the coefficients a i b j . Each coefficient can be interpreted as the area of a rectangle. Accordingly, for fixed itemî ∈ {1, . . . , n}, all rectangles corresponding to coefficients aˆı b j , j = 1, . . . , n, have the same width, and all rectangles corresponding to coefficients a j bˆı , j = 1, . . . , n, have the same height. Note that the objective function can be rewritten as which can be interpreted as choosing a subset S ⊂ {1, . . . , n} of items such that the area of the rectangle with width i∈S a i and height i∈S b i is maximized.

Example 1
We consider the following instance of RKP: max (7, 12, 2, 5, 4) x · (6, 3, 8, 5, 10) x The corresponding rectangles are plotted in Fig. 1. Each rectangle has the same position in the overall rectangle as the corresponding coefficient p i j = a i b j in the profit matrix P. The optimal solution x = (0, 1, 0, 0, 1) generates an objective function value that corresponds to the highlighted area in the figure.

Bounds
The structure of the profit matrix P implies an easy computation of bounds for RKP. In the following, we assume that all instances are defined or reordered such that Let S n denote the symmetric group of order n and π ∈ S n denote a permutation of {1, . . . , n}. More specifically, consider π such that Using the sorted coefficients a i , b π( j) of the objective function, one can compute an upper bound for RKP in a straightforward way.
Lemma 1 For every feasible solution x ∈ {0, 1} n of RKP, the following inequality holds: Note that, in general, this upper bound does not correspond to a solution of RKP since the value of a variable x i may be differently defined w. r. t. the respective sorting of the coefficients. As soon as Eq. (1) holds, the upper bound U corresponds to a feasible solution of RKP and this solution is optimal.
Proof We consider the objective function of RKP: The cardinality constraint restricts the number of selected items to k. Due to the ordering of coefficients a i , it is for every feasible solution x of RKP. Analogously, due to the definition of the permutation π , we know that for every feasible solution x of RKP. Thus, The solution x is feasible and realizes U. Hence, x is optimal and U is a tight upper bound.
A lower bound on RKP can also be obtained by using the sorting of the coefficients. Letx andx ∈ {0, 1} n be defined as follows: For notational convenience, let κ:= k 2 , κ:= k 2 , and κ:= k 2 . If k is even, the equality κ = κ = k 2 holds, i.e.,x andx are identical. Remark 1 The definition ofx guarantees that at least the product is realized in the objective function. Due to the ordering of the coefficients a i and b π( j) , this is the maximal possible value that a product of κ coefficients a i and κ coefficients b j can achieve. The same holds analogously forx. This property is important to prove an approximation quality in the following, see the proof of Theorem 1.
We define L:= f (x) and L:= f (x). The following example shows a connection between the bound computation and the visualization of RKP as a selection of a subset of rectangular areas.
In this context, we show that the following inequality holds: Referring to the description of Example 2, the right-hand side of this inequality corresponds to the area of the κ · κ largest rectangles in the left upper part of the overall rectangle (see also Remark 1). We partition the area corresponding to the lower bound into four distinct areas to show that the inequality holds: In the context of Example 2, the four terms resulting from this partition correspond, in this order, to the four areas (I to IV) in Fig. 3.
Analogously, using the definition ofx, it holds that:

Approximation algorithms
The results of Sect. 3.2 naturally motivate an approximation algorithm, see Algorithm 1. It computes the solutionsx andx and outputs the better alternative as an approximate solution.
The computation ofx andx and of their objective function values L and L can be realized in time O(n). Therefore, with a time complexity of O(n log n), the sorting of the coefficients determines the time complexity of Algorithm 1.
Theorem 1 Algorithm 1 is a polynomial time 4.5-approximation algorithm for the rectangular knapsack problem.
Proof Algorithm 1 returns a feasible solution in polynomial time O(n log n).
In summary, this yields the approximation factor ρ: As presented in the proof of Theorem 1, we can guarantee better results for even values of k. Also, if k is odd the quality of the approximation increases for increasing values of k. However, in the worst case the approximation ratio is tight as is shown in the following example.

Example 3
Consider an instance of RKP with n ≥ 3k, M ∈ R, and coefficients a 1 = . . . = a k = M a k+1 = . . . = a n−k = M − 1 a n−k+1 = . . . = a n = 1 Algorithm 1 computes a lower bound solution with L even = (κ · M + κ · 1) 2 = k 2 4 (M + 1) 2 for even values of k and for odd values of k, respectively. As one can easily see, one optimal solution is given by Thus, for increasing values of M the approximation ratio tends towards for even values of k and for odd values of k ≥ 3, respectively. Note that, for fixed values of k, ρ odd exactly matches the approximation ratios given in Remark 3.

Improvements of the approximation algorithm
In  (1), . . . , π(κ)} are not disjoint, i.e., if n i=1x i < k. Hence, it is possible to increase the lower bound value by including further items. Algorithm 2 demonstrates a possible procedure to compute an improved lower bound L impr that takes this into account.
An additional parameter k , which we call adaptive capacity, is introduced to increase the sets {1, . . . , κ} and {π(1), . . . , π(κ)}, and, therefore, increase the number of selected items, without violating the constraint. In the beginning, k is set to k. After computing the lower bound solutionx as defined in (2), the algorithm tests whether k items are selected or not. In the latter case, the adaptive capacity k is increased by the difference k − n i=1x i . A re-computation ofx, using k as capacity, allows to include more items in accordance with the ordering of the respective coefficients a i or b π(i) which compensates for the fact that the original sets are not disjoint. Subsequently, it is tested again if the constraint is satisfied with equality. If not, the adaptive capacity k is further increased. Otherwise, the algorithm continues by computingx using the current value of the parameter k as capacity and testing which of the lower bound values is larger.

Lemma 3
If the solutionx allows to increase the adaptive capacity k to k + (k − n i=1x i ) (in Step 10 of Algorithm 2), then this increase is also feasible for the computation ofx. Proof For ease of notation, we assume that we are examining the iteration where the adaptive capacity k is increased for the first time from the capacity k to k = k + (k − n i=1x i ). The following discussion can be applied in an analogous manner to all further iterations by adapting the notation accordingly.
Lemma 4 Let n be the number of items and let k be the capacity of RKP. Algorithm 2 terminates, has a worst case time complexity of O(n log n), and a worst case approximation ratio of 4.5.

Proof
The while-loop for computing the solutionx with n i=1x i = k is critical for the termination of Algorithm 2. In the first iteration, at least κ variables are set to 1. The parameter k is increased by at least 1 in each consecutive iteration and, thus, in at least every second iteration an additional entry ofx is set to 1. Hence, after at most 2 · κ + 1 = k iterations k variables have been selected forx and the loop terminates.
We take advantage of the ordering of the coefficients to set only new variables to 1 if the adaptive capacity is increased. Thus, the execution of the while loop requires O(k). The complexity of Algorithm 2 is determined by the sorting algorithm (cf. Algorithm 1), i.e., Algorithm 2 has a worst case time complexity of O(n log n).
The approximation ratio is at most 4.5, since the heuristic solution of Algorithm 1 gives a lower bound on the heuristic solution of Algorithm 2. In the worst case, the approximation ratio is tight, since Algorithm 2 computes the same heuristic solution for the RKP instance of Example 3 as Algorithm 1.

Example 4
We apply the improved approximation algorithm, Algorithm 2, on the instance of RKP of Example 2. The solutionx is defined by the set with |Ĩ| = n i=1x i = 4 < 5. Thus, the adaptive capacity can be set to k :=5 + (5 − 4) = 6. The re-computation ofx leads tõ with |Ĩ| = n i=1x i = 5. Hence, the cardinality constraint is tight and, since k is even, the solution x = (1, 1, 1, 1, 0, 1, 0, 0) generates the improved lower bound L = f (x) = 920. The optimal solution of the instance x * is identical to the improved lower bound solution x * = x.
As proven above, setting the adaptive capacity to k :=6 is also feasible forx. The solutionx defined by k = 5 corresponds to the set with |Î| = n i=1x i = 4 < 5. Hence, one additional item can be included, resulting again in the same lower bound solution x = (1, 1, 1, 1, 0, 1, 0, 0) (with k = 6). The areas of rectangles corresponding to the lower bounds L and L based on the first computations ofx andx (Algorithm 1), respectively, and the improved lower bound L (Algorithm 2) are shown in Fig. 4.
A second improvement of Algorithm 1 is motivated differently: Without an analysis of the input data, the distribution of the entries of the coefficient vectors a and b is unknown, and, thus, there might be better selections then deciding equally according to both of the orderings. One possible approach is to compute various alternative solutions, still based on the sorting of the coefficients, and select the best solution.The alternatives can be defined by setting the variables only corresponding to the sorting of a, to the sorting of b, and by all alternatives in between, i.e., x 1 = . . . = x j = 1, x π(1) = . . . = x π(k− j) = 1 and x i = 0 for all remaining indices, for 0 ≤ j ≤ k. We call this approach shifted selection.
The combination of both improvements, shifted selection and adaptive capacity, is formalized in Algorithm 3. If the adaptive capacity is updated, further items are included according to the ordering of the coefficients b π(i) . Other strategies are possible: include further items according to the ordering of the coefficients a i , alternate between both orderings or include arbitrarily chosen items. The quality of those strategies strongly depends on the given problem instance.
Algorithm 3 Improved approximation algorithm for RKP with shifted selection wrt. a and b and adaptive capacity Input: coefficients a = (a 1 , . . . , a n ) sorted in non-increasing order, b = (b 1 , . . . , b n ) , capacity k 1: x:=0 n and L:=0 2: compute permutation π ∈ S n such that In practice, it is quite intuitive to assume that Algorithm 3 leads to better approximation results than that of the basic version of Algorithm 1. However, the theoretical approximation ratio is the same. The approximation ratio is at most 4.5, since the solutionsx andx of Algorithm 1 are included in the set of alternatives. In the worst case, the approximation ratio is tight, since Algorithm 3 computes the same heuristic solution for the RKP instance of Example 3 as Algorithm 1.

Computational experiments
In this section, the quality of all presented variants of the approximation algorithm is evaluated experimentally on a wide range of RKP instances. We implemented the basic variant of the approximation algorithm (Algorithm 1), the improved variants with adaptive capacity (Algorithm 2), the shifted selection, and the combined version of these two improvements (Algorithm 3); we will refer to the solution quality returned by these variants as L basic , L impr , L shift , L comb , respectively. All algorithm variants were implemented in C. The QKP solver by Caprara et al. (1999) was used to compute the optimal solutions of RKP (see Pisinger 2016). The computational experiments were performed on an Intel Quadcore 3.2 GHz with 4 GB RAM running Linux compiled with gcc 4.8.
Three different classes of instances were generated to test the algorithms: For each type of instances, four different constraint slacknesses c k , with k = c k ·n , were chosen: c k = 0.1, c k = 0.25, c k = 0.5, and c k = 0.75. The instance sizes were n = 100, 200, 300, 400, except for the negatively correlated instances, where problems with n = 25, 50, 75 were generated, additionaly. For the latter instance class, the QKP solver was not able to solve instances with n ≥ 75 and k ≥ 14 within one hour of CPU-time. For each combination of instance class, size and constraint slackness, 10 instances were generated. Noteworthy, all approximation algorithms required at most 0.01 seconds for all instances tested.
Tables 1, 2, 3 and 4 present the average results obtained for the three classes of instances where columns z * /L • refer to the average approximation ratios obtained by the four algorithm variants and column U/L * gives an upper bound on the approximation ratio, with L * = min{L basic , L impr , L shift }. In general, the results indicate that the approximation quality of all algorithm variants is much better than the guaranteed approximation ratio of 4.5 and that the improved versions yield even better results than the basic variant, except on negatively correlated instances, for which all versions presented a similar performance. Moreover, the instance size does not play a strong role on the approximation ratio. However, the performance of the four variants seem to be affected by the instance type. In the following, we discuss the results in more detail for each instance type. Uncorrelated instances The experimental results in Table 1 suggest that the improved variant performs better as the constraint slackness increases and that a larger capacity value k improves the approximation (see Remark 3). Differently, the basic variant does not seem to be affected by c k and presents the worst approximation ratio in all cases. The shifted and combined variants present the best approximation ratio for small c k . Both variants also improve the approximation for larger c k but not as much as for the improved variant, which gives the best approximation ratio.
Positively correlated instances Table 2 shows that the basic variant has the worst approximation ratio, although still far from the theoretical bound. For this variant, many items seem to be selected twice due to both orderings, which can be expected for positive correlated instances since the orderings of the coefficients to both objective functions should be rather similar. Noteworthy, all the three improved variants are close to an approximation ratio of 1.0. We observed that as soon as the equality in the constraint is ensured, all improved variants solve the problem to optimality.
Negatively correlated instances For this instance type, all variants present a similar approximation ratio; see Tables 3 and 4. In fact, the basic variant generates very good approximation results with approximation ratios close to 1.0 for the tested instances (see Table 3). This is especially good, since the exact algorithm could not solve larger instances in less then one hour of CPU-time whereas the approximation algorithm can compute a high quality approximation in less then 0.01 seconds. For these larger instances the upper bounds on the approximation ratio U/L • show a similar behavior as for the small instances. For small c k the bound takes values around 3.5 and improves for larger values of c k up to around 1.5.

Conclusion
We presented a geometric interpretation of the rectangular knapsack problem. Upper and lower bounds for the problem can be computed directly by sorting the coefficients of the objective function. Based on these bound computations, we introduced a polynomial time approximation algorithm for RKP that provides an approximation ratio of 4.5. In practice, however, the algorithm can be further improved by selecting additional items if the cardinality constraint is not met with equality. Furthermore, the selection strategy for items can be modified We tested all algorithm variants on knapsack instances with three different correlation structures, up to 400 items, and four different constraint slacknesses. The approximations were computed in 0.01 seconds or less per instance. We observed that in practice the approximation ratios of all algorithms are much better than the theoretical ratio of 4.5. Thus, our approximation algorithms are an efficient tool to compute approximations of good quality for RKP.
In the future it would be interesting to integrate the bound computations in a branchand-bound procedure to formulate an exact algorithm for RKP. Furthermore, the results seem to be transferable to higher dimensions, where we think of problems of the form The bound computations and algorithm formulations should be convertible without problems, whereas the proof of an approximation ratio may become more complicated due to more possible cases that may occur. We also suggested a field of application for RKP. Finding a representative solution of the bi-objective cardinality constrained knapsack problem that maximizes the hypervolume with the origin as reference point is modeled by the rectangular knapsack problem. It is, therefore, very interesting for future research.