Computer Science
Computer Science-2021.gz:143339
Computer Science
761197dd2e38636cc7ad57592812f410671e9d56
Faster Goal-Oriented Shortest Path Search for Bulk and Incremental Detailed Routing

We develop new algorithmic techniques for VLSI detailed routing. First, we improve the goal-oriented version of Dijkstra's algorithm to find shortest paths in huge incomplete grid graphs with edge costs depending on the direction and the layer, and possibly on rectangular regions. We devise estimates of the distance to the targets that offer better trade-offs between running time and quality than previously known methods, leading to an overall speed-up. Second, we combine the advantages of the two classical detailed routing approaches - global shortest path search and track assignment with local corrections - by treating input wires (such as the output of track assignment) as reservations that can be used at a discount by the respective net. We show how to implement this new approach efficiently.


Introduction
The task of VLSI routing [3,19] is to connect the set of pins of every net on a chip by wires so that wires of different nets are sufficiently far apart and various other constraints are met. Typically, one first computes a global routing, a rough packing of wires that ignores all local constraints but guarantees that the wires in certain areas do not require more space than available. Global routing allows for globally optimizing objectives such as power consumption and timing constraints [13,8].
The output of global routing then restricts the search space for every net in detailed routing, where many complicated rules need to be obeyed and one essentially routes one net at a time. While the detailed routing graph formed by routing tracks on an entire chip can contain about 10 13 vertices on 10-20 layers, the restricted area corresponding to the global routing solution for a net results in a much smaller detailed routing graph, with rarely more than 10 8 vertices. Nevertheless, these subgraphs are still huge, and there are millions of nets to connect. Two general strategies have been proposed (cf. [19]).
The first approach is based on a fast subroutine to find a shortest path that connects two metal components, each of which can consist of a pin or a set of previously computed wires connecting a subset of the pins of that net. The subgraph is given by the global routing solution, excluding vertices and edges that would result in a conflict to previously routed wires. To allow for an efficient packing of wires and to model various aspects such as signal delays, one uses different costs for horizontal and for vertical edges on each layer as well as for vias connecting two adjacent layers.
The second approach first considers the layers one after the other and assigns wires to routing tracks so that the most important detailed routing rules are satisfied, at least for most wires. This is often called track assignment [18,4]. Then detailed routing tries to correct violations locally. A very similar detailed routing problem occurs when a detailed routing has already been computed, but a few changes to the input have been made (for example corrections of the logical behavior or to speed up signals that arrived too late). In both cases, one asks for an incremental detailed routing, largely following the input but deviating where necessary. However, local corrections are often not possible if the routing is very dense.
One classical speed-up technique of Dijkstra's shortest path algorithm [5] is to use reduced costs, based on a feasible potential that estimates the distance to the targets [7,16,11]. Instead of the undirected graph with the original edge cost c(e), we orient each edge in both ways and run Dijkstra's algorithm with the reduced cost c π (e) := c(e) − π(v) + π(w) for every edge e directed from v to w, where the vertex potentials π are chosen so that c π is nonnegative and π(t) = 0 for every target t. These conditions imply that π(v) is a lower bound on the distance between v and the closest target. The better this lower bound is, the fewer vertices this goal-oriented version of Dijkstra's algorithm must label before it knows a shortest path to a target.
Hence, there is a trade-off between a possible preprocessing time, the query time to compute the potential of a vertex, and the quality of the lower bound. For example, in subgraphs of unweighted grid graphs, the ℓ 1 -distance to the nearest target can be a reasonable choice for π [9]. A better estimate, which however requires substantial preprocessing, was suggested by [14]. In this paper, we propose new methods with better trade-offs than were known before.
Moreover, we combine the advantages of the two classical detailed routing approaches mentioned above. Our new, more global approach treats given input wires as so-called reservations and encourages, but not forces, the detailed router to follow the reservations where feasible. This is achieved by finding a shortest path where reservations can be used at a discount.
However, this does not work well together with the classical goal-oriented techniques. For example, if there are some reservations (edges) that can be used at a 50 % discount, the ℓ 1 -distance would have to be divided by 2 in order to induce a feasible potential. This would often be a very inaccurate estimate, leading to an increased number of labels in Dijkstra's algorithm and hence larger running time. We show that our better potentials make goal-oriented Dijkstra not only as fast as without reservations, but in fact faster. Overall, this yields a new efficient incremental detailed routing algorithm.

