Flower snark experiment

experiments/flower-snarks.py

+#!/usr/bin/python
+
+DOC = """
+Calculate that the number of CDCs of Flower snark $J_k$ is $16^{(1-o(1))k}$.
+
+We calculate the step matrix $M$, show that 16 is its eigenvalue by finding
+its eigenvector. The eigenvector is guessed for the ratios of coefficients
+of $M^10 x$ and $M^9 x$ where $x$ is the initial vector. This eigenvector
+shares one non-zero coordinate with the final vector.
+"""
+
+import sys, os
+sys.path.append(os.path.dirname(__file__) + "/..")
+
+from graph_tools.all import *
+from sage.all import identity_matrix, QQ, DiGraph, matrix
+
+def test_eigenvalue(M, x):
+  H = M - x*identity_matrix(QQ, M.ncols(), sparse=True)
+  if H.is_singular():
+    print("  Value %i is an eigenvalue of M with geometric multiplicity %i" %
+      (x, H.ncols() - H.rank()))
+
+def matrix_to_rev_graph(M):
+  assert M.ncols() == M.nrows()
+  return DiGraph([range(M.ncols()), M.dict().keys()],
+    format="vertices_and_edges", loops=True)
+
+
+if __name__ == "__main__":
+  print(DOC)
+  print("Stabilizing...")
+  S = FlowerSnark.stabilize(CircuitDoubleCover, jordan=False)
+  print("Done")
+  M = S["step_matrix"]
+
+  print("\nThe rank of the step matrix M: %s" % M.rank())
+  upto = 16
+  print("Testing if integers up to %i are eigenvalues:" % upto)
+  for x in range(upto + 1):
+    test_eigenvalue(M, x)
+
+  fin = S["finalize"]
+
+  print("\nCheking correctness of the matrix")
+  N = 16
+  A = S["initial_vector"]
+  for i in range(3, N):
+    nu1 = FlowerSnark.graph(i).eval(CircuitDoubleCover)
+    A = M * A
+    nu2 = fin * A
+    assert nu1 == nu2
+    print("  ν(J_%i) = %i" % (i, nu1))
+  print("  Done")
+
+  print("\nCalculating the eigenvector for eigenvalue 16")
+  to_dict = lambda g: { b.boundary: b.value for b in g }
+  gadget_to_dict = lambda i: to_dict(FlowerSnark.gadget(i).eval_gadget(CircuitDoubleCover))
+
+  X10 = gadget_to_dict(10)
+  X9 = gadget_to_dict(9)
+  B = [ b for b in X10.keys() if b in X9 and X10[b] / X9[b] > 15 ]
+  m = min( X10[b] for b in B )
+  vec = { b: int(round(X10[b] / m)) for b in B }
+  print("  Guessed eigenvector (with %i non-zero coordinates):" % len(vec))
+  for b, v in sorted(vec.items()):
+    print("    %i at %s" % (v, b))
+  
+  print("  Checking it...")
+  nex = to_dict(FlowerSnark._next_gadget(FakeGadget(6,
+    [ Boundary(b, v) for b, v in vec.items() ])).eval_gadget(CircuitDoubleCover))
+
+  assert nex.keys() == vec.keys()
+  for b in vec:
+    assert vec[b] * 16 == nex[b], "Oops at %b" % b
+
+  assert set(S["variables"]).issuperset(vec.keys())
+  index = { b: i for i, b in enumerate(S["variables"]) }
+  vec_ = matrix(QQ, M.ncols(), 1, { (index[b], 0): v for b, v in vec.items() })
+  assert M * vec_ == 16 * vec_
+
+  print("  Done")
+
+  print("  Boundaries of with the non-zero value in both the eigenvector and the final vector:")
+  for i in range(vec_.nrows()):
+    if vec_[i, 0] > 0 and fin[0, i] > 0:
+      print("    %s" % (S["variables"][i],))
+