# coding: utf-8
""" Base Blades Module """
from itertools import combinations, permutations
import itertools



class GAbase(object):
	def __init__(self, n, metric):
		self.n = n
		self.metric = metric
		self.N = 2**n


def b_sort(xs):
	""" Sorts list stable, in-place. Counts the number of swaps 
	beetween adjacent elements. Corresponds to organizing the 
	vectors in the outer product using its anti-symmetry.
	
	b^a = -a^b   (outer product anti-symmetry on vectors)

	xs -- list (base blade representation)
	return
		xs -- sorted xs (in-place)
		s -- sign, (-1)**num_of_swaps
	"""
	s = 1
	for i in range(len(xs)):
		j = i
		while j > 0 and xs[j-1] > xs[j]:
			xs[j-1],xs[j] = xs[j],xs[j-1]
			s *= -1
			j -= 1

	return (xs, s)





def b_reduce(metric, xs):
	""" Reduces pairs of identical elements from sorted list.
	The pair (i,i) will be extracted as the number given by the metric[i].
	Corresponds to the geometric product on the base blades because 
	of the following:
	
	e*e = e|e + e^e    (only on vectors)
	e|e = metric[e]    (only in orthonormal basis?)	
	e^e = 0

	also:	
	e1*(f*e2) = f*(e1*e2)

	ex: [1 2 2 2 3 4 4] -> ([1,2,3], metric[2]*metric[4])

	metric -- list with norms [|_1|, |e1|, |e2|, ..., |en|]
	xs -- base blade repr (list)
	return -- pair (rs, f)
		rs -- reduced list
		f -- iloczyn wszystkich wartości redukcji
	"""

	rs = []
	f = 1.0
	
	l = len(xs)
	if l < 1:
		pass
	else:
		i = 1
		c = xs[0]
		j = 1
		while i < l:
			x = xs[i]
			if c == x:
				j += 1
			else:
				if j == 1:
					rs.append(c)
				elif j == 2:
					f *= metric[c]
				c = x
				j = 1
			i += 1

		if j == 1:
			rs.append(c)
		elif j == 2:
			f *= metric[c]
		
	return (rs, f)





		
def b_mul(metric, r1, r2):
	""" Geometric multiplication on the base blades representations 
	metirc
	r1, r2 -- tuples
	return
		r -- tuple
		f -- scalar coefficient
	"""
	r12 = list()    # need list
	r12.extend(r1)
	r12.extend(r2)
	r12,s = b_sort(r12)   # in-place operation
	r,f = b_reduce(metric, r12)
	return (tuple(r), s*f)





def construct_name(n, r):
	"""
	n -- wymiar przestrzeni
	r -- reprezentacja (lista indexow elementow w bazie)
	"""
	if len(r) == n:
		return 'I'
	if len(r) == 0:
		return '1'
	prefix = 'e'
	return prefix+''.join(map(str, r))

def construct_base(n):
	""" GA is generated by vector space.
	n -- wymiar przestrzeni
	"""
	index_by_repr = {}
	repr_by_index = {}
	index = 0
	vectors = range(1,n+1) # [1 2 3]
	for subspace_dim in range(n+1): # [0 1 2 3]
		for r in combinations(vectors, subspace_dim):
			index_by_repr[r] = index
			repr_by_index[index] = r
			index += 1
	return index_by_repr, repr_by_index


def construct_bst(ga, repr_by_index):			
	""" base symbol table """
	name_by_index = {}
	for i,r in repr_by_index.items():
		name_by_index[i] = construct_name(ga, r)
	return name_by_index



def construct_brt(repr_by_index):
	""" Construct base reverse table
	"""
	# base size
	N = len(repr_by_index)
	
	brt = [None]*N
	for i,r in repr_by_index.items():
		f = bb_reverse_coef(r)
		brt[i] = i,f
	return brt



def bb_reverse_coef(r):
	l = b_gr(r)
	s = (-1) ** (((l-1)*l)/2)
	return s



def construct_bmt(index_by_repr, repr_by_index, op):
	""" Construct multiplication table over base elements.
	index_by_repr, repr_by_index -- 2way index:repr mapping of base elements
	op -- operation taking two representations and returning (repr, float)
	"""
	assert len(index_by_repr) == len(repr_by_index)	
	
	# base size
	N = len(index_by_repr)
	# base multiplication table (2D array)
	bmt = [([None]*N) for _ in range(N)]
	for i1,r1 in repr_by_index.items():
		for i2,r2 in repr_by_index.items():
			r,f = op(r1, r2)
			i = index_by_repr[r]
			bmt[i1][i2] = (i,f) 

	return bmt


def format_bmt(bmt, bst):
	""" Construct human-readable bmt in textual form 
	bst -- numer elementu w bazie na symbol
	"""
	N = len(bmt)  # wielkosc bazy
	rows = []
	
	for i in range(N):
		pmt = bmt[i] # partial mul table (one line)
		line = []
		for j in range(N):
			k,f = pmt[j]  
			if abs(f) == 1.:
				name = bst[k]
				line.append('{0:>5}'.format(('-' if f<0 else ' ')+name))
			elif abs(f) == 0.:
				line.append('    0')
			else:
				assert 0
			
		rows.append(' '.join(line))
	
	return '\n'.join(rows)


def b_gr(r):
	""" Grade of the blade """
	return len(r)

def b_part(bl, k):
	""" k-graded part of the blade """
	if len(bl) == k:
		return bl, 1.0
	else:
		return (), 0.




def demo(n):
	""" Prints bmt-s of n-dimensional GA """
	N = 2**n
	#ga = GAbase(n = n, metric=[1,1,1,1,1,1,1,1,1])
	index_by_repr, repr_by_index = construct_base(n)
	bst = construct_bst(n, repr_by_index)
	m = [None,1,1,1]
	
	def b_geo(r1,r2): 
		return b_mul(m, r1, r2)

	def b_out(r1,r2):
		r3,f3 = b_mul(m, r1, r2)
		r5,f5 = b_part(r3, b_gr(r1)+b_gr(r2))
		return r5, f5*f3
	
	def b_inn(r1,r2):
		r3,f3 = b_mul(m, r1, r2)
		r5,f5 = b_part(r3, abs(b_gr(r1)-b_gr(r2)))
		return r5, f3*f5
	
	def b_left(r1,r2):
		r3,f3 = b_mul(m, r1, r2)
		r5,f5 = b_part(r3, b_gr(r2)-b_gr(r1))
		return r5, f3*f5

	def b_right(r1,r2):
		r3,f3 = b_mul(m, r1, r2)
		r5,f5 = b_part(r3, b_gr(r1)-b_gr(r2))
		return r5, f3*f5

	bmt_geo = construct_bmt(index_by_repr, repr_by_index, b_geo)
	bmt_out = construct_bmt(index_by_repr, repr_by_index, b_out)
	bmt_inn = construct_bmt(index_by_repr, repr_by_index, b_inn)
	bmt_left = construct_bmt(index_by_repr, repr_by_index, b_left)

	print 'geometric:'	
	print format_bmt(bmt_geo, bst)
	print 'outer:'
	print format_bmt(bmt_out, bst)
	print 'inner:'
	print format_bmt(bmt_inn, bst)
	print 'left:'
	print format_bmt(bmt_left, bst)

	brt = construct_brt(repr_by_index)		
	print 'brt'
	print brt


if __name__ == '__main__':
	demo(3)
