An exact algorithm for biobjective integer programming problems

# An exact algorithm for biobjective integer programming problems

Saliha Ferda Doğan111All authors contributed equally to this work. All are from Bilkent University, Department of Industrial Engineering, Ankara, 06800 Turkey ferda.dogan@bilkent.edu.tr    Özlem Karsu ozlemkarsu@bilkent.edu.tr    Firdevs Ulus firdevs@bilkent.com.tr
July 2, 2019
###### Abstract

We propose an exact algorithm for solving biobjective integer programming problems, which arise in various applications of operations research. The algorithm is based on solving Pascoletti-Serafini scalarizations to search specified regions (boxes) in the objective space and returns the set of nondominated points. We develop variants of the algorithm, where the choice of the scalarization model parameters differ; and demonstrate their performance through computational experiments both as exact algorithms and as solution approaches under time restriction. The results of our experiments show the satisfactory behaviour of our algorithm, especially with respect to the number of mixed integer programming problems solved compared to an existing approach. The experiments also demonstrate that different variants have advantages in different aspects: while some variants are quicker in finding the whole set of nondominated solutions, other variants return good-quality solutions in terms of representativeness when run under time restriction.

Keywords: Biobjective integer programming, Pascoletti-Serafini scalarization, Algorithms.

## 1 Introduction

In many operations research applications such as scheduling, task assignment and transportation, the underlying problem is an integer programming problem. Moreover, a vast amount of these problems require more than a single criterion to be considered, leading to biobjective (multiobjective) integer programming problems.

In this study, we focus on biobjective integer programming problems (BOIP) and propose an algorithm that returns the whole set of nondominated points of these problems. There are a number of solution approaches that have been designed for BOIP in the literature, most of which explore the objective (criterion) space by repetitively solving single objective optimization problems related to the BOIP, called scalarization problems (or simply, scalarizations). A scalarization is formulated by means of a real-valued scalarizing function of the objective functions of the BOIP, auxiliary scalar or vector variables and/or parameters ([7]).

There are several scalarizations proposed in the literature. The widely-used ones are the weighted sum scalarization ([9, 17]), the -constraint scalarization ([8]) and the (weighted) Chebyshev scalarization ([2, 15]). Most of the current algorithms in the literature solve these scalarizations or their modifications repetitively to find the set of nondominated solutions. Commonly used ones are the perpendicular search and the -constraint algorithm, which are based on weighted sum scalarization and -constraint scalarization, respectively ([4, 5, 8]). Examples of algorithms using weighted Chebychev scalarizations are proposed in [12] and [14], where a modified version of the scalarization is used. There are also two-phase algorithms, which generate supported nondominated points in the first phase and find the unsupported nondominated points by exploring the triangles defined by two consecutive supported nondominated points in the second phase ([11, 16]). Recently, the balanced box algorithm is proposed in [1] and a two-stage algorithm which combines the balanced box and -constraint algorithms is discussed in [6].

We propose an exact solution algorithm that finds the whole set of nondominated solutions to biobjective integer programming problems. The algorithm is based on Pascoletti-Serafini scalarization ([10]). We adapt this scalarization model for biobjective integer programming settings and develop different variants of the algorithm. We compare these variants with respect to number of (mixed) integer programming problems solved and solution time. We also test the performances of the variants under time limit and report on the representativeness of the obtained solution sets using the (scaled) coverage error ([3, 13]).

The structure of the paper is as follows. In Section 2 we give the preliminaries and the problem definition. In Sections 3 and 4 we explain the base algorithm and its variants, respectively. We test the performances of the algorithm and report the results of our experiments in Section 5. We conclude our discussion in Section 6.

## 2 Preliminaries and problem definition

A general bi-objective integer programming problem is formulated as

 min"{z(x)=(z1(x),z2(x))T|x∈X⊆Zn}, (P)

where , are integer-valued objective functions. The set represents the feasible set in the decision space and the set represents the feasible set in the objective space.

Throughout the paper we will use the following notation for vector inequalities:

 z(x′)≤z(x) :⟺zi(x′)≤zi(x)~{}for~{}i∈{1,2}; z(x′)⪇z(x) :⟺z(x′)≤z(x)~{}and~{}z(x′)≠z(x); z(x′)
###### Definition 2.1.

dominates (strictly dominates) if (). If there exists no such that dominates (strictly dominates) , then is nondominated (weakly nondominated) and is efficient (weakly efficient).

The set of all nondominated vectors is denoted by . The ideal and nadir points of problem () are as follows:

 s0:=(minx∈Xz1(x),minx∈Xz2(x))T,u0:=(maxz∈Nz1,maxz∈Nz2)T.

