Analyzing tradeoffs between zonal constraints and accessibility in facility location

Computers Ops Res. Vol. 21, No. 1, pp. 79-99, 1994 Printed in Great Britain. All rights reserved

0305-0548/94 $6.00 + 0.00 Copyright 0 1993 Pergamon Press Ltd

ANALYZING TRADEOFFS BETWEEN ZONAL CONSTRAINTS AND ACCESSIBILITY IN FACILITY LOCATION Ross. A. GERRARD~ and RICHARD L. CHURCH$§ Department of Geography and National Center for Geographic Information and Analysis, University of California at Santa Barbara, Santa Barbara, CA 93106, U.S.A. (Received September 1991; accepted December 1992)

Scope and Purpose-Over the last three decades, many researchers in geography, regional science, operations research, and other disciplines have become interested in modeling to locate facilities. Examples of such facilities include warehouses, sand piles for urban snow and ice removal, public libraries, fire stations, and centers for many other kinds of activities. One of the most useful location models, the p-median problem (PMP), finds the optimal locations of a number of facilities on a network such that the total distance traveled by all users is minimized, assuming that each user is assigned to or served by his closest facility. In applying a location model such as the PMP, it is not uncommon for decision makers to focus attention on the distribution of facilities. To address this concern, location models have been extended to include constraints on distributional requirements. Distributional constraints are added in order to address equity. We propose a new PMP model which is a multiobjective model which can be used to analyze the tradeoff between accessibility and equity. An efficient solution approach is presented and computational results are given.

Abstract-One recent extension of the PMP is the zonally constrained median problem. This model recognizes that site selection often is influenced by the desire to distribute equitably the impacts or benefits of facilities by locating them among multiple regions, districts, or zones. Zonal constraints can be used in one form to ensure a minimum number of facilities in any zone and in another form to prevent too many facilities in any zone. However, a planner’s desire to meet zonal constraints can conflict with the desire to maximize system-wide public accessibility (minimize total distance traveled). Non-inferior compromise solutions which partially enforce zonal constraints could be most helpful to decision-makers, especially in a sensitive political climate. This paper presents a constrained multiobjective model (denoted the extended zonally constrained median problem, or EZCOMP) which can identify both supported non-dominated solutions and unsupported non-dominated solutions (which would be missed using the weighting approach to multiobjectives). A special Lagrangian relaxation is exploited in the proposed solution methodology. This is a first attempt at using a Lagrangian based approach to identify unsupported non-dominated solutions to a location model. Results on two data sets with different types of zones show the Lagrangian approach to be. efficient compared to linear-integer programming and a vertex substitution heuristic, even in the solution of problems of over 23,000 variables and 23,000 constraints.

1. INTRODUCTION

