Virtual interpolation of discrete multi-objective programming solutions with probabilistic operation

Interpolação virtual de soluções de programação multi-objetivo discreta com operação probabilística

OTIMIZAÇÃO EM CONTROLE E AUTOMAÇÃO

Virtual interpolation of discrete multi-objective programming solutions with probabilistic operation

Interpolação virtual de soluções de programação multi-objetivo discreta com operação probabilística

Ricardo C. Silva^I; Edilson F. Arruda^II; Fabrício O. Ourique^III

^IInstitute of Science and Technology, Federal University of São Paulo, Rua Talim, 330, São José dos Campos, SP, 12231-280, Brazil. [e-mail]

^IIDepartment of Production Enginnering, Alberto Luiz Coimbra Institute for Graduate Studies and Research in Engineering, Federal University of Rio de Janeiro, P.O. Box 68507, Rio de Janeiro, RJ, 21941-972, Brazil. [e-mail]

^IIIComputer Engineering, Federal University of Pampa, Campus Bagé, Travessa 45, 1650, Bagé, RS, 96413-170, Brazil. [e-mail]

ABSTRACT

This work presents a novel framework to address the long term operation of a class of multi-objective programming problems. The proposed approach considers a stochastic operation and evaluates the long term average operating costs/profits. To illustrate the approach, a two-phase method is proposed which solves a prescribed number of K mono-objective problems to identify a set of K points in the Pareto-optimal region. In the second phase, one searches for a set of non-dominated probability distributions that define the probability that the system operates at each point selected in the first phase, at any given operation period. Each probability distribution generates a vector of average long-term objectives and one solves for the Pareto-optimal set with respect to the average objectives. The proposed approach can generate virtual operating points with average objectives that need not have a feasible solution with an equal vector of objectives. A few numerical examples are presented to illustrate the proposed method.

Keywords: Pareto-optimality, Dynamic operation, Discrete optimization

RESUMO

Apresenta-se, neste trabalho, um novo método para se tratar a operação em horizonte infinito de uma classe de problemas de programação multi-objetivo. A abordagem proposta considera uma operação estocástica e avalia o custo/lucro médio em horizonte infinito. Para ilustrar a abordagem proposta, um método de duas fases é proposto para resolver um número pré-determinado de K problemas multi-objetivos, de modo a se identificar um conjunto com K pontos pertencentes à região Pareto-ótima. A segunda fase consiste em se buscar um conjunto não-dominado de distribuições de probabilidade no domínio dos K pontos de operação selecionados na primeira fase, que define a probabilidade de se escolher cada um dos pontos de operação em todo instante de decisão. Cada distribuição de probabilidade gera um vetor de funções objetivo no longo prazo e determina-se o conjunto Pareto-ótimo com respeito às médias dos objetivos. A abordagem proposta pode gerar pontos de operação virtuais, com médias de função objetivo que não necessariamente correspondem a um ponto factível de operação. Alguns experimentos numéricos são desenvolvidos para ilustrar a abordagem proposta.

Palavras-chave: Otimização de Pareto, Operação dinâmica, Otimização discreta

1 INTRODUCTION

Humans have to choose among a finite number of alternatives to make decisions according to some criterion or objective. However, many real-world problems do not involve a single objective. Rather, they have multiple, possibly conflicting objectives to be optimized simultaneously. The wide range of problems in this category motivated the development of the field of Multiple Criteria Decision Making (MCDM), also known as Multi-Objective Programming (MOP). This field dates back to the early works by PARETO (1897a) and PARETO (1897b) that introduced the concept of corporative equilibrium. A survey about the history and social motivations for the development of this equilibrium concept is presented by STADLER (1979); it covers the interval between 1776 and 1960. CHANKONG and HAIMES (1983) present a formulation of the MOP problem in the space of the objectives, as well as some definitions of Pareto-optimality concepts that are used in most of the existing MOP methods.

As the term suggests, MCDM is a branch of Operational Research (OR) primarily concerned with decision problems involving multiple objectives, e.g. (ZIONTS, 1992). The fundamental guidelines regarding the decision making in multiobjective problems can be found in the work by KEENEY and RAIFFA (1976). An effort is conducted by REID and CITRON (1971) to establish some basic properties of the noninferior surface in general multi-objective problems. In another line of research, conditions for the existence of a single optimal solution to multi-objective problems were investigated in (ATHANS and GEERING, 1973; SALUKVADZE, 1974). A study on sensitivity analysis tools, designed to investigate the influence in the output of any given change in the input, was presented by INSUA (1990).

