Art of Unlocking  ††footnotetext: The work of the first author (O. O.) was
supported by the grant PhyMath ANR-19-CE40-00. Part of the work of the second
author (S.S.) has been carried out at Skoltech and at IITP RAS. The support of
Russian Science Foundation 20-41-09009 is gratefully acknowledged.

Oleg Ogievetsky1,2,3 and Senya Shlosman1,3,4

1Aix Marseille Université, Université de Toulon, CNRS,

CPT UMR 7332, 13288, Marseille, France

2I.E.Tamm Department of Theoretical Physics, Lebedev Physical Institute,
Leninsky prospect 53, 119991, Moscow, Russia

3Inst. of the Information Transmission Problems, RAS, Moscow, Russia

4Skolkovo Institute of Science and Technology, Moscow, Russia

###### Abstract

We describe our recent results concerning the rigidity/unlockability
properties of clusters of rigid bodies sliding over the unit sphere.

##  1 Introduction: Balls

The most known question in the area of movable/rigid configurations deals with
the FCC and HCP configurations of 12 unit balls touching the central unit ball
$\mathbb{B}$ (shown in red):

Figure 1: FCC and HCP configurations

The question is: can one roll the 12 balls over $\mathbb{B}$ in such a way
that they cease to touch (or kiss) each other? (The number of kissing points
of the outer spheres is 24 for both configurations.) When the answer is
positive, we say that the configuration can be unlocked.

The fact that the configuration FCC can be unlocked is well known. One of
unlocking moves is described in [CS]. The situation with HCP cluster is less
known. In [T] it is stated that the configuration HCP is rigid. But in fact it
is not the case, and the unlocking moves for both FCC and HCP are presented in
[KKLS].

We present quite simple moves which unlock FCC and HCP.

### FCC configuration.

The balls of the FCC configuration are centered at the vertices of the
cuboctahedron.

During the move the top three balls and the bottom three balls remain fixed.
Consider the three “triangles of balls” sharing one ball with the top
triangle. (One of them – the triangle $BCD$ – is visible in Fig. 2.) Now
rotate each of these three triangles with the same velocities around their top
vertices, e.g. the triple $C,B,D$ is rolled over the central ball, as a solid,
with the ball $C$ fixed (Fig. 2).

Figure 2: FCC unlocking

To see that the balls do not hit each other during the move, consider the ball
$B$, which kisses four balls $A$, $C$, $D$ and $E$. The triple $BCD$ moves as
a solid, so $B$ keeps kissing $C$ and $D$, and no conflict happens here. The
ball $B$, being rotated around $C$, goes below the equator, while $A$, being
rotated around $F$, goes above the equator, with the same equatorial
projection of speeds, so the distance between them increases. The balls $C,B$
and $E$ lie on the same great circle, and both $C$ and $E$ do not move, so the
distance $BE$ also increases.

These considerations suffice due to the rotation symmetry (by $2\pi/3$ angle
around the vertical axis) of the cuboctahedron.

### HCP configuration.

The balls of the HCP configuration are centered at the vertices of the Johnson
solid called the triangular orthobicupola. They can be seen as three ‘rhombic’
configurations, one of them being $ABCD$ on the picture:

Figure 3: HCP unlocking

Here $A,B,C,D$ are the centers of the corresponding four balls. The move
consists in rolling each of three rhombi around their centers, with equal
velocities. That is, one applies to the four balls of each rhombus the
rotation around an axis connecting the center of the rhombus (the mid of the
segment $BD$ for the rhombus $ABCD$) and the origin.

To see the absence of conflicts here, we observe that the balls $E$ and $C$
move downward, with the same horizontal projections of velocities, so their
mutual distance increases (since their common distance to the north pole
increases). The ball $B$ goes up, while the ball $F$ goes down, so they avoid
each other. The rest again follows from the same rotational $\mathbb{Z}_{3}$
symmetry of the orthobicupola.

In both cases, FCC and HCP, after unlocking, the twelve balls can be moved
apart and be positioned at the vertices of the regular icosahedron. Then some
free space appears between the unit balls; their radii can be blown up to the
value $r=\left(\sqrt{\frac{5+\sqrt{5}}{2}}-1\right)^{-1}\approx 1.10851$, when
they finally start to kiss, making 30 kissing points. Yet this extra space is
not large enough to incorporate the 13-th unit ball, as was shown in [SW],
settling the famous discussion started by Newton and Gregory.

##  2 Cylinders

Our interest in the cylinder case started from the question of Kuperberg, [K]:
one can easily put six unit cylinders around a central unit ball $B:$

Figure 4: Configuration $C_{6}$

The question is whether one can arrange seven such cylinders around $B$ in a
non-intersecting way. The question looks as an insult to intuition, but to
this day we do not have a rigorous answer to it. However one cannot arrange
eight unit cylinders [BW]).

The configuration shown looks quite solid; though it is not rigid –

Figure 5: Non-rigidity of $C_{6}$

– it looks pretty tight. So it came as a big surprise that there are
configurations of six unit cylinders around the unit ball, where the cylinders
do not touch each other. Such configurations were found by M. Firsching, [F].
In his example he was able to position six cylinders of the radius
$r=1.049659$ around the unit ball. This example was obtained by a numerical
exploration of the corresponding 18-dimensional configuration manifold.

In our paper [OS1] we found a way to unlock the configuration $C_{6}$ in a
symmetric (with respect to the dihedral group $\mathbb{D}_{3}$) manner, which
enabled us to improve the Firsching’s $1.049659$ to the value

$r_{\mathfrak{m}}=\frac{1}{8}\left(3+\sqrt{33}\right)\approx 1.093070331...$
(1)

The corresponding $\mathbb{D}_{3}$-symmetric configuration $C_{\mathfrak{m}}$
is shown below:

Figure 6: The initial configuration $C_{6}$ and the record configuration
$C_{\mathfrak{m}}$

The unlocking move looks as follows. Assume that the cylinders on the left
picture of Fig. 6 (the configuration $C_{6}$) point at the North. We first
describe a three-dimensional family of moves. The triple of pink cylinders go
upward by $\varphi,$ ‘horizontally’ by $-\varkappa,$ and the three vectors
$\uparrow$ are rotated (around the axis joining the origin and the tangent
point, counterclockwise if viewed from its tip) by $\delta,$ while the three
remaining cylinders go downward by $\varphi$, ‘horizontally’ by $\varkappa,$
and the three vectors $\uparrow$ are rotated, in the same way as the upper
ones, by $\delta$. Now, optimising $\varkappa$ and $\delta$ for each
$\varphi$, we obtain a curve $\gamma$ in the moduli space of six cylinders.
For each point on the trajectory $\gamma$ the configuration is
$\mathbb{D}_{3}$-symmetric: its symmetries are 120o rotation around the axis
North-South and 180o rotation around a perpendicular axis. The common radius
of the cylinders grows as the configuration moves along $\gamma$ up to a
certain point which corresponds to the configuration $C_{\mathfrak{m}}$ and
then decreases. It turns out that the trajectory $\gamma$ is an algebraic
curve (given by the relations (25)-(27) in [OS1]), and this is the reason for
the exact equality $\left(\ref{01}\right)$. Y. Kallus made a movie featuring
our path $\gamma$, [Ka].

It is interesting to note that all the angles describing the configuration
$C_{\mathfrak{m}}$ are pure geodetic, in the sense of [CRS]: an angle $\alpha$
is pure geodetic if the square of its sine is rational.

We believe that our configuration $C_{\mathfrak{m}}$ is the record
configuration, and in any configuration the radii of six equal non-
intersecting cylinders tangent to the unit ball $B$ are $\leq
r_{\mathfrak{m}}.$ But we have no proof of this maximality. We can prove,
however, a local version of this statement: if $C$ is a small proper
perturbation of the configuration $C_{\mathfrak{m}},$ then the radii $r$ of
its cylinders satisfy: $r<r_{\mathfrak{m}}.$ Here proper means that $C$ is not
a rotation of $C_{\mathfrak{m}},$ and small means that the six tangent points
and six directions of the (equal) cylinders comprising $C$ are close to those
of $C_{\mathfrak{m}}.$ The proof of this statement is a content of our paper
[OS2]. The proof is quite involved, since the function $r$ on our
18-dimensional configuration manifold is not smooth. It turns out that the
differentials of the distances between the tangent (to the unit central
sphere) generatrices of the cylinders obey a convex linear dependence
$\Lambda$ at the point $C_{\mathfrak{m}}$ on the curve $\gamma$. Let $E$ be
the linear subspace of the tangent space on which all differentials vanish.
Then we have proved that a sufficient condition for a local maximum is as
follows: the same linear combination $\Lambda$ of the second differentials is
negatively defined on $E$. We were able to check that this sufficient
condition holds at the point $C_{\mathfrak{m}}$ [OS1].

Remark. Let $F$ denote the function $\max(F_{1},\ldots,F_{k})$ where
$F_{1},\ldots,F_{k}$ are analytic functions on a vector space, and
$F_{1}(x_{0})=\ldots=F_{k}(x_{0})$ at some point $x_{0}$. The condition “the
restriction of $F$ to any straight line passing through the point $x_{0}$ has
a local maximum at $x_{0}$” is not sufficient to guarantee that the function
$F$ has a local maximum at the point $x_{0}$, see example in [OS2], section 5.
However we believe that if the restriction of $F$ to any analytic curve
passing through the point $x_{0}$ has a local maximum at $x_{0}$ then the
function $F$ does have a local maximum at the point $x_{0}$.

Without the analyticity assumption on $F_{i}$-s the above conjecture does not
hold. We present an example of just one $C^{\infty}$ function $\Phi$ whose
restriction to any analytic path passing through
$\left(0,0\right)\in\mathbb{R}^{2}$ has a local minimum at the origin, while
the origin is not a local minimum of $\Phi$. Consider the function

$\psi\left(x\right)=\left\\{\begin{tabular}[c]{ll}$\exp\left(-\frac{1}{x}\right)$&if
$x>0$\\\ $0$&if $x\leq 0$\end{tabular}\ \ \right.,$

and let $\eta\left(x\right)$ be any $C^{\infty}$ function on $\mathbb{R}^{1}$
with support on the segment $\left[-\frac{1}{2},\frac{1}{2}\right],$ which is
i) even, ii) decaying on $[0,\infty),$ iii) $\eta\left(0\right)=1.$ We define
the $C^{\infty}$ function $\Phi$ on $\mathbb{R}^{2}$ by

