# Robust Model of Discrete Competitive Facility Location Problem with Partially Proportional Rule

Authors: Wuyang Yu
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

When consumers faced with the choice of competitive chain facilities that offer exclusive services, current rules cannot describe these customers’ behaviors very well. So we propose a partially proportional rule to represent this kind of customer behavior. In addition, the exact demands of customers in many real-world environments are often difficult to determine. This is contradicting to the assumption in most studies of the competitive facility location problem. For the competitive facility location problem with the partially proportional rule, we establish a robust optimization model to handle the uncertainty of customers’ demands. We propose two methods to solve the robust model by studying the properties of the counterpart problem. The first methodMIP is presented by solving a mixed-integer optimization model of the counterpart problem directly. The second method SAS is given by embedding a sorting subalgorithm into the simulated annealing framework, in which the sorting subalgorithm can easily solve the subproblem. The effects of the budget and the robust control parameter to the location scheme are analyzed in a quasi-real example. The result shows that changes in the robust control parameter can affect the customer demands that were captured by the new entrants, thereby changing the optimal solution for facility location. In addition, there is a threshold of the robust control parameter for any given budget. Only when the robust control parameter is larger than this threshold, the market share captured by the new entering firm increases with the increases of this parameter. Finally, numerical experiments show the superiority of the algorithm SAS in large-scare competitive facility location problems.


---

## Body

## 1. Introduction

Competitive facility location (CFL) is the problem of locating new facilities in competitive markets with the aim of maximizing market share [1, 2]. A lot of competitive facility location models have been proposed in recent years because these models are very useful in many competitive situations [3]. There are various competitive facility location models depending on the components to be considered. These components usually include competitive type, location space, and customer behavior (see survey paper [4]).Suppose that there are different facilities offering similar goods or services, customer behavior refers to the way that a customer how to spend his buying power on these facilities. It is obvious that customer behavior is a key ingredient in the competitive facility location problem [5]. The two most common customer behavior rules employed in literature are the binary rule and the proportional rule [6]. The binary rule dates back to the duopoly model proposed by Hotelling [7]. It is traditionally assumed that customers patronize the nearest facility to be served, but other characteristics of the facilities can also be taken into account [8–10]. In fact, if we use the parameter (attraction) to uniformly represent the different characteristics considered by competitive facilities, then the essence of the binary rule is that the customers always patronize the most attractive facility. The proportional rule is first proposed by Huff [11], which assumes that the customers patronize all facilities in proportion to their attractions [12–16]. Many researchers have studied different customer behavior rules, but most of them can be regarded as variants of the above two rules. Such as in [17], the authors studied the customer behavior with minimal attraction requirement. Drezner et al. [18] investigated the competitive facility location problem based on the concept of “radius of influence.” In their model, although the calculation of the attraction is based on the facility’s cover distance, the division of the market share still relies on the proportional rule. Qi et al. [19] observed that customers usually patronize facilities within a range that they feel is convenient. They proposed a proportional rule with a limited range to study the competitive facility location problem. Fernández et al. [20] studied the customer behavior that customer only patronizes those facilities that the attraction greater than or equal to a threshold value.According to the classification of Suáez-Vega et al. [21], there is another basic rule of customer behavior called the partially binary rule. Following this rule, the customer first selects the most attractive facility from each firm and then splits his demand among those facilities proportionally to their attractions [22]. The partially binary rule is usually studied together with the binary rule and/or the proportional rule. For example, Suárez-Vega et al. [21] and Biesinger et al. [23] studied the competitive facility location problems in the network space and the discrete space, respectively. Six scenarios have been considered in both papers, they are combinations of two service types (essential and unessential) and three customer behavior rules (binary, proportional, and partially binary). Fernández et al. [24] proposed two new heuristic algorithms to solve competitive facility location problems with the binary rule and the partially binary rule, respectively. Fernández et al. [25] studied a continuous competitive single-facility location and design problem with the partially binary rule. The comparison with the proportional rule reveals their location results are quite different.Most of the customer behaviors in competitive facility location problems can be described by these three rules. However, when consumers faced with the choice of competitive chain facilities that offer exclusive services, we find that these three rules cannot describe customers’ behaviors well. For example, when somebody wants to apply for a bank account for deposit/withdrew or other services, because each bank has several business service points (including manual service points and ATMs), a customer usually uses the total attraction of all business service points as the evaluation value of the bank. The typical customer behavior is as follows: the customer chooses the most attractive bank to apply a bank account, and then he will patronize the business service points of this bank in proportion to their attractions. From the typical behavior of how customers choose competitive chain facilities that provide exclusive services, we can find that neither of the binary rule, the proportional rule, nor the partially binary rule adequately describes this kind of customer behavior. In this paper, we propose a new kind of customer behavior rule to describe this kind of customer behavior, and we call this rule as the partially proportional rule.In most studies of the competitive facility location problem, it is usually assumed that the demands of customers are known. However, due to different reasons such as lack of historical data of customers or the high-cost to conduct a comprehensive market survey, the exact demands of customers in many real-world environments are often difficult to determine. To the best of our knowledge, only two literatures involved the uncertainty of the demands of customers. Shiode and Drezner [26] assumed that the weights of demand points are stochastic, and then they studied the competitive facility location problem in the Stackelberg game framework. Beresnev and Melnikov [27] studied the competitive facility location problem in the situation of there are several alternative demand scenarios. The uncertainty of demands in the above two papers is actually expressed by random variables, but it is quite difficult to determine the corresponding random distribution in many cases. Therefore, we use interval data to express the uncertainty of demands and then adopt a robust optimization method to deal with this uncertainty.The remainder of this paper is organized as follows: Section2 is devoted to present the definition of the partially proportional rule and the robust optimization method. Section 3 consists the robust model of competitive facility location problem with the partially proportional rule. Section 4 includes two solution methods MIP and SAS to solve the presented model. The effects of the budget and the robust control parameter in the location scheme are presented in Section 5, and the superiority of SAS is also shown in large-scale numerical tests. Finally, some conclusions are presented in Section 6.


## 2. Problem Definition

### 2.1. Partially Proportional Rule

In this section, we first introduce a suitable rule to describe the behavior of how customers choose competitive chain facilities that provide exclusive services. Suppose there are several firms in the competitive market, each firm has several chain facilities that offer exclusive services. A customer only chooses one of these services due to the exclusive of these services. If the customer chooses one of these firms, any of the firm’s chain facilities can provide the corresponding exclusive service to the customer. In order to describe well this kind of customer behavior, thepartially proportional rule is defined as follows: a customer chooses the most attractive firm from all firms at first and then splits his buy powers on the selected firm’s facilities in proportion to their attractions.In order to illustrate the difference between the binary rule, the proportional rule, the partially binary rule, and the partially proportional rule, we give an example to show how customer will patronize facilities under different behavior rules. Suppose that there are three firms A, B, and C in the competing market. Each firm has several chain facilities, firm A has facilities 1, 4, 5, and 9, firm B has facilities 2, 6, and 7, and firm C has facilities 3, 8, and 10. For simplicity, suppose that for a customer, the attractionak of the facility k is equal to its serial number, i.e., ak=k. The patronizing patterns under different customer behavior rules are shown in Figure 1.Figure 1
Patronizing patterns under different customer behavior rules.If the customer follows the binary rule, he only patronizes the facility 10 because this facility has the largest attraction among all facilities. If the customer follows the proportional rule, he patronizes all the facilities and splits his demands to facilityk with the proportion k/55, where 55 is the total attraction of all facilities. If the customer follows the partially binary rule, he selects the most attractive facility from each firm at first, that is, facility 9 of firm A, facility 7 of firm B, and facility 10 of firm C. Then he patronizes facilities 9, 7, and 10 and splits his demands in proportion to 9/26, 7/26, and 10/26, respectively. If the customer follows the partially proportional rule, he compares the total attraction of each firm at first. In this example, the total attraction of firms A, B, and C is 19, 15, and 21, respectively. The most attractive firm is C; hence, he patronizes facilities 3, 8, and 10 of firm C and splits his demands in proportion to 3/21, 8/21, and 10/21, respectively.Note that for a customer, if there are multiple firms that are most attractive, then the customer patronizes all the facilities of these firms and splits his demands in proportion to their attractions.


