Optimal Location of Fast Charging Stations for Mixed Tra ﬃ c of Electric Vehicles and Gasoline Vehicles Subject to Elastic Demands

: With the rapid development of electric vehicles (EVs), one of the urgent issues is how to deploy limited charging facilities to provide services for as many EVs as possible. This paper proposes a bilevel model to depict the interaction between tra ﬃ c ﬂow distribution and the location of charging stations (CSs) in the EVs and gasoline vehicles (GVs) hybrid network. The upper level model is a maximum ﬂow-covering model where the CSs are deployed on links with higher demands. The lower level model is a stochastic user equilibrium model under elastic demands (SUE-ED) that considers both demands uncertainty and perceived path constraints, which have a signiﬁcant inﬂuence on the distribution of link ﬂow. Besides the path travel cost, the utility of charging facilities, charging speed, and waiting time at CSs due to space capacity restraint are also considered for the EVs when making a path assignment in the lower level model. A mixed-integer nonlinear program is constructed, and the equivalence of SUE-ED is proven, where a heuristic algorithm is used to solve the model. Finally, the network trial and sensitivity analysis are carried out to illustrate the feasibility and e ﬀ ectiveness of the proposed model.


Introduction
The market share of new energy vehicles without fossil fuels has been increasing rapidly in recent years, especially electric vehicles (EVs), which provide better performance, higher efficiency, and zero emissions [1,2]. To promote the usage of EVs, governments around the world have successively issued a series of active policies to provide subsidies for EV purchases and to deploy public charging infrastructures at convenient locations [2,3]. However, the market acceptance of EVs seems to fluctuate with the decreasing subsidy. On the one hand, limited driving range, insufficient public charging infrastructures, and longer charging time are the main reasons for at-home charging. On the other hand, the deployment of EV charging stations may induce more trips [4,5], that is, demands between each origin to destination (OD) pair that respond to the available fast charging are elastic, so does the route choice. Therefore, it is especially important to study how to effectively and reasonably deploy EV public charging stations within the hybrid network of rapidly increasing EVs and dominated gasoline vehicles (GVs), which should help reduce the range anxiety of EV users and maximize the coverage rate of EVs.
According to previous studies, a variety of factors affect the location of charging stations, including, among others, preference and travel choice behavior of users [6,7], travel demand of users [8], information of the en-route energy consumption of EVs [9][10][11][12], information of the remained

Literature Review
The conventional and dominant three types of location optimization models include the point demand model [19,20], the flow demand model [20][21][22][23], and the multi-objective optimization model [24][25][26][27]. The P-Median model assumes that "charging demand is generated in the road network node" and is widely used as one of the point demand models. The flow demand model is formulated on the basis of the Flow Capturing Location Model (FCLM). For the first time, some researches proposed the Flow Refueling Location Model (FRLM), where mileage limitations were explicitly considered in facility location issues [28][29][30][31]. One branch of FRLM aims to maximize the demand coverage by locating a fixed number of charging facilities, which is called the maximum coverage location problem. Although the multi-objective optimization models have the advantage of addressing more complicated experimental requirements, they are not good at dealing with the uncertainty planning problem [32].
In a hybrid network with both EVs and GVs, the layout of EV charging stations and the distribution of EV flow affect mutually. Elastic demands in a hybrid network result in many more uncertainties. Most of the previous studies used the User Equilibrium Model (UE) to locate charging facilities. In these studies, Xu et al. dealt with the user equilibrium problem in a hybrid transport network with battery switching stations and road grade constraints [33]. Jing et al. gave a comprehensive discussion of the equilibrium network model [34]. Jiang et al. introduced the path distance constraints into the UE model [13]. Zheng et al. proposed a bilevel model where the upper layer minimized travel costs, and the lower layer aimed to find the path-constrained EVs equalization flow [35]. The bounded rationality of EV users led to much more complicated energy consumption [7,36,37], and, therefore, route choice behavior reflected more unobserved heterogeneity, which resulted in various elastic demands.
The network design problem with elastic demands, where the induced or transferred OD demand is the subject of responses in traveler itinerary choices to enjoy the improvement of new infrastructure, have several formulations with various motivations. Ge et al. [38] is one of the early attempts to consider both the proportion of EVs and the charge rate of EVs when determining the elastic charging demands from the total number of vehicles on road connections. The elastic demand was formulated based either on feedback of congested travel and congested station on route choices [39] or on the assumption that charge demand between OD pairs follows a nonlinear inverse cost function without considering the pre-generating paths and charging combinations [40].
It should be noted that it is hard to consider all these constraints simultaneously, say, the elastic demand of the road network, the capacity of the charging stations, and the range limit of EVs by using the SUE-ED model in location problems of public EV charging stations. In this study, a novel bilevel public charging station location model that combines a flow-capturing location model and a multi-objective optimization model is proposed.