One way to tackle MOP problems involves defining a mapping from the set of feasible solutions to function space. Such a mapping is often called utility function and allows a complete ordering of the set of feasible solutions. Therefore, the knowledge of the utility function makes it possible to reformulate a MOP as an ordinary mono-objective optimization problem with the utility function acting as the objective. However, in many situations it is not possible to identify or assess a utility function. In those cases, a partial ordering of the feasible solutions in terms of the values of the individual objective functions can be used. Under such an ordering, one can identify particular solutions that cannot be thoroughly improved, i.e. solutions such that no other solution exists that improves all objective functions simultaneously. Any solution of this type is called nondominated, and multi-objective programming methods are developed to identify the set of nondominated solutions, also known as the Pareto-optimal set or Pareto front (or subsets of it). After a set of nondominated solutions is identified, it is handed in to the decision maker and she is assigned the task of choosing a particular solution from this set.

KIM and de WECK (2006) try to map the Pareto front by using a linear combination of neighbor points, which are obtained by a traditional weighted sum method for multiobjective optimization. This method is called adaptive weighted sum because it does not allow two or more different weight vectors to achieve the same point in the Pareto front. The optimality conditions to some dynamic multiobjective programming problems cannot be applied, and ESOGBUE et al. (2006) define alternative conditions for this kind of problems. For problems with non-mensurable objective functions, evolutionary computation has arisen as a good alternative to obtain a satisfactory Pareto front within a resonable time window. Some evolutionary algorithms for MOP problems can be found in (DEB, 2001; COELLO et al., 2002). For a survey about concepts and specific methods to solve multiobjective optimization problems, we refer to MARLER and ARORA (2004) and (CHINCHULUUN and PARDALOS, 2007). These methods, which can be divided in three main categories: (i) a-Priori methods, (ii) a-Posteriori methods, and (iii) interactive methods, can be exact or heuristic. Exact methods use concepts such as differentiability to reach the optimal solution, whereas heuristic methods, such as evolutionary algorithms, sweep a selected part of the domain in their search for a satisfactory solution. A disadvantage of heuristic methods is that they are not guaranteed to find the optimal solution. Nevertheless, heuristics typically run faster and are less computationally demanding. A wide range of difficult real-world multiobjective problems have been dealt with by means of evolutionary computation. CARDOSO et al. (2009) develop an evolutionary approach to control biological pests in soybeans crop. EVERSON and FIELDSEND (2006) apply an evolutionary strategy in safety related systems, which monitor the location of airplanes by radars and send alert messages to the flight controllers. ABIDO (2006) applies evolutionary algorithms in electrical power dispatch problems. AMORIM et al. (2009) propose an evolutionary algorithm called AEMO and apply it to a study on the decentralization of the Brazilian electrical energy system. In (DIAS et al., 2008) a memetic algorithm is applied in multiobjective dynamic location problems, whereas in (TRALHÃO et al., 2010) a variation of the facility location problem is presented whose goal is to determine a set of waste collection points for different sorts of waste in a urban region.

Even though the MOP formulation of any given system is static, it can be related to the dynamic operation of the system. Think, for instance, of a multi-objective allocation problem in which clients are allocated to existing facilities first thing in the morning on every business day. If the decision maker criterion does not change during a given stretch of time, the same allocation is selected every single day during that time interval. In such a problem, a vector comprised of the averages of all objective functions during the selected interval (which in this case coincides with the vector of objective functions at the selected operating point) can be used as an alternative objective function. It is worth pointing out that, although the MOP literature is rich and vast, there are not, to the best of our knowledge, works dealing with the dynamic or probabilistic operation of systems and multi-objective criterion. It is the probabilistic selection of operating points for the dynamic operation of MOP systems that we address in this work.

The contribution of this paper is to introduce a novel framework for dynamic multiple criteria decision making (MCDM). The paper proposes a probabilistic selection of the operating point for systems with static or slowly changing decision making criteria. The proposed framework can be viewed as a generalization of the classical static MOP framework. Under the former, one is concerned with selecting an operating point for the system under multiple (possibly conflicting) optimization criteria. The latter, on the other hand, deals with the dynamic operation of the system. This involves selecting different operating points and establishing rules that govern the switch between any two of the selected operating points.

The paper is organized as follows. Section 2 presents the mathematical formulation and the basic concepts of multiobjective programming. Section 3 describes and motivates the proposed probabilistic approach to multi-objective programming. Section 4 features a procedure for selecting of operating points and establishing transition probabilities among these points. Section 5 presents a numerical simulation for a facility location problem and an analysis of the obtained outcomes. Finally, concluding remarks can be found in Section 6.

2 TRADITIONAL MATHEMATICAL FORMULATION

A classical multi-objective problem can be formulated as follows:

where F = [ƒ₁, ƒ₂, ..., ƒ_m]^T, (m > 2) is a vector of objectives and Ω is the set of feasible solutions.