### 2.2. Robust Optimization

A robust model with interval data uncertainty was first proposed by Soyster [28], but his method usually leads to a solution that is too conservative because of always considering the worst case. Ben-Tal and Nemirovsky [29] addressed the over-conservatism by applying an ellipsoid uncertainty set, but the robust counterpart model cannot conserve the complexity of the original problem. Bertsims and Sim [30] presented a method with the flexibility to adjust the degree of conservativeness while maintaining the linear robust counterparts. Before establishing a robust model, we first briefly describe the robust model introduced in [30]. The authors considered the following problem:(1)maxc′x,s.t.Ax≤b,l≤x≤u,where x is the n-dimensional decision variable, and xj is the j-th component of x. Some parameters of the coefficient matrix A are uncertain. For a particular row i of A, let Ji denote the set of coefficients aij,j∈Ji that are uncertain. Suppose each entry aij takes values in a¯ij−a^ij,a¯ij+a^ij. For every i, a parameter Γi∈0,Ji is introduced as the robust control parameter. The goal is to be protected against all cases that up to Γi of these coefficients are allowed to change, and one coefficient ait changes by Γi−Γia^it. The role of the parameter Γi is to adjust the robustness against the level of conservatism of the solution. Then the robust model of the original problem is(2)maxc′x,s.t.∑j∈Ja¯ijxj+βx,Γi≤bi,∀i,−yj≤xj≤yj,∀j,l≤x≤u,y≥0,where(3)βx,Γi=maxSi∪tiSi⊆Ji,Si=Γ,ti∈Ji\Si∑j∈Sja^ijyj+Γi−Γia^itiyti.


## 2.1. Partially Proportional Rule

In this section, we first introduce a suitable rule to describe the behavior of how customers choose competitive chain facilities that provide exclusive services. Suppose there are several firms in the competitive market, each firm has several chain facilities that offer exclusive services. A customer only chooses one of these services due to the exclusive of these services. If the customer chooses one of these firms, any of the firm’s chain facilities can provide the corresponding exclusive service to the customer. In order to describe well this kind of customer behavior, thepartially proportional rule is defined as follows: a customer chooses the most attractive firm from all firms at first and then splits his buy powers on the selected firm’s facilities in proportion to their attractions.In order to illustrate the difference between the binary rule, the proportional rule, the partially binary rule, and the partially proportional rule, we give an example to show how customer will patronize facilities under different behavior rules. Suppose that there are three firms A, B, and C in the competing market. Each firm has several chain facilities, firm A has facilities 1, 4, 5, and 9, firm B has facilities 2, 6, and 7, and firm C has facilities 3, 8, and 10. For simplicity, suppose that for a customer, the attractionak of the facility k is equal to its serial number, i.e., ak=k. The patronizing patterns under different customer behavior rules are shown in Figure 1.Figure 1
Patronizing patterns under different customer behavior rules.If the customer follows the binary rule, he only patronizes the facility 10 because this facility has the largest attraction among all facilities. If the customer follows the proportional rule, he patronizes all the facilities and splits his demands to facilityk with the proportion k/55, where 55 is the total attraction of all facilities. If the customer follows the partially binary rule, he selects the most attractive facility from each firm at first, that is, facility 9 of firm A, facility 7 of firm B, and facility 10 of firm C. Then he patronizes facilities 9, 7, and 10 and splits his demands in proportion to 9/26, 7/26, and 10/26, respectively. If the customer follows the partially proportional rule, he compares the total attraction of each firm at first. In this example, the total attraction of firms A, B, and C is 19, 15, and 21, respectively. The most attractive firm is C; hence, he patronizes facilities 3, 8, and 10 of firm C and splits his demands in proportion to 3/21, 8/21, and 10/21, respectively.Note that for a customer, if there are multiple firms that are most attractive, then the customer patronizes all the facilities of these firms and splits his demands in proportion to their attractions.


## 2.2. Robust Optimization

A robust model with interval data uncertainty was first proposed by Soyster [28], but his method usually leads to a solution that is too conservative because of always considering the worst case. Ben-Tal and Nemirovsky [29] addressed the over-conservatism by applying an ellipsoid uncertainty set, but the robust counterpart model cannot conserve the complexity of the original problem. Bertsims and Sim [30] presented a method with the flexibility to adjust the degree of conservativeness while maintaining the linear robust counterparts. Before establishing a robust model, we first briefly describe the robust model introduced in [30]. The authors considered the following problem:(1)maxc′x,s.t.Ax≤b,l≤x≤u,where x is the n-dimensional decision variable, and xj is the j-th component of x. Some parameters of the coefficient matrix A are uncertain. For a particular row i of A, let Ji denote the set of coefficients aij,j∈Ji that are uncertain. Suppose each entry aij takes values in a¯ij−a^ij,a¯ij+a^ij. For every i, a parameter Γi∈0,Ji is introduced as the robust control parameter. The goal is to be protected against all cases that up to Γi of these coefficients are allowed to change, and one coefficient ait changes by Γi−Γia^it. The role of the parameter Γi is to adjust the robustness against the level of conservatism of the solution. Then the robust model of the original problem is(2)maxc′x,s.t.∑j∈Ja¯ijxj+βx,Γi≤bi,∀i,−yj≤xj≤yj,∀j,l≤x≤u,y≥0,where(3)βx,Γi=maxSi∪tiSi⊆Ji,Si=Γ,ti∈Ji\Si∑j∈Sja^ijyj+Γi−Γia^itiyti.


## 3. Proposed Model

Suppose there are several firms in a region, each providing exclusive services through their chain facilities. A new entering firm wants to set up several chain facilities to offer its exclusive service. The problem of the new entering firm is how to select the locations of its chain facilities to gain the largest market share under a given budget. There is a set of discrete customers that can aggregate to demand points. The demand of each customer is uncertain, but the range of its value is known. The following notations are used:K: set of existing firms, indexed by kI: set of demand points (customers), indexed by iJk: set of firm k’s chain facilities, k is an existing firmJ: set of existing facilities, i.e., J=∪k∈KJkJn: set of potential facility locations, indexed by jwi: demand of customer i∈I, it is uncertain within the interval ci,ci+difj: opening cost at potential facility location j∈Jndij: distance between customer i and facility jaij: attraction that customer i feels from facility jG: budget of the new entering firmyj,j∈Jn: binary variable, yj=1 if facility j is select to open and 0 otherwiseThe goal of the new entrant is to maximize the market share. Hence, the competitive facility location problem for the new entering firm is as follows:(4)z=max∑i∈Iwia˜iy,s.t.∑j∈Jnfjyj≤G,yj∈0,1,j∈Jn.where a˜iy indicates the proportion of demand that the new entering firm captured from customer i. ∀i∈I, let Ri=maxk∈K∑j∈Jkaij, and let Ki=k′∈K∑j∈Jk′aij=Ri. If customers follow the partially proportional rule, then a˜iy can be expressed as(5)a˜iy=1,∑j∈Jnaijyj>Ri,1Ki+1,∑j∈Jnaijyj=Ri,0,∑j∈Jnaijyj<Ri.By adding binary variablesxi+,xi−, and let ϵ be a positive real number that is small enough, condition (5) can be expressed with some linear constraints. Hence, the above competitive facility location model can be transformed into an integer optimization model as follows:(6)z=max∑i∈IwiKixi++xi−Ki+1,(7)s.t.∑j∈Jnfjyj≤G,xi+Ri+ε≤∑j∈Jnaijyj,i∈I,xi−Ri−ε≤∑j∈Jnaijyj,i∈I,xi+,xi−,yj∈0,1,i∈I,j∈Jn.Note that for anyyj, if Ri+ε≤∑j∈Jnaijyj, then the optimal solution is xi+=1,xi−=1. In this case, the new entering firm gets the total demand of customer i, i.e., a˜iy=1. If Ri−ε≤∑j∈Jnaijyj<Ri+ε, then the optimal solution is xi+=0,xi−=1, and the proportion of demand that the new entering firm get from the customer i is a˜iy=1/Ki+1. If ∑j∈Jnaijyj<Ri−ε, then the optimal solution is xi+=0,xi−=0, that is, a˜iy=0. The role of ϵ is to compare the size of Ri and ∑j∈Jnaijyj; hence, any positive number that is less than mini∈Iminy∈yj∑j∈Jnyj=sRi−∑j∈Jnaijyj can be used as ϵ.LetDx,y denote the solution set of (7), and DΓ=μ∑i∈Iμi≥Γ,0≤μi≤1,i∈I, where 0≤Γ≤I is a parameter which is used to control the total uncertainty of wi. Naturally, the uncertainty of wi can be expressed as wi=ci+μidi with μ∈DΓ. Similar to the robust model that is proposed in [30], the robust model of the above optimization problem can be expressed as follows:(8)RPZ=maxx,y∈Dx,y∑i∈IciKixi++xi−Ki+1+minμ∈DΓ∑i∈IμidiKixi++xi−Ki+1.


