import random

tests = [
  'example',
  'STprod',
  'STprod1',
  'truncate',
]
coverage = {t:0 for t in tests}

def divstep(delta,f,g):
  kx = parent(f)
  x = kx.gen()
  if delta>0 and g[0]!=0:
    return 1-delta,g,kx((g[0]*f-f[0]*g)/x)
  return 1+delta,f,kx((f[0]*g-g[0]*f)/x)

def T(delta,f,g):
  kx = parent(f)
  x = kx.gen()
  M2kx = MatrixSpace(kx.fraction_field(),2)
  if delta>0 and g[0]!=0:
    return M2kx((0,1,g[0]/x,-f[0]/x))
  return M2kx((1,0,-g[0]/x,f[0]/x))

def S(delta,f,g):
  kx = parent(f)
  x = kx.gen()
  M2Z = MatrixSpace(ZZ,2)
  if delta>0 and g[0]!=0:
    return M2Z((1,0,1,-1))
  return M2Z((1,0,1,1))

k = GF(7)
kx.<x> = k[]
delta = 1
f = 2+7*x+1*x^2+8*x^3+2*x^4+8*x^5+1*x^6+8*x^7
g = 3+1*x+4*x^2+1*x^3+5*x^4+9*x^5+2*x^6

expected = [
  (1,[2,0,1,1,2,1,1,1],[3,1,4,1,5,2,2]),
  (0,[3,1,4,1,5,2,2],[5,2,1,3,6,6,3]),
  (1,[3,1,4,1,5,2,2],[1,4,4,0,1,6]),
  (0,[1,4,4,0,1,6],[3,6,1,2,5,2]),
  (1,[1,4,4,0,1,6],[1,3,2,2,5]),
  (0,[1,3,2,2,5],[1,2,5,3,6]),
  (1,[1,3,2,2,5],[6,3,1,1]),
  (0,[6,3,1,1],[1,4,4,2]),
  (1,[6,3,1,1],[0,2,4]),
  (2,[6,3,1,1],[5,3]),
  (-1,[5,3],[4,5,5]),
  (0,[5,3],[6,4]),
  (1,[5,3],[2]),
  (0,[2],[6]),
  (1,[2],[]),
  (2,[2],[]),
  (3,[2],[]),
  (4,[2],[]),
  (5,[2],[]),
  (6,[2],[]),
]

for n in range(20):
  assert (delta,list(f),list(g)) == expected[n]
  delta,f,g = divstep(delta,f,g)
coverage['example'] += 1

for loop in range(1000):
  q = random.choice([2,3,5,7,11,13,17,19,23,29])
  k = GF(q)
  kx.<x> = k[]
  coeffs1 = randrange(100)
  coeffs2 = randrange(100)
  delta = randrange(-9,10)
  h = kx(sum(k.random_element()*x^i for i in range(coeffs1)))
  if h[0] == 0: continue
  f = h*kx(sum(k.random_element()*x^i for i in range(coeffs2)))
  if f[0] == 0: continue
  g = h*kx(sum(k.random_element()*x^i for i in range(coeffs2)))

  deltan = {}
  fn = {}
  gn = {}
  Tn = {}
  Sn = {}

  for n in range(20):
    deltan[n] = delta
    fn[n] = f
    gn[n] = g
    Tn[n] = T(delta,f,g)
    Sn[n] = S(delta,f,g)

    delta,f,g = divstep(delta,f,g)

  M2kx = MatrixSpace(kx.fraction_field(),2)
  M2Z = MatrixSpace(ZZ,2)

  for n in range(-10,10):
    for m in range(-10,10):
      if n >= m and m >= 0:
        Tprod = M2kx(prod(Tn[i] for i in reversed(range(m,n))))
        Sprod = M2Z(prod(Sn[i] for i in reversed(range(m,n))))
        assert vector((fn[n],gn[n])) == Tprod * vector((fn[m],gn[m]))
        assert vector((1,deltan[n])) == Sprod * vector((1,deltan[m]))
        if n == m:
          assert Tprod == 1
          assert Sprod == 1
          coverage['STprod1'] += 1
        if n > m:
          assert kx(x^(n-m-1)*Tprod[0][0]).degree() <= n-m-1
          assert kx(x^(n-m-1)*Tprod[0][1]).degree() <= n-m-1
          assert kx(x^(n-m)*Tprod[1][0]).degree() <= n-m-1
          assert kx(x^(n-m)*Tprod[1][1]).degree() <= n-m-1
          assert Sprod[0][0] == 1
          assert Sprod[0][1] == 0
          assert Sprod[1][0] >= 2-(n-m)
          assert Sprod[1][0] <= n-m
          assert (Sprod[1][0] - (n-m)) % 2 == 0
          assert Sprod[1][1] in [1,-1]
          coverage['STprod'] += 1

  primedeltan = {}
  primefn = {}
  primegn = {}
  primeTn = {}
  primeSn = {}

  for m in range(20):
    t = randrange(20)
    delta = deltan[m]
    f = fn[m].truncate(t)
    g = gn[m].truncate(t)
    for n in range(m,20):
      primedeltan[n] = delta
      primefn[n] = f
      primegn[n] = g
      primeTn[n] = T(delta,f,g)
      primeSn[n] = S(delta,f,g)

      delta,f,g = divstep(delta,f,g)

    for n in range(20):
      if m < n and n <= m + t:
        assert primeTn[n-1] == Tn[n-1]
        assert primeSn[n-1] == Sn[n-1]
        assert primedeltan[n] == deltan[n]
        assert primefn[n].truncate(t-(n-m-1)) == fn[n].truncate(t-(n-m-1))
        assert primegn[n].truncate(t-(n-m)) == gn[n].truncate(t-(n-m))
        coverage['truncate'] += 1

for t in tests:
  print coverage[t],t