Problem Assumptions and EV Paths Analysis
The proposed location model has two optimization objectives. The upper layer employs an improved maximum coverage model to maximize the coverage of the total EV flow by deploying a given number of charging stations on the links where the EV flow is the largest. The charging stations can serve more EV users by deploying charging facilities on the links where most EV drivers pass, which is an effective way to improve the utilization rate of public charging facilities and alleviate range anxieties of electric vehicle drivers [4]. The lower layer uses an SUE-ED model with incomplete information. The assignment result of SUE is a decisive factor for the placement of the charging stations. There is an inherent difference in the driving behavior of each user (EV users and GV users), in the hybrid network. In particular, range limit, charging time, and location of charging facilities have an important impact on the path choose behaviors of EV users.

Notations
i: The set of traveler types in the network. i presents an element, mainly including EVs and GVs. When specific instructions are needed, subscripts g and e are used to indicate variables or parameters related to GVs and EVs, respectively; W: The set of OD pairs. w is an element of the set; R w : The set of paths between all OD pairs; p: The number of charging stations subject to the budget; U: The utility value of the charging stations, resulting in a reduction in travel cost; K·U: The waiting cost in charging stations, where K is a constant coefficient; l w k : Length of path k between OD pair w; l ri k : Length of path k between node r and I; l ij k : Length of path k between node i and j; l js k : Length of path k between node j and s; R e : The range limit of EVs; ε: Charging power. Unit: hour/mile; t w ck : The minimum charging time; v a : Traffic flow on link a, including EVs and GVs, so v a = v ae + v ag ; t a (v ag , v ae ): Actual travel time through link a; t 0 a : Free travel time through link a; H a : Capacity of link a; x a : Binary parameter. If the charging facility is on link a, x a = 1, otherwise x a = 0; δ w ak : Binary parameter. If the path k crosses the link a, δ w ak = 1, otherwise δ w ak = 0; q w i : Elastic travel demand of type-i traveler between OD pair w; q w i : Maximum potential demand of type-i traveler between OD pair w; f w ki : Traffic flow of type-i traveler on path k between OD pair w; p w ki : The probability that type-i traveler chooses the path k between OD pair w; c w k : The actual travel time of the path k between OD pair w, and c w k = t a ( v a )δ w ak ; c w ki : The generalized travel cost of type-i traveler on path k between OD pair w; C wi : The expected perceived travel cost of the type-i traveler between OD pair w; D wi (·): The traffic demand function of the type-i traveler between OD pair w; c w ( → x ): Vector of the actual travel cost of all paths between OD pair w; θ i : a non-negative parameter that characterizes the uncertainty of type-i traveler's understanding of the path travel time.

