- Email: [email protected]
Semi-obnoxious single facility location in Euclidean space
Available online at www.sciencedirect.com
Computers & Operations Research 30 (2003) 2191 – 2209
www.elsevier.com/locate/dsw
Semi-obnoxious single facility location in Euclidean space Emanuel Melachrinoudisa;∗ , Zaharias Xanthopulosb a
Department of Mechanical, Industrial and Manufacturing Engineering, Northeastern University, 334 Snell Engineering Center, Boston, MA 02115 5000, USA b KOHLER, Supply Chain Management, Kohler, WI 53044, USA Received 1 June 2001; received in revised form 1 March 2002
Abstract This paper deals with the problem of determining within a bounded region the location for a new facility that serves certain demand points. For that purpose, the facility planners have two objectives. First, they attempt to minimize the undesirable e2ects introduced by the new facility by maximizing its minimum Euclidean distance with respect to all demand points (maximin). Secondly, they want to minimize the total transportation cost from the new facility to the demand points (minisum). Typical examples for such “semi-obnoxious” facilities are power plants that, as polluting agents, are undesirable and should be located far away from demand points, while cost considerations force planners to have the facility in close proximity to the customers. We describe the set of e6cient solutions of this bi-criterion problem and propose an e6cient algorithm for its solution.
Scope and purpose It is becoming increasingly di6cult to site necessary but potentially polluting (semi-obnoxious) facilities such as power plants, chemical plants, waste dumps, airports or train stations. More systematic decision-aid tools are needed to generate several options that balance the public’s concerns with the interests of the developer or location planner. In this paper, a model is presented that generates the best possible sites (e6cient solutions) with respect to two con9icting criteria: maximize distance from population centers and minimize total transportation costs. Having all e6cient solutions at hand, the two sides can select one that best compromises their criteria. A very interesting property found is that most of these e6cient solutions are on edges of a Voronoi diagram. An algorithm is developed that constructs the complete trajectory of e6cient solutions. ? 2003 Elsevier Ltd. All rights reserved. Keywords: Location; Semi-obnoxious facility; Multiple criteria; E6cient set
∗
Corresponding author. Tel.: +1-617-373-4850; fax: +1-617-373-2921. E-mail address: [email protected] (E. Melachrinoudis).
0305-0548/03/$ - see front matter ? 2003 Elsevier Ltd. All rights reserved. PII: S 0 3 0 5 - 0 5 4 8 ( 0 2 ) 0 0 1 4 0 - 5
2192
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
1. Introduction With the advent of more stringent environmental standards, the resurgence of environmental interest groups, and a greater awareness of the public on the potential dangers of pollution, it becomes increasingly important that the placement of polluting (semi-obnoxious) facilities undergo stringent investigation and analysis. Typical examples are the selection of the location for a nuclear plant, a chemical plant, a waste dump or an airport. Ideally, these facilities should be located as far as possible from residential areas, ecosystems, water reservoirs, etc., in order to minimize the undesirable e2ects introduced into the environment by their presence. However, such a location might not be acceptable due to the high transportation costs that increase with distance from the existing facilities. This means that the farthest location found may never be built due to the unreasonably high transportation costs. Therefore, these two con9icting objectives need to be considered simultaneously. This paper addresses the problem of determining within a bounded region the location for a single semi-obnoxious facility using two objectives. The Hrst objective is to minimize the unpleasant e2ects introduced by the new facility, expressed as maximizing the minimum Euclidean distance between the new and existing facilities (maximin). The second objective is to minimize the total transportation cost between the new facility and the existing facilities, expressed as a sum of weighted Euclidean distances (minisum). For certain applications the most suitable distance metric appears to be the Euclidean one. For example, in the case of a power plant, power lines usually follow Euclidean paths and pollution is carried on a straight line in the wind’s direction away from the power plant. Assuming that there are no prevailing wind directions, the risk imposed on a demand point decreases with its Euclidean distance from the power plant. Other polluting agents, such as radiation leaking from a nuclear plant, deHnitely follows Euclidean paths that are omni-directional rays emitted from the facility. Another popular distance metric in the location literature is the rectilinear one. This rectilinear distance metric is more appropriate to use in settings where travel between facilities occurs on network of streets or aisles, such that each street, or aisle, is parallel to two rectangular directions. As it was justiHed above, the Euclidean is the distance metric of choice in this paper. For this bi-objective problem, properties of e6cient solutions are established and an algorithm is developed to construct the complete set of e6cient solutions. Location planners can then choose a solution that best compromises the two con9icting objectives. The generation of many “good” choices that include the ones, which optimize each individual objective, is deHnitely a strength of the bi-criterion model over the single criteria ones. In addition, the maximin criterion alone, which has been traditionally used for undesirable facility location, generates optimal solutions that are very often at the boundary of the feasible region (see Hansen [1]). The latter is certainly not desirable when on the other side of the boundary there are population centers belonging to a di2erent administrative region. Another location in this case which is a little closer to the demand points but which o2ers a reduction in transportation costs may be more appealing to a location planner. When multiple agents exhibiting con9icting objectives participate in the locational decision process, such as in the case of a utility company and the public, multiobjective approaches are more appropriate than single-objective ones for analyzing location problems. See also Cohon [2], Current et al. [3], Avella et al. [4] and Xanthopulos et al. [5] for support of this argument. The rest of this paper is organized as follows: The bi-criterion problem is formulated in Section 2. In Section 3, methods for solving each single objective problem are reviewed and a method for
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2193
generating e6cient solutions of the bi-criterion problem is described. Properties of the e6cient set are also developed in Section 3. Section 4 includes the algorithm for constructing the e6cient set, its computational complexity and an example problem that illustrates the steps of the algorithm and the trajectory shape of the e6cient set. The paper ends with a summary and conclusion in Section 5. 2. Problem formulation Prior to the formulation of the bi-objective model, the two individual objective problems are being reviewed. The maximin objective seeks a location X for an undesirable facility as far away as possible from the closest existing facility: Q1 : max F1 (X ) = min X − Pi ; X ∈S
16i6n
where S is a bounded convex set in R2 , Pi is a demand point in S; i = 1; : : : ; n, and X − Pi is the Euclidean distance between X and Pi . S will be approximated by a convex polygon having m edges, Bj Bj+1 ; j = 1; : : : ; m, where Bm+1 ≡ B1 . Equivalently, S = {X ∈ R2 |hj (X ) = Aj X − cj 6 0; j = 1; : : : ; m}, where hj (X ) = Aj · X − cj 6 0 is the linear constraint associated with the jth edge. Shamos and Hoey [6] solved Problem Q1 for S being the convex hull of the demand points Pi by utilizing the Voronoi diagram (see deHnition of the Voronoi diagram at the beginning of the next section). Dasarathy and White [7] extended their algorithm for a region S that is a convex polygon, i.e., Problem Q1 . Melachrinoudis and Cullinane [8] used a combinatorial algorithm to solve the weighted version of Q1 and Melachrinoudis and Smith [9] improved that algorithm by utilizing a Voronoi structure. The second objective that is considered in this paper is the Euclidean distance minisum objective, commonly known as the Weber problem. According to this problem, a point X is to be found on the plane that minimizes the sum of weighted Euclidean distances from itself to n Hxed points Pi . This problem can be used to Hnd the location of a production facility whose products are to be distributed to n customers. The objective is to minimize the total transportation cost for delivering the products to the customers (demand points): n Q2 : min F2 (X ) = vi X − Pi ; X ∈S
i=1
where vi represents the transportation cost per unit distance between the new location X and the demand point Pi ; i = 1; : : : ; n. There are many early notable contributions in both solution procedures and theory for Q2 , such as Weiszfeld [10], Miehle [11], Kuhn and Kuenne [12], Witzgall [13], Love [14] and Katz [15]. An optimal solution of Q2 can be found by the well-known iterative procedure of Weiszfeld [10]. This procedure however has a 9aw. Since F2 (X ) is nondi2erentiable at demand points, Kuhn [16] showed that the procedure fails if the iteration falls on a demand point. Several variations of the procedure have been proposed to correct this problem. For example, Ostresh [17] deHned a slightly di2erent step when iteration falls on a demand point and Wesolowsky and Love [18] redeHned the Euclidean
2194
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
distance by introducing a hyperbolic approximation. Other variations of this procedure were targeted in speeding up convergence, which is generally slow in the vicinity of demand points. Finally, many other methods di2erent from the Weiszfeld procedure have been proposed in the last 20 years. For a comprehensive review of the history, theory and solution procedures of the minisum problem, see Wesolowsky [19]. The problem under consideration is to Hnd the set of all e6cient solutions of the bi-objective Problem B: n B: max F1 (X ) = min X − Pi ; min F2 (X ) = vi X − Pi : X ∈S
16i6n
X ∈S
i=1
There are only three papers that consider multiple criteria in conjunction with a continuous solution space; two of them deal with rectilinear distances. Mehrez et al. [20] considered a weighted objective function, which combines both minimax and maximin criteria. Melachrinoudis [21] used Fourier– Motzkin elimination to solve the bi-objective minisum/maximin model, taking advantage of the small size of the resulting linear bi-objective problems. Such problems are easier to solve due to the piecewise linearity of the objective functions, while problems involving Euclidean distances such as the one undertaken in this paper cannot be decomposed into linear programs. Therefore, the solution methodology developed in this paper is completely di2erent from the ones developed in the above two papers. It is based on the Karush–Kuhn–Tucker (K–K–T) conditions for nonlinear programming and on Voronoi diagrams. The third paper by Brimberg and Juel [22] deals with Euclidean distances. The Hrst criterion (minisum) is the same as the one considered in this study. The second criterion however is di2erent, being a minisum of negative powers of distances, representing an aggregate undesirable e2ect. Not only the combination of criteria but the solution procedure is also di2erent as well from the one considered in this paper. E6cient points are found using a numerical procedure that solves a parametric system of two di2erential equations. The authors note the di6culties encountered when the e6cient trajectory is discontinuous. In this paper, we deal directly with discontinuities in the e6cient set by determining the endpoints of e6cient segments based on their mathematical properties. In the discrete solution space, there has been an increasing e2ort dealing with the multiobjective location problem of siting multiple undesirable facilities, such as in Ratick and White [23], Erkut and Neuman [24], Current and Ratick [25] and Melachrinoudis et al. [26]. 3. Model development Dasarathy and White [7] solved Problem Q1 using Voronoi polygons. Given a set of n demand points, the Voronoi polygon VPi associated with the demand point Pi is the set of all points X ∈ R2 such that X − Pi 6 min16l6n X − Pl . The Voronoi polygons VPi ; i = 1; : : : ; n, constitute a partition of R2 . The interior points of a Voronoi polygon VPi are closer to Pi than to any other demand point Pj = Pi . The Voronoi edge associated with Voronoi polygons VPi , VPj is on the perpendicular bisector of Pi Pj , i.e., all its points are equidistant from Pi and Pj . A Voronoi vertex is a common point of at least three adjacent Voronoi polygons VPi , VPj , VPk and is therefore equidistant from at least three demand points Pi , Pj and Pk . A Voronoi diagram is deHned as the set of Voronoi edges and Voronoi vertices. For a detailed coverage of Voronoi diagrams see Okabe et al. [27].
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2195
The optimal solution of Problem Q1 can be found by constructing the Voronoi diagram and enumerating a set of points as follows: Theorem 1. Consider the set U consisting of vertices of the Voronoi diagram, intersection points of the edges of the Voronoi diagram with the edges of S, and vertices of S. The optimal solution of Problem Q1 is in the set U. Proof of this theorem as well as a procedure for enumerating the elements of U to obtain the optimal solution of Problem Q1 can be found in Dasarathy and White [7]. Problem Q1 , being nonconvex, is likely to contain many local maxima. As it will be shown later in this paper, local maxima of Problem Q1 are important in identifying discontinuities in the e6cient set of Problem B. The following corollary screens the elements of the set U for local maxima. Corollary 1. (a) A Voronoi vertex belonging to Voronoi polygons VPi , VPj , and VPk is a local maximum of Problem Q1 if it is in the interior of the triangle Pi Pj Pk . (b) An intersection point of a Voronoi edge (equidistant from Pi ; Pj ) with an edge Bj Bj+1 of S is a local maximum of Problem Q1 if it is in the interior of the line segment whose endpoints are the perpendicular projections of Pi and Pj onto the line Bj Bj+1 . (c) Let Pi be the closest demand point to a vertex Bj of S and let Dj and Dj+1 be the perpendicular projections of Pi on the lines that pass through the edges of S that are adjacent to Bj . Bj is a local maximum of Problem Q1 , if Dj and Dj+1 do not both lie on the same one of the two closed half-spaces in R2 , generated by the line Pi Bj . For proof, see Melachrinoudis and Cullinane [8]. Based on Theorem 1 and Corollary 1, the set of all local maxima of Q1 can be constructed. Problem Q2 is a convex problem and its optimum occurs within the convex hull of Pi ; i = 1; : : : ; n, which is a subset of S. Therefore, Q2 can be treated as the classical unconstrained minisum problem. Its solution is found by setting the gradient of the objective function equal to the zero vector: n vi (X − Pi ) C0 : ∇F2 (X ) = = 0: X − Pi i=1 As it was pointed out in the previous section, C0 is solved by one of the variations of the Weiszfeld’s iterative procedures that resolve the issue of nonexistence of ∇F2 (X ) at demand points Pi ; i = 1; : : : ; n. Before we proceed with the generation of e6cient solutions we give the deHnition of e6ciency, nondominance and nondominated tradeo2 rate from Steuer [28] after being adapted for Problem B. Whereas S denotes the feasible region in decision space, C denotes the feasible region in criterion space, i.e., C = {Z = (z; w) | z = F1 (X ); w = F2 (X ); X ∈ S}: De nition 1. A point Xe ∈ S is e6cient if and only if there does not exist another X ∈ S such that F1 (X ) ¿ F1 (Xe ), F2 (X ) 6 F2 (Xe ) and (F1 (X ); F2 (X )) = (F1 (Xe ); F2 (Xe )).
2196
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
De nition 2. A criterion vector Ze = (ze ; we ); Ze ∈ C, is nondominated if and only if there does not exist another criterion vector, Z = (z; w); Z ∈ C, such that z ¿ ze , w 6 we and Z = Ze . Let Ef denote the set of all e6cient points and Nd the set of all nondominated criterion vectors for Problem B. Given the e6cient set Ef , the nondominated tradeo2 rate at each e6cient point Xe provides useful information to the decision maker in making his selection of the most-preferred or best-compromise point from Ef . De nition 3. The nondominated tradeo2 rate between objective z = F1 (X ) and w = F2 (X ) at an e6cient point Xe is d z : R(Xe ) = dw Xe ∈Ef The above rate represents the amount of improvement gained in the maximin objective by selecting another e6cient solution X over Xe in which the transportation cost, represented by the minisum objective, increases by a unit. Equivalently, if the location planner is willing to sacriHce $1000 (money unit) in total transportation costs, the next good choice for location (e6cient point) will be further away from the closest existing facility, as compared to the reference point Xe , by a distance equal to the nondominated tradeo2 rate at Xe . The best-compromise solution is the e6cient solution for which the decision maker’s desirable tradeo2 rate between objectives is equal to the nondominated tradeo2 rate for each pair of objectives (see Cohon [2]). By optimizing each one of the two objectives, we can obtain two e6cient solutions of Problem B. If XU is a unique optimal solution of Problem Q1 , it is an e6cient solution of Problem B. Otherwise, denote by T1 the optimal set and select as XU the one that solves minX ∈T1 F2 (X ). Similarly, if XL is a unique optimal solution of Problem Q2 , it is an e6cient solution of Problem B. If the optimal solution of Problem Q2 is not unique, denote by T2 the optimal set and select as XL the one that solves maxX ∈T2 F1 (X ). Clearly, e6cient points XU and XL deHne the range of nondominance, [F1 (XL ); F1 (XU )] for the maximin criterion and [F2 (XL ); F2 (XU )] for the minisum criterion. It is easy to see that any solution X , whose criterion vector components (one of the two or both) are outside the above ranges, is not e6cient for Problem B. Conversely, any e6cient solution of Problem B, Xe , has both criterion vector components within the above two ranges. The only point of the e6cient set that could lie on a demand point is its starting point XL . An e6cient point Xe , other than the starting point, cannot lie on a demand point. Otherwise, its e6ciency is contradicted since F1 (Xe ) = 0, F2 (Xe ) ¿ F2 (XL ) while F1 (XL ) ¿ 0. Therefore, ∇F2 (Xe ) is properly deHned over Ef , except possibly at Xe = XL , whenever XL lies on a demand point. Of course, demand points can be tested a priori for minisum optimality (see Wesolowsky, [19]). A variation of Weiszfeld’s procedure is needed only when none of the demand points is minisum optimum. Consider the constrained problem Q3 (b), F2 (XL ) 6 b 6 F2 (XU ), formed from Problem B by converting the minisum objective to a minisum constraint and introducing a new variable z to simplify
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2197
the objective function of the maximin objective: Q3 (b):
max
F1 (X ) = z
Subject to
F2 (X ) 6 b; gi (X; z) = z 2 − X − Pi 2 6 0
i = 1; : : : ; n;
hj (X ) = Aj X − cj 6 0
j = 1; : : : ; m:
In the presence of alternative optima, Problem Q3 (b) returns the one that optimizes the second objective. A point Xe is an e6cient solution of Problem B if and only if Xe is an optimal solution of the parametric Problem Q3 (b) for some b, F2 (XL ) 6 b 6 F2 (XU ). As b increases within that range, Q3 (b) generates e6cient points whose criterion vectors increase in both components. The constrained approach for generating e6cient solutions in a multiobjective problem was formally introduced by Haimes et al. [29] and was Hrst applied to water resources problems by Cohon and Marks [30] and Haimes [31]. An alternative approach for generating e6cient solutions is the well-known weighting technique that Hrst appeared in Zadeh [32]. However, as Cohon and Marks [30] point out, when the e6cient set is nonconvex the weighting technique may fail. Q3 (b) is a mathematical programming problem with similar properties to Problem Q1 , i.e., it is nonlinear and nonconvex, complicated further by the nonlinear constraint F2 (X ) 6 b. An optimal solution always exists because the last set of constraints deHnes a bounded feasible region for X and X − Pi is therefore bounded from above, for all Pi . Problem Q3 (b) may have many local maxima, out of which the global maximum has to be identiHed. Instead of directly solving Q3 (b) for several values of b to approximate the e6cient set by a discrete subset, an algorithm is developed in this paper that generates the complete e6cient set. But Hrst the loci of e6cient points are identiHed and their properties are established based on the following lemma and corollaries. Thus Q3 (b) is used in developing the properties of the efHcient segments that constitute the e6cient set. In addition to the feasibility conditions of Q3 (b), an optimal point Xe should satisfy the following optimality and complementary slackness K–K–T conditions: C1 (b):
∇F1 (Xe ) − u · ∇F2 (Xe ) −
n
!i ∇gi (Xe ; ze ) −
i=1
u · [F2 (Xe ) − b] = 0; !i gi (Xe ; ze ) = 0; "j hj (Xe ) = 0;
m
"j ∇hj (Xe ) = 0;
j=1
u ¿ 0;
!i ¿ 0; i = 1; : : : ; n; "j ¿ 0; j = 1; : : : ; m:
Lemma 1. At an optimal point Xe of Problem Q3 (b), (a) gi (Xe ; ze )=0 for at least one i; i=1; : : : ; n, and (b) Either: (i) F2 (Xe ) = b, and hj (Xe ) = 0, j = 1; : : : ; m, or
2198
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
(ii) F2 (Xe ) = b, and hj (Xe ) = 0 for at least one j, or (iii) F2 (Xe ) ¡ b, which further implies that Xe is a local maximum of Problem Q1 . See proof of Lemma 1 in Appendix A. 3.1. Discontinuities in the nondominated set In a linear multiobjective problem, the e6cient set Ef is always connected. In a nonlinear multiobjective problem, however, Ef may not be connected (Steuer [28]). The same holds for the nondominated set Nd . We expect therefore discontinuities in the nondominated and e6cient sets of Problem B. We start with discontinuities in Nd . To elaborate, observe the last case (iii) of Lemma 1, i.e., if an e6cient point Xe solves Q3 (b) and the minisum constraint is not binding at Xe , then Xe is a local maximin point. In this case, a very small increase in b will still make Xe optimal for Q3 (b). This implies that the nondominated set will be discontinuous at the criterion vector of Xe . The following two corollaries deal with discontinuity in the nondominated set. Corollary 2. Let Xe be an optimal solution of Problem Q3 (b) for some b, F2 (XL ) 6 b 6 F2 (XU ), and Ze be its criterion vector. (a) If F2 (Xe ) ¡ b, then there exist no e>cient solutions X such that F2 (Xe ) ¡ F2 (X ) 6 b, i.e., the nondominated set Nd is discontinuous at Ze . (b) If there are no nondominated criterion vectors Z(z1 ; z2 ) such that F2 (Xe ) = b1 ¡ z2 6 b, then Xe is a local maximin point. See proof of Corollary 2 in Appendix A. Let bl and bu be the lowest and highest value of b, respectively, for which the optimal solution of Q3 (b) is a local maximin point, Xe1 . Then according to Corollary 2, the nondominated set Nd is discontinuous at Ze1 = (F1 (Xe1 ); bl ) and there are no nondominated criterion vectors Z = (F1 (Xe1 ); b) for bl ¡ b 6 bu . The next corollary complements Corollary 2 by establishing continuity in the nondominated set in terms of the binding status of the minisum constraint of Q3 (b). Corollary 3. If for every b in the range [b1 ; b2 ] the optimal solution of Problem Q3 (b), Xe , is such that F2 (Xe ) = b, the subset of nondominated criterion vectors Z(z1 ; z2 ) for z2 ∈ [b1 ; b2 ] comprises a connected subset of Nd (continuous nondominated frontier). See proof of Corollary 3 in Appendix A. A geometrical interpretation of the role local maximin points play in generating discontinuities follows: The minisum constraint in Q3 (b) restricts X to be within the minisum contour F2 (X ) = b. Assume that for some value of b; Q3 (b) returns a local maximin point Xe . When b increases further, the minisum contour expands providing more choices for X . The new choices of X however are not better since Xe is the best point in its local neighborhood. 3.2. Characterization and properties of e>cient segments We investigate now the consequences of case (i) and (ii) of part b of Lemma 1. Let us start with case (i) by assuming that F2 (Xe ) = b and gk (Xe ; ze ) = 0, the latter implying that Xe is closest to Pk
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2199
than to any other demand point. In terms of the K–K–T conditions, C1 (b), this means that u ¿ 0, !i = 0 i = k, "j = 0 for all j and the optimality condition for its X component can be simpliHed as, −u∇F2 (Xe ) + 2!k (Xe − Pk ) = 0. Setting !k = 2!k =u, the complete set of K–K–T conditions becomes C2 (b):
!k (Xe − Pk ) = ∇F2 (Xe ); F2 (Xe ) = b; gk (Xe ; ze ) = 0; !k ¿ 0:
With the exception of XL , in the case it falls on a demand point Pk , all points of the Hrst e6cient segment satisfy conditions C2 (b). This Hrst e6cient segment starts at XL and terminates at some point Y which is either on a Voronoi edge or on the boundary of S. A method for determining Y is developed Hrst and a trajectory analysis, as an alternative to conditions C2 (b), is presented next to identify segment XL Y . Let VPk be the Voronoi polygon containing the minisum point XL and let Vk = VPk ∩ S be its feasible subset. Vk is bounded by a set of edges Dk , which are either Voronoi edges or line segments of the boundary of S. Let XL Y = Vk ∩ Ef be the Hrst e6cient segment of Ef . This segment will contain e6cient points that are closer to Pk than to any other demand point. Since Y lies either on a Voronoi edge or on the boundary of S, it can be found as follows: For each edge Ek E‘ ∈ Dk , the following system of equations, based on conditions C2 (b), are solved: C3 :
!k (X − Pk ) = ∇F2 (X ); X ∈ E k E‘ ; !k ¿ 0:
Among all the solutions obtained, the one that minimizes F2 (X ) is the endpoint Y of the Hrst e6cient segment. Given a value for b, conditions C2 (b) provide a system of three nonlinear equations in three unknowns to solve for Xe . However, an easier way to Hnd Xe is described next, based on the convex combination of the two criteria: minX F3 (X ; !) = −!F1 (X ) + (1 − !)F2 (X ). For X ∈ Vk , this problem can be simpliHed substantially by setting F1 (X ) = X − Pk . It subsequently becomes similar to the minisum Problem Q2 , except that its weights have been parameterized using !: Q4 (!): min F3 (X ; !) = (1 − !)vi X − Pi + [(1 − !)vk − !]X − Pk : X
i=k
Q4 (!) has to be solved for ! ∈ (0; !∗ ] to generate points of the Hrst e6cient segment, where is found as the solution of anyone of the two linear equations from (C0 ), after setting X = Y (already computed) and using the parameterized weights from Q4 (!). For practical reasons, only a few dispersed e6cient points on segment XL Y can be generated and presented to the decision maker. Consider now case (i) of Lemma 1 in which at least two of the constraints gi (X; z) 6 0 are binding. It is obvious that if exactly two of these constraints are binding Xe is on a Voronoi edge. If three or more are binding, Xe is a Voronoi vertex. Case (ii) of Lemma 1 indicates that e6cient points may also occur on the edges and vertices of S, whenever one or two of the constraints hj (X ) 6 0 are binding at Xe . !∗
2200
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
Concluding the discussion of Lemma 1 and the subsequent analysis of K–K–T, we can characterize the components of the e6cient set as either (a) points satisfying conditions C2 (b), or equivalently solutions of Q4 (!), or (b) line segments on Voronoi edges, or (c) line segments of the boundary of S. The above two corollaries state that discontinuities in the nondominated set occur only at criterion vectors corresponding to some local maximin points. Discontinuities in the e6cient set however may occur at other e6cient points as well, that are not local maximin points. Corollary 4 establishes properties of e6cient points at which discontinuities occur. Corollary 4. Let Xe1 be the endpoint of an e>cient segment of Ef , and X 2 be the starting point of the next e>cient segment as the two objectives increase. Then, either (i) Xe1 is a local maximin point or (ii) Xe1 and X 2 are alternative optima of Problem Q3 (b) for some b. See proof of Corollary 4 in Appendix A. Based on Corollary 4 we develop procedures next for identifying endpoints of e6cient segments. For example, given the starting point of the current e6cient segment, and the line on which it lies, such as a Voronoi edge, what we need to determine is the endpoint of that e6cient segment and the beginning of the next one as b increases. Then we can repeat this step until the global maximin point is encountered. These procedures will be useful in constructing e6cient segments. Without loss of generality, let us assume that the current e6cient segment is on a Voronoi edge Ek E‘ starting at Y, and the direction of increase for both objectives is from Y toward Ek . According to Corollary 4 case (ii), if there is a discontinuity in the e6cient set at some point Xk , Xk ∈ YEk ; Xk = Ek; and the starting point of the next e6cient segment is Xq , Xq ∈ Eq Es , then Xk and Xq are alternative optima of Q3 (b) for some b; F2 (Y ) 6 b ¡ F2 (Ek ). In order to determine points Xk and Xq , therefore, the following system of equations should be solved for every Voronoi edge or edge of S, Eq Es , where Eq Es ∈ Vq for q = 1; : : : ; n: C4 :
Xk − Pk = Xq − Pq ; F2 (Xk ) = F2 (Xq ); Xk ∈ YEk ; X q ∈ E q Es :
Consider two points Xk ∈ YEk and Xq ∈ Eq Es such that F2 (Xk )=F2 (Xk )−Ub and F2 (Xq )=F2 (Xq )−Ub for a small Ub. Since Xk is e6cient, F1 (Xk ) ¿ F1 (Xq ). F1 (Xk ) = F1 (Xq ) implies that, F1 (Xk ) − F1 (Xk )=Ub ¡ F1 (Xq ) − F1 (Xq )=Ub, i.e., the tradeo2 rate between objective F1 (X ) and F2 (X ) on Eq Es at Xq is greater than the tradeo2 rate on Ek E‘ at Xk . In addition, since the tradeo2 rate is continuous on line segments (see Appendix A) the optimal solution of Q3 (b + Ub) will be on Eq Es . If there is no solution to C4 for any edge Eq Es , then the complete segment YEk is e6cient. If in addition, Ek is a local maximin point, then Ek is the endpoint of the e6cient segment YEk according to Corollary 4 case (i). Theorem 1 in conjunction with Corollary 1 can be used to determine if Ek is a local maximin point. In order to determine the starting point of the next e6cient segment, the following system of equations should be solved for every Voronoi edge or edge of S; Eq Es , where
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2201
Eq Es ∈ Vq for q = 1; : : : ; n: C5 : X − Pq = F1 (Ek ); X ∈ Eq E s ; F2 (X ) ¿ F2 (Ek ): If G is the set of solutions of C5 , an optimal solution Y from G is selected, such that F2 (Y ) = minX ∈G F2 (X ). Although Y is dominated by Ek , a point at an inHnitesimal distance from Y lying on Eq Es in a direction of increase for both objectives will be e6cient, and the new e6cient segment lies on Eq Es . 4. Construction of the e#cient set Rather than solving Q3 (b) for several values of b to generate a discrete subset of the e6cient set, the following algorithm generates a complete trajectory of the e6cient set. Based on the characterization of e6cient segments and the properties of their endpoints, established in the previous section, the algorithm identiHes each e6cient segment by its endpoints and the equations of the line on which it lies, as parameter b scans the range [F2 (XL ); F2 (XU )]. Step 1 of the algorithm determines the optimal solutions XU and XL for problems Q1 , and Q2 , respectively. XL is the starting point of the trajectory of the e6cient set. Step 2 constructs the Voronoi diagram and partitions S into n polygons. Step 3 identiHes the polygon Vk containing XL and determines the smallest value of b for which the curve F2 (X ) = b is tangent to an edge of the polygon Vk . The tangent point, Y, and point XL are the endpoints of the Hrst e6cient segment of Ef . The points of that segment are solutions of the system of Eqs. C2 (b) or equivalently solutions of Q4 (!). Step 4 Hnds e6cient segments on Voronoi edges and edges of S by identifying endpoints at which discontinuities in the e6cient set occur that are caused by alternative optima in Problem Q3 (b). Step 5 Hnds e6cient segments by identifying endpoints, which are the local maximin points. The algorithm terminates when the endpoint of the e6cient set is encountered, XU . 4.1. Algorithm Step 1. Determine the maximin (XU ) and minisum (XL ) points by solving Problem Q1 and Q2 , respectively. Compute the nondominance range for each objective, i.e., [F1 (XL ); F1 (XU )] and [F2 (XL ); F2 (XU )]. Step 2. Construct the Voronoi diagram and Hnd its intersections with the edges of S. The Voronoi diagram partitions S into polygons. Each polygon Vi contains a point Pi ; i = 1; : : : ; n and is bounded by a set of edges Di ; i = 1; : : : ; n which are either Voronoi edges or edge segments of the boundary of S. Step 3. Identify the demand point Pk corresponding to the polygon Vk , containing XL . Find the smallest value of b for which the solution of Q3 (b) is a point of the boundary of Vk by solving the system of equations C3 for all edges Ek E‘ ∈ Dk . Denote such a point by Y and the corresponding edge of Vk by Ek E‘ . Points of the e6cient segment XL Y are found by solving Problem Q4 (!) for various values of !. Without loss of generality let Ek be the endpoint of Ek E‘ such that movement from Y toward Ek increases the distance from Pk .
2202
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
Step 4.1. For every edge and boundary side Eq Es ∈ Di ; i = 1; : : : ; n, Hnd Xk ∈ YEk and Xq ∈ Eq Es such that Xk − Pk = Xq − Pq
and
F2 (Xk ) = F2 (Xq ):
If there is no solution, the segment YEk is e6cient; go to Step 4.2. Otherwise, if there is a solution (Xk ; Xq ) a new e6cient segment on edge Eq Es starts at Xq . Set Y = Xq , Ek = Eq and E‘ = Es and repeat Step 4.1. Step 4.2. If F1 (Ek ) = F1 (XU ) STOP. Otherwise, go to Step 5. Step 5. If Ek is a local maximin point: Denote by G the set of points X such that X − Pi = F1 (Ek ); X ∈ Eq Es ; F2 (X ) ¿ F2 (Ek );
for Eq Es ∈ Di ; i = 1; : : : ; n:
Let Y ∈ G be the point and Eq Es the corresponding edge such that F2 (Y ) = min F2 (X ): X ∈G
A new e6cient segment starts at Y on edge Eq Es . Set Ek = Eq ; E‘ = Es , and go to Step 4.1. Otherwise, an e6cient segment will start at Ek along the edge emanating from Ek on which X − Pk increases. Let that edge be Ek Er . Set Y = Ek , Ek = Er , and go to Step 4.1. Although unlikely, it is possible that ties occur in certain steps of the above algorithm. In particular, if more than one point Xq is found in Step 4.1 then the point that yields the largest tradeo2 rate on the respective edge Eq Es should be selected. If there is a tie in the tradeo2 rate as well, then multiple e6cient segments start at points Xq and the algorithm should be executed completely for each new segment. A similar tie breaking procedure should be used if more than one point Y is found in Step 5. 4.2. Nondominated tradeo@ rate The nondominated tradeo2 rate between the maximin and minisum objective at an e6cient point, Xe , can be determined using its deHnition given earlier in this paper. The nondominated tradeo2 rate can be reduced to a formula (see Appendix A) for e6cient line segments, lying on edges Eq Es ∈ Vk , i.e., for points on Voronoi edges or edges of S, cos ’ X e ∈ E q Es ; R(Xe ) = ∇F2 (Xe ) · u where u is a unit vector on Eq Es , and ’ is the angle between u and the direction Pk Xe . For points belonging to nonlinear e6cient segments, the nondominated tradeo2 rate can be estimated by generating pairs of e6cient points which are very close to each other and directly computing the ratio of change in the two objective functions, Uz1 =Uz2 . As the e6cient set moves from one segment to another the nondominated tradeo2 rate changes value, being discontinuous at all endpoints of e6cient segments, except for the beginning and end of the e6cient set. At these points, where discontinuity occurs, a left and a right value have to be determined.
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2203
4.3. Complexity analysis of the algorithm Certain steps of the algorithm are bounded by polynomial times while others, such as solving the minisum problem or a system of equations, require numerical methods that converge within certain predetermined accuracy. Steps 1 and 2 are executed once and amount to solving for XL and XU . XL can be found by a variant of Weiszfeld’s procedure. The construction of the Voronoi diagram, its intersection with the edges of S and subsequently determining XU requires O(m log2 n + n log n + n log m) time (see Dasarathy and White, [7]). Since there exist O(n) Voronoi edges and m edges in S; O(m + n) systems C3 of three equations in three unknowns have to be solved in Step 3 to Hnd point Y . Step 4.1 requires the solution of O(n + m) systems of two equations in two unknowns. Finally, the number of executions of iterative Steps 4 and 5 combined is equal to the number of e6cient segments, k, lying on Voronoi edges and edges of S. Although some edges may be repeated, many others may not be included in the trajectory of e6cient set. Our conjecture is that k =O(n+m). The following example has been developed to apply the steps of the algorithm and illustrate the trajectory of the e6cient set, the nondominated set, and the nondominated tradeo2 rate. 4.4. An example problem A 13 by 13-square, S, is shown in Fig. 1. Inside S there are four demand points Pi ; i = 1; 2; 3; 4 for which information is given in Table 1. The Voronoi diagram partitions S into four polygons Vi ; i = 1; 2; 3; 4. Voronoi edges are shown as dashed lines. The set R of all candidates for local maximin points according to Theorem 1 is constructed as R = {Y2 ; Y4 ; Y7 ; Y8 ; Y9 ; Y10 ; Y11 ; Y12 ; Y13 ; Y14 }. Among these candidates, Y2 and Y4 are Voronoi vertices, Y7 ; Y8 ; Y10 and Y11 are intersections of Voronoi edges with edges of S and Y9 ; Y12 ; Y13 and Y14 are the vertices of S. The elements of the set R after being screened for local maxima according to Corollary 1 yield the set of local maximin points {Y2 ; Y9 ; Y10 ; Y11 ; Y12 ; Y14 }. Point Y9 is the global maximin point, XU =Y9 . The solution to the minisum Problem Q2 yields P2 as the minisum point, XL = P2 . The nondominance range is [F1 (XL ); F1 (XU )] = [0; 9:85] and [F2 (XL ); F2 (XU )] = [34:39; 125:53] for the maximin and minisum objectives, respectively. As the minisum objective function increases within its nondominance range, the algorithm provides all endpoints of e6cient segments and their criterion vectors as displayed in Table 2. The trajectory of the e6cient set, Ef , is shown in Fig. 1 with solid line segments. Starting at P2 and ending at Y9 , Ef consists of seven e6cient segments, i.e., Ef = {P2 Y1 ; Y1 Y2 ; Y3 Y4 ; Y4 Y5 ; Y6 Y7 ; Y7 Y8 ; Y8 Y9 }. The Hrst e6cient segment, P2 Y1 , is the set of e6cient points which are the solutions to Q4 (!), for 0 6 ! 6 0:778. Several values of ! in that range were considered and the solution points to Q4 (!) were connected in Fig. 1. The next four segments are parts of Voronoi edges, and the last two segments are parts of edges of S. Similarly, the nondominated set, Nd , has been constructed in Fig. 2, Nd = {ZL Z1 ; Z1 Z2 ; Z3 Z4 ; Z4 Z5 ; Z6 Z7 ; Z7 Z8 ; Z8 ZU }. The e6cient set exhibits both types of discontinuities as deHned in Corollary 4. The Hrst discontinuity occurs at Y2 , a local maximin point. The starting point of the next e6cient segment, Y3 , is determined in step 5 of the algorithm. It should be noted that Y3 is not an e6cient point, i.e., the e6cient segment Y3 Y4 is open at Y3 . Similarly, its criterion vector Z3 is dominated. To exclude these points from their respective e6cient segments, Y3 and Z3 are marked by solid circles in Figs. 1 and 2, respectively. The second discontinuity of Ef occurs at Y5 . Y5 and Y6 , determined in
2204
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209 Y7
Y 9 Y8
Y12
Y6
V3
V2 P3
•
Y1
V1 P2
Y2
•
P1
Y3
•
Y10
•
P4
Y4 Y5
S
V4
Y14
Y11
Y13
Fig. 1. Feasible region and e6cient set. Table 1 Data for the example problem i
1
2
3
4
Pi vi
(4,4) 2
(8,7) 4
(11,10) 3
(13,4) 2
Step 4.1 of the algorithm, are alternative optima of Problem Q3 (b) for b = 79:46. As it is illustrated in Fig. 2, there is a single discontinuity in the nondominated set, occurring at Z2 , which is the criterion vector of the already mentioned local maximin point Y2 . The nondominated vectors of the alternative optima Y5 and Y6 are identical, Z5 ≡ Z6 . The nondominated tradeo2 rate R between the maximin and minisum objectives is shown in Fig. 3 as a function of z1 . The segment RL R1 was constructed using several values of the nondominated tradeo2 rate calculated as the ratio of change in the two objectives for selected pairs of e6cient points. All other segments were constructed by determining the nondominated tradeo2 rate for several points along the e6cient edges using the formula provided in Appendix A. The nondominated tradeo2
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2205
Table 2 Endpoints of e6cient segments and their criterion vectors Endpoints of e6cient segments
Criterion vector
Nondominated tradeo2 rate
P2 (8; 7) = XL Y1 (10:21; 7:79) Y2 (11:25; 6:75) Y3 (9:75; 4:25) Y4 (8:5; 2:17) Y5 (8:5; 0:38) Y6 (5:70; 12:30) Y7 (5; 13) Y8 (0:38; 13) Y9 (0; 13) = XU
ZL (0:0; 34:39) Z1 (2:35; 40:40) Z2 (3:26; 44:85) Z3 (3:26; 48:72) Z4 (4:86; 63:54) Z5 (5:77; 79:46) Z6 (5:77; 79:46) Z7 (6:71; 89:15) Z8 (9:70; 122:35) Z9 (9:85; 125:53)
0.567 0.2937,0.2875 0.1815 0.1011 0.1079,0.1192 0.106094 0.0968 0.0960,0.0801 0.0935,0.0444 0.0480
130
ZU
120
Z8
110 100 90
Z7
z2
80
Z5 = Z6
70 60
Z3
50 40 30
Z2
Z1
ZL
Z4
20 10 0 0
1
2
3
4
5
6
7
z1
8
9
10
11
12
13
Fig. 2. Nondominated set.
rate, as shown in Fig. 3, is discontinuous at endpoints of e6cient segments. For these endpoints the nondominated tradeo2 rates are displayed in the third column of Table 2. Two numbers are listed, a left and a right value separated by a comma, at the endpoints of the e6cient segments where discontinuities occur. The solution of the system of nonlinear equations required in Steps 3, 4.1, and 5 of the algorithm were obtained using the Mathcad software package (Mathcad 3.1 [33]).
2206
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
0.6 RL
0.5
R
0.4
R1
0.3
0.2
0.1 RU 0 0
1
2
3
4
5
z1
6
7
8
9
10
Fig. 3. Nondominated tradeo2 rate.
5. Conclusion In many instances, the simultaneous consideration of more than one objective improves the realism of the mathematical model and provides the decision maker with more solution options in the form of the e6cient set. This is especially true for the location of a semi-obnoxious facility such as a nuclear plant, a chemical plant, or a waste dump. For these location decisions transportation cost is very important to be ignored and a single maximin objective may be inadequate. In addition, since very often the maximin optimal solution is found at the boundary of the feasible region, the generation of additional solution options in the interior of the feasible region is very desirable. The Euclidean distance maximin–minisum location model is a nonlinear and nonconvex bi-objective mathematical programming problem. The constrained approach is used to convert it to a single objective parametric nonlinear problem. Properties of e6cient solutions are established using the K–K–T conditions. It was shown that the e6cient set consists of segments not necessarily connected which are either (a) loci of points that maximize distance from a demand point to minisum contours, or (b) parts of Voronoi edges, or (c) parts of the edges of the feasible region. A complete trajectory of the e6cient set is constructed by determining the endpoints of e6cient segments. Despite the above merits, the presented model for semi-obnoxious facility location points to a number of directions for future work. The computational complexity of the developed algorithm has been delineated; however, extensive computational experiments can be performed to determine the algorithm’s average complexity as the number of existing facilities increase. In addition, future
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2207
research should point to the direction of multifacility semi-obnoxious facility location under multiple objectives. Such a problem is expected to be much more complex than its single facility counterpart. Since the (multiple) semi-obnoxious facilities should not be crowded in a single location, additional (separation) constraints have to be imposed. This in turn, will spoil the elegant Voronoi structure of the problem and will not enable the given algorithm to be generalized for multiple facilities. Acknowledgements The authors wish to thank Professor H.A. Eiselt for his valuable comments on a previous draft of this paper as well as the two anonymous referees who provided valuable suggestions for improving the paper. Appendix A. Proof of Lemma 1. (a) The optimality K–K–T condition, solved for its z-component, yields 1 − 2 · n z · !i = 0. Since Problem Q3 (b) has a bounded objective value, at least one of the !i ’s has to be i=1
positive, i.e., gi (Xe ; ze ) = 0 for at least one i. (b) Constraint F2 (X ) 6 b is either binding or nonbinding at Xe . If it is binding, we have case (i) and (ii). If it is nonbinding, then u = 0 in C1 (b) and the above K–K–T conditions reduce to the K–K–T conditions satisHed by a local maximum of Problem Q1 . Therefore, F2 (Xe ) ¡ b implies that Xe is a local maximum of problem Q1 (local maximin point). Proof of Corollary 2. (a) By Lemma 1, F2 (Xe ) ¡ b implies that Xe is a local maximin point. Assume that there exists an e6cient solution X of Problem B such that F2 (Xe ) ¡ F2 (X ) 6 b, i.e., X is a feasible solution of Problem Q3 (b). Since both X and Xe are e6cient, F1 (X ) ¿ F1 (Xe ), thus contradicting the optimality of Xe to Problem Q3 (b). (b) If Xe is not a local maximin point, then there exists a direction of improvement Ux, where Ux is very small and F1 (Xe +Ux)¿F1 (Xe ). In addition, since Xe is optimal for Q3 (b1 ), b1 ¡F2 (Xe + Ux) = b1 + Ub1 ¡ b for a small positive Ub1 and the solution of Q3 (b1 + Ub1 ) will yield e6cient point Xe + Ux. The latter contradicts the nonexistence of nondominated criterion vectors Z(z1 ; z2 ) such that F2 (Xe ) = b1 ¡ z2 = b1 + Ub1 6 b. Proof of Corollary 3. Let Xe and Xe be the optimal solutions of problems Q3 (b) and Q3 (b + -), respectively, for a very small positive - and b1 6 b ¡ b + - 6 b2 . Since both Xe and Xe are e6cient, b + - = F2 (Xe ) ¿ F2 (Xe ) = b implies that F1 (Xe ) ¿ F1 (Xe ). Since the maximum distance of any point Pi to the contour F2 (X ) = b is a continuous function of b, a very small increase - in the RHS of constraint F2 (X ) 6 b will result in a very small increase in the optimal objective function value of Problem Q3 (b). Therefore, for z2 ∈ [b1 ; b2 ] the nondominated set is connected. Proof of Corollary 4. (i) F2 (X ) is either continuous or not over the e6cient set at Xe1 . If it is not continuous, then according to Corollary 2, Xe1 is a local maximin point.
2208
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
(ii) If F2 (X ) is continuous over the e6cient set at Xe1 , then F2 (X 2 ) = F2 (Xe1 ) = b. Since both points are e6cient, this further implies that, F1 (X 2 ) = F1 (Xe1 ), i.e., Xe1 and X 2 are alternative optima of Problem Q3 (b). A.1. Nondominated tradeo@ rate on a line segment Eq Es Let Xe be an e6cient point such that Xe ∈ Eq Es and let Ze (ze ; we ) be its criterion vector, such that ze = Xe − Pk . Let u(cos .; sin .) be the unit vector on Eq Es having a direction from Eq to Es , where . is the angle between vector u and the x-axis. Consider another e6cient point on the line segment at an inHnitesimal distance from Xe , i.e., Xe + u dt having criterion vector Z(z; w). Then after ignoring dt 2 terms we have d(z 2 ) = 2z d z = Xe + u dt − Pk 2 − Xe − Pk 2 = 2(Xe − Pk ) · u dt; or d z (Xe − Pk ) · u = cos ’; where ’ is the angle between u and Pk Xe : = dt Xe Xe − Pk The directional derivative of w = F2 (X ), on the line segment Eq Es at Xe is dw = ∇F2 (Xe ) · u: dt Xe
By dividing the two derivatives we obtain d z cos ’ : R(Xe ) = = dw Xe ∇F2 (Xe )u References [1] Hansen P, Peeters D, Thisse JF. On the location of an obnoxious facility. Sistemi Urbani 1981;3:299–317. [2] Cohon JL. Multiobjective programming and planning. New York: Academic Press, 1978. [3] Current J, ReVelle C, Cohon J. The median shortest path problem: a multiobjective approach to analyze cost vs. accessibility in the design of a transportation network. Transportation Science 1987;21:188–97. [4] Avella P, Benati S, Martinez LC, Dalby K, Girolamo DD, Dimitrijevic B, Ghiani G, Giannikos I, Guttmann N, Hultberg TH, Fliege J, Marin A, Marquez MM, Ndiaye MM, Nickel S, Peeters P, Brito DP, Policastro S, Saldanha de Gama FA, Zidda P. Some personal views on the current state and the future of locational analysis. European Journal of Operational Research, 1998;104:269 –87. [5] Xanthopulos Z, Melachrinoudis E, Solomon M. Interactive multiobjective group decision making with interval parameters. Management Science 2000;46:1585–601. [6] Shamos MI, Hoey D. Closest-point problems. In: Proceedings of the 16th Annual IEEE Symposium on Foundations of Computer Science, 1975. p. 151– 62. [7] Dasarathy B, White LJ. A maximin location problem. Operations Research 1980;28/6:1385–401. [8] Melachrinoudis E, Cullinane TP. Locating an undesirable facility within a geographical region using the maximin criterion. Journal of Regional Science 1985;25/1:115–27. [9] Melachrinoudis E, MacGregor Smith J. An O (mn2 ) algorithm for the maximin problem in E 2 . Operations Research Letters 1995;18:25–30. [10] Weiszfeld E. Sur le Point Lequel la Somme des Distances de n Points Donnes Est Minimum. Tohoku Mathematics Journal 1937;43:355–86. [11] Miehle W. Link-length minimization in networks. Operations Research 1958;6:232–43.
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2209
[12] Kuhn HW, Kuenne E. An e6cient algorithm for the numerical solution of the generalized weber problem in spatial economics. Journal of Regional Science 1962;4(2):21–33. [13] Witzgall C. Optimal location of a central facility: mathematical models and concepts. National Bureau of Standards Report 8388. Gaithersburg, MD: NBS, 1965. [14] Love RF. A note on the convexity of siting deposits. The International Journal of Production Research 1967;6: 153–4. [15] Katz NI. Local convergence in fermat’s problem. Mathematical Programming 1974;6:89–104. [16] Kuhn HW. A note on fermat’s problem. Mathematical Programming 1973;4:98–107. [17] Ostresh Jr. L. On the convergence of a class of iterative methods for solving the weber location problem. Operations Research 1978;26:597–609. [18] Wesolowsky GO, Love RF. A nonlinear approximation method for solving a general rectangular distance weber problem. Management Science 1972;18:656–63. [19] Wesolowsky GO. The weber problem: history and perspectives. Location Science 1993;1/1:5–23. [20] Mehrez A, Sinuany-Stern Z, Stulman A. A single facility location problem with a weighted maximin–minimax rectilinear distance. Computers and Operations Research 1985;12(1):51–60. [21] Melachrinoudis E. Bicriteria location of a semi-obnoxious facility. Computers & Industrial Engineering 1999;37: 581–593. [22] Brimberg J, Juel H. A bicriteria model for locating a semi-desirable facility in the plane. European Journal of Operational Research, 1998;106:144 –51. [23] Ratick SJ, White AC. A risk-sharing model for locating noxious facilities. Environment and Planning B: Planning and Design 1988;15:165–79. [24] Erkut E, Neuman S. A multiobjective model for locating undesirable facilities. Annals of Operational Research 1993;40:209–27. [25] Current J, Ratick S. A model to assess risk, equity and e6ciency in facility location and transportation of hazardous materials. Location Science 1995;3(3):187–201. [26] Melachrinoudis E, Min H, Wu X. A multiobjective model for the dynamic location of landHlls. Location Science 1995;3(3):143–66. [27] Okabe A, Boots B, Sugihara K. Spatial tesselations, concepts and applications of voronoi diagrams. Chichester, England: Wiley, 1992. [28] Steuer ER. Multiple criteria optimization: theory, computation, and application. Malabar, FL: Krieger, 1989. [29] Haimes YY, Wismer DA, Lasdon LS. On Bicriterion formulation of the integrated system identiHcation and system optimization. IEEE Transactions on Systems, Man and Cybernetics 1971;1:296–7. [30] Cohon JL, Marks DH. Multiobjective screening models and water resources investment. Water Resources Research 1973;9(4):826–36. [31] Haimes YY. Multiple objectives in water and related land resources planning. Proceedings Seminar. Colorado River Basin, Modeling Studies, Utah State University, Logan, 1975. p. 55 –70. [32] Zadeh L. Optimality and non–scalar–valued performance criteria. IEEE Transactions on Automatic Control 1963;AC-8/59. [33] Mathcad 3.1. User’s guide windows version. MathSoft Inc., Cambridge, MA, 1992. Emanuel Melachrinoudis is an Associate Professor of Industrial Engineering and Operations Research at Northeastern University. He received his Ph.D. in Operations Research from the University of Massachusetts. His research interests are in the areas of Location Analysis, Routing, Supply Chain Management and Multiple Criteria Decision Making. He has published in journals such as Management Science, European Journal of Operational Research, Naval Research Logistics, IIE Transactions, Networks, International Transactions of Operational Research, and Transportation Science. Zaharias Xanthopulos is the Manager of Forecasting at KOHLER. He received his Ph.D. in Industrial Engineering at Northeastern University. His research interests are in the areas of Multiple Criteria Decision Making, Location Analysis, Stochastic Optimization, Group Decision Making and Data Mining. He has published in journals such as Management Science and Journal of Mathematical Analysis and Applications.
Computers & Operations Research 30 (2003) 2191 – 2209
www.elsevier.com/locate/dsw
Semi-obnoxious single facility location in Euclidean space Emanuel Melachrinoudisa;∗ , Zaharias Xanthopulosb a
Department of Mechanical, Industrial and Manufacturing Engineering, Northeastern University, 334 Snell Engineering Center, Boston, MA 02115 5000, USA b KOHLER, Supply Chain Management, Kohler, WI 53044, USA Received 1 June 2001; received in revised form 1 March 2002
Abstract This paper deals with the problem of determining within a bounded region the location for a new facility that serves certain demand points. For that purpose, the facility planners have two objectives. First, they attempt to minimize the undesirable e2ects introduced by the new facility by maximizing its minimum Euclidean distance with respect to all demand points (maximin). Secondly, they want to minimize the total transportation cost from the new facility to the demand points (minisum). Typical examples for such “semi-obnoxious” facilities are power plants that, as polluting agents, are undesirable and should be located far away from demand points, while cost considerations force planners to have the facility in close proximity to the customers. We describe the set of e6cient solutions of this bi-criterion problem and propose an e6cient algorithm for its solution.
Scope and purpose It is becoming increasingly di6cult to site necessary but potentially polluting (semi-obnoxious) facilities such as power plants, chemical plants, waste dumps, airports or train stations. More systematic decision-aid tools are needed to generate several options that balance the public’s concerns with the interests of the developer or location planner. In this paper, a model is presented that generates the best possible sites (e6cient solutions) with respect to two con9icting criteria: maximize distance from population centers and minimize total transportation costs. Having all e6cient solutions at hand, the two sides can select one that best compromises their criteria. A very interesting property found is that most of these e6cient solutions are on edges of a Voronoi diagram. An algorithm is developed that constructs the complete trajectory of e6cient solutions. ? 2003 Elsevier Ltd. All rights reserved. Keywords: Location; Semi-obnoxious facility; Multiple criteria; E6cient set
∗
Corresponding author. Tel.: +1-617-373-4850; fax: +1-617-373-2921. E-mail address: [email protected] (E. Melachrinoudis).
0305-0548/03/$ - see front matter ? 2003 Elsevier Ltd. All rights reserved. PII: S 0 3 0 5 - 0 5 4 8 ( 0 2 ) 0 0 1 4 0 - 5
2192
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
1. Introduction With the advent of more stringent environmental standards, the resurgence of environmental interest groups, and a greater awareness of the public on the potential dangers of pollution, it becomes increasingly important that the placement of polluting (semi-obnoxious) facilities undergo stringent investigation and analysis. Typical examples are the selection of the location for a nuclear plant, a chemical plant, a waste dump or an airport. Ideally, these facilities should be located as far as possible from residential areas, ecosystems, water reservoirs, etc., in order to minimize the undesirable e2ects introduced into the environment by their presence. However, such a location might not be acceptable due to the high transportation costs that increase with distance from the existing facilities. This means that the farthest location found may never be built due to the unreasonably high transportation costs. Therefore, these two con9icting objectives need to be considered simultaneously. This paper addresses the problem of determining within a bounded region the location for a single semi-obnoxious facility using two objectives. The Hrst objective is to minimize the unpleasant e2ects introduced by the new facility, expressed as maximizing the minimum Euclidean distance between the new and existing facilities (maximin). The second objective is to minimize the total transportation cost between the new facility and the existing facilities, expressed as a sum of weighted Euclidean distances (minisum). For certain applications the most suitable distance metric appears to be the Euclidean one. For example, in the case of a power plant, power lines usually follow Euclidean paths and pollution is carried on a straight line in the wind’s direction away from the power plant. Assuming that there are no prevailing wind directions, the risk imposed on a demand point decreases with its Euclidean distance from the power plant. Other polluting agents, such as radiation leaking from a nuclear plant, deHnitely follows Euclidean paths that are omni-directional rays emitted from the facility. Another popular distance metric in the location literature is the rectilinear one. This rectilinear distance metric is more appropriate to use in settings where travel between facilities occurs on network of streets or aisles, such that each street, or aisle, is parallel to two rectangular directions. As it was justiHed above, the Euclidean is the distance metric of choice in this paper. For this bi-objective problem, properties of e6cient solutions are established and an algorithm is developed to construct the complete set of e6cient solutions. Location planners can then choose a solution that best compromises the two con9icting objectives. The generation of many “good” choices that include the ones, which optimize each individual objective, is deHnitely a strength of the bi-criterion model over the single criteria ones. In addition, the maximin criterion alone, which has been traditionally used for undesirable facility location, generates optimal solutions that are very often at the boundary of the feasible region (see Hansen [1]). The latter is certainly not desirable when on the other side of the boundary there are population centers belonging to a di2erent administrative region. Another location in this case which is a little closer to the demand points but which o2ers a reduction in transportation costs may be more appealing to a location planner. When multiple agents exhibiting con9icting objectives participate in the locational decision process, such as in the case of a utility company and the public, multiobjective approaches are more appropriate than single-objective ones for analyzing location problems. See also Cohon [2], Current et al. [3], Avella et al. [4] and Xanthopulos et al. [5] for support of this argument. The rest of this paper is organized as follows: The bi-criterion problem is formulated in Section 2. In Section 3, methods for solving each single objective problem are reviewed and a method for
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2193
generating e6cient solutions of the bi-criterion problem is described. Properties of the e6cient set are also developed in Section 3. Section 4 includes the algorithm for constructing the e6cient set, its computational complexity and an example problem that illustrates the steps of the algorithm and the trajectory shape of the e6cient set. The paper ends with a summary and conclusion in Section 5. 2. Problem formulation Prior to the formulation of the bi-objective model, the two individual objective problems are being reviewed. The maximin objective seeks a location X for an undesirable facility as far away as possible from the closest existing facility: Q1 : max F1 (X ) = min X − Pi ; X ∈S
16i6n
where S is a bounded convex set in R2 , Pi is a demand point in S; i = 1; : : : ; n, and X − Pi is the Euclidean distance between X and Pi . S will be approximated by a convex polygon having m edges, Bj Bj+1 ; j = 1; : : : ; m, where Bm+1 ≡ B1 . Equivalently, S = {X ∈ R2 |hj (X ) = Aj X − cj 6 0; j = 1; : : : ; m}, where hj (X ) = Aj · X − cj 6 0 is the linear constraint associated with the jth edge. Shamos and Hoey [6] solved Problem Q1 for S being the convex hull of the demand points Pi by utilizing the Voronoi diagram (see deHnition of the Voronoi diagram at the beginning of the next section). Dasarathy and White [7] extended their algorithm for a region S that is a convex polygon, i.e., Problem Q1 . Melachrinoudis and Cullinane [8] used a combinatorial algorithm to solve the weighted version of Q1 and Melachrinoudis and Smith [9] improved that algorithm by utilizing a Voronoi structure. The second objective that is considered in this paper is the Euclidean distance minisum objective, commonly known as the Weber problem. According to this problem, a point X is to be found on the plane that minimizes the sum of weighted Euclidean distances from itself to n Hxed points Pi . This problem can be used to Hnd the location of a production facility whose products are to be distributed to n customers. The objective is to minimize the total transportation cost for delivering the products to the customers (demand points): n Q2 : min F2 (X ) = vi X − Pi ; X ∈S
i=1
where vi represents the transportation cost per unit distance between the new location X and the demand point Pi ; i = 1; : : : ; n. There are many early notable contributions in both solution procedures and theory for Q2 , such as Weiszfeld [10], Miehle [11], Kuhn and Kuenne [12], Witzgall [13], Love [14] and Katz [15]. An optimal solution of Q2 can be found by the well-known iterative procedure of Weiszfeld [10]. This procedure however has a 9aw. Since F2 (X ) is nondi2erentiable at demand points, Kuhn [16] showed that the procedure fails if the iteration falls on a demand point. Several variations of the procedure have been proposed to correct this problem. For example, Ostresh [17] deHned a slightly di2erent step when iteration falls on a demand point and Wesolowsky and Love [18] redeHned the Euclidean
2194
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
distance by introducing a hyperbolic approximation. Other variations of this procedure were targeted in speeding up convergence, which is generally slow in the vicinity of demand points. Finally, many other methods di2erent from the Weiszfeld procedure have been proposed in the last 20 years. For a comprehensive review of the history, theory and solution procedures of the minisum problem, see Wesolowsky [19]. The problem under consideration is to Hnd the set of all e6cient solutions of the bi-objective Problem B: n B: max F1 (X ) = min X − Pi ; min F2 (X ) = vi X − Pi : X ∈S
16i6n
X ∈S
i=1
There are only three papers that consider multiple criteria in conjunction with a continuous solution space; two of them deal with rectilinear distances. Mehrez et al. [20] considered a weighted objective function, which combines both minimax and maximin criteria. Melachrinoudis [21] used Fourier– Motzkin elimination to solve the bi-objective minisum/maximin model, taking advantage of the small size of the resulting linear bi-objective problems. Such problems are easier to solve due to the piecewise linearity of the objective functions, while problems involving Euclidean distances such as the one undertaken in this paper cannot be decomposed into linear programs. Therefore, the solution methodology developed in this paper is completely di2erent from the ones developed in the above two papers. It is based on the Karush–Kuhn–Tucker (K–K–T) conditions for nonlinear programming and on Voronoi diagrams. The third paper by Brimberg and Juel [22] deals with Euclidean distances. The Hrst criterion (minisum) is the same as the one considered in this study. The second criterion however is di2erent, being a minisum of negative powers of distances, representing an aggregate undesirable e2ect. Not only the combination of criteria but the solution procedure is also di2erent as well from the one considered in this paper. E6cient points are found using a numerical procedure that solves a parametric system of two di2erential equations. The authors note the di6culties encountered when the e6cient trajectory is discontinuous. In this paper, we deal directly with discontinuities in the e6cient set by determining the endpoints of e6cient segments based on their mathematical properties. In the discrete solution space, there has been an increasing e2ort dealing with the multiobjective location problem of siting multiple undesirable facilities, such as in Ratick and White [23], Erkut and Neuman [24], Current and Ratick [25] and Melachrinoudis et al. [26]. 3. Model development Dasarathy and White [7] solved Problem Q1 using Voronoi polygons. Given a set of n demand points, the Voronoi polygon VPi associated with the demand point Pi is the set of all points X ∈ R2 such that X − Pi 6 min16l6n X − Pl . The Voronoi polygons VPi ; i = 1; : : : ; n, constitute a partition of R2 . The interior points of a Voronoi polygon VPi are closer to Pi than to any other demand point Pj = Pi . The Voronoi edge associated with Voronoi polygons VPi , VPj is on the perpendicular bisector of Pi Pj , i.e., all its points are equidistant from Pi and Pj . A Voronoi vertex is a common point of at least three adjacent Voronoi polygons VPi , VPj , VPk and is therefore equidistant from at least three demand points Pi , Pj and Pk . A Voronoi diagram is deHned as the set of Voronoi edges and Voronoi vertices. For a detailed coverage of Voronoi diagrams see Okabe et al. [27].
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2195
The optimal solution of Problem Q1 can be found by constructing the Voronoi diagram and enumerating a set of points as follows: Theorem 1. Consider the set U consisting of vertices of the Voronoi diagram, intersection points of the edges of the Voronoi diagram with the edges of S, and vertices of S. The optimal solution of Problem Q1 is in the set U. Proof of this theorem as well as a procedure for enumerating the elements of U to obtain the optimal solution of Problem Q1 can be found in Dasarathy and White [7]. Problem Q1 , being nonconvex, is likely to contain many local maxima. As it will be shown later in this paper, local maxima of Problem Q1 are important in identifying discontinuities in the e6cient set of Problem B. The following corollary screens the elements of the set U for local maxima. Corollary 1. (a) A Voronoi vertex belonging to Voronoi polygons VPi , VPj , and VPk is a local maximum of Problem Q1 if it is in the interior of the triangle Pi Pj Pk . (b) An intersection point of a Voronoi edge (equidistant from Pi ; Pj ) with an edge Bj Bj+1 of S is a local maximum of Problem Q1 if it is in the interior of the line segment whose endpoints are the perpendicular projections of Pi and Pj onto the line Bj Bj+1 . (c) Let Pi be the closest demand point to a vertex Bj of S and let Dj and Dj+1 be the perpendicular projections of Pi on the lines that pass through the edges of S that are adjacent to Bj . Bj is a local maximum of Problem Q1 , if Dj and Dj+1 do not both lie on the same one of the two closed half-spaces in R2 , generated by the line Pi Bj . For proof, see Melachrinoudis and Cullinane [8]. Based on Theorem 1 and Corollary 1, the set of all local maxima of Q1 can be constructed. Problem Q2 is a convex problem and its optimum occurs within the convex hull of Pi ; i = 1; : : : ; n, which is a subset of S. Therefore, Q2 can be treated as the classical unconstrained minisum problem. Its solution is found by setting the gradient of the objective function equal to the zero vector: n vi (X − Pi ) C0 : ∇F2 (X ) = = 0: X − Pi i=1 As it was pointed out in the previous section, C0 is solved by one of the variations of the Weiszfeld’s iterative procedures that resolve the issue of nonexistence of ∇F2 (X ) at demand points Pi ; i = 1; : : : ; n. Before we proceed with the generation of e6cient solutions we give the deHnition of e6ciency, nondominance and nondominated tradeo2 rate from Steuer [28] after being adapted for Problem B. Whereas S denotes the feasible region in decision space, C denotes the feasible region in criterion space, i.e., C = {Z = (z; w) | z = F1 (X ); w = F2 (X ); X ∈ S}: De nition 1. A point Xe ∈ S is e6cient if and only if there does not exist another X ∈ S such that F1 (X ) ¿ F1 (Xe ), F2 (X ) 6 F2 (Xe ) and (F1 (X ); F2 (X )) = (F1 (Xe ); F2 (Xe )).
2196
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
De nition 2. A criterion vector Ze = (ze ; we ); Ze ∈ C, is nondominated if and only if there does not exist another criterion vector, Z = (z; w); Z ∈ C, such that z ¿ ze , w 6 we and Z = Ze . Let Ef denote the set of all e6cient points and Nd the set of all nondominated criterion vectors for Problem B. Given the e6cient set Ef , the nondominated tradeo2 rate at each e6cient point Xe provides useful information to the decision maker in making his selection of the most-preferred or best-compromise point from Ef . De nition 3. The nondominated tradeo2 rate between objective z = F1 (X ) and w = F2 (X ) at an e6cient point Xe is d z : R(Xe ) = dw Xe ∈Ef The above rate represents the amount of improvement gained in the maximin objective by selecting another e6cient solution X over Xe in which the transportation cost, represented by the minisum objective, increases by a unit. Equivalently, if the location planner is willing to sacriHce $1000 (money unit) in total transportation costs, the next good choice for location (e6cient point) will be further away from the closest existing facility, as compared to the reference point Xe , by a distance equal to the nondominated tradeo2 rate at Xe . The best-compromise solution is the e6cient solution for which the decision maker’s desirable tradeo2 rate between objectives is equal to the nondominated tradeo2 rate for each pair of objectives (see Cohon [2]). By optimizing each one of the two objectives, we can obtain two e6cient solutions of Problem B. If XU is a unique optimal solution of Problem Q1 , it is an e6cient solution of Problem B. Otherwise, denote by T1 the optimal set and select as XU the one that solves minX ∈T1 F2 (X ). Similarly, if XL is a unique optimal solution of Problem Q2 , it is an e6cient solution of Problem B. If the optimal solution of Problem Q2 is not unique, denote by T2 the optimal set and select as XL the one that solves maxX ∈T2 F1 (X ). Clearly, e6cient points XU and XL deHne the range of nondominance, [F1 (XL ); F1 (XU )] for the maximin criterion and [F2 (XL ); F2 (XU )] for the minisum criterion. It is easy to see that any solution X , whose criterion vector components (one of the two or both) are outside the above ranges, is not e6cient for Problem B. Conversely, any e6cient solution of Problem B, Xe , has both criterion vector components within the above two ranges. The only point of the e6cient set that could lie on a demand point is its starting point XL . An e6cient point Xe , other than the starting point, cannot lie on a demand point. Otherwise, its e6ciency is contradicted since F1 (Xe ) = 0, F2 (Xe ) ¿ F2 (XL ) while F1 (XL ) ¿ 0. Therefore, ∇F2 (Xe ) is properly deHned over Ef , except possibly at Xe = XL , whenever XL lies on a demand point. Of course, demand points can be tested a priori for minisum optimality (see Wesolowsky, [19]). A variation of Weiszfeld’s procedure is needed only when none of the demand points is minisum optimum. Consider the constrained problem Q3 (b), F2 (XL ) 6 b 6 F2 (XU ), formed from Problem B by converting the minisum objective to a minisum constraint and introducing a new variable z to simplify
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2197
the objective function of the maximin objective: Q3 (b):
max
F1 (X ) = z
Subject to
F2 (X ) 6 b; gi (X; z) = z 2 − X − Pi 2 6 0
i = 1; : : : ; n;
hj (X ) = Aj X − cj 6 0
j = 1; : : : ; m:
In the presence of alternative optima, Problem Q3 (b) returns the one that optimizes the second objective. A point Xe is an e6cient solution of Problem B if and only if Xe is an optimal solution of the parametric Problem Q3 (b) for some b, F2 (XL ) 6 b 6 F2 (XU ). As b increases within that range, Q3 (b) generates e6cient points whose criterion vectors increase in both components. The constrained approach for generating e6cient solutions in a multiobjective problem was formally introduced by Haimes et al. [29] and was Hrst applied to water resources problems by Cohon and Marks [30] and Haimes [31]. An alternative approach for generating e6cient solutions is the well-known weighting technique that Hrst appeared in Zadeh [32]. However, as Cohon and Marks [30] point out, when the e6cient set is nonconvex the weighting technique may fail. Q3 (b) is a mathematical programming problem with similar properties to Problem Q1 , i.e., it is nonlinear and nonconvex, complicated further by the nonlinear constraint F2 (X ) 6 b. An optimal solution always exists because the last set of constraints deHnes a bounded feasible region for X and X − Pi is therefore bounded from above, for all Pi . Problem Q3 (b) may have many local maxima, out of which the global maximum has to be identiHed. Instead of directly solving Q3 (b) for several values of b to approximate the e6cient set by a discrete subset, an algorithm is developed in this paper that generates the complete e6cient set. But Hrst the loci of e6cient points are identiHed and their properties are established based on the following lemma and corollaries. Thus Q3 (b) is used in developing the properties of the efHcient segments that constitute the e6cient set. In addition to the feasibility conditions of Q3 (b), an optimal point Xe should satisfy the following optimality and complementary slackness K–K–T conditions: C1 (b):
∇F1 (Xe ) − u · ∇F2 (Xe ) −
n
!i ∇gi (Xe ; ze ) −
i=1
u · [F2 (Xe ) − b] = 0; !i gi (Xe ; ze ) = 0; "j hj (Xe ) = 0;
m
"j ∇hj (Xe ) = 0;
j=1
u ¿ 0;
!i ¿ 0; i = 1; : : : ; n; "j ¿ 0; j = 1; : : : ; m:
Lemma 1. At an optimal point Xe of Problem Q3 (b), (a) gi (Xe ; ze )=0 for at least one i; i=1; : : : ; n, and (b) Either: (i) F2 (Xe ) = b, and hj (Xe ) = 0, j = 1; : : : ; m, or
2198
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
(ii) F2 (Xe ) = b, and hj (Xe ) = 0 for at least one j, or (iii) F2 (Xe ) ¡ b, which further implies that Xe is a local maximum of Problem Q1 . See proof of Lemma 1 in Appendix A. 3.1. Discontinuities in the nondominated set In a linear multiobjective problem, the e6cient set Ef is always connected. In a nonlinear multiobjective problem, however, Ef may not be connected (Steuer [28]). The same holds for the nondominated set Nd . We expect therefore discontinuities in the nondominated and e6cient sets of Problem B. We start with discontinuities in Nd . To elaborate, observe the last case (iii) of Lemma 1, i.e., if an e6cient point Xe solves Q3 (b) and the minisum constraint is not binding at Xe , then Xe is a local maximin point. In this case, a very small increase in b will still make Xe optimal for Q3 (b). This implies that the nondominated set will be discontinuous at the criterion vector of Xe . The following two corollaries deal with discontinuity in the nondominated set. Corollary 2. Let Xe be an optimal solution of Problem Q3 (b) for some b, F2 (XL ) 6 b 6 F2 (XU ), and Ze be its criterion vector. (a) If F2 (Xe ) ¡ b, then there exist no e>cient solutions X such that F2 (Xe ) ¡ F2 (X ) 6 b, i.e., the nondominated set Nd is discontinuous at Ze . (b) If there are no nondominated criterion vectors Z(z1 ; z2 ) such that F2 (Xe ) = b1 ¡ z2 6 b, then Xe is a local maximin point. See proof of Corollary 2 in Appendix A. Let bl and bu be the lowest and highest value of b, respectively, for which the optimal solution of Q3 (b) is a local maximin point, Xe1 . Then according to Corollary 2, the nondominated set Nd is discontinuous at Ze1 = (F1 (Xe1 ); bl ) and there are no nondominated criterion vectors Z = (F1 (Xe1 ); b) for bl ¡ b 6 bu . The next corollary complements Corollary 2 by establishing continuity in the nondominated set in terms of the binding status of the minisum constraint of Q3 (b). Corollary 3. If for every b in the range [b1 ; b2 ] the optimal solution of Problem Q3 (b), Xe , is such that F2 (Xe ) = b, the subset of nondominated criterion vectors Z(z1 ; z2 ) for z2 ∈ [b1 ; b2 ] comprises a connected subset of Nd (continuous nondominated frontier). See proof of Corollary 3 in Appendix A. A geometrical interpretation of the role local maximin points play in generating discontinuities follows: The minisum constraint in Q3 (b) restricts X to be within the minisum contour F2 (X ) = b. Assume that for some value of b; Q3 (b) returns a local maximin point Xe . When b increases further, the minisum contour expands providing more choices for X . The new choices of X however are not better since Xe is the best point in its local neighborhood. 3.2. Characterization and properties of e>cient segments We investigate now the consequences of case (i) and (ii) of part b of Lemma 1. Let us start with case (i) by assuming that F2 (Xe ) = b and gk (Xe ; ze ) = 0, the latter implying that Xe is closest to Pk
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2199
than to any other demand point. In terms of the K–K–T conditions, C1 (b), this means that u ¿ 0, !i = 0 i = k, "j = 0 for all j and the optimality condition for its X component can be simpliHed as, −u∇F2 (Xe ) + 2!k (Xe − Pk ) = 0. Setting !k = 2!k =u, the complete set of K–K–T conditions becomes C2 (b):
!k (Xe − Pk ) = ∇F2 (Xe ); F2 (Xe ) = b; gk (Xe ; ze ) = 0; !k ¿ 0:
With the exception of XL , in the case it falls on a demand point Pk , all points of the Hrst e6cient segment satisfy conditions C2 (b). This Hrst e6cient segment starts at XL and terminates at some point Y which is either on a Voronoi edge or on the boundary of S. A method for determining Y is developed Hrst and a trajectory analysis, as an alternative to conditions C2 (b), is presented next to identify segment XL Y . Let VPk be the Voronoi polygon containing the minisum point XL and let Vk = VPk ∩ S be its feasible subset. Vk is bounded by a set of edges Dk , which are either Voronoi edges or line segments of the boundary of S. Let XL Y = Vk ∩ Ef be the Hrst e6cient segment of Ef . This segment will contain e6cient points that are closer to Pk than to any other demand point. Since Y lies either on a Voronoi edge or on the boundary of S, it can be found as follows: For each edge Ek E‘ ∈ Dk , the following system of equations, based on conditions C2 (b), are solved: C3 :
!k (X − Pk ) = ∇F2 (X ); X ∈ E k E‘ ; !k ¿ 0:
Among all the solutions obtained, the one that minimizes F2 (X ) is the endpoint Y of the Hrst e6cient segment. Given a value for b, conditions C2 (b) provide a system of three nonlinear equations in three unknowns to solve for Xe . However, an easier way to Hnd Xe is described next, based on the convex combination of the two criteria: minX F3 (X ; !) = −!F1 (X ) + (1 − !)F2 (X ). For X ∈ Vk , this problem can be simpliHed substantially by setting F1 (X ) = X − Pk . It subsequently becomes similar to the minisum Problem Q2 , except that its weights have been parameterized using !: Q4 (!): min F3 (X ; !) = (1 − !)vi X − Pi + [(1 − !)vk − !]X − Pk : X
i=k
Q4 (!) has to be solved for ! ∈ (0; !∗ ] to generate points of the Hrst e6cient segment, where is found as the solution of anyone of the two linear equations from (C0 ), after setting X = Y (already computed) and using the parameterized weights from Q4 (!). For practical reasons, only a few dispersed e6cient points on segment XL Y can be generated and presented to the decision maker. Consider now case (i) of Lemma 1 in which at least two of the constraints gi (X; z) 6 0 are binding. It is obvious that if exactly two of these constraints are binding Xe is on a Voronoi edge. If three or more are binding, Xe is a Voronoi vertex. Case (ii) of Lemma 1 indicates that e6cient points may also occur on the edges and vertices of S, whenever one or two of the constraints hj (X ) 6 0 are binding at Xe . !∗
2200
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
Concluding the discussion of Lemma 1 and the subsequent analysis of K–K–T, we can characterize the components of the e6cient set as either (a) points satisfying conditions C2 (b), or equivalently solutions of Q4 (!), or (b) line segments on Voronoi edges, or (c) line segments of the boundary of S. The above two corollaries state that discontinuities in the nondominated set occur only at criterion vectors corresponding to some local maximin points. Discontinuities in the e6cient set however may occur at other e6cient points as well, that are not local maximin points. Corollary 4 establishes properties of e6cient points at which discontinuities occur. Corollary 4. Let Xe1 be the endpoint of an e>cient segment of Ef , and X 2 be the starting point of the next e>cient segment as the two objectives increase. Then, either (i) Xe1 is a local maximin point or (ii) Xe1 and X 2 are alternative optima of Problem Q3 (b) for some b. See proof of Corollary 4 in Appendix A. Based on Corollary 4 we develop procedures next for identifying endpoints of e6cient segments. For example, given the starting point of the current e6cient segment, and the line on which it lies, such as a Voronoi edge, what we need to determine is the endpoint of that e6cient segment and the beginning of the next one as b increases. Then we can repeat this step until the global maximin point is encountered. These procedures will be useful in constructing e6cient segments. Without loss of generality, let us assume that the current e6cient segment is on a Voronoi edge Ek E‘ starting at Y, and the direction of increase for both objectives is from Y toward Ek . According to Corollary 4 case (ii), if there is a discontinuity in the e6cient set at some point Xk , Xk ∈ YEk ; Xk = Ek; and the starting point of the next e6cient segment is Xq , Xq ∈ Eq Es , then Xk and Xq are alternative optima of Q3 (b) for some b; F2 (Y ) 6 b ¡ F2 (Ek ). In order to determine points Xk and Xq , therefore, the following system of equations should be solved for every Voronoi edge or edge of S, Eq Es , where Eq Es ∈ Vq for q = 1; : : : ; n: C4 :
Xk − Pk = Xq − Pq ; F2 (Xk ) = F2 (Xq ); Xk ∈ YEk ; X q ∈ E q Es :
Consider two points Xk ∈ YEk and Xq ∈ Eq Es such that F2 (Xk )=F2 (Xk )−Ub and F2 (Xq )=F2 (Xq )−Ub for a small Ub. Since Xk is e6cient, F1 (Xk ) ¿ F1 (Xq ). F1 (Xk ) = F1 (Xq ) implies that, F1 (Xk ) − F1 (Xk )=Ub ¡ F1 (Xq ) − F1 (Xq )=Ub, i.e., the tradeo2 rate between objective F1 (X ) and F2 (X ) on Eq Es at Xq is greater than the tradeo2 rate on Ek E‘ at Xk . In addition, since the tradeo2 rate is continuous on line segments (see Appendix A) the optimal solution of Q3 (b + Ub) will be on Eq Es . If there is no solution to C4 for any edge Eq Es , then the complete segment YEk is e6cient. If in addition, Ek is a local maximin point, then Ek is the endpoint of the e6cient segment YEk according to Corollary 4 case (i). Theorem 1 in conjunction with Corollary 1 can be used to determine if Ek is a local maximin point. In order to determine the starting point of the next e6cient segment, the following system of equations should be solved for every Voronoi edge or edge of S; Eq Es , where
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2201
Eq Es ∈ Vq for q = 1; : : : ; n: C5 : X − Pq = F1 (Ek ); X ∈ Eq E s ; F2 (X ) ¿ F2 (Ek ): If G is the set of solutions of C5 , an optimal solution Y from G is selected, such that F2 (Y ) = minX ∈G F2 (X ). Although Y is dominated by Ek , a point at an inHnitesimal distance from Y lying on Eq Es in a direction of increase for both objectives will be e6cient, and the new e6cient segment lies on Eq Es . 4. Construction of the e#cient set Rather than solving Q3 (b) for several values of b to generate a discrete subset of the e6cient set, the following algorithm generates a complete trajectory of the e6cient set. Based on the characterization of e6cient segments and the properties of their endpoints, established in the previous section, the algorithm identiHes each e6cient segment by its endpoints and the equations of the line on which it lies, as parameter b scans the range [F2 (XL ); F2 (XU )]. Step 1 of the algorithm determines the optimal solutions XU and XL for problems Q1 , and Q2 , respectively. XL is the starting point of the trajectory of the e6cient set. Step 2 constructs the Voronoi diagram and partitions S into n polygons. Step 3 identiHes the polygon Vk containing XL and determines the smallest value of b for which the curve F2 (X ) = b is tangent to an edge of the polygon Vk . The tangent point, Y, and point XL are the endpoints of the Hrst e6cient segment of Ef . The points of that segment are solutions of the system of Eqs. C2 (b) or equivalently solutions of Q4 (!). Step 4 Hnds e6cient segments on Voronoi edges and edges of S by identifying endpoints at which discontinuities in the e6cient set occur that are caused by alternative optima in Problem Q3 (b). Step 5 Hnds e6cient segments by identifying endpoints, which are the local maximin points. The algorithm terminates when the endpoint of the e6cient set is encountered, XU . 4.1. Algorithm Step 1. Determine the maximin (XU ) and minisum (XL ) points by solving Problem Q1 and Q2 , respectively. Compute the nondominance range for each objective, i.e., [F1 (XL ); F1 (XU )] and [F2 (XL ); F2 (XU )]. Step 2. Construct the Voronoi diagram and Hnd its intersections with the edges of S. The Voronoi diagram partitions S into polygons. Each polygon Vi contains a point Pi ; i = 1; : : : ; n and is bounded by a set of edges Di ; i = 1; : : : ; n which are either Voronoi edges or edge segments of the boundary of S. Step 3. Identify the demand point Pk corresponding to the polygon Vk , containing XL . Find the smallest value of b for which the solution of Q3 (b) is a point of the boundary of Vk by solving the system of equations C3 for all edges Ek E‘ ∈ Dk . Denote such a point by Y and the corresponding edge of Vk by Ek E‘ . Points of the e6cient segment XL Y are found by solving Problem Q4 (!) for various values of !. Without loss of generality let Ek be the endpoint of Ek E‘ such that movement from Y toward Ek increases the distance from Pk .
2202
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
Step 4.1. For every edge and boundary side Eq Es ∈ Di ; i = 1; : : : ; n, Hnd Xk ∈ YEk and Xq ∈ Eq Es such that Xk − Pk = Xq − Pq
and
F2 (Xk ) = F2 (Xq ):
If there is no solution, the segment YEk is e6cient; go to Step 4.2. Otherwise, if there is a solution (Xk ; Xq ) a new e6cient segment on edge Eq Es starts at Xq . Set Y = Xq , Ek = Eq and E‘ = Es and repeat Step 4.1. Step 4.2. If F1 (Ek ) = F1 (XU ) STOP. Otherwise, go to Step 5. Step 5. If Ek is a local maximin point: Denote by G the set of points X such that X − Pi = F1 (Ek ); X ∈ Eq Es ; F2 (X ) ¿ F2 (Ek );
for Eq Es ∈ Di ; i = 1; : : : ; n:
Let Y ∈ G be the point and Eq Es the corresponding edge such that F2 (Y ) = min F2 (X ): X ∈G
A new e6cient segment starts at Y on edge Eq Es . Set Ek = Eq ; E‘ = Es , and go to Step 4.1. Otherwise, an e6cient segment will start at Ek along the edge emanating from Ek on which X − Pk increases. Let that edge be Ek Er . Set Y = Ek , Ek = Er , and go to Step 4.1. Although unlikely, it is possible that ties occur in certain steps of the above algorithm. In particular, if more than one point Xq is found in Step 4.1 then the point that yields the largest tradeo2 rate on the respective edge Eq Es should be selected. If there is a tie in the tradeo2 rate as well, then multiple e6cient segments start at points Xq and the algorithm should be executed completely for each new segment. A similar tie breaking procedure should be used if more than one point Y is found in Step 5. 4.2. Nondominated tradeo@ rate The nondominated tradeo2 rate between the maximin and minisum objective at an e6cient point, Xe , can be determined using its deHnition given earlier in this paper. The nondominated tradeo2 rate can be reduced to a formula (see Appendix A) for e6cient line segments, lying on edges Eq Es ∈ Vk , i.e., for points on Voronoi edges or edges of S, cos ’ X e ∈ E q Es ; R(Xe ) = ∇F2 (Xe ) · u where u is a unit vector on Eq Es , and ’ is the angle between u and the direction Pk Xe . For points belonging to nonlinear e6cient segments, the nondominated tradeo2 rate can be estimated by generating pairs of e6cient points which are very close to each other and directly computing the ratio of change in the two objective functions, Uz1 =Uz2 . As the e6cient set moves from one segment to another the nondominated tradeo2 rate changes value, being discontinuous at all endpoints of e6cient segments, except for the beginning and end of the e6cient set. At these points, where discontinuity occurs, a left and a right value have to be determined.
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2203
4.3. Complexity analysis of the algorithm Certain steps of the algorithm are bounded by polynomial times while others, such as solving the minisum problem or a system of equations, require numerical methods that converge within certain predetermined accuracy. Steps 1 and 2 are executed once and amount to solving for XL and XU . XL can be found by a variant of Weiszfeld’s procedure. The construction of the Voronoi diagram, its intersection with the edges of S and subsequently determining XU requires O(m log2 n + n log n + n log m) time (see Dasarathy and White, [7]). Since there exist O(n) Voronoi edges and m edges in S; O(m + n) systems C3 of three equations in three unknowns have to be solved in Step 3 to Hnd point Y . Step 4.1 requires the solution of O(n + m) systems of two equations in two unknowns. Finally, the number of executions of iterative Steps 4 and 5 combined is equal to the number of e6cient segments, k, lying on Voronoi edges and edges of S. Although some edges may be repeated, many others may not be included in the trajectory of e6cient set. Our conjecture is that k =O(n+m). The following example has been developed to apply the steps of the algorithm and illustrate the trajectory of the e6cient set, the nondominated set, and the nondominated tradeo2 rate. 4.4. An example problem A 13 by 13-square, S, is shown in Fig. 1. Inside S there are four demand points Pi ; i = 1; 2; 3; 4 for which information is given in Table 1. The Voronoi diagram partitions S into four polygons Vi ; i = 1; 2; 3; 4. Voronoi edges are shown as dashed lines. The set R of all candidates for local maximin points according to Theorem 1 is constructed as R = {Y2 ; Y4 ; Y7 ; Y8 ; Y9 ; Y10 ; Y11 ; Y12 ; Y13 ; Y14 }. Among these candidates, Y2 and Y4 are Voronoi vertices, Y7 ; Y8 ; Y10 and Y11 are intersections of Voronoi edges with edges of S and Y9 ; Y12 ; Y13 and Y14 are the vertices of S. The elements of the set R after being screened for local maxima according to Corollary 1 yield the set of local maximin points {Y2 ; Y9 ; Y10 ; Y11 ; Y12 ; Y14 }. Point Y9 is the global maximin point, XU =Y9 . The solution to the minisum Problem Q2 yields P2 as the minisum point, XL = P2 . The nondominance range is [F1 (XL ); F1 (XU )] = [0; 9:85] and [F2 (XL ); F2 (XU )] = [34:39; 125:53] for the maximin and minisum objectives, respectively. As the minisum objective function increases within its nondominance range, the algorithm provides all endpoints of e6cient segments and their criterion vectors as displayed in Table 2. The trajectory of the e6cient set, Ef , is shown in Fig. 1 with solid line segments. Starting at P2 and ending at Y9 , Ef consists of seven e6cient segments, i.e., Ef = {P2 Y1 ; Y1 Y2 ; Y3 Y4 ; Y4 Y5 ; Y6 Y7 ; Y7 Y8 ; Y8 Y9 }. The Hrst e6cient segment, P2 Y1 , is the set of e6cient points which are the solutions to Q4 (!), for 0 6 ! 6 0:778. Several values of ! in that range were considered and the solution points to Q4 (!) were connected in Fig. 1. The next four segments are parts of Voronoi edges, and the last two segments are parts of edges of S. Similarly, the nondominated set, Nd , has been constructed in Fig. 2, Nd = {ZL Z1 ; Z1 Z2 ; Z3 Z4 ; Z4 Z5 ; Z6 Z7 ; Z7 Z8 ; Z8 ZU }. The e6cient set exhibits both types of discontinuities as deHned in Corollary 4. The Hrst discontinuity occurs at Y2 , a local maximin point. The starting point of the next e6cient segment, Y3 , is determined in step 5 of the algorithm. It should be noted that Y3 is not an e6cient point, i.e., the e6cient segment Y3 Y4 is open at Y3 . Similarly, its criterion vector Z3 is dominated. To exclude these points from their respective e6cient segments, Y3 and Z3 are marked by solid circles in Figs. 1 and 2, respectively. The second discontinuity of Ef occurs at Y5 . Y5 and Y6 , determined in
2204
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209 Y7
Y 9 Y8
Y12
Y6
V3
V2 P3
•
Y1
V1 P2
Y2
•
P1
Y3
•
Y10
•
P4
Y4 Y5
S
V4
Y14
Y11
Y13
Fig. 1. Feasible region and e6cient set. Table 1 Data for the example problem i
1
2
3
4
Pi vi
(4,4) 2
(8,7) 4
(11,10) 3
(13,4) 2
Step 4.1 of the algorithm, are alternative optima of Problem Q3 (b) for b = 79:46. As it is illustrated in Fig. 2, there is a single discontinuity in the nondominated set, occurring at Z2 , which is the criterion vector of the already mentioned local maximin point Y2 . The nondominated vectors of the alternative optima Y5 and Y6 are identical, Z5 ≡ Z6 . The nondominated tradeo2 rate R between the maximin and minisum objectives is shown in Fig. 3 as a function of z1 . The segment RL R1 was constructed using several values of the nondominated tradeo2 rate calculated as the ratio of change in the two objectives for selected pairs of e6cient points. All other segments were constructed by determining the nondominated tradeo2 rate for several points along the e6cient edges using the formula provided in Appendix A. The nondominated tradeo2
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2205
Table 2 Endpoints of e6cient segments and their criterion vectors Endpoints of e6cient segments
Criterion vector
Nondominated tradeo2 rate
P2 (8; 7) = XL Y1 (10:21; 7:79) Y2 (11:25; 6:75) Y3 (9:75; 4:25) Y4 (8:5; 2:17) Y5 (8:5; 0:38) Y6 (5:70; 12:30) Y7 (5; 13) Y8 (0:38; 13) Y9 (0; 13) = XU
ZL (0:0; 34:39) Z1 (2:35; 40:40) Z2 (3:26; 44:85) Z3 (3:26; 48:72) Z4 (4:86; 63:54) Z5 (5:77; 79:46) Z6 (5:77; 79:46) Z7 (6:71; 89:15) Z8 (9:70; 122:35) Z9 (9:85; 125:53)
0.567 0.2937,0.2875 0.1815 0.1011 0.1079,0.1192 0.106094 0.0968 0.0960,0.0801 0.0935,0.0444 0.0480
130
ZU
120
Z8
110 100 90
Z7
z2
80
Z5 = Z6
70 60
Z3
50 40 30
Z2
Z1
ZL
Z4
20 10 0 0
1
2
3
4
5
6
7
z1
8
9
10
11
12
13
Fig. 2. Nondominated set.
rate, as shown in Fig. 3, is discontinuous at endpoints of e6cient segments. For these endpoints the nondominated tradeo2 rates are displayed in the third column of Table 2. Two numbers are listed, a left and a right value separated by a comma, at the endpoints of the e6cient segments where discontinuities occur. The solution of the system of nonlinear equations required in Steps 3, 4.1, and 5 of the algorithm were obtained using the Mathcad software package (Mathcad 3.1 [33]).
2206
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
0.6 RL
0.5
R
0.4
R1
0.3
0.2
0.1 RU 0 0
1
2
3
4
5
z1
6
7
8
9
10
Fig. 3. Nondominated tradeo2 rate.
5. Conclusion In many instances, the simultaneous consideration of more than one objective improves the realism of the mathematical model and provides the decision maker with more solution options in the form of the e6cient set. This is especially true for the location of a semi-obnoxious facility such as a nuclear plant, a chemical plant, or a waste dump. For these location decisions transportation cost is very important to be ignored and a single maximin objective may be inadequate. In addition, since very often the maximin optimal solution is found at the boundary of the feasible region, the generation of additional solution options in the interior of the feasible region is very desirable. The Euclidean distance maximin–minisum location model is a nonlinear and nonconvex bi-objective mathematical programming problem. The constrained approach is used to convert it to a single objective parametric nonlinear problem. Properties of e6cient solutions are established using the K–K–T conditions. It was shown that the e6cient set consists of segments not necessarily connected which are either (a) loci of points that maximize distance from a demand point to minisum contours, or (b) parts of Voronoi edges, or (c) parts of the edges of the feasible region. A complete trajectory of the e6cient set is constructed by determining the endpoints of e6cient segments. Despite the above merits, the presented model for semi-obnoxious facility location points to a number of directions for future work. The computational complexity of the developed algorithm has been delineated; however, extensive computational experiments can be performed to determine the algorithm’s average complexity as the number of existing facilities increase. In addition, future
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2207
research should point to the direction of multifacility semi-obnoxious facility location under multiple objectives. Such a problem is expected to be much more complex than its single facility counterpart. Since the (multiple) semi-obnoxious facilities should not be crowded in a single location, additional (separation) constraints have to be imposed. This in turn, will spoil the elegant Voronoi structure of the problem and will not enable the given algorithm to be generalized for multiple facilities. Acknowledgements The authors wish to thank Professor H.A. Eiselt for his valuable comments on a previous draft of this paper as well as the two anonymous referees who provided valuable suggestions for improving the paper. Appendix A. Proof of Lemma 1. (a) The optimality K–K–T condition, solved for its z-component, yields 1 − 2 · n z · !i = 0. Since Problem Q3 (b) has a bounded objective value, at least one of the !i ’s has to be i=1
positive, i.e., gi (Xe ; ze ) = 0 for at least one i. (b) Constraint F2 (X ) 6 b is either binding or nonbinding at Xe . If it is binding, we have case (i) and (ii). If it is nonbinding, then u = 0 in C1 (b) and the above K–K–T conditions reduce to the K–K–T conditions satisHed by a local maximum of Problem Q1 . Therefore, F2 (Xe ) ¡ b implies that Xe is a local maximum of problem Q1 (local maximin point). Proof of Corollary 2. (a) By Lemma 1, F2 (Xe ) ¡ b implies that Xe is a local maximin point. Assume that there exists an e6cient solution X of Problem B such that F2 (Xe ) ¡ F2 (X ) 6 b, i.e., X is a feasible solution of Problem Q3 (b). Since both X and Xe are e6cient, F1 (X ) ¿ F1 (Xe ), thus contradicting the optimality of Xe to Problem Q3 (b). (b) If Xe is not a local maximin point, then there exists a direction of improvement Ux, where Ux is very small and F1 (Xe +Ux)¿F1 (Xe ). In addition, since Xe is optimal for Q3 (b1 ), b1 ¡F2 (Xe + Ux) = b1 + Ub1 ¡ b for a small positive Ub1 and the solution of Q3 (b1 + Ub1 ) will yield e6cient point Xe + Ux. The latter contradicts the nonexistence of nondominated criterion vectors Z(z1 ; z2 ) such that F2 (Xe ) = b1 ¡ z2 = b1 + Ub1 6 b. Proof of Corollary 3. Let Xe and Xe be the optimal solutions of problems Q3 (b) and Q3 (b + -), respectively, for a very small positive - and b1 6 b ¡ b + - 6 b2 . Since both Xe and Xe are e6cient, b + - = F2 (Xe ) ¿ F2 (Xe ) = b implies that F1 (Xe ) ¿ F1 (Xe ). Since the maximum distance of any point Pi to the contour F2 (X ) = b is a continuous function of b, a very small increase - in the RHS of constraint F2 (X ) 6 b will result in a very small increase in the optimal objective function value of Problem Q3 (b). Therefore, for z2 ∈ [b1 ; b2 ] the nondominated set is connected. Proof of Corollary 4. (i) F2 (X ) is either continuous or not over the e6cient set at Xe1 . If it is not continuous, then according to Corollary 2, Xe1 is a local maximin point.
2208
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
(ii) If F2 (X ) is continuous over the e6cient set at Xe1 , then F2 (X 2 ) = F2 (Xe1 ) = b. Since both points are e6cient, this further implies that, F1 (X 2 ) = F1 (Xe1 ), i.e., Xe1 and X 2 are alternative optima of Problem Q3 (b). A.1. Nondominated tradeo@ rate on a line segment Eq Es Let Xe be an e6cient point such that Xe ∈ Eq Es and let Ze (ze ; we ) be its criterion vector, such that ze = Xe − Pk . Let u(cos .; sin .) be the unit vector on Eq Es having a direction from Eq to Es , where . is the angle between vector u and the x-axis. Consider another e6cient point on the line segment at an inHnitesimal distance from Xe , i.e., Xe + u dt having criterion vector Z(z; w). Then after ignoring dt 2 terms we have d(z 2 ) = 2z d z = Xe + u dt − Pk 2 − Xe − Pk 2 = 2(Xe − Pk ) · u dt; or d z (Xe − Pk ) · u = cos ’; where ’ is the angle between u and Pk Xe : = dt Xe Xe − Pk The directional derivative of w = F2 (X ), on the line segment Eq Es at Xe is dw = ∇F2 (Xe ) · u: dt Xe
By dividing the two derivatives we obtain d z cos ’ : R(Xe ) = = dw Xe ∇F2 (Xe )u References [1] Hansen P, Peeters D, Thisse JF. On the location of an obnoxious facility. Sistemi Urbani 1981;3:299–317. [2] Cohon JL. Multiobjective programming and planning. New York: Academic Press, 1978. [3] Current J, ReVelle C, Cohon J. The median shortest path problem: a multiobjective approach to analyze cost vs. accessibility in the design of a transportation network. Transportation Science 1987;21:188–97. [4] Avella P, Benati S, Martinez LC, Dalby K, Girolamo DD, Dimitrijevic B, Ghiani G, Giannikos I, Guttmann N, Hultberg TH, Fliege J, Marin A, Marquez MM, Ndiaye MM, Nickel S, Peeters P, Brito DP, Policastro S, Saldanha de Gama FA, Zidda P. Some personal views on the current state and the future of locational analysis. European Journal of Operational Research, 1998;104:269 –87. [5] Xanthopulos Z, Melachrinoudis E, Solomon M. Interactive multiobjective group decision making with interval parameters. Management Science 2000;46:1585–601. [6] Shamos MI, Hoey D. Closest-point problems. In: Proceedings of the 16th Annual IEEE Symposium on Foundations of Computer Science, 1975. p. 151– 62. [7] Dasarathy B, White LJ. A maximin location problem. Operations Research 1980;28/6:1385–401. [8] Melachrinoudis E, Cullinane TP. Locating an undesirable facility within a geographical region using the maximin criterion. Journal of Regional Science 1985;25/1:115–27. [9] Melachrinoudis E, MacGregor Smith J. An O (mn2 ) algorithm for the maximin problem in E 2 . Operations Research Letters 1995;18:25–30. [10] Weiszfeld E. Sur le Point Lequel la Somme des Distances de n Points Donnes Est Minimum. Tohoku Mathematics Journal 1937;43:355–86. [11] Miehle W. Link-length minimization in networks. Operations Research 1958;6:232–43.
E. Melachrinoudis, Z. Xanthopulos / Computers & Operations Research 30 (2003) 2191 – 2209
2209
[12] Kuhn HW, Kuenne E. An e6cient algorithm for the numerical solution of the generalized weber problem in spatial economics. Journal of Regional Science 1962;4(2):21–33. [13] Witzgall C. Optimal location of a central facility: mathematical models and concepts. National Bureau of Standards Report 8388. Gaithersburg, MD: NBS, 1965. [14] Love RF. A note on the convexity of siting deposits. The International Journal of Production Research 1967;6: 153–4. [15] Katz NI. Local convergence in fermat’s problem. Mathematical Programming 1974;6:89–104. [16] Kuhn HW. A note on fermat’s problem. Mathematical Programming 1973;4:98–107. [17] Ostresh Jr. L. On the convergence of a class of iterative methods for solving the weber location problem. Operations Research 1978;26:597–609. [18] Wesolowsky GO, Love RF. A nonlinear approximation method for solving a general rectangular distance weber problem. Management Science 1972;18:656–63. [19] Wesolowsky GO. The weber problem: history and perspectives. Location Science 1993;1/1:5–23. [20] Mehrez A, Sinuany-Stern Z, Stulman A. A single facility location problem with a weighted maximin–minimax rectilinear distance. Computers and Operations Research 1985;12(1):51–60. [21] Melachrinoudis E. Bicriteria location of a semi-obnoxious facility. Computers & Industrial Engineering 1999;37: 581–593. [22] Brimberg J, Juel H. A bicriteria model for locating a semi-desirable facility in the plane. European Journal of Operational Research, 1998;106:144 –51. [23] Ratick SJ, White AC. A risk-sharing model for locating noxious facilities. Environment and Planning B: Planning and Design 1988;15:165–79. [24] Erkut E, Neuman S. A multiobjective model for locating undesirable facilities. Annals of Operational Research 1993;40:209–27. [25] Current J, Ratick S. A model to assess risk, equity and e6ciency in facility location and transportation of hazardous materials. Location Science 1995;3(3):187–201. [26] Melachrinoudis E, Min H, Wu X. A multiobjective model for the dynamic location of landHlls. Location Science 1995;3(3):143–66. [27] Okabe A, Boots B, Sugihara K. Spatial tesselations, concepts and applications of voronoi diagrams. Chichester, England: Wiley, 1992. [28] Steuer ER. Multiple criteria optimization: theory, computation, and application. Malabar, FL: Krieger, 1989. [29] Haimes YY, Wismer DA, Lasdon LS. On Bicriterion formulation of the integrated system identiHcation and system optimization. IEEE Transactions on Systems, Man and Cybernetics 1971;1:296–7. [30] Cohon JL, Marks DH. Multiobjective screening models and water resources investment. Water Resources Research 1973;9(4):826–36. [31] Haimes YY. Multiple objectives in water and related land resources planning. Proceedings Seminar. Colorado River Basin, Modeling Studies, Utah State University, Logan, 1975. p. 55 –70. [32] Zadeh L. Optimality and non–scalar–valued performance criteria. IEEE Transactions on Automatic Control 1963;AC-8/59. [33] Mathcad 3.1. User’s guide windows version. MathSoft Inc., Cambridge, MA, 1992. Emanuel Melachrinoudis is an Associate Professor of Industrial Engineering and Operations Research at Northeastern University. He received his Ph.D. in Operations Research from the University of Massachusetts. His research interests are in the areas of Location Analysis, Routing, Supply Chain Management and Multiple Criteria Decision Making. He has published in journals such as Management Science, European Journal of Operational Research, Naval Research Logistics, IIE Transactions, Networks, International Transactions of Operational Research, and Transportation Science. Zaharias Xanthopulos is the Manager of Forecasting at KOHLER. He received his Ph.D. in Industrial Engineering at Northeastern University. His research interests are in the areas of Multiple Criteria Decision Making, Location Analysis, Stochastic Optimization, Group Decision Making and Data Mining. He has published in journals such as Management Science and Journal of Mathematical Analysis and Applications.