Often the objectives are conflicting and there exists no feasible solution that minimizes all objectives simultaneously. To tackle this problem, one uses the concept of non-dominated solutions. This concept was introduced in Vilfredo Pareto's classical works (PARETO, 1897a; PARETO, 1897b). A solution x* ∈ Ω is said to be non-dominated (or Pareto-optimal, or efficient) if there exists no alternative solution in Ω that is best or equal in all objectives simultaneously and that strictly improves at least one objective.

For any pair x, y ∈ Ω, x ≠ y, define the following partial ordering in terms of the objectives in Problem (1):

When F(x) < F(y) we say that x dominates y. On the other hand, if F(x) > F(y), we say that x is dominated by y. Note that, when a feasible solution dominates another one, it is certainly a better solution to (1) than its counterpart; when it is dominated by the counterpart, it is certainly a worse solution to (1). Based on the relations in (2) one can construct the following sets for any point x⁰∈ Ω Observe that:

1. Ω = Ω _< (x⁰) ∪Ω_>(x⁰) ∪Ω _~ (x⁰).

2. Ω_< (x⁰) and Ω_>(x⁰) denote, respectively, the set of points in Ω that dominate and are dominated by x⁰.

3. Ω _~ (x⁰) lists every point in Ω which cannot be compared to x⁰ according to the partial ordering employed.

4. For any global (efficient) solution to x^*∈ Ω to (1), Ω∩Ω _< (x^*) = Ø.

A local (locally efficient) solution to (1) can also be defined in terms of the sets introduced above. We say that x^* is a local solution to the problem in question whenever

(x*, ε) ∩Ω<(x*) = Ø,

for some scalar ε > 0, where (x^*, ε) = {x ∈ Ω : ║x - x^*║ < ε}; ║·║ denotes the Euclidean norm.

Specific methods to solve multi-objective programming problems were addressed in (CHANKONG and HAIMES,1983; MARLER and ARORA, 2004; CHINCHULUUN and PARDALOS, 2007). These methods are classified according to the instant the decision maker applies their criteria. Three methods are proposed: (i) a-Priori method, where the decision maker assigns preferences to the objective functions a priori, thus obtaining a single mono-objective criterion; (ii) a-Posteriori method, where the decision maker strives to create a set of (locally) efficient solutions to choose from a posteriori; and (iii) interactive methods, where the decision maker informs their preferences during the search for an efficient solution. The first category is used when the decision maker has preferences concerning the objectives. In this category, one can highlight the lexicographic method, which requires a list of objectives in decreasing order of preference and the bounded objective function method, that selects a single most preferred objective function and transforms the remaining objectives into constraints with fixed bounds chosen by decision maker. The second category aims at generating a set of non-dominated solutions to be chosen from at a later stage by the decision maker. This category features the weighted method, which obtains each point in the non-dominated set as a solution to a transformed single objective problem, where the objective is a convex combination of the original objective functions, and the ε-constrained method, which is similar to the bounded objective function, but where the bounds are not flexible. The last category obtains the decision maker's preferences during the run of the algorithm, with a view of generating a desirable operating point. Under this framework, one can highlight the STEP method, which is exact and the evolutionary algorithm GENOCOP, e.g, (SAKAWA and YAUCHI, 1998). It is worth pointing out that evolutionary algorithms are used to seek solutions in all three categories. However, they are typically classified as a-Posteriori methods.

The traditional formulation described above will help us introduce in the next section the motivation of this work. We shall work with a generalized dynamic approach that results in a multi-objective linear formulation, regardless of the structure of the objective functions in (1). The classical formulation serves as a reference for comparison, and can also be used to generate the operation points to be used in the dynamic operation of the system being optimized.

3 THE PROPOSED PROBABILISTIC FORMULATION

For systems with a dynamic, repetitive operation, one can divide the time horizon in discrete periods and redefine Problem (1) in terms of the long term average cost of the objective functions. Let X_k, k > 0, denote the operating point in period k and define the following problem in

Note that problem (1) can be viewed as a specialized version of the proposed problem (Eq. (4)), with X_k = x, ∀k > 0.

In the proposed approach, we assume that at each discrete period k > 0, the system operates in a given point x ∈ Ω with probability µ(x) > [0, 1]. In such a setting, X_k, k > 0 becomes an homogeneous Markov chain (BRÉMAUD, 1999) with transition matrix P = p_xy, with p_xy = µ(y), ∀x, y ∈ Ω and problem (4) becomes

where the last equality is a consequence of the strong law of large numbers (BRÉMAUD, 1999). The last equality also follows by noting that µ: Ω → [0, 1] is a stationary distribution for process X_k, k > 0, e.g. (BRÉMAUD, 1999, Theorem 3.4.1).

Note that, although Problem (5) can be viewed as a generalization of Problem (1), once the operating points are set, it becomes a linear programming problem, regardless of the definitions of the objective functions ƒ_i, i ∈ {1, ..., m}. In contrast, the latter (static) problem is only linear when all objective functions are linear. Hence, the solution to the generalized Problem (5) can be, in fact, easier to find than the solution to Problem (1).