A lexicographic optimization problem with two objective functions is given by

 lexmin{zi(x),zj(x)|x∈X}, (1)

where and . Solving (1) means solving the following two (single objective) optimization problems: First, ; and given an optimal solution of the first model, . Solving a lexicographic optimization yields an efficient solution.

In general, scalarization models are solved in order to find (weakly) efficient solutions. Throughout, Pascoletti-Serafini scalarization, which employs two parameters, a reference point and a direction is employed. The model is as follows:

 min{ρ | x∈X, z(x)≤s+ρd,ρ∈R}. (2)
###### Lemma 2.2 ([10]).

If is an optimal solution of (2) for some and , then is weakly efficient.

###### Lemma 2.3.

If is an optimal solution of (2) for some and then and are equal in at least one component.

###### Proof.

Proof Assume to the contrary that and are different in both components. That is, . Hence there exists such that , which contradicts the optimality of . ∎

When a subset of is found through an algorithm or a procedure, in order to measure how well represents the set of all nondominated points (), it is possible to use the ‘coverage error’ that is introduced by [13]. Similar measures are used in the literature to measure representativeness, an example is the coverage gap measure used recently in [3]. Here we provide the definition of coverage error for the special case where Chebyshev metric is used. We also introduce the scaled version as in [3].

###### Definition 2.4.

Let be a representative subset. The coverage error of with respect to is

 CE(¯N,n):=min¯n∈¯N(max{|n1−¯n1|,|n2−¯n2|}),

the coverage error of is

 CE(¯N)=maxn∈NCE(¯N,n)

and the scaled coverage error of is

 SCE(¯N)=CE(¯N)max{u01−s01,u02−s02},

where and are the nadir and the ideal points, respectively.

## 3 The algorithm

Throughout the algorithm the search regions in the objective space are referred to as boxes. A box is defined by three points in the criterion space, namely the starting point , the nondominated point which defines the first component of the starting point and the nondominated point which defines the second component of the starting point, and denoted as follows

 b(s,p,t):={y∈R2∣s1≤y1≤p1, s2≤y2≤t2}.

Note that it is possible to define the box using only and . However, we keep the starting point in the definition as it is used in the scalarization models.

The general idea of the algorithm can be described as follows. At the beginning, two sets namely and , are defined to denote the set of nondominated points and boxes to be investigated, respectively. For initialization, two corner points of the nondominated set are found by solving and . Let the optimal objective function vectors of these models be , , respectively. We initialize as and as , where is the ideal point. Clearly, the initial box includes all nondominated points. See Figure 1 for the illustration of the initial region.

At each iteration, the algorithm searches one box from set to find a (weakly) nondominated point by solving a Pascoletti-Serafini scalarization. In order to ensure finding a nondominated point, an (two) extra model(s) is (are) solved and the obtained nondominated point(s) is (are) added to . Then, the explored box is discarded and if at least a new nondominated point is found, two new boxes are added to to be searched in the next iterations. The algorithm continues until there are no boxes to explore. The pseudocode of the algorithm is given by Algorithm 1.

At an arbitrary iteration, a box from set is selected and the following optimization problem is solved to search the box

 min{ρ|x∈X,z(x)≤sb+ρd,z1(x)≤ pb1−ϵ,z2(x)≤ tb2−ϵ}, (R(b,d))

where a is direction vector set to and . This is a slightly modified Pascoletti-Serafini model. The last two constraints are added to prevent finding the nondominated points and , which are already found in the previous iterations. If this problem is infeasible, then there is no nondominated point other than and in the box. Otherwise, let the optimal solution of () be and the corresponding (weakly) nondominated point be , see Lemma 2.2. Note that is the step size and defines the point , which has at least one common component with , see Lemma 2.3. Since the scalarization only guarantees that is weakly nondominated, the following problem(s) is (are) solved to ensure that a nondominated point is found. If the first components of and are equal () then,

 min{z2(x)|x∈X, z1(x)=z1(xb)} (P1(xb))

is solved and, if the second components are equal () then,

 min{z1(x)|x∈X, z2(x)=z2(xb)} (P2(xb))

is solved. Notice that it is possible to have and in this case, both problems are solved. Let the solutions of () and () be and , respectively and and be the corresponding points in the criterion space. If only () is solved, then is set to and symmetrically, if only () is solved is set to (to be used in partitioning) (see lines 8-17 in Algorithm 1). Then, is updated accordingly (lines 18-23).

If both () and () are solved, it is possible to find two nondominated points and in the same iteration. In this case, both and are added to . See Figures 4-4 for illustrations of these cases.