## 4. Solution Method

In this section, we propose two methods to solve the robust optimization problemRP. The first method MIP is based on directly solving an integer model that is equivalent to the origin robust model RP.


### 4.1.MIP Method

Givenxi+,xi−,i∈I, let us consider the suboptimization problem SP:(9)min∑i∈IμidiKixi++xi−Ki+1,SPs.t.∑i∈Iμi≥Γ,0≤μi≤1,i∈I.The dual problem ofSP is as follows:(10)maxΓθ−∑i∈Ivi,s.t.θ−vi≤diKixi++xi−Ki+1,i∈I,θ≥0,vi≥0,i∈I.So the robust optimization modelRP can be converted to RP1 as follows:(11)maxΓθ−∑i∈Ivi+∑i∈IciKixi++xi−Ki+1,RP1s.t.∑j∈Jnfjyj≤Gxi+Ri+ε≤∑j∈Jnaijyj,i∈I,xi−Ri−ε≤∑j∈Jnaijyj,i∈I,θ−vi≤diKixi++xi−Ki+1,i∈I,θ≥0,vi≥0,xi+,xi−,yj∈0,1,i∈I,j∈Jn.Since the equivalence modelRP1 is a mixed integer optimization model, the first method is to solve it directly using some commerce software such as CPLEX and GUROBI. We denote this method as Method MIP.


### 4.2.SAS Method

The second method to solveRP is based on the simulated annealing framework. The key point is that we found the structure of the optimal solution to the subproblem. Now we first present a Lemma 1.Lemma 1.
For0<a1≤a2≤⋯≤an,0≤g≤n, let Dg=x∑i=1nxi≥g,0≤xi≤1,i=1,2,…,n. Consider the optimization problem OP: minx∈Dgz=∑i=1naixi. Let g¯ be the largest integer less than or equal to g. If x∗ is defined as follows, it is the optimal solution of OP:(12)xi∗=1,i≤g¯,g−g¯,i=g¯+1,0,i>g¯+1.Proof.
Firstly, we prove that∑i=1nxi∗=g is a necessary condition for the optimal solution x∗ of problem OP. Assuming this conclusion is wrong, i.e., there is ∑i=1nxi∗−g=ε>0. Let xi′∗=minxi∗>0,i=1,2,…,n. We can construct a new solution x˜ as follows: set x˜i≔xi∗,∀i≠i′, and x˜i′≔xi′∗−minxi′∗,ε. It is easy to know that ∑i=1nx˜i=∑i=1nxi∗−minxi′∗,ε≥g, but zx∗−zx˜=ai′minxi′∗,ε>0. This is a contradiction.
Secondly, we prove that there exists an optimal solution, where at most one componentxi takes the value 0<xi<1, and the other components are 0 or 1. Assuming this conclusion is wrong, i.e., for every optimal solution, there are at least two components xi1 and xi2 such that 0<xi1,xi2<1. Without loss of generality, we can assume that i1<i2. Now we can construct a new solution. Set x˜i≔xi,∀i≠i1,i2, and x˜i1,x˜i2 as follows:(1)
Ifxi1+xi2≤1, then x˜i1≔xi1+xi2,x˜i2≔0(2)
Ifxi1+xi2>1, then x˜i1≔1,x˜i2≔xi1+xi2−1
Obviously,x˜∈Dg, for case (1) we know that zx−zx˜=ai1xi1+ai2xi2−ai1xi1+xi2=ai2−ai1xi2≥0. For case (2), we know that zx−zx˜=ai1xi1+ai2xi2−ai1+ai2xi1+xi2−1=1−xi1ai2−ai1≥0. So x˜ is also an optimal solution. This is a contradiction.
Lastly, we prove that the solution defined in (11) is an optimal solution of OP. Therefore, we know that the optimal solution of (OP) must satisfy ∑i=1nxi=g, and at most one component xi takes the value 0<xi<1. So we can conclude that there are g¯ components xi equal to 1, and one component xi equals g−g¯ and all other components are 0. Therefore, according to Chebyshev’s sum inequality [31], we know that x∗ defined by (11) is the optimal solution of the problem OP.
By applying the above Lemma1 to the problem SP, we can get the following Theorem 1.Theorem 1.
For givenxi+,xi−, let I1=i∈Ixi+=1,xi−=1,I2=i∈Ixi+=0,xi−=1,I3=i∈Ixi+=0,xi−=0. Let d¯i=di,i∈I1,d¯i=di/Ki+1,i∈I2. Then μ∗ is an optimal solution of the problem SP if it is defined as follows:(1)
IfΓ≤I3, μi∗=1,i∈I3, and μi∗=0,i∈I2∪I3(2)
IfΓ>I3, record the sequence of d¯ii∈I1∪I2 in the ascending order as i1,i2,…,is, where s=I1∪I2, then μj∗=1,j∈I3, and for j∈I1∪I2, μj∗ defined as follows:(13)μj∗=1,j<ik′+1,Γ−Γ¯,j=ik′+1,0,j>ik′+1,where k′=Γ¯−I3.Proof.
The problemSP can be rewritten as follows:(14)SP1zx=minμ∈DΓ∑i∈I1μidi+∑i∈I2μidiKi+1.(1)
IfΓ≤I3, from μi∗=1,i∈I3, and μi∗=0,i∈I1∪I2, we can get that μ∗∈DΓ. Obviously, the objective function corresponding to μ∗ is 0. Hence, we can conclude that μ∗ is an optimal solution of SP in this case.(2)
IfΓ>I3, firstly, we prove that the optimal solution of SP1 satisfies μi=1,∀i∈I3. Assuming this conclusion is wrong, i.e., there is a i′∈I3 such that μi′<1. Let σ=1−μi′>0. It is easy to construct a new solution as follows:
Setμ¯i≔1,∀i∈I3. Find any one i˜∈I1∪I2 such that μi˜>0, let μ¯i˜≔μi˜−minσ,μi˜ and μ¯i≔μi for other i≠i˜. We can get that zμ>zμ¯. This is a contradiction. Hence, there must be μi=1, ∀i∈I3. We can conclude that μi∗ is an optimal solution of SP1 according to Lemma 1.
Given a facility location solutiony, we propose a sorting-based algorithm to solve SP as follows, this algorithm can be used to calculate the market share that was captured by the new entering firm:SortA:(1)
Defined¯i,i∈I, as follows:(15)d¯i=di,ifRi<∑j∈Jnaijyj,diKi+1,ifRi=∑j∈Jnaijyj,0,ifRi>∑j∈Jnaijyj.And definec¯i,i∈I, in a similar way (just substitute di by ci).(2)
Record the sequence ofd¯ in the ascending order as d˜. Take into account the uncertainty of wi, the market share of the new entering firm is(16)Zy=∑i∈Ic¯i+∑i≤Γ¯d˜i+Γ−Γ¯d˜Γ¯+1.
The simulated annealing algorithm was introduced by Kirkpatrick et al. [32]. By embedding the sorting-based algorithm (SortA) into the simulated annealing framework, we present the solution method in Algorithm 1.
The process of CreateInitialSolution (ProblemSize) is as follows:
Calculatesj=∑i∈IAij/fj,∀j∈Jn, by sorting sj in the descend order. Record the rearranged index sequence as j1,j2,…,jN. Find the maximal index jt such that fj1+fj2+⋯+fjt≤G, then let yj′=1 for j=j1,j2,…,jt and yj′=0 for j=jt+1,…,jN to generate an initial feasible solution y′∈D, where D=y∑j∈Jnfjyj≤G,yj∈0,1,j∈Jn.
The process of CreateNeighborSolution (Ycurrent) is as follows:Lety′≔Ycurrent, calculate the remaining budget after opening new facilities Rbudget≔G−∑j∈Jnfjyj′. Let I1y′≔j∈Jnyj′=1,I0y′≔j∈Jnyj′=0. Randomly choose two subsets SI1⊂I1y′ and SI0⊂I0y′ such that ∑j∈SI0fj−∑j∈SI1fj≤Rbudget, and then let yj″≔0,j∈SI1, yj″≔1,j∈SI0, and yj″≔yj′ for other j.Algorithm 1: SAS method.
Input: ProblemSize,T0, Tend, r, L
Output:Ybest
Y
current ⟵ CreateInitialSolution (ProblemSize)
Obj
Y
current ⟵ algorithm SortA
Y
best ⟵ Ycurrent
While (T0≥Tend)
For (i=1 To L)
Ytemp ⟵ CreateNeighborSolution (Ycurrent)
ObjYtemp ⟵ algorithm SortA
If (ObjYtemp ≥ ObjYcurrent)
Ycurrent ⟵ Ytemp
If (ObjYtemp ≥ ObjYbest)
Ybest ⟵ Ytemp
End
Else If (Exp(ObjYtemp−ObjYcurrent/T0) > Rand())
Ycurrent ⟵ Ytemp
End
End
T0 ⟵ r⋅T0
End
Return (Ybest)


