#!/usr/bin/python

"""this script takes an integer linear combination of words in Z/nZ * Z/mZ, 
generated by x and y, and computes scl of it using scallop by lifting the 
chain to the finite index free normal subgroup which is the kernel of the 
map to Z/nZ x Z/mZ.  It's generated by strings of the form 
x^iy^jx^{n-i}y^{m-j}, so each generator is identified with a pair (i,j), and 
generator (i,j) is given letter (i-1)*(m-1) + (j-1)

the generators can be denoted (i,j,k), where k is 0 or 1 depending on whether 
it's an inverse (1=inverse)

It takes whatever power of the chain is necessary so that it's in the kernel

the syntax is 
./scallopFC.py [-q] [-G] g1 g2 n m C
where g1 and g2 are generators and n and m are the orders and C is the chain
-q makes it be quiet (just print the answer)
-G makes it use Gurobi
"""

import fractions
import subprocess
import time
import sys

alphabet = list('abcdefghijklmnopqrstuvwxyz')

SCALLOPDIR = './'
FNULL = open('/dev/null', 'w')

def scl(C, do5=False, timeout=1e8, rational=False, m=None):
  runString = ['scallop', '-local'] if not do5 else ['scallop', '-local', '-p5']
  if m != None:
    runString.append('-m' + m)
  sca = subprocess.Popen(runString+[str(C[0][i]) + C[1][i] for i in xrange(len(C[0]))],  \
                         executable=SCALLOPDIR+'scallop',         \
                         stdout=subprocess.PIPE,                  \
                         stderr=FNULL,                            \
                         cwd=SCALLOPDIR)
  startTime = time.time()
  while sca.poll() == None:
    if time.time() - startTime > timeout:
      try:
        sca.terminate()
        return -3.0
      except:
        return -4.0
    time.sleep(0.01)      
  scaout, scaerr = sca.communicate()
  s = scaout
  return fractions.Fraction(s.split('=')[1])




def multiply_words(w1, w2=None):
  if type(w1) == str:
    i = 0
    w1L = len(w1)
    w2L = len(w2)
    while i < w1L and i < w2L and w1[w1L-i-1] == w2[i].swapcase():
      i += 1
    return w1[:w1L-i] + w2[i:]
  elif type(w1) == list:
    return reduce(multiply_words, w1)

def word_power(w1, n):
  ans = ''
  for i in xrange(n):
    ans = multiply_words(ans, w1)
  return ans

def cyc_red(w):
  LW = len(w)
  if len(w) == 0 or len(w) == 1:
    return w
  else:
    i = 0
    while w[i] == w[LW-i-1].swapcase():
      i+=1
    return w[i:LW-i]

  
def inverse(w):
  return w[::-1].swapcase()


def gcd(a,b=None):
  if type(a) == int:
    t1 = max(a,b)
    t2 = min(a,b)
    while t2 != 0:
      t3 = t1 % t2
      t1 = t2
      t2 = t3
    return t1;
  elif type(a) == list:
    return reduce(gcd, a)


def lcm(a,b=None):
  if type(a) == int:
    if a==0 or b==0:
      return 1
    return a*b/gcd(a,b)
  elif type(a) == list:
    return reduce(lcm, a)


def write_gen(g, n, m):
  """write a gen in the form (i,j,k) as a letter"""
  i,j,k = g
  if k == 0:
    return alphabet[(i-1)*(m-1) + (j-1)]
  else:
    return alphabet[(i-1)*(m-1) + (j-1)].swapcase()
    