Many location models have been defined on networks where potential locations are restricted to a set of discrete choices on the network. Such location problems often address the issue of placing facilities where they can most efficiently serve dispersed demand. One network-based model which has commanded much attention over the last 25 years is the p-median problem (PMP). This model seeks to find facility locations which minimize the total distance incurred by all users as each assigns to his closest facility. The lower the total distance, the greater the accessibility of the locations. Hakimi [l and 21 introduced the problem on a network and provided a proof which allows the tR. A. Gerrard is a graduate student in the Department of Geography at the University of California at Santa Barbara. He received his B.A. and M.A. degrees in Geography at U.C.S.B. His interests include facility location modeling and multiobjective techniques. SAuthor for correspondence. §R. L. Church is a Professor of Geography at the University of California at Santa Barbara. He received his Ph.D. in Environmental Systems Engineering at The Johns Hopkins University. His research interests include systems modeling applied to a variety of economic, environmental, urban, and transportation related topics. His papers have appeared in Transportation Science, Geographical Analysis, Papers in Regional Science, and Computers dr Operations Research, as well as other journals. 79

Ross A. GERRARD and RICHARDL. CHURCH

80

search for an optimal solution to be restricted to discrete locational decisions. ReVelle and Swain [3] followed with a mathematical programming formulation whose solution under linear programming usually terminates with an all-integer result. The PMP has been applied to the location of various centers such as emergency service vehicle garages, parks, libraries, support facilities for urban snow removal equipment [4], post offices [S], and sea-based crude oil terminals [6]. In response to specific applications, the PMP has been the subject of a wide-ranging set of

14

27 12

ZONE 4

28

16

17

40

ZONE 2

24

3.5 50 37

21

51 Fig. 1. Four zones, Swain S-node data.

Analyzing tradeoffs between zonal constraints and accessibility in facility location

81

mathematical augmentations. Several extended median problems have been designed to account for real-world circumstances not addressed in the original model. These extensions have included (among many others): Rojeski and ReVelle [‘/I, Swain [S], Mirchandani and Odoni [9], Narula and Ogbu [lo], Weaver and Church [I 1 and 121. An extension of the median model which has recently been of interest invoives constraining facility placement in predefined districts or zones. This is especially appropriate for making applications more sensitive to political constraints, such as those which are geographical in character. For example, Fig. 1 depicts an urban area aggregated into a set of 55 nodes [13], with the whole area divided into four zones. Suppose that funding has been made available for the expansion of the city’s recreational parks and playgrounds. It makes sense that the parks should be located to maximize user accessibility over the entire area. However, suppose also that each of the four zones is a city district or ward with representation in the local government. In such an eventuality, one would foresee that spreading the funding out over the respective districts could be a political imperative. Each representative would want to ensure ‘a piece of the pie’ for his district in the form of at least one center ‘close’ to his constituents. As another example, suppose that undesirable activities were to be located (e.g. trash incinerators, drug rehabilitation centers). No representative would be willing to ‘bear all of the burden’ of such activities concentrated in his district alone. Such political constraints are spatially distributional in character. These two examples illustrate how the spatial distribution of benefits and burdens can be a major limitation on the location of certain public facilities and services. While equitable distribution among zones may be a political imperative, global efficiency (for our purposes minimum weighted distance) remains a legitimate system objective. As such, identifying the potential ~~u~eo~ between the two objectives could be a most useful tool for effective planning of certain facilities. This paper presents a generalized model (denoted EZCOMP, for extended zonally constrained median problem) which addresses the issue of finding the non-inferior tradeoff between locational solutions which maximize public accessibility (efficiency) and locational solutions which meet all requirements of spatial distribution (equity). Section 2 provides a brief review of relevant literature on zonally constrained median models. In Section 3, we introduce the general multiobjective model which describes our tradeoff problem. Section 4 presents a solution procedure based on a Lagrangian relaxation of a constrained multiobjective program. Computational results on two data sets, one a large (152-node) network, are reported in Section 5, including comparisons of Lagrangian relaxation, a substitution heuristic, and linear programming. Concluding remarks are offered in Section 6.

2. RECENT

RELEVANT

LITERATURE

A median model which incorporates regional distribution requirements is the zonally constrained median problem (ZCOMP), formulated by Church [14]. This model retains the minimum weighted distance objective but requires the dist~bution of facilities according to a number of zonal constraints which augment the original formulation of ReVelle and Swain [3]. Specifically, a minimum number (Pr’“) of servers (facilities) is mandated for each zone z. A maximum number (Py) is also required for each zone, where 0 < ,Fi” G Py. Py” and Py must be fixed in accordance with the number of defined zones and with p for the problem to be feasible. The upper limit constraint (P,“““) would be especially useful in the case of undesirable centers or ‘obnoxious’ facilities. In ZCOMP, assignments can take place across zonal boundaries, but zones do not overlap. Note also that the ‘zones’ need not be geographical in character but may be based on essentially any common characteristic of a bundle of potential facility sites. This means that the ‘zonal’ constraints may be applied to modifying the distribution of facilities in ways other than to achieve spatial equity. ReVelle and Elzinga [ 153 have structured and solved a special case of ZCOMP with PF” = 1 for each zone, with no demand assignments allowed across zonal bounda~es (and no overlap of zones), and where an upper limit fPy) can be specified for a zone. Berman et al. [f6] have structured and solved a related problem. in their work, PFi”= 1 for each region/zone, but zones are allowed to overlap, and assignments can cross zonal boundaries, although no upper limit (Py) can be specified for a zone. Berman et al. also solved a p-center problem under similar constraints. Recent developments of Church [14], ReVelle and Elzinga [15], and Berman et al. [i63 clearly

Ross A. GERRARD and RICHARD L. CHURCH

82

indicate recent interest regarding zonal constraints in location-allocation problems. Of these three recent problems, the model in this paper (EZCOMP) most closely resembles ZCOMP [14] but adds a multiobjective capability not seen in any of the published works on zonally constrained median problems. Table 1 presents a comparison of EZCOMP with the three recent papers.

3. EXTENDED

ZONALLY

CONSTRAINED

MEDIAN

PROBLEM

(EZCOMP)

Many planning decisions are made in an environment of convicting objectives. The advent of zonal constraints raises the issue of the potential tradeoff between maximizing accessibility and fairly distributing services among the various regions [17). In this context, decision-makers might ask if there are any compromises between absolute enforcement of zonal constraints (the ZCOMP solution) and the complete absence of distributional concerns (the PMP solution). Decision-makers could use such information to identify what they consider to be a ‘fair’ solution which appropriately trades-off efficiency and equity. The following model is proposed as a technique which can be used to identify such compromises, and the PMP and ZCOMP solutions, within a unified framework. The formulation is:

minimize minimize

subject to

CXij>l

for each i

(1)

Xij ~ Yj

for each i and j

(2)

C

(3)

yj=P

c yjppzmin-en

for each z

(41

C Yj~P~+~”

for each

(5)

.WL

z

ML

Xij=O or 1

for each i and j

(6)

5-O or 1 F”&O T>O

for each j for each z for each z

(7) (8) (9)

where : i refers to a demand node j refers to a potential facility site 2 refers to a particular zone or region al =demand (population) at node i d, =shortest distance on the network between demand i and facility j = the fixed number of facilities to be located P Z, = the set of potential facility nodes in zone z Pf”” = the minimum number of facilities which may be located in zone z Py = the maximum number of facilities which may be located in zone z 1 if demand i assigns to facility x,

=

0 if not

j

83

Analyzing tradeoffs between zonal constraints and accessibility in facility location

y, f S,, S Gin y

= 1 if node j receives a facility Oifnot = weighted distance objective = total zonal constraint violations objective = the number of facilities in zone z less than J’Fi” = the number of facilities in zone z greater than Pyx.

The fo~ulation of EZCOMP retains almost all of the structure of ZCOMP and the PMP. There are two objectives to EZCOMP. One objective is to minimize total weighted distance between all patrons and their closest facilities. The other is to minimize total viotations of the zonal constraints where the violations are represented by the variables Fin and y, which appear in constraint sets (4) and (5). Constraint set (4) establishes that the variable I’Fi” will represent the number of facilities less than PF” located in zone z. The index z represents specially defined sets of potential facility sites. These sets most intuitively would be geographic zones or districts but could stand for sites which have some other common element. Constraint set (5) establishes that the variable r” will represent the number of facilities greater than fy which have been located in zone z. Together cz (VT’“+ y) 1s ’ th e total number of zonal constraint violations, and this expression is minimized in the objective along with the weighted distance objective. To review the rest of the problem, constraint set (1) requires that each demand node be assigned; constraint set (2) ensures that any assignment is limited to those locations where a facility has been placed. Constraint (3) specifies that p facilities are to be located. The constraints in (6) and (7) ensure the ‘all or nothing’ assignment of demand nodes and that a center is either located at a node or not. Constraints (8) and (9) are non-negativity constraints. For our purposes we will assume that each of n nodes is both a point of demand and a potential point for a facility location. EZCOMP is quite general in that it maintains many of the properties addressed by Church [14], ReVelle and Elzinga [lSJ. and Berman et al. [16] while expanding on those models in certain respects. EZCOMP expands on ReVelle and Elzinga in that PFin can equal zero, expands on Berman et al. in that a maximum number of facilities (Pr”) can be specified for a zone, and expands on both ReVelle and Elzinga and Berman et al. in that Py” can be greater than one. EZCOMP is structurally very similar to Church [14], but differs in its allowance for zonal constraint violation in a multiobjective context. Another look at Table 1 shows that EZCOMP retains all of the noted properties of previous work except for zonal overlap. It is unique in that it is multiobjective and can identify the non-inferior tradeoff between optimal accessibility and zonal constraint enforcement. It is important to note that since the p-median problem (known to be NP-complete) is a special case of ZCOMP and EZCOMP, both ZCOMP and EZCOMP are NP-complete problems. A multiobjective problem such as EZCOMP may be approached using either of the two classical methods that have been used to generate tradeoffs between conflicting objectives: the weighting method and the constraint method [ 181. For EZCOMP, the weighting method involves combining the two objectives for a given set of weights (W, and W,) in the following way: minimize

W,S,, f WJ,,

subject to

EZCOMP constraints (l)-(9)

For a given set of weights WI >O and W, 20, the weighted model may be solved to yield a non-inferior combination of public accessibility and adherence to zonal location requirements. The relative values of the weights determine which of the two objectives is favored in any particular Table 1. Comparison of recent zonally constrained median models ReVelle and Elzinga [ 151 NO

Zones can overlap? Assignments across zonal boundaries? P,i, can be greater than one? PFin can equal zero?

NO NO NO

Py can be specified? Tradeoff analysis possible?

YES NO

Berman et ol. [Ifi]

Church [ 141

This paper

YES YES NO YES (with slight modification) NO NO

NO YES YES YES

NO YES YES YES

YES NO

YES YES

Ross A. GERRARDand RICHARDL. CHURCH

84

solution. A high relative weight on total distance will cause the solution to tend towards that of the PMP. A high relative weight on the zonal conditions will cause the solution to approach that of ZCOMP. Depending on the values assigned to the weights, solutions between the two extremes can be generated [18]. The general structure of the constraint method is to optimize only one objective while requiring a specific performance for the other objective by means of a constraint. This will be our approach. Accordingly, our approach to EZCOMP is: minimize

S,,

subject to

EZCOMPconstraints

(l)-(9)

WV

&GQ

Constraint (10) limits the number of zonal constraint transgressions to a number Q, which is specified by the user. If Q=O, then the optimal solution to the problem is clearly the zonally constrained median solution. If Q equals or exceeds the number of violations in the p-median solution, then the p-median solution is optimal. A compromise solution can be generated between these two extremes by specifying an appropriate level of Q. The constraint method is especially attractive because it may identify some noninferior solutions which are missed by the weighting method. A study of Fig. 2 will clarify this point. For purposes of illustration, the point of interest in Fig. 2 will always be designated point ‘C’. For a value of Q between 0 and the maximum possible number of zonal constraint violations, there are three possibilities for the corresponding weighted

1

low

low (b)

(a)

1 weighted distance

weighted distance

high

high 1

I

I

hi; low # zonal constraint violations (Q)

low bib # zonal constraint violations (Q)

weighted distance

:A . -. . **.

*. ..

*. : . . . . . . . . . . .. : . B

Of-

.c I

bib low # zonal constraint violations (Q) Fig. 2. (a) Point C is a supported nondominated solution. Point C can be identified either by constrained or weighting methods. (b) Point C is an unsupported nondominated solution. Point C can be identified by constrained method only. (c) Point C cannot be identified by either constrained or weighting methods (point C is an inferior solution).

85

Analyzing tradeoffs between zonal constraints and accessibility in facility location

distance value on a noninferior frontier. The first possibility is a value which lies above the line segment connecting the two neighboring noninferior points [see Fig. 2(a), point C]. Such a solution (point C) is a ‘supported nondominated’ solution and can be identified either by the weighting method or by the constraint method [19]. The second possibility is a value [Fig. 2(b), point C] beneath the aforementioned line segment but still above (less in terms of weighted distance) the more stringent solution which specifies a lower number of zonal violations (point B). Point C in Fig. 2(b) cannot be identified by a weighted sums method because it is dominated by a linear combination of points A and B. However, it is still a noninferior solution and thus should be of interest to both the modeler and the decision maker. Such a point is an unsupported nondominated solution [ 191, or a ‘gap point’. Unlike the weighting method, the constraint method does not suffer the deficiency of being ‘blind’ to unsupported nondominated solutions. The third possibility is a value [Fig. 2(c), point C] higher in weighted distance than the solution with a lower number of zonal violations (point B). Such a solution is inferior to solution B, which offers both lower weighted distance and fewer regional constraint violations. The generation of the complete noninferior set of solution to bi-objective problems has been the subject of increasing attention. Solanki [20] has introduced a general method for estimating the noninferior set of points and has applied it to a location problem. His method should be applicable to other mixed integer programming problems as well. Rather than finding the exact noninferior set, Solanki finds an approximation within a specified error tolerance. The distance between feasible points in objective space and a noninferior point is measured by a weighted Tchebycheff distance metric. Solanki’s method has the advantage of identifying unsupported nondominated solutions such as in Fig. 2(b). Solanki used linear programming with branch and bound to solve the small sized covering problems given in his example. Solanki’s approach cannot be easily extended to the much larger p-median model. Due to the need to detect unsupported nondominated alternatives to assess the efficiency/equity tradeoff between the PMP solution and the ZCOMP solution, we have chosen to explore the possibility of using a constraint method to solve the EZCOMP model. In contrast to Solanki, our method will be specific to one model (EZCOMP) and will be an exact (not an estimation) procedure. The remainder of this paper presents the details of our approach along with an application. This is a unique application since virtually all multiobjective p-median models have been solved using a weighting approach.

4. SOLVING

EZCOMP

USING

LAGRANGIAN

RELAXATION

In the previous section, we discussed the possible application of both the weighting and the constrained methods for solving the extended zonally constrained median problem (EZCOMP). This section presents details of the use of Lagrangian relaxation to solve the constrained version of EZCOMP. In the early 197Os, the method of Lagrange multipliers began to be widely recognized as a powerful solution approach to various linear-integer programming problems, often in company with a branch-and-bound procedure. Much of this surge of interest was based on the work of Rockafellar [21, 221. In the p-median literature, Lagrangian relaxation was introduced by Narula et al. [23] and used and extended by Christofides and Beasley [24], Hanjoul and Peeters [25], as well as others. Many problems have benefited from this approach: optimizing bank float [26]; the traveling salesman problem [27], [28]; and the capacitated plant location problem [29]. A good review paper has been provided by [30]. Fisher has also authored a more basic pedagogical review [31]. The latter work states the three major questions of any attempt at using Lagrangian relaxation: (1) Which constraints should be relaxed to form the Lagrangian dual? (2) How can good multipliers be computed? (3) How can a feasible primal solution be calculated from a Lagrangian solution?

dual

The ideal Lagrangian dual is a problem which retains most of the structure of the primal but whose solution is easier to obtain than the primal. Because the dual is less constrained, solutions

Ross A. GERRARD and RICHARDL. CHURCH

86

to the dual will always (in the case of any minimization problem such as a minimum weighted distance model) be as small or smaller than any solution to the primal. This property of the dual is useful for establishing lower bounds on the optimal primal solution and thus for evaluating the performance of any given primal solution. The valid lower bound on the optimal primal is the highest dual solution which has been found. If we can identify a high dual bound, sufficiently close to the lowest known primal solution, then that primal solution can be verified as optimal. Finding such a tight lower bound usually requires multiple solutions of the dual problem. In the interest of rapid computation, the dual must be easily solvable, and a ‘good’ set of multipliers must be provided if we are to get a good bound. Finally, a feasible primal must be easily calculated from the dual. Like Narula, Ogbu, and Samuelsson [23] did for the PMP, we have chosen to relax the assignment constraints (1). This yields the following dual: Lagrangian Dual of EZCOMP:

LR

G, = 1 c a,d,X, + 1 ui 1 -c

minimize

i

subject to

j

i(

Xii j

>

Constraints (2)-(10)

where : Ui is the Lagrange multiplier associated with the ith assignment constraint G, = Lagrangian dual objective function. By rearranging terms: G, = C Ui+ min 11 I

(aidij - ui)Xij.

1 i

Since the expression Ci cj (aidij-ui)Xij is to be minimized, it is clear that i should assign to j (Xi,= 1 and thus rj= 1) only in company with negative values of (aidij-ui), starting with the very smallest. Based on this reasoning, one may define Ej=Ci min(&j-ui, 0). Ej is the sum of the negative contributions to the objective as a result of allowing assignment to j (i.e. locating at j). To optimize G,,, the problem is reduced to minimizing the sum of the Ej values as long as p sites are selected and the other constraints are met. The dual formulation can then be expressed as:

Equivalent Lagrangian Dual for EZCOMP

minimize

c EjYj+c i

subject to

ui I

c Yj=p ; Yj~ PiW. c

Yj
en

for each z

(4)

foreach z

(5)

jeZ=

~(~“+~)
(10)

1

for each j

(7)

F”>O

for each z

(8)

V>O

for each z

(9).

This equivalent dual formulation is an extended version of that developed for the PMP by Narula et al. [23]. However, in the case of EZCOMP, we are faced not merely with locating p facilities

Analyzing tradeoffs between zonal constraints and accessibility in facility location

87

as in the PMP, but we are also faced with the zonal constraints and the limit that must be maintained on their violation. These constraints add a significant degree of complexity to that of solving the PMP. The equivalent dual for the classical PMP, lacking constraints (4), (S), and (lo), can be solved by picking greedily the p best sites in terms of the most negative Ej values. Unfortunately, constraints (4), (5), and (10) do add a degree of complexity that at first inspection might cause one to drop this relaxation as a fruitful approach. However, consider the possibility that we select the p sites with the most negative Ej values. Then four cases can arise. Case I

For the p sites chosen, cz (1/1”‘”+ r) < Q. Constraints (4), (5), and (10) have been met, and an optimal solution to the dual is identified. Case 2

For the p sites chosen, cz (Tin+ v)=Q +0 and no zone has been allocated more than P$,, facilities. All violations arise from the allocation of fewer than PTi” facilities to some zones. Denote such zones as zones (c). Other zones may have been allocated between PFi” and Py facilities. Denote such zones as zones (b). After selecting the p lowest Ej values, a few additional steps suffice to find the optimal relaxed dual with Q or fewer violations. First of all, consider only zones (b). Group all the sites belonging to all zones (b) and order them from lowest to highest in terms of Ej values. The 8 highest, or most marginal, Ej values are selected and their respective sites eliminated from the solution, as long as each zone (b) retains at least Z’yi” sites [and thus remains a zone (b)]. In the place of these 8 most marginal sites from zones (b) come the 8 best sites (lowest E, values) that were not selected from zones (c), as long as no such zone goes above Py”. By introducing 0 sites in zones of type (c), 8 violations of zonal constraints are eliminated. We are guaranteed the optimal solution because the 8 best sites which meet our limit on zonal constraint violations are chosen. To calculate the final objective function value of the Lagrangian dual, add the Ej values of the selected sites to ‘& ui. Case 3

For the p sites chosen, cz (yin + y) = Q + 0 and no zone has been allocated less than ,Tin facilities. All violations arise from the allocation of too many facilities to some zones. Denote zones with more than Py sites as zones (a). Denote zones (b) as the remaining zones, where the number of locations is between Plfi” and PFax. First group sites belonging to zones (a) and remove the 0 most marginal sites as long as each such zone does not fall below Py”. Then select the 0 best additional sites in zones (b) as long as no more than Py sites are placed in any zone. The result is an optimal dual solution which has no more than Q violations. Case 4

For the p sites chosen, cz (c” + y) = Q + 0 and some zones have been allocated more than Pyx facilities, some less than Pyi”, and some in between Pri” and PFx. Zones (a)-(c) are defined as before. The procedure for Case 4 is somewhat more complicated than for Cases 2 and 3. For whatever value that 8 may be, one does the following: of the Ej values of all sites selected from zones (a), take the highest or least negative). Call the site am l of the Ej values of all sites selected from zones (b), take the highest or least negative). Call the site b” l of the Ej values of all sites not selected from zones (c), take or most negative). Cali the site c’ a of the Ej values of all sites not selected from zones (b), take or most negative). Call the site 6’. l

most marginal (i.e. most marginal (i.e. the best (i.e. lowest the best (i.e. lowest

Now compare the three following possibilities for modifying the solution: (1) remove am and replace it with b’. (Reduces violations by one) (2) remove am and replace it with c’. (Reduces violations by two) (3) remove b” and replace it with c’. (Reduces violations by one.)

88

Ross A. GERRARD

and RICHARD L. CHURCH

Of these three possibilities, take the one that yields the least increase in Cj EjYj and make the switch. Return to recompute the a”, b”, cl, and b’ values. Continue the process until cz (y+ yy;nax)
(3)

(4) (5)

problem with a given set of multipliers (using one of the four cases given above) if the new dual solution is higher than the currently highest known dual solution (dual bound), install the new dual solution as the lower bound on the primal from the dual solution (which consists of p sites), find the associated primal solution by assigning each demand node to its closest facility and calculating total weighted distance. Note that the dual solution does not violate the constraint cr (en+ r)
In order to calculate new Uis, we employed the well-known method of subgradient optimization, which has been widely applied in procedures involving Lagrange multipliers. A general explanation of subgradient optimization may be found in Held et al. [32] and Fisher [30]. Narula et al. [23J described the technique for their approach to the PMP, and we use the same basic approach here. The first iteration of any problem is always solved using an arbitrary set of multipliers (say ui=minj, i aidij for every i). At each iteration k+ 1 the multipliers are updated as follows: ui(k + 1) = max (0, u,(k) + stepsize (k)*sk(i)) A(zphat - zhat(k)) stepwise (k) = 1 CS”(W where : /I zphat zhat(k) Sk(i)

=alternatively 0.5 or 1.0 (alternates every five iterations) = best known primal solution = dual value solved at iteration k = (1 - Cj X,j) = ith component of the subgrad~ent at iteration k.

5. COMPUTATIONAL

RESULTS

The previous section outlined a solution procedure for EZCOMP utilizing the relaxation of the assignment constraints and the subgradient optimization of the Lagrangian dual bound. Other logical solution procedures for the EZCOMP linear program exist. One would be the use of linear programming (LP) with branch and bound. Another candidate would be the use of the Teitz and Bart [33] (TAB) heuristic modified for use with zonal constraints. All three of these techniques have been used to solve EZCOMP, and comparative results are presented below. The use of mixed integer programming was carried out by writing a Fortran 77 program to create an MPS file which was submitted to the IBM Optimization Subroutine Library (OSL) version 1.2 for solution. A branch and bound procedure to resolve fractional solutions was

89

Analyzing tradeoffs between zonal constraints and accessibility in facility location Table 2. EZCOMP

results (4--S Sacihties) on 4-zone swain data (%-nodes), comparison linear-integer-pfo~amming

of Lagrangian

relaxation,

TAB, and

No. of facilities (zonal constraint set)

Pattern (figure label)

No. of zonal constraint violations

Total weighted distance

Lagrangian time [No. of iterations, No. of BB nodes]

TAB time (20 runs) [No. of substitutions, No. of TAB cycles] (No. of optima, No. of runs)

4(A)

1,15,20,22

0-t

3558.268

0.04 s [32,0]

0.60 s [185,48]

30.64 s [1606,0]

(1 l/20) 0.53 s [202,52] (2Oi20)

24.81 s [1378,0]

4(A)

1,3,22.41

18

3283.862

0.02 s C22.01

403)

15,30,41,53 (iv)

ot

5641.191

0.13s [81-O]

4(B) 4(B)

15,21,24,30 12,15,21,30 (N4)

5660.076 5438.989

0.10 s [64,0]

4(B)

1,15,20,22 (ni3f

3558.268

0.04 s [33.0]

4(B)

1,22,31,41 (N’)

3489.188

0.07 s 153.01

4(B)

1.3,22,41 (N’)

3283.862

0.03 s [30.0]

1,3.20,22,36

3036.809

0.09 s [73.0]

1,3,16,17,41

2989.424

0.07 s [S3.0]

1,3,10,22,36

25

2950.413

0.04 s [30.0]

6(A1

1,3,16,17,20,36 (M”)

at

2742.571

0.12 s [86,0]

6b4)

1,3,16,17,24,41 (M*)

1

2732.940

0.85 s [601.13]

6641

1,3.10,16,17,36 (W’)

4

2656.174

0.06 s [41 .O]

1,3,16,17,21,24,41

0t

2531.754

0.07 s [41,0]

OS9 s [227,53] (0120) 0.46 s [215,60] (20/20) 0.46 s [ 174.4513 (7120) 0.56 s [205,48] (14/20) 0.58 s [202,52] (20/20)

LP time [No. ofiterations, No. of BB nodes]

14.06 s [795.0] 33.31 MIN [51006,2572] 14.54 s [775.0] 21.83 s 11229.81 15.84 s [X03.0]

0.57 s [204.44] (20120) 0.47 s [185,43] (2!20) 0.58 s 1204,531 (19~20)

26.73 s [1480.0-J

0.48 s [199,45] (20/20) 0.46 s C196.451 (2i20) 0.67 s 121I.561 (20’20)

22.71 s [1265,0]

0.76 s [242,47]

10.21 s [584,0]

18.09 s Cl iOS.O] 21.95 s [1244,0]

24.01 s [1166,8] 14.35 s [916.0]

(6i20)

1,3,9,16,17,21,36

11'

2440.860

0.05 s [34.0]

1.2,3,10,16,17.36

3§

2427.411

0.07 s [52.0]

1,3,16,17,21,24,41,50

ot

2441.766

0.08 s C44.01

2299.065

0.93 s C542.171

1,3.9,15,16,17,21,24

0.75 s [263,&l] (20/~20) 0.60 s 1228,521 (10:20)

17.11 s [991,0]

0.71 s [301,50] (3:20) 0.80 s 1265,521

8.37 s [563,0]

11.55 s [747,0]

9.47 s [577.2]

(1 l/20)

1,2,3,6,16,17,21,36

2224.469

Zonal constraint set A Zone

P”“’ 2

pin E

1 2 3 4

2 2 2 2

I 1 1 1

0.12 s C80.01

0.72 s [250,48] (20/20)

15.96 s [944,0]

Zonal constraint set B Zone

P”“” I

pn’” I

1 2 3 4

2 2 0 2

1 2 0 I

tZCClMP solution. :Best TAB solution. $p-median solution. 1iNo solution with seven facilities and two violations was identified. When Q was specified as 2, the one-violation solution was identified, indicating that all two-violation solutions are dominated. The results given are for Q = 1. The patterns labeled M’ (Fig. 3) and N’ and N’ (Fig. 4) are unsupported nondominated (‘gap point’) solutions.

automatically invoked by OSL if necessary. None of the heuristic preprocessors available on OSL for integer programming was used. The TAB heuristic was programmed in Fortran 77. This heuristic has been extensively tested and found very effective for the PMP [34]. Our modified method begins with an initial solution of p facilities. The list of n nodes is then processed sequentialiy. Each node which is not currently a facility is evaluated as a replacement for each of the current sites. Desirable substitutions are those which (1) reduce total weighted distance and (2) help to reduce zonal constraint violations to a number less than or equal to Q. All potential substitutions are evaluated

90 Table 3. EZCOMP

No. of facilities (zonal constrained 12 12 12 12 12 12 12 12 12 17 17 17 17 17 17 17

Ross results

(12 and

17 facilities) linear-integer

Pattern 1,7,9,12,14,20, 24,29,33,34,53,55 1,3,9,12,14,15, 21,24,29,33,50,55 1,3,9,12,14,15, 20,24,29,33,53,55 1,3,4,6,12,14,15, 21,24,29,33,55 1,2,3,6,12,14,15, 21,24,29,33,55 1,3,4,6,12,14,15, 16,21,24,29,55 1,2,3,4,6,12,14, 21,24,29,33,55 1,2 3 33 4> 6 ,12,15, 16,21,24,29,55 1,2,3,4,6,12,14, 16,21,24,29,55 1,3,4,6,7,12,14,16,17, 21,24,25,29,30,39,46,50 1,2,3,4,6,7,12,14,16, 21,24,25,29,33,34,50,55 1,2,3,4,6,7,12,14,16, 21,24,25,29,30,39,46,50 1,2,3,4,6,7,12,14,16, 17,21,24,25,29,32,50,55 1,2,3,4,5,6,7,12,14, 16,21,24,25,29,33,50,55 1,2,3,4,5,6,7,12,14, 16,17,21,24,25,29,50,55 1,2.3,4,5,6,7,10,12, 14,16,21,24,25,29,50,55

tZCOMP solution. $Best TAB solution. $p-median solution. “Unsupported nondominated

A. GERRARDand RICHARDL. CHURCH on ltzone swain data (55-nodes), comparison of Lagrangian programming (P,mz= 2 and P:“’ = 1 for each zone)

No. of zonal constraint violations

Total weighted distance

Lagrangian time No. of iterations, No. of BB nodes]

ot

2091.729

0.16 s [52,0]

1

1856.513

0.06 s [22,0]

11

1863.452

-

2

1784.719

0.16s [58,0]

39

1745.383

0.11 s [48,0]

3f

1755.454

-

4

1688.440

0.08 s [34,0]

4:

1693.001

-

5§

1659.175

ot

TAB time (20 runs) @IO. of substitutions, No. of TAB cycles] (No. of optima, No. of runs) 1.24 s [464,69] (9/20) 1.14 s c437.593 (O/20)

0.97 s [375,63] (2/20) 1.15 s [352,57] (O/20)

relaxation,

TAB, and

LP time [No. of iterations, No. of BB nodes] 7.91 s [535,0] 6.25 s [381,0]

7.18 s [496,0] 18.23 s [600,23]

1.22 s [344,58] (O/20)

6.73 s [476,0]

0.05 [23,0]

1.03 s [350,56] (11/20)

6.05 s [432,0]

1377.906

1.20 s [346,0]

6.27 s [376,0]

1268.151

2.21 s [697,61]

1.40 s [453,55] (19/20) 1.36 s [442,60] (O/20)

1272.887

-

1231.865

0.22 s [76,0]

6.87 s [311,0]

1203.107

0.08 s [24,0]

1190.019

0.09 s [33,0]

1186.000

0.50 s [186,0]

1.39 s [447,59] (9/20) 1.28 s [430&O] (3/20) 1.26s [411,55] (19/20) 1.56 s [413,69] (20/20)

7.99 s [337.2]

5.26 s [352,0] 4.57 s [315,0] 4.82 s [303,0]

(‘gap point’) solution.

based on a combination of both criteria. Substitutions which affect the number of zonal constraint violations are appropriately rewarded or punished. Any potential exchange from a solution with more than Q zonal constraint violations to a solution with fewer violations is assigned a very largely negative (i.e. favorable) value; an exchange which would increase violations from a level less than or equal to Q to a level greater than Q is assigned a huge positive (i.e. unfavorable) value; potential substitutions which would begin from a solution of fewer than Q violations and would not result in a solution of any more than Q violations are considered exclusively on improvement in weighted distance. Of the p possible substitutions, the one which would most improve the objective is made. Another node not currently in the solution is then similarly evaluated and the best of the p possible exchanges is made. Once the list of nodes has been exhausted, a TAB cycle has been completed. Substitutions continue in this manner until a complete TAB cycle is undertaken where no possible substitutions are found which would improve the objective. At this point the process terminates. The results presented in this paper are the best result out of twenty runs of TAB, each with a different initial solution. The final of the three solution procedures, the Lagrangian relaxation method, was initiated using multipliers defined as suggested by Narula et al.1 ui=minj+i aidij. The Lagrangian procedure for EZCOMP was programmed in Fortran 77. All the solution procedures were implemented on an IBM POWERserver RS/6000 model 930 running AIX 3.1. The Lagrangian relaxation method was carried out if necessary as part of a branch and bound procedure designed to resolve any ‘duality gaps’ which might occur after the initial attempt at convergence. For the results reported here, a duality gap of no more than 0.10% (for some problems 0.0010%) of the dual’s value was allowed. If a gap of this size or smaller could not be achieved

Analyzing tradeoffs between zonal constraints and accessibility in facility location

91

after an appropriate computational investment, the branch and bound routine was invoked. The process builds a tree of branches and nodes as it takes sites in the best identified primal and successively forces them in and out of newly calculated solution to create new tree nodes. The best primal (and its facility pattern) and best dual originally calculated identify the first node of the tree. From this node (and any subsequent node) two additional nodes representing two new problems may be formed: one node which forces a site from the preceding node out of its solution and one node which forces the site in. At the node which forces a site out, a new Lagrangian problem with the necessary restriction (and with any other restrictions resulting from previous such decisions in the tree) is solved. Tree nodes which force a site from the preceding node in involve only transferring the predecessor node’s solution to the new node, as obviously it is the same solution as the predecessor node. The choice of the node on which to branch is determined by which active node has the lowest dual bound. The branch-and-bound strategy is basically an attempt to drive the dual bound of any branch of the tree to within the specified tolerance of the best known primal (which commonly is found at the first node but may be identified at a subsequent node) and thus eliminate that branch from further consideration as a source of an improved primal. Any node in the tree is either an interior node (and not eligible for further branching), an active terminal node (and thus eligible for branching), or a fathomed terminal node. A node is fathomed if p sites h&e been fixed into the solution, n-p sites have been fixed out, or the dual bound at that node has risen to within 0.10% (or whatever tolerance has been specified) of the best known primal in the tree. Once branch-and-bound has been invoked, the quickest problem resolution occurs if each of the tree nodes created quickly results in a dual bound high enough for the node to be fathomed. This would typically occur if the first node of the tree identifies the global optimal solution, but fails to achieve a dual bound high enough to verify optimality. In this case, each branching (creation of two nodes, one fixing a site out and other fixing a site in) often leads quickly to the fathoming of the ‘out’ node, and successively only the ‘in’ node is available for branching, until p sites have been fixed in, at which point the process would terminate with a total of 2p+ 1 total nodes in the tree. To achieve this desirable outcome (if branch-and-bound must be used), it is necessary to make a sufficient investment in the first tree node. Investing extra iterations at the first node aids the efficiency of the branch and bound tree in two ways: (1) it increases the chances of finding the optimal primal (thus making it easier for subsequent nodes to be fathomed due to excessively high dual bounds) and (2) the Lagrange multipliers at each tree node are saved and transferred as the starting multipliers of offspring nodes. This second property means a better start for Lagrangian problems well down in the tree from the first node. Many examples of this advantage have been confirmed in computational testing of EZCOMP. Overall, the choice of how many iterations to invest at the first tree node and at subsequent tree nodes must be determined through computational experimentation. The EZCOMP problems solved here are based on the following data sets: the %-node data of Swain [13] divided into 4 zones (Fig. 1) and 12 zones (Fig. 5), and the 152-node data of Migereko [35] divided into 17 zones (Fig. 6). Each node is both a demand point and a potential facility site, and n is defined as the number of nodes. The 55-node, 4-zone problem contains 3090 constraints and 3088 variables; the 55-node, 1Zzone problem contains 3106 constraints and 3104 variables. The 152-node, 17-zone problem contains 23,292 constraints and 23,290 variables, a relatively large problem compared to most p-median problems solved in the literature. (For the linear programs run under OSL, reduce the totals above in each case by n, as X, was used as the site variable rather than Yj.) The following tables contain problem results for solution by Lagrangian relaxation, TAB, and LP: Table 2 for 55 nodes and 4 zones, Table 3 for 55 nodes and 12 zones, and Tables 4 and 5 for 152 nodes and 17 zones. The tables offer the problem parameters, solutions, and time comparisons for the three methods. For the Lagrangian problems in Tables 2-4, the maximum number of iterations allowed at the first branch-and-bound tree node is 500; at any subsequent tree node, 150. For the results of the larger problems in Table 5, the respective iterations were increased to 15,000 and 1200. Despite having a more complicated solution procedure (as compared to the classical PMP) for each problem, EZCOMP under Lagrangian relaxation is as computationally efficient as a weighted sums model in comparisons on the 55-node data with 4 zones and 12 zones [36] while identifying

92

Ross A.

GERRARD and RICHARDL. CHURCH

alternatives that would be missed by the weighting method. Tables 2 and 3 contain results from the S-node data. For all S-node results, a very tight dual bound ((primal-dual)/dual
2600

2650

2700 wt. dist. = 2732.940 non-dominated solution

weighted distance

2750 wt. dist. = 2742.571 (ZCOMP solution)

2800

2850

2

1

0

# zonal constraint violations (Q) Fig. 3. Four zone, six facility noninferior tradeoff, 55-node data, zonal constraint set A (refer to Table 2, six-facility results).

93

Analyzing tradeoffs between zonal constraints and accessibility in facility location

wt. dist. = 3489.188

weighted distance

5000

5500 wt. dist. = 5641.191 (ZCOMP solution)

I

I

4

I

I

2 1 3 # zonal constraint violations (Q)

I

0

Fig. 4. Four zone, four facility noninferior tradeoff, S-node data, zonal constraint set B (refer to Table 2, four

Table 4. EZCOMP

No. of facilities

results

(5 and 7 facilities)

on 17-zone Uganda linear-integer

data (152~nodes), comparison programming

Pattern

No. of zonal constraint violations

Total weighed distance

27,42,64,106,130

ot

225,75 I .5

2.03 s [294,0]

22,28,63,106,130

1

187,984.8

1.61 s [262,0]

22,28,64,110.133

2

175.650.6

0.47 s [75,0]

22,28,61,114,119 22,28.37,64,110

2: 3

176,356.5 165,958.2

5.37,55,76,110

4

159,220.7

0.67 s [Ill,01

5.32,55,76,1 IO

55

157,072.3

1.27 s [215,0]

22,28,42,61,106,114,129

0t

166,776.3

4.39 s [633,15]

22,28,37,64,106,114,129

I

152.264.4

5.35 s [802,27]

Lagrangian time [No. of iterations, No. of BB nodes]

__ 3.53 s [567,1 I]

5,28,37,55,76,110,129

2

145,941.a

11.44 s [1745,61]

4.5,37,55.76,110,129

3

140.368.9

4.80 s [722,15]

5,22,37,57,64-l 10,114

4

l37v235.2

3.78 s [600,15]

5,22,32,57.64,110,114

5

135605.3

5.11 s [810,29]

134.020.9

4.10 s [665,15]

132.500.3

3.37 s [562,15]

8.22,32,57.64,110,114 8,23,32,59&t,]

TZCOMP solution. fBest TAB solution @p-median solution.

10,114

6 78

of Lagrangian

TAB time (20 runs) [No. of substitutions, No. of TAB cycles] (No. of optima, No. of runs) 3.17 s [675,60] (20/20) 3.16 s [658,60] (20120) 2.83 s [618,59] (0120)

facility

relaxation,

results).

TAB, and

LP time [No. of iterations, No. of BB nodes] 35.9min C13889.01 38.4min [14038.0] 29.4min[I1478.0]

2.95 s [569,62] (t4/20) 3.20 s [564,69] (g/20) 3.07 s [531,66] (8/20)

46.7min[l8115.2]

3.23 s [720,57] (7/20) 3.14 s [733,58] (20/20) 3.33 s [701,65] (5120) 3.37 s [678,63] (9120) 3.11 s [649,61] (5/20) 3.44 s [655,65] (l/20) 3.56 s [643,70] (6/20) 3.78 s [645,73] (12120)

37.3min[l3912.0]

19.5 min [8069,0] 20.8 min [9056.0]

43.5min[15361,0] 44.3 min [15950,6] 35.6min Cl4255.41 27.7 min [12347,0] 29.0min[ll719,6] 40.5 min [ 15460.61 27.6 min [ 10400.0]

Ross A. GERRARD and REHARD L. CHURCH

94

on the desires of the decision maker, alternative NZ or alternative N4 might be preferred. However, if a weighted sums approach (rather than our constrained method) were used, only the convex envelope defined by N’, N3, and N5 would be revealed. Potentially desirable nondominated solutions would remain obscured. To test the model further on the 55”node data, EZCOMP was applied to problems of locating 12 and 17 facilities on the 12-zone example of Fig. 5. Each zone was allotted a maximum of two sites and a minimum of one. See Table 3 for a list of the solutions. Again, the Lagrangian procedure is very competitive with TAB and with LP, even with the relatively advanced OSL software. TAB shows several instances in Table 3 of failure to identify the optimum even with twenty runs from different starting solutions.

ZONE 11

ZONE7

26

18

\

27

19

7

ZONE 5

\

2

1 ._

44

Fig. 5. Twelve zones, Swain S-node

data.

95

Analyzing tradeoffs between zonal constraints and accessibility in facility location

35

44

39

43

ZONE 2 [0,2] 38 36

46

34 37

47

ZONE 3 (0,2] 4s

30

48

32 31

ZONE 5 [O,O]

1a ZONE 10 (1,2] 127

128 126 &30

122

ZONE 12 [O,O]

ZONE 13 (0,2]

KEY: ZONE # [Pzmi”,Ppar]

N

/

Fig. 6. Seventeen zones, Uganda 152~node data.

125

123

124

Ross A. GERRARD

96

and RICHARD L. CHURCH

The Uganda data (152 nodes) presents a significantly greater computational challenge, creating an EZCOMP program of over 23,000 variables and 23,000 constraints. Table 4 shows results forthe 17-zone data (Fig. 6) for p = 5 and p = 7. The Lagrangian relaxation procedure usually takes a cpu time comparable to twenty TAB runs and is much faster than the LP procedure. The specified tolerance was a bound within 0.10% of optimality. This is a less tight bound than for the smaller 55-node problems, but if a 0.10% bound can be verified, it still results in either the optimal solution or a solution very close to optimal. A more challenging problem results when p= 15 on the 152~node, 17-zone problem (Table 5). The lower and upper bounds for each zone are given in Fig. 6. The problem is to optimize weighted distance while placing 15 facilities over 17 zones where 5 zones are required to have at least one facility, and three zones are allowed no facilities at all. In these circumstances, the zonal constraints become more difficult to meet. The performance of the Lagrangian method noticeably degrades compared to the p = 5 and p = 7 results in Table 4. The TAB and LP times do not similarly lengthen.

Table 5. EZCOMP results (15 facilities) on 17-zone Uganda data (15Xnodes), comparison of Lagrangian relaxation, TAB, and linear-integer programming TAB time (20 runs) [No. of substitutions, No. of TAB cycles] (No. of optima, No. of runs)

LP time [No. of iterations, No. of BB nodes]

Total weighted distance

Lagrangian time [No. of iterations, No. of FIB nodes]

ot

98982.2

2.06 mitt [12270.0]

5.23 s [1003,77] (13/20)

31.4min[l0532,0]

3,5,10,25,28, 37,54,57,76,90, 106,110,114,129,140

1

95590.0

2.59 min [16435,31]

5.78 s [1013,7771 (12/20)

37.7min [12914,2]

15

4,8,17,25,28, 37,54,57,69,90, 106,110,114,129,140

2

92,473.3

2.38 min [15103.31]

5.48 s [1023,77] (3/20)

37.7min[l3155.0]

15

3,5,10,25,28, 37,45,54,69,76, 90,110,114,129,137

3

90,086.4

1.55 min [99SS,O]

5.21 s [1046,81] (5120)

12.3 min [5162,0]

15

4,8,17,25,28,37, 45,54,69,90,110, 114,129,137,140 (R2) 3,5,10,25,28,37, 50,54.69,76,90, 110,114,129,137

4

88,372.5

5.64 s [1012,78] (l/20)

9.1 min [3843,6]

471

88,37&S

2.53 min [16397,31]

15

4,8,17,25,28, 33,37,45,55,59, 90,110,114,129,137

5

86494.4

2.39 min [15733,31]

4.99 s [1010,82] (l/20)

32.2min[l2338,0]

I5

4,8,17,25,28, 33,37,50,55,69, 90,l IO,114,129,137

6

84.95 I .6

1.12 min [7569,0]

6.11 s [1028,90] (S/20)

27.9 min [10179,0]

15

3,8,16,25,28,33, 37,50,55,6490, 110,114,129.137 (RI)

7

84,713.6

1.54 min [1~70,0]

6.00 s [1~7,~]

17.2 min [6317,2]

4,8,17,25,28, 33,37,50,55,69, 90.110,114,137,141

8

3,8,16,25,28, 33,37,50.55,64, 90.1 IO,114,137,141

9t

3,8,16,25,31, 33,37,50,55,64, 90,110,114,137,141

115

No. of facilities

Pattern (fiaure label)

15

3,5,25,28,37, 42,54,57,76,90

15

15

15

15

15

No. of zonal constraint violations

114/20)

84.154.7

1.29 min [8989,0]

5.67s [1007,87]

23.5 min [9013,0]

(S/20) 83,916,7

1.22 min [8547,0]

5.39 s [1010,88]

IO.7 min [5109,0]

u4/201

83,798.9

1.61 min [12609,0]

4.79s [993,76] (9/20)

19.0 min [7464,0]

tZCOMP solution. :No solution with 15 facilities and 10 violations was identified. When Q was specified as 10, the nine-violation sohition was identified, indicating that ail ten-violation solution are dominated. The results given are for Q=9. &median solution. 1lIndicates non-optimal Lagrangian solution, although within 0.10% bound. The patterns labeled R’ and RZ are both unsupported nondominated solutions and are graphed in Fig. 7.

97

Analyzing tradeoffs between zonal constraints and accessibility in facility location

80000

(PMP solution) .

non-dominated

85000

solution

R2 un.wqGgd

non-dominated

\

9oooo

weighted distance 95000

“\u (ZCOMP

100000

105000

L-‘_-11

10

9

8

7

6

5

4

3

2

solulion)

1

0

# zonal constraint violations (Q) Fig. 7. Seventeen zone, fifteen facility noninferior tradeoff, 152-node data (refer to Table 5).

Extensive testing of the Lagrangian problems was done to find which parameters would result in the most favorable solution times. Lagrangian relaxation still retains the virtue of verification of the optimum (or very close to optimum) in considerably less time than that required by LP. For the results in Table 5, a single run of TAB was used to provide an initial primal solution for the Lagrangian procedure. When branch and bound was invoked, the resulting tree always contained 2p+ 1 (31) nodes, indicative of the type of branch and bound tree discussed above. Finally, note that for p=4, the Lagrangian solution terminates at slightly worse than optimum (but within the required 0.10% bound). Figure 7 graphs the points shown in Table 5 as a tradeoff curve. Note the presence of two unsupported nondominated solutions in Fig. 7. It should be restated that few multiobjective location models are capable of identifying these ‘gap point’ solutions, as other models are almost exclusively solved by using the weighting method. 6. FINAL

REMARKS

The zonally constrained median problem is a recent development in median location modeling. Special constraints are used to make distributional issues endogenous to the PMP. Recent papers by ReVelle and Elzinga [lS]; Berman et al. [16]; and Church L-143demonstrate interest in this approach. Facility location can be a controversial issue, whether the services to be provided are desirable or undesirable. ZCOMP [14] has the capability to mitigate such controversy by distributing facilities among sets of sites. These sets may be defined on the basis of essentially any criterion, but they are most easily interpreted as geographical zones, jurisdictions, or neighborhoods. The zonal constraints are maintained while retaining an objective of overall geographic efficiency. Unfortunately, imposing constraints to achieve regional equity generally leads to decreases in system-wide efficiency. More flexibility can be provided for locational decisions and increases in efficiency can result if one is willing to relinquish some of the requirements for spatial distribution. This paper has introduced a general multiobjective model (EZCOMP) which can be used to analyze the conflict between zonal constraints and accessibility in facility location. Of the recently published

Ross A. GERRARDand RICHARDL. CHURCH

98

works, it is closest to Church [14] in structure but is unique in that it can analyze the problem in the multiobjective context just described. Compared to ZCOMP, EZCOMP can identify a far greater variety of potential solutions (noninferior set) for a given set of zones and zonal constraints. The constraint method applied to EZCOMP was used to generate noninferior combinations of overall accessibility (efficiency) and satisfaction of regional location demands (distribution). The model has proven capable of identifying the complete noninferior set, including unsupported nondominated solutions. Lagrangian relaxation, the Teitz and Bart heuristic (in modified form) and linear-integer programming were shown to be viable solution procedures. However, LP/IP on OSL was almost always at least ten times slower than Lagrangian relaxation, was commonly 500 times slower, and at times was over 1000 times slower. In addition, a unique solution procedure for a Lagrangian dual problem was demonstrated. This technique has expanded the application of Lagrangian relaxation to specially constrained cases of the p-median problem. This process may offer a basis for further application of the constraint method for solving other multiobjective facility location problems. Acknowledgements-We would like to thank the IBM Corporation for making possible the use of the OSL software and an RS/6000 workstation through generous support to the National Center for Geographic Information and Analysis at the University of California at Santa Barbara.

REFERENCES 1. S. L. Hakimi, Optimum locations of switching centers and the absolute centers and medians of a graph. Ops Res. 11, 450-459 (1964). 2. S. L. Hakimi, Optimum distribution of switching centers in a communication network and some related graph theoretic problems. Ops Res. 13, 462-475 (1965). 3. C. S. ReVelle and R. W. Swain, Central facilities location. Geogr. Analysis 2, 30-42 (1970). 4. K. Reinert, T. Miller and H. Dickerson, A location-assignment model for urban snow and ice control operations. J. Urb. analysis public Mgmf 8, 175-191 (1985). 5. C. Partoune and D. Peeters, Un modele de localisation des bureaux de poste. Annalspublic coop. econ. 51,53-68 (1980). 6. P. R. Odell, K. E. Rosing and H. Beke-Vogelaar, Optimising the oil pipeline system in the U.K. sector of the North Sea. Energy Policy 4, 50-55 (1976). 7. P. Rojeski and C. S. ReVelle, Central facilities location under an investment constraint. Geogr. Analysis2,343-360 (1970). 8. R. W. Swain, Central facilities location with multiple planning periods. Department of Industrial and Systems

Engineering, Ohio State University, Columbus, Ohio (1976). 9. P. B. Mirchandani and A. R. Odoni, Locations of medians on stochastic networks. Transport. Sri. 13 86-97 (1979). 10. S. C. Narula and U. I. Ogbu, An hierarchical location-allocation problem. Omega 7, 137-143 (1979). 11. J. R. Weaver and R. L. Church, A median location model with nonclosest facility service. Transport. Sci. 19,58-74 (1985). 12. J. R. Weaver and R. L. Church, The formal and computational relationship of the supporting median problem to the p-median problem. Transport. Res. 21B, 323-329 (1987). 13. R. W. Swain, A decomposition algorithm for a class of facility location problems. Ph.D. thesis, Cornell University, Ithaca, N.Y. (1971). 14. R. L. Church, The regionally constrained p-median problem Geogr. Analysis 22, 22-32 (1990). 15. C. S. ReVelle and D. J. Elzinga, An algorithm for facility location in a districted region. Environ. P/an. B 16,41-50 (1989). 16. 0. Berman, D. Einav and G. Handler, The zone-constrained location problem on a network. Europ. J. opl. Res. 53, 14-24 (1991). 17. M. Mandell, Modelling effectiveness-equity tradeoffs in public service delivery systems. Mgmt Sci. 37,467-482 (1991). 18. J. Cohon, Multiobjective Programming and Planning, pp. 100-105, 115. Academic Press, New York, N.Y. (1978). 19. R. E. Steuer, Multiple Criteria Optimization: Theory, Computation, and Application, pp. 431-433. Wiley, New York,

N.Y. (1986). 20. R. Solanki, Generating the noninferior set in mixed integer biobjective linear programs: an application to a location problem. Computers Ops. Res. 18, 1-15 (1991). 21. R. T. Rockafellar, Convex Analysis. Princeton University Press, Princeton, N.J. (1970). 22. R. T. Rockafellar, The Theory of Subgradients and Its Applications to Problems of Optimization: Convex and Nonconvex Functions. Heldermann, Berlin (1981). 23. S. C. Narula, U. I. Ogbu and H. M. Samuelsson, An algorithm for the p-median problem. Ops. Res. 25,709-713 (1977). 24. N. Christofides and J. E. Beasley, A tree search algorithm for the p-median problem. E. J. opl Res.10,196-204 (1982). 25. P. Hanjoul and D. Peeters, A comparison of two dual-based procedures for solving the p-median problem. E. J. opl Res. 20, 387-396 (1985). 26. G. Comuejols, M. L. Fisher and G. L. Nemhauser, Location of bank accounts to optimize float: an analytic study of exact and approximate algorithms. Afgmr Sci. 23, 789-810 (1977). 27. M. Held and R. M. Karp,The traveling salesman problem and minimum spanning trees. Ops Res. 18,1138- 1162 (1970). 28. E. Balas and N. Christofides, A restricted Lagrangian approach to the traveling salesman problem. Math/ Program. 21, 19-46 (1981).

29. A. Geoffiion and R. McBride, Lagrangian relaxation applied to capacitated facility location problems. AIIE Trans.

10,4O-47 (1978).

AnaIyzing tradeoffs between zonal constraints and accessibility in facility location

99

30. M. L. Fisher, The Lagrangian relaxation method for solving integer programming problems. Mgmt Sci. 27, 1-18 (1981). 31. M. L. Fisher, An applications oriented guide to Lagrangian relaxation. Interfaces 15, lo-21 (1985). 32. M. Held, P. Wolfe and H. Crowder, Validation of subgradient optimization. Mathl Program. 6, 62-88 (1974). 33. M. B. Teitz and P. Bart, Heuristic methods for estimating the generalized vertex median of a weighted graph. Ops Res. 16,955961 (1968). 34. K. E. Rosing, E. L. Hi&man and H. Rosing-Vogelaar, A note comparing optimal and heuristic solutions to the p-median problem. Geogr. AnaJysis 11,86-89 (1979). 3.5. D. Migereko, An analysis of the coffee cooperative marketing system in Busoga, Uganda: transportation and facilities location. M.A. thesis in Geography, University of California, Santa Barbara, Calif. (1983). 36. R. A. Gerrard, Accounting for political imperatives in regional facility location. M.A. thesis in Geography, University of California, Santa Barbara, Calif. (1989).