Problem statement
We number the layers 1, . . . , l and let V = Z × Z × {1, . . . , l} and be the vertex set and edge set of an infinite grid with l layers. Edges connecting adjacent layers are called vias, edges in x-direction are horizontal and edges in y-direction vertical. We will consider finite subgraphs of G = (V, E). These subgraphs correspond to the area defined by the global routing solution and to the restriction to the routing tracks that can be used for the current net. Often, many vertices of these subgraphs will have degree 2 and will not be considered explicitly, but we ignore this here for the sake of a simpler exposition. Every layer has a preference direction (↔ or ); edges in the other direction are more expensive. Horizontal and vertical layers alternate. Moreover, the layers have very different electrical properties which is reflected by appropriate edge costs. In the simplest model, the cost of an edge depends only on its direction and the layer: let c ↔ z , c z > 0 for z ∈ {1, . . . , l} and c z,z+1 > 0 for z ∈ {1, . . . , l − 1}; then In a more general model, a rectilinear grid induces rectangular regions, called tiles, and the cost also depends on the tile. Let be integer coordinates that define the rectangular tiles and set la y e r 1 ξ1 ξ2 ξ3 υ1 υ2 υ3 la y e r 2 in this example. The cost of the green via edges on ξ 2 is min{c 11 1,2 , c 21 1,2 }, which is c 11 1,2 here. The edges that are not drawn have cost infinity because, in this example, they lie outside the area corresponding to the global routing solution. Now we have costs c ij↔ z , c ij z , c ij z,z+1 > 0 that also depend on the tile and define the edge costs accordingly. If an edge belongs to more than one tile, the minimum cost applies. See Figure 1 for an example. We allow that two (but not three) consecutive coordinates are identical, i.e., ξ i = ξ i+1 or υ j = υ j+1 , in order to model a cheap cost at one x-or y-coordinate only.
With this more general model, one can, for example, punish wires on low layers near the electrical source of a net (which would lead to poor delays) or implement a discount on reservations as we will describe in detail in Section 5.2. Moreover, we can set edge costs to infinity outside the area corresponding to the global routing solution so that the distances in (G, c) reflect necessary detours that are implied by routing in this subgraph.
Given a finite subgraph G ′ = (V ′ , E ′ ) of G and sets S, T ⊆ V ′ , we look for a shortest path from S to T in G ′ with respect to the length function c. For a goal-oriented path search, we define a potential π(v) for every vertex v ∈ V ′ by the distance to T in G (instead of G ′ ): The idea is that distances in G are much easier to compute than in the subgraph G ′ (we will see how fast), but often still give a good lower bound. The reason is that (G, c) has a simple structure, given by the tiles, la y e r 1 s1 s2 s3 la y e r 2 Figure 2: A target set T consisting of three target rectangles on two layers. The coarsest partition into tiles that is consistent with T is shown by the dotted lines. The figure also shows three possible query locations s 1 , s 2 , and s 3 , for which we might be interested in the distance to the closest target. The shortest s 1 -T -path, the shortest s 2 -T -path and the shortest s 3 -T -path for a cost function depending only on direction and layer are shown. In this example, the preference directions of layer 1 and 2 are ↔ and , respectively. Nevertheless it is cheapest for the s 3 -T -path to stay on layer 1 since its vertical segment is very short.
while G ′ can be very complicated as it does not contain vertices or edges whose use would result in a conflict to nets routed previously. This allows us to use Dijkstra's algorithm with the reduced costs c π in the digraph resulting from G ′ by orienting every edge in both ways, since the reduced costs are nonnegative. Some of our algorithms work best for simple targets. We often assume that T is represented as the union of t rectangles, where a rectangle is a vertex set of the form {(x, y, z) ∈ V | ξ − ≤ x ≤ ξ + , υ − ≤ y ≤ υ + , z = ζ} for some ξ − , ξ + , υ − , υ + ∈ Z and ζ ∈ {1, . . . , l}. Often t is small in practice. Sometimes it will be useful to assume that this representation is consistent with the partition of V into tiles in the following sense: each of the t rectangles representing T fits into the grid, i.e., is of the form This can be achieved by adding at most 2t new x-coordinates ξ i and at most 2t new y-coordinates υ i . We call this procedure refining the grid with respect to the targets. See Figure 2 for an example of the empty grid (i.e., p = q = 0) refined with respect to several target rectangles.

Previous work and our results
In the simple model without regions (i.e., for p = q = 0), one can query π easily in O(tl 2 ) time without preprocessing; see Proposition 2. We show that this can be reduced to O(tl); see Theorem 3. With a preprocessing time polynomial in t and l, we obtain a query time of O(log(t + l)); see Theorem 7. These results will be presented in Sections 2 and 3.
For the more general model, which is the subject of Section 4, Peyer et al. [14] refined the grid with respect to the targets and showed that then the restriction of π : V → R ≥0 to V ij z is the minimum of k 2

Model Preprocessing time
Query time Reference Table 1: Various methods to compute π(v), possibly after preprocessing. The running times depend on the number t of target rectangles, the number l of layers, and in the general model on the numbers p and q of coordinates that define the (p + 1)(q + 1) regions, and on the number k of different horizontal and vertical edge costs (cf. (1)). Note that k ≤ 2(p + 1)(q + 1)l.
affine functions for any i, j, z, where k is the number of different horizontal and vertical edge costs, i.e., They also showed that all these functions can be computed in O((p + t)(q + t)lk 4 log(p + q + l + t)) time, allowing O(k 2 ) time queries after this preprocessing (plus O(log(p + q + t)) to find the region containing the given vertex by binary search; here and henceforth p and q refer to the original number of rows and columns, before refining the grid). We make multiple improvements over the approach of Peyer et al. [14]. By considering domination between affine functions with different slopes, we reduce the number of affine functions which are needed to describe the minimum. By first computing the distances from the edges on the boundaries of the tiles to the targets, we can compute these affine functions faster. Finally, we use a regional query data structure to reduce query time. For any 0 < ε ≤ 1, we can obtain an algorithm with preprocessing time O((p + t)(q + t) min{k, (p + q + 1)l}l 1+ε 1 ε log(p + q + l + t)) and query time O(log(p + q + t) + 1 ε log(k + l)). See Table 1 for an overview.
Our second contribution is a new approach to incremental routing. Rather than trying to correct a given infeasible input routing with local transformations only, we compute a new routing from scratch, at least for all nets for which the input routing does not obey all rules. However, to compute a solution similar to the input where reasonable, we reserve the space occupied by legal input wires for the respective net and allow to use edges corresponding to input wires at a discount. By letting each input wire be a separate tile V ij z , we can model the discount in the cost function c and work with reduced costs efficiently. When most of the input routing can be used, we can find a shortest path much faster than without a discount.
This makes this new approach not only useful for incremental routing, but also for bulk routing. Treating the output of a track assignment as reservations (wherever it is legal) and then pursuing our new incremental routing approach can combine the advantages of the two classical bulk routing approaches, successive shortest paths and track assignment with local corrections. We explain our new approach in detail in Section 5, where we also show experimental results.

