Single-facility huff location problems on networks

Ann Oper Res (2014) 222:175–195 DOI 10.1007/s10479-013-1445-x Single-facility huff location problems on networks Rafael Blanquero ·Emilio Carrizosa · Amaya Nogales-Gómez ·Frank Plastria Published online: 17 September 2013 © Springer Science+Business Media New York 2013 Abstract Huff location problems have been extensively analyzed within the field of competitive continuous location. In this work, two Huff location models on networks are addressed, by considering that users go directly to the facility or they visit the facility in their way to a destination. Since the problems are multimodal, a branch and bound algorithm is proposed, in which two different bounding strategies, based on Interval Analysis and DC optimization, are used and compared. Computational results are given for the two bounding procedures, showing that problems of rather realistic size can be solved in reasonable time. Keywords Huff location models ·Location on networks ·DC optimization ·Interval analysis 1 Introduction In this paper we address two competitive location models (Blanquero and Carrizosa 2009a; Drezner and Drezner 2004; Fernández et al. 2007;Huff1964,1966;Plastria2001)ona This work has been partially supported by IMUS (Mathematics Institute of University of Seville) and by projects MTM2012-36163, Ministerio de Economía y Competitividad, Spain, P11FQM7603, P10TIC6064 and FQM-329, Junta de Andalucía, Spain, all with EU ERF Funds. R. Blanquero ·E. Carrizosa ·A. Nogales-Gómez (B) Departamento de Estadística e Investigación Operativa Facultad de Matemáticas, Universidad de Sevilla, Seville, Spain e-mail: [email protected] R. Blanquero e-mail: [email protected] E. Carrizosa e-mail: [email protected] F. Plastria MOSI, Vrije Universiteit Brussel, Brussels, Belgium e-mail: [email protected] 176 Ann Oper Res (2014) 222:175–195 network (Dearing and Shier 1983; Labbé et al. 1995). Let N=(V,E) beanetwork,with node set Vand edge set E. The length of the edge e∈Eis given, and it is denoted by le.The distance between two nodes ai,ajis calculated as the shortest path (Labbé et al. 1995) from aito aj. For each e∈E, with end nodes ai,aj, we identify each x∈[0,le]with the point in the edge eat distance xfrom aiand distance le−xfrom aj. This way, we obtain that, for any vertex ak∈V, the distance d(x,ak)from xto akis, as a function of x,aconcave piecewise linear function, given by: d(x, ak)=minfiak(x),fjak(x) fiak(x) =x+d(ai,ak)f jak(x) =(le−x)+d(aj,ak) (1) We refer the reader to Labbé et al. (1995) for a comprehensive introduction to location models on networks. The remainder of this paper is structured as follows. Two different location models on networks are considered. The first model (Sect. 2) is the classic Huff location model, as addressedinBermanetal.(2011), in which customers perceive the facility attractiveness in terms of their distance to the facility, while the second model (Sect. 3) is new, called the Huff origin-destination (OD) trip model. In the OD trip model the facility attraction is a function of the length of the shortest path from the origin to the destination through the facility. In Sect. 4a branch and bound algorithm is designed to solve both problems and two different procedures to obtain bounds are presented. Computational results are given in Sect. 5, comparing the two bounding strategies implemented. Some concluding remarks are presented in Sect. 6. 2 Huff location model In this model, the finite set Vof vertices of the network represents users, asking for a certain service. Each user a∈Vhas demand ωa≥0. Such demand is being patronized by different existing facilities, located at points x1,...,xron the network, so that the demand captured by facility at xifrom user ais inversely proportional to a positive nondecreasing function of the distance d(a,xi)from the user at ato the facility at xi. In other words, the demand captured by the facility at xifrom the user at ais given by ωa 1/ϕa(d(a, xi)) r j=11/ϕa(d(a, xj)).(2) The usual choice for each ϕahas the form ϕa(d) =dλa.Whenλa=2foralla∈V,wehave the so-called gravitational model. A new firm is entering the market, by locating one facility at some point on the network. This perturbs market share, since the new facility at xwill capture a demand from a∈V equal to ωa 1/ϕa(d(a, x)) 1/ϕa(d(a, x)) +r j=11/ϕa(d(a, xj)).(3) Here ϕais assumed to be non-negative, non-decreasing and twice continuously differentiable in R+. Expression (3) must be carefully considered when ϕa(0)=0. Indeed, in such case, if some xjexists with xj=a∈V, then the full demand of awill be captured by xj. Ann Oper Res (2014) 222:175–195 177 Hence, such awill not be taken into account, and thus we assume in what follows without loss of generality that xj/∈V,j=1,...,r. The goal of the entering firm is the maximization of its market share. This is written as the following optimization problem: max x∈[0,le],e∈E a∈V ωa 1/ϕa(d(a, x)) 1/ϕa(d(a, x)) +r j=11/ϕa(d(a, xj)).(4) Defining for each a∈Vthe positive constant βa, βa= r  j=1 1 ϕa(d(a, xj)),(5) it follows that problem (4) can be rewritten as max x∈[0,le],e∈EF(x) (6) with Fdefined as F(x)= a∈V ωa 1 1+βaϕa(d(a, x)).(7) This problem has been addressed in Berman et al. (2011), in which a branch and bound algorithm with interval analysis bounds is proposed. 3 Huff OD trip model In this model, we have the set {{a,b},a ∈V,b ∈V}, of origin-destination pairs. The trip from origin ato destination bwill imply a demand ωab >0. Consumers in their way from origin ato destination bwill stop at one facility; they choose among the rexisting facilities, x1,...,xr, and the new facility at xas in the model described in Sect. 2: the demand from a user in his way from ato bcaptured by each facility is inversely proportional to a positive nondecreasing function of the length of the path from the origin to the destination via the facility. In other words, the demand captured by xfrom users in their trip from origin ato destination bis given by ωab 1/ϕab(d(a,x) +d(x,b)) 1/ϕab(d(a,x) +d(x,b))+r j=11/ϕab(d(a,xj)+d(xj,b)).(8) As in the other model, the goal of the entering firm is the maximization of its market share. This is written as the following optimization problem: max x∈[0,le],e∈E a,b∈V ωab 1/ϕab(d(a,x) +d(x,b)) 1/ϕab(d(a,x) +d(x,b))+r j=11/ϕab(d(a,xj)+d(xj,b)).(9) Defining for each a,b ∈Vthe positive constant βab, βab = r  j=1 1 ϕab(d(a, xj)+d(xj,b)), 178 Ann Oper Res (2014) 222:175–195 it follows that problem (9) can be rewritten as max x∈[0,le],e∈EF(x) with Fdefined as F(x)= a,b∈V ωab 1 1+βabϕab(d(a,x) +d(x,b)).(10) We see that both the location and OD trip models yield an optimization problem of the same form, namely, (6), with a rather similar function F,asgivenby(7)and(10) respectively. Both functions Fcan be written in the form F(x)= δ∈Δ ωδ 1 1+βδϕδ(dδ(x)).(11) where Δ=Vand,for δ∈Δ, dδ(x) =d(a,x) for the location model, and Δ={{a,b},a,b∈V}and,for δ∈Δ, dδ(x) =d(a,x)+d(x,b) for the OD trip model. 4 Solving the models 4.1 Multimodality The optimization problems described in Sect. 2and 3are, in general, multimodal, and standard optimization methods get stuck at local optima. This can be seen in simple examples, even when the network is a segment and is illustrated in the following examples for the Huff location problem introduced in Sect. 2. First, data for one problem with two users and r=1 facility on a segment were randomly generated. The objective function Fof such instance is plotted in Fig. 1(left). One can see that the problem is bimodal. 100 runs of a local search procedure starting with a random point were performed. In Fig. 1(right) one can see the histogram of the objective values provided by the optimizer: below a 50 % of the runs yielded the global optimum, whereas the remaining runs stopped at the local not globally optimal solution. Another instance, with 500 user and 100 facilities was also generated. In Fig. 2the problem is shown to be multimodal, and just below a 14 % of the runs solved with the optimizer yielded the global optimum. In the right part of Figs. 1–2the xaxis is normalized, so that 1 corresponds to the best found objective value, z∗,andabaratx∈[0,1]indicates that an objective value x·z∗was found. Since multimodality appears even in the simplest cases, and local search procedures can get stuck at very bad local optima, global optimization tools are needed if the global optimum is sought. Ann Oper Res (2014) 222:175–195 179 Fig. 1 Users =2, facilities =1, λa=1∀a Fig. 2 Users =500, facilities =100, λa=20 ∀a 4.2 Algorithm We will use a branch and bound algorithm (Drezner and Suzuki 2004;Hansenetal.1985; Plastria 1992) to solve the problems under consideration. We first outline the algorithm that finds an optimal solution within a relative accuracy of ε, and later give the details on the bounding process, which will exploit monotonicity properties or the fact that the objective function is DC on each edge (Blanquero and Carrizosa 2009b; Horst and Thoai 1999;Tuy 1995; Tuy et al. 1995). We remind the reader that a function his DC if it can be written as h=h+−h−,whereh+,h−are both convex; the expression h+−h−is called a DC decomposition of h. In the description of the algorithm, we have Δ=Vfor the Huff location model and Δ={{a,b},a,b∈V}for the Huff OD trip model. Phase 1: Initialization – Fix the required accuracy ε>0. –SetLB =0. – Compute the all-pairs distance matrix. – Calculate βδ,∀δ∈Δ. 180 Ann Oper Res (2014) 222:175–195 – Set the list Hof remaining segments as empty. Phase 2: Prepare the list of segments – Consider the edge eas segment with its nodes as the segment vertices. – The value of the objective function is evaluated at the segment midpoint. If the value is greater than LB,thenLB is updated to such value and the midpoint stored as incumbent. – Calculate an upper bound for the segment e,UB(e). – In case UB(e) ≥LB ·(1+ε) insert einto H. Phase 3: Branch and bound process. Repeat as long as no stop is reached: – Select from Hthe highest upper bound segment, UBmax, with LB as ε-optimal value and the incumbent as an ε-optimal solution. If UBmax ≤LB ·(1+ε), stop the algorithm with LB as the solution. – The highest upper bound segment, UBmax, is selected for a split at its midpoint into two smaller segments. – The value of the objective function at the midpoint of the two small segments is calculated. If any of these values is greater than LB,thenLB is updated and all segments from H whose upper bound is lower than LB are discarded. – An upper bound for each small segment is calculated. – All segments whose upper bound is greater than LB ·(1+ε) are added to H. In the algorithm above, the segment midpoint yielding the best upper bound is given as ε-optimal solution. Observe that the full list of segments in Hcontain all ε-optimal solutions, and can thus be used in a two-phase process, as suggested in the GBSSS algorithm (Plastria 1992). Let us detail the algorithm steps previously outlined. To calculate the all-pairs distance matrix we use the Floyd algorithm (Bertsekas 1998), which uses the edge length matrix to build recursively the distance matrix. The computation of the bounds requires more detail. The algorithm needs the calculation of an upper bound, UB(s), for each segment s. We present two procedures, one based on Interval Analysis, and the other in properties of DC functions. 4.3 Interval analysis bound When the distance dδ(x) decreases, the market share as given by the objective function (11) increases. Hence we obtain an upper bound on the market share F(x)for any location xon asegments=[x0,x1]⊂e∈Eby replacing dδ(x) by the lowest possible of these distances on the segment. Defining therefore such lowest distance for the location model as in Berman et al. (2011)dδ(s) =min{d(a,x0), d(a, x1)}for δ=a, and for the OD trip model as dδ(s) = min{d(a,x0)+d(b,x0), d(a, x1)+d(b,x1)}for δ={a,b}, we obtain the Interval Analysis bound UBIA(s) = δ∈Δ ωδ 1 1+βδϕδ(dδ(s)) (12) 4.4 DC bound An upper bound obtained making use of the fact that the objective function is DC on each edge exploits the following properties: Ann Oper Res (2014) 222:175–195 181 Proposition 4.1 Let I⊂Rbe an interval.Let d:I→Rbe a concave function on I,and let g:R→Rbe DC,with a DC decomposition given by g(x) =g+(x) −g−(x),with both g+ and g−non-increasing functions.Then,the function f:I→Rdefined as f(x)=g(d(x)) is DC on Iand a DC decomposition is given by f(x)=f+(x) −f−(x),where f+(x) = g+(d(x)) and f−(x) =g−(d(x)). Proof The proof follows directly from the fact that the compositions g+(d(x)) and g−(d(x)) are also convex functions.  Remark 4.1 Proposition 4.1 makes use of a function gwhich can be written as the difference of two convex functions, g+and g−. Since such convex functions are also non-increasing, it turns out that gbelongs to a subclass of DC functions, namely DCM functions, as introduced in Blanquero and Carrizosa (2009a), which are those functions expressed as the difference of two convex monotonic functions, as g+,g−are. See Blanquero and Carrizosa (2009a)for further properties. We are now in position to give a bound for Fon an edge exploiting the fact that Fis DC. Let dδ(x) be the concave function given in (1). Assuming ϕδ(d) =dλ δ,(11) can be rewritten as F(x)= δ∈Δ ωδ 1 1+βδdλ δ(x) (13) Let us define a simpler function: g(t) =1 1+βtλ(14) The following DC decomposition is known for function g(Bello et al. 2011): c=λ−1 (λ +1)β 1/λ g+(t) =g(c) +g(c)(t −c) if t≤c g(t) if t>c g−(t) =g(c) +g(c)(t −c) −g(t) if t≤c 0ift>c That means that g(t) =g+(t) −g−(t). Applying Proposition 4.1 we have that a DC decomposition for Fas defined in (13)is: F+ δ(x) =ωδg+dδ(x) F− δ(x) =ωδg−dδ(x) F(x)= δ∈Δ F+ δ(x) − δ∈Δ F− δ(x) = δ∈ΔF+ δ(x) −F− δ(x)(15) 182 Ann Oper Res (2014) 222:175–195 To construct an upper bound, UBDC, first we obtain a convex minorant of F− δ(x) as in Blanquero and Carrizosa (2009a): F− δ(x) ≥F− δ(x0)+ξδ(x −x0) F(x)≤ δ∈ΔF+ δ(x) −F− δ(x0)−ξδ(x −x0)for ξδ∈∂F− δ(x0), (16) where ∂F− δ(x0)denotes the set of subgradients of F− δat x0(Tuy 1995). Define for each δ∈Δthe function from (16): U(x)= δ∈ΔF+ δ(x) −F− δ(x0)−ξδ(x −x0)for ξδ∈∂F− δ(x0)(17) An upper bound for a segment sis obtained: UBDC(s) =maxU(v1), U(v2) v1,v2being vertices of s. 5 Computational results The algorithm described in Sect. 4was programmed in an Intel Fortran Compiler XE 12.0 and executed with an Intel Core i7 computer with 8.00 Gb of RAM memory at 2.8 Ghz. The solutions were found within an accuracy of 10−10. The Huff location model was first tested on the 55-node and 134-edge Swain’s (1971) network (Serra et al. 1999; Marianov and Serra 2002), see the Appendix,Table6, with both bounding strategies. Several instances of the problem were generated using different values for the number rof existing facilities, ranging from low-saturated markets (r=10 % of the number of edges of the network, |E|) to high-saturated markets (r≤90 %|E|). For each value of r,10 different problems were solved. Each problem is obtained by randomly and independently generating the demands (each vertex of the network is assumed to have a demand uniformly distributed in the interval (0,1)) and the location of the existing facilities. To generate the r facility locations, redges are randomly chosen with replacement; on each edge, the facility location is generated following a uniform distribution. The results are shown in Table 1. The percentage rof existing facilities is shown in the first column. Then, it is reported the minimum, maximum, mean and standard deviation (std) for the number of iterations, i.e., the number of executions of Phase 3, the branch and bound (B& B) list size and the CPU time. The branch and bound list size is the maximum size of the data structure used for storage reached during the algorithm execution. In all cases, the best solution is found in less than 0.02 seconds for DC bounds, and 0.75 seconds for Interval Analysis bounds. It is remarkable that DC bounds lead to a very stable procedure, as can be seen for all values of rin memory requirements as well as computational time. On the other hand, when using Interval Analysis bounds one can see very extreme cases. There is always a huge difference between the minimum and maximum for the number of iterations, branch and bound list size and time. This means that in the ten runs that are solved for each value of r, Interval Analysis bounds are quite erratic for problems of the same difficulty. Therefore, we have an algorithm that when using DC bounds Ann Oper Res (2014) 222:175–195 183 Table 1 Results of the 55-node and 134-edge Swain’s network for the Huff location model r(%|E|)Iterations B&B list size Time min max mean±std min max mean±std min max mean±std DC 10 111 196 144.40±34.68 9 61 31.50±17.25 0.00 0.02 0.00±0.00 20 80 193 140.10±34.32 15 83 42.70±22.76 0.00 0.00 0.00±0.00 30 51 209 117.80±47.86 25 96 66.70±22.42 0.00 0.00 0.00±0.00 40 50 256 176.40±55.33 15 107 84.80±30.51 0.00 0.02 0.00±0.00 50 81 290 172.20±67.40 74 126 107.50±14.92 0.00 0.02 0.00±0.00 60 121 317 166.50±60.54 72 121 110.10±16.09 0.00 0.02 0.00±0.00 70 94 323 190.00±64.79 113 131 123.80±6.36 0.00 0.02 0.00±0.00 80 145 332 209.70±58.43 117 133 123.70±5.62 0.00 0.02 0.00±0.01 90 137 296 211.40±55.51 111 129 121.40±5.66 0.00 0.02 0.00±0.00 Interval 10 158 318 223.70±50.81 26 120 64.20±31.72 0.00 0.02 0.00±0.00 20 185 178915 18100.40±56504.50 29 58170 5890.10±18369.30 0.00 0.59 0.06±0.19 30 203 70911 17779.60±25278.44 75 25465 6219.20±8904.76 0.00 0.20 0.05±0.07 40 243 243922 24639.40±77048.06 96 84425 8551.40±26659.27 0.00 0.75 0.08±0.24 50 187 37594 9132.40±13971.49 118 13306 3143.50±4794.67 0.00 0.11 0.03±0.05 60 186 29111 7734.70±12130.03 118 9927 2703.30±4174.80 0.00 0.09 0.03±0.04 70 212 10748 3029.00±4148.73 138 3489 1038.00±1349.47 0.00 0.03 0.01±0.02 80 172 1950 866.60±677.23 127 648 310.10±192.99 0.00 0.02 0.00±0.00 90 224 22520 4243.00±7225.49 147 7248 1442.20±2338.18 0.00 0.09 0.02±0.03 190 Ann Oper Res (2014) 222:175–195 Table 4 (Continued) Network Nodes Edges DC bounds Interval analysis bounds Iterations B&B list size Time Iterations B&B list size Time mean±std mean±std mean±std mean±std mean±std mean±std UR752 580 1735 232.10±33.86 160.20±11.69 17.04±0.79 330.80±70.44 160.20±11.69 16.15±1.04 UR762 593 2089 266.80±15.37 202.00±13.52 21.13±0.40 489.40±32.67 202.00±13.52 21.10±0.35 UR132 605 1122 120.50±54.90 24.60±2.17 11.48±1.48 186.00±28.42 24.60±2.17 10.93±0.41 UR735 662 1200 78.90±13.30 40.00±1.49 13.37±0.42 120.20±4.89 40.00±1.49 12.79±0.15 UR142 709 1815 115.00±16.57 94.70±13.98 23.59±0.65 179.30±30.44 94.70±13.98 22.73±0.62 UR745 713 1616 146.90±31.56 76.30±2.91 23.02±1.27 172.90±26.59 76.30±2.91 20.52±0.47 UR755 724 1966 213.60±28.51 252.60±24.49 30.11±1.20 290.40±53.88 252.60±24.49 27.37±1.11 UR765 741 2278 213.70±9.41 111.20±3.77 35.85±0.29 250.60±11.48 111.20±3.77 31.92±0.33 UR737 744 1315 107.60±17.78 49.40±7.18 20.64±0.81 137.70±8.43 49.40±7.18 19.15±0.22 UR747 745 1659 122.60±31.26 56.00±2.31 25.11±1.31 163.90±15.77 56.00±2.31 23.22±0.33 UR757 748 1969 155.90±28.57 165.10±16.54 30.83±1.30 256.00±22.63 165.10±16.54 29.70±0.50 UR767 749 2314 274.60±58.70 236.60±21.26 39.65±2.59 333.50±22.87 236.60±21.26 34.84±0.51 UR152 766 2390 208.30±39.61 155.40±11.81 39.98±1.88 268.30±16.64 155.40±11.81 36.43±0.38 UR162 802 2897 407.20±107.72 288.50±31.36 62.94±5.80 504.60±39.82 288.50±31.36 55.42±1.04 UR135 892 1619 49.80±9.35 33.00±0.82 34.40±0.69 132.50±16.64 33.00±0.82 36.79±0.61 UR145 929 2117 136.40±21.10 115.20±6.73 57.26±1.59 243.80±18.41 116.30±6.11 57.94±0.74 UR155 975 2680 174.70±31.59 195.80±21.34 83.70±2.69 242.80±25.37 195.80±21.34 85.42±4.06 UR137 980 1744 94.60±13.78 48.90±5.99 53.14±1.17 120.30±3.09 48.90±5.99 54.90±0.28 UR165 980 3068 155.10±26.17 181.40±18.46 93.71±2.34 236.00±37.40 181.40±18.46 98.32±1.94 UR147 996 2254 89.00±13.47 132.30±7.53 69.48±1.34 169.50±10.26 132.30±7.53 71.75±2.24 UR157 1000 2690 154.00±21.55 139.70±12.14 83.01±1.90 256.30±13.62 139.70±12.14 82.30±0.99 UR167 1000 3083 229.10±33.85 254.70±24.43 100.38±3.03 313.70±42.59 254.70±24.43 98.70±2.65 Ann Oper Res (2014) 222:175–195 191 Table 5 Results of test instances for the Huff OD trip model with 90 %|E|facilities Network Nodes Edges DC bounds Interval analysis bounds Iterations B&B list size Time Iterations B&B list size Time mean±std mean±std mean±std mean±std mean±std mean±std KROB150G 150 296 88.50±3.66 62.70±4.30 0.22±0.00 137.80±3.71 62.70±4.30 0.21±0.01 KROA150G 150 297 92.70±19.97 41.70±6.00 0.22±0.02 167.20±10.12 41.70±6.00 0.22±0.01 PR152G 152 296 71.90±3.00 58.50±9.69 0.20±0.00 122.60±2.01 58.50±9.69 0.20±0.01 RAT195G 195 336 55.90±4.82 93.90±6.42 0.36±0.01 119.10±11.12 93.90±6.42 0.38±0.02 KROB200G 200 386 90.90±22.34 34.10±3.96 0.47±0.04 144.90±5.63 34.10±3.96 0.45±0.01 KROA200G 200 392 78.60±3.27 37.30±1.89 0.45±0.01 129.60±2.17 37.30±1.89 0.44±0.01 TS225G 225 306 65.80±2.49 57.10±3.28 0.46±0.02 120.20±0.42 57.10±3.28 0.46±0.01 UR532 298 597 97.40±4.01 59.10±5.22 1.53±0.02 155.60±4.38 59.10±5.22 1.41±0.01 UR542 343 862 161.40±8.77 150.30±6.55 3.11±0.04 240.00±14.00 150.30±6.55 2.83±0.05 UR552 388 1135 113.80±4.32 99.40±4.48 4.74±0.06 211.60±3.75 99.40±4.48 4.60±0.03 UR562 416 1403 204.00±10.14 195.20±8.47 7.69±0.09 356.40±13.94 199.10±12.56 7.29±0.09 UR732 452 915 79.90±3.38 45.10±2.42 5.44±0.08 130.20±2.15 45.10±2.42 5.22±0.04 UR535 458 812 81.60±4.53 44.40±4.62 5.19±0.08 137.80±9.73 44.40±4.62 4.98±0.09 UR545 476 1104 86.80±14.46 91.40±4.35 7.42±0.23 148.20±8.11 91.40±4.35 7.21±0.10 UR555 490 1305 156.00±29.30 133.00±23.67 10.30±0.48 223.20±18.90 133.00±23.67 9.49±0.17 UR537 493 868 66.60±2.63 39.10±6.40 6.44±0.15 127.20±3.33 39.10±6.40 6.25±0.03 UR565 496 1513 153.50±9.62 128.70±2.98 11.96±0.19 281.70±28.41 136.30±5.58 11.43±0.24 UR547 498 1112 79.40±3.50 77.50±3.92 8.31±0.09 139.90±4.36 77.50±3.92 7.90±0.04 UR557 498 1310 132.10±10.12 122.90±13.92 10.40±0.18 270.80±18.26 122.90±13.92 10.26±0.17 UR567 499 1426 112.00±2.67 95.10±2.77 10.84±0.03 187.50±5.25 95.10±2.77 10.26±0.04 UR742 538 1325 162.20±13.79 162.30±9.13 13.17±0.29 246.00±16.07 162.30±9.13 12.39±0.21 192 Ann Oper Res (2014) 222:175–195 Table 5 (Continued) Network Nodes Edges DC bounds Interval analysis bounds Iterations B&B list size Time Iterations B&B list size Time mean±std mean±std mean±std mean±std mean±std mean±std UR752 580 1735 194.70±11.72 160.10±7.37 20.15±0.36 361.80±73.69 160.10±7.37 19.55±0.96 UR762 593 2089 303.80±77.63 191.50±15.41 27.38±2.02 445.70±38.44 191.50±15.41 24.96±0.57 UR132 605 1122 73.50±13.33 26.10±0.32 12.95±0.37 185.20±28.54 26.10±0.32 13.41±0.60 UR735 662 1200 72.80±3.55 43.60±1.51 16.70±0.15 123.50±4.17 43.60±1.51 16.46±0.81 UR142 709 1815 105.50±3.89 94.60±6.74 29.68±0.23 157.00±6.02 94.60±6.74 27.89±0.32 UR745 713 1616 127.60±44.94 81.50±2.59 28.40±1.88 188.00±4.69 81.50±2.59 25.80±0.18 UR755 724 1966 211.70±16.17 259.10±6.17 38.18±0.69 308.20±28.83 259.10±6.17 34.46±0.62 UR765 741 2278 135.80±4.73 111.90±3.96 41.83±0.35 258.40±15.54 111.90±3.96 40.24±0.58 UR737 744 1315 69.50±2.22 54.00±2.67 24.32±0.25 142.20±2.57 54.00±2.67 24.75±0.48 UR747 745 1659 86.50±8.45 56.90±2.08 30.30±0.31 162.50±14.03 56.90±2.08 29.43±1.03 UR757 748 1969 145.80±8.12 172.90±4.18 38.24±0.35 288.10±35.09 172.90±4.18 37.58±1.21 UR767 749 2314 195.00±13.23 241.50±10.74 46.02±0.75 327.40±22.74 241.50±10.74 44.82±1.23 UR152 766 2390 178.80±26.39 160.30±3.50 49.72±1.31 267.00±5.29 160.30±3.50 46.27±1.14 UR162 802 2897 328.10±22.26 324.70±14.35 74.00±1.49 550.80±38.55 324.70±14.35 68.54±1.19 UR135 892 1619 50.50±10.22 33.50±0.53 42.73±1.64 136.60±10.36 33.50±0.53 45.95±0.72 UR145 929 2117 133.10±5.65 127.50±3.69 67.30±0.42 266.10±24.08 127.50±3.69 71.57±1.10 UR155 975 2680 171.90±7.31 204.00±16.49 99.39±0.66 266.40±15.06 204.00±16.49 98.18±0.90 UR137 980 1744 72.80±3.22 51.00±4.08 61.46±0.28 124.60±4.20 51.00±4.08 62.38±0.52 UR165 980 3068 161.20±8.42 199.10±7.92 111.06±0.75 235.40±4.81 199.10±7.92 111.63±0.68 UR147 996 2254 91.40±7.89 135.60±4.48 82.76±0.73 171.10±3.21 135.60±4.48 86.18±0.51 UR157 1000 2690 136.50±4.28 136.10±2.18 98.28±0.41 247.50±7.53 136.10±2.18 99.93±0.80 UR167 1000 3083 216.80±7.25 262.90±9.26 118.59±0.77 327.70±19.31 262.90±9.26 116.10±1.34 Ann Oper Res (2014) 222:175–195 193 The computational experience reported shows that large networks can be successfully handled with both bounding procedures where the DC procedure seems to be more stable in both time and memory requirements. Extensions of these problems to the multifacility case deserve further study. Appendix Table 6 The Swain data set (Marianov and Serra 2002; Serra et al. 1999) Initial node Final node Edge length Initial node Final node Edge length Initial node Final node Edge length 123.1623 13 47 4.4721 25 49 7.0000 152.0000 14 16 10.2956 26 36 8.6023 184.4721 14 22 12.0416 27 38 9.2195 113 2.2361 14 27 9.8489 27 54 7.2111 14310.0000 15 18 6.4031 29 31 3.6056 144 4.1231 15 25 7.2801 30 33 4.0000 234.4721 15 31 7.0000 30 45 5.0000 243.0000 15 36 6.7082 32 33 5.0000 283.1623 15 41 6.3246 32 38 4.4721 242 2.8284 15 42 6.3246 32 45 4.0000 374.2426 16 22 6.7082 34 43 6.0828 383.1623 16 27 6.4031 34 45 3.0000 319 6.4031 16 32 6.3246 35 36 7.6158 330 5.0000 16 33 6.7082 35 48 7.6158 331 6.0828 16 38 7.2111 36 48 6.0000 334 5.3852 17 22 5.0990 37 47 10.0499 453.0000 17 23 5.3852 37 50 10.4403 492.0000 17 28 7.0000 37 53 7.8102 442 2.2361 18 23 9.0000 38 43 8.2462 511 1.4142 18 26 7.0711 38 45 7.2111 513 2.2361 18 29 5.3852 38 54 6.4031 693.6056 18 31 4.4721 38 55 7.8102 610 5.0000 18 36 8.2462 39 40 11.7047 615 5.8310 19 22 8.2462 39 54 6.0828 641 4.2426 19 23 5.3852 39 55 9.8489 642 5.0990 19 29 3.6056 40 46 8.9443 715 5.8310 19 30 5.0990 40 55 8.2462 731 3.6056 19 31 5.0990 43 44 7.2801 742 4.2426 20 21 4.4721 43 45 6.3246 834 3.6056 20 25 9.2195 43 46 5.0000 910 6.0000 20 41 8.2462 43 55 5.0000 911 4.1231 20 49 6.0000 44 46 7.2111 942 3.6056 20 51 9.2195 44 47 6.0000 10 20 8.0623 21 37 7.2111 46 47 11.6619 194 Ann Oper Res (2014) 222:175–195 Table 6 (Continued) Initial node Final node Edge length Initial node Final node Edge length Initial node Final node Edge length 10 21 9.0000 21 51 7.2801 46 52 15.6205 10 37 7.8102 22 33 4.2426 46 55 6.0000 10 41 6.0828 23 26 12.2066 47 52 16.1245 10 47 8.6023 23 28 7.0711 47 53 5.6569 11 13 2.2361 23 29 4.4721 48 49 6.4031 11 47 3.6056 24 26 7.2111 49 51 12.0416 12 14 10.6301 24 35 4.4721 50 52 8.6023 12 17 6.3246 24 36 9.0554 50 53 5.8310 12 22 8.6023 25 36 5.6569 52 53 12.1655 12 28 7.8102 25 41 4.1231 54 55 10.1980 13 44 2.8284 25 48 4.4721 References Bello, L., Blanquero, R., & Carrizosa, E. (2011). On minimax-regret Huff location models. Computers and Operations Research,38, 90–97. Berman, O., Drezner, Z., & Krass, D. (2011). Big segment small segment global optimization algorithm on networks. Networks,58(1), 1–11. Bertsekas, D. P. (1998). Network optimization: continuous and discrete models. Belmont: Athenas Scientific. Blanquero, R., & Carrizosa, E. (2009a). Continuous location problems and big triangle small triangle: constructing better bounds. Journal of Global Optimization,45, 389–402. Blanquero, R., & Carrizosa, E. (2009b). On covering methods for d.c. optimization. Journal of Global Optimization,18, 265–274. Corberán, A., & Sanchis, J. M. (2007). A branch & cut algorithm for the windy general routing problem and special cases. Networks,49, 245–257. Dearing, P. M., & Shier, D. R. (1983). Optimal locations for a class of nonlinear single-facility location problems on a network. Operations Research,31(2), 292–303. Drezner, T., & Drezner, Z. (2004). Finding the optimal solution to the Huff competitive location model. Computational Management Science,1, 193–208. Drezner, Z., & Suzuki, A. (2004). The big triangle small triangle method for the solution of non-convex facility location problems. Operations Research,52, 128–135. Fernández, J., Pelegrín, B., Plastria, F., & Tóth, B. (2007). Solving a Huff-like competitive location and design model for profit maximization in the plane. European Journal of Operational Research,170, 1274–1287. Hansen, P., Peeters, D., Richard, D., & Thisse, J. F. (1985). The minisum and minimax location problems revisited. Operations Research,33(6), 1251–1265. Horst, R., & Thoai, N. V. (1999). Dc programming: overview. Journal of Optimization Theory and Applications,103, 1–43. Huff, D. L. (1964). Defining and estimating a trading area. Journal of Marketing,8, 28–34. Huff, D. L. (1966). A programmed solution for approximating an optimum retail location. Land Economics, 8(42), 293–303. Labbé, M., Peeters, D., & Thisse, J. F. (1995). Location on networks. In M. O. Ball et al. (Eds.), Handbooks in operations research and management science (Vol. 8, pp. 551–624). Amsterdam: Elsevier. Marianov, V., & Serra, D. (2002). Location-allocation of multiple-server service centers with constrained queues or waiting times. Annals of Operations Research,111, 35–50. Plastria, F. (1992). GBSSS: the generalized big square small square method for planar single-facility location. European Journal of Operational Research,62(2), 163–174. Plastria, F. (2001). Static competitive facility location: an overview of optimisation approaches. European Journal of Operational Research,129(3), 461–470. Reinelt, G. (1991). Tsplib—a traveling salesman problem library. ORSA Journal on Computing,3(4), 376– 384. Ann Oper Res (2014) 222:175–195 195 Serra, D., ReVelle, C., & Rosing, K. (1999). Surviving in a competitive spatial market: the threshold capture model. Journal of Regional Science,39(4), 637–652. Tuy, H. (1995). DC optimization: theory, methods and algorithms, handbook of global optimization.Dordrecht: Kluwer Academic. Tuy, H., Al-Khayyal, F., & Zhou, F. (1995). A d.c. optimization method for single facility location problems. Journal of Global Optimization,7, 209–227.