import random

tests = [
  'example',
  'gcddegree',
  'gcdx',
  'recipx',
  'recipxnone',
  'Gdelta',
  'Gf',
  'Vv',
  'zerotopright',
  'deggamma',
  'gamma',
  'gammav',
  '1overx',
  'x2',
  'x4',
]
coverage = {t:0 for t in tests}

def divstepsx(n,t,delta,f,g):
  assert t >= n and n >= 0
  f,g = f.truncate(t),g.truncate(t)
  kx = f.parent()
  x = kx.gen()
  u,v,q,r = kx(1),kx(0),kx(0),kx(1)

  while n > 0:
    f = f.truncate(t)
    if delta > 0 and g[0] != 0: delta,f,g,u,v,q,r = -delta,g,f,q,r,u,v
    f0,g0 = f[0],g[0]
    delta,g,q,r = 1+delta,(f0*g-g0*f)/x,(f0*q-g0*u)/x,(f0*r-g0*v)/x
    n,t = n-1,t-1
    g = kx(g).truncate(t)

  M2kx = MatrixSpace(kx.fraction_field(),2)
  return delta,f,g,M2kx((u,v,q,r))

def gcddegree(R0,R1):
  d = R0.degree()
  assert d > 0 and d > R1.degree()
  f,g = R0.reverse(d),R1.reverse(d-1)
  delta,f,g,P = divstepsx(2*d-1,2*d-1,1,f,g)
  return delta//2

def gcdx(R0,R1):
  d = R0.degree()
  assert d > 0 and d > R1.degree()
  f,g = R0.reverse(d),R1.reverse(d-1)
  delta,f,g,P = divstepsx(2*d-1,3*d-1,1,f,g)
  return f.reverse(delta//2)/f[0]

def recipx(R0,R1):
  d = R0.degree()
  assert d > 0 and d > R1.degree()
  f,g = R0.reverse(d),R1.reverse(d-1)
  delta,f,g,P = divstepsx(2*d-1,2*d-1,1,f,g)
  if delta != 0: return
  kx = f.parent()
  x = kx.gen()
  return kx(x^(2*d-2)*P[0][1]/f[0]).reverse(d-1)

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))

def bezout(R0,R1):
  # sage gcd and xgcd fail to return monic for, e.g., (2*x,0)
  if R0 == 0 and R1 == 0: return 0,0,0
  if R1 == 0:
    return R0/R0.leading_coefficient(),1/R0.leading_coefficient(),0
  if R0 == 0:
    return R1/R1.leading_coefficient(),0,1/R1.leading_coefficient()
  return xgcd(R0,R1)

k = GF(7)
kx.<x> = k[]
R0 = 2*x^7+7*x^6+1*x^5+8*x^4+2*x^3+8*x^2+1*x+8
R1 = 3*x^6+1*x^5+4*x^4+1*x^3+5*x^2+9*x+2
d = R0.degree()
f = kx(x^d*R0(1/x))
g = kx(x^(d-1)*R1(1/x))
assert f == 2+7*x+1*x^2+8*x^3+2*x^4+8*x^5+1*x^6+8*x^7
assert g == 3+1*x+4*x^2+1*x^3+5*x^4+9*x^5+2*x^6
delta = 1
for n in range(13):
  delta,f,g = divstep(delta,f,g)
assert (delta,f,g) == (0,2,6)
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(20)
  coeffs2 = randrange(1,100)
  coeffs3 = randrange(coeffs2)
  h = kx(sum(k.random_element()*x^i for i in range(coeffs1)))
  if h[0] == 0: continue
  R0 = h*kx(sum(k.random_element()*x^i for i in range(coeffs2)))
  if R0[0] == 0: continue
  R1 = h*kx(sum(k.random_element()*x^i for i in range(coeffs3)))
  
  M2kx = MatrixSpace(kx.fraction_field(),2)

  if R0.degree() > 0:
    if R0.degree() > R1.degree():
      G,U,V = bezout(R0,R1)
      assert G == U*R0+V*R1

      assert gcddegree(R0,R1) == G.degree()
      coverage['gcddegree'] += 1

      assert gcdx(R0,R1) == G
      coverage['gcdx'] += 1

      if G == 1:
        assert recipx(R0,R1) == V
        coverage['recipx'] += 1
      else:
        assert recipx(R0,R1) == None
        coverage['recipxnone'] += 1

      d = R0.degree()

      delta = 1
      f = kx(x^d*R0(1/x))
      g = kx(x^(d-1)*R1(1/x))

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

      for n in range(2*d):
        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)

      assert G.degree() == deltan[2*d-1]/2
      coverage['Gdelta'] += 1

      denom = fn[2*d-1](0)
      assert G == x^G.degree()*fn[2*d-1](1/x)/denom
      coverage['Gf'] += 1

      v = M2kx(prod(Tn[i] for i in reversed(range(0,2*d-1))))[0][1]
      assert V == x^(-d+1+G.degree())*v(1/x)/denom
      coverage['Vv'] += 1

      if G == 1:
        assert kx(v*x^(2*d-2)).degree() < d
        coverage['zerotopright'] += 1

  delta = 1
  f = R0
  g = R1
  d = R0.degree() + randrange(10)

  if f(0) and f.degree() <= d and g.degree() < d:
    gamma = bezout(f,g)[0]

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

    for n in range(2*d):
      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)

    assert gamma.degree() == fn[2*d-1].degree()
    coverage['deggamma'] += 1
    
    l = fn[2*d-1].leading_coefficient()
    assert gamma == fn[2*d-1]/l
    coverage['gamma'] += 1

    v = M2kx(prod(Tn[i] for i in reversed(range(0,2*d-1))))[0][1]
    left = (x^(2*d-2)*gamma) % f
    right = (kx(x^(2*d-2)*v)*g/l) % f
    assert left == right
    coverage['gammav'] += 1

    if f.degree() > 0:
      assert (kx((1-f/f(0))/x)*x)%f == 1
      coverage['1overx'] += 1

  if d >= 3:
    f = x^d-1
    assert bezout(x^(2*d-2),f)[:2] == (1,x^2)
    coverage['x2'] += 1

  if d >= 5:
    f = sum(x^i for i in range(d+1))
    assert bezout(x^(2*d-2),f)[:2] == (1,x^4)
    coverage['x4'] += 1

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