Propose Assumptions
This paper focus on a hybrid network, denoted as G = (N,A), where N is the set of nodes, A is the set of links. W is used to represent a set of OD pairs, and w is an element, w ∈ W. R w represents the set of paths between w. The main differences between GVs and EVs are the limitation of travel distance and the composition of travel costs. Without loss of generality, the following assumptions are made: (a) There are two types of users in the mixed network: GVs and EVs, where EVs have an identical range limit; (b) The travel demand of GVs and EVs between each OD pair is elastic. And two types of vehicles have incomplete information about the travel cost; (c) Each EV is fully charged at its origin; (d) The level of anxiety and risk-taking behaviors of EV drivers are not considered in this model; (e) The charging time of electric vehicles is linearly related to mileage; (f) The charging facilities will be placed in the middle of links whose EV flow ranks p; (g) Deploying the charging facility on the road/path will increase the attractiveness of the route, which is called the utility of the charging facilities U, reducing the path travel cost; (h) Charging facilities have a fixed charging capacity. It is assumed that the waiting cost is proportional to the attraction value, expressed by K·U.

Analysis of EV Paths
When performing SUE-ED traffic assignment, the GVs can choose any route because of no range limit. However, the selection of feasible path sets is required before the distribution of EV flow. To describe the three types of EV paths more clearly, the concept of sub-path is described as follows.
Suppose that the origin is r and the destination is s in a pair OD w, and there are two charging stations i, j, which are also regarded as nodes, located in the arc of the path k ( Figure 1). When k ri , k ij , k js no longer includes other charging facilities except i, j, they are called sub-paths [41]. f) The charging facilities will be placed in the middle of links whose EV flow ranks p; 155 g) Deploying the charging facility on the road/path will increase the attractiveness of the route, 156 which is called the utility of the charging facilities U, reducing the path travel cost; 157 h) Charging facilities have a fixed charging capacity. It is assumed that the waiting cost is 158 proportional to the attraction value, expressed by · . 159

Analysis of EV Paths 160
When performing SUE-ED traffic assignment, the GVs can choose any route because of no 161 range limit. However, the selection of feasible path sets is required before the distribution of EV 162 flow. To describe the three types of EV paths more clearly, the concept of sub-path is described as 163 follows. 164 Suppose that the origin is r and the destination is s in a pair OD w, and there are two charging 165 stations i, j, which are also regarded as nodes, located in the arc of the path k ( Figure 1). When 166 , , no longer includes other charging facilities except i, j, they are called sub-paths [41]. 167 168 Figure 1. Description of the sub-path.

169
Based on the relationship between the travel mileage limit and the sub-path distance, the 170 three travel paths of EVs are discussed, as shown in  Based on the relationship between the travel mileage limit R e and the sub-path distance, the three travel paths of EVs are discussed, as shown in Table 1: Scenario 1. A travel path that can be completed without charging. When l w k ≤ R e , the EV drivers can reach the destination without charging. The generalized travel cost is c w ke = c w k (no charging stations) or c w ke = c w k − U (with charging stations). Scenario 2. A travel path that cannot be completed after charging. When l ri k ≥ R e or l ij k ≥ R e or l js k ≥ R e , the path k is an infeasible path. And the generalized path travel cost becomes infinite (c w ke = ∞) and the probability that the drivers select the path is zero. Scenario 3. A travel path that can be completed after charging. The path k is defined as a feasible path, if l w k ≥ R e and l ri k ≤ R e and l ij k ≤ R e and l js k ≤ R e . The drivers can complete the trip by charging at least once and the minimum charging time is t w ck = ε . . . (l w k − R e ). The generalized travel cost is expressed as c w ke = c w k + t w ck − U + K·U = c w k + t w ck + (K − 1)U, consisting of four parts: travel time, charging time, charging facilities utility, and charging stations' waiting cost.
For GVs, the general path travel cost is c w kg = c w k .

2.
A path that cannot be completed after charging l ri k ≥ R e or l i j k R e or l js k R e c w ke = ∞ 3. A path that can be completed after charging l w k ≥ R e and l ri k ≤ R e and l i j k ≤ R e and l js k ≤ R e c w ke = c w k + t w ck + (K − 1)U

