A Colored Traveling Salesman Problem with Varying City Colors

Anhui Provincial Key Laboratory of Power Electronics, and Motion Control, Anhui University of Techonlogy, Maanshan 243032, China Ministry of Education Key Laboratory of Measurement and Control of Complex Engineering Systems, Southeast University, Nanjing 210096, China Southeast University Shenzhen Research Institute, Shenzhen 518063, China Department of Electrical and Computer Engineering, New Jersey Institute of Technology, Newark 07102, NJ, USA 4e Institute of Systems Engineering and Collaborative Laboratory for Intelligent Science and Systems, Macau University of Science and Technology, Macau 999078, China


Introduction
A colored traveling salesman problem (CTSP) generalizes the well-known multiple traveling salesman problem [1]. Each salesman in it has a single color, and each city has one to multiple salesmen's colors. A city allows only one salesman of the same color to visit exactly once. Colors are used to characterize the cities' accessibility for salesmen and can be applied to some real-world scheduling problems in which the tasks to be performed own different accessibilities for movable processing devices [2]. Li and Meng apply the CTSP to schedule multibridge machining systems for avoiding collision [3,4]. In addition, He and Hao [5] propose a two-phase local search for solving the CTSP. Recently, three CTSP variants have been proposed, that is, large-scale general CTSP [6], biobjective CTSP [7], and precedence-constrained CTSP [8]. [22] analyzes the effects of the edge weight changes on the fitness landscapes of TSP. Recently, Mavrovouniotis et al. [23] investigate a memetic ant colony optimization algorithm with local search to address symmetric and asymmetric DTSPs, and the experimental results show the efficiency of their proposed algorithm. Chen et al. [16] investigate the value of choosing the next stop to visit in a multistop trip based on current traffic conditions to minimize the expected total travel time of the tour, and they model this problem as a Markov decision process. In dynamic scenarios, datasets (city coordinate information) are generally incomplete and variable. Baykasoglu and Durmuşoglu [24] address a DTSP in which cities may be added to/ deleted from the city domain with a multiagent-based approach.
In logistics and distribution, traffic network topology and customer demands are changing over time. If some practical constraints in logistics distribution are introduced into DTSP, it becomes a dynamic vehicle routing problem (DVRP) [25]. DVRP requires solving not only a dynamic path planning problem with varying edge weights but also the problem with some uncertain constraints in real delivery. It is more complex than the DTSP. DVRP can be divided into two categories according to the customer's location information. One is those with known customer locations, but a traffic network topology varies and customer demands are random. Kim et al. [26] propose a DVRP model with nonstationary stochastic travel time under traffic congestion, and a Markov decision process model is adopted to solve this problem. A vehicle routing problem with stochastic demands (VRPSDs) is the one in which customer locations are known, but the demands of customers are random. For example, a flexible solution method is adopted to transform VRPSD into a small set of capacitated VRP, where Monte Carlo simulation with statistical distributions is used to estimate customer demands [27]. Mavrovouniotis et al. [23] present a DVRP with stochastic demands where customers arrive dynamically and randomly, and thereby their requests are not completely known. To solve it, an ant algorithm with an immigrant scheme is developed. Customer service time may be uncertain in a DVRP. Vonolfen and Affenzeller [28] resolve a dynamic pickup and delivery problem with time windows by using a novel waiting for heuristic strategy based on historical request information. Juan et al. [27] propose a method to estimate truck arrival time at each customer location for truck route planning in nonstationary stochastic networks.
Another DVRP category is the real-time or online VRP, that is, new customer orders are randomly generated, and both customer locations and demands are unknown. In this situation, it is necessary to reschedule vehicles based on previous schedules, that is, redistribute these vehicles for all not-served customers. Many scholars have studied the rule of accepting or rejecting new service requirements in a DVRP [29][30][31]. Bopardikar et al. [32] investigate a DVRP with the Poisson distribution orders and design routing policies for a service vehicle to maximize the number of demands at the steady state. Pavone et al. [33] propose a distributed and adaptive segmentation algorithm to solve a DVRP and validate its effectiveness. e above studies focus on path optimization over traffic networks with varying edge weights and all the cities with the same accessibility. So far, there is no research on the optimization problem with dynamic city accessibilities, which can well describe a supply and distribution problem faced by online retailers. erefore, this work focuses on delivery route optimization for online retail with multiwarehouse and dynamic task accessibility.
For example, aiming at delivering goods faster and better, online retailers tend to build their own logistics system consisting of widely located warehouses in different countries or regions containing goods of different types and stocks. Amazon owns a total of 90 warehouses in the world, JD Company has built 143 warehouses in 43 cities of China, and Suning Logistics Group has 12 large-scale warehouses in China.
Once random customer orders are generated, an online retailer has to provide the customers with their ordered goods. e stocks of the goods differ at distinct warehouses and are constantly changing with sales and supplies. erefore, decision-making is needed to decide which warehouse should provide the goods to a specific customer, and optimize the delivery routes of the different goods to different customers to meet customer requirements and reduce the supply cost. So far, there is no satisfying decision-making mechanism for this kind of online retail.
Suppose that each warehouse has one unique color. It can be considered as a salesman with the color to deliver goods for customers by its own transport system. ose online ordering customers, located in different areas, can be provided with the ordered goods stocked at one of the different warehouses. In other words, such a customer carries the warehouses' colors, and thereby the customers can be regarded as the cities with one to multiple warehouses' colors. us, the warehouses' supply problem can be described as a CTSP. With the customers' online purchases, each warehouse's stock changes over time. When a commodity ordered by a customer is out of stock in a warehouse, the customer loses his/her color. is means that a customer whose warehouses' colors will change in the sales process. erefore, we propose a new CTSP with varying city colors (CTSP-VCC). It can be applied to supply decision and delivery route optimization for online retail with multiwarehouse and varying task accessibility.
Considering the outstanding performance of Variable Neighborhood Search (VNS) [2] in solving static CTSP, we adopt it to solve CTSP-VCC. To overcome its instability under a dynamic circumstance, greedy initialization and best population immigrant are incorporated into it to maintain its fast convergence and solution diversity.
is work intends to make the following contributions: (1) It formulates a CTSP-VCC in which city colors are changing over time. CTSP-VCC can be applied to supply decision and delivery route optimization for a online retail with multiwarehouse. (2) It constructs a rigorous mathematical optimization problem for CTSP-VCC and designs a dynamic environment simulator to mimic the change of city colors.
2 Discrete Dynamics in Nature and Society (3) It develops an improved VNS to solve CTSP-VCC in which a greedy operation and a novel population immigrant scheme are adopted to achieve better convergence and solution quality. e remainder of this paper is organized as follows: Section 2 introduces CTSP-VCC. Section 3 presents a greedy-initialization VNS with a best population immigrant scheme. Section 4 shows simulation results. Section 5 concludes the work.

Problem Formulation
CTSP uses the colors to characterize city accessibility to salesmen [1,2]. It owns m salesmen and n cities, where m < n, m, n ∈ Z � 1, 2, . . . { }. It can be formulated over a digraph G � (V, E) with vertex set V � 1, 2, . . . , n { } representing n cities; each edge (i, j) ∈ E and i ≠ j represents a visit cost between cities i and j, associated with weight ω ij . Salesman k∀k ∈ Z m � 1, 2, . . . , m { } is allocated with color c k . A city color matrix (CCM) M is defined to indicate the accessibility of cities to salesmen.
where a row vector represents the color set allocated to city i, Clearly,1 ≤ |c(i)| ≤ m, namely, city v i allows at least one salesman and at most m salesmen to visit. e j-th column vector in M expresses which cities in the city set can be visited by salesman j. e accessible city set with respect to salesman j is denoted by In addition, we define the set of single-colored cities T j � i|a ij � 1, j a ij � 1, i ∈ V to collect the cities that carry the same color c j , j ∈ Z m . Here, we define a dynamic CTSP whose city colors are changing over time. CCM in CTSP-VCC is not a constant matrix, that is, a ij , i ∈ V, j ∈ Z m , becomes a function a ij (t) that varies over time. Let a ij (t) � 1, if city i has color c j at time t; otherwise a ij (t) � 0.
Next, a 0-1 integer programming model of CTSP-VCC is given as follows: the binary access variables Subject to: Each salesman k is required to start from and return to depot d k ∈ V. e equation is as follows: Salesman k cannot visit a city whose color set does not includec k : Each city can be visited exactly once by one and only properly colored salesman: A salesman may visit a multicolor city, and a pair of entry and exit is required if it is the case: Any solution consisting of several disconnected subtours for a salesman must be forbidden: where u ik is the number of nodes visited in the tour of salesman k from d k to the node i.

Lemma 1. CTSP is a special case of CTSP-VCC.
Proof. In CCM, a ij (t) denotes the variable color of city in CTSP-VCC. If a ij (t) is an invariant constant (i.e., a ij (t) � a ij (t 0 )), it means that the model turns to be a stationary one, i.e., CTSP.

Lemma 2.
e solution space size of CTSP-VCC is m k�1 |R k |!, where |R k | denotes the number of cities visited by salesman k, and m k�1 |R k | � n.
Proof. Let R k denote the route of salesman k. |R k | represents the number of cities in route k, and it is a variable, Proof. CTSP is NP-hard [1], and it has been proved to be a special case of CTSP-VCC by Lemma 1. e solution space size of CTSP-VCC is m k�1 |R k |! given by Lemma 2, and its optimal solution varies over time, that is, with the variable city color a ij (t), the optimal solution varies in its solution space. erefore, CTSP-VCC is also NP-hard.
For this type of dynamic optimization problem, it is challenging to develop algorithms to track the optimal solution changes. In particular, the solution feasibility changes represent changes in city colors, (i.e., a feasible original Discrete Dynamics in Nature and Society solution path may become no longer feasible because of changes in city colors).

Enhanced VNS for CTSP-VCC
VNS has been proved to be an efficient method to address CTSP [2], so we adopt it to resolve CTSP-VCC. Considering the optimal solution of CTSP-VCC is no longer a fixed solution and varies over time, we introduce a best-population immigrant scheme and greedy-initialization method into VNS to enhance global searching ability. Furthermore, given that the varying city colors can change the solution feasibility, we introduce a verification operation into VNS to judge the solution feasibility. Some solutions must be deleted if it is infeasible after changes in city colors.

Solution Encoding.
e direct-route encoding is adopted to encode each salesman's route directly for CTSP-VCC. A coding example of CTSP-VCC with n � 10 and m � 3 (see Figure 1), namely, Cities 2-8-7-1 and 10-4-3, 9-5-6 are visited by salesmen 1-3, respectively. Let R k denote the route of salesman k, excluding his depot. e solution is denoted by R � (R 1 , R 2 , . . ., R m ). e permutation of cities in the search space size is

Basic VNS.
e core thought of VNS is to search a better solution in the neighborhood of a high-quality solution. First, to ensure the population diversity of VNS, a random initialization method is used in VNS to generate initial individuals, that is, all salesmen select randomly the cities carrying the same color to establish their routes. en, shaking and local search operations are performed in each generation of VNS, and the current best solution is updated if a better one is obtained. On one hand, Interchange and Relocation are adopted to perturb the current incumbent solution to generate its neighborhood solution set in the shaking operation (see . e detailed operations of them are given in [2]. On the other hand, a local search explores the neighborhood of the solution selected by shaking to obtain a better one. Particularly, two operations, that is, city-removal and reinsertion operations, carry the deep neighborhood search. Once a better solution is obtained in the local search, it is saved for the next iteration.
In addition, if the environment changes, VNS needs to analyze the feasibility of the original solution. If the solution is feasible and better than the current solutions, it is used to update the best solution; otherwise, it is deleted. e flowchart of basic VNS is given in Figure S1 in the Supplementary Material.

Improved VNS.
A greedy-initialization scheme and a best-population immigrant scheme are introduced to improve the performance of basic VNS. Greedy-initialization means that a two-stage greedy is adopted to generate a better initial solution, that is, the greedy search with the minimum insertion for single color cities first and the multicolor cities next. An immigrant strategy with the best population is used to prevent the premature convergence of VNS. First, the current best individual is selected to perform immigrants after VNS's local search. After the immigrant operation, several new individuals are generated, and they are used to replace those original individuals which are of poor quality. e detailed operations of the best-population immigrant scheme are shown in Algorithm 1.
On one hand, the greedy operation can be used to generate an initial solution of high quality during the initialization of VNS, and a good initial solution can benefit the quick search for a satisfactory solution. On the other hand, the best-population immigrant scheme can prevent premature convergence in a single search direction, which can help VNS to jump out of a local optimum and search more space.
e flowchart of greedy ability VNS with a bestpopulation immigrant scheme is illustrated in Figure 5.

Time Complexity of Improved VNS.
We name the improved VNSs as VNS, VNSB, GVNS, and GVNSB. ey represent in turn the random-initialized VNS, VNS with the best population immigrant, greedy-initialization VNS, and greedy-initialization VNS with a best population immigrant scheme.
Let g and p denote the generations and population size, respectively, and the time complexity of VNSs are investigated as follows: the greedy-initialization takes O(|R 1 | 2 + · · · + |R m | 2 ) � O( m k�1 |R k | 2 ) time in constructing the initial routs. In VNS, the shaking operation needs O(n) time. e least-cost insertion is adopted in the local search of VNS, and it requires O(n) time. In addition, in the stage of Best-Population Immigrant, the rank operation and 2-opt operation take O((p)log(p)) and O( m k�1 |R k | 2 ), respectively. erefore, we can obtain that GVNSB takes

Dynamic Environmental
Simulations. e dynamics of CTSP-VCC are reflected by its city colors' change over time. First, let M change over time whereM � [a ij (t)] n×m . en, to describe the change of city colors, each a ij (t) is changed every τ generations during the run of algorithms with the probability ρ, where a ij (t) � 1 , 0 < ρ < 1, 0 < i < n, 0 < j < m, and a ij (t 0 ) is the initial value. e environmental dynamics can be easily tuned by two parameters, τ and ρ. e former controls the change speed, while the latter manipulates the severity each time the environment changes. e dynamic environment of CTSP-VCC is generated by Algorithm 2.

Experiment Design.
We have implemented all the algorithms in C++ with Microsoft Visual Studio 2012 and executed them on Dell Computer Optiplex 390 running Windows 7 with CPU Intel Core i 3 at 3.40 GHz and RAM 4 GB. e offline performance F BOG is adopted to evaluate the overall performance of an algorithm: Let N � 50 denote the total number of runs, and G � 10 × τrepresents the total number of generations for a run. F BOGij is the best-of-generation fitness of generation i at run j, and F BOG is the average fitness of 50 runs.
To test the performances of the algorithms, our experiments are designed as follows: first, test cases, Eil101, A280, and Rat575 are derived from the datasets in TSPLIB [34] while their CCMs are generated randomly according to the number of cities and salesmen. en, the two parameters for dynamic environment control τ is set to be 30 or 50, and ρ is set at 4 different levels, i.e., 0.1, 0.2, 0.5 and 1.0. In addition, each VNS is executed independently 50 times for each test case. One-tailed t-test with 98 degrees of freedom at a 0.05 level of significance is adopted to compare the statistical results of VNSs in dynamic environments.

Comparison among the Proposed Algorithms.
To compare the performances among the proposed algorithms, the statistical analysis of the experimental results is given as follows. Figures 6 and 7 show the offline performances, and Figures 8 and 9 illustrates the dynamic performances for each test. In the dynamic environments with different change speeds and intensities, the statistical results of compared algorithms by one-tailed t-test with 98 degrees of freedom at a 0.05 level of significance are given in Tables 1 and 2, where "s+", "s-", and "∼" denote that one (1) Begin (2) Input: x 1 , x 2 , . . . , x ps , k � 0; (3) Rank the ps solutions from small to large, where x 1 is the smallest one and x ps is the largest one.
Discrete Dynamics in Nature and Society algorithm is significantly better than, significantly worse than, and similar to another algorithm, respectively. e standard deviations of 50 results of the algorithms for each test are given in Tables 3 and 4.

Effectiveness of Greedy-Initialization.
e offline performances illustrated in Figures 6 and 7 show that GVNS outperforms VNS.
For example, the result from GVNS in Eil101-6, at τ � 50 and ρ � 0.1, is 2515.8 and it is 15.1% less than that from VNS. In A280-15, the result from GVNS is 5836.6 and it is only 22.4% of that from VNS.
In Figure 9, the dynamic performances show that, the greedy-initialization can effectively improve the initial solution quality of VNS. For example, the results obtained by GVNS in the environment at τ � 50and ρ � 0.2 are 2357.6, (1) Begin (2) Input: gen; End If (13) End If (14) End If (15) Output: a ij (t).  Discrete Dynamics in Nature and Society 8029.4, and 22478.7, that is, 82%, 34%, and 28% of those by VNS in Eil101-5, A280-8, and Rat575-10, respectively. Similarly, the t-test results illustrated in Tables 1 and 2 show that GVNSB outperforms VNS, VNSB, and GVNS significantly. erefore, we can draw a conclusion that the greedyinitialization is very effective for VNS to solve CTSP-VCC. e reason is that a better initial solution is critically important for population evolution, and it can guide individuals to perform their search in a more promising solution space of its neighborhood. Indeed, VNS can find a better solution quickly and effectively at the initial stage by means of greedy-initialization, and it also can obtain other diverse solutions by random initialization. e greedy initialization can guarantee VNS not only an initial solution quality but also the diversity of the population. It is why GVNS and GVNSB perform better than VNS and VNSB.

Effectiveness of a Best Population Immigrant Scheme.
e results of algorithms in a dynamic environment with τ � 50 are illustrated in Figure 9. Specially, in Eil101-4 at ρ � 0.1, the results obtained by VNS and VNSB are 2670.2 and 1398.9, respectively. VNSB's result is improved about 47.6% after introducing a best population immigrant scheme. Furthermore, the result of GVNSB is 1293.8, just 62.8% of that with GVNS. In addition, the result in A280-8 shows that VNSB obtains 9783.9 and is 39.8% of that from Discrete Dynamics in Nature and Society VNS, and GVNSB obtains 5328.5 and it is 58.7% of that from GVNS. Besides, the results in three cases of Rat 575 show that we obtain better results after introducing the best population immigrant scheme into VNS. All the curves in Figure 6 show that the fitness decreases as ρ increases. e reason is that the increase of city colors with ρ enlarges the solution space of a problem. us, the algorithms may obtain better solutions in the new solution spaces. Furthermore, Figure 7 shows the results at τ � 30, like Figure 6, where VNS and GVNS with a best population immigrant scheme can obtain better results than those of VNS and GVNS. Figures 8 and 9 illustrate the dynamic performances in some test cases. e algorithms have different convergences.
Obviously, the algorithms converge well when the best population immigrant scheme is applied. For example, the fitness value of the first generation in VNS is 34764 and 112337 in A280-15 and Rat575-15, at dynamic environment τ � 50 and ρ � 0.5, respectively. After one generation, they are 33745.6 and 108974, respectively. Although the fitness of value in the first generation in VNSB is similar to that in VNS, after one-generation evolution, they are 56.4% and 60.2% of those in VNS, respectively. Table 3 shows the standard deviation of 50 runs for each test case. e obtained results illustrate that the standard deviations of each algorithm are markedly different at the beginning. But they decrease significantly by introducing the best population immigrant scheme into the algorithms. For  Tables 1 and 2 illustrate the t-test results of the algorithms in dynamic environments with ρ increasing from 0.1 to 1.0, while τ � 50 and 30. Concretely, the t-test results of Eil101, A280, and Rat575 show that VNSB is significantly better than VNS, and GVNSB is significantly better than GVNS. e results are sufficient for us to validate that the best population immigrant scheme is very effective for VNS to solve CTSP-VCC. Tables 1 and 2, the performance of GVNS is significantly worse than that of VNSB in three cases of Eil101 with ρ � 0.1, 0.2, and 0.5. is result is also illustrated more clearly in Figure 6. In particular, in dynamic environment at τ � 50 and ρ � 0.5, the results in GVNS in Eil101-4, Eil101-5, and Eil101-6 are 2342, 2015, and 2132.1, respectively. ose results obtained by VNSB are 1494.5, 1561.8, and 1734.6, respectively, all less than their corresponding results by GVNS.

Performance Comparison between Greedy-Initialization and Best-Population Immigrant. In
In addition, we find that the pros and cons of GVNS and VNSB swap with the increase of the number of cities. For example, in Figure 7, the results obtained by GVNS are better than those by VNSB in A280-8 and A280-15, but in A280-5 with ρ � 0.1 and 0.2, the result from GVNS is worse than that from VNSB. While in the dynamic environment at ρ � 0.5 and 1.0, the results have changed, that is, the result    Figures 8 and 9 show that GVNSB always achieves the best results among the algorithms. Furthermore, Tables 1 and 2 suggest that the introduction of the best population immigrant scheme into GVNS can guide the algorithm to track the change of the optimal solution rapidly in a dynamic environment. Hence, GVNSB outperforms significantly other three algorithms in terms of the offline performance, dynamic performance, and standard deviation.
In summary, both the best population immigrant scheme and greedy-initialization scheme can help VNS converge rapidly and obtain much better solutions, and they are beneficial for VNS to track the dynamic environment of CTSP-VCC.

Conclusions
is work proposes a dynamic CTSP in which city colors change with time. It can be applied to a supply decision problem of online retailers with multiwarehouses and varying supply availability. en, a variable neighborhood search algorithm is proposed to address this new CTSP. It uses directroute encoding and a two-stage greedy-initialization algorithm to obtain an optimal/near-optimal solution. In addition, a best population immigrant scheme is adopted to track the change of the optimal solution rapidly in a dynamic CTSP. Furthermore, to mimic the dynamic environment of CTSP-VCC, we construct a dynamic environment simulator that can be easily tuned via two parameters, i.e., change speed and severity.
To test the performances of the proposed algorithm, i.e., greedy-initialization VNS with the best population  immigrant, extensive experiments are conducted in terms of offline performance, dynamic performance, and standard deviation. e results show that VNS can obtain a desired solution in its initialization by means of a greedy-initialization scheme. Solution quality and stability can be significantly improved with the help of the best population immigrant scheme. Specially, because of the two schemes, the results of the baseline VNS can be improved about 30%-50%. In summary, the proposed best population immigrant scheme and greedy-initialization scheme are both effective for VNS to track the change of optimal solutions in the new CTSP and help the algorithm obtain better solutions. Our ongoing work is to apply CTSP-VCC and VNSs to some real-life problems [35,36].

Data Availability
e data used to support the findings of this study have been deposited in the Github (Link:https://github.com/optmlearing/vehicle-routing-problems).

Conflicts of Interest
e authors declare that they have no conflicts of interest.