As for any (weakly) nondominated point , the dominated region and the dominating region () can not contain any nondominated points, the current box is split into two boxes using and . More specifically, the first box is formed as , where and the second box is formed as , where . See Figures 4-4 for the illustrations of newly formed boxes for different cases.

Finally, we avoid searching regions which can not have any new nondominated points, by taking the advantage of the integrality of the problem , and the structure of a box. The boxes which do not satisfy and are eliminated since they can not include any nondominated points other than and . After new boxes are defined and their sizes are checked to make sure that they can include nondominated points, they are added to set to be searched in the next iterations. Then, the searched box is removed from the set (lines 24-29). The algorithm repeats the steps which are introduced above until there is no box in .

The algorithm works correctly and returns the set of all nondominated points after finitely many iterations. These are shown by the following two propositions.

###### Proposition 3.1.

Algorithm 1 works correctly: It returns the set of all nondominated points.

###### Proof.

Proof The points that are added to set are guaranteed to be nondominated. Indeed, () is a Pascoletti-Serafini scalarization with box contraints and by Lemma 2.2, it returns a weakly efficient solution. By solving () and/or (), finding an efficient solution is guaranteed. Moreover, by the structure of defining the new boxes, it is guaranteed that the set of all boxes to be searched () include all the remaining (if any) nondominated points at any time through the algorithm. ∎

###### Proposition 3.2.

Algorithm 1 solves integer programs, where is the number of nondominated points, is the number of cases where , is the number of sub-cases that two nondominated points are found and is the number of eliminated boxes using the elimination rule.

###### Proof.

Proof The following expression, parts of which will be explained in detail, shows the number of models solved:

 (4)(a)+(1)(b)+(2C2)(c)+2(N−2−2C2)(d)+(N−2)(e)+(C−C2)(f)−E(g)

At the beginning of Algorithm 1, two lexicographical minimization problems are solved to find and (a) and one () problem is solved to search the initial box (b). points are found in number of cases ( and two solutions are found), each of these points lead to a new box, hence a new () model (c). For the rest of the nondominated points, (), each point results in two new boxes (and hence two () models to be solved) (d). As for the () models: -2 points are found by solving a single second stage model (either () or ()) (e). Moreover, when and only a single nondominated point is found (in number of cases), we solve an extra () or (), which does not yield a new point (f). Finally, boxes are eliminated, avoiding the () models that would otherwise have been solved (g). ∎

## 4 The variants

The algorithm described in Section 3 can be modified in different ways. First of all, in each iteration the direction parameter can be chosen according to the current box instead of fixing it to (Fixed). We consider two alternatives: to choose as the diagonal direction of the current box, (Changing) and to choose as the direction starting from towards the nadir point, (Nadir). Secondly, the new search regions added to in each iteration can be chosen differently. In addition to the base version that is described in Section 3, we consider employing in defining the new regions. Accordingly, we use instead of in order to define a new box if it yields a smaller region than the base version, see Figures 6 and 6. This is done by replacing lines 12 and 17 of Algorithm 1 with and , respectively. Notice that since the corners of the newly formed boxes are not necessarily integer valued, one also needs to change the elimination rule slightly. More specifically, for this variant the strict inequalities on lines 26 and 28 of Algorithm 1 are replaced by greater than or equal to signs.

The six variants can be seen in Table 1. The base version that is described in Section 3 is denoted by FDN.

## 5 Computational results

We examine the efficiency of the algorithms by solving knapsack and assignment problem instances which are used in [1]222The instances are available at http:// hdl.handle.net/1959.13/1036183. Both problem types contain four classes, A, B, C, D each with five instances. The first set consists of biobjective Knapsack Problem (KP) instances with 375, 500, 625 and 750 variables. The second set consists of biobjective Assignment Problem (AP) instances with 200 200 and binary variables.

The algorithms are coded in C++ and all mixed integer programming models are solved using CPLEX 12.6. Only a single thread is used. All of the instances are run in a computer with Intel Xeon CPU E5-1650 3.6 GHz processor and 32 GB RAM. Computation times are given in central processing unit (CPU) seconds.

We first conduct preliminary experiments on type A knapsack and assignment instances in order to compare the performances of the algorithm variants. In Table 2, we report the average values for the number of nondominated points (), the number of all (mixed) integer programming problems solved, the solution time (in CPU seconds), the number of () models solved, average time for solving one () model, average time for solving one () model, , and , see Proposition 3.2.