Distances without preprocessing in the simple model
In the simple model, there is always a shortest path with a very simple structure: Lemma 1. Let c : E → R >0 depend only on direction and layer, and let r, s ∈ V . Then there is a shortest path P between r and s in (G, c) that consists of at most one sequence of horizontal edges, at most one sequence of vertical edges, and hence at most three sequences of vias.
Proof. Let P be a shortest path, and let P [v,w] and P [v ′ ,w ′ ] be two maximal subpaths of P in the same direction (all-horizontal or all-vertical), say from v to w and from v ′ to w ′ , respectively, and let P [w,v ′ ] be the subpath in between. Suppose, without loss of generality, that P [v,w] and P [v ′ ,w ′ ] are horizontal paths and that the cost of an edge of P [v,w] is not more expensive than the cost of an edge of P [v ′ ,w ′ ] (note that these paths may be on different layers). Then translating P [v ′ ,w ′ ] by adding w − v ′ to all its vertices, translating P [w,v ′ ] by adding w ′ − v ′ to all its vertices, and swapping these two paths in P yields a walk from r to s with one maximal horizontal subpath less and at most the same number of maximal vertical subpaths. Moreover, the cost does not increase. If the walk is not a path, we can shortcut it to a path. By induction, the assertion follows.
Hence, in order to compute a shortest path, we can enumerate the layers on which the horizontal sequence and the vertical sequence are, and which of the two comes first: Proposition 2. Let c : E → R >0 depend only on direction and layer. Then, without preprocessing, one can compute dist (G,c) (s, T ) for any given s ∈ V and given T ⊆ V consisting of t rectangles in O(tl 2 ) time.
Proof. Enumerate over all t rectangles that T is composed of. For each such rectangle R, we can determine the vertex r ∈ R that is closest to s (geometrically) in constant time. Then, for each pair of layers z ↔ , z ∈ {1, . . . , l}, we consider two paths that connect r and s. The first one is composed of the path of vias that goes from r to layer z ↔ , followed by the horizontal path that goes to the x-coordinate of s, followed by the path of vias that goes to layer z , followed by the vertical path that goes to the y-coordinate of s, followed by the path of vias that goes to s. Some of these subpaths can be empty. The second one is constructed analogously, swapping the roles of r and s. By Lemma 1, one of these 2l 2 paths must be optimal.
We now show how to improve on this, obtaining a linear dependence on the number of layers: Proof. Again we enumerate over all t rectangles that T is composed of, and for each such rectangle R, we determine the vertex r ∈ R that is closest to s (geometrically) in constant time. First compute the total cost c z1,z2 of a path of vias between layer z 1 and layer z 2 for all z 1 , where z r and z s denote the layers of r and s, respectively. This can easily be done in O(l) time.
Now we compute the minimum length of a path from r to s that (when traversed from r to s) consists of a path of vias, then a horizontal path, then a path of vias, then a vertical path, then a path of vias. We will then do the same with exchanging the roles of r and s, and the smaller of the two is the distance between r and s by Lemma 1.
For each layer z ∈ {1, . . . , l}, consider the vertex v z on layer z whose y-coordinate is the one of r and whose x-coordinate is the one of s. We first compute the distanced z between r and v z in the subgraph of G that contains no horizontal edges on the layers 1, . . . , z − 1. By Lemma 1, a shortest path from r to v z consists of a path of vias, a single horizontal path, and another path of vias. Henced z is either the sum of c zr ,z (which we have precomputed) and c ↔ z times the difference of the x-coordinates of r and s, or d z+1 + c z,z+1 (if z < l), whichever is smaller. This implies thatd 1 , . . . ,d l can be computed in reverse order in total time O(l). Figure 3 for an illustration.
Next we compute the length of a path from r to s that consists of a shortest path from r to v z , a vertical path on layer z, and a path of vias. It is given by d z + c z,zs plus c z times the difference of the y-coordinates of r and s, and thus can now be computed in constant time. Taking the minimum over all layers z yields a shortest path from r to s by Lemma 1 and thus completes the proof.

Logarithmic query time in the simple model
In [1], the first author showed how to achieve O(log(t + l)) query time with polynomial preprocessing time, but did not care about the degree of the polynomial. We will show how to do the same with only O(t 2 l 3 log l) preprocessing time. We will need two ingredients before we can prove this. The first ingredient is an efficient algorithm for computing the intersection of three-dimensional half-spaces: Figure 3: Illustration of the algorithm described in the proof of Theorem 3. The vertices correspond to values computed during the algorithm, the edges correspond to the paths used to derive them. The two vertices in the same gray box are at the same geometric location v z = (x s , y r , z). The algorithm starts by computingd l and then derivesd z for z ∈ {1, . . . , l − 1} as the minimum ofd z+1 + c z,z+1 and the length of the path using layer z for the horizontal edges. Then it propagates the values from bottom to top to compute d 1 , . . . , d l . For each z ∈ {1, . . . , l} we combine the resulting r-v z -path with the paths indicated on the right to get an r-s-path. Finally, we take the cheapest of the l paths to obtain a shortest r-s-path.