## 4.1.MIP Method

Givenxi+,xi−,i∈I, let us consider the suboptimization problem SP:(9)min∑i∈IμidiKixi++xi−Ki+1,SPs.t.∑i∈Iμi≥Γ,0≤μi≤1,i∈I.The dual problem ofSP is as follows:(10)maxΓθ−∑i∈Ivi,s.t.θ−vi≤diKixi++xi−Ki+1,i∈I,θ≥0,vi≥0,i∈I.So the robust optimization modelRP can be converted to RP1 as follows:(11)maxΓθ−∑i∈Ivi+∑i∈IciKixi++xi−Ki+1,RP1s.t.∑j∈Jnfjyj≤Gxi+Ri+ε≤∑j∈Jnaijyj,i∈I,xi−Ri−ε≤∑j∈Jnaijyj,i∈I,θ−vi≤diKixi++xi−Ki+1,i∈I,θ≥0,vi≥0,xi+,xi−,yj∈0,1,i∈I,j∈Jn.Since the equivalence modelRP1 is a mixed integer optimization model, the first method is to solve it directly using some commerce software such as CPLEX and GUROBI. We denote this method as Method MIP.


## 4.2.SAS Method

The second method to solveRP is based on the simulated annealing framework. The key point is that we found the structure of the optimal solution to the subproblem. Now we first present a Lemma 1.Lemma 1.
For0<a1≤a2≤⋯≤an,0≤g≤n, let Dg=x∑i=1nxi≥g,0≤xi≤1,i=1,2,…,n. Consider the optimization problem OP: minx∈Dgz=∑i=1naixi. Let g¯ be the largest integer less than or equal to g. If x∗ is defined as follows, it is the optimal solution of OP:(12)xi∗=1,i≤g¯,g−g¯,i=g¯+1,0,i>g¯+1.Proof.
Firstly, we prove that∑i=1nxi∗=g is a necessary condition for the optimal solution x∗ of problem OP. Assuming this conclusion is wrong, i.e., there is ∑i=1nxi∗−g=ε>0. Let xi′∗=minxi∗>0,i=1,2,…,n. We can construct a new solution x˜ as follows: set x˜i≔xi∗,∀i≠i′, and x˜i′≔xi′∗−minxi′∗,ε. It is easy to know that ∑i=1nx˜i=∑i=1nxi∗−minxi′∗,ε≥g, but zx∗−zx˜=ai′minxi′∗,ε>0. This is a contradiction.
Secondly, we prove that there exists an optimal solution, where at most one componentxi takes the value 0<xi<1, and the other components are 0 or 1. Assuming this conclusion is wrong, i.e., for every optimal solution, there are at least two components xi1 and xi2 such that 0<xi1,xi2<1. Without loss of generality, we can assume that i1<i2. Now we can construct a new solution. Set x˜i≔xi,∀i≠i1,i2, and x˜i1,x˜i2 as follows:(1)
Ifxi1+xi2≤1, then x˜i1≔xi1+xi2,x˜i2≔0(2)
Ifxi1+xi2>1, then x˜i1≔1,x˜i2≔xi1+xi2−1
Obviously,x˜∈Dg, for case (1) we know that zx−zx˜=ai1xi1+ai2xi2−ai1xi1+xi2=ai2−ai1xi2≥0. For case (2), we know that zx−zx˜=ai1xi1+ai2xi2−ai1+ai2xi1+xi2−1=1−xi1ai2−ai1≥0. So x˜ is also an optimal solution. This is a contradiction.
Lastly, we prove that the solution defined in (11) is an optimal solution of OP. Therefore, we know that the optimal solution of (OP) must satisfy ∑i=1nxi=g, and at most one component xi takes the value 0<xi<1. So we can conclude that there are g¯ components xi equal to 1, and one component xi equals g−g¯ and all other components are 0. Therefore, according to Chebyshev’s sum inequality [31], we know that x∗ defined by (11) is the optimal solution of the problem OP.
By applying the above Lemma1 to the problem SP, we can get the following Theorem 1.Theorem 1.
For givenxi+,xi−, let I1=i∈Ixi+=1,xi−=1,I2=i∈Ixi+=0,xi−=1,I3=i∈Ixi+=0,xi−=0. Let d¯i=di,i∈I1,d¯i=di/Ki+1,i∈I2. Then μ∗ is an optimal solution of the problem SP if it is defined as follows:(1)
IfΓ≤I3, μi∗=1,i∈I3, and μi∗=0,i∈I2∪I3(2)
IfΓ>I3, record the sequence of d¯ii∈I1∪I2 in the ascending order as i1,i2,…,is, where s=I1∪I2, then μj∗=1,j∈I3, and for j∈I1∪I2, μj∗ defined as follows:(13)μj∗=1,j<ik′+1,Γ−Γ¯,j=ik′+1,0,j>ik′+1,where k′=Γ¯−I3.Proof.
The problemSP can be rewritten as follows:(14)SP1zx=minμ∈DΓ∑i∈I1μidi+∑i∈I2μidiKi+1.(1)
IfΓ≤I3, from μi∗=1,i∈I3, and μi∗=0,i∈I1∪I2, we can get that μ∗∈DΓ. Obviously, the objective function corresponding to μ∗ is 0. Hence, we can conclude that μ∗ is an optimal solution of SP in this case.(2)
IfΓ>I3, firstly, we prove that the optimal solution of SP1 satisfies μi=1,∀i∈I3. Assuming this conclusion is wrong, i.e., there is a i′∈I3 such that μi′<1. Let σ=1−μi′>0. It is easy to construct a new solution as follows:
Setμ¯i≔1,∀i∈I3. Find any one i˜∈I1∪I2 such that μi˜>0, let μ¯i˜≔μi˜−minσ,μi˜ and μ¯i≔μi for other i≠i˜. We can get that zμ>zμ¯. This is a contradiction. Hence, there must be μi=1, ∀i∈I3. We can conclude that μi∗ is an optimal solution of SP1 according to Lemma 1.
Given a facility location solutiony, we propose a sorting-based algorithm to solve SP as follows, this algorithm can be used to calculate the market share that was captured by the new entering firm:SortA:(1)
Defined¯i,i∈I, as follows:(15)d¯i=di,ifRi<∑j∈Jnaijyj,diKi+1,ifRi=∑j∈Jnaijyj,0,ifRi>∑j∈Jnaijyj.And definec¯i,i∈I, in a similar way (just substitute di by ci).(2)
Record the sequence ofd¯ in the ascending order as d˜. Take into account the uncertainty of wi, the market share of the new entering firm is(16)Zy=∑i∈Ic¯i+∑i≤Γ¯d˜i+Γ−Γ¯d˜Γ¯+1.
The simulated annealing algorithm was introduced by Kirkpatrick et al. [32]. By embedding the sorting-based algorithm (SortA) into the simulated annealing framework, we present the solution method in Algorithm 1.
The process of CreateInitialSolution (ProblemSize) is as follows:
Calculatesj=∑i∈IAij/fj,∀j∈Jn, by sorting sj in the descend order. Record the rearranged index sequence as j1,j2,…,jN. Find the maximal index jt such that fj1+fj2+⋯+fjt≤G, then let yj′=1 for j=j1,j2,…,jt and yj′=0 for j=jt+1,…,jN to generate an initial feasible solution y′∈D, where D=y∑j∈Jnfjyj≤G,yj∈0,1,j∈Jn.
The process of CreateNeighborSolution (Ycurrent) is as follows:Lety′≔Ycurrent, calculate the remaining budget after opening new facilities Rbudget≔G−∑j∈Jnfjyj′. Let I1y′≔j∈Jnyj′=1,I0y′≔j∈Jnyj′=0. Randomly choose two subsets SI1⊂I1y′ and SI0⊂I0y′ such that ∑j∈SI0fj−∑j∈SI1fj≤Rbudget, and then let yj″≔0,j∈SI1, yj″≔1,j∈SI0, and yj″≔yj′ for other j.Algorithm 1: SAS method.
Input: ProblemSize,T0, Tend, r, L
Output:Ybest
Y
current ⟵ CreateInitialSolution (ProblemSize)
Obj
Y
current ⟵ algorithm SortA
Y
best ⟵ Ycurrent
While (T0≥Tend)
For (i=1 To L)
Ytemp ⟵ CreateNeighborSolution (Ycurrent)
ObjYtemp ⟵ algorithm SortA
If (ObjYtemp ≥ ObjYcurrent)
Ycurrent ⟵ Ytemp
If (ObjYtemp ≥ ObjYbest)
Ybest ⟵ Ytemp
End
Else If (Exp(ObjYtemp−ObjYcurrent/T0) > Rand())
Ycurrent ⟵ Ytemp
End
End
T0 ⟵ r⋅T0
End
Return (Ybest)