This paper addresses systems with static or slowly changing decision making criteria. Such systems yield long time intervals of operation under the same criterion. In that scenario, we can assume that the length of the operating window under a given decision criterion, denoted here as K, tends to infinity. At each window, one is faced with an static problem of selecting a probability distribution µ in the set Π of feasible distributions from Ω to the interval [0, 1] and the objective is to find the set of non-dominated distributions with respect to the vector of average objectives

where F = [ƒ₁ ƒ₂ ... ƒ_n]^T is the vector of objective functions; the superscript T denotes the transpose of a vector.

Note that, regardless of the functional form of the objectives ƒ_i: Ω → , the expected objective function F: Π → is linear. Therefore, the proposed probabilistic approach can significantly simplify the original problem depending on the functional forms of the objective functions ƒ_i.

3.1 Selection of a Finite Set of Operating Points

Within the general proposed framework, one applies Eq. (6) to find a probability distribution that shall determine the probability of selecting each x ∈ Ω as the operating point at a given period. An obvious simplification is to select a finite set of operating points in O_p ⊂ Ω and make the system operate exclusively within this set. This is equivalent to setting µ(x) = 0, ∀x ∉ O_p. It opens the possibility of selecting operating points within the static Pareto-optimal set. It is this approach that we will further exploit in the next section.

Note that the proposed approach is founded in the repetitive operation of a given system. The operation of the system can be divided in discrete periods and, at the onset of each period, the decision maker chooses the operation point to be utilized in that period. Any operating point x ∈ O_p can be selected with probability µ(x) > 0. At each new period this process is repeated, and the long term average operating cost of the system can be assessed with the help of Eq. (6), where µ(x) = 0, if x ∉ O_p.

3.2 Motivations and Discussions

To illustrate the proposed approach, consider the following hypothetic two-objective optimization problem with only two feasible solutions, as depicted in Fig. 1.

Using the traditional static MOP approach, the decision maker would choose either a or b as the operating point for the MOP problem and the long term average cost would be either F(a) or F(b). It is worth noting that F(a) = [ƒ₁(a) ƒ₂(a)]^T and F(b) = [ƒ₁(b) ƒ₂(b)]^T are two-dimensional vectors comprising both objective functions.

By setting a probability distribution µ: {a, b} → [0, 1] and applying formulation (5), one can emulate a relaxed version of the problem with an expanded operating region. While the static formulation must choose to operate either at point a, paying a long term average cost F(a) or at point b, paying F(b), the dynamic formulation allows one to select any convex combination of F(a) and F(b) as the long time operating cost of the system. It can be argued that each probability distribution generates a virtual operating point c, with average long term cost

F(c) = µ(a)F(a) + µ(b)F(b)

It is worth noting that a real operating point d ∈ Ω, with F(d) = F(c) need not exist. Hence, the proposed approach can, in fact, expand the operating region by generating virtual operating points whose objective function cannot be emulated by any point within the feasible set Ω.

Observe in Fig. 2 that virtual operating points with average operating costs in the line connecting F(a) and F(b) can be generated by setting appropriate probability distributions. The virtual operating point c is only an example corresponding to an specific probability distribution . In Fig. 2, the values ƒ_i(c) = (a) ƒ_i(a) + (b) ƒ_i(b), i = 1, 2. Observe also that the virtual operating region expands the operating region from the discrete set Ω = {a, b} to a continuous set, whose objective function vectors are depicted in Fig. 2.

As illustrated in Fig. 2, one possible application of the proposed approach is to expand the operating region of multi-objective problems with sparse feasible solutions, providing the decision maker with more alternatives to choose from. One obvious example is Integer Multi-Objective Programming (IMOP). By defining a set of feasible operating points (vertices) in Ω and defining a probability distribution in that set, the proposed approach can emulate any long term average objective vector within an hypercube whose vertices comprise the objective vectors at each selected operating point. That may generate a larger operating region than the classical Pareto-optimal set, as illustrated in the example. Observe that the generation of a virtual probabilistic operating point that does not have an equivalent solution in Ω is possible and it does not involve operation outside of Ω.

Another possible application of the proposed formulation is to facilitate the management of the system operation in the long run. With system operation being constrained to a small number of vertex points - such as a few points selected from the Pareto-optimal set, also known as the Pareto front - the setup (transition) costs and times tend to decrease as time elapses, creating a very specialized environment. Moreover, to change virtual operating points one only needs to alter the probability distributions µ: Ω → [0, 1], which does not disturb the operation of the system by adding unexperienced operating points. This can be advantageous, since the Pareto-optimal region often contains an infinite number of points and any change in the operating point within the original MOP framework will very likely involve a transition to an operating point never experienced before. The lack of experience with such a transition may render the operation troublesome and expensive in terms of setup time and cost.