Theorem 4 ([15]). The intersection of a set of n half-spaces in three-dimensional space can be computed in O(n log n) time. If the intersection is nonempty, it is a convex polyhedron that is presented as a minimal set of inequalities, the cycle of edges surrounding each face, and the coordinates of their endpoints.
The second ingredient is an algorithm for solving the planar point location problem, a well-studied problem in computational geometry. In the following, we regard any connected, closed part of a line in R 2 as a line segment and any region whose boundary consists of a finite number of line segments as a polygon. Theorem 5 ([12]). Let L be a finite set of line segments that intersect only at their endpoints and let (P i ) i∈I be the (open polygonal) connected components of R 2 \ L. Then there is a data structure that requires O(|L| log|L|) preprocessing time and, given any query point p ∈ R 2 , can then determine an index i ∈ I such that p lies in the closure of P i in O(log|L|) query time.
Point location algorithms which attain the same theoretic guarantees as [12], but successively improve practical performance and ease of implementation have been described in [10,6,17].
Theorems 4 and 5 can be combined in order to obtain a data structure for storing affine functions efficiently such that their pointwise minimum can be queried in logarithmic time: Then there is a data structure that requires O(|F | log|F |) preprocessing time and, given any query point p ∈ R, can then determine the value After projecting the result into the plane, we obtain a subdivision of the rectangle R into at most |F | convex polygons and a minimizing function f ∈ F for each polygon.
Second, we initialize a data structure to solve the point location problem within that subdivision. By Euler's formula, since each vertex (with possible exception of the four corners of R) in this subdivision is incident to at least three edges, the subdivision contains at most 3|F | + 1 line segments. Hence, this preprocessing can be implemented to run in O(|F | log|F |) time by Theorems 4 and 5.
When given a query location p ∈ R, we look up a polygon containing the query location in O(log|F |) time by Theorem 5, and evaluate the function attaining the minimum on that polygon in constant time. If the query point is on the boundary of multiple polygons, it suffices to evaluate the minimizing function of any one of these polygons.
We can now prove the main result of this section: Theorem 7. Let c : E → R >0 depend only on direction and layer, and let T ⊆ V consist of t rectangles. Then there is a data structure that requires O(t 2 l 3 log l) preprocessing time and, for any given s ∈ V , can Proof. We first interpret our instance as an instance of the general model by choosing p = q = 0, and then refine the grid with respect to the targets, which yields O(t 2 l) tiles V ij z . We compute an independent data structure for each tile.
By Lemma 1, a shortest path from any s ∈ V to any r ∈ T contains at most one sequence of horizontal and at most one sequence of vertical edges. For a fixed set of directions and layers of these sequences, for example north on layer z n and west on layer z w , there are at most two such paths, depending on whether we first go north and then west or vice versa, having possibly different via costs. Given a fixed tile V ij z , for every s = (x s , y s , z) ∈ V ij z , we are interested in the cost of a shortest such path to T , which is given by an affine function . To obtain the constant, one takes the minimum over all r = (x r , y r , z r ) ∈ T that lie northwest of V ij z , considering the sum of y r c z n − x r c ↔ z w and the via cost in the cheaper of the two cases regarding the order of the two directions. If there is no r ∈ T northwest of V ij z , no path of the required structure exists, and the constant can be considered to be ∞.
Similarly, an affine function for each of the 1 + 4l + 4l 2 combinations can be defined, and f (i,j,z) : (x s , y s ) → dist(s, T ) is the pointwise minimum of these. Here (−) means that we go to the target only by vias (or we are already at a target). We next show how to compute all O(t 2 l 3 ) affine functions in total time O(t 2 l 3 ). We describe this for one combination (nw, z n , z w ); it works analogously for the other combinations. Here we need to compute the constants c ij z := c (i,j,z,nw,z n ,z w ) for all tiles V ij z . For each layer z, we do this from northwest to southeast. We start with c ij z = ∞ whenever i = 0 or j = q and then set where c z1,z2,z3,z4 is the total cost of the vias to go from layer z 1 to layer z 2 to layer z 3 to layer z 4 ; see Figure 4 for an illustration. Since each of the t target rectangles shows up O(l 3 ) times (once for each combination and each z), all these O(t 2 l 3 ) affine functions can be computed in O(t 2 l 3 ) time.
Rather than storing just a list of these 1 + 4l + 4l 2 affine functions for each V ij z , we now build a data structure for each tile, using Lemma 6, in order to obtain a logarithmic query time. The total preprocessing time required to build these data structures is O(t 2 l 3 log l).
When given a query location s ∈ V , we find the correct tile V ij z and hence the correct data structure in O(log t) time by performing two binary searches. Then we look up the distance from s to T in O(log l) time by Lemma 6. If the query point is on the boundary of multiple tiles, it suffices to evaluate the minimizing function of any one of these tiles.

The general model
In this section, we develop an algorithm to compute the potential π(v) = dist (G,c) (v, T ) for any v in the general model efficiently after preprocessing. We will assume T to be consistent with the grid, i.e., we have already refined the grid if it was not. Our preprocessing will work on the horizontal and vertical line segments of the grid, i.e., the sets Hor ij The exposition will focus on the horizontal line segments; vertical segments can be handled analogously. Our algorithm consists of two preprocessing steps and a query step. The first preprocessing step is a variant of Dijkstra's algorithm. For its correctness, the following observation about the structure of shortest paths, similar to Lemma 1 for the simple model, is essential: Proof. By the same argument as in Lemma 1, there is a shortest path P such that every subpath of P that is entirely in the interior of one column {(x, y, z) ∈ V | ξ i ≤ x ≤ ξ i+1 } contains at most one horizontal sequence of edges and every subpath of P that is entirely in the interior of one row {(x, y, z) ∈ V | υ j ≤ y ≤ υ j+1 } contains at most one vertical sequence of edges. It is easy to see that such a path satisfies the claim.
Our algorithm will first compute dist (G,c) (s, T ) for all s lying in horizontal segments of the grid. By applying Lemma 8 inductively, we may assume that the corresponding shortest paths never use a horizontal edge or a via in the interior of a tile.
In the variant of Dijkstra's algorithm that we apply, we do not mark all vertices whose labels are guaranteed to be permanent explicitly, but the set of these vertices at some point in the algorithm is implicitly given: it consists of all vertices that have a smaller label than the one considered in the current iteration. In every iteration, several vertices might be added to this set and several updates of neighboring vertices are performed. This will be useful in order to decrease the required number of iterations significantly compared to the regular Dijkstra's algorithm, while the additional updates in every iteration can be performed with negligible overhead. This is possible because all horizontal edges along some horizontal segment have the same cost, and hence, the corresponding labels can be represented by affine functions that can be stored and updated very efficiently.
More precisely, for each horizontal segment Hor ij z , the algorithm maintains a set F ij z of affine functions f : [ξ i , ξ i+1 ] → R ≥0 such that each value f (x) corresponds to the length of a path between (x, υ j , z) and T . At any point during the algorithm, for every vertex (x, υ j , z), the minimum value min{f (x) | f ∈ F ij z } can be considered to be its current label. If ξ i = ξ i+1 , we simply store that value. Otherwise we keep only those functions which are not dominated, i.e., attain the pointwise minimum in more than one point. For storing the sets F ij z , we will use the following result. See Figure 5 for an illustration.
x f i (x), g(x) x Figure 5: Update of the data structure described in Lemma 9 when adding an affine function g to a set here f 3 and f 4 are deleted, x + f2 decreases, key(f 2 ) increases, and x − f5 increases.