## 5. Numerical Example

### 5.1. A Quasi-Real Example

Consider the 49-node data set described in Daskin [33] consisting of the capitals of the continental United States plus Washington, DC. The lower and upper bounds of the customer demands (ci,c¯i) are proportional to the population of the state (i) in 1890 and 1990, respectively. Hence di≔c¯i−ci. Now, there are two firms (A and B) that offer similar exclusive services through their chain facilities. Firm A has 5 facilities located in the states of Alabama, Arizona, Iowa, Maine, and Nebraska, respectively. Firm B has 5 facilities located in the states of New Hampshire, Ohio, South Dakota, Utah, and Wyoming, respectively. The attraction aij that the customer i feels from facility j is calculated as 1/dij+1,∀i∈I,j∈J∪Jn. The opening cost at the potential facility locations j∈Jn is also taken from the data set of Daskin. Before the new firm C entered the competitive market, the market was completely divided into two parts by firm A and firm B. The division pattern of the market share is shown in Figure 2.Figure 2
Division pattern of the market (blue circles: customers; red triangles: firm A’s facilities; black rectangles: firm B’s facilities; plum red regions: customers captured by firm A; white regions: customers captured by firm B).Now firm C wants to enter this market by opening several chain facilities. Any customer point without an existing facility can be regarded as a potential location for new facilities. The goal of firm C is to maximize the market share captured by itself under the constraint of a limited budgetG. The robust control parameter Γ can take any value in the interval 0,49. For simplicity, Γ takes integer values in our experiments. If the value of Γ is 0 or 49, then wi is equal to the lower limit ci or the upper limit ci+di, respectively. This means that the demand of each customer has been determined. Hence, we limit the parameter Γ to the integers in the interval 1,48.ForG=100,000, when the parameter Γ≤42, the optimal decision of firm C is to open two facilities in Pennsylvania and Mississippi, respectively. The opening cost of these two facilities is 93000. The location of facilities and the customers that captured by each firm are shown in Figure 3.Figure 3
Division pattern of the market withG=100000 and Γ≤42 (blue circles: customers; red triangles: firm A’s facilities; yellow circles: firm C’s facilities; black rectangles: firm B’s facilities; plum red regions: customers captured by firm A; white regions: customers captured by firm B; orange regions: customers captured by firm C).However, whenΓ≥43, the optimal locations of firm C’s facilities are Pennsylvania and Michigan, respectively. The opening cost of these two facilities is 86,800. Similarly, all information is shown in Figure 4.Figure 4
Division pattern of the market withG=100,000 and Γ≥43 (blue circles: customers; red triangles: firm A’s facilities; yellow circles: firm C’s facilities; black rectangles: firm B’s facilities; plum red regions: customers captured by firm A; white regions: customers captured by firm B; orange regions: customers captured by firm C).For simplicity, we refer to location plans corresponding toΓ=42 and Γ=43 as Plan 1 and Plan 2, respectively. Why does firm C’s location plan change when Γ changes from 42 to 43? In order to find the cause of this change, we list the customer demand ranges captured by firm C in Table 1. The overall contents of Table 1 are divided into two parts: Plan 1 and Plan 2. Take the part of Plan 1 as an example, the first column is the customers captured by firm C, the second column is the demand intervals of customers. The column μ42∗ denotes the optimal solution of μ in SP when Γ=42, and the definition of μ43∗ is similar.Table 1
Comparison of optimal solutions in Plan 1 and Plan 2.
Customer
Demand interval
μ
42
∗
μ
43
∗
Plan 1
Delaware
[1.68493, 6.66168]
1
1
Illinois
[38.26352, 114.30602]
0
0
Louisiana
[11.18588, 42.19973]
1
1
Maryland
[10.42390, 47.81468]
0
1
Mississippi
[12.89600, 25.73216]
1
1
New Jersey
[14.44933, 77.30188]
0
0
North Carolina
[16.17949, 66.28637]
0
0
Pennsylvania
[52.58113, 118.81643]
0
0
Virginia
[16.55980, 61.87358]
0
0
Wisconsin
[16.93330, 48.91769]
0
0
Washington DC
[2.30392, 6.06900]
1
1
Plan 2
Delaware
[1.68493, 6.66168]
1
1
Florida
[3.91422, 129.37926]
0
0
Georgia
[18.37353, 64.78216]
0
0
Kentucky
[18.58635, 36.85296]
1
1
Maryland
[10.42390, 47.81468]
1
1
Michigan
[20.93890, 92.95297]
0
0
New Jersey
[14.44933, 77.30188]
0
0
North Carolina
[16.17949, 66.28637]
0
0
Pennsylvania
[52.58113, 118.81643]
0
0
Virginia
[16.55980, 61.87358]
0
1
Washington DC
[2.30392, 6.06900]
1
1Given a facility location plany, denote the set of customers captured by firm C as Cuy. Then the minimal demands captured by firm C can be expressed as ∑i∈Cuyci+μi∗di. According to the Lemma 1 proposed in Section 4.2, the values of μ42∗ and μ43∗ are easy to obtain. When Γ=42, the demands captured by firm C with Plan 1 and Plan 2 are 246.0530 and 240.3947, respectively. So in this case, Plan 1 is better than Plan 2. However, when Γ=43, the demands captured by firm C with Plan 1 and Plan 2 are 283.4438 and 285.7085, respectively. Obviously, Plan 2 is better than Plan 1 in this case. Therefore, we can get a conclusion that changes in the robust control parameter can affect the customer demands that were captured by the new entrants, thereby changing the optimal solution for facility location.ForG=100,000,150,000,200,000,250,000,300,000, the relationships between the demands captured by firm C and the values of Γ are shown in Figure 5. From Figure 5, it is obvious that a larger budget makes the entering firm capture more demands. This is in line with our intuition, as more budgets expand the feasible range and the value of the objective function increases accordingly. For a given budget, there is a threshold for the parameter Γ. Only when the parameter Γ is larger than this threshold, the market share that was captured by the new entering firm increases with the increase of Γ. In addition, the threshold of the robust control parameter Γ decreases as the budget increases. Therefore, uncertainty will have a greater impact on the location of competitive facilities when the budget is larger.Figure 5
Curves between firm C’s market share andΓ (G=100,000,150,000,200,000,250,000,300,000).LetG=100000. For each Γ∈5,10,15,20,25,30,35,40,45, generate scenarios ξii=1,2,…,10 as follows: randomly select Γ number of customers, let their demands to be the corresponding upper bounds, and the demands of the remaining customers to be the corresponding lower bounds. For each ξi, we can solve the deterministic model of the competitive facility location problem to get the optimal location scheme for firm C. Denote the demands captured by firm C in the worst case as ZDξi. Let ZD be the average of ZDξi and record the customer demands obtained by firm C based on the robust model as ZR. For each value of Γ, the comparison of ZD and ZR is listed in Table 2.Table 2
Comparison of the robust model and deterministic model (G=100,000).
Γ
5
10
15
20
25
30
35
40
45
Z
D
167.4710
161.7401
154.7905
156.5370
145.0612
160.5284
153.0439
158.2890
352.7972
Z
R
193.4612
193.4612
193.4612
193.4612
193.4612
193.4612
193.4612
202.2030
382.2240From Table2, ZD<ZR for all values of Γ. In fact, for the uncertain demands expressed by interval data, the solution obtained by the robust model is better than the solution obtained by other models in the worst case.The optimal solutions of firm C’s competitive facility location problem are presented in Table3. The first column denotes the budget G of the firm C, the second column denotes the range of the parameter Γ, the third column denotes the facilities’ opening cost of the firm C, the fourth column denotes the optimal facility locations, and the last column denotes the customers captured by the firm C.Table 3
Locations for different budget G and parameterΓ.
G
Γ
Cost
Locations
Captured customers
100,000
[1, 42]
93,000
Pennsylvania
Delaware, Illinois, Louisiana, Maryland, Mississippi, New Jersey
Mississippi
North Carolina, Pennsylvania, Virginia, Wisconsin, Washington DC
[43, 48]
86,800
Pennsylvania
Delaware, Florida, Georgia, Kentucky, Maryland, Michigan
Michigan
New Jersey, North Carolina, Pennsylvania, Virginia, Washington DC
150,000
[1, 48]
141,400
Pennsylvania
Delaware, Florida, Georgia, Illinois, Kentucky, Louisiana
Mississippi
Maryland, Michigan, Mississippi, New Jersey, North Carolina
Michigan
Pennsylvania, Virginia, Wisconsin, Washington DC
200,000
[1, 33]
200,000
Indiana
Delaware, Illinois, Indiana, Kentucky, Louisiana
Louisiana
Maryland, Michigan, Mississippi, Missouri, New Jersey
New Jersey
Pennsylvania, Virginia, Wisconsin, Washington DC
[34, 40]
197,100
Kansas
Arkansas, Delaware, Florida, Georgia, Illinois
Michigan
Kansas, Kentucky, Louisiana, Maryland, Michigan
Missouri
Mississippi, Missouri, New Jersey, North Carolina
Pennsylvania
Pennsylvania, Virginia, Wisconsin, Washington DC
[41, 48]
174,300
Michigan
Delaware, Florida, Georgia, Illinois, Kentucky, Louisiana, Maryland
Missouri
Michigan, Mississippi, New Jersey, North Carolina, New York
New Jersey
North Carolina, Pennsylvania, Virginia, Wisconsin, Washington DC
250,000
[1, 42]
242,000
Indiana
Arkansas, Delaware, Florida, Georgia, Illinois, Indiana
Michigan
Kentucky, Louisiana, Maryland, Michigan, Mississippi, Missouri
Missouri
New Jersey, New York, North Carolina, Pennsylvania
New Jersey
South Carolina, Tennessee, Virginia, Wisconsin, Washington DC
[43, 48]
247,900
Louisiana
California, Delaware, Florida, Georgia, Illinois, Kentucky
Michigan
Louisiana, Maryland, Michigan, Mississippi, Missouri, Nevada
New Jersey
New Jersey, New York, North Carolina, Oregon, Pennsylvania
Oregon
South Carolina, Virginia, Washington, Wisconsin, Washington DC
300,000
[1, 37]
299,600
Arkansas
Arkansas, Delaware, Florida, Georgia, Illinois
Indiana
Indiana, Kentucky, Louisiana, Maryland, Michigan
Michigan
Mississippi, Missouri, New Jersey, New York, North Carolina
New Jersey
Oklahoma, Pennsylvania, South Carolina, Tennessee
Oklahoma
Texas, Virginia, Wisconsin, Washington DC
[38, 45]
292,700
Idaho
Arkansas, California, Delaware, Florida, Georgia, Idaho
Indiana
Illinois, Indiana, Kentucky, Louisiana, Maryland, Michigan
Minnesota
Minnesota, Mississippi, Missouri, Nevada, New Jersey
Oklahoma
North Carolina, Oklahoma, Oregon, Pennsylvania, Texas
Pennsylvania
Virginia, Washington, Wisconsin, Washington DC
[46, 48]
296,400
Michigan
Arkansas, California, Delaware, Florida, Georgia
Missouri
Illinois, Kentucky, Louisiana, Maryland, Michigan
New Jersey
Mississippi, Missouri, Nevada, New Jersey, New York
Oklahoma
North Carolina, Oklahoma, Oregon, Pennsylvania, Texas
Oregon
Virginia, Washington, Washington DC