The proposed probabilistic approach can also be employed as a method to manage customer dissatisfaction. It may be used to add flexibility to the operation of the system in the long run and help managing marginal costs which may be difficult to model precisely, such as the cost of customer dissatisfaction. You may alternate solutions, so that the time operating under unattractive solutions is mitigated for selected classes of clients. Such a measure may help maintain clients with conflicting preferences reasonably satisfied with the long term operation of some given service that they share.

As presented, the probabilistic approach we introduce, with the long run average cost as the performance function, generalizes the classical MOP framework and provides new degrees of freedom in the management of the system. Moreover, the original problem becomes a linear MOP problem regardless of the structure of the objective functions. For illustration purposes, we suggest in the next section a specialized version of the approach in which a finite set of operating points is selected a priori.

4 PROPOSED METHOD

The proposed method uses two phases to solve the proposed dynamic multi-objective programming problem. In the first phase, we make use of a conventional technique aimed at dealing with Problem in (1). This technique, which is called weighted global criterion method (MARLER and ARORA, 2004; CHINCHULUUN and PARDALOS, 2007), makes use of the following problem

where w is an m-dimensional vector of weights satisfying w(i) = 1; and F is the m-dimensional objective vector. It builds the Pareto-optimal region as the set of the outcomes of Problem (7) for all possible values of w ∈ ⁿ that satisfy the constraints in (7).

The first phase of the method here proposed solves Problem (7) for K possible values of w ∈ ⁿ, where K > m. This results in K points in the Pareto-optimal set which are then selected as the possible operating points of the dynamic probabilistic MOP formulation (5). These K points add up to form the operating point set O_p, defined in Section 3.1. Hence, the probability of selecting any operating point {x ∈ Ω: x ∉ O_p} is nill.

In order to solve Problem (7) for a given value of w, one needs to find a point x ∈ Ω that satisfies the Karush-Kun-Tucker optimality conditions (BAZARAA et al., 1993; NEMHAUSER and WOSLEY, 1999). The obtained solutions satisfy

Before we proceed to the second phase, we point out that each solution in (8) is an efficient solution to Problem (1). To check that, it suffices to verify that Ω∩Ω _< (x^*) = Ø, where Ω _< (x^*) is defined in (3).

In the second phase, the rationale is to search for solutions whose long term average cost lie in the polytope generated by the set of operating points

where vector F is defined in (1). One way to accomplish that is to solve the problem

where µ_j is the probability that the system is operating in at any given time instant. Note that the set of operating points is S_O = . It is clear to see that Problem (10) is linear and we can use the multiobjective simplex method to solve it, as described in (EHRGOTT et al., 2007).

Remark 1 Note that Problem (10) is always linear and, as a result, convex. Hence, whereas the form of the Pareto front of the static formulation depends on the of the structure of the objective functions and the set Ω in (1), the Pareto front of the proposed formulation is always convex.

The proposed application of the probabilistic approach introduced in Section 3 is related to the classical weighting method to solve MOP problems, e.g. (CHANKONG and HAIMES, 1983; KEENEY and RAIFFA, 1976). The difference is that, whereas the latter method fixes the weights a priori and generates a single mono-objective solution to be used as the operating point (for each possible convex combination of weights, a single mono-objective solution can be generated), the proposed formulation fixes the operating points in Ω a priori and solves for the non-dominated set of probability distributions D = {µ_x, x ∈ O_p}. Note that each probability distribution can be viewed as a convex combination of weights. Thus, any non-dominated solution to Problem (10) assigns weights to each mono-objective solution in (7).

The goal is to minimize the summation

µ_j1(x_j), ƒ₂(x_j), ..., ƒ_n(x_j)] by using a probability distribution function. We call the response virtual solution, because the expected values of the objectives may not correspond to the objective vector of any feasible point x ∈ Ω. That happens because of the well know fact that the expected value of a random variable need not be part of the sample space (ROSS, 2009).

Note that the proposed method is merely one possible implementation of the probabilistic approach introduced in Section 3. It serves to illustrate the benefits and shortcomings of approach. In principle, the choice of the set of operating points O_pcan be extended to the whole set of feasible operating points (Ω) or any subset of it. Hence, the probabilistic approach in Section 3 has the potential to completely emulate (and extend) the classical Pareto-optimal region.

In the numerical examples in the Section 5, the set of vertices were generated as in (9). The experiments illustrate the application of the proposed probabilistic approach in an integer MOP problem.

4.1 Computational Issues

The proposed approach is founded in selecting a finite set of operating points O_p, to be employed for solving (10). The first phase of the method comprising of selecting the set O_p, while the second phase involves finding the Pareto front of Problem (10).