Lemma 9
. Let x − , x + ∈ R be given with x − < x + . We can maintain a data structure that • stores a set F of affine functions g :

points where the function attains the minimum is an interval
• stores for every function g the interval [x − g , x + g ] and the minimum value key(g) := min{g(x − g ), g(x + g )} of g on this interval, and Moreover, each update increases the key of at most one function that remains in F , and never decreases a key.
Proof. We store the set of affine functions as a binary search tree in which they are sorted by their slopes. When a function g is added, we check whether g is dominated and otherwise insert it into the search tree in time O(log |F |). Then, starting from g, we iterate forward and backward, deleting all functions f that are identical to g or dominated by g on the whole interval [x − f , x + f ]. Each of these operations requires time O(log |F |). If a function f is partly dominated by g, then its interval and, if necessary, its key are updated. It can be required to shorten two intervals, left and right of the interval of g. However, only one key can increase due to the concavity of the pointwise minimum.
In addition to storing each set F ij z as specified in Lemma 9, we maintain a binary heap representing all functions in {F ij z | i ∈ {0, . . . , p}, j ∈ {1, . . . , q}, z ∈ {1, . . . , l}} that have not been processed yet, using the keys defined in Lemma 9. The functions that are added to or removed from some F ij z must be added to or removed from the heap at the same time, and whenever a key changes, it must be updated also in the heap.
The algorithm starts by initializing F ij z := ∅ for all i ∈ {0, . . . , p}, j ∈ {1, . . . , q}, and z ∈ {1, . . . , l}. If Hor ij z ⊆ T , we add the constant function x → 0 to the corresponding set F ij z . If not the whole segment, but one (or both) of its endpoints lies in T , we add the affine function describing the distance to this endpoint, i.e., x → min{c ij↔ z , c In every iteration, a function f with minimum value key(f ) is chosen and removed from the heap. The function f describes the labels of a subset of some horizontal segment Hor ij z , corresponding to the interval [x − f , x + f ]. We now propagate the labels from these vertices to the neighboring horizontal segments by computing at most six new affine functions: Using Lemma 10, we derive that g is not strictly dominated in x ′ . Hence, key(g) ≤ g(x ′ ) = dist (G,c) (s ′ , T ) < dist (G,c) (s, T ) < K. Since keys never decrease during the algorithm by Lemma 9, g was chosen and removed from the heap in a prior iteration. Now note that it is possible to choose i, j, i ′ , and j ′ such that s ∈ Hor ij z and the propagation of g adds a function f satisfying f (x) = c(P ) to F ij z . This concludes the first part of the inductive step.
For the second part of the inductive step, consider all i and j with s ∈ Hor ij z . Note that the possible choices of i and j are independent and consecutive. We know from the first part that there are some i 1 and j 1 such that there is a function f 1 ∈ F i1j1 z with f 1 (x) = dist (G,c) (s, T ) < K. Since this function is not dominated in x, we conclude key(f 1 ) < K. This means that this function has already been propagated and has added functions with the same value at x to all adjacent F ij z for which s ∈ Hor ij z . Successive application concludes the proof.
Note that both steps are necessary in the above proof. Figure 6 contains one such situation: we are looking for a function h ∈ F 12 1 which gives us the correct distance from the top right corner of the left tile to T , i.e., with h(4) = dist (G,c) ((4, 1, 1), T ) = 18. Since the vertical cost is lower in the right tile, the first step gives us the function g : x → 4 + x ∈ F 21 1 , which is propagated to f : x → 14 + x ∈ F 22 1 , which in turn is propagated to the desired function h : x → 26 − 2x ∈ F 12 1 . To ensure that the algorithm terminates and has the desired running time, we first show: Lemma 12. Let k ′ := min{k, (q + 1)l}. The number of slopes of functions ever added to F ij z is at most 2k ′ + 1. In particular, this also bounds the cardinality of F ij z at any stage.
Proof. By induction, the slope of every affine function added to F ij z during the algorithm is either zero or an element of Hence, if there were more than 2k ′ + 1 functions in F ij z , then there would be two different functions having the same slope. But then one of them is strictly dominated by the other one and would be removed from F ij z , contradicting the specification in Lemma 9.
In order to achieve the desired running time, we bound the number of iterations by (p+1)q(2k ′ +1)l. Each iteration chooses and removes a function f from the heap. Let (x, υ j , z) ∈ Hor ij z such that f (x) = key(f ). By Lemma 11, f (x) = dist((x, υ j , z), T ) and, by Lemma 10, each function g added to F ij z later on satisfies g(x) ≥ dist G,c ((x, υ j , z), T ). This means no other function with the same slope as f can ever be an element of F ij z in the future. Since the number of slopes of functions in F i ′ j ′ z ′ was bounded by 2k ′ + 1 in Lemma 12, we obtain the claimed bound on the number of iterations.
Finally, we need to implement each iteration in amortized time O(log(p + q + l)). Since each iteration generates at most six new functions and the number of iterations is bounded by (p + 1)q(2k ′ + 1)l, at most 2(p + 1)ql + 6(p + 1)q(2k ′ + 1)l functions are added in total, including the initialization. This also bounds the size of the binary heap to O(pq 2 l 2 ) such that each heap operation can be performed in time O(log(p + q + l)). By Lemma 9, each added function causes at most one increase-key operation. The total number of deletions is clearly bounded by the total number of insertions. Thus, the total time needed for the heap operations is given by O(pqk ′ l log(p + q + l)). The same holds for the time required for the updates of the sets F ij z by Lemma 9. This proves: Theorem 13. There is an algorithm that computes for each horizontal segment Hor ij z a set F ij z of at most , υ j , z), T ) for all (x, υ j , z) ∈ Hor ij z . The algorithm can be implemented to run in O(pqk ′ l log(p + q + l)) time, where again k ′ = min{k, (q + 1)l}.
We will now consider how to implement queries. We do so for all pairs of i and j independently. If we wanted to execute a query right after executing the algorithm described by Theorem 13 without any further preprocessing, the best we could do is O(log(p + q) + l log k) query time: we first compute i and j by binary search. In each of the 4l segments Hor ij z ′ , Ver , and Ver ij z ′ (z ′ ∈ {1, . . . , l}), we then compute the closest point r to the query location s. The distance from s to r can be computed in amortized constant time and the distance from r to T can be looked up in O(log k) time in the data structure storing the affine functions of that segment.
The above can be seen as an evaluation of the minimum of O(min{k, (p+q +1)l}l) affine functions. Hence it might be worthwhile building up the data structure described by Lemma 6 in an additional preprocessing step. This would speed up our queries to O(log(p + q + l)) time, however at the cost of an additional O(min{k, (p + q + 1)l}l log(p + q + l)) preprocessing time per tile.
The following result will obtain a trade-off between these two alternatives. By choosing the trade-off factor 0 < ε ≤ 1 to be a small constant, we obtain a query time of O(log(p + q + l)) after a preprocessing time which is arbitrarly close to O(pq min{k, (p + q + 1)l}l log(p + q + l)).