### 5.2. Performances of Two Methods

To evaluate the performances of the two methods proposed in this paper, we generated the following 12 experimental examples. There are already 3 firms in the competitive market, each with 5 chain facilities.M is the number of customers, N is the number of potential facility locations, and G is the budget of the new entering firm. The attraction aij are generated randomly within 1,100. The lower bound of demand ci and the ranges of uncertain demand di are generated randomly within 10,100. The opening costs fj are generated randomly within 500,1000. The robust control parameter Γ takes value as N/2. The method MIP is solved by CPLEX 12.5 to get the exact solutions, and the method SAS is coded on the platform of Matlab R2015b. Parameters in the method SAS are T0=1000,Tend=1,r=0.99,andL=200. The settings of these parameters are determined based on the results of multiple calculations, and they ensure that the algorithm SAS can find an approximate optimal solution in a short time. Due to the randomness of the algorithm SAS, each example is calculated 20 times. The results include the objective values, CPU times for both methods, and the gaps between objective values and CPU times are presented in Table 4. Here the gaps of two indicators between two methods are defined as follows: Gap (Obj, Time) = [(Obj, Time)(MIP) − (Ave, Time)(SAS)]/(Ave, Time)(SAS).Table 4
Comparison of the methodMIP and the method SAS.
(M, N, G)
MIP
SAS
GAP
Time (s)
Obj
Time (s)
Max
Ave
Min
Time
Obj
(100, 30, 2000)
4.6613
255
8.6972
255
255
255
−0.46
0
(100, 30, 3000)
147.9407
938
5.6810
938
938
938
25.04
0
(100, 30, 4000)
486.3300
4435
6.0592
4435
4435
4435
79.26
0
(100, 50, 2000)
8.9027
306
7.0733
306
306
306
0.26
0
(100, 50, 3000)
1437.0
2331
8.0774
2331
2327.7
2265
176.90
0.0014
(100, 50, 4000)
4080.4
4618
8.2981
4618
4618
4618
490.73
0
(200, 30, 2000)
21.4404
413
7.2421
413
413
413
1.96
0
(200, 30, 3000)
330.1236
255
7.2915
255
255
255
44.28
0
(200, 30, 4000)
1422.3
7807
8.3080
7807
7807
7807
170.20
0
(200, 50, 2000)
147.4972
508
8.8188
508
508
508
15.73
0
(200, 50, 3000)
1843.5
3941
10.2389
3926
3926
3926
179.05
0.0038
(200, 50, 4000)
44326
8514
9.3528
8493
8492.7
8490
4738.4
0.0025From Table4, it is obvious that only for small-scale problems, method MIP is better than method SAS both in objective value and CPU time. When the scale of problems increases, the CPU time required for the method MIP increases dramatically, while the CPU time required for the method SAS remains stable. At the same time, the objective values obtained by the method SAS are only slightly worse than that obtained by the method MIP. Therefore, the method SAS is a good choice to solve large-scale competitive facility location problems.