$\Phi\left(x,y\right)=\left\\{\begin{tabular}[c]{ll}$\exp\left(-\frac{1}{x^{2}}\right)\eta\left(\frac{y-\psi\left(x\right)}{\psi\left(x\right)}\right)$&for
$x>0$\\\ $0$&for $x\leq 0$\end{tabular}\ \ \right..$

The function $\Phi$ is indeed of $C^{\infty}$ class. The only problematic
point is the origin $\left(0,0\right)$. Any partial derivative (as
$x\rightarrow 0^{+},y\rightarrow 0$) of $\Phi$ belongs to the vector space
(invariant with respect to the partial derivatives) of finite linear
combinations of functions
$e^{-x^{-2}+a_{1}x^{-1}}y^{a_{2}}x^{-a_{3}}\eta^{(a_{4})}\left(ye^{x^{-1}}-1\right)$
for some $a_{1},a_{2},a_{3},a_{4}\in\mathbb{Z}_{\geq 0}$ and hence vanishes.

Let us show that for any analytic path
$\gamma:\left[0,1\right]\rightarrow\mathbb{R}^{2}$,
$\gamma\left(0\right)=\left(0,0\right),$ the function
$\Phi\left(\gamma\left(t\right)\right)$ vanishes on some segment $0\leq
t<u_{\gamma},$ so $0$ is a local minimum of the function
$\Phi\left(\gamma\left(t\right)\right)$ on $\left[0,1\right],$ while the point
$\left(0,0\right)\in\mathbb{R}^{2}$ is evidently not a local minimum of the
function $\Phi.$ To see this, note that the support of the function $\Phi$
lies inside the ‘beak’ $\beta=\left\\{\left(x,y\right):x\geq
0,\frac{1}{2}\psi\left(x\right)\leq y\leq
2\psi\left(x\right)\right\\}\subset\mathbb{R}^{2}.$ Let
$\gamma\left(t\right)=\left(x\left(t\right),y\left(t\right)\right)$ be an
analytic path, defined by two analytic functions
$x\left(t\right),y\left(t\right).$ Then
$x\left(t\right)=a^{\prime}t^{k^{\prime}}+O\left(t^{k^{\prime}+1}\right),$
$y\left(t\right)=a^{\prime\prime}t^{k^{{}^{\prime\prime}}}+O\left(t^{k^{{}^{\prime\prime}}+1}\right)$
for some real $a^{\prime},a^{\prime\prime}$ and integer
$k^{\prime},k^{\prime\prime},$ and therefore the path $\gamma$ in the vicinity
of the origin is a graph of a function
$y_{\gamma}\left(x\right)=bx^{k^{{}^{\prime\prime}}/k^{\prime}}+o\left(x^{k^{{}^{\prime\prime}}/k^{\prime}}\right)$
for some real $b$. (Without loss of generality we can assume that
$k^{\prime}>0.)$ Therefore there exists a value $\varepsilon_{\gamma}>0$ such
that the graph of the function $y_{\gamma}$ does not intersect the beak
$\beta$ for $0<x<\varepsilon_{\gamma},$ because the function
$\exp\left(-\frac{1}{x}\right)$ is smaller than any power of $x$ in the
appropriate vicinity of $0.$ Therefore the function $\Phi$ vanishes on that
piece of $\gamma.$

Apparently, to probe some $C^{\infty}$ function, one needs all $C^{\infty}$
paths, and not just analytic paths.

The local maximality of the configuration $C_{\mathfrak{m}}$ implies that it
is rigid, i.e. cannot be unlocked. It seems to play the role of the
icosahedral configuration of 12 kissing balls above, while the six unit
cylinders can be rolled away from each other and create some free space
between them on the sphere. Yet whether this space is sufficient for the
seventh unit cylinder to be squeezed in ( $\equiv$ Kuperberg question) is
unknown, in contrast to the 13 unit balls problem. The thesis of O. Yardimci,
[Y], contains the theorem, proven together with A. Bezdek, that one cannot put
7 cylinders of the Firsching radius $1.049659$ in contact with the unit
central ball, in a non-intersecting way.

In addition to configurations on Figs. 4 and 5, W. Kuperberg pointed out yet
another configuration of six unit cylinders around the unit ball:

Figure 7: Configuration $O_{6}$

It looks rigid, and in our paper [OS3] we have proven this rigidity: for every
small proper perturbation of the configuration $O_{6}$, the radii $r$ of its
cylinders satisfy: $r<1.$ In particular, the configuration $O_{6}$ also cannot
be unlocked. The proof turned out to be even more involved. The differentials
of the distances between the tangent (to the unit central sphere) generatrices
of the cylinders obey three convex linear dependencies $\Lambda_{1}$,
$\Lambda_{2}$ and $\Lambda_{3}$. Let again $E$ be the linear subspace of the
tangent space on which all differentials vanish and denote by
$q_{1},q_{2},q_{3}$ the same linear combinations $\Lambda_{1}$, $\Lambda_{2}$
and $\Lambda_{3}$ of the second derivatives. We have checked that the system
$q_{1}(x)>0,q_{2}(x)>0,q_{3}(x)>0$ of inequalities has no solution on $E$ and
proved that this is a sufficient condition for a local maximum [OS3, OS5].
Interestingly, no convex linear combination of the forms $q_{1},q_{2},q_{3}$
is negatively defined on $E$ so we could not just apply the technics like
Sylvester criterion. We think that it is a beautiful program – to develop the
theory (examples, criteria, classification) of tuples $q_{1},\ldots,q_{j}$ of
quadratic forms on a vector space for which the system
$q_{1}(x)>0,\ldots,q_{j}(x)>0$ of inequalities has no solution.

Sometimes this configuration $O_{6}$ is called ‘octahedral’ because, probably,
the tangency points lie at the vertices of the regular octahedron. We present
an interpretation of the configuration $O_{6}$ which relates it to the
configuration of rotated edges of the regular tetrahedron. This
interpretation, which shows that it rather deserves to be called the
tetrahedral configuration is as follows. Consider the configuration of the
tangent to the unit sphere lines which are continuations of the edges of the
regular tetrahedron. The points of the sphere at which tangent lines pass are
the edge middles of the tetrahedron:

Figure 8: $\delta$-process

Then each edge is rotated around the diameter of the unit sphere, passing
through the middle of the edge, by an angle $\delta$. On our picture the point
$A$ (in green) is the middle of the edge $UV$. The line, passing through the
point A and rotated by the angle $\delta$, is shown in red. The lines passing
through other five middles of edges are rotated by the same angle $\delta,$ in
accordance with the group $\mathbb{A}_{4}$ of the proper symmetries of the
tetrahedron. This motion we call $\delta$-process.

Let is replace each rotated line by a cylinder of radius $r,$ tangent to the
sphere, in such a way that the line is its generatrix, and (some of) the
cylinders are kissing, which uniquely defines $r$ as a function
$r\left(\delta\right)$ of the angle $\delta$. Then for $\delta=0,$ $\pi/2$ the
radius $r=0,$ while it is maximal at $\delta=\pi/4,$ with
$r\left(\pi/4\right)=1.$ This is precisely the configuration $O_{6}.$

Following this interpretation we have introduced in [OS4] configurations of
tangent cylinders for two remaining pairs of dual Platonic bodies, that is,
for the octahedron/cube $\left(\mathcal{O}/\mathcal{C}\right)$ and
icosahedron/dodecahedron $\left(\mathcal{I}/\mathcal{D}\right)$; for the pair
tetrahedron/tetrahedron this is precisely the configuration $O_{6}$. The
number of cylinders in our configurations is equal to the number of edges of
either of Platonic bodies in the pair, that is, twelve for the pair
octahedron/cube and thirty for the pair icosahedron/dodecahedron. The
corresponding radii $r_{\mathcal{O}/\mathcal{C}},\
r_{\mathcal{I}/\mathcal{D}}$ of the cylinders are obtained as a result of the
maximisation of the functions
$r_{\mathcal{O}/\mathcal{C}}\left(\delta\right),\
r_{\mathcal{I}/\mathcal{D}}\left(\delta\right)$ over the rotation angle
$\delta.$

Pair $\mathcal{O}/\mathcal{C}$. For the pair $\mathcal{O}/\mathcal{C}$ the
optimal configuration is the following:

Figure 9: $\mathcal{O}/\mathcal{C}$ maximal configuration of cylinders, view
from a vertex of the cube

Each cylinder touches four other cylinders. The corresponding value
$\delta_{\mathcal{O}/\mathcal{C}}$ satisfies
$\tan(\delta_{\mathcal{O}/\mathcal{C}})=\frac{3^{1/4}}{\sqrt{2}},$
(approximately, $\delta_{\mathcal{O}/\mathcal{C}}\simeq 0.23856\pi\simeq
0.74946$). The corresponding radius of touching cylinders is

$r_{\mathcal{O}/\mathcal{C}}=\frac{\sqrt{3}-1}{1+2\sqrt{2}-\sqrt{3}}\approx
0.3492\ .$

Pair $\mathcal{I}/\mathcal{D}$. For the pair $\mathcal{I}/\mathcal{D}$ the
optimal configuration is the following:

Figure 10: Maximal $\mathcal{I}/\mathcal{D}$ configuration, view from the tip
of a 3-fold axis

Here each cylinder touches eight other cylinders. The optimal angle
$\delta_{\mathcal{I}/\mathcal{D}}$ is given by

$\delta_{\mathcal{I}/\mathcal{D}}=\arctan\left(\sqrt{t_{0}}\right)\ ,$

where $t_{0}\simeq 0.694356$ and is a root of the polynomial

$5t^{6}-80t^{5}+190t^{3}-4t^{2}-84t+9\ .$

Approximately,

$\delta_{\mathcal{I}/\mathcal{D}}\simeq 0.694707\ .$

The corresponding radius of cylinders is approximately

$r_{\mathcal{I}/\mathcal{D}}\simeq 0.115558\ .$

Another image of this configuration:

Figure 11: Maximal $\mathcal{I}/\mathcal{D}$ configuration, view from the tip
of a 5-fold axis

is shown on the cover.

We conjecture that both record configurations – $\mathcal{O}/\mathcal{C}$ and
$\mathcal{I}/\mathcal{D}$ – are rigid. This probably can be checked by using
our machinery developed in [OS2] and [OS3], but the computations are too
lengthy, and we have not attempted to perform them.

Another interesting phenomena happens during the $\delta$-process at the angle
values when the function $r_{\mathcal{I}/\mathcal{D}}\left(\delta\right)$
vanishes. That means that some of the 30 lines are intersecting. There are
three values of $\delta$ at which this happens. The patterns emerging are: ten
linked triangles, five linked tetrahedra and six linked pentagonal stars.

Figure 12: $\mathcal{I}/\mathcal{D}$ minima

## References

  * [BW] P. Brass and C. Wenk, _On the number of cylinders touching ball_ , Geom. Dedicata 81 (2000), 281–284. DOI 10.1023/A:1005290809501
  * [CRS] J. H. Conway, C. Radin and L. Sadun, _On angles whose squared trigonometric functions are rational_ , Discrete & Computational Geometry 22 (3) (1999), 321-332. DOI 10.1007/PL00009463
  * [CS] Conway, J.H. and Sloane, N.J.A., 2013. _Sphere packings, lattices and groups_ (Vol. 290). Springer Science & Business Media.
  * [T] L. Fejes-Tóth, _Lagerungen in der Ebene auf der Kugel und im Raum_ , Springer-Verlag, 2nd ed. (1972).
  * [F] M. Firsching, _Optimization Methods in Discrete Geometry_ , Berlin (2016).
  * [Ka] Y. Kallus page, https://ykallus.github.io/demo/cyl2.html
  * [K] W. Kuperberg, _How many unit cylinders can touch a unit ball?_ Problem 3.3, in: DIMACS Workshop on Polytopes and Convex Sets, Rutgers University (1990).
  * [KKLS] R. Kusner, W. Kusner, J. C. Lagarias, and S. Shlosman, _Configuration Spaces of Equal Spheres Touching a Given Sphere: The Twelve Spheres Problem_ , in: Ambrus G., Barany I., Boroczky K., Fejes Toth G., Pach J. (eds) New Trends in Intuitive Geometry, Bolyai Society Mathematical Studies, vol 27. Springer, Berlin, Heidelberg (2018); arXiv:1611.10297
  * [OS1] O. Ogievetsky and S. Shlosman, _The Six Cylinders Problem:_ $\mathbb{D}_{3}$_-symmetry Approach_ , Discrete Comput. Geom. 65 (2021) no. 2, 385–404. DOI 10.1007/s00454-019-00064-3
  * [OS2] O. Ogievetsky and S. Shlosman, _Extremal Cylinder Configurations I: Configuration_ $C_{\mathfrak{m}}$; Discrete Comput. Geom. 66 (2021) no. 1, 140–164. DOI 10.1007/s00454-020-00244-6
  * [OS3] O. Ogievetsky and S. Shlosman, _Extremal Cylinder Configurations II: Configuration_ $O_{6}$; Experimental Mathematics (2019). DOI 10.1080/10586458.2019.1641768
  * [OS4] O. Ogievetsky and S. Shlosman, _Platonic compounds of cylinders_ , in: Integrability, Quantization, and Geometry, Proceedings of Symposia in Pure Mathematics, vol. 103.2, Amer. Math. Soc., Providence, RI, 2021, 447–480. DOI 10.1090/pspum/103.2/01864. arXiv:1904.02043
  * [OS5] O. Ogievetsky and S. Shlosman, _Critical configurations of solid bodies and the Morse theory of MIN functions_ ; (Russian). Uspekhi Mat. Nauk, 74 (2019) no. 4, 59–86. DOI 10.4213/rm9899 Translation in Russian Mathematical Surveys 74 (2019) no. 4, 631–657. DOI 10.1070/RM9899
  * [SW] K. Schutte and B. L. van der Waerden, Das Problem der dreizehn Kugeln, Math. Ann. 125 (1953), 325-334.
  * [Y] O. Yardimci, _On The Number of Cylinders Touching a Sphere_ , PhD thesis, Auburn University, 2019.


[Bill]bred [Keld_Stefan]kblue

# Local elimination in the traveling salesman problem

William Cook
Combinatorics and Optimization
University of Waterloo
Supported by a Humboldt Research Award held at the Research Institute for
Discrete Mathematics, University of Bonn. Keld Helsgaun
Department of People and Technology
Roskilde University  Stefan Hougardy
Research Institute for Discrete Mathematics
University of Bonn  Rasmus T. Schroeder
DB Systel GmbH, Frankfurt am Main

###### Abstract

Hougardy and Schroeder (WG 2014) proposed a combinatorial technique for
pruning the search space in the traveling salesman problem, establishing that,
for a given instance, certain edges cannot be present in any optimal tour. We
describe an implementation of their technique, employing an exact TSP solver
to locate $k$-opt moves in the elimination process. In our computational
study, we combine LP reduced-cost elimination together with the new
combinatorial algorithm. We report results on a set of geometric instances,
with the number of points $n$ ranging from 3,038 up to 115,475. The test set
includes all TSPLIB instances having at least 3,000 points, together with 250
randomly generated instances, each with 10,000 points, and three currently
unsolved instances having 100,000 or more points. In all but two of the test
instances, the complete-graph edge sets were reduced to under 3$n$ edges. For
the three large unsolved instances, repeated runs of the elimination process
reduced the graphs to under 2.5$n$ edges.

##

The classic Dantzig-Fulkerson-Johnson paper [5] on the traveling salesman
problem (TSP) introduces the cutting-plane method as an effective means for
obtaining linear-programming (LP) relaxations for problems in discrete
optimization and integer programming. Dantzig et al. also describe briefly, in
the case of the TSP, how the relaxation can be used to initiate a
combinatorial procedure, aimed at reducing the search space for an optimal
solution. We study an algorithm based on this idea, proposed by Hougardy and
Schroeder [12] in 2014.

In our work, we consider the symmetric form of the TSP, where the distance to
travel from point $a$ to point $b$ is the same as the distance to travel from
$b$ to $a$. An instance of the problem can be described as a complete
(undirected) graph $K_{n}=(V,E(K_{n}))$ with $n$ nodes $V$, edges $E(K_{n})$,
and edge lengths $(d_{e}:e\in E(K_{n}))$. The nodes of $K_{n}$ represent the
points to visit and the edge lengths indicate the distance to travel between
the ends of each edge. The problem asks to find a tour of minimum total edge
length, where a tour is a circuit of $K_{n}$ visiting each node exactly once,
also known as a Hamiltonian circuit.

In the class of LP relaxations considered by Dantzig et al., a tour $T$ with
edge set $E(T)\subseteq E(K_{n})$ is represented in variables $(x_{e}:e\in
E(K_{n}))$ by setting $x_{e}=1$ for each $e\in E(T$) and $x_{e}=0$ for each
$e\notin E(T)$; constraints are linear equations and inequalities satisfied by
all such tour vectors; the objective function is to minimize
$\sum(d_{e}x_{e}:e\in E(K_{n}))$. A dual solution to such an LP relaxation
provides a lower bound on the length of any tour for the given TSP instance.
Letting $L$ denote the value of the lower bound and $U$ the length of a known
TSP tour, then $\Delta\equiv U-L\geq 0$ is called the integrality gap with
respect to the relaxation and tour. An important observation of Dantzig et al.
[5] is that if a variable $x_{e}$ has reduced cost in the LP dual solution
greater than $\Delta$, then $x_{e}$ can be eliminated from the problem, since
the corresponding edge $e$ cannot be included in any optimal TSP tour.
Similarly, if the reduced cost of $x_{e}$ is less than $-\Delta$, then the
value of $x_{e}$ can be fixed to 1, since edge $e$ must be present in every
optimal tour. Concerning this, Dantzig et al. [5] write the following, where
$E$ is used in place of $\Delta$.

> During the early stages of the computation, $E$ may be quite large and very
> few links can be dropped by this rule; however, in the latter stages often
> so many links are eliminated that one can list all possible tours that use
> the remaining admissible links.

In 1959, they published a second paper [6] on the topic, writing “it is our
belief that a linear-programming, combinatorial approach affords a practical
way of solving traveling-salesman problems” but note in conclusion that they
did not “indicate how one could make the combinatorial analysis a routine
procedure.”

Reduced-cost elimination and variable fixing are now standard in TSP
computation, as well as in general integer programming and discrete
optimization. In comparison, little attention has been given to further
restricting the solution space via combinatorial analysis of the type proposed
by Dantzig et al. The goal of this paper is to demonstrate, within the context
of the TSP, the potential power of the combination of the two techniques.

The combinatorial algorithm we study is an edge-elimination/fixing method
developed by Hougardy and Schroeder [12]. The central idea was proposed by
Jonker and Volgenant [13], in the context of a spanning-tree-based branch-and-
bound algorithm for the TSP. Jonker and Volgenant write the following.

> The simple identification rules are based on the fact that a solution can be
> optimal only if it is 2-optimal.

That is, they eliminate certain edges $e$ by proving any tour $T$ containing
$e$ can be improved by exchanging two edges. Hougardy and Schroeder extend
this to general local-improvement moves, possibly exchanging large sets of
edges in a non-optimality proof.

In Section 1, we present Hougardy-Schroeder’s method as a 2-person game,
between an “edge eliminator” and a “tour builder.” Edge $e$ is
eliminated/fixed by demonstrating a winning strategy for the eliminator.
Details of our engineering steps to implement this local elimination process
are given in Section 2. In Section 3, we report on a computational study of
geometric instances, with the number of points $n$ ranging from 3,038 up to
115,475, including all TSPLIB instances having at least 3,000 points. In 14 of
the 16 examples, the complete-graph edge sets were reduced to under 3$n$
edges, and also in 14 of the 16 examples we fixed at least 0.1$n$ edges. In
Section 4, we discuss a data structure to record the elimination process,
certifying the correctness of the computations. We present such certificates
for reduced edge sets of the three open TSP instances in our test set, each
having 100,000 or more points. In each of the three cases, repeated runs of
the elimination process reduced the edge sets to under 2.5$n$ edges. The
certificates allow the elimination results to potentially be employed in an
exact solver, as we discuss in Section 5.

##  1 Local elimination

Suppose we have an instance of the TSP specified by $K_{n}$ and edge lengths
$(d_{e}:e\in E(K_{n}))$. Rather than working with the full edge set, we
consider a graph $G=(V,E)$ where $E\subseteq E(K_{n})$ and $E$ contains all
optimal tours. For example, as in our computational study, $E$ could be the
set of edges that remain after LP reduced-cost elimination is applied to the
full instance. For $u,v\in V$, we denote by $uv$ the edge having ends
$\\{u,v\\}$ and we say $u$ is a neighbor of $v$ if $uv\in E$.

We will discuss a technique for proving an edge $e=ab\in E$ is not contained
in any optimal TSP tour.

###  1.1 Witnesses

Let ${\cal Q}=\\{F_{1},F_{2},\ldots,F_{q}\\}$ be a family of edge sets such
that, for each $i=1,\ldots,q$, $F_{i}\subseteq E$ and $e\in F_{i}$. We call
${\cal Q}$ a witness family for edge $e$, or an $e$-witness family, if it has
the property that every tour $T$ of $G$ having $e\in E(T)$ also has
$F_{i}\subseteq E(T)$ for some $i\in\\{1,\ldots,q\\}$.

  * •

Example 1: A simple $e$-witness family is to take, for a positive integer
$t<n$, the edge sets of all paths of length $t$ containing edge $e$ in the
“center” of the path. To be precise, we say a path $P$ with ordered edges
$(e_{1},e_{2},\ldots,e_{t})$ is $e$-centered if $e_{\lceil t/2\rceil}=e$,
where $\lceil t/2\rceil$ denotes $t/2$ rounded up to the nearest integer. A
tour $T$ with $e\in E(T)$ must contain an $e$-centered path for any specified
$t$ such that $1\leq t\leq n.$

  * •

Example 2: Every tour $T$ of $G$ with $e=ab\in E(T)$ must contain a second
edge meeting node $a$. Therefore, another simple $e$-witness family consists
of the sets $\\{e,ax\\}$ for each neighbor $x$ of $a$ such that $x\neq b$.

  * •

Example 3: Similarly, for any $v\in V\setminus\\{a,b\\}$, we obtain an
$e$-witness family by taking $\\{e,vx,vy\\}$ all for pairs $x,y$ of neighbors
of $v$.

We will use these three examples as building blocks to construct complex
witness families. The goal of the construction process is to produce (possibly
large) edge sets $F$ that are incompatible with TSP optimality, that is, there
does not exist an optimal tour $T$ with $F\subseteq E(T)$. If edge $e$ has a
witness family ${\cal Q}$ such that each $F\in{\cal Q}$ is incompatible with
TSP optimality, then $e$ can be eliminated from the TSP instance.

###  1.2 Nowhere $k$-optimality

A $k$-opt move, for a tour $T$ and integer $k\geq 2$, consists of replacing
$k$ edges in $T$ to obtain a shorter tour $T^{\prime}$. That is, we delete $k$
edges and reconnect the resulting paths into a tour, where the set of $k$
edges added to reconnect the paths has total length less than that of the
deleted edges. For convenience, we permit the case where a deleted edge is one
of those added back, so a $k$-opt move may also be an $l$-opt move for some
$l<k$.

A $k$-opt move certifies tour $T$ is not optimal. We make repeated use of this
simple fact in a certificate showing that $F\subseteq E$ is incompatible with
TSP optimality.

In this discussion, it is convenient to express $F$ as the union of edge sets
of node-disjoint paths $p_{1},\ldots,p_{m}$. We call
$P=\\{p_{1},p_{2},\ldots,p_{m}\\}$ a path system. For each $i=1,\ldots,m$, let
$\vv{p}_{i}$ denote path $p_{i}$ oriented from one end node $s_{i}$ to the
other end node $t_{i}$, and let $\reflectbox{$\vv{\reflectbox{$p$}}$}_{i}$
denote the reverse path, from $t_{i}$ to $s_{i}$.

Letting ${\cal T}_{F}$ denote the set of tours of $G$ containing $F$, then an
orientation of any $T\in{\cal T}_{F}$ will visit consecutively the edges in
each path $p_{i}\in P$, either in the $\vv{p}_{i}$ or
$\reflectbox{$\vv{\reflectbox{$p$}}$}_{i}$ direction. To fix the orientation
of $T$, we select $s_{1}$ as the tour’s start node and orient $T$ to begin
with a traversal of $\vv{p}_{1}$. The tour then traverses $P$ in the path
ordering

$[\vv{p}_{1},\overline{p}_{\sigma(2)},\ldots,\overline{p}_{\sigma(m)}]$

where $\sigma$ is a permutation of $\\{1,\ldots,m\\}$ with $\sigma(1)=1$ and
$\overline{p}_{i}$ denotes either $\vv{p}_{i}$ or
$\reflectbox{$\vv{\reflectbox{$p$}}$}_{i}$, for each $i=2,\ldots,m.$ An
illustration of such a traversal is given in Figure 1.

Figure 1: Tour $T$ traversing the ordered path system
$[\vv{p}_{1},\vv{p}_{3},\reflectbox{$\vv{\reflectbox{$p$}}$}_{2}]$.

We partition ${\cal T}_{F}$ according to each tour’s traversal of the path
system $P$. Thus there are $2^{m-1}(m-1)!$ classes in the partition, one for
each possible path ordering. Let $\cal{C}$ be a class in the partition and
consider a tour $T\in{\cal C}$. Observe that a $k$-opt move for $T$ such that
all deleted edges belong to $F$ is also a $k$-opt move for all other tours in
${\cal C}$, since the validity for the $k$-opt move does not depend on the
segments of $T$ that connect the end nodes of the path ordering. See the
illustration in Figure 2.

Figure 2: 4-opt move, deleting $f_{1},f_{2},f_{3},f_{4}\in F$ and adding edges
$e_{1},e_{2},e_{3},e_{4}$.

We call $F\subseteq E$ nowhere $k$-optimal if each class of tours in the path-
ordering partition has a $k$-opt move with all deleted edges in $F$. Finding
such a collection of $k$-opt moves is our mechanism for showing $F$ is
incompatible with TSP optimality. The fact that the edges in the $k$-opt moves
are entirely among the nodes in $F$ allows the elimination algorithm to work
locally in the graph $G$, permitting the application of the process on large-
scale instances of the TSP.

###  1.3 Hamilton versus Tutte

The complexity of proving a set $F$ is nowhere $k$-optimal grows rapidly with
$|F|$ and with the number of paths in the corresponding path system $P_{F}$.
Thus the elimination procedure aims to build a witness family ${\cal Q}$ for a
target edge $e\in E$ such that each $F\in{\cal Q}$ has small values for these
two measures.

The process is akin to a game, where one player, Hamilton, reveals edges
purported to be in an optimal tour and the second player, Tutte, attempts to
show the revealed edge set is nowhere $k$-optimal. As the first move, Tutte
selects an integer $t\geq 1$ and requests Hamilton reveal the edges $F$ in an
$e$-centered path of length $t$ containing $e$. In each subsequent move, Tutte
selects a node $v\in V$ and requests Hamilton reveal two edges meeting $v$. If
$v$ is an end node of a path in $P_{F}$, then Hamilton’s reveal adds a single
edge to $F$. Otherwise, Hamilton’s reveal adds two edges to $F$. The game ends
in one of two ways. Tutte wins if $F$ is shown to be nowhere $k$-optimal,
while Hamilton wins if the complexity of $F$, measured as function of $|F|$
and $|P_{F}|$, exceeds a prescribed threshold.

To eliminate edge $e$, we show that Tutte has a winning strategy. That is, we
prepare a response to each possible Hamilton reveal.

  * •

In the first move, we consider every extension of $e$ to an $e$-centered path
of length $t$ in graph $G$.

  * •

In a Tutte move where the selected $v$ is an end node in $P_{F}$, we consider
every extension of the path containing $v$ by a single edge $(v,w)\in E$.

  * •

In a Tutte move where the selected $v$ is not part of $P_{F}$, we consider
each pair of edges $(u,v),(v,w)\in E$.

Tying this back to our earlier discussion, note that the Hamilton reveals
follow the witness building blocks described in Section 1.1.

The steps of a game strategy can be recorded in a tree, where non-leaf tree
nodes are associated with Tutte moves and tree edges are associated with the
possible Hamilton reveals. Thus a non-leaf tree node has a child for every
Hamilton reveal that can be made in response to the specified Tutte move. Each
node $x$ of the tree has a corresponding edge set $F_{x}$, obtained as the
union of the Hamilton reveals associated with the edges on the path from $x$
up to the root of the tree. In a winning strategy, each $F_{x}$ is below the
complexity threshold and for leaf nodes $x$, the set $F_{x}$ is shown to be
nowhere $k$-optimal. The corresponding witness family is

${\cal Q}=\\{F_{x}:x\mbox{\ is a leaf node of the game tree}\\}.$

An example of a Hamilton-Tutte game tree is illustrated in Figure 3. In this
toy instance, the target edge $e$ is in the center of the graph and the
initial Tutte move requests only the path of length $t=1$. Subsequent Tutte
moves are indicated by the grey nodes, numbered 1, 2, and 3. At each node of
the tree, the union of the Hamilton reveals (the sets $F_{x}$) are indicated
by the thick (red) edges. There are five leaf nodes where the indicated edge
sets would need to be shown to be nowhere $k$-optimal.

Figure 3: Hamilton-Tutte game tree with initial path length $t=1$.

In general, to find a winning game tree, we adopt a backtracking search
procedure, where we consider a small number of candidate Tutte moves at each
internal tree node $x$. While processing $x$, if a candidate Tutte move leads
to a Hamilton reveal exceeding the complexity threshold, we consider $x$’s
next Tutte candidate. If we exhaust the supply of candidates at $x$, we
backtrack to $x$’s parent node in the game tree, and consider the next Tutte
candidate for the parent.

The tree-building process is outlined in pseudocode in Algorithm 1. The main
function, HamiltonTutte, takes as input $e\in E$ and a length $t$ for the
initial $e$-centered path to be revealed by Hamilton. Backtracking is handled
by recursive calls to a function TutteMove, that considers one by one the list
of Tutte candidates for an edge set $F$.

Algorithm 1 Building a Hamilton-Tutte tree.

HamiltonTutte$e$, $t$ $\triangleright$ $e$ an edge, $t$ an integer each
$e$-centered $t$-length path $P$ containing edge $e$ TutteMove$E(P)$ returns
Failure

Failure

Success

TutteMove$F$ $\triangleright$ $F$ an edge set $F$ exceeds the complexity
threshold

Failure $F$ is nowhere $k$-optimal

Success each candidate Tutte move $v$ HamiltonReveal$F,v$ returns Success

Success

Failure

HamiltonReveal$F,v$ $\triangleright$ $F$ an edge set, $v$ a node each Hamilton
reveal $H$ for node $v$ TutteMove$F\cup H$ returns Failure

Failure

Success

This is a rather simple strategy for producing the witness family ${\cal Q}$,
but it has flexibility in its possible implementation. For example, short
searches can be used in a pre-processing step to potentially improve the
practical performance of TSP algorithms for modest-sized instances, such as in
the original work of Jonker and Volgenant [13]. While long searches give the
possibility to reduce the search space to such a degree that new techniques
can be considered in the exact solution of very large instances.

##  2 Implementation

In this section, we present general ideas for improving the time to search for
a winning Hamilton-Tutte tree, together with specific choices made in our
implementation of these ideas.

###  2.1 Verifying nowhere $k$-optimality

We begin with techniques for showing an edge set $F\subseteq E$ is nowhere
$k$-optimal, and thus incompatible with TSP optimality.

####  2.1.1 Reusing $k$-opt moves

As in the previous section, we express $F$ as a path system
$P_{F}=\\{p_{1},p_{2},\ldots,p_{m}\\}$, where for each $i=1,\ldots,m$ we let
$s_{i},t_{i}$ denote the end nodes of path $p_{i}$. To show $F$ is nowhere
$k$-optimal, we must consider all path orderings

$P_{F}(\sigma,{\bf
b})\equiv[\vv{p}_{1},\overline{p}_{\sigma(2)},\ldots,\overline{p}_{\sigma(m)}]$

where $\sigma$ is a permutation of $\\{1,\ldots,m\\}$ with $\sigma(1)=1$,
${\bf b}$ is an m-dimensional binary vector with ${\bf b}[1]=0$, and for
$i=2,\ldots,m$

$\overline{p}_{i}\equiv\begin{cases}\vv{p}_{i}&\text{if ${\bf b}[i]=0$}\\\
\reflectbox{$\vv{\reflectbox{$p$}}$}_{i}&\text{if ${\bf b}[i]=1.$}\end{cases}$

Notice that $P_{F}(\sigma,{\bf b})$ can be specified by listing the edges
$\\{o_{1},\ldots,o_{m}\\}$ used to join the paths into a circuit following the
prescribed ordering. To be precise, for $i=1,\ldots,m$, we let

$o_{i}\equiv\begin{cases}(t_{\sigma(i)},s_{\sigma(i+1)})&\text{if ${\bf
b}[\sigma(i)]=0$ and ${\bf b}[\sigma(i+1)]=0$}\\\
(t_{\sigma(i)},t_{\sigma(i+1)})&\text{if ${\bf b}[\sigma(i)]=0$ and ${\bf
b}[\sigma(i+1)]=1$}\\\ (s_{\sigma(i)},s_{\sigma(i+1)})&\text{if ${\bf
b}[\sigma(i)]=1$ and ${\bf b}[\sigma(i+1)]=0$}\\\
(s_{\sigma(i)},t_{\sigma(i+1)})&\text{if ${\bf b}[\sigma(i)]=1$ and ${\bf
b}[\sigma(i+1)]=1$}\end{cases}$

where $\sigma(m+1)=1$. We refer to this as the outside matching for
$P_{F}(\sigma,{\bf b})$ and denote the edge set by $O_{F}(\sigma,{\bf b})$. In
Figure 4 we illustrate the outside matching for the tour displayed in Figure
1; in this example $\sigma=(1,3,2)$ and ${\bf{b}}=[0,1,0]$.

Figure 4: Outside matching $\\{o_{1},o_{2},o_{3}\\}$ and ordered path system
$[\vv{p}_{1},\vv{p}_{3},\reflectbox{$\vv{\reflectbox{$p$}}$}_{2}]$.

By construction, adding $O_{F}(\sigma,{\bf b})$ to $F$ gives a TSP tour
through the union of the ends of $F$; we denote this tour by
$T_{F}(\sigma,{\bf b})$. For all choices of $\sigma$ and ${\bf b}$, we must
find a $k$-opt move for $T_{F}(\sigma,{\bf b})$ that deletes only edges in
$F$, leaving the outside matching in the tour. We discuss below the techniques
we use to find such a $k$-opt move. Here we observe that, once found, a single
$k$-opt move can possibly be applied to handle more than one choice of
$\sigma$ and ${\bf b}$.

Indeed, consider a $k$-opt move for the tour $T_{F}(\sigma,{\bf b})$ that
deletes only edges in $F$ and let $T^{\prime}$ denote the shorter tour
produced by the move. Deleting the outside matching $O_{F}(\sigma,{\bf b})$
from $T^{\prime}$ leaves node-disjoint paths $q_{1},\ldots,q_{m}$ that
together contain all nodes in $T_{F}(\sigma,{\bf b})$. Each path $q_{j}$
starts and ends at nodes $v_{j_{1}}$ and $v_{j_{2}}$ from the set
$\\{s_{1},t_{1},\ldots,s_{m},t_{m}\\}$. Letting $i_{j}$ denote the edge
$v_{j_{1}}v_{j_{2}}$ for each $j=1,\ldots,m$, we call
$\\{i_{1},\ldots,i_{m}\\}$ the inside matching for the $k$-opt move. This
construction is illustrated in Figure 5 for a 4-opt move applied to the
example displayed in Figure 4.

Figure 5: $4$-opt move and inside matching $\\{i_{1},i_{2},i_{3}\\}$.

The union of two matchings $\\{o_{1},\ldots,o_{m}\\}$ and
$\\{i_{1},\ldots,i_{m}\\}$ forms a circuit through the nodes
$\\{s_{1},t_{1},\ldots,s_{m},t_{m}\\}$. Moreover, the $k$-opt move that
produced $T^{\prime}$ can be applied to any tour
$T_{F}(\hat{\sigma},{\bf\hat{b}})$ such that the ordered path system
$P_{F}(\hat{\sigma},\hat{\bf b})$ corresponds to an outside matching
$\\{\hat{o}_{1},\ldots,\hat{o}_{m}\\}$ that forms a circuit with
$\\{i_{1},\ldots,i_{m}\\}$. Such a second outside matching is illustrated in
Figure 6, corresponding to the ordered path system
$[\vv{p}_{1},\reflectbox{$\vv{\reflectbox{$p$}}$}_{2},\vv{p}_{3}]$.

Figure 6: Outside matching $\\{\hat{o}_{1},\hat{o}_{2},\hat{o}_{3}\\}$
compatible with $\\{i_{1},i_{2},i_{3}\\}$.

We make direct use of this observation, creating a list of the $2^{m-1}(m-1)!$
possible outside matchings to process, corresponding to all choices of
$\sigma$ and ${\bf b}$. When a $k$-opt move is found, we extract the inside
matching $\\{i_{1},\ldots,i_{m}\\}$ and remove from the list all remaining
outside matchings that form a circuit with $\\{i_{1},\ldots,i_{m}\\}$.

####  2.1.2 Small incompatible sets

Complex paths systems $P_{F}$ may require many $k$-opt moves to prove $F$ is
incompatible with TSP optimality, even when reusing inside matchings. It is
therefore useful to have simple sufficient conditions for incompatibility,
allowing us to quickly handle a portion of the Hamilton reveals.

To begin, in one of the first research papers on the TSP, Flood [9] discusses
the case where $F$ consists of two edges.

> There is one useful general theorem, which is quickly discovered by each one
> who considers the traveling-salesman problem. In the euclidean plane it
> states simply that the minimal tour does not intersect itself, and this
> intersection condition generalizes easily for arbitrary [edge lengths]
> $a_{\alpha\beta}$.

The general condition is $F=\\{ab,xy\\}$ is incompatible with TSP optimality
if

$\max\\{d_{ax}+d_{by},d_{ay}+d_{bx}\\}<d_{ab}+d_{xy}.$ (1)

This follows from the fact that deleting $ab$ and $xy$ from any tour
containing both edges yields two paths that can be reconnected into a tour by
adding either $\\{ax,by\\}$ or $\\{ay,bx\\}$. Condition (1) implies that, in
either case, we improve the length of the tour with the 2-opt move. We make
frequent use of this simple test, checking for incompatible pairs of edges as
we consider new Hamilton reveals.

A second simple incompatibility test arises when $F$ consists of three edges
$ab$, $xy$, and $yz$. Here we assume the five nodes $a,b,x,y,z$ are distinct.
Thus $P_{F}$ consists of paths $a-b$ and $x-y-z$. In this case, the condition

$d_{ay}+d_{by}+d_{xz}<d_{ab}+d_{xy}+d_{yz}$ (2)

is sufficient to show $F$ is incompatible with TSP optimality. Indeed, the
single 3-opt move, replacing $F$ by $ay,by,$ and $xz$ is valid for any tour
containing $F$. The move is illustrated in Figure 7.

Figure 7: A 3-opt move valid for any tour containing edges $ab$, $xy$, and
$yz$.

Its usefulness comes from the fact that Hamilton reveals, following an
isolated-node Tutte move, are paths of length two that can be immediately
tested for incompatibility with the edge $ab$ we are trying to eliminate.

####  2.1.3 Brute force for small $k$-opt

For general sets $F$, we search for $k$-opt moves applied to a single tour,
rather than aiming for universal conditions such as (1) and (2). In this
search, as long as $|F|$ is not too large, small $k$-opt moves can be handled
by brute force, considering every subset of $k$ edges to delete from $F$ and
examining all possible ways to reconnect the paths into a tour.

In our implementation, we start with $k=3$ and proceed to $k=4$ and $k=5$ if
we do not find a $k$-opt move. To improve the speed of these computations, at
the expense of possibly missing a move, we consider only $k$-opt moves that
include deleting the edge $ab$ we are trying to eliminate.

####  2.1.4 Adopting a TSP solver

Moving up to larger values of $k$, we build a TSP instance consisting of nodes

$V_{F}=\\{v:v\mbox{\ is an end of an edge in}\ F\\}$

and edge lengths obtained from the full TSP instance. For a specified
permutation $\sigma$ and binary vector ${\bf b}$, we create the tour
$T_{F}(\sigma,{\bf b})$ consisting of the paths $P_{F}$ and the outer matching
$O_{F}(\sigma,{\bf b})$.

Working with this TSP instance, any tour $T^{\prime}$ that includes all of the
edges in $O_{F}(\sigma,{\bf b})$ and has length less than that of
$T_{F}(\sigma,{\bf b})$ corresponds to a $k$-opt move. Indeed, we delete the
edges $T_{F}(\sigma,{\bf b})\setminus T^{\prime}$ and add the edges
$T^{\prime}\setminus T_{F}(\sigma,{\bf b})$.

Any heuristic or exact TSP code can be used to search for $T^{\prime}$. If the
code does not directly support the option of fixing edges to be in the tour, a
large constant $M$ can be added to the cost of each edge not in
$O_{F}(\sigma,{\bf b})$. Note that the validity of the elimination procedure
is not dependent on the correctness of the TSP solver that is adopted, since
the elimination code can easily verify that the $k$-opt move obtained does
indeed produce a shorter tour.

In our implementation, we adopt the Held-Karp 1-tree code [10] included with
the Concorde TSP solver [1]. The branch-and-bound search is terminated early
if either a tour shorter than $T_{F}(\sigma,{\bf b})$ is found or a specified
maximum number of search nodes is reached.

###  2.2 Tutte moves

Our method for selecting Tutte moves is based on simple heuristic ideas.
First, we prefer nodes $v$ that are in close proximity to the edge $ab$ we are
trying to eliminate. The measure of closeness is either geometric distance to
the center of $ab$ or, for non-geometric instances, the minimum number of
edges in a path to $a$ or $b$ in the graph $G$. Secondly, we take a greedy
approach, preferring moves that generate a fewer number of Hamilton reveals
requiring further Tutte moves to process. This greedy measure is obtained with
a look-ahead routine, selecting a small number of candidate nodes $v$ and
counting the number of Hamilton reveals that cannot be discarded using the
2-opt condition (1), the 3-opt condition (2), and brute-force tests for 3-opt,
4-opt, and 5-opt moves. Thirdly, to have variety in the search, for each
Hamilton reveal we consider both a small number of nodes $v$ that are not
contained in the current set $F$ and, if each of these moves returns a
failure, a small number of nodes $v$ that are end nodes of paths in $P_{F}$.

When the initial Tutte move requests only a path of length one, we find it
useful to consider a pair of follow-up Tutte moves, selecting as the second
and third moves nodes $c$ and $d$, distinct from $a$ and $b$, such that edge
$cd$ is either not in the graph $G$ or $cd$ is incompatible with $ab$. The
fact that $cd$ is not available as an edge in an optimal tour adds structure
to the corresponding pair of Hamilton reveals and appears to benefit our
search process, when compared to making independent choices for $c$ and $d$.
If the small number of pairs of Tutte moves $c$ and $d$ we try each return a
failure, then we carry out a second search, requesting an initial Tutte move
consisting of an $e$-centered path of length two.

###  2.3 Thresholds

The elimination process requires complexity thresholds on the edge sets $F$
and path systems $P_{F}$ created by Hamilton reveals. In our code, this is
accomplished by setting $|P_{F}|\leq 5$ and giving an upper bound on $|F|$
controlled by command-line parameters for the length of the path requested in
the initial Tutte move and for the cumulative number of edges added by
subsequent Hamilton reveals.

To further restrict the overall complexity of the Hamilton-Tutte trees, for
each Tutte move, other than the initial path request, we also bound the
maximum number of Hamilton reveals that are not discarded by condition (1),
condition (2), and brute-force tests for 3-opt, 4-opt, and 5-opt moves. In our
code this bound is set to 5, and switched to 1 for tree nodes for which 10 or
more edges have been added to $F$, beyond the edges in the initial path. This
practice of switching the bound to 1 permits the code to carry out deep
searches in some circumstances.

###  2.4 Fixing edges

The discussion has focused on eliminating edges from the TSP search space, but
the techniques can also be applied to fixing edges in the tour, that is,
proving an edge $e=ab$ must be present is every optimal tour. This can be
done, for example, by showing that for all pairs of edges $au,av$ meeting node
$a$ but not meeting node $b$ and all pairs of edges $bx,by$ meeting node $b$
but not meeting node $a$, the set $\\{au,av,bx,by\\}$ is incompatible with TSP
optimality. Ruling out all such pairs of two-edge paths through $a$ and $b$,
we can conclude that $ab$ must be in every optimal tour. See the illustration
in Figure 8.

Figure 8: Pair of two-edge paths for fixing edge $ab$.

###  2.5 Non-pairs

In long computations, it can be beneficial to list pairs of edges $xy,yz$ that
are known to be incompatible with TSP optimality. Hamilton reveals containing
such pairs can then be discarded during edge elimination or edge fixing. To
build this non-pair list, for all nodes $y$ and edges $xy,yz$ meeting $y$, we
search for a winning Hamilton-Tutte tree starting with the two-edge path
$x-y-z$. In this setting, the leaf nodes of the game tree correspond to a
$\\{xy,yz\\}$-witness family, defined analogously to the single-edge witness
families introduced in Section 1.1.

In our implementation of this pair-elimination process, if the initial Tutte
move requests only the path of length two, then as the second Tutte move we
consider neighbors $v$ of $y$ such that $v\neq x$ and $v\neq z$. This allows
the algorithm to take immediate advantage of the fact that $yv$ cannot be in a
tour containing $x-y-z$.

##  3 Computation

The elimination process requires additional choices and settings to produce a
working code suitable for large-scale instances. In this section we present
further details, together with test results using our implementation.

###  3.1 Test instances

Our study focuses on the potential use of edge elimination/fixing to enhance
exact solution algorithms for the TSP. We thus consider as a test set
instances of size beyond the scope of simple optimization methods. The main
source is the standard TSPLIB [18, 19] collection, which contains 13 examples
having 3,000 or more nodes, the largest having 85,900.

Each of the large TSPLIB instances is defined by a set of points in the
Euclidean plane. Travel distances are determined either by the EUC_2D norm,
defined as the Euclidean distance rounded to the nearest integer, or the
CEIL_2D norm, defined as the Euclidean distance rounded up to the next
integer.

In addition to the 13 TSPLIB instances, we include in our test set three
larger examples, E100k.0, mona-lisa100k, and usa115475, all adopting the
EUC_2D norm.

  * •

E100k.0 contains 100,000 points with integer coordinates drawn uniformly from
the 1,000,000 $\times$ 1,000,000 square. It was created in 2000 by David S.
Johnson for the DIMACS TSP Challenge. Code for generating the point set is
available at http://dimacs.rutgers.edu/archive/Challenges/TSP/download.html.

  * •

mona-lisa100k was created in 2009 by Robert Bosch using his TSP Art methods
[2] to render da Vinci’s Mona Lisa. The 100,000-point instance can be
downloaded at https://www.math.uwaterloo.ca/tsp/data/ml/monalisa.html.

  * •

usa115475 consists of the latitude and longitude of 115,475 cities, towns, and
villages in the United States. The point set, created in 2012, is available at
https://www.math.uwaterloo.ca/tsp/data/usa/index.html.

These three examples have not yet been solved to proven optimality, continuing
to provide a challenge for new optimization techniques.

###  3.2 Sparse edge sets

The 16 TSP instances in our test set are defined on complete graphs. To obtain
sparse graphs for the elimination process, we use as a pre-processing routine
reduced-cost elimination. To carry this out, we adopt the LKH and Concorde TSP
solvers. LKH is a TSP heuristic code developed by Helsgaun [11] and Concorde
is the exact solver of Applegate et al. [1], based on the cutting-plane
method.

For the 11 instances having fewer than 20,000 nodes, we execute the pre-
processing in a “live” fashion, simulating the potential use of edge
elimination by a TSP solver. In these cases, we make a run of LKH to produce a
high-quality tour and a run of Concorde to produce an LP relaxation.111LKH
version 2.09 was run with settings CANDIDATE_SET_TYPE = POPMUSIC,
MAX_CANDIDATES = 7, RUNS = 1. Concorde was run with command-line option -mC48
(to adopt the local cuts procedure up to chunksize 48) and with an option to
continue the cutting-plane method until no violated cuts have been added to
the LP relaxation. Due to an excessively long run time, the fl3795 instance
was run with Concorde -mC28. Reduced-cost elimination is then carried out
using the LP relaxation and an upper bound of 1 greater than the length of the
LKH tour, to ensure that only edges not in any optimal tour are eliminated.
The results are reported in Table 1.

Table 1: LP Reduced-Cost Elimination (LKH+Concorde) Name | LKH | Concorde | Opt Ratio | LP Edges
---|---|---|---|---
pcb3038 | 163.7s | 5460.2s | 1.00006 | 6883 | 2.3$n$
fl3795 | 4236.2s | 100004.6s | 1.00053 | 61609 | 16.2$n$
fnl4461 | 179.4s | 5074.1s | 1.00004 | 10343 | 2.3$n$
rl5915 | 457.9s | 40867.6s | 1.00011 | 29143 | 4.9$n$
rl5934 | 490.3s | 9072.4s | 1.00011 | 28678 | 4.8$n$
pla7397 | 9102.7s | 80363.4s | 1.00011 | 42180 | 5.7$n$
rl11849 | 7935.2s | 43422.1s | 1.00010 | 78993 | 6.7$n$
usa13509 | 7906.9s | 28402.2s | 1.00009 | 136732 | 10.1$n$
brd14051 | 10510.4s | 61880.5s | 1.00008 | 95959 | 6.8$n$
d15112 | 11117.5s | 38843.6s | 1.00011 | 166499 | 11.0$n$
d18512 | 17449.5s | 56388.8s | 1.00008 | 178825 | 9.7$n$

The TSPLIB identifier in the ‘Name’ column includes the number of nodes in the
graph. The running times for LKH and Concorde are reported in seconds on a
single core of a linux workstation equipped with two 6-core Intel Xeon CPU
E5-2620 2.00GHz processors. The ‘Opt Ratio’ column reports the length of the
LKH tour divided by the value of the LP relaxation. We remark that in 7 of the
11 test instances, LKH actually delivered an optimal tour. Finally, the ‘LP
Edges’ column reports the number of edges that were not eliminated via the
reduced-cost technique. For instances over 5,000 nodes, this number ranges
from 4.8$n$ up to 11.0$n$, where $n$ denotes the number of nodes in the
instance.

For the 5 largest instances in our test set, reduced-cost elimination based on
single runs of LKH and Concorde produce very large edge sets that are
difficult to manage, both in the algorithm and in the storage of the data. For
these instances, that are either unsolved or tractable only with very long
runs of Concorde, we used the best available tour (found by repeated runs of
LKH) and the best available LP relaxation (found by runs of Concorde). The
reduced-cost elimination results are reported in Table 2.

Table 2: LP Reduced-Cost Elimination (Best Available Tour+LP) Name | Opt Ratio | Edges
---|---|---
pla33810 | 1.000164 | 1,602,936 | 47.4$n$
pla85900 | 1.000007 | 266,273 | 3.1$n$
E100k.0 | 1.000084 | 8,280,715 | 82.8$n$
mona-lisa100k | 1.000019 | 1,509,670 | 15.1$n$
usa115475 | 1.000104 | 25,009,702 | 216.6$n$

###  3.3 Computer code

The elimination process is implemented in 7,781 lines of code written in the C
programming language. The full code includes also 7,364 lines of code from the
Concorde library, handling various utility functions, as well as the Held-Karp
1-tree solver used to check for possible $k$-opt moves. The code links only
the standard C library, allowing it to be run on a wide variety of computing
platforms. The source code is available at
https://github.com/bicobico2/ElimTSP.

The main executable function, called elim, takes as input a full TSP instance
in TSPLIB file format and a sparse graph $G$ expressed as a list of edges.
Optionally, a list of non-pairs, a list of fixed edges, and an input tour can
be specified. These three optional data sets can improve the computation by
ruling out certain Hamilton reveals or by reducing the number of edges to
process. In the case of the input tour, the code will not attempt to eliminate
any edge that is in the tour and it will not attempt to fix any edge that is
not in the tour.

The code is set up to operate in a boss-worker parallel mode. Workers receive
from the boss a small subset of edges $S$ to process, they carry out the
elimination algorithm for each edge in $S$, and return to the boss the
elimination status for each edge in $S$. The boss is in a listening mode,
responding to worker requests for new subsets of edges to process or to
receive the result of worker computations. Communication between the boss and
workers is enabled by socket-based message-passing routines from the Concorde
library.

###  3.4 Bootstrapping

The number of Hamilton reveals in response to a given Tutte move depends
directly on the sparsity of the input graph $G$. This dependency makes it
possible to run the elimination process in a bootstrapping fashion, quickly
eliminating edges with small Hamilton-Tutte trees in an initial pass, then
applying the algorithm again starting with the reduced graph. This can be done
repeatedly, slowly increasing the code’s parameters to allow more complex
searches for Hamilton-Tutte trees to eliminate additional edges.

In our bootstrapping run, we control the complexity of the Hamilton-Tutte
search via the following three settings.

  * –

neighbors: the number of nodes, in the proximity of the edge to be eliminated,
that are considered as possible candidates for Tutte moves.

  * –

depth: the maximum depth of a Tutte node in the Hamilton-Tutte tree.

  * –

fast: in fast mode, we test only one candidate $(c,d)$ pair for the first two
Tutte moves and we do not employ the TSP solver to search for $k$-opt moves
with $k>5$.

In addition, after a sequence of fast elimination rounds, we run a computation
to build a non-pair list to improve the speed of further rounds of edge
elimination. Finally, using the reduced edge set and the non-pair list, we run
a computation to determine fixed edges, that is, edges that must be in every
optimal tour.

The bootstrapping loop proceeds in levels, where the settings are increased to
the next level if less than 5% of the remaining edges are eliminated or less
than 25% of the remaining pairs are eliminated or less than 5% of the non-
fixed edges in the input tour have been fixed. The level settings we use are
displayed in Table 3.

Table 3: Settings in Bootstrap Run Type | Level | Depth | Neighbors
---|---|---|---
Fast edges | 1 | 2 | 5
Fast edges | 2 | 2 | 10
Fast edges | 3 | 3 | 5
Fast edges | 4 | 3 | 10
Fast edges | 5 | 4 | 10
Fast edges | 6 | 4 | 15
Fast edges | 7 | 5 | 25
Fast edges | 8 | 6 | 25
Non-pairs | 1 | 2 | 5
Non-pairs | 2 | 3 | 10
Non-pairs | 3 | 4 | 25
Edges | 1 | 4 | 25
Edges | 2 | 4 | 50
Edges | 3 | 5 | 25
Edges | 4 | 6 | 25
Edges | 5 | 6 | 50
Fixed | 1 | 6 | 25

A full run adopts eight levels of fast edge elimination, followed by three
levels of pair elimination, then five levels of standard edge elimination, and
finally one level of edge fixing.

Computational results for the loop are reported in Table 4.

Table 4: Edge Elimination (48-core network)  Name | Fast edges | Non-pairs | Edges | Fixed | Clock
---|---|---|---|---|---
pcb3038 | 6521 | 2.1$n$ | 49.4% | 5548 | 1.8$n$ | 934 | 0.31$n$ | 497s
fl3795 | 21028 | 5.4$n$ | 20.2% | 15921 | 4.2$n$ | 628 | 0.17$n$ | 27203s
fnl4461 | 9928 | 2.2$n$ | 54.8% | 9012 | 2.0$n$ | 973 | 0.22$n$ | 462s
rl5915 | 17582 | 3.0$n$ | 29.0% | 11504 | 1.9$n$ | 3503 | 0.59$n$ | 7461s
rl5934 | 17114 | 2.9$n$ | 33.3% | 12007 | 2.0$n$ | 3258 | 0.55$n$ | 6255s
pla7397 | 21779 | 2.9$n$ | 72.6% | 19989 | 2.7$n$ | 1484 | 0.20$n$ | 4339s
rl11849 | 38471 | 3.2$n$ | 35.3% | 32239 | 2.7$n$ | 4042 | 0.34$n$ | 7882s
usa13509 | 39686 | 2.9$n$ | 45.1% | 35211 | 2.6$n$ | 2950 | 0.22$n$ | 3844s
brd14051 | 43400 | 3.1$n$ | 50.8% | 39488 | 2.8$n$ | 1577 | 0.11$n$ | 4177s
d15112 | 46662 | 3.1$n$ | 50.7% | 42271 | 2.8$n$ | 1528 | 0.10$n$ | 3893s
d18512 | 58537 | 3.2$n$ | 50.9% | 58537 | 2.9$n$ | 1814 | 0.10$n$ | 5196s
pla33810 | 143197 | 4.2$n$ | 38.8% | 124331 | 3.7$n$ | 2580 | 0.08$n$ | 88483s
pla85900 | 232768 | 2.7$n$ | 74.3% | 214910 | 2.5$n$ | 17372 | 0.20$n$ | 29680s
E100k.0 | 324292 | 3.2$n$ | 43.9% | 297603 | 3.0$n$ | 16879 | 0.17$n$ | 39065s
mona-lisa100k | 299185 | 3.0$n$ | 80.4% | 278341 | 2.8$n$ | 225 | 0.00$n$ | 13762s
usa115475 | 374029 | 3.2$n$ | 45.9% | 336642 | 2.9$n$ | 13517 | 0.12$n$ | 67042s

The runs were carried out on a 48-core network of 4 linux compute servers,
each equipped with two 6-core Intel Xeon CPU E5-2620 2.00GHz processors. The
‘Clock’ column reports the wall-clock time in seconds for each test run. The
‘Fast edges’ and ‘Edges’ columns report the remaining edges after the
corresponding rounds of elimination, expressed both as the total number and as
a multiple of $n$, the number of nodes in the given instance. The ‘Non-pairs’
column reports the number of non-pairs, expressed as a percentage of the total
number of two-edge paths in the reduced graph. Finally, the ‘Fixed’ column
reports the number of fixed edges, expressed as the total number and as a
multiple of $n$. Averaged over the 16 test instances, the proportion of wall-
clock time used in each division was 9% for fast edges, 5% for non-pairs, 74%
for edges, and 12% for fixed.

The elimination runs give only minor reductions for pcb3038, fnl4461, and
pla859000, where reduced-cost elimination had already reduced the edge counts
to 2.3$n$, 2.3$n$, and 3.1$n$, respectively. In all other cases the edge count
was reduced by at least a factor of two, resulting in under 3$n$ edges in 11
of these instances, including the three examples having 100,000 or more nodes.
Also, in 14 of the test instances the runs fixed at least 0.1$n$ of the edges.

###  3.5 Randomly-generated Euclidean instances

As a stress-test for the code, we generated 250,000 Euclidean instances, each
with $n=100$ points having integer coordinates drawn uniformly from the 10,000
$\times$ 10,000 square. The integers were generated as a single stream using
the RngStream package written by Pierre L’Ecuyer; see L’Ecuyer [14] and
L’Ecuyer et al. [15]. Starting edge sets for the instances were produced with
reduced-cost elimination, using the subtour LP relaxation and upper bound set
to 1.001 times the length of a tour found by Concorde’s Lin-Kernighan
heuristic; these choices create larger edge sets than those produced by
Concorde’s LP relaxation and an LKH tour, and avoid cases that give a zero
integrality gap.

In this test, we did not provide input tours. This allowed us to check for
possible errors in the implementation of the algorithm, running Concorde to
find an optimal tour and then verifying that no tour edge was eliminated and
only tour edges were fixed. For elimination and fixing, we again used the
bootstrap loop described in Table 3; this loop is far too time consuming for
such small instances, but it allowed us to test the main components of our
elimination code.

The results of the test are summarized in Table 5. Each run was carried out on
a single hyperthread of an Intel Xeon CPU E5-2620 2.00GHz processor. Columns
‘LP Edges’, ‘Reduced’, and ‘Fixed’ report the number of edges after reduced-
cost elimination, the number edges after the elimination process, and the
number of fixed edges, respectively. The column labeled ‘Elim Time’ reports
the CPU time in seconds for the elimination process. The rows report mean,
median, minimum, and maximum values, taken over the 250,000 instances. We note
that in 1,830 instances, the elimination process was able to fix all 100 edges
in the unique optimal tour.

Table 5: 250,000 Random Euclidean Instances with $n$ = 100 Value | LP Edges | Reduced | Fixed | Elim Time
---|---|---|---|---
Mean | 430.5 | 177.8 | 49.9 | 1014.0s
Median | 403 | 183 | 46 | 1015.7s
Minimum | 159 | 100 | 22 | 186.5s
Maximum | 1673 | 247 | 100 | 3208.5s

In a second test, we generated 250 Euclidean instances with $n=10,000$ and
integer coordinates drawn uniformly from the 10,000,000 $\times$ 10,000,000
square. For these larger instances, the starting edge set was computed via
reduced-cost elimination using a Concorde LP relaxation and an upper bound of
1 greater than the length of an optimal tour. Again, we did not provide input
tours and verified that in each case the elimination results were consistent
with an optimal tour computed by Concorde. The results are reported in Table
6, for runs carried out on a single hyperthread of the Intel Xeon processor.
In each instance, the reduced edge set has under 2.72$n$ edges.

Table 6: 250 Random Euclidean Instances with $n$ = 10,000 Value | LP Edges | Reduced | Fixed | Elim Time
---|---|---|---|---
Mean | 40017.7 | 23305.2 | 2497.7 | 142745.1s
Median | 39809 | 23456 | 2438 | 141126.4s
Minimum | 22393 | 17482 | 1905 | 52195.2s
Maximum | 66318 | 27182 | 4195 | 253582.2s

###  3.6 Fast elimination

The bootstrap loop produces strong results, at the expense of considerable
computation time. In potential applications where shorter times are desired,
an elimination process can be run with Hamilton-Tutte trees having at most two
Tutte moves, taking advantage of the simple path systems that arise in this
case. This is implemented in the KH-elim code, available at
https://github.com/bicobico2/ElimTSP.

It proved useful to begin the process with a pass through the full edge set,
adopting the Jonker-Volgenant [13] strategy of using 2-opt moves to
potentially show very quickly that an edge $(a,b)$ can be eliminated. In this
phase, as a possible single Tutte move, we consider a small number of nodes
$y$ in a neighborhood of $(a,b)$. Rather than considering each pair of edges
meeting node $y$ (that is, every Hamilton reply), the Jonker-Volgenant
strategy is to make the much quicker check that each individual edge $(x,y)$
meeting $y$ satisfies the incompatibility condition (1) with $(a,b)$, thus
showing $(a,b)$ can be eliminated. Following this, we make another pass,
considering for each remaining edge $(a,b)$ up to ten candidate $(c,d)$ pairs
of Tutte moves, where $c$ and $d$ belong to the neighborhood of $(a,b)$. In
this second phase, we proceed as in the full elimination code, considering all
possible Hamilton replies to the $(c,d)$ pair.

Computational results are reported in Table 7 for all test instances having
fewer than 20,000 nodes.

Table 7: Single $(c,d)$ Pair of Tutte Moves Name | LP Edges | Reduced (1 core) | Reduced (44 cores)
---|---|---|---
pcb3038 | 6883 | 2.3$n$ | 6466 | 2.1$n$ | 0.0s | 6466 | 2.1$n$ | 0.1s
fl3795 | 61609 | 16.2$n$ | 31017 | 8.2$n$ | 50.7s | 31043 | 8.2$n$ | 2.2s
fnl4461 | 10343 | 2.3$n$ | 9865 | 2.2$n$ | 0.0s | 9864 | 2.2$n$ | 0.1s
rl5915 | 29143 | 4.9$n$ | 24596 | 4.2$n$ | 1.5s | 25107 | 4.2$n$ | 0.2s
rl5934 | 28678 | 4.8$n$ | 23426 | 3.9$n$ | 1.3s | 23448 | 4.0$n$ | 0.2s
pla7397 | 42180 | 5.7$n$ | 23935 | 3.2$n$ | 2.9s | 24385 | 3.3$n$ | 0.4s
rl11849 | 78993 | 6.7$n$ | 61792 | 5.2$n$ | 6.0s | 61792 | 5.2$n$ | 0.6s
usa13509 | 136732 | 10.1$n$ | 60448 | 4.5$n$ | 33.5s | 60415 | 4.5$n$ | 1.4s
brd14051 | 95959 | 6.8$n$ | 60894 | 4.3$n$ | 9.5s | 60939 | 4.3$n$ | 0.6s
d15112 | 166499 | 11.0$n$ | 73840 | 4.9$n$ | 39.6s | 73851 | 4.9$n$ | 1.9s
d18512 | 178825 | 9.7$n$ | 92220 | 5.0$n$ | 35.3s | 92279 | 5.0$n$ | 1.7s

The computations were carried out on a linux server equipped with two 22-core
Intel Xeon Gold 6238 CPU @ 2.10GHz. The column labeled “Reduced (1 core)”
lists results for runs using a single core of the server and the column
labeled “Reduced (44 cores)” lists results obtained using 44 computation
threads, where inter-thread communication is handled by the OpenMP application
programming interface. For each instance we report the number of remaining
edges, the number of remaining edges as a multiple of the number of nodes $n$,
and the wall-clock time in seconds. The run eliminated on average 43% of the
LP edges, in each case in under 3 seconds of wall-clock time when using the
full 44 cores of the server.

##  4 Certifying results

A long-term goal of this line of research is to build edge elimination into a
tool that can enhance the TSP cutting-plane method, following general ideas
outlined in Dantzig et al. [5, 6]. The target is to extend the range of
current exact solvers, such as the Concorde code.

In this exact-solution setting, it can be helpful to record the Hamilton-Tutte
tree that eliminates a specified edge. The recorded tree can be used in a
shorter verification run to certify the edge is not in any optimal tour. By
this means, results can be gathered from remote (possibly non-trusted)
computing platforms, using a range of parameter settings or variations of the
basic elimination code. The collected results can then be independently
verified, before edges are removed from the cutting-plane computation.

###  4.1 Hamilton-Tutte tree storage

When recording a Hamilton-Tutte tree, we do not store its leaf nodes, since
path systems associated with leaf nodes can be shown to be incompatible with
TSP optimality without need of additional information. With each non-leaf tree
node $x$, other than the root of the tree, we store the corresponding Tutte
move together with the Hamilton reveal associated with the tree edge joining
$x$ to $x$’s parent. Child nodes of $x$ are stored as a list.

Internally, the tree is stored using two types of data structures: httree and
htnode. An httree contains the target type of the tree (edge elimination, edge
fixing, or pair elimination) and the algorithm settings and parameters that
were used in the tree’s construction. An httree contains also a pointer to the
htnode representing the root of the tree.

An htnode for a non-leaf tree node contains the graph node indices for the
corresponding Tutte move and, except for the root of the tree, the graph node
indices for the edges in the Hamilton reveal that led to the Tutte move. An
htnode contains also two (possibly null) pointers, child and next. These
pointers allow us to store the children of a tree node as a linked list: child
points to the head of the list, child$\rightarrow$next points to the second
node of the list, and so on.

Together with the internal data structures, we need a portable file format for
storing and transferring the trees. For this, the node indices and algorithm
parameters are recorded as integers, and are thus easy to store in either text
or binary format in a file. To handle the htnode pointers, we assign each node
of the tree an integer index and use this integer value in place of pointers
to the node. In the file format, we include as header information the edge or
pair we are targeting, together with the algorithm settings. Then, for each
tree node $x$, we have a line containing the index of the parent of $x$,
followed by $x$’s index, the Hamilton reveal information, the Tutte move
information, and finally the indices of the children of $x$.

###  4.2 Verification run

The verification process works through the nodes of the Hamilton-Tutte tree,
starting at the root. The tree nodes are processed in a depth-first order,
that is, while processing a node $x$ we recursively process all children of
$x$ and then return to continue work at $x$’s parent node. The process does
not depend on the correctness of the computations used to build the tree, but
rather uses the stored information to simplify a new elimination run.

To process a tree node $x$, we first list all possible Hamilton reveals in
response to $x$’s Tutte move. This is done as in the full elimination code,
using the input graph $G$ and the path system $P_{x}$ associated with $x$
(that we have by adding $x$’s Hamilton reveal to the path system associated
with the parent of $x$). For each Hamilton reveal $H$ in our list, we first
check if there is a child node $c$ with $H$ as its Hamilton reveal. If there
is such a child, we process node $c$ before continuing with the next reveal in
our list. If there is no such child, then we verify that the path system
obtained by adding $H$ to $P_{x}$ is incompatible with TSP optimality. This
last step is done with the techniques used in the main elimination code.

The verification process is simpler and faster than the initial elimination
runs, since we avoid the backtracking search for Tutte moves and also avoid
tests for incompatibility at non-leaf nodes of the Hamilton-Tutte tree.

###  4.3 Sparse edge sets for the 100k-point instances

We used the verification process in a study of the three largest instances in
our test collection, both as a test of the limits of our elimination code and
to provide edge sets that can potentially be adopted in exact solution attacks
on these currently unsolved examples of the TSP.

The elimination runs made on these instances involved ad hoc changes to the
default settings and minor modifications to the source code, exploring
alternatives to our standard search directions. The computations were carried
out on a 288-core network of linux compute servers, equipped with two 6-core
Intel Xeon CPU E5-2620 2.00GHz processors. The total computational time
(adding together the times on the individual cores) was over 10 core-years in
each instance.

The above details make clear this is not a reproducible experiment. But with
the stored Hamilton-Tutte tree for each eliminated edge, for each fixed edge,
and for each non-pair, we can verify the final outcome of the many ad hoc
runs.

The results are reported in Table 8. The ‘Edges’ and ‘Fixed’ columns list the
number of edges that were not eliminated and the number of fixed edges,
respectively. The ‘Verification’ column reports the running time in seconds to
verify the results and the size in gigabytes for the (text) Hamilton-Tutte
tree files. Each of the sparse graphs has average degree under 5.0, improving
the average degrees of 6.4, 6.0, and 6.4 in the corresponding results reported
in Table 4.

Table 8: Sparse edge sets for 100k+ instances Name | Edges | Fixed | Verification
---|---|---|---
E100k.0 | 241682 | 2.42$n$ | 26106 | 0.26$n$ | 319,731s | 4.6GB
mona-lisa100k | 248972 | 2.49$n$ | 1635 | 0.02$n$ | 189,629s | 5.5GB
usa115475 | 286655 | 2.48$n$ | 22394 | 0.19$n$ | 1,950,878s | 6.1GB

##  5 Remarks

We conclude with several observations on the elimination code and its possible
application.

###  5.1 Code enhancements

Several ideas, described below, showed promise in our computational tests, but
were ultimately removed from the code in an effort to reduce its complexity,
particularly in regards to the verification procedure.

####  5.1.1 Full-witness families

In Section 1.1, we defined the notion of an $e$-witness family, where $e$ is a
target edge to eliminate from a TSP instance specified by a graph $G=(V,E)$
and edge lengths $(d_{f}:f\in E)$. We can similarly define a full-witness
family as a family ${\cal Q}=\\{F_{1},F_{2},\ldots,F_{q}\\}$ of edge sets such
that, for each $i=1,\ldots,q$, $F_{i}\subseteq E$ and every optimal tour $T$
of $G$ has $F_{i}\subseteq E(T)$ for some $i\in\\{1,\ldots,q\\}$. A full-
witness family ${\cal Q}$ can be used as a tool for showing multiple edges can
be eliminated. Indeed, to eliminate an edge $f$, it suffices to show that for
each $i=1,\ldots,q$ the set $F_{i}\cup\\{f\\}$ is incompatible with TSP
optimality. Similarly, during the construction of a Hamilton-Tutte tree, an
edge set $F$ can be shown to be incompatible with TSP optimality by verifying
that for each $i=1,\ldots,q$ the set $F_{i}\cup F$ is incompatible.

To construct a full-witness family, a Hamilton-Tutte tree (with no starting
edge or path) can be grown in a breadth-first manner, working to a prescribed
upper bound on the number of leaf nodes. Such trees can be grown in local
regions of the graph $G$, for example, by selecting a random set of nodes $S$
and growing a tree for each $s\in S$ with Tutte nodes restricted to a neighbor
set of $s$. The process can be bootstraped, in the sense that previously
constructed trees can be used to reduce the size of a new tree, using the
incompatibility test described above.

####  5.1.2 Order-specific Tutte moves

Our code for building a Hamilton-Tutte tree follows the outline given in
Algorithm 1 in Section 1.2. In this algorithm, a Hamilton reveal, producing a
path system $P=\\{p_{1},p_{2},\ldots,p_{m}\\}$, leads to a failure if one or
more of the $2^{m-1}(m-1)!$ path orderings of $P$ cannot be shown to have a
$k$-opt move. An alternative approach is to make a deeper search for an
individual failed path ordering, considering an additional Tutte move to
handle only this specific case.

An advantage of using such order-specific Tutte moves is that there are at
most $2m$ possibilities for inserting each of the Hamilton reveals into the
failed path ordering. This has the potential for a large savings in
computation, if there are few failed path orderings.

####  5.1.3 Metric excess

Hougardy and Schroeder [12] describe a notion of “metric excess” that can be
used to show a path of three edges is incompatible with TSP optimality. The
idea is to consider an additional Tutte move $z$, together with a fast
mechanism to check that each Hamilton reply to $z$ leads to an improving 3-opt
move. For details see [12].

####  5.1.4 LP reduced costs

We used LP reduced-cost elimination to create the test sets used in our study,
but our elimination code itself does not make use of LP results. It is natural
to consider integrating the two approaches. For example, during the
construction of a Hamilton-Tutte tree, an edge set $F$ is incompatible with
TSP optimality if the sum of the reduced costs of the edges in $F$ exceeds
$\Delta$, the integrality gap with respect to the LP relaxation and a known
TSP tour. This gives the possibility of avoiding a time-consuming search for a
proof that the set is nowhere $k$-optimal. Or, as a second example, one can
target certain edges for elimination or fixing based on reduced costs (large
positive values suggest edges for elimination and large negative values
suggest edges for fixing).

###  5.2 Using elimination in TSP algorithms

We describe below several possibilities for exploiting the elimination
routine.

####  5.2.1 LP relaxations over sparse edge sets

Unlike reduced-cost elimination, our combinatorial routines can remove edges
that take on positive values in an optimal solution to an LP relaxation,
potentially improving both the LP bound and the computational performance of
an LP solver. To illustrate this, in Table 9 we present LP information for the
three 100k+ instances. In each of the examples, the “(full)” values refer to
optimizing over the complete graph and the “(sparse)” values refer to
optimizing over the sparse edge sets with fixed edges set to one, using the
results described in Table 8.

Table 9: Improved LP optimality ratios/gaps with sparse edge sets LP Relaxation | E100k.0 | mona-lisa100k | usa115475
---|---|---|---
Subtour (full) | 1.0064790 | 1.0009764 | 1.0073171
Subtour (sparse) | 1.0058662 | 1.0009670 | 1.0066773
Gap Closed | 9.4% | 0.9% | 8.7%
Concorde (full) | 1.0000848 | 1.0000193 | 1.0001036
Concorde (sparse) | 1.0000826 | 1.0000193 | 1.0001022
Gap Closed | 1.6% | 0.0% | 1.4%
Concorde+Cuts (full) | 1.0000837 | 1.0000193 | 1.0001028
Concorde+Cuts (sparse) | 1.0000796 | 1.0000193 | 1.0000951
Gap Closed | 4.9% | 0.0% | 8.2%
Current Best | 1.0000486 | 1.0000127 | 1.0000803

The first block of rows in the table reports results for the standard subtour
relaxation of the TSP. This LP is the starting point for TSP cutting-plane
methods, going back to the work of Dantzig, Fulkerson, and Johnson [5]. The
first two rows list the optimality ratio provided by the relaxation, that is,
the value of the best known tour divided by the value of the dual LP solution.
The “Gap Closed” row reports the percentage of the integrality gap closed
using the sparse edge sets. The second block of rows reports the corresponding
results for the best LP relaxations produced by the Concorde code working with
the full edge sets. These relaxations were obtained using multiple runs of
Concorde, following the strategies outlined in Section 16.4 of Applegate et
al. [1]. The dual LP values were again improved for E100k.0 and usa115475, but
the sparse edge set did not improve the bound for the Mona Lisa instance.

In the third block of rows, we applied Concorde’s cutting-plane routines to
the LP relaxations reported in the second block. In these tests, we used the
command-line parameter “-mC100” to repeatedly run Concorde’s cutting-plane
loop, in each iteration increasing by 4 the value of $t_{max}$, which controls
the size of chunks considered in Concorde’s local-cuts routine. (See Section
16.1 of Applegate et al. [1].) The runs for E100k.0 and mona-lisa100k finished
in under 70 hours, while the two runs for usa115475 were terminated after 320
hours (with the full instance reaching $t_{max}=44$ and the spare instance
reaching $t_{max}=40$).

Finally, the “Current Best” row reports the optimality ratio for the current
best LP relaxations, found using the reduced edge sets and multiple additional
runs of Concorde. We must note that, whereas the first three blocks
demonstrate improvements due only to edge elimination, the better bounds
reported in this final row are due also in part to new cutting-plane search
techniques that were applied to the sparse LP relaxations. Indeed, it is
difficult to separate the contribution of the sparse graphs. Nonetheless, the
large improvement over the best full-graph results are suggestive of the
potential for application of the elimination routine.

An important observation is that, when comparing the full-graph runs and the
sparse-graph runs from the third block of rows in Table 9, there was a
significant reduction in the percentage of CPU time spent in solving the LP
relaxations arising in the cutting-plane process, as reported in Table 10.

Table 10: Percentage of computation time used by the LP solver LP Relaxation | E100k.0 | mona-lisa100k | usa115475
---|---|---|---
Concorde+Cuts (full) | 29.4% | 69.9% | 97.5%
Concorde+Cuts (sparse) | 25.1% | 57.1% | 93.0%

This would contribute to the improved current bounds obtained in computations
carried out with the reduced edge sets.

####  5.2.2 General MIP techniques

With a sufficiently sparse edge set, it is possible to apply general MIP
techniques to improve an LP relaxation for a TSP instance. This makes it
possible to combine the combinatorial cutting planes studied in the TSP
literature and the algebraic methods that have been very successfully applied
in MIP solvers. Note, however, that for exact TSP solutions it is important to
adopt numerically safe MIP cutting planes and bounds, such as those described
in [3] and [7].

####  5.2.3 Planar graphs

Another direction is to use the elimination process to produce graphs having
structural properties that can be exploited in TSP solution methods. For
example, the reduced graph for mona-lisa100k is planar, opening up several
possibilities for improving the execution of the cutting-plane method.

First, for planar graphs, Letchford [16] described a polynomial-time
separation algorithm for an important class of TSP inequalities called domino-
parity constraints. The implementation of his algorithm by Cook et al. [4]
required a heuristic method for extracting a planar graph from an LP solution,
but this step can now be avoided in the Mona Lisa example, resulting in exact
separation for this class of constraints.

Second, the solution set of the subtour relaxation of the TSP is a polytope in
$\mathbb{R}^{|E|}$, called the subtour polytope, defined by exponentially many
linear constraints. For planar graphs, Rivin [20] gave an extended formulation
of the subtour polytope, showing it is a projection of a polytope in
$\mathbb{R}^{3|E|}$ having no more than $|V|+6|E|$ constraints; a nice
description of Rivin’s formulation can be found in Pashkovich [17]. For Mona
Lisa, we created an integer-programming model consisting of the Riven system
together with the constraints from the current best LP relaxation, reported in
Table 9. Using the IBM CPLEX 22.1.0.0 mixed-integer programming solver, a
branch-and-bound run (after 7,342 search nodes) reported a lower bound of
5,757,132, potentially closing 19.2% of the integrality gap if the computation
could be carried out with numerically safe methods.

####  5.2.4 Edge sets in local-search heuristics

A reduced edge set can also be adopted in heuristic-search methods, serving as
a candidate set of edges for inclusion in a TSP tour. Indeed, the best-known
tour for the E100k.0 instance was found by Keld Helsgaun with LKH in 2013,
using results of an early version of our elimination code. This computation
reduced the tour length from 225,786,958 to 225,784,127. The reduction of
2,830 units improved the gap between the best-known tour length and current LP
bound by over 20.5%. Note that in this context it is also possible to employ
non-exact methods for reducing an edge set, such as in the study by Fischer
and Merz [8].

## References

  * [1] D. L. Applegate, R. E. Bixby, V. Chvátal, W. Cook. 2006. The Traveling Salesman Problem: A Computational Study. Princeton University Press. ISBN: 9780691129938.
  * [2] R. Bosch. 2019. Opt Art: From Mathematical Optimization to Visual Design. Princeton University Press. ISBN: 9780691164069.
  * [3] W. Cook, S. Dash, R. Fukasawa, M. Goycoolea. 2009. Numerically safe Gomory mixed-integer cuts. INFORMS Journal on Computing 21, 641–649. https://doi.org/10.1287/ijoc.1090.0324.
  * [4] W. Cook, D. G. Espinoza, M. Goycoolea. 2007. Computing with domino-parity inequalities for the traveling salesman problem (TSP). INFORMS Journal on Computing 19, 356–365. https://doi.org/10.1287/ijoc.1060.0204.
  * [5] G. Dantzig, R. Fulkerson, S. Johnson. 1954. Solution of a large-scale traveling-salesman problem. Journal of the Operations Research Society of America 2, 393–410. https://doi.org/10.1287/opre.2.4.393.
  * [6] G. B. Dantzig, D. R. Fulkerson, S. M. Johnson. 1959. On a linear-programming, combinatorial approach to the traveling-salesman problem. Operations Research 7, 58–66. https://doi.org/10.1287/opre.7.1.58.
  * [7] L. Eifler, A. Gleixner. 2023. A computational status update for exact rational mixed integer programming. Mathematical Programming 197, 793–812. https://doi.org/10.1007/s10107-021-01749-5.
  * [8] T. Fischer, P. Merz. 2007. Reducing the size of traveling salesman problem instances by fixing edges. In: C. Cotta, J. van Hemert (editors), Evolutionary Computation in Combinatorial Optimization, EvoCOP 2007. Lecture Notes in Computer Science 4446, Springer, 72–83. https://doi.org/10.1007/978-3-540-71615-0_7.
  * [9] M. M. Flood. 1956. The traveling-salesman problem. Operations Research 4, 61–75. https://doi.org/10.1287/opre.4.1.61.
  * [10] M. Held, R. M. Karp. 1971. The traveling-salesman problem and minimum spanning trees: Part II. Mathematical Programming 1, 6–25. https://doi.org/10.1007/BF01584070.
  * [11] K. Helsgaun. 2000. An effective implementation of the Lin-Kernighan traveling salesman heuristic. European Journal of Operational Research 126, 106–130. https://doi.org/10.1016/S0377-2217(99)00284-2.
  * [12] S. Hougardy, R. T. Schroeder. 2014. Edge elimination in TSP instances. In: Graph-Theoretic Concepts in Computer Science, WG 2014, D. Kratsch, I. Todinca (editors), Lecture Notes in Computer Science 8747, Springer, 275–286. https://doi.org/10.1007/978-3-319-12340-0_23.
  * [13] R. Jonker, T. Volgenant. 1984. Nonoptimal edges for the symmetric traveling salesman problem. Operations Research 32, 837–846. https://doi.org/10.1287/opre.32.4.837.
  * [14] P. L’Ecuyer. 1999. Good parameters and implementations for combined multiple recursive random number generators. Operations Research 7, 159–164. https://doi.org/10.1287/opre.47.1.159.
  * [15] P. L’Ecuyer, R. Simard, E. J. Chen, W. D. Kelton. 2002. An objected-oriented random-number package with many long streams and substreams. Operations Research 50, 1073–1075. https://doi.org/10.1287/opre.50.6.1073.358.
  * [16] A. N. Letchford. 2000. Separating a superclass of comb inequalities in planar graphs. Mathematics of Operations Research 25, 443–454. https://doi.org/10.1287/moor.25.3.443.12213.
  * [17] K. Pashkovich. 2012. Extended Formulations for Combinatorial Polytopes. Ph.D. Thesis, Otto-von-Guericke-Universität Magdeburg, Germany. http://dx.doi.org/10.25673/3850.
  * [18] G. Reinelt. 1991. TSPLIB—A traveling salesman problem library. ORSA Journal on Computing 3, 376–384. https://doi.org/10.1287/ijoc.3.4.376.
  * [19] G. Reinelt. 1995. TSPLIB95. http://comopt.ifi.uni-heidelberg.de/software/TSPLIB95/.
  * [20] I. Rivin. 2003. Combinatorial optimization in geometry. Advances in Applied Mathematics 31, 242–271. https://doi.org/10.1016/S0196-8858(03)00093-9.
  * [21] X. Zhong. 2022. Probabilistic analysis of edge elimination for Euclidean TSP. https://doi.org/10.48550/arXiv.1809.10469.