Theorem 14.
Let 0 < ε ≤ 1, let c : E → R >0 depend on tile and direction, and let T ⊆ V be consistent with the grid. Then there is a data structure that requires O(pq min{k, (p + q + 1)l}l 1+ε 1 ε log(p + q + l)) preprocessing time and, for any given s ∈ V , can then determine dist (G,c) (s, T ) in O(log(p + q) + 1 ε log(k + l)) query time.
Proof. We first apply Theorem 13 to compute dist (G,c) (v, T ) for all v in some Hor ij z and (analogously) for all v in some Ver ij z . Now our preprocessing will consider all combinations of i ∈ {0, . . . , p} and j ∈ {0, . . . , q} separately. Given a query location (x, y, z) ∈ V ij z , i and j can be determined in O(log(p + q)) time by binary search. Hence we fix i and j from now on.
We refer to the union of the up to 4l segments Hor ij z ′ , Ver Now suppose that a shortest s-T -path touches the boundary at least once. Consider a (without loss of generality, horizontal) segment H of the boundary which is touched first by one such path. By Lemma 8, we may pick our shortest path such that it starts with a sequence of vias to the layer of H, followed by a sequence of vertical edges to H. Hence, we can compute the cost of the path until it first touches H as an affine function in y, where s = (x, y, z). The rest of the path is a shortest path from some point in H to T , so we already computed its cost as a minimum of O(min{k, (q + 1)l}) affine functions in x, by Theorem 13. By iterating over each of these affine functions for every horizontal and vertical segment in the boundary, we can express the cost of a shortest s-T -path as a minimum of O(min{k, (p + q + 1)l}l) affine functions in x and y. Note that the only dependence of these functions on z is the cost of the initial via stack. We will exploit this now.
Instead of building a separate data structure for each layer, each involving all these affine functions, we distinguish between the cases whether the shortest path from s to T begins without vias, with vias up to a higher level, or with vias down to a lower layer. For the first case, we build up the point location data structure just as in Lemma 6, but each involving only the O(min{k, (p + q + 1)l}) affine functions on the boundary segments on that layer. We call these data structures D =z for z ∈ {1, . . . , l}.
For the other two cases, we build data structures D ↑  + 1 − a) min{k, (p + q + 1)l} log(p + q + l)) time and then allows for queries in O(log(k + l)) time.
The main advantage is that we can use the same data structure D ↑[a,b] for all layers z ∈ {1, . . . , a} because the cost of the via stack from layer z to layer a is a constant term that depends only on z. We can design the data structures so that each affine function on the boundary shows up only in at most 2l ε (instead of l) of these data structures.
To this end, we consider a balanced arborescence A whose leaves are the layers 1, . . . , l, such that A has maximum out-degree d = ⌈l ε ⌉ and depth ⌈ 1 ε ⌉, and, for every vertex v of A, the set of leaves reachable from v in A is a consecutive range L(v). Let parent(v) denote the parent of v (unless v is the root). For every vertex v except for the root, let L > (v) := {z ∈ L(parent(v)) | z > z ′ for all z ′ ∈ L(v)} and L < (v) := {z ∈ L(parent(v)) | z < z ′ for all z ′ ∈ L(v)}. Then we store D ↑L > (v) and D ↓L < (v) , unless this layer range is empty; see Figure 7. To answer a query for a point on layer z, we ask D =z and then traverse the path from z to the root in the arborescence A, and for each vertex v on that path (except for the root), we ask D ↑L > (v) and D ↓L < (v) . We have to query at most 2⌈ 1 ε ⌉ + 1 data structures, and each of these queries takes O(log(k + l)) time. To bound the preprocessing time, we see that each layer z appears in the layer range of only d data structures on each level of the arborescence (one per vertex whose parent's layer range contains z), and hence O(l ε 1 ε ) overall. Hence, the total preprocessing time is O(min{k, (p + q + 1)l}l 1+ε 1 ε log(p + q + l)). Since we do this for all i and j, the theorem follows.
Corollary 15. Let 0 < ε ≤ 1, let c : E → R >0 depend on tile and direction, and let T ⊆ V , not necessarily consistent with the grid. Then there is a data structure that requires O((p + t)(q + t) min{k, (p + q + 1)l}l 1+ε 1 ε log(p + q + l + t)) preprocessing time and, for any given s ∈ V , can then determine dist (G,c) (s, T ) in O(log(p + q + t) + 1 ε log(k + l)) query time. Proof. Refine the grid with respect to the targets, then apply Theorem 14 to ε and the refined instance. Note that this refinement does not increase the number of different costs.
We remark that the same ideas used in this section could be applied in the simple model to obtain a faster preprocessing time of O(t 2 l 2+ε 1 ε log l) at the cost of a slower query time of O(log t + 1 ε log l).

Implementation
With some modifications that we will describe below, we implemented the algorithms presented in the previous sections as part of BonnRoute [1,2], a detailed router developed at the University of Bonn in joint work with IBM. BonnRoute is the main detailed routing tool used by IBM for the design of its processor chips. Up to parallelization and conflict resolution, BonnRoute routes one net after the other. Each net is routed by iteratively connecting two of its components by a path until the net is fully connected, i.e., one component remains. The path search is the algorithmic core of BonnRoute and requires approximately 80-90 % of the total runtime.
Using shortest paths frequently leads to wirings that violate the design rules which need to be satisfied to ensure that the layout can be manufactured. Respecting even simple design rules is NP-hard [1]. BonnRoute uses a framework consisting of multiple components for avoiding violations. First, every computed path is  Table 2: Testbed consisting of nine real-world instances. Tech refers to the naming of the technology nodes by the foundry. The vertex set V of the detailed routing graph is the set of all locations which are on track with respect to at least one wire type on the current layer and would be on track when projected to at least one adjacent layer. Both (total length of) Wires and (total number of) Vias refer to the detailed routing computed by BonnRoute. Calls is the number of calls to our Dijkstra implementation. Since BonnRoute may try out different side constraints to find a path with as few design rule violations as possible, this does not match the number of computed paths. t and pql are the arithmetic mean over all calls to the preprocessing.
handed to a post-processing routine, which attempts to resolve violations locally. Second, we apply multilabeling, i.e., we search for shortest paths in a modified graph that can have multiple copies of each vertex (and different edges). The modifications are done in such a way that certain design rule violations are avoided. Finally, we impose restrictions to avoid violations at the start and end of a path. For further implementation details, see [1,2]. All experiments were performed on the same AMD EPYC 7601 machine with 64 CPUs and 1024 GB main memory using 64 threads. Table 2 gives an overview of our testbed. It consists of nine real-world instances from three recent IBM processor chips in 7 nm and 5 nm technology nodes. We started all experiments on the same instance from the same snapshot, which was taken right before the detailed routing. At this point, a (three-dimensional) global routing and possibly an allowed layer range were already computed for each net.
We use edge costs as they have been developed for many years in real design practice. They have three main components. The base cost does not depend on the net and models the amount of routing resources consumed by a path. Wiring against the preferred direction of a layer is ten times as expensive as wiring in the preferred direction. Apart from that, wires in x-and y-direction have the same base cost on all layers. The base cost of vias is chosen such that a via bridge, i.e., a path consisting of two vias on the same layer and a single segment of wiring in preferred direction between them, is cheaper than the direct connection between the two endpoints if and only if it blocks strictly fewer additional tracks. This means that the precise via costs depend heavily on the precise design rules and track patterns of the technology. On the highest layers, the base cost of vias can be more than ten times more expensive than on the lowest layers. This is because the thicker wires on the high layers require a track pattern with larger spacing. The second component of the cost function is an additive penalty, increasing the cost of wires outside of the assigned layer range. This is a heuristic approach to avoid timing failures due to wires on lower layers having more resistance. Finally, we restrict our path search to vertices which are inside of the area corresponding to the global routing solution. In the context of our general cost model, this can be seen as setting the cost to infinity outside of this area. Table 3 compares the performance of path searches with the original edge costs and with the reduced costs using three feasible potentials. Each of these three potentials is the distance to T in the same supergraph G of G ′ , but with respect to different edge costs c. For the simple and the general model, our implementation differs from the description in the previous sections as described below.
• Since t is usually small (3.3 on average), we iterate over all target rectangles in T for every query.
When computing the distance from the query location to one of the target rectangles in the simple model, we know the start and end layer. If we guess the lowest and highest layer used by a shortest  Table 3: Performance of the following four different feasible potentials on our testbed. In the rows without potential, each query returns 0 in constant time. When using ℓ 1 -distance, the distances in x-and y-direction are scaled by the minimal c ↔ z and c z , respectively, over all z. An O(l) preprocessing computes these two numbers and the total cost c 1,z of vias from layer 1 to each layer z; for the distance of two points on layers z < z ′ , we then use c 1,z ′ − c 1,z . The query returns the ℓ 1 -distance to T in O(t) time by iterating over the targets. In the simple and general rows, the shortest distance to T in the respective models is returned. Here the only difference between the two is that the general model restricts to the area corresponding to the global routing solution (outside of it, the costs are infinite). Implementation details for the last two potentials are described in the text. Runtimes are summed over all 64 threads except for the last column, which shows the total wall time of the overall BonnRoute run. path, the distance can be computed in constant time (Proposition 2). Due to the special structure of our cost function, at most l combinations need to be considered. In a preprocessing step, we compute those combinations that can be optimal for some query locations.
• The implementation of the general model uses a modified version of the algorithm used in Theorem 13.
Since the distance from a single tile to T can be expressed as a minimum of very few affine functions (on every instance, the average is below 1.2), it is more efficient to compute these functions instead of the distance from the horizontal and vertical segments to T . On the other hand, instead of Lemma 9, we need to use a more complicated data structure to maintain the set of non-dominated functions. We store the convex polygon of points on which each function attains the minimum. This way, insertion can be implemented to run in O(|F | log |F |) time. During each query, we find the correct tile in O(log(p + q + t)) time using binary search and evaluate all non-dominated functions of that tile on the query location to compute the minimum.
The results show that the general potential performs significantly better than the simple potential, which already performs much better than the ℓ 1 -distance potential. Both the number of labels and the runtime improve on every instance, even when considering the additional preprocessing time. The relative improvement differs a lot between different instances. Most of this difference can be explained by some situations in which we get only minor improvements by our potentials: • If no path is found, all reachable vertices in the graph are labeled. None of the potentials show any improvement on these instances. In fact, the path searches without potential are the fastest since they do not need any query time, with a total of 285 hours (summed over all instances and all 64 threads).
With the three potentials, these path searches take a total of 329, 337, and 331 hours, respectively.
• After a path search failed, BonnRoute may perform a backup path search which allows routing through existing wires at high cost (and then would rip-up such wires and try to re-route them). Since these rip-up costs are not modeled in any of our potentials, a large portion of the graph may be labeled regardless of which potential is used. On such instances, the order of the potentials regarding their performance is the same as when looking at all instances, but the relative improvements are much smaller.
The column Standard Dijkstra calls in Table 3 excludes these situations and hence shows an even larger gain than the column All Dijkstra calls. The question how to model rip-up costs efficiently when computing potentials remains for future research.

Reservations and discounts for incremental routing
In chip design practice, there are two main scenarios where a detailed routing is not computed from scratch, using just a global routing as input, but in an incremental way, using an approximate detailed routing as input. The first scenario is when a detailed routing has already been computed, but now a few changes have been made, for example in order to correct the logical function of the chip or to improve its timing behavior. The second scenario is when a step in between global and detailed routing is used, typically called track assignment, that maps the global wires to routing tracks in a way that obeys most -but not all -design rules.
In both scenarios of incremental routing, we get an almost feasible detailed routing as part of the input, and the task is to compute a completely feasible detailed routing by doing only few changes. While it is not exactly specified what "few" means, the motivation is that the input routing has already been optimized, for example with respect to the timing behavior of the chip; moreover, one aims at saving runtime.
The traditional approach to incremental routing is to check for violations of design rules (e.g., wires of different nets overlapping) and to try to repair such violations locally, in a relatively small area around that violation. While this can be parallelized very well, many violations cannot be repaired locally, and then the overall approach may fail or resort to global path searches as backup. Moreover, if the wiring of a net needs to be repaired in multiple places, the final result can be quite bad, for example with too many detours to meet timing constraints.  : Example of re-routing a net using reservations after changes to the input have been made. The green and red wires connected pins p 1 , p 2 , and the old position of p 3 . Now suppose p 3 has been moved and is no longer connected. Moreover, pin q from a different net has been moved and now makes the red piece of wire illegal. Now we want to connect p 1 , p 2 , and the new position of p 3 and use much of the old wiring. We convert all of the old wiring (green and red) to global wires and add a global wire (not shown) that connects to the new position of p 3 . Next we create reservations. Even though all of the green wires are legal, we may choose to create reservations only for the thick green wires, e.g., if we expect the harm of blocking other nets to outweigh the benefit of keeping them usable for this net. When this net is being routed, we may end up discarding the dashed part of the reservations, using the rest of the reservations and adding the blue wires. The solid wires then constitute the new routing for this net.
We suggest to repair violations globally but with a preference of using the initial solution. To this end, we convert any detailed wire in the input to a global wire and possibly a reservation. A reservation reserves that space for the particular net. When other nets are routed earlier, this space is blocked. Therefore, reservations are created only for (parts of) detailed wires that do not conflict with other detailed wires in the input. For an example situation, see Figure 8.
Once a net is routed for which we have created reservations, we would like to encourage, but not force, the net to use the reserved space. We do this by defining a discount factor 0 < δ < 1 and multiplying all edge costs on a reservation of that net by δ. There are two reasons for using this incentive to route nets similarly as in the input: first, during detailed routing, we cannot do a complete timing analysis (this would be too slow), and the input routing has already been analyzed carefully. Second, we would hope for a speedup if the reservation serves as a useful guide how to route a net.
However, this speedup does not come automatically. In fact, with the traditional goal-oriented search techniques, reservations would lead to a slow-down. For example, if we define the potential π to be the ℓ 1 -distance to the nearest target, multiplied by the minimum edge weight in that direction, then we have to multiply it by δ if there are any reservations in that direction (no matter how useful). Our generalized framework, however, allows us to refine the grid not only with respect to the targets, but also with respect to the reservations, and define individual (discounted) costs on the edges corresponding to reservations.
If we use the output of track assignment as input, this often is a good solution on higher layers, where  we have mostly longer wires, but much less so on lower layers, which are primarily used for pin access.
(This is because long wires on low layers have a high resistance and thus poor delay.) In this case we may define reservations only on high layers and let the pin access and the short wires be freely determined by the detailed router.
To evaluate the effect of reservations on incremental routing, we compare two different algorithms. Both of them replace the same subset of the input wires by global wires and compute a new solution. What subset is chosen depends on the scenario and is described below. The algorithm no reservations is our standard bulk routing algorithm, starting from scratch using these global wires without any additional information. The algorithm reservations creates reservations for all input wires that are legal, except for short dangling wires. Both algorithms use the distance in the general model as their potential. In the algorithm reservations, the cost of reservations is multiplied by a discount factor of δ = 3 4 (i.e., 25% discount). The general model takes this discount into account.
We compare these two algorithms on 15 instances belonging to two different scenarios. In the first scenario, we start from an input in which all nets are connected and almost all wires are legal. More precisely, our input is a snapshot taken right after those bulk routing runs in Table 3 which used the general potential. We do not keep any part of our old solutions fixed, but replace all the detailed wires in nets we connected by global wires. The results of this experiment can be seen in Table 4.
The second scenario in which we evaluate the effect of reservations is a detailed routing which is no longer legal due to timing optimization. Table 5 gives an overview of this part of the testbed, consisting of six snapshots taken in a production flow just before incremental detailed routing. Unlike in the previous scenario, we keep pin-to-pin paths fixed if they consist only of detailed wires. Any wire not in such a path will be replaced by a global wire. See Figure 8 for an example. The performance of both algorithms on these instances is shown in Table 6.