## 5.1. A Quasi-Real Example

Consider the 49-node data set described in Daskin [33] consisting of the capitals of the continental United States plus Washington, DC. The lower and upper bounds of the customer demands (ci,c¯i) are proportional to the population of the state (i) in 1890 and 1990, respectively. Hence di≔c¯i−ci. Now, there are two firms (A and B) that offer similar exclusive services through their chain facilities. Firm A has 5 facilities located in the states of Alabama, Arizona, Iowa, Maine, and Nebraska, respectively. Firm B has 5 facilities located in the states of New Hampshire, Ohio, South Dakota, Utah, and Wyoming, respectively. The attraction aij that the customer i feels from facility j is calculated as 1/dij+1,∀i∈I,j∈J∪Jn. The opening cost at the potential facility locations j∈Jn is also taken from the data set of Daskin. Before the new firm C entered the competitive market, the market was completely divided into two parts by firm A and firm B. The division pattern of the market share is shown in Figure 2.Figure 2
Division pattern of the market (blue circles: customers; red triangles: firm A’s facilities; black rectangles: firm B’s facilities; plum red regions: customers captured by firm A; white regions: customers captured by firm B).Now firm C wants to enter this market by opening several chain facilities. Any customer point without an existing facility can be regarded as a potential location for new facilities. The goal of firm C is to maximize the market share captured by itself under the constraint of a limited budgetG. The robust control parameter Γ can take any value in the interval 0,49. For simplicity, Γ takes integer values in our experiments. If the value of Γ is 0 or 49, then wi is equal to the lower limit ci or the upper limit ci+di, respectively. This means that the demand of each customer has been determined. Hence, we limit the parameter Γ to the integers in the interval 1,48.ForG=100,000, when the parameter Γ≤42, the optimal decision of firm C is to open two facilities in Pennsylvania and Mississippi, respectively. The opening cost of these two facilities is 93000. The location of facilities and the customers that captured by each firm are shown in Figure 3.Figure 3
Division pattern of the market withG=100000 and Γ≤42 (blue circles: customers; red triangles: firm A’s facilities; yellow circles: firm C’s facilities; black rectangles: firm B’s facilities; plum red regions: customers captured by firm A; white regions: customers captured by firm B; orange regions: customers captured by firm C).However, whenΓ≥43, the optimal locations of firm C’s facilities are Pennsylvania and Michigan, respectively. The opening cost of these two facilities is 86,800. Similarly, all information is shown in Figure 4.Figure 4
Division pattern of the market withG=100,000 and Γ≥43 (blue circles: customers; red triangles: firm A’s facilities; yellow circles: firm C’s facilities; black rectangles: firm B’s facilities; plum red regions: customers captured by firm A; white regions: customers captured by firm B; orange regions: customers captured by firm C).For simplicity, we refer to location plans corresponding toΓ=42 and Γ=43 as Plan 1 and Plan 2, respectively. Why does firm C’s location plan change when Γ changes from 42 to 43? In order to find the cause of this change, we list the customer demand ranges captured by firm C in Table 1. The overall contents of Table 1 are divided into two parts: Plan 1 and Plan 2. Take the part of Plan 1 as an example, the first column is the customers captured by firm C, the second column is the demand intervals of customers. The column μ42∗ denotes the optimal solution of μ in SP when Γ=42, and the definition of μ43∗ is similar.Table 1
Comparison of optimal solutions in Plan 1 and Plan 2.
Customer
Demand interval
μ
42
∗
μ
43
∗
Plan 1
Delaware
[1.68493, 6.66168]
1
1
Illinois
[38.26352, 114.30602]
0
0
Louisiana
[11.18588, 42.19973]
1
1
Maryland
[10.42390, 47.81468]
0
1
Mississippi
[12.89600, 25.73216]
1
1
New Jersey
[14.44933, 77.30188]
0
0
North Carolina
[16.17949, 66.28637]
0
0
Pennsylvania
[52.58113, 118.81643]
0
0
Virginia
[16.55980, 61.87358]
0
0
Wisconsin
[16.93330, 48.91769]
0
0
Washington DC
[2.30392, 6.06900]
1
1
Plan 2
Delaware
[1.68493, 6.66168]
1
1
Florida
[3.91422, 129.37926]
0
0
Georgia
[18.37353, 64.78216]
0
0
Kentucky
[18.58635, 36.85296]
1
1
Maryland
[10.42390, 47.81468]
1
1
Michigan
[20.93890, 92.95297]
0
0
New Jersey
[14.44933, 77.30188]
0
0
North Carolina
[16.17949, 66.28637]
0
0
Pennsylvania
[52.58113, 118.81643]
0
0
Virginia
[16.55980, 61.87358]
0
1
Washington DC
[2.30392, 6.06900]
1
1Given a facility location plany, denote the set of customers captured by firm C as Cuy. Then the minimal demands captured by firm C can be expressed as ∑i∈Cuyci+μi∗di. According to the Lemma 1 proposed in Section 4.2, the values of μ42∗ and μ43∗ are easy to obtain. When Γ=42, the demands captured by firm C with Plan 1 and Plan 2 are 246.0530 and 240.3947, respectively. So in this case, Plan 1 is better than Plan 2. However, when Γ=43, the demands captured by firm C with Plan 1 and Plan 2 are 283.4438 and 285.7085, respectively. Obviously, Plan 2 is better than Plan 1 in this case. Therefore, we can get a conclusion that changes in the robust control parameter can affect the customer demands that were captured by the new entrants, thereby changing the optimal solution for facility location.ForG=100,000,150,000,200,000,250,000,300,000, the relationships between the demands captured by firm C and the values of Γ are shown in Figure 5. From Figure 5, it is obvious that a larger budget makes the entering firm capture more demands. This is in line with our intuition, as more budgets expand the feasible range and the value of the objective function increases accordingly. For a given budget, there is a threshold for the parameter Γ. Only when the parameter Γ is larger than this threshold, the market share that was captured by the new entering firm increases with the increase of Γ. In addition, the threshold of the robust control parameter Γ decreases as the budget increases. Therefore, uncertainty will have a greater impact on the location of competitive facilities when the budget is larger.Figure 5
Curves between firm C’s market share andΓ (G=100,000,150,000,200,000,250,000,300,000).LetG=100000. For each Γ∈5,10,15,20,25,30,35,40,45, generate scenarios ξii=1,2,…,10 as follows: randomly select Γ number of customers, let their demands to be the corresponding upper bounds, and the demands of the remaining customers to be the corresponding lower bounds. For each ξi, we can solve the deterministic model of the competitive facility location problem to get the optimal location scheme for firm C. Denote the demands captured by firm C in the worst case as ZDξi. Let ZD be the average of ZDξi and record the customer demands obtained by firm C based on the robust model as ZR. For each value of Γ, the comparison of ZD and ZR is listed in Table 2.Table 2
Comparison of the robust model and deterministic model (G=100,000).
Γ
5
10
15
20
25
30
35
40
45
Z
D
167.4710
161.7401
154.7905
156.5370
145.0612
160.5284
153.0439
158.2890
352.7972
Z
R
193.4612
193.4612
193.4612
193.4612
193.4612
193.4612
193.4612
202.2030
382.2240From Table2, ZD<ZR for all values of Γ. In fact, for the uncertain demands expressed by interval data, the solution obtained by the robust model is better than the solution obtained by other models in the worst case.The optimal solutions of firm C’s competitive facility location problem are presented in Table3. The first column denotes the budget G of the firm C, the second column denotes the range of the parameter Γ, the third column denotes the facilities’ opening cost of the firm C, the fourth column denotes the optimal facility locations, and the last column denotes the customers captured by the firm C.Table 3
Locations for different budget G and parameterΓ.
G
Γ
Cost
Locations
Captured customers
100,000
[1, 42]
93,000
Pennsylvania
Delaware, Illinois, Louisiana, Maryland, Mississippi, New Jersey
Mississippi
North Carolina, Pennsylvania, Virginia, Wisconsin, Washington DC
[43, 48]
86,800
Pennsylvania
Delaware, Florida, Georgia, Kentucky, Maryland, Michigan
Michigan
New Jersey, North Carolina, Pennsylvania, Virginia, Washington DC
150,000
[1, 48]
141,400
Pennsylvania
Delaware, Florida, Georgia, Illinois, Kentucky, Louisiana
Mississippi
Maryland, Michigan, Mississippi, New Jersey, North Carolina
Michigan
Pennsylvania, Virginia, Wisconsin, Washington DC
200,000
[1, 33]
200,000
Indiana
Delaware, Illinois, Indiana, Kentucky, Louisiana
Louisiana
Maryland, Michigan, Mississippi, Missouri, New Jersey
New Jersey
Pennsylvania, Virginia, Wisconsin, Washington DC
[34, 40]
197,100
Kansas
Arkansas, Delaware, Florida, Georgia, Illinois
Michigan
Kansas, Kentucky, Louisiana, Maryland, Michigan
Missouri
Mississippi, Missouri, New Jersey, North Carolina
Pennsylvania
Pennsylvania, Virginia, Wisconsin, Washington DC
[41, 48]
174,300
Michigan
Delaware, Florida, Georgia, Illinois, Kentucky, Louisiana, Maryland
Missouri
Michigan, Mississippi, New Jersey, North Carolina, New York
New Jersey
North Carolina, Pennsylvania, Virginia, Wisconsin, Washington DC
250,000
[1, 42]
242,000
Indiana
Arkansas, Delaware, Florida, Georgia, Illinois, Indiana
Michigan
Kentucky, Louisiana, Maryland, Michigan, Mississippi, Missouri
Missouri
New Jersey, New York, North Carolina, Pennsylvania
New Jersey
South Carolina, Tennessee, Virginia, Wisconsin, Washington DC
[43, 48]
247,900
Louisiana
California, Delaware, Florida, Georgia, Illinois, Kentucky
Michigan
Louisiana, Maryland, Michigan, Mississippi, Missouri, Nevada
New Jersey
New Jersey, New York, North Carolina, Oregon, Pennsylvania
Oregon
South Carolina, Virginia, Washington, Wisconsin, Washington DC
300,000
[1, 37]
299,600
Arkansas
Arkansas, Delaware, Florida, Georgia, Illinois
Indiana
Indiana, Kentucky, Louisiana, Maryland, Michigan
Michigan
Mississippi, Missouri, New Jersey, New York, North Carolina
New Jersey
Oklahoma, Pennsylvania, South Carolina, Tennessee
Oklahoma
Texas, Virginia, Wisconsin, Washington DC
[38, 45]
292,700
Idaho
Arkansas, California, Delaware, Florida, Georgia, Idaho
Indiana
Illinois, Indiana, Kentucky, Louisiana, Maryland, Michigan
Minnesota
Minnesota, Mississippi, Missouri, Nevada, New Jersey
Oklahoma
North Carolina, Oklahoma, Oregon, Pennsylvania, Texas
Pennsylvania
Virginia, Washington, Wisconsin, Washington DC
[46, 48]
296,400
Michigan
Arkansas, California, Delaware, Florida, Georgia
Missouri
Illinois, Kentucky, Louisiana, Maryland, Michigan
New Jersey
Mississippi, Missouri, Nevada, New Jersey, New York
Oklahoma
North Carolina, Oklahoma, Oregon, Pennsylvania, Texas
Oregon
Virginia, Washington, Washington DC


