Add code to generate supernecklaces and -bracelets

Makefile

@@ -7,3 +7,12 @@ tsp/%.tsp: atsp/%.atsp
 
 atsp/%.atsp:
 	bin/mkatsp.py $* > "$@"
+
+necklace/%.gtsp:
+	bin/mkatsp.py --necklace $* > "$@"
+
+bracelet/%.gtsp:
+	bin/mkatsp.py --bracelet $* > "$@"
+
+alternating/%.gtsp:
+	bin/mkatsp.py --alternating $* > "$@"

README.md

@@ -106,3 +106,6 @@ concorde -o "$SUPERPERM_DIR"/concorde/out/6.out "$SUPERPERM_DIR"/tsp/6.tsp
 * The directory `superpermutations` contains various known superpermutations that are interesting for some reason.
 
 * Want to see what's in the superpermutations? Check out the `demutator` directory. (Be wary of the C code.)
+
+* To search for "supernecklaces" and "superbracelets" (strings which contain
+  each permutation only up to rotation and/or reflection), run `make necklace/6.gtsp` or `make bracelet/6.gtsp`. The resulting input files can be solved using the GLKH solver, available from http://webhotel4.ruc.dk/~keld/research/GLKH/ .

bin/mkatsp.py

@@ -2,6 +2,7 @@
 # -*- encoding: utf-8 -*-
 from __future__ import division
 
+from collections import defaultdict
 import math
 from itertools import permutations
 import optparse
@@ -10,6 +11,9 @@
 parser.add_option("-b", "--bound", type="int", help="only include edges up to this weight")
 parser.add_option("-n", "--no-cyclic", action="store_true", help="no edges between non-adjacent cyclic permutations")
 parser.add_option("-s", "--simple", action="store_true", help="only include edges from a.b -> b.a^r")
+parser.add_option("--necklace", action="store_true", help="consider permutations equivalent under rotation")
+parser.add_option("--bracelet", action="store_true", help="consider permutations equivalent under rotation and reflection")
+parser.add_option("--alternating", action="store_true", help="show weights and 1-cycles in the alternating graph")
 
 (options, args) = parser.parse_args()
 if len(args) != 1: parser.error("Wrong number of arguments")
@@ -27,6 +31,18 @@
 n_perms = len(perms)
 ordered = tuple(range(N))
 
+def cyclically_equivalent(p, q):
+    sp = "".join([ SYMBOLS[i] for i in p ])
+    sq = "".join([ SYMBOLS[i] for i in q ])
+    return sp in sq+sq
+
+def alternating_distance(p, q):
+    if p == q:
+        return 0
+    if cyclically_equivalent(p, q):
+        return INF
+    return max(distance((p[-1],) + p[:-1], q) - 1, 0)
+
 def distance(p, q):
     if q == ordered: return 0
 
@@ -42,15 +58,49 @@ def distance(p, q):
                 weight = n
                 break
 
-    if weight > 1 and options.no_cyclic:
-        sp = "".join([ SYMBOLS[i] for i in p ])
-        sq = "".join([ SYMBOLS[i] for i in q ])
-        if sp in sq+sq: return INF
+    if weight > 1 and options.no_cyclic and cyclically_equivalent(p, q):
+        return INF
 
     return weight
 
-print "NAME : superperm %d" % (N,)
-print "TYPE : ATSP"
+def normalise_necklace(p):
+    n = list(p)
+    i = n.index(0)
+    return n[i:] + n[:i]
+
+def normalise_bracelet(p):
+    n = normalise_necklace(p)
+    i = n.index(1)
+    j = n.index(2)
+    if i < j:
+        n.reverse()
+        return normalise_necklace(n)
+    return n
+
+def print_classes(normalise):
+    classes = defaultdict(set)
+    for (i, p) in enumerate(permutations(range(N))):
+        key = "".join(map(str, normalise(p)))
+        classes[key].add(i + 1)
+    num_classes = len(classes)
+    print "GTSP_SETS: %d" % (len(classes),)
+    print "GTSP_SET_SECTION"
+    for (i, c) in enumerate(classes.values()):
+        print "%d %s -1" % (i + 1, " ".join(map(str, c)))
+    print "EOF"
+
+if options.necklace:
+    print "NAME: supernecklace %d" % (N,)
+    print "TYPE: AGTSP"
+elif options.bracelet:
+    print "NAME: superbracelet %d" % (N,)
+    print "TYPE: AGTSP"
+elif options.alternating:
+    print "NAME: alternating %d" % (N,)
+    print "TYPE: AGTSP"
+else:
+    print "NAME : superperm %d" % (N,)
+    print "TYPE : ATSP"
 print "DIMENSION : %d" % (n_perms,)
 print "EDGE_WEIGHT_TYPE : EXPLICIT"
 print "EDGE_WEIGHT_FORMAT : FULL_MATRIX"
@@ -59,5 +109,14 @@ def distance(p, q):
 
 print "EDGE_WEIGHT_SECTION :"
 
-for p in perms:
-    print " ".join([ str(distance(p, q)) for q in perms ])
+if options.alternating:
+    for p in perms:
+        print " ".join([ str(alternating_distance(p, q)) for q in perms ])
+else:
+    for p in perms:
+        print " ".join([ str(distance(p, q)) for q in perms ])
+
+if options.necklace or options.alternating:
+    print_classes(normalise_necklace)
+if options.bracelet:
+    print_classes(normalise_bracelet)

bin/printtour.py

@@ -3,12 +3,17 @@
 from __future__ import division
 
 from itertools import permutations
-import sys
+import optparse
+
+parser = optparse.OptionParser(usage="%prog [options] N tour_filename")
+parser.add_option("-r", "--rotations", action="store_true",
+        help="fill in rotations of visited permutations")
+(options, args) = parser.parse_args()
 
 SYMBOLS = "123456789ABCDEFG"
 
-n = int(sys.argv[1])
-tour_filename = sys.argv[2]
+n = int(args[0])
+tour_filename = args[1]
 
 perms = list(permutations(range(n)))
 n_perms = len(perms)
@@ -24,6 +29,11 @@ def read_tour(tour_filename):
                 if ix == -1: break
                 yield "".join([ SYMBOLS[i] for i in perms[ix - 1] ])
 
+def rotations(xs):
+    for x in xs:
+        for i in range(n):
+            yield x[i:] + x[:i]
+
 def squash(xs):
     return reduce(lambda x, y: x + y[overlap(x, y):], xs, "")
 
@@ -33,4 +43,7 @@ def overlap(x, y):
         if x[len(x)-i:] == y[:i]:
             return i
 
-print squash(read_tour(tour_filename))
+tour = read_tour(tour_filename)
+if options.rotations:
+    tour = rotations(tour)
+print squash(tour)