A Gröbner bases methodology for solving multiobjective polynomial integer programs

0.23 MB·

by  Víctor Blanco

A GRO¨BNER BASES METHODOLOGY FOR SOLVING MULTIOBJECTIVE POLYNOMIAL
INTEGER PROGRAMS
V.BLANCOANDJ.PUERTO
Abstract. Multiobjectivediscreteprogrammingisawell-knownfamilyofoptimizationproblemswithalargespectrum
ofapplications. Thelinearcasehasbeentackledbymanyauthorsduringthelastyears. However,thepolynomialcase
has not been deeply studied due to its theoretical and computational difficulties. This paper presents an algebraic
approach for solving these problems. We propose a methodology based on transforming the polynomial optimization
problem in the problem of solving one or more systems of polynomial equations and we use certain Gro¨bner bases to
solve these systems. Different transformations give different methodologies that are analyzed and compared from a
9 theoretical pointofviewandbysomecomputational experiments viathealgorithmsthattheyinduce.
0
0
2

b
e
F 1. Introduction

8 A multiobjective polynomial program consists of a finite set of polynomial objective functions and a finite set of

polynomialconstraints(in inequality orequationform), andsolvingthatproblemmeans obtaining the setofminimal
]
C elements in the feasible region defined by the constraints with respect to the partial order induced by the objective
O functions.
Polynomialprogramshave a wide spectrum of applications. Examples of them are capital budgeting [17], capacity
.
h planning [6], optimization problems in graph theory [3], portfolio selection models with discrete features [2, 16] or
t
a chemical engineering [20], among many others. The reader is referred to [18] for further applications.
m Polynomialprogramminggeneralizeslinearandquadraticprogrammingandcanserveasatooltomodelengineering
[ applications that are expressedby polynomialequations. Eventhose problems with transcendentalterms such as sin,

log,andradicalscanbe reformulatedbymeansofTaylorseriesasapolynomialprogram. Asurveyofthe publications
1
on general nonlinear integer programming can be found in [9].
v
4 We study here multiobjective polynomial integer programs (MOPIP). Thus, we assume that the feasible vectors
0 have integer components and that there are more than one objective function to be optimized. This change makes
3 single-objective and multiobjective problems to be treated in a totally different manner, since the concept of solution
1 is not the same.
.
2 In this paper, we introduce a new methodology for solving general MOPIP based on the construction of reduced
0 Gr¨obnerbasesofcertainidealsrelatedtotheproblemandonsolvingtriangularsystemsofpolynomialequationsgiven
9
by those bases. Gr¨obner bases were introduced by Bruno Buchberger in 1965 in his PhD Thesis [7]. He named it
0
Gr¨obnerbasis payingtribute to his advisorWolfgang Gr¨obner. This theoryemergedas a generalization,fromthe one
:
v variable case to the multivariate polynomial case, of the Euclidean algorithm,Gaussian elimination and the Sylvester
i
X resultant. Oneofthe outcomesofGr¨obnerBasesTheorywasits applicationto linearintegerprogramming[8,14,22].
Later, Blanco and Puerto [5] introduces a new notion of partial Gr¨obner basis for toric ideals in order to solve
r
a multiobjective linear integer programs. A different approach for solving linear integer programs was developed by
Bertsimas et al. [4] based on the application of Gr¨obner bases for solving system of polynomial equations. This
alternative use of Gr¨obner bases is also used in the paper by Ha¨gglof et al. [13] for solving continuous polynomial
optimization problems. Further details about Gr¨obner bases can be found in [10, 11].
We describedifferentapproachesforsolvingMOPIPusingGro¨bnerbaseswhicharebasedonreducingthe problem
toseveraloptimalityconditions: thenecessaryKarush-Kuhn-Tucker,theFritz-JohnandthemultiobjectiveFritz-John
optimality conditions.
The paper is structured as follows. In the next section we give some preliminaries in multiobjective polynomial
integer optimization. We present in Section 3 our first algorithm for solving MOPIP using only the triangularization
property of lexicographic Gr¨obner bases. Section 4 is devoted to two different algorithms for solving MOPIP using a
ChebyshevlikescalarizationandtheKarush-Kuhn-TuckerortheFritz-Johnoptimalityconditions. Thelastalgorithm,
basedonthemultiobjectiveFritz-Johnoptimalitycondition,isdescribedinSection5. Finally,inSection6,wecompare
the algorithms with the results of some computational experiments and its analysis.
Date:December10,2008.
ThisresearchwaspartiallysupportedbySpanishresearchgrantnumbersMTM2007-67433-C02-01andP06-FQM-01366.
1
2. The multiobjective integer polynomial problem
The goal of this paper is to the solve multiobjective polynomial integer programs (MOPIP):
min (f (x),...,f (x))
1 k
s.t. g (x) ≤0 j =1,...,m
(MOPIP ) j
f,g,h h (x) =0 r =1,...,s
r
x ∈Zn
+
with f ,...,f ,g ,...,g ,h ,...,h polynomials in R[x ,...,x ] and the constraints defining a bounded feasible
1 k 1 m 1 s 1 n
region. Therefore, from now on we deal with MOPIP and we denote f = (f ,...,f ), g = (g ,...,g ) and
f,g,h 1 k 1 m
h = (h ,...,h ). If the problem had no equality (resp. inequality) constraints, we would denote it by MOPIP
1 r f,g
(resp. MOPIP ), avoiding the nonexistent term.
f,h
However, MOPIP can be transformed to an equivalent multiobjective polynomial binary problem. Since the
f,g,h
feasible region {x∈Rn :g (x)≤0, h (x)=0, j =1,...,m,r =1,...,s} is assumedto be bounded, it can be always
+ j r
n
embedded in a hypercube [0,u ]n. Then, every component in x, x , has an additional, but redundant, constraint
i i
iY=1
x ≤ u . We write x in binary form, introducing new binary variables z with values in {0,1}, x = ⌊logui⌋2jz ,
i i i ij i j=0 ij
substituting every xi in MOPIPf,g,h we obtain an equivalent 0−1 problem. P
Then, from now on, without loss of generality, we restrict ourselves to multiobjective polynomial binary programs
(MOPBP) in the form:
min (f (x),...,f (x))
1 k
s.t. g (x) ≤0 j =1,...,m
(MOPBP ) j
f,g,h h (x) =0 r =1,...,s
r
x ∈{0,1}n
If the problem had no equality (resp. inequality) constraints, we would denote the problem by MOPBP (resp.
f,g
MOPBP ), avoiding the nonexistent term.
f,h
The number of solutions of the above problem is finite, since the decision space is finite. Thus, the number of
feasible solutions is, at most |{0,1}n|=2n.
It is clear that MOPBP is not a standard optimization problem since the objective function is a k-coordinate
f,g,h
vector, thus inducing a partial order among its feasible solutions. Hence, solving the above problem requires an
alternative concept of solution, namely the set of nondominated (or Pareto-optimal) points.
A feasible vector x ∈ Rn is said to be a nondominated (or Pareto optimal) solution of MOPIP if there is no
f,g
other feasible vector y such that
b
f (y)≤f (x) ∀j =1,...,k
j j
withatleastonestrictinequalityforsomej. Ifxisanondominatedsolution,thevectorf(x)=(f (x),...,f (x))∈Rk
b 1 k
is called efficient.
We say that a feasible solution, y, is dominated by a feasible solution, x, if f (x) ≤ f (y) for all i = 1,...,k and
i i
f(x)6=f(y). We denote by X the set of all nondominated solutions for MOPBP and by Y the image under the
E f,g,h E
objective functions of X , that is, Y ={f(x):x∈X }. Note that X is a subset of Rn (decision space) and Y is
E E E E E
a subset of Rk (space of objectives).
From the objective functions f =(f ,...,f ), we obtain a partial order on Zn as follows:
1 k
x≺ y :⇐⇒f(x)(cid:0)f(y) or x=y.
f
Note that since f :Rn →Rk, the above relation is not complete. Hence, there may exist incomparable vectors.
In the following sections we describe some algorithms for solving MOPIP using tools from algebraic geometry. In
particular, in each of these methods, we transform our problem in a certain system of polynomial equations, and we
use Gr¨obner bases to solve it.
3. Obtaining nondominated solutions solving systems of polynomial equations
In this section we present the first approachfor solving multiobjective polynomial integer programs using Gr¨obner
bases. Forthismethod,wetransformtheprograminasystemofpolynomialequationsthatencodesthesetoffeasible
solutionsandits objectivevalues. Solvingthatsysteminthe objectivevalues,andthen, selecting the minimalonesin
the partial componentwise order, allows us to obtain the associate feasible vectors, thus, the nondominated solutions.
ThroughthissectionwesolveMOPBP . Withoutlossofgenerality,wereducethegeneralproblemtotheproblem
f,h
without inequality constraints since we can transform inequality constraints to equality constrains as follows:
(1) g(x)≤0⇐⇒g(x)+z2 =0,z ∈R.
where the quadratic term, z2, assures the nonnegativity of the slack variable and then, less than or equal to type
inequality. Initially, we suppose that all the variables are binary. In Remark 3.1 we describe how to modify the
algorithm to incorporate the above slack variables.
This approach consists of transforming MOPBP to an equivalent problem such that the objective functions are
f,h
part of the constraints. For this transformation, we add k new variables, y ,...,y to the problem, encoding the
1 k
objective values for all feasible solutions. The modified problem is:
min (y ,...,y )
1 k
s.t. h (x) =0 r=1,...,s
r
(2) y −f (x) =0 j =1,...,k
j j
x (x −1) =0 i=1,...,n
i i
y ∈Rk x∈Rn
where integralityconstraintsare codified as quadraticconstraints so,MOPBP is a polynomialcontinuous problem.
f,h
The algorithm consists of, first, obtaining the set of feasible solutions of Problem (2) in the y variables; then,
selecting from that set those solutions that are minimal with respect to the componentwise order, obtaining the set
of efficient solutions of MOPBP . The feasible solutions in the x-variables associated to those efficient solutions
f,h
correspond with the nondominated solutions of MOPBP .
f,h
Then, first, we concentrate in describing a procedure for solving the system of polynomial equations that encodes
the feasible region of Problem (2), i.e. the solutions of
h (x) =0 for all r =1,...,s
r
(3) y −f (x) =0 for all j =1,...,k
j j
x (x −1) =0 for all i=1,...,n.
i i
For analyzing the system (3) we use Gr¨obner bases as a tool for solving systems of polynomial equations. Further
details can be found in the book by Sturmfels [21].
Thesetofsolutionsof(3)coincideswiththeaffinevarietyofthefollowingpolynomialidealinR[y ,...,y ,x ,...,x ]:
1 k 1 n
I =hh (x),...,h (x),y −f (x),...,y −f (x),x (x −1),...,x (x −1)i.
1 m 1 1 k k 1 1 n n
Note that I is a zero dimensional ideal since the number of solutions of the equations that define I is finite. Let V(I)
denote the affine variety of I. If we restrict I to the family of variables x (resp. y) the variety V(I ∩R[x ,...,x ])
1 n
(resp. V(I ∩R[y ,...,y ])) encodes the set of feasible solutions (resp. the set of possible objective values) for that
1 k
problem.
Applying the elimination property,the reduced Gr¨obnerbasis forI, G, with respectto the lexicographicalordering
with y ≺ ··· ≺ y ≺ x ≺ ··· ≺ x gives us a method for solving system (3) sequantially, i.e., solving in one
k 1 n 1
indeterminate at a time. Explicitly, the shape of G is:
1) G contains one polynomial in R[y ]: p (y )
k k k
2) G contains one or severalpolynomials in R[y ,y ] : p1 (y ,y ),...,pmk−1(y ,y ).
k−1 k k−1 k−1 k k−1 k−1 k
.
.
.
k+1) G contains one or severalpolynomials in R[x ,y ,...,y ] : q1(x ,y),...,qsn(x ,y).
n 1 k n n n n
.
.
.
k+n) G contains one or severalpolynomials in R[x ,y ,...,y ] : q1(x ,...,x ,y),...,qs1(x ,...,x ,y).
n 1 k 1 1 n n 1 n
Then, with this structure of G, we can solve, in a first step, the system in the y variables using those polynomial in
G that only involve this family of variables as follows: we first solve for y in p (y ) = 0, obtaining the solutions:
k k k
y1,y2,.... Then, for fixed yr, we find the common roots of p1 ,p2 ,... getting solutions y1 ,y2 ,... and so
k k k k−1 k−1 k−1,r k−1,r
on, until we have obtained the roots for p (y ,...,y ). Note that at each step we only solve one-variable polynomial
1 1 k
equations.
We denote by Ω the above set of solutions in vector form
Ω={(yˆ ,...,yˆ ):p (yˆ )=0,p1 (yˆ ,yˆ )=0,...,pmk−1(yˆ ,yˆ )=0,...
1 k k k k−1 k−1 k k−1 k−1 k
p1(yˆ ,yˆ ,...,yˆ )=0,...,pm1(yˆ ,yˆ ,...,yˆ )=0}.
1 1 2 k 1 1 2 k
As westatedabove,Ω isthesetofallpossiblevaluesofthe objectivefunctions atthe feasiblesolutionsofMOPBP .
f,h
We are looking for the nondominated solutions that are associated with the efficient solutions. From Ω, we can select
the efficient solutions as those that are minimal with respect to the componentwise order in Rk. So, we can extract
from Ω the set of efficient solutions, Y :
E
Y ={(y∗,...,y∗)∈Ω:6∃(y′,...,y′)∈Ω with y′ ≤y∗ for j =1,...,k and (y′,...,y′)6=(y∗,...,y∗)}
E 1 k 1 k j j 1 k 1 k
Once wehaveobtainedthe solutionsinthey variablesthatareefficientsolutionsforMOPBP , wecompute with an
f,h
analogousprocedurethe nondominatedsolutionsassociatedtothe y-valuesinY . It consistsofsolvingthe triangular
E
system given by G for the polynomial where the x-variables appear once the values for the y-variables are fixed to be
each of the vectors in Y .
E
A pseudocode for this procedure is described in Algorithm 1.
Theorem 3.1. Algorithm 1 either provides all nondominated and efficient solutions or provides a certificate of infea-
sibility whenever G={1}.
Algorithm 1: Solving MOPIP by solving systems of polynomial equations
Input : f ,...,f , h ,...h ∈R[x ,...,x ]
1 k 1 s 1 n
Initialization: I =hf −y ,...,f −y ,h ,...,h ,x (x −1),...,x (x −1)i.
1 1 k k 1 s 1 1 n n
Algorithm:
Step 1.: Compute a Gr¨obner basis, G, for I with respect to a lexicographic order with
y ≺···≺y ≺x ≺···≺x .
k 1 n 1
Step 2.: Let Gy =G∩R[y ,...,y ] be a Gr¨obner basis for Iy =I∩R[y ,...,y ], for l=0,...,k−1. (By
l l+1 k l l+1 k
the Elimination Property).
(1) Find all yˆ ∈V(Gy ).
k k−1
(2) Extend every yˆ to (yˆ ,yˆ )∈V(Gy ).
k k−1 k k−2
.
.
.
(k−1) Extend every (yˆ ,...,yˆ ) to (yˆ ,yˆ ,...,yˆ )∈V(Gy).
3 k 2 3 k 1
(k) Find all yˆ such that (yˆ ,...,yˆ )∈V(Gy).
1 1 k 0
Step 3.: Select from V(Gy) the minimal vectors with respect to the usual componentwise order in Rk. Set Y
0 E
this subset.
Step 4.: Let G =G∩R[y ,...,y ,x ,...,x ] a Gr¨obner basis for I ∩R[y ,...,y ,x ,...,x ], for
l 1 k l+1 n l 1 k l+1 n
l=0,...,n−1. (By the Elimination Property). Denote by
S ={(yˆ ,...,yˆ ,xˆ ,...,xˆ ):(yˆ ,...,yˆ )∈Y , and ∃(x ,...,x ) such that (x ,...,x ) is feasible } for
l 1 k l+1 n 1 k E 1 l 1 n
l=0,...,n−1.
(1) Find all xˆ such that (yˆ ,...,yˆ ,xˆ )∈V(G )∩S .
n 1 k n n−1 n−1
(2) Extend every xˆ to ((yˆ ,...,yˆ ,xˆ ,xˆ )∈V(G )∩S .
n 1 k n−1 n n−2 n−2
.
.
.
(n−1) Extend every (yˆ ,...,yˆ ,xˆ ,...,xˆ ) to (yˆ ,...,yˆ ,xˆ ,xˆ ,...,x )∈V(G )∩S .
1 k 3 n 1 k 2 3 n 1 1
(n) Find all xˆ such that (yˆ ,...,yˆ ,xˆ ,...,xˆ )∈V(G )∩S .
1 1 k 1 n 0 0
Set X =π (V(G )∩S ), where π denotes the projection over the x-variables.
E x 0 0 x
Output: Y the set of efficient solutions and X the set of nondominated solutions for MOPBP .
E E f,h
Proof. Suppose that G6={1}. Then, Gy has exactly one element, namely p(y ). This follows from the observation
k−1 k
that I ∩R[y ] is a polynomial ideal in one variable, and therefore, needs only one generator.
k
Solvingp(y )=0we obtaineveryyˆ ∈V(Gy ). Sequentiallywe obtainyˆ extending yˆ to the partialsolutions
k k k−1 k−1 k
(yˆ ,yˆ ) in V(Gy ) and so on.
k−1 k k−1
By the Extension Theorem, this is always possible in our case.
Continuing in this way and applying the Extension Theorem, we can obtain all solutions (yˆ ,...,yˆ ) in V(G∩
1 k
R[y ,...,y ]. These vectors are all the possible objective values for all feasible solutions of the problem. Selecting
1 k
from V(G∩R[y ,...,y ]) those solutions that are not dominated in the componentwise order in Rk, we obtain Y .
1 k E
Following a similar scheme in the x- variables, we have the set V(G )∩S∗ encoding all efficient (in the first k
0 0
coordinates) and nondominated (in the last n coordinates) solutions.
Finally, if G = {1}, then, the ideal I coincides with R[y ,...,y ,x ,...,x ], indicating that V(I) is empty (it
1 k 1 n
is the set of the common roots of all polynomials in R[y ,...,y ,x ,...,x ]). Then, we have an infeasible integer
1 k 1 n
problem. (cid:3)
Remark 3.1. In the case when we have added slack variables, as explained in (1), we slightly modify the above
algorithm solving first in the slack variables and selecting those solutions that are real numbers. Then continue with
the procedure described in Algorithm 1.
Remark 3.2. The Gro¨bner basis, G, computed for solving the system of polynomial equations can be computed with
respect to any other elimination ordering. The only conditions that are required for that ordering is that it allows to
separate the family of x−variables from the family of y-variables and such that the system of polynomials given by that
basis allows solving first for the y-variables and then for the x-variables sequentially.
4. Obtaining nondominated solutions by the Chebyshev norm approach
In this sectionwe describe two more methods for solving MOPIPbased ona different rationale,namely scalarizing
the multiobjective problem andsolvingit as a parametricsingle-objectiveproblem. We proposea methodology based
ontheapplicationofoptimalityconditionstoafamilyofsingle-objectiveproblemsrelatedtoouroriginalmultiobjective
problem. Themethodsconsistoftwomainsteps: afirststepwherethemultiobjectiveproblemisscalarizedtoafamily
ofsingle-objectiveproblemssuchthateachnondominatedsolutionis anoptimalsolutionforatleastoneofthe single-
objective problems in that family; and a second step that consists of applying necessary optimality conditions to
each one of the problems in the family, to obtain their optimal solutions. Those solutions are only candidates to be
nondominated solutions of the multiobjective problem since we just use necessary conditions.
For the first step, the scalarization,we use a weightedChebyshev normapproach. Other weightedsum approaches
could be used to transform the multiobjective problem in a family of single-objective problems whose set of solutions
contains the set of nondominated solutions of our problem. However, the Chebyshev approach seems to be rather
adequate since it does not require to impose extra hypothesis to the problem. This approach can be improved for
problems satisfying convexity conditions, where alternative well-known results can be applied (see [15] for further
details).
Forthe secondstep, weusethe Fritz-JohnandKarush-Kuhn-Tuckernecessaryoptimality conditions,givingustwo
different approaches. In this section we describe both methodologies since each of them has its own advantages over
the other.
ForapplyingtheChebyshevnormscalarization,weusethefollowingresultthatstateshowtotransformourproblem
toafamily ofsingle-objectiveproblems,andhowtoobtainnondominatedsolutionsfromthe optimalsolutionofthose
single-objective problems. Further details and proofs of this result can be found in [15].
Theorem 4.1 (Corollary 11.21 in [15]). Let MOPBP be feasible. x∗ is a nondominated solution of MOPBP
f,g,h f,g,h
if and only if there are positive real numbers ω ,...,ω >0 such that x∗ is an image unique solution of the following
1 k
weighted Chebyshev approximation problem:
min γ
s.t. ω (f (x)−yˆ)−γ ≤0 i=1,...,k
i i i
g (x) ≤0 j =1,...,m
(P ) j
ω h (x) =0 r =1,...,s
r
x (x −1) =0 i=1,...,n
i i
γ ∈R x∈Rn
where yˆ = (yˆ ,...,yˆ ) ∈ Rk is a lower bound of f = (f ,...,f ), i.e., yˆ ≤ f (x) for all feasible solution x and
1 k 1 k i i
i=1,...,k.
According to the above result, every nondominated solution of MOPBP is the unique solution of (P ) for
f,g,h ω
some ω > 0. We apply, in the second step, necessary optimality conditions for obtaining the optimal solutions for
those problems(takingω asparameters). Thesesolutionsarecandidatestobe nondominatedsolutionsofouroriginal
problem. Actually, every nondominated solution is among those candidates.
Inthefollowingsubsectionswedescribetheabove-mentionedtwomethodologiesforobtainingtheoptimalsolutions
for the scalarized problems (P ) for each ω.
ω
4.1. The Chebyshev-Karush-Kuhn-Tucker approach: The first optimality conditions that we apply are the
Karush-Kuhn-Tucker (KKT) necessary optimality conditions, that were stated, for the general case, as follows (see
e.g. [1] for further details):
Theorem 4.2 (KKT necessary conditions). Consider the problem:
min f(x)
s.t. g (x) ≤0 =1,...,m
(4) j
h (x) =0 r =1,...,s
r
x ∈Rn
Let x∗ be a feasible solution, and let J ={j :g (x∗)=0}. Suppose that f and g , for j =1,...,m, are differentiable
j j
at x∗, that g , for j 6∈J, is continuous at x∗,and that h , for r=1,...,s, is continuously differentiable at x∗. Further
j r
suppose that ∇g , for j ∈ I, and ∇h , for r = 1,...,s, are linearly independent (regularity conditions). If x∗ solves
j r
Problem 4 locally, then there exist scalars λ , for j =1,...,m, and µ , for r =1,...,s, such that
j r
m s
∇f(x∗)+ λ ∇g (x∗)+ µ ∇h (x∗) =0
j j r r
(KKT) Xj=1 Xr=1
λ g (x∗) = 0 for j =1,...,m
j j
λ ≥ 0 for j =1,...,m
j
From the above theorem the candidates to be optimal solutions for Problem (4) are those that either satisfy the
KKT conditions (in the case where all the functions involved in Problem (4) are polynomials, this is a system of
polynomialequations)ordonotsatisfytheregularityconditions. Notethatthesetwosetsare,ingeneral,notdisjoint.
Regularity conditions can also be formulated as a system of polynomial equations when the involved functions are
all polynomials. Let x∗ be a feasible solution for Problem (4), x∗ does not verify the regularity conditions if there
exist scalars λ , for j ∈J, and µ , for r =1,...,s, not all equal to zero, such that:
j r
s
(Non-Regularity) λ ∇g + µ ∇h =0
j j r r
Xj∈I Xr=1
The above discussion justifies the following result.
Corollary 4.1. Let x∗ be a nondominated solution for MOPBP . Then, x∗ is a solution of the systems of
f,g,h
polynomial equations (5) or (6), for some ω >0.
k
1− νi =0 
(5) k ν ω ∇f (x)+ m λ ∇g (x)+ s µ ∇h (x)+ n β δ (2xXi=−11)=0 
i i i j j j r l il i
Xi=1 Xj=1 Xr=1 Xl=1
for some λ g (x)=0, for j =1,...,m, λ ≥0, for j =1,..ν.i,ωmi,(fain(dx)ν−≥yˆi0),−foγr=i=0,1i,=...1,,k....,k
j j j i
k
νi =0 
(6) k ν ω ∇f (x)+ m λ ∇g (x)+ s µ ∇h (x)+ n β δ (2xXi=−11)=0 
i i i j j j r l il i
Xi=1 Xj=1 Xr=1 Xl=1
with x∈Rn such that g (x)≤0, for j =1,...,m, h (x)=0, fωoir(rfi=(x1),−..yˆ.i,)s−γ ≤0,i=1,...,k
j r
λ ≥0, for j =1,...,m, and ν ≥0, for i=1,...,k.
j i
Let XKKT denote the set of solutions,in the x-variables,of system (5) andlet XNR denote the set of solutions,in
E E
the x-variables, of system (6) (the problem is solved avoiding inequality constraints, then every solution is evaluated
to check if it satisfies the inequality constraints).
For solving these systems (Chebyshev-KKT and Non-Regulariry), we use a Gr¨obner bases approach. Let I be the
ideal generated by the involved equations.
Let us consider a lexicographical order over the monomials in R[x,γ,λ,ν,µ,β] such that x ≺ γ ≺ λ≺ ν ≺ µ ≺ β.
Then, the Gr¨obner basis, G, for I with this order has the following triangular shape:
• G contains one polynomial in R[x ]: p (x )
n n n
• G contains one or severalpolynomials in R[x ,x ] : p1 (x ,x ),...,pm1 (x ,x ).
n−1 n n−1 n−1 n n−1 n−1 n
.
.
.
• G contains one or severalpolynomials in R[x] : p1(x ,...,x ),...,pmn(x ,...,x ).
1 1 n 1 1 n
• The remaining polynomials involve variables x and at least one γ, λ, µ, ν or β.
We are interested in finding only the values for the x-variables, so, we avoid the polynomials in G that involve any
oftheotherauxiliaryvariables. Ingeneral,wearenotabletodiscussaboutthevaluesoftheparametersγ,λ,µ,ν and
β. Needless to say that in those cases when we can do it, some values of x may be discarded simplifying the process.
We denote by Gx the subset of G that contains only polynomials in the x-variables. By the ExtensionTheorem, Gx is
a Gr¨obner basis for I∩R[x ,...,x ].
1 n
Solvingthe systemgivenby Gx andcheckingfeasibility ofthose solutions,weobtainas solutionsthose ofourKKT
or Non-Regularity original systems.
ItisclearthatthesetofnondominatedsolutionsofourproblemisasubsetofXKKT∪XNR,sinceeitherasolution
E E
is regular, and then, KKT conditions are applicable or it satisfies the non regularity conditions. However, the set
XKKT ∪XNR may contain dominated solutions, so, at the end we must remove the dominated ones to get only X .
E E E
The steps to solve Problem MOPBP using the Chebyshev-KKT approach are summarized in Algorithm 2.
f,g,h
Algorithm 2: Summary of the procedure for solving MOPBP using Chebyshev-KKT approach.
Input : f ,...,f , g ,...g ,h ,...,h ∈R[x ,...,x ]
1 k 1 m 1 r 1 n
Algorithm:
Step 1: Formulate the Chebyshev scalarization of MOPBP . (Problem (P ))
f,g,h ω
Step 2: Solve System (5) in the x-variables: XKKT.
E
Step 3: Solve System (6) in the x-variables: XNR.
E
Step 4: Remove from XKKT ∪XNR the subset of dominated solutions: X .
E E E
Output: X the set of nondominated solutions for MOPBP
E f,g,h
Theorem 4.3. Algorithm 2 solves Problem MOPBP in a finite number of steps.
f,g,h
4.2. The Chebyshev-Fritz-John approach. Analogously to the previous approach, once we have scalarized the
original multiobjective problem to a family of single-objective problems, in this section we apply the Fritz-John (FJ)
conditions to all the problems in this family. The following well-known result justifies the use of FJ conditions to
obtain candidates to optimal solutions for single-objective problems. Proofs and further details can be found in [1].
Theorem 4.4 (FJ necessary conditions). Consider the problem:
min f(x)
s.t. g (x) ≤0 =1,...,m
(7) j
h (x) =0 r =1,...,s
r
x ∈Rn
Let x∗ be a feasible solution, and let J ={j :g (x∗)=0}. Suppose that f and g , for j =1,...,m, are differentiable
j j
at x∗, and that h , for r = 1,...,s, is continuously differentiable at x∗. If x∗ locally solves Problem (7), then there
r
exist scalars λ , for j =1,...,m, and µ , for r =1,...,s, such that
j r
m s
λ ∇f(x∗)+ λ ∇g (x∗)+ µ ∇h (x∗) =0
0 j j r r
Xj=1 Xr=1
(FJ) λ g (x∗) = 0 for j =1,...,m
j j
λ ≥ 0 for j =1,...,m
j
(λ ,λ,µ) 6=(0,0,0)
0
Note that, in the FJ conditions, regularity conditions are not required to set the result.
Corollary4.2. Letx∗ beanondominatedsolutionforMOPBP . Then, x∗ isasolutionofthesystemofpolynomial
f,g,h
equations (8) for some ν ,λ ,µ ,β , for i=1,...,k,l=1,...,n, j =1,...,m, r =1,...,s and ω >0.
i j r l
k
λ − ν =0
0 i 
(8) k νiωi∇fi(x)+ m λj∇gj(x)+ s µj∇hr(x)+ n βlδil(2xiXi=−11)=0 
Xi=1 Xj=1 Xr=1 Xl=1
where λ ≥ 0, for j = 1,...,m, ν ≥ 0, for i = 1,...ν,ik(ωain(dfiδ(x,)−forλyˆijl),g−ij(=xγ))1==,.00.,,.ij,=n=,10d,,.e..n..o.,t,kems the Kronecker delta
j i li
1 if l=i
function, δ = .
li (cid:26) 0 otherwise
Let XFJ denote the set of solutions, in the x-variables, that are feasible solutions of MOPBP and solutions of
E f,g,h
system (8).
The set of nondominated solutions of our problem is a subset of XFJ, since every nondominated solution is an
E
optimal solution for some problem in the form of (8), and every solution of this single-objective problem is a solution
of the FJ system.
However, dominated solutions may appear in the set of solutions of (8), so, a final elimination process is to be
performed to select only the nondominated solutions.
The steps to solve MOPBP using the Chebyshev-FJ approachare summarized in Algorithm 3.
f,g,h
Algorithm 3: Summary of the procedure for solving MOPBP using the Chebyshev-FJ approach.
Input : f ,...,f , g ,...g ,h ,...,h ∈R[x ,...,x ]
1 k 1 m 1 r 1 n
Algorithm:
Step 1: Formulate the Chebyshev scalarization of MOPBP . (Problem (P ))
f,g,h ω
Step 2: Solve system (8) in the x-variables for any value of ω >0: XFJ.
E
Step 3: Remove from XFJ the set of dominated solutions: X .
E E
Output: X the set of nondominated solutions for Problem (MOPBP )
E f,g,h
Theorem 4.5. Algorithm 3 solves MOPBP in a finite number of steps.
f,g,h
The last part of the section is devoted to show how to solve the Chebyshev-FJ system using Gr¨obner bases.
Consider the following polynomial ideal
k k m s n
I =hλ0−Xνi,Xνiωi∇fi(x)+Xλj∇gj(x)+Xµj∇hr(x)+Xβlδil(2xi−1),ν1(ω1(f1(x)−yˆ1)−γ),...,νk(ωk(fk(x)−
i=1 i=1 j=1 r=1 l=1
yˆk)−γ),λ1g1(x),...,λmgm(x)i
in the polynomial ring R[x,γ,λ,ν,µ,β].
Let us consider a lexicographical order over the monomials in R[x,γ,λ,ν,µ,β] such that x ≺ γ ≺ λ≺ ν ≺ µ ≺ β.
Then, the Gr¨obner basis, G, for I with this order has the following triangular shape:
• G contains one polynomial in R[x ]: p (x )
n n n
• G contains one or severalpolynomials in R[x ,x ] : p1 (x ,x ),...,pm1 (x ,x )
n−1 n n−1 n−1 n n−1 n−1 n
···
• G contains one or severalpolynomials in R[x] : p1(x ,...,x ),...,pmn(x ,...,x )
1 1 n 1 1 n
• The remainder polynomials involve variables x and at least one of γ, λ, µ, ν or β.
We are interested in finding only the values for the x-variables, so, we avoid the polynomials in G that involve
any of the other auxiliary variables. We denote by Gx the subset of G that contains only all the polynomials in the
x-variables. By the Extension Theorem, Gx is a Gr¨obner basis for I∩R[x ,...,x ].
1 n
Solving the system given by Gx, we obtain as solutions, those of our FJ original system.
Remark 4.1 (Convex Case). In the special case where both objective functions and constraints are convex, sufficient
KKT conditions can be applied. If the feasible solution x∗ satisfies KKT conditions, and all objective and constraints
functions are convex, then x∗ is a nondominated solution. As a particular case, this situation is applicable to linear
problems.
In this case, we may choose a linear scalarization instead of the Chebyshev scalarization. With this alternative
approach, the scalarized problem is
k
min t f (x)
s s
Xs=1
s.t. g (x) ≤0 j =1,...,m
j
h (x) =0 r =1,...,s
r
x (x −1) =0 i=1,...,n
i i
for t ,...,t >0.
1 k
Then, by Corollary 11.19 in [15], and denoting by S the feasible region, if f(S)+Rk is convex, then each x∗ is a
+
nondominated solution if and only if x∗ is a solution of Problem 4.1 for some t ,...,t >0.
1 k
Using both results, necessary and sufficient conditions are given for that problem and the removing step is avoided.
Remark 4.2 (Single-Objective Case). The same approach can be applied to solve single-objective problems. In this
case, KKT (or FJ) conditions can be applied directly to the original problem, without scalarizations.
5. Obtaining nondominated solutions by multiobjective optimality conditions
Inthis section,we addressthe solutionofMOPBP by directly applyingnecessaryconditionsfor multiobjective
f,g,h
problems. With theseconditionswedonotneedtoscalarizethe problem,asinthe abovesection,avoidingsomesteps
in the process followed in the previous sections.
The following result states the Fritz-John necessary optimality conditions for multiobjective problems.
Theorem 5.1 (Multiobjective FJ necessary conditions, Theorem 3.1.1. in [19]). Consider the problem:
min (f (x),...,f (x))
1 k
s.t. g (x) ≤0 =1,...,m
(9) j
h (x) =0 r =1,...,s
r
x ∈Rn
Let x∗ a feasible solution. Suppose that f , for i = 1,...,k, g , for j = 1,...,m and h , for r = 1,...,s, are
i j r
continuously differentiable at x∗. If x∗ is a nondominated solution for Problem 9, then there exist scalars ν , for
i
i=1,...,k, λ , for j =1,...,m, and µ , for r =1,...,s, such that
j r
k m s
ν ∇f (x∗)+ λ ∇g (x∗)+ µ ∇h (x∗) =0
i i j j r r
Xi=1 Xj=1 Xr=1
(MO-FJ) λjgj(x∗) = 0 for j =1,...,m
λ ≥ 0 for j =1,...,m
j
ν ≥ 0 for i=1,...,k
i
(ν,λ,µ) 6=(0,0,0)
Withthisresult,onecansolvethesystemgivenbythenecessaryconditionstoobtaincandidatestobenondominated
solutions for our problem. For solving this system, we use lexicographicalGr¨obner bases as in the abovesections. We
summarize the algorithm for solving the multiobjective polynomial problem in Algorithm 4.
Algorithm 4: Summary ofthe procedureforsolvingMOPBPusing the multiobjective FJoptimality conditions.
Input : f ,...,f , g ,...g ,h ,...,h ∈R[x ,...,x ]
1 k 1 m 1 r 1 n
Algorithm:
Step 1: Solve system (MO-FJ): XMOFJ.
E
Step 2: Remove from XMOFJ the subset of dominated solutions: X .
E E
Output: X the set of nondominated solutions for Problem MOPBP
E f,g,h
Remark5.1. Inthespecialcasewherebothobjectivefunctionsandconstraintsareconvex,Theorem5.1givessufficient
nondominance conditions for MOPBP requiring that ν >0 (see Theorem 3.1.8 in [19]).
f,g,h i
6. Computational Experiments
A series of computational experiments have been performed in order to evaluate the behavior of the proposed
solution methods. Programs have been coded in MAPLE 11 and executed in a PC with an Intel Core 2 Quad
processor at 2x 2.50 Ghz and 4 GB of RAM. The implementation has been done in that symbolic programming
language, available upon request, in order to make the access easy to both optimizers and algebraic geometers.
Werunthealgorithmsforthreefamiliesofbinarybiobjectiveandtriobjectiveknapsackproblems: linear,quadratic
and cubic, and for a biobjective portfolio selection model. For each problem, we obtain the set of nondominated
solutions as well as the CPU times for computing the corresponding Gr¨obner bases associated to the problems, and
the total CPU times for obtaining the set of solutions.
We give a short description of the problems where we test the algorithms. In all cases, we use binary variables
x , j = 1,...,n, where x = 1 means that the item (resp. security) j is selected for the knapsack (resp. portfolio)
j j
problem.
(1) Biobjective (linear) knapsack problem (biobj linkn): Assume that n items are given. Item j has associated
costsq1,q2 fortwodifferenttargets,andaunitprofita ,j =1,...,n. Thebiobjectiveknapsackproblemcalls
j j j
for selecting the item subsets whose overall profit ensures a knapsack with value at least b, so as to minimize
(in the nondominance sense) the overall costs. The problem may be formulated:
n n
min( q1x , q2x )
j j j j
Xj=1 Xj=1
n
s.t. a x ≥b, x∈{0,1}n
i i
Xi=1
(2) Biobjective cubicknapsack problem(biobj cubkn): Assumethatnitemsaregivenwhereitemj hasaninteger
profita . Inadditionwearegiventwon×n×nmatricesP1 =(p1 )andP2 =(p2 ),wherep1 andp2 are
j ijk ijk ijk ijk
the costs for each of the targets if the combinationof items i,j,k is selected for i<j <k; and two additional
n×n matrices Q1 =(q1) and Q2 =(q2), where q1 and q2 are the costs for the two different targets if both
ij ij jj ij
items iandj areselectedfori<j Thebiobjectivecubic knapsackproblemcallsforselectingthe itemsubsets
whose overall profit exceeds the purpose of the knapsack b, so as to minimize the overall costs. The problem
may be formulated:
n n n−2 n−1 n n n n−2 n−1 n
min q1 x x + p1 x x x , q2 x x + p2 x x x
ij i j ijl i j l ij i j ijl i j l
(cid:0)Xi=1Xj=i Xi=1j=Xi+1l=Xj+1 Xi=1Xj=i Xi=1j=Xi+1l=Xj+1 (cid:1)
n
s.t. a x ≥b, x∈{0,1}n
i i
Xi=1
(3) Biobjective quadratic knapsack problem (biobj qkn): This problem may be seen as a special case of the
biobjective cubic knapsack problem when there are no cost correlations between triplets.
(4) Triobjective (linear) knapsack problem (triobj linkn): Assume that n items are given where item j has an
integer profit a . In addition, we are given three vectors q1 = (q1), q2 = (q2) and q3 = (q3), where q1,q2 and
j j j j j j
q3 arethe costsforthreedifferenttargetsifj isselected. The triobjectiveknapsackproblemcallsforselecting
j
the item subsets whose overall profit ensures a profit for the knapsack at least b, so as to minimize (in the
nondominance sense) the overall costs. The problem is:
n n n
min( q1x , q2x , q3x )
j j j j j j
Xj=1 Xj=1 Xj=1
n
s.t. a x ≥b, x∈{0,1}n
i i
Xi=1
(5) Triobjective cubic knapsack problem (triobj cubkn): We aregivenn items whereitem j hasanintegerprofit
a . Inaddition,wearegiventhreen×n×nmatricesP1 =(p1 ),P2 =(p2 )andP3 =(p3 ),wherep1 ,p2
j ijk ijk ijk ijk ijk
and p3 are the costs for each of the targets if the combination of items i,j and k is selected for i < j < k;
ijk
and three additional n×n matrices Q1 =(q1),Q2 =(q2) and Q3 =(q2), where q1 ,q2 and q3 are the costs
ij ij ij jj ij ij
for three different targets if both items i and j are selected for i<j The triobjective cubic knapsack problem
callsforselectingthe itemsubsetswhoseoverallprofitensuresavalueofb,soasto minimizethe overallcosts.
The problem is:
n n n−2 n−1 n n n n−2 n−1 n
min q1 x x + p1 x x x , q2 x x + p2 x x x ,
ij i j ijl i j l ij i j ijl i j l
(cid:0)Xi=1Xj=i Xi=1j=Xi+1l=Xj+1 Xi=1Xj=i Xi=1j=Xi+1l=Xj+1
n n n−2 n−1 n
q3 x x + p3 x x x
ij i j ijl i j l
Xi=1Xj=i Xi=1j=Xi+1l=Xj+1 (cid:1)
n
s.t. a x ≥b, x∈{0,1}n
i i
Xi=1
(6) Triobjective quadratic knapsack problem (triobj qkn): This problem may be seen as a special case of the
triobjective cubic knapsack problem when there are no cost correlations between triplets.
(7) Biobjective portfolio selection (portfolio): Consider a market with n securities. An investor with initial
wealth b seeks to improve his wealth status by investing it into these n risky securities. Let X be the
i
random return per a lot of the i-th secutiry (i = 1,...,n). The mean, µ = E[X ], and the covariance,
i i
σ = Cov(X ,X ),i,j = 1,...,n, of the returns are assumed to be known. Let x be a decision variable
ij i j i
that takes value 1 if the decision-maker invests in the i-th security and 0 otherwise. Denote the decision
vectorbyx=(x ,...,x ). Then, the randomreturnfora inversionvectorx fromthe securitiesis x X
1 n i=1 i i
n
P
and the mean and variance of this random variable are E[ x X ] = µ x and Var( x X ) =
i=1 i i i i i=1 i i
P Xi=1 P
n n
x x σ .
i j ij
Xi=1Xj=1
Leta be thecurrentpriceofthei-thsecurity. Then,ifaninvestorlooksforminimizinghisinvestmentrisk
i
and simultaneously maximizing the expected return with that investment, the problem can be formulated as:
n n n n
min(Var( x X ),−E[ x X ]) min( σ x x ,− µ x )
i i i i ij i j i i
 
s.t. n a xXi=1≤b, x∈{0Xi,=11}n  ⇒ s.t. nXi=1a x ≤b, xXi=∈1{0,1}n 
i i i i
For each of the above 7Xic=la1sses of problems, we considerinstances raXin=d1omly generated as follows: a is randomly
i
drawn in [−10,10] and the coefficients of the objective functions, qk, pk , σ and µ , range in [−10,10]. Once the
ij ijl ij i
constraint vector, (a ,...,a ), is generated, the right hand side, b, is randomly generated in [1,| n a |]. For each
1 n i=1 i
type of instances and each value of n in [2,13] we generated 5 instances. P
Tables 1 and 3 contain a summary of the average results over the different instances generated for the above
problems. Each algorithm is labeled conveniently: alg1 corresponds with Algorithm 1, kkt is Algorithm 2, kkt sl
is Algorithm 2 where the inequality is transformed to an equation using a slack variable, fj is Algorithm 3, fj sl is
Algorithm 3 where the inequality is transformed to an equation using a slack variable and mofj stands for Algorithm
4. For each of these algorithms we present the CPU time for computing the corresponding Gr¨obner basis (tgb), the
total CPU time for obtaining the set of nondominated solutions (ttot), the number of nondominated solutions (#nd)
and the number of variables involved in the resolution of the problem (#vars).
Fromthosetables,thereadermaynotethatAlgorithm1isfasterthantheothersforthesmallestinstances,although
the CPU times for this algorithm increase faster than for the others and it is not able to obtain solutions when the
size of the problem is around 12 variables. The algorithms based on Chebyshev scalarization (kkt, kkt sl, fj and
fj sl) are better than alg1 for the largest instances. The differences between these four methods are meaningful,
but the algorithms based on the KKT conditions are, in almost all the instances, faster than those based on the FJ
conditions. Note that considering slack variablesto avoidthe inequality constraintis not better, since the CPU times
when the slackvariableis consideredarelarger. Finally, the bestalgorithm,inCPUtime, is mofjsinceexceptfor the
small instances is the fastest and it was able to solve larger instances.
Onemaythinkthatthelaststepofourmethods,i.e. removingdominatedsolutions,shouldbemoretimeconsuming
inalg1thanintheremainingmethodssincealg1doesnotuseoptimalityconditions. However,fromourexperiments
this conclusion is not clearly supported. Actually, although this process, in time consuming, when the dimension of
the problem increases this time is rather small compared with the effort necessary to obtain the Gr¨obner bases.
Table 2 shows some information about each of the presented algorithms. For a multiobjective problem with n
variables, m polynomial inequality constraints given by g = (g ,...,g ), s polynomial equality constraints given by
1 m
h = (h ,...,h ) and k objectives functions given by f = (f ,...,f ), Table 3 shows the number of variables (#var),
1 s 1 k
the numberofgenerators(#gen)andthe maximaldegrees(maxdeg)ofthe initialpolynomialidealsrelatedtotheeach
of the algorithms. These numbers inform us about the theoretical complexity of the algorithms. The computation
of Gr¨obner bases depends of the number of variables (in general, double exponential) and of the size of the initial
systemofgenerators(degreesandnumberofpolynomials). Actually,itisknownthatcomputingaGr¨obnerbasisusing
Buchberger Algorithm is doubly exponential in the number of variables. Some complexity bounds for this algorithm
involving #var, #gen and maxdeg can be found in [12].