## 5.2. Performances of Two Methods

To evaluate the performances of the two methods proposed in this paper, we generated the following 12 experimental examples. There are already 3 firms in the competitive market, each with 5 chain facilities.M is the number of customers, N is the number of potential facility locations, and G is the budget of the new entering firm. The attraction aij are generated randomly within 1,100. The lower bound of demand ci and the ranges of uncertain demand di are generated randomly within 10,100. The opening costs fj are generated randomly within 500,1000. The robust control parameter Γ takes value as N/2. The method MIP is solved by CPLEX 12.5 to get the exact solutions, and the method SAS is coded on the platform of Matlab R2015b. Parameters in the method SAS are T0=1000,Tend=1,r=0.99,andL=200. The settings of these parameters are determined based on the results of multiple calculations, and they ensure that the algorithm SAS can find an approximate optimal solution in a short time. Due to the randomness of the algorithm SAS, each example is calculated 20 times. The results include the objective values, CPU times for both methods, and the gaps between objective values and CPU times are presented in Table 4. Here the gaps of two indicators between two methods are defined as follows: Gap (Obj, Time) = [(Obj, Time)(MIP) − (Ave, Time)(SAS)]/(Ave, Time)(SAS).Table 4
Comparison of the methodMIP and the method SAS.
(M, N, G)
MIP
SAS
GAP
Time (s)
Obj
Time (s)
Max
Ave
Min
Time
Obj
(100, 30, 2000)
4.6613
255
8.6972
255
255
255
−0.46
0
(100, 30, 3000)
147.9407
938
5.6810
938
938
938
25.04
0
(100, 30, 4000)
486.3300
4435
6.0592
4435
4435
4435
79.26
0
(100, 50, 2000)
8.9027
306
7.0733
306
306
306
0.26
0
(100, 50, 3000)
1437.0
2331
8.0774
2331
2327.7
2265
176.90
0.0014
(100, 50, 4000)
4080.4
4618
8.2981
4618
4618
4618
490.73
0
(200, 30, 2000)
21.4404
413
7.2421
413
413
413
1.96
0
(200, 30, 3000)
330.1236
255
7.2915
255
255
255
44.28
0
(200, 30, 4000)
1422.3
7807
8.3080
7807
7807
7807
170.20
0
(200, 50, 2000)
147.4972
508
8.8188
508
508
508
15.73
0
(200, 50, 3000)
1843.5
3941
10.2389
3926
3926
3926
179.05
0.0038
(200, 50, 4000)
44326
8514
9.3528
8493
8492.7
8490
4738.4
0.0025From Table4, it is obvious that only for small-scale problems, method MIP is better than method SAS both in objective value and CPU time. When the scale of problems increases, the CPU time required for the method MIP increases dramatically, while the CPU time required for the method SAS remains stable. At the same time, the objective values obtained by the method SAS are only slightly worse than that obtained by the method MIP. Therefore, the method SAS is a good choice to solve large-scale competitive facility location problems.


## 6. Conclusion

When consumers faced with the choice of competitive chain facilities that offer exclusive services, all current rules cannot describe these customers’ behaviors very well. In this paper, the partially proportional rule is proposed to describe this kind of customer behavior. After expressing the uncertainty of customer demand with interval data, we use robust optimization method to study the competitive facility location problem with uncertain demand. Two methods are proposed to solve the presented robust optimization model. One is an exact methodMIP, and it is based on directly solving an integer model equivalent to the original robust model. The other is a heuristic method SAS, which is a combination of the simulated annealing framework and a sorting-based algorithm. A quasi-real example is analyzed in detail. One result that is found through this example is the budget of the new entering firm is an important factor for capturing more market share. Another result is for any given budget, there is a threshold of the robust control parameter Γ. If the parameter Γ is larger than this threshold, the market share captured by the new entering firm will increase as this parameter increases. Random example tests show that for small scale problems, the method MIP is a good choice. But for large-scale problems, the method SAS is better than MIP because the former only requires very little computation time, while maintaining a fairly good calculation result.For customer behavior of choosing competitive facilities that provide exclusive services, the partially proportional rule is more suitable to describe it than other rules. The proposed model can be used to solve the competitive facility location problem with different customer behavior rules, as long as there is uncertainty in customer demand. We are looking forward to studying the competitive facility problem with more uncertain parameters in the future.



---

*Source: 3107431-2020-01-21.xml*