Overall, we see that partitioning a box using nondominated points is a better box defining strategy leading to smaller number of problems solved, hence smaller solution times, except the KP case with changing direction according to nadir. This good performance is mostly due to the increase in the number of boxes that are eliminated (E) with our elimination rule (see lines 26 and 28 of Algorithm 1).

We observe that the original algorithm (FDN) consistently performs good in terms of solution time over all test instances, being the fastest algorithm for KP and the second fastest for AP.

Based on these results, we conduct further preliminary experiments with FDN, CDN and NDN variants. Since finding the whole set of nondominated points might be computationally demanding for most biobjective integer programming problems, early termination performances of the algorithms are also worth considering. Therefore, we run these variants with predetermined time limits and report the quality of the set of nondominated points obtained. We calculate the coverage error measure, given by Definition 2.4.

Table 4 shows the performance results for the three algorithm variants when they are run with time limits for class A of KP and AP. The time limit is set as 300 and 700 seconds for KP and AP, respectively. This corresponds to approximately 30% of time required to find the whole set of nondominated points. The table shows the number of nondominated points found (), the coverage error (CE) and the scaled coverage error (SCE) values for each variant. It is seen that CDN significantly outperforms the other variants with respect to coverage. This result is expected as setting the direction as the diagonal vector of the box to be searched encourages the algorithm to find scattered solutions across the Pareto frontier and provides a highly representative set even at the early stages of the algorithm. In Figures 9-9, we provide the solution sets found when KP instances are solved with time limited versions of FDN, CDN and NDN, respectively.

In the first set of preliminary experiments, we eliminated the variants defining the box using , and concluded that FDN, CDN, NDN are worth further consideration, FDN being the most computationally efficient one. In the second set of experiments with time limits we have observed that CDN is the top-performer. Based on these preliminary observations, we decided to perform the main experiments with FDN as it is computationally more promising, and CDN as it outperforms the other variants under time restriction.

Tables 5 and 6 show the results of our main experiments, in which we compare FDN and CDN over all instances of KP and AP. We report the average values for the number of nondominated solutions (), the number of models solved, the total solution times, the number of () models solved and average time spent to solve () and () models as well as , and . The results verify the observations made at the preliminary experiments: although FDN solves more integer programming problems in total, it solves less of the more difficult () models, hence it works faster than CDN. Moreover, when FDN is used, the number of cases where () is significantly larger than that of CDN. This is due to the nature of the direction vector used; moving along the search region in a fixed direction of , the algorithm visits integer corners more often compared to a diagonal direction vector. This increases the cases where () and () are both solved within a box. Note that, in only a small portion of these cases two new nondominated points are found, implying that one () is solved unnecessarily. However, since these models are much easier to solve compared to (), solving more of these does not significantly affect the computational performance of FDN.

We also provide a comparison with the balanced box algorithm ([1]). Since the algorithms are coded and run on different platforms, we cannot compare the solution times with those of the balanced box algorithm. Note that the balanced box algorithm solves exactly problems. Therefore, it will solve more models for all of the problem instances considered; indeed it solves 25.5%, 36.5% more problems than our best algorithm variant on average for KP and AP, respectively.

CDN works slower compared to FDN but our preliminary experiments show that it is promising when used with time limits. We verified this observation by performing experiments for the whole KP and AP sets with time limit, the results of which are provided in Table 7.

Overall, one can conclude that both variants are powerful in different aspects. When used to find the complete set of nondominated points, FDN works better since it runs faster. However, CDN is a very promising variant when run with a time limit since it quickly provides a highly representative subset of solutions.

### 5.1 An extension of CDN

When we examine the results of average time spend for an () model, we observe that there is significant difference between FDN and CDN for class C of AP, see Table 6. In these instances, average time spent per () model in CDN is more than twice of the time spent in FDN. To investigate this further, we check the solution times of each individual () model solved in CDN for these instances. We see that the majority of the total time is occupied by only few models. To overcome this issue of extreme solution times, we modify CDN and solve each () model under a time limit. If the model is aborted due to the time limit, we slightly modify the direction and solve the model with the new direction parameter. That is, we change line 4 of Algorithm 1 as follows: \@floatalgocf[htbp]     \end@float

We refer to this extension of CDN with time limited () models as TL-CDN. We compare the performance of TL-CDN (where a time limit of 50 seconds is used for each () model) with those of FDN and CDN in class C of AP. The results are presented in Table 8. When we compare the number of integer programming problems solved by the algorithms, we observe that CDN is the best algorithm and it is closely followed by TL-CDN, as expected. When we analyse the run times and average () solution times of TL-CDN and CDN, we observe that there is a significant improvement when TL-CDN is used, indicating that the extension is successful.

