Comparing Greedy Constructive Heuristic Subtour Elimination Methods for the Traveling Salesman Problem

This paper further defines the class of fragment constructive heuristics used to compute feasible solutions for the Traveling Salesman Problem into arc-greedy and node-greedy subclasses. Since these subclasses of heuristics can create subtours, two known methodologies for subtour elimination on symmetric instances are reviewed and are expanded to cover asymmetric problem instances. This paper introduces a third novel methodology, the Greedy Tracker, and compares it to both known methodologies. Computational results are generated across multiple symmetric and asymmetric instances. The results demonstrate the Greedy Tracker is the fastest method for preventing subtours for instances below 400 nodes. A distinction between fragment constructive heuristics and the subtour elimination methodology used to ensure the feasibility of resulting solutions enables the introduction of a new node-greedy fragment heuristic called Ordered Greedy.


Introduction
Applegate et al [1] describes the traveling salesman problem as, "Given a set of cites along with the cost of travel between each pair of them, the traveling salesman problem, or TSP for short, is the problem of finding the cheapest way of visiting all the cities and returning to the starting point." It can also be mathematically defined as, given a complete graph G = (V, E), cities are represented via the graph vertices g ∈ G , and edges e ∈ E represent the paths between the cities where the edge weights are the distances between each city: What is the shortest tour that visits all vertices once and returns to the starting vertex? When traditional linear programming methods were applied to the TSP, intractability issues arose [1]. It has since been shown this is because the TSP falls into a class of known computationally 'hard' problems called NP-Complete [2]. No one has yet developed an efficient method for universally solving large instances of NP-complete problems to optimality [3]. The inclusion of the TSP in the set of NP-complete problems motivates the usage of other solving techniques such as heuristics.
One class of heuristics utilize a greedy-type methodology, where the best choice, according to a predefined parameter, is selected at each step of the method. An example of a greedy-type method for the TSP is the Multiple-Fragment heuristic [4], where the shortest available arc is iteratively added to form a tour. However, this greedy heuristic runs the risk of generating subtours, or disconnected tours of less than size N (where N is the number of nodes present in the graph) that prevent a single continuous tour from being formed. Some research has been completed to develop methodologies that avoid subtours when utilizing the arc-greedy heuristic [4,5].
This paper further defines the class of Fragment Heuristics defined by Bentley [4] and discusses the importance of subtour elimination methodologies within this class. A novel subtour elimination methodology for the arc-greedy heuristic is introduced and compared to two known subtour elimination methodologies. Extensions to all three subtour elimination methodologies to handle asymmetric instances are introduced. Computational results are generated across multiple TSP instances for each method to compare run times by size and geometry.
2 Literature Review

Types of TSP Construction Heuristics
Bentley's [4] 1990s paper, in addition to many other contributions, established classifications for TSP constructive heuristics. Bentley breaks TSP construction heuristics into three general categories: heuristics that grow fragments, heuristics that grow tours, and heuristics based on trees. The latter two classes each ensure viable TSP tours utilizing a consistent methodology for each heuristic in its class. For example, the class "heuristics that grow tours" contains heuristics that follow the Insertion or Addition expansion rules. These heuristics each start with a partial tour consisting of a single node which is then gradually grown by adding new nodes. This is accomplished by deleting an edge and adding two new edges, according to the specified rule-set, to connect the new node to the current tour, ensuring a valid TSP tour is constructed. However, the "heuristics that grow fragments" class does not have a consistent methodology to avoid premature partial circuits, or subtours.

Heuristics that Grow Fragments
In this section we will briefly introduce the heuristics described in the "Heuristics that Grow Fragments" section of Bentley's [4] paper.

Nearest Neighbor
The nearest neighbor (NN) heuristic was first applied to the TSP in a 1954 paper by Flood [6] but was introduced as the "next closest city method." The process was later refined by Dacey [7] and coined with its eventual name. The NN starts at an arbitrary city, and successively visits the closest unvisited city. Note that the nearest neighbor heuristic maintains a single path fragment that originates at the predetermined starting city, and cannot be closed into a cycle until every node has been visited. Therefore the decision of "which arc to add" is limited to only those arcs that leave the current end node of the fragment, this yields an algorithm run time of O(N 2 ). Future work by Bentley [4] allowed this heuristic to perform in O(N log N). This methodology allows NN to quickly create an initial tour which avoids sub tours. However, NN is extremely sensitive to the choice of starting node, especially in larger instances. Bentley also introduces a second variant of NN called the Double-Ended Nearest Neighbor (DENN) which allows the fragment to grow from both ends.

Multiple-Fragment
The Multiple-Fragment heuristic (MF) was first introduced by Papadimitiou and Steiglitz [8] as a modification of a process first seen in a 1968 paper by Steiglitz and Weiner [9]. The heuristic is a more complex greedy-type TSP heuristic where all edges of the graph are sorted from shortest to longest. Edges are then added to the tour starting with the shortest arc as long as the addition of this arc will not make it impossible to complete a tour. Specifically, this means avoiding adding edges that make early cycles, and also avoiding creation of vertices of degree three. This process, as originally proposed, required O(N 2 log N) time. However, Bentley was able to speed up this process to O(N log N) [4] in a paper introducing his MF version. This yields a similar run time to NN while maintaining a similar worst case solution quality. MF's tour construction methodology causes the heuristic to only produce a single solution for each instance where NN can arrive at different solutions based on a different starting point. When compared to the average NN solution over all starting points, MF tends to outperform NN on an instance-to-instance basis [4,10].

Methodology/Discussion
At present it is entirely unclear when referencing the "Greedy Heuristic" for the TSP if the heuristic under discussion is NN or MF. This is acknowledged in [11]. Much of this confusions seems to draw from the poor naming conventions used with respect to how a greedy methodology works. In accordance with the framework established by Talibi [12] constructive heuristics involve two choices. First, determine a set of elements, S j = {e 1,j , e 2,j , ..., e p,j } which comprise the neighborhood of the current solution, that is the set of possible choices at each iteration j. Second, define a methodology to choose an element e i,j from this set S j . Thus the framework of how this set S j is defined determines in large part how the heuristic constructs the solution. In the case of these TSP fragment heuristics we propose a more comprehensive classification that splits this class into two greedy-type methodologies; node-greedy and arc-greedy heuristics based on how the set S j is defined.

Node-Greedy
For the first proposed fragment greedy type, node-greedy, the set of arcs available at each iteration is defined by a given node. The key defining characteristic of this sub-class is that the set S j is limited to the arcs incident to a given node. Depending on the measure of merit being greedily optimized this sub-class results in different heuristics. NN is an example of this node-greedy sub-class, where the starting node is user defined, from that point on the elements of S j are defined as the arcs incident to the node at the head of the fragment. The choice of element from S j associated with NN is the shortest arc.

Arc-Greedy
The arc-greedy type of fragment greedy TSP heuristic considers the set of all available arcs at each iteration, which will not make it impossible to complete a tour. The first criteria to narrow the scope of this set is to identify all arcs that will not cause a node to have a degree of more than 2 when added. The second criteria, which is much more difficult, is to identify what arcs will cause a tour of less than size N, or subtour, to form. Consider the following tour construction utilizing the multiple-fragment arc-greedy constructive heuristic methodology on a 5 node TSP instance with only the degree criteria. After adding the first two shortest arcs A-B and B-C, we can see from the distance matrix that arc A-C is the next shortest and still ensures that all nodes in the graph do not exceed a degree of 2. However, adding this arc creates a subtour, which would prevent the heuristic from ever constructing a feasible TSP tour (see Figure 1). Part of the confusion with regards to the Figure 1 Multiple-Fragment arc-greedy heuristic's can generate subtours "Greedy Heuristic" revolves around Bentley's MF being viewed as what we describe as the arc-greedy heuristic. However, we argue the true significance of Bentley's MF is the development of a subtour tracking and elimination methodology, which can then be applied to any member of the fragment class of TSP constructive heuristics.

Subtour Tracking and Elimination Methodologies
In addition to Bentley's MF methodology [4], a well known, though to the best of our knowledge, undocumented methodology we dub exhaustive loop also can be used to track and eliminate subtours. This paper introduces a third novel method for eliminating subtours while using an arc-greedy constructive heuristic we call the Greedy Tracker (GT).

Directional vs. Non-Directional
TSP instances can be either symmetric or asymmetric, similarly subtour tracking methodologies can be either directional or non-directional. Non-directional subtour tracking methodologies construct a tour with no regard to the direction of travel for each arc while ensuring no node has a degree of more than two. This methodology can only be used with symmetric TSP instances where the distance to travel from node to node is equal in both directions. This poses some computational advantages as only n * (n + 1)/2 arcs need to be initially sorted. Directional subtour tracking methodologies can be used on either symmetric or asymmetric instances when the direction of arc travel is either of importance to the final solution and/or takes different costs to travel based on direction. In a directional scenario, all arcs of each direction n 2 − n, are sorted from shortest to longest and rather than tracking the total degree of each node, each node can only be entered and left once, ensuring a continuous direction throughout the tour.

Multiple-Fragment
The Non-Directional variant of MF is well documented throughout TSP literature [4,5,10]. However, a formalization of the Directional variant is missing from this literature. This process appears to have been utilized in a 1999 paper by Glover et al [13]. However, neither pseudocode nor an explicit description of altering the methodology for directional instances is available. The primary alterations are to consider all arcs of each direction and the Degree array is split into 'To' and 'From' arrays. The process then continues as described by Bentley [4], where after each arc is added the tails of the associated fragment are updated to ensure no subtours are formed. Pseudocode for this modified methodology is included as Algorithm 1.

Exhaustive Loop
The Exhaustive Loop (EL) is not well documented in academic literature, often simply referenced as "the standard way." To the best of our knowledge no formal coverage of this method exists in literature. EL cycles through every edge in the fragment containing the most recently added edge. Once an edge e ij is added to the partial tour, node i is identified as the "start node" and node j will be set as "current node." A trace along the current partial tour then begins. At each step of the trace the "current node," node j, is checked to see if it is connected to another node k via edge e jk in the partial tour. If it is, then node k becomes the new "current node." If the trace returns back to the "start node" in under N steps, where N = |V| (the number of nodes in the instance), then the added edge e ij has created a subtour and is an illegal edge. If no edge leaves the "current node" the addition of edge e ij is valid and the current portion of the tour is still a fragment. Each time an edge is added, a count is incremented and the process continues until N-1 edges have been added upon which the last two endpoints are connected.
When applied to the earlier example, after adding edge A-C the heuristic identifies node A as the starting node and Node C as the current node (seen in Figure 2). The heuristic finds Node C has a degree of 2 and finds the other connected arc C-B. Node B becomes the current node and the heuristic verifies that the current node is not the same as the start node (seen in Figure 3). Once again, the heuristic finds the new current node, Node B, has a degree of 2, finds the other connected arc B-A, and updates the current node to Node A. This time when the heuristic checks the current and start node, it realizes they are the same (Figure 4).    found in Algorithm 2. The EL can also be modified to handle a directional methodology by splitting the 'Degree' array into a 'To' and 'From' array and rather than checking if current node has a degree of 2, check if the current node has a value of 1 in 'From' array

Greedy Tracker
We now introduce a novel way to track the progress of the arc-greedy construction heuristic, and ensure subtours are not created. This new method is called the "greedy tracker" (GT). The GT tracks a node's connection with other nodes when constructing a TSP tour. While if Both nodes of current edge have degree <2 then 5: Set Start = Tail of current edge 6: Set Current = Head of current edge 7: while Continue = True do 8: if Current is Tail to Another Edge then 9: Set Next Node = Head of found edge 10: if Next Node = Start then Next Edge in List 24: end while 25: Connect Hamilton Path we define both non-directional and directional GT variants, it is conceptually easier to visualize the GT using its directional variant on a symmetric instance and then generalize the process for use on asymmetric instances or to the non-directional variant. Thus, the following introduction to the GT utilizes the directional variant on a symmetric matrix and is accomplished using the following structures: • X = binary n by n matrix of x ij , • F = binary n by 1 array of f i , • T = binary n by 1 array of t i , • x ij = 0 if arc from i to j is eligible, greater than 0 if not eligible, • f i = binary for whether node i has been left, • t i = binary for whether node i has been entered.
These structures track each move to prevent sub tours. Given our prior example the initial condition of GT, and associated structures, can be seen in Figure 5: The X (identity), F (From), and T (To) structures can be seen above on the left and, for ease of reference, the associated distance matrix from the TSP instance can be seen on the right. The 1s loaded on the diagonal of the X matrix (where i=j) signal ineligible moves. Note that the diagonal on the distance matrix has been colored grey, and set to zero, to correspondingly show these ineligible arcs. The distance matrix indicates that the shortest arc is either from A to B or vice versa, thus arc A to B is selected. The X, F, and T matrices are updated with 1s to indicate this move, as shown in Figure 6. Next, the column of the X matrix associated with the new arc is processed. Every row where a 1 appears is combined with the 'From' row of the created arc. Figure 7 illustrates this operation. As seen in Figure 7, since Row 2 has a 1 in the same column as our new arc, the two rows were combined so that any 1s that were in the Row 1 are now also in Row 2. Note that for the example we only show values of 1 so as to not detract from their purpose of referring to an ineligible move. However, in the code the values in each row will actually be added and values of greater than 1 will appear. For ease of reference in this example ineligible values in the distance matrix are turned grey. As can be seen in Figure 8, distances that correspond with a 1 in the X matrix are ineligible moves. Note that any row or column that has a 1 in the T or F arrays is also marked as an ineligible move. This information is utilized in the first step of the next iteration where the shortest available arc is identified. As seen in Figure 9, the shortest available arc is B-C and once again the X matrix, From, and To arrays are updated with 1s to indicate the move. The column of X associated with the "To" node of the new arc is processed and every row where a 1 appears is combined with the "From" row of the created arc which can be seen in Figure 10. All the distances that correspond with a 1 in the X matrix are marked as ineligible moves in the distance matrix, as well as any distances associated with a 1 in the T and F arrays. The resulting step can be seen in Figure 11. The grey numbers in the distance matrix indicates that adding arc A-C is no longer possible because node C already has an edge entering it. This process prevents the formation of the subtour. The process shown above continues until all nodes have been visited which creates a Hamiltonian Path. The final connection to complete the tour is made using the T and F arrays as each will have one index that is still empty.
Real-world implementation of GT requires adjustments to achieve a more linear computational growth for larger instances by reducing the total number of operations that occur within each iteration. This is achieved by removing the addition of values with respect to nodes that have been exited and entered. This decreases the dimensionality of the GT as the tour is constructed and is possible because once a node has been entered (resp. exited) no more arcs may enter (resp. exit) that node. Therefore, it is unnecessary to track what arcs could produce a subtour by entering (resp. exiting) that node. Consider for example the same 5-by-5 instance. After completing the row additions after adding arc A-B, column B can be deleted. Figure 12 shows the resulting GT and distance matrix. This process can also be applied to rows to further reduce dimensionality. When working with a non-directional instance a column and row would be deleted after the node had a degree of 2.
Breaking down the process of the row and column delete methodology leads to the realization that fundamentally after an arc is added the GT is only identifying the arc(s) which would form a subtour for the fragment including that new arc. Bentleys MF on the other hand tracks the tails of every fragment and prevents the head and tail of a fragment from connecting. There appears to be some advantage in GTs approach as only one number is added to the X matrix at each step. However, the GT must search for and identify the exact arc that needs to be moved, whereas while MF references more nodes and tails, there is no searching necessary.
Pseudocode for the GT subtour elimination method is found in Algorithm 3.

Results and Analysis
To compare the run times for the MF, EL, and GT subtour elimination methods all three methods were run across 14 symmetric and 9 asymmetric TSP instances.

TSP Instances
TSP instances are available in an online library, TSPLIB, maintained by Ruprecht-Karls-Universitat Heidelberg located in Baden-Wurttemberg, Germany [14]. For the purposes of this research testing was performed on the instances seen in Table 1, where the alpha prefix is an identifier and the numerical suffix indicates the instance size (in number of nodes).

Testing
Initial tests verified that each subtour method (MF, EL, and GT) resulted in the same tour for all TSP instances. These tests were conducted with both directional and non-directional versions of codes on symmetric TSP instances. In addition, the directional code versions were run on asymmetric TSP instances.
Once testing verified each method produced identical greedy tours; that is all directional code variants produced identical tours, and all non-directional variants produced identical tours, the remaining testing focused on computational run-time comparisons. Each methodology was placed in the same arc-greedy heuristic shell so that testing would fairly compare the speed of the three subtour tracking and elimination methodologies. Bentley [4] and Wang [5] each utilized advanced computer techniques (k-d trees) and additional data structures to speed up the process of finding the next shortest arc available. However, since neither of these effect the speed of the subtour tracking and elimination methodologies they were not utilized.
Speed tests were conducted utilizing the R package "microbenchmark." 100 iterations of each code were run to create summary statistics on 13 different symmetric TSP instances and 9 asymmetric instances. Both symmetric and asymmetric instances were tested to determine if symmetry effected run time.

Symmetric Instance Results
Mean run times for a variety of symmetric TSP instances utilizing each of the methodologies can be seen in Table 2. When looking at the directional variants, the GT tends to be the fastest methodology on small instances followed by EL and MF. Once instance size reaches around 442, MF takes over as the fastest method for eliminating subtours. This largely is due to it's linear growth in operation count as instance size grows. For larger instances, the heuristic is conducting the same number of operations at each step. While the operations are slower for small instances, once the problem becomes larger it proves to be the most efficient.
As expected all Non-Directional subtour elimination variants run faster on symmetric instances versus their Directional counterparts. GT remains the fastest for instances up to around 442 nodes in size. EL is second fastest up till 343 nodes. At 442 nodes MF is fastest method. 2.

Asymmetric Instance Results
The directional variants of each subtour elimination codes were also run on Asymmetric TSP instances to compare runtimes to determine if any trends changed. The mean runtimes are in Table 3.
For the asymmetrical instances we note that MF became the fastest methodology when the instance size reached 358 nodes, but was beaten by the GT at an instance size of 403 nodes. This variation leads us to believe that some of the subtour methodologies may have computational advantages for specific TSP instances dependent on how the tour is constructed. For example, specific instances may have nodes spaced in such a fashion that for a majority of the tour construction the tracking methodology is maintaining a small number of large fragments. It is possible that one of the methodologies is computationally more efficient for these instances and may be computationally slower in TSP instances were nodes force an arc-greedy tour construction that results in many small fragments to be tracked. Future analysis of instance geometry and resulting runtimes should be conducted to test whether this hypothesis holds any merit. With the exception of the variation seen for the 403 node instance, prior overall trends from the Symmetric instances remain, where GT is competitive for small to medium sized asymmetric instances, but MF is fastest for larger instances.

Future Improvements/Analysis
The portion of the GT most susceptible to computational growth is the search to identify what tail of the path to move. If this search process growth can be limited, it is possible that the modified GT could outperform MF for larger instances as well. Some possible methodologies to limit computational growth include a better implementation of the row and column delete methodology in conjunction with a new row and column generation methodology. Size of the search operations could also be reduced drastically especially during the early iterations of the arc-greedy heuristic by only generating nodes and tails as needed.
The significance of separating the subtour tracking methodologies from the underlying greedy heuristic methodology is the ability to then use these subtour elimination methodologies to develop new fragment heuristics to build viable TSP tours. For example consider the constructive heuristic called the Ordered-Greedy (OG) heuristic. The OG heuristic is a node-greedy heuristic that takes as input a complete ordered list of nodes. Starting at the top of the list, each node is considered in turn and the available set of choices S j at each step is limited to the feasible arcs originating at that node. What differentiates the OG from NN, another node-greedy heuristic, is that multiple fragments may exist during the tour construction.
The motivation for the OG heuristic is to apply a more structured approach in which nodes are given priority in connecting to their nearest neighbors. Nodes higher in the list have maximum flexibility with minimal concern for node degree or subtours and thus typically choose better arcs than nodes later in the list which experience significantly less degrees of freedom in their legal arc choices. The quality of the solution found is thus heavily dependent on the order of the list.
To introduce the methodology of the OG heuristic, consider the following example. In this example an ordered list of D,E,C,B,A has been, through some unspecified fashion, predetermined. This ordering of the nodes list is reflected, for ease of reference, in the X matrix and To and From arrays on the left-side, and in the distance matrix on the rightside of Figure 13 whose rows are now sorted according to this list order. The constructive heuristic now makes greedy decisions starting at the top of this list and working down. The first greedy decision is made with respect to node D. The greediest, or shortest arc, from node D is arc D-E as indicated above. This arc and its associated node is tracked via the GT so the next decision can be made. The next decision is made with respect to node E. This is not due to node E being the head of the previous arc added, but rather because it is second in the provided ordered list: D,E,C,B,A. Looking at the row in the distance matrix associated with node E along with the GT output that captures ineligible moves (as seen in Figure 14), it can be seen that the shortest legal arc available is arc E-B. Figure 14 OG example -Iteration 2. OG picks greediest legal arc incident to node E This process continues row by row until the final row is reached which is where the To and From arrays are scanned to find the final legal arc as seen in Figure 15.
After adding arc A-D, the resulting tour becomes A-D-E-B-C-A which is also the optimal tour for this TSP instance. Pseudocode for the OG is found in Algorithm 4.

Conclusion
As an NP-hard combinatorial optimization problem, the TSP is often solved via heuristic methodologies. One of the biggest considerations when constructing solutions is avoiding  Nodes.Visited = Nodes.Visited+1 13: end while 14: Connect Hamilton Path subtours, or a loop of interconnected nodes that prevents a single continuous tour amongst all cities within the instance. This paper introduced a novel subtour elimination methodology for the arc-greedy heuristic that is compared to two known subtour elimination methodologies. Computational results were generated across multiple TSP instances for each method.
When utilizing an arc-greedy type heuristic, additional steps must be taken to ensure that subtours are avoided and resulting tour is a valid TSP solution. This paper recognized two accepted arc-greedy subtour elimination methodologies, the Exhaustive loop and Bentley's MF, and compared them to our Greedy Tracker. The comparison utilized both directional and non-directional variants of each code on 14 symmetric TSP instances and the directional variants on 9 asymmetric instances.
The results of the comparison between each of these arc-greedy subtour elimination methodologies showed that the GT was the fastest tracking methodology for small to medium sized instances. However, Bentley's MF maintains the computational advantage for larger instances.
However, these results also indicated that given a more efficient coding implementation of methodology used for the X Matrix, the GT could become the preferred methodology for all instance sizes. For future research, the GT should be modified to handle a new row/column generation and delete technique to minimize the computational time utilized in the searching portions of the GT.