As mentioned in Section 4, each point in the finite set O_pcan be found as a solution to Problem (7): a mono-objective version of the original static formulation in (1). Hence, the computational cost of finding the operating points in O_p is a function of both the cardinality of the set O_p and the computational cost of Problem (7), which depends on the nature of the problem being solved. For example, if Problem (1) happens to be NP-hard, such as a multi-objective traveling salesman problem, then finding O_p also becomes NP-hard. Otherwise, the first phase of the proposed approach can be concluded in polynomial time.

A typical approach to solving (10) in the second phase, here employed, is to find a finite set of points belonging to an arbitrarily fine discretized grid of the Pareto front. Hence, the computational burden depends on the number of points in the selected grid. Observe that, since Problem (10) is linear, each point in its Pareto front can be found in polynomial time, regardless of the nature of the original static problem. Hence, the second phase of the proposed approach can always be concluded in polynomial time.

5 NUMERICAL EXPERIMENTS

Following (FERNÁNDEZ and PUERTO, 2003), we used combinations of mono-objective uncapacitated plant location problems to illustrate the proposed method. The individual problems, which are defined on three different groups , are extracted from the web source (http://www-eio.upc.es/~elena/sscplp/index.html). These groups are divided by number of clients, M, and plants, N. The capacity constraints of the source problems are ignored. We solved three two-objective problems generated as combinations of problems p₁, ..., p₆, each with 20 clients and 10 plants, and two three-objective problems generated as combinations of p₂₆, p₂₇, p₂₈, each with 50 clients and 20 plants, and p₅₀, p₅₁, p₅₂, each with 90 clients and 30 plants.

The points that comprise the Pareto front of each problem solved in this work were obtained by means of the a-posteriori weighted method. Each point corresponds to the solution of a mono-objective problem whose objective is a convex combination of the original objective functions. Every point was found by means of the BINTPROG function, which solves a binary integer programming problem, and is part of the optimization toolbox built in the MATLAB^® 7.8.0 program.

The computational experiments were all performed in a PC with 2.26GHZ Intel^® Core^TM 2 Duo processor, 4GB RAM running Ubuntu 10.10 operational system.

5.1 Formulation of the plant location problem

The multicriteria uncapacitated plant location problem (MU-PLP) is among the most visited binary integer programming problems. Its binary decision variables are

The MUPLP can be formulated in the following way:

where M = {1, ..., m} represents the set of indices for clients and N = {1, ..., n} denotes the set of indices for plants. and , respectively, are the set-up cost and allocation cost for any r ∈ {1, ..., p}, i ∈ M, j ∈ N, and p is the dimension of the objective vector. The constraints are to insure that each client is assigned to only one plant and to guarantee that no client is unattended.

Some elements of the Pareto front of Problem (11) can be selected and the proposed probabilistic formulation becomes:

where µ_j is the probability that the system is operating at some given element , j = 1, ..., K, within the Pareto front. The K elements from the Pareto front can be obtained by means of the a-posteriori weighted method.

5.2 Experimental results

We applied the probabilistic approach to three example problems with two objectives, labeled Example 1 to Example 5. Example 1 is the combination of the instances p₁ and p₂ in the web source; Example 2 combines p₃ and p₄ and Example 3 combines p₅ and p₆ Example 4 combines p₂₆, p₂₇, and p₂₈; Example 5 combines p₅₀, p₅₁, and p₅₂. We applied the method in Section 4 and selected two solutions belonging to the Pareto front, i.e. K = 2, in the first example as the possible operating points for the probabilistic operation. For the second example, three solutions have been selected that belong to the Pareto front, while five solutions from the Pareto front have been selected for the third example.. Finally, we present two problems with three objective functions ant ten solutions have been selected from the Pareto front in both. Figures 3 to ⁷ contain the numerical results.

Note that the extreme points in each figure correspond to the individual solutions to Problems (7) and belong to the original Pareto-optimal region. Observe that, the more evenly distributed the probability distribution is, the closer the corresponding virtual operating point gets to the center of the virtual operating region, which is formed by the convex combination of the K extreme points. Observe also that the virtual operating points between the extreme points are always dominated by some point in the original Pareto-optimal region. However, as we increase K, the Pareto front of the proposed probabilistic approach gets closer to that of the traditional static formulation. One can see in Fig. 5 that the Pareto front of the proposed approach, with K = 5 operating points, almost coincides with the Pareto front of the original formulation.

Finding a function that maps the Pareto front can be very involving. Hence, the classical methods that solve multiobjective problems obtain a limited Pareto-optimal set and it is very sparse sometimes. In such problems, a probabilistic operation can be expand the operating region, providing an infinite number of solutions to choose from, and thus adding flexibility and additional trade-offs for the decision maker.

In Figure 3, two extreme points were selected and the probabilistic operating region is obtained as the set of all possible convex combinations of the extreme points. One can easily see that, in this case, the virtual interpolation front, i.e. the Pareto-optimal region for the probabilistic operation, coincides with the probabilistic operating region itself. In other words, each point in the probabilistic operating region is Pareto-optimal within the proposed formulation. The computational time required to find the static Pareto front was 11.6070s, while the two points of the virtual interpolation front were found in 0.1930s.

For Example 2, we select a middle point in addition to the extreme points in the static Pareto front to compose the operating set O_p. The resulting probabilistic operating region is depicted in Fig. 4. Note that the Pareto front of the proposed formulation gets closer to that of the static formulation. This illustrates that even a small subset of operating points can render the performance of the probabilistic method almost equivalent to the performance of the classical approach. Moreover, as observed above, the range of possibilities for the decision maker to choose from increases dramatically by applying the probabilistic approach. The computational time required to find the static Pareto front was 15.6160s, while the three points of the virtual interpolation front were found in 0.3480s.

For Example 3, we select five points to generate the probabilistic operating region. In this example, the virtual interpolation front is even closer the Pareto front than in Example 2. The results are depicted in Figure 5. Observe that, the more we increase the number of operating points in O_p, the closer the probabilistic Pareto front gets to the static front. Naturally, there is a compromise between the cardinality of the selected set of extreme points and the computational effort needed to identify those points. It is worth noting that the proposed examples suggest that a modest number of operating points can result in a performance which closely approximates that of the original static formulation, while providing extra flexibility and additional trade-offs for the decision making process. The computational time required to the static Pareto front was 9.4340s, while the five points of the virtual interpolation front were found in 0.4830s.

Observe also that one can assign to the probabilistic operating set O_p the whole classical Pareto-optimal set. Hence, it is not difficult to see that the proposed probabilistic approach is a generalization of the static MOP formulation.

For Example 4, we select ten points to generate the probabilistic operating region. The results are depicted in Figure 6. The blue points form Pareto front and they were divided into layers to better represent the frontier uncover the virtual interpolation front, which is green, while the operating region is yellow. Note that the Pareto front of the proposed formulation gets very close to that of the static formulation. The computational time required to find the static Pareto front was 28.7020s, while the ten points of the virtual interpolation were found in 4.4490s.

For Example 5, we select ten points to generate the probabilistic operating region. The results are depicted in ^{Figure 7}. The computational time required to reach find the static Pareto front was 101.8310s, while the ten points of the virtual interpolation were found in17.3220s.

6 CONCLUDING REMARKS

This work presents a novel dynamic approach to a class of static MOP problems. The approach involves a probabilistic selection of the operating point at each discrete operation period. The objective is to minimize the long term average multi-dimensional objective vector. The approach is applied in the case where a set of points in the original MOP Pareto-optimal region is selected, and one is concerned to selecting the Pareto-optimal set of probability distributions over the set of operating points.

The proposed approach can add flexibility to the decision maker by creating a new, expanded set of virtual operating points, each corresponding to a given probability distribution. Moreover, a careful selection of a finite set of operating points can render the probabilistic Pareto-optimal region almost equivalent to the original MOP Pareto front. This limits the number of operating points, rendering the long term operation more efficient and specialized, while also maintaining the performance requirements.

One possible drawback of the method is that, although it can limit the number of operating points to be found in the static Pareto front, the cost of finding each such point still depends on the problem at hand, and one may be forced to seek approximate solutions in the case of NP-hard and/or large-scale problems. Future research directions involve contrasting the proposed method with purely stochastic formulations, such as multi-objective Markov decision processes, as well as studying the effect of the proposed formulation in non-convex MOP problems. A conjecture that remains to be verified is that the probabilistic Pareto front on the latter problems can, in fact, dominate some solutions in the static Pareto front.

Artigo submetido em 21/05/2010 (Id.: 01152)

Revisado em 26/09/2010, 07/02/2011

Aceito sob recomendação do Editor Associado Prof. Ivan Nunes Da Silva

• ABIDO, M. A. (2006). Multiobjective evolutionary algorithms for electric power dispatch problem, IEEE Transactions on Evolutionary Computation 10(3): 315329.
• AMORIM, E. A., ROMERO, R. & MANTOVANI, J. R. S. (2009). Fluxo de potência ótimo descentralizado utilizando algoritmos evolutivos multiobjetivo, Controle & Automação 20: 217232.
• ATHANS, A. & GEERING, H. P. (1973). Necessary and sufficient conditions for differenciable nonscalar-valued functions to attain extrema, IEEE Transaction on Automatic Control 18(2): 132139.
• BAZARAA, M. S., SHERALI, H. D. & SHETTY, C. M. (1993). Nonlinear programming - theory and algorithms, second edn, Jonh Wiley & Sons, New York, USA.
• BRÉMAUD, P. (1999). Gibbs fields, monte carlo simulation, and queues, Springer-Verlag, New York, USA.
• CARDOSO, R. T. N., da CRUZ, A. R., WANNER, E. F. & TAKAHASHI, R. H. C. (2009). Multi-objective evolutionary optimization of biological pest control with impulsive dynamic in soybeans crops, Bulletin of Mathematical Biology 71: 14631481.
• CHANKONG, V. & HAIMES, Y. Y. (1983). Multiobjective decision making: Theory and Methodology, Vol. 8 of North Hollando series in system science and engineering, North Holland, New York, USA.
• CHINCHULUUN, A. & PARDALOS, P. M. (2007). A survey of recent developments in multiobjective optimization, Annals of Operational Research 154: 2950.
• COELLO, C. A. C., VAN VELDHUIZEN, D. A. & LAMONT, G. B. (2002). Evolutionary algorithms for solving multi-objective problems, Kluwer Academic Publishers, New York. ISBN 0-3064-6762-3.
• DEB, K. (2001). Multi-objective optimization using evolutionary algorithms, John Wiley & Sons, LTD, Chichester, UK.
• DIAS, J., CAPTIVO, M. E. & CLÍMACO, J. (2008). A memetic algorithm for multi-objective dynamic location problems, Journal of Global Optimization 42: 221253.
• EHRGOTT, M., PUERTO, J. & RODRÍGUEZ-CHÍA, A. M. (2007). Primal-dual simplex method for multiobjective linear programming, Journal of Optimization Theory and Applications 134(3): 483497.
• ESOGBUE, A. O., SONG, Q. & YOUNG, D. (2006). Non-euler-lagrangian pareto-optimality conditions for dynamic multiple criterion decision problems, Mathematical Methods of Operations Research 63(3): 525542.
• EVERSON, R. M. & FIELDSEND, J. E. (2006). Multiobjective optimization of safety related systems: An application to short-term conflict alert, IEEE Transactions on Evolutionary Computation 10(2): 187198.
• FERNÁNDEZ, E. & PUERTO, J. (2003). Multiobjective solution of the uncapacitated plant location problem, European Journal of Operational Research 145: 509529.
• INSUA, D. R. (1990). Sensitivity analysis in multiobjective decision making, in M. BECKMANN & W. KRELLE (eds), Lecture Notes in Economics and Mathematical Systems, Vol. 347, Springer-Verlag, Berlin, Germany.
• KEENEY, R. L. & RAIFFA, H. (1976). Decisions with mutiples objectives: preferences and calue trade-off, Wiley, New York, USA.
• KIM, I. Y. & de WECK, O. L. (2006). Adaptive weighted sum method for multiobjective optimization: a new method for pareto front generation, Structural and Multidisciplinary Optimization 31(2): 105116.
• MARLER, R. T. & ARORA, J. S. (2004). Survey of multi-objective optimization methods for engineering, Structural and Multidisciplinary Optimization 26(6): 369395.
• NEMHAUSER, G. L. & WOSLEY, L. A. (1999). Integer and combinatorial optimization, Wiley, New York, USA.
• PARETO, V. (1897a). Cours d'economique politique, Vol. I, Macmillan, Paris, FR.
• PARETO, V. (1897b). Le cours d'economique politique, Vol. II, Macmillan, London, UK.
• REID, R. W. & CITRON, S. J. (1971). On noninferior performance index vectors, Journal of Optimization Theory and Applications 7(1): 1128.
• ROSS, S. M. (2009). Introduction to probability and statistics for engineers and scientists, fourth edn, Academic Press, New York, USA.
• SAKAWA, M. & YAUCHI, K. (1998). Coevolutionary genetic algorithms for nonconvex nonlinear programming problems: Revised genocop iii, Cybernetics and Systems 29(8): 885899.
• SALUKVADZE, M. (1974). On the existence of solutions in problems of optimization under vector-valued criteri, Journal of Optimization Theory and Applications 13(2): 203217.
• STADLER, W. (1979). a survey of multicriteria optimization ot the vector maximum problem, part i: 1776-1960, Journal of Optimization Theory and Applications 29(1): 152.
• TRALHÃO, L., COUTINHO-RODRIGUES, J. & ALÇADA-ALMEIDA, L. (2010). A multiobjective modeling approach to locate multi-compartment containers for urban-sorted waste, Waste Management 30: 24182429.
• ZIONTS, S. (1992). Some thoughts on research in mutiple criteria decision making, Computers Operations Research 19(7): 567570.

Publication Dates

• Publication in this collection
  26 Sept 2011
• Date of issue
  Aug 2011

History

• Reviewed
  07 Feb 2011
• Received
  21 May 2010