We also run TL-CDN with predetermined time limits for class C of AP and report the quality of the solution set (using coverage error) in Table 9, by comparing it with FDN and CDN. It is seen that TL-CDN outperforms CDN in terms of coverage error in most instances. The results show that this modification (TL-CDN) is successful in significantly reducing run time without sacrificing from performance in representativeness.

## 6 Conclusion

We propose an exact solution approach for biobjective integer programming problems based on solving Pascoletti-Serafini scalarizations to search for nondominated points within boxes in the objective space. We generate different variants of the algorithm based on how the boxes are defined and how the direction vector in the scalarization problem is set.

We compare the performances of the algorithm variants both with and without time limits and determine the ones that outperform the others. Our results indicate that the variants using nondominated points to define the boxes are better. Moreover, although the variant using a fixed direction vector of leads to more (mixed) integer programming problems solved, it requires less computational time since it solves less of the more difficult scalarization models. We, however observe that the variant setting direction with respect to the diagonal of the box to be searched is still promising since it returns a highly representative subset (measured using coverage error) of the set of nondominated points when it is run with a time limit. We suggest an extension to this variant, which has lower solution times and better coverage error results.

We prove that the algorithm terminates and show through computational experiments that the best variant solves less problems than a recent algorithm, the balanced box method.

## References

• [1] Natashia Boland, Hadi Charkhgard, and Martin Savelsbergh. A criterion space search algorithm for biobjective integer programming: The balanced box method. INFORMS Journal on Computing, 27(4):735–754, 2015.
• [2] V. Joseph Bowman. On the relationship of the Tchebycheff norm and the efficient frontier of multiple-criteria objectives. In Hervé Thiriez and Stanley Zionts, editors, Multiple Criteria Decision Making, pages 76–86, Berlin, Heidelberg, 1976. Springer Berlin Heidelberg.
• [3] GÃ¶khan Ceyhan, Murat KÃ¶ksalan, and Banu Lokman. Finding a representative nondominated set for multi-objective mixed integer programs. European Journal of Operational Research, 272(1):61 – 77, 2019.
• [4] L.G. Chalmet, L. Lemonidis, and D.J. Elzinga. An algorithm for the bi-criterion integer programming problem. European Journal of Operational Research, 25(2):292 – 300, 1986.
• [5] V. Chankong and Y.Y. Haimes. Multiobjective Decision Making Theory and Methodology. Elsevier Science, New York, 1983.
• [6] Rui Dai and Hadi Charkhgard. A two-stage approach for bi-objective integer linear programming. Operations Research Letters, 46(1):81 – 87, 2018.
• [7] Matthias Ehrgott. Multicriteria Optimization (2. ed.). Springer, 2005.
• [8] Y Haimes, L Lasdon, and D Wismer. On a bicriterion formulation of the problems of integrated system identification and system optimization. IEEE Transactions on Systems, Man, and Cybernetics, SMC-1(3):296–297, 1971.
• [9] Juhani Koski. Multicriteria Truss Optimization, pages 263–307. Springer US, Boston, MA, 1988.
• [10] A. Pascoletti and P. Serafini. Scalarizing vector optimization problems. Journal of Optimization Theory and Applications, 42(4):499–524, Apr 1984.
• [11] Anthony Przybylski, Xavier Gandibleux, and Matthias Ehrgott. Two phase algorithms for the bi-objective assignment problem. European Journal of Operational Research, 185(2):509 – 533, 2008.
• [12] Ted Ralphs, Matthew Saltzman, and Margaret Wiecek. An improved algorithm for solving biobjective integer programs. Annals OR, 147:43–70, 10 2006.
• [13] Serpil Sayın. Measuring the quality of discrete representations of efficient sets in multiple objective mathematical programming. Mathematical Programming, 87(3):543–560, May 2000.
• [14] Serpil Sayın and Panos Kouvelis. The multiobjective discrete optimization problem: A weighted min-max two-stage optimization approach and a bicriteria algorithm. Management Science, 51(10):1572 – 1581, 2005.
• [15] Ralph E. Steuer and Eng-Ung Choo. An interactive weighted Tchebycheff procedure for multiple objective programming. Mathematical Programming, 26:326–344, 10 1983.
• [16] E.L. Ulungu and J. Teghem. The two phases method: An efficient procedure to solve bi-objective combinatorial optimization problems. Foundations of Computing and Decision Sciences, 20(2):149–165, 1995.
• [17] L. Zadeh. Optimality and non-scalar-valued performance criteria. IEEE Transactions on Automatic Control, 8(1):59–60, 1963.
You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters