On matchings, T -joins, and arc routing in road networks

Matchings and T - joins are fundamental and much-studied concepts in graph theory and combinatorial optimization. One important application of matchings and T -joins is in the computation of strong lower bounds for arc routing problems (ARPs). An ARP is a special kind of vehicle routing problem, in which the demands are located along edges or arcs, rather than at nodes. We point out that the literature on applying matchings and T -joins to ARPs does not fully exploit the structure of real-life road networks. We propose some ways to exploit this structure. Computational results show significant running time improvements, without deteriorating the quality of the lower bounds.


INTRODUCTION
Matchings are a fundamental concept in graph theory and combinatorial optimization, with a wide array of applications (see, e.g., [20,26,30]). Given an undirected graph G = (V, E), with |V| even, a perfect matching is a set of edges that meets each vertex exactly once. If we are also given a weight w e for each e ∈ E, the minimum-weight perfect matching problem (WPM for short) calls for a perfect matching of minimum total weight. Edmonds [13] showed that this problem can be solved in polynomial time, and a variety of efficient algorithms are now available (e.g., [12,19,26]).
A closely related concept is that of a T-join. Given a graph G as before, and a set T ⊆ V with |T| even, a T-join is a set of edges that meets each vertex in T an odd number of times, and each other vertex an even number of times. Given edge weights as before, the minimum-weight T-join problem (WTJ) calls for a T-join of minimum total weight. Figure 1 illustrates the concept of a T-join. On the left is a graph G, with nodes in T represented by hollow circles, and weights on the edges. On the right is the same graph, with the MWJ solution indicated with thick lines. (In the online version of this article, the hollow circles and thick lines are in red. ) Edmonds and Johnson [14] showed that any WTJ instance can be transformed into a WPM instance, and thereby solved efficiently. (For details, and alternative algorithms, see Section 2.1.) They also pointed out that T-joins can be used to solve a problem that they called the Chinese postman problem (CPP), in honor of Mei-Gu Guan [23].
In the CPP, G represents a road network, and the weight of an edge represents the time taken to traverse the corresponding road. A postman wishes to traverse every road at least once, as quickly as possible, starting and finishing at the same node. Any CPP instance can be reduced to a WTJ instance by setting T to the set of vertices that are incident on an odd number of edges. The edges in the optimal T-join then represent roads that need to be traversed twice.
The CPP is an example of an arc routing problem (ARP). An ARP is a special kind of vehicle routing problem, in which the demands are located along edges or arcs, rather than at nodes (e.g., [10]). Whereas the CPP can be solved in polynomial time, most ARPs of interest are  -hard. Matchings have been used to compute useful lower bounds for such ARPs (e.g., [2, 3, 7, 22, 27, 32-35, 37, 38]). Graph G with nodes in T indicated In a recent project with real-life instances, we encountered large-scale ARPs, with over ten thousand edges. For these particular ARPs, matching techniques gave acceptable lower bounds, but used an excessive amount of both time and memory. The purpose of this article is to show that one can dramatically reduce the amount of computational effort needed by matching techniques, by exploiting the structure of real-life road networks. For ease of presentation, we focus on the CPP and another prominent ARP, the so-called capacitated arc routing problem or CARP [22]. This article has a simple structure. The key literature is reviewed in Section 2. Some observations concerning road networks are given in Section 3. In Section 4, we show how to solve the CPP more quickly. In Section 5, we show how to compute lower bounds for the CARP more quickly. In Section 6, we make some remarks about possible extensions of our technique to other ARPs.
Throughout this article, all graphs are undirected, simple, and loopless, unless otherwise specified. We let n and m denote |V| and |E|, respectively. We also let deg(i) denote the degree of node i in G. We also assume that, in any WTJ instance, all weights are non-negative.

LITERATURE REVIEW
We now briefly review the relevant literature. Section 2.1 covers algorithms for WPM and WTJ. Section 2.2 deals with the application of matchings and T-joins to ARPs. Section 2.3 deals with planar graphs.

Matchings and T-joins
For surveys of WPM algorithms, see [12,19,20,26,30]. At present, the fastest strongly polynomial algorithm is that of Gabow [18,19], which runs in O (n(m + n log n)) time. Among the algorithms that are only weakly polynomial, we mention the one of Duan et al. [12]. It assumes that all weights are integers, and it runs in O Software for matching has lagged behind the theory to some extent. At present, the most effective available routines for WPM are those of Mehlhorn and Schäfer [31] and Kolmogorov [25]. They perform very well in practice, but run in O(nm log n) and O(n 2 m) time, respectively.
As for WTJ, Edmonds and Johnson [14] showed that it can be reduced to WPM as follows. First, compute shortest paths between all pairs of nodes in T. Then construct a complete undirected "auxiliary" graph, say G + , with T as its node set. Set the weight of each edge in G + to the length of the corresponding shortest path in G. Solve WPM in G + , and then replace each edge in the matching with the corresponding shortest path in G. We call this approach "Ed-Jo." Figure 2 shows how Ed-Jo works for the WTJ instance on the left of Figure 1. We remark that Dijkstra's single-source shortest path algorithm can be implemented to run in O(m + n log n) time [17], and therefore the shortest-path phase in Ed-Jo takes O (|T| (m + n log n)) time. The matching phase takes O(|T| 3 ) time, since G + is complete.
Korte and Vygen [26] presented an alternative approach to WTJ, which we will call "Ko-Vy." Like Ed-Jo, Ko-Vy involves the solution of WPM in an auxiliary graph. However, the auxiliary graph is different. Consider a given node v ∈ V. If deg(v) is odd and v ∈ T, or deg(v) is even and v ∉ T, then v is replaced with a clique on deg(v) nodes. We call the new nodes "clones".
If v does not satisfy the stated condition, it is replaced by a clique on deg(v) + 1 nodes. We call the additional node a "parity correction" node. The new edges are given a weight of zero. Figure 3 shows how Ko-Vy works for the WTJ instance in Figure 1. The auxiliary graph is shown on the left, with parity correction nodes represented as green diamonds. An optimal WPM solution is shown on the right.
Ko-Vy is very easy to implement, with no need for shortest-path calculations. On the other hand, the auxiliary graph has O(m) nodes and O(nm) edges. As a result, constructing the auxiliary graph takes O(nm) time and solving the WPM takes O(n 2 m 2 ) time. More efficient reductions from WTJ to WPM, which use auxiliary graphs with only O(m) nodes and edges, appear in [4,8]. These reductions enable one to solve WTJ in only O ( m 2 log n ) time. They are however trickier to implement. At the time of writing, the fastest known WTJ algorithm is that of Gabow [19], which runs in O (|T|(m + n log n)) time. It is based on a conversion of WTJ into a capacitated b-matching problem. We omit details for brevity.

Applications to arc routing
We mentioned above that the CPP can be converted to WPM or WTJ. The authors have also used WPM and/or WTJ to compute lower bounds for  -hard ARPs (e.g., [2,3,27,33,34]). For brevity, we focus here on the CARP [22], which has received the most attention. In the CARP, we are given a graph G = (V, E) and a set E R ⊆ E of required edges. Node 0 represents the depot. For each e ∈ E, we are given a positive cost c e . For each e ∈ E R , we are given a positive demand q e . An unlimited fleet of identical vehicles, each of capacity Q, is located at the depot. The task is to find a minimum-cost set of routes, each starting and ending at the depot, such that each required edge is serviced on exactly one route, and the total demand on each route does not exceed Q.
Many lower bounds for the CARP have been proposed which use WPM as a sub-routine [3,7,22,32,35,37,38]. For brevity, we review only one in detail: the bound called "LB1" in [7]. To describe it, we need some more notation. We let V R and V O denote the set of nodes in G that are incident on at least one required edge, and on an odd number of required edges, respectively. We also let K denote ⌈ ∑ e∈E R q e ∕Q ⌉ , which is a lower bound on the number of routes. Finally, for notational simplicity, we assume that no required edges are incident on the depot. (If this is not the case, then we can make it so by adding a dummy node to G, along with a dummy edge of zero cost.) LB1 works as follows. First, compute shortest paths from the depot to each node in V R . Let v 1 be the closest node to the depot, v 2 be the second closest, and so on. Let r i be the number of required edges incident on node v i . Let s be the smallest integer such that r 1 + · · · + r s ≥ 2K, and let V c denote {v 1 , … , v s }. (The subscript c is to remind us that the nodes in V c are "close" to the depot.) Let A consist of 2K "dummy" nodes, representing copies of the depot. Let B contain r i "clones" of vertex v i , for i = 2, … , s. Let S contain all nodes in V O \ V c , except the depot. Construct a complete graph, say G, with node set A ∪ B ∪ S. The weight of each edge in G is set to the cost of the shortest path between the corresponding nodes in G, with one exception: edges between dummy nodes are given infinite weight. Finally, solve the WPM in G. Then LB1 is the weight of the WPM solution plus the cost of the required edges.
The above-mentioned procedure is illustrated in Figure 4. On the left, we see a graph G, with required and non-required edges represented by red lines and gray dashed lines, respectively. For each edge e ∈ E, the cost c e and demand q e are indicated. We suppose that Q = 4 and that the depot is node v 0 . Note that K = 2 and s = 3 for this instance. On the right, we show the auxiliary graph G, with an optimal WPM solution indicated by thick red lines. The weight of the matching is 37, which yields a lower bound of 37 + 18 = 55.
Among the other matching-based bounds for the CARP, we mention only the node duplication lower bound (NDLB) from [35]. The procedure is very similar to the one for LB1, but the auxiliary graph is much larger than G, with up to 2(|E R |+K) nodes. Ahr [1] showed that NDLB is slightly stronger than LB1, though this comes at the cost of a significantly increased running time.

Planar graphs
Now we recall some facts concerned with planar graphs. The first is Euler's theorem, which states that, in a planar graph, m ≤ 3n − 6. This implies that, when G is planar, one can solve WPM in O ( n 2 log n ) time. Lipton and Tarjan [29] showed that, in fact, it is possible to solve planar WPM in only O(n 3∕2 log n) time. Their algorithm exploits a famous result in their earlier paper [28], which states that, in any planar graph, there exists a set S ⊂ V such that (a) and (b) if S is removed from G, each connected component in the resulting graph contains no more than 2n/3 nodes. A suitable S, called a separator, can be found in linear time [28]. Henzinger et al. [24] gave an algorithm for single-source shortest-paths in planar graphs, which runs in only O(n) time. This implies that the shortest-path phase of Ed-Jo algorithm can be conducted in only O(n|T|) time in the planar case. The matching phase, however, still takes O(|T| 3 ) time.
It is also worth noting that the reduction from WTJ to WPM in Barahona [4] preserves planarity. Together with the Lipton-Tarjan result, this implies that planar WTJ can be solved in O(n 3∕2 log n) time. A more direct algorithm, with the same running time, is given in Barahona [5].
To close this section, we mention that road networks are often planar and, even if not, they invariably contain small separators. Fast algorithms for finding small separators in road networks, together with encouraging computational results, are given in [11] and [36].

ON ROAD NETWORKS
In this section, we make some observations regarding real-life road networks. Our first observation is that road networks are not just sparse; they also have bounded degree. Although this is obvious intuitively, we decided to compute the degree distribution for a sample of 12 cities from across the world, using the Python package OSMnx [9]. For each city, we considered the closest 1000 nodes to a central landmark, and computed their degrees. The results are shown in Table 1.
As one might expect, the vast majority of nodes have degree 3 or 4. We expected the maximum degree to be 5, but we see instead that it is 7. On the other hand, only 3 nodes (out of 12 000) have degree 7. We remark that the average degree ranged from 2.71 (Hanoi) to 3.71 (New York).
The reason that this is relevant is that, in Ko-Vy, a node of degree d in G becomes a clique of size d in the auxiliary graph (or d + 1 if a parity correction node is needed). If G has bounded degree, the size of each clique in the auxiliary graph will be O(1), and the auxiliary graph will have only O(n) nodes and edges. Not only that, but Ko-Vy is much easier to implement than either Ed-Jo or the methods in [4] and [8]. Thus, for road networks, Ko-Vy may be an attractive alternative to those methods. We will see in the next section that this is indeed the case.
Next, we considered the issue of planarity. For each city, we computed the graph induced by the given 1000 nodes. We then tested planarity using a package called Python Planarity 1 . We were surprised to find that none of the graphs were planar. Closer inspection revealed that all of them could be made planar by deleting a very small number of edges (corresponding to over-or under-passes).
Finally, we decided to test the claim made in [11] and [36] that road networks tend to have small separators. We found that, indeed, all 12 graphs contained separators of size O ( √ n ) . Moreover, such separators can be found easily. One way is as follows. Let G = (V, E) be the graph in question, and let G − be a planar subgraph obtained by deleting some edges corresponding to over-or under-passes. Compute a separator S ⊂ V in G − using the linear-time algorithm in [28]. By definition, removing S from V causes G − to become disconnected. Let F be the set of edges in E that have end-nodes in different connected components. If F = ∅, we are done. Otherwise, use a greedy algorithm to construct a minimal set of nodes T ⊂ V \ S that covers the edges in F. By definition, S ∪ T is a separator in G.
We conjecture that the planar WTJ algorithm of Barahona [5] could be modified, using the small separators mentioned, to solve the CPP in road networks in O(n 3∕2 log n) time.

THE CHINESE POSTMAN PROBLEM IN ROAD NETWORKS
Armed with the facts mentioned in the previous section, we can now analyze the theoretical running times of various approaches to the CPP in road networks:

3.
The approaches in [4] and [8]. These also yield auxiliary graphs with O(n) nodes and edges. Thus, they take the same time (asymptotically) as Ko-Vy.
The above analysis suggests that Ko-Vy should be faster than Ed-Jo when applied to large-scale CPP instances on road networks. To test this, we created 18 CPP instances, with cities selected from Paris, London and Moscow, and with n selected from 1000, 2000, 5000, 10 000, 20 000, and 50 000. The nodes were selected as in the previous subsection, and all edges having both end-nodes in the given set of nodes were put into the graph. Note that none of the resulting 18 graphs were planar. The cost of each edge was set to the length of the corresponding road, rounded to the nearest meter. The data for each instance, along with the optimal solution values, is made available at the Lancaster University Data Repository 2 .
We remark that some of the graphs contained loops (edges that connect a node with itself). Most of these loops represented small areas, near entrances of hotels, where taxis can change direction. We removed them, for simplicity of implementation.
For Ed-Jo, we implemented our own version of Dijkstra's single-source shortest path algorithm. We used a binary heap, since it is much easier to code than a Fibonacci heap, yet its running time is O(n log n) for sparse graphs. To solve the WPM instances, we used the open-source software package Blossom V [25]. Although it runs in O(n 2 m) time in theory, it performs extremely well in practice, as we will see below. For all experiments, we used a Lenovo ThinkPad laptop with an i5-8250U processor, running under Windows 10 at 1.6 GHz with 16 GB of RAM. Table 2 shows the results obtained when transforming each of the 18 CPP instances into WPM instances. For each instance, we show the name of the city, the number of nodes (n), and the number of edges (m). Then, for each approach, we show the number of nodes (ñ) and edges (m) in the resulting auxiliary graph. It is clear that Ed-Jo leads to auxiliary graphs with far fewer nodes than Ko-Vy. On the other hand, the number of edges is drastically bigger. Table 3 shows, for each instance and each approach, the time taken to construct the auxiliary graph (T1) and the time taken to solve the WPM instance (T2). All times are reported in seconds. A dash indicates that Blossom V had to be aborted due to memory limitations.
It is clear that the traditional approach, Ed-Jo, is to slow to be used for large-scale instances. In particular, the shortest-path phase is time-consuming. In principle, one could obtain a speed-up by computing the shortest paths in parallel. Nevertheless, the matching phase of Ed-Jo is itself rather slow. Moreover, we found that the matching phase consumed a great deal of memory for these instances.
We were surprised to find that the other approach, Ko-Vy, is orders of magnitude faster. Indeed, it takes less than one second for every instance. This is no doubt due to the fact that road networks have bounded degree.
For interest, we also implemented a third reduction from WTJ to WPM, due to Barahona [4]. The resulting auxiliary graphs were a bit larger than the ones that we obtained with Ko-Vy. As a result, the WPM phase took about twice as long to solve. We omit details for brevity. The main conclusion from this section is that, for large-scale CPP instances on road networks, Ko-Vy is preferable to Ed-Jo.

THE CAPACITATED ARC ROUTING PROBLEM IN ROAD NETWORKS
We now move on to consider the CARP. We start by analyzing the time taken to compute LB1 in the case of road networks. The procedure begins by solving |V O | shortest-path problems in G. Since G is sparse, this takes O (|V O | n log n) time. The next step  We now present a procedure that yields lower bounds of comparable quality to those of LB1, but which exploits the properties of road networks. The procedure incorporates concepts from both LB1 and Ko-Vy.
The first steps in our algorithm are identical to that of LB1. That is, we compute K (the lower bound on the number of vehicles), the shortest paths in G from the depot to each node in V R , the values r 1 , … , r s , and s (the smallest integer such that r 1 + · · · + r s ≥ 2K). We then let V c denote {v 1 , … , v s }. The next step is different. We take G and apply the Korte-Vygen procedure with T set to V O , the set of nodes that are incident on an odd number of required edges. The resulting graph will be called G + . Figure 5 shows the graph G + for the CARP instance on the left of Figure 4. We remind the reader that K = 2 and s = 3 for this example.
Graph G + for the CARP instance in Figure 4 As before, we call the copies of the original nodes "clones," and any additional nodes "parity correction" nodes. If a clone is incident on a required edge, it will be called an "R-clone." Note that, if a node in G is incident on t required edges, then t of its clones will be R-clones. In particular, there are r i R-clones of node v i . Note also that, since G has bounded degree, G + has only O(n) nodes and edges.
For a given node v ∈ V, we will call the corresponding set of nodes in G + a "gadget," and denote it by g(v). In Figure 5, each gadget is enclosed in a circle. We also let g ′ (v) denote the corresponding set of nodes inG, the auxiliary graph that is used for

the same parity if and only if there is no parity correction node in g(v).)
We now create an even bigger graph, called G ++ . We begin by taking G + and adding 2K extra "dummy" nodes, representing copies of the depot. We will call the dummy nodes D 1 , D 2 , … , D 2K . We then add several sets of edges: • We connect D 1 to the first R-clone of v 1 , D 2 to the second R-clone of v 1 , and so on. In this way, we connect each dummy node to its own R-clone. The weight of each additional edge is the length of the shortest path from the depot to the corresponding member of V c . • For each i ∈ V O \ V c , we pick one R-clone as a "representative" of i. We then add an edge between each dummy node and each representative. The weight of each edge is the length of the shortest path from the depot to the given node i.
If p > 0, we have p R-clones of v s that are not connected to any dummy nodes. We pick one of those as a "representative," and connect each dummy node to it. The weight of each additional edge is the length of the shortest path from the depot to v s . Figure 6 shows the graph G ++ for the same CARP instance as before. Since K = 2, there are four dummy nodes. We have s = 3 and V c = {v 1 , v 2 , v 3 }. Thus, two dummy nodes are connected to R-clones of v 1 , one is connected to the R-clone of v 2 , and the other is connected to an R-clone of v 3 . We also have V O \ V c = {v 5 }. Thus, all dummy nodes are connected to a representative of v 5 . Finally, we have p = 5 − 4 = 1. Since p is positive and s = 3, all dummy nodes are connected to a representative of v 3 (namely, the unique R-clone of v 3 that was not already connected to a dummy node).
Note that, when G is sparse and has bounded degree, G ++ has only O(n) nodes and O(n + K|V O |) edges. We will show that the solution to the WPM in G ++ yields a valid lower bound for the CARP. We will need the following lemma. Proof. Let M be a perfect matching inG that matches k dummy nodes with R-clones of v s , where k > p ′ . We will show that we can construct a perfect matching M ′ , of no larger cost, such that k − 2 dummy nodes are matched with R-clones of v s . Let j 1 and j 2 be two R-clones of v s that are currently matched with dummy nodes. We assume w.l.o.g. that those dummy nodes are D 1 and D 2 . Also let I denote the set of R-clones of nodes in V c \ {v s } that are currently not matched with dummy nodes, and note that |I| ≥ k ≥ 2. We consider two cases: Case 1: there exist two R-clones in I, say i 1 , i 2 , that were matched in M. We obtain a matching of no larger cost by deleting the edges {D 1 , j 1 }, {D 2 , j 2 } and {i 1 , i 2 }, and adding the edges {j 1 , j 2 }, {D 1 , i 1 } and {D 2 , i 2 } (see Figure 7). To see why, note that the nodes in V c represented by i 1 and i 2 are no further from the depot than v s .
Graph G ++ for the same CARP instance Proof.
Recall that |g(v) | − |g ′ (v) |≥0 and even for every v ∈ V. Let M be a perfect matching inG. We construct a perfect matching M ′ in G ++ , as follows: • We first match the dummy nodes. If k R-clones of a given node v are matched with dummy nodes inG, we match k corresponding R-clones of v with dummy nodes in G ++ . We can apply this procedure until all dummy nodes are matched since, by Lemma 1, there are no more than 2K − ∑ s−1 i=1 r i R-clones of v s matched with dummy nodes inG. • We then take the other edges in M, if any, and construct the corresponding edges in M ′ . For a matched edge {i, j} inG, we match in G ++ pairs of clones of the nodes in the shortest path between the corresponding nodes in G. We can repeat this until all matched edges inG are considered since, for any node v ∈ V, the gadget g(v) contains at least deg(v) clones. • Finally, we check to see if there are any unmatched nodes in G ++ . Consider a node v ∈ V. By construction, the number of unmatched copies of v in G ++ , if any, must be even. Thus, if such copies of v exist, we can match pairs of them using edges in G ++ of zero cost.
By construction, M and M ′ have the same cost. ▪ Figure 8 illustrates Theorem 1 for our CARP example. The graph G ++ is shown once more, but with an optimal WPM solution indicated by thick red lines. The weight of the matching is 37, so the resulting lower bound is 37 + 18 = 55. We remark that LB1 has the same value for this example.
To test the new bounding procedure, we converted each of our 18 CPP instances into CARP instances, as follows. The depot is located near the center. Edges incident on the depot were made non-required. Each other edge was made required with probability 1/2. The demand of each required edge was set to a random integer between 1 and 10. To make the fleet size Graph G ++ and an optimal matching (K) realistic, we set the vehicle capacity (Q) in such a way that K is |E R |/500, rounded up to the nearest integer. The cost of traversing each edge is set to be the length of the corresponding road, in meters, rounded up to the nearest integer. The data for each instance is made at the Lancaster University Data Repository as well.
As before, we used Blossom V to solve the matching problems. For the largest instances, with 50 000 nodes, our PC ran into memory difficulties. Accordingly, we used a workstation with Intel Xeon E5-2640 v3 processor at 2.6 GHz and 32 GB of RAM for those instances. Table 4 gives information concerning the size of the auxiliary graphs. It is clear that our transformation leads to auxiliary graphs with more nodes, but the number of edges is dramatically smaller. Table 5 shows, for each instance and each approach, the time taken to construct the auxiliary graph (T1), the time taken to solve the WPM instance (T2), and the resulting lower bound. One can see that, for these instances, our procedure is much faster than the one in [7]. Moreover, our lower bound matches LB1 in 16 out of 18 cases, and is only slightly weaker in the remaining two cases.
For interest, we also implemented NDLB from [35]. The resulting auxiliary graphs were much larger than the ones from either LB1 or our procedure, and they took much longer to compute. The lower bounds were slightly better than LB1, by around 5% on average, but the WPM phase took around ten times longer than it did for LB1. We also experienced even more problems with memory. We omit details for brevity. We close this section with a remark. In our approach, we use 2K dummy nodes. It would be much more efficient to use just one dummy node, and require it to have degree 2K. One would then have to solve a minimum-weight f -factor problem in G ++ , instead of a WPM. Unfortunately, at the time of writing, no efficient open-source f -factor code was available.

CONCLUDING REMARKS
We have shown that, by exploiting the special structure of road networks, we can solve large-scale CPP instances in less than a second, and compute strong lower bounds for large-scale CARP instances within a few seconds. Our procedures are orders of magnitude faster than previous approaches in the literature. Moreover, any future improvements in matching software will make our procedures even better.
We believe that our approach could be fairly easily adapted to other ARPs, such as ARPs with multiple depots [15], time deadlines, [16,27] and intermediate facilities [21], or ARPs on mixed graphs [6].