chiark - git - mdw - ocb-tv/blob - find-stretch.sage

   1 #! /usr/local/bin/sage
   2 ### -*-python-*-
   3
   4 k = GF(2)
   5
   6 def memoize(f):
   7   memo = dict()
   8   def _(*args):
   9     try:
  10       return memo[args]
  11     except KeyError:
  12       r = f(*args)
  13       memo[args] = r
  14       return r
  15   return _
  16
  17 @memoize
  18 def Abase(w, c):
  19   return matrix(k, [[j == i for j in xrange(w)]
  20                     for i in xrange(w)] +
  21                    [[j == i or j == i + c for j in xrange(w)]
  22                     for i in xrange(w)])
  23
  24 @memoize
  25 def A(w, c, i): return Abase(w, c)[i:i + w, 0:w]
  26
  27 def B(w, c, i, j): return A(w, c, i) + A(w, c, j)
  28
  29 def urank(w, c, i): return A(w, c, i).rank()
  30
  31 def domlen(w, c):
  32   for i in xrange(w):
  33     if urank(w, c, i) < w: return ZZ(i)
  34     for j in xrange(i):
  35       if B(w, c, i, j).rank() < w: return ZZ(i)
  36   return ZZ(w) #unlikely
  37
  38 def shifts(w):
  39   for c1 in xrange(0, 8):
  40     for c8 in xrange(w - 8, -8, -8):
  41       yield c8 + c1
  42
  43 def stretch_shift(w):
  44   want_bits = (w - 1).nbits()
  45   best_bits, best_dom, best_c = 0, 0, -1
  46   for c in shifts(w):
  47     d = domlen(w, c)
  48     bits = d.nbits()
  49     if bits == want_bits: return c, d
  50     elif bits > best_bits: best_bits, best_dom, best_c = bits, d, c
  51   return best_c, best_dom
  52
  53 for w in [64, 96, 128, 192, 256]:
  54   c, dom = stretch_shift(w)
  55   print '%3d: %3d [%d]' % (w, c, dom)