def rewrite_word(w, g1, g2, n, m):
  """rewrite a word in the kernel in terms of gens.  no simplification."""
  if len(w) == 0:
    return ''
  #first, convert the word to a list of powers
  pows = [[w[0].lower(), 0]]
  wLen = len(w)
  i=0
  while i < wLen:
    if w[i].lower() == pows[-1][0]:
      if w[i].islower():
        pows[-1][1] += 1
      else:
        pows[-1][1] -= 1
    else:
      if w[i].islower():
        pows.append( [w[i], 1] )
      else:
        pows.append( [w[i].lower(), -1] )
    i += 1
  
  #print pows
  #make all the powers positive
  for i in xrange(len(pows)):
    if pows[i][0] == g1:
      pows[i][1] %= n
    else:
      pows[i][1] %= m
  
  #scan through and remove any errant zero powers
  while True:
    changed_something = False
    for i in xrange(len(pows)):
      if pows[i][0] == g1:
        if pows[i][1] % n == 0:
          changed_something = True
          del pows[i]
          break
        else:
          pows[i][1] = pows[i][1] % n
      if pows[i][0] == g2:
        if pows[i][1] % m == 0:
          changed_something = True
          del pows[i]
          break
        else:
          pows[i][1] = pows[i][1] % m
    if changed_something:
      i=0
      while i<len(pows)-1:
        if pows[i][0] == pows[i+1][0]:
          pows[i][1] += pows[i+1][1]
          del pows[i+1]
        else:
          i += 1
    else:
      break
          
  
  if len(pows) == 0:
    return ''
  
  #now go through from the beginning and pull off the generators
  W = []  #denotes the gens, (i,j,I), where I=0 means no inverse, I=1 means inverse
  while len(pows) > 4:
    #print pows
    i,j,k,l = pows[0][1], pows[1][1], pows[2][1], pows[3][1]
    #the next gen is simply the first two characters
    if pows[0][0] == g1:
      W.append( (i, j, 0) )
    else:
      W.append( (j, i, 1) )
    #rewrite the beginning of the word so everything is still in the kernel
    if pows[0][0] == g1:  #this gen is a reguler gen
      if (i+k)%n == 0:    #the a powers are already ok
        if (j+l)%m == 0:  #and the b powers are too!
          pows = pows[4:]
        else:
          pows = [ [g2, ( l-(m-j) )%m] ] + pows[4:]
      else:
        pows = [ [g2, j], [g1, (k-(n-i))%n] ] + pows[3:]
    
    else:                 #this gen is an inverse
      if (i+k)%m == 0:    #the b powers are ok
        if (j+l)%n == 0:  #the a powers are too!
          pows = pows[4:]
        else:
          pows = [ [g1, (l-(n-j))%n] ] + pows[4:]
      else:
        pows = [ [g1, j], [g2, (k-(m-i))%m] ] + pows[3:]
  #end of the while loop
  #there's one remaining; throw it on
  #print pows
  i,j = pows[0][1], pows[1][1]
  if pows[0][0] == g1:
    W.append( (i, j, 0) )
  else:
    W.append( (j, i, 1) )
  if len(pows) != 4:
    print "Word not in the kernel?"
  #now write it in terms of actual gens
  W = ''.join([write_gen(c, n, m) for c in W])
  return W
  
  
#the canonical form is cyclically reduced and lexicographically first
#Just do this brute force for now, because whatever
def cyclic_canonical_form(w):
  """return the canonical form of the cyclic word w (i.e. [multiplier, word])"""
  W = cyc_red(w)
  wLen = len(W)
  if wLen == 0 or wLen == 1:
    return [1,w]
  its_bad = True
  for possible_period in xrange(1, wLen):
    if wLen % possible_period != 0:
      continue
    its_bad = False
    for j in xrange(possible_period, wLen, possible_period):
      if W[j:j+possible_period] != W[:possible_period]:
        its_bad = True
        break
    if its_bad:
      continue
    else:
      break
  if its_bad:
    possible_period = wLen
  #now possible_period is definitely the period
  multiplier = wLen/possible_period
  W = W[:possible_period]
  wLen = possible_period
  best_index = 0
  for i in xrange(1, wLen):
    if W[i:] + W[:i] < W[best_index:] + W[:best_index]:
      best_index = i
  return [multiplier, W[best_index:] + W[:best_index]]
  
  
def powers_needed_to_kernelize(C, g1, g2, n, m):
  """return a list of the powers needed to bring each word into the kernel."""
  ans = []
  for w in C:
    letter_counts = [w.count(g1) - w.count(inverse(g1)), \
                     w.count(g2) - w.count(inverse(g2)) ]
    while letter_counts[0] <= 0:
      letter_counts[0] += n
    while letter_counts[1] <= 0:
      letter_counts[1] += m
    min_pows = ( lcm(n, letter_counts[0])/letter_counts[0], \
                 lcm(m, letter_counts[1])/letter_counts[1] )
    ans.append(lcm(*min_pows))
  return ans


def reduce_chain(C):
  """returns [p/q, D], where C = p/q * D, and D is in canonical form."""
  new_chain = [ [], [] ]
  for i in xrange(len(C[0])):
    weight, w = C[0][i], C[1][i]
    if w == '':
      continue
    cfw = cyclic_canonical_form(w)
    if cfw[1] in new_chain[1]:
      new_chain[0][ new_chain[1].index(cfw[1]) ] += cfw[0]*weight
    else:
      new_chain[0].append( weight*cfw[0] )
      new_chain[1].append( cfw[1] )
    #print new_chain, cfw

  if len(new_chain[0]) == 0:
    return [1, [[],[]]]
  
  multiplier = gcd(new_chain[0])
  if multiplier != 1:
    new_chain[0] = [x/multiplier for x in new_chain[0]]
  
  return [multiplier, new_chain]
  
  
  
#if C is a list of words, then that's the chain.  it also accepts 
#[weights, words]
def lift_chain(C, g1, g2, n, m):
  """lift a chain to one in the kernel.  Returns [p/q,c], where 
  c is the chain, and p/q  is the multiplier (i.e. C is p/q times the 
  projection of c)"""
  if type(C) == str:
    new_C = [[1], [C]]
  elif type(C[0]) == str:
    new_C = [[1 for i in xrange(len(C))], C]
  else:
    new_C = C
  
  #first, get a list of the powers needed to get a word into the kernel
  powers_needed = powers_needed_to_kernelize(new_C[1], g1, g2, n, m)
  P = lcm(powers_needed)
  multiplier = fractions.Fraction(1, P)
  c = [ new_C[0], [word_power(w, P) for w in new_C[1]] ]
  
  #now, we need to conjugate each word by everything in G/H = Z/n x Z/m
  new_c = [ [], [] ]
  for i in xrange(len(c[0])):
    for j in xrange(n):
      for k in xrange(m):
        new_c[0].append(c[0][i])
        new_c[1].append( multiply_words( [ j*g1 + k*g2, c[1][i], ((m-k)%m)*g2 + ((n-j)%n)*g1] ) )
  c = new_c
  #and update the multiplier
  multiplier *= fractions.Fraction(1, n*m)
  
  #now we rewrite the words in the subgroup
  c[1] = [rewrite_word(w, g1, g2, n, m) for w in c[1]]
  
  #and reduce the chain to simplify
  new_multiplier, c = reduce_chain(c)
  multiplier *= new_multiplier
  
  return [multiplier, c]
  

def scl_FC(C, g1, g2, n, m, quiet=False, method=None):
  if type(C) == list:
    if type(C[0]) == list:
      multiplier, c = lift_chain(C, g1, g2, n, m)
    else:
      multiplier, c = lift_chain([ [1 for i in xrange(len(C))], C ], g1, g2, n, m)
  else:
    multiplier, c = lift_chain( [ [1], [C] ], g1, g2, n, m)
  if len(c[0]) == 0:
    return 0
  if not quiet:
    print "Running: scallop " + ' + '.join([ str(c[0][i]) + c[1][i]  for i in xrange(len(c[0]))] )
    sys.stdout.flush()
  s = scl(c, m=method)
  return multiplier * s

if len(sys.argv) > 2:
  quiet = False
  G = False
  while sys.argv[1][0] == '-':
    if sys.argv[1] == '-q':
      quiet = True
      del sys.argv[1]
    if sys.argv[1] == '-G':
      G = True
      del sys.argv[1]
  if sys.argv < 5:
    print "Read scallopFC.py for syntax"
    sys.exit(0)
  g1 = sys.argv[1]
  g2 = sys.argv[2]
  n = int(sys.argv[3])
  m = int(sys.argv[4])
  C = [[],[]]
  i = 5
  while i < len(sys.argv):
    j=0
    while not sys.argv[i][j].isalpha():
      j += 1
    if j==0:
      C[0].append(1)
    else:
      C[0].append( int(sys.argv[i][:j]) )
    C[1].append( sys.argv[i][j:] )
    i += 1
  
  print scl_FC(C, g1, g2, n, m, quiet, method=(None if G==False else 'GUROBI'))
  