Establish a Double-Layer Model
In this section, a bilevel optimization model for the location of charging facilities is proposed. The upper level model determines the location of the charging facilities by selecting the top p links. The lower level model calculates users' generalized travel cost, and randomly allocates the flow demand between the OD pairs to the filtered paths through known charging facilities location.

Upper Level Problem
The upper model is designed to maximize the total covered EV link flow by deploying a given number of charging facilities. The EV flow is covered when the charging facility is on the link. That is Equation (1) is the objective function of the upper model, and Equation (2) is the budget constraint indicating the number of charging stations p in a given network.

Lower Level Problem
The Bureau of Public Road (BPR) function in the lower model is employed [42].
Since U, t w ck , k are parameters that are independent of the equilibrium flow of the link, the Jacobian matrix [43] in the lower layer problem is as follows: It proves that the Jacobian matrix is symmetrical, and the lower layer model can be established as a convex function problem. It is assumed that all types of travelers make path selection in a random manner. According to the random utility theory, the probability that type-i traveler chooses the path k between OD pair w is: Assume that elastic travel demand q w i of type-i traveler is a strictly monotonically decreasing function of the expected minimum travel time between OD pair w, and with an upper bound.
According to the discrete selection theory [44,45], define C wi as: Energies 2020, 13,1964 6 of 16 For the type-i traveler, its SUE condition can be expressed as: In the mixed network, SUE-ED problem can be described by the equivalent mathematical programming model: Subject to : Equation (10) is the flow conservation constraint. Equation (11) is the non-negative constraint of path flow for type-i travelers. Equation (12) is the non-negative constraint of type-i traveler's OD flow demand. Equation (13) is the correlation between link flow and path flow. The novelty of this problem is that the introduction of sub-paths in the Equation (14) can generate a feasible set of paths in advance from the finite paths between each OD pair. The superscripts i and j include the origin r and the destination s of all OD pairs in Equation (14).
It is necessary to prove the equivalence between the solution of the proposed program Equation (9) and the solution of the SUE-ED model. The generalized Lagrangian function [46] of the mathematical model is constructed as follows: According to the Kuhn Tucker conditions [47], Equation (15) must satisfy the following conditions at the extreme point: Energies 2020, 13, 1964 7 of 16 The partial derivative of f w ki for Equation (15) is derived as follows: Note that if f w ki = 0, then ∂L ∂ f w ki does not exist. The above formula is only valid when f w ki > 0. So µ w ki = 0.
Equation (22) shows that type-i traveler follows the Logit model to select the travel path, which satisfies the SUE-ED condition described in Equation (8).
Calculate the partial derivative of q w i for Equation (15) as follows: Because ln q w i should exist, so ν wi = 0. Then derive from Equation (21): Comparing Equation (23) and Equation (24), the inverse function of the flow demand function is as follows: Then q w i = D wi (C wi ). This shows that Equations (9)- (14) can be used to represent a multi-user SUE problem under elastic demand. The BPR function t a (v a ) is a strictly monotonically increasing function of the link flow v a . The objective function is a strict convex function about the link flow vector v and the path flow vector f. At the same time, the constraints of Equations (9)-(14) are linear equality constraints and non-negative constraints, so its solution space is a convex set. According to the optimization theory, the strict convex function defined on the convex set only has one optimal solution.

Solution Method
The proposed bilevel model is an NP-hard (non-deterministic polynomial hard) problem because the lower level model is a mixed-integer nonlinear program with nonconvex path choosing, subject to uncertainties in elastic demands, where it is difficult to find a polynomial time complexity algorithm [48]. Thus, a heuristic method is employed to solve the problem of SUE-ED and MFC in an iterative manner. The detailed processes in Figure 2 are as follows:

249
Step 1: Input the road network parameters, such as the length of the link and the capacity of the 250 links, and find all the paths between the OD pair w, and record them as the initial path set. 251 Step 2: Initialize. There is no charging facility in the road network and relax the range limit of 252 EVs. Based on the travel time of the zero-flow link (0) = , calculate the effective path impedance 253 , and obtain the initial link flows (1)and (1). Let the upper generation counter z = 1. 254 Step 3: Sort all the links flow of EVs in ascending order to find the top p ranked links, and arrange 255 the charging facilities in the middle of the p links. Then increase the upper iteration counter by 1. 256 Step 4: Perform SUE-ED traffic distribution on the network with the location of the charging 257 facilities obtained by step 3. The detailed steps are as follows: 258 Step 4.1: Check feasible paths and update the path set. If the length of any sub-path is greater 259 than , remove the path from the initial set of paths. The program stops if there is no feasible path 260 between any OD pair. If there is at least one feasible path between each OD pair, proceed to the next 261 step. 262 Step 4.2: Carry out random loading of traffic flow in the network with charging facilities. 263 Calculate the generalized path travel cost ̅ , elastic flow demand between OD pairs, the probability 264 that the path is selected to obtain updated road link flows (2) and (2). Then set the iteration 265 counter n = 1. 266 Step 4.3: Repeat the random loading of step 4.2 to obtain additional link flow { }, { }. 267 Step 4.4: Use the predetermined step size sequence { }:  = 1/ , n= 1,2,...,∞.

268
Step respectively. Otherwise, set n=n+ 1 and go to step 4.3. 273 Step 5: Repeat step 3 and update charging facilities' location. The program stops until the 274 location of the charging facilities is no longer changed; that is, the maximum coverage flow remains 275 stable; otherwise, go to step 4. 276 Step 1: Input the road network parameters, such as the length of the link and the capacity of the links, and find all the paths between the OD pair w, and record them as the initial path set.

Numerical Analysis 277
Step 2: Initialize. There is no charging facility in the road network and relax the range limit of EVs. Based on the travel time of the zero-flow link t a (0) = t 0 a , calculate the effective path impedance c w k , and obtain the initial link flows V ag (1) and V ae (1). Let the upper generation counter z = 1.
Step 3: Sort all the links flow of EVs in ascending order to find the top p ranked links, and arrange the charging facilities in the middle of the p links. Then increase the upper iteration counter by 1.
Step 4: Perform SUE-ED traffic distribution on the network with the location of the charging facilities obtained by step 3. The detailed steps are as follows: Step 4.1: Check feasible paths and update the path set. If the length of any sub-path is greater than R e , remove the path from the initial set of paths. The program stops if there is no feasible path between any OD pair. If there is at least one feasible path between each OD pair, proceed to the next step.
Step 4.2: Carry out random loading of traffic flow in the network with charging facilities. Calculate the generalized path travel cost c w ke , elastic flow demand between OD pairs, the probability that the path is selected to obtain updated road link flows V ag (2) and V ae (2). Then set the iteration counter n = 1.
Step 4.3: Repeat the random loading of step 4.2 to obtain additional link flow {y ag }, {y ae }.
Step 4.5: Calculate the current flow of links by Method of Successive Average [49]. V n+1 ag = V n ag + ( 1 n )(y n ag − V n ag ); V n+1 ae = V n ae + ( 1 n )(y n ae − V n ae ).
Step 4.6: If the convergence condition is met where ω indicates the convergence accuracy, then proceed to step 5. {V n+1 ag }, {V n+1 ae } are the sets of balanced link flow for the GVs and EVs, respectively. Otherwise, set n = n + 1 and go to step 4.3.
Step 5: Repeat step 3 and update charging facilities' location. The program stops until the location of the charging facilities is no longer changed; that is, the maximum coverage flow remains stable; otherwise, go to step 4.

Numerical Analysis
The model is applied to the Nguyen-Dupuis network [32,45] for a case study in this section, which has been widely used in transportation network researches in the past decades. Due to the small scale of the Nguyen-Dupuis network, the paths can be enumerated to better analyze the relationship between the location of the charging station and the feasible path and the traffic volume of the section. At the same time, the elastic demand between OD pairs, the effects of charging speed, range limitation, charging facilities utility, and waiting cost on the location of charging facilities were evaluated. The test network (shown in Figure 3) consisted of 13 nodes, 19 segments, 25 paths (shown in Table 2), and 4 OD pairs (1,3), (1,2), (4,3), (4,2). The traveler's OD demand function uses a linear function. The familiarity of the two types of travelers to the network was measured by the parameters θ e and θ g , respectively. The parameters were set as follows during the trial calculation. q w e ( C we ) = 400 − 7 C we , q w g ( C wg ) = 400 − 7 C wg , θ e = 0.1, θ g = 0.1, U = 5, R e = 20, ε = 1, K = 0.5, p = 3. The trial results are shown in Table 3.
Energies 2020, 13, x FOR PEER REVIEW 9 of 16 Energies 2020, 13, x; doi: FOR PEER REVIEW www.mdpi.com/journal/energies The model is applied to the Nguyen-Dupuis network [32,45] for a case study in this section, 280 which has been widely used in transportation network researches in the past decades. Due to the 281 small scale of the Nguyen-Dupuis network, the paths can be enumerated to better analyze the 282 relationship between the location of the charging station and the feasible path and the traffic volume 283 of the section. At the same time, the elastic demand between OD pairs, the effects of charging speed, 284 range limitation, charging facilities utility, and waiting cost on the location of charging facilities were 285 evaluated. The test network (shown in Figure 3) consisted of 13 nodes, 19 segments, 25 paths (shown 286 in Table 2 Table 3.
Energies 2020, 13, x; doi: FOR PEER REVIEW www.mdpi.com/journal/energies link flow. For EVs, the feasible paths were less than GVs because of the range limit, and the travel 312 cost was not only affected by the link flow as GVs, but also affected by the utility of charging facilities, 313 charging speed, and waiting cost at charging stations. Therefore, the travel demand of EVs was much 314 lower than GVs, which inspired us to increase the range of EVs and the accessibility of charging 315 facilities to promote the travel demand and development of EVs. The elastic fluctuation of flow 316 demands affect the distribution of link flow and further affect the travel cost and the location of 317 charging facilities, which will react to the travel demand. In the continuous iteration and mutual 318 influence, the final equilibrium was reached. The change trend of elastic demand in other OD pairs 319 was consistent with the OD pair (1,3), which was matched with the actual travel demand. 320 321

322
The results of sensitivity analysis about range limits, charging speed, charging facilities utility, 323 and charging stations waiting cost are shown in Figure 5 -8. 324 In Figure 5, the impact of EV range limit on the locations of charging facilities is revealed. When 325 the range limits change, the charging facilities were located differently. As the range limit increased, 326 there were more feasible paths within a certain range, and the EV covered flow first increases and 327 then decreased, and finally increases when = 45. 328 The results of sensitivity analysis about range limits, charging speed, charging facilities utility, and charging stations waiting cost are shown in Figures 5-8.
In Figure 5, the impact of EV range limit on the locations of charging facilities is revealed. When the range limits change, the charging facilities were located differently. As the range limit increased, there were more feasible paths within a certain range, and the EV covered flow first increases and then decreased, and finally increases when R e = 45.
Energies 2020, 13, x FOR PEER REVIEW 12 of 16 329 Figure 5. The sensitivity analysis of the range limit.

330
Then the effect of charging speed on the location of the charging facilities is examined in Figure  331 6. Formula  Then the effect of charging speed on the location of the charging facilities is examined in Figure 6. Formula t w ck = ε·(l w k − R e ) indicates that the charging speed directly affects the charging time and general travel cost, which results in different charging facility layouts. Different values of ε can present different charging methods. ε = 0.1 means fast charging and ε = 10 means slow charging. EV drivers tend to choose a path with a long travel distance without charging facilities, rather than a path with a very slow charging speed. Then the effect of charging speed on the location of the charging facilities is examined in Figure  331 6. Formula =  · ( − ) indicates that the charging speed directly affects the charging time and 332 general travel cost, which results in different charging facility layouts. Different values of ε can 333 present different charging methods. ε = 0.1 means fast charging and ε = 10 means slow charging. EV 334 drivers tend to choose a path with a long travel distance without charging facilities, rather than a 335 path with a very slow charging speed. 336 337 Figure 6. The sensitivity analysis of charging speed.

338
In Figure 7, the utility of the charging facilities reflects the attractiveness of the link for EV users. 339 The greater the utility value, the higher the level of anxiety for EV users, so they are more likely to 340 choose a path with charging facilities. If there are many types of EVs in the network, the utility value 341 may not be very large for EVs with large battery capacity. 342 In Figure 7, the utility of the charging facilities reflects the attractiveness of the link for EV users. The greater the utility value, the higher the level of anxiety for EV users, so they are more likely to choose a path with charging facilities. If there are many types of EVs in the network, the utility value may not be very large for EVs with large battery capacity.

344
The relationship between the charging stations waiting cost coefficient K and the charging 345 station position is shown in Figure 8. When 0 < K < 1, the charging cost is reduced. When K > 1, the 346 charging cost increases, and the total coverage of the electric vehicle begins to decrease. When 347 deploying a fast charging station, a slightly larger capacity should be considered to reduce the 348 occurrence of K > 1. 349 350 Figure 7. The sensitivity analysis of charging facilities' utility.
The relationship between the charging stations waiting cost coefficient K and the charging station position is shown in Figure 8. When 0 < K < 1, the charging cost is reduced. When K > 1, the charging cost increases, and the total coverage of the electric vehicle begins to decrease. When deploying a fast charging station, a slightly larger capacity should be considered to reduce the occurrence of K > 1.

344
The relationship between the charging stations waiting cost coefficient K and the charging 345 station position is shown in Figure 8. When 0 < K < 1, the charging cost is reduced. When K > 1, the 346 charging cost increases, and the total coverage of the electric vehicle begins to decrease. When 347 deploying a fast charging station, a slightly larger capacity should be considered to reduce the 348 occurrence of K > 1. interact with each other and continuously iterate to reach equilibrium, which was more consistent 361 with the actual travel situation. A hybrid integer nonlinear programming method based on the 362 method of successive average (MSA) was constructed to prove the equivalence and the uniqueness 363 of the SUE-ED model with range constraints. Finally, a network trial was conducted to examine the 364

Conclusions
This paper studies the deployment of public charging stations in a hybrid network to maximize the service efficiency of the charging facilities. A bilevel model was proposed to depict the interaction between the mixed link flow and the location of the charging stations. The link flow obtained by the lower SUE-ED model is the key to determine the location of the CSs. In the upper level model, the CSs were arranged on the link flow ranking p to achieve the maximum coverage. Four important factors were taken into accounts in the lower level SUE model, including the range limit of EVs that affected the path choice, the elastic travel demand closely related to the distribution of link flow, the road congestion effect, and the capacity of charging facilities in travel costs. These four elements interact with each other and continuously iterate to reach equilibrium, which was more consistent with the actual travel situation. A hybrid integer nonlinear programming method based on the method of successive average (MSA) was constructed to prove the equivalence and the uniqueness of the SUE-ED model with range constraints. Finally, a network trial was conducted to examine the impact of elastic demand between OD pairs, the range limit, charging speed, the charging facilities' utility, and waiting cost on the location problem.
It should be noted that the actual road network systems are very diverse, especially in urban areas, and significantly differ from the example of the Nguyen-Dupuis Network. More realistic factors are required to be considered in the future when designing the location of charging stations, not only by technical and operational factors but also by social factors. And assumptions can be relaxed appropriately in future work. The nonlinearity of the charging time, the uncertainties in EVs energy consumption, as well as the bounded rationality